# «1 Introduction This paper is devoted to the numerical solution of the scalar convection– diﬀusion equation −ε ∆u + b · u=f in Ω, u = ub ...»

On discontinuity–capturing methods

for convection–diﬀusion equations

Volker John1 and Petr Knobloch2

Universit¨t des Saarlandes, Fachbereich 6.1 – Mathematik,

a

Postfach 15 11 50, 66041 Saarbr¨cken, Germany, john@math.uni-sb.de

u

Charles University, Faculty of Mathematics and Physics,

Sokolovsk´ 83, 186 75 Praha 8, Czech Republic, knobloch@karlin.mff.cuni.cz

a

Summary. This paper is devoted to the numerical solution of two–dimensional

steady scalar convection–diﬀusion equations using the ﬁnite element method. If the popular streamline upwind/Petrov–Galerkin (SUPG) method is used, spurious oscillations usually arise in the discrete solution along interior and boundary layers. We review various ﬁnite element discretizations designed to diminish these oscillations and we compare them computationally.

1 Introduction This paper is devoted to the numerical solution of the scalar convection– diﬀusion equation −ε ∆u + b · u=f in Ω, u = ub on ∂Ω, (1) where Ω ⊂ R2 is a bounded domain with a polygonal boundary ∂Ω, ε 0 is constant and b, f and ub are given functions.

If convection strongly dominates diﬀusion, the solution of (1) typically contains interior and boundary layers and solutions of Galerkin ﬁnite element discretizations are usually globally polluted by spurious oscillations. To enhance the stability and accuracy of these discretizations, various stabilization strategies have been developed during the past three decades. One of the most eﬃcient procedures is the SUPG method developed by Brooks and Hughes [2]. Unfortunately, the SUPG method does not preclude spurious oscillations localized in narrow regions along sharp layers and hence various terms introducing artiﬁcial crosswind diﬀusion in the neighbourhood of layers have been proposed to be added to the SUPG formulation. This procedure is often referred to as discontinuity capturing (or shock capturing). The literature on discontinuity–capturing methods is rather extended and numerical tests published in the literature do not allow to draw conclusions concerning their advantages and drawbacks. Therefore, the aim of this paper is to 2 Volker John and Petr Knobloch provide a review of various discontinuity–capturing methods and to compare these methods computationally.

The plan of the paper is as follows. In the next section, we recall the Galerkin discretization of (1) and, in Section 3, we formulate the SUPG method.

Section 4 contains a review and a computational comparison of discontinuity– capturing methods and, in Section 5, we present our conclusions.

2 Galerkin’s ﬁnite element discretization We introduce a triangulation Th of the domain Ω consisting of a ﬁnite number of open polygonal elements K. We assume

4 Methods diminishing spurious oscillations in layers In this section, we present a review and a computational comparison of most of the methods introduced during the last two decades to diminish the oscillations arising in discrete solutions of the problem (1). These methods can be divided into upwinding techniques and into methods adding additional artiﬁcial diﬀusion to the SUPG discretization (2). The artiﬁcial diﬀusion may be either isotropic, or orthogonal to streamlines, or based on edge stabilizations.

These four classes of methods will be discussed in the following subsections.

It is not possible to describe here thoroughly the ideas on which the design of the methods relies, see [11] for a more comprehesive description. Generally, one can say that the methods are based either on convergence analyses or on investigations of the discrete maximum principle (called DMP in the following) or on heuristic arguments. As we shall see, most of the methods will be nonlinear. The computational comparison of the methods will be performed

**by means of two test problems speciﬁed by the following data of (1):**

** Example 1. Ω = (0, 1)2, ε = 10−7, b = (cos(−π/3), sin(−π/3))T, f = 0, ub (x, y) = 0 for x = 1 or y ≤ 0.**

7, ub (x, y) = 1 otherwise.

** Example 2. Ω = (0, 1)2, ε = 10−7, b = (1, 0)T, f = 1, ub = 0.**

The solution of Ex. 1 possesses an interior layer and exponential boundary layers whereas the solution of Ex. 2 possesses parabolic and exponential boundary layers but no interior layers. All results were computed on uniform N × N triangulations of the type depicted in Fig. 1. Unless stated otherwise, we used the conforming linear ﬁnite element P1, N = 20 for Ex. 1 and N = 10 for Ex. 2. The SUPG solutions of Ex. 1 and 2 are shown in Fig. 2 and 5, respectively. It is important that the parameter τ is optimal for the P1 element in the sense that the SUPG method approximates the boundary layers at y = 0 in Ex. 1 and at x = 1 in Ex. 2 sharply and without oscillations.

**4.1 Upwinding techniques**

Initially, stabilizations of the Galerkin discretization of (1) imitated upwind ﬁnite diﬀerence techniques. However, like in the ﬁnite diﬀerence method, the upwind ﬁnite element discretizations remove the unwanted oscillations but the accuracy attained is often poor since too much numerical diﬀusion is introduced. According to our experiences, one of the most successful up

winding techniques is the improved Mizukami–Hughes (IMH) method, see Knobloch [14]. It is a nonlinear Petrov–Galerkin method for P1 elements which satisﬁes the DMP on weakly acute meshes. In contrast with many other upwinding methods for P1 elements satisfying the DMP, the IMH method adds much less numerical diﬀusion and provides rather accurate solutions, cf. Knobloch [15]. The IMH solution for Ex. 1 is depicted in Fig. 3. For Ex. 2, it is even nodally exact.

4.2 Methods adding isotropic artiﬁcial diﬀusion Hughes et al. [10] came with the idea to change the upwind direction in the SUPG term of (2) by adding a multiple of the function bh which is the projection of b into the direction of uh. This leads to the additional term

on the left–hand side of (2), where σ is a nonnegative stabilization parameter.

Since bh depends on uh, the resulting method is nonlinear. Hughes et al. [10] proposed to set σ = max{0, τ (bh ) − τ (b)} where we use the notation τ (b ) for τ deﬁned by (3) with b replaced by b. Other deﬁnitions of σ in (4) were proposed by Tezduyar and Park [17]. Since the term (4) equals to

which reduces the amount of artiﬁcial diﬀusion along the z h direction.

Do Carmo and Gale˜o [6] also introduced a feedback function which should a minimize the inﬂuence of the term (6) in regions where the solution u of (1) is smooth. Since this approach was rather involved, do Carmo and Alvarez [5] introduced another procedure (still deﬁned using several formulas) suppressing the addition of the artiﬁcial diﬀusion in regions where u is smooth.

Again, the term (6) can be written in the form (5), now with ε = σ |Rh (uh )|2 /| uh |2. To prove error estimates, Knopp et al. [16] proposed to replace this ε, on any K ∈ Th, by

4.3 Methods adding artiﬁcial diﬀusion orthogonally to streamlines Since the streamline diﬀusion introduced by the SUPG method seems to be enough along the streamlines, an alternative approach to the above methods is to modify the SUPG discretization (2) by adding artiﬁcial diﬀusion in the crosswind direction only as considered by Johnson et al. [13]. A straightforward generalization of their approach leads to the additional term

4.4 Edge stabilizations Another stabilization strategy for linear simplicial ﬁnite elements was introduced by Burman and Hansbo [4]. The term to be added to the left–hand side of (2) is deﬁned by

where t∂K is a unit tangent vector to the boundary ∂K of K, ΨK (uh ) = diam(K) (C1 ε + C2 diam(K)) maxE⊂∂K | [|nE · uh |]E |, nE are normal vectors to edges E of K, [|v|]E denotes the jump of a function v across the edge E and C1, C2 are appropriate constants. Burman and Hansbo proved that, using an edge stabilization instead of the SUPG term, the DMP is satisﬁed. Other choices of ΨK (uh ) based on investigations of the DMP were recently proposed by Burman and Ern. However, all these edge stabilizations add much more artiﬁcial diﬀusion than the best methods of the previous subsections.

5 Conclusions Our computations indicate that, among the methods mentioned in this paper, the best ones are: the IMH method [14], the method of do Carmo, Gale˜o [6] a deﬁned by (6), (7), the method of Almeida and Silva [1] deﬁned by (6), (8), the MC method introduced below (12) and the MBE method deﬁned by (10), (13). The IMH method can be used for the P1 element only but gives best results in this case. The other methods can be successfully also applied to other ﬁnite elements as Figs. 7 and 8 show (for the conforming quadratic nc element P2 and the nonconforming Crouzeix–Raviart element P1 ). However, much more comprehensive numerical studies are still necessary to obtain clear conclusions of the advantages and drawbacks of the discontinuity–capturing methods.

Acknowledgements The work of Petr Knobloch is a part of the research project MSM 0021620839 ﬁnanced by MSMT and it was partly supported by the Grant Agency of the Charles University in Prague under the grant No. 343/2005/B–MAT/MFF.

8 Volker John and Petr Knobloch References

1. Almeida, R.C., Silva, R.S.: A stable Petrov–Galerkin method for convection– dominated problems. Comput. Methods Appl. Mech. Eng., 140, 291–304 (1997)

2. Brooks, A.N., Hughes, T.J.R.: Streamline upwind/Petrov–Galerkin formulations for convection dominated ﬂows with particular emphasis on the incompressible Navier–Stokes equations. Comput. Methods Appl. Mech. Eng., 32, 199–259 (1982)

3. Burman, E., Ern, A.: Nonlinear diﬀusion and discrete maximum principle for stabilized Galerkin approximations of the convection–diﬀusion–reaction equation.

Comput. Methods Appl. Mech. Eng., 191, 3833–3855 (2002)

4. Burman, E., Hansbo, P.: Edge stabilization for Galerkin approximations of convection–diﬀusion–reaction problems. Comput. Methods Appl. Mech. Eng., 193, 1437–1453 (2004)

5. do Carmo, E.G.D., Alvarez, G.B.: A new stabilized ﬁnite element formulation for scalar convection–diﬀusion problems: The streamline and approximate upwind/Petrov–Galerkin method. Comput. Methods Appl. Mech. Eng., 192, 3379–3396 (2003)

6. do Carmo, E.G.D., Gale˜o, A.C.: Feedback Petrov–Galerkin methods for conveca tion–dominated problems. Comput. Methods Appl. Mech. Eng., 88, 1–16 (1991)

7. Codina, R.: A discontinuity–capturing crosswind–dissipation for the ﬁnite element solution of the convection–diﬀusion equation. Comput. Methods Appl.

Mech. Eng., 110, 325–342 (1993)

8. Gale˜o, A.C., Almeida, R.C., Malta, S.M.C., Loula, A.F.D.: Finite element anala ysis of convection dominated reaction–diﬀusion problems. Appl. Numer. Math., 48, 205–222 (2004)

9. Gale˜o, A.C., do Carmo, E.G.D.: A consistent approximate upwind Petrov– a Galerkin method for convection–dominated problems. Comput. Methods Appl.

Mech. Eng., 68, 83–95 (1988)

10. Hughes, T.J.R., Mallet, M., Mizukami, A.: A new ﬁnite element formulation for computational ﬂuid dynamics: II. Beyond SUPG. Comput. Methods Appl.

Mech. Eng., 54, 341–355 (1986)

11. John, V., Knobloch, P.: A comparison of spurious oscillations at layers diminishing (SOLD) methods for convection–diﬀusion equations: Part I. Preprint Nr. 156, FR 6.1 – Mathematik, Universit¨t des Saarlandes, Saarbr¨cken (2005) a u

12. Johnson, C.: Adaptive ﬁnite element methods for diﬀusion and convection problems. Comput. Methods Appl. Mech. Eng., 82, 301–322 (1990)

13. Johnson, C., Schatz, A.H., Wahlbin, L.B.: Crosswind smear and pointwise errors in streamline diﬀusion ﬁnite element methods. Math. Comput., 49, 25–38 (1987)

14. Knobloch, P.: Improvements of the Mizukami–Hughes method for convection– diﬀusion equations. Preprint No. MATH–knm–2005/6, Faculty of Mathematics and Physics, Charles University, Prague (2005)

15. Knobloch, P.: Numerical solution of convection–diﬀusion equations using upwinding techniques satisfying the discrete maximum principle. Submitted to the Proceedings of the Czech–Japanese Seminar in Applied Mathematics 2005

16. Knopp, T., Lube, G., Rapin, G.: Stabilized ﬁnite element methods with shock capturing for advection–diﬀusion problems. Comput. Methods Appl. Mech. Eng., 191, 2997–3013 (2002)

17. Tezduyar, T.E., Park, Y.J.: Discontinuity–capturing ﬁnite element formulations for nonlinear convection–diﬀusion–reaction equations. Comput. Methods Appl.

Mech. Eng., 59, 307–325 (1986)