- Email: [email protected]
Interpolation methods for the construction of the shape function space of nonconforming finite elements
Computer methods in applied mechanics and enginaaring
Computer Methods in Applied Mechanics and Engineering 122 ( 1995) 93-103
Interpolation methods for the construction of the shape function space of nonconforming finite elements J.B. Gao”***‘, T.M. Shih b a Department of Mathematics, Huazhong University of Science and Technology, Wuhan 430074, People’s Republic of China h Department of Applied Mathematics, The Hong Kong Polytechnic University, Kowloon, Hong Kong
Received 1 December 1993; revised 1 April 1994
Abstract Based on the multivariate polynomial interpolation theory, a method is derived to construct shape function spaces which can pass relevant convergent tests. As applications of this method, the shape function spaces of de Veubeke’s, Zienkiewicz’s and Specht’s elements are obtained. It is also shown that this method can be used to construct new finite elements.
1. Introduction Galerkin’s finite element method is a simple and efficient way for numerically approximating the solution of elliptic problems, especially for solving problems of continuum mechanics possessing infinite degree-offreedom. In order to solve a fourth-order boundary value problem such as the plate bending problem with the clamped boundary conditions, we have to choose higher degree polynomial finite elements in virtue of C’ compatibility conditions. The difficulty caused by higher degree elements can be surmounted by the socalled microelement method with splitting procedure, such as the Hsiech-Clough-Tocher splitting triangle and the Powell-Sabin triangle, which result in more complicated computation. It is well known that the nonconforming elements were introduced, such as Zienkiewicz’s triangle (Co class) (cf. [ l] ), de Veubeke’s triangle (cf. [ 2]>, the TRUNC element and the quasi-conforming plate element (cf. [ 3-51). because the requirement of conforming conditions is very harsh. Unfortunately, the nonconforming element methods cannot always ensure the convergence of approximated solutions. Irons introduced the “Patch Test” in continuum mechanics to conclude the convergence, which was one of the greatest contributions to the ‘theory’ and ‘practice’ of the finite element method (cf. [ 6,7] ) . Since then some other techniques have been developed to check sufficient requirements for convergence such as the “Interpolation Test”, the “Generalized Patch Test” (cf. [ 81) (sufficient and necessary test), and the “F-E-M-Test” (cf. [ 91) . In general, interpolation test and patch test can be used to develop useful, accurate and convergent nonconforming finite elements (cf. 16,101). Since test conditions may be viewed as interpolation conditions, the construction of nonconforming finite element is actually an interpolation problem. The problem is how to construct the associated shape function space. For example, Specht (cf. [ lo] ) wrote ‘The required three higher terms are assumed to be linear combinations of the following cubic and quartic terms: (f&2, $$$s, &,$i, ,t$&&, * Corresponding author. ’ Currently visiting The Hong Kong Polytechnic University. 0045-7825/95/$09.50 @ 1995 Ekevier Science S.A. All rights reserved SSDIOO45-7825(94)00728-7
94
J.B. Gao, TM. ShihlCotnput.
Methods Appl. Mech. Engrg. 122 (1995) 93-103
and St&&. This assumption is successful.“. However, it is not at all obvious how Specht gave the shape function space for his element. Another example is the famous de Veubeke’s triangle which passes the patch and interpolation tests. The nodal parameters of de Veubeke’s element are the function values at vertices of the triangle and the means of function values and first derivatives in outward normal direction on each edge of the triangle. How do we choose a subspace of 7rs as shape function space? We should find a method to determine some polynomial subspace (shape function space) from which the interpolation is possible and unique in terms of the interpolation conditions (such as the patch or interpolation tests or F-E-M-Test conditions). The recent work on multivariate polynomial interpolation is a useful for approaching this problem (cf. [ 11,121) . The following problem is thus interesting and important: how do you construct a subspace of polynomials so that the shape function spaces for given nodal interpolation parameters and the corresponding finite elements can pass the relevant tests? Here we aim to derive a method which can answer the above question. As examples, de Veubeke’s element, Zienkiewicz’s element and Specht’s element, in Sections 3, 4 and 5, respectively, are constructed by the new interpolation method. A new triangular element, refered to as the modified Specht’s element, is obtained. (,&$3,
2. The fundament
of interpolation
methods
The problem of multivariate polynomial interpolation is a complicated and classical problem. The generalization of univariate polynomial interpolation to the multivariate context is very difficult. Most efforts went into identifying sets of points @ C RS on which one can interpolate uniquely from rk (polynomial space of total degree < k) for some k. Our method comes from another direction. We wish to determine a certain corresponding polynomial space Q directly for any given interpolation conditions. Although such a space Q is not unique (cf. [ 1 l] ), it is enough for us to choose just one of them, with some criterions such as a lower degree about Q. An interpolation problem can be presented in the following abstract way. Let H be a finite dimensional linear subspace of a linear space X and n be a finite dimensional linear subspace of the dual of X (over C or R). The interpolation problem on H can be described as follows: Given f E X, find an h E H such that Ah = Af,
V’A E A.
(2.1)
If there is, for each f E X, exactly one solution h E H such that (2.1) is valid, the interpolation problem (H, A) is said to be correct, where A and H are the (space of) interpolation conditions and the space of interpolation for the problem (H, A), respectively. Recall (without proof) the following well-known characterization of correctness. LEMMA 2.1. Let H and A be$nite dimensional linear subspaces of a linear space X and its dual, respectively. The following three statements are equivalent: (I) The interpolation problem (H, A) is correct. (2) The linear map TH : h + (Aj h)? is one-one and onto, where ( Aj); is any basis of A. (3) A % H’ in the sense that the map F* : h --f AI* is one-one and onto. In the rest of this paper let X = Ao, the space of all functions analytic at the origin, the collection of all polynomials on RS. Our main interest is the polynomial interpolation interpolation conditions in linear functional form as follows:
[BIP(D,) : f + (p(D)f)(N
and denote by r problem with all
Vf E Ao,
where[B]f=f(B)foranyf,D=(Di ,..., D,),cr=(a(l) ,..., a(s)),andD*=D;Y(‘) ,..., D?“). Associate each point 0 in RS with a polynomial space Pe and define the space of interpolation conditions A := I{]%
(2.2)
: P E PO},
Vq E 7~ and f E Ao, define a linear functional
as
in the following
form:
95
J.B. Gao, TM. Shih/Comput. Method.? Appl. Mech. Engrg. 122 (199.5) 93-103 q*(f)
= (s(D>f)(O)
Then the following
FfW~!
= p%(O).
(2.3)
equality holds,
bbq E rt
[dp(D)q = q*(eep),
(2.4)
where ee : x E RS -+ ee.x. The map q --+ q* is linear and linear embedding of rr in the dual of Ao. A polynomial space Q is to be determined for which H = C es ~0 c Ao, then the interpolation problems (Q, A) From the problem (H, Q*), Q* can be used to represent Lemma 2.1. This means that it is enough to find H’. THEOREM 2.2. [I 1,121 H’ Z HI with HI =
( fl
one-one
in the sense of (2.3),
hence it provides
a
(Q, A) is correct. Let Q* = {q* : q E Q} and and (H, Q*) are equivalent in the sense of (2.4). H’, the dual of H, and conversely according to
: Vf E H}, where f 1 is the homogeneous polynomial
of
largest degree j for which f(x)
=fl(x)+o(IIxllj),
asx+O.
EXAMPLE 2.3. Let f (x, y) = x2 + xy + y3 + x*y + xy2, then f 1(x, y) = x2 + xy. We have the following
procedure
f =clfl
+c2f2+...+c,f,
=
L(cl,C2,...,C,)~xn
c
of HI. For any f E H, let { fj};
for the construction
be a basis of H, then
a
where La(ct,c2,. . . , c,) are linear forms of ct, c2,. . . .c,. Denote Ia\ = a( 1) + a(2) + . . . + a(s). Then we get the following THEOREM 2.4. [II,121
result.
Assume that above conditions are satisfied, then there exists an integer mg for which
~La(cl,c2,...,cn)x”
: L,(c1,c2,...,cn)=0,
EXQMPLE 2.5. For interpolation conditions n = span{f(xu), find a polynomial subspace Q such that (Q, A) is correct. Define a finite dimensional
space associated
lal
f’(xl),
f(x2))
(x0
< x1 6 x2
and x0 $x2),
with A as
H = {cl eXO*+c2x exlx +cg exzx : ci E R}. For any f E H, let f = ci exox +c2x
e”” +c3 exzx = cl + c3 + x( xocl + c2 + X2C3) +
x2
+ (;x;c,
XIC2 +
;x;c, )+2 (
1
3
,x&l
+
;x;c,
+
If ct + cs $0, then fl = cl + c3 = c (constant); hence 1 E HI. If cl + cg = 0 but xnct + c:! + ~2~3 SO, then fl = ((x2 - xo)c3 + c2)x; hence x E HI. If CI + ~3 = 0, and XOCI + ~2 + ~2~3 = 0, then f = ~3 (~2 - ~0) (T - ~1) x2 + (; (x;
Ifxt
+ x2x0 + xi, fl
$T,then
If xi = T, Dehne
-
fx?)
x3 +
=cg(x2-x0)(
. . .. T-xt)x2;hencex2
span{l,x,x2},
x1 #-
{ span{l,x,x3},
x1 =
Q=
x0 + x2
x03
x2 ’
2
’
>
c3( ~2 -
~0)
E H~.Hencema=2andH~=span{l,x,x2}.
and hence x3 E Hi. Consequently,
then fl =c3 -x3
+.
..
ixg3
mo = 3 and HI = span{ 1, x, x3}.
96
J.B. Gao. T.M. Shih/Comput. Methods Appl. Mech. Engrg. I22 (1995) 93-103
then (Q, LI) is correct.
3. Interpolation
construction
of de Veubeke’s element
K with vertices at 9 = (xi, yi), i = 1,2,3,
Given a triangle 51 =x2--3,
(2=x3--1,
53=x1
--X2,
vl=Y2-Y3,
r)2=Y3-Y1~
773=Yl
-Y2,
~:=s:+&
z;=s;+r1;,
denote by A the area of K and put
(3.1)
I;=&+&
(3.2)
Let ri and IZ~ (i = 1, 2, 3) be the unit tangent and outward normal on the edge Pi+tPi+Z opposite vertex pi (i = 1, 2, 3), respectively, and Wi = Xi+tyi+z - xi+zyi+t (i, i + 1, i + 2 E 23). NOW let Ai denote the area coordinates related to the vertices Pi, i.e., x = XI At + X2A2 + Y =
Yl Al +
to the
X3A3,
(3.3)
Y2A2 + Y3A3,
1 = A1 + A2 + A3, such that the triangle K is transformed into the standard simplex K* = { (AI, AZ, As) 1 At +A2+A3 = 1, Ai > 0). By (3.3) any function u(x,y) defined on K can be associated with a function w(At, A2,As) defined on K* such that u(x, Y> =
44,
(3.4)
A29 A3).
Now it is not difficult to prove that Zig=-: 1 +--I
at4
tig+ri-l (
aw
dAi_1
s+ri+l-
I
lfl
aw
dAi_ 1>
9 (3.5)
aAi+l'
where i,il,i+ 1 E Zs. De Veubeke’s element (cf. [ 21) is triangular, which satisfies the F-E-M-Test (cf. [ 91) (also passes the patch test (cf. [ 61) ) that ensures convergence of the element. The nine nodal parameters are the function values at the vertices Pi and at the middle points Mi = i (Pi-1 + Pi+1), and the mean values of first order derivatives in outward normal direction on each side of K. There are nine interpolation conditions and in general the one has to give the other linearly-independent interpolation condition so that the interpolation problem (7r3, LI) is correct; II is the interpolation condition space of ten dimensions mentioned above. Some authors had given the methods adding the tenth interpolation condition. For example, Shi showed in [ 131 that the condition given in [ 141 is linearly dependent with de Veubeke’s parameters (or interpolations) on r. Now a shape function space B of de Veubeke’s element can easily be obtained with the new interpolation approach without the tenth condition. With the transformation (3.3) the nine interpolation conditions about u (defined on K) can be transformed as, with respect to w on K*, U(Pl)
= w(l,O,O),
u(P2)
= w(O,l,O),
u(P3)
= w(O,O, l),
(3.6)
J.B. Gao, TM. Shih/Compur. Methods Appl. Mech. Engrg. 122 (1995) 93-103
1
pJ
au
-ds=-& II s pz ant
1
tt [( +
$0,
r3
(
pl I ps
G
au
-ds=-&-
2
an2
t2 I( + rl
+ r3
2
1,O) +4g(o,
;) + $(O,O,
1)
2
1)
1) + -$(l,O,O) 2
g?o,o, 3
1) +4$+,0,
-g(O,O,
1) +4$$0,$
3
(3.7) >
>
(3.8)
;) + $(l,O,O) 3
>
+ -$l,O,O)
1
1
p2at4
1 13
>
;, 5 + -$o,o,
I> +4g(;.o,
1
1)
1
+4$0&
g(o,l,o) 2
g(o,o,
(
;, ;) + $(O,O, 1
-g(O,
+ 12
1
I,O) +4$0,
I
97
Jp,
-d,=-&an3
$(l,O,O) 3
3 Ht3
+4$+,
i,O) + -$O,l,O) 3
3
>
(3.9) + 12 (
1 -$(l,O,O)
+4$$
+
E(l,O,O)
+4g(i,
rl
1
1 i,O) + -$O,l,O)
i,O) + $(O, 2
>
l,O)
I
3.1. (3.7)-(3.9) are Simpson’s formulae for the integrals on the left hand side. It is exact because one can expect ZJto be a polynomial of degree 3. In fact, the integral can be expanded as a series in derivatives, yet the final results will be consistent.
REMARK
From Section 2, the space of interpolation conditions with respect to w( At, AZ, A3), (3.6)-( 3.9) is corresponding to a space H of functions analytic at the origin (0, 0,O) in which each member is of the following form, where ci (i = 1, 2,. . . ,9) are constants, Ct en1 +c2 eA2+c3 eAJ+qet(Al+Az) +cs ef(A2+h)+cg ez'(h+Al)
f=
+c7(r2A1 +rlA2 +t3A3)(eAl+4et(Al+A2) +e*2) + cS(r3A2 + r2A3 + tlAl)(eA2 +C9(QA3
=
(Cl
+
+4ei(A2+ha)
+e*3)
+ rgAl +t2A2)(eA3 +4ef(A3+Al)+e"l) C2 + c3 + c4 + c5 + cg)
+ (Cl + k C4 -f- f C6 •k 6r2c7 •k 6tlc8+6r3c9)A1
+ (C2 + iC4
+ kc5
+ 6rlc7 + 6r3c8+6t2C9)A2
+(C3+~c5+~c6+6t3C7+6r2c8+6rIc9)A3 + ( ; Cl + $ C4 + f C,5 + 3r2C7 + 3qcg) A:
+ (i
C2 -k $ C4 +
1 1 + ( - C3 + - C5 2 8
f C5 + jr-1 C7 + 3ryZg)
1 8
hi
+ - C6 + 3r2C8 + 3rlcg)A$
98
J.B. Gao. TM. Shih/Comput. Methods Appl. Mech. Engrg. 122 (1995) 93-103 + (;
c4 - 3t3C7 + 3tlc8
+ 3t2Cg)AlA2
+ (;
cg + 3t3C7 - 3tlc8
+ 3t$9)&&
+ (;
Cfj + 3t3C7 + 3tlcs
- 3t$g)&hl
+ (;
Cl + &
c4 +
f
c6
+
r2c7
+
r3c9)A;
+(;c2+-&4+-&S+r3Cg+rlc7)$
+(~c3+-&+-$6+rlC9+r2c~)~~
+(~C4+(rl+~2)c7+t2C9)A:A2+(~CS+(12+r3)C8+t3C7)A~A3 +(kC6f(rl
+r3)c9+tlc8)&1
lkr2)C7$tlCg)hlh;
+($Cq+(rl
(3.10)
+(~c~+(r2+~3)c~+t2c9)~2~:+(~c6+(r~+r3)c9+t3c7)h3~~ +(t3C7+tlcg+t2Cg)hlA2A3+..*. Let
the
coefficients
of the quadratic
CI = -6r2c7
+ 6tlcg
-
6r3c9,
c2 = -6r3cg
+
6t2c9
-
6rlc7,
c3 = -6rlc9
+
6t3c7
-
6r2Cg,
c4 = 12t3c7 -
12t1cfj -
C5 =
12tlc8 - l&C9
-
12t3C7,
12t2c9
-
12t,c8.
Substituting +8+2C9)A,A;
-
12t3c7
(3.11)
be zero, then we have the following
relations
(3.11)
12t2c9,
c(j
=
items in (3.10)
into (3.10)) we get +sA:+4c9A$++9A;
+(-$C7-~t*C*+~Cg)A;A2
+(~C7-~C~-~t2C9)A;A3+(-~t3C7-~c*+&C9)A2A~
+(-$t3q+$
+(-?c7+
Q - &
C9)&
+(~C7-~t,C8+9)~3h~+(t3c7+t~c~+t2c9)~,A2~3+”’.
Hence by Theorem
2.4,
p=HI
=7T2(Al,A2,A3)
${Cl(2A;
- A;A2 -
A,A;
+ A;A3 - 3A2A; - 3A;A,
+ A3A; +4AlA2Aj)
+
- 3A3A; +4AlA2A3)
+
- AsA; +4AlA2A3)
:
C2(2A:
- 3A;A2 + AlA;
- A;A3 - A2A; + A;Al
c3(2Az
+ AfA2 - 3AlA;
- 3AiA3 + A2A: -
AiA,
ci E R}. Thus, F =
HI is the shape function space of de Veubeke’s element with Al + A2 + A3 = 1. It is not difficult
to prove that P=r2(h,A2,h3)
@{clA,(2A;+
A2A3) +
czA2(2hi
+
AlA3)
C3h3(?&
+
hlh2)
+
:
ci E I?}. Shi [ 131 defined the shape function
of de Veubeke’s element as
J.B. Gao, T.M. Shih/Comput. Methods Appl. Mech. Engrg. 122 (1995) 93-103 p =r2(h,h2,A3)
@{Cd&2
-
&A;
+ CAlA2A3)
+
C2( &3
-
h24
+ Chl A2A3)
+
C3(AiAl
-
h3Af
+ CAlA2h3)
99
:
Ci E R}
with a certain constant c $0. It can be proved that 7 $ P. This element, (K, P, (3.6)-( 3.9) ), passes the patch test because it passes the strong Fl and strong F2 tests (cf. [91>.
4. Interpolation
construction
of Zienkiewicz’s
element
Assume that all conditions in the previous section are satisfied. Zienkiewicz’s nine interpolation conditions are U(Pi), Let u(x,y)
u,(Pi),
Uy(Pi),
E w(At,A2,As)
d&y)
= &
t$(X,Y)
= -&
i= 1,293.
by (3.3),
v,$
1
srg 1
u(P2)
U(PI) = w(l,O,O), G(P1)
= &
[
ux(P2)
= &
+(P2)
= -&
Kc(P3)
= &
uy(P3)
= -&
[
[
I
conditions
1 u(P3)
+f2~(lr0,0)
190)
1
= w(O,O,
+s,~(l,o,o)
+1)3-gy0r1,0)
+52-$0,1,0> 2
+t,-$)o,
+772g(o,o,l) 2
+93-g(o,o,1)
~,~(O.O.
1) +~,-$(o,o,l)
+t3$o,o,
1
I
What is the shape function W(pl23)
2
space P of Zienkiewicz’s -2&W(Pi)
i=l
3
element?
-kkDW(Pi)(Pi-P123)
1 1 I I I I ,
,
l,O>
111-$(0,0,1) 1
l),
,
+I?~(l.o.o)
+7?2~(o,1,0) 2
I
(4.2)
(4.1) can be changed into the following:
= w(O,l,O),
1
&-$(O,
,
.
+t3E
~,~(l,O,O)
&$(O,l.O)
[
3I
+772-$(1.010) 2
W$(1,0,0)
u!.(f?) =-A [
then
+RE
In terms of (4.2)) the interpolation
The
(4.1)
+72$+7)3-$ 2
1
element is also triangular.
(4.3)
,
,
1)
.
Ciarlet (cf. [ 141) defined P as follows
=O
(4.4)
i=l
where Pi23 is the interior centroid of the triangle with vertices at PI, P2 and Ps. But this constraint is complicated. Now, by the new interpolation method, we can give a new shape function space ?r depending only on (4.4). As in the previous section, associated with interpolation conditions (4.4)) an analytical function subspace with elements of the following form can be constructed. h(At,A2,As)
=CteA’+C2eA2+CseAi + C4(51A1 + &A2 +&As)
e** +Cs(771At + 72A2 + 773A3)e*l
+ G(&AI
+ &A2 + &As) eAz+C7(771At + 772A2+ 73A3) eA2
+ W&AI
+ &A2 +&As)
cAs +G(vIAI
+ 712A2+ 773A3)eAS .
(4.5)
100
J.B. Gao, T.M. Shih/Comput.
Methods Appl. Mech. Engrg. 122 (1995) 93-103
Expand h( At, AZ, A3) at (0, O,O> as a power series and let the coefficients and without loss of generality, assume 53 SO, then we have c,
= -c2
c2
=
2
-
of items with degree < 2 be zero,
c3,
h’3 - %d3
c7
+
24
53 c3 = -253cs
c8
53
I
252r12 ,n9
63
’
- 2’313%
(4.6)
Now the cubic items of h( Al, AZ, A3) become
where A= $(&Cs+~lCs),B Denote
=
and C = ~(&C~-!-TIC~).
;(&2C8+772C9),
~=T~@{A(;A;+A,A;-;A~-A~A~)+ B( ;A; + A2A; - -A ;:
- A2A;) +
C(;A;+A,A;--A
-A2A:)
A, B,C =7r2 @{Cl
;:
:
E R}
(A, - A3j3 + C2(A2
-
AlI3
+
C3(A3 - A2)3 : Ci E R). Then the interpolation problem ( (4.1) , p) is correct, but the finite element (K, p, (4.1) ) is just the Bergan’s modified element, hence it can not pass the patch test unless the net lines of the triangulation are parallel to three different slopes. In [ 151, one of the authors used the method proposed here to give some other new finite elements for plate bending problem.
J.B.
5. Interpolation
Gao, TM.
construction
Shih/Comput.
Methods Appl.
101
Mech. Engrg. 122 (1995) 93-103
of Specht’s element
As Zienkiewicz’s element fails in passing the patch test, Specht (cf. [ lo] ) proposed a modified displacement basis based on interpolation test which is equivalent to the strong Fl test and F2 test in [ 91, i.e., finite elements satisfying SF, Auds = 0 and SF, A$ ds = 0 (cf. [ 10, p. 7091) for interior boundaries of the element. The interpolation conditions (or nodal parameters) of Specht’s element are the parameters of Zienkiewicz’s element, i.e. (4.1). Thus, if, on each triangle K = Pt P2P3, s?’ u ds and J?’ & ds are determined only by the nodal parameters propose the interpolation U(Pi>,
Ux(fi),
on PiPi+ then the interpolation test will be satisfied. Therefore, problem on K with the following interpolation conditions
Uy(pi>,
ly’uds,
f’
Let n-4 be the subspace of polynomials THEOREM
5. I. The interpolation
Consider
the following
computed,
respectively,
sub-interpolation
p, ,udS=+[u(P~) IP,
(5.1) we have:
(or4, (5.1) ) is correct. problem, i.e., interpolation
formulae:
I
I
+U(Pi+l)]
i= 1,2,3.
of order 4, then, obviously,
problem
by the following
&ds,
what we need is to
+ %
s(s)
condition
L 1,
- g(Pj+,)
$‘+I u ds and Jp’
&
ds are
(5.2)
and (5.3) By Theorem 5.1 a nine-dimensional subspace P of 7~4 according to conditions U(Pi), u,( Pi), uy( Pi) (i = 1,2,3), (5.2) and (5.3) can be determined. It is not difficult to prove that the subspace P thus obtained is the Specht’s modified shape function space. In the above discussion we may note that the computation of J; uds is not important while u( Pi) are proposed. For the finite element function u, we readily proved that [‘I Au ds = o( hllu112) if u is continuous at Pi and Pi+, , and U( Pi) and u( Pi+1) are element parameters, thus this element passes the Fl test (cf. [ 91) . Now let us consider the following interpolation conditions, P,, I U(Pi)*
Ux(Pi),
s P, One has to find a subspace Sections 3 and 4, we have THEOREM Q(K)
@{C3(r2hI CI (r3A2
i= 1,2,3.
(5.4)
Q of 7r4 such that (Q, (5.4) ) is correct. By the methods similar to those used in
5.2. The interpolation =T3
au
EdS, I I
Uy(Pi),
+
problem rlh2
+ r2A3
, (5.4) ) is correct, where
(Q(K)
+
t3h3)($
+
+
tlAl)(Az
+ Ai + 3AlA2A3) +
Cz(rlA3 + rnAl + tzA2)(Ai
Ai + 3A1A2A3)
+ Ai +3AtA:!As)
+
:
where t1C1 + t2C2 + t3C3 = 0 and Ci E R}. At last by (5.3) computing the interpolation conditions J:+’ & ds, a nine-parameter subspace F of Q is obtained while this P associated with the nine interpolation conditions A = {u( Pi), u,( pi), uy ( pi>, i = 1,2,3} is unisohble. The finite element (K, P, A) will pass the Fl test and the strong F2 test, and hence it converges for fourth order problems.
J.B. Gao, TM. Shih/Comput.
102
The dimension D(u)
of Q(K)
Methods Appl. Mech. Engrg. 122 (1995) 93-103
is 12. Let
(~(~>~~,(~l>~~,(~l)~~(~2)~~,(~2)~~y(~2),~(~3~,~,(~3),~~t~3)r~lO(~)rd,,(U),d,2(U))t
=
=(dl(u),...,dl2(U))I,
where
40(u) = By (5.3)
I
p2au
p,
-dds,
dt,(u)
an
=
I
4 au
-dds, p2 an
d,2(u)
=
I
pt au
- ds. p, afl
we have
D(u)
= Ct2x9D(u),
(5.5)
where C is a 12 x 9 matrix and b(u) Let U= /%A1 +
P2A2
+PII(P~~I
+P3A3
+P4hh2
+/%A2A3
. . . ,dg(u))‘. +&A:Az
+/%A3&
+&A;As
+/%A;A,
+~,oAtA2As
+0A2+t3A3)(A:+A:+3AtA2A3)
+/%2(r3A2+rzAs +Pts(rtA3
= (dt(u),
+t~At)(Ai+
+ r3At + tzAz)(Ai
A: f3A1A2A3)
(5.6)
+ A: + ~AIAzA~),
where +
@II
tlPl2 +
Let X= (PI,/32,.
t2P13 = 0.
(5.7)
. . , /?ts)‘. Substituting
where G is a 13 x 13 nonsingular
(5.6) and (5.7) into (5.5) one gets
matrix. Consequently,
n={u(Pi),u,(Pi),u,(Pi),
a new finite element can be defined as follows:
i= 1,2,3},
and defined by (5.6)
and (5.7) and X = G-’
(C;(n))}.
The finite element (K, P(K) , A) thus obtained passes the Fl test and strong F2 test (cf. [ 9]), convergent. We will refer to this element as the modified Specht’s element.
hence it is
Acknowledgement The authors to the referee Committee of Foundation of by HUST.
gratefully thank Professor Hughes of Stanford University for his many kind helps. Thanks also for his good ideas and comments. This project has been partially supported by the Research The Hong Kong Polytechnic University. J.B. Gao is supported by the Postdoctoral Science China, the Natural Science Foundation of China and the Science Foundation for Youths provided
References 11J 121 131 [41
Z.C. Shi, The generalized patch test for Zienkiewicz’s triangles, J. Comput. Math. 1 ( 1984) 279-286. f? Lascaux and P. Lesaint, Some nonconforming finite elements for plate bending problems, RAIRO Anal. Numer. 9 ( 1975) 9-53. S.C. Chen and Z.C. Shi, Double set parameter method of constructing stiffness matrices, Numer. Math. Sinica 13 ( 1991) 486-496. Z.C. Shi and SC. Chen, A direct analysis of the 9-parameter quasi-conforming plate element, Numer. Math. Sinica 12 ( 1990) 76-84. [ 51 H.Q. Zhang, The generalized patch and 9-parameter quasi-conforming element, in: Feng Kang and J.L. Lions, eds., Proceedings of the China-France symposium on finite element methods (Science Press, Gordan and Breach, 1983) 566-583.
J.B. Gao, TM. Shih/Compul. Methods Appl. Mech. Engrg. 122 (1995) 93-103
103
161 G.F? Bazeley, Y.K. Cheung, B.M. Irons and O.C. Zienkiewicz, Triangular elements in plate bending, conforming and nonconforming
solutions, Proceedings 1st Conference on Matrix Methods in Structures Mechanics, AFFDL-TR-CC-80 (Wrig Patterson Air Force Base, Ohio, 1986) 547-576. [7] R.L. Taylor, J.C. Simo, O.C. Zienkiewicz and A.C.H. Chan, The patch test - a condition for assessing FEM convergence, Int. J. Numer. Methods Engrg. 22 (1986) 39-62. [S] E Stummel, The generalized patch test, SIAM J. Numer. Anal. 16 ( 1979) 449-471. ]9] Z.C. Shi, The F-E-M-Test for convergence of nonconforming finite element, Math. Comput. 49 (1987) 391-405. [ IO] B. Specht, Modified shape functions for the three-node plate bending element passing the patch test, bit. J. Numer. Methods Engrg. 26 (1988) 705-715. ] I I ] C. deBoor and A. Ron, On multivariate polynomial interpolation, Con&. Approx. 6 (1990) 287-302. [ 121 C. deBoor and A. Ron, Computational aspects of polynomial interpolation in several variables, Math. Comput. 58 ( 1992) 705-727. [ 131 Z.C. Shi, A unified formation of shape functions for two kinds of nonconforming plate bending elements, Numer. Math. Sinica 8 ( 1986) 428-496. 1I4 ] PG. Ciarlet, The finite element method for elliptic problems (North-Holland, Amsterdam, 1978). [ 15 1 J.B. Gao, A new class of nonconforming finite elements, Numer. Math., J. Chin. Univ. 2 ( 1993) 186-194.
Computer Methods in Applied Mechanics and Engineering 122 ( 1995) 93-103
Interpolation methods for the construction of the shape function space of nonconforming finite elements J.B. Gao”***‘, T.M. Shih b a Department of Mathematics, Huazhong University of Science and Technology, Wuhan 430074, People’s Republic of China h Department of Applied Mathematics, The Hong Kong Polytechnic University, Kowloon, Hong Kong
Received 1 December 1993; revised 1 April 1994
Abstract Based on the multivariate polynomial interpolation theory, a method is derived to construct shape function spaces which can pass relevant convergent tests. As applications of this method, the shape function spaces of de Veubeke’s, Zienkiewicz’s and Specht’s elements are obtained. It is also shown that this method can be used to construct new finite elements.
1. Introduction Galerkin’s finite element method is a simple and efficient way for numerically approximating the solution of elliptic problems, especially for solving problems of continuum mechanics possessing infinite degree-offreedom. In order to solve a fourth-order boundary value problem such as the plate bending problem with the clamped boundary conditions, we have to choose higher degree polynomial finite elements in virtue of C’ compatibility conditions. The difficulty caused by higher degree elements can be surmounted by the socalled microelement method with splitting procedure, such as the Hsiech-Clough-Tocher splitting triangle and the Powell-Sabin triangle, which result in more complicated computation. It is well known that the nonconforming elements were introduced, such as Zienkiewicz’s triangle (Co class) (cf. [ l] ), de Veubeke’s triangle (cf. [ 2]>, the TRUNC element and the quasi-conforming plate element (cf. [ 3-51). because the requirement of conforming conditions is very harsh. Unfortunately, the nonconforming element methods cannot always ensure the convergence of approximated solutions. Irons introduced the “Patch Test” in continuum mechanics to conclude the convergence, which was one of the greatest contributions to the ‘theory’ and ‘practice’ of the finite element method (cf. [ 6,7] ) . Since then some other techniques have been developed to check sufficient requirements for convergence such as the “Interpolation Test”, the “Generalized Patch Test” (cf. [ 81) (sufficient and necessary test), and the “F-E-M-Test” (cf. [ 91) . In general, interpolation test and patch test can be used to develop useful, accurate and convergent nonconforming finite elements (cf. 16,101). Since test conditions may be viewed as interpolation conditions, the construction of nonconforming finite element is actually an interpolation problem. The problem is how to construct the associated shape function space. For example, Specht (cf. [ lo] ) wrote ‘The required three higher terms are assumed to be linear combinations of the following cubic and quartic terms: (f&2, $$$s, &,$i, ,t$&&, * Corresponding author. ’ Currently visiting The Hong Kong Polytechnic University. 0045-7825/95/$09.50 @ 1995 Ekevier Science S.A. All rights reserved SSDIOO45-7825(94)00728-7
94
J.B. Gao, TM. ShihlCotnput.
Methods Appl. Mech. Engrg. 122 (1995) 93-103
and St&&. This assumption is successful.“. However, it is not at all obvious how Specht gave the shape function space for his element. Another example is the famous de Veubeke’s triangle which passes the patch and interpolation tests. The nodal parameters of de Veubeke’s element are the function values at vertices of the triangle and the means of function values and first derivatives in outward normal direction on each edge of the triangle. How do we choose a subspace of 7rs as shape function space? We should find a method to determine some polynomial subspace (shape function space) from which the interpolation is possible and unique in terms of the interpolation conditions (such as the patch or interpolation tests or F-E-M-Test conditions). The recent work on multivariate polynomial interpolation is a useful for approaching this problem (cf. [ 11,121) . The following problem is thus interesting and important: how do you construct a subspace of polynomials so that the shape function spaces for given nodal interpolation parameters and the corresponding finite elements can pass the relevant tests? Here we aim to derive a method which can answer the above question. As examples, de Veubeke’s element, Zienkiewicz’s element and Specht’s element, in Sections 3, 4 and 5, respectively, are constructed by the new interpolation method. A new triangular element, refered to as the modified Specht’s element, is obtained. (,&$3,
2. The fundament
of interpolation
methods
The problem of multivariate polynomial interpolation is a complicated and classical problem. The generalization of univariate polynomial interpolation to the multivariate context is very difficult. Most efforts went into identifying sets of points @ C RS on which one can interpolate uniquely from rk (polynomial space of total degree < k) for some k. Our method comes from another direction. We wish to determine a certain corresponding polynomial space Q directly for any given interpolation conditions. Although such a space Q is not unique (cf. [ 1 l] ), it is enough for us to choose just one of them, with some criterions such as a lower degree about Q. An interpolation problem can be presented in the following abstract way. Let H be a finite dimensional linear subspace of a linear space X and n be a finite dimensional linear subspace of the dual of X (over C or R). The interpolation problem on H can be described as follows: Given f E X, find an h E H such that Ah = Af,
V’A E A.
(2.1)
If there is, for each f E X, exactly one solution h E H such that (2.1) is valid, the interpolation problem (H, A) is said to be correct, where A and H are the (space of) interpolation conditions and the space of interpolation for the problem (H, A), respectively. Recall (without proof) the following well-known characterization of correctness. LEMMA 2.1. Let H and A be$nite dimensional linear subspaces of a linear space X and its dual, respectively. The following three statements are equivalent: (I) The interpolation problem (H, A) is correct. (2) The linear map TH : h + (Aj h)? is one-one and onto, where ( Aj); is any basis of A. (3) A % H’ in the sense that the map F* : h --f AI* is one-one and onto. In the rest of this paper let X = Ao, the space of all functions analytic at the origin, the collection of all polynomials on RS. Our main interest is the polynomial interpolation interpolation conditions in linear functional form as follows:
[BIP(D,) : f + (p(D)f)(N
and denote by r problem with all
Vf E Ao,
where[B]f=f(B)foranyf,D=(Di ,..., D,),cr=(a(l) ,..., a(s)),andD*=D;Y(‘) ,..., D?“). Associate each point 0 in RS with a polynomial space Pe and define the space of interpolation conditions A := I{]%
(2.2)
: P E PO},
Vq E 7~ and f E Ao, define a linear functional
as
in the following
form:
95
J.B. Gao, TM. Shih/Comput. Method.? Appl. Mech. Engrg. 122 (199.5) 93-103 q*(f)
= (s(D>f)(O)
Then the following
FfW~!
= p%(O).
(2.3)
equality holds,
bbq E rt
[dp(D)q = q*(eep),
(2.4)
where ee : x E RS -+ ee.x. The map q --+ q* is linear and linear embedding of rr in the dual of Ao. A polynomial space Q is to be determined for which H = C es ~0 c Ao, then the interpolation problems (Q, A) From the problem (H, Q*), Q* can be used to represent Lemma 2.1. This means that it is enough to find H’. THEOREM 2.2. [I 1,121 H’ Z HI with HI =
( fl
one-one
in the sense of (2.3),
hence it provides
a
(Q, A) is correct. Let Q* = {q* : q E Q} and and (H, Q*) are equivalent in the sense of (2.4). H’, the dual of H, and conversely according to
: Vf E H}, where f 1 is the homogeneous polynomial
of
largest degree j for which f(x)
=fl(x)+o(IIxllj),
asx+O.
EXAMPLE 2.3. Let f (x, y) = x2 + xy + y3 + x*y + xy2, then f 1(x, y) = x2 + xy. We have the following
procedure
f =clfl
+c2f2+...+c,f,
=
L(cl,C2,...,C,)~xn
c
of HI. For any f E H, let { fj};
for the construction
be a basis of H, then
a
where La(ct,c2,. . . , c,) are linear forms of ct, c2,. . . .c,. Denote Ia\ = a( 1) + a(2) + . . . + a(s). Then we get the following THEOREM 2.4. [II,121
result.
Assume that above conditions are satisfied, then there exists an integer mg for which
~La(cl,c2,...,cn)x”
: L,(c1,c2,...,cn)=0,
EXQMPLE 2.5. For interpolation conditions n = span{f(xu), find a polynomial subspace Q such that (Q, A) is correct. Define a finite dimensional
space associated
lal
f’(xl),
f(x2))
(x0
< x1 6 x2
and x0 $x2),
with A as
H = {cl eXO*+c2x exlx +cg exzx : ci E R}. For any f E H, let f = ci exox +c2x
e”” +c3 exzx = cl + c3 + x( xocl + c2 + X2C3) +
x2
+ (;x;c,
XIC2 +
;x;c, )+2 (
1
3
,x&l
+
;x;c,
+
If ct + cs $0, then fl = cl + c3 = c (constant); hence 1 E HI. If cl + cg = 0 but xnct + c:! + ~2~3 SO, then fl = ((x2 - xo)c3 + c2)x; hence x E HI. If CI + ~3 = 0, and XOCI + ~2 + ~2~3 = 0, then f = ~3 (~2 - ~0) (T - ~1) x2 + (; (x;
Ifxt
+ x2x0 + xi, fl
$T,then
If xi = T, Dehne
-
fx?)
x3 +
=cg(x2-x0)(
. . .. T-xt)x2;hencex2
span{l,x,x2},
x1 #-
{ span{l,x,x3},
x1 =
Q=
x0 + x2
x03
x2 ’
2
’
>
c3( ~2 -
~0)
E H~.Hencema=2andH~=span{l,x,x2}.
and hence x3 E Hi. Consequently,
then fl =c3 -x3
+.
..
ixg3
mo = 3 and HI = span{ 1, x, x3}.
96
J.B. Gao. T.M. Shih/Comput. Methods Appl. Mech. Engrg. I22 (1995) 93-103
then (Q, LI) is correct.
3. Interpolation
construction
of de Veubeke’s element
K with vertices at 9 = (xi, yi), i = 1,2,3,
Given a triangle 51 =x2--3,
(2=x3--1,
53=x1
--X2,
vl=Y2-Y3,
r)2=Y3-Y1~
773=Yl
-Y2,
~:=s:+&
z;=s;+r1;,
denote by A the area of K and put
(3.1)
I;=&+&
(3.2)
Let ri and IZ~ (i = 1, 2, 3) be the unit tangent and outward normal on the edge Pi+tPi+Z opposite vertex pi (i = 1, 2, 3), respectively, and Wi = Xi+tyi+z - xi+zyi+t (i, i + 1, i + 2 E 23). NOW let Ai denote the area coordinates related to the vertices Pi, i.e., x = XI At + X2A2 + Y =
Yl Al +
to the
X3A3,
(3.3)
Y2A2 + Y3A3,
1 = A1 + A2 + A3, such that the triangle K is transformed into the standard simplex K* = { (AI, AZ, As) 1 At +A2+A3 = 1, Ai > 0). By (3.3) any function u(x,y) defined on K can be associated with a function w(At, A2,As) defined on K* such that u(x, Y> =
44,
(3.4)
A29 A3).
Now it is not difficult to prove that Zig=-: 1 +--I
at4
tig+ri-l (
aw
dAi_1
s+ri+l-
I
lfl
aw
dAi_ 1>
9 (3.5)
aAi+l'
where i,il,i+ 1 E Zs. De Veubeke’s element (cf. [ 21) is triangular, which satisfies the F-E-M-Test (cf. [ 91) (also passes the patch test (cf. [ 61) ) that ensures convergence of the element. The nine nodal parameters are the function values at the vertices Pi and at the middle points Mi = i (Pi-1 + Pi+1), and the mean values of first order derivatives in outward normal direction on each side of K. There are nine interpolation conditions and in general the one has to give the other linearly-independent interpolation condition so that the interpolation problem (7r3, LI) is correct; II is the interpolation condition space of ten dimensions mentioned above. Some authors had given the methods adding the tenth interpolation condition. For example, Shi showed in [ 131 that the condition given in [ 141 is linearly dependent with de Veubeke’s parameters (or interpolations) on r. Now a shape function space B of de Veubeke’s element can easily be obtained with the new interpolation approach without the tenth condition. With the transformation (3.3) the nine interpolation conditions about u (defined on K) can be transformed as, with respect to w on K*, U(Pl)
= w(l,O,O),
u(P2)
= w(O,l,O),
u(P3)
= w(O,O, l),
(3.6)
J.B. Gao, TM. Shih/Compur. Methods Appl. Mech. Engrg. 122 (1995) 93-103
1
pJ
au
-ds=-& II s pz ant
1
tt [( +
$0,
r3
(
pl I ps
G
au
-ds=-&-
2
an2
t2 I( + rl
+ r3
2
1,O) +4g(o,
;) + $(O,O,
1)
2
1)
1) + -$(l,O,O) 2
g?o,o, 3
1) +4$+,0,
-g(O,O,
1) +4$$0,$
3
(3.7) >
>
(3.8)
;) + $(l,O,O) 3
>
+ -$l,O,O)
1
1
p2at4
1 13
>
;, 5 + -$o,o,
I> +4g(;.o,
1
1)
1
+4$0&
g(o,l,o) 2
g(o,o,
(
;, ;) + $(O,O, 1
-g(O,
+ 12
1
I,O) +4$0,
I
97
Jp,
-d,=-&an3
$(l,O,O) 3
3 Ht3
+4$+,
i,O) + -$O,l,O) 3
3
>
(3.9) + 12 (
1 -$(l,O,O)
+4$$
+
E(l,O,O)
+4g(i,
rl
1
1 i,O) + -$O,l,O)
i,O) + $(O, 2
>
l,O)
I
3.1. (3.7)-(3.9) are Simpson’s formulae for the integrals on the left hand side. It is exact because one can expect ZJto be a polynomial of degree 3. In fact, the integral can be expanded as a series in derivatives, yet the final results will be consistent.
REMARK
From Section 2, the space of interpolation conditions with respect to w( At, AZ, A3), (3.6)-( 3.9) is corresponding to a space H of functions analytic at the origin (0, 0,O) in which each member is of the following form, where ci (i = 1, 2,. . . ,9) are constants, Ct en1 +c2 eA2+c3 eAJ+qet(Al+Az) +cs ef(A2+h)+cg ez'(h+Al)
f=
+c7(r2A1 +rlA2 +t3A3)(eAl+4et(Al+A2) +e*2) + cS(r3A2 + r2A3 + tlAl)(eA2 +C9(QA3
=
(Cl
+
+4ei(A2+ha)
+e*3)
+ rgAl +t2A2)(eA3 +4ef(A3+Al)+e"l) C2 + c3 + c4 + c5 + cg)
+ (Cl + k C4 -f- f C6 •k 6r2c7 •k 6tlc8+6r3c9)A1
+ (C2 + iC4
+ kc5
+ 6rlc7 + 6r3c8+6t2C9)A2
+(C3+~c5+~c6+6t3C7+6r2c8+6rIc9)A3 + ( ; Cl + $ C4 + f C,5 + 3r2C7 + 3qcg) A:
+ (i
C2 -k $ C4 +
1 1 + ( - C3 + - C5 2 8
f C5 + jr-1 C7 + 3ryZg)
1 8
hi
+ - C6 + 3r2C8 + 3rlcg)A$
98
J.B. Gao. TM. Shih/Comput. Methods Appl. Mech. Engrg. 122 (1995) 93-103 + (;
c4 - 3t3C7 + 3tlc8
+ 3t2Cg)AlA2
+ (;
cg + 3t3C7 - 3tlc8
+ 3t$9)&&
+ (;
Cfj + 3t3C7 + 3tlcs
- 3t$g)&hl
+ (;
Cl + &
c4 +
f
c6
+
r2c7
+
r3c9)A;
+(;c2+-&4+-&S+r3Cg+rlc7)$
+(~c3+-&+-$6+rlC9+r2c~)~~
+(~C4+(rl+~2)c7+t2C9)A:A2+(~CS+(12+r3)C8+t3C7)A~A3 +(kC6f(rl
+r3)c9+tlc8)&1
lkr2)C7$tlCg)hlh;
+($Cq+(rl
(3.10)
+(~c~+(r2+~3)c~+t2c9)~2~:+(~c6+(r~+r3)c9+t3c7)h3~~ +(t3C7+tlcg+t2Cg)hlA2A3+..*. Let
the
coefficients
of the quadratic
CI = -6r2c7
+ 6tlcg
-
6r3c9,
c2 = -6r3cg
+
6t2c9
-
6rlc7,
c3 = -6rlc9
+
6t3c7
-
6r2Cg,
c4 = 12t3c7 -
12t1cfj -
C5 =
12tlc8 - l&C9
-
12t3C7,
12t2c9
-
12t,c8.
Substituting +8+2C9)A,A;
-
12t3c7
(3.11)
be zero, then we have the following
relations
(3.11)
12t2c9,
c(j
=
items in (3.10)
into (3.10)) we get +sA:+4c9A$++9A;
+(-$C7-~t*C*+~Cg)A;A2
+(~C7-~C~-~t2C9)A;A3+(-~t3C7-~c*+&C9)A2A~
+(-$t3q+$
+(-?c7+
Q - &
C9)&
+(~C7-~t,C8+9)~3h~+(t3c7+t~c~+t2c9)~,A2~3+”’.
Hence by Theorem
2.4,
p=HI
=7T2(Al,A2,A3)
${Cl(2A;
- A;A2 -
A,A;
+ A;A3 - 3A2A; - 3A;A,
+ A3A; +4AlA2Aj)
+
- 3A3A; +4AlA2A3)
+
- AsA; +4AlA2A3)
:
C2(2A:
- 3A;A2 + AlA;
- A;A3 - A2A; + A;Al
c3(2Az
+ AfA2 - 3AlA;
- 3AiA3 + A2A: -
AiA,
ci E R}. Thus, F =
HI is the shape function space of de Veubeke’s element with Al + A2 + A3 = 1. It is not difficult
to prove that P=r2(h,A2,h3)
@{clA,(2A;+
A2A3) +
czA2(2hi
+
AlA3)
C3h3(?&
+
hlh2)
+
:
ci E I?}. Shi [ 131 defined the shape function
of de Veubeke’s element as
J.B. Gao, T.M. Shih/Comput. Methods Appl. Mech. Engrg. 122 (1995) 93-103 p =r2(h,h2,A3)
@{Cd&2
-
&A;
+ CAlA2A3)
+
C2( &3
-
h24
+ Chl A2A3)
+
C3(AiAl
-
h3Af
+ CAlA2h3)
99
:
Ci E R}
with a certain constant c $0. It can be proved that 7 $ P. This element, (K, P, (3.6)-( 3.9) ), passes the patch test because it passes the strong Fl and strong F2 tests (cf. [91>.
4. Interpolation
construction
of Zienkiewicz’s
element
Assume that all conditions in the previous section are satisfied. Zienkiewicz’s nine interpolation conditions are U(Pi), Let u(x,y)
u,(Pi),
Uy(Pi),
E w(At,A2,As)
d&y)
= &
t$(X,Y)
= -&
i= 1,293.
by (3.3),
v,$
1
srg 1
u(P2)
U(PI) = w(l,O,O), G(P1)
= &
[
ux(P2)
= &
+(P2)
= -&
Kc(P3)
= &
uy(P3)
= -&
[
[
I
conditions
1 u(P3)
+f2~(lr0,0)
190)
1
= w(O,O,
+s,~(l,o,o)
+1)3-gy0r1,0)
+52-$0,1,0> 2
+t,-$)o,
+772g(o,o,l) 2
+93-g(o,o,1)
~,~(O.O.
1) +~,-$(o,o,l)
+t3$o,o,
1
I
What is the shape function W(pl23)
2
space P of Zienkiewicz’s -2&W(Pi)
i=l
3
element?
-kkDW(Pi)(Pi-P123)
1 1 I I I I ,
,
l,O>
111-$(0,0,1) 1
l),
,
+I?~(l.o.o)
+7?2~(o,1,0) 2
I
(4.2)
(4.1) can be changed into the following:
= w(O,l,O),
1
&-$(O,
,
.
+t3E
~,~(l,O,O)
&$(O,l.O)
[
3I
+772-$(1.010) 2
W$(1,0,0)
u!.(f?) =-A [
then
+RE
In terms of (4.2)) the interpolation
The
(4.1)
+72$+7)3-$ 2
1
element is also triangular.
(4.3)
,
,
1)
.
Ciarlet (cf. [ 141) defined P as follows
=O
(4.4)
i=l
where Pi23 is the interior centroid of the triangle with vertices at PI, P2 and Ps. But this constraint is complicated. Now, by the new interpolation method, we can give a new shape function space ?r depending only on (4.4). As in the previous section, associated with interpolation conditions (4.4)) an analytical function subspace with elements of the following form can be constructed. h(At,A2,As)
=CteA’+C2eA2+CseAi + C4(51A1 + &A2 +&As)
e** +Cs(771At + 72A2 + 773A3)e*l
+ G(&AI
+ &A2 + &As) eAz+C7(771At + 772A2+ 73A3) eA2
+ W&AI
+ &A2 +&As)
cAs +G(vIAI
+ 712A2+ 773A3)eAS .
(4.5)
100
J.B. Gao, T.M. Shih/Comput.
Methods Appl. Mech. Engrg. 122 (1995) 93-103
Expand h( At, AZ, A3) at (0, O,O> as a power series and let the coefficients and without loss of generality, assume 53 SO, then we have c,
= -c2
c2
=
2
-
of items with degree < 2 be zero,
c3,
h’3 - %d3
c7
+
24
53 c3 = -253cs
c8
53
I
252r12 ,n9
63
’
- 2’313%
(4.6)
Now the cubic items of h( Al, AZ, A3) become
where A= $(&Cs+~lCs),B Denote
=
and C = ~(&C~-!-TIC~).
;(&2C8+772C9),
~=T~@{A(;A;+A,A;-;A~-A~A~)+ B( ;A; + A2A; - -A ;:
- A2A;) +
C(;A;+A,A;--A
-A2A:)
A, B,C =7r2 @{Cl
;:
:
E R}
(A, - A3j3 + C2(A2
-
AlI3
+
C3(A3 - A2)3 : Ci E R). Then the interpolation problem ( (4.1) , p) is correct, but the finite element (K, p, (4.1) ) is just the Bergan’s modified element, hence it can not pass the patch test unless the net lines of the triangulation are parallel to three different slopes. In [ 151, one of the authors used the method proposed here to give some other new finite elements for plate bending problem.
J.B.
5. Interpolation
Gao, TM.
construction
Shih/Comput.
Methods Appl.
101
Mech. Engrg. 122 (1995) 93-103
of Specht’s element
As Zienkiewicz’s element fails in passing the patch test, Specht (cf. [ lo] ) proposed a modified displacement basis based on interpolation test which is equivalent to the strong Fl test and F2 test in [ 91, i.e., finite elements satisfying SF, Auds = 0 and SF, A$ ds = 0 (cf. [ 10, p. 7091) for interior boundaries of the element. The interpolation conditions (or nodal parameters) of Specht’s element are the parameters of Zienkiewicz’s element, i.e. (4.1). Thus, if, on each triangle K = Pt P2P3, s?’ u ds and J?’ & ds are determined only by the nodal parameters propose the interpolation U(Pi>,
Ux(fi),
on PiPi+ then the interpolation test will be satisfied. Therefore, problem on K with the following interpolation conditions
Uy(pi>,
ly’uds,
f’
Let n-4 be the subspace of polynomials THEOREM
5. I. The interpolation
Consider
the following
computed,
respectively,
sub-interpolation
p, ,udS=+[u(P~) IP,
(5.1) we have:
(or4, (5.1) ) is correct. problem, i.e., interpolation
formulae:
I
I
+U(Pi+l)]
i= 1,2,3.
of order 4, then, obviously,
problem
by the following
&ds,
what we need is to
+ %
s(s)
condition
L 1,
- g(Pj+,)
$‘+I u ds and Jp’
&
ds are
(5.2)
and (5.3) By Theorem 5.1 a nine-dimensional subspace P of 7~4 according to conditions U(Pi), u,( Pi), uy( Pi) (i = 1,2,3), (5.2) and (5.3) can be determined. It is not difficult to prove that the subspace P thus obtained is the Specht’s modified shape function space. In the above discussion we may note that the computation of J; uds is not important while u( Pi) are proposed. For the finite element function u, we readily proved that [‘I Au ds = o( hllu112) if u is continuous at Pi and Pi+, , and U( Pi) and u( Pi+1) are element parameters, thus this element passes the Fl test (cf. [ 91) . Now let us consider the following interpolation conditions, P,, I U(Pi)*
Ux(Pi),
s P, One has to find a subspace Sections 3 and 4, we have THEOREM Q(K)
@{C3(r2hI CI (r3A2
i= 1,2,3.
(5.4)
Q of 7r4 such that (Q, (5.4) ) is correct. By the methods similar to those used in
5.2. The interpolation =T3
au
EdS, I I
Uy(Pi),
+
problem rlh2
+ r2A3
, (5.4) ) is correct, where
(Q(K)
+
t3h3)($
+
+
tlAl)(Az
+ Ai + 3AlA2A3) +
Cz(rlA3 + rnAl + tzA2)(Ai
Ai + 3A1A2A3)
+ Ai +3AtA:!As)
+
:
where t1C1 + t2C2 + t3C3 = 0 and Ci E R}. At last by (5.3) computing the interpolation conditions J:+’ & ds, a nine-parameter subspace F of Q is obtained while this P associated with the nine interpolation conditions A = {u( Pi), u,( pi), uy ( pi>, i = 1,2,3} is unisohble. The finite element (K, P, A) will pass the Fl test and the strong F2 test, and hence it converges for fourth order problems.
J.B. Gao, TM. Shih/Comput.
102
The dimension D(u)
of Q(K)
Methods Appl. Mech. Engrg. 122 (1995) 93-103
is 12. Let
(~(~>~~,(~l>~~,(~l)~~(~2)~~,(~2)~~y(~2),~(~3~,~,(~3),~~t~3)r~lO(~)rd,,(U),d,2(U))t
=
=(dl(u),...,dl2(U))I,
where
40(u) = By (5.3)
I
p2au
p,
-dds,
dt,(u)
an
=
I
4 au
-dds, p2 an
d,2(u)
=
I
pt au
- ds. p, afl
we have
D(u)
= Ct2x9D(u),
(5.5)
where C is a 12 x 9 matrix and b(u) Let U= /%A1 +
P2A2
+PII(P~~I
+P3A3
+P4hh2
+/%A2A3
. . . ,dg(u))‘. +&A:Az
+/%A3&
+&A;As
+/%A;A,
+~,oAtA2As
+0A2+t3A3)(A:+A:+3AtA2A3)
+/%2(r3A2+rzAs +Pts(rtA3
= (dt(u),
+t~At)(Ai+
+ r3At + tzAz)(Ai
A: f3A1A2A3)
(5.6)
+ A: + ~AIAzA~),
where +
@II
tlPl2 +
Let X= (PI,/32,.
t2P13 = 0.
(5.7)
. . , /?ts)‘. Substituting
where G is a 13 x 13 nonsingular
(5.6) and (5.7) into (5.5) one gets
matrix. Consequently,
n={u(Pi),u,(Pi),u,(Pi),
a new finite element can be defined as follows:
i= 1,2,3},
and defined by (5.6)
and (5.7) and X = G-’
(C;(n))}.
The finite element (K, P(K) , A) thus obtained passes the Fl test and strong F2 test (cf. [ 9]), convergent. We will refer to this element as the modified Specht’s element.
hence it is
Acknowledgement The authors to the referee Committee of Foundation of by HUST.
gratefully thank Professor Hughes of Stanford University for his many kind helps. Thanks also for his good ideas and comments. This project has been partially supported by the Research The Hong Kong Polytechnic University. J.B. Gao is supported by the Postdoctoral Science China, the Natural Science Foundation of China and the Science Foundation for Youths provided
References 11J 121 131 [41
Z.C. Shi, The generalized patch test for Zienkiewicz’s triangles, J. Comput. Math. 1 ( 1984) 279-286. f? Lascaux and P. Lesaint, Some nonconforming finite elements for plate bending problems, RAIRO Anal. Numer. 9 ( 1975) 9-53. S.C. Chen and Z.C. Shi, Double set parameter method of constructing stiffness matrices, Numer. Math. Sinica 13 ( 1991) 486-496. Z.C. Shi and SC. Chen, A direct analysis of the 9-parameter quasi-conforming plate element, Numer. Math. Sinica 12 ( 1990) 76-84. [ 51 H.Q. Zhang, The generalized patch and 9-parameter quasi-conforming element, in: Feng Kang and J.L. Lions, eds., Proceedings of the China-France symposium on finite element methods (Science Press, Gordan and Breach, 1983) 566-583.
J.B. Gao, TM. Shih/Compul. Methods Appl. Mech. Engrg. 122 (1995) 93-103
103
161 G.F? Bazeley, Y.K. Cheung, B.M. Irons and O.C. Zienkiewicz, Triangular elements in plate bending, conforming and nonconforming
solutions, Proceedings 1st Conference on Matrix Methods in Structures Mechanics, AFFDL-TR-CC-80 (Wrig Patterson Air Force Base, Ohio, 1986) 547-576. [7] R.L. Taylor, J.C. Simo, O.C. Zienkiewicz and A.C.H. Chan, The patch test - a condition for assessing FEM convergence, Int. J. Numer. Methods Engrg. 22 (1986) 39-62. [S] E Stummel, The generalized patch test, SIAM J. Numer. Anal. 16 ( 1979) 449-471. ]9] Z.C. Shi, The F-E-M-Test for convergence of nonconforming finite element, Math. Comput. 49 (1987) 391-405. [ IO] B. Specht, Modified shape functions for the three-node plate bending element passing the patch test, bit. J. Numer. Methods Engrg. 26 (1988) 705-715. ] I I ] C. deBoor and A. Ron, On multivariate polynomial interpolation, Con&. Approx. 6 (1990) 287-302. [ 121 C. deBoor and A. Ron, Computational aspects of polynomial interpolation in several variables, Math. Comput. 58 ( 1992) 705-727. [ 131 Z.C. Shi, A unified formation of shape functions for two kinds of nonconforming plate bending elements, Numer. Math. Sinica 8 ( 1986) 428-496. 1I4 ] PG. Ciarlet, The finite element method for elliptic problems (North-Holland, Amsterdam, 1978). [ 15 1 J.B. Gao, A new class of nonconforming finite elements, Numer. Math., J. Chin. Univ. 2 ( 1993) 186-194.