turbo_moc.core.classes module#
- class turbo_moc.core.classes.FluidManager(fluid_name: str, backend: str, P0, T0, Q0)[source]#
Bases:
objectManages the single instance of the jaxprop fluid handler and stores the persistent stagnation properties (s0, h0) and critical state.
Methods
V_of_nu(nu_val)Inverse of nu_of_V: velocity at a given real-gas Prandtl-Meyer angle (radians).
a_of_V(V)Speed of sound at velocity V along the same s=s0 isentrope used for nu_of_V (cheap table lookup, no fluid call).
calculate_static_properties(V, s0, h0)Calculates static properties for a single scalar velocity V.
nu_of_V(V)Real-gas Prandtl-Meyer angle (radians) at velocity V, referenced to the critical/jump state.
rho_of_V(V)Density at velocity V along the same s=s0 isentrope used for nu_of_V (cheap table lookup, no fluid call).
calculate_stagnation_state
- V_of_nu(nu_val)[source]#
Inverse of nu_of_V: velocity at a given real-gas Prandtl-Meyer angle (radians).
- a_of_V(V)[source]#
Speed of sound at velocity V along the same s=s0 isentrope used for nu_of_V (cheap table lookup, no fluid call).
- calculate_static_properties(V, s0, h0)[source]#
Calculates static properties for a single scalar velocity V.
Calls the raw (non-JAX) Fluid object directly instead of routing through jitted_calculate_static_properties + FluidJAX.get_state’s nested JIT/pure_callback layers. This is the sole caller of that path and V is always a single scalar here (ThermoFlowProp, called once per point throughout Stage 2/3 marching – never batched), so the JIT layer was paying a full XLA trace+compile tax on nearly every call for zero reuse benefit (same diagnosis and fix as FluidManager._build_nu_table’s per-point loop).
- class turbo_moc.core.classes.InternalPoint(fluid_manager: FluidManager, x=None, y=None)[source]#
Bases:
objectClass representing a single point in a Method of Characteristics (MoC) grid. Stores geometry and flow properties. Includes a method to compute the sonic condition (Sauer real-gas approach) using a 1D root finder and the main Sauer function for the sonic line.
Methods
Inv_Wall_pnt_vertical(pwall, p33, p11, wall_ang)Wall-point unit process for an initial data line that is VERTICAL (constant x -- e.g. the flat flashing throat front, x=0 for every y).
ReflexLine(n, l)Reflex Line position calculator (ReflexZone).
compute_from_characteristics(paramP, paramM)Computes a new MOC point.
Calculates the critical (sonic, M=1) state.
corner_point_prop(x_corner, y_corner, theta_rad)Direct (search-free) sharp-corner unit process for a centred Prandtl-Meyer expansion fan.
Recompute thermo/flow properties only if u or v changed since the last ThermoFlowProp call -- see _thermo_fresh above.
euler_predictor_corrector(paramP, paramM[, ...])Euler Predictor-Corrector optimized with Object Reuse and Thermophysical Properties (a, theta) calculated ONLY ONCE before the iteration loop.
interior_point_prop(paramP, paramM)Exact planar (delta=0) interior-point unit process.
solve_centerline_point(paramP, paramM)Calculates the properties of a new point on the Axis of Symmetry (y=0).
to_dict()Inv_Wall_pnt
copy_from
- Inv_Wall_pnt_vertical(pwall, p33, p11, wall_ang)[source]#
Wall-point unit process for an initial data line that is VERTICAL (constant x – e.g. the flat flashing throat front, x=0 for every y). Inv_Wall_pnt interpolates position along the “line 3-1” by x, which is degenerate when x33 == x11; here the same shooting method is redone parametrizing by y instead (x is fixed along the line, so any y between p33.y and p11.y is a valid position). Otherwise identical to Inv_Wall_pnt.
- ReflexLine(n, l)[source]#
Reflex Line position calculator (ReflexZone). Calculates the position of the new wall point (or reflection point) a distance ‘l’ away, projected along the Mach angle plus flow angle.
- Parameters:
n (float): Step multiplier (1 in the case of wall projection, n_step in others). l (float): Calculated distance based on mass flow ratio.
- compute_from_characteristics(paramP, paramM, EPC=0)[source]#
Computes a new MOC point.
EPC defaults to 0 (no corrector): interior_point_prop / solve_centerline_point now solve the EXACT planar (theta, nu) compatibility relations directly (see their docstrings), so there is no thermodynamic linearisation error left for a corrector step to reduce. euler_predictor_corrector still uses the old linearised Q/R1/S algebra and would only reintroduce that error if applied here, so it is skipped unless explicitly requested.
- compute_sonic_state()[source]#
Calculates the critical (sonic, M=1) state. FIX: Reads the pre-calculated critical state properties from the Manager.
- corner_point_prop(x_corner, y_corner, theta_rad)[source]#
Direct (search-free) sharp-corner unit process for a centred Prandtl-Meyer expansion fan.
Within a centred simple wave emanating from a sharp corner, theta and nu are related EXACTLY by theta - theta_ref = nu(V) - nu_ref (Zebbiche & Youbi 2007, eqs. 1-2 restricted to a simple wave). Since the corner’s reference state (the throat/post-jump flow, before any turning) has theta_ref=0 and is exactly where this class’s nu(V) table is zeroed (see FluidManager._build_nu_table), this collapses to theta = nu(V) at the corner itself for ANY prescribed turning angle theta – no dependence on neighbouring front data at all.
This replaces Inv_Wall_pnt / Inv_Wall_pnt_vertical’s iterative shooting search (which interpolates against bracketing front points to locate a wall intersection) for the MLN sharp-corner case, where the wall position is fixed and the flow angle at each fan angle is prescribed by construction – there is nothing to search for. Used by moc_solver_mln.py’s run_stage_3_kernel_mln.
- ensure_fresh()[source]#
Recompute thermo/flow properties only if u or v changed since the last ThermoFlowProp call – see _thermo_fresh above.
- euler_predictor_corrector(paramP, paramM, max_iter=1000, tol_x=1e-06, tol_v=1e-06)[source]#
Euler Predictor-Corrector optimized with Object Reuse and Thermophysical Properties (a, theta) calculated ONLY ONCE before the iteration loop.
- interior_point_prop(paramP, paramM)[source]#
Exact planar (delta=0) interior-point unit process.
Zebbiche & Youbi (2007), eqs. 1-2, 9-23: for planar flow the characteristic compatibility relations are the exact invariants:
theta + nu(V) = constalong C- (paramM’s characteristic)theta - nu(V) = constalong C+ (paramP’s characteristic)
so theta3 and nu3 (hence V3) at the new point follow directly from the two foot points’ theta and nu, with no linearisation at all in the thermodynamic part. This replaces the general (axisymmetric- capable) linearised Q/R1/S solve, which is only approximate for this delta=0 problem and was the source of accumulating error in earlier attempts. Geometry (the (x, y) intersection) still needs the standard linearised treatment, unchanged from before.
- solve_centerline_point(paramP, paramM)[source]#
Calculates the properties of a new point on the Axis of Symmetry (y=0).
Exact planar (delta=0) compatibility: theta + nu(V) = const along the C- characteristic through paramM (Zebbiche & Youbi 2007, eqs. 1, 9, 27-28). With the axis reflection boundary condition theta=0, this gives nu(V_axis) = paramM.theta + nu(paramM.V) directly – no linearised Q/R1/S algebra needed. Geometry (intersection) is unchanged: paramP supplies the mirror-image C+ line used to locate the axis point.
- turbo_moc.core.classes.MakeVec(x1, y1, x2, y2)[source]#
Converts coordinates into a normalized unit vector along the segment P2 -> P1.
- turbo_moc.core.classes.MassFlowCalc(IO)[source]#
Calculates the mass flow rate based on the original formula, FIXED to use ThermoFlowProp() and JAX math functions (jnp).
- turbo_moc.core.classes.RotVec(vec, angle_deg)[source]#
Rotates a vector (u, v) by an angle (degrees). Returns [u_rot, v_rot].
- turbo_moc.core.classes.build_convergent_inlet(y_t, rho_d, L_conv, n_points=30)[source]#
Purely geometric (no MOC/flow physics) upstream extension of the wall, from the inlet plane (x=-L_conv) back to the throat (x=0, y=y_t) – an “S-curve” convergent contour made of two tangent circular arcs of EQUAL radius rho_d, meeting at a single inflection point where concavity flips (a reversed-curvature/ogive convergent inlet, the same construction long-radius flow nozzles use).
Arc 1 (throat -> inflection) is literally the SAME circle as the kernel region’s own throat rounding (
turbo_moc.core.classes.rotation_around_center, center(0, rho_d+y_t), same radiusrho_d) – continued backward (negative angle) instead of forward, so the wall is curvature-continuous across the throat, not just position/slope- continuous.Arc 2 (inflection -> inlet) is arc 1’s own point reflection through the inflection point P1 (Q -> 2*P1 - Q for every point Q on arc 1). A point reflection through a point ON a curve preserves the tangent LINE there while flipping curvature sign – exactly “changes concavity” – and, since the throat point itself (ang=0) has zero slope by the throat’s own left-right symmetry, its reflection (the inlet-plane end point) also has EXACTLY zero slope: the wall meets the inlet duct already parallel to the axis, not just “almost” parallel, with no separate solve needed for where that happens.
Splitting L_conv exactly in half between the two arcs is what makes arc 2’s radius come out equal to arc1’s (rho_d) automatically – algebraically, R2 = L2/sin(phi1) = L2/(L1/rho_d) = rho_d when L1 == L2 – rather than requiring a free second radius.
- Parameters:
- y_tfloat
Throat half-height (m).
- rho_dfloat
Throat/kernel radius of curvature (m) – same value already used for the divergent side’s own throat rounding (
rotation_around_center’srho_d).- L_convfloat
Total axial length of the convergent inlet section (m), split evenly between the two arcs (L1 = L2 = L_conv/2). Must satisfy L_conv/2 < rho_d (arc 1 cannot reach further upstream than its own quarter-circle) – raises ValueError otherwise.
- n_pointsint
Number of points per arc (2*n_points total, excluding the duplicated throat point).
- Returns:
- dict
{“x”, “y”}, ordered from the inlet plane (x=-L_conv) to just before the throat (x=0 itself is NOT included, so this can be directly prepended to a wall array that already starts there).
- turbo_moc.core.classes.get_fluid_manager(fluid_name: str, backend: str, P0, T0, Q0) FluidManager[source]#
Returns a cached FluidManager for (fluid_name, backend, P0, T0, Q0), building it only on first use for that combination.
FluidManager depends only on the fluid and inlet condition, not on nozzle geometry – building it (notably its Prandtl-Meyer table, see _build_nu_table) is the dominant one-time cost per solver instantiation. A design/optimization loop that sweeps geometry at a FIXED operating point should call this once and pass the result to MOCSolver/ MOCSolverMLN’s fluid_manager= argument, instead of letting each one rebuild it from scratch. Changing the fluid or the operating point is a different cache key – it does NOT reuse a stale manager.
Q0 = float(‘nan’) marks a subcooled-liquid/single-phase inlet (use P0,T0 instead) – this is the common case here, so NaN needs explicit handling: two separately-constructed float(‘nan’) objects are neither is nor == to each other (NaN != NaN by IEEE 754), so using Q0 directly as a dict-key component would silently never hit the cache for this case (every call would look like a distinct, never-before-seen key). Collapsed to a fixed sentinel string instead so all NaN-Q0 calls at the same (fluid, P0, T0) correctly share one cached manager.
- turbo_moc.core.classes.jax_solve_critical_state(fluid, h0, s0, mach_tol=0.05)[source]#
Critical-state solver with a flashing fallback.
Tries to find the smooth sonic point (Ma=1) via bracketed root-finding.
If no smooth root exists — or the “root” brentq returns actually sits on a discontinuous jump in sound speed rather than a true zero-crossing, which happens when a subcooled-liquid isentrope crosses the saturation line: a(P) jumps from the high single-phase liquid branch to the low collapsed two-phase branch over an arbitrarily small pressure step, so bisection can satisfy its P-tolerance while landing on either side of the jump with a residual nowhere near zero — falls back to locating the saturation pressure and nudging slightly into the two-phase dome to pick up the collapsed (HEM) speed of sound. Mirrors the approach in moc_tool_flashing/srcMOC/moc_classes.py, with an added post-hoc validation of the smooth root (that file trusts brentq unconditionally, which silently returns a bogus state when the isentrope crosses a jump like this one).
- turbo_moc.core.classes.jitted_calculate_static_properties(fluid, V: Array, s0: float, h0: float)[source]#
- turbo_moc.core.classes.jitted_corrector_algebra(xM, yM, uM, vM, aM, thetaM, xP, yP, uP, vP, aP, thetaP, u_P_ref, v_P_ref, u_M_ref, v_M_ref, delta)[source]#
- turbo_moc.core.classes.jitted_moc_predictor_algebra(xM, yM, lamM, SM, QM, R1M, uM, vM, xP, yP, lamP, SP, QP, R1P, uP, vP)[source]#
JIT-compiled function to solve MOC geometry and linear compatibility system (Predictor Step).
Parameters are properties of the C- foot point (M) and the C+ foot point (P).
- turbo_moc.core.classes.rotation_around_center(y_t, ang, rho_d)[source]#
Rotate the point (0, y_t) about the center located at (0, rho_d + y_t) by an angle ang (in degrees), counter-clockwise.
This is used in the kernel-region construction of the Method of Characteristics to generate rotated characteristic boundaries starting from the throat radius.
- Parameters:
- y_tfloat
The initial y-coordinate of the point to be rotated (the x-coordinate is fixed to 0).
- angfloat
Rotation angle in degrees. Positive = counter-clockwise rotation.
- rho_dfloat
Kernels offset distance; the rotation center is at (0, rho_d + y_t).
- Returns:
- x_rotfloat
Rotated x-coordinate of the point.
- y_rotfloat
Rotated y-coordinate of the point.
Notes
- The transformation performed is:
Translate the point so the center becomes the origin
Apply rotation matrix R(ang)
Translate back to the original center