Document Type: Research Paper
Authors
Department of Solid Mechanics, Faculty of Mechanical Engineering, University of Kashan, Kashan, I.R. Iran
Abstract
Keywords
Thermoelastic Analysis of a Functionally Graded Simple Blade Using FirstOrder Shear Deformation Theory
M.M. Hosseini Mirzaei, M. Arefi^{*}, A. Loghman
Department of Solid Mechanics, Faculty of Mechanical Engineering, University of Kashan, Kashan, I.R. Iran
KEYWORDS 

ABSTRACT 
Simple blade Functionally graded materials Firstorder shear deformation theory (FSDT)

In this article, the thermoelastic behavior of a functionally graded simple blade subjected to the mechanical and thermal loadings is presented, applying a semianalytical method and a variable thickness cantilever beam model. A specific temperature gradient is employed between the root and the edges of the beam. It is assumed that the mechanical and thermal properties are longitudinal direction dependent pursuant to volume percent of reinforcement. The approach is composed of several steps, including adoption of firstorder shear deformation theory, applying beam division accompanying the longitudinal direction, imposing global boundary conditions, and deliberating the continuity conditions. As a result, longitudinal and transverse displacements, and consequently longitudinal, shear and effective stresses are acquired. The analysis is performed for three different distributions of reinforcement particles and pure matrix. Minimum effective and shear stresses distribution belong to the blade with 0% reinforcement at root and 40% reinforcement at tip surface. It has also been discovered that application of reinforcement particles have reasonable effect on the longitudinal and transverse deflections. 
Nowadays, the application of functionally graded materials in the structures, devices, and parts of components is developing. Most of these materials have desirable resistance to heat and consequently are employed in heatexposed components. On the other hand, beams are one of the most famous and wellknown elements whose analysis can be beneficial in different ways. In this study, a semianalytical solution is presented for a functionally graded simple blade. A simple blade modelled by variable thickness cantilever beam is contemplated in order to gain the idea of deformations and stresses. Subsequently, applying the present method, thermoelasticanalysis of the functionally graded blade made of ZrO2/Ti6Al4V is proposed, and the results are extracted and compared with each other.
Kapania and Raciti reviewed the developments in the analysis of laminated beams and plates [1]. In a special section of the review of the previous work, they discussed the various sheardeformation theories for plates and beams and pointed to several methods of analysis pursuant to the finite element method. They implied to one of the finite element methods that are based on the highorder shear deformation theory. The application of shear deformation theory was also deliberated by researchers in order to analyze beams, especially composite beams. In 1998, Shi et al. developed a model in consonance to finite element for the analysis of composite beam elements. Interestingly, the basis of this model is pursuant to the shear deformation theory [2]. In 2001, Sankar presented an elasticity solution for functionally graded beams [3]. In that research, mechanical properties such as Young’s modulus and shear modulus are a function of the thickness of the beam, and the analytical method does not respond to the application of functionally graded materials that change their properties along the beam. Sankar and Tzeng inspected the thermal stresses in functionally graded material beams [4]. Here, properties are also variable in the thickness of beam. In another study, Chakrabortya et al. examined the threelayer beam by the method of finite element relations [5]. This beam consists of two layers of the dissimilar material at the bottom and lowers and a core of a functionally graded material. In this study, the properties are the function of the core thickness. Nirmala et al. indicated an analytical solution for thermoelastic stresses in a composite beam with functionally graded material core [6]. Moreover, Kadoli et al. inspected a static analysis of functionally graded beams applying higher order shear deformation theory [7]. In one part of their research, they compared the results of the use of firstorder and higherorder shear deformation theory and concluded that the dissimilarity in results was very small and insignificant. In another similar study in 2008, X.F. Li provided a formulation for the static and dynamic analysis of functionally graded Timoshenko and Euler–Bernoulli beams [8]. In this research, the natural frequencies and mode shapes of a simply supported beam were obtained and reported.
Investigations on the thermal buckling of functionally graded beams were performed by Kiani and Eslami, where the gradation in properties was through the thickness [9]. This research surveyed different boundary conditions under various heat distributions. Loghman et al. examined timedependent thermoelastic creep analysis of a rotating disk made of Al–SiC composite [10]. They divide the disk thickness into a finite section, and by applying the principle of continuity and boundary conditions solved the differential equations for each section, and subsequently computed deflections and, consequently, thermoelastic stresses and strains. Furthermore, by this analysis, a creep analysis was performed for the disk. Xu and Zhou provided a twodimensional thermoelastic analysis of homogeneous beams with variable thicknesses under thermomechanical forces for simple support [11]. They initially assumed the harmonic series with unknown coefficients for temperature and stress distribution, and then acquired the coefficients applying the upper and lower boundary conditions and edges boundary conditions. They validated their analysis with the finite element method. In another article, Nguyen et al. proposed a method for static analysis and free vibration of the beams made of axially loaded functionally graded materials based on the firstorder shear deformation theory [12]. They explored the influence of the nonhomogeneous index which governs the volume fraction gradation of functionally graded materials on displacements, natural frequencies, buckling loads, and load–frequency curves as well as corresponding mode shapes. Hien Ta Duy et al. inspected the free vibration of FGM beams with variable crosssection [13] . They assumed material properties to vary in thickness following the exponential law. They developed numerical analyses based on finite element method by using a beam element and verified their analytical analysis by this numerical method. Niknam et al. investigated nonlinear bending of a tapered beam consisted of functionally graded material subjected to thermal and mechanical load with arbitrary boundary conditions [14]. As in the previous studies described here, they assumed variable properties in the thickness direction, and provide an analytical closedform solution for which there is no axial force. For a more general state with axial force in the longitudinal direction of the beam, they applied a Galerkin technique and a Generalized Differential Quadrature (GDQ) method for satisfaction of equilibrium. Finally, the results of the two methods were compared.
Applying firstorder shear deformation theory to analyze the famous elements made of functionally graded material is a useful and promising function. Arefi et al. provided an analytical solution for an FG cylindrical shell subjected to mechanical and thermal loading by firstorder shear deformation theory [15]. Heshmati and Daneshmand inquired the effect of different profile variations on vibrational properties of nonuniform beams made of graded porous materials by applying Timoshenko beam theory [16]. Recently, an ameliorated mathematical model is presented in order to examine the free vibration behavior of postbuckled tapered functionally graded material beam, subjected to the thermal load by Amlan Paul and Debabrata Das [17]. They revealed reasonable effects of shear nonlinearly and temperaturedependent thermal conductivity on dynamics of the beam.
In this study, the variable cantilever beam is assumed as a simple blade. It is clear that application of FEM and numerical methods are capable for analysis of exact geometry of the beam and other complicate geometries. For example, Akbari et al. developed a meshless method known as the Meshless Local PetrovGalerkin (MLPG) method for analysis of thermoelastic waves in a twodimensional FGM domain [18]. Rectangular domains have been analyzed under transient thermal loading by this method. Areias and Rabczuk developed a beneficial algorithm for FEMbased computational fracture of plates and shells made of both brittle and ductile materials [19]. This algorithm is simple and effective. Thai et al. proposed a computational method based on the modified strain gradient theory and higherorder shear deformation theory for static bending and free vibration analyses of FG carbon nanotubereinforced composite microplates [20]. A novel algorithm based on combining a NURBSbased inverse analysis with both kinematic and constitutive nonlinearities employed by VuBac et al. [21]. This new method employed widely for analyzing thin shell structures.
In this study, we present a thermoelastic analysis of a functionally graded blade made of ZrO2/Ti6Al4V composite applying firstorder shear deformation theory. Unlike previous studies, not only the material properties are variable longitudinally, but the rotating beam geometry is varying functionally along the beam as well.
Mechanical and thermal properties of the material are linearly variable in the longitudinal direction based on the volume percent of the constituent. Pursuant to this assumption, properties can be described as below.
(1) 
As variable x is the coordinate in conformity with blade length and L is the length of the blade, is the volume percent of reinforcement particles at point x, and are volume percent of reinforcement at the root and tip of the beam, respectively.
Here is the xdependent property, and are the pure matrix and reinforcement property, respectively. This beam is exposed to distributed forces representing the aerodynamical force. The root of the beam is clamped, and the tip of the beam is free. Thermal loading is a temperature field as a result of temperature difference between root and tip of the beam with arbitrary function in terms of the beam length. The sketch of the beam is depicted in Fig. 1. The length of the beam is L. The variation of the beam thickness is expressed by the h(x). In fact, this function exhibits the profile of the beam.
In this section, by the FSDT, an elastic solution is presented. The displacements in the X and Z directions are in accordance with the theory of elasticity as:
(2) 
As u and w are longitudinal and transverse displacement respectively indicates the longitudinal displacement of the neutral axis represents the transverse displacement and is the rotation component. Applying the linear strain–displacement relations, and also considering thermal strain, longitudinal and shear strain components are written asfollow.
(3) 
As is the normal strain in the longitudinal direction and is the shearing strain in xz plane. Generalized Hooke’s law by considering the thermal strains, in this case, is written as below.

(4) 
As and are normal and shear stress in the xz plane. Substituting the strain components from Eq.(3) into the strain energy density equation and integrating over the volume of the structure, the total strain energy of the beam is acquired. The variation of the strain energy in conformity with the principle of virtual work is as below.
(5) 
As , and are the force and moment per unit length defined as follows.
(6) 
The external work due to mechanical loads is represented by:
Fig. 1. sketch of a beam subjected to distributed forces and thermal loading 

(7) 
where W is the total external work, _{ }is the centrifugal inertia body force, is the external torque and is the transverse load. The loading matrix is written as below.
(8) 
In the absence of external torque, F_{2}(x) is set to zero in Eq. (8). According to Eq. (8), the variation of external work is expressed as
(9) 
Pursuant to the principle of virtual work, the strain energy variations are equal to the variation of external work, namely:

(10) 
Deliberating Eq. (5), Eq. (10) and Eq (9) one can obtain:
(11) 
By substituting from relation (6) into Eq. (11) the following constitutive differential equations of the problem are acquired.
(12) 
As A_{i}(x) (i=1, 2, 3, 4, 5, 6) is given in the Appendix A. A_{i}(x) are unknown constants. In order to solve these dissimilar equations with variable coefficients, a semi analytical division method has been applied method [22].
In this method, the solution domain is divided into a large number of subdomains. For each division, properties of the midpoint are assumed to be constant for the same division. Consequently, in the set of constitutive differential equations, variable coefficients are changed into constant coefficients. With respect to the local and global boundary conditions between every two adjacent subdomains, constitutive differential equations yield a set of linear algebraic equations.
Here, a large number of divisions in the longitudinal direction are contemplated (Fig. 2). For each division, the Eqs. (12) are solved and, applying the local boundary conditions introduced in relation (13), and global boundary conditions (14), a set of algebraic equations containing constant coefficients ofeach division are acquired. Solving these equations, constants are obtained, then displacements, stresses, and strains are achieved.
The local boundary conditions that are as a result of the continuity condition are: as below.
(13) 
The tip of the blade is free, and root of the blade is fixed. On that account, the root and tip global boundary conditions of the cantilever FGM blade are as follows. According to the division method, Eq. (12) for Kth division is rewritten as:
Fig. 2. dividing xdirection of cantilever FGM beam 

(14) 
(15) 
As the constants A_{i} (x^{k}) (i = 1, 2, 3, 4, 5, 6) are given in the Appendix B.
An analysis is performed for a cantilever FGM beam whose nonhomogeneity is as a result of variable property in the longitudinal direction, using the method described in the previous section. Upper and lower surface profiles are mathematically represented by the following relations (16) and (17). In these equations x is in millimeters.
(16) 

(17) 
With the above mentioned profiles thickness as a function of x is represented by Eq. (18).
(18) 
In this composite, Poisson’s ratio is constantly supposed, and ν = 0.29. Also, it is assumed that L = 200 mm, r_{o}=300 mm, T_{root}=300 K, T_{tip}=325 K, ω=7000 rpm. The mechanical and thermal material properties are presented in Table 1.
For four different distributions of reinforcement particle A: V_{root}=0 & V_{tip}=40, B: V_{root}=20 & V_{tip}=20, C: V_{root}=0 & V_{tip}=0 and D: V_{root}=40 & V_{tip}=0,the distribution of longitudinal and transverse displacements and shear stresses and effective stresses is computed by the present method.
Figs. 36 depicted respectively longitudinal displacement, transverse displacement, shear, and effective stress of all four distributions nondimensionally. It is clear that the longitudinal tension in the neutral axis is zero.
Table 1. Mechanical and thermal properties
Property (Unit) 
ZrO2 
Ti6Al4V 
E (GPa) 
175 
114.5 
G (GPa) 
69.9 
47.5 
5575 
4470 

Fig. 3. distribution of longitudinal displacement
Fig. 4. distribution of transverse displacement
Fig. 5. distribution of longitudinal stress
Fig. 6. distribution of shear stress
The distributed transverse forces P(x) and the temperature distribution T(x) are mathematically represented in relations (19) and (20) as follows.
(19) 


(20) 
Since the longitudinal tip deflection of the blade is essential here, effect of temperature rise on the maximum longitudinal tip deflection of the blade is inspected. Case A is selected for this study because the maximum effective stress belongs to this case. In Table 2, effect of temperature rise on the maximum longitudinal tip deflection is reported
It is clear that the tip deflection is increased due to temperature rise as expected.
In order to demonstrate the effect of transverse distributed load on the maximum deflections and stresses, analysis for three different transverse loads P(x), 1.2P(x), and 1.5P(x) is carried out, and results are compared with each other. In Table 3, the effect of increasing transverse load on the maximum deflections and stresses is reported. Pursuant to the Table 3, increasing transverse loads increases the maximum transverse deflection, maximum shear stress, and effective stress and decreases the maximum longitudinal deflection
In order to operate the precision and validity of the current study, applying Abaqus software, displacements distribution and effective stress distribution are computed for the case A and compared with the same case obtained from the proposed method. The CPS8R element in the commercial software Abaqus was used for the analysis (ABAQUS Documentation User‘s Manual).
Table 2. Effect of temperature on the maximum longitudinal tip deflection

T(x) 
1.2T(x) 
1.5T(x) 
uE_{Ti6Al4V}/ ρ_{Ti6Al4V}L^{3}ω^{2} 
1.0669 
1.1151 
1.1873 
Table 3. Effect of transverse load on the maximum deflections and stresses

P(x) 
1.2P(x) 
1.5P(x) 
uE_{Ti6Al4V}/ ρ_{Ti6Al4V}L^{3}ω^{2 } 
1.0669 
1.0668 
1.0666 
wE_{Ti6Al4V}/ ρ_{Ti6Al4V}L^{3}ω^{2 } 
0.004342 
0.00521 
0.006513 
σ_{x}E_{Ti6Al4V}/ ρ_{Ti6Al4V}L^{2}ω^{2 } 
2.308 
2.308 
2.308 
σ_{xz}E_{Ti6Al4V}/ ρ_{Ti6Al4V}L^{2}ω^{2 } 
1.107 
1.329 
1.661 
σ_{eff}E_{Ti6Al4V}/ ρ_{Ti6Al4V}L^{2}ω^{2} 
3.00982 
3.2600 
3.689 
Sensitivity analysis indicated the influence of uncertain input parameters on the variance of the outputs [23]. For example, Sensitivity analysis made possible to find the maximum allowable principal stress, and Young’s modulus of the epoxy matrix has the most significant effect on the variance in the fracture energy [24]. Recently, this analysis was applied to extract the key input parameters influencing the energy conversion factor of flexoelectric materials [25].
Here, distribution A is selected for validation. Also, point S=x/L=0.5 is used to sensitivity to mesh size. In order to find the number of elements suitable for accurate analysis, the distribution of displacements and stresses at S is plotted versus the number of elements. The result of this examination is revealed in Fig. 7.
Fig. 7. distribution of effective stress
In this section, in order to validate the method, distribution of longitudinal and transverse displacements and effective stresses in finite element solutions and current research are compared in a beam made of ZrO2/Ti6Al4V composite with the arbitrary distribution. Here, distribution A is assumed for the reinforcement. The comparison is illustrated in Figs. 811.
Fig. 8. Illustrates the fact that if the number of elements is more than 195, the convergence of the results will be achieved for the finite element analysis
Fig. 9. Comparison of the distribution of longitudinal displacement in the finite element analysis and the current
Fig. 10. Comparison of the distribution of transverse displacement in the finite element analysis and the current research
Fig. 11. comparison of effective stress in the finite element analysis and the current research
As it is observed from Figs. 911 that there is a decent agreement between the method presented in this study and the finite element method.
In this study, a semianalytical analysis for thermoelastic behavior of the blade made of ZrO2/Ti6Al4V composite is presented. In fact, this beam is a simplified model of the turbine’s blade. The distributed force represents the aerodynamic force on the blade. Loading is an arbitrary transverse distributed force and a distributed temperature field. The novelty of the current study is that not only the material properties are variable longitudinally but also the beam geometry is varying functionally in conformity with the beam as well.
The governing equations were derived by applying the firstorder shear deformation theory. The solution to the differential equation is carried out through employing a semianalytical method, the continuity conditions at the interfaces, and the global boundary conditions at the root and tip of the beam for each subdomain.
The distribution of reinforcement particles in the matrix of pure Ti6Al4V is assumed to change linearly. However, regarding the proposed method, there is no limitation for the nonhomogeneous distribution, i.e. the model is generic with respect to the distribution. For five different distributions of reinforcement particles A: V_{root}=0 & V_{tip}=40, B: V_{root}=20 & V_{tip}=20, C: V_{root}=0 & V_{tip}=0 and D: V_{root}=40 & V_{tip}=0,the distribution of longitudinal and transverse displacements, longitudinal stresses, shear stresses, and effective stresses are computed applying the proposed method. It was found that in distribution A, the maximum longitudinal displacement and stresses are the minimum rather than all cases. To validate the results, the longitudinal and transverse displacements and effective stress for distribution A were compared to the finite element model. It is found that there is good agreement between the method presented in this study and the finite element method.
Nomenclature
Temperature field (K°) 
T(x) 
Distributed transverse load (MPa) 
P(x) 
Beam length (mm) 
L 
value of the property at x point 
M(x) 
Modulus of Elasticity (GPa) 
E 
Shear Modulus (GPa) 
G 
the thermal expansion coefficient (K^{1}) 

Density (kg/m^{3}) 
ρ 
Acknowledgements
This work was supported by the University of Kashan under Grant number: 468797.
Appendices
Appendix A
Appendix B
[1] Kapania R K, Raciti S. Recent advances in analysis of laminated beams and plates. Part I  Sheareffects and buckling. AIAA Journal 1989; 27(7): 923935.
[2] Shi G, Lam K Y, Tay T E. On efficient finite element modeling of composite beams and plates using higherorder theories and an accurate composite beam element. Composite Structures 1998; 41(2): 159165.
[3] Sankar B V. An elasticity solution for functionally graded beams. Composites Science and Technology 2001; 61(5): 689696.
[4] Sankar B V, Tzeng J T. Thermal Stresses in Functionally Graded Beams. AIAA Journal 2002; 40(6): 12281232.
[5] Chakraborty A, Gopalakrishnan S, Reddy J N. A new beam finite element for the analysis of functionally graded materials. International Journal of Mechanical Sciences 2003; 45(3): 519539.
[6] Nirmala K, Upadhyay P C, Prucz J, Lyons D. Thermoelastic Stresses in Composite Beams with Functionally Graded Layer. Journal of Reinforced Plastics and Composites 2006; 25(12): 12411254.
[7] Kadoli R, Akhtar K, Ganesan N. Static analysis of functionally graded beams using higher order shear deformation theory. Applied Mathematical Modelling 2008; 32(12): 25092525.
[8] Li X F. A unified approach for analyzing static and dynamic behaviors of functionally graded Timoshenko and Euler–Bernoulli beams. Journal of Sound and Vibration 2008; 318(4): 12101229.
[9] Kiani Y, Eslami M R. Thermal buckling analysis of functionally graded material beams. International Journal of Mechanics and Materials in Design 2010; 6(3): 229238.
[10] Loghman A, Ghorbanpour Arani A, Shajari A R, Amir S. Timedependent thermoelastic creep analysis of rotating disk made of Al–SiC composite. Archive of Applied Mechanics 2011; 81(12): 18531864.
[11] Xu Y, Zhou D. Twodimensional thermoelastic analysis of beams with variable thickness subjected to thermomechanical loads. Applied Mathematical Modelling 2012; 36(12): 58185829.
[12] Nguyen TK, Vo T P, Thai HT. Static and free vibration of axially loaded functionally graded beams based on the firstorder shear deformation theory. Composites Part B: Engineering 2013; 55 147157.
[13] Duy H, Van T, Noh HC. Eigen analysis of functionally graded beams with variable crosssection resting on elastic supports and elastic foundation. 2014; 52 10331049.
[14] Niknam H, Fallah A, Aghdam M M. Nonlinear bending of functionally graded tapered beams subjected to thermal and mechanical loading. International Journal of NonLinear Mechanics 2014; 65(Supplement C): 141147.
[15] Arefi M, Faegh R K, Loghman A. The effect of axially variable thermal and mechanical loads on the 2D thermoelastic response of FG cylindrical shell. Journal of Thermal Stresses 2016; 39(12): 15391559.
[16] Heshmati M, Daneshmand F. Vibration analysis of nonuniform porous beams with functionally graded porosity distribution. Proceedings of the Institution of Mechanical Engineers, Part L: Journal of Materials: Design and Applications 2018; 1464420718780902.
[17] Paul A, Das D. Free vibration behavior of tapered functionally graded material beam in thermal environment considering geometric nonlinearity, shear deformability and temperaturedependent thermal conductivity. Proceedings of the Institution of Mechanical Engineers, Part L: Journal of Materials: Design and Applications 2018; 1464420718759376.
[18] Akbari R A, Bagri A, Bordas S, Rabczuk T. Analysis of Thermoelastic Waves in a TwoDimensional Functionally Graded Materials Domain by the Meshless Local PetrovGalerkin (MLPG) Method. 2010; 65 2774.
[19] Areias P, Rabczuk T. Finite strain fracture of plates and shells with configurational forces and edge rotations. International Journal for Numerical Methods in Engineering 2013; 94(12): 10991122.
[20] Thai C H, Ferreira A J M, Rabczuk T, NguyenXuan H. Sizedependent analysis of FGCNTRC microplates based on modified strain gradient elasticity theory. European Journal of Mechanics  A/Solids 2018; 72 521538.
[21] VuBac N, Duong T X, Lahmer T, Zhuang X, Sauer R A, Park H S, Rabczuk T. A NURBSbased inverse analysis for reconstruction of nonlinear deformations of thin shell structures. Computer Methods in Applied Mechanics and Engineering 2018; 331 427455.
[22] Kordkheili S A H, Naghdabadi R. Thermoelastic analysis of a functionally graded rotating disk. Composite Structures 2007; 79(4): 508516.
[23] VuBac N, Lahmer T, Zhuang X, NguyenThoi T, Rabczuk T. A software framework for probabilistic sensitivity analysis for computationally expensive models. Advances in Engineering Software 2016; 100 1931.
[24] Hamdia K M, Silani M, Zhuang X, He P, Rabczuk T. Stochastic analysis of the fracture toughness of polymeric nanoparticle composites using polynomial chaos expansions. International Journal of Fracture 2017; 206(2): 215227.
[25] Hamdia K M, Ghasemi H, Zhuang X, Alajlan N, Rabczuk T. Sensitivity and uncertainty analysis for flexoelectric nanostructures. Computer Methods in Applied Mechanics and Engineering 2018; 337 95109.