Finite element formulations for large strain, large displacement problems are formulated using a kinematic description based on the corotational components of the velocity strain. The corotational components are defined in terms of a system that rotates with each element and approximates the rotation of the material. To account for rotations of the material relative to this element system, extra terms are introduced in the velocity strain equations. Although this formulation is incremental, in explicitly integrated transient problems it compares very well with formulations that are not

Belytschko, T.

Stolarski, H.

Tsay, C. S.

English

NASA Technical Reports Server

COROTATIONAL VELOCITY STRAIN  FORMULATIONS FOR NONLINEAR 
ANALYSIS OF BEAMS AND AXISYMMETRIC  SHELLS* 
Ted Belytschko 
Department of Civil Engineering 
Northwestern University, Evanston, Illinois 
H. Stolarski  
I n s t i t u t e  of Fundamental Technological Research 
Polish Academy of Sciences, Warsaw, Poland 
C.S. Tsay 
Northwestern University 
SUMMARY 
Finite element formulations for large strain, large displacement problems 
are formulated using a kinematic description based on the corotat ional  compo- 
nents of the velocity strain.  The corotational components are defined in terms 
of a system that rotates with each element and approximates the rotat ion of  the 
material. To account for rotations of the  material r e l a t i v e  t o  t h i s  element 
system, extra terms are introduced in the velocity strain equations.  Although 
this formulation i s  incremental, in expl ic i t ly  in tegra ted  t rans ien t  problems it 
conpares very w e l l  with formulations that are not. Its s implici ty ,  and i t s  com- 
pa t ib i l i ty  wi th .  cons t i tu t ive  equat ions  found i n  "hydro" codes make it very 
a t t r a c t i v e  f o r  t h i s  class of problems. 
INTRODUCTION 
Nonlinear structures are conventionally treated by kinematic descriptions 
t h a t  are e s sen t i a l ly  Lagrangian in  na tu re ,  i n  t ha t  t he  measure of deformation 
is direct ly  re la ted to  the total  displacements .  Several  types of  Lagrangian 
formulations are frequently used: formulations based on the  Green s t r a i n  o r  
. Almansi s t r a i n  [1 ,2]  and formulations  based on coro ta t iona l  s t re tch  [3,4]. 
Although velocity strain formulations have been used extensively for  non- 
l inear  so l ids ,  as exemplified in the work of Key [ S I ,  little study has been made 
of the application of these formulations to structures.  Hughes and Liu  [61 
have presented a formulation based on the global components of the  ve loc i ty  
s t r a in .  
In this  paper ,  a corotat ional  veloci ty  s t ra in  formulat ion w i l l  be presented 
i n  which the  components of the  ve loc i ty  strain are expressed in a framework 
that  rotates  with the mater ia l ;  formulat ions of  this  type have been studied by 
Green and Naghdi [ 7 ] .  The formulation is then special ized to  f ini te  e lements  
*This work was supported by the   E lec t r ic  Power Research Ins t i t u t e .  
26 3 
https://ntrs.nasa.gov/search.jsp?R=19800024267 2020-03-21T15:57:12+00:00Z
by assuming that the rotation within an element i s  e i ther  cons tan t  or t h a t   t h e  
va r i a t ion  in  the  ro t a t ion  f i e ld  i s  small  or  moderate. 
The potent ia l  benefi ts  of  these methods are s igni f icant .  The basic equa- 
t ions  are simpler than Green s t r a i n   o r  Almansi equations, which endows the  re- 
su l t i ng  computer programs with both simplicity and speed. The stress conjugate 
to  the  coro ta t iona l  ve loc i ty  s t ra in  is the  Cauchy stress tensor expressed in 
the corotational system. Any constitutive equations based on  Cauchy stress and 
ve loc i ty  s t r a in  can therefore be used. Furthermore, the corotational stress 
and s t ress-s t ra in  matr ix  are both materially objective,  so no Jaumann type cor- 
rections need be made fo r  t he  stress s t a t e  and the formulation is  d i r ec t ly  
.appl icable  to  anisotropic  materials, which i s  not true of the formulations 
given in [5 ]  o r  [81 . 
In the next section, the fundamental equations for the corotational velo- 
c i ty  s t ra in  formulat ion are presented. Next, the general  equations for a f i -  
n i t e  element application of this formulation i s  given. In order to i l lustrate  
the  s implici ty  of  the method, we then give the formulation for a beam element 
assuming a constant rotation in the element.  More complex re la t ions  which ac- 
count for the variation of rotation in an element are then given. The last 
section gives some examples of the appl icat ion of  this  method to  nonl inear  
t rans ien t  problems. 
BASIC EQUATIONS 
we w i l l  use a kinematical and stress description by Green and Naghdi [7] . 
Let us denote the material coordinates of the s t ructure  by X . ,  t he  spa t i a l  co- 
ordinates by x i ,  t h e  displacements by u i  and the  ve loc i t i e s  hy vi. Then 
u =  x - xi i i (1) 
and the  deformation  gradient F i s  given by i j  
From the polar decomposition theorem (see [SI) it follows that the deformation 
gradient can be expressed as a pure deformation, which i s  expressed by a SF- 
metric  matrix U and a r i g i d  body rotat ion % in the  form kl' 
26 4 
The rotation matrix R i l  is  orthogonal, so that 
where 6 is  the  Kronecker de l ta .  l m  
We w i l l  denote the cyrdinate system which i s  rotated by the  r ig id  body 
motion of the material by x i  and call it a corotational coordinate system. 
This system is re la ted  t o  xi by 
x = R x  
A 
i j 
and i t s  or ientat ion var ies  from point  to  point  in  the material. 
The veloci ty  strain ( r a t e  of deformation) tensor i s  given by 
and the  coro ta t iona l  ve loc i ty  s t ra in ,  which i s  simply the same tensor with i t s  
components expressed in the corotational coordinates,  i s  given by 
A - 
Dkl - ?ik " j R  Dij (71 
The state of stress w i l l  be represented by the corotat ional  stress fjj 
which are the corotat ional  components of the Cauchy (physical) stress T 
the  t w o  are re la ted  by 
+j 
i j '  
A 
Tkl = R R i k  j R  Tij 
26 5 
The corotat ional  components of the stress are frame invariant ,  SO t h a t  the  
veloci ty  strain is re la ted  t o  t h e  rate of  corotat ional  stress by 
h A A 
'i j k l  Dkl 
where the matrix C fo r  a .material depends on the  state of stress and 
s ta te  var iab les  su&dk&s the  y ie ld  stress but is independent of material rota- 
tion, regardless of whether the material is  isotropic  or  anisotropic .  This  is  
a key advantage of corotational formulations.  If  the velocity strain and 
Cauchy stress are expressed in a fixed coordinate system, a Jaumann ra%e i s  
required  to  provide frame invariance, bu t  more importantly;  the  matrix c,i jklqmust 
a l so  be modified t o  account for the r i g i d  body rotation. Furthermore, such 
formulations are quite awkward in  s t ruc tura l  theor ies  where it i s  of ten conv- 
venient  to  dis t inguish veloci t ies  tangent  and normal to  the current  configu-  
ra t ion.  
A 
For a material i n  the domain n, t h e  rate of i n t e rna l  work i s  given by 
W = jQ Dij Tij dR = Dij Tij dR A A  
FINITE EIXMENT EQUATIONS 
We consider an'element which curre?tly occupies a volume R . Its nodal e 
displacements are utI, nodal velocities V i 1  and nodal forces fg". We represent 
the velocit ies withln the element by shape functions 
A A A  
V .  
1 
= N (x) viI 
_. 
I 
where NI are the shape functions which are expressed in ternls of the corota- 
t ional coordinates.  Throughout this paper,  upper case subscripts w i l l  refer t o  
nodal values, as exemplified in Eq. ( 1 2 )  , and the  ind ic ia1  summation convention 
w i l l  a l so  app ly  to  these subscripts.  
The pr inc ip le  of  v i r tua l  work gives 
266 
int A A  
viI fiI  = Je D~~ T~~ dn 
n 
For el.ements other than the simplest ,  i.e. those with l inear shape functions, 
the orientation of the corotational coordinates w i l l  vary within the element as 
shown i n  Fig. 1. Several  alternatives are then available for handling the right 
hand side  of Eq. (13) : 
i. a single corotational coordinate system G i  can be chosen for  the  ele- 
ment as shown i n  Fig. 1 and the  relative rotations ignored; 
ii. a single corotational coordinate system G can  be  chosen fo r  t he  ele- 
ment and the  ro t a t ions  r e l a t ive  to  xi can be accounted for by modify- 
ing the veloci ty  strain equations; 
i 
iii. the relat ive rotat ions can be accounted f o r  by using the transforma- 
t ions  (7)  and (9)  a t  each point of the element. 
For t h e  f i r s t  a l t e r n a t i v e ,  t h e  use of Eq. (12) gives 
V i 
so the use of the transformation (5) and the arbitrariness of of vi= gives 
It should be observed t h a t   t h e  stress is expressed in terms of a single coro- 
tational coordinate system throughout the element. Therefore, i f  w e  consider 
a b e a m  with a constant  axial  stress T, it follows that the only nonzero 
nodal forces l i e  along the Z axis regardless of the curvature of the beam. This 
anomaly can yield spurious results whenever the  f l exura l  s t i f fnes s  is s m a l l ,  
s ince it introduces parasi t ic  bending in  s ta tes  of  pure membrane stress, cf 
no ] .  
The second a l t e rna t ive  is t o  introduce velocity strain re la t ions  which 
account for the  var ia t ion in  rotat ions of  the element  but  to  express their compo- 
nents in the element system. If we represent these relations by 
26 7 
then Eq. (15) becomes 
i n t  - 
fiI = %i DmnkI mn T dR 
Re 
In a subsequent section, forms of Eq. (16) for  beams and axisymmetric she l l s  
w i l l  be presented. Higher order formulations as exemplified by Eq. (17) do 
provide better accuracy, particularly for relatively coarse meshes, but they do 
not eliminate parasit ic bending. 
A SIMPLE BEAM FORMULATION 
In order  to  i l lus t ra te  the  appl ica t ion  of a corotat ional  veloci ty  s t ra in ,  
we w i l l  f i r s t  consider a beam element with the simplest corotational formula- 
t ion where the nodal forces are evaluated by Eq.(15). The notation used i s  
shown i n  Fig. 2 .  We w i l l  embed the element corotational coordinate within t h e  
element so tha t  the z - axis always connects nodes 1 and 2.  Euler-Bernoulli 
beam theory w i l l  be used, so tha t  the  ve loc i t ies  through the depth are com- 
pletely defined by veloci t ies  of the  middle surface V; and V- Y 
so tha t  
where commas denote different ia t ion with respect  to  the subsequent variables. 
The veloci ty  f ie ld  7- w i l l  be approximated by l inear  shape functions, and 
the transverse velocity f l e l d  of the midline by cubic shape functions so  that  x 
268 
where t h e   r i g i d  body p a r t  of the  ve loc i ty  f ie ld  has  been omitted Since it 
causes no s t r a in .  The nodal velocities associated with deformation are thus 
and the conjugate nodal forces 
where m are  the  nodal moments.  Combining Eqs.(l8) and (19),  w e  obtain I - 
0. - 
DX 
= -  "x2 - 5 [(65-4)01 + (65-2)(b21 R 
Equation 115) then gives 
where A is the cross-sectional area of the beam. 
The remaining nodal forces can be obtained from equilibrium 
- - 
fxl - -?x2 
(24) 
269 
HIGHER ORDER VELOCITY 
STRAIN EXPRESSIONS 
The v e l o c i t y  s t r a i n  e x p r e s s i o n ,  E q . ( z 2 ) ,  is exac t  for  a beam only when t h e  
element 's  midline is  coinc ident  w i t h  the  x-axis ,  which corresponds t o  the chord 
between t h e  t w o  nodes. I f  the beam element is  i n i t i a l l y  c u r v e d  or i f  it de- 
forms enough so tha t  the  mid l ine  d i sp laces  subs t an t i a l ly  from the x-axis, these  
equat ions w i l l  be i n  error because of the following effects: 
1. The x and axes are no longer   coincident  w i t h  the  midline  and i t s  
normal , r e spec t ive ly ,  so Eqs . (18) are i n  error because vx and v are 
not  a long the midl ine and i t s  normal. Y 
2. The  volume in t eg ra t ion  neg lec t s  t h e  deformation of the beam relative 
t o  t h e  x ax i s .  
In  o rde r  t o  accoun t  fo r  t he  f i r s t  e f f ec t  w i thou t  t r ans fo rming  between co- 
ord ina te  sys tems wi th in  the  e lement ,  h igher  order  ve loc i ty  s t ra in  re la t ions  a re  
developed, using t h e  same b a s i c  ideas employed i n  [41 f o r  c o r o t a t i o n a l  stretch 
theories .  For  t h i s  purpose, the displacement of t h e  m i d l i n e  r e l a t i v e  t o  the 
chord, or x a x i s  , i s  denoted by Y. It w i l l  be assumed t h a t  
- = O ( e )  Y R U,- = o(e) X y,-- = o(e2) xx 
and accuracy of order  e2 i s  assumed; terms of h igher  order  and  cer ta in  o ther  
terms are neglected.  Because of space   l imi ta t ions ,  w e  omit t h e  d e r i v a t i o n ,  
and  g ive  on ly  the  f ina l  r e su l t  
K 1 + 
For an axisymmetric shell element, w e  l e t  x cor respond to  the  radial  coor- 
d ina t e ,  r ,  and €I be the  c i rcumferent ia l  coord ina te .  The corresponding relat ions 
are 
270 
A 
DX = d - Y K ~  
- 
1 
where is the  hoop v e l o c i t y   s t r a i n .  The terms dl  and K~ are i d e n t i c a l  t o  
those f o r   t h e  beam, E q s  . (26) , while  e 
vX 
d 2 = - - -  
- (V cos a + v s i n  a) 
r r X Y 
1 
Y I X  V - COS a +  - [ r ( c o s  a - Y I X  s i n  a) - Y s inacosdV - (28b) - 1 2 X , X  r K = -  Y ? X  
L 1 
where r is  the  cur ren t  rad ia l  coord ina te  of  a p o i n t ,  a the current angle between x and x ,  a s  shown i n  F i g .  2.  For a l l  app l i ca t ions  w e  have considered so f a r  , 
t he  second terms i n  t h e  e x p r e s s i o n s  f o r  K and tc2 have been insignif icant .  1 
EXAMPLES 
Resul ts  are g iven  for  a clamped r i n g  shown in  F ig .  3 ,  f o r  which experimen- 
t a l  r e s u l t s  are r e p o r t e d  i n  [ll]; numerical  resul ts  have been reported in  [3] .  
E x p l i c i t  time i n t e g r a t i o n  w i t h  t h e  c e n t r a l  d i f f e r e n c e  method and a lumped m a s s  
matr ix  w a s  used. The material model i s  e l a s t i c - p l a s t i c  w i t h  a Mises y i e l d  con- 
d i t i o n  and isotropic hardening. Because the width of the ring i s  l a rge  compared 
t o  i t s  th i ckness ,  p l ane  s t r a in  w a s  assumed in  the z-direct ion.  Furthermore,  the 
compressibi l i ty  of  the elastic s t r a i n s  w a s  cons idered  negl ib le ,  so t h e  h e i g h t  h 
of  the cross-sect ion w a s  modified by 
h = hoRo/R 
for  bo th  the  ve loc i ty  s t ra in  computa t ions ,  Eq . ( l8b) ,  and the  nodal  force  cornput- 
t a t ions ,  Eq . (23 ) .  The nodal  forces  w e r e  obtained by numerical  quadrature using 
f ive  poin ts  th rough the  th ickness ,  and two poin ts  a long  the  length  wi th  a trape- 
zo ida l  ru le .  
The displacements for the midpoint of t h e  clamped r i n g  are compared f o r  
p lane  stress and plane s t ra in ,  with and without  the area cor rec t ion  of  Eq . (29) ,  
i n  T a b l e  1. As can  be  seen ,  the  e f fec ts  of the  p lane  s t ra in  assumpt ion  are 
very significant,  changing the r e s u l t  by  20%. The e f f ec t  o f  t he  th i ckness  cor- 
r e c t i o n  is  less pronounced bu t  neve r the l e s s  no t  neg l ib l e .  
2 71 
TABLE 1 
EFFECTS OF ASSUMPTIONS ON MIDPOINT 
DEFLECTION  OF CLAMPED R I N G  
Assumptions Midpo in t   de f l ec t ion ( in )  
P l a n e  s t r a i n ,  v a r i a b l e  t h i c k n e s s  
P l a n e  s t r a i n ,  c o n s t a n t  t h i c k n e s s  
Plane stress, va r i ab le  th i ckness  
Plane stress, cons t an t  t h i ckness  
2 -99 
3.06 
3 -43 
3.53 
The p lane  s t r a in  so lu t ion  wi th  th i ckness  co r rec t ion  compares  best wi th  
the experiment ,  so w e  w i l l  restrict  a l l  subsequen t  compar i sons  to  th i s  case. 
The time h i s to ry  o f  t he  midpo in t  i s  compared t o  t h e  e x p e r i m e n t  i n  F i g ;  4.  The 
repor ted  exper imenta l  resu l t s  exhib i t  cons iderable  snapback ,  which  are absen t  
i n  the  computa t ions .  F igu res  5 and 6 compare t h e  deformed  shape  and s t r a i n  
t ime  h i s to r i e s  w i th  the  expe r imen t .  Overa l l ,  t he  ag reemen t  is  q u i t e  good  and 
comparable t o   t h a t  of t h e  n u m e r i c a l  r e s u l t s  r e p o r t e d  i n  [ 3 ] .  
272 
I 
REFE-NCES 
1. K . J .  Bathe  and H. Ozdemir, "Elastic-Plastic Large Deformation Static and 
Dynamic Analysis," Computers  and S t ruc tu res ,  V o l .  6, 1976, pp. 81-92. 
2. S. Yaghmai and E.P. Popov, "Incremental Analysis of Large Deflections of 
3. T. Belytschko  and B.J .  Hsieh, "Nonlinear Transient Finite Element Analysis 
with Convected Coordinates," International Journal for  Numerical Methods 
in  Engineering, Vol. 7,  1973,  pp. 255-271. 
4. T. Belytschko  and L.W. G l a u m ,  "Applications  of  Higher Order Corota t iona l  
S t re tch  Theor ies  t o  Nonlinear Finite Element Analysis," Computers and Struc- 
t u r e s ,  Vol. 10,  1979,  pp. 175-182. 
5. S.W. Key, "HONDO - A F i n i t e  Element Computer  Program for  the Large Defor- 
mation Dynamic Response of Axisymmetric Solids," Report  SLA-74-QO39, Sandia 
Laboratories,  April  1974. 
6. T . J . R .  Hughes,  and Y.K. Liu, "Nonlinear Finite Element Analysis of She l l s  - 
Par t  I ,  Three Dimensional Shells," to be published Computer Methods i n  
Applied Mechanics and Engineering. 
7. A.E. Green  and P.M. Naghdi, "A General Theory of an Elastic-Plastic Contin- 
uum," Archives of Rational  Mechanics  and Analysis, Vol. 18, 1965, pp. 
251-281. 
8. M.L. Wilkins,   "Calculation  of  Elastic-Plastic  Flow,"  Report  UCRL-7322- 
Revision 1, Lawrence Radiat ion Laboratory,  Universi ty  of  Cal i fornia ,  
Livermore,  1969. 
9.  L. Malvern, Introduction t o  t h e  Mechanics  of a Continuous Medium, Prent ice-  
H a l l  Incorporated,  Englewood C l i f f s ,  New Je r sey ,  1969. 
10. T. Belytschko, "Nonlinear Analysis - Descript ions and Numerical  Stabi l i ty ,"  
i n  Shock and Vibration Computer Programs, e d i t e d  by W. Pilkey and B. P i lkey ,  
The Shock and Vibration Information Center, 1975, pp. 537-562. 
11. H.A. Balmer and E.A. Witmer, "Theoretical-Experimental Correlation of Large 
Dynamic and Permanent Deformation of Impulsively Loaded Simple Structures ,"  
Report FDL-TDR-64-108, A i r  Force Fl ight  Dyn. L a b . ,  Wright-Patterson A i r  
Force Base, Ohio,  1964. 
27 3 
Y - 
/ 
DEFORMED 
UNDEFORMED 
I X 
Fig. 1. Deformed  and undeformed beam element showing 
element corotational coorgiGates x,? and local  
corotational coordinates x,y. 
Y I  
Fig. 2. Nomenclature 
for  nodal forces 
and veloci t ies .  
I node 1 x ,  ( r )  
27 4 
P’ 
E =  
Y ’  
2.61 X I bs-secz/ in4 
1.04 X lo7 psi J!n (!l/l0) 
0.3 
EXPERIMENT MODEL 
Fig. 3. Experiment [ l l ]  and model f o r  clamped r ing  problem. 
- 1  X EXPERIMENT 
4 t  - COMPUTATION 
0 
TIME (MILLISECONDS) 
Fig. 4. Comparison of computed centerpoint  (node 1 i n  
Fig.3) ver t ical  displacement  with experimental  
results i l l ] .  
275 
- 3 j : : : : : : : : : : : I  
0 i 2 3 4 5 6 
X - A X I S  
Fig. 5 .  Deformed configuration of clamped ring 
compared t o  experimental r e s u l t s  [ l l ] .  
7, I 
I I 
0 1 2 3 0 2 3 
TIME (MILLISECONDS) TIME (MILLISECONDS) 
Fig. 6. Comparison of computed strains against experi- 
mental resu l t s  [ll]; angle i n  parenthesis i s  for 
the experimental record. 
276 


Corotational velocity strain formulations for nonlinear analysis of beams and axisymmetric shells

Abstract

Similar works

Full text

Available Versions

NASA Technical Reports Server