Blood Flow & Poiseuille's Law — The Fluid Mechanics of Circulation

1. Poiseuille Flow Derivation

For a viscous, incompressible, steady, laminar flow in a cylindrical tube, the Navier-Stokes equations reduce to a radial ODE:

2. The Fourth-Power Law

3. Blood as a Non-Newtonian Fluid

Whole blood does not obey Newton's law of viscosity (η = constant) — it is a shear-thinning suspension:

4. Vascular Tree and Murray's Law

5. Windkessel Model of Arterial Compliance

6. Turbulence and Heart Sounds

7. JavaScript Vascular Network Simulator

// Vascular network modelled as resistor network
// Nodes: pressure values; Edges: vessel segments with Poiseuille resistance

function poiseuille_resistance(radius, length, viscosity = 0.003) {
  // η_blood ≈ 3 mPa·s, returns resistance [Pa/(m³/s)]
  return (8 * viscosity * length) / (Math.PI * radius ** 4);
}

// Simple branching tree: heart → aorta → arteries[2] → arterioles[4]
function buildVascularTree(r0 = 0.015, L0 = 0.3, generations = 4) {
  const vessels = [];
  function branch(r, L, parentPressure, depth) {
    const R_h = poiseuille_resistance(r, L);
    vessels.push({ r, L, R_h, depth });
    if (depth < generations) {
      const r_child = r / Math.pow(2, 1/3); // Murray's law, bifurcation
      branch(r_child, L * 0.7, parentPressure, depth + 1);
      branch(r_child, L * 0.7, parentPressure, depth + 1);
    }
  }
  branch(r0, L0, 13330, 0); // 100 mmHg aortic pressure in Pa
  return vessels;
}

// Windkessel 2-element model
function windkessel({
  C = 7.5e-9,    // arterial compliance [m³/Pa]
  R = 1e8,       // peripheral resistance [Pa·s/m³]
  P0 = 10600,    // diastolic pressure [Pa] ≈ 80 mmHg
  HR = 70,       // heart rate [bpm]
  SV = 70e-6     // stroke volume [m³] = 70 mL
}, steps = 2000) {
  const dt = 0.001; // 1 ms
  const T = 60 / HR; // cardiac period [s]
  const t_sys = T * 0.35; // systolic duration ≈ 35% of cycle
  const Q_peak = SV / (t_sys * 0.5); // triangular ejection

  let P = P0, t = 0;
  const record = [];

  for (let i = 0; i < steps; i++) {
    const t_in_cycle = t % T;
    let Q_in = 0;
    if (t_in_cycle < t_sys) {
      // Triangular pulse: rise to Q_peak then fall
      Q_in = Q_peak * (1 - Math.abs(2*t_in_cycle/t_sys - 1));
    }
    const dP = (Q_in - P/R) / C;
    P += dP * dt;
    t += dt;
    if (i % 10 === 0) record.push({t: t.toFixed(3), P_mmHg: (P/133.3).toFixed(1)});
  }
  return record;
}

const pulse = windkessel({});
console.log(`Peak P: ${Math.max(...pulse.map(p => +p.P_mmHg)).toFixed(0)} mmHg (systolic)`);

8. Clinical Pathology