74 Fernando Flores and Eugenio Onate N x oN t3= loxx (19) From these expressions it is also possible to compute in the original configuration the element area A,the outer normals (n,n2)'at each side and the side lengths I.Equation (19)also allows to evaluate the thickness ratio A in the deformed configuration and the actual normal t3. Direction t3 can be seen as a reference direction.If a different direction than that given by (19)is chosen,at an angle 0 with the former,this has an influence of order 02 in the computation of the curvatures (see (23)below).This justifies (19) for the definition of ts as a function exclusively of the three nodes of the central triangle,instead of using the 6-node isoparametric interpolation. The numerical evaluation of the line integral in(16)results in a sum over the integration points at the element boundary which are,in fact,the same points used for evaluating the gradients when computing the membrane strains.As one integration point is used over each side,it is not necessary to distinguish between sides(i)and integration points(Gi). The explicit form of the gradient evaluated at each side Gi(3)from the quadratic interpolation is P2 Ni.2N2N.2N+3.2 (20) Pi+3 We note again that the gradient at each mid side point Gi depends only on the coordinates of the three nodes of the central triangle and on those of an additional node in the patch,associated to the side i where the gradient is computed. In this way the curvatures can be computed by LM 0 0 L (21) P2 An alternative form to express the curvatures,which is useful when their varia- tions are needed.is to define the vectors 3 hij (22) k=1 This gives K=h订·t3 (23) The variation of the curvatures can be obtained as 0 =空 0 LM LM +[a] (L出h+L) (LM22+L22) t3·dui) (24) (LMi2+L20i2)
74 Fernando Flores and Eugenio O˜nate t3 = ϕM 1 × ϕM 2 |ϕM 1 × ϕM 2 | = λ ϕM 1 × ϕM 2 . (19) From these expressions it is also possible to compute in the original configuration the element area A0 M , the outer normals (n1, n2) i at each side and the side lengths l M i . Equation (19) also allows to evaluate the thickness ratio λ in the deformed configuration and the actual normal t3. Direction t3 can be seen as a reference direction. If a different direction than that given by (19) is chosen, at an angle θ with the former, this has an influence of order θ2 in the computation of the curvatures (see (23) below). This justifies (19) for the definition of t3 as a function exclusively of the three nodes of the central triangle, instead of using the 6-node isoparametric interpolation. The numerical evaluation of the line integral in (16) results in a sum over the integration points at the element boundary which are, in fact, the same points used for evaluating the gradients when computing the membrane strains. As one integration point is used over each side, it is not necessary to distinguish between sides (i) and integration points (Gi). The explicit form of the gradient evaluated at each side Gi (3) from the quadratic interpolation is ϕi 1 ϕi 2 = Ni 1,1 Ni 2,1 Ni 3,1 Ni i+3,1 Ni 1,2 Ni 2,2 Ni 3,2 Ni i+3,2 ⎡ ⎢ ⎣ ϕ1 ϕ2 ϕ3 ϕi+3 ⎤ ⎥ ⎦ . (20) We note again that the gradient at each mid side point Gi depends only on the coordinates of the three nodes of the central triangle and on those of an additional node in the patch, associated to the side i where the gradient is computed. In this way the curvatures can be computed by κ = 2 *3 i=1 ⎡ ⎣ LM i,1 0 0 LM i,2 LM i,2 LM i,1 ⎤ ⎦ t3 · ϕi 1 t3 · ϕi 2 (21) An alternative form to express the curvatures, which is useful when their variations are needed, is to define the vectors hij = *3 k=1 LM k,iϕk j + LM k,jϕk i (22) This gives κij = hij · t3 (23) The variation of the curvatures can be obtained as δκ = 2 *3 i=1 ⎡ ⎣ LM i,1 0 0 LM i,2 LM i,2 LM i,1 ⎤ ⎦ '*3 i=1 Ni j,1(t3 · δuj ) Ni j,2(t3 · δuj ) + Ni i+3,1(t3 · δui+3) Ni i+3,2(t3 · δui+3) ( − *3 i=1 ⎡ ⎣ (LM i,11 11 + LM i,22 11) (LM i,11 22 + LM i,22 22) (LM i,11 12 + LM i,22 12) ⎤ ⎦(t3 · δui) (24)
Applications of a Rotation-Free Triangular Shell Element 75 where the projections of the vectors hi;over the contravariant base vectors have been included e%=hp0,a,i,j=1,2 (25) with p州=入p当xts (26) p9=-入pH×tg (27) In above expressions superindexes in L and u refer to element numbers whereas subscripts denote node numbers.As before the superindex M denotes values in the central triangle (Fig.1.a).Note that as expected the curvatures (and their variations)in the central element are a function of the nodal displacements of the six nodes in the four elements patch. Details of the derivation of (12)and(24)can be found in [10].The explicit expressions of the membrane and curvature matrices can be found in [12].The derivation of the element stiffnes matrix is described in [10,12].Also in [10,12] details of the quasi-static formulation and the fully explicit dynamic formulation are given. It must be noted that while the membrane strains are linear the curvature strains are constant.A full numerical integration of the stiffness matrix terms requires three points for the membrane part and one point for the bending part.Numerical experiments show that: when using one or three integration points the element is free of spurious energy modes and passes the patch test for initial curved surfaces the element with full (three point)integration leads to some membrane locking.This defect dissapears if one integration point is used for the membrane stiffness term. It can also be observed that: for large strain elastic or elastic-plastic problems membrane and bending parts can not be integrated separately,and a numerical integration trought the thick- ness must be performed for explicit integrators (hydro codes)is much more effective to use only one integration point for both the membrane and bending parts. Above arguments lead to reccomended the use of one integration point for both membrane and bending parts.This element is termed EBSTI to distinguish from the fully integrated one. 2.3 Boundary Conditions Elements at the domain boundary,where an adjacent element does not exist,de- serve a special attention.The treatment of essential boundary conditions associated to translational constraints is straightforward,as they are the degrees of freedom of the element.The conditions associated to the normal vector are crucial in this for- mulation for bending.For clamped sides or symmetry planes,the normal vector t3 must be kept fixed(clamped case),or constrained to move in the plane of symmetry (symmetry case).The former case can be seen as a special case of the latter,so we
Applications of a Rotation-Free Triangular Shell Element 75 where the projections of the vectors hij over the contravariant base vectors ϕ˜ Mα have been included α ij = hij · ϕ˜Mα , α, i, j = 1, 2 (25) with ϕ˜M1 = λ ϕM2 × t3 (26) ϕ˜M2 = −λ ϕM1 × t3 (27) In above expressions superindexes in LM j and δuk j refer to element numbers whereas subscripts denote node numbers. As before the superindex M denotes values in the central triangle (Fig. 1.a). Note that as expected the curvatures (and their variations) in the central element are a function of the nodal displacements of the six nodes in the four elements patch. Details of the derivation of (12) and (24) can be found in [10]. The explicit expressions of the membrane and curvature matrices can be found in [12]. The derivation of the element stiffnes matrix is described in [10, 12]. Also in [10, 12] details of the quasi-static formulation and the fully explicit dynamic formulation are given. It must be noted that while the membrane strains are linear the curvature strains are constant. A full numerical integration of the stiffness matrix terms requires three points for the membrane part and one point for the bending part. Numerical experiments show that: • when using one or three integration points the element is free of spurious energy modes and passes the patch test • for initial curved surfaces the element with full (three point) integration leads to some membrane locking. This defect dissapears if one integration point is used for the membrane stiffness term. It can also be observed that: • for large strain elastic or elastic-plastic problems membrane and bending parts can not be integrated separately, and a numerical integration trought the thickness must be performed • for explicit integrators (hydro codes) is much more effective to use only one integration point for both the membrane and bending parts. Above arguments lead to reccomended the use of one integration point for both membrane and bending parts. This element is termed EBST1 to distinguish from the fully integrated one. 2.3 Boundary Conditions Elements at the domain boundary, where an adjacent element does not exist, deserve a special attention. The treatment of essential boundary conditions associated to translational constraints is straightforward, as they are the degrees of freedom of the element. The conditions associated to the normal vector are crucial in this formulation for bending. For clamped sides or symmetry planes, the normal vector t3 must be kept fixed (clamped case), or constrained to move in the plane of symmetry (symmetry case). The former case can be seen as a special case of the latter, so we