Ghostscript wrapper for D:\Digitalizacja\MTS83_t21z1_4_PDF_artyku³y\mts83_t21z1.pdf M ECH AN I KA TEORETYCZNA I S TOS OWANA 1, 21 (1983) EF F ECTIVE  N ON LIN EAR AN ALYSIS  OF   ARBITRARY  TH IN   SH ELLS  BY  SI M P LE F IN ITE ELEM EN TS M I C H AŁ  K L E I B E R AN D R Z E J  Z A C H A R S K I IPPT   PAN W arszawa 1. Introduction Available  analytical  solutions  to  structural  problems  of  general  shells  are  limited  in scope and in general do not apply  to arbitrary  shapes, load conditions, irregular  stiffening and  support  conditions, cut  outs, and  many  ether  aspects  of  practical design.  In the case of  a  simultaneous  action  of  inelastic  material  properties  and  large  deformations  in  free- form  shells  the situation is  particularly  difficult  and successful  analytical  approaches'even for  simple  geometries  and  loadings  can  hardly  be  expected.  The finite  element method has  consequently  come to the fore  as  an approach to structural analysis  of  shells  because of  its  facility  to  deal  with  these complications. A  comprehensive review of  the history  of the finite  element developments for  shell  analysis  can be found in  [1]. This history evolved from  very  simple,  flat  elements  to  extremely  complicated  double  curved  elements  and covered  essentially  the  whole  range  of  different  approximations  applied  to  the  classical and  nonclassical  shell  equations  considered  at the  element level. U nfortunately  it is still fair  to say that the final  finite  element formulation for nonlinear thin  shell  analysis  is  far  from  being  settled.  F rom  the engineering  point of  view the pre- vailing problem is the lack  of elements which can show both good accuracy and  efficiency. The  apparantly  letter  deeply  curved  elements  quickly  because  so  complicated  that  its- attendant  high  computational cost  have  pretty well precluded the acceptance and general use of these elements. Their use seems  to  be particularly  ruled out in the case of  extensive nonlinearities  of  the  problem  which  requires  the  element  striffnesses  to  be  recalculated a  large  number  of  times.  It  is  therefore  understandable  that  the  researches  in  the  field of  nonlinear shell  analysis  turned back  to more simple elements. In this way,  as concluded in  [2], the history  of the development of  general  shell finite  elements has come  full  circle.. The  simplest  possible  geometrical  representation  of  a  doubly  curved  shell  surface  is- a  facet  approximation  by  flat  elements.  The  extremely  simple  and  efficient  formulation that.can be achieved by flat  elements make them very well suited for nonlinear applications, in particular when used in the framework  of the updated Lagrangian description of motion. The wide  range  of  numerical examples  studied in  [3 - 11]  indicate that flat  finite  elements may  be  very  useful  in  the  analysis  of  nonlinear  shell  problems.  However,  this approach has also  its well- known  deficiences  to mention only 26  M.  KLEIBER,  A.  ZACHARSKI (a)  the  exclusion  of  the  coupling  of  stretching  and  bending  within  the  elements,  (b)  the difficulty  of treating junctions where  all elements  are co­planar  and  (c) the  presence  of „discontinuity"  bending  moments,  which  do  not  appear  in  the  continuosly­curved actual  structure,  at  the  element  juncture  lines.  Fortunately,  many  of  these  can  be dealt  with  through  various  devices and  additional  computational  effort. The  use  of  planar  elements,  in  which  the  membrane  and  plate  bending  striffness  are derived  from  displacement  patterns  of  different  forms  cannot  insure  complete  compa­ tibility  of the assemblage which  is needed  for  convergence  of the  sequence  of finite  element solutions  to  a  true  solution.  (The  second  criterion  necessary  for  convergence  claims  the folloving:  the  displacement  functions  have  to  be  of  such  a  form  that  if  nodal  generalized displacements  are  compatible  with  a  constant  strain  condition  such  constant  strain  will in  fact  be  obtained.  Note  that  this  criterion  incorporates  in  fact  the  commonly  quoted requirement  of  rigid  body  displacements  as  these  are  a  particular  case  of  constant  strain displacement.  In  our  approach  the  second  criterion  will  always  be  satisfied). The effect  of the kinematic incompatibility  can be  expected  to  diminish with  decreasing mesh  size.  In  the  present  development  an  extensive  study  of  this  phenomenon  has  been made  showing the  essentially monotonie  convergence  to  the  true  solution  for  a  wide  range of  shell  geometries  and  external  loading  patterns. Formulating  a  nonlinear  shell problem,  one  of the  basic  decisions implied  by  efficiency considerations  is  the  selection  of  a  frame  of  reference.  For  large  displacement  analysis of  thin  shells  the  total  Lagrangian  formulation  has  been  adopted  by  most  authors.  As  it is  demonstrated  in  the  following,  however,  the  updated  Lagrangian  description  offers remarkable  simplifications  in  the  formulation  and  this  approach  is used  in  this  paper. In  the  present  study  both  the  geometrical  and  material  nonlinearities  are  taken  into account.  The  former  make  it  possible  to  solve  the  linearized  and  nonlinearized  (solved by  means  of  the  step­by­step  procedure)  structural  stability  problems.  The  latter  lead to  the  inelastic  analysis  performed  basing  upon  an  elasto­viscoplastic  material  model  as proposed  in  [12]  and  subsequently  discussed  in  [13 ­17].  This  approach  seems  to  have some  significant  advantages  over the  classical  rate­independent  elastic­plastic  formulation. First  it  produces  and  additional  numerical  effect  which  stabilizes  the  iteration  procedures used  in  the  program.  Second,  it  allows more  rational  generalizations  towards  the  inclusion of  the  dynamic  effects  into  the  solution  process.  And  third,  it  is  in  a  certain  sense  more general  approach  as  the  classical  elastoplastic  solutions  can  be  recovered  in  the  limit  as stationary  non­viscous  solutions  of  the  viscous  problem. The  presentation  is  necessarily  brief  and  no  explicit  forms  of  the  stiffness  matrices are given. More details on both the theoretical and numerical parts of the study are aVailable in  [10, 11]. 2.  Coordinate  systems The  geometry  of  the  shell  is  replaced  by  an  assembly  of  flat  elements  of  triangular and/or  quadrilateral  shape  (the  latter  elements  being  composed  of  four  flat  triangles not  necessarily  forming  one  plane),  of  Sec. 3. To  describe  the  geometry  and  stiffness  pro­ perties  of the  idealized  structure we use the  following  coordinate  systems,  Fig.  1: E F F E C T I VE  N O N LI N E AR  AN ALYSIS 27 { x y z j  global  coordinates local  coordinates { x  y ź}  local coordinates for triangular  element {x  y z]   local  coordinates for quadrilateral  element F ig.  1 (a)  global  cartesian  coordinate system  {xyz}, (b)  at  each  nodal  point:  moving  (e.g.  stepwise  updated)  surface  coordinate  system {£x £ 2 f 3} where  f3  — axis is always taken normal to the current shell surface  while the ix  — and £ 2 —   a x e s  have the tangential directions at the particular nodal point, (c)  for  each triangular  element: local  cartesian  coordinates  {xyz}, (d) for  each triangular  sub- element  of  a  quadrilateral  element: local  cartesian coordinates {xj a z a },  a  =   1, 2,  3, 4, (e)  for  each  quadrilateral  element:  local  cartesian  coordinates  {xyz}  with  the  axes  x and  y  lying  in an „ averaged"  plane „ tangential" to the element. This plane is formed by  minimizing  the  sum  of  the  squares  of  the  normal  distances  from  the  plane  to the exterior  nodes  of the  quadrilateral. The  common  Reference  frame  to  which  all  element  matrices  are  transformed  prior to  the  assembly  of  stiffnesses  is  herein  called  base  coordinates.  There  are  two  different possibilities  to choose the base  coordinates. We  describe  them brifly  in Sec. 3. 3.  F inite elements The  detailed  description  of  the  nonlinear  formulation  for  shell  elements  used  in  the present analysis  can be found  in  [11] and will be summarized below.  Since shell  behaviour is  characterized  by  both membrane  action  and bending  action it  is  essential  to  recognize 28  M.  KXEIBER,  A.  ZACHARSKI both of these in evaluating the element stiffness  properties. Two finite  elements are  available in the program  both  being based  on a flat  triangular  shell  element  as described  in  [10]: (a)  plane  triangular  element  with  eighteen  degrees  of  freedom1'  composed  of: (al)  linear  displacement  membrane  element  with  six  d.o.f.  (three  corner  nodes with  two in­plane  displacement  components  at  each  of  them), (a2)  fully  compatible  Kirchhoff  plate  element  with  twelve  d.o.f.  (one  transverse displacement  and two plate­type  rotations  at  each  of the  corners  and  at the mid­point  of the triangle).  The element  is in  fact  a  super­element  composed of  three  triangular  elements  allowing  the  C1 — transverse  displacement  (full compatibility!)  to  vary  as  cubic  polynominal  within  each  triangular  sub­ element.  In  other  words  the  transverse  displacement  approximation  for the element  is  formed  from  the  polynomial  spline  of  the  degree  3  and  smoot­ hness  1. We note  also  that  the  above  properties  imply  the linear  curvature variation  within  each  sub­triangle. Since the three  mid­point  plate­type  d.o.f.  are local  to the element  they  are  condensed (by means  of the  inverse  Gauss  elimination)  prior  to the assembly  procedure.  This  results in the total  number  of fifteen  d.o.f.  for the triangular  shell  element  (two membrane­type and  three  plate­type  d.o.f.  at  each  exterior,  corner  nodes), (b)  non­plane  quadrilateral  element  with  forty  one d.o.f.  composed  of four  triangular elements  each of them  based  upon, (bl)  quadratic  displacement  membrane  element  with  twelve  d.o.f.  (three  corner and  three  mid­side  nodes with  two in­plane  displacement  components  at  each of them), (b2)  plate  element  described  above, cf. (a2), The  external  (with  regard  to the quadrilateral)  boundaries  of the four  triangular ele­ ments  are  additionally  constrained  to  deform  linearly.  These  constraints  eliminate  the exterior  mid­side  nodes  of  the  element,  thereby  reducing  the  connectivity  (band  width) which  must be considered in the direct  solution  of the  nodal  point  equilibrium  equations. In  this way the total  number  of d.o.f.  is reduced  to thirty  three  (forty  one less  two d.o.f. at each of the four  mid­side nodes). Moreover,  since there exist in the quadrilateral  element thirteen d.o.f. which are local to the element, the corresponding  static condensation  reduces finally  the global  number  of d.o.f. to twenty. In  this  way we consistently  end up with  two finite  element:  triangular  (three  nodes) and  quadrilateral  (four  nodes)  with five  d.o.f.  at each  node. As we already mentioned in Sec. 2 the  generalized  displacements have to be  transformed from  local  to  a  common  frame  of reference  so that  the assembly  procedure  could  be  ef­ fectively  performed.  It  is  here  assumed  that  the  rotatiaral  d.o.f.  are  always  referred  to the surface  coordinates  | j  and | 3 . For the translational  d.o.f. it is left  to the user to  choose between  two common  coordinate  systems:  the global  system  {xyz}  or the surface  system {f 112 £3 }• In this way we can practically  perform  the  assembly  calculation  either in a  „mi­ I }  The  term  „degrees  of freedom"  will be further  refered  to as d.o.f. EFFECTIVE  NONLINEAR  ANALYSIS 29 xed"  coordinate  system  (the rotational  d.o.f. in the  {^  f 2 £ 3  } coordinates, the translational d.o.f.  in the  {x, y,z)  coordinates)  or in the surface  coordinate  system  alone  (all the d.o.f. referred  to the  { £ 1 1 2 1 3 } coordinates). 4.  Incremental  description  of  motion  and the solution  procedure It  is  generally  accepted  that  the  incremental  approach  is  the  most  effective  way of handling  nonlinear  structural  problems.  In the case  of elastic  structures  it is only  an  alter­ native  to  other  solution  algorithms  while  for  the inelastic  structures  the step­by­step ap­ proach is^n general unavoidable because  of the incremental  nature of the material response. The solution  algorithm accepted in the present paper relies entirely on the stepwise linearized solutions which  enable us to trace the characteristic  load — displacement  curves  describing the  nonlinear  behaviour  of the  shell  structures  analysed. At each  solution  step an iteration algorithm  over the  residual  out­of­balance  forces is planned  to be additionally  implemented to improve the solution by restoring  exact  equilibrium. The solution  algorithm is controlled by  parameters  that  are input  to the  computer  program. Fig.  2 A  typical  triangular  element  with  the  nodes  1 ­ 2 ­ 3  will  be referred  to  the local car­ tesian  coordinates  {x~yz},  Fig. 2. The x­axis  is taken  to  coincide  with  the middle  line of the  1 ­ 2  side  of  the triangle,  the  j­axis  lies  in  the  element  middle  plane  and is  directed towards  the node  3 while the i­axis is  chosen  so that  the  {xyz}  —system  be" right­handed. The  unit  vectors  of the  system  are built  a new at each  incremental  step  basing  upon new nodal  coordinates  and in  accordance  with  the above  definition. We  assume  that  the  solution  for  the  kinematic  and  static  variables  for  all time  steps from  a  time  t0  to  the current  time  f,­inclusive,  is known,  and that  the solution  for  time t+At  is required  next.  According  to  the  concept  of the updated  Lagrangian  description we  take  the  configuration  'C  at  time  t  as  a  reference  state  to  describe  the  incremental motion  *C*^  C,  Fig.  3.  Referning  to  this  configuration  we  introduce  in  the  plane  of the  element  the  second  Piola­Kirchhoff  stress  tensor  SA/l  A,A~  1,2  which  describes 30 M .  K L E I BE R ,  A.  Z AC H AR SK I F ig.  3 the  curren t stress  at  tim e  t,  an d  th e  G reen  strain  ten sor  E M ,  A,  A  =  1,2  which  describe the  curren t  strain  in  th e  element.  The  stresses  S AA   are  assumed  t o  be  equilibrated  by  t h e given  external  forces/ ^  k  — 1, 2, 3  acting  upon  the  shell. T h e increm ental changes  of  th e load  are  den oted by  Ap k ;  they  give rise  t o in crem en tal stresses  AS AA ,  incremental  strains  AE M   and  increm ental  displacements  AU A ,  AW ,  th e latter  referred  to  th e  x,  y,  ź  directions, respectively.  Th e  new  total  stresses  S' AA +AS AA   are assumed  t o  be  in  th e equilibrium  with  th e new  total  external  forces  p k +Ap k . T h e incremental  strain  displacement  relation ship  is  taken  in  t h e  form (1) which  means  th at  the  only  geometric  n on lin ear  effect  included  in  the  form ulation  is  th e influence  of  the  transverse  displacement  upon  th e  m em bran e  strain s.  I t  is  broadly kn own  t h at  the  results  obtain ed  within  th is  approxim ation  are  sufficiently  accurate  for majority  of  mildly  n on lin ear practical  problem s. We  in troduce next  th e shape functions  for  th e trian gular  plate  element  as (2) (3) where  th e  vector  A t / ( W)  collects  th e  generalized  n odal  displacements  for  th e  elem en t  con- sidered.  Its  conjugate  in tern al n odal force  vector  V[%  i  satisfies  th e virtual  work  equation of  the  form S A jE M dV  = (4 ) E F F E C T I VE  N ON LI N E AR  AN ALYSIS  31 where  F is the element volume,  dulN )  is an arbitrary  (virtual) variation of the  displacement) vector  and  6E AA   is the  corresponding variation  of E AA   taken for Au A   =  A W  =  0, e.g.  at the  beginning of the step considered  (in the configuration 'C) &E AA  =  - ^ {Su AtA   + 5u A , A )- z5w tliA .  (5) F or  the configuration  t+atC  we  write the similar  virtual  work  quation as /   (S AA +AS AA )dE AA dV  =   (U^ +AU™)^  du[f xl ,  ' (6) where  now  dE AA   is the variation of E AA   taken in the new  configuration  tJcMC, e.g. Y  (dw >A w iA +w tA dw tA )ź dw iM   (7) By using  eqs.  (4) -  (7) we immediately get A U^T  du< N> = J \ sAA 1  (3w, A wtA +  wtA dw,A)  + y   •   ( 8 ) +  AS A A  — {du A , A   + du Ai A)+   2   (dw tA w iA   +  \ v iA dw, A )- zdw, AA \ \ clV N oting  the  symmetry  of S AA   and  AS AA   and neglecting  the third- order terms  eq.  (8) can be conveniently written as A U mT  (5«w  -   {A up+A  UP)  dup,  (9) where p runs  over  the sequence  1,2, ...,  15; the summations with  respect to is implicitly assumed  on the right - hand side of eq.  (9), tyfy  )i (10)(/ v in  the nonlinear contribution to the nodal incremental forces  while -   f  AS AA   [ 1  C ^  +  < Z > i, i)- z£ fU   dV,  (11) v  L  J describes  the corresponding linear  contribution.  D efining  the resultant forces by W2  1,12  7i/ 2 N * -   j  S xx dz,  N y  =  J  S yy dz,  N xy   =  f  S xy dz,  '  (12) - Ml  - h\ l  - ft/2 and  the resultant moments by / ;/ 2  A/2  -   hj2 M x   =   /   S xx zdz,  M y =  f  Sy y zdz,  M xy   m  f  S xy zdz  (13) - A/3  - H\ l  - hl2 eqs.  (10), (11) can be transformed  to the form ,  (14) (15) 32 M.  KLEIBER,  A.  ZACHARSKI "where A  is the  triangle  area  and Bntisn3 — (16) (17) (18) (19) (20) The incremental  forces  and  moments  in  eq.  (15) are  defined  by  replacing in  the  definitions (12),  (13)  the  stress  components  S^  by  the  corresponding  incremental  stress  components ASM. Due  to  the  arbitrariness  of  the  variations  dum  eq.  (9)  yields  the  increment  of  the  in­ ternal  nodal  generalized  force  vector  as  > AVW  — AU(N>+AUm-  (22) In  order  to  maintain  the  equilibrium  at  a  given  node  of  the  discetization  mesh  the  sum ^U(int)  o j  ikg  i n t e r n a i  generalized  incremental  forces  coming  from  all the  neighbouring elements  must  be  equal  to  the  external  load  Zli?at\Ar)  = AR^  (24) where Ar  is  the  vector  of  the  generalized  displacements  of  the  whole  assemblage  of  the shell  finite  elements.  The  explicit  form  of  eq.  (24)  is  regarded  as  the  fundamental  rela­ tionship  for  the  static  incremental  analysis  of  shells.  Such  an  explicit  evaluation  of  all the matrices used  above is disregarded here; the reader is referred  to  [10,11] for  the  details. The  elemental  geometric  stiffness  matrix  is  defined  as '.A  (25) The  geometric  stiffness  matrix  was  implemental  in  the  program  in  a  slightly  modified version.  This  so­called  inconsistent  formulation  corresponds  to  the  fact  that  the  shape EFFECTIVE  NONLINEAR  ANALYSIS 33 ASXX ASxy E 1-v2 1 V 0 V 1 0 0  " 0 1-v 2 AEXX AEyy AV which  more  compactly  reads AS3xl  = D3x2AE3xL where E is the  Young  modulus  and v the  Poisson  ratio. We write  next  eq.  (1) as functions  #  used  in  the  derivation  of  A(a)  was  different  (simpler)  than  that  used  to  obtain the  constitutive  stiffness  matrix  to be discussed  below.  As stated  in  [18], for instance, such  an  approximation  leads  to  practically  acceptable  results  allowing  the  significant computer  time  savings  as compared  to the  consistent  geometric  stiffness  evaluation. We  note  that  in order  to calculate  the  force  vector A f7  (or  the  geometric  stiffness matrix  £(cr)) the  state  of  stress  at  the  beginning  of step  is the  essential information  required. In contrast, the evaluation  of the vector A t/(iV) requires the material properties of the element to  the  known. Let  us  start with the  assumption  of the linear  elastic behaviour  of the shell material,  e.q. (26) (27) AE3xl = J E S x l + J E 3 x l , (28) where AE  and A~E are the corresponding linear and  nonlinear  (with respect to the incremen­ tal displacements)  parts  of the  incremental  strain. The  vector AE  is defined  as ZlE3 x l  = [BZ3xl5-zBJzxis]Au[Vxi  (29) with  the  matrices Bm, Bb  given  in  eqs.  (20),  (21). Without  recalling  explicitly the  definition of  the  nonlinear  part  we  note  only  that  the  appropriate  expression  is independent  of the coordinate  z, cf.eq.  (1).  By  using  eqs.  (26) ­ (29)  we  arrive at AS = DAE  = D[(Bi -zBZ)AuW+/ll)]  (30) or,  performing  the  linearization, at AS = D[{Bl-zBj)Au^].  (31) The  generalization  of the  above  approach  to include  the  inelastic  analysis  capability is  achieved  here  by  specifying  an  inelastic  constitutive  law  to  be used  instead  of  the  elastic law  given  by  eq.  (26).  The  analysis  will  be  based  upon  the  elasto­viscoplastic  constitutive assumptions first  proposed in  [12] and later explored numerically by many authors,  [13­ 17]. The  elastic­viscoplastic  material  is defined  by  the  following  constitutive  relation AS  = D(AE-AE\  =  yAt 1 [iSxx- Syyf a%2  i/^C2  i  CC2  _ 1 C  O  I_a/;C2xy 2 > J X X  — S 2Syy — S yy yy — Sxx 6S xy (39) The  square  root  in  the  denominator  on  the  right­hand  side  of  (39)  is  introduced  to  nor­ rnalize the vector  ­—­. dS For  hardening  materials  the  yield  limit  ,  (42) /IE&"> = y l / 3 P E ^ ) 3 +  (AE$»y+ AE™AE™  +  (AE%ń )2.  (43) (the  latter1 definition  takes  into  account  the  inelastic  incompressibility  of  the material), and assuming that a 0  =  50( E U ) ) , we define the hardening modulus as Using eqs.  (31) and  (32) we arrive at AS  =  D [B) n  - zBj)Au^ - AE^ ].  (45) Performing  the integrations  (12), (13), for the value of AS given by eq.  (45) we  get incre- mental  membrane  forces AN 3  x ,  =  hDi,  3 flff,.M A Ufti!  -  AN *  (46) and the incremental bending moments 1 L ( 4 7 ) where  the „ initial" generalized  forces  AN *  and AM are defined by ZliV* »> D  J  AE^ Mz,  (48) - A/2 ft/ 2 zIM* m,  D  j  AW ^ zdz.  (49) - A/2 Similarly  as  before,  cf.eq.  (15)  the expression  for  Z lZ / ^ ! can now be obtained  in  the form where  / c(e) is the elastic elemental stiffness  matrix used in the linear version  of the program, [10]  and  based  now  on the  current  configuration  at the  beginning  of the  load  step,  and AU iN ) *  forms  a vector  of the „ additional" incremental nodal forces  and moments  defined by =   1 [B m AN *- B b AM*]dA.  (51) A The  last  two  expressions  suggest  an iterative  „initial  force"  procedure to be performed at each incremental step. The  application  of  the  direct  assembly  procedure  leads  to the fundamental  matrix equation  describing  the large  displacement  inelastic  shell  problem in the  form [K^ +K m ]Ar  =  ARW +AR*  (52) 3 *  • - . . ' • , 36  M .  K L E I BE R ,  A.  Z AC H AR SK I in  which  th e matrices K{e),  K^   an d th e vectors  Ar, AR* are t h e  global  coun terparts of th e elemental  matrices  &(e), k<- a) an d th e vector  AuS N\   AUiW*  while  zJ i?( e xt )  is  th e  external load  vector,  ef.  eq. (23).  Because  in th e case  of  inelastic  m aterial  properties  th ere  is n o linear  stress  distribution  across  the  thickness  of th e  shell  (and, in fact,  n o functional ap- proximation  for  such  a  stress  distribution  can be  rationally  assumed  a  priori),  th e shell element  is  considered  as  composed  of  layers.  To effectively  find  the inelastic  elemental nodal  forces  we proceed as  follows: 1.  F rom  eq. (39) the incremental  inelastic  strain  is  calculated  for  th e given  stress  S in all  layers  at each nodal poin t. 2.  U sing  th e trapezoidal  rule  the across- thickness  in tegration s  are performed  accordin g to  eqs. (48),  (49)  resulting  in th e  initial  m em brane  forces  AN *  an d ben din g  m om en ts AM*. 3.  Eq. (41) is  used  to  find  t h e  initial  n odal  generalized  forces  AU(N )*.  N um erical  area integration is carrried  ou t  at this  stage  by simlpy  assuming 3 A #<">* =  i -  A J T (B mi AN ?  - B bl AMf)  (53) / = i where  the index  „ / ", i =   1, 2, 3,  refers  to th e n odal values  of th e trian gular  element. 4.  Eq. (52) is  obtained  as a result  of th e  direct  assembly  procedure. I t is  th en  solved  for AT ,  and this  calculation is followed  by: —  evaluation  of the new  increm ental stresses  by eq.  (45), —  evaluation  of  th e  new total  stresses  according  to th e kn own  stress  accum ulation procedure, —  evaluation  of the new  „ in itial  lo ad "  vector —  solution  of  eq.  (52) for  the improved  value  of  the in crem en tal  displacem ent AT . T h e  iteration  process  is  continued  un til  convergence  is  achieved  u p t o a  desired ac- curacy.  I n the  present  study  th e  convergence  is m on itored by using  alternatively  th e  con- dition s: t d ( r ) where  k stan ds for the  k- th.  iteration . 5.  Linearized  stability  analysis F o r  one- parameter (proportion al) loadings  th e fundam ental  m atrix  equation  describing the  static  problem  of  elastic  th in  shells  can be  presented  in  a  convenient,  approxim ate form  as, cf.eq. (52) =  A«Łext>  (54) EFFECTIVE  NONLINEAR  ANALYSIS  37 where  K(a)(o*)  signifies  symbolically  dependence  of  the  initial  stress  matrix  on  the  stress state  a*  which  corresponds  linearly  to  a  reference  external  load  2?(JX<) while  X  is the  scalar load  multiplier.  The  approximation  in  eq.  (54)  consists  essentially  in  using  the  initial coordinates  for  setting  up  the  stiffness  matrices  K(c)  and  jfiT(cr),  which  is rigourously  valid for  small  deformation  problems  only. The  purpose  of  the  linearized  analysis  of  the  shell  stability  is  to  check  the  uniqueness of  the  solution  A r  of  eq.  (54)  for  each  given  value  of  the  parameter  X.  The  ponits  of  such a  non­uniqueness  are  called  bifurcation  (or  branching)  points  on  the  primary  equilibrium path  in  the  load­displacement  space  X—Ar.  According  to  the  definition  at  the  bifurcation point  1  =  Xcr the  relations  hold [#«> +  Ac rK«>*)] Ay,  =  XcrR + lcrK<°\=0.3 clumped  t  = 12.7mm e d 9 e  L=  508 mm R=   2540 mm 6  =  0.2  r ad PtKNl 2 - 1  - F ig.  8 1401 1 1, present 121] i solution p m o x=2.a 122) 1 KN  Ę ,ax=2.23KN \ 10  - ZAimml EFFECTIVE  NONLINEAR  ANALYSIS 41 completely  free.  The  analysis  was  continued  up  to  the  point  at  which  singularity of  the  total  stiffness  matrix  appeared.  One­quarter  of  the  panel  was  discretized by  a  6 x 6 finite  element  mesh.  The  comparison  of  the  present  results  with  those discussed in  [21], [22] is shown in Fig. 8. Good  agreement  of the results is observed. V.  Geometrically  nonlinear  analysis of another  cylindrical  shell, Fig. 9. i  i  r h=3.175mm E=3102.75 N/mm2 ^=0.3 present solution [23] [241 -2 -8 -10-4 -6 ZA[mm] Fig.  9 The  circular  cylindrical  shell  portrayed  in  Fig.  9.  is  clamped  along  all  four  boun­ daries and  subjected  to uniform inward radial loading. One­quarter  of the shell was discretized  by  an  uniform  8 x 8  finite  element  mesh  and  32 equal  load  steps were applied.  Fig.  9  shows  a  very  good  agreement  of  the  present  results  as  compared against those given in  [23], [24]. Almost the same results were obtained in the present study by using an uniform  6 x 6  mesh and 40 load increments. VI.  Geometrically  nonlinear  analysis  of  a  spherical  shell,  Fig. 10. The  shell is  subjected  to  a  concentrated  load  at  the  apex;  all  edges are hinged and immoveable. One­quarter  of the  shell was discretized by a 6 x 6 finite  element mesh and  40 load  increments were used. As a matter  of fact  the problem was considered under the apex displacement  control rather than under the force control. This made it possible to  get through the limit point  on the load­displacement  diagram without any difficulties.  The results were foundjto  be sufficiently  accurate, cf. Fig. 10 cf. [23], their  further  improvement  is possible by  simply using more load increments. VII.  Limit load  analysis  of a  quadratic plate, under uniform  loading, Fig. 11. The upper  and lower limit load  estimates  are given as, cf.  [25] 24M0 Pk =  —FT  V. Mo  = —  aoh 2. 1/3 42 M.  KLEIBER,  A.  ZACHARSKI 1.2 1.0 0.8 0.6 OM 02 0 a=508rnm t  =2.54 mm E = 6 8 . 9 5 * i 0 3 N / m m 2 \>=0.3 30=248.22 N/mm* w o =mippoint  deflection w j = f o r  first  yielding i / / , 1 I I ? 2 x 2 4x4 5x5 i 3 [25] [251 present 1 4 !  1  1 11 ! solution i 5  6 W0/W() Fig. 11 The 5 x 5 finite  element mesh was used  each element  consisting  of 11 layers  assumed for  the observation  of the across­thickness  plastic  zone  development.  The  results of  the  analysis are shown in Fig.  11,  12. VIII.  Inelastic,  large  displacement  analysis  of a spherical  cap,  Fig.  13. The  shell  is  lunged  at the boundary  and subjected  to  a  concentrated  force  acting at  the appex.  The material  and geometric  data  are the same  as in  Example  II, cf. Fig.  5. The material  is assumed  to be ideally  plastic  with  a0  =  137.9  N/mm 2 . The displacement  control  of the process was used.  The first  five  incremental  apex  displa­ cements  were  assummed  to be  equal  to  ­0.0254  mm. which  was followed  by the m_ MMm ^ yielded  at  8  layers yielded  at  6  layers yielded  at  4  layers p/p,=1 Fig.  12 1 2  5  8  11  14  17  20  x Fig.  13 1.0 2.0  3.0"  4.0 ­zlmml  . Fig.  14 [43] 44 M.  KXEIBER,  A.  ZACHARSKI Fig.  15 [N/mm2] 0.003 0.002 0.001 geometrically  nonlinearj solution linear  solution [26] present solution X  E=210Q0N/mm; 0=0.0 100  wB[mm] 25 steps of  —0.127 mm. The present  results  are  shown in  Fig.  14 and  15 and  compa­ red  with those  given  in  [7]. IX.  Inelastic,  large  displacement  analysis  of  a  cylindrical  panel,  Fig.  16. The  shell  is  subjected  to  the  uniformly  distributed  loading  acting  in  the  negative z­direction. The  load  was  assumed  to  increase  from  0 to  0.00315  N/rnm3  in  the  24  equal  incre­ ments.  The  results  are  given  in  Fig.  16  cf.  [26]'. X.  Linearized  stability  of  elastic  quadratic  plates,  Tabl.  3. The linearized  stability  formulation  described  in  Sec.  5 was  the  base  for  calculating the  buckling  loads  for  the  two  differently  supported  plates.  One­quarter  of  the EFFECTIVE  NONLINEAR  ANALYSIS 45 plate was analysed.  It is seen that the accuracy of the  solution is strongly dependent on  the number of finite  elements used in the calculation. This can be partly attri- buted to the simplified  assumptions used in deriving  the geometric  stiffness  matrix. Table  3 JT£ E  [ _h_] 3 simply suported plate a - I I L- a- rT 7 a ) clumped plate a .mroi a - t iłiil a ł Mesh 6x6 9x9 14x14 6x6 9x9 14x14 19x19 Coefficient K 4.192 4.044 4.016 5.947 5.506 5.374 5.337 error 4.8 1.1 0.4 12.0 3.7 1.2 0.5 References 1.  R. H .  G ALLAG H ER,  Shell elements,  P roc. World Congress on F inite Element  Method  in Struct.  M ech. Bournemouth,  England,  1975. 2.  M . D .  OLSON ,  T. W.  BEARD EN ,  A  simple flat  triangular shell  element revisted,  Int. J.  N um. M eths. Eng.  14,  51  -   68,  1979. 3.  R. W.  CLOU OH ,  C. P.  JOH N SON ,  A  finite  element approximation  for  the  analysis of  thin shells, Inst. Y.  Sol.  Struct.  4,  43 -  60,  1968. 4.  H . WEN N ERSTRÓM, N onlinear shell analysis performed with flat  elements, Proc. I n t. Conf. F inite Elements in  N onlinear  Mechanics,  G eilo, N orway,  August  1977, pp. 285 -  302. 5.  J. H . ARQYRIS, P. C. D U N N E , G . A.  MALEJAN N AKIS,  E. SCH ELKE, A  simple triangular facet  shell element with applications to  linear and non- linear equilibrium and elastic stability problems,  Comp. M eths. Appl. M ech.  Eng.  10,  371- 403  and  11, 97- 131,  1977. 6.  M. KLEIBER, A  triangular finite  element for  large deformation elasto- plastic  analysis of  arbitrary shells, Bull.  Acad. Polon.  Sci. Ser. Sci. Techn. 26, N o 2, 61- 71,  1978. 7.  J.  H .  ARG YRIS, H .  BALMER,  M .  KLEIBER, U .  H IN D EN LAN G ,  N atural formulation of  large inelastic de- formations for  shells of arbitrary shape- Application of the T RUMP  element,  Comp. M eths. Appl. Mech. Eng.  22,  361  -  389,  1980. 8.  G .  H ORRIG MOE, P . G.  BERG AN ,  N onlinear  analysis of free- form  shells by flat  finite  elements, Comp. M eths.  Appl.  Mech.  Eng.  16,  11- 35,  1978. 9.  M .  KLEIBER,  H .  STOLARSKI, N umerical analysis  of  elastic- plastic  shells in the range  of  large  static de- formations,  (in  Polish),  Institute  of  F undamental  Technological  Research,  Warsaw  1976, Results obtained  in  the course  of  research  sponsored  by  the  Centrum  Techniki  Okrę towej,  contract  EU / B/ 246/ 74. 10.  M . KLEIBER, A. ZACH ARSKI, SHEL 1N —  linear finite  element analysis of thin free- form shells,  (in Polish), Institute  of  F undamental  Technological  Research,  Report  N o 51/ 1978,  Warsaw. 11.  A.  ZACH ARSKI,  N onlinear static  analysis of  thin shells of arbitrary shape,  (in Polish),  P h. D . thesis,  In- stitute  of  F undamental  Technological  Research,  Report  N o 15/ 1982, Warsaw. 12.  P .  PERZYN A,  T he constitutive  equations for  rate  sensitive plastic  materials, Quart.  Appl.  M ath.,  20, 321  -  332,  1963. 46  M .  KLE I BE R ,  A-.  ZACH ARSKI 13.  O. C.  Z I E N KI E WI C Z ,  I . C.  CORMEAU ,-  Visco- plasticity  and creep  in  elastic  schkls- a unified  numerical  so- lution approach,  I n t .  J. N u m .  M eths. En g. 8, 821 -  845, 1974. 14.  S.  N AG ARAJAN ,  E. P .  P OP OV,  Plastic  and viscoplastic analysis of  axisymmetric  shells, I n t . J.  Sol.  Struct. 11,  1 -   19,  1975. 15.  B.  KRAKELAN D ,  L arge  displacement  analysis  of  shells  considering  elasto- plastic  and  elasto- viscoplastie materials,  Report N o 77 -  6, D ec. 1977, D iv.  of  Struct. M ech., The U niv. of  Trondheim,  N orway. 16.  M .  KLEIH ER, N atural finite  elements  and  large  deformation  elasto- viscoplasticity,  Bull.  Acad.  P olon . Sci.  Set- . Sxci.  T ech n ,  26,  73 -   81,  1978. 17.  J. H .  AR G YR I S, J. St. D OLTSI N I S, M . KLE I BE R , Incremental discretized formulation  in non- linear mechanics and  finite  strain  elasto- plasticity- natural approach  P art .  I I , C om p.  M eth s.  M ech.  E n g. 14, N o 2,  259 - 294,  1978. 18.  R. W.  C LÓU G H ,  C.  A.  F E LI P P A,  A  refined  quadrilateral  element  for  analysis  of  plate  bending,  P ro c. Conf. Matrix  Methods  in  Struct.  Mech., Wright—P atterson ,  Ohio,  1968. 19.  R.  BARES,  T ables for  the  analysis of  plates,  slubs  and diaphragms  based  on  the  elastic  theory,  2nd ed, , Bauverlay, G mbh., Wiesbaden,.  1971,  G erman- English edition. 20.  W.  KAN OK- N U KU LCH AI,  R. L.  TAYLOR, T. I.  H U G H ES, A  large deformation formulation for  shell analysis by  the  finite  element method,  Comp. Struct., Vol.  13, 1981. 21.  K. J.  BATHER,  S.  BOLOURCHI, A geometric and material nonlinear plate and shell element,  Comp. Struct., Vol.  11,  1980. 22.  H .  PARISH ,  L arge displacements  of  shells including  material nonlimarities, Comp.  M eth. Appl.  M ech. Engng.,  Vol.  27,  1981. 23.  G .  D H ATT,  Instability  of  thin shells  by  the finite  element method,  Symp.  of  I n t.  Assoc.  of  Shell  Struct., Vienna,  1970. 24.  C. TAH I AN I , L.  LACH AN CE, L inear and non- linear analysis of  thin shallow  shells  by  mixed finite  elements, Comp.  Struct.,  Vol.  5,  1975. 25.  P . BERG AN , N onlinear analysis of plates considering geometric and material effects, D iv. of  Struct. M ech. the  N orwegian  Inst. of  Techn., Report N o 72 -  1,  May  1972. •  26.  P . BERG AN ,  G . H ORRIG MÓE, B.  KRXKELAN D,  T.  SPREID E, Solution techniques for  non- linear finite element problems,  Int. J. N um. M eth. Engng., Vol.  12,  1978. P  e 3  io  M  e 9*< I > E K T H BH Ł ia  H E JI H H E flH Ł lK  AH AJI H 3  T O H K H X  OE OJI O^E K  n P O H 3BO J I bH O fł ctOP M BI B  paSoTe  oScyjKfleHbi  OCHOBHŁIC acnei n o jiyq e H B ix  pe3yjibT aT O B.  y ^ r e H b i n p o 6jieiwbi  6ojiBniH X  n ep eM em eH H H   O G O J I O ^ B K  H   H e y n p y r n x  C B O H C T B  M a i e p H a n a .  O 6cyjK fleH M   T o i o i t e P a6oia S t r e s z c z e n i e W  pracy  przedstawiono  podstawy  statycznej,  nieliniowej  analizy  powł ok  dowolnego  kształ tu  metodą elementów  skoń czonych.  Zastosowano  wzglę dnie  proste  elementy  skoń czone,  upraszczają c  w  ten  sposób wprowadzanie  danych  wejś ciowych  oraz  interpretację   otrzymywanych  wyników.  U wzglę dniono  proble- matykę   duż ych  przemieszczeń  powł ok  oraz  niesprę ż yste  wł asnoś ci  materiał u. Omówiono  również  zagad- nienie analizy  zlinearyzowanej  statecznoś ci  powł ok. Praca  zilustrowana  jest  licznymi  przykł adami.