Chapter 5 Other M e t h o d s in Space D i s c r e t i z a t i o n Continuous problems have several fine structures and it is desirable to realize them in the discretizing process. In this chapter, we describe three techniques of FEM; lumping of mass, upwind difference, and mixed elements. The first two are concerned with the dissipative structure, and the last one is related to the variational structure. Those methods are efficient particularly in nonlinear problems as we shall describe in the following chapter. In the last two sections, we describe the method of boundary elements and that of charge simulations. Both of them are associated with fundamental solutions for specific operators and reduce the quantities of computation in some cases.
5.1
L u m p i n g of M a s s
We suppose that f~ C R 2 is a polygon. In Chapter 2, we studied the finite element approximation of the initial boundary value problem
ut-Au=O
i n f ~ x ( 0 , T)
with
Ulo~ = 0
ult= o = Uo(X).
and
This problem, however, is provided with the dissipative structure such as the maximum principle and the L 1 contraction. Those properties remain valid even in nonlinearly perturbed problems and control basic features of the solutions. Numerical schemes kept with those properties, therefore, are desirable from both theoretical and practical points of view. Method of lumping of mass realizes them.
S c h e m e of
Lumping
Recall that {%} is a family of regular triangulations of gt, and Z-h the totality of vertices subject to Th locating in Ft. The set of continuous functions, linear on each T E Th and taking the value 0 on the boundary, is denoted by Vh. For a E Zh, the function Wa E Vh is given by 171
5. Other Methods in Space Discretization
172
10 (x = a)
Wa(X) =
(X e Zh \ {a}).
The family {Wa J a C 27h} forms a basis of Vh and the interpolation operator rrh" W --+ Vh is defined by
7ChV= E
v(a)Wa,
aEZh
where W = {v e C(U) Jv = 0 on Of/}. Each a C Zh takes the barycentric domain denoted by D~. It is nothing but the closed convex hull of barycenters of the triangles in rh containing a as their vertices. Let 1 0
~,(x) =
(xeDa) (xe~\D~).
and denote by Vh the vector space generated by {~a I a C 2"h}. The linear transformation Mh " Vh ---+Vh, sometimes referred to as the lumping operator, is defined through Wa ~-+ Wa. Let M;: Vh 9 -+ Vh be the dual operator associated with the L 2 inner product, and
I ( h = M;: Mh" V~ --~ Vh. Given u0 E W, we can show that the scheme d -~KhUh + AhUh =
0
with
Uh(0) = 7rhuo
(5.1)
in Vh has the desired properties, where Ah indicates the finite element approximation of - A . In the weak form, it is written as d
dt (gh, wa) + (VUh, Vwa) = 0
with
(Uh.(0), Wa)
=
(7rhUO,Wa)
for a E Z-h, where 'uh = MhUh. The operator Mh is an isomorphism and Kt71 = M,71 (M,7:)* is well-defined. Problem (5.1) is equivalent to d'It h
d--7 + K~'A,,~,,,
:
0
with
u,,(O) :
~o.
There are works by T. Ushijima concerning L 2 and L ~176 convergences of approximate solutions. Here, we study L 1 and L ~ structures of (5.1) in details. In this section, Xh and Xh denote the spaces Vh, and Vh equipped with the topology induced by Ll(f~), respectively.
5.1. L u m p i n g of Mass
173
/2-stability Before doing so, we show L p stability and convergence. First, we have a constant C >_ 1 independent of h satisfying
c-'
IIM~x,~ll,.,,,(~) <
xh
I~,(~)_< c IIMhxhll.,(r)
(s.2)
for Xh E X h , T C rh, and p C [1, co]. This implies
(s.a) In fact, take the canonical reference element T of which vertices are P0(0, 0), PI(1, 0), and /:'2(0, 1). Inequality (5.2) is then reduced to the case T = T from the regularity of {%}. Linear functions on T, provided with suitable boundary conditions, form a finite (up to three) dimensional vector space and the desired estimate follows because any two norms are equivalent there. Furthermore, because the extremal cases p = 1 and p = oc are valid, that inequality is uniform in p E [1, oc] by Riesz-Thorin's interpolation theorem. This gives (5.2). It is possible to evaluate C _> 1 explicitly by elementary calculations in use of the least angle of the vertices of T for instance, and is left to the reader. As a consequence, lumping operator is shown to have the following property, where p E [1, cxD]: Uh e X h ---+U
in LP(t2)
=:::a
To see this, take g > 0 and w E W satisfying
M,~
Jlu-
MhUh ---+ U
(5.4)
in LP(t2)
wJls,(a) < e. In use of
- ~ = M~ [~, - ~ , ~ ] + ( M , ~ , ~ - ~) + (~ - ~).
we have
IIMhUh -- UlIL~(~) < c (llu,~ - ~ I,..,,(~) + c + I1'~ - ",,-~'~ll~oo(~)) + IIMu,,,,.~ - wll,r.,o~(~)
9
Here, w is uniformly continuous on f~ and we obtain lim sup I l M h u h - Ullr,(a ) _< Ca. hi0
Because e > 0 is arbitrary, we get the conclusion. Now we study the uniform boundedness of Kh " X h --+ X h with respect to t h e / 2 norm. Taking -Xh E X h and v E W, we have
_< II~hllLp'(~)" IIMh~vlILp(~)<---II~IILp'(~)IIP, vlIL~(~)
5. Other Methods in Space Discretization
174
with 1 / p ' + 1/p = 1, where Ph'L2(f~) ~ Xh denotes the L 2 projection. If {7-h} satisfies the inverse assumption, {Ph} is regarded as a family of uniformly bounded operators on LP(gt) by Theorem 1.5. In this case we have
IIPhvll~.(~) _< c IIv
(5.5)
~,(~)
for p E [1, cx~]. The inequality
I(M;N,,..'v)I _< c ~hll~,,(r~)IIv
L,(~)
follows and hence
IIM;~,.II~p,(~) _< C I ~,.11...(~) is obtained. Combining this with (5.5), we have
I K.~.II~(~) _< C l ~ l ~(n)
(xh e x~).
(5.6)
If p = 2, the inverse assumption is not necessary for this inequality to derive, because IPh r2(n),r2(n) -< 1 always holds true. The L 2 projection P h " L2(f~) --~ Xh has the same property and we have I](Kh) -1X~hIIL~(~)-< C IX~hllL~(a)
(X~h C X~)
(5.7)
if {Th} satisfies the inverse assumption, and without it if p = 2.
Maximum
Principle
In the following theorem, the inverse assumption is not necessary. 5.1. If any T C Th is an acute or right triangle, scheme (5.1) is provided with the maximum principle so that
Theorem
0 <_ fh E Xh
==~
e-tg':'n"fh >_ 0
(5.8)
holds for the associated semigroup {e-tK;i ' Ah }t>_O on Xh. To see this, let Bh be the set of vertices subject to rh locating on Cgf~. Even for a E Bh, we can define a continuous function Wa on f~, linear on each T E ~-t~by I
Wa(X) -'~
0
(x = a) (:1;, E •h U ~h \ {a}).
Then we have
E aEZhUBh
wa(a)= 1
(z E-~).
(5.9)
5.1. L u m p i n g of Mass
175
We can also take the barycentric domain Da even for a
Wa(X) =
9•h" Taking Wa as
(xeDa) (X 9 \ Da)
1 0
similarly, we have
Z Wa(X)---1
(x
9-~).
(5.10)
aE:ThUI3h Let a, b E Zh t2 Bh with a =fl b. It is obvious t h a t
(Wa, Wb) "--0.
(6.11)
(V~a, VWb) <_0
(5.12)
The inequality
is a consequence of the following l e m m a proven similarly to L e m m a 1.9. Lemma
5.2.
We have VWa" VWb ~ 0
(5.13)
on each T E ~-h if a, b C ffh U Bh and a sC b. Proof: I f a E Z h U B h is not a v e r t e x o f T E ~-h, then Vwa = 0 o n T and (5.13) is obvious. Let a be a vertex of T and Ta be the face of T not containing a. Because the left-hand side of (5.13) is invariant under the r o t a t i o n of the variable x, we may assume t h a t T~ is located on the line x2 = 0 and a lies in {x2 > 0}. In this case OWa/OXl = 0 and hence 0w~ 0Wb V w ~ . Vwb = Ox2 " Ox2 follows. Let {Pj [ 0 < _ i < _ 2 } be the vertices o f T with P0 = a. Because T i s a r i g h t or acute triangle, the vector e = ( 0 , - 1 ) is expressed as 2
e = ~_s aj eoj j=l
with some aj >_ 0, where e0j = P 0 ~ . j = 1, 2. It holds t h a t 0W a
Ox2
Here e 0 j . Vw~ is a negative c o n s t a n t on T for
2
= e . VWa = ~ a j e o j j=l
VWa 9 <_ O.
5. Other Methods in Space Discretization
176
On the other hand we have e0j Vwb 9 _> 0 (j = 1, 2) if b -r a, and
Owb > 0 O3g 2
--
holds similarly on T. Therefore, inequality (5.13) follows and the proof is complete. We are now ready to give the following.
Proof of Theorem 5.1: Note that (5.8) is reduced to A > 0, fh E Xh
===> m_ax(Ih + AK,71Ah) -1 fh _< max{0, m_ax fh}. fl
(5.14)
f2
Here, the well-definedness of (Ih + AKs
-1 is also included. In fact, (5.14)implies
fh <_ o
~
(I,, + ~ K s
f. <_ 0
fh > 0
==~
(Ih + AKs
-1 fh >_ 0.
or equivalently,
Then, tile fundamental relation
e-tKs assures (5.8). L e t [ . ]+ = m a x { . , 0 } .
Ah
=
lim (Ih + tn- 1Ks lAb) -n 71,'----+00
(5.1,5)
Property (5.14)follows if max_ Uh _< max_ rrh [fh]+
(.5.16)
is proven for Uh C= Xh and fh = (Ih + AI(s In fact, then fl + = 0 implies u h+ = 0 and hence fh = 0 gives 'uh = 0 is obtained. This means the well-definedness of (Ih + AI(,71Ah)-I " X,, ---+Xh and also (5.14). For (5.16) to prove we may suppose that the maximum of the left-hand side is attained at some a E :2h with a non-negative value, because the right-hand side is non-negative. Suppose those situations. It holds that
(a~, ~h) + A (v,,,h, Vx~) = (L, ~,~) for Xh E Xh. Let Xh = IDa"
(a~, tea)+ A (v,,,,,,, W*,a) = (/,,, r&) We have (gZh, g&) = uh(a)Isupp Wal
and
(-fs ~ )
= fh(a) supp ~ l .
(5.17)
5.1. Lumping of Mass
177
Furthermore, the equality (VUh' VWa)
---
E Uh(b) (VWb, VWa) bEZh
=
~
(~(b) - ~ ( a ) ) ( W ~ , W a ) + ~.(a) ~
bEZh
--
(W~, W a )
bEZh
~
bEZh\{a}
(~th(b) -- ,/th(a))(VWb, VWa) -- "/th(a) ~
(V~Ub, Y e a )
bEI3h
holds by (5.9). Because a E/-h attains a non-negative maximum of uh we have (Vuh, VWa) ~ 0 by (5.12). Those relations imply
m~x_ ~ = ~(~) <_A(a) <_m_~xfZ gt
or equivalently (5.16). The proof is complete.
Ll-contraction Now, we proceed to the following. T h e o r e m 5.3. Under the assumption of Theorem 5.1, scheme (5.1) has the property of
L 1 contraction expressed as O <_ fh E Xh
=:~
f Mhe-tK[lAhfhdx
(5.18)
where {e-tK: 1mh}t>_O denotes the associated semigroup on Xh. Proof: Similarly, property (5.18) is reduced to ~ M~ (Zh + )~KhlAh)-l fh dx < faM~fh dx
(5.19)
for 0 <_ fi~ E Xh and )~ > 0. We have only show it for fi~ = Wa with a E Ih. Writing
Uh = (In + AKhlLh) -1 Wa = E
~bWb~
bEIh
we have ~b ~ 0 by (5.14). Relation (5.17) with A -- Wa and )(~h -- wc is expressed as E
~b { (Wb, We) -~ /~ (VWb, VWc)} : (We, We),
(5.20)
bEIh
where c E :Z-h. We shall show that
E = ~ M h ( I h + ) ~ K h l L h ) -1 Wa d x -
= ~
~ bEIh,CEBh
c~ (w~, w ~ ) .
MhWa dx
(5.2~)
178
5. Other Methods in Space Discretization
This implies (5.19) by (5.12) and ~b _> 0. In fact, we have
b6 Ih
b6 Ih ,c6 Ih U Bh
b61 h ,c6 Bh
c6 Ih U Bh
c6 B h
b,cE I h
by (5.10) and (5.20). Applying (5.9) for the last term, we get
bE Ih
cE Bh
cE Bh
This means (5.21)by (5.11).Property (5.19) has been proven. Ll-estimate
In connection with the L 1 structure, elliptic estimate (1.69) is also worth noting. We shall study its discretization here. We emphasize that the convexity of f~ is not necessary in the following theorem. Furthermore, not all T E rh must be acute or right triangles. Remember, however, that {rh} is supposed to be regular. T h e o r e m 5.4. If {7h} satisfies the inverse assumption, then q E [1,2) admits C > 0 satisfying
IIA;~K,,PhlI~I(.),..,,~(~) <
(s.22)
c.
If any corner point of [2 does not have the angle 27r furthermore, the convergence lira hi0
IlAhlI(,~Pj~v -
A-lv[Iwl,q(n)
= 0
(5.23)
holds for v E LI(f~), where Ph denotes the L 2 projection. Proof: Because {Th} satisfies the inverse assumption, {Ph} is regarded as a family of uniformly bounded operators on L2(f~) and we have (5.5). If we apply (5.6) for p = 1, inequality (5.22) is a consequence of (1.69). Taking v E Ll(f~), we get it as
Turning to the proof of (5.23), we note that it is reduced to limhlo IA-h:I(hPhv -- A -lvllH:(a ) - 0
(v E L2(~)) .
(5.24)
5.1. Lumping of Mass
179
To see this, let us admit (5.24) for the moment. Given v E Ll(f~) and e > 0, we can take w E L2(f~) satisfying IIw - VlILI<~) < c. In use of (5.22), (1.102), and (1.69), we have
IIA-;lKhPh~- A-iv
wl,q(f~)
]A-s
_<
(v - w)[lWl,~(~) +
IIA-~ ( w - ~)1 w~,~<~)
<_ C (c + IlA~'KhPh~ -- A - i ~ l l / r ( a ) ) . Therefore, letting h I 0 and e ; 0, we get (5.23). We now prove (5.24). Taking v E L~(a), we have [ A h l K h P h - A -1] v
=
[ A h l p h - A - I ] W + A h l p h ( I ( h - y) Phv
=
( R h -- I ) A - i v
+ R h A -1 (Kh - I ) Phv,
where Rh" V --~ Xh denotes the Ritz operator and relation (1.24) is m a d e use of:
Rh A-1 = Ah 1Ph Because v E L2(f~), it holds t h a t A - i v E V = H~(f~). We have lim hlO
II(R
- i) A - i v
V
= 0
by (1.39). On the other hand we have IIRhllv, v ~ c , and therefore the relation lim IIR~A -1 (K~ - I ) P , ~ l l v hi0
= 0
follows from lim IIA -1 (Kh - I ) P ~ l l ~
= 0
hi0
(5.25)
Thus, we have only to derive (5.25) for ( 5 . 2 4 ) t o prove. We shall show w-lim (Kh -- I) Pt~v = 0 hlO
in L2 (f~).
(5.26)
T h a n k s to the assumption to f~, the operator L -1 9L2(f~) + U~(f~) = V is compact. Relation (5.26) implies (5.25), and in this way, (5.24) is reduced to (5.26). For the proof of (5.26), we note t h a t (5.6) with p = 2 implies
II(Kh -- I)PhlIL=<~),L=<~) ~ C. Therefore, (5.26) is reduced to the case of v E W similarly. In fact, assume (5.26) for v = w E W. Given v E L2(t2) and c > 0, we take w C W satisfying
5. Other M e t h o d s in Space Discretization
180 We put Jh =
(Kh -
I) Ph for simplicity and take u E
L2(f~).
T h e n we have
I(&v,~,)l < II& ( v - ~)IIL=(~)II~llL=(a) + I(&w,,~)l _< g e Ilull~=(~) + I(Jhw,~)l. Letting h + 0 and g + 0, we get (5.26) for v E L2(f2). a similar argument now reduces (5.26) (for v E W ) to
(5.27)
lim ((Kh -- I) Phv, w) = 0 h,t O
for w C W. We have only to verify (5.27) for v, w E W. For this purpose we note t h a t
((ir
- ~) e~,,, ~)
=
(M;M,,P,,v, ~) - (e~v, w)
=
(MsPhv, MhPhw) - (Phv, w)
=
(M,,Phv, MhP,, (I - re,,)w) + (Mhre,,v, M h r e , , w ) - (Phv, w) + (MhPh (I - re,,)v, Mhrehw).
For the first term of the right-hand side, we make use of (5.3) and (5.6) as
IIM,~PhvlIL,(~) < c IIP,,v gl(~'2) ~ C II,, L1(~2) and
IIMhPh (I
-
reh)w
IL=(a) < IIPh (I -- ~h) WlIL=<~> ~ C I1(1 -- ~h)W L=(~),
respectively. It follows that
I(MhPhV, MhPh (I - reh) W)[ < C l v l l z x ( a ) [ ( I - reh) ~IIL~(~) ~ 0 from the uniform continuity of w c W . Similarly, we have
I(MhPh (I -- re,,)v, Mhrehw)l
,0.
For the rest terms we make use of (1.40): lim liP,,,,- vl e=(~) = o hlO
and also lim []Mhrehv- V[[f~o(a) = lim [[MhrehW -- W [L~(a) = 0, h~0
h~O
which is a consequence of the uniform continuity of v and w. We get that
I(A6,ro, v, MhrehW) -- (Phv, w)l - - ~ 0 as h i 0, and relation (5.27) has been proven.
[]
An analogy of the above theorem is verified for schemes without lumping more easily.
5.2. Upwind Finite Elements
5.2
181
U p w i n d Finite E l e m e n t s
Continue to suppose that Ft C R 2 is a polygon. Let b(x) be a Lipschitz continuous vector field, and f, g scalar fields on t2, respectively. The problem -Au+b.
Vu=f
inf,,
u=g
oncgf~
(5.28)
describes the phenomenon of convection-diffusion, so that u(x) is subject to the flow generated by b(x) and diffuses by itself. In this sense, the system keeps basic features of the problems in fluid dynamics. If the convection term is stronger than that of diffusion, then the standardized m e t h o d of w is not efficient. Here, one sees t h a t the value of u at x is under the influence of its values upwind with respect to b(x). Upwind finite element method is based on this observation and we describe the simplest case in this section. We make use of the following notations of the previous section. Let {Th} be a family of regular triangulations of acute type so that any T E Th is an acute or right triangle. The set of piecewise linear continuous functions on f~ is denoted by Wh, and Vh = Wh M V, where V = H0X(f~). The sets of interior and boundary vertices are denoted by 2-h and Bh, respectively. Each a E :Z-hU Bh admits We E Wh with the value 1 at a and 0 at each b C :Z-hU Bh \ {a}. Interpolation operator 7rh'C(gt) ~ Wh is defined by
rrhf= E
f(a)Wa"
aCZhUBh
We also make use of the lumping operator Mh " Wh -+ Wh and set MhXh = Xh in short. Given a C 2-h, we take T(a) E rh which meets the vector b(a) with the origin a, and set
Bauh = ~
(b(a).
VUhlT(a) ) Wa.
aEIh
We study the scheme to find Uh E Wh satisfying Uh -- r~hg E Vh and ( v u , . Vx ) +
=
(5.29)
for Xh E ~ . We can show the following Theorem
5.5. Suppose one of the following conditions:
1) Each a E Ih admits a sequence a l , . . . , am C Ih, am+t E 13h with ai+l is a vertex of T(ai) for 1 <_ i < m. 2) Each T E Th is an acute triangle.
Then, (5.29) is uniquely solvable and the discrete maximum principle f, g < 0
follows. (f, Xh).
~
uh < 0
(5.30)
The right-hand side of (5.29) may be replaced by the term without lumping,
5. Other Methods in Space Discretization
182
Discrete maximum principle (5.30) can imply the L ~176 bound of the Ritz operator (1.59) and also IIu -uhliL~(a) = O(h), although details are not described here. If p (~h, Xh) is added to the left-hand side of (5.29), the conclusion of Theorem 5.5 holds without the assumption, where p > 0.
Proof of Theorem 5.5: From the Fredholm alternative, the solvability follows from the uniqueness in (5.29). The latter is a consequence of the discrete maximum principle. We have only to derive (5.30), supposing that 'uh E Wh is a solution of (5.29). It is proven similarly to (5.14) and we shall describe the sketch. Because uh = 7rhg on oqf~, we have only to show urn(a) 5 O, supposing that m a x g u h is attained at a E Zh. Let Xh = w~ in (5.29). We have
(w,,,, v x , ) > 0
(5.31)
(BhUh, Xh) = supp Wa] (b(a). Vuh)r(a) --> 0.
(5.32)
and
Because the right-hand side is non-negative, the left-hand sides of (5.31) and (5.32) are equal to 0. If the first assumption holds, then a~ attains maxg'uh for i = 1, 2 , . . . , m + 1. This implies m a x ~ u h _< 0 by 'uh = rrhg _< 0 on oqf~. If the second assumption holds, we have (VWa, Vwb) < 0 for a : / b . Then, the inequality (5.31) implies that 'uh is a constant on f~. Then m a x g u h _< 0 follows similarly. To prove the latter part, we take A > 0 and u~, E Wh satisfying u h - rrh9 E Vh and
(~h, ~,~) + a (v,,,,h, vxh) + (~h,,,,, ~,~) = (~-Tf, ~,~)
(5.aa)
for )C~ E Vh. Without those assumptions, we can show m_ax a uh _< Inax { m_ax a rrhf, I ~ X rrh } .
(5.34)
In fact, if a E 5~ attains max~u,, we have (B,,u,,, ~,,) >_ O. Then, relation (5.33) implies (5.34). The proof is complete. []
5.3
Mixed Finite Elements
The m e t h o d in consideration is concerned with the variational problem with constraints. To describe the idea, we take the Dirichlct problem for the Poisson equation Ap = f
in [2
with
p = 0
on cgf~,
(5.35)
where f~ C R 2 is a polygon. Putting 'u = Vp, we llav(; u-
Vp = 0
and
V . ',, = .f.
(5.36)
183
5.3. Mixed Finite Elements
Let V = L2(f~) 2, W = H~(~), and B ' = V : W ~ V. The dual operator B = ( B ' ) ' : V "~ V' --, W _~ W' is defined through the representation theorem of Riesz. Similarly, f 6 L2(ft) determines 9 e W', regarded as an element of W, by (g, q)w = (f, q)L2(fl) for q 6 W. Under those preparations, (5.35) has an abstract form, to find (u,p) 6 V x W satisfying
u + B'p = 0
B u = 9.
and
(5.37)
As' we shall see, this problem is realized as an Euler equation on u for a variational problem with a constraint, where p acts as a Lagrangian multiplier. Then, it is transformed into another variational problem without constraints, where the solution (u, p) is to be a saddle point. Method of mixed (hybrid) finite elements arises naturally as a discretization of those structures. Although it had been applied to many problems in engineering related to solid mechanics, hydrodynamics, electro-magnetic theory, and so forth, the discovery of such a fine structure made it much more reliable practically. This section describe the fundamental ideas. Abstract
Theory
We take the abstract theory first. Let V and W be real Hilbert spaces, and ,4 and B be bounded bilinear forms on V x V and V x W, respectively. We have
1"t4(u, ~)1 ~ C1 for u , v E V, where C1 > 0 is a constant. variational problem with constraints,
II~lIv IIvll~
Given g c W and F E V', consider the
(5.39)
inf ,7,
u
where U = {v E V l B ( v , q ) =
(5.38)
(9, q)w for any q E W} and
y(~)
= ~
A(v, ~) - F(~)
The bounded linear operator B : V ~ W is defined by
(B~, q)~ = z(~, q)
(q ~ w ) .
We have U = {v E V I B y = g} and hence U # 0 if and only if g 6 Ran B, which means the existence of Uo E V satisfying Buo = g. Then it holds that U = {u0} + Ker B. It is a closed atiine space in V and the s t a n d a r d argument described in w the condition
A(v, ~ ) ~ ~ I vll~
(v
e
Ker B)
applies. If
(5.40)
5. Other M e t h o d s in Space Discretization
184
holds with a constant a > 0, variational problem (5.39) is uniquely solvable. The solution u C V is characterized by u-
u0 E U
and
.A(u, v) = F(v)
(v C Ker B ) .
(5.41)
Let c~ : V' + V be the canonical isomorphism defined through Riesz' theorem, and A : V + V the bounded linear operator associated with A:
A ( ~ , ~) = (Au, ~ ) ~
(~, ~
9v ) .
Tile second equality of (5.41) is expressed as
(oF-
Au, v)v = 0
(v
9Ker B ) ,
or equivalently, a F - A u E (Ker B) • in V. Here, the dual operator of B is realized as B' : W + V through Riesz' theorem: V' " V, W' _~ W, and then (Ker B) • = Ran B' follows. Under the assumption that Ran B' : closed
C V,
(5.42)
u solves (5.41) if and only if it admits p E W satisfying
A u + B'p = o F
B u = g.
and
(5.43)
Note that p acts as a Lagrangian multiplier here, and (5.37) obeys a form of (5.43). If one prefers to tile weak formulation, (5.43) may be replaced by
.A(u, v) + B(v, p) = F(v)
B(u, q) = (g, q)w
and
(5.44)
for a n y v C V a n d q C W . Problem (5.43) is reduced to another variational problem without constraints. In fact, (u, p) solves (5.43) if and only if it is a stationary point of 1
,]('u, q) = -~A(v, v) - F ( v ) + (q, B v - g)w
(5.45)
on V • W. Actually (u, p) is a saddle point,. To see this, take q E W and v C V. If (u, p) solves (5.43), then we have 2(,,,, p) - 2(,,,, q) = (~, - q, B,,, - 0 ) ~
= 0
and
! A ( . . , .~,)- ~A(.,,,, . . ) 2
F(.v-
u)+
(~, B ( v -
~))~ 1
= A(.,~, v - u) - r(,,, - ..,) + (p, B ( v - ..,))~ + ~ A ( v - .,,,, ,,, - u) 1 A ( , v - ,u, v - u) a F , v - "u)v + -~
=
(Au + B'p-
-
_1 2 A (v - u, v _ u) > O.
185
5.3. Mixed Finite Elements
In particular, the relation 2 ( u , q) _< 2 ( u , p) _< ST(v, p)
(v e V, q E W)
(5.46)
follows as is expected. Relation (5.42) is equivalent to the existence of k > 0 satisfying
IlB'qllv =
sup
~o~
~li~
> kllq
-
for q C (Ker B') -L. Because (Ker B') • = Ran B, this means sup
B(v, q) > k Ilq w
v~.\{0~ II~ll. -
(q E Ran B).
(5.47)
Here, the closed range theorem says that Ran B' C V is closed if and only if Ran B C W is so. Because (Ker B) -c = Ran B', similarly it is equivalent to the existence of k > 0 satisfying
sup
qeW\{0}
15(v, q) > k II~llv q Iw --
(v e R a n / 3 ' ) .
(5.48)
We have the following. T h e o r e m 5.6. Let V and W be Hilbert spaces over R, and ~4 and 13 be bounded bilinear forms on V x V and V x W , respectively. Suppose (5.38), (5.~0), and (5.47) (or equivalently (5.~8)). Let A : V ---+ V and B : V ---+ W be the bounded linear operators associated with A and 13, respectively, and take g E Ran B. Then, the solution (u, p) r V x W of (5.~3) (or ( 5 . ~ ) ) exists, and is characterized as a saddle point of 2 defined by (5.~3). Here, u C V is unique, while p E W is unique up to an additive element in K e r B ' . Actually, p is taken uniquely in Ran B, and then we have Ilullv + IIPlIw < C (IIFIIv, + Ilgllf),
(5.49)
where C = C (5, C1, k) > 0 is a constant. Proof: From the above arguments follow the equivalence between (5.46) and (5.43), the existence of the solution (u, p), and the uniqueness of u. The orthogonal decomposition u = Ul + u2 with Ul C Ker B and u2 C Ran B' gives Bu2 = g, and assumption (5.48) assures
I1~11~ _< k'-I Ilgll~. Furthermore, A u l = ~ F -
Au2 - B'p and hence
A(ul, v) = f(?3) -- A(u2, v) follows for v E Ker B. Putting v = ul, we obtain
II~lll. _< ~_1( FII., + 61 lu~ll.) _< ~-1 (llrl ., + c~k-' I gll~).
186
5. O t h e r M e t h o d s in Space Discretization
Similarly, the orthogonal decomposition p = Pl + P2 with Pl E Ker B ~ and p2 E Ran B gives B~p2 = o F - A u and Ilp211w <_ k - i (IiFlIv , -}-Ca IlulIv) by (5.47). Inequality (5.49) follows and p is taken uniquely in Ran B except for an additive element in Ker B'. T h e proof is complete. [] Poisson
Equation
Above t h e o r e m justifies ( 5 . 3 6 ) i n the following way. Let V = L2(f~) 2, W = H~(f~), A(u, v) = (u, v)c2(n)2, and B(v, q) = - (v, Vq)r2(n)2. The bounded linear operator t9' : I/V --~ V is realized as B ' = V. From Poincar6's inequality, IIB'qllv -IIVqlIL~(n)~ provides a n o r m on W. Hence d a n B = w , Ran B ' = v , (5.47), and (5.48) follow. Given f E L2(f~), we set F = 0 and define g E W ' ~_ W by ( g , q ) w = (f,q)L2(n) for q E W. The abstract theorem assures a saddle point (u, p) E V x W of 2(v,q)=
1
~11 v 2L2(n)2 -- (v, Vq)c:(a)2 -- (f, q)L:(n)
satisfying (5.43). We have (U, V)L2(f~)2 ~- ( V p , V)L2(12)2
and -
('a, Vq)L2(fl)2 : (f, q)L2(n)
for any v E L2(gt) 2 and q E H~(f~). This implies p E H~(f~), Vp = u E H I ( ~ ) , and V u : f E L2(ft) 2. In particular, p E H2(f~)VI H~(ft) is a solution to (5.35). P r o b l e m (5.35) is formulated differently in this framework. Let V : H(div, f~):
H(div, ~) : {~ c r~(~)~ I v . ~ c r~(~)} It is a Hilbert space with the inner product
(~, V)v : (u, ~)~(~) + ( v . u , v . ,v)L~(~) 9 Taking W = L2(ft), we set A(u, v) = (u, 'V)L2(fl)2 for u, v E V, and B(v, q) = ( V . v , q)r2(a) for v E V, q E W. This time tile b o u n d e d linear operator B : V + W is realized as B = V., and we have Ker B = {v E V I V . v = 0}. Condition (5.40) holds as A(v, v ) =
IIv ~/
(v E Keg B ) .
We shall show t h a t inequality (5.47) holds for any q E W = L2(~). This implies Ran B = W in particular. In fact, take (b E Hg(f~) satisfying V - V c b = q. We have [IV011L=(a) _< C IlqlIL=(a) and hence v = Vq5 E V satisfies v 2V -- I V r IL2(~) 2 2 2 -[- ]lV " V 4) 2g2(fi) _< ( C 2 + 1) Iq IL2(fl) 9
187
5.3. M i x e d F i n i t e E l e m e n t s
T h e n (5.47) follows as
z3(~ q) ,
=
I~lIg
Let F = 0 and g = f
llqll ~
w
1
>
IlVllv -
v/C 2 + 1
Ilqllw.
9 W = L~(f~) in (5.45). There exists a saddle point (u, p)
9VxW
of d(v,q)
1
2
: ~ Ilvll L 2 (a)2
-t-
(q, V " v - f)L2(a)
Relation (5.43) reads;
- (~, ~ ) ~ ( ~ ) : (p, v
~9) ~ ( ~ )
and ( V . u, q)L2(a) = (f, g)c2(a), where v C H(div, f~) and q C L2(f~) are arbitrary. This implies f = V - u , p E H2(f~), and Ploa = 0. In particular, p E H2(ft) A HI(f~) solves (5.35).
Vp = u,
The Stokes System In the above problem, f~ may be an arbitrary b o u n d e d Lipschitz domain in I~n. We have an i m p o r t a n t example formulated as above for this case; the stationary Stokes system described as -Au+Vp=f
and
in f2
V.u=0
(5.50)
with u = 0
on 0~.
(5.51)
Here, u and p are vector and scalar fields on f~, respectively. For V = HI(f~) n and W = L2(f~), its weak form is introduced; to find (u, p) E V x W satisfying
~(
V u <9 Vv) d x - (p, V . v)r2(n) : ( f , V)L2(a)n
(5.52)
and
(v.~,
q)~(~) = 0
(5.53)
for any v E V and q E W. Here, OU i z,3
OV i
5. O t h e r M e t h o d s in Space Discretization
188
for u = (u i) and v = (vi). We introduce the b o u n d e d bilinear forms
A(~, ~) = s ( w | w )
dx
and
and take F ( v ) = (f, V)L2(n). and 9 = 0. Condition (5.40)is a consequence of Poincar6's inequality. To examine the closedness of Ran B or Ran B', let q E Ker B'. This means (q, V V)L2(a 9 )= 0 for any v E V. The space Ker B' is composed of constant functions and hence Ran B = (Ker B') • = {q E L2(~) I (q, 1)L2(n) = 0} follows. If 0f~ is Lipschitz continuous, it is known t h a t any q E Ran B admits v E V satisfying V-v = q
Ilvll~ <
and
C llqlIw,
(5.54)
where C > 0 is a constant determined by f~. In this way condition (5.47) holds true. In use of
IVvll~= ~
~
g~<~)
1,3
for v = (vi), problem (5.50) with (5.51) is formulated as (5.44), or finding a saddle point of d ( v , q) -- ~1
Vvll~-
(f, V)L2(n), 4- (q, V ' V ) g 2 ( n )
oil H~(f2) n x L2(f~), where T h e o r e m 5.6 is applicable. The variational structure induces the discretization. Under the assumptions of Theorem 5.6, we prepare a family of finite dimensional subspaces of V x W denoted by {Vh x Wh}h>0. Assume the property that any (v, q) E V x W admits a sequence {(vh, Wh)} in (Vh, Wh) E ~ X Wh satisfying lira (ll v - vhlIv 4-IIq - qhlIw) = o. td.O
The linear operator Bh : Vh---* Wh is defined by
13(vh, qh) = (S,,vh, qh)w for vh E Vh and qh E Wh.
5.3. Mixed Finite Elements
189
It is reasonable to assume (5.40) and (5.47) with V, W, and B replaced by Vh, Wh, and gh, respectively, with constants ~ > 0 and k > 0 independent of h: inf ~or~\~o~
r
inf
(5.55)
~h) ~___5,
I1~11~ sup
B(Vh, qh)
> k.
(5.56)
Then, Ritz-Galerkin method (5.44) is realized as to find (Uh, Ph) E Vh • Wh satisfying
.A(~h, vh)+ Z(vh,ph) = P(vh)
(vh ~ Vh)
(5.57)
and
s ( ~ , q~) : (g, q ~ ) ~
(q~
9w ~ ) .
(5.s8)
Those relations are equivalent to
AhUh + B~hPh = PhcrF
and
BhUh = Qh9,
respectively, where Ah : Vh --~ Vh denotes the linear operator associated with Alyhxvh, B~ " Wh ~ Ws --~ l/h ~ V/~ the dual operator of Bh " Vh ~ l/Vh, and Ph" V ~ Vh and Qh : W ~ Wh the orthogonal projections. This means that Qhg E Ran Bh is necessary for this discretized problem to have a solution. Then Theorem 5.6 applies and there exists a unique (uh, Ph) E Vh x Ran Bh satisfying (5.57) and (5.58). This solution is also a saddle point of ,7 restricted to Vh x Wh, with ,~ defined by (5.45). Based on those considerations, it is natural to assign the requirement QhRan B c Ran Bh. Taking orthogonal complements in Wh, this is equivalent to assuming Ker B~ C Ker B'. Under those circumstances the following theorem holds, where (u, p) denotes the unique saddle point of J on V x W belonging to V x Ran B. T h e o r e m 5.7. We have Uh]ly + l i P - Phllw < C ( inf Ilu - Vhl]y + inf liP-- qh]lw~ \ VhE Vh qh E Vh / with a constant C > 0 independent of h. Proof: Equalities (5.43) and (5.57)-(5.58)imply A ( u - Uh, Wh) + B(Wh, p - Ph) = 0
and
13(u - Uh, rh) = 0
for any wh E V~ and rh E Wh. Taking vt~ E Vh and qh C Wh arbitrarily, we have A(v~ - uh, wh) + S(w,~, q~ - p~) = A(vh -- u, wh) + S(w~, q~ -- p) and B(~
-
~,. <~) : s ( ~ . - ~. <~).
(5.59)
5. Other M e t h o d s in Space Discretization
190 or equivalently,
(5.6o)
Ah(Vh -- Uh) + B'h(qh -- Ph) = PhA(vh - u) + PhB'(qh - p)
and s,,(v,
- u,,) = Q,S(v,,-
(5.61)
u).
Theorem 5.6 is applicable to (5.60) and (5.61) with V, W, A, and B replaced by Vh, Wh, Ah, and Bh, respectively. We have I~h - ~hl v +
Iqh - p , llw
<_ c tlPhA(v,,. - ~)llv + C ( PhB'(qh -- P)I v +
QhB(v,,.- ~) ~) _< C ((tlmll + IIBII)Ilvh- **llv + IIBII IIq,~ - PlIw),
and hence II?t
~
'/t h v + I p
-
PhlIw -----
I1~, - v,,. v + I p
-
qh
_< c ' (11,,- vhllw + l i p -
w
+
Iv,,. - ",,,,. Iv + Ilqh
-
p,,.
w
qhllw)
follows. Thus, inequality (5.59) has been proven. Although it is not easy to construct such a family {Vh x Wh}h>0 with the required properties in T h e o r e m 5.7, many ways are known now.
5.4
Boundary
Element
Method
(BEM)
If one makes use of the fundamental solution, then tim boundary value problem is reduced to an integral equation on the boundary. Boundary clernertt method arises, associated with that structure. Let f~ C R n (n = 2, 3) b c a bounded domain with Lipschitz t)oundary 0f~, and consider tile boundary value problem -A'~J = f
in [2,
'u-- g
on ~)f~.
The fundamental solution F(:r) =
1
~
1
log T~
1
47r [:,:[
(rz-- 2)
(rz = 3)
satisfies - A F = ~(0), where d (lcnot(.'s Dira("s delta function. If on(: takes w(:r) =
.[~ F(:r-
:,j).f(y) dy,
(5.62)
5.4. Boundaly Element Method
191
problem (5.62) is reduced to the case f = 0: -Av=0
int2,
v=g
on0t2
(5.63)
The fundamental solution induces also single and double layer potentials
F(x) = foa V(y)F(x - y) dS v and
o
~ r ( ~ - v) ds~
respectively. For ~ E 0t2, we put A+(~)=
lira
xCgt+~
A(x),
where ~+ = ~c and ~ _ = ~. If V is continuous, we have [F] + =_ F + -
F_ = 0
on 0~.
(5.64)
Given x c 0t2, we can take the principal value for F(x) to define. On the other hand (O/Onv)F(x- .) is regarded as an Ll(0t2)-valued continuous function of x E 0t2. Letting
Ho(~) = foa V(Y) ~---~vF(~- Y) ClSy for ~ E 0t2, we have 2H+ -- i V + 2H0
on 0~.
(5.65)
In the Fredholm theory, property (5.65) is applied to reduce (5.63) to an integral equation on the boundary. In fact, we have v(x)=-2foa#(y)o~v['(x-y)
dSy
for x E ~, with # satisfying
#(~) - 2 / o # ( y ) o ~ v F ( ~ - y) dS v = 9(~) for ~ C 0t2. Linear operator K :
C(OFt) ~ C(Ot2) defined by f
K,
(5.66)
=
0
jo ,(v)-~ r(. - v) dS v
is compact by Ascoli-Arz~la's theorem. Unique solvability of (5.66) then follows from the Fredholm alternative; in fact 1/2 is not an eigenvalue of K.
5. Other Methods in Space Discretization
192 However, the simple iteration scheme #n-I-1 -- g -~- 2 K p n
(n = 0, 1, 2 , ' ' " )
does not always converge because - 1 / 2 is now an eigenvalue of K. In particular, }--~n(2K)nr does not converge for the corresponding eigenfunction r - 1. To see this, let us apply Green's formula to -Av=0 We get
inf,,
Ov On =g
on Oft.
(5.67)
~(*)=fo (-~(v)~ r(*-y)+v(v)r(*-.~J))ds~
for x C f~. Therefore, from (5.65), it follows that
/o
o
/o
for ~ E 0f~. Because v = 1 is a solution to (5.67) for g = 0, this implies the assertion that - 1 / 2 is an eigenvalue of K. The boundary element method avoids the difficulty in the following way. First, one derives from (5.63) that v(x)=fo n
O--~F(z - Y) + ~---~v(y) F(x (-g(Y) On 9 _y ) ) dS~
(5.68)
for x E ~. Then, (5.65) gives 1
-~v(~) + f o n v ( y ) ~ F ( ~ -
~.])dS~ = fon J~v(y) . F(~- y) dS~.
If the integral equation ~o q(~)F(.(-:t j ) d S y =
~f(~)+
Jo.f(:~j)~F(~-y)dSy
(5.69)
in ~ C oqf~ is solved for q = (O/On)v, then formula (5.68) gives the solution v(x). This section is devoted to the wellposedness of the scheme; we show the unique solvability of (5.69). It was founded by several authors with the efficiency of the application of Ritz-Galerkin's or the collocation method. Let
and
A(p,q)= fo~.s ,, v(e)q(e)r(e- ,),~,s~,s,, The following theorem indicates that the bilinear form A( , ) is realized as a bounded coercive bilinear form on H-1/2(0f~) x H-l/'2(Of~). Then, the unique solvability of (5.69) is a direct consequence.
5.4. Boundaw Element Method Theorem
193
5.8. There exist constants C1, C2 > 0 determined by gt such that
c, Ilqll~H-1/2(0~"~) __< A(q, q) _< 02
Ilqll 2u_,,,~(o~)
(5.7o)
forqEX. Proof: Fourier transformation of a rapidly decreasing function f C S ( R ~) is given by .T" = f (~) = f
e-'Z~ f (x) dx,
n
where z = v/L~. Given q c X, we take Tq E D ' ( R ~) defined by (r
Tq) ---
fOa qr dS
for r E D ( R n) = C~~ Its support is contained in O~. It is also regarded as a compact measure on R ~ and Tq is taken as a C ~ function on R n. Now we show
A(v, q)= (Z~)-nJ~ 1___~(~)tq(~) d~
(5.7].)
for p, q E X, with the right-hand side converging absolutely. In fact, if n = 2, we have f'({) = p.v. ~
1
+ c(~({)
(5.72)
with a constant c. From the assumption, we have Tq(0) = 0. Regarding P C 8'(R2), we can define P Tq 9 C S'(R2). By (5.64), it is regarded as a continuous function on R2:
r 9T~(.) = f
JO
r ( ~ - y)q(y)
dSy .
We also have
f (r 9T~) = t . ~q =
p.v. ~
1
+ ~a(r
}
Tq(~)
t~(r 9 = ir
Let ~b E $ ( R 2) be real-valued. The inverse Fourier transformation is given as
~(~) = (2~)-~ a~f~ ~'~ ~(~) d~ and hence
f~ f0r(~-y)q(y)dS~r
= (r,
Tq,r (5.73)
5. Other Methods in Space Discretization
194 follows. We take p C C ~ ( R 2) satisfying 0_
supppc{lxl<_l},
and
[
JR
2
p dx = l.
Let p~(x) = e-2p(x/e) for e > 0. Given p c X, we put G = P~ * Tp. The left-hand side of (5.73) for r = ~ converges to (F Tq, 9 Tp) = A(p, q) as e I 0. To handle with the right-hand side, we shall show that 72 1 2
B~c~us~ Ir
IT,
I~1 9 -IT,I, th~n the right-hand sid~ of (5.73) convcr~r to
i~l= q(~)" ~-~p(~) d~ as e $ 0 by the dominated convergence theorem, and hence equality (5.71) follows. It is easy to see that (5.74) is reduced to the case p = q. Then we have
=
I_<~
I~l ~
d~ +
,>
1
1~
12
d~ = I + II.
The second term of the right-hand side is estimated as
II<_2
~1+1~12 d~-CllYqll2_,(~=).
Here, we have
[[Tqll/_/_I(IR2) -- sup { <2/2,Tq) I w ~ HI(I~2), [Iwtti/l(iR2) ~__1 } .-~ s u p { ~ n q v d S =
v E H'/2(Of~), ,.wl H,/2(on) <_ 1} (5.75)
]qlln-'/~(on)"
On the other hand, we have ]7~q(()- ]bq(0)] -< ].~a (e .... ~ - 1 ) q ( x ) d S x _< [e -'x~ - 1[ ,,,/2(oa)[IqllH-1/:(on) <- Cl~l Those relations imply
<-
ll ll "
Ilqll.-,/~(0n).
5.5. Charge Simulation Method
195
Because ] = I for n = 2, the right-hand side of (5.70) is proven. The left-hand side of (5.70) follows similarly from (5.75)as
I~
dg > ~f -
~ 1+1~
i5 dg ~
Ilq
I~-~/=(o~) 9
If n = 3, we have 1
f'(~) = ~-~ E L~oc(IR3).
(5.76)
It also holds t h a t [Tq(O)[- 1(1, Tq)l
~ O q [H_1/2(0~).
Noting those relations, we get (5.70) and (5.71) similarly. The proof is complete.
5.5
Charge
Simulation
Method
(CSM)
The m e t h o d in consideration is also concerned with (5.63) and make use of the fundamental solution F(x). We take N points y ~ , - . . , YN outside f~ and the solution v(x) is approximated by
N
(5.rr)
j=l
Here, the coefficients Q 1 , " " , QN are determined by
v N(xj) = f(xy)
(5.78)
for j = 1 , . . . , N , where X l , ' ' " ,XN are taken on 0f~. Sometimes {yl," "" ,YN} and { x l , . . - , X N } are called the charge points and the collocation points, respectively. In this section we study the case that gt is a two dimensional disc, ft = {x E R 2 I Ixl < P} with p > 0, and the charge and the collocation points are so located as yj = Re ~(j-1)~ and xj = pe ~(j-1)~ respectively, where z = x/~-l, R > p and 0 = 2rc/N.
Wellposedness Equation (5.78) is expressed as G Q = f , where G = (gij) with gij = F ( x i - yj), Q = T ( Q 1 , ' " , QN), and f = T ( f ( x l ) , . . . , f(XN)). In the above case, it follows t h a t 1
gij = -~--~ log p - Rez(J-i)~ I ~ Lj-i.
5. Other Methods in Space Discretization
196 Therefore, G is a cyclic m a t r i x given as
Lo LN-1
G
L1 Lo
9
9 9 LN-1 9 9 9 LN-2 9
.
L1
.
L2
..9
Lo
and hence :
p=O \ k=O
follows for w = e *~ L e t t i n g N-1
~;(z) = ~
~r
(z- P~),
k=O
we o b t a i n N-1
det G = I - I ~P(P)"
(5.79)
p=0
We can show the following. 5.9. It holds that
Lemma
- ~ log Iz~ - ~ 1 ~p (z) = I 1 N
for
Izl
1 ( ~z~lR)
(p =- 0 mod N) Iml
Z TM
(otherwise)
(5.80)
< R.
This implies cpp(p) r 0 for p = 1 , . . . , N holds by (5.79). Theorem
5.10.
1, in particular. Hence tile following t h e o r e m
The matrix G is non-singular if and only if RN
--
pN
?s 1,
(5.81)
in which case scheme (5.78) is uniquely solvable. Proof of Lemma 5.9: In the following, the equivalence class associated with m o d u l o N is d e n o t e d by "-". If p = 0, then w p = 1. In this case, we have N-1
qpp(z) -
1 27r Z
k=O
=
1
N-1
1 log Iz - Rwkl -- --~-~ log 17I Iz -- Rwkl
-2--; log I~~ - ~ 1 .
k=O
5.5. Charge Simulation M e t h o d
197
In the other cases, p = 1 , . . . , N -
r(~-R~
1, we have 1
~) =
-2---~ log ]z - Re,,,k ]
1 log ]Ra,,k] + log 1---R-2~1 (logR+Re ]og(l w - k z ) )
2~ 1
- -R-7-
log R -
Re
27r
-
n--1 n
Writing z = rP ~ we obtain
F (z-
R w k)
=
n=l n
1 27r
log R - Re
1 27r
log/7-E~nn n=l
1
w - n k e mO
~
R
(~,~o~-~
+
~-,~0~)
'
and hence
~(~) = - ~ 1 ~ ~
~1
log R - ~
k=0
~r
(e,,Ow_nk + e_~nOwnk )
n=l
follows. Because of N-1
Z J~
{ N (e - o)
=
k=0
0
(5.s2)
(otherwise)
we get
~(~)
1
_1
=
4~
_--
1 Z 47r
n=l
n>l
r
_1 n
p-n=_O
4 7r n GZ
n
-~
~
-~l
(eznOod(p_n) k +
e_mOcd(p+n)k )
k=0
r
n P n ~N + T~ 1
1 -n
r
e-'~~
n>l
(5),m, ,too e
p+n=O
,
m----p
and the proof is complete.
[]
5. Other Methods in Space Discretization
198 Convergence
Now we study the problem of convergence. Let f~ (n C Z) be Fourier coeMcients of f so that
f(x) = E fnem~ ncZ
I,fnl
holds for x = pe ~~E O~t. We suppose ~
+cr
<
(5.83) First, we note tim following.
nEZ
L e m m a 5.11.
We have vN(x) = ~ f,~n(p)-l~,~(x) n,EZ
(5.84)
.for x E f~ under the assumption of (5.81). Proof." We take 1 N-~
gij = "-~ ~
~)p(P) -102p(i-j)
p=0
and set d = (gij). Then, the (i, j) component of G G is given as N
E
1
Lk-i-N
k=l
~gP(P)-lcup(k-J)
=
p=0
=
1 N-1
N
-- E ~gP(D)-I~ N p=0 1
N-1
N E cdP(i-J) : 5ij p=O
by (5.82). We have (~ = G -1 and hence
D(I
QJ
=
~=
=
ff
~ k=l
1/
N ~=0
~?p(p)-la2p(j-1) p=0
~-p(k-i)f
(t00)k-1)
k=l
follows from Q = G - i f . Plugging this into (5.77), we get
~p(p)-lcoP(J-1) ~"~ (.x.)-P(k-1)f (/9(,,0k-l)
uN(x) = ~ -~ j=l
p=0
1 ~
=
N
N
j=l
k=l
N z_., ~pp(p)-l ~--~ ~p0-,) F ( x - R J - ' ) E p=O
1 g-1
--" N E
p=O
r (x-/~(.,o j-l)
k=l
~-pef (owe)~p(X)
(,pp(p)-i e=O
~-p(k-i)f
(pcdk-1)
5.5. Charge Simulation Method
199
by the definition of ~p(z). This implies
(;)
:
vN (x)
p,j =0
1 N
N-1
N-1
E (~gP(R)-I(cgP(X)E fn E CO(n-p)J p=0
nEZ
j =0
N-1
: E ~gP(P)-199p(X)E fn
(5.85)
n=p
p:0
by (5.83) and (5.82). The function Cp(Z) is harmonic and hence
I~p(p) -l~p(x) l _~ 1 holds. The right-hand side of (5.85) converges absolutely and we obtain (5.11) by (~gq= (~gp for q _= p. The proof is complete. [] Equality (5.83) implies r ) bl znO e nEZ
for x = rP ~ C Ft. We obtain eznO
~(x) } ~(p)
The function eN(X) is also harmonic and hence
II~N .o(~) _
max
Ixl=p
follows with an=
I~,(~)1 ~ ~ IAI ~
(5.86)
nEZ
I
sup Iezn~ oe[o,2~-) I
~,~(pe ~e) ~n(p)
Now we show the following. L e m m a 5.12. We have an = a _ n < m a x { 2 ,
1+
log log
(j~N __ fiN)
(587)
for n E Z and
an ~_
.for N >> 1, respectively.
NllogRI 8n ( p ) N - 2 n N-n-R
(n=0) (5.88) (n= 1,
,N-~)
5. Other Methods in Space Discretization
200
Proof:
T h e relation ~n(z) = ~_n(z) implies a - n = N
1
an.
If n ~ O, we have
Iml
m~n
and hence / an _< sup | 1 + oe[o,2~-) \
~..
(pe~~
<2
~,,(p)
follows. If n - 0, we have 1 I~,~(p~'~ <- ~I log lp~e'N~ RN] I < ~1 flog (pN + RN) I.
This implies log
(R N -k- pN)
log
(~N
an_<1+
__
pN)
and inequality (5.87) is proven. To show (5.88), we note
a~ _< I~.(p)1-1 sup ]e'~~
~ (pe'~
oe[o,,)
If n = 1 , . . . , N, we have N
~lml
N~ ~ 47rn '
where/3 =
p/R.
We also have
l~,,,o~,,(p) _ ~,, (p~,o) [
N =
_
7;7
(~ .. ~ -
Here, the r i g h t - h a n d side is e s t i m a t e d from above by N 47r
S ~§
~{
(~ ....0
_ ~,(~+.)0) + 9~N-" (e ....o _ ~_ ,(~N-,~)o} )
e=l gN + n < N k(flgN+n --
2rre= 1
=
27r(N-n)
gN+n +gN-n =
N/3N-n -< 2 ~ - ( N - n )
gN- n fleN-n)
g N + n t3(e~2n
(1+
oo
)Z/3(eg=l
+g.N-n
1)N __ N[3N-n - 27r(N-n)
1 +/3 2n 1 - / 3 N"
5.5. Charge Simulation Method Because
/~N <
201
1/2 for N large, then we get 2N
~(N-~)
~N-n "
Inequality (5.88) for the cases of n = 1, 2 , . . - , N - 1 is proven in this way. The case n = 0 follows similarly in use of the first equality of (5.80). The proof is complete. [] We are ready to prove the results on convergence. In the first theorem, the b o u n d a r y value f is supposed to be real-analytic. Therefore, the solution v of (5.63) has a harmonic extension to {x C R~ I Ixl < to}, where r0 > p. 5.13. Suppose R r 1 and (5.81). Then, there exists a constant C = C~,n,ro >
Theorem
0 such that
I vN -,11~(~)
-~ c,-~ sup I~(~)1, Ixl =To
wh ere
~-=
{
~T;;o (~o < R ~)
plR
(pro > t~~)
so that r C (0, 1). If R = 1, then v N(O) = 0 always holds, regardress with the values yy and Qj. Therefore, assumption R 7~ 1 cannot be eliminated. The following theorem deals with the case of more rough data. 5.14. Under the assumptions R 5r 1 and (5.81), the following results holds, where fn (n C Z) denote the Fourier coefficients of f indicated by (5.83):
Theorem
(i) If (ii)
~-~n IAI < §
If fn =
O
t h ~ limN_o~ II~ - ~11~(~) = 0.
(1~1-~),
then
I ~ - ~11~(~) = O (N-~+'),
wh~
~ > 1.
Proof of Theorem 5.13: The function f is supposed to be real-analytic, and therefore, Cauchy's inequality indicated as
p) i~i sup ~(~)1
IAI _< ~0
(5.s9)
~,=~o
holds. Letting m = [N/2], we split the right-hand side of (5.86) as
I1~11,_~(~)
m
_~ ISola0 +
~( n=l
SI -~- S2--~- S 3 .
oo
AI + IS-~l)an
+
E n=rn+l
(IAI + If-hi)an
(5.90)
5. Other Methods in Space Discretization
202
T h e first a n d t h e s e c o n d t e r m s of t h e r i g h t - h a n d side are e s t i m a t e d by (5.89) a n d (5.88), w h e r e a = R ~/(pro): $1
$2
8
<
( P ) N sup Iv(x)l ,
-
NJlogRI
<
Z
--
n=l
2
R
Ixi=ro
P
9
8n
N -
sup
16 ~
Ixl=ro
p N-2n n
v(z)l
sup Iv(~)l Ixl=ro
ct n n=l
N
<
ol- 1 \ro/ 16 1 -- oz
Izl=~o
--(P)x sup I~(~)1 ~
(~
Izl=ro
< 1).
O n t h e o t h e r h a n d , t h e last t e r m is e s t i m a t e d by (5.87):
s3 _< 2M~ ~
n=m+l
=
TO
sup I~(*)1
Ixl=ro
2Mn <~00) P m+l 1 1 - -~ sup I~(~)1 ro Izl=ro
<
-
1 - -p-
sup
ro
Ixl=ro
v(x)l.
If a > 1, we have p2/R2 < p/ro a n d hence
N/2
II<~ll~o~(r~) < _
~ c
~ --'
Vo
+o
s
TO
+ c
p
sup Iv(*)l
Ixl=ro
sup Iv(~,)l
Izl=,'o
follows. If a < 1, on t h e o t h e r h a n d , we have p2/R2 > p/ro a n d hence o b t a i n
Tile p r o o f is c o m p l e t e .
Proof of Theorem 5.14: We have I.f-I - IlfllL~o(o~)-
[] T a k i n g m = [N/3], we
estimate
Commentaly to Chapter 5
203
each term of the right-hand side of (5.90) as follows: S1 < -
(fl)N
8 NllogR
R
s~ _
S3_<(I+MR) ~
Ilflc~(Oa),
(I/~l+fAn).
n=[N/3]
The last term is o(1) and O ( N -a+l) according to the first and the second assumptions on fn, respectively, and the proof is complete. []
Commentary to Chapter 5 5.1. An important feature of hyperbolic or dispersive systems is the law of conservation. M. Mori and others tried to realize them in discretized schemes. The precise definition of the barycentric domain is as follows. Let T be the twodimensional simplex with the vertices P1, P2, P3. Then, its barycentric coordinate
{A(x", T, Pi)} 3i = l of a point P = P(x) C R 2 is defined in the following way: Put Pi = (x~, x~, x~) and set
(~) 1 A(x;T,P~)=~det
1 1 x I xl
1 x~
for x = (Xl, X2) and A = det
1 1 1 x I x21 x31 x I x~ x~
Now, given a E Zh [.J Bh, we put Tta = {T E 7 h l a ~ T} and
rer~ Then, the barycentric domain Da with respect to a is defined as
Da = U {z E T I A(x;T,b)_< A(x;T,a) TeT~
(b E I~ \ {a})}.
5. Other Methods in Space Discretization
204
T. Ushijima's works on scheme (5.1) were published in Ushijima [396] and [395]. RieszThorin's interpolation theorem is described in [127] and so forth. As is described in the commentary to Chapter 1, see Brezis and Strauss [65] for the L ~ estimate of elliptic operators. For the fact that A -~ : L2(t2) ~ H(~(t2) is compact if t2 has no corner with angle 2re, see Ne~as [286] and Grisvard [164]. Results obtained in this section keeps to hold for the general dimensional case under the following condition for the triangulation. Given T E rh, let V(T) be the set of vertices of T. Furthermore, for a E V(T), let Ta be the n - 1 dimensional simplex (i.e. face) of which vertices are composed of the elements in V(T) \ {a}. Write f(a,T) for the line containing a and perpendicular to Ta. Then the condition is described as
If n - 1, this assumption always holds. If n - 2, it is equivalent to saying that any T E % is an acute or right triangle. Here, n denotes the space dimension of the domain in consideration. See w for related arguments. Relation (5.15) is obvious because dim Xh < +o~. But it provides a key identity for E. Hille to construct (C0) semigroups on Banach spaces. K. Yosida made use of another aI)proximation for t h a t purpose, now referred to as the Yosida approximation. Some variants are made use of in w167 and 6.3. See also Hille and Phillips [179], Yosida [410], Kato [205], and Tanabe [378]. 5.2. The properties of L ~ stability and convergence were established by Tabata [373], [374] for scheme (5.29) with # (Eh, Xh) added to the left-hand side and without lumping in the right-hand side , where # -> #o = C(7)Ilbll 2L ~ ( n ) " Here, C(7) > 0 is a constant determined by the parameter 7 > 0 arising in connection with the regularity of {%}. See also [375]. 5.3. The inverse operator of the canonical injection a : V' --, V is sometimes referred to as the duality map. Note that the bounded linear operator A : V --, V is different from that of the previous sections defined through the Gel'fand triple V C H C V'. Variational problems of saddle point type arise in many areas. See Aubin and Ekeland [15], Rabinowitz [320], and Suzuki [369]. In connection with H(div, f~), the function space H(rot, f~) is also introduced when f~ C IR3. It is composed of the vector field 'v E L2(t2, R 3) satisfying V x v E L2(t2, Ra). See Girault and Raviart [160] for those spaces and applications. The existence of v satisfying (5.54) is also proven there. See also L e m m m a 7.32 in w Instead of (5.52) with (5.53), the Stokes system may be formulated as to find (u, p) E H(~(f~) n x L2(t2) satisfying 2~
(e~j('u), e,j(v))L~(n ) - (p, V v)L:(n 9 ) = (.f, v),:(n),~
i,j
for ~ c H ] ( ~ ) n ~nd ( V . ~,, q)~(~) = 0 for q c L : ( ~ ) , w h ~
1 (O'u ~ Ou~) r
= ~ \o.~, + ~
Colnmelatary to Chapter 5
205
This way is more natural from the physical point of view. Then, Korn's inequality takes place for Poincar~'s one to confirm the assumptions of Theorem 5.6. The continuous Stokes system and its finite element discretization were studied by many people, including Cattabriga [70], Kellogg and Osborn [217], and Crouzeix and Raviart [101], Bercovier and Pironneau [30], Glowinski and Pironneau [163], LeTallec [237], Stenberg [357], respectively. Mixed finite elements for magnetostatic and electrostatic problems were studied by Kikuchi [219], [220], [221], [222], [223]. Concerning the scalar equation, the Stokes system, and the equation of linear elasticity, see Brezzi, Douglas, Jr., and Marini [67], Raviart and Thomas [325], Crouzeix and Falk [99], Crouzeix and Raviart [101], Falk [122], and the references therein. Inequality (5.56) is referred to as Babu~ka-Brezzi-Kikuchi's inf-sup condition. It was introduced by Babu~ka [18], Kikuchi [218], and Brezzi [66] independently. Concerning actual examples of Vh and Wh satisfying the assumptions of Theorem 5.7 in use of finite elements, see Giraut and Raviart [160] and Brezzi and Fortin [68]. Under this condition, the discretized nonstationary Stokes system is treated similarly by the semigroup theory. See Okamoto [307] for details. The Navier-Stokes system is the fundamental equation in fluid mechanics. Operator theoretical approach for nonstationary problems was initiated by Kato and Fujita [207] and Fujita and Kato [143]. For other references, see the commentary to w Its finite element discretization was studied by Okamoto [308] and Heywood and Rannacher [175], independently. Error estimates given by the former are a priori; they hold only under reasonable assumptions on the initial value. The latter's are, on the other hand, a posteriori so that are valid under some additional estimates of the solution. Up to now, fundamental theory for this system is unsatisfactory, especially for the case of three space dimensions. A posteriori approaches are intended to compensate such a situation by numerical computations. See also Bernardi and Raugel [34], Heywood and Rannacher [176], [177], [178], and Rannacher [322]. 5.4. For layer potentials, see Courant and Hilbert [94], Kellog [216], and Garabedian [155]. Theorem 5.8 was proven by LeRoux [233] for n = 2 and Nedelec and Planchard [288] for n = 3. We have followed the method of Okamoto [309] for the proof, where q E X is taken in L2(cgf~). In the actual computation, the method of collocation is more realistic than that of Galerkin's. Iso [191] and Hayakawa and Iso [168] made error analysis for this method applied to the Neumann boundary value problem. Related works were done by Arnold and Wendland [13], [14]. Several monographs and surveys were published concerning the boundary element method. See Sloan [350] and the references therein. For the theory of distributions, particularly equalities (5.72) and (5.76), see Yosida [410] and Schwartz [348]. 5.5 CSM is also called the fundamental solutions method. Since the pioneering work by Bogomolny [41], the efficiency of the scheme has been clarified rigorously. Theorems 5.10 and 5.13 were obtained by Katsurada and Okamoto [213]. Theorem 5.14 was proven by Katsurada [210]. For other charge and collocation points, it can happen that the
206
5. O t h e r M e t h o d s in Space Discretization
coefficient matrix G becomes singular. Some examples are given in [210]. Problem (5.63) is invariant under the following transformation: (1) x H c~x, yj H c~yj, where c~ is a constant (2) f (x) ~ f (x) + c and v ( x ) ~ v ( x ) + c, where c is a constant. Scheme (5.77) has lack of such properties. K. Murota proposed N
j=l
where Q o , ' " ,
Q N are determined by VN(Xj) = f ( x j )
(j=
1,...,N)
and N
k=l
This scheme is provided with such properties. In this case, tile conditions R g - oN ~ 1 and R ~ 1 are not necessary for the wellposedness and the error estimate to hold. Katsurada [211], [212] studied the case that c%2 is a Jordan curve in use of the conformal transformation. An application to the free boundary problem was made by Shoji [351]. Furthermore CSM is a powerful method to obtain numerical conformal mappings. See, for this topic, Amano [6], [7], and a m a n o et. al [8].