diff --git a/motors/synrm_3f_36s.py b/motors/synrm_3f_36s.py index dd3a840..93af0c0 100644 --- a/motors/synrm_3f_36s.py +++ b/motors/synrm_3f_36s.py @@ -6,6 +6,7 @@ from machine_design.designs.computation import ComputationBase from machine_design.designs.geometry import GeometryBase +from machine_design.generic.transforms import to_dq class Geometry(GeometryBase): @@ -174,12 +175,13 @@ def set_oper_params(self): f = 50 # [Hz] RotSpeed = 60 * f / self.geometry.PolePairs # [rpm] w = 2 * np.pi * f + self.InitPos = -30 # deg self.oper_params = { "Id": "0.0A", "Iq": "0.0A", "epsI": "atan2(Iq,Id)", # current angle "Im": "sqrt(Id^2+Iq^2)", - "InitPos": "-30deg", + "InitPos": f"{self.InitPos}deg", "w": f"{w}Hz", "RotSpeed": f"{RotSpeed}rpm", "Nper": "1/6", # number of included periods @@ -187,33 +189,17 @@ def set_oper_params(self): } def set_solution_expressions(self): - self.solution_expressions = "Moving1.Torque" + self.solution_expressions = [ + "Moving1.Position", + "Moving1.Torque", + *[f"FluxLinkage(Phase{p})" for p in "ABC"], + *[f"InducedVoltage(Phase{p})" for p in "ABC"], + *[f"InputCurrent(Phase{p})" for p in "ABC"], + *[f"L(Phase{x},Phase{y})" for x in "ABC" for y in "ABC"], + ] def set_output_vars(self): - self.output_vars = { - "pos": "(Moving1.Position -InitPos) * Poles/2", - "cos0": "cos(pos)", - "cos1": "cos(pos-2*PI/3)", - "cos2": "cos(pos-4*PI/3)", - "sin0": "sin(pos)", - "sin1": "sin(pos-2*PI/3)", - "sin2": "sin(pos-4*PI/3)", - "Lad": "L(PhaseA,PhaseA)*cos0 + L(PhaseA,PhaseB)*cos1 + L(PhaseA,PhaseC)*cos2", - "Laq": "L(PhaseA,PhaseA)*sin0 + L(PhaseA,PhaseB)*sin1 + L(PhaseA,PhaseC)*sin2", - "Lbd": "L(PhaseB,PhaseA)*cos0 + L(PhaseB,PhaseB)*cos1 + L(PhaseB,PhaseC)*cos2", - "Lbq": "L(PhaseB,PhaseA)*sin0 + L(PhaseB,PhaseB)*sin1 + L(PhaseB,PhaseC)*sin2", - "Lcd": "L(PhaseC,PhaseA)*cos0 + L(PhaseC,PhaseB)*cos1 + L(PhaseC,PhaseC)*cos2", - "Lcq": "L(PhaseC,PhaseA)*sin0 + L(PhaseC,PhaseB)*sin1 + L(PhaseC,PhaseC)*sin2", - "L_d": "(Lad*cos0 + Lbd*cos1 + Lcd*cos2) * 2/3", - "L_q": "(Laq*sin0 + Lbq*sin1 + Lcq*sin2) * 2/3", - "Flux_d": "(FluxLinkage(PhaseA)*cos0+FluxLinkage(PhaseB)*cos1+FluxLinkage(PhaseC)*cos2)*2/3", - "Flux_q": "-(FluxLinkage(PhaseA)*sin0+FluxLinkage(PhaseB)*sin1+FluxLinkage(PhaseC)*sin2)*2/3", - "Ui_d": "(InducedVoltage(PhaseA)*cos0+InducedVoltage(PhaseB)*cos1+InducedVoltage(PhaseC)*cos2)*2/3", - "Ui_q": "-(InducedVoltage(PhaseA)*sin0+InducedVoltage(PhaseB)*sin1+InducedVoltage(PhaseC)*sin2)*2/3", - "I_d": "(InputCurrent(PhaseA)*cos0 + InputCurrent(PhaseB)*cos1 + InputCurrent(PhaseC)*cos2)*2/3", - "I_q": "-(InputCurrent(PhaseA)*sin0 + InputCurrent(PhaseB)*sin1 + InputCurrent(PhaseC)*sin2)*2/3", - "Irms": "sqrt(I_d^2+I_q^2)/sqrt(2)", - } + self.output_vars = {} def set_post_params(self): self.post_params = { # reports @@ -266,8 +252,49 @@ def inductance_computation(self, m2d: Maxwell2d) -> None: m2d.change_inductance_computation(compute_transient_inductance=True, incremental_matrix=False) def set_variables(self, m2d: Maxwell2d, Id, Iq): + self.Id, self.Iq = Id, Iq m2d.variable_manager["Id"] = f"{Id}A" m2d.variable_manager["Iq"] = f"{Iq}A" def extract_results(self, solutions): - return solutions.data_magnitude() + position = np.array(solutions.data_real("Moving1.Position")) + torque = np.array(solutions.data_real("Moving1.Torque")) + + theta_el = np.deg2rad(position - self.InitPos) * self.geometry.PolePairs + flux_phases = np.stack([np.array(solutions.data_real(f"FluxLinkage(Phase{p})")) for p in "ABC"], axis=-1) + vind_phases = np.stack([np.array(solutions.data_real(f"InducedVoltage(Phase{p})")) for p in "ABC"], axis=-1) + current_phases = np.stack([np.array(solutions.data_real(f"InputCurrent(Phase{p})")) for p in "ABC"], axis=-1) + + Flux_d, Flux_q = to_dq(flux_phases, theta_el, harmonic=1) + Ui_d, Ui_q = to_dq(vind_phases, theta_el, harmonic=1) + I_d, I_q = (v / 1e3 for v in to_dq(current_phases, theta_el, harmonic=1)) + # mA to A + + L_raw = [np.stack([np.array(solutions.data_real(f"L(Phase{x},Phase{y})")) for y in "ABC"], axis=-1) / 1e9 for x in "ABC"] + # nH to H + + L_d_row = np.zeros((len(position), 3)) + L_q_row = np.zeros((len(position), 3)) + for i, L_row in enumerate(L_raw): + L_d_row[:, i], L_q_row[:, i] = (v * m for v, m in zip(to_dq(L_row, theta_el, harmonic=1), (3 / 2, -3 / 2))) + + L_d, _ = to_dq(L_d_row, theta_el, harmonic=1) + _, L_q_raw = to_dq(L_q_row, theta_el, harmonic=1) + L_q = -L_q_raw + + Irms = np.sqrt(I_d**2 + I_q**2) / np.sqrt(2) + + out = { + "Flux_d": Flux_d, + "Flux_q": Flux_q, + "Ui_d": Ui_d, + "Ui_q": Ui_q, + "I_d": I_d, + "I_q": I_q, + "L_d": L_d, + "L_q": L_q, + "Irms": Irms, + "Moving1.Torque": torque, + } + + return out diff --git a/tests/test_ansys_regression.py b/tests/test_ansys_regression.py index e910202..5bf8111 100644 --- a/tests/test_ansys_regression.py +++ b/tests/test_ansys_regression.py @@ -30,6 +30,18 @@ SEEDS = (42, 43) NUM_CORES = 1 current_setpoint = (1.5, 1.5) +check_expressions = [ + "Flux_d", + "Flux_q", + "Ui_d", + "Ui_q", + "I_d", + "I_q", + "L_d", + "L_q", + "Irms", +] +# we should verify all the above 9 vars + Moving1.Torque (a total of 10 vars with different seeds) def _load_legacy_design_class(): @@ -90,6 +102,20 @@ def test_legacy_and_live_design_match_across_rotor_rebuilds(tmp_path): legacy_design, live_design = None, None try: legacy_design = _create_design(LegacyDesign, "RegressionTest_legacy", str(tmp_path / "legacy.aedt")) + legacy_design.solution_expressions = ["Moving1.Torque", *check_expressions] + + def _legacy_extract_results(solutions): + out = {} + for expr in legacy_design.solution_expressions: + val = np.array(solutions.data_real(expr)) + if expr in ("I_d", "I_q"): + val = val / 1e3 # mA to A + elif expr in ("L_d", "L_q"): + val = val / 1e9 # nH to H + out[expr] = val + return out + + legacy_design.extract_results = _legacy_extract_results live_geometry = Geometry() live_computation = Computation(live_geometry) @@ -106,7 +132,10 @@ def test_legacy_and_live_design_match_across_rotor_rebuilds(tmp_path): assert torque_legacy is not None, f"legacy torque is None for seed {seed}" assert torque_live is not None, f"live torque is None for seed {seed}" - np.testing.assert_allclose(torque_live, torque_legacy, rtol=1e-6, err_msg=f"mismatch for seed {seed}") + np.testing.assert_allclose(torque_live["Moving1.Torque"], torque_legacy["Moving1.Torque"], rtol=1e-3, err_msg=f"mismatch for Moving1.Torque, seed {seed}") + for expr in check_expressions: + np.testing.assert_allclose(torque_live[expr], torque_legacy[expr], rtol=1e-3, err_msg=f"mismatch for {expr}, seed {seed}") + finally: if live_design is not None: live_design.close_project() diff --git a/tests/test_synrm_3f_36s.py b/tests/test_synrm_3f_36s.py index a05fa99..3f5fa68 100644 --- a/tests/test_synrm_3f_36s.py +++ b/tests/test_synrm_3f_36s.py @@ -1,9 +1,9 @@ -import re from types import SimpleNamespace import numpy as np import pytest +from machine_design.generic.transforms import to_dq from motors.synrm_3f_36s import Computation, Geometry @@ -95,8 +95,12 @@ def test_derived_params(geometry): assert geometry.rotor_r_max == pytest.approx(39.5 - 0.225) -def test_solution_expressions(computation): - assert computation.solution_expressions == "Moving1.Torque" +def test_solution_expressions_is_a_list_ending_in_torque(computation): + assert isinstance(computation.solution_expressions, list) + assert computation.solution_expressions[0] == "Moving1.Position" + assert computation.solution_expressions[1] == "Moving1.Torque" + assert computation.solution_expressions[-1] == "L(PhaseC,PhaseC)" + assert len(computation.solution_expressions) == 20 def test_udp_par_list_stator(geometry): @@ -124,20 +128,8 @@ def test_udp_par_list_stator(geometry): ] -def test_output_vars_only_reference_earlier_keys(computation): - # Ansys evaluates output variables in insertion order, so a formula referencing - # another output variable must not reference one that is defined later. - defined = set() - all_keys = set(computation.output_vars.keys()) - for key, formula in computation.output_vars.items(): - later_keys = all_keys - defined - {key} - referenced_later = [k for k in later_keys if re.search(rf"\b{re.escape(k)}\b", formula)] - assert not referenced_later, f"'{key}' formula references not-yet-defined {referenced_later}" - defined.add(key) - - -def test_output_vars_known_formula(computation): - assert computation.output_vars["Irms"] == "sqrt(I_d^2+I_q^2)/sqrt(2)" +def test_output_vars_is_empty(computation): + assert computation.output_vars == {} def test_post_params_plot_names_unique(computation): @@ -154,6 +146,20 @@ def test_mm_to_str_raises_without_mm_suffix(geometry): geometry.mm_to_str("geom_params", "SlotNumber") -def test_extract_results_returns_data_magnitude(computation): - solutions = SimpleNamespace(data_magnitude=lambda: [1.0, 2.0, 3.0]) - assert computation.extract_results(solutions) == [1.0, 2.0, 3.0] +def test_extract_results_converts_units(computation): + computation.set_variables(SimpleNamespace(variable_manager={}), Id=1.0, Iq=2.0) + + def data_real(expr): + letters = ["ABC".index(part[0]) + 1 for part in expr.split("Phase")[1:]] + value = float(np.prod(letters) if letters else 7.0) + return np.full(2, value * 1e9 if expr.startswith("L(") else value) + + solutions = SimpleNamespace(data_real=data_real, primary_sweep_values=[0.0, 1.0]) + out = computation.extract_results(solutions) + + theta_el = np.deg2rad(np.full(2, 7.0) - computation.InitPos) * computation.geometry.PolePairs + current_phases = np.stack([np.full(2, i + 1.0) for i in range(3)], axis=-1) + expected_I_d, _ = to_dq(current_phases, theta_el, harmonic=1) + + assert out["I_d"] == pytest.approx(expected_I_d / 1e3) + assert out["Moving1.Torque"] == pytest.approx(np.full(2, 7.0))