- Email: [email protected]
Modeling the stress field in laminated composite plates near discontinuities
Composite Structures 3 (1985) 145-166
Modeling the Stress Field in Laminated Composite Plates Near Discontinuities
Eric R. Johnson and Brian L. K e m p Aerospace and Ocean Engineering Department, Virginia Polytechnic Institute and State University, Blacksburg, Virginia 24061, USA
ABSTRACT A structural model is developed using Reissner's variational principle to analyze stress gradient problems in laminated composite plates. A statically admissible stress field, which has explicit dependence on the thickness coordinate, is selected for each layer. Both interlaminar traction and displacement continuity are enforced between layers. This model is similar to an earlier one developed by Pagano, but has four dependent variables less per layer than Pagano's model. It is shown that the model can duplicate Pagano's results for the same number of subdivisions through the thickness. The uniform axial extension of a [0/90Islaminate is used as the example problem for comparison to Pagano' s theory and to illustrate the methodology.
1 INTRODUCTION Failure prediction capability in laminated composite structures depends upon a reasonably accurate analysis of the stress response. For example, to predict the initiation of delamination, or the stress redistribution in the vicinity of damage, suggests that the theory used to determine stresses includes all stress components. Continuity of displacements and stresses at interfaces between laminae should also be required of the theory. Elasticity theory would satisfy these criteria, but currently threedimensional elasticity analysis is intractable for practical laminated composite structures. On the other hand, structural theories for 145 Composite Structures 0263-8223/85/$03-30 O Elsevier Applied Science Publishers Ltd, England, 1985. Printed in Great Britain
146
Eric R. Johnson, Brian L. Kemp
laminates which are based on assumptions for the displacements in the thickness coordinate through the entire laminate do not satisfy these criteria, but are more amenable to solution. In general, structural theories based on assumed displacement fields (this includes finite element formulations of this type) cannot satisfy continuity of stresses along laminae interfaces. In particular, classical structural theory based on the Kirchhoff-Love hypotheses neglect thickness stress components entirely. The deficiency of displacement-based structural theories for the analysis of stress concentration problems in laminated composites is discussed in some detail by Pagano.' In his paper, Pagano establishes the criteria for a new laminate structural theory; i.e. all stress components are nonzero, continuity of displacements and stresses is satisfied at interfaces, and that individual layers satisfy force and moment equilibrium. Pagano develops a structural theory which satisfies these criteria and which is based on assumed stress fields within each layer. Reissner's variational principle is employed to obtain a consistent mathematical formulation. The purpose of this paper is to present a simplification of Pagano's original theory. Pagano's theory has 23N dependent variables, where N is the number of layers. It is shown that with a different assumption on the thickness normal stress (o-:) in each layer the number of dependent variables in the formulation can be reduced by four per layer. Thus, in the general case, the model developed herein would have 19N dependent variables. It is also shown that this simpler model can duplicate Pagano's results for the same number of layers. This is illustrated by solving the problem of a [0/90]~ laminate subject to uniform axial extension.
2 L A M I N A T E V A R I A T I O N A L PRINCIPLE The mechanism for deriving the structural theory is Reissner's variational principle. The representation of this principle for laminated bodies of the type considered in this paper is given by Pagano. ~ The body is a flat rectangular plate consisting of an integral number (N) of layers which are flat rectangular sheets having identical dimensions. Each layer is assumed to be homogeneous, anisotropic, and linear elastic. Also displacements and strains are assumed to be small when the laminate is subjected to prescribed surface tractions and/or displacements. Beginning at one of the two surfaces of the body which are parallel to the interfacial planes,
Modelingthestressfieldinlaminatedplates
147
the layers are n u m b e r e d consecutively from one to N. For the kth layer (k = 1, 2, • . . , N) let Vk designate its volume, ~_~k~ j (i, j = 1,2,3) denote the stresses, u]k~ the displacements, and let W ~k~ represent the complementary strain energy density. Also let S' represent the portion of the laminate's boundary surface area S on which tractions are prescribed, and let S " represent the portion of S on which displacements are prescribed. The traction components ri acting on a surface are related to the stresses by ri = njo-j~, in which nj are the components of the unit outward normal vector to the surface, and repeated indices are summed from one to three. Prescribed surface tractions are denoted by ~-~, and prescribed displacements by fi~. Let a comma designate partial differentiation. Neglecting body forces, the Reissner functional is
J = k~ ' =
uj.i) - W
dV k
!/k
-- f s, ~'iuidS- ~s,, "t"i(l'li--(ti)dS
(1)
The equations of laminate elasticity follow from the vanishing of the first variation of J, which is written as 8J = 0. Both stresses and displacements are varied in this process. After applying Gauss's integral t h e o r e m and using symmetry of the stress tensor, the results of this procedure are
k~=,f vk { [2 (uiJ+ujJ)- ~] 80"~i--o'J"'SU'} 'k)dV k + L, ('~-r,)Su~dS- ~, ~r,(u,- fi,)dS +
[r}k)~Ulkl+'~k+l)SUlk+l)]dI~ = 0 k = I
(2)
Ik'"
In eqn (2) Ik designates the area of the interracial plane c o m m o n to layers k and k + 1, ,?/k) are the interracial traction components acting on the upper surface of the kth layer, and ~.]k+~) are the interracial traction c o m p o n e n t s acting on the lower surface of layer k + 1. Let Ik,' It," and It'", designate the portions of the interracial planes on which tractions are prescribed, displacements are prescribed, and neither tractions nor displacements are prescribed, respectively. Regions I~, belong to S' and I~' belong to S", such that integrals over the interfacial plane in eqn (2) only
148
Eric R. Johnson, Brian L. Kemp j
6X
Ex
F g . 1. Laminate configuration and nonzero stresses.
involve portion I~". In regions I;," of the interfacial planes continuity conditions are enforced. The vanishing of the volume integrals in eqn (2) leads to constitutive equations and equilibrium equations of linear elasticity. Vanishing of the surface integrals on S requires that one term in each of the products zl u t, "/'2U2, and r3 u3 be prescribed at each point on S. T h e general formulation given above is specialized to the problem of uniform axial extension of a symmetric orthotropic laminate, which is shown in Fig. 1. The length of the laminate is 2L, the width is 2b, and its thickness is 2H. With the origin of a right-handed Cartesian coordinate system (x,y,z) at the center of the laminate, the x-, y-, and z- axes are taken parallel to the length, width, and thickness, respectively. Each layer is considered to be orthotropic with material axes parallel to the Cartesian coordinate axes. Consequently the strain-stress equations for a layer are
Ey Ez
Yyz
=
k ' Y xyJ
where
Sij(i,j
=
$21
$22
531
$32
S33
0
0
0
S~
0
0
0
0
,2 . . . . .
SYM.
t
O"y or z
Cryz
0
866
6) are the compliances.
O'xy
(3)
Modeling the stressfield in laminatedplates
149
Let u, v, and w be the displacement components in the x-, y-, and zdirections, respectively. The laminate is subjected to a uniform axial strain ~x by prescribed displacements u = ---~L atx --- - L . The laminate is assumed sufficiently long such that any end effects associated with the applied end surface displacements are negligible in the x = 0 plane. In addition the laminate may be subjected to prescribed surface tractions which are spatially uniform in x and do not have an x-direction c o m p o n e n t (~-~= 0). In this situation the stresses and strains are i n d e p e n d e n t of x, and the nonzero stresses components are o-~, o-~, O-z,and o'~ (see Fig. 1). Within this class of boundary value problems in laminate elasticity are the free edge tensile coupon problem and the problem of cylindrical bending. For this class of problems eqn (2) becomes
2L
{ [Ex-(sllO'x t-sl2Ory-sl3Orz](kjaO'(xk}
~,
+ [v,y - (S,20-x+ S,,o'y + S2~o':)]~kl~-~)~
+ [w,z- (S,30"x+ S,~rv + $330"~)]~k~~o-I:kl 4- [ V , z 4-
W , y - - r,,3 440ryz]l(k) O0"vz e, (k)
- - [Ory.y -~
O'zy.z](kJ~v (k)
[o-vz.y+ o'~.z]Sw(~ } dAk
+2L f
÷z)w]dC
[(ry-~'r)~Sv+ C'
[&ry(v-O)+brz(w-Cv)]dC
+2L I
(4)
C"
N
I
+b
[Uzy ov (k)__ -ro'~k)~w(k) -__(k)~.
+2Lk~= I I
(k+l)~ OV (k+l)
O'zy
b
_ O.(zk+l~ aW~k+ ,i]j~,,dy =
0
In eqn (4) A k denotes the area of the kth layer in the x = 0 plane, and C is the edge curve of the laminate area in this plane.
150
Eric R. Johnson, Brian L. Kemp
3 STRUCTURAL MODEL The structural theory follows from explicit assumptions on the stresses in thickness coordinate. Let h represent the thickness of the layers, and let be a dimensionless thickness coordinate in a typical layer such that (k~ - 1 -< g -< I. It is assumed that the in-plane normal stress components O-x (k) ana or,, within each layer vary linearly in ~ by the familiar expressions --
o.Ik ) __
N!kl ''t
6M~k~
h
(5)
h 2 ~ (i = x,y)
in which the layer force and moment resultants are (N,, M o~l~l =
o-iIk~ 1, - ~ ~
(i = x,y)
(6)
I
In the following the superscript (k) indicating a layer variable will be omitted unless it is needed for clarity. To obtain a stress field which is statically admissible the stress components must satisfy the elasticity equilibrium equations. These are
The assumed stress component o-y (5) is substituted into eqn (7) and indefinite integration on ~ is performed. This results in an expression for the shear stress o-~ywhich is quadratic in ~, contains an arbitrary constant, and also contains derivatives of the resultants (6). To eliminate these derivatives, eqn (7) is first integrated with respect to Z2from - 1 to + 1, and second is multiplied by ~ and integrated from - 1 to + 1 again. The results of this procedure are the layer equilibrium equations
Ny,y+s2-s1 ---
0
h My.y+ V y - 5 ( s l + s 2 ) = 0
(9) (10)
Modelingthestressfield in laminatedplates
151
where sl = tr~yat ~ = - 1
(11)
s2 = o-~yat ~ = +1
(12)
and V r = (h/2)
j-+'~rzyd~
(13)
-1
U s i n g eqns (9)-(13) the expression for o-~yo b t a i n e d from indefinite integration of (7) can be s h o w n to be °'~r = 1 (sl + Sz) (3~ 2 - 1) + This e q u a t i o n for the shear stress is substituted into eqn (8) and indefinite integration on ~ of the resulting equation gives o'~ as a cubic polynomial in ~. T h e derivative of the shear resultant V r appears in this expression for trz as well as a second arbitrary constant. The derivative of V r is eliminated by the third layer equilibrium equation obtained by integrating eqn (8) with respect to g from - 1 to + 1 and using eqn (13). This third e q u a t i o n is
Vy,y+p2-pl
= 0
(15)
pl = o'zat ~ = - 1
(16)
p2 = o-~at ~ = +1
(17)
where
Using eqns (15)-(17) the final form for o'z is trz = ~ ( p , + p 2 ) +
( p 2 - p l ) ( ~ - ~3/3)
+ (~ )(s2-sO,y(1-~2)+ (~)(s,+sz),y(~-~ 3)
(18)
E q u a t i o n s (5), (14), and (18) are the assumptions for the stress comp o n e n t s o-i (i -- x, y). o'~y, and trz, respectively. These constitute a
152
Eric R. Johnson, Brian L. Kemp
statically admissible field if layer equilibrium eqns (9), (10) and (15) are satisfied. The assumptions for o'x, o-,,, and o-..vare identical to Pagano's ~ assumptions. However, the assumption for o-: (18) is different from Pagano's. In eqn (18) derivatives of the interlaminar shear stress comp o n e n t s in y appear. Pagano's assumption does not contain derivatives of the interlaminar shear stresses. Instead, his expression contains two mathematical stress quantities N: and M. which are given by definition (6) if i is set equal to z. Assumptions (5), (14), and (18) for the stress components are substituted into the variational principle (4). Explicit integration on /~ is carried out. Although no assumptions on the displacements have been m a d e , weighted integrals in ~;on the displacements occur as a result of this explicit integration. Weighted integrals of a displacement function f are defined by (f,?,f) =
If'
f ( l , {, ~Z)d(
(lV)
Also in this integration process displacements evaluated at the bottom and top of a layer occur. Let vl k~ and wl kl denote displacements at the b o t t o m of the kth layer, and v~kl and w!k~denote the displacements at the top. After completing the explicit integration on ~, eqn (4) reduces to a single integral in y. In this reduced form of eqn (4), integration by parts is p e r f o r m e d on derivatives (in y) of the variations in the interlaminar shear stresses. These derivatives originate from the 8cr~ ~term in eqn (4). Also in the process of manipulating the variational principle, it is assumed that the prescribed tractions at y = -+b have the same functional dependence on ~ as do o-,. and or:,., and interlaminar traction continuity is assumed. The resulting structural model is as follows.
3.1 Layer displacement--stress equations he~ = S I,N~ + S,eN~. + SI3N:
0 = -SItM~h
S12Mv + SI3M:
g ~,,. = S p N ~ + S,~2Nv + $23N:
(20) (21) (22)
Modeling the stressfield in laminatedplates
153
h2, --~v,y = -Si2Mx- S22My + 823Mz
(23)
2v + (h/2) ( ~ - ff),y = S~[(4/5) V , - (h/lS)(s, + s2)]
(24)
In eqns (20)-(23) the quantities Nz and Mz are defined by N: = (h/2)(p, +P2) + (h2/12)(s2-sl),r
(25)
Mz = (h2/lO)(p2- p,) + (h3/lZO)(s,
(26)
q'- S2),y
3.2 Layer equilibrium = 0
(27)
My,y+ Vv- (h/2)(s, + s2) = 0
(28)
Vyy+p2-pl =0
(29)
Ny,y+S2-Si
3.3 Interlaminar continuity (k = 1,2 . . . . .
N - 1)
s?' = s? +'~
(30)
p?) = plk+l)
(31)
+ 3~ )
-
h2
(s, - 4s2) ] h2
] (k)
+ -~ S13X2.y + ]~$23 Y2.y+ h~3:~Z2.y -
( v,-~(7-3~')-lSaa[Vr+~(s2-4sO] h2
h2
2 q-'-~S13Xl,y-~-~S23Yl,r-l- h S33Zl,y
-~w~ -
J
~ (~, - ~,) - s , ~ 2 ~ -
s~'~-
{2 3 -~w,-~(~-~)+s,~Y(,+
I (k+l)
(32)
= 0
s~2~
s~',+s~2,
} (k+,)
=o
(33)
Eric R. J o h n s o n ,
154
Brian
L. Kemp
Equations (32) and (33) result from interlaminar traction continuity. In these equations it is consistent with the variational principle (4) to speci~ the interlaminar displacement continuity conditions v~k ' = vl k+''
(34)
w~ k) = wll k+tt
(k = 1,2
.....
(35)
1)
N-
Thus displacements at the interfaces do not appear explicitly in the continuity conditions (32) and (33). Also in eqns (32) and (33) the following stress quantities at the bottom and top of a layer are defined
Xl ] _ Nx
YI ] Y21 -
6 Mx
(36)
Ny 6 M~ h +-5 h i
{Z1} = [11/210 Z2 [13/420
(37)
13/420-t//105 11/210 -h/140
h/140]
t:it
h/105]
~S,,v[
(39)
= h __--h-5 T ~'1 ]
f'2
J
N.,,
(38)
12 My
(40)
Dh + - -5h ~ 9
fg2 = fO
26
- 13h/6
llh/3]
I Sl.y ~S2.y
(41)
Modeling the stress field in laminated plates
155
3.4 Boundary conditions A t the b o t t o m surface (z = - / 4 ) and the top surface of the laminate (z = + H) one term in the following products is prescribed as a function of Y
SI'IVl 1), plt)wl '1, S~V~ N), p~mw[~; --b < y < b
(42)
If the displacements are prescribed, then they are related to the solution variables by v, = ( v - 3~)/2 + (I/lO)S44[Vv + (h/3)(s2- 4s~)]
- (h2/12)[S13X, + S23YI + 12S33Z,],v
(43)
v2 = (V+ 3~)/2 - (I/IO)S~4[Vv+ (h/3)(s~- 4s2)] - (h2/12)[S,3X2 + S,.3Y2 + 12S33Z2].y
(44)
2 ~w, = 3(~ - v~)/(2h) - S,3-~',- $2~ Y, - S3~Z.,
(45)
2 ~w2 = 3(~ - ff)/(2h) + S~3X2 + $23~'2 + $33Z...
(46)
A t the edges y = _+b one term in the following products is prescribed for each layer
NyV, MRS,, Vy(ff~- if)
(47)
A t the b o t t o m , at laminae interfaces, and at the top of the laminate one t e r m in the following products is prescribed at the edges y -- +-b. s l " { ~ + (~ - 3ff)/2}
(48)
s~k){ - [w + (~ - 3ff)/2] t~+ ~)+ [~ - (~ - 3ff)/2] tk'} (k = 1,2 . . . . . N - 1)(49) s~m { ~, - (W - 3if)/2 }tin
(50)
Eric R. Johnson, Brian L. Kemp
156
If the displacement quantities at y -- ± b are prescribed in the edge conditions (48)-(50), then they are related to the solution variables at y = ±bby (3/h)[~ + ( ~ - 3ff)/2] ~1 = [S,3X, + $23Y,
+
12S33Z1] Ill
(51)
- (3/h)[~ + (~ - 3ff)/2] 'k +" + (3/h)[~, - (ff - 3v~)/2]~k' =
_[Si3S -~-
I + $23Y I +
12S33ZI](k
~1)
[513X2 ~t_$23 Y2 + 12S33Z2] Ik) (k = 1,2 . . . . .
N - 1)
(3/h)[~, - (~ - 3ff)/2] ~u~ = [S,3X2 + $23Y2 + 12S33Z2]~m
(52)
(53)
3.5 Summary of model An accounting of the dependent variables and governing equations is shown in Table 1. There are 12N variables and equations. At y = - b there are 4 N + 1 edge conditions; 3N from eqns (47) and N + 1 from eqns (48)-(50). For this class of boundary value problems Pagano's general theory reduces to 16N dependent variables and equations, and at y = - b there would be 5N edge conditions. The reason the present model has 4 N variables less than Pagano's theory is that the additional stress quantities Nz and Mz are not introduced into o-~. The difference in the number of edge conditions at y = ± b is due to Pagano's more general formulation. H e does not initially assume continuity of stresses at the interfaces, whereas the derivation of the present model does. Consequently the continuity of the interlaminar shear stress at y = ± b is equivalent to N - 1 edge conditions. Adding these N - 1 edge conditions to the 4N + 1 edge conditions which are explicitly stated equal the 5N edge conditions of Pagano's theory.
4 SOLUTION A numerical solution to the mathematical model developed in the previous section was obtained using the computer code P A S V A R T . 2This code can solve a system of first order linear or nonlinear ordinary differ-
Total
Layer resultants Layer displacements Interlaminar tractions Upper and lower surface tractions ~, ~, ~ - ff P h Sh p2, S2 PII), st1), p~U), s~N)
Nx, Mx, Ny, My, Vy
Variables
12N
5N 3N 4 ( N - 1) 4
Number
TABLE 1 Summary of the Mathematical Model
Constitutive (20)-(24) Equilibrium (27)--(29) Continuity (30)-(35) Prescribed tractions or displacements (42)-(46)
Equations
12N
4
5N 3N 4 ( N - 1)
Number
Eric R. Johnson, Brian L. Kemp
158
ential equations. The method of solution is by finite differences with variable step size. To utilize this solver the field equations in Sections 3.1 to 3.3 were put into first order form. This was accomplished by eliminating those variables that appear algebraically. Equations (20) and (21) were used to eliminate N(~k~ and M~k~ (k = 1,2 . . . . . N) and eqns (31) and (33) were used to eliminate pl k~ and p~k~ (k = 1,2 . . . . . N - 1). This reduced the n u m b e r of variables and equations from 12N to 8N + 2. The first order system has the form dY dy
-
JY + Cf
(54)
where J is a (8N + 2)-"coefficient matrix, C is a known 8N + 2 by 9 matrix, f(y) is a 9 by 1 vector function, and the 8N + 2 solution vector is
y(y)V ---s(:+ I.v,
/oI~ ~1,'.'
+
s(i I )
,
s(Z+, +<+~ . . . . .
V(I)~,(i)
-v- v
,
(~_
~)(i)~/(i)
-.'+
(x> , s~,+',i?+ •~ ;+( ~ / ) } _.
/~f(i) V(i) , +--'v' • ++' ,
(55)
The vector f(y) contains the applied strain ~x and the prescribed surface quantities and their derivatives, see eqn (42). The system (54) is solved subject to the 4N + 1 edge conditions at y = b and 4N + 1 conditions at y = - b , see eqns (47)-(50).
5 RESULTS Numerical results are presented for the [0/90]~ laminate subject to uniform axial extension, traction free surfaces at z = +-H, and traction free edges at y - ~ _+b. This free-edge problem has been extensively studied by many researchers. Salamon 3surveys most of this work prior to 1980. Recent work significant to this study has been done by Pagano and Soni +and Wang and Choi. ' Pagano and Soni have incorporated Pagano's t original theory into a global-local modeling scheme in order to analyze laminates of many layers, and hence avoid the computer overflow/ underflow limitations encountered by Pagano' when the number of layers exceeded six. Using linear elasticity theory and assuming each lamina is homogeneous, Wang and Choi have determined the strength of the
Modeling the stress field in laminated plates
159
singularity of the interlaminar normal stress o'~ at the free edge between the 0 ° ply and 90 ° ply of a graphite--epoxy laminate. Pagano's ~theory and the simplification of it presented herein do not exhibit this singularity. H o w e v e r , the finite interlaminar normal stress at the free edge grows monotonically with an increasing number (N) of mathematical subdivisions through the thickness. Since the purpose of this paper is to simplify Pagano's original theory, the results presented are compared to his results only. Because of symmetry the quarter section 0 <- z <- H and 0 <- y <- b in the x = 0 plane is analyzed. The lower surface (z = 0) and upper surface (z -= H) conditions prescribed are slt) = wl,) = s~U) = p~m = 0,
0
(56)
The conditions prescribed at the edge y = 0 are -~lk) = ~ k ) =
v~k) = 0
(k = 1 , 2 . . . .
N)
(57)
and sl 1' = s~ 2 ) =
...
= s~ N~ = o
(58)
The free-edge conditions prescribed at y = b are N(k) = ill(k) =. V(vk) = 0 y ---y .
(k = 1,2 . . . . . N)
(59)
and sl '1 = s~ 2 1 =
...
= s~ ~ = 0
(60)
The material properties are taken to be those used by Pagano. ~They are Ell = 138 GPa, Eez
=
E33
=
14-5 GPa
Gi2 = G13 = G23 = 5.86GPa 1)12 = 1.'13= 1/23~-- 0.21 where subscripts 1, 2, and 3 refer to the fiber, transverse, and thickness directions, respectively. Figures 2-4 show the results of the present model for N = 2 and N = 6. The solid lines are the solutions for N = 2 and the dashed lines are the
160
Eric R. Johnson, Brian L. Kemp
20
- - - - - - T
. . . . . . .
i
F
- "•
T
z ©o
16 o
12
90 °
0
b=4H
0 08
S
04
TH/2 ~ y
.....
N 2
Ref I o
6
0
[ j
v
C
~ y / f t
0
-0.4~__~ 0
, 02
__
04
,
~
06
08
ty ";'
J 10
y/b Fig. 2. M i d p l a n e normal stress distributions for N = 2 and N = 6.
20 Z J6
o°
no
90 ° 1.2
0
0 0 . 8
t
-
b:4H N 2 .... 6
2 b Ref I o 0 !
04b" 0
-0.4~ 0
0.2
04
O6
08
0
y/b Fig. 3. Interlaminar normal stress distributions for N = 2 and N = 6
Modeling the stress field in laminated plates
161
]
2.0 Z
0° 19_
L
90 ° 0
O
~-y
b:4H N 2
12
v . . . .
6
b Ref I
t
/
o
D
08 T >. N
b
04.
I
00
0'2
0 4
016
018
ylb
Fig. 4. I n t e r l a m i n a r s h e a r stress distributions for N = 2 a n d N = 6.
solutions for N = 6. For N = 2 each lamina is treated as a single layer, and for N = 6 each lamina is subdivided into three layers. Increasing values of N result in a more refined analysis. The circles in Figs 2-4 represent data points from Pagano's solution with N = 2, and the squares represent his solutions for N = 6. The interlaminar normal stress along the midplane (z = 0) is shown in Fig. 2, the interlaminar normal stress along the 0/90 interface in Fig. 3, and the interlaminar shear stress along the 0/90 interface is shown in Fig. 4. These results demonstrate that the present model can duplicate Pagano's solution. Also Table 2 presents numerical values for the interlaminar normal stresses at the free edge with increasing values of N for the present model and Pagano's theory. It is apparent from this table that the 0/90 interface value of the interlaminar normal stress at the free edge is increasing monotonically with increasing N. Also the free-edge value of the interlaminar normal stress at the midplane is apparently increasing asymptotically to a finite value with increasing N. In their study on the suitability of the displacement-formulated finite element method to compute interlaminar stresses for the free-edge problem, Raju, et a l . , 6 point out that the stress tensor may not be symmetric at a stress discontinuity or singularity. As stated earlier, Wang
Eric R. Johnson, Brian L. Kemp
162
TABLE 2 l n t e r l a m i n a r N o r m a l S t r e s s e s at t h e F r e e E d g e w i t h I n c r e a s i n g N
N
~rz(b, O)/(~xlO +6) (kPa)
o-=(b, H/2)/(Ex 10 +6) (kPa)
Present model
Pagano' s model
Present model
Pagano's model
2 4
1.316 1.872
1-317 --
1.453 1-834
1-453 --
6 8
1.974 1-996
1.975 --
2.057 2.220
2.057 --
10
2.001
--
2.347
--
and Choi 5 have recently determined the strength in the singularity in cr~ at the 0/90 interface and free edge for the [0/90]~ graphite--epoxy laminate. Thus it is very likely that o-~y4: o-vz at (y, z) = (b, H/2) for the linear elasticity solution. This situation can be approximated by the present model even though the stress tensor is symmetric in the interior and at the boundary. Instead of prescribing s~k~to vanish at the 0/90 interface and y = b in the edge conditions (49), we prescribe for the layers adjacent to the 0/90 interface and at y = b - [~ + ( ~ - 3 f f ) / 2 ] (k+')+[~ - ( # - 3 f f ) / 2 ] (k) = 0
(61)
This condition is related to the dependent variables in the layers adjacent to the 0/90 interface by eqn (52). This set of boundary conditions is called Case 2, whereas the boundary conditions given in eqns (56)-(60) are called Case 1. The only difference between the two cases is that the shear stress vanishes at the 0/90 interface at y -- b in Case 1, whereas it does not necessarily vanish in Case 2. It is recognized that prescribing the layer displacements adjacent to the 0/90 interface at y = b as given by eqn (61) is arbitrary, since the actual free-edge displacement w(b, z) is unknown. Prescribing condition (61) places a restriction on the displacement response between these two adjacent layers. However, this restriction represents conditions on the displacements beyond the degree to which they are represented in the interior. That is, the thickness displacement w in the interior of a layer is known only to the extent that the integrated displacement quantity ~ - ~ is known. The integrated displacements and ~ - 3~ appearing in eqn (61) do not occur as dependent variables in the model.
Modeling the stress field in laminated plates 5.0
"~
l
I
I
I
163
I
I
2.5 O°
(3
o..
I
/
90°
20 0
o
,
b=4H
I ! I
b
N=IO 1.5
CASE I
f
,ff B.C.
//1 J
J
i
£,J
t
1.0 N
b I
0.5
'k,
0
I
075
I
080
[
I
I
085
090
0.95
I 0
y/b
Fig. 5. Interlaminar shear stress distributions for Case 1 and Case 2 free-edge interlaminar boundary conditions.
Figures 5 and 6 show the interlaminar shear and normal stresses, respectively, along the 0/90 interface for Cases 1 and 2 and for N = 10. The results for Case 1 are given by the solid line and the results for Case 2 are given by the dashed line. As shown in Fig. 5 the interlaminar shear stresses in the interior are the same for both cases, but they differ significantly near the free edge. In Case 1 the shear stress distribution must adjust rapidly to satisfy the prescribed value (zero) at the free edge, whereas in Case 2 the shear stress increases monotonically to a maximum value at the free edge. Figure 6 shows that the interlaminar normal stress varies only slightly between the two cases in the vicinity of the free edge. The normal stress at the free edge is greater in Case 2 than in Case 1. Figure 7 is a plot of the shear stress O-yzat y = b in Case 2 as a function of z for N = 10. It is zero everywhere except in the two layers adjacent to the interface. (However, in the two layers adjacent to the 0/90 interface o-y: has an average value of zero since Vy is prescribed to vanish at y = b.) The maximum magnitudes for the interlaminar shear stresses at the 0/90 interface are given in Table 3 for Cases 1 and 2 and for various values of N. In Case 1 the location of the maximum is also given. In Case 2 the m a x i m u m always occurs at y = b. Maximum shear stresses are signifi-
164
Eric R. Johnson, Brian L. Kemp
,z 0
0
I 5
. . . .
7
9
b
b:4H
I
,,
/
t:!O I 0
-
-
/ CASE
I B C
I
/ /
I
0
~
l
0
I
075
~ ....
080
'085
I
....
090
J.J
±
I0
095
y/b
Fig. 6. Interlaminar normal stress distributions for Case 1 and Case 2 free-edge interlaminar boundary conditions.
i0
~
•
~
I
I
1]
~z
0 °
] 1 , O- Z
075 0
z/H
b=4H N:IO
b
±
I __i
0 5
i 025t-
o I-i
I
3o
J--
I
2o
1
I
-Io
1
I
o
J
Io
210
30
O-yz(b,z)/(ExlO 8) kPa Fig. 7. Free-edge shear stress distribution for Case 2 interlaminar boundary condition.
165
Modeling the stress field in laminated plates TABLE 3
Maximum Magnitudes of the Interlaminar Shear Stress for Cases 1 and 2 and for Various N N
2 4 6 8 10
- ~rzy(y, HI2) I(~x 10 +6) (kPa) Case 1 a t y / b
Case2aty/b = I
1-339 at 0.9607 1-500 at 0.9792 1.583 at 0-9861 1.635 at 0.9883 1.671 at 0.9905
2.026 2.446 2.586 2.667 2.721
TABLE 4
Maximum Interlaminar Normal Stress for Cases 1 and 2 and for Various N N
2 4 6 8 10
orz (b, H/2)/(~xlO +6) (kPa) Case 1
Case 2
1.453 1.834 2-058 2.220 2.347
1.451 1.879 2.100 2-262 2.389
cantly greater in Case 2 than for Case 1, and the differences are increasing with increasing N. However, the rate of increase diminishes with increasing N. Table 4 gives the maximum values of the interlaminar normal stress for Cases 1 and 2 and for various values of N. The differences are small, and remain nearly constant for N > 2.
6 CONCLUDING REMARKS The purpose of this paper is to present a simplified version of Pagano's ~ theory for the analysis of interlaminar stresses in laminated composite plates. It has been shown that the present model can duplicate Pagano's results for the [0/90]s laminate subjected to uniform axial extension. The model developed in this paper, however, has four dependent variables
166
Eric R. Johnson, Brian L. Kemp
less per layer than Pagano's theory. This reduction in the number of d e p e n d e n t variables is significant and should reduce computational costs. Two boundary conditions were examined for the free edge problem. In Case 1 the interlaminar shear stress at the free edge of the 0/90 interface was prescribed to vanish. In Case 2 this free edge shear stress was not prescribed and the alternate condition (61) was prescribed to vanish. Thus in Case 2 the interlaminar shear stress at the free edge is not necessarily zero. Case 2 is an attempt to model the singular behavior present in the linear elasticity solution at the free edge. It was found that the maximum interlaminar shear stresses are significantly greater for Case 2 boundary conditions than for Case 1 (see Table 3). For the type of model developed here, the interlaminar normal stresses at the free edge are finite and do not differ by very much between the two boundary condition cases (see Table 4).
ACKNOWLEDGEMENT Research was sponsored by the Air Force Office of Scientific Research, Air Force Systems Command, USAF, under Contract F49620-82-C-0035, and under Grant AFOSR-83-0191. The US Government is authorized to reproduce and distribute reprints for Governmental purposes notwithstanding any copyright notation thereon.
REFERENCES 1. Pagano, N. J., Stress fields in composite laminates, Int. J. Solids Structures, 14 (1978) 385-400. 2. Lentini, M. and Pereyra, V., An adaptive finite difference solver for nonlinear two point boundary value problems with mild boundary layers, S I A M J. Numerical Analysis, 14 (1977) 91-111. 3. Salamon, Nicholas J., An assessment of the interlaminar stress problem in laminated composites, J. Composite Matl. Supplement, 14 (1980) 177-94. 4. Pagano, N. J. and Soni, S. R., Global-local laminate variational model, Int. J. Solids Structures, 19 (1983) 207-28. 5. Wang, S. S. and Choi, I., Boundary-layer effects in composite laminates: Part I--Free edge stress singularities, J. Appl. Mech., 49 (1982) 541-8. 6. Raju, I. S., Whitcomb, J. D. and Goree, J. G., A new look at numerical analyses of free-edge stresses in composite laminates, NASA Tech. Paper 1751 (1980).
Modeling the Stress Field in Laminated Composite Plates Near Discontinuities
Eric R. Johnson and Brian L. K e m p Aerospace and Ocean Engineering Department, Virginia Polytechnic Institute and State University, Blacksburg, Virginia 24061, USA
ABSTRACT A structural model is developed using Reissner's variational principle to analyze stress gradient problems in laminated composite plates. A statically admissible stress field, which has explicit dependence on the thickness coordinate, is selected for each layer. Both interlaminar traction and displacement continuity are enforced between layers. This model is similar to an earlier one developed by Pagano, but has four dependent variables less per layer than Pagano's model. It is shown that the model can duplicate Pagano's results for the same number of subdivisions through the thickness. The uniform axial extension of a [0/90Islaminate is used as the example problem for comparison to Pagano' s theory and to illustrate the methodology.
1 INTRODUCTION Failure prediction capability in laminated composite structures depends upon a reasonably accurate analysis of the stress response. For example, to predict the initiation of delamination, or the stress redistribution in the vicinity of damage, suggests that the theory used to determine stresses includes all stress components. Continuity of displacements and stresses at interfaces between laminae should also be required of the theory. Elasticity theory would satisfy these criteria, but currently threedimensional elasticity analysis is intractable for practical laminated composite structures. On the other hand, structural theories for 145 Composite Structures 0263-8223/85/$03-30 O Elsevier Applied Science Publishers Ltd, England, 1985. Printed in Great Britain
146
Eric R. Johnson, Brian L. Kemp
laminates which are based on assumptions for the displacements in the thickness coordinate through the entire laminate do not satisfy these criteria, but are more amenable to solution. In general, structural theories based on assumed displacement fields (this includes finite element formulations of this type) cannot satisfy continuity of stresses along laminae interfaces. In particular, classical structural theory based on the Kirchhoff-Love hypotheses neglect thickness stress components entirely. The deficiency of displacement-based structural theories for the analysis of stress concentration problems in laminated composites is discussed in some detail by Pagano.' In his paper, Pagano establishes the criteria for a new laminate structural theory; i.e. all stress components are nonzero, continuity of displacements and stresses is satisfied at interfaces, and that individual layers satisfy force and moment equilibrium. Pagano develops a structural theory which satisfies these criteria and which is based on assumed stress fields within each layer. Reissner's variational principle is employed to obtain a consistent mathematical formulation. The purpose of this paper is to present a simplification of Pagano's original theory. Pagano's theory has 23N dependent variables, where N is the number of layers. It is shown that with a different assumption on the thickness normal stress (o-:) in each layer the number of dependent variables in the formulation can be reduced by four per layer. Thus, in the general case, the model developed herein would have 19N dependent variables. It is also shown that this simpler model can duplicate Pagano's results for the same number of layers. This is illustrated by solving the problem of a [0/90]~ laminate subject to uniform axial extension.
2 L A M I N A T E V A R I A T I O N A L PRINCIPLE The mechanism for deriving the structural theory is Reissner's variational principle. The representation of this principle for laminated bodies of the type considered in this paper is given by Pagano. ~ The body is a flat rectangular plate consisting of an integral number (N) of layers which are flat rectangular sheets having identical dimensions. Each layer is assumed to be homogeneous, anisotropic, and linear elastic. Also displacements and strains are assumed to be small when the laminate is subjected to prescribed surface tractions and/or displacements. Beginning at one of the two surfaces of the body which are parallel to the interfacial planes,
Modelingthestressfieldinlaminatedplates
147
the layers are n u m b e r e d consecutively from one to N. For the kth layer (k = 1, 2, • . . , N) let Vk designate its volume, ~_~k~ j (i, j = 1,2,3) denote the stresses, u]k~ the displacements, and let W ~k~ represent the complementary strain energy density. Also let S' represent the portion of the laminate's boundary surface area S on which tractions are prescribed, and let S " represent the portion of S on which displacements are prescribed. The traction components ri acting on a surface are related to the stresses by ri = njo-j~, in which nj are the components of the unit outward normal vector to the surface, and repeated indices are summed from one to three. Prescribed surface tractions are denoted by ~-~, and prescribed displacements by fi~. Let a comma designate partial differentiation. Neglecting body forces, the Reissner functional is
J = k~ ' =
uj.i) - W
dV k
!/k
-- f s, ~'iuidS- ~s,, "t"i(l'li--(ti)dS
(1)
The equations of laminate elasticity follow from the vanishing of the first variation of J, which is written as 8J = 0. Both stresses and displacements are varied in this process. After applying Gauss's integral t h e o r e m and using symmetry of the stress tensor, the results of this procedure are
k~=,f vk { [2 (uiJ+ujJ)- ~] 80"~i--o'J"'SU'} 'k)dV k + L, ('~-r,)Su~dS- ~, ~r,(u,- fi,)dS +
[r}k)~Ulkl+'~k+l)SUlk+l)]dI~ = 0 k = I
(2)
Ik'"
In eqn (2) Ik designates the area of the interracial plane c o m m o n to layers k and k + 1, ,?/k) are the interracial traction components acting on the upper surface of the kth layer, and ~.]k+~) are the interracial traction c o m p o n e n t s acting on the lower surface of layer k + 1. Let Ik,' It," and It'", designate the portions of the interracial planes on which tractions are prescribed, displacements are prescribed, and neither tractions nor displacements are prescribed, respectively. Regions I~, belong to S' and I~' belong to S", such that integrals over the interfacial plane in eqn (2) only
148
Eric R. Johnson, Brian L. Kemp j
6X
Ex
F g . 1. Laminate configuration and nonzero stresses.
involve portion I~". In regions I;," of the interfacial planes continuity conditions are enforced. The vanishing of the volume integrals in eqn (2) leads to constitutive equations and equilibrium equations of linear elasticity. Vanishing of the surface integrals on S requires that one term in each of the products zl u t, "/'2U2, and r3 u3 be prescribed at each point on S. T h e general formulation given above is specialized to the problem of uniform axial extension of a symmetric orthotropic laminate, which is shown in Fig. 1. The length of the laminate is 2L, the width is 2b, and its thickness is 2H. With the origin of a right-handed Cartesian coordinate system (x,y,z) at the center of the laminate, the x-, y-, and z- axes are taken parallel to the length, width, and thickness, respectively. Each layer is considered to be orthotropic with material axes parallel to the Cartesian coordinate axes. Consequently the strain-stress equations for a layer are
Ey Ez
Yyz
=
k ' Y xyJ
where
Sij(i,j
=
$21
$22
531
$32
S33
0
0
0
S~
0
0
0
0
,2 . . . . .
SYM.
t
O"y or z
Cryz
0
866
6) are the compliances.
O'xy
(3)
Modeling the stressfield in laminatedplates
149
Let u, v, and w be the displacement components in the x-, y-, and zdirections, respectively. The laminate is subjected to a uniform axial strain ~x by prescribed displacements u = ---~L atx --- - L . The laminate is assumed sufficiently long such that any end effects associated with the applied end surface displacements are negligible in the x = 0 plane. In addition the laminate may be subjected to prescribed surface tractions which are spatially uniform in x and do not have an x-direction c o m p o n e n t (~-~= 0). In this situation the stresses and strains are i n d e p e n d e n t of x, and the nonzero stresses components are o-~, o-~, O-z,and o'~ (see Fig. 1). Within this class of boundary value problems in laminate elasticity are the free edge tensile coupon problem and the problem of cylindrical bending. For this class of problems eqn (2) becomes
2L
{ [Ex-(sllO'x t-sl2Ory-sl3Orz](kjaO'(xk}
~,
+ [v,y - (S,20-x+ S,,o'y + S2~o':)]~kl~-~)~
+ [w,z- (S,30"x+ S,~rv + $330"~)]~k~~o-I:kl 4- [ V , z 4-
W , y - - r,,3 440ryz]l(k) O0"vz e, (k)
- - [Ory.y -~
O'zy.z](kJ~v (k)
[o-vz.y+ o'~.z]Sw(~ } dAk
+2L f
÷z)w]dC
[(ry-~'r)~Sv+ C'
[&ry(v-O)+brz(w-Cv)]dC
+2L I
(4)
C"
N
I
+b
[Uzy ov (k)__ -ro'~k)~w(k) -__(k)~.
+2Lk~= I I
(k+l)~ OV (k+l)
O'zy
b
_ O.(zk+l~ aW~k+ ,i]j~,,dy =
0
In eqn (4) A k denotes the area of the kth layer in the x = 0 plane, and C is the edge curve of the laminate area in this plane.
150
Eric R. Johnson, Brian L. Kemp
3 STRUCTURAL MODEL The structural theory follows from explicit assumptions on the stresses in thickness coordinate. Let h represent the thickness of the layers, and let be a dimensionless thickness coordinate in a typical layer such that (k~ - 1 -< g -< I. It is assumed that the in-plane normal stress components O-x (k) ana or,, within each layer vary linearly in ~ by the familiar expressions --
o.Ik ) __
N!kl ''t
6M~k~
h
(5)
h 2 ~ (i = x,y)
in which the layer force and moment resultants are (N,, M o~l~l =
o-iIk~ 1, - ~ ~
(i = x,y)
(6)
I
In the following the superscript (k) indicating a layer variable will be omitted unless it is needed for clarity. To obtain a stress field which is statically admissible the stress components must satisfy the elasticity equilibrium equations. These are
The assumed stress component o-y (5) is substituted into eqn (7) and indefinite integration on ~ is performed. This results in an expression for the shear stress o-~ywhich is quadratic in ~, contains an arbitrary constant, and also contains derivatives of the resultants (6). To eliminate these derivatives, eqn (7) is first integrated with respect to Z2from - 1 to + 1, and second is multiplied by ~ and integrated from - 1 to + 1 again. The results of this procedure are the layer equilibrium equations
Ny,y+s2-s1 ---
0
h My.y+ V y - 5 ( s l + s 2 ) = 0
(9) (10)
Modelingthestressfield in laminatedplates
151
where sl = tr~yat ~ = - 1
(11)
s2 = o-~yat ~ = +1
(12)
and V r = (h/2)
j-+'~rzyd~
(13)
-1
U s i n g eqns (9)-(13) the expression for o-~yo b t a i n e d from indefinite integration of (7) can be s h o w n to be °'~r = 1 (sl + Sz) (3~ 2 - 1) + This e q u a t i o n for the shear stress is substituted into eqn (8) and indefinite integration on ~ of the resulting equation gives o'~ as a cubic polynomial in ~. T h e derivative of the shear resultant V r appears in this expression for trz as well as a second arbitrary constant. The derivative of V r is eliminated by the third layer equilibrium equation obtained by integrating eqn (8) with respect to g from - 1 to + 1 and using eqn (13). This third e q u a t i o n is
Vy,y+p2-pl
= 0
(15)
pl = o'zat ~ = - 1
(16)
p2 = o-~at ~ = +1
(17)
where
Using eqns (15)-(17) the final form for o'z is trz = ~ ( p , + p 2 ) +
( p 2 - p l ) ( ~ - ~3/3)
+ (~ )(s2-sO,y(1-~2)+ (~)(s,+sz),y(~-~ 3)
(18)
E q u a t i o n s (5), (14), and (18) are the assumptions for the stress comp o n e n t s o-i (i -- x, y). o'~y, and trz, respectively. These constitute a
152
Eric R. Johnson, Brian L. Kemp
statically admissible field if layer equilibrium eqns (9), (10) and (15) are satisfied. The assumptions for o'x, o-,,, and o-..vare identical to Pagano's ~ assumptions. However, the assumption for o-: (18) is different from Pagano's. In eqn (18) derivatives of the interlaminar shear stress comp o n e n t s in y appear. Pagano's assumption does not contain derivatives of the interlaminar shear stresses. Instead, his expression contains two mathematical stress quantities N: and M. which are given by definition (6) if i is set equal to z. Assumptions (5), (14), and (18) for the stress components are substituted into the variational principle (4). Explicit integration on /~ is carried out. Although no assumptions on the displacements have been m a d e , weighted integrals in ~;on the displacements occur as a result of this explicit integration. Weighted integrals of a displacement function f are defined by (f,?,f) =
If'
f ( l , {, ~Z)d(
(lV)
Also in this integration process displacements evaluated at the bottom and top of a layer occur. Let vl k~ and wl kl denote displacements at the b o t t o m of the kth layer, and v~kl and w!k~denote the displacements at the top. After completing the explicit integration on ~, eqn (4) reduces to a single integral in y. In this reduced form of eqn (4), integration by parts is p e r f o r m e d on derivatives (in y) of the variations in the interlaminar shear stresses. These derivatives originate from the 8cr~ ~term in eqn (4). Also in the process of manipulating the variational principle, it is assumed that the prescribed tractions at y = -+b have the same functional dependence on ~ as do o-,. and or:,., and interlaminar traction continuity is assumed. The resulting structural model is as follows.
3.1 Layer displacement--stress equations he~ = S I,N~ + S,eN~. + SI3N:
0 = -SItM~h
S12Mv + SI3M:
g ~,,. = S p N ~ + S,~2Nv + $23N:
(20) (21) (22)
Modeling the stressfield in laminatedplates
153
h2, --~v,y = -Si2Mx- S22My + 823Mz
(23)
2v + (h/2) ( ~ - ff),y = S~[(4/5) V , - (h/lS)(s, + s2)]
(24)
In eqns (20)-(23) the quantities Nz and Mz are defined by N: = (h/2)(p, +P2) + (h2/12)(s2-sl),r
(25)
Mz = (h2/lO)(p2- p,) + (h3/lZO)(s,
(26)
q'- S2),y
3.2 Layer equilibrium = 0
(27)
My,y+ Vv- (h/2)(s, + s2) = 0
(28)
Vyy+p2-pl =0
(29)
Ny,y+S2-Si
3.3 Interlaminar continuity (k = 1,2 . . . . .
N - 1)
s?' = s? +'~
(30)
p?) = plk+l)
(31)
+ 3~ )
-
h2
(s, - 4s2) ] h2
] (k)
+ -~ S13X2.y + ]~$23 Y2.y+ h~3:~Z2.y -
( v,-~(7-3~')-lSaa[Vr+~(s2-4sO] h2
h2
2 q-'-~S13Xl,y-~-~S23Yl,r-l- h S33Zl,y
-~w~ -
J
~ (~, - ~,) - s , ~ 2 ~ -
s~'~-
{2 3 -~w,-~(~-~)+s,~Y(,+
I (k+l)
(32)
= 0
s~2~
s~',+s~2,
} (k+,)
=o
(33)
Eric R. J o h n s o n ,
154
Brian
L. Kemp
Equations (32) and (33) result from interlaminar traction continuity. In these equations it is consistent with the variational principle (4) to speci~ the interlaminar displacement continuity conditions v~k ' = vl k+''
(34)
w~ k) = wll k+tt
(k = 1,2
.....
(35)
1)
N-
Thus displacements at the interfaces do not appear explicitly in the continuity conditions (32) and (33). Also in eqns (32) and (33) the following stress quantities at the bottom and top of a layer are defined
Xl ] _ Nx
YI ] Y21 -
6 Mx
(36)
Ny 6 M~ h +-5 h i
{Z1} = [11/210 Z2 [13/420
(37)
13/420-t//105 11/210 -h/140
h/140]
t:it
h/105]
~S,,v[
(39)
= h __--h-5 T ~'1 ]
f'2
J
N.,,
(38)
12 My
(40)
Dh + - -5h ~ 9
fg2 = fO
26
- 13h/6
llh/3]
I Sl.y ~S2.y
(41)
Modeling the stress field in laminated plates
155
3.4 Boundary conditions A t the b o t t o m surface (z = - / 4 ) and the top surface of the laminate (z = + H) one term in the following products is prescribed as a function of Y
SI'IVl 1), plt)wl '1, S~V~ N), p~mw[~; --b < y < b
(42)
If the displacements are prescribed, then they are related to the solution variables by v, = ( v - 3~)/2 + (I/lO)S44[Vv + (h/3)(s2- 4s~)]
- (h2/12)[S13X, + S23YI + 12S33Z,],v
(43)
v2 = (V+ 3~)/2 - (I/IO)S~4[Vv+ (h/3)(s~- 4s2)] - (h2/12)[S,3X2 + S,.3Y2 + 12S33Z2].y
(44)
2 ~w, = 3(~ - v~)/(2h) - S,3-~',- $2~ Y, - S3~Z.,
(45)
2 ~w2 = 3(~ - ff)/(2h) + S~3X2 + $23~'2 + $33Z...
(46)
A t the edges y = _+b one term in the following products is prescribed for each layer
NyV, MRS,, Vy(ff~- if)
(47)
A t the b o t t o m , at laminae interfaces, and at the top of the laminate one t e r m in the following products is prescribed at the edges y -- +-b. s l " { ~ + (~ - 3ff)/2}
(48)
s~k){ - [w + (~ - 3ff)/2] t~+ ~)+ [~ - (~ - 3ff)/2] tk'} (k = 1,2 . . . . . N - 1)(49) s~m { ~, - (W - 3if)/2 }tin
(50)
Eric R. Johnson, Brian L. Kemp
156
If the displacement quantities at y -- ± b are prescribed in the edge conditions (48)-(50), then they are related to the solution variables at y = ±bby (3/h)[~ + ( ~ - 3ff)/2] ~1 = [S,3X, + $23Y,
+
12S33Z1] Ill
(51)
- (3/h)[~ + (~ - 3ff)/2] 'k +" + (3/h)[~, - (ff - 3v~)/2]~k' =
_[Si3S -~-
I + $23Y I +
12S33ZI](k
~1)
[513X2 ~t_$23 Y2 + 12S33Z2] Ik) (k = 1,2 . . . . .
N - 1)
(3/h)[~, - (~ - 3ff)/2] ~u~ = [S,3X2 + $23Y2 + 12S33Z2]~m
(52)
(53)
3.5 Summary of model An accounting of the dependent variables and governing equations is shown in Table 1. There are 12N variables and equations. At y = - b there are 4 N + 1 edge conditions; 3N from eqns (47) and N + 1 from eqns (48)-(50). For this class of boundary value problems Pagano's general theory reduces to 16N dependent variables and equations, and at y = - b there would be 5N edge conditions. The reason the present model has 4 N variables less than Pagano's theory is that the additional stress quantities Nz and Mz are not introduced into o-~. The difference in the number of edge conditions at y = ± b is due to Pagano's more general formulation. H e does not initially assume continuity of stresses at the interfaces, whereas the derivation of the present model does. Consequently the continuity of the interlaminar shear stress at y = ± b is equivalent to N - 1 edge conditions. Adding these N - 1 edge conditions to the 4N + 1 edge conditions which are explicitly stated equal the 5N edge conditions of Pagano's theory.
4 SOLUTION A numerical solution to the mathematical model developed in the previous section was obtained using the computer code P A S V A R T . 2This code can solve a system of first order linear or nonlinear ordinary differ-
Total
Layer resultants Layer displacements Interlaminar tractions Upper and lower surface tractions ~, ~, ~ - ff P h Sh p2, S2 PII), st1), p~U), s~N)
Nx, Mx, Ny, My, Vy
Variables
12N
5N 3N 4 ( N - 1) 4
Number
TABLE 1 Summary of the Mathematical Model
Constitutive (20)-(24) Equilibrium (27)--(29) Continuity (30)-(35) Prescribed tractions or displacements (42)-(46)
Equations
12N
4
5N 3N 4 ( N - 1)
Number
Eric R. Johnson, Brian L. Kemp
158
ential equations. The method of solution is by finite differences with variable step size. To utilize this solver the field equations in Sections 3.1 to 3.3 were put into first order form. This was accomplished by eliminating those variables that appear algebraically. Equations (20) and (21) were used to eliminate N(~k~ and M~k~ (k = 1,2 . . . . . N) and eqns (31) and (33) were used to eliminate pl k~ and p~k~ (k = 1,2 . . . . . N - 1). This reduced the n u m b e r of variables and equations from 12N to 8N + 2. The first order system has the form dY dy
-
JY + Cf
(54)
where J is a (8N + 2)-"coefficient matrix, C is a known 8N + 2 by 9 matrix, f(y) is a 9 by 1 vector function, and the 8N + 2 solution vector is
y(y)V ---s(:+ I.v,
/oI~ ~1,'.'
+
s(i I )
,
s(Z+, +<+~ . . . . .
V(I)~,(i)
-v- v
,
(~_
~)(i)~/(i)
-.'+
(x> , s~,+',i?+ •~ ;+( ~ / ) } _.
/~f(i) V(i) , +--'v' • ++' ,
(55)
The vector f(y) contains the applied strain ~x and the prescribed surface quantities and their derivatives, see eqn (42). The system (54) is solved subject to the 4N + 1 edge conditions at y = b and 4N + 1 conditions at y = - b , see eqns (47)-(50).
5 RESULTS Numerical results are presented for the [0/90]~ laminate subject to uniform axial extension, traction free surfaces at z = +-H, and traction free edges at y - ~ _+b. This free-edge problem has been extensively studied by many researchers. Salamon 3surveys most of this work prior to 1980. Recent work significant to this study has been done by Pagano and Soni +and Wang and Choi. ' Pagano and Soni have incorporated Pagano's t original theory into a global-local modeling scheme in order to analyze laminates of many layers, and hence avoid the computer overflow/ underflow limitations encountered by Pagano' when the number of layers exceeded six. Using linear elasticity theory and assuming each lamina is homogeneous, Wang and Choi have determined the strength of the
Modeling the stress field in laminated plates
159
singularity of the interlaminar normal stress o'~ at the free edge between the 0 ° ply and 90 ° ply of a graphite--epoxy laminate. Pagano's ~theory and the simplification of it presented herein do not exhibit this singularity. H o w e v e r , the finite interlaminar normal stress at the free edge grows monotonically with an increasing number (N) of mathematical subdivisions through the thickness. Since the purpose of this paper is to simplify Pagano's original theory, the results presented are compared to his results only. Because of symmetry the quarter section 0 <- z <- H and 0 <- y <- b in the x = 0 plane is analyzed. The lower surface (z = 0) and upper surface (z -= H) conditions prescribed are slt) = wl,) = s~U) = p~m = 0,
0
(56)
The conditions prescribed at the edge y = 0 are -~lk) = ~ k ) =
v~k) = 0
(k = 1 , 2 . . . .
N)
(57)
and sl 1' = s~ 2 ) =
...
= s~ N~ = o
(58)
The free-edge conditions prescribed at y = b are N(k) = ill(k) =. V(vk) = 0 y ---y .
(k = 1,2 . . . . . N)
(59)
and sl '1 = s~ 2 1 =
...
= s~ ~ = 0
(60)
The material properties are taken to be those used by Pagano. ~They are Ell = 138 GPa, Eez
=
E33
=
14-5 GPa
Gi2 = G13 = G23 = 5.86GPa 1)12 = 1.'13= 1/23~-- 0.21 where subscripts 1, 2, and 3 refer to the fiber, transverse, and thickness directions, respectively. Figures 2-4 show the results of the present model for N = 2 and N = 6. The solid lines are the solutions for N = 2 and the dashed lines are the
160
Eric R. Johnson, Brian L. Kemp
20
- - - - - - T
. . . . . . .
i
F
- "•
T
z ©o
16 o
12
90 °
0
b=4H
0 08
S
04
TH/2 ~ y
.....
N 2
Ref I o
6
0
[ j
v
C
~ y / f t
0
-0.4~__~ 0
, 02
__
04
,
~
06
08
ty ";'
J 10
y/b Fig. 2. M i d p l a n e normal stress distributions for N = 2 and N = 6.
20 Z J6
o°
no
90 ° 1.2
0
0 0 . 8
t
-
b:4H N 2 .... 6
2 b Ref I o 0 !
04b" 0
-0.4~ 0
0.2
04
O6
08
0
y/b Fig. 3. Interlaminar normal stress distributions for N = 2 and N = 6
Modeling the stress field in laminated plates
161
]
2.0 Z
0° 19_
L
90 ° 0
O
~-y
b:4H N 2
12
v . . . .
6
b Ref I
t
/
o
D
08 T >. N
b
04.
I
00
0'2
0 4
016
018
ylb
Fig. 4. I n t e r l a m i n a r s h e a r stress distributions for N = 2 a n d N = 6.
solutions for N = 6. For N = 2 each lamina is treated as a single layer, and for N = 6 each lamina is subdivided into three layers. Increasing values of N result in a more refined analysis. The circles in Figs 2-4 represent data points from Pagano's solution with N = 2, and the squares represent his solutions for N = 6. The interlaminar normal stress along the midplane (z = 0) is shown in Fig. 2, the interlaminar normal stress along the 0/90 interface in Fig. 3, and the interlaminar shear stress along the 0/90 interface is shown in Fig. 4. These results demonstrate that the present model can duplicate Pagano's solution. Also Table 2 presents numerical values for the interlaminar normal stresses at the free edge with increasing values of N for the present model and Pagano's theory. It is apparent from this table that the 0/90 interface value of the interlaminar normal stress at the free edge is increasing monotonically with increasing N. Also the free-edge value of the interlaminar normal stress at the midplane is apparently increasing asymptotically to a finite value with increasing N. In their study on the suitability of the displacement-formulated finite element method to compute interlaminar stresses for the free-edge problem, Raju, et a l . , 6 point out that the stress tensor may not be symmetric at a stress discontinuity or singularity. As stated earlier, Wang
Eric R. Johnson, Brian L. Kemp
162
TABLE 2 l n t e r l a m i n a r N o r m a l S t r e s s e s at t h e F r e e E d g e w i t h I n c r e a s i n g N
N
~rz(b, O)/(~xlO +6) (kPa)
o-=(b, H/2)/(Ex 10 +6) (kPa)
Present model
Pagano' s model
Present model
Pagano's model
2 4
1.316 1.872
1-317 --
1.453 1-834
1-453 --
6 8
1.974 1-996
1.975 --
2.057 2.220
2.057 --
10
2.001
--
2.347
--
and Choi 5 have recently determined the strength in the singularity in cr~ at the 0/90 interface and free edge for the [0/90]~ graphite--epoxy laminate. Thus it is very likely that o-~y4: o-vz at (y, z) = (b, H/2) for the linear elasticity solution. This situation can be approximated by the present model even though the stress tensor is symmetric in the interior and at the boundary. Instead of prescribing s~k~to vanish at the 0/90 interface and y = b in the edge conditions (49), we prescribe for the layers adjacent to the 0/90 interface and at y = b - [~ + ( ~ - 3 f f ) / 2 ] (k+')+[~ - ( # - 3 f f ) / 2 ] (k) = 0
(61)
This condition is related to the dependent variables in the layers adjacent to the 0/90 interface by eqn (52). This set of boundary conditions is called Case 2, whereas the boundary conditions given in eqns (56)-(60) are called Case 1. The only difference between the two cases is that the shear stress vanishes at the 0/90 interface at y -- b in Case 1, whereas it does not necessarily vanish in Case 2. It is recognized that prescribing the layer displacements adjacent to the 0/90 interface at y = b as given by eqn (61) is arbitrary, since the actual free-edge displacement w(b, z) is unknown. Prescribing condition (61) places a restriction on the displacement response between these two adjacent layers. However, this restriction represents conditions on the displacements beyond the degree to which they are represented in the interior. That is, the thickness displacement w in the interior of a layer is known only to the extent that the integrated displacement quantity ~ - ~ is known. The integrated displacements and ~ - 3~ appearing in eqn (61) do not occur as dependent variables in the model.
Modeling the stress field in laminated plates 5.0
"~
l
I
I
I
163
I
I
2.5 O°
(3
o..
I
/
90°
20 0
o
,
b=4H
I ! I
b
N=IO 1.5
CASE I
f
,ff B.C.
//1 J
J
i
£,J
t
1.0 N
b I
0.5
'k,
0
I
075
I
080
[
I
I
085
090
0.95
I 0
y/b
Fig. 5. Interlaminar shear stress distributions for Case 1 and Case 2 free-edge interlaminar boundary conditions.
Figures 5 and 6 show the interlaminar shear and normal stresses, respectively, along the 0/90 interface for Cases 1 and 2 and for N = 10. The results for Case 1 are given by the solid line and the results for Case 2 are given by the dashed line. As shown in Fig. 5 the interlaminar shear stresses in the interior are the same for both cases, but they differ significantly near the free edge. In Case 1 the shear stress distribution must adjust rapidly to satisfy the prescribed value (zero) at the free edge, whereas in Case 2 the shear stress increases monotonically to a maximum value at the free edge. Figure 6 shows that the interlaminar normal stress varies only slightly between the two cases in the vicinity of the free edge. The normal stress at the free edge is greater in Case 2 than in Case 1. Figure 7 is a plot of the shear stress O-yzat y = b in Case 2 as a function of z for N = 10. It is zero everywhere except in the two layers adjacent to the interface. (However, in the two layers adjacent to the 0/90 interface o-y: has an average value of zero since Vy is prescribed to vanish at y = b.) The maximum magnitudes for the interlaminar shear stresses at the 0/90 interface are given in Table 3 for Cases 1 and 2 and for various values of N. In Case 1 the location of the maximum is also given. In Case 2 the m a x i m u m always occurs at y = b. Maximum shear stresses are signifi-
164
Eric R. Johnson, Brian L. Kemp
,z 0
0
I 5
. . . .
7
9
b
b:4H
I
,,
/
t:!O I 0
-
-
/ CASE
I B C
I
/ /
I
0
~
l
0
I
075
~ ....
080
'085
I
....
090
J.J
±
I0
095
y/b
Fig. 6. Interlaminar normal stress distributions for Case 1 and Case 2 free-edge interlaminar boundary conditions.
i0
~
•
~
I
I
1]
~z
0 °
] 1 , O- Z
075 0
z/H
b=4H N:IO
b
±
I __i
0 5
i 025t-
o I-i
I
3o
J--
I
2o
1
I
-Io
1
I
o
J
Io
210
30
O-yz(b,z)/(ExlO 8) kPa Fig. 7. Free-edge shear stress distribution for Case 2 interlaminar boundary condition.
165
Modeling the stress field in laminated plates TABLE 3
Maximum Magnitudes of the Interlaminar Shear Stress for Cases 1 and 2 and for Various N N
2 4 6 8 10
- ~rzy(y, HI2) I(~x 10 +6) (kPa) Case 1 a t y / b
Case2aty/b = I
1-339 at 0.9607 1-500 at 0.9792 1.583 at 0-9861 1.635 at 0.9883 1.671 at 0.9905
2.026 2.446 2.586 2.667 2.721
TABLE 4
Maximum Interlaminar Normal Stress for Cases 1 and 2 and for Various N N
2 4 6 8 10
orz (b, H/2)/(~xlO +6) (kPa) Case 1
Case 2
1.453 1.834 2-058 2.220 2.347
1.451 1.879 2.100 2-262 2.389
cantly greater in Case 2 than for Case 1, and the differences are increasing with increasing N. However, the rate of increase diminishes with increasing N. Table 4 gives the maximum values of the interlaminar normal stress for Cases 1 and 2 and for various values of N. The differences are small, and remain nearly constant for N > 2.
6 CONCLUDING REMARKS The purpose of this paper is to present a simplified version of Pagano's ~ theory for the analysis of interlaminar stresses in laminated composite plates. It has been shown that the present model can duplicate Pagano's results for the [0/90]s laminate subjected to uniform axial extension. The model developed in this paper, however, has four dependent variables
166
Eric R. Johnson, Brian L. Kemp
less per layer than Pagano's theory. This reduction in the number of d e p e n d e n t variables is significant and should reduce computational costs. Two boundary conditions were examined for the free edge problem. In Case 1 the interlaminar shear stress at the free edge of the 0/90 interface was prescribed to vanish. In Case 2 this free edge shear stress was not prescribed and the alternate condition (61) was prescribed to vanish. Thus in Case 2 the interlaminar shear stress at the free edge is not necessarily zero. Case 2 is an attempt to model the singular behavior present in the linear elasticity solution at the free edge. It was found that the maximum interlaminar shear stresses are significantly greater for Case 2 boundary conditions than for Case 1 (see Table 3). For the type of model developed here, the interlaminar normal stresses at the free edge are finite and do not differ by very much between the two boundary condition cases (see Table 4).
ACKNOWLEDGEMENT Research was sponsored by the Air Force Office of Scientific Research, Air Force Systems Command, USAF, under Contract F49620-82-C-0035, and under Grant AFOSR-83-0191. The US Government is authorized to reproduce and distribute reprints for Governmental purposes notwithstanding any copyright notation thereon.
REFERENCES 1. Pagano, N. J., Stress fields in composite laminates, Int. J. Solids Structures, 14 (1978) 385-400. 2. Lentini, M. and Pereyra, V., An adaptive finite difference solver for nonlinear two point boundary value problems with mild boundary layers, S I A M J. Numerical Analysis, 14 (1977) 91-111. 3. Salamon, Nicholas J., An assessment of the interlaminar stress problem in laminated composites, J. Composite Matl. Supplement, 14 (1980) 177-94. 4. Pagano, N. J. and Soni, S. R., Global-local laminate variational model, Int. J. Solids Structures, 19 (1983) 207-28. 5. Wang, S. S. and Choi, I., Boundary-layer effects in composite laminates: Part I--Free edge stress singularities, J. Appl. Mech., 49 (1982) 541-8. 6. Raju, I. S., Whitcomb, J. D. and Goree, J. G., A new look at numerical analyses of free-edge stresses in composite laminates, NASA Tech. Paper 1751 (1980).