Investigation of Through-Thickness Stresses in Composite Laminates Using Layerwise Theory

In this study, an analyticalmethod is developed to exactly obtain the interlaminar stresses near the free edges of laminated composite plates under the bending moment based on the reduced form of elasticity displacement field for a long laminate. The analytical and numerical studies were performed based on the Reddy’s layerwise theory for the boundary layer stresses within cross-ply, symmetric, angle-ply, and general composite laminates. Finally, a variety of numerical results are presented for the interlaminar normal and shear stresses along the interfaces and through thickness of laminates near the free edges. The results showed high stress gradient of interlaminar normal and shear stresses near the edges of laminates.


Introduction
Due to high specific strength and stiffness of fiber reinforced polymer composite, laminated composites have found augmented use in many industrial applications.In the boundarylayer regions, due to geometry and material discontinuities, the interlaminar stresses exhibit much higher values than those predicted by the classical lamination theory (CLT).These highly concentrated stresses cause delaminate failure in the laminates.However, no exact solution is found for elasticity equations because of inherent complexities involved in the problem of finding exact stress values in the edges.Hence, different analytical and numerical methods for finding the interlaminar stresses are used to describe the interlaminar stresses at the free edges of composite laminates.Complete literature surveys on this subject are available in review articles of Kant and Swaminathan [1] which obviously show the detailed path of development of methods.The first approximate analysis of interlaminar stresses was presented by Puppo and Evensen [2].They studied interlaminar shear stresses in an idealized laminate consisting of orthotropic layers separated by isotropic shear layers with interlaminar normal stress being neglected through the laminate.Other approximate analytical methods utilized to examine the problem consist of the use of the higher-order plate theory by Pagano [3], the perturbation technique by Hsu and Herakovich [4], the boundary layer theory by Tang and Levy [5], and the approximate elasticity solutions by Pipes and Pagano [6].An approximate theory is also employed by Pagano [7] based on the assumed inplane stresses and the use of Reissner's variational principle.Wang and Choi [8] utilized Lekhnistskii's stress potential and the theory of anisotropic elasticity for examining the free edge singularities.A variational approach concerning Lekhnitskii's stress functions is used by Yin [9] for the evaluation of free-edge stresses in laminates under uniaxial tension, bending, and torsion.Wang and Crossman [10] developed a quasi-three-dimensional finite element solution to determine the free-edge stresses in a symmetric balanced composite laminate under uniaxial tension and uniform thermal loading.Whitcomb et al. [11] studied the differences in numerical results for interlaminar stresses obtained by various methods (finite difference methods, finite element methods, and perturbation techniques).Boundary element method and the integral equation theory were used by Davì [12] to study the stresses in a general laminate under uniform axial strain.Carrera and Demasi [13,14]  by using the Reissner mixed variational theorem (RMVT).They compared the numerical results for interlaminar stresses in several finite-element models and elasticity theory within composite laminates and sandwich plates.Nguyen and Caron [15] employed a multiparticle finite element method to study the interlaminar stresses near the free edges of general composite laminates under mechanical and thermal loading.Robbins and Reddy [16] used a displacement-based variable kinematic global-local finite element method.Mittelstedt and Becker [17] utilized Reddy's layerwise laminate plate theory to find the closed-form analysis of free-edge effects in layered plates of arbitrary nonorthotropic layups.The approach consists of the subdivision of the physical laminate layers into an arbitrary number of mathematical layers through the plate thickness.Na [18] used a finite element model based on the layerwise theory.He employed the von Kármán type nonlinear strains to analyze damage in laminated composite beams.In his formulation, the Heaviside step function is employed to express the discontinuous interlaminar displacement field at the delaminated interfaces.Plagianakos and Saravanos [19] presented a higher-order layerwise theoretical framework, which enables prediction of the static response of thick composite and sandwich composite plates.Ullah et al. [20] carried out some experimental tests to characterize the behavior of a woven CFRP material under large-deflection bending.Twodimensional finite element (FE) models were implemented in the commercial code Abaqus.They performed series of simulations to study the deformation behavior and damage in CFRP for cases of high-deflection bending.Hélénon et al. [21] presented an experimental and numerical investigation into failure of T-shaped laminated composite structures.Three out-of-plane bending cases are studied.They found that very high free-edge maximum principal transverse tensile stresses perpendicular to the fiber direction occur at the failure locations.Thai et al. [22] indicated an isogeometric finite element formulation for static, free vibration and buckling analysis of laminated composite and sandwich plates.Their method allows removing shear correction factors and improves the accuracy of transverse shear stresses.Thai et al. [23] investigated the behavior of laminated composites using several high order or layerwise finite element calculations.A layerwise model and its dedicated  ∘ eight-node finite element were specifically developed for interlaminar stresses analysis in free edge problem.Malekzadeh [24] developed a high accuracy and rapid convergence hybrid approach for two-dimensional static analyses of circular arches with different boundary conditions.The method essentially consists of a layerwise theory used for the thickness direction in conjunction with differential quadrature method in the axial direction.
There have been very limited works to study the interlaminar stresses subjected to the bending moment.In the present paper, by the use of Reddy's LWT, an analytical solution is presented to evaluate interlaminar stresses in cross-ply, symmetric, angle-ply, and general composite laminates under the bending moment.To commence with, based on physical arguments regarding the deformations of a long generally and other laminated composite plates, an appropriate reduced elasticity displacement field is established.The boundarylayer stresses within the laminate are obtained analytically based on Reddy's LWT.

Elasticity Displacement
Field.An Nth-layered composite plate (with arbitrary lamination) under the bending moment is considered as shown in Figure 1.The coordinate system (, , ) is located at the middle plane of the laminate, that is of thickness ℎ, width 2, and length 2 and is assumed to be long in the  direction so that the strains away from the ends ( = ±) of the laminate are functions of only  and .
The integrations of the three-dimensional elasticity strain-displacement relations [25] inside the kth layer of the laminate will generate the most general form of displacement field which can be shown to be [26]  ()  1 (, , ) =  () 4  +  () 6  +  () 2  +  () (, ) where  () 1 ,  () 2 , and  () 3 represent the displacement components of the material point (, , ) in the , , and  International Journal of Engineering Mathematics 3 directions, respectively.To satisfy the displacement continuity conditions at the interfaces of the adjoining layers, the integration constants in (1) must be the same for all layers within the laminate.Therefore, the most general displacement field in the kth layer is shown: If loading conditions at  = − and  are similar, based on the physical grounds, the following conditions hold: Upon imposing these restrictions on (1) it is readily seen that the constants  4 and  5 must vanish and the displacement field in ( 1) is reduced to the following equations: Also, by replacing  () (, ) by − 3  +  () (, ) in (4a), it becomes apparent that terms involving  3 in (4a)-(4c) can be neglected since no strains are generated by such terms (these terms will correspond to an infinitesimal rigid-body rotation of the laminate about the z-axis in Figure 1).Thus, the most general form of the displacement field of an arbitrary laminate is given by The displacement field in (5a), (5b), and (5c) may be used, in principle, for computing the stress field in any composite laminate subjected to arbitrary combinations of self-equilibrating mechanical and uniform hygrothermal loads.The terms involving  1 ,  2 , and  6 in (5a), (5b), and (5c) demonstrate certain global deformations that occurred in the laminate.On the other hand, the unknown functions appearing in (5a), (5b), and (5c) illustrate the local deformations.They will be evaluated by LWT of Reddy that happened in a laminate.
For cross-ply laminates subjected to mechanical loadings, based on physical grounds, the following restrictions hold (see Figure 1): Upon imposing these conditions on (5a), it is readily seen that  () (, ) = 0. Also it is concluded from (6b) and (5b) that  1 = 0. Thus, for cross-ply laminates the most general form of the displacement field is given: Also, for symmetric laminates based on physical grounds the following conditions must hold (see Figure 1): Upon imposing these conditions on (5a), (5b), and (5c), it is found that  1 =  6 = 0. Thus, for symmetric laminates the most general form of displacement field is expressed by Finally, for antisymmetric angle-ply laminates with  layers subjected to mechanical loadings, the following condition must hold: ) . ( Concluding from (10) and (5a),  6 = 0 and, therefore, the most general displacement field for such laminates are written as International Journal of Engineering Mathematics The displacement fields in ( 7), (9), and ( 11) can be represented in one place as where   =   = 0 and   = 1 for general cross-ply laminates,   =   = 0 and   = 1 for symmetric laminates,   =   = 1 and   = 0 for angle-ply laminates, and   =   =   = 1 for general laminates.

Layerwise Laminated Plate Theory of Reddy
In this section, LWT is used for calculation of interlaminar stresses in arbitrary laminated composite plates with free edges.In LWT, it is possible to replace the actual physical layers to many mathematical layers [27].This can be done by subdividing each physical layer through the thickness, to get the proper accuracy.The theory assumes that the displacement components of a generic point in the laminate are given by [27] It should be emphasized that a repeated index indicates summation from 1 to  + 1.The functions   (, ), V  (, ), and   (, ) represent the displacements of the points initially located at the kth plane of the laminate in the , , and  directions, respectively.The variable  in (13) is the total number of numerical layers introduced in any laminate.Also the functions   's are the global approximation functions of the thickness coordinate which are assumed here to be linear [27].This function expressed as where the local Lagrangian interpolation functions    () ( = 1, 2) related with the kth surface in the laminate are defined as where ℎ  is the thickness of the kth mathematical layer.Based on the displacement field in (12), the displacement field of LWT in (13) takes the following form: It is noted that, by the means of through-the-thickness linear interpolation functions, the continuity of displacement components through the thickness of laminate is identically satisfied.On the other hand, the transverse strain components remain discontinuous at the interfaces.They will ultimately emphasize the feasibility of having continuous interlaminar stresses at the interfaces of adjoining layers by increasing the number of numerical layers through the physical laminate.Substituting (15) into the linear straindisplacement relations of elasticity [25] yields where a prime in ( 16) denotes ordinary differentiation with respect to a suitable independent variable (i.e., either  or ).
The three-dimensional constitutive law within the kth layer (with arbitrary fiber orientation) of the laminate may also be shown as [28] ( ( where [] () is the transformed (i.e., off axis) stiffness matrix of the kth layer.By using the strain-displacement relations International Journal of Engineering Mathematics 5 (16) in the principle of the minimum total potential energy [25], the equilibrium equations within LWT are found as Also the traction-free boundary conditions at the free edges of the laminate are given as In the same way, upon substitution of ( 16) into (17) where the extra laminate rigidities appearing in (24a), (24b), and (24c) are defined in the Appendix.The system of equations in (23a), (23b), and (23c) shows 3( + 1) coupled ordinary differential equations with constant coefficients which may be indicated in a matrix form as where The coefficient matrices [], [], and [] in (25) are displayed in the Appendix.Also from the conditions in (3) and the displacement field in (23a), (23b), and (23c) it is observed that the functions   ,   , and   are all odd functions of the independent variable .It can be confirmed that the general solution of (25) may be presented as where [sinh()] is a 3( + 1) × 3( + 1) diagonal matrix, and Also [] and ( 2 1 ,  2 2 , . . .,  2 3(+1) ) are the model matrix and eigenvalues of (−[] −1 []), respectively.{} is an unknown vector representing 3( + 1) integration constants.The constants   ( = 1, 2, 6) must be calculated within LWT analysis.So, the boundary conditions in (20) are first imposed to found the vector {} in terms of the unknown parameters   ( = 1, 2, 6).These constants are next calculated in terms of the specific bending moment  0 by satisfaction of the global equilibrium conditions in (19a), (19b), and (19c).

Numerical Results and Discussions
In this part, several numerical examples are shown for the distribution of interlaminar stress.All physical layers are assumed to be of equal thickness (=0.5 mm).The mechanical properties of unidirectional graphite/epoxy T300/5208 used in this study are shown in Table 1 [28].Also, each physical ply is modeled as being made up of 12 numerical layers within LWT (i.e.,  = 12).In addition, the width thickness ratio (i.e., 2/ℎ) is assumed to be equal to 10. Furthermore, the stress components are normalized as where  0 =  0 /ℎ 2 .Figure 2 shows the distribution of interlaminar normal and shear stresses along the width of [45 ∘ /−45 ∘ ]  laminate.Good correspondence was found between the present results and the results obtained by Goodsell and Pipes [29] which confirms the accuracy of the present solution.Figure 3 shows the distribution of interlaminar stresses   and   along the −45 ∘ /0 ∘ and   along the 0 ∘ /90 ∘ interfaces of the general [45 ∘ /− 45 ∘ /0 ∘ /90 ∘ /30 ∘ /− 30 ∘ ] laminate.It is seen that the interlaminar stresses show high stress gradient near the free edge.It is observed that the interlaminar stresses   and   grow quickly near the free edge, while being zero in the inner region of the laminate.Also the magnitude of the transverse shear stress   is greater than that of transverse normal stress   .On the other hand,   rises toward the free edge and decreases suddenly to zero at free edge.
The interlaminar normal stress   along the 90 ∘ /0 ∘ interface of the unsymmetric cross-ply [90 ∘ /0 ∘ /90 ∘ /0 ∘ /90 ∘ /0 ∘ ] laminate for various width-to-thickness ratios is demonstrated in Figure 4.It is observed that by decreasing the width-to-thickness ratio, the boundary-layer region can be expanded towards the internal region of the laminate with its width being almost equal to the thickness of the laminate.It is seen that the magnitude of the interlaminar stress at the free edge is not varied as the width-to-thickness ratio of the laminate is changed.In addition, the highly localized nature of interlaminar stresses near and exactly at the free edges of the laminate is apparent.
The distribution of interlaminar stresses along the 0 ∘ /90 ∘ and 90 ∘ /0 ∘ interfaces of the unsymmetric cross-ply [0 ∘ /90 ∘ /0 ∘ /90 ∘ ] and [90 ∘ /0 ∘ /90 ∘ /0 ∘ ] laminates, respectively, are compared in Figure 5.It is seen that the maximum interfacial value   at the 0 ∘ /90 ∘ and 90 ∘ /0 ∘ interfaces take place at the free edge.Figure 6 shows the variation of interlaminar normal stress   on through-the-thickness and near the free edge of the angle-ply [45 ∘ /− 45 ∘ /45 ∘ /− 45 ∘ ] laminate as the free edge is approached.It is observed that  the maximum negative and positive values of   occur within the top −80 ∘ layer and the bottom 80 ∘ layer at the free edge (i.e.,  = ), respectively.It is also noticed that   diminishes away from the free edge as the interior region of the laminate is approached.Moreover,   reduces by moving slightly away from the free edge and its distribution becomes smoother.
The variations of the interlaminar shear stress   at  =  through the thickness of the [0 ∘ /60 ∘ /−60 ∘ ]  laminate are displayed in Figure 7.By increasing the number of numerical layers in each lamina,   becomes slightly closer to zero but International Journal of Engineering Mathematics    it never becomes zero.It should be noted that increasing the number of subdivisions results in no convergence for   and   at interface-edge junction of two dissimilar layers.The numerical values of these components (  and   ) continue to grow as the number of sublayers increased.Since the generalized stress resultants    instead of   are forced to disappear at the free edge in LWT, the numerical value of   may never become zero at the interface-edge intersection (even by increasing the number of sublayers in each physical layer).
The distribution of interlaminar normal stress   on through thickness of the general [−15 ∘ /90 ∘ /90 ∘ /−15 ∘ ] laminate is displayed in Figure 8.It is noted that increasing the number of numerical layers () in each lamina has no significant effect on the values of the interlaminar stress   within the boundary-layer region of the laminate, expect exactly at the free edge (i.e.,  = ).Furthermore, the interlaminar normal stress   grows rapidly in the vicinity of the free edge, while being zero in the interior region of the laminate.

Conclusions
In the present paper, analytical solutions are developed to analyze the interlaminar stresses for different lay-up configurations in laminated composite plates subjected to the bending moment.The solutions were obtained based on the reduced elasticity displacement field for long laminated plates.It is found that the components of the displacements field are composed of two distinct parts, signifying the global and local deformations within laminates.The local displacement functions as well as the interlaminar stresses through the layers of the laminate were found by Reddy's layerwise theory (LWT).The numerical results based on LWT were shown for the free-edge interlaminar stresses through the thickness and across the interfaces of various cross-ply, symmetric, angle-ply, and general composite laminates.

4 −
z at z =h/4  z at z  yz  yz at z  yz at z = h/