Three-dimensional semi-analytical solutions for the transient response of functionally graded material cylindrical panels with various boundary conditions

Abstract

In this paper, a 3D semi-analytical method is proposed by introducing the Durbin’s Laplace transform, as well as its numerical inversion method, state space approach and differential quadrature method to analyse the transient behaviour of functionally graded material cylindrical panels. Moreover, to investigate the effectiveness of the proposed semi-analytical solution, four boundary conditions are used to undertake the analyses. Comparing the proposed approach with other theoretical methods from the literatures, we see better agreements in the natural frequencies. Besides, the semi-analytical solution acquires nearly the same transient response as those obtained by ANSYS. Convergence studies indicate that the proposed method has a quick convergence rate with growing sample point numbers along the length direction, so do layer numbers increase along the radial direction. The effects of thickness/outer radius ratio, length/outer radius ratio and functionally graded indexes are also studied. When carbon nanotube is added to functionally graded material cylindrical panel, the composite structures have been reinforced greatly. The proposed 3D semi-analytical method has high accuracy for the analysis of composite structures. This study can serve as a foundation for solving more complicated environments such as fluid–structure interaction of flexible pipe or thermal effect analysis of functionally graded material in aerospace field.

Introduction

Functionally graded materials (FGMs), whose property gradient is caused by chemical composition, atomic order or microstructure, have attracted much attention of many groups.^1–4 Because of their advantageous stiffness-to-weight ratio and strength-to-weight ratio as well as their tendency for high performance, FGM structures play important roles in ocean engineering, fuselage and submarine components. Therefore, the research on the dynamic response of FGM structure is of great significance to the development of the frontier in aerospace engineering, civil engineering, ocean engineering, etc.

Over the course of history, various plates/shells theories have been used to analyse FGM structures. Using the framework of the non-local strain gradient theory and Hamilton principle, Li et al. ⁵ deduced the governing equations and analysed the natural frequencies of FGM beam. Parida and Mohanty ⁶ investigated the free vibration of FGM plates by the use of higher order shear deformation plate theory (HSPT) on the foundation of Winkler–Pasternak. To draw a more valid conclusion, Cheng and Batra ⁷ derived equations for an FG plate whose responses were retrieved by either the FSDT or third-order shear deformation theory (TSDT). Based on two kinds of shell theories and Higher-order shear theory, Frikha and Dammak ⁸ studied the nonlinear mechanical response of FG shells. Similarly, Yaghoubshahi et al. ⁹ deduced the governing equations of laminated plate problems under static load by using HSDT and virtual work principle. Quan and Duc ¹⁰ investigated the dynamic response and nonlinear vibration of FGM thick shells through the use of third-order shear deformation shell theory. Mantari et al. ¹¹ developed a new HSDT for elastic composite/sandwich plates and shells. In the case of tangential stress-free boundary conditions, the theory provides a holistic explanation for a full distribution of the transverse shear strains. Consequently, a shear correction factor was not needed. Generally speaking, the CPT is appropriate for the analysis of thin structures, while HSDT and TSDT go well with medium thick and thick structures. However, such theories always ignore some stress variables or displacements and consequently lead to mathematical mistakes. ¹²

To get minimize potential problems, state space method (SSM) is utilized throughout the experiment. The advantages of SSM are that all the fundamental equations of three-dimensional (3D) elasticity are exactly satisfied and all nine elastic variables are taken into account. ¹³ Without any assumption of stress and displacement, SSM is efficient at dealing with the static and dynamic problems of FGM structures. Zeng et al. ¹⁴ investigated the natural frequencies of FG circular arches by employing the Fourier series expansion and state space formulation. Xu ¹⁵ studied the fundamental response of annular, circular and sectorial plates, and established new state space formulations by adding two stress functions and extra displacement functions. The numerical results were very parallel with those of FEM. Furthermore, to consider various boundary conditions, the differential quadrature method (DQM) is effectively used to discretize the governing equations. ¹⁶ Based on the 3D elasticity theory, Alibeigloo et al.^17,18 analysed free vibration of FG cylindrical plates and shells embedded in piezoelectric layers by using SSM and DQM. Nie and Zhong ¹⁹ proposed a semi-analytical method which integrated the SSM and DQM to analyse the free vibration of FG annular sectorial plates. Besides, many analytical methods are also used to solve the governing equations of FGM structures such as variational iteration method^20,21 and perturbation method.^22,23 Odibat and Momani ²⁴ compared these two methods when solving different types of differential equations of fractional order. Although the mentioned theories and methods can obtain good results, the above work focused on the statics or free vibration of composite structures.

Hence, some researchers have been investigating the dynamic vibration of FGM structures by other kinds of methods. Based on the hybrid numerical method (HNM) combining with a reduced-basis method (RBM), Huang and Huang ²⁵ investigated the real-time transient vibration of FGM plates. Zhou et al. ²⁶ studied the transient thermoelastic response of FG rectangular plates by utilising state space approach, numerical Laplace transformations and shooting methods. Nezhadi et al. ²⁷ studied the response of FG shells with impulse loads by integrating the Hamilton’s principle and Rayleigh–Ritz method. Frikha et al. ²⁸ analysed the dynamic behaviour of FG carbon nanotubes-reinforced composite shell structures. Based on the elasticity theory and Hamilton’s principle, the transient response in thermal environment of multi-layered FG shells was presented by Malekzadeh et al. ²⁹ Selahi et al. ³⁰ developed a hybrid method by using 3D elasticity theory for dynamic behaviour of FG truncated conical shell, with DQM discretizing the governing equations in both spatial and time domains. Liang et al.^{31t a} investigated the 3D transient response of FGM rectangular plates, annular sector plates and cylindrical shells under various boundary conditions. The responses were calculated by a proposed approach combining the SSM with DQM and Durbin’s ³⁴ numerical method.

As mentioned above, the dynamics of FGM structures have garnered a lot of attention from many researchers. But to the author’s best knowledge, the dynamic response of the FGM cylindrical panels using this 3D semi-analytical method have not been found yet, especially studies considering various boundary conditions. Besides, there is no reliable comparative verification of dynamics in other researchers’ work. On those circumstances, this paper gives trustworthier comparative solutions with nature frequencies in other literatures and transient behaviours by ANSYS. In the recent years, many analytical and semi-analytical methods, such as integral transform ³⁵ and unified Jacobi–Ritz method, ³⁶ have been applied for solving the vibration problems. The work in this paper can also serve as a benchmark of transient vibration analysis of FGM structures.

The paper is organised in the following manner. The backgrounds for the static and dynamic researches of FGM structures are provided. Then, the problem description is introduced by the three-dimensional linear elasticity theory. A concise introduction of the relevant fundamental methods is given next. Governing equations for four different boundary conditions are presented based on the SSM, and the semi-analytical solutions are derived based on the DQM and Durbin’s approach. In contrast to the finite element analysis (FEA), the proposed method is validated, and the influences of boundary conditions, geometric, material, and computational parameters are also studied. The research of carbon nanotube (CNT) reinforced FGM composite structures is investigated. The conclusions are then presented in the last section.

Problem description

A linear elastic FGM cylindrical panel is established, and its length l, central angle α, outer radius a and inner radius b are all considered. The model and the coordinate system are depicted in Figure 1. The mechanical material properties of the panel vary gradually along the radial r direction in an arbitrary manner. By extending the laminated plate model to FGM panel, the FGM panel consists of K-layers of graded materials along the radial directions.

Figure 1.

FGM cylindrical panel and coordinate system.

Based on the small deformation assumption, the second order of strain is negligible. Assuming the material of each layer is orthotropic among the coordinate planes, the linear stress–displacement relationships for an arbitrary layer can be expressed as the following matrix form

$$\left[\begin{array}{c}{\sigma }_r\\ {\sigma }_{\theta }\\ {\sigma }_z\\ {\tau }_{z\theta }\\ {\tau }_{\mathit{rz}}\\ {\tau }_{r\theta }\end{array}\right]=\left[\begin{array}{cccccc}{C}_{11}& C_{12}& C_{13}& 0& 0& 0\\ C_{12}& C_{22}& C_{23}& 0& 0& 0\\ C_{13}& C_{23}& C_{33}& 0& 0& 0\\ 0& 0& 0& C_{44}& 0& 0\\ 0& 0& 0& 0& C_{55}& 0\\ 0& 0& 0& 0& 0& C_{66}\end{array}\right]\left[\begin{array}{l}{\partial }_r{u}_r\hfill \\ (u_r+{\partial }_{\theta }{u}_{\theta })/r\hfill \\ {\partial }_z{u}_z\hfill \\ {\partial }_z{u}_{\theta }+{\partial }_{\theta }{u}_z/r\hfill \\ {\partial }_z{u}_r+{\partial }_r{u}_z\hfill \\ ({\partial }_{\theta }{u}_r-u_{\theta })/r+{\partial }_r{u}_{\theta }\hfill \end{array}\right]$$

(1)

where C₁₁, C₁₂, C₁₃, C₂₂, C₂₃, C₃₃, C₄₄, C₅₅, C₆₆ are the elastic stiffness coefficients; σ_r, σ_θ, σ_z are the radial stress components; τ_zθ, τ_rz, τ_rθ are the shear stress components; u_r, u_θ, u_z are the displacements components.

ρ is the mass density. With the absence of body force, the equilibrium equation is given as

$$\begin{array}{l}{\partial }_z{\sigma }_z+\frac{{\tau }_{\mathit{rz}}}{r}+{\partial }_r{\tau }_{\mathit{rz}}+\frac{{\partial }_{\theta }{\tau }_{z\theta }}{r}-\rho {\partial }_t^2{u}_z=0,\\ \frac{{\partial }_{\theta }{\sigma }_{\theta }}{r}+\frac{2{\tau }_{r\theta }}{r}+{\partial }_r{\tau }_{r\theta }+{\partial }_z{\tau }_{z\theta }-\rho {\partial }_t^2{u}_{\theta }=0,\\ {\partial }_r{\sigma }_r+\frac{{\sigma }_r-{\sigma }_{\theta }}{r}+{\partial }_z{\tau }_{\mathit{rz}}+\frac{{\partial }_{\theta }{\tau }_{r\theta }}{r}-\rho {\partial }_t^2{u}_r=0\end{array}$$

(2)

Four different boundaries are set in z = 0 or l. All of them can be written as

Clamped (z = 0)–Clamped (z = l)

$$\text{at}z=0,u_r=u_{\theta }=u_z=0$$

(3)

$$\text{at}z=l,u_r=u_{\theta }=u_z=0$$

(4)

Clamped (z = 0)–Simply supported (z = l)

$$\text{at}z=0,u_r=u_{\theta }=u_z=0$$

(5)

$$\text{at}z=l,{\sigma }_z=u_r=u_{\theta }=0$$

(6)

Clamped (z = 0)–Free (z = l)

$$\text{at}z=0,u_r=u_{\theta }=u_z=0$$

(7)

$$\text{at}z=l,{\sigma }_z={\tau }_{z\theta }={\tau }_{\mathit{rz}}=0$$

(8)

Simply supported (z = 0)–Simply supported (z = l)

$$\text{at}z=0,{\sigma }_z=u_r=u_{\theta }=0$$

(9)

$$\text{at}z=l,{\sigma }_z=u_r=u_{\theta }=0$$

(10)

The boundary conditions at inner (r = b) and outer (r = a) surfaces are given as

$$\text{at}\text{inner}\text{surface}\left(r=b\right),{\sigma }_r=f_{\mathit{rb}},{\tau }_{r\theta }=f_{\theta b},\text{ and }{\tau }_{\mathit{rz}}=f_{\mathit{zb}}$$

(11)

$$\text{at}\text{outer}\text{surface}\left(r=a\right),{\sigma }_r=f_{\mathit{ra}},{\tau }_{r\theta }=f_{\theta a},\text{ and }{\tau }_{\mathit{rz}}=f_{\mathit{za}}$$

(12)

Numerical method

An efficient numerical inversion for Laplace transform

The Laplace transform, which is very efficient in dealing with complex differential equations, is widely used in structural dynamics.^37,38 To seek a solution to the complicated dynamic problems, Durbin’s ³⁴ numerical inversion method was proposed which always generates almost the same results as that given by analytical inversion methods in a short time range. The numerical inversion method for Laplace transform can be expressed by the following formula. ³⁴

$$f(\stackrel{\rightharpoonup }{t})=\frac{\mathrm{2\hspace{0.17em}}\exp (\alpha \stackrel{\rightharpoonup }{t})}{T}\left\{\frac{\tilde{f}(\alpha )}{2}+{\displaystyle \sum _{k=1}^K\mathrm{Re}\left[\tilde{f}(\alpha +k\text{πi/}T)\right]\cos (k\text{π}\stackrel{\rightharpoonup }{t}\text{/}T)}\right\}$$

(13)

where α = 5/T, T = 2×T_d, K is a sufficient number, and T_d is the observation period.

Differential quadrature method

The function of DQM is to discretise the fundamental equations by turning the partial derivative to polynomials along the z direction. Based on the high-order polynomials, the weighted coefficients of DQM are computed by a series of mathematical formulation. ³⁹ Now, assuming $\tilde{\overline{f}}\left(\stackrel{\rightharpoonup }{r},\stackrel{\rightharpoonup }{z},\stackrel{\rightharpoonup }{s}\right)$ is a continuous function, the i-th order partial derivative along the z direction can be expressed by a sum of values in the spatial domain

$$\frac{{\partial }^i\tilde{\overline{f}}\left(\stackrel{\rightharpoonup }{r},{\stackrel{\rightharpoonup }{z}}_m,\stackrel{\rightharpoonup }{s}\right)}{\partial {\stackrel{\rightharpoonup }{z}}^i}={\displaystyle \sum _{n=1}^M{A}_{\mathit{mn}}^{(i)}\hspace{0.17em}\tilde{\overline{f}}\left(\stackrel{\rightharpoonup }{r},{\stackrel{\rightharpoonup }{z}}_n,\stackrel{\rightharpoonup }{s}\right)}$$

(14)

where m is from 1 to sampling points M; n is from 1 to M−1; $A_{\mathit{mn}}^{(i)}$ is the weighting coefficients.

Semi-analytical solution

According to the stress–displacement relationships and equilibrium differential equations, a semi-analytical method which requires employing the SSM and DQM is proposed to solve the governing equations under various boundary conditions. To obtain the natural frequencies and transient response of the cylindrical panel, the SSM and one-dimensional differential quadrature rule are used to establish a linear eigenvalue system, firstly along the radial direction and then established along the length direction.

Normalization

To simplify the calculation, the variables can be normalised in below

$$\left\{\begin{array}{c}\stackrel{\rightharpoonup }{C}ij\\ {\stackrel{\rightharpoonup }{\sigma }}_r\\ {\stackrel{\rightharpoonup }{\sigma }}_{\theta }\\ {\stackrel{\rightharpoonup }{\sigma }}_z\\ \stackrel{\rightharpoonup }{\tau }r\theta \\ \stackrel{\rightharpoonup }{\tau }z\theta \\ \stackrel{\rightharpoonup }{\tau }rz\end{array}\right\}=\frac{1}{C_{33}}\left\{\begin{array}{c}{\mathit{C}}_{\mathit{ij}}\\ {\sigma }_r\\ {\sigma }_{\theta }\\ {\sigma }_z\\ \tau r\theta \\ \tau z\theta \\ {\tau }_{\mathit{rz}}\end{array}\right\},\left\{\begin{array}{c}\stackrel{\rightharpoonup }{r}\\ \stackrel{\rightharpoonup }{z}\\ \stackrel{\rightharpoonup }{t}\\ \stackrel{\rightharpoonup }{s}\\ {\stackrel{\rightharpoonup }{u}}_r\\ {\stackrel{\rightharpoonup }{u}}_{\theta }\\ {\stackrel{\rightharpoonup }{u}}_z\end{array}\right\}={\left\{\begin{array}{ccccccc}\frac{1}{a}& \frac{1}{l}& \frac{c}{a}& \frac{a}{c}& \frac{1}{a}& \frac{1}{a}& \frac{1}{a}\end{array}\right\}}^T{\left\{\begin{array}{c}r\\ z\\ t\\ s\\ \mathit{ur}\\ u\theta \\ \mathit{uz}\end{array}\right\}}^T$$

(15)

where i, j = 1,2,…,6, and c is the wave velocity.

$$c={\left(C_{33}/\rho \right)}^{1/2}$$

(16)

To separate the variable θ, Fourier series are employed to expand the displacement and stress components into a form of trigonometric function. ³² The wave number is j, the panel’s sector angle is α.

$$\left\{\begin{array}{c}{\stackrel{\rightharpoonup }{\sigma }}_r(\stackrel{\rightharpoonup }{r},\theta ,\stackrel{\rightharpoonup }{z},t)\\ {\stackrel{\rightharpoonup }{u}}_r(\stackrel{\rightharpoonup }{r},\theta ,\stackrel{\rightharpoonup }{z},t)\\ {\stackrel{\rightharpoonup }{u}}_{\theta }(\stackrel{\rightharpoonup }{r},\theta ,\stackrel{\rightharpoonup }{z},t)\\ {\stackrel{\rightharpoonup }{u}}_z(\stackrel{\rightharpoonup }{r},\theta ,\stackrel{\rightharpoonup }{z},t)\\ {\stackrel{\rightharpoonup }{\tau }}_{\mathit{rz}}(\stackrel{\rightharpoonup }{r},\theta ,\stackrel{\rightharpoonup }{z},t)\\ {\stackrel{\rightharpoonup }{\tau }}_{r\theta }(\stackrel{\rightharpoonup }{r},\theta ,\stackrel{\rightharpoonup }{z},t)\\ {\stackrel{\rightharpoonup }{\sigma }}_{\theta }(\stackrel{\rightharpoonup }{r},\theta ,\stackrel{\rightharpoonup }{z},t)\\ {\stackrel{\rightharpoonup }{\sigma }}_z(\stackrel{\rightharpoonup }{r},\theta ,\stackrel{\rightharpoonup }{z},t)\\ {\stackrel{\rightharpoonup }{\tau }}_{z\theta }(\stackrel{\rightharpoonup }{r},\theta ,\stackrel{\rightharpoonup }{z},t)\end{array}\right\}={\displaystyle \sum _{j=0}^{\infty }\left\{\begin{array}{c}{\overline{\sigma }}_r(\stackrel{\rightharpoonup }{r},\stackrel{\rightharpoonup }{z},t)\sin \left(j\text{π}\theta /\alpha \right)\\ {\overline{u}}_r(\stackrel{\rightharpoonup }{r},\stackrel{\rightharpoonup }{z},t)\sin \left(j\text{π}\theta /\alpha \right)\\ {\overline{u}}_{\theta }(\stackrel{\rightharpoonup }{r},\stackrel{\rightharpoonup }{z},t)\cos \left(j\text{π}\theta /\alpha \right)\\ {\overline{u}}_z(\stackrel{\rightharpoonup }{r},\stackrel{\rightharpoonup }{z},t)\sin \left(j\text{π}\theta /\alpha \right)\\ {\overline{\tau }}_{\mathit{rz}}(\stackrel{\rightharpoonup }{r},\stackrel{\rightharpoonup }{z},t)\sin \left(j\text{π}\theta /\alpha \right)\\ {\overline{\tau }}_{r\theta }(\stackrel{\rightharpoonup }{r},\stackrel{\rightharpoonup }{z},t)\cos \left(j\text{π}\theta /\alpha \right)\\ {\overline{\sigma }}_{\theta }(\stackrel{\rightharpoonup }{r},\stackrel{\rightharpoonup }{z},t)\sin \left(j\text{π}\theta /\alpha \right)\\ {\overline{\sigma }}_z(\stackrel{\rightharpoonup }{r},\stackrel{\rightharpoonup }{z},t)\sin \left(j\text{π}\theta /\alpha \right)\\ {\overline{\tau }}_{z\theta }(\stackrel{\rightharpoonup }{r},\stackrel{\rightharpoonup }{z},t)\cos \left(j\text{π}\theta /\alpha \right)\end{array}\right\}}$$

(17)

By substituting equations (1) and (2), governing equation can be obtained in the form of state space approach. After employing the Fourier expansion and Laplace transform, the fundamental equations are derived as

$$\begin{array}{l}\frac{d{\tilde{\overline{\sigma }}}_r}{d\overrightarrow{r}}=\left(\overrightarrow{\rho }{\overrightarrow{s}}^2+\frac{{\eta }_2}{{\overrightarrow{r}}^2}\right)\tilde{\overline{u}}r+\frac{a{\eta }_3}{\mathit{lr}}\frac{\partial \tilde{\overline{u}}}{\partial \overrightarrow{z}}-\frac{j\pi {\eta }_2\tilde{\overline{u}}\theta }{\alpha {\overrightarrow{r}}^2}-\frac{{\eta }_1{\overline{\sigma }}_r}{\overrightarrow{r}}-\frac{a\partial \tilde{\overline{\tau }}rz}{l\partial \overrightarrow{z}}+\frac{j\pi \tilde{\overline{\tau }}r\theta }{\alpha \overrightarrow{r}}\\ \frac{d{\tilde{\overline{u}}}_r}{d\overrightarrow{r}}=\frac{{\tilde{\overline{\sigma }}}_r}{{\overrightarrow{C}}_{11}}-\frac{{\overrightarrow{C}}_{12}}{{\overrightarrow{C}}_{11}\overrightarrow{r}}{\tilde{\overline{u}}}_r+\frac{j\pi {\overrightarrow{C}}_{12}}{\alpha \overrightarrow{r}{\overrightarrow{C}}_{11}}{\tilde{\overline{u}}}_{\theta }-\frac{a{\overrightarrow{C}}_{13}}{l{\overrightarrow{C}}_{11}}\frac{\partial {\tilde{\overline{u}}}_z}{\partial \overrightarrow{z}}\end{array}$$

$$\begin{array}{l}\frac{d{\tilde{\overline{u}}}_{\theta }}{d\overrightarrow{r}}=\frac{-j\pi {\tilde{\overline{u}}}_r}{\alpha \overrightarrow{r}}+\frac{{\tilde{\overline{u}}}_{\theta }}{\overrightarrow{r}}+\frac{\tilde{\overline{\tau }}r\theta }{{\overrightarrow{C}}_{66}}\\ \frac{d{\tilde{\overline{u}}}_z}{d\overrightarrow{r}}=-\frac{a\partial {\tilde{\overline{u}}}_r}{l\partial \overrightarrow{z}}+\frac{\tilde{\overline{\tau }}rz}{{\overrightarrow{C}}_{55}}\\ \frac{d\tilde{\overline{\tau }}rz}{d\overrightarrow{r}}=-\frac{a{\eta }_3\partial {\tilde{\overline{u}}}_r}{l\overrightarrow{r}\partial \overrightarrow{z}}+\left(\overrightarrow{\rho }{\overrightarrow{s}}^2-\frac{a^2{\eta }_4{\partial }^2}{l^2\partial {\overrightarrow{z}}^2}+\frac{j^2{\pi }^2{\overrightarrow{C}}_{44}}{{\alpha }^2{\overrightarrow{r}}^2}\right){\tilde{\overline{u}}}_z+\frac{\mathit{aj}\pi {\eta }_5\partial {\tilde{\overline{u}}}_{\theta }}{l\alpha \overrightarrow{r}}-\frac{a{\overrightarrow{C}}_{13}\partial {\tilde{\overline{\sigma }}}_r}{l{\overrightarrow{C}}_{11}\partial \overrightarrow{z}}-\frac{\tilde{\overline{\tau }}rz}{\overrightarrow{r}}\\ \frac{d\tilde{\overline{\tau }}r\theta }{d\overrightarrow{r}}=-\frac{j\pi {\eta }_2{\tilde{\overline{u}}}_r}{\alpha {\overrightarrow{r}}^2}+\frac{\mathit{aj}\pi {\eta }_5\partial {\tilde{\overline{u}}}_z}{l\alpha \overrightarrow{r}\partial \overrightarrow{z}}+\left(\overrightarrow{\rho }{\overrightarrow{s}}^2-\frac{a^2{\overrightarrow{C}}_{44}{\partial }^2}{l^2\partial {\overrightarrow{z}}^2}+\frac{j^2{\pi }^2{\overrightarrow{\eta }}_2}{{\alpha }^2{\overrightarrow{r}}^2}\right){\tilde{\overline{u}}}_{\theta }-\frac{j\pi {\overrightarrow{C}}_{12}{\tilde{\overline{\sigma }}}_r}{\alpha \overrightarrow{r}{\overrightarrow{C}}_{11}}-\frac{2\tilde{\overline{\tau }}r\theta }{\overrightarrow{r}}\end{array}$$

(18)

where

$$\begin{array}{l}{\eta }_1\text{=1}-\frac{{\stackrel{\rightharpoonup }{C}}_{12}}{{\stackrel{\rightharpoonup }{C}}_{11}}\text{,}{\eta }_2={\stackrel{\rightharpoonup }{C}}_{22}-\frac{{\stackrel{\rightharpoonup }{C}}_{12}^2}{{\stackrel{\rightharpoonup }{C}}_{11}},{\eta }_3={\stackrel{\rightharpoonup }{C}}_{23}-\frac{{\stackrel{\rightharpoonup }{C}}_{12}{\stackrel{\rightharpoonup }{C}}_{13}}{{\stackrel{\rightharpoonup }{C}}_{11}},{\eta }_4={\stackrel{\rightharpoonup }{C}}_{33}-\frac{{\stackrel{\rightharpoonup }{C}}_{13}^2}{{\stackrel{\rightharpoonup }{C}}_{11}},\\ {\eta }_5=\left[-{\stackrel{\rightharpoonup }{C}}_{12}{\stackrel{\rightharpoonup }{C}}_{13}+{\stackrel{\rightharpoonup }{C}}_{11}\left({\stackrel{\rightharpoonup }{C}}_{23}+{\stackrel{\rightharpoonup }{C}}_{44}\right)\right]/{\stackrel{\rightharpoonup }{C}}_{11}\end{array}$$

(19)

The other three variables can be determined by

$$\begin{array}{l}{\tilde{\overline{\sigma }}}_{\theta }=\frac{{\eta }_2}{\stackrel{\rightharpoonup }{r}}{\tilde{\overline{u}}}_r+\frac{a{\eta }_3}{l}\frac{\partial {\tilde{\overline{u}}}_z}{\partial \stackrel{\rightharpoonup }{z}}-\frac{j\pi {\eta }_2}{\alpha \stackrel{\rightharpoonup }{r}}{\tilde{\overline{u}}}_{\theta }+\frac{{\stackrel{\rightharpoonup }{C}}_{12}}{{\stackrel{\rightharpoonup }{C}}_{11}}{\tilde{\overline{\sigma }}}_r,\\ {\tilde{\overline{\sigma }}}_z=\frac{{\eta }_3}{\stackrel{\rightharpoonup }{r}}{\tilde{\overline{u}}}_r+\frac{a{\eta }_4}{l}\frac{\partial {\tilde{\overline{u}}}_z}{\partial \stackrel{\rightharpoonup }{z}}-\frac{j\pi {\eta }_3}{\alpha \stackrel{\rightharpoonup }{r}}{\tilde{\overline{u}}}_{\theta }+\frac{{\stackrel{\rightharpoonup }{C}}_{13}}{{\stackrel{\rightharpoonup }{C}}_{11}}{\tilde{\overline{\sigma }}}_r,\\ {\tilde{\overline{\tau }}}_{z\theta }=\frac{j\pi {\stackrel{\rightharpoonup }{C}}_{44}}{\alpha \stackrel{\rightharpoonup }{r}}{\tilde{\overline{u}}}_z+\frac{a{\stackrel{\rightharpoonup }{C}}_{44}}{l}\frac{\partial {\tilde{\overline{u}}}_{\theta }}{\partial \stackrel{\rightharpoonup }{z}}\end{array}$$

(20)

As M is the sum of sample points, m-th is the location of calculating the sample point. By applying DQM on equations (18) and (20), the new state space equations can be rewritten as shown below

$$\begin{array}{l}\frac{\text{d}{\tilde{\overline{\sigma }}}_{\mathit{rm}}}{\text{d}\stackrel{\rightharpoonup }{r}}=-\frac{a}{l}{\displaystyle \sum }_{n=1}^M{A}_{\mathit{mn}}^{\left(1\right)}{\tilde{\overline{\tau }}}_{\mathit{rzn}}+\frac{a{\eta }_3}{l\stackrel{\rightharpoonup }{r}}{\displaystyle \sum }_{n=1}^M{A}_{\mathit{mn}}^{\left(1\right)}{\tilde{\overline{u}}}_{\mathit{zn}}+\left({\stackrel{\rightharpoonup }{s}}^2\stackrel{\rightharpoonup }{\rho }+\frac{{\eta }_2}{{\stackrel{\rightharpoonup }{r}}^2}\right){\tilde{\overline{u}}}_{\mathit{rm}}\\ -\frac{j\text{π}{\eta }_2{\tilde{\overline{u}}}_{\theta m}}{\alpha {\stackrel{\rightharpoonup }{r}}^2}-\frac{{\eta }_1{\tilde{\overline{\sigma }}}_{\mathit{rm}}}{\stackrel{\rightharpoonup }{r}}+\frac{j\text{π}{\tilde{\overline{\tau }}}_{r\theta m}}{\alpha \stackrel{\rightharpoonup }{r}},\\ \frac{\text{d}{\tilde{\overline{u}}}_{\mathit{rm}}}{\text{d}\stackrel{\rightharpoonup }{r}}=\frac{{\tilde{\overline{\sigma }}}_{\mathit{rm}}}{{\stackrel{\rightharpoonup }{C}}_{11}}-\frac{{\stackrel{\rightharpoonup }{C}}_{12}}{\stackrel{\rightharpoonup }{r}{\stackrel{\rightharpoonup }{C}}_{11}}{\tilde{\overline{u}}}_{\mathit{rm}}+\frac{j\text{π}{\stackrel{\rightharpoonup }{C}}_{12}}{\alpha \stackrel{\rightharpoonup }{r}{\stackrel{\rightharpoonup }{C}}_{11}}{\tilde{\overline{u}}}_{\theta m}-\frac{a{\stackrel{\rightharpoonup }{C}}_{13}}{l{\stackrel{\rightharpoonup }{C}}_{11}}{\displaystyle \sum }_{n=1}^M{A}_{\mathit{mn}}^{\left(1\right)}{\tilde{\overline{u}}}_{\mathit{zn}},\\ \frac{\text{d}{\tilde{\overline{u}}}_{\theta m}}{\text{d}\stackrel{\rightharpoonup }{r}}=-\frac{j\text{π}{\tilde{\overline{u}}}_{\mathit{rm}}}{\alpha \stackrel{\rightharpoonup }{r}}+\frac{{\tilde{\overline{u}}}_{\theta m}}{\stackrel{\rightharpoonup }{r}}+\frac{{\tilde{\overline{\tau }}}_{r\theta m}}{{\stackrel{\rightharpoonup }{C}}_{66}},\end{array}$$

$$\begin{array}{l}\frac{\text{d}{\tilde{\overline{u}}}_{\mathit{zm}}}{\text{d}\stackrel{\rightharpoonup }{r}}=-\frac{a}{l}{\displaystyle \sum }_{n=1}^M{A}_{\mathit{mn}}^{\left(1\right)}{\tilde{\overline{u}}}_{\mathit{rn}}+\frac{{\tilde{\overline{\tau }}}_{\mathit{rzm}}}{{\stackrel{\rightharpoonup }{C}}_{55}},\\ \frac{\text{d}{\tilde{\overline{\tau }}}_{\mathit{rzm}}}{\text{d}\stackrel{\rightharpoonup }{r}}=-\frac{a{\eta }_3}{l\stackrel{\rightharpoonup }{r}}{\displaystyle \sum }_{n=1}^M{A}_{\mathit{mn}}^{\left(1\right)}{\tilde{\overline{u}}}_{\mathit{rn}}-\frac{a^2{\eta }_4}{l^2}{\displaystyle \sum }_{n=1}^M{A}_{\mathit{mn}}^{\left(2\right)}{\tilde{\overline{u}}}_{\mathit{zn}}+\frac{\mathit{aj}\text{π}{\eta }_5}{l\alpha \stackrel{\rightharpoonup }{r}}{\displaystyle \sum }_{n=1}^M{A}_{\mathit{mn}}^{\left(1\right)}{\tilde{\overline{u}}}_{\theta n}\\ -\frac{{\tilde{\overline{\tau }}}_{\mathit{rzm}}}{\stackrel{\rightharpoonup }{r}}-\frac{a{\stackrel{\rightharpoonup }{C}}_{13}}{l{\stackrel{\rightharpoonup }{C}}_{11}}{\displaystyle \sum }_{n=1}^M{A}_{\mathit{mn}}^{\left(1\right)}{\tilde{\overline{\sigma }}}_{\mathit{rn}}+\left({\stackrel{\rightharpoonup }{s}}^2\stackrel{\rightharpoonup }{\rho }+\frac{j^2{\text{π}}^2{\stackrel{\rightharpoonup }{C}}_{44}}{{\alpha }^2{\stackrel{\rightharpoonup }{r}}^2}\right){\tilde{\overline{u}}}_{\mathit{zm}},\\ \frac{\text{d}{\tilde{\overline{\tau }}}_{r\theta m}}{\text{d}\stackrel{\rightharpoonup }{r}}=\frac{\mathit{aj}\text{π}{\eta }_5}{l\alpha \stackrel{\rightharpoonup }{r}}{\displaystyle \sum }_{n=1}^M{A}_{\mathit{mn}}^{\left(1\right)}{\tilde{\overline{u}}}_{\mathit{zn}}-\frac{j\text{π}{\eta }_2{\tilde{\overline{u}}}_{\mathit{rm}}}{\alpha {\stackrel{\rightharpoonup }{r}}^2}+\left({\stackrel{\rightharpoonup }{s}}^2\stackrel{\rightharpoonup }{\rho }+\frac{j^2{\text{π}}^2{\eta }_2}{{\alpha }^2{\stackrel{\rightharpoonup }{r}}^2}\right){\tilde{\overline{u}}}_{\theta m}\\ -\frac{2{\tilde{\overline{\tau }}}_{r\theta m}}{\stackrel{\rightharpoonup }{r}}-\frac{j\text{π}{\stackrel{\rightharpoonup }{C}}_{12}}{\alpha \stackrel{\rightharpoonup }{r}{\stackrel{\rightharpoonup }{C}}_{11}}{\tilde{\overline{\sigma }}}_{\mathit{rm}}-\frac{a^2{\stackrel{\rightharpoonup }{C}}_{44}}{l^2}{\displaystyle \sum }_{n=1}^M{A}_{\mathit{mn}}^{\left(2\right)}{\tilde{\overline{u}}}_{\theta n}\end{array}$$

(21)

For the m-th sample point, the other three variables can be determined by

$$\begin{array}{c}{\tilde{\overline{\sigma }}}_{\theta m}=\frac{a{\eta }_3}{l}{\displaystyle \sum }_{n=1}^M{A}_{\mathit{mn}}^{\left(1\right)}{\tilde{\overline{u}}}_{\mathit{zn}}+\frac{{\eta }_2{\tilde{\overline{u}}}_{\mathit{rm}}}{\stackrel{\rightharpoonup }{r}}-\frac{j\text{π}{\eta }_2{\tilde{\overline{u}}}_{\theta m}}{\alpha \stackrel{\rightharpoonup }{r}}+\frac{{\stackrel{\rightharpoonup }{C}}_{12}}{{\stackrel{\rightharpoonup }{C}}_{11}}{\tilde{\overline{\sigma }}}_{\mathit{rm}},\\ {\tilde{\overline{\sigma }}}_{\mathit{zm}}=\frac{a{\eta }_4}{l}{\displaystyle \sum }_{n=1}^M{A}_{\mathit{mn}}^{\left(1\right)}{\tilde{\overline{u}}}_{\mathit{zn}}+\frac{{\eta }_3{\tilde{\overline{u}}}_{\mathit{rm}}}{\stackrel{\rightharpoonup }{r}}-\frac{j\text{π}{\eta }_3{\tilde{\overline{u}}}_{\theta m}}{\alpha \stackrel{\rightharpoonup }{r}}+\frac{{\stackrel{\rightharpoonup }{C}}_{13}}{{\stackrel{\rightharpoonup }{C}}_{11}}{\tilde{\overline{\sigma }}}_{\mathit{rm}},\\ {\tilde{\overline{\tau }}}_{z\theta m}=\frac{a{\stackrel{\rightharpoonup }{C}}_{44}}{l}{\displaystyle \sum }_{n=1}^M{A}_{\mathit{mn}}^{\left(1\right)}{\tilde{\overline{u}}}_{\theta n}+\frac{j\text{π}{\stackrel{\rightharpoonup }{C}}_{44}}{\alpha \stackrel{\rightharpoonup }{r}}{\tilde{\overline{u}}}_{\mathit{zm}}\end{array}$$

(22)

To apply the DQM, Laplace transform and Fourier series on equations (3) to (10), the cylindrical panel’s boundary conditions are given as

C-C

$$z=0,{\tilde{\overline{u}}}_{r1}={\tilde{\overline{u}}}_{\theta 1}={\tilde{\overline{u}}}_{z1}=0$$

(23)

$$z=l,{\tilde{\overline{u}}}_{\mathit{rM}}={\tilde{\overline{u}}}_{\theta M}={\tilde{\overline{u}}}_{\mathit{zM}}=0$$

(24)

C-S

$$z=0,{\tilde{\overline{u}}}_{r1}={\tilde{\overline{u}}}_{\theta 1}={\tilde{\overline{u}}}_{z1}=0$$

(25)

$$z=l,{\tilde{\overline{\sigma }}}_{\mathit{zM}}={\tilde{\overline{u}}}_{\theta M}={\tilde{\overline{u}}}_{\mathit{rM}}=0$$

(26)

C-F

$$z=0,{\tilde{\overline{u}}}_{r1}={\tilde{\overline{u}}}_{\theta 1}={\tilde{\overline{u}}}_{z1}=0$$

(27)

$$z=l,{\tilde{\overline{\sigma }}}_{\mathit{zM}}={\tilde{\overline{\tau }}}_{z\theta M}={\tilde{\overline{\tau }}}_{\mathit{rzM}}=0$$

(28)

S-S

$$z=0,{\tilde{\overline{\sigma }}}_{z1}={\tilde{\overline{u}}}_{\theta 1}={\tilde{\overline{u}}}_{r1}=0$$

(29)

$$z=l,{\tilde{\overline{\sigma }}}_{\mathit{zM}}={\tilde{\overline{u}}}_{\theta M}={\tilde{\overline{u}}}_{\mathit{rM}}=0$$

(30)

Similarly, for the inner and outer surface, the transformed boundaries are derived as

$$r=b,{\widetilde{\overline{\mathbf{\sigma }}}}_{\mathit{rb}}={\widetilde{\overline{\mathbf{f}}}}_{\mathit{rb}},{\widetilde{\overline{\mathbf{\tau }}}}_{r\theta b}={\widetilde{\overline{\mathbf{f}}}}_{\theta b},\text{\hspace{1em}and\hspace{1em}}{\widetilde{\overline{\mathbf{\tau }}}}_{\mathit{rzb}}={\widetilde{\overline{\mathbf{f}}}}_{\mathit{zb}}$$

(31)

$$r=a,{\widetilde{\overline{\mathbf{\sigma }}}}_{\mathit{ra}}={\widetilde{\overline{\mathbf{f}}}}_{\mathit{ra}},{\widetilde{\overline{\mathbf{\tau }}}}_{r\theta a}={\widetilde{\overline{\mathbf{f}}}}_{\theta a},\text{\hspace{1em}and\hspace{1em}}{\widetilde{\overline{\mathbf{\tau }}}}_{\mathit{rza}}={\widetilde{\overline{\mathbf{f}}}}_{\mathit{za}}$$

(32)

where ${\widetilde{\overline{\mathbf{\sigma }}}}_{\mathit{rb}}={\{{\tilde{\overline{\sigma }}}_{r1},\dots {\tilde{\overline{\sigma }}}_{\mathit{rM}}\}}_{r=b}{\hspace{0.17em}}^{\text{T}}$ , ${\widetilde{\overline{\mathbf{f}}}}_{\theta a}={\{{\tilde{\overline{f}}}_{\theta 1},\dots {\tilde{\overline{f}}}_{\theta M}\}}_{r=a}{\hspace{0.17em}}^{\text{T}}$ ,…

Semi-analytical solution

As shown in equation (21), the differential equation with so many variables is too complex to be solved by analytical approaches. Hence, for the sake of simplification, an approximate method is performed by introducing a radial local coordinate ξ_k, an approximation ξ_k/R_k≪1 is considered as a thin cylindrical panel in each layer. ξ_k is located at the centre of k-th layer. ⁴⁰

$$\frac{1}{\stackrel{\rightharpoonup }{r}}=\frac{\left(1-{\lambda }_k\right)}{R_k}\text{,\hspace{1em}}\frac{1}{\stackrel{\rightharpoonup }{r}}=\frac{\left(1-2{\lambda }_k\right)}{R_k^2}$$

(33)

where ξ_k = r−R_k and λ_k = ξ_k/R_k. According to the assumption, λ_k can be neglected. Hence, equation (21) is supposed to be derived by the form of matrix

$$\frac{\text{d}\tilde{\overline{Q}}(\lambda )}{\text{d}\lambda }=\mathbf{H}\cdot \tilde{\overline{Q}}(\lambda )$$

(34)

where $\tilde{\overline{Q}}=\{{\widetilde{\overline{\mathbf{\sigma }}}}_r,\hspace{0.17em}{\widetilde{\overline{\mathbf{u}}}}_r,\hspace{0.17em}{\widetilde{\overline{\mathbf{u}}}}_{\theta },\hspace{0.17em}{\widetilde{\overline{\mathbf{u}}}}_z,\hspace{0.17em}{\widetilde{\overline{\mathbf{\tau }}}}_{\mathit{rz}},\hspace{0.17em}{\widetilde{\overline{\mathbf{\tau }}}}_{r\theta }\}$ , ${\widetilde{\overline{\mathbf{\sigma }}}}_r\hspace{0.17em}=\hspace{0.17em}{\{{\tilde{\overline{\sigma }}}_{r1},\dots {\tilde{\overline{\sigma }}}_{\mathit{rM}}\}}^{\text{T}}$ , ${\widetilde{\overline{\mathbf{u}}}}_r\hspace{0.17em}=\hspace{0.17em}{\{{\tilde{\overline{u}}}_{r1},\dots {\tilde{\overline{u}}}_{\mathit{rM}}\}}^{\text{T}}$ ,… The constitution of matrix H can be found in Appendix 1. After applying different boundary conditions to the equations (23) to (30), the various H is given in Appendix 2.

The expression of equation (34) is supposed to be

$$\tilde{\overline{Q}}(\lambda )=\exp \left[{\mathbf{H}}_k\left(\lambda -{\lambda }_{k-1}\right)\right]\cdot \tilde{\overline{Q}}({\lambda }_{k-1}).\left(\frac{-h}{2R_k}\le \lambda \le \frac{h}{2R_k}\right)$$

(35)

Equation (35) at $\lambda ={\lambda }_k$ yields

$$\tilde{\overline{Q}}({\lambda }_k)=\exp \left(\mathbf{H}{h}_k\right)\cdot \tilde{\overline{Q}}({\lambda }_{k-1})$$

(36)

where h_k is the thickness of k-th layer.

Subsequently

$$\tilde{\overline{Q}}({\lambda }_{k+1})=\exp \left({\mathbf{H}}_{k+1}{h}_{k+1}\right)\cdot \tilde{\overline{Q}}({\lambda }_k)=\exp \left({\mathbf{H}}_{k+1}{h}_{k+1}\right)\cdot \exp \left({\mathbf{H}}_k{h}_k\right)\cdot \tilde{\overline{Q}}({\lambda }_{k-1})$$

(37)

For all K layers have the same behaviour, the FGM cylindrical panel’s state vectors at the outer surface can be determined by the inner surface

$$\tilde{\overline{Q}}\left(b\right)=\mathbf{T}\left(a-b\right)\cdot \tilde{\overline{Q}}\left(a\right)$$

(38)

where

$$\mathbf{T}\left(h\right)={\displaystyle \prod _{k=K}^1\exp \left({\mathbf{H}}_k{h}_k\right)}$$

(39)

$$\mathbf{T}\left(h\right)\text{=}\underset{6\times 6}{\mathit{Tij}}(i,j=1,2,\dots ,6)$$

(40)

At the outer and inner surfaces, the boundaries equations (31) and (32) are applied into equation (38)

$${\left\{\begin{array}{c}{\widetilde{\overline{\mathbf{f}}}}_r\\ {\widetilde{\overline{\mathbf{f}}}}_z\\ {\widetilde{\overline{\mathbf{f}}}}_{\theta }\end{array}\right\}}_{r=b}=\left[\begin{array}{cc}\begin{array}{ccc}{\mathbf{T}}_{11}& {\mathbf{T}}_{12}& {\mathbf{T}}_{13}\\ {\mathbf{T}}_{51}& {\mathbf{T}}_{52}& {\mathbf{T}}_{53}\\ {\mathbf{T}}_{61}& {\mathbf{T}}_{62}& {\mathbf{T}}_{63}\end{array}& \begin{array}{ccc}{\mathbf{T}}_{14}& {\mathbf{T}}_{15}& {\mathbf{T}}_{16}\\ {\mathbf{T}}_{54}& {\mathbf{T}}_{55}& {\mathbf{T}}_{56}\\ {\mathbf{T}}_{64}& {\mathbf{T}}_{65}& {\mathbf{T}}_{66}\end{array}\end{array}\right]{\left\{\begin{array}{c}{\widetilde{\overline{\mathbf{f}}}}_r\\ {\widetilde{\overline{\mathbf{u}}}}_r\\ {\widetilde{\overline{\mathbf{u}}}}_{\theta }\\ {\widetilde{\overline{\mathbf{u}}}}_z\\ {\widetilde{\overline{\mathbf{f}}}}_z\\ {\widetilde{\overline{\mathbf{f}}}}_{\theta }\end{array}\right\}}_{r=a}$$

(41)

where

$${\widetilde{\overline{\mathbf{f}}}}_r={\left\{{\tilde{\overline{f}}}_{r1},\dots ,{\tilde{\overline{f}}}_{\mathit{rM}}\right\}}^{\text{T}},\hspace{1em}{\widetilde{\overline{\mathbf{f}}}}_z={\left\{{\tilde{\overline{f}}}_{z1},\dots ,{\tilde{\overline{f}}}_{\mathit{zM}}\right\}}^{\text{T}},\hspace{1em}{\widetilde{\overline{\mathbf{f}}}}_{\theta }={\left\{{\tilde{\overline{f}}}_{\theta 1},\dots ,{\tilde{\overline{f}}}_{\theta M}\right\}}^{\text{T}}$$

Numerical results

Firstly, compared with data obtained by other theoretical methods, the natural frequencies of a clamped FGM cylindrical panel are investigated for different functionally graded index γ. Subsequently, the displacements along r direction obtained by the proposed method are compared with FEM to validate the proposed method. And then, the research regarding convergence are carried out for different cases. In conclusion, this allows the effects of cylindrical panel’s geometric size and functionally graded index γ to be investigated.

Nature frequencies study

In comparison with previous results in literature,^41,42 the natural frequencies of the FGM cylindrical panel are investigated (Figure 2). The geometric size of the cylindrical panel is L/h = 10, L/b = 0.1 with clamped boundary condition varies in the functionally graded index γ. And x is its chord length. The material properties of this FGM panel are stainless steel (SUS304) and silicon nitride (Si₃N₄). The constituents of the two kind of material are described in Table 1.

Figure 2.

Geometry of FGM cylindrical panel.

Table 1.

FGM cylindrical panels’ parameters. ⁴¹

Constituent	Elastic modulusE (GPa)	Poisson’s ratioμ	Mass densityρ (kg/m3)
SUS304 (metal)	207.7877	0.31776	8166
Si3N4 (ceramic)	322.2715	0.24	2370

In Table 2, non-dimensional frequency parameter is given by ${\omega }^{*}=\omega l^2\sqrt{{\rho }_{\mathit{m}}\mathit{h}/D_m^{*}}$ in which $D_m^{*}={\mathit{E}}_{\mathit{m}}{h}^3/12(1-{\mu }_m^2)$ .

Table 2.

Natural frequencies of the first mode for FGM cylindrical panels.

γ	Yang and Shen 42	Pradyumna and Bandyopadhyay 41	Proposed method
0	74.518	72.9613	73.8433
0.2	57.479	60.0269	59.6093
2	40.750	39.1457	39.0261

The purpose of free vibration analysis is to validate the reliability of the proposed approach. As we know, the boundary conditions, wave number and geometric size are the main factors to decide the values of natural frequency. In this section, the boundary condition is simply supported and wave number i is equal to 1. The results of the first mode frequencies are shown in Table 2 for different functionally graded index γ along with those of Yang and Shen ⁴² and Pradyumna and Bandyopadhyay. ⁴¹ It is indicated that the values obtained by the proposed method agree with Yang and Shen’s ⁴² and Pradyumna and Bandyopadhyay’s ⁴¹ methods.

Transient responses study

In order to further demonstrate the accuracy of the proposed approach, its displacements along the r direction are compared with the results by commercial software ANSYS. The cases with various variables are given in Table 3.

Table 3.

Cases with different variables and boundaries.

Case number	Variation law of material properties	central angle α	FG index γ	Load parameter j	Lengthl (m)	OuterRadiusa (m)	Widthb (m)	Boundary conditionsat z = 0, l
1	Exponential	30°	5	3	1	1	0.98	C-C
2	Exponential	45°	0.2	1	2	1	0.98	C-F
3	Power	90°	0.5	3	0.7	1	0.98	C-S
4	Power	60°	5	1	1.2	1	0.98	S-S

For the orthotropic properties, the stiffness coefficient and mass density are supposed to change following the variation law along the r direction. By applying exponential law, the material properties can be expressed as

$$\{C_{\mathit{ij}},\rho \}=\{{\widehat{C}}_{\mathit{ij}},\widehat{\rho }\}\times \exp \left[\gamma \left(a-r\right)/\left(a-b\right)\right]$$

(42)

In our studies and experiment, the material properties of this FGM panel are made up with zirconia and aluminium. Moreover, the distribution of volume fractions along the radial direction obeyed the power laws can be defined as

$$\begin{array}{c}{V}_C\left(r\right)={\left[\left(a-r\right)/\left(a-b\right)\right]}^{\gamma },\\ V_A\left(r\right)=1-{\left[\left(a-r\right)/\left(a-b\right)\right]}^{\gamma }\end{array}$$

(43)

Here, A means aluminium and C represents ceramic zirconia. Therefore, the material properties of this FGM cylindrical panel are derived as

$$\{C_{\mathit{ij}},\rho \}=\{C_{\mathit{ij},C},{\rho }_C\}{V}_C\left(r\right)+\{C_{\mathit{ij},A},{\rho }_A\}{V}_A\left(r\right)$$

(44)

The fundamental material properties in this paper are shown in Table 4.

Table 4.

Fundamental material properties. ⁴³

Constituent	Elastic modulusE (GPa)	Poisson’s ratioμ	Mass densityρ (kg/m3)
Aluminium (Al)	70	0.3	2700
Zirconia (ZrO2)	200	0.3	5700

Then the distribution properties of each layer can be defined as

$${\{C_{\mathit{ij}},\rho \}}_k=\frac{{\displaystyle {\int }_{z_{k-1}}^{z_k}\{C_{\mathit{ij}},\rho \}}\text{(}r\text{)d}r}{r_k-r_{k-1}}$$

(45)

where k = 1,2,…,K.

The force applying to the outer surface (r=a) of the cylindrical panel is given as ³¹

$$\mathrm{at}r=a,{\sigma }_r={10}^{-3}{C}_{33,A}\sin \left(\frac{j\pi }{\alpha }\theta \right)\exp \left(-t{c}_A/l\right)$$

(46)

where $c_A={\left(C_{33,A}/{\rho }_A\right)}^{1/2}$ is the longitudinal wave velocity of aluminium. ⁴⁴

At first, the layer number K and sampling point M are 4 and 21, respectively. The finite element models established by software ANSYS are given in Figure 3. Solid 64 of the structure elements are used in the finite element analysis. The central position (r=(a + b)/2, θ = 0, z = l/2) is chosen to calculate the results. In Figure 4, the numerical results between the proposed approach and FEM are compared with each other. It is shown that the results predicted by the two methods agree with each other. Compared with ANSYS, the method has a higher computational efficiency in calculating the transient response of the cylindrical panel, and relevant data can be found in Table 5. Moreover, the numerical solutions are independent from conditions determined by geometry, FG index and other boundary conditions.

Figure 3.

FE models of (a) case 1, (b) case 2, (c) case 3, and (d) case 4.

Figure 4.

Deflection history: (a) case 1; (b) case 2; (c) case 3; (d) case 4.

Table 5.

Time required for calculating the dynamic response of the FGM cylindrical panel.

Methods	NSP	1	2	3	Average time (s)
Proposed method	5	68.234	69.375	69.156	68.922
	11	162.154	139.295	179.751	160.4
	21	402.031	423.969	475.75	433.917
	31	879.189	967.252	974.125	940.189
ANSYS	/	1323.000	1293.688	1295.516	1304.068

NSP: number of sampling points.

Convergence studies

The importance of convergence studies lies in their purpose of clarifying the effects of different factors. And the condition which is used to carry out the research is the same as the conditions put in place for case 4 in Table 3.

Firstly, a series of sampling points are used to analyse the convergence rate along the z- direction. The values of M are 5, 11, 17, 21 and 31, respectively. The relation between the normalised deflection u_r/a and time history at central position are shown in Figure 5. It is indicated that, the more sampling point numbers, the more accurate the results will be. When the numbers of sampling points are at a relatively large level, the results do not change much. It also displays that the proposed method converges fast with the increasing sampling points.

Figure 5.

Deflection history with different numbers of sampling points.

Secondly, a series of layer numbers are employed in order to analyse the convergence rate along the r direction. The values of K are 2, 4, 8 and 12, respectively. The relation between the normalised deflection u_r/a and time history at central position are shown in Figure 6. It is indicated that the cases with K are larger than 4 and have nearly the same results. Therefore, the proposed approach has a high convergence rate with the change of layer numbers.

Figure 6.

Deflection history with different layer numbers.

Effects of l/a and (a−b)/a

To study the effects of l/a and (a-b)/a, the deflection along the radial direction of FGM cylindrical panels are examined here. The sampling point number and the layer number are 21 and 4 respectively. Also, the condition which is used to carry out the researches is the same as case 4 in Table 3.

Firstly, the effect of the length/outer radius ratio (l/a) is investigated. A series of l/a = 1, 1.2, 1.6, and 2 are used here. The relation between the normalised deflection u_r/a and time history at the central position are shown in Figure 7. It shows that the deflection of FGM cylindrical panels increases as length/outer radius ratio l/a increases. The reason is that the flexibility of the whole structure becomes greater when the ratio increases. Thus, the cylindrical panel deforms easier under the same external force, and the vibration period will also be longer.

Figure 7.

Deflection history with different length/outer radius ratio.

Secondly, the effect of the thickness/outer radius ratio (a-b)/a is studied. A series of (a-b)/a = 0.02, 0.04, 0.08 and 0.16 is employed here. Both the length and outer radius are fixed to be 1 m. The relation between the normalised deflection u_r/a and time history at central position are plotted in Figure 8. As indicated in the figure, the deflection shows a trend that is actually opposite to what is observed with the thickness/outer radius ratio. When the ratio increases, the deflection actually decreases. Because, on the one hand, the structural flexibility will become less with the increase of thickness, and on the other hand, the bearing capacity of the structure is strengthened with the increase of thickness. So that it can resist external force more effectively, and the vibration period will also be shortened.

Figure 8.

Deflection history with different thickness/outer radius ratio.

Effect of FG index

To study the effects of FG index γ, the deflection along the radial direction of FGM cylindrical panel is examined here. The sampling point number and the layer number are 21 and 4, respectively. Furthermore, the condition which is used to carry out the research is the same as case 4 in Table 3.

In this section, the paper has studied two different variation laws. On the one hand, the power variation law means that the fundamental properties distribute in a power law. On the other hand, the exponential variation law assumes the fundamental properties distribute in an exponential law. And both of them employed various FG parameters γ = 0.2, 0.5, 2 and 5 in this FGM cylindrical panel.

Two different variation laws show distinct differences on the transient response of FGM cylindrical panels at central position in Figure 9. In terms of exponential variation law, the deflection of the panel decreases as functionally graded index γ increases. But for the FGM panel with a power variation law, the deflection of the panel increases as functionally graded index γ increases. It is because the smaller value of FG index is, the more uniform the volume fraction becomes, which leading to the higher structural strength and weaker vibration response. However, the result of the exponential law is on the contrary.

Figure 9.

Deflection history with different variation law: (a) power law; (b) exponential law.

CNT-reinforced FGM

The discovery of CNT has been a great breakthrough in many applications.^45,46 Due to their superior mechanical, electrical and chemical properties, CNT can be utilized to strengthen composite materials. In this section, the uniform distribution of CNT is taken into consideration to improve the ability of resisting the FGM cylindrical panel’s vibration. The condition which is used to carry out the research is the same as case 3 in Table 3, but the time domain is different. The single-wall CNT (SWCNT) material properties are assumed to be ρ_CNT = 1400 kg/m³, μ_CNT = 0.175, E_CNT = 5.6466 TPa. The volume fraction of CNT for the FGM cylindrical panel is assumed to be V_CNT. The SWCNT uniformly distributed along the thickness direction of the composite cylindrical panel can be depicted in Figure 10. The relationship between V_CNT and V_FGM should satisfy with

$$V_{\mathit{CNT}}+V_{\mathit{FGM}}=1$$

(47)

Figure 10.

Uniform distribution of CNT.

The comparison of FGM and CNT-reinforced FGM at the central position is shown in Figure 11. Besides, four boundary conditions in this circumstance are also depicted.

Figure 11.

Deflection history with FGM and CNT reinforced FGM.

Figure 11 indicates that SWCNT greatly enhances the mechanical characteristics of FGM and strengthens the ability of composite structure to resist external force. As shown in Figure 12, clamp-clamp condition gives the smallest deflection in these four boundary conditions.

Figure 12.

Deflection history with four different boundary conditions.

Conclusions

To analyse the dynamic vibration of FGM cylindrical panels under different boundaries, a 3D semi-analytical method is proposed by introducing the Durbin’s numerical inversion method, state space approach, and differential quadrature method. Compared with many numerical methods, the semi-analytical results in this paper have higher accuracy, because there are fewer errors and limitations caused by the introduction of stress and displacement assumptions.

To compare the proposed approach with other theoretical methods from the literature, the natural frequencies show consistent results. Especially, the reliable comparative verification of transient response produced by the proposed method and FEM on different boundary conditions is presented. Convergence studies indicate that the proposed method has a quick convergence rate with growing sample point numbers along the length direction and layer numbers along the radial direction. For the study of different length/outer radius ratio, the flexibility of the whole structure becomes greater when the ratio increases. Thus, the cylindrical panel deforms easier under the same external force, and the vibration period will also become longer. The effect of thickness/outer radius ratio is also studied. The increase of thickness can resist external force more effectively, and the vibration period will also be shortened. The results compared with the power law and exponential law show completely different. To combine this characteristic with engineering, FGM can be designed to meet complex circumstances and optimize the structural performance. When the uniform distribution of CNT is added to FGM cylindrical panel, the mechanical properties of composite structure have been reinforced greatly. The proposed 3D semi-analytical method has high accuracy for the analysis of composite structures and can be used to study the behaviour of functionally graded porous material in the further research. This study can serve as a foundation for solving more complicated problems such as fluid–structure interaction of flexible pipe or thermal effect analysis of FGM in aerospace field.