Static Behaviour of Functionally Graded Rotating Cantilever Beams Using B-Spline Collocation Technique

The present work reports the static behaviour of functionally graded rotating beam based on Timoshenko beam theory, which includes the effect of shear deformation. The principle of virtual displacement is applied to derive a governing equation for the functionally graded (FG) rotating beam, considering the centrifugal stiffening effect. The B-spline collocation technique is employed to solve the differential equation and material properties are function of power law distribution. The effect of material power index and rotational speed on the static characteristics of functionally graded rotating cantilever beams has been investigated for two different slenderness ratios. The obtained results demonstrate the beam deflection and the beam normal and shear stresses for the power gradient and rotational speed values. The results help to conclude that material gradation can be utilized to improve the functioning of the rotating structures such as helicopter rotor


INTRODUCTION
Functionally graded materials are a new class under composites where the material constituent is tailored in a way to give a gradual variation of material properties, thus, avoiding failure due to interfacial stress concentration, delamination, and cracking of the matrix as commonly found in composites. Due to superior properties, a functionally graded material finds wide range of applications in structures. Rotating beams are often used to model structural elements such as helicopter blades, airplane propellers, windmill rotors, and other blades of turbomachinery. The authors in paper [1] have described the unpredictable failure of the composite material of helicopter blade under influence of various force. The authors [2], have studied helicopter blades to investigate damage growth in material that resulted in matrix cracking, debonding/delamination, and fibre breakage using the finite element method.
In their paper [3], authors have conducted a case study of a military helicopter crash that indicated failure in blade due to the fatigue mechanism caused by an incision made during spar maintenance. Therefore, it becomes necessary to evaluate the material of the rotating structure to avoid accumulation. of stresses, cracking and debonding of material caused due to damage of the component.
The study of Euler-Bernoulli beams was done by the rational interpolation functions to solve governing differential equation for Euler-Bernoulli rotating beam to analyse static and dynamic behaviour for homo-genous material. A modified Gazelle helicopter blade has been analysed to study the natural frequency and mode shape by [5] using Lanczos method to extract decomposition matrix and eigenvalue. [6]also developed new beam element to analysis the behaviour of functionally graded material for a beam structure, where static and free vibration analysis of functionally graded material (FGM) beam under thermal loading was done. Early work reported by [7] and [8] focussed on formulating a decoupled equation of motion based on Timoshenko beam theory. Where using a power series solution to evaluate the dynamic characteristic of a rotating isotropic beam.
The work by [9] has considered the effect of Coriolis force on natural frequency of rotating beam, [10] has investigated flap-wise bending vibration along with coupled lag-wise bending and axial vibration to study the effect of Coriolis force, slenderness ratio at high angular velocity using power series solution. [11] approached through the finite element method to analysis natural frequency, time responses, and distributions of deformation and stresses for functionally graded rotating beam. Similarly, [12] also used the finite element method to study a nonlinear rotating beam model which considers the effect of geometric stiffness, Coriolis force applied to isotropic, composite, and functionally graded beam. Early studies on rotating beam have ignored the term ρIΩ 2 θ, [13] has identified the significance of this term at higher speed. In paper [14] the authors have studied free vibration of the Timoshenko beam using a differential transform method to examine various geometrical parameters on vibrational characteristics of rotating FG beam. The study by [15] analysed free vibration of a tapered rotating beam using the frobenius method of series solution. In another article by [16], introduced NURBS-based finite element formulation to investigate flap-wise vibration. [17], in their work, have analysed stress distribution and deformation of a functionally graded beam using the B-spline collocation method. In this paper, [18] the authors have investigated the effect of tapering on a rotating beam. [19] solved the differential equation using spectral finite element method to analysis mode shape and modal frequency for uniformly tapered rotating beam.
In their widely acclaimed work by [20], the authors have derived a single fourth order governing partial differential equation for static and dynamic analysis of non-rotating beam. Euler-Bernoulli and Timoshenko beam theory was used to study deflection and stress distribution in functionally graded material. Wave propagation and free vibration was studied for dynamic analysis. Results reported by [21] suggest that the derived flexural stiffness function of a rotating beam whose eigen pair is equivalent to a non-rotating beam can be employed to check the code and develop the actual beam to counter the stiffening effect due to the centrifugal force. [22] as compared two methods, that is, the differential transformation method (DTM) and the differential quadrature element method, to analyse the free vibration of FG tapered rotating beam. For various boundary conditions, [23] studied free vibration of axially functionally graded nonuniform cross-section beam using Fredholm integral equations. In paper [24] the authors have used the ANSYS workbench to study the static and dynamic characteristics of functionally graded material. A recent analysis by [25] highlights stress and vibration analysis for a functionally graded rotating beam where the effect of material gradation and rotation speed has been studied. [26] has investigated bending deformation and vibration behaviour of functionally graded rotating beam using novel floating frame reference (FFR) formulation. [27] has investigated non-linear behaviour of rotating beam at large amplitude where hardening /softening effect at various rotating speed has been analysed. Axially graded nonuniform beam was investigated by [28] for wrinkle type of elastic foundation, where using energy method dynamic analysis was done on FG beam. Using the Rayleigh-Ritz method coupled with Kane's method [29] has studied vibration characteristics of rotating pretwisted tapered FG beam for various geometric and material parameters.
An experimental investigation of helicopter rotor blade was done in paper [30], where vibratory testing was done to test the material properties of composite helicopter blades. In their another paper [31] the author has done fatigue testing of helicopter rotor blade for composite laminated material. Static and dynamic experimental tests were performed by [32] for the rotor tail blade made up of composite laminated materials. The paper [33] discuss about the design of rotor blade of helicopter to satisfy the cross-sectional properties of composite material.
Other structural models where rotating effect is considered for static analysis is disc where authors like [34] has studied stresses and deformation of rotating FG disc for various boundary condition using finite element method. [35] has reported limit elastic speed analysis of rotating FG beam based on von mises criteria. [36] has also done analysis radial, tangential stress and strain for FG rotating disc. [37] and [38] has derived equation of motion for axial, flap wise and chordwise motion using Hamilton's principle and has analysed natural frequencies and time response for isotropic and FGM beam. [39] has studied vibrational behaviour and stress of perforated rotating beam using finite element analysis. The study shows the effect of perforated configuration and rotating speed on stress distribution and vibration response of the model. [40] has presented static response of porous rotating thick truncated cone for functionally graded material.
The present work lays the foundation of the mathematical formulation of rotating functionally graded beam using the B-spline collocation method and discusses the static behaviour of such structures under the effect of different parameters. The effect of rotational speed on the deflection and stress distribution on the FG beam for various values of power index has been investigated for two different slenderness ratios. Based on the literature review, it is observed that there is hardly any literature present that directly deals with static analysis of functionally graded rotating beam. The structure of the paper is as follows: The subsequent section, i.e., section 2 presents the formulation of rotating beam using principle of virtual displacement where the property of materials used is varied in vertical direction. This section also introduces to the numerical method used to solve the governing equation i.e., Bspline collocation method. Verification of proposed numerical method is done by validating the preliminary results with those obtained using ANSYS and is reported in section 3. Section 4 reports and discusses the obtained results. Section 5 concludes and summarises the work reported in this article.

METHODOLOGY
A functionally graded cantilever beam of length L, width b, height h rotating at Ω rad/sec is considered. The beam is fixed to a hub of radius, r and the hub is free to rotate about the vertical axis, as shown in figure 1. Material grading of the beam is assumed to vary in the transverse direction across the height of the beam according to the power law (1). The FG beam is made of metal on the top surface and ceramic on the bottom surface. The effective properties generalized as P, i.e., young's modulus (E), density(ρ), Shear modulus(G) and Poisson's ratio (µ) of the functionally graded beam, are evaluated using the power law function as follows: where k s is the shear correction factor [20], λ is the power index, T and B denote the top and bottom surfaces of the beam, (where the top and bottom surface comprises of metal and ceramic with young's modulus  The governing equations are derived using the principal of virtual work in the framework of Timoshenko beam theory to incorporate the effect of shear deformation. Displacement variables (u, v, w, and θ) are considered in the x, y, z-direction, and rotation of the cross section, respectively. The displacement field is taken as: Here, u 0 and θ denotes axial displacement and rotation of the midplane. Therefore, normal strain ε xx and shear deformation γ xz for moderate deformation is given by Using the principle of virtual work displacement, one can state that work done by the actual force is zero if and only if the body is in equilibrium. Therefore, the sum of the virtual work done by the internal and external forces is equal to zero.
 , is the strain energy and is given by: Here, , is the external work done by centrifugal force on the beam.
Centrifugal force N R is given by.
In the absence of applied axial force xx N , substituting (6), (9), (11) in (4) , the governing equation for rotating beam is given by.
To simplify the calculation, (12) is converted into a single equation with a single unknown. Therefore, from [20]: where F is a new auxiliary function of length dimension. On substituting the above equation for 'w' and 'θ', a single fourth order equation is obtained. 4 2 In the absence of axial force, bending moment and shear force in the beam is given by: The authors, in their paper [20], have formulated the governing equations for a non-rotating beam and simplified them to a single equation whose solution is approached using direct integration method. However, in the present case it is difficult to obtain solutions using direct techniques available for non-linear cases. Hence the solutions to the governing equation for rotating FG beam is approached using numerical techniques. On solving (15) by using an appropriate numerical method, the auxiliary function 'F' can be determined which can be used to find other depended variables. In the present study solution of (15) is obtained using the collocation technique with b-spline basis as approximating function.
The B-spline collocation method provides a closedform, piecewise continuous solution where the differential equation is satisfied at a finite number of points called collocation points. The parametric coordinate 't' is defined in terms of order of the resulting polynomial, also known as the knot vector. It is assumed that the knot vector is of open type given by: where, n+1 are the number of control points and 'k' is the order of the polynomial spline. Once the knot vector is selected, B-spline basis function N i,k (t) is defined by a recursive relationship that is known as the Cox-de Boor recursion formula.
The resulting B-spline curve is defined using the (n+1) control points B 0 , B 1 , B 2 .... B n and the basis function as: The above function 'F' is substituted in the governing equation as approximating polynomial such that it satisfies the given boundary conditions. As the number of unknowns in approximating polynomial ex-ceeds the number of boundary conditions, hence we need to force the solutions at discrete points called collocation points. In this study, the collocation points or abscissa control points are selected using the Grevillea abscissa approach. The position vector of the control point 'B i ' used to define the B-spline curve with the parametric coordinate 't' is calculated as follows: (20) In the current analysis, a polynomial of degree 4 or above will be suitable to solve (15). A sixth-order polynomial (fifth degree) has been selected using openuniform knot vector with no intermediate points to act as an approximating polynomial. Considering the case of material and geometric nonlinearity related to problem and to investigate wider levels of research, a higher order of b-spline functions is taken. A MATLAB code is developed indigenously to numerically investigate the problem of rotating FG beams.

VERIFICATION OF METHOD WITH ANSYS
To verify the correctness and accuracy of the present method, the results obtained by the B-spline collocation method are compared with the results obtained through commercial finite element package of ANSYS. A rotating cantilever beam of uniform cross-section with length of beam = 0.5 m and height = 0.125 is considered. The materials used for the FG beam are of ceramic metal composition, which specifically have the following properties: Metal part is aluminium with young's modulus= 70GPa, density = 2700Kg/m 3 , poison's ratio =0.3 while the ceramic material is silicon carbide with young's modulus = 413.68GPa, density = 3100kg/m 3 and poison's ratio =0.3. The beam rotates at a uniform speed of Ω=100 rad/sec.
The response of rotating beam is investigated using the MATLAB code developed as per the mathematical formulation discussed in previous section and is verified using the ANSYS software. Comparison between the obtained result for deflection and normal stresses in Figure 3 and Table 1 shows that a good agreement exists between results obtained using MATLAB and ANSYS.

RESULTS AND DISCUSSION
In the present study, the static deflection characteristic and stress distribution in a rotating functionally graded beam have been investigated. The rotating beam is assumed to be composed of metal and ceramic; the material properties smoothly vary according to the power law across the height of the beam; properties of material have been listed in Table 2. Two cases have been taken with slenderness ratio L/h=10 and L/h=20 of constant crosssection i.e., height and width =0.125 m. The proposed computational model is used to inves-tigate the effect of centrifugal force on stress distri-bution as well as deflection throughout the beam under self-body weight.

Effect of Rotation
The effect of rotation on functionally graded beams at different speeds and at specific power gradient has been investigated. The centrifugal force generated due to rotation shows different behavior for the slenderness ratio L/h=10 and L/h=20 on static deflection and stresses. Transverse deflection for different rotating speeds at a specific power gradient are shown in figure  4 and 5 at slenderness-ratio (L/h) 10, 20 respectively. From the figure it can be depicted that, under nonrotating conditions, i.e., Ω=0 rad/sec, the beam deflects under self-weight and as soon as rotation is given, this deflection starts to increase. The deflection increases further as the rotational speed increases. In the case of L/h = 10 ( Figure 4) typical behavior in deflection is observed where higher speeds give higher deflections. The gap amongst deflection at different speeds continues to decrease with increasing power gradient, which is the lowest for ceramics and at λ=5 whereas for metal the gap is the highest. Stiffness of the structure depends on the two factors, firstly on stiffness of the material and secondly on the geometry of the structure. Higher value of power index increases the stiffness whereas higher value of slenderness ratio reduces the stiffness. For the slenderness ratio, L/h = 20, Figure 5 shows an inconsistent behavior, where the beam deflects downwards till Ω=100 rad/sec and thereafter the beam deflects upward. This is due to the increase in length of the beam, which in turn decreases the flexural rigidity. This flexural rigidity increases with increase in the material gradient as it can be seen from the figure that the beam axis is moving toward the neutral axis with increasing power gradient. The comparison of the induced normal stresses with the slenderness ratio (L/h) 10, 20 is presented in Figures  6 and 7, showing the normal stresses for different speeds for specific material gradation. For L/h = 10, the graph shows smooth variation of normal stresses along the height of the beam for functionally graded material, where there is an increase in the normal stress with increase in the rotational speed at all the value of power gradient. However, the effect of an increase in centrifugal force with rotational speed is minimal in the case of L/h=10. For L/h=20, a typical variation of normal stresses is observed due to decrease in flexural rigidity. The normal stresses for Ω = 0, 50 rad/sec shows the tensile stress on the upper surface and the compressive stress on the lower surface. At Ω = 100 rad/sec a sudden change in the direction of normal stress is observed due to fluctuation of beam in upward direction. On further increase in rotational speed, the normal stress is reported to increase.
Similarly, for shear stresses, the effect of rotation is least for of L/h = 10 as shown in Figure 8. The maximum shear stress increases with increase in the rotational speed, and this increase in the maximum stresses is reduced with increase in power gradient. For a slenderness ratio of L/h = 20, Figure 9 shows that shear stress increases with increase in rotational speed. However, the fluctuation in the shear stress is observed at the rotational speed of Ω = 100, 200 rad/sec. For metal, the maximum shear stress is quite high for Ω = 200 rad/sec. Hence the alteration in shear stress is relatively high in case of L/h =20 and least in case of L/h = 10.

Effect of Power Gradation
Figures 10 and 11 present the influence of the power gradient on the deflection and stresses of the rotating beam as the centrifugal stiffness increases with the increase in the rotation speed. Figure 10 shows that in the slenderness ratio L/h = 10, the increasing material gradient increases the stiffness of the material; as a result, there is a decrease in the deflection. Maximum deflection is for full metal at non-rotating and rotating condition, while least for ceramic and λ = 5. For a nonrotating beam under self-load, the beam deflection keeps decreasing with increasing power gradient.
As soon as the beam starts rotating due to centrifugal force, the deflection in the beam increases. At 100 rad / sec, the metal with the lowest value for flexural rigidity is overshadowed by centrifugal stiffness; with further increase in rotational speed, the beam with a lower value of the material gradient, i.e., λ =0.2 shows a higher deflection. At higher speeds, the stiffness due to the centrifugal force overcomes the flexural stiffness due to the material gradient and the beam deflects in the upward direction. For pure metal due to the least flexible stiffness, the material deflects first in the upward direc-tion followed by the beam at λ=0.2, 0.5, 1, 5 and ceramic with increasing rotational speed.  Figure 12 shows a change in normal stress with rotational speed where the stress varies linearly for isotropic material and non-linearly for functionally graded material. For rotating as well as non-rotating beam the tensile stress on top surface and compressive stress on bottom surface decreases with increase in power gradient for L/h =10. Although for the slenderness ratio of L/h = 20 (Figure 13), there is a decrease in normal stress, however, the distribution fluctuates with an increase in the power gradient at higher rotational speed. At Ω =200 rad/sec, minimum normal stresses are ob-tained at λ =0.2 while for lower speed minimum normal stresses is obtained at λ = 5.0. After λ = 5, the neutral axis moves back towards the centroidal axis. However, there is an increase in shear stress with increase in power gradient for non-rotating beam at λ =full metal, 0.2, 0.5, 1, 5. On further increase in power gradient, the maximum value of shear stress decreases. In the case of a rotating beam, the minimum value of maximum shear stress is obtained at power gradient λ =0.2. Maximum value of shear stresses for L/h=20 is higher as compared to L/h =10. An enormous increase in the shear stress is observed for full metal case at the rotational speed of 200 rad/sec, while for other power index, not much increase in the value of shear stress is reported.

Effect of rotation along the length of the beam
Normal stresses along the beam length for slenderness ratio of L/h =10, 20 have been compared in Figure 16 for λ= 0.2, 1, 5 at 200 rad/sec. In both cases the normal stress at the free end gives linear vertical curve whose value is zero across the height of the beam. As one moves towards the fixed end, the curvature of nonlinear curves increases along the length at L/4, L/2, 3L/4, L distance from fixed end of the beam.
The beam deflection for L/h =20 in figure 16 (c), (d) is in opposite direction to that of L/h=10 in figure 16 (a), (b) at 200 rad / sec, therefore, the direction of the deflection curve for normal stresses is also in the opposite direction. Similar is the case for shear stresses shown in Figure 17, where due to fluctuation of beam in upward direction at 200 rad/sec, there is a shift in the shear stress distribution for slenderness ratio L/h= 20 as compared to the beam with L/h=1.

CONCLUSION
The B-spline collocation method is a numerical technique suitable for use with standard geometries, like rectangles, that provides piecewise-continuous, closed-form solution. A MATLAB code was developed to solve the mathematical model of a rotating functionally graded beam. The solution obtained is verified using the ANSYS software. In this work, the static behaviour of a transversally functionally graded rotating beam is presented. This study appears to be the first to analyse the effect of functional grading on deflection and the stresses developed when the beam is under the influence of centrifugal force such as those found in the applications of rotor blades and slender rotating structures. The results also highlights the crucial part played by slenderness ratio in stability of the rotating beam and effect of functional grading on slenderness ratio.