对齐test_mql保存时刻与PNVO事件边界

This commit is contained in:
huojiarong committed 2026-07-21 07:56:25 +00:00
1 parent 7c68caed27
commit a0855018d5
5 files changed
+345 -35

No files matched your search

@@ -1,7 +1,9 @@
from __future__ import annotations
import argparse
from dataclasses import dataclass, field
from datetime import UTC, datetime
from math import nextafter
from pathlib import Path
from PythonModels.core.solver import SolveIVPConfig
@@ -35,6 +37,8 @@ DEFAULT_FULL_STATE_COMPARISON_DATA_PATHS = (
"v1@mass_friction_endstops_18",
"acc1@mass_friction_endstops_18",
"dm1@pneumatic_69",
"xv@pn_morifice_1",
"dm2@pn_morifice_1",
)
@@ -68,6 +72,22 @@ class TestMqlFullStateFlowDiagnostic:
final_canonical_abs_error_kg_s: float
@dataclass(frozen=True)
class TestMqlFullStatePnvoDiagnostic:
alias: str
initial_time_s: float
final_time_s: float
initial_python_opening: float
initial_amesim_opening: float
final_python_opening: float
final_amesim_opening: float
initial_python_mass_flow_kg_s: float
initial_amesim_mass_flow_kg_s: float
final_python_mass_flow_kg_s: float
final_amesim_mass_flow_kg_s: float
final_mass_flow_abs_error_kg_s: float
@dataclass(frozen=True)
class TestMqlFullStateComparisonRun:
system: TestMqlSystem
@@ -159,6 +179,33 @@ class TestMqlFullStateComparisonRun:
),
)
def pnvo_diagnostic(
self,
alias: str = "pn_morifice_1",
) -> TestMqlFullStatePnvoDiagnostic:
opening = self.signal_diagnostic(f"xv@{alias}")
mass_flow = self.signal_diagnostic(f"dm2@{alias}")
initial_python_mass_flow = mass_flow.initial_python_value * 1.0e-3
initial_amesim_mass_flow = mass_flow.initial_amesim_value * 1.0e-3
final_python_mass_flow = mass_flow.final_python_value * 1.0e-3
final_amesim_mass_flow = mass_flow.final_amesim_value * 1.0e-3
return TestMqlFullStatePnvoDiagnostic(
alias=alias,
initial_time_s=opening.initial_time_s,
final_time_s=opening.final_time_s,
initial_python_opening=opening.initial_python_value,
initial_amesim_opening=opening.initial_amesim_value,
final_python_opening=opening.final_python_value,
final_amesim_opening=opening.final_amesim_value,
initial_python_mass_flow_kg_s=initial_python_mass_flow,
initial_amesim_mass_flow_kg_s=initial_amesim_mass_flow,
final_python_mass_flow_kg_s=final_python_mass_flow,
final_amesim_mass_flow_kg_s=final_amesim_mass_flow,
final_mass_flow_abs_error_kg_s=abs(
final_python_mass_flow - final_amesim_mass_flow
),
)
def diagnostics_by_final_abs_error(self):
return tuple(
sorted(
@@ -185,6 +232,22 @@ class TestMqlFullStateComparisonRun:
)
@dataclass(frozen=True)
class TestMqlPnvoEventBoundaryDiagnostic:
orifice_alias: str
event_time_s: float
integration_stop_time_s: float
data_paths: tuple[str, ...]
python_values_by_data_path: dict[str, float]
amesim_values_by_data_path: dict[str, float]
def abs_error(self, data_path: str) -> float:
return abs(
self.python_values_by_data_path[data_path]
- self.amesim_values_by_data_path[data_path]
)
@dataclass(frozen=True)
class TestMqlFullStateComparisonPathConfig:
archive_path: Path = field(
@@ -201,9 +264,9 @@ class TestMqlFullStateComparisonExecutionConfig:
write_comparison_csv: bool = True
data_paths: tuple[str, ...] | None = DEFAULT_FULL_STATE_COMPARISON_DATA_PATHS
solver: SolveIVPConfig = field(
default_factory=lambda: SolveIVPConfig(t_stop=1.0e-5, max_step=1.0e-6)
default_factory=lambda: SolveIVPConfig(t_stop=1.0e-2, max_step=1.0e-3)
)
t_eval: tuple[float, ...] | None = (0.0, 1.0e-5)
t_eval: tuple[float, ...] | None = (0.0, 1.0e-2)
inlet_node_pressure_pa: float = 15.31e6
resistance_boundary_pressure_pa: float = 15.29e6
inlet_node_temperature_k: float = 293.15
@@ -306,6 +369,94 @@ def run_test_mql_full_state_comparison(
return run, output_dir
def run_test_mql_pnvo_event_boundary_diagnostic(
config: TestMqlFullStateComparisonScriptConfig | None = None,
*,
orifice_alias: str = "pn_morifice_1",
) -> TestMqlPnvoEventBoundaryDiagnostic:
config = config or TestMqlFullStateComparisonScriptConfig()
system = TestMqlSystem(archive_path=config.paths.archive_path)
control = system.pneumatic_assembly.variable_orifice_controls[orifice_alias]
event_time_s = control.step.step_time_s
integration_stop_time_s = nextafter(event_time_s, 0.0)
solver_template = config.execution.solver
solver = SolveIVPConfig(
t_start=solver_template.t_start,
t_stop=integration_stop_time_s,
method=solver_template.method,
rtol=solver_template.rtol,
atol=solver_template.atol,
max_step=solver_template.max_step,
)
spec = system.discover_pneumatic_branch_topology().chamber_segment_specs[0]
closure_kwargs = {
"inlet_node_pressure_pa": config.execution.inlet_node_pressure_pa,
"resistance_boundary_pressure_pa": (
config.execution.resistance_boundary_pressure_pa
),
"inlet_node_temperature_k": config.execution.inlet_node_temperature_k,
"resistance_boundary_temperature_k": (
config.execution.resistance_boundary_temperature_k
),
}
solution = system.simulate_full_state_from_spec(
spec,
config=solver,
t_eval=[solver.t_start, integration_stop_time_s],
**closure_kwargs,
)
state_vector = [row[-1] for row in solution.y]
closure = system.full_state_closure_from_spec(spec, **closure_kwargs)
data_paths = (
"press@pn_c1_8",
"dm1@pneumatic_69",
f"xv@{orifice_alias}",
f"dm2@{orifice_alias}",
)
python_values = closure.data_path_values(
time_s=event_time_s,
state_vector=state_vector,
data_paths=data_paths,
)
amesim_results = load_test_mql_amesim_results(config.paths.archive_path)
amesim_values = {
data_path: interpolate_series_value(
amesim_results.times,
amesim_results.series(data_path),
event_time_s,
)
for data_path in data_paths
}
return TestMqlPnvoEventBoundaryDiagnostic(
orifice_alias=orifice_alias,
event_time_s=event_time_s,
integration_stop_time_s=integration_stop_time_s,
data_paths=data_paths,
python_values_by_data_path=python_values,
amesim_values_by_data_path=amesim_values,
)
def format_test_mql_pnvo_event_boundary_summary(
diagnostic: TestMqlPnvoEventBoundaryDiagnostic,
) -> str:
lines = [
"Model: test_mql",
f"Mode: PNVO event boundary diagnostic ({diagnostic.orifice_alias})",
f"Event time: {diagnostic.event_time_s}",
f"Integrated left limit: {diagnostic.integration_stop_time_s}",
"Observation side: right-continuous STEP0 opening",
]
for data_path in diagnostic.data_paths:
lines.append(
f" - {data_path}: "
f"python={diagnostic.python_values_by_data_path[data_path]}, "
f"amesim={diagnostic.amesim_values_by_data_path[data_path]}, "
f"abs_error={diagnostic.abs_error(data_path)}"
)
return "\n".join(lines) + "\n"
def format_test_mql_full_state_comparison_summary(
run: TestMqlFullStateComparisonRun,
) -> str:
@@ -315,6 +466,10 @@ def format_test_mql_full_state_comparison_summary(
f"Samples: {run.sample_count}",
f"Compared signals: {run.signal_count}",
f"Output schema signals: {run.output_schema.signal_count}",
f"AMESim first saved interval: "
f"{run.amesim_results.times[1] - run.amesim_results.times[0]}",
f"Final time aligns with AMESim sample: "
f"{_matches_amesim_sample(run, run.output.times[-1])}",
f"Max absolute error: {run.comparison.max_abs_error}",
f"Max relative error: {run.comparison.max_rel_error}",
]
@@ -368,6 +523,27 @@ def format_test_mql_full_state_comparison_summary(
f"{flow_diagnostic.final_canonical_abs_error_kg_s}",
]
)
pnvo_diagnostic = _pn_morifice_1_diagnostic(run)
if pnvo_diagnostic is not None:
lines.append(f"PNVO diagnostic: {pnvo_diagnostic.alias}")
lines.extend(
[
f" - initial_python_opening="
f"{pnvo_diagnostic.initial_python_opening}",
f" - initial_amesim_opening="
f"{pnvo_diagnostic.initial_amesim_opening}",
f" - final_python_opening="
f"{pnvo_diagnostic.final_python_opening}",
f" - final_amesim_opening="
f"{pnvo_diagnostic.final_amesim_opening}",
f" - final_python_mass_flow_kg_s="
f"{pnvo_diagnostic.final_python_mass_flow_kg_s}",
f" - final_amesim_mass_flow_kg_s="
f"{pnvo_diagnostic.final_amesim_mass_flow_kg_s}",
f" - final_mass_flow_abs_error_kg_s="
f"{pnvo_diagnostic.final_mass_flow_abs_error_kg_s}",
]
)
chamber_diagnostic = _largest_chamber_rhs_diagnostic(run)
if chamber_diagnostic is not None:
lines.append(
@@ -389,6 +565,16 @@ def format_test_mql_full_state_comparison_summary(
return "\n".join(lines) + "\n"
def _matches_amesim_sample(
run: TestMqlFullStateComparisonRun,
time_s: float,
) -> bool:
return any(
abs(float(sample_time) - float(time_s)) <= 1.0e-12
for sample_time in run.amesim_results.times
)
def _pneumatic_69_flow_diagnostic(run: TestMqlFullStateComparisonRun):
try:
return run.pnl0001_mass_flow_diagnostic()
@@ -396,6 +582,13 @@ def _pneumatic_69_flow_diagnostic(run: TestMqlFullStateComparisonRun):
return None
def _pn_morifice_1_diagnostic(run: TestMqlFullStateComparisonRun):
try:
return run.pnvo_diagnostic()
except KeyError:
return None
def _largest_chamber_rhs_diagnostic(run: TestMqlFullStateComparisonRun):
diagnostic = run.largest_final_abs_error_diagnostic
if diagnostic is None or "@" not in diagnostic.data_path:
@@ -408,6 +601,18 @@ def _largest_chamber_rhs_diagnostic(run: TestMqlFullStateComparisonRun):
def main() -> None:
parser = argparse.ArgumentParser()
parser.add_argument(
"--pnvo-event-boundary",
action="store_true",
help="compare the STEP0 left-state/right-opening boundary at t=0.04 s",
)
args = parser.parse_args()
if args.pnvo_event_boundary:
diagnostic = run_test_mql_pnvo_event_boundary_diagnostic()
print(format_test_mql_pnvo_event_boundary_summary(diagnostic), end="")
return
run, output_dir = run_test_mql_full_state_comparison()
print(format_test_mql_full_state_comparison_summary(run), end="")
print(f"Output directory: {output_dir}")
+27 -1
View File
@@ -3909,7 +3909,14 @@ def _default_full_state_data_paths() -> tuple[str, ...]:
for signal in ("x1", "v1", "acc1")
]
return tuple(
[*chamber_paths, *piston_paths, *mass_paths, "dm1@pneumatic_69"]
[
*chamber_paths,
*piston_paths,
*mass_paths,
"dm1@pneumatic_69",
"xv@pn_morifice_1",
"dm2@pn_morifice_1",
]
)
@@ -3973,6 +3980,10 @@ _PNL0001_MASS_FLOW_DIAGNOSTIC_FIELDS = {
"pneumatic_69": "p4_port3_remote_chamber_to_line_flow",
}
_VARIABLE_ORIFICE_MASS_FLOW_DIAGNOSTIC_FIELDS = {
"pn_morifice_1": "p4_port3_remote_orifice_to_node_flow",
}
class TestMqlFullStateClosure:
def __init__(
@@ -4162,6 +4173,21 @@ class TestMqlFullStateClosure:
_PNL0001_MASS_FLOW_DIAGNOSTIC_FIELDS[alias],
)
values_by_data_path[data_path] = -canonical_flow_kg_s * 1.0e3
elif alias in self.pneumatic_assembly.variable_orifices:
orifice = self.pneumatic_assembly.variable_orifices[alias]
if signal == "xv":
values_by_data_path[data_path] = orifice.opening
elif (
signal == "dm2"
and alias in _VARIABLE_ORIFICE_MASS_FLOW_DIAGNOSTIC_FIELDS
):
flow_kg_s = getattr(
snapshot.pneumatic,
_VARIABLE_ORIFICE_MASS_FLOW_DIAGNOSTIC_FIELDS[alias],
)
values_by_data_path[data_path] = flow_kg_s * 1.0e3
else:
raise KeyError(data_path)
elif alias in chamber_properties_by_alias:
chamber_properties = chamber_properties_by_alias[alias]
chamber = chamber_by_alias[alias]