# Forced Vibration Analysis of Laminated Composite Shells Reinforced with Graphene Nanoplatelets Using Finite Element Method

Authors: Trung Thanh Tran; Van Ke Tran; Pham Binh Le; Van Minh Phung; Van Thom Do; Hoang Nam Nguyen
Journal: Advances in Civil Engineering
 (2020)

---

## Abstract

This paper carries out forced vibration analysis of graphene nanoplatelet-reinforced composite laminated shells in thermal environments by employing the finite element method (FEM). Material properties including elastic modulus, specific gravity, and Poisson’s ratio are determined according to the Halpin–Tsai model. The first-order shear deformation theory (FSDT), which is based on the 8-node isoparametric element to establish the oscillation equation of shell structure, is employed in this work. We then code the computing program in the MATLAB application and examine the verification of convergence rate and reliability of the program by comparing the data of present work with those of other exact solutions. The effects of both geometric parameters and mechanical properties of materials on the forced vibration of the structure are investigated.


---

## Body

## 1. Introduction

Nowadays, due to the development of science and technology, a variety of new materials have been applied widely to engineering applications such as composite materials, functionally graded materials (FGM), and piezoelectric materials (PZT). To enhance the strengths of structures, there are many common ways like adding stiffeners, using folded structures, or reinforcing the structures with smart materials. In these techniques, the carbon nanotube (CNT)-reinforced composite shell is one of the most modern structures [1]. However, in comparison with nanographene or graphene nanoplatelet (GPL) structures, carbon nanotube (CNT) structures are not only inflexible but also have a high cost [2]. In recent years, GPL reinforcement for composite materials (GPLRC) has been manufactured by Stankovich et al. [3–5]. The investigations found that only a small amount of GPL in the material could significantly improve contemporaneously its physical characteristics [6, 7]. In comparison with conventional composite materials containing more than 60% carbon fiber by the volume, GPLRC contains only a low graphene nanoplatelet (GPL) ratio (0.5%–20% weight) [8–10]. Therefore, GPLRC is concerned by many scientists all over the world. In civil engineering, cylindrical shells (CYL) and spherical shell (SPH) reinforced with GPL are considered with many different approaches.In [11], nonlocal continuum theory and molecular dynamics simulations were used to discover the free vibration of single and double-layer graphene structures. More works dealt with nonlocal methods. Ghannadpour et al. [12] performed bending, buckling, and vibration analysis of Euler beams, and the same mechanical response analysis of their work [13] for Timoshenko beam was also carried out. Chandra et al. [14] investigated natural frequencies and oscillation mode shapes of GRP composite plates with several types of boundary supports using numerical approaches. The linear buckling and nonlinear buckling of multilayer beam structures reinforced by GPL with the Halpin–Tsai model were studied by Yang and his co-workers [15], where effects of material characteristics of each component on the vibration of the structure are examined. Feng and his colleagues [16] developed a relatively complete investigation of nonlinear static bending of multilayer beams reinforced with graphene. Then, he [17] extended his study to carry out the nonlinear free vibration problems. The material characteristics calculated by the Halpin–Tsai model and the von Kármán nonlinear deformation-displacement were used in both studies. Kitipornchai and colleagues [18] studied free vibration and stabilization of FGM beams with reinforced holes in the graphene. Lin et al. [19] carried out the investigation of mechanical characteristics of grapheme-reinforced components. In [20], Zhao and his co-workers performed the static bending and free vibration analysis of FGM plates reinforced with nanographene by using the newly proposed Halpin–Tsai model. Mirzaei and Kiani [21] examined the thermal buckling of composite structures (plate structures) reinforced with graphene based on FSDT where the volume of graphene in each component layer was different from one another.Shen and his co-workers [22] carried out the nonlinear bending, buckling, and postbuckling analysis of composite beams reinforcing graphene on an elastic foundation or in a thermal environment. Shen et al. [23, 24] introduced the results of nonlinear bending and vibration analysis of FGM plate-reinforced graphene under the effect of thermal load. Then, they [25, 26] developed the linear and nonlinear thermostats of composite structures reinforcing graphene resting on an elastic foundation employing high-order shear deformation theory (HSDT). In [27], the author used both the Von Kármán-type nonlinearity and the third-order shear deformation theory (TSDT) to carry out the vibration analysis of graphene-reinforced cylinders resting on the elastic foundation (in the thermal environment) [28]. For buckling problems, Wang and colleagues [29] introduced the buckling analysis of cylindrical composite shells using FEM based on the Halpin–Tsai model and mixed rule in order to determine the mechanical characteristics of each component layer. Then, the buckling problem of composite CYL-reinforced graphene with holes and FGM CYL reinforcing graphene was presented in [30, 31]. They [32] also studied buckling of dynamic multilayer nanocomposite beams under the thermal load. Wu et al. [33] carried out the parameter buckling of the grapheme-reinforced multilayer nanocomposites under both thermal and mechanical loads. Kiani [34] examined the large deformation of the composite plate with reinforced graphene in which the masses of graphene in each component were divergent. Rout et al. [35] studied on nonlinear vibration of composite shell types (CYL, SPH, HYP, and so on) reinforced with graphene using higher-order shear deformation theory (HSDT). Recently, Duan et al. [36] introduced a novel formulation by employing the discrete singular convolution (DSC) for free vibration analysis of circular thin plates with uniform and stepped thickness. Civalek [37] investigated the free oscillation of CNT and FGM structures. Nguyen et al. [38, 39] computed dynamic responses of the composite shell with shear connectors using FEM.In the above studies, to our knowledge, most plate and shell structures are reinforced with graphenes, where a few of which dealt with graphene-reinforced composite shells. In addition, the forced vibration problems of graphene-reinforced composite shells with the different mass fractions of graphene of each component in thermal environments are still limitations. The main goals of this work are to analyze the forced vibration of the shell structures (Figure1) with reinforced graphene performing in a thermal environment using FEM. Herein, we consider four models with different distributions of graphene including UD-type, X-type, V-type, and O-type to study in this paper. For UD structure, the dispensation rule of graphene shells is homogeneous in each component layer. On the other hand, for X structure, the graphene plates are found to be minimal in the middle surface and gradually increase to the top and bottom layers. For V-type structure, the volume of each part of graphene is detected as maximum in the top layer and minimum in the bottom layer. For O structure, the volume fraction of graphene is maximal in middle surface of the structure and minimum in the top layer and the bottom layer. The temperature-dependent material characteristics of the multilayer structure are evaluated by using the Halpin–Tsai extension method. In this paper, the authors investigated the influences of both geometrical and physical properties on the forced vibration with two cases of shells (cylindrical shell (CYL) and spherical shell (SPH)) in the thermal environment.Figure 1
Composite shell with 4 graphene reinforcement forms.The organization of this paper is divided into 4 main sections. Section2 presents briefly the equation of motion of the shell element by using finite element method (FEM). Numerical results and discussions are shown in Section 3. Some highlighted results are concluded in Section 4.


## 2. Finite Element Formulations

### 2.1. Equation of Motion of the Shell Element

Consider a curvilinear shell model reinforced graphene nanoplatelets with the following dimensions: lengtha, width b, and thickness h, as shown in Figure 1.R
x
and
R
y denote the radii of the shells:(i)
CYL:Rx≠0,Ry=∞ (Case 1)(ii)
SPH:Rx;Ry≠0 (Case 2)The displacement field is defined through FSDT as follows:(1)q=ux,y,z,t=u0x,y,0,t+zφxx,y,0,t,vx,y,z,t=v0x,y,0,t+zφyx,y,0,t,wx,y,z,t=w0x,y,0,t,where u,v, and w are the displacement components at any point x,y,z of the shell; u0,v0, and w0 are the displacement components at the midsurface of the shell; and φx and φy are the rotation angles of the cross section of the structure around the y-axis and x-axis.The stress-strain relationship is expressed as follows:(2)εx=∂u∂x=∂u0∂x+w0Rx+z∂φx∂x;εy=∂v∂y=∂v0∂y+w0Ry+z∂φy∂y;γxy=∂v∂x+∂u∂y=∂u0∂y+∂v0∂x+zk∂φx∂y+∂φy∂x+121Ry−1Rx∂v0∂x−∂u0∂y;γxz=∂w0∂x+∂u∂z=∂w0∂x+φx−u0Rx;γyz=∂w0∂y+∂v∂z=∂w0∂y+φy−u0Ry;Equation (2) can be shortened in the matrix form as follows:(3)ε=εxεyγxy=ε0+zκ;γ=γyzγzx,in which(4)ε0=εx0εy0γxy0=∂u0∂x+w0Rx∂v0∂y+w0Ry∂u0∂y+∂v0∂x+2w0Rxy;κ=κxκyκxy=∂φx∂x∂φy∂y∂φx∂y+∂φy∂x+121Ry−1Rx∂v0∂x−∂u0∂y,γ=γxzγyz=−u0Rx−v0Rxy+∂w0∂x+φx−u0Rxy−v0Ry+∂w0∂y+φy.The stress and strain fields of layerk has a relationship as follows:(5)σk=Dkεk−εkT;τk=56Gkγk,where DkandGk are the bending rigidity matrix and shear rigidity matrix of layer k:(6)Qk=Q¯11Q¯12Q¯16Q¯12Q¯22Q¯26Q¯16Q¯26Q¯66k=m2n2−2mnn2m22mnmn−mnm2−n2kQ11Q120Q12Q22000Q66km2n2−2mnn2m22mnmn−mnm2−n2kT,Gk=Q¯44Q¯45Q¯45Q¯55k=m−nnmkQ5500Q44km−nnmkT,εkT=T−T0α=ΔTαxαyαxy=ΔTm2α1+n2α2n2α1+m2α22mnα1−α2,with(7)m=cosθ;n=sinθ,Q11=E11−ν12ν21,Q22=E21−ν12ν21,Q12=ν12E21−ν12ν21,Q66=G12,Q55=G23,Q44=G13,where θ is the fiber angle of the k layer; α1andα2 are coefficients of thermal expansion of the k layer in the x and y directions. Herein, we employ the extended Halpin–Tsai method to estimate the effective elastic modulus as well as shear modulus of the structures which are reinforced by using graphene nanoplatelets. In the present considerations, we choose two options of reinforcement of graphene nanoplatelets, which are zigzag (θ = 0°) and armchair (θ = 90°). The operative material characteristics of composite structures reinforced with graphene are presented as [27](8)E1=η11+2aG/hGλ111−λ11VGEm,E2=η21+2aG/hGλ221−λ22VGEm,G12=η311−λ12VGGm,where ηi (i = 1–3) represent graphene proficiency parameters calculated by comparing the results of the work in [27]. Em and Gm are, respectively, elastic modulus and shear modulus. The dimensions of the graphene sheet include length aG, width bG, and effective thickness hG. The volume component of graphene is VG. The coefficients λij are defined as(9)λ11=η1E11G/Em−1E11G/Em+2aG/hG,λ22=η1E22G/Em−1E22G/Em+2bG/hG,λ12=G12G/Em−1G12G/Em,where E11G, E22G, and G12G are Young’s moduli and shear modulus of the graphene structures. The density of the mass and Poisson’s ratio of the composite structures are defined as(10)ν12=VGν12G+Vmνm,ρ=VGρG+Vmρm,where(i)
Vm is the matrix volume fraction(ii)
G12 and νm are, respectively, Poisson’s ratio of the grapheme and matrix(iii)
ρm are ρG are, respectively, the mass density of the graphene and the matrixThe parameters of the thermal development along the longitudinal and transverse directions of composite structures reinforced with graphene can be defined as follows:(11)α1=VGE11Gα11G+VmEmαGVGE11G+VmEm,α2=1+ν12GVGα22G+1+νmVmαm−ν12α1.In this work, we employ the 8-node isoparametric element, where each node consists of five degrees of freedom (DOF); the degree of freedom of nodei is qei, and the total DOFs of the shell element qe is expressed as follows:(12)qei=u0v0w0φxφyT;i=1,8¯,qe=qe1qe2qe3qe4qe5qe6qe7qe8T,u0=∑i=18Niξ,ηu0i,v0=∑i=18Niξ,ηv0i,w0=∑i=18Niζ,ηw0i,φx=∑i=18Niξ,ηφxi,φy=∑i=18Niξ,ηφyi,in which Nii=1,8¯ can be found in Ref. [39].By substituting in the expression for verifying displacement of element, we have(13)ε=B1+zB2qe,γ=B3qe,where B1,B2,andB3 are defined as follows:(14)B1=B11B21B31B41B51B61B71B81;B2=B12B22B32B42B52B62B72B82;B3=B13B23B33B43B53B63B73B83;where Bi1, Bi2, and Bi3 can be found in Appendix A.To obtain the dynamic equation, we employ the weak form for each element, and then we obtain(15)∫Vδq˙eTρq˙edV+∫Vδε−εTTσdV+∫VδγTτdV−δqeT∫ANwptdA=0.By substituting equations (1) and (13) into equation (15), we have the dynamic equation of the shell element as follows:(16)Meq¨e+Keqe=Fet+FeTt,with(17)Ke=∫AB1TD0B1dA+∫AB1TD1B2dA+∫AB2TD1B1dA+∫AB2TD2B2dA+∫AB3TD3B3dA,where(18)D0;D1;D2=∫−h/2h/21;z;z2Qdz;D3=56∫−h/2h/2Gdz,Me=∫ANTρ∫−h/2h/2LTLdzNdA,where NandL can be found in Appendix B.(19)Fet=∫AptNwTdA,FeTt=ΔT∫AB1TA+B2TBαdA,in which(20)Nw=Nw1Nw2Nw3Nw4Nw5Nw6Nw7Nw8,Nwi=00Ni00.For the case of taking into account the structural damping ratio, we obtain the forced oscillation connection of the structure element as follows:(21)Meq¨e+Ceq˙e+Keqe=Fet+FeTt,in which Ce=αMe+βKe and αandβ are Rayleigh drag coefficients, which are calculated as in the literature [40].


### 2.2. The Differential Equation of Vibration

From the differential equation of forced vibration of the shell element (21), we obtain the differential equation of forced vibration of entire composite shell structure as follows:(22)Mq¨+Cq+Kq=Ft+FTt,in which M, C, K, Ft, and FTt are, respectively, the global mass matrix, the global structural damping matrix, the global stiffness matrix, the global load vector, and the global thermal load vector. These matrices and vectors are assembled from the element matrices and vectors, correspondingly. These are linear differential equations, where the right-hand side of the equation varies over time. To analyze these equations, we employ the Newmark-beta approach [40]. The program is coded in the MATLAB (MathWorks, Natick, MA, USA) environment with the algorithm flowchart of Newmark as shown in Figure 2.Figure 2
Flowchart of Newmark-beta method for solving the dynamic response problem of the shell.WhenC=0 and Ft+FTt=0, from equation (24), natural frequencies can be obtained by figuring out the roots of the following equation:(23)Mq¨+Kq=0.By assuming thatq=q0sinωt≠0, we obtain(24)K−ω2Mq0sinωt=0⟶K−ω2M=0,where ω is the natural frequency.Now, we have the flowchart of the Newmark-beta method with symbols “n” is the total number of integral steps and “j” is the number of integral step.Step 1.
Determine the initial conditions(25)q0=q0;q˙0=q˙0.
From the initial conditions, we obtain(26)q¨0=M0−1F0−K0q0−C0q˙0.Step 2.
By approximatingq¨t+Δt,q˙t+Δt by qt+Δt, we have(27)q¨t+Δt=a0qt+Δt−qt−a2q˙t−a3q¨t,q˙t+Δt=q˙t+a6q¨t+a7q¨t+Δt,where(28)a0=2γΔt2;a1=2αγΔt;a2=2γΔt;a3=1γ−1;a4=2αγ−1;a5=αγ−1Δt;a6=1−αΔt;a7=αΔt,in which αandγ are defined by the assumption that the acceleration varies in each calculating step; the authors select the linear law for the varying acceleration:(29)q¨τ=q¨t+τΔtq¨t+Δt−q¨twitht≤τ≤t+Δtthenα=12;γ=13.
The condition to stabilize the roots is expressed as follows:(30)Δt≤12ϖmax1α−γorΔtTmin≤12π21α−γ.Step 3.
Compute the efficiency stiffness matrix and the efficiency nodal force vector:(31)K∗=K+a0M+a1C,F∗=Ft+Δt+Ma0qt+a2q˙t+a3q¨t+Ca1qt+a4q˙t+a5q¨t.Step 4.
Determine the nodal displacement vectorqt+Δt:(32)Kt+Δt∗qt+Δt=Ft+Δt∗⟹qt+Δt=Kt+Δt∗−1Ft+Δt∗,and then repeat the loop until the time runs out.


## 2.1. Equation of Motion of the Shell Element

Consider a curvilinear shell model reinforced graphene nanoplatelets with the following dimensions: lengtha, width b, and thickness h, as shown in Figure 1.R
x
and
R
y denote the radii of the shells:(i)
CYL:Rx≠0,Ry=∞ (Case 1)(ii)
SPH:Rx;Ry≠0 (Case 2)The displacement field is defined through FSDT as follows:(1)q=ux,y,z,t=u0x,y,0,t+zφxx,y,0,t,vx,y,z,t=v0x,y,0,t+zφyx,y,0,t,wx,y,z,t=w0x,y,0,t,where u,v, and w are the displacement components at any point x,y,z of the shell; u0,v0, and w0 are the displacement components at the midsurface of the shell; and φx and φy are the rotation angles of the cross section of the structure around the y-axis and x-axis.The stress-strain relationship is expressed as follows:(2)εx=∂u∂x=∂u0∂x+w0Rx+z∂φx∂x;εy=∂v∂y=∂v0∂y+w0Ry+z∂φy∂y;γxy=∂v∂x+∂u∂y=∂u0∂y+∂v0∂x+zk∂φx∂y+∂φy∂x+121Ry−1Rx∂v0∂x−∂u0∂y;γxz=∂w0∂x+∂u∂z=∂w0∂x+φx−u0Rx;γyz=∂w0∂y+∂v∂z=∂w0∂y+φy−u0Ry;Equation (2) can be shortened in the matrix form as follows:(3)ε=εxεyγxy=ε0+zκ;γ=γyzγzx,in which(4)ε0=εx0εy0γxy0=∂u0∂x+w0Rx∂v0∂y+w0Ry∂u0∂y+∂v0∂x+2w0Rxy;κ=κxκyκxy=∂φx∂x∂φy∂y∂φx∂y+∂φy∂x+121Ry−1Rx∂v0∂x−∂u0∂y,γ=γxzγyz=−u0Rx−v0Rxy+∂w0∂x+φx−u0Rxy−v0Ry+∂w0∂y+φy.The stress and strain fields of layerk has a relationship as follows:(5)σk=Dkεk−εkT;τk=56Gkγk,where DkandGk are the bending rigidity matrix and shear rigidity matrix of layer k:(6)Qk=Q¯11Q¯12Q¯16Q¯12Q¯22Q¯26Q¯16Q¯26Q¯66k=m2n2−2mnn2m22mnmn−mnm2−n2kQ11Q120Q12Q22000Q66km2n2−2mnn2m22mnmn−mnm2−n2kT,Gk=Q¯44Q¯45Q¯45Q¯55k=m−nnmkQ5500Q44km−nnmkT,εkT=T−T0α=ΔTαxαyαxy=ΔTm2α1+n2α2n2α1+m2α22mnα1−α2,with(7)m=cosθ;n=sinθ,Q11=E11−ν12ν21,Q22=E21−ν12ν21,Q12=ν12E21−ν12ν21,Q66=G12,Q55=G23,Q44=G13,where θ is the fiber angle of the k layer; α1andα2 are coefficients of thermal expansion of the k layer in the x and y directions. Herein, we employ the extended Halpin–Tsai method to estimate the effective elastic modulus as well as shear modulus of the structures which are reinforced by using graphene nanoplatelets. In the present considerations, we choose two options of reinforcement of graphene nanoplatelets, which are zigzag (θ = 0°) and armchair (θ = 90°). The operative material characteristics of composite structures reinforced with graphene are presented as [27](8)E1=η11+2aG/hGλ111−λ11VGEm,E2=η21+2aG/hGλ221−λ22VGEm,G12=η311−λ12VGGm,where ηi (i = 1–3) represent graphene proficiency parameters calculated by comparing the results of the work in [27]. Em and Gm are, respectively, elastic modulus and shear modulus. The dimensions of the graphene sheet include length aG, width bG, and effective thickness hG. The volume component of graphene is VG. The coefficients λij are defined as(9)λ11=η1E11G/Em−1E11G/Em+2aG/hG,λ22=η1E22G/Em−1E22G/Em+2bG/hG,λ12=G12G/Em−1G12G/Em,where E11G, E22G, and G12G are Young’s moduli and shear modulus of the graphene structures. The density of the mass and Poisson’s ratio of the composite structures are defined as(10)ν12=VGν12G+Vmνm,ρ=VGρG+Vmρm,where(i)
Vm is the matrix volume fraction(ii)
G12 and νm are, respectively, Poisson’s ratio of the grapheme and matrix(iii)
ρm are ρG are, respectively, the mass density of the graphene and the matrixThe parameters of the thermal development along the longitudinal and transverse directions of composite structures reinforced with graphene can be defined as follows:(11)α1=VGE11Gα11G+VmEmαGVGE11G+VmEm,α2=1+ν12GVGα22G+1+νmVmαm−ν12α1.In this work, we employ the 8-node isoparametric element, where each node consists of five degrees of freedom (DOF); the degree of freedom of nodei is qei, and the total DOFs of the shell element qe is expressed as follows:(12)qei=u0v0w0φxφyT;i=1,8¯,qe=qe1qe2qe3qe4qe5qe6qe7qe8T,u0=∑i=18Niξ,ηu0i,v0=∑i=18Niξ,ηv0i,w0=∑i=18Niζ,ηw0i,φx=∑i=18Niξ,ηφxi,φy=∑i=18Niξ,ηφyi,in which Nii=1,8¯ can be found in Ref. [39].By substituting in the expression for verifying displacement of element, we have(13)ε=B1+zB2qe,γ=B3qe,where B1,B2,andB3 are defined as follows:(14)B1=B11B21B31B41B51B61B71B81;B2=B12B22B32B42B52B62B72B82;B3=B13B23B33B43B53B63B73B83;where Bi1, Bi2, and Bi3 can be found in Appendix A.To obtain the dynamic equation, we employ the weak form for each element, and then we obtain(15)∫Vδq˙eTρq˙edV+∫Vδε−εTTσdV+∫VδγTτdV−δqeT∫ANwptdA=0.By substituting equations (1) and (13) into equation (15), we have the dynamic equation of the shell element as follows:(16)Meq¨e+Keqe=Fet+FeTt,with(17)Ke=∫AB1TD0B1dA+∫AB1TD1B2dA+∫AB2TD1B1dA+∫AB2TD2B2dA+∫AB3TD3B3dA,where(18)D0;D1;D2=∫−h/2h/21;z;z2Qdz;D3=56∫−h/2h/2Gdz,Me=∫ANTρ∫−h/2h/2LTLdzNdA,where NandL can be found in Appendix B.(19)Fet=∫AptNwTdA,FeTt=ΔT∫AB1TA+B2TBαdA,in which(20)Nw=Nw1Nw2Nw3Nw4Nw5Nw6Nw7Nw8,Nwi=00Ni00.For the case of taking into account the structural damping ratio, we obtain the forced oscillation connection of the structure element as follows:(21)Meq¨e+Ceq˙e+Keqe=Fet+FeTt,in which Ce=αMe+βKe and αandβ are Rayleigh drag coefficients, which are calculated as in the literature [40].


## 2.2. The Differential Equation of Vibration

From the differential equation of forced vibration of the shell element (21), we obtain the differential equation of forced vibration of entire composite shell structure as follows:(22)Mq¨+Cq+Kq=Ft+FTt,in which M, C, K, Ft, and FTt are, respectively, the global mass matrix, the global structural damping matrix, the global stiffness matrix, the global load vector, and the global thermal load vector. These matrices and vectors are assembled from the element matrices and vectors, correspondingly. These are linear differential equations, where the right-hand side of the equation varies over time. To analyze these equations, we employ the Newmark-beta approach [40]. The program is coded in the MATLAB (MathWorks, Natick, MA, USA) environment with the algorithm flowchart of Newmark as shown in Figure 2.Figure 2
Flowchart of Newmark-beta method for solving the dynamic response problem of the shell.WhenC=0 and Ft+FTt=0, from equation (24), natural frequencies can be obtained by figuring out the roots of the following equation:(23)Mq¨+Kq=0.By assuming thatq=q0sinωt≠0, we obtain(24)K−ω2Mq0sinωt=0⟶K−ω2M=0,where ω is the natural frequency.Now, we have the flowchart of the Newmark-beta method with symbols “n” is the total number of integral steps and “j” is the number of integral step.Step 1.
Determine the initial conditions(25)q0=q0;q˙0=q˙0.
From the initial conditions, we obtain(26)q¨0=M0−1F0−K0q0−C0q˙0.Step 2.
By approximatingq¨t+Δt,q˙t+Δt by qt+Δt, we have(27)q¨t+Δt=a0qt+Δt−qt−a2q˙t−a3q¨t,q˙t+Δt=q˙t+a6q¨t+a7q¨t+Δt,where(28)a0=2γΔt2;a1=2αγΔt;a2=2γΔt;a3=1γ−1;a4=2αγ−1;a5=αγ−1Δt;a6=1−αΔt;a7=αΔt,in which αandγ are defined by the assumption that the acceleration varies in each calculating step; the authors select the linear law for the varying acceleration:(29)q¨τ=q¨t+τΔtq¨t+Δt−q¨twitht≤τ≤t+Δtthenα=12;γ=13.
The condition to stabilize the roots is expressed as follows:(30)Δt≤12ϖmax1α−γorΔtTmin≤12π21α−γ.Step 3.
Compute the efficiency stiffness matrix and the efficiency nodal force vector:(31)K∗=K+a0M+a1C,F∗=Ft+Δt+Ma0qt+a2q˙t+a3q¨t+Ca1qt+a4q˙t+a5q¨t.Step 4.
Determine the nodal displacement vectorqt+Δt:(32)Kt+Δt∗qt+Δt=Ft+Δt∗⟹qt+Δt=Kt+Δt∗−1Ft+Δt∗,and then repeat the loop until the time runs out.


## 3. Numerical Examples and Discussion

Based on the finite element equations, the authors then code the calculation program in the MATLAB environment to examine the influence of some geometric parameters and mechanical properties on the forced vibration of the shell. In this study, we consider four different cases of graphene reinforcement: UD-type, X-type, O-type, and V-type with ten composite layers, which depend on the mass fraction of graphene. For the whole analysis, the volume component of graphene in the all considered cases of ten composite layers is expressed as follows [35]:(i)
Case 1: UD-type: [0.07/0.07/0.07/0.07/0.07]s(ii)
Case 2: X-type: [0.11/0.09/0.07/0.05/0.03]s(iii)
Case 3: O-type: [0.03/0.05/0.07/0.09/0.11]s(iv)
Case 4: V-type: [(0.11)2/(0.09)2/(0.07)2/(0.05)2/(0.03)2]Herein, we choose polymethyl methacrylate (PMMA) as the matrix, which has temperature-dependent material properties asEm = (3.52–0.0034 T) GPa, ρm=1150kg/m3, αm=45+0.0225T−T0×10−6K−1, and νm=0.34, in which T0 is the temperature inside the room, and it is equal to 300 K. Temperature-dependent material characteristics of monolayer graphene with aG = 14.76 nm, bG = 14.77 nm, hG = 0.188 nm, and ν12G=0.177 are shown in Table 1. PMMA (C5O2H8)n is a transparent thermoplastic, which is often used in the plate and shell structures due to its advantages as a light weight, shatter-resistant alternative to glass. The same material can be used as inks, coatings, a casting resin, and so on. Although PMMA is not considered as a type of familiar silica-based glass, the substance, similar to many thermoplastics, is often technically categorized as a type of glass (in that it is called as a no-crystalline vitreous substance); therefore, it is named as acrylic glass. PMMA is an economical alternative to polycarbonate (PC) with high tensile and flexural strength, transparency, polishability, and chemical and heat resistance; especially, it can be protected from UV tolerance. Nonmodified PMMA behaves in a brittle manner under load, typically when suffering impact load, and is more prone to scratching than other common inorganic glass, but modified PMMA sometimes could gain high ability in this property.Table 1
Temperature-dependent material characteristics of monolayer graphene.
T (K)
E
11
G (TPa)
E
22
G (TPa)
G
12
G (TPa)
α
12
G (10−6 K−1)
α
22
G (10−6 K−1)
300
1.812
1.807
0.683
−0.90
−0.95
400
1.769
1.763
0.691
−0.35
−0.40In Table2, we illustrate the efficiency parameters which are related to the volume component of graphene.Table 2
Efficiency parameters related to volume component of graphene.
T (K)
V
G
η
1
η
2
η
3
300
0.03
2.929
2.855
11.842
0.05
3.068
2.962
15.944
0.07
3.013
2.966
23.575
0.09
2.647
2.609
32.816
0.11
2.311
2.260
33.125
400
0.03
2.977
2.896
13.928
0.05
3.128
3.023
15.229
0.07
3.060
3.027
22.588
0.09
2.701
2.603
28.868
0.11
2.405
2.337
29.527


### 3.1. Verification Examples

Example 1.
Let us consider a (0°/90°/0°/90°/0°)s composite shell with geometrical parameters:a = b, thickness h = a/10, radii Rx = Ry = 50 h (for elliptic paraboloid-EPS), Rx = 50 h, Ry = ∞ (for cylindrical-CLY). Material properties of the shell are presented in Tables 1 and 2. The comparative data of the first nondimensional free vibration ω¯=ωa2/hρ0/E0 with those of Rout et al. [35] are shown in Table 3 (ρ0 and E0 are mass density and elastic modulus of shell at the room temperature T = T0 = 300 K, respectively).
In ref. [35], Rout et al. employed the finite element method (FEM) based on high-order shear deformation theory (HSDT) with the 8-node isoparametric element, where each node consists of 7 degrees of freedom (DOF). For shell structures with medium thicknesses (h/a = 1/100 ÷ 1/10) and the radius (R/a > 50), using the 8-node isoparametric element with 5 DOF per one element to analyze free vibration of the shell is very convenient. From Table 3, we can see clearly that the computed results between this work and the analytical method [35] are in good agreement; therefore, it demonstrates that our proposed theory and program are verified for free vibration problems and the convergence rate is guaranteed with the mesh 8 × 8. However, in order to improve accuracy, we used 10 × 10 mesh.Table 3
The first nondimensional fundamental frequencies with different meshes.
Cases of GPL
UD-type
X-type
O-type
V-type
BC
SSSS
CCCC
SSSS
CCCC
SSSS
CCCC
SSSS
CCCC
SPH (spherical shell)
Mesh
4 × 4
28.9466
45.7319
29.7675
44.5987
23.9914
39.9593
25.9421
41.9856
6 × 6
28.9305
45.6370
29.7532
44.5122
23.9789
39.8547
25.9291
41.8878
8 × 8
28.9279
45.6238
29.7508
44.4995
23.9769
39.8417
25.9270
41.8752
10 × 10
28.9273
45.6202
29.7502
44.4959
23.9763
39.8384
25.9291
41.8719
12 × 12
28.9270
45.6189
29.7500
44.4946
23.9761
39.8372
25.9264
41.8706
14 × 14
28.9269
45.6183
29.7499
44.4940
23.9760
39.8367
25.9262
41.8701
[35]
29.5920
44.0300
30.8370
43.2450
25.3060
39.4860
27.4230
40.5730
CLY (cylindrical shell)
Mesh
4 × 4
27.0804
43.7319
28.4903
41.3521
22.7309
38.5900
24.5928
40.6802
6 × 6
27.0659
43.6467
28.4761
41.2750
22.7179
38.4899
24.5792
40.5859
8 × 8
27.0635
43.6348
28.4737
41.2635
22.7158
38.4780
24.5770
40.5741
10 × 10
27.0629
43.6316
28.4731
41.2603
22.7152
38.4750
24.5764
40.5711
12 × 12
27.0627
43.6304
28.4729
41.2590
22.7150
38.4739
24.5762
40.5700
14 × 14
27.0626
43.6299
28.4728
41.2585
22.7149
38.4734
24.5761
40.5695
[35]
28.4430
42.4540
29.7630
41.7740
23.9830
37.9270
25.9830
39.0500Example 2.
Consider a fully clamped square plate with parameters that can be found in [41]: a = b = 1 m and h/a = 10. Material characteristics are Young’s modulus E = 30 GPa, Poisson’s ratio ν = 0.3, and the density ρ = 2800 kg/m3. The structure is subjected to distribution blast load p0 = 104 Pa. The nondimensional displacement is calculated by the formula w∗=100Eh3/12p0a41−ν2w0. The calculation program is used with the mesh 10 × 10 (let the radii of the shell be very large Rx=Ry=∞, the elastic modulus be E1 = E2 = E, and the fiber angle of the k layer be θ = 0°). The deflection of the center point of the structure and those of [41] are shown in Figure 3 (integral time is 5 ms, and acting time of load is 2 ms).
We can see from Figure3 that the results are similar to each other (both shape and value). This proves that our program is verified.Figure 3
The deflection of the center point of the structure over time.


### 3.2. Effects of Some Parameters on Forced Vibration of the Shell

Next, we study the effects of some parameters on forced vibration of the shell structure; we consider a fully clamped composite shell consisting of ten composite layers with the following geometrical parameters: the lengtha = 1 m, the width b, the thickness h, and radii of the shell Rx = Ry = R; material parameters are shown in Tables 1 and 2. All cases of graphene reinforcement are UD-type, X-type, O-type, and V-type as mentioned. The external uniform load pt acts on the upper surface of the shell in the z-direction.(33)pt=ΔPΦ⋅Ft;Ft=1−tτhd,0≤t≤τhd,0,otherwise,withΔPΦ=0.20679.106N/m2,τhd=0.028s.Nondimensional displacement and velocity of the center point over time are given as follows:(34)w∗=100h03EtΔPΦa4wa2,b2,v∗=tmaxh03EtΔPΦa4va2,b2,σx∗=1010h04ΔPΦa2σxa2,b2,0,σxy∗=102h04ΔPΦa2σxya2,b2,0,withh0=a25,tmax=0.15s,where wa/2,b/2 and va/2,b/2 are the deflection and velocity of the centroid of the shell, respectively.The nondimensional frequencies are calculated asw∗=ωa3/h0ρm/Em0 (where Em0 is value of Em at T = 300 K). For fully clamped shell, the first five nondimensional frequencies are presented in Table 4. From numerical results, it can be seen that with the same boundary condition and the temperature in 4 cases of graphene reinforcement, X-type gives the largest nondimensional frequencies, followed by UD-type, V-type, and O-type, respectively. This demonstrates that the composite shell reinforced with graphene X-type is the strongest one and the smallest one is O-type. For all distributions of graphene, the temperature decreases the nondimensional fundamental frequencies. Therefore, “stiffness” of the shell would be reduced.Table 4
The first five nondimensional fundamental frequencies of fully clamped shell.
UD-type
X-type
O-type
V-type
T = 300 K
ω
1
∗
55.0585
57.9374
46.0600
49.0584
ω
2
∗
105.4688
111.5374
87.5778
93.7707
ω
3
∗
105.7189
111.9037
87.6186
93.9600
ω
4
∗
151.0943
160.8222
125.7792
135.2560
ω
5
∗
184.6815
194.3923
154.0573
164.2965
T = 400 K
ω
1
∗
51.9456
54.5550
43.5648
46.7749
ω
2
∗
99.3956
104.9268
82.7213
89.1410
ω
3
∗
99.6000
105.2008
82.7251
89.4610
ω
4
∗
141.8134
150.0833
118.5870
128.0876
ω
5
∗
174.0027
182.9746
145.3327
156.2593First four vibration mode shapes of the shell structure are presented in Figure4. For the case of clamping all edges of the shell, the second and the third vibration mode shapes are similar to each other (the second fundamental frequency is approximated the third fundamental frequency). This phenomenon is suitable to the actual symmetrical shell structures under the same supporting conditions.Figure 4
First four vibration mode shapes of fully clamped SPH: (a) mode 1; (b) mode 2; (c) mode 3; (d) mode 4.
(a)
(b)
(c)
(d)


### 3.3. Influence of Four Different Types of Graphene-Reinforced Structures

In this section, we investigate the effect of four different types of graphene reinforcement. Let us consider a fully clamped (CCCC) composite shell with geometrical parametersa = b = 1 m, h = a/25, and R/a = 10. Nondimensional deflection, velocity, and the stress (when the deflection obtains the maximum value) of the centroid of the shell at the room temperature T = T0 = 300 K are presented in Figure 5. In the case of T = 400 K, nondimensional results (including velocity, stress and deflection) of the center point of the shell are presented in Figure 6 and the maximum value is introduced in Table 5.Figure 5
Nondimensional deflection, stress, and velocity of the center point atT = 300 K: (a) deflection w∗; (b) velocity v∗; (c) stress σx; (d) stress σxy.
(a)
(b)
(c)
(d)Figure 6
Nondimensional deflection, stress, and velocity of the center point atT = 400 K: (a) deflection w∗; (b) velocity v∗; (c) stress σx; (d) stress σxy.
(a)
(b)
(c)
(d)Table 5
Effect of four different types of graphene reinforcement on nondimensional velocity and deflection of the center point.
Temperature
T
=
300
K
T
=
400
K
Maximum value
w
∗
v
∗
w
∗
v
∗
UD-type
0.0705
0.1464
0.0873
0.1655
X-type
0.0616
0.1383
0.0768
0.1493
O-type
0.1068
0.2062
0.1268
0.2251
V-type
0.0931
0.1846
0.1097
0.2009From Figures5 and 6 and Table 5, we can see that in four cases of graphene reinforcement, X-type gives the smallest displacement, velocity, and stress responses of the centroid of shell structure, followed by the UD-type, V-type, and O-type, respectively. This demonstrates that with the same geometrical parameters, the composite shell reinforced with graphene X-type is the strongest one and the smallest one is O-type. Due to the amount of graphene nanoplatelets added in the original material, the structure will become stiffer. The distribution law of graphene nanoplatelets also affects strongly on the stiffness of structures. For example, X-type has more nanographene particles focused on the surface of the shell while O-type has more nanographene particles focused on the neutral plane of the shell; therefore, the stiffness of X-type is higher than that of O-type. Similarly, in the cases of nanographene particles distributing uniformly in layers of the structure (UD-type) and reducing smoothly from the top surface to the bottom surface, the stiffeners are smaller than those of X-type. It can be confirmed that nanographene particles are distributed more on the surface of the shell and less gradually into the core (X-type); we can obtain the structure, where the surface layer is stiffer than the core layer. This characteristic will help the structure can absorb energy better. One more thing is that in the thermal environment, all types of reinforced shell structures give the smaller displacement, velocity, and stress responses than those of the structures in the room temperature. This phenomenon shows that the temperature reduces the stiffness of the structure, in which the thermal stress appears; after the acting process, the displacement response tends to move away from position 0.


#### 3.3.1. Influence ofa/h Ratio

Next, to investigate the influence ofa/h ratio, we consider the shell structure with geometrical parameters b = a (while a is fixed), the value of a/h ratio gets 10, 25, 50 and 75, R/a = 10. The nondimensional deflection, velocity, and the stress (when the deflection reaches the maximum value) of the centroid of the shell are presented in Figures 7 and 8, and the maximum value is listed in Table 6.Figure 7
Influence ofa/h ratio on nondimensional deflection, velocity, and stress of the center point at T = 300 K of V-type shell: (a) deflection w∗; (b) velocity v∗; (c) stress σx; (d) stress σxy.
(a)
(b)
(c)
(d)Figure 8
Influence ofa/h ratio on nondimensional deflection, velocity, and stress of the center point at T = 400 K of V-type shell: (a) deflection w∗; (b) velocity v∗; (c) stress σx; (d) stress σxy.
(a)
(b)
(c)
(d)Table 6
Influence ofa/h ratio on the nondimensional deflection, velocity, and stress of the centroid point
a/h
10
25
50
75
Maximum value
w
∗
v
∗
w
∗
v
∗
w
∗
v
∗
w
∗
v
∗
TemperatureT=300K
UD-type
0.0052
0.0182
0.0705
0.1464
0.4652
0.6832
1.1324
1.4665
X-type
0.0046
0.0171
0.0616
0.1383
0.4311
0.6469
1.0903
1.4413
O-type
0.0074
0.0212
0.1068
0.2062
0.6239
0.8364
1.4910
1.6342
V-type
0.0065
0.0180
0.0931
0.1846
0.5753
0.7604
1.4027
1.6145
TemperatureT=400K
UD-type
0.0072
0.0176
0.0873
0.1655
0.5418
0.6793
1.3206
1.5670
X-type
0.0064
0.0164
0.0768
0.1493
0.5058
0.6837
1.2767
1.5529
O-type
0.0100
0.0237
0.1268
0.2251
0.7014
0.9191
1.6752
1.7314
V-type
0.0086
0.0210
0.1097
0.2009
0.6417
0.8327
0.9758
1.6705We can see in Figures7 and 8 and Table 6, for all cases of graphene-reinforced shells (UD-type, X-type, O-type, and V-type), when the a/h ratio increases, nondimensional deflection, velocity, and stress (change over time) of the center point also increase rapidly. For the cases of nanographene distribution, X-type gives the smallest displacement, velocity, and stress responses, followed by UD-type and V-type; the largest one is O-type. From Figures 8(c) and 8(d), we can see that the stress (by the thickness direction) of the center point of the shell (when the displacement obtains the maximum value) is suitable to the nanographene distribution. It is appropriate in this situation because with the same geometric dimensions, materials, boundary conditions, and thinner shells have smaller stiffness. And the temperature is also the reason of increasing displacement, velocity, and stress of the centroid point.


#### 3.3.2. Influence ofa/b Ratio

Next, we examine the influence of the length-to-width ratio on the nondimensional deflection and velocity of the center point of the shell witha = 1 m; a/b gets the values of 0.5, 1, 1.5, and 2. Geometrical parameters are h = a/25 and R/a = 10. The numerical results of the nondimensional deflection, velocity, and stress (when the deflection reaches the maximum value) of the centroid of the shell are shown in Figures 9 and 10, and the maximum values are listed in Table 7.Figure 9
Influence of length-to-width ratio on nondimensional deflection, velocity, and stress of the centroid atT = 300 K of O-type shell: (a) deflection w∗; (b) velocity v∗; (c) stress σx; (d) stress σxy.
(a)
(b)
(c)
(d)Figure 10
Influence of length-to-width ratio on nondimensional deflection, velocity, and stress of the center point atT = 400 K of O-type shell: (a) deflection w∗; (b) velocity v∗; (c) stress σx; (d) stress σxy.
(a)
(b)
(c)
(d)Table 7
Influence of length-to-width ratio on the nondimensional deflection and velocity of the center point.
a/b
0.5
1.0
2.0
4.0
Maximum value
w
∗
v
∗
w
∗
v
∗
w
∗
v
∗
w
∗
v
∗
TemperatureT=300K
UD-type
0.1258
0.2321
0.0705
0.1464
0.0086
0.0340
0.0006
0.0036
X-type
0.1186
0.2465
0.0616
0.1383
0.0076
0.0311
0.0005
0.0035
O-type
0.1892
0.2915
0.1068
0.2062
0.0129
0.0475
0.0008
0.0049
V-type
0.1671
0.2994
0.0931
0.1846
0.0112
0.0452
0.0007
0.0045
TemperatureT=400K
UD-type
0.1531
0.2732
0.0873
0.1655
0.0107
0.0402
0.0007
0.0042
X-type
0.1377
0.2283
0.0768
0.1493
0.0095
0.0346
0.0007
0.0039
O-type
0.2250
0.3244
0.1268
0.2251
0.0158
0.0415
0.0010
0.0055
V-type
0.1963
0.2949
0.1097
0.2009
0.0135
0.0473
0.0009
0.0051Now, we can see in Figures9 and 10 and Table 7 that when increasing a/b ratio, the nondimensional deflection and velocity of the center point overtime decrease for all cases of graphene-reinforced shells (UD-type, X-type, O-type, and V-type). This demonstrates that the stiffness of the shell gets larger, especially when a/b ratio equals to 4. The stress (by the thickness direction) of the center point of the shell is suitable to the nanographene distribution (see Figures 10(c) and 10(d)). This can be understood obviously that when the shape of structure gets smaller, while boundary condition and other parameters are not changed, the structure will become stronger.


## 3.1. Verification Examples

Example 1.
Let us consider a (0°/90°/0°/90°/0°)s composite shell with geometrical parameters:a = b, thickness h = a/10, radii Rx = Ry = 50 h (for elliptic paraboloid-EPS), Rx = 50 h, Ry = ∞ (for cylindrical-CLY). Material properties of the shell are presented in Tables 1 and 2. The comparative data of the first nondimensional free vibration ω¯=ωa2/hρ0/E0 with those of Rout et al. [35] are shown in Table 3 (ρ0 and E0 are mass density and elastic modulus of shell at the room temperature T = T0 = 300 K, respectively).
In ref. [35], Rout et al. employed the finite element method (FEM) based on high-order shear deformation theory (HSDT) with the 8-node isoparametric element, where each node consists of 7 degrees of freedom (DOF). For shell structures with medium thicknesses (h/a = 1/100 ÷ 1/10) and the radius (R/a > 50), using the 8-node isoparametric element with 5 DOF per one element to analyze free vibration of the shell is very convenient. From Table 3, we can see clearly that the computed results between this work and the analytical method [35] are in good agreement; therefore, it demonstrates that our proposed theory and program are verified for free vibration problems and the convergence rate is guaranteed with the mesh 8 × 8. However, in order to improve accuracy, we used 10 × 10 mesh.Table 3
The first nondimensional fundamental frequencies with different meshes.
Cases of GPL
UD-type
X-type
O-type
V-type
BC
SSSS
CCCC
SSSS
CCCC
SSSS
CCCC
SSSS
CCCC
SPH (spherical shell)
Mesh
4 × 4
28.9466
45.7319
29.7675
44.5987
23.9914
39.9593
25.9421
41.9856
6 × 6
28.9305
45.6370
29.7532
44.5122
23.9789
39.8547
25.9291
41.8878
8 × 8
28.9279
45.6238
29.7508
44.4995
23.9769
39.8417
25.9270
41.8752
10 × 10
28.9273
45.6202
29.7502
44.4959
23.9763
39.8384
25.9291
41.8719
12 × 12
28.9270
45.6189
29.7500
44.4946
23.9761
39.8372
25.9264
41.8706
14 × 14
28.9269
45.6183
29.7499
44.4940
23.9760
39.8367
25.9262
41.8701
[35]
29.5920
44.0300
30.8370
43.2450
25.3060
39.4860
27.4230
40.5730
CLY (cylindrical shell)
Mesh
4 × 4
27.0804
43.7319
28.4903
41.3521
22.7309
38.5900
24.5928
40.6802
6 × 6
27.0659
43.6467
28.4761
41.2750
22.7179
38.4899
24.5792
40.5859
8 × 8
27.0635
43.6348
28.4737
41.2635
22.7158
38.4780
24.5770
40.5741
10 × 10
27.0629
43.6316
28.4731
41.2603
22.7152
38.4750
24.5764
40.5711
12 × 12
27.0627
43.6304
28.4729
41.2590
22.7150
38.4739
24.5762
40.5700
14 × 14
27.0626
43.6299
28.4728
41.2585
22.7149
38.4734
24.5761
40.5695
[35]
28.4430
42.4540
29.7630
41.7740
23.9830
37.9270
25.9830
39.0500Example 2.
Consider a fully clamped square plate with parameters that can be found in [41]: a = b = 1 m and h/a = 10. Material characteristics are Young’s modulus E = 30 GPa, Poisson’s ratio ν = 0.3, and the density ρ = 2800 kg/m3. The structure is subjected to distribution blast load p0 = 104 Pa. The nondimensional displacement is calculated by the formula w∗=100Eh3/12p0a41−ν2w0. The calculation program is used with the mesh 10 × 10 (let the radii of the shell be very large Rx=Ry=∞, the elastic modulus be E1 = E2 = E, and the fiber angle of the k layer be θ = 0°). The deflection of the center point of the structure and those of [41] are shown in Figure 3 (integral time is 5 ms, and acting time of load is 2 ms).
We can see from Figure3 that the results are similar to each other (both shape and value). This proves that our program is verified.Figure 3
The deflection of the center point of the structure over time.


## 3.2. Effects of Some Parameters on Forced Vibration of the Shell

Next, we study the effects of some parameters on forced vibration of the shell structure; we consider a fully clamped composite shell consisting of ten composite layers with the following geometrical parameters: the lengtha = 1 m, the width b, the thickness h, and radii of the shell Rx = Ry = R; material parameters are shown in Tables 1 and 2. All cases of graphene reinforcement are UD-type, X-type, O-type, and V-type as mentioned. The external uniform load pt acts on the upper surface of the shell in the z-direction.(33)pt=ΔPΦ⋅Ft;Ft=1−tτhd,0≤t≤τhd,0,otherwise,withΔPΦ=0.20679.106N/m2,τhd=0.028s.Nondimensional displacement and velocity of the center point over time are given as follows:(34)w∗=100h03EtΔPΦa4wa2,b2,v∗=tmaxh03EtΔPΦa4va2,b2,σx∗=1010h04ΔPΦa2σxa2,b2,0,σxy∗=102h04ΔPΦa2σxya2,b2,0,withh0=a25,tmax=0.15s,where wa/2,b/2 and va/2,b/2 are the deflection and velocity of the centroid of the shell, respectively.The nondimensional frequencies are calculated asw∗=ωa3/h0ρm/Em0 (where Em0 is value of Em at T = 300 K). For fully clamped shell, the first five nondimensional frequencies are presented in Table 4. From numerical results, it can be seen that with the same boundary condition and the temperature in 4 cases of graphene reinforcement, X-type gives the largest nondimensional frequencies, followed by UD-type, V-type, and O-type, respectively. This demonstrates that the composite shell reinforced with graphene X-type is the strongest one and the smallest one is O-type. For all distributions of graphene, the temperature decreases the nondimensional fundamental frequencies. Therefore, “stiffness” of the shell would be reduced.Table 4
The first five nondimensional fundamental frequencies of fully clamped shell.
UD-type
X-type
O-type
V-type
T = 300 K
ω
1
∗
55.0585
57.9374
46.0600
49.0584
ω
2
∗
105.4688
111.5374
87.5778
93.7707
ω
3
∗
105.7189
111.9037
87.6186
93.9600
ω
4
∗
151.0943
160.8222
125.7792
135.2560
ω
5
∗
184.6815
194.3923
154.0573
164.2965
T = 400 K
ω
1
∗
51.9456
54.5550
43.5648
46.7749
ω
2
∗
99.3956
104.9268
82.7213
89.1410
ω
3
∗
99.6000
105.2008
82.7251
89.4610
ω
4
∗
141.8134
150.0833
118.5870
128.0876
ω
5
∗
174.0027
182.9746
145.3327
156.2593First four vibration mode shapes of the shell structure are presented in Figure4. For the case of clamping all edges of the shell, the second and the third vibration mode shapes are similar to each other (the second fundamental frequency is approximated the third fundamental frequency). This phenomenon is suitable to the actual symmetrical shell structures under the same supporting conditions.Figure 4
First four vibration mode shapes of fully clamped SPH: (a) mode 1; (b) mode 2; (c) mode 3; (d) mode 4.
(a)
(b)
(c)
(d)


## 3.3. Influence of Four Different Types of Graphene-Reinforced Structures

In this section, we investigate the effect of four different types of graphene reinforcement. Let us consider a fully clamped (CCCC) composite shell with geometrical parametersa = b = 1 m, h = a/25, and R/a = 10. Nondimensional deflection, velocity, and the stress (when the deflection obtains the maximum value) of the centroid of the shell at the room temperature T = T0 = 300 K are presented in Figure 5. In the case of T = 400 K, nondimensional results (including velocity, stress and deflection) of the center point of the shell are presented in Figure 6 and the maximum value is introduced in Table 5.Figure 5
Nondimensional deflection, stress, and velocity of the center point atT = 300 K: (a) deflection w∗; (b) velocity v∗; (c) stress σx; (d) stress σxy.
(a)
(b)
(c)
(d)Figure 6
Nondimensional deflection, stress, and velocity of the center point atT = 400 K: (a) deflection w∗; (b) velocity v∗; (c) stress σx; (d) stress σxy.
(a)
(b)
(c)
(d)Table 5
Effect of four different types of graphene reinforcement on nondimensional velocity and deflection of the center point.
Temperature
T
=
300
K
T
=
400
K
Maximum value
w
∗
v
∗
w
∗
v
∗
UD-type
0.0705
0.1464
0.0873
0.1655
X-type
0.0616
0.1383
0.0768
0.1493
O-type
0.1068
0.2062
0.1268
0.2251
V-type
0.0931
0.1846
0.1097
0.2009From Figures5 and 6 and Table 5, we can see that in four cases of graphene reinforcement, X-type gives the smallest displacement, velocity, and stress responses of the centroid of shell structure, followed by the UD-type, V-type, and O-type, respectively. This demonstrates that with the same geometrical parameters, the composite shell reinforced with graphene X-type is the strongest one and the smallest one is O-type. Due to the amount of graphene nanoplatelets added in the original material, the structure will become stiffer. The distribution law of graphene nanoplatelets also affects strongly on the stiffness of structures. For example, X-type has more nanographene particles focused on the surface of the shell while O-type has more nanographene particles focused on the neutral plane of the shell; therefore, the stiffness of X-type is higher than that of O-type. Similarly, in the cases of nanographene particles distributing uniformly in layers of the structure (UD-type) and reducing smoothly from the top surface to the bottom surface, the stiffeners are smaller than those of X-type. It can be confirmed that nanographene particles are distributed more on the surface of the shell and less gradually into the core (X-type); we can obtain the structure, where the surface layer is stiffer than the core layer. This characteristic will help the structure can absorb energy better. One more thing is that in the thermal environment, all types of reinforced shell structures give the smaller displacement, velocity, and stress responses than those of the structures in the room temperature. This phenomenon shows that the temperature reduces the stiffness of the structure, in which the thermal stress appears; after the acting process, the displacement response tends to move away from position 0.


### 3.3.1. Influence ofa/h Ratio

Next, to investigate the influence ofa/h ratio, we consider the shell structure with geometrical parameters b = a (while a is fixed), the value of a/h ratio gets 10, 25, 50 and 75, R/a = 10. The nondimensional deflection, velocity, and the stress (when the deflection reaches the maximum value) of the centroid of the shell are presented in Figures 7 and 8, and the maximum value is listed in Table 6.Figure 7
Influence ofa/h ratio on nondimensional deflection, velocity, and stress of the center point at T = 300 K of V-type shell: (a) deflection w∗; (b) velocity v∗; (c) stress σx; (d) stress σxy.
(a)
(b)
(c)
(d)Figure 8
Influence ofa/h ratio on nondimensional deflection, velocity, and stress of the center point at T = 400 K of V-type shell: (a) deflection w∗; (b) velocity v∗; (c) stress σx; (d) stress σxy.
(a)
(b)
(c)
(d)Table 6
Influence ofa/h ratio on the nondimensional deflection, velocity, and stress of the centroid point
a/h
10
25
50
75
Maximum value
w
∗
v
∗
w
∗
v
∗
w
∗
v
∗
w
∗
v
∗
TemperatureT=300K
UD-type
0.0052
0.0182
0.0705
0.1464
0.4652
0.6832
1.1324
1.4665
X-type
0.0046
0.0171
0.0616
0.1383
0.4311
0.6469
1.0903
1.4413
O-type
0.0074
0.0212
0.1068
0.2062
0.6239
0.8364
1.4910
1.6342
V-type
0.0065
0.0180
0.0931
0.1846
0.5753
0.7604
1.4027
1.6145
TemperatureT=400K
UD-type
0.0072
0.0176
0.0873
0.1655
0.5418
0.6793
1.3206
1.5670
X-type
0.0064
0.0164
0.0768
0.1493
0.5058
0.6837
1.2767
1.5529
O-type
0.0100
0.0237
0.1268
0.2251
0.7014
0.9191
1.6752
1.7314
V-type
0.0086
0.0210
0.1097
0.2009
0.6417
0.8327
0.9758
1.6705We can see in Figures7 and 8 and Table 6, for all cases of graphene-reinforced shells (UD-type, X-type, O-type, and V-type), when the a/h ratio increases, nondimensional deflection, velocity, and stress (change over time) of the center point also increase rapidly. For the cases of nanographene distribution, X-type gives the smallest displacement, velocity, and stress responses, followed by UD-type and V-type; the largest one is O-type. From Figures 8(c) and 8(d), we can see that the stress (by the thickness direction) of the center point of the shell (when the displacement obtains the maximum value) is suitable to the nanographene distribution. It is appropriate in this situation because with the same geometric dimensions, materials, boundary conditions, and thinner shells have smaller stiffness. And the temperature is also the reason of increasing displacement, velocity, and stress of the centroid point.


### 3.3.2. Influence ofa/b Ratio

Next, we examine the influence of the length-to-width ratio on the nondimensional deflection and velocity of the center point of the shell witha = 1 m; a/b gets the values of 0.5, 1, 1.5, and 2. Geometrical parameters are h = a/25 and R/a = 10. The numerical results of the nondimensional deflection, velocity, and stress (when the deflection reaches the maximum value) of the centroid of the shell are shown in Figures 9 and 10, and the maximum values are listed in Table 7.Figure 9
Influence of length-to-width ratio on nondimensional deflection, velocity, and stress of the centroid atT = 300 K of O-type shell: (a) deflection w∗; (b) velocity v∗; (c) stress σx; (d) stress σxy.
(a)
(b)
(c)
(d)Figure 10
Influence of length-to-width ratio on nondimensional deflection, velocity, and stress of the center point atT = 400 K of O-type shell: (a) deflection w∗; (b) velocity v∗; (c) stress σx; (d) stress σxy.
(a)
(b)
(c)
(d)Table 7
Influence of length-to-width ratio on the nondimensional deflection and velocity of the center point.
a/b
0.5
1.0
2.0
4.0
Maximum value
w
∗
v
∗
w
∗
v
∗
w
∗
v
∗
w
∗
v
∗
TemperatureT=300K
UD-type
0.1258
0.2321
0.0705
0.1464
0.0086
0.0340
0.0006
0.0036
X-type
0.1186
0.2465
0.0616
0.1383
0.0076
0.0311
0.0005
0.0035
O-type
0.1892
0.2915
0.1068
0.2062
0.0129
0.0475
0.0008
0.0049
V-type
0.1671
0.2994
0.0931
0.1846
0.0112
0.0452
0.0007
0.0045
TemperatureT=400K
UD-type
0.1531
0.2732
0.0873
0.1655
0.0107
0.0402
0.0007
0.0042
X-type
0.1377
0.2283
0.0768
0.1493
0.0095
0.0346
0.0007
0.0039
O-type
0.2250
0.3244
0.1268
0.2251
0.0158
0.0415
0.0010
0.0055
V-type
0.1963
0.2949
0.1097
0.2009
0.0135
0.0473
0.0009
0.0051Now, we can see in Figures9 and 10 and Table 7 that when increasing a/b ratio, the nondimensional deflection and velocity of the center point overtime decrease for all cases of graphene-reinforced shells (UD-type, X-type, O-type, and V-type). This demonstrates that the stiffness of the shell gets larger, especially when a/b ratio equals to 4. The stress (by the thickness direction) of the center point of the shell is suitable to the nanographene distribution (see Figures 10(c) and 10(d)). This can be understood obviously that when the shape of structure gets smaller, while boundary condition and other parameters are not changed, the structure will become stronger.


## 3.3.1. Influence ofa/h Ratio

Next, to investigate the influence ofa/h ratio, we consider the shell structure with geometrical parameters b = a (while a is fixed), the value of a/h ratio gets 10, 25, 50 and 75, R/a = 10. The nondimensional deflection, velocity, and the stress (when the deflection reaches the maximum value) of the centroid of the shell are presented in Figures 7 and 8, and the maximum value is listed in Table 6.Figure 7
Influence ofa/h ratio on nondimensional deflection, velocity, and stress of the center point at T = 300 K of V-type shell: (a) deflection w∗; (b) velocity v∗; (c) stress σx; (d) stress σxy.
(a)
(b)
(c)
(d)Figure 8
Influence ofa/h ratio on nondimensional deflection, velocity, and stress of the center point at T = 400 K of V-type shell: (a) deflection w∗; (b) velocity v∗; (c) stress σx; (d) stress σxy.
(a)
(b)
(c)
(d)Table 6
Influence ofa/h ratio on the nondimensional deflection, velocity, and stress of the centroid point
a/h
10
25
50
75
Maximum value
w
∗
v
∗
w
∗
v
∗
w
∗
v
∗
w
∗
v
∗
TemperatureT=300K
UD-type
0.0052
0.0182
0.0705
0.1464
0.4652
0.6832
1.1324
1.4665
X-type
0.0046
0.0171
0.0616
0.1383
0.4311
0.6469
1.0903
1.4413
O-type
0.0074
0.0212
0.1068
0.2062
0.6239
0.8364
1.4910
1.6342
V-type
0.0065
0.0180
0.0931
0.1846
0.5753
0.7604
1.4027
1.6145
TemperatureT=400K
UD-type
0.0072
0.0176
0.0873
0.1655
0.5418
0.6793
1.3206
1.5670
X-type
0.0064
0.0164
0.0768
0.1493
0.5058
0.6837
1.2767
1.5529
O-type
0.0100
0.0237
0.1268
0.2251
0.7014
0.9191
1.6752
1.7314
V-type
0.0086
0.0210
0.1097
0.2009
0.6417
0.8327
0.9758
1.6705We can see in Figures7 and 8 and Table 6, for all cases of graphene-reinforced shells (UD-type, X-type, O-type, and V-type), when the a/h ratio increases, nondimensional deflection, velocity, and stress (change over time) of the center point also increase rapidly. For the cases of nanographene distribution, X-type gives the smallest displacement, velocity, and stress responses, followed by UD-type and V-type; the largest one is O-type. From Figures 8(c) and 8(d), we can see that the stress (by the thickness direction) of the center point of the shell (when the displacement obtains the maximum value) is suitable to the nanographene distribution. It is appropriate in this situation because with the same geometric dimensions, materials, boundary conditions, and thinner shells have smaller stiffness. And the temperature is also the reason of increasing displacement, velocity, and stress of the centroid point.


## 3.3.2. Influence ofa/b Ratio

Next, we examine the influence of the length-to-width ratio on the nondimensional deflection and velocity of the center point of the shell witha = 1 m; a/b gets the values of 0.5, 1, 1.5, and 2. Geometrical parameters are h = a/25 and R/a = 10. The numerical results of the nondimensional deflection, velocity, and stress (when the deflection reaches the maximum value) of the centroid of the shell are shown in Figures 9 and 10, and the maximum values are listed in Table 7.Figure 9
Influence of length-to-width ratio on nondimensional deflection, velocity, and stress of the centroid atT = 300 K of O-type shell: (a) deflection w∗; (b) velocity v∗; (c) stress σx; (d) stress σxy.
(a)
(b)
(c)
(d)Figure 10
Influence of length-to-width ratio on nondimensional deflection, velocity, and stress of the center point atT = 400 K of O-type shell: (a) deflection w∗; (b) velocity v∗; (c) stress σx; (d) stress σxy.
(a)
(b)
(c)
(d)Table 7
Influence of length-to-width ratio on the nondimensional deflection and velocity of the center point.
a/b
0.5
1.0
2.0
4.0
Maximum value
w
∗
v
∗
w
∗
v
∗
w
∗
v
∗
w
∗
v
∗
TemperatureT=300K
UD-type
0.1258
0.2321
0.0705
0.1464
0.0086
0.0340
0.0006
0.0036
X-type
0.1186
0.2465
0.0616
0.1383
0.0076
0.0311
0.0005
0.0035
O-type
0.1892
0.2915
0.1068
0.2062
0.0129
0.0475
0.0008
0.0049
V-type
0.1671
0.2994
0.0931
0.1846
0.0112
0.0452
0.0007
0.0045
TemperatureT=400K
UD-type
0.1531
0.2732
0.0873
0.1655
0.0107
0.0402
0.0007
0.0042
X-type
0.1377
0.2283
0.0768
0.1493
0.0095
0.0346
0.0007
0.0039
O-type
0.2250
0.3244
0.1268
0.2251
0.0158
0.0415
0.0010
0.0055
V-type
0.1963
0.2949
0.1097
0.2009
0.0135
0.0473
0.0009
0.0051Now, we can see in Figures9 and 10 and Table 7 that when increasing a/b ratio, the nondimensional deflection and velocity of the center point overtime decrease for all cases of graphene-reinforced shells (UD-type, X-type, O-type, and V-type). This demonstrates that the stiffness of the shell gets larger, especially when a/b ratio equals to 4. The stress (by the thickness direction) of the center point of the shell is suitable to the nanographene distribution (see Figures 10(c) and 10(d)). This can be understood obviously that when the shape of structure gets smaller, while boundary condition and other parameters are not changed, the structure will become stronger.


## 4. Conclusions

In this paper, the numerical results of mechanical responses of laminated composite shells reinforced with graphene nanoplatelets are explored. The authors use FEM based on FSDT that has the following advantages:(i)
Simple formulations for a theoretical representation(ii)
The computed results of free vibration and forced vibration obtained by this approach compared to other solutions show a good agreementFrom new computed results, several conclusions may be drawn as follows:(i)
With the same geometric parameters, the X-type has the largest hardness and the smallest one is the O-type structure. In the temperature environment, the shells with all reinforcement types are greater in the displacement response, velocity, and stress than those of the shells in the room temperature. This demonstrates that heat reduces the stiffness of the shell, and thermal stresses appear in the structure. Thus, due to the effect of the external load, the displacement response tends to move away from position 0.(ii)
When increasinga/b ratio, the stiffness of the shell gets larger, especially when a/b ratio equals 4. This can be understood obviously that when the shape of structure gets smaller, while boundary conditions and other parameters are not changed, the structure will become stronger.(iii)
We also understand that in all cases, the maximum displacement of the centroid of the structures performing in the temperature environmentT = 400 K is approximately 10–20% larger than that in case of T = 300 K. Hence, we would like to express that the temperature has a small influence on the deflections of the structure under the impact loads in the thermal environment. These remarkable results are a reliable suggestion for engineers when using these laminated composite shells reinforced with graphene nanoplatelets in the thermal environment.Finally, based on achieved numerical results, the proposed program is able to analyze the static bending, dynamic response, nonlinear problems, etc., with complicated structures, which is not convenient to solve by analytical and other methods.



---

*Source: 1471037-2020-01-03.xml*