Free Transverse Vibration Analysis of Axially Functionally Graded Tapered Euler-Bernoulli Beams through Spline Finite Point Method

A new model for the free transverse vibration of axially functionally graded (FG) tapered Euler-Bernoulli beams is developed through the spline finite pointmethod (SFPM) by investigating the effects of the variation of cross-sectional andmaterial properties along the longitudinal directions. In the proposed method, the beam is discretized with a set of uniformly scattered spline nodes along the beam axis instead of meshes, and the displacement field is approximated by the particularly constructed cubic B-spline interpolation functions with good adaptability for various boundary conditions. Unlike traditional discretization and modeling methods, the global structural stiffness and mass matrices for beams of the proposed model are directly generated after spline discretizationwithout needing elementmeshes, generation, and assembling.Theproposedmethod shows the distinguished features of high modeling efficiency, low computational cost, and convenience for boundary condition treatment. The performance of the proposed method is verified through numerical examples available in the published literature. All results demonstrate that the proposed method can analyze the free vibration of axially FG tapered Euler-Bernoulli beams with various boundary conditions. Moreover, high accuracy and efficiency can be achieved.


Introduction
Functionally graded (FG) materials were initially proposed by material scientists from Japan in 1984 as thermal barrier materials; since then, the development of engineering application and theoretical research on structures with FG materials has made a great breakthrough because of the particular properties of these materials, such as high stiffness and thermal resistance [1][2][3]. FG materials are multiphase composites with the volume fraction of phase varying through any desired direction [4]. FG material gradation along an arbitrary direction can be expressed with power (polynomial) laws and exponential laws, which are the most renowned theories for FG material description. The gradient variation for FG beams may also be oriented in the thickness direction, longitudinal direction, and bidirection (oriented in both the thickness and lengthwise directions).
During the past decade, the problems of the stability, dynamic response, and linear/nonlinear free vibration behavior of the FG structures with material properties varying along the thickness direction have been thoroughly investigated by researchers around the world, and lots of numerical and analytical solutions have been obtained [5][6][7][8][9][10][11]. However, given the existence of variable coefficients on the governing equation, the free vibration analysis of beams with nonuniform cross sections and variable material properties along the longitudinal direction becomes more complicated than the beams with uniform cross sections and homogenous materials [3]. And so far, few analytical solutions are found for the free vibration of the axially FG tapered beams due to the difficulty of mathematical treatment of the problem; for example, Elishakoff et al. [12,13] proposed the close-form solutions for axially FG tapered beam with material gradation that follows the polynomial law using the semi-inverse method; however, the proposed semi-inverse method cannot be applied for graded beams of any axial nonhomogeneity. Therefore, a large number of numerical methods have been proposed in the reported literatures by other researchers, according to which lots of numerical models, based on different beam theories, namely, the Euler-Bernoulli beam theory (EBT) [2,3,[14][15][16][17][18][19][20] and the Timoshenko beam theory (TBT) [21][22][23][24][25], have been constructed to investigate the free vibration of the tapered beams considering the material property homogeneity or nonhomogeneity. However, it is well understood that the effects of the shear deformation and rotary inertia of the cross section for the long slender beams can be neglected; thus the EBT models are capable of analyzing the vibration problems of the long slender beams and the reliable results with high accuracy can be achieved. In this paper, the authors mainly engage in establishing a new EBT model for the free vibration of axially tapered beams (long slender beams).
Based on EBT, Shahba et al. [2] firstly investigated the stability and free vibration of the axially FG Euler beams considering cross-sectional breadth and height both linearly tapered along the longitudinal direction by using the finite element method (FEM); then, in order to improve computational efficiency and modify the convergence conditions, Shahba and Rajasekaran [3] developed other models by using the differential transform element method (DTEM) and differential quadrature element method of the lowest order (DQEL) to investigate the same problems. Results show that the DTEM and DQEL models obtained better convergence and had fewer elements than the FEM model. However, with the increase of the taper ratio, all models proposed by Shahba et al. needed to solve a larger scale of eigenvalue equations in order to achieve a stable accuracy. Huang and Li [14] studied free vibration problems of axially FG tapered Euler beams through a new approach that transformed the governing equation into the Fredholm integral equation with consideration of the case of material properties in polynomial and exponential variations along longitudinal direction. Li et al. [15] developed the exact frequency equation for axially exponentially FG Euler beams, according to which the results showed that the natural frequencies noticeably depend on exponential gradient parameter and boundary conditions. Akgöz and Civalek [16] developed a model combined with the EBT and modified coupled stress theory to study the free vibration of the axially FG micro-Euler beams, in which the Rayleigh-Ritz solution method was utilized to obtain the natural transverse frequencies. Mao and Pietrzko [17] studied the free vibration of tapered Euler beams with a continuously exponential variation of breadth and a constant height along the length through Adomian decomposition method. Hsu et al. [18] developed a model to investigate the free vibration of homogenous tapered Euler beams with elastic supports through Adomian modified decomposition method. Rajasekaran [4] andÖzdemir and Kaya [19] applied differential transform method (DTM) to study the free bending vibration of rotating tapered Euler beams by considering material nonhomogeneity and homogeneity along the longitudinal direction. Banerjee et al. [20] applied the dynamic stiffness method (DSM) to analyze the free bending vibration of the rotating homogeneous tapered Euler beams; though the high accuracy was obtained, the modeling process was not convenient. Other models have been investigated based on TBT; Shahba et al. [21] developed the FEM model to study the free vibration of the axially FG taper Timoshenko beams. Rajasekaran and Norouzzadeh Tochaei [22,23] developed different models through DTM and DQEL to study axially FG Timoshenko beams by considering the effects of rotating and nonrotating. Ozgumus and Kaya [24] investigated model through DTM for the free bending vibration of rotating homogeneous tapered Timoshenko beams by considering the effects of taper ratio, rotating speed, and slenderness ratio. Yardimoglu [25] investigated a new FEM model with displacement functions derived from the coupled displacement field of transverse displacement and cross-sectional rotation to study the free bending vibration of the rotating Timoshenko beams.
Among the aforementioned numerical models for the axially tapered beams, such as FEM, DTM, DTEM, and DSM models, the disadvantages are obvious. (1) For the FEM models, the element generation and assembling are both needed in the modeling process; the displacement interpolation functions of FEM models were unable to fully capture the variation of the cross section and material properties along the longitudinal direction due to the homogenous uniform beams elements were utilized, and with the increase of the taper ratio and gradient parameter, the results' error and elements numbers increase noticeably [2,21,25]. (2) For the DTM and DTEM models, although both models obtained high accuracy, these two numerical models needed to solve a large scale of eigenvalue equations, thereby leading to a high computational cost; what is more, the DTM model encountered convergence problems in some cases [3,4,19].
(3) For the DSM models, the most common disadvantages were their complicated modeling process [20]. Therefore, the aforementioned numerical models need to be improved, and a new model with convenience modeling process, high accuracy, and low computational cost will benefit the analysis of the axially tapered Euler beams.
Spline finite point method (SFPM) was firstly proposed by Professor Qin in 1979 [26,27], who proposed a cubic B-spline function-based method for structural analysis. Given that the interpolation of structural displacement is constructed through a set of scattered spline nodes along one direction instead of meshes, the element generations and assembling are both unnecessary for the modeling process of SFPM; thus, the SFPM contains the excellent feature of no mesh distortion. The process of modeling and computer programming can be achieved easily, which leads to high computation accuracy and efficiency because the interpolation functions used in SFPM have explicit forms. The performances of SFPM on the structures with uniform cross section and material property that keeps constant along the longitudinal direction were validated by reported works on the FGM beams/plates, laminates, and other structures [28][29][30][31][32]. Li et al. [28] investigated the free vibration of piezoelectric smart FGM plates using SFPM. Then, Li et al. [29] investigated a new model through bidirectional B-spline finite point method (B-SFPM) to analyze the displacement and stress distribution of a piezoelectric laminate composite plate combined with the problems of material parameter identification; Zhou and Li [30] studied the free vibrations of asymmetrical sandwich plates with laminated faces and orthotropic core through SFPM. Qin et al. [31,32] developed different models through SFPM to study the nonlinear behaviors of composite plates and the seismic behavior of a box arch bridge. Other researchers applied the spline strip method to the free vibration of rectangular composite laminates [33,34], which also validated the advantages of application of the cubic B-spline functions. However, to the best of the authors' knowledge, the SFPM have not been applied to analyze the free vibration of the axially FG tapered beams before, and similar works or models have been rarely found in the reported literatures as well. In this study, the authors introduce the SFPM to analyze the vibration behavior of the axially FG tapered Euler-Bernoulli beams for the first time, which provide a novelty approach for solving such problems.
The main objective of this study is to develop new simple models to study the free vibrations of axially FG tapered Euler-Bernoulli beams by using SFPM. The new models are able to achieve faster modeling process, higher accuracy, and lower computational cost compared with those of reported numerical models. On the basis of SFPM, the axially FG tapered beam is discretized with a set of uniformly scattered spline nodes along the beam axis instead of meshes, and the displacement field is approximated by particularly constructed cubic B-spline interpolation functions, which possess good adaptability for various boundary conditions. Both the global structural stiffness and mass matrices for beams are directly generated after spline discretization, thereby making the modeling convenient. Because of the employment of presented displacement interpolation functions, the derived structural stiffness and mass matrices being in the explicit forms are capable of directly taking account of any variations of the material properties and cross-section profile for the axially FG beams directly. A suitable spline node number is discussed and selected to maintain high computational accuracy based on the convergence analysis. The processes of boundary condition treatment and derivation of the structural matrices are introduced in detail. The effects of the variation of cross-sectional and material properties along the longitudinal direction are investigated. The performance of the proposed method is verified with numerical examples available in the published literature. All results demonstrate that the proposed method can analyze the free vibration of axially FG tapered Euler-Bernoulli beams with various boundary conditions with quick convergence and high accuracy. Figure 1 shows the axially FG tapered beams in length . The beam with cross-sectional parameters, such as the crosssectional area ( ) and the second moment of inertia ( ), and with material properties, such as Young's modulus ( ) and the mass density ( ), continuously vary along the lengthwise direction, which are the functions of longitudinal coordinate . According to EBT, the following assumptions are made: (1) the cross section of FG beams remains plane after deformation; (2) the slenderness ratio ( /ℎ 0 ) is big enough that the effects of rotary inertia and shear deformation inplane can be ignored; (3) the transverse deformation is small enough that only the linear vibration response is considered. Given that the shear deformation and rotary inertia are not considered, the axial and transverse displacement fields for axially FG tapered beams are presented, respectively, as follows:

Structural Model
where , , and are the spatial coordinates and is the transverse displacement. The axial and shear strains are, respectively, obtained as follows by using (1): The stresses can be given as where is the axial stress. The shear stress is equal to zero because shear deformation is not considered in the model. The strain and kinetic energies of the axially FG tapered Euler beams can be written, respectively, as The following can be obtained by applying the Hamilton principle: where Π is the total potential energy function for the axially FG tapered beam, which can be derived by some manipulations from (4) and (5). The derivation of Π can be seen in a reported work [19] and can be expressed by The governing equations of motion for an axially FG tapered Euler beam undergoing transvers free vibration are derived as follows through the variation principle:

Definition of Cubic B-Spline Function.
The cubic B-spline function 3 ( ) is adopted in this study and is defined by the following expressions [27]: According to the characteristics of 3 ( ), its derivative and integration can be expressed by the following:    Figure 3 illustrates that the beam is modeled via SFPM considering an axially FG tapered Euler beam described in Figure 1. The beam is discretized along the beam axis through uniformly scattered spline nodes on the entire beam. The tapered beam can be divided by the interval ℎ in direction, and the local coordinates of the corresponding nodes in the axis of the beam can be presented by

SFPM Modeling.
where stands for the number of subdivisions of beams in direction, which is also called the maximum number of spline nodes; 0 denotes the coordinate of left end of beam; and ℎ is the length between two spline nodes.
In accordance with the modeling process for SFPM, the linear combination of cubic B-spline functions defined in (8), also called spline interpolation functions ( ), is applied to approximate the transverse displacement ( , ) of the beam, which can be expressed by where [ ] denotes the spline interpolation functions formed by the linear combination of cubic B-spline functions and { } represents the vector of unknown generalized displacement parameters.
In this paper, the spline interpolation functions [ ], which are given in the literature [27], can be constructed by where [ ] is ( +3)×( +3) transformation matrix and [ 3, ] are a group of cubic B-spline functions, which is expressed by Shock and Vibration ] .
In this study, spline interpolation functions [ ] are characterized as the following boundary conditions: where the subscript " " denotes the derivative of with respect to .
With (14) substituted into (6), the total potential energy function Π for the free transverse vibration of the axially FG tapered Euler beams can be expressed by in which where the subscript double dot of {} denotes the differentiation with respect to time; [ ] and [ ] are, respectively, the ( + 3) × ( + 3) stiffness and mass matrices; and {} stands for the ( +3)×1 vector of unknown generalized acceleration parameters.
The following is obtained through the variation principle: Moreover, the free vibration equation for axially FG tapered beam can be expressed as The following eigenvalue equation must be solved to perform free vibration analysis: where is the natural frequency and { } is the mode shape parameter vectors. However, the boundary conditions must be manipulated first before solving (27). Figure 1, the axially FG tapered beam with a rectangular cross section, whose breath and height are tapered linearly, is considered. Thus, the cross-sectional area ( ) and second moment of inertia ( ) vary along the beam axis:

Cross-Sectional Variation. As shown in
In (28), and ℎ are the breadth and height taper ratios, respectively; and 0 and 0 are the cross-sectional area and second moment of inertia, respectively, at = 0. The tapered beam becomes uniform when = ℎ = 0, and the tapered beam theoretically tapers to be a point at = when = ℎ = 1.

Material Property Variation.
In this study, two laws are considered for the variation of material property along the lengthwise direction of the beam in accordance with the reported literature [3,14].
(i) Material Type I: Power Law or Polynomial Law. The material properties of the axially FG beam are considered to vary along the lengthwise direction following the power law. The graded functions are simplified into a generalized form [3]: wherẽ( = 1, . . . , 1) and̃( = 1, . . . , 2) stand for the constant coefficients derived from the polynomial expansion of ( ) and ( ) and 1 and 2 represent the maximum power number for the expansion of ( ) and ( ). For homogenous beams, only the items̃0 = 0 and̃0 = 0 are not zero in (31).
(ii) Material Type II: Exponential Law. The material properties of the axially FG beam are considered to vary along the lengthwise direction following the exponential law, which can be simplified as wherẽ1 and̃2 correspond to the coefficient of the simplified form of the exponential functions, such as Young's modulus ( ) and mass density ( ), respectively. is the gradient parameter that describes the volume fraction change of both constituents involved. For homogenous beams, only the item̃1 is not zero in (32).

Stiffness and Mass
Matrices. According to (24), the calculation of [ ] and [ ] depends on the variation of the cross section and the graded functions of the material property.
Then, (33) is substituted into (24). The stiffness and mass matrices are obtained: in which (ii) The cross sections are considered to be defined by (28), and the material properties are defined by (32) (exponential law). Thus, the following is obtained: where represents the gradient index of the exponentially FG materials. Moreover̃1 1, ,̃1 2, ( = 0, . . . , 4) and 11, ,̃1 2, ( = 0, 1, 2) stand for the constant coefficients of the expanded form for the product of ( ) ( ) and ( ) ( ), respectively, with which the polynomial items and exponential items should be distinguished.
Then, (36) is substituted into (24), and the stiffness and mass matrices are obtained: Shock and Vibration 7 in which The matrices defined as (35) and (38) are presented in detail in Appendix A, which are the basic matrices for the derivation of stiffness and mass matrices. The integration methods are also proposed.

Boundary Conditions. Boundary conditions should be imposed on [ ]
and [ ] matrices before solving the eigenvalue equation of (27) to obtain the natural transverse frequencies for the beams. According to (14) and (22), the end constraint of the transverse displacement = ( , ) for beams in various boundary conditions can be expressed by { } in the proposed model as follows: (i) The following is obtained for pined-pined (P-P) beams: (ii) The following is obtained for the clamped-clamped (C-C) beams: (iii) The following is obtained for clamped-pined (C-P) beams: (iv) The following is obtained for clamped-free (C-F) beams: (ii) The following is obtained for C-C beams: (iii) The following is obtained for C-P beams: (iv) The following is obtained for C-F beams: where   After applying the boundary conditions, the stiffness and mass matrices are substituted into (27). The natural frequencies and the corresponding generalized mode shape parameters { } are determined by solving (27) with the standard computer subroutines. Then, the mode displacement at the arbitrary location along the beam axis can be calculated by substituting { } into (14).

Flowchart for Analysis of the Axially FG Tapered Beam via SFPM.
The flowchart for free vibration of axially FG tapered beam via SFPM is shown in Figure 4 to state the proposed method visually. It is obvious that the proposed SFPM model is convenient for studying the vibration problems of the axially tapered Euler beams.

Numerical Examples
The nondimensional natural transverse frequencies are  (12) to check the convergence of the presented method (i.e., SFPM). The first four nondimensional natural frequencies, , are tabulated in Table 1 and compared with the exact solutions [35] and other numerical methods [14]. The comparison indicated that the SFPM solution has rapid convergence. As increases from 5 to 25, the first three nondimensional natural frequencies are identical to the exact solutions. This finding indicates that the present method is efficient except for the case of = 5, in which the fourth and fifth nondimensional natural frequencies do not match the exact solution because of the insufficient order of the eigenvalue equation, that is, (27). However, the proposed model does not have any convergence problems. As increases from 15 to 25, computational time and cost increase, but the results of the natural frequencies via SFPM vary slightly. Thus, we take spline node number = 15 in the following examples. Moreover, the first five nondimensional natural transverse frequencies for C-C, C-P, and P-P beams are calculated and tabulated in Table 2 with = 15. This finding verifies that SFPM can be applied to analyze the free vibration problem of beams with various boundary conditions.

Effect of the Variable Cross Section.
The second example is devoted to the Euler-Bernoulli beam with a rectangular cross section and has a breadth and height that taper linearly, as shown in Figure 2. Thus, the cross-sectional area ( ) and the second moment of inertia ( ) that vary along axis are defined as (28), as mentioned earlier.
Firstly, the natural transverse frequencies of the tapered homogenous beams with various boundary conditions  (C-F, C-C, and P-P) are calculated by using SFPM, for all beams, = 15 is selected, and two cases of cross-sectional variation, namely, the case of only the height tapered linearly ( ℎ = 0 ) and the case of both breadth and height tapered linearly ( = ℎ = 0 ), are both investigated. In order to validate the performance of the proposed method, FEM models (3D), which are implemented in ANSYS (see in Figure 5), are employed to investigate the free vibration of the tapered homogenous beams. In the FEM models, 3D solid element SOLID64 is employed, and all beams are meshed into 1600 elements to ensure the computational accuracy for fundamental and high-order frequencies and to avoid producing the abnormal element due to the increase of taper ratios. The nondimensional natural transverse frequencies for C-F, C-C, and P-P beams are tabulated in Tables 3-5, which indicate that the results of SFPM model agree very well with those of FEM models and DTM model reported by Rajasekaran [4]. The taper ratio ( 0 ) plays an important role in the tapered homogenous beams. For all beams, the natural frequencies decrease with the increase of the taper ratio ( 0 ) except for the fundamental frequencies of C-F beams, which increase with 0 ; for C-F beams, it is obvious that the frequencies of the case of height tapered linearly ( ℎ = 0 ) are smaller than those of the case of both height and breadth tapered linearly ( = ℎ = 0 ), while, for C-C and P-P beams, the frequencies of the case of height tapered linearly ( ℎ = 0 ) are equal to those of the case of both height and breadth tapered linearly ( = ℎ = 0 ). However, in this study, given = 15 is employed, the SFPM model shows a considerable decrease of the dimensions of the eigenvalue equations compared with the FEM and DTM models [4], leading to remarkable decrease in computational cost. Furthermore, FEM models implemented in ANSYS are unable to analyze the free vibration of tapered where 0 and 0 are Young's modulus and mass density, respectively, at = 0.
The material properties defined by (47) follow the polynomial law distribution, and (47) can be simplified as (31), with 1 = 1 and 2 = 2:  Figure 6: Convergence of natural frequencies for axially FG tapered cantilever beams with respect to spline node number N ( = ℎ = 0.8).
in which̃= 0 ( = 0, 1) ; Using (48), the structural stiffness and mass matrices can be calculated by (34). In Figure 6, the convergence of natural frequencies for axially FG tapered cantilever beams ( = ℎ = 0.8) (material gradation following the polynomial law) with respect to the spline node number is discussed; it is noticeable that the natural frequencies achieve a convergence value when increase from 5 to 15 (about 4-digit precision), and when increase from 15 to 25, the frequencies vary a little but the computational cost increased remarkably. Therefore, = 15 is a suitable choice for the proposed SFPM models to ensure both the accuracy and efficiency for the free vibration of the axially FG tapered beams with different taper ratios.
The natural transverse frequencies for axially FG tapered beams with different breadth taper ratio ( ) and height taper ratio ( ℎ ) are calculated by using the proposed method, and the results are tabulated in Tables 6-8. The results of SFPM agree with the results of DTEM reported by Shahba and Rajasekaran [3]. This observation indicates that the proposed method can be used to investigate the free vibration of axially FG tapered beams with high accuracy. The variations of the first two nondimensional natural frequencies of the axially FG beams with various boundary conditions for different values of taper ratios ( , ℎ ) are illustrated in Figures 7-9. Figures 7-9 illustrate that the varying trend of the natural frequencies is distinguished from beams with different boundary conditions. In general, the natural transverse frequencies are considerably affected by the variation of height taper ratio ℎ . Figure 7 and Table 6 indicate that the increase of taper ratios ℎ and leads to the increase of the fundamental natural frequencies for C-F beams, whose second frequencies decrease with ℎ and slightly increase with . For C-C and P-P beams, Figures  8 and 9 and Tables 7 and 8 indicate that all frequencies remarkably decrease with ℎ , whereas the variations with the increase of are not monotonic. Thus, predicting the natural frequencies for axially FG tapered beams of different   boundary conditions with respect to the taper ratios is impossible.
Notice that the proposed model needs considerably less computational cost and modeling time compared with the DTM and DTEM models developed by Rajasekaran [4] and Shahba et al. [2], respectively. In the DTM and DTEM models, 132 simple points and 12 elements are used, and the eigenvalue equation in the dimension of 578 × 578 needs to be solved. However, in the proposed model, as = 15 is selected based on the convergence study, the eigenvalue equation is only in the dimension of 18 × 18, which shows that the proposed method has higher remarkable ascending efficiency than the reported methods. Furthermore, the natural transverse frequencies for tapered beams with consideration of material nonhomogeneity and homogeneity are distinguished by comparing the results in Table 3 with those in Table 4. The material nonhomogeneity defined by (47) decreases the natural frequencies for tapered C-F beams compared with homogenous beams. Thus, the material properties obviously affect the dynamic characteristics of the tapered beams.

Effect of Axial Gradient
Parameter. Finally, the free vibration for tapered rectangular beams with axially exponentially graded material is investigated to verify the adaptability of the proposed method to different gradient functions. Similar works have been conducted by Huang and Li [14]. The material properties are defined as follows: in which 1 and 2 correspond to the material properties, such as Young's modulus ( ) and mass density ( ) at = 0 and = , respectively, and is the gradient parameter that describes the volume fraction change of both constituents involved. The variation of ( ) against is shown in Figure 10 for 1 = 3 2 . The material properties defined by (50) can be classified into two types of distribution that depends on the value of gradient index ; one is the polynomial type that is defined by (31) ( = 0) with 1 = 2 = 1 and the other is exponential type that is defined by (32) ( ̸ = 0), which can be simplified into the following equations: Using (51), the stiffness and mass matrices can be calculated by (34) and (37) for = 0 and ̸ = 0, respectively. In the following calculations, the two materials chosen are aluminum (Al) and zirconia (ZrO 2 ). The material properties are presented by [14] ZrO 2 : The nondimensional frequencies are defined by  In consideration of the two cases in this paper, one is Al rich near the end = 0 and ZrO 2 rich near the end at = (Case 1), whereas the other is Al rich near the end at = and ZrO 2 rich near the end at = 0 (Case 2). The free vibration of axially FG beams with uniform cross-sectional ( = ℎ = 0) and material properties defined as (50) is calculated to verify the proposed method (i.e., SFPM). The results for the fundamental nondimensional natural transverse frequencies are tabulated in Table 9, and the high-order frequencies calculated with SFPM are tabulated in Table 10. Figure 11 gives the convergence of the natural frequencies for axially exponentially FG cantilever beams (C-F) ( = −10) with respect to spline node number ; the variations of the first three nondimensional natural frequencies of axially FG beams with respect to gradient parameter are depicted in Figures 12-15. Figure 11 indicates that convergence solutions of natural frequencies with high accuracy will be obtained for the proposed model of the axially exponentially FG beams at = 15. Therefore, for axially FG beams with different gradient index investigated in this section, let = 15 be the suitable choice for all models to achieve the demands for both the high accuracy and low computational cost. Table 9 indicates that the results of fundamental frequencies of the proposed model  Figure 11: Convergence of natural frequencies for axially exponentially FG cantilever beams with respect to spline node number N ( = −10). agree with those reported by Huang and Li [14] very well. The comparison of the results for Cases 1 and 2 listed in Tables 9  and 10 indicates that the gradient index is important for the natural frequencies.
As shown in Figures 12 and 13, the natural frequencies of P-P beams and C-C beams for Case 1 with index are identical to those for Case 2 with index − . However, the varying tendency of the natural frequencies for C-F and C-P beams with respect to gradient index is different because of the asymmetry of the boundary conditions, as shown in Figures 14 and 15. The fundamental natural transverse frequencies for Case 2 are larger than those for Case 1 with similar index . However, the high-order natural frequencies, which show a complex variation tendency, do not always follow this law. For all beams, the natural frequencies for Case 2 show an ascending tendency with the increase of gradient index, whereas the natural frequencies for Case 1 show an opposite trend.
The natural transverse frequencies of the tapered beams ( = ℎ = 0 ) with material nonhomogeneity defined by (50) ( = −10), which belongs to exponential law, are calculated by the proposed model to investigate the effects of the taper ratio, and the results are tabulated in Table 11, which indicate that all natural transverse frequencies for the beams with different boundary conditions generally decrease with the taper ratio, except for the beams with C-F boundary, whose fundamental natural frequencies (as shown in italic bold fonts in Table 11) increase with the taper ratio. By contrast, the high-order frequencies show an opposite trend, and the same varying trends have been found for the tapered beams with material gradation following the polynomial law (Section 4.2) with respect to the taper ratios.
The results obtained by the numerical examples in this section indicate that the proposed model shows good adaptability to the variation of cross sections, boundary conditions, and material gradient functions in either a polynomial or an exponential gradient function. And for the free vibration of tapered beams with different type of material nonhomogeneity, the results of higher accuracy and lower computational cost are obtained for the proposed SFPM model compared with other reported models. What is more, the natural frequencies of the tapered beams of the same boundary conditions with homogenous and nonhomogeneous materials vary in a similar trend toward the variation of the taper

Conclusion
In this study, a new model for the analysis of the free vibration of axially FG tapered Euler-Bernoulli beams is developed with SFPM. The effects of the variation in cross-sectional and material properties along the longitudinal directions are investigated. In the proposed approach, the beam is discretized with a set of uniformly scattered spline nodes along the beam axis direction instead of meshes, thereby resulting in convenience for modeling and programming. The displacement field of the beam is approximated with cubic B-spline interpolation functions, which shows good adaptability to different boundary conditions. Given the special features of SFPM, the global structural stiffness and mass matrices for axially FG tapered beams are directly generated after spline discretization instead of element meshes, generation, and assembling. The effectiveness of the proposed method is confirmed by comparing our numerical results with those available in the literature, including the case of rectangular beams with linearly varying cross-sectional breadth or height and axially FG tapered beams with polynomial or exponential function gradient. All results show that the proposed method can analyze the free vibration of Euler-Bernoulli beams with arbitrary axial gradient and varying cross sections in various boundary conditions. Moreover, a stable solution with high accuracy can be achieved with a small spline node number.     (34) and (37). The structural stiffness and mass matrices for axially FG tapered beams can be easily calculated.