Simulate Microfluidic Cell Sorting in Matter.js

Simulating microfluidic cell sorting via hydrodynamic channel splitting in Matter.js requires approximating continuous fluid-particle interactions within a 2D rigid-body physics engine. Because Matter.js is designed for rigid bodies rather than computational fluid dynamics (CFD), hydrodynamic phenomena—such as laminar velocity profiles, Stokes' drag, and dividing streamlines at a channel bifurcation—must be modeled by applying programmatic forces to cell-like bodies at each engine update. This guide explains how to construct channel boundaries, define distinct cell properties, simulate laminar drag forces, and achieve particle separation at a fluidic junction.

1. Constructing the Channel Geometry

The microfluidic device consists of an inlet channel, a constriction or alignment zone, and a bifurcation dividing into two or more outlet channels. In Matter.js, channel walls are constructed using static rigid bodies.

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

const engine = Engine.create({ gravity: { x: 0, y: 0 } }); // Zero gravity
const world = engine.world;

// Channel dimensions
const channelWidth = 60;
const wallThickness = 10;

// Example: Static walls creating a main channel that branches into two symmetric outlets
const walls = [
  // Main channel walls
  Bodies.rectangle(150, 100 - channelWidth / 2, 300, wallThickness, { isStatic: true }),
  Bodies.rectangle(150, 100 + channelWidth / 2, 300, wallThickness, { isStatic: true }),
  
  // Bifurcation splitter (central wedge)
  Bodies.polygon(350, 100, 3, 30, { isStatic: true, angle: Math.PI / 6 }),
  
  // Upper outlet outer wall
  Bodies.rectangle(450, 50, 200, wallThickness, { isStatic: true, angle: -Math.PI / 12 }),
  // Lower outlet outer wall
  Bodies.rectangle(450, 150, 200, wallThickness, { isStatic: true, angle: Math.PI / 12 })
];

Composite.add(world, walls);

2. Defining Cell Particles

Cells are represented as circular dynamic bodies. Sorting in hydrodynamic channel splitting typically relies on size differences (pinched flow fractionation) or deformability. Differentiate cell types by radius and density:

function createCell(x, y, type) {
  const isLarge = type === 'large';
  const radius = isLarge ? 8 : 4; // Distinct hydrodynamic radii
  const mass = isLarge ? 2.0 : 0.5;

  return Bodies.circle(x, y, radius, {
    mass: mass,
    friction: 0,
    frictionAir: 0, // Handled manually via Stokes' drag
    restitution: 0.1,
    label: type
  });
}

3. Hydrodynamic Flow Field Approximation

In a microfluidic channel, Reynolds numbers are low (\(Re \ll 1\)), resulting in laminar, parabolic Poiseuille flow.

For a channel of width \(W\) centered at \(y_c\), the fluid velocity \(u(y)\) parallel to the channel axis is:

\[u(y) = u_{max} \left( 1 - \left( \frac{y - y_c}{W / 2} \right)^2 \right)\]

At the bifurcation, the fluid flux splits based on the hydraulic resistance of the outlet channels. Define an analytical vector field function \(\vec{v}_{\text{fluid}}(x, y)\) representing the flow velocity at any coordinate.

function getFluidVelocity(x, y) {
  const yc = 100;
  const halfWidth = channelWidth / 2;
  
  // Main straight channel
  if (x < 300) {
    const distFromCenter = Math.abs(y - yc);
    if (distFromCenter >= halfWidth) return { x: 0, y: 0 };
    const uMax = 4;
    const vx = uMax * (1 - Math.pow(distFromCenter / halfWidth, 2));
    return { x: Math.max(0, vx), y: 0 };
  }
  
  // Splitting region: streamline deflection
  const uOutlet = 3;
  if (y < yc) {
    // Upper branch vector (upward and rightward)
    return Vector.normalise({ x: 1, y: -0.5 });
  } else {
    // Lower branch vector (downward and rightward)
    return Vector.normalise({ x: 1, y: 0.5 });
  }
}

4. Applying Stokes' Drag Force

Instead of native physics engine damping, apply Stokes' drag force (\(\vec{F}_d = 6 \pi \mu R (\vec{v}_{\text{fluid}} - \vec{v}_{\text{cell}})\)) on every simulation tick via the beforeUpdate event.

const dynamicViscosity = 0.05;

Events.on(engine, 'beforeUpdate', () => {
  const cells = Composite.allBodies(world).filter(b => !b.isStatic);

  cells.forEach(cell => {
    const vFluid = getFluidVelocity(cell.position.x, cell.position.y);
    const radius = cell.circleRadius || 4;

    // Relative velocity between fluid and cell
    const vRelX = vFluid.x - cell.velocity.x;
    const vRelY = vFluid.y - cell.velocity.y;

    // Stokes' drag: F = 6 * pi * dynamicViscosity * R * vRel
    const dragCoeff = 6 * Math.PI * dynamicViscosity * radius;
    const force = {
      x: dragCoeff * vRelX * 0.001, // Scaled for Matter.js engine units
      y: dragCoeff * vRelY * 0.001
    };

    Body.applyForce(cell, cell.position, force);
  });
});

5. Achieving Sorting via Hydrodynamic Alignment

To separate particles at the junction:

  1. Focusing/Pinching: Align all incoming cells against one wall using a sheath flow (simulated via an asymmetric velocity field or an initial offset position).
  2. Radial Offset: Because the rigid boundary prevents a cell's center of mass from getting closer to the wall than its radius \(R\), larger cells occupy streamlines further toward the channel center than smaller cells.
  3. Streamline Splitting: Position the bifurcation splitter so that the dividing streamline passes between the center-of-mass positions of the aligned large and small cells. Large cells are carried across the dividing line into the outer stream, while small cells remain in the wall-adjacent stream.