diff --git a/docs/figures/optimization/composite_maps.png b/docs/figures/optimization/composite_maps.png new file mode 100644 index 0000000..028e089 Binary files /dev/null and b/docs/figures/optimization/composite_maps.png differ diff --git a/docs/figures/optimization/curves.png b/docs/figures/optimization/curves.png new file mode 100644 index 0000000..c7eea7b Binary files /dev/null and b/docs/figures/optimization/curves.png differ diff --git a/docs/figures/optimization/path_animation.gif b/docs/figures/optimization/path_animation.gif new file mode 100644 index 0000000..f4c682b Binary files /dev/null and b/docs/figures/optimization/path_animation.gif differ diff --git a/docs/figures/optimization_case_study.py b/docs/figures/optimization_case_study.py new file mode 100644 index 0000000..e826e67 --- /dev/null +++ b/docs/figures/optimization_case_study.py @@ -0,0 +1,332 @@ +"""Reproduce the reference pair's before/after ISO-contact optimization ratings. + +Run from the repository root with the docs extra installed:: + + python docs/figures/optimization_case_study.py --stage ratings + +This stage invokes no FEM solver. It exports inputs and every returned ISO/AGMA factor; +AGMA safety factors remain null when qualified AGMA material stress numbers are absent. +""" + +from __future__ import annotations + +import argparse +import hashlib +import json +import os +import shutil +from dataclasses import asdict, replace +from enum import Enum +from pathlib import Path + +from iso6336_figures import readme_pair + +from gearbox import Transmission, __version__ +from gearbox.fem.post import LineOfAction, contact_sweep +from gearbox.optimization import ( + feasible_pinion_profile_shift_range, + minimum_contact_stress_profile_shifts, +) +from gearbox.standards import agma2101, iso6336 + +POSITIONS = 35 # sampled contact positions per design +MARGIN = 0.15 # mm kept from the exact path ends, where the contact opens or closes + + +def designs() -> dict[str, Transmission]: + """Return the unchanged reference pair and the optimum of each standard's contact stress.""" + original = readme_pair() + pairs = {"baseline": original} + for standard in ("ISO", "AGMA"): + shifts = minimum_contact_stress_profile_shifts(original, standard=standard) + pairs[f"optimized_{standard.lower()}"] = replace( + original, + pinion=replace(original.pinion, profile_shift_coefficient=shifts.pinion), + wheel=replace(original.wheel, profile_shift_coefficient=shifts.wheel), + ) + return pairs + + +def ratings(pairs: dict[str, Transmission]) -> dict: + """Collect inputs, geometry and public rating results without material-data substitutions.""" + before = pairs["baseline"] + output = { + "schema_version": 1, + "software_version": __version__, + "objectives": { + "optimized_iso": "minimize_governing_iso_contact_stress", + "optimized_agma": "minimize_agma_contact_stress", + }, + "nominal_force_n": before.nominal_tangential_force, + "feasible_pinion_shift_range": feasible_pinion_profile_shift_range(before), + "designs": {}, + } + for name, pair in pairs.items(): + output["designs"][name] = { + "inputs": asdict(pair), + "geometry": { + "x1": pair.pinion.profile_shift_coefficient, + "x2": pair.wheel.profile_shift_coefficient, + "center_distance_mm": pair.working_center_distance, + "contact_ratio": pair.transverse_contact_ratio, + }, + "ISO": { + "pitting": asdict(iso6336.calculate_pitting(pair)), + "bending": asdict(iso6336.calculate_bending(pair)), + }, + "AGMA": { + "pitting": asdict(agma2101.calculate_pitting(pair)), + "bending": asdict(agma2101.calculate_bending(pair)), + }, + } + return output + + +def _positions(pair: Transmission) -> list[float]: + """Sample positions on the line of action, ``MARGIN`` inside the exact path ends.""" + line = LineOfAction.of(pair) + start, end = -line.approach + MARGIN, line.recess - MARGIN + return [start + (end - start) * i / (POSITIONS - 1) for i in range(POSITIONS)] + + +def sweep(pairs: dict[str, Transmission], cache: Path) -> dict: + """Solve both designs at the same nominal force over their full paths of contact. + + Resumable: solved positions are read back from ``cache``. Root stress is the envelope over + each sector's fillet nodes, not the history of one followed root. + """ + force = pairs["baseline"].nominal_tangential_force + output = {"schema_version": 1, "software_version": __version__, "nominal_force_n": force, + "positions_per_design": POSITIONS, "margin_mm": MARGIN, "designs": {}} # fmt: skip + for name, pair in pairs.items(): + line = LineOfAction.of(pair) + points = contact_sweep(pair, cache / name, _positions(pair), force) + output["designs"][name] = { + "path_mm": {"approach": line.approach, "recess": line.recess}, + "samples": [ + { + "position_mm": p.position, + "normalized_position": (p.position + line.approach) + / (line.approach + line.recess), + "max_pressure_mpa": p.followed.pressure, + "followed_tooth": p.followed.tooth, + "root_stress_pinion_mpa": p.root_stress_pinion, + "root_stress_wheel_mpa": p.root_stress_wheel, + "requested_torque_nmm": p.requested_torque, + } + for p in points + ], + } + return output + + +LABELS = { + "baseline": ("original (x1 = 0)", "#555555"), + "optimized_iso": ("ISO optimum", "#1f77b4"), + "optimized_agma": ("AGMA optimum", "#d95f02"), +} + + +EDGE_PRESSURE = 2000.0 # MPa; a second pair above this is in tip-edge contact at its transition + + +def _edge_contact(point) -> bool: + """True when a pair other than the followed one carries a tip-edge pressure peak.""" + return any(p is not point.followed and p.pressure > EDGE_PRESSURE for p in point.patches) + + +def _solved(pairs: dict[str, Transmission], cache: Path) -> dict: + """Read the solved sweep of every design from the cache; stop if any position is missing.""" + force = pairs["baseline"].nominal_tangential_force + points = {} + for name, pair in pairs.items(): + positions = _positions(pair) + key = hashlib.sha1(repr((pair, force, [])).encode()).hexdigest()[:12] + missing = [ + x for x in positions if not (cache / name / f"sweep_{key}_{x:+.4f}.npz").is_file() + ] + if missing: + raise SystemExit(f"{name}: {len(missing)} positions not solved yet; not solving") + points[name] = contact_sweep(pair, cache / name, positions, force) + return points + + +def animation(pairs: dict[str, Transmission], cache: Path, out: Path) -> None: + """Animate the composite von Mises state of the three designs along the path of contact. + + Frame ``i`` shows sample ``i`` of every design, i.e. the same normalized position. + """ + import matplotlib.pyplot as plt + import numpy as np + from matplotlib.animation import PillowWriter + + from gearbox.fem.post import plot_field + + points = _solved(pairs, cache) + out.mkdir(parents=True, exist_ok=True) + top = 800.0 # MPa, fixed colour scale for every frame and design + fig, axes = plt.subplots(1, 3, figsize=(15, 5.4)) + fig.subplots_adjust(left=0.04, right=0.88, bottom=0.1, top=0.82, wspace=0.12) + colour_axis = fig.add_axes([0.9, 0.15, 0.015, 0.62]) + writer = PillowWriter(fps=4) + with writer.saving(fig, str(out / "path_animation.gif"), dpi=90): + for i in range(POSITIONS): + contour = None + for ax, (name, series) in zip(axes, points.items(), strict=True): + ax.clear() + p = series[i] + centre = np.asarray(p.followed.point) + contour = plot_field(ax, p.stress("von_mises"), vmax=top, edges=False) + ax.plot(*centre, "o", mfc="white", mec="k", ms=5) + ax.set(xlim=(centre[0] - 7, centre[0] + 7), ylim=(centre[1] - 7, centre[1] + 7)) + ax.set_xticks([]) + ax.set_yticks([]) + note = "\nsecond pair in tip-edge contact" if _edge_contact(p) else "" + ax.set_title( + f"{LABELS[name][0]}\ns = {p.position:+.2f} mm, " + f"contact pressure {p.followed.pressure:.0f} MPa{note}", + fontsize=9, + ) + colour_axis.clear() + fig.colorbar(contour, cax=colour_axis, label="von Mises stress [MPa]") + fig.suptitle( + f"Composite stress state along the path of contact, sample {i + 1}/{POSITIONS} " + f"(normalized position {i / (POSITIONS - 1):.2f} of the sampled path)" + ) + writer.grab_frame() + print("Wrote path_animation.gif") + + +def figures(pairs: dict[str, Transmission], cache: Path, out: Path) -> None: + """Draw the before/after curves and the composite stress maps from the solved sweeps. + + Reads only the sweep cache: no solver is started, and a missing position stops the stage. + """ + import matplotlib.pyplot as plt + import numpy as np + + from gearbox.fem.post import plot_field + + points = _solved(pairs, cache) + force = pairs["baseline"].nominal_tangential_force + out.mkdir(parents=True, exist_ok=True) + + fig, axes = plt.subplots(1, 3, figsize=(16, 4.8)) + for name, series in points.items(): + line = LineOfAction.of(pairs[name]) + label, colour = LABELS[name] + u = [(p.position + line.approach) / (line.approach + line.recess) for p in series] + for ax, values in zip( + axes, + ( + [p.followed.pressure for p in series], + [p.root_stress_pinion for p in series], + [p.root_stress_wheel for p in series], + ), + strict=True, + ): + ax.plot(u, values, "o-", ms=3, color=colour, label=label) + flagged = [i for i, p in enumerate(series) if _edge_contact(p)] + ax.plot( + [u[i] for i in flagged], + [values[i] for i in flagged], + "o", + ms=6, + mfc="white", + mec=colour, + label=None, + ) + for ax, title in zip( + axes, + ("contact pressure of the followed tooth", "pinion root envelope", "wheel root envelope"), + strict=True, + ): + ax.set(title=title, xlabel="normalized position on the full path of contact (0 = A, 1 = E)") + ax.set_ylabel("pressure [MPa]" if "pressure" in title else "max principal stress [MPa]") + ax.spines[["top", "right"]].set_visible(False) + ax.grid(alpha=0.25) + axes[0].legend(frameon=False) + axes[0].plot([], [], "o", mfc="white", mec="k", label="second pair in tip-edge contact") + axes[0].legend(frameon=False, fontsize=8) + fig.suptitle(f"FEM, 2D plane strain, common Ft = {force:.0f} N, 35 positions per design") + fig.tight_layout() + fig.savefig(out / "curves.png", dpi=170) + plt.close(fig) + + # The middle of the path lies in the single-pair zone of all three designs; the pitch point C + # is not a common position (the AGMA optimum engages only after it). + peaks = {} + for name, series in points.items(): + line = LineOfAction.of(pairs[name]) + middle = (line.recess - line.approach) / 2 + peaks[name] = min(series, key=lambda p: abs(p.position - middle)) + windows = {} + for name, p in peaks.items(): + centre = np.asarray(p.followed.point) + windows[name] = (centre, np.linalg.norm(p.nodes - centre, axis=1) < 7.0) + vmax = max(float(np.nanmax(p.von_mises[windows[n][1]])) for n, p in peaks.items()) + fig, axes = plt.subplots(1, 3, figsize=(16, 5.6)) + contour = None + for ax, (name, p) in zip(axes, peaks.items(), strict=True): + centre = windows[name][0] + contour = plot_field(ax, p.stress("von_mises"), vmax=vmax, edge_width=0.1) + ax.plot(*centre, "o", mfc="white", mec="k", ms=6) + ax.set(xlim=(centre[0] - 7, centre[0] + 7), ylim=(centre[1] - 7, centre[1] + 7)) + ax.set_title( + f"{LABELS[name][0]}\ncontact pressure {p.followed.pressure:.0f} MPa at s = {p.position:+.2f} mm", + fontsize=10, + ) + ax.set_xlabel("x [mm]") + fig.colorbar(contour, ax=list(axes), shrink=0.8, label="von Mises stress [MPa]") + fig.suptitle( + "Composite stress state: contact and bending in the same solution at the middle of the path of contact (same colour scale)" + ) + fig.savefig(out / "composite_maps.png", dpi=170) + plt.close(fig) + print("Wrote curves.png and composite_maps.png") + + +def _json_value(value): + if isinstance(value, Enum): + return value.value + raise TypeError(f"Cannot encode {type(value).__name__} as study data") + + +def main() -> None: + """Run the selected reproducible study stage.""" + parser = argparse.ArgumentParser(description=__doc__) + parser.add_argument( + "--stage", + choices=["ratings", "fem", "figures", "animation"], + default="ratings", + help="ratings: no solver; fem: resumable CalculiX sweeps (needs ccx); figures/animation: plots from the solved cache, no solver", + ) + parser.add_argument( + "--output-dir", + type=Path, + default=Path("examples/output/optimization-case-study"), + help="data output folder (default: examples/output/optimization-case-study)", + ) + args = parser.parse_args() + args.output_dir.mkdir(parents=True, exist_ok=True) + if args.stage == "animation": + animation(designs(), args.output_dir / "cache", args.output_dir / "figures") + return + if args.stage == "figures": + figures(designs(), args.output_dir / "cache", args.output_dir / "figures") + return + if args.stage == "fem": + if not (os.environ.get("GEARBOX_CCX") or shutil.which("ccx")): + raise SystemExit("CalculiX not found: put ccx on the PATH or set GEARBOX_CCX") + data = sweep(designs(), args.output_dir / "cache") + name = "fem.json" + else: + data = ratings(designs()) + name = "ratings.json" + text = json.dumps(data, indent=2, default=_json_value, allow_nan=False) + "\n" + (args.output_dir / name).write_text(text) + print(f"Wrote {name} for the baseline and both optimized designs.") + + +if __name__ == "__main__": + main() diff --git a/docs/optimization.md b/docs/optimization.md index df8beb8..13b2425 100644 --- a/docs/optimization.md +++ b/docs/optimization.md @@ -152,6 +152,64 @@ pinion root (430 → 428 MPa) while the wheel root gets 40 % worse. The contact around the start point, so the optimizer finds only a modest gain. Optimization tells you where the feasible region ends; whether you want to be there is an engineering decision. +## Case study: one pair, two optima, one FEM check + +The same 22/40-tooth pair as above (module 2.5 mm, 40 kW, 1450 rpm, application factor 1.3) is +optimized twice with the contact-stress objective, once per standard, and every result is then +checked with the 2D FEM model at the same nominal tangential force (9579 N). The shift sum stays +at zero, so the centre distance never changes. + +| Design | $x_1$ | ISO $\sigma_{H}$ (MPa) | ISO min. $S_H$ | ISO min. $S_F$ | AGMA $\sigma_c$ (MPa) | FEM peak pressure (MPa) | FEM wheel root (MPa) | +|---|---|---|---|---|---|---|---| +| Original | 0 | 1492 | 0.85 | 1.76 | 1846 | 1392 | 294 | +| ISO optimum | +0.326 | 1439 | 0.88 | 1.69 | 1762 | 1318 | 328 | +| AGMA optimum | +1.023 | 1496 | 0.85 | 1.16 | 1689 | 1305 | 466 | + +Three things stand out: + +* **Each standard finds its own optimum.** The ISO optimum equalizes the pinion and wheel contact + stresses; the AGMA optimum goes all the way to the largest feasible shift. Judged by the *other* + standard, each one looks worse: the AGMA optimum does not improve the ISO contact stress at all. +* **Lower contact stress costs wheel-root margin.** The wheel root gets thinner as the pinion gets + stronger. The FEM agrees on the direction and size of the trend: wheel-root stress rises 12 % + with the ISO optimum and 58 % with the AGMA one, while the peak contact pressure falls only 5-6 %. +* **The three designs engage differently.** The path of contact moves along the line of action + as the shifts change, so the curves below use the normalized position on the whole path + (0 = A, 1 = E). + +```{figure} figures/optimization/curves.png +:alt: Contact pressure and pinion and wheel root stress along the path of contact for the original design and the ISO and AGMA optima. + +Contact pressure and root-stress envelopes along the full path of contact, 35 FEM positions per +design. Open markers flag positions where a second pair is in tip-edge contact and its pressure +peak distorts the neighbouring values. +``` + +The FEM solves contact and bending in the same simulation, so each result is a complete stress +state, not a local number. The first figure shows it at the middle of each path, the animation +follows it along the whole path. + +```{figure} figures/optimization/composite_maps.png +:alt: Von Mises stress maps of the three designs at the middle of the path of contact. + +Von Mises stress at the middle of each path of contact, same colour scale. +``` + +```{figure} figures/optimization/path_animation.gif +:alt: Animation of the von Mises stress of the three designs along the path of contact. + +Von Mises stress along the path of contact, one frame per FEM position. +``` + +This is a demonstration of the tools, not a research study. The FEM model is 2D plane strain +without friction, no mesh-convergence study was done, and the ISO and AGMA values are factored +rating stresses while the FEM values are nominal: they show trends, they do not validate the +standards. Tip-edge contact where the next pair engages produces local peaks that are artifacts +of the model. The ISO safety factors below 1 reflect the example loading, not a recommendation. + +To reproduce it, run `python docs/figures/optimization_case_study.py --stage ratings`, then +`--stage fem` (CalculiX needed, about three hours), then `--stage figures` and `--stage animation`. + ## Limits * **Internal gears** raise `ValueError`: the geometric limits used here are those of external gears. diff --git a/odd/tasks/optimization-case-study.md b/odd/tasks/optimization-case-study.md new file mode 100644 index 0000000..f88083e --- /dev/null +++ b/odd/tasks/optimization-case-study.md @@ -0,0 +1,34 @@ +# Optimization before-and-after case study + +Extend the guide with a reproducible ISO-contact optimization study and real FEM curves. +Keep numerical-library algorithms unchanged; joint safety optimization is a separate feature. + +## Specs +S1: "sabes me gustaria hacer una cosa en el apartado de la optimizacion, poner una transmision, resolverla con las normas, luego optimizarla, resolverla con las normas y FEM, y ver todos los resultados, esa curva de contacto completo antes y despues de la optimizacion" +S2: "para contacto y flexion" +S3: "Minimizar contacto ISO — optimizador existente" (explicit objective selection). +S4: "pero puedes agregar la otra idea como un issue de feature" + +## Tasks +T1 [done] S4 | route: repository feature form, duplicate search, one create and exact read-back | reference: github.com/efirvida/python-gearbox#15 +T2 [done] S1-S3 | route: public-API case generator and deterministic tests; ISO/AGMA full factor outputs, original and optimized geometry | commit: pending +T3 [done] S1-S3 | route: pilot then resumable 35-position CalculiX sweeps per design at equal nominal force | commit: pending +T4 [done] S1-S3 | route: composite-state von Mises maps per design, original before/after pressure and root-envelope figures, tables and narrative, links from optimization | commit: pending +T5 [done] S1-S3 | route: public example tests, data/provenance checks, strict Sphinx, lint/format and visual inspection | commit: pending + +## Log +L1: "sabes me gustaria hacer una cosa en el apartado de la optimizacion, poner una transmision, resolverla con las normas, luego optimizarla, resolverla con las normas y FEM, y ver todos los resultados, esa curva de contacto completo antes y despues de la optimizacion" +L2: "para contacto y flexion" +L3: User explicitly selected "Minimizar contacto ISO — optimizador existente" rather than the proposed joint minimum-safety objective. +L4: "pero puedes agregar la otra idea como un issue de feature" +L5: Created and read-back confirmed https://github.com/efirvida/python-gearbox/issues/15, enhancement, no duplicate. It proposes maximizing min(SH1,SH2,SF1,SF2) with fixed shift sum; not implemented in this study. +L6: Implementation plan: use README pair22/40, normal module2.5mm, width34mm,40kW,1450rpm,10000h, application factor1.3, initial shifts0/0. Optimize only the shift distribution at fixed sum. Use both standards before/after; do not invent AGMA stress numbers or safety factors. +L7: FEM uses the actual common nominal tangential force (about9579N), not the older8000N illustration. Nominal FEM and factored standard stresses are different quantities; no claimed normative validation. Use35 sample positions per design, normalize the sampled path for comparison, retain physical positions and disclose excluded exact endpoints. +L8: API inspection: contact_sweep returns followed-tooth patch pressure and largest principal stress over each sector's fillet nodes. Label root curves as sector envelopes, not the history of one followed root. Results are2D plane strain, frictionless, with common solver/mesh settings; no3D run or mesh-convergence claim. +L9: First documentation version is complete on master582d61b. Docs39987991175 and CI39987991426 passed after a formatting-only test correction. New study is isolated on docs/optimization-case-study; do not push/merge it as part of the already closed first-version delivery without a new publication decision. +L10: Pilot with conda env gearbox-fem ccx 2.23 (GEARBOX_CCX=~/miniconda3/envs/gearbox-fem/bin/ccx; Doctorado ccx_2.21 lacks libs): ~3 min/position, C-point baseline p=1338 MPa vs optimized 1082 MPa. Full 70-solve sweep launched in background (--stage fem, resumable cache in examples/output/optimization-case-study/cache), ~3.5 h. +L11: User asked to optimize with BOTH standards and compare against FEM (comparison matrix). Added optimized_agma (x1=1.023, minimize AGMA contact); designs now baseline/optimized_iso/optimized_agma, 105 solves (~5 h). Matrix: 3 designs x ISO/AGMA ratings x FEM. +L12: User approved adding von Mises/principal maps of the composite (contact+bending, same solve) stress state at the max-pressure position per design, from cached nodal fields; no new solves. 2D plane strain, frictionless. +L13: FEM sweep complete (3 designs x 35). C is outside the AGMA-optimum path (starts s=+0.17 mm), so composite maps use the middle of each path (single-pair zone in all three). Open: unexplained pinion-root spike (~298 MPa, u~0.68) in AGMA design; no mesh-convergence claim. +L14: Pinion-root spike (AGMA, s=+7.04/+7.33) coincides with a second pair in tip-edge contact (4000-10000 MPa transition peaks in ALL designs, EDGE_PRESSURE=2000 flags them). Data kept; flagged with open markers. Added --stage animation (path_animation.gif, fixed 800 MPa scale). User declined a table. +L15: User authorized publish, merge and deploy (no scientific rigor needed). Docs section written, strict Sphinx, ruff, 318 tests pass; PR to master. diff --git a/tests/test_optimization_case_study.py b/tests/test_optimization_case_study.py new file mode 100644 index 0000000..d6dfe6c --- /dev/null +++ b/tests/test_optimization_case_study.py @@ -0,0 +1,80 @@ +"""Reproducible public commands for the documented optimization case.""" + +import json +import subprocess +import sys +from copy import deepcopy +from pathlib import Path + +import pytest + +SCRIPT = Path(__file__).resolve().parents[1] / "docs/figures/optimization_case_study.py" + + +def test_ratings_compare_both_standards_without_changing_operating_conditions(tmp_path): + result = subprocess.run( + [sys.executable, str(SCRIPT), "--stage", "ratings", "--output-dir", str(tmp_path)], + capture_output=True, + text=True, + check=False, + ) + assert result.returncode == 0, result.stdout + result.stderr + assert result.stdout == "Wrote ratings.json for the baseline and both optimized designs.\n" + assert result.stderr == "" + data = json.loads((tmp_path / "ratings.json").read_text()) + assert sorted(p.name for p in tmp_path.iterdir()) == ["ratings.json"] + assert data["objectives"] == { + "optimized_iso": "minimize_governing_iso_contact_stress", + "optimized_agma": "minimize_agma_contact_stress", + } + assert data["nominal_force_n"] == pytest.approx(9579.2317, rel=1e-7) + before, after, agma = ( + data["designs"][name] for name in ("baseline", "optimized_iso", "optimized_agma") + ) + assert before["geometry"]["x1"] == before["geometry"]["x2"] == 0 + assert after["geometry"]["x1"] == pytest.approx(0.325917, abs=1e-5) + assert after["geometry"]["x1"] + after["geometry"]["x2"] == pytest.approx(0) + assert after["geometry"]["center_distance_mm"] == before["geometry"]["center_distance_mm"] + assert after["geometry"]["contact_ratio"] >= 1.2 + restored = deepcopy(after["inputs"]) + for member in ("pinion", "wheel"): + restored[member]["profile_shift_coefficient"] = 0 + assert restored == before["inputs"] + assert before["ISO"]["pitting"]["contact_stress_pinion"] == pytest.approx(1492.234, rel=1e-6) + assert before["ISO"]["bending"]["tooth_root_stress_pinion"] == pytest.approx(456.901, rel=1e-5) + assert before["AGMA"]["pitting"]["contact_stress"] == pytest.approx(1845.574, rel=1e-6) + assert before["AGMA"]["bending"]["bending_stress_wheel"] == pytest.approx(462.138, rel=1e-5) + assert agma["geometry"]["x1"] == pytest.approx(1.02303, abs=1e-4) + assert agma["geometry"]["x1"] + agma["geometry"]["x2"] == pytest.approx(0) + assert agma["geometry"]["contact_ratio"] >= 1.2 + assert agma["AGMA"]["pitting"]["contact_stress"] < before["AGMA"]["pitting"]["contact_stress"] + for design in (before, after, agma): + for standard in ("ISO", "AGMA"): + for mode in ("pitting", "bending"): + assert design[standard][mode]["dynamic_factor"] >= 1 + if standard == "AGMA": + assert design[standard][mode]["safety_factor_pinion"] is None + assert design[standard][mode]["safety_factor_wheel"] is None + assert max(after["ISO"]["pitting"][f"contact_stress_{g}"] for g in ("pinion", "wheel")) < max( + before["ISO"]["pitting"][f"contact_stress_{g}"] for g in ("pinion", "wheel") + ) + # Contact optimization must not hide the wheel-root trade-off. + assert ( + after["ISO"]["bending"]["tooth_root_stress_wheel"] + > before["ISO"]["bending"]["tooth_root_stress_wheel"] + ) + + +def test_help_describes_the_ratings_stage_without_starting_a_solver(tmp_path): + result = subprocess.run( + [sys.executable, str(SCRIPT), "--help"], + cwd=tmp_path, + capture_output=True, + text=True, + check=False, + ) + assert result.returncode == 0, result.stdout + result.stderr + assert "--stage" in result.stdout and "ratings" in result.stdout + assert "--output-dir" in result.stdout + assert result.stderr == "" + assert not list(tmp_path.iterdir())