# Study of Dynamic Behavior of Multilayered Clamped Composite Skewed Hypar Shell Roofs under Impact Load

Authors: Sanjoy Das Neogi; Amit Karmakar; Dipankar Chakravorty
Journal: Journal of Engineering
 (2013)

---

## Abstract

With advancement in the field of structural engineering, hunt for smarter materials has channelised the research towards the application of composite material. It is the high specific weight and specific stiffness of this material that have drawn the interest of different industrial sectors. Civil engineers also picked up composites to use it as a roofing material. Laminated composite shells, which can cover large column-free area and reduces dead weight of structure, show vulnerability under sudden impact due to their low transverse shear resistances. This study utilises finite element tool to investigate the dynamic response of a multilayered laminated composite hypar shells for fully clamped boundary condition. This class of shells is unique in a sense that the curvature has only the radius of cross curvature and these shells do not admit easy closed form solution particularly when the boundary conditions are complicated. Contact behavior of impactor and impacted mass has been modeled by modified Hertzian contact law and time-dependent equations are solved using Newmark’s time integration technique. Basic aim is to analyse the shell for symmetrically placed multilayered angle and cross ply lamination under different impact velocities.


---

## Body

## 1. Introduction

Composite material has gained popularity in different industrial sectors. It is the high specific weight and stiffness of this material that have gained the popularity for composite. With other sectors, civil engineering industries also started using the laminated composites as roofing material in the situation where large column-free space is required. The high specific weight also causes reduction in dead weight of structure and thus reduces seismic weight and foundation cost. In spite of these advantages, low transverse shear strength makes them vulnerable under sudden impact. It is obvious that an accurate modeling of contact behavior is one of the most important steps to analyze an impact response problem. The classical contact law between elastic solids derived by Hertz [1] was employed by different authors from time to time. The problem gets more complicated in case of composite materials and the Hertzian contact law which was derived for homogeneous isotropic materials may not be adequate. Tan and Sun [2] undertook an experimental program on a graphite/epoxy laminated plate to establish an empirical indentation law. The theoretical basis of their experiment was checked by a nine-noded plate finite element. Time histories of contact force and displacement were reported by Sun and Chen [3] for simply supported initially stressed plate under impact. The contact law as proposed by Tan and Sun [2] was utilized. They used a steel ball as an impactor and a ten-layered simply supported composite plate as the impacted mass. Impact analysis of shell structure was first reported by Toh et al. [4] for an orthotropic laminated cylindrical shell under low velocity impact generated by a solid striker. Shim et al. [5] studied the elastic response of glass/epoxy laminated composite ogival shells subjected to low velocity impact at any arbitrary location by a solid striker.They reported an analytic biharmonic polynomial solution. A finite element model, with and without geometric nonlinearity, was presented by Kistler and Waas [6] for a laminated composite cylindrical shell subjected to transverse central impact. Chun and Lam [7] proposed a numerical formulation to calculate the dynamic response of fully clamped laminated composite plates subjected to low velocity impact. They utilised Lagrange’s principle to derive the nonlinear, second order governing differential equation and Hertzian contact law. A centrally impacted three-layered laminated graphite-epoxy composite plate, of different ply orientation was investigated. Their study excluded any problem on shell forms. Karmakar and Kishimoto [8] undertook a transient dynamic finite element analysis to study the response of centrally impacted delaminated composite pretwisted rotating cylindrical shells under low velocity impact. Effects of transverse shear deformation and rotary inertia were included in their study. A transient dynamic response of rotating shallow shells subjected to low velocity impact was reported by Karmakar and Kishimoto [9].A look through the literature reveals the fact that impact response of civil engineering shell structures has not received due attention though such shells may often be subjected to impact loads. Impact forces on shell roofs are encountered due to snow fall, due to air borne debris in cyclone prone areas and in other similar situations. So the necessity to study impact response of composite shell roofs is felt. A parallel review in the area of composite shells reveals that the industrially important hypar shell needs a lot more attention although some important aspects of these shell forms were reported recently by Sahoo and Chakravorty [10–12]. The only report on impact response of composite hypar shell was due to Das Neogi et al. [13] where they discussed some aspects of such two-layered simply supported shells. Hence, this paper aims to carry out impact response of composite hypar shells for fully clamped boundary condition.


## 2. Mathematical Formulation

The shell surface considered for the present study is hypar shell (Figure1) with doubly curved and anticlastic surface having cross curvature 1/Rxy only and 1/Rx=1/Ry=0. A doubly curved thin shallow shell of uniform thickness h made of homogeneous, laminated composite, linearly elastic material is considered. A shell is characterized as shallow if any infinitesimal line element of its middle surface is approximated by the length of its projection on the XY plane. This implies that
(1)(∂z∂x)2≪1,(∂z∂y)2≪1,(∂z∂x)(∂z∂y)≪1.
Moreover, the lateral boundary of a shallow shell is approximated by its projection on the XY plane with regard to its boundary conditions. According to Vlasov [14], the above conditions are practically satisfied for shells with a rise to span ratio less than 1/5 and cross curvature is approximately represented as
(2)1Rxy=∂2z∂x∂y.Figure 1
Surface of a skewed hypar shell and degrees of freedom.The generalised displacement vector of an element is expressed in terms of the shape functions and nodal degrees of freedom as(3)[u]=[N]{de},
that is,
(4){u}={uvwαβ}=[[Ni][Ni][Ni][Ni][Ni]]{uiviwiαiβi}.The element stiffness and mass matrices are derived by using the minimum energy principle. The element stiffness matrix is(5)[Kshe]=∬[B]T[D][B]dxdy.The generalised inertia matrix per unit area includes the translatory and rotatory inertia terms. Incorporating both the translatory and rotatory inertia terms, the generalised inertia matrix takes the following form:(6)[M]=∬[N]T[ρ][N]dxdy.The dynamic equilibrium equation of the target shell for low velocity impact is given by the following equation:(7)[M]{δ..}+k{δ}={F},
where [M] and [K] are global mass and elastic stiffness matrices, respectively. {δ} is the  global displacement vector. For the impact problem, {F} is given as
(8){F}={000⋯Fc⋯000}T.HereFc  is the contact force given by the indentation law and the equation of motion of the rigid impactor is given  as
(9)mi..ωi..+Fc=0,
where mi.. and ωi..  are  the mass and acceleration of the impactor, respectively. A power law was proposed by Yang and Sun [15] based on static indentation tests using steel ball as an indentor. This contact law accounted for the permanent indentation after unloading cycle; that is, collisions upon the rebound of the target structure after the first period of contact were considered. The modified version of the above mentioned contact law was proposed by Tan and Sun [2] and was utilised by Sun and Chen [3]. The contact force model following Sun and Chen [3] is incorporated in the present finite element formulation. If k is the contact stiffness and αm is the maximum local indentation, the contact force Fc during loading is given by
(10)Fc=kα1.5,0<α≤αm.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. At an instant, the maximum indentation takes place and as a result maximum contact force is also obtained. At this instant, displacement of the impactor also attains the maximum value [16]. There after the displacement of the impactor gradually decreases, but the target point displacement keeps on changing and finally increases to a maximum value and some point of time these two displacements become equal [16]. This leads to zero value of indentation and eventually the contact force becomes zero. At this instant, the impactor loses the 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 [3]. 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. This process is known as reloading. If Fm is the maximum contact force at the beginning of unloading and αm is the maximum indentation during loading, the contact force Fc  for unloading and reloading is expressed as [3]unloading phase:(11)Fc=Fm[α-α0αm-α0]2.5;reloading phase:(12)Fc=Fm[α-α0αm-α0]1.5,where α0  denotes the permanent indentation in a loading-unloading cycle.
(13)α0=βc(αm-αp)ifαm>αcrα0=0ifαm<αcr,
where βcis a material dependent constant and  αcr  is the critical indentation beyond which permanent indentation occurs, and  the values are 0.094 and .01667 cm, respectively for graphite-epoxy composite [3]. Equations (7) and (9) are solved using Newmark constant-acceleration time integration algorithm in the present analysis. Equation (7) may be expressed in iteration form at each time step as follows:.
(14)[K-]{Δ}t+Δti+1=Δt24{F}t+Δti+[M]{b}i,
where
(15)[K-]=Δt24[K]+[M](16){b}t={Δ}t+Δt{Δ.}t+Δt24{Δ..}t.In the above equations (14) and (15), [K] and [M] are the global stiffness and mass matrices obtained after assembling the matrices obtained at element level and imposing the boundary conditions. The degrees of freedom that are locked are selected through an input matrix and the corresponding terms in the global stiffness and mass matrices and load vectors are eliminated. The approach enables imposition of arbitrary boundary conditions in the problem although the present results are obtained only for clamped edges.  The same solution scheme is also utilized for solving the equation of motion of the impactor, that is, (7). In (14), i is the number of iterations within a time step. It is to be noted that a modified contact force Ft+Δti obtained from the previous iteration is used to solve the current response  {Δ}t+Δti+1. The iteration procedure is continued until the equilibrium criterion is met.


## 3. Numerical Examples

Numerical problems are solved with two different objectives. The present approach is applied to solve natural frequencies of graphite-epoxy twisted plates which are structurally similar to skewed hypar shells. This is to validate both the stiffness and mass matrix formulations in present finite element code. The problem details are furnished with Table1. Another problem solved earlier by Chun and Lam [7] related to the impact response of clamped composite plate, is taken up as the second benchmark to validate the impact formulation and the results are furnished as Figure 2.Table 1
Nondimensional natural frequencies(ϖ) for a three-layer graphite epoxy twisted plate.
Angle of twistθ  (deg)
0
15
30
45
60
75
90
ϕ
=
15
°
Qatu and Leissa [17]
1.0035
0.9296
0.7465
0.5286
0.3545
0.2723
0.2555
Present formulation
0.9990
0.9257
0.7445
0.5279
0.3542
0.2720
0.2551
ϕ
=
30
°
Qatu and Leissa [17]
0.9566
0.8914
0.7205
0.5149
0.3443
0.2606
0.2436
Present formulation
0.9490
0.8842
0.7181
0.5142
0.3447
0.2613
0.2444
E
11
=
138 GPa, E22=8.96 GPa, G12=7.1 GPa, ν12=0.3, a/b=1, a/h=100.Figure 2
Contact force history of clamped plate.E11=142.73 GPa, E22=13.79 GPa, G12=4.64 GPa, ν12=0.30, ρ=1.61×103 kg/m3, a=b=0.14 m, h=3.81×10-3 m, mass of striker =0.014175 kg, velocity of striker  =22.6 m/s, contact stiffness (kc)=1×108 N/m1.5.Apart from the problems mentioned above, impact response of skewed hypar shells being impacted at the central point is also studied for clamped boundary condition, different laminations, and impact velocities. The details of the problems which are the authors’ own are given below.(i)
Boundary condition: Clamped (CL).(ii)
Lamination:00/900/00(CP).(iii)
Velocity of impact in (m/s): 1, 2, 3, 5, 7, 10.(iv)
Details of shell geometry:a=1.0 m, b=1.0 m, h=0.02 m, c=0.2 m.(v)
Material details:E11=120 GPa, E22=7.9 GPa, G12=G23=G13=5.5 GPa ν12=0.30, ρ=1.58×10-5 N-sec2/cm4.The converged value (Table2) for a time step Δt=3 μs is adopted for the present analysis.Table 2
Convergence study for time step.
Time step
Maximum contact force (N)
(µs)
0°/90°/0°
2
291.38795
2.5
301.58181
2.75
311.0179
3
318.5276
3.25
324.95893
3.5
329.38339


## 4. Results and Discussions

The results of Table1 show that the fundamental frequency values of the twisted plates obtained by the present formulation agree closely to those reported by Qatu and Lessia [17]. Thus the correct incorporation of stiffness and mass matrix formulation in the present code is established. Figure 2 shows the time variation of the contact force induced in a composite plate under low velocity impact as obtained by Chun and Lam [7]. The values obtained by the present formulation are also presented graphically in the same figure in a different style. Here again excellent agreement of results is observed which establishes the correctness of impact formulation. The displacement history and dynamic magnification of several points (Figure 3) are studied along with the node of impact. For studying the impact response of clamped (CL) cross ply (CP) shell Figures 4, 5, 6, 7, 8, 9, and 10 and Table 3 are to be studied. All the results of contact force and displacement that are presented, in either graphical or tabular form, are arrived at after a study of time step convergence.Table 3
Maximum contact force, maximum dynamic displacement, equivalent static load, dynamic magnification factor for different velocities.
Boundary condition and ply orientation
Velocity (m/s)
Maximum impact load (N)
Maximum displacement (m)
Equivalent static load  (N)
Dynamic magnification factor
Clamped0°/90°/0°
1
318.5276
0.0000534
1880
3.75
2
736.1553
0.000111
3.910
3.46
3
1204.515
0.00017
5.980
2.90
5
2242.325
0.000292
10300
2.75
7
3373.736
0.000417
14700
2.50
10
5198.173
0.000612
21500
2.24Figure 3
Points at which deflection is measured.Impact response of clamped cross ply (CL/CP) composite hypar shells for impact velocity 1 m/s.
(a)
(b)
(c)
(d)Impact response of clamped cross ply (CL/CP) composite hypar shells for impact velocity 2 m/s.
(a)
(b)
(c)
(d)Impact response of clamped cross ply (CL/CP) composite hypar shells for impact velocity 3 m/s.
(a)
(b)
(c)
(d)Impact response of clamped cross ply (CL/CP) composite hypar shells for impact velocity 5 m/s.
(a)
(b)
(c)
(d)Impact response of clamped cross ply (CL/CP) composite hypar shells for impact velocity 7 m/s.
(a)
(b)
(c)
(d)Impact response of clamped cross ply (CL/CP) composite hypar shells for impact velocity 10 m/s.
(a)
(b)
(c)
(d)Variation of maximum impact load, maximum displacement, equivalent static load, and dynamic magnification factor with velocity for simply supported angle ply (CL/CP) composite hypar shells.
(a)
(b)
(c)
(d)The finite element mesh adopted is also based on force and displacement convergence criteria. When low velocity normal impact response of clamped cross ply shell is studied being struck by the spherical impactor centrally it is observed that the contact force shows a sort of parabolic variation with a single peak. After a given time span of 100μs or less, the contact force dies down to a null value. It is interesting to note that higher the impactor velocity higher is the contact force as expected, but the force dies down to a null value earlier. This behavior may be attributed to the fact that higher the velocity, more rapid is the elastic rebound of the impactor followed by detachment which causes the contact force to decay out. It is also very interesting to observe that the time instant corresponding to peak contact force and that for peak displacement do not match. This is because the resultant displacement at any time instant is a cumulative effect of the instantaneous contact force value and the inertia effect of the previous instant. The figures showing the transient displacements reflect the fact that vibration continues even after the force dies down with successively occurring peaks though, the peak values are less in magnitude than the highest peak which occurs a bit after the instant of maximum contact force but before the full decay of it.To estimate the equivalent static load (ESL) corresponding to a particular impactor velocity, a concentrated load at the centre (point of impact) is applied and adjusted to yield a central displacement equal to the maximum dynamic displacement. The magnitude of the central displacement is calculated with the peak contact force applied at the point of impact as a static concentrated load. The central displacement obtained under such a load when divides the maximum dynamic displacements yields dynamic magnification factor (DMF). The variations of maximum contact force, maximum dynamic displacement, and equivalent static load (ESL) with impact velocity are almost linear and all of the three above mentioned values are increasing functions of the impact velocity. However the dynamic magnification-factor (DMF) and the impact velocity show a logarithmic relation and the DMF is a decreasing function of velocity of impact.It may also be noted that at some of the points, for example, point-2, the displacement at some instant may exceed the same occurring at the point of impact. But such displacements will never govern the design as they are always less than the absolute maximum displacement occurring at the point of contact.


## 5. Conclusions

The following conclusions may be derived from the present study.(1)
The close agreement of the results obtained by the present method with those available in the published literature establishes the correctness of the approach used here.(2)
Under the influence of normal low velocity impact, the contact force shows a parabolic combined loading and unloading curve with a single peak for the practical class of shells considered here. Higher magnitude of impact velocity results in higher value of the peak contact force but due to a sharp elastic rebound, the total duration of contact force is less for higher velocity of impactor.(3)
The time instants at which the maximum contact force and the maximum dynamic displacement occur show a phase difference and interestingly in some cases, the maximum displacement and hence stresses may occur even after the contact force dies down totally. Thus it is concluded that the study may be stopped only after when the major peaks of the dynamic displacement die down and not after the full decay of the contact force only.(4)
The maximum contact force, the peak dynamic displacement, and the equivalent static load are all increasing functions of impactor velocity, the relations being almost linear. However, the dynamic magnification factor shows a logarithmically decreasing tendency with increase of the velocity of impact.(5)
At some of the points, the displacement at some instant may exceed the same occurring at the point of impact. But such displacements will never govern the design as they are always less than the absolute maximum displacement occurring at the point of contact.



---

*Source: 192176-2013-03-31.xml*