• 検索結果がありません。

Formulation of the problem

N/A
N/A
Protected

Academic year: 2022

シェア "Formulation of the problem"

Copied!
10
0
0

読み込み中.... (全文を見る)

全文

(1)

ON A POSTERIORI ERROR ESTIMATORS IN THE FINITE ELEMENT METHOD

ON ANISOTROPIC MESHES

MANFRED DOBROWOLSKIy, STEFFEN GR ¨AFz,ANDCHRISTOPH PFLAUMx

Abstract. On anisotropic finite element meshes, the standard residual based error indicator is derived and it is proved that it is not efficient if the aspect ratio deteriorates. For a nonlocal error indicator it is proved that it is reliable and efficient independent of the aspect ratio. This is also confirmed by some numerical calculations.

Key words. finite elements, a posteriori estimators, anisotropic meshes.

AMS subject classifications. 65N30, 65N15.

1. Formulation of the problem. Let

IR

2be a bounded polygonal domain. Con- sider Poisson’s equation

,

u

=

f

in

; u

=0 on

@

:

(1.1)

For approximating problem (1.1) we use the standard conforming finite element method on an anisotropic triangulation=fgwhich is defined by the following conditions:

a) The intersection of two triangles is void or coincides with a common side or vertex.

b) The interior angles of the triangles are bounded from above by

< :

c) Let

U

()denote the union of the triangles adjacent to

:

It is assumed that each

U

()can be rotated such that it can be represented as the im- age of an isotropic reference configuration

U

^()^ of size O(1) under the mapping

x

i=

h

i

x

^i.

(1.2)

The last condition ensures that the direction of the anisotropic mesh does not change too rapidly. Since condition (1.2b) guarantees that in the anisotropic case two of the sides ofare long and nearly perpendicular to the small side, there exists a local orthogonal coordinate sys- tem(

e

1

;e

2)=(

e

1()

;e

2())where

e

1can be chosen to be the direction of one of the larger sides. Similarly, the long and short local step sizes are denoted by(

h

1

;h

2)=(

h

1()

;h

2())

:

The sets of long and small sides are denoted by,land,s

;

respectively.

For

k

=1

;

2we define the standard finite element spaces consisting of continuous and piecewise linear or quadratic shape functions,

S

k=

v

2

C

():

v

j 2

IP

kfor all2and

v

j@=0

:

Denoting by

I

k :

C

0()\

H

01;2() !

S

k the standard interpolation operators we obtain from Theorem 2 in [1] the estimates

jj

D

1(

u

,

I

k

u

)jj2;

ch

,11 X

jj=k +1

h

jj

D

u

jj2;

;

(1.3)

Received October 19, 1998. Accepted December 22, 1998. Recommended for publication by M. Eiermann.

yDepartment of Applied Mathematics and Statistics, University of W ¨urzburg, Am Hubland, D-97074 W ¨urzburg, Germany,([email protected])

zDepartment of Applied Mathematics and Statistics, University of W ¨urzburg, Am Hubland, D-97074 W ¨urzburg, Germany,([email protected])

xDepartment of Applied Mathematics and Statistics, University of W ¨urzburg, Am Hubland, D-97074 W ¨urzburg, Germany,([email protected])

36

(2)

jj

D

2(

u

,

I

k

u

)jj2;

c

X

jj=k

h

jj

D

Du

jj2;

;

(1.4)

where

D

i =@ei@ and

h

=

h

11

h

22

:

The finite element approximation

P

k

u

2

S

kis defined by

(

DP

k

u;Dv

)=(

f;v

) 8

v

2

S

k

:

(1.5)

We are interested in a posteriori estimates for the error

e

=

u

,

P

1

u

by local error indicators

satisfying

m

1jj

De

jj2,

T

(

h;f

)X

2

m

2jj

De

jj2+

T

(

h;f

)

;

(1.6)

where

T

(

h;f

)is usually a small term depending on

f

and converging to0for

h

!0

:

For earlier work on a posteriori error estimators on isotropic meshes we refer to Babu˘ska and Rheinboldt [2] and to the survey Verf¨urth [10]. Of special interest are the residual based indicator of Verf¨urth [9] and the indicators of Bank and Weiser [4]. The crucial point of anisotropic a posteriori estimating is the fact that all classical estimators deteriorate if the aspect ratio

a

() =

h

1()

=h

2()tends to infinity. Siebert [8] solves this problem by lo- cally balancing the directional errors avoiding anisotropic overrefinement. On the other hand, overrefinement occurs in elliptic systems where one equation is singularly perturbed and the others are not.

The outline of the paper is as follows. In section 2, we show that the standard error estimator based on the residual does not satisfy (1.6) with constants

m

1

;m

2independent of the aspect ratio

a:

Section 3 is devoted to the study of a nonlocal error estimator inspired by the third indicator of Bank and Weiser [4]. Despite the fact that the estimator is nonlo- cal, it is proved that it can be computed economically on isotropic and anisotropic meshes.

On anisotropic meshes, the estimator shows a significant propagation of local errors along the small mesh direction

e

2

;

which clearly indicates that local a posteriori error estimation is impossible as long as the standard normjj

D

jjis used. Some numerical computations demonstrate that the nonlocal indicator behaves exactly as predicted by the theory.

2. An error estimator based on local residuals. By

R

1 :

H

01;2 !

S

1 we denote

the local approximation operator constructed by Scott and Zhang [7] which satisfies, on an isotropic mesh with mesh parameter

h

=1,

jj

v

,

R

1

v

jj2^

c

jj

Dv

jj2U(

^

)

;

(2.1)

jj

v

,

R

1

v

jj2,^

c

jj

Dv

jj2U(

^

,)

;

(2.2)

where

U

(,)^ consists of the union of the triangles adjacent to,^

:

In view of condition (1.2) c) we can transform (2.1) , (2.2) by

x

i=

h

i

x

^iand obtain

jj

v

,

R

1

v

jj2

c

n

h

21jj

D

1

v

jj2U()+

h

22jj

D

2

v

jj2U()o

;

(2.3)

jj

v

,

R

1

v

jj2,

ch

,12 n

h

21jj

D

1

v

jj2U(,)+

h

22jj

D

2

v

jj2U(,)o

;

,2,l

;

(2.4)

jj

v

,

R

1

v

jj2,

ch

,11 n

h

21jj

D

1

v

jj2U(,)+

h

22jj

D

2

v

jj2U(,)o

;

,2,s

:

(2.5)

(3)

38 A posteriori estimators

Using the orthogonality relation(

De;Dv

) = 0for all

v

2

S

1 and integration by parts, it follows that

jj

De

jj2=(

De;D

(

e

,

R

1

e

))

= X

Z

(,

u

)(

e

,

R

1

e

)

dx

+Z

@

D

n

e

(

e

,

R

1

e

)

d

= X

Z

f

(

e

,

R

1

e

)

dx

+X

, Z

,

[

D

n

P

1

u

]J(

e

,

R

1

e

)

d;

where[

D

n

v

]Jdenotes the ”jump” of the normal derivative

D

n

v

across,

:

The right hand side can be bounded by Cauchy’s inequality, and (2.3) - (2.5) which gives

jj

De

jj2

c

X

jj

f

jjn

h

21jj

D

1

e

jj2U()+

h

22jj

D

2

e

jj2U()o1=2

+

c

X

,l

jj[

D

n

P

1

u

]Jjj,n

h

,12

h

21jj

D

1

e

jj2U(,)+

h

2jj

D

2

e

jj2U(,)o1=2

+

c

X

,

s

jj[

D

n

P

1

u

]Jjj,n

h

1jj

D

1

e

jj2U(,)+

h

,11

h

22jj

D

2

e

jj2U(,)o1=2

:

In view of the fact that

h

1

h

2we have found the a posteriori bound

jj

De

jj2

c

X

h

21jj

f

jj2+

c

X

,

l

h

,12

h

21jj[

D

n

P

1

u

]Jjj2,+

c

X

,

s

h

1jj[

D

n

P

1

u

]Jjj2,

(2.6)

with local mesh sizes

h

i()

:

Denoting by,()the set of the sides ofwe define the local estimator

by

2= X

,

l

\,()

h

,12

h

21jj[

D

n

P

1

u

]Jjj2,+ X

,s\,()

h

1jj[

D

n

P

1

u

]Jjj2,

:

(2.7)

Though the right hand side definitely deteriorates for

h

2

<< h

1, one can argue that the corresponding jump[

D

n

P

1

u

]Jbecomes smaller in this case. But in the sequel, we will prove that

P

2leads to an arbitrarily large overestimation of the true errorjj

De

jj2if the aspect ratio tends to infinity.

As an example, we consider the orthogonal subdivision of the unit square with mesh sizes

h

^1

; h

^2

; h

^1 =

a

^

h

2

;a >

1

:

In order to transform this mesh to an isotropic mesh we use

x

2=

a x

^2

;x

1=

x

^1and get the operator

Lu

=,

D

112

u

,

a

2

D

222

u

on the rectangle[0

;

1](0

;a

)

:

Denoting by

S

1bthe space of continuous and piecewise bilinear functions satisfying a0–boundary condition the finite element method is defined by

(

D

1

P

1

u;D

1

v

)+

a

2(

D

2

P

1

u;D

2

v

)=(

f;v

) 8

v

2

S

1b

:

(2.8)

Let,ibe the set of edges with direction

e

i

;i

=1

;

2

:

By a similar analysis as before we obtain for the error

e

=

u

,

P

1

u;

jj

D

1

e

jj2+

a

2jj

D

2

e

jj2=X

Z

Lu

(

e

,

R

1

e

)

dx

+X

,2 Z

,

[

D

1

P

1

u

]J(

e

,

R

1

e

)

dx

2

(4)

+

a

2X

,1 Z

,

[

D

2

P

1

u

]J(

e

,

Re

)

dx

1

ch

X

jj

f

jjjj

De

jjU()+

ch

1=2X

,2

jj[

D

1

P

1

u

]Jjj,jj

De

jjU(,)

+

ca

2

h

1=2X

,1

jj[

D

2

P

1

u

]Jjj,jj

De

jjU(,)

:

By Young’s inequality it follows that(

a

1)

jj

D

1

e

jj2+

a

2jj

D

2

e

jj2

ch

2jj

f

jj2+

ch

X

,2

jj[

D

1

P

1

u

]Jjj2,+

ca

4

h

X

,1

jj[

D

2

P

1

u

]Jjj2,

:

The local error indicator is now defined by

2 =

a

4

h

X

,1\,()

jj[

D

2

P

1

u

]Jjj2,+

h

X

,2\,()

jj[

D

1

P

1

u

]Jjj2,

:

(2.9)

We remark that we obtain the same error indicator by simply transforming (2.6) using

x

1 =

x

^1

;x

2=

a x

^2

;h

=^

h

1=

a h

^2

;

jj

D

^

e

jj2!

a

,1,jj

D

1

e

jj2+

a

2jj

D

2

e

jj2

; h

^21jj

f

^jj2!

a

,1

h

2jj

f

jj2

;

X

,

l

^

h

,12 ^

h

21jj[

D

n

P

1

u

^]Jjj2^,

! X

,1

a

3

h

jj[

D

2

P

1

u

]Jjj2,

;

X

,s

^

h

1jj[

D

n

P

1

u

^]Jjj2^,

! X

,2

a

,1

h

jj[

D

1

P

1

u

]Jjj2,

:

LEMMA 2.1. For the finite element approximation

P

1

u

in (2.8) we have the error esti- mate

jj

D

1

e

jj2+

a

2jj

D

2

e

jj2

ca

2

h

2jj

D

2

u

jj2

:

Proof. From the interpolation estimates (1.3), (1.4) it follows that

jj

D

1

e

jj2+

a

2jj

D

2

e

jj2=(

D

1

e;D

1(

u

,

I

1

u

))+

a

2(

D

2

e;D

2(

u

,

I

1

u

))

1

2

jj

D

1

e

jj2+

a

2

2

jj

D

2

e

jj2+

ca

2

h

2k

D

2

u

k2

:

Let

D

+2 be the forward finite difference operator approximating

D

2

;

i.e.

D

2+

v

(

x

)=

h

1(

v

(

x

+

he

2),

v

(

x

))

:

LEMMA2.2. For the finite element approximation

P

1

u

in (2.8) we have the estimate

jj

D

+2

D

1

e

jj20+

a

2jj

D

2+

D

2

e

jj20

ca

2

h

jj

u

jj23;2

(5)

40 A posteriori estimators

for every0

; >

0

;

and0

< h

h

0(0)

:

Proof. Since the subdivision is uniform we have

(

D

+2

D

1

e;D

1

v

)+

a

2(

D

+2

D

2

e;D

2

v

)=0

(2.10)

for all

v

2

S

1b with dist(supp (

v

)

;@

) 2

h:

Let 0 1 be domains that are sufficiently far away from

@

;

let

be a cut–off function with respect tof0

;

1g

;

i.e.

2

C

01(1)with

=1in0

:

For

we have the estimatesj

D

k

j

c; k

=1

;

2

;

with a

constant c depending on0

;

1

:

Using (2.10) we obtain that

Z

0

j

D

2+

D

1

e

j2+

a

2j

D

+2

D

2

e

j2

dx

Z

j

D

+2

D

1

e

j2+

a

2j

D

+2

D

2

e

j2

dx

=

,

D

2+

D

1

e;D

1(

D

2+

e

,

v

)

+

a

2,

D

+2

D

2

e;D

2(

D

+2

e

,

v

)

,(

D

+2

D

1

e;D

2+

eD

1

),

a

2(

D

+2

D

2

e;D

+2

eD

2

)

:

We choose

v

=

I

1(

D

+2

e

)and obtain from (1.3), (1.4) that

jj

D

2+

D

1

e

jj20+

a

2jj

D

+2

D

2

e

jj20

ch

jj

D

2+

D

1

e

jj1jj

D

2(

D

+2

e

)jj1

+

ca

2

h

jj

D

2+

D

2

e

jj1jj

D

2(

D

+2

e

)jj1

(2.11)

+

c

jj

D

2+

D

1

e

jj1jj

D

+2

e

jj1

(2.12)

+

a

2jj

D

2+

D

2

e

jj1jj

D

+2

e

jj1

:

Let2be a slightly larger domain than1

;

such that

jj

D

+2

e

jj1 jj

D

2

e

jj2

(see [6] p.161). By Lemma 2.1, this term can then be bounded by

jj

D

2+

e

jj1

ch

jj

u

jj2;2

:

(2.13)

Moreover, we have the simple inequality

jj

D

2(

D

+2

e

)jj1jj

D

2

D

2+

u

jj1+

c

jj

DD

+2

e

jj1+

c

jj

D

+2

e

jj1

:

Applying the interpolation estimates (1.3), (1.4) and the usual inverse inequality to the second term on the right hand side of the last expression, we obtain

jj

DD

2+

e

jj1jj

DD

2+(

u

,

I

1

u

)jj1+jj

DD

2+(

I

1

u

,

P

1

u

)jj1

ch

jj

u

jj2;2+

ch

,1jj

D

+2(

I

1

u

,

u

)jj1+

ch

,1jj

D

+2(

u

,

P

1

u

)jj1

(2.14)

c

jj

u

jj2;2

;

from which it follows that

jj

D

2(

D

2+

e

)jj1

c

jj

u

jj3;2

:

Inserting the last inequality and (2.13), (2.14) into (2.11), we obtain

jj

D

2+

D

1

e

jj20+

a

2jj

D

+2

D

2

e

jj20

ch

jj

D

+2

D

1

e

jj1jj

u

jj3;2

(2.15)

+

ca

2

h

jj

D

2+

D

2

e

jj1jj

u

jj3;2

:

(6)

In view of the property thatjj

D

j

D

2+

e

jj1 =jj

D

+2

D

j

e

jj1for

j

=1

;

2and (2.14) we have

jj

D

+2

D

j

e

jj1

c

jj

u

jj2;2

c

jj

u

jj3;2

which leads to

jj

D

+2

D

1

e

jj20+

a

2jj

D

2+

D

2

e

jj20

ch

jj

u

jj23;2+

ca

2

h

jj

u

jj23;2

ca

2

h

jj

u

jj23;2

:

Remark: Note that one can get a better estimation than stated in the Lemma by iterating (2.15) for a sequence of domains0 1

:::

k

:::

. Arguing in this way we would get the inequality

jj

D

2+

D

1

e

jj20+

a

2jj

D

+2

D

2

e

jj20

ca

2

h

2,jj

u

jj23;2

for all

>

0.

Now letbe a square with upper neighboring square~and common side,

:

For

v

2

S

0b

we have

Z

j

D

+2

D

2

v

j2

dx

1

dx

2=

h

12

Z

j

D

2

v

(

x

+

he

2),

D

2

v

(

x

)j2

dx

1

dx

2

:

In view of the fact that

D

2

v

depends only on the variable

x

1we get

jj

D

+2

D

2

v

jj2= 1

h

Z

,

[

D

2

v

]2J

dx

1

:

Let0

be a fixed domain. Denoting the set of sides of0in direction

e

1by,0we

obtain from Lemma 2.2 that

ha

4X

,0

jj[

D

2

P

1

u

]Jjj2,

a

4

h

2jj

D

+2

D

2

P

1

u

jj20

(2.16)

1

2

a

4

h

2jj

D

+2

D

2

u

jj20,

a

4

h

2jj

D

+2

D

2

e

jj20

1

2

a

4

h

2jj

D

+2

D

2

u

jj20,

ca

4

h

3jj

u

jj23;2;

:

Choosing a smooth function

u

which behaves likesin

x

2in0such that

jj

D

2+

D

2

u

jj0

c

1

;

jj

u

jj3;2;

c

2

;

with constants

c

1

;c

2

>

0, then (2.16) shows, that for

h

sufficiently small

ha

4X

,

0

jj[

D

2

P

1

u

]Jjj2, 1

4

a

4

h

2

c

21

:

From Lemma 2.1 we conclude that the error estimator gives an overestimation with factor

a

2

for functions of this type.

3. A nonlocal error indicator. In this section we return to Poisson’s equation (1.1) and its finite element approximation (1.5). Recall that

I

k

;k

=1

;

2

;

are the standard interpolation operators into the spaces

S

k, and define the space

S

0=f

v

2

S

2:

I

1

v

=0g

;

(7)

42 A posteriori estimators

which means that the elements of

S

0vanish in the nodal points of the triangulation.

The nonlocal version of the third error indicator of Bank–Weiser [4] is given byo

e

2

S

0

such that

(

D

o

e;Dv

)=(

De;Dv

) 8

v

2

S

0

:

(3.1) In view of

(

De;Dv

)=(

f;v

),(

DP

1

u;Dv

)

the right-hand side of (3.1) is known if

P

1

u

is known.

Let us compare

e

owith the original third indicator of [4] which also gives some insight into error propagation on anisotropic meshes. Let

S

~be the discontinuous version of

S

0

;

i.e.

S

~

consists of all piecewise quadratic functions vanishing in the nodal points of the triangulation.

The indicator~

e

2

S

~is then defined by

(

D

~

e;Dv

)=

F

(

v

) 8

v

2

S;

~

(3.2) where

F

(

v

)=(

f;v

)+1

2 X

,2,() Z

,

[

D

n

P

1

u

]J[

v

]A

d;

and where[

v

]Ais the average of

v

on the neighboring triangles of,

:

Summing (3.2) over and using integration by parts yield

X

(

D

~

e;Dv

)=X

(

De;Dv

),X

, Z

,

[

D

n

e

]J[

v

]A 8

v

2

S:

~

Comparing this with (3.1) shows thato

e

is the continuous and nonlocal counterpart of~

e:

For the

actual computation of~

e;

a33linear system has to be solved on each trianglein contrast to the large system required for the computation of

e :

o On the other hand, a complicated computation using the symbolic program ”mathematica” shows that the system corresponding to (3.1) is strictly diagonally dominant if the largest interior angle is bounded by

<

2

:

Let

0

<

<

2

:

Since the triangles with interior angles between

and

are compactly parametrized we obtain that the system in (3.1) is uniformly strictly diagonally dominant in this class of triangles and can efficiently be solved by the simple Gauß–Seidel method.

Furthermore, we conclude that local errors decrease exponentially on such meshes. This is the reason why local error estimation is possible. For isotropic triangles with angles bounded by

<

the reasoning is similar. Since local error estimation is also possible in this case, we conclude that the system in (3.1) must be strictly dominant “in the mean” and can again be solved by the Gauß–Seidel method.

Now we turn to the anisotropic case and consider the orthogonal mesh with parameters

h

2

<< h

1shown in Fig. 1. The entries of the matrix in (3.1) can be represented by stencils.

For the midpoints of the larger sides we obtain

S

l=

0

@

0 ,1 0

0 2 0

0 ,1 0 1

A

+

O

h

2

h

1

;

whereas the entries of the stencil for the shorter sides are of the type

S

s=

0

@

0 0 0

0 2 0 1

A

+

O

h

2

h

1

:

(8)

h1

h2

FIG. 3.1. An anisotropic mesh.

Thus, local errors propagate in the direction of the small sides with a stencil approximating

,

h

22

D

yy2

:

Since the indicatoro

e

is reliable and efficient in this case, we believe that local error estimation is inherently inaccurate on anisotropic meshes.

On the anisotropic mesh shown in Fig. 1, the indicator

e

ocan efficiently be determined with the aid of a Block–Gauß–Seidel method since there is nearly no coupling normal to the smaller sides. On general anisotropic meshes, we use a mesh–dependent Block–Gauß–Seidel method where points coupled by smaller sides are updated simultaneously.

Since the local error indicator is equivalent to the residual based indicator on anisotropic meshes, the counterexample given in section 2 can also be applied. It remains to prove that the nonlocal indicator is reliable and efficient. We start with a preparatory result.

LEMMA3.1. There is a constant

c

0depending only on the angle

in condition (1.2) b), but not on the shape of the trianglesuch that

jj

DI

1

v

jj

c

0jj

Dv

jj 8

v

2

S

2

:

Proof. On isotropic triangles, this estimate can simply be proved by using a reference element and transforming the corresponding estimate to

:

In the anisotropic case, our proof requires some notations, but is simpler and shows that

c

01for moderate

:

Let

P

1

;P

2

;P

3

be the nodes ofnumbered counterclockwise. The edge opposite to

P

iis denoted by

e

iwith

midpoint

a

i

:

The derivative in direction

e

iis denoted by

D

i

:

For

w

2

IP

2()we have

D

i

w

(

a

i)= 1

j

e

ij(

w

(

P

i+1),

w

(

P

i,1))

;

(3.3)

Z

w

(

x

)

dx

=

()

3 3

X

i=1

w

(

a

i)

:

(3.4)

Assuming that(

e

1

;e

2)is a pair of a larger and a smaller side we obtain

jj

DI

1

v

jj2

c

jj

D

1

I

1

v

jj2+jj

D

2

I

1

v

jj2

;

where

c

depends on the angle

in condition (1.2) b) . Using (3.3), (3.4) and the fact that

D

i

I

1

v

is constant in

;

we obtain

jj

DI

1

v

jj2

c

()j

D

1

I

1

v

(

a

1)j2+j

D

2

I

1

v

(

a

2)j2

c

jj

Dv

jj2

:

参照

関連したドキュメント