---
title: Closest Point Method for Surface PDEs
url: https://www.emergentmind.com/topics/closest-point-method-cpm
type: topic
---

# Closest Point Method for Surface PDEs

The Closest Point Method (CPM) is an embedding-based numerical framework for the discretization and solution of partial differential equations (PDEs) defined on smooth, possibly open or closed, manifolds immersed in Euclidean space. CPM transfers the complexity of intrinsic surface operators to the embedding space by exploiting the closest-point extension, allowing surface derivatives to be replaced by standard Cartesian derivatives operating on closest-point-extended data. This enables the application of classical finite-difference methods, interpolation, and time-stepping schemes on narrow uniform Cartesian bands surrounding the surface, thus facilitating high-order, efficient, and robust solvers for a broad class of PDEs including elliptic, parabolic, reaction–diffusion, and variational problems. CPM supports both explicit and implicit algorithms, is straightforward to implement, parallelizes efficiently, and accommodates nonlinear, variable-coefficient, and higher-order surface operators while maintaining rigorous convergence and stability guarantees [1307.5657].

## 1. Mathematical Foundations and Closest-Point Extension

Let $\Gamma \subset \mathbb{R}^n$ be a smooth manifold (e.g., surface or curve). The method is constructed around the **closest-point mapping**
$$
\operatorname{cp}: B(\Gamma) \to \Gamma, \quad \operatorname{cp}(x) = \arg\min_{y \in \Gamma} \|x - y\|
$$
where $B(\Gamma)$ is a narrow tubular neighborhood ("band") around $\Gamma$. For any function $u$ defined on $\Gamma$, the **closest-point extension operator $E$** is defined as
$$
(Eu)(x) := u(\operatorname{cp}(x)),\quad x \in B(\Gamma).
$$
By construction, $(E^2 = E)$ and $Eu$ is constant in the normal direction to $\Gamma$.

Critical analytic identities underlie the CPM:
- **Gradient/Divergence/Laplacian principles**: For $v = Eu$, and $x$ on the surface $y = \operatorname{cp}(x)$,
  - $\nabla (Eu)(y) = \nabla_\Gamma u(y)$,
  - $\operatorname{div}(Eg)(y) = \operatorname{div}_\Gamma g(y)$,
  - For Euclidean cp: $\Delta (Eu)(y) = \Delta_\Gamma u(y)$.
  
These equivalence properties justify replacing surface differential operators with Cartesian derivatives of extended functions for all intrinsic surface operators $L_\Gamma$, including nonlinear or higher-order cases [1307.5657].

## 2. Embedding PDEs and Penalized Modified Equations

Given a time-dependent surface PDE $u_t(y, t) = L_\Gamma u(y,t)$ on $\Gamma$, CPM constructs an embedded formulation in the band:
- Extend to $v(x, t)= Eu(x, t)$. On $\Gamma$, $v_t = L_\Gamma u = [L (Eu)]|_\Gamma$. Extension via $E$ to $B(\Gamma)$ gives:
  $$
  v_t = E L v, \qquad v = E v.
  $$
- Rather than enforce $v=Ev$ as a constraint, CPM introduces a **penalized single PDE**:
  $$
  v_t = E L v - \sigma (v - Ev), \quad x \in B(\Gamma)
  $$
  for penalty parameter $\sigma > 0$. The side-by-side or constraint formulation and this penalized system are equivalent as $v - Ev \to 0$ exponentially with time [1307.5657].

## 3. Discretization: Method of Lines and Semi-Discrete System

To discretize, one introduces a uniform Cartesian grid of spacing $h$ covering the active band:
- Let $v \in \mathbb{R}^N$ denote the vector of unknowns $v_i \approx v(x_i)$ at grid nodes $x_i$.
- **Cartesian finite-difference matrices** approximate operators (e.g., $L$ for Laplacian with 5-point or 7-point stencil).
- **Discrete extension/interpolation operator $E$** at each grid node implements $(Ev)(x_i) = \sum_j w_{ij} v_j$, where $w_{ij}$ are polynomial interpolation weights for $cp(x_i)$ over a local Cartesian stencil.

The **semi-discrete method-of-lines (MOL) system** is:
$$
v'(t) = E L v - \sigma (I - E) v.
$$

Time-stepping can involve explicit Euler:
$$
v^{n+1} = (I - \Delta t\,\sigma(I - E)) v^n + \Delta t\, E L v^n
$$
or implicit schemes (e.g., implicit Euler, IMEX, BDF2, Runge–Kutta), with appropriate treatment of stiff/nonlinear/penalty terms.

The solution on $\Gamma$ at any given time is obtained by restricting $v^n$ back to the surface grid, either at grid points $x_{i(j)}$ with $\operatorname{cp}(x_{i(j)}) = y_j$ or via further interpolation [1307.5657].

## 4. Closest-Point Construction, Interpolation, and Implementation

For surfaces given as level sets ($\phi(x)=0$), the closest point map is:
$$
\operatorname{cp}(x) = x - \phi(x)\, \frac{\nabla\phi(x)}{|\nabla\phi(x)|}.
$$
For parameterized or triangulated $\Gamma$, local search or iterative projection (e.g., Newton's method) is applied to find $cp(x)$ onto the mesh [1307.5657].

The interpolation stencil for $E$ typically uses polynomial interpolation of degree $p$ (commonly $p=3$ for second-order methods). The band width must cover all finite-difference and interpolation stencils.
- $E$ is a sparse matrix (with per-row sparsity determined by stencil size).
- No explicit "reinitialization" is required: the penalty term $\sigma (I-E)$ enforces normal constancy.

## 5. Consistency, Stability, and Convergence

### 5.1 Consistency
- Finite-difference truncation in $L$ is $O(h^2)$ for the 5-point Laplacian.
- Interpolation error in $E$ is $O(h^{p+1})$.
- Penalty term error is $O(h^{p+1-\alpha})$ if $\sigma = O(h^{-\alpha})$.
- To retain second order, one typically sets $p \geq 2$ and $\alpha \leq 2$, e.g., $\sigma = 2d/h^2$ for embedding dimension $d$.

### 5.2 Stability
- The penalty $-\sigma (v-Ev)$ drives $v-Ev$ to zero exponentially: $v_t = -\sigma (v-Ev)$, so $w_t = -\sigma w$, ensuring decay for $\sigma > 0$.
- Stability restrictions for Forward Euler: $|1-\Delta t\,\sigma| \leq 1 \implies \Delta t \leq 2/\sigma$.
- Diffusion (from $E L$) imposes $\Delta t \leq h^2/(2d)$; with $\sigma = 2d/h^2$, both restrictions coincide.

### 5.3 Convergence
- For Laplace–Beltrami and canonical problems (on the circle, sphere, etc.), second-order convergence in the maximum norm is observed with $p=3, \sigma=2d/h^2, \Delta t=h^2/(4d)$, e.g., for the heat equation.
- For higher-order PDEs (biharmonic), implicit time-stepping and higher-order interpolation ($p \geq 4$) are employed, still yielding second-order convergence [1307.5657].

## 6. Model Problems and Representative Applications

**Canonical equations and their embedded CPM forms include:**
- **Heat equation**: $u_t = \Delta_\Gamma u$ becomes $v_t = E \Delta v - \sigma (v-Ev)$.
- **Surface biharmonic**: $u_t = -\Delta_\Gamma^2 u \to v_t = -E E \Delta \Delta v - \sigma (v - E v)$.
- **Reaction–diffusion systems** (e.g., Gray–Scott): Coupled CPM-embedded PDEs for each field, with IMEX time-stepping, directly applied to triangulated and genus-3 surfaces.
- **Curvature-dependent diffusion**: $u_t = \operatorname{div}_\Gamma(a(\kappa)\nabla_\Gamma u)$, with $\kappa$ (mean curvature) computed as the two-norm of the Cartesian Laplacian applied to each cp-component; standard variable-coefficient stencils used.

The approach uniformly applies to both closed and open geometries, as well as to non-linear or variable-coefficient surface differential operators [1307.5657].

## 7. Comparison With Alternative CPM Formulations

CPM as a method-of-lines unifies several earlier strands:
- **Ruuth–Merriman (2008)**: Explicit two-step alternating evolution and extension, restricted to explicit time-stepping and limited higher-order capabilities.
- **Merriman–Ruuth implicit CPM (2009)**: Introduces diagonal stabilization and Laplacian-specific penalty, less general.
- **RBF-based and mesh-free variants** (Piret, Flyer–Wright): Avoid regular grids in favor of scattered-data interpolation to increase order, at the cost of dense matrix operations.

The method-of-lines CPM:
- Permits explicit, implicit, and high-order time-integration schemes in a unified PDE framework.
- Supports nonlinear, higher-order, and variable-coefficient surface operators.
- Only requires standard sparse finite-difference stencils and sparse interpolation.
- Achieves $O(N)$ per-time-step complexity for $N$ band nodes.
- Admits rigorous stability and convergence analysis paralleling classical MOL theory [1307.5657].

**Summary table comparing CPM variants:**

| Variant               | Time-Stepping         | Surface Operators      | Implementation Features    |
|-----------------------|----------------------|-----------------------|---------------------------|
| Ruuth–Merriman (2008) | Explicit only        | Laplace-Beltrami      | Two-step extension        |
| MOL-CPM (this method) | Explicit/implicit/HO | General ($L_\Gamma$)  | Penalized single PDE, O(N)|
| RBF (Piret et al.)    | Varied               | High-order, all $L_\Gamma$ | Mesh-free, dense ops     |

The CPM in the method-of-lines framework provides a robust, parallelizable, and theoretically sound strategy for solving a broad spectrum of PDEs on arbitrary smooth surfaces using only local Cartesian operations and standard ODE solvers [1307.5657].

Source: https://www.emergentmind.com/topics/closest-point-method-cpm