PythonModels/scripts/run_test_mql_full_state_comparison.py

@@ -498,6 +498,22 @@ class TestMqlPnl0003EnergyDiagnostic:
     amesim_port_2_pn2vol_dtemp_k_s: float
     amesim_port_1_pn2vol_dtemp_residual_k_s: float
     amesim_port_2_pn2vol_dtemp_residual_k_s: float
+    amesim_port_1_pressure_derivative_fd_pa_s: float
+    amesim_port_2_pressure_derivative_fd_pa_s: float
+    amesim_port_1_pressure_derivative_mass_component_pa_s: float
+    amesim_port_1_pressure_derivative_temperature_component_pa_s: float
+    amesim_port_1_pressure_derivative_eos_sum_pa_s: float
+    amesim_port_1_pressure_derivative_eos_residual_pa_s: float
+    amesim_port_2_pressure_derivative_mass_component_pa_s: float
+    amesim_port_2_pressure_derivative_temperature_component_pa_s: float
+    amesim_port_2_pressure_derivative_eos_sum_pa_s: float
+    amesim_port_2_pressure_derivative_eos_residual_pa_s: float
+    python_port_1_pressure_derivative_mass_component_pa_s: float
+    python_port_1_pressure_derivative_temperature_component_pa_s: float
+    python_port_1_pressure_derivative_eos_sum_pa_s: float
+    python_port_2_pressure_derivative_mass_component_pa_s: float
+    python_port_2_pressure_derivative_temperature_component_pa_s: float
+    python_port_2_pressure_derivative_eos_sum_pa_s: float
     python_port_1_temperature_k: float
     python_port_2_temperature_k: float
     python_port_1_pressure_pa: float
@@ -1634,6 +1650,31 @@ def _actual_reference_enthalpy_flow(
     )
 
 
+def _pressure_derivative_components_from_density_temperature(
+    *,
+    gas: AmesimPneumaticGas,
+    density_kg_m3: float,
+    temperature_k: float,
+    density_derivative_kg_m3_s: float,
+    temperature_derivative_k_s: float,
+) -> tuple[float, float, float]:
+    density_step = max(abs(density_kg_m3) * 1.0e-6, 1.0e-9)
+    if density_kg_m3 <= density_step:
+        density_step = density_kg_m3 * 0.5
+    temperature_step = max(abs(temperature_k) * 1.0e-6, 1.0e-4)
+    dp_drho = (
+        gas.pressure(density_kg_m3 + density_step, temperature_k)
+        - gas.pressure(density_kg_m3 - density_step, temperature_k)
+    ) / (2.0 * density_step)
+    dp_dtemp = (
+        gas.pressure(density_kg_m3, temperature_k + temperature_step)
+        - gas.pressure(density_kg_m3, temperature_k - temperature_step)
+    ) / (2.0 * temperature_step)
+    mass_component = dp_drho * density_derivative_kg_m3_s
+    temperature_component = dp_dtemp * temperature_derivative_k_s
+    return mass_component, temperature_component, mass_component + temperature_component
+
+
 def _pnl0003_energy_diagnostic(
     *,
     closure: object,
@@ -1811,6 +1852,40 @@ def _pnl0003_energy_diagnostic(
             time_s=time_s,
         )
     )
+    amesim_p1_series = amesim_results.series(f"p1@{line_alias}")
+    amesim_p2_series = amesim_results.series(f"p2@{line_alias}")
+    amesim_fd_dpressure_1 = _series_finite_difference_at(
+        times=amesim_results.times,
+        values=amesim_p1_series,
+        time_s=time_s,
+    )
+    amesim_fd_dpressure_2 = _series_finite_difference_at(
+        times=amesim_results.times,
+        values=amesim_p2_series,
+        time_s=time_s,
+    )
+    (
+        amesim_dpressure_1_mass,
+        amesim_dpressure_1_temperature,
+        amesim_dpressure_1_sum,
+    ) = _pressure_derivative_components_from_density_temperature(
+        gas=line.gas,
+        density_kg_m3=amesim_mass_1 / line.compliance_volume,
+        temperature_k=amesim_t1,
+        density_derivative_kg_m3_s=amesim_sdm1_kg_s / line.compliance_volume,
+        temperature_derivative_k_s=amesim_fd_dtemp_1,
+    )
+    (
+        amesim_dpressure_2_mass,
+        amesim_dpressure_2_temperature,
+        amesim_dpressure_2_sum,
+    ) = _pressure_derivative_components_from_density_temperature(
+        gas=line.gas,
+        density_kg_m3=amesim_mass_2 / line.compliance_volume,
+        temperature_k=amesim_t2,
+        density_derivative_kg_m3_s=amesim_sdm2_kg_s / line.compliance_volume,
+        temperature_derivative_k_s=amesim_fd_dtemp_2,
+    )
 
     current_derivative_1, current_derivative_2 = line.derivatives_from_connections(
         port_1_m_flow=port_1_flow,
@@ -1824,6 +1899,28 @@ def _pnl0003_energy_diagnostic(
     python_current_dtemp_2 = (
         current_derivative_2.U - port_2.u * current_derivative_2.m
     ) / (line.state_2.m * line.gas.cv)
+    (
+        python_dpressure_1_mass,
+        python_dpressure_1_temperature,
+        python_dpressure_1_sum,
+    ) = _pressure_derivative_components_from_density_temperature(
+        gas=line.gas,
+        density_kg_m3=line.state_1.m / line.compliance_volume,
+        temperature_k=port_1.T,
+        density_derivative_kg_m3_s=current_derivative_1.m / line.compliance_volume,
+        temperature_derivative_k_s=python_current_dtemp_1,
+    )
+    (
+        python_dpressure_2_mass,
+        python_dpressure_2_temperature,
+        python_dpressure_2_sum,
+    ) = _pressure_derivative_components_from_density_temperature(
+        gas=line.gas,
+        density_kg_m3=line.state_2.m / line.compliance_volume,
+        temperature_k=port_2.T,
+        density_derivative_kg_m3_s=current_derivative_2.m / line.compliance_volume,
+        temperature_derivative_k_s=python_current_dtemp_2,
+    )
 
     python_center_dh_1 = _actual_reference_enthalpy_flow(
         gas=line.gas,
@@ -2186,6 +2283,42 @@ def _pnl0003_energy_diagnostic(
         amesim_port_2_pn2vol_dtemp_residual_k_s=(
             amesim_pn2vol_dtemp_2 - amesim_fd_dtemp_2
         ),
+        amesim_port_1_pressure_derivative_fd_pa_s=amesim_fd_dpressure_1,
+        amesim_port_2_pressure_derivative_fd_pa_s=amesim_fd_dpressure_2,
+        amesim_port_1_pressure_derivative_mass_component_pa_s=(
+            amesim_dpressure_1_mass
+        ),
+        amesim_port_1_pressure_derivative_temperature_component_pa_s=(
+            amesim_dpressure_1_temperature
+        ),
+        amesim_port_1_pressure_derivative_eos_sum_pa_s=amesim_dpressure_1_sum,
+        amesim_port_1_pressure_derivative_eos_residual_pa_s=(
+            amesim_dpressure_1_sum - amesim_fd_dpressure_1
+        ),
+        amesim_port_2_pressure_derivative_mass_component_pa_s=(
+            amesim_dpressure_2_mass
+        ),
+        amesim_port_2_pressure_derivative_temperature_component_pa_s=(
+            amesim_dpressure_2_temperature
+        ),
+        amesim_port_2_pressure_derivative_eos_sum_pa_s=amesim_dpressure_2_sum,
+        amesim_port_2_pressure_derivative_eos_residual_pa_s=(
+            amesim_dpressure_2_sum - amesim_fd_dpressure_2
+        ),
+        python_port_1_pressure_derivative_mass_component_pa_s=(
+            python_dpressure_1_mass
+        ),
+        python_port_1_pressure_derivative_temperature_component_pa_s=(
+            python_dpressure_1_temperature
+        ),
+        python_port_1_pressure_derivative_eos_sum_pa_s=python_dpressure_1_sum,
+        python_port_2_pressure_derivative_mass_component_pa_s=(
+            python_dpressure_2_mass
+        ),
+        python_port_2_pressure_derivative_temperature_component_pa_s=(
+            python_dpressure_2_temperature
+        ),
+        python_port_2_pressure_derivative_eos_sum_pa_s=python_dpressure_2_sum,
         python_port_1_temperature_k=port_1.T,
         python_port_2_temperature_k=port_2.T,
         python_port_1_pressure_pa=port_1.p,
@@ -2916,6 +3049,22 @@ def format_test_mql_pnvo_event_window_summary(
                     f"amesim_pn2vol_dT2={pnl0003_energy.amesim_port_2_pn2vol_dtemp_k_s}, "
                     f"amesim_pn2vol_dT1_residual={pnl0003_energy.amesim_port_1_pn2vol_dtemp_residual_k_s}, "
                     f"amesim_pn2vol_dT2_residual={pnl0003_energy.amesim_port_2_pn2vol_dtemp_residual_k_s}, "
+                    f"amesim_fd_dP1={pnl0003_energy.amesim_port_1_pressure_derivative_fd_pa_s}, "
+                    f"amesim_fd_dP2={pnl0003_energy.amesim_port_2_pressure_derivative_fd_pa_s}, "
+                    f"amesim_dP1_m={pnl0003_energy.amesim_port_1_pressure_derivative_mass_component_pa_s}, "
+                    f"amesim_dP1_T={pnl0003_energy.amesim_port_1_pressure_derivative_temperature_component_pa_s}, "
+                    f"amesim_dP1_sum={pnl0003_energy.amesim_port_1_pressure_derivative_eos_sum_pa_s}, "
+                    f"amesim_dP1_residual={pnl0003_energy.amesim_port_1_pressure_derivative_eos_residual_pa_s}, "
+                    f"amesim_dP2_m={pnl0003_energy.amesim_port_2_pressure_derivative_mass_component_pa_s}, "
+                    f"amesim_dP2_T={pnl0003_energy.amesim_port_2_pressure_derivative_temperature_component_pa_s}, "
+                    f"amesim_dP2_sum={pnl0003_energy.amesim_port_2_pressure_derivative_eos_sum_pa_s}, "
+                    f"amesim_dP2_residual={pnl0003_energy.amesim_port_2_pressure_derivative_eos_residual_pa_s}, "
+                    f"python_dP1_m={pnl0003_energy.python_port_1_pressure_derivative_mass_component_pa_s}, "
+                    f"python_dP1_T={pnl0003_energy.python_port_1_pressure_derivative_temperature_component_pa_s}, "
+                    f"python_dP1_sum={pnl0003_energy.python_port_1_pressure_derivative_eos_sum_pa_s}, "
+                    f"python_dP2_m={pnl0003_energy.python_port_2_pressure_derivative_mass_component_pa_s}, "
+                    f"python_dP2_T={pnl0003_energy.python_port_2_pressure_derivative_temperature_component_pa_s}, "
+                    f"python_dP2_sum={pnl0003_energy.python_port_2_pressure_derivative_eos_sum_pa_s}, "
                     f"python_t1={pnl0003_energy.python_port_1_temperature_k}, "
                     f"python_t2={pnl0003_energy.python_port_2_temperature_k}, "
                     f"python_p1_abs={pnl0003_energy.python_port_1_pressure_pa}, "

tests/test_run_test_mql_full_state_comparison.py

@@ -443,6 +443,22 @@ class RunTestMqlFullStateComparisonScriptTests(unittest.TestCase):
                     amesim_port_2_pn2vol_dtemp_k_s=-530.0,
                     amesim_port_1_pn2vol_dtemp_residual_k_s=-10.0,
                     amesim_port_2_pn2vol_dtemp_residual_k_s=10.0,
+                    amesim_port_1_pressure_derivative_fd_pa_s=-2100000.0,
+                    amesim_port_2_pressure_derivative_fd_pa_s=-5400000.0,
+                    amesim_port_1_pressure_derivative_mass_component_pa_s=-3100000.0,
+                    amesim_port_1_pressure_derivative_temperature_component_pa_s=900000.0,
+                    amesim_port_1_pressure_derivative_eos_sum_pa_s=-2200000.0,
+                    amesim_port_1_pressure_derivative_eos_residual_pa_s=-100000.0,
+                    amesim_port_2_pressure_derivative_mass_component_pa_s=-6100000.0,
+                    amesim_port_2_pressure_derivative_temperature_component_pa_s=800000.0,
+                    amesim_port_2_pressure_derivative_eos_sum_pa_s=-5300000.0,
+                    amesim_port_2_pressure_derivative_eos_residual_pa_s=100000.0,
+                    python_port_1_pressure_derivative_mass_component_pa_s=-3000000.0,
+                    python_port_1_pressure_derivative_temperature_component_pa_s=950000.0,
+                    python_port_1_pressure_derivative_eos_sum_pa_s=-2050000.0,
+                    python_port_2_pressure_derivative_mass_component_pa_s=-6000000.0,
+                    python_port_2_pressure_derivative_temperature_component_pa_s=850000.0,
+                    python_port_2_pressure_derivative_eos_sum_pa_s=-5150000.0,
                     python_port_1_temperature_k=293.15,
                     python_port_2_temperature_k=282.2,
                     python_port_1_pressure_pa=15.3e6,
@@ -659,6 +675,13 @@ class RunTestMqlFullStateComparisonScriptTests(unittest.TestCase):
         self.assertIn("amesim_back_dT2=-600.0", summary)
         self.assertIn("amesim_fwd_dT2=120.0", summary)
         self.assertIn("amesim_pn2vol_dT2=-530.0", summary)
+        self.assertIn("amesim_fd_dP1=-2100000.0", summary)
+        self.assertIn("amesim_dP1_m=-3100000.0", summary)
+        self.assertIn("amesim_dP1_T=900000.0", summary)
+        self.assertIn("amesim_dP1_residual=-100000.0", summary)
+        self.assertIn("amesim_dP2_sum=-5300000.0", summary)
+        self.assertIn("python_dP1_sum=-2050000.0", summary)
+        self.assertIn("python_dP2_T=850000.0", summary)
         self.assertIn("python_real_ref_dT2=-650.0", summary)
         self.assertIn("python_amesim_center_dT2=-610.0", summary)
         self.assertIn("python_amesim_ext_dT2=-620.0", summary)