Abstract

We propose the characteristics-mixed method for approximating the solution to a convection-dominated transport problem. The new method is a combination of characteristic approximation to handle the convection part in time and a mixed finite element spatial approximation to deal with the diffusion part. Boundary conditions are incorporated in a natural fashion. The scheme is locally conservative since fluid is transported along the approximate characteristics on the discrete level and the test function can be piecewise constant. Analysis shows that the scheme has much smaller time-truncation errors than those of standard methods.

1. Introduction

Given an open bounded domain with a smooth boundary and a time interval , we consider in this paper the following convection-diffusion equation: where(1) denotes, for example, the concentration of a possible substance; (2) represents the velocity of the flow;(3) and denote the gradient and the divergence operators, respectively;(4) is sufficiently smooth and there exist constants and such that ;(5) denotes a source term.

This equation governs such phenomena as the flow of heat within a moving fluid, the transport of dissolved nutrients or contaminants within the groundwater, and the transport of a surfactant or tracer within an incompressible oil in a petroleum reservoir.

For convenience, we assume (1) is -periodic, this is physically reasonable since no-flow boundaries are generally treated by reflection, and because in general interior flow patterns are much more important than boundary effects. Thus, the no-flow boundary condition above can be dropped. Because of molecular diffusion, is uniformly positive. Although this implies that the equation is uniformly parabolic, in many applications the Peclet number is quite high. Thus, convection dominates diffusion and the equation is nearly hyperbolic in nature. The concentration often develops sharp fronts that are nearly shocking. It is well known that strictly parabolic discretization schemes applied to the problem do not work well when it is convection dominated. It is especially difficult to approximate well the sharp fronts and conserve the material or mass in the system.

Effective discretization schemes should recognize to some extent the hyperbolic nature of the equation. Many such schemes have been developed, such as the explicit method of characteristics, upstream-weighted finite difference schemes [1], higher-order Godunov schemes [2, 3], the streamline diffusion method [4], the least-squares-mixed finite element method [5, 6], the modified method of characteristics-Galerkin finite element procedure (MMOC-Galerkin) [79], and the Eulerian-Lagrangian localized adjoint method (ELLAM) [10]. Each method has its advantages and disadvantages. Explicit characteristic and Godunov schemes require that a Courant-Friedrich-Lewy (CFL) time step constraint is to be imposed. Upstream weighting tends to introduce into the solution an excessive amount of numerical diffusion near the sharp fronts. Compared with upstream weighting, the streamline diffusion method and the least-squares-mixed finite element method reduce the amount of numerical diffusion. It adds a user-defined amount of bias in the direction of streamline. In ELLAM, it can be difficult to evaluate the resulting integrals.

We concentrate on MMOC-Galerkin. It is an implicit scheme, so reasonably large time steps may be used, and it does not numerically diffuse the fronts to a particularly excessive degree. Unfortunately, it has certain inherent difficulties, especially with regard to local mass balance.

Mixed finite element method has been found to be very useful for solving flow equations not only for its optimal approximation of the unknown function and the vector flux, but also for its mass conservation due to the use of piecewise test function. Its theories and applications have been extensively discussed, for example, see [1114].

Therefore, to purse the high performance of an algorithm, one should take a nice combination of MMOC and the mixed finite element method. In this way, characteristics-mixed finite element method was considered. However, most integrals in this method involve the transformation images of the test function, so it is difficult to compute in practice.

In this paper, we propose a new characteristics-mixed finite element scheme. It is similar to MMOC-Galerkin in that we approximate the hyperbolic part of the equation along the characteristics. We use, however, a mixed finite element spatial discretization of the equation. Since piecewise constant functions can be used as test functions, the scheme conserves mass locally and globally in the discrete level. One of our goals in the paper is to yield in time truncation. Since the solution changes much less rapidly in the characteristic direction, the bound in time truncation will permit the use of larger time steps, with corresponding improvements in efficiency, at no cost in accuracy. For we use the mixed finite element, we can simultaneously approximate the scalar unknown and the vector flux. In error analysis, we introduce some new techniques to overcome difficulties in dealing with the combination of MMOC and the mixed finite element. The optimal error estimates for the unknown function and the vector flux are given.

An outline of the paper is as follows. In Section 2, we define an approximation to the characteristics and the characteristics-mixed finite element formulation of the problem. We give the proof of the existence and uniqueness of the discrete problem in Section 3. Several important lemmas are given in Section 4, the critical argument in this section is the proof that for any , is bounded by a constant multiple of . The optimal order estimates for in are presented in Section 5.

2. The Characteristics-Mixed Finite Element Formulation

We begin this section by introducing some notations. We denote by the standard Sobolev space of -differential functions in . Let be its norm and let be the norm of or , where we omit if . When , we let denote the corresponding space defined on and let denote its norm.

We also use the following spaces that incorporate time dependence. If is a Sobolev space and suitably smooths on , we let where if , the integral is replaced by the essential supremum.

We list the assumptions about the coefficients as follows: here and throughout this paper, denotes different constants in different places.

We also assume that the solution of (1) satisfies the following: Under the above assumptions, we begin to discretize the problem. Let and let the characteristic direction associated with the operator be denoted by , where then, (1)(a) can be put in the form the following:

Define the spaces: , and let , then the mixed variational problem corresponding to (1) is to find , such that where and it is this form that will be discretized below.

Remark 1. The restriction on the space can be removed if the Dirichlet boundary condition is imposed, that is, the results presented in the following section also hold in the case of the first or the second type boundary value problems.
Let be a quasi-regular polygonaization of and let be the associated Raviart-Thomas-Nedelec space [15] of index . In the procedure to be used, we will consider a time step and approximate the solution at times . The characteristic derivative will be approximated basically in the following manner.
Let and note that so, our characteristics-mixed finite element method is the determination of the map satisfying the relations following: where , and , are defined in (18) for .

Remark 2. This scheme is an implicit scheme and it is difficult to compute in practice. For simplicity, we can use instead of . Then and (10) becomes an explicit scheme, but the following results still hold.

3. The Existence and Uniqueness of the Solution of the Discrete Problem

In this section, we give the proof of the existence and uniqueness of the solution of the discrete problem (10).

Theorem 3. There exists a unique solution to the characteristics-mixed finite element procedure (10).

Proof. Let and be two sets of bases, respectively, and let then (10) may be written in the following matrix form.
Find such that where Since and are bases, respectively, and , and are two symmetric positive definite matrices. Solve from (13)(b), then substitute it into (13)(a), we have Then, where Since and are globally Lipschitz continuous about , is globally Lipschitz continuous about too. With given by (10)(c) it is clear that this nonlinear system of algebraic equations may be solved for small , so that (10) defines a unique discrete solution for .

In order to analysis the convergence of the method, it is convenient to introduce the mixed elliptic projection associated with our equations.

Let be given by the following relations

Let , and . By [11] we know that for or and , consider It follows easily that We also have By (8), (10), and (18), we obtain the error equation in the following form

4. Several Lemmas

To obtain the optimal error estimates, we introduce the following three basic lemmas. These three lemmas are crucial to our main arguments.

Lemma 4. For any function , there exists a function , such that where is a constant.

Lemma 5. For any , if is satisfying then, there exists a constant such that

Lemma 6. Let , . If satisfying Then,

Proof. Let . By Lemma 5 we have To obtain the desired estimate, what we should deduce is to demonstrate the following: Thus, consider About the transformation , because of the smoothness and periodicity of is a differentiable homeomorphism of onto itself for sufficiently small (see [9]). By using the Jacobian of this transformation and a change of variables, we have It is easy to obtain Next, let be the unit vector in the direction of ; then
We now obtain (30) by combining (33), (34), and (36). A similar change of variables demonstrates that Then it is easy to obtain

5. Error Estimates

Under the above assumptions about , , , and , we can derive the optimal order estimates of and .

Theorem 7. Let denote the solution of (10) and (8), respectively. Suppose that , we have for sufficiently small the following:

Proof. Lets first show the estimate of (38). Taking and in (23), we obtain
For the first term on the right-hand side of the above equation, following the treatment manner in [8], we will bound through the representation below involving an integral in the parameter along the tangent to the characteristics from to . Denote the coordinates of the point on the segment by . The standard backward difference quotient error equation is given by analogously, along the tangent to the characteristics following: Taking the norm of this error term, we obtain to relate this to a standard norm of , consider the following transformation: By the assumption of (3), is invertible for sufficiently small and its determinant is . Since our problem is periodic, for any fixed is a homeomorphism of onto itself, so the same is true for . It follows that and the first term on the right-hand side of (40) is bounded by Concerning the second term on the right-hand side of (40), write as the sum , then, by Lemma 6 and (23)(b), under the assumption , we have Now we come to the term , by (4) we have the following: About the last term in the right-hand side of (40), by (4), we have and this completes the treatment of the right-hand side of (40).
By using (33), the left-hand side of (40) is bounded by , combining (46)–(49) with (40) to give the following recursion relation: If (50) is multiplied by and summed in time, and if it is noted that (10)(c) implies that , then, it follows that For sufficiently small , by Gronwall's lemma, we obtain Recall that the error , (19), (21), and (53) together imply that This completes the first part of our conclusion.
In order to show the estimate (39), by (23)(b), we get Choosing and , respectively, in (23)(a) and (54) and adding them to obtain The left-hand side satisfies the inequality: Substitute it into (55), we have Bound on the right-hand side of (56) one by one, we have Concerning and , we can easily get Substituting (58)-(59) into (57) to obtain Noting that , we have from (23)(b). To multiply (60) by and add it in time, we obtain By (23)(b) and Lemma 6, we have where the assumption is used.
Similarly, we get Substituting (62) and (63) into (61), we get
We see that an application of Gronwall's lemma would complete the second part of our argument if we did not have the term on the right hand side of (64). We will use an induction argument to prove the hypothesis as follows: Since , by the inverse property and (19), it can be easily verified that . Using (19), (21), (52), and (64) we see that (65) clearly holds for if is sufficiently small. Assume that (65) holds for , using (64), we obtain Then, an application of Gronwall's lemma gives us an estimate which together with (19), (21), and (52) and the inverse property of , show that the induction hypothesis (65) holds for for sufficiently small . This completes the induction argument.
Finally by Gronwall's lemma, it follows from (64) that Note that , by (52) and (19)–(21), we obtain This also shows that

Note that our results in this paper do not cover the cases of nonlinear convection-dominated systems, which are of importance in some applications, particularly in reservoir simulations.

Acknowledgment

This research is partially supported by the China Postdoctoral Science Foundation funded project (2011M501149), the Humanity and Social Science Foundation of the Ministry of Education of China (12YJCZH303), the Special Fund Project for Post Doctoral Innovation of Shandong Province (201103061), and Independent Innovation Foundation of Shandong University, IIFSDU (IFW12109).