新增test_mql气动组件原语

AmesimModels/test_mql/README.md

@@ -42,6 +42,11 @@
   - 当前覆盖压缩因子、摩尔体积、密度和由密度反算压力。
   - 已内置 `HELIUM_PR`，供 `test_mql` 默认使用。
 
+- `PythonModels/components/amesim_pneumatic.py`
+  - 提供 `test_mql` 后续会用到的 AMESim 气动组件原语。
+  - 当前包含氦气 PR 压力闭合的气动容腔、标准可压缩孔口流量、单位换算 helper。
+  - 这层仍是首版近似原语，后续必须通过 AMESim CSV 对齐修正。
+
 ## 物性约定
 
 AMESim 模型中 `test_mql` 使用氦气，Python 侧当前通过 `HELIUM_PR` 使用 Peng-Robinson 状态方程计算气体压缩因子和密度。当前物性层先覆盖状态方程相关量，完整焓/内能偏差函数后续在接气室能量方程时再补。
@@ -88,7 +93,7 @@ AMESim 模型中 `test_mql` 使用氦气，Python 侧当前通过 `HELIUM_PR`
 ## 验证方式
 
 ```bash
-python3 -m py_compile PythonModels/core/peng_robinson.py PythonModels/systems/test_mql.py PythonModels/systems/test_mql_config.py PythonModels/scripts/run_test_mql.py
+python3 -m py_compile PythonModels/components/amesim_pneumatic.py PythonModels/core/peng_robinson.py PythonModels/systems/test_mql.py PythonModels/systems/test_mql_config.py PythonModels/scripts/run_test_mql.py
 python3 -m PythonModels.scripts.run_test_mql
 python3 -m unittest discover -s tests -t .
 ```

PythonModels/components/amesim_pneumatic.py

@@ -0,0 +1,216 @@
+from __future__ import annotations
+
+from dataclasses import dataclass
+from math import pi, sqrt
+
+from PythonModels.core.base import AlgebraicComponent, DynamicComponent
+from PythonModels.core.medium import ThermodynamicProperties
+from PythonModels.core.peng_robinson import HELIUM_PR, PengRobinsonFluid
+from PythonModels.core.ports import PortState
+from PythonModels.core.state import VolumeState
+
+
+@dataclass(frozen=True)
+class AmesimPneumaticGas:
+    """Caloric constants plus Peng-Robinson EOS for AMESim pneumatic components."""
+
+    fluid: PengRobinsonFluid = HELIUM_PR
+    cp: float = 5193.0
+    cv: float = 3116.0
+
+    @property
+    def gamma(self) -> float:
+        return self.cp / self.cv
+
+    @property
+    def R_gas(self) -> float:
+        return self.fluid.specific_gas_constant
+
+    def density(self, pressure: float, temperature: float) -> float:
+        return self.fluid.density(pressure, temperature)
+
+    def pressure(self, density: float, temperature: float) -> float:
+        return self.fluid.pressure_from_density(temperature, density)
+
+    def specific_internal_energy(self, temperature: float) -> float:
+        return self.cv * temperature
+
+    def specific_enthalpy(self, temperature: float) -> float:
+        return self.cp * temperature
+
+    def temperature_from_internal_energy(self, specific_internal_energy: float) -> float:
+        if self.cv <= 0.0:
+            raise ValueError("cv must be positive.")
+        return specific_internal_energy / self.cv
+
+
+HELIUM_PNEUMATIC_GAS = AmesimPneumaticGas()
+
+
+def liters_to_m3(value: float) -> float:
+    return value * 1.0e-3
+
+
+def mm2_to_m2(value: float) -> float:
+    return value * 1.0e-6
+
+
+def diameter_mm_to_area_m2(diameter_mm: float) -> float:
+    diameter_m = diameter_mm * 1.0e-3
+    return pi * diameter_m * diameter_m / 4.0
+
+
+class AmesimPneumaticVolume(DynamicComponent):
+    """First-pass AMESim pneumatic control volume using helium PR pressure closure."""
+
+    def __init__(
+        self,
+        name: str,
+        volume: float,
+        gas: AmesimPneumaticGas = HELIUM_PNEUMATIC_GAS,
+        p0: float = 101_325.0,
+        T0: float = 293.15,
+    ) -> None:
+        if volume <= 0.0:
+            raise ValueError("volume must be positive.")
+        super().__init__(name=name)
+        self.volume = volume
+        self.gas = gas
+        rho0 = gas.density(p0, T0)
+        m0 = rho0 * volume
+        U0 = m0 * gas.specific_internal_energy(T0)
+        self.state = VolumeState(m=m0, U=U0)
+        self.port_a = PortState()
+
+    @classmethod
+    def from_liters(
+        cls,
+        name: str,
+        volume_liters: float,
+        gas: AmesimPneumaticGas = HELIUM_PNEUMATIC_GAS,
+        p0: float = 101_325.0,
+        T0: float = 293.15,
+    ) -> "AmesimPneumaticVolume":
+        return cls(name=name, volume=liters_to_m3(volume_liters), gas=gas, p0=p0, T0=T0)
+
+    def get_state_vector(self) -> list[float]:
+        return self.state.as_vector()
+
+    def set_state_vector(self, values: list[float]) -> None:
+        self.state = VolumeState.from_vector(values)
+
+    def properties(self) -> ThermodynamicProperties:
+        if self.state.m <= 0.0:
+            raise ValueError("volume mass must stay positive.")
+        T = self.gas.temperature_from_internal_energy(self.state.U / self.state.m)
+        rho = self.state.m / self.volume
+        p = self.gas.pressure(rho, T)
+        u = self.state.U / self.state.m
+        h = self.gas.specific_enthalpy(T)
+        self.port_a.p = p
+        self.port_a.h_outflow = h
+        return ThermodynamicProperties(p=p, T=T, rho=rho, u=u, h=h)
+
+    def derivatives(self, inlet_h: float, m_flow: float) -> VolumeState:
+        return VolumeState(m=m_flow, U=m_flow * inlet_h)
+
+
+class AmesimPneumaticOrifice(AlgebraicComponent):
+    """First-pass PNOR001/PNVO001-style compressible helium orifice.
+
+    This is a calibrated placeholder boundary for the Python port. It preserves
+    AMESim-style area and coefficient inputs, but final parity must be checked
+    against AMESim CSV results before treating it as numerically equivalent.
+    """
+
+    def __init__(
+        self,
+        name: str,
+        area: float,
+        flow_coefficient: float = 1.0,
+        gas: AmesimPneumaticGas = HELIUM_PNEUMATIC_GAS,
+        opening: float = 1.0,
+    ) -> None:
+        if area < 0.0:
+            raise ValueError("area must be non-negative.")
+        if flow_coefficient < 0.0:
+            raise ValueError("flow_coefficient must be non-negative.")
+        super().__init__(name=name)
+        self.area = area
+        self.flow_coefficient = flow_coefficient
+        self.gas = gas
+        self.opening = opening
+        self.port_a = PortState()
+        self.port_b = PortState()
+
+    @classmethod
+    def from_mm2(
+        cls,
+        name: str,
+        area_mm2: float,
+        flow_coefficient: float = 1.0,
+        gas: AmesimPneumaticGas = HELIUM_PNEUMATIC_GAS,
+        opening: float = 1.0,
+    ) -> "AmesimPneumaticOrifice":
+        return cls(
+            name=name,
+            area=mm2_to_m2(area_mm2),
+            flow_coefficient=flow_coefficient,
+            gas=gas,
+            opening=opening,
+        )
+
+    @property
+    def effective_area(self) -> float:
+        return self.area * max(self.opening, 0.0)
+
+    def mass_flow(self, p_a: float, p_b: float, upstream_temperature: float) -> float:
+        if p_a == p_b or self.effective_area == 0.0 or self.flow_coefficient == 0.0:
+            return 0.0
+        if p_a > p_b:
+            return compressible_orifice_mass_flow(
+                upstream_pressure=p_a,
+                downstream_pressure=p_b,
+                upstream_temperature=upstream_temperature,
+                area=self.effective_area,
+                flow_coefficient=self.flow_coefficient,
+                gas=self.gas,
+            )
+        return -compressible_orifice_mass_flow(
+            upstream_pressure=p_b,
+            downstream_pressure=p_a,
+            upstream_temperature=upstream_temperature,
+            area=self.effective_area,
+            flow_coefficient=self.flow_coefficient,
+            gas=self.gas,
+        )
+
+
+def compressible_orifice_mass_flow(
+    *,
+    upstream_pressure: float,
+    downstream_pressure: float,
+    upstream_temperature: float,
+    area: float,
+    flow_coefficient: float,
+    gas: AmesimPneumaticGas = HELIUM_PNEUMATIC_GAS,
+) -> float:
+    if upstream_pressure <= 0.0 or downstream_pressure < 0.0:
+        raise ValueError("pressures must be non-negative and upstream pressure must be positive.")
+    if upstream_temperature <= 0.0:
+        raise ValueError("upstream_temperature must be positive.")
+    if area < 0.0 or flow_coefficient < 0.0:
+        raise ValueError("area and flow_coefficient must be non-negative.")
+    if downstream_pressure >= upstream_pressure or area == 0.0 or flow_coefficient == 0.0:
+        return 0.0
+
+    gamma = gas.gamma
+    pressure_ratio = max(downstream_pressure / upstream_pressure, 0.0)
+    critical_ratio = (2.0 / (gamma + 1.0)) ** (gamma / (gamma - 1.0))
+    coefficient = flow_coefficient * area * upstream_pressure / sqrt(gas.R_gas * upstream_temperature)
+    if pressure_ratio <= critical_ratio:
+        flow_function = sqrt(gamma) * (2.0 / (gamma + 1.0)) ** ((gamma + 1.0) / (2.0 * (gamma - 1.0)))
+    else:
+        term = pressure_ratio ** (2.0 / gamma) - pressure_ratio ** ((gamma + 1.0) / gamma)
+        flow_function = sqrt((2.0 * gamma / (gamma - 1.0)) * max(term, 0.0))
+    return coefficient * flow_function

tests/test_amesim_pneumatic_components.py

@@ -0,0 +1,92 @@
+from __future__ import annotations
+
+import unittest
+
+from PythonModels.components.amesim_pneumatic import (
+    HELIUM_PNEUMATIC_GAS,
+    AmesimPneumaticOrifice,
+    AmesimPneumaticVolume,
+    compressible_orifice_mass_flow,
+    diameter_mm_to_area_m2,
+    liters_to_m3,
+    mm2_to_m2,
+)
+
+
+class AmesimPneumaticComponentsTest(unittest.TestCase):
+    def test_unit_conversions(self) -> None:
+        self.assertAlmostEqual(liters_to_m3(15.0), 0.015)
+        self.assertAlmostEqual(mm2_to_m2(78.5), 78.5e-6)
+        self.assertAlmostEqual(diameter_mm_to_area_m2(10.0), 7.853981633974483e-5)
+
+    def test_volume_initial_state_matches_requested_pressure_temperature(self) -> None:
+        volume = AmesimPneumaticVolume.from_liters(
+            name="pn_general_chamber",
+            volume_liters=57.0,
+            p0=15.3e6,
+            T0=293.15,
+        )
+        props = volume.properties()
+
+        self.assertAlmostEqual(props.p, 15.3e6, delta=15.3e6 * 1.0e-12)
+        self.assertAlmostEqual(props.T, 293.15)
+        self.assertGreater(props.rho, 20.0)
+        self.assertLess(props.rho, 30.0)
+
+    def test_orifice_returns_signed_mass_flow(self) -> None:
+        orifice = AmesimPneumaticOrifice.from_mm2(
+            name="pn_orifice_18",
+            area_mm2=78.5,
+            flow_coefficient=0.9,
+        )
+
+        forward = orifice.mass_flow(15.3e6, 1.0e6, 293.15)
+        reverse = orifice.mass_flow(1.0e6, 15.3e6, 293.15)
+
+        self.assertGreater(forward, 0.0)
+        self.assertAlmostEqual(reverse, -forward)
+        self.assertEqual(orifice.mass_flow(1.0e6, 1.0e6, 293.15), 0.0)
+
+    def test_choked_flow_is_independent_of_lower_downstream_pressure(self) -> None:
+        base = compressible_orifice_mass_flow(
+            upstream_pressure=15.3e6,
+            downstream_pressure=1.0e6,
+            upstream_temperature=293.15,
+            area=78.5e-6,
+            flow_coefficient=0.9,
+            gas=HELIUM_PNEUMATIC_GAS,
+        )
+        lower_back_pressure = compressible_orifice_mass_flow(
+            upstream_pressure=15.3e6,
+            downstream_pressure=0.1e6,
+            upstream_temperature=293.15,
+            area=78.5e-6,
+            flow_coefficient=0.9,
+            gas=HELIUM_PNEUMATIC_GAS,
+        )
+
+        self.assertGreater(base, 0.0)
+        self.assertAlmostEqual(lower_back_pressure, base)
+
+    def test_subcritical_flow_decreases_as_back_pressure_rises(self) -> None:
+        low_back_pressure = compressible_orifice_mass_flow(
+            upstream_pressure=1.0e6,
+            downstream_pressure=0.6e6,
+            upstream_temperature=293.15,
+            area=78.5e-6,
+            flow_coefficient=0.9,
+        )
+        high_back_pressure = compressible_orifice_mass_flow(
+            upstream_pressure=1.0e6,
+            downstream_pressure=0.9e6,
+            upstream_temperature=293.15,
+            area=78.5e-6,
+            flow_coefficient=0.9,
+        )
+
+        self.assertGreater(low_back_pressure, high_back_pressure)
+        self.assertGreater(high_back_pressure, 0.0)
+
+
+if __name__ == "__main__":
+    unittest.main()