v2.0 · FastAPI + WebGL
Regenerative Cooling Design Tool

A rocket chamber wall is one bad channel geometry away from burning through.

First-principles 1D thermal / hydraulic / structural solver for regeneratively-cooled liquid rocket engine chambers — real-gas combustion, real-fluid coolant properties, coupled Bartz–Gnielinski–Chen heat transfer, and seven independent design margins that either pass or don't. No curve-fit shortcuts, no black box.

7
design margins / run
120
axial stations
54.4
MW/m² peak q, demo case
4113
geometries screened, 1 pass
solver.py — run_analysis_stream()
localhost:8000 — ANALYSIS
Live analysis view: 3D nozzle heat-flux map, axial heat flux curve, and the seven-row margins table
Live run, RP-1/LOX, 45 kN, 90 rectangular channels, CuCrZr liner. Peak heat flux 54.4 MW/m² at the throat; the margins table on the right is computed from the same march, station by station — four of seven checks fail for this particular geometry, which is the point: the tool reports what's actually true, not what looks good in a demo.
Why this exists

Regen cooling is where student rocket designs quietly fail

Most of the tools available to a university propulsion team are one of two things: a CEA-only performance spreadsheet that never touches the wall, or a single-point 1D model that assumes constant coolant properties and reports one number for the whole chamber. Neither catches a throat that dry-boils, a rib that fatigues in 40 cycles, or a channel that seizes 6× over its pressure budget.

What breaks first, in order

Hot-wall temperature and thermal stress peak at the throat because heat flux tracks At/A to the 0.9 power. Coolant-side film coefficient is lowest exactly where the tube is narrowest for structural reasons. Two-phase boiling, if it happens, happens downstream of the throat where bulk temperature has climbed — and dryout there is a wall-burnthrough event, not a performance hit.

What this tool actually checks

Every run produces: max hot-wall temperature vs. melt/service limit, max thermal hoop stress vs. 0.67σy(T), coolant ΔP vs. 25% of inlet pressure, max coolant velocity vs. an erosion limit, coolant-wall and bulk temperature vs. coking limits, critical-heat-flux margin, and estimated low-cycle-fatigue life via Coffin–Manson. All seven, every station, every run.

01 · System

One causal march, three fidelity tiers

Inputs flow forward through combustion → gas dynamics → a per-station coupled solve → structural margins. The same physics core is exposed at three speeds, because a 4000-point optimizer sweep and a single certified design case have very different runtime budgets.

UserInputs Pc · O/F · thrust · geometry · wall alloy · coolant │ ▼ combustion.py Cantera HP-flash equilibrium (real product composition)validated table fallback if Cantera/mechanism unavailable │ T0, γ, cp, Pr, MW, x_CO2, x_H2O ▼ materials.py gasdynamics.py area–Mach march, Brent's method per station k(T), σ_y(T), E(T), │ M(x), T_static(x), T_recovery(x) α(T), T_melt, LCF coeffs ▼ ┌──────────────────────────────────────────────────────────────┐ │ solver.py — per-station fixed-point march │ │ │ │ Bartz h_g(x) ⇄ wall series-resistance + rib fin ⇄ h_c(x) │ │ (heat_transfer.py, pressure.py Δp) │ │ │ │ + implicit tridiagonal axial-conduction pass (one-shot, │ │ scipy.linalg.solve_banded — unconditionally stable) │ └──────────────────────────────────────────────────────────────┘ │ T_wg(x), T_wc(x), q(x), ΔP, σ_thermal(x) ▼ structures.py 7 design margins · Coffin–Manson LCF · Kutateladze CHF │ ├──▶ run_batch_analysis_fast() vectorised NumPy grid screen (optimizer) ├──▶ differential_evolution() continuous refinement past the grid (SciPy) └──▶ run_monte_carlo() Gaussian tolerance-stack → P(all margins pass)

Three solver tiers, same physics

PathUsed byNotes
run_analysis_streamLive UI (WS)Full-fidelity, streams per-station
run_analysis / _march/api/analyze, testsFull-fidelity, single shot
run_batch_analysis_fastOptimizer sweepVectorised, documented simplifications

Coolant flow direction

Counter-current cooling (the default, and the harder case) integrates the coolant energy balance backward from the nozzle exit while the gas-side march runs forward from the chamber — the two loops meet through the fixed-point iteration at each station, not through a global relaxation sweep.

localhost:8000 — input panel
Full parameter panel: engine, geometry, channel, and coolant inputs
Every physical input is exposed — propellant pair, chamber/throat/exit geometry, channel cross-section and packing, 20 wall alloys, 10+ coolants with real CoolProp backends. Nothing is hardcoded to one engine class.
02 · Correlations

The correlation stack, in march order

Numbered the way the solver actually walks through them — gas side, through the wall, into the coolant, and out as a pressure drop and a stress. Every correlation is cited; none of the constants are guessed.

01
Combustion equilibrium Cantera HP-flash
Real chemical-equilibrium solve at the actual injector mixture ratio and pressure — not a fixed T₀/γ lookup. Falls back to a literature-validated static table (18 propellant pairs) with a logged warning if Cantera or the mechanism is unavailable.
c* = √(R·T₀) / [√γ · (2/(γ+1))(γ+1)/(2(γ-1))]
gri30 / h2o2 / nDodecane mechanismsmixture-averaged transport18-pair table fallback
02
Nozzle gas dynamics area–Mach, isentropic
Mach number solved per station from the local area ratio via Brent's method (bounded root-find, not a lookup table), frozen composition downstream of the throat.
1/M · [(2/(γ+1))(1 + ((γ-1)/2)M²)](γ+1)/(2(γ-1)) = A/A*
03
Hot-gas film coefficient Bartz, Sutton & Biblarz 9th ed. Eq. 8-24
Classical Bartz correlation with the property-variation correction σ and a throat radius-of-curvature term — the curvature exponent matters: a sharper throat fillet measurably raises local hg.
hg = (0.026/Dt0.2) · (μ0.2cp/Pr0.6) · (Pc/c*)0.8 · (Dt/Rcurv)0.1 · (At/A)0.9 · σ
04
Gas-side radiation gray-gas exponential absorption
Radiative flux driven by the actual CO₂/H₂O mole fractions read out of the Cantera equilibrium state, not an assumed composition, with the beam length scaled to the local chamber/nozzle diameter so radiation correctly tapers through the expansion.
εg = 1 − exp[−kabs·(pCO2+pH2O)·Lbeam]
05
Wall conduction series resistance + 2D rib-fin + axial coupling
Radial path is a series resistance through a temperature-dependent-k liner with a 2D shape-factor correction for rib fin efficiency (NASA SP-8089 style). A second pass then couples axial (streamwise) conduction between stations — the throat's sharp temperature peak physically bleeds heat into its cooler neighbours. Solved as one implicit tridiagonal system, not an under-relaxed outer loop that can diverge.
qgen,i + KL(Ti-1−Ti) + KR(Ti+1−Ti) = heff,iAi(Ti−Tcool,i)
scipy.linalg.solve_bandedunconditionally stableO(n) direct solve
06
Coolant-side heat transfer Gnielinski + Chen + Krasnoshchekov–Protopopov
Single-phase turbulent flow uses Gnielinski (coupled to the same Colebrook friction factor used for ΔP, so wall roughness affects both consistently), corrected for entrance length (Boelter, turbulent tube entry), channel curvature (Dean number, Bergles–Kraus), and supercritical property variation (Krasnoshchekov–Protopopov, both the ρ and cp ratio terms). Subcritical boiling uses the full Chen (1963/66) superposition model — forced-convective and nucleate-boiling contributions summed, with the Forster–Zuber nucleate term and a Reynolds-based suppression factor.
htp = F·hl + S·hnb  Chen (1963)
Nusc = NuDB · (ρwb)0.3 · (cp,film/cp,b)0.4  Krasnoshchekov–Protopopov
07
Pressure drop Colebrook + Lockhart–Martinelli + Idelchik
Darcy friction factor from a Colebrook solve (Brent's method against absolute channel roughness, not a Moody chart read), with a Chisholm (1973, C=20) two-phase multiplier applied automatically once local quality crosses the boiling threshold, plus Idelchik entry/exit minor losses keyed to actual inlet geometry (sharp vs. rounded).
φlo² = 1 + C/Xtt + 1/Xtt²  Chisholm (1973)
08
Structural margins hoop stress + Coffin–Manson + Kutateladze
Thermal hoop stress from mid-wall temperature-dependent E and α, checked against 0.67σy(T). Low-cycle-fatigue life estimated via Coffin–Manson from the same thermal strain range. Critical heat flux checked against the Kutateladze/Zuber pool-boiling bound — called out explicitly as a conservative screening check, not a forced-convection-validated CHF correlation (see Roadmap).
σ = E(Tm)·α(Tm)·ΔT / [2(1−ν)]
09
Design-space search vectorised sweep + differential evolution + Monte Carlo
A vectorised NumPy pass screens a discrete grid in seconds, SciPy's differential evolution then refines continuously past the grid, and a separate Monte Carlo pass perturbs manufacturing/process tolerances (Gaussian, 1σ) through the full-fidelity solver to report reliability as P(all 7 margins pass), not a single point safety factor.
03 · Live output

Every station, every quantity, live

The WebGL nozzle view is colour-mapped directly from solver arrays (heat flux, Twg, Twc, thermal stress) — the same numbers driving the Plotly traces beside it, streamed station-by-station over a WebSocket as the march runs.

04 · Honesty

The demo case fails four of seven margins — on purpose

The shipped default (45 kN, RP-1/LOX, 90 channels, CuCrZr) is under-cooled. That's deliberate: a tool that always shows green checkmarks in its own screenshots is a tool you can't trust. Below is the actual margins table from the run captured above.

CheckValueLimitUnitsRatioStatus
Max hot-wall temp Twg1.18k1.35kK1.14PASS
Max thermal stress575.10M93.26MPa0.16FAIL
Coolant ΔP (w/ SF)7.69M2.00MPa0.26FAIL
Max coolant velocity67.1780.00m/s1.19PASS
Coolant-wall temp Twc820.86700.00K0.85FAIL
Coolant bulk temp645.26560.00K0.87FAIL
Estimated LCF life387.56100.00cycles3.88PASS

Regression coverage (test_margin_agreement.py) sweeps a grid of channel geometries through both the vectorised optimizer surrogate and the full march, and asserts the two never flip a pass/fail decision and agree on total ΔP within 15% — the fast path used for a 4000-point sweep is checked against the slow path used for the final number.

05 · Search

4113 geometries in, one passing configuration out

Sweep channel count, height, width, rib width, wall thickness and coolant flow rate; pre-filter anything that fails packing (channels + ribs exceeding the throat circumference); vectorised-batch solve every survivor; then let differential evolution refine continuously past the grid.

06 · Implementation

Stack

Python 3.11
FastAPI
WebSocket streaming
NumPy / SciPy
Cantera (real-gas equilibrium)
CoolProp (real-fluid properties)
Three.js WebGL
Plotly.js
Vanilla JS — no framework
Single-file frontend
07 · What's not in here

Parked, not hidden

Two known gaps, left out on purpose rather than papered over with a guessed constant.

PARKED

Katto–Ohno critical heat flux

The Kutateladze/Zuber pool-boiling bound used today is a conservative stand-in — it has no mass-flux (G) term, so it doesn't capture forced-convection CHF enhancement. The correct upgrade is Katto–Ohno (1984), a five-regime correlation whose constants sit behind paywalled primary literature. Rather than guess coefficients for a safety-critical wall-burnthrough margin, it's left out until a verifiable source (textbook table or the original paper) is in hand.

PARKED

Shifting-equilibrium nozzle expansion

Combustion products are equilibrated once, at chamber conditions, and expanded with frozen composition (constant γ, T₀) — standard practice, and conservative for Isp. True shifting equilibrium would re-run a Cantera SP-flash at every station, which changes Isp, which changes required ṁ, which cascades through the entire coupled solve. Real, but high blast-radius; not added silently.