Error Estimate for a Finite Difference Crank-Nicolson-ADI Scheme for a Class of Nolinear Parabolic Isotropic Systems

1. Introduction

Phase transition processes like solidification are often understood using phase field models. These models couple the energy (heat) equation for temperature with a nonlinear parabolic partial differential equation (p.d.e.) that includes another unknown, the phase, which takes values such as zero in the solid phase and one in the liquid phase. Phase field models provide a framework for describing microstructures and capturing the dynamics of interfaces and complex patterns during solidification. The phase field distinguishes between states of matter, with the phase variable set to zero in the solid phase and one in the liquid phase. This allows the simulation of phenomena like dendritic growth and intricate patterns at the solid-liquid interface.

Phase field models [6] are versatile and can be adapted to study various physical processes beyond solidification, including melting, precipitation, and alloy formation. They incorporate various physical parameters and boundary conditions, making them suitable for theoretical and practical applications in materials science. The coupling of the heat equation with the phase field equation enables the study of thermal gradients, cooling rates, and other factors on microstructural evolution during phase transitions, providing insights into the mechanisms driving these processes and aiding in material design.

These models have been widely used to describe phase transition phenomena such as solidification. In these models, the energy equation for temperature is coupled with a nonlinear parabolic p.d.e. involving the phase, which takes values like zero in the solid and one in the liquid phase, showing a large gradient near the interface. The solution from initial conditions should show a sharp moving front representing the interface and describe complex geometric patterns (dendrites) in solidification problems. Accurate numerical solutions require fine spatial discretizations around the interface and small time steps, which are time-consuming even in two dimensions.

In this note, we consider a specific system in two dimensions, simplifying a model by McFadden, Wheeler, Sekerka, Wang et al. [13]−[16]. This “isotropic” model consists of a system of two p.d.e’s of the form

$$\begin{array}{c}{\displaystyle l{\varphi }_t=\nabla ·\left(A\left(\varphi \right)\nabla \varphi \right)+f(\varphi ,u),}\hfill \\ {\displaystyle u_t=\Delta u+{\left[p\left(\varphi \right)\right]}_t,}\hfill \end{array}$$

(1)

where $\varphi =\varphi (x,y,t)$ is the phase indicator function and $u=u(x,y,t)$ is the temperature, both defined on a rectangle $\Omega$ in the $x,y$ plane for $t\ge 0$. The functions f and p are given smooth scalar functions of their arguments and A is a $2\times 2$ matrix given by

$$A\left(\varphi \right)=\left(\begin{array}{cc}\hfill a\left(\varphi \right)\hfill & 0\\ 0\hfill & b\left(\varphi \right)\end{array}\right),$$

where $a,b$ are smooth function of $\varphi .$ The system (1) is supplemented by given initial conditions $\varphi (x,y,0)={\varphi }_0(x,y)$, $u(x,y,0)=u_0(x,y)$, $(x,y)\in \Omega$, and boundary conditions of Neumann or Dirichlet type for $\varphi$ and u on the boundary $\partial \Omega$ of $\Omega$ for $t\ge 0$. In its fully nonlinear “anisotropic” version, where A is a function of $\nabla \varphi ,$ the system has been solved numerically by Wang [15], and Wang and Sekerka [16], using an ’explicit-implicit’ finite difference scheme. This scheme employs the explicit Euler method to advance the phase field over a temporal step via the first p.d.e., and then uses an ADI (Alternating Direction Implicit) scheme to solve for the temperature field in the second p.d.e. These references include numerous interesting numerical computations and measurements of the efficiency of the underlying numerical technique. In related work, Rappaz and his collaborators [2]–[4] considered similar systems and proved the existence and uniqueness of weak solutions. They also constructed and implemented fully discrete, adaptive finite element methods, using them to simulate dendritic growth in the anisotropic case. Furthermore, in [9]–[11], we constructed and implemented a Crank-Nicolson-ADI finite difference scheme for both equations in the general anisotropic case, demonstrating how it can be used in a parallel algorithm to efficiently solve the dendrite generation problem. In this note, we solve the isotropic system numerically using a finite difference discretization of the Crank-Nicolson-ADI type for both equations. We consider the initial-boundary-value problem with homogeneous Dirichlet boundary conditions for $\varphi$ and u on the boundary of $\Omega .$ (The Neumann problem may be similarly analyzed.) In Section 2, we state the p.d.e. problem and the assumptions on its coefficients. Section 3 presents preliminary results on various finite difference approximations. In Section 4, we state and prove the main result of the paper, establishing that the finite difference approximations to u and $\varphi$ converge in the discrete $L_2$ and $H_1$ norms, respectively, with bounds of order $\Delta t^2+h^2$, where $\Delta t$ is the time step and h the uniform spatial discretization mesh length, under a stability condition of the form ${\displaystyle \frac{\Delta t}{h}}\le \sigma$. We conclude by presenting simple numerical experiments that verify the orders of convergence. This comprehensive approach not only provides valuable insights into the numerical methods but also validates their efficiency and accuracy through practical implementation and testing. The results demonstrate the robustness of the Crank-Nicolson-ADI scheme in handling complex phase transition problems and pave the way for future research and applications in the field of materials science and computational physics.

2. The Problem

This “isotropic” model consists of a system (1) of two p.d.e’s is supplemented by given initial conditions $\varphi (x,y,0)={\varphi }_0(x,y)$, $u(x,y,0)=u_0(x,y)$, $(x,y)\in \Omega$, and boundary conditions of Neumann or Dirichlet type for $\varphi$ and u on the boundary $\partial \Omega$ of $\Omega$ for $t\ge 0$. Let $\Omega =(\alpha ,\beta )\times (\alpha ,\beta ).$ For $(x,y,t)\in \overline{\Omega }\times [0,T],$ where $T>0,$ we consider the initial-boundary-value problem

$$\left\{\begin{array}{cc}{\displaystyle ll\frac{\partial \varphi }{\partial t}-\frac{\partial }{\partial x}\left(a\left(\varphi \right)\frac{\partial \varphi }{\partial x}\right)-\frac{\partial }{\partial y}\left(b\left(\varphi \right)\frac{\partial \varphi }{\partial y}\right)=f(\varphi ,u),}\hfill & (x,y,t)\in \overline{\Omega }\times [0,T],\\ {\displaystyle \frac{\partial u}{\partial t}-\frac{{\partial }^{2}u}{\partial x^2}-\frac{{\partial }^{2}u}{\partial y^2}=\frac{\partial }{\partial t}p\left(\varphi \right),}\hfill & (x,y,t)\in \overline{\Omega }\times [0,T],\\ \varphi (x,y,t)=0,u(x,y,t)=0\hfill & (x,y,t)\in \partial \Omega \times [0,T],\\ \varphi (x,y,0)=g^{\varphi }(x,y)\hfill & (x,y)\in \overline{\Omega },\\ u(x,y,0)=g^u(x,y)\hfill & (x,y)\in \overline{\Omega },\end{array}\right.$$

(2)

where a and $b,$ are smooth functions of $\varphi ,$ such that $0<c_a^0\le a\left(\xi \right)\le c_a^1,$$0<c_b^0\le b\left(\xi \right)\le c_b^1,$ for $\xi \in \mathrm{I}\mathrm{R},$ and f and p are smooth functions defined on $\mathrm{I}{\mathrm{R}}^2$ and $\mathrm{I}\mathrm{R}$ respectively, with $f(0,0)=0$ and $p\left(0\right)=0.$ We suppose that $g^{\varphi }=g^u=0$ on $\partial \Omega .$ Furthermore, we make the assumption (in order to simplify the proofs of the error estimates) that the functions $f,a,b$, and p are globally Lipschitz functions of their arguments, i.e. that there exists a constant C such that:

$$|f(z_1,z_2)-f(Z_1,Z_2)|\le C(|z_1-Z_1|+|z_2-Z_2\left|\right),$$

(3)

$$|a\left(z_1\right)-a\left(Z_1\right)|\le C|z_1-Z_1|,$$

(4)

$$|b\left(z_1\right)-b\left(Z_1\right)|\le C|z_1-Z_1|,$$

(5)

$$|p\left(z_1\right)-p\left(Z_1\right)|\le C|z_1-Z_1|,$$

(6)

for all $z_1,z_2,Z_1,Z_2\in \mathrm{I}\mathrm{R}.$ We will assume that the initial-boundary-value problem (83) has a unique solution $(u,\varphi )$, smooth enough for the purpose of the numerical approximations.

3. Notation and Preliminary Results

We discretize the initial-boundary-value problem (83) as follows. Let $h=(\beta -\alpha )/(J+1),$ where J is a positive integer and $x_i=\alpha +ih,$$y_j=\alpha +jh,$$0\le i,j\le J+1$. We define ${\Omega }_h:=\left\{(x_i,y_j),i,j=1,\dots ,J\right\}$ and $\partial {\Omega }_h:=\left\{(x_i,y_j),i=0ori=J+1orj=0orj=J+1\right\}$. Let N be positive integer and $t^n=n\Delta t,$$n=0,\dots ,N,$ where $\Delta t=T/N,$ and define $t^{n+{\displaystyle \frac{1}{2}}}=t^n+\Delta t/2,$$x_{i\pm {\displaystyle \frac{1}{2}}}=x_i\pm h/2,$$y_{j\pm {\displaystyle \frac{1}{2}}}=y_j\pm h/2,$ and

${\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}:=\left\{U={(U_{00},\dots ,U_{J+1,J+1})}^T\in \mathrm{I}{\mathrm{R}}^{(J+2)\times (J+2)}:\right.$$\left.U_{ij}=0on\partial {\Omega }_h\right\}.$ We approximate the solution $u,\varphi$ of (83) by mesh functions $U^n,{\Phi }^n\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}$ as follows: For $n=0,1,2,\dots ,N,$ we approximate the vectors ${(u(x_0,y_0,t^n),\dots ,u(x_{J+1},y_{J+1},t^n))}^T,$${(\varphi (x_0,y_0,t^n),\dots ,\varphi (x_{J+1},y_{J+1},t^n))}^T$ by $U^n,{\Phi }^n\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}$ satisfying the following finite difference scheme:

$$\left.\begin{array}{cc}\left(\mathrm{i}\right){\Phi }_{ij}^0=g^{\varphi }(x_i,y_j),\hfill & (x_i,y_j)\in {\Omega }_h\bigcup \partial {\Omega }_h,\\ \left(\mathrm{ii}\right)U_{ij}^0=g^u(x_i,y_j),\hfill & (x_i,y_j)\in {\Omega }_h\bigcup \partial {\Omega }_h,\\ \mathrm{For} n=0,1,\dots ,N-1:\hfill & \\ {\displaystyle \left(\mathrm{iii}\right)\frac{{\Phi }_{ij}^{n+\frac{1}{2}}-{\Phi }_{ij}^n}{\frac{1}{2}\Delta t}-{\delta }_x\left(a_{i-\frac{1}{2},j}^n{\delta }_{\overline{x}}{\Phi }_{ij}^{n+\frac{1}{2}}\right)}\hfill & \\ {\displaystyle -{\delta }_y\left(b_{i,j-\frac{1}{2}}^n{\delta }_{\overline{y}}{\Phi }_{ij}^n\right)=F_{ij}^n,}\hfill & (x_i,y_j)\in {\Omega }_h,\\ {\displaystyle \left(\mathrm{iv}\right)\frac{{\Phi }_{ij}^{n+1}-{\Phi }_{ij}^{n+\frac{1}{2}}}{\frac{1}{2}\Delta t}-{\delta }_x\left(a_{i-\frac{1}{2},j}^n{\delta }_{\overline{x}}{\Phi }_{ij}^{n+\frac{1}{2}}\right)}\hfill & \\ {\displaystyle -{\delta }_y\left(b_{i,j-\frac{1}{2}}^n{\delta }_{\overline{y}}{\Phi }_{ij}^{n+1}\right)=F_{ij}^n,}\hfill & (x_i,y_j)\in {\Omega }_h,\\ \left(\mathrm{v}\right){\Phi }_{ij}^{n+{\displaystyle \frac{1}{2}}}={\Phi }_{ij}^{n+1}=0,\hfill & (x_i,y_j)\in \partial {\Omega }_h,\\ {\displaystyle \left(\mathrm{vi}\right)\frac{U_{ij}^{n+\frac{1}{2}}-U_{ij}^n}{\frac{1}{2}\Delta t}-{\delta }_{x\overline{x}}{U}_{ij}^{n+\frac{1}{2}}}\hfill & \\ {\displaystyle -{\delta }_{y\overline{y}}{U}_{ij}^n=\frac{P_{ij}^{n+1}-P_{ij}^n}{\Delta t},}\hfill & (x_i,y_j)\in {\Omega }_h,\\ {\displaystyle \left(\mathrm{vii}\right)\frac{U_{ij}^{n+1}-U_{ij}^{n+\frac{1}{2}}}{\frac{1}{2}\Delta t}-{\delta }_{x\overline{x}}{U}_{ij}^{n+\frac{1}{2}}}\hfill & \\ {\displaystyle -{\delta }_{y\overline{y}}{U}_{ij}^{n+1}=\frac{P_{ij}^{n+1}-P_{ij}^n}{\Delta t},}\hfill & (x_i,y_j)\in {\Omega }_h,\\ \left(\mathrm{viii}\right)U_{ij}^{n+{\displaystyle \frac{1}{2}}}=U_{ij}^{n+1}=0,\hfill & (x_i,y_j)\in \partial {\Omega }_h,\end{array}\right\}$$

(7)

where $P_{i,j}^n:=p\left({\Phi }_{i,j}^n\right),$ ${\delta }_x{v}_{ij}:={\displaystyle \frac{v_{i+1,j}-v_{i,j}}{h}},$ ${\delta }_{\overline{x}}{v}_{ij}:={\displaystyle \frac{v_{i,j}-v_{i-1,j}}{h}},$

${\delta }_y{v}_{ij}:={\displaystyle \frac{v_{i,j+1}-v_{i,j}}{h}},{\delta }_{\overline{y}}{v}_{ij}:={\displaystyle \frac{v_{i,j}-v_{i,j-1}}{h}},$${\delta }_{x\overline{x}}{v}_{ij}:={\displaystyle \frac{v_{i+1,j}-2v_{i,j}+v_{i-1,j}}{h^2}},$

${\delta }_{y\overline{y}}{v}_{ij}:={\displaystyle \frac{v_{i,j+1}-2v_{i,j}+v_{i,j-1}}{h^2}},$${\Delta }_h{v}_{ij}:={\delta }_{x\overline{x}}{v}_{ij}+{\delta }_{y\overline{y}}{v}_{ij}.$ We also let for $(x_i,y_j)\in {\Omega }_h$ $a_{i+{\displaystyle \frac{1}{2}},j}^n:=a\left({\displaystyle \frac{{\Phi }_{i+1,j}^n+{\Phi }_{i,j}^n}{2}}\right),$ $a_{i-{\displaystyle \frac{1}{2}},j}^n:=a\left({\displaystyle \frac{{\Phi }_{i,j}^n+{\Phi }_{i-1,j}^n}{2}}\right),$ $a_{i,j+{\displaystyle \frac{1}{2}}}^n:=a\left({\displaystyle \frac{{\Phi }_{i,j+1}^n+{\Phi }_{i,j}^n}{2}}\right),$ $a_{i,j-{\displaystyle \frac{1}{2}}}^n:=a\left({\displaystyle \frac{{\Phi }_{i,j}^n+{\Phi }_{i,j-1}^n}{2}}\right).$ (The quantities $b_{i\pm {\displaystyle \frac{1}{2}},j}^n,$ $b_{i,j\pm {\displaystyle \frac{1}{2}}}^n$ are defined analogously.) In addition, we let

$$F_{ij}^n:=\left\{\begin{array}{cc}f\mathrm{(}{\check{\Phi }}_{ij}^{{\displaystyle \frac{1}{2}}},{\check{U}}_{ij}^{{\displaystyle \frac{1}{2}}}\mathrm{)},\hfill & \mathrm{if}n=0,\\ f({\displaystyle \frac{3}{2}}{\Phi }_{ij}^n-{\displaystyle \frac{1}{2}}{\Phi }_{ij}^{n-1},{\displaystyle \frac{3}{2}}{U}_{ij}^n-{\displaystyle \frac{1}{2}}{U}_{ij}^{n-1}),\hfill & \mathrm{if}n\ge 1,\end{array}\right.$$

where

$${\check{\Phi }}_{ij}^{{\displaystyle \frac{1}{2}}}:=g^{\varphi }(x_i,y_j)+{\displaystyle \frac{\Delta t}{2}}\left\{a\left(g^{\varphi }(x_i,y_j)\right)g_{xx}^{\varphi }(x_i,y_j)+b\left(g^{\varphi }(x_i,y_j)\right)g_{yy}^{\varphi }(x_i,y_j)\right.$$

$$\left.+a^{\prime }\left(g^{\varphi }(x_i,y_j)\right){\left(g_x^{\varphi }(x_i,y_j)\right)}^2+b^{\prime }\left(g^{\varphi }(x_i,y_j)\right){\left(g_y^{\varphi }(x_i,y_j)\right)}^2+f(g^{\varphi }(x_i,y_j),g^u(x_i,y_j)\right\}$$

and

$${\check{U}}_{ij}^{{\displaystyle \frac{1}{2}}}:=g^u(x_i,y_j)+{\displaystyle \frac{\Delta t}{2}}\left\{\Delta g^u(x_i,y_j)+{\displaystyle \frac{2}{\Delta t}}{p}^{\prime }\left(g^{\varphi }(x_i,y_j)\right)[{\check{\Phi }}_{ij}^{{\displaystyle \frac{1}{2}}}-g^{\varphi }(x_i,y_j)]\right\}.$$

Subtracting (48.(vi)) from (48.(vii)) we have for $(x_i,y_j)\in {\Omega }_h$

$${\displaystyle \frac{U_{ij}^{n+1}-2U_{ij}^{n+\frac{1}{2}}+U_{ij}^n}{\frac{1}{2}\Delta t}}-{\delta }_{y\overline{y}}(U_{ij}^{n+1}-U_{ij}^n)=0.$$

Solving for $U_{ij}^{n+{\displaystyle \frac{1}{2}}}$ we get

$$U_{ij}^{n+{\displaystyle \frac{1}{2}}}={\displaystyle \frac{U_{ij}^{n+1}+U_{ij}^n}{2}}-{\displaystyle \frac{\Delta t}{4}}{\delta }_{y\overline{y}}(U_{ij}^{n+1}-U_{ij}^n).$$

(8)

Replacing the above relation in (48.(vi)) we get for $0\le n\le N-1$

$${\displaystyle \frac{U_{ij}^{n+1}-U_{ij}^n}{\Delta t}}-{\Delta }_h{\displaystyle \frac{U_{ij}^{n+1}+U_{ij}^n}{2}}+{\displaystyle \frac{\Delta t}{4}}{\delta }_{x\overline{x}}{\delta }_{y\overline{y}}(U_{ij}^{n+1}-U_{ij}^n)={\displaystyle \frac{P_{ij}^{n+1}-P_{ij}^n}{\Delta t}},(x_i,y_j)\in {\Omega }_h.$$

(9)

Conversely, suppose that (50) holds. Then, defining $U^{n+{\displaystyle \frac{1}{2}}}$from (49) we see that (48.(vi)) and (48.(vii)) are valid. Hence the relations (48.(vi)-(vii)) are equivalent to (49)-(50). Similarly, subtracting (48.(iii)) from (48.(iv)), we see, for $(x_i,y_j)\in {\Omega }_h$

$${\displaystyle \frac{{\Phi }_{ij}^{n+1}-2{\Phi }_{ij}^{n+\frac{1}{2}}+{\Phi }_{ij}^n}{\frac{1}{2}\Delta t}}-{\delta }_y\left(b_{i,j-{\displaystyle \frac{1}{2}}}^n{\delta }_{\overline{y}}({\Phi }_{ij}^{n+1}-{\Phi }_{ij}^n)\right)=0,$$

and therefore

$${\Phi }_{ij}^{n+{\displaystyle \frac{1}{2}}}={\displaystyle \frac{{\Phi }_{ij}^{n+1}+{\Phi }_{ij}^n}{2}}-{\displaystyle \frac{\Delta t}{4}}{\delta }_y\left(b_{i,j-{\displaystyle \frac{1}{2}}}^n{\delta }_{\overline{y}}({\Phi }_{ij}^{n+1}-{\Phi }_{ij}^n)\right).$$

(10)

Finally, substituting in (48.(iii)) we get for $0\le n\le N-1,(x_i,y_j)\in {\Omega }_h:$

$${\displaystyle \frac{{\Phi }_{ij}^{n+1}-{\Phi }_{ij}^n}{\Delta t}}-{\delta }_x\left(a_{i-{\displaystyle \frac{1}{2}},j}^n{\delta }_{\overline{x}}\left({\displaystyle \frac{{\Phi }_{ij}^{n+1}+{\Phi }_{ij}^n}{2}}\right)\right)-{\delta }_y\left(b_{i,j-{\displaystyle \frac{1}{2}}}^n{\delta }_{\overline{y}}\left({\displaystyle \frac{{\Phi }_{ij}^{n+1}+{\Phi }_{ij}^n}{2}}\right)\right)$$

$$+{\displaystyle \frac{\Delta t}{4}}{\delta }_x\left(a_{i-{\displaystyle \frac{1}{2}},j}^n\left({\delta }_{\overline{x}}\left({\delta }_y\left(b_{i,j-{\displaystyle \frac{1}{2}}}^n\left({\delta }_{\overline{y}}({\Phi }_{ij}^{n+1}-{\Phi }_{ij}^n)\right)\right)\right)\right)\right)=F_{ij}^n.$$

(11)

Conversely, let (52) hold. Then, if we define ${\Phi }^{n+{\displaystyle \frac{1}{2}}}$ from (51), we see that (48.(iii)) and (48.(iv)) are equivalent to (52) and (51).

While the solution algorithm is given by the equations (48), the error analysis will be based on (50) and (52). It is clear that solving in (48) for ${\Phi }^{n+{\displaystyle \frac{1}{2}}},$${\Phi }^{n+1},$$U^{n+{\displaystyle \frac{1}{2}}},$$U^{n+1}$ requires the solution of $J,$$J\times J$ symmetric tridiagonal systems of equations in each case. There systems are clearly invertible in view of our assumptions on the functions a and $b.$

In what follows we list a few auxiliary helpful results whose proofs are standard and may be found in [10].

Lemma 1.

a) 

Let $\varphi \in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}},\psi \in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}.$ Then

$$\sum _{i=1}^J({\varphi }_{i+1,j}-{\varphi }_{i,j}){\psi }_{ij}=-\sum _{i=1}^{J+1}{\varphi }_{ij}({\psi }_{i,j}-{\psi }_{i-1,j}).$$

b) 

Let $\varphi \in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}},\psi \in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}.$ Then

$$\sum _{j=1}^J({\varphi }_{i,j+1}-{\varphi }_{i,j}){\psi }_{ij}=-\sum _{j=1}^{J+1}{\varphi }_{ij}({\psi }_{i,j}-{\psi }_{i,j-1}).$$

c) 

If $\varphi ,v\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}},$ then for $(x_i,y_j)\in {\Omega }_h:$

i) 

${\delta }_{x\overline{x}}{v}_{ij}:={\delta }_x{\delta }_{\overline{x}}{v}_{ij}={\delta }_{\overline{x}}{\delta }_x{v}_{ij},{\delta }_{y\overline{y}}{v}_{ij}:={\delta }_y{\delta }_{\overline{y}}{v}_{ij}={\delta }_{\overline{y}}{\delta }_y{v}_{ij}.$

ii) 

${\delta }_x{\delta }_y{v}_{ij}={\delta }_y{\delta }_x{v}_{ij},$${\delta }_{\overline{x}}{\delta }_y{v}_{ij}={\delta }_y{\delta }_{\overline{x}}{v}_{ij},$${\delta }_{\overline{y}}{\delta }_x{v}_{ij}={\delta }_x{\delta }_{\overline{y}}{v}_{ij}.$

iii) 

${\delta }_{\overline{x}\overline{y}}{v}_{ij}:={\delta }_{\overline{x}}{\delta }_{\overline{y}}{v}_{ij}={\delta }_{\overline{y}}{\delta }_{\overline{x}}{v}_{ij}.$

iv) 

${\delta }_{\overline{x}}\left({\varphi }_{ij}\left({\delta }_{\overline{y}}{v}_{ij}\right)\right)=\left[{\delta }_{\overline{x}}{\varphi }_{ij}{\delta }_{\overline{y}}{v}_{ij}+{\varphi }_{ij}{\delta }_{\overline{x}\overline{y}}{v}_{ij}\right]-h{\delta }_{\overline{x}}{\varphi }_{i,j}{\delta }_{\overline{x}\overline{y}}{v}_{i,j},$

${\delta }_{\overline{y}}\left({\varphi }_{ij}\left({\delta }_{\overline{x}}{v}_{ij}\right)\right)=\left[{\delta }_{\overline{y}}{\varphi }_{ij}{\delta }_{\overline{x}}{v}_{ij}+{\varphi }_{ij}{\delta }_{\overline{x}\overline{y}}{v}_{ij}\right]-h{\delta }_{\overline{y}}{\varphi }_{i,j}{\delta }_{\overline{x}\overline{y}}{v}_{i,j}.$

v) 

${\delta }_x\left({\varphi }_{ij}{v}_{ij}\right)=v_{ij}{\delta }_x{\varphi }_{ij}+{\varphi }_{ij}{\delta }_x{v}_{ij}+h\left({\delta }_x{\varphi }_{i,j}\right)\left({\delta }_x{v}_{i,j}\right).$

For $v,w\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}},$ define the ${\ell }_2$ inner product on ${\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}$ by ${\displaystyle {(v,w)}_h:=h^2\sum _{i=1}^J\sum _{j=1}^J{v}_{ij}{w}_{ij},}$ with corresponding norm ${\displaystyle {\mathrm{\parallel }v\mathrm{\parallel }}_h:=h{\left(\sum _{i=1}^J\sum _{j=1}^J{v}_{ij}^2\right)}^{\frac{1}{2}}.}$

Lemma 2.

Define ${\Delta }_h:{\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}\to {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}$ by

$${\left({\Delta }_{h}v\right)}_{i,j}:=\left\{\begin{array}{cc}{\displaystyle ll\frac{v_{i+1,j}+v_{i-1,j}+v_{i,j+1}+v_{i,j-1}-4v_{ij}}{h^2},}\hfill & (x_i,y_j)\in {\Omega }_h,\\ 0,\hfill & (x_i,y_j)\in \partial {\Omega }_h\end{array}\right.$$

Then, there hold that

$${({\Delta }_{h}v,w)}_h={(v,{\Delta }_{h}w)}_h\forall v,w\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}},$$

(12)

and

$${({\Delta }_{h}v,v)}_h=-{\left|v\right|}_{1,h}^2\forall v\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}},$$

(13)

where ${\left|v\right|}_{1,h}^2:=h^2{\sum }_{i=1}^{J+1}{\sum }_{j=1}^{J+1}[{\left({\delta }_{\overline{x}}{v}_{ij}\right)}^2+{\left({\delta }_{\overline{y}}{v}_{ij}\right)}^2]$ is a discrete $H_1$ norm.

Lemma 3.

Define $L_h^n:{\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}\to {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}$ by

$${\left(L_h^{n}v\right)}_{i,j}:=\left\{\begin{array}{cc}{\displaystyle ll{\delta }_x\left(a_{i-\frac{1}{2},j}^n{\delta }_{\overline{x}}{v}_{ij}\right)+{\delta }_y\left(b_{i,j-\frac{1}{2}}^n{\delta }_{\overline{y}}{v}_{ij}\right)}\hfill & (x_i,y_j)\in {\Omega }_h,\\ 0\hfill & (x_i,y_j)\in \partial {\Omega }_h.\end{array}\right.$$

Then, there holds that

$${(L_h^{n}v,w)}_h={(v,L_h^{n}w)}_h\forall v,w\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}$$

(14)

and

$${(L_h^{n}v,v)}_h=-h^2\sum _{j=1}^{J+1}\sum _{i=1}^{J+1}\left\{a_{i+{\displaystyle \frac{1}{2}},j}^n{\left({\delta }_{\overline{x}}{v}_{ij}\right)}^2+b_{i,j+{\displaystyle \frac{1}{2}}}^n{\left({\delta }_{\overline{y}}{v}_{ij}\right)}^2\right\}\le 0\forall v\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}.$$

(15)

Lemma 4.

Define, for $v,w\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}},$ the bilinear form ${\beta }_h(v,w)$ as

$${\displaystyle {\beta }_h^n(v,w):={\left({\delta }_x\left(a^n\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left({\delta }_{\overline{y}}v\right)\right)\right)\right)\right),w\right)}_h:=}$$

$${\displaystyle h^2\sum _{i=1}^J\sum _{j=1}^J{\delta }_x\left(a_{i-\frac{1}{2},j}^n\left({\delta }_{\overline{x}}\left({\delta }_y\left(b_{i,j-\frac{1}{2}}^n\left({\delta }_{\overline{y}}{v}_{ij}\right)\right)\right)\right)\right)w_{ij}.}$$

Then, for some constants $\lambda ,\mu$ independent of $h,$ there holds that

$${\beta }_h^n(v,v)\ge {\lambda \left|v\right|}_{\overline{x}\overline{y},h}^2-\mu {\left|v\right|}_{1,h}^2\forall v\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}},$$

(16)

where ${\displaystyle {\left|v\right|}_{\overline{x}\overline{y},h}^2:=h^2\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{\left({\delta }_{\overline{x}\overline{y}}v\right)}^2.}$

4. Error Estimation

4.1. Local Error

We define first the local error of the scheme (48). Let $(u,\varphi )$ be the solution of (43). For $1\le n\le N-1$ we let

$${\epsilon }_{ij}^n:={\displaystyle \frac{u(x_i,y_j,t^{n+1})-u(x_i,y_j,t^n)}{\Delta t}}-{\Delta }_h{\displaystyle \frac{u(x_i,y_j,t^{n+1})+u(x_i,y_j,t^n)}{2}}$$

$$+{\displaystyle \frac{\Delta t}{4}}{\delta }_{x\overline{x}}{\delta }_{y\overline{y}}(u(x_i,y_j,t^{n+1})-u(x_i,y_j,t^n))$$

$$-{\displaystyle \frac{P\left(\varphi (x_i,y_j,t^{n+1})\right)-P\left(\varphi (x_i,y_j,t^n)\right)}{\Delta t}},$$

(17)

and

$${\zeta }_{ij}^n:={\displaystyle \frac{\varphi (x_i,y_j,t^{n+1})-\varphi (x_i,y_j,t^n)}{\Delta t}}$$

$$+{\displaystyle \frac{1}{h^2}}\left\{-a\left({\displaystyle \frac{\tilde{\varphi }(x_{i+1},y_j,t^n)+\tilde{\varphi }(x_i,y_j,t^n)}{2}}\right)\right.$$

$$*\left[\left({\displaystyle \frac{\varphi (x_{i+1},y_j,t^{n+1})+\varphi (x_{i+1},y_j,t^n)}{2}}\right)-\left({\displaystyle \frac{\varphi (x_i,y_j,t^{n+1})+\varphi (x_i,y_j,t^n)}{2}}\right)\right]$$

$$+a\left({\displaystyle \frac{\tilde{\varphi }(x_i,y_j,t^n)+\tilde{\varphi }(x_{i-1},y_j,t^n)}{2}}\right)$$

$$*\left[\left({\displaystyle \frac{\varphi (x_i,y_j,t^{n+1})+\varphi (x_i,y_j,t^n)}{2}}\right)-\left({\displaystyle \frac{\varphi (x_{i-1},y_j,t^{n+1})+\varphi (x_{i-1},y_j,t^n)}{2}}\right)\right]$$

$$-b\left({\displaystyle \frac{\tilde{\varphi }(x_i,y_{j+1},t^n)+\tilde{\varphi }(x_i,y_j,t^n)}{2}}\right)$$

$$*\left[\left({\displaystyle \frac{\varphi (x_i,y_{j+1},t^{n+1})+\varphi (x_i,y_{j+1},t^n)}{2}}\right)-\left({\displaystyle \frac{\varphi (x_i,y_j,t^{n+1})+\varphi (x_i,y_j,t^n)}{2}}\right)\right]$$

$$+b\left({\displaystyle \frac{\tilde{\varphi }(x_i,y_j,t^n)+\tilde{\varphi }(x_i,y_{j-1},t^n)}{2}}\right)$$

$$\left.*\left[\left({\displaystyle \frac{\varphi (x_i,y_j,t^{n+1})+\varphi (x_i,y_j,t^n)}{2}}\right)-\left({\displaystyle \frac{\varphi (x_i,y_{j-1},t^{n+1})+\varphi (x_i,y_{j-1},t^n)}{2}}\right)\right]\right\}$$

$$+{\displaystyle \frac{\Delta t}{4}}{\delta }_x\left(a_{i-{\displaystyle \frac{1}{2}},j}^{n+{\displaystyle \frac{1}{2}}}\left({\delta }_{\overline{x}}\left({\delta }_y\left(b_{i,j-{\displaystyle \frac{1}{2}}}^{n+{\displaystyle \frac{1}{2}}}\left({\delta }_{\overline{y}}\left(\varphi (x_i,y_j,t^{n+1})-\varphi (x_i,y_j,t^n)\right)\right)\right)\right)\right)\right)$$

$$-f(\tilde{\varphi }(x_i,y_j,t^n),\tilde{u}(x_i,y_j,t^n)),$$

(18)

where $\tilde{\varphi }(x_i,y_{j-1},t^n):={\displaystyle \frac{3}{2}}\varphi (x_i,y_{j-1},t^n)-{\displaystyle \frac{1}{2}}\varphi (x_i,y_{j-1},t^{n-1})$ and $\tilde{u}(x_i,y_{j-1},t^n):={\displaystyle \frac{3}{2}}u(x_i,y_{j-1},t^n)-{\displaystyle \frac{1}{2}}u(x_i,y_{j-1},t^{n-1}).$ Then by Taylor’s theorem we have the following local error estimate, whose proof, long but straightforward, is omitted.

Lemma 5.

Suppose that the solution $(u,\varphi )$ of (43) is sufficient smooth. Then, there exist constants $C_{\epsilon },C_{\zeta }$ independent of $\Delta t,$h such that

$$\underset{\begin{array}{c}c1\le n\le N-1\hfill \\ (x_i,y_j)\in {\Omega }_h\hfill \end{array}}{max}\left|{\epsilon }_{ij}^n\right|\le C_{\epsilon }(\Delta t^2+h^2),$$

(19)

$$\underset{\begin{array}{c}c1\le n\le N-1\hfill \\ (x_i,y_j)\in {\Omega }_h\hfill \end{array}}{max}\left|{\zeta }_{ij}^n\right|\le C_{\zeta }(\Delta t^2+h^2).$$

(20)

4.2. Error Estimate

We now prove the main result of paper.

Theorem 1.

Suppose that the solution $(u,\varphi )$ of (43) is sufficient smooth. Let $(U^n,{\Phi }^n)$ be the solution of (48). Suppose that ${\displaystyle \frac{\Delta t}{h}}\le \sigma ,$ where σ is sufficiently small and that $\mathrm{\parallel }{z}^1{\mathrm{\parallel }}_h=O(\Delta t^2+h^2),$ where $z^1={\varphi }^1-{\Phi }^1.$ Then, there exists a constant C independent of $\Delta t$ and h such that

$$\underset{0\le n\le N}{max}{\mathrm{\parallel }{u}^n-U^n\mathrm{\parallel }}_h\le C(\Delta t^2+h^2),$$

(21)

$$\underset{0\le n\le N}{max}{|{\varphi }^n-{\Phi }^n|}_{1,h}\le C(\Delta t^2+h^2).$$

(22)

Proof:

In the sequel C will denote generic constants independent of $\Delta t$ and $h.$ At various points in the proof we follow some techniques of Dendy, [5], extending them from the case of a single parabolic equation considered in [5] to the case of the parabolic system considered herein.

Let $e^n:=u^n-U^n,z^n:={\varphi }^n-{\Phi }^n.$ From (59) and (52) we obtain, for $1\le n\le N-1,$ and $1\le i,j\le J:$

$$\left.\begin{array}{c}{\displaystyle \frac{z_{ij}^{n+1}-z_{ij}^n}{\Delta t}}\hfill \\ -{\delta }_x\left(a(x_{i-{\displaystyle \frac{1}{2}}},y_j,t^{n+{\displaystyle \frac{1}{2}}},{\displaystyle \frac{\tilde{\varphi }(x_i,y_j,t^n)+\tilde{\varphi }(x_{i-1},y_j,t^n)}{2}}){\delta }_{\overline{x}}\left({\displaystyle \frac{{\varphi }_{ij}^{n+1}+{\varphi }_{ij}^n}{2}}\right)\right)\hfill \\ -{\delta }_y\left(b(x_i,y_{j-{\displaystyle \frac{1}{2}}},t^{n+{\displaystyle \frac{1}{2}}},{\displaystyle \frac{\tilde{\varphi }(x_i,y_j,t^n)+\tilde{\varphi }(x_i,y_{j-1},t^n)}{2}}){\delta }_{\overline{y}}\left({\displaystyle \frac{{\varphi }_{ij}^{n+1}+{\varphi }_{ij}^n}{2}}\right)\right)\hfill \\ +{\delta }_x\left(a(x_{i-{\displaystyle \frac{1}{2}}},y_j,t^{n+{\displaystyle \frac{1}{2}}},{\displaystyle \frac{\tilde{\Phi }(x_i,y_j,t^n)+\tilde{\Phi }(x_{i-1},y_j,t^n)}{2}}){\delta }_{\overline{x}}\left({\displaystyle \frac{{\Phi }_{ij}^{n+1}+{\Phi }_{ij}^n}{2}}\right)\right)\hfill \\ +{\delta }_y\left(b(x_i,y_{j-{\displaystyle \frac{1}{2}}},t^{n+{\displaystyle \frac{1}{2}}},{\displaystyle \frac{\tilde{\Phi }(x_i,y_j,t^n)+\tilde{\Phi }(x_i,y_{j-1},t^n)}{2}}){\delta }_{\overline{y}}\left({\displaystyle \frac{{\Phi }_{ij}^{n+1}+{\Phi }_{ij}^n}{2}}\right)\right)\hfill \\ +{\displaystyle \frac{\Delta t}{4}}{\delta }_x\left(a(x_{i-{\displaystyle \frac{1}{2}}},y_j,t^{n+{\displaystyle \frac{1}{2}}},{\displaystyle \frac{\tilde{\varphi }(x_i,y_j,t^n)+\tilde{\varphi }(x_{i-1},y_j,t^n)}{2}})\right.\hfill \\ \left.*\left({\delta }_{\overline{x}}\left({\delta }_y\left(b\left({\displaystyle \frac{\tilde{\varphi }(x_i,y_j,t^n)+\tilde{\varphi }(x_i,y_{j-1},t^n)}{2}}\right)\left({\delta }_{\overline{y}}({\varphi }_{ij}^{n+1}-{\varphi }_{ij}^n)\right)\right)\right)\right)\right)\hfill \\ -{\displaystyle \frac{\Delta t}{4}}{\delta }_x\left(a(x_{i-{\displaystyle \frac{1}{2}}},y_j,t^{n+{\displaystyle \frac{1}{2}}},{\displaystyle \frac{\tilde{\Phi }(x_i,y_j,t^n)+\tilde{\Phi }(x_{i-1},y_j,t^n)}{2}})\right.\hfill \\ \left.*\left({\delta }_{\overline{x}}\left({\delta }_y\left(b\left({\displaystyle \frac{\tilde{\Phi }(x_i,y_j,t^n)+\tilde{\Phi }(x_i,y_{j-1},t^n)}{2}}\right)\left({\delta }_{\overline{y}}({\Phi }_{ij}^{n+1}-{\Phi }_{ij}^n)\right)\right)\right)\right)\right)\hfill \\ =f(x_i,y_j,t^{n+{\displaystyle \frac{1}{2}}},\tilde{\varphi }(x_i,y_j,t^n),\tilde{u}(x_i,y_j,t^n))\hfill \\ -f(x_i,y_j,t^{n+{\displaystyle \frac{1}{2}}},\tilde{\Phi }(x_i,y_j,t^n),\tilde{U}(x_i,y_j,t^n))+{\zeta }_{ij}^n\hfill \end{array}\right\}$$

(23)

For $n\ge 1$ we define ${\varphi }_{i,j}^n:=\varphi (x_i,y_j,t^n)$, ${\tilde{\varphi }}_{ij}^n:={\displaystyle \frac{3}{2}}{\varphi }_{ij}^n-{\displaystyle \frac{1}{2}}{\varphi }_{ij}^{n-1},$ $a_{i-{\displaystyle \frac{1}{2}},j}^n\left(\widehat{\varphi }\right):=a\left({\displaystyle \frac{{\tilde{\varphi }}_{ij}^n+{\tilde{\varphi }}_{i-1,j}^n}{2}}\right)=a\left({\widehat{\varphi }}_{i-{\displaystyle \frac{1}{2}},j}^n\right),$ $i=1,\dots ,J+1,j=0,\dots ,J+1.$ $b_{i,j-{\displaystyle \frac{1}{2}}}^n\left(\widehat{\varphi }\right):=b\left({\displaystyle \frac{{\tilde{\varphi }}_{ij}^n+{\tilde{\varphi }}_{i,j-1}^n}{2}}\right)=b\left({\widehat{\varphi }}_{i,j-{\displaystyle \frac{1}{2}}}^n\right),$ $i=0,\dots ,J+1,j=1,\dots ,J+1,$ where

$${\widehat{\varphi }}_{i-{\displaystyle \frac{1}{2}},j}^n={\displaystyle \frac{1}{2}}\left({\displaystyle \frac{3}{2}}{\varphi }_{ij}^n-{\displaystyle \frac{1}{2}}{\varphi }_{ij}^{n-1}+{\displaystyle \frac{3}{2}}{\varphi }_{i-1,j}^n-{\displaystyle \frac{1}{2}}{\varphi }_{i-1,j}^{n-1}\right)$$

$$={\displaystyle \frac{3}{4}}\left({\varphi }_{ij}^n+{\varphi }_{i-1,j}^n\right)-{\displaystyle \frac{1}{4}}\left({\varphi }_{ij}^{n-1}+{\varphi }_{i-1,j}^{n-1}\right),$$

$${\widehat{\varphi }}_{i,j-{\displaystyle \frac{1}{2}}}^n={\displaystyle \frac{1}{2}}\left({\displaystyle \frac{3}{2}}{\varphi }_{ij}^n-{\displaystyle \frac{1}{2}}{\varphi }_{ij}^{n-1}+{\displaystyle \frac{3}{2}}{\varphi }_{i,j-1}^n-{\displaystyle \frac{1}{2}}{\varphi }_{i,j-1}^{n-1}\right)$$

$$={\displaystyle \frac{3}{4}}\left({\varphi }_{ij}^n+{\varphi }_{i,j-1}^n\right)-{\displaystyle \frac{1}{4}}\left({\varphi }_{ij}^{n-1}+{\varphi }_{i,j-1}^{n-1}\right),$$

We define similarly the discrete quantities $\tilde{\Phi },$$\widehat{\Phi },$$a_{i-{\displaystyle \frac{1}{2}},j}^n\left(\widehat{\Phi }\right)$ and $b_{i,j-{\displaystyle \frac{1}{2}}}^n\left(\widehat{\Phi }\right)$ and we use the notation

$$F_{ij}^n(\tilde{\varphi },\tilde{u})-F_{ij}^n(\tilde{\Phi },\tilde{U}):=f({\tilde{\varphi }}_{ij}^n,{\tilde{u}}_{ij}^n)-f({\tilde{\Phi }}_{ij}^n,{\tilde{U}}_{ij}^n).$$

We finally define for $1\le \nu \le N$

$$\mathrm{\parallel }{z}^{\nu }{\mathrm{\parallel }}_n^2:=h^2\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}\left(a_{i-{\displaystyle \frac{1}{2}},j}^n\left({\widehat{\Phi }}^n\right){\left({\delta }_{\overline{x}}{z}_{ij}^{\nu }\right)}^2+b_{i,j-{\displaystyle \frac{1}{2}}}^n\left({\widehat{\Phi }}^n\right){\left({\delta }_{\overline{y}}{z}_{ij}^{\nu }\right)}^2\right).$$

(24)

The proof proceeds inductively. Our induction hypothesis is: Suppose that the following inequality holds for $1\le m\le M:$

$${\displaystyle \frac{1}{\Delta t}}\sum _{n=1}^{m-1}\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2+\mathrm{\parallel }{z}^m{\mathrm{\parallel }}_{m-1}^2+\Delta t{\displaystyle \frac{\lambda }{4}}\sum _{n=1}^{m-1}{|z^{n+1}-z^n|}_{\overline{x}\overline{y},h}^2$$

$$\le C_2\Delta t\sum _{n=1}^{m-1}{\mathrm{\parallel }{z}^n\mathrm{\parallel }}_{n-1}^2+C_2(\Delta t^4+h^4).\left(I\right)$$

Obviously, $\left(I\right)$ holds for $m=1.$

Taking the ${\ell }^2$ inner product of both sides of (64) with $z^{n+1}-z^n$ we have for $n\ge 2$

$$\left.\begin{array}{c}\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2\hfill \\ -\Delta t{\left({\delta }_x\left(a^n\left(\widehat{\Phi }\right){\delta }_{\overline{x}}\left({\displaystyle \frac{z^{n+1}+z^n}{2}}\right)\right)+{\delta }_y\left(b^n\left(\widehat{\Phi }\right){\delta }_{\overline{y}}\left({\displaystyle \frac{z^{n+1}+z^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ +{\displaystyle \frac{\Delta t^2}{4}}{\left({\delta }_x\left(a^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{y}}(z^{n+1}-z^n)\right)\right)\right)\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ =\Delta t{\left(F^n(\tilde{\varphi },\tilde{u})-F^n(\tilde{\Phi },\tilde{U}),z^{n+1}-z^n\right)}_h+\Delta t{\left({\zeta }^n,z^{n+1}-z^n\right)}_h\hfill \\ +\Delta t{\left({\delta }_x\left([a^n\left(\widehat{\varphi }\right)-a^n\left(\widehat{\Phi }\right)]{\delta }_{\overline{x}}\left({\displaystyle \frac{{\varphi }^{n+1}+{\varphi }^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ +\Delta t{\left({\delta }_y\left([b^n\left(\widehat{\varphi }\right)-b^n\left(\widehat{\Phi }\right)]{\delta }_{\overline{y}}\left({\displaystyle \frac{{\varphi }^{n+1}+{\varphi }^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ -{\displaystyle \frac{\Delta t^2}{4}}{\left({\delta }_x\left(a^n\left(\widehat{\varphi }\right)\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left(\widehat{\varphi }\right)\left({\delta }_{\overline{y}}({\varphi }^{n+1}-{\varphi }^n)\right)\right)\right)\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ +{\displaystyle \frac{\Delta t^2}{4}}{\left({\delta }_x\left(a^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{y}}({\varphi }^{n+1}-{\varphi }^n)\right)\right)\right)\right)\right),z^{n+1}-z^n\right)}_h.\hfill \end{array}\right\}$$

(25)

By Lemma 8, the symmetry of $L_h^n,$ and our hypotheses on $a,b,$ we have, using the equivalence of the norms ${\mathrm{\parallel }·\mathrm{\parallel }}_n$ and ${|·|}_{1,h}$ with constants independent of h and $\Delta t,$ and the notation ${\widehat{z}}^n:={\widehat{\varphi }}^n-{\widehat{\Phi }}^n:$

$$\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_n^2:=h^2\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}\left(a_{i-{\displaystyle \frac{1}{2}},j}^n\left({\widehat{\Phi }}^n\right){\left({\delta }_{\overline{x}}{z}_{ij}^n\right)}^2+b_{i,j-{\displaystyle \frac{1}{2}}}^n\left({\widehat{\Phi }}^n\right){\left({\delta }_{\overline{y}}{z}_{ij}^n\right)}^2\right)\le$$

$$\le \mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+C{\left(|{\widehat{z}}^n-{\widehat{z}}^{n-1}|{\delta }_{\overline{x}}{z}^n,{\delta }_{\overline{x}}{z}^n\right)}_h$$

$$+C{\left(|{\widehat{z}}^n-{\widehat{z}}^{n-1}|{\delta }_{\overline{y}}{z}^n,{\delta }_{\overline{y}}{z}^n\right)}_h+C\Delta t{\mathrm{\parallel }{z}^n\mathrm{\parallel }}_{n-1}^2.$$

Therefore

$$-{\left({\delta }_x\left(a^n\left(\widehat{\Phi }\right){\delta }_{\overline{x}}\left({\displaystyle \frac{z^{n+1}+z^n}{2}}\right)\right)+{\delta }_y\left(b^n\left(\widehat{\Phi }\right){\delta }_{\overline{y}}\left({\displaystyle \frac{z^{n+1}+z^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h\ge$$

$$\ge {\displaystyle \frac{1}{2}}\mathrm{\parallel }{z}^{n+1}{\mathrm{\parallel }}_n^2-{\displaystyle \frac{1}{2}}[\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+C{\left(|{\widehat{z}}^n-{\widehat{z}}^{n-1}|{\delta }_{\overline{x}}{z}^n,{\delta }_{\overline{x}}{z}^n\right)}_h$$

$$+C{\left(|{\widehat{z}}^n-{\widehat{z}}^{n-1}|{\delta }_{\overline{y}}{z}^n,{\delta }_{\overline{y}}{z}^n\right)}_h+C^{\prime }\Delta t\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2].$$

(26)

In addition, there holds that

$${\left(|{\widehat{z}}^n-{\widehat{z}}^{n-1}|{\delta }_{\overline{x}}{z}^n,{\delta }_{\overline{x}}{z}^n\right)}_h\le {\displaystyle \frac{C}{h^{1/3}}}\left(\right|{\widehat{z}}^n{|}_{1,h}+|{\widehat{z}}^{n-1}{|}_{1,h})\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2.$$

(27)

To see this consider Hölder’s inequality:

$$\left|\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{f}_{ij}{g}_{ij}\right|\le {\left(\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{\left|f_{ij}\right|}^p\right)}^{{\displaystyle \frac{1}{p}}}{\left(\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{\left|g_{ij}\right|}^q\right)}^{{\displaystyle \frac{1}{q}}}for{\displaystyle \frac{1}{p}}+{\displaystyle \frac{1}{q}}=1,$$

(28)

and

$$\left|\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{f}_{ij}{g}_{ij}{r}_{ij}\right|\le {\left(\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{\left|f_{ij}\right|}^p\right)}^{{\displaystyle \frac{1}{p}}}{\left(\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{\left|g_{ij}\right|}^q\right)}^{{\displaystyle \frac{1}{q}}}{\left(\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{\left|r_{ij}\right|}^s\right)}^{{\displaystyle \frac{1}{s}}}$$

(29)

$for{\displaystyle \frac{1}{p}}+{\displaystyle \frac{1}{q}}+{\displaystyle \frac{1}{s}}=1.$

Applying (70) for $p=6,q={\displaystyle \frac{12}{5}},ands={\displaystyle \frac{12}{5}},$ we have,

$${\left(|{\widehat{z}}^n-{\widehat{z}}^{n-1}|{\delta }_{\overline{x}}{z}^n,{\delta }_{\overline{x}}{z}^n\right)}_h\le \mathrm{\parallel }{\widehat{z}}^n-{\widehat{z}}^{n-1}{\mathrm{\parallel }}_{6,h}\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n{\mathrm{\parallel }}_{{\displaystyle \frac{12}{5}},h}{\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n\mathrm{\parallel }}_{{\displaystyle \frac{12}{5}},h},$$

where for $v\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}$ we define ${\displaystyle {\mathrm{\parallel }v\mathrm{\parallel }}_{p,h}:=(h^2\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}|v_{ij}{{|}^p)}^{1/p}}$ (with ${\mathrm{\parallel }v\mathrm{\parallel }}_{2,h}\equiv {\mathrm{\parallel }v\mathrm{\parallel }}_h$).

Therefore

$${\left(|{\widehat{z}}^n-{\widehat{z}}^{n-1}|{\delta }_{\overline{x}}{z}^n,{\delta }_{\overline{x}}{z}^n\right)}_h\le \mathrm{\parallel }{\widehat{z}}^n-{\widehat{z}}^{n-1}{\mathrm{\parallel }}_{6,h}{\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n\mathrm{\parallel }}_{{\displaystyle \frac{12}{5}},h}^2.$$

From (69) we have

$$\left|\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{f}_{ij}^r{g}_{ij}^s\right|\le {\left(\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{\left|f_{ij}\right|}^{rp}\right)}^{{\displaystyle \frac{1}{p}}}{\left(\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{\left|g_{ij}\right|}^{sq}\right)}^{{\displaystyle \frac{1}{q}}}for{\displaystyle \frac{1}{p}}+{\displaystyle \frac{1}{q}}=1.$$

(30)

Using (71) with $r={\displaystyle \frac{9}{5}},s={\displaystyle \frac{3}{5}},p={\displaystyle \frac{10}{9}},q=10$ we see that

$$\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n{\mathrm{\parallel }}_{{\displaystyle \frac{12}{5}},h}\le \mathrm{\parallel }{\delta }_{\overline{x}}{z}^n{\mathrm{\parallel }}^{3/4}{\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n\mathrm{\parallel }}_{6,h}^{1/4}.$$

Note now that the inverse inequality ${\mathrm{\parallel }v\mathrm{\parallel }}_{p,h}\le C{h}^{-N({\displaystyle \frac{1}{2}}-{\displaystyle \frac{1}{p}})}{\mathrm{\parallel }v\mathrm{\parallel }}_h$ holds on ${\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}.$ For $p=6$ we have

$$\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n{\mathrm{\parallel }}_{6,h}\le C{h}^{-2({\displaystyle \frac{1}{2}}-{\displaystyle \frac{1}{6}})}\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n{\mathrm{\parallel }}_h=C{h}^{-{\displaystyle \frac{2}{3}}}{\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n\mathrm{\parallel }}_h.$$

Hence

$$\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n{\mathrm{\parallel }}_{{\displaystyle \frac{12}{5}},h}\le C\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n{\mathrm{\parallel }}_h^{3/4}(h^{-{\displaystyle \frac{2}{3}}}\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n{\mathrm{\parallel }}_h{)}^{1/4}=C{h}^{-{\displaystyle \frac{2}{12}}}\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n{\mathrm{\parallel }}_h=C{h}^{-{\displaystyle \frac{1}{6}}}{\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n\mathrm{\parallel }}_h.$$

By the above we condude that

$$\left(|{\widehat{z}}^n-{\widehat{z}}^{n-1}|{\delta }_{\overline{x}}{z}^n,{\delta }_{\overline{x}}{z}^n\right)\le \mathrm{\parallel }{\widehat{z}}^n-{\widehat{z}}^{n-1}{\mathrm{\parallel }}_{6,h}{\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n\mathrm{\parallel }}_{{\displaystyle \frac{12}{5}},h}^2$$

$$\le C|{\widehat{z}}^n-{\widehat{z}}^{n-1}{|}_{1,h}{\left(C{h}^{-{\displaystyle \frac{1}{6}}}{\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n\mathrm{\parallel }}_h\right)}^2=C{h}^{-{\displaystyle \frac{1}{3}}}|{\widehat{z}}^n-{\widehat{z}}^{n-1}{|}_{1,h}{\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n\mathrm{\parallel }}_h^2$$

$$\le C{h}^{-{\displaystyle \frac{1}{3}}}|{\widehat{z}}^n-{\widehat{z}}^{n-1}{|}_{1,h}{\mathrm{\parallel }{z}^n\mathrm{\parallel }}_{n-1}^2,$$

from which (68) follows. Now, using the induction hypotheses $\left(I\right)$ we will show that $|{\widehat{z}}^n{|}_{1,h},{\left|{\widehat{z}}^{n-1}\right|}_{1,h}$ and

${\displaystyle max\left\{\underset{i,j}{max}\right|{\delta }_{\overline{x}}{\widehat{\Phi }}_{ij}^n|,\underset{i,j}{max}|{\delta }_{\overline{y}}{\widehat{\Phi }}_{ij}^n\left|\right\}}$ are bounded for $n\le M.$

Indeed, from $\left(I\right)$ we have for a small enagh $ϵ>0$ (that we consider fixed from now on) that

$$\mathrm{\parallel }{z}^m{\mathrm{\parallel }}_{m-1}^2\le C\Delta t\sum _{n=1}^m{\mathrm{\parallel }{z}^n\mathrm{\parallel }}_{n-1}^2+C\Delta t^4,m\le M.$$

Hence, from Gronwall’s lemma we obtain

$$|z^m{|}_{1,h}^2\le e^{Cm\Delta t}\Delta t^{4}form\le M.$$

(31)

Hence, using the stability condition ${\displaystyle \frac{\Delta t}{h}}\le \sigma ),$ we see that

$$\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}|\left[{\left({\delta }_{\overline{x}}{\widehat{z}}_{ij}^m\right)}^2+{\left({\delta }_{\overline{y}}{\widehat{z}}_{ij}^m\right)}^2\right]\le C{\displaystyle \frac{\Delta t^4}{h^2}}\le C\Delta t^2.$$

Hence

$$max\left(\underset{i,j}{max}|{\delta }_{\overline{x}}{\widehat{z}}_{ij}|,\underset{i,j}{max}\left|{\delta }_{\overline{y}}{\widehat{z}}_{ij}\right|\right)\le \sqrt{{\sum }_{i=1}^{J+1}{\sum }_{j=1}^{J+1}\left[{\left({\delta }_{\overline{x}}{\widehat{z}}_{ij}^m\right)}^2+{\left({\delta }_{\overline{y}}{\widehat{z}}_{ij}^m\right)}^2\right]}\le C\Delta t,$$

and therefore

$$|{\delta }_{\overline{x}}{\widehat{\Phi }}_{ij}^n|\le |{\delta }_{\overline{x}}{\widehat{z}}_{ij}^n|+|{\delta }_{\overline{x}}{\widehat{\varphi }}_{ij}^n|\le C,$$

(32)

since $|{\delta }_{\overline{x}}{\widehat{\varphi }}_{ij}^n|$ is bounded, as it is easily seen. Similarly we see that $|{\delta }_{\overline{y}}{\widehat{\Phi }}_{ij}^n|\le C.$

Now, from (72), supposing that $n\le M,$ we have, by the stability conttion that

$${\displaystyle |{\widehat{z}}^n{|}_{1,h}+{\left|{\widehat{z}}^{n-1}\right|}_{1,h}\le e^{CT}\Delta t^2\le C{e}^{CT}\Delta th\le C{e}^{CT}{h}^{\frac{2}{3}}\Delta t{h}^{\frac{1}{3}}\le C\Delta t{h}^{\frac{1}{3}}.}$$

Therefore from (68) we have, for $n\le M$

$$\left(|{\widehat{z}}^n-{\widehat{z}}^{n-1}|{\delta }_{\overline{x}}{z}^n,{\delta }_{\overline{x}}{z}^n\right)\le C\Delta t{\mathrm{\parallel }{z}^n\mathrm{\parallel }}_{n-1}^2,$$

(33)

and a similar bound for the term $\left(|{\widehat{z}}^n-{\widehat{z}}^{n-1}|{\delta }_{\overline{y}}{z}^n,{\delta }_{\overline{y}}{z}^n\right).$ Substituting (74) in (67) we have, for $n\le M$

$$-{\left({\delta }_x\left(a^n\left(\widehat{\Phi }\right){\delta }_{\overline{x}}\left({\displaystyle \frac{z^{n+1}+z^n}{2}}\right)\right)+{\delta }_y\left(b^n\left(\widehat{\Phi }\right){\delta }_{\overline{y}}\left({\displaystyle \frac{z^{n+1}+z^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h$$

$$\ge {\displaystyle \frac{1}{2}}\mathrm{\parallel }{z}^{n+1}{\mathrm{\parallel }}_n^2-{\displaystyle \frac{1}{2}}(\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+C\Delta t\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2)$$

(34)

Hence, from (66) and (76) we obtain, for $n\le M$

$\begin{array}{c}\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2\hfill \\ +\Delta t\left[{\displaystyle \frac{1}{2}}\mathrm{\parallel }{z}^{n+1}{\mathrm{\parallel }}_n^2-{\displaystyle \frac{1}{2}}(\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+C\Delta t\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2)\right]\hfill \\ +{\displaystyle \frac{\Delta t^2}{4}}{\left({\delta }_x\left(a^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{y}}(z^{n+1}-z^n)\right)\right)\right)\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ \le \Delta t{\left(F^n(\tilde{\varphi },\tilde{u})-F^n(\tilde{\Phi },\tilde{U}),z^{n+1}-z^n\right)}_h+\Delta t{\left({\zeta }^n,z^{n+1}-z^n\right)}_h\hfill \\ +\Delta t{\left({\delta }_x\left([a^n\left(\widehat{\varphi }\right)-a^n\left(\widehat{\Phi }\right)]{\delta }_{\overline{x}}\left({\displaystyle \frac{{\varphi }^{n+1}+{\varphi }^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ +\Delta t{\left({\delta }_y\left([b^n\left(\widehat{\varphi }\right)-b^n\left(\widehat{\Phi }\right)]{\delta }_{\overline{y}}\left({\displaystyle \frac{{\varphi }^{n+1}+{\varphi }^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ -{\displaystyle \frac{\Delta t^2}{4}}{\left({\delta }_x\left(a^n\left(\widehat{\varphi }\right)\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left(\widehat{\varphi }\right)\left({\delta }_{\overline{y}}({\varphi }^{n+1}-{\varphi }^n)\right)\right)\right)\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ +{\displaystyle \frac{\Delta t^2}{4}}{\left({\delta }_x\left(a^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{y}}({\varphi }^{n+1}-{\varphi }^n)\right)\right)\right)\right)\right),z^{n+1}-z^n\right)}_h.\hfill \end{array}$

Therefore, using Lemma 9, we conclude for $n\le M$ that

$$\left.\begin{array}{c}\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2+\Delta t\left[{\displaystyle \frac{1}{2}}\mathrm{\parallel }{z}^{n+1}{\mathrm{\parallel }}_n^2-{\displaystyle \frac{1}{2}}(\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+C\Delta t\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2)\right]\hfill \\ +{\displaystyle \frac{\Delta t^2}{4}}\left(\lambda \right|z^{n+1}-z^n{|}_{\overline{x}\overline{y},h}^2-\mu |z^{n+1}-z^n{|}_{1,h}^2)\le \hfill \\ \le \Delta t{\left(F^n(\tilde{\varphi },\tilde{u})-F^n(\tilde{\Phi },\tilde{U}),z^{n+1}-z^n\right)}_h+\Delta t{\left({\zeta }^n,z^{n+1}-z^n\right)}_h\hfill \\ +\Delta t{\left({\delta }_x\left([a^n\left(\widehat{\varphi }\right)-a^n\left(\widehat{\Phi }\right)]{\delta }_{\overline{x}}\left({\displaystyle \frac{{\varphi }^{n+1}+{\varphi }^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ +\Delta t{\left({\delta }_y\left([b^n\left(\widehat{\varphi }\right)-b^n\left(\widehat{\Phi }\right)]{\delta }_{\overline{y}}\left({\displaystyle \frac{{\varphi }^{n+1}+{\varphi }^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ -{\displaystyle \frac{\Delta t^2}{4}}{\left({\delta }_x\left(a^n\left(\widehat{\varphi }\right)\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left(\widehat{\varphi }\right)\left({\delta }_{\overline{y}}({\varphi }^{n+1}-{\varphi }^n)\right)\right)\right)\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ +{\displaystyle \frac{\Delta t^2}{4}}{\left({\delta }_x\left(a^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{y}}({\varphi }^{n+1}-{\varphi }^n)\right)\right)\right)\right)\right),z^{n+1}-z^n\right)}_h.\hfill \end{array}\right\}$$

(35)

We now define the following quantities

$\begin{array}{ccc}I& :=& {\left({\delta }_x\left([a^n\left(\widehat{\varphi }\right)-a^n\left(\widehat{\Phi }\right)]{\delta }_{\overline{x}}\left({\displaystyle \frac{{\varphi }^{n+1}+{\varphi }^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ & & +{\left({\delta }_y\left([b^n\left(\widehat{\varphi }\right)-b^n\left(\widehat{\Phi }\right)]{\delta }_{\overline{y}}\left({\displaystyle \frac{{\varphi }^{n+1}+{\varphi }^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h,\hfill \\ II& :=& {\left({\delta }_x\left(a^n\left(\widehat{\varphi }\right)\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left(\widehat{\varphi }\right)\left({\delta }_{\overline{y}}({\varphi }^{n+1}-{\varphi }^n)\right)\right)\right)\right)\right),z^{n+1}-z^n\right)}_h,\hfill \\ III& :=& {\left({\delta }_x\left(a^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{y}}({\varphi }^{n+1}-{\varphi }^n)\right)\right)\right)\right)\right),z^{n+1}-z^n\right)}_h.\hfill \end{array}$

From (72), since $n\le M,$ we have $\mathrm{\parallel }{\widehat{z}}^n{\mathrm{\parallel }}_h\le C{\left|{\widehat{z}}^n\right|}_{1,h}\le e^{CT}\Delta t^2.$

Long but straightforward calculations yield now for $\epsilon >0$ that

$$\Delta t\left|I\right|\le C_{\epsilon }\Delta t^6+\epsilon {\mathrm{\parallel }{z}^{n+1}-z^n\mathrm{\parallel }}_h^2,$$

(36)

and

$$\Delta t^2|III-II|\le C_{\epsilon }\Delta t^6+\epsilon \Delta t^2{|z^{n+1}-z^n|}_{\overline{x},\overline{y},h}^2.$$

(37)

Substituting in (76) we have

$$\begin{array}{c}\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2+\Delta t\left[{\displaystyle \frac{1}{2}}\mathrm{\parallel }{z}^{n+1}{\mathrm{\parallel }}_n^2-{\displaystyle \frac{1}{2}}(\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+C\Delta t\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2)\right]\hfill \\ +{\displaystyle \frac{\Delta t^2}{4}}\left(\lambda \right|z^{n+1}-z^n{|}_{\overline{x}\overline{y},h}^2-\mu |z^{n+1}-z^n{|}_{1,h}^2)\le \hfill \\ \le \Delta t{\left(F^n(\tilde{\varphi },\tilde{u})-F^n(\tilde{\Phi },\tilde{U}),z^{n+1}-z^n\right)}_h+\Delta t{\left({\zeta }^n,z^{n+1}-z^n\right)}_h\hfill \\ +C_{{\epsilon }^{\prime }}\Delta t^6+{\epsilon }^{\prime }{\mathrm{\parallel }{z}^{n+1}-z^n\mathrm{\parallel }}_h^2\hfill \\ +C_{{\epsilon }^{\prime \prime }}\Delta t^6+{\epsilon }^{\prime \prime }\Delta t^2{|z^{n+1}-z^n|}_{\overline{x},\overline{y},h}^2,\hfill \end{array}$$

where ${\epsilon }^{\prime },{\epsilon }^{\prime \prime }>0$ arbitrary.

Hence, using the Cauchy-Schwarz inequality and the arithmetic geometric mean inequality we have

$$\begin{array}{c}\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2+\Delta t\left[{\displaystyle \frac{1}{2}}\mathrm{\parallel }{z}^{n+1}{\mathrm{\parallel }}_n^2-{\displaystyle \frac{1}{2}}(\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+C\Delta t\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2)\right]\hfill \\ +{\displaystyle \frac{\Delta t^2}{4}}\left(\lambda \right|z^{n+1}-z^n{|}_{\overline{x}\overline{y},h}^2-\mu |z^{n+1}-z^n{|}_{1,h}^2)\le \hfill \\ \le \left({\displaystyle \frac{\Delta t^2}{4{\epsilon }^{\prime \prime \prime }}}\mathrm{\parallel }{F}^n(\tilde{\varphi },\tilde{u})-F^n(\tilde{\Phi },\tilde{U}){\mathrm{\parallel }}_h^2+{\epsilon }^{\prime \prime \prime }{\mathrm{\parallel }{z}^{n+1}-z^n\mathrm{\parallel }}_h^2\right)\hfill \\ +\left({\displaystyle \frac{\Delta t^2}{4{\epsilon }^{\prime \prime \prime }}}\mathrm{\parallel }{\zeta }^n{\mathrm{\parallel }}_h^2+{\epsilon }^{\prime \prime \prime }{\mathrm{\parallel }{z}^{n+1}-z^n\mathrm{\parallel }}_h^2\right)\hfill \\ +C_{{\epsilon }^{\prime },{\epsilon }^{\prime \prime }}\Delta t^6+{\epsilon }^{\prime }\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2+{\epsilon }^{\prime \prime }\Delta t^2{|z^{n+1}-z^n|}_{\overline{x},\overline{y},h}^2.\hfill \end{array}$$

In addition since $|F_{\varphi }|\le c,|F_u|\le c$ we have from the mean value theorem

$$\mathrm{\parallel }{F}^n(\tilde{\varphi },\tilde{u})-F^n(\tilde{\Phi },\tilde{U}){\mathrm{\parallel }}_h\le L_1\mathrm{\parallel }{\tilde{z}}^n{\mathrm{\parallel }}_h+L_2{\mathrm{\parallel }{\tilde{e}}^n\mathrm{\parallel }}_h$$

$$\le L_1^{\prime }(\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_h+\mathrm{\parallel }{z}^{n-1}{\mathrm{\parallel }}_h)+L_2^{\prime }(\mathrm{\parallel }{e}^n{\mathrm{\parallel }}_h+\mathrm{\parallel }{e}^{n-1}{\mathrm{\parallel }}_h).$$

Since for $n\le M$ (72) gives $\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_h\le C{\left|z^n\right|}_{1,h}\le e^{CT}\Delta t^2$, we have for ${\epsilon }^{\prime \prime \prime }>0$

$$\mathrm{\parallel }{F}^n(\tilde{\varphi },\tilde{u})-F^n(\tilde{\Phi },\tilde{U}){\mathrm{\parallel }}_h^2\le C\Delta t^4+C(\mathrm{\parallel }{e}^n{\mathrm{\parallel }}_h^2+\mathrm{\parallel }{e}^{n-1}{\mathrm{\parallel }}_h^2),$$

$$\begin{array}{c}\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2\hfill \\ +\Delta t\left[{\displaystyle \frac{1}{2}}\mathrm{\parallel }{z}^{n+1}{\mathrm{\parallel }}_n^2-{\displaystyle \frac{1}{2}}(\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+C\Delta t\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2)\right]\hfill \\ +{\displaystyle \frac{\Delta t^2}{4}}\left(\lambda \right|z^{n+1}-z^n{|}_{\overline{x}\overline{y},h}^2-\mu |z^{n+1}-z^n{|}_{1,h}^2)\hfill \\ \le \left({\displaystyle \frac{\Delta t^2}{4{\epsilon }^{\prime \prime \prime }}}(C\Delta t^4+C(\mathrm{\parallel }{e}^n{\mathrm{\parallel }}_h^2+\mathrm{\parallel }{e}^{n-1}{\mathrm{\parallel }}_h^2\left)\right)+{\epsilon }^{\prime \prime \prime }{\mathrm{\parallel }{z}^{n+1}-z^n\mathrm{\parallel }}_h^2\right)\hfill \\ +\left({\displaystyle \frac{\Delta t^2}{4{\epsilon }^{\prime \prime \prime }}}\mathrm{\parallel }{\zeta }^n{\mathrm{\parallel }}_h^2+{\epsilon }^{\prime \prime \prime }{\mathrm{\parallel }{z}^{n+1}-z^n\mathrm{\parallel }}_h^2\right)\hfill \\ +C_{{\epsilon }^{\prime },{\epsilon }^{\prime \prime }}\Delta t^6+{\epsilon }^{\prime }\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2+{\epsilon }^{\prime \prime }\Delta t^2{|z^{n+1}-z^n|}_{\overline{x},\overline{y},h}^2.\hfill \end{array}$$

Then, for suitable ${\epsilon }^{\prime },{\epsilon }^{\prime \prime },{\epsilon }^{\prime \prime \prime }=O\left(\epsilon \right),$ we see that

$$\begin{array}{c}(1-\epsilon )\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2\hfill \\ +{\displaystyle \frac{\Delta t}{2}}\left[\mathrm{\parallel }{z}^{n+1}{\mathrm{\parallel }}_n^2-{\mathrm{\parallel }{z}^n\mathrm{\parallel }}_{n-1}^2\right]+\Delta t^2({\displaystyle \frac{\lambda }{4}}-\epsilon ){|z^{n+1}-z^n|}_{\overline{x}\overline{y},h}^2\hfill \\ \le C\Delta t^2\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+\Delta t^2{C}_{\epsilon }(\mathrm{\parallel }{e}^n{\mathrm{\parallel }}_h^2+\mathrm{\parallel }{e}^{n-1}{\mathrm{\parallel }}_h^2)\hfill \\ +\Delta t^2{C}_{\epsilon }\mathrm{\parallel }{\zeta }^n{\mathrm{\parallel }}_h^2+C_{\epsilon }\Delta t^6+{\displaystyle \frac{\Delta t^2}{4}}\mu {|z^{n+1}-z^n|}_{1,h}^2.\hfill \end{array}$$

Hence, since $|z^{n+1}+z^n{|}_{1,h}\le {\displaystyle \frac{C_{*}}{h}}{\mathrm{\parallel }{z}^{n+1}+z^n\mathrm{\parallel }}_h,$ at ${\displaystyle \frac{\Delta t}{h}}\le \sigma ,$ we have that ${\displaystyle \frac{\Delta t^2}{4}}\mu |z^{n+1}+z^n{|}_{1,h}^2\le C^{\prime }{\sigma }^2{\mathrm{\parallel }{z}^{n+1}+z^n\mathrm{\parallel }}_h$ where $C^{\prime }={\displaystyle \frac{C_{*}\mu }{4}}$

Therefore

$$\begin{array}{c}(1-\epsilon -C^{\prime }{\sigma }^2){\mathrm{\parallel }{z}^{n+1}-z^n\mathrm{\parallel }}_h^2\hfill \\ +{\displaystyle \frac{\Delta t}{2}}\left[\mathrm{\parallel }{z}^{n+1}{\mathrm{\parallel }}_n^2-{\mathrm{\parallel }{z}^n\mathrm{\parallel }}_{n-1}^2\right]+\Delta t^2({\displaystyle \frac{\lambda }{4}}-\epsilon ){|z^{n+1}-z^n|}_{\overline{x}\overline{y},h}^2\hfill \\ \le C\Delta t^2\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+\Delta t^2{C}_{\epsilon }(\mathrm{\parallel }{e}^n{\mathrm{\parallel }}_h^2\hfill \\ +\mathrm{\parallel }{e}^{n-1}{\mathrm{\parallel }}_h^2)+\Delta t^2{C}_{\epsilon }{\mathrm{\parallel }{\zeta }^n\mathrm{\parallel }}_h^2+C_{\epsilon }\Delta t^6,n\le M.\hfill \end{array}$$

Summing from $n=1tom$ we have

$${\displaystyle \frac{2(1-\epsilon -C^{\prime }{\sigma }^2)}{\Delta t}}\sum _{n=1}^m\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2+{\mathrm{\parallel }{z}^{m+1}\mathrm{\parallel }}_m^2$$

$$+2\Delta t({\displaystyle \frac{\lambda }{4}}-\epsilon )\sum _{n=1}^m{|z^{n+1}-z^n|}_{\overline{x}\overline{y},h}^2$$

$$\le C\Delta t\sum _{n=1}^m\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+C_{\epsilon }\Delta t\sum _{n=1}^m{\mathrm{\parallel }{e}^n\mathrm{\parallel }}_h^2$$

$$+C_{\epsilon }\Delta t\sum _{n=1}^m{\mathrm{\parallel }{\zeta }^n\mathrm{\parallel }}_h^2+C_{\epsilon }\Delta t^{4}1\le m\le M,\forall \epsilon >0.$$

(38)

Now, from (58)and (50) we have for $0\le n\le N-1,1\le i,j\le J:$

$${\displaystyle \frac{e_{ij}^{n+1}-e_{ij}^n}{\Delta t}}-{\Delta }_h\left({\displaystyle \frac{e_{ij}^{n+1}+e_{ij}^n}{2}}\right)+{\displaystyle \frac{\Delta t}{4}}{\delta }_{x\overline{x}}{\delta }_{y\overline{y}}(e_{ij}^{n+1}-e_{ij}^n)$$

$$={\displaystyle \frac{1}{\Delta t}}\left\{[P(x_i,y_j,t^{n+1},\varphi (x_i,y_j,t^{n+1}))-P(x_i,y_j,t^n,\varphi (x_i,y_j,t^n))]\right\}$$

$$-{\displaystyle \frac{1}{\Delta t}}\left\{[P(x_i,y_j,t^{n+1},\Phi (x_i,y_j,t^{n+1}))-P(x_i,y_j,t^n,\Phi (x_i,y_j,t^n))]\right\}+{\epsilon }_{ij}^n.$$

We define ${\widehat{P}}_{ij}{\left(v\right)}^n:=P(x_i,y_j,t^{n+1},v(x_i,y_j,t^{n+1}))-P(x_i,y_j,t^n,v(x_i,y_j,t^n)).$ Taking the ${\ell }^2$ inner product of both sides of the above relation with $e^{n+1}+e^n$ we have

$${\left(e^{n+1}-e^n,e^{n+1}+e^n\right)}_h-\Delta t{\left({\Delta }_h\left({\displaystyle \frac{e^{n+1}+e^n}{2}}\right),e^{n+1}+e^n\right)}_h$$

$$+{\displaystyle \frac{\Delta t^2}{4}}{\left({\delta }_{x\overline{x}}{\delta }_{y\overline{y}}(e^{n+1}-e^n),e^{n+1}+e^n\right)}_h$$

$$={\left(\widehat{P}{\left(\varphi \right)}^n-\widehat{P}{(\Phi )}^n,e^{n+1}+e^n\right)}_h+\Delta t{\left({\epsilon }^n,e^{n+1}+e^n\right)}_h.$$

(39)

Using the symmetry of ${\beta }_h(·,·),$ Lemma 4 (with $a\equiv 1,b\equiv 1),$ and Lemma 7 and summing the above relation from $n=0$ to $m^{\prime }-1$ we have, since $e^0=0$

$$\mathrm{\parallel }{e}^{m^{\prime }}{\mathrm{\parallel }}_h^2+{\displaystyle \frac{\Delta t}{2}}\sum _{n=0}^{m^{\prime }-1}{|e^{n+1}+e^n|}_{1,h}^2$$

$$\le \sum _{n=0}^{m^{\prime }-1}{\left(\widehat{P}{\left(\varphi \right)}^n-\widehat{P}{(\Phi )}^n,e^{n+1}+e^n\right)}_h+\Delta t\sum _{n=0}^{m^{\prime }-1}{\left({\epsilon }^n,e^{n+1}+e^n\right)}_h$$

$$\le \sum _{n=0}^{m^{\prime }-1}\mathrm{\parallel }\widehat{P}{\left(\varphi \right)}^n-\widehat{P}{(\Phi )}^n{\mathrm{\parallel }}_h\mathrm{\parallel }{e}^{n+1}+e^n{\mathrm{\parallel }}_h+\Delta t\sum _{n=0}^{m^{\prime }-1}\mathrm{\parallel }{\epsilon }^n{\mathrm{\parallel }}_h{\mathrm{\parallel }{e}^{n+1}+e^n\mathrm{\parallel }}_h$$

$$\le {\displaystyle \frac{1}{4\epsilon \Delta t}}\sum _{n=0}^{m^{\prime }-1}\mathrm{\parallel }\widehat{P}{\left(\varphi \right)}^n-\widehat{P}{(\Phi )}^n{\mathrm{\parallel }}_h^2+\epsilon \Delta t\sum _{n=0}^{m^{\prime }-1}{\mathrm{\parallel }{e}^{n+1}+e^n\mathrm{\parallel }}_h^2$$

$$+{\displaystyle \frac{\Delta t}{4\epsilon }}\sum _{n=0}^{m^{\prime }-1}\mathrm{\parallel }{\epsilon }^n{\mathrm{\parallel }}_h^2+\Delta t\epsilon \sum _{n=0}^{m^{\prime }-1}{\mathrm{\parallel }{e}^{n+1}+e^n\mathrm{\parallel }}_h^2$$

$$\le {\displaystyle \frac{1}{4\epsilon \Delta t}}\sum _{n=0}^{m^{\prime }-1}\mathrm{\parallel }\widehat{P}{\left(\varphi \right)}^n-\widehat{P}{(\Phi )}^n{\mathrm{\parallel }}_h^2+C_{\epsilon }\Delta t\sum _{n=0}^{m^{\prime }-1}\mathrm{\parallel }{\epsilon }^n{\mathrm{\parallel }}_h^2+\Delta t\epsilon \sum _{n=0}^{N-1}{\mathrm{\parallel }{e}^{n+1}+e^n\mathrm{\parallel }}_h^2.$$

In view of (79) we have

$$\mathrm{\parallel }{e}^{m^{\prime }}{\mathrm{\parallel }}_h^2\le C\left(C_{\epsilon }\Delta t\sum _{n=0}^{m^{\prime }-1}\mathrm{\parallel }{e}^n{\mathrm{\parallel }}_h^2+C_{\epsilon }\Delta t\sum _{n=1}^{m^{\prime }-1}\mathrm{\parallel }{\zeta }^n{\mathrm{\parallel }}_h^2+C_{\epsilon }\Delta t^4+{\mathrm{\parallel }{z}^1\mathrm{\parallel }}_0^2\right)$$

$$+c\Delta t\sum _{n=0}^{m^{\prime }-1}{\mathrm{\parallel }{\epsilon }^n\mathrm{\parallel }}_h^2.$$

Hence, by Gronwall lemma we infer that

$$\underset{0\le n\le m^{\prime }}{max}{\mathrm{\parallel }{e}^n\mathrm{\parallel }}_h^2\le C\left[\Delta t\sum _{n=1}^{m-1}\mathrm{\parallel }{\zeta }^n{\mathrm{\parallel }}_h^2+\Delta t\sum _{n=0}^{m-1}\mathrm{\parallel }{\epsilon }^n{\mathrm{\parallel }}_h^2+\Delta t^4+{\mathrm{\parallel }{z}^1\mathrm{\parallel }}_0^2\right].$$

Therefore, by (60) and (61)

$$\underset{0\le n\le m}{max}{\mathrm{\parallel }{e}^n\mathrm{\parallel }}_h^2\le C\left[\Delta t\sum _{n=1}^{m-1}{(\Delta t^2+h^2)}^2+\Delta t^4+{\mathrm{\parallel }{z}^1\mathrm{\parallel }}_0^2\right]$$

$$\le C\left(N\Delta t{(\Delta t^2+h^2)}^2+{\mathrm{\parallel }{z}^1\mathrm{\parallel }}_0^2\right)$$

$$\le C\left({(\Delta t^2+h^2)}^2+{\mathrm{\parallel }{z}^1\mathrm{\parallel }}_0^2\right).$$

Hence, using $\mathrm{\parallel }{z}^1{\mathrm{\parallel }}_0=O(\Delta t^2+h^2)$ we have

$$\underset{0\le n\le m}{max}{\mathrm{\parallel }{e}^n\mathrm{\parallel }}_h\le C(\Delta t^2+h^2).$$

(40)

Since

$$\Delta t\sum _{n=0}^m{\mathrm{\parallel }{e}^n\mathrm{\parallel }}_h^2\le C{(\Delta t^2+h^2)}^2,$$

(41)

substituting (82) and (61) in (79) we have

$${\displaystyle \frac{2(1-\epsilon -C^{\prime }{\sigma }^2)}{\Delta t}}\sum _{n=1}^m\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2+\mathrm{\parallel }{z}^{m+1}{\mathrm{\parallel }}_m^2+2\Delta t({\displaystyle \frac{\lambda }{4}}-\epsilon )\sum _{n=1}^m{|z^{n+1}-z^n|}_{\overline{x}\overline{y},h}^2$$

$$\le C\Delta t\sum _{n=1}^m{\mathrm{\parallel }{z}^n\mathrm{\parallel }}_{n-1}^2+C_{\epsilon }(\Delta t^4+h^4)$$

(42)

Take now $\epsilon >0$ so that $2(1-\epsilon -C^{\prime }{\sigma }^2)\ge 1and2({\displaystyle \frac{\lambda }{4}}-\epsilon )\ge {\displaystyle \frac{\lambda }{4}}$ so (83) gives

$${\displaystyle \frac{1}{\Delta t}}\sum _{n=1}^m\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2+\mathrm{\parallel }{z}^{m+1}{\mathrm{\parallel }}_m^2+\Delta t{\displaystyle \frac{\lambda }{4}}\sum _{n=1}^m{|z^{n+1}-z^n|}_{\overline{x}\overline{y},h}^2$$

$$\le C\Delta t\sum _{n=1}^m{\mathrm{\parallel }{z}^n\mathrm{\parallel }}_{n-1}^2+C(\Delta t^4+h^4),$$

which is $\left(I\right)$ for m instead of $m-1.$ Hence, the inductive step is complete and $\left(I\right)$ holds for $1\le m\le N.$ We conclude from (81) that (62) holds.

Finally, it follows from $\left(I\right)$ that

$$\mathrm{\parallel }{z}^m{\mathrm{\parallel }}_{m-1}^2\le C\Delta t\sum _{n=1}^{m-1}{\mathrm{\parallel }{z}^n\mathrm{\parallel }}_{n-1}^2+C(\Delta t^4+h^4),1\le m\le N$$

from which, by Gronwall’s lemma it follows that

$$\underset{1\le m\le N}{max}{\mathrm{\parallel }{z}^m\mathrm{\parallel }}_{m-1}\le C(\Delta t^4+h^4),$$

which implies (63).

▪

We may easily obtain from Taylor’s theorem an approximation at the first step satisfying $\mathrm{\parallel }{z}^1{\mathrm{\parallel }}_h=O(h^2+\Delta t^2).$

5. Numerical Verification of the Order of Convergence

In the sequel, we solved the system (43) under Neumann boundary conditions on $\partial \Omega .$ (We obtained similar rates for the problems with Dirichlet boundary conditions.) We took $\Omega =[0,1]\times [0,1],$ $T=0.1,$$a=b=1,$ $f=\varphi (1-\varphi )u/(1+0.25u),$$p={\varphi }^3(10-15\varphi +6{\varphi }^2).$ We added a suitable nonhomogeneous term so that the solution of the associated initial-boundary-value problem with Neumann boundary conditions on $\partial \Omega$ was $(\varphi ,u)=(e^{-t}cos\left(x(x-1)\right)cos\left(y(y-1)\right),e^{-t}cos\left(\pi x\right)cos\left(\pi y\right)).$

Table 1. Errors and order of convergence of the ADI-ADI scheme, T=0.1.

	${\mathrm{\parallel }·\mathrm{\parallel }}_{\mathit{\infty }}$		${\mathrm{\parallel }·\mathrm{\parallel }}_{\mathit{h}}$
J	errors	order	errors	order
50	$\varphi :1.920e-04$	$--$	$8.664e-05$	$--$
	$u:1.047e-04$	$--$	$4.762e-05$	$--$
75	$\varphi :8.708e-05$	$1.95$	$3.851e-05$	$1.99$
	$u:4.725e-05$	$1.96$	$2.119e-05$	$1.99$
112	$\varphi :3.979e-05$	$1.93$	$1.731e-05$	$1.97$
	$u:2.143e-05$	$1.94$	$9.513e-06$	$1.97$
168	$\varphi :1.806e-05$	$1.94$	$7.747e-06$	$1.98$
	$u:9.622e-06$	$1.97$	$4.231e-06$	$1.99$
252	$\varphi :8.254e-06$	$1.93$	$3.486e-06$	$1.96$
	$u:4.320e-06$	$1.97$	$1.883e-06$	$1.99$
378	$\varphi :3.814e-06$	$1.90$	$1.582e-06$	$1.94$
	$u:1.944e-06$	$1.96$	$8.383e-07$	$1.99$
567	$\varphi :1.795e-06$	$1.85$	$7.300e-07$	$1.90$
	$u:8.783e-07$	$1.95$	$3.736e-07$	$1.99$

Table 1. Errors and order of convergence of the ADI-ADI scheme, T=0.1.

	${\mathrm{\parallel }·\mathrm{\parallel }}_{\mathit{\infty }}$		${\mathrm{\parallel }·\mathrm{\parallel }}_{\mathit{h}}$
J	errors	order	errors	order
50	$\varphi :1.920e-04$	$--$	$8.664e-05$	$--$
	$u:1.047e-04$	$--$	$4.762e-05$	$--$
75	$\varphi :8.708e-05$	$1.95$	$3.851e-05$	$1.99$
	$u:4.725e-05$	$1.96$	$2.119e-05$	$1.99$
112	$\varphi :3.979e-05$	$1.93$	$1.731e-05$	$1.97$
	$u:2.143e-05$	$1.94$	$9.513e-06$	$1.97$
168	$\varphi :1.806e-05$	$1.94$	$7.747e-06$	$1.98$
	$u:9.622e-06$	$1.97$	$4.231e-06$	$1.99$
252	$\varphi :8.254e-06$	$1.93$	$3.486e-06$	$1.96$
	$u:4.320e-06$	$1.97$	$1.883e-06$	$1.99$
378	$\varphi :3.814e-06$	$1.90$	$1.582e-06$	$1.94$
	$u:1.944e-06$	$1.96$	$8.383e-07$	$1.99$
567	$\varphi :1.795e-06$	$1.85$	$7.300e-07$	$1.90$
	$u:8.783e-07$	$1.95$	$3.736e-07$	$1.99$

Table 2 shows the errors and the experimental orders of convergence of the numerical solution computed by the ADI-ADI scheme for the this problem, in the discrete maximum norm ${\mathrm{\parallel }·\mathrm{\parallel }}_{\infty }$ and also in the ${\ell }_2$ norm ${\mathrm{\parallel }·\mathrm{\parallel }}_h.$ We computed with $h=1/(J+1),$ where $J=50,75,\dots ,567,$ and $\Delta t=0.01h.$ The observed order of convergence is approximately equal to 2 for both variables in both norms. (For two consecutive runs with values $h_1,h_2$ of $h,$ giving errors $e_1,e_2,$ respectively, the rate of convergence was computed as $log(e_2/e_1)/log(h_2/h_1).$)

With this scheme we performed two numerical experiments where we took (cf. Wang[15]) on $\Omega =[0,1]\times [0,1],$ ${\varphi }_0(x,y)={\displaystyle \frac{1}{2}}\left[1+tanh{\displaystyle \frac{r-R_0}{2\sqrt{2}\epsilon }}\right],$ $u_0(x,y)=\left\{\begin{array}{cc}ll{t}_s\hfill & r<R_0\\ {\displaystyle \frac{\left[ln\left(\frac{r}{R_0}\right)+t_{s}ln\left(\frac{r}{R_{00}}\right)\right]}{ln\left(\frac{R_0}{R_{00}}\right)}}\hfill & R_0\le R_{00}\\ -1\hfill & r\ge R_{00}\end{array}\right.$, $r=\sqrt{x^2+y^2}$,

$a=b=1$, $q=1/m(1+{\delta }_{\mu }cos\left(4\theta \right))$, where $\theta =arctan({\varphi }_y/{\varphi }_x),$$f={\displaystyle \frac{1}{{\epsilon }^2}}(\varphi (1-\varphi )(\varphi -{\displaystyle \frac{1}{2}}+30\epsilon \alpha S{\displaystyle \frac{u}{1+0.25u}}\varphi (1-\varphi )),$$p={\varphi }^3(10-15\varphi +6{\varphi }^2)/S$, S=dimensionless supercooling of the melt, m=ratio of the capillary to the kinetic length and $\epsilon =$ratio of the average interface thickness parameter.

In the isotropic case, with values of physical parameters given by $t_s=0.01,$$R_0=0.1,$$R_{00}=2R_0,$$S=0.8,$ $m=0.1,$$\alpha =70,$${\delta }_{\mu }=0,$ and $\epsilon =1/400,$ the evolution of the solution is’ given in Figures 1, $2\left(a\right)$(we took $h={10}^{-3}$, $\Delta t=2·{10}^{-5}$).

6. The Problem

This “isotropic” model consists of a system (1) of two p.d.e’s is supplemented by given initial conditions $\varphi (x,y,0)={\varphi }_0(x,y)$, $u(x,y,0)=u_0(x,y)$, $(x,y)\in \Omega$, and boundary conditions of Neumann or Dirichlet type for $\varphi$ and u on the boundary $\partial \Omega$ of $\Omega$ for $t\ge 0$. Let $\Omega =(\alpha ,\beta )\times (\alpha ,\beta ).$ For $(x,y,t)\in \overline{\Omega }\times [0,T],$ where $T>0,$ we consider the initial-boundary-value problem

$$\left\{\begin{array}{cc}{\displaystyle \frac{\partial \varphi }{\partial t}-\frac{\partial }{\partial x}\left(a\left(\varphi \right)\frac{\partial \varphi }{\partial x}\right)-\frac{\partial }{\partial y}\left(b\left(\varphi \right)\frac{\partial \varphi }{\partial y}\right)=f(\varphi ,u),}\hfill & (x,y,t)\in \overline{\Omega }\times [0,T],\\ {\displaystyle \frac{\partial u}{\partial t}-\frac{{\partial }^{2}u}{\partial x^2}-\frac{{\partial }^{2}u}{\partial y^2}=\frac{\partial }{\partial t}p\left(\varphi \right),}\hfill & (x,y,t)\in \overline{\Omega }\times [0,T],\\ \varphi (x,y,t)=0,u(x,y,t)=0\hfill & (x,y,t)\in \partial \Omega \times [0,T],\\ \varphi (x,y,0)=g^{\varphi }(x,y)\hfill & (x,y)\in \overline{\Omega },\\ u(x,y,0)=g^u(x,y)\hfill & (x,y)\in \overline{\Omega },\end{array}\right.$$

(43)

where a and $b,$ are smooth functions of $\varphi ,$ such that $0<c_a^0\le a\left(\xi \right)\le c_a^1,$$0<c_b^0\le b\left(\xi \right)\le c_b^1,$ for $\xi \in \mathrm{I}\mathrm{R},$ and f and p are smooth functions defined on $\mathrm{I}{\mathrm{R}}^2$ and $\mathrm{I}\mathrm{R}$ respectively, with $f(0,0)=0$ and $p\left(0\right)=0.$ We suppose that $g^{\varphi }=g^u=0$ on $\partial \Omega .$ Furthermore, we make the assumption (in order to simplify the proofs of the error estimates) that the functions $f,a,b$, and p are globally Lipschitz functions of their arguments, i.e. that there exists a constant C such that:

$$|f(z_1,z_2)-f(Z_1,Z_2)|\le C(|z_1-Z_1|+|z_2-Z_2\left|\right),$$

(44)

$$|a\left(z_1\right)-a\left(Z_1\right)|\le C|z_1-Z_1|,$$

(45)

$$|b\left(z_1\right)-b\left(Z_1\right)|\le C|z_1-Z_1|,$$

(46)

$$|p\left(z_1\right)-p\left(Z_1\right)|\le C|z_1-Z_1|,$$

(47)

for all $z_1,z_2,Z_1,Z_2\in \mathrm{I}\mathrm{R}.$ We will assume that the initial-boundary-value problem (43) has a unique solution $(u,\varphi )$, smooth enough for the purpose of the numerical approximations.

7. Notation and Preliminary Results

We discretize the initial-boundary-value problem (43) as follows. Let $h=(\beta -\alpha )/(J+1),$ where J is a positive integer and $x_i=\alpha +ih,$$y_j=\alpha +jh,$$0\le i,j\le J+1$. We define ${\Omega }_h:=\left\{(x_i,y_j),i,j=1,\dots ,J\right\}$ and $\partial {\Omega }_h:=\left\{(x_i,y_j),i=0ori=J+1orj=0orj=J+1\right\}$. Let N be positive integer and $t^n=n\Delta t,$$n=0,\dots ,N,$ where $\Delta t=T/N,$ and define $t^{n+{\displaystyle \frac{1}{2}}}=t^n+\Delta t/2,$$x_{i\pm {\displaystyle \frac{1}{2}}}=x_i\pm h/2,$$y_{j\pm {\displaystyle \frac{1}{2}}}=y_j\pm h/2,$ and

${\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}:=\left\{U={(U_{00},\dots ,U_{J+1,J+1})}^T\in \mathrm{I}{\mathrm{R}}^{(J+2)\times (J+2)}:\right.$$\left.U_{ij}=0on\partial {\Omega }_h\right\}.$ We approximate the solution $u,\varphi$ of (43) by mesh functions $U^n,{\Phi }^n\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}$ as follows: For $n=0,1,2,\dots ,N,$ we approximate the vectors ${(u(x_0,y_0,t^n),\dots ,u(x_{J+1},y_{J+1},t^n))}^T,$${(\varphi (x_0,y_0,t^n),\dots ,\varphi (x_{J+1},y_{J+1},t^n))}^T$ by $U^n,{\Phi }^n\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}$ satisfying the following finite difference scheme:

$$\left.\begin{array}{cc}\left(\mathrm{i}\right){\Phi }_{ij}^0=g^{\varphi }(x_i,y_j),\hfill & (x_i,y_j)\in {\Omega }_h\bigcup \partial {\Omega }_h,\\ \left(\mathrm{ii}\right)U_{ij}^0=g^u(x_i,y_j),\hfill & (x_i,y_j)\in {\Omega }_h\bigcup \partial {\Omega }_h,\\ \mathrm{For} n=0,1,\dots ,N-1:\hfill & \\ {\displaystyle \left(\mathrm{iii}\right)\frac{{\Phi }_{ij}^{n+\frac{1}{2}}-{\Phi }_{ij}^n}{\frac{1}{2}\Delta t}-{\delta }_x\left(a_{i-\frac{1}{2},j}^n{\delta }_{\overline{x}}{\Phi }_{ij}^{n+\frac{1}{2}}\right)}\hfill & \\ {\displaystyle -{\delta }_y\left(b_{i,j-\frac{1}{2}}^n{\delta }_{\overline{y}}{\Phi }_{ij}^n\right)=F_{ij}^n,}\hfill & (x_i,y_j)\in {\Omega }_h,\\ {\displaystyle \left(\mathrm{iv}\right)\frac{{\Phi }_{ij}^{n+1}-{\Phi }_{ij}^{n+\frac{1}{2}}}{\frac{1}{2}\Delta t}-{\delta }_x\left(a_{i-\frac{1}{2},j}^n{\delta }_{\overline{x}}{\Phi }_{ij}^{n+\frac{1}{2}}\right)}\hfill & \\ {\displaystyle -{\delta }_y\left(b_{i,j-\frac{1}{2}}^n{\delta }_{\overline{y}}{\Phi }_{ij}^{n+1}\right)=F_{ij}^n,}\hfill & (x_i,y_j)\in {\Omega }_h,\\ \left(\mathrm{v}\right){\Phi }_{ij}^{n+{\displaystyle \frac{1}{2}}}={\Phi }_{ij}^{n+1}=0,\hfill & (x_i,y_j)\in \partial {\Omega }_h,\\ {\displaystyle \left(\mathrm{vi}\right)\frac{U_{ij}^{n+\frac{1}{2}}-U_{ij}^n}{\frac{1}{2}\Delta t}-{\delta }_{x\overline{x}}{U}_{ij}^{n+\frac{1}{2}}}\hfill & \\ {\displaystyle -{\delta }_{y\overline{y}}{U}_{ij}^n=\frac{P_{ij}^{n+1}-P_{ij}^n}{\Delta t},}\hfill & (x_i,y_j)\in {\Omega }_h,\\ {\displaystyle \left(\mathrm{vii}\right)\frac{U_{ij}^{n+1}-U_{ij}^{n+\frac{1}{2}}}{\frac{1}{2}\Delta t}-{\delta }_{x\overline{x}}{U}_{ij}^{n+\frac{1}{2}}}\hfill & \\ {\displaystyle -{\delta }_{y\overline{y}}{U}_{ij}^{n+1}=\frac{P_{ij}^{n+1}-P_{ij}^n}{\Delta t},}\hfill & (x_i,y_j)\in {\Omega }_h,\\ \left(\mathrm{viii}\right)U_{ij}^{n+{\displaystyle \frac{1}{2}}}=U_{ij}^{n+1}=0,\hfill & (x_i,y_j)\in \partial {\Omega }_h,\end{array}\right\}$$

(48)

where $P_{i,j}^n:=p\left({\Phi }_{i,j}^n\right),$${\delta }_x{v}_{ij}:={\displaystyle \frac{v_{i+1,j}-v_{i,j}}{h}},$${\delta }_{\overline{x}}{v}_{ij}:={\displaystyle \frac{v_{i,j}-v_{i-1,j}}{h}},$

${\delta }_y{v}_{ij}:={\displaystyle \frac{v_{i,j+1}-v_{i,j}}{h}},{\delta }_{\overline{y}}{v}_{ij}:={\displaystyle \frac{v_{i,j}-v_{i,j-1}}{h}},$${\delta }_{x\overline{x}}{v}_{ij}:={\displaystyle \frac{v_{i+1,j}-2v_{i,j}+v_{i-1,j}}{h^2}},$

${\delta }_{y\overline{y}}{v}_{ij}:={\displaystyle \frac{v_{i,j+1}-2v_{i,j}+v_{i,j-1}}{h^2}},$${\Delta }_h{v}_{ij}:={\delta }_{x\overline{x}}{v}_{ij}+{\delta }_{y\overline{y}}{v}_{ij}.$ We also let for $(x_i,y_j)\in {\Omega }_h$ $a_{i+{\displaystyle \frac{1}{2}},j}^n:=a\left({\displaystyle \frac{{\Phi }_{i+1,j}^n+{\Phi }_{i,j}^n}{2}}\right),$$a_{i-{\displaystyle \frac{1}{2}},j}^n:=a\left({\displaystyle \frac{{\Phi }_{i,j}^n+{\Phi }_{i-1,j}^n}{2}}\right),$ $a_{i,j+{\displaystyle \frac{1}{2}}}^n:=a\left({\displaystyle \frac{{\Phi }_{i,j+1}^n+{\Phi }_{i,j}^n}{2}}\right),$$a_{i,j-{\displaystyle \frac{1}{2}}}^n:=a\left({\displaystyle \frac{{\Phi }_{i,j}^n+{\Phi }_{i,j-1}^n}{2}}\right).$ (The quantities $b_{i\pm {\displaystyle \frac{1}{2}},j}^n,$$b_{i,j\pm {\displaystyle \frac{1}{2}}}^n$ are defined analogously.) In addition, we let

$$F_{ij}^n:=\left\{\begin{array}{cc}f({\check{\Phi }}_{ij}^{{\displaystyle \frac{1}{2}}},{\check{U}}_{ij}^{{\displaystyle \frac{1}{2}}}),\hfill & \mathrm{if}n=0,\\ f({\displaystyle \frac{3}{2}}{\Phi }_{ij}^n-{\displaystyle \frac{1}{2}}{\Phi }_{ij}^{n-1},{\displaystyle \frac{3}{2}}{U}_{ij}^n-{\displaystyle \frac{1}{2}}{U}_{ij}^{n-1}),\hfill & \mathrm{if}n\ge 1,\end{array}\right.$$

where

$${\check{\Phi }}_{ij}^{{\displaystyle \frac{1}{2}}}:=g^{\varphi }(x_i,y_j)+{\displaystyle \frac{\Delta t}{2}}\left\{a\left(g^{\varphi }(x_i,y_j)\right)g_{xx}^{\varphi }(x_i,y_j)+b\left(g^{\varphi }(x_i,y_j)\right)g_{yy}^{\varphi }(x_i,y_j)\right.$$

$$\left.+a^{\prime }\left(g^{\varphi }(x_i,y_j)\right){\left(g_x^{\varphi }(x_i,y_j)\right)}^2+b^{\prime }\left(g^{\varphi }(x_i,y_j)\right){\left(g_y^{\varphi }(x_i,y_j)\right)}^2+f(g^{\varphi }(x_i,y_j),g^u(x_i,y_j)\right\}$$

and

$${\check{U}}_{ij}^{{\displaystyle \frac{1}{2}}}:=g^u(x_i,y_j)+{\displaystyle \frac{\Delta t}{2}}\left\{\Delta g^u(x_i,y_j)+{\displaystyle \frac{2}{\Delta t}}{p}^{\prime }\left(g^{\varphi }(x_i,y_j)\right)[{\check{\Phi }}_{ij}^{{\displaystyle \frac{1}{2}}}-g^{\varphi }(x_i,y_j)]\right\}.$$

Subtracting (48.(vi)) from (48.(vii)) we have for $(x_i,y_j)\in {\Omega }_h$

$${\displaystyle \frac{U_{ij}^{n+1}-2U_{ij}^{n+\frac{1}{2}}+U_{ij}^n}{\frac{1}{2}\Delta t}}-{\delta }_{y\overline{y}}(U_{ij}^{n+1}-U_{ij}^n)=0.$$

Solving for $U_{ij}^{n+{\displaystyle \frac{1}{2}}}$ we get

$$U_{ij}^{n+{\displaystyle \frac{1}{2}}}={\displaystyle \frac{U_{ij}^{n+1}+U_{ij}^n}{2}}-{\displaystyle \frac{\Delta t}{4}}{\delta }_{y\overline{y}}(U_{ij}^{n+1}-U_{ij}^n).$$

(49)

Replacing the above relation in (48.(vi)) we get for $0\le n\le N-1$

$${\displaystyle \frac{U_{ij}^{n+1}-U_{ij}^n}{\Delta t}}-{\Delta }_h{\displaystyle \frac{U_{ij}^{n+1}+U_{ij}^n}{2}}+{\displaystyle \frac{\Delta t}{4}}{\delta }_{x\overline{x}}{\delta }_{y\overline{y}}(U_{ij}^{n+1}-U_{ij}^n)={\displaystyle \frac{P_{ij}^{n+1}-P_{ij}^n}{\Delta t}},(x_i,y_j)\in {\Omega }_h.$$

(50)

Conversely, suppose that (50) holds. Then, defining $U^{n+{\displaystyle \frac{1}{2}}}$from (49) we see that (48.(vi)) and (48.(vii)) are valid. Hence the relations (48.(vi)-(vii)) are equivalent to (49)-(50). Similarly, subtracting (48.(iii)) from (48.(iv)), we see, for $(x_i,y_j)\in {\Omega }_h$

$${\displaystyle \frac{{\Phi }_{ij}^{n+1}-2{\Phi }_{ij}^{n+\frac{1}{2}}+{\Phi }_{ij}^n}{\frac{1}{2}\Delta t}}-{\delta }_y\left(b_{i,j-{\displaystyle \frac{1}{2}}}^n{\delta }_{\overline{y}}({\Phi }_{ij}^{n+1}-{\Phi }_{ij}^n)\right)=0,$$

and therefore

$${\Phi }_{ij}^{n+{\displaystyle \frac{1}{2}}}={\displaystyle \frac{{\Phi }_{ij}^{n+1}+{\Phi }_{ij}^n}{2}}-{\displaystyle \frac{\Delta t}{4}}{\delta }_y\left(b_{i,j-{\displaystyle \frac{1}{2}}}^n{\delta }_{\overline{y}}({\Phi }_{ij}^{n+1}-{\Phi }_{ij}^n)\right).$$

(51)

Finally, substituting in (48.(iii)) we get for $0\le n\le N-1,(x_i,y_j)\in {\Omega }_h:$

$${\displaystyle \frac{{\Phi }_{ij}^{n+1}-{\Phi }_{ij}^n}{\Delta t}}-{\delta }_x\left(a_{i-{\displaystyle \frac{1}{2}},j}^n{\delta }_{\overline{x}}\left({\displaystyle \frac{{\Phi }_{ij}^{n+1}+{\Phi }_{ij}^n}{2}}\right)\right)-{\delta }_y\left(b_{i,j-{\displaystyle \frac{1}{2}}}^n{\delta }_{\overline{y}}\left({\displaystyle \frac{{\Phi }_{ij}^{n+1}+{\Phi }_{ij}^n}{2}}\right)\right)$$

$$+{\displaystyle \frac{\Delta t}{4}}{\delta }_x\left(a_{i-{\displaystyle \frac{1}{2}},j}^n\left({\delta }_{\overline{x}}\left({\delta }_y\left(b_{i,j-{\displaystyle \frac{1}{2}}}^n\left({\delta }_{\overline{y}}({\Phi }_{ij}^{n+1}-{\Phi }_{ij}^n)\right)\right)\right)\right)\right)=F_{ij}^n.$$

(52)

Conversely, let (52) hold. Then, if we define ${\Phi }^{n+{\displaystyle \frac{1}{2}}}$ from (51), we see that (48.(iii)) and (48.(iv)) are equivalent to (52) and (51).

While the solution algorithm is given by the equations (48), the error analysis will be based on (50) and (52). It is clear that solving in (48) for ${\Phi }^{n+{\displaystyle \frac{1}{2}}},$${\Phi }^{n+1},$$U^{n+{\displaystyle \frac{1}{2}}},$$U^{n+1}$ requires the solution of $J,$$J\times J$ symmetric tridiagonal systems of equations in each case. There systems are clearly invertible in view of our assumptions on the functions a and $b.$

In what follows we list a few auxiliary helpful results whose proofs are standard and may be found in [10].

Lemma 6.

a) 

Let $\varphi \in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}},\psi \in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}.$ Then

$$\sum _{i=1}^J({\varphi }_{i+1,j}-{\varphi }_{i,j}){\psi }_{ij}=-\sum _{i=1}^{J+1}{\varphi }_{ij}({\psi }_{i,j}-{\psi }_{i-1,j}).$$

b) 

Let $\varphi \in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}},\psi \in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}.$ Then

$$\sum _{j=1}^J({\varphi }_{i,j+1}-{\varphi }_{i,j}){\psi }_{ij}=-\sum _{j=1}^{J+1}{\varphi }_{ij}({\psi }_{i,j}-{\psi }_{i,j-1}).$$

c) 

If $\varphi ,v\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}},$ then for $(x_i,y_j)\in {\Omega }_h:$

i) 

${\delta }_{x\overline{x}}{v}_{ij}:={\delta }_x{\delta }_{\overline{x}}{v}_{ij}={\delta }_{\overline{x}}{\delta }_x{v}_{ij},{\delta }_{y\overline{y}}{v}_{ij}:={\delta }_y{\delta }_{\overline{y}}{v}_{ij}={\delta }_{\overline{y}}{\delta }_y{v}_{ij}.$

ii) 

${\delta }_x{\delta }_y{v}_{ij}={\delta }_y{\delta }_x{v}_{ij},$${\delta }_{\overline{x}}{\delta }_y{v}_{ij}={\delta }_y{\delta }_{\overline{x}}{v}_{ij},$${\delta }_{\overline{y}}{\delta }_x{v}_{ij}={\delta }_x{\delta }_{\overline{y}}{v}_{ij}.$

iii) 

${\delta }_{\overline{x}\overline{y}}{v}_{ij}:={\delta }_{\overline{x}}{\delta }_{\overline{y}}{v}_{ij}={\delta }_{\overline{y}}{\delta }_{\overline{x}}{v}_{ij}.$

iv) 

${\delta }_{\overline{x}}\left({\varphi }_{ij}\left({\delta }_{\overline{y}}{v}_{ij}\right)\right)=\left[{\delta }_{\overline{x}}{\varphi }_{ij}{\delta }_{\overline{y}}{v}_{ij}+{\varphi }_{ij}{\delta }_{\overline{x}\overline{y}}{v}_{ij}\right]-h{\delta }_{\overline{x}}{\varphi }_{i,j}{\delta }_{\overline{x}\overline{y}}{v}_{i,j},$

${\delta }_{\overline{y}}\left({\varphi }_{ij}\left({\delta }_{\overline{x}}{v}_{ij}\right)\right)=\left[{\delta }_{\overline{y}}{\varphi }_{ij}{\delta }_{\overline{x}}{v}_{ij}+{\varphi }_{ij}{\delta }_{\overline{x}\overline{y}}{v}_{ij}\right]-h{\delta }_{\overline{y}}{\varphi }_{i,j}{\delta }_{\overline{x}\overline{y}}{v}_{i,j}.$

v) 

${\delta }_x\left({\varphi }_{ij}{v}_{ij}\right)=v_{ij}{\delta }_x{\varphi }_{ij}+{\varphi }_{ij}{\delta }_x{v}_{ij}+h\left({\delta }_x{\varphi }_{i,j}\right)\left({\delta }_x{v}_{i,j}\right).$

For $v,w\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}},$ define the ${\ell }_2$ inner product on ${\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}$ by ${\displaystyle {(v,w)}_h:=h^2\sum _{i=1}^J\sum _{j=1}^J{v}_{ij}{w}_{ij},}$ with corresponding norm ${\displaystyle {\mathrm{\parallel }v\mathrm{\parallel }}_h:=h{\left(\sum _{i=1}^J\sum _{j=1}^J{v}_{ij}^2\right)}^{\frac{1}{2}}.}$

Lemma 7.

Define ${\Delta }_h:{\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}\to {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}$ by

$${\left({\Delta }_{h}v\right)}_{i,j}:=\left\{\begin{array}{cc}{\displaystyle \frac{v_{i+1,j}+v_{i-1,j}+v_{i,j+1}+v_{i,j-1}-4v_{ij}}{h^2},}\hfill & (x_i,y_j)\in {\Omega }_h,\\ 0,\hfill & (x_i,y_j)\in \partial {\Omega }_h\end{array}\right.$$

Then, there hold that

$${({\Delta }_{h}v,w)}_h={(v,{\Delta }_{h}w)}_h\forall v,w\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}},$$

(53)

and

$${({\Delta }_{h}v,v)}_h=-{\left|v\right|}_{1,h}^2\forall v\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}},$$

(54)

where ${\left|v\right|}_{1,h}^2:=h^2{\sum }_{i=1}^{J+1}{\sum }_{j=1}^{J+1}[{\left({\delta }_{\overline{x}}{v}_{ij}\right)}^2+{\left({\delta }_{\overline{y}}{v}_{ij}\right)}^2]$ is a discrete $H_1$ norm.

Lemma 8.

Define $L_h^n:{\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}\to {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}$ by

$${\left(L_h^{n}v\right)}_{i,j}:=\left\{\begin{array}{cc}{\displaystyle {\delta }_x\left(a_{i-\frac{1}{2},j}^n{\delta }_{\overline{x}}{v}_{ij}\right)+{\delta }_y\left(b_{i,j-\frac{1}{2}}^n{\delta }_{\overline{y}}{v}_{ij}\right)}\hfill & (x_i,y_j)\in {\Omega }_h,\\ 0\hfill & (x_i,y_j)\in \partial {\Omega }_h.\end{array}\right.$$

Then, there holds that

$${(L_h^{n}v,w)}_h={(v,L_h^{n}w)}_h\forall v,w\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}$$

(55)

and

$${(L_h^{n}v,v)}_h=-h^2\sum _{j=1}^{J+1}\sum _{i=1}^{J+1}\left\{a_{i+{\displaystyle \frac{1}{2}},j}^n{\left({\delta }_{\overline{x}}{v}_{ij}\right)}^2+b_{i,j+{\displaystyle \frac{1}{2}}}^n{\left({\delta }_{\overline{y}}{v}_{ij}\right)}^2\right\}\le 0\forall v\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}.$$

(56)

Lemma 9.

Define, for $v,w\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}},$ the bilinear form ${\beta }_h(v,w)$ as

$${\displaystyle {\beta }_h^n(v,w):={\left({\delta }_x\left(a^n\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left({\delta }_{\overline{y}}v\right)\right)\right)\right)\right),w\right)}_h:=}$$

$${\displaystyle h^2\sum _{i=1}^J\sum _{j=1}^J{\delta }_x\left(a_{i-\frac{1}{2},j}^n\left({\delta }_{\overline{x}}\left({\delta }_y\left(b_{i,j-\frac{1}{2}}^n\left({\delta }_{\overline{y}}{v}_{ij}\right)\right)\right)\right)\right)w_{ij}.}$$

Then, for some constants $\lambda ,\mu$ independent of $h,$ there holds that

$${\beta }_h^n(v,v)\ge {\lambda \left|v\right|}_{\overline{x}\overline{y},h}^2-\mu {\left|v\right|}_{1,h}^2\forall v\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}},$$

(57)

where ${\displaystyle {\left|v\right|}_{\overline{x}\overline{y},h}^2:=h^2\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{\left({\delta }_{\overline{x}\overline{y}}v\right)}^2.}$

8. Error Estimation

8.1. Local Error

We define first the local error of the scheme (48). Let $(u,\varphi )$ be the solution of (43). For $1\le n\le N-1$ we let

$${\epsilon }_{ij}^n:={\displaystyle \frac{u(x_i,y_j,t^{n+1})-u(x_i,y_j,t^n)}{\Delta t}}-{\Delta }_h{\displaystyle \frac{u(x_i,y_j,t^{n+1})+u(x_i,y_j,t^n)}{2}}$$

$$+{\displaystyle \frac{\Delta t}{4}}{\delta }_{x\overline{x}}{\delta }_{y\overline{y}}(u(x_i,y_j,t^{n+1})-u(x_i,y_j,t^n))$$

$$-{\displaystyle \frac{P\left(\varphi (x_i,y_j,t^{n+1})\right)-P\left(\varphi (x_i,y_j,t^n)\right)}{\Delta t}},$$

(58)

and

$${\zeta }_{ij}^n:={\displaystyle \frac{\varphi (x_i,y_j,t^{n+1})-\varphi (x_i,y_j,t^n)}{\Delta t}}$$

$$+{\displaystyle \frac{1}{h^2}}\left\{-a\left({\displaystyle \frac{\tilde{\varphi }(x_{i+1},y_j,t^n)+\tilde{\varphi }(x_i,y_j,t^n)}{2}}\right)\right.$$

$$*\left[\left({\displaystyle \frac{\varphi (x_{i+1},y_j,t^{n+1})+\varphi (x_{i+1},y_j,t^n)}{2}}\right)-\left({\displaystyle \frac{\varphi (x_i,y_j,t^{n+1})+\varphi (x_i,y_j,t^n)}{2}}\right)\right]$$

$$+a\left({\displaystyle \frac{\tilde{\varphi }(x_i,y_j,t^n)+\tilde{\varphi }(x_{i-1},y_j,t^n)}{2}}\right)$$

$$*\left[\left({\displaystyle \frac{\varphi (x_i,y_j,t^{n+1})+\varphi (x_i,y_j,t^n)}{2}}\right)-\left({\displaystyle \frac{\varphi (x_{i-1},y_j,t^{n+1})+\varphi (x_{i-1},y_j,t^n)}{2}}\right)\right]$$

$$-b\left({\displaystyle \frac{\tilde{\varphi }(x_i,y_{j+1},t^n)+\tilde{\varphi }(x_i,y_j,t^n)}{2}}\right)$$

$$*\left[\left({\displaystyle \frac{\varphi (x_i,y_{j+1},t^{n+1})+\varphi (x_i,y_{j+1},t^n)}{2}}\right)-\left({\displaystyle \frac{\varphi (x_i,y_j,t^{n+1})+\varphi (x_i,y_j,t^n)}{2}}\right)\right]$$

$$+b\left({\displaystyle \frac{\tilde{\varphi }(x_i,y_j,t^n)+\tilde{\varphi }(x_i,y_{j-1},t^n)}{2}}\right)$$

$$\left.*\left[\left({\displaystyle \frac{\varphi (x_i,y_j,t^{n+1})+\varphi (x_i,y_j,t^n)}{2}}\right)-\left({\displaystyle \frac{\varphi (x_i,y_{j-1},t^{n+1})+\varphi (x_i,y_{j-1},t^n)}{2}}\right)\right]\right\}$$

$$+{\displaystyle \frac{\Delta t}{4}}{\delta }_x\left(a_{i-{\displaystyle \frac{1}{2}},j}^{n+{\displaystyle \frac{1}{2}}}\left({\delta }_{\overline{x}}\left({\delta }_y\left(b_{i,j-{\displaystyle \frac{1}{2}}}^{n+{\displaystyle \frac{1}{2}}}\left({\delta }_{\overline{y}}\left(\varphi (x_i,y_j,t^{n+1})-\varphi (x_i,y_j,t^n)\right)\right)\right)\right)\right)\right)$$

$$-f(\tilde{\varphi }(x_i,y_j,t^n),\tilde{u}(x_i,y_j,t^n)),$$

(59)

where $\tilde{\varphi }(x_i,y_{j-1},t^n):={\displaystyle \frac{3}{2}}\varphi (x_i,y_{j-1},t^n)-{\displaystyle \frac{1}{2}}\varphi (x_i,y_{j-1},t^{n-1})$ and $\tilde{u}(x_i,y_{j-1},t^n):={\displaystyle \frac{3}{2}}u(x_i,y_{j-1},t^n)-{\displaystyle \frac{1}{2}}u(x_i,y_{j-1},t^{n-1}).$ Then by Taylor’s theorem we have the following local error estimate, whose proof, long but straightforward, is omitted.

Lemma 10.

Suppose that the solution $(u,\varphi )$ of (43) is sufficient smooth. Then, there exist constants $C_{\epsilon },C_{\zeta }$ independent of $\Delta t,$h such that

$$\underset{\begin{array}{c}c1\le n\le N-1\hfill \\ (x_i,y_j)\in {\Omega }_h\hfill \end{array}}{max}\left|{\epsilon }_{ij}^n\right|\le C_{\epsilon }(\Delta t^2+h^2),$$

(60)

$$\underset{\begin{array}{c}c1\le n\le N-1\hfill \\ (x_i,y_j)\in {\Omega }_h\hfill \end{array}}{max}\left|{\zeta }_{ij}^n\right|\le C_{\zeta }(\Delta t^2+h^2).$$

(61)

8.2. Error Estimate

We now prove the main result of paper.

Theorem 2.

Suppose that the solution $(u,\varphi )$ of (43) is sufficient smooth. Let $(U^n,{\Phi }^n)$ be the solution of (48). Suppose that ${\displaystyle \frac{\Delta t}{h}}\le \sigma ,$ where σ is sufficiently small and that $\mathrm{\parallel }{z}^1{\mathrm{\parallel }}_h=O(\Delta t^2+h^2),$ where $z^1={\varphi }^1-{\Phi }^1.$ Then, there exists a constant C independent of $\Delta t$ and h such that

$$\underset{0\le n\le N}{max}{\mathrm{\parallel }{u}^n-U^n\mathrm{\parallel }}_h\le C(\Delta t^2+h^2),$$

(62)

$$\underset{0\le n\le N}{max}{|{\varphi }^n-{\Phi }^n|}_{1,h}\le C(\Delta t^2+h^2).$$

(63)

Proof:

In the sequel C will denote generic constants independent of $\Delta t$ and $h.$ At various points in the proof we follow some techniques of Dendy, [5], extending them from the case of a single parabolic equation considered in [5] to the case of the parabolic system considered herein.

Let $e^n:=u^n-U^n,z^n:={\varphi }^n-{\Phi }^n.$ From (59) and (52) we obtain, for $1\le n\le N-1,$ and $1\le i,j\le J:$

$$\left.\begin{array}{c}{\displaystyle \frac{z_{ij}^{n+1}-z_{ij}^n}{\Delta t}}\hfill \\ -{\delta }_x\left(a(x_{i-{\displaystyle \frac{1}{2}}},y_j,t^{n+{\displaystyle \frac{1}{2}}},{\displaystyle \frac{\tilde{\varphi }(x_i,y_j,t^n)+\tilde{\varphi }(x_{i-1},y_j,t^n)}{2}}){\delta }_{\overline{x}}\left({\displaystyle \frac{{\varphi }_{ij}^{n+1}+{\varphi }_{ij}^n}{2}}\right)\right)\hfill \\ -{\delta }_y\left(b(x_i,y_{j-{\displaystyle \frac{1}{2}}},t^{n+{\displaystyle \frac{1}{2}}},{\displaystyle \frac{\tilde{\varphi }(x_i,y_j,t^n)+\tilde{\varphi }(x_i,y_{j-1},t^n)}{2}}){\delta }_{\overline{y}}\left({\displaystyle \frac{{\varphi }_{ij}^{n+1}+{\varphi }_{ij}^n}{2}}\right)\right)\hfill \\ +{\delta }_x\left(a(x_{i-{\displaystyle \frac{1}{2}}},y_j,t^{n+{\displaystyle \frac{1}{2}}},{\displaystyle \frac{\tilde{\Phi }(x_i,y_j,t^n)+\tilde{\Phi }(x_{i-1},y_j,t^n)}{2}}){\delta }_{\overline{x}}\left({\displaystyle \frac{{\Phi }_{ij}^{n+1}+{\Phi }_{ij}^n}{2}}\right)\right)\hfill \\ +{\delta }_y\left(b(x_i,y_{j-{\displaystyle \frac{1}{2}}},t^{n+{\displaystyle \frac{1}{2}}},{\displaystyle \frac{\tilde{\Phi }(x_i,y_j,t^n)+\tilde{\Phi }(x_i,y_{j-1},t^n)}{2}}){\delta }_{\overline{y}}\left({\displaystyle \frac{{\Phi }_{ij}^{n+1}+{\Phi }_{ij}^n}{2}}\right)\right)\hfill \\ +{\displaystyle \frac{\Delta t}{4}}{\delta }_x\left(a(x_{i-{\displaystyle \frac{1}{2}}},y_j,t^{n+{\displaystyle \frac{1}{2}}},{\displaystyle \frac{\tilde{\varphi }(x_i,y_j,t^n)+\tilde{\varphi }(x_{i-1},y_j,t^n)}{2}})\right.\hfill \\ \left.*\left({\delta }_{\overline{x}}\left({\delta }_y\left(b\left({\displaystyle \frac{\tilde{\varphi }(x_i,y_j,t^n)+\tilde{\varphi }(x_i,y_{j-1},t^n)}{2}}\right)\left({\delta }_{\overline{y}}({\varphi }_{ij}^{n+1}-{\varphi }_{ij}^n)\right)\right)\right)\right)\right)\hfill \\ -{\displaystyle \frac{\Delta t}{4}}{\delta }_x\left(a(x_{i-{\displaystyle \frac{1}{2}}},y_j,t^{n+{\displaystyle \frac{1}{2}}},{\displaystyle \frac{\tilde{\Phi }(x_i,y_j,t^n)+\tilde{\Phi }(x_{i-1},y_j,t^n)}{2}})\right.\hfill \\ \left.*\left({\delta }_{\overline{x}}\left({\delta }_y\left(b\left({\displaystyle \frac{\tilde{\Phi }(x_i,y_j,t^n)+\tilde{\Phi }(x_i,y_{j-1},t^n)}{2}}\right)\left({\delta }_{\overline{y}}({\Phi }_{ij}^{n+1}-{\Phi }_{ij}^n)\right)\right)\right)\right)\right)\hfill \\ =f(x_i,y_j,t^{n+{\displaystyle \frac{1}{2}}},\tilde{\varphi }(x_i,y_j,t^n),\tilde{u}(x_i,y_j,t^n))\hfill \\ -f(x_i,y_j,t^{n+{\displaystyle \frac{1}{2}}},\tilde{\Phi }(x_i,y_j,t^n),\tilde{U}(x_i,y_j,t^n))+{\zeta }_{ij}^n\hfill \end{array}\right\}$$

(64)

For $n\ge 1$ we define ${\varphi }_{i,j}^n:=\varphi (x_i,y_j,t^n)$, ${\tilde{\varphi }}_{ij}^n:={\displaystyle \frac{3}{2}}{\varphi }_{ij}^n-{\displaystyle \frac{1}{2}}{\varphi }_{ij}^{n-1},$ $a_{i-{\displaystyle \frac{1}{2}},j}^n\left(\widehat{\varphi }\right):=a\left({\displaystyle \frac{{\tilde{\varphi }}_{ij}^n+{\tilde{\varphi }}_{i-1,j}^n}{2}}\right)=a\left({\widehat{\varphi }}_{i-{\displaystyle \frac{1}{2}},j}^n\right),$ $i=1,\dots ,J+1,j=0,\dots ,J+1.$ $b_{i,j-{\displaystyle \frac{1}{2}}}^n\left(\widehat{\varphi }\right):=b\left({\displaystyle \frac{{\tilde{\varphi }}_{ij}^n+{\tilde{\varphi }}_{i,j-1}^n}{2}}\right)=b\left({\widehat{\varphi }}_{i,j-{\displaystyle \frac{1}{2}}}^n\right),$ $i=0,\dots ,J+1,j=1,\dots ,J+1,$ where

$${\widehat{\varphi }}_{i-{\displaystyle \frac{1}{2}},j}^n={\displaystyle \frac{1}{2}}\left({\displaystyle \frac{3}{2}}{\varphi }_{ij}^n-{\displaystyle \frac{1}{2}}{\varphi }_{ij}^{n-1}+{\displaystyle \frac{3}{2}}{\varphi }_{i-1,j}^n-{\displaystyle \frac{1}{2}}{\varphi }_{i-1,j}^{n-1}\right)$$

$$={\displaystyle \frac{3}{4}}\left({\varphi }_{ij}^n+{\varphi }_{i-1,j}^n\right)-{\displaystyle \frac{1}{4}}\left({\varphi }_{ij}^{n-1}+{\varphi }_{i-1,j}^{n-1}\right),$$

$${\widehat{\varphi }}_{i,j-{\displaystyle \frac{1}{2}}}^n={\displaystyle \frac{1}{2}}\left({\displaystyle \frac{3}{2}}{\varphi }_{ij}^n-{\displaystyle \frac{1}{2}}{\varphi }_{ij}^{n-1}+{\displaystyle \frac{3}{2}}{\varphi }_{i,j-1}^n-{\displaystyle \frac{1}{2}}{\varphi }_{i,j-1}^{n-1}\right)$$

$$={\displaystyle \frac{3}{4}}\left({\varphi }_{ij}^n+{\varphi }_{i,j-1}^n\right)-{\displaystyle \frac{1}{4}}\left({\varphi }_{ij}^{n-1}+{\varphi }_{i,j-1}^{n-1}\right),$$

We define similarly the discrete quantities $\tilde{\Phi },$$\widehat{\Phi },$$a_{i-{\displaystyle \frac{1}{2}},j}^n\left(\widehat{\Phi }\right)$ and $b_{i,j-{\displaystyle \frac{1}{2}}}^n\left(\widehat{\Phi }\right)$ and we use the notation

$$F_{ij}^n(\tilde{\varphi },\tilde{u})-F_{ij}^n(\tilde{\Phi },\tilde{U}):=f({\tilde{\varphi }}_{ij}^n,{\tilde{u}}_{ij}^n)-f({\tilde{\Phi }}_{ij}^n,{\tilde{U}}_{ij}^n).$$

We finally define for $1\le \nu \le N$

$$\mathrm{\parallel }{z}^{\nu }{\mathrm{\parallel }}_n^2:=h^2\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}\left(a_{i-{\displaystyle \frac{1}{2}},j}^n\left({\widehat{\Phi }}^n\right){\left({\delta }_{\overline{x}}{z}_{ij}^{\nu }\right)}^2+b_{i,j-{\displaystyle \frac{1}{2}}}^n\left({\widehat{\Phi }}^n\right){\left({\delta }_{\overline{y}}{z}_{ij}^{\nu }\right)}^2\right).$$

(65)

The proof proceeds inductively. Our induction hypothesis is: Suppose that the following inequality holds for $1\le m\le M:$

$${\displaystyle \frac{1}{\Delta t}}\sum _{n=1}^{m-1}\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2+\mathrm{\parallel }{z}^m{\mathrm{\parallel }}_{m-1}^2+\Delta t{\displaystyle \frac{\lambda }{4}}\sum _{n=1}^{m-1}{|z^{n+1}-z^n|}_{\overline{x}\overline{y},h}^2$$

$$\le C_2\Delta t\sum _{n=1}^{m-1}{\mathrm{\parallel }{z}^n\mathrm{\parallel }}_{n-1}^2+C_2(\Delta t^4+h^4).\left(\mathrm{I}\right)$$

Obviously, $\left(I\right)$ holds for $m=1.$

Taking the ${\ell }^2$ inner product of both sides of (64) with $z^{n+1}-z^n$ we have for $n\ge 2$

$$\left.\begin{array}{c}\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2\hfill \\ -\Delta t{\left({\delta }_x\left(a^n\left(\widehat{\Phi }\right){\delta }_{\overline{x}}\left({\displaystyle \frac{z^{n+1}+z^n}{2}}\right)\right)+{\delta }_y\left(b^n\left(\widehat{\Phi }\right){\delta }_{\overline{y}}\left({\displaystyle \frac{z^{n+1}+z^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ +{\displaystyle \frac{\Delta t^2}{4}}{\left({\delta }_x\left(a^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{y}}(z^{n+1}-z^n)\right)\right)\right)\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ =\Delta t{\left(F^n(\tilde{\varphi },\tilde{u})-F^n(\tilde{\Phi },\tilde{U}),z^{n+1}-z^n\right)}_h+\Delta t{\left({\zeta }^n,z^{n+1}-z^n\right)}_h\hfill \\ +\Delta t{\left({\delta }_x\left([a^n\left(\widehat{\varphi }\right)-a^n\left(\widehat{\Phi }\right)]{\delta }_{\overline{x}}\left({\displaystyle \frac{{\varphi }^{n+1}+{\varphi }^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ +\Delta t{\left({\delta }_y\left([b^n\left(\widehat{\varphi }\right)-b^n\left(\widehat{\Phi }\right)]{\delta }_{\overline{y}}\left({\displaystyle \frac{{\varphi }^{n+1}+{\varphi }^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ -{\displaystyle \frac{\Delta t^2}{4}}{\left({\delta }_x\left(a^n\left(\widehat{\varphi }\right)\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left(\widehat{\varphi }\right)\left({\delta }_{\overline{y}}({\varphi }^{n+1}-{\varphi }^n)\right)\right)\right)\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ +{\displaystyle \frac{\Delta t^2}{4}}{\left({\delta }_x\left(a^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{y}}({\varphi }^{n+1}-{\varphi }^n)\right)\right)\right)\right)\right),z^{n+1}-z^n\right)}_h.\hfill \end{array}\right\}$$

(66)

By Lemma 8, the symmetry of $L_h^n,$ and our hypotheses on $a,b,$ we have, using the equivalence of the norms ${\mathrm{\parallel }·\mathrm{\parallel }}_n$ and ${|·|}_{1,h}$ with constants independent of h and $\Delta t,$ and the notation ${\widehat{z}}^n:={\widehat{\varphi }}^n-{\widehat{\Phi }}^n:$

$$\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_n^2:=h^2\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}\left(a_{i-{\displaystyle \frac{1}{2}},j}^n\left({\widehat{\Phi }}^n\right){\left({\delta }_{\overline{x}}{z}_{ij}^n\right)}^2+b_{i,j-{\displaystyle \frac{1}{2}}}^n\left({\widehat{\Phi }}^n\right){\left({\delta }_{\overline{y}}{z}_{ij}^n\right)}^2\right)\le$$

$$\le \mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+C{\left(|{\widehat{z}}^n-{\widehat{z}}^{n-1}|{\delta }_{\overline{x}}{z}^n,{\delta }_{\overline{x}}{z}^n\right)}_h$$

$$+C{\left(|{\widehat{z}}^n-{\widehat{z}}^{n-1}|{\delta }_{\overline{y}}{z}^n,{\delta }_{\overline{y}}{z}^n\right)}_h+C\Delta t{\mathrm{\parallel }{z}^n\mathrm{\parallel }}_{n-1}^2.$$

Therefore

$$-{\left({\delta }_x\left(a^n\left(\widehat{\Phi }\right){\delta }_{\overline{x}}\left({\displaystyle \frac{z^{n+1}+z^n}{2}}\right)\right)+{\delta }_y\left(b^n\left(\widehat{\Phi }\right){\delta }_{\overline{y}}\left({\displaystyle \frac{z^{n+1}+z^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h\ge$$

$$\ge {\displaystyle \frac{1}{2}}\mathrm{\parallel }{z}^{n+1}{\mathrm{\parallel }}_n^2-{\displaystyle \frac{1}{2}}[\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+C{\left(|{\widehat{z}}^n-{\widehat{z}}^{n-1}|{\delta }_{\overline{x}}{z}^n,{\delta }_{\overline{x}}{z}^n\right)}_h$$

$$+C{\left(|{\widehat{z}}^n-{\widehat{z}}^{n-1}|{\delta }_{\overline{y}}{z}^n,{\delta }_{\overline{y}}{z}^n\right)}_h+C^{\prime }\Delta t\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2].$$

(67)

In addition, there holds that

$${\left(|{\widehat{z}}^n-{\widehat{z}}^{n-1}|{\delta }_{\overline{x}}{z}^n,{\delta }_{\overline{x}}{z}^n\right)}_h\le {\displaystyle \frac{C}{h^{1/3}}}\left(\right|{\widehat{z}}^n{|}_{1,h}+|{\widehat{z}}^{n-1}{|}_{1,h})\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2.$$

(68)

To see this consider Hölder’s inequality:

$$\left|\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{f}_{ij}{g}_{ij}\right|\le {\left(\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{\left|f_{ij}\right|}^p\right)}^{{\displaystyle \frac{1}{p}}}{\left(\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{\left|g_{ij}\right|}^q\right)}^{{\displaystyle \frac{1}{q}}}for{\displaystyle \frac{1}{p}}+{\displaystyle \frac{1}{q}}=1,$$

(69)

and

$$\left|\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{f}_{ij}{g}_{ij}{r}_{ij}\right|\le {\left(\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{\left|f_{ij}\right|}^p\right)}^{{\displaystyle \frac{1}{p}}}{\left(\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{\left|g_{ij}\right|}^q\right)}^{{\displaystyle \frac{1}{q}}}{\left(\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{\left|r_{ij}\right|}^s\right)}^{{\displaystyle \frac{1}{s}}}$$

(70)

$for{\displaystyle \frac{1}{p}}+{\displaystyle \frac{1}{q}}+{\displaystyle \frac{1}{s}}=1.$

Applying (70) for $p=6,q={\displaystyle \frac{12}{5}},ands={\displaystyle \frac{12}{5}},$ we have,

$${\left(|{\widehat{z}}^n-{\widehat{z}}^{n-1}|{\delta }_{\overline{x}}{z}^n,{\delta }_{\overline{x}}{z}^n\right)}_h\le \mathrm{\parallel }{\widehat{z}}^n-{\widehat{z}}^{n-1}{\mathrm{\parallel }}_{6,h}\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n{\mathrm{\parallel }}_{{\displaystyle \frac{12}{5}},h}{\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n\mathrm{\parallel }}_{{\displaystyle \frac{12}{5}},h},$$

where for $v\in {\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}$ we define ${\displaystyle {\mathrm{\parallel }v\mathrm{\parallel }}_{p,h}:=(h^2\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}|v_{ij}{{|}^p)}^{1/p}}$ (with ${\mathrm{\parallel }v\mathrm{\parallel }}_{2,h}\equiv {\mathrm{\parallel }v\mathrm{\parallel }}_h$).

Therefore

$${\left(|{\widehat{z}}^n-{\widehat{z}}^{n-1}|{\delta }_{\overline{x}}{z}^n,{\delta }_{\overline{x}}{z}^n\right)}_h\le \mathrm{\parallel }{\widehat{z}}^n-{\widehat{z}}^{n-1}{\mathrm{\parallel }}_{6,h}{\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n\mathrm{\parallel }}_{{\displaystyle \frac{12}{5}},h}^2.$$

From (69) we have

$$\left|\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{f}_{ij}^r{g}_{ij}^s\right|\le {\left(\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{\left|f_{ij}\right|}^{rp}\right)}^{{\displaystyle \frac{1}{p}}}{\left(\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}{\left|g_{ij}\right|}^{sq}\right)}^{{\displaystyle \frac{1}{q}}}for{\displaystyle \frac{1}{p}}+{\displaystyle \frac{1}{q}}=1.$$

(71)

Using (71) with $r={\displaystyle \frac{9}{5}},s={\displaystyle \frac{3}{5}},p={\displaystyle \frac{10}{9}},q=10$ we see that

$$\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n{\mathrm{\parallel }}_{{\displaystyle \frac{12}{5}},h}\le \mathrm{\parallel }{\delta }_{\overline{x}}{z}^n{\mathrm{\parallel }}^{3/4}{\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n\mathrm{\parallel }}_{6,h}^{1/4}.$$

Note now that the inverse inequality ${\mathrm{\parallel }v\mathrm{\parallel }}_{p,h}\le C{h}^{-N({\displaystyle \frac{1}{2}}-{\displaystyle \frac{1}{p}})}{\mathrm{\parallel }v\mathrm{\parallel }}_h$ holds on ${\mathbf{X}}_{\mathbf{h}}^{\mathbf{0}}.$ For $p=6$ we have

$$\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n{\mathrm{\parallel }}_{6,h}\le C{h}^{-2({\displaystyle \frac{1}{2}}-{\displaystyle \frac{1}{6}})}\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n{\mathrm{\parallel }}_h=C{h}^{-{\displaystyle \frac{2}{3}}}{\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n\mathrm{\parallel }}_h.$$

Hence

$$\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n{\mathrm{\parallel }}_{{\displaystyle \frac{12}{5}},h}\le C\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n{\mathrm{\parallel }}_h^{3/4}(h^{-{\displaystyle \frac{2}{3}}}\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n{\mathrm{\parallel }}_h{)}^{1/4}=C{h}^{-{\displaystyle \frac{2}{12}}}\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n{\mathrm{\parallel }}_h=C{h}^{-{\displaystyle \frac{1}{6}}}{\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n\mathrm{\parallel }}_h.$$

By the above we condude that

$$\left(|{\widehat{z}}^n-{\widehat{z}}^{n-1}|{\delta }_{\overline{x}}{z}^n,{\delta }_{\overline{x}}{z}^n\right)\le \mathrm{\parallel }{\widehat{z}}^n-{\widehat{z}}^{n-1}{\mathrm{\parallel }}_{6,h}{\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n\mathrm{\parallel }}_{{\displaystyle \frac{12}{5}},h}^2$$

$$\le C|{\widehat{z}}^n-{\widehat{z}}^{n-1}{|}_{1,h}{\left(C{h}^{-{\displaystyle \frac{1}{6}}}{\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n\mathrm{\parallel }}_h\right)}^2=C{h}^{-{\displaystyle \frac{1}{3}}}|{\widehat{z}}^n-{\widehat{z}}^{n-1}{|}_{1,h}{\mathrm{\parallel }{\delta }_{\overline{x}}{z}^n\mathrm{\parallel }}_h^2$$

$$\le C{h}^{-{\displaystyle \frac{1}{3}}}|{\widehat{z}}^n-{\widehat{z}}^{n-1}{|}_{1,h}{\mathrm{\parallel }{z}^n\mathrm{\parallel }}_{n-1}^2,$$

from which (68) follows. Now, using the induction hypotheses $\left(I\right)$ we will show that $|{\widehat{z}}^n{|}_{1,h},{\left|{\widehat{z}}^{n-1}\right|}_{1,h}$ and

${\displaystyle max\left\{\underset{i,j}{max}\right|{\delta }_{\overline{x}}{\widehat{\Phi }}_{ij}^n|,\underset{i,j}{max}|{\delta }_{\overline{y}}{\widehat{\Phi }}_{ij}^n\left|\right\}}$ are bounded for $n\le M.$

Indeed, from $\left(I\right)$ we have for a small enagh $ϵ>0$ (that we consider fixed from now on) that

$$\mathrm{\parallel }{z}^m{\mathrm{\parallel }}_{m-1}^2\le C\Delta t\sum _{n=1}^m{\mathrm{\parallel }{z}^n\mathrm{\parallel }}_{n-1}^2+C\Delta t^4,m\le M.$$

Hence, from Gronwall’s lemma we obtain

$$|z^m{|}_{1,h}^2\le e^{Cm\Delta t}\Delta t^{4}form\le M.$$

(72)

Hence, using the stability condition ${\displaystyle \frac{\Delta t}{h}}\le \sigma ),$ we see that

$$\sum _{i=1}^{J+1}\sum _{j=1}^{J+1}|\left[{\left({\delta }_{\overline{x}}{\widehat{z}}_{ij}^m\right)}^2+{\left({\delta }_{\overline{y}}{\widehat{z}}_{ij}^m\right)}^2\right]\le C{\displaystyle \frac{\Delta t^4}{h^2}}\le C\Delta t^2.$$

Hence

$$max\left(\underset{i,j}{max}|{\delta }_{\overline{x}}{\widehat{z}}_{ij}|,\underset{i,j}{max}\left|{\delta }_{\overline{y}}{\widehat{z}}_{ij}\right|\right)\le \sqrt{{\sum }_{i=1}^{J+1}{\sum }_{j=1}^{J+1}\left[{\left({\delta }_{\overline{x}}{\widehat{z}}_{ij}^m\right)}^2+{\left({\delta }_{\overline{y}}{\widehat{z}}_{ij}^m\right)}^2\right]}\le C\Delta t,$$

and therefore

$$|{\delta }_{\overline{x}}{\widehat{\Phi }}_{ij}^n|\le |{\delta }_{\overline{x}}{\widehat{z}}_{ij}^n|+|{\delta }_{\overline{x}}{\widehat{\varphi }}_{ij}^n|\le C,$$

(73)

since $|{\delta }_{\overline{x}}{\widehat{\varphi }}_{ij}^n|$ is bounded, as it is easily seen. Similarly we see that $|{\delta }_{\overline{y}}{\widehat{\Phi }}_{ij}^n|\le C.$

Now, from (72), supposing that $n\le M,$ we have, by the stability conttion that

$${\displaystyle |{\widehat{z}}^n{|}_{1,h}+{\left|{\widehat{z}}^{n-1}\right|}_{1,h}\le e^{CT}\Delta t^2\le C{e}^{CT}\Delta th\le C{e}^{CT}{h}^{\frac{2}{3}}\Delta t{h}^{\frac{1}{3}}\le C\Delta t{h}^{\frac{1}{3}}.}$$

Therefore from (68) we have, for $n\le M$

$$\left(|{\widehat{z}}^n-{\widehat{z}}^{n-1}|{\delta }_{\overline{x}}{z}^n,{\delta }_{\overline{x}}{z}^n\right)\le C\Delta t{\mathrm{\parallel }{z}^n\mathrm{\parallel }}_{n-1}^2,$$

(74)

and a similar bound for the term $\left(|{\widehat{z}}^n-{\widehat{z}}^{n-1}|{\delta }_{\overline{y}}{z}^n,{\delta }_{\overline{y}}{z}^n\right).$ Substituting (74) in (67) we have, for $n\le M$

$$-{\left({\delta }_x\left(a^n\left(\widehat{\Phi }\right){\delta }_{\overline{x}}\left({\displaystyle \frac{z^{n+1}+z^n}{2}}\right)\right)+{\delta }_y\left(b^n\left(\widehat{\Phi }\right){\delta }_{\overline{y}}\left({\displaystyle \frac{z^{n+1}+z^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h$$

$$\ge {\displaystyle \frac{1}{2}}\mathrm{\parallel }{z}^{n+1}{\mathrm{\parallel }}_n^2-{\displaystyle \frac{1}{2}}(\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+C\Delta t\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2)$$

(75)

Hence, from (66) and (76) we obtain, for $n\le M$

$\begin{array}{c}\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2\hfill \\ +\Delta t\left[{\displaystyle \frac{1}{2}}\mathrm{\parallel }{z}^{n+1}{\mathrm{\parallel }}_n^2-{\displaystyle \frac{1}{2}}(\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+C\Delta t\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2)\right]\hfill \\ +{\displaystyle \frac{\Delta t^2}{4}}{\left({\delta }_x\left(a^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{y}}(z^{n+1}-z^n)\right)\right)\right)\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ \le \Delta t{\left(F^n(\tilde{\varphi },\tilde{u})-F^n(\tilde{\Phi },\tilde{U}),z^{n+1}-z^n\right)}_h+\Delta t{\left({\zeta }^n,z^{n+1}-z^n\right)}_h\hfill \\ +\Delta t{\left({\delta }_x\left([a^n\left(\widehat{\varphi }\right)-a^n\left(\widehat{\Phi }\right)]{\delta }_{\overline{x}}\left({\displaystyle \frac{{\varphi }^{n+1}+{\varphi }^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ +\Delta t{\left({\delta }_y\left([b^n\left(\widehat{\varphi }\right)-b^n\left(\widehat{\Phi }\right)]{\delta }_{\overline{y}}\left({\displaystyle \frac{{\varphi }^{n+1}+{\varphi }^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ -{\displaystyle \frac{\Delta t^2}{4}}{\left({\delta }_x\left(a^n\left(\widehat{\varphi }\right)\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left(\widehat{\varphi }\right)\left({\delta }_{\overline{y}}({\varphi }^{n+1}-{\varphi }^n)\right)\right)\right)\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ +{\displaystyle \frac{\Delta t^2}{4}}{\left({\delta }_x\left(a^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{y}}({\varphi }^{n+1}-{\varphi }^n)\right)\right)\right)\right)\right),z^{n+1}-z^n\right)}_h.\hfill \end{array}$

Therefore, using Lemma 9, we conclude for $n\le M$ that

$$\left.\begin{array}{c}\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2+\Delta t\left[{\displaystyle \frac{1}{2}}\mathrm{\parallel }{z}^{n+1}{\mathrm{\parallel }}_n^2-{\displaystyle \frac{1}{2}}(\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+C\Delta t\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2)\right]\hfill \\ +{\displaystyle \frac{\Delta t^2}{4}}\left(\lambda \right|z^{n+1}-z^n{|}_{\overline{x}\overline{y},h}^2-\mu |z^{n+1}-z^n{|}_{1,h}^2)\le \hfill \\ \le \Delta t{\left(F^n(\tilde{\varphi },\tilde{u})-F^n(\tilde{\Phi },\tilde{U}),z^{n+1}-z^n\right)}_h+\Delta t{\left({\zeta }^n,z^{n+1}-z^n\right)}_h\hfill \\ +\Delta t{\left({\delta }_x\left([a^n\left(\widehat{\varphi }\right)-a^n\left(\widehat{\Phi }\right)]{\delta }_{\overline{x}}\left({\displaystyle \frac{{\varphi }^{n+1}+{\varphi }^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ +\Delta t{\left({\delta }_y\left([b^n\left(\widehat{\varphi }\right)-b^n\left(\widehat{\Phi }\right)]{\delta }_{\overline{y}}\left({\displaystyle \frac{{\varphi }^{n+1}+{\varphi }^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ -{\displaystyle \frac{\Delta t^2}{4}}{\left({\delta }_x\left(a^n\left(\widehat{\varphi }\right)\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left(\widehat{\varphi }\right)\left({\delta }_{\overline{y}}({\varphi }^{n+1}-{\varphi }^n)\right)\right)\right)\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ +{\displaystyle \frac{\Delta t^2}{4}}{\left({\delta }_x\left(a^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{y}}({\varphi }^{n+1}-{\varphi }^n)\right)\right)\right)\right)\right),z^{n+1}-z^n\right)}_h.\hfill \end{array}\right\}$$

(76)

We now define the following quantities

$\begin{array}{ccc}I& :=& {\left({\delta }_x\left([a^n\left(\widehat{\varphi }\right)-a^n\left(\widehat{\Phi }\right)]{\delta }_{\overline{x}}\left({\displaystyle \frac{{\varphi }^{n+1}+{\varphi }^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h\hfill \\ & & +{\left({\delta }_y\left([b^n\left(\widehat{\varphi }\right)-b^n\left(\widehat{\Phi }\right)]{\delta }_{\overline{y}}\left({\displaystyle \frac{{\varphi }^{n+1}+{\varphi }^n}{2}}\right)\right),z^{n+1}-z^n\right)}_h,\hfill \\ II& :=& {\left({\delta }_x\left(a^n\left(\widehat{\varphi }\right)\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left(\widehat{\varphi }\right)\left({\delta }_{\overline{y}}({\varphi }^{n+1}-{\varphi }^n)\right)\right)\right)\right)\right),z^{n+1}-z^n\right)}_h,\hfill \\ III& :=& {\left({\delta }_x\left(a^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{x}}\left({\delta }_y\left(b^n\left(\widehat{\Phi }\right)\left({\delta }_{\overline{y}}({\varphi }^{n+1}-{\varphi }^n)\right)\right)\right)\right)\right),z^{n+1}-z^n\right)}_h.\hfill \end{array}$

From (72), since $n\le M,$ we have $\mathrm{\parallel }{\widehat{z}}^n{\mathrm{\parallel }}_h\le C{\left|{\widehat{z}}^n\right|}_{1,h}\le e^{CT}\Delta t^2.$

Long but straightforward calculations yield now for $\epsilon >0$ that

$$\Delta t\left|I\right|\le C_{\epsilon }\Delta t^6+\epsilon {\mathrm{\parallel }{z}^{n+1}-z^n\mathrm{\parallel }}_h^2,$$

(77)

and

$$\Delta t^2|III-II|\le C_{\epsilon }\Delta t^6+\epsilon \Delta t^2{|z^{n+1}-z^n|}_{\overline{x},\overline{y},h}^2.$$

(78)

Substituting in (76) we have

$$\begin{array}{c}\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2+\Delta t\left[{\displaystyle \frac{1}{2}}\mathrm{\parallel }{z}^{n+1}{\mathrm{\parallel }}_n^2-{\displaystyle \frac{1}{2}}(\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+C\Delta t\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2)\right]\hfill \\ +{\displaystyle \frac{\Delta t^2}{4}}\left(\lambda \right|z^{n+1}-z^n{|}_{\overline{x}\overline{y},h}^2-\mu |z^{n+1}-z^n{|}_{1,h}^2)\le \hfill \\ \le \Delta t{\left(F^n(\tilde{\varphi },\tilde{u})-F^n(\tilde{\Phi },\tilde{U}),z^{n+1}-z^n\right)}_h+\Delta t{\left({\zeta }^n,z^{n+1}-z^n\right)}_h\hfill \\ +C_{{\epsilon }^{\prime }}\Delta t^6+{\epsilon }^{\prime }{\mathrm{\parallel }{z}^{n+1}-z^n\mathrm{\parallel }}_h^2\hfill \\ +C_{{\epsilon }^{\prime \prime }}\Delta t^6+{\epsilon }^{\prime \prime }\Delta t^2{|z^{n+1}-z^n|}_{\overline{x},\overline{y},h}^2,\hfill \end{array}$$

where ${\epsilon }^{\prime },{\epsilon }^{\prime \prime }>0$ arbitrary.

Hence, using the Cauchy-Schwarz inequality and the arithmetic geometric mean inequality we have

$$\begin{array}{c}\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2+\Delta t\left[{\displaystyle \frac{1}{2}}\mathrm{\parallel }{z}^{n+1}{\mathrm{\parallel }}_n^2-{\displaystyle \frac{1}{2}}(\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+C\Delta t\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2)\right]\hfill \\ +{\displaystyle \frac{\Delta t^2}{4}}\left(\lambda \right|z^{n+1}-z^n{|}_{\overline{x}\overline{y},h}^2-\mu |z^{n+1}-z^n{|}_{1,h}^2)\le \hfill \\ \le \left({\displaystyle \frac{\Delta t^2}{4{\epsilon }^{\prime \prime \prime }}}\mathrm{\parallel }{F}^n(\tilde{\varphi },\tilde{u})-F^n(\tilde{\Phi },\tilde{U}){\mathrm{\parallel }}_h^2+{\epsilon }^{\prime \prime \prime }{\mathrm{\parallel }{z}^{n+1}-z^n\mathrm{\parallel }}_h^2\right)\hfill \\ +\left({\displaystyle \frac{\Delta t^2}{4{\epsilon }^{\prime \prime \prime }}}\mathrm{\parallel }{\zeta }^n{\mathrm{\parallel }}_h^2+{\epsilon }^{\prime \prime \prime }{\mathrm{\parallel }{z}^{n+1}-z^n\mathrm{\parallel }}_h^2\right)\hfill \\ +C_{{\epsilon }^{\prime },{\epsilon }^{\prime \prime }}\Delta t^6+{\epsilon }^{\prime }\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2+{\epsilon }^{\prime \prime }\Delta t^2{|z^{n+1}-z^n|}_{\overline{x},\overline{y},h}^2.\hfill \end{array}$$

In addition since $|F_{\varphi }|\le c,|F_u|\le c$ we have from the mean value theorem

$$\mathrm{\parallel }{F}^n(\tilde{\varphi },\tilde{u})-F^n(\tilde{\Phi },\tilde{U}){\mathrm{\parallel }}_h\le L_1\mathrm{\parallel }{\tilde{z}}^n{\mathrm{\parallel }}_h+L_2{\mathrm{\parallel }{\tilde{e}}^n\mathrm{\parallel }}_h$$

$$\le L_1^{\prime }(\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_h+\mathrm{\parallel }{z}^{n-1}{\mathrm{\parallel }}_h)+L_2^{\prime }(\mathrm{\parallel }{e}^n{\mathrm{\parallel }}_h+\mathrm{\parallel }{e}^{n-1}{\mathrm{\parallel }}_h).$$

Since for $n\le M$ (72) gives $\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_h\le C{\left|z^n\right|}_{1,h}\le e^{CT}\Delta t^2$, we have for ${\epsilon }^{\prime \prime \prime }>0$

$$\mathrm{\parallel }{F}^n(\tilde{\varphi },\tilde{u})-F^n(\tilde{\Phi },\tilde{U}){\mathrm{\parallel }}_h^2\le C\Delta t^4+C(\mathrm{\parallel }{e}^n{\mathrm{\parallel }}_h^2+\mathrm{\parallel }{e}^{n-1}{\mathrm{\parallel }}_h^2),$$

$$\begin{array}{c}\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2\hfill \\ +\Delta t\left[{\displaystyle \frac{1}{2}}\mathrm{\parallel }{z}^{n+1}{\mathrm{\parallel }}_n^2-{\displaystyle \frac{1}{2}}(\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+C\Delta t\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2)\right]\hfill \\ +{\displaystyle \frac{\Delta t^2}{4}}\left(\lambda \right|z^{n+1}-z^n{|}_{\overline{x}\overline{y},h}^2-\mu |z^{n+1}-z^n{|}_{1,h}^2)\hfill \\ \le \left({\displaystyle \frac{\Delta t^2}{4{\epsilon }^{\prime \prime \prime }}}(C\Delta t^4+C(\mathrm{\parallel }{e}^n{\mathrm{\parallel }}_h^2+\mathrm{\parallel }{e}^{n-1}{\mathrm{\parallel }}_h^2\left)\right)+{\epsilon }^{\prime \prime \prime }{\mathrm{\parallel }{z}^{n+1}-z^n\mathrm{\parallel }}_h^2\right)\hfill \\ +\left({\displaystyle \frac{\Delta t^2}{4{\epsilon }^{\prime \prime \prime }}}\mathrm{\parallel }{\zeta }^n{\mathrm{\parallel }}_h^2+{\epsilon }^{\prime \prime \prime }{\mathrm{\parallel }{z}^{n+1}-z^n\mathrm{\parallel }}_h^2\right)\hfill \\ +C_{{\epsilon }^{\prime },{\epsilon }^{\prime \prime }}\Delta t^6+{\epsilon }^{\prime }\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2+{\epsilon }^{\prime \prime }\Delta t^2{|z^{n+1}-z^n|}_{\overline{x},\overline{y},h}^2.\hfill \end{array}$$

Then, for suitable ${\epsilon }^{\prime },{\epsilon }^{\prime \prime },{\epsilon }^{\prime \prime \prime }=O\left(\epsilon \right),$ we see that

$$\begin{array}{c}(1-\epsilon )\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2\hfill \\ +{\displaystyle \frac{\Delta t}{2}}\left[\mathrm{\parallel }{z}^{n+1}{\mathrm{\parallel }}_n^2-{\mathrm{\parallel }{z}^n\mathrm{\parallel }}_{n-1}^2\right]+\Delta t^2({\displaystyle \frac{\lambda }{4}}-\epsilon ){|z^{n+1}-z^n|}_{\overline{x}\overline{y},h}^2\hfill \\ \le C\Delta t^2\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+\Delta t^2{C}_{\epsilon }(\mathrm{\parallel }{e}^n{\mathrm{\parallel }}_h^2+\mathrm{\parallel }{e}^{n-1}{\mathrm{\parallel }}_h^2)\hfill \\ +\Delta t^2{C}_{\epsilon }\mathrm{\parallel }{\zeta }^n{\mathrm{\parallel }}_h^2+C_{\epsilon }\Delta t^6+{\displaystyle \frac{\Delta t^2}{4}}\mu {|z^{n+1}-z^n|}_{1,h}^2.\hfill \end{array}$$

Hence, since $|z^{n+1}+z^n{|}_{1,h}\le {\displaystyle \frac{C_{*}}{h}}{\mathrm{\parallel }{z}^{n+1}+z^n\mathrm{\parallel }}_h,$ at ${\displaystyle \frac{\Delta t}{h}}\le \sigma ,$ we have that ${\displaystyle \frac{\Delta t^2}{4}}\mu |z^{n+1}+z^n{|}_{1,h}^2\le C^{\prime }{\sigma }^2{\mathrm{\parallel }{z}^{n+1}+z^n\mathrm{\parallel }}_h$ where $C^{\prime }={\displaystyle \frac{C_{*}\mu }{4}}$

Therefore

$$\begin{array}{c}(1-\epsilon -C^{\prime }{\sigma }^2){\mathrm{\parallel }{z}^{n+1}-z^n\mathrm{\parallel }}_h^2\hfill \\ +{\displaystyle \frac{\Delta t}{2}}\left[\mathrm{\parallel }{z}^{n+1}{\mathrm{\parallel }}_n^2-{\mathrm{\parallel }{z}^n\mathrm{\parallel }}_{n-1}^2\right]+\Delta t^2({\displaystyle \frac{\lambda }{4}}-\epsilon ){|z^{n+1}-z^n|}_{\overline{x}\overline{y},h}^2\hfill \\ \le C\Delta t^2\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+\Delta t^2{C}_{\epsilon }(\mathrm{\parallel }{e}^n{\mathrm{\parallel }}_h^2\hfill \\ +\mathrm{\parallel }{e}^{n-1}{\mathrm{\parallel }}_h^2)+\Delta t^2{C}_{\epsilon }{\mathrm{\parallel }{\zeta }^n\mathrm{\parallel }}_h^2+C_{\epsilon }\Delta t^6,n\le M.\hfill \end{array}$$

Summing from $n=1tom$ we have

$${\displaystyle \frac{2(1-\epsilon -C^{\prime }{\sigma }^2)}{\Delta t}}\sum _{n=1}^m\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2+{\mathrm{\parallel }{z}^{m+1}\mathrm{\parallel }}_m^2$$

$$+2\Delta t({\displaystyle \frac{\lambda }{4}}-\epsilon )\sum _{n=1}^m{|z^{n+1}-z^n|}_{\overline{x}\overline{y},h}^2$$

$$\le C\Delta t\sum _{n=1}^m\mathrm{\parallel }{z}^n{\mathrm{\parallel }}_{n-1}^2+C_{\epsilon }\Delta t\sum _{n=1}^m{\mathrm{\parallel }{e}^n\mathrm{\parallel }}_h^2$$

$$+C_{\epsilon }\Delta t\sum _{n=1}^m{\mathrm{\parallel }{\zeta }^n\mathrm{\parallel }}_h^2+C_{\epsilon }\Delta t^{4}1\le m\le M,\forall \epsilon >0.$$

(79)

Now, from (58)and (50) we have for $0\le n\le N-1,1\le i,j\le J:$

$${\displaystyle \frac{e_{ij}^{n+1}-e_{ij}^n}{\Delta t}}-{\Delta }_h\left({\displaystyle \frac{e_{ij}^{n+1}+e_{ij}^n}{2}}\right)+{\displaystyle \frac{\Delta t}{4}}{\delta }_{x\overline{x}}{\delta }_{y\overline{y}}(e_{ij}^{n+1}-e_{ij}^n)$$

$$={\displaystyle \frac{1}{\Delta t}}\left\{[P(x_i,y_j,t^{n+1},\varphi (x_i,y_j,t^{n+1}))-P(x_i,y_j,t^n,\varphi (x_i,y_j,t^n))]\right\}$$

$$-{\displaystyle \frac{1}{\Delta t}}\left\{[P(x_i,y_j,t^{n+1},\Phi (x_i,y_j,t^{n+1}))-P(x_i,y_j,t^n,\Phi (x_i,y_j,t^n))]\right\}+{\epsilon }_{ij}^n.$$

We define ${\widehat{P}}_{ij}{\left(v\right)}^n:=P(x_i,y_j,t^{n+1},v(x_i,y_j,t^{n+1}))-P(x_i,y_j,t^n,v(x_i,y_j,t^n)).$ Taking the ${\ell }^2$ inner product of both sides of the above relation with $e^{n+1}+e^n$ we have

$${\left(e^{n+1}-e^n,e^{n+1}+e^n\right)}_h-\Delta t{\left({\Delta }_h\left({\displaystyle \frac{e^{n+1}+e^n}{2}}\right),e^{n+1}+e^n\right)}_h$$

$$+{\displaystyle \frac{\Delta t^2}{4}}{\left({\delta }_{x\overline{x}}{\delta }_{y\overline{y}}(e^{n+1}-e^n),e^{n+1}+e^n\right)}_h$$

$$={\left(\widehat{P}{\left(\varphi \right)}^n-\widehat{P}{(\Phi )}^n,e^{n+1}+e^n\right)}_h+\Delta t{\left({\epsilon }^n,e^{n+1}+e^n\right)}_h.$$

(80)

Using the symmetry of ${\beta }_h(·,·),$ Lemma 9 (with $a\equiv 1,b\equiv 1),$ and Lemma 7 and summing the above relation from $n=0$ to $m^{\prime }-1$ we have, since $e^0=0$

$$\mathrm{\parallel }{e}^{m^{\prime }}{\mathrm{\parallel }}_h^2+{\displaystyle \frac{\Delta t}{2}}\sum _{n=0}^{m^{\prime }-1}{|e^{n+1}+e^n|}_{1,h}^2$$

$$\le \sum _{n=0}^{m^{\prime }-1}{\left(\widehat{P}{\left(\varphi \right)}^n-\widehat{P}{(\Phi )}^n,e^{n+1}+e^n\right)}_h+\Delta t\sum _{n=0}^{m^{\prime }-1}{\left({\epsilon }^n,e^{n+1}+e^n\right)}_h$$

$$\le \sum _{n=0}^{m^{\prime }-1}\mathrm{\parallel }\widehat{P}{\left(\varphi \right)}^n-\widehat{P}{(\Phi )}^n{\mathrm{\parallel }}_h\mathrm{\parallel }{e}^{n+1}+e^n{\mathrm{\parallel }}_h+\Delta t\sum _{n=0}^{m^{\prime }-1}\mathrm{\parallel }{\epsilon }^n{\mathrm{\parallel }}_h{\mathrm{\parallel }{e}^{n+1}+e^n\mathrm{\parallel }}_h$$

$$\le {\displaystyle \frac{1}{4\epsilon \Delta t}}\sum _{n=0}^{m^{\prime }-1}\mathrm{\parallel }\widehat{P}{\left(\varphi \right)}^n-\widehat{P}{(\Phi )}^n{\mathrm{\parallel }}_h^2+\epsilon \Delta t\sum _{n=0}^{m^{\prime }-1}{\mathrm{\parallel }{e}^{n+1}+e^n\mathrm{\parallel }}_h^2$$

$$+{\displaystyle \frac{\Delta t}{4\epsilon }}\sum _{n=0}^{m^{\prime }-1}\mathrm{\parallel }{\epsilon }^n{\mathrm{\parallel }}_h^2+\Delta t\epsilon \sum _{n=0}^{m^{\prime }-1}{\mathrm{\parallel }{e}^{n+1}+e^n\mathrm{\parallel }}_h^2$$

$$\le {\displaystyle \frac{1}{4\epsilon \Delta t}}\sum _{n=0}^{m^{\prime }-1}\mathrm{\parallel }\widehat{P}{\left(\varphi \right)}^n-\widehat{P}{(\Phi )}^n{\mathrm{\parallel }}_h^2+C_{\epsilon }\Delta t\sum _{n=0}^{m^{\prime }-1}\mathrm{\parallel }{\epsilon }^n{\mathrm{\parallel }}_h^2+\Delta t\epsilon \sum _{n=0}^{N-1}{\mathrm{\parallel }{e}^{n+1}+e^n\mathrm{\parallel }}_h^2.$$

In view of (79) we have

$$\mathrm{\parallel }{e}^{m^{\prime }}{\mathrm{\parallel }}_h^2\le C\left(C_{\epsilon }\Delta t\sum _{n=0}^{m^{\prime }-1}\mathrm{\parallel }{e}^n{\mathrm{\parallel }}_h^2+C_{\epsilon }\Delta t\sum _{n=1}^{m^{\prime }-1}\mathrm{\parallel }{\zeta }^n{\mathrm{\parallel }}_h^2+C_{\epsilon }\Delta t^4+{\mathrm{\parallel }{z}^1\mathrm{\parallel }}_0^2\right)$$

$$+c\Delta t\sum _{n=0}^{m^{\prime }-1}{\mathrm{\parallel }{\epsilon }^n\mathrm{\parallel }}_h^2.$$

Hence, by Gronwall lemma we infer that

$$\underset{0\le n\le m^{\prime }}{max}{\mathrm{\parallel }{e}^n\mathrm{\parallel }}_h^2\le C\left[\Delta t\sum _{n=1}^{m-1}\mathrm{\parallel }{\zeta }^n{\mathrm{\parallel }}_h^2+\Delta t\sum _{n=0}^{m-1}\mathrm{\parallel }{\epsilon }^n{\mathrm{\parallel }}_h^2+\Delta t^4+{\mathrm{\parallel }{z}^1\mathrm{\parallel }}_0^2\right].$$

Therefore, by (60) and (61)

$$\underset{0\le n\le m}{max}{\mathrm{\parallel }{e}^n\mathrm{\parallel }}_h^2\le C\left[\Delta t\sum _{n=1}^{m-1}{(\Delta t^2+h^2)}^2+\Delta t^4+{\mathrm{\parallel }{z}^1\mathrm{\parallel }}_0^2\right]$$

$$\le C\left(N\Delta t{(\Delta t^2+h^2)}^2+{\mathrm{\parallel }{z}^1\mathrm{\parallel }}_0^2\right)$$

$$\le C\left({(\Delta t^2+h^2)}^2+{\mathrm{\parallel }{z}^1\mathrm{\parallel }}_0^2\right).$$

Hence, using $\mathrm{\parallel }{z}^1{\mathrm{\parallel }}_0=O(\Delta t^2+h^2)$ we have

$$\underset{0\le n\le m}{max}{\mathrm{\parallel }{e}^n\mathrm{\parallel }}_h\le C(\Delta t^2+h^2).$$

(81)

Since

$$\Delta t\sum _{n=0}^m{\mathrm{\parallel }{e}^n\mathrm{\parallel }}_h^2\le C{(\Delta t^2+h^2)}^2,$$

(82)

substituting (82) and (61) in (79) we have

$${\displaystyle \frac{2(1-\epsilon -C^{\prime }{\sigma }^2)}{\Delta t}}\sum _{n=1}^m\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2+\mathrm{\parallel }{z}^{m+1}{\mathrm{\parallel }}_m^2+2\Delta t({\displaystyle \frac{\lambda }{4}}-\epsilon )\sum _{n=1}^m{|z^{n+1}-z^n|}_{\overline{x}\overline{y},h}^2$$

$$\le C\Delta t\sum _{n=1}^m{\mathrm{\parallel }{z}^n\mathrm{\parallel }}_{n-1}^2+C_{\epsilon }(\Delta t^4+h^4)$$

(83)

Take now $\epsilon >0$ so that $2(1-\epsilon -C^{\prime }{\sigma }^2)\ge 1and2({\displaystyle \frac{\lambda }{4}}-\epsilon )\ge {\displaystyle \frac{\lambda }{4}}$ so (83) gives

$${\displaystyle \frac{1}{\Delta t}}\sum _{n=1}^m\mathrm{\parallel }{z}^{n+1}-z^n{\mathrm{\parallel }}_h^2+\mathrm{\parallel }{z}^{m+1}{\mathrm{\parallel }}_m^2+\Delta t{\displaystyle \frac{\lambda }{4}}\sum _{n=1}^m{|z^{n+1}-z^n|}_{\overline{x}\overline{y},h}^2$$

$$\le C\Delta t\sum _{n=1}^m{\mathrm{\parallel }{z}^n\mathrm{\parallel }}_{n-1}^2+C(\Delta t^4+h^4),$$

which is $\left(I\right)$ for m instead of $m-1.$ Hence, the inductive step is complete and $\left(I\right)$ holds for $1\le m\le N.$ We conclude from (81) that (62) holds.

Finally, it follows from $\left(I\right)$ that

$$\mathrm{\parallel }{z}^m{\mathrm{\parallel }}_{m-1}^2\le C\Delta t\sum _{n=1}^{m-1}{\mathrm{\parallel }{z}^n\mathrm{\parallel }}_{n-1}^2+C(\Delta t^4+h^4),1\le m\le N$$

from which, by Gronwall’s lemma it follows that

$$\underset{1\le m\le N}{max}{\mathrm{\parallel }{z}^m\mathrm{\parallel }}_{m-1}\le C(\Delta t^4+h^4),$$

which implies (63).

▪

We may easily obtain from Taylor’s theorem an approximation at the first step satisfying $\mathrm{\parallel }{z}^1{\mathrm{\parallel }}_h=O(h^2+\Delta t^2).$

9. Numerical Verification of the Order of Convergence

In the sequel, we solved the system (43) under Neumann boundary conditions on $\partial \Omega .$ (We obtained similar rates for the problems with Dirichlet boundary conditions.) We took $\Omega =[0,1]\times [0,1],$ $T=0.1,$$a=b=1,$ $f=\varphi (1-\varphi )u/(1+0.25u),$$p={\varphi }^3(10-15\varphi +6{\varphi }^2).$ We added a suitable nonhomogeneous term so that the solution of the associated initial-boundary-value problem with Neumann boundary conditions on $\partial \Omega$ was $(\varphi ,u)=(e^{-t}cos\left(x(x-1)\right)cos\left(y(y-1)\right),e^{-t}cos\left(\pi x\right)cos\left(\pi y\right)).$

Table 2. Errors and order of convergence of the ADI-ADI scheme, T=0.1.

	${\mathrm{\parallel }·\mathrm{\parallel }}_{\mathit{\infty }}$		${\mathrm{\parallel }·\mathrm{\parallel }}_{\mathit{h}}$
J	errors	order	errors	order
50	$\varphi :1.920e-04$	$--$	$8.664e-05$	$--$
	$u:1.047e-04$	$--$	$4.762e-05$	$--$
75	$\varphi :8.708e-05$	$1.95$	$3.851e-05$	$1.99$
	$u:4.725e-05$	$1.96$	$2.119e-05$	$1.99$
112	$\varphi :3.979e-05$	$1.93$	$1.731e-05$	$1.97$
	$u:2.143e-05$	$1.94$	$9.513e-06$	$1.97$
168	$\varphi :1.806e-05$	$1.94$	$7.747e-06$	$1.98$
	$u:9.622e-06$	$1.97$	$4.231e-06$	$1.99$
252	$\varphi :8.254e-06$	$1.93$	$3.486e-06$	$1.96$
	$u:4.320e-06$	$1.97$	$1.883e-06$	$1.99$
378	$\varphi :3.814e-06$	$1.90$	$1.582e-06$	$1.94$
	$u:1.944e-06$	$1.96$	$8.383e-07$	$1.99$
567	$\varphi :1.795e-06$	$1.85$	$7.300e-07$	$1.90$
	$u:8.783e-07$	$1.95$	$3.736e-07$	$1.99$

Table 2. Errors and order of convergence of the ADI-ADI scheme, T=0.1.

	${\mathrm{\parallel }·\mathrm{\parallel }}_{\mathit{\infty }}$		${\mathrm{\parallel }·\mathrm{\parallel }}_{\mathit{h}}$
J	errors	order	errors	order
50	$\varphi :1.920e-04$	$--$	$8.664e-05$	$--$
	$u:1.047e-04$	$--$	$4.762e-05$	$--$
75	$\varphi :8.708e-05$	$1.95$	$3.851e-05$	$1.99$
	$u:4.725e-05$	$1.96$	$2.119e-05$	$1.99$
112	$\varphi :3.979e-05$	$1.93$	$1.731e-05$	$1.97$
	$u:2.143e-05$	$1.94$	$9.513e-06$	$1.97$
168	$\varphi :1.806e-05$	$1.94$	$7.747e-06$	$1.98$
	$u:9.622e-06$	$1.97$	$4.231e-06$	$1.99$
252	$\varphi :8.254e-06$	$1.93$	$3.486e-06$	$1.96$
	$u:4.320e-06$	$1.97$	$1.883e-06$	$1.99$
378	$\varphi :3.814e-06$	$1.90$	$1.582e-06$	$1.94$
	$u:1.944e-06$	$1.96$	$8.383e-07$	$1.99$
567	$\varphi :1.795e-06$	$1.85$	$7.300e-07$	$1.90$
	$u:8.783e-07$	$1.95$	$3.736e-07$	$1.99$

Table 2 shows the errors and the experimental orders of convergence of the numerical solution computed by the ADI-ADI scheme for the this problem, in the discrete maximum norm ${\mathrm{\parallel }·\mathrm{\parallel }}_{\infty }$ and also in the ${\ell }_2$ norm ${\mathrm{\parallel }·\mathrm{\parallel }}_h.$ We computed with $h=1/(J+1),$ where $J=50,75,\dots ,567,$ and $\Delta t=0.01h.$ The observed order of convergence is approximately equal to 2 for both variables in both norms. (For two consecutive runs with values $h_1,h_2$ of $h,$ giving errors $e_1,e_2,$ respectively, the rate of convergence was computed as $log(e_2/e_1)/log(h_2/h_1).$)

With this scheme we performed two numerical experiments where we took (cf. Wang[15]) on $\Omega =[0,1]\times [0,1],$ ${\varphi }_0(x,y)={\displaystyle \frac{1}{2}}\left[1+tanh{\displaystyle \frac{r-R_0}{2\sqrt{2}\epsilon }}\right],$ $u_0(x,y)=\left\{\begin{array}{cc}ll{t}_s\hfill & r<R_0\\ {\displaystyle \frac{\left[ln\left(\frac{r}{R_0}\right)+t_{s}ln\left(\frac{r}{R_{00}}\right)\right]}{ln\left(\frac{R_0}{R_{00}}\right)}}\hfill & R_0\le R_{00}\\ -1\hfill & r\ge R_{00}\end{array}\right.$, $r=\sqrt{x^2+y^2}$,

$a=b=1$, $q=1/m(1+{\delta }_{\mu }cos\left(4\theta \right))$, where $\theta =arctan({\varphi }_y/{\varphi }_x),$$f={\displaystyle \frac{1}{{\epsilon }^2}}(\varphi (1-\varphi )(\varphi -{\displaystyle \frac{1}{2}}+30\epsilon \alpha S{\displaystyle \frac{u}{1+0.25u}}\varphi (1-\varphi )),$$p={\varphi }^3(10-15\varphi +6{\varphi }^2)/S$, S=dimensionless supercooling of the melt, m=ratio of the capillary to the kinetic length and $\epsilon =$ratio of the average interface thickness parameter.

In the isotropic case, with values of physical parameters given by $t_s=0.01,$$R_0=0.1,$$R_{00}=2R_0,$$S=0.8,$ $m=0.1,$$\alpha =70,$${\delta }_{\mu }=0,$ and $\epsilon =1/400,$ the evolution of the solution is’ given in Figures 1, $2\left(a\right)$(we took $h={10}^{-3}$, $\Delta t=2·{10}^{-5}$).

Figure 3. $\varphi ,u$ the isotropic case. (a):3D plot of $\varphi ,$ (b):trace at $x=0.1$ of $\varphi ,$ (c):3D plot of $u,$ (d): trace at $x=0.1$ of u (all at $t=2·{10}^{-3}$).

Figure 3. $\varphi ,u$ the isotropic case. (a):3D plot of $\varphi ,$ (b):trace at $x=0.1$ of $\varphi ,$ (c):3D plot of $u,$ (d): trace at $x=0.1$ of u (all at $t=2·{10}^{-3}$).

(a)	(b)
(c)	(d)

Figure 4. Succesive contours of the interface of every 100 time steps for the isotropic case.

Figure 4. Succesive contours of the interface of every 100 time steps for the isotropic case.