# On the Finite Volume Element Method for Self-Adjoint Parabolic Integrodifferential Equations

Authors: Mohamed Bahaj; Anas Rachid
Journal: Journal of Mathematics
 (2013)

---

## Abstract

Finite volume element schemes for non-self-adjoint parabolic integrodifferential equations are derived and stated. For the spatially discrete scheme, optimal-order error estimates inL2, H1, and Lp, W1,p norms for 2≤p<∞, are obtained. In this paper, we also study the lumped mass modification. Based on the Crank-Nicolson method, a time discretization scheme is discussed and related error estimates are derived.


---

## Body

## 1. Introduction

The main purpose of this paper is to study semidiscrete and full discrete finite volume element method (FVE) for parabolic integrodifferential equation of the form(1)ut-∇·(A(x,t)∇u)-∫0t‍∇·(B(x,t,s)∇u(s))ds=f(x,t),inΩ×(0,T],u=0,on∂Ω×(0,T],u(·,0)=u0,inΩ,
where Ω is a bounded domain in ℝd,  d=2,3, with smooth boundary ∂Ω, and T<∞. Here A(t), a non-self-adjoint second-order strongly elliptic, and B(t,s), an arbitrary second-order linear partial differential operator, both with coefficients depending smoothly on x and t,  f=f(x,t) and u0(x) are known functions, which are assumed to be smooth and satisfy certain compatibility conditions for x∈Ω and t=0, so that (1) has a unique solution in certain Sobolev space. Problem (1) occurs in nonlocal reactive flows in porous media, viscoelasticity, and heat conduction through materials with memory.Finite volume method is an important numerical tool for solving partial differential equations. It has been widely used in several engineering fields, such as fluid mechanics, heat and mass transfer, and petroleum engineering. The method can be formulated in the finite difference framework or in the Petrov-Galerkin framework. Usually, the former one is called finite volume method [1], marker and cell (MAC) method [2], or cell-centered method [3], and the latter one is called finite volume element method (FVE) [4–9], covolume method [10], or vertex-centered method [11, 12]. We refer to the monographs [13, 14] for general presentation of these methods. The most important property of FVE method is that it can preserve the conservation laws (mass, momentum, and heat flux) on each control volume. This important property, combined with adequate accuracy and ease of implementation, has attracted more people to do research in this field.Recently, the authors in [8, 15] studied FVE method for general self-adjoint elliptic problems. The authors in [16] presented and analyzed the semidiscrete and full discrete symmetric finite volume schemes for a class of parabolic problems. In [6, 7] the authors have studied FVE for one- and two-dimensional parabolic integrodifferential equations and have obtained an optimal-order estimate in the L2-norm. The regularity required on the exact solution u is W3,p for p>1 which is higher when compared to that for finite element methods.The aim of this paper is to study the convergence of FVE discretization for a nonself-adjoint parabolic integrodifferential problem (1). Both spatially discrete scheme and discrete-in-time scheme are analyzed, and optimal error estimates in L2 and H1 norms are proved using only energy method. We also explore and generalize that idea to develop the lumped mass modification and Lp estimates, 2≤p<∞. Our analysis avoids the use of semigroup theory, and the regularity requirement on the solution is the same of that of finite element method. Furthermore, based on the Crank-Nicolson method the fully discrete scheme is analyzed and the related optimal error estimates are established.This paper is organized as follows. In Section2, we introduce some notations and present some preliminary materials to be used later. The Ritz-Volterra projection to finite volume element spaces is introduced and related estimates are carried out in Section 3. In Section 4, we estimate the error of the finite volume element approximations derived in the previous section. In Section 5, the lumped mass is presented and optimal estimates in L2 and H1 norms are obtained Finally, the Crank-Nicolson scheme is studied in Section 6.


## 2. Finite Volume Element Scheme

In this section, we introduce some material which will be used repeatedly hareafter. Throughout this paper,C (with or without index) denotes a generic positive constant which does not depend on the spatial and time discretization parameters h and k,  respectively.


### 2.1. Notations

We will use∥·∥m and |·|m (resp., ∥·∥m,p and |·|m,p) to denote the norm and seminorm of the Sobolev space Hm(Ω) (resp., Wm,p(Ω)). The scalar product and norm in L2(Ω) are denoted by (·,·) and ∥·∥, respectively. Let H10(Ω) be the standard Sobolev subspace of H1(Ω) of functions vanishing on ∂Ω.The weak form of (1) is used to find u(·,t):[0,T]→H01(Ω),  such that
(2)(ut,v)+A(t;u,v)+∫0t‍B(t,s;u(s),v)ds=(f,v),∀v∈H01(Ω),u(0)=u0,
where
(3)A(t;u,v)=∫Ω‍A(x,t)∇u·∇v,B(t,s;u(s),v)=∫Ω‍B(x,t,s)∇u(s)·∇v.Let𝒯h be a decomposition of Ω into triangles (for the 2D case) or tetrahedral (for the 3D case) with h=maxhK, where hK is the diameter of the element K∈𝒯h.In order to describe the FVE method for solving problem (1), we will introduce a dual partition 𝒯h* based upon the original partition 𝒯h whose elements are called control volumes. We construct the control volumes in the same way as in [7, 17]. Let zK be a point of K∈𝒯h. In the 2D case, on each edge e of K a point qe is selected; then we connect zK with line segments to qe; thus, partitioning K into three quadrilaterals Kz,  z∈Zh(K), where Zh(K) are the vertices of K. Then with each vertex z∈Zh=∪K∈𝒯hZh(K) we associate a control volume Vz, which consists of the union of the subregions Kz, sharing the vertex  z  (see Figure 1).(a) A sample region with blue lines indicating the corresponding control volumeVz. (b) A triangle K partitioned into three subregions Kz.
(a)
(b)Similarly, in the 3D case, on each of the four facesSi,i=1,…,4, a point qSi,i=1,…,4, is selected, and on each of the six edges e a point qe is selected. On each of the two faces S1 and S2 of K sharing an edge e, we connect qSi,i=1,2,  with qe and with zK by line segments, thus, partitioning K into twelve tetrahedron Kz,  z∈Zh(K) (see Figure 2). Then for z∈Zh the control volume Vz consists of the union of the subregions Kz sharing the vertex z. Thus, we finally obtain a group of control volumes covering the domain Ω, which is called the dual partition 𝒯h* of the triangulation 𝒯h. We denote by Zh0 the set of interior vertices and Nh=#Zh0. For a vertex zi∈Zh0,  let Π(i) be the index set of those vertices that along with zi are in some element of Th (Figure 2).Figure 2
A tetrahedronK partitioned into twelve subregions Kz.There are various ways to introduce a regular dual partition𝒯h*. In this paper, we will also use the construction of the control volumes in which we let zK be the barycenter of K∈𝒯h. In the 2D case, we choose qe to be the midpoint of the edge e (see Figure 3).z
K is the barycenter of K, and qe is to be the midpoint of the edge e.
(a)
(b)In the 3D case, we chooseqe to be the midpoint of the edge e and qSi to be the barycenter of the face Si  (Figure 4).Figure 4
q
e is the midpoint of the edge e, and  qSi is the barycenter of the face Si.We call the partition𝒯h* regular or quasiuniform, if there exists a positive C>0 such that
(4)C-1h2≤meas(Vz)≤Ch2,∀Vz∈𝒯h*.If the finite element triangulation𝒯h is quasiuniform, that is, there exists a positive C>0 such that
(5)C-1h2≤meas(K)≤Ch2,∀K∈𝒯h,
then the dual partition𝒯h*is also quasiuniform.Based on the triangulationTh, let Sh be the standard conforming finite element space of piecewise linear functions, defined on the triangulation Th as follows:
(6)Sh={v∈𝒞(Ω):v|Kislinear∀K∈Th,andv|Γ=0}.LetIh:𝒞(Ω)→Sh be the standard interpolation operators, such that
(7)Ihu=∑z∈Zh0‍vz(t)φz(x),∀v∈Sh,
where {φz}z∈Zh0 are the standard basis functions of Sh and vz(t)=v(t;z).


### 2.2. Construction of the FVE Scheme

We formulate the FVE method for the problem (1) as follows: Given a z∈Zh0, integrating (1)1 over the associated control volume Vz and applying Green’s formula, we obtain an integral conservation as follows form:
(8)∫Vz‍ut-∫∂Vz‍A(x,t)∇u·nds-∫∂Vz‍B(x,t,s)∇u·nds=∫Vz‍f(x,t),
wherendenotes the unit outer normal vector to ∂Vz.LetIh*:𝒞(Ω)→Sh* be the transfer operator defined by
(9)Ih*v=∑z∈Zh0‍v(z)χz,∀v∈Sh,
where
(10)Sh*={v∈L2(Ω):vi|Vzisconstant,∀z∈Zh0},
and χz is the characteristic function of the control volume Vz.Now fort>0 and for an arbitrary Ih*v, we multiply (8) by v(z), and sum over all z∈Zh0. Then the semidiscrete FVE approximation uh of (1) is a solution to the following problem: find uh(t)∈Sh for t>0 such that
(11)(uht,vh)+A(t;uh,vh)+∫0t‍B(t,s;uh(s),vh)ds=(f,vh)vh∈Sh*,uh(0)=u0h∈Sh.Here the bilinear formsA(t;u,v) and B(t,s;u,v) are defined by
(12)A(t;u,v)={-∑z∈Zh0vi∫∂Vz‍A(x,t)∇u·nds,(u,v)∈((H01∩H2)∪Sh)×Sh*,∫Ω‍A(x)∇u·∇vdx,(u,v)∈H01×H01,B(t,s;u,v)={-∑z∈Zh0vz∫∂Vz‍B(x,t,s)∇u·nds,(u,v)∈((H01∩H2)∪Sh)×Sh*,∫Ω‍B(x,t,s)∇u·∇vdx,(u,v)∈H01×H01.Let(13)uh=∑‍j=1Nhαz(t)φz(x),α(t)=(α1(t),α2(t),…,αNh(t))T.
Then, we can rewrite scheme (11)1 as systems of ordinary differential equations as follows:
(14)Mhα′(t)+Ah(t)α(t)+∫0t‍Bh(t,s)α(s)ds=Fh(t).
Here Fh(t)=(f1(t),f2(t),…,fNh(t))T,  the mass matrix Mh={Mhij}={(φi,χj)} is tridiagonal, and both Ah(t)={A(t;φi,χj)} and Bh(t,s)={B(t,s;φi,χj)} are positive definites.In order to describe features of the bilinear forms defined in (11), we introduce some discrete norms on Sh in the same way as in [7]:
(15)∥vh∥0,h2=(vh,vh)0,h=(Ih*vh,Ih*vh),|vh|1,h2=∑xi∈Zh0∑xj∈Π(i)meas(Vi)(vi-vjdij)2,∥vh∥1,h2=∥vh∥0,h2+|vh|1,h2,|||vh|||2=(vh,Ih*vh),
where dij=d(xi,xj),  the distance between xi and xj. Obviously, these norms are well defined for vh∈Sh* as well and ∥vh∥0,h=|||vh|||.Hereafter we state the equivalence of discrete norms∥·∥0,h and ∥·∥1,h with usual norms ∥·∥ and ∥·∥1 on Sh, respectively.Lemma 1 (see [7]).
There exist two positive constantsC0 and C1 such that for all vh∈Sh, we have
(16)C0∥vh∥0,h≤∥vh∥≤C1∥vh∥0,h,∀vh∈Sh,C0|||vh|||≤∥vh∥≤C1|||vh|||,∀vh∈Sh,C0∥vh∥1,h≤∥vh∥1≤C1∥vh∥1,h,∀vh∈Sh.Next we recall some properties of the bilinear forms (see [7, 18])Lemma 2 (see [7]).
There exist two positive constantsC and C0 such that for all uh,vh∈Sh,  we have
(17)A(uh,Ih*vh)≤C∥uh∥1∥vh∥1,∀uh,vh∈Sh,A(vh,Ih*vh)≥C0∥vh∥12,∀vh∈Sh.The following lemmas are proved in [3, 7], which give the key feature of the bilinear forms in the FVE method.Lemma 3 (see [3]).
Assume thatφ∈W01,p. Then, one has
(18)A(t;φ,vh)-A(t;φ,Ih*vh)=∑‍K∈τh∫∂K‍(A(t)∇φ·n)(vh-Ih*vh)ds-∑‍K∈τh∫K‍(∇·A(t)∇φ)(vh-Ih*vh)ds,∀vh∈Sh.
The aforementioned identity holds true when A(·,·) is replaced by B(t,s;·,·).Lemma 4 (see [3]).
Assume thatφ∈Sh. Then, one has
(19)A(t;φ,χ)-A(t;φ,Ih*χ)≤Ch|φ|1,p|χ|1,q.
Furthermore, for φ∈W01,p∩W2,p,  we have
(20)A(t;φ,χ)-A(t;φ,Ih*χ)≤Ch∥φ∥2,p∥χ∥1,q.


## 2.1. Notations

We will use∥·∥m and |·|m (resp., ∥·∥m,p and |·|m,p) to denote the norm and seminorm of the Sobolev space Hm(Ω) (resp., Wm,p(Ω)). The scalar product and norm in L2(Ω) are denoted by (·,·) and ∥·∥, respectively. Let H10(Ω) be the standard Sobolev subspace of H1(Ω) of functions vanishing on ∂Ω.The weak form of (1) is used to find u(·,t):[0,T]→H01(Ω),  such that
(2)(ut,v)+A(t;u,v)+∫0t‍B(t,s;u(s),v)ds=(f,v),∀v∈H01(Ω),u(0)=u0,
where
(3)A(t;u,v)=∫Ω‍A(x,t)∇u·∇v,B(t,s;u(s),v)=∫Ω‍B(x,t,s)∇u(s)·∇v.Let𝒯h be a decomposition of Ω into triangles (for the 2D case) or tetrahedral (for the 3D case) with h=maxhK, where hK is the diameter of the element K∈𝒯h.In order to describe the FVE method for solving problem (1), we will introduce a dual partition 𝒯h* based upon the original partition 𝒯h whose elements are called control volumes. We construct the control volumes in the same way as in [7, 17]. Let zK be a point of K∈𝒯h. In the 2D case, on each edge e of K a point qe is selected; then we connect zK with line segments to qe; thus, partitioning K into three quadrilaterals Kz,  z∈Zh(K), where Zh(K) are the vertices of K. Then with each vertex z∈Zh=∪K∈𝒯hZh(K) we associate a control volume Vz, which consists of the union of the subregions Kz, sharing the vertex  z  (see Figure 1).(a) A sample region with blue lines indicating the corresponding control volumeVz. (b) A triangle K partitioned into three subregions Kz.
(a)
(b)Similarly, in the 3D case, on each of the four facesSi,i=1,…,4, a point qSi,i=1,…,4, is selected, and on each of the six edges e a point qe is selected. On each of the two faces S1 and S2 of K sharing an edge e, we connect qSi,i=1,2,  with qe and with zK by line segments, thus, partitioning K into twelve tetrahedron Kz,  z∈Zh(K) (see Figure 2). Then for z∈Zh the control volume Vz consists of the union of the subregions Kz sharing the vertex z. Thus, we finally obtain a group of control volumes covering the domain Ω, which is called the dual partition 𝒯h* of the triangulation 𝒯h. We denote by Zh0 the set of interior vertices and Nh=#Zh0. For a vertex zi∈Zh0,  let Π(i) be the index set of those vertices that along with zi are in some element of Th (Figure 2).Figure 2
A tetrahedronK partitioned into twelve subregions Kz.There are various ways to introduce a regular dual partition𝒯h*. In this paper, we will also use the construction of the control volumes in which we let zK be the barycenter of K∈𝒯h. In the 2D case, we choose qe to be the midpoint of the edge e (see Figure 3).z
K is the barycenter of K, and qe is to be the midpoint of the edge e.
(a)
(b)In the 3D case, we chooseqe to be the midpoint of the edge e and qSi to be the barycenter of the face Si  (Figure 4).Figure 4
q
e is the midpoint of the edge e, and  qSi is the barycenter of the face Si.We call the partition𝒯h* regular or quasiuniform, if there exists a positive C>0 such that
(4)C-1h2≤meas(Vz)≤Ch2,∀Vz∈𝒯h*.If the finite element triangulation𝒯h is quasiuniform, that is, there exists a positive C>0 such that
(5)C-1h2≤meas(K)≤Ch2,∀K∈𝒯h,
then the dual partition𝒯h*is also quasiuniform.Based on the triangulationTh, let Sh be the standard conforming finite element space of piecewise linear functions, defined on the triangulation Th as follows:
(6)Sh={v∈𝒞(Ω):v|Kislinear∀K∈Th,andv|Γ=0}.LetIh:𝒞(Ω)→Sh be the standard interpolation operators, such that
(7)Ihu=∑z∈Zh0‍vz(t)φz(x),∀v∈Sh,
where {φz}z∈Zh0 are the standard basis functions of Sh and vz(t)=v(t;z).


## 2.2. Construction of the FVE Scheme

We formulate the FVE method for the problem (1) as follows: Given a z∈Zh0, integrating (1)1 over the associated control volume Vz and applying Green’s formula, we obtain an integral conservation as follows form:
(8)∫Vz‍ut-∫∂Vz‍A(x,t)∇u·nds-∫∂Vz‍B(x,t,s)∇u·nds=∫Vz‍f(x,t),
wherendenotes the unit outer normal vector to ∂Vz.LetIh*:𝒞(Ω)→Sh* be the transfer operator defined by
(9)Ih*v=∑z∈Zh0‍v(z)χz,∀v∈Sh,
where
(10)Sh*={v∈L2(Ω):vi|Vzisconstant,∀z∈Zh0},
and χz is the characteristic function of the control volume Vz.Now fort>0 and for an arbitrary Ih*v, we multiply (8) by v(z), and sum over all z∈Zh0. Then the semidiscrete FVE approximation uh of (1) is a solution to the following problem: find uh(t)∈Sh for t>0 such that
(11)(uht,vh)+A(t;uh,vh)+∫0t‍B(t,s;uh(s),vh)ds=(f,vh)vh∈Sh*,uh(0)=u0h∈Sh.Here the bilinear formsA(t;u,v) and B(t,s;u,v) are defined by
(12)A(t;u,v)={-∑z∈Zh0vi∫∂Vz‍A(x,t)∇u·nds,(u,v)∈((H01∩H2)∪Sh)×Sh*,∫Ω‍A(x)∇u·∇vdx,(u,v)∈H01×H01,B(t,s;u,v)={-∑z∈Zh0vz∫∂Vz‍B(x,t,s)∇u·nds,(u,v)∈((H01∩H2)∪Sh)×Sh*,∫Ω‍B(x,t,s)∇u·∇vdx,(u,v)∈H01×H01.Let(13)uh=∑‍j=1Nhαz(t)φz(x),α(t)=(α1(t),α2(t),…,αNh(t))T.
Then, we can rewrite scheme (11)1 as systems of ordinary differential equations as follows:
(14)Mhα′(t)+Ah(t)α(t)+∫0t‍Bh(t,s)α(s)ds=Fh(t).
Here Fh(t)=(f1(t),f2(t),…,fNh(t))T,  the mass matrix Mh={Mhij}={(φi,χj)} is tridiagonal, and both Ah(t)={A(t;φi,χj)} and Bh(t,s)={B(t,s;φi,χj)} are positive definites.In order to describe features of the bilinear forms defined in (11), we introduce some discrete norms on Sh in the same way as in [7]:
(15)∥vh∥0,h2=(vh,vh)0,h=(Ih*vh,Ih*vh),|vh|1,h2=∑xi∈Zh0∑xj∈Π(i)meas(Vi)(vi-vjdij)2,∥vh∥1,h2=∥vh∥0,h2+|vh|1,h2,|||vh|||2=(vh,Ih*vh),
where dij=d(xi,xj),  the distance between xi and xj. Obviously, these norms are well defined for vh∈Sh* as well and ∥vh∥0,h=|||vh|||.Hereafter we state the equivalence of discrete norms∥·∥0,h and ∥·∥1,h with usual norms ∥·∥ and ∥·∥1 on Sh, respectively.Lemma 1 (see [7]).
There exist two positive constantsC0 and C1 such that for all vh∈Sh, we have
(16)C0∥vh∥0,h≤∥vh∥≤C1∥vh∥0,h,∀vh∈Sh,C0|||vh|||≤∥vh∥≤C1|||vh|||,∀vh∈Sh,C0∥vh∥1,h≤∥vh∥1≤C1∥vh∥1,h,∀vh∈Sh.Next we recall some properties of the bilinear forms (see [7, 18])Lemma 2 (see [7]).
There exist two positive constantsC and C0 such that for all uh,vh∈Sh,  we have
(17)A(uh,Ih*vh)≤C∥uh∥1∥vh∥1,∀uh,vh∈Sh,A(vh,Ih*vh)≥C0∥vh∥12,∀vh∈Sh.The following lemmas are proved in [3, 7], which give the key feature of the bilinear forms in the FVE method.Lemma 3 (see [3]).
Assume thatφ∈W01,p. Then, one has
(18)A(t;φ,vh)-A(t;φ,Ih*vh)=∑‍K∈τh∫∂K‍(A(t)∇φ·n)(vh-Ih*vh)ds-∑‍K∈τh∫K‍(∇·A(t)∇φ)(vh-Ih*vh)ds,∀vh∈Sh.
The aforementioned identity holds true when A(·,·) is replaced by B(t,s;·,·).Lemma 4 (see [3]).
Assume thatφ∈Sh. Then, one has
(19)A(t;φ,χ)-A(t;φ,Ih*χ)≤Ch|φ|1,p|χ|1,q.
Furthermore, for φ∈W01,p∩W2,p,  we have
(20)A(t;φ,χ)-A(t;φ,Ih*χ)≤Ch∥φ∥2,p∥χ∥1,q.


## 3. Ritz-Volterra Projection and Related Estimates

Following [7, 19, 20], we define the Ritz-Volterra projection Vh(t):H01→Sh as follows:
(21)A(t;u-Vhu,Ih*vh)+∫0t‍B(t,s;u(s)-Vhu(s),Ih*vh)ds=0,t>0,∀vh∈Sh.ThisVh(t) is an elliptic projection with memory of u into Sh*. It is easy to see that (21) is actually a system of integral equations of Volterra type. In fact, if Vh(t)u=∑j=1Nhαj(t)φj(x), then (21) can be rewritten as
(22)Ah(t)α(t)+∫0t‍Bh(t,s)α(s)ds=Fh(t),
where Ah(t),  Bh(t,s) are matrices and α(t),  Fh(t) are vectors, defined via
(23)α(t)=(α1(t),α2(t),…,αNh(t))T,Fhk(t)=A(t;u,χk)+∫0t‍B(t,s;u(s),χk)ds,k=1,2,…,Nh,Ah(t)=A(t,φk(x),χl),Bh(t,s)=B(t,s;φk(x),χl).From the positivity ofA (Lemma 2) and the linearity of (22), we see that the system (22) possesses a unique solution α(t). Consequently, Vh(t)u in (21) is well defined.Setρ=u-Vh(t)u. The following lemma was proved in [7], which shows the H1 error estimate for ρ and its temporal derivative.Lemma 5 (see [7]).
Assume thatDtnu∈L∞(H01∩H2) for all 0≤n≤k, for some integer k≥0.  Then, for T>0 fixed there is a constant C=C(T;k)>0, independent of h and u, such that for all 0≤n≤k and 0<t<T,
(24)∥ρ(t)∥1≤Ch(∥u∥2+∫0t‍∥u∥2ds),∥Dtnρ(t)∥1≤Ch(∑‍i=0n∥Dtiu∥2+∫0t‍∥u∥2ds).Now we establishL2 error estimate for ρ and its temporal derivative which improves Theorem  2.2 in [7]. This estimate is optimal with respect to the order.Lemma 6.
Assume that, for some integerk≥0,  Dtnu∈L∞(H01∩H2)  for all 0≤n≤k. Then, for T>0 fixed there is a constant C=C(T;k)>0, independent of h and u, such that for all 0≤n≤k and 0<t<T,
(25)∥ρ(t)∥≤Ch2(∥u∥2+∫0t‍∥u∥2ds),∥Dtnρ(t)∥≤Ch2(∑‍i=0n∥Dtiu∥2+∫0t‍∥u∥2ds).Proof.
The proof will proceed by duality argument. Letψ∈H2(Ω)∩H01(Ω) be the solution of
(26)A*(t)ψ=ρ,inΩψ=0,in∂Ω.
The solution ψ∈H2(Ω)∩H01(Ω) satisfies the following regularity estimate:
(27)∥ψ∥2≤C∥ρ∥.
Multiplying this equation by ρ and then taking L2 innerproduct over Ω, we obtain the following:
(28)∥ρ∥2=A(t;ρ,ψ)=A(t;ρ,ψ-Rhψ)+A(t;ρ,Rhψ-Ih*(Rhψ))-∫0t‍B(t,s;ρ(s),Ih*Rhψ-Rhψ)ds-∫0t‍B(t,s;ρ(s),Rhψ-ψ)ds-∫0t‍B(t,s;ρ(s),ψ)ds=I1+I2+I3+I4+I5.
We have
(29)|I1|+|I4|≤Ch2(∥u∥2+∫0t‍∥u∥2ds)∥ψ∥2.
Applying Lemma 4, we obtain
(30)|I2|+|I3|≤Ch2(∥u∥2+∫0t‍∥u∥2ds)∥ψ∥2.
Finally, we have
(31)|I5|≤∫0t‍(ρ(s),B*(t,s)ψ)ds≤C(∫0t‍∥ρ∥ds)∥ψ∥2;
then, we have
(32)∥ρ∥≤Ch2(∥u∥2+∫0t‍∥u∥2ds)+C(∫0t‍∥ρ∥ds).
Finally, an application of Gronwall’s lemma yields the first estimate.
The second inequality follows in a similar fashion.Lemma 7.
There exists a constantC independent of h such that
(33)∥ρ∥0,p+h∥ρ∥1,p≤Ch2(∥u∥2,p+∫0t‍∥u∥2,pds).Proof.
Letρx be an arbitrary component of ∇ρ, with p and q conjugate indices; we have ∥ρx∥p=sup{(ρx,φ);φ∈𝒞0∞(Ω),∥φ∥q=1}.
For any suchφ, let ψ be the solution of
(34)A*(t;ψ,v)=-(φx,v)∀v∈H01(Ω),ψ=0,on∂Ω.
It follows from the regularity theory for the elliptic problem that
(35)∥ψ∥1,q≤Cp∥φ∥q=Cp.
We then have by application of (21) that
(36)(ρx,φ)=A(t;ρ,ψ)=A(t;ρ,ψ-Rhψ)+A(t;ρ,Rhψ-Ih*(Rhψ))+∫0t‍B(t,s;ρ(s),Ih*(Rhψ))ds=I1+I2+I3,A(t;ρ,ψ-Rhψ)=A(t;Rhu-u,ψ)=-((Rhu-u)x,φ)≤Ch∥u∥2,p.
Applying Lemma 4, we have
(37)I2=A(t;u,Rhψ-Ih*(Rhψ))-A(t;Vhu,Rhψ-Ih*(Rhψ))≤Ch∥u∥2,p.
Finally, I3 is estimated as follows:
(38)I3=∫0t‍B(t,s;ρ(s),Ih*(Rhψ))ds≤Cp∫0t‍∥ρ∥1,pds.
Combining these estimates, we get
(39)∥ρ∥1,p≤Ch∥u∥2,p+Cp∫0t‍∥ρ∥1,pds,
hence by Gronwall’s lemma
(40)∥ρ∥1,p≤Ch(∥u∥2,p+∫0t‍∥u∥2,pds).
The derivation of the error estimate in Lp is similar to the case when p=2.


## 4. Error Estimates for Semidiscrete Approximations

We split the errore(t)=u(t)-uh(t) as follows:
(41)e(t)=(u(t)-Vhu(t))+(Vhu(t)-uh(t))=ρ+θ.It is easy to see thatθ=Vhu(t)-uh(t)∈Sh satisfies an error equation of the form
(42)(θt,Ih*vh)+A(t;θ,Ih*vh)+∫0t‍B(t,s;θ(s),Ih*vh)ds=-(ρt,Ih*vh),vh∈Sh.Since the estimates ofρ are already known, it is enough to have estimates for θ.We will prove a sequence of lemmas which lead to the following result.Lemma 8.
There is a positive constantC independent of h such that
(43)|||θ(t)|||≤C(|||θ(0)|||2+∫0t‍∥ρt∥ds).Proof.
Sinceθ∈Sh we may take vh=θ in (42) to obtain
(44)12ddt|||θ(t)|||2+c∥θ∥12≤∥ρt∥∥θ∥+C∫0t‍∥θ∥1ds∥θ∥1≤∥ρt∥∥θ∥+12c∥θ∥12+C∫0t‍∥θ∥12ds,
and hence by integration and Lemma 1, we have
(45)||θ(t)||2+∫0t‍∥θ∥12ds≤C(|||θ(0)|||2+∫0t‍∥ρt∥∥θ∥ds+∫0t‍∫0s‍∥θ(τ)∥12dτds).
Gronwall’s lemma now implies the following:
(46)|||θ(t)|||2+∫0t‍∥θ∥12ds≤C(|||θ(0)|||2+∫0t‍∥ρt∥∥θ∥ds)≤C|||θ(0)|||2+12sups≤t∥θ(s)∥2+(∫0t‍∥ρt∥ds)2.
Since this holds for all ∈J, we may conclude that
(47)||θ(t)||≤C(|||θ(0)|||+∫0t‍∥ρt∥ds).Remark 9.
If the initial value was chosen so that∥u0h-u0∥≤Ch2∥u0∥2,  then ∥θ(0)∥≤∥u0h-u0∥+∥Vhu0-u0∥≤Ch2∥u0∥2. One can derive
(48)|||θ(t)|||≤Ch2(∥u0∥2+∫0t‍∥ut∥2ds).Lemma 10.
There is a positive constantC independent of h such that
(49)∫0t‍∥θt∥2ds+∥θ∥12≤C(∥θ(0)∥12+∫0t‍∥ρt∥2ds).Proof.
Setvh=θt in (42) to get
(50)∥θt∥2+12ddtA(t;θ,Ih*θ)=-(ρt,Ih*θt)-∫0t‍B(t,s;θ(s),Ih*θt(t))ds+12At(t;θ,Ih*θ)+12[A(t;θt,Ih*θ)-A(t;θ,Ih*θt)]≤12∥ρt∥2+12|||θt|||2+At(t;θ,Ih*θ)+12[A(t;θt,Ih*θ)-A(t;θ,Ih*θt)]-ddt∫0t‍B(t,s;θ(s),Ih*θ(t))ds+B(t,t;θ(t),Ih*θ(t))+∫0t‍Bt(t,s;θ(s),Ih*θ(t))ds.
Then(51)∥θt∥2+ddtA(t;θ,Ih*θ)≤∥ρt∥2-ddt∫0t‍B(t,s;θ(s),Ih*θ)ds+C(∥θ∥12+∫0t‍∥θ(s)∥12ds)+12[A(t;θt,Ih*θ)-A(t;θ,Ih*θt)].
In addition, recall that(52)A(t;uh,Ih*vh)-A(t;vh,Ih*uh)≤Ch∥uh∥1∥vh∥1,∀uh,vh∈Sh;
then applying an inverse inequality and using kickback argument, we obtain
(53)[A(t;θt,Ih*θ)-A(t;θ,Ih*θt)]≤Ch∥θt∥1∥θ∥1≤C∥θt∥∥θ∥1≤ε∥θt∥2+C∥θ∥12.
Combining these estimates, we derive
(54)∥θt∥2+ddtA(t;θ,Ih*θ)≤∥ρt∥2-ddt∫0t‍B(t,s;θ(s),Ih*θ)ds+C(∥θ∥12+∫0t‍∥θ(s)∥12ds).
So after integration in time and using the weak coercivity of A(t;θ,Ih*θ), we get
(55)∫0t‍∥θt∥2ds+c0∥θ∥12≤c0∥θ(0)∥12+∫0t‍∥ρt∥2ds+∫0t‍B(t,s;θ(s),Ih*θ)ds+C∫0t‍∥θ(s)∥12ds≤c0∥θ(0)∥12+c2∥θ∥12+C(∫0t‍∥ρt∥2+∥θ(s)∥12ds),
and by Gronwall’s lemma,
(56)∫0t‍∥θt∥2ds+c∥θ∥12≤C(∥θ(0)∥12+∫0t‍∥ρt∥2ds).Remark 11.
Ifθ(0)=0, then
(57)∫0t‍∥θt∥2ds+c∥θ∥12≤Ch2(∫0t‍∥ut∥22ds).Theorem 12 (error estimates inL2 and H1-norms).
Letu,  uh be the solutions of (2) and (11), respectively. Assume that u,ut∈L∞(H01∩H2).(a)
Letu0h be chosen so that ∥u0h-u0∥≤Ch2∥u0∥2. Then for T>0 fixed there is a constant C=C(T) independent of h, such that for all 0<t<T,
(58)∥uh(t)-u(t)∥≤Ch2(∥u0∥2+∫0t‍∥ut∥2ds).(b)
Letu0h be chosen so that ∥u0h-u0∥1≤Ch∥u0∥2. Then for T>0 fixed there is a constant C=C(T) independent of h, such that for all 0<t<T,
(59)∥uh(t)-u(t))∥1≤Ch(∥u0∥2+∫0t‍∥ut∥2ds).We now prove error estimates for FVE approximations inLp and W1,p-norms.Theorem 13 (error estimates inLpandW1,p-norms).
Letu,  uh be the solutions of (2) and (11), respectively and u0h=Vhu0. Assume that u,ut∈L∞(H01∩W2,p). For h sufficiently small, we have
(60)∥u-uh∥0,p≤Ch2(∥u0∥2+∫0t‍∥ut∥2ds),∥u-uh∥1,p≤Ch(∥u0∥2+∥u∥2,p+∫0t‍∥ut∥2ds).Proof.
If2≤p<∞, by the following Sobolev embedding inequality
(61)∥θ∥0,p≤C∥θ∥1,
then the first desired estimate follows from Lemmas 7 and 10.
Givenφ∈𝒞0∞(Ω), find ψ∈H01(Ω) such that
(62)A(t)*ψ=-φx,inΩ,ψ=0,on∂Ω,∥ψ∥1,q≤∥φ∥0,q.
We have
(63)((u-uh)x,φ)=A(t;u-uh,ψ)=A(t;u-uh,ψ-Rhψ)+A(t;u-uh,Rhψ-Ih*Rhψ)-∫0t‍B(t,s;(u-uh)(s),Ih*Rhψ)ds-((u-uh)t,Ih*Rhψ)=I1+I2+I3+I4,|I1|≤|A(t;u-Rhu,ψ)|≤C∥u-Rhu∥1,p∥ψ∥1,q≤Ch∥u∥2,p∥ψ∥1,q.
By Lemma 4(64)|I2|≤A(t;u-uh,Rhψ-Ih*Rhψ)≤Ch(|u-uh|1,p+|u|2,p)∥ψ∥1,q,|I3|≤∫0t‍∥u-uh∥1,pds∥ψ∥1,q,|I4|≤(∥u-uh∥)∥ψ∥≤Ch2(∥u0∥2+∫0t‍∥ut∥2ds)∥ψ∥1,q,
where we have used the fact ∥ψ∥≤∥ψ∥1,r,  r>1. Combining these estimates, we get
(65)|((u-uh)x,φ)|≤Ch(∥u0∥2+∥u∥2,p+∫0t‍∥ut∥2ds)∥ψ∥1,q,∥(u-uh)x∥0,p=sup((u-uh)x,φ)∥φ∥0,q≤Ch|u-uh|1,p+Ch(∥u0∥2+∥u∥2,p+∫0t‍∥ut∥2ds).
Hence using the Poincaré inequality, we have forhsufficiently small
(66)∥u-uh∥1,p≤Ch(∥u0∥2+∥u∥2,p+∫0t‍∥ut∥2ds).We compare the relationship between covolume solution and the Galerkin finite element solution.Corollary 14.
Letu~hbe the finite element solution to (2); that is,
(67)(u~ht,vh)+A(t;u~h,vh)+∫0t‍B(t,s;u~h(s),vh)ds=(f,vh),vh∈Sh,u~h(0)=Rhu0.
For h sufficiently small, we have
(68)∥(u~h-uh)∥1,p≤C(h∥u-uh∥1,p+∥(u-uh)t∥+∥(u~h-u)t∥+∫0t‍(∥(u-uh)(s)∥1,p+∥(u-u~h)(s)∥1,p)ds)≤C(u)h.Proof.
By (2) and (67),
(69)((u~h-u)t,vh)+A(t;u~h-u,vh)+∫0t‍B(t,s;(u~h-u)(s),vh)ds=0,vh∈Sh.
Consider the following auxiliary problem. For any such φ, let ψ be the solution of the following:
(70)A(t)*ψ=-φx,inΩ,ψ=0,on∂Ω,
with
(71)∥ψ∥1,q≤∥φ∥0,q,((u~h-uh)x,φ)=A(t;u~h-uh,ψ)=A(t;u~h-uh,ψ-Rhψ)+A(t;u-uh,Rhψ)-A(t;u-uh,Ih*Rhψ)-((u-uh)t,Ih*Rhψ)-∫0t‍B(t,s;(u-uh)(s),Ih*Rhψ)ds+A(t;u~h-u,Rhψ)=[A(t;u-uh,Rhψ)-A(t;u-uh,Ih*Rhψ)]-((u-uh)t,Ih*Rhψ)-((u~h-u)t,Rhψ)-∫0t‍B(t,s;(u-uh)(s),Ih*Rhψ)ds-∫0t‍B(t,s;(u~h-u)(s),Rhψ)ds=I1+I2+I3.
On the other hand,
(72)|I1|≤Ch∥u-uh∥1,p∥ψ∥1,q,|I2|≤C(∥(u-uh)t∥+∥(u~h-u)t∥)∥ψ∥≤C(∥(u-uh)t∥+∥(u~h-u)t∥)∥ψ∥1,q,
where we have used the fact ∥ψ∥≤∥ψ∥1,r,  r>1(73)|I3|≤∫0t‍(∥(u-uh)(s)∥1,p+∥(u-u~h)(s)∥1,p)ds∥ψ∥1,q∥(u~h-uh)x∥0,p=supφ∈𝒞0∞((u~h-uh)x,φ)∥φ∥0,q≤C(h∥u-uh∥1,p+∥(u-uh)t∥+∥(u~h-u)t∥+∫0t‍(∥(u-uh)(s)∥1,p+∥(u-u~h)(s)∥1,p)ds).
We deduce the result from the known finite element estimates.Remark 15.
In order to estimate∥(u-uh)t∥, by differentiating (42) with respect to t, we obtain
(74)(θtt,Ih*vh)+A(t;θt,Ih*vh)+At(t;θt,Ih*vh)+B(t,t;θ,Ih*vh)+∫0t‍Bt(t,s;θ(s),Ih*vh)ds=-(ρtt,Ih*vh).
Setting vh=θt, we obtain
(75)12ddt|||θt|||2+c∥θt∥12≤∥ρtt∥∥θt∥+12c∥θt∥12+C∥θ∥12+∫0t‍∥θ∥12ds≤∥ρtt∥∥θt∥+12c∥θt∥12+C∫0t‍∥θt∥12ds.
Using kickback argument, integrating and applying Gronwall’s lemma, we deduce
(76)∥θt∥≤C(∥θt(0)∥+∫0t‍∥ρtt∥1ds)≤Ch2(∥u0∥2+∥u1∥2+∫0t‍∥utt∥2ds).


## 5. The Lumped Mass Finite Volume Element Method

In this section, we restrict our study to the 2D case. A simple way to define the lumped mass scheme [21] is to replace the mass matrix Mh in (14) by the diagonal matrix M¯h obtained by taking for its diagonal elements the numbers M¯hii=∑j=1NhMhij or by lumping all masses in one row into the diagonal entry. This makes the inversion of the matrix in front of α′(t) a triviality. We will therefore study the matrix problem
(77)M¯hα′(t)+Ah(t)α(t)+∫0t‍Bh(t,s)α(s)ds=Fh(t).We know that the lumped mass method defined by (77) above is equivalent to
(78)(Ih*uht,Ih*vh)+A(t;uh,Ih*vh)+∫0t‍B(t,s;uh(s),Ih*vh)ds=(f,Ih*vh),vh∈Sh.Our alternative interpretation of this procedure will be to think of (77) as being obtained by evaluating the first term in (78) by numerical quadrature. Let K be a triangle of the triangulation Th, let xj,j=1,2,3,  be its vertices, and consider the quadrature formula
(79)QK,h(f)=13areaK∑‍j=13f(xj)≃∫K‍fdx.We may then define the associated bilinear form inSh×Sh*, using the quadrature scheme, by the following:
(80)(vh,ηh)h=∑‍K∈ThQK,h(vhηh)=∑‍xi∈Nha`vh(xi)ηh(xi)|Vxi|,∀vh∈Sh,ηh∈Sh*.We note that∥vh∥h2=(vh,Ih*vh)h is a norm in Sh which is equivalent to the L2-norm uniformly in h; that is, there exist two positive constants C1 and C2 such that for all vh∈Sh, we have
(81)C0∥vh∥≤∥vh∥h≤C1∥vh∥,∀vh∈Sh.We note that the aforementioned definition(vh,ηh)h may be used also for ηh∈Sh  and that (vh,wh)h=(vh,Ih*wh)h for vh,wh∈Sh.The lumped mass method defined by (78) is equivalent to
(82)(uht,Ih*vh)h+A(t;uh,Ih*vh)+∫0t‍B(t,s;uh(s),Ih*vh)ds=(f,Ih*vh),vh∈Sh.We introduce the quadrature error(83)εh(vh,wh)=(vh,wh)h-(vh,wh),Lemma 16 (see [21]).
Letvh,wh∈Sh. Then
(84)|εh(vh,wh)|≤Ch2∥∇vh∥∥∇wh∥.Theorem 17.
Letuh and u be the solutions of (82) and (2), respectively, and assume uh(0)=Rhu0. Then we have for the error in the lumped mass semidiscrete method, for t≥0, the following:
(85)∥uh(t)-u(t)∥≤Ch2(∥u0∥2+∥u∥2+∫0t‍∥ut∥2ds).Proof.
In order to estimate∥θ∥, we write
(86)(θt,Ih*vh)h+A(t;θ,Ih*vh)+∫0t‍B(t,s;θ(s),Ih*vh)ds=(uht,Ih*vh)h+A(t;uh,Ih*vh)+∫0t‍B(t,s;uh(s),Ih*vh)ds-((Vhu)t,Ih*vh)h-A(t;Vhu,Ih*vh)-∫0t‍B(t,s;Vhu(s),Ih*vh)ds=(f,Ih*vh)-((Vhu)t,Ih*vh)h-A(t;u,Ih*vh)-∫0t‍B(t,s;u(s),Ih*vh)=(ut,Ih*vh)-((Vhu)t,Ih*vh)h=-(ρt,Ih*vh)-((Vhu)t,Ih*vh)h+((Vhu)t,Ih*vh).
We rewrite
(87)((Vhu)t,Ih*vh)h-((Vhu)t,Ih*vh)=((Vhu)t,Ih*vh)h-((Vhu)t,vh)+((Vhu)t,vh)-((Vhu)t,Ih*vh)=εh((Vhu)t,vh)+((Vhu)t,vh)-((Vhu)t,Ih*vh).(θt,Ih*vh)h+A(t;θ,Ih*vh)+∫0t‍B(t,s;θ(s),Ih*vh)ds=-(ρt,Ih*vh)+εh((Vhu)t,vh)+((Vhu)t,vh)-((Vhu)t,Ih*vh).
Setting vh=θ in (87), we obtain
(88)12ddt∥θ∥h2+c0∥θ∥12≤∥ρt∥∥θ∥+12c0∥θ∥12+C∫0t‍∥θ∥12ds+εh((Vhu)t,θ)+((Vhu)t,θ)-((Vhu)t,Ih*θ).
Using Lemma 16 and the inverse estimate, we get
(89)|εh(Vhut,θ)|≤Ch2∥∇(Vhu)t∥∥∇θ∥≤Ch2∥∇ut∥∥∇θ∥≤Ch∥∇ut∥∥θ∥,
we have
(90)|((Vhu)t,θ)-((Vhu)t,Ih*θ)|≤Ch∥∇ut∥∥θ∥.
Using Young’s inequality and Gronwall’s lemma to eliminate ∥θ∥1 on the right-hand side and using integration in t, we get the result. 							(91)12ddt∥θ∥h2+c0∥θ∥≤∥ρt∥∥θ∥+Ch∥∇ut∥∥θ∥.
Using Young’s inequality to eliminate ∥θ∥ on the right hand side it becomes.
Using integration int, we get the result.We will now show that theH1-norm error bound of theorem remains valid for the lumped mass method (82).Theorem 18.
Letuh and u be the solutions of (82) and (2), respectively, and assume
(92)uh(0)=Rhu0,∥u1h(0)-u1∥≤Ch2∥u1∥2.
Then, we have for the error in the lumped mass semidiscrete method, for t≥0, the following:
(93)∥uh(t)-u(t)∥1≤Ch2(∥u0∥2+∥u1∥2+∫0t‍∥utt∥2ds).Proof.
Settingvh=θt in (87), we obtain
(94)∥θt∥h2+12ddtA(t;θ,Ih*θ)=12At(t;θ,Ih*θ)+12[A(t;θt,Ih*θ)-A(t;θ,Ih*θt)]B(t,t;θ(t),Ih*θ(t))+∫0t‍Bt(t,s;θ(s),Ih*θ(t))ds-ddt∫0t‍B(t,s;θ(s),Ih*θ(t))ds-(ρt,Ih*θt)-εh((Vhu)t,θt)+((Vhu)t,θ)-((Vhu)t,Ih*θ).
It follows, thus, that using integration in t and Gronwall’s lemma, we have
(95)∫0t‍∥θt∥h2+∥θ∥12≤C∥∇θ(0)∥2+C∫0t‍∥ρt∥∥θt∥ds+Ch2∫0t‍∥ut∥12ds.


## 6. Full Discretization

Let∂¯Un=(Un-Un-1)/k be the backward difference quotient of Un; assume that Ah=PhA is a discrete analogue of A (similarly Bh=PhB), where Ph:L2(Ω)→Sh* the L2 projection is defined by
(96)(Phv,Ih*vh)=(v,Ih*vh),v∈L2(Ω),vh∈Sh.In order to define fully discrete approximation of (11), we discretize the time by taking tn=nk,  k>0,  n=1,2,…   and use the numerical quadrature
(97)∫0tn-1/2‍g(s)ds≈∑‍k=1nωn,kg(tk-1/2),tn-1/2=(n-12)k.Here{ωn,k} are the integration weights and we assume that the following error estimate is valid:
(98)qn(g)=∫0tn-1/2‍g(s)ds-∑‍k=1nωn,kg(tk-1/2)≤Ck2∫0tn‍(|g′|+|g′′|)ds.Now, define our complete discrete FVE approximation of (11) by the following: find Un∈Sh for n=1,2,…, such that for all vh∈Sh(99)(∂¯Un,Ih*vh)+A(tn-1/2;Un-1/2,Ih*vh)+∑k=1nωn,kB(tn-1/2,tk-1/2,Uk-1/2,Ih*vh)=(fn-1/2,Ih*vh),U0inSh,
where Un-1/2=(Un+Un-1)/2.Theorem 19.
Letu(t) and Un be the solutions of problem (2) and its complete discrete scheme (99), respectively. Then for any T>0 there exists a positive constant C=C(T)>0, independent of h, such that for 0<tn≤T(100)∥u(tn)-Un∥≤Ch2(∥u0∥2+∫0tn‍∥ut∥2ds)+Ck2(∫0tn‍(∥u∥2+∥ut∥2+∥utt∥2+∥uttt∥)ds).Proof.
Let us split the error into two parts:u(tn)-Un=ρn+θn,  where  ρn=u(tn)-Vhu(tn) and θn=Vhu(tn)-Un, and let W=Vhu(t)∈Sh be the Ritz-Volterra projection of u. Then from (2) and (99) we have for all vh∈Sh the following:
(101)(∂¯θn,Ih*vh)+A(tn-1/2;θn-1/2,Ih*vh)+∑‍k=1nωn,kB(tn-1/2,tk-1/2,θk-1/2,Ih*vh)=-(rn,Ih*vh),∀vh∈Sh,
where
(102)rn=rn1+rn2+rn3+rn4,rn1=∂¯ρn,rn2=∂¯u(tn)-ut(tn-1/2),rn3=A(tn-1/2;(u(tn)+u(tn-1))2-u(tn-1/2)),rn4=qn(BhW)=∑‍k=1nωn,kBh(tn-1/2,tk-1/2,Wk-1/2,Ih*vh)-∫0tn-1/2‍B(tn,s,W(s),Ih*vh)ds.
In fact, by Taylor expansion,(103)un+1=un+ku′(tn)+∫tntn+1‍u′′(s)(tn+1-s)ds=un+ku′(tn)+k22u′′(tn)+k36u(3)(tn)+16∫tntn+1‍u(4)(s)(tn+1-s)3ds,
we have
(104)∥rn1∥=∥∂¯ρn∥≤1k∫tn-1tn‍∥ρt∥ds≤Ch2k∫tn-1tn‍∥utt∥2ds,∥rn2∥=∥∂¯u(tn)-ut(tn-1/2)∥=1k∥∫tn-1tn‍(ut(s)-ut(tn-1/2))ds∥≤Ck∫tn-1tn‍∥u(3)(s)∥ds,∥rn3∥=∥A(tn-1/2;u(tn)+u(tn-1)2-u(tn-1/2),Ih*vh)∥≤Ck∫tn-1tn‍∥Autt(s)∥ds≤Ck∫tn-1tn‍∥utt∥2ds.
In addition, the quadrature error satisfies
(105)∥rn4∥=qn-1/2(BhW)=∑‍k=1nωn,kB(tn-1/2,tk-1/2,Wk-1/2,Ih*vh)-∫0tn-1/2‍B(tn,s,W(s),Ih*vh)ds≤Ck2∫0tn‍∥(BhW)ss∥ds≤Ck2∫0tn‍(∥u∥2+∥ut∥2+∥utt∥2)ds.k∑‍n=1N∥rn∥≤Ch2∫0tn‍∥utt∥2ds+Ck2∫0tn‍(∥u∥2+∥ut∥2+∥utt∥2+∥u(3)∥)ds.
Taking vh=θn-1/2 in (101) and noting that (∂¯θn,Ih*θn-1/2)=(1/2)∂¯|||θn|||2, there is
(106)|||θn|||2-|||θn-1|||2+2kc|||θn-1/2|||12≤Ck2∑‍k=1n∥θk-1/2∥1∥θn-1/2∥1+Ck∥rn∥∥θn-1/2∥≤kc|||θn-1/2|||12+Ck2∑‍k=1n∥θk-1/2∥12+Ck∥rn∥∥θn-1/2∥.
Summing from n=1 to N and then, after cancelling the common factor and using Gronwall’s lemma, we obtain
(107)|||θN|||2≤C|||θ0|||2+Ck∑‍k=1N∥rn∥(∥θk∥+∥θk-1/2∥),
and then
(108)|||θN|||≤C|||θ0|||+Ck∑‍n=1N∥rn∥,
the theorem follows from the estimates of ρn and rn.



---

*Source: 464893-2013-08-06.xml*