Finite element methods for the stochastic mean curvature flow

# Finite element approximations of the stochastic mean curvature flow of planar curves of graphs

Xiaobing Feng Department of Mathematics
The University of Tennessee
Knoxville, TN 37996, U.S.A.
Yukun Li Department of Mathematics
The University of Tennessee
Knoxville, TN 37996, U.S.A.
and  Andreas Prohl Mathematisches Institut, Universität Tübingen, Auf der Morgenstelle 10, D-72076, Tübingen, Germany
###### Abstract.

This paper develops and analyzes a semi-discrete and a fully discrete finite element method for a one-dimensional quasilinear parabolic stochastic partial differential equation (SPDE) which describes the stochastic mean curvature flow for planar curves of graphs. To circumvent the difficulty caused by the low spatial regularity of the SPDE solution, a regularization procedure is first proposed to approximate the SPDE, and an error estimate for the regularized problem is derived. A semi-discrete finite element method, and a space-time fully discrete method are then proposed to approximate the solution of the regularized SPDE problem. Strong convergence with rates are established for both, semi- and fully discrete methods. Computational experiments are provided to study the interplay of the geometric evolution and gradient type-noises.

###### Key words and phrases:
Stochastic mean curvature flow, level set method, finite element method, error analysis.
###### 1991 Mathematics Subject Classification:
65M12, 65M15, 65M60,
The work of the first two authors was partially supported by the NSF grant DMS-1016173.

## 1. Introduction

The mean curvature flow (MCF) refers to a one-parameter family of hypersurfaces which starts from a given initial surface and evolves according to the geometric law

 Vn(t,⋅)=H(t,⋅),

where and denote respectively the normal velocity and the mean curvature of the hypersurface at time . The MCF is the best known curvature-driven geometric flow which finds many applications in differential geometry, geometric measure theory, image processing and materials science and have been extensively studied both analytically and numerically (cf. [10, 16, 24, 29, 32] and the references therein).

As a geometric problem, the MCF can be described using different formulations. Among them, we mention the classical parametric formulation [18], Brakke’s varifold formulation [2], De Giorgi’s barrier function formulation [9, 3, 4], the variational formulation [1], the level set formulation [26, 14, 8], and the phase field formulation [13, 19]. We remark that different formulations often lead to different solution concepts and also lead to developing different analytical (and numerical) concepts and techniques to analyze and approximate the MCF. However, all these formulations of the MCF give rise to difficult but interesting nonlinear geometric partial differential equations (PDEs), and the resolution of the MCF then depends on the solutions of these nonlinear geometric PDEs. One interesting feature of the MCF is the development of singularities, in particular singularities which may occur in finite time, even when the initial hypersurface is smooth. The singularities may appear in different forms such as self-intersection, pinch-off, merging, and fattening. To understand and characterize these singularities have been the focus of the analytical and numerical research on the MCF (cf. [8, 10, 14, 15, 24, 29, 32], and the references therein).

For application problems, there is a great deal of interest to include stochastic effects, and to study the impact of special noises on regularities of solutions, as well as their long-time behaviors. The uncertainty may arise from various sources such as thermal fluctuation, impurities of the materials, and the intrinsic instabilities of the deterministic evolutions. In this paper we consider the following form of a stochastically perturbed mean curvature flow:

 (1.1) Vn=H(t,⋅)+ϵ˙Wt,

where denotes a white in time noise, and is a constant. It is easy to check that (cf. [31, 12]) the level set formulation of (1.1) is given by the following nonlinear parabolic stochastic partial differential equation (SPDE):

 (1.2) df=|∇x′f|divx′(∇x′f|∇x′f|)dt+ϵ|∇x′f|∘dWt,

where with denotes the level set function so that is represented by the zero level set of , and ‘’ refers to the Stratonovich interpretation of the stochastic integral. Again, stochastic effects are modeled by a standard -valued Wiener process which is defined on a given filtered probability space .

In the case that is a -dimensional graph, that is, , equation (1.2) reduces to

 (1.3) du=√1+|∇xu|2divx(∇xu√1+|∇xu|2)dt+ϵ√1+|∇xu|2∘dWt.

To the best of our knowledge, a comprehensive PDE theory for the SPDE (1.3) is still missing in the literature. For the case , (1.3) reduces to the following one-dimensional nonlinear parabolic SPDE:

 (1.4) du =∂2xu√1+|∂xu|2dt+ϵ√1+|∂xu|2∘dWt =∂x(arctan(∂xu))dt+ϵ√1+|∂xu|2∘dWt.

Here stands for the derivative of with respect to . This Stratonovich SPDE can be equivalently converted into the following Itô SPDE:

 (1.5) du =[ϵ22∂2xu+(1−ϵ22)∂2xu√1+|∂xu|2]dt+ϵ√1+|∂xu|2dWt =∂x(ε22∂xu+(1−ε22)arctan(∂xu))dt+ϵ√1+|∂xu|2dWt.

As is evident from (1.4), (1.5), the stochastic mean curvature flow (1.3) for may be interpreted as a gradient flow with multiplicative noise. Recently, Es-Sarhir and von Renesse [12] proved existence and uniqueness of (stochastically) strong solutions for (1.4) by a variational method, based on the Lyapunov structure of the problem (cf. [12, property (H3)]) which replaces the standard coercivity assumption (cf. [12, property (A)]). As is pointed out in [12], mild solutions for (1.4) may not be expected due to its quasilinear character.

The primary goal of this paper is to develop and analyze by a variational method some semi-discrete and fully discrete finite element methods for approximating (with rates) the strong solution of the Itô form (1.5) of the stochastic MCF. The error analysis presented in this paper differs from most existing works on the numerical analysis of SPDEs, where mild solutions are mostly approximated with the help of corresponding discrete semi-groups (see  [20] and the references therein). We also note that the error estimates derived in [17] which hold for general quasilinear SPDEs do not apply to (1.5) because the structural assumptions, such as the coercivity assumption [17, cf. Assumption 2.1, (ii)] and the strong monotonicity assumption [17, cf. Assumption 2.2, (i)] fail to hold for (1.5), and also the regularity assumptions [17, cf. Assumption 2.3] are not known to hold in the present case. In this paper, we use a variational approach similar to [17, 6, 7] to analyze the convergence of our finite element methods. One main difficulty for approximating the strong solution of (1.5) with certain rates is caused by the low regularity of the solution. To circumvent this difficulty, we first regularize the SPDE (1.5) by adding an additional linear diffusion term to the drift coefficient of (1.5); as a consequence the related drift operator in (2.3) becomes strongly monotone, and the corresponding solution process is then -valued in space. However, it is due to the ‘gradient-type’ noise that a relevant Hölder estimate in the -norm for the solution seems not available, which is necessary to properly control time-discretization errors. In order to circumvent this problematic issue, we proceed first with the spatial discretization (3.1)-(3.2); we may then use an inverse finite element estimate, and the weaker Hölder estimate (3.17) for the process to control time-discretization errors. We remark that addressing space discretization errors first requires to efficiently cope with the limited regularity of Lagrange finite element functions in the context of required higher norm estimates, which is overcome by a perturbation argument (cf. Proposition 3.2).

The remainder of this paper consists of three additional sections. In section 2 we first recall some relevant facts about the solution of (1.5) from [12]; we then present an analysis for the regularized problem. The main result of this section is to prove an error bound for in powers of . In section 3 we propose a semi-discrete (in space) and a fully discrete finite element method for the regularized equation (2.3) of the SPDE (1.5). The main result of this section is the strong -error estimate for the finite element solution. Finally, in section 4 we present several computational results to validate the theoretical error estimate, and to study relative effects due to geometric evolution and gradient-type noises.

## 2. Preliminaries and error estimates for a PDE regularization

The standard function and space notation will be adapted in this paper. For example, denotes the Sobolev space on the interval , and . Let denote the -inner product on . The triple stands for a given probability space. For a random variable , we denote by the expected value of .

We first quote the following existence and uniqueness result from [12] for the SPDE (1.5) with periodic boundary conditions.

###### Theorem 2.1.

Suppose that and fix . Let . There exists a unique strong solution to SPDE (1.4) with periodic boundary conditions and attaining the initial condition , that is, there exists a unique -valued -adapted process such that -almost surely

 (u(t),φ)I = (u0,φ)I−ϵ22∫t0(∂xu,∂xφ)Ids −(1−ϵ22)∫t0(arctan(∂xu),∂xφ)I]ds +ϵ∫t0(√1+|∂xu|2,φ)IdWs∀φ∈H1(I)∀t∈[0,T].

Moreover, satisfies for some independent of ,

 (2.2) supt∈[0,T]E[∥u(t)∥2H1(I)]≤C.

It is not clear if such a regularity can be improved from the analysis of [12] because of the difficulty caused by the gradient-type noise. In particular, -regularity in space, which would be desirable in order to derive some rates of convergence for finite element methods, seems not clear. To overcome this difficulty, we introduce the following simple regularization of (1.5):

 (2.3) duδ=[(δ+ϵ22)∂2xuδ+(1−ϵ22)∂2xuδ√1+|∂xuδ|2]dt+ϵ√1+|∂xuδ|2dWt.

To make this indirect approach successful, we need to address the well-posedness and regularity issues for (2.3) and to estimate the difference between the strong solutions of (2.3) and of (1.5).

###### Theorem 2.2.

Suppose that and , where is independent of . Let . Then there exists a unique strong solution to SPDE (2.3) with periodic boundary conditions and initial condition , that is, there exists a unique -valued -adapted process such that there holds -almost surely

 (2.4) (uδ(t),φ)I =(uδ0,φ)I−(δ+ϵ22)∫t0(∂xuδ,∂xφ)Ids +ϵ∫t0(√1+|∂xuδ|2,φ)IdWs∀φ∈H1(I)∀t∈[0,T].

Moreover, satisfies

 (2.5) supt∈[0,T]E[12∥∂xuδ(t)∥2L2(I)]+δE[∫T0∥∂2xuδ(s)∥2L2(I)ds]≤E[12∥∂xuδ0∥2L2(I)].
###### Proof.

Existence of can be shown in the same way as done in Theorem 2.1 (cf.  [12]). To verify (2.5), we proceed formally and apply Ito’s formula with to (a Galerkin approximation of) the solution to get

 =12∥∂xuδ0∥2L2(I)+ε22∫t0∥∥∂x√1+|∂xuδ|2∥∥2L2ds+Mt

where

 Mt:=ϵ∫t0(∂x√1+|∂xuδ(s)|2,∂xuδ)IdWs

is a martingale. Taking expectation yields

 E[12∥∂xuδ(t)∥2L2(I)+∫t0[δ∥∂2xuδ∥2L2 +(1−ε22)∥∥∂2xuδ√1+|∂xuδ|2∥∥2L2]ds] ≤E[12∥∂xuδ0∥2L2].

Hence, (2.5) hold. The proof is complete. ∎

Next, we shall derive an upper bound for the error as a low order power function of .

###### Theorem 2.3.

Suppose that . Let and denote respectively the strong solutions of the initial-boundary value problems (1.5) and (2.3) as stated in Theorems 2.1 and 2.2. Then there holds the following error estimate:

 (2.6) supt∈[0,T]E[∥uδ(t)−u(t)∥2L2(I)]+δE[∫T0∥∂x(uδ(s)−u(s))∥2L2(I)ds]≤CTδ.
###### Proof.

Let . Subtracting (2.1) from (2.4) we get that -a.s.

 (eδ(t),φ)I =−∫t0[δ(∂xu,∂xφ)I+(δ+ϵ22)(∂xeδ,∂xφ)I +(1−ϵ22)(arctan(∂xuδ)−arctan(∂xu),∂xφ)I]ds+Mt

for all and , with the martingale

 Mt:=ϵ∫t0(√1+|∂xuδ|2−√1+|∂xu|2,φ)IdWs.

By Itô’s formula we get

 (2.7) ∥eδ(t)∥2L2(I) =−2∫t0[δ(∂xu,∂xeδ)I+(δ+ϵ22)∥∂xeδ∥∥2L2(I) +(1−ϵ22)(arctan(∂xuδ)−arctan(∂xu),∂xeδ)I]ds +ϵ2∫t0∥∥√1+|∂xuδ|2−√1+|∂xu|2∥∥2L2(I)ds +2ϵ∫t0(√1+|∂xuδ|2−√1+|∂xu|2,eδ)IdWs.

Taking expectations on both sides, and using the monotonicity property of the function and the inequality yield

 E[∥eδ(t)∥2L2(I)] ≤δE[∫T0[∥∂xu∥2L2(I)+∥∂xeδ∥2L2(I)]ds],

which and (2.2) imply that

 E[∥eδ(t)∥2L2(I)]+δE[∫t0∥∂xeδ∥2L2(I)ds] ≤δE[∫T0∥∂xu∥2L2(I)ds] ≤(CT)δ.

The desired estimate (2.6) follows immediately. The proof is complete. ∎

## 3. Finite element methods

In this section we propose a fully discrete finite element method to solve the regularized SPDE (2.3) and to derive an error estimate for the finite element solution. This goal will be achieved in two steps. We first present and study a semi-discrete in space finite element method and then discretize it in time to obtain our fully discrete finite element method.

### 3.1. Semi-discretization in space

Let be a quasiuniform partition of . Define the finite element spaces

 Vhr:={vh∈C0(¯¯¯I);vh|[xj,xj+1]∈Pr([xj,xj+1]),j=0,1,⋯,J},

where denotes the space of all polynomials of degree not exceeding ) on . Our semi-discrete finite element method for SPDE (2.3) is defined by seeking such that -almost surely

 (3.1) (uδh(t),vh)I =(uδh(0),vh)I−(δ+ϵ22)∫t0(∂xuδh,∂xvh)Ids +ϵ∫t0(√1+|∂xuδh|2,vh)IdWs∀vh∈Vhr∀t∈[0,T], (3.2) uδh(t,0) =uδh(t,1)∀t∈[0,T],

where , and denotes the -projection operator from to .

To rewrite the above weak form in the equation form, we introduce the discrete (nonlinear) operator by

 (3.3) (Aδhwh,vh)I :=(δ+ϵ22)(∂xwh,∂xvh)I +(1−ϵ22)(arctan(∂xwh),∂xvh)I∀wh,vh∈Vhr.

Then (3.1) can be equivalently written as

 (3.4) duδh(t)=−Aδhuh(t)dt+ϵPh(√1+|∂xuδh(t)|2)dWt.
###### Proposition 3.1.

For , there is a unique solution to scheme (3.1). Moreover, there holds

 (3.5) sup0≤t≤TE[12∥∥uδh(t)∥∥2L2(I)]+δE[∫T0∥∥∂xuδh(s)∥∥2L2(I)ds]
###### Proof.

Well-posedness of (3.4) follows from the standard theory for stochastic ODEs with Lipschitz drift and diffusion. To verify (3.5), applying Itô’s formula to and using (3.4) we get

 (3.6) ∥uδh(t)∥2L2(I) =∥uδh(0)∥2L2(I)−2∫t0(Aδhuδh(s),uδh(s))Ids +ϵ2∫t0∥∥Prh√1+|∂xuδh(s)|2∥∥2L2(I)ds +2ϵ∫t0(Prh√1+|∂xuδh(s)|2,uδh)IdWs.

It follows from the definitions of and that

 (3.7) ∥uδh(t)∥2L2(I) ≤∥uδh(0)∥2L2(I)−(2δ+ϵ2)∫t0∥∥∂xuδh(s)∥∥2L2(I)ds −(2−ϵ2)∫t0(arctan(∂xuδh(s)),∂xuδh(s))Ids +ϵ2∫t0[1+∥∥∂xuδh(s)∥∥2L2(I)]ds +2ϵ∫t0(√1+|∂xuδh(s)|2,uδh)IdWs.

Then (3.5) follows from applying expectation to (3.7), and using the coercivity of arctan. The proof is complete. ∎

An a priori estimate for in stronger norms is more difficult to obtain, which is due to low global smoothness and local nature of finite element functions. We shall derive some of these estimates in Proposition 3.2 using a perturbation argument after establishing error estimates for .

To derive error estimates for , we introduce the elliptic -projection , i.e., for any , is defined by

 (3.8) (∂x[Rrhw−w],∂xvh)I+(Rrhw−w,vh)I=0∀vh∈Vhr.

The following error bounds are well-known (cf. [5]),

 (3.9) ∥∥w−Rrhw∥∥L2(I)+h∥∥w−Rrhw∥∥H1(I)≤Ch2∥w∥H2(I).
###### Theorem 3.1.

Let . Then there holds

 (3.10) supt∈[0,T]E[∥∥uδ(t)−uδh(t)∥∥2L2(I)] +δE[∫T0∥∥∂x[uδ(s)−uδh(s)]∥∥2L2(I)ds] ≤Ch2(1+δ−2).
###### Proof.

Let

 eδ(t):=uδ(t)−uδh(t),ηδ:=uδ(t)−Rrhuδ(t),ξδ:=Rrhuδ(t)−uδh(t).

Then . Subtracting (3.1) from (2.4) we obtain the following error equation which holds -almost surely:

 (3.11) (eδ(t),vh)I+(δ+ϵ22)∫t0(∂xeδ(s),∂xvh)Ids =−(1−ϵ22)∫t0(arctan(∂xuδ(s))−arctan(∂xuδh(s)),∂xvh)Ids +ϵ∫t0(√1+|∂xuδ(s)|2−√1+|∂xuδh(s)|2,vh)IdWs+(eδ(0),vh)I

for all . Substituting and rearranging terms leads to

 (3.12) (ξδ(t),vh)I+(δ+ϵ22)∫t0(∂xξδ(s),∂xvh)Ids =ϵ∫t0(√1+|∂xuδ(s)|2−√1+|∂xuδh(s)|2,vh)IdWs −(δ+ϵ22)∫t0(ηδ(s),vh)Ids−(ηδ(t),vh)I+(eδ(0),vh)I.

Applying Itô’s formula with , and using (3.12) and (3.3) we obtain

 (3.13) ∥∥ξδ(t)∥∥2L2(I)+(2δ+ϵ2)∫t0∥∥∂xξδ(s)∥∥2L2(I)ds =−(2−ϵ2)∫t0(arctan(∂xuδ(s))−arctan(∂xRrhuδ(s)),∂xξδ(s))Ids +ϵ2∫t0∥∥√1+|∂xuδ(s)|2−√1+|∂xuδh(s)|2∥∥2L2(I)ds +2ϵ∫t0(√1+|∂xuδ(s)|2−√1+|∂xuδh(s)|2,ξδ(s))IdWs

By the monotonicity of arctan, (3.9), (2.5), and the inequality we have

 E[∫t0(arctan(∂xRrhuδ(s))−arctan(∂xuδh(s)),∂xξδ(s))Ids]≥0, (2−ϵ2)E[∫t0(arctan(∂xuδ(s))−arctan(∂xRrhuδ(s)),∂xξδ(s))Ids] ≤E[∫t0(δ4∥∂xξδ(s)∥2L2(I)+4δ−1∥∥∂xuδ(s)−∂xRrhuδ(s)∥∥2L2(I))ds] ≤δ4E[∫t0∥∂xξδ(s)∥2L2(I)ds]+Ch2δ−2, E[ϵ2∫t0∥∥√1+|∂xuδ(s)|2−√1+|∂xuδh(s)|2∥∥2L2(I)ds] ≤E[(ϵ2+δ4)∫t0∥∂xξδ(s)∥2L2(I)ds]+Cδ−1E[∫t0∥∂xηδ∥2L2(I)ds] ≤(ϵ2+δ4)E[∫t0∥∂xξδ(s)∥2L2(I)ds]+Ch2δ−2, E[(eδ(0),ξδ(t))I]≤14E[∥ξδ(t)∥2L2(I)]+E[∥eδ(0)∥2L2(I)]≤14∥ξδ(t)∥2L2(I)+Ch2.

Taking the expectation in (3.13) and using the above estimates then yields

 (3.14) +3δE[∫T0∥∥∂xξδ(s)∥∥2L2(I)ds]≤Ch2(1+δ−2).

Finally, (3.10) follows from the triangle inequality, (3.9), and (3.14). The proof is complete. ∎

###### Remark 3.1.

(a) Estimate (3.10) is optimal in the -norm, but suboptimal in the -norm. The suboptimal rate for the -error is caused by the stochastic effect, i.e., the second term on the right-hand side of (3.13), and it is also caused by the lack of the space-time regularity in for .

(b) The proof still holds if the elliptic projection is replaced by the -projection .

We now use estimate (3.14) to derive some stronger norm estimates for . To this end, we define the discrete Laplacian by

 (3.15) (∂2hwh,vh)I=−(∂xwh,∂xvh)I∀wh,vh∈Vhr,

and the -projection by

 (Prhw,vh)I=(w,vh)I∀vh∈Vhr.
###### Proposition 3.2.

For there hold the following estimates for the solution of scheme (3.1):

 (3.16) sup0≤t≤TE[∥∂xuδh(t)∥2L2(I)]+δE[∫T0∥∂2huδh(s)∥2L2(I)ds]≤C(1+δ−2), (3.17) E[∥uδh(t)−uδh(s)∥2L2(I)+δ2∫ts∥∂x[uδh(ζ)−uδh(s)]∥2L2(I)dζ]

Notice that