close

Вход

Забыли?

вход по аккаунту

?

512

код для вставкиСкачать
INTERNATIONAL JOURNAL FOR NUMERICAL METHODS IN ENGINEERING
Int. J. Numer. Meth. Engng. 47, 427– 475 (2000)
Guaranteed computable bounds for the exact error in the
nite element solution—Part II: bounds for the energy norm
of the error in two dimensions†
T. Strouboulis1 , I. Babuska2;∗ and S. K. Gangaraj1
1
Department of Aerospace Engineering; Texas A&M University; College Station; TX 77843-3141; U.S.A.
2 Texas Institute for Computational and Applied Mathematics; The University of Texas at Austin;
Austin; TX 78712; U.S.A.
SUMMARY
This paper addresses the computation of guaranteed upper and lower bounds for the energy norm of the exact
error in the nite element solution. These bounds are constructed in terms of the solutions of local residual
problems with equilibrated residual loads and are rather sharp, even for coarse meshes. The sharpness of the
bounds can be further improved by employing a few iterations of a relatively inexpensive iterative scheme.
The main result is that the bounds are guaranteed for the energy norm of the exact error, unlike the bounds
which have been proposed in [13, 14] which are guaranteed only for the energy norm of the error with respect
to an enriched (truth-mesh) nite element solution. Copyright ? 2000 John Wiley & Sons, Ltd.
KEY WORDS: a posteriori error estimation; nite element method; guaranteed error bounds; upper and lower error
bounds
1. INTRODUCTION
In this paper we address the problem of computation of upper and lower bounds EuUFE , and EuLFE ,
for the energy norm of the exact error |||uEX − uFE |||
, namely
EuLFE 6|||uEX − uFE |||
6EuUFE
(1)
where uEX and uFE are the exact and nite element solutions of the Laplacian in two dimensions,
respectively. The main tools for computing the bounds are solutions of local residual problems of
the form:
∗ Correspondence to: I. Babuska, TICAM, The University of Texas at Austin, TAY 2-400, Austin, TX 78712, U.S.A.
E-mail: babuska@brahma.ticam.utexas.edu
† Dedicated to the memory of Dr. Richard H. Gallagher
Contract=grant
Contract=grant
Contract=grant
Contract=grant
Contract=grant
sponsor:
sponsor:
sponsor:
sponsor:
sponsor:
U.S. Army Research Oce; contract=grant number: DAAL03-G-028
National Science Foundation; contract=grant number: MSS-9025110
U.S. Oce of Naval Research; contract=grant numbers: N00014-96-1-0021 and N00014-96-1-1015
U.S. Oce of Naval Research; contract=grant number: N00014-90-J-1030
National Science Foundation; contract=grant numbers: DMS-91-20877 and DMS-95-01841
CCC 0029-5981/2000/020427–49$17.50
Copyright ? 2000 John Wiley & Sons, Ltd.
Received 17 February 1999
428
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
!
Find ê Q
!; uFE ∈ Q! such that
EQ
!
B! ê Q
!; uFE ; v = R!; uFE (v) ∀v ∈ Q!
(2)
where ! is the typical subdomain in which the residual problem is posed, Q! is the solution
EQ
space, B! is the bilinear form over !, and R!;
uFE is an equilibrated local residuum. For the
one-dimensional version of this paper, see [1] and for previous work on the equilibrated residual
problems, see [2–11]. Letting Q! = U! , the energy space in the subdomain, we get the exact
!
indicators |||ê U
!; uFE |||! , and we can use them to get the guaranteed upper bound
rP
!
2
|||uEX − uFE |||
6E uUFE =
|||ê U
(3)
!; uFE |||!
!
def p
where ||| · |||
= B
(·; ·), is the energy norm over . Here E uUFE is the exact upper bound or
upper estimator, which is not computable and is often approximated by employing Q! = S! ⊂ U! ,
where S! is a discrete subdomain space obtained by the h, p or hp restriction of the nite element
space in the subdomain, and we have
rP
EuSFE =
|||ê S!;! uFE |||2! 6E uUFE
(4)
!
We will call EuSFE! the computed estimate. It can be proved that EuSFE! , is a guaranteed upper estimate
for the error with respect to an enriched or truth mesh solution uFE
e
S
|||uFE
e − uFE |||
6EuFE
(5)
Note however that in many cases we have
EuSFE 6|||uEX − uFE |||
6EuUFE
(6)
and hence EuSFE cannot be relied upon as an upper bound for |||uEX − uFE |||
, which is the quantity
S
of interest and not |||uFE
e − uFE |||
. In [1] and in this paper we show that EuFE can be corrected by
realizing that
r
P
S!
2
!
(7)
E uU = (EuS )2 + |||ê U
!; u − ê !; u |||!
FE
FE
FE
!
FE
and by constructing a computable upper estimate
S!
!
|||ê U
!; uFE − ê !; uFE |||! 6Eeˆ S!;!u
(8)
FE
to get the computable upper bound
|||uEX −
uFE |||
6EuUFE
s
=
(EuSFE )2 +
P
!
Eeˆ S!;!u
2
(9)
FE
A lower bound for the energy norm can be obtained from the observation that, for any function
’ in the energy space U
R(’)
R(v)
6 sup
= |||R|||−1; = |||uEX − uFE |||
|||’|||
v∈U |||v|||
Copyright ? 2000 John Wiley & Sons, Ltd.
(10)
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
429
where R is the global residuum. However, unless ’ is close to the exact error uEX −uFE , this bound
could be too small and hence useless. To construct an admissible ’, we rst construct the global
Q
Q!
indicator ê Q
uFE patching together the subdomain error indicator functions, namely ê uFE |! = ê !; uFE .
Q
Q
However ê uFE 6∈ U, is not admissible because ê uFE is discontinuous at the vertices and the edges
across the subdomains. Nevertheless, we can construct a gap function eˆ Qu , such that êQ
uFE +
FE
eˆ Qu ∈ U, is admissible and we have the lower bound
FE
EQ
uFE =
R ê Q
uFE + eˆ Q
u
FE
|||ê Q
uFE + eˆ Q |||
u
6|||uEX − uFE |||
(11)
FE
Another lower bound may be obtained by exploiting the orthogonality properties of the error and
the error indicators to obtain the exact lower estimator
r
P
=
(E uUFE )2 − |||eˆ Uu |||2! 6|||uEX − uFE |||
(12)
EU
uFE
!
FE
and its computable version
r
P
E SuFE = (EuSFE )2 − |||eˆ Su |||2! 6E U
uFE 6|||uEX − uFE |||
(13)
The lower estimate (12) was obtained from the identity
r
P
|||uEX − uFE |||
= (EuSFE )2 − |||eˆ Su |||2! + |||qEX |||2
(14)
where qEX ∈ U, is the exact equilibrium correction which is obtained from
P
∀v ∈ U
B
(qEX ; v) = − B! eˆ Su ; v
(15)
!
FE
!
!
FE
FE
Employing now nite element approximation qFE of qEX , which is computed from the same space
as uFE , we have
|||qEX |||2
= |||qFE |||2
+ |||qEX − qFE |||2
(16)
Note that qFE can be computed by utilizing the already available factorized global stiness matrix.
We can now repeat the above steps (9) and (13) to construct computable bounds for |||qEX −qFE |||
,
namely
E SqFE 6|||qEX − qFE |||
6EqUFE
Using (16) and (17) into (14) we get the enhanced bounds
r
P
ẼuU;FE(1) = (EuSFE )2 − |||eˆ Su |||2! + |||qFE |||2
+ (EqUFE )2
!
and
ẼuL;FE(1) =
FE
r
P
(EuSFE )2 − |||eˆ Su |||2! + |||qFE |||2
+ (E SqFE )2
Copyright ? 2000 John Wiley & Sons, Ltd.
!
FE
(17)
(18)
(19)
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
430
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
Figure 1. The polygonal domain for the model problem:
(a) The domain , its boundary @
which consists of the
Dirichlet parts D1 ; D2 ; D3 (shown with thick gray line),
and the Neumann part N (shown with thin line); (b) the
coarse regular mesh of curvilinear quadrilaterals which is
used for the computation of the nite element solution
Figure 2. Examples of partition of the mesh into patches:
(a) The patches coincide with the elements; (b) the
patches consist of the vertex patches for the re-entrant corners, and the elements that remain after the vertex patches
for the re-entrant corners have been formed
Repeating the above procedure in a recursive way we obtain an iteration for improving the bounds.
Each iteration employs the factorization of the global stiness matrix which was obtained during
the computation of the nite element solution and the factorization of the local stiness matrices
which are obtained during the solution of the local residual problems. At the end of each iteration
we can compute the reliability interval
Ru(m)
=
FE
m)
ẼuU;FE(m) − ẼuL;(
FE
|||uFE |||
(20)
As we will see below, a few iterations are sucient to produce dramatic improvement in the
accuracy of the bounds and the size of the reliability interval.
Following this Introduction, in Section 2 we give the description of the boundary value problem,
and the local residual problems. In Section 3 we give the construction of the computable guaranteed
upper bound, in Section 4 we give the computable guaranteed lower bound and in Section 5 we
give the recursion for improving the bounds.
2. THE MODEL PROBLEM, ITS FINITE ELEMENT APPROXIMATION, THE RESIDUAL
EQUATION AND THE EQUILIBRATED LOCAL RESIDUAL PROBLEMS
We will employ the mixed boundary-value problem for the Laplacian as our model problem. Let be the polygonal domain shown in Figure 1(a) with boundary ≡ @
; = D ∪ N ; D ∩ N = ∅,
with D = D1 ∪ D2 ∪ D3 , as shown in Figure 1(a). The model problem reads:
Find uEX which satises
− uEX = 0
uEX |
1
D
∪
@uEX @n N Copyright ? 2000 John Wiley & Sons, Ltd.
3
D
= 0;
=0
in (21)
uEX | D2 = 1
(22)
(23)
N
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
The variational formulation
of the general mixed boundary-value
problem reads
n
o
def
def p
Find uEX ∈ U
; u = u|kukU = B
(u; u)¡∞; u| D = u such that
Z
Z
Z
def
def
B
(uEX ; v) =
∇uEX · ∇v = L
(v) =
fv +
gv ∀v ∈ U
; 0
431
(24)
N
In the example employed throughout the paper we used zero loads f ≡ 0; g ≡ 0, zero Dirichlet
1D = 0; u|
3D = 0, and constant non-zero Dirichlet boundary
boundary conditions on 1D and 3D , u|
2
1
2
2D = 1, where , D , D , 3D , and N , are shown in Figure 1(a).
condition on D , u|
def
Let Th
= {j }nelem
j = 1 , denote a nite element mesh, which is a subdivision of into nonnelem
S
=
overlapping curvilinear quadrilaterals, namely j ; j ∩ k = ∅ for j 6= k. Here we will
j=1
employ the coarse regular mesh shown in Figure 1(b) and the meshes which are obtained from
the nested subdivision this mesh. By regular mesh we mean a mesh for which the intersection of
the boundary of any two elements, @j ∩ @k ; j 6= k is either empty, or a vertex, or an entire edge
for both the elements j , and k . Given a mesh Th
we will construct the corresponding set of
mapped bi-p nite element (trial) solutions
n
o
p
def
(25)
Sh;p u (
) = u ∈ C 0 (
) | u| ◦ F ∈ Ŝ ∀ ∈ Th
; u| D = u ⊆ U
; u
and the corresponding space of test-functions Sh;p 0 (
) ⊆ U
; 0 . Here F : ˆ 7→ is the transformation
def
which maps the master square ˆ = (−1; 1)2 into the curvilinear quadrilateral (in the results of
p
this paper F was constructed using the blending function approach given in [12]), Ŝ is the bi-p
element space over ˆ (here we employed p = 2). We will denote the nite element solution by
uS hp and we will obtain it by solving the following discrete problem:
Find uS hp ∈ Sh;p u (
) such that
(26)
B
uS hp ; v = L
(v) ∀vh ∈ Sh;p 0 (
)
In general we have uS hp 6= uEX . This can be veried by computing the residuals corresponding
to uS hp , namely, the interior residuals
def
r;uS p = f + uS hp |
(27)
h
and the jumps in the normal derivative
T

(∇uS hp |k − ∇uS hp |j )(x)·n (x); x ∈ = @j @k ; j 6= k



def
x∈⊆ N
(28)
J; uS p (x) = 2(g(x) − ∇uS hp (x)·n N (x));

h


0;
x∈⊆ D
T
where n is the unit normal on = @j @k , which points into k , and n N is the exterior unit
normal on N .
We will denote the exact error in uS hp by
def
p = uEX − u p
eSEX
Sh
h
Copyright ? 2000 John Wiley & Sons, Ltd.
(29)
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
432
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
We will also employ the error relative to an enriched nite element solution uS p+kn , namely
h=2
e
p+k
Sh=2
n
p
Sh
def
= uS p+kn − uS hp
(30)
h=2
p+k
Here uSp+k
is the nite element solution from Sh=2
n ; u (
), the set of trial solutions of degree p+k on
h=2n
the mesh Th=2
which
is
obtained
from
T
,
using
n uniform nested subdivisions. The enriched nite
n
h
S p+kn
element solution uS p+kn is often called the truth-mesh solution, and the error e Sh=2
p , the truth-mesh
h=2
h
error; see for example [13,14].
S p+kn
h=2
. Substituting uEX = eSEX
Let us now give the variational problems satised by eSEX
p , and e p
p +uS p ,
S
h
h
h
h
into the variational statement B
(uEX ; v) = L
(v), for all v ∈ U
; 0 , we get that the exact error eSEX
p ,
h
is the solution of the variational problem:
Find eSEX
p ∈ U
; 0 such that
h
def
p ; v) = R
; u p (v) = L
(v) − B
(u p ; v)
B
(eSEX
Sh
S
h
∀v ∈ U
; 0
h
(31)
This variational problem is called the residual equation, andp+k the linear functional R
;uS p :
S
n
h
U
; 0 7→ R, the global residuum. Similarly, letting uS p+kn = uS hp + e Sh=2
p , into the variational statement
h=2
for uS p+kn , namely
h
h=2
p+k
∀v ∈ Sh=2
n ; 0 (
)
B
(uS p+kn ; v) = L
(v)
h=2
(32)
S p+kn
is the solution of the variational problem:
it can be seen that e Sh=2
p
p+k
Sh=2
n
Find e S p
h
h
p+k
∈ Sh=2
n ; 0 (
)
such that
S p+kn ; v = R
; uS p (v)
B
e Sh=2
p
Noting that B
(uS hp ; v) =
P
∈Th
h
p+k
∀v ∈ Sh=2
n ; 0 (
)
h
B (uS hp ; v), and integrating by parts the term B (uS hp ; v) we get
!
Z
P 1
r; uS p v +
J; uS p v
R
; uS p (v) =
h
h
h
h
∈Th
∈Th
⊆ @ 2 n
o
def p
def
B (w; w)¡∞; w|@∩ D = 0 7→ R, given by
Here R;UNEQ
uS p : U; 0 = w| ||w||U =
h
Z
Z
P 1
def
(v)
=
r
v
+
J; uS p v
R;UNEQ
; uS p
uS p
h
h
h
⊆ @ 2 P
(33)
R;UNEQ
uS p (v| ) =
P
Z
(34)
(35)
is the unequilibrated element residuum, namely the linear functional which corresponds to the
loads r; uS p and 12 J; uS p , ⊆ @. We call R;UNEQ
uS p unequilibrated because, in general, we have
h
h
h
Z
R;UNEQ
u p (1) =
S
h
Copyright ? 2000 John Wiley & Sons, Ltd.
r; uS p +
h
P 1
⊆ @ 2
Z
J; uS p 6= 0
(36)
h
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
433
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
Next we construct a new set of element residua
Z
Z
Z
P
P
def
EQ; p
UNEQ
J; u p v = R
(v) +
sg(; ); uS p v
R; u p (v) = r; uS p v +
S
h
S
⊆ @
h
⊆ @
h
(37)
h
where
J; u p = 12 J; uS p + sg(; ); uS p
S
h
h
(38)
h
and sg(; ) = ± 1, such that sg(; ) + sg(; ∗ ) = 0, for every interior edge = @ ∩ @∗ , and
sg(; ) = 0, for ⊆ @
, and ; uS p ∈ L2 (), is the edge correction for the edge , which will be
h
specied below. Note that
J; u
∗
S
p
h
+ J; u p = J; uS p
S
∀ = @ ∩ @∗
h
h
and hence we have an exact splitting of the global residuum, namely
P EQ; p
R; u p (v| ) ∀v ∈ U
; 0
R
; uS p (v) =
S
∈Th
h
(39)
h
Following Ladeveze [4] the edge corrections ; uS p will be constructed such that the orthogonality
h
condition
p
p
R;EQ;
u p (v) = 0 ∀v ∈ Sh; 0 () =
def
S
h
n p
v v | ◦ F ∈ Ŝ ; v|@∩
o
=0
D
(40)
is satised in each element (see also [4–10] and Appendix I for the construction). Note that for
each element with @ ∩ D = ∅, we have 1 ∈ Sh;p 0 (), and (40) implies the equilibrium condition
Z
Z
P
(1)
=
r
+
J; u p = 0
(41)
R;EQ;p
p
;
u
u p
S
S
h
⊆ @
h
S
h
We will now introduce the local residual problems which will be used in the a posteriori error
def
estimation. Let h = {!j }npatches
denote a partition of the mesh Th
into non-overlapping patches
j=1
of elements
=
npatches
S
j=1
! j ;
! j =
nel(!
S j)
k=1
!
k j
T
!j !k = ∅; j 6= k
(42)
!
where nel(!j ), is the number of elements in the patch !j , and k j , is the kth element in !j .
Figures 2(a) and 2(b) show two examples of partition of the mesh Th
of Figure 1(b) into two
sets of patches. Figure 2(a) shows a partition for which the patches are the same as the elements.
Figure 2(b) shows a partition which employs the vertex patches for the re-entrant corners (each
of these patches consists of the elements connected to the vertex of each re-entrant corner), and
that rest of the patches coincide with the elements that remain after the vertex patches have been
formed.
The a posteriori error estimates
will be based on the following patch
n o residual problems:
def p
def
EX
Find ê !j ; u p ∈ U!j ;0 = v ||v||U!j = B!j (v; v)¡ ∞; v|@!j ∩ D = 0 such that
S
h
def P
EQ; p
p
R;EQ;
B!j ê EX
!j ; u p ; v = R!j ; u p (v) =
u p (v| )
S
h
Copyright ? 2000 John Wiley & Sons, Ltd.
S
h
∈Th
⊆ !j
S
∀v ∈ U!j ;0
(43)
h
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
434
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
In the case that @!j ∩ D = , (43) is a mixed boundary value problem and it has a unique solution.
In the cases that @!j ∩ D = ∅, (43) is a Neumann problem; in this case however R!EQ;p
j ; uS p satises
h
P EQ;p
EQ;p
EX
R; uS p (1) = 0, and ê !j ; u p exists and is unique up to
the consistency condition R!j ; uS p (1) =
h
S
h
∈Th
⊆ !j
h
a constant, which will be xed by imposing the condition
Z
ê EX
!j ; u p = 0
S
!j
(44)
h
We will refer to (43) as the patch residual problem in !j ; in the case that !j coincides with an
element, we will call (43) the element residual problem. The functions ê EX
!j ; u p will be called the
S
h
patch residual error indicator functions, and the quantities
EX
EX
!j ;u p = ||ê !j ; u p ||U!j
S
S
h
(45)
h
the exact patch error indicators. Summing the squares of the indicators, we get the exact patch
residual estimator
s
2 r P
P
def
EX
2
E
h ; u p =
=
||ê EX
(46)
EX
!j ; u p ||U!
!j ;u p
S
S
!j ∈
h
h
S
!j ∈
h
h
h
j
EX
EX
Below we will prove that ||eSEX
p ||U 6E
h ; u p , which means that E
h ; u
S
h
S
h
estimator for ||eSEX
p ||U , the energy norm of the exact error.
p
h
is a guaranteed upper
h
EX
The exact indicators EX
!j ; u p , and the exact estimator E
h ; u
S
!
Th=2j n
S
h
p
h
are not computable and must be
approximated. Let
denote a nite element mesh in !j which is obtained by employing n
uniform nested subdivisions of the restriction of Th
in !j , and
n
o
p+k
def
!
p+k
v ∈ C 0 (!j ) v |!j ◦ F!j ∈ Ŝ
∀!j ∈ Th=2j n ; v|@!j ∩ D = 0
(47)
Sh=2
n ;0 (!j ) =
p+k
the restriction of the test space Sh=2
n ; 0 (
) in !j . We will approximate the exact error indicator
S p+kn
h=2
function ê EX
!j ; u p , by the computed error indicator function ê !j ; uS p , which is the solution of the
S
h
h
following discrete problem:
S p+kn
p+k
Find ê !h=2
j ; uS p ∈ Sh=2n ; 0 (!j ) such that
h
S p+kn
EQ;p
B!j ê !h=2
j ; uS p ; v = R!j ; u p (v)
S
h
In the cases that @!j ∩
h
p+k
∀v ∈ Sh=2
n ;0 (!j )
(48)
S p+kn
D
= ∅, we will x ê !h=2
j ; uS p by imposing the additional condition
h
Z
p+k
Sh=2
n
!j
Copyright ? 2000 John Wiley & Sons, Ltd.
ê !j ; uS p = 0
(49)
h
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
435
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
We will also employ the computed patch residual error indicators
S p+kn
S p+kn
h=2
!h=2
j ; uS p = ||ê !j ; uS p ||U!j
h
(50)
h
and the computed patch residual estimator
s
S p+kn
def
=
E
h=2
h ;uS p
P
!j ∈
h
h
S p+kn
S p+kn 2
!h=2
j ; uS p
s
P
=
!j ∈
h
h
S p+kn
S p+kn
2
||ê !h=2
j ; uS p ||U!
h
S p+kn
||U 6E
h=2
, and hence E
h=2
Below we will prove that ||e Sh=2
p
h ;u p
h ;u
S
h
for ||e
p+k
Sh=2
n
p
Sh
h
S
(51)
j
is a guaranteed upper estimator
p
h
||U , the energy norm of the truth-mesh error. We will also show that in many cases
S p+kn
S p+kn
EX
p ||U 6E
||e Sh=2
||U 6E
h=2
6||eSEX
p
h ; u
h; u p
S
h
S
h
h
(52)
p
h
S p+kn
, to be a guaranteed upper estimator for ||eSEX
and hence we cannot rely E
h=2
p ||U .
h ;u p
S
h
h
3. THE UPPER BOUND FOR THE ENERGY NORM OF THE EXACT ERROR
EX
We will now prove the main results about the estimates E
h; u
p
S
h
omit the subscript uS hp from the various quantities.
S p+kn
and E
h=2
. Below we will often
h; u p
S
h
Theorem 3.1 (Upper bounds). Consider the setting and the notations given in Section 1. We
have
EX
p ||U 6E
||eSEX
h ; u
p
h
S
h
(53)
and
S p+kn
S p+kn
||e Sh=2
||U 6E
h=2
p
h; u
S
h
p
h
(54)
Further; we have
S p+kn
EX
!h=2
j ;uS p 6!j ; u
h
S
p
h
(55)
and hence
S p+kn
E
h=2
6E
EX
h; u
h; u p
S
Copyright ? 2000 John Wiley & Sons, Ltd.
h
S
p
h
(56)
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
436
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
Proof. Employing the denition of the energy norm using duality, the residual equation and the
exact splitting of the global residuum, we get
B
eSEX
p ;v
p ||U =
sup
||eSEX
h
P
R
(v)
= sup
6 sup
v∈U
; 0 ||v||U
v∈U
; 0
h
||v||U
v∈U
; 0
!j ∈
h
|R!EQ;p
(v|!j )|
j
(57)
||v||U
It follows that
v
EQ;p
2
u
(v|
)
R
!j
!j ||v||U!
u P
P
(v)
R!EQ;p
j
j
EX
t
6
sup
||eS p ||U 6 sup
h
||v||U!j
||v||U
||v||2U!
v∈U
; 0 !j ∈
h
!j ∈
h v∈U!j ;0
(58)
j
where we have used the Cauchy–Schwarz inequality. From the denition of ê EX
!j , as the solution
EQ;p
of B!j (ê EX
!j ; v) = R!j (v), for all v ∈ U!j ;0 , we get
def
||L(U!j ;0 ;R) = sup
||R!EQ;p
j
v∈U!j ;0
B!j (êEX
R!EQ;p
(v)
!j ; v)
j
= sup
= ||êEX
!j ||U!j
||v||U!j
||v||U!j
v∈U!j ;0
and substituting into (58) we get (53).
S p+kn
S p+kn
||U 6E
h=2
The proof of the bound ||e Sh=2
p
h ;u
h
S
p
h
(59)
follows exactly the same steps (57)–(59) where the
p+k
p+k
supremums are taken over Sh=2
n ;0 (
) and Sh=2n ;0 (!j ) instead of U
; 0 and U!j ;0 respectively. In
summary, we have
v
2
u
EQ;p
p+k
u P
(v)
R!EQ;p
Sh=2
R
(v)
n
j
t
6
supv∈S p+kn (!j )
(60)
||e S p ||U = sup
h=2 ;0
h
||v||U
||v||2U!
!j ∈
h
v∈S p+k (
)
h=2n ;0
j
S p+kn
we get
and using the denition of ê !h=2
j
S p+kn B!j ê !h=2
;v
S p+kn
j
sup
= ||ê !h=2
||U!j
j
||v||U!j
v∈S p+k (!j )
R!EQ;p
(v)
j
=
sup
||v||
p+k
U!j
v∈S
(!j )
h=2n ;0
(61)
h=2n ;0
S p+kn
p+k
EX
Let us now prove that !h=2
j ;uS p 6!j ;u p . Noting that Sh=2n ;0 (!j ) ⊂ U!j ;0 , and subtracting (48) from
S
h
h
(43) we get the orthogonality
B!j ê EX
!j ; u
S
p
h
S p+kn
− ê !h=2
j ; uS p ; v = 0
h
p+k
∀v ∈ Sh=2
n ;0 (!j )
(62)
and hence
S p+kn
h=2
2
2
EX
||ê EX
!j ; u p ||U!j = ||ê !j ; uS p ||U!j + ||ê !j ; u
S
h
Copyright ? 2000 John Wiley & Sons, Ltd.
h
S
S p+kn
p
h
2
− ê !h=2
j ; uS p ||U!
j
(63)
h
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
437
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
p+k
p+k
p+k
Sh=2
Sh=2
Sh=2
n
n
n
p+k
where we have used that B!j ê EX
!j ; u p − ê !j ; uS p ; ê !j ; uS p = 0, because ê !j ; uS p ∈ Sh=2n ;0 (!j ). This can
S
h
h
h
h
also be written as
S p+kn
!h=2
j ;uS p
s
EX
!j ;u
=
h
2
p
S
h
S p+kn
h=2
2
− ||ê EX
!j ; u p − ê !j ; uS p ||U!
S
(64)
j
h
h
S p+kn
EX
and we have !h=2
j ;uS p 6!j ;u p . Further, summing (63) over all the patches in h , we get
S
h
P
!j ∈
h
h
2
||ê EX
!j ; u p ||U! =
S
j
h
P
P
S p+kn
!j ∈
h
2
||ê !h=2
j ; uS p ||U! +
j
h
!j ∈
h
S p+kn
||ê EX
!j ; u
p
S
h
2
− ê !h=2
j ; uS p ||U!
h
j
(65)
which can also be written as
S p+kn
E
h=2
h ;uS p
s
2
EX
E
h ;u
=
S
h
p
h
−
P
S p+kn
!j ∈
h
h=2
2
||ê EX
!j ;u p − ê !j ;uS p ||U!
S
j
h
h
(66)
S p+kn
EX
6E
.
and we get that E
h=2
h ;u p
h ;u p
S
S
h
h
S p+kn
h=2
It is rather easy to prove that the estimate ||eSEX
p ||U 6E
h ;u
S
h
p
h
does not hold in general.
Counterexample 3.1. Consider the case of the example problem where Th
is the mesh shown
in Figure 1(b), and let k = 1; n = 3. For this example
S3
EX
p ||U ¡ ||e p ||U
≈ 0·894||eSEX
ETh=8
S
h ;u p
S
S p+kn
Hence we may have E
h=2
h ;u
S
for ||eSEX
p ||U .
h
h
h
S p+kn
p
h
h=2
¡ ||eSEX
p ||U , and E
h ;u
S
h
p
h
cannot be relied to be an upper estimator
h
We will now get a computable upper bound for ||eSEX
p ||U by employing the identity
h
s
EX
=
E
h ;u p
S
h
S p+kn
E
h=2
h ;u
p
S
h
2
+
P
S p+kn
!j ∈
h
and constructing an upper estimate for ||ê EX
!j ; u
p
S
h
h=2
2
||ê EX
!j ;u p − ê !j ;uS p ||U!
S
h
h
j
(67)
S p+kn
− ê !h=2
j ; uS p ||U!j , in each patch !j . We will rst
h
consider the case that h = Th
, in which the patches are the elements. We have the following
result:
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
438
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
S p+k
h
Theorem 3.2 (Computable explicit upper estimator for ||ê EX
; u p −ê ; uS p ||U ). Consider the
S
h
h
setting and the notations given in Section 2. Then
S p+k
||ê EX
; u
p
S
h
− ê ;huS p ||U 6EEXPL
p+k
S
; eˆ;uh
h
(68)
p
h
S
where
def
; eˆ;uh
P
S p+k
h
EEXPL
p+k = C1 (p + k)||r; u p + ê ; u p ||L2 () +
S
S
S
h
h
p
h
S
The constants C1 and C2 are given by
s
max |J |
C1 (q)6
Ĉ1 (q);
K
⊂ @
C2 (p + k)||J; u
s
C2 (q)6
p
S
h
−
@ Shp+k
ê ; uS p ||L2 ()
@n
h
max |J |
Ĉ2 (q)
K
(69)
(70)
where
K = min(|J | min (J−T J−1 ))
(71)
J (resp. Jj ) is the Jacobian of the transformation of the element (resp. edge); and min (J−T J−1 )
is the minimum eigenvalue of the matrix J−T J−1 ; and
Ĉ1 (q) = sup
v∈V
ˆ ˆ
with
V qˆ =
q
||v̂||L2 ( ˆ)
;
||v̂||U ˆ
Ĉ2 (q) = sup
q
v∈V
ˆ ˆ
||v̂||L2 ()ˆ
||v̂||U ˆ
Z
q
v ∈ U ˆ B ˆ (v; w) = 0 ∀w ∈ Ŝ ; v = 0
(72)
def
(73)
ˆ
Proof. Letting the patches !j coincide with the elements of the mesh, we get the exact element
error indicator functions ê EX
; u p , as the solutions of the variational problems:
S
Find
ê EX
; uS p
h
h
∈ U; 0 such that
EQ;p
B ê EX
; u p ; v = R; u p (v)
S
S
h
∀v ∈ U; 0
(74)
h
S p+k
Similarly, the computed element error indicator functions ê ;huS p , are the solutions of the discrete
h
problems:
Find ê
Shp+k
; uS p
h
∈ Sh;p+k
0 () such that
S p+k
B ê ;huS p ; v = R;EQ;p
u p (v)
S
h
h
∀v ∈ Sh;p+k
0 ()
(75)
Subtracting (75) from (74) for v ∈ Sh;p+k
0 (), and using the notation
S p+k
def
e S p+k = ê EX
; u
eˆ;uh
Copyright ? 2000 John Wiley & Sons, Ltd.
S
p
h
p
h
− ê ;huS p
(76)
h
S
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
439
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
S p+k
for the exact error in ê ;huS p , we get
h
B e S p+k ; v = REQ;p
eˆ;uh
p+k
S
;uS p ;eˆ h
h
p
h
S
S p+k
def
h
(v) = R;EQ;p
u p (v) − B (ê; uS p ; v)
S
∀v ∈ U; 0
h
h
(77)
and the orthogonality condition
B e S p+k ; v = REQ;p
eˆ;uh
S
p+k
;uS p ;eˆ h
p
S
h
(v) = 0
∀v ∈ Sh;p+k
0 ()
(78)
h
Further, for each element with @ ∩ D = ∅, from (44) and (49), we have
Z
Z p+k
Z
S
EX
e S p+k = ê ; u p − ê ;huS p = 0
eˆ;uh
S
p
h
S
h
(79)
h
S p+k
It follows that e S p+k , the error in ê ;huS p , is the solution of the variational problem:
eˆ;uh
h
p
h
S
p+k
Find e S p+k ∈ V
such that
; 0
eˆ;uh
p
h
S
B e S p+k ; v = REQ;p
p+k
S
;uS p ;eˆ h
h
eˆ;uh p
S
h
where
p+k def
=
V
; 0
p+k
∀v ∈ V
; 0
(v)
(80)
Z
p+k
v ∈ U B (v; w) = 0 ∀w ∈ S0 ();
v = 0 if @ ∩
D
=∅
(81)
From the duality, we get
REQ;p S p+k (v)
||e S p+k ||U =
eˆ;uh
p
S
h
sup
p+k
v∈V
; 0
;uS p ;eˆ h
h
||v||U
= ||REQ;p
S
p+k
;uS p ;eˆ h
||L(U; 0 ;R)
(82)
h
S p+k
Further, integrating by parts the term B (ê;huS p ; v) in (77) and using (37) we get
h
Z REQ;p
p+k
S
;uS p ;eˆ h
h
(v) =
r; uS p + ê
h
Shp+k
; uS p
h
v+
P
⊂ @
Z J; u
S
p
h
−
@ Shp+k
ê; uS p
@n
h
v
(83)
From the triangle inequality in (83) we have
Z Z @ Shp+k
Shp+k
P EQ;p
r; uS p + ê; uS p v+
J; u p −
ê; uS p v:
|R
p+k (v)|6
S
S
h
@n
h
h
h
h
;uS p ;eˆ
⊂ @
h
and employing the Cauchy–Schwarz inequality we get
P @ Shp+k S p+k ê; uS p ||v||L2 ()
|REQ;p S p+k (v)|6r; uS p + ê h 2 ||v||L2 () +
J; uS p −
@n
L ()
L2 ()
h
h
h
;uS p ;eˆ h
⊂ @
(84)
h
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
440
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
and from (82) we get
||e S p+k ||U 6 sup
eˆ;uh
p+k
v∈V
; 0
p
h
S
S p+k r; uS p + ê h h
P ||v||L2 ()
||v||L2 ()
@ Shp+k +
−
ê
J
; uS p p
;
u
S
@n
L2 () ||v||U
L2 () ||v||U
h
h
⊂ @
(85)
and letting
C1 (p + k) = sup
p+k
v∈V
; 0
||v||L2 ()
;
||v||U
C2 (p + k) = sup
p+k
v∈V
; 0
||v||L2 ()
||v||U
(86)
and we obtain
S p+k ||e S p+k ||U 6C1 (p + k)r; uS p + ê;huS p eˆ;uh
h
h
p
h
S
L2 ()
+
P
⊂ @
C2 (p + k)J; u
p
S
h
−
@ Shp+k ê; uS p @n
L2 ()
h
(87)
The constants C1 , and C2 , are given in terms of the corresponding constants Ĉ1 , and Ĉ2 , in the
master element ,
ˆ as shown below. Employing that
q
||v||L2 () 6 max |J | || v̂||L2 ( ˆ) ;
q
||v||L2 () 6 max |J | || v̂||L2 ()ˆ
(88)
where v ∈ U , J = @x=@x̂, is the Jacobian matrix and |J |, is the determinant of the Jacobian
matrix, and
p
||v||U ¿ K || v̂||U ˆ
(89)
K = min(|J |min (J−T J−1 ))
(90)
where
and (88) and (89) are proven in Appendix II; we get
s
C1 (q)6
max |J |
Ĉ1 (q);
K
Copyright ? 2000 John Wiley & Sons, Ltd.
s
C2 (q)6
max |J |
Ĉ2 (q)
K
(91)
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
441
where
Ĉ1 (q) = sup
q
v∈V
ˆ ˆ
||v̂||L2 ( ˆ)
;
||v̂||U ˆ
Ĉ2 (q) = sup
q
v∈V
ˆ ˆ
||v̂||L2 ()ˆ
||v̂||U ˆ
(92)
Let us now give some details about how the constants Ĉ1 (q), and Ĉ2 (q) are computed.
Computation of Ĉ1 and Ĉ2 : The computation of Ĉ1 (q), and Ĉ2 (q) will be done by noting that
p
p
Ĉ2 (q) = lim Ĉ2 (q)
Ĉ1 (q) = lim Ĉ1 (q);
p→∞
p→∞
(93)
where
p
Ĉ1 (q) = sup
v̂∈Vq=ˆ p
and
Vˆ q=p =
def
|| v̂||L2 ( ˆ)
;
|| v̂||U ˆ
p
Ĉ2 (q) = sup
v̂∈Vq=ˆ p
|| v̂||L2 ()ˆ
|| v̂||U ˆ
Z
p q
v ∈ Ŝ B ˆ (v; w) = 0 ∀w ∈ Ŝ ;
v=0
(94)
(95)
ˆ
p
p
In order to compute Ĉ1 (q), and Ĉ2 (q) we need to construct a basis for Vˆ q=p. Let i ; i = 1; 2; : : : ;
(p + 1)2 , be the bi-p element shape functions and let
Z
i+1
ˆ
˜
Z
1 ; i = 1; 2; : : : ; (p + 1)2 − 1
(96)
i = i+1 −
1
ˆ
We then have
Z
ˆ
˜i = 0
(97)
Group the above functions ˜i , into lower order
Li = ˜i ;
i = 1; 2; : : : ; (q + 1)2 − 1
(98)
j = 1; 2; : : : ; (p + 1)2 − (q + 1)2
(99)
and higher order
˜
H
j = i+(q+1)2 −1 ;
and compute the matrices D and E with entries
L
Dij = − B ˆ (H
i ; j );
Ejk = B ˆ (Lj ; Lk )
(100)
def
for i = 1; 2; : : : ; N q=p = (p + 1)2 − (q + 1)2 , and j; k = 1; 2; : : : ; (q + 1)2 − 1, and determine F = E−1 .
Now let
!
2
2
(q+1)
P (q+1)
P
H
˜
Dij Fjk Lk ; i = 1; 2; : : : ; N q=p;
(101)
i = i +
k=1
Copyright ? 2000 John Wiley & Sons, Ltd.
j=1
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
442
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
and note that the ˜ i , are orthogonal to any Lm , for m = 1; 2; : : : ; (q + 1)2 − 1, because
!
2
2
(q+1)
P (q+1)
P
L
H
L
˜
Dij Fjk Bˆ (L ; Lm )
Bˆ (i ; m ) = Bˆ (i ; m ) +
j=1
k=1
= −Dim +
= −Dim +
2
(q+1)
P
2
(q+1)
P
k=1
j=1
2
(q+1)
P
j=1
Dij Fjk
2
(q+1)
P
Di;j
k=1
|
k
!
Ekm
!
Fjk Ekm
{z
=0
}
jm
R
where jm is the Kronecker delta. Further, we also have ˆ ˜ i = 0, from (97). Hence ˜ i ; i = 1; 2; : : : ;
N q=p; form the basis of Vˆq=p. Therefore, any v ∈ Vˆq=p, can be expressed in the form
v̂ =
q= p
NP
i=1
ai ˜ i
(102)
and we have
||v̂||2L2 ()
ˆ
=a A
T
q= p
p
Aq=
ij =
a;
and
||v̂||2Uˆ = aT Bq=pa;
Bijq=p =
Z
ˆ
ˆ
˜i ˜j
(103)
@˜i @˜j
@˜i @˜j
+
@x̂1 @x̂1
@x̂2 @x̂2
Similarly, we have
||v̂||2L2 ()ˆ = aT Cq=pa;
Z
Cijq=p =
(104)
Z
ˆ
˜i ˜j
(105)
It follows that
p
aT Aq=pa
(Cˆ1 )2 = max T q=p ;
q= p a B
a
a∈RN
and hence
p
Cˆ1 =
r
max kq=p;
k=1; :::; N q=p
p
aT Cq=pa
(Cˆ2 )2 = max T q=p
q= p a B
a
a∈RN
p
Cˆ2 =
r
max
k=1; :::; N q=p
kq=p
(106)
(107)
where kq=p and kq=p are the eigenvalues of the generalized eigenvalue problems
Aq=py = kq=pBq=py;
Cq=py = kq=pBq=py
(108)
p
p
We computed Cˆ1 ; and Cˆ2 using the Harwell Library routine EA12AD (see [15]). For each q,
p
p
we employed a suciently high p to obtain convergence of Cˆ1 (q); and Cˆ1 (q) to 12 signicant
digits. Table I shows the converged values of Cˆ1 (q) and Cˆ2 (q) for q = 0; 1; 2; : : : ; 8, and the values
of p for which the convergence was achieved.
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
443
Table I. Guaranteed upper bound for the energy norm of the error: The
constants Cˆ1 (q) and Cˆ2 (q) for q = 1; 2; : : : ; 8, and the corresponding
values of p which was used to compute them
q
Cˆ1 (q)
p
Cˆ2 (q)
p
0
1
2
3
4
5
6
7
8
0·636619772367
0·630301952842
0·609347491519
0·568383672094
0·530296680623
0·500349587345
0·464059359874
0·430324985456
0·400678456146
8
9
12
14
16
18
20
22
24
0·816496580927
0·810315743719
0·781457077283
0·756306743339
0·690556011657
0·621895745689
0·558985734593
0·498567556786
0·426746908655
11
14
15
17
19
21
24
26
28
ˆ
ˆ
Note that the values
R of C1 (0) and C2 (0) can also be computed analytically. Any function v̂ ∈ Uˆ
with zero average ˆ v̂ = 0 can be expanded in a cosine series in the form
∞
P
v̂ =
m; n=0
m+n6=0
am; n cos
m
n
(1 + x̂1 ) cos
(1 + x̂2 )
2
2
(109)
We then have
||v̂||2Uˆ =
2
4
∞
P
m; n=0
m+n6=0
(m2 + n2 )a2m; n ;
||v̂||2L2 ()
ˆ =
∞
P
m; n=0
m+n6=0
a2m; n
(110)
and hence
∞
P
2
Cˆ1 (0) = sup
0
v∈
ˆ V̂ˆ
||v̂||2L2 ()
ˆ
||v̂||2Uˆ
=
4
sup
2 am; n
m; n=0
m+n6=0
∞
P
m; n=0
m+n6=0
a2m; n
(111)
(m2 + n2 )a2m; n
The above supremum is achieved for the cases m = 0 and n = 1, or m = 1 and n = 0, and we get
2
Cˆ1 (0) = ≈ 0·636619772367
(112)
which agrees with the value given in Table I. Further, in order to compute Cˆ2 (0), let us consider
the edge ˆ1 which joins the vertices (−1; −1) and (1; −1) in the master element. We have
∞
2 ∞ ∞ 2
∞
P
P
P P 1
2
nam; n
am; n =
||v̂||L2 (ˆ1 ) =
m=1
n=1
m=1
n=1 n
∞
∞
∞
∞ P
∞
P
P 1
P 2 2
2 P
=
6
n
a
n2 a2m; n
(113)
m; n
2
6 m=1 n=1
m=1
n=1 n
n=1
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
444
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
∞
P
where we used Cauchy–Schwarz inequality, and the identity
1=n2 = 2 =6. Therefore, we have
n=1
∞
P
||v̂||2L2 (ˆ1 )
2
Cˆ2 (0) = sup
||v̂||2Uˆ
0
v∈
ˆ V̂ˆ
m; n=0
m+n6=0
∞
P
2
2
6 sup
3 am; n
m; n=0
m+n6=0
n2 a2m; n
(114)
(m + n2 )a2m; n
The above supremum is achieved for the case m = 0, and for all n and we get
r
2
ˆ
≈ 0·816496580927
C2 (0) =
3
which agrees with the value given in Table I.
Above we have constructed the guaranteed upper estimate
v
2
2
u
P
U; Shp+k def u
S p+k
EX
EX
EEXPL
||eS p ||U 6ETh ; u p 6ETh ; u p = t EThh ; u p +
p+k
S
S
h
S
h
S
h
∈Th
h
; eˆ;uh
(115)
(116)
p
S
h
This estimate is valid for the special case that h ≡ Th
, and that the solutions of the element
S p+k
residual problems ê ;EXu p , are approximated by the p method, namely ê ;huS p .
S
h
h
Let us now consider the more general case where the element residual problems are solved by
S p+kn
, is computed. As before, we have
the h method, namely the case for which the estimate ETh=2
h; u p
S
s
S p+kn
ETh=2
h; u
=
ETEX
h ;u p
S
2
p
S
h
h
+
P
∈Th
h
S p+kn
||ê EX
; u
− ê ;h=2uS p ||2U
p
S
h
(117)
h
and
||ê EX
; u
S
P
S p+kn
p
h
− ê ;h=2uS p ||2U =
h
0 ∈Th=2
n
||ê EX
; u
S
S p+kn
p
h
− ê ;h=2uS p ||2U0
(118)
h
where Th=2
n is the mesh obtained from n nested subdivisions of the element .
S p+kn
h=2
We will now construct an upper estimate for ||ê EX
; u p − ê; uS p ||U , exactly as we constructed the
S
h
h
upper estimate for ||uEX − uS hp ||U in (116) above. We have
v
u p+k+‘ 2
p+k
p+k+‘
P
S
Sh=2
u S n
n
n
h=2
EX
EX
2
h=2
+
||
ê
−
ê
||ê EX
p+k
0
p+k = t E
0
; u p − ê ; uS p ||U 6E
S n
;eˆ
;eˆ ||U0
S
S
h
h
n
; eˆ h=2
Th=2
n ; eˆ h=2
Th=2
n 0 ∈Th=2
n
(119)
Using the explicit upper estimate
S p+k+‘
n
h=2
EXPL
||êEX
0 ; eˆ − ê 0
p+k+‘
; eˆ ||U0 6E
S
n
(120)
0 ; eˆ;uh=2p
S
h
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
S p+k+‘
n
S p+kn
h=2
where we have used the subscript ê in êEX
0 ; eˆ , and ê0 ; eˆ
(119) we get
S p+kn
h=2
EX
||êEX
; u p − ê; uS p ||U 6E
S
h
h
p+k
h=2n
S
; eˆ
Th=2
n 445
indicate ê h=2 , and employing (120) into
6E
p+k+‘
U; Sh=2
n
S
(121)
p+k
n
; eˆ h=2
Th=2
n where
E
p+k+‘
U; Sh=2
n
S
p+k
n
; eˆ h=2
Th=2
n v
u p+k+‘ 2
P EXPL 2
S
E S p+k+‘
=u
t E h=2n S p+kn +
def
; eˆ h=2
Th=2
n n
0 ∈Th=2
n
(122)
0 ; eˆ;uh=2p
S
h
The guaranteed upper bound for the case that h 6= Th
follows similar steps as depicted in
Figure 3.
We have thus constructed the guaranteed upper estimator for ||e SEX
p ||U .
h
Theorem 3.3 (Guaranteed upper bound for ||eSEX
p ||U ). Consider the above setting and the notah
tions given above. We have
p+k p+k+‘
U; Sh=2
n =Sh=2n
def
EX
p ||U 6E
6E
=
||eSEX
p
h ; uS
h ; uS p
h
h
h
v
u p+k 2
2
P U; S p+k+‘
u Sh=2n
h=2n p+k
t E
; u p +
E
S n
h
!
S
Th=2j n ;eˆ!h=2
j
!j ∈
h
h
(123)
Let us now give an example which illustrates the above constructions.
Example 1 (The upper bound for the energy norm of the exact error). We considered the
model problem and the mesh Th
shown in Figure 1 and biquadratic elements p = 2. We computed
the energy norm of the exact error ||eSEX
2 ||U , by employing an overkill mesh Th;ovk , shown in
h
Figure 4 which is constructed such that the relative error ||eSEX
2 ||U =||uS 2 ||U , is 0·01 per cent accurate.
h
h
This mesh was obtained by rening the original mesh Th
, uniformly three times then by rening
the elements adjacent to the re-entrant corners ve more times.
The overkill error was computed using the GFEM (see [16]) space of degree p = 8 and by
including the rst term of the asymptotic expansion of the solution of the Laplacian at each
reentrant corner.
, by solving the element residual problems using
We also computed the exact upper bound ETEX
h
the rened meshes shown in Figure 5 and p = 8; these meshes are the restriction of the mesh
shown in Figure 4 which was used in the overkill.
We have
||eSEX
2 ||U
h
||uS 2 ||U
= 29·17 per cent6
h
ETEX
h
= 48·22 per cent
||uS 2 ||U
h
Next, we obtained the element error indicators by employing the p and the h methods to solve
the element residual problems.
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
446
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
Figure 3. Guaranteed upper bound for the energy norm of the error: Schematic depiction of the steps involved in the
construction of the guaranteed upper bound for ||eEX
p ||U : (a) The mesh Th which is employed in the computation of the
S
h
!
nite element solution uS p ; (b) The partition of Th
into the patches !j and the meshes Th=2j n in the computation of the
h
S
p+k
n
error indicator functions ê !h=2
j; u
S
p
h
!
!
; (c) The mesh Th=21 n employed in the patch !1 ; (d) The elements 0 ∈ Th=21 n which are
used for the solution of the element residual problems
S p+k
(a) Element residual problems solved using the p method: We computed EThh , for p + k = 3;
S p+k
4; : : : ; 8. Figure 6 shows the graphs of the eectivity indices of the estimators EThh =||eSEX
2 ||U ,
S p+k
S p+k
h
S p+k
h
||U¡EThh ¡||eSEX
and ||eSh2 ||U =||eSEX
2 ||U , for p + k = 3; 4; : : : ; 8. Note that, we have ||e 2
2 ||U .
S
h
h
h
S3
h
n
(b) Element residual problems solved using the h method: We computed EThh=2 , for n = 0; 1; 2; 3.
S3
S3
n
n
h=2
||U =
Figure 7 shows the graphs of the eectivity indices EThh=2 =||eSEX
2 ||U , and the ratios ||e 2
S
h
S3
h
n
h=2
||eSEX
||U ¡
2 ||U , for n = 0; 1; 2; 3. Note that as expected from Theorem 3.1, we have, ||e 2
S
S3
h
h
n
EThh=2 ¡||eSEX
2 ||U .
h
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
447
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
Figure 4. The overkill mesh Th;ovk
used for the computation of ||eEX
p ||U : This mesh was obtained by rening the original
S
h
mesh Th
uniformly three times, then by rening the elements adjacent to the re-entrant corners ve more times
Figure 5. The overkill meshes employed for the computation of the exact element error indicators EX
!j : These meshes are
the restrictions of the overkill mesh of Figure 4 in the elements of the original coarse mesh
Let us also construct the computable guaranteed upper bounds by employing the solutions of
the element residual problems obtained by the p and the h methods.
(a) Element residual problems solved using the p method: Figure 8 compares the graphs of the
U; Shp+k
eectivity indices ETh
as expected we have
EX
EX
=||eSEX
2 ||U , and ETh =||eS 2 ||U , for p + k = 3; 4; : : : ; 8. We note that
h
h
U; Shp+k
EX
||eSEX
2 ||U ¡ET ¡ET
h
h
h
(b) Element residual problems solved using the h method: In this case we employed k = 1
and ‘ = 1 and computed the upper bound E
p+k p+k+‘
U; Sh=2
n =Sh=2n
Figure 9 compares the graphs of the eectivity indices
ETEX
=||eSEX
2 ||U ,
h
h
3
U; Sh=2
n
, which we will denote by ETh
.
EX
||eSEX
2 ||U =||eS 2 ||U
h
h
and
U; S 3 n
ETh h=2 =||eSEX
2 ||U ,
for n = 0; 1; 2; 3 and as expected, we have
h
3
U; Sh=2
n
EX
||eSEX
2 ||U ¡ET ¡ET
h
h
h
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
448
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
Figure 6. The upper bound for energy norm of the error: Approximate solution of the element residual problems using
S
p+k
h
the p method. Graphs of the eectivity index of the exact estimator ETEX=||eEX
2 ||U , the computed estimator ET
h
S
h
h
p+k
h
S2
h
and the ratio of the energy norm of the truth-mesh error to the energy norm of the exact error ||e
S
p+k
p + k = 3; 4; : : : ; 8. Note that ET h
h
S
p+k
S
=||eEX
2 ||U ,
||U =||eEX
|| ,
S2 U
S
h
versus
h
is an upper estimator for ||e h2 ||U , but it underestimates the energy norm of the exact
S
h
error ||eEX
2 ||U
S
h
S p+kn
p+k
U; Sh=2
n
Let us also mention that the bounds E
EX
, E
h=2
, and E
h
h
h
S p+kn
E
h=2
(Xint ),
h
depend on the particular choice of
U; S p+kn
the edge corrections, i.e. we have E
EX
= E
EX
(Xint ),
and E
h h=2 (Xint ) where Xint is the
h
h
vector of edge-corrections for the interior edges.
In Appendix I we give the Ladeveze’s method for constructing the edge corrections which are
unique up to a set of free constants, one for each interior vertex of the mesh. The constants can
be chosen equal to zero and we call this choice as the Bank–Weiser equilibration [10] or the
constants can be chosen such that they minimize a weighted norm at each vertex (see [4–6]) and
we call this choice as Ladeveze equilibration. Another factor which aects these bounds is the
order of equilibration. As we have seen in Appendix I, we can construct the edge corrections Xint ,
such that the orthogonality condition
p
∀v ∈ Sh;0
()
EQ;p
R;u
p (v) = 0
S
h
(124)
holds for each element. However, a weaker condition is sucient for the existence of the estimate,
namely the consistency condition
EQ;0
R;u
p (1) = 0
S
(125)
h
In general, we can employ a qth-order equilibration for which
EQ; q
R;u
p (v) = 0
S
h
Copyright ? 2000 John Wiley & Sons, Ltd.
q
∀v ∈ Sh;0
()
(126)
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
449
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
Figure 7. The upper bound for energy norm of the error: Approximate solution of the element residual problems using
S3 n
h=2
the h method. Graphs of the eectivity index of the exact estimator ETEX=||eEX
=||eEX
2 ||U , the computed estimator ET
2 ||U ,
h
S
3
Sh=2
n
and the ratio of the energy norm of the truth-mesh error to the energy norm of the exact error ||e
p+k
h
S2
h
n = 0; 1; 2; 3; 4. Note that we have ||e
S
S
h
h
S2
p+k
S
||U ¡ET h
h
h
h
||U =||eEX
|| ,
S2 U
for
h
¡||eEX
2 ||U
S
h
with 06q6p. We will express the dependence of the bounds on the order and the type of
BW(q)
Lad(q)
(resp. E
EX;
) is the
equilibration by introducing additional indices, for example E
EX;
h
h
exact upper estimator based on Bank–Weiser (resp. Ladeveze) equilibration of order q.
Example 2 (The inuence of the choice of equilibration on the upper bound). Figure 10(a)
S p+k ; EQ(q)
EQ(q)
h
(resp. Figure 10(b)) shows the graphs of ETEX;
=||eSEX
=||eSEX
resp.
2 ||U , and ETh
2 ||U ,
h
h
h
U; Shp+k ;EQ(q)
=||eSEX
for EQ = BW, and EQ = Lad, and for q = 1; 2. From these results it can
ETh
2 ||U
h
EQ(q)
be seen that for the example problem considered here ETEX;
is practically independent of
h
S p+k ; EQ(q)
U; Shp+k ; EQ(q)
, and the upper bounds ETh
the choice of equilibration and EThh
signicantly by the choice of the equilibration.
Let us now summarize the conclusions in this section:
are not inuenced
1. We have the guaranteed upper estimator
EX
p ||U 6E
||eSEX
h =
r P
h
!j ⊂ h
2
||êEX
;u p ||U!
S
h
j
(127)
given by
2. However, the computable version of E
EX
h
S p+kn
E
h=2
h
Copyright ? 2000 John Wiley & Sons, Ltd.
s
=
P
!j ⊂ h
S p+kn
||ê;uh=2S p ||2U!
h
j
(128)
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
450
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
Figure 8. The upper bound for energy norm of the error: Approximate solutions of the element residual problems using
p+k
h
U; S
the p method. Graphs of the eectivity indices ET
S
p+k
U; S
ET h
h
→
p+k
h
U; S
=||eEX
2 ||U for p + k = 3; 4; : : : ; 8. Note that ET
h
ETEX ,
h
h
h
¿ETEX and
h
as p + k → ∞
Figure 9. The upper bound for energy norm of the error: Approximate solutions of the element residual problems
3
U; Sh=2
n
using the h method. Graphs of the eectivity indices ET
3
U; Sh=2
n
ET
h
→
3
U; Sh=2
n
=||eEX
2 ||U for n = 0; 1; 2; 3. Note that ET
S
h
ETEX ,
h
as
h
h=2n
h
¿ETEX and
h
→0
is not a guaranteed upper bound for the energy norm of the exact error ||eSEX
p ||U and we
h
could have
S p+kn
EX
p ||U 6E
¡||eSEX
E
h=2
h
h
h
Copyright ? 2000 John Wiley & Sons, Ltd.
(129)
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
451
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
Figure 10. The inuence of the choice of equilibration on the upper bound: Eect of the type and the order of
EX; EQ(q)
the equilibration. (a) Graphs of ET
h
p+k
U; S
ET h
h
; EQ(q)
=||eEX
||
S2 U
S
p+k
h
=||eEX
2 ||U , and ET
S
h
h
; EQ(q)
EX; EQ(q)
=||eEX
2 ||U ; (b) graphs of ET
S
h
h
=||eEX
2 ||U and
S
h
for EQ = BW, and EQ = Lad, and for q = 1; 2. Note that for the considered example, the eect
h
of the order and the type of equilibration is practically negligible
3. Employing the identity
v
u p+k 2
p+k
u S n
P
Sh=2
n
EX
EX
2
+
||
ê
−
ê
E
h = t E
h=2
;uS p ||U!
;u p
h
!j ⊂ h
Copyright ? 2000 John Wiley & Sons, Ltd.
S
h
h
j
(130)
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
452
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
and constructing a computable upper estimate for
S p+kn
p+k+‘
U; Sh=2
n
h
Th=2j n ;eˆ!h=2
j
h=2
||êEX
;u p − ê;uS p ||U!j 6E
S
h
S
!
(131)
p+k
n
we get the guaranteed computable upper estimator
v
2
u p+k 2
p+k
p+k+‘
u Sh=2n
P
U;
S
U; Sh=2
n
h=2n
EX
t
p ||U 6E
||eSEX
E
E
6E
=
+
p+k
h
h
h
S
h
!j ⊂ h
!
n
(132)
Th=2j n ;eˆ!h=2
j
4. There is no signicant dierence between the values of the estimators for various types and
orders of equilibration.
4. THE LOWER BOUND FOR THE ENERGY NORM OF THE EXACT ERROR
The residual error indicator functions êEX
!j ; u p , can be put together into a function dened over ,
S
which we will denote by
êEX
uS p
h
∈ L (
), namely
2
h
def
EX
êEX
u p |!j = ê!j ; u p ;
S
S
h
!j ∈ h
(133)
h
Note that êEX
u p is smooth in the interior of elements but it may have jump discontinuities at the
S
h
vertices and across the interior edges of the mesh. We will denote these discontinuities by
eˆEX
u p
G
S
h
def
EX
EX
= <êEX
u p = = êu p (x + 0n ) − êu p (x − 0n )
S
S
h
h
S
(134)
h
eˆEX
u p
S
and we will call G h , the gap of êEX
u p on the edge . Here , is an interior edge and n , is the unit
S
h
eˆEX
u p
S
normal assigned to . For the edges on the boundary ⊂ @
, we let G h ≡ 0. Let ETh denote the
set of all edges in the mesh and let us dene the broken energy space with prescribed gaps on
all the edges, namely
(
def
UBroken
=
Ex
G eˆu p
eˆEX
u p
S
2
v ∈ L (
) v| ∈ U;0 ; <v= = − G h ; ∈ ETh
)
(135)
T s
h h
Broken
EX
eu p
and we have êEX
u p ∈ U eˆ Ex . Note also that the exact error eS p does not have gaps, namely G S ≡
EX
S
h
0, for all ∈ ETh and
G
h
p
Tus
h h
p ∈ U
; 0
eSEX
h
Copyright ? 2000 John Wiley & Sons, Ltd.
h
(136)
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
453
The dierence between the exact error and the exact error indicator
def
EX
EX
p
eEX
ˆEX = eS − êu p
u p
S
h
S
h
(137)
h
∈ UBroken
and eEX
satises the following orthogonality condition:
Then we have eEX
eˆ Ex
ˆEX
ˆEX
G
u p
u p
p
Tus
h h
S
h
S
h
Lemma 4.1. Consider the above setting and the notations. Then
P
∈Th
B eEX
;
v
=0
EX
ˆ
∀v ∈ U
; 0
u p
(138)
S
h
Proof. For any v ∈ U
; 0 , we have
P
P
p ;v
B eEX
; v = B
eSEX
−
B!j êEX
!j ; u p ; v
ˆEX
u p
∈Th
h
S
S
!j ∈
h
h
= R
; uS p (v) −
h
P
!j ∈
h
h
R!EQ;p
j ; uS p (v) = 0
(139)
h
From the above orthogonality it follows that eEX
has the minimal energy norm as shown in
ˆEX
u p
S
h
the following result.
Lemma 4.2. Consider the above setting and notations. Then we have
rP
∈Th
rP
||eEX
||2 ¡
ˆEX U
u p
S
h
∈Th
||eˆEX
||2U
u
∀eˆEX
∈ UBroken
Ex
u
G eˆ
p
h
p
h
S
S
p
Tus
h h
(140)
Proof. Note that
∈ U
; 0
= eEX
ˆEX − eˆEX
u
u p
S
h
(141)
p
h
S
and hence
P
∈ Th
||eˆEX
||2U =
u
p
S
h
P
∈ Th
||eEX
− ||2U =
ˆEX
u p
S
h
P
∈ Th
where we have used the orthogonality condition (4:6) that
and omitting ||||2U , (140) follows.
||eEX
||2 + ||||2U
ˆEX U
u p
(142)
S
h
P
∈ Th
B eEX
;
v
= 0, for all v ∈ U
; 0 ,
ˆEX
u p
S
h
Lemma 4.3. Consider the above setting and the notations. Then
p ||U =
||eSEX
h
Copyright ? 2000 John Wiley & Sons, Ltd.
r
P
EX )2 −
(E
||eEX
||2
ˆEX U
h
∈ Th
u p
(143)
S
h
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
454
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
EX
Proof. We have eEX
= eSEX
p − êu p , and hence
ˆEX
u p
S
h
S
h
P
∈ Th
h
2
||eEX
||2 = ||eSEX
p ||
U−2
ˆEX U
u p
h
S
h
P
∈Th
P
EX
2
+
B êEX
||êEX
u p ; eS p
u p ||U
S
h
h
S
∈Th
(144)
h
and noting that
P
P
EX
EX
EX
EX
EX 2
=
B êEX
;
e
REQ;p (eSEX
p
p ) = R
(e p ) = B
(e p ; e p ) = ||e p ||
u p
U
S
S
S
S
S
S
∈Th
h
h
h
∈ Th
h
h
h
(145)
h
we get (143).
Combining (140) and (143) we get the following guaranteed lower bound for ||eSEX
p ||U .
h
Theorem 4.1. Consider the above setting and the notations. Then
r
P
def
L
EX
p ||U
EX ) =
(
(E
)2 −
||eˆEX
||2U 6 ||eSEX
∀eˆEX
∈ UBroken
E
Ex
p
e
ˆ
;
u
h
h; u p
u
u
u
G eˆ
S
h
S
p
h
S
h
p
h
∈Th
h
S
p
Tus
h h
p
h
S
(146)
and
L
L
eˆEX
6 E
E
h
h; u
u
p
S
h
p
S
h
EX
p
eEX
ˆEX = ||eS ||U
u p
(147)
h
S
h
Let us now give one particular method of construction of the gap function eˆEX
. Let XTh denote
u
p
h
S
the set of vertices in the mesh. For X ∈ XTh , note that
êEX
uS p (X )
is well dened because
h
EX
1+
(!j )
êEX
u p |!j = ê!j ; u p ∈ H
S
where ¿0, and êEX
!j ; u
S
p
h
S
h
(148)
h
is continuous at the vertices of !j . Figure 11(b) depicts the error indicator
functions along the common edges of the shaded elements shown in Figure 11(a) and illustrates
eˆEX
u p
the gap at the vertex node X . Therefore, we dene the vertex gap GX; , at the vertex X of ,
S
h
def
= êEX
GX; eˆEX
u p | (X ) −
u
p
S
h
S
h
nX
1 P
êEX | X (X );
nX j=1 uS hp j
X ∈ XTh
(149)
X
is the set of elements which are connected to X . Let e(1)
be given by piecewise
where {Xj }nj=1
ˆEX
u p
S
h
mapped bilinear functions
|
e(1)
ˆEX
u p
S
h
=−
4
P
k=1
eˆEX
u p
GX ; ’Xk
S
h
k
(150)
where ’Xk ∈ S01 (), is the mapped bilinear basis function corresponding to the vertex Xk . Figure
, along the common edges of the shaded elements shown in
11(c) depicts the gap function e(1)
ˆEX
u p
S
h
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
455
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
Figure 11. The lower bound for energy norm of the error: Construction of the gap functions for the lower bound.
(a) The mesh indicating the four shaded elements for which the error indicator functions and the gap functions are
(1)
shown; (b) the error indicator functions along the edges AX and BX ; (c) the gap function EX , along the common edges
eˆu
(2)
of the shaded elements; (d) the gap function eˆEX
u p;
S
h
1,
2
p
h
S
along the common edges of the shaded elements
(1)
Figure 11(a). The function êEX
u p + eˆEX is continuous at the vertices but has gaps across the edges.
S
u p
h
S
h
Let
(
|
e(2)
ˆEX
u p
S
h
def
∈ U; <eˆEX +(1) = =
u
eˆEX
p
h
S
u p
(1)
v ∈ U;0 v| = − 12 <êEX
uS p + eˆEX = ; ⊂ @
h
S
u p
)
(151)
S
h
h
then we have
= e(1)
+ e(2)
∈ UBroken
e(1)+(2)
Ex
ˆEX
ˆEX
ˆEX
G eˆ
u p
u p
u p
h
h
h
S
Copyright ? 2000 John Wiley & Sons, Ltd.
S
S
p
Tus
h h
(152)
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
456
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
L
and we can now compute the lower bound E
h; u
S
L
E
h; u
S
e(1)+(2)
=
ˆEX ;OPT
p
h
u p
S
h
(2)
eˆEX
u p
S
h
e(1)+(2)
. We can also optimize this bound by
ˆEX
p
h
u p
S
h
L
max
(1)+(2)
E
eˆEX
h ; uS p
u p
h
∈ U EX
S
(1)
h
; <eˆu +
=
eˆEX
p
h
S
(153)
u p
S
h
1=2
(),
This maximum can also be expressed in terms of a set of unknown edge-functions ∈ H0;0
(2)
∗
by letting eˆ EX ∈ U; <eˆEX +(1) = , be such that for every = @ ∩ @ , we have
u
u ;E
eˆ EX
p
h
p
h
S
S
u p
S
h
(1)
| = − <êEX
(2)
u + eˆEX = ;
eˆEX ; p
h
u p
S
(1)
(2)
| ∗ = − (1 − )<êEX
u + eˆEX =
eˆEX ; u p
S
p
h
u p
S
h
u p
S
S
h
(154)
S
h
h
and
; v = − B (1)
; v ∀v ∈ U;0
B (2)
eˆEX ;
eˆEX
u p
u p
h
h
S
(155)
S
and we have
L
E
h; u
p
S
h
L
(1)+(2)
= max E
h; u
eˆEX ;OPT
u p
p
S
h
∈ C 0 ()
∈E
Th
S
h
+ (2)
(1)
eˆEX
eˆEX ;
u p
u p
h
h
S
(156)
S
(1)
This optimum value is, in general, less than ||eSEX
p ||U . Note that we could have also dened EX
eˆ
u p
h
S
h
in terms of the unknown weights at the vertices of the mesh and then the optimal value would
(2)
, along the common edges
coincide with ||eSEX
p ||U . Figure 11(d) depicts the gap functions EX
eˆ ;1=2
u p
h
S
h
of the shaded elements shown in Figure 11(a).
Remark 4.1. The solution of the Dirichlet problem (155) exists and is unique. This is veried by
(1)
(2)
1=2
noting that <êEX
u +eˆEX = ∈ H0;0 (); ∈ @ and for ∈ Th . Therefore, eˆEX ;1=2 can be obtained by the
p
h
S
u p
u p
h
h
S
S
superposition of the solutions of four Dirichlet boundary value problems for which homogeneous
(1)
boundary conditions are imposed on three edges and <êEX
u +eˆEX = are imposed on the fourth edge.
p
h
S
u p
S
h
Let us note that the lower bound
L
((1)+(2)
)
E
eˆEX
h ; uS p
u p
h
ê EX
!j ; u
S
p
h
is based on the exact error indicator functions
S
h
and it is not computable. Let us now consider the case that the computed error indicator
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
457
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
S p+kn
Broken
functions ê !h=2
j ; uS p are employed. We also note that for S p+k ∈ U S p+k , we have
eˆ;uh
h
h
eˆ;u
p
p
h
S
S
GT
h
h
S p+kn
|!j = ê EX
S p+k |!j − eˆEX
!j ; u
u
eˆ;uh
S
p
h
S
p
h
S
− ê !h=2
j ; uS p
p
h
(157)
h
and hence
|| S p+k ||2U!j = || S p+k ||2U!j + ||ê EX
!j ; u
eˆ;uh
eˆ;uh
p
h
S
S p+kn
S
p
h
p
h
2
− ê !h=2
j ; uS p ||U!
j
where we have used the orthogonality (3:10) that B!j ê EX
!j ; u
S
p+k
Sh=2
n ;0 (!j )
and that
(158)
h
S
p+k
S p+k |!j ∈ Sh=2
n ;0 (!j ).
eˆ;uh p
p
h
S p+kn
− ê !h=2
j ; uS p ; v = 0, for all v ∈
h
Using (158) in the result in Theorem 4.1, we get the
S
h
computable lower bound
L
E
h ; uS p
h
v
u
P
u S p+kn 2
p ||U
|| S p+k ||2U 6 ||eSEX
S p+k = t E
h=2h ; u p −
eˆ;uh
S
p
h
S
eˆ;uh
∈ Th
h
(159)
h
p
h
S
Let us now compute the lower bound for the example problem.
Example 3. (The lower bound for the energy norm of the error). Figure 12 (resp. Figure. 13)
L;Shp+k
3
L; Sh=2
n
L
EX
=||eSEX
2 ||U ) with ETh =||eS 2 ||U for p + k = 3; 4; : : : ; 8 (resp.
h
h
p+k
L;Sh=2
n
L
to indicate the lower bound E
S p+k : We
n = 0; 1; 2; 3). Here used the notation ETh
h; u p
compares ETh
=||eSEX
2 ||U (resp. ETh
h
S
h
eˆ;uh
p
h
S
note that we have
L; Shp+k
E Th
¡ ETLh ;
3
L; Sh=2
n
ETh
¡ETLh
i.e. the lower bounds obtained from the approximate solutions of the element residual problems
are smaller than ETLh . Figures. 12 and 13 show that the lower bounds can be very small relative to
the energy norm of the exact error and hence it is necessary to compute enhanced lower bounds.
Let us now summarize the conclusions in this section:
1. We have the guaranteed lower bound
L
E
h; u
S
p
h
eˆEX
u
p
h
S
def r EX 2
P
p ||U
=
||eˆEX
||2U 6 ||eSEX
E
h ; u p −
u
S
h
∈Th
p
h
S
(160)
h
∈ UBroken
is the gap function which is constructed such that êEX
∈ U
; 0 .
where eˆEX
u p + eˆEX
eˆ Ex
u
u
p
h
S
G
p
Tus
h h
Copyright ? 2000 John Wiley & Sons, Ltd.
S
h
p
h
S
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
458
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
Figure 12. The lower bound for energy norm of the error: Approximate solutions of the element residual problems by
p+k
h
L; S
the p method. Graphs of the eectivity indices of the lower bounds ET
have
L;S 3
ET h
h
h
=||eEX
2 ||U for p + k = 2; 3; : : : ; 8. Note that we
S
h
= 0 and for p + k65 the lower bound for the energy norm of the exact error is practically useless
Figure 13. The lower bound for energy norm of the error: Approximate solutions of the element residual problems by
3
L; Sh=2
n
the h method. Graphs of the eectivity indices of the lower bounds ET
h
L; S 3
ET h
h
=||eEX
2 ||U for n = 0; 1; 2; 3. Note that we have
S
h
= 0 and for n62 the lower bound for the energy norm of the exact error is practically useless
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
459
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
2. The computable version of the lower bound is given by
L
E
h; u
S
S p+k
p
h
eˆ;uh
v
u S p+k P
u
2
n
p ||U
|| S p+k ||2U 6 ||eSEX
= t E
h=2h ; u p −
S
p
h
S
eˆ;uh
∈Th
h
(161)
h
p
h
S
S p+kn
such that êuh=2
where S p+k ∈ UBroken
p + S p+k ∈ U
; 0 .
eˆ Ex
S
eˆ;uh
G
p
Tus
h h
p
h
S
eˆ;uh
h
p
h
S
3. The computed lower bound can be very small relative to ||eSEX
p ||U , in some cases and we
h
L
L
.
have E
h ; u p S p+k 6 E
h ; u p eˆEX
u
S
eˆ;uh
h
S
p
h
S
p
h
S
h
5. ITERATIVE IMPROVEMENT OF THE BOUNDS FOR THE ENERGY NORM
OF THE EXACT ERROR
In this section we give the construction of the improved upper and lower bounds. Let us rst
consider the case that the element residual problems are solved exactly. Let
def
EX
ẽEX
u p | = êu p + eˆEX
u
S
S
h
(162)
p
h
h
S
denote enhanced error indicator function and let
def
EX
EX
qeEX
˜EX = eS p − ẽu p
u p
S
h
S
h
(163)
h
be its dierence with the exact error eSEX
p . For any v ∈ U
; 0 , we have
h
EX
EX
B
qeEX
˜EX ; v = B
eS p ; v − B
ẽ ; v
u p
S
h
h
P
p ;v −
B êEX
;v
= B
eSEX
u p + eˆEX
u
h
∈Th
p ;v −
= B
eSEX
|
h
P
∈Th
{z
S
p
h
h
S
REQ;p (v) −
}
P
∈Th
B eˆEX
;v
u
p
S
h
(164)
0
is the solution of the variational problem:
It follows that qeEX
˜EX
u p
S
h
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
460
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
Find qeEX
∈ U
; 0 such that
˜EX
u p
S
h
P
B eˆEX
; v ∀v ∈ U
; 0
B
qeEX
˜EX ; v = −
u
u p
p
h
∈Th
S
h
(165)
S
We will call qeEX
˜EX the exact equilibrium correction. We now have:
u p
S
h
Lemma 5.1. Consider the above setting and the notations. Then
2
EX 2
p ||
||eSEX
U = (E
h ) −
h
P
∈Th
||eˆEX
||2U + ||qeEX
||2
˜EX U
u
u p
p
h
S
(166)
S
h
EX
EX
Proof. Substituting (162) into (163) and using eˆEX
EX = eS p − êu p , we get
u p
S
h
S
h
h
EX
eˆEX
= eˆEX
EX − qe˜EX
u
p
h
S
u p
u p
h
h
(167)
S
S
and therefore
P
∈Th
||eˆEX
||2U =
u
p
h
S
P
∈Th
2
||eˆEX
EX ||U +
u p
S
h
P
∈Th
||qeEX
||2
˜EX U
u p
(168)
S
h
EX
= 0, from Lemma 4.1 in which we proved the orthogonality
where we have used B eˆEX
EX ; qe˜EX
u
u
P
that ∈Th B eˆEX
EX ; v = 0, for all v ∈ U
; 0 . Using (168) in Lemma 4.3 where we showed that
u
r
P
EX )2 −
2
(E
||eˆEX
||eSEX
p ||U =
EX ||U , we get the proof of (166).
h
p
h
p
h
S
S
p
h
S
h
∈Th
u p
S
h
Sp
h
Let qe˜EX
be the nite element approximation of qeEX
, namely, the solution of the following
˜EX
u p
u p
h
h
S
S
problem:
Find q
Shp
e˜EX
u p
∈ Shp such that
S
h
P
Shp
;v = −
B eˆEX
;v
B
qe˜EX
u
u p
p
h
∈Th
S
h
S
∀v ∈ Shp
(169)
By subtracting (169) from (165) we get the orthogonality condition
Shp
B
qeEX
˜EX − qe˜EX ; v = 0
u p
u p
h
h
S
Copyright ? 2000 John Wiley & Sons, Ltd.
∀v ∈ Shp
(170)
S
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
461
and it follows that
Sp
Sp
2
2
EX
2
h
h
||qeEX
˜EX ||U = ||qe˜EX ||U + ||qe˜EX − qe˜EX ||U
u p
u p
u p
u p
h
h
h
h
S
S
S
(171)
S
Let
Sp
def
h
ũ(0)
= uShp + ẽEX
u p + qe˜EX
Sp
S
h
h
(172)
u p
S
h
denote the enhanced nite element solution. From (163) and (172) we have
Sp
(0)
h
qeEX
˜EX − qe˜EX = uEX − ũS p
u p
u p
h
h
S
(173)
h
S
satises the orthogonality condition
and ũ(0)
Sp
h
;v =0
B
uEX − ũ(0)
Sp
h
∀v ∈ Shp
(174)
We then have
L
EX
E
6||uEX − ũ(0)
|| 6E
Sp U
; ũ(0)
; ũ(0)
h
L
resp. E
; ũ(0)
EX
where E
; ũ(0)
h
h
p
S
h
p
S
h
p
h
h
h
S
(175)
p
h
S
is the upper bound (resp. lower bound) which is computed by
substituting ũ(0)
in place of uShp in Section 3 (resp. Section 4). Using this in (171), we get
Sp
h
2
2
Shp
Shp
L
EX 2
EX
||qe˜EX
||2U + E
6||q
||
6||q
||2 + E
EX
(0)
U
e˜
e˜EX U
; ũ
; ũ(0)
h
u p
S
p
h
S
h
u p
u p
h
h
S
h
S
(176)
p
h
S
Substituting (176) into (166), we obtain
(1)
U; (1)
p
ẼL;
h ; u p 6||uEX − uSh ||U 6Ẽ
h ; u p
S
S
h
(177)
h
where
s
(1) def
ẼU;
h ; uS p =
h
EX
E
h; u
S
2
p
h
Copyright ? 2000 John Wiley & Sons, Ltd.
−
P
∈Th
Sp
EX
h
||eˆEX
||2U + ||qe˜EX
||2U + E
u
; ũ(0)
p
S
h
u p
S
h
h
p
h
2
(178)
S
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
462
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
is the enhanced upper bound and
s
2
2
P
Shp
L; (1) def
2
2
L
EX
EX
E
−
||
||
+
||q
||
+
E
Ẽ
h ; u p =
(0)
EX
eˆu
U
U
e˜
h; u p
; ũ
S
S
h
h
u p
p
h
∈Th
S
h
S
h
(179)
p
h
S
is the enhanced lower bound.
We note that
2 2 2 2
P
(1)
(1)
EX
L
− ẼL;
− E
=
||eˆEX
||2U
ẼU;
(0)
h ; u p
h ; u p = E
;ũ(0)
u
;ũ
S
S
h
h
h
h
p
h
S
p
h
(180)
p
h
∈Th
S
S
is the gap function corresponding to ũ(0)
. It follows that the dierence between the
where eˆEX
u p
Sp
h
S
h
upper and the lower bounds can be decreased iteratively by decreasing the gaps in the enhanced
(1)
U; (1)
nite element solution. We also note that the enhanced bounds ẼL;
h ; u p and Ẽ
h ; u p , are based on
S
h
S
h
the exact solutions of the residual problems êEX
u p and hence they are not computable. Following
S
h
similar steps as in [1] we get the computable enhanced bounds and using the symbol ũ(m) , to
denote ũ(m)
we obtain following recursive relation:
Sp
h
U; S p+kn def
Ẽ
h ; ũh=2(m) =
v
u S p+k 2
U; S p+k 2
P
u
Shp
n
n
2
2
t E
h=2
−
||
+
||q
||
+
E
h ; ũh=2(m+1)
p+k ||
(m)
p+k
U
U
S n
S
h ;ũ
eˆ h=2
(m)
∈Th
n
(181)
e˜ h=2
(m)
u˜
u˜
for m = 0; 1; : : : ; M − 1 where
S p+kn
= ũ(m)
+ ẽũh=2
ũ(m+1)
(m) + q
Sp
Sp
h
h
Shp
S
(182)
p+k
n
e˜ h=2
(m)
u˜
and for m = M , we let
U; S p+kn
U
Ẽ
h ;ũ(M ) = E
h ;ũh=2
(M )
(183)
Similarly a recursive relation can be derived for the lower bound,
v
u S p+k 2 P
L; S p+k 2
p+k
L; Sh=2
Shp
n def u
n
n
2
2
Ẽ
h ; ũ(m) = t E
h=2
−
||
+
||q
||
+
Ẽ
h ;ũh=2(m+1)
p+k ||
(m)
p+k
U
U
S
S
h ;ũ
n
∈Th
eˆ h=2
(m)
u˜
n
(184)
e˜ h=2
(m)
u˜
for m = 0; 1; 2; : : : ; M − 1 and for m = M , we let
v
u p+k 2
P
u S n
L
L
−
|| S p+kn ||2U
Ẽ
h ; ũ(M ) = Ẽ
h ;ũ(M ) = t E
h=2
(M )
h ;ũ
∈Th
eˆ h=2
(M )
(185)
u˜
Let us now employ these constructions for the example problem.
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
463
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
Figure 14. Iterative improvement of the bounds for the energy norm of the error: Exact solutions of the element residU; (m)
L; (m)
ual problems. Graphs of the eectivity indices ẼTh =||eEX
=||eEX
2 ||U , and ẼTh
2 ||U , for m = 0; 1; 2; 3. Note that both
S
S
h
h
the bounds practically converge in three iterations to the energy norm of the exact error
Example 4 (Iterative improvement of the upper and the lower bounds). Let us rst consider
the case that the element residual problems are solved exactly. Figure 14 shows the graphs of
(m)
L; (m)
=||eSEX
=||eSEX
the eectivity indices ẼU;
2 ||U , and ẼTh
2 ||U , for m = 0; 1; 2; 3. We note that both the
Th
h
h
bounds practically converge in three iterations to the energy norm of the exact error.
Next we will report the computable improved bounds. As before, we will denote the estimate
S p+k
obtained by the p-method (resp. h-method) of solution of the element residual problems by EThh
S3 n resp. ETh=2
; recall that these estimates are not upper bounds for the energy norm of the exh
U; Sh8 ; (m)
act error. Let us denote the corresponding enhanced guaranteed upper bounds by ẼTh
U; S 3 ; (m)
ẼTh h=8 ,
respectively and the lower bounds by
L; S 8 ; (m)
ẼTh h ,
and
Figure 15 compares the graphs of the eectivity indices of
L; Sh8 ; (m)
and ẼTh
h
S8
S 8 ; (m)
h
EX
h
=||eSEX
2 ||U with ||e 2 ||U =||eS 2 ||U . We note that after three iterations ẼTh
S
h
h
h
derestimates
and
U; S 8 ; (m)
ẼTh h
S 3 ; (m)
S3
always un-
is the guaranteed upper bound for the energy norm of the exact
error. Figure 16 compares the graphs of the eectivity indices ẼTh=8
h
L; S 3 ; (m)
h
L; Sh8 ; (m)
, and ẼTh
S 8 ; (m)
practically converge to the truth-mesh error estimate, for m = 0; 1; 2; 3. Further, ẼThh
||eSEX
p ||U
h
, and
L; S 3 ; (m)
ẼTh h=8 , respectively.
S 8 ; (m)
U; Sh8 ; (m)
ẼThh =||eSEX
=||eSEX
2 ||U , ẼTh
2 ||U ,
3
U; Sh=8
; (m)
=||eSEX
2 ||U , ẼTh
h
=||eSEX
2 ||U ,
h
n
h=2
||U =||eSEX
and ẼTh h=8 =||eSEX
2 ||U with ||e 2
2 ||U . As before, the estimates converge to the energy
Sh
h
h
norm of the exact error after three iterations.
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
464
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
Figure 15. Iterative improvement of the bounds for the energy norm of the error: Approximate solutions of the elS 8 ; (m)
ement residual problems using the p method. Graphs of the eectivity indices ẼTh
h
and
L; S 8 ; (m)
ẼT h
=||eEX
|| ,
h
S2 U
S 8 ; (m)
ẼTh
,
h
for m = 0; 1; 2; 3. Note that after three iterations
h
S8
U; Sh8 ; (m)
verged to the truth-mesh error estimate ||e h2 ||U , and ẼT
S
h
h
and
U; Sh8 ; (m)
=||eEX
2 ||U , ẼT
S
h
L; S 8 ; (m)
ẼT h
,
h
h
=||eEX
2 ||U ,
S
h
have practically con-
, is a guaranteed upper bound for the energy norm
of the exact error
From Figures 14 and 15 it can be seen that the enhanced bounds converge to constant values
which are dierent from ||eSEX
p ||U . The converged values of the bounds depend on values of p + k
h
S p+kn
and n used for computing the error indicator functions êuh=2
p .
S
h
U; Shp+k ; (m)
Figure 17 shows the graphs of the eectivity indices ẼTh
||eSEX
2 ||U , for p +
h
U; S 3 n ; (m)
=||eSEX
ẼTh h=2
2 ||U ,
h
L; Shp+k ; (m)
=||eSEX
2 ||U and ẼTh
h
=
k = 3; 4; 6; 8; while Figure 18 shows the graphs of the eectivity indices
L; S n ; (m)
and E˜Th h=2
=||eSEX
2 ||U , for n = 1; 2; 3. It can be seen that
3
h
p+k
U; S
lim E˜Th h
; (m)
p+k→∞
p+k
=
L; S
lim E˜Th h
p+k→∞
; (m)
p ||U
= ||eSEX
h
(186)
and
U; S n ; (m)
L; S n ; (m)
p ||U
E˜Th h=2
E˜Th h=2
= lim
= ||eSEX
lim
n
n
3
h=2 →0
3
h=2 →0
h
(187)
The relative dierence between the upper bound and the lower bound will be called the reliability
interval
S p+kn ; (m) def
=
RThh=2
p+k
U; Sh=2
n ; (m)
E˜Th
L; S p+k
; (m)
h=2n
− E˜Th
||uS 2 ||U
(188)
h
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
465
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
Figure 16. Iterative improvement of the bounds for the energy norm of the error: Approximate solutions of the eleS 3 ; (m)
ment residual problems using the h method. Graphs of the eectivity indices E˜T h=8
h
3
L; Sh=8
; (m)
and E˜T
h
=||eEX
|| ,
S2 U
h
3
˜U; Sh=8 ; (m) =||eEX
=||eEX
2 ||U , E
2 ||U ,
S
h
S3
3
L; Sh=8
h
h
for m = 0; 1; 2; 3. Note that after three iterations E˜T h=8 , and E˜T
to the truth-mesh error estimate ||e
3
Sh=2
n
s2
h
3
U; Sh=8
||U , and E˜T
h
S
h
, practically converge
is the guaranteed upper bound for the energy norm of the
exact error
The reliability intervals for the bounds shown in Figures 17 and 18 are shown in Figures 19
and 20 respectively. From Figure 19 (resp. Figure 20) we see that the reliability interval is less
than 0·1 for p + k ¿ 6 (resp. n ¿ 2).
From the above examples we conclude that
1. The improved upper and lower bounds converge to the energy norm of the exact error only
in the case that the element residual problems are solved exactly.
2. In the case that the approximate solutions of the element residual problems are employed, the
improved estimates converge to the the truth-mesh error estimate and not the energy norm
of the exact error.
3. The guaranteed upper bound for the energy norm of the exact error is obtained by taking
into account the explicit estimate correction term as shown above.
4. The improved upper and lower bounds converge after a few iterations, i.e. the reliability
interval converges to a constant value. The reliability interval can be decreased by solving
the element residual problems accurately.
Let us now give examples to show the inuence of the approximate solutions of the residual
problems on the reliability interval.
Example 5 (Bounds for the energy norm of the error using the patch residual method). Let
us group the elements which are connected to the re-entrant corners into patches, !1 ; !2 ; : : : ; !Npat
(see Figure 21(a)) and solve the patch-residual problems. We will denote the enhanced upper and
lower bounds which are computed by solving patch residual problems shown in Figure 21 by
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
466
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
Figure 17. Iterative improvement of the bounds for the energy norm of the error: Approximate solutions of the element
p+k
; (m)
h
U; S
residual problems using the p method. Graphs of the eectivity indices E˜T
h
p+k
; (m)
h
L; S
˜
=||eEX
2 ||U , and ET
S
for p + k = 3; 4; 6; 8, and m = 0; 1; 2; 3
h
h
=||eEX
2 ||U ,
S
h
Figure 18. Iterative improvement of the bounds for the energy norm of the error: Approximate solutions of the element
3
U; Sh=2
n ; (m)
residual problems using the h method. Graphs of the eectivity indices E˜T
h
n = 1; 2; 3, and m = 0; 1; 2; 3
3
L; Sh=2
n ; (m)
˜
=||eEX
2 ||U , and ET
S
h
=||eEX
2 ||U , for
S
h
h
E U; Sh=2n ; (m) and EL; Sh=2n ; (m) respectively. Figure 22 shows the graphs of the eectivity indices of the
3
U; S 3 n ; (m)
˜L; Sh=2n ; (m) =||eEX
=||eEX
improved upper bounds E˜ h=2
2 ||U and the improved lower bounds E
2 ||U
3
3
Th
Sh
Th
Sh
for n = 0; 1; 2; 3. By comparing Figure 22 with Figure 18 we note that the reliability interval
obtained by the patch residual method is smaller than the reliability interval obtained by the
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
467
Figure 19. Iterative improvement of the bounds for the energy norm of the error: Approximate solutions of the element
S
p+k
; (m)
, for p + k = 3; 4; 6; 8, and m = 0; 1; 2; 3.
residual problems using the p method. Graphs of the reliability intervals RT h
h
Note that the reliability interval is less than 0·1 for p + k ¿ 6
Figure 20. Iterative improvement of the bounds for the energy norm of the error: Approximate solutions of the element
S 3 n ; (m)
, for n = 1; 2; 3, and m = 0; 1; 2; 3.
residual problems using the h method. Graphs of the reliability intervals RT h=2
h
Note that the reliability interval is less than 0·1 for n ¿ 2
element residual method, i.e.
S3
n ; (m)
R
h=2
h
S3
n ; (m)
¡RThh=2
(189)
The above example shows that reliability interval can be improved by solving the patch residual
problems in the elements which are connected to the re-entrant corners.
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
468
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
Figure 21. Bounds for the energy norm of the error using the patch residual method: Approximate solution of the patch
residual problems. The meshes of curvilinear quadrilaterals of degree p = 3 constructed within the elements to solve the
p+k
; (m)
h=2n
U; S
patch residual problems and obtain the estimate E
h
: (a) n = 0; (b) n = 1; (c) n = 2
Figure 22. Iterative improvement of the bounds for the energy norm of the error: Approximate solutions of the patch
residual problems as shown in Figures 21(a)–(c). Graphs of the eectivity indices of the improved upper bounds
3
U; Sh=2
n ; (m)
E˜
3
L; Sh=2
n ; (m)
˜
=||eEX
2 ||U and the improved lower bounds E
S
h
h
h
=||eEX
2 ||U for n = 0; 1; 2; 3, and m = 0; 1; 2; 3
S
h
Let us now consider dierent methods of constructions of the edge corrections which were
mentioned in Section 1.
Example 6 (Comparison of the bounds based on dierent methods of construction of the
edge corrections). In this example we show that type of splitting of the jumps for the construction
of the equilibrated residuals in (37), practically, does not inuence the converged values of the
reliability intervals. We considered Ladeveze and Bank–Weiser methods of construction of the
edge corrections as in Example 2 and computed the values of the upper and the lower bounds.
Figure 23 compares the graphs of the eectivity indices of the improved upper bounds
3
U;
S 3 ; (m)
˜L; Sh ; (m) =||eEX
E˜ h =||eEX
2 ||U and the improved lower bounds E
2 ||U for the above types of
Th
Sh
Copyright ? 2000 John Wiley & Sons, Ltd.
Th
Sh
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
469
Figure 23. Comparison of the bounds based on dierent methods of constructions of the edge corrections: Approximate
solutions of the element residual problems as shown in Figure 2(a). Graphs of the eectivity indices of the improved
3
3
U; S ; (m)
˜L; Sh ; (m) =||eEX
upper bounds E˜T h =||eEX
2 ||U and the improved lower bounds ET
2 ||U for q = 1; 2 for Ladeveze and the
h
S
S
h
h
h
Bank–Weiser methods of construction of the edge corrections
splittings. We note that for m¿3 the eectivity indices of the bounds obtained from the four
types of the edge corrections are practically close.
In this section we have shown the iterative improvement of the bounds for the energy norm of
the error and we note the following:
1. If the local residual problems are solved exactly, the upper and the lower bounds converge to
the energy norm of the exact error. If the local residual problems are solved approximately,
the lower bound converges to the energy norm of the error with respect to the enhanced
solution (truth-mesh error). The guaranteed upper bound is obtained by taking into account
the error in the approximate solutions of the residual problems.
2. The converged value of the reliability interval decreases if the local residual problems are
solved accurately. We also note that the converged value of the reliability interval obtained
using the patch residual problems is smaller than the reliability interval obtained using the
element residual problems.
3. The converged values of the upper and the lower bounds are practically independent of the
method of construction of the equilibrated element residuals.
6. CONCLUSIONS
In this paper we addressed the problem of a posteriori estimation of the bounds for the energy
norm of the exact error in the energy norm. The major conclusions of this paper are:
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
470
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
1. The upper bound for the global energy norm of the error is obtained by employing the exact
solutions of the local residual problems.
2. In order to construct guaranteed upper bound based on the computed error indicators, we
must construct an upper estimate for the element energy norm of the error approximate solutions
of the local residual problems. We constructed an explicit estimate which is suciently accurate
for computing the guaranteed upper bound.
3. A lower bound for the global energy norm of the error can be obtained by constructing the
‘gap function’ for the solutions of element residual problems.
4. The improved upper and lower bounds for the global energy norm of the error can be
constructed by employing the ‘equilibrium corrections’. The improved upper and lower bounds,
the global energy norm of the error, can be improved by employing an iterative procedure
which employs the factorized global stiness matrix and the solutions to the element residual
problems.
The objective of this paper was to illustrate the concepts and the procedure to compute the
bounds for the energy norm of the exact error in the nite element solution given in Reference
[1], for two-dimensional heat-conduction problems. In forthcoming papers, we will address to
construction of bounds for the outputs (temperature, heat-ux at a point, ux-intensity factors,
etc.) computed either from the nite element solution or a recovered solution obtained by local
averaging (e.g. the ZZ-SPR). We will also address the case of the elasticity problem.
APPENDIX I
I.1. Construction of the edge corrections for the equilibrated element residuals
Let us describe the detailed construction of the edge corrections ; uS p by the Ladeveze method
h
[2–5] and the Bank–Weiser method [9]. For simplicity we will omit the subscript uS hp . We let
=
p+1
P
i=1
(i) (i)
(190)
where (i) ; i = 1; 2, are the mapped linear functions on and (i+1) ; i = 2; : : : ; p. are the mapped
polynomial functions of degree i. We will construct ; ⊂ @, such that REQ (v) = 0, for all v ∈
p
(). This is achieved by determining (i) such that
Sh;0
REQ (; j ) = 0
(191)
where ; j ; j = 1; 2; : : : ; p + 1, are the mapped element shape functions. Below will employ a
special form of (i) such that they satisfy the following orthogonality condition
Z
0 if j¡i
(i)
; j =
(192)
1 if j = i
where i = 3; : : : ; p + 1, and j = 1; 2; : : : ; p + 1. This enables us to determine (i) ; i = 1; 2, by solving
a small system of equations and (i) ; i = 3; : : : ; p + 1, can be determined explicitly, as shown
below.
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
471
Figure 24. Construction of the edge corrections for the equilibrated element residuals: An example which depicts the
enumeration of the elements jX , and the edges kX in the patch associated with the vertex X
nel(X )
Let X denote an interior vertex node, jX j=1 , a set of elements connected to X , and
X nedge(X )
the set of edges in the mesh patch connected X (see for example Figure 24). Let
k k=1
1
X ∈ Sh be the hat function associated with the vertex X . In the rst step we construct by
employing the orthogonality condition that R
(v) = 0, for v ∈ Sh;p 0 (
):
R
(X ) =
nel(X
P)
j=1
REQ;p
(X |jX ) = 0
X
and we get the following set of equations:
Z
P
(
)
+
sg(kX ; jX )kX X = 0;
RUNEQ
X
X
j
kX
⊂ @jX
j = 1; 2; : : : ; nel(X )
and using (190) together with the orthogonalities in (192) we get
Z
P
(1)
(2) (2)
UNEQ
sg(kX ; jX ) (1)
R X (X ) +
X = 0;
X X + X X
j
⊂@jX
k
kX
(193)
j
k
k
k
j = 1; 2; : : : ; nel(X )
(194)
(195)
(2)
Further, we will also construct (1)
X and X such that the following orthogonalities are satised:
k
k
Z
Z
kX
(1)
X X = 1;
k
kX
(2)
X X = 0
k
(196)
Therefore, we can now write the system of equations (195) as follows:
MX XX = bX
Copyright ? 2000 John Wiley & Sons, Ltd.
(197)
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
472
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
where b X =
n
onel(!X )
−RUNEQ
(X )
X
j
j=1
; M ∈ Rnel; nedge and kX = (1)
for k = 1; 2; : : : ; nedge. In the
X
k
example shown in Figure 24 the above system of equations is




(X ) 
−RUNEQ
X



  1X  
1





1
0
0
−1 








UNEQ




X 
(
)
−R


 −1


X
X
1
0
0
2
2


=



 0
−1
1
0 
3X 




(X ) 
−RUNEQ
X








3








0
0
−1
1
X


 −RUNEQ ( ) 
4
4X
(198)
X
We note that the solution to (197) exists but it is not unique (see [2–5] for details).
The Bank–Weiser equilibration involves choosing a particular solution
X
nedge
= 0;
X
nedge−j
=−
nedge
P
k=nedge−j+1
RUNEQ
(X );
X
k
j = 1; 2; : : : ; nedge − 1
(199)
The Ladeveze-type of splitting involves choosing a particular solution for (197), which minimizes
the ‘2 norm,
s
X 2
nedge
P
k
X
(200)
||X ||‘2 ;{| X |−2 }nedge =
X
k
k=1
|
k=1
k |
After the linear edge corrections have been determined, the higher-order corrections can be computed from
Z
i−1
P (j)
(j) ; i ; i = 2; : : : ; p + 1; ⊆ @
(201)
(i) = − RUNEQ (; i ) −
j=1
Note that the higher-order corrections can be determined explicitly and are dened uniquely for
each edge.
APPENDIX II
II.1. Bounds for the norms of the functions in an element
def
Lemma I.1. Let ˆ = (−1; 1) × (−1; 1) denote the master element and let x = x( x̂1 ; x̂2 ) denote
the transformation of the master element co-ordinates ( x̂1 ; x̂2 ) to x and let


@x1 @x1
@ x̂2 
def  @ x̂ 1

(202)
J = 
 @x2 @x2 
@ x̂1
@ x̂2
denote the Jacobian matrix of the transformation. Then we have
r
max |J | ||v̂||L2 (ˆ)
||v||L2 () 6
Copyright ? 2000 John Wiley & Sons, Ltd.
(203)
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
and
473
p
||v||U ¿ K ||v̂||Uˆ
(204)
K = min(|J |min ((J )−T (J )−1 ))
(205)
where
def
|J |¿0 is the determinant of the the Jacobian matrix in (202) and min ((J )−T (J )−1 ) is the
minimum eigenvalue of (J )−T (J )−1 .
Proof. The proof of (203) is obtained by noting that
Z
Z
v(x1 ; x2 ) dx1 dx2 =
v̂( x̂1 ; x̂2 )|J | d x̂1 d x̂2
(206)
ˆ
and hence
Z
||v||2L2 () =
Z
v2 =
ˆ
Z
v̂2 |J |6 max |J |
v̂2
(207)
ˆ
Further
1 @x2 @ x̂1
1 @x1
@ x̂1
=
;
=−
@x1 |J | @ x̂2 @x2
|J | @ x̂2
1 @x2 @ x̂2
1 @x1
@ x̂2
=−
;
=
@x1
|J | @ x̂1 @x2 |J | @ x̂1
and we obtain


@v 


 @x 


@ x̂1

@x1
def
1
=
∇v =



@v 
@ x̂1



@x2
@x2
T
@v̂ @v̂
def
;
and
where ∇x̂ v̂ =
@ x̂1 @ x̂2


@ x̂2  @v̂ 




@x1 
 @ x̂1 = (J )−1 ∇x̂ v̂
 @v̂ 
@ x̂2  



@x2
@ x̂2
(208)

@ x̂2
@x1 

@ x̂1 
@x1
(209)

@ x̂2

1
 @x2
(J )−1 =
|J |  @ x̂1
−
@x2
−
is the inverse of the Jacobian matrix (202). By employing (206) and (207), we get
Z
Z
2
T
||v||U = (∇v) (∇v) = ((J )−1 ∇x̂ v̂)T ((J )−1 ∇x̂ v̂)|J |
Z
=
ˆ
ˆ
(∇x̂ v̂)T (J )−T (J )−1 (∇x̂ v̂)|J |
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
474
T. STROUBOULIS, I. BABUSKA
AND S. K. GANGARAJ
Z
¿
ˆ
min ((J )−T (J )−1 )|J |(∇x̂ v̂)T (∇x̂ v̂)
¿ K
Z
ˆ
(∇x̂ v̂)T (∇x̂ v̂)
(210)
where min ((J )−T (J )−1 ) is the minimum eigenvalue of the matrix
(J )−T (J )−1

2 2

@x1
@x2
@x2
@x2
@x1
@x1
−
+
−
@ x̂2
@ x̂2
@ x̂1
@ x̂2
@ x̂1
@ x̂2 
1 


=


2 2
|J |2  @x @x @x @x 
@x1
@x2
1
1
2
2
−
+
−
@ x̂1
@ x̂2
@ x̂1
@ x̂2
@ x̂1
@ x̂1
(211)
ACKNOWLEDGEMENTS
The work of T. Strouboulis and S. K. Gangaraj was supported by the U.S. Army Research Oce
under Grant DAAL03-G-028, by the National Science Foundation under Grant MSS-9025110 and
by the U.S. Oce of Naval Research under Grants N00014-96-1-0021 and Grant N00014-96-11015. The support of Dr. Kailasam Iyer of the U.S. Army Research Oce and Dr. Dick Lau of
the Oce of Naval Research is greatly appreciated. The work of I. Babuska, was supported by
the U.S. Oce of Naval Research under Grant N00014-90-J-1030 and by the National Science
Foundation under Grants DMS-91-20877 and DMS-95-01841.
REFERENCES
1. Babuska I, Strouboulis T, Gangaraj SK. Guaranteed computable bounds for the exact error in the nite element
solution—Part I: one dimensional model problem. Computer Methods in Applied Mechanics and Engineering 1999;
176:51–79.
2. Babuska I, Strouboulis T, Upadhyay CS, Gangaraj SK, Copps K. Validation of a posteriori error estimators by numerical
approach. International Journal for Numerical Methods in Engineering 1994; 37:1073–1123.
3. Babuska I, Strouboulis T, Upadhyay CS, Gangaraj SK. A model study of element residual estimators for linear
elliptic problems: the quality of the estimators in the interior of meshes of triangles and quadrilaterals. Computers
and Structures 1995; 57:1009–1028.
4. Ladeveze P. Comparaison de modeles de milieux continus, These de Doctorat d’Etat es Science, Mathematiques,
L’Universite P. et M. Curie, Paris, 1975.
5. Ladeveze P, Leguillon D. Error estimate procedure in the nite element method and applications. SIAM Journal on
Numerical Analysis 1983; 20(3):485–509.
6. Ladeveze P, Maunder EAW. A general method for recovering equilibrating element tractions. Computer Methods in
Applied Mechanics and Engineering 1996; 137:111–151.
7. Ainsworth M, Oden JT. A-posteriori error estimation in nite element analysis. Computer Methods in Applied
Mechanics and Engineering: Computational Mechanics Advances 1997; 142:1–88.
8. Ainsworth M, Oden JT. A-posteriori error estimators for second order elliptic systems, Part 2. An optimal order process
for calculating self-equilibrating uxes. Computational Mathematics and Applications 1993; 26:75–87.
9. Stein E, Ahamad R. An equilibrium method for stress calculation using nite element displacement models. Computer
Methods in Applied Mechanics and Engineering 1977; 10:175–198.
10. Bank RE, Weiser A. Some a-posteriori error estimators for elliptic partial dierential equations. Mathematics of
Computation 1985; 44:283–301.
11. Babuska I, Strouboulis T. The Finite Element Method and its Reliability. Oxford University Press: Oxford, 2000,
to appear.
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
GUARANTEED COMPUTABLE BOUNDS FOR THE EXACT ERROR—PART II
475
12. Szabo BA. Babuska I. Finite Element Analysis. Wiley: New York, 1991.
13. Paraschivoiu M, Peraire J, Patera AT. A-posteriori nite element bounds for linear functional outputs of elliptic partial
dierential equations. Computer Methods in Applied Mechanics and Engineering 1997; 150:289–312.
14. Paraschivoiu M, Patera AT. A hierarchical duality approach to bounds for the bounds of partial dierential equations.
Computer Methods in Applied Mechanics and Engineering 1998; 158:389–407.
15. Harwell Subroutine Library, AEA Technology, Oxfordshire, England.
16. Babuska I, Strouboulis T, Copps K. The design and analysis of the generalized nite element method. Computer
Methods in Applied Mechanics and Engineering, 1999, in press.
Copyright ? 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 47, 427–475 (2000)
Документ
Категория
Без категории
Просмотров
7
Размер файла
551 Кб
Теги
512
1/--страниц
Пожаловаться на содержимое документа