Nodal Pressure Network Analysis

The Mathematical Engine Behind SootSafe

Abstract

Fire-safety ventilation systems must maintain precise pressure differentials across compartment boundaries during a fire event. Verifying that a designed system meets those requirements requires solving a system of non-linear equations that couples fluid resistance, air density, and fan performance across the full building network. This paper describes the mathematical model and numerical method used by SootSafe to perform that calculation.

1. The Physical System

A ventilated building can be modelled as a directed graph. Rooms, shafts, and outdoor reference points become nodes (junctions). Ducts, leakage paths, fans, and dampers become edges (branches). At steady state, the pressure at every junction is uniquely determined by two physical laws operating simultaneously:

Mass conservation at each junction — mass flow in equals mass flow out.
Resistance law on each branch — the pressure difference across a component drives a flow through it.

Because air density varies with temperature (a 600 °C fire room and a 20 °C corridor have air that is roughly 3× different in density), and because most resistance laws are quadratic in flow, the resulting system is non-linear. Simple matrix inversion is insufficient; an iterative method is required.

2. Mathematical Formulation

2.1 Nodes as Conservation Laws

Let the network contain $N$ junctions with unknown gauge pressures $\mathbf{p} = (p_1, \ldots, p_N)^\top$ . A subset of nodes may have fixed (boundary) pressures — typically the outdoor reference at $0\,\text{Pa (gauge)}$ .

At every variable junction $i$ , mass must balance:

$F_i(\mathbf{p}) = \sum_{j \in \mathcal{N}(i)} \rho_{ij}^{\uparrow} \, Q_{ij}(\mathbf{p}) \;+\; \dot{m}_i^{\text{target}} = 0$

where

$\mathcal{N}(i)$ — neighbours connected to junction $i$ ,
$\rho_{ij}^{\uparrow}$ — air density at the upstream end of branch $ij$ (determined by flow direction),
$Q_{ij}$ — volumetric flow on branch $ij$ (m³/s, positive toward $i$ ),
$\dot{m}_i^{\text{target}}$ — externally imposed mass flow, e.g. a mechanical supply or exhaust (kg/s).

The full system is $\mathbf{F}(\mathbf{p}) = \mathbf{0}$ , a vector of $N$ non-linear equations.

2.2 Air Density and the Ideal Gas Law

Air density at each junction is pre-computed once from the ideal gas law:

$\rho_i = \frac{P_{\text{atm}}}{R_{\text{air}} \, T_i}$

Symbol	Value	Unit
$P_{\text{atm}}$	101 325	Pa
$R_{\text{air}}$	287.05	J / (kg·K)
$T_i$	node temperature	K

Reference: ISO 10294-3.

Using upstream density on each branch (rather than an average or fixed value) ensures that mass is conserved exactly even when adjacent compartments are at very different temperatures. Without this correction, a single branch connecting a 20 °C corridor to a 600 °C fire room would violate $\sum \dot{m} = 0$ by a factor proportional to the density ratio.

3. Branch Resistance Laws

3.1 The Orifice Flow Equation

For all passive components — leakage gaps, passive vents, orifice plates — the governing relationship is:

$Q = \operatorname{sign}(\Delta P) \cdot \sqrt{\frac{|\Delta P|}{k}}$

where $k\,[\text{Pa·s}^2/\text{m}^6]$ is the component's hydraulic resistance and $\Delta P = p_i - p_j$ is the pressure difference driving flow from junction $i$ to junction $j$ . Equivalently:

$|\Delta P| = k \, Q^2$

For a passive vent with flow coefficient $C$ (m³/s per Pa $^{1/2}$ , per EN 13141 / ISO 10294-3):

$k = \frac{1}{C^2}$

3.2 Circular Ducts — Darcy–Weisbach

For a circular duct of length $L$ , diameter $D$ , and roughness $\varepsilon$ :

$k = \frac{8 \, f \, \rho \, L}{\pi^2 \, D^5}$

This follows directly from the Darcy–Weisbach pressure-loss formula $\Delta P = f \, \frac{L}{D} \, \frac{\rho v^2}{2}$ after substituting $v = Q/A$ with $A = \pi D^2/4$ .

The Darcy friction factor $f$ is evaluated at the current Reynolds number $Re = \rho Q D / (A \mu)$ , using $\mu = 1.81 \times 10^{-5}\,\text{Pa·s}$ for air at 20 °C.

3.3 The Haaland Friction Factor

Rather than iterating the implicit Colebrook–White equation, SootSafe uses the Haaland (1983) explicit formula for turbulent flow ( $Re \geq 3000$ ):

$\frac{1}{\sqrt{f}} = -1.8 \log_{10}\!\left[\left(\frac{\varepsilon / D}{3.7}\right)^{1.11} + \frac{6.9}{Re}\right]$

Reference: Haaland, S. E. (1983). Simple and Explicit Formulas for the Friction Factor in Turbulent Pipe Flow. Journal of Fluids Engineering, 105(1), 89–90.

For laminar flow ( $Re < 2300$ ), the Hagen–Poiseuille result applies:

$f = \frac{64}{Re}$

In the transition zone ( $2300 \leq Re < 3000$ ), $f$ is linearly interpolated between the laminar and turbulent values. This smooth blending is essential: a discontinuous friction factor produces a discontinuous Jacobian, which causes the Newton–Raphson solver to oscillate rather than converge.

3.4 Other Passive Components

Component	Resistance formula
Passive vent (flow coeff. $C$ )	$k = 1/C^2$
Local loss (loss coeff. $\zeta$ , diameter $D$ )	$k = 8\zeta\rho / (\pi^2 D^4)$
Damper (open fraction $f \in (0,1]$ )	$k = 1 / f^2$
Check valve	$k$ from underlying component; $Q = 0$ clamped if flow is reversed

3.5 Fans — Active Pressure Sources

A fan adds pressure rather than consuming it. Its characteristic curve is a quadratic polynomial (IEC 60050-415):

$\Delta P_{\text{fan}}(Q) = a Q^2 + b Q + c \qquad (Q \geq 0)$

where $a < 0$ (the curve droops at high flow), $b$ and $c$ are coefficients fitted to manufacturer data. For reverse flow ( $Q < 0$ ) the fan contribution is set to zero (stall model). The net branch resistance for a branch containing both a duct and a fan is:

$\Delta P_{\text{branch}} = k_{\text{duct}} \cdot Q^2 \; - \; \Delta P_{\text{fan}}(Q)$

Series components within a branch are summed: total resistance $k = \sum_m k_m$ and total fan boost $\Delta P_{\text{fan}} = \sum_m \Delta P_{\text{fan},m}$ .

4. Numerical Solution — Newton–Raphson

4.1 Why Newton–Raphson?

The quadratic resistance law $|\Delta P| = k Q^2$ makes $\mathbf{F}(\mathbf{p})$ non-linear. Newton–Raphson (NR) is the standard choice: it achieves quadratic convergence near the solution, meaning the number of correct digits roughly doubles each iteration, typically reaching machine precision in 3–8 steps.

The iteration is:

$\mathbf{p}^{(k+1)} = \mathbf{p}^{(k)} + \alpha \, \boldsymbol{\delta}$

where $\boldsymbol{\delta}$ solves the linear system

$\mathbf{J}(\mathbf{p}^{(k)}) \, \boldsymbol{\delta} = -\mathbf{F}(\mathbf{p}^{(k)})$

and $\mathbf{J}_{ij} = \partial F_i / \partial p_j$ is the Jacobian matrix, and $\alpha \in (0, 1]$ is a step-size chosen by a line search.

4.2 Jacobian by Finite Differences

The Jacobian is computed column-by-column via forward differences:

$J_{ij} = \frac{F_i(\mathbf{p} + h \, \mathbf{e}_j) - F_i(\mathbf{p})}{h}$

with an adaptive step size:

$h = \max\!\left(10^{-3},\; 10^{-6} \cdot |p_j|\right) \quad \text{[Pa]}$

The absolute floor of $10^{-3}$ Pa prevents underflow at low-pressure nodes; the relative component keeps accuracy at high pressures (e.g. 1 500 Pa in a pressurised stairwell).

4.3 The "Linear Bridge" — Regularisation at Zero Pressure Difference

The orifice flow equation $Q = \sqrt{|\Delta P|/k}$ has an infinite slope at $\Delta P = 0$ . If the solver ever proposes equal pressures on both sides of a branch, the Jacobian entry diverges and the iteration collapses.

The fix is a linear approximation for small pressure differences:

$Q(\Delta P) = \begin{cases} \dfrac{\Delta P \sqrt{1/k}}{\Delta P_{\text{lam}}} & |\Delta P| < \Delta P_{\text{lam}} \\[8pt] \operatorname{sign}(\Delta P) \sqrt{|\Delta P|/k} & |\Delta P| \geq \Delta P_{\text{lam}} \end{cases}$

with $\Delta P_{\text{lam}} = 1\,\text{Pa}$ . In the linear regime the slope is $\sqrt{1/k}$ , which is twice the slope of the orifice tangent at $\Delta P_{\text{lam}}$ — a modest over-estimate, but one that keeps the Jacobian bounded and well-conditioned everywhere.

At the operating pressures relevant to fire safety ( $50\text{–}1\,500\,\text{Pa}$ ), the $1\,\text{Pa}$ regularisation zone is negligible.

4.4 Convergence Criterion

At each iteration the infinity-norm of the residual is tested:

$\|\mathbf{F}(\mathbf{p}^{(k)})\|_\infty < \tau, \qquad \tau = \max\!\left(10^{-6},\; 10^{-8} \sum_i |\dot{m}_i^{\text{target}}|\right) \quad \text{[kg/s]}$

The absolute floor $10^{-6}\,\text{kg/s}$ prevents declaring convergence on a trivially empty network; the relative term scales naturally with large forced flows.

4.5 Backtracking Line Search

A full Newton step ( $\alpha = 1$ ) can overshoot the solution, especially at the first few iterations when the linearisation is coarse. A simple backtracking scheme halves $\alpha$ until the residual norm decreases:

$\|\mathbf{F}(\mathbf{p}^{(k)} + \alpha\boldsymbol{\delta})\|_2 < \|\mathbf{F}(\mathbf{p}^{(k)})\|_2$

This guarantees descent at every step and provides robustness for networks with strongly non-linear components (check valves, near-stall fans).

4.6 Jacobian Recycling

The dominant cost in each iteration is building and factorising the $N \times N$ Jacobian. When the residual norm is contracting rapidly (ratio $< 0.5$ ), the current Jacobian is still a good local approximation and can be reused for the next solve without rebuilding. The factorisation is discarded only when:

the norm ratio exceeds 1.0 (the solution moved to a different regime), or
the line search accepted a reduced step $\alpha < 0.9$ (the operating point shifted substantially), or
two consecutive iterations showed slow contraction.

In practice this reduces Jacobian builds from $O(N_{\text{iter}})$ to 2–3 per solve, cutting run time by 50–70 % on large networks.

5. Assumptions and Limitations

Assumption	Consequence	Direction
Steady state	Transient effects (door opening, fan spin-up) are not modelled	Neutral — the pseudo-transient simulation (see Pseudo-Transient Fire Simulation) calls this solver once per second to reconstruct time dependence
Incompressible flow	Valid when duct velocities are well below the speed of sound, as in all building HVAC applications	Neutral — accurate to within < 0.1 % at typical HVAC velocities
Single-phase, dry air	Moisture and smoke particle transport are not included	Mildly non-conservative — smoke aerosols slightly increase effective viscosity and duct resistance; ignoring them marginally over-predicts airflow
Adiabatic compartments	Heat transfer between rooms is neglected; temperatures are inputs, not outputs	Conservative if input temperatures come from an adiabatic fire model — adiabatic models over-predict room temperature, which lowers air density and amplifies expansion-driven pressure rise
Scalar density per node	Vertical temperature stratification within a compartment is not resolved	Mildly non-conservative — in a real fire room the hot buoyant layer near the ceiling is substantially hotter than the volume-average temperature; using the average under-predicts the density difference that drives pressure build-up
No flow momentum	Inertial effects (e.g. bends with dynamic pressure recovery) are captured only through loss coefficients, not the full momentum equation	Mildly non-conservative — loss coefficients approximate local pressure losses but do not recover dynamic pressure at expansions; actual resistance may be slightly higher than modelled

The conservatism of the overall calculation depends primarily on the temperature boundary conditions supplied to the solver. The solver's own numerical assumptions are accurate or, at worst, mildly non-conservative; the principal conservative margin in fire-safety applications is introduced by the thermal model upstream (see Pseudo-Transient Fire Simulation, §4.1).

6. Standards Referenced

Standard	Scope
ISO 10294-3	Air leakage measurement; ideal gas air density reference conditions
EN 13141	Orifice flow coefficient convention for ventilation components
IEC 60050-415	Fan characteristic curve conventions
Haaland (1983)	Explicit friction factor formula for turbulent pipe flow

7. Summary

SootSafe models a ventilated building as a pressure-driven flow network and solves for the equilibrium pressure field using Newton–Raphson iteration. The key technical contributions are:

Upstream density propagation — correctly carries temperature-dependent air density from source to sink, ensuring mass conservation in multi-temperature networks.
Haaland friction — an explicit friction model that eliminates the iterative inner loop of Colebrook–White, covering laminar, transition, and fully turbulent regimes in a single continuous formula.
Linear bridge regularisation — prevents solver collapse when branch pressure differences pass through zero, without introducing artificial resistance.
Jacobian recycling with backtracking — achieves robust convergence in 3–8 iterations for typical building networks.

These properties together allow the engine to handle fire scenarios — where a single building simultaneously contains ambient, hot-smoke, and pressurised-stairwell zones — with the same reliability as a routine HVAC calculation.