Fifth Order IIM Discretization Scheme

• Fifth Order IIM discretization scheme is a high-order finite-difference method that adapts its stencil near interfaces to preserve fifth-order spatial accuracy.
• It employs interface-aware strategies such as ghost-value extrapolation and specialized upwind or irregular stencils tailored to Navier–Stokes, Fokker–Planck, and elliptic problems.
• The method demonstrates enhanced stability and efficiency by rigorously enforcing jump conditions and subdomain-local reconstructions, achieving robust convergence even with discontinuous coefficients.

A fifth-order immersed interface method (IIM) discretization scheme is a high-order finite-difference construction designed to preserve fifth-order spatial accuracy in the presence of immersed boundaries, material interfaces, or coefficient discontinuities. In the literature provided here, the term denotes several closely related interface-aware discretizations: a fifth-order upwind IIM treatment of the advection term for incompressible Navier–Stokes equations on a collocated Cartesian grid (Ji et al., 20 Aug 2025), a fifth-order finite-difference / immersed-interface-type scheme for Fokker–Planck equations with drift-admitting jumps (Chen et al., 2019), and a fifth-order irregular-point stencil embedded in a hybrid sixth-/fifth-order method for elliptic interface problems with discontinuous and high-contrast coefficients (Feng et al., 2022). Across these formulations, the defining objective is the same: retain high-order convergence without permitting interface crossings or naive ghost closures to contaminate consistency or stability.

1. Definition and common numerical structure

In immersed-interface discretization, the underlying PDE is posed on a Cartesian mesh, while the geometric boundary or interface is allowed to cut arbitrarily through the grid. Fifth-order IIM schemes therefore modify only the interface-adjacent part of the discretization, rather than remeshing the domain.

The three formulations considered here differ in PDE type and local algebra, but they share a common structural pattern. Each scheme separates smooth subdomains from interface-affected regions, uses interface conditions explicitly, and replaces standard cross-interface differencing by interface-aware reconstruction. In the Navier–Stokes setting, the fifth-order component is the advection operator $\mathbf{u}\cdot\nabla f$, closed by sixth-order ghost-value extrapolation near immersed inflow boundaries (Ji et al., 20 Aug 2025). In the Fokker–Planck setting, fifth-order accuracy is recovered by setting each drift discontinuity as both a solution point and a flux point, and then applying interpolation and differentiation only within each smooth subdomain (Chen et al., 2019). In the elliptic setting, irregular interface-crossing points use a 13-point fifth-order stencil derived from Taylor expansion, jump conditions, and transmission formulas (Feng et al., 2022).

A concise comparison is given below.

Setting	Fifth-order component	Interface mechanism
Incompressible Navier–Stokes	Upwind advection discretization	Mixed normal/extreme stencil with sixth-order ghost extrapolation
Fokker–Planck with drift jumps	Flux and flux-derivative discretization	Interface points placed on solution and flux grids; no cross-jump stencils
Elliptic interface problems	Irregular-point finite-difference stencil	13-point interface stencil with transmission-based coefficient matching

This suggests that “fifth-order IIM discretization scheme” is not a single universal stencil, but a class of interface-corrected high-order discretizations whose fifth-order property is preserved by local interface algebra rather than by free-space finite differences alone.

2. Fifth-order advection IIM for incompressible Navier–Stokes

The most explicit use of the phrase in the supplied material is the immersed advection discretization developed for the incompressible Navier–Stokes equations

$\frac{\partial \vb{u}}{\partial t} + \vb{u}\cdot\nabla \vb{u} = -\nabla p + \nu \Delta \vb{u}, \qquad \nabla\cdot \vb{u}=0,$

with immersed boundary or interface $\Gamma$ and no-slip condition

$\vb{u}(\vb{x},t)=\vb{u}_b(\vb{x},t), \qquad \vb{x}\in\Gamma.$

The distinctive claim is that fourth-order spatial accuracy for the full solver is enabled by a novel fifth-order IIM discretization of the advection term, while diffusion, gradients, divergence, and Poisson operators use existing high-order interface-corrected centered finite differences (Ji et al., 20 Aug 2025).

In free space, the method starts from the standard fifth-order upwind finite-difference formula, applied dimension-split in each coordinate direction. In one dimension,

$\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{-2f_{k-3}+15f_{k-2}-60f_{k-1}+20f_{k}+30f_{k+1}-3f_{k+2}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{3f_{k-2}-30f_{k-1}-20f_{k}+60f_{k+1}-15f_{k+2}+2f_{k+3}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$

However, the paper shows by GKS stability analysis that simply combining this “normal” fifth-order stencil with previous high-order IIM ghost-point construction can be unconditionally unstable near immersed inflow boundaries (Ji et al., 20 Aug 2025).

The remedy is a mixed fifth-order IIM advection scheme. Within a near-boundary strip of thickness $3h$, the normal stencil is replaced by an “extreme” fifth-order upwind stencil,

$\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{3f_{k-4} - 20f_{k-3} + 60f_{k-2} - 120f_{k-1} + 65f_k + 12f_{k+1}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{-12f_{k-1} - 65f_{k} + 120f_{k+1} - 60f_{k+2} + 20f_{k+3} - 3f_{k+4}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$

This stencil remains formally fifth-order accurate, but is more dissipative and more robust near the interface. Its local truncation error constant is reported as roughly twice that of the normal stencil, but only the $3h$ near-boundary strip uses it (Ji et al., 20 Aug 2025).

The resulting Navier–Stokes spatial pairing is denoted $(5,4)$: fifth-order advection, fourth-order centered discretizations for all other operators. The formal order separation is central. The advection term is high enough not to limit the global order, whereas the fourth-order pressure, viscous, and divergence operators set the overall spatial accuracy to fourth order (Ji et al., 20 Aug 2025).

3. Interface closure, ghost values, and moving-boundary embedding

The immersed-interface closure for the Navier–Stokes advection operator is based on local polynomial extrapolation. For a boundary control point $\mathbf{x}_c$, a least-squares polynomial $\frac{\partial \vb{u}}{\partial t} + \vb{u}\cdot\nabla \vb{u} = -\nabla p + \nu \Delta \vb{u}, \qquad \nabla\cdot \vb{u}=0,$0 is constructed on a half-elliptical neighborhood of fluid points, with decaying weights away from the wall for advection in order to enhance stability. The ghost extension is then defined by

$\frac{\partial \vb{u}}{\partial t} + \vb{u}\cdot\nabla \vb{u} = -\nabla p + \nu \Delta \vb{u}, \qquad \nabla\cdot \vb{u}=0,$1

For Dirichlet data, $\frac{\partial \vb{u}}{\partial t} + \vb{u}\cdot\nabla \vb{u} = -\nabla p + \nu \Delta \vb{u}, \qquad \nabla\cdot \vb{u}=0,$2 is supplied by the boundary condition. For Neumann data, the paper uses a discrete approximation to the normal derivative,

$\frac{\partial \vb{u}}{\partial t} + \vb{u}\cdot\nabla \vb{u} = -\nabla p + \nu \Delta \vb{u}, \qquad \nabla\cdot \vb{u}=0,$3

which can be inverted to reconstruct $\frac{\partial \vb{u}}{\partial t} + \vb{u}\cdot\nabla \vb{u} = -\nabla p + \nu \Delta \vb{u}, \qquad \nabla\cdot \vb{u}=0,$4 when $\frac{\partial \vb{u}}{\partial t} + \vb{u}\cdot\nabla \vb{u} = -\nabla p + \nu \Delta \vb{u}, \qquad \nabla\cdot \vb{u}=0,$5 is prescribed (Ji et al., 20 Aug 2025).

For advection, the ghost values are built with a sixth-order extrapolation polynomial so that the fifth-order stencil can be applied across immersed inflow boundaries without loss of formal order. This sixth-order interface closure is the reason the paper characterizes the construction as a novel fifth-order IIM discretization: the finite-difference operator is fifth order, and the ghost data are sufficiently accurate to support that order near the interface (Ji et al., 20 Aug 2025).

The same local machinery extends to immersed interfaces with jumps

$\frac{\partial \vb{u}}{\partial t} + \vb{u}\cdot\nabla \vb{u} = -\nabla p + \nu \Delta \vb{u}, \qquad \nabla\cdot \vb{u}=0,$6

using separate least-squares polynomials $\frac{\partial \vb{u}}{\partial t} + \vb{u}\cdot\nabla \vb{u} = -\nabla p + \nu \Delta \vb{u}, \qquad \nabla\cdot \vb{u}=0,$7 and $\frac{\partial \vb{u}}{\partial t} + \vb{u}\cdot\nabla \vb{u} = -\nabla p + \nu \Delta \vb{u}, \qquad \nabla\cdot \vb{u}=0,$8 and a coupled system for the interface values $\frac{\partial \vb{u}}{\partial t} + \vb{u}\cdot\nabla \vb{u} = -\nabla p + \nu \Delta \vb{u}, \qquad \nabla\cdot \vb{u}=0,$9. Although this jump formulation is used generally for scalar advection-diffusion and conjugate heat transfer, the specific fifth-order result concerns stability and fifth-order accuracy near inflow boundaries (Ji et al., 20 Aug 2025).

For moving boundaries, the geometric configuration changes during the Runge–Kutta stages. The method therefore uses stagewise extrapolation and zeroing operators:

• $\Gamma$0: extrapolation with Dirichlet/no-slip information into the body,
• $\Gamma$1: extrapolation without boundary constraints for terms such as the advection right-hand side,
• $\Gamma$2: zero out values inside the body after the stage update.

The paper states that the body-CFL condition $\Gamma$3 makes the mixed space-time error terms remain high order. In the full solver, the advection operator $\Gamma$4 is embedded in a low-storage Runge–Kutta projection method through the stage update

$\Gamma$5

followed by pressure projection and velocity correction (Ji et al., 20 Aug 2025).

4. Fifth-order interface discretization for drift-admitting jumps in Fokker–Planck equations

A distinct fifth-order immersed-interface-type discretization appears in the one-dimensional Fokker–Planck problem

$\Gamma$6

where the drift $\Gamma$7 may be discontinuous at points $\Gamma$8. The physically relevant interface conditions are continuity of the propagator $\Gamma$9 and continuity of the probability current

$\vb{u}(\vb{x},t)=\vb{u}_b(\vb{x},t), \qquad \vb{x}\in\Gamma.$0

The paper’s starting observation is that an earlier fifth-order scheme for smooth drifts dropped to second order when the drift contained jumps, because some derivative approximations of the flux still used stencils that crossed the discontinuity (Chen et al., 2019).

The improved fifth-order design eliminates that failure mode by a geometric-grid reformulation: each drift discontinuity is both a solution point and a flux point, and interpolation or differencing is performed only within each smooth subdomain. For two jumps $\vb{u}(\vb{x},t)=\vb{u}_b(\vb{x},t), \qquad \vb{x}\in\Gamma.$1, the domain is decomposed into

$\vb{u}(\vb{x},t)=\vb{u}_b(\vb{x},t), \qquad \vb{x}\in\Gamma.$2

The solution grids $\vb{u}(\vb{x},t)=\vb{u}_b(\vb{x},t), \qquad \vb{x}\in\Gamma.$3 and flux grids $\vb{u}(\vb{x},t)=\vb{u}_b(\vb{x},t), \qquad \vb{x}\in\Gamma.$4 are arranged so that $\vb{u}(\vb{x},t)=\vb{u}_b(\vb{x},t), \qquad \vb{x}\in\Gamma.$5 and $\vb{u}(\vb{x},t)=\vb{u}_b(\vb{x},t), \qquad \vb{x}\in\Gamma.$6 (Chen et al., 2019).

The flux evaluation has a strictly subdomain-local form. First, $\vb{u}(\vb{x},t)=\vb{u}_b(\vb{x},t), \qquad \vb{x}\in\Gamma.$7 at flux points is computed by upwind splitting,

$\vb{u}(\vb{x},t)=\vb{u}_b(\vb{x},t), \qquad \vb{x}\in\Gamma.$8

The left and right reconstructed values $\vb{u}(\vb{x},t)=\vb{u}_b(\vb{x},t), \qquad \vb{x}\in\Gamma.$9 are obtained from fifth-order interpolation matrices $\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{-2f_{k-3}+15f_{k-2}-60f_{k-1}+20f_{k}+30f_{k+1}-3f_{k+2}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{3f_{k-2}-30f_{k-1}-20f_{k}+60f_{k+1}-15f_{k+2}+2f_{k+3}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$0, with endpoint data modified by averages at the jumps. The derivative $\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{-2f_{k-3}+15f_{k-2}-60f_{k-1}+20f_{k}+30f_{k+1}-3f_{k+2}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{3f_{k-2}-30f_{k-1}-20f_{k}+60f_{k+1}-15f_{k+2}+2f_{k+3}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$1 is then computed by fifth-order difference matrices $\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{-2f_{k-3}+15f_{k-2}-60f_{k-1}+20f_{k}+30f_{k+1}-3f_{k+2}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{3f_{k-2}-30f_{k-1}-20f_{k}+60f_{k+1}-15f_{k+2}+2f_{k+3}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$2, again subdomain by subdomain. The flux itself is

$\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{-2f_{k-3}+15f_{k-2}-60f_{k-1}+20f_{k}+30f_{k+1}-3f_{k+2}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{3f_{k-2}-30f_{k-1}-20f_{k}+60f_{k+1}-15f_{k+2}+2f_{k+3}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$3

and the interface fluxes are averaged to enforce continuity: $\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{-2f_{k-3}+15f_{k-2}-60f_{k-1}+20f_{k}+30f_{k+1}-3f_{k+2}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{3f_{k-2}-30f_{k-1}-20f_{k}+60f_{k+1}-15f_{k+2}+2f_{k+3}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$4 Finally, $\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{-2f_{k-3}+15f_{k-2}-60f_{k-1}+20f_{k}+30f_{k+1}-3f_{k+2}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{3f_{k-2}-30f_{k-1}-20f_{k}+60f_{k+1}-15f_{k+2}+2f_{k+3}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$5 is computed within each subdomain by fifth-order matrices $\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{-2f_{k-3}+15f_{k-2}-60f_{k-1}+20f_{k}+30f_{k+1}-3f_{k+2}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{3f_{k-2}-30f_{k-1}-20f_{k}+60f_{k+1}-15f_{k+2}+2f_{k+3}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$6, yielding the semi-discrete ODE system

$\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{-2f_{k-3}+15f_{k-2}-60f_{k-1}+20f_{k}+30f_{k+1}-3f_{k+2}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{3f_{k-2}-30f_{k-1}-20f_{k}+60f_{k+1}-15f_{k+2}+2f_{k+3}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$7

Time stepping uses the same third-order Runge–Kutta method as the earlier 2018 work (Chen et al., 2019).

The paper explicitly interprets the method as enforcing

$\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{-2f_{k-3}+15f_{k-2}-60f_{k-1}+20f_{k}+30f_{k+1}-3f_{k+2}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{3f_{k-2}-30f_{k-1}-20f_{k}+60f_{k+1}-15f_{k+2}+2f_{k+3}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$8

through grid placement, point-value averaging, flux averaging, and one-sided within-subdomain interpolation and differentiation. The decisive claim is that this structural change restores fifth-order convergence for discontinuous drifts (Chen et al., 2019).

5. Fifth-order irregular-point stencil in elliptic interface problems

A third fifth-order IIM formulation arises in elliptic interface problems with discontinuous and high-contrast variable coefficients: $\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{-2f_{k-3}+15f_{k-2}-60f_{k-1}+20f_{k}+30f_{k+1}-3f_{k+2}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{3f_{k-2}-30f_{k-1}-20f_{k}+60f_{k+1}-15f_{k+2}+2f_{k+3}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$9 with interface conditions

$3h$0

Here the fifth-order scheme is not the sole discretization, but the irregular-point component of a hybrid method: regular interior points use a compact 9-point stencil of sixth-order accuracy, whereas interface-cut irregular interior points use a 13-point stencil of fifth-order accuracy (Feng et al., 2022).

For an irregular point $3h$1, the local stencil is

$3h$2

The associated discrete operator is

$3h$3

Its right-hand side contains source terms and jump terms derived from interface data and transmission formulas. The paper states that this irregular stencil achieves fifth-order accuracy at irregular interior points, expressed as

$3h$4

so that the truncation error is $3h$5 in the sense of the scheme construction (Feng et al., 2022).

The coefficient construction is based on a base point $3h$6 on the interface, Taylor expansions, and a transmission formula that writes derivatives of $3h$7 in terms of derivatives of $3h$8, $3h$9, $\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{3f_{k-4} - 20f_{k-3} + 60f_{k-2} - 120f_{k-1} + 65f_k + 12f_{k+1}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{-12f_{k-1} - 65f_{k} + 120f_{k+1} - 60f_{k+2} + 20f_{k+3} - 3f_{k+4}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$0, and $\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{3f_{k-4} - 20f_{k-3} + 60f_{k-2} - 120f_{k-1} + 65f_k + 12f_{k+1}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{-12f_{k-1} - 65f_{k} + 120f_{k+1} - 60f_{k+2} + 20f_{k+3} - 3f_{k+4}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$1. The coefficients are represented as

$\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{3f_{k-4} - 20f_{k-3} + 60f_{k-2} - 120f_{k-1} + 65f_k + 12f_{k+1}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{-12f_{k-1} - 65f_{k} + 120f_{k+1} - 60f_{k+2} + 20f_{k+3} - 3f_{k+4}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$2

with $\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{3f_{k-4} - 20f_{k-3} + 60f_{k-2} - 120f_{k-1} + 65f_k + 12f_{k+1}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{-12f_{k-1} - 65f_{k} + 120f_{k+1} - 60f_{k+2} + 20f_{k+3} - 3f_{k+4}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$3 for the irregular-point fifth-order construction, and the moment-matching conditions produce a linear system $\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{3f_{k-4} - 20f_{k-3} + 60f_{k-2} - 120f_{k-1} + 65f_k + 12f_{k+1}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{-12f_{k-1} - 65f_{k} + 120f_{k+1} - 60f_{k+2} + 20f_{k+3} - 3f_{k+4}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$4. For the 13-point scheme, the paper reports that $\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{3f_{k-4} - 20f_{k-3} + 60f_{k-2} - 120f_{k-1} + 65f_k + 12f_{k+1}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{-12f_{k-1} - 65f_{k} + 120f_{k+1} - 60f_{k+2} + 20f_{k+3} - 3f_{k+4}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$5 is initially a $\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{3f_{k-4} - 20f_{k-3} + 60f_{k-2} - 120f_{k-1} + 65f_k + 12f_{k+1}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{-12f_{k-1} - 65f_{k} + 120f_{k+1} - 60f_{k+2} + 20f_{k+3} - 3f_{k+4}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$6 matrix and is reduced to $\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{3f_{k-4} - 20f_{k-3} + 60f_{k-2} - 120f_{k-1} + 65f_k + 12f_{k+1}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{-12f_{k-1} - 65f_{k} + 120f_{k+1} - 60f_{k+2} + 20f_{k+3} - 3f_{k+4}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$7 using an identity among coefficient sums (Feng et al., 2022).

Within this hybrid framework, the fifth-order IIM role is sharply localized: it resolves the interface-cut points where compact sixth-order regular stencils are not available. This suggests that fifth-order IIM discretization can function either as a global interface method, as in the Fokker–Planck scheme, or as a localized irregular-point correction within a higher-order composite method, as in the elliptic scheme.

6. Accuracy, stability, and computational implications

The supplied papers emphasize that fifth-order immersed-interface accuracy is meaningful only if the interface closure does not degrade order or destabilize the discretization.

For the Navier–Stokes advection scheme, the principal stability issue is that the naive combination of the standard fifth-order upwind stencil with high-order IIM ghosting is unconditionally unstable near immersed inflow boundaries by GKS analysis. The mixed strategy addresses this by using the standard fifth-order stencil in the interior and the extreme fifth-order stencil only within $\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{3f_{k-4} - 20f_{k-3} + 60f_{k-2} - 120f_{k-1} + 65f_k + 12f_{k+1}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{-12f_{k-1} - 65f_{k} + 120f_{k+1} - 60f_{k+2} + 20f_{k+3} - 3f_{k+4}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$8 of immersed inflow boundaries. The paper states that the resulting IIM advection scheme achieves comparable accuracy and stability to a standard free-space fifth-order upwind finite-difference stencil (Ji et al., 20 Aug 2025). In advection-only tests with a star-shaped obstacle, fifth-order convergence in both $\left(\frac{\partial f}{\partial x}\right)_k = \begin{dcases} \frac{3f_{k-4} - 20f_{k-3} + 60f_{k-2} - 120f_{k-1} + 65f_k + 12f_{k+1}}{60h} + \mathcal{O}(h^5), & u_k \ge 0,\[6pt] \frac{-12f_{k-1} - 65f_{k} + 120f_{k+1} - 60f_{k+2} + 20f_{k+3} - 3f_{k+4}}{60h} + \mathcal{O}(h^5), & u_k < 0. \end{dcases}$9 and $3h$0 norms is reported for uniform translation in a periodic domain and for rigid rotation inside a closed star-shaped domain. In the full Navier–Stokes tests, the $3h$1 scheme gives fifth-order velocity convergence in free space when advection dominates, fourth-order pressure convergence, $3h$2-order convergence of divergence under fixed CFL, and fourth-order accuracy on stationary and moving immersed boundaries. The paper also reports an efficiency comparison: in a pitching-plate example, the fourth-order scheme reaches convergence at roughly half to two-thirds of the linear resolution of a second-order scheme, translating to roughly an order-of-magnitude savings in 2D spatio-temporal degrees of freedom under fixed CFL (Ji et al., 20 Aug 2025).

For the Fokker–Planck scheme, the main issue is not a separate stability theorem but the consistency failure created by stencils that straddle drift discontinuities. Once those crossings are removed, the numerical evidence shows fifth-order behavior for both smooth and discontinuous drifts (Chen et al., 2019). In the pure dry friction example with drift $3h$3, domain $3h$4, jump at $3h$5, $3h$6, $3h$7, zero current boundary conditions, and time step

$3h$8

the reported $3h$9 errors at $(5,4)$0 decrease from $(5,4)$1 for $(5,4)$2 to $(5,4)$3 for $(5,4)$4, with rates approaching $(5,4)$5. The corresponding $(5,4)$6 errors decrease from $(5,4)$7 to $(5,4)$8, approximately fifth order as well (Chen et al., 2019).

For the elliptic hybrid method, the paper does not provide a full stability theorem, but reports that the method remains robust for coefficient ratios such as $(5,4)$9 and $\mathbf{x}_c$0 (Feng et al., 2022). Numerical experiments confirm sixth-order global behavior for the hybrid assembly while using the fifth-order irregular-point stencil as designed. In Example 3 with contrast ratio $\mathbf{x}_c$1, the reported $\mathbf{x}_c$2-error orders are around $\mathbf{x}_c$3–$\mathbf{x}_c$4, and the $\mathbf{x}_c$5-error orders are around $\mathbf{x}_c$6–$\mathbf{x}_c$7. The paper also notes that replacing the 13-point irregular stencil by a 9-point irregular stencil significantly increases the error, while only slightly reducing the condition number (Feng et al., 2022).

Taken together, these results indicate that fifth-order IIM discretization is most effective when the interface treatment is made structurally compatible with the underlying PDE operator: ghost extrapolation and mixed upwinding for advective immersed boundaries, subdomain-local solution/flux grids for discontinuous drifts, and transmission-matched irregular stencils for elliptic jumps.

7. Position within high-order interface methods

Within the supplied literature, fifth-order IIM discretization occupies an intermediate but consequential position between lower-order interface methods and globally sixth-order regular-grid schemes. In the Navier–Stokes paper, fifth-order advection is deliberately paired with fourth-order operators so that the advection term does not become the accuracy bottleneck (Ji et al., 20 Aug 2025). In the elliptic paper, fifth-order irregular-point discretization is paired with sixth-order regular and boundary stencils because compact high-order formulas at interface-cut points are more restrictive (Feng et al., 2022). In the Fokker–Planck paper, fifth-order accuracy is recovered not by enlarging the order of the time integrator, but by reorganizing the interface geometry of the discretization so that no interpolation or differentiation crosses a jump (Chen et al., 2019).

A common misconception is that free-space fifth-order finite differences automatically retain fifth-order performance after an immersed-interface closure is added. The supplied papers show otherwise. In one case, naive ghost-point coupling is unconditionally unstable near immersed inflow boundaries (Ji et al., 20 Aug 2025); in another, flux differentiation across a drift jump reduces a formally fifth-order method to second order (Chen et al., 2019). The elliptic work likewise motivates the use of a dedicated irregular-point stencil because lower-complexity interface discretizations suffer from larger errors at cut cells (Feng et al., 2022).

Accordingly, a fifth-order IIM discretization scheme is best understood as a stabilized, interface-consistent high-order closure rather than merely a fifth-order free-space stencil placed next to an interface. Its distinguishing feature is the preservation of fifth-order local accuracy through explicit handling of jump conditions, ghost values, interface geometry, or one-sided reconstruction, depending on the governing PDE and the nature of the embedded discontinuity.