Free Vibration of Axially Functionally Graded Beam

Abstract

Axially functionally graded (AFG) beam is a special kind of nonhomogeneous functionally gradient material structure, whose material properties vary continuously along the axial direction of the beam by a given distribution form. There are several numerical methods that have been used to analyze the vibration characteristics of AFG beams, but it is difficult to obtain precise solutions for AFG beams because of the variable coefficients of the governing equation. In this topic, the free vibration of AFG beam using analytical method based on the perturbation theory and Meijer G-Function are studied, respectively. First, a detailed review of the existing literatures is summarized. Then, based on the governing equation of the AFG Euler-Bernoulli beam, the detailed analytic equations are derived on basis of the perturbation theory and Meijer G-function, where the nature frequencies are demonstrated. Subsequently, the numerical results are calculated and compared, meanwhile, the analytical results are also confirmed by finite element method and the published references. The results show that the proposed two analytical methods are simple and efficient and can be used to conveniently analyze free vibration of AFG beam.

1. Introduction

Functionally gradient materials (FGMs) make a composite material by varying the microstructure from one material to another material with a specific gradient. It can be designed for specific function and applications. If it is for thermal or corrosive resistance or malleability and toughness, both strengths of the material may be used to avoid corrosion, fatigue, fracture, and stress corrosion cracking. FGMs are usually made into several structures, such as beams [1, 2, 3, 4], plates [5, 6, 7, 8], and shells [9, 10, 11, 12]. In this area, the variation of material properties in functionally graded beams may be oriented in transverse (thickness) direction or/and longitudinal/axial (length) direction.

For functionally graded beams with thickness-wise gradient variation, there have been many studies devoted to this topic. Lee et al. [13] establish an accurate transfer matrix method to analyze the free vibration characteristics of FGM beams whose Young’s modulus and density vary continuously with the height of the beam section through power law distribution. Su et al. [14] developed the dynamic stiffness method to investigate the free vibration behavior of FGM beams. Jing et al. [15] introduced a new approach by combining the cell-centered finite volume method and Timoshenko beam theory to analyze static and free vibration of FGM beams. Ait Atmane et al. [16] investigated the free vibration of a nonuniform FGM beams with exponentially varying width and material properties. Sina et al. [17] studied the free vibration of FGM beams by analytical method based on the traditional first-order shear deformation theory. Sharma [18] investigated the computational characteristics of harmonic differential quadrature method for free vibration of functionally graded piezoelectric material beam, which the material properties are assumed to have a power law or sigmoid law variation across the depth. Li et al. [19] proposed a high-order shear theory for free vibration of FGM beams with continuously varying material properties under different boundary conditions. Celebi et al. [20] employed the complementary function method to investigate the free vibration analysis of simply supported FGM beams, which the material properties change arbitrarily in the thickness direction. Chen et al. [21] studied the nonlinear free vibration behavior of shear deformable sandwich porous FGM beam based on the von kármán type geometric nonlinearity and Ritz method. Nazemnezhad and Hosseini-Hashemi [22] examined the nonlinear free vibration of FGM nanobeams with immovable ends using the multiple scale method.

As the FGMs are good for severe conditions, thermal-mechanical effect on FGM structures has attracted broad attention. In this field, Farzad Ebrahimi and Erfan Salari obtained outstanding achievements. Considering the thermal-mechanical effect and size-dependent thermo-electric effect, the buckling and vibration behavior of FGM nanobeams are studied [23, 24, 25, 26]. Considering the concept of neutral axis, they [27] studied the free buckling and vibration of FGM nanobeams using semi-analytical differential transformation method. To discuss the effect of the shear stress, Reddy’s higher-order shear deformation beam theory is introduced to study the vibration of the FGM structures [28, 29, 30]. Ebrahimi et al. [31, 32, 33] also studied vibration characteristics of FGM beams with porosities. Based on nonlocal elasticity theory, the nonlocal temperature-dependent vibration of FGM nanobeams were studied in thermal environment [34, 35, 36].

Another significant class of functionally graded beams is those with lengthwise varying material properties. It is difficult to obtain precise solutions for axially functionally graded (AFG) beams because of the variable coefficients of the governing equation. To solve this problem, a great deal of methods has been used to analyze the vibration characteristics of AFG beams. By assuming that the material constituents vary throughout the longitudinal directions according to a simple power law, Alshorbagy et al. [37] developed a two-node, six-degree-of-freedom finite element method (FEM) in conjunction with Euler-Bernoulli beam theory to detect the free vibration characteristics of a functionally graded beam. Shahba et al. [38, 39] used the FEM to study the free vibration of an AFG-tapered beam based on Euler-Bernoulli and Timoshenko beam theory. Shahba and Rajasekaran [40] studied the free vibration analysis of AFG-tapered Euler-Bernoulli beams employing the differential transform element method. Liu et al. [41] applied the spline finite point method to investigate the same problems. Rajasekaran [42] researched the free bending vibration of rotating AFG-tapered Euler-Bernoulli beams with different boundary conditions using the differential transformation method and differential quadrature element method. Rajasekaran and Tochaei [43] carried out the free vibration analysis of AFG Timoshenko beams using the same method. Huang and Li [44] studied the free vibration of variable cross-sectional AFG beams. The differential equation with variable coefficients is combined with the boundary conditions and transformed into Fredholm integral equation. By solving Fredholm integral equation, the natural frequencies of AFG beams can be obtained. Huang et al. [45] proposed a new approach for investigating the vibration behaviors of AFG Timoshenko beams with nonuniform cross section by introducing an auxiliary function. Huang and Rong [46] introduced a simple approach to deal with the free vibration of nonuniform AFG Euler-Bernoulli beams based on the polynomial expansion and integral technique. Hein and Feklistova [47] solved the vibration problems of AFG beams with various boundary conditions and varying cross sections via the Haar wavelet series. Xie et al. [48] presented a spectral collocation approach based on integrated polynomials combined with the domain decomposition technique for free vibration analyses of beams with axially variable cross sections, moduli of elasticity, and mass densities. Kukla and Rychlewska [49] proposed a new approach to study the free vibration analysis of an AFG beam; the approach relies on replacing functions characterizing functionally graded beams with piecewise exponential functions. Zhao et al. [50] introduced a new approach based on Chebyshev polynomial theory to investigate the free vibration of AFG Euler-Bernoulli and Timoshenko beams with nonuniform cross sections. Fang and Zhou [51, 52] researched the modal analysis of rotating AFG-tapered Euler-Bernoulli and Timoshenko beams with various boundary conditions employing the Chebyshev-Ritz method. Li et al. [53, 54] obtained the exact solutions for the free vibration of FGM beams with material profiles and cross-sectional parameters varying exponentially in the axial direction, where assumptions of Euler-Bernoulli and Timoshenko beam theories have been applied, respectively. Sarkar and Ganguli [55] studied the free vibration of AFG Timoshenko beams with different boundary conditions and uniform cross sections. Akgöz and Civalek [56] examined the free vibrations of AFG-tapered Euler-Bernoulli microbeams based on Bernoulli-Euler beam and modified couple stress theory. Yuan et al. [57] proposed a novel method to simplify the governing equations for the free vibration of Timoshenko beams with both geometrical nonuniformity and material inhomogeneity along the beam axis, and a series of exact analytical solutions are derived from the reduced equations for the first time. Yilmaz and Evran [58] investigated the free vibration of axially layered FGM short beams using experimental and FEM simulation, which the beams are manufactured by using the powder metallurgy technique using different weight fractions of aluminum and silicon carbide powders.

Till now, there also are plenty of literatures devoted to the free vibration for nonuniform beams. Boiangiu et al. [59] obtained the exact solutions for free bending vibrations of straight beams with variable cross section using Bessel’s functions and proposed a transfer matrix method to determine the natural frequencies of a complex structure of conical and cylindrical beams. Garijo [60] analyzed the free vibration of Euler-Bernoulli beams of variable cross section employing a collocation technique based on Bernstein polynomials. Arndt et al. [61] presented an adaptive generalized FEM to determine the natural frequencies of nonuniform Euler-Bernoulli beams. The spline-method of degree 5 defect 1 is proposed by Zhernakov et al. [62] to determine the natural frequencies of beam with variable cross section. Wang [63] studied the vibration of a cantilever beam with constant thickness and linearly tapered sides by means of a novel accurate, efficient initial value numerical method. Silva and Daqaq [64] solved the linear eigenvalue problem exactly of a slender cantilever beam of constant thickness and linearly varying width using the Meijer G-function approach. Rajasekaran and Khaniki [65] applied the FEM to research the vibration behavior of nonuniform small-scale beams in the framework of nonlocal strain gradient theory. Çalım [66] investigated the dynamic behavior of nonuniform composite beams employing an efficient method of analysis in the Laplace domain. Yang et al. [67] applied the power series method to investigate the natural frequencies and the corresponding complex mode functions of a rotating tapered cantilever Timoshenko beam. Clementi [68] analytically determined the frequency response curves of a nonuniform beam undergoing nonlinear oscillations by the multiple time scale method. Wang [69] proposed the differential quadrature element method for the natural frequencies of multiple-stepped beams with an aligned neutral axis. Abdelghany [70] utilized the differential transformation method to examine the free vibration of nonuniform circular beam.

The asymptotic development method, which is a kind of perturbation analysis method, is always used to solve nonlinear vibration equations. For example, Chen et al. [71, 72] studied the nonlinear dynamic behavior of axially accelerated viscoelastic beams and strings based on the asymptotic perturbation method. Ding et al. [73, 74] studied the influence of natural frequency of transverse vibration of axially moving viscoelastic beams and the steady-state periodic response of forced vibration of dynamic viscoelastic beams based on the multi-scale method. Chen [75] used the asymptotic perturbation method to analyze the finite deformation of prestressing hyperelastic compression plates. Hao et al. [76] employed the asymptotic perturbation method to obtain the nonlinear dynamic responses of a cantilever FGM rectangular plate subjected to the transversal excitation in thermal environment. Andrianov and Danishevs’kyy [77] proposed an asymptotic method for solving periodic solutions of nonlinear vibration problems of continuous structures. Based on the asymptotic expansion method of Poincaré-Lindstedt version [78], the longitudinal vibration of a bar and the transverse vibration of a beam under the action of a nonlinear restoring force are studied. The asymptotic development method is applied to obtain an approximate analytical expression of the natural frequencies of nonuniform cables and beams [79, 80]. Cao et al. [81, 82] applied the asymptotic development method to analyze the free vibration of nonuniform axially functionally graded (AFG) beams. Tarnopolskaya et al. [83] gave the first comprehensive study of the mode transition phenomenon in vibration of beams with arbitrarily varying curvature and cross section on the basis of asymptotic analysis.

The present topic focus on the free vibration of AFG beams with uniform or nonuniform cross sections using analytical method: the asymptotic perturbation method (APM) and Meijer G-function. First, the governing differential equation for free vibration of nonuniform AFG beam is summarized and rewritten in a form of a dimensionless equation based on Euler-Bernoulli beam theory. Second, the analytic equations are then derived in detail in Sections 3 and 4, respectively, where the nature frequencies are obtained and compared with the results of the finite element method and the published references. Finally, the conclusions are presented.

2. Governing equation of the AFG beam

This studied free vibration of the axially functionally graded beam, which is a nonuniform and nonhomogeneous structure because of the variable width and height, as shown in Figure 1. Based on Euler-Bernoulli beam theory, the governing differential equation of the beam can be written as

Figure 1.

The geometry and coordinate system of an AFG beam.

$$\frac{{\partial }^2}{\partial x^2}\left[E\left(x\right)I\left(x\right)\frac{{\partial }^{2}w\left(x,t\right)}{\partial x^2}\right]+\rho \left(x\right)A\left(x\right)\frac{{\partial }^{2}w\left(x,t\right)}{\partial t^2}=0,0\le x\le L$$ E1

where $w\left(x,t\right)$ is the transverse deflection at position x and time t; $L$ is the length of the beams; $E\left(x\right)I\left(x\right)$ is the bending stiffness, which is determined by Young’s modulus $E\left(x\right)$ and the area moment of inertia $I\left(x\right)$ ; and $\rho \left(x\right)A\left(x\right)$ is the unit mass length of beam, which is determined by volume mass density $\rho \left(x\right)$ and variable cross-sectional area $A\left(x\right)$ .

Because of the particularity of AFG beam, bending stiffness $E\left(x\right)I\left(x\right)$ and unit mass $\rho \left(x\right)A\left(x\right)$ will change with the axis coordinates, which makes the original constant coefficient differential equation become variable coefficient differential equation and to some extent increases the difficulty of solving. In order to facilitate calculation, we simplify the calculation process by introducing dimensionless parameters. Reference flexural stiffness ${\mathit{EI}}_0$ and reference mass $\rho A_0$ are introduced, and the above two dimensionless parameters are invariant. Suppose $E\left(x\right)I\left(x\right)={\mathit{EI}}_0+\overline{E\left(x\right)I\left(x\right)}$ and $\rho \left(x\right)A\left(x\right)=\rho A_0+\overline{\rho \left(x\right)A\left(x\right)}$ , where ${\mathit{EI}}_0$ and $\rho A_0$ are the invariant parts and $\overline{E\left(x\right)I\left(x\right)}$ and $\overline{\rho \left(x\right)A\left(x\right)}$ represent flexural stiffness and mass per unit length, respectively, and vary with the axial coordinates. Here, we introduce a dimensionless space variable $\xi =x/L$ and a dimensionless time variable $\tau =\frac{t}{L^2}\sqrt{\frac{{\mathit{EI}}_0}{\rho A_0}}$ ; Eq. (1) can be rewritten in the dimensionless form:

$$\frac{{\partial }^2}{\partial {\xi }^2}\left\{\left[1+f_1\left(\xi \right)\right]\frac{{\partial }^{2}w\left(\xi ,\tau \right)}{\partial {\xi }^2}\right\}+\left[1+f_2\left(\xi \right)\right]\frac{{\partial }^{2}w\left(\xi ,\tau \right)}{\partial {\tau }^2}=0,0\le \xi \le 1$$ E2

where

$$f_1\left(\xi \right)=\frac{\overline{E\left(\xi \right)I\left(\xi \right)}}{{\mathit{EI}}_0}\text{and}{f}_2\left(\xi \right)=\frac{\overline{\rho \left(\xi \right)A\left(\xi \right)}}{\rho A_0}$$ E3

are the nondimensional varying parts of the flexural stiffness and of the mass per unit length, respectively.

3. Asymptotic perturbation method

3.1 Equation deriving

In this section, the APM is introduced to obtain a simple proximate formula for the nature frequency of the AFG beam. Firstly, we assume that

$$w\left(\xi ,\tau \right)=W\left(\xi \right)sin\left(\mathit{\omega \tau }\right)$$ E4

where $W\left(\xi \right)$ is the amplitude of vibration and $\omega$ is the circular frequency of vibration. We obtain the following equation by substituting Eq. (4) with Eq. (2):

$$\frac{d^2}{d{\xi }^2}\left\{\left[1+f_1\left(\xi \right)\right]\frac{d^{2}W}{d{\xi }^2}\right\}-{\omega }^2\left[1+f_2\left(\xi \right)\right]W=0,0\le \xi \le 1$$ E5

To use the APM, a small perturbation parameter $\epsilon$ is introduced:

$$f_1\left(\xi \right)\to \epsilon f_1\left(\xi \right),f_2\left(\xi \right)\to \epsilon f_2\left(\xi \right)$$ E6

According to the Poincaré-Lindstedt method [78, 79, 80, 81, 82], we assume the expansion for $\omega$ and $W\left(\xi \right)$ as

$$\begin{array}{l}\omega ={\omega }_0+{\mathit{\epsilon \omega }}_1+{\epsilon }^2{\omega }_2+\dots ,\\ W\left(\xi \right)=W_0\left(\xi \right)+\epsilon W_1\left(\xi \right)+{\epsilon }^2{W}_2\left(\xi \right)+\dots .\end{array}$$ E7

Substituting these expressions with governing Eq. (5) and then expanding the expressions into a $\epsilon$ -series, Eqs. (8) and (9) are obtained by equating the coefficients of ${\epsilon }^0$ and ${\epsilon }^1$ to zero, yielding a sequence of problems for the unknowns ${\omega }_i$ and $W_i\left(\xi \right)$ :

$$\frac{d^4{W}_0}{d{\xi }^4}-{\omega }_0^2{W}_0=0$$ E8

$$\frac{d^4{W}_1}{d{\xi }^4}-{\omega }_0^2{W}_1+h_1\left(\xi \right)-2{\omega }_1{\omega }_0{W}_0=0$$ E9

where

$$h_1\left(\xi \right)=2\frac{{\mathit{df}}_1\left(\xi \right)}{\mathit{d\xi }}\frac{d^3{W}_0}{d{\xi }^3}+\frac{d^2{f}_1\left(\xi \right)}{d{\xi }^2}\frac{d^2{W}_0}{d{\xi }^2}+{\omega }_0^2\left[f_1\left(\xi \right)-f_2\left(\xi \right)\right]W_0$$ E10

For Eq. (8), the following general solution can be obtained:

$$W_0=Asin\left(\mathit{k\xi }\right)+Bcos\left(\mathit{k\xi }\right)+Csinh\left(\mathit{k\xi }\right)+Dcosh\left(\mathit{k\xi }\right)$$ E11

where

$$k=\sqrt{{\omega }_0}$$ E12

For simplicity, we consider clamped-free (C-F) beams, and the corresponding boundary conditions are

$$W_0=\frac{{\mathit{dW}}_0}{\mathit{d\xi }}=0,\xi =0$$ E13

$$\frac{d^2{W}_0}{d{\xi }^2}=\frac{d^3{W}_0}{d{\xi }^3}=0,\xi =1$$ E14

Then, the following equations from equation can be obtained:

$$\begin{array}{l}A+C=0\\ B+D=0\\ \frac{C}{D}=\frac{sink-sinhk}{cosk+coshk}\end{array}$$ E15

and the frequency equation is

$$coskcoshk+1=0$$ E16

The spatial mode shape can be obtained as

$$W_0=cosh\left(\mathit{k\xi }\right)-cos\left(\mathit{k\xi }\right)+\frac{C}{D}\left[sinh\left(\mathit{k\xi }\right)-sin\left(\mathit{k\xi }\right)\right]$$ E17

Now, the solution of the first-order equation is analyzed. In Eq. (9), both $h_1\left(\xi \right)$ and $W_1$ are linearly correlated with $W_0$ . Based on the theory of ordinary differential equations [84], the solvability conditions of linear differential equations can be expressed by the orthogonality of solutions of homogeneous systems of equations. At the same time, according to the orthogonality of modal vibration theory, the solution of Eq. (9) exists under the condition of the solvability of differential equation:

$${\int }_0^1\left[h_1\left(\xi \right)-2{\omega }_1{\omega }_0{W}_0\right]W_0\mathit{d\xi }=0$$ E18

is satisfied. As a result,

$${\omega }_1=\frac{{\int }_0^1{h}_1\left(\xi \right)W_0\mathit{d\xi }}{2{\omega }_0{\int }_0^1{W}_0^2\mathit{d\xi }}$$ E19

Because $h_1\left(\xi \right)$ is linearly correlated with $W_0$ , the former equations indicate that the arbitrary amplitude of $W_0$ does not impact ${\omega }_1$ . This finding yields the first-order correction of the natural frequency ${\omega }_0$ corresponding to a nonuniform and homogeneous beam.

Integrating by parts, we obtain

$$\begin{array}{l}{\int }_0^1{h}_1\left(\xi \right)W_0\mathit{d\xi }={\left.\left(\frac{{\mathit{df}}_1}{\mathit{d\xi }}\frac{d^2{W}_0}{d{\xi }^2}{W}_0+f_1\frac{d^3{W}_0}{d{\xi }^3}{W}_0-f_1\frac{d^2{W}_0}{d{\xi }^2}\frac{{\mathit{dW}}_0}{\mathit{d\xi }}\right)\right|}_0^1\\ +{\int }_0^1\left[f_1{\left(\frac{d^2{W}_0}{d{\xi }^2}\right)}^2-{\omega }_0^2{f}_2{W}_0^2\right]\mathit{d\xi }\end{array}$$ E20

By definition we have

$$f_1\left(\xi \right)=\frac{\overline{E\left(\xi \right)I\left(\xi \right)}}{E_{0}I}=\frac{E\left(\xi \right)I\left(\xi \right)-{\mathit{EI}}_0}{{\mathit{EI}}_0}$$ E21

so that

$$\begin{array}{l}{\left.\left(\frac{{\mathit{df}}_1}{\mathit{d\xi }}\frac{d^2{W}_0}{d{\xi }^2}{W}_0+f_1\frac{d^3{W}_0}{d{\xi }^3}{W}_0-f_1\frac{d^2{W}_0}{d{\xi }^2}\frac{{\mathit{dW}}_0}{\mathit{d\xi }}\right)\right|}_0^1+{\int }_0^1{f}_1{\left(\frac{d^2{W}_0}{d{\xi }^2}\right)}^2\mathit{d\xi }\hfill \\ =\left\{\frac{d\left[E\left(\xi \right)I\left(\xi \right)\right]}{{\mathit{EI}}_0\mathit{d\xi }}\frac{d^2{W}_0}{d{\xi }^2}{W}_0+\frac{E\left(\xi \right)I\left(\xi \right)}{{\mathit{EI}}_0}\frac{d^3{W}_0}{d{\xi }^3}{W}_0-\frac{E\left(\xi \right)I\left(\xi \right)}{{\mathit{EI}}_0}\frac{d^2{W}_0}{d{\xi }^2}\frac{{\mathit{dW}}_0}{\mathit{d\xi }}\right.\hfill \\ {\left.\left.+\frac{d^2{W}_0}{d{\xi }^2}\frac{{\mathit{dW}}_0}{\mathit{d\xi }}-\frac{d^3{W}_0}{d{\xi }^3}{W}_0\right\}\right|}_0^1+\frac{1}{{\mathit{EI}}_0}{\int }_0^{1}E\left(\xi \right)I\left(\xi \right){\left(\frac{d^2{W}_0}{d{\xi }^2}\right)}^2\mathit{d\xi }-{\int }_0^1{\left(\frac{d^2{W}_0}{d{\xi }^2}\right)}^2\mathit{d\xi }\end{array}$$ E22

Choosing the reference bending stiffness

$${\mathit{EI}}_0=\frac{{\left.\left\{\frac{d\left[E\left(\xi \right)I\left(\xi \right)\right]}{\mathit{d\xi }}\frac{d^2{W}_0}{d{\xi }^2}{W}_0+E\left(\xi \right)I\left(\xi \right)\frac{d^3{W}_0}{d{\xi }^3}{W}_0-E\left(\xi \right)I\left(\xi \right)\frac{d^2{W}_0}{d{\xi }^2}\frac{{\mathit{dW}}_0}{\mathit{d\xi }}\right\}\right|}_0^1+{\int }_0^{1}E\left(\xi \right)I\left(\xi \right){\left(\frac{d^2{W}_0}{d{\xi }^2}\right)}^2\mathit{d\xi }}{{\left.\left(\frac{d^3{W}_0}{d{\xi }^3}{W}_0-\frac{d^2{W}_0}{d{\xi }^2}\frac{{\mathit{dW}}_0}{\mathit{d\xi }}\right)\right|}_0^1+{\int }_0^1{\left(\frac{d^2{W}_0}{d{\xi }^2}\right)}^2\mathit{d\xi }}$$ E23

we have ${\left.\left(\frac{{\mathit{df}}_1}{\mathit{d\xi }}\frac{d^2{W}_0}{d{\xi }^2}{W}_0+f_1\frac{d^3{W}_0}{d{\xi }^3}{W}_0-f_1\frac{d^2{W}_0}{d{\xi }^2}\frac{{\mathit{dW}}_0}{\mathit{d\xi }}\right)\right|}_0^1+{\int }_0^1{f}_1{\left(\frac{d^2{W}_0}{d{\xi }^2}\right)}^2\mathit{d\xi }=0$ .

Analogously, we choose

$$\rho A_0=\frac{{\int }_0^1\rho \left(\xi \right)A\left(\xi \right)W_0^2\mathit{d\xi }}{{\int }_0^1{W}_0^2\mathit{d\xi }}$$ E24

giving ${\int }_0^1{f}_2{W}_0^2\mathit{d\xi }=0$ . Then, we obtain ${\omega }_1=0$ . These values are the properties of the equivalent homogeneous beam having the same frequency (at least up to the first order) as the given nonuniform AFG beam.

The nth natural circular frequency of the AFG beam can be derived as

$${\lambda }_n=\frac{1}{L^2}\sqrt{\frac{{\mathit{EI}}_0}{\rho A_0}}{\omega }_0$$ E25

Each order of frequency of ${\omega }_0$ can be determined by Eq. (16) (in turn, positive numbers from small to large). The required variables have been computed by the above expression. Eq. (25) is an approximate formula for the natural frequencies of variable cross-sectional AFG beams.

In order to show the applicability of this method, we study other supporting conditions, and we can easily get the corresponding boundary conditions of Eqs. (13) and (14). Due to the limited space, the detailed derivation process is omitted, and the final results are shown in Table 1.

Boundary conditions	Frequency equation	Mode shape
Simply supported (S-S)	$sink=0$	$W_0=sin\left(\mathit{k\xi }\right)$
Clamped-pinned (C-P)	$tank-tanhk=0$	$W_0=cosh\left(\mathit{k\xi }\right)-cos\left(\mathit{k\xi }\right)-\frac{coshk-cosk}{sinhk-sink}\left[sinh\left(\mathit{k\xi }\right)-sin\left(\mathit{k\xi }\right)\right]$
Clamped-clamped (C-C)	$coskcoshk-1=0$	$W_0=cosh\left(\mathit{k\xi }\right)-cos\left(\mathit{k\xi }\right)+\frac{sink+sinhk}{cosk-coshk}\left[sinh\left(\mathit{k\xi }\right)-sin\left(\mathit{k\xi }\right)\right]$

Table 1.

Frequency equations and mode shapes for various beams.

3.2 Numerical results and discussion

Based on the above analysis, four kinds of AFG beams with various taper ratios are considered, as shown in Figure 2. The numerical simulations are carried out, and the results are compared with the published literature results to verify the validity of the proposed method.

Figure 2.

Geometry and coordinate system of an AFG beam for different taper ratios: (a) case 1, $c_b=c_h=0$ ; (b) case 2, $c_h=0,c_b\ne 0$ ; (c) case 3, $c_b=0,c_h\ne 0$ ; and (d) case 4, $c_b\ne 0,c_h\ne 0$ .

In Figure 2, $B_L$ and $B_R$ are the width of the left and right ends of the beams, respectively, and $H_L$ and $H_R$ are the height of the left and right ends of the beams, respectively. Here, we assume that the geometric characteristics of AFG beams vary linearly along the longitudinal direction. Therefore, the variation of cross-sectional area $A\left(x\right)$ and moment of inertia $I\left(x\right)$ along the beam axis can be clearly expressed as follows:

$$A\left(x\right)=A_L\left(1-c_b\frac{x}{L}\right)\left(1-c_h\frac{x}{L}\right),I\left(x\right)=I_L\left(1-c_b\frac{x}{L}\right){\left(1-c_h\frac{x}{L}\right)}^3$$ E26

where $c_b=1-\frac{B_R}{B_L}$ and $c_h=1-\frac{H_R}{H_L}$ are the breadth and height taper ratios, respectively. $A_L$ and $I_L$ are cross-sectional area and area moment of inertia of the beam left sides, respectively. It is instructive to remember that if $c_b=c_h=0$ , the beam would be uniform; if $c_h=0,c_b\ne 0$ , the beam would be tapered with constant height; if $c_b=0,c_h\ne 0$ , the beam would be tapered with constant width; and if $c_b\ne 0,c_h\ne 0$ , the beam would be double tapered. These four cases are corresponding to Figure 2(a)–(d), respectively. Moreover, the material properties such as Young’s modulus $E\left(x\right)$ and mass density $\rho \left(x\right)$ along the beam axis are assumed as

$$E\left(x\right)=E_L\left(1+\frac{x}{L}\right),\rho \left(x\right)={\rho }_L\left[1+\frac{x}{L}+{\left(\frac{x}{L}\right)}^2\right]$$ E27

where $E_L$ and ${\rho }_L$ are Young’s modulus and mass density of the beam left sides, respectively.

Based on the introduced analytical equation, the first third-order nondimensional natural frequencies ( ${\mathrm{\Omega }}_n={\lambda }_n{L}^2\sqrt{{\rho }_L{A}_L/E_L{I}_L}$ ) of the four cases of nonuniform AFG beams with different boundary configurations were obtained. The results were listed in Tables 2–7, respectively, and it also was compared with those of published work by Shahba et al. [38].

Boundary condition		First mode	Second mode	Third mode
C-F	Present	2.439	18.437	54.339
	Ref. [38]	2.426	18.604	55.180
S-S	Present	9.053	35.834	80.470
	Ref. [38]	9.029	36.372	81.732
C-C	Present	20.585	56.251	109.869
	Ref. [38]	20.472	56.549	110.947

Table 2.

Nondimensional natural frequencies of the AFG uniform beam (case 1) with different boundary conditions.

$c_b$		C-F			S-S			C-C
		First mode	Second mode	Third mode	First mode	Second mode	Third mode	First mode	Second mode	Third mode
0.2	Present	2.613	18.887	54.951	9.068	35.957	80.772	20.457	56.196	110.003
	Ref. [38]	2.605	19.004	55.534	9.060	36.342	81.685	20.415	56.472	110.862
0.4	Present	2.854	19.483	55.753	9.088	36.117	81.165	20.294	56.124	110.177
	Ref. [38]	2.851	19.530	56.023	9.087	36.315	81.645	20.288	56.298	110.671
0.6	Present	3.214	20.311	56.853	9.113	36.332	81.697	20.079	56.028	110.411
	Ref. [38]	3.214	20.296	56.800	9.099	36.297	81.624	20.019	55.921	110.250
0.8	Present	3.832	21.542	58.453	9.147	36.638	82.456	19.783	55.892	110.743
	Ref. [38]	3.831	21.676	58.435	9.069	36.277	81.639	19.385	54.971	109.142

Table 3.

Nondimensional natural frequencies of the AFG-tapered beam with constant height (case 2) and different boundary conditions.

$c_h$		C-F			S-S			C-C
		First mode	Second mode	Third mode	First mode	Second mode	Third mode	First mode	Second mode	Third mode
0.2	Present	2.5054	17.2596	49.4982	8.1416	32.1888	72.2680	18.2420	50.1851	98.2992
	Ref. [38]	2.5051	17.3802	50.0491	8.1341	32.5236	73.1138	18.2170	50.4801	99.1734
0.4	Present	2.6293	16.2995	45.3519	7.2793	28.9717	65.1203	16.3027	45.0600	88.4345
	Ref. [38]	2.6155	16.0705	44.6181	7.1531	28.4747	63.9942	15.8282	44.0246	86.6272
0.6	Present	2.8535	15.6697	42.2358	6.4872	26.3694	59.4850	14.9152	41.2502	80.9747
	Ref. [38]	2.7835	14.6508	38.7446	6.0357	24.1101	54.0921	13.2293	36.9653	72.8740
0.8	Present	3.2889	15.5662	40.6554	5.7966	24.6371	55.9734	14.2233	39.1823	76.7690
	Ref. [38]	3.0871	13.1142	32.1309	4.6520	19.1314	42.6954	10.2235	28.7492	56.8109

Table 4.

Nondimensional natural frequencies of the AFG-tapered beam with constant width (case 3) and different boundary conditions.

$c_b$		0.2	0.4	0.6	0.8
$c_h$	First mode
0.2	Present	2.6873	2.9380	3.3113	3.9455
	Ref. [38]	2.6863	2.9336	3.2993	3.9219
0.4	Present	2.8226	3.0877	3.4796	4.1377
	Ref. [38]	2.7987	3.0486	3.4181	4.0471
0.6	Present	3.0640	3.3506	3.7700	4.4625
	Ref. [38]	2.9699	3.2237	3.5985	4.2355
0.8	Present	3.5271	3.8475	4.3081	5.0458
	Ref. [38]	3.2794	3.5401	3.9232	4.5695
$c_h$	Second mode
0.2	Present	17.7225	18.3289	19.1598	20.3725
	Ref. [38]	17.7501	18.2379	18.9501	20.2432
0.4	Present	16.7822	17.4061	18.2458	19.4418
	Ref. [38]	16.4092	16.8571	17.5139	18.7164
0.6	Present	16.1771	16.8214	17.6687	18.8380
	Ref. [38]	14.9567	15.3627	15.9616	17.0694
0.8	Present	16.0947	16.7493	17.5836	18.6877
	Ref. [38]	13.3850	13.7466	14.2848	15.2955
$c_h$	Third mode
0.2	Present	50.2194	51.1534	52.4114	54.1995
	Ref. [38]	50.3934	50.8645	51.6029	53.1332
0.4	Present	46.1970	47.2734	48.6925	50.6520
	Ref. [38]	44.9504	45.4003	46.0957	47.5129
0.6	Present	43.2042	44.4117	45.9613	48.0269
	Ref. [38]	39.0605	39.4844	40.1304	41.4236
0.8	Present	41.7065	42.9817	44.5636	46.5828
	Ref. [38]	32.4229	32.8123	33.3986	34.5521

Table 5.

Nondimensional natural frequencies of the AFG double-tapered beam (case 4); boundary conditions: C-F.

$c_b$		0.2	0.4	0.6	0.8
$c_h$	First mode
0.2	Present	8.1682	8.2018	8.2456	8.3051
	Ref. [38]	8.1462	8.1498	8.1336	8.0646
0.4	Present	7.3172	7.3647	7.4262	7.5089
	Ref. [38]	7.1455	7.1254	7.0794	6.9703
0.6	Present	6.5357	6.5960	6.6732	6.7754
	Ref. [38]	6.0082	5.9638	5.8868	5.7351
0.8	Present	5.8537	5.9240	6.0128	6.1283
	Ref. [38]	4.6046	4.5355	4.4264	4.2283
$c_h$	Second mode
0.2	Present	32.4133	32.7007	33.0819	33.6118
	Ref. [38]	32.5123	32.5079	32.5164	32.5326
0.4	Present	29.2971	29.7076	30.2419	30.9665
	Ref. [38]	28.4822	28.5003	28.5370	28.5928
0.6	Present	26.7834	27.2965	27.9493	28.8091
	Ref. [38]	24.1371	24.1791	24.2469	24.3497
0.8	Present	25.1032	25.6683	26.3683	27.2590
	Ref. [38]	19.1803	19.2509	19.3590	19.5300
$c_h$	Third mode
0.2	Present	72.8179	73.5237	74.4625	75.7732
	Ref. [38]	73.0959	73.0903	73.1116	73.1855
0.4	Present	65.9158	66.9202	68.2291	70.0069
	Ref. [38]	64.0054	64.0350	64.1007	64.2374
0.6	Present	60.4922	61.7392	63.3243	65.4089
	Ref. [38]	54.1330	54.1992	54.3126	54.5207
0.8	Present	57.0969	58.4547	60.1303	62.2530
	Ref. [38]	42.7677	42.8742	43.0436	43.3451

Table 6.

Nondimensional natural frequencies of the AFG double-tapered beam (case 4); boundary conditions: S-S.

$c_b$		0.2	0.4	0.6	0.8
$c_h$	First mode
0.2	Present	18.2779	18.3231	18.3818	18.4612
	Ref. [38]	18.1996	18.1286	17.9437	17.4566
0.4	Present	16.4975	16.7396	17.0484	17.4563
	Ref. [38]	15.8498	15.8350	15.7367	15.4025
0.6	Present	15.2512	15.6622	16.1771	16.8423
	Ref. [38]	13.2896	13.3319	13.3238	13.1529
0.8	Present	14.6662	15.2004	15.8587	16.6925
	Ref. [38]	10.3229	10.4255	10.5168	10.5339
$c_h$	Second mode
0.2	Present	50.4430	50.7713	51.2035	51.7981
	Ref. [38]	50.4565	50.3599	50.1017	49.3728
0.4	Present	45.6257	46.3346	47.2495	48.4763
	Ref. [38]	44.0553	44.0370	43.9027	43.4066
0.6	Present	42.0890	43.1214	44.4245	46.1234
	Ref. [38]	37.0509	37.1137	37.1104	36.8678
0.8	Present	40.2151	41.4614	42.9975	44.9420
	Ref. [38]	28.8912	29.0409	29.1842	29.2402
$c_h$	Third mode
0.2	Present	98.2992	99.7466	100.8219	102.313
	Ref. [38]	99.1474	99.0414	98.7543	97.9046
0.4	Present	88.4345	90.9806	92.8200	95.3023
	Ref. [38]	86.6608	86.6414	86.4932	85.9176
0.6	Present	80.9747	84.4598	86.8967	90.0855
	Ref. [38]	72.9681	73.0382	73.0375	72.7615
0.8	Present	76.7690	80.8426	83.5867	87.0576
	Ref. [38]	56.9674	57.1341	57.2991	57.3787

Table 7.

Nondimensional natural frequencies of the AFG double-tapered beam (case 4); boundary conditions: C-C.

Table 2 shows the first third-order natural frequencies of the AFG beam, case of Figure 2(a), which is uniform but nonhomogeneous. It can be clearly seen that the analytical results obtained from the asymptotic development method are in good agreement with those given by Ref. [38].

As can be seen from Tables 3 and 4, the first third-order dimensionless natural frequencies of AFG conical beams with only varying width or height are studied, respectively. It is easy to find the following conclusions. This method has higher accuracy on the equal height AFG-tapered beam. When the height changes, there is a certain fractional error in the AFG-tapered beam.

According to Figure 2(d), when the height and width of AFG beams change simultaneously, we can see that AFG beams are not uniform. The natural frequencies of three boundary conditions (free clamping, simply supported, and clamping) are studied in Tables 5, 6, 7. From the data in the table, it can be clearly found that the natural frequencies of AFG beams at low order are in good agreement with Ref. [38], while at high order, there are some errors in the natural frequencies.

4. The method of Meijer G-function

4.1 Equation deriving

In this section, the Meijer G-function is introduced to obtain the formula of the nature frequency of the AFG beam. Here, a special case of AFG beam is considered, where the cross section is uniform. Thus, in Eq. (1), Young’s modulus $E\left(x\right)$ and material mass density $\rho \left(x\right)$ are variable parameters, and the area moment of inertia $I$ and the cross-sectional area $A$ are invariant. To solve the governing equation, two parameters are firstly introduced to depict the functional gradient parameter equation:

$$E\left(x\right)=E_L\left(1-f_E\frac{x}{L}\right),\rho \left(x\right)={\rho }_L\left(1-f_{\rho }\frac{x}{L}\right)$$ E28

where $f_E=1-\frac{E_R}{E_L}$ , ${\rho }_E=1-\frac{{\rho }_R}{{\rho }_L}$ . $E_L$ and $E_R$ are Young’s modulus at the left/right end of the beam, and ${\rho }_L$ and ${\rho }_R$ are the mass density at the left/right end of the beam. Eq. (2) is then rewritten as

$${\left[\left(1-f_{E}x\right)w^{″}\right]}^{\prime \prime }+\left(1-f_{\rho }x\right)\ddot{w}=0$$ E29

Based on the vibration theory, we assume $w\left(x,t\right)=\varphi \left(x\right)q\left(t\right)$ , where $q_n\left(\mathrm{t}\right)=A_{n}cos{{\beta }_n}^{2}t+B_{n}sin{{\beta }_n}^{2}t$ and ${{\beta }_n}^2$ is the modal frequency for dimensionless. The governing equation is then derived as

$${\left[\left(1-f_{E}x\right){{\varphi }_n}^{\prime \prime }\right]}^{\prime \prime }-{{\beta }_n}^4\left(1-f_{\rho }x\right){\varphi }_n=0$$ E30

Next, Meijer G-function will be used to solve the linear partial differential equation. The general expression of Meijer G-function differential equation is written as

$$\left[{\left(-1\right)}^{\left(p-m-n\right)}z\prod _{l=1}^p\left(\mathrm{z}\frac{d}{\mathit{dz}}+1-a_l\right)-\prod _{k=1}^q\left(\eta \frac{d}{\mathit{dz}}-b_k\right)\right]G\left(\mathrm{z}\right)=0$$ E31

where $m,n,p$ and $q$ are integers satisfying $0\le m\le q,0\le \mathrm{n}\le p$ , G is the dependent variable also known as the Meijer G-function, z is the independent variable, and $a_l$ and $b_k$ are real numbers.

A definition of the Meijer G-function is given by the following path integral in the complex plane, called the Mellin-Barnes type:

$$G_{p,q}^{m,n}⟨\begin{array}{l}{a}_1\dots a_n,a_{n+1}\dots a_p\\ b_1\dots b_m,b_{m+1}\dots b_q\end{array}|z⟩=\frac{1}{2\mathit{\pi i}}{\int }_{\tau }\frac{{\prod }_{k=1}^m\mathrm{\Gamma }\left(\xi -b_k\right){\prod }_{k=1}^n\mathrm{\Gamma }\left(1-a_k+\xi \right)}{{\prod }_{k=1}^p\mathrm{\Gamma }\left(\xi -a_k\right){\prod }_{k=m+1}^q\mathrm{\Gamma }\left(1-b_k+\xi \right)}{z}^{-\xi }\mathit{d\xi }$$ E32

where an empty product is interpreted as 1, $0\le m\le q,0\le \mathrm{n}\le p$ , and the parameters are such that none of the poles of $\mathrm{\Gamma }\left({\mathrm{b}}_j-\xi \right),\left(j=1\dots m\right)$ coincides with the poles of $\mathrm{\Gamma }\left(1-a_j+\xi \right),\left(j=1\dots n\right)$ . Where $i$ is a complex number such that i² $=$ −1.

A special case of Eq. (31) can be expanded by assuming $n=p=0$ and $q=4$ . We can get that

$$\begin{array}{l}{z}^4\frac{d^{4}G}{{\mathit{dz}}^4}+\left(6-\sum _{k=1}^4{b}_k\right){\mathrm{z}}^3\frac{d^{3}G}{{\mathit{dz}}^4}+\left(7-3\sum _{n=1}^4{b}_k+\sum _{k,l=1}^4{b}_k{b}_l\right){\mathrm{z}}^2\frac{d^{2}G}{{\mathit{dz}}^2}\\ +\left[1-\sum _{k=1}^4{b}_k+\sum _{k,l=1}^4{b}_k{b}_l-\left(b_1{b}_2+b_1{b}_3+b_2{b}_3\right)b_4\right]z\frac{\mathit{dG}}{\mathit{dz}}\\ -\left[{\left(-1\right)}^{-m}z-\prod _{k=1}^4{b}_k\right]G=0\end{array}$$ E33

where $k\ne l$ . Although Eq. (30) is not similar to Eq. (33), the two equations can be similar by introducing some transformations:

$${\varphi }_n\left(x\right)=G\left({\mathrm{z}}_n\left(x\right)\right),z_n\left(x\right)={\left(\frac{{\beta }_n}{4f_E}\right)}^4{\left(1-f_{E}x\right)}^4$$ E34

Eq. (30) is transformed into

$${{\eta }_n}^3{G}^{\prime \prime \prime \prime }+5{{\eta }_n}^2{G}^{\prime \prime \prime }+\frac{69}{16}{\eta }_n{G}^{\prime \prime }+\frac{9}{32}{G}^{\prime }-\frac{1-f_{\rho }x}{1-f_{E}x}G=0$$ E35

Because of the difficulty of solving the differential equation with variable coefficients, we can simplify Eq. (35). Let $1-f_{E}x=1-f_{\rho }x=1-\mathit{Fx}$ ; it can be rewritten as

$${{\eta }_n}^3{G}^{\prime \prime \prime \prime }+5{{\eta }_n}^2{G}^{\prime \prime \prime }+\frac{69}{16}{\eta }_n{G}^{\prime \prime }+\frac{9}{32}{G}^{\prime }-G=0$$ E36

In order to solve the above equation, the coefficients of the ordinary differential Eqs. (33) and (36) are the same, so we can calculate the corresponding values, as shown in Table 8.

Case	$b_1$	$b_2$	$b_3$	$b_4$
1	1/2	1/4	0	1/4
2	1/4	1/4	0	1/2
3	0	1/4	1/2	1/4
4	1/2	0	1/4	1/4
5	0	1/2	1/4	1/4
6	0	1/4	1/4	1/2
7	1/4	0	1/4	1/2
8	1/4	1/2	0	1/4
9	1/4	0	1/2	1/4

Table 8.

Possible case of unknown constant $b_k$ of Meijer G-function equation.

One set of data can be selected from Table 8 and expressed in the form of closed solutions of Meijer G-function:

$$\begin{array}{l}{\phi }_{1n}\left(x\right)=G_{0,4}^{2,0}\left(\frac{1}{4},\frac{1}{4},0,\frac{1}{2}|\frac{{{\beta }_n}^4{\left(1-\mathit{Fx}\right)}^4}{256F^4}\right)\\ \end{array}$$ E37

$$\begin{array}{l}{\phi }_{2n}\left(x\right)=G_{0,4}^{1,0}\left(\frac{1}{4},0,\frac{1}{2},\frac{1}{4}|-\frac{{{\beta }_n}^4{\left(1-\mathit{Fx}\right)}^4}{256F^4}\right)\\ \end{array}$$ E38

$$\begin{array}{l}{\phi }_{3n}\left(x\right)=G_{0,4}^{1,0}\left(0,\frac{1}{2},\frac{1}{4},\frac{1}{4}|-\frac{{{\beta }_n}^4{\left(1-\mathit{Fx}\right)}^4}{256F^4}\right)\\ \end{array}$$ E39

$${\phi }_{4n}\left(x\right)=G_{0,4}^{1,0}\left(\frac{1}{2},\frac{1}{4},\frac{1}{4},0|-\frac{{{\beta }_n}^4{\left(1-\mathit{Fx}\right)}^4}{256F^4}\right)$$ E40

Modal modes of beams:

$${\varphi }_n\left(x\right)=C_{1n}{\phi }_{1n}\left(x\right)+C_{2n}{\phi }_{2n}\left(x\right)+C_{3n}{\phi }_{3n}\left(x\right)+C_{4n}{\phi }_{4n}\left(x\right),\mathrm{n}=1,2,3,\dots$$ E41

In order to determine the undetermined coefficients $C_i$ and ${\beta }_n$ , the boundary conditions of beams need to be considered:

1. C-F:

$$\left(\begin{array}{llll}{\phi }_{1n}\left(0\right)& {\phi }_{2n}\left(0\right)& {\phi }_{3n}\left(0\right)& {\phi }_{4n}\left(0\right)\\ {\phi }_{1n}^{\prime }\left(0\right)& {\phi }_{2n}^{\prime }\left(0\right)& {\phi }_{3n}^{\prime }\left(0\right)& {\phi }_{4n}^{\prime }\left(0\right)\\ {\phi }_{1n}^{″}\left(1\right)& {\phi }_{2n}^{″}\left(1\right)& {\phi }_{3n}^{″}\left(1\right)& {\phi }_{4n}^{″}\left(1\right)\\ {\phi }_{1n}^{‴}\left(1\right)& {\phi }_{2n}^{‴}\left(1\right)& {\phi }_{3n}^{‴}\left(1\right)& {\phi }_{4n}^{‴}\left(1\right)\end{array}\right)\left\{\begin{array}{l}{C}_{1n}\\ C_{2n}\\ C_{3n}\\ C_{4n}\end{array}\right\}=\left\{\begin{array}{l}0\\ 0\\ 0\\ 0\end{array}\right\}$$ E42

1. C-P:

$$\left(\begin{array}{llll}{\phi }_{1n}\left(0\right)& {\phi }_{2n}\left(0\right)& {\phi }_{3n}\left(0\right)& {\phi }_{4n}\left(0\right)\\ {\phi }_{1n}^{\prime }\left(0\right)& {\phi }_{2n}^{\prime }\left(0\right)& {\phi }_{3n}^{\prime }\left(0\right)& {\phi }_{4n}^{\prime }\left(0\right)\\ {\phi }_{1n}\left(1\right)& {\phi }_{2n}\left(1\right)& {\phi }_{3n}\left(1\right)& {\phi }_{4n}\left(1\right)\\ {\phi }_{1n}^{″}\left(1\right)& {\phi }_{2n}^{″}\left(1\right)& {\phi }_{3n}^{″}\left(1\right)& {\phi }_{4n}^{″}\left(1\right)\end{array}\right)\left\{\begin{array}{l}{C}_{1n}\\ C_{2n}\\ C_{3n}\\ C_{4n}\end{array}\right\}=\left\{\begin{array}{l}0\\ 0\\ 0\\ 0\end{array}\right\}$$ E43

1. S-S:

$$\left(\begin{array}{llll}{\phi }_{1n}\left(0\right)& {\phi }_{2n}\left(0\right)& {\phi }_{3n}\left(0\right)& {\phi }_{4n}\left(0\right)\\ {\phi }_{1n}^{″}\left(0\right)& {\phi }_{2n}^{″}\left(0\right)& {\phi }_{3n}^{″}\left(0\right)& {\phi }_{4n}^{″}\left(0\right)\\ {\phi }_{1n}\left(1\right)& {\phi }_{2n}\left(1\right)& {\phi }_{3n}\left(1\right)& {\phi }_{4n}\left(1\right)\\ {\phi }_{1n}^{″}\left(1\right)& {\phi }_{2n}^{″}\left(1\right)& {\phi }_{3n}^{″}\left(1\right)& {\phi }_{4n}^{″}\left(1\right)\end{array}\right)\left\{\begin{array}{l}{C}_{1n}\\ C_{2n}\\ C_{3n}\\ C_{4n}\end{array}\right\}=\left\{\begin{array}{l}0\\ 0\\ 0\\ 0\end{array}\right\}$$ E44

1. C-C:

$$\left(\begin{array}{llll}{\phi }_{1n}\left(0\right)& {\phi }_{2n}\left(0\right)& {\phi }_{3n}\left(0\right)& {\phi }_{4n}\left(0\right)\\ {\phi }_{1n}^{″}\left(0\right)& {\phi }_{2n}^{″}\left(0\right)& {\phi }_{3n}^{″}\left(0\right)& {\phi }_{4n}^{″}\left(0\right)\\ {\phi }_{1n}\left(1\right)& {\phi }_{2n}\left(1\right)& {\phi }_{3n}\left(1\right)& {\phi }_{4n}\left(1\right)\\ {\phi }_{1n}^{\prime }\left(1\right)& {\phi }_{2n}^{\prime }\left(1\right)& {\phi }_{3n}^{\prime }\left(1\right)& {\phi }_{4n}^{\prime }\left(1\right)\end{array}\right)\left\{\begin{array}{l}{C}_{1n}\\ C_{2n}\\ C_{3n}\\ C_{4n}\end{array}\right\}=\left\{\begin{array}{l}0\\ 0\\ 0\\ 0\end{array}\right\}$$ E45

4.2 Numerical results and discussion

Based on the above analysis, the natural frequencies of beams under different boundary conditions can be solved. Meanwhile, the results of finite element method are also conducted to verify the accuracy of the analytical results. Here, we use the power law gradient of the existing AFG beams [44], and the material properties of AFG beams change continuously along the axial direction. Therefore, the expressions of Young’s modulus $E\left(x\right)$ and mass density $\rho \left(x\right)$ are listed in detail:

$$Y\left(x\right)=\left\{\begin{array}{l}{Y}_L\left(1-\frac{e^{\mathit{\alpha x}/L}-1}{e^{\alpha }-1}\right)+Y_R\frac{e^{\mathit{\alpha x}/L}-1}{e^{\alpha }-1},\alpha \ne 0,\\ Y_L\left(1-\frac{x}{L}\right)+Y_R\frac{x}{L},\alpha =0.\end{array}\right.$$ E46

where $Y_L$ and $Y_R$ denote the corresponding material properties of the left and right sides of the beam, respectively. $\alpha$ is the gradient parameter describing the volume fraction change of both constituents involved. When gradient parameter $\alpha$ is equal and less than zero, Young’s modulus and mass density at the left end are less than those at the right end. When $\alpha$ equals zero, the beam is equivalent to a uniform Euler-Bernoulli beam, and Young’s modulus and mass density of the beam do not change with the length direction of the beam.

The variation of $Y\left(x\right)$ along the axis direction of the beam can be shown in Figure 3 for $Y_R=3Y_L$ . In order to show the practicability of this method, we choose the existing materials to study. The materials of AFG beams are composed of aluminum ( $\mathit{Al}$ ) and zirconia ( ${\mathit{ZrO}}_2$ ).The left and right ends of the beam are pure aluminum and pure zirconia, respectively. The material properties of AFG beams are given in detail in Table 9. We choose the sizes of commonly used beams which are L = 0.2 m, B = 0.02 m, and H = 0.001 m.

Figure 3.

Variation of the material properties defined by Eq. (46) with $Y_R=3Y_L$ .

Properties	Unit	Aluminum	Zirconia
E	GPa	70	200
$\rho$	Kg/m^3	2702	5700

Table 9.

Material properties of the AFG beam [44].

In order to verify the correctness of this method, some finite element simulation software is used to verify its correctness. In this paper, we analyze the natural frequencies of uniform AFG beams under different boundary conditions. In the process of finite element analysis, the AFG beam is transformed into a finite length model by using the delamination method [85]. At the same time, the AFG beam is delaminated along the axial direction. As shown in Figure 4, the material properties change along the axial direction, and the material properties of the adjacent layers are different. In order to analyze the performance of the beam, the uniform element is used to mesh each layer. In order to make the natural frequencies of AFG beams more precise, we can increase the number of layers and refine the finite element meshes.

Figure 4.

Finite segment model of the AFG beam.

In the Meijer G-function method, in order to solve the linear natural frequencies of beams under different boundary conditions, the determinant of the coefficient matrix of Eqs. (42)–(45) is equal to zero. Finally, linear natural frequencies of beams with different boundary conditions of the first four orders are listed in Table 10.

F	Order	C-F		C-P		S-S		C-C
		Present	FEM	Present	FEM	Present	FEM	Present	FEM
$\alpha =0.9$	1	2.4641	2.4651	4.0787	4.0789	3.0888	3.0891	4.5585	4.5579
	2	5.2251	5.2265	7.1520	7.1516	6.2895	6.2893	7.6920	7.6908
	3	8.2540	8.2560	10.2762	10.2778	9.4410	9.4420	10.8549	10.8560
	4	11.3209	11.324	13.4075	13.4092	12.5854	12.5869	14.0136	14.0148
$\alpha =0.5$	1	2.0774	2.0772	4.0055	4.0056	3.1344	3.1344	4.7098	4.7096
	2	4.8497	4.8491	7.1104	7.1100	6.2859	6.2861	7.8364	7.8360
	3	7.9496	7.9501	10.2396	10.2409	9.4278	9.4294	10.9827	10.9836
	4	11.0455	11.0652	13.3748	13.3760	12.5668	12.5707	14.1256	14.1273
$\alpha =1.7$	1	1.6098	1.6104	3.7738	3.7738	3.1279	3.1277	4.6896	4.6897
	2	4.4786	4.4792	6.9816	6.9816	6.2887	6.2884	7.8189	7.8187
	3	7.7326	7.7332	10.1477	10.1490	9.4315	9.4325	10.9679	10.9696
	4	10.9067	10.9079	13.3045	13.3042	12.5726	12.5737	14.1139	14.1158
$\alpha =2.7$	1	1.5370	1.5377	3.7214	3.7217	3.1188	3.1212	4.6629	4.6630
	2	4.4142	4.4142	6.9470	6.9471	6.2907	6.2904	7.7947	7.7948
	3	7.6920	7.6931	10.1207	10.1222	9.4350	9.4360	10.9482	10.9497
	4	10.8759	108,772	13.2807	13.2823	12.5762	12.5774	14.0972	14.0988

Table 10.

Comparisons between FEM and numerical calculation of linear dimensionless natural frequencies of AFG beams with different boundary conditions.

From Table 10, we can see that the results of finite element method are similar to those of Meijer G-function and the error is small. This can prove the accuracy of the method in frequency calculation on the one hand. In Figure 5, we can find that the first third-order dimensionless natural frequencies of C-F beams are in good agreement with FEM and numerical calculation. With the gradual increase of gradient parameter F, the dimensionless natural frequency of C-F beam increases gradually, and the change speed is accelerated. At the same time, the FEM and numerical simulation errors are very small, so the precise linear natural frequencies can be obtained.

Figure 5.

Dimensionless natural frequencies of C-F beams vary with parameter F: (a) fundamental frequencies, (b) second-order frequencies, and (c) third-order frequencies.

5. Conclusions

FGMs are innovative materials and are very important in engineering and other applications. Despite the variety of methods and approaches for numerical and analytical investigation of nonuniform FG beams, no simple and fast analytical method applicable for such beams with different boundary conditions and varying cross-sectional area was proposed. In this topic, two analytical approaches, the asymptotic perturbation and the Meijer G-function method, were described to analyze the free vibration of the AFG beams.

Based on the Euler-Bernoulli beam theory, the governing differential equations and related boundary conditions are described, which is more complicated because of the partial differential equation with variable coefficients. For both the asymptotic perturbation and the Meijer G-function method, the variable flexural rigidity and mass density are divided into invariant parts and variable parts firstly. Different analytical processes are then carried out to deal with the variable parts applying perturbation theory and the Meijer G-function, respectively. Finally, the simple formulas are derived for solving the nature frequencies of the AFG beams with C-F boundary conditions followed with C-C, C-S, and C-P conditions, respectively. It is observed that natural frequency increases gradually with the increase of the gradient parameter.

Accuracy of the results is also examined using the available data in the published literature and the finite element method. In fact, it can be clearly found that result of the APM is more accurate in low-order mode, which is caused by the defect of the perturbation theory. However, the APM is simple and easily comprehensible, while the Meijer G-function method is more complex and unintelligible for engineers. In general, the results show that the proposed two analytical methods are efficient and can be used to analyze the free vibration of AFG beams.

Acknowledgments

The authors gratefully acknowledge the support of the National Natural Science Foundation of China (Grant Nos. 11672008, 11702188, and 11272016).