- Email: [email protected]
On mixed methods for fourth-order problems
COMPUTER METHODS IN APPLIED @ NORTH-HOLLAND PUBLISHING
MECHANICS COMPANY
AND ENGINEERING
24 (1980) 13-34
ON MIXED METHODS FOR FOURTH-ORDER PROBLEMS A. QUARTERONI Zstitutodi Analisi Numerica de1 C.N.R., Pavia, Italy
Revised
A mixed finite element conditions is described. The problems. The computation methods for block matrices: couple of problems for the Some numerical results for
Received 24 September 1979 manuscript received 15 October
1979
method for the problem v + aZA2v =x with different types of boundary method converges, and it is well suited for the analysis of various evolution of the discrete solution is made by applying a sequence of iterative correspondingly to each iteration either a couple of Poisson problems or a identity operator are solved, according to the value of the parameter (T. two model examples are presented.
1. Introduction
The following problem is considered: 2,
E
H*(LQ, (I + a*A*)u = x
in 0,
(1.1)
where 0 is an open bounded domain of fR*, I and A* are respectively the identity and the biharmonic operators, (T is a positive number, and x is a given function of H-‘(O) (for the Hilbert and Sobolev spaces used in the paper we refer to [7]). The boundary conditions for u are given by either u
=*=o an
over r
(1.2)
over r,
(1.3)
or u=Au=O where r is the boundary operator. Setting
of 0, n is the outward normal direction to r, and A is the Laplace
(1.4) and integrating
by parts, we have (formally)
(1.5)
14
A. Quarteroni, On mixed methods for fourth-order problems
Let 4 be the function of HA(n) satisfying: -vAq
=x. Then with
(1.6)
w=uAu+q
it follows that u and w satisfy
(24w}E Hh(fq x
w,
=0
v+uAw
in g’(n),
w-uAv=q
(1.7)
in g”(0),
where two cases are considered: (i) if the boundary &/an = 0 comes (ii) if the boundary follows that w =
conditions (1.2) are given, then W = If’&?). (the second condition out to be a natural boundary condition); conditions (1.3) are prescribed, then W = Hb(R) since from (1.6) it 0 over r.
The aim of the paper is to analyze finite element approximations (1.7), namely:
{u,4
E
~b(~) x
of problems which generalize
w
v+uAw=p
in g’(0),
w-aAv=q
in B’(0),
(1 .S)
where p and q are two given functions of L2(fl). Remark that the problem (1.1) can be written in infinitely many ways in the form (1.8) since there are infinitely many possibilities of writing x as p - aAq with p E L2(L?) and q E H’(R). On the other hand, systems like (1.8) arise from some time-dependent boundary value problems such as the vibrations of an elastic plate or the Schrodinger equation (cf. [3]). In the quoted examples boundary conditions like (1.2) and (1.3) apply respectively, while (T behaves like the discretization time step At. Using the Green formula (1.5) the weak variational formulation of (1.8) becomes:
{u,WIE Hll(fq v4 EHX(fq
x
w (1.9)
VIJEW Hereafter (e, *) denotes the inner product of L2(R). A displacement finite element approximation of (1.9) is analyzed in the paper; clearly such an approximation also defines a mixed-type finite element scheme for the problem (1.1). Convergence results, and asymptotic error estimates in the energy norm and in the L2-norm are presented. In section 2 the case W = Hi is considered; a lot of simplifications apply since in this case the problem (1.9) is coercive. The harder situation related to the case W = H’(fl) is analyzed in section 3, while
A. Qua~er5ni, On mixed methods forf5u~~-order problems
IS
section 4 is devoted to discussion about the computational aspects of the problem. If W = H~(fIln),an iterative method for block matrices is used to solve the approximate problem. Correspondingly to each step either a couple of Poisson problems or a couple of problems for the identity operator are solved according to the value of the parameter (r. When W = H’(k!), i.e. when boundary conditions are prescribed for zt and its normal derivative, generalizing to our case the techniques developed by Ciarlet and Glowinski and by Glowinski and Pironneau (cf. [2], [4], [5]), then it is possible to reduce the approximate problem to a problem quite similar to the one corresponding to the case W = Hk(ft). Some numerical results for a model problem are also presented. An application of these results to time-dependent problem can be found in [8]. For any real number s the symbol ]].]lsdenotes the norm of I$“&?), I.1denotes the norm of L2(J2). Moreover C denotes a generic constant independent of the mesh size.
2. Analysis of problem (1.9) when W = ~~~(~) In this section the notation foflows:
A is continuous;
moreover,
v = HA(O) x HA(n) is used. Define the bilinear form A over v as
it is coercive since
As a matter of fact, from the Poincare inequality it follows that
Note that since the test functions written as
belong to the same space, the probiem
(1.9) can also be
where
and (,> denotes the duality between H-‘(0) and ~~(~). From the Lax-Milgram lemma the problem (2.3) admits a unique solution for any p, 4 E H-‘(a). Moreover, from (2.2) and (2.3) the following a priori estimate holds:
(2.5)
A. Quarteroni, On mixed methods for fourth-order problems
16
If the data
p, q
belong to L2(R), then with 4 = -V and $ = w from (2.1) it follows that
(2.6)
If-4+ I4 5 IPI + Id
Moreover, using (1.5) in (2.1) one obtains UAW = p - v and uAu = w - q; so the estimate (2.6) leads to 11412 + Ml2 5 CaPI
+
ISI).
(2.7)
2.1. Discretization of (2.3) and error bounds Assume that 0 is a convex polygonal decompositions of J2 into polygons T. Then
6=
h = m;x(diam
UT
TET,,
domain
and let {T,,} be a family of regular
T).
(2.8)
For any integer k zz 0 let P,(T) denote the set of polynomials of degree not greater than k over T, and define
VI, = (4 E
c"(fi))($IT
l&‘-)
VT E W,
(2.9)
& = v, f-UG(R),
(2.10)
v, = & x;,.
(2.11)
The following approximation
vv E I&l(n)
flH’(f2)
VW E H’(0)
properties
are fulfilled by the above spaces (cf. e.g. [l]):
Inf ]]V- 4(]i 5 Ch’-i]]V]]r, 4& I$Jh llw - tilli 9 Ch’-illWJl~,
for i=O, 1 and irr~k+l. The following problem is a finite element displacement
By the same arguments be proved:
used for the continuous
approximation
(2.12) (2.13)
of (2.3):
problem the following stability properties
can
l~hl+lwhl4Pl+l~l
(2.15)
]]Vh]]l + IIWhlll zs $=g (IlPII-1 + 1l~ll-d~
(2.16)
17
A. Quarteroni, On mixed methods for fourth-order problems
As far as the error bound is concerned, the following result applies: THEOREM 2.1. Assume that the solution V, w of problem (2.3) belongs to the space H’+‘(fi). Then the following asymptotic convergence estimates hold:
llu- Uh((l + Ilw- wz(II5
$
h"(llullk+l
]?J- &,I + ]w - wJ,(-= C63’2hk+‘(llv((k+,
+IIWllk+dr
(2.17)
+ IIWllk+,)
(2.18)
Proof. From (2.3) and (2.14)
A({Q-
vh,
W -
WI,},
bf’,
$1)
=
0
%#,
+I) E
(2.19)
vh.
Then (2.17) follows from (2.12) and (2.19) using standard arguments (cf. e.g. [l]). An application of the Aubin-Nitsche trick to problems (2.3) and (2.14) leads to: lu - Uhl+ (w - WI 5 Cb-‘((Iv
- uhlll+ llw - %ll1),
and so the estimate (2.18) follows thanks to (2.17).
cl
3. Analysis of (1.9) for W = H’(R) In this section problem (1.9) with W = H’(LI) is considered. This time the coupling (2.1) for the equations of (1.1) cannot be applied. Let us define therefore a new bilinear and continuous form over the space w = H:(R) x H’(O):
B(b, WI,{A 44)= (VT4) + 6%@I) - 4w
4) + aa(v $ci).
(3.1)
The problem (1.1) can be written as follows: {II, w} E W
WV w], I+, $]) = ](P, q]> (4, +]lo
W,
$1 E w,
(3.2)
where HP, 4], {+Ytw”
= (P,d4+ (4,#cl>.
(3.3)
Remark 3.1. Let (Ybe a positive integer and consider
The operator
I + aA* is an isomorphism
between Hi(O) nH3(0)
and H-‘(Q);
1144135 WI-1
for a positive constant C depending
moreover: (3.5)
on Q (cf. [7]). Also, following the arguments
of [6], there
18
A. Quarteroni, On mixed methods for fourth-order problems
exists a real number v E IO, l] depending IId43+”
5
ca
on the shape of 0 such that
IfI,
(3.6)
assuming f E L*(0), The estimates (3.5) and (3.6) (for some Y 2 0) still hold if 4 satisfies other boundary conditions, for instance if +=a~$=0
over f.
or A4
=@!b=o
•1
over f.
an
THEOREM 3.1. For any u > 0 and for any p, q EL’(LJ) problem solution. Proof. (A) Uniqueness. Setting 4 = z, and t,,Q = W, then from (3.2)
I4 +
I4 5 IPI+ 141.
(3.2) has a unique
(3.7)
Let us emphasize that the bound (3.7) does not depend on u. If p = q = 0, then from (3.7) one obtains u = w = 0. (I?) Existence. Consider the problem (I+u*A*)c#J =p
in 0,
+=0
over r.
Setting u1 = 4 and w1 = aA+, then from the arguments such that {v,, wl} is the unique solution of ~1E H&Ll) n H3+“(0)),
WIE H’+“(0),
v,+aAwl=p
in L*(o)),
w1-uAv,
=0
(3.8)
of remark 3.1 there exists v E IO, l]
(3.9)
in L’(0).
Consider now (I + (~*A*)t,b= q
in 0,
A$=?=
0
Setting wz = $, uz = -ad+, then by the same arguments such that {v*, w2} is the unique solution of 02 E H:(fl) nH’+*((n),
w2 E H3+” (ill),
u2+aAw;l=O
in L*(a)),
wz-oAvz=q
in L*(a).
over f.
(3.10)
it follows that there exists p E [0, l]
(3.11)
A. Quarteroni, On mixed methods for fourth-order problems Now
19
let v = v1 + v2 and w = w1 + w2; then problems (3.9) and (3.11) lead to v
E
H;(n)
v +aAw
f-W+*(n),
w E H’+*(n), (3.12)
=p,
w-aAv=q for a suitable ,$ E [0, 11. Then {v, w} E w and solves (3.2). q THEOREM 3.2. Let .$ E [0, l] be the parameter appearing in (3.12) and let {v, w} be the solution of (3.2). There exists a positive constant C independent of u such that
llvllr+ lb4 5
3
(IPI+ Id),
Proof. The Green formula for the biharmonic (A2u,w)=(Av,Aw)+
(3.13)
r E [O, 1+ 61.
operator gives (formally)
($gw)-((A%$$
(3.14)
where (-, -) and ((e, -)) indicate dualities over r in some suitable spaces. Let us refer to the constructive proof given in step (B) of the previous theorem. From (3.8) and (3.14)
Id’+ a2b412 = (P,4),
(3.15)
(At$12+ u~IA~~#J(~ = (p, A2$).
(3.16)
The equality (3.15) gives
141~IPL
bPl +PI.
(3.17)
From the second inequality of (3.17) one obtains, for a constant C independent
of (T: (3.18)
Moreover,
from (3.16) and (3.17)
IA’41~$2 IPI. Interpolating
(3.19)
2 and I-1(see [7]), then from (3.18) and (3.19) the norms I]*]]_ (3.20)
llA2c#+~+s5 Ca-‘3+6”2jpI. Still using an interpolation
theorem,
from (3.19) and (3.20) the following estimates
can be
20
A. Quarteroni, On mixed methods for fourth-order problems
obtained:
IlVllls= 1144 5 c~-“‘*lPl~ llw,lls-2= llv,ll, +
~~d#4-2
5
&lis
5
C~1-s’2bi,
s E 1%3 + !$I,
(3.21)
51,
(3.22)
r E [O,1+ 61.
(3.23)
s
Ilwllr5 cfr’21PI>
E
+
[0,3
In a similar way, starting from (3.10) and using (3.14) and the interpolation it follows that IIw*lls
llu2l(s-2
(12)2(lr
11~11s 5 c~-%l7
=
=
+
()d+I(,-2
5
~lb~ls
5
C~‘ps’2~ql,
Ilwzllr 22c~-“*ld
theorems (see [7]),
s E [O,3 + SIT 3 + 51,
(3.25)
r E [O,1+
(3.26)
s
E
(3.24)
[O,
61.
Since v = ul + v2 and w = wl + w2, the result (3.13) applies thanks to (3.23) (3.26) and the triangular inequality. cl 3.1. Discretization of (3.2) and error bounds The,,notations used in subsection 2.1 are preserved hereafter. Moreover, we set S,, = V,, 6’h = S,, x S,, and assume in addition to (2.8) that there exists a positive constant 6 independent of h such that h 5 6(diam T)
A decomposition
VT ETh.
(3.27)
{T,,} satisfying (3.27) is said quasiuniform.
For such a decomposition
the
inverse inequality applies, i.e.
for a suitable constant C (cf. e.g. [l]). M oreover, inf-sup conditions hold:
it can be easily verified that the following
where p is a positive constant independent of h. Consider the following displacement finite element approximation
of (3.2):
21
A. Quarteroni, On mixed methods for fourth-order problems
Setting 4 = uh and J, = wh by (3.31) one gets
If p = q = 0,then &, = wh = 0; therefore, (3.31) admits a unique solution. THEOREM 3.3.Assume a = 1, and let {v, w} and (uh, wh} be respectively the corresponding solutions of (3.2) and (3.31). Assuming 21 and w regular, the following asymptotic error estimate holds: Ilu - t&II, + Iw - w,,l 5 Chk-*>
(3.33)
where the constant C depends on Ijullk+land Ilw&+1. Proof. From (3.2) and (3.31) it follows that (V-Uh,+&+wh,&)=O
v4
E
ih,
(3.34a)
(W-Wh,$)+a@-&,,$)=O
vt,b
E
sh.
(3.34b)
In the proof 4, $ shall denote generic functions of g,, and Sh, respectively. (1) EstimateOf IW
-
W,,(
~W-wh~2=(W-Wh,W-~)+(W--h,1(1--h)
(from (3.34b))
=(W-wh,W-$)--a(t?t?h,t,b-W)=-a(Z)-uh,W-Wh)
=(W-Wh,W-(CI)+a(~-~h,W-(Ir)+a(~-~,Wh-W)+~(~h-~,W-Wh),
(3.35) (i) bound for u(v -
uh,
w
-
$):
Then thanks to (3.30) there exists 6 E Sh such that U(U - u,,, w -
$) 5 Ilw- ~~~1{~~~ - d& + P-‘llill%~
= ((W - I,b~~~{~~TJ - (6111 + P-‘(a(C$5
-
Uhr
II;)}
u, 6) + u(u - uh, 6)))
IlW- $ll,{llu- 4111 + p-‘(llu - +I\1+ lw- Whl)}
(from (3.34b)), (3.36)
(ii) bound for u(v - 4, wh - w): U(U - 4, w,, - w) 5 5
Ilu- $lll(llw- $111 + llwh- @IId c((U- ~~~,{~~~ - #(I,+ h-‘(lwh- w1+ lw-
$1))
(from (3.28)), (3.37)
A. Quarteroni, On mixed methods for fourth-order problems
22
(iii) bound for a(~$, - f#h W -
wh):
(from (3.34a))
&h-$,W-Wh)=-(u-Uh,+Uh)
= -(ti - v,,, v - vh) - (0 - vh, (a - v) 5 -(v - v,,, 4 - 0)s (v - v,,[lv - +[ lv-~h(+-~I+I~-%t
%
Iv
-
41
+
P-‘i&‘(d6
-
2)~ +)+
a(u
-
oh,
5 ~-~I+P~‘(ll~-4lll+Iw-wtl)
(from (3.30))
9))
(from (3.34b)).
Hence (3.38) Substituting
(3.36), (3.37), (3.38) into (3.35) it follows that
Iw - &I2
d~‘<~~w - till’ + h-‘(l&t’ - whl+ Iw - $I))+ Iv - 41<1v - $1+/Iv- dl + Iw- Whl)h (3.39) 5
c{lw - wh((w- +I + IIw- $lll(llu- +(I’+ lw - whl>+lb -
(2) Estimate of 11~1 - uhj)l Using (3.30) once more, one gets I/v - Vhlll
5
/Iv - (j$lj’+
114- vh)),s I/v- 41)’+~-‘tl&‘“(~
-
vh’
‘I’
From (3.34b) I(V -
Z&II,
5 [IV - &‘([I +@-‘(I&‘(@ 5 (1 +p-‘)llV
- VT i)+
- +I(, +p-‘lw
-
whj
(wh - w, 6)) 5
c(((u
-
411’
+
iw
-
wd)’
(3.40)
Going back to (3.39) and using the results (2.12) (2.13) (w - wh125 c{lw - wh(hk+’+ hk(hk + (w +
hk(hk
+
whl)
hk + h-‘/w - w,,l)+ hk+‘(hk+’ + hk + (w -
(note that q5 and rj~are arbitrary functions of &, and Sh respectively). IW -
w,,l~ C’hkp’,
cl
=
C’(llv(lk+‘,
Whl)j
Therefore
(IW((k+l).
(3.41)
Finally, thanks to (3.40) (2.12) and (3.41) the following estimate applies:
[Iv - v,,)I,5
C,hk-',
c2 =
c2(((v((k+‘,
((Wllk+l)*
(3.42)
A. Quarteroni,
On mixed methods
for fourth-orderproblems
23
Clearly (3.41) and (3.42) give (3.33). ~~~ORE~
0
3.4. Let B be an arbitrary positive number; if the assumptions
of theorem 3.3
apply, then Ilk>- Uh((,I C(hk-’ + Q -‘hk + (Z--“h’+’+ f_-“W),
(3.43)
Iw- WJ I C(hk + ,P’
(344)
+ a-~hk+’ + k%zk).
Proof. This is quite similar to the proof of the previous theorem.
q
3.2. L’-error bounds for linear elements When piecewise linear elements are used, i.e. if k = 1, the estimates (3,33), (3.43) and (3.44) do not provide an asymptotic convergence for 11~- r.+,((,and Iw - whl. The aim of this subsection is to provide an asymptotic error estimate in the L2-norms even for this case. Let us define the Ritz operators by ~,:H:,(.n)~~~:a(t:-i;ht’,~)=O
V& E &_
(3.45)
VJ, E Sh,
(3.46)
vu E H;,(n) n w(a)
lu _ ii,Vl + hll1,- &VII,5 chqlull,.
(3.47)
VW E H’(.n)
jw - R,wf + hllw - RhwI(, 5 C!z’llw&
(3.48)
Rh : H’(L+
Sk : a(w - RhW, J/) = 0
(w - R,,w, 1) = 0. The following results apply (see e.g. 191):
for 1 I r zz 2. Generalizing (see [lo]), we can prove
the techniques
developed
by Scholz for the biharmonic
problem
THEOREM 3.5. Let CT= 1, k = 1 and (u, w}, {uh, wh} be respectively the solutions of (3.2) and (3.31). Then for some constant C depending on I/r&, ~~A~~l~,~~~~~, Ilwli2the following result applies:
Iu- vhl+ h”2(ln
hllw - wh( 5 Chlln h12.
(3.49)
FTOOf. Setting el = v - v4 and e2 = w - wh, (3.34) gives (er.ct,)--a(e2,#)=0
V$ E &,
(3SOa)
tc,, 4)+ &I? $) = 0
vJ/ E Sk.
(3.50b)
24
A. ~uu~ero~i, On mined methods for fourth-order problems
(1) Estimate of le,j Let 5 = &,u,
q = &W ;
IWh- r112 = (Wb-
then (3.51)
w,wb--r))+(w-q,wh--v)
(from (3SOb))
(Wk- w, Wk -~)=-(e2,Wr*-r))=~~~i,Wb-r))
=a(V-StWh-rl)f~(Z‘-~k,Wk-77), a(~-Ub.Wh-rl)=~(5-~h,Wh-W) =
-f&,
r - %uh)
(3 S2)
(from (3.46)) (from {3.5Oa))
= -(er, e,)+ fez, 2,- 8)” (6, 0 - 5).
(3.53)
Let us note the result (cf. [lOI)
where C is a positive constant independent state
of U, b, and $J. An application
of (3.54) allows to
(3.55)
la(v - 6, wh - r))l s CIP211n hl~~AullL~~~~~w~ - r)l.
(3.56) (3.57)
(3.58)
Let 6 EW@?il) satisfy 6 -t-A26 = el in 0. Assuming estimate applies for a suitable Y E IO, l] (cf. remark 3.1):
0
regular,
the following
a priori
(3.59)
A. Qwarteroni, On mixed methods for fourth-order problems
25
Moreover, (3.60)
)e112= (0 + A*& el) = (0, el)-- u(A0, e,) = (0, et) - u(A0 - R,,AO, e,) - a(R,,AO, el) - a(RhAO, el) = (e2, R,Ae)
(from (3.50b)) (from (1.5))
= (e*, RhAe - AtI) - a(e2, e)
= (e2, RhAe - de) - a(e2, 8 - k,O) - (e,, k8)
(from (3.50a))
(3.61)
Applying (3.47), (3.48) and (3.54), the following estimates are deduced:
la(Ae - R,,Ae, e,)l = [a(Ae - R,,Ae, u - t)ls
C2h’+Y~l~~~2~~e~(3+v,
I(e2,R,,Ae - AtI)1 5 C3h1+Y~~8~~3+v~e2~, (u(e,, 8 - rZ,e)ls \U(W - 7),8 - &e)l
+ Iu(~ - wh,8 - iLe)l,
5 C4h211wl12(l1412 + IIA~llL~~nJ+ Gh1’211n 4 - llAf4hde21.
Since from Sobolev’s imbedding theorems H3+V(0) C W22m(L?), then IIAellL-cn,5 C611e((3+v, therefore, from (3.59), (3.60) (3.61) and the previous inequalities one obtains (e112 I C{h21e112+ h1+“((1412+ llwl12>leIl+ (II’+” + h”*lln
Wle2ll4l.
For h sufficiently small the term 1 - Ch is positive, then IelII C{h1'"(llu(12 +
Ilwl12)+ (h’+”
+ h1’2(ln hl)(e&.
(3.62)
Using (3.58) from (3.62) it follows that (ell 5 C{h”“(l(u(12+ jl~11~)+(h3’*+“(lnhj +
(h3’2\1nhi+
+ hlln h12]lAz&-~n~
h2’y)~fl}
(3.63)
From (3.63) it easily follows leI( I Chlln
hi2
for a positive constant depending by (3.64), it follows that
(3.64) Going back to the estimate (3.58) on Ilul12,IlAo IIL=(R)and (Jw(12.
[es1I Ch1’2(ln hj. Now the result (3.49) applies thanks to (3.64) and (3.65).
(3.65) 0
26
A. Quarteroni, On mixed methods for fourth-order problems
4. Numerical algorithms to solve (2.3) and (3.2) The aim of this section is to propose a method for solving problem (2.3) and a general method for reducing (3.2) to a problem whose numerical solution is quite similar to the one of (2.3). To do this, some notations are needed. Assume that the known functions p and 4 are accurately approximated by two piecewise polynomials ph and @, of Sh. For fixed k, let fi and N (fi< N) be the dimensions of Sh and Sh, respectively;
moreover
of N independent
let {&} be a set of N independent
functions of &,, and let {$j} be a set
functions of S,,, Set:
A = [A,],
Aij = (4iy4j);
c = [Cj],
Cij
D=
= ($iy $j);
toi,],
Dij
B = lBij 19 Bij
=
=
a($i,
4j),
a (41)$j).
4.1. Solution of problem (2.3) The discrete problem (2.3) is equivalent
to the linear system
;;u=f
(4.1)
where k is a symmetric, positive definite matrix defined by i;=
[
u;
-iD .
(4.2)
1
moreover, u =(%BY,
f =(PdI)
and 2
ai+i
=
vh,
2
Pi6
=
why
pi
=
(ph,
+i),
4i
=
(qh,
6i).
The system (4.1) can be solved by the following iterative scheme: For an arbitrarily fixed /!I” define the sequence cy” and p” by solving: aDa”“’ +A@" Aan+'
= q,
- uD@"+' = p.
Solving for the unknowns,
(4.3) n 20.
(4.4)
eqs. (4.3) and (4.4) become (4.5)
A. Quarteroni, On mixed methods for fourth-order problems
27
The scheme (4.5) defines a block Gauss-Seidel iterative procedure for approximating (Y and 6. A sufficient condition for convergence is that the maximum norm of the iteration matrix be less than 1; this is equivalent to require that the maximum norm of the matrix a-‘D-IA be less than 1. Note that by definition the elements of D-‘A have the upper bound Ch2 for some positive constant C; hence, setting C1 = C’, a practical condition for convergence is h2 <
C,a.
(4.6)
If the parameters h, u do not fulfill the condition (4.6), then instead of (4.3) (4.4) the following scheme can be used: For any given vector /3” compute the sequence a”, p” by solving
AP“+I
+ g&p+
zz
q,
Aa*+’ -uDjY’=p,
(4.7) It 20.
(4.8)
The explicit form of such a scheme is (4.9) By the same arguments used to achieve (4.6) it can be seen that a sufficient condition in order to get convergence of the iterative scheme (4.9) is a
c C2h2,
(4.10)
for a suitable constant C2. Let us remark that for evolution problems u behaves like the discretization time-step At; then for h sufficiently small one of the conditions (4.6), (4.10) holds except for At - h2. However, it should be emphasized that the best situation is that in (4.6) at least for implicit schemes. When At - h2 or (4.10) applies, then an explicit scheme seems more suitable for the computations. Remark 4.1. It can be verified that at any step n to solve (4.5) means to solve two Poisson problems: -uAv ;+l = W ;: - q,,,
(4.11)
-UAW;+’ = -v”h+’+Ph.
(4.12)
On the other side, to solve (4.9) means to solve the following two problems for the identity operator: V ;+’ W;+’
= =
-UAW “h+ ph,
(4.13)
uAv;+’ + q,,,
(4.14)
in both cases a starting function wi = Xip’l4l is needed.
28
A. ~~a?ie~on~, On mixed meiho~ for fourth-order problems
4.2. Solution of problem (3.31) Setting
u =bm,
f =cPd
and c y&i =
Uhr
2
&$j
=
wit,
pi
=
(ph,
(bill
qj
=
(qh,
$j),
i
the
problem (3.31) is equivalent
to the linear system (4.15)
Hu=f,
where the matrix II is given by (4.16)
H =
Since the matrix B is not square (%< N), systems like (4.15) cannot be solved separately as two sets of equations, and so the method proposed in the previous subsection for the system (4.1) does not apply. On the other hand H is not well featured (it is undefinite). In this section, following [5] closely, we introduce a method for solving (3.31) which avoids the direct solution of the system (4.15). Let (u, w} be the solution of (3.2) and let A be the (unknown) trace of AU over f. Clearly from (1.8) it follows that the solution of (3.2) satisfies v+uAw
=p
/ w-uAv=q
in 0, in fz,
1V=O iw=ai\+q
(4.17) over IY
In order to write the last equality over F, the trace of q over 1” must exist; that happens if, for instance, q E H’(f2) or even if q belongs to the space N(A, fz) that we are going to define. Let H(A, a) = {# E L*(O)]A& E L2(Q)], and for any p E H-“2(r) ‘,y+aAo
=O
let A (cc) = (x, o} E H*(J2)
H(A, L?) be the solution of
in 0,
w-crdx=O
in 0,
x=0
over r,
.o =a#.&
X
over IY
(4.18)
A. Quarteroni. On mixed merhods for fourth-order problems
Moreover,
let {v,,, w,,} be the solution of vu -t UAW,!= p
in 0,
wo - udvo = q
in L?,
V()= 0
over I:
WI& =q
over f.
(4.19)
Setting: (6. @}= A(A),
(4.20)
it follows from (4.17). . . . , (4.20) that v = VI,+ 6,
w =
w,,+ G.
(4.21)
Consider now the operator
/.L+u”~
,:H-“*(f)+Hii2(f);
over jT,
(4.22)
where by definition (x, o} = A @); E is an isomorphism since A is an isomorphism, and the application which maps 4 into its first-order trace ~~~~~ over r is an isomo~hism between HZ(a) and H”‘(f) (cf. e.g. f71). Finally set e : H -l/‘(f) x H’/2(f)-+
IR;
e(A, p) = (EA,p).
(4.23)
Here the duality is between H”‘(f) and H “‘(I-‘). The form e is continuous moreover by (4.20) and the de~nition of e the following relation holds:
e(A, p) = u2 -
Ir$$dy=u -
dy
(using the Green formula (3.14))
(u2A2,y + x, 6) + h, i2)-t-a’(Ax, AC) = (,y, a) + (w, @)
Hence, for any 4 = (4,, &) and $ = (#,, (cIz)belonging to L2(L!) x L’(R), [#%J/l= @I> $1) +
(42.
$2)
and bilinear;
(from (4.18)) (4.24) setting (4.25)
we get from (4.24) e(k/J)=[Nq,NP)l; hence, e is symmetric
(4.26)
and coercive over H ^ ‘l”(r). Since a~/a~ = 0 over f (cf. (1 .l)) from
A. Quarteroni, On mixed methods for fourth-order problems
30
(4.20), . . . , (4.22), it follows that the trace A of Au over T is the unique solution of: (4.27) Setting once more h, w} = A(p), from (4.18) (4.19) (4.27) and the Green formula (1.5) for the Laplace operator it follows that
I dvo
-CT2
I‘ an
/A
d y = -o@ut,, VW)- u(&,,
w) v/AE w-1’2(r),
= -a(Vuo, VW)- (wo, 0) -+(43 @)
and so (4.27) becomes
vp
E
Er-“yr).
(4.28)
Therefore, in order to obtain the datum h of (4.17) we must solve the variational problem (4.28) (which in turns requires, for the determination of the right-hand side, the solution of problem (4.19)). Passing to the numerical approximation of (u, w), let {oh, wh) be the solution of (3.2). Let A& be a finite-dimensional space contained in S,, such that Sk =&+MJ,.
(4.29)
We shall think of Mh as the complementary of $, into S,. Using Lagrange interpolation, any function of Mh can be characterized boundary nodes. Let Ah be the component of wh into Mb. Setting Wh
by its values at the
=/i,+(w,-rc,)
and following [2], the triplet {v,,, wh, h;r} can be characterized
(4.30) as the unique solution of
where for any 5 belonging to Mh we denote by Sh(e) the following subset of Sh:
Since qh E Sh let us set
A. Quarteroni, On mixed methods for fourth-order problems
The following problem is an approximation
31
of (4.19):
(4.35) Moreover,
for any j.&hE A& let (xh, ah}
=
be the solution of
Ah&)
and set
{i?k,'c3h) = Ah@k)v
(4.37)
where hh EM,, is defined by:
ik = CrAh+
(4.38)
ijk.
From (4.31), (4.35), (4.36) and (4.37) it follows (in analogy with (4.21)) that vh
=
f,+jk +
ek,
Wk =
WOk+
ek(hh,
pk)=
w-k-
(4.39)
Now, following (4.26), we set ek
:Mk
x Mh
+R
;
The form eh is bilinear, continuous follows that hk is the solution of eh(hk,
pk)=
-u(bOk,
[Ak@k),
Akbk)].
(4.40)
and coercive. From (4.32), (4.35), (4.36), (4.38) and (4.39) it
VWk)-(WOk7
ok)+
(% %)
(in analogy with (4.28)). Therefore, in order to compute the solution {t& wk} of (3.21, we must (i) compute the solution {&,h,WOk}of (4.35), (ii) solve the yariational problem (4.41), (iii) COmpUte hh by (4.38), (iv) solve the problem (4.31). Now let us make some remarks about steps (i), . . . , (iv). Step (i), (it,). We can easily check that the problem
(4.41)
32
A. Quarteroni, On mixed methods for fourth-order problems
(with the right-hand
side eventually zero) is equivalent
to the linear system
iiz
(4.43)
=b
whose matrix is defined by (4.2). Then (4.43) can be solved by the method proposed in subsection 4.1. The nonvanishing boundary condition for wh only changes the right-hand side of the system. Therefore, both problems (4.31) and (4.35) can be reduced to a linear system like (4.43). Step (ii). Let us denote by {pi}:, a set of basis functions for Mh (N is the number of boundary nodes). Setting
Eij = et,bi, hi,
E
mi}
=
=
[Eij],
Ah
bi
P,
1,
Ci = -(T(hOhy
), A
=
[A’],
it easily follows that (4.41) is equivalent
c
=
vwi)-
[c-i],
(WOh,
@i)+
(qh,
wi),
1 Ii,jlN,
to the linear system
EA = C.
(4.44)
The matrix E is symmetric and positive definite; so (4.44) can be solved by direct methods (for instance by the Cholesky method) or by iterative methods. (Using for instance the conjugate gradient method, an explicit construction of E can be avoided provided two Dirichlet problems for -A are solved at each iteration (see [5, section 5.41 for details). Moreover, E does not depend on the terms p and q but only on the geometry of 0. This is an important feature since for problems that require a recursive solution of (4.15) (this is the case of the evolution problems), the solution of (4.44) must be done only once. Since the functions U()h,W()h,qh and wi are piecewise polynomials, the right-hand side of (4.44) can be correctly computed without using quadrature formulas. Moreover, in order to compute xi and wi, we must solve the problem (4.36) with ,& = pi. That in turn implies that a systems like (4.43) must be solved (cf. the discussion about steps (i), (iv)), where b consists of the contributions given by the boundary values of wi. Such a system can be solved by the method proposed in subsection 4.1. In order to compute C, the solution of N systems like (4.43) is required. Step (iii). Once Ah is known, the computation of & only requires the computation of &,, 0 and that is trivial since, if Mh is the orthogonal space of sh with respect to S,, for the H’ inner product, we get (4.45) the summation being restricted to the boundary nodes (7,). Finally we can conclude that the computation of the solution {zlh,wh}
of (3.2) requires the
33
A. Quarteroni, On mixed methods for fourth -order problems
w(x,,x,f=oAv+i v(~&f=$D21x,)P”o?%~ t ? ! \
i i i ! ! \
i i i i i i
i i
_____ “___(I(J_1., \I rel.error
_--._
--___ i
\
i 0 ii+4
i i i
--------_.-.-,-
1
-----
! ! ! \
6
8
10
(J= 1. t re, ,?rror U=.l I . (i=
6_
1.
.,
\
j rel. error
* l/h
_.
u- 1. \ ~_ .11 rel.
error
for V for W
i i i i
i i i
tar V for W
i
i i
--r,, Fig. 1. Results for problem (2.3) with data in (4.46).
4
6
8
f0
c
l/h
Fig. 2. Results for problem (3.2) with data in (4.47).
solution of N + 2 systems like (4.43) (by the method proposed in the previous subsection) plus cl the solution of the system (4.44) whose matrix is well featured. Rebury 4.2. We refer to [ll) and [12] for the analysis of methods similar to the ones of section 4.2 related to steady and unsteady Stokes problems using the velocity-pressure formulation via also a mixed finite element method. q In figs. 1 and 2 we show some results related to problems (2.3) and (3.2) respectively, solved by the methods proposed in section 4. We assume 0 = (0,l) x (0, l), we consider decompositions of fi by squares and piecewise bilinear polynomials. Fig. 1 is related to a problem like (2.3) solved by the algorithm of section 4.1; the model problem corresponds to the following data: p(x,, x*) = #+M(xz)-
2~7M&)X(&),
q(xl, x*) = X(~I)X(=z) - 2G(#(G) f tfi(X2))T
(4.46)
34
A. Quarteroni, On mixed methods for fourth-order problems
i = 1,2. where 4 (Xi) = Xi(1 - Xi), x(x~) = sin TXi, Fig. 2 is related to a problem like (3.2) solved by the method proposed corresponding to the data
&I,
x2)= 1.
in section 4.2 and
(4.47)
The automatic procedures terminated when the distance between two successive iterates was less than 1O-4in the maximum norm. The relative errors are evaluated in the L2-discrete norm. References [l] P.G. Ciarlet, The finite element method for elliptic problems (North-Holland, 1978). [2] P.G. Ciarlet and R. Glowinski, Dual iterative techniques for solving a finite element approximation of the biharmonic equation, Comp. Meths. Appl. Mech. Eng. 5 (1975) 277-295. [3] J.E. Dendy, Jr., An alternating direction method for Schrodinger’s equation, SIAM J. Numer. Anal. 14 (1977) 1028-1032. [4] R. Glowinski and 0. Pironneau, Sur la resolution, via une approximation par elements finis mixtes du problemes de Dirichlet pour l’operateur biharmonique par une methode “quasi-directe” et diverses methodes iteratives, Rapport 197 (IRIA, Le Chesnay, 1976). [S] R. Glowinski and 0. Pironneau, Numerical methods for the biharmonic equation and for the two-dimensional Stokes problem, SIAM Rev. 21 (1979) 167-212. [6] P. Grisvard, Singularites des solutions du probleme de Stokes dans un polygone. To appear. [7] J.L. Lions and E. Magenes, Non homogeneous boundary value problems and applications I, II (Springer, 1972). [8] A. Quarteroni, Mixed approximation of evolution problems, Comp. Meths. Appl. Mech. Eng. 24 (1980) to appear. [9] R. Scholz, Approximation von Sattelpunkten mit finiten Elementen, Bonner Math. Schriften 89 (1976) 53-66. [lo] R. Scholz, A mixed method for 4th order problems using linear finite elements, R.A.I.R.O. Analyse Numerique 12 (1978) 85-90. [ll] M.O. Bristeau, R. Glowinski, J. Periaux, P. Perrier and 0. Pironneau, Application of optimal control and finite element methods to the calculation of transonic flows and incompressible viscous flows, Rapport de recherche n. 294, Laboria, Apr. 1978 (also to appear in: B. Hunt (ed.), Numerical methods in applied fluid dynamics (Academic Press, London). [12] R. Glowinski and 0. Pironneau, On numerical methods for the Stokes problem, in: R. Glowinski, E.Y. Rodin and O.C. Zienkiewicz (eds.), Energy methods in finite element analysis (Wiley, Chichester, 1979) ch. 13, 243-264.
MECHANICS COMPANY
AND ENGINEERING
24 (1980) 13-34
ON MIXED METHODS FOR FOURTH-ORDER PROBLEMS A. QUARTERONI Zstitutodi Analisi Numerica de1 C.N.R., Pavia, Italy
Revised
A mixed finite element conditions is described. The problems. The computation methods for block matrices: couple of problems for the Some numerical results for
Received 24 September 1979 manuscript received 15 October
1979
method for the problem v + aZA2v =x with different types of boundary method converges, and it is well suited for the analysis of various evolution of the discrete solution is made by applying a sequence of iterative correspondingly to each iteration either a couple of Poisson problems or a identity operator are solved, according to the value of the parameter (T. two model examples are presented.
1. Introduction
The following problem is considered: 2,
E
H*(LQ, (I + a*A*)u = x
in 0,
(1.1)
where 0 is an open bounded domain of fR*, I and A* are respectively the identity and the biharmonic operators, (T is a positive number, and x is a given function of H-‘(O) (for the Hilbert and Sobolev spaces used in the paper we refer to [7]). The boundary conditions for u are given by either u
=*=o an
over r
(1.2)
over r,
(1.3)
or u=Au=O where r is the boundary operator. Setting
of 0, n is the outward normal direction to r, and A is the Laplace
(1.4) and integrating
by parts, we have (formally)
(1.5)
14
A. Quarteroni, On mixed methods for fourth-order problems
Let 4 be the function of HA(n) satisfying: -vAq
=x. Then with
(1.6)
w=uAu+q
it follows that u and w satisfy
(24w}E Hh(fq x
w,
=0
v+uAw
in g’(n),
w-uAv=q
(1.7)
in g”(0),
where two cases are considered: (i) if the boundary &/an = 0 comes (ii) if the boundary follows that w =
conditions (1.2) are given, then W = If’&?). (the second condition out to be a natural boundary condition); conditions (1.3) are prescribed, then W = Hb(R) since from (1.6) it 0 over r.
The aim of the paper is to analyze finite element approximations (1.7), namely:
{u,4
E
~b(~) x
of problems which generalize
w
v+uAw=p
in g’(0),
w-aAv=q
in B’(0),
(1 .S)
where p and q are two given functions of L2(fl). Remark that the problem (1.1) can be written in infinitely many ways in the form (1.8) since there are infinitely many possibilities of writing x as p - aAq with p E L2(L?) and q E H’(R). On the other hand, systems like (1.8) arise from some time-dependent boundary value problems such as the vibrations of an elastic plate or the Schrodinger equation (cf. [3]). In the quoted examples boundary conditions like (1.2) and (1.3) apply respectively, while (T behaves like the discretization time step At. Using the Green formula (1.5) the weak variational formulation of (1.8) becomes:
{u,WIE Hll(fq v4 EHX(fq
x
w (1.9)
VIJEW Hereafter (e, *) denotes the inner product of L2(R). A displacement finite element approximation of (1.9) is analyzed in the paper; clearly such an approximation also defines a mixed-type finite element scheme for the problem (1.1). Convergence results, and asymptotic error estimates in the energy norm and in the L2-norm are presented. In section 2 the case W = Hi is considered; a lot of simplifications apply since in this case the problem (1.9) is coercive. The harder situation related to the case W = H’(fl) is analyzed in section 3, while
A. Qua~er5ni, On mixed methods forf5u~~-order problems
IS
section 4 is devoted to discussion about the computational aspects of the problem. If W = H~(fIln),an iterative method for block matrices is used to solve the approximate problem. Correspondingly to each step either a couple of Poisson problems or a couple of problems for the identity operator are solved according to the value of the parameter (r. When W = H’(k!), i.e. when boundary conditions are prescribed for zt and its normal derivative, generalizing to our case the techniques developed by Ciarlet and Glowinski and by Glowinski and Pironneau (cf. [2], [4], [5]), then it is possible to reduce the approximate problem to a problem quite similar to the one corresponding to the case W = Hk(ft). Some numerical results for a model problem are also presented. An application of these results to time-dependent problem can be found in [8]. For any real number s the symbol ]].]lsdenotes the norm of I$“&?), I.1denotes the norm of L2(J2). Moreover C denotes a generic constant independent of the mesh size.
2. Analysis of problem (1.9) when W = ~~~(~) In this section the notation foflows:
A is continuous;
moreover,
v = HA(O) x HA(n) is used. Define the bilinear form A over v as
it is coercive since
As a matter of fact, from the Poincare inequality it follows that
Note that since the test functions written as
belong to the same space, the probiem
(1.9) can also be
where
and (,> denotes the duality between H-‘(0) and ~~(~). From the Lax-Milgram lemma the problem (2.3) admits a unique solution for any p, 4 E H-‘(a). Moreover, from (2.2) and (2.3) the following a priori estimate holds:
(2.5)
A. Quarteroni, On mixed methods for fourth-order problems
16
If the data
p, q
belong to L2(R), then with 4 = -V and $ = w from (2.1) it follows that
(2.6)
If-4+ I4 5 IPI + Id
Moreover, using (1.5) in (2.1) one obtains UAW = p - v and uAu = w - q; so the estimate (2.6) leads to 11412 + Ml2 5 CaPI
+
ISI).
(2.7)
2.1. Discretization of (2.3) and error bounds Assume that 0 is a convex polygonal decompositions of J2 into polygons T. Then
6=
h = m;x(diam
UT
TET,,
domain
and let {T,,} be a family of regular
T).
(2.8)
For any integer k zz 0 let P,(T) denote the set of polynomials of degree not greater than k over T, and define
VI, = (4 E
c"(fi))($IT
l&‘-)
VT E W,
(2.9)
& = v, f-UG(R),
(2.10)
v, = & x;,.
(2.11)
The following approximation
vv E I&l(n)
flH’(f2)
VW E H’(0)
properties
are fulfilled by the above spaces (cf. e.g. [l]):
Inf ]]V- 4(]i 5 Ch’-i]]V]]r, 4& I$Jh llw - tilli 9 Ch’-illWJl~,
for i=O, 1 and irr~k+l. The following problem is a finite element displacement
By the same arguments be proved:
used for the continuous
approximation
(2.12) (2.13)
of (2.3):
problem the following stability properties
can
l~hl+lwhl4Pl+l~l
(2.15)
]]Vh]]l + IIWhlll zs $=g (IlPII-1 + 1l~ll-d~
(2.16)
17
A. Quarteroni, On mixed methods for fourth-order problems
As far as the error bound is concerned, the following result applies: THEOREM 2.1. Assume that the solution V, w of problem (2.3) belongs to the space H’+‘(fi). Then the following asymptotic convergence estimates hold:
llu- Uh((l + Ilw- wz(II5
$
h"(llullk+l
]?J- &,I + ]w - wJ,(-= C63’2hk+‘(llv((k+,
+IIWllk+dr
(2.17)
+ IIWllk+,)
(2.18)
Proof. From (2.3) and (2.14)
A({Q-
vh,
W -
WI,},
bf’,
$1)
=
0
%#,
+I) E
(2.19)
vh.
Then (2.17) follows from (2.12) and (2.19) using standard arguments (cf. e.g. [l]). An application of the Aubin-Nitsche trick to problems (2.3) and (2.14) leads to: lu - Uhl+ (w - WI 5 Cb-‘((Iv
- uhlll+ llw - %ll1),
and so the estimate (2.18) follows thanks to (2.17).
cl
3. Analysis of (1.9) for W = H’(R) In this section problem (1.9) with W = H’(LI) is considered. This time the coupling (2.1) for the equations of (1.1) cannot be applied. Let us define therefore a new bilinear and continuous form over the space w = H:(R) x H’(O):
B(b, WI,{A 44)= (VT4) + 6%@I) - 4w
4) + aa(v $ci).
(3.1)
The problem (1.1) can be written as follows: {II, w} E W
WV w], I+, $]) = ](P, q]> (4, +]lo
W,
$1 E w,
(3.2)
where HP, 4], {+Ytw”
= (P,d4+ (4,#cl>.
(3.3)
Remark 3.1. Let (Ybe a positive integer and consider
The operator
I + aA* is an isomorphism
between Hi(O) nH3(0)
and H-‘(Q);
1144135 WI-1
for a positive constant C depending
moreover: (3.5)
on Q (cf. [7]). Also, following the arguments
of [6], there
18
A. Quarteroni, On mixed methods for fourth-order problems
exists a real number v E IO, l] depending IId43+”
5
ca
on the shape of 0 such that
IfI,
(3.6)
assuming f E L*(0), The estimates (3.5) and (3.6) (for some Y 2 0) still hold if 4 satisfies other boundary conditions, for instance if +=a~$=0
over f.
or A4
=@!b=o
•1
over f.
an
THEOREM 3.1. For any u > 0 and for any p, q EL’(LJ) problem solution. Proof. (A) Uniqueness. Setting 4 = z, and t,,Q = W, then from (3.2)
I4 +
I4 5 IPI+ 141.
(3.2) has a unique
(3.7)
Let us emphasize that the bound (3.7) does not depend on u. If p = q = 0, then from (3.7) one obtains u = w = 0. (I?) Existence. Consider the problem (I+u*A*)c#J =p
in 0,
+=0
over r.
Setting u1 = 4 and w1 = aA+, then from the arguments such that {v,, wl} is the unique solution of ~1E H&Ll) n H3+“(0)),
WIE H’+“(0),
v,+aAwl=p
in L*(o)),
w1-uAv,
=0
(3.8)
of remark 3.1 there exists v E IO, l]
(3.9)
in L’(0).
Consider now (I + (~*A*)t,b= q
in 0,
A$=?=
0
Setting wz = $, uz = -ad+, then by the same arguments such that {v*, w2} is the unique solution of 02 E H:(fl) nH’+*((n),
w2 E H3+” (ill),
u2+aAw;l=O
in L*(a)),
wz-oAvz=q
in L*(a).
over f.
(3.10)
it follows that there exists p E [0, l]
(3.11)
A. Quarteroni, On mixed methods for fourth-order problems Now
19
let v = v1 + v2 and w = w1 + w2; then problems (3.9) and (3.11) lead to v
E
H;(n)
v +aAw
f-W+*(n),
w E H’+*(n), (3.12)
=p,
w-aAv=q for a suitable ,$ E [0, 11. Then {v, w} E w and solves (3.2). q THEOREM 3.2. Let .$ E [0, l] be the parameter appearing in (3.12) and let {v, w} be the solution of (3.2). There exists a positive constant C independent of u such that
llvllr+ lb4 5
3
(IPI+ Id),
Proof. The Green formula for the biharmonic (A2u,w)=(Av,Aw)+
(3.13)
r E [O, 1+ 61.
operator gives (formally)
($gw)-((A%$$
(3.14)
where (-, -) and ((e, -)) indicate dualities over r in some suitable spaces. Let us refer to the constructive proof given in step (B) of the previous theorem. From (3.8) and (3.14)
Id’+ a2b412 = (P,4),
(3.15)
(At$12+ u~IA~~#J(~ = (p, A2$).
(3.16)
The equality (3.15) gives
141~IPL
bPl +PI.
(3.17)
From the second inequality of (3.17) one obtains, for a constant C independent
of (T: (3.18)
Moreover,
from (3.16) and (3.17)
IA’41~$2 IPI. Interpolating
(3.19)
2 and I-1(see [7]), then from (3.18) and (3.19) the norms I]*]]_ (3.20)
llA2c#+~+s5 Ca-‘3+6”2jpI. Still using an interpolation
theorem,
from (3.19) and (3.20) the following estimates
can be
20
A. Quarteroni, On mixed methods for fourth-order problems
obtained:
IlVllls= 1144 5 c~-“‘*lPl~ llw,lls-2= llv,ll, +
~~d#4-2
5
&lis
5
C~1-s’2bi,
s E 1%3 + !$I,
(3.21)
51,
(3.22)
r E [O,1+ 61.
(3.23)
s
Ilwllr5 cfr’21PI>
E
+
[0,3
In a similar way, starting from (3.10) and using (3.14) and the interpolation it follows that IIw*lls
llu2l(s-2
(12)2(lr
11~11s 5 c~-%l7
=
=
+
()d+I(,-2
5
~lb~ls
5
C~‘ps’2~ql,
Ilwzllr 22c~-“*ld
theorems (see [7]),
s E [O,3 + SIT 3 + 51,
(3.25)
r E [O,1+
(3.26)
s
E
(3.24)
[O,
61.
Since v = ul + v2 and w = wl + w2, the result (3.13) applies thanks to (3.23) (3.26) and the triangular inequality. cl 3.1. Discretization of (3.2) and error bounds The,,notations used in subsection 2.1 are preserved hereafter. Moreover, we set S,, = V,, 6’h = S,, x S,, and assume in addition to (2.8) that there exists a positive constant 6 independent of h such that h 5 6(diam T)
A decomposition
VT ETh.
(3.27)
{T,,} satisfying (3.27) is said quasiuniform.
For such a decomposition
the
inverse inequality applies, i.e.
for a suitable constant C (cf. e.g. [l]). M oreover, inf-sup conditions hold:
it can be easily verified that the following
where p is a positive constant independent of h. Consider the following displacement finite element approximation
of (3.2):
21
A. Quarteroni, On mixed methods for fourth-order problems
Setting 4 = uh and J, = wh by (3.31) one gets
If p = q = 0,then &, = wh = 0; therefore, (3.31) admits a unique solution. THEOREM 3.3.Assume a = 1, and let {v, w} and (uh, wh} be respectively the corresponding solutions of (3.2) and (3.31). Assuming 21 and w regular, the following asymptotic error estimate holds: Ilu - t&II, + Iw - w,,l 5 Chk-*>
(3.33)
where the constant C depends on Ijullk+land Ilw&+1. Proof. From (3.2) and (3.31) it follows that (V-Uh,+&+wh,&)=O
v4
E
ih,
(3.34a)
(W-Wh,$)+a@-&,,$)=O
vt,b
E
sh.
(3.34b)
In the proof 4, $ shall denote generic functions of g,, and Sh, respectively. (1) EstimateOf IW
-
W,,(
~W-wh~2=(W-Wh,W-~)+(W--h,1(1--h)
(from (3.34b))
=(W-wh,W-$)--a(t?t?h,t,b-W)=-a(Z)-uh,W-Wh)
=(W-Wh,W-(CI)+a(~-~h,W-(Ir)+a(~-~,Wh-W)+~(~h-~,W-Wh),
(3.35) (i) bound for u(v -
uh,
w
-
$):
Then thanks to (3.30) there exists 6 E Sh such that U(U - u,,, w -
$) 5 Ilw- ~~~1{~~~ - d& + P-‘llill%~
= ((W - I,b~~~{~~TJ - (6111 + P-‘(a(C$5
-
Uhr
II;)}
u, 6) + u(u - uh, 6)))
IlW- $ll,{llu- 4111 + p-‘(llu - +I\1+ lw- Whl)}
(from (3.34b)), (3.36)
(ii) bound for u(v - 4, wh - w): U(U - 4, w,, - w) 5 5
Ilu- $lll(llw- $111 + llwh- @IId c((U- ~~~,{~~~ - #(I,+ h-‘(lwh- w1+ lw-
$1))
(from (3.28)), (3.37)
A. Quarteroni, On mixed methods for fourth-order problems
22
(iii) bound for a(~$, - f#h W -
wh):
(from (3.34a))
&h-$,W-Wh)=-(u-Uh,+Uh)
= -(ti - v,,, v - vh) - (0 - vh, (a - v) 5 -(v - v,,, 4 - 0)s (v - v,,[lv - +[ lv-~h(+-~I+I~-%t
%
Iv
-
41
+
P-‘i&‘(d6
-
2)~ +)+
a(u
-
oh,
5 ~-~I+P~‘(ll~-4lll+Iw-wtl)
(from (3.30))
9))
(from (3.34b)).
Hence (3.38) Substituting
(3.36), (3.37), (3.38) into (3.35) it follows that
Iw - &I2
d~‘<~~w - till’ + h-‘(l&t’ - whl+ Iw - $I))+ Iv - 41<1v - $1+/Iv- dl + Iw- Whl)h (3.39) 5
c{lw - wh((w- +I + IIw- $lll(llu- +(I’+ lw - whl>+lb -
(2) Estimate of 11~1 - uhj)l Using (3.30) once more, one gets I/v - Vhlll
5
/Iv - (j$lj’+
114- vh)),s I/v- 41)’+~-‘tl&‘“(~
-
vh’
‘I’
From (3.34b) I(V -
Z&II,
5 [IV - &‘([I +@-‘(I&‘(@ 5 (1 +p-‘)llV
- VT i)+
- +I(, +p-‘lw
-
whj
(wh - w, 6)) 5
c(((u
-
411’
+
iw
-
wd)’
(3.40)
Going back to (3.39) and using the results (2.12) (2.13) (w - wh125 c{lw - wh(hk+’+ hk(hk + (w +
hk(hk
+
whl)
hk + h-‘/w - w,,l)+ hk+‘(hk+’ + hk + (w -
(note that q5 and rj~are arbitrary functions of &, and Sh respectively). IW -
w,,l~ C’hkp’,
cl
=
C’(llv(lk+‘,
Whl)j
Therefore
(IW((k+l).
(3.41)
Finally, thanks to (3.40) (2.12) and (3.41) the following estimate applies:
[Iv - v,,)I,5
C,hk-',
c2 =
c2(((v((k+‘,
((Wllk+l)*
(3.42)
A. Quarteroni,
On mixed methods
for fourth-orderproblems
23
Clearly (3.41) and (3.42) give (3.33). ~~~ORE~
0
3.4. Let B be an arbitrary positive number; if the assumptions
of theorem 3.3
apply, then Ilk>- Uh((,I C(hk-’ + Q -‘hk + (Z--“h’+’+ f_-“W),
(3.43)
Iw- WJ I C(hk + ,P’
(344)
+ a-~hk+’ + k%zk).
Proof. This is quite similar to the proof of the previous theorem.
q
3.2. L’-error bounds for linear elements When piecewise linear elements are used, i.e. if k = 1, the estimates (3,33), (3.43) and (3.44) do not provide an asymptotic convergence for 11~- r.+,((,and Iw - whl. The aim of this subsection is to provide an asymptotic error estimate in the L2-norms even for this case. Let us define the Ritz operators by ~,:H:,(.n)~~~:a(t:-i;ht’,~)=O
V& E &_
(3.45)
VJ, E Sh,
(3.46)
vu E H;,(n) n w(a)
lu _ ii,Vl + hll1,- &VII,5 chqlull,.
(3.47)
VW E H’(.n)
jw - R,wf + hllw - RhwI(, 5 C!z’llw&
(3.48)
Rh : H’(L+
Sk : a(w - RhW, J/) = 0
(w - R,,w, 1) = 0. The following results apply (see e.g. 191):
for 1 I r zz 2. Generalizing (see [lo]), we can prove
the techniques
developed
by Scholz for the biharmonic
problem
THEOREM 3.5. Let CT= 1, k = 1 and (u, w}, {uh, wh} be respectively the solutions of (3.2) and (3.31). Then for some constant C depending on I/r&, ~~A~~l~,~~~~~, Ilwli2the following result applies:
Iu- vhl+ h”2(ln
hllw - wh( 5 Chlln h12.
(3.49)
FTOOf. Setting el = v - v4 and e2 = w - wh, (3.34) gives (er.ct,)--a(e2,#)=0
V$ E &,
(3SOa)
tc,, 4)+ &I? $) = 0
vJ/ E Sk.
(3.50b)
24
A. ~uu~ero~i, On mined methods for fourth-order problems
(1) Estimate of le,j Let 5 = &,u,
q = &W ;
IWh- r112 = (Wb-
then (3.51)
w,wb--r))+(w-q,wh--v)
(from (3SOb))
(Wk- w, Wk -~)=-(e2,Wr*-r))=~~~i,Wb-r))
=a(V-StWh-rl)f~(Z‘-~k,Wk-77), a(~-Ub.Wh-rl)=~(5-~h,Wh-W) =
-f&,
r - %uh)
(3 S2)
(from (3.46)) (from {3.5Oa))
= -(er, e,)+ fez, 2,- 8)” (6, 0 - 5).
(3.53)
Let us note the result (cf. [lOI)
where C is a positive constant independent state
of U, b, and $J. An application
of (3.54) allows to
(3.55)
la(v - 6, wh - r))l s CIP211n hl~~AullL~~~~~w~ - r)l.
(3.56) (3.57)
(3.58)
Let 6 EW@?il) satisfy 6 -t-A26 = el in 0. Assuming estimate applies for a suitable Y E IO, l] (cf. remark 3.1):
0
regular,
the following
a priori
(3.59)
A. Qwarteroni, On mixed methods for fourth-order problems
25
Moreover, (3.60)
)e112= (0 + A*& el) = (0, el)-- u(A0, e,) = (0, et) - u(A0 - R,,AO, e,) - a(R,,AO, el) - a(RhAO, el) = (e2, R,Ae)
(from (3.50b)) (from (1.5))
= (e*, RhAe - AtI) - a(e2, e)
= (e2, RhAe - de) - a(e2, 8 - k,O) - (e,, k8)
(from (3.50a))
(3.61)
Applying (3.47), (3.48) and (3.54), the following estimates are deduced:
la(Ae - R,,Ae, e,)l = [a(Ae - R,,Ae, u - t)ls
C2h’+Y~l~~~2~~e~(3+v,
I(e2,R,,Ae - AtI)1 5 C3h1+Y~~8~~3+v~e2~, (u(e,, 8 - rZ,e)ls \U(W - 7),8 - &e)l
+ Iu(~ - wh,8 - iLe)l,
5 C4h211wl12(l1412 + IIA~llL~~nJ+ Gh1’211n 4 - llAf4hde21.
Since from Sobolev’s imbedding theorems H3+V(0) C W22m(L?), then IIAellL-cn,5 C611e((3+v, therefore, from (3.59), (3.60) (3.61) and the previous inequalities one obtains (e112 I C{h21e112+ h1+“((1412+ llwl12>leIl+ (II’+” + h”*lln
Wle2ll4l.
For h sufficiently small the term 1 - Ch is positive, then IelII C{h1'"(llu(12 +
Ilwl12)+ (h’+”
+ h1’2(ln hl)(e&.
(3.62)
Using (3.58) from (3.62) it follows that (ell 5 C{h”“(l(u(12+ jl~11~)+(h3’*+“(lnhj +
(h3’2\1nhi+
+ hlln h12]lAz&-~n~
h2’y)~fl}
(3.63)
From (3.63) it easily follows leI( I Chlln
hi2
for a positive constant depending by (3.64), it follows that
(3.64) Going back to the estimate (3.58) on Ilul12,IlAo IIL=(R)and (Jw(12.
[es1I Ch1’2(ln hj. Now the result (3.49) applies thanks to (3.64) and (3.65).
(3.65) 0
26
A. Quarteroni, On mixed methods for fourth-order problems
4. Numerical algorithms to solve (2.3) and (3.2) The aim of this section is to propose a method for solving problem (2.3) and a general method for reducing (3.2) to a problem whose numerical solution is quite similar to the one of (2.3). To do this, some notations are needed. Assume that the known functions p and 4 are accurately approximated by two piecewise polynomials ph and @, of Sh. For fixed k, let fi and N (fi< N) be the dimensions of Sh and Sh, respectively;
moreover
of N independent
let {&} be a set of N independent
functions of &,, and let {$j} be a set
functions of S,,, Set:
A = [A,],
Aij = (4iy4j);
c = [Cj],
Cij
D=
= ($iy $j);
toi,],
Dij
B = lBij 19 Bij
=
=
a($i,
4j),
a (41)$j).
4.1. Solution of problem (2.3) The discrete problem (2.3) is equivalent
to the linear system
;;u=f
(4.1)
where k is a symmetric, positive definite matrix defined by i;=
[
u;
-iD .
(4.2)
1
moreover, u =(%BY,
f =(PdI)
and 2
ai+i
=
vh,
2
Pi6
=
why
pi
=
(ph,
+i),
4i
=
(qh,
6i).
The system (4.1) can be solved by the following iterative scheme: For an arbitrarily fixed /!I” define the sequence cy” and p” by solving: aDa”“’ +A@" Aan+'
= q,
- uD@"+' = p.
Solving for the unknowns,
(4.3) n 20.
(4.4)
eqs. (4.3) and (4.4) become (4.5)
A. Quarteroni, On mixed methods for fourth-order problems
27
The scheme (4.5) defines a block Gauss-Seidel iterative procedure for approximating (Y and 6. A sufficient condition for convergence is that the maximum norm of the iteration matrix be less than 1; this is equivalent to require that the maximum norm of the matrix a-‘D-IA be less than 1. Note that by definition the elements of D-‘A have the upper bound Ch2 for some positive constant C; hence, setting C1 = C’, a practical condition for convergence is h2 <
C,a.
(4.6)
If the parameters h, u do not fulfill the condition (4.6), then instead of (4.3) (4.4) the following scheme can be used: For any given vector /3” compute the sequence a”, p” by solving
AP“+I
+ g&p+
zz
q,
Aa*+’ -uDjY’=p,
(4.7) It 20.
(4.8)
The explicit form of such a scheme is (4.9) By the same arguments used to achieve (4.6) it can be seen that a sufficient condition in order to get convergence of the iterative scheme (4.9) is a
c C2h2,
(4.10)
for a suitable constant C2. Let us remark that for evolution problems u behaves like the discretization time-step At; then for h sufficiently small one of the conditions (4.6), (4.10) holds except for At - h2. However, it should be emphasized that the best situation is that in (4.6) at least for implicit schemes. When At - h2 or (4.10) applies, then an explicit scheme seems more suitable for the computations. Remark 4.1. It can be verified that at any step n to solve (4.5) means to solve two Poisson problems: -uAv ;+l = W ;: - q,,,
(4.11)
-UAW;+’ = -v”h+’+Ph.
(4.12)
On the other side, to solve (4.9) means to solve the following two problems for the identity operator: V ;+’ W;+’
= =
-UAW “h+ ph,
(4.13)
uAv;+’ + q,,,
(4.14)
in both cases a starting function wi = Xip’l4l is needed.
28
A. ~~a?ie~on~, On mixed meiho~ for fourth-order problems
4.2. Solution of problem (3.31) Setting
u =bm,
f =cPd
and c y&i =
Uhr
2
&$j
=
wit,
pi
=
(ph,
(bill
qj
=
(qh,
$j),
i
the
problem (3.31) is equivalent
to the linear system (4.15)
Hu=f,
where the matrix II is given by (4.16)
H =
Since the matrix B is not square (%< N), systems like (4.15) cannot be solved separately as two sets of equations, and so the method proposed in the previous subsection for the system (4.1) does not apply. On the other hand H is not well featured (it is undefinite). In this section, following [5] closely, we introduce a method for solving (3.31) which avoids the direct solution of the system (4.15). Let (u, w} be the solution of (3.2) and let A be the (unknown) trace of AU over f. Clearly from (1.8) it follows that the solution of (3.2) satisfies v+uAw
=p
/ w-uAv=q
in 0, in fz,
1V=O iw=ai\+q
(4.17) over IY
In order to write the last equality over F, the trace of q over 1” must exist; that happens if, for instance, q E H’(f2) or even if q belongs to the space N(A, fz) that we are going to define. Let H(A, a) = {# E L*(O)]A& E L2(Q)], and for any p E H-“2(r) ‘,y+aAo
=O
let A (cc) = (x, o} E H*(J2)
H(A, L?) be the solution of
in 0,
w-crdx=O
in 0,
x=0
over r,
.o =a#.&
X
over IY
(4.18)
A. Quarteroni. On mixed merhods for fourth-order problems
Moreover,
let {v,,, w,,} be the solution of vu -t UAW,!= p
in 0,
wo - udvo = q
in L?,
V()= 0
over I:
WI& =q
over f.
(4.19)
Setting: (6. @}= A(A),
(4.20)
it follows from (4.17). . . . , (4.20) that v = VI,+ 6,
w =
w,,+ G.
(4.21)
Consider now the operator
/.L+u”~
,:H-“*(f)+Hii2(f);
over jT,
(4.22)
where by definition (x, o} = A @); E is an isomorphism since A is an isomorphism, and the application which maps 4 into its first-order trace ~~~~~ over r is an isomo~hism between HZ(a) and H”‘(f) (cf. e.g. f71). Finally set e : H -l/‘(f) x H’/2(f)-+
IR;
e(A, p) = (EA,p).
(4.23)
Here the duality is between H”‘(f) and H “‘(I-‘). The form e is continuous moreover by (4.20) and the de~nition of e the following relation holds:
e(A, p) = u2 -
Ir$$dy=u -
dy
(using the Green formula (3.14))
(u2A2,y + x, 6) + h, i2)-t-a’(Ax, AC) = (,y, a) + (w, @)
Hence, for any 4 = (4,, &) and $ = (#,, (cIz)belonging to L2(L!) x L’(R), [#%J/l= @I> $1) +
(42.
$2)
and bilinear;
(from (4.18)) (4.24) setting (4.25)
we get from (4.24) e(k/J)=[Nq,NP)l; hence, e is symmetric
(4.26)
and coercive over H ^ ‘l”(r). Since a~/a~ = 0 over f (cf. (1 .l)) from
A. Quarteroni, On mixed methods for fourth-order problems
30
(4.20), . . . , (4.22), it follows that the trace A of Au over T is the unique solution of: (4.27) Setting once more h, w} = A(p), from (4.18) (4.19) (4.27) and the Green formula (1.5) for the Laplace operator it follows that
I dvo
-CT2
I‘ an
/A
d y = -o@ut,, VW)- u(&,,
w) v/AE w-1’2(r),
= -a(Vuo, VW)- (wo, 0) -+(43 @)
and so (4.27) becomes
vp
E
Er-“yr).
(4.28)
Therefore, in order to obtain the datum h of (4.17) we must solve the variational problem (4.28) (which in turns requires, for the determination of the right-hand side, the solution of problem (4.19)). Passing to the numerical approximation of (u, w), let {oh, wh) be the solution of (3.2). Let A& be a finite-dimensional space contained in S,, such that Sk =&+MJ,.
(4.29)
We shall think of Mh as the complementary of $, into S,. Using Lagrange interpolation, any function of Mh can be characterized boundary nodes. Let Ah be the component of wh into Mb. Setting Wh
by its values at the
=/i,+(w,-rc,)
and following [2], the triplet {v,,, wh, h;r} can be characterized
(4.30) as the unique solution of
where for any 5 belonging to Mh we denote by Sh(e) the following subset of Sh:
Since qh E Sh let us set
A. Quarteroni, On mixed methods for fourth-order problems
The following problem is an approximation
31
of (4.19):
(4.35) Moreover,
for any j.&hE A& let (xh, ah}
=
be the solution of
Ah&)
and set
{i?k,'c3h) = Ah@k)v
(4.37)
where hh EM,, is defined by:
ik = CrAh+
(4.38)
ijk.
From (4.31), (4.35), (4.36) and (4.37) it follows (in analogy with (4.21)) that vh
=
f,+jk +
ek,
Wk =
WOk+
ek(hh,
pk)=
w-k-
(4.39)
Now, following (4.26), we set ek
:Mk
x Mh
+R
;
The form eh is bilinear, continuous follows that hk is the solution of eh(hk,
pk)=
-u(bOk,
[Ak@k),
Akbk)].
(4.40)
and coercive. From (4.32), (4.35), (4.36), (4.38) and (4.39) it
VWk)-(WOk7
ok)+
(% %)
(in analogy with (4.28)). Therefore, in order to compute the solution {t& wk} of (3.21, we must (i) compute the solution {&,h,WOk}of (4.35), (ii) solve the yariational problem (4.41), (iii) COmpUte hh by (4.38), (iv) solve the problem (4.31). Now let us make some remarks about steps (i), . . . , (iv). Step (i), (it,). We can easily check that the problem
(4.41)
32
A. Quarteroni, On mixed methods for fourth-order problems
(with the right-hand
side eventually zero) is equivalent
to the linear system
iiz
(4.43)
=b
whose matrix is defined by (4.2). Then (4.43) can be solved by the method proposed in subsection 4.1. The nonvanishing boundary condition for wh only changes the right-hand side of the system. Therefore, both problems (4.31) and (4.35) can be reduced to a linear system like (4.43). Step (ii). Let us denote by {pi}:, a set of basis functions for Mh (N is the number of boundary nodes). Setting
Eij = et,bi, hi,
E
mi}
=
=
[Eij],
Ah
bi
P,
1,
Ci = -(T(hOhy
), A
=
[A’],
it easily follows that (4.41) is equivalent
c
=
vwi)-
[c-i],
(WOh,
@i)+
(qh,
wi),
1 Ii,jlN,
to the linear system
EA = C.
(4.44)
The matrix E is symmetric and positive definite; so (4.44) can be solved by direct methods (for instance by the Cholesky method) or by iterative methods. (Using for instance the conjugate gradient method, an explicit construction of E can be avoided provided two Dirichlet problems for -A are solved at each iteration (see [5, section 5.41 for details). Moreover, E does not depend on the terms p and q but only on the geometry of 0. This is an important feature since for problems that require a recursive solution of (4.15) (this is the case of the evolution problems), the solution of (4.44) must be done only once. Since the functions U()h,W()h,qh and wi are piecewise polynomials, the right-hand side of (4.44) can be correctly computed without using quadrature formulas. Moreover, in order to compute xi and wi, we must solve the problem (4.36) with ,& = pi. That in turn implies that a systems like (4.43) must be solved (cf. the discussion about steps (i), (iv)), where b consists of the contributions given by the boundary values of wi. Such a system can be solved by the method proposed in subsection 4.1. In order to compute C, the solution of N systems like (4.43) is required. Step (iii). Once Ah is known, the computation of & only requires the computation of &,, 0 and that is trivial since, if Mh is the orthogonal space of sh with respect to S,, for the H’ inner product, we get (4.45) the summation being restricted to the boundary nodes (7,). Finally we can conclude that the computation of the solution {zlh,wh}
of (3.2) requires the
33
A. Quarteroni, On mixed methods for fourth -order problems
w(x,,x,f=oAv+i v(~&f=$D21x,)P”o?%~ t ? ! \
i i i ! ! \
i i i i i i
i i
_____ “___(I(J_1., \I rel.error
_--._
--___ i
\
i 0 ii+4
i i i
--------_.-.-,-
1
-----
! ! ! \
6
8
10
(J= 1. t re, ,?rror U=.l I . (i=
6_
1.
.,
\
j rel. error
* l/h
_.
u- 1. \ ~_ .11 rel.
error
for V for W
i i i i
i i i
tar V for W
i
i i
--r,, Fig. 1. Results for problem (2.3) with data in (4.46).
4
6
8
f0
c
l/h
Fig. 2. Results for problem (3.2) with data in (4.47).
solution of N + 2 systems like (4.43) (by the method proposed in the previous subsection) plus cl the solution of the system (4.44) whose matrix is well featured. Rebury 4.2. We refer to [ll) and [12] for the analysis of methods similar to the ones of section 4.2 related to steady and unsteady Stokes problems using the velocity-pressure formulation via also a mixed finite element method. q In figs. 1 and 2 we show some results related to problems (2.3) and (3.2) respectively, solved by the methods proposed in section 4. We assume 0 = (0,l) x (0, l), we consider decompositions of fi by squares and piecewise bilinear polynomials. Fig. 1 is related to a problem like (2.3) solved by the algorithm of section 4.1; the model problem corresponds to the following data: p(x,, x*) = #+M(xz)-
2~7M&)X(&),
q(xl, x*) = X(~I)X(=z) - 2G(#(G) f tfi(X2))T
(4.46)
34
A. Quarteroni, On mixed methods for fourth-order problems
i = 1,2. where 4 (Xi) = Xi(1 - Xi), x(x~) = sin TXi, Fig. 2 is related to a problem like (3.2) solved by the method proposed corresponding to the data
&I,
x2)= 1.
in section 4.2 and
(4.47)
The automatic procedures terminated when the distance between two successive iterates was less than 1O-4in the maximum norm. The relative errors are evaluated in the L2-discrete norm. References [l] P.G. Ciarlet, The finite element method for elliptic problems (North-Holland, 1978). [2] P.G. Ciarlet and R. Glowinski, Dual iterative techniques for solving a finite element approximation of the biharmonic equation, Comp. Meths. Appl. Mech. Eng. 5 (1975) 277-295. [3] J.E. Dendy, Jr., An alternating direction method for Schrodinger’s equation, SIAM J. Numer. Anal. 14 (1977) 1028-1032. [4] R. Glowinski and 0. Pironneau, Sur la resolution, via une approximation par elements finis mixtes du problemes de Dirichlet pour l’operateur biharmonique par une methode “quasi-directe” et diverses methodes iteratives, Rapport 197 (IRIA, Le Chesnay, 1976). [S] R. Glowinski and 0. Pironneau, Numerical methods for the biharmonic equation and for the two-dimensional Stokes problem, SIAM Rev. 21 (1979) 167-212. [6] P. Grisvard, Singularites des solutions du probleme de Stokes dans un polygone. To appear. [7] J.L. Lions and E. Magenes, Non homogeneous boundary value problems and applications I, II (Springer, 1972). [8] A. Quarteroni, Mixed approximation of evolution problems, Comp. Meths. Appl. Mech. Eng. 24 (1980) to appear. [9] R. Scholz, Approximation von Sattelpunkten mit finiten Elementen, Bonner Math. Schriften 89 (1976) 53-66. [lo] R. Scholz, A mixed method for 4th order problems using linear finite elements, R.A.I.R.O. Analyse Numerique 12 (1978) 85-90. [ll] M.O. Bristeau, R. Glowinski, J. Periaux, P. Perrier and 0. Pironneau, Application of optimal control and finite element methods to the calculation of transonic flows and incompressible viscous flows, Rapport de recherche n. 294, Laboria, Apr. 1978 (also to appear in: B. Hunt (ed.), Numerical methods in applied fluid dynamics (Academic Press, London). [12] R. Glowinski and 0. Pironneau, On numerical methods for the Stokes problem, in: R. Glowinski, E.Y. Rodin and O.C. Zienkiewicz (eds.), Energy methods in finite element analysis (Wiley, Chichester, 1979) ch. 13, 243-264.