# Optimal Rate of Convergence for a Nonstandard Finite Difference Galerkin Method Applied to Wave Equation Problems

Authors: Pius W. M. Chin
Journal: Journal of Applied Mathematics
 (2013)

---

## Abstract

The optimal rate of convergence of the wave equation in both the energy and theL2-norms using continuous Galerkin method is well known. We exploit this
technique and design a fully discrete scheme consisting of coupling the nonstandard
finite difference method in the time and the continuous Galerkin method in the space
variables. We show that, for sufficiently smooth solution, the maximal error in the
L2-norm possesses the optimal rate of convergence O(h2+(Δt)2) where h is the mesh size and Δt is the time step size. Furthermore, we show that this scheme replicates
the properties of the exact solution of the wave equation. Some numerical experiments
should be performed to support our theoretical analysis.


---

## Body

## 1. Introduction

Most physical phenomena such as the acoustics, electromagnetic, and elastic problems are modeled by the wave equation. The qualitative solution of the model in the spacetime domain is always a very delicate but a fundamental issue that needs careful study. Our point of departure of this paper is to consider the following model of the wave equation: findu(x,t) such that
(1)∂2u∂t2-Δu=fin(x,t)∈Ω×(0,T),(2)u(x,0)=u0(x),(3)∂u(x,0)∂t=u1(x)x∈Ω,
where Ω is a smooth bounded domain in ℝ2 with smooth boundary ∂Ω. Problem (1)–(3) consists of constant coefficients and f(x,t) the source term, u0(x) and u1(x) are prescribed as the initial data. Furthermore, (0,T) is taken to be a finite time interval and the boundary conditions satisfied by u(x,t) are given by
(4)u(x,t)=0,(x,t)∈∂Ω×(0,T).
The methods which have been heavily used for the study of the wave equation (1)–(4) in physical life are the continuous as well as the discontinuous Galerkin methods; see [1, 2] for more details. The reason for these methods may be due to the way they deal with heterogeneous media and arbitrary shaped geometric objects, represented by unstructured grids.The advantages of the continuous Galerkin method are enormous. Firstly, the convergence theory of this method is based on lower regularity (differentiability) requirements than the finite difference and the spectral methods. Secondly, the method retains the important energy conservation properties provided by the discrete version of the initial/boundary valued problem such as the one under consideration. Thirdly, the computation and the analysis from the Galerkin method could be extended in the approximation of the nonlinear wave equation. Furthermore, the method could be applicable to problems of any desired order of accuracy. For more details on these advantages, see in [3–5].The apriori error estimate for continuous Galerkin approximation of the wave equation (1)–(4) was first derived by Dupont [6] and later improved by Baker [7], both for continuous and discrete time schemes. Gekeler [8] analyzed general multistep methods for the time discretization of second-order hyperbolic equation, when a Galerkin procedure is used in space. The nonclassical finite element treatment of the wave equation can be seen in Johnson [9] and Richter [10].In this paper, we exploit and present a reliable technique consisting of coupling the nonstandard finite difference (NSFD) method in time and the continuous Galerkin (CG) method in the space variables. A similar approach was done for the first time using parabolic problems more specifically the diffusion equations in the nonsmooth domain as seen in [11]. The NSFD method was initiated by Mickens in [12] and major contributions to the foundation of the NSFD method could be seen in [13, 14]. Since its initiation, the NSFD method has been extensively applied to many concrete problems in engineering and science; see [12, 15, 16] for an overview. This paper compliments the technique used in [11]. The technique is geared toward obtaining a sufficiently smooth solution, the maximal error in the L2-norm, and to show that the error across the entire interval convergences optimally as O(h2+Δt2) where h is the mesh size and Δt is the time step size. The reliability of this technique comes from the fact that the NSFD-CG method preserves both the energy features and the hyperbolicity of the exact solution of the wave equation (1)–(4).The organization of the paper is as follows. Under Section2, we review some of the useful spaces and their notations needed in the paper. In Section 3, we gather essential tools necessary to prove the main result of the paper. We present, in Section 4, a reliable scheme NSFD-CG and show that the numerical solution obtained from this scheme attains the optimal convergence rate in the energy as well as in the L2-norms. Furthermore, we show that the scheme under consideration replicates the properties of the exact solution. Section 5 is devoted to some numerical experiments using a numerical example which confirms the optimal rate of convergence of the solution proved analytically in Section 4. The concluding remarks are given in Section 6 and these underline how the work fits in the existing literature and also how it can be extended for further work.


## 2. Notations

In this section, we will review some of the spaces which we will be using in the paper together with their notations and possibly properties. Fors≥0, Hs(Ω) will denote the Sobolev space of real-valued functions on Ω, and the norm on Hs(Ω) will be denoted by ∥·∥s. See [17] for the definitions and the relevant properties of these spaces. In a particular case, where s=0 the space H0(Ω)=L2(Ω) and its inner product together with the norm will be denoted, respectively, by
(5)(u,v)=∫Ω‍uvdx,u,v∈L2(Ω),∥u∥L2(Ω)={(u,u)}1/2,u∈L2(Ω).
In addition, C0∞(Ω) will denote the space of infinitely differentiable functions with support compactly contained in Ω. The space H01(Ω) will denote the subspace of H1(Ω) obtained by completing C0∞(Ω) with respect to the norm ∥·∥1. Following [17], for X a Hilbert space, we will more generally use the Sobolev space Hs[(0,T);X], where s≥0 and in the case where s=0 we will have H0[(0,T);X]≡L2[(0,T);X] with norm
(6)∥v∥L2[(0,T);X]=(∫0T‍∥v(·,t)∥Xdt)1/2.
In practice, X will be the Sobolev space Hm(Ω) or H0m(Ω). Associated with (1) is the bilinear form
(7)a(u,v)=∫Ω‍∇u∇vdx,u,v∈H1(Ω).a(·,·) will be symmetric and positive definite; that is,
(8)a(u,v)=a(v,u),a(u,u)≥0.


## 3. The Continuous Galerkin Method

Having the previously mentioned notations in place, we proceed under this section to gather essential tools necessary to prove the main result of our paper. We begin by stating the following weak problem of (1)–(4).Findu∈L2[(0,T);H01(Ω)] given such that
(9)(∂2u(·,t)∂t2,v)-a(u(·,t),v)=(f(·,t),v)∀v∈H01(Ω),t≥0,(u(·,0),v)=(u0,v),(∂u(·,0)∂t,v)=(u1,v).
See [8] for the existence and the uniqueness of a solution u of (9). Hence forth, in appropriate places to follow, additional conditions on the regularity of u which guarantee the convergence results will be imposed. We continue next by providing the framework for stating the discrete version of (9). To this end, we let 𝒯h be a regular family of triangulations of Ω- consisting of compatible triangles T of diameter hT≤h; see [17] for more. For each mesh size 𝒯h, we associate the finite element space Vh of continuous piecewise linear function that are zero on the boundary
(10)Vh:={vh∈C0(Ω-);vh|∂Ω=0,vh|T∈P1,∀T∈𝒯h},
where P1 is the space of polynomials of degree less than or equal to 1 and Vh is a finite dimensional subspace of the Sobolev space H01(Ω). It is well known that Vh parametrized by h∈(0,1) possesses the following approximation properties: there exists a constant C such that if v∈H01(Ω)∩H2(Ω) we have
(11)infχ∈Vh{∥v-χ∥+h∥v-χ∥1}≤Ch2∥v∥2.
By the use of the energy method and Gronwall’s Lemma, there exists a discrete Galerkin solution uh∈Vh such that
(12)(∂2uh∂t2,vh)-a(uh,vh)=(f,vh),∀vh∈Vh,t∈[0,T],(uh,vh)=(Phu0,vh),(∂uh∂t,vh)=(Phu1,vh),
where Ph denote the L2-projection onto Vh. Furthermore, we let u∈H2(Ω) and we define the Galerkin projection πhu of u in Vh by requiring that
(13)a(πhu,v)=a(u,v),∀v∈Vh,(14)∂k(πhu)∂tk=πh(∂ku∂tk),k=0,1,2,t∈[0,T].
The previous essential tools lead to the immediate consequence of the approximation properties (11).Lemma 1.
Letu be the solution of (9). Then, there exists a unique mapping πhu∈L2[(0,T);Vh] which satisfies (13). Furthermore, if for some integer k≥0, ∂ku/∂tk∈Lp[(0,T);H2(Ω)], then
(15)∂k(πhu)∂tk∈Lp[(0,T);Vh],∥(∂∂t)k[u-πhu]∥Lp[(0,T);L2(Ω)]≤Ch2∥(∂∂t)ku∥Lp[(0,T);H2(Ω)]
for some constant C independent of the mesh size.


## 4. Coupled Nonstandard Finite Difference and Continuous Galerkin Methods

Instead of the continuous Galerkin method summarized previously, we present in this section a reliable scheme NSFD-CG consisting of the nonstandard finite difference method in time and the continuous Galerkin method in the space variable. We show that the numerical solution obtained from the scheme NSFD-CG attained the optimal convergence in both the energy andL2-norms. We proceed, in this regard, by letting the step size tn=nΔt for n=0,1,2,…,N. For a sufficiently smooth function v(x,t), we set
(16)(∂∂t)kvn=(∂∂t)kv(·,tn),vn=v(·,tn),k≥0.
With this, we proceed to find the NSFD-CG approximation {Uhn} such that Uhn≈uhn at discrete time tn. That is, find a sequence {Uhn}n=0N in Vh such that
(17)(δ2Uhn,vh)+a(Uhn,vh)=(fn,vh),∀vh∈Vh,n=1,2,…,N-1,
where
(18)δ2Uhn=Uhn+1-2Uhn+Uhn-1(ψ(Δt))2,n=1,2,…,N-1
and ψ(Δt)=2sin(Δt/2) restricted between 0<ψ(Δt)<1. The initial conditions Uh0∈Vh and Uh1∈Vh are given by
(19)(Uh0,vh)=(Phu0,vh),(20)(Uh1,vh)=(Uh0+ψ(Δt)Phu1+(ψ(Δt))22U~h0,vh),∀vh∈Vh,
where U~h0∈Vh is defined by
(21)(U~h0,vh)=(f0,vh)-a(u0,vh),∀vh∈Vh.
If f=0 in (1), we will have in view of (17) an exact scheme
(22)(Uhn+1-2Uhn+Uhn-14sin2(Δt/2),vh)+(∇Uhn,∇vh)=0,
which according to Mickens [12] replicates both the energy preserving features and the properties of the exact solution (1)–(4). In order to state and prove the main result, we need a framework on which this result is based. To this end, we proceed by decomposing the error denoted by en=un-Uhn into the following error equation:
(23)en=un-whn+whn-Uhn=ηn-ϕn,
where whn=πhun∈Vh is the Galerkin projector of un. Due to the regularity assumption mentioned earlier, the exact solution u of (1)–(4) satisfies
(24)(∂2un∂t2,vh)+a(un,vh)=(fn,vh),∀vh∈Vh,n=1,2,…,N.
Subtracting (17) from (24) and using some properties of the Galerkin projection in the space we have
(25)(δ2whn-δ2Uhn,vh)+a(whn-un,vh)=(δ2whn-∂2un∂t2,vh),
from where we obtain in view of (23)
(26)(δ2ϕn,vh)+(ϕn,vh)=(rn,vh),∀vh∈Vh,n=1,2,…,N-1,
where rn=δ2whn-∂2un/∂t2. In view of the necessity for error bound, we set
(27)rn={δ2whn-∂2un∂t2,forn=1,2,…,N-1ϕ1-ϕ0(ψ(Δt))2forn=0,
from where we define
(28)Rn=ψ(Δt)∑m=0n‍rm.
The previously mentioned framework leads to the following main result.Theorem 2.
Letu be the solution of (9) and {Uhn}n=0N in Vh a sequence defined by (17)–(21). Suppose that u∈L2[(0,T);H2(Ω)], ∂u/∂t∈L2[(0,T);H2(Ω)], and ∂ku/∂tk∈L2[(0,T);L2(Ω)] for k=3,4. Then, there exists a constant C>0 independent of Δt and the mesh refinement size h:
(29)max0≤n≤N∥u(·,Δt)-Uhn∥0≤C[h2+(Δt)2].
Furthermore, the discrete solution replicates all the properties of solution of the hyperbolic equation in the limiting case of the space independent equation.We prove the previous Theorem2 thanks to the following series of results.Proposition 3.
The following result holds for a constantC>0 that is independent of Δt and the mesh refinement size h:
(30)maxn=0,…,N∥en∥0≤C(∥e0∥0+maxn=0,…,N∥ηn∥0+ψ(Δt)maxn=0,…,N-1∥Rn∥0).Proof.
In view of (23), we have the following relation using triangular inequality:
(31)maxn=0,…,N∥en∥0≤maxn=0,…,N∥ϕn∥0+maxn=0,…,N∥ηn∥0.
We bound ϕn in (31) by first taking partial sums of the first term of (26), seconded by adding the remaining terms from n=1 to m for 1≤m≤N-1, and multiplying both sides by ψ(Δt) yields
(32)(ϕm+1-ϕmψ(Δt),vh)-(ϕ1-ϕ0ψ(Δt),vh)+ψ(Δt)∑n=1m‍a(ϕm,vh)=ψ(Δt)∑n=1m‍(rn,vh).
In view of (28), we can define
(33)Φm=ψ(Δt)∑n=1m‍ϕm,Φ0=0
and using this in (32) we have
(34)(ϕm+1-ϕmψ(Δt),vh)+a(Φm,vh)=(Rm,vh),∀vh∈Vh.
If we choose vh=ϕm+1+ϕm∈Vh in (34) and multiply the result by ψ(Δt), this yields
(35)∥ϕm+1∥02-∥ϕm∥02+ψ(Δt)a(Φm,ϕm+1+ϕm)=ψ(Δt)∑m=1n-1‍(Rm,ϕm+1+ϕm),for0≤m≤N-1.
Summing this from m=0 to m=n-1 for 1≤n≤N-1 gives
(36)∥ϕn∥02-∥ϕ0∥02+ψ(Δt)∑m=0n-1‍a(Φm,ϕm+1+ϕm)=ψ(Δt)∑m=0n-1‍(Rm,ϕm+1+ϕm).
Since a is symmetric, Φ0=0, and
(37)Φm+1-Φm-1=ψ(Δt)(ϕm+ϕm+1),m=1,2,…,N-1,
then, we deduce in view of the third term of (36)
(38)ψ(Δt)∑m=0n-1‍a(Φm,ϕm+1+ϕm)=∑m=0n-1‍a(Φm,Φm+1+Φm-1)=∑m=0n-1‍a(Φm,Φm+1)-∑m=0n-1‍a(Φm+1,Φm-1)=a(Φn-1,Φn).
By symmetry and coercivity properties of a and the fact that
(39)Φn-Φn-1=ψ(Δt)ϕnforn=1,2,…,N,
we have in view of (38)
(40)(ψ(Δt))22a(ϕn,ϕn)=12a(Φn-Φn-i,Φn-Φn-1)=12a(Φn,Φn)-a(Φn,Φn-1)+12a(Φn-1,Φn-1)
and so
(41)(ψ(Δt))22a(ϕn,ϕn)≥-a(Φn,Φn-1)
and this together with (36) yields
(42)∥ϕn∥02-(ψ(Δt))22a(ϕn,ϕn)≤∥ϕ0∥02+ψ(Δt)∑m=0n-1‍(Rn,ϕm+ϕm+1),for1≤n≤N.
By Poincare inequality together with the continuity and coercivity of a, we have the following inequality:
(43)∥ϕn∥02-(ψ(Δt))22∥ϕn∥02≤∥ϕ0∥02+ψ(Δt)∑m=0n-1‍(Rn,ϕm+ϕm+1).
We suppose at this stage that h and Δt satisfy the CFL condition such that ψ(Δt)/h<1. With this condition, the previous inequality becomes
(44)CΔt∥ϕn∥02≤∥ϕ0∥02+ψ(Δt)∑m=0n-1‍(Rn,ϕm+ϕm+1),1≤n≤N
from where we have 1/2<CΔt=1-(ψ(Δt))2/2h2<1 and 1/2<CΔt<1. By the use of Cauchy-Schwarz inequality on (44) and some algebraic manipulations, we have
(45)CΔt∥ϕn∥02≤∥ϕ0∥02+ψ(Δt)∑m=0n-1‍(Rm,ϕm+ϕm+1)≤∥ϕ0∥02+ψ(Δt)∑m=0n-1‍∥Rm∥0(∥ϕm∥0+∥ϕm+1∥0)≤∥ϕ0∥02+2ψ(Δt)∑m=0n-1‍∥Rm∥0maxm=0,…,n∥ϕm∥0≤∥ϕ0∥02+CΔt2maxm=1,…,n∥ϕm∥02+2CΔt(ψ(Δt)∑m=0n-1‍∥Rn∥0)2for0≤m≤n
and since the right-hand side is independent of n, then
(46)(CΔt-CΔt2)∥ϕm∥02≤∥ϕ0∥02+2CΔt(ψ(Δt)∑n=0N-1‍∥Rn∥0)2
which then implies that
(47)CΔt2∥ϕm∥02≤∥ϕ0∥02+2CΔt(ψ(Δt)∑n=0N-1‍∥Rn∥0)2.
Furthermore, we have
(48)∥ϕn∥02≤2CΔt∥ϕn∥02+4(ψ(Δt))2CΔt(∑n=0N-1‍∥Rn∥0)2
and hence the following result:
(49)maxn=0,…,N∥ϕn∥0≤2CΔt∥ϕ0∥0+2ψ(Δt)CΔt∑n=1N-1‍∥Rn∥0.
With (49) and the fact that
(50)∥ϕ0∥0≤∥e0∥0+∥η0∥0
we have in view of (31) the desired estimate of the proposition.We now proceed to bound the term containingRn on the right-hand side of Proposition 3. We achieve this by estimating in the L2-norm the function rn for the cases n=0 and follow by the case when n≥1.Lemma 4.
There holds(51)∥r0∥0≤C[(ψ(Δt))-1h2∥∂u∂t∥L2[(0,T);H2(Ω)]+ψ(Δt)∥∂2u∂t2∥L2[(0,T);L2(Ω)]],
with a constant C>0 that is independent of h and the mesh size.Proof.
In view of (27), we have r0=(ψ(Δt))-2(ϕ1-ϕ0). We estimate ∥ϕ1-ϕ0∥0 by taking vh∈Vh arbitrary as follows:
(52)(ϕ1-ϕ0,vh)=(wh1-Uh1,vh)-(wh0-Uh0,vh)=((πh-I)(u1-u0),vh)+(u1-U1,vh),
where we have used (u0-Uh0,vh)=(u0-Phu0,v)=0 in view of (19). It now follows from (52) that
(53)|((πh-I)(u1-u0),vh)|≤∫0t1‍|(∂∂t(πh-I)u(·,s),vh)|dsByLemma1≤∫0t1‍|((πh-I)∂u(·,s)∂t,vh)|dsBy(18)≤Cψ(Δt)h2∥∂u∂t∥L2[(0,T);H2(Ω)]∥vh∥0.
We now proceed to estimate the term (u1-Uh1,v) in (52) as follows: by Taylor’s formula and the fact u(·,0)=u0 and ∂u(·,0)/∂t=u1 in (1)–(3), we have
(54)u1=u0+ψ(Δt)u1+(ψ(Δt))22∂2u(·,0)∂t2+12∫0t1‍(ψ(Δt)-s)2∂3u(·,s)∂t3ds.
In view of the definition of Uh1 in (20) and the fact that
(55)(u0-Phu0,v)=0,(u1-Phu1,v)=0
we have
(56)Uh1=Uh0+ψ(Δt)Phu1+(ψ(Δ))22U~h0.
We then deduce from (54) and (56) that
(57)(u1-Uh1,vh)=(ψ(Δt))22(∂2u(·,0)∂t2-U~h0,vh)+12∫0t1‍(ψ(Δt)-s)2(∂3u(·,s)∂t3,vh)ds.
But by the definition of U~h0 in (20) we have in view of (21) that
(58)(∂2u(·,s)∂t2-U~h0,vh)=(f0,vh)-a(u0,vh)-(f0,vh)+a(u0,vh)=0
implying that
(59)|(u1-Uh1,v)|≤12∫0t1‍(ψ(Δt)-s)2|(∂3u(·,s)∂t3,v)|ds≤C(ψ(Δt))3∥∂3u∂t3∥L2[(0,T);L2(Ω)]∥v∥0.
Since ϕ1-ϕ0∈Vh, then (52) together with (53) and (59) yields
(60)∥ϕ1-ϕ0∥0≤C(ψ(Δt)h2∥∂u∂t∥L2[(0,T);H2(Ω)]+(ψ(Δt))3∥∂3u∂t3∥L2[(0,T);L2(Ω)])
and dividing throughout by (ψ(Δt))2 yields
(61)∥r0∥0≤C((ψ(Δt))-1h2∥∂u∂t∥L2[(0,T);H2(Ω)]+ψ(Δt)∥∂3u∂t3∥L2[(0,T);L2(Ω)])
as required.Lemma 5.
For1≤n≤N-1, there holds
(62)∥rn∥0≤C(h2ψ(Δt)∫tn-1tn+1‍∥∂2u(·,s)∂t2∥2ds+ψ(Δt)∫tn-1tn+1‍∥∂4u(·,s)∂t4∥0ds)
with a constant C>0 that is independent of h and mesh size.Proof.
In view of (27), we have
(63)rn=δ2whn-∂2un∂t2,forn=1,2,…,N-1.
By the triangular inequality, we have
(64)∥rn∥0=∥δ2whn-∂2un∂t2∥0≤∥δ2(πh-I)un∥0+∥δ2whn-∂2un∂t2∥0
and using the following identity
(65)v(·,tn+1)-2v(·,tn)+v(·,tn-1)=Δt∫tn-1tn+1‍(1-|s-tn|Δt)∂2v(·,s)∂t2ds
on the first term of the right-hand side of (64) we have
(66)∥δ2(πh-I)un∥0≤1ψ(Δt)∫tn-1tn+1‍(1-|s-tn|Δt)∥∂2(πh-I)u∂t2∥0ds≤Ch2ψ(Δt)∫tn-1tn+1‍∥∂2u(·,s)∂t2∥2ds
after the use of Lemma 1 and (14).
The second term of (64) can be estimated by the use of the identity
(67)δ2whn-∂2un∂t2=16(Δt)2∫tn-1tn+1‍(Δt-|s-tn|)3∂4u(·,s)∂t4ds
which is obtained from Taylor’s formulae with integral remainder. This is deduced as follows:
(68)∥δ2whn-∂2un∂t2∥0≤16(Δt)2∫tn-1tn+1‍(Δt-|s-tn|)3∥∂4u(·,s)∂t4∥0ds≤16(Δt)2∫tn-1tn+1‍(ψ(Δt))3∥∂4u(·,s)∂t4∥0ds≤ψ(Δt)6∫tn-1tn+1‍∥∂4u(·,s)∂t4∥0ds
using the relation that ψ(Δt)-|s-tn|≤ψ(Δt) in (68).In view of (64) using (66) and (68), we have the desired result for the bound of ∥rn∥0. We now assemble Lemmas 4 and 5 in the next proposition to obtain the bound Rn as follows.Proposition 6.
For0≤n≤N-1, there holds
(69)∥Rn∥0≤C(ψ(Δt))2(∥∂3u∂t3∥L2[(0,T);L2(Ω)]+∥∂4u∂t4∥L2[(0,T);L2(Ω)])+Ch2(∥∂u∂t∥L2[(0,T);H2(Ω)]+∥∂2u∂t2∥L2[(0,T);H2(Ω)]).Proof.
Using the bounds on∥rn∥0 in Lemmas 4 and 5, we have
(70)∥Rn∥0≤ψ(Δt)∥r0∥0+ψ(Δt)∑m=1N-1‍rm≤C(ψ(Δt))2(∥∂3u∂t3∥L2[(0,T);L2(Ω)]+∥∂4u∂t4∥L2[(0,T);L2(Ω)])+Ch2(∥∂u∂t∥L2[(0,T);H2(Ω)]+∥∂2u∂t2∥L2[(0,T);H2(Ω)])
and the proof is completed.With all these results, we are now in the position to prove the main result as follows.Proof of Theorem2.
Using Proposition3 and the fact that
(71)ψ(Δt)∑n=0N-1‍∥Rn∥0≤Tmaxn=0,…,N∥Rn∥0
we have
(72)maxn=0,…,N∥en∥≤C(∥e0∥0+maxn=0,…,N∥ηn∥0+Tmaxn=0,…,N-1∥Rn∥0).
By the use of Lemma 1, we can bound the second term on the right-hand side as follows.
(73)maxn=0,…,N∥ηn∥0≤Ch2∥u∥L2[(0,T);H2(Ω)].
From the approximation property of the L2-projection, we have
(74)∥e0∥0≤Ch2∥u0∥2≤Ch2∥u∥L2[(0,T);H2(Ω)].
Finally, by the bound of maxn=0,…,N-1∥Rn∥0 obtained via Proposition 6, we have the result
(75)maxn=0,…,N∥en∥0≤Ch2∥u∥L2[(0,T);H2(Ω)]+C(ψ(Δt))2(∥∂3u∂t3∥L2[(0,T);L2(Ω)]+∥∂4u∂t4∥L2[(0,T);L2(Ω)]).
In view of the relationship
(76)ψ(Δt)=2sin(Δt2)≈Δt+0[(Δt)2]
as proposed for such schemes in Mickens [12], we have the required estimate as Δt→0. Furthermore, using the fact that its uniform convergence results imply pointwise convergence for x∈Ω′ completes the proof.


## 5. Numerical Experiments

In this section, we present the numerical experiments on problem (1) using both the standard finite difference (SFD) and NSFD-CG methods. These experiments are performed in Ω=(0,1)2×(0,T) where Ω was discretized using regular meshes of sizes h=1/M in the space and Δt=T/N in the time. The M and N in such a discretization denote the number of nodes and time respectively. The initial data was considered to be u0(x)=x1x2(1-x1)(1-x2) and u1(x)=0 where these data were deduced from the exact solution
(77)u(x,t)=x1x2(1-x1)(1-x2)cos(πt).
Using the previous exact solution we obtained the right-hand side f of (1). In the computation, (1) together with (17)–(20) led to a system of equations
(78)AX=b.
In solving for X in the above system, we took the following values of M=10, 15, 20, 25, 50, and 100. For a fix M=20, N=10, and T=1, we had Figures 1–3 illustrating various solutions corresponding to their respective schemes.Figure 1
The Exact solution.Figure 2
Approximate solution for SFD-CG scheme.Figure 3
Approximate solution for NSFD-CG scheme.Figure1 shows the exact solution, Figure 2 the solution from the SFD-CG, and Figure 3 the solution from the NSFD-CG schemes, followed by Tables 1 and 2 which demonstrate various optimal rates of convergence in both H1 and L2-norms of these schemes.Table 1
NSFD method error in bothL2 and H1-norms.
M
L
2-error
H
1-error
RateL2
RateH1
10
3.800E-3
1.67E-2
15
1.700E-3
1.110E-2
1.983
1.0073
20
1.000E-3
8.200E-3
1.8445
0.9688
25
6.591E-4
6.700E-3
1.8678
1.0134
50
2.012E-4
3.400E-3
1.7119
0.9786
100
5.098E-5
1.800E-3
1.9806
0.9175Table 2
SFD method error in bothL2 and H1-norms.
M
L
2-error
H
1-error
RateL2
RateH1
10
3.400E-3
1.700E-2
15
1.700E-3
1.130E-2
1.709
1.007
20
1.050E-3
8.500E-3
1.674
0.989
25
1.600E-3
6.800E-3
1.887
1.000
50
5.170E-4
3.500E-3
1.629
0.958
100
1.799E-4
1.900E-3
1.522
0.881These optimal rates of convergence were calculated by using the formula(79)Rate=ln(e2/e1)ln(h2/h1),
where h1 and h2 together with e1 and e2 are successive triangle diameters and errors, respectively. These results are self-explanatory and we could conclude that the results as shown by these experiments exhibit the desired results as expected from our theoretical analysis.


## 6. Conclusion

We presented a reliable scheme of the wave equation consisting of the nonstandard finite difference method in time and the continuous Galerkin method in the space variable (NSFD-CG). We proved theoretically that the numerical solution obtained from this scheme attains the optimal rate of convergence in both the energy and theL2-norms. Furthermore, we showed that the scheme under investigation replicates the properties of the exact solution of the wave equation. We proceeded by the help of a numerical example and showed that the optimal rate of convergence as proved theoretically is guaranteed in both the energy and the L2-norms. This convergence results hold for any fully discrete NSFD-CG method where the scheme under consideration has a bilinear form which is symmetric, continuous, and coercive.The method presented in this paper can be extended to the nonlinear hyperbolic or parabolic problems with either smooth or nonsmooth domain if at all these cases followed the procedure as proposed by Mickens [12]. We will also exploit another form of nonstandard finite differential method as proposed in [18, 19].



---

*Source: 520219-2013-12-30.xml*