Free longitudinal vibrations of functionally graded tapered axial bars by pseudospectral method

In this work, the problem of free longitudinal vibration of rods with variable cross-sectional area and material properties is investigated using the pseudospectral method. With the gradation of material properties like modulus of elasticity and mass density in the axial direction, the results corresponding to a functionally graded axial bar are obtained using the proposed pseudospectral formulation. The pseudospectral formulation used is relatively easy to implement and powerful in analyzing vibration problems. With the help of several numerical examples, the non-dimensional natural frequencies of rods obtained using the pseudospectral method are compared with those obtained by the analytical solution, generalized finite element method, the discrete singular convolution method and differential transformation method. The numerical results obtained show that the proposed technique allows boundary conditions to be incorporated easily and yields results with good accuracy and faster convergence rates than other methods.


Introduction
Functionally graded (FG) material structures find extensive use in modern engineering as their material properties can be tailored to meet the requirements of different applications [1].The FG axial bars are a class of inhomogeneous rods/bars with material properties varying continuously in desired spatial directions.The inhomogeneous rods are known to provide a suitable distribution of strength and weight for engineering structures [2].A study of the vibration characteristics of these rods is a subject of considerable scientific interest that has wide applications in aerospace, civil and mechanical engineering [3][4][5].
The governing differential equations of motion for longitudinal vibration of non-uniform rods have variable coefficients introduced by variable cross-sectional area.In addition, for functionally graded tapered axial bars, the varying material properties add up to the previous variable coefficients in the governing equation increasing the complexity in vibration analysis of these bars/rods.In general, methods of vibration analysis are classified as analytical and numerical methods while some are semi-analytical that use a combination of both the approaches.It is generally agreed that obtaining analytical/closed-form solutions of this vibration problem is possible only for some specific geometric and material functions and are often cumbersome.Application of numerical methods becomes essential to obtain the solution of the problem for general cases.There are many different methods for the numerical solution of differential equations which include Finite Element Method (FEM), Differential Quadrature Method (DQM), Differential Transformation Method (DTM), Discrete Singular Convolution (DSC) method and Spectral/Pseudo Spectral methods.Several studies have been dedicated to the problem of exact/analytic solutions for longitudinal vibrations of non-uniform rods [6][7][8][9][10][11][12][13] while some researchers [3,[14][15][16][17] have developed numerical methods.
Eisenberger [6] obtained the exact longitudinal vibration of a rod with polynomial variation in the cross-sectional area and mass distribution along the member employing an exact element method.Bapat [7] introduced an efficient approach to solve the vibration problem of a rod composed of  uniformly tapered sections and studied the vibrations of a rod composed of uniform, conical, exponential and catenoidal tapers.Abrate [8] has shown that there is a class of non-uniform rods for which the equations of motion can be transformed into the wave equation and considered area of variation of the form  (1 +   ⁄ ) .Kumar and Sujith [9] obtained exact analytical solutions for the longitudinal vibration of rods with polynomial and sinusoidal area variations of form ( + ) ,  sin ( + ).Li [10] carried out a functional transformation of the governing differential equation and presented the exact solutions for certain functional forms of an involved parameter.The differential transformation method was used in [11] to obtain exact solutions of the free axial vibration of a tapered bar with linearly varying mass and area.Raj and Sujith [12] derived a general analytical technique that helps in obtaining cross-sectional area variations which give specific, closed-form solutions for the longitudinal vibration of rods and have presented numerical solution for area variation of the form  exp( ),  exp().Yardimoglu and Aydin [13] obtained the exact longitudinal vibration characteristics of rods with sinusoidal area cross-sections of the form sin ( + ) , in terms of associated Legendre functions.In [14], a domain partition power series method for the vibration analysis of variable cross-section rods is presented and the longitudinal vibrations of a rod with linearly varying area is studied.The differential quadrature method was used in [15] to analyze the free vibrations of a general non-uniform rod.In [16], the free vibration analysis of elastic rods was investigated using a novel global collocation method.An adaptive generalized finite element (GFEM) method was proposed in [17] to study the free longitudinal vibration analysis of straight bars and trusses.The results obtained using GFEM for sinusoidal and polynomial area variations of rod were compared with those obtained using finite element and composite elements methods.Several closed form solutions have been derived by the semi-inverse method for the problem of determining eigenvalues of inhomogeneous structures in [18].It is also shown in [19] than an inhomogeneous rod can possess an exponential mode shape and derived closed-form solutions that can be utilized as model solutions for verification purposes.Recently, [3] investigated the longitudinal free vibrations of non-uniform rods by using the Discrete Singular Convolution (DSC) approach.The non-dimensional natural frequencies, obtained were compared with those obtained using exact methods and the differential quadrature method.
In most of the aforementioned works, the focus was on vibration of rods having variable area cross-sections.Recently, the longitudinal vibration of a functionally graded tapered axial bar was considered in [20] wherein, the material gradation was accounted for by varying the modulus of elasticity and mass density.Two new approaches based on Differential Transformation Method (DTM), namely Differential Transform Element Method (DTEM) and Differential Quadrature Element method of Lowest order (DQEL) were introduced by them to numerically study the free longitudinal vibration of tapered functionally graded (FG) rods.Further, in [21], the longitudinal vibrations of a  linearly tapered axial bar with exponential variation in mass density and modulus of elasticity was investigated using Finite Element Method (FEM) by proposing a new element that could be derived by using the basic principles of structural mechanics.
The purpose of this paper is to implement the pseudospectral method for the longitudinal vibration analysis of tapered  axial bars.Pseudospectral methods (PSM) are powerful tools for solving differential equations with high accuracy in a simple domain.The method has been successfully applied to the free vibration analysis of curved Timoshenko beams [22], Timoshenko beams and axisymmetric Mindlin plates [23], linear and non-linear beams [24], cylindrical helical springs [25] and Functionally graded Timoshenko beams [26].There have been different approaches in implementing the pseudospectral formulations based on their treatment of boundary conditions.In this work, an efficient formulation of the PSM involving a change of independent variable that greatly simplifies computer programs for solving differential equations is used to transform the governing differential equation together with the boundary conditions into an SRI HARIKRISHNA PILLUTLA, SUDHEER GOPINATHAN, VASUDEVA RAO YERIKALAPUDY eigenvalue problem, in which the eigenvalue represents the physical property of the vibration problem, namely the natural frequency.The problem is then solved using standard eigensolvers.The remainder of the paper is organized as follows.In Section 2, we introduce the pseudospectral method.In Section 3, we outline the problem and the solution methodology.In Section 4, the method is applied to examples that outline the convergence behavior, perform validations and demonstrate the advantages of the method.

Pseudospectral method
The spectral methods arise from the fundamental problem of approximation of a function by interpolation on an interval and are very much successful for the numerical solution of differential equations [27].The PS method can be considered as a spectral method that performs a collocation process and delivers exponential (spectral) accuracy under the condition of smoothness.The starting point for the Pseudospectral method is the representation of the approximate solution ( ()) of the solution (()) of the given differential equation by a linear combination of basis functions ( ()), ( = 0,1,2, … ., ): The best choice for  () are the eigenfunctions of a singular Sturm-Liouville problem such as the Chebyshev polynomials.Although the Legendre, Hermite and Laguerre polynomials are eigenfunctions of Sturm-Liouville problem, Chebyshev polynomials are more effective in handling the boundaries [28].The Chebyshev polynomials are recursive orthogonal polynomials that can be defined as: where  is an integer [29].Although defined for all  , Chebyshev polynomials are a stable representation only on (−1, 1).
In PS method, one is usually content with obtaining an approximation to the solution on a discrete set of grid points.The PS method, demands that the differential equation be exactly satisfied at these set of points known as collocation points [30].The accuracy, stability and rate of convergence of the numerical solutions depend on the choice of the collocation point [31] and the best choice for the collocation points is to use nodes that are clustered near the edges of the interval with an asymptotic density proportional to (1 −  ) ⁄ as  → ∞ [32].The Chebyshev polynomial extreme points that are the Chebyshev-Gauss-Lobatto sampling grid points on −1, 1 are known to have such density properties [33].Therefore, these points are chosen as collocation points in this work.The internal points which are the extrema of  () are given by: with  = 0,  denoting the end points and  = 1, 2, …,  − 1 the internal points.These points correspond to: The pseudospectral method associates the grid of collocation points with the basis set and the coefficients  in the expansion (1) are found by requiring the residual function  −  where  is the operator of the differential equation  = () to vanish at the collocation points.There are several ways of implementing the PS method of which the differentiation based formulation [34,35], integration based formulation [36,37] and straight forward implementation [22,23] are normally used.The differentiation based formulation is widely used in which one can construct the differentiation matrices using the recurrence finite-differences formulae on arbitrarily spaced grids [34,38] or utilizing the approach followed in [35,39].However, incorporation of boundary conditions when using differentiation matrices is complicated.The present work utilizes the straight forward implementation utilizing the Chebyshev polynomial-to-cosine change of variable together with a proper selection of internal collocation points based on the boundary conditions.

Mathematical formulation and solution methodology
The differential equation that governs the free longitudinal vibration of a axially graded bar reads as [40]: Here  denotes the axial coordinate, (, ) denotes longitudinal displacement at any position  and time .(), () and () are the modulus of elasticity, cross-sectional area and mass density respectively, which are assumed to vary along the axial coordinate  with  and  being the end coordinates of the rod.The classical boundary conditions considered are: Fixed-Free bar: The longitudinal displacement is assumed to be separable in space and time where the time dependence is harmonic with angular frequency  i.e., (, ) = () where () represents the mode shape.Denote ()() = () where () is the axial rigidity that varies along the axial coordinate .
Eq. ( 5) becomes () = − ()(): where: The governing equation becomes: We now outline the solution methodology for Eq. ( 8) using PSM.The natural interval  ∈ −1,1 in which the collocation points are located may be adapted to the interval of interest  ∈  ,  by the transformation: where ℎ =  −  .Due to change of variable: Applying the transformation, the solution function () is now a function of  and is assumed to be of the form: Here: It is to be noticed that the trigonometric form of the Chebyshev polynomial which is a change of coordinate is used for the representation of the approximate solution.This is coupled with analytic differentiation and collocation at internal points followed by enforcement of boundary conditions at end points.
The transformed governing differential equation is: If we collocate Eq. ( 17) at the points given by Eq. ( 4) and enforce the boundary conditions at  = 0,  =  which are given by: Fixed-Fixed bar: ( = 0) = ( = ) = 0, (18a) Fixed-Free bar: ( = 0) =  ( = ) = 0. (18b) Then we get a system of equations that is expressed as a matrix eigen value equation.The resulting equation is then solved using a standard eigensolver.

Numerical results and discussion
In this section we examine the convergence behavior of the PS method, validate the accuracy by comparing the results obtained with other analytical/numerical methods.Finally, we apply the method to obtain the non-dimensional frequencies of vibration of a functionally graded axial bar.In all the examples considered, the rods are of unit length.The details of the work carried out in this section are presented in Table 1.A diagrammatic view of typical rods considered herein is given in Fig. 1.

Convergence behaviour of PS method
To study the convergence behaviour of the PS method in obtaining the frequency values of the vibration of non-uniform rods, some examples have been examined.As a first example, we consider the longitudinal free vibration of a fixed-free non-uniform bar with polynomial variation of cross-section area, constant elasticity modulus  and constant density  and unit length.The cross-section area varies as () =  (1 + 2) where  is a reference cross-section area.Abrate [8] has presented exact analytical solution for longitudinal free vibration of bars with polynomial area variations.The exact non-dimensional frequency parameters of rods given [8] are set as basis values and the relative error of the PS method with  values varying from 10 to 16 is analyzed.The relative errors of the numerical values of the first six non-dimensional natural frequencies given by Ω =    ⁄ ; computed using the PS method are presented in Table 2.It is observed that the proposed method requires less number of grid points to obtain a relative error of the order of 10 -5 in the computation of the first six non-dimensional frequency values of a non-uniform rod with polynomial variation of cross section area.
As a second example we examine the convergence of the PS method by considering the free axial vibration of a tapered fixed-free bar.The area of cross-section () and mass density are assumed to vary linearly i.e. () = 2  and () = ()() = 2  where  is the characteristic area and  the characteristic mass.The computations are carried out for a unit length rod.The exact first two non-dimensional frequency values given in [11] are taken as basis values and the relative error of the PS method with  varying from 6 to 9 are presented in Table 3.It is observed that the present method achieves a relative error of order 10 -6 with fewer grid points.We further study the convergence of the frequency values obtained for fixed-fixed and fixed-free non-uniform rods with cross-section area given by () =  .exp (−0.5 ).We examine our results with the non-dimensional frequency values presented in [3] as basis values.Using the four digit DSC results that were computed in [3] with 500 grid points for a unit length rod, the relative error of the PS method for fixed-fixed and fixed-free rods with  taking values  = 20, 22, 25 and 28 is presented in Table 4.The relative errors of the non-dimensional frequency parameter given by Ω =    ⁄ for the first 13 modes are given in the table.It is observed that precision of the order of 10 -5 is obtained for all the 13 modes with  = 28 (29 collocation points).

Comparison with GFEM
In this example, the longitudinal free vibration of a fixed-fixed non-uniform rod with area variation given by () =  ( + 1) where  is a reference cross-section area is analyzed and the non-dimensional frequencies obtained using PS method are compared with the corresponding results obtained using GFEM method.The GFEM method is an iterative approach whose main goal is to increase the accuracy of eigenvalues related to a chosen vibration mode.In the method [17] each tangent frequency was obtained by different iterative analysis.The results presented in [17] brought out the narrow precision of GFEM than the c-version of composite element method and the ℎ-version of FEM in free longitudinal vibration analysis of uniform and non-uniform straight bars for the same degree of freedom.For the considered example, the results of the PS method have been obtained using  = 16 (17 collocation points).Table 5 presents the first four non-dimensional eigenvalues    ⁄ obtained by these methods.The results show that the PS method which is easier and simple to implement achieves the same accuracy as that of the adaptive GFEM with relatively few grid points.

Rods with sinusoidal variation of cross-section area
A numerical study of the longitudinal vibrations of non-uniform rods with area variations of the form () =  sin ( + 1) is carried out.The non-dimensional eigenvalues are computed with  = 29 (i.e. 30 collocation points) for various parameter (, , ) values under different boundary conditions and the lowest six eigenvalues for each end conditions are listed in the tables.The non-dimensional frequencies of fixed-fixed and fixed-free non-uniform rods with variation in area given by () =  sin ( + ),  =  − 2,  = 1 for  = 1, 2, 3, 4 is computed.The computations have been carried out with  = 29 (30 collocation points) to obtain a precision of 6 digits.The results obtained are given in Tables 6 and 7 and are found to be in good agreement with the exact results obtained by [13].

Rods with exponential variation of cross-section area
In this section, the variation in the natural frequencies of rod with exponential area variation of the forms () =  exp ( ) and () =  exp () are analyzed.The length of the rod is considered to be one unit with  = 0.1 units and  = 1.1 units as in [12].The fundamental frequencies for fixed-fixed and fixed-free rods are computed using PS method for  = -2.5, -2.0, -1.5, -0.5, 0.5, 1.0, 1.5, 2.0, 2.5 and  = -3, -1, 1 and 3 in the first area variation.A comparison of the values obtained is carried out with the corresponding the exact values in [12].The frequency values with 6 significant digits for the three end conditions are given in Tables 8, 9 and have been obtained using  = 29.The results bring out the accuracy of the method.Now, consider a fixed-fixed non-uniform rod with area variation given by () =  exp () The natural frequencies (Ω) for different values of the parameters are given in Table 10.The fundamental, first overtone and second overtone frequency values obtained using the PS method are presented in the tables for  = -2.5, -1.5, -0.5, 0.5, 1.5, 2.5 and  = -3, -1, 1 and 3.These values are presented along with the corresponding exact values given in [12].Though exact values are available for 3 significant digits, the approximate values are given for 6 significant digits for future comparison.It is observed from the results that the proposed PS formulation achieves good accuracy/convergence with less computational effort.-3 -1 1 3 Raj and Sujith [12] PS method Raj and Sujith [12] PS method Raj and Sujith [12] PS method Raj and Sujith [12] PS method -2.5

Functionally graded tapered axial bar
In the first example, we consider a tapered bar with a rectangular cross-section whose breath taper ratio is  , height taper ratio is  and whose cross-section is given by () =  (1 −  )(1 −  ).The modulus of elasticity and mass density vary as () =  (1 + ) and () =  (1 +  +  ) respectively.For such an FG axial bar, Shahba and Rajasekharan [20] computed the non-dimensional frequencies using DQEL and DTEM methods.Both the methods are known to considerably improve the convergence rate of conventional Differential Transform Method (DTM).PSM is used in the present work to compute the non-dimensional longitudinal frequencies given by  =    ⁄ , ( = 1,2,3) with  = 10 (11 collocation point) under Fixed-Free and Fixed-Fixed boundary conditions.The first three frequency values computed using the PSM method are compared with those obtained using DQEL/DTEM methods in Tables 11,12.It is observed that the frequencies are significantly dependent on the cross-section area which is in turn dependent on both the height and breadth taper ratios.The accuracy of the proposed method is evident from the results presented in the Tables 11, 12.The ease of implementation is the added advantage.In the second example we consider a functionally graded tapered axial bar with () =  (1 − ); () =   and () =   where  is the taper ratio.For such an FG axial bar Shahba et al [21] computed the nondimensional frequencies  =    ⁄ , ( = 1,2,3) using Finite element method by introducing the concept of basic displacement functions.The pseudospectral method is used in the present work to compute the first three non-dimensional frequencies of the axial bar with  = 16 under fixed-fixed, fixed-free boundary conditions.The results obtained using the present method are compared with those presented in [21] in Table 13.It is observed from the table that the results obtained using PSM is slightly lower than that of [21].The results presented in [21] were calculated utilizing 20 elements to model the bar for a taper ratio of  = 0.1 and higher number of elements for other taper ratios.Moreover, it is observed from literature [21], that the stiffness method generally provides results that are higher than the exact ones.We therefore, feel that the results obtained here are of high precision with fewer collocation points.

Conclusions
The main contribution of the present study consists in applying a pseudospectral formulation for solving the governing differential equation of motion for tapered functionally graded axial bars.The method utilizes a Chebyshev polynomial-to-cosine change of variable that simplifies the computer implementation and achieves good precision.Numerical experiments that include variation in geometrical properties (area of cross-section) and/or material properties (mass density, modulus of elasticity) have been considered and the results obtained with less computational effort exhibits good accuracy and flexibility for longitudinal vibration analysis of tapered axial bars.

FREE
LONGITUDINAL VIBRATIONS OF FUNCTIONALLY GRADED TAPERED AXIAL BARS BY PSEUDOSPECTRAL METHOD.

Table 1 .
A layout of the numerical work carried out in Section 4

Table 5 .
Comparison of non-dimensional frequencies of fixed-fixed non-uniform rods

Table 6 .
Non-dimensional natural frequencies of fixed-fixed rods with

Table 7 .
Non-dimensional natural frequencies of clamped-free (fixed-free) rods with

Table 10 .
Fundamental, first and second overtone frequency values for fixed-fixed rods with () =  exp()

Table 11 .
Non-dimensional frequencies of a fixed-free axially FG bar

Table 12 .
Non-dimensional frequencies of a fixed-fixed axially FG bar

Table 13 .
Non-dimensional frequencies of a FG tapered axial bar