Simulating the circulation, then closing the loop on it

The base model divides the body into thirteen compartments — lungs, brain, heart, liver and GI, kidneys, muscle, skin, bone, fat, thyroid, sex glands, adrenals, and a catch-all — each with its own volume, resistance, compliance, and inertance, plus arteries and veins. A four-chamber heart with time-varying elastance and pressure-driven valves supplies the flow.

The first question the model answers is whether pulsatility matters. One unit of solute is injected and tracked as it circulates, run twice: once with the heart actually beating, and once replaced by steady flow.

The useful result is the trade. Pulsatility costs roughly nineteen times the compute and buys almost nothing on the slow measures: mean transit time through the lungs differs by two hundredths of a second, and both conditions converge on the same well-mixed final concentration, which checks out against one unit distributed across total blood volume. Where it does matter is recirculation timing, which arrives about five seconds later with a beating heart, and peak amplitude, which is far higher. So if the question is steady-state distribution, constant flow is the right model; if it is peak exposure or arrival timing, it is not.

Watching the animation, the fast compartments are the ones you would expect from perfusion per unit mass — heart, brain, kidneys, thyroid, adrenals. Bone, fat, and the catch-all are slowest in both conditions.

Part two turns the open transport model into a closed-loop physiological controller. Carbon dioxide is generated by every tissue in proportion to its blood supply, carried by the same transport machinery, and cleared only at the lungs. The clearance rate is not fixed — it is set by the controller, which makes the whole thing a negative feedback loop.

Mapped onto control-system blocks: the chemoreceptors are the sensor and comparator, the CNS is the controller, and the diaphragm and lungs are the plant. Partial pressure of CO₂ goes in, pH comes out of the Henderson-Hasselbalch equation, and minute ventilation comes out of a parametric function fitted to the physiological pH–PPCO₂–MV relationship. Ventilation then drives depletion at the lungs, D = PPCO₂ × 0.005 MV, which feeds back into the CO₂ the tissues are still producing.

The ventilation function itself is a three-point parametric fit through the clinical relationship, anchored at 30 mmHg / pH 7.52 / 11 L·min⁻¹ and 50 mmHg / pH 7.15 / 43 L·min⁻¹. Homeostatic targets are a PPCO₂ of 25 to 45 mmHg, a pH of 7.4, and minute ventilation between 3 and 50 L·min⁻¹.

A model of a controller is only worth as much as its response to being knocked off target, so the loop is disturbed at t = 300 s with a simulated vomiting episode: acid is lost from the liver and GI compartment, producing metabolic alkalosis and a step rise in CO₂ in that tissue.

Transport delay from the disturbance to the measured signal is 10 seconds, and to the command signal another 10 — the time it physically takes disturbed blood to reach the chemoreceptors and for the response to reach the lungs. The loop then settles back to its pre-disturbance values over a homeostatic time constant of roughly 100 seconds, with no sustained offset and no oscillation. The controller is stable and the delays are the ones circulation time would predict.

Around 1,900 lines of MATLAB across six files. Compartment properties load from a CSV rather than being hard-coded, so the body can be re-partitioned without touching the solver.

The transport routine is the core of it. For each connection it reads the flow, checks the sign to decide which compartment is upstream, computes the flux as concentration times flow over volume, and applies it as a loss to one compartment and a gain to the other. Handling the sign explicitly is what lets the same routine serve both the systemic tissues and the reversed pulmonary path.