Skip to content

core.mass_balance

Right-hand side of the chamber-pressure ODE solved during internal-ballistics simulation. A single control-volume mass balance covers both choked and sub-critical nozzle flow (Seidel 1965, Eq. 35), and is shared by the SRM and LRE state integrations — the caller supplies the appropriate \(\dot{m}_{in}\) (propellant regression for SRM, summed injector flows for LRE).

machwave.core.mass_balance

compute_chamber_pressure_mass_balance(chamber_pressure, external_pressure, mass_flow_in, free_chamber_volume, throat_area, k, R, flame_temperature, discharge_coefficient=1.0)

Right-hand side of the chamber pressure ODE from a control-volume mass balance.

Handles both choked and sub-critical nozzle flow (Seidel 1965, Eq. 35).

Parameters:

Name Type Description Default
chamber_pressure float

Chamber pressure [Pa].

required
external_pressure float

External pressure [Pa].

required
mass_flow_in float

Mass flow rate into the chamber [kg/s].

required
free_chamber_volume float

Chamber free volume [m^3].

required
throat_area float

Nozzle throat area [m^2].

required
k float

Isentropic exponent of the mix.

required
R float

Gas constant per molecular weight [J/(kg·K)].

required
flame_temperature float

Flame temperature [K].

required
discharge_coefficient float

Discharge coefficient.

1.0

Returns:

Type Description
tuple[float]

Derivative of chamber pressure with respect to time, as a one-tuple.

Source code in machwave/core/mass_balance.py
def compute_chamber_pressure_mass_balance(
    chamber_pressure: float,
    external_pressure: float,
    mass_flow_in: float,
    free_chamber_volume: float,
    throat_area: float,
    k: float,
    R: float,
    flame_temperature: float,
    discharge_coefficient: float = 1.0,
) -> tuple[float]:
    """
    Right-hand side of the chamber pressure ODE from a control-volume mass balance.

    Handles both choked and sub-critical nozzle flow (Seidel 1965, Eq. 35).

    Args:
        chamber_pressure: Chamber pressure [Pa].
        external_pressure: External pressure [Pa].
        mass_flow_in: Mass flow rate into the chamber [kg/s].
        free_chamber_volume: Chamber free volume [m^3].
        throat_area: Nozzle throat area [m^2].
        k: Isentropic exponent of the mix.
        R: Gas constant per molecular weight [J/(kg·K)].
        flame_temperature: Flame temperature [K].
        discharge_coefficient: Discharge coefficient.

    Returns:
        Derivative of chamber pressure with respect to time, as a one-tuple.
    """
    critical_pressure_ratio = get_critical_pressure_ratio(k=k)
    pressure_ratio = external_pressure / chamber_pressure

    if pressure_ratio <= critical_pressure_ratio:  # choked
        H = (k**0.5) * (2 / (k + 1)) ** ((k + 1) / (2 * (k - 1)))
    else:
        H = (
            ((2 * k / (k - 1)) ** 0.5)
            * pressure_ratio ** (1 / k)
            * (1 - pressure_ratio ** ((k - 1) / k)) ** 0.5
        )

    mass_flow_out = (
        discharge_coefficient
        * chamber_pressure
        * throat_area
        * H
        / (R * flame_temperature) ** 0.5
    )

    chamber_pressure_derivative = (R * flame_temperature / free_chamber_volume) * (
        mass_flow_in - mass_flow_out
    )
    return (chamber_pressure_derivative,)