A monotone finite element method for reaction-drift-diffusion equations with discontinuous reaction coefficients

TL;DR

This work addresses the challenge of discretizing reaction-drift-diffusion equations with discontinuous reaction coefficients while preserving monotonicity. It develops a dual-mesh, finite-volume-type quadrature that generalizes mass lumping within a multilinear tensor-product finite-element framework, yielding a nonnegative diagonal reaction operator and $O(h^2)$ consistency even at jumps. The authors prove local and global error bounds, refine interface estimates using uniform trace inequalities, and extend the approach to a fully monotone tensor-product scheme (EAFE), with rigorous $H^1$ and $L^2$ convergence results and a demonstrated supercloseness to the Galerkin solution. Numerical experiments in PBSRD contexts validate the theory, showing global $L^2$ convergence rates of $O(h^2)$ and supporting the practical robustness of the method for discontinuous coefficients and evolving interfaces.

Abstract

We prove an abstract convergence result for a family of dual-mesh based quadrature rules on tensor products of simplical meshes. In the context of the multilinear tensor-product finite element discretization of reaction-drift-diffusion equations, our quadrature rule generalizes the mass-lump rule, retaining its most useful properties; for a nonnegative reaction coefficient, it gives an $O(h^2)$-accurate, nonnegative diagonalization of the reaction operator. The major advantage of our scheme in comparison with the standard mass lumping scheme is that, under mild conditions, it produces an $O(h^2)$ consistency error even when the integrand has a jump discontinuity. The finite-volume-type quadrature rule has been stated in a less general form and applied to systems of reaction-diffusion equations related to particle-based stochastic reaction-diffusion simulations (PBSRD); in this context, the reaction operator is \textit{required} to be an $M$-matrix and a standard model for bimolecular reactions has a discontinuous reaction coefficient. We apply our convergence results to a finite element discretization of scalar drift-diffusion-reaction model problem related to PBSRD systems, and provide new numerical convergence studies confirming the theory.

Figures (4)

• Figure 1: Left: A finite element primal mesh on a disk (black) and its dual (blue, dashed). Right: the supporting primal and dual mesh elements for an interior mesh node.
• Figure 2: Convergence rate of the reaction-drift-diffusion finite element solution computed using the lumping method (red line, circular) markers, and the averaging method (blue line, square markers), compared with an $O(h^2)$ rate of convergence, with $\overline{\kappa} = 1$. The parameter $\overline{\kappa}$, represents the strength of the potential force and also influences the magnitude of the gradient for the smooth part of the reaction kernel.
• Figure 3: Convergence rate of the reaction-drift-diffusion finite element solution computed using the lumping method (red line, circular) markers, and the averaging method (blue line, square markers), compared with an $O(h^2)$ rate of convergence, with $\overline{\kappa} = 5$. The parameter $\overline{\kappa}$, represents the strength of the potential force and also influences the magnitude of the gradient for the smooth part of the reaction kernel.
• Figure 4: $L^2$-norms of $u_h - u_h^G$ (red line, circular markers) and $\nabla(u_h - u_h^G)$ (blue line, square markers) for the 2D (left) and 3D (right) problem, where $u_h^G$ is the Galerkin solution of the drift-diffusion equation and $u_h$ is computed using the finite-volume quadrature rule. The $H^1$ errors are $O(h^{\frac{3}{2}})$, a half-order reduction in the convergence rate compared to when the reaction coefficient is smooth.

Theorems & Definitions (35)

• Remark 2.1
• Definition 2.2
• Lemma 2.3: Bramble-Hilbert for bilinear forms
• Lemma 2.4
• proof
• Lemma 2.5
• proof
• Remark 2.6
• Remark 2.7
• Remark 2.8