Loading [MathJax]/jax/output/SVG/jax.js
Research article

An adaptive finite element method based on Superconvergent Cluster Recovery for the Cahn-Hilliard equation

  • Received: 22 October 2022 Revised: 12 December 2022 Accepted: 26 December 2022 Published: 09 January 2023
  • In this study, we construct an error estimate for a fully discrete finite element scheme that satisfies the criteria of unconditional energy stability, as suggested in [1]. Our theoretical findings, in more detail, demonstrate that this system has second-order accuracy in both space and time. Additionally, we offer a powerful space and time adaptable approach for solving the Cahn-Hilliard problem numerically based on the posterior error estimation. The major goal of this technique is to successfully lower the calculated cost by controlling the mesh size using a Superconvergent Cluster Recovery (SCR) approach in accordance with the error estimation. To demonstrate the effectiveness and stability of the suggested SCR-based algorithm, numerical results are provided.

    Citation: Wenyan Tian, Yaoyao Chen, Zhaoxia Meng, Hongen Jia. An adaptive finite element method based on Superconvergent Cluster Recovery for the Cahn-Hilliard equation[J]. Electronic Research Archive, 2023, 31(3): 1323-1343. doi: 10.3934/era.2023068

    Related Papers:

    [1] Hao Wu . A review on the Cahn–Hilliard equation: classical results and recent advances in dynamic boundary conditions. Electronic Research Archive, 2022, 30(8): 2788-2832. doi: 10.3934/era.2022143
    [2] Liupeng Wang, Yunqing Huang . Error estimates for second-order SAV finite element method to phase field crystal model. Electronic Research Archive, 2021, 29(1): 1735-1752. doi: 10.3934/era.2020089
    [3] Hao Wang, Wei Yang, Yunqing Huang . An adaptive edge finite element method for the Maxwell's equations in metamaterials. Electronic Research Archive, 2020, 28(2): 961-976. doi: 10.3934/era.2020051
    [4] Shan Jiang, Li Liang, Meiling Sun, Fang Su . Uniform high-order convergence of multiscale finite element computation on a graded recursion for singular perturbation. Electronic Research Archive, 2020, 28(2): 935-949. doi: 10.3934/era.2020049
    [5] Zengyan Zhang, Yuezheng Gong, Jia Zhao . A remark on the invariant energy quadratization (IEQ) method for preserving the original energy dissipation laws. Electronic Research Archive, 2022, 30(2): 701-714. doi: 10.3934/era.2022037
    [6] Suayip Toprakseven, Seza Dinibutun . A weak Galerkin finite element method for parabolic singularly perturbed convection-diffusion equations on layer-adapted meshes. Electronic Research Archive, 2024, 32(8): 5033-5066. doi: 10.3934/era.2024232
    [7] Changling Xu, Huilai Li . Two-grid methods of finite element approximation for parabolic integro-differential optimal control problems. Electronic Research Archive, 2023, 31(8): 4818-4842. doi: 10.3934/era.2023247
    [8] Zuliang Lu, Fei Huang, Xiankui Wu, Lin Li, Shang Liu . Convergence and quasi-optimality of $ L^2- $norms based an adaptive finite element method for nonlinear optimal control problems. Electronic Research Archive, 2020, 28(4): 1459-1486. doi: 10.3934/era.2020077
    [9] Victor Ginting . An adjoint-based a posteriori analysis of numerical approximation of Richards equation. Electronic Research Archive, 2021, 29(5): 3405-3427. doi: 10.3934/era.2021045
    [10] Hsueh-Chen Lee, Hyesuk Lee . An a posteriori error estimator based on least-squares finite element solutions for viscoelastic fluid flows. Electronic Research Archive, 2021, 29(4): 2755-2770. doi: 10.3934/era.2021012
  • In this study, we construct an error estimate for a fully discrete finite element scheme that satisfies the criteria of unconditional energy stability, as suggested in [1]. Our theoretical findings, in more detail, demonstrate that this system has second-order accuracy in both space and time. Additionally, we offer a powerful space and time adaptable approach for solving the Cahn-Hilliard problem numerically based on the posterior error estimation. The major goal of this technique is to successfully lower the calculated cost by controlling the mesh size using a Superconvergent Cluster Recovery (SCR) approach in accordance with the error estimation. To demonstrate the effectiveness and stability of the suggested SCR-based algorithm, numerical results are provided.



    The Cahn-Hilliard equation is solved in this study using the adaptive finite element approach.

    {ut=εΔ2u+1εΔf(u)in Ω×(0,T],nu=n(εΔu+1εf(u))=0,on Ω×(0,T],u(x,0)=u0in Ω. (1.1)

    In the equation above, we use to Ω signify a limited domain in Rd(d=2 or 3) with Lipschitz boundary Ω, and n stands for the unit outward normal to Ω, u also stands for the phase-field variable that will be solved. f(u) is a nonlinear function. ε is a constant. T is the ultimate time.

    The Cahn-Hilliard (CH) equation, often known as the common phase field model, was first published in 1958 in a key paper by Cahn and Hillard's [2], which examined the thermodynamic phenomenon of mutual diffusion of two substances (such as alloys and polymers). In general, this equation can be used to express the intricate phase separation and coarsening phenomena in solids, especially in materials science and fluid dynamics, for examples, see [3,4,5,6] and the references therein. The Cahn-Hilliard equation, which holds multiple time scales and spatial scales, is a rigid and nonlinear fourth-order PDE. So it is challenging to identify the precise solution.

    The majority of studies in recent years have focused on finite difference techniques or the Fourier-spectral approach (periodic boundary), and many authors have devoted their time to examining different Cahn-Hilliard equation versions, including the viscous Cahn-Hilliard equation, the surface Cahn-Hilliard equation, the CHNS equation, etc. and research issues from space and time that are related. Barrett and Elliott presented to study the Cahn-Hilliard system by using conforming and nonconforming finite element methods [7,8,9]. Shen and Xu described in [10] the treatment of the Cahn-Hilliard equation by the scalar auxiliary variable (SAV) method. Zhao and Xiao combined the surface finite element method (SFEM) with stabilized semi-implicit model to solve the surface Cahn-Hilliard model [11]. There are many parallels between the Allen-Cahn equation and the Cahn-Hilliard equation since they are phase field flows of the same energy in different regions, such as energy decreasing. In [12], Shen and Yang analyzed the discrete versions of the Allen-Cahn and Cahn-Hilliard equations. Numerous research has been conducted on the Allen-Cahn equation [13,14,15], and Chen and Huang applied the Superconvergent Cluster Recovery (SCR) approach for the Allen-Cahn equation [16]. The method was also used in [17] to solve the CHNS equation, but the Cahn-Hilliard equation has not yet seen this treatment. In [17], the Cahn-Hilliard component is treated as a first-order format combing with a time-space adaptive method. In our paper, we analyze the Cahn-Hilliard equation using a second-order time discrete method and convex splitting for the nonlinear factor f, combined with the SCR method.

    There are two reasons why we use adaptive methods in this study. First, due to the presence of the tiny parameter ε, it produces the phenomenon of the interfacial layer. When meshing on this interfacial layer, different from other parts, we need fine encryption, and the adaptive strategy can effectively solve this problem. Second, it takes a very long time for the numerical simulation of the Cahn-Hilliard model to reach the steady state, therefore using adaptive mesh creation makes sense in terms of time and money savings. At present, numerous writers have examined the Cahn-Hilliard equation's adaptive technique, including [11], where the discrete chemical potential is employed as an error estimator to evaluate the fluctuation in numerical energy. Based on the a posteriori error estimations, the standard mesh refinement procedures were used in [18] for automatic mesh refining. In [17], the SCR based posteriori error estimators were constructed for space discretization of phase field variable and velocity function. The adaptive time-stepping method was applied to change the time step in [19,20].

    In this work, we first provide an error estimate for the second-order fully discrete technique suggested in [1]. The fully discrete scheme is then introduced with an effective time-space adaptive technique called SCR. We need understand that the SCR approach has the same or higher accuracy while being more effective and affordable than other adaptive techniques. The SCR strategy was originally proposed by Huang and Yi [21], and they applied it to the Allen-Cahn system [16]. Moreover, for the definition of error estimator, the time discretization error estimator is given by the difference of numerical approximation between adjacent time steps, and the spatial discretization error estimator is defined by gradient, that is, the difference between the reconstructed gradient and the numerical gradient.

    The remainder of the essay is structured as follows. For the Cahn-Hilliard equation, we offer a second-order fully discrete scheme in Section 2, after which we examine the error estimation and unconditional energy stability. In Section 3, we introduce the characteristics of the superconvergent cluster recovery operator, and define the error estimation operator in time and space. Several numerical tests illustrate the effectiveness of our adaptive procedure in Section 4.

    First of all, we introduce a new variable μ=:εΔu+1εf(u) and rewrite Eq (1.1) as follows.

    {ut=Δμin Ω×(0,T],μ=εΔu+1εf(u)in Ω×(0,T],nu=nμ=0on Ω×(0,T],u(x,0)=u0in Ω, (2.1)

    where ε>0 is a small parameter, which represents the thickness of the transition interface between materials, and u means the phase field variable, u0 H1(Ω) as the initial value at t=0. The solution is driven into the two pure states u=±1 by the nonlinear part f(u)=F(u), and F(u) is a double-well potential function, which is given by F(u)=14(u21)2.

    Actually, the Cahn-Hilliard equation could be viewed as the H1gradient flow of the free energy E(u):

    E(u)=Ω(1εF(u)+ε2|u|2)dx, (2.2)

    Furthermore, by a simple mathematical derivation, for the Cahn-Hilliard equation, the following energy law can be derived:

    ddtE(u)=Ω|(εΔu+1εf(u))|2dx0, (2.3)

    that the energy-decreasing property is satisfied by the free energy E(u), over time,

    E(u)(tn+1)E(u)(tn),nN. (2.4)

    Besides, the Cahn-Hilliard equation satisfies the law of conservation of mass [2].

    ddtΩu(x,t)dx=0,0tT. (2.5)

    Next, we will present a few notations that will be used all across the paper before we get started. The standard Sobolev space will be indicated by the notation Hm=Wm,2(Ω), and its associated norm of Hm can be represented by m. As well as and (,) respectively stand in for norm and the inner product of L2. Besides, we use the symbol to signify the L norm.

    Let size h be the largest element diameter for every eh and Th={eh} be the set of regularly shaped triangles. Let Vh be the corresponding space in finite dimensions for piecewise linear continuous functions:

    Vh={vH1(Ω):v|ehP1(eh),ehTh},

    and its basis functions are the standard Lagrange basis functions ϕz(zNh, where the Nh represents the set of the triangular vertices). Naturally, we denote V0h=VhH10(Ω).

    In addition, for the real number p0 and satisfying the condition v(t), we define the following norm:

    vLp(0,T;X)=(T0v(t)pXdt)12.

    Then, by using the above symbol, we first introduce the weak form of Eq (2.1), which is as follow: find uL(0,T;H1(Ω)), wL2(0,T;H1(Ω))

    {(ut,w)+(μ,w)=0wH1(Ω),(μ,v)=ε(u,v)+1ε(f(u),v)vH1(Ω),(u(0)u0,q)=0qH1(Ω). (2.6)

    To solve Eq (2.6), the semi-discrete system needs to find uhC1(0,T;Vh)

    {(uh,t,wh)+(μh,wh)=0whVh,(μh,vh)=ε(uh,vh)+1ε(f(uh),vh)vhVh,(uh(0)u0h0,qh)=0qhVh, (2.7)

    where u0h is an approximation of u0 in Vh.

    Divide [0,T] into Nst subintervals In:=(tn1,tn], uniformly, 0=t0<t1<...<tN=T. Let t represent the time step size so T=Nt, n=1,...,N. The fully discrete scheme for the nonlinear approximation, which is developed from a modified Crank-Nicolson approach, is as follows: find unhVnh,n=1,2,...,N

    {(unhΠnun1ht,wh)+(μn12h,wh)=0whVnh,(μn12h,vh)=ε(^unh,vh)+1ε(f(unh,un1h),vh)vhVnh,(uh(0)u0h0,qh)=0qhV0h, (2.8)

    where ^unh, and f(unh,un1h) are defined as

    ^unh=unh+Πnun1h2,f(unh,un1h)=(unh)3+(unh)2Πnun1h+unh(Πnun1h)2+(Πnun1h)34unh+Πnun1h2, (2.9)

    as well as μn12 is the approximation at the midpoint tn12=(tn1+tn)/2 (directly computed), and Πn represents the interpolation into the finite element space Vh. Next, we consider the mass conservative, taking wh=1 in (2.8), the following results can be easily obtained

    Ωun=Ωun1=Ωu0.

    Obviously, the fully discrete scheme is mass conservative.

    In addition, We use the convex partitioning method to divide f into two parts, i.e., the linear convex part and the nonlinear concave part. The linear convex part is then handled by the modified Crank-Nicolson method, while the linear convex part is handled by the Midpoint approximation (MP) method. So far, the modified Crank-Nicolson method has been applied several times to the Cahn-Hilliard equation to deal with the nonlinear concave terms, such as [22,23,24,25], and the linear convex terms in these works of literature are treated by BDF2.

    Lemma 2.1. (Discrete Gronwall lemma) [26]. Let C and t be non-negative constants, and ak,bk,ck,dk be non-negative sequences satisfying

    ak+tnk=0bktn1k=0dkak+tn1k=0ck+C0,n1,

    then

    an+tnk=0bkexp(tn1k=0dk)(tn1k=0ck+C0),n1.

    Theorem 2.1. The solution of the fully discrete scheme Eqs (2.8) and (2.9), which is unconditionally energy stable, satisfies

    E(unh)E(un1h),n1.

    Proof. Taking wh=tμn12h and vh=unhun1h in Eq (2.8) we obtain

    {(unhun1h,μn12h)+tμn12h2=0,(μn12h,unhun1h)=ε2(unh2un1h2)+1ε(f(unh,un1h),unhun1h),

    where

    1ε(f(unh,un1h),unhun1h)=1ε((unh)3+(unh)2un1h+unh(un1h)2+(un1h)34unh+un1h2,unhun1h)=1ε14(unh)414(un1h)412(unh)2+12(un1h)2dx,

    by integrating the upper formula, we obtain

    0=tμn12h2+ε2(unh2un1h2)+1ε14(unh)414(un1h)412(unh)2+12(un1h)2dx.

    Also for energy we have

    E(unh)E(un1h)=1εF(unh)+ε2|unh|21εF(un1h)+ε2|un1h|2dx=1ε14(unh)414(un1h)412(unh)2+12(un1h)2dx+ε2(unh2un1h2),

    the equation is finally obtained

    tμn12h2+E(unh)E(un1h)=0,

    which means that the energy-decreasing feature is maintained by the provided fully discrete scheme Eqs (2.8) and (2.9). We omit Πn since it implies interpolation into the finite element space Vnh and has no bearing on the theoretical derivation procedure.

    To derive the error estimates, we define the elliptic projection operator as Ph:H1Vh, which satisfies

    ((Phvv),αh)=0,αhVh. (2.10)

    From the literature [27], we know the following two inequalities hold,

    vPhvCh2v2,(vPhv)tCh2vt2. (2.11)

    Then, we denote the error functions as,

    ˜en=Phu(tn)unh,ˆen=u(tn)Phu(tn),ˉen12=Phμ(tn12)μn12h,ˇen12=μ(tn12)Phμ(tn12).

    We may quickly derive by using the Taylor expansion with integral residuals and Young's inequality [12].

    Rn12st3tntn1uttt(t)2sdt,s=1,0, (2.12)
    Rn22st3tntn1utt(t)2s+2dt,s=1,0, (2.13)

    where

    Rn1=u(tn)u(tn1)tut(tn12),
    Rn2=(u(tn)+u(tn1)2u(tn12)).

    Theorem 2.2. Let unh and u(tn) represent, respectively, the solutions of Eqs (2.8), (2.9) and (2.1). After that, for uC(0,T;H2(Ω)), utL2(0,T;H2(Ω))L2(0,T;L4(Ω)), uttL2(0,T;H1(Ω)), and utttL2(0,T;H1(Ω)) we get

    ukhu(tk)+(tεkn=0μn12hμ(tn12)2)12C(ε,T)(K1(ε,u)t2+K2(ε,u)h2),

    where

    C(ε,T)exp(Tε),K1(ε,u)=ε(utttL2(0,T;H1)+uttL2(0,T;H1))+1ε(uttL2(0,T;L2)+utL2(0,T;L4)),K2(ε,u)=u02+εutL2(0,T;H2)+1εμC(0,T;H2)+uC(0,T;H2).

    Proof. Subtracting Eq (2.8) from the weak form of Eq (2.6) on Ω at tn12, we obtain

    1t(˜en˜en1,wh)+(ˉen12,wh)=(Rn11t(IPh)(u(tn)u(tn1)),wh),(ˉen12+ˇen12,vh)=ε2(˜en+˜en1,vh)+ε(Rn2,vh)+1ε(f(u(tn12))f(unh,un1h),vh),

    then taking wh=12εt(˜en+˜en1) and vh=tˉen12 respectively, and summing up the two identities above, we gain

    ˜en2˜en12+2tεˉen122=t(Rn1,˜en+˜en1)(IPh)(u(tn)u(tn1),˜en+˜en1)2t(Rn2,ˉen12)+2tε2(f(u(tn12))f(unh,unn1),ˉen12)2tε(ˇen12,ˉen12):=I+II+III+IV+V,

    where

    I:=t(Rn1,˜en+˜en1),II:=(IPh)(u(tn)u(tn1),˜en+˜en1),III:=2t(Rn2,ˉen12),IV:=2tε2(f(u(tn12))f(unh,un1h),ˉen12),V:=2tε(ˇen12,ˉen12).

    the Cauchy inequality with ε and Young's inequality were used to estimate I, II,III, and V as following:

    ItRn1˜en+˜en1t2(εRn12+1ε˜en+˜en12)εt42tntn1uttt(t)21dt+t2ε˜en2+t2ε˜en12,II(IPh)(u(tn)u(tn1))˜en+˜en1ε2tntn1(IPh)ut(t)2dt+12ε˜en2+12ε˜en12,III3tεRn22+t3εˉen1223εt4tntn1utt(t)21dt+t3εˉen122,V2tε(3ˇen12)(13ˉen12)3tεˇen122+t3εˉen122.

    For simplicity, we denote u(tn)=un and u(tn12)=un12. We consider the fourth term Ⅳ as follows:

    IV=2tε2(f(un12)f(unhun1h),ˉen12)=2tε2(unhun1h2un12+(un12)3(unh)3+(unh)2un1hh+unh(un1h)2+(un1h)34)=2tε2(unhun1h2un12+(un12)3(unun12)3+(unun12)3gn+gn(unh)3+(unh)2un1h+unh(un1h)2+(un1h)34,ˉen12):=2tε2(IV1+IV2+IV3+IV4,ˉen12),

    where

    gn=(un)3+(un)2un1+un(un1)2+(un1)34,IV1:=unh+un1h2,IV2:=(un12)3(un+un12)3,IV3:=(un+un12)3gn,IV4:=gn(unh)3+(unh)2un1h+unh(un1h)2+(un1h)34.

    We continue our analysis of these four terms using the Taylor expansion with integral remainders, Cauchy-Schwarz inequality, and Eqs (2.12) and (2.13), we arrive at:

    IV1=˜enˆen˜en1ˆen1+un+un12un12˜en+˜en12+ˆen+ˆen12+un+un12un12˜en+˜en12+ˆen+ˆen12+t3tntn1utt(t)2dt,
    IV2=(un12)3(un+un12)33ξ21(un+un12un12)3ξ21t3tntn1utt(t)2dt,
    IV3=(un)3(un)2un1un(un1)2+(un1)38=((un)2(un1)2)(unun1)8ξ24(unun1)2ξ24t3tntn1u2t(t)2dt,

    where ξ1 lies between (un+un1)/2 and un12, ξ2 lies between un and un1, ζ1 and ζ2 are located between un1 and un.

    IV4=1423(un)3+13(un+un1)3+23(un1)323(unh)3+13(unh+un1h)3+23(un1h)3=16((un)3(unh)3)+16((un1)3(un1h)3)+112((un+un1)3(unh+un1h)3)12ξ23(ununh)+12ξ24(un1un1h)+14ξ25(ununh+un1un1h)C(un1un1h+ununh)C(ˆen1+˜en1+ˆen+˜en)C(ˆen1+˜en1+ˆen+˜en),

    where ξ3 lies between un and un1, ξ4 lies between un1 and un1h, and ξ5 lies between un+un1 and unh+un1h. To guarantee the existence of ξ3,ξ4,ξ5, we need the condition [28] (Suppose C, a positive constant that is unaffected by t or h, exists).

    unhC,1nN.

    Further, combining the above inequalities IVi,i=1,2,3,4 into IV, we arrive at

    IV=2tε2(IV1+IV2+IV3+IV4,ˉen12)2tε2IV1+IV2+IV3+IV4ˉen122tε2(IV1+IV2+IV3+IV4)ˉen122tε2(˜en+˜en12+ˆen+ˆen12+t3tntn1utt(t)2dt+3ξ21t3tntn1utt(t)2dt+ξ24t3tntn1u2t(t)2dt+C(ˆen1+˜en1+ˆen+˜en))ˉen123tε(˜en+˜en12+ˆen+ˆen12+t3tntn1utt(t)2dt+3ξ21t3tntn1utt(t)2dt+ξ24t3tntn1u2t(t)2dt+C(ˆen1+˜en1+ˆen+˜en))2+t3εˉen122Ct4εtntn1utt(t)2dt+Ct4εtntn1u2t(t)2dt+C(ˆen12+˜en12+ˆen2+˜en2)+t3εˉen122.

    Moreover, when we combine the terms I,II,III,IV and V, we arrive at the following estimates,

    ˜en2˜en12+2tεˉen122εt42tntn1uttt(t)21dt+ε2tntn1(IPh)utt(t)2dt+3tεˇen122+t3εˉen122+Ct4εtntn1utt(t)2dt+Ct4εtntn1u2t(t)2dt+3εt4tntn1utt(t)21dt+t3εˉen122+C(ˆen12+˜en12+ˆen2+˜en2)+t3εˉen122.

    Using Eq (2.11), we may add up the aforementioned inequality for n=1,...,k(kT/t), and we find,

    ˜ek2˜e02+tεkn=1ˉen12212t4εuttt(t)2L2(0,T;H1)+ε2Ch4ut(t)2L2(0,T;H2)+3tεkn=1ˇen122+Ct4εutt(t)2L2(0,T;L2)+Ct4εut(t)L2(0,T;L4)+3εt4utt(t)2L2(0,T;H1)+Ckn=1(ˆen12+˜en12+ˆen2+˜en2).

    Lastly, by using the triangular inequality and the discrete Gronwall lemma to the above inequality, we may conclude the proof.

    The SCR technique takes derivatives to determine the recovered gradient at recovered locations after fitting a linear polynomial to solution values in a set of suitable sampling points around the vertex. In comparison to other methods, the SCR recovery technique implements adaptive algorithms more simply while saving on computing costs. The crucial component of the SCR technique is the introduction of a posterior error estimate operator. The process for getting the recovered gradient, which is broken down into three parts, is described in the sections that follow.

    Step1: We plan to restore the gradient of an interior vertex a=a0Nh in Ω (Nh represents the mesh nodes), try to choose some points symmetrically distributed around z, and then form a new point set A={ai=(xi,yi)}, 1in (n4). Specially, if select mesh nodes, we will get better results.

    Step2: The recovery gradient operator Gh:VhVh×Vh is defined as follows

    (Ghuh)(a)=pa(x,y).

    Let Kz represent the sample points of convex polygons, the linear polynomial pa(x,y) can be found

    pa(x,y):=argminp1P1ni=0(p1(ai)uh(ai))2,

    In order to overcome the instability caused by small parameter h, we denote F by

    F:(x,y)(ψ,φ)=(x,y)(x0,y0)h,

    where h:=max{|xix0|,|yiy0|}, i=1,2,...,n. Then the fitting polynomial can be written as

    pa(x,y)=PTm=ˆPTˆm,

    with

    PT=(1,x,y),ˆPT=(2,ψ,φ),mT=(m1,m2,m3),ˆmT=(ˆm1,ˆm2,ˆm3).

    This ˆm is the coefficient vector, it satisfies the linear systems

    ATAˆm=ATu,

    where

    A=(1ψ0φ01ψ1φ11ψnφn)andu=(u(a0)u(a1)u(an)).

    Finally, the recovered gradient can be derived

    Ghu=pa=(m2m3)=1h(ˆm2ˆm3)

    Step3: the recovered gradient Ghu on Ω is obtained by the interpolation

    Ghu=aNhGhu(a)ϕa.

    Next, the following are the properties of the operator Gh was proved in [21].

    a. For ehTh, there is a constant C which is independent of the value of h, satisfied

    Ghv0,ehC|v|1,KvVh,

    where K=3i=1Ki with Ki indicates the components containing the sample points coming from the ith vertex of eh.

    b. The recovering point a=a0=(x0,y0) is the center of the circle, and the sampling points ai=(xi,yi),i=1,2,...,n(n4) equally distributed around it. we derive

    |v(a)Ghv(a)|Ch2,vW3(eh).

    We use ηnt and ηh,eh as discretization error indicators in time and space respectively, and then adjust the time step and mesh adaption in the algorithm.

    ηnt:=unhΠnun1h,ηh,eh:=Ghuhuh0,eh,η2h=ehThη2h,eh.

    During SCR discretization, the time step changes with the error estimate, and the grid size changes according to the time step, so each tn constructs a new grid, and continues to construct a new space Vh on that grid and Πn represents the interpolation into the finite element space Vh.

    Remark 3.1. If the gradient is on the boundary Ω, we can treat it in some methods such as a higher order extension, extrapolation, and take the average value.

    Remark 3.2. The sample points we selected are placed with maximum symmetry around a. On the basis of this, we talk about the sample points in the next two situations. One is that the approximation order can be raised if they are also nodes. The other is that it will only improve recovery accuracy if they are not nodes.

    Remark 3.3. The SCR strategy can be used not only for the Crank-Nicolson approximation of time-discrete the Cahn-Hilliard equations, but also for other time-discrete formats of the Cahn-Hilliard equations, for example, we can apply to the BDF approximation mentioned in the literature [29][30], and its stability and unique solvability have been demonstrated. Unlike the format in this paper, the format has three levels of time, so we need unh and un1h to approximate un+1h.

    {(3un+1h4Πnunh+Πnun1h2t,vh)+(wn+1h,vh)=0,vhVh,(wn+1h,φh)=ε2(un+1h,φh)+((un+1h)32Πnunh+Πnun1h,φh)+At((un+1hΠnunh),φh),φhVh,

    To demonstrate the viability of our strategy, some numerical findings are offered in this section. We discuss a few characteristics of the Cahn-Hilliard Eqs (2.8) and (2.9), which are fully discrete and have various beginning value requirements.

    Example 4.1. In the first test, we use the initial condition listed below to investigate Eq (2.1).

    u0(x,y)=tanh(((x0.3)2+y20.22)/ε)tanh(((x+0.3)2+y20.22)/ε)×tanh((x2+(y0.3)20.22)/ε)tanh((x2+(y+0.3)20.22)/ε), (4.1)

    where the parameter Ω=[1,1]×[1,1],andε=0.01.

    Sequences of the mesh adaption results and the numerical solutions for six distinct time steps are shown in Figure 1. We can observe that when the node count decreases, the mesh refinement moves in accordance with the level-set zeros. Figure 2 shows the energy level dropping over time. Furthermore, it is also possible to perceive intuitively how the number of nodes and time steps change over time. We see that during the stage of energy's rapid evolution, a small time-step is chosen to capture the change in the numerical solution, and as the system settles into a stable state, the time steps increase.

    Figure 1.  Example 4.1 (CH), results of mesh adaptation and numerical solutions.
    Figure 2.  Example 4.1 (CH), from left to right: the distribution of energy, node, and time step.

    Example 4.2. In terms of the Cahn-Hilliard Eq (2.1). Allowing for the following initial condition, let Ω=[1,1]×[1,1].

    u0(x,y)=tanh((x2+y20.152)/ε)×tanh(((x0.31)2+y20.152)/ε)×tanh(((x+0.31)2+y20.152)/ε)×tanh((x2+(y0.31)20.152)/ε)×tanh((x2+(y+0.31)20.152)/ε)×tanh(((x0.31)2+(y0.31)20.152)/ε)×tanh(((x0.31)2+(y+0.31)20.152)/ε)×tanh(((x+0.31)2+(y0.31)20.152)/ε)×tanh(((x+0.31)2+(y+0.31)20.152)/ε), (4.2)

    with the parameter ε=0.01.

    In this illustration, Figure 3 shows a succession of mesh adaptation results and their numerical solutions. It is therefore evident via six different time steps that the mesh transformation adapts the rule of the zeros level-set and the number of nodes is decreasing. We infer from Figure 4 that the energy that decays over time, the number of nodes, and the number of time steps all change with time. This numerical result therefore shows that the time steps will grow greater than before when the system stabilizes.

    Figure 3.  Example 4.2 (CH), results of mesh adaptation and numerical solutions.
    Figure 4.  Example 4.2 (CH), from left to right: the distribution of energy, node, and time step.

    Example 4.3. In the third example, we apply the following initial value to the Cahn-Hilliard Eq (2.1)

    u0(x,y)=tanh(((x0.3)2+y20.252)/ε)tanh(((x+0.3)2+y20.32)/ε), (4.3)

    where the parameter Ω=[1,1]×[1,1],andε=0.01.

    Figure 5 shows the various adaptive meshes and the associated numerical solutions. It is evident that the mesh transformation applies the level-set rule for zeros. The discrete-time history in Figure 6 also shows that the energy degrades over time. Additionally, Figure 6 shows the connection between time and the quantity of nodes (time-step). In the early stages of the simulation, small steps are used to record changes in the numerical solution, and when the system stabilizes in the latter stages, bigger time steps are used.

    Figure 5.  Example 4.3 (CH), results of mesh adaptation and numerical solutions.
    Figure 6.  Example 4.3 (CH), from left to right: the distribution of energy, node, and time step.

    Example 4.4. In the last example, Eq (2.1) with the following initial value was analyzed

    u0(x,y)=0.01rand(x,y), (4.4)

    where the parameter Ω=[1,1]×[1,1],andε=0.01.

    In Figure 7, we plot the different adaptive meshes and their contour plots of the numerical solutions. Figure 8, it displays energy is getting more and more small over time and the distribution of the time steps as well as nodes. We come to the conclusion quite similar to the previous examples.

    Figure 7.  Example 4.4 (CH), results of mesh adaptation and numerical solutions.
    Figure 8.  Example 4.4 (CH), from left to right: the distribution of energy, node, and time step.

    Remark 4.1. When ε takes a value less than 0.01, the change of node number throughout the process still follows the law of zero level set, and the energy is still declining. But the time step will become smaller and the computational cost of the program will become larger, so finally, in the paper we choose ε=0.01 as the experimental result.

    Table 1 displays the whole discrete format's spatial and temporal precision Eqs (2.8) and (2.9). We take the initial value u0(x,y)=0.01rand(x,y), the parameter ε=0.01, and T=1. Since the precise solution to the initial equation is unknown, the numerical error ehu:=uh(x,T)uh2(x,T) will be determined by the difference between the coarse and fine grids, and the convergence order log2(ehu/eh2u) will be determined by the ratio of the errors. The table shows that the accuracy converges to 2, which is consistent with our theoretical findings for both spatial and temporal accuracy.

    Table 1.  The spatial convergence order and the time convergence order when ε=0.01.
    h umhumh20 Rate t Rate
    116 0.052470400 116
    132 0.018185700 1.528700 132 2.18632
    164 0.003449260 2.398450 164 2.52708
    1128 0.000945365 1.867340 1128 1.78946

     | Show Table
    DownLoad: CSV

    The fundamental idea of this paper is to obtain the spatial discretization operator by SCR method, and then use it as the index of spatial discretization of the Cahn-Hilliard equation. Additionally, the index of temporal discretization is determined by the variation in approximate solutions between adjacent time steps. First of all, in order to solve the Cahn-Hilliard equation, we first develop a second-order scheme that is unconditionally energy stable. Then, the fully discrete system has second-order accuracy both time and space, as determined by error estimation. Finally, these numerical findings are used to demonstrate the efficacy of this approach.

    Tian's research was supported by Shanxi Scholarship Council of China (No. 2021029) and the 2021 Shanxi Science and Technology Cooperation and Exchange Special Program (No. 202104041101019). Chen's reserach was partially supported by NSFC Project (12201010).

    No potential conflict of interest was reported by the authors.



    [1] F. Guillen-Gonzalez, G. Tierra, Second order schemes and time-step adaptivity for Allen-Cahn and Cahn-Hilliard models, Comput. Math. Appl., 68 (2014), 821–846. https://doi.org/10.1016/j.camwa.2014.07.014 doi: 10.1016/j.camwa.2014.07.014
    [2] J. Cahn, J. Hilliard, Free energy of a nonuniform system. Ⅰ. interfacial free energy, J. Chem. Phys., 28 (1958), 258–267. https://doi.org/10.1063/1.1744102 doi: 10.1063/1.1744102
    [3] A. Karma, W. Rappel, Quantitative phase-field modeling of dendritic growth in two and three dimensions, Phys. Rev. E, 57 (1998), 4323–4349. https://doi.org/10.1103/PhysRevE.57.4323 doi: 10.1103/PhysRevE.57.4323
    [4] S. Allen, J. Cahn, A microscopic theory for antiphase boundary motion and its application to antiphase domain coarsening, Acta. Matall., 27 (1979), 1085–1095. https://doi.org/10.1016/0001-6160(79)90196-2 doi: 10.1016/0001-6160(79)90196-2
    [5] R. Kobayashi, Modeling and numerical simulations of dendritic crystal growth, Physica D Nonlinear Phenom., 63 (1993), 410–423. https://doi.org/10.1016/0167-2789(93)90120-P doi: 10.1016/0167-2789(93)90120-P
    [6] M. Gurtin, D. Polignone, J. Vinals, Two-phase binary fluids and immiscible fluids described by an order parameter, Math. Models Methods Appl. Sci., 6 (1996), 815–831. https://doi.org/10.1142/S0218202596000341 doi: 10.1142/S0218202596000341
    [7] J. Barret, J. Blowey, H. Garcke, Finite element approximation of the Cahn-Hilliard equation with degenerate mobility, SIAM J. Numer. Anal., 37 (1999), 286–318. https://doi.org/10.1137/S0036142997331669 doi: 10.1137/S0036142997331669
    [8] C. Elliott, D. French, A nonconforming finite-element method for the two-dimensional Cahn-Hilliard equation, SIAM J. Numer. Anal., 26 (1989), 884–903. https://doi.org/10.1137/0726049 doi: 10.1137/0726049
    [9] C. Elliott, D. French, Numerical studies of the Cahn-Hilliard equation for phase separation, IMA J. Appl. Math., 38 (1987), 97–128. https://doi.org/10.1093/imamat/38.2.97 doi: 10.1093/imamat/38.2.97
    [10] J. Shen, J. Xu, J. Yang, The scalar auxiliary variable (SAV) approach for gradient flows, J. Comput. Phys., 353 (2018), 407–416. https://doi.org/10.1016/j.jcp.2017.10.021 doi: 10.1016/j.jcp.2017.10.021
    [11] S. Zhao, X. Xiao, X. Feng, An efficient time adaptivity based on chemical potential for surface Cahn-Hilliard equation using finite element approximation, Appl. Math. Comput., 369 (2020), 124901. https://doi.org/10.1016/j.amc.2019.124901 doi: 10.1016/j.amc.2019.124901
    [12] J. Shen, X. Yang, Numerical approximations of Allen-Cahn and Cahn-Hilliard equations, Discrete Contin. Dyn. Syst., 28 (2010), 1669–1691. https://doi.org/10.3934/dcds.2010.28.1669 doi: 10.3934/dcds.2010.28.1669
    [13] Y. Huang, W. Yang, H. Wang, J. Cui, Adaptive operator splitting finite element method for Allen-Cahn equation, Numer. Methods Partial Differ. Equations, 35 (2019), 1290–1300. https://doi.org/10.1002/num.22350 doi: 10.1002/num.22350
    [14] D. Kay, A. Tomasi, Color image segmentation by the Vector-Valued Allen-Cahn Phase-Field Model: A multigrid solution, IEEE. Trans. Image Process., 18 (2009), 2330–2339. https://doi.org/10.1109/TIP.2009.2026678 doi: 10.1109/TIP.2009.2026678
    [15] X. Feng, Y. Li, Analysis of interior penalty discontinuous Galerkin methods for the Allen-Cahn equation and the mean curvature flow, IMA J. Numer. Anal., 35 (2015), 1622–1651. https://doi.org/10.1093/imanum/dru058 doi: 10.1093/imanum/dru058
    [16] Y. Chen, Y. Huang, N. Yi, A SCR-based error estimation and adaptive finite element method for the Allen-Cahn equation, Comput. Math. Appl., 78 (2019), 204–223. https://doi.org/10.1016/j.camwa.2019.02.022 doi: 10.1016/j.camwa.2019.02.022
    [17] Y. Chen, Y. Huang, N. Yi, A decoupled energy stable adaptive finite element method for Cahn-Hilliard-Navier-Stokes equations, Commun. Comput. Phys., 29 (2021), 1186–1212. https://doi.org/10.4208/cicp.OA-2020-0032 doi: 10.4208/cicp.OA-2020-0032
    [18] D. Mao, L. Shen, A. Zhou, Adaptive finite element algorithms for eigenvalue problems based on local averaging type a posteriori error estimates, Adv. Comput. Math., 25 (2006), 135–160. https://doi.org/10.1007/s10444-004-7617-0 doi: 10.1007/s10444-004-7617-0
    [19] Z. Zhang, Z. Qiao, An adaptive time-stepping Strategy for the Cahn-Hilliard Equation, Commun. Comput. Phys., 11 (2012), 1261–1278. https://doi.org/10.4208/cicp.300810.140411s doi: 10.4208/cicp.300810.140411s
    [20] Y. Li, Y. Choi, J. Kim, Computationally efficient adaptive time step method for the Cahn-Hilliard equation, Comput. Math. Appl., 73 (2017), 1855–1864. https://doi.org/10.1016/j.camwa.2017.02.021 doi: 10.1016/j.camwa.2017.02.021
    [21] Y. Huang, N. Yi, The superconvergent cluster recovery method, J. Sci. Comput., 44 (2010), 301–322. https://doi.org/10.1007/s10915-010-9379-9 doi: 10.1007/s10915-010-9379-9
    [22] A. Diegel, C. Wang, S. Wise, Stability and convergence of a second order mixed finite element method for the Cahn-Hilliard equation, arXiv preprint, 2016, arXiv: 1411.5248. https://doi.org/10.48550/arXiv.1411.5248
    [23] J. Guo, C. Wang, S. Wise, X. Yue, An H2 convergence of a second-order convex-splitting, finite difference scheme for the three-dimensional Cahn–Hilliard equation, Commun. Math. Sci., 14 (2016), 489–515. https://dx.doi.org/10.4310/CMS.2016.v14.n2.a8 doi: 10.4310/CMS.2016.v14.n2.a8
    [24] K. Cheng, C. Wang, S. Wise, X. Yue, A second-order, weakly energy-stable pseudo-spectral scheme for the Cahn–Hilliard equation and its solution by the homogeneous linear iteration method, J. Sci. Comput., 69 (2016), 1083–1114. https://doi.org/10.1007/s10915-016-0228-3 doi: 10.1007/s10915-016-0228-3
    [25] J. Guo, C. Wang, S. Wise, X. Yue, An improved error analysis for a second-order numerical scheme for the Cahn-Hilliard equation, J. Comput. Appl. Math., 388 (2021), 113300. https://doi.org/10.1016/j.cam.2020.113300 doi: 10.1016/j.cam.2020.113300
    [26] J. Shen, On error estimates of the projection methods for the Navier-Stokes equations: first-order schemes, Math. Comput., 65 (1996), 1039–1065.
    [27] X. Feng, A. Prohl, Error analysis of a mixed finite element method for the Cahn-Hilliard equation, Numer. Math., 99 (2004), 47–84. https://doi.org/10.1007/s00211-004-0546-5 doi: 10.1007/s00211-004-0546-5
    [28] C. Li, Y. Huang, N. Yi, An unconditionally energy stable second order finite element method for solving the Allen–Cahn equation, J. Comput. Appl. Math., 353 (2019), 38–48. https://doi.org/10.1016/j.cam.2018.12.024 doi: 10.1016/j.cam.2018.12.024
    [29] Y. Yan, W. Chen, C. Wang, S. Wise, A second-order energy stable BDF numerical scheme for the Cahn-Hilliard equation, Commun. Comput. Phys., 23 (2018), 572–602.
    [30] K. Cheng, W. Feng, C. Wang, S. Wise, An energy stable fourth order finite difference scheme for the Cahn-Hilliard equation, J. Comput. Appl. Math., 362 (2019), 574–595. https://doi.org/10.1016/j.cam.2018.05.039 doi: 10.1016/j.cam.2018.05.039
  • This article has been cited by:

    1. Youngjin Hwang, Seokjun Ham, Chaeyoung Lee, Gyeonggyu Lee, Seungyoon Kang, Junseok Kim, A simple and efficient numerical method for the Allen–Cahn equation on effective symmetric triangular meshes, 2023, 31, 2688-1594, 4557, 10.3934/era.2023233
  • Reader Comments
  • © 2023 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)
通讯作者: 陈斌, bchen63@163.com
  • 1. 

    沈阳化工大学材料科学与工程学院 沈阳 110142

  1. 本站搜索
  2. 百度学术搜索
  3. 万方数据库搜索
  4. CNKI搜索

Metrics

Article views(1953) PDF downloads(170) Cited by(1)

Figures and Tables

Figures(8)  /  Tables(1)

Other Articles By Authors

/

DownLoad:  Full-Size Img  PowerPoint
Return
Return

Catalog