eduzhai > Physical Sciences > Materials Sciences >

Geometrically nonlinear analysis of composite laminates under low velocity impact

  • sky
  • (0) Download
  • 20211030
  • Save International Journal of Composite M aterials 2013, 3(6B): 23-39 DOI: 10.5923/s.cmaterials.201310.04 Geometrically Non-Linear Analysis of Composite Laminated Plates Subjected to Low-Velocity Impact Simon Wang1,*, Xiuqin Zhang2, Yingshun Zhang3, Chenyu Zhang1 1Dep artment of Aeronautical and Automotive Engineer in g, Lou ghborou gh University, Loughborou gh LE11 3TU, UK 2Taiyuan University of Technology, 79 West Yingze Street, Taiyuan 030024, Shanxi, China 3Doosan Babcock Ener gy Limited, Porterfield Road, Renfr ew PA4 8DJ, UK Abstract A B-spline fin ite strip model is developed in the context of a layer-wise plate theory for analysing the geometrically non-linear transient response of laminated composite plates subjected to transverse low-velocity impact. To simp lify the comp licated contact analysis, a Hert z-type contact law has been incorporated into the finite strip (FS) model for accounting for the contact behaviour. The model includes the geometrical non-linearity through use of von Karman's non-linear strain-displacement relationship. The resulting non-linear dynamic problem is solved using the Newmark time-stepping scheme together with Newton-Raphson iteration. Several nu merical applications are described and a close comparison is found between the results calculated through the present model and the existing analytical and experimental res u lts . Keywords Laminated Co mposite, Lo w-velocity Impact, Geo metrical non-linearity, Transient, B-spline Finite St rip 1. Introduction Laminated composite structures have played an important role in aircraft, vehicles and many other demanding applications because of their high strength-to-weight and high stiffness-to-weight ratios. As it is well known, however, these structures are very susceptible to low-velocity impact damage, which can be introduced during manufacture and in service. Impact induced damage can have a significant effect on the strength, stability and reliability of the structures. Therefore, great concern has been received on the low-velocity impact of the structures[1]. For understanding low-velocity impact response of composite structures, several analytical approaches have been used by a number of researchers. These approaches vary from simp le mathematical models, such as spring-mass models and energy-balance models, to more co mp licated dynamic analysis of the impact events. In spring-mass models, the structure is represented by an assemblage of springs and masses. The low-velocity impact event is simu lated by a discrete system with a few degrees of freedom. These kinds of model are generally used to estimate the impact force h istory for the impact events in which the structures behave quasi-statically. Caprino et al[2] used a single degree of freedo m model to p redict elastic impact response of a small g lass cloth-polyester panel due to a drop * Corresponding author: (Simon Wang) Published online at Copyright © 2013 Scientific & Academic Publishing. All Rights Reserved weight impact. In their model, the panel is represented as a linear spring of stiffness k, which corresponds to the static rig idity of the panel at the impact point. Similar models are also used by others[3,4,5]. Sh ivaku mar et al[6] suggested a two degree of freedo m model, in wh ich the impactor and the plate are treated as two masses, and four springs are used to represent the contact stiffness and the bending, shear and memb rane rig idities of the plate. Bucinell et al[7] emp loyed the same model to study the response of composite plates to impacts. In the two-degrees of freedo m model of Sjoblo m et al[8], the plate and contact rigidity is modelled separately using two springs. Toh et al[9,10] also used a two-degrees of freedom model to predict the impact force. Gong and Lam[11] proposed an imp roved two-degrees of freedom spring-mass model by imp lementing the structural damping to determine the contact force between the target and striker during impact. Another kind of frequently used simp le models for the impact events in which the structures behave quasi-statically are energy-balance models[6,12]. In these models, it is assumed that the velocity of the impactor becomes zero when the structure reaches its maximu m deflection. Energy-balance models allo w d irect estimate of the maximu m contact force, considering the conservation of energy without having to compute the entire time history. Whilst they are useful in understanding the main features of low-velocity impact events, the above-mentioned simp le models are inadequate in accounting for the dynamic nature of the composite structures under the low-velocity impact. For understanding the initiation and propagation of the low-velocity impact damage as well as the interaction between the damage and p late dynamic response, dynamic 24 Simon Wang et al.: Geometrically Non-Linear Analysis of Composite Laminated Plates Subjected to Low-Velocity Impact analysis is often needed for accurate predict ion of low-velocity impact behaviour o f the co mposite laminates. In addition, it is the case that for many situations the responses of the structures in low-velocity impact events cannot be viewed as quasi-static ones. Therefore, it is very important to have an insight into dynamic response of the impacted structures. Dynamic analysis of low-velocity impact response of la minated composites generally involves the global dynamic response of the laminates and the local contact between the impactor and the structure. To model the contact, one of the possibilit ies is based on a formu lation of the joint contact problem for the system of impactor-target. This is tightly connected with the part icular nu merical method to be applied, for instance, finite element method (FEM), fin ite difference method, or other method based on some variat ional principles. In FEM analysis of the impact contact problem, contact elements are used to model the laminate and the impactor[13,14,15,16,17,18,19]. The solution of the contact problem may be achieved using several contact algorithms such as the Lag range mult iplier method[13] and the penalty method[19]. The Lagrange multip lier method has the advantage of enforcing the exact constraints, but induces additional para meters which substantially enla rge the overall size of the equations to be solved. The penalty method is relatively simple and does not require additional equations. With contact algorithms, the contact state is detected at each time step and the contact constraints are imposed to the contacting nodes/elements once contact occurs. Since contact is a non-linear problem, fu ll analysis of the impact contact for the system of impactor-target is a co mplex and time-consuming procedure. Hence, simp lificat ion of the contact problem is often made through use of the contact laws in impact analysis of composite laminates. Yang and Sun[20] have proposed a contact law for contact between a sphere and a co mposite laminate based on static indentation tests. This contact law accounts for permanent indentation after unloading cycles and uses different relat ionships between the contact force and the indentation for loading, unloading and reloading processes. Tan and Sun[21] further studied the unloading and reloading process and proposed a modified version of the contact law. Experimental observations[22,23] confirm that loading rate effects during low-velocity impact of co mposites are insignificant. This suggests that statically determined contact laws may be used for impact analysis of composite structures. The above-mentioned contact laws developed by Sun and his co-workers have gained extensive use in analysing the dynamic response of co mposite laminates to impact. It should be noted that these contact laws do not give the contact stress distribution under the indentor. Determination of the contact stress distribution has to appeal to analytical method. In addition, these contact laws have not accounted for damage effect on fo rce-indentation relationship. Wu and Shyu[24] showed that the contact phenomenon is different in small and large indentation stages due to occurrence of laminate damage. In the small indentation stage where the plate is intact, the change of laminate stacking sequence has an insignificant effect on the force-indentation relationship. Beyond the small indentation stage, damage occurs and the indentation spring is stiffened as a result of matrix crack and delamination damage of the laminates. The contact behaviour during low-velocity impact is very similar to that in a static test. Because of the complexity o f the impact prob lem, closed-form exact solutions exist only for simple cases. In most situations, approximate analytical and/or numerical methods have to be adopted. In conjunction with contact laws, dynamic analyses of low-velocity impact of co mposite laminates have been carried out analyt ically by a nu mber of researchers. Sun and Chattopadhyay[25] studied a simp ly-supported orthotropic plate subjected to central impact using the first order shear deformable p late theory (FSDPT) developed by Whitney and Pagano[26]. Dobyns[27] used the same method but replaced the concentrated impact load by a uniform patch pressure to avoid transverse shear force singularity at the contact point. Qian and Swanson[28] examined two solution techniques which were based on series expansion. One of them was based on a Rayleigh-Ritz approach with nu merical integration in time, and the other was an analytical approach using Laplace transformation of the governing differential equations, but requiring a linearisation of the contact law. Carvalho and Soares[29] studied a simp ly-supported composite plate subjected to an impact load utilising the technique of Fourier series expansion for the solution of the dynamic p late equations. Co mparison was made between the numerical solution technique based on Newmark integration method and the analytical formu lation using the Lap lace transform technique. Pierson and Vaziri[30] presented a Fourier series solution that retains the frequencies associated with rotary inert ia effects. A double Fourier series expansion and the Timoshenko small-increment method were used by Ambur et al[31,32] for determining the transient response of simp ly supported, rectangular laminated plates subjected to impact loads. Heit zer[33] studied the interaction of an impactor and a clamped, an isotropic plate at lo w-velocities by assuming a series expansion for the plate deflection. Large deflection was taken into account and the contact law was used. It is well known that FEM is one of the most powerful tools of solution in structural analysis. Besides 2D or 3D fin ite element analyses of the joint impactor-target system[13-19,34-36], many researchers have employed the FEM in conjunction with contact laws for the analyses of composite laminates due to low-velocity impact. Sun and his collaborators[21,37,38] used static indentation laws and FEM based on FSDPT to analyse the impact responses of composite laminates. The init ially stressed composite laminates were studied in[37] and dynamic large deflection was considered in[38]. Wu and Chang[39] developed a 3D fin ite element code for t ransient dynamic analysis of laminated composite plates due to transverse foreign object impact in wh ich a contact law is incorporated. Co mbined with fa ilure criteria this 3D finite ele ment code has been used International Journal of Composite M aterials 2013, 3(6B): 23-39 25 by Wu and Springer[40], Choi and Chang[41], Finn and Springer[42,43], and Scarponi et al[44] for failure analyses of impacted laminates. Choi and Chang[45-47] developed a 2D FEM for studying the impact damage of laminated composite beam resulting fro m the line-loading impact. Davies and Zhang et al[48-51] investigated impact-induced damage using FEM based on Mindlin’s plate theory combined with a contact law[22]. The FEM based on FSDPT and a Hert zian-type indentation law was emp loyed by Hu et al[36] for analysing the transient response of composite laminates with mu ltiple delaminations subjected to low-velocity impact. It appears that the majo rity of p revious analytical and numerical analyses are limited to s mall-deflection behaviour. Although such linear analyses are practical and useful for many impact problems of composite structures, it is often the case that the geometric non-linearity (GNL) has very significant effects on the impact response. In this paper, a layer-wise B-spline finite strip model developed in[52] is extended through introduction of time dimension for the GNL transient analysis of laminated composite plates subjected to transverse low-velocity impact of spherical object within the context of a layer-wise plate theory[53]. It is noted that in[54,55] a B-spline FSM has been used for transient analysis of laminated composite plates subjected to dynamic loads based on the first order shear deformable plate theory (FSDPT). In the present study, the geometrical non-linearity is taken into consideration by use of von Karman 's non-linear strain-displacement relationship. The dynamic problem is solved using the Newmark t ime-stepping scheme in the time do main and the solution of the resulting non-linear equations is sought with Newton-Raphson iteration. The impact problem in the small-deflection regime is briefly discussed as well. Several numerical applicat ions are presented. It is noted that material damage and delamination are not dealt with in the current work. 2. Problem Description and Fundamentals m, v0 A x1 h B x2 x3 Figure 1. T ransverse impact on a rectangular plate Consider a laminated rectangular plate with arbitrary lay-ups subjected to transverse low-velocity impact of a spherical object as shown in Figure 1. The impactor is of radius rs, mass m and initial velocity v0. It is assumed that the vibrations of the elastic impactor can be neglected. An orthogonal Cartesian co-ordinate system xi (i=1,2,3) is used in this paper. 2.1. Basic Plate Equati ons In the present study, the layer-wise plate theory proposed by Reddy[53] is adopted for representation of displacement behaviour in the plate. Through the thickness direction, the laminated plate is divided into a desired number, N, of numerical layers, which can be less than, equal to, or greater than the number of physical layers. The assumed displacements ( µα , u3 ) at a general point of the laminate, with t as the time d imension, take the form as what follows. uα ( xβ , x3 , t) = uα ( xβ , t) + Ψ i (x3 )u i α ( xβ , t ) (1) u3 (xβ , x3 ,t) = u3 (xβ ,t) where the usual Cartesian indicial notation is adopted. The Greek subscripts, which take values 1 and 2, and subscript 3 refer to x1, x2, x3 directions, respectively. Repeated indices imply the summation convention. Superscript i, ranged from 1 to N, is related to the nodes through the thickness. uα and u3 denote the three displacement co mponents of the reference plane (x3=0). Ψ i are piecewise continuous functions and defined as Ψ i = Rik~ Ψ~ k~ (2) where i is ranged from 1 to N and k~ fro m 1 to N+1. Ψ~ k~ is defined in terms of linear Lagrange interpolation functions   x3 − x3k −1 Ψk (x3) =  xxx33k3kk++−11x−−3kxx−33k1   0 x3 ∈[x3k−1, x3k ]∈[x13, x3N +1] (3) x3 ∈[x3k , x3k+1]∈[x13, x3N +1] x3 ∉[x13, x3N +1] and Rik~ is defined as 1 Rik=   x3m / x3m+1   0 i = k < m or i = k −1 ≥ m i= k= m (4) others Here m is the number of the layer at wh ich the reference plane is. The resulting displacement configuration is shown in Figure.2. 26 Simon Wang et al.: Geometrically Non-Linear Analysis of Composite Laminated Plates Subjected to Low-Velocity Impact uα = uα + Ψiuαi u3 = u3 x3 u uαN uαN −1 uαm x3N +1 N-th layer x3N h/2 x m+1 3 x3m m-th layer xα h/2 x32 1st layer x31 uα2 uα1 ( x3m x m+1 3 uαm ) Figure 2. In-plane displacement configuration of layerwise theory The expressions of the strain-displacement relat ionship can be given by substitution of the displacement field Eq.1 in the Green's expressions for in-plane non-linear strains and neglecting higher order terms in a manner consistent with von Karman's assumptions. They are: ε= αβ εαβ + Ψiεαi β (5) γ= α 3 γα 3 + Ψi,3γ i α 3 where εαβ = 1 2 (uα ,β + uβ ,α ) + 1 2 u3,α u3,β = εαi β 1 2 (uαi ,β + uβi ,α ) (6) = γα 3 u= 3,α γαi 3 uαi The strain energy of the laminate considered can be obtained as an integral over the reference plane area S. The result is: ∫∫ = U 1 2 S (Nαβ εαβ + Nαi β εαi β + Qα 3γα 3 + Qαi 3γαi 3)dx1dx2 (7) where  Nαβ   Nαi β    =   Aαβγω Bαi βγω Bαkβγω Dαikβγω    ε ε γω k γω    Qα Qαi 3 3    =   Aα 3β Bαi 3β 3 3 Bαk 3β 3 Dαik3β 3    γ  γ β k β 3 3    (8) Nαβ , Nαi β , Qα 3 , and Qαi 3 are generalised stresses per unit length. Aαβγω , Bαi βγω , etc, are the stiffness coefficients of the laminate wh ich are defined as ∫ ( Aαβγω , Bαi = βγω , Dαikβγω ) h2 −h 2 Qαβγω (1, Ψi , Ψ iΨ k )dx3 ∫ ( Aα 3β 3, Bαi 3= β 3, Dαik3β 3) h2 −h 2 Qα 3β 3 (1, Ψ ,i3 , Ψ,i3Ψ,k3 )dx3 (9) where Qαβγω and Qα 3β 3 denote the transformed reduced stiffness coefficients. Qαβγω possess symmetry in the indices α and β, γ and ω, and the pairs of αβ and γω. The similar sy mmetrical property also applies to Qα 3β 3 . 2.2. Fi nite Stri p Approxi mation The whole plate is modelled with a number o f finite strips along the crosswise x2-direction. Each strip is further partitioned longitudinally into q sections. A typical individual finite strip element (quadratic strip) is shown in Figure 3. Over the strip each of the fundamental International Journal of Composite M aterials 2013, 3(6B): 23-39 27 displacement quantities uα , u3 , and uαi can be approximated as a function of mult iplicative type, in which q+k B-spline functions of degree k are used in the longitudinal x1-direction and simp le polynomia ls of degree n in the crosswise x2-direction. Mathematically, the fundamental displace ment quantities are given as uα  Φ    u3 u1i u i 2     = N I      Φ Φ −1     Φ d  α dd13i d i 2 I I I I       (10) Here dα , d 3 and d i α are colu mn matrices of generalised displacement parameters associated with uα , u3 and uαi , respectively. The row matrices Φ and Φ −1 are modified B-spline function bases of order k and k-1 in the longitudinal x1 -direction. NI are polyno mial functions of degree n in the crosswise x2 -direction and here Lagrangian shape functions are used. The superscript i, ranged fro m 1 to N, is related to the nodes through the thickness. The capital superscript I, ranging fro m 0 to n, denotes the number of a reference line in the crosswise x2 -direction of the strip. b A h x1 uαN uα ,u3 x2 u 1 α uα = uα + Ψiuαi x3 u3 = u3 Figure 3. A typical layerwise finite strip element f are colu mn matrices of generalised displacement parameters and generalised force due to applied loading, F is impact contact force, c is a column matrix, resulting fro m the product of the polynomials in the x2 -direction and the B-spline functions in the x1-direction at the impacting point of the plate, and cTd defines the deflection of the reference plane of the plate at impact point. The relationship between the contact force F and the indentation depth α is assumed as[20] F = kcα3/ 2 when loading (12a) F = Fm [(α − α0 ) /(α m − α0 )]5 / 2 when unloading(12b) F = Fm [(α − α 0 ) /(α m − α 0 )]3/ 2 when reloading (12c) where kc is a constant. Fm denotes the maximu m contact force just before unloading, αm is the indentation corresponding to Fm, and α0 is the permanent indentation during the loading/unloading cycle. kc can be determined by experiment or simply calculated fro m the modified Hertz contact coefficient proposed by Sun[20] as kc = 4 3 rs [(1 − ν 2 s ) / Es + 1/ E2 ]−1 (13) here rs, νs and Es are the radius, the Poisson's ratio and the Young's modulus of the impactor, respectively, and E2 is the transverse modulus of elasticity of the plate. Eq.12 can be expressed in a general form as F = ks (α − α0 )q (14) where ks and q are the contact coefficient and a constant, respectively. It is noticed that for the loading, un loading and reloading process α0, ks and q may be different. The indentation depth α is the difference of the displacement of the impactor and the deflection o f the reference plane of the plate at impact point, so that (F / ks )1/ q + α0 = v0t − 1 m s − cTd (= α) (15) where tt s = ∫0 ∫0 Fdtdt (16) Eqs.11 and 15 together with appropriate initial conditions define the present impact dynamic problem. It is here assumed that d , d and d are known at t ime t = 0 and F , s , s are equal to zero at time t = 0 . 3. Governing Equations of Motion In the absence of damping, the governing equations of the plate motion can be given through use of Hamilton's principle in terms of strain energy U, kinetic energy T and potential energy W as Md + Kd = Fc + f (11) where M and K denote the mass matrix and the effective stiffness matrix, respectively, and K are functions of d, d and 4. Problem Solution In the time do main the impact problem defined in Eqs.11 and 15 can be solved with a few of different approaches. In the present investigation the solution of the problem is sought using the popular Newmark time integral algorith m. Suppose that the time period concerned is divided into a number of equal time intervals Δt and consider the satisfaction of Eqs.11 and 15 at time (n+1)∆t and write Md n+1 + Kd n+1 = Fn+1c + f n+1 (17a) 28 Simon Wang et al.: Geometrically Non-Linear Analysis of Composite Laminated Plates Subjected to Low-Velocity Impact (Fn+1 / ks )1/ q + α0 = v0 (n + 1)Δt − 1 m sn+1 − cT d n+1 (17b) with appro ximat ions for d n+1 d n+1 , sn+1 and sn+1 as d n+1 = dn + Δtd n + Δt 2 2 (1 − β2 )dn + Δt 2 2 β2d n+1 ( 18a ) d n+1 = d n + Δt(1 − β1 )d n + Δtβ1d n+1 s n +1 = sn + Δtsn + Δt 2 2 (1 − β2 )Fn + Δt 2 2 β2 Fn+1 (18b) sn+1 = sn + Δt(1 − β1 )Fn + Δtβ1Fn+1 where the parameters β1 and β2 can be chosen as such that good approximation properties to the algorithm are yielded. In this study the constant-average-acceleration version of Newmark method is used, and consequently, β1=β2=1/2. The algorithms y ield unconditional stability for linear problems at least. Belytschko and Schoeberle[56] use the Newmark-β method for a nonlinear structural dynamics problem. Their proof of unconditional stability with β = 0.25 using an energy method applies only to nonlinear material properties but they state that numerical results show unconditional stability with geometric nonlinearity also. Substitution of Eq.18 into Eq.17 leads to a set of non-linear equations in wh ich the unknowns are dn+1 and Fn+1. The non-linear equations can be written in terms of residual forces as R = Mdn+1 + K(dn + Δtd n + Δt 2 2 (1 − β 2 )dn + Δt 2 2 β2dn+1) (19a) − Fn+1c − fn+1 = 0 r = cT d n+1+ (Fn+1 / k s )1/ q + 1 m (sn + Δtsn + Δt 2 2 (1 − β 2 )Fn + Δt 2 2 β 2 Fn+1 ) (19b) + α 0 − v0 (n + 1)Δt = 0 To seek the solution of the non-linear problem defined in Eq.19, the Ne wton-Raphson iteration scheme is here adopted. It is noted that, to the first order, Eq .19 can be appro ximated as R i+1 ≈ Ri + (M + Δt 2 2 β 2 K i T )δd i n+1 − δFni+1c = 0 ( 20a ) r i+1 ≈ r i + Δt 2 2 β 2c T δd i n+1 + k t δFni+1 = 0 (20b ) where K i T = ∂(Kd n+1 ) ∂d n+1 d n+1 = d i n +1 (21) k i t = ( Fi n +1 )1 / q −1 qk 1 s / q + 1 m Δt 2 2 β2 δd i n+1 = d i +1 n+1 − d i n+1 δFni+1 = F i+1 n+1 − Fi n+1 (2 2) Here s u p ers crip t i denotes the iteration counter. K i T is the symmetric tangent stiffness matrix evaluated at iterat ion i. k i t is defined with the condition of Fni+1 > 0 . δd i n+1 and δFni+1 are increments, or iterative corrections, of the unknowns d n+1 and Fn+1 at iteration i, respectively, and they can be given from Eq.20 as δ din+1 = −(M + Δt 2 2 β2K i T + 1 kti Δt 2 2 β2ccT )−1(Ri + ri kti c) δ Fni+1 = − k1ti (ri + Δt 2 2 β2cTδ din+1) (23) Through the use of Eq.22 with the corrections given in Eq.23, the iteration procedure results in improved solution for the unknowns d n+1 and Fn+1. The procedure is repeated until the following convergence criteria a re satisfied δd iT n+1 δd i n+1 <ε d d i T i n+1 n+1 and δFni+1 Fni+1 <ε (24) where ε = 0.0005 is used in this work. Up to now the impact dynamic transients are concerned with geomet ric non-linearity and the solution procedure is described in details as above. For comparison, the impact problem is also investigated within the context of small deflection in this paper. In fact there exists no difference between the geometric non-linear and linear analyses of the impact problem in the temporal appro ximation. W ithin a time step, however, the solution procedure used for the linear transient response is different fro m that described above for the non-linear one, and is here discussed below. In the case of geometric linear analyses of the impact problem, the non-linear strain terms are ignored in the strain-displacement relationship. As a result, the stiffness matrix K beco mes independent of the plate displacement. Fro m Eq.19a the following equations can be given (M + Δt 2 2 β 2 K )(d c n + 1 , d f n+1 ) = (25) (c , f n+1 − K(d n + Δtd n + Δt 2 2 (1 − β 2 )dn )) and d n+1 = Fn +1d c n +1 + d f n +1 (26) Eqs.25 and 26 together with Eq.19b define the present impact problem in terms of small deflection. It is obvious that Eq.25 defines two sets of linear algebraic equations with respect to unknowns d  c n +1 and d f n +1 . These two s ets of linear equations can be solved simultaneously since they share the same coefficient matrix. The solution procedures involve triangular decomposition of the coefficient matrix and then forward elimination and back substitution. Most of the computing efforts are expended in the procedure of triangular deco mposition of the coefficient matrix during the solution of the equations. Noting the coefficient matrix is International Journal of Composite M aterials 2013, 3(6B): 23-39 29 independent of time increments, the triangular decomposition of the coefficient matrix is performed ahead of the time iteration in the present methodology and in each time increment only the forward elimination and back substitution are involved. With this strategy, the computer efforts can be greatly reduced. Substitution of Eq.26 into Eq.19b leads to a non-linear equation with respect to Fn+1, and again, the Newton-Raphson iteration technique is adopted for the solution of the non-linear equation. The iterative corrections for Fn+1 can be written as δ Fni+1 =− k i t ri + c T d c n +1 (27) where k i t and r i are the same as those given before and d c n +1 are evaluated fro m Eq.19. 5. Numerical Applications The above analysis procedure has been imp lemented in an in-house computer software. It can be used for impact problems and, neglecting impact-related terms, for plate transient response under dynamic loading. Nu merical tests have been conducted to verify the valid ity of the present computational model, which include a series of examp les of composite laminate plates under low-velocity impact. The selected examples of the impact problem refer to two typical kinds of low-velocity impact events, i.e. s mall impact mass with relatively high speed and large mass with low speed. In all the applications, the quadratic Lagrangian shape 350 functions are applied in the crosswise x2-direction. Unless otherwise specified, the simply supported boundary conditions are defined such that the plate is simply supported for out-of-plane behaviour whilst the in-plane movement of the sides is constrained in the tangent direction but allowed in the norma l direction. For the impact proble ms considered, the contact law defined in Eq.12a is used for both the loading and unloading processes. The details are presented in the following subsections. 5.1. Impact of a[0/90/0/90/0]s Plate Consider the impact proble m involving central impact of a simp ly supported square T300/934 graphite/epoxy [0/90/0/90/0]s plate by a steel object with a spherical impact tip of diameter 12.7mm. The plate has dimensions of 200 × 200 × 2.69 mm. The material properties of a T300/934 graphite/epoxy lamina are given as E1 = 120 GPa, E2 = 7.9 GPa, ν12 = 0.3 G12 =G13 = G23 = 5.5GPa, ρ = 1580 kg/m3 Two cases of impact events are considered here, i.e., impact by a 7.5g mass with a velocity of 3m/s and by a 19.6g mass with a velocity o f 10m/s, respectively. In the present study, the contact law defined in Eq.12a is simp ly emp loyed for the whole contact processing. Since the laminate is thin one numerical layer is used with shear correction factor k 2 55 = 0.88075 and k 2 44 = 0.73173 , evaluated by the method in Reference[26] Case 1 Impact by a 7.5g mass of 3m/s ∆t=20×10-6 Sec 300 ∆t=10×10-6 Sec ∆t=2×10-6 Sec 250 ∆t=1×10-6 Sec Half plate, 200 5 quadratic strips, 10 sections ks=25644 N/mm3/2 150 Contact force (N) 100 50 0 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 Time (Sec) x 10 -4 (a) Contact force 30 Simon Wang et al.: Geometrically Non-Linear Analysis of Composite Laminated Plates Subjected to Low-Velocity Impact 0.5 0.4 0.3 Central displacement (mm) 0.2 0.1 ∆t=20×10-6 Sec 0 ∆t=10×10-6 Sec ∆t=2×10-6 Sec -0.1 ∆t=1×10-6 Sec Half plate, -0.2 5 quadratic strips, 10 sections ks=25644 N/mm3/2 -0.3 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 Time (Sec) x 10 -3 (b) Central deflection Figure 4. Convergence with respect to time step for the problem of a[0/90/0/90/0]s plate impacted by a 7.5g mass of 3m/s 400 ne=2, q=4 350 ne=3, q=6 ne=3, q=10 300 ne=5, q=10 ne=10, q=10 ne=10, q=20 250 ∆t=2×10-6 Sec 200 Contact force (N) 150 100 50 0 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 Time (Sec) x 10 -4 (a) Contact force International Journal of Composite M aterials 2013, 3(6B): 23-39 31 0.5 0.4 0.3 Central displacement (mm) 0.2 0.1 0 -0.1 -0.2 -0.3 0 ne=2, q=4 ne=3, q=6 ne=3, q=10 ne=5, q=10 ne=10, q=10 ne=10, q=20 ∆t=2×10-6 Sec 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 Time (Sec) x 10 -3 (b) Central deflection Figure 5. Convergence with respect to plate mesh for the problem of a[0/90/0/90/0]s plate impacted by a 7.5g mass of 3m/s Nu merical tests are carried out on the impact problem of the[0/90/0/90/0]s laminated plate by a 7.5g mass with a velocity of 3m/s. The examp le chosen here is the one studied by Sun and Chen[37], and Qian and Swanson[28]. The convergence of time integration is first studied. Taking advantage of symmetry in this thin p late geometry, half of the plate (crosswise direction of the strip) is modelled using five quadratic strips in the crosswise direction. Each strip has ten sections in the longitudinal direction and one numerical layer through the thickness. The contact coefficient kc is evaluated using Eq.13, wh ich is equal to 25644 N/mm1.5. The time increments are taken as 20 µs, 10 µs, 2 µs and 1 µs. Figures 4a and 4b show the contact force history and the central deflection of the plate, respectively. It can be seen that all the chosen time increments give stably converged solutions. As the time increment decreases, convergence behaviour can be observed. It is noted that the results for contact force and deflection histories obtained using a time increment of 20 µs are fairly reasonable. The results obtained by using a time increment of 2 µs are almost identical to those using a time increment of 1 µs, as well as 0.1 µs (not presented in the figure), within p lotting accuracy. The predicted maximu m contact forces are 303.2 N, 291.0 N, 290.2 N, and 290.4 N, respectively, when using time steps of 20µs, 10µs 2µs and 1µs, respectively. Next, the convergence with respect to the plate mesh is investigated. Half of the plate is modelled using a nu mber, ne, of quadratic strips with d ifferent section number q as (1) ne=2, q=4, (2) ne=3, q=6, (3) ne=3, q=10, (4) ne=5, q =10, (5) ne=10, q=10, and (6) ne=10, q=20. The time increment is taken as 2 µs. The contact coefficient kc is equal to 25644 N/mm1.5. Figures 5a and 5b show the convergence of contact force and the central deflection histories with respect to the plate mesh. A converged trend can be observed in the figures as the number of the strip and sections increases. It is found that use of coarse mesh yields reasonable deflection prediction but unacceptable estimat ion of the contact force. Both the contact force and deflect ion histories calculated by using 5 strips with 10 sections are nearly as good as those using 10 strips with 20 sections. The follo wing set of numerical tests is carried out to study the effect of the contact coefficients. For the same impact problem, Sun and Chen[37] meshed one quarter of the plate with 8×8 n ine-node quadrilateral fin ite elements and used a contact law as Eq.12 with k c = 44683 N/ mm1.5, whereas Qian[28] and Swanson applied the Rayleigh-Rit z technique (100×100 modes) and incorporated the contact law as Eq.12a with the contact area changing at each time step and kc = 26540 N/mm1.5. In the present study, the contact law defined in Eq.12a is used with a series of contact coefficients, kc = 44683 N/mm1.5 (the same as Sun and Chen for the loading processing), kc = 36266 N/ mm1.5 and kc = 25644 N/mm1.5 (evaluated fro m Eq.13). Figure 6 shows the comparison of the contact force and central deflection histories. In the present FSM analysis, 20 quadratic strips with 40 sections are used to mesh the whole plate and the time increment is taken as 1 µs. The pred icted maximu m contact forces are 303.8 N, 297.0 N, and 287.0 N, respectively, when using contact coefficients of 44638 N/mm1.5, 36266 N/ mm1.5 and 25644 N/ mm1.5, respectively. It is noted that the contact force and deflection histories predicted by the present FSM using contact coefficient of 25644 N/mm1.5 (close to that of Qian 32 Simon Wang et al.: Geometrically Non-Linear Analysis of Composite Laminated Plates Subjected to Low-Velocity Impact and Swanson) are basically the same as those of Qian and Swanson using the Rayleigh-Rit z method. Use of contact coefficients 44638 N/mm1.5, the same value as used by Sun and Chen, gives an agreed prediction for the deflection history to that of Sun and Chen, but a slightly h igher maximu m contact force than that of Sun and Chen. Since spline functions are the smoothest piecewise polynomials, they possess the blended advantages of smooth analytical functions and versatile polynomials. As fewer degrees of freedom are required for the same accuracy, the present FSM is much more co mputationally efficient in comparison with the conventional FEM in which conventional polynomials are used. More details can be found in authors’ earlier work[52, 54, 55]. 350 FEM, Sun and Chen Rayleigh-Ritz, Qian and Swanson 300 Present FSM, ne=20, q=40, ∆t=1µ, whole plate ks=44683 N/mm3/2 250 ks=36266 N/mm3/2 ks=25644 N/mm3/2 200 Contact force (N) 150 100 50 0 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 Time (Sec) x 10 -4 (a) Contact force 0.45 0.4 0.35 Central displacement (mm) 0.3 0.25 0.2 0.15 0.1 0.05 0 0 FEM, Sun and Chen Rayleigh-Ritz, Qian and Swanson Present FSM, ne=20, q=40, ∆t=1µ, whole plate ks=44683 N/mm3/2 ks=25644 N/mm3/2 0.2 0.4 0.6 0.8 1 Time (Sec) 1.2 1.4 1.6 x 10 -3 (b) Central deflection Figure 6. Comparison of predicted contact force and central deflection histories of a[0/90/0/90/0]s plate impacted by a 7.5g mass of velocity 3m/s

... pages left unread,continue reading

Document pages: 17 pages

Please select stars to rate!


0 comments Sign in to leave a comment.

    Data loading, please wait...