Robust Equilibrated Residual Error Estimator for ... - Semantic Scholar

Report 2 Downloads 54 Views
Robust Equilibrated Residual Error Estimator for Diffusion Problems: Conforming Elements∗ Zhiqiang Cai†

Shun Zhang‡

November 11, 2011 Abstract. This paper analyzes an equilibrated residual a posteriori error estimator for the diffusion problem. The estimator, which is a modification of those in [9, 30], is based on the Prager-Synge identity and on a local recovery of an equilibrated flux. Numerical results for an interface test problem show that the modification is necessary for the robustness of the estimator. When the distribution of diffusion coefficients is local quasi-monotone, it is shown theoretically that the estimator is robust with respect to the size of jumps.

1

Introduction

Let Ω be a bounded polygonal/polyhedronal domain in 0 independent of α and h: sup τ

k−1 ∈RT−1,z

bz (τ , (v, µ)) ≥ β|||(v, µ)|||α,h,z kτ kα−1 ,h,z

∀ (v, µ) ∈ Pzk−1 × Mzk−1 .

(5.9)

k−1 Proof. Choose a τ˜ ∈ RT−1,z such that

Z K

(˜ τ + α∇v) · ∇qdx = 0 ∀ q ∈ Pk−2 (K)

∀ K ∈ Tz

and that

τ˜ · nK |F =

 α ,h α  sgn(K, F ) F [[v]] + K (µ + {v}w )   h 2 hF  F

α

K (µ + v)   2 h  F 

F ∈ EK ∩ EI ,z , F ∈ EK ∩ EN ,z , if z ∈ ND ∪ NN ,

(5.10)

F ∈ EK ∩ E0,z ,

0

k−1 where sgn(K, F ) = nK · nF . (See [12] for existence of τ˜ ∈ RT−1,z satisfying the above conditions.) Obviously, (5.10) implies

[[˜ τ · nF ]] =

αF ,a (µ + {v}w ) hF

∀ F ∈ Ej,z

and {˜ τ · nF } w =

12

αF ,h [[v]] hF

∀ F ∈ EI ,z ,

which, together with (5.3), gives b(˜ τ , (v, µ)) = |||(v, µ)|||2α,h,z .

(5.11)

For every K ∈ Tz , by the standard scaling argument, there exists a constant C > 0 independent of α and the mesh size such that 



X

k˜ τ k20,K ≤ C kαK ∇vk20,K + hK

F ∈EK ∩EI ,z

α ,h α ksgn(K, F ) F [[v]] + K (µ + {v}w )k20,F + αK Λ(v, µ, N ) , hF 2 hF

which, together with (5.2), gives  −1/2

kαK

τ˜ k20,K

≤ C kα

1/2

∇vk20,K

X

+

F ∈EK ∩EI ,z



αF ,h α k[[v]]k20,F + K kµ + {v}w k20,F hF hF





+ Λ(v, µ, N ) .

Hence, there exists a constant C˜ > 0 independent of α and h such that ˜ kτ kα−1 ,h,z ≤ C|||(v, µ)|||α,h,z . ˜ This completes the proof of the which, together with (5.11), leads to (5.9) with β = 1/C. lemma. k−1 Theorem 5.1. Problem (3.2) has a unique solution (σ ∆ , wz , λz ) ∈ RT−1,z × Pzk−1 × Mzk−1 that Tz satisfies the following a priori estimate:

kσ ∆ k −1 ,h,z + |||(wz , λz )|||α,h,z ≤ C(ac , β) Tz α

sup (v,µ)∈Pzk−1 ×Mzk−1

Rz (v) + Jz (µ) , |||(v, µ)|||α,h,z

(5.12)

where the constant C(ac , β) > 0 is independent of the mesh size and jumps. Proof. The theorem follows from the abstract theory of saddle point problem (see, e.g., [12]) and Lemmas 5.4, 5.5, and 5.6.

5.2

Local Efficiency Bound

For any z ∈ N , let ω ˆ z = {K ∈ ωz : αK = max αK 0 }. 0 K

∈ωz

Assume that the distribution of the coefficients αK for all K ∈ T is locally quasi-monotone [24] which is slightly weaker than Hypothesis 2.7 in [7]. For convenience of readers, we restate it here. Definition 5.7. Given a vertex z ∈ N , the distribution of the coefficients αK , K ∈ ωz , is said to be quasi-monotone with respect to the vertex z if there exists a subset ω ˜ K ,z,qm of ωz such that the union of elements in ω ˜ K ,z,qm is a Lipschitz domain and that • if z ∈ N \ND , then {K} ∪ ω ˆz ⊂ ω ˜ K ,z,qm and αK ≤ αK 0 ∀K 0 ∈ ω ˜ K ,z,qm ; • if z ∈ ND , then K ∈ ω ˜ K ,z,qm , ∂ ω ˜ K ,z,qm ∩ ΓD 6= ∅, and αK ≤ αK 0 ∀K 0 ∈ ω ˜ K ,z,qm . The distribution of the coefficients αK , K ∈ T , is said to be locally quasi-monotone if it is quasimonotone with respect to every vertex z ∈ N . 13

For an element K ∈ Tz and any v ∈ H 1 (K), let v¯K =

1 |K|

Z K

v dx and v¯K ,F =

1 |F |

Z F

v|K dx

∀ F ∈ EK

be the averages of v|K over K and face F ∈ EK , respectively. It is well known (see, e.g, (4.7) of [31]) that the following Poincaré-Fridrichs inequalities hold on element K with diameter hK ; i.e., there exists a positive constant C independent of hK such that kv − v¯K k0,K ≤ ChK k∇vk0,K

and kv − v¯K ,F k0,K ≤ ChK k∇vk0,K

∀ v ∈ H 1 (K).

(5.13)

Define the following piecewise H 1 function spaces Vz =

 2 1  {v ∈ L (ωz ) : v|K ∈ H (K)} 

z ∈ N \ND ,

{v ∈ L2 (ωz ) : v|K ∈ H 1 (K) and v|ΓD = 0}

z ∈ ND .

ˆ z , define v¯z = v¯K 0 . Obviously, Pzk ⊂ Vz . Let K 0 ∈ ω Lemma 5.8. Suppose the distribution of α is quasi-monotone with respect to the vertex z ∈ N . For any v ∈ Vz , there exists a constant C independent of the mesh size and α such that  X

1/2 (v − v¯z )k20,K ≤ C  h−2 K kα

 2 X α1/2 Z F [[v]]ds  kα1/2 ∇vk20,K + hd

(5.14)

 2 X α1/2 Z F [[v]]ds  kα1/2 ∇vk20,K + hd

(5.15)

X

F ∈EI ,z

K∈Tz

K∈Tz

F

F

when z ∈ N \ND and that  1/2 2 vk0,K ≤ C  h−2 K kα

X K∈Tz

X K∈Tz

F ∈EI ,z

F

F

when z ∈ ND . Proof. We will only show the validity of (5.14) since (5.15) may be proved in a similar fashion. By the definition of v¯z and (5.13), there exists a K 0 ∈ ω ˆ z such that v¯z = v¯K 0 and that kα1/2 (v − v¯z )k0,K 0 = kα1/2 (v − v¯K 0 )k0,K 0 ≤ ChK 0 kα1/2 ∇vk0,K 0 . For any K ∈ Tz , the triangle inequality gives 1/2 1/2 1/2 kαK (v − v¯z )k0,K ≤ kαK (v − v¯K )k0,K + kαK (¯ vK 0 − v¯K )k0,K .

Since α is quasi-monotone with respect to z, there is a connected path, {K = K0 , K1 , · · · , Kl = K 0 } with Ki ∈ Tz , from K to K 0 such that αK ≤ αK1 ≤ · · · ≤ αK 0 . Denote by Fi the common face between Ki−1 and Ki , then αK ≤ αF1 ≤ · · · ≤ αFl .

14

Using the triangle inequality, we have 1/2 1/2 1/2 kαK (¯ vK 0 − v¯K )k0,K = αK |K|1/2 |¯ vK0 − v¯Kl | ≤ αK |K|1/2

l−1 X

|¯ vKi − v¯Ki+1 |

i=0 1/2 ≤ αK |K|1/2

l−1  X



|¯ vKi − v¯Ki ,Fi+1 | + |¯ vKi ,Fi+1 − v¯Ki+1 ,Fi+1 | + |¯ vKi+1 ,Fi+1 − v¯Ki+1 | .

i=0

Since Tz is quasi-uniform, it then follows from the fact that αK ≤ αKi , (5.13), and the triangle inequality that 1/2 1/2 αK |K|1/2 |¯ vKi − v¯Ki ,Fi+1 | ≤ CkαK (¯ vKi − v¯Ki ,Fi+1 )k0,Ki i





1/2 1/2 1/2 ≤ C kαK (v − v¯Ki )k0,Ki + kαK (v − v¯Ki ,Fi+1 )k0,Ki ≤ ChKi kαK ∇vk0,Ki . i

i

i

Similarly, we have 1/2 1/2 αK |K|1/2 |¯ vKi+1 ,Fi+1 − v¯Ki+1 | ≤ ChKi+1 kαK ∇vk0,Ki+1 , i+1

1/2

αK |K|

1/2

1/2

|¯ vKi ,Fi+1 − v¯Ki+1 ,Fi+1 | ≤ αF

i+1

Z [[v]]ds . Fi+1

1−d/2

hF

i+1

Combining the above inequalities gives 1/2 h−1 K kαK (v

l X

− v¯z )k0,K ≤ C

1/2

kαK ∇vk0,Ki + i

i=0

l X

1/2

αF

i=1

i+1

! Z [[v]]ds . Fi+1

−d/2

hF

i+1

Squaring both sides, using the Cauchy-Schwarz inequality, and summing up over all K ∈ Tz imply (5.14). This completes the proof of the lemma. Corollary 5.9. Under the assumptions of Lemma 5.8, for any v ∈ Vz , there exists a constant C independent of the mesh size and α such that  X

h−2 kα1/2 (v − v¯z )k20,K ≤ C  K

K∈Tz

X

kα1/2 ∇vk20,K +

X α1/2 F F ∈EI ,z

K∈Tz

hF



k[[v]]k20,F 

(5.16)

when z ∈ N \ND and  X

h−2 kα1/2 vk20,K ≤ C  K

K∈Tz

X

kα1/2 ∇vk20,K +

X α1/2 F F ∈EI ,z

K∈Tz

hF



k[[v]]k20,F 

when z ∈ ND . Proof. The corollary is an immediate consequence of Lemma 5.8 and the fact that Z

d−1

[[v]]ds ≤ hF 2 k[[v]]k0,F .

F

15

(5.17)

Denote the vertex based local residual error indicator by 

ηˆz = 

X h2 K

αK

K∈Tz

1/2

X

krK k20,K +

F ∈Ej,z

hF kj k2  αF ,a F 0,F

,

(5.18)

its local efficiency was proved in [7, 24]; i.e., there exists a constant C > 0 is independent of α and the mesh size such that ηˆz ≤ Ckα1/2 ∇(u − uT )k0,˜ωz , (5.19) where ω ˜ z = ωz ∪ {K : K and ∂ωz shares an edge/face }. Note that there are no oscillation terms in (5.19) due to the assumption that f is a piecewise polynomial of degree k − 1. It is well known (see, e.g., [29]) that the residual functional has the following L2 -representation: X

R(v) = f (v) − (A∇uT , ∇v) =

X

(rK , v)K −

(jF , v)F

1 ∀ v ∈ HD (Ω).

F ∈EI ∪EN

K∈T

Define the local residual on patch ωz by X

Rz (v) := R(vφz ) =

X

(rK ,z , v)K −

(jF ,z , v)F

1 ∀ v ∈ HD (Ω)

F ∈Ej,z

K∈Tz

with rK ,z = φz rK and jF ,z = φz jF . By the Galerkin orthorgnality, i.e., Rz (1) = R(φz ) = 0, we then have the following local orthogonality property X

(rK ,z , 1)K −

X

(jF ,z , 1)F = 0.

(5.20)

F ∈Ej,z

K∈Tz

Theorem 5.10. (Efficiency) Under the assumptions of Lemma 5.8, the local indicators ηz and ηK are efficient; i.e., there exists a constant C > 0 independent of α and the mesh size such that ηz ≤ Ckα1/2 ∇ek0,˜ωz

and

ηK ≤ C

X

kα1/2 ∇ek0,˜ωz .

(5.21)

z∈NK

Proof. Squaring both sides of the first inequality in (5.21) and summing up over all z ∈ NK imply the second inequality in (5.21). To prove the validity of the first inequality in (5.21), by Theorem 5.1 and (5.19), it suffices to show that sup (v,µ)∈Pzk−1 ×Mzk−1

Rz (v) + Jz (µ) ≤ C ηˆz |||(v, µ)|||α,h,z

or, equivalently, I ≡ Rz (v) + Jz (µ) ≤ C ηˆz |||(v, µ)|||α,h,z

∀ (v, µ) ∈ Pzk−1 × Mzk−1 .

To do so, from the local orthogonality property in (5.20) and the facts (rK ,z , 1)K = (¯ rK ,z , 1)K

∀ K ∈ Tz

and (jF ,z , 1)F = (¯jF ,z , 1)F

∀ F ∈ Ej,z ,

we have that for an arbitrary constant c Rz (c) − Jz (c) =

X

(¯ rK ,z , c)K −

X F ∈EI ,z ∪EN ,z

K∈Tz

16

(¯jF ,z , c)F = 0,

(5.22)

which implies I = Rz (v − c) + Jz (µ + c) = Rz (v − c) − Jz ({v − c}w ) + Jz (µ + {v}w )

(5.23)

for any (v, µ) ∈ Pzk−1 × Mzk−1 . It follow from the triangle inequality, the facts that k¯ rK ,z k0,K ≤ krK k0,K and k¯jF ,z k0,F ≤ kjF k0,F , and the Cauchy-Schwarz inequality that I ≤

X

X

k¯ rK ,z k0,K kv − v¯z k0,K +

F ∈Ej,z

K∈Tz



≤ ηˆz 

k¯jF ,z k0,F (k{v − v¯z }w k0,F + kµ + {v}w k0,F )

X α K

K∈Tz



≤ C ηˆz 

h2K

kv −

K∈Tz

X α ,a  F

+

F ∈Ej,z

X α K

h2K

v¯z k20,K

kv − v¯z k20,K +

hF

X α ,a F F ∈Ej,z



X α ,h F

≤ C ηˆz kα1/2 ∇h vk20,ωz +

F ∈EI ,z

hF

hF

1/2

k{v −

v¯z }w k20,F

+ kµ +

{v}w k20,F

 

1/2

kµ + {v}w k20,F 

X α ,a F

k[[v]]k20,F +

F ∈Ej,z

hF

1/2

kµ + {v}w k20,F 

= C ηˆz |||(v, µ)|||α,h,z . This proves the validity of (5.22) and, hence, the theorem.

6

Numerical Experiments

In this section, we report some numerical results for an interface problem with intersecting interfaces used by many authors, e.g., [19, 14, 15, 16], which is considered as a benchmark test problem. For this test problem, we show numerically that the local minimization procedure is essential. That is, estimators based on the σ ∆ as in [8, 9] will produce a non-optimal mesh and, Tz ,e hence, are not robust. To this end, let Ω = (−1, 1)2 and u(r, θ) = rγ µ(θ) in the polar coordinates at the origin with µ(θ) being a smooth function of θ (see, e.g., [14]). The function u(r, θ) satisfies the interface equation with A = αI, ΓN = ∅, f = 0, and (

α(x) =

R

in (0, 1)2 ∪ (−1, 0)2 ,

1

in Ω \ ([0, 1]2 ∪ [−1, 0]2 ).

The γ depends on the size of the jump. In our test problem, γ = 0.1 is chosen and is corresponding to R ≈ 161.4476387975881. Note that the solution u(r, θ) is only in H 1+γ− (Ω) for any  > 0 and, hence, it is very singular for small γ at the origin. This suggests that refinement is centered around the origin. Let uT ∈ S 1 be the linear finite element approximation. We start with the coarsest triangulation T0 obtained from halving 16 congruent squares by connecting the bottom left and upper right corners. A sequence of meshes is generated by using standard adaptive meshing algorithm that 17

adopts the maximum marking strategy: mark those elements such that ηK ≥ 0.5 max η0 . Marked K 0 ∈T K η triangles are refined by bisection. Define the effectivity index: eff-index := , kα1/2 ∇(u − uT )k0,Ω kα1/2 ∇(u − uT )k0,Ω ≤ tol. We report numerical and use the following stopping criteria: rel-err := kα1/2 ∇uk0,Ω results with tol = 0.05. Broken RT0 elements are used to approximate the error flux. Denote by ξ the error estimator based on the σ ∆ similar to that in [8, 9] and by η = ηK the error estimator defined in (4.2), i.e., Tz ,e ∆ ∆ based on σ Tz = σ Tz ,e + σ ∆ . Here, the correction σ ∆ is given by (see Remark 3.2) Tz ,d Tz ,d 

σ∆ =− Tz ,d

A−1 σ ∆ , ∇⊥ φz Tz ,e



(A−1 ∇⊥ φz , ∇⊥ φz )

∇⊥ φz ,

where φz is the nodal basis function associated with the vertex z ∈ N . Figures 1 and 2 clearly show that adaptive mesh refinement using ξ as the indicator introduces unnecessary refinements along the interfaces and that the convergence rate is not optimal. 1

||A

1/2

(¢ u ï¢ uh) ||0 and j

10

||A1/2(¢ u ï¢ uh) ||0 0

10

j

ï1

10

ï2

10

1

10

Figure 1: mesh generated by ξ

2

10

3

10 number of nodes

4

10

5

10

Figure 2: error and extimator ξ

Mesh generated by η is shown in Figure 3. The refinement is centered at origin, The mesh is optimal with no over refinement along the interface. Similar meshes for this test problem generated by other error estimators can be found in [14, 15]. The comparison of the error and the η is shown in Figure 4. The effectivity index is close to 1.4. Moreover, the slope of the log(dof)log(relative error) for η is −1/2, which indicates the optimal decay of the error with respect to the number of unknowns.

References [1] M. Ainsworth, Robust a posteriori error estimation for nonconforming finite element approximation, SIAM J. Numer. Anal., 42:6 (2005), 2320-2341. 1 [2] M. Ainsworth, A posteriori error estimation for lowest order Raviart-Thomas mixed finite elements, SIAM J. Sci. Comput., 30 (2007), 189-204. 1 18

1

10

||A

1/2

(¢ u ï¢ uh) ||0 and d

||A1/2(¢ u ï¢ u ) || h

0

10

0

d

ï1

10

ï2

10

1

10

Figure 3: mesh generated by η

2

10

3

10 number of nodes

4

10

5

10

Figure 4: error and extimator η

[3] M. Ainsworth and J. T. Oden, A unified approach to a posteriori error estimation using element residual methods, Numer. Math., 65 (1993), 23–50. 1 [4] M. Ainsworth and J. T. Oden, A Posteriori Error Estimation in Finite Element Analysis, John Wiley & Sons, Inc., 2000. 1 [5] R. Bank and A. Weiser, Some a posteriori error estimators for elliptic partial differential equations, Math. Comp., 44 (1985), 283–301. 1 [6] I. Babuška, J. Osborn, and J. Pitkäranta, Analysis of mixed methods using mesh dependent norms, Math. Comp., 35 (1980), pp. 1039–1062. 5.1 [7] C. Bernardi and R. Verfürth, Adaptive finite element methods for elliptic equations with non-smooth coefficients, Numer. Math., 85:4 (2000), 579-608. 1, 5.2, 5.2 [8] D. Braess, Finite Elements: Theory, Fast Solvers and Applications in Solid Mechanics, 3rd edition, Cambridge University Press, Cambridge, UK, 2007. 1, 3, 6 [9] D. Braess and J. Schöberl, Equilibrated residual error estimator for edge elements, Math. Comp., 77 (2008), 651–672. (document), 1, 1, 6 [10] D. Braess, V. Pillwein, and J. Schöberl, Equilibrated residual error estimates are p-robust, Comput. Methods Appl. Mech. Engrg. 198 (2009) 1189–1197. 1, 3, 5 [11] D. Braess and R. Verfürth, A posteriori error estimators for the Raviart-Thomas element, SIAM J. Numer. Anal., 33 (1996), 2431–2444. 5.1 [12] F. Brezzi and M. Fortin, Mixed and Hybrid Finite Element Methods, Springer-Verlag, New York, 1991. 1, 3, 3, 3, 5, 5.1, 5.1 [13] Z. Cai, X. Ye, and S. Zhang, Discontinuous Galerkin finite element methods for interface problems: A priori and a posteriori error estimations, SIAM J. Numer. Anal., 49:5 (2011), 1761-1787. 1

19

[14] Z. Cai and S. Zhang, Recovery-based error estimator for interface problems: Conforming linear elements, SIAM J. Numer. Anal., 47:3 (2009), 2132–2156. 1, 6, 6 [15] Z. Cai and S. Zhang, Recovery-based error estimator for interface problems: Mixed and nonconforming elements, SIAM J. Numer. Anal., 48:1 (2010), 30–52. 1, 6, 6 [16] Z. Cai and S. Zhang, Flux recovery and a posteriori error estimators: Conforming elements for scalar elliptic equations, SIAM J. Numer. Anal., 48:2 (2010), 578–602. 1, 6 [17] P. Destuynder and B. Métivet, Explicit error bounds for a nonconforming finite element method, SIAM J. Numer. Anal., 35:5 (1998), 2099–2115. 1 [18] P. Destuynder and B. Métivet, Explicit error bounds in a conforming finite element method, Math. Comp. 68 (1999), 1379-1396. 1 [19] K-Y. Kim, A posteriori error analysis for locally conservative mixed methods, Math. Comp., 76 (2007), 43-66. 1, 6 [20] P. Ladevèze and D. Leguillon, Error estimate procedure in the finite element method and applications, SIAM J. Numer. Anal., 20:3 (1983), 485–509.. 1 [21] R. Luce and B. I. Wohlmuth, A local a posteriori error estimator based on equilibrated fluxes, SIAM J. Numer. Anal., 42:4 (2004), 1394–1414. 1 [22] J. C. Nédélec, Mixed finite elements in