Copyrighted Materials 0p UyPress rm CHAPTER NINE Finite Element Analysis The finite element method offers a practical means of calculating the deformations of,and stresses and strains in,complex structures.A detailed description of the finite element method is beyond the scope of this book.Instead,we focus on those features specific to composite materials. Finite element analysis consists of the following major steps: 1.A mesh encompassing the structure is generated (Fig.9.1). 2.The stiffness matrix [k]of each element is determined. 3.The stiffness matrix [K]of the structure is determined by assembling the ele- ment stiffness matrices. 4.The loads applied to the structure are replaced by an equivalent force system such that the forces act at the nodal points. 5.The displacements of the nodal points d are calculated by k]d f. (9.1) where f is the force vector representing the equivalent applied nodal forces (Fig.9.1). 6.The vector d is subdivided into subvectors 6,each 6 representing the displace- ments of the nodal points of a particular element. 7.The displacements at a point inside the element are calculated by u=[W6, (9.2) where the vector u represents the displacements and [N]is the matrix of the shape vectors. 8.The strains at a point inside the elements are calculated by e=[B6, (9.3) where [B]is the strain-displacement matrix. 395
CHAPTER NINE Finite Element Analysis The finite element method offers a practical means of calculating the deformations of, and stresses and strains in, complex structures. A detailed description of the finite element method is beyond the scope of this book. Instead, we focus on those features specific to composite materials. Finite element analysis consists of the following major steps: 1. A mesh encompassing the structure is generated (Fig. 9.1). 2. The stiffness matrix [k] of each element is determined. 3. The stiffness matrix [K] of the structure is determined by assembling the element stiffness matrices. 4. The loads applied to the structure are replaced by an equivalent force system such that the forces act at the nodal points. 5. The displacements of the nodal points d are calculated by [K] d = f, (9.1) where f is the force vector representing the equivalent applied nodal forces (Fig. 9.1). 6. The vector d is subdivided into subvectors δ, each δ representing the displacements of the nodal points of a particular element. 7. The displacements at a point inside the element are calculated by u = [N] δ, (9.2) where the vector u represents the displacements and [N] is the matrix of the shape vectors. 8. The strains at a point inside the elements are calculated by ε = [B] δ, (9.3) where [B] is the strain–displacement matrix. 395
396 FINITE ELEMENT ANALYSIS nodal forces nodal points elements Figure 9.1:Structure and its finite element mesh. 9.The stresses at a point inside the element are calculated by =[E]E, (9.4) where [E]is the stiffness matrix characterizing the material. The element stiffness matrix,referred to in Step 2,is defined as [k6=e, (9.5) where fe represents the forces acting at the nodal points of the element.The element stiffness matrix is1 内=[A'A[Ba业 (9.6 where V is the volume of the element. The preceding steps apply to structures made of either isotropic or compos- ite materials.The only difference between isotropic and composite structures is in the material stiffness matrix [E].In the following we present expressions for [E]. 9.1 Three-Dimensional Element The stress-strain relationships for a three-dimensional element are(Eq.2.20) C11 C12 C13 C14 C15 C16 Ex 6 C21 C22 C23 C24 C25 C26 31 C32 C33 C34 C35 C36 y C C42 C43 C44 C45 C46 (9.7) C51 Cs2 C53 C54 Cs6 C61 C62 C63 C64 C65 C66 [E] where [E]is the stiffness matrix for a three-dimensional element. 1R.D.Cook,D.S.Malkus,and M.E.Plesha,Concepts and Applications of Finite Element Analysis 3rd edition.John Wiley and Sons,New York,1989.p.110
396 FINITE ELEMENT ANALYSIS nodal forces nodal points elements Figure 9.1: Structure and its finite element mesh. 9. The stresses at a point inside the element are calculated by σ = [E] ε, (9.4) where [E] is the stiffness matrix characterizing the material. The element stiffness matrix, referred to in Step 2, is defined as [k] δ = fe, (9.5) where fe represents the forces acting at the nodal points of the element. The element stiffness matrix is1 [k] = ) (V) [B] T [E] [B] dV, (9.6) where V is the volume of the element. The preceding steps apply to structures made of either isotropic or composite materials. The only difference between isotropic and composite structures is in the material stiffness matrix [E]. In the following we present expressions for [E]. 9.1 Three-Dimensional Element The stress–strain relationships for a three-dimensional element are (Eq. 2.20) σx σy σz τyz τxz τxy = C11 C12 C13 C14 C15 C16 C21 C22 C23 C24 C25 C26 C31 C32 C33 C34 C35 C36 C41 C42 C43 C44 C45 C46 C51 C52 C53 C54 C55 C56 C61 C62 C63 C64 C65 C66 % &' ( [E] x y z γyz γxz γxy , (9.7) where [E] is the stiffness matrix for a three-dimensional element. 1 R. D. Cook, D. S. Malkus, and M. E. Plesha, Concepts and Applications of Finite Element Analysis. 3rd edition. John Wiley and Sons, New York, 1989, p. 110
9.3 BEAM ELEMENT 397 9.2 Plate Element In the absence of shear deformation,the force-strain relationships for a thin-plate element are (Eq.3.21) A A1 A6 B11 B12 B16 A12 2 A26 B12 B22 B N A6 A26 A66 B16 B26 B66 Bu (9.8) B12 B16 Du D12 D16 B12 B22 B26 D12 D2 D26 Ky B16 B26 B66 D16 D26 D66 Kxy [E where [E]is the stiffness matrix for a plate element without shear deformations.In the presence of shear deformation,the force-strain relationships are (Eqs.5.13- 5.15) N A A2 A16 B11 B12 B16 0 0 A2 An 6 B12 B22 B26 0 A6 6 A66 B6 B26 B6 0 0 M Bu1 B12 B16 D11 D12 D 0 0 Kx M B12 B22 B26 D12 D22 D26 0 B16 B26 B66 D16 D26 D66 0 0 Kxy 0 0 0 0 0 0 0 0 0 0 0 0 s Yyz [目 (9.9) where [E]is the stiffness matrix for a plate element with shear deformations. Frequently,the behavior of shells can be described by replacing the curved surface with small,flat elements.The stiffness matrices above are applicable to such flat shell elements. 9.3 Beam Element For a beam element,the force-strain relationships are arbitrary layup,no shear deformation,no restrained warping (Eq.6.2) N P P2 P3 P4 P2 P2 P2s P24 Pi3 P2s P33 (9.10) 4 P2a P P 11 [可
9.3 BEAM ELEMENT 397 9.2 Plate Element In the absence of shear deformation, the force–strain relationships for a thin-plate element are (Eq. 3.21) Nx Ny Nxy Mx My Mxy = A11 A12 A16 B11 B12 B16 A12 A22 A26 B12 B22 B26 A16 A26 A66 B16 B26 B66 B11 B12 B16 D11 D12 D16 B12 B22 B26 D12 D22 D26 B16 B26 B66 D16 D26 D66 % &' ( [E] o x o y γ o xy κx κy κxy , (9.8) where [E] is the stiffness matrix for a plate element without shear deformations. In the presence of shear deformation, the force–strain relationships are (Eqs. 5.13– 5.15) Nx Ny Nxy Mx My Mxy Vx Vy = A11 A12 A16 B11 B12 B16 0 0 A12 A22 A26 B12 B22 B26 0 0 A16 A26 A66 B16 B26 B66 0 0 B11 B12 B16 D11 D12 D16 0 0 B12 B22 B26 D12 D22 D26 0 0 B16 B26 B66 D16 D26 D66 0 0 000000 S11 S12 000000 S12 S22 % &' ( [E] o x o y γ o xy κx κy κxy γxz γyz , (9.9) where [E] is the stiffness matrix for a plate element with shear deformations. Frequently, the behavior of shells can be described by replacing the curved surface with small, flat elements. The stiffness matrices above are applicable to such flat shell elements. 9.3 Beam Element For a beam element, the force–strain relationships are arbitrary layup, no shear deformation, no restrained warping (Eq. 6.2) N My Mz T = P11 P12 P13 P14 P12 P22 P23 P24 P13 P23 P33 P34 P14 P24 P34 P44 % &' ( [E] o x 1 ρy 1 ρz ϑ (9.10)
398 FINITE ELEMENT ANALYSIS orthotropic,no shear deformation,no restrained warping (Eq.6.8) N EA 0 0 0 就介 0 Elyy 0 0 (9.11) 0 0 0 0 [间 orthotropic,no shear deformation,restrained warping(Eqs.6.8 and 6.233) N EA 0 0 0 0 0 Elyy Ely: 0 0 0 Elyz El. 0 0 1101 (9.12) 0 0 0 0 熙 0 0 0 0 29 [可 orthotropic,shear deformation,restrained warping (Eqs.7.30,7.32,7.34,7.36) N EA 0 0 0 0 0 0 0 e 0 Elyy Elyz 0 0 0 0 0 dx: d 0 Elye El 0 0 0 0 0 dxy d 0 0 EL 0 0 0 0 doB dx 0 0 GI 0 0 0 29 0 0 0 0 Snt Yy 0 0 0 0 0 0 0 0 0 0 5 (E] (9.13) where [E]is the stiffness matrix for a beam element. 9.4 Sublaminate A laminate consisting of several plies may be analyzed by either plate(flat shell) or three-dimensional elements (Fig.9.2).For thick laminates neither of these Figure 9.2:Thick laminate (top),analysis with plate elements (left).analysis with three-dimensional elements(right)
398 FINITE ELEMENT ANALYSIS orthotropic, no shear deformation, no restrained warping (Eq. 6.8) N My Mz T = EA 0 00 0 EI yy EI yz 0 0 EI yz EI zz 0 00 0 GI t % &' ( [E] o x 1 ρy 1 ρz ϑ (9.11) orthotropic, no shear deformation, restrained warping (Eqs. 6.8 and 6.233) N My Mz Mω T sv = EA 0 0 00 0 EI yy EI yz 0 0 0 EI yz EI zz 0 0 00 0 EI ω 0 00 00 GI t % &' ( [E] o x 1 ρy 1 ρz −dϑ dx ϑ (9.12) orthotropic, shear deformation, restrained warping (Eqs. 7.30, 7.32, 7.34, 7.36) N My Mz Mω T sv V y V z T ω = EA 0 0 0 0000 0 EI yy EI yz 0 0000 0 EI yz EI zz 0 0000 00 0 EI ω 0000 00 00 GI t 000 00 0 00 Syy Syz Syω 00 0 00 Syz Szz Szω 00 0 00 Syω Szω Sωω % &' ( [E] o x −dχz dx −dχy dx −dϑB dx ϑ γy γz ϑS , (9.13) where [E] is the stiffness matrix for a beam element. 9.4 Sublaminate A laminate consisting of several plies may be analyzed by either plate (flat shell) or three-dimensional elements (Fig. 9.2). For thick laminates neither of these Figure 9.2: Thick laminate (top), analysis with plate elements (left), analysis with three-dimensional elements (right).
9.4 SUBLAMINATE 399 Laminate Sublaminates FE Mesh plies sublaminates Figure 9.3:Thick laminate(left),sublaminates(middle),and the finite element mesh(right). is practical.Plate elements give inaccurate results.Three-dimensional elements require that the material be uniform throughout the element,and,hence,an ele- ment must contain a single layer or adjacent identical layers.This may result in a very large number of elements,making the numerical computation difficult and often infeasible. We can overcome these difficulties by dividing the laminate into sublaminates (Fig.9.3).Each layer in the sublaminate may be monoclinic,orthotropic,trans- versely isotropic,or isotropic.The thickness of each element is the same as the thickness of the corresponding sublaminate.The stiffness matrix [E]of such a sublaminate is defined by the relationship Ex =[ (9.14) 7y2 气xy The bar denotes average stresses and strains.It is convenient to represent this expression in terms of the compliance matrix [ J11 J12 J13 J14 J15 J16] x Jy J24 J25 J26 J31 J32 /33 J54 J35 J36 J41 J J 5 J4 (9.15) J51 J52 J53 J54 J55 56 J62 163 J64 165 J66 元y where [=]-1 (9.16)
9.4 SUBLAMINATE 399 Laminate Sublaminates FE Mesh plies sublaminates Figure 9.3: Thick laminate (left), sublaminates (middle), and the finite element mesh (right). is practical. Plate elements give inaccurate results. Three-dimensional elements require that the material be uniform throughout the element, and, hence, an element must contain a single layer or adjacent identical layers. This may result in a very large number of elements, making the numerical computation difficult and often infeasible. We can overcome these difficulties by dividing the laminate into sublaminates (Fig. 9.3). Each layer in the sublaminate may be monoclinic, orthotropic, transversely isotropic, or isotropic. The thickness of each element is the same as the thickness of the corresponding sublaminate. The stiffness matrix [E] of such a sublaminate is defined by the relationship σ x σ y σ z τ yz τ xz τ xy = [E] x y z γ yz γ xz γ xy . (9.14) The bar denotes average stresses and strains. It is convenient to represent this expression in terms of the compliance matrix [J ] x y z γ yz γ xz γ xy = J11 J12 J13 J14 J15 J16 J21 J22 J23 J24 J25 J26 J31 J32 J33 J34 J35 J36 J41 J42 J43 J44 J45 J46 J51 J52 J53 J54 J55 J56 J61 J62 J63 J64 J65 J66 % &' ( [J ] σ x σ y σ z τ yz τ xz τ xy , (9.15) where [E] = [J ] −1 . (9.16)
400 FINITE ELEMENT ANALYSIS K 尽： Figure 9.4:Illustration of a sublaminate. The x,y,z coordinate system is shown in Figure 9.4.The average stresses are defined as 1 石xhsJh) xdz= 石：=02 1 1 yhsJh Ny Tyz=Tyz (9.17) 1 元xy=hsJh) xydi=hs Nxy Txz Txz- The second equalities in the left-hand column are written by virtue of Eq.(3.9). The average strains are 1「 e.dz Ex -Ex 1 7= Ey=Ey (9.18) 7=元小 Yxadz Txy=Yxy where hs is the thickness of the sublaminate(Fig.9.4).The terms in the right-hand columns of Eqs.(9.17)and(9.18)show that the stresses o,tyz,txz and the strains er,Ey,Yry do not vary across the thickness. In the following we derive the elements of the compliance matrix. 9.4.1 Step 1.Elements of [due to In-Plane Stresses In this step we determine the elements in the first,second,and sixth columns of the matrix [To this end,we impose the average in-plane stresses y,and Try on the sublaminate(Fig.9.5,top).Since,yr are zero,the stress-strain
400 FINITE ELEMENT ANALYSIS hs z y x z y x hs 1 … ks Ks zk −1 zk Figure 9.4: Illustration of a sublaminate. The x,y,z coordinate system is shown in Figure 9.4. The average stresses are defined as σ x = 1 hs ) (hs) σxdz= 1 hs Nx σ z = σz σ y = 1 hs ) (hs) σydz= 1 hs Ny τ yz = τyz τ xy = 1 hs ) (hs) τxydz= 1 hs Nxy τ xz = τxz. (9.17) The second equalities in the left-hand column are written by virtue of Eq. (3.9). The average strains are z = 1 hs ) (hs) zdz x = x γ yz = 1 hs ) (hs) γyzdz y = y γ xz = 1 hs ) (hs) γxzdz γ xy = γxy, (9.18) where hs is the thickness of the sublaminate (Fig. 9.4). The terms in the right-hand columns of Eqs. (9.17) and (9.18) show that the stresses σz, τyz, τxz and the strains x, y, γxy do not vary across the thickness. In the following we derive the elements of the compliance matrix. 9.4.1 Step 1. Elements of [J ] due to In-Plane Stresses In this step we determine the elements in the first, second, and sixth columns of the matrix [J ]. To this end, we impose the average in-plane stresses σ x, σ y, and τ xy on the sublaminate (Fig. 9.5, top). Since σ z, τ yz, τ xz are zero, the stress–strain
9.4 SUBLAMINATE 401 0. L(1+E) h(1+e) Figure 9.5:Illustration of Step 1.The ply stress and the corresponding average stress on the sublaminate. relationship (Eq.9.15)may be written as J12 J16 J21 J J26 (9.19) 16 J62 e}=[J31 J32 J36] (9.20) 「J4 J46 6 (9.21) J56 69 The strains are uniform across the thickness(Eq.9.18).Under these conditions Kx,Ky,Kxy are zero,and we have (see Egs.3.21 and 3.7) (9.22) where [A]is the tensile stiffness matrix of the sublaminate.(The summation in Eq.3.20 is performed from 1 to Ks,where Ks is the number of layers or ply groups in the sublaminate;see Fig.9.4).For the sublaminate Egs.(9.17),(9.18),and(9.22) yield [A (9.23) J11 J12 J16 J21 J26 J61 J62 J6小
9.4 SUBLAMINATE 401 σx σx L(1+ )x hs (1+ )z σx σx hs L L σx σy τxy Figure 9.5: Illustration of Step 1. The ply stress and the corresponding average stress on the sublaminate. relationship (Eq. 9.15) may be written as x y γ xy = J11 J12 J16 J21 J22 J26 J61 J62 J66 σ x σ y τ xy (9.19) {z} = [J31 J32 J36] σ x σ y τ xy (9.20) 1 γ yz γ xz6 = J41 J42 J46 J51 J52 J56! σ x σ y τ xy . (9.21) The strains are uniform across the thickness (Eq. 9.18). Under these conditions κx, κy, κxy are zero, and we have (see Eqs. 3.21 and 3.7) Nx Ny Nxy = [A] o x o y γ o xy = [A] x y γxy , (9.22) where [A] is the tensile stiffness matrix of the sublaminate. (The summation in Eq. 3.20 is performed from 1 to Ks, where Ks is the number of layers or ply groups in the sublaminate; see Fig. 9.4). For the sublaminate Eqs. (9.17), (9.18), and (9.22) yield x y γ xy = hs [A] −1 % &' ( J11 J12 J16 J21 J22 J26 J61 J62 J66 σ x σ y τ xy . (9.23)
402 FINITE ELEMENT ANALYSIS By comparing Eqs.(9.19)and(9.23),we have Ju Jv 67 J26 =h[4-1 (9.24) J61 J62J66 Owing to the Poisson effect,the in-plane forces introduce a normal strain ea in the z direction.For a single layer this strain is (see Eq.2.133) Ox e2=3 23 S36 (9.25) Txy where ox,oy,and txy are the stresses in the plies.The average normal stress x across the laminate and the normal ply stress ox are illustrated in Figure 9.5 (bottom).The stresses in a single layer are (Eq.3.13) =[@ (9.26) Equations(9.25),(9.26),(9.23),and(9.18)give e2=323536]ohs[4 (9.27) By combining Eqs.(9.18)and(9.27),and by replacing the integral with a sum- mation,for the sublaminate we obtain Ks 石x 5g5356](k-2-1)@k)[4 (9.28) k=1 [U31J32J36] By comparing Eqs.(9.20)and (9.28),we have 1g=之Sg5k(a-a-@[4, (9.29)
402 FINITE ELEMENT ANALYSIS By comparing Eqs. (9.19) and (9.23), we have J11 J12 J16 J21 J22 J26 J61 J62 J66 = hs [A] −1 . (9.24) Owing to the Poisson effect, the in-plane forces introduce a normal strain z in the z direction. For a single layer this strain is (see Eq. 2.133) z = [S13 S23 S36] σx σy τxy , (9.25) where σx, σy, and τxy are the stresses in the plies. The average normal stress σ x across the laminate and the normal ply stress σx are illustrated in Figure 9.5 (bottom). The stresses in a single layer are (Eq. 3.13) σx σy τxy = [Q] x y γxy . (9.26) Equations (9.25), (9.26), (9.23), and (9.18) give z = [S13 S23 S36][Q]hs [A] −1 σ x σ y τ xy . (9.27) By combining Eqs. (9.18) and (9.27), and by replacing the integral with a summation, for the sublaminate we obtain z = * Ks k=1 ([S13 S23 S36](zk − zk−1)[Q]k)[A] −1 % &' ( [J31 J32 J36] σ x σ y τ xy . (9.28) By comparing Eqs. (9.20) and (9.28), we have [J31 J32 J36] = * Ks k=1 ([S13 S23 S36]k (zk − zk−1)[Q]k)[A] −1 , (9.29)
9.4 SUBLAMINATE 403 where zk,Z-1 are illustrated in Figure 9.4.Since Tyz=Txz=0,we have that yyz and yxz are zero (Eqs.2.26 and 2.21),and Eq.(9.21)now becomes 63 (9.30) To satisfy this equation,the preceding elements of the compliance matrix must be zero: J J42 J46 001 51J52 J56」 0 00 (9.31) 9.4.2 Step 2.Elements of [due to Out-of-Plane Normal Stresses In this step we determine the elements in the third column of the matrix []To accomplish this we consider a sublaminate in which there are only o:stresses.To form such a sublaminate,first we consider a sublaminate restrained along its edges (ex =ey=yxy =0)and apply a uniform stress o:=:(Stage 1,Fig.9.6),which introduces in-plane stresses,,y.Second,we apply in-plane stressesx, -y,-Txy(Stage 2,Fig.9.6).Third,we superimpose Stages 1 and 2 and arrive at the sublaminate inside which the in-plane average stresses oy,Try and the transverse shear stresses y,x are zero and the stress-strain relationships (Eq.9.15)are 「J3 J23 (9.32) Yxy J63 e2=J33可： (9.33) - (9.34) Stage 1.Following the outline above,we apply o:to the sublaminate,which is restrained along the edges(x=y=y=0).For a single layer,Eqs.(2.27)and 个个 个个个个个 6 ↓↓↓↓↓ ↓↓↓↓ Stage 1 Stage 2 Stage 3 Figure 9.6:Illustration of Step 2.Sublaminate restrained along its edges subjected to o:and the resulting stressx(Stage 1);unrestrained sublaminate subjected tox(Stage 2);unrestrained sublaminate subjected to o:(Stage 3);(y and Txy are not shown)
9.4 SUBLAMINATE 403 where zk, zk−1 are illustrated in Figure 9.4. Since τ yz = τ xz = 0, we have that γyz and γxz are zero (Eqs. 2.26 and 2.21), and Eq. (9.21) now becomes 0 0 = J41 J42 J46 J51 J52 J56! σ x σ y τ xy . (9.30) To satisfy this equation, the preceding elements of the compliance matrix must be zero: J41 J42 J46 J51 J52 J56! = 000 000! . (9.31) 9.4.2 Step 2. Elements of [J ] due to Out-of-Plane Normal Stresses In this step we determine the elements in the third column of the matrix [J ]. To accomplish this we consider a sublaminate in which there are only σz stresses. To form such a sublaminate, first we consider a sublaminate restrained along its edges (x = y = γxy = 0) and apply a uniform stress σz = σ z (Stage 1, Fig. 9.6), which introduces in-plane stresses σ x, σ y, τ xy. Second, we apply in-plane stresses −σ x, −σ y, −τ xy (Stage 2, Fig. 9.6). Third, we superimpose Stages 1 and 2 and arrive at the sublaminate inside which the in-plane average stresses σ x, σ y, τ xy and the transverse shear stresses τ yz, τ xz are zero and the stress–strain relationships (Eq. 9.15) are x y γ xy = J13 J23 J63 σ z (9.32) z = J33σ z (9.33) γ yz γ xz = J43 J53! σ z. (9.34) Stage 1. Following the outline above, we apply σz to the sublaminate, which is restrained along the edges (x = y = γxy = 0). For a single layer, Eqs. (2.27) and σx σx + = σx σx σz σz Stage 1 Stage 2 Stage 3 Figure 9.6: Illustration of Step 2. Sublaminate restrained along its edges subjected to σz and the resulting stress σ x (Stage 1); unrestrained sublaminate subjected to σ x (Stage 2); unrestrained sublaminate subjected to σz (Stage 3); (σ y and τ xy are not shown)
404 FINITE ELEMENT ANALYSIS (2.20)give 0:=C33E2 (9.35) C23 (9.36) Txy These equations may be rearranged to yield e= (9.37) 1 ℃3 (9.38) By combining Egs.(9.37)and(9.18),and by replacing the integrals by sum- mations,for the sublaminate we obtain (9.39) (C33)k Equations (9.38)and(9.18)yield 卧 -2 (C3) 0 (9.40) This equation gives the in-plane stresses in the restrained sublaminate Stage 2.We apply the equal and opposite of the stresses,calculated by Eq.(9.40), to the sublaminate.The corresponding strains are obtained from Eqs.(9.28),(9.23), and (9.40).Equations (9.28)and (9.40)result in E2=-[31 J32 J36] =-[J31 Zk Zk-1 (9.41) (C33)k Equations(9.23)and(9.40)give Ex J11 J12 J16 J21 J61 J16 J26 k-k-1 (9.42) 16 s (C3)k
404 FINITE ELEMENT ANALYSIS (2.20) give σz = C33z (9.35) σx σy τxy = C13 C23 C63 z. (9.36) These equations may be rearranged to yield z = 1 C33 σz (9.37) σx σy τxy = C13 C23 C63 1 C33 σz. (9.38) By combining Eqs. (9.37) and (9.18), and by replacing the integrals by summations, for the sublaminate we obtain z = 1 hs * Ks k=1 zk − zk−1 (C33)k σz. (9.39) Equations (9.38) and (9.18) yield σ x σ y τ xy = 1 hs * Ks k=1 C13 C23 C63 k zk − zk−1 (C33)k σz. (9.40) This equation gives the in-plane stresses in the restrained sublaminate. Stage 2. We apply the equal and opposite of the stresses, calculated by Eq. (9.40), to the sublaminate. The corresponding strains are obtained from Eqs. (9.28), (9.23), and (9.40). Equations (9.28) and (9.40) result in z = −[J31 J32 J36] σ x σ y τ xy = −[J31 J32 J36] 1 hs * Ks k=1 C13 C23 C63 k zk − zk−1 (C33)k σz. (9.41) Equations (9.23) and (9.40) give x y γ xy = − J11 J12 J16 J21 J22 J26 J61 J62 J66 σ x σ y τ xy = − J11 J12 J16 J21 J22 J26 J61 J62 J66 1 hs * Ks k=1 C13 C23 C63 k zk − zk−1 (C33)k σz. (9.42)