Uncertainty propagation in dynamics of composite plates: A semi-analytical non-sampling-based approach

Abstract

In this study, the influences of spatially varying stochastic properties on free vibration analysis of composite plates were investigated via development of a new approach named the deterministic-stochastic Galerkin-based semi-analytical method. The material properties including tensile modulus, shear modulus, and density of the plate were assumed to be spatially varying and uncertain. Gaussian fields with first-order Markov kernels were utilized to define the aforementioned material properties. The stochastic fields were decomposed via application of the Karhunen-Loeve theorem. A first-order shear deformation theory was assumed, following which the displacement field was defined using admissible trigonometric modes to derive the potential and kinetic energies. The stochastic equations of motion of the plate were obtained using the variational principle. The deterministic-stochastic Galerkin-based method was utilized to find the probability space of natural frequencies, and the corresponding mode shapes of the plate were determined using a polynomial chaos approach. The proposed method significantly reduced the size of the mathematical models of the structure, which is very useful for enhancing the computational efficiency of stochastic simulations. The methodology was verified using a stochastic finite element method and the available results in literature. The sensitivity of natural frequencies and corresponding mode shapes due to the uncertainty of material properties was investigated, and the results indicated that the higher-order modes are more sensitive to uncertainty propagation in spatially varying properties.

Keywords

composite plate / spatially varying stochastic properties / Galerkin method / polynomial chaos approach / semi-analytical approach

Introduction

Composite structures are widely used in different industries. Uncertainties in the mechanical properties of these structures are common because the fabrication processes of the materials used have many parameters that are not controllable. In addition, such uncertainties can also be caused by any type of mechanical or metallurgical damage. Uncertainty quantification of structures in free vibration is essential for structural reliability assessment. It has a very important application in vibration-based structural health monitoring, wherein natural frequency variation is required. In addition, vibrational modes are needed to find the location of damages. It is worth mentioning that modal space is also necessary in post-processing simulations such as in seismic and harmonic analyses.

Stochastic simulations and sensitivity analysis are performed for different purposes. Hamdia and Rabczuk [¹] investigated the sensitivity of fracture toughness of polymer nanocomposites with respect to variations in system parameters by employing different sensitivity analysis techniques. Vu-Bac et al. [²] studied the stochastic bulk properties of amorphous polyethylene based on molecular dynamics simulations. The stochastic interfacial characteristics of polymeric nanocomposites were evaluated by using polynomial regression, moving least squares and hybrid of quadratic polynomial, and moving least squares regression methods [³]. A unified framework for the stochastic prediction of mechanical properties of polymeric nanocomposites was proposed [⁴]. Vu-Bac et al. [⁵] carried out uncertainty quantification for the multiscale modeling of polymer nanocomposites with correlated parameters. The mechanical properties of polyethylene based on a fully atomistic model were investigated by using stochastic frameworks [⁶]. A useful software framework for probabilistic sensitivity analysis of computationally expensive models was presented in Ref. [⁷].

Uncertainties in the properties of composite structures are considered in the literature as random variables or random fields. In practice, experimental data show spatially varying Gaussian or non-Gaussian properties [⁸,⁹]; thus, it is very important to simulate the stochasticity of properties with random fields. The number of studies that considered a spatially varying stochasticity of composite structure properties in their simulations is limited in comparison to studies based on a random variable assumption. Both static [¹⁰–¹³] and dynamic responses of composites [¹⁴,¹⁵] are addressed considering the spatial stochasticity of properties.

Another aspect that distinguishes the studies investigating stochasticity in composite structures is the method used for the stochastic simulation. The basic method is the Monte Carlo approach, which is widely applied in stochastic responses of composite structures with uncertain material properties because of its flexibility. For example, the static response of composite plates with uncertain material properties was studied by Navaneetha Raj et al. [¹⁶] using the Monte Carlo approach. The nonlinear free vibration of composite plates was examined in Ref. [¹⁷] by assuming the mechanical properties as random variables. The dynamics of thin-walled composite beams with stochastic material properties were studied in Ref. [¹⁸]. Sasikumar et al. [¹²] investigated the failure of composite beams with spatially varying uncertain material properties. The uncertainty of composite ply thickness and its effect on the reliability of the structure was analyzed by Zhang et al. [¹⁹]. The multivariate Monte Carlo approach was utilized for the uncertainty quantification of composite structures with defects [²⁰]. Despite its popularity, the Monte Carlo method is based on stochastic sampling, which leads to high computational cost.

Perturbation-based approaches as well as the first-order and second-order reliability methods have been introduced to enhance the computational efficiency of sampling-based approaches. For research using this category of approaches, we can refer to Lal and Singh [²¹], who focused on the stochastic free vibration of composite plates in a thermal environment. Other studies based on these methods dealt with free and forced vibrations of uncertain composite plates [²²,²³] or the static response and stability of uncertain structures employing perturbation-based approaches [²⁴–²⁶]. The deflection statistics of uncertain sandwich plates were evaluated by Grover et al. [²⁷]. A stochastic fracture analysis of composite plates was performed using a perturbation-based approach [²⁸]. A multiscale uncertainty propagation analysis in composite laminated plates has been carried out by utilizing the aforementioned approach [²⁹]. The accuracy of these methods decreases when the uncertainty level increases.

Spectral approaches based on polynomial chaos methods have been developed to enhance the computational efficiency and applicability of the aforementioned methods. For studies using these approaches, we can refer to the bending analysis of composite plates with uncertain properties using the polynomial chaos approach [³⁰]. A possibilistic approach based on polynomial chaos expansion was utilized to quantify uncertainties in the dynamics of composite plates [³¹]. Sepahvand et al. [³²] investigated the free vibration of orthotropic plates using a generalized polynomial chaos expansion. Umesh and Ganguli [³³] utilized the aforementioned approach to control the vibration of uncertain composite plates. In addition, this method was applied in acoustoelastic problems of uncertain composite plates [³⁴]. Another application of this method is uncertainty propagation analysis in aeroelastic stability of uncertain composite wings [³⁵]. Recently, a novel polynomial chaos-based approach was proposed to quantify uncertainties in composite structures with random and interval uncertainties [³⁶]. This approach was used to predict the fracture toughness of nanocomposites [³⁷].

Recent developments in stochastic simulations focused on deriving surrogate models such as response surface methods or metamodel-based approaches, which dramatically enhance the computational efficiency of stochastic simulations. For this category, we can refer to the generalized high-dimensional model representation (GHDMR) approach. The thermal uncertainty propagation of composite plates in free vibration was investigated by utilizing this method [³⁸]. GHDMR was used for analyzing uncertainties in ply orientation and material properties and their effects on natural frequencies of composite plates [³⁹,⁴⁰]. Another research evaluated the delamination effects on natural frequencies of uncertain composite plates by applying GHDMR [⁴¹]. The artificial neural network (ANN) is another metamodel-based approach that was utilized for uncertainty quantification of composite materials [⁴²] and dynamics of composite structures [⁴³,⁴⁴]. Another efficient surrogate model is the Kriging model. This method was applied for uncertainty quantification of the natural frequency of composite structures [⁴⁵,⁴⁶] and the mechanical properties of composite materials [⁴⁷]. The efficient radial basis function method is another metamodel-based approach that was utilized in stochastic dynamic analysis of composite beams [¹⁵]. A useful comparative review of recent metamodel-based approaches is presented in Ref. [⁴⁸].

It can be found from the literature that the most recent studies have focused on developing metamodels or surrogate models to decrease the computational cost of stochastic simulations based on the finite element method (FEM). However, it is worth mentioning that when there are many random variables or high uncertainty levels in the structure, the metamodel-based approaches have some limitations in presenting the total stochastic modal space [⁴⁹] (natural frequencies and corresponding mode shapes).

In high-frequency vibration evaluations, such as seismic analysis and vibration-based health monitoring, high-dimensional mathematical models (mass, damping, and stiffness matrices) lead to very high computational cost in stochastic simulations when the modal space is required.

In the present research, a stochastic Galerkin method is used to find the stochastic modal space by utilizing the polynomial chaos approach for plates with spatially varying Gaussian properties. The other issues that will be addressed in the present research are the spatial stochasticity of properties and their effects on natural frequencies and corresponding mode shapes. A semi-analytical approach is proposed by using admissible trigonometric modes for decreasing the size of the mathematical models of the structure. The Kirchhoff-Love theory [⁵⁰,⁵¹] and the first-order shear deformation plate theory (FSDPT) can be applied to define the displacement fields; however, the FSDPT is utilized in the present research because different side-to-thickness ratios are considered. For thin plates, the Kirchhoff-Love theory has the highest efficiency. It should be noted here that in the case of plates with material properties as random variables, the admissible modes method is applied in the literature. The present research shows that in the case of plates with spatially varying random properties, the admissible modes can be used as the span set for the stochastic modal space with good accuracy and a relatively smaller mathematical model, in comparison to the FEM-based approach. Finally, the combination of the aforementioned method and the polynomial chaos approach leads to a new semi-analytical non-sampling-based approach to clarify the stochastic dynamics of composite plates with spatially varying random properties. Owing to the orthogonality of admissible modes and the smaller mathematical model size (in comparison to that of the stochastic FEM (SFEM)), the computational cost of the stochastic simulations is significantly reduced. In addition, the stochastic modal space is achieved for possible post-processing. This is extremely important particularly in the case of high-frequency vibration problems where a very small element size is required in the SFEM, which leads to a very high-dimensional mathematical model. This methodology is not limited to analyzing vibration of composite plates. It can be used for dynamic analysis of similar plate and beam structures with spatially varying stochastic properties.

Problem statement

The free vibration problem will be investigated for composite plates with different side to thickness ratios (a/h) (see Fig. 1(a)). The plate is assumed to be square and simply supported in all edges, as shown Fig. 1(b).

The mechanical properties and density of laminas are modeled as stochastic fields [⁹] in the main coordinates of the plate (x, y, z), as displayed in Fig. 1(a). It should be noted here that spatial stochasticity could be caused by uncertainties in the fiber volume ratio or mechanical damages. The first-order Markov process with exponential kernels is utilized to define the stochastic fields. Exponential kernels are suitable to define spatially varying stochastic properties based on experimental studies in Ref. [⁹].

The off-axis mechanical properties of each lamina will be obtained by applying a transformation matrix.

Theoretical background

In this part, the definition of stochastic fields, problem formulation, and solution methodology are presented.

Stochastic properties

To define stochastic properties as stochastic fields, an exponential kernel as a first-order Markov kernel, as indicated in Eq. (1), is utilized:

(1)C_{\mathrm{RP}}={\sigma }_{\mathrm{RP}}^{\mathrm{2}}\exp (-\left|(x-x^{\prime })/l_x\right|-\left|(y-y^{\prime })/l_y\right|),

where {\sigma }_{\mathrm{RP}}^{\mathrm{2}}is the variance of the random properties, i.e., E_{\mathrm{11}}, E_{\mathrm{22}}=\left({\upsilon }_{\mathrm{21}}/{\upsilon }_{\mathrm{12}}\right)E_{\mathrm{11}}as Young’s moduli along the fiber direction and normal to it, respectively, G_{\mathrm{12}}=G_{\mathrm{13}}, G_{\mathrm{23}}as the in-plane and out-of-plane shear moduli, respectively, and \rhoas the density of the plate. Hence, there are four independent stochastic variables. l_xand l_yare the correlation lengths in the x and y directions. Under the assumption of deterministic and constant minor and major Poisson’s ratios, four independent stochastic fields are defined to simulate the aforementioned independent stochastic properties. By means of defining the stochastic fields in the area of x\in [-a,a],y\in [-b,b] using the aforementioned kernels and employing the Karhunen-Loeve (KL) theorem, the stochastic fields will be decomposed in the form of Eq. (2) [⁵²].

P(x,y,\mathrm{\Omega })

=P_{\mathrm{m}}\{\mathrm{1}+COV(P\left){\displaystyle \sum _{i=\mathrm{1}}^n\sqrt{{\lambda }_{xi}{\lambda }_{yi}}}{\zeta }_i^P\right(\mathrm{\Omega }\left)\mathrm{\hspace{0.17em}}{\phi }_i\right(x\left)\mathrm{\hspace{0.17em}}{\phi }_i\right(y))\},

(2)\{\begin{array}{l}E\left({\zeta }_i^{}\right)=\mathrm{0},\hfill \\ \mathrm{var}\left({\zeta }_i^{}\right)=\mathrm{1},\hfill \end{array}

where P(x,y,\omega ) is the spatially varying stochastic property (E₁₁, G₁₂, G_23, and \rho), \Omega is a probabilistic vector, P_{\mathrm{m}} is the mean value of the stochastic property, COV\left(P\right) is the coefficient of variation (equal to the standard deviation per mean value) of the stochastic property, {\zeta }_i^{\mathrm{P}}are standard independent normal random variables of the corresponding random property (P), and \lambda and {\phi }_i are eigenvalues and eigenfunctions, respectively, in the Fredholm eigenvalue problem, as defined by Eqs. (3)–(5) [⁵²].

\{\begin{array}{ll}c-{\omega }_i\tan \left({\omega }_{i}d\right)=\mathrm{0},\hfill & i:\mathrm{o}\mathrm{d}\mathrm{d},\hfill \\ {\omega }_i+c\tan \left({\omega }_{i}d\right)=\mathrm{0},\hfill & i:\mathrm{e}\mathrm{v}\mathrm{e}\mathrm{n},\hfill \end{array}

(3)\begin{array}{l}c=\mathrm{1}/l_x\mathrm{\hspace{0.17em}},\mathrm{1}/l_y,\hfill \\ d=a,b,\hfill \end{array}

(4){\lambda }_z={\displaystyle \frac{\mathrm{2}c}{{\omega }_z^{\mathrm{2}}+c^{\mathrm{2}}}},{\lambda }_z={\lambda }_{x_i},{\lambda }_{y_i},

(5){\phi }_i\left(v\right)=\{\begin{array}{ll}{\displaystyle \frac{\cos \left({\omega }_{i}v\right)}{\sqrt{a+\sin \left(\mathrm{2}{\omega }_{i}d\right)/\mathrm{2}{\omega }_i}}},\hfill & i:\mathrm{odd,}\hfill \\ {\displaystyle \frac{\sin \left({\omega }_{i}v\right)}{\sqrt{a-\sin \left(\mathrm{2}{\omega }_{i}d\right)/\mathrm{2}{\omega }_i}}},\hfill & i:\mathrm{even,}\hfill \end{array}

v=x,y,d=a,b.

The number of terms n in Eq. (2) can be extracted by a criterion (Eq. (6)):

(6){\displaystyle \frac{{\displaystyle \sum _{i=\mathrm{1}}^n{\lambda }_i}}{{\displaystyle \sum _{i=\mathrm{1}}^{\infty }{\lambda }_i}}}\ge \tau .

In this equation, \tauis usually equal to 90% [⁵³].

Strain definition

The stochastic displacement field can be defined using the first-order shear deformation theory, as shown in Eq. (7):

(7)\begin{array}{l}u(x,y,\mathrm{\Omega })=u_{\mathrm{0}}(x,y,\mathrm{\Omega })+z{\varphi }_x(x,y,\mathrm{\Omega }),\hfill \\ v(x,y,\mathrm{\Omega })=v_{\mathrm{0}}(x,y,\mathrm{\Omega })+z{\varphi }_y(x,y,\mathrm{\Omega }),\hfill \\ w(x,y,\mathrm{\Omega })=w_{\mathrm{0}}(x,y,\mathrm{\Omega }),\hfill \end{array}

where u is the displacement along the x-axis, u_{\mathrm{0}} is the pure displacement of the mid surface (z = 0) of the plate along the x-axis, v is the displacement along the y-axis, v_{\mathrm{0}} is the pure displacement of the mid surface of the plate along the y-axis, and w is the displacement in the z direction. The terms {\varphi }_xand {\varphi }_yare the pure rotations about the y-axis and x-axis, respectively.

The stochastic linear strains can be extracted as the in-plane and out-of-plane strains, as indicated by Eqs. (8)–(9), respectively:

(8)\left\{\begin{array}{l}{\epsilon }_{xx}\left(\mathrm{\Omega }\right)\hfill \\ {\epsilon }_{yy}\left(\mathrm{\Omega }\right)\hfill \\ {\epsilon }_{xy}\left(\mathrm{\Omega }\right)\hfill \end{array}\right\}=\left\{\begin{array}{l}{\epsilon }_{xx}^{\left(\mathrm{0}\right)}\left(\Omega \right)\hfill \\ {\epsilon }_{yy}^{\left(\mathrm{0}\right)}\left(\mathrm{\Omega }\right)\hfill \\ {\epsilon }_{xy}^{\left(\mathrm{0}\right)}\left(\mathrm{\Omega }\right)\hfill \end{array}\right\}+z\left\{\begin{array}{l}{\epsilon }_{xx}^{\left(\mathrm{1}\right)}\left(\mathrm{\Omega }\right)\hfill \\ {\epsilon }_{yy}^{\left(\mathrm{1}\right)}\left(\mathrm{\Omega }\right)\hfill \\ {\epsilon }_{xy}^{\left(\mathrm{1}\right)}\left(\mathrm{\Omega }\right)\hfill \end{array}\right\},

(9)\left\{\begin{array}{l}{\gamma }_{yz}\left(\mathrm{\Omega }\right)\hfill \\ {\gamma }_{xz}\left(\mathrm{\Omega }\right)\hfill \end{array}\right\}=\left\{\begin{array}{l}{\gamma }_{yz}^{\left(\mathrm{0}\right)}\left(\mathrm{\Omega }\right)\hfill \\ {\gamma }_{xz}^{\left(\mathrm{0}\right)}\left(\mathrm{\Omega }\right)\hfill \end{array}\right\},

where the terms with indices are defined in Eqs. (10)–(12).

(10)\left\{\begin{array}{lll}{\epsilon }_{xx}^{\left(\mathrm{0}\right)}\hfill & {\epsilon }_{yy}^{\left(\mathrm{0}\right)}\hfill & {\epsilon }_{xy}^{\left(\mathrm{0}\right)}\hfill \end{array}\right\}=\left\{\begin{array}{lll}{u}_{\mathrm{0},x}\hfill & v_{\mathrm{0},y}\hfill & u_{\mathrm{0},y}+\mathrm{\hspace{0.17em}}{v}_{\mathrm{0},x}\hfill \end{array}\right\},

(11)\left\{\begin{array}{lll}{\epsilon }_{xx}^{\left(\mathrm{1}\right)}\hfill & {\epsilon }_{yy}^{\left(\mathrm{1}\right)}\hfill & {\epsilon }_{xy}^{\left(\mathrm{1}\right)}\hfill \end{array}\right\}=\left\{\begin{array}{lll}{\varphi }_{x,x}\hfill & {\varphi }_{y,y}\hfill & {\varphi }_{x,y}+\mathrm{\hspace{0.17em}}{\varphi }_{y,x}\hfill \end{array}\right\},

(12)\left\{\begin{array}{ll}{\gamma }_{xz}^{\left(\mathrm{0}\right)}\hfill & {\gamma }_{yz}^{\left(\mathrm{0}\right)}\hfill \end{array}\right\}=\left\{\begin{array}{ll}{\omega }_{\mathrm{0},x}+{\varphi }_x\hfill & {\omega }_{\mathrm{0},y}+{\varphi }_y\hfill \end{array}\right\}.

Constitutive equations

The spatially varying on-axis stochastic stiffness matrix is defined by Eq. (13):

Q_{\mathrm{11}}={\displaystyle \frac{E_{\mathrm{1}}(x,y,\mathrm{\Omega })}{\mathrm{1}-{\nu }_{\mathrm{12}}{\nu }_{\mathrm{21}}}},\mathrm{\hspace{0.17em}}{Q}_{\mathrm{12}}={\displaystyle \frac{{\nu }_{\mathrm{12}}{E}_{\mathrm{2}}(x,y,\mathrm{\Omega })}{\mathrm{1}-{\nu }_{\mathrm{12}}{\nu }_{\mathrm{21}}}}={\displaystyle \frac{{\nu }_{\mathrm{21}}{E}_{\mathrm{1}}(x,y,\mathrm{\Omega })}{\mathrm{1}-{\nu }_{\mathrm{12}}{\nu }_{\mathrm{21}}}},

Q_{\mathrm{22}}={\displaystyle \frac{E_{\mathrm{2}}(x,y,\mathrm{\Omega })}{\mathrm{1}-{\nu }_{\mathrm{12}}{\nu }_{\mathrm{21}}}},\mathrm{\hspace{0.17em}}{Q}_{\mathrm{66}}=G_{\mathrm{12}}(x,y,\mathrm{\Omega }),\mathrm{\hspace{0.17em}}{Q}_{\mathrm{44}}=G_{\mathrm{23}}(x,y,\mathrm{\Omega }),

(13)Q_{\mathrm{55}}=G_{\mathrm{13}}(x,y,\mathrm{\Omega })\mathrm{.}

The constitutive equation in the off-axis coordinates for a single lamina with index k is defined as

{\left\{\begin{array}{l}{\sigma }_{xx}\hfill \\ {\sigma }_{yy}\hfill \\ {\sigma }_{xy}\hfill \end{array}\right\}}^k={\left[\begin{array}{lll}{\overline{Q}}_{\mathrm{11}}\hfill & {\overline{Q}}_{\mathrm{12}}\hfill & {\overline{Q}}_{\mathrm{16}}\hfill \\ {\overline{Q}}_{\mathrm{21}}\hfill & {\overline{Q}}_{\mathrm{22}}\hfill & {\overline{Q}}_{\mathrm{26}}\hfill \\ {\overline{Q}}_{\mathrm{16}}\hfill & {\overline{Q}}_{\mathrm{26}}\hfill & {\overline{Q}}_{\mathrm{66}}\hfill \end{array}\right]}^k\left\{\begin{array}{l}{\epsilon }_{xx}\hfill \\ {\epsilon }_{yy}\hfill \\ {\epsilon }_{xy}\hfill \end{array}\right\},

(14){\left\{\begin{array}{l}{\sigma }_{yz}\hfill \\ {\sigma }_{xz}\hfill \end{array}\right\}}_{}^k={\left[\begin{array}{ll}{\overline{Q}}_{\mathrm{44}}\hfill & {\overline{Q}}_{\mathrm{45}}\hfill \\ {\overline{Q}}_{\mathrm{45}}\hfill & {\overline{Q}}_{\mathrm{55}}\hfill \end{array}\right]}^k\left\{\begin{array}{l}{\gamma }^{\mathrm{0}}{}_{yz}\hfill \\ {\gamma }^{\mathrm{0}}{}_{xz}\hfill \end{array}\right\},

where \overline{Q} is the off-axis stiffness matrix. It should be noted that the off-axis stiffness matrix can be obtained by applying the rotation matrix on the on-axis stiffness matrix for each lamina.

Potential and kinetic energies

Assuming linear elastic strains, the stochastic potential energy is written as follows:

U_{\mathrm{P}}\left(\mathrm{\Omega }\right)={\displaystyle \sum _{i=\mathrm{1}}^{Nlayer}{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}{\displaystyle \underset{V}{\int }({\sigma }_{xx}{\epsilon }_{xx}+{\sigma }_{yy}{\epsilon }_{yy}+{\sigma }_{xy}{\epsilon }_{xy}+K_{\mathrm{s}}{\sigma }_{yz}{\epsilon }_{yz}}}

(15){\displaystyle {\displaystyle +K_{\mathrm{s}}{\sigma }_{xz}{\epsilon }_{xz})}\mathrm{d}V},

where Nlayer is the total number of laminas in the composite plate, and K_{\mathrm{s}}is the shear correction factor in the FSDPT, which is taken to be 5/6 here.

Considering a spatially varying stochastic plate density \rho (x,y,\mathrm{\Omega }), the stochastic kinetic energy can be written as Eq. (16):

T_{\mathrm{P}}\left(\mathrm{\Omega }\right)={\displaystyle \sum _{i=\mathrm{1}}^{Nlayer}{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}{\displaystyle \underset{V}{\int }\rho (x,y,\mathrm{\Omega })}[{\dot{u}}_{}^{\mathrm{2}}(x,y,\mathrm{\Omega })+{\dot{v}}_{}^{\mathrm{2}}(x,y,\mathrm{\Omega })}

(16){\displaystyle +{\dot{w}}_{}^{\mathrm{2}}(x,y,\mathrm{\Omega })]\mathrm{d}V},

Substituting the definitions of u, v, and w from Eq. (7) into Eq. (16) leads to another expression of the kinetic energy:

(17)T_{\mathrm{P}}\left(\mathrm{\Omega }\right)={\displaystyle \sum _{i=\mathrm{1}}^{Nlayer}{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}{\displaystyle \underset{V}{\int }(I_{\mathrm{0}}}\left[{\dot{u}}_{\mathrm{0}}^{\mathrm{2}}+{\dot{v}}_{\mathrm{0}}^{\mathrm{2}}+{\dot{w}}_{\mathrm{0}}^{\mathrm{2}}\right]+I_{\mathrm{1}}\left[\mathrm{2}{\dot{\varphi }}_x^{}{\dot{u}}_{\mathrm{0}}^{}+\mathrm{2}{\dot{v}}_{\mathrm{0}}^{}{\dot{\varphi }}_y^{}\right]+I_{\mathrm{2}}\left[{\dot{\varphi }}_x^{\mathrm{2}}+{\dot{\varphi }}_y^{\mathrm{2}}\right])\mathrm{d}x\mathrm{d}y},

where I_j(x,y,\mathrm{\Omega })are the spatially varying stochastic rotary moments of inertia in the plate, which are defined by Eq. (18).

(18)\begin{array}{l}{I}_0(x,y,\mathrm{\Omega })={\displaystyle {\int }_{-h/\mathrm{2}}^{h/\mathrm{2}}\rho }\mathrm{\hspace{0.17em}}(x,y,\mathrm{\Omega })\mathrm{d}z,\hfill \\ I_{\mathrm{1}}(x,y,\mathrm{\Omega })={\displaystyle {\int }_{-h/\mathrm{2}}^{h/\mathrm{2}}z}\mathrm{\hspace{0.17em}}\rho \mathrm{\hspace{0.17em}}(x,y,\mathrm{\Omega })\mathrm{d}z,\mathrm{\hspace{0.17em}}\hfill \\ I_{\mathrm{2}}(x,y,\mathrm{\Omega })={\displaystyle {\int }_{-h/\mathrm{2}}^{h/\mathrm{2}}{z}^{\mathrm{2}}}\rho \mathrm{\hspace{0.17em}}(x,y,\mathrm{\Omega })\mathrm{d}z.\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}\hfill \end{array}

Stochastic equations of motion

The stochastic equations of motion of free vibration problems can be derived by applying Hamilton’s principle, as given by Eq. (19):

(19)\delta {\displaystyle {\int }_{t_{\mathrm{1}}}^{t_{\mathrm{2}}}(T}-U_{\mathrm{P}})\mathrm{d}t=\mathrm{0.}

The stochastic equations of motion are derived as follows:

(20){\displaystyle \frac{\partial N_{xx}}{\partial x}}+{\displaystyle \frac{\partial N_{xy}}{\partial y}}=I_0{\ddot{u}}_0+I_{\mathrm{1}}{\ddot{\varphi }}_x,

(21){\displaystyle \frac{\partial N_{xy}}{\partial x}}+{\displaystyle \frac{\partial N_{yy}}{\partial y}}=I_0{\ddot{v}}_0+I_{\mathrm{1}}{\ddot{\varphi }}_y,

(22){\displaystyle \frac{\partial Q_x}{\partial x}}+{\displaystyle \frac{\partial Q_y}{\partial y}}=I_0\ddot{w},

(23){\displaystyle \frac{\partial M_{xx}}{\partial x}}+{\displaystyle \frac{\partial M_{xy}}{\partial y}}-Q_x=I_{\mathrm{1}}{\ddot{u}}_0+I_{\mathrm{2}}{\ddot{\varphi }}_x,

(24){\displaystyle \frac{\partial M_{xy}}{\partial x}}+{\displaystyle \frac{\partial M_{yy}}{\partial y}}-Q_y=I_{\mathrm{1}}{\ddot{v}}_0+I_{\mathrm{2}}{\ddot{\varphi }}_y,

where N_{ij}, M_{ij}, and Q_i are stochastic stress resultants and can be written as Eqs. (25)–(27):

(25)\left\{\begin{array}{l}{N}_{xx}(x,y,\mathrm{\Omega })\hfill \\ N_{yy}(x,y,\mathrm{\Omega })\hfill \\ N_{xy}(x,y,\mathrm{\Omega })\hfill \end{array}\right\}={\displaystyle {\int }_{-h/\mathrm{2}}^{h/\mathrm{2}}\left\{\begin{array}{l}{\sigma }_{xx}(x,y,\mathrm{\Omega })\hfill \\ {\sigma }_{yy}(x,y,\mathrm{\Omega })\hfill \\ {\sigma }_{xy}(x,y,\mathrm{\Omega })\hfill \end{array}\right\}}\mathrm{d}z,

(26)\left\{\begin{array}{l}{M}_{xx}(x,y,\mathrm{\Omega })\hfill \\ M_{yy}(x,y,\mathrm{\Omega })\hfill \\ M_{xy}(x,y,\mathrm{\Omega })\hfill \end{array}\right\}={\displaystyle {\int }_{-h/\mathrm{2}}^{h/\mathrm{2}}\left\{\begin{array}{l}{\sigma }_{xx}(x,y,\mathrm{\Omega })\hfill \\ {\sigma }_{yy}(x,y,\mathrm{\Omega })\hfill \\ {\sigma }_{xy}(x,y,\mathrm{\Omega })\hfill \end{array}\right\}}z\mathrm{d}z,

(27)\left\{\begin{array}{l}{Q}_x(x,y,\Omega )\hfill \\ Q_y(x,y,\mathrm{\Omega })\hfill \end{array}\right\}=K_{\mathrm{s}}{\displaystyle {\int }_{-h/\mathrm{2}}^{h/\mathrm{2}}\left\{\begin{array}{l}{\sigma }_{xz}(x,y,\mathrm{\Omega })\hfill \\ {\sigma }_{yz}(x,y,\mathrm{\Omega })\hfill \end{array}\right\}}\mathrm{d}z.

Equations (20)–(24) are the same deterministic equations of motion [⁵⁴] with stochastic terms. The above equations can be discretized using the deterministic Galerkin approach. Equations (20)–(24) can be written in weak forms. The weight functions in the Galerkin method can be assumed equal to the admissible modes or shape functions in the FEM.

For proposing suitable admissible modes, the boundary conditions are specified. The plate is assumed to be simply supported in all edges, as shown in Eq. (28).

(28)\begin{array}{l}x=\mathrm{0},a,\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}W=U={\varphi }_y=\mathrm{0},\hfill \\ y=\mathrm{0},b,\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}W=V={\varphi }_x=\mathrm{0.}\hfill \end{array}

The admissible modes are defined for each degree of freedom considering essential boundary conditions as follows:

(29)\begin{array}{c}\hfill {\left\{u_{\mathrm{0}}\right\}}_{m,n}=U_{mn}\left(\mathrm{\Omega }\right)\sin \left({\displaystyle \frac{n\mathrm{\pi }}{l}}x\right)\cos \left({\displaystyle \frac{m\mathrm{\pi }}{l}}y\right)\cos {\omega }_{mn}t,\hfill \\ \hfill n,m=\mathrm{1,...,}N,\hfill \\ {\left\{v_{\mathrm{0}}\right\}}_{m,n}=V_{mn}\left(\mathrm{\Omega }\right)\cos \left({\displaystyle \frac{n\mathrm{\pi }}{l}}x\right)\sin \left({\displaystyle \frac{m\mathrm{\pi }}{l}}y\right)\cos {\omega }_{mn}t\mathrm{,}\hfill \\ n,m=\mathrm{1,...,}N,\\ {\left\{w_{\mathrm{0}}\right\}}_{m,n}=W_{mn}\left(\mathrm{\Omega }\right)\sin \left({\displaystyle \frac{n\mathrm{\pi }}{l}}x\right)\sin \left({\displaystyle \frac{m\mathrm{\pi }}{l}}y\right)\cos {\omega }_{mn}t\mathrm{,}\hfill \\ n,m=\mathrm{1,...,}N,\\ {\left\{{\varphi }_x\right\}}_{m,n}=X_{mn}\left(\mathrm{\Omega }\right)\cos \left({\displaystyle \frac{n\mathrm{\pi }}{l}}x\right)\sin \left({\displaystyle \frac{m\mathrm{\pi }}{l}}y\right)\cos {\omega }_{mn}t\mathrm{,}\hfill \\ n,m=\mathrm{1,...,}N,\\ {\left\{{\varphi }_y\right\}}_{m,n}=Y_{mn}\left(\mathrm{\Omega }\right)\sin \left({\displaystyle \frac{n\mathrm{\pi }}{l}}x\right)\cos \left({\displaystyle \frac{m\mathrm{\pi }}{l}}y\right)\cos {\omega }_{mn}t\mathrm{,}\hfill \\ n,m=\mathrm{1,...,}N,\end{array}

where U_{mn}, V_{mn}, W_{mn}, X_{mn}, and Y_{mn} are unknown stochastic coefficients in generalized coordinates, \mathrm{\Omega } is the probabilistic vector in Eq. (2), and {\mathrm{\Omega }}_{mn} is the frequency of harmonic response. Since there is no semi-analytical solution for free vibration of plates with spatially varying stochastic properties in the literature, the SFEM is employed for verifying the proposed semi-analytical method to ensure that the set of admissible modes (Eq. (29)) is a suitable spanning set for the stochastic modal space. Utilizing the iso-parametric 8-node serendipity element for this purpose, the displacement field in the SFEM can be written as Eq. (30):

(30)\left\{\mathrm{\Lambda }\right\}={\displaystyle \sum _{i=\mathrm{1}}^{NN}{N}_i}(\eta ,\mathit{\xi }){\left\{\mathrm{\Lambda }\right\}}_i\left(\mathrm{\Omega }\right),

where \left\{\mathrm{\Lambda }\right\}is a vector with five degrees of freedom in the FSDPT, {\left\{\mathrm{\Lambda }\right\}}_i\left(\mathrm{\Omega }\right)is the degree of freedom of node i, N_iis the shape function with respect to natural coordinates (\eta ,\mathit{\xi }), and NN is the number of nodes, which is equal to 8 for the present SFEM. The standard shape functions are written as Eq. (31):

(31)\begin{array}{c}\hfill N_j={\displaystyle \frac{\mathrm{1}}{\mathrm{4}}}(\mathrm{1}+{\mathit{\xi }}_j\mathit{\xi })(\mathrm{1}+{\eta }_j\eta )({\mathit{\xi }}_j\mathit{\xi }+{\eta }_j\eta -\mathrm{1}),\hfill \\ \hfill \mathrm{for}\mathrm{\hspace{0.17em}}\mathrm{corner}\mathrm{\hspace{0.17em}}\mathrm{nodes}\mathrm{,}\mathrm{\hspace{0.17em}}j=\mathrm{1},\mathrm{2},\mathrm{3},\mathrm{4}\hfill \\ \hfill N_j={\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}(\mathrm{1}-{\mathit{\xi }}^{\mathrm{2}})(\mathrm{1}+{\eta }_j\eta ),\hfill \\ \hfill \mathrm{for}\mathrm{\hspace{0.17em}}\mathrm{mid-side}\mathrm{\hspace{0.17em}}\mathrm{nodes}\mathrm{,}\mathrm{\hspace{0.17em}}j=\mathrm{5},\mathrm{7}\hfill \\ \hfill N_j={\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}(\mathrm{1}+{\mathit{\xi }}_j\mathit{\xi })(\mathrm{1}-{\eta }^{\mathrm{2}}),\mathrm{\hspace{0.17em}}\hfill \\ \hfill \mathrm{for}\mathrm{\hspace{0.17em}}\mathrm{mid-side}\mathrm{\hspace{0.17em}}\mathrm{nodes}\mathrm{,}\mathrm{\hspace{0.17em}}j=\mathrm{6},\mathrm{8},\hfill \end{array}

where {\mathit{\xi }}_jand {\eta }_jare initialized natural coordinates of nodes in the standard 8-node serendipity element.

Applying the deterministic Galerkin method on the weak forms of Eqs. (20)–(24) by utilizing the two proposed approaches leads to the stochastic discretized equations of motion in the generalized coordinate space, as indicated by Eq. (32).

(\overline{M}+{\displaystyle \sum _{i=\mathrm{1}}^n{M}_i{\mathit{\xi }}_i^{\rho }})\left\{\ddot{X}\right\}\left(\mathrm{\Omega }\right)

(32)+(\overline{K}+{\displaystyle \sum _{i=\mathrm{1}}^n{K}_i{\mathit{\xi }}_i^{E_{\mathrm{11}}}+{\displaystyle \sum _{i=\mathrm{1}}^n{K}_i{\mathit{\xi }}_i^{G_{\mathrm{12}}}+{\displaystyle \sum _{i=\mathrm{1}}^n{K}_i{\mathit{\xi }}_i^{G_{\mathrm{23}}}}}})\left\{X\right\}\left(\mathrm{\Omega }\right)=\mathrm{0},

where \overline{M}and \overline{K}are the averages of the mass and stiffness matrices, respectively; M_iand K_iare the mass and elastic stiffness matrices corresponding to the terms in the KL expansion, respectively; and {\mathit{\xi }}_i^{\rho }are the standard random variables in Eq. (2). Assuming a pure harmonic response, Eq. (32) can be written in the form of a standard eigenvalue problem as follows:

(\overline{K}+{\displaystyle \sum _{i=\mathrm{1}}^n{K}_i{\mathit{\xi }}_i^{E_{\mathrm{11}}}+{\displaystyle \sum _{i=\mathrm{1}}^n{K}_i{\mathit{\xi }}_i^{G_{\mathrm{12}}}+{\displaystyle \sum _{i=\mathrm{1}}^n{K}_i{\mathit{\xi }}_i^{G_{\mathrm{23}}}}}})\left\{X\left(\mathrm{\Omega }\right)\right\}

(33)=\lambda (\mathrm{\Omega })(\overline{M}+{\displaystyle \sum _{i=\mathrm{1}}^n{M}_i{\mathit{\xi }}_i^{\rho }})\left\{X\left(\mathrm{\Omega }\right)\right\},

where \lambda \left(\mathrm{\Omega }\right) is the stochastic eigenvalue (square of natural frequency) and X\left(\mathrm{\Omega }\right)is the stochastic eigenvector. This equation can be solved using various methods: direct Monte Carlo (which is computationally expensive), metamodel-based approaches, etc. In the present research, the intrusive polynomial chaos method is utilized to solve Eq. (33) by applying the stochastic Galerkin projection. This method is advantageous because the stochastic modal space is obtained and the computational cost is much lower than in the Monte Carlo simulation. The intrusive stochastic Galerkin projection method is currently applied in the SFEM-based approach [⁵⁵]. Because of the relatively large mathematical model size in the FEM, combining it with the intrusive stochastic Galerkin approach is not computationally efficient. In the present research, the stochastic Galerkin approach is utilized in combination with the admissible modes method. This combination has a much smaller mathematical model size, which leads to an efficient semi-analytical non-sampling-based approach. Equation (33) is also solved by the FEM-based approach for verification purposes. In the case of the SFEM, stochastic simulations are carried out using the Monte Carlo method to verify the proposed semi-analytical non-sampling approach.

Stochastic Galerkin approach

In this section, the stochastic Galerkin projection is described and the stochastic eigenvalues and eigenvectors will be determined using a non-sampling-based approach. Equation (33) can be written in a simple form as Eq. (34).

(34)K\left(\mathrm{\Omega }\right)\left\{X\right(\mathrm{\Omega }\left)\right\}=\lambda \left(\mathrm{\Omega }\right)M\left(\mathrm{\Omega }\right)\left\{X\right(\mathrm{\Omega }\left)\right\}.

Equation (34) can be derived for shifted eigenvalues by employing the eigenvalue shift theorem, as given by Eq. (35) [⁵⁵].

(35)\begin{array}{l}\left(K\right(\mathrm{\Omega })-\eta (\mathrm{\Omega }\left)M\right(\mathrm{\Omega }\left)\right)X\left(\mathrm{\Omega }\right)=\widehat{\lambda }\left(\mathrm{\Omega }\right)M\left(\mathrm{\Omega }\right)X\left(\mathrm{\Omega }\right),\hfill \\ \widehat{\lambda }\left(\mathrm{\Omega }\right)=\lambda \left(\mathrm{\Omega }\right)-\eta \left(\mathrm{\Omega }\right),\hfill \end{array}

where \widehat{\lambda } is the shifted stochastic eigenvalue, and \etais the stochastic shift parameter. It should be noted here that the first eigenvalue \eta is equal to zero.

The stochastic eigenvalues \lambda \left(\mathrm{\Omega }\right)and \widehat{\lambda }\left(\mathrm{\Omega }\right), eigenvector \left\{X\right(\mathrm{\Omega }\left)\right\}, and stochastic shift parameter \eta \left(\mathrm{\Omega }\right) can be expanded in the probability space using the polynomial chaos expansion as follows:

(36)\mathrm{\hspace{0.17em}}\lambda \left(\Omega \right)={\displaystyle \sum _{i=\mathrm{0}}^P{\lambda }_i}{\psi }_i\left(\mathrm{\Omega }\right),

(37)\mathrm{\hspace{0.17em}}\widehat{\lambda }\left(\mathrm{\Omega }\right)={\displaystyle \sum _{i=\mathrm{0}}^P{\widehat{\lambda }}_i}{\psi }_i\left(\mathrm{\Omega }\right),\mathrm{\hspace{0.17em}}\mathrm{\hspace{0.17em}}

(38)X\left(\mathrm{\Omega }\right)={\displaystyle \sum _{i=\mathrm{0}}^P{x}_i}{\psi }_i\left(\mathrm{\Omega }\right),

(39)\eta \left(\mathrm{\Omega }\right)={\displaystyle \sum _{i=\mathrm{0}}^P{\eta }_i}{\psi }_i\left(\mathrm{\Omega }\right),

where {\lambda }_i,{\widehat{\lambda }}_i,{\eta }_i, and x_iare unknown coefficients, which should be found; {\psi }_i({\mathit{\xi }}_{\mathrm{1}}^{},{\mathit{\xi }}_{\mathrm{2}}^{},\mathrm{...},{\mathit{\xi }}_n^{})are orthogonal probabilists’ Hermite polynomials of order t, which can be extracted as follows [⁵²]:

(40)\begin{array}{l}{\psi }_t({\mathit{\xi }}_{i_1},\mathrm{...},{\mathit{\xi }}_{i_t})=e^{\frac{1}{2}{\mathit{\xi }}^{\mathrm{T}}\mathit{\xi }}{(-1)}^t\frac{{\partial }^t}{\partial {\mathit{\xi }}_{i_1}\mathrm{\cdots }\partial {\mathit{\xi }}_{i_t}}{e}^{-\frac{1}{2}{\mathit{\xi }}^{\mathrm{T}}\mathit{\xi }},\hfill \\ \mathit{\xi }={[{\mathit{\xi }}_1,\mathrm{...},{\mathit{\xi }}_n]}^{\mathrm{T}}.\hfill \end{array}

The number of terms in Eqs. (36–39) is equal to P + 1 and P can be calculated as follows:

(41)P={\displaystyle \frac{(d+t)!}{d!t!}}-\mathrm{1},

where d is the dimensionality (probabilistic space) or number of stochastic terms in the KL expansion (n) multiplied by the number of independent stochastic fields (equal to 4 in the present research).

Based on the stochastic Galerkin projection, both sides of Eq. (35) should be multiplied by weight functions that are applied in expanding the solution. Here, the weight functions are the same, i.e., {\psi }_i({\mathit{\xi }}_{\mathrm{1}}^{},{\mathit{\xi }}_{\mathrm{2}}^{},\mathrm{...},{\mathit{\xi }}_n^{}), which are used in the definition of eigenvalues and eigenvectors. By substituting Eqs. (36)–(39) into Eq. (35) and multiplying {\psi }_i({\mathit{\xi }}_{\mathrm{1}}^{},{\mathit{\xi }}_{\mathrm{2}}^{},\mathrm{...},{\mathit{\xi }}_n^{}) on both sides of Eq. (35), Eq. (42) is obtained by applying the probabilistic expected value operator.

(42)k=\mathrm{0},\mathrm{...},P\mathrm{.}

The above equations can be transformed into a simple form as follows:

(43)(\mathbb{K}-{\mathbb{K}}_r)\mathbb{X}=\widehat{\lambda }\mathbb{M}\mathbb{X}.

Equation (43) is the same equation derived by Salim et al. for SFEM-based approaches [⁵⁵]. An iterative procedure is applied for solving Eq. (43). By determining {\lambda }_i,{\widehat{\lambda }}_i,{\eta }_i, and x_i, the stochastic eigenvalues and eigenvectors are found. Assuming m degrees of freedom for the mass and stiffness matrices in Eq. (35), the total degree of freedom for generalized mass matrix \mathbb{M}and generalized stiffness matrix \mathbb{K}in Eq. (43) is equal to m(P+1). Hence, the computational efficiency of the proposed method is strongly dependent on the value of m and P.

Owing to the orthogonality of Hermite polynomials, the mean and standard deviations of the eigenvalues and eigenvectors can be found easily.

Verification

Two problems were solved to verify the presented polynomial chaos-based semi-analytical approach (PC-SAA). In the first problem, the two proposed approaches (SFEM and PC-SAA) were utilized to solve the free vibration of a composite plate with stochastic properties as random variables. Subsequently, the results were compared with available results in the literature. A square composite plate with a/h = 10 and stochastic tensile modulus (E₁₁) as a random variable was considered. The other mechanical and physical properties were considered to be deterministic and are presented in Table 1.

To simulate the stochastic property as a stochastic variable, it was assumed that n = 1 in the KL expansion and the eigenvalues and eigenvectors were equal to 1 in Eq. (2). In consideration of 16 modes for each degree of freedom in the PC-SAA, a 10 × 10 mesh was sufficient in the SFEM to obtain agreement between the results of the COV of the first eigenvalue for different COVs of tensile modulus. Fourth-order Hermite polynomials (t = 4) were applied to achieve suitable accuracy in the results. The COV of the first eigenvalue extracted from the two presented approaches and the results from Ref. [⁵⁶] are plotted in Fig. 2.

Figure 2 shows good agreement between the results. Hence, the formulation and computer codes were verified for the plate with stochastic properties as random variables. The number of degrees of freedom m in the above simulation is 80 and 1457 for the PC-SAA and SFEM, respectively, considering the boundary condition (Eq. (28)).

For the plates with spatially varying properties, there is no similar case in the literature. The SFEM based on the Monte Carlo simulation was applied to ensure that the proposed set of admissible modes in Eq. (29) is a suitable span set for the stochastic modal space. The free vibration of a square composite plate with unit length a/h equal to 10 and 100 and [0,90,90,0] stacking sequence was investigated using the two presented methods. In this step, the average values of the applied material properties are listed in Table 2.

All material properties including tensile modulus, shear modulus, and density of the plate were considered to be spatially varying and modeled as random fields of Eq. (2) with exponential kernels of Eq. (1). The correlation lengths in both directions were assumed to be 2 m. Three terms in the KL expansion (n = 3) were enough to satisfy the criterion (Eq. (6)). The SFEM based on Monte Carlo simulation was applied with 10000 samples to reach convergence of the COV of the first eigenvalue (variation of less than 1%). The number of applied modes in the PC-SAA and mesh size in the SFEM were the same as in the previous problem. Fourth-order Hermite polynomials (t = 4) were applied to achieve suitable accuracy in the results. The COV of the first eigenvalue extracted from the two presented methods for different uncertainty levels in random properties are displayed in Fig. 3.

Good agreement in the results implies the efficiency of the proposed semi-analytical approach. It should be noted here that the number of degrees of freedom m is equal to 80 and 1457 for the PC-SAA and SFEM, respectively. Considering four independent fields and n = 3, the probabilistic dimensionality is equal to 12. Assuming t = 4, the number of terms in the polynomial chaos expansion in Eqs. (36)–(39), P+1, is equal to 140. Hence, the total degree of freedom of the generalized mass and stiffness matrices (M, K) in Eq. (43) is equal to 140 × 80 for the PC-SAA. If we apply the stochastic Galerkin projection method to the FEM-based approach, this value is equal to 140 × 1457. This demonstrates the computational efficiency of the PC-SAA in comparison to the SFEM-based approach.

Case studies

Plates with properties listed in Table 2, unit length, and different side to thickness ratios (a/h) with [0,90,90,0] stacking sequence were considered as case studies. The correlation length in both directions was equal to 2 m. Three terms in the KL expansion were employed. At first, the sensitivity of the COV of the first natural frequency of the composite plate due to uncertainty propagation in the spatially varying mechanical properties was evaluated individually. The results of the uncertainty propagation in the first natural frequency due to uncertainties in the tensile modulus and density of the plates are illustrated in Fig. 4.

From Fig. 4, it can be observed that the relation between the COV of natural frequencies and COV of random properties remains linear considering spatial stochasticity of the mechanical properties. This is the same as the observed relation considering mechanical properties as random variables (see Fig. 2). Thinner plates are more sensitive owing to the spatial stochasticity of tensile modulus. However, in the case of density, the aforementioned sensitivity remains constant for different (a/h) values. The results reveal that the uncertainty propagation in the density of plates typically has greater effects on the stochastic natural frequencies than that in the tensile modulus. It is important to consider this issue when mechanical experiments are required to assess the structural reliability of structures.

Similar results are depicted in Fig. 5 for the in-plane and out-of-plane shear moduli.

As in the previous results, the relation remains linear but the stochastic natural frequency of thicker plates is more sensitive to the uncertainty of shear modulus. The results indicate that the in-plane shear modulus has greater effect on the stochastic natural frequency than the out-of-plane shear modulus. To show the importance of the uncertainties in each random property and compare their effects on the stochastic natural frequency individually, Fig. 6 is presented.

It can be observed that the spatial stochasticities of the plate density and tensile modulus have the greatest effects on the stochastic natural frequency. The out-of-plane shear modulus has the least effect on the uncertainty of natural frequency. It is important to consider this sensitivity analysis when mechanical experiments are required to find the spatial stochasticity of mechanical properties for possible post-processing such as in reliability estimation or health assessment of structures.

The effects of the combination of uncertainty propagation in all aforementioned random properties on the natural frequencies were determined, as displayed in Fig. 7.

It can be seen from the results that the relation between the COV of all spatially varying random properties, which were considered simultaneously, and the COV of the first and second natural frequencies remains linear. The order of uncertainty propagation in the natural frequencies is similar in both the first and second natural frequencies.

As mentioned before, uncertainty in vibration modes is important especially when these modes are required for assessing the structural health or finding the location of possible damages. The proposed method has good capability of providing the stochastic modes for any type of vibration post-processing. For example, vibration modes are necessary to extract modal strain energy when structural health monitoring is based on the vibration response of structures. For plates with the previously defined geometries and material properties and a/h = 50, the first five transversal vibration modes will be presented. It should be noted that the presented modes were normalized with respect to the stochastic mass matrix. The contour of deterministic modes and three typical stochastic samples of the first five modes (transversal) are listed in Table 3.

The results show that spatially varying stochastic properties lead to stochastic vibrational modes that can have completely different shapes in comparison to the deterministic mode shapes. The sensitivity of the first two modes with respect to uncertainties in properties is lower than that of higher-order modes. Hence, the sensitivity of higher-order modes should be considered when the spatial stochasticity is under investigation with a vibration-based approach. The uncertainty level and correlation length of stochastic fields can have significant effects on the stochastic mode shapes. Higher-order vibrational modes are required to track spatially varying stochastic mechanical properties.

Conclusions

A stochastic free vibration analysis of composite plates with spatially varying stochastic properties was conducted. A deterministic-stochastic Galerkin-based semi-analytical approach was developed to solve the stochastic equations of motion. In comparison with the SFEM, the size of the mathematical model of the structure is reduced significantly by applying the proposed semi-analytical approach. The stochastic Galerkin projection based on the intrusive polynomial chaos approach was utilized to extract the stochastic modal space. The reduction in the size of mathematical models led to an enhancement in the computational efficiency of the stochastic simulations compared to the SFEM. The proposed approach can present the stochastic modal space. The sensitivity of stochastic natural frequencies with respect to uncertainty propagation in the spatially varying mechanical properties including tensile modulus, shear modulus, and density of the plates was studied. Uncertainties in the tensile modulus and density of the plates have the greatest effect on the stochastic natural frequencies, while the out-of-plane shear modulus has the least effect. Considering Gaussian field properties, there is still a linear relation between the COV of natural frequencies and COV of random properties. The results indicate that spatially varying stochastic properties can significantly affect the vibrational mode shapes, and higher-order mode shapes are more sensitive to uncertainty propagation in random field material properties.

References

[1] Hamdia K M, Rabczuk T. Key parameters for fracture toughness of particle/polymer nanocomposites: Sensitivity analysis via XFEM modeling approach. In: Proceedings of the 7th International Conference on Fracture Fatigue and Wear. FFW 2018. Lecture Notes in Mechanical Engineering. Singapore: Springer, 2019

[2] Vu-Bac N, Lahmer T, Keitel H, Zhao J, Zhuang X, Rabczuk T. Stochastic predictions of bulk properties of amorphous polyethylene based on molecular dynamics simulations. Mechanics of Materials, 2014, 68: 70–84

[3] Vu-Bac N, Lahmer T, Zhang Y, Zhuang X, Rabczuk T. Stochastic predictions of interfacial characteristic of polymeric nanocomposites (PNCs). Composites. Part B, Engineering, 2014, 59: 80–95

[4] Vu-Bac N, Silani M, Lahmer T, Zhuang X, Rabczuk T. A unified framework for stochastic predictions of mechanical properties of polymeric nanocomposites. Computational Materials Science, 2015, 96: 520–535

[5] Vu-Bac N, Rafiee R, Zhuang X, Lahmer T, Rabczuk T. Uncertainty quantification for multiscale modeling of polymer nanocomposites with correlated parameters. Composites. Part B, Engineering, 2015, 68: 446–464

[6] Vu-Bac N, Zhuang X, Rabczuk T. Uncertainty quantification for mechanical properties of polyethylene based on fully atomistic model. Materials (Basel), 2019, 12(21): 3613

[7] Vu-Bac N, Lahmer T, Zhuang X, Nguyen-Thoi T, Rabczuk T. A software framework for probabilistic sensitivity analysis for computationally expensive models. Advances in Engineering Software, 2016, 100: 19–31

[8] Sasikumar P, Suresh R, Vijayaghosh P K, Gupta S. Experimental characterisation of random field models for CFRP composite panels. Composite Structures, 2015, 120: 451–471

[9] Sriramula S, Chryssanthopoulos M K. An experimental characterisation of spatial variability in GFRP composite panels. Structural Safety, 2013, 42: 1–11

[10] Umesh K, Ganguli R. Matrix crack detection in composite plate with spatially random material properties using fractal dimension. Computers, Materials & Continua, 2014, 41(3): 215–239

[11] Sasikumar P, Suresh R, Gupta S. Analysis of CFRP laminated plates with spatially varying non-Gaussian inhomogeneities using SFEM. Composite Structures, 2014, 112: 308–326

[12] Sasikumar P, Suresh R, Gupta S. Stochastic finite element analysis of layered composite beams with spatially varying non-Gaussian inhomogeneities. Acta Mechanica, 2014, 225(6): 1503–1522

[13] Sriramula S, Chryssanthopoulos M K. Quantification of uncertainty modelling in stochastic analysis of FRP composites. Composites. Part A, Applied Science and Manufacturing, 2009, 40(11): 1673–1684

[14] Naskar S, Mukhopadhyay T, Sriramula S. Probabilistic micromechanical spatial variability quantification in laminated composites. Composites. Part B, Engineering, 2018, 151: 291–325

[15] Naskar S, Mukhopadhyay T, Sriramula S, Adhikari S. Stochastic natural frequency analysis of damaged thin-walled laminated composite beams with uncertainty in micromechanical properties. Composite Structures, 2017, 160: 312–334

[16] Navaneetha Raj B, Iyengar N G R, Yadav D. Response of composite plates with random material properties using FEM and Monte Carlo simulation. Advanced Composite Materials, 1998, 7(3): 219–237

[17] Singh B N, Bisht A K S, Pandit M K, Shukla K K. Nonlinear free vibration analysis of composite plates with material uncertainties: A Monte Carlo simulation approach. Journal of Sound and Vibration, 2009, 324(1–2): 126–138

[18] Piovan M T, Ramirez J M, Sampaio R. Dynamics of thin-walled composite beams: Analysis of parametric uncertainties. Composite Structures, 2013, 105: 14–28

[19] Zhang S, Zhang L, Wang Y, Tao J, Chen X. Effect of ply level thickness uncertainty on reliability of laminated composite panels. Journal of Reinforced Plastics and Composites, 2016, 35(19): 1387–1400

[20] Butler R, Dodwell T J, Haftka R T, Kim N H, Kim T, Kynaston S, Scheichl R. Uncertainty quantification of composite structures with defects using multilevel Monte Carlo simulations. In: The 17th AIAA Non-Deterministic Approaches Conference. Kissimmee, FL: AIAA, 2015

[21] Lal A, Singh B N. Stochastic free vibration of laminated composite plates in thermal environments. Journal of Thermoplastic Composite Materials, 2010, 23(1): 57–77

[22] Onkar A K, Yadav D. Non-linear free vibration of laminated composite plate with random material properties. Journal of Sound and Vibration, 2004, 272(3–5): 627–641

[23] Onkar A K, Yadav D J C S. Forced nonlinear vibration of laminated composite plates with random material properties. Composite Structures, 2005, 70(3): 334–342

[24] Onkar A K, Yadav D. Non-linear response statistics of composite laminates with random material properties under random loading. Composite Structures, 2003, 60(4): 375–383

[25] Singh B N, Iyengar N G R, Yadav D. Effects of random material properties on buckling of composite plates. Journal of Engineering Mechanics, 2001, 127(9): 873–879

[26] Chaudhuri A, Chakraborty S. Reliability of linear structures with parameter uncertainty under non-stationary earthquake. Structural Safety, 2006, 28(3): 231–246

[27] Grover N, Sahoo R, Singh B N, Maiti D K. Influence of parametric uncertainties on the deflection statistics of general laminated composite and sandwich plates. Composite Structures, 2017, 171: 158–169

[28] Lal A, Palekar S P. Stochastic fracture analysis of laminated composite plate with arbitrary cracks using X-FEM. International Journal of Mechanics and Materials in Design, 2017, 13(2): 195–228

[29] Sakata S, Okuda K, Ikeda K. Stochastic analysis of laminated composite plate considering stochastic homogenization problem. Frontiers of Structural and Civil Engineering, 2015, 9(2): 141–153

[30] da Silva C R Á Jr, Beck A T. Bending of stochastic Kirchhoff plates on Winkler foundations via the Galerkin method and the Askey-Wiener scheme. Probabilistic Engineering Mechanics, 2010, 25(2): 172–182

[31] Dey S, Mukhopadhyay T, Khodaparast H H, Adhikari S. Fuzzy uncertainty propagation in composites using Gram-Schmidt polynomial chaos expansion. Applied Mathematical Modelling, 2016, 40(7–8): 4412–4428

[32] Sepahvand K, Marburg S, Hardtke H J. Stochastic free vibration of orthotropic plates using generalized polynomial chaos expansion. Journal of Sound and Vibration, 2012, 331(1): 167–179

[33] Umesh K, Ganguli R. Material uncertainty effect on vibration control of smart composite plate using polynomial chaos expansion. Mechanics of Advanced Materials and Structures, 2013, 20(7): 580–591

[34] Sepahvand K, Scheffler M, Marburg S J A A. Uncertainty quantification in natural frequencies and radiated acoustic power of composite plates: Analytical and experimental investigation. Applied Acoustics, 2015, 87: 23–29

[35] Scarth C, Cooper J E, Weaver P M, Silva G H. Uncertainty quantification of aeroelastic stability of composite plate wings using lamination parameters. Composite Structures, 2014, 116: 84–93

[36] Chen X, Qiu Z. A novel uncertainty analysis method for composite structures with mixed uncertainties including random and interval variables. Composite Structures, 2018, 184: 400–410

[37] Hamdia K M, Silani M, Zhuang X, He P, Rabczuk T. Stochastic analysis of the fracture toughness of polymeric nanoparticle composites using polynomial chaos expansions. International Journal of Fracture, 2017, 206(2): 215–227

[38] Dey S, Mukhopadhyay T, Sahu S K, Li G, Rabitz H, Adhikari S. Thermal uncertainty quantification in frequency responses of laminated composite plates. Composites. Part B, Engineering, 2015, 80: 186–197

[39] Dey S, Mukhopadhyay T, Spickenheuer A, Adhikari S, Heinrich G. Bottom up surrogate based approach for stochastic frequency response analysis of laminated composite plates. Composite Structures, 2016, 140: 712–727

[40] Dey S, Mukhopadhyay T, Adhikari S. Stochastic free vibration analysis of angle-ply composite plates—A RS-HDMR approach. Composite Structures, 2015, 122: 526–536

[41] Mukhopadhyay T, Naskar S, Karsh P K, Dey S, You Z. Effect of delamination on the stochastic natural frequencies of composite laminates. Composites. Part B, Engineering, 2018, 154: 242–256

[42] Hamdia K M, Lahmer T, Nguyen-Thoi T, Rabczuk T. Predicting the fracture toughness of PNCs: A stochastic approach based on ANN and ANFIS. Computational Materials Science, 2015, 102: 304–313

[43] Dey S, Naskar S, Mukhopadhyay T, Gohs U, Spickenheuer A, Bittrich L, Sriramula S, Adhikari S, Heinrich G. Uncertain natural frequency analysis of composite plates including effect of noise—A polynomial neural network approach. Composite Structures, 2016, 143: 130–142

[44] Dey S, Mukhopadhyay T, Spickenheuer A, Gohs U, Adhikari S. Uncertainty quantification in natural frequency of composite plates-An Artificial neural network based approach. Advanced Composites Letters, 2016, 25(2): 43–48

[45] Mukhopadhyay T, Naskar S, Dey S, Adhikari S. On quantifying the effect of noise in surrogate based stochastic free vibration analysis of laminated composite shallow shells. Composite Structures, 2016, 140: 798–805

[46] Dey S, Mukhopadhyay T, Adhikari S. Stochastic free vibration analyses of composite shallow doubly curved shells—A Kriging model approach. Composites. Part B, Engineering, 2015, 70: 99–112

[47] Ghasemi H, Rafiee R, Zhuang X, Muthu J, Rabczuk T. Uncertainties propagation in metamodel-based probabilistic optimization of CNT/polymer composite structure using stochastic multi-scale modeling. Computational Materials Science, 2014, 85: 295–305

[48] Dey S, Mukhopadhyay T, Adhikari S. Metamodel based high-fidelity stochastic analysis of composite laminates: A concise review with critical comparative assessment. Composite Structures, 2017, 171: 227–250

[49] Pichler L, Pradlwarter H J, Schuëller G I. A mode-based meta-model for the frequency response functions of uncertain structural systems. Computers & Structures, 2009, 87(5–6): 332–341

[50] Vu-Bac N, Duong T X, Lahmer T, Zhuang X, Sauer R A, Park H S, Rabczuk T. A NURBS-based inverse analysis for reconstruction of nonlinear deformations of thin shell structures. Computer Methods in Applied Mechanics and Engineering, 2018, 331: 427–455

[51] Vu-Bac N, Duong T X, Lahmer T, Areias P, Sauer R A, Park H S, Rabczuk T. A NURBS-based inverse analysis of thermal expansion induced morphing of thin shells. Computer Methods in Applied Mechanics and Engineering, 2019, 350: 480–510

[52] Ghanem R G, Spanos P D. Spectral stochastic finite-element formulation for reliability analysis. Journal of Engineering Mechanics, 1991, 117(10): 2351–2372

[53] Fish J, Wu W. A nonintrusive stochastic multiscale solver. International Journal for Numerical Methods in Engineering, 2011, 88(9): 862–879

[54] Reddy J N. Mechanics of Laminated Composite Plates and Shells: Theory and Analysis. New York: CRC press, 2003

[55] Mulani S, Kapania R, Walters R. Stochastic eigenvalue problem with polynomial chaos. In: The 47th AIAA/ASME/ASCE/AHS/ASC Structures, Structural Dynamics, and Materials Conference. Newport, RI: AIAA, 2006

[56] Singh B N, Yadav D, Iyengar N G R. A C° element for free vibration of composite plates with uncertain material properties. Advanced Composite Materials, 2002, 11(4): 331–350

RIGHTS & PERMISSIONS

Higher Education Press