A Correction Function Method for Poisson Problems with Interface Jump Conditions

In this paper we present a method to treat interface jump conditions for constant coefficients Poisson problems that allows the use of standard “black box” solvers, without compromising accuracy. The basic idea of the new approach is similar to the Ghost Fluid Method (GFM). The GFM relies on corrections applied on nodes located across the interface for discretization stencils that straddle the interface. If the corrections are solution-independent, they can be moved to the right-hand-side (RHS) of the equations, producing a problem with the same linear system as if there were no jumps, only with a different RHS. However, achieving high accuracy is very hard (if not impossible) with the “standard” approaches used to compute the GFM correction terms. In this paper we generalize the GFM correction terms to a correction function, defined on a band around the interface. This function is then shown to be characterized as the solution to a PDE, with appropriate boundary conditions. This PDE can, in principle, be solved to any desired order of accuracy. As an example, we apply this new method to devise a 4th order accurate scheme for the constant coefficients Poisson equation with discontinuities in 2D. This scheme is based on (i) the standard 9-point stencil discretization of the Poisson equation, (ii) a representation of the correction function in terms of bicubics, and (iii) a solution of the correction function PDE by a least squares minimization. Several applications of the method are presented to illustrate its robustness dealing with a variety of interface geometries, its capability to capture sharp discontinuities, and its high convergence rate.

💡 Research Summary

The paper introduces a novel “Correction Function Method” (CFM) for solving constant‑coefficient Poisson equations that contain discontinuities across an interface. Traditional approaches such as the Immersed Interface Method (IIM) or the Ghost Fluid Method (GFM) either modify the finite‑difference stencil (leading to non‑symmetric, more complex linear systems) or use pre‑computed correction terms that are typically only first‑order accurate. The authors generalize the GFM idea by replacing pointwise correction terms with a continuous correction function D(x) defined on a thin band ΩΓ surrounding the interface Γ.

Mathematically, D(x) = u⁺(x) – u⁻(x) satisfies a Poisson‑type PDE inside the band: ∇²D = f⁺ – f⁻, together with Dirichlet and Neumann jump conditions D = a and ∂ₙD = b on Γ, where a and b are the prescribed jumps in the solution and its normal derivative, and f⁺, f⁻ are smooth extensions of the source term. Although this is an elliptic Cauchy problem (generally ill‑posed), the authors argue that in the numerical context the high‑frequency components of the data are filtered and the solution is only required in a distance O(h) from the interface, which mitigates the instability.

Implementation proceeds in three stages. First, the interface is represented implicitly by a level‑set function φ(x) and its gradient is reconstructed using a Gradient‑Augmented Level Set (GA‑LS) technique, providing high‑order curvature and normal information without particle tracking. Second, within ΩΓ the correction function is approximated by bicubic polynomials (16 coefficients per cell), ensuring C¹ continuity across cell boundaries. Third, the coefficients are obtained by a least‑squares minimization of the residuals of the PDE (∇²D – (f⁺–f⁻)) and the interface conditions, weighted according to the relative position of grid points and the interface.

Because D(x) can be evaluated at any point near Γ, the standard 9‑point second‑order central difference stencil for the Laplacian can be used unchanged throughout the domain. The only modification is the addition of D‑dependent terms to the right‑hand side of the discrete Poisson equation. Consequently, the global linear system retains the symmetric, banded structure of the original Poisson problem, allowing the use of existing fast solvers (multigrid, FFT‑based, etc.).

Two concrete schemes are built to demonstrate the framework: a second‑order and a fourth‑order accurate method. The fourth‑order scheme employs the same bicubic representation but solves the correction‑function PDE with fourth‑order accuracy and uses a fourth‑order nine‑point Laplacian (including cross‑derivative terms). Numerical experiments on a variety of interface geometries—circles, ellipses, highly curved shapes, and multiple intersecting interfaces—show optimal convergence rates in both L₂ and maximum norms. The method remains robust when the interface is nearly aligned with the grid or when it cuts through grid cells at very small angles, situations that typically degrade the accuracy of dimension‑by‑dimension Taylor‑expansion based GFM corrections.

Compared with IIM, CFM does not require stencil modification or the assembly of non‑symmetric matrices, yet it achieves comparable or higher order accuracy and sharp interface representation. The approach is also dimension‑agnostic; extending to three dimensions merely involves using trivariate bicubic (or triquartic) basis functions and solving the analogous correction‑function PDE, making CFM a flexible tool for multiphysics problems where Poisson‑type equations with jump conditions arise (e.g., electrostatics with dielectric interfaces, heat conduction across composite materials).

In summary, the Correction Function Method provides a systematic way to embed interface jump conditions into standard Poisson solvers, preserving the efficiency of black‑box linear algebra while delivering high‑order accuracy and robust handling of complex interface geometries.