Mass-lumped finite elements with maximum-bound preservation and original-energy dissipation for the degenerate-mobility Allen-Cahn equation

Introduction

The Allen-Cahn equation is a fundamental phase-field model used to describe interfacial dynamics and phase separation. It was originally introduced to describe antiphase boundary motion in crystalline solids [1]. A classic numerical analysis of the model and its relation to mean-curvature flow can be found in [5]. In many applications, the relaxation rate depends on the local phase state and is therefore modeled through as a nonlinear mobility. In this paper, we investigate the generalized Allen-Cahn model described by the following system:

Here, $\Omega\subset\mathbb{R}^d$ ($d=2,3$) is a bounded polygonal or polyhedral domain, $\mathbf{n}$ represents the outward-pointing unit normal vector, and $\varepsilon>0$ is the interfacial width. Hereafter, $(\cdot,\cdot)$ and $\|\cdot\|_0$ denote the $L^2(\Omega)$ inner product and norm, respectively. Moreover,

is the bistable potential with $f(u)=F'(u)=u^3-u$.

The mobility $M(u)$ is assumed to be nonnegative for $u\in[-1,1]$ and is allowed to degenerate at the pure phases. One example is $M(u)=(1-u^2)^m$, where $m$ is a positive integer. The considered model has two important structural properties. If $\|u_0\|_{L^\infty(\Omega)}\leq1$, the solution remains bounded by the pure phases, while the Ginzburg-Landau energy

formally satisfies

Maintaining both properties at the discrete level is even more crucial when $M$ is degenerate because violating the maximum-bound can easily cause the mobility to become negative and lose its physical interpretation.

Considerable work has been devoted to structure-preserving discretizations of Allen-Cahn-type equations. For nonlinear and possibly degenerate mobility, Shen, Tang and Yang [12] constructed maximum-principle-preserving finite-difference (FD) schemes for the generalized Allen-Cahn equation. Meanwhile, Du and Hou [3] developed a second-order, maximum-bound-preserving FD method for a generalized model with nonlinear mobility and convection. Second-order BDF and stabilized Crank-Nicolson schemes for general mobility were subsequently analyzed in [7]. More recently, Hou, Liu and Ju [9] obtained a second-order FD method that simultaneously preserves the maximum-bound principle and dissipates the original energy. Hou, Zhang and Zhu [8] developed an unconditionally maximum-bound-preserving second-order scheme. Fu et al. [6] addressed degenerate mobilities through dynamic stabilization. Their first-order method unconditionally preserves both the maximum-bound principle and energy dissipation, while the second-order formulation satisfies a conditional maximum-bound principle and a uniform original energy bound.

The finite element method (FEM) provides greater geometric flexibility, but preserving these structures simultaneously is more challenging. Barrett and Blowey [2] developed an early finite element analysis for degenerate Allen-Cahn/Cahn-Hilliard systems. For the classical Allen-Cahn equation, energy-stable FEMs have been studied in [11]. Mass-lumped, maximum-principle-preserving formulations have also been developed for surface Allen-Cahn equations [13]. Subsequently, Yang, Yi, and Zhang [14] proposed high-order mass-lumped FEM schemes that unconditionally preserve the maximum principle for the classical Allen-Cahn equation. However, FEM theory for nonlinear mobility, which can vanish at pure phases, is comparatively underdeveloped. This is especially true when a rigorous nodal maximum-bound and dissipation of the original Ginzburg-Landau energy are both required within the same fully discrete formulation.

Motivated by this observation, we construct a mass-lumped $P_1$ finite element scheme for (1) in conjunction with convex-splitting for the double-well potential. The diagonal property of mass lumping, combined with the sign structure of the stiffness matrix on non-obtuse simplicial meshes, enables us to establish a maximum-bound principle for the scheme using a direct nodal argument. Meanwhile, the convex splitting results in the original discrete energy law. We prove that, for arbitrary time step sizes, the method satisfies

In particular, the analysis remains valid when $M(\pm1)=0$. We conducted numerical experiments to assess the accuracy of the proposed scheme, confirming that it maintains both the maximum-bound principle and the energy-dissipation property.

Mass-lumped finite element approximation

Let $\mathcal{T}_h$ be a conforming and shape-regular simplicial partition of $\Omega$, consisting of triangles for $d=2$ and tetrahedra for $d=3$, and let $h=\max_{K\in\mathcal{T}_h}\operatorname{diam}(K)$. We consider the continuous piecewise linear finite element space

We denote by $\{x_i\}_{i=1}^{N_h}$ the vertices of $\mathcal{T}_h$ and by $\{\phi_i\}_{i=1}^{N_h}$ the corresponding nodal basis satisfying $\phi_i(x_j)=\delta_{ij}$. For any nodal functions $v$ and $w$, we introduce the mass-lumped inner product

and introduce the induced norm $\lVert v_h\rVert_h^2=(v_h,v_h)_h$. Thus, the associated mass matrix is diagonal with strictly positive diagonal entries. We recall that we also use the standard inner product (without the subscript $h$). Mass lumping is particularly useful in maximum-principle analysis since it preserves the nodal character of the discrete evolution [10, 14].

For the stiffness matrix $A = (a_{ij})_{i,j=1}^{N_h}$, let $a_{ij} \coloneqq(\nabla\phi_j, \nabla\phi_i)$. Throughout the analysis, we impose the standard non-obtuse mesh condition

This condition (4) ensures that the finite element stiffness matrix $A$ possesses the monotonicity structure required by a discrete maximum principle [10]. Since the nodal basis forms a partition of unity and no essential boundary degrees of freedom are eliminated under the homogeneous Neumann condition, we have

To obtain an energy-stable temporal approximation, we use a convex splitting of the double-well potential

where both $F_c$ and $F_e$ are convex. Consequently, $F_c'(s) = s^3$, and $F_e'(s) = s$. This decomposition enables us to treat the convex contribution implicitly and the concave contribution explicitly. This is the key mechanism behind the discrete energy law [4, 11].

Next, let $0 = t_0 < t_1 < \cdots< t_N = T$ be a uniform temporal partition with time step $\tau= T/N$, and define

Assuming $u_0 \in C(\overline{\Omega})$, we set $U_h^0 = I_hu_0$, where $I_h$ denotes the nodal interpolation operator. Given $U_h^n \in V_h$, the proposed mass-lumped FEM seeks to find $(U_h^{n+1}, \mu_h^{n+1}) \in V_h \times V_h$ such that, for all $v_h, q_h \in V_h$,

Here and below, nonlinear quantities appearing in $(\cdot,\cdot)_h$ are evaluated nodally. In particular, $M_i^n \coloneqq M(U_h^n(x_i))$. No strictly positive lower bound on the mobility is imposed. We only require

so that genuinely degenerate mobilities, such as

are included. Notice that if $M_i^n = 0$, then equation (7) simply implies $U_h^{n+1}(x_i) = U_h^n(x_i)$ at the corresponding node. Therefore, the discrete formulation remains well-defined at pure phases.

For later use, we record the nodal representation of (7). Choosing $v_h = q_h = \phi_i$ yields the following system:

where $U_i^n \coloneqq U_h^n(x_i)$, and $\mu_i^n \coloneqq\mu_h^n(x_i)$. Eliminating $\mu_i^{n+1}$ from (9), we obtain

Together with (4)–(5), the representation (12) will provide a direct nodal mechanism for establishing the maximum-bound principle.

Finally, corresponding to the Ginzburg–Landau energy in (2), we define the mass-lumped discrete energy

The combination of the diagonal lumped mass matrix and the stiffness-matrix sign structure yields the discrete maximum-bound. Meanwhile, the convex splitting in (8) leads directly to dissipation of the original energy (13). These two essential properties are established in the next section.

Maximum-bound and energy-dissipation properties

We will now establish the two main structural properties of the scheme (7). Throughout this section, we will assume that the mesh condition (4) and the mobility assumption (9) hold.

Well-posedness

Suppose that $M_i^n := M(U_i^n) \ge0$ for all $i$. If $M_i^n = 0$, then (10) immediately yields $U_i^{n+1}=U_i^n$. On the remaining set $\mathcal{I}_+^n := \{i:M_i^n>0\}$, eliminating $\mu_i^{n+1}$ from (9) gives

The system (14) is the Euler equation of the functional

defined on the affine subspace satisfying $V_i=U_i^n$ for $i\notin\mathcal{I}_+^n$. Since the first term in (15) is strictly convex on the active nodes and the remaining quadratic and quartic terms are convex, $\mathcal{J}_n$ is coercive and strictly convex. Therefore, it admits a unique minimizer. The chemical potential is subsequently uniquely determined from (11). Therefore, the fully discrete problem is well-posed.

Maximum-bound analysis

The following elementary property of the stiffness operator is central to the maximum-bound analysis.

Proof. Using (5), we have

Since $a_{ij}\le0$ for $j\ne i$, and at a nodal maximum, $V_j-V_i\le0$, the first assertion is proven. The argument at a nodal minimum is identical, but with the inequalities reversed.

We are now ready to establish the discrete maximum-bound principle.

Proof. We proceed by induction. The assertion is true for $n=0$. Suppose that $|U_i^n| \le1$ for all nodes $x_i$. By (9), we have $M_i^n=M(U_i^n)\ge0$. Assume, to the contrary, that

By Lemma 3.1,

Hence, from the nodal equation (12),

On the other hand, the induction hypothesis gives

which is a contradiction. Therefore, $\max_i U_i^{n+1}\le1$. Similarly, suppose that

From Lemma 3.1 we get the following

Using (12), we obtain

whereas the induction hypothesis yields

This is again a contradiction. Thus $-1\le U_i^{n+1}\le1$ at every node. Since $U_h^{n+1}$ is piecewise linear, its extrema on every simplex are attained at the vertices, and hence (16) follows. The induction is complete.

A consequence of Theorem 3.2 is that $M(U_h^n)\ge0$ at every time level, including the degenerate case $M(\pm1)=0$.

Energy-dissipation properties

We next show that the same scheme dissipates the original discrete Ginzburg-Landau energy (13).

Proof. First, we set $\delta U_h^{n+1}:=U_h^{n+1}-U_h^n$. Taking $v_h=\mu_h^{n+1}$ in (7a) and multiplying by $\tau$ gives

Taking $q_h=\delta U_h^{n+1}$ in (7b), we obtain

The gradient contribution satisfies the exact identity

Moreover, for arbitrary $a,b\in\mathbb{R}$ and $F(s)=\frac{1}{4}(s^2-1)^2$, it is straightforward that

Since

the last term in (23) is nonnegative. Applying (23) nodewise in the lumped inner product gives

Substitution of (22) and (24) into (21), followed by (20), yields (15). Finally, Theorem 3.2 ensures $M(U_i^n)\geq0$ for every node, and therefore

Then, the claimed energy dissipation property follows, along with (16), which implies (17).

Thus, the proposed mass-lumped FEM preserves the intrinsic structures of the degenerate-mobility Allen-Cahn dynamics without any restriction on the time-step size. In particular, the maximum-bound keeps the numerical solution within the physically admissible interval $[-1,1]$, while (15) shows that the original discrete Ginzburg-Landau energy, rather than a modified one, is dissipated at every time step.

Numerical experiments

In this section, we present two numerical examples that demonstrate the accuracy and structure-preserving properties of the proposed method. The computational domain is triangulated by dividing each uniform rectangular cell into two right triangles. Thus, the resulting mesh satisfies the non-obtuse condition (4). All nonlinear systems are solved by Newton iteration using $U_h^n$ as the initial guess. The iterative procedure stops once the relative residual decreases below the prescribed tolerance of $10^{-10}$. Unless otherwise stated, we use the degenerate mobility $M(u)=1-u^2$.

Example 1. (Accuracy test) First, we consider the forced problem

on $\Omega=(0,1)^2$, supplemented with homogeneous Neumann boundary conditions. The source term is chosen such that

is the exact solution. In particular, since $\Delta u=-2\pi^2u$, we take

The exact solution satisfies the homogeneous Neumann condition and remains strictly within the admissible interval $[-1,1]$. We set $\varepsilon=0.1$, $T=1$, and choose $\tau=h^2$ so that the first-order temporal error does not dominate the spatial discretization error. At the final time, we measure

The corresponding convergence rates are computed by

The numerical errors are reported in Table 1. as expected, the results exhibit the expected second-order convergence in the $L^2$ norm and first-order convergence in the $H^1$ seminorm for the piecewise linear finite element approximation.

$h$	$\tau$	$|u(T)-U_h^N|_0$	Rate	$|\nabla(u(T)-U_h^N)|_0$	Rate
$1/8$	$1.5625 \times10^{-2}$	$3.36204 \times10^{-3}$	-	$8.36695 \times10^{-2}$	-
$1/16$	$3.9062 \times10^{-3}$	$8.69831 \times10^{-4}$	$1.951$	$4.05445 \times10^{-2}$	$1.045$
$1/32$	$9.7656 \times10^{-4}$	$2.19438 \times10^{-4}$	$1.987$	$2.01082 \times10^{-2}$	$1.012$
$1/64$	$2.4414 \times10^{-4}$	$5.49916 \times10^{-5}$	$1.997$	$1.00347 \times10^{-2}$	$1.003$

Example 2. (Phase separation with degenerate mobility) Next, we examine the long-time phase evolution of the unforced Allen-Cahn problem (1) on $\Omega=(0,1)^2$. We set $\varepsilon=0.02$ with $M(u)=1-u^2$ and initialize the nodal values with a small random perturbation $U_h^0(x_i)=0.05\,\xi_i$, where $\xi_i\sim\mathcal{U}[-1,1]$. Hence,

The computation is performed up to time $T=10$ on a uniform triangulation with $h=1/128$ and time step $\tau=10^{-2}$. To examine the two structural properties established in Theorems 3.2 and 3.3, we monitor the following:

Representative phase-field configurations, along with the histories of $U_{\max}^n$, $U_{\min}^n$, and $\mathcal{E}^n$ are displayed in Figure 1. The nodal extrema remain in the admissible interval, $-1\leq U_{\min}^n\leq U_{\max}^n\leq1$, while the original discrete Ginzburg-Landau energy decreases monotonically, $\mathcal{E}^{n+1}\leq\mathcal{E}^n$, in agreement with the theoretical results. Furthermore, the mobility $M(u)=1-u^2$ decreases as $|u|$ approaches one and vanishes at the pure phases. This illustrates the effect of mobility degeneracy on the interfacial evolution.

Conclusion

We have developed a mass-lumped $P_1$ FEM for the Allen-Cahn problem involving nonlinear, potentially degenerate mobility. Combining mass lumping with a convex-splitting treatment of the double-well potential preserves the nodal maximum-bound principle in the proposed scheme, which dissipates the original discrete Ginzburg-Landau energy without any restriction on the time step size. The analysis remains valid when the mobility vanishes at the pure phases, and the discrete energy law further reveals the associated numerical dissipation mechanisms. Numerical experiments confirmed the expected finite element accuracy and verified the maximum-bound and energy-dissipation properties. These results provide a simple, robust, structure-preserving finite element framework for Allen-Cahn dynamics with degenerate mobility.