Investigation of the Free Vibrations of Radial Functionally Graded Circular Cylindrical Beams Based on Differential Quadrature Method

Xiaojun Huang，Liaojun Zhang，Renyu Ge，Hanbo Cui and Zhedong Xu

1College of Water Conservancy and Hydropower Engineering，Hohai University，Nanjing，210098，China

2School of Architecture and Civil Engineering，Anhui Polytechnic University，Wuhu，241000，China

ABSTRACT In the current research，an effective differential quadrature method(DQM)has been developed to solve natural frequency and vibration modal functions of circular section beams along radial functional gradient.Based on the high-order theory of transverse vibration of circular cross-section beams，lateral displacement equation was reconstructed neglecting circumferential shear stress.Two equations coupled with deflection and rotation angles were derived based on elastic mechanics theory and further simplified into a constant coefficient differential equation with natural frequency as eigenvalue.Then，differential quadrature method was applied to transform the eigenvalue problem of the derived differential equation into a set of algebraic equation eigenvalue problems.Natural frequencies of the free vibrations of cylindrical beams with circular cross-sections were calculated at one time，and corresponding modal functions were solved together.The obtained numerical results indicated that the natural frequencies of functionally graded(FG)circular cylindrical beams obtained using differential quadrature method agreed with the results reported in related literatures.In addition，influences of varying gradient parameters on the modal shapes of circular cylindrical beams were found to be strongly consistent with previous studies.Numerical results further validated the feasibility and accuracy of the developed differential quadrature method in solving the transverse vibration of FG circular cross-section beams.

KEYWORDS Functionally graded materials;circular cylindrical beams;natural frequency;modal function;differential quadrature method

1 Introduction

Functionally graded materials(FGMs)are a new class of advanced composite materials and are extensively employed in many industries such as civil engineering，mechanical engineering，space engineering，etc.So far，research works have generally focused on the free vibrations and buckling of radially graded rectangular cross-section beams.Currently，the free vibration of functionally graded(FG)beams with varying material properties along the direction of one axis have attracted great attention among researchers.Little attention has been paid in recent years to the free vibrations of FG beams with varying material properties along thickness direction.Xiang et al.[1]investigated the free and forced vibrations of beams by varying material elastic modulus along thickness direction under given initial thermal stress.Chakraborty et al.[2] developed a new beam element for the evaluation of the thermoelastic behavior of FG beams taking into account shear deformation.Ching et al.[3]applied meshless local Petrov-Galerkin method to numerically investigate simply supported FGM beams.Sankar et al.[4-6]obtained the elastic solution of FGM beams with shear deformation under transverse or thermal loadings.

As important structural elements，solid or hollow elastic circular cross-section beams are extensively employed in military，mechanical，aviation and civil engineering applications [7] in different forms such as steel ropes of bridges，aircraft landing gears and circular cross-section beams.According to continuum theory，solid or hollow elastic circular cross-section beams could be assumed as circular cross-section beam and bar models with slight error[8-10].Thus，studying the vibrations of circular cross-section beams is of great theoretical and practical significance.Extensive research has been conducted on circular cylindrical beams and a series of beam theories such as Euler-Bernoulli，Rayleigh and Timoshenko beam theories have been proposed.In the early 1980s，Levinson [11]derived an equation for the calculation of high-order beam theory for bending and free vibrations of rectangular section beams taking into account the influences of inertia term and shear deformation with no need for modified shear coefficient.However，because of the lack of a higher-order theory for FG circular cylindrical beams，few research works have been performed on the mechanical properties of FG beams with circular sections.Loy et al.[12]investigated the free vibrations of FG cylindrical shells.Oh et al.[13] studied the vibration and stability of thin-walls with FG circular cross-sections using a simplified one-dimensional beam model.Law et al.[14]employed wave propagation method to explore the vibrations of cylindrical tubes.Based on nonlocal strain gradient theory，the nonlinear vibrations of beams under different FG distributions were studied by Gao et al.[15-17].A thirdorder shear deformation beam model was developed by Ma et al.[18-20]to investigate the dynamic behaviors of straight hollow cylinders with annular cross-section，in which shear stress vanished on the inner and outer surfaces of the pipe.

Although several methods have been developed for the analysis of the free vibrations of FG beams with gradient along axial，radial or thickness directions，they generally suffer from complex computation methods and tedious solving processes.Differential quadrature method(DQM)was firstly introduced by Bellman and Casti in 1971 and presented the advantages of simple formulation，convenient programming，less computation and high precision [21].Based on the high-order theory of transverse vibrations of circular-section beams.DQM was applied to calculate the free vibration natural frequencies of radial FG circular cylindrical beams at one time and corresponding modal functions could be solved together.According to the basic theory of DQM，the eigenvalue problem of ordinary differential equation for the free vibrations of radial FG circular cylindrical beams was transformed into a set of eigenvalue problems of algebraic equations.Then，the free vibration natural frequencies of radial FG circular cylindrical beams were calculated at one time and corresponding modal functions could be solved together.A circular cylindrical tube was taken as a special case of a bi-layered cylinder where material properties were vanished on the bottom surface.As a special case of an FG beam with circular section characterized by bi-layer material structure，the natural frequencies of thin cylindrical tubes with different boundary conditions can also be calculated by DQM.

2 Theoretical Formulations

A circular cylindrical beam with radial non-homogeneity was considered in this section.For such a structure，a Cartesian coordinate system(x，y，z)was assumed such thatx-axis was considered as the neutral axis of cylindrical beam，positivez-axis was directed upward and perpendicular tox-axis，and corresponding elastic displacement components were written as(u，v，w).Meanwhile，polar cylindrical coordinate system(x，r，θ)and corresponding displacement vector(u，ur，uθ)were also employed.When a circular cross-section cylindrical beam is subjected to transverse loads inxoz-plane，distinctu-axial andw-transverse displacements are obvious.When the beam is restricted to transversely move inxoz-plane，the elastic displacement componentvalongy-direction can be approximately assumed independent fromx.In Fig.1，Ris the radius of circular cylindrical beam andqis distributed force applied along a principal axis of beam corresponding toz-direction.It was concluded from this figure that physical quantities in the two coordinate systems had the following geometric relations，as seen the Eq.(1).

Figure 1:Schematic diagram of cylindrical beams with corresponding coordinates

The relationship between shear strainγxron circumference and radial displacementurwas stated as Eq.(2a).

Substitution ofurin Eq.(1)into Eq.(2a)gave Eq.(2b).

When a circular section was exposed to transverse load，axial and transverse displacements were dominated forms.When only cylindrical beam deformation inxozplane was considered，displacement alongydirection was assumed to be constant along the axis，i.e.，v=v(y，z).In order to make the shear stress on any circumference boundary zero，letγxrtake the following form as Eq.(2c).

whereζ(x)is an undetermined function.To calculateζ(x)，axial displacement was written in the following form as Eq.(2d).

Substitution of Eq.(2d)into Eq.(2c)gave Eq.(2e).

On both sides of the above Eq.(2e)，thezintegral with respect to the variable resulted in Eq.(2f).

The rotation of cross-sectionψwas introduced along centroidal axis as Eq.(2g).

Simultaneous Eqs.(2d)，(2f)and(2g)，elimination functionζ(x)，then the expression of the axial displacement is thus derived as follows Eq.(3):

Based on Saint-Venant principle，the stress components of the beam could be neglected;σyy=σzz=τyz= 0[22].Consequently，τxr= 0，and the condition of freedom on circumferential boundary was automatically satisfied from Eq.(4a).

Equations for other stress componentsσxx，τxyandτxzwere demonstrated as Eq.(4b).

Due to the equilibrium of internal forces，bending momentMand shear forceQat any positionxof beam could be expressed as Eq.(5).

where

whereE(r)andG(r)are Young’s modulus and shear modulus of cylindrical beams depending onr，respectively，whenr= 0，E(r)=EcandG(r)=Gccorrespond to material properties in cylindrical beam axis center andr=R，E(r)=EsandG(r)=Gscorrespond to material properties on cylindrical beam surface.Because cylindrical beam was deflected alongz-direction，noty-direction，under applied bending moment，only two motion equations were involved for flexural deformation inxoz-plane as Eq.(7).

Substitution of Eq.(4)into Eq.(7)yielded Eq.(8).

where

whereρ(r)is beam mass density depending onr，whenr=0，ρ(r)=ρcis material density at cylindrical beam center whenr=R，andρ(r)=ρsis material density on circular cross-section beam surface.Introducing auxiliary functionF(x)to decouple the formulations presented in Eq.(8)into a single governing equation gave deflectionwand rotationψas Eqs.(10a)and(10b).

By substituting Eq.(10)into Eq.(8)，the first equation was naturally satisfied and substituting the second equation into Eq.(8)yielded Eq.(11).

In order to analyze the free vibrations of radial FG circular cylindrical beams，the introduced auxiliary functionFcould take the following form Eq.(12).

whereωis the free vibrations of circular frequency.Substitute Eq.(12)into Eq.(11)supposingq=0 gave the governing equation of radial FG circular cylindrical beams as Eq.(13).

The corresponding boundary conditions of simply-supported(SS)，free-free(FF)，clampedclamped(CC)and clamped-free(CF)beams were Eq.(14).

Thus，the solution of the natural frequency of FG cylindrical beams converted to the calculation of Eq.(13)under the boundary condition Eq.(14).For the convenience of programming，was assumed，the governing equation of the free vibrations of circular cross-section beams in Eq.(13)was rewritten into a dimensionless equation as Eq.(15).

To express generality，a CC beam was taken as an example.From Eqs.(10)and(12)，the followings could be inferred Eqs.(16)and(17).

3 Calculation of the Natural Frequencies of Circular Section Beams by DQM

The vibration modal functionf(ξ)of the beams was considered to be differentiable along beam length interval [0，1] and interval [0，1] was divided intondiscrete elements，assumingξ0= 0 andξn=1.In this research，non-uniform grid in Eq.(18a)(root of Chebyshev polynomial)was applied to represent node coordinatesξiand equal size uniform grid in Eq.(18b)represented node coordinatesξi:

In the interval[0，1]，the function and its derivatives at each node were represented by a weighted linear sum of the values of discreten+ 1 node functions;the functionf(ξ)was described using Lagrange interpolation function as Eq.(19).

wheref(ξ)is Lagrange polynomial and is stated as Eq.(20).

From Eq.(19)，the first derivative of functionf(ξ)was obtained as Eq.(21).

It was divided intonsegments along beam length interval [0，1].Therefore，Eq.(21)was transformed as Eq.(22).

Eq.(22)was reduced to vector form as Eq.(23).

4 Numerical Examples and Discussion

4.1 Free Vibrations of Radial FG Cylindrical Beams

To investigate the effects of gradient parameters on the natural frequencies of circular cross-section beams，for simplification，the material properties of circular cylindrical beams were assumed to follow power-law functions as Eq.(36).

wheremis gradient parameter reflecting volume fraction change effect on material properties.In the following calculations(presented in Tables 1-4)，a circular cylindrical beam with the radius to length ratio ofR/L= 0.1 was applied andnrefers to the number of discrete elements along beam length interval(0，1).Zirconia and aluminum were chosen as outer and the inner surface materials，respectively，and the material properties of which were[24]:

Al:Ec=70GPa，Gc=27GPa，=2702kg/m3

ZrO2:Es=200GPa，Gs=77GPa，ρs=5700kg/m3

Table 1:Relationships between the dimensionless natural frequencyof a circular FF along radial functional gradient(n=20)

Table 1:Relationships between the dimensionless natural frequencyof a circular FF along radial functional gradient(n=20)

m The stationing way The first The second The third The fourth The fifth The sixth 1 Uniform grid 0.1027264 0.1579136 0.4048480 0.7746017 1.223927 1.244700 Non-uniform grid 0.0327436 0.1749075 0.4240046 0.7231154 1.041864 1.379669 Literature[25] 0.0325 0.1759 0.4238 0.7165 1.0358 1.6936 5 Uniform grid 0.1824912 0.3995242 0.6674286 1.143029 1.349056 1.782809 Non-uniform grid 0.0321374 0.1716784 0.4104955 0.690696 0.996193 1.609794 Literature[25] 0.0320 0.1710 0.4095 0.6888 0.9919 1.6085?

Table 2:Relationships between dimensionless natural frequencyand the number n of a circular SS along radial functional gradient(m=5)

Table 2:Relationships between dimensionless natural frequencyand the number n of a circular SS along radial functional gradient(m=5)

Number of discrete elements n The first The second The third The fourth The fifth The sixth n=6 0.08826002 0.3300874 1.019301 1.800898 1.906845 2.206325 n=8 0.08836415 0.3104719 0.6096987 0.9886772 1.539500 1.800898 n=12 0.08834945 0.3094812 0.5980237 0.9222708 1.280307 1.693504

Table 2(continued)Number of discrete elements n The first The second The third The fourth The fifth The sixth n=16 0.08834393 0.3091891 0.5952669 0.9089863 1.235276 1.570097 n=20 0.08834288 0.3091361 0.5947968 0.9068756 1.228627 1.553210 Literature[25]0.0885 0.3095 0.5954 0.9075 1.2291 1.5531

Table 3:Dimensionless natural frequencyof a circular CF along radial functional gradient(n=20)

Table 3:Dimensionless natural frequencyof a circular CF along radial functional gradient(n=20)

Note:“/”in the table indicates that calculated values were not given in reference[25].The following is similar.

m The first The second The third The fourth The fifth The sixth m=0.2 Literature[25]0.1921 0.4616 0.7837 1.1234 1.4657 1.7929 DQM solution 0.1919242 0.4610495 0.7827532 1.122112 1.464038 1.790692 m=1 Literature[25]0.1943 0.4626 0.7795 1.1102 1.4404 1.7478 DQM solution 0.1940521 0.4619524 0.7785301 1.108935 1.438912 1.746147 m=5 Literature[25]0.1902 0.4489 0.7514 1.0639 1.3705 1.6520 DQM solution 0.1899352 0.4483912 0.7506514 1.062987 1.372359 1.650890 m=7 Literature[25]//////DQM solution 0.1871886 0.4423731 0.7411917 1.050353 1.356883 1.634063

Table 4:Dimensionless natural frequencyof a circular CC along radial functional gradient(n=20)

Table 4:Dimensionless natural frequencyof a circular CC along radial functional gradient(n=20)

m The first The second The third The fourth The fifth The sixth m=0.2 Literature[25] 0.1681 0.3915 0.6674 0.9772 1.3106 1.6582 DQM solution 0.1724183 0.4029486 0.6839595 0.9945545 1.325356 1.668221 m=1 Literature[25] 0.1684 0.3892 0.6606 0.9642 1.2895 1.6273 DQM solution 0.1721427 0.3986689 0.6738341 0.9774612 1.300228 1.633864 m=5 Literature[25] 0.1631 0.3742 0.6331 0.9219 1.2301 1.5492

Note:“/”in the table indicates that calculated values were not given in reference[25].The following is similar.

Table 1 summarizes the first six dimensionless natural frequencies of CF in radial FG calculated by DQM(n= 20，m= 1 andm= 5).It was found that the numerical results obtained from non-uniform grid distribution method were identical to those reported in literature [25];however，large error or even distortions were observed between the numerical results obtained from uniform and equal step size distribution method compared to those reported in literature[25].The numerical results given in Table 1 show that Eq.(18b)could be used for discrete elements along beam length.The numerical results calculated by DQM were unstable or even distorted，while higher accuracy could be achieved by Eq.(18a).Table 2 summarizes the first six dimensionless natural frequencies of CC in radial functional gradient calculated by DQM(m=5)，where Eq.(18a)is applied as the distribution mode of discrete nodes along beam length.As was seen in the table，the numerical results obtained for the first six natural frequencies were close to those reported in reference[25]except for the first and second natural frequencies whenn=6 and 8.Other numerical results obtained for the fourth to sixth natural frequencies had large errors and even distortions.However，the numerical results obtained for the first six natural frequencies were close to those reported in reference[25]through increasingn.In other words，with the increase of the number of discrete elementsn，numerical results showed higher accuracy.

Tables 3 and 4 present the first six dimensionless natural frequencies of FF and CC in radial functional gradient calculated by DQM，respectively(n= 20，m= 0.2，1，5 andm= 7).As seen from Tables 3 and 4，numerical results obtained from non-uniform grid distribution method showed excellent consistency with those reported in reference[25]，which indicated that the high-order theory of transverse vibrations of circular cross-section beams were correct and DQM was effective and accurate for the calculation of the natural frequencies of radial functional gradient circular crosssection beams.The calculation results presented in Tables 3 and 4 also showed that the dimensionless natural frequencies of beams were decreased with the increase of gradient materialm.According to Eq.(36)，this was because the average stiffness of beams was monotonically decreased with the increase of gradient materialm.

Fig.2 presents the first 15 non-dimensional frequency parameters of five radius-to-length ratios calculated by DQM(n=20，m= 3).For smaller values ofR/L，i.e.，R/L=0.002 and 0.005，relations between non-dimensional frequency parameter and mode number were almost linear.The relations became non-linear with the increase of radius-to-length ratios.

Figure 2:Non-dimensional frequency parameter Ω of beams under different boundary conditions

In this research，the natural frequencies of the free vibrations of radial FG circular cylindrical beams were calculated using DQM as well as the first and second modal curves of various gradient parameters under different boundary conditions were obtained.Figs.3 to 6 indicate that the variations of gradient parameters had almost negligible effects on the modal curves of CF and SS，but significant effects on the first and second modal curves of FF and CC.The results were consistent with FG rectangular CF along thickness on effect of material parameters on the modal behavior[26].

Figure 3:Modal diagram of radial FG circular CF(a)First modal diagram(b)Second modal diagram

Figure 4:Modal diagram of radial FG circular SS(a)First modal diagram(b)Second modal diagram

Figure 5:Modal diagram of radial FG circular FF(a)First modal diagram(b)Second modal diagram

Figure 6:Modal diagram of radial FG circular CC(a)First modal diagram(b)Second modal diagram

4.2 Free Vibrations of Isotropic Homogeneous Material Beams

Based on Eq.(36)，whenm= 0，supposingE(r)=Es=E，G(r)=Gs=G，ρ(r)=ρs=ρ，radial FG circular cylindrical beams were simplified to isotropic homogeneous materials beams，whereE，Gandρare，respectively，the elastic modulus，shear elastic modulus and density of the uniform material.SinceR/L= 0.1 andμ= 0.3，based on high-order theory of transverse vibration of circular cross-section beams，the first eight dimensionless natural frequenciesfor the circular cross-section beams of isotropic homogeneous materials were calculated by DQM，wheren= 20 and the distribution mode of discrete elements was non-uniform grid.The obtained results are listed in Tables 5 and 6 and compared to those reported in reference [27].It was found that numerical results obtained from Euler-Bernoulli beam theory were generally too high because shear deformation was not considered.However，the results showed excellent consistency with those obtained by DQM in literature [27] which were very close to the results obtained from Timoshenko beam theory.It was proved that DQM was feasible and accurate in calculating the natural frequencies of radial functional gradient circular cross-section beams with homogeneous materials.Unlike Euler-Bernoulli and Timoshenko beam theories，a higher-order theoretical beam model was developed which did not require shear correction factor and planar cross-section assumption after deformation.

Table 5:The dimensionless natural frequencyof CF

Table 5:The dimensionless natural frequencyof CF

The frequency Euler-Bernoulli beam[27]Timoshenko beam[27]Literature[27] DQM solution 1 0.02835 0.02772 0.02755 0.0275142 2 0.1776 0.1545 0.1522 0.1518469 3 0.4974 0.3774 0.3717 0.3702360 4 0.9747 0.6400 0.6312 0.6291503 5 1.6113 0.9327 0.9143 0.9117888 6 2.04070 1.2156 1.2075 1.2052101 7 3.3619 1.5079 1.5017 1.5009302 8 4.4759 1.7808 1.7792 1.7806117

Table 6:The dimensionless natural frequencyof CC

Table 6:The dimensionless natural frequencyof CC

The frequency Euler-Bernoulli beam[27]Timoshenko beam[27]Literature[27] DQM solution 1 0.1804 0.1524 0.1483 0.1479448 2 0.4972 0.3596 0.3485 0.3470239 3 0.9747 0.6098 0.5932 0.5901438 4 1.6113 0.8817 0.8634 0.8591617 5 2.4070 1.1673 1.1502 1.1460887 6 3.3619 1.4601 1.4467 1.4441774 7 4.4759 1.7565 1.7488 1.7494903 8

4.3 Free Vibrations of Isotropic Cylindrical Tubes

Isotropic cylindrical tubes could be regarded as special FG materials with bi-layer material structures.The inner of bi-layer had vanishing material properties(i.e.，material constant was zero)，while the outer of bi-layer was homogeneous.Therefore，the material properties of cylindrical tubes could be rewritten as Eq.(37).

wherehis tube thickness andRis cylindrical tube section radius.Substituting Eq.(37)into Eq.(13)gave the natural frequency governing equation of a FG circular cylindrical tube.R/L= 0.05，h/R=0.002，and material Poisson’s ratioμ=0.3 were assumed.The distribution mode of discrete elements was non-uniform grid and the number of discrete elements wasn= 20.The dimensionless natural frequenciesof cylindrical tubes under two boundary conditions could be calculated.The results obtained from DQM and those reported in reference [12] are listed in Tables 7 and 8，respectively.As seen from these tables，our results were completely consistent with those in existing literatures，the accuracy of DQM in calculating the natural frequencies of the free vibrations of FG circular cylindrical tubes were validated.

Table 7:The dimensionless natural frequencyof SS cylindrical tubes

Table 7:The dimensionless natural frequencyof SS cylindrical tubes

The frequency order DQM solution Literature[27] Literature[28] Literature[29]1 0.01603352 0.0160 0.0161 0.0168 2 0.05845778 0.0583 / /3 0.11672024 0.1166 / /4 0.18304303 0.1827 / /

Table 8:Dimensionless natural frequencyof clamped cylindrical tubes

Table 8:Dimensionless natural frequencyof clamped cylindrical tubes

The frequency order DQM solution Literature[27] Literature[30] Literature[31]1 0.03435762 0.0343 0.0344 0.0349 2 0.08407543 0.0839 0.0848 0.0874 3 0.14421041 0.1445 / /4 0.20891862 0.2097 / /

5 Conclusion

In this research，based on the high-order theory of transverse vibration of circular cross-section beams，the calculation equation of natural frequency was converted into a differential equation with natural frequency as eigenvalue.From the above analyses and derivations，the differential equation was transformed into a standard generalized algebraic eigenvalue equation by differential quadrature method theory and the obtained algebraic equations were solved byQRmethod.Free vibration natural frequencies of a cylindrical beam with circular cross-section was calculated and corresponding modal function curve was achieved.The following main conclusions were drawn:

(1)In this paper，the natural frequencies of the transverse vibrations of FG circular cross-section beams under different boundary conditions were calculated by DQM.The obtained results were consistent with those reported related literature，which indicated that DQM was highly efficient and precise in calculating the natural frequencies of FG circular cross-section beams.

(2)Numerical results showed that uniform grid was used for discrete elements along the direction of beam length，the numerical results calculated by DQM were unstable or even distorted，but non-uniform grid was applied to obtain calculation results which were consistent with the results of reported in the existing literatures.

(3)The circular cylindrical beams of homogeneous materials and cylindrical tubes can be regarded as the special cases of the circular cross-section beams of FGMs.In this paper，based on highorder theory for the transverse vibrations of radial FG circular cylindrical beams，the natural frequencies of circular cross-section beams with homogeneous materials and cylindrical tubes were calculated by DQM.The consistence of numerical results with those reported in existing literature further verified the computational accuracy of DQM.

(4)The natural frequencies of circular cross-section beams with FG materials and first and second modal curves of different gradient parameters under various boundary conditions could be calculated by DQM.The obtained numerical results showed almost negligible effect on the modals of CF and SS when changing gradient parameters，but these effects were especially obvious for FF and CC.

Acknowledgement:The authors would like to thank the anonymous reviewers for carefully reading the article and this research was financially supported by the National key Research and Development Plan of Ministry of Science and Technology of the People’s Republic of China(2017YFC0404903).

Funding Statement:The National key Research and Development Plan of Ministry of Science and Technology of the People’s Republic of China(2017YFC0404903).

Conflicts of Interest:The authors declare that they have no conflicts of interest to report regarding the present study.