🔢 mathematics

Convergence of the Immersed Interface Method in Linear Elasticity

This paper establishes and numerically validates that the $L^2$ -norm error between the exact solution of a linear elasticity problem with a continuous interface force and the solution with a quadrature-approximated force is of the same order as the quadrature error, utilizing fundamental solutions, singularity removal, and the Extended Trace Theorem.

Original authors: Sabia Asghar, Qiyao Peng, Etelvina Javierre, Fred J. Vermolen

Published 2026-03-13

📖 5 min read🧠 Deep dive

CC BY 4.0

Original authors: Sabia Asghar, Qiyao Peng, Etelvina Javierre, Fred J. Vermolen

Original paper licensed under CC BY 4.0 (http://creativecommons.org/licenses/by/4.0/). ✨ This is an AI-generated explanation of the paper below. It is not written or endorsed by the authors. For technical accuracy, refer to the original paper. Read full disclaimer

Imagine you are trying to understand how a soft, jelly-like tissue (like skin or a tumor) reacts when a tiny cell inside it pulls on its surroundings. This is a bit like a person standing in the middle of a trampoline and pulling the fabric down with their hands.

This paper is about a specific mathematical tool called the Immersed Interface Method. Think of this method as a way to simulate that "pulling" force without having to rebuild the entire trampoline mesh every time the person moves their hands. Instead of moving the trampoline's structure, we just tell the math, "Hey, there's a force here, and here, and here."

Here is the breakdown of what the researchers did, using simple analogies:

1. The Problem: The "Infinite" Pull

In the real world, a cell doesn't pull at just one single mathematical point; it pulls along its entire surface. Mathematically, this is described as a continuous integral—a smooth, unbroken sum of forces all along the cell's edge.

However, computers are bad at handling smooth, unbroken things. They are good at counting discrete steps. So, to make the computer understand, the researchers had to chop that smooth surface into tiny little pieces (like cutting a pizza into slices) and apply the force at the center of each slice. This is called numerical quadrature.

The Big Question: If we chop the smooth pull into tiny discrete chunks, how much does the answer change? Is the "pizza slice" version of the force close enough to the "smooth whole pizza" version?

2. The Solution: The "Ghost" Helper

The forces in this problem are tricky because they act like Dirac delta distributions. In plain English, these are forces that are infinitely strong at a single point but zero everywhere else. It's like poking a balloon with a needle so sharp that the math breaks down right at the tip.

Because of this "needle," the standard math tools (which expect smooth, gentle curves) break. The solution gets "singular" (it blows up) right where the force is applied.

To fix this, the authors used a clever trick called Singularity Removal.

The Analogy: Imagine you are trying to measure the temperature of a room, but there is a tiny, super-hot flame in the corner that makes the thermometer explode if you get too close.
The Trick: Instead of trying to measure the whole room at once, you split the problem into two parts:
1. The Flame: You calculate the effect of the hot flame using a special "fundamental solution" (a pre-made map of how heat spreads from a single point).
2. The Room: You calculate the rest of the room (the smooth part) using standard math.
3. The Sum: You add them together.

This allows them to handle the "explosive" math separately from the "smooth" math.

3. The Discovery: The "Pizza Slice" Accuracy

The main goal of the paper was to prove that chopping the force into slices (the quadrature method) doesn't ruin the answer.

They proved mathematically that:

If you use a standard "Midpoint Rule" (taking the force from the center of each slice), the error in your final answer is exactly the same size as the error in your slicing method.
If you cut the pizza into 100 slices, your error is small. If you cut it into 10,000 slices, your error gets much, much smaller.
Crucially, this holds true even in 3D (like a sphere) and even if the slice-cutting method isn't perfect.

The Metaphor: It's like trying to draw a perfect circle using straight lines. If you use a few long lines, the circle looks blocky. If you use thousands of tiny lines, it looks smooth. The researchers proved that the "blockiness" of your final drawing is directly tied to how small your lines are. You don't get any extra weird errors just because you used the Immersed Interface Method.

4. Why Does This Matter?

This might sound like abstract math, but it's vital for biology and medicine.

Cancer Research: Tumors grow by pushing and pulling on the tissue around them. To simulate how a tumor invades the body, scientists need to model these forces accurately.
Wound Healing: Cells pull on skin to close a wound.
Efficiency: The "Immersed Interface Method" is much faster than traditional methods because you don't have to rebuild the computer mesh every time a cell moves.

The Bottom Line

The authors showed that this fast, efficient way of simulating cell forces is mathematically sound. They proved that the "approximation" (cutting the force into slices) doesn't introduce any hidden, mysterious errors. The error you see is exactly what you expect from the slicing method itself.

This gives scientists the confidence to use these faster computer models to study complex biological processes, like how cancer spreads or how tissues heal, without worrying that their math is "broken" by the way they simplified the forces.

1. Problem Statement

The paper addresses the numerical accuracy of the Immersed Interface Method (IIM) (also known as the Immersed Boundary Method) when applied to linear elasticity problems involving forces exerted on internal interfaces (e.g., cell membranes or tumor boundaries) within a computational domain.

Physical Context: In biological modeling (e.g., cell migration, tumor growth), forces are often modeled as being distributed over a closed curve (2D) or surface (3D) $\Gamma$ inside a domain $\Omega$ .
Mathematical Challenge: Theoretically, these forces are represented by an integral of a Dirac delta distribution over the interface:
$-\nabla \cdot \sigma(u) = \int_{\Gamma} Q(x')n(x')\delta(x-x') d\Gamma(x')$
However, in numerical simulations, this integral must be approximated by a discrete summation (quadrature rule):
$-\nabla \cdot \sigma(u_h) = \sum_{j=1}^{m} w_j Q(x_j)n(x_j)\delta(x-x_j)$
Core Question: The authors investigate the convergence behavior between the solution of the continuous integral formulation ( $u$ ) and the discrete quadrature formulation ( $u_h$ ). Specifically, they aim to prove that the error introduced by the quadrature approximation is of the same order as the error of the quadrature rule itself, independent of the finite element discretization error.

2. Methodology

The authors employ a rigorous analytical approach combined with numerical validation, utilizing the following key techniques:

Fundamental Solutions (Green's Functions): Unlike traditional Finite Element Method (FEM) analyses that rely solely on variational formulations, this study utilizes the fundamental solutions of linear elasticity in unbounded domains ( $\mathbb{R}^2$ and $\mathbb{R}^3$ ). These solutions exhibit singular behavior at the point of force application and do not belong to the standard Sobolev space $H^1$ .
Singularity Removal Technique: To handle the non-smoothness of the solution caused by the Dirac delta distribution, the authors decompose the total displacement field $u$ $u$ into two parts:
$u = u^* + v$
- $u^*$ : The particular solution obtained by convolving the fundamental solution with the force distribution (handling the singularity).
- $v$ : An auxiliary solution that satisfies the homogeneous elasticity equation with boundary conditions derived from $u^*$ . Since $u^*$ is smooth away from the interface, $v$ lies in $H^1(\Omega)$ , allowing for standard variational analysis.
Extended Trace Theorem: The proof relies on the Extended Trace Theorem to establish the convergence of the auxiliary solutions ( $v$ and $v_h$ ) in the $L^2$ -norm, linking the boundary error (on $\partial\Omega$ ) to the interior error.
Quadrature Error Analysis: The authors analyze the error of the Midpoint Rule (and general quadrature rules) applied to the integral over the polygonal/polyhedral interface $\Gamma$ . They prove that for smooth functions, the error is bounded by $O(h^r)$ , where $h$ is the mesh size and $r$ is the order of the quadrature.
Numerical Validation: Simulations were conducted in 2D and 3D using Python (FEniCS) for FEM and custom codes for fundamental solutions. They compared:
1. The exact integral solution ( $u^*$ ).
2. The quadrature-based fundamental solution ( $u^*_h$ ).
3. The full FEM solution with singularity removal ( $u_h$ ).
4. The classical FEM solution ( $u^{FEM}_h$ ).

3. Key Contributions

Theoretical Convergence Proof: The paper establishes a rigorous proof that the $L^2$ -norm of the difference between the integral-based solution and the quadrature-based solution is of the same order as the quadrature error ( $O(h^r)$ ). This holds for both bounded and unbounded domains and for dimensions $d=2, 3$ .
Decoupling Errors: The analysis explicitly separates the error caused by the quadrature approximation of the interface force from the error caused by the spatial discretization (FEM mesh). This clarifies that the convergence rate of the immersed interface method is fundamentally limited by the quadrature rule used on the interface, not just the mesh size of the background grid.
Extension of Trace Theorems: The work successfully applies the Extended Trace Theorem to elasticity problems involving Dirac delta distributions, proving convergence in regions away from the singularity ( $\Omega_\delta$ ).
Comparison with Boundary Element Method (BEM): The authors highlight that while their method shares the use of fundamental solutions with BEM, it avoids the computational complexity of solving dense linear systems and iterative boundary condition adjustments, offering a more efficient formulation for moving geometries.

4. Results

Theoretical Results:
- Proposition 1 & 2: Proved that the difference between the exact field solution ( $u^*$ ) and the quadrature approximation ( $u^*_h$ ) converges as $O(h^r)$ in the $L^2$ -norm on any surface $S$ strictly away from the interface $\Gamma$ .
- Theorem 1 & 2: Proved that the auxiliary solutions ( $v, v_h$ ) and the total solutions ( $u, u_h$ ) converge in the $L^2$ -norm with the same order as the quadrature rule ( $O(h^r)$ ).
Numerical Results:
- 2D Simulations: Using a Midpoint rule ( $r=2$ ), the observed order of convergence for the $L^2$ -norm was approximately 2.0, matching the theoretical prediction.
- 3D Simulations: Similarly, the 3D simulations demonstrated a quadratic convergence rate ( $O(h^2)$ ) for both cases where the evaluation surface intersected the cell and where it did not.
- Singularity Removal vs. Classical FEM: The solutions utilizing the singularity removal technique ( $u_h$ ) showed superior convergence compared to classical FEM approaches ( $u^{FEM}_h$ ), which suffered from logarithmic degradation ( $O(h^2 \ln h)$ ) near singularities.

5. Significance

Validation of IIM: The study provides a mathematical foundation confirming that the Immersed Interface Method is a robust and accurate approach for modeling forces on internal interfaces in linear elasticity, provided a sufficiently accurate quadrature rule is used.
Biological Modeling: The results are directly applicable to biomedical simulations (e.g., tumor growth, wound healing) where cells exert forces on the extracellular matrix. It assures researchers that the discretization of cell boundaries does not introduce uncontrolled errors if the quadrature is handled correctly.
Future Directions: The authors note that while this work focuses on linear elasticity, the framework sets the stage for more complex models (viscoelasticity, morpho-poro-viscoelasticity). However, they acknowledge that future work will face challenges when fundamental solutions are unavailable or when governing equations become nonlinear, requiring advanced numerical strategies.

In summary, this paper bridges the gap between the theoretical formulation of immersed forces (integrals) and their numerical implementation (summations), proving that the numerical error is strictly controlled by the quadrature rule's order, thereby validating the method's precision for complex biological simulations.