Standard versus Asymptotic Preserving Time Discretizations for the Poisson-Nernst-Planck System in the Quasi-Neutral Limit

TL;DR

This work tackles the stiffness of the Poisson–Nernst–Planck system in the quasi-neutral limit by comparing time discretizations that are asymptotic-preserving with respect to the Debye length $\varepsilon$. It introduces a reformulation in terms of $\mathcal{C}$ and $\mathcal{Q}$ that yields an AP/AA time discretization when coupled with IMEX Runge–Kutta schemes, and contrasts it with a standard $($c_+,c_-)$ formulation and a standard implicit Euler method. The study provides a detailed spatial discretization via ghost finite elements on level-set domains and systematically analyzes stability and accuracy across multiple IMEX schemes, identifying the $$(\mathcal{C},\mathcal{Q})$$ formulation as the most robust for small $\varepsilon$. These results enable efficient and reliable simulations of ion transport in quasi-neutral regimes and lay groundwork for coupling with fluid dynamics in dynamic boundary settings such as oscillating bubble interfaces.

Abstract

In this paper, we investigate the correlated diffusion of two ion species governed by a Poisson-Nernst-Planck (PNP) system. Here we further validate the numerical scheme recently proposed in \cite{astuto2025asymptotic}, where a time discretization method was shown to be Asymptotic-Preserving (AP) with respect to the Debye length. For vanishingly Debye lengths, the so called Quasi-Neutral limit can be adopted, reducing the system to a single diffusion equation with an effective diffusion coefficient \cite{CiCP-31-707}. Choosing small, but not negligible, Debye lengths, standard numerical methods suffer from severe stability restrictions and difficulties in handling initial conditions. IMEX schemes, on the other hand, are proved to be asymptotically stable for all Debye lengths, and do not require any assumption on the initial conditions. In this work, we compare different time discretizations to show their asymptotic behaviors.

Figures (13)

• Figure 1: Scheme of the experimental setup. The top-left inset shows the behavior of anions and cations at the bubble air–water interface: cations (blue) have hydrophilic heads, while anions (red), with their hydrophobic tails, result to be inside the bubble and their hydrophilic heads at the surface.
• Figure 2: Discretization of the computational domain. $\Omega$ is the light blue region inside the unit square $R$. (a): classification of the grid points: the blue points are the internal ones while the red circles denote the ghost points. (b): points of intersection between the grid and the circular boundary $\Gamma_\mathcal{B}$.
• Figure 3: The AP diagram. $\mathcal{P}^\varepsilon$ is the original problem and $\mathcal{P}^\varepsilon_h$ its numerical approximation characterized by a discretization parameter $h$. The AP property corresponds to the request that $\mathcal{P}^0_h$ is consistent with $\mathcal{P}^0$ as $\varepsilon \to 0$, independently of $h$.
• Figure 4: Time accuracy orders of the system (\ref{['eq:serial']}-\ref{['eq:serial1']}) at final time $t = 0.1$, for different values of $\varepsilon$.
• Figure 5: Time accuracy orders of the system (\ref{['eq:vectorial']}-\ref{['eq:vectorial1']}) at final time $t = 0.1$, for different values of $\varepsilon$.

Theorems & Definitions (4)

• Proposition 1
• Proposition 2
• Proposition 3
• Proposition 4