Received 17 July 2016; accepted 9 August 2016; published 12 August 2016
Functionally Graded Materials (FGMs) are new advanced heat resistant materials that are used in modern technologies as smart structures. Besides excellent thermal properties, they are erosion and corrosion resistant and have high fracture resistance. The fundamental theory is to mix metal and ceramic such that the material properties continuously vary from one basic material to the other. Actually, the governing equations for the temperature and stress distribution are coordinate dependent. As well, other material properties are functions of position.
There is massive amount of research in both of functionally graded beams and the thermal behaviors of beams. An elastic solution is given for a functionally graded beam by Snaker  . A beam theory similar to the Euler- Bernoulli beam theory is developed, and the beam theory results are compared with elastic solutions. Shahsiah et al.  studied the thermal buckling of functionally graded beams. The normalized functions proportional to the thermal buckling loads for thin beams made of functionally graded material are derived when beam is under uniform temperature rise and axial temperature difference. On the other hand Davoodinik et al.  have derived the fully non-linear equations governing the non-linear deformation of the flexible cantilever Tapered Functionally Graded Beam (TFGB) under combined inclined end loading and intermediate loading and could be solved by the semi-analytical method. In addition Davoodinik et al.  had derived an analysis of thermal behavior of Functionally Graded Beam (FGB). They assumed that the distribution of material properties followed exponential function, while for thermal loading the steady state of heat conduction with exponentially and hyperbolic variations through the thickness of FGB.
Esfahani et al.  had studied the thermal buckling and post-buckling analysis of Functionally Graded Material (FGM) Timoshenko beams resting on a non-linear elastic foundation. Different types of boundary conditions such as clamped, simply supported, and rolled edges are assumed for edge supports. Giunta et al.  had investigated the mechanical behavior of three-dimensional beams subjected to thermal stresses. Conjunction, the temperature field was obtained by exactly solving Fourier’s heat conduction equation and, as classically done by a staggered solution approach, it was considered as an external load within the mechanical analysis. Arbind et al.  developed the governing equations for a general third order theory for functionally graded beam considering the modified couple stress theory and von Kármán nonlinearity using Hamilton’s principle.
Cheng et al.  have obtained a closed form solution for the thermo-mechanical deformations of an isotropic linear thermo-elastic functionally graded elliptic plate rigidly clamped at the edges. Javaheri et al.  had derived the equilibrium and stability equations for the rectangular functionally graded plates on the basis of classical plate theory. Resulting equations were employed to obtain the closed-form solution for the critical buckling temperature. Kadoli et al.  had studied the static behavior of functionally graded metal-ceramic (FGM) beams under ambient temperature. The FGB with variation of volume fraction of metal or ceramic based on power law exponent was considered. Using the principle of stationary potential energy, the finite element form of static equilibrium equation for FGB was presented. Sankar et al.  have introduced a solution for thermo- elastic equilibrium equations for a functionally graded beam in a closed-form to obtain the axial stress distribution. The thermo-elastic constants of the beam and the temperature were assumed to vary exponentially through the thickness. The Poisson ratio was held constant. The exponential variation of the elastic constants and the temperature allowed exact solution for the plane thermo-elastic equations. Based on the nonlinear first-order shear deformation beam theory (FBT) and the physical neutral surface concept, Ma et al.  have derived the governing equations for both the static behavior and the dynamic response of FGB subjected to uniform in-plane thermal loading. Moreover Zhong et al.  had studied the cantilever functionally graded beam subjected to different loads. In terms of Airy stress function a general two-dimensional solution was presented for a cantilever functionally graded beam, assuming that all the elastic moduli of the material had the same variations along the beam-thickness direction. Chi et al.  had studied the Mechanical behavior of functionally graded material plates under transverse load. Based on the classical plate theory and Fourier series expansion, the series solutions of power-law FGM, sigmoid FGM, and exponential FGM plates were obtained. The analytical solutions of the previous distributions FGM plates were proved by the numerical results of finite element method  . While in  they evaluated the numerical solutions directly from theoretical formulations and calculated by finite element method using MARC program. On the other hand, Sun et al.  had numerically investigated the buckling and post-buckling thermo-mechanical deformations of a Functionally Graded Material (FGM) Timoshenko beam resting on a two-parameter nonlinear elastic foundation and subjected to only a temperature rise with the shooting method. Meanwhile, Akbas  had investigated the free vibration and static bending analysis of Functionally Graded (FG) beams resting on Winkler foundation within Euler Bernoulli beam theory and Timoshenko beam theory. The material properties of the beam changed in the thickness direction according to power-law distributions
This study is intended to introduce a mathematical analysis to the functionally graded beam under thermal loading. In this study the material properties are assumed to vary according to a power law while the Poisson’s ratio is considered to be constant along the beam height. The heat distribution along the beam height is assumed to have two behaviors: power distribution (P-state) and exponential distribution (E-state).
2. Theory and Formulation
2.1. Material Gradient of FGB
The concept of functionally gradient materials is to make a composite material by varying the microstructure from one material to another one with a specific gradient. This enables the material to have the best interaction of both materials. The transition between the two materials can usually be approximated by means of a power- law function (P-FGB), exponential function (E-FGB), or sigmoid function (S-FGB). In this study, the material properties distribution is defined by power function through the thickness direction of the beam (i.e. z-direction) along the height (h) from z = 0 until z = h (Figure 1). While, the Poisson’s ratio is considered constant along the beam height. Consider a P-FGB with different variations of modulus of elasticity, thermal expansion and thermal conductivity through the thickness direction of beam, which gives:
where E, a and k are the modulus of elasticity, the thermal expansion and the thermal conductivity of the beam, respectively. The constants Eo, ao, ko, ω, b and l can be obtained by boundary conditions.
2.2. The Governing Equations
In accordance with the first order shear deformation theory, a point with a distance z from lower surface (z = 0) will be moved, after deformation with certain displacement  . The axial displacement at this point with distance z can be presented by
where uo is the displacement at the lower surface (i.e. z = 0) and w is the transverse deformation, ∂w/∂x and fy are the rotations and all of them are independent of z-direction  . Thus the strain in the x-direction could be described as follows,
The stress-strain relations for the FGB which subjected to thermal load, and plane strain condition are,
where, ν and T(z) are the axial stress on the surface with distance z from zero surface, the Poisson’s ratio and the temperature distribution along the z-direction of the beam, respectively. It is required to evaluate the axial stress distribution along the height of the beam taking into consideration the behavior of material which is described in Equation (1). Thus Equation (4-b) is the axial stress distribution along the z-direction, furthermore it involves three unknowns which are independent of z-direction. These unknown terms could be derived from
Figure 1. Rectangular cross section FGB.
the boundary conditions and the equilibrium equations.
2.3. The Axial Stress
When the beam is in equilibrium, the axial resultant forces in the x-direction must be zero value, which yields; (assuming the thermal loading only along distribution in z-direction);
where, and b are the length and width of the beam, respectively. While the equilibrium equation for the resultant moment is, due to the resultant moments about Y-axis appears by thermal effect only (MTy);
where Iy is the moment of inertia about Y-axis, MMy is the mechanical moment about the Y-axis and MTy is the moment due to the thermal affect which could expressed as:
At this instant, Equation (5) and Equation (6) give two of the three required equations for determination the three unknowns in Equation (4-b). The third one could be provided through the boundary conditions which are, for the simply supported beam:
where is the total moment about Y-axis. Thus using Equations (5), (6), (10), the axial stress in Equation (4-b) could be obtained. Regard as conventional definitions for simplicity as follows:
Substituting from Equation (11) in Equation (5), Equation (6), Equation (10); yields, the following form;
It is clear that A2 is directly obtained, while the coefficient A1 and A3 could be given as follows:
Therefore, using Equation (13), the axial stress could be obtained as:
Moreover, Equation (14) could introduce the axial stress distribution taking into consideration the temperature distribution along the beam height, which will be discussed next.
2.4. The Transverse Shear Stress
The thermal bending occurring in the x-z plane could be presented as follows; 
where txz is the transverse shear stress. Rearranging Equation (15) gives:
Using Equation (14) into Equation (16), yields;
For the simply supported beam it is obvious that the shear stress is going to zero value because of absent of second derivatives of uo, fy and the third derivatives of w.
2.5. Steady State Temperature Distribution
Firstly starting with the power distribution of temperature (P-State). The governing equation of heat conduction; assuming one-dimensional transient heat conduction is given as; 
Thus the heat distribution, T(z) through the thickness h, in the direction z is, according to Equation (1-c);
where constants B1 and B2 could be obtained from the thermal boundary conditions, which are:
Therefore, the thermal stress distribution σxx(z) along the height of the beam (h) could be obtained.
Whereas for the exponential distribution of temperature (E-State); the governing equation of heat conduction could be obtained; referring to Equation (18) but with the exponential distribution of temperature. The heat distribution, T(z) through the thickness h, in the direction z is,
where D1 and D2 could be obtained from the thermal boundary conditions; assume again
Moreover, the thermal stress distribution σxx(z) along the height of the beam (h) could be obtained with the exponential distribution of temperature.
3. Case Study
Consider an elastic square cross section beam (Figure 1) with simply supported ends and with height and width equal to 0.05 (m). First assuming the P-state for temperature distribution, which yields the temperature distribution as in Equation (19) with temperature boundary conditions:
where T1 and T2 present the temperature of stress free state i.e. temperatures at surface z = 0 and surface z = h, respectively.
The effect of the parameters (ω, b and l) on the normalized thermal stress (which is σxx(z) divided by σn) distributions are illustrated Figures 2-6, where
The final expression of the axial stress, Equation (14), as a function of A1, A2 and A3 can be obtained by using a MATLAB program to perform the integrations of the expressions of Equation (11-b).
Figure 2. The normalized axial stress distribution () along the thickness of the beam (z/h) for different values of constants;, , , , , b = h = 0.05 m.
Figure 3. The axial stress distribution () along the thickness of the beam (z/h) for different values of constants;, , , , (−ve),.
Figure 4. The axial stress distribution () along the thickness of the beam (z/h) for different values of constants, , , , , ,.
Figure 5. The axial stress distribution () along the thickness of the beam (z/h) for different values of constants, , , , , ,.
Figure 6. The normalized axial stress distribution () along the thickness of the beam (z/h) with different values of constants, , , , , ,.
The effect of the positive values of parameters (ω, b and l) on the normalized axial stress () is shown in Figure 2. The values of the normalized axial stresses increase with increasing of the value the beam coefficients until they reach higher values with ω, b and l more than 3. While The values of the normalized stress are almost constant along the beam height with small values of ω, b and l. On the other hand, the effects of the negative values of parameters (ω, b and l) on the normalized stress are shown in Figure 3. As the negative value increases of (ω, b and l), the normalized stress, also increases.
Figure 4 illustrates the effect of the parameter (ω) on the normalized axial stress with constant values of b and l (b = l = 3). The normalized stress increases as the parameter (ω) increases until it reaches almost 3 ´ 104 MPa with (ω) equals to 4 at h = z. Moreover, the effect of the parameter (b) on the normalized stress with constant values of ω and l (ω = l = 3) is shown in Figure 5. The thermal expansion parameter (b) has a positive increasing effect on the normalized stress along the beam height (the normalized stress reaches 2.5 ´ 104 at b = 4).
Whereas the effect of the thermal conductivity parameter (l) on the normalized stress with constant values of b and ω (b = ω = 3) is shown in Figure 6. It is noticed that the increases of the parameter (l) has a remarkable effect on increasing the values of normalized stress especially at the outer fiber of the beam (z/h = 1.0).
4. Numerical Study
After Consider the elastic square cross section beam (Figure 1) with simply supported ends and with height and width equal to 0.05 (m), with and. Where T1 and T2 are indicated the temperature of stress free state i.e. temperatures at surface z = 0 and surface at z = h, respectively. Taking a FGB with material properties given in Table 1. Assume at surface z = 0 (unmixed Ni) and at surface z = h (unmixed Al2O3).
With dividing the thermal responses of beams to two parts, one part for explaining the thermal behavior of beam without any end supports (referred as part A) and another part for representing the thermal behavior of beam with considering the end support influence (referred as part B) when the beam is exposed to P-state or E-state, for the temperature distribution. Thus for these two parts and considering Equation (14), the following can be obtained;
where the σA, σB and eA, eB are the stress and strain resultants of part A and part B, respectively. Then the total stresses and strains are as follows;
Furthermore, comparing the beam in Table 1 with a traditional beam, made of one material, (such as Ni only) will be done (could be named also as Timoshenko’s beam). Referring to the beam in Table 1 as beam (F) (related to FGB); while the other beam will be referred as beam (N) (related to Ni). Figure 7 and Figure 8 show the thermal resulting distribution for stress and strain, respectively for the beam in Table 1 with the P-state temperature distribution. It could be noticed from Figure 7 and Figure 8 that both of stress (A) and stress (B) are almost symmetrically distributed along the beam height. On the other hand, Figure 9 shows the normalized axial stress distribution along the beam height for both beam (F) and beam (N) for P-state. Indeed both of the curves in Figure 9 starts at the same point which is obviously true, while the stress for beam (F) increasing gradually whereas the stress for beam (N) is almost a sinusoids distribution.
Alternatively, Figure 10 and Figure 11 show the thermal resulting distribution for stress and strain, respectively for the beam in Table 1 with the E-state temperature distribution. It could be recognized from Figure 10 and Figure 11 that both of stress (A) and stress (B) are less symmetrically distributed along the beam height, that was because of the great effect of the exponential temperature distribution.
Moreover, Figure 12 shows the stress distribution by Rahimi and Davoodinik  which assumed the material properties distribution follows exponential function whereas in this work the material properties distribution follows power function. It is noticed that the results obtained in the present work (Figure 10) are in good agreement with which obtained by Ref.  (Figure 12).
Table 1. Materils properties of Al2O3-Ni beam  .
Figure 7. The total stress distribution along the beam height with P-state.
Figure 8. The total strain distribution along the beam height with P-state.
Figure 9. The normalized axial stress along the beam height for beam (F) and beam (N) (P-state).
Figure 10. The total stress distribution along the beam height for E-state.
In addition, Figure 13 shows the normalized axial stress distribution along the beam height for both beam (F) and beam (N) for E-state. Again both of the curves in Figure 13 starts at the same point while, the stress for beam (N) has values higher than the values of stress for beam (F). This behavior is due to the noticeable effect of the E-state distribution of the temperature on the beam (N) which is almost totally elastic.
Figure 11. The total strain distribution along the beam height with E-state.
Figure 12. The total stress distribution along the beam height with E-state dotted right stress (A), dotted left stress (B), solid in middle total stress  .
Figure 13. The normalized axial stress along the beam height (z/h) for beam (F) and beam (N) (E-state).
A mathematical analysis was introduced to study the functionally graded beam under thermal loading. In this study the material properties were assumed to change according to a power law while the Poisson’s ratio is considered constant along the beam height. The heat distribution along the beam height is assumed for two cases: P-state and E-state. The form of temperature distribution has an enormous effect on the resultant thermal stress distributions. It was found that the exponential distribution of thermal loading has more effect than the power distribution on the thermal stress distribution. The functionally graded beam was compared with a traditional beam which also could be called as a Timoshenko’s beam. The analysis introduced in the present work could be extended to any functionally graded beam with power function material property and different beams dimensions.