- Email: [email protected]
Basis for isoparametric stress elements
COMPUTER METHODS IN APPLIED MECHANICS 0 NORTH-HOLLAND PUBLISHING COMPANY
AND ENGINEERING
BASIS FOR ISOPARAMETRIC
2 (1973)
43-63
STRESS ELEMENTS
John ROBINSON Consultant, Vicarage Road, Verwood, Dorset, England Consultant to Atkins Research and Development, EPSOM, Surrey, England
Revised
Received 8 December 197 1 manuscript received 2 August
197 2
The family of the so-called ‘isoparametric strain (displacement) elements’ is restricted to membranes and solids. The reason for this restriction has led to the development of a new family based on stress assumptions; these elements will be referred to as ‘isoparametric stress elements’. This family contains plates and solids but no membranes. The omission of a particular element in each family is consistent with the plate-membrane.analogies. The basic flexibility matrix of an isoparametric stress element is singular since the zero stress state is directly included. The rank technique is adopted to automatically extract the zero stress modes such that the element can be completely interchangeable between any finite element system. The theory for stress assumed isoparametric “quadrilateral” plate bending elements with curved boundaries is given. A brief presentation of the theory for isoparametric stress solid elements is also included.
0. Introduction Some of the most interesting developments in recent years in the field of finite element methods have been the family of ‘Isoparametric Elements’, introduced by Zienkiewicz and Company [ 1, 23, and the family of ‘Lumina Elements’, introduced by Argyris [3,41, based on a more general concept. These elements are based on strain (displacement) assumptions and will now be referred to more generally as ‘isoparametric strain elements’. It is interesting to note that this family contains only membranes and solids. The reason for this lies in the original definition of an isoparametric element, that is, the displacement interpolation functions are assumed to be the same as the shape interpolation functions. Therefore, for the existence of an isoparametric strain element the number of displacement coordinates must equal the number of shape or modal coordinates. For membrane elements there are two shape or modal coordinates (x,~) and two displacement coordinates (u, IJ); for solid elements there are three shape coordinates (x, y, z) and three displacement coordinates (u, u, w). However, for plate bending elements there are two shape coordinates (x,~) but only one displacement coordinate (w), hence, there are no isoparametric strain plate bending elements. In this paper a new family of isoparametric elements is presented in which the elements are based on stress assumptions; these will be referred to as ‘isoparametric stress elements’. In this family the stress interpolation functions (Airy-type stress functions) are assumed to be the same as the shape or modal interpolation functions. For plate bending elements there are two shape coordinates (x, y) and two stress coordinates (@i, (a,); for solid elements there are three shape coordinates (x,y, z) and three stress coordinates ((a,, @*, aa). However, for membrane elements
J. Robinson,
44
Basis for isoparametric stress elements
there are two shape coordinates (x,~) but only one stress coordinate (a,,), hence, there are no isoparametric stress membrane elements. These new developments lead to a more general definition of an isoparametric element; an isoparametric element is one whose elastic (stress or strain) interpolation functions are assumed to be the same as the shape interpolation functions. Such an element exists when the number of elastic (stress or strain) and shape coordinates are equal. The definition and existence of isoparametric forms are shown clearly in table 1. Table 1 Isoparametric elements and their existence. Analogies
Isoparametric form for stress elements Stress elements
Strain elements
Membrane elements
Membrane elements
aZQ1
au
OX = y--
@l
____
EX = TX
u,u
9
X>Y
azQl O =’ a,*
au EY= ay
X,1’
No
a*Ql TXY= - axay
au au YXY=G+G
Yes
Plate bending elements Ql,
Plate bending elements a*, 1 _=_ Rx ax*
@2
W
i _a*~
X,Y
R~
X>Y
ay*
No
Yes
Solid elements
Solid elements a2Q3 OX=-+ay* Ql>Qz>Q3
X,Y>Z
0
Y
a*Q* a,*
a*Q,
a*Q3
az*
ax*
=-+-
a2Q2 oz=-+ax*
au eY = G
a*Q,
u,u,w
aw
X,YJ
ay2
01
or x, y,z
a2Q3 7xY = - axay
Yes
a*Q, Tyz = - aya;
Tzx = ~-
a*Q* azax
--I_
Isoparametric form for strain elements
au au %Xy= G + ay
x, y,z Yes
au + 2 YZX
=a,
ax
J. Robinson,
Basis for isoparametric stress elements
45
The isoparametric stress element theory leads to a singular flexibility matrix which corresponds to element stress degrees of freedom (nodal values of the assumed stress functions). In this form the element can only be used in a force method or combined method I (Lagrangian Multiplier). To enable complete interchangeability of the element between any finite element system it is necessary to extract the zero stress state to give a natural (non-singular) flexibility matrix corresponding to independent force variables. It is shown that this can be carried out automatically using the rank technique. This technique was originally developed for the automatic selection of redundancies in the finite element force method. The detailed theory for stress assumed isoparametric “quadrilateral” plate bending elements with curved boundaries is developed but has been restricted to elements with no distributed applied loading. The theory for isoparametric stress solid elements is briefly discussed. This paper concentrates on elements whose shape and elastic (displacement or stress) functions are the same (isoparametric). However, a more general family of displacement elements (Hermes family) was developed by Argyris and Company 151 in which the shape functions are the same as previously but the displacement functions are different. The latter functions contain displacement derivatives as degrees of freedom. In ref. [ 61 Sander introduced derivatives of stress functions for equilibrium elements. Hence, the Hermes concept can also be extended to stress elements. The work in this paper has also been influenced by the work of Clough [71, Morley [ 81, Timoshenko and Woinowsky-Krieger [ 91, Fraeijs de Veubeke [10,11 I.
1. Isoparametric stress plate bending elements The theory for isoparametric “quadrilateral” will now be developed in detail. 1.1. Shape functions The shape of an element coordinates, that is,
can be expressed
plate bending elements
in terms of interpolation
with curved boundaries
functions
and its nodal
N
x=
c Pj& 77>Xj
(1)
i=l
and
(2) where P,(t, Q)are interpolation functions in the curvilinear coordinates g and 77,xi and yi are the local x,; coordinates if node i, and N is the total number of nodes for the element. See fig. l(a). 1.2. Stress functions Two stress functions, that is,
a’1 and Cp,, are now assumed in the same form as the shape functions,
J. Robinson,
46
Basis for isoparametric stress elements
(3) and N
@P, = C Pitt, rl)+,j
(4)
i=l
where the @ri and Qzj terms are the values of the stress functions tively. The nodal values will be denoted by the vector,
and referred
to as the ‘element
@, and Q2 at node i, respec-
stress degrees of freedom’. Corner Nodes
a,o,@,@
Oxy
Local Axes
(a)
(b) Fig. 1. Plate element
with curved
boundaries.
1.3. Moment functions The case considered here will be for no distributed applied load (p,=O). The moments and M,, can be expressed in terms of the stress functions [ 8, 101 as,
M,, M,
(6)
(7)
J. Robinson,
Using these relationships librium equation,
the moment
a=M,,
a=M, __-_2_
axay
ax2
Basis for isoparametric stress elements
functions
will automatically
satisfy the incremental
47
equi-
a=M, _
+
--0.
(9)
ay2
The plate shears, Q, and Q,, are given by, Q,2!$-~
(10)
a,=!$-!!!$
Substituting that is,
from eqs. (6), (7) and (8) into (IO) gives the shears in terms of the stress functions,
Substituting
eqs. (3) and (4) into (6), (7) and (8) results in, -
_
aa
Mx
-
ay
M=’MY M XY 2
a@,
=
(12)
ax
1 -__ a% 5 ( ay
a@2 + ax
)_
where TMa, and
=
[T,(t, ~1,T2(t, 17),.... T&t, v), .... TN&
~11
(13)
48
J. Robinson,
Basis for isoparametric
stress elements
(14)
The sign convention
and shears is shown in fig. 2.
for the moments
aMyx Myx+ -ddy
Y
ay
t
My’ aM’dy
ay-
.
T ----------
CJ,+ sdy (I)
ay
dy
a
L_.
Qxl
x
Fig. 2. Incremental plate loading.
The interpolation functions Pi are defined in terms of the curvilinear coordinates t and 1); a change in derivatives is therefore required in eq. (14). The Cartesian and curvilinear derivatives are related by,
(15)
where J is the Jacobian
J=
matrix,
that is,
r-1
73x
ay
-ax a7)
-ay
agaE aa
(16)
J. Robinson, Basisfor isoparametricstress elements
Differentiating
49
eqs. ( 1) and (2) in accordance with ( 16) gives the Jacobian matrix as,
(17)
The derivatives of eq. (14) are now obtained using eqs. (15) and (17), that is, aPi =J-’
ar -aPi ’
_~ aq_.I 1.4. Stress and strain fields
The stress and strain fields are respectively, (19) and E= Ee-f-Eo
G?O)
uG= (0, cry Q-_} = stress field
(21)
E = i E, e,, 3;,,) = total strain field
(22)
where
(23) = stress related strains zO= stress-free initial strains
(24
50
J. Robinson,
Basis for isoparamettic stress elements
TI =
(25)
t = plate thickness
E = modulus
of elasticity
v = Poisson’s ratio. The sign convention
for stresses si shown in fig. 3.
Fig. 3. Incremental
Substituting
stresses.
from eq. (12) into (19>, (20, (23) gives the stress field UC7 = LA?lCl
(26)
T o@ =z
(27)
where
and the total strain field E=
T,,@ ma + E,,
(28)
J. Robinson,
Basis for isoparametric stress elements
51
where
(F )We, &&A
Tea= z
(29)
.
1.5. Distributed loads on curved boundaries The normal bending moment (M,), twisting moment (AI,,) and shear (Q,) distributions curved boundary ([ = f 1 or 77= f 1) are obtained using the relations [ 91 -
Mn=
Kt :
Q, -
where for a constant
sine
=
e
~0s~
sin8
cost9
c0se -
J&
c0se
sine -?ay
ML
cos2 e - sin2 8
4
-
(
c0se a+sinO ay
&
111
-
1 (30)
MXY _
-1 % 1ic32ayc 1 w
i M2ax
c0se =
-sine
- 2 sine cosf3
[ boundary,
+
and for a constant
sin’e
on a
(31)
(sign W)
aJ’*% Caq ) 1
E
q boundary,
+
ay
))
at
(-sign (7)))
2%
0
The boundary moments and shear have a right hand vector system (n, t, z), see fig. l(c). The angle is positive as shown in fig. I(b), again following a vector system.
(32)
52
J. Robinson,
Basis for isoparametric stress elements
1.6. Kirchhoff shear and corner forces The Kirchhoff shear (VW) distribution
a”nt
(
V, = Q, - a~
!J
on a boundary
(g = + 1 or TJ= +_1) is given by,
1
(33)
where S is measured along a boundary and positive as shown in fig. l(a). The directional tives a/as and a/an are obtained using the relations ‘a
-!
as
cosf3
=
a
an
Substituting
1 -aax
- sin8
deriva-
-aax
1_I =c
cos0
a
sin8 1
.ay-
(34)
-a
ay
from eq. (15) into (34) yields
-a- as
a
=cj-’
an -1 a
The Kirchhoff
H 2
(35)
a
G
concentrated
corner
force at node j ( Wj> is given by
(36) where M,, is the twisting moment distribution on boundary ij and M,, is the twisting moment distribution on boundary jk, both of which are evaluated at ,$ = tj and 7) = r)j. Noting the arrangement of subscripts the remaining corner forces are given by
w, = Wjkl> The boundary
w,=‘,,i>
wi = w,ij .
loadings and corner forces have a right-hand
(37) vector sign convention
(n, t, z).
1.7. Element singular jlexibility matrix and initial deformations The force-deformation relationship for the element, in terms of the stress degrees of freedom {Gm,}, is derived by applying the principle of virtual forces. For the jth value, (38) where * moi = jth stress degree of freedom
J. Robinson, (~~1,.
=
correspondingjth
Basis for isoparametric stress elements
stress field
“WW = deformation corresponding V = volume of element indicates virtual quantities. Substituting
53
to @‘moi
from eqs. (26) and (28) into (38) gives
(39) Therefore,
for all virtual values,
Gmoa,,,o = G,,,,, s T&U&@,,
+Q
dl’
(40)
V
where & m. is a diagonal matrix. Integrating over the volume and solving for the deformations
a
mo
1
t2
=
s
ss
-t/2
--1
The force-deformation a
ma=f
ma
gives
1
T,&,U’eQ@m,
+ co> Det
U)
4
&
dz
.
(41)
-1
relationship
in terms of stress degrees of freedom
is therefore
+%nao
@ ma
(42)
where f,,
= element t/2
singular flexibility 1
matrix based on stress assumptions
1
T,f, T,, Det (J) dt dq dz = -tL
and
ama0
-1J
= vector of corresponding
-1
(44)
-1
from eqs. (27) and (29) into (43) yields
L7=#
1 ’
Substituting
1
T&,EODetCJ)dijdqdz.
-t/2
Substituting
initial deformation
sss 1
t2
=
(43)
-1J
-1
j TAG TEaTMa Det W) dt -1
from eq. (27) into (44) gives
drl
(45)
J. Robinson,
54
amoo = _tY/z
(j$
Basis for isoparametric stress elements
TAQ TIE0 Det (J) d&idq dz
The flexibility matrix f,, is singular since the element independent set of variables.
(46)
stress degrees of freedom,
Qmo, are not an
1.8. Element boundary force systems The normal bending moment (M,) and Kirchhoff shear (V,) distributions on the element boundaries, together with the Kirchhoff corner forces, can be expressed in terms of the stress degrees of freedom using the relationships given by eqs. (30) to (37) and the moment functions of eq. ( 12). The generalised forces on the ith boundary are denoted by Hi=
[
IcIrzi ‘ni
1
(47)
which in turn can be expressed
in terms of the stress degrees of freedom
as
Hi = THQi(Dmo
(48)
For all boundaries H,
(49)
= TH@@mo
where I
H, = W,
I -_- i Hi i ___}
(50)
and r-
T l____.__.l H01
The Kirchhoff
corner forces can also be expressed
*rn = ?v,@,, Equations
(51)
in terms of the stress degrees of freedom
(52)
.
(49) and (52) can now be combined
as
together
to give
NnI = ~mNA?Io
(53)
N,
(54)
where = {II,,, i W,}
J. Robinson,
Basis for isoparametric stress elements
55
and E mN@ =
(55)
Matrix E mNQgives the element boundary force systems for unit values of the stress degrees of freedom, each system being in equilibrium. 1.9. Element stress output The stress field within an isoparametric plate element is given by eq. (26) that is
=a= T&4,,
.
The moment functions are given by eq. (12), that is M = T&Dmo
.
The stress output for a plate element would be, say, moments and stresses at selected points within the element. This can be expressed in terms of the stress degrees of freedom as u ma
(56)
= SmoQ(Dmcr
where u mo
=
vector of required output data for an isoparametric stress element
Sma@
=
stress output matrix which is derived using matrices T,, specified locations within the element.
and TMQ)at
1.10. Characteristic element matrices The force-deformation relation (eq. (42)), boundary force systems (eq. (53)), and stress output (eq. (56)) is basically all the information that is required to use the isoparametric stress element in a force method or combined method I. However, to obtain complete interchangeability of any type of element between the various analysis systems [ 121 (displacement, force or combined) requires, for elastostatic problems, four characteristic matrices, namely; 1. The element natural elastic matrix (stiffness or flexibility) 2. The initial stress-free deformation vector 3. The element assembly matrix 4. The element output matrix. These matrices define the elastic behaviour, the spatial assembly into the structure, and the required output information for each element, they have been described and demonstrated in considerable detail by Robinson and Haggenmacher [ 121. 1.10.1. Independent element stress degrees of freedom The force-deformation relationship of eq. (42) includes the zero stress state, in other words, the stress degrees of freedom {am,} are not an independent set of variables, this leads to a singular
J. Robinson,
56
Basis for isoparametric stress elements
flexibility matrix. To use an isoparametric stress element in a displacement method or combined method II programme requires the extraction of the zero stress state to give the element natural flexibilitv matrix and corresponding initial deformation vector [ 121. This can be carried out automatically by using “The Rank Technique” [ 13, 141 which was originally developed for the automatic selection of redundancies in the force method. Consider now eq. (42), that is u M0 =fmA?3a+ Q?n,o ’ This equation a
can be written
MO-
as
a ma0 = am,e
=f*9m7
(57)
corresponding where a,,_ are the elastic deformations In the augmented form, eq. (57) becomes
Applying form
the rank technique
[:;:_+?]
to Oma.
to eq. (58) will result, after rearranging
{am0
(y,,,)
= [_j
The lower part of eq. (59) is now extracted
equations,
.
(59)
to give
The number of equations in this system corresponds to the number of dependent of freedom. Application of the rank technique to eq. (60) will result in a rncie
=
in the equivalent
stress degrees
(61)
Tado4?me
selected independent where d,,,,, are the automatically sides of eq. (61) by T& and solve for d,,,. Hence d moe
=
Tdare
= (T:doT&)T:do
elastic deformations.
Premultiply
both
(62)
Tdaa%
where
Substituting
(63)
.
from eq. (57) into (62) gives
d mae
= Tdwr(ama
-a,,O)
.
(64)
J. Robinson,
Basis for isoparametric stress elements
57
Equation (64) expresses the independent elastic deformations d,,, in terms of the original set of and the corresponding initial deformations, amoO. This non-independent total deformations, urn,, equation can also be written as
d moe
=
(65)
dmu- dmoo
where
d rnlJ = ~dcioQ?no
(66)
= total independent
deformations
and dma0
= TdaoamoO (67) =
Rearranging
initial stress-free
independent
deformations.
eq. (65) gives dmu
= dmue
+ dmoO
(68)
’
The relationship between the non-independent stress degrees of freedom a,, and the indepencan now be established by applying the principle of virtual dent stress degrees of freedom F,, displacements for virtual a,, . For the jth value (69) Now, from eq. (66), (70) Substituting
eq. (70) into (69) gives (71)
Therefore,
for all virtual values (72)
where Grn, is a diagonal matrix.
a mu = T@FOFmO
Hence (73)
where (74)
58
J. Robinson,
Basis for isoparametric stress elements
1.10.2. First and second characteristic matrices The elements being considered are based on stress assumptions, in this case, the elastic matrix which is derived first is the element natural flexibility matrix D,,. The force-deformation relationship for the element is now obtained in two steps. Firstly, substitute eq. (42) into (66), noting eq. (74), this gives
Secondly,
substitute d ma =
eq. (73) into (75) to give
GFavmoGFoFma +amacJ .
Hence, the force-deformation F,, is d
m(~
=
DmoFma
relationship
+
(76)
for the element
in terms of independent
force variables
dmoO
(77)
where D mo
= TiFofmoT@Fo
= element
natural flexibility
(78)
matrix
and d mu0
= GT%,O
(79) = initial stress-free
(independent)
deformation
vector.
1.10.3. Third characteristic matrix The element force assembly matrix for a stress element can be given in various forms. This matrix transforms the independent force variables F,, into chosen external element forces (boundary distributions and/or nodal forces). Substituting eq. (73) into (53) gives
(80)
Nm= EmNFFmo where E mNF
= EmN@
T@Fo
.
(81)
In eq. (80), matrix Em, is the element force assembly matrix in the local system, this matrix transforms the independent force variables into boundary force distributions and corner forces. It should be again noted that {IV,} contains the static supports. This form of assembly matrix is typical of equilibrium elements [ 1 11. In the traditional finite element methods and programmes the assembly matrix for a stress element would transform the independent force variables into equivalent nodal forces resolved into the global system, that is
J. Robinson,
Qma
=E
mQF
Basis for isoparametric stress elements
59
F ma
(82)
is the vector of element nodal forces in the global system (including static supports) and EmQF is the element force assembly matrix in the global system. Equation (82) is obtained by applying the principle of virtual displacements to find the equivalent nodal forces which are in equilibrium with the boundary forces as given by eq. (80). This procedure was described in detail in ref. [ 121. Various types of assembly matrices have just been discussed; in this light a node must now be considered not only as a discrete point but also as an interface. The element force assembly matrix E,, (third characteristic matrix) can now take either of the discussed forms depending on the method of forming the nodal equilibrium equations, that is, nodal in the most general sense. where Q,,
1.10.4. Fourth characteristic matrix The stress field within an isoparametric plate element can be expressed dent force variables by substituting eq. (73) into (26), that is,
in terms of the indepen-
=a = =oFFmo
(83)
T oF
(84)
where
The moment
= =o@=@Fa
functions
are obtained
by substituting
eq. (73) into (12), hence
M =
‘MFFmu
(85)
TMF
= =M@=@Fo
(86)
where
The stress output for the element is established by mutual agreement between the developer and user, this may be stresses and moments at specific locations within the element. The stress output is expressed as 0
ma = ‘maF
(87)
Fmo
where 0
mo
= vector of required
output
data for a stress element
S maE’ -- stress element
output matrix,for the force method and combined method I. This matrix is derived using matrices ToF and TMF at specified locations within the element.
1.1 1. Interchangeability of elements To achieve complete interchangeability of the characteristic matrices of any type of element between the various systems for finite element analysis it is essential to work with independent
J. Robinson,
60
stress or strain variables. Therefore, the characteristic matrices are;
Basis for isoparametric stress elements
for the interchangeability
of an isoparametric
(1) D,,
- -- Natural
flexibility
matrix (eqs. (78) and (45))
(2)
k?o
-- - Initial deformation
vector (eqs. (79) and (46))
(3)
E,,
- - - Assembly
(4)
%lcrF - -- Output
matrix (eq. (8 1) or (82) depending
on the definition
stress element
of a node)
matrix (eq. (87)).
These matrices are those required for a force method or combined method I (fig. 4). To use an isoparametric stress element in a displacement method or combined method II (fig. 4) the required characteristic matrices are = 0;:
---------,
(1)
k,,
(2)
Fmao = D;‘,d,,,
(3)
u,,
= EL, ---------,
(4)
Loll
= LoF D;‘, - - - -, Output
These interchangeability ity is shown in fig. 4.
Natural
stiffness matrix
(88)
- - - - -, Initial force vector Assembly
relations
(89)
matrix
matrix
(90)
.
(91)
were developed
in ref. [ 121 and the process of interchangeabil-
lsoparametric Stress
Element
I Dmti
Natural
d rnEO
Initial
E rnG
Assembly
s
Output
rnGF
Flexibility
Matrix
Deformation
Vector
Matrix Matrix
1 k m.5 = 02, F moo = %h alno 5
7,
dm,,
= EA,
moD = S,&c
I
Force
Combined
Method
Method
Fig. 4. Interchangeability
Combined
I
Method
of isoparametric
Displacement II
Method
stress element.
J. Robinson,
Basis for isoparametric stress elements
61
2. Isoparame tric stress solid elements In this paper the theoretical developments for stress assumed isoparametric solid elements will be restricted to the essential steps which define an isoparametric element, namely, the shape functions, the stress functions, and the expressions for global stresses in terms of stress functions.
2.1. Shape functions The shape of a solid element can be expressed coordinates of the nodes, that is,
in terms of interpolation
functions
and the global
N x = c
@L
(92)
77, rwj
i=l
N
Y=Z pjcb?,w,
(93)
and N
z =2
‘i(t,
(94)
1), !t>‘i
where P,(& r], 5) are interpolation functions in the curvilinear coordinates g, 7) and 5, Xi, Yi and Zi are the global coordinates of node i, and N is the total number of nodes for the element. See fig. 5.
rr’
Y
Zp--------_ Fig. 5. Isoparametric
2.2. Stress functions Three stress functions, tions, that is,
X solid element.
Cpr, a, and ap,, are now assumed in the same form as the shape func-
J. Robinson,
62
Basis for isoparametric stress elements
N @p, =
’i=1pj(t, rl,
(95)
3‘)@h,j
N (96)
and
N (97)
where the @,j, Qzj and Qp,j terms are the values of the stress functions Qp,, Cp, and a13 at node i, respectively. The nodal values, referred to as element stress degrees of freedom, are denoted by (98) 2.3. Stress field The global stress field
=a= ((5xx can be expressed
UYY
0z.z TXY TYZ
in terms of the stress functions
a*@ a*@ +__-A *XZ-2 aY* az* a*+, a;cp, (JY= --- + az* dX2 a*+ a% (JZ=__-_2+ _~’ ax* ar* a**3
rxy
= -
rzx1
yi>y,
[ 151 as
(99)
(100)
(101)
a*@, ryz
= - rjaz
,
rzx
a*+*
= - __ azax
.
(102)
References L’l 1. Ergatoudis,
B.M. Irons and O.C. Zienkiewicz, “Curved, isoparametric, ‘quadrilateral’ elements for finite element analysis”, Int. J. Solids and Structures 31 (1968) 31-42. Symp. on I21 I. Ergatoudis, B.M. Irons and O.C. Zienkiewicz, “Three dimensional analysis of arch dams and their foundations”, arch dams, Institution of Civil Engineers (1968) analysis of two arch dams by a finite element method, Part II”, Symp. [31 J.H. Argyris and S.C. Redshaw, “Thre-_ dimensional on Arch Dams, Institution of Civil Engineers (1968) method”, Aeronatical Journal, Tech. Note 11, 72 (1968). 141 J.H. Argyris, “The LUMINA element for the matrix displacement
J. Robinson,
Basis for isoparametric stress elements
63
[S] J.H. Argyris, I. Fried and D.W. Scharpf, “The Hermes 8 element for the matrix displacement method”, Aeronautical Journal, Tech. Note 12, 72 (1968). [6] G. Sander, “Application of the dual analysis principle”, Proc. IUTAM Colloquium on High Speed Computing of Elastic Structures, Liege, Belgium (1970). [ 71 R.W. Clough, “Comparison of three dimensional finite elements”, Proc. Symp. on Application of Finite Element Methods in Civil Engineering, School of Engineering, Vanderbilt University, Nashville, Tennessee (1969) pp. l-26. (81 L.S.D. Morley, “A triangular equilibrium element with linearly varying bending moments for plate bending problems”, J. Roy. Aero. Sot. 71 (1967) 715. [9] S.P. Timoshenko and S. Woinowsky-Krieger, Theory of plates and shells (McGraw-Hill Book Company, Inc., New York, 1959). [lo] B.M. Fraeijs de Veubeke, “Basis of a well conditioned force program for equilibrium models via the Southwell slab analogies”, AFFDL-TR-67-10 (1967). [ 1 l] B.M. Fraeijs de Veubeke and G. Sander, “An equilibrium model for plate bending”, Int. J. Solids and Structures 4 (1968) 447-468. (121 J. Robinson and G.W. Haggenmacher, “Basis for element interchangeability in finite element programs”, Proceedings of the Third Conference on Matrix Methods in Structural Mechanics, Wright-Patterson Air Force Base, Ohio (1971). [ 131 J. Robinson, Structural matrix analysis for the engineer (John Wiley and Sons, Inc., New York, 1966). [ 141 J. Robinson and G.W. Haggenmacher, “Optimization of redundancy selection in the finite element force method”, AIAA Journal 8, No. 8 (1970). [15] E.E. Sechler, Elasticity in engineering (John Wiley and Sons, Inc., New York, 1952).
AND ENGINEERING
BASIS FOR ISOPARAMETRIC
2 (1973)
43-63
STRESS ELEMENTS
John ROBINSON Consultant, Vicarage Road, Verwood, Dorset, England Consultant to Atkins Research and Development, EPSOM, Surrey, England
Revised
Received 8 December 197 1 manuscript received 2 August
197 2
The family of the so-called ‘isoparametric strain (displacement) elements’ is restricted to membranes and solids. The reason for this restriction has led to the development of a new family based on stress assumptions; these elements will be referred to as ‘isoparametric stress elements’. This family contains plates and solids but no membranes. The omission of a particular element in each family is consistent with the plate-membrane.analogies. The basic flexibility matrix of an isoparametric stress element is singular since the zero stress state is directly included. The rank technique is adopted to automatically extract the zero stress modes such that the element can be completely interchangeable between any finite element system. The theory for stress assumed isoparametric “quadrilateral” plate bending elements with curved boundaries is given. A brief presentation of the theory for isoparametric stress solid elements is also included.
0. Introduction Some of the most interesting developments in recent years in the field of finite element methods have been the family of ‘Isoparametric Elements’, introduced by Zienkiewicz and Company [ 1, 23, and the family of ‘Lumina Elements’, introduced by Argyris [3,41, based on a more general concept. These elements are based on strain (displacement) assumptions and will now be referred to more generally as ‘isoparametric strain elements’. It is interesting to note that this family contains only membranes and solids. The reason for this lies in the original definition of an isoparametric element, that is, the displacement interpolation functions are assumed to be the same as the shape interpolation functions. Therefore, for the existence of an isoparametric strain element the number of displacement coordinates must equal the number of shape or modal coordinates. For membrane elements there are two shape or modal coordinates (x,~) and two displacement coordinates (u, IJ); for solid elements there are three shape coordinates (x, y, z) and three displacement coordinates (u, u, w). However, for plate bending elements there are two shape coordinates (x,~) but only one displacement coordinate (w), hence, there are no isoparametric strain plate bending elements. In this paper a new family of isoparametric elements is presented in which the elements are based on stress assumptions; these will be referred to as ‘isoparametric stress elements’. In this family the stress interpolation functions (Airy-type stress functions) are assumed to be the same as the shape or modal interpolation functions. For plate bending elements there are two shape coordinates (x, y) and two stress coordinates (@i, (a,); for solid elements there are three shape coordinates (x,y, z) and three stress coordinates ((a,, @*, aa). However, for membrane elements
J. Robinson,
44
Basis for isoparametric stress elements
there are two shape coordinates (x,~) but only one stress coordinate (a,,), hence, there are no isoparametric stress membrane elements. These new developments lead to a more general definition of an isoparametric element; an isoparametric element is one whose elastic (stress or strain) interpolation functions are assumed to be the same as the shape interpolation functions. Such an element exists when the number of elastic (stress or strain) and shape coordinates are equal. The definition and existence of isoparametric forms are shown clearly in table 1. Table 1 Isoparametric elements and their existence. Analogies
Isoparametric form for stress elements Stress elements
Strain elements
Membrane elements
Membrane elements
aZQ1
au
OX = y--
@l
____
EX = TX
u,u
9
X>Y
azQl O =’ a,*
au EY= ay
X,1’
No
a*Ql TXY= - axay
au au YXY=G+G
Yes
Plate bending elements Ql,
Plate bending elements a*, 1 _=_ Rx ax*
@2
W
i _a*~
X,Y
R~
X>Y
ay*
No
Yes
Solid elements
Solid elements a2Q3 OX=-+ay* Ql>Qz>Q3
X,Y>Z
0
Y
a*Q* a,*
a*Q,
a*Q3
az*
ax*
=-+-
a2Q2 oz=-+ax*
au eY = G
a*Q,
u,u,w
aw
X,YJ
ay2
01
or x, y,z
a2Q3 7xY = - axay
Yes
a*Q, Tyz = - aya;
Tzx = ~-
a*Q* azax
--I_
Isoparametric form for strain elements
au au %Xy= G + ay
x, y,z Yes
au + 2 YZX
=a,
ax
J. Robinson,
Basis for isoparametric stress elements
45
The isoparametric stress element theory leads to a singular flexibility matrix which corresponds to element stress degrees of freedom (nodal values of the assumed stress functions). In this form the element can only be used in a force method or combined method I (Lagrangian Multiplier). To enable complete interchangeability of the element between any finite element system it is necessary to extract the zero stress state to give a natural (non-singular) flexibility matrix corresponding to independent force variables. It is shown that this can be carried out automatically using the rank technique. This technique was originally developed for the automatic selection of redundancies in the finite element force method. The detailed theory for stress assumed isoparametric “quadrilateral” plate bending elements with curved boundaries is developed but has been restricted to elements with no distributed applied loading. The theory for isoparametric stress solid elements is briefly discussed. This paper concentrates on elements whose shape and elastic (displacement or stress) functions are the same (isoparametric). However, a more general family of displacement elements (Hermes family) was developed by Argyris and Company 151 in which the shape functions are the same as previously but the displacement functions are different. The latter functions contain displacement derivatives as degrees of freedom. In ref. [ 61 Sander introduced derivatives of stress functions for equilibrium elements. Hence, the Hermes concept can also be extended to stress elements. The work in this paper has also been influenced by the work of Clough [71, Morley [ 81, Timoshenko and Woinowsky-Krieger [ 91, Fraeijs de Veubeke [10,11 I.
1. Isoparametric stress plate bending elements The theory for isoparametric “quadrilateral” will now be developed in detail. 1.1. Shape functions The shape of an element coordinates, that is,
can be expressed
plate bending elements
in terms of interpolation
with curved boundaries
functions
and its nodal
N
x=
c Pj& 77>Xj
(1)
i=l
and
(2) where P,(t, Q)are interpolation functions in the curvilinear coordinates g and 77,xi and yi are the local x,; coordinates if node i, and N is the total number of nodes for the element. See fig. l(a). 1.2. Stress functions Two stress functions, that is,
a’1 and Cp,, are now assumed in the same form as the shape functions,
J. Robinson,
46
Basis for isoparametric stress elements
(3) and N
@P, = C Pitt, rl)+,j
(4)
i=l
where the @ri and Qzj terms are the values of the stress functions tively. The nodal values will be denoted by the vector,
and referred
to as the ‘element
@, and Q2 at node i, respec-
stress degrees of freedom’. Corner Nodes
a,o,@,@
Oxy
Local Axes
(a)
(b) Fig. 1. Plate element
with curved
boundaries.
1.3. Moment functions The case considered here will be for no distributed applied load (p,=O). The moments and M,, can be expressed in terms of the stress functions [ 8, 101 as,
M,, M,
(6)
(7)
J. Robinson,
Using these relationships librium equation,
the moment
a=M,,
a=M, __-_2_
axay
ax2
Basis for isoparametric stress elements
functions
will automatically
satisfy the incremental
47
equi-
a=M, _
+
--0.
(9)
ay2
The plate shears, Q, and Q,, are given by, Q,2!$-~
(10)
a,=!$-!!!$
Substituting that is,
from eqs. (6), (7) and (8) into (IO) gives the shears in terms of the stress functions,
Substituting
eqs. (3) and (4) into (6), (7) and (8) results in, -
_
aa
Mx
-
ay
M=’MY M XY 2
a@,
=
(12)
ax
1 -__ a% 5 ( ay
a@2 + ax
)_
where TMa, and
=
[T,(t, ~1,T2(t, 17),.... T&t, v), .... TN&
~11
(13)
48
J. Robinson,
Basis for isoparametric
stress elements
(14)
The sign convention
and shears is shown in fig. 2.
for the moments
aMyx Myx+ -ddy
Y
ay
t
My’ aM’dy
ay-
.
T ----------
CJ,+ sdy (I)
ay
dy
a
L_.
Qxl
x
Fig. 2. Incremental plate loading.
The interpolation functions Pi are defined in terms of the curvilinear coordinates t and 1); a change in derivatives is therefore required in eq. (14). The Cartesian and curvilinear derivatives are related by,
(15)
where J is the Jacobian
J=
matrix,
that is,
r-1
73x
ay
-ax a7)
-ay
agaE aa
(16)
J. Robinson, Basisfor isoparametricstress elements
Differentiating
49
eqs. ( 1) and (2) in accordance with ( 16) gives the Jacobian matrix as,
(17)
The derivatives of eq. (14) are now obtained using eqs. (15) and (17), that is, aPi =J-’
ar -aPi ’
_~ aq_.I 1.4. Stress and strain fields
The stress and strain fields are respectively, (19) and E= Ee-f-Eo
G?O)
uG= (0, cry Q-_} = stress field
(21)
E = i E, e,, 3;,,) = total strain field
(22)
where
(23) = stress related strains zO= stress-free initial strains
(24
50
J. Robinson,
Basis for isoparamettic stress elements
TI =
(25)
t = plate thickness
E = modulus
of elasticity
v = Poisson’s ratio. The sign convention
for stresses si shown in fig. 3.
Fig. 3. Incremental
Substituting
stresses.
from eq. (12) into (19>, (20, (23) gives the stress field UC7 = LA?lCl
(26)
T o@ =z
(27)
where
and the total strain field E=
T,,@ ma + E,,
(28)
J. Robinson,
Basis for isoparametric stress elements
51
where
(F )We, &&A
Tea= z
(29)
.
1.5. Distributed loads on curved boundaries The normal bending moment (M,), twisting moment (AI,,) and shear (Q,) distributions curved boundary ([ = f 1 or 77= f 1) are obtained using the relations [ 91 -
Mn=
Kt :
Q, -
where for a constant
sine
=
e
~0s~
sin8
cost9
c0se -
J&
c0se
sine -?ay
ML
cos2 e - sin2 8
4
-
(
c0se a+sinO ay
&
111
-
1 (30)
MXY _
-1 % 1ic32ayc 1 w
i M2ax
c0se =
-sine
- 2 sine cosf3
[ boundary,
+
and for a constant
sin’e
on a
(31)
(sign W)
aJ’*% Caq ) 1
E
q boundary,
+
ay
))
at
(-sign (7)))
2%
0
The boundary moments and shear have a right hand vector system (n, t, z), see fig. l(c). The angle is positive as shown in fig. I(b), again following a vector system.
(32)
52
J. Robinson,
Basis for isoparametric stress elements
1.6. Kirchhoff shear and corner forces The Kirchhoff shear (VW) distribution
a”nt
(
V, = Q, - a~
!J
on a boundary
(g = + 1 or TJ= +_1) is given by,
1
(33)
where S is measured along a boundary and positive as shown in fig. l(a). The directional tives a/as and a/an are obtained using the relations ‘a
-!
as
cosf3
=
a
an
Substituting
1 -aax
- sin8
deriva-
-aax
1_I =c
cos0
a
sin8 1
.ay-
(34)
-a
ay
from eq. (15) into (34) yields
-a- as
a
=cj-’
an -1 a
The Kirchhoff
H 2
(35)
a
G
concentrated
corner
force at node j ( Wj> is given by
(36) where M,, is the twisting moment distribution on boundary ij and M,, is the twisting moment distribution on boundary jk, both of which are evaluated at ,$ = tj and 7) = r)j. Noting the arrangement of subscripts the remaining corner forces are given by
w, = Wjkl> The boundary
w,=‘,,i>
wi = w,ij .
loadings and corner forces have a right-hand
(37) vector sign convention
(n, t, z).
1.7. Element singular jlexibility matrix and initial deformations The force-deformation relationship for the element, in terms of the stress degrees of freedom {Gm,}, is derived by applying the principle of virtual forces. For the jth value, (38) where * moi = jth stress degree of freedom
J. Robinson, (~~1,.
=
correspondingjth
Basis for isoparametric stress elements
stress field
“WW = deformation corresponding V = volume of element indicates virtual quantities. Substituting
53
to @‘moi
from eqs. (26) and (28) into (38) gives
(39) Therefore,
for all virtual values,
Gmoa,,,o = G,,,,, s T&U&@,,
+Q
dl’
(40)
V
where & m. is a diagonal matrix. Integrating over the volume and solving for the deformations
a
mo
1
t2
=
s
ss
-t/2
--1
The force-deformation a
ma=f
ma
gives
1
T,&,U’eQ@m,
+ co> Det
U)
4
&
dz
.
(41)
-1
relationship
in terms of stress degrees of freedom
is therefore
+%nao
@ ma
(42)
where f,,
= element t/2
singular flexibility 1
matrix based on stress assumptions
1
T,f, T,, Det (J) dt dq dz = -tL
and
ama0
-1J
= vector of corresponding
-1
(44)
-1
from eqs. (27) and (29) into (43) yields
L7=#
1 ’
Substituting
1
T&,EODetCJ)dijdqdz.
-t/2
Substituting
initial deformation
sss 1
t2
=
(43)
-1J
-1
j TAG TEaTMa Det W) dt -1
from eq. (27) into (44) gives
drl
(45)
J. Robinson,
54
amoo = _tY/z
(j$
Basis for isoparametric stress elements
TAQ TIE0 Det (J) d&idq dz
The flexibility matrix f,, is singular since the element independent set of variables.
(46)
stress degrees of freedom,
Qmo, are not an
1.8. Element boundary force systems The normal bending moment (M,) and Kirchhoff shear (V,) distributions on the element boundaries, together with the Kirchhoff corner forces, can be expressed in terms of the stress degrees of freedom using the relationships given by eqs. (30) to (37) and the moment functions of eq. ( 12). The generalised forces on the ith boundary are denoted by Hi=
[
IcIrzi ‘ni
1
(47)
which in turn can be expressed
in terms of the stress degrees of freedom
as
Hi = THQi(Dmo
(48)
For all boundaries H,
(49)
= TH@@mo
where I
H, = W,
I -_- i Hi i ___}
(50)
and r-
T l____.__.l H01
The Kirchhoff
corner forces can also be expressed
*rn = ?v,@,, Equations
(51)
in terms of the stress degrees of freedom
(52)
.
(49) and (52) can now be combined
as
together
to give
NnI = ~mNA?Io
(53)
N,
(54)
where = {II,,, i W,}
J. Robinson,
Basis for isoparametric stress elements
55
and E mN@ =
(55)
Matrix E mNQgives the element boundary force systems for unit values of the stress degrees of freedom, each system being in equilibrium. 1.9. Element stress output The stress field within an isoparametric plate element is given by eq. (26) that is
=a= T&4,,
.
The moment functions are given by eq. (12), that is M = T&Dmo
.
The stress output for a plate element would be, say, moments and stresses at selected points within the element. This can be expressed in terms of the stress degrees of freedom as u ma
(56)
= SmoQ(Dmcr
where u mo
=
vector of required output data for an isoparametric stress element
Sma@
=
stress output matrix which is derived using matrices T,, specified locations within the element.
and TMQ)at
1.10. Characteristic element matrices The force-deformation relation (eq. (42)), boundary force systems (eq. (53)), and stress output (eq. (56)) is basically all the information that is required to use the isoparametric stress element in a force method or combined method I. However, to obtain complete interchangeability of any type of element between the various analysis systems [ 121 (displacement, force or combined) requires, for elastostatic problems, four characteristic matrices, namely; 1. The element natural elastic matrix (stiffness or flexibility) 2. The initial stress-free deformation vector 3. The element assembly matrix 4. The element output matrix. These matrices define the elastic behaviour, the spatial assembly into the structure, and the required output information for each element, they have been described and demonstrated in considerable detail by Robinson and Haggenmacher [ 121. 1.10.1. Independent element stress degrees of freedom The force-deformation relationship of eq. (42) includes the zero stress state, in other words, the stress degrees of freedom {am,} are not an independent set of variables, this leads to a singular
J. Robinson,
56
Basis for isoparametric stress elements
flexibility matrix. To use an isoparametric stress element in a displacement method or combined method II programme requires the extraction of the zero stress state to give the element natural flexibilitv matrix and corresponding initial deformation vector [ 121. This can be carried out automatically by using “The Rank Technique” [ 13, 141 which was originally developed for the automatic selection of redundancies in the force method. Consider now eq. (42), that is u M0 =fmA?3a+ Q?n,o ’ This equation a
can be written
MO-
as
a ma0 = am,e
=f*9m7
(57)
corresponding where a,,_ are the elastic deformations In the augmented form, eq. (57) becomes
Applying form
the rank technique
[:;:_+?]
to Oma.
to eq. (58) will result, after rearranging
{am0
(y,,,)
= [_j
The lower part of eq. (59) is now extracted
equations,
.
(59)
to give
The number of equations in this system corresponds to the number of dependent of freedom. Application of the rank technique to eq. (60) will result in a rncie
=
in the equivalent
stress degrees
(61)
Tado4?me
selected independent where d,,,,, are the automatically sides of eq. (61) by T& and solve for d,,,. Hence d moe
=
Tdare
= (T:doT&)T:do
elastic deformations.
Premultiply
both
(62)
Tdaa%
where
Substituting
(63)
.
from eq. (57) into (62) gives
d mae
= Tdwr(ama
-a,,O)
.
(64)
J. Robinson,
Basis for isoparametric stress elements
57
Equation (64) expresses the independent elastic deformations d,,, in terms of the original set of and the corresponding initial deformations, amoO. This non-independent total deformations, urn,, equation can also be written as
d moe
=
(65)
dmu- dmoo
where
d rnlJ = ~dcioQ?no
(66)
= total independent
deformations
and dma0
= TdaoamoO (67) =
Rearranging
initial stress-free
independent
deformations.
eq. (65) gives dmu
= dmue
+ dmoO
(68)
’
The relationship between the non-independent stress degrees of freedom a,, and the indepencan now be established by applying the principle of virtual dent stress degrees of freedom F,, displacements for virtual a,, . For the jth value (69) Now, from eq. (66), (70) Substituting
eq. (70) into (69) gives (71)
Therefore,
for all virtual values (72)
where Grn, is a diagonal matrix.
a mu = T@FOFmO
Hence (73)
where (74)
58
J. Robinson,
Basis for isoparametric stress elements
1.10.2. First and second characteristic matrices The elements being considered are based on stress assumptions, in this case, the elastic matrix which is derived first is the element natural flexibility matrix D,,. The force-deformation relationship for the element is now obtained in two steps. Firstly, substitute eq. (42) into (66), noting eq. (74), this gives
Secondly,
substitute d ma =
eq. (73) into (75) to give
GFavmoGFoFma +amacJ .
Hence, the force-deformation F,, is d
m(~
=
DmoFma
relationship
+
(76)
for the element
in terms of independent
force variables
dmoO
(77)
where D mo
= TiFofmoT@Fo
= element
natural flexibility
(78)
matrix
and d mu0
= GT%,O
(79) = initial stress-free
(independent)
deformation
vector.
1.10.3. Third characteristic matrix The element force assembly matrix for a stress element can be given in various forms. This matrix transforms the independent force variables F,, into chosen external element forces (boundary distributions and/or nodal forces). Substituting eq. (73) into (53) gives
(80)
Nm= EmNFFmo where E mNF
= EmN@
T@Fo
.
(81)
In eq. (80), matrix Em, is the element force assembly matrix in the local system, this matrix transforms the independent force variables into boundary force distributions and corner forces. It should be again noted that {IV,} contains the static supports. This form of assembly matrix is typical of equilibrium elements [ 1 11. In the traditional finite element methods and programmes the assembly matrix for a stress element would transform the independent force variables into equivalent nodal forces resolved into the global system, that is
J. Robinson,
Qma
=E
mQF
Basis for isoparametric stress elements
59
F ma
(82)
is the vector of element nodal forces in the global system (including static supports) and EmQF is the element force assembly matrix in the global system. Equation (82) is obtained by applying the principle of virtual displacements to find the equivalent nodal forces which are in equilibrium with the boundary forces as given by eq. (80). This procedure was described in detail in ref. [ 121. Various types of assembly matrices have just been discussed; in this light a node must now be considered not only as a discrete point but also as an interface. The element force assembly matrix E,, (third characteristic matrix) can now take either of the discussed forms depending on the method of forming the nodal equilibrium equations, that is, nodal in the most general sense. where Q,,
1.10.4. Fourth characteristic matrix The stress field within an isoparametric plate element can be expressed dent force variables by substituting eq. (73) into (26), that is,
in terms of the indepen-
=a = =oFFmo
(83)
T oF
(84)
where
The moment
= =o@=@Fa
functions
are obtained
by substituting
eq. (73) into (12), hence
M =
‘MFFmu
(85)
TMF
= =M@=@Fo
(86)
where
The stress output for the element is established by mutual agreement between the developer and user, this may be stresses and moments at specific locations within the element. The stress output is expressed as 0
ma = ‘maF
(87)
Fmo
where 0
mo
= vector of required
output
data for a stress element
S maE’ -- stress element
output matrix,for the force method and combined method I. This matrix is derived using matrices ToF and TMF at specified locations within the element.
1.1 1. Interchangeability of elements To achieve complete interchangeability of the characteristic matrices of any type of element between the various systems for finite element analysis it is essential to work with independent
J. Robinson,
60
stress or strain variables. Therefore, the characteristic matrices are;
Basis for isoparametric stress elements
for the interchangeability
of an isoparametric
(1) D,,
- -- Natural
flexibility
matrix (eqs. (78) and (45))
(2)
k?o
-- - Initial deformation
vector (eqs. (79) and (46))
(3)
E,,
- - - Assembly
(4)
%lcrF - -- Output
matrix (eq. (8 1) or (82) depending
on the definition
stress element
of a node)
matrix (eq. (87)).
These matrices are those required for a force method or combined method I (fig. 4). To use an isoparametric stress element in a displacement method or combined method II (fig. 4) the required characteristic matrices are = 0;:
---------,
(1)
k,,
(2)
Fmao = D;‘,d,,,
(3)
u,,
= EL, ---------,
(4)
Loll
= LoF D;‘, - - - -, Output
These interchangeability ity is shown in fig. 4.
Natural
stiffness matrix
(88)
- - - - -, Initial force vector Assembly
relations
(89)
matrix
matrix
(90)
.
(91)
were developed
in ref. [ 121 and the process of interchangeabil-
lsoparametric Stress
Element
I Dmti
Natural
d rnEO
Initial
E rnG
Assembly
s
Output
rnGF
Flexibility
Matrix
Deformation
Vector
Matrix Matrix
1 k m.5 = 02, F moo = %h alno 5
7,
dm,,
= EA,
moD = S,&c
I
Force
Combined
Method
Method
Fig. 4. Interchangeability
Combined
I
Method
of isoparametric
Displacement II
Method
stress element.
J. Robinson,
Basis for isoparametric stress elements
61
2. Isoparame tric stress solid elements In this paper the theoretical developments for stress assumed isoparametric solid elements will be restricted to the essential steps which define an isoparametric element, namely, the shape functions, the stress functions, and the expressions for global stresses in terms of stress functions.
2.1. Shape functions The shape of a solid element can be expressed coordinates of the nodes, that is,
in terms of interpolation
functions
and the global
N x = c
@L
(92)
77, rwj
i=l
N
Y=Z pjcb?,w,
(93)
and N
z =2
‘i(t,
(94)
1), !t>‘i
where P,(& r], 5) are interpolation functions in the curvilinear coordinates g, 7) and 5, Xi, Yi and Zi are the global coordinates of node i, and N is the total number of nodes for the element. See fig. 5.
rr’
Y
Zp--------_ Fig. 5. Isoparametric
2.2. Stress functions Three stress functions, tions, that is,
X solid element.
Cpr, a, and ap,, are now assumed in the same form as the shape func-
J. Robinson,
62
Basis for isoparametric stress elements
N @p, =
’i=1pj(t, rl,
(95)
3‘)@h,j
N (96)
and
N (97)
where the @,j, Qzj and Qp,j terms are the values of the stress functions Qp,, Cp, and a13 at node i, respectively. The nodal values, referred to as element stress degrees of freedom, are denoted by (98) 2.3. Stress field The global stress field
=a= ((5xx can be expressed
UYY
0z.z TXY TYZ
in terms of the stress functions
a*@ a*@ +__-A *XZ-2 aY* az* a*+, a;cp, (JY= --- + az* dX2 a*+ a% (JZ=__-_2+ _~’ ax* ar* a**3
rxy
= -
rzx1
yi>y,
[ 151 as
(99)
(100)
(101)
a*@, ryz
= - rjaz
,
rzx
a*+*
= - __ azax
.
(102)
References L’l 1. Ergatoudis,
B.M. Irons and O.C. Zienkiewicz, “Curved, isoparametric, ‘quadrilateral’ elements for finite element analysis”, Int. J. Solids and Structures 31 (1968) 31-42. Symp. on I21 I. Ergatoudis, B.M. Irons and O.C. Zienkiewicz, “Three dimensional analysis of arch dams and their foundations”, arch dams, Institution of Civil Engineers (1968) analysis of two arch dams by a finite element method, Part II”, Symp. [31 J.H. Argyris and S.C. Redshaw, “Thre-_ dimensional on Arch Dams, Institution of Civil Engineers (1968) method”, Aeronatical Journal, Tech. Note 11, 72 (1968). 141 J.H. Argyris, “The LUMINA element for the matrix displacement
J. Robinson,
Basis for isoparametric stress elements
63
[S] J.H. Argyris, I. Fried and D.W. Scharpf, “The Hermes 8 element for the matrix displacement method”, Aeronautical Journal, Tech. Note 12, 72 (1968). [6] G. Sander, “Application of the dual analysis principle”, Proc. IUTAM Colloquium on High Speed Computing of Elastic Structures, Liege, Belgium (1970). [ 71 R.W. Clough, “Comparison of three dimensional finite elements”, Proc. Symp. on Application of Finite Element Methods in Civil Engineering, School of Engineering, Vanderbilt University, Nashville, Tennessee (1969) pp. l-26. (81 L.S.D. Morley, “A triangular equilibrium element with linearly varying bending moments for plate bending problems”, J. Roy. Aero. Sot. 71 (1967) 715. [9] S.P. Timoshenko and S. Woinowsky-Krieger, Theory of plates and shells (McGraw-Hill Book Company, Inc., New York, 1959). [lo] B.M. Fraeijs de Veubeke, “Basis of a well conditioned force program for equilibrium models via the Southwell slab analogies”, AFFDL-TR-67-10 (1967). [ 1 l] B.M. Fraeijs de Veubeke and G. Sander, “An equilibrium model for plate bending”, Int. J. Solids and Structures 4 (1968) 447-468. (121 J. Robinson and G.W. Haggenmacher, “Basis for element interchangeability in finite element programs”, Proceedings of the Third Conference on Matrix Methods in Structural Mechanics, Wright-Patterson Air Force Base, Ohio (1971). [ 131 J. Robinson, Structural matrix analysis for the engineer (John Wiley and Sons, Inc., New York, 1966). [ 141 J. Robinson and G.W. Haggenmacher, “Optimization of redundancy selection in the finite element force method”, AIAA Journal 8, No. 8 (1970). [15] E.E. Sechler, Elasticity in engineering (John Wiley and Sons, Inc., New York, 1952).