Microsoft Word - 1.doc Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 LARGE DISPLACEMENT ELASTIC-PLASTIC ANALYSIS OF LAMINATED COMPOSITE PLATES UNDER IN-PLANE COMPRESSIVE LOADS Dr. Husain M. Husain Professor in Civil Engineering, University of Tekrit Dr. Nameer A. Alwash Professor in Civil Engineering, University of Babylon Dr.Haider K. Ammash Lecturer. in Civil Engineering, University of Al-Qadisiya Abstract A nonlinear finite element method is adopted for the large displacement elastic-plastic static analysis of anisotropic plates under in-plane compressive loads. The analysis is based on the two- dimensional layered approach with classical and higher order shear deformation theory with five, seven, and nine degrees of freedom per node. Nine-node Lagrangian isoparametric quadrilateral elements are used for the discretization of the laminated plates. Effects of orthotropy of individual layers, through-thickness shear deformation, fiber’s orientation angle, and fiber waviness on the large displacement elastic-plastic static analysis are considered. The plate is analyzed with a range of number of sequences (k) of sine wave fiber (1-12) and with a range of the amplitude of fiber path (∆) of sine wave fiber (0.0-0.5). The conclusion it is shown that the behavior of the laminated plate is very sensitive to the shape of fibers (straight or sine wave), also the behavior of the plate with sine wave fiber depends on the amplitude of the fiber and the number of sequences of the fiber, and so the capacity of the laminated plate with sine wave fiber and under in-plane compressive load in the direction of waviness is higher than the capacity of the plate with sine wave fiber and under in- plane compressive load orthogonal to the direction of waviness by approximate value (42%). Keyword: Large Displacement, elastic-plastic analysis, finite element method, composite plate, in-plane compressive loads. حمل ضغط في تأثيراللدن للصفائح غير متماثلة الخواص تحت- المرنةالكبير اإلزاحةليل حت المستوي حيدر آاظم عماش.د القادسيةجامعة/مدرس نمير عبد األمير علوش. د جامعة بابل/أستاذ مساعد حسينحسين محمد .د جامعة تكريت/أستاذ -:الخالصة اللدن االستاتيكي للصفائح غير متماثلة الخواص تحت - المرن ةالكبير اإلزاحة تم تقديم طريقة العناصر المحددة الالخطية لتحليل واعتمدت ) two-dimensional layered approach( الطبقة ثنائي البعد الدراسة طريقة هتبنت هذ . حمل ضغط في المستوي مع ) classical and higher order shear deformation theory(نظرية التشوهات القصية الكالسيكية وذات المرتبة العليا . ةي ذي العقد التسع لتمثيل الصفائح الطبق )Lagrangian ( الآرانجت حرية لكل عقدة، تم توظيف عنصراسبع وتسع درجوخمس، الكبيره اإلزاحة على تحليل األلياف تموج واأللياف زاوية تدوير الخواص للطبقات المنفردة و تعامدخذ بنظر االعتبار تأثير ُأ قدو من ارتفاع قمة االلياف جال، وم)12-1( من تسلسل االلياف المتموجة مجالتم تحليل الصفيحة مع . اللدن االستاتيكي-المرن Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 2 مستقيم او ( لشكل االلياف ة حساس جداً ي طبقرف الصفائح الُ تص ن من النتائج المستحصلة، يمكن مالحظة ا ). 0.5-0.0(المتموجة قابلية وآذلك ان ’عتمد على ارتفاع قمة االلياف وعدد تسلسل التموجات مف الصفائح ذات االلياف المتموجة تصّروأيضا ). متموج ة والمتموجة ي التحمل للصفائح الطبقبقة والمتموجة االلياف وتحت حمل ضغط في اتجاه التموج اآثر من قابلية ّطالتحمل للصفائح المُ .)%42(االلياف وتحت حمل ضغط عمودي على التموج بحوالي Notations Symbol Description a, b Plate dimensions in x and y-directions, respectively. [B] Strain-nodal displacement matrix. D Flexural rigidity = ( )23 112 vEt − . Ei Modulus of elasticity in i-direction. Ef Modulus of elasticity of fiber. Em Modulus of elasticity of matrix. Fi, Fij Strength tensors of the first and second order, respectively. {F} External load vector. F Yield function. G Shear modulus. hL-hL-1 Distance from plate middle surface to the upper and lower surface of Lth lamina. h Plate thickness. [ ]oK Constant linear elastic stiffness matrix [ ]LK Initial or large displacement matrix [ ]σK Initial stress stiffness matrix [ ]0TK Tangent stiffness matrix. xyyx MMM ,, Bending and twisting moments (per unit width) (on yz, xz, and both yz and xz- sections). *** ,, xyyx MMM Higher order bending and twisting moments (per unit width) (on yz, xz, and both yz and xz-sections). xyyx NNN ,, In-plane stress resultants (per unit width) (on yz, xz, and both yz and xz-sections). *** ,, xyyx NNN Higher order in-plane stress resultants (per unit width) (on yz, xz, and both yz and xz-sections). Px In-plane applied load in x-direction. Qx, Qy Transverse shearing forces (per unit width) (on yz and xz-sections). wo Amplitude of initial imperfection. x,y,z Coordinates. o ijγ Shear strain in ij-plane at middle surface. *o ijγ Higher order shear strain in ij-plane at middle surface. { }ε Strain vector. { }oε Middle surface strain vector. iε Normal strain in i-direction. o iε Normal strain in i-direction at middle surface. *o iε Higher order normal strain in i-direction at middle surface. Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 3 ηξ , Curvilinear coordinates system. θ Fiber’s orientation angle. θx,θy Rotations of transverse normals in the (xz) and (yz) planes. ** , yx θθ Higher order rotations of transverse normals in the (xz) and (yz) planes. o iκ Bending curvature in i-plane at middle surface. o ijκ Bending curvature in ij-plane at middle surface. *o iκ Higher order bending curvature in i-plane at middle surface. *o ijκ Higher order bending curvature in ij-plane at middle surface. vi Poisson’s ratio in i-direction. oσ Yield stress of steel Introduction The composite material is made of two or more macro-constituent's materials essentially soluble or mixable into each other, usually a reinforcing material supported in a compatible matrix, where the sum of properties of each constituent taken separately, are assembled in prescribed amounts to achieve specific physical and chemical properties [Jones,1999](9). The reinforcement materials may form as continuous or discontinuous fibers, flakes, fillers or particles embedded in the matrix material. Fibers in various forms (mat yarn, woven roving, and chopped strands) are inherently much stiffer and stronger than the same material in bulk form. The matrix material works as a binder material giving the composite a protection and supports its bulk form and stress transfer when the fiber is broken. Typically, the matrix is of considerably lower density, stiffness and strength than those of the fibers. Layered composite material plates are extensively used in the construction of aerospace, civil, marine, automotive and other high performance structures, and during the operation of this structure, it is subjected to static and dynamic loads, as shown in Figure (1). Therefore, there exists a need for investigating the response of layered (laminated) composite material plates subjected to such types of loading. The finite element method has been applied with great success to geometrical and material nonlinearities in continuum and structural problems. The geometric nonlinearity is modeled by well-known formulations, the total or updated Lagrangian coordinates, while success in modeling the material nonlinearities depends on the validity of the constitutive models to be used. Elseifi [1998](7) presented nonlinear finite element method for the post-buckling analysis of stiffened composite panels with geometric imperfections. A four node, six degrees of freedom per node, rectangular, conformal element was used. Transverse shear effects were neglected since the width to thickness ratio of the panels under consideration is over 500. A maximum stress failure criterion was added to the finite element code to predict the panel post-buckling failure load. A new integration technique that mixes symbolic closed-form function manipulation and Gaussian quadrature numerical integration had been introduced in order to reduce the required computation time for each analysis. His study did not take the effect of fiber orientation, number of fibers, and type of initial imperfection. Shukla and Nath [2000](15) presented a post-buckling analysis of shear deformable cross-ply laminated rectangular plate subjected to the combination of in-plane edge compressive mechanical loading and thermal loads due to a linearly varying temperature across the thickness. Their formulation was based on the first order shear deformation theory and Von- Karman type nonlinearity. They observed that their present method was quite efficient in obtaining the buckling and post-buckling response of a laminated composite plate under thermomechanical Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 4 loading. More recently, Zou and Qiao [2002](16) presented a higher order finite strip method for post-buckling behavior of imperfect composite plates subjected to progressive end shortening. The arbitrary nature of initial geometric imperfection induced during manufacturing was accounted for in the analysis. The nonlinear equilibrium equations were solved by Newton-Raphson method. That study showed that the post-buckling behavior of an imperfect composite plate depends not only on the material lay up, snap-to-thickness and anisotropy of the laminate, but also on the direction of induced out-of-plane imperfection. From the preceding review of literature, it is clear that there is no study which considers the nonlinear static analysis of isolated laminated plate under axial compression load by taking into account the effect of type of fiber (straight or wavy). There is also a little amount of literature that takes into account the higher order displacement model of nine degrees of freedom per node with different types of lamination. Laminated Plate Theories A laminated plate is a series of laminas bonded together to act as an integral structural element. Thus, a laminate is not a material but instead a structural element with essential features of both material properties and geometry. The stiffness and strength of such a composite material with structural configuration are obtained from the properties of the constituent laminas, and thus the macromechanical behavior of a laminate is the main topic of this section. The lamination so described can be considered as a single layer with "rule of mixtures" representation of the interaction between the multiple laminas in a plate or shell [Jones, 1999](9). There are two categories of theories, equivalent single layer and three dimensional elasticity theories. In the first category, the material properties of the constituent layer are smeared to form a hypothetical single layer whose properties are equivalent to through thickness integrated sum of its constituents, and this category contains classical lamination theory, first order shear deformation theory, and higher order shear deformation theory. The higher order shear deformation theories are more efficient to represent the transverse shear deformation, through-thickness displacement and strains. The assumption of a higher order plate theory can also be used within the equivalent layer formulation [Jones, 1999](9). The strain expressions derived from the displacement field were considered by [Ali, 2004](2) with nine degrees of freedom per node as follows: ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( )tyxwtzyxw tyxztyxvztyxztyxvtzyxv tyxztyxuztyxztyxutzyxu o yoyo xoxo ,,,,, ,,,,,,,,,,, ,,,,,,,,,,, ** ** = θ++θ+= θ++θ+= 32 32 (1) in which the parameters (u, v, w, θx, θy, * xθ , and * yθ ) are defined previously, * ou , and * ov are the corresponding higher order terms in Taylor's series expression and they are also defined at the middle plane. The strain-displacement relations after differentiating Equation (1) are: * * ** ** ** y o yzyyz x o xzxxz xy o xyxy o xyxy y o yy o yy x o xx o xx zz y w z v zz x w z u zzz x v y u zzz y v zzz x u ϕ+γ+ϕ= ∂ ∂ + ∂ ∂ =γ ϕ+γ+ϕ= ∂ ∂ + ∂ ∂ =γ κ+γ+κ+γ= ∂ ∂ + ∂ ∂ =γ κ+ε+κ+ε= ∂ ∂ =ε κ+ε+κ+ε= ∂ ∂ =ε 2 2 32 32 32 (2) Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 5 Where y w x w xyyx x v y u y v x u o yy o xx yx xy y y x x ooo xy oo y oo x ∂ ∂ += ∂ ∂ += ∂ ∂ + ∂ ∂ = ∂ ∂ = ∂ ∂ = ∂ ∂ + ∂ ∂ = ∂ ∂ = ∂ ∂ = θϕ θϕ θθ κ θ κ θ κ γεε ,, ,, * o *o yz * o *o xz * o * o*o xy * o*o y * o*o x v u x v y u , y v , x u 2 2 =γ =γ ∂ ∂ + ∂ ∂ =γ ∂ ∂ =ε ∂ ∂ =ε (3) Also, all the strains above are defined in the middle-plane of the laminate. By substitution from Equation (3) into the stress-strain relations given by the following Equation: ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ = ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ yz xz xy y x yz xz xy y x QQ QQ QQQ QQQ QQQ γ γ γ ε ε τ τ τ σ σ 4445 4555 662616 262212 161211 000 000 00 00 00 (4) After complete integration, the stress-resultant/strain relations of the laminate are as follows(4): ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ κ κ κ κ κ κ γ ε ε γ ε ε ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ = ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ * * * * * * * * * * * * xy y x xy y x o xy o y o x o xy o y o x xy y x xy y x xy y x xy y x HHHFFFGGGEEE HHHFFFGGGEEE HHHFFFGGGEEE FFFDDDEEEBBB FFFDDDEEEBBB FFFDDDEEEBBB GGGEEEFFFDDD GGGEEEFFFDDD GGGEEEFFFDDD EEEBBBDDDAAA EEEBBBDDDAAA EEEBBBDDDAAA M M M M M M N N N N N N 662616662616662616662616 262212262212262212262212 161211161211161211161211 662616662616662616662616 262212262212262212262212 161211161211161211161212 662616662616662616662616 262212262212262212262212 161211161211161211161211 662616662616662616662616 262212262212262212262212 161211161211161211161211 (5) and, ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ ϕ ϕ γ γ ϕ ϕ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ = ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ * y * x *o yz *o xz y x * y * x y x y x FFEEDD FFEEDD EEDDBB EEDDBB DDBBAA DDBBAA Q Q S S Q Q 444544454445 455545554555 444544454445 455545554555 444544454445 455545554555 (6) all coefficients in A, B, D, E, F, G, and H groups are defined as follows: Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 6 ∑ = −−= NL L LLijij hhQA 1 1 )( i, j = 1, 2, 6 or i, j= 4, 5 (7 a) ( )∑ = −−= NL L LLijij hhQB 1 2 1 221 )( i, j = 1, 2, 6 or i, j= 4, 5 (7 b) ( )∑ = −−= NL L LLijij hhQD 1 3 1 331 )( i, j = 1, 2, 6 or i, j= 4, 5 (7 c) ( )∑ = −−= NL L LLijij hhQE 1 4 1 441 )( i, j = 1, 2, 6 or i, j= 4, 5 (7 d) ( )∑ = −−= NL L LLijij hhQF 1 5 1 551 )( i, j = 1, 2, 6 or i, j= 4, 5 (7 e) ( )∑ = −−= NL L LLijij hhQG 1 6 1 661 )( i, j = 1, 2, 6 (7 f) ( )∑ = −−= NL L LLijij hhQH 1 7 1 771 )( i, j = 1, 2, 6 (7 g) The present study explores the idea of tailoring the profile of reinforcing fibers to improve the buckling strength of composite plates. This study investigates the effect of waviness of fibers on the post buckling curves, as shown in Figure (4), and this waviness is of the form: ( ) ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ π α= a xk xy sin (8) such that the angle of fiber orientation θ varies along the longitudinal x-axis as: ( ) ( )xkk a xk a k dx dy ππ∆=⎟ ⎠ ⎞ ⎜ ⎝ ⎛ ππα ==θ cos.cos.tan (9) where a = plate length; k= number of half sine waves; and α= wave amplitude. Two normalized variables, ∆=α/a and axx /= , are introduced. The main objective is to study the effect of fiber waviness, characterized by k and ∆, on the static and dynamic buckling behavior of composite laminates. The fiber can also be rotated in any direction with the x-axis, as shown in Figure (4), by using the following expression: ( ) ( )β+β= sincos yxxn (10) where xn represent the x-coordinate for a rotated fiber, and β is the angle of the waviness fiber. The angle of fiber orientation in Equation (9) is variable with x-coordinate and instead of the constant angle used for straight fibers. Figure (5) shows the principal material directions aligned with the lamina axes by angle (β ). Geometrical Nonlinearity The components of the Green-Lagrangian strain vector are known in terms of local derivatives of the displacements for the plate element as [Pica, et al., 1979](13): [ ] ⎭ ⎬ ⎫ ⎩ ⎨ ⎧ε + ⎭ ⎬ ⎫ ⎩ ⎨ ⎧ε + ⎭ ⎬ ⎫ ⎩ ⎨ ⎧ ε ε + ⎭ ⎬ ⎫ ⎩ ⎨ ⎧ε = ⎪ ⎪ ⎪ ⎭ ⎪⎪ ⎪ ⎬ ⎫ ⎪ ⎪ ⎪ ⎩ ⎪⎪ ⎪ ⎨ ⎧ γ γ γ ε ε =ε 000 p I P L s o b o p o yz xz xy y x z (11) where the linear mid-plane strains are: Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 7 [ ] ⎪ ⎪ ⎪ ⎭ ⎪ ⎪ ⎪ ⎬ ⎫ ⎪ ⎪ ⎪ ⎩ ⎪ ⎪ ⎪ ⎨ ⎧ ∂ ∂ + ∂ ∂ ∂ ∂ ∂ ∂ =ε x v y u y v x u p o (12) Equation (12) represents the in-plane strain. Also: [ ] ⎪ ⎪ ⎪ ⎭ ⎪ ⎪ ⎪ ⎬ ⎫ ⎪ ⎪ ⎪ ⎩ ⎪ ⎪ ⎪ ⎨ ⎧ ∂ θ∂ + ∂ θ∂ ∂ θ∂ ∂ θ∂ =ε xy y x yx y x b o (13) Equation (13) represents the bending strain. Also: [ ] ⎭ ⎬ ⎫ ⎩ ⎨ ⎧ ϕ ϕ =ε y xs o (14) Equation (14) represents the shear strain. Moreover: [ ] ⎪ ⎪ ⎪ ⎪ ⎭ ⎪⎪ ⎪ ⎪ ⎬ ⎫ ⎪ ⎪ ⎪ ⎪ ⎩ ⎪⎪ ⎪ ⎪ ⎨ ⎧ ∂ ∂ ∂ ∂ ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ ∂ ∂ ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ ∂ ∂ =ε y w x w y w x w p L 2 2 2 1 2 1 (15) Equation (15) represents the nonlinear component of in-plane strain. Finally: [ ] ⎪ ⎪ ⎪ ⎭ ⎪ ⎪ ⎪ ⎬ ⎫ ⎪ ⎪ ⎪ ⎩ ⎪ ⎪ ⎪ ⎨ ⎧ ∂ ∂ ∂ ∂ + ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ =ε y w x w y w x w y w y w x w x w oo o o p I (16) Equation (16) represents the initial strain due to initial deflection. The vector components of Equation (11) represent the generalized strains. It can be noted that the vector ( pL p L o p ε+ε+ε ) reproduce the Marguerre strain expression for plate. Variational Equation of Equilibrium The variation of strain dε due to the virtual displacements du, generally dε is given as the sum of the variation of the linear and nonlinear generalized strains as: Lo ddd ε+ε=ε (17) Since oε is a linear function of displacement, [ ]duBd oo =ε (18) Also, [ ]duBd LL =ε (19) Thus, the total strain-nodal displacement matrix for total strains is: Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 8 [ ] [ ] [ ]Lo BBB += (20) in which [ ]oB is the same matrix as in the linear strain analysis and [ ]LB depends on the displacements. [ ] [ ] ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ =⎥ ⎦ ⎤ ⎢ ⎣ ⎡ = 00 0 0 0 bL Lb o p o o B B B B B , (21) where [ ]LB can be found by taking the variation of the nonlinear strain components { }pLε with respect to the displacements. This nonlinear strain components of Equation (11) can be written in a more convenient form as: [ ] [ ]{ }θ= ⎪ ⎪ ⎭ ⎪⎪ ⎬ ⎫ ⎪ ⎪ ⎩ ⎪⎪ ⎨ ⎧ ∂ ∂ ∂ ∂ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ = ⎪ ⎪ ⎪ ⎪ ⎭ ⎪⎪ ⎪ ⎪ ⎬ ⎫ ⎪ ⎪ ⎪ ⎪ ⎩ ⎪⎪ ⎪ ⎪ ⎨ ⎧ ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ ∂ ∂ ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ ∂ ∂ ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ ∂ ∂ ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ ∂ ∂ =ε θA y w x w x w y w y w x w y w x w y w x w p L 2 1 0 0 2 1 2 1 2 1 2 2 (22) where the displacement gradients with respect to the lateral displacements (w) are: { } ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ ∂ ∂ ∂ ∂ =θ y w x w (23) Then, the variation of the nonlinear component of the in-plane strain is obtained from Equation (22) in terms of the virtual gradients dθ as: θ=ε θ dAd p L (24) where ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ =θ x w y w y w x w A 0 0 (25) in which, Aθ represents the gradients of total displacements and, ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ ∂ ∂ ∂ ∂ =θ y w x w dd (26) represents the gradients of incremental displacements. The displacement gradients (θ) of Equation (23) may now be written in terms of the nodal displacements (u) and Cartesian derivatives of the shape functions as: { } [ ]{ }uG=θ (27) where [ ] [ ]nGGGG .......21= (28) and, Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 9 ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ ∂ ∂ ∂ ∂ = 0000 0000 y N x N G i i i (29) The above equation represents the gradient for five degrees of freedom per node. Taking the variation of Equation (22) as follows: { } [ ]{ } [ ] { } [ ] { } [ ][ ] { }udGAdAdAAdd pL θθθθ =θ=θ+θ=ε 2 1 2 1 (30) hence immediately, by definition [ ] [ ][ ]GAB obL θ= (31) In order to incorporate the imperfections in the formulation, the strain due to imperfections as given in Equation (11), thus can be written the imperfection strain components as follows: { } [ ]{ }θ= ⎪ ⎪ ⎭ ⎪⎪ ⎬ ⎫ ⎪ ⎪ ⎩ ⎪⎪ ⎨ ⎧ ∂ ∂ ∂ ∂ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ =ε 12 1 22 20 02 2 1 A y w x w x w y w y w x w oo o o I p (32) Following the same analysis as for the nonlinear strains, the combination of the effects of the two strains into one [ ]bLB is defined as: [ ] [ ][ ]GAB bL 2= (33) where [ ] [ ] [ ]12 AAA += θ [ ] ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ ∂ ∂ + ∂ ∂ ∂ ∂ + ∂ ∂ ∂ ∂ + ∂ ∂ ∂ ∂ + ∂ ∂ = x w x w y w y w y w y w x w x w A oo o o 22 20 02 2 (34) In the present study, the imperfection is assumed to be of sinusoidal function over the plate as: ( ) ⎟⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ π ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ π = yx oo L yn L xn wyxw sinsin, (35) where Lx, and Ly are the dimensions of the plate in the x, and y-direction, respectively. wo is the maximum value of the initial imperfection at the plate center. The variation in the potential energy of deformation for a plate element with large deflection can be written as: [ ] dVBdU T V σ= ∫ (36) Substituting Equations (36) can give the equilibrium equations written as: ( ) [ ] 0=−σ=Ψ ∫ dWdVBu V T (37) where Ψ represents the sum of external and internal generalized forces. Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 10 Clearly, the solution of Equation (37) will have to be approached iteratively. In order to use an incremental solution procedure, the relation between du and dΨ must be found. Thus, taking appropriate variation of Equation (37) with respect to du: [ ] [ ] dV du d BdV du Bd du d T V T V σ σ ∫∫ += Ψ (38) where the variation of the external load with respect to displacements is equal to zero, and thus Equation (38) can be written at another form: [ ] [ ]KdV du Bd du d T V L += ∫ σ Ψ (39) where [ ] [ ] [ ][ ] [ ] [ ]Lo T A KKdABDBK +== ∫ (40) The first term of Equation (39) can be written as: [ ] [ ]σσ KdVdu Bd T V L =∫ (41) where [ ]σK is a symmetric matrix dependent on the stress level. This matrix is known as initial stress matrix or geometric matrix. Thus, [ ] [ ] [ ]( ) [ ] duKduKKKd TLo =++=Ψ σ (42) with [ ]TK being the total, or tangent stiffness matrix. Tangent Stiffness Matrix The tangent stiffness matrix can be written as: [ ] [ ] [ ] [ ]σ++= KKKK LoT (43) where [ ]oK is the constant linear elastic stiffness matrix and can be written as: [ ] [ ] [ ][ ]dABDBK o T A oo ∫= (44) [ ]LK is the initial or large displacement matrix which is quadratically dependent upon displacement u, and can be written as: [ ] [ ] [ ][ ] [ ] [ ][ ] [ ] [ ][ ]dABDBdABDBdABDBK A o A T LL A T LL T oL ∫ ∫∫ ++= (45) Finally [ ]σK is the initial stress stiffness matrix which has to be found by using the definition of Equation (41). By taking the variation of Equation (20) then: [ ] [ ] ⎥⎦ ⎤ ⎢ ⎣ ⎡ = 0 00 Tb L T L Bd Bd (46) This on substitution into Equations (41) and (34) gives: [ ] [ ] [ ]∫ ⎪ ⎪ ⎪ ⎪ ⎭ ⎪⎪ ⎪ ⎪ ⎬ ⎫ ⎪ ⎪ ⎪ ⎪ ⎩ ⎪⎪ ⎪ ⎪ ⎨ ⎧ ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ =σ A xy y x xy y x TT M M M N N N AdG K 0 00 (47) Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 11 However, using the mathematical properties of the matrix [A], this matrix can be written as: [ ] [ ]daG NN NN N N N Ad yxy xyx xy y x T ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ = ⎪ ⎭ ⎪ ⎬ ⎫ ⎪ ⎩ ⎪ ⎨ ⎧ (48) and finally one can obtain [ ] [ ]⎥⎦ ⎤ ⎢ ⎣ ⎡ = σ σ bK K 0 00 (49) Thus, [ ] [ ] [ ]daG NN NN GK A yxy xyxT ∫ ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ =σ (50) Failure criteria for laminated plate structure The stresses in an individual lamina are fundamental to control the failure initiation and progression in the laminate. The strength of each individual lamina is assessed separately by considering the stresses acting on it along the material axes. The initial failure of a lamina is governed by exceeding the maximum limit prescribed by a failure criterion. The determination of failure load is very essential in understanding the failure process as well as the reliability and safety of structures. The ultimate load that makes the plate fail is calculated based on Tsai-Wu criterion for general composite materials and on Hashin criterion for fiber composite materials as follows [Jones, 1999](9): 611 ,......,,; ==σσ+σ jiFF jiijii (51) where Fi and Fij are strength tensors of the second and fourth order respectively and the usual contracted tensor notation is used except that 235134 τστσ == , ,and 126 τσ = . Equation (51) is obviously very complicated thus it will restrict the above attention to the reduction of above equation for an orthotropic lamina under plane stress conditions: 12666 2 555 2 44432233113 2112 2 333 2 222 2 111332211 =σ+σ+σ+σσ+σσ+ σσ+σ+σ++σ+σ+σ+σ FFFFF FFFFFFF (52) The terms that are linear in the stresses are useful in representing different strengths in tension and in compression. The terms that are quadratic in the stresses are the more or less usual terms to represent an ellipsoid in stress space, where 03 =F indicates that to the shear strength of a material in compression and in tension is similar, and 03 =σ in z-direction. Also, the shear strength of a material is equal in three dimensions and equal to S. Thus, the terms of Fi is: 50 23 50 13 50 12 266255244 332211 321 1 2 1 1 2 11 2 1 111 111 111111 . .. ... ... , .. ,, . , . , . ,, ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛− = ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛− =⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛− = ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ =⎟ ⎠ ⎞ ⎜ ⎝ ⎛ =⎟ ⎠ ⎞ ⎜ ⎝ ⎛ = ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ =⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ =⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ = ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ −=⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ −=⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ −= tcct tcctctct ctctct ctctct ZZYY F ZZXX F YYXX F T F S F R F ZZ F YY F XX F ZZ F YY F XX F (53) Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 12 where, =tc XX , The axial or longitudinal strength in compression and tension. =tc YY , The transverse strength in compression and tension. =tc ZZ , The transverse strength in compression and tension. =STR ,, Shear strength of the material. Equation (52) becomes as: 12666 2 555 2 4442112 2 222 2 1112211 =σ+σ+σ+σσ+σ++σ+σ+σ FFFFFFFF (54) Equation (54) is suitable for the elastic-plastic analysis of anisotropic materials. For matrix cracking failure, two different failure criteria are used depending on whether the transverse normal stress, 22σ , is in tension or in compression. The failure index, 2 me , is defined in terms of transverse tensile strength, Yt , transverse compressive strength, Yc, and in-plane shear strength, R, and is expressed as: 2 12 2 22 2 222 2 1 2 ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ +⎟ ⎠ ⎞ ⎜ ⎝ ⎛ + ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ −⎟ ⎠ ⎞ ⎜ ⎝ ⎛ = RRR Y Y e c c m τσσ for 022 <σ (55) and, 2 12 2 222 ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ +⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ = RY e t m τσ for 022 >σ (56) where (em) is the failure index for matrix cracking. Matrix cracking is assumed to occur when the failure index (em) exceeds unity. Fiber-matrix shear failure is assumed to be dependent on a combination of axial stress, 11σ , and shear stress,τ12, and is expressed as follows: 2 12 2 112 ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ +⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ = RX e t s τσ for 011 >σ and 2 12 2 11 ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ <⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ RX t τσ (57) and, 2 12 2 112 ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ +⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ = RX e c s τσ for 011 <σ and 2 12 2 11 ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ <⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ RX t τσ (58) where (es) is the failure index for fiber-matrix shearing, Xt is the tensile strength along the fiber direction and Xc is the compressive strength along the fiber direction. Equations (57) and (58) predict that when the failure (es) exceeds unity, fiber-matrix shearing dominated failure occurs. Fiber breakage failure occurs in tension due to the combination of axial stress and shear stress while the failure in compression is governed by buckling as expressed in terms of only axial stress. The criterion for breakage failure is expressed as follows: 2 12 2 112 ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ +⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ = RX e t f τσ for 011 >σ (59) and, 22 112 ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ = c f X e σ for 011 <σ (60) where (es) is the failure index for fiber breakage. The fiber breakage failure occurs when (es) exceeds unity. Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 13 Applications and Discussions In order to verify the reliability of the adopted numerical method, some case studies reported by other researchers are utilized. The study of the composite plate will be introduced. Comparison with experimental investigation of composite plate The post-buckling and failure characteristics of flat, rectangular graphite-epoxy panels with and without holes that are loaded in axial compression have been examined in an experimental study by Starnes and Rouse [1981] and in a theoretical study by Elseifi [1998](7). The panels were fabricated from commercially available unidirectional Thornel 300 graphite-fiber tapes preimpregnated with 450 K cure Narmco 5208 thermosetting epoxy resin. Typical lamina properties for this graphite-epoxy system are 131.0 GPa for the longitudinal Young’s modulus, 13.0 GPa for the transverse Young’s modulus, 6.4 GPa for the in-plane shear modulus, 0.38 for the major Poisson’s ratio (v12 ), and 0.14 mm for the lamina thickness. The loaded ends of the panel were clamped by fixtures during testing, and the unloaded edges were simply supported by knife-edge restrains to prevent the panels from buckling. In the present study, the plate is analyzed by using nine-node isoparametric Lagrangian finite elements with nine degrees of freedom per node. The panel is 0.508 m long, 0.178 m wide, and 24-ply orthotropic laminate with [(45o/0o /-45o)2, (45o /0o /-45o)2, (45o /0o /90o)4] stacking sequence. The modeling approach of the quarter plate was based on using two elements in the short direction, and three elements in the long direction. The finite element mesh used is shown in Figure (7). In order to efficiently proceed beyond the critical buckling point in the post-buckling analysis of the panel, an initial geometric imperfection in the same shape as the first buckling mode was assumed. The amplitude of this initial imperfection is (1%) of the total laminate thickness. Figure (8) shows the out-of-plane deflection (w) near a point of maximum deflection (node i ) normalized by the panel thickness h as a function of the normalized load. From this figure, it can be noticed that good agreement exists with the experimental results with a difference not more than (6%). On the other hand, the present results are closer to the experimental investigation. Comparison with theoretical investigation of composite plate A square cross-ply laminated plate with simply supported edges and initial imperfection was analyzed by Zou and Qiao [2002](16). Lateral in-plane expansion is allowed at the loaded ends and the unloaded edges can be moved in the plane but remain straight. The layer material and geometry properties are presented in Figure (9). The slenderness ratio is set as (b/h=20), and it represents a moderately thick laminate. The laminated plate contains eight equal-thickness layers in [0o/90o]4 layup. The initial imperfection (wo/h) is given by (0.0 and 0.1) by which the shape is considered to be a sinusoidal curve. Zou and Qiao used higher order finite strip method and solved the nonlinear equations by Newton-Raphson method. In the present study, a quarter of the laminate is modeled with (2×2) mesh of nine-node isoparametric Lagrangian element with nine degrees of freedom per node. Numerical results and response comparisons with Zou and Qiao [2002](16) are shown in Figures (10) and (11) for axial load versus total deflection and axial load versus end shortening strain (ε). The present results are really close to those of Zou and Qiao [2002](16) with a difference of not more than (15%). Parametric Study A parametric study is performed to assess the influence of several important parameters on the elastic-plastic large displacement analysis of a composite laminated plate subjected to in-plane compressive load. The selected parametric studies are summarized as follows: 1. The effect of through-thickness shear deformation. 2. The effect of fiber’s orientation angle. 3. The effect of degree of orthotropy of individual layers. 4. The effect of fiber waviness. Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 14 Each one of the above parameters was studied individually by analyzing a type of laminated composite plate. In all cases, a nine-node element was used and also one quadrant of the plate was analyzed due to symmetry and (2×2) mesh is used in the cross-ply and straight fiber plates while for angle-ply and sine wave fiber plates, were analyzed by considering the full plates with (4×4) element mesh. Lateral in-plane expansion is allowed at the loaded ends and the unloaded edges can be moved in the in-plane direction but remain straight. The initial imperfection shape is considered to be a sinusoidal curve. The following geometry and layer material properties of high graphite epoxy are used in the analysis: (E1=172.5 GPa; E2=7.08 GPa; G12=G13=3.45 GPa, G23=1.38 GPa; Ef=341.42 GPa; Em=3.58 GPa; Vf=0.5; Vm=0.5; v12=v13=v23=0.25, Xt=Xc=1450 MPa, Yt=36 MPa, Yc=230 MPa, S=62 MPa)[ Parhi, et al., 2001](83). The geometry properties are (a=1.0 m, a/b=1). 1. Effect of through-thickness shear deformation To show the effect of transverse shear deformation on the large displacement elastic-plastic analysis of laminated composite plate under in-plane compressive load, a simply supported plate with a range of slenderness ratio (b/h) from (20) to (120), with symmetric cross-ply and antisymmetric cross-ply arrangement and with six layers was analyzed. The initial imperfection is (wo/h= 0.1) by which the shape is considered to be a sinusoidal curve. Figures (12) and (13) present the load-deflection curves of the symmetric cross-ply, and the antisymmetric cross-ply laminated composite plate under in-plane loading and with slenderness ratio (b/h=20) by taking the effects of transverse shear deformation through the degrees of freedom per node of the element. Figures (14) and (15) show the effect of shear deformation of symmetric cross-ply, and antisymmetric cross-ply laminated composite plate under in-plane loading with range of slenderness ratio(b/h) (20-120). 2. Effect of fiber’s orientation angle To study the effect of fiber’s orientation angle on the large displacement elastic-plastic analysis of laminated composite plates under in-plane compressive load, a square simply supported laminated plate with two layers was analyzed. The initial imperfection is neglected in the present study. Figure (16) shows the effect of fiber’s orientation on the nonlinear analysis of composite laminated plate under in-plane compressive load. From this figure, it could be noticed that the ultimate strength of the plate with (0o/90o) gives ultimate load (651.2 kN/m). This orientation’s fiber means that it is the optimum for a plate under in-plane compressive load. 3. Effects of fiber waviness To show the effect of fiber waviness on the large displacement elastic-plastic analysis of laminated composite plate, a square simply supported plate, with six layers was considered. The shape of fiber was considered to follow a sinusoidal curve. The effects of this type are (number of sequence (k), amplitude of wave (∆), and fiber’s orientation. The initial imperfection (wo/h= 0.1) by which the shape is considered to be sinusoidal curve. The value of amplitude of sine wave fiber is varying (∆=0.05-0.5) and the number of sequences of the sine wave fiber (k) was considered changeable in the range of (1-12). Figures (17) and (18) present the load-deflection curves of laminated with symmetric cross- ply composite plate under in-plane compressive load and with sine wave fibers with a range of amplitude (0.05-0.5). Figures (19)-(22) present the ultimate strength-fiber path amplitude (∆) curves for the laminated, with symmetric cross-ply, and antisymmetric cross-ply composite plates under in-plane compressive load and with a range of number of sequences (k) (1-12), respectively. Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 15 Conclusions A nonlinear finite element method is adopted for the large displacement elastic-plastic static analysis of anisotropic plates under in-plane compressive load. Many of effects were considered in the present study such as the effect of orthotropy of individual layers, through-thickness shear deformation, fiber’s orientation angle, and fiber waviness on the large displacement elastic-plastic static analysis. The following conclusions are drawn with regard to the results obtained for the anisotropic plates under in-plane static loading as follows, 1) The capacity of a laminated plate with sine wave fiber under in-plane compressive load in the direction of waviness is more than the capacity of the plate under in-plane compressive load orthogonal to the direction of waviness by (42%) for a plate with sine wave fiber (k=12, ∆=0.4). 2) The ultimate strength of a cross-ply laminated plate is greater than the ultimate strength of an angle-ply laminated plate with the same laminas. 3) The ultimate strength of a symmetric cross-ply laminated plate is greater than the ultimate strength of an antisymmetric cross-ply laminated plate with the same laminas. 4) A symmetric cross-ply laminated plate with sine wave fiber (k=12, ∆=0.4) gives in-plane loading (610 kN/m) and this represents the peak capacity and this case is the best. References 1. Abid-Ali, A. K., “Stress Analysis of Laminated Fiber Reinforced Composite Cylinder”, M.Sc. Thesis, University of Babylon, Hilla, Iraq, 2000. 2. Ali, N. H., “Finite Element Dynamic Analysis of Laminated Composite Plates Including Damping Effect”, M.Sc. Thesis, University of Babylon, Hilla, Iraq, 2004. 3. Amash, H. K., “Post-buckling and Post-Yielding Analysis of Imperfect Thin Plate by Finite Difference Method”, M.Sc. Thesis, University of Babylon, Hilla, Iraq, 2003. 4. Amash, H. K., “Nonlinear Static and Dynamic Analysis of Laminated Plates Under In-plane Forces”, Ph.D. Thesis, University of Babylon, Hilla, Iraq, 2008. 5. Chia, C.Y., “Nonlinear Analysis of Plates”, McGraw-Hill International Book Company, 1980. 6. Cook, R.D., “Finite Element Modeling for Stress Analysis”, John Wiley & Sons, Inc., 1995. 7. Elseifi, M., “A New Scheme for the Optimum Design of Stiffened Composite Panels with Geometric Imperfections”, Ph.D. Thesis, Virginia Polytechnic Institute, U.S.A., 1998. 8. Hashin, Z. “Failure Criteria for Unidirectional Fiber Composites”, ASME, J. Appl. Mech., Vol.47, 1980, pp.329-334. 9. Jones, R.M., “Mechanics of Composite Materials”, Second Edition, Taylor and Francis Inc., U.S.A., 1999. 10. Kaw, A., “Mechanics of Composite Materials”, Second Edition, Taylor and Francis Group, LLC, 2006. 11. Kommineni, J. R., and Kant, T. “Geometrically Non-linear Transient Co Finite Element Analysis of Composite and Sandwich Plates with a Refined Theory.” Struct. Eng. And Mech., Vol.1, No.1, 1993, pp87-102. 12. Pandey, M. D., “Effect of Fiber Waviness on Buckling Strength of Composite Plates.” ASCE, J. Eng. Mech., Vol.125, No.10, 1999, pp.1173-1179. 13. Pica, A., Wood, R.D., and Hinton, E. “Finite Element Analysis of Geometrically Nonlinear Plate Behavior Using a Mindlin Formulation.” Comp. & Struct., Vol.11, 1979, pp.203-215. 14. Pica, A., Wood, R.D “Post-Buckling Behavior of Plates and Shells Using a Mindlin Shallow Shell Formulation.” Comp. & Struct., Vol.12, 1980, pp.759-768. 15. Shukla, K.K., and Nath, Y., “Thermomechanical Post-Buckling of Cross-Ply Laminated Rectangular Plates”, ASCE, J. Eng. Mech., Vol.128, No.1, 2002, pp.93-101. 16. Zou, G., and Qiao, P., “Higher Order Finite Strip Method for Post-Buckling Analysis of Imperfect Composite Plates" J. Eng. Mech., Vol.128, No.9, Sep., 2002, pp.1008-1015. Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 16 Max. Deflection ((w+wo)/h) Present study Theoretical results(31) Experimental results(*) Load (Px/Pcr) 0.010 0.010 0.010 0.00 0.013 0.020 0.021 0.10 0.017 0.030 0.042 0.30 0.035 0.060 0.063 0.60 0.064 0.110 0.125 0.75 0.210 0.280 0.310 0.90 0.770 0.645 0.850 1.10 1.122 1.021 1.176 1.28 1.495 1.380 1.580 1.50 1.645 1.500 1.720 1.60 1.930 1.720 2.020 1.80 2.068 1.810 2.166 1.90 2.203 1.900 2.333 2.00 * : Given by Starnes and Rouse [1981] as mentioned in Elseifi [1998](7). Figure (1): Use of composites in the space shuttle [Kaw,2006](10) Table (1): Comparison of results with experimental and theoretical studies of composite laminated plate under in-plane compressive load in x-direction, (a/b=2.854, wo/h=0.1) Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 17 Figure (2): Laminated plate with several lamina orientations[Jones, 1999](9) LAMINAS LAMINATE RESPONSE? CONDITIONS? h2 h1 h0 hn-1 hn-2 hn-3 hn h Middle surface n Layer number 1 2 Figure (3): Geometry of an NL-layered laminate [Jones,1999](9) Fiber (4): (a) Lamina with variable fiber orientation; (b) Geometry of sinusoidal fiber path [Pandey,1999](12) Lamina Fibers y x b a α k=2 Sinusoidal fiber path. (a) (b) Figure (5): Laminate plate with sine wave fibers aligned with x-axis β=0 β Figure (6): Idealizated stress-strain relationship of uniaxial loading behavior for orthrotropic plate[Jones, 1999](9) Ultimate point uε ε fσ σ Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 18 0.0 0.5 1.0 1.5 2.0 2.5 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 Experiment by Starnes Elseifi Present study (30) In -p la ne lo ad (P x/ P cr ) Deflection ratio ((w+wo)/h) Figure (8): Post-buckling curve for a rectangular thin composite laminated plate under in-plane compressive load in x- direction Px * (31) Figure (7): Finite element mesh and boundary condition for the quarter of the composite plate under in-plane compressive load in x-direction x u=0 v=0 θx=0 θy≠0 u=0 v=0 θx=0 θy≠0 u≠0;v=0 θy=0;θx=0 u=0;v=0 θx=0;θy=0 y 89 mm 25 4 m m i Px Px Figure (9): Details of a square laminated composite plate under in-plane compressive load and material properties a=200 (mm) b=200 (mm) h=10 mm E1/E2=20 G12/E2=G13/E2=0.6 G23/E2=0.5 E2=10 GPa v12=v13=v23=0.25 Px x y Figure (10): Load-deflection curve of antisymmetric cross-ply imperfect rectangular thick composite laminated plate under compressive load,(b/h=20, a/b=1) End shortening strain ( ) Figure (11): Load-end shortening strain curve of antisymmetric cross-ply imperfect rectangular thick composite laminated plate under compressive load,(b/h=20, a/b=1) Px 0.00 0.25 0.50 0.75 1.00 1.25 0 500 1000 1500 2000 2500 3000 3500 4000 4500 5000 5500 6000 6500 7000 7500 Zou and Qiao Present study (124) In -p la ne lo ad (k N /m ) Deflection ratio ((w+wo)/h) Px 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 0 500 1000 1500 2000 2500 3000 3500 4000 4500 5000 5500 6000 6500 7000 7500 Zou and Qiao Present study (124) In -p la ne lo ad (k N /m ) (16) (16) Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 19 PxPx Figure (12): Load-deflection curve of symmetric cross-ply imperfect square composite laminated plate under in-plane compressive load with (b/h=20;a/b=1;wo /h=0.1) Figure (13): Load-deflection curve of antisymmetric cross-ply imperfect square composite laminated plate under in-plane compressive load with (b/h=20;a/b=1; wo /h=0.1) Deflection ratio ((w+wo)/h) Deflection ratio ((w+wo)/h) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 0 2500 5000 7500 10000 12500 15000 17500 Antisymmetric cross-ply 5 DOF (no shear effect) 7 DOF (first order shear effect) 9 DOF (higher order shear effect) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0 2500 5000 7500 10000 12500 15000 17500 20000 Symmetric cross-ply 5 DOF (no shear effect) 7 DOF (first order shear effect) 9 DOF (higher order shear order) In -p la ne lo ad (k N /m ) In -p la ne lo ad (k N /m ) 5 DOF (first order shear effect) 7 DOF (higher order shear effect) 9 DOF (higher order shear effect) 5 DOF (first order shear effect) 7 DOF (higher order shear effect) 9 DOF (higher order shear effect) Deflection ratio ((w+wo)/h) Figure (16): Effect of orientation’s fiber on the large displacement elastic-plastic analysis of composite laminated plate under in-plane compressive load (b/h=100; wo /h=0.0;a/b=1) 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 0 50 100 150 200 250 300 350 400 450 500 550 600 650 700 Orientation of layers 0/90 30/-30 45/-45 60/-60 75/-75 Px In -p la ne lo ad (k N /m ) Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 20 Slenderness ratio ((b/h) Figure (14): Effect of number of degrees of freedom on the large elastic-plastic analysis of symmetric cross-ply composite laminated plate under in-plane compressive load with a range of slenderness ratio (b/h) Slenderness ratio ((b/h) Figure (15): Effect of number of degrees of freedom on the large elastic-plastic analysis of antisymmetric cross-ply composite laminated plate under in-plane compressive load with a range of slenderness ratio (b/h) U lt im at e lo ad (k N /m ) U lt im at e lo ad (k N /m ) 20 40 60 80 100 120 0 2500 5000 7500 10000 12500 15000 17500 20000 20 40 60 80 100 120 0 2500 5000 7500 10000 12500 15000 17500 5 DOF (no shear effect) 7 DOF (first order shear effect) 9 DOF (higher order shear effect) 5 DOF (first order shear effect) 7 DOF (higher order shear effect) 9 DOF (higher order shear effect) 5 DOF (first order shear effect) 7 DOF (higher order shear effect) 9 DOF (higher order shear effect) Figure (17): Load-deflection curve of simply supported square laminated composite plate with sine wave fiber with a range of fiber amplitude (0.05-0.50) and under in-plane compressive load, (wo/h =0.1;b/h=100;a/b=1.0;k=4) Deflection ratio ((w+wo)/h) Px 0 1 2 3 0 50 100 150 200 250 300 350 400 450 500 550 600 Symmetric cross-ply 0/90/0/0/90/0 k=4, D=0.05 k=4, D=0.10 k=4, D=0.20 k=4, D=0.30 k=4, D=0.40 k=4, D=0.50 In -p la ne lo ad (k N /m ) k=4, ∆=0.20 k=4, ∆=0.30 k=4, ∆=0.40 Symmetric cross-ply 0/90/0/0/90/0 k=4, ∆=0.50 k=4, ∆=0.05 k=4, ∆=0.10 Laminated plate (0)6 k=4, ∆=0.05 k=4, ∆=0.10 k=4, ∆=0.20 k=4, ∆=0.30 k=4, ∆=0.40 k=4, ∆=0.50 Figure (18): Load-deflection curve of simply supported square symmetric cross-ply composite plate with sine wave fiber with a range of fiber amplitude (0.05-0.50) and under in-plane compressive load, (wo/h =0.1;b/h=100;a/b=1.0;k=4) Deflection ratio ((w+wo)/h) Px 0 1 2 3 4 5 0 100 200 300 400 500 600 Laminated plate k=4,D=0.05 k=4,D=0.10 k=4,D=0.20 k=4,D=0.30 k=4,D=0.40 k=4,D=0.50 In -p la ne lo ad (k N /m ) Al-Qadisiya Journal For Engineering Sciences Vol. 2 No. 1 Year 2009 21 0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40 0.45 0.50 250 300 350 400 450 500 550 Laminated plate k=1 k=2 k=4 k=8 k=12 Fiber path amplitude (∆) Figure (19): Ultimate load-fiber path amplitude curve of simply supported square laminated composite plate under in-plane compressive load in x-direction and with a range of number of sequences (1-12),(wo/h= 0.1, b/h=100, a/b=1.0) k=1 k=2 k=4 k=8 k=12 Laminated plate (0)6 U lt im at e lo ad (k N /m ) Px Fiber path amplitude (∆) U lt im at e lo ad (k N /m ) Figure (20): Ultimate load-fiber path amplitude curve of simply supported square laminated composite plate under in-plane compressive load in y-direction and with a range of number of sequences (1-12),(wo/h=0.1, b/h=100, a/b=1.0) Py 0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40 0.45 0.50 300 350 400 450 500 550 600 650 Symmetric cross-ply k=1 k=2 k=4 k=8 k=12 Laminated plate (0)6 k=1 k=2 k=4 k=8 k=12 0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40 0.45 0.50 250 300 350 400 450 500 550 600 650 Laminated plate k=1 k=2 k=4 k=8 k=12 U lt im at e lo ad (k N /m ) Fiber path amplitude (∆) Figure (21): Ultimate load-fiber path amplitude curve of simply supported square symmetric cross-ply composite plate under in-plane compressive load and with a range of number of sequences (1-12),(wo/h= 0.1, b/h=100, a/b=1.0) k=1 k=2 k=4 k=8 k=12 Symmetric cross-ply 0/90/0/0/90/0 Px Fiber path amplitude (∆) U lt im at e lo ad (k N /m ) Figure (22): Ultimate load-fiber path amplitude curve of simply supported square antisymmetric cross-ply composite plate under in-plane compressive load and with a range of number of sequences (1-12),(wo/h= 0.1, b/h=100, a/b=1.0) Px 0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40 0.45 0.50 300 350 400 450 500 550 600 650 Antisymmetric cross-ply k=1 k=2 k=4 k=8 k=12 k=1 k=2 k=4 k=8 k=12 Antisymmetric cross-ply 0/90/0/90/0/90