Simulating Spinal Compression in Matter.js

This article explains how to model biomechanical spinal disc compression under heavy weightlifting loads using the Matter.js 2D physics engine. By representing vertebrae as rigid rectangular bodies and intervertebral discs as tuned, damped spring constraints, you can accurately simulate the axial deformation, shear stress, and viscoelastic behavior of a human spine under dynamic mechanical load.

1. Conceptual Physics Model

A realistic spinal simulation in a 2D rigid-body engine requires breaking the spine into individual functional spinal units (FSUs). Each unit consists of:

2. Setting Up Vertebrae Bodies

Define vertebrae with realistic proportional dimensions and mass. A heavier mass simulates the cumulative mass of the torso, while an external mass represents a barbell.

const { Engine, Render, Runner, Bodies, Composite, Constraint, Body, Vector } = Matter;

const engine = Engine.create();
const world = engine.world;

// Vertebra properties
const vertebraWidth = 40;
const vertebraHeight = 15;
const vertebraMass = 2; // kg representation

const vertebrae = [];
const numVertebrae = 5; // e.g., Lumbar spine (L1 to L5)
const startY = 300;
const startX = 200;
const discRestLength = 8;

for (let i = 0; i < numVertebrae; i++) {
  const vertebra = Bodies.rectangle(
    startX,
    startY - i * (vertebraHeight + discRestLength),
    vertebraWidth,
    vertebraHeight,
    {
      mass: vertebraMass,
      frictionAir: 0.02,
      render: { fillStyle: '#e0e0e0' }
    }
  );
  vertebrae.push(vertebra);
  Composite.add(world, vertebra);
}

// Anchor the base (Sacrum)
Body.setStatic(vertebrae[0], true);

3. Implementing Viscoelastic Discs

A single constraint at the center creates a hinge, failing to model resistance to bending. To properly simulate disc compression and rotational stiffness, connect each adjacent pair of vertebrae with two parallel constraints (anterior and posterior):

const discStiffness = 0.08; // Modulates axial resistance
const discDamping = 0.1;     // Viscous damping to dissipate energy

for (let i = 0; i < numVertebrae - 1; i++) {
  const lower = vertebrae[i];
  const upper = vertebrae[i + 1];
  const halfWidth = vertebraWidth / 2 - 4;

  // Left constraint (Anterior longitudinal support)
  const leftConstraint = Constraint.create({
    bodyA: lower,
    pointA: { x: -halfWidth, y: -vertebraHeight / 2 },
    bodyB: upper,
    pointB: { x: -halfWidth, y: vertebraHeight / 2 },
    stiffness: discStiffness,
    damping: discDamping,
    render: { strokeStyle: '#ff4d4d', lineWidth: 3 }
  });

  // Right constraint (Posterior longitudinal support)
  const rightConstraint = Constraint.create({
    bodyA: lower,
    pointA: { x: halfWidth, y: -vertebraHeight / 2 },
    bodyB: upper,
    pointB: { x: halfWidth, y: vertebraHeight / 2 },
    stiffness: discStiffness,
    damping: discDamping,
    render: { strokeStyle: '#ff4d4d', lineWidth: 3 }
  });

  Composite.add(world, [leftConstraint, rightConstraint]);
}

4. Applying Heavy Weightlifting Loads

To simulate actions like squats or deadlifts, load the topmost vertebra axially with additional mass or apply a downward vector force.

// Adding a heavy barbell load to the top vertebra
const topVertebra = vertebrae[vertebrae.length - 1];

const barbell = Bodies.circle(startX, topVertebra.position.y - 30, 20, {
  mass: 100, // Heavy external load
  render: { fillStyle: '#333333' }
});

// Rigidly couple barbell to top vertebra
const loadConstraint = Constraint.create({
  bodyA: topVertebra,
  pointA: { x: 0, y: -vertebraHeight / 2 },
  bodyB: barbell,
  pointB: { x: 0, y: 0 },
  stiffness: 0.9
});

Composite.add(world, [barbell, loadConstraint]);

Alternatively, simulate dynamic ground-reaction acceleration or axial trunk muscle co-contraction using Body.applyForce:

// Simulate sudden dynamic load phase
Events.on(engine, 'beforeUpdate', () => {
  Body.applyForce(topVertebra, topVertebra.position, { x: 0, y: 0.05 });
});

5. Calculating Disc Compression and Strain Metrics

Measure compression quantitatively by monitoring the relative Euclidean distance between neighboring vertebrae on every engine tick.

function getDiscCompression(lowerBody, upperBody, initialSpacing) {
  const currentDistance = Vector.magnitude(
    Vector.sub(upperBody.position, lowerBody.position)
  );
  const currentDiscHeight = currentDistance - vertebraHeight;
  const deformation = initialSpacing - currentDiscHeight;
  const strain = deformation / initialSpacing; // Normal strain: ΔL / L0

  return { deformation, strain };
}

// Read compression of the L5-S1 equivalent joint
Events.on(engine, 'afterUpdate', () => {
  const { deformation, strain } = getDiscCompression(vertebrae[0], vertebrae[1], discRestLength);
  
  if (strain > 0.4) {
    // Disc failure or excessive compressive strain threshold
    console.warn(`Critical disc compression: ${(strain * 100).toFixed(1)}%`);
  }
});

6. Tuning Engine Solver Iterations

Matter.js uses iterative impulse resolution. Heavy masses resting on spring constraints cause joint separation and numerical instability unless engine iterations are increased:

engine.positionIterations = 12;
engine.velocityIterations = 8;
engine.constraintIterations = 10;

Increasing constraint iterations ensures the paired springs properly transfer compressive forces down through the entire kinetic chain without jittering or unrealistically flattening.