- Email: [email protected]
Finite elements satisfying all governing equations inside the element
Ei&ierscirnce
Copyright @J 1994 Ltd Printed in Gnat B&in. All ri&ts wcrwd a0457949p4 sd.00 + 0.00
Pergamon
FINITE ELEMENTS SATISFYING ALL GOVERNING EQUATIONS INSIDE THE ELEMENT V. KOMPIS Institute of Mechanics. Technical Unjversjty of Transport and ~ommunj~tions,
&%a,
Slovak Republic (Received 11 August 1992) introduce a formulation of the finite element method for the solution of elasticity problems in which all the governing equations will be satisfied inside the element. The displacements inside the elements are approximated by polynomials, whose coefficients are given from equilibrium equations. Another independent displacement field, which is common to adjacent elements, is defined on the element boundaries only. The element stiffness equations are derived using two modified variational principles. The element can have arbitrary shape and can be connected to isoparametric displacement elements and/or can be used as a transition element. Although the solution is shown for two-dimensional elasticity, the extension on three-dimensional and/or other field problems is straightforward. AlMract-We
1. INTRODUCTION
differentiation in corresponding direction and the Einstein’s summation is used for repeated indices. Assuming a homogeneous, isotropic material, the constitutive relationships are
Whereas in the standard isoparametric finite element method (FEM) the approximate solution is designed to satisfy interelement compatibility of displacements, attempts have been made to derive the formulations satisfying all governing equations inside the elements. Trefftz presented first a variational method[l] in which the approximate solution was designed to satisfy a priori the differential equations of the problem. From these principles the boundary element method [2] was formulated. For the FEM incorporation of the Trefftz principles new variational principles for interelement connection derived from more general principles [3] were employed [4,5]. The aim of the present work is to present a new method of presentation of displacement fields satisfying all governing equations inside the elements in polynomial form [6]. Two variational formulations for interelement connections are used and corresponding stiffness matrices for two-dimensional (2D) elastic solution are presented. The extension on three-dimensional (3D) and/or other field problems is straightfo~ard.
CFIl= A$, + f&z 022= BE,, + A~22 612- a2, = 2Gt,*,
(2)
where A = E/(1 - v*) B=vA
G=fE/(l
+v)
(3)
and E is Young’s modulus and v is Poisson’s ratio. The strain tensor tij may be expressed in terms of the displacement vector u, Cii= $(U$. j + uj, J.
2. GOVERNING EQUATIONS
The internal equilibrium is then expressed (from (2x4)) by the Lam&-Navier equations
The governing equations for 2D elasticity probequations
lems are the equiiibrium
Oii,$+ & = 0
(i*j = 1,2),
Au,, II +
(11
where (T,~is the stress tensor in Cartesian coordinates xi, 6, is the body force vector acting on a unit volume. The index after a comma denotes the
(4)
Fu2.12
+
G4.22
=
F~,.~+&~z+%I,=
-4
-b2,
(3
where F=BfG=fE/(l 273
-v).
(6)
274
KOMPIS
V. 1
Table
{Ai) = {a)‘), a!*‘, ui3’, . . . }‘,
(10)
1
*;
where T denotes transposed vector/matrix. Carrying out the second derivatives of the polynomial approximation of the displacements (8) and setting it into the Lam&-Navier equations (5) and comparing the coefficients of equal polynomial terms, some coefficients of {Ai} can be expressed by means of the others. We can see, that in this way all coefficients corresponding to inner terms in the Pascal triangle in Table 1 (e.g. the terms with the mixed product of coordinate components) can be expressed by means of coefficients of the outer (single component) terms. We can see also that some coefficients of the u, component of the displacement field are expressed by coefficients of its second component u2 and vice versa. Derivation of the dependent coefficients in each row in Table 1 from the independent ones could seem to decrease the efficiency of the procedure. However, we need to do this only once for each material, even for different degrees of interpolation polynomials, because this dependence is given by material constants only. Moreover an efficient procedure can be found if we realize the form of the matrices (some submatrices are diagonal the others are three-diagonal). If we chose independent coefficients of the third order to be the terms xix2 and x,x: and these of the fourth order x:x2 and x,x: (the solution for these terms was in this particular case simpler), we are able to write the polynomial satisfying all governing equations (body forces are not included) in the form
XI x2 _ XIX2 *f
x: XI-4 x:x* x; x: x,x: x:x: 4x2 d 3. FINITE ELEMENT SOLUTION
3.1. Displacement fields satisfying the governing equations inside the elements Let us suppose the interpolation form
polynomial in the
P(x) = [a”)]o + [amx, + a(3)x2],
+ [a’“‘x: + a”‘x, x2 + a(6)xi]2 + [a(‘)x: + a(*)xix2 + at9)x, xi + a(‘O)x&
+ [fJ”“xf + u”“x~x2 + a(‘3)xix: + ~2(‘~)x, xi + a(15)xj]4+
.. ,
(7)
where the lower index after square brackets denotes the order of corresponding polynomial terms. Thus the full polynomial of n th order will include the terms OfO-, l-,..., and n th order. In the matrix notation the ith component of the displacement field will be approximated as
4 = [Pl{4),
(8)
where [p] is a row matrix of polynomial terms [PI =[[110, [x,rx*l,r [x:,xIx2rx:L..
. .I
(9)
{a 1 = ]Bl{c, 19
and {A,} is a column matrix of its coefficients
where the coefficients of the matrix [B] = [B, , B,] are Table 2
By’ = x, B\4’
=
B’P’=
x;
By = x* /q5’
-x:A/(3G)
X:X,
By’ =x3*I 2 @y’ = /y/I) = qw
=
*;
B\” = x,x;
-x;G/(3,4)
By’ = x x3 = 0
B\“‘=
:2:,x,G,F
B\14’= -2x,x,A/F
B\15’= -x;
B\‘@= -xi F/(3B)
B’j”=(H,x:+x;A/~)/4-
B\‘*’ = (x:G/F
- l.Sx;x;A/F
+ H,x;)/4
B$” = Bi2’ = B$j’ = 0
F/(3A)
Bi6’ = -x;F/(3G)
B$“= -x;F/(3A)
Bf’=(x:H,+x;G/F)/4-
B$‘j’
zz x;
Bi15’
=
X:X,
8s”’ = x:.x2
+x:A/F)/4B$“’ =x, -
x;G/(3A)
l.Sx;x;G/F
Bi4’ = -2x,x,A/F
Bi5’= -2x,x,G/F BI”=(x:H, B$‘O’= 1
(11)
l.Sx:x;G/F B”2’ *
_ -x2
B\‘4’
=
x;
By@ = x,x;
- x;A/(3B)
B$‘*’ = x, x;,
where H, = (G2- F’)/(AF)
H, = (AZ - F2)/(GF)
l.Sxfx;A/F
FE satisfying all goveming equations inside the element given in Table 2 and {u) = {u,, u2) is the vector of displacements. The vector {c, 1 is a vector of unknown coefficients for an element. The number of unknown coefficients in this vector is equal to 2(1 + 2n), where n is the order of interpolation polynomial Ip]. The components of the deformation tensor will be given by
will be defined on all element boundaries, variation
I
where
{E)= (q,, 622,62fT [cl
=I&,,
(4,,+B2.1)/21~ (14)
4.2,
and substituting this into (2) we obtain the components of the stress tensor fcf = (G,,, oz2, ctzf by
_ ri
by PI =
tr; 1 -v2
--L 0
1
0
0
1
0
I
(1 -v)/2
StT(v, - ~7)dT = 0
(19)
will enforce the continuity of displacements over the whole element boundary r, (it must be satisfied for each element separately). Moreover the variation &*{t,-I)dT I r,
=0
(20)
>
will enforce the continuity of tractions. r, and T, are the parts of the boundary where displacements and tractions are given respectively. r, is the interelement boundary. So the interelement equilib~~ (continuity of tractions) between the adjacent elements is given by the first integral (20) and thus this variation also establishes that the equilbrium on the whole will be satisfied. If the independent displacement field fii is approximated by means of shape functions NW and displacements U@ in nodal points, the functions &will be C, continuous on continuous parts of the element boundaries, i.e.
(16)
and from it also the tractions fi = aijnj,
and the
f
Wt c dr + or in the vector form
275
(17)
where nj is a component of a unit normal vector of the surface. In vectorial form the tractions are {t) = {f,* t,JT
for each boundary segment of integration. In this notation each boundary segment of an element has the same meaning as an element by the boundary element method (BEM). W is the vector of displacements in thejth nodal point of the element boundary. Defining now
3.2. The element formulation The displacement field (11) will be generally incompatible between the elements and the tractions which will be defined from these fields. Variational principles for general discontinuous fields were given by Prayer [3]. If we define a new boundary displacement field I$, which will be completely independent of ui (and consequently of the related boundary vector vi), two successful approaches can be found by variational formulations defined by Jirousek [S]. The two independent fields are imposed to fuifil the compatibility of incompatible fields in an integral sense. By the first approach (A) the independent field ci
From eqn (19) we get the dependence between the
V. KOMPIS
276
coefficients of interpolating polynomials {c,} and of the nodal degrees of freedom {U?} in the form =
(26)
W,l{c,~ = [G,l{ve)
s r,
and eqn (20) gives the relation
+
KKl%J = Fe 1
qce 1’
(27) =
and substituting
--B{cyJT [B]qdT
s ru
(bIDI KITfi dr
g{cyJT_(pc*)
(35)
{c,} from (26) into (27) we get
(28)
KlVW’~~A~~,) = {pr)
GdTt,dT = ,(,}’ i
I r,
r,
which can be written in the form
[NIT([nl[Dl[CI)dT(c,J
=8{~cWc*F{ce).
(36)
We obtain from (31) where
]Hfl{c,j
&I= [~elr[fW’[~,l
(30)
is the stiffness matrix of the element. The second approach (B) will be based on the definition of the independent displacement field 6, on the interelement boundaries only. The variation -
=
[Woe,) + ICI
and from (32) [GflT{c,} = 0. Substituting
[GflTW:I~‘[Gtl~~~~= -[G:lTWbl-‘{P:
6tr(o,-6)dT
1 (39)
I r. +
(38)
(37) into (38) we get
St:@, - a) dI-
6o,T(t, -7) dT +
(37)
=0
(31)
Kl{ V, )= If?)? will enforce the continuity of boundary conditions (if the boundaries T, and TU are parts of the boundary of the element f,) and the continuity of displacements on interelement boundaries. Further the variation
(40)
where
[&I= [Gfl~fW’[G:I
(41)
is the element stiffness matrix and Srt, dT = 0 will enforce the continuity Let us denote now 6t:v,
s r.+r’, = S{cy}’
of tractions.
GvJt,
dT -
dl-
s r,
s
r +r (blPl[WPld~~c,~ Y f
(32)
{P:} = -[G:lrIHfl~‘{P:}
(42)
is the vector of static equivalent loads. We can see, that both A and B formulations give the same element matrices for elements not containing boundaries I’, and r, The second formulation, B, does not contain nodal points on the rU and T, boundaries and so the matrices are built similarly as classical elements with reduced degrees of freedom on these boundaries. The number of nodes in an element cannot be greater than 2n + 1, where n is the order of the interpolation polynomial (according to Sec. 3.1). The element can be connected to compatible displacement elements which have the same shape functions along the adjacent boundaries. The shape of the element can be arbitrary and also the degree of shape functions along each boundary segment can be arbitrarily chosen. So it can be used as a transition element for the different meshes.
FE satisfying all governing equations inside the element
277
points. We can also compare the normal stresses on the loaded edge with the prescribed boundary tractions. The calculated values (related to the boundary tractions) varied from 0.75 to 1.08, while by the displacement FE fo~ulation from 0.24 to 1.89. A curved beam by the pure bending (Fig. 2)
Fig. I. FE model of the cantilever.
Formulation for 3D or other types of fields is straightforward.
The exact solution [7] containing non-polynomial terms is 10.364 on the inner radius and -9.658 on the outer radius. The B formulation gave stresses varying from 10.312 to 10.420 and from -9.612 to -9.768, respectively. A hole in a band (Fig. 3)
4. EXAMPLES
The following three examples will show some properties of the described elements. A five-point Gauss formula was used for integration along the element boundaries and the stresses were calculated in the Gauss nodal points everywhere. A cantilever loaded by a uniform pressure (Fig. 1)
If the boundary tractions were given to get the solution described by the polynomial of the fourth order, the fo~ulation B gives the exact solution for any shape of interelement boundaries. If the clamped end boundary conditions are given by a, = 0 and tr is quadratic as before, the exact solution is not known and it is not described by the polynomial of the fourth order. The maximal stress on the element boundary (in the closest Gauss nodal point of the clamped boundary) was 0.9468 of that obtained for the former boundary conditions. If the boundary was fully clamped (u, = 0, ur = 0), the stress on the boundary with unifo~ pressure was 0.9704 (related to that of the former solution) and 0.9186 on the opposite edge. Supposing the linear distribution of normal stress on the clamped boundary, we get the stresses 0.9288 in corresponding
(a)
Oniy one quarter of the region was modelled by FE, because of the double symmetry of the problem. The convergency of the fo~ulation was examined on this example. The theoretical stresses [7] are 0.75 and 4.30 in the points A and B, respectively (for the traction on the edge E-F equal to 1.0). Models consisting of two, three and eight elements are shown in Fig. 3 (the two-element model is the same as the three-node model, only the region A-B-D-C is modelled by one element). The stresses calculated in the point nearest to B were: 2.24, 3.65 and 4.38, respectively, for those models. The displacements u, were 0.00986, 0.00888 and 0.00876 in the point D and 0.00868, 0.00808 and 0.00811 in the point C for corresponding models. These results were obtained for the formulation A. A similar behaviour was observed for the formulation B.
5. CONCLUSIONS
The given FE formulation satisfying all the governing equations inside the element is very attractive, it can be incorporated into the existing FE codes. The
W
Fig. 2. FE models of the curved beam. (a) The formulation A. (b) The formulation ff.
V. KOMPIS
278
C
B
D (a)
Fig. 3. FE model of a hole in a band. (a) Three-element model. (b) Eight-element model, the formulation A. (c) Eight-element model. the formulation B.
elements can be connected to the classical displacement elements having the same nodal points on the adjacent interelement boundaries. The shape of the
element can be more general than that by the classical formulations. The transition between elements of different orders is very simply realized. There is a possibility to define special functions which enable us to model various local effects [4, $81. In this case, the special functions which satisfy the equilibrium equations will be included directly into the matrix [B] of the eqn (11). The B fo~ulation is especially attractive, because it introduces the elimination of all inner degrees of freedom on the element level and enables better satisfaction of boundary conditions and local boundary effects. The transition to the 3D problems and/or other
field problems is straightforward. Acknowledgement--The author appreciates the cooperation of the Slovenske Energeticke Strojame TImaCe for the loan of a PC computer for this research work.
REFERENCES I. E. Trefftz, Ein Gegenstueck zum Ritzschen Verfahren. Proc. 2nd Congr. Appf. Mech., pp. 131-13’7, Zurich (1926). 2. M. A. Jaswon, Integral equation method in potential theory. Proc. Roy. Sot. 275, 23-32 (1963). 3. W. Prager, Variational principles of linear elastostatics for discontinuous displacements, strains and stresses. In Recent Progress in Applied Mechanics, pp. 463474. The Folke-Odqvist Volume, Aimqvist and Wiksell, Stockholm f 1967). 4. J. Jirousek and N. Leon, A powerful finite element for plate bending. Comput. Meth. appl. Mech. Engng 12, 77-96 (1977). 5. J. Jirousek, Basis for development
of large finite elements locally satisfying all held equations. Comput. Meth. anal. Mech. Ennna 14. 65-92 (1978). 6. V. KompiS, Appro~i~a~on of displ~~m~nts in FEM satisfying all field equations inside the elements. Proc. 1st ConJ on Mech., Vol. 2, pp. 88-91, Prague (1987). 7. S. Timoshenko and J. N. Goodier, Theory of Elasticity, 2nd Bdn. McGraw-Hill, New York (1951). 8. J. Jirousek and A. Venkatesh, Hybrid-Trefftz plane elasticity elements with p-method capabilities. Internal report LSC 91/7, Swiss Federal Inst. of Technology, Lausanne (1991).
Copyright @J 1994 Ltd Printed in Gnat B&in. All ri&ts wcrwd a0457949p4 sd.00 + 0.00
Pergamon
FINITE ELEMENTS SATISFYING ALL GOVERNING EQUATIONS INSIDE THE ELEMENT V. KOMPIS Institute of Mechanics. Technical Unjversjty of Transport and ~ommunj~tions,
&%a,
Slovak Republic (Received 11 August 1992) introduce a formulation of the finite element method for the solution of elasticity problems in which all the governing equations will be satisfied inside the element. The displacements inside the elements are approximated by polynomials, whose coefficients are given from equilibrium equations. Another independent displacement field, which is common to adjacent elements, is defined on the element boundaries only. The element stiffness equations are derived using two modified variational principles. The element can have arbitrary shape and can be connected to isoparametric displacement elements and/or can be used as a transition element. Although the solution is shown for two-dimensional elasticity, the extension on three-dimensional and/or other field problems is straightforward. AlMract-We
1. INTRODUCTION
differentiation in corresponding direction and the Einstein’s summation is used for repeated indices. Assuming a homogeneous, isotropic material, the constitutive relationships are
Whereas in the standard isoparametric finite element method (FEM) the approximate solution is designed to satisfy interelement compatibility of displacements, attempts have been made to derive the formulations satisfying all governing equations inside the elements. Trefftz presented first a variational method[l] in which the approximate solution was designed to satisfy a priori the differential equations of the problem. From these principles the boundary element method [2] was formulated. For the FEM incorporation of the Trefftz principles new variational principles for interelement connection derived from more general principles [3] were employed [4,5]. The aim of the present work is to present a new method of presentation of displacement fields satisfying all governing equations inside the elements in polynomial form [6]. Two variational formulations for interelement connections are used and corresponding stiffness matrices for two-dimensional (2D) elastic solution are presented. The extension on three-dimensional (3D) and/or other field problems is straightfo~ard.
CFIl= A$, + f&z 022= BE,, + A~22 612- a2, = 2Gt,*,
(2)
where A = E/(1 - v*) B=vA
G=fE/(l
+v)
(3)
and E is Young’s modulus and v is Poisson’s ratio. The strain tensor tij may be expressed in terms of the displacement vector u, Cii= $(U$. j + uj, J.
2. GOVERNING EQUATIONS
The internal equilibrium is then expressed (from (2x4)) by the Lam&-Navier equations
The governing equations for 2D elasticity probequations
lems are the equiiibrium
Oii,$+ & = 0
(i*j = 1,2),
Au,, II +
(11
where (T,~is the stress tensor in Cartesian coordinates xi, 6, is the body force vector acting on a unit volume. The index after a comma denotes the
(4)
Fu2.12
+
G4.22
=
F~,.~+&~z+%I,=
-4
-b2,
(3
where F=BfG=fE/(l 273
-v).
(6)
274
KOMPIS
V. 1
Table
{Ai) = {a)‘), a!*‘, ui3’, . . . }‘,
(10)
1
*;
where T denotes transposed vector/matrix. Carrying out the second derivatives of the polynomial approximation of the displacements (8) and setting it into the Lam&-Navier equations (5) and comparing the coefficients of equal polynomial terms, some coefficients of {Ai} can be expressed by means of the others. We can see, that in this way all coefficients corresponding to inner terms in the Pascal triangle in Table 1 (e.g. the terms with the mixed product of coordinate components) can be expressed by means of coefficients of the outer (single component) terms. We can see also that some coefficients of the u, component of the displacement field are expressed by coefficients of its second component u2 and vice versa. Derivation of the dependent coefficients in each row in Table 1 from the independent ones could seem to decrease the efficiency of the procedure. However, we need to do this only once for each material, even for different degrees of interpolation polynomials, because this dependence is given by material constants only. Moreover an efficient procedure can be found if we realize the form of the matrices (some submatrices are diagonal the others are three-diagonal). If we chose independent coefficients of the third order to be the terms xix2 and x,x: and these of the fourth order x:x2 and x,x: (the solution for these terms was in this particular case simpler), we are able to write the polynomial satisfying all governing equations (body forces are not included) in the form
XI x2 _ XIX2 *f
x: XI-4 x:x* x; x: x,x: x:x: 4x2 d 3. FINITE ELEMENT SOLUTION
3.1. Displacement fields satisfying the governing equations inside the elements Let us suppose the interpolation form
polynomial in the
P(x) = [a”)]o + [amx, + a(3)x2],
+ [a’“‘x: + a”‘x, x2 + a(6)xi]2 + [a(‘)x: + a(*)xix2 + at9)x, xi + a(‘O)x&
+ [fJ”“xf + u”“x~x2 + a(‘3)xix: + ~2(‘~)x, xi + a(15)xj]4+
.. ,
(7)
where the lower index after square brackets denotes the order of corresponding polynomial terms. Thus the full polynomial of n th order will include the terms OfO-, l-,..., and n th order. In the matrix notation the ith component of the displacement field will be approximated as
4 = [Pl{4),
(8)
where [p] is a row matrix of polynomial terms [PI =[[110, [x,rx*l,r [x:,xIx2rx:L..
. .I
(9)
{a 1 = ]Bl{c, 19
and {A,} is a column matrix of its coefficients
where the coefficients of the matrix [B] = [B, , B,] are Table 2
By’ = x, B\4’
=
B’P’=
x;
By = x* /q5’
-x:A/(3G)
X:X,
By’ =x3*I 2 @y’ = /y/I) = qw
=
*;
B\” = x,x;
-x;G/(3,4)
By’ = x x3 = 0
B\“‘=
:2:,x,G,F
B\14’= -2x,x,A/F
B\15’= -x;
B\‘@= -xi F/(3B)
B’j”=(H,x:+x;A/~)/4-
B\‘*’ = (x:G/F
- l.Sx;x;A/F
+ H,x;)/4
B$” = Bi2’ = B$j’ = 0
F/(3A)
Bi6’ = -x;F/(3G)
B$“= -x;F/(3A)
Bf’=(x:H,+x;G/F)/4-
B$‘j’
zz x;
Bi15’
=
X:X,
8s”’ = x:.x2
+x:A/F)/4B$“’ =x, -
x;G/(3A)
l.Sx;x;G/F
Bi4’ = -2x,x,A/F
Bi5’= -2x,x,G/F BI”=(x:H, B$‘O’= 1
(11)
l.Sx:x;G/F B”2’ *
_ -x2
B\‘4’
=
x;
By@ = x,x;
- x;A/(3B)
B$‘*’ = x, x;,
where H, = (G2- F’)/(AF)
H, = (AZ - F2)/(GF)
l.Sxfx;A/F
FE satisfying all goveming equations inside the element given in Table 2 and {u) = {u,, u2) is the vector of displacements. The vector {c, 1 is a vector of unknown coefficients for an element. The number of unknown coefficients in this vector is equal to 2(1 + 2n), where n is the order of interpolation polynomial Ip]. The components of the deformation tensor will be given by
will be defined on all element boundaries, variation
I
where
{E)= (q,, 622,62fT [cl
=I&,,
(4,,+B2.1)/21~ (14)
4.2,
and substituting this into (2) we obtain the components of the stress tensor fcf = (G,,, oz2, ctzf by
_ ri
by PI =
tr; 1 -v2
--L 0
1
0
0
1
0
I
(1 -v)/2
StT(v, - ~7)dT = 0
(19)
will enforce the continuity of displacements over the whole element boundary r, (it must be satisfied for each element separately). Moreover the variation &*{t,-I)dT I r,
=0
(20)
>
will enforce the continuity of tractions. r, and T, are the parts of the boundary where displacements and tractions are given respectively. r, is the interelement boundary. So the interelement equilib~~ (continuity of tractions) between the adjacent elements is given by the first integral (20) and thus this variation also establishes that the equilbrium on the whole will be satisfied. If the independent displacement field fii is approximated by means of shape functions NW and displacements U@ in nodal points, the functions &will be C, continuous on continuous parts of the element boundaries, i.e.
(16)
and from it also the tractions fi = aijnj,
and the
f
Wt c dr + or in the vector form
275
(17)
where nj is a component of a unit normal vector of the surface. In vectorial form the tractions are {t) = {f,* t,JT
for each boundary segment of integration. In this notation each boundary segment of an element has the same meaning as an element by the boundary element method (BEM). W is the vector of displacements in thejth nodal point of the element boundary. Defining now
3.2. The element formulation The displacement field (11) will be generally incompatible between the elements and the tractions which will be defined from these fields. Variational principles for general discontinuous fields were given by Prayer [3]. If we define a new boundary displacement field I$, which will be completely independent of ui (and consequently of the related boundary vector vi), two successful approaches can be found by variational formulations defined by Jirousek [S]. The two independent fields are imposed to fuifil the compatibility of incompatible fields in an integral sense. By the first approach (A) the independent field ci
From eqn (19) we get the dependence between the
V. KOMPIS
276
coefficients of interpolating polynomials {c,} and of the nodal degrees of freedom {U?} in the form =
(26)
W,l{c,~ = [G,l{ve)
s r,
and eqn (20) gives the relation
+
KKl%J = Fe 1
qce 1’
(27) =
and substituting
--B{cyJT [B]qdT
s ru
(bIDI KITfi dr
g{cyJT_(pc*)
(35)
{c,} from (26) into (27) we get
(28)
KlVW’~~A~~,) = {pr)
GdTt,dT = ,(,}’ i
I r,
r,
which can be written in the form
[NIT([nl[Dl[CI)dT(c,J
=8{~cWc*F{ce).
(36)
We obtain from (31) where
]Hfl{c,j
&I= [~elr[fW’[~,l
(30)
is the stiffness matrix of the element. The second approach (B) will be based on the definition of the independent displacement field 6, on the interelement boundaries only. The variation -
=
[Woe,) + ICI
and from (32) [GflT{c,} = 0. Substituting
[GflTW:I~‘[Gtl~~~~= -[G:lTWbl-‘{P:
6tr(o,-6)dT
1 (39)
I r. +
(38)
(37) into (38) we get
St:@, - a) dI-
6o,T(t, -7) dT +
(37)
=0
(31)
Kl{ V, )= If?)? will enforce the continuity of boundary conditions (if the boundaries T, and TU are parts of the boundary of the element f,) and the continuity of displacements on interelement boundaries. Further the variation
(40)
where
[&I= [Gfl~fW’[G:I
(41)
is the element stiffness matrix and Srt, dT = 0 will enforce the continuity Let us denote now 6t:v,
s r.+r’, = S{cy}’
of tractions.
GvJt,
dT -
dl-
s r,
s
r +r (blPl[WPld~~c,~ Y f
(32)
{P:} = -[G:lrIHfl~‘{P:}
(42)
is the vector of static equivalent loads. We can see, that both A and B formulations give the same element matrices for elements not containing boundaries I’, and r, The second formulation, B, does not contain nodal points on the rU and T, boundaries and so the matrices are built similarly as classical elements with reduced degrees of freedom on these boundaries. The number of nodes in an element cannot be greater than 2n + 1, where n is the order of the interpolation polynomial (according to Sec. 3.1). The element can be connected to compatible displacement elements which have the same shape functions along the adjacent boundaries. The shape of the element can be arbitrary and also the degree of shape functions along each boundary segment can be arbitrarily chosen. So it can be used as a transition element for the different meshes.
FE satisfying all governing equations inside the element
277
points. We can also compare the normal stresses on the loaded edge with the prescribed boundary tractions. The calculated values (related to the boundary tractions) varied from 0.75 to 1.08, while by the displacement FE fo~ulation from 0.24 to 1.89. A curved beam by the pure bending (Fig. 2)
Fig. I. FE model of the cantilever.
Formulation for 3D or other types of fields is straightforward.
The exact solution [7] containing non-polynomial terms is 10.364 on the inner radius and -9.658 on the outer radius. The B formulation gave stresses varying from 10.312 to 10.420 and from -9.612 to -9.768, respectively. A hole in a band (Fig. 3)
4. EXAMPLES
The following three examples will show some properties of the described elements. A five-point Gauss formula was used for integration along the element boundaries and the stresses were calculated in the Gauss nodal points everywhere. A cantilever loaded by a uniform pressure (Fig. 1)
If the boundary tractions were given to get the solution described by the polynomial of the fourth order, the fo~ulation B gives the exact solution for any shape of interelement boundaries. If the clamped end boundary conditions are given by a, = 0 and tr is quadratic as before, the exact solution is not known and it is not described by the polynomial of the fourth order. The maximal stress on the element boundary (in the closest Gauss nodal point of the clamped boundary) was 0.9468 of that obtained for the former boundary conditions. If the boundary was fully clamped (u, = 0, ur = 0), the stress on the boundary with unifo~ pressure was 0.9704 (related to that of the former solution) and 0.9186 on the opposite edge. Supposing the linear distribution of normal stress on the clamped boundary, we get the stresses 0.9288 in corresponding
(a)
Oniy one quarter of the region was modelled by FE, because of the double symmetry of the problem. The convergency of the fo~ulation was examined on this example. The theoretical stresses [7] are 0.75 and 4.30 in the points A and B, respectively (for the traction on the edge E-F equal to 1.0). Models consisting of two, three and eight elements are shown in Fig. 3 (the two-element model is the same as the three-node model, only the region A-B-D-C is modelled by one element). The stresses calculated in the point nearest to B were: 2.24, 3.65 and 4.38, respectively, for those models. The displacements u, were 0.00986, 0.00888 and 0.00876 in the point D and 0.00868, 0.00808 and 0.00811 in the point C for corresponding models. These results were obtained for the formulation A. A similar behaviour was observed for the formulation B.
5. CONCLUSIONS
The given FE formulation satisfying all the governing equations inside the element is very attractive, it can be incorporated into the existing FE codes. The
W
Fig. 2. FE models of the curved beam. (a) The formulation A. (b) The formulation ff.
V. KOMPIS
278
C
B
D (a)
Fig. 3. FE model of a hole in a band. (a) Three-element model. (b) Eight-element model, the formulation A. (c) Eight-element model. the formulation B.
elements can be connected to the classical displacement elements having the same nodal points on the adjacent interelement boundaries. The shape of the
element can be more general than that by the classical formulations. The transition between elements of different orders is very simply realized. There is a possibility to define special functions which enable us to model various local effects [4, $81. In this case, the special functions which satisfy the equilibrium equations will be included directly into the matrix [B] of the eqn (11). The B fo~ulation is especially attractive, because it introduces the elimination of all inner degrees of freedom on the element level and enables better satisfaction of boundary conditions and local boundary effects. The transition to the 3D problems and/or other
field problems is straightforward. Acknowledgement--The author appreciates the cooperation of the Slovenske Energeticke Strojame TImaCe for the loan of a PC computer for this research work.
REFERENCES I. E. Trefftz, Ein Gegenstueck zum Ritzschen Verfahren. Proc. 2nd Congr. Appf. Mech., pp. 131-13’7, Zurich (1926). 2. M. A. Jaswon, Integral equation method in potential theory. Proc. Roy. Sot. 275, 23-32 (1963). 3. W. Prager, Variational principles of linear elastostatics for discontinuous displacements, strains and stresses. In Recent Progress in Applied Mechanics, pp. 463474. The Folke-Odqvist Volume, Aimqvist and Wiksell, Stockholm f 1967). 4. J. Jirousek and N. Leon, A powerful finite element for plate bending. Comput. Meth. appl. Mech. Engng 12, 77-96 (1977). 5. J. Jirousek, Basis for development
of large finite elements locally satisfying all held equations. Comput. Meth. anal. Mech. Ennna 14. 65-92 (1978). 6. V. KompiS, Appro~i~a~on of displ~~m~nts in FEM satisfying all field equations inside the elements. Proc. 1st ConJ on Mech., Vol. 2, pp. 88-91, Prague (1987). 7. S. Timoshenko and J. N. Goodier, Theory of Elasticity, 2nd Bdn. McGraw-Hill, New York (1951). 8. J. Jirousek and A. Venkatesh, Hybrid-Trefftz plane elasticity elements with p-method capabilities. Internal report LSC 91/7, Swiss Federal Inst. of Technology, Lausanne (1991).