Expanded Mixed Finite Element Method for the Two-Dimensional Sobolev Equation

Abstract

Expanded mixed finite element method is introduced to approximate the two-dimensional Sobolev equation. This formulation expands the standard mixed formulation in the sense that three unknown variables are explicitly treated. Existence and uniqueness of the numerical solution are demonstrated. Optimal order error estimates for both the scalar and two vector functions are established.

1. Introduction

In this paper, we consider the following Sobolev equation: 𝑢_𝑡(𝐱,𝑡)−∇⋅𝑎(𝐱,𝑡)∇𝑢_𝑡+𝑏₁()((),()(],𝑢()()[],𝑢()𝐱,𝑡∇𝑢𝐱,𝑡)=𝑓𝐱,𝑡𝐱,𝑡∈Ω×0,𝑇𝐱,𝑡=0,𝐱,𝑡∈𝜕Ω×0,𝑇𝐱,0=𝑢₀(𝐱),𝐱∈Ω, (1) in a bounded domain Ω⊂𝑅² with Lipschitz continuous boundary 𝜕Ω. Equation (1) has a wide range of applications in many mathematical and physical problems [1, 2], for example, the percolation theory when the fluid flows through the cracks, the transfer problem of the moisture in the soil, and the heat conduction problem in different materials. So there exists great and actual significance to research Sobolev equation. Up till now, there are some different schemes studied to solve this kind of equation (see [1–4] for instance).

The mixed finite element method, which is a finite element method [5] with constrained conditions, plays an important role in the research of the numerical solution for partial differential equations. Its general theory was proposed by Babuska [6] and Brezzi [7]. Falk and Osborn [8] improved their theory and expanded the adaptability of the mixed finite element method. The mixed finite element method (see [9, 10] for instance) is wildly used for the modeling of fluid flow and transport, as it provides accurate and locally mass conservative velocities.

The main motivation of the expanded mixed finite element method [11–15] is to introduce three (or more) auxiliary variables for practical problems, while the traditional finite element method and mixed finite element method can only approximate one and two variables, respectively. The expanded mixed finite element method also has some other advantages except introducing more variables. It can treat individual boundary conditions. Also, this method is suitable for differential equation with small coefficient (close to zero) which does not need to be inverted. Consequently, this method works for the problems with small diffusion or low permeability terms in fluid problems. Using this method, we can get optimal order error estimates for certain nonlinear problems, while standard mixed formulation sometimes gives only suboptimal error estimates.

The object of this paper is to present the expanded mixed finite element method for the Sobolev equation. We conduct theoretical analysis to study the existence and uniqueness and obtain optimal order error estimates. The rest of this paper is organized as follows. In Section 2, the mixed weak formulation and its mixed element approximation are considered. In Section 3, we prove the existence and uniqueness of approximation form. In Section 4, some lemmas are given. In Section 5, optimal order semidiscrete error estimates are established.

Throughout the paper, we will use 𝐶, with or without subscript, to denote a generic positive constant which does not depend on the discretization parameter ℎ. Vectors will be expressed in boldface. At the same time, we show a useful 𝜀-Cauchy inequality 𝑎𝑏≤𝜀𝑎²+1𝑏4𝜀²,𝜀>0. (2)

2. Mixed Weak Form and Mixed Element Approximation

For 1≤𝑞≤+∞ and 𝑚 any nonnegative integer, let 𝑊^𝑚,𝑞(Ω)=𝑓∈𝐿^𝑞(Ω)∣𝐷^𝛼𝑓∈𝐿^𝑞(Ω),|𝛼|≤𝑚 (3) denote the Sobolev spaces [16] endowed with the norm 𝑓_𝑚,𝑞,Ω=_{|𝛼|≤𝑚}𝐷^𝛼𝑓^𝑞_𝐿𝑞(Ω)^1/𝑞 (4) (the subscript Ω will always be omitted).

Let 𝐻^𝑚(Ω)=𝑊^𝑚,2(Ω) with the norm ‖⋅‖_𝑚=‖⋅‖_𝑚,2. The notation ‖⋅‖ will mean ‖⋅‖_𝐿²_(Ω) or ‖⋅‖_𝐿2(Ω)². We denote by (,) the inner product in either 𝐿²(Ω) or 𝐿²(Ω)²; that is, 𝜉,𝜂=_Ω𝜉𝜂𝑑𝐱,(𝝃,𝜼)=_Ω𝝃⋅𝜼𝑑𝐱. (5) The notation ⟨,⟩ denotes the 𝐿²-inner product on the boundary of Ω𝜉,𝜂=_𝜕Ω𝜉𝜂𝑑𝑠. (6)

To formulate the weak form, let 𝑏=𝑏₁/𝑎; we introduce two vector variables 𝝀=−∇𝑢_𝑡−𝑏∇𝑢,𝝈=𝑎𝝀=−𝑎∇𝑢_𝑡−𝑏₁∇𝑢. (7) Letting 𝐜=−∇𝑏, we rewrite (1) as the following system: 𝑢_𝑡+∇⋅𝝈=𝑓,(𝐱,𝑡)∈Ω×(0,𝑇],𝝀+∇𝑢_𝑡()()(],()(],𝑢()()(],𝑢()−∇𝑏𝑢+𝐜𝑢=0,𝐱,𝑡∈Ω×0,𝑇𝝈−𝑎𝝀=0,𝐱,𝑡∈Ω×0,𝑇𝐱,𝑡=0,𝐱,𝑡∈𝜕Ω×0,𝑇𝐱,0=𝑢₀(𝐱),𝐱∈Ω. (8) Let 𝐧 be the unit exterior normal vector to the boundary of Ω. Then (8) is formulated in the following expanded mixed weak form: find (𝑢,𝝀,𝝈)∈𝑊×Λ×𝐕, such that 𝑢_𝑡+()=,()−𝑢,𝑤div𝝈,𝑤𝑓,𝑤∀𝑤∈𝑊,0<𝑡≤𝑇,𝝀,𝐯_𝑡+()()−()(𝑢(0))=𝑢+𝑏𝑢,div𝐯𝐜𝑢,𝐯=0,∀𝐯∈𝐕,0<𝑡≤𝑇,𝝈,𝝉𝑎𝝀,𝝉=0,∀𝝉∈Λ,0<𝑡≤𝑇,,𝑤₀,𝑤,∀𝑤∈𝑊, (9) where 𝑊=𝐿²(Ω),Λ=𝐿²(Ω)²,𝐕=𝐻(div,Ω)=𝐯∈𝐿²(Ω)²∣∇⋅𝐯∈𝐿²(Ω). (10)

Let 𝑇_ℎ be a quasiregular polygonalization of Ω (by triangles or rectangles), with ℎ being the maximum diameter of the elements of the polygonalization. Let 𝑊_ℎ×Λ_ℎ×𝐕_ℎ be a conforming mixed element space with index 𝑘 and discretization parameter ℎ. 𝑊_ℎ×Λ_ℎ×𝐕_ℎ is an approximation to 𝑊×Λ×𝐕. There are many conforming (or compatible) mixed element function spaces such as Raviart-Thomas elements [17], BDFM elements [18, 19]. Some RT type mixed elements are listed in Table 1. Here, 𝑃_𝑘(𝑇) is the polynomial up to 𝑘 order in two-dimensional domain used in triangle, while 𝑄_𝑘,𝑙(𝑇) is the polynomial up to 𝑘,𝑙 in each dimension used in rectangle.

Table 1 

Some RT type mixed elements.

Replacing the three variables by their approximation, we get the expanded mixed finite element approximation problem: find (𝑢_ℎ,𝝀_ℎ,𝝈_ℎ)∈𝑊_ℎ×Λ_ℎ×𝐕_ℎ, such that 𝑢_ℎ,𝑡,𝑤+div𝝈_ℎ,𝑤=𝑓,𝑤,∀𝑤∈𝑊_ℎ𝝀,0<𝑡≤𝑇,_ℎ,𝐯−𝑢_ℎ,𝑡+𝑏𝑢_ℎ,div𝐯+𝐜𝑢_ℎ,𝐯=0,∀𝐯∈𝐕_ℎ𝝈,0<𝑡≤𝑇,_ℎ,𝝉−𝑎𝝀_ℎ,𝝉=0,∀𝝉∈Λ_ℎ(𝑢(0))=𝑢,0<𝑡≤𝑇,,𝑤₀,𝑤,∀𝑤∈𝑊_ℎ, (11) where 𝑊_ℎ(𝐸)=𝜔∈𝑊∶𝜔|_𝐸∈𝑊_ℎ(𝐸),∀𝐸∈𝑇_ℎ,Λ_ℎ(𝐸)=𝝉∈Λ∶𝝉|_𝐸∈𝐕_ℎ(𝐸),∀𝐸∈𝑇_ℎ,𝐕_ℎ(𝐸)=𝐯∈𝐕∶𝐯|_𝐸∈𝐕_ℎ(𝐸),∀𝐸∈𝑇_ℎ. (12)

The error analysis next makes use of three projection operators. The first operator is the Raviart-Thomas projection (or Brezzi-Douglas-Marini projection) Π_ℎ∶𝐕→𝐕_ℎ; Π_ℎ satisfies ∇⋅𝐯−Π_ℎ𝐯,𝜔_ℎ=0,∀𝜔_ℎ∈𝑊_ℎ. (13) The following approximation properties are well known. 𝐯−Π_ℎ𝐯≤𝐶ℎ^𝑟‖𝐯‖_𝑟,12≤𝑟≤𝑘+1,∇⋅(𝐯−Π_ℎ𝐯)_𝑟≤𝐶ℎ^𝑟‖∇⋅𝐯‖_𝑟,0≤𝑟≤𝑘+1. (14)(15) The other two operators are the standard 𝐿²-projection [5] 𝑃_ℎ and 𝑅_ℎ onto 𝑊_ℎ and Λ_ℎ, respectively, 𝜔−𝑃_ℎ𝜔,∇⋅𝐯_ℎ=0,∀𝜔∈𝑊,𝐯_ℎ∈𝐕_ℎ,𝝉−𝑅_ℎ𝝉,𝝉_ℎ=0,∀𝝉∈Λ,𝝉_ℎ∈Λ_ℎ. (16) They have the approximation properties 𝜔−𝑃_ℎ𝜔≤𝐶ℎ^𝑟‖𝜔‖_𝑟,0≤𝑟≤𝑘+1,𝝉−𝑅_ℎ𝝉≤𝐶ℎ^𝑟‖𝝁‖_𝑟,0≤𝑟≤𝑘+1. (17) The two projections Π_ℎ and 𝑃_ℎ preserve the commuting property div∘Π_ℎ=𝑃_ℎ∘div∶𝐻¹(Ω)²⟶𝑊_ℎ. (18)

3. Existence and Uniqueness of Approximation Form

In this section, we consider the existence and uniqueness of the solution of (11).

Lemma 1. Equation (11) has a unique solution.

Proof. In fact, this equation is linear; it suffices to show that the associated homogeneous system 𝑢_ℎ,𝑡,𝑤+div𝝈_ℎ,𝑤=0,∀𝑤∈𝑊_ℎ𝝀,0<𝑡≤𝑇,_ℎ,𝐯−𝑢_ℎ,𝑡+𝑏𝑢_ℎ,div𝐯+𝐜𝑢_ℎ,𝐯=0,∀𝐯∈𝐕_ℎ𝝈,0<𝑡≤𝑇,_ℎ,𝝉−𝑎𝝀_ℎ,𝝉=0,∀𝝉∈Λ_ℎ(𝑢(0)),0<𝑡≤𝑇,,𝑤=0,∀𝑤∈𝑊_ℎ, (19) has only the trivial solution. In the first equation of (19), if we take 𝜔=𝑢_ℎ,𝑡 and 𝜔=div𝝈_ℎ, respectively, then we have that 𝑢_ℎ,𝑡≤div𝝈_ℎ,div𝝈_ℎ≤𝑢_ℎ,𝑡. (20) By (20), it is easy to see that div𝝈_ℎ=𝑢_ℎ,𝑡. (21) Choosing 𝜔=div𝝈_ℎ, 𝐯=𝝈_ℎ, and 𝝉=𝝀_ℎ in (19), we get 𝑎𝝀_ℎ,𝝀_ℎ+div𝝈_ℎ,div𝝈_ℎ=𝑏𝑢_ℎ,div𝝈_ℎ−𝐜𝑢_ℎ,𝝈_ℎ. (22) In the third equation of (19), letting 𝝉=𝝈_ℎ and 𝝉=𝝀_ℎ, respectively, then 𝝈_ℎ≤𝐶𝝀_ℎ,𝝀_ℎ≤𝐶𝝈_ℎ. (23)(24) Using the 𝜀-Cauchy inequality to (22), we have 𝑎₀𝝀_ℎ²+div𝝈_ℎ²≤𝐶𝑢_ℎ²+12div𝝈_ℎ²+𝑎₀2𝝀_ℎ². (25) Note that (21), we deserve 𝝀_ℎ≤𝐶𝑢_ℎ,𝑢_ℎ,𝑡≤𝐶𝑢_ℎ. (26)(27) From (27), we know 𝑢_ℎ=^𝑡₀𝑢_ℎ,𝑡𝑑𝑡≤𝐶^𝑡₀𝑢_ℎ,𝑡𝑑𝑡≤𝐶^𝑡₀𝑢_ℎ𝑑𝑡. (28) Using Gronwall’s inequality to (28), we can prove 𝑢_ℎ=0. Further, from (26) and (23), we know that 𝝀_ℎ=𝟎, 𝝈_ℎ=𝟎. The proof is completed.

4. Some Lemmas

In the study of parabolic equations, we usually introduce a mixed elliptic projection associated with our equations. Define a map: find (𝑢_ℎ,𝝀_ℎ,𝝈_ℎ)∈𝑊_ℎ×Λ_ℎ×𝐕_ℎ, such that div𝝈−𝝈_ℎ,𝑤=0,∀𝑤∈𝑊_ℎ𝝀,0<𝑡≤𝑇,𝝀−_ℎ,𝐯−𝑢_𝑡−𝑢_ℎ,𝑡+𝑏𝑢−𝑢_ℎ,div𝐯+𝐜𝑢−𝑢_ℎ,𝐯=0,∀𝐯∈𝐕_ℎ𝝈,0<𝑡≤𝑇,𝜎−_ℎ,𝝉−𝑎𝝀−𝝀_ℎ,𝝉=0,∀𝝉∈Λ_ℎ𝑢(0)−𝑢,0<𝑡≤𝑇,_ℎ(0),𝑤=0,∀𝑤∈𝑊_ℎ. (29) Similarly to Lemma 1, we can prove that system (29) has a unique solution.

Now we give some error estimates of (𝑢_ℎ,𝝀_ℎ,𝝈_ℎ). Define 𝑢−𝑢_ℎ=𝑝_ℎ𝑢−𝑢_ℎ+𝑢−𝑝_ℎ𝑢=𝑢₁+𝑢₂,𝝀₁=𝝀−𝝀_ℎ,𝝈₁=𝝈−𝝈_ℎ. (30) System (29) can be rewritten as follows: div𝝈₁,𝑤=0,∀𝑤∈𝑊_ℎ,𝝀₁,𝐯−𝑢_1,𝑡+𝑏𝑢₁+𝑢₂,div𝐯+𝐜𝑢₁+𝑢₂,𝐯=0,∀𝐯∈𝐕_ℎ,𝝈₁,𝝉−𝑎𝝀₁,𝝉=0,∀𝝉∈Λ_ℎ.𝑢(0)−𝑢_ℎ(0),𝑤=0,∀𝑤∈𝑊_ℎ. (31) Now we consider the estimates of 𝑢₁ and 𝑢_1,𝑡.

Lemma 2. Let (𝑢₁,𝝀₁,𝝈₁) be the solution of system (31). If 𝑢, 𝝀, and 𝝈 are sufficiently smooth, then there exist positive constants 𝐶 such that 𝑢₁≤𝐶^𝑡₀ℎ𝝀₁+ℎ^{1−𝛿𝑘0}𝝀₁+ℎ𝑢₂+𝑢₂₋₁𝑢𝑑𝑡,_1,𝑡≤𝐶ℎ𝝀₁+ℎ^{1−𝛿𝑘0}𝝀₁+ℎ𝑢₂+𝑢₂₋₁+𝑢₁, (32) where 𝛿_𝑘0=0 for 𝑘=0, and 𝛿_𝑘0=0 for 𝑘≥1.

Proof. Let 𝜙∈𝐻²(Ω)⋂𝐻¹₀(Ω) be the solution of the following problem: ∇⋅𝑎∇𝜙=𝜓,∀𝑥∈Ω,𝜙|_𝜕Ω=0. (33) Then we know 𝜙₂≤𝐶𝜓. (34) For 0<𝑡≤𝑇, from the second equation of (31), we have that 𝑢_1,𝑡,𝜓=𝑢_1,𝑡=𝑢,div𝑎∇𝜙_1,𝑡,divΠ_ℎ=𝜆𝑎∇𝜙₁,Π_ℎ𝑎∇𝜙+𝐜𝑢₁+𝑢₂,Π_ℎ−𝑏𝑢𝑎∇𝜙₁+𝑢₂,divΠ_ℎ𝑎∇𝜙=𝐼₁+𝐼₂+𝐼₃. (35) It is easy to see that 𝐼₁=𝝀₁,Π_ℎ=𝝀𝑎∇𝜙₁,Π_ℎ𝑎∇𝜙−𝑎∇𝜙+𝝀₁=𝝀,𝑎∇𝜙₁,Π_ℎ𝑎∇𝜙−𝑎∇𝜙+𝑎𝝀₁,∇𝜙=𝐸₁+𝐸₂,𝐸₁≤𝐶ℎ𝝀₁𝜙₂≤𝐶ℎ𝝀₁𝜓,𝐸₂=𝑎𝝀₁,∇𝜙−𝑝_ℎ𝜙+𝑎𝝀₁,∇𝑝_ℎ𝜙=𝑎𝝀₁,∇𝜙−𝑝_ℎ𝜙+𝝈₁,∇𝑝_ℎ𝜙=𝑎𝝀₁,∇𝜙−𝑝_ℎ𝜙−div𝝈₁,𝑝_ℎ𝜙=𝑎𝝀₁,∇𝜙−𝑝_ℎ𝜙≤𝐶𝝀₁∇𝜙−𝑝_ℎ𝜙≤𝐶ℎ^{1−𝛿𝑘𝑜}𝜙₂𝝀₁≤𝐶ℎ^{1−𝛿𝑘𝑜}𝜓𝝀₁. (36) Now, we estimate 𝐼₂, 𝐜𝑢₁,Π_ℎ𝑢𝑎∇𝜙≤𝐶₁Π_ℎ𝑢𝑎∇𝜙−𝑎∇𝜙+𝑎∇𝜙≤𝐶₁ℎ+𝑢∇𝜙∇𝜙≤𝐶₁𝜓,𝐜𝑢₂,Π_ℎ=𝑎∇𝜙𝐜𝑢₂,Π_ℎ𝑎∇𝜙−𝑎∇𝜙+𝐜𝑢₂ℎ𝑢,𝑎∇𝜙≤𝐶₂+𝑢₂₋₁𝜓. (37) We turn to consider 𝐼₃, −𝑏𝑢₁,divΠ_ℎ𝑢𝑎∇𝜙≤𝐶₁∇𝜙≤𝐶𝑢₁𝜓,−𝑏𝑢₂,divΠ_ℎ𝑎∇𝜙≤−𝑏−𝑝⁰_ℎ𝑏𝑢₂,divΠ_ℎ𝑎∇𝜙−_{𝑒∈𝑇ℎ}𝑝⁰_ℎ𝑏𝑢₂,divΠ_ℎ𝑎∇𝜙_𝑒=−𝑏−𝑝⁰_ℎ𝑏𝑢₂,divΠ_ℎ𝑢𝑎∇𝜙≤𝐶ℎ₂𝜓, (38) where 𝑝⁰_ℎ𝑏 is the piecewise constant interpolation of function 𝑏. Using the previous estimates, we obtain 𝑢_1,𝑡≤𝐶ℎ𝝀₁+ℎ^{1−𝛿𝑘0}𝝀₁+ℎ𝑢₂+𝑢₂₋₁+𝑢₁. (39) So we deserve 𝑢₁=^𝑡₀𝑢_1,𝑡𝑑𝑡≤𝐶^𝑡₀𝑢_1,𝑡𝑑𝑡≤𝐶^𝑡₀ℎ𝝀₁+ℎ^{1−𝛿𝑘0}𝝀₁+ℎ𝑢₂+𝑢₂₋₁+𝑢₁𝑑𝑡. (40) Applying Gronwall’s inequality to (40), we have 𝑢₁≤𝐶^𝑡₀ℎ𝝀₁+ℎ^{1−𝛿𝑘0}𝝀₁+ℎ𝑢₂+𝑢₂₋₁𝑑𝑡. (41) Noting the estimates of 𝐼₁, 𝐼₂, and 𝐼₃, the proof is completed.

Define ‖𝑢,𝝀,𝝈‖_𝑟=‖𝑢‖_𝑟+‖𝝀‖_𝑟+‖𝝈‖_𝑟,|‖𝑢,𝝀,𝝈‖|_𝑟=^𝑡₀‖𝑢‖_𝑟+‖𝝀‖_𝑟+‖𝝈‖_𝑟𝑑𝑡. (42)

Lemma 3. Let (𝑢,𝝀,𝝈) and (𝑢_ℎ,𝝀_ℎ,𝝈_ℎ) be the solution of (19) and (29), respectively. If 𝑢, 𝝀, and 𝝈 are sufficiently smooth, then there exist positive constants 𝐶 such that(I)if 𝑘≥0, then div𝝈−𝝈_ℎ≤𝐶ℎ^𝑟‖𝝈‖_𝑟+1𝑝,0≤𝑟≤𝑘+1,_ℎ𝑢−𝑢_ℎ≤𝐶ℎ^{𝑟+1−𝛿𝑘0}|‖𝑢,𝝀,𝝈‖|_𝑟𝑝,1≤𝑟≤𝑘+1,_ℎ𝑢_𝑡−𝑢_ℎ,𝑡≤𝐶ℎ^{𝑟+1−𝛿𝑘0}‖𝑢,𝝀,𝝈‖_𝑟,1≤𝑟≤𝑘+1, (43)(44)(45)(II)if 𝑘=0, then 𝑢−𝑢_ℎ≤𝐶ℎ‖𝑢‖₁+|‖𝑢,𝝀,𝝈‖|₁,𝑢_𝑡−𝑢_ℎ,𝑡≤𝐶ℎ𝑢_𝑡1+‖𝑢,𝝀,𝝈‖₁+|‖𝑢,𝝀,𝝈‖|₁,𝝀−𝜆_ℎ≤𝐶ℎ‖𝑢,𝝀,𝝈‖₁+|‖𝑢,𝝀,𝝈‖|₁𝝈,𝝈−_ℎ≤𝐶ℎ‖𝑢,𝝀,𝝈‖₁+|‖𝑢,𝝀,𝝈‖|₁, (46)(47)(48)(49)(III)if 𝑘≥1,2≤𝑟≤𝑘+1, then 𝑢−𝑢_ℎ≤𝐶ℎ^𝑟‖𝑢‖_𝑟+|‖𝑢,𝝀,𝝈‖|_𝑟−1,𝑢_𝑡−𝑢_ℎ,𝑡≤𝐶ℎ^𝑟𝑢_𝑡𝑟+‖𝑢,𝝀,𝝈‖_𝑟+|‖𝑢,𝝀,𝝈‖|_𝑟−1,𝝀−𝝀_ℎ≤𝐶ℎ^𝑟‖𝑢,𝝀,𝝈‖_𝑟+|‖𝑢,𝝀,𝝈‖|_𝑟−1𝝈,𝝈−_ℎ≤𝐶ℎ^𝑟‖𝑢,𝝀,𝝈‖_𝑟+|‖𝑢,𝝀,𝝈‖|_𝑟−1. (50)(51)(52)(53)

Proof. For 0<𝑡≤𝑇, the proof proceeds in three steps as follows.
(I) In the first equation of (29), taking 𝜔=div𝝈₁, we know div𝝈₁,div𝝈₁=div𝝈₁,div𝝈−Π_ℎ𝝈≤div𝝈₁div𝝈−Π_ℎ𝝈, (54) which implies that div𝝈₁≤div𝝈−Π_ℎ𝝈≤𝐶ℎ^𝑟‖𝝈‖_𝑟+1,0≤𝑟≤𝑘+1. (55) Choosing 𝜔∈𝑊_ℎ, note that (div(𝝈−𝝈_ℎ),𝜔)=0 and (div(𝝈−Π_ℎ𝝈),𝜔)=0, we get divΠ_ℎ𝝈−𝝈_ℎ=0. (56) Combing (15) with (56), we can prove (43).
(II) In (29), letting 𝐯=Π_ℎ𝝈−𝝈_ℎ, together with (43), we obtain 𝝀−𝝀_ℎ,Π_ℎ𝝈−𝝈_ℎ−𝐜𝑢−𝑢_ℎ,Π_ℎ𝝈−𝝈_ℎ=0. (57) Hence 𝑅_ℎ𝝀−𝝀_ℎ,Π_ℎ𝝈−𝝈_ℎ−𝐜𝑢−𝑢_ℎ,Π_ℎ𝝈−𝝈_ℎ=0. (58) In the third equation of (29), choosing 𝝉=𝑅_ℎ𝝀−𝝀_ℎ, we get 𝝈−𝝈_ℎ,𝑅_ℎ𝝀−𝝀_ℎ=𝑎𝝀−𝝀_ℎ,𝑅_ℎ𝝀−𝝀_ℎ. (59) So we know 𝝈−Π_ℎ𝝈,𝑅_ℎ𝝀−𝝀_ℎ+Π_ℎ𝝈−𝝈_ℎ,𝑅_ℎ𝝀−𝝀_ℎ=𝑎𝝀−𝝀_ℎ,𝑅_ℎ𝝀−𝝀_ℎ. (60) Note that (44) and (45), we obtain Π_ℎ𝝈−𝝈_ℎ,𝑅_ℎ𝝀−𝝀_ℎ=𝐜𝑢_ℎ−𝑢_ℎ,Π_ℎ𝝈−𝝈_ℎ=𝑎𝝀−𝝀_ℎ,𝑅_ℎ𝝀−𝝀_ℎ−𝝈−Π_ℎ𝝈,𝑅_ℎ𝝀−𝝀_ℎ, (61) which implies that Π_ℎ𝝈−𝝈_ℎ,𝑅_ℎ𝝀−𝝀_ℎ≤𝐶𝑢−𝑢_ℎΠ_ℎ𝝈−𝝈_ℎ≤𝐶𝑢₁+𝑢₂Π_ℎ𝝈−𝝈_ℎ. (62) Now we estimate ‖Π_ℎ𝝈−𝝈_ℎ‖, 𝝈−𝝈_ℎ,Π_ℎ𝝈−𝝈_ℎ=𝑎𝝀−𝝀_ℎ,Π_ℎ𝝈−𝝈_ℎ. (63) It is easy to see that Π_ℎ𝝈−𝝈_ℎ²=Π_ℎ𝝈−𝝈,Π_ℎ𝝈−𝝈_ℎ+𝑎𝝀−𝝀_ℎ,Π_ℎ𝝈−𝝈_ℎ=Π_ℎ𝝈−𝝈,Π_ℎ𝝈−𝝈_ℎ+𝑎𝝀−𝑅_ℎ𝝀,Π_ℎ𝝈−𝝈_ℎ+𝑎𝑅_ℎ𝝀−𝝀_ℎ,Π_ℎ𝝈−𝝈_ℎ≤𝐶Π_ℎ𝝈−𝝈+𝝀−𝑅_ℎ𝝀+𝑢₁+𝑢₂×Π_ℎ𝝈−𝝈_ℎ. (64) By the properties of projection Π_ℎ, we get Π_ℎ𝝈−𝝈_ℎ≤𝐶ℎ^𝑟‖𝑢‖_𝑟+‖𝝀‖_𝑟+‖𝝈‖_𝑟+𝑢₁,0≤𝑟≤𝑘+1, (65) which proves (49) and (53).
Now we estimate ‖𝑅_ℎ𝝀−𝝀_ℎ‖, 𝑎𝑅_ℎ𝝀−𝝀_ℎ,𝑅_ℎ𝝀−𝝀_ℎ=𝑎𝝀−𝝀_ℎ,𝑅_ℎ𝝀−𝝀_ℎ+𝑎𝑅_ℎ𝝀−𝝀,𝑅_ℎ𝝀−𝝀_ℎ=𝝈−Π_ℎ𝝈,𝑅_ℎ𝝀−𝝀_ℎ+Π_ℎ𝝈−𝝈_ℎ,𝑅_ℎ𝝀−𝝀_ℎ+𝑎𝑅_ℎ𝝀−𝝀,𝑅_ℎ𝝀−𝝀_ℎ=𝝈−Π_ℎ𝝈,𝑅_ℎ𝝀−𝝀_ℎ+𝐜𝑢−𝑢_ℎ,Π_ℎ𝝈−𝝈_ℎ+𝑎𝑅_ℎ𝝀−𝝀,𝑅_ℎ𝝀−𝝀_ℎ≤𝝈−Π_ℎ𝝈+𝐶𝝀−𝑅_ℎ𝝀𝑅_ℎ𝝀−𝝀_ℎ+𝐶𝑢₁+𝑢₂Π_ℎ𝝈−𝝈_ℎ≤𝐶𝝈−Π_ℎ𝝈²+𝑢₁²+𝑢₂²+𝑅_ℎ𝝀−𝝀²+𝑐₀2𝑅_ℎ𝝀−𝝀_ℎ², (66) where 𝑎𝑅_ℎ𝝀−𝝀_ℎ,𝑅_ℎ𝝀−𝝀_ℎ≥𝑐₀𝑅_ℎ𝝀−𝝀_ℎ². (67) From (66) and (67), we know 𝑅_ℎ𝝀−𝝀_ℎ≤𝐶ℎ^𝑟‖𝑢‖_𝑟+‖𝝀‖_𝑟+‖𝝈‖_𝑟+𝑢₁,𝜆−𝜆_ℎ≤𝐶ℎ^𝑟‖𝑢‖_𝑟+‖𝜆‖_𝑟+‖𝝈‖_𝑟+𝑢₁,0≤𝑟≤𝑘+1. (68) By Lemma 2, if 𝑘=0, we have that 𝝀−𝝀_ℎ=𝝀₁≤𝐶ℎ‖𝑢‖₁+‖𝝀‖₁+‖𝝈‖₁+^𝑡₀‖𝑢‖₁+‖𝝀‖₁+‖𝝈‖₁𝑑𝜏. (69) If 𝑘≥1,2≤𝑟≤𝑘+1, we have that 𝝀−𝝀_ℎ≤𝐶ℎ^𝑟‖𝑢‖_𝑟+‖𝝀‖_𝑟+‖𝝈‖_𝑟+^𝑡₀‖𝑢‖_𝑟−1+‖𝝀‖_𝑟−1‖𝝈‖_𝑟−1𝑑𝑡. (70) Using the estimate of ‖Π_ℎ𝝈−𝝈_ℎ‖, we obtain 𝝈−𝝈_ℎ≤𝐶ℎ^𝑟‖𝑢‖_𝑟+‖𝝀‖_𝑟+‖𝝈‖_𝑟+𝑢₁,0≤𝑟≤𝑘+1. (71) If 𝑘=0, we get 𝝈−𝝈_ℎ≤𝐶ℎ‖𝑢‖₁+‖𝝀‖₁+‖𝝈‖₁+^𝑡₀‖𝑢‖₁+‖𝝀‖₁+‖𝝈‖₁𝑑𝑡. (72) If 𝑘≥1,2≤𝑟≤𝑘+1, we get 𝝈−𝝈_ℎ≤𝐶ℎ^𝑟‖𝑢‖_𝑟+‖𝜆‖_𝑟+‖𝝈‖_𝑟+^𝑡₀‖𝑢‖_𝑟−1+‖𝝀‖_𝑟−1+‖𝝈‖_𝑟−1𝑑𝑡. (73)
(III) By Lemma 2, we have that 𝑝_ℎ𝑢−𝑢_ℎ=𝑢₁≤𝐶ℎ^{𝑟+1−𝛿𝑘0}×^𝑡₀‖𝑢‖_𝑟+‖𝝀‖_𝑟+‖𝝈‖_𝑟𝑝𝑑𝜏,0≤𝑟≤𝑘+1,_ℎ𝑢_𝑡−𝑢_ℎ,𝑡≤𝐶ℎ^{𝑟+1−𝛿𝑘0}‖𝑢‖_𝑟+‖𝝀‖_𝑟+‖𝝈‖_𝑟+^𝑡₀‖𝑢‖_𝑟+‖𝝀‖_𝑟+‖𝝈‖_𝑟𝑑𝑡,0≤𝑟≤𝑘+1. (74) If 𝑘=0, we know 𝑢−𝑢_ℎ≤𝑢−𝑝_ℎ𝑢+𝑢₁≤𝐶ℎ‖𝑢‖₁+^𝑡₀‖𝑢‖₁+‖𝝀‖₁+‖𝝈‖₁,𝑢𝑑𝑡_𝑡−𝑢_ℎ,𝑡≤𝐶ℎ𝑢_𝑡1+‖𝑢‖₁+‖𝝀‖₁+‖𝝈‖₁+^𝑡₀‖𝑢‖₁+‖𝝀‖₁+‖𝝈‖₁𝑑𝑡. (75) If 𝑘≥1,2≤𝑟≤𝑘+1, we know 𝑢−𝑢_ℎ≤𝐶ℎ^𝑟‖𝑢‖_𝑟+^𝑡₀‖𝑢‖_𝑟−1+‖𝝀‖_𝑟−1+‖𝝈‖_𝑟−1,𝑢𝑑𝑡_𝑡−𝑢_ℎ,𝑡≤𝐶ℎ^𝑟𝑢_𝑡𝑟+‖𝑢‖_𝑟−1+‖𝝀‖_𝑟−1+^𝑡₀‖𝑢‖_𝑟−1+‖𝝀‖_𝑟−1+‖𝝈‖_𝑟−1𝑑𝑡. (76) The proof is completed.

Remark 4. The estimate results of (46)–(53) are optimal. If 𝑘≥1, the estimate results of (43)–(45) are superconvergent.

5. Main Result

In this section, we consider error estimates for the continuous-in-time mixed finite element approximation. Define 𝑢−𝑢_ℎ=𝑢−𝑢_ℎ+𝑢_ℎ−𝑢_ℎ=𝑢₃+𝑢₄,𝝀−𝝀_ℎ=𝝀−𝝀_ℎ+𝝀_ℎ−𝝀_ℎ=𝝀₁+𝝀₂,𝝈−𝝈_ℎ=𝝈−𝝈_ℎ+𝝈_ℎ−𝝈_ℎ=𝝈₁+𝝈₂,𝑢,𝑢_𝑡,𝝀,𝝈_𝑟=‖𝑢‖_𝑟+𝑢_𝑡𝑟+‖𝝀‖_𝑟+‖𝝈‖_𝑟,𝑢,𝑢_𝑡,𝝀,𝝈_𝑟=^𝑡₀‖𝑢‖_𝑟+𝑢_𝑡𝑟+‖𝝀‖_𝑟+‖𝝈‖_𝑟𝑑𝑡. (77)

Theorem 5. Let (𝑢,𝝀,𝝈) and (𝑢_ℎ,𝝀_ℎ,𝝈_ℎ) be the solution of (9) and (11), respectively. If 𝑢, 𝝀, and 𝝈 are sufficiently smooth, then there exist positive constants 𝐶 such that(I)if 𝑘=0, then 𝑢−𝑢_ℎ≤𝐶ℎ‖𝑢‖₁+𝑢,𝑢_𝑡,𝝀,𝝈₁,𝑢_𝑡−𝑢_ℎ,𝑡≤𝐶ℎ𝑢,𝑢_𝑡,𝝀,𝝈₁+|‖𝑢,𝝀,𝝈‖|₁,𝝀−𝝀_ℎ≤𝐶ℎ𝑢,𝑢_𝑡,𝝀,𝝈₁+𝑢,𝑢_𝑡,𝝀,𝝈₁,𝝈−𝝈_ℎ≤𝐶ℎ𝑢,𝑢_𝑡,𝝀,𝝈₁+𝑢,𝑢_𝑡,𝝀,𝝈₁,div𝝈−𝝈_ℎ≤𝐶ℎ‖𝑢‖₁+𝑢_𝑡1+‖𝝈‖₂+𝑢,𝑢_𝑡,𝝀,𝝈₁, (78)(79)(80)(81)(82)(II)if 𝑘≥1 and 2≤𝑟≤𝑘+1, then 𝑢−𝑢_ℎ≤𝐶ℎ^𝑟‖𝑢‖_𝑟+𝑢,𝑢_𝑡,𝝀,𝝈_𝑟−1,𝑢_𝑡−𝑢_ℎ,𝑡≤𝐶ℎ^𝑟𝑢_𝑡𝑟+‖𝑢‖_𝑟−1+‖𝝀‖_𝑟−1+‖𝝈‖_𝑟−1+^𝑡₀‖𝑢‖_𝑟−1+𝑢_𝑡𝑟+‖𝝀‖_𝑟−1+‖𝝈‖_𝑟−1𝑑𝑡,𝝀−𝝀_ℎ≤𝐶ℎ^𝑟‖𝑢,𝝀,𝝈‖_𝑟+^𝑡₀𝑢_𝑡𝑟+‖𝑢‖_𝑟−1+‖𝝀‖_𝑟−1+‖𝝈‖_𝑟−1𝑑𝑡,𝝈−𝝈_ℎ≤𝐶ℎ^𝑟‖𝑢,𝝀,𝝈‖_𝑟+^𝑡₀𝑢_𝑡𝑟+‖𝑢‖_𝑟−1+‖𝝀‖_𝑟−1+‖𝝈‖_𝑟−1𝑑𝑡,div𝝈−𝝈_ℎ≤𝐶ℎ^𝑟𝑢_𝑡𝑟+‖𝝀‖_𝑟+‖𝝈‖_𝑟+1+‖𝑢‖_𝑟−1+^𝑡₀𝑢_𝑡𝑟+‖𝑢‖_𝑟−1+‖𝝀‖_𝑟−1+‖𝝈‖_𝑟−1𝑑𝑡. (83)(84)

Proof. From (9) and (19), we have the following error equation: 𝑢_4,𝑡+𝑢_3,𝑡,𝑤+div𝝈₂,𝑤=0,∀𝑤∈𝑊_ℎ𝝈,0<𝑡≤𝑇,₂,𝝉−𝑎𝝀₂,𝝉=0,∀𝝉∈Λ_ℎ𝝀,0<𝑡≤𝑇,₂,𝐯−𝑢_4,𝑡+𝑏𝑢₄,div𝐯+𝐜𝑢₄,𝐯=0,∀𝐯∈𝐕_ℎ𝑢,0<𝑡≤𝑇,₄(0)=0. (85) In (85), choosing 𝜔=𝑢_4,𝑡+𝑝_ℎ(𝑏𝑢₄), 𝝉=𝝀₂, and 𝐯=𝝈₂, we have 𝝀₂,𝝈₂+𝑢_4,𝑡,𝑢_4,𝑡=−𝑢_3,𝑡,𝑢_4,𝑡+𝑝_ℎ𝑏𝑢₄−𝑢_4,𝑡,𝑝_ℎ𝑏𝑢₄−𝐜𝑢₄,𝝈₂. (86) In the second equation of (85), letting 𝝉=𝝈₂, we know 𝝈₂≤𝐶𝝀₂. (87) Further, using 𝜀-inequality to (86), we get 𝑢_4,𝑡²+𝝀₂²≤𝐶𝑢₄²+𝑢_3,𝑡². (88) Hence, 𝑢₄=^𝑡₀𝑢_4,𝑡𝑑𝜏≤𝐶^𝑡₀𝑢_4,𝑡𝑑𝜏≤𝐶^𝑡₀𝑢₄+𝑢_3,𝑡𝑑𝑡. (89) Using Gronwall’s inequality to (89), we obtain 𝑢₄≤𝐶^𝑡₀𝑢_3,𝑡𝑑𝑡. (90) Further, we know 𝝈₂+𝑢_4,𝑡+𝝀₂≤𝐶𝑢_3,𝑡+^𝑡₀𝑢_3,𝑡𝑑𝑡. (91) Noting the first equation of (78), we get div𝝈₂≤𝑢_3,𝑡+𝑢_4,𝑡≤𝐶𝑢_3,𝑡+^𝑡₀𝑢_3,𝑡𝑑𝑡. (92) So we have that 𝑢−𝑢_ℎ≤𝑢₃+𝐶^𝑡₀𝑢_3,𝑡𝑢𝑑𝑡,_𝑡−𝑢_ℎ,𝑡≤𝐶𝑢_3,𝑡+^𝑡₀𝑢_3,𝑡𝑑𝑡,𝝀−𝝀_ℎ≤𝝀₁+𝐶𝑢_3,𝑡+^𝑡₀𝑢_3,𝑡𝑑𝑡,𝝈−𝝈_ℎ≤𝝈₁+𝐶𝑢_3,𝑡+^𝑡₀𝑢_3,𝑡𝑑𝑡,div𝝈−𝝈_ℎ≤div𝝈₁+div𝝈₂. (93) Together with the results of Lemmas 2 and 3, the proof is completed.

Remark 6. The estimate results of (78)–(84) are optimal.

Acknowledgments

This work is supported by the National Natural Science Foundation of China (nos. 11171190, 11101246, and 11101247) and the Natural Science Foundation of Shandong Province (ZR2011AQ020). The authors thank the anonymous referees for their constructive comments and suggestions, which led to improvements in the presentation.