Transient dynamic response of delaminated composite rotating shallow shells subjected to impact

In this paper a transient dynamic finite element analysis is presented to study the response of delaminated composite pretwisted rotating shallow shells subjected to low velocity normal impact. Lagrange’s equation of motion is used to derive the dynamic equilibrium equation and moderate rotational speeds are considered wherein the Coriolis effect is negligible. An eight noded isoparametric plate bending element is employed in the finite element formulation incorporating rotary inertia and effects of transverse shear deformation based on Mindlin’s theory. To satisfy the compatibility of deformation and equilibrium of resultant forces and moments at the delamination crack front a multipoint constraint algorithm is incorporated which leads to unsymmetric stiffness matrices. The modified Hertzian contact law which accounts for permanent indentation is utilized to compute the contact force, and the time dependent equations are solved by Newmark’s time integration algorithm. Parametric studies are performed in respect of location of delamination, angle of twist and rotational speed for centrally impacted graphite-epoxy composite cylindrical shells.


Introduction
Delamination in fibre-reinforced composites resulting from interlaminar debonding of constituting laminae causes strength degradation and can promote instability.Due to high specific strength composite materials are advantageous in a weight sensitive application such as turbomachinery blades which with low aspect ratio could be idealized as twisted rotating cylindrical shells or plates (Fig. 1).But the laminated curved structures namely, shallow shells have relatively low through-the-thickness strength and are susceptible to high strain rate loading caused by impact of foreign objects.On the other hand, the delaminated structures exhibit new deformation characteristics depending on size and location of the delamination, and the presence of invisible delamination can be detected with the help of prior knowledge of these characteristics.Moreover, the initial stress system in a rotating shell due to centrifugal body forces may aggravate the delamination damage due to impact.Hence, insight into the transient dynamic response of delaminated composite pretwisted rotating shells subjected to localized contact loading is of significant concern as a precursor to the application of twisted composite shells in the critical parts of aero-engines in order to ensure operational safety.The works on pretwisted composite plates were first carried out by Qatu and Leissa [13].They determined the non-dimensional natural frequencies of stationary plates using the laminated shallow shell theory in conjunction with the Ritz method.The investigations by Wang et al. [10] and Shaw et al. [7] provided the first known results on natural frequencies of rotating composite plates.Pagano and Soni [15], Seshu et al. [17] and Bhumbla et al. [18] analyzed rotating composite blades using plate model.The analyses focussed the derivation of two approximate analytical models for determining the stress-strain field [15], the details of fabrication and experimentation of glass-epoxy composite blades for determining natural frequencies and steady state strain along with the effect of rotational speed on it [17] and the vibration characteristics including geometric nonlinearity [18].Later on, McGee and Chu [16] used the Ritz method to carry out the three-dimensional continuum vibration analysis including full geometric nonlinearities and centrifugal accelerations in the blade kinematics.
The combined effect of rotation and pretwist on composite plates was demonstrated by Bhumbla and Kosmatka [19].A nonlinear finite element technique was developed to study the nonlinear static deflection and vibration behaviour of spinning pretwisted composite plates.Of late, Kee and Kim [24] used twisted cylindrical shell model to carry out the free vibration analyses of rotating composite blades.The impact response of initially stressed composite plates received attention in the investigations by Sun and Chattopadhyay [5], and Sun and Chen [6].In both the works simply supported boundary condition was considered, but the former work involved the solution of nonlinear integral equation and the later one dealt with numerical method employing the shear flexible nine noded isoparametric plate finite element.Although delamination is one of the most feared damage modes in laminated composites the impact behaviour of delaminated structures has not been addressed until the investigations by Sekine et al. [9] and Hu et al. [14] wherein simply supported plates with single and multiple delamination were considered for the analyses.
So it is clear that no attention has been paid so far to the initially stressed delaminated composite plates or shells especially under impact load.To the best of the authors' knowledge there is no literature available, which deals with delaminated composite pretwisted rotating shells subjected to low velocity impact and considers its transient dynamic analysis by finite element method.The present work is aimed at investigating the effects of delamination on low velocity normal impact response of composite pretwisted rotating shallow shells.The study herein concentrates upon long cylindrical shells as defined by Aas-Jakobsen's parameters.The finite element model is based on Lagrange's equation of motion and the investigation is carried out for moderate rotational speeds for which the Coriolis effect is negligible.An eight noded isoparametric plate bending element is employed in the finite element formulation.Effects of transverse shear deformation based on Mindlin's theory and rotary inertia are included.To satisfy the compatibility of deformation and equilibrium of resultant forces and moments at the delamination crack front a multipoint constraint algorithm is incorporated which leads to unsymmetric stiffness matrices.The modified Hertzian contact law [6,21,23] which accounts for permanent indentation is utilized to compute the contact force, and the time dependent equations are solved by Newmark's time integration algorithm (constant-average-acceleration method).Parametric studies are performed in respect of location of delamination, angle of twist and rotational speed for graphite-epoxy composite cylindrical shells subjected to impact at the centre.Numerical solutions obtained for symmetric laminate illustrate the information necessary for failure analysis and predicting impact damage.

Governing equations
A shallow shell is characterized by its middle surface which is defined by the equation [3] where R x , R y and R xy denote the radii of curvature in the x and y directions and the radius of twist, respectively.The dynamic equilibrium equation of the target shell for moderate rotational speeds is derived employing Lagrange's equation of motion and neglecting Coriolis effect the equation in global form is expressed as [1] [ where [M ], [K] and [K σ ] are global mass, elastic stiffness and geometric stiffness matrices, respectively.{F (Ω 2 )} is the nodal equivalent centrifugal forces and {δ} is the global displacement vector.[K σ ] depends on the initial stress distribution and is obtained by the iterative procedure [11,20,22] upon solving For the impact problem, {F } is given as where F C is the contact force and the equation of motion of the rigid impactor is obtained as [6] Note that m i and ẅi are mass and acceleration of impactor, respectively.

Contact force
The evaluation of the contact force depends on a contact law which relates the contact force with the indentation.The present analyses concentrate upon rotating cylindrical shallow shells having a large ratio of the radius of curvature to its thickness as well as a high value of width to thickness ratio.Yang and Sun [21] proposed a power law based on static indentation tests using steel balls as indentors.This contact law accounted for permanent indentation after unloading cycles i.e. collisions upon the rebound of the target structure after the first period of contact were considered.The modified version of the above contact law obtained by Tan and Sun [23] was used in the analyses of Sun and Chen [6].Centrifugal forces arising out of rotation generate initial stresses and therefore the contact force model following Sun and Chen [6] has been incorporated into the present finite element formulation as the study focuses numerical simulation of delaminated graphite-epoxy cylindrical shells.If k is the contact stiffness and α m is the maximum local indentation, the contact force F C during loading is evaluated as [6] The contact stiffness for cylindrical shell as proposed by Yang and Sun [21] is determined from the following relation where R i , E i and υ i are the radius, modulus of elasticity and Poisson's ratio of the impactor, and R s and E 2 are the radius of curvature and transverse modulus of elasticity of the composite cylindrical shell target.
The indentation parameter α depends on the difference of the displacements of the impactor and the target structure at any instant of time, and so also the contact force.The values of α are changing with time because of time varying displacements of both the rigid impactor and the target structure.So at an instant the maximum indentation takes place and as a result the maximum contact force is also obtained.At this instant the displacement of the impactor also attains the maximum value [25].Thereafter, the displacement of the impactor gradually decreases, but the target point displacement keeps on changing and finally increases to a maximum value and at some point of time these two displacements become equal [25].This leads to zero value of indentation and eventually the contact force becomes zero.At this instant the impactor looses contact with the target.The process after attaining the maximum contact force till the reduction of contact force to zero value is essentially referred to as unloading [6].If the mass of the impactor is not very small, a second impact may occur upon the rebound of the target structure leading to a same phenomenon of contact deformation and attainment of maximum contact force.There is always a time gap between the end of the first contact period and the beginning of the second impact.This process is known as reloading.If F m is the maximum contact force at the beginning of unloading and α m is the maximum indentation during loading, the contact force F C for unloading and reloading are simulated as [6] and where α 0 denotes the permanent indentation in a loading-unloading cycle and is determined as [6] α 0 = 0 when α m < α cr , and where β c is a constant and α cr is the critical indentation beyond which permanent indentation occurs, and the respective values are 0.094 and 1.667 × 10 −2 cm for graphite-epoxy composite [6].
Neglecting the contribution of plate displacements along global x and y directions, the indentation α is given as [2,6,25] where w i and w p are displacement of impactor and target plate displacement along global z direction at the impact point (x c , y c ), respectively and φ is the angle of twist.The components of force at the impact point in global directions are given by [2] F ix = 0, F iy = F C sin φ and

Multipoint constraints
Figure 2 shows the plate elements at a delamination crack front.The nodal displacements of the delaminated elements 2 and 3 at the crack tip are expressed as (Gim [4]) ) where ūj , vj and wj are the mid-plane displacements, zj is the z coordinate of the mid-plane of element j and θ x , θ y are the rotations about x and y axes, respectively.The transverse displacements and rotations of elements 1, 2 and 3 at a common node have the same values expressed as [4]  The in-plane displacements of all the three elements at the crack tip are equal and consequently, the mid-plane displacements of elements 1, 2 and 3 are related as [4] ū2 = ū1 − z2 θ x and v2 = v1 − z2 θ y (14) ū3 = ū1 − z3 θ x and v3 = v1 − z3 θ y (15) Equations ( 13), ( 14) and ( 15) relating the nodal displacements and rotations of elements 1, 2 and 3 at the delamination crack tip, are the multipoint constraint equations used in the finite element formulation to satisfy the compatibility of displacements and rotations.
Mid-plane strains between elements 2, 3 and 1 are related as [4] {ε} j = {ε} 1 + zj {κ} (j = 2, 3) where {ε} represents the strain vector and {κ} is the curvature vector being identical at the crack tip for elements 1, 2 and 3.The in-plane stress resultants, {N } and the moment resultants, {M } of elements 2 and 3 can be expressed as [4] {N where [A], [B] and [D] are the extension, bending-extension coupling and bending stiffness coefficients of the composite laminate, respectively.It is quite obvious from Eqs ( 17) and ( 18) that the elastic stiffness matrices are unsymmetric in nature.

Results and discussion
The results obtained from the computer code developed on the basis of present finite element modelling are validated with those in the literature both in respect of impact and delamination model.Figure 3(a) shows time histories of contact force for symmetrically laminated cross ply ([0/90/0/90/0] S ) composite plate under simply supported boundary condition as analysed by Sun and Chen [6] using finite element technique, while Fig. 3(b) presents variation of natural frequencies of composite cantilever beam with relative position of the delamination as computed by Krawczuk et al. [12].The discrepancies in the results with those of Sun and Chen could be attributed to the fact that the present analyses consider the full laminated plate whereas Sun and Chen used quarter plate  along with symmetry boundary conditions.Moreover, the number of nodes used in the finite elements, the order of numerical integration and the converged values of time step are different in two cases.
The study is carried out to investigate the effects of relative position of delamination, angle of twist and rotational speed on the normal impact response of graphite-epoxy composite rotating shallow shells.Long cylindrical shells as defined by Aas-Jackobsen's parameters (Γ < 7, Λ < 0.12) with bending stiff ([0 2 /± 30] S ) laminate configuration [8] are considered for the analyses and Aas-Jackobsen's parameters are given as Accordingly the dimensions of length (L), width (b), radius of curvature (R) and thickness (h) of the shells are adopted as 0.34 m, 0.2 m, 0.2 m and 0.002 m, respectively.The above dimensions satisfy the limit of the values of Γ and Λ for long shell which are computed as 5.86 and 0.099, respectively.Considering the complete planform of the shell a uniform mesh division of 8 × 8 has been used for the analyses.For all the cases shells are centrally impacted by a spherical steel ball of 0.0127 m diameter with an initial velocity of 3.0 m per second.The values of the contact stiffness coefficient, (k) and mass density of the impactor used in this study are 0.805 × 10 9 N/m 1.5 and 7960 kg/m 3 , respectively [6].The following material properties [6] are adopted for computation:

Transient response
Deflections at the impact point (L/2, b/2) are obtained and the optimum value of time step is chosen after performing convergence study.Parametric studies are conducted for symmetrically laminated bending stiff ([0 2 /± 30] S ) graphite-epoxy composite shells with respect to location of delamination, rotational speed and angle of twist.Transient response is obtained for cantilever shell and twisted shell of angle 30 • corresponding to non-dimensional speed of rotation, Ω = Ω/ω 0 = 0,1 using the converged value of time step of 1.0 µs, where Ω and ω 0 are the speed of rotation and fundamental natural frequency of non-rotating shell, respectively.In each case 25% delamination is considered irrespective of the location of delamination.The effect of rotational speed on contact force and central deflection for cantilever and twisted shells are shown in Figs 4a-5b wherein mid-plane delamination is considered with relative position of its centre line being 0.25 from the fixed end.Fundamental natural frequencies of cantilever and twisted shells are computed to be 87.3Hz and 40.0 Hz, respectively.It is found that contact force history does    not change with rotational speed during the loading cycle (Fig. 4a and 4b).The change in centrifugal stiffening (geometric stiffness) of the shell with increase of speed does not significantly affect the maximum contact force, despite the fact that higher initial stresses result at higher speed of rotation.Higher value of pretwist angle is found to intensify the contact force but no significant change with speed is observed even for unloading cycle.Rotational speed has a noticeable effect on contact force for cantilever shell during the unloading cycle.Contact duration of the unloading cycle for stationary cantilever shell is longer compared to that of rotating one with more fluctuations in character.This finding corroborates the fact that stiffening because of higher initial stresses reduces the duration of contact period [6].Deflection histories (Fig. 5a and 5b) for both the shells have no variation with speed during the loading cycle.The speed of rotation for cantilever shell is more and the displacement is due to the combined effect of steady state    centrifugal force and contact force.Hence, the peak value of normalized deflection for cantilever shell is found to be more in comparison with that of twisted shell.Of course the change of displacement with speed is observed during the unloading cycle, although variation is less for twisted shell with higher value in case of stationary condition as expected due to the softening effect compared to rotating condition.Earlier completion of unloading cycle for rotating cantilever shell leads to persistence of centrifugal effect in the later stages and as a result the central deflection shoots up considerably, but the stationary shell shows more or less steady values of displacement during the unloading cyle.The effect of relative position of delamination on contact force and central deflection for cantilever and twisted rotating shells (Ω = 1) are shown in Figs 6a-7b wherein delamination is considered with relative position its centre line being 0.25 from the free end, and the relative locations of delamination across thickness are 0.125, 0.5 and 0.75 considered from top surface.Relative location across the thickness has no sharp effect on the contact force history (Fig. 6a and 6b) as in contrast to natural frequency [12] which is the maximum when the delamination is located near the top or bottom surface and decreases gradually when the location moves towards the mid-plane which leads to the minimum value.In other words, it can be said that elastic stiffness also changes in the same manner as that of natural frequency with respect to location of delamination across the thickness.In conformity with this the peak values of the contact force in the later part of the unloading cycle of twisted shells are also found to occur, the maximum being for relative position at 0.125 and the minimum for relative location at 0.5.It is to be noted that for rotating cantilever shell contact period is completed earlier for the case of 0.75 compared to that of 0.125.Although the contact force is invariant irrespective of location of delamination for both twisted and cantilever rotating shells, but the absolute peak value of contact force is more in case of twisted shell.Like the contact force central deflection (Fig. 7a and 7b) is also influenced by the relative position of delamination in the later stages of unloading cycle, but no significant variation is noticed during the loading cycle.Twisted shells show decreasing trend of values of displacements with respect to relative location of delamination across thickness in order of 0.5, 0.75 and 0.125 as expected as the elastic stiffness also increases in this order.The similar trend is also found in case of rotating cantilever shells.Compared to the case of twisted shell in case of the cantilever one deflection builds up quickly when the contact period is over which is also less for cantilever shell.For twisted shells deflections shoot up earlier to attain the maximum value and thereafter die out gradually the absolute peak value during the span of analyses being smaller.In conformity with the contact force history in the later stages of the unloading cycle for the twisted shells it can be said that deflections will drop down at a faster rate which is not exhibited for cantilever shells.Natural frequencies of both cantilever and twisted delaminated shells have higher values when the delamination is located near the free surface (top or bottom) compared to its location at the mid-plane and accordingly the elastic stiffness also changes implying lower value corresponding to mid-plane delamination.This is also reflected in the displacement history during the unloading process and the maximum values of displacement are obtained for the mid-plane delamination.

Conclusions
A transient dynamic finite element method incorporating multipoint constraint algorithm is developed to study the normal impact response of delaminated composite rotating shallow shells.The finite element method is formulated based on Lagrange's equation of motion and employing an eight noded isoparametric plate bending element.Contact force is invariant with rotational speed during the loading cycle for both cantilever and twisted delaminated composite shells, but rotating cantilever shell because of higher initial stresses leads to reduction of unloading period.Rotational effect is pronounced in case of cantilever shell only and the centrifugal effect is considerably manifested in sharp increase of central deflection during the unloading cycle.Relative position of delamination across thickness has significant effect for both cantilever and twisted rotating shells during the unloading cycle only.Contact force is the maximum when the delamination is located near the top or bottom surface of the laminate and the minimum value of contact force is obtained for mid-plane delamination.Consequently, a reverse trend in displacement is also found.

Fig. 3 .
Fig. 3. (a) The contact force of centrally impacted cross ply ([0/90/0/90/0] S ) composite plate under simply supported boundary condition.(b) The influence of the relative position of the delamination on the first natural frequency of the composite cantilever beam.

Fig. 4 .
Fig. 4. (a) The effect of rotational speed on contact force for bending stiff ([0 2 /±30] S ) graphite-epoxy composite cantilever cylindrical shell with 25% mid-plane delamination (relative position of centre line at 0.25 from the fixed end).(b) The effect of rotational speed on contact force for bending stiff ([0 2 /±30] S ) graphite-epoxy composite twisted cylindrical shell with 25% mid-plane delamination (relative position of centre line at 0.25 from the fixed end).

Fig. 5 .
Fig. 5. (a) The effect of rotational speed on central deflection for bending stiff ([0 2 /±30] S ) graphite-epoxy composite cantilever cylindrical shell with 25% mid-plane delamination (relative position of centre line at 0.25 from the fixed end).(b) The effect of rotational speed on central deflection for bending stiff ([0 2 /±30] S ) graphite-epoxy composite twisted cylindrical shell with 25% mid-plane delamination (relative position of centre line at 0.25 from the fixed end).

Fig. 6 .
Fig. 6.(a) The effect of relative position of delamination across thickness on contact force for bending stiff ([0 2 /±30] S ) graphite-epoxy composite cantilever cylindrical shell with 25% delamination (relative position of centre line at 0.25 from the free end).(b) The effect of relative position of delamination across thickness on contact force for bending stiff ([0 2 /±30] S ) graphite-epoxy composite twisted cylindrical shell with 25% delamination (relative position of centre line at 0.25 from the free end).

Fig. 7 .
Fig. 7. (a) The effect of relative position of delamination across thickness on central deflection for bending stiff ([0 2 /±30] S ) graphite-epoxy composite cantilever cylindrical shell with 25% delamination (relative position of centre line at 0.25 from the free end).(b) The effect of relative position of delamination across thickness on central deflection for bending stiff ([0 2 /±30] S ) graphite-epoxy composite twisted cylindrical shell with 25% delamination (relative position of centre line at 0.25 from the free end).