## Services on Demand

## Journal

## Article

## Indicators

- Cited by SciELO
- Access statistics

## Related links

- Cited by Google
- Similars in SciELO
- Similars in Google

## Share

## Computational & Applied Mathematics

##
*On-line version* ISSN 1807-0302

### Comput. Appl. Math. vol.31 no.3 São Carlos 2012

#### http://dx.doi.org/10.1590/S1807-03022012000300004

**A weighted mass explicit scheme for convection-diffusion equations **

**Vitoriano Ruas ^{*}**

Research associate, Departamento de Engenharia Mecânica, PUC-Rio - Catholic University of Rio de Janeiro, Brazil. Research collaborator, Institut Jean le Rond d'Alembert, UPMC - Université Pierre et Marie Curie, Paris 6, France. Visiting professor, Instituto de Ciências Matemáticas e Computação, USP - University of São Paulo, São Carlos, Brazil. E-mails: vitoriano.ruas@pq.cnpq.br / vitoriano.ruas@upmc.fr / vitoriano@icmc.usp.br

**ABSTRACT**

An explicit scheme based on a weighted mass matrix, for solving time-dependent convection-diffusion problems was recently proposed by the author and collaborators. Convenient bounds for the time step, in terms of both the method's weights and the mesh step size, ensure its stability in space and time, for piecewise linear finite element discretisations in any space dimension. In this work we study some techniques for choosing the weights that guarantee the convergence of the scheme with optimal order in the space-time maximum norm, as both discretisation parameters tend to zero.

**Mathematical subject classification:** Primary: 65M60; Secondary: 76Rxx.

**Key words:** Convection-diffusion, convergence, explicit scheme, finite elements, piecewise linear, time-dependent, weighted mass.

**1 Introduction**

This work deals with an explicit scheme introduced in [6] , for the numerical time integration of the convection-diffusion equations, discretised in space by techniques based on variational formulations such as the finite element method. In this framework, since the mid-eighties, the most widespread manner to deal with dominant convection has been the use of stabilizing procedures based on the space mesh parameter, among which the streamline upwind Petrov-Galerkin (SUPG) technique introduced by Hughes & Brooks (cf. [2] ) is one of the most popular.

The author and collaborators studied in [7] a contribution in this direction, based on a standard Galerkin approach, and a space discretisation of theconvection-diffusion equations with piecewise linear finite elements, combined with a non standard explicit forward Euler scheme for the time integration. The main theoretical result in that work, states that the numerical solution is stable in the maximum norm in both space and time (and even convergent with order *h*|*lnh*| if the mesh is of the acute type [10] ), provided that roughly the time step is bounded by the space mesh parameter *h* multiplied by a mesh-independent constant, for a high Péclet number. As it should be clarified, the scheme under consideration follows similar principles to the one long exploited by Kawahara and collaborators, for simulating convection dominated phenomena (see e.g. [4] , among several other papers published by them before and later on). The originality of our contribution relies on the fact that it not only introduces a reliable scheme for any space dimension, but also exhibits rigorous conditions for it to provide converging sequences of approximations in the sense of the space-time maximum norm.

The main purpose of this work, is to specify procedures for the determination of sets of weights that characterize our explicit scheme. As for the theoretical contribution, we prove that for such choices of the weights the method converges with optimal order in the maximum space-time norm.

An outine of the paper is as follows: In Section 2 we recall the problem to solve together with the type of discretisation corresponding to our explicit scheme; more particularly the weighted manner to deal with the mass matrices on both sides of the discrete equations is described. In Section 3 we further recall the stability results that hold for the method being considered, in the sense of the space and time maximum norm, together with the conditions to be fulfilled by the sequences of meshes and time steps in order to ensure convergence. We proceed in Section 4 by studying in detail some particular choices of the weights associated with the scheme, that satisfy the conditions leading to both stability and convergence. Finally in Section 5 we consider some implementation aspects of the method and give corresponding numerical results.

**2 The problem to solve and its discretisation**

Let us consider a time-dependent convection-diffusion problem described as follows:

Find a scalar valued function *u*(**x**, *t*) defined in × [0, ∞), Ω being a bounded open subset of ℜ^{N} with boundary ∂Ω, N = 1, 2 or 3, such that,

where *u _{t}* represents the first order derivative of

*u*with respect to

*t*,

*ν*is a positive constant and

**a**is a given solenoidal convective velocity at every time

*t*, assumed to be uniformly bounded in Ω × (0, ∞). The data

*f*and

*g*are respectively, a given forcing function belonging to

*L*

^{∞}[Ω × (0, ∞)] , and a prescribed value in L

^{∞}[∂Ω × (0, ∞)] . We further assume that

*u*

^{0}∈

*L*

^{∞}(Ω) and that for every

**x**∈ Ω

*g*(

**x**, ·) is of bounded variation in (0, ∞). In (1)

*ν*represents the inverse of the Péclet number.

Without loss of essential aspects, in all the sequel we assume that Ω is an interval if *N* = 1, a polygon if *N* = 2 or a polyhedron if *N* = 3. In so doing we consider a partition _{h} of Ω into *N*- simplices, with maximum edge length equal to *h*. We assume that _{h} satisfies the usual compatibility conditions for finite element meshes, and that it belongs to a quasi-uniform family of partitions. We further define a second mesh parameter *h*_{min} as the minimum height of all the elements of _{h} if *N* = 2 or 3 and the minimum length of *K* ∈ _{h} if *N* = 1.

Let *N _{h}* be the number of nodes of

_{h}, denoted by

*P*,

_{j}*j*= 1, 2,...,

*N*. We assume that these are numbered in such a manner that the first

_{h}*I*nodes are located in the interior of Ω and the remaining

_{h}*N*-

_{h}*I*nodes are located on ∂Ω. Now for every

_{h}*K*∈

_{h}we denote by

*P*

_{1}(

*K*) the space of polynomials of degree less than or equal to one defined in

*K*. In so doing we introduce the following spaces or manifolds associated with

_{h}:

We further introduce for any function *ϕ* defined in *C*^{0}(∂Ω) the following manifold of *V _{h}*:

Now let be the field of satisfying (*P _{j}*) =

*u*

^{0}(

*P*) for every vertex

_{j}*P*of

_{j}_{h}, and Δ

*t*> 0 be a given time step. Defining

*g*on ∂Ω by

^{n}*g*(·) =

^{n}*g*(·,

*n*Δ

*t*),

*f*in Ω by

^{n}*f*(·) =

^{n}*f*(·,

*n*Δ

*t*) and

**a**

*in Ω by*

^{n}**a**

^{n}(·) =

**a**(·,

*n*Δ

*t*), for

*n*= 1, 2,..., idealistically we wish to determine approximations (·) of u(·,

*n*Δ

*t*) for

*n*∈

^{*}, by solving the following finite element discrete problem described below, corresponding to a modification of the first order forward Euler scheme.

For *n* successively equal to 1, 2,..., we wish to determine ∈ of the form

where *φ** _{j}* is the canonical basis function of

*V*associated with the

_{h}*j*-th node of

_{h}, i.e.

*P*, ∈ ℜ being the value of at

_{j}*P*. We denote by

_{j}*S*the support of

_{j}*φ*

*and by Π*

_{j}_{j}its measure.

The unknowns for *n* = 1, 2,..., are recursively determined by the following expressions:

where, setting := **a*** ^{n}*(

*P*) and :=

_{i}*f*(

^{n}*P*), the coefficients and are given by:

_{i}Coefficient is the well-known lumped mass diagonal matrix (cf. [9] ) given by . The mass matrix coefficients on the right hand side of (2) in turn are defined by a weighted quadrature formula described as follows.

Let *M _{i}* be the number of nodes different from

*P*lying in the closure of

_{i}*S*, i.e.

_{i}_{i}, and be such nodes for

*j*= 1, 2,...,

*M*with 1

_{i}__<__

*k*

_{j}__<__

*N*. Let also be the measure fractions associated with given by:

_{h}and be corresponding strictly positive weights satisfying,

Notice that since each *N*-simplex in *S _{i}* appears in exactly

*N*measure fractions , we necessarily have:

Now selecting the nodes in _{i} different from *P _{i}*, we define,

together with

where is the *i-th* diagonal coefficient of the standard consistent massmatrix ([9] ) given by

Naturally enough, by definition, = 0 if *P _{j}* does not lie in

_{i}.

Typically we may choose = for every *j* and for every node *P _{i}*, thereby generating a weighted combination of the lumped mass and the consistent mass matrix (cf. [9] ) on the right hand side of (2), with weights equal to

respectively. However, except for the case of uniform meshes, in principle this is not the choice to make, if one wishes to reach the best results in terms of accuracy.

**3 Stability and convergence of the scheme**

In this Section we recall the stability and convergence results proven in [7] that hold for the above defined weighted mass scheme. In short they state that, provided Δ*t* is chosen conveniently small with respect to the spatial mesh parameter, the scheme (2) is stable in the sense of the maximum norm. Moreover under a suitable condition on the mesh it converges in both space and time in the same sense, as *h* and Δ*t* go to zero.

First we have to define the following quantities:

• A = |

a_{i}(t)|;•

ω=

The following theorem proved in [7] states the stability result that holds for scheme (2).

**Theorem 1. ** *If *Δ*t fulfills the condition *

*then the finite element solution sequence ** given by ** generated by *(2)* for n *= 1, 2,...* satisfies the following stability result for every m *∈* , whereby ||F*||_{0, }_{∞}_{, D}* denotes the L*^{∞}*-norm of a function F defined in an open set D of *ℜ^{N}*, and BV*[*G*] * represents the standard norm of a function G*(*t*)* having bounded variation for t *∈ (0, ∞):

The above stability result can be refined as follows, in the particular casewhere the partition _{h} is of the *acute type* (see e.g. [10] ).

**Theorem 2. (cf. [7] ). ***Assume that the partition *_{h}* is of the acute type (cf. [10] ). Then if *Δ*t satisfies the condition*

*the finite element solution sequence ** given by ** generated by *(2)* for n *= 1, 2,...* satisfies the stability condtion *(10).

In [7] we give error estimates for the approximations of the solution of (1) generated by (2) under condition (9). In particular we recall here that, provided the weights are suitably chosen, this scheme provides convergent approximations in the maximum norm, as both *h* and Δ*t* go to zero, under the assumption (11) of Theorem 2. For this purpose we use Sobolev spaces *W ^{m}*

^{,}

^{∞}(

*D*) equipped with the standand norm and seminorm denoted respectively by ||·||

_{m, }

_{∞}

_{, D}and |·|

_{m, }

_{∞}

_{, D}, where

*m*is a non negative integer and

*D*is a subset of ℜ

*(cf [1] ).*

^{N}As usual a suitable consistency result is needed, which together with the stability results given in this Section leads to convergence. Actually the consistency of our scheme is a consequence of the following lemma:

**Lemma 1. (cf. [7] ). ***Let P _{i} be a node of _{h}, for i *∈ {1, 2,...,

*I*}

_{h}*, and*

*be the vector leading from P*

_{i}to its neighbor*, that is, the j-th node belonging to*= 1, 2,...

_{i}, j*, M*

_{i}. Then there exists strictly positive weights*satisfying*(5)

*such that*

Then we can establish the validity of the following convergence result for scheme (2):

**Theorem 2. (cf. [7] ). ***Let the strictly positive weights **, *∀*i *∈ {1, 2,..., *I _{h}*}

*and*∀

*j*∈ {1, 2,...,

*M*}

_{i}*, satisfy*(12)-(5)

*. Assume that for a given finite time T*> 0

*both the solution u of*(1)

*and u*

_{t}belong to W^{2, }^{∞}(Ω)

*. Assume also that*∀

*t*∈ [0,

*T*] ,

**a**(·,

*t*) ∈ [

*W*

^{1, }

^{∞}(Ω)]

^{N}

*and f*(·,

*t*) ∈

*W*

^{1, }

^{∞}(Ω)

*, and*(

*u*)

_{t}

_{t}belongs to L^{∞}(Ω)

*. Finally let a strictly positive integer k*Δ

_{T}be defined as the minimum of all integers k such that the quantity*t : =*

*fulfills the condition*(11)

*. Let us further assume that, besides belonging to a quasiuniform family, all the partitions*ω

_{h}in use are of the acute type, and the quantity*is bounded below away from zero independently of h. Then there exists a constant C independent of u, h and*Δ

*t, such that the following estimate applies:*

** 4 Consistent choices of the weights**

In this section we describe two coherent strategies to determine a set of *M _{i}* strictly positive weights, for each mesh inner node

*P*, that are proven to be bounded below away from zero, independently of the mesh parameter

_{i}*h*.

Let us first consider the one-dimensional case. From (12) and (5) it is trivially seen that the pair of weights (, ) associated with inner node *P _{i}* is uniquely defined by the equations:

where and are the lengths of the intervals of _{h} having *P _{i}* as the right and the left end, respectively. This yields = and = .

Next we switch to the case *N* > 1. In principle for *N* = 2 or *N* = 3 there are infinitely many solutions, except for the case of the least possible value of *M _{i}*, i.e.

*M*=

_{i}*N*+ 1, in which the solution is necessarily unique. The two constructions described below allow for the unique determination of a set of weights satisfying (12) and (5), and incidentally they apply even to the particular case where this set is unique.

*4.1 A first set of weights*

Let := be the set of unit vectors corresponding to the vectors , that is, = /, where is the modulus of . Let also := , and *E _{i}* = ∪ . Setting

*K*=

_{i}*card*(

*E*), we clearly have

_{i}*K*

_{i}__>__

*M*. We number the

_{i}*M*vectors of

_{i}*E*∩ in the same manner as the vectors of .

_{i}Let us first consider the case where = = *E _{i}*. Then for every

*j*

_{1}∈ {1,...,

*M*} there is necessarily another

_{i}*j*

_{2}∈ {1,...,

*M*} such that + = 0. This implies in particular that

_{i}*M*must be an even number. Hence we may choose the weights in pairs, say (, ), in the same way as in the one-dimensional case (cf. (14)). More specifically, we number the vectors in

_{i}*E*in such a manner that the first

_{i}*M*/2 ones form a subset of

_{i}*E*whose vectors do not have any vector opposite to it in this subset, and from

_{i}*M*/2 + 1 up to

_{i}*M*the vectors in the complementary subset consisting of corresponding opposite vectors. In so doing the weights satisfy:

_{i}for all pairs (*j*_{1}, *j*_{2}) ∈ {1,..., *M _{i}*/2} ×{

*M*/2 + 1,...,

_{i}*M*}, such that + =

_{i}**0**. Notice that the set of weights determined by solving (15) trivially satisfy both (12) and (5).

Next we assume that ≠ (or equivalently, *E _{i}*≠ ). In this case we have

*K*>

_{i}*M*, and we number the

_{i}*K*-

_{i }*M*vectors in

_{i}*E*that do not belong to from

_{i}*M*+ 1 up to

_{i }*K*in an arbitrary order, say ,

_{i}*k*=

*M*+ 1,...,

_{i }*K*. By assumption the intersection of

_{i}_{i}with the half straight line with origin at

*P*and oriented in the sense and direction of for

_{i}*k*>

*M*, cannot be an edge

_{i}*P*with

_{i}*j*

__<__

*M*. Letting the segment

_{i}*P*be such intersection, where is necessarily a point of the boundary of

_{i}*S*, it follows that it is contained in either a single

_{i}*N*-simplex of

*S*for

_{i}*N*= 2 or

*N*= 3, or in a common face of exactly two neighboring tetrahedra for

*N*= 3. In any case, the vector leading from

*P*to , for

_{i}*k*=

*M*+ 1,...,

_{i}*K*, still denoted by , is a non trivial convex combination of either exactly

_{i}*N*edge vectors with 1

__<__

*j*

__<__

*M*pertaining to the same

_{i}*N*-simplex of

*S*, or of two such edge vectors pertaining to a common face of two neighboring tetrahedra of

_{i}*S*. Let

_{i}*J*be the number of such edge vectors, i.e., either

*J*=

*N*for

*N*= 2 or

*N*= 3, or

*J*= 2 for

*N*= 3 only, and , for

*l*= 1,...,

*J*with 1

__<__

*m*

_{l}__<__

*M*, be the vertices of

_{i}*S*whose convex combination yields for

_{i}*k*=

*M*+ 1,...,

_{i }*K*. Let us denote the modulus of by , for

_{i}*k*=

*M*+ 1,...,

_{i }*K*too. In so doing we have:

_{i}where the 's for *l* = 1,..., *J*, are coefficients of a non trivial convexcombination, that is, 0 < < 1, l = 1,..., *J* and *k* = *M _{i }*+ 1,...,

*K*, with

_{i}Now we momentarily assign to each node of *S _{i}* different from

*P*, and to each point ,

_{i}*j*=

*M*+ 1,...,

_{i }*K*, the same measure fraction, say

_{i}and the weight

Thanks to the fact that by construction, we necessarily have:

together with

as one can easily check.

Now we replace in (17) the 's for *k* = *M _{i }*+ 1,...,

*K*, by the expressiongiven by (16). Then rearranging the terms in the resulting expression, we establish that relation (12) holds for weights defined in the following manner:

_{i}In (19) *C _{i}* is a normalizing constant allowing (5) to hold. Notice that, provided the weight increments are all non-negative, the value of

*C*is strictly positive. The values of in turn, simply account for the sum of the contributions of for a given

_{i}*j*

__<__

*M*, to the vectors for

_{i}*k*>

*M*, expressed by (16), respectively multiplied by , and the corresponding convex combination coefficient. More specifically we have,

_{i}where = 0 if does not belong to *S _{i}* ∩ and = for the pertaining convex combination coefficients in (16), that is, all the 's such that

*m*=

_{l}*j*. Notice that by construction (5) holds for the weights defined by (19)-(20).

An important result that holds in connection with the above construction is

**Theorem 4.1.** *The weights generated in the way prescribed in this sub-section are all bounded below by a strictly positive constant **ω** independent of the mesh step size, whatever the mesh one might consider in the quasi-uniform family of meshes in use.*

**Proof. ** The case where *K _{i}* =

*M*is trivial and hence we give a detailed proof only for the case

_{i}*K*>

_{i}*M*.

_{i}First of all we note that

Next from (20) and (5) we have,

Since

and from (18) it holds that , taking into account (21) we obtain,

It follows from (22) that for a suitable mesh independent constant *c*_{0} we have *C _{i}*

__>__

*c*

_{0}∀

*i*. Indeed

*K*

_{i}__<__2

*M*and for no mesh of the quasi-uniform family of meshes {

_{i}_{h}}

_{h}under consideration,

*M*exceeds the value

_{i}*c*

^{- N}, where c is a mesh independent constant such that

*ρ*

__>__

*ch*for every

*h*,

*ρ*being the minimum over the elements of

_{h}of the radii of the largest inscribed balls in the elements of

_{h}for a given

*h*(cf. [3] ).

Finally, since > 0 for all *i*, *j* we have,

It immediately follows from (23) that there exists a suitable mesh independent strictly positive constant *c _{N}* such that

**Remark 4.2. ** (24) clearly indicates that *c _{N}* may play the role of the parameter

*ω*in the relations (9) and (11). However for very distorted meshes such a value of

*ω*may be largely under-evaluated, and for this reason in practical computations it is advisable to determine this parameter simply as the minimum of all the 's for a given mesh, keeping in mind that, according to (24) such a value is necessarily bounded below away from zero independently of the mesh size.

*4.2 A second set of weights*

Here we consider again the sets and . The case where = is to be treated in the same manner as above. Hence we may assume that ≠ . Distinguishing eventually coincident vectors in and , we put together all the vectors in both sets. Otherwise stated we are dealing with a set of exactly 2*M _{i}* vectors not necessarily distinct, say := , where the first

*M*vectors are those of and the last

_{i}*M*vectors are given by , for

_{i}*j*= 1,...,

*M*. Now let the segment

_{i}*P*be the intersection of

_{i}_{i}with the half straight line with origin at

*P*and oriented in the sense and direction of for k >

_{i}*M*, being necessarily a point of the boundary of

_{i}*S*. It follows that

_{i}*P*either coincides with an edge

_{i}*P*, where is a vertex of

_{i}*S*, or is contained in either a single

_{i}*N*-simplex of

*S*for

_{i}*N*= 2 or

*N*= 3, or in a common face of exactly two neighboring tetrahedra for

*N*= 3. In any case, the vector leading from

*P*to , for

_{i}*k*=

*M*+ 1,..., 2

_{i }*M*, still denoted by , is a convex combination of at most

_{i}*N*edge vectors with 1

__<__

*j*

__<__

*M*pertaining to the same

_{i}*N*-simplex of

*S*. Let

_{i}*J*be the number of such edge vectors, i.e.,

*J*= 1 if coincides with a given vertex of

*S*,

_{i}*J*= 2 for

*N*= 3 only if is contained in a common face of two neighboring tetrahedra of

*S*, or

_{i}*J*=

*N*for

*N*= 2 or

*N*= 3 otherwise. Let , for

*l*= 1,...,

*J*with 1

__<__

*m*

_{l}__<__

*M*, be the vertices of

_{i}*S*whose convex combination yields for

_{i}*k*=

*M*+ 1,..., 2

_{i }*M*. Here again we denote the modulus of by , for

_{i}*k*=

*M*+ 1,..., 2

_{i }*M*too. In so doing we have:

_{i}where the 's for *l* = 1,..., *J*, are coefficients of a convex combination, that is, 0 __<__ __<__ 1, *l* = 1,..., *J* and *k* = *M _{i }*+ 1,..., 2

*M*, with

_{i}Now we assign to each point , *k* = *M _{i }*+ 1,..., 2

*M*the measure fraction , and the weight = , where are provisional weights respectively associated with the vertices of

_{i}*S*different from

_{i}*P*satisfying = . Then similarly to the case of (15) we necessarily have:

_{i}together with

as one can easily check.

Now we replace in (26) the 's for *k* = *M _{i }*+ 1,..., 2

*M*, by the expression given by (25). Then rearranging the terms in the resulting expression, we establish that relation (12) holds for weights defined in the following manner:

_{i}In (28) *C _{i}* is a normalizing constant allowing (5) to hold. Notice that, provided the weight increments are all non-negative, the value of

*C*is strictly positive. The values of in turn, simply account for the sum of the contributions of for a given

_{i}*j*

__<__

*M*, to the vectors for

_{i}*k*>

*M*, expressed by (25), respectively multiplied by , and the corresponding convex combination coefficient. More specifically we have,

_{i}where = 0 if does not belong to *S _{i}* ∩ and = for the pertaining convex combination coefficients in (25), that is, all the 's such that

*m*=

_{l}*j*. Notice that by construction (5) holds for the weights defined by (28)-(29).

By arguments in all similar to those in the proof of Theorem 4.1 we can prove,

**Theirem 4.3.** *The weights generated in the way prescribed in this sub-section are all bounded below by a strictly positive constant **ω** independent of the mesh step size, whatever the mesh one might consider in the quasi-uniform family of meshes in use.*

**Remark 4.4. ** The pairs of expressions (19)-(20) and (28)-(29) defining two a priori different sets of weights attempt to take into account three main factors playing a role in the influence of the vector in equation (12). First of all its modulus, since the smaller it is, the larger must be the corresponding weight. Next the associated measure fraction , which goes in the same sense as . Finally the increments , which introduce the necessary adjustment in the weight values, in order to take into account the eventual existence of vectors in *L _{i}* roughly opposite to , that is, making angles with it close to

*π*. In this sense we may assert that the above described procedure for determining the 's, can be viewed as both close to optimal choices of the weights of our method, and straightforward manners to compute them.

**5 Implementation and Numerical Aspects**

It is not difficult to figure out from the construction of the sets of weights described in the last Section, that both choices give rise to a straightforward implementation of the method. Nevertheless we would like to point out that it is wiser to assemble node by node both the weighted mass matrix and the matrix generated by the discretisation of the convection-diffusion operator, instead of the usual element-by-element procedure. All that is needed for this purpose is a table of node numbers and an associated integer pointer vector having as many components as there are nodes. In the table the series of numbers of the neighbors of each inner node are successively stored, while in the pointer vector the *i - th* component contains the position in the table corresponding to the first neighbor of the *i - th* inner node, thereby allowing to uniquely identify all the series of neighbors in the table.

Next we check the performance of the method studied in this paper as compared to a well-established technique to deal with dominant convection in convection-diffusion problems, namely, the least squares formulation. Here the latter is implemented in connection with the Crank-Nicholson scheme for the time-integration, as described in [5] . This comparative study is illustrated by means of some results extracted from [8] , in the framework of a test-problem described below, where uniform meshes are used, thereby allowing the use of weights all equal to 1/(*N* + 2).

Take Ω to be the unit square (0, 1) ×(0, 1) and a space discretisation based on a uniform *L* ×*L* mesh in which every square cell is subdivided into two triangles by means of the diagonal parallel to the line *x*_{1} = *x*_{2}. We take *T* = 0.1 and *ν* = 10^{-k}, for *k* = 2 and *k* = 5, and **a** = (1/*π*; 1/*π*). For this choice the Péclet number equals 10* ^{k}*/

*π*. An exact solution is considered to be a function with a double boundary layer in the neighborhood of the edges given by

*x*

_{1}= 1 and

*x*

_{2}= 1. More specifically we take

where for *z* ∈ [0, 1] ,

Equation (1) has effectively the above solution, provided the right hand side is given by *f*(*x*_{1}, *x*_{2}, *t*) = *e ^{-t}* [

*s*(

*x*

_{1}) cos(

*π*

*x*

_{2}) +

*s*(

*x*

_{2}) cos(

*π*

*x*

_{1}) ] , and the prescribed initial and boundary values are those of

*u*.

In self explanatory Tables 1 through 4 we summarize the results generated by both methods in the solution of the above test problem, by showing both *L*^{2} and *L*^{∞} relative errors.

Observation of Tables 1 and 2 leads to the conclusion that both methods being compared simulate correctly thicker boundary layers, that is, those corresponding to a moderate Péclet number close to 10^{2}. Moreover, as one can infer from the above results, our scheme is much less accurate than the least-squares formulation in this case. This is quite natural for the latter is a second order method in the *L*^{2}-norm, and empirically in the *L*^{∞}-norm too, whereas the former is a quasi first order method in the L^{∞}-norm (cf. [7] ). Nevertheless the least-squares approach failed completely in the case of a thin boundary layer corresponding to Pé roughly equal to 10^{5}, as clearly indicated in Tables 3 and 4, while our explicit scheme showed a much better behavior. However even so the errors in the maximum norm are not so small. In this respect it is worthwhile commenting that such results are to be expected. Indeed the maximum error values tend to occur precisely in the interior of the thin boundary layer, which cannot be reached by any numerical method, at least not for the degree of mesh refinement used in the above test.

**Acknowledgement. ** The author is thankful to Carlos Tomei from the Department of Mathematics of PUC-Rio for fruitful discussions.

**References**

[1] R.A. Adams, *Sobolev Spaces*, Academic Press, N.Y., (1975). [ Links ]

[2] A.N. Brooks and T.J.R. Hughes, *The streamline upwind/Petrov-Galerkin formulation for convection dominated flows with particular emphasis on the Navier-Stokes equations*, Computer Methods in Applied Mechanics and Engineering, **32** (1982), 199-259. [ Links ]

[3] P.G. Ciarlet, *The Finite Element Method for Elliptic Problems*, Noth-Holland,Amsterdam (1978). [ Links ]

[4] M. Kawahara, N. Takeuchi and Y. Yoshida, *Two step finite element methods for Tsunami wave propagation analysis*, International Journal of Numerical Methods in Engineering, **12** (1978), 331-351. [ Links ]

[5] R. Leal Toledo and V. Ruas, *Numerical analysis of a least-squares finite element method for the time-dependent advection-diffusion equation*, Journal of Computational and Applied Mathematics, **235** (2011), 3615-3631. [ Links ]

[6] V. Ruas and A.P. Brasil Jr., *A stable explicit method for time-dependent convection-diffusion equations*, Proc. ICNAAM, Greece, (2007), 480-483. [ Links ]

[7] V. Ruas, A.P. Brasil Jr. and P. Trales, *An explicit scheme for solving convection-diffusion equations*, Japan Journal of Industrial and Applied Mathematics, **26** (2009), 65-91. [ Links ]

[8] V. Ruas, M. Kischinhevsky and R. Leal Toledo, *Elementos finitos em formulação mista de mínimos quadrados para a simulação da convecção-difusão em regime transiente*, to appear in Revista Internacional de Métodos Numéricos para Cálculo y Diseño en Ingeniería (January 2013). [ Links ]

[9] G. Strang and G.J. Fix, *An Analysis of the Finite Element Method*, Prentice-Hall, Englewood Cliffs (1973). [ Links ]

[10] M. Tabata, *Uniform convergence of the upwind finite element approximation for semilinear parabolic problems*, Journal of Mathematics of the Kyoto University, **18-2** (1978), 327-351. [ Links ]

Received: 22/IX/12.

Accepted: 28/IX/12.

#CAM-698/12.

*Supported by CNPq grant 307996/2008-5.