A posteriori error estimates for the Allen-Cahn problem

TL;DR

This paper establishes a posteriori error estimates for fully-discrete Galerkin methods solving the Allen-Cahn equation, providing bounds that depend polynomially on the inverse interface length and improving existing error bounds.

Contribution

It introduces new a posteriori error estimates for the Allen-Cahn problem using spectral estimates and elliptic reconstruction, enhancing previous bounds in various norms.

Findings

Error estimates depend polynomially on inverse interface length

Spectral estimate for linearized operator is crucial

Improved bounds in $L_2(H^1)$ and $L_ extinfty(L_2)$ norms

Abstract

This work is concerned with the proof of \emph{a posteriori} error estimates for fully-discrete Galerkin approximations of the Allen-Cahn equation in two and three spatial dimensions. The numerical method comprises of the backward Euler method combined with conforming finite elements in space. For this method, we prove conditional type \emph{a posteriori} error estimates in the $L^{}_4(0,T;L^{}_4(\Omega))$-norm that depend polynomially upon the inverse of the interface length $\epsilon$. The derivation relies crucially on the availability of a spectral estimate for the linearized Allen-Cahn operator about the approximating solution in conjunction with a continuation argument and a variant of the elliptic reconstruction. The new analysis also appears to improve variants of known \emph{a posteriori} error bounds in $L_2(H^1)$, $L_\infty^{}(L_2^{})$-norms in certain regimes.

Equations

$$u_t-\Delta u+\frac{1}{{ϵ}^2}(u^3-u)$$

$$u$$

$$u(\cdot ,0)$$

$$\displaystyle\|v\|_{L_{p}(0,T;V)}=\Big{(}\int_{0}^{T}\|v\|^{p}_{V}\mathop{}\!\mathrm{d}t\Big{)}^{1/p},\ p\in[1,+\infty),\quad\|v\|_{L_{\infty}(0,T;V)}=\operatorname*{ess.sup}_{t\in[0,T]}\|v\|_{V},$$

$$\Vert v{\Vert }_{L_4(\Omega )}$$

$$ab\le \delta a^p+C(p,q){\delta }^{-\frac{q}{p}}{b}^q,\text{where}1/p+1/q=1,$$

$$⟨u_t(t),v⟩+\left(\nabla u(t),\nabla v\right)+{ϵ}^{-2}\left(u^3(t)-u(t),v\right)$$

$$(u(0),v)$$

$$\displaystyle(u(T),v(T))+\int_{0}^{T}\Big{(}-\langle u,v_{t}\rangle+(\nabla u,\nabla v)+\epsilon^{-2}(u^{3}-u,v)\Big{)}\mathop{}\!\mathrm{d}t$$

$$=$$

$$V_h^n:=\{\chi \in C(\bar{\Omega });\chi {|}_{\tau }\in {\mathbb{P}}_{\kappa }(\tau ),\forall \tau \in {\mathcal{T}}_h\},$$

$$\displaystyle V^{n}_{hk}:=\left\{X:[0,T]\rightarrow V^{n}_{h};X\in L_{2}(0,T;H^{1}_{0}(\Omega));\ X|_{J_{n}}\in\mathbb{P}_{0}\big{[}J_{n};V^{n}_{h}\big{]}\right\},$$

$$k_n^{-1}\left(U_h^n-U_h^{n-1},X\right)+\left(\nabla U_h^n,\nabla X\right)+{ϵ}^{-2}\left(F(U_h^n),X\right)$$

$$U_h^0$$

$$k_n^{-1}\left(U_h^n-{\mathcal{P}}_h^n{U}_h^{n-1}\right)-{\Delta }_h^n{U}_h^n+{ϵ}^{-2}{\mathcal{P}}_h^{n}F(U_h^n)={\mathcal{P}}_h^n{f}^n.$$

$$\left(\nabla {\omega }^n,\nabla v\right)=⟨g_h^n,v⟩,\text{ for all}v\in H_0^1(\Omega ),$$

$$g_h^n:$$

$$-k_n^{-1}\left({\mathcal{P}}_h^n{U}_h^{n-1}-U_h^{n-1}\right).$$

$$\left(\nabla ({\omega }^n-U_h^n\right),\nabla X)=0,\text{ for all}X\in V_{hk}^n.$$

$$U_h(t):={\ell }_{n-1}(t)U_h^{n-1}+{\ell }_n(t)U_h^n,\text{and}\omega (t):={\ell }_{n-1}(t){\omega }^{n-1}+{\ell }_n(t){\omega }^n,$$

$$U_{h,t}(t)=\frac{U_h^n-U_h^{n-1}}{k_n},$$

$$e:=u-U_h=\theta -\rho ,\text{where}\theta :=\omega -U_h,\rho :=\omega -u.$$

$$⟨{\rho }_t,v⟩+\left(\nabla \rho ,\nabla v\right)+{ϵ}^{-2}\left(F(U_h)-F(u),v\right)$$

$$=$$

$$\left(\nabla (\omega -{\omega }^n),\nabla v\right)$$

$$\displaystyle\quad\quad\ +\|f^{n}-f^{n-1}\|_{L_{2}(\Omega)}\Big{)}\|v\|_{L_{2}(\Omega)}$$

$$(\nabla (\omega -{\omega }^n),\nabla v)={\ell }_{n-1}(t)\left(\nabla ({\omega }^{n-1}-{\omega }^n),\nabla v\right)$$

$$=$$

$$g_h^n$$

$$=k_n^{-1}\left(U_h^{n-1}-U_h^n\right)-{ϵ}^{-2}F(U_h^n)+f^n,$$

$${\mathcal{L}}_1:$$

$${\mathcal{L}}_2:$$

$${\Theta }_1\equiv {\Theta }_1(t):$$

$${\Theta }_2\equiv {\Theta }_2(t):$$

$$A(t):$$

$$\frac{1}{4}{\int }_0^{\text{uptau}}\Vert \rho {\Vert }_{L_4(\Omega )}^4\mathrm{d}t+\frac{1}{8}\Vert {\int }_0^{\text{uptau}}\nabla {\rho }^2(s)\mathrm{d}s{\Vert }_{L_2(\Omega )}^2+\frac{1}{2}\Vert \rho (\text{uptau}){\Vert }_{L_2(\Omega )}^2$$

Full text

1. Introduction

The Allen-Cahn problem comprises of a singularly perturbed parabolic semi-linear parabolic partial differential equation (PDE) together with suitable initial and boundary conditions, viz.,

[TABLE]

we assume that $\Omega\subset\mathbb{R}^{d}$ is a convex, polygonal ($d=2$) or polyhedral ($d=3$) domain of the Euclidean space $\mathbb{R}^{d}$, $T\in\mathbb{R}^{+}$, $0<\epsilon\ll 1$, for sufficiently smooth initial condition $u_{0}$ and forcing function $f$ (precise regularity statements will be given below).

The problem (1) belongs to the class of the so-called phase field PDEs models for solidification of a pure material, originally introduced by Allen & Cahn [3] to describe the phase separation process of a binary alloy at a fixed temperature. The nonlinear function $F(u):=u^{3}-u$ is the derivative of the classical double-well potential $\int F(u){\rm d}u$. Due to the nature of the non-linearity, the solution $u$ develops time-dependent interfaces $\Gamma_{t}:=\{x\in\Omega:u(x,t)=0\}$, separating regions for which $u\approx 1$ from regions where $u\approx-1$. The solution moves from one region to another within the, so-called, diffuse interfaces of length $O(\epsilon)$. For a recent comprehensive review of phase field models and their relationship to geometric flows, we refer to [12].

Realistically, $\epsilon$ should be orders of magnitude smaller than the physical domain of simulation. Therefore, the accurate and efficient numerical solution of such phase field models requires the resolution of the dynamic diffuse interfaces. This means that the discretisation parameters of any numerical method used should provide sufficient numerical resolution to approximate the interface evolution accurately. In the context of finite element methods, this is typically achieved via the use of very fine meshes in the vicinity of the interface region. In an effort to simulate at a tractable computational cost, especially for $d=3$, it is essential to design adaptive algorithms which are able to dynamically modify the local mesh size.

A standard error analysis of finite element approximations of (1) leads to a priori estimates with unfavourable exponential dependence on $\epsilon^{-1}$. This is impractical even for moderately small interface length $\epsilon$. The celebrated works [8, 10, 2] showed that uniform bounds for the principal eigenvalue of the linearized Allen-Cahn spatial operator about the solution $u$ are possible as long as the evolving interface is smooth (cf., (20) below). Such spectral estimates are used in the seminal work [14] whereby a priori error bounds with only polynomial dependence on $\epsilon^{-1}$ for finite element methods have been proven, enabling also the proof of convergence to the sharp-interface limit. Moreover, assuming the validity of a spectral estimate about the exact solution $u$, allowed the proof of the first conditional-type a posteriori error bounds for finite element methods approximating the Allen-Cahn problem in $L_{2}(H^{1})$-norm, for which the condition depends only polynomially on $\epsilon$; this was presented in the influential works [18, 13].

This direction of research has taken a further leap forward with the seminal works [4, 6, 7], whereby the principle eigenvalue of the linearized spatial Allen-Cahn operator about the numerical solution $U_{h}$ is used instead, in an effort to arrive to fully computable a posteriori error estimates in $L_{2}(H^{1})$- and $L_{\infty}(L_{2})$-norms, the latter using the elliptic reconstruction framework [20, 19]. We also mention [16] whereby a posteriori error bounds in the $L_{\infty}(L_{r})$-norms, $r\in[2,\infty]$ are proven.

When the interface $\Gamma_{t}$ undergoes topological changes, however, e.g., when an interface collapses, unbounded velocities occur and the all-important principal eigenvalue $\lambda$ can scale like $\lambda\sim\epsilon^{-2}$ on a time interval of length comparable to $\epsilon^{2}$. This crucial observation, made in [7], showed that the principal eigenvalue can be assumed to be $L_{1}$-integrable with respect to the time variable allowing, in turn, for robust conditional a posteriori error analysis under topological changes in $L_{2}(H^{1})$- and $L_{\infty}(L_{2})$-norms.

In a recent work [9], a priori bounds for the $L_{4}(L_{4})$-norm error have been proved, which appear to deliver a rather favourable $\epsilon^{-1}$-polynomial dependence on the respective constant, noting that $L_{4}(L_{4})$-norm is present in the stability of the spatial Allen-Cahn operator upon multiplication of (1) by $u$ and integration with respect to space and to time. An immediate question is whether proving conditional a posteriori error bounds in $L_{4}(L_{4})$-norm norm can also improve the dependence of the condition on the interface length $\epsilon$. Motivated by this, in this work, we prove conditional a posteriori error bounds for the $L_{4}(L_{4})$-norm for a backward Euler-finite element method. The proof is valid under the hypothesis of the existence of a spectral estimate under topological changes in the spirit of [7]. The argument uses a carefully constructed test function, in conjunction with a continuation argument and a new variant of the elliptic reconstruction introduced in [15]. As a result of the method of proof, the new a posteriori error analysis provides also new $L_{\infty}(L_{2})$- and $L_{2}(H^{1})$-norm a posteriori error bounds which appear to, at least formally, be valid under less stringent smallness condition compared to results from the literature.

The remainder of this work is structured as follows. The model problem is introduced In Section 2. Section 3 include the definition of the numerical method along with the elliptic and time reconstructions needed for the proof of the main results. The key estimates and the main result are stated and proven in Section 4. Section 5 completes the derivation of fully computable error bounds by estimating the terms appearing in the residuals of the main results.

2. Model problem

We denote by $L_{p}(\Omega)$, $1\leq p\leq\infty$ the standard Lebesgue spaces with corresponding norms $\|\cdot\|_{L_{p}(\Omega)}$. Let also $W^{k,p}(\Omega)$ is the $k$th order of Sobolev space based on $L_{p}(\Omega)$ and $H^{k}(\Omega):=W^{k,2}(\Omega)$, $k\geq 0$, along with the corresponding norms $\|\cdot\|_{W^{k,p}(\Omega)}$ and $\|\cdot\|_{H^{k}(\Omega)}$, respectively. Set $H_{0}^{1}(\Omega):=\{v\in H^{1}(\Omega):v|_{\partial\Omega}=0\}$. We shall denote by $\langle\cdot,\cdot\rangle$ the duality pairing between $H^{-1}(\Omega)$ and $H^{1}_{0}(\Omega)$, which becomes the standard $L_{2}(\Omega)$ inner product $\left(\cdot,\cdot\right)$ when the arguments are sufficiently smooth. The respective Bochner spaces are denoted by $L_{p}(0,T;V)$, endowed with the norms:

[TABLE]

with $V$ a Banach space with norm $\|\cdot\|_{V}$.

We shall make extensive use of the classical Gagliardo-Nirenberg-Ladyzhenskaya inequalities (GNL) reading:

[TABLE]

for all $v\in H^{1}_{0}(\Omega)$ with $\tilde{c}>0$, independent of $v$. For later use, we also recall a basic algebraic estimate, often referred to as the Young’s inequality: for any $\delta>0$, we have

[TABLE]

for any $a,b\geq 0$ and $p,q>1$, for some $C(p,q)>0$ independent of $a,b$.

Let $f\in L_{\infty}(0,T;L_{4}(\Omega))$ and $u_{0}\in L_{\infty}(\Omega)$. Then, for a.e. $t\in(0,T]$, we seek $u\in L_{2}(0,T;H_{0}^{1}(\Omega))\cap H^{1}(0,T;H^{-1}(\Omega))$, such that

[TABLE]

for all $v\in H_{0}^{1}(\Omega).$ Integrating for $t\in(0,T]$, and integrating by parts the above becomes: find $u\in L_{2}(0,T;H_{0}^{1}(\Omega))\cap L_{\infty}(0,T;L_{2}(\Omega))$, such that

[TABLE]

for all $v\in L_{2}(0,T;H_{0}^{1}(\Omega))\cap H^{1}(0,T;H^{-1}(\Omega))$.

3. The fully discrete scheme and reconstructions

We shall first present a fully discrete scheme for the Allen-Cahn problem (1) by combining the lowest order discontinuous Galerkin time-stepping method with conforming finite elements in space. Further, we shall define suitable space and time reconstructions of the fully discrete scheme, which will be crucial for the proof of the a posteriori error bounds below.

3.1. Discretisation

Let $0=t_{0}<t_{1}<\ldots<t_{N}=T$. We partition the time interval $[0,T]$ into subintervals $J_{n}:=(t_{n-1},t_{n}]$ and we denote by $k_{n}:=t_{n}-t_{n-1}$, $n=1,\ldots,N$ each time step.

Let also $\{\mathcal{T}^{n}_{h}\}_{n=0}^{N}$ be a sequence of conforming and shape-regular triangulations of the domain $\Omega$, that are allowed to be modified between time steps. We define the meshsize function, $h_{n}:\Omega\rightarrow\mathbb{R}$, by $h_{n}(x):=\text{diam}(\tau)$, $x\in\tau$ for $\tau\in\mathcal{T}^{n}_{h}$. To each $\mathcal{T}^{n}_{h}$ we associate the finite element space:

[TABLE]

with $\mathbb{P}_{\kappa}$ denoting the $d$-variate space of polynomials of degree at most $\kappa\in\mathbb{N}$. The whole theory presented below remains valid if box-type elements are used and respective polynomial spaces of degree $\kappa$ on each variable.

We say that a set of triangulations is compatible when they are constructed by different refinements of the same (coarser) triangulation. Given two compatible triangulations $\mathcal{T}^{n-1}_{h}$ and $\mathcal{T}^{n}_{h}$, we consider their finest common coarsening $\hat{\mathcal{T}}^{n}_{h}:=\mathcal{T}^{n}_{h}\wedge\mathcal{T}^{n-1}_{h}$ and set $\hat{h}_{n}:=\max(h_{n},h_{n-1})$. The partial order relation among the triangulations implies that $V^{n-1}_{h}\subset V^{n}_{h}$. Furthermore, we denote by $\mathcal{S}^{n}_{h}$ the interior mesh skeleton of $\mathcal{T}^{n}_{h}$, and we define the sets $\hat{\mathcal{S}}^{n}_{h}:=\mathcal{S}^{n}_{h}\cap\mathcal{S}^{n-1}_{h}$ and $\check{\mathcal{S}}^{n}_{h}:=\mathcal{S}^{n}_{h}\cup\mathcal{S}^{n-1}_{h}.$

Approximations will be subordinate to the time partition. A finite element space $V^{n}_{h}\subset H^{1}_{0}(\Omega)$ is specified on each time interval $J_{n}$, $n=1,\ldots,N$. Then, we seek approximate solutions from the space

[TABLE]

with \mathbb{P}_{0}\big{[}J_{n};V_{h}\big{]} denoting the space of constant polynomials over $J_{n}$, having values in $V^{n}_{h}$; these functions are allowed to be discontinuous at the nodal points, but are taken to be continuous from the left.

3.2. Fully discrete scheme

For brevity, we set $F(v):=v^{3}-v$. The backward Euler-finite element method reads: for each $n=1,\ldots,N$, find $U^{n}_{h}\in V^{n}_{hk}$, such that

[TABLE]

for every $X\in V^{n}_{hk}$, with $f^{n}:=f(t_{n})$ and $\mathcal{P}^{n}_{h}$ denoting the orthogonal $L_{2}$-projection operator onto $V^{n}_{hk}$.

Let now $\Delta^{n}_{h}:V^{n}_{hk}\rightarrow V^{n}_{hk}$ defined by $\left(-\Delta^{n}_{h}V,X\right)=\left(\nabla V,\nabla X\right)$, for all $V,X\in V^{n}_{hk}$, i.e., the discrete Laplacian. This allows for the strong representation of (6) as

[TABLE]

We now introduce a variant of the elliptic reconstruction [20, 19, 15], which will be instrumental in the proof of the a posteriori error bounds below.

Definition 3.1** (elliptic reconstruction).**

For each $n=0,1,\ldots,N,$ we define the elliptic reconstruction $\omega^{n}\in H^{1}_{0}(\Omega)$ to be the solution of the elliptic problem

[TABLE]

where

[TABLE]

Remark 3.2* (Galerkin orthogonality).*

We observe that $\omega^{n}$ satisfies

[TABLE]

This relation implies that $\omega^{n}-U^{n}_{h}$ is orthogonal to $V^{n}_{hk}$ with respect to the Dirichlet inner product, a crucial property that allows to use a posteriori error bounds for elliptic problems to estimate various norms of $\omega^{n}-U^{n}_{h}$ from above; we refer to Section 5 for a detailed discussion.

Definition 3.3** (time reconstruction).**

For $t\in J_{n}$, $n=1,\ldots,N$, we set

[TABLE]

where $\ell_{n}$ the piecewise linear Lagrange basis function with $\ell_{n}(t_{k})=\delta_{kn}$.

The above definition implies that the time derivative of $U_{h}$,

[TABLE]

is the discrete backward difference at $t_{n}$.

4. A posteriori error estimates

We shall now use the reconstructions defined above, together with non-standard energy and continuation arguments and a spectral estimate for the linearized steady-state problem about the approximate solution $U_{h}$, to arrive at a posteriori error bounds in the $L_{4}(L_{4})$-, $L_{2}(H^{1})$- and $L_{\infty}(L_{2})$-norms.

4.1. Error relation

We begin by splitting the total error as follows:

[TABLE]

In view of Remark 3.2, $\theta$ can be estimated by a posteriori error bounds for elliptic problems in various norms.

Also, $\rho$ satisfies an equation of the form (4) with a fully computable right-hand side that consists of $\theta$ and the problem data. To see this, (4) along with Definitions 3.1 and 3.3 and elementary manipulations lead to the following result.

Lemma 4.1** (error equation).**

On $J_{n}$, $n=1,\ldots,N$ and for all $v\in H_{0}^{1}(\Omega)$, we have

[TABLE]

Therefore, norms of $\rho$ can be estimated through PDE stability arguments; this will be performed below. Before doing so, however, we further estimate the term involving the elliptic reconstructions on the right-hand side from (14).

Lemma 4.2**.**

On $J_{n}$, $n=1,\dots,N$, we have

[TABLE]

for all $v\in H^{1}_{0}(\Omega)$.

Proof.

From (11) and Definition 3.1, we can write

[TABLE]

Then, using (9) in conjunction with (7), we obtain

[TABLE]

and correspondingly for $g_{h}^{n-1}$. Combining the above, the result already follows. ∎

4.2. Energy argument

We begin by introducing some notation. We define

[TABLE]

on each $J_{n}$, $n=1,\dots,N$, noting that $\mathcal{L}_{2}\equiv\mathcal{L}_{2}(t)$; for $n=1$ we adopt the convention that $U^{-1}_{h}=U^{0}_{h}$.

Moreover, for brevity, we also set

[TABLE]

[TABLE]

where $C_{0}:=(C_{PF}\tilde{c}^{2}+1)/2$, $C_{1}:=9+9C_{PF}\tilde{c}^{2}+6^{4}11^{2}C^{2}_{PF}\tilde{c}^{4}$, $C_{2}:=2\cdot 3^{7}C^{2}_{PF}\tilde{c}^{4}$, where $C_{PF}$ is the constant of the Poincaré-Friedrichs inequality $\|v\|\leq C_{PF}\|\nabla v\|$ and $\tilde{c}$ as in (2).

Lemma 4.3** ($d=2$).**

Let $d=2$ and $u$ be the solution of (4) and $\omega$ as in (11). Assume that $\rho(t)\in W^{1,4}_{0}(\Omega)$ for a.e. $t\in(0,T]$. Then, for any $\uptau\in(0,T]$, we have

[TABLE]

where

[TABLE]

Proof.

Using Taylor’s theorem, we immediately deduce

[TABLE]

Let $\phi:[0,\uptau]\times\Omega\rightarrow\mathbb{R}$, φορ $0<\uptau\leq T$, such that

[TABLE]

Hypothesis $\rho\in W^{1,4}_{0}(\Omega)$ implies that $\phi\in H_{0}^{1}(\Omega)$. Setting $v=\phi$ in (14), we have

[TABLE]

Observing now the identities

[TABLE]

elementary calculations yield

[TABLE]

We shall further estimate each $I_{j}$. We begin by splitting $I_{1}$ into

[TABLE]

Applying Hölder, GNL for $d=2$, Poincaré-Friedrichs and Young inequalities gives, respectively,

[TABLE]

The Cauchy-Schwarz and Young inequalities also yield $I^{2}_{1}\leq\frac{1}{2}\mathcal{L}_{2}+\frac{1}{2}\|\rho\|^{2}_{L_{2}(\Omega)}.$ Likewise, we split $I_{3}$ as follows:

[TABLE]

yielding the following bounds

[TABLE]

From Lemma 4.2 and working as before, we have

[TABLE]

[TABLE]

Next, we split $I_{2}$ as follows:

[TABLE]

and, using Hölder, Poincaré-Friedrichs and Young inequalities, we deduce

[TABLE]

Next, we split

[TABLE]

which can be further bounded as follows:

[TABLE]

In the same spirit, we also have

[TABLE]

and, thus,

[TABLE]

Next, we consider the splitting

[TABLE]

and we have the following bounds:

[TABLE]

Next, we set

[TABLE]

and we further estimate as follows:

[TABLE]

For $I_{6}$ and $I_{9}$, we work collectively as follows:

[TABLE]

and estimate:

[TABLE]

Finally for the last term on the right-hand side of (18), we have

[TABLE]

Applying the above estimates into (18) and integrating with respect to $t\in(0,\uptau)$ and observing the identities

[TABLE]

along with elementary manipulations, the result already follows. ∎

The use of the dimension-dependent GNL inequalities (2) necessitates certain modifications in the above argument when $d=3$, which we now provide. For brevity, we shall only provide the terms which are handled differently to the proof of the two-dimensional case from Lemma 4.3.

Lemma 4.4** ($d=3$).**

Let $d=3$, $u$ the solution of (4) and $\omega$ as in (11). Assume that $\rho(t)\in W^{1,4}_{0}(\Omega)$ for a.e. $t\in(0,T]$. Then, for any $\uptau\in(0,T]$, we have

[TABLE]

where

[TABLE]

with $\tilde{C}_{0}:=(C^{1/2}_{PF}\tilde{c}^{2}+1)/2$, $\tilde{C}_{1}:=9+9C^{1/2}_{PF}\tilde{c}^{2}+6^{4}11^{2}C_{PF}\tilde{c}^{4}$, $\tilde{C}_{2}:=3^{7}C_{PF}\tilde{c}^{4}$.

Proof.

Starting from (18), we discuss only the different treatment of the terms $I_{j}$, $j=6,9,11$; the estimation of the remaining terms is identical to the proof of Lemma 4.3 and is, therefore, omitted. To that end, we begin by setting $\zeta(\theta,U_{h}):=\|\theta\|_{L_{\infty}(\Omega)}+\|U_{h}\|_{L_{\infty}(\Omega)}$. Then, we have

[TABLE]

using (2) for $d=3$. Similarly, we have

[TABLE]

Likewise, using completely analogous arguments, we have

[TABLE]

The estimation of the remaining $I_{j}$ on the right-hand of (18) are completely analogous to the two-dimensional case with the difference that one applies (2) for $d=3$. Collecting all the estimates, we arrive at the desirable result. ∎

Remark 4.5*.*

In the a posteriori error estimation literature for evolution problems, $\mathcal{L}_{1}$ and the term $\epsilon^{-4}\|F(U_{h})-F(U^{n}_{h})\|^{2}_{L_{2}(\Omega)}$ of $\mathcal{L}_{2}$ are often referred to as the time error estimates, while $\|f^{n}-f\|^{2}_{L_{2}(\Omega)}$ is the data approximation. $\Theta_{1}$ represents the mesh change and $\Theta_{2}$ (or $\tilde{\Theta}_{2}$, respectively) is often termed as the spatial error estimate. These will be presented in detail in Section 5.

Remark 4.6*.*

We stress that the above result remains valid for the case of Neumann boundary conditions, upon modifying slightly the definition of the elliptic reconstruction (8) to eliminate the undetermined mode. Moreover, this can be done in such a way to recover (2) for terms involving $\rho$. This is not done here in the interest of simplicity of the presentation only.

4.3. Spectral estimates

To ensure polynomial dependence of the resulting estimates on $\epsilon^{-1}$, a widely used idea is to employ spectral estimates of the principal eigenvalue of the linearized Allen-Cahn operator:

[TABLE]

The celebrated works [8, 10, 2] showed that $\lambda$ can be bounded independently of $\epsilon$ for the case of smooth, evolved interfaces. This idea was used in the seminal works [14] for the proof of a priori and [18, 13] for a posteriori error bounds for finite element methods in various norms with constants depending upon $\epsilon^{-1}$ only in a polynomial fashion. The a priori nature of the spectral estimate (20) is somewhat at odds, however, with the presence of $\lambda$ in a posteriori error bounds. This difficulty was overcome in the seminal work [4] by first linearizing about the numerical solution $U_{h}$, viz.,

[TABLE]

and by then proving verifiable eigenvalue approximation error bounds. The latter ensure that it is possible to compute principle eigenvalue approximations $\Lambda_{h}>0$, such that $\Lambda_{h}\geq\lambda_{h}$; we refer to [4, Section 5] for the detailed construction. In short, it has been shown that for linear conforming finite element spaces, ($\kappa=1$,) it is possible to construct $\Lambda_{h}(t)\geq\lambda_{h}(t)$ for almost all $t\in(0,T]$ upon assuming that $\|U_{h}\|_{L_{\infty}(\Omega)}$ remains bounded independently of $\epsilon^{-1}$.

The $\epsilon$-independence $\lambda$, (resp. $\lambda_{h}$, $\Lambda_{h}$,) however, is not guaranteed when the evolving interfaces are subjected to topological changes. This is an important challenge, since phase-field approaches are preferred over sharp-interface models exactly due to their ability evolve interfaces past topological changes. To address this, in [7] (cf., also [5, 6]) a crucial observation on the temporal integrability of $\lambda$ under topological changes was given: during topological changes we have $\lambda\sim\epsilon^{-2}$, but only for time periods of length $\epsilon^{2}$. Therefore, it has been postulated that there exists an $m>0$, such that

[TABLE]

holds for some constant $C>0$ independent of $\epsilon$, for some $m\geq 0$; notice that for $m=0$, we return to the earlier case of no topological changes. A number of numerically validated scenarios justifying (22) for the scalar Allen-Cahn and its vectorial counterpart, the Ginzburg-Landau equation, can be found in [7]. Moreover, a construction for a $\Lambda_{h}\in L_{1}(0,T)$ such that

[TABLE]

has been provided in [7, Proposition 3.8].

The above motivate the following assumption on the behaviour of the principal eigenvalue $\lambda_{h}$, which we shall henceofrth adopt.

Assumption 4.7**.**

We postulate the validity of one of the following options:

• (I)

we assume that the zero level set $\Gamma_{t}=\{x\in\Omega:u(x,t)=0\}$ is sufficiently smooth. Then, for almost every $t\in(0,T]$, there exists a computable bound $\Lambda_{h}(t)\geq\lambda_{h}(t)$ which is independent of $\epsilon$.

• (II)

there exists an $m>0$, such that $\int_{0}^{T}\lambda_{h}(t)\mathop{}\!\mathrm{d}{t}\leq C+\log{(\epsilon^{-m})}$ for some constant $C>0$ independent of $\epsilon$ and we can construct a $\Lambda_{h}\in L_{1}(0,T)$ such that (23) holds.

Of course, Assumption 4.7(I) is a special case of Assumption 4.7(II), arising when $m=0$. Nonetheless, when Assumption 4.7(I) is valid, the resulting a posteriori error estimates will have more favourable dependence on the final time $T$ than the estimates that are possible under the more general Assumption 4.7(II).

We shall prove a posteriori error estimates under the more general Assumption 4.7(II), commenting, nevertheless, on the differences that would arise in the proof under 4.7(I) instead.

4.4. Continuation argument

We begin by noting that, compared to the state-of-the-art estimates of [7, 6], there are three additional terms on the right hand side of (16), (19), due to the use of the special test function (17): $\|\theta\|_{L_{4}(0,T;L_{4}(\Omega))}$ and $\|\theta_{t}\|_{L_{4}(0,T;L_{4}(\Omega))}$ which arise naturally and are symmetric with respect to the $\|.\|_{L_{4}(0,T;L_{4}(\Omega))}$ norm that is to be estimated, while the additional term $\|.\|_{L_{6}(0,T;L_{6}(\Omega))}$ can be compensated by the presence of the additional terms $A(t)$ (weighted norms) appearing on the left-hand side. Since the $L_{6}(0,T;L_{6}(\Omega))$-norm does not arise naturally in the Allen-Cahn energy functions, we have opted in dropping the $A(t)$ terms in the analysis below.

Assuming that $\Lambda_{h}$ is available, we set $v=\rho\in H_{0}^{1}(\Omega)$ in (21), to deduce

[TABLE]

For $d=2$, we work as follows. Upon setting

[TABLE]

$\mathcal{D}_{2}:=\max\{4,\alpha(U_{h})+2{\Lambda}_{h}(t)(1-\epsilon^{2})+2\},$ and $\mathcal{B}_{2}:=\max\{16\beta(\theta,U_{h}),\gamma(\theta,U_{h})\},$ we use (24) on the left-hand side of (16), we note that $-F^{\prime}(U_{h})\leq 1$, and ignore $\int_{0}^{\uptau}A(t)\mathop{}\!\mathrm{d}t$, to arrive at

[TABLE]

where $\bar{\mathcal{B}}_{2}:=\sup_{t\in[0,T]}\mathcal{B}_{2}(t)$.

Now, we set $E_{2}:=\exp\left(\int_{0}^{T}\mathcal{D}_{2}(t)\mathop{}\!\mathrm{d}{t}\right)$ and, for $d=2,3$, we use the abbreviation

[TABLE]

for the collection of semi-norms on the left-hand side of the last estimate. With this notation, we define the set

[TABLE]

The set $I_{2}$ is non-empty because $0\in I_{2}$ and the left-hand side depends continuously on $\uptau$. We set $\uptau^{*}:=\max I_{2}$, and we assume that $\uptau^{*}<T$; we aim to arrive at a contradiction. Hence, using the definition of the set $I_{2}$, we deduce

[TABLE]

If the last term on the right-hand side of the last estimate is bounded above by $\eta^{4}_{2}$, or, equivalently, if it holds

[TABLE]

then for all $0\leq\uptau\leq\uptau^{*}$ we have

[TABLE]

Since $\frac{1}{8}\|\int_{0}^{\uptau}\nabla\rho^{2}(s)\mathop{}\!\mathrm{d}s\|^{2}_{L_{2}(\Omega)}+\frac{1}{2}\|\rho\|^{2}_{L_{2}(\Omega)}\leq\mathcal{N}_{[0,\uptau],2}(\rho)$, Grönwall’s Lemma implies

[TABLE]

upon setting $\uptau=\uptau^{*}$. This contradicts the hypothesis $\uptau^{*}<T$ and, therefore, proves that $I_{2}=[0,T]$.

Likewise for $d=3$, we insert the spectral estimate (24) into (19), and we work as for $d=2$. Setting

[TABLE]

$\mathcal{D}_{3}:=\max\{4,\alpha(U_{h})+2{\Lambda}_{h}(t)(1-\epsilon^{2})+3\},$ $\mathcal{B}_{3}:=\max\{16\tilde{\beta}(\theta,U_{h}),\tilde{\gamma}(\theta,U_{h})\},$ and $\bar{\mathcal{B}}_{3}:=\sup_{t\in[0,T]}\mathcal{B}_{3}(t)$, through the same argumentation, we conclude that now the set $I_{3}:=\{\uptau\in[0,T]:\mathcal{N}_{[0,\uptau],3}(\rho)\leq 4\eta^{4}_{3}E_{3}\}$ equals $[0,T]$ upon assuming the condition

[TABLE]

The above argument has already confirmed the validity of the following result.

Lemma 4.8**.**

Assume that (25) holds when when $d=2$ or (26) holds when $d=3$. Then, we have the bound

[TABLE]

4.5. Main results

Now we are ready to present the main error estimate in the $L_{4}(0,T;L_{4}(\Omega))$-norm, from which we can easily arrived at a fully computable a posteriori estimate in Section 5.

Theorem 4.9**.**

Let $u_{0}\in L_{\infty}(\Omega)$ and $f\in L_{\infty}(0,T;L_{4}(\Omega))$, $\Omega\subset\mathbb{R}^{d}$, $d=2,3$. Let $u$ be the solution of (4) and $U_{h}$ is its approximation (6). Then, under Assumption (4.7)(II)* and the condition*

[TABLE]

the following error bound holds

[TABLE]

Proof.

Ignoring nonnegative terms on the left-hand side of (27), we have

[TABLE]

the proof follows by a triangle inequality. ∎

Remark 4.10*.*

Under the more restrictive Assumption 4.7(I), the continuation argument presented in Section 4.4 remains analogous with minor alterations. Specifically, if we set $m=0$ and we replace ${\Lambda}_{h}(t)$ by $\lambda$ and $E_{d}=\exp\left(\int_{0}^{T}\mathcal{D}_{d}(t)\mathop{}\!\mathrm{d}(t)\right)$ by $E_{d}=\exp\left(\bar{\mathcal{D}}_{d}T\right)$, with $\bar{\mathcal{D}}_{d}:=\sup_{t\in[0,T]}\max\{4,\alpha(U_{h})+2\lambda(1-\epsilon^{2})+d\},$ $d=2,3$, Theorem 4.9 remains valid.

Remark 4.11*.*

We stress that Theorem 4.9 holds also in cases whereby it is not possible to assume that $\|U_{h}\|_{L_{\infty}(0,T;L_{\infty}(\Omega))}$ is bounded independently of $\epsilon$. We note, however, that $\|U_{h}\|_{L_{\infty}(0,T;L_{\infty}(\Omega))}$ remains uniformly bounded with respect to $\epsilon$ and the mesh parameters in all scenarios of practical interest we are aware of and it is typically required in scenarios ensuring the validity of Assumption 4.7.

It is instructive to discuss in detail the dependence of the various terms appearing in (28) and (29) to assess the practicality of the resulting a posteriori error bound below. The computational challenge for $\epsilon\ll 1$ is manifested by the satisfaction of the condition (28). Indeed as $\epsilon\to 0$ the condition (28) becomes increasingly more stringent to be satisfied, necessitating meshes to be increasingly locally fine enough so as to reduce the estimator $\eta_{d}$; this results to proliferation of the numerical degrees of freedom. Once $\eta_{d}$ is small enough, an adaptive algorithm could make use of Theorem 4.9 for further estimation, which requires (28) to be valid.

Assume for argument’s sake that $\|U_{h}^{n}\|_{L_{\infty}(\Omega)}\leq C^{\prime}$ for all $n=1,\ldots,N$ for some $\epsilon$-independent constant $C^{\prime}>0$. Also, we have

[TABLE]

The $L_{\infty}(\Omega)$-norm of each $\theta^{n}$ will be further estimated in Section 5. For the moment, if also assume that $\|\theta^{n}\|_{L_{\infty}(\Omega)}\leq C^{\prime}$ uniformly with respect to $\epsilon$, then we can conclude that $6\leq\bar{\mathcal{B}}_{d}\leq CC^{\prime}$, $d=2,3$ and, therefore,

[TABLE]

for some generic constants $C>0$, independent of $\epsilon$, upon noting that $\sqrt[4]{6}>1.5$.

Moreover, in the case of smooth developed interfaces (Assumption 4.7(I)), one expects that $E_{d}\sim 1$ as highlighted in the classical works [8, 10]. When topological changes take place, we can follow [7] and postulate that $E_{d}\sim\epsilon^{-m}$, $m>0$. With the above convention, we find that (28) becomes

[TABLE]

for some constant $G_{d}\geq 1$ for all $m\geq 0$, thus encapsulating simultaneously both cases of Assumption 4.7.

Hence, the $\epsilon$-dependence for the condition (28) appears to be less stringent than in the respective conditional a posteriori in $L_{\infty}(L_{2})$- and $L_{2}(H^{1})$-norms from [4, 7, 6], which reads, roughly speaking, $\tilde{\eta}\leq c\epsilon^{4+3m}$ for the corresponding estimator $\tilde{\eta}$ and some constant $c>0$. Therefore, seeking to prove a posteriori error estimates for the $L_{4}(L_{4})$-norm error is, in our view, justified, as they can be potentially used to drive space-time adaptive algorithms without excessive numerical degree of freedom proliferation. This is an significant undertaking in its own right and will be considered in detail elsewhere.

The new a posteriori error analysis appears to also improve the $\epsilon$-dependence on the condition for $L_{2}(H^{1})$- and $L_{\infty}(L_{2})$-norm bounds compared to [13, 4, 7, 6] in certain cases. Of course, the different method of proof above results to different terms appearing in $\eta_{d}$ above compared to the respective conditional a posteriori error bounds from [13, 4, 7, 6]. Therefore, the performance of the proposed estimates above has to be assessed numerically before any conclusive statements can made. In particular, we have the following result.

Proposition 4.12** ($L_{2}(H^{1})$- and $L_{\infty}(L_{2})$-norm estimates).**

With the hypotheses of Theorem 4.9 and, assuming condition (28), we have the bounds

[TABLE]

Therefore, in the same setting as before, we have (28) implies

[TABLE]

If we accept that $\eta^{2}_{d}\sim\tilde{\eta}$ from [4, 7, 6], for the sake of the argument, at least at the level of the conditional estimate, (28) gives formally favourable dependence on $\epsilon$ when $d=2$ and $m\geq 0$ and also when $d=3$ and $m\geq 1/2$, compared to the respective dependence $\tilde{\eta}\leq c\epsilon^{4+3m}$ from [7, 6].

5. Fully computable upper bound

The bound in Theorem 4.9 is still not fully computable, due various terms involving $\theta$ and $\rho(0)$, which we shall now further estimate by computable quantities.

5.1. Initial condition estimates

For the terms involving $\rho(0)$, we have

[TABLE]

The Sobolev norms of $\theta$ appearing on $\eta_{d}$ can be further estimated by a posteriori bounds for elliptic problems; see, e.g., [22, 1]. We focus, therefore, in the derivation of $L_{p}$-norm a posteriori error bounds for elliptic problems for $\theta$ and for $\theta_{t}$ via suitable duality arguments. Although the derivation is somewhat standard, we prefer to present it here with some level of detail to highlight the regularity assumptions required. Specifically, consider the dual problem:

[TABLE]

on an $\Omega\subset\mathbb{R}^{d}$ convex domain. Then, there exists a constant $C_{\Omega}>0$, depending on the domain $\Omega$, such that

[TABLE]

we refer to [17] for details.

5.2. Spatial error estimates

We shall estimate $\Theta_{2}$ by residual-type estimators due to the presence of non-Hilbertian norms. In view of Remark 3.2 above, $\theta^{n}=w^{n}-U_{h}^{n}$ is the error of the elliptic problem (8), so we can further estimate norms of $\theta$ once we have estimators of the form

[TABLE]

at our disposal for $p=2,4,6$. Therefore, from (11) we have

[TABLE]

giving

[TABLE]

for $\hat{c}>0$ an algebraic constant.

Let $2\leq p<+\infty$. To determine the estimator $\mathcal{E}$ precisely, we set $\psi=\theta^{n}$ on (30) and we have

[TABLE]

from Remark 3.2, with $\mathcal{I}^{n}_{h}:W^{1,1}(\Omega),\rightarrow V^{n}_{hk}$ denoting the standard Scott-Zhang interpolation operator that satisfies optimal approximation properties [21]. Continuing in standard fashion, we have

[TABLE]

with $r_{n}:=\mathit{g}^{n}_{h}+\Delta U_{h}^{n}$, for some constant $C_{SZ}>0$ independent of $h_{n}$ and of the functions involves, using the approximation properties of $\mathcal{I}_{h}^{n}$; here $\llbracket\nabla U^{n}_{h}\rrbracket$ is the jump across the internal edge $\varepsilon$. Then, the elliptic regularity estimate (31) implies that

[TABLE]

the element residual at time $t_{n}$.

For the limiting case $p=+\infty$, we can take

[TABLE]

with $\ell_{h,d}=\left(\ln{(1/h_{n})}\right)^{\alpha_{d}}$, where $\alpha_{2}=2$ and $\alpha_{3}=1$; we refer to [11] for details.

5.3. Mesh change estimates

The general strategy of time extensions in (11), (12) consists in decomposing $\theta_{t}$ as follows

[TABLE]

with $U_{h}^{m}\in V^{m}$, $m=1,\ldots,N$. Since $V^{n}_{hk}\neq V^{n-1}_{hk}$in general, we define the Scott-Zang interpolation operator $\hat{\mathcal{I}}_{h}^{n}:H^{1}_{0}(\Omega)\rightarrow V^{n}_{hk}\cap V^{n-1}_{hk}$ relative to the finest common coarsening $\hat{\mathcal{T}}^{n}_{h}$ of $\mathcal{T}^{n}_{h}$ and $\mathcal{T}_{h}^{n-1}$. The latter allows to apply the Galerkin orthogonality property of the elliptic reconstruction in $V^{n}_{hk}\cap V^{n-1}_{hk}$. Moreover, we have the following approximation result: for all $\varepsilon\in\check{\mathcal{S}}_{h}^{n}\backslash\hat{\mathcal{S}}_{h}^{n}$ and $1\leq p<\infty$ it holds that

[TABLE]

where $\hat{h}_{n}:=\max(h_{n},h_{n-1})$, with $\omega(\varepsilon)$ denoting the neighbourhood of elements sharing the face $\varepsilon$, where, as before, the positive constant $C_{SZ}$ depends only on the shape regularity of the triangulation.

Setting $\psi=\theta_{t}$ on (30), we work as before to deduce

[TABLE]

with $\check{\mathcal{S}}_{h}^{n}$ denoting the finite element space subordinate to the coarsest common refinement $\check{\mathcal{T}}_{h}^{n}$ of $\mathcal{T}_{h}^{n}$ and $\mathcal{T}_{h}^{n-1}$ and $\partial X_{n}:=(X_{n}-X_{n-1})/k_{n}$ for some sequence $\{X_{n}\}_{n\in\mathbb{N}\cup\{0\}}$. Standard estimation via Hölder’s inequality and (34) give, in turn,

[TABLE]

Finally, the assumed elliptic regularity (31), gives the a posteriori error estimator

[TABLE]

for which we have $\|\theta_{t}\|^{p}_{L_{p}(\Omega)}\leq\hat{\mathcal{E}}^{p}(U_{h,t},g_{h,t};L_{p}(\Omega))$.

Aknowledgments

EHG acknowledges the financial support of The Leverhulme Trust via a research project grant (grant no. RPG-2015-306). DP acknowledges the financial support of the Stavros Niarchos Foundation within the framework of project ARCHERS.