Skip to content

Latest commit

 

History

History
1466 lines (961 loc) · 92.5 KB

File metadata and controls

1466 lines (961 loc) · 92.5 KB

MIME Node Taxonomy

Scientific and engineering rationale for each category of physics node in the MIcrorobotics Multiphysics Engine

Purpose: This document is handed to a literature research agent that will identify specific governing equations, discretisation schemes, and validated parameter regimes for each node type. It provides the physical context, input/output relationships, and known challenges — enough for a researcher to know what to look for in the literature and why it matters for microrobotics simulation.

Scope: This document covers the categories of physics nodes that MIME needs. It does not name specific node classes, define APIs, or specify implementation details. Those are downstream of the literature research.

Key differentiators of the MIME approach: MIME's node taxonomy reflects two commitments that distinguish it from general robotics simulation frameworks. First, the physics nodes are designed for regulatory traceability: every node carries validated regime bounds (what conditions it has been verified in), hazard hints (what can go wrong), and biocompatibility metadata (what the robot is made of and what that means biologically). This is not documentation overhead — it is the evidence chain required for IEC 62304 SOUP assessment and ISO 14971 hazard identification when MIME is incorporated into a regulated medical device. Second, because MIME is built on MADDENING's differentiable physics graph, every node's update function is JAX-traceable. This means the simulation is differentiable where the underlying physics permits — nodes modelling well-behaved continuous physics (magnetic response, rigid body dynamics in the Stokes regime, diffusion, drug release kinetics) are fully differentiable, while nodes modelling discontinuous or stiff physics are explicitly flagged as differentiability-limited throughout this document. Gradient-based controller training, sensitivity analysis, and uncertainty propagation are first-class capabilities, not afterthoughts. A trajectory dataset (the LeRobot model) cannot support this; a static FEM simulation (the COMSOL model) cannot support this. The combination of medical-grade regulatory infrastructure and differentiable physics is MIME's core scientific contribution.


Table of Contents

  1. External Apparatus Nodes
  2. Robot Body Physics Nodes
  3. Fluid Environment Nodes
  4. Biological Tissue Nodes
  5. Therapeutic Payload Nodes
  6. Sensing / Imaging Physics Nodes
  7. Cross-Cutting Concerns
  8. Phase Assignment Summary

1. External Apparatus Nodes

These nodes model hardware that sits outside the patient's body and generates the fields or energy that actuate the microrobot. They are the physical interface between the control system and the in-vivo environment. In the MIME architecture, external apparatus nodes are what the ControlPolicy commands — they translate control signals into physical fields that propagate into the body.

1.1 Magnetic Field Sources

What it models

The spatial and temporal structure of magnetic fields generated by hardware outside the body. Three main hardware configurations are relevant:

Rotating permanent magnets: A permanent magnet (typically NdFeB, diameter 10–50 mm) mounted on a motor at a controlled distance from the patient. Generates a rotating dipole field. The field strength decays as 1/r^3. The rotation frequency is the primary control input. This is the simplest and most common actuation hardware for helical microrobots.

Helmholtz/Maxwell coil arrays: Electromagnetic coil pairs (or arrays of 3–8 coils) that generate a uniform or near-uniform field over a working volume. Helmholtz pairs produce uniform fields; Maxwell pairs produce uniform gradients. Combined arrays can independently control field direction, magnitude, and gradient — enabling 5-DOF (or higher) control of the microrobot. Current in each coil is the control input.

Gradient systems: Dedicated gradient coils (similar to MRI gradient coils) that produce strong, spatially varying fields for pulling forces. Distinguished from Helmholtz/Maxwell arrays by higher gradient strength and potentially faster slew rates. Important for gradient-based steering where the robot has no onboard propulsion and moves by magnetic gradient pulling.

Why it is needed

  • Magnetic actuation is the primary propulsion mechanism for most microrobot designs in CSF. Without an accurate external field model, the torque and force on the robot cannot be computed.
  • The field spatial structure matters: a dipole field from a rotating magnet is highly non-uniform, and the robot's behaviour depends strongly on where in the field it sits. A Helmholtz coil produces a uniform field where the robot's position within the working volume doesn't matter (to first order).
  • Separating the external apparatus from the robot response enables: calibration simulations (characterise the field without the robot), multi-robot simulations (N robots share one field), and uncertainty injection at the correct physical boundary.

Inputs and outputs

  • Inputs (from ControlPolicy via ExternalInputSpec): rotation frequency, field magnitude, coil currents, gradient strengths, field direction angles
  • Outputs (to robot body via edges): magnetic field vector B(x,t) at the robot's position, field gradient tensor dB/dx at the robot's position
  • Parameters: magnet geometry, coil geometry, magnet-to-patient distance, magnet dipole moment, coil turns/radius

Key physical parameters

  • Dipole moment of the permanent magnet (A·m^2)
  • Coil current amplitudes and phases (A)
  • Working distance (mm) — strongly affects field strength due to 1/r^3 decay
  • Working volume dimensions (mm) — region over which field is acceptably uniform
  • Maximum field strength (mT) and gradient (T/m)
  • Maximum rotation frequency (Hz) — limited by motor or power amplifier
  • Field homogeneity (%) over the working volume

Known challenges

  • Near-field vs. far-field: close to the magnet, the dipole approximation breaks down and full multipole expansion or numerical (FEM) field solutions are needed. Literature should identify the distance-to-size ratio at which the dipole approximation is valid.
  • Eddy currents and shielding: conductive biological tissue and metallic implants distort the applied field. Not typically modelled in initial simulations but relevant for clinical translation.
  • Temporal bandwidth: coil arrays have inductance-limited bandwidth. The field cannot change instantaneously — there is a slew rate limit. This matters for high-frequency actuation and fast gradient steering.
  • Mutual inductance in coil arrays: currents in one coil affect the field of others. The inverse problem (desired field → required currents) requires the full mutual inductance matrix.

Differentiability status: [CONDITIONALLY DIFFERENTIABLE] — differentiability depends on the specific discretisation and parameter regime. Assess per implementation.

Relationship to other node categories

  • Robot body (magnetic response): the external field is the primary input to the magnetic response calculation
  • Sensing (MRI): the external actuation field may interfere with MRI imaging fields — this cross-talk is a known challenge for MRI-guided microrobotics
  • Control: the external apparatus is the actuator in the control loop

Phase assignment

  • Rotating permanent magnet model: Phase 1 (essential for basic helical robot simulation)
  • Helmholtz/Maxwell coil array model: Phase 1 (essential for controlled-field experiments)
  • Gradient steering model: Phase 2 (needed for gradient-based actuation modalities)

Existing implementations to examine

  • OctoMag system (ETH Zurich) — 8-coil electromagnetic manipulation system with published field models
  • MagnetoSuture (Martel group) — rotating permanent magnet models
  • OpenMagnetics — open-source magnetic field computation library

Implemented nodes (current)

Algorithm ID Node Role Notes
MIME-NODE-001 ExternalMagneticFieldNode Helmholtz / far-field Uniform 2D rotating field; appropriate when field_gradient = 0 is acceptable AND the uniform-field assumption holds in the workspace
MIME-NODE-100 MotorNode Magnet rotor stage Single-axis rotary motor (DC brushed); torque / voltage / velocity command modes; semi-implicit Euler
MIME-NODE-101 PermanentMagnetNode Bar-magnet field producer point_dipole / current_loop / coulombian_poles field models; jax.jacrev ∇B; Earth-field superposition
MIME-NODE-102 RobotArmNode URDF-driven articulated arm FK + CRBA + RNEA + forward dynamics, all pure JAX

The Motor + PermanentMagnet (+ optional RobotArm) chain is the right choice when (a) field gradient matters, (b) the magnet has a finite or tracked physical pose, (c) misalignment / wobble effects are under study, or (d) the demo scene needs a rendered apparatus. ExternalMagneticFieldNode remains the right choice when none of those apply — it has lower setup overhead and is a faithful model of a Helmholtz coil pair driven in quadrature.


1.2 Acoustic Sources

What it models

The acoustic pressure field generated by ultrasound transducers or transducer arrays outside the body. Two configurations:

Focused ultrasound transducers: Single-element or phased-array transducers that generate a focused acoustic beam at a target point inside the body. The pressure field has a focal zone where the intensity is highest. Used for acoustic streaming (steady flow generated by sound absorption) and radiation force.

Acoustic tweezers: Arrays of transducers that create structured pressure fields (standing waves, vortex beams) capable of trapping and manipulating microparticles. Use interference patterns to create potential wells.

Why it is needed

  • Acoustic actuation is the primary alternative to magnetic actuation, especially for robots without magnetic components
  • Ultrasound can penetrate deep tissue (unlike optical methods) and does not require line-of-sight
  • Acoustic streaming generates bulk fluid flow that can propel robots even without onboard response elements
  • Radiation force from focused ultrasound can push or trap microrobots
  • Acoustic actuation can trigger drug release from microbubble-coated robots

Inputs and outputs

  • Inputs (from ControlPolicy): transducer drive voltage/frequency, phased-array element phases, focal point coordinates
  • Outputs (to robot body via edges): acoustic pressure p(x,t) at robot position, radiation force vector, streaming velocity field
  • Parameters: transducer geometry, operating frequency (MHz), focal length, aperture diameter, medium acoustic impedance

Key physical parameters

  • Operating frequency (typically 1–20 MHz for therapeutic ultrasound)
  • Acoustic pressure at focus (MPa) — determines radiation force and streaming
  • Mechanical index (MI) — safety parameter, MI < 1.9 required for diagnostic use
  • Thermal index (TI) — tissue heating concern
  • Focal zone dimensions (mm) — axial and lateral resolution
  • Attenuation coefficient of tissue/fluid (dB/cm/MHz)

Known challenges

  • Skull/bone transmission: acoustic propagation through bone (relevant for transcranial CSF access) causes severe aberration and attenuation. Phase correction using CT-derived skull models is an active research area.
  • Nonlinear propagation: at high pressures, the acoustic wave distorts nonlinearly (shock formation). Linear models break down.
  • Scattering from microstructures: the microrobot itself, blood cells, and tissue microstructure scatter the acoustic field. Multiple scattering is computationally expensive.
  • Standing wave artefacts: reflections from tissue interfaces create unintended standing waves that may trap the robot at nodes rather than the intended focal point.

Differentiability status: [CONDITIONALLY DIFFERENTIABLE] — differentiability depends on the specific discretisation and parameter regime. Assess per implementation.

Relationship to other node categories

  • Robot body (microbubble dynamics): the acoustic pressure field drives bubble oscillation
  • Robot body (rigid body): radiation force contributes to robot motion
  • Environment (fluid flow): acoustic streaming generates fluid flow that interacts with the ambient CSF flow
  • Therapeutic (drug release): ultrasound can trigger release from acoustically-responsive carriers

Phase assignment

  • Phase 2 (important for non-magnetic actuation modalities)

Existing implementations to examine

  • k-Wave — open-source MATLAB/C++ acoustic simulation toolbox (widely used in focused ultrasound research)
  • FOCUS (Fast Object-oriented C++ Ultrasound Simulator)
  • BabelBrain — transcranial focused ultrasound simulation

1.3 Optical Sources

What it models

Optical trapping forces and photoacoustic excitation from laser sources. Two mechanisms:

Optical tweezers: Tightly focused laser beams that trap dielectric microparticles via gradient force. The trap stiffness depends on the laser power, beam waist, and particle refractive index contrast.

Photoacoustic excitation: Pulsed laser light absorbed by the microrobot or a coating generates thermoelastic expansion that produces an acoustic pulse, which in turn generates a propulsive force or triggers drug release.

Why it is needed

  • Optical tweezers provide extremely precise (nm-resolution) positioning and force measurement for microrobots in vitro
  • Photoacoustic mechanisms enable light-triggered drug release and propulsion
  • Essential for in-vitro validation and characterisation, even if not applicable in vivo at depth

Inputs and outputs

  • Inputs: laser power, wavelength, focal position, pulse parameters (for photoacoustic)
  • Outputs: trapping force at robot position, photoacoustic pressure pulse
  • Parameters: numerical aperture, beam waist, wavelength, pulse duration, absorption coefficient

Key physical parameters

  • Laser wavelength (typically 532 nm, 808 nm, or 1064 nm)
  • Laser power (mW to W)
  • Trap stiffness (pN/nm) — depends on power, NA, and particle properties
  • Penetration depth in tissue (mm) — severely limited in vivo (< 1 mm at visible wavelengths, up to 5 mm in the near-IR window)

Known challenges

  • Penetration depth: optical methods cannot reach targets deeper than a few millimetres in tissue, severely limiting in-vivo applicability for CSF targets
  • Phototoxicity: prolonged laser exposure damages cells
  • Scattering: tissue is highly scattering at optical wavelengths, defocusing the beam

Differentiability status: [CONDITIONALLY DIFFERENTIABLE] — differentiability depends on the specific discretisation and parameter regime. Assess per implementation.

Relationship to other node categories

  • Robot body: optical trapping force contributes to motion
  • Therapeutic: photoacoustic excitation can trigger release

Phase assignment

  • Advanced (relevant for in-vitro validation setups; limited in-vivo applicability for CSF targets)

Existing implementations to examine

  • Optical Tweezers Toolbox (Alexander et al.) — MATLAB/Python T-matrix optical force computation
  • Lumerical FDTD — commercial optical simulation (for reference)

2. Robot Body Physics Nodes

These nodes model the physics of the microrobot itself — its motion, its response to external fields, its structural mechanics, and its interaction with boundaries. The microrobot is the central physical object in the simulation.

2.1 Rigid Body Dynamics

What it models

The translational and rotational motion of a rigid microrobot in a viscous fluid, under the combined effects of all applied forces and torques. At the microscale (characteristic length 10–500 um, Reynolds number 10^-4 to 10^-1), inertia is negligible compared to viscous drag — the equations of motion are overdamped.

The state is the robot's position (3D) and orientation (quaternion or rotation matrix). In the low-Re regime, the velocity is instantaneously proportional to the applied force (no acceleration phase). The dynamics are:

  • Translation: F_total = F_magnetic + F_drag + F_gravity + F_buoyancy + F_acoustic + F_contact + F_Brownian
  • Rotation: T_total = T_magnetic + T_drag + T_contact + T_Brownian

At low Re, F_drag = -R_T · v and T_drag = -R_R · omega, where R_T and R_R are the translational and rotational resistance tensors (shape-dependent).

Why it is needed

This is the most fundamental node in the entire simulation. Every microrobot simulation requires tracking the robot's position and orientation. All other physics (fluid coupling, drug delivery, sensing) reference the robot's pose.

  • Answers: "Where is the robot? How fast is it moving? What is its orientation?"
  • Enables: trajectory prediction, navigation planning, closed-loop control evaluation

Inputs and outputs

  • Inputs (from other nodes via edges): magnetic force/torque, drag force/torque, contact force/torque, gravity/buoyancy, Brownian force/torque, acoustic radiation force
  • Outputs (to other nodes via edges): position, velocity, orientation, angular velocity
  • Parameters: robot mass, moment of inertia tensor, characteristic length, shape (for resistance tensor), density

Key physical parameters

  • Characteristic length (um) — determines Re and the relative importance of forces
  • Mass and density (kg, kg/m^3) — often close to fluid density for neutrally buoyant designs
  • Shape descriptor — sphere, prolate ellipsoid, helical, cylindrical — determines resistance tensor
  • Resistance tensor components (N·s/m for translation, N·m·s for rotation) — shape-dependent, from Stokes flow solutions or slender body theory

Known challenges

  • Resistance tensor for complex shapes: for a helical robot, the resistance tensor couples translation and rotation (propulsion matrix). Analytical solutions exist for simple helices (Lighthill slender body theory, resistive force theory) but not for arbitrary geometries. Pre-computed tensors from Stokes flow BEM solvers may be needed.
  • Orientation representation: quaternions avoid gimbal lock but require normalisation at each timestep. Rotation matrices are redundant (9 params for 3 DOF). Rodrigues vectors have singularities. Literature should identify the best representation for JAX-traceable rigid body dynamics.
  • Near-wall hydrodynamics: the resistance tensor changes near boundaries (wall correction factors). The free-space Stokes resistance is only valid far from walls.
  • Brownian motion: at the microscale, thermal fluctuations are significant. Brownian forces scale as sqrt(k_B T / dt) and must be included for accurate trajectory statistics, but they are stochastic and complicate JAX traceability (requires explicit RNG key management).
  • Overdamped vs. underdamped: most microrobot simulations assume zero Reynolds number (overdamped). But at the upper end of the size range (hundreds of um), finite-Re effects may matter. Need to identify the crossover.

Differentiability status: [DIFFERENTIABLE] — gradients are reliable for policy training, sensitivity analysis, and uncertainty propagation within the validated physical regime. Caveat: near-wall corrections involving position-dependent resistance tensors require smooth interpolation of the correction factors to avoid gradient discontinuities at the wall-correction onset distance.

Relationship to other node categories

  • Magnetic response: receives torque and force
  • Fluid environment: bidirectional — sends position/velocity to flow solver, receives drag force/torque
  • Surface contact: receives contact forces at boundaries
  • Sensing: position and orientation are the primary quantities observed by the imaging system
  • Phase tracking: orientation is read by the phase observer
  • Therapeutic: position determines where drug is released

Phase assignment

  • Phase 1 (essential — every simulation needs this)

Existing implementations to examine

  • Resistive Force Theory (Lighthill 1976, Gray & Hancock 1955) — analytical drag for slender helices
  • Regularised Stokeslet method (Cortez 2001) — numerical Stokes flow for arbitrary shapes
  • MagneBot (Nelson group) — published microrobot dynamics models
  • PyStokes (Singh & Adhikari) — Python library for Stokes flow and microswimmer dynamics

2.2 Magnetic Response

What it models

The interaction between the microrobot's onboard magnetic element (permanent magnet, soft-magnetic material, or superparamagnetic particles) and the externally applied magnetic field. This interaction produces a torque and/or force on the robot.

Permanent magnet response: The robot has a fixed magnetic dipole moment m. In an external field B, the torque is T = m x B (tends to align m with B). In a field gradient, the force is F = nabla(m · B). For a helical robot, the torque drives rotation about the helix axis, which generates propulsion via the shape-coupling in the resistance tensor.

Soft-magnetic response: The robot's magnetisation depends on the applied field (susceptibility chi). The induced moment is m = chi · V · B / mu_0. The torque and force depend on the field and its gradient, but the moment itself changes as the field changes.

Superparamagnetic response: Similar to soft-magnetic but with Langevin-function saturation at high fields and zero remanence (no permanent moment).

Why it is needed

  • Converts the external field (from the apparatus node) into the forces and torques that drive robot motion
  • The step-out phenomenon (when the robot can't keep up with the rotating field) is fundamentally a magnetic response effect — the viscous drag torque exceeds the maximum magnetic torque
  • Different magnetic materials have qualitatively different responses (permanent vs. soft vs. superparamagnetic) and the simulation must capture this

Inputs and outputs

  • Inputs: external field B(x,t) and gradient dB/dx from apparatus node; robot orientation from rigid body node
  • Outputs: magnetic torque and force vectors to rigid body node; magnetisation state (for sensing/imaging interaction)
  • Parameters: magnetic dipole moment magnitude and direction (body frame), susceptibility tensor, saturation magnetisation, Langevin parameter, remanent magnetisation

Key physical parameters

  • Permanent magnet moment (A·m^2) — typically 10^-12 to 10^-9 A·m^2 for microrobots
  • Susceptibility chi (dimensionless) — for soft-magnetic materials
  • Saturation magnetisation M_s (A/m) — maximum achievable magnetisation
  • Curie temperature (if thermal effects matter)
  • Demagnetisation factors (shape-dependent)

Known challenges

  • Demagnetisation: the internal field in the robot differs from the external field due to shape-dependent demagnetisation. For non-ellipsoidal shapes, the demagnetisation tensor is non-uniform.
  • Hysteresis: some magnetic materials exhibit hysteresis (the magnetisation depends on field history). This matters for soft-magnetic robots that undergo repeated cycling.
  • Dipole-dipole interaction: for multi-robot systems, the magnetic field of one robot affects another. This is an N^2 interaction that becomes important for swarm configurations.
  • Thermal demagnetisation: body temperature (37°C) may reduce the magnetic moment of some materials

Differentiability status: [DIFFERENTIABLE] — gradients are reliable for policy training, sensitivity analysis, and uncertainty propagation within the validated physical regime.

Relationship to other node categories

  • External apparatus: receives field and gradient
  • Rigid body: sends torque and force, receives orientation
  • Phase tracking: magnetisation direction determines phase relative to external field
  • Sensing (MRI): magnetic moment creates susceptibility artefacts in MRI images

Phase assignment

  • Permanent magnet response: Phase 1 (essential for basic helical robot)
  • Soft-magnetic response: Phase 2
  • Superparamagnetic response: Phase 2

Existing implementations to examine

  • Magpylib — Python library for magnetic field calculation from permanent magnets and coils
  • OOMMF/MuMax3 — micromagnetic simulation (overkill for rigid-body level, but reference for material models)

2.3 Elastic/Flexible Body Mechanics

What it models

The deformation of a flexible microrobot under applied loads. Relevant for:

Flagellar robots: artificial bacterial flagella (ABFs) are helical filaments that are not perfectly rigid. At high rotation frequencies or under strong fluid loads, they deform — changing the pitch, amplitude, and propulsion characteristics.

Compliant robots: some designs use flexible joints or compliant mechanisms to convert magnetic torque into locomotion (e.g., a multi-link swimmer where each link is a permanent magnet connected by flexible hinges).

Elastic capsules: drug-carrying microrobots with elastic shells that deform under flow shear or acoustic pressure.

Why it is needed

  • Rigid body assumption fails when the robot's shape changes under load
  • Flagellar deformation affects propulsion efficiency: a deformed helix generates different thrust than a rigid one
  • Compliant mechanisms are an active design space for microrobots
  • Drug release from elastic capsules depends on the deformation state (stress-triggered release)

Inputs and outputs

  • Inputs: applied forces/torques along the body, fluid shear stress, internal magnetic torques at joints
  • Outputs: deformed shape (centreline or surface mesh), effective resistance tensor (if feeding back to rigid body dynamics), stress distribution (for release triggering)
  • Parameters: Young's modulus, bending stiffness (EI), torsional stiffness (GJ), material density, Poisson's ratio, joint geometry

Key physical parameters

  • Bending stiffness EI (N·m^2) — determines how much the helix bends under fluid load
  • Sperm number Sp = L(omega * zeta_perp / EI)^(1/4) — dimensionless parameter governing the transition from rigid to flexible behaviour. Sp << 1 means rigid-like; Sp >> 1 means highly flexible.
  • Number of links/joints (for multi-link swimmers)
  • Joint stiffness (N·m/rad) for compliant mechanisms

Known challenges

  • Large deformations: microrobot flagella can undergo large (finite) deformations, requiring geometrically nonlinear elasticity rather than small-strain theory
  • Fluid-structure interaction: the deformation depends on the flow, which depends on the deformed shape — a two-way coupled problem. This requires coupling group iteration in MADDENING's graph framework.
  • Discretisation of slender bodies: Kirchhoff rod theory or Cosserat rod theory for the flagellum centreline; shell theory for capsules. Different discretisation approaches (finite elements, discrete elastic rods, Euler-Bernoulli beam elements) have different tradeoffs.
  • Stability: explicit time integration of flexible bodies is conditionally stable; the timestep is limited by the stiffest bending mode. Implicit integration may be needed.

Differentiability status: [CONDITIONALLY DIFFERENTIABLE] — differentiability depends on the specific discretisation and parameter regime. Assess per implementation.

Relationship to other node categories

  • Rigid body: may replace or couple with rigid body dynamics — a deformable robot has different effective resistance than a rigid one
  • Fluid environment: bidirectional FSI coupling
  • Therapeutic: stress-triggered drug release depends on deformation state
  • Magnetic response: distributed magnetisation along a flexible body means the torque varies along the length

Phase assignment

  • Phase 2 (important for flagellar and compliant robot designs)

Existing implementations to examine

  • Elastica (Gazzola group) — Python/C++ Cosserat rod simulator for soft robots
  • SOFA Framework — soft body physics simulation (medical/surgical focus)
  • Discrete Elastic Rods (Bergou et al. 2008, 2010) — widely used for slender bodies

2.4 Acoustic Bubble Dynamics

What it models

The oscillation dynamics of a gas microbubble that is part of the microrobot's body or coating. When excited by an external ultrasound field, the bubble undergoes radial oscillation described by the Rayleigh-Plesset equation (or its extensions). This oscillation generates:

  • Microstreaming: steady fluid flow around the oscillating bubble, which propels the robot
  • Radiation force: the bubble scatters the incident acoustic field, experiencing a net force
  • Sonochemical effects: at high oscillation amplitudes, inertial cavitation can generate extreme local temperatures and pressures (relevant for drug release)

Why it is needed

  • Acoustic bubble-based propulsion is a major alternative to magnetic actuation
  • Microbubble coatings are used for ultrasound contrast (imaging) and as drug carriers
  • The bubble response is highly nonlinear and frequency-dependent — the propulsion characteristics depend sensitively on the driving frequency relative to the bubble's resonance frequency

Inputs and outputs

  • Inputs: acoustic pressure p(t) from external acoustic node; ambient pressure from fluid environment; temperature
  • Outputs: bubble radius R(t), radial velocity dR/dt, microstreaming force, radiation force
  • Parameters: equilibrium bubble radius R_0, gas type (air, SF6, C3F8), shell properties (viscosity, elasticity for coated bubbles), surrounding fluid properties

Key physical parameters

  • Equilibrium radius R_0 (um) — determines resonance frequency (Minnaert frequency ~ 3/R_0 MHz for um-scale bubbles in water)
  • Shell elasticity chi_s (N/m) and viscosity kappa_s (Pa·s·m) — for lipid/polymer-coated bubbles
  • Driving pressure amplitude (kPa to MPa)
  • Polytropic exponent of the gas
  • Surface tension (N/m)

Known challenges

  • Stiff ODE: the Rayleigh-Plesset equation is stiff, especially near collapse. Explicit integration requires extremely small timesteps. Implicit integration (via MADDENING's implicit_residual()) is strongly recommended.
  • Nonlinear resonance: the bubble's resonance frequency shifts with amplitude (nonlinear softening/hardening). Linear models underpredict the response at moderate driving pressures.
  • Coated bubble models: lipid and polymer coatings add viscoelastic shell terms that modify the dynamics significantly. Multiple competing shell models exist in the literature (Marmottant, de Jong, Church-Hoff).
  • Collapse and rebound: inertial cavitation involves violent bubble collapse followed by rebound. The physics near collapse involves extreme conditions (T > 5000 K, p > 100 MPa inside the bubble) that push models to their limits.
  • Coupling to flow: the oscillating bubble generates a time-varying flow field around it. Computing this flow and feeding it back to the rigid body motion requires careful multi-timescale handling (bubble oscillates at MHz, robot translates at Hz).

Differentiability status: [DIFFERENTIABILITY-LIMITED] — Rayleigh-Plesset dynamics near bubble collapse requires implicit solvers; differentiating through the implicit Newton loop is computationally punishing and often numerically unstable at collapse. Gradients through this node should not be used for policy training or sensitivity analysis. Use in inference/forward simulation only, or apply a smoothed surrogate if differentiability is required.

Relationship to other node categories

  • External acoustic: receives driving pressure
  • Rigid body: streaming force and radiation force contribute to robot motion
  • Fluid environment: bubble oscillation disturbs the local flow field
  • Therapeutic: cavitation can trigger drug release; bubble dynamics determines release timing

Phase assignment

  • Phase 2 (important for acoustic actuation modality)

Existing implementations to examine

  • BubbleDynamics (Prosperetti) — reference implementations of Rayleigh-Plesset variants
  • COMSOL bubble dynamics module — commercial FEM reference
  • Marmottant shell model (2005) — widely used for coated microbubbles

2.5 Surface Interaction / Contact Mechanics

What it models

The mechanical interaction between the microrobot and solid boundaries — vessel walls, tissue surfaces, channel walls, other robots. At the microscale, surface interactions include:

  • Hydrodynamic wall effects: the drag on a sphere near a wall increases (lubrication theory). The resistance tensor becomes position-dependent near boundaries.
  • Contact/collision: when the robot touches a wall, a repulsive contact force prevents penetration. Models range from simple penalty methods (spring-like repulsion) to adhesion models (JKR, DMT).
  • Adhesion: van der Waals, electrostatic, and specific (ligand-receptor) adhesion forces can cause the robot to stick to walls. Rolling and sliding adhesion affect navigation along vessel walls.
  • Friction: at the microscale, friction is dominated by adhesion rather than gravity. The Amontons-Coulomb model may not apply; adhesion-based friction models are more appropriate.

Why it is needed

  • Microrobots in CSF channels must navigate without getting stuck to walls
  • Contact with vessel walls is inevitable in narrow channels (ventricles, aqueducts)
  • Wall effects change the drag coefficient by factors of 2–10× at separations comparable to the robot size
  • Adhesion is a major failure mode: a robot that adheres to the vessel wall is lost
  • Navigation strategies (wall-following, channel centring) depend on accurate wall-interaction models

Inputs and outputs

  • Inputs: robot position/velocity (from rigid body), wall geometry (from environment mesh), material properties
  • Outputs: contact force, friction force, adhesion force (to rigid body as additive boundary inputs)
  • Parameters: gap distance, wall material, surface energy, Hamaker constant, surface charge density, adhesion ligand density

Key physical parameters

  • Gap distance h (um) — determines the strength of wall effects
  • Hamaker constant A_H (J) — van der Waals interaction strength (typically 10^-20 to 10^-19 J for biological surfaces)
  • Surface charge / zeta potential (mV) — electrostatic interaction (DLVO theory)
  • Contact stiffness k (N/m) — for penalty-based contact
  • Adhesion energy (J/m^2)

Known challenges

  • Contact detection: determining when and where the robot contacts the wall requires distance computation between the robot geometry and the wall mesh. For complex geometries (helical robots in curved channels), this is non-trivial.
  • Stiffness of contact: penalty-based contact introduces stiff forces that require small timesteps or implicit integration. The penalty parameter must be tuned to balance stiffness (preventing penetration) against numerical stability.
  • Adhesion hysteresis: JKR and DMT models involve snap-in/snap-off transitions that are discontinuous — problematic for JAX traceability. Smoothed approximations are needed.
  • Surface roughness: at the microscale, surface roughness affects both hydrodynamic lubrication and contact mechanics. Perfectly smooth surfaces are unrealistic.
  • Specific adhesion: ligand-receptor bonds (relevant for targeted drug delivery) are stochastic and involve binding/unbinding kinetics

Differentiability status: [DIFFERENTIABILITY-LIMITED] — penalty-based contact forces are discontinuous; jax.grad through a collision event produces unreliable gradients unless expensive smooth approximations are used. Gradients through this node should not be used for policy training or sensitivity analysis. Use in inference/forward simulation only, or apply a smoothed surrogate if differentiability is required.

Relationship to other node categories

  • Rigid body: receives position, sends contact/adhesion forces
  • Fluid environment: near-wall hydrodynamics may be handled by the fluid solver or by analytical corrections in this node
  • Therapeutic: targeted adhesion to specific tissue types is a drug delivery strategy

Phase assignment

  • Hydrodynamic wall corrections: Phase 1 (needed for any channel simulation)
  • Penalty-based contact: Phase 2
  • Adhesion models: Phase 2
  • Specific adhesion kinetics: Advanced

Existing implementations to examine

  • Brenner (1961) — classical wall correction factors for Stokes drag
  • Goldman, Cox, Brenner (1967) — translation and rotation near a plane wall
  • JKR/DMT adhesion models (Johnson et al. 1971, Derjaguin et al. 1975)
  • LIGGGHTS/DEM — open-source discrete element method with contact models

2.6 Magnetic Navigation in Confined Geometry

This is not a separate node but a distinct physics regime that sits at the intersection of rigid body dynamics, magnetic response, and near-wall hydrodynamics. It deserves its own treatment because it is the core operational regime of a helical robot navigating CSF channels — where most clinical navigation happens.

What it models

The phase diagram of helical robot behaviour in a rotating magnetic field near a wall or inside a channel. In free space, the robot synchronises with the rotating field and propels forward. As confinement increases or frequency rises, the robot transitions between distinct locomotion modes:

  • Synchronised precession: robot tumbles in a cone-shaped precession pattern, generating net propulsion along the field rotation axis. This is the normal operating mode.
  • Step-out: driving frequency exceeds the maximum synchronised frequency. The robot loses phase lock with the field, rotation becomes irregular, and propulsion drops abruptly. The step-out frequency is the single most important characterisation parameter for a helical microrobot.
  • Wobbling: intermediate regime near step-out, where the robot oscillates between synchronised and asynchronous states. The frequency-averaged velocity is reduced but non-zero. This regime has distinct dynamics from both synchronised and fully stepped-out states.
  • Wall-walking: in a confined channel, the robot rolls along the channel wall rather than swimming freely through the centre. This is a distinct locomotion mode that can be deliberately exploited for navigation along vessel walls — it provides a surface to push against and can be more efficient than free swimming in narrow channels.

The boundaries between these regimes depend on:

  • The ratio of magnetic torque to viscous drag torque (Mason number, Mn)
  • The driving frequency relative to the step-out frequency (frequency ratio)
  • The channel-to-robot size ratio (confinement ratio, Lambda)
  • The robot's proximity to the wall (gap distance)

Why it is needed

A robot navigating the aqueduct of Sylvius (channel diameter ~1.5 mm, robot diameter ~0.5 mm, confinement ratio Lambda ~ 3) is heavily confined. The confinement fundamentally changes the physics:

  • The step-out frequency is confinement-dependent — a robot that synchronises at 20 Hz in free fluid may step out at 12 Hz in a confined channel due to increased wall drag. Ignoring this means B1 (step-out detection) produces physically incorrect results for channel geometries.
  • Wall-walking mode does not exist in free space. It emerges only when the robot is near a boundary.
  • The propulsion velocity in a channel differs from the free-space velocity by a factor that depends on confinement ratio and gap distance.
  • The transition between locomotion modes involves bifurcations (sudden changes in behaviour) that are difficult to represent smoothly in a JAX-traceable simulation (requires jnp.where branching, not Python if).

Key dimensionless parameters

The literature search should focus on the following dimensionless groups:

  • Mason number Mn = (viscous drag torque) / (magnetic torque) = (8 * pi * mu * omega * a^3) / (m * B) — governs the synchronisation-to-step-out transition. Step-out occurs when Mn exceeds a critical value that depends on robot geometry.
  • Confinement ratio Lambda = D_channel / L_robot — governs wall effects. Lambda >> 1 means free swimming; Lambda ~ 1–5 means significant confinement; Lambda < 1 means the robot is larger than the channel.
  • Frequency ratio Omega/Omega_c — ratio of driving frequency to critical (step-out) frequency. Omega_c itself depends on Lambda.
  • Gap Reynolds number Re_gap = rho * U * h / mu — Reynolds number based on the gap between robot and wall. Even when the bulk Re is very low, the gap Re can be non-negligible in tight confinement.

Known challenges

  • Confinement-dependent step-out frequency: the relationship between Lambda and the critical frequency is not a simple analytical formula for realistic robot geometries. It requires either pre-computed lookup tables from detailed Stokes flow simulations or empirical correlations from the literature. The literature search should identify published correlations for helical swimmers in cylindrical channels.
  • Mode bifurcations: the transitions between synchronised, wobbling, wall-walking, and stepped-out states involve discontinuities in propulsion velocity. In a JAX-traceable simulation, these must be smoothed or handled via jnp.where branching. The literature should identify the mathematical structure of these bifurcations.
  • Coupled wall hydrodynamics: the wall correction to the resistance tensor is position-dependent and geometry-dependent. Near a curved wall (inside a cylindrical channel), the corrections differ from those near a flat wall.
  • Benchmark implications: B1 (step-out detection) and B5 (step-out recovery) must be specified with a confinement ratio to be physically meaningful. A step-out frequency measured in free space is not the same as the one in a channel.

Relationship to other nodes

This is an emergent behaviour of the interaction between ExternalMagneticFieldNode, MagneticResponseNode, RigidBodyNode, CSFFlowNode, and SurfaceContactNode. No separate node is needed, but:

  • The rigid body node must accept confinement-dependent resistance tensor corrections
  • The phase tracking node may need to track which locomotion mode the robot is in (extending its responsibilities beyond simple phase error tracking — the literature search should identify whether this requires a separate mode classification node or can be inferred from phase error dynamics)
  • The surface contact node must provide the wall-proximity corrections that enter the resistance tensor

Phase assignment

  • Phase 1 (essential — the confinement correction to step-out frequency is needed for B1 to be physically meaningful in a channel geometry)

Existing implementations to examine

  • Purcell, E.M. (1977). "Life at Low Reynolds Number." American Journal of Physics 45(1):3–11 — foundational
  • Lauga, E. & Powers, T.R. (2009). "The hydrodynamics of swimming microorganisms." Reports on Progress in Physics 72(9):096601 — review of locomotion mode transitions
  • Peyer, K.E., Zhang, L., Nelson, B.J. (2013). "Bio-inspired magnetic swimming microrobots for biomedical applications." Nanoscale 5(4):1259–1272 — confinement effects
  • Abbott, J.J. et al. (2009). "How Should Microrobots Swim?" International Journal of Robotics Research 28(11–12):1434–1447 — step-out analysis
  • Fischer, P. et al. — multiple publications on helical microswimmer dynamics, step-out, and confinement (Max Planck Institute for Medical Research / Physical Intelligence department)
  • Zöttl, A. & Stark, H. (2012). "Nonlinear dynamics of a microswimmer in Poiseuille flow." Physical Review Letters 108:218104 — confinement effects on microswimmer dynamics

2.7 Cross-cutting: Geometry Dependency of Spatial Nodes

What it is

Several node categories in §3 (fluid environment), §4 (tissue), and §5 (therapeutic) require a spatial domain — a geometry within which the physics is solved. This geometry ranges from a simple parametric description (a straight cylinder of given diameter and length, for bench-top validation) to a patient-specific mesh derived from MRI segmentation (for clinically realistic simulation). The relationship between MIME physics nodes and their geometry source is a first-class architectural concern, not an implementation detail.

Why it needs a defined interface

Without a defined GeometrySource interface, each spatial node would either hard-code a geometry type (making it impossible to swap parametric for mesh geometries) or defer geometry entirely to MICROBOTICA (making the node untestable without the full simulator). Neither is acceptable. The GeometrySource interface in mime/core/geometry.py defines the contract that spatial nodes depend on, decoupling the physics from the geometry loading mechanism.

Geometry types required by phase

Benchmark Geometry type Provider
B0 (experimental validation) Parametric cylinder matching the experimental paper's channel Author-specified
B1, B2 (physics benchmarks) Parametric cylinder or free-space (no walls) Author-specified
B4-T1 Parametric cylinder (D=2mm, L=50mm) MIME built-in
B4-T2 Neurobotika-derived ventricular mesh Neurobotika pipeline
B4-T3 Pathological anatomy mesh variant Neurobotika pipeline (pathological variants)

Key design requirements for GeometrySource

  • Must be usable without MICROBOTICA (enables standalone MIME testing and CI benchmarks)
  • Must support versioning so that a MimeAssetSchema remains linked to the exact geometry used during benchmarking
  • Must provide at minimum: domain bounds, boundary surface representation, and coordinate frame
  • Parametric subtypes (cylinder, sphere, torus) must be serialisable to/from USD as typed prims, and as JSON for standalone MIME testing without MICROBOTICA
  • Mesh subtypes must carry a provenance reference (Neurobotika version, segmentation parameters, source MRI metadata) as USD metadata attributes on the mesh prim

Phase assignment

  • Parametric GeometrySource (cylinder, sphere): Phase 0 (needed before any spatial node can be implemented)
  • Mesh GeometrySource (Neurobotika-derived): Phase 2 (needed for B4-T2 and realistic drug diffusion)
  • Pathological anatomy variants: Phase 3 / Advanced

3. Fluid Environment Nodes

These nodes model the physiological fluid medium in which the microrobot operates. The dominant fluid environments for medical microrobotics are cerebrospinal fluid (CSF), blood, and interstitial fluid.

3.1 Stokes / Creeping Flow

What it models

The steady or quasi-steady flow of a viscous fluid at very low Reynolds number (Re << 1), where inertial effects are negligible. This is the Stokes equation:

-nabla p + mu nabla^2 u = f, nabla · u = 0

At the microrobot scale (L ~ 10-500 um, U ~ 10-1000 um/s, nu ~ 0.7e-6 m^2/s for CSF), Re ranges from 10^-4 to 10^-1 — firmly in the creeping flow regime.

The Stokes equation is linear, time-reversible (Purcell's scallop theorem), and allows superposition of solutions. This simplifies the computation: the flow field around the robot can be decomposed into contributions from translation, rotation, and external forcing.

Why it is needed

  • The ambient flow field determines the drag force and torque on the robot
  • Stokes flow is the correct physical model for microscale dynamics in biological fluids
  • The linearity of Stokes flow enables efficient computation via Green's functions (Stokeslets) and boundary element methods
  • Understanding the flow field around the robot is essential for drug transport modelling

Inputs and outputs

  • Inputs: robot position and velocity (boundary condition), external forcing (gravity, pressure gradient), domain geometry (from mesh or analytical description)
  • Outputs: velocity field u(x), pressure field p(x), drag force and torque on the robot
  • Parameters: fluid viscosity mu, density rho, domain geometry, boundary conditions

Key physical parameters

  • Dynamic viscosity of CSF: mu ~ 0.7-1.0 mPa·s (close to water at 37°C)
  • Density of CSF: rho ~ 1005 kg/m^3
  • Characteristic velocity: 10 um/s to 1 mm/s (robot) + ambient CSF flow
  • Channel dimensions: lateral ventricles (several cm), aqueduct of Sylvius (1-2 mm diameter), subarachnoid space (variable)

Known challenges

  • Complex geometry: the CSF spaces have highly complex, patient-specific geometry. Solving Stokes flow in realistic ventricular anatomy requires mesh generation from MRI/CT data and finite element or boundary element methods.
  • Moving boundary: the robot moves through the fluid, changing the boundary conditions at each timestep. This requires either re-meshing (expensive), immersed boundary methods, or regularised Stokeslet approaches.
  • Efficiency: full 3D Stokes flow solutions are expensive. For large domains (full ventricular system), lattice Boltzmann or singularity methods may be more practical than FEM.
  • Flow-structure coupling: bidirectional coupling between robot motion and fluid flow requires iterative solution at each timestep (or coupling group iteration in MADDENING).

Differentiability status: [CONDITIONALLY DIFFERENTIABLE] — differentiability depends on the specific discretisation and parameter regime. Assess per implementation.

Relationship to other node categories

  • Rigid body: bidirectional — receives robot position/velocity as boundary conditions, sends drag force/torque
  • Diffusion: the flow field advects drug concentration (advection-diffusion equation)
  • Tissue: CSF flow interacts with compliant tissue boundaries (choroid plexus, ependyma)

Phase assignment

  • Phase 1 (essential — every simulation with fluid coupling needs this)

Existing implementations to examine

  • MADDENING's LBMPipeNode — lattice Boltzmann method already implemented in MADDENING (potential Mode 1 wrapping for CSF)
  • PyStokes — regularised Stokeslet method for microswimmers
  • Firedrake/FEniCS — finite element Stokes solvers
  • OpenFOAM — general CFD (overkill but useful for validation reference data)

3.2 Pulsatile Flow

What it models

The time-varying flow of CSF (or blood) driven by cardiac and respiratory cycles. CSF flow is not steady — it oscillates with each heartbeat (cardiac pulsation, ~1 Hz) and with respiration (~0.25 Hz). The flow pattern depends on the anatomical location:

  • Aqueduct of Sylvius: oscillatory, nearly sinusoidal, peak velocity ~5-10 cm/s
  • Lateral ventricles: slow circulation driven by choroid plexus secretion (~0.35 mL/min total production), with superimposed cardiac pulsation
  • Subarachnoid space: complex oscillatory flow patterns driven by brain pulsation

In blood vessels, pulsatile flow is driven by the cardiac cycle (Womersley flow). The Womersley number Wo = R * sqrt(omega/nu) characterises the ratio of pulsatile inertia to viscous effects.

Why it is needed

  • CSF flow is not steady — ignoring pulsatility means ignoring a dominant component of the flow that the robot must navigate through
  • Pulsatile flow creates time-varying forces on the robot that affect navigation accuracy
  • The phase of the cardiac cycle relative to the robot's actuation determines whether the flow helps or hinders navigation
  • Drug transport in CSF is significantly affected by pulsatile mixing (Taylor dispersion)
  • Benchmark B4 (closed-loop navigation) requires realistic flow conditions to be meaningful

Inputs and outputs

  • Inputs: flow rate waveform Q(t) (from physiological data or analytical model), anatomy geometry
  • Outputs: time-varying velocity field u(x,t), time-varying pressure gradient, pulsatile drag on robot
  • Parameters: mean flow rate, pulsation amplitude, cardiac frequency, respiratory frequency, Womersley number

Key physical parameters

  • Cardiac pulsation frequency: ~1 Hz (60-100 bpm)
  • Respiratory frequency: ~0.2-0.3 Hz (12-18 breaths/min)
  • Peak CSF velocity in aqueduct: ~5-10 cm/s
  • Mean CSF production rate: ~0.35 mL/min (total)
  • Womersley number in aqueduct: Wo ~ 5-15 (significant inertial effects during pulsation)

Known challenges

  • Multi-frequency drive: CSF pulsation has both cardiac and respiratory components. The superposition creates complex, quasi-periodic flow patterns.
  • Patient variability: pulsatile flow patterns vary significantly between patients (age, pathology, intracranial compliance). Parametric models must cover a range of physiological conditions.
  • Computational cost: time-dependent 3D flow in complex geometry is expensive. Reduced-order models (1D network models, modal decomposition) may be needed for real-time use.
  • Phase-locked actuation: the control system may want to synchronise actuation with the cardiac cycle to exploit favourable flow phases. This requires the flow model to be driven by a realistic cardiac waveform.

Differentiability status: [CONDITIONALLY DIFFERENTIABLE] — differentiability depends on the specific discretisation and parameter regime. Assess per implementation.

Relationship to other node categories

  • Stokes flow: pulsatile flow extends the Stokes flow node with time-dependent forcing
  • Rigid body: time-varying drag affects robot trajectory
  • Diffusion: pulsatile flow drives Taylor dispersion of drug concentration
  • Sensing: flow velocity is a quantity that Doppler ultrasound can measure

Phase assignment

  • Phase 1 (important for realistic CSF simulation — can start with analytical Womersley profiles in simple geometries before moving to full 3D)

Existing implementations to examine

  • Womersley analytical solutions (1955) — pulsatile flow in circular tubes
  • Phase-contrast MRI data — provides in-vivo CSF flow velocity measurements for validation
  • OpenFOAM pimpleFoam — transient incompressible flow solver

3.3 Non-Newtonian Rheology

What it models

The deviation of biological fluids from Newtonian viscosity behaviour. While CSF is approximately Newtonian (viscosity close to water), two important biological fluids are not:

Blood: exhibits shear-thinning behaviour (viscosity decreases with shear rate) due to red blood cell aggregation at low shear rates and deformation/alignment at high shear rates. Models: Carreau-Yasuda, Cross, Casson, power-law.

Mucus: highly viscoelastic (both viscous and elastic response) due to its polymer gel structure. Relevant if the microrobot must traverse mucous barriers (e.g., respiratory mucosa, cervical mucus). Models: Maxwell, Oldroyd-B, Giesekus, FENE-P.

Why it is needed

  • If the microrobot operates in blood (e.g., vascular microrobotics), Newtonian models significantly underpredict the viscosity at the low shear rates characteristic of microrobot-scale flows
  • Viscoelastic fluids can generate elastic "swimming" forces that don't exist in Newtonian fluids — some microrobot designs exploit this
  • Drug transport through mucus barriers is governed by viscoelastic rheology
  • For CSF-only simulations, non-Newtonian effects are minimal — but the framework should support them for broader applicability

Inputs and outputs

  • Inputs: velocity field (or velocity gradient tensor) from flow solver
  • Outputs: stress tensor (to flow solver as modified viscosity or extra stress), effective viscosity field
  • Parameters: zero-shear viscosity, infinite-shear viscosity, relaxation time, power-law exponent, Carreau parameters

Key physical parameters

  • Zero-shear viscosity of blood: ~100 mPa·s (at hematocrit 45%)
  • Infinite-shear viscosity of blood: ~3-4 mPa·s
  • Characteristic shear rate for transition: ~1-10 s^-1
  • Relaxation time for mucus: ~1-100 s (highly variable)
  • Storage modulus G' and loss modulus G'' of mucus (Pa)

Known challenges

  • Constitutive equation selection: many competing models exist for blood and mucus. The choice depends on the shear rate regime and the specific application.
  • High Weissenberg number problem: viscoelastic flow simulations become numerically unstable at high Weissenberg numbers (Wi = relaxation time × shear rate). This is a classical challenge in computational rheology.
  • Multi-scale: red blood cell effects occur at the cellular scale (~8 um), comparable to the microrobot size. Continuum rheology models may not be appropriate when the robot is the same size as the blood cells.
  • Temperature dependence: viscosity is temperature-dependent (body temperature is 37°C, not room temperature)

Differentiability status: [CONDITIONALLY DIFFERENTIABLE] — differentiability depends on the specific discretisation and parameter regime. Assess per implementation.

Relationship to other node categories

  • Stokes/pulsatile flow: modifies the viscosity used in the flow equation
  • Rigid body: non-Newtonian drag differs from Newtonian drag (shear-thinning reduces drag at high velocities)
  • Diffusion: viscoelastic fluids hinder diffusion (anomalous diffusion in polymer gels)

Phase assignment

  • Shear-thinning blood models (Carreau-Yasuda): Phase 2 (needed for vascular applications)
  • Viscoelastic mucus models: Advanced (needed for specific barrier-crossing scenarios)

Existing implementations to examine

  • RheoTool (OpenFOAM) — viscoelastic flow solver
  • FEniCS with dolfin-adjoint — supports non-Newtonian constitutive models
  • Oldroyd-B implementation in various open-source CFD codes

3.4 Diffusion / Species Transport

What it models

The spreading of a dissolved substance (drug, tracer, nutrient) in the fluid medium by molecular diffusion and advection. Governed by the advection-diffusion equation:

dc/dt + u · nabla c = D nabla^2 c + S

where c is concentration, u is fluid velocity, D is diffusion coefficient, and S is a source/sink term (e.g., drug release from the robot, uptake by tissue).

Why it is needed

  • Drug delivery is the primary clinical application. After the robot releases its payload, the drug must reach the target tissue at a therapeutic concentration. The concentration field evolves by diffusion and advection.
  • Answers: "What concentration does the drug reach at the target? How long after release? What fraction is lost to the bulk CSF?"
  • Needed for benchmark B3 (drug release kinetics)
  • Tracer transport enables validation against experimental dye studies

Inputs and outputs

  • Inputs: velocity field (from flow solver) for advection; source terms (from drug release node); boundary conditions (tissue uptake, clearance)
  • Outputs: concentration field c(x,t); concentration at target location; total drug mass in domain (mass conservation check)
  • Parameters: diffusion coefficient D, Peclet number Pe = UL/D, initial concentration, boundary condition type

Key physical parameters

  • Diffusion coefficient of typical drugs in CSF: D ~ 10^-10 to 10^-9 m^2/s (small molecules: ~10^-9; proteins: ~10^-10; nanoparticles: ~10^-11)
  • Peclet number for microrobot-released drugs: Pe ~ 1-1000 (regime where both advection and diffusion matter)
  • CSF clearance rate: ~0.35 mL/min (bulk flow removal)
  • Tissue uptake rate (varies by drug and target)

Known challenges

  • Numerical diffusion: advection-dominated transport (high Pe) suffers from numerical diffusion in first-order schemes, smearing sharp concentration fronts. Higher-order schemes (MUSCL, WENO) or lattice Boltzmann methods mitigate this.
  • Point source: drug released from a microrobot is essentially a moving point source. Resolving the sharp concentration gradient near the source requires either very fine mesh resolution or analytical inner solutions matched to numerical outer solutions.
  • Long time scales: drug diffusion operates on much longer timescales than robot motion. A 100 um drug cloud takes ~10 s to diffuse 100 um in CSF. Multi-rate timestepping (MADDENING feature) is essential.
  • Complex geometry: diffusion in realistic CSF spaces with narrow passages (aqueduct) and large reservoirs (ventricles) requires 3D mesh resolution

Differentiability status: [DIFFERENTIABLE] — gradients are reliable for policy training, sensitivity analysis, and uncertainty propagation within the validated physical regime.

Relationship to other node categories

  • Flow: velocity field drives advection
  • Therapeutic (drug release): source term comes from the release node
  • Tissue: boundary conditions for tissue uptake and barrier transport

Phase assignment

  • Phase 2 (needed for drug delivery simulation; needed for B3)

Existing implementations to examine

  • MADDENING's HeatNode — solves the same equation (heat diffusion ~ mass diffusion by analogy); potential Mode 1 wrapping with appropriate parameters
  • LBM passive scalar transport (already in MADDENING's LBMPipeNode)
  • FEniCS advection-diffusion solver

4. Biological Tissue Nodes

These nodes model the interaction between the simulation and the surrounding biological tissue. They represent the boundary between the fluid domain and the living tissue.

4.1 Tissue Deformation / Contact

What it models

The mechanical response of biological tissue (vessel walls, brain parenchyma, ependymal lining) to forces applied by the microrobot or by fluid pressure. Biological tissues are soft, viscoelastic, and often anisotropic.

Two main scenarios:

Vessel wall compliance: the walls of CSF channels and blood vessels are not rigid. They deform under fluid pressure (pulsatile distension) and under contact forces from the microrobot. The wall compliance affects the flow field (fluid-structure interaction) and the contact mechanics.

Soft tissue indentation: when the microrobot contacts or pushes against tissue (e.g., brain surface for targeted delivery), the tissue deforms. The indentation depth and contact area depend on the robot's force and the tissue's mechanical properties.

Why it is needed

  • Rigid wall assumptions are inaccurate for biological channels — vessel walls distend with each cardiac pulse, changing the channel diameter and flow velocity
  • Contact between the robot and tissue may cause damage. Predicting tissue stress is essential for safety — the hazard_hints for these nodes feed directly into ISO 14971 risk assessment
  • Compliant walls change the near-wall hydrodynamics that affect robot navigation
  • Drug delivery to tissue depends on the contact area and pressure at the delivery site

Inputs and outputs

  • Inputs: fluid pressure on wall, contact force from robot, cardiac pressure waveform
  • Outputs: wall displacement, wall velocity (for FSI boundary condition), contact area, tissue stress
  • Parameters: Young's modulus, Poisson's ratio, wall thickness, viscoelastic relaxation time, ultimate tensile stress (damage threshold)

Key physical parameters

  • Brain tissue Young's modulus: ~1-10 kPa (extremely soft)
  • Blood vessel wall Young's modulus: ~100 kPa to 1 MPa (varies by vessel type and age)
  • Ependymal lining thickness: ~10-50 um
  • CSF channel wall compliance (distensibility): varies with intracranial pressure
  • Tissue damage threshold stress: varies by tissue type

Known challenges

  • Constitutive modelling: brain tissue is nonlinear, viscoelastic, and anisotropic. Simple linear elastic models are inadequate for large deformations. Hyperelastic models (Ogden, Mooney-Rivlin) or poroelastic models (brain as a porous medium with interstitial fluid) may be needed.
  • Patient variability: tissue properties vary significantly with age, pathology, and anatomical location
  • Fluid-structure interaction: compliant walls require coupling between the flow solver and the structural solver. This is computationally expensive and may require coupling group iteration in MADDENING.
  • Damage modelling: predicting when tissue damage occurs requires stress/strain failure criteria. This is safety-critical — a hazard_hint for these nodes.

Differentiability status: [CONDITIONALLY DIFFERENTIABLE] — differentiability depends on the specific discretisation and parameter regime. Assess per implementation.

Relationship to other node categories

  • Flow: bidirectional FSI — wall motion changes the flow domain
  • Rigid body / surface contact: wall deformation affects contact mechanics
  • Therapeutic: drug delivery to tissue depends on contact conditions

Phase assignment

  • Rigid wall model (baseline, no deformation): implicit in Phase 1 (boundary conditions for flow)
  • Compliant wall model: Phase 2
  • Tissue damage modelling: Advanced (safety-critical, requires careful validation)

Existing implementations to examine

  • FEBio — finite element software for biomechanics (widely used for brain tissue)
  • NiftySim — GPU-accelerated soft tissue simulation
  • Brain tissue constitutive models: Ogden, Mooney-Rivlin, Bilston (2001)

4.2 Biological Barrier Transport

What it models

The transport of drugs (or the microrobot itself) across biological barriers — membranes, cell layers, and tissue structures that separate compartments. The most important barriers for microrobotics drug delivery:

Blood-brain barrier (BBB): the endothelial cell layer lining cerebral blood vessels, with tight junctions that restrict paracellular transport. Most drugs cannot cross the BBB. The microrobot may carry drugs across the BBB by disrupting the tight junctions (using focused ultrasound), by transcytosis (nanoparticle-mediated), or by direct intrathecal injection (bypassing the BBB entirely via CSF).

Ependymal barrier: the cell layer lining the ventricles. Less restrictive than the BBB but still a barrier to drug penetration into brain parenchyma.

Cell membrane: if the microrobot is designed to enter cells (nano-scale), membrane crossing involves endocytosis or membrane disruption.

Why it is needed

  • The therapeutic efficacy of microrobot-delivered drugs depends on whether the drug reaches its target in the brain parenchyma
  • BBB disruption by focused ultrasound is a major research area — modelling the disruption and subsequent drug transport is essential for treatment planning
  • Even in CSF delivery (bypassing the BBB), the drug must cross the ependymal barrier to reach deep brain tissue

Inputs and outputs

  • Inputs: drug concentration at barrier surface (from diffusion node), barrier permeability (may be time-varying if disrupted), ultrasound parameters (for BBB disruption)
  • Outputs: drug flux across barrier, drug concentration on the other side, barrier integrity state
  • Parameters: permeability coefficient P (cm/s), barrier thickness, tight junction opening rate (for BBB disruption), recovery time constant

Key physical parameters

  • BBB permeability to small molecules: P ~ 10^-7 to 10^-5 cm/s (varies dramatically by drug)
  • Ependymal permeability: higher than BBB (~10^-5 to 10^-3 cm/s for small molecules)
  • BBB disruption duration after focused ultrasound: ~4-6 hours
  • BBB recovery time constant: ~hours to days

Known challenges

  • Multi-scale: BBB transport involves molecular-scale phenomena (tight junction opening, receptor-mediated transcytosis) that must be represented as effective permeability coefficients in the continuum model
  • Patient and region variability: BBB permeability varies by brain region and with pathology (stroke, tumor, inflammation)
  • Safety: BBB disruption creates a window for both therapeutic drug entry and harmful substance entry. This is a critical hazard_hint.
  • Limited in-vivo data: barrier permeability measurements are difficult and uncertain

Differentiability status: [CONDITIONALLY DIFFERENTIABLE] — differentiability depends on the specific discretisation and parameter regime. Assess per implementation.

Relationship to other node categories

  • Diffusion: receives concentration at barrier surface, sends flux as boundary condition
  • External acoustic: focused ultrasound drives BBB disruption
  • Therapeutic (pharmacokinetics): barrier transport determines drug bioavailability

Phase assignment

  • Advanced (requires detailed pharmacokinetic modelling; important for treatment planning but not for basic robot navigation simulation)

Existing implementations to examine

  • Kety-Schmidt model — classical BBB permeability model
  • PBPK (Physiologically Based Pharmacokinetic) models — compartmental barrier transport
  • Focused ultrasound BBB disruption models (McDannold, Hynynen groups)

5. Therapeutic Payload Nodes

These nodes model the drug delivery mechanics — how the drug leaves the microrobot and reaches the target tissue.

5.1 Drug Release Kinetics

What it models

The rate at which drug is released from the microrobot's payload compartment. Release can be:

  • Passive diffusion: drug slowly leaks from a porous coating or matrix
  • pH-triggered: the coating dissolves at a specific pH (e.g., acidic tumor microenvironment)
  • Magnetically triggered: an alternating magnetic field heats magnetic nanoparticles, melting a thermosensitive coating
  • Ultrasonically triggered: acoustic cavitation from the microbubble coating disrupts the drug carrier
  • Enzymatically triggered: specific enzymes in the target tissue degrade the coating

The release profile follows a kinetic model (first-order, Higuchi, Korsmeyer-Peppas, or mechanistic diffusion).

Why it is needed

  • The timing and rate of drug release directly determine therapeutic efficacy
  • If the drug is released too early (before reaching the target), it is wasted in the bulk CSF
  • If the drug is released too slowly, the therapeutic window may not be reached
  • Answers: "How much drug has been released at time t? What is the release rate? When is the payload exhausted?"
  • Needed for benchmark B3

Inputs and outputs

  • Inputs: trigger signal (pH, temperature, ultrasound pressure, enzyme concentration), elapsed time since trigger
  • Outputs: drug release rate (mass/time), remaining payload fraction, drug source term for diffusion node
  • Parameters: total payload mass, release half-life, trigger threshold, coating thickness, pore size, matrix diffusivity

Key physical parameters

  • Total drug payload (pg to ng for microrobots)
  • Release half-life (minutes to hours, depending on mechanism)
  • Trigger threshold (pH 6.5 for tumor targeting, temperature 42°C for magnetic heating)
  • Coating degradation rate

Known challenges

  • Multi-mechanism release: real drug carriers may involve multiple simultaneous release mechanisms (e.g., passive diffusion + triggered burst release)
  • Carrier degradation: the carrier structure changes as drug is released (pore opening, shell thinning), changing the release kinetics over time
  • Environmental sensitivity: release kinetics may depend on local pH, temperature, ionic strength — quantities that vary in space and time in the body
  • Stochastic release: at the single-robot level, release events (coating rupture, pore opening) may be stochastic. Ensemble statistics are needed.

Differentiability status: [DIFFERENTIABLE] — gradients are reliable for policy training, sensitivity analysis, and uncertainty propagation within the validated physical regime.

Relationship to other node categories

  • Diffusion: the released drug enters the fluid as a source term
  • External acoustic/magnetic: provides trigger signal for stimuli-responsive release
  • Rigid body: position determines where the drug is released
  • Biological barriers: released drug must cross barriers to reach the target

Phase assignment

  • First-order and diffusion-based release: Phase 2
  • Stimuli-responsive (triggered) release: Phase 2
  • Mechanistic multi-mechanism models: Advanced

Existing implementations to examine

  • Higuchi model (1963) — classical matrix diffusion release
  • Korsmeyer-Peppas model (1983) — power-law release kinetics
  • COMSOL Drug Delivery module — commercial reference

5.2 Pharmacokinetic Transport

What it models

The fate of the released drug in the body — uptake by target tissue, clearance by bulk CSF flow, metabolism, and redistribution. This is the compartmental pharmacokinetics layer:

  • Target tissue uptake: rate of drug absorption by the intended target (tumor, neural tissue)
  • CSF clearance: drug removal by CSF flow into the venous system (arachnoid granulations)
  • Metabolism: drug degradation by enzymes in CSF or tissue
  • Redistribution: drug spread to non-target tissues

Why it is needed

  • Drug release is only the first step. The therapeutic outcome depends on the drug reaching the target at a sufficient concentration for a sufficient duration.
  • Answers: "What is the drug concentration at the target tissue over time? Is the therapeutic window achieved? What fraction of the dose reaches the target vs. being cleared?"
  • Essential for treatment planning and dose optimisation

Inputs and outputs

  • Inputs: drug concentration field (from diffusion node), tissue properties, clearance rates
  • Outputs: tissue drug concentration, total target exposure (AUC), clearance rate, therapeutic window status
  • Parameters: tissue partition coefficient, clearance rate constant, metabolic half-life, target binding affinity

Key physical parameters

  • CSF turnover time: ~6-8 hours (total volume ~150 mL, production ~0.35 mL/min)
  • Drug half-life in CSF (highly drug-dependent)
  • Tissue-to-CSF partition coefficient
  • Target receptor binding affinity Kd

Known challenges

  • Compartmental vs. distributed: simple compartmental PK models (1-compartment, 2-compartment) may not capture the spatial heterogeneity of drug distribution in complex CSF anatomy. Spatially-distributed PK (coupling with the diffusion node) is more accurate but more expensive.
  • Patient variability: PK parameters vary significantly between patients
  • Nonlinear binding: at high local concentrations (near the robot), receptor saturation and nonlinear binding kinetics may apply
  • Blood-CSF exchange: drugs in CSF can be absorbed into the bloodstream and vice versa

Differentiability status: [CONDITIONALLY DIFFERENTIABLE] — differentiability depends on the specific discretisation and parameter regime. Assess per implementation.

Relationship to other node categories

  • Diffusion: provides the concentration field
  • Drug release: provides the source term
  • Biological barriers: determines transcellular drug transport rates
  • Tissue: tissue properties affect uptake

Phase assignment

  • Compartmental PK: Phase 2
  • Spatially-distributed PK (coupled to diffusion): Advanced

Existing implementations to examine

  • SimBiology (MATLAB) — compartmental PK modelling
  • PBPK models for intrathecal drug delivery (Linninger group)
  • CSF pharmacokinetics models (Hladky & Barrand 2014)

6. Sensing / Imaging Physics Nodes

These nodes model the physics of the localisation system — how the microrobot's position and state are measured. In the MIME architecture, sensing nodes produce raw measurement data that is consumed by the UncertaintyModel to generate realistic observation noise.

6.1 MRI Signal Formation

What it models

The formation of magnetic resonance imaging signals that are used to localise and track the microrobot. MRI is the most promising in-vivo tracking modality for microrobots in CSF because it provides 3D imaging without ionising radiation.

Key phenomena:

  • Susceptibility artefacts: the microrobot's magnetic material creates local field inhomogeneities that distort the MRI signal, producing characteristic artefacts. These artefacts can be larger than the robot itself and are used for tracking (the artefact is the "signature" of the robot).
  • Contrast mechanisms: T1, T2, T2* contrast — the robot's presence changes the local relaxation rates
  • k-space sampling: MRI acquires data in spatial frequency domain (k-space). The acquisition trajectory (Cartesian, radial, spiral) determines temporal resolution, spatial resolution, and artefact characteristics.
  • Motion artefacts: robot motion during MRI acquisition causes ghosting and blurring

Why it is needed

  • MRI is the clinical imaging modality for CSF spaces. Real-time MRI tracking of microrobots is an active research area.
  • Understanding MRI artefacts from the microrobot is essential for: (a) estimating position from the artefact, (b) predicting when tracking will fail, (c) designing robots with optimal MRI signatures
  • The sensing physics determines the ultimate achievable closed-loop control performance (you can't control what you can't see)
  • Populates SensingMeta in the asset schema

Inputs and outputs

  • Inputs: robot position, orientation, magnetisation state (from robot body nodes); MRI scanner parameters (field strength, sequence, resolution); surrounding tissue properties (T1, T2, proton density)
  • Outputs: simulated MRI image or k-space data, estimated robot position from artefact, position uncertainty (feeds into UncertaintyModel), tracking confidence
  • Parameters: B0 field strength (T), voxel size (mm), TR/TE timing (ms), susceptibility of robot material (ppm), sequence type

Key physical parameters

  • Scanner field strength: 1.5T or 3T (clinical); 7T+ (research)
  • Voxel size: 1-3 mm isotropic (clinical MRI) — the robot (< 1 mm) is sub-voxel
  • Temporal resolution: 10-500 ms per frame (depending on sequence)
  • Susceptibility of NdFeB: ~1000 ppm — creates very large artefacts
  • Susceptibility of iron oxide nanoparticles: ~10-100 ppm — smaller, more manageable artefacts

Known challenges

  • Artefact modelling: the susceptibility artefact depends on the robot's shape, magnetisation, and orientation relative to B0. Computing the artefact requires solving the magnetostatic field perturbation and then simulating the MRI signal formation. This is computationally intensive.
  • Sub-voxel localisation: the robot is smaller than the MRI voxel. Position estimation from the artefact requires model-based inversion (fitting an artefact model to the image data).
  • Temporal resolution vs. spatial resolution tradeoff: fast imaging (for real-time tracking) sacrifices spatial resolution and SNR
  • Interference with actuation field: the external magnetic field used for actuation interferes with the MRI field. MRI-compatible actuation hardware is required, and the simulation must model this interference.
  • Partial volume effects: at the interface between CSF and tissue, voxels contain both substances, affecting the signal

Differentiability status: [CONDITIONALLY DIFFERENTIABLE] — differentiability depends on the specific discretisation and parameter regime. Assess per implementation.

Relationship to other node categories

  • Robot body (magnetic response): magnetisation state affects susceptibility artefact
  • External apparatus: actuation field interferes with imaging field
  • Uncertainty model: MRI physics determines position noise, dropout probability, and tracking confidence

Phase assignment

  • Phase 2 (important for realistic closed-loop simulation)

Existing implementations to examine

  • JEMRIS — open-source MRI simulator (Bloch equation based)
  • MRiLab — MATLAB MRI simulation
  • Susceptibility artefact models from Martel group, Keenan/Bhatt publications

6.2 Ultrasound Imaging

What it models

The formation of ultrasound images used to track the microrobot. Ultrasound provides real-time imaging with high temporal resolution (> 1000 fps for ultrafast imaging) but limited spatial resolution and penetration in some anatomical contexts.

Key modes:

  • B-mode imaging: standard brightness-mode imaging showing tissue echogenicity. The microrobot appears as a bright spot (if it scatters strongly) or a shadow (if it absorbs).
  • Doppler imaging: detects motion by measuring frequency shifts. Can track the microrobot's velocity directly.
  • Pulse-echo localisation: measuring the time-of-flight of reflected pulses to determine distance to the robot.
  • Plane-wave ultrafast imaging: transmit unfocused plane waves and reconstruct images computationally. Enables very high frame rates (> 10 kHz) for real-time tracking.

Why it is needed

  • Ultrasound provides much higher temporal resolution than MRI (ms vs. hundreds of ms), enabling faster control loops
  • Ultrasound is portable and inexpensive compared to MRI
  • For acoustically-actuated robots, the same transducer array may be used for both actuation and imaging
  • Doppler ultrasound can measure both robot velocity and ambient flow velocity (useful for flow-aware control)

Inputs and outputs

  • Inputs: robot position, velocity, size, acoustic impedance; medium properties; transducer parameters
  • Outputs: simulated B-mode image, Doppler velocity estimate, position estimate, tracking confidence
  • Parameters: transducer frequency, imaging depth, frame rate, transmit sequence, medium speed of sound

Key physical parameters

  • Imaging frequency: 5-50 MHz (higher frequency = better resolution, lower penetration)
  • Spatial resolution: lambda/2 ~ 30-300 um (frequency-dependent)
  • Frame rate: 30-10000 fps
  • Penetration depth: 1-10 cm (frequency-dependent)
  • Speed of sound in tissue: ~1540 m/s

Known challenges

  • Acoustic access: the skull blocks most ultrasound. Transcranial ultrasound imaging for CSF spaces requires specialised low-frequency transducers and the bone creates severe aberration.
  • Speckle noise: ultrasound images are inherently speckled (coherent interference), making small object detection difficult
  • Shadowing and reverberation: the robot may shadow structures behind it or create reverberations
  • Registration with anatomy: ultrasound lacks the 3D anatomical context of MRI. Co-registration with pre-operative MRI may be needed.

Differentiability status: [CONDITIONALLY DIFFERENTIABLE] — differentiability depends on the specific discretisation and parameter regime. Assess per implementation.

Relationship to other node categories

  • External acoustic: same hardware may serve both actuation and imaging (time-division multiplexing)
  • Robot body: acoustic impedance mismatch determines echo strength
  • Uncertainty model: US physics determines position noise and tracking characteristics

Phase assignment

  • Phase 2 (important for high-temporal-resolution tracking)

Existing implementations to examine

  • k-Wave — also supports pulse-echo ultrasound simulation
  • Field II — ultrasound transducer simulation program
  • MUST (Matlab UltraSound Toolbox)

6.3 Fluorescence / Optical Imaging

What it models

Fluorescence imaging of microrobots labelled with fluorescent markers. The robot (or its payload) emits light when excited by a specific wavelength. A camera or fiber-optic system detects the emitted light.

Why it is needed

  • Fluorescence is the standard tracking method for in-vitro experiments (microfluidic channels, ex-vivo tissue)
  • Provides excellent spatial resolution (diffraction-limited, ~250 nm) and temporal resolution (video rate or faster)
  • Essential for validating simulation predictions against bench-top experiments
  • Not applicable for in-vivo use at depth (penetration < 1 mm in tissue)

Inputs and outputs

  • Inputs: robot position, fluorophore concentration, excitation light intensity, tissue optical properties
  • Outputs: fluorescence intensity at detector, apparent position (centroid), SNR
  • Parameters: excitation/emission wavelengths, quantum yield, photobleaching rate, camera exposure, NA of objective

Key physical parameters

  • Excitation/emission wavelength pair (nm)
  • Quantum yield (dimensionless)
  • Photobleaching rate (depends on fluorophore and illumination)
  • Tissue autofluorescence background

Known challenges

  • Depth limitation: fluorescence cannot image through more than ~1 mm of tissue. Limited to superficial or ex-vivo applications.
  • Photobleaching: fluorophore intensity decreases with cumulative exposure
  • Scattering: tissue scatters fluorescence, blurring the image of deeper structures
  • Autofluorescence: tissue has intrinsic fluorescence that creates background noise

Differentiability status: [CONDITIONALLY DIFFERENTIABLE] — differentiability depends on the specific discretisation and parameter regime. Assess per implementation.

Relationship to other node categories

  • Robot body: fluorophore is attached to the robot
  • Uncertainty model: fluorescence physics determines in-vitro tracking performance

Phase assignment

  • Advanced (useful for in-vitro validation scenarios)

6.4 Electromagnetic Localisation

What it models

Non-imaging electromagnetic methods for detecting the microrobot's position, including:

  • Inductive sensing: the robot's magnetic material is detected by pickup coils (similar to metal detector). The signal strength and phase indicate position.
  • Impedance-based sensing: the robot's presence changes the local electrical impedance. Measured via implanted electrodes or external coils.
  • Magnetic particle imaging (MPI): a specialised technique that directly images the spatial distribution of magnetic nanoparticles by exploiting their nonlinear magnetisation response.

Why it is needed

  • MPI provides direct, quantitative imaging of magnetic material with no tissue background signal — potentially the ideal modality for magnetic microrobot tracking
  • Inductive sensing is simpler and cheaper than MRI, potentially suitable for point-of-care tracking
  • These modalities complement MRI and ultrasound with different tradeoffs

Inputs and outputs

  • Inputs: robot position, magnetic moment, gradient field configuration (for MPI)
  • Outputs: sensor signal, estimated position, spatial resolution estimate
  • Parameters: coil geometry, drive field amplitude/frequency, receiver sensitivity

Key physical parameters

  • MPI drive field frequency: typically 25 kHz
  • MPI spatial resolution: ~1 mm (current state of the art)
  • Inductive sensing range: limited to a few cm
  • Sensitivity: depends on magnetic moment and distance

Known challenges

  • MPI availability: MPI scanners are rare and expensive; limited clinical adoption so far
  • Spatial resolution: current MPI resolution (~1 mm) is adequate for tracking but not for detailed imaging
  • Multi-robot disambiguation: distinguishing multiple robots requires different magnetic signatures

Differentiability status: [CONDITIONALLY DIFFERENTIABLE] — differentiability depends on the specific discretisation and parameter regime. Assess per implementation.

Relationship to other node categories

  • Robot body (magnetic response): magnetic moment determines signal strength
  • External apparatus: the MPI drive field interacts with the actuation field
  • Uncertainty model: determines tracking accuracy for electromagnetic localisation

Phase assignment

  • Advanced (specialised modalities, relevant when MPI becomes more accessible)

Existing implementations to examine

  • OpenMPIData — open-source MPI reconstruction
  • MPIReco.jl — Julia MPI reconstruction framework

7. Cross-Cutting Concerns

7.1 Brownian Motion / Thermal Fluctuations

Not a separate node category, but a force contribution that must be included in the rigid body dynamics node. At the microscale, thermal fluctuations generate random forces and torques:

  • F_Brownian ~ sqrt(2 * k_B * T * gamma_T / dt) * N(0,1) (translational)
  • T_Brownian ~ sqrt(2 * k_B * T * gamma_R / dt) * N(0,1) (rotational)

where gamma_T and gamma_R are translational and rotational friction coefficients. These scale inversely with robot size — significant for robots < 10 um, negligible for robots > 100 um.

Brownian motion requires explicit RNG key management in JAX. It is implemented as a stochastic forcing term in the rigid body node, not as a separate node.

7.2 Gravity and Buoyancy

Simple but important: gravitational settling and buoyancy depend on the density difference between the robot and the fluid. For robots near neutral buoyancy (density ratio ~1.0), settling is slow but non-negligible over long timescales.

Implemented as a constant force term in the rigid body node, not a separate node.

7.3 Multi-Timescale Coupling

Different physics operate on vastly different timescales:

Physics Timescale
Acoustic bubble oscillation ~1 us (MHz)
Magnetic field rotation ~10 ms (10-100 Hz)
Robot translation ~100 ms
CSF pulsation ~1 s (1 Hz cardiac)
Drug diffusion ~10-100 s
Drug release ~minutes to hours
Pharmacokinetics ~hours

MADDENING's multi-rate timestepping handles this: each node declares its own delta_t, and the graph manager derives a base timestep (GCD) with conditional application. Coupling groups handle bidirectional coupling between nodes at different rates.

7.4 Geometry / Mesh Management

Patient-specific anatomy (ventricle geometry, vessel diameters, brain surface) is an input to multiple nodes (flow, contact, tissue). Geometry handling is not a physics node but a shared concern:

  • Robot geometry: parametric (sphere, ellipsoid, helix with pitch/radius/length) or mesh (STL/OBJ)
  • Anatomy geometry: surface mesh from MRI segmentation (Neurobotika atlas or patient-specific)
  • Both are passed via geometry_source parameter in SimulationNode.__init__

Geometry management is a MICROBOTICA responsibility (USD scene format, mesh library). MIME nodes consume geometry as a constructor parameter, not as runtime state.

7.5 Out of Scope: Electroosmotic and Electrophoretic Effects

Electroosmotic flow (fluid motion driven by an applied electric field in a charged channel) and electrophoresis (motion of charged particles in an electric field) are significant in microfluidic lab-on-chip devices. They are explicitly out of scope for MIME Phase 1–3 for the following reasons:

  1. CSF channels are not electroosmotically driven — the dominant flows are pressure-driven (choroid plexus secretion) and pulsatile (cardiac and respiratory cycles).
  2. External electric fields strong enough to drive electroosmosis are not used in the target clinical application (magnetic and acoustic actuation dominate).
  3. The microrobot designs of interest use magnetic or acoustic actuation, not electrostatic.

This is a deliberate scope exclusion, not an oversight. If future MIME node designs involve electrically-actuated robots (e.g., self-electrophoretic Janus particles) or electrophoretic drug transport, this note should be revisited and the corresponding node categories added to the taxonomy.


8. Phase Assignment Summary

Phase 0 — Foundation (blocking for Phase 1)

Category Node Type Rationale
Cross-cutting Parametric GeometrySource interface defined Blocking for all spatial node implementation
Cross-cutting B0 experimental dataset identified and cited Must be confirmed before Phase 1 implementation begins

Phase 1 — Essential (basic magnetic helical robot in CSF)

Category Node Type Rationale
External apparatus Rotating magnet / Helmholtz coil field model Primary actuation hardware
Robot body Rigid body dynamics (6-DOF, overdamped) Every simulation needs this
Robot body Permanent magnet response Converts field to torque/force
Robot body Phase tracking (observer) Needed for B1, B5
Environment Stokes / creeping flow Drag computation
Environment Pulsatile flow (analytical Womersley) CSF is pulsatile
Robot body Hydrodynamic wall corrections Needed for channel navigation
Robot body Confinement-dependent step-out model (wall corrections to step-out frequency) Needed for B1 to be meaningful in channel geometry; core clinical navigation regime

Phase 2 — Important (drug delivery, realistic environment)

Category Node Type Rationale
Robot body Soft-magnetic response Broader material support
Robot body Flexible body mechanics Flagellar/compliant robots
Robot body Acoustic bubble dynamics Acoustic actuation modality
Robot body Surface contact / adhesion Wall interaction
External apparatus Gradient field model Gradient-based steering
External apparatus Focused ultrasound Acoustic actuation
Environment Non-Newtonian rheology Blood applications
Environment Diffusion / species transport Drug transport (B3)
Tissue Compliant wall model FSI
Therapeutic Drug release kinetics B3
Therapeutic Compartmental pharmacokinetics Drug fate
Sensing MRI signal formation Realistic imaging
Sensing Ultrasound imaging High-speed tracking

Advanced (specific modalities, edge cases)

Category Node Type Rationale
External apparatus Optical sources Limited in-vivo applicability
Tissue Biological barrier transport BBB modelling
Tissue Tissue damage modelling Safety-critical
Therapeutic Spatially-distributed PK Full PK modelling
Sensing Fluorescence imaging In-vitro only
Sensing Electromagnetic localisation / MPI Specialised modalities
Environment Viscoelastic (mucus) rheology Specific barriers