Free Vibration Analysis of Reinforced Composite Functionally Graded Plates with Steady State Thermal Conditions

The present paper deals with free vibration of functionally graded fiber reinforced rectangular plates subjected to thermal loads. The rectangular plates are assumed orthotropic. The continuous grading fiber reinforced plates have a smooth variation in matrix volume fraction in the thickness direction. Two different types of volume fraction profiles through the thickness of plate are proposed: classic and symmetric. As the plate is thick, the equations of motion are derived based on three dimensional theory of elasticity. Inevitably, 3D General differential quadrature method is used instead of regular solving methods in order to discretize equations of motion equations as linear set of algebraic equations. The effects of temperature, volume fraction profiles, and boundary conditions are investigated. Some interesting conclusions obtained when the material properties were assumed to be temperature-dependent. It has been observed that temperature and functionality of FG plate have significant effect on the natural frequencies of the plate.


INTRODUCTION
Functionally Graded Materials (FGMs) was introduced in order to avoid problems associated with material mismatch at the layer interfaces in compare with composite materials.The FGM structures are characterized by variable mechanical properties due to the through-the-thickness variation of the volume fraction distribution of the two constituents and the arbitrary thickness profile.By assigning a continuous gradual variation of the mechanical properties along a specified direction, Latin American Journal of Solids and Structures 14 (2017) 886-905 these composites do not show discontinuities in the material.As a consequence, the residual stresses and the stress concentrations that commonly affect a laminated structure can be reduced by mixing two or more constituents according to a specific graded distribution of the volume fraction [Tornabene et al. (2017)].
Functionally graded fiber reinforced composites materials have been widely used in the aerospace marine, and other engineering industries, recently.These new type of FGMs are often proposed to utilize in high temperature as well as corrosive environments.The effect of environment, temperature distribution and other aspects on the material properties of FGM composites and nonhomogeneous materials were studied by many researchers; some focused on reinforcement by fiber or carbon nanotube (Suresh et al. (1998); Nejati et al. (2016); Nejati et al. (2016); Yas et al. (2016), Yas et al. (2010)).Studies depicted that higher temperature and especially corrosive conditional work decrease the elasticity constants and degrade the strength of composites.Using FGMs have been based on this hypothesis that it can produce a material with good resistance on high temperature and corrosive environment.
Some studies have been done focusing on thermal analysis and temperature effects on FGM plates by Yang et al. (2001).Khalili et al. (2012) focused on this effect on free vibration analysis of the FG properties which were assumed temperature dependent.
Within these analyses, some research are based on the classical plate theory i.e. neglecting the effect of transverse shear deformation whereas some others have used elasticity approach to obtain free frequency of laminated and functionally graded composites plates and shells by Matsunaga (2009).
There are limited works about the mechanical stress analysis, free and forced vibration analysis subjected to thermos mechanical loading.A closed-form solution for vibration frequencies of simply supported thick plates was presented by Xiang et al. (1996).Hosseini-Hashemi et al. (2010) presented analytical solution for free vibration of moderately thick rectangular plates which are composed of particle reinforced FGM and supported on elastic foundation.The analysis was performed based on first-order shear deformation plate theory.
In the recent years, analytical and numerical methods have been applied by the research workers in order to study the plate vibration for different shapes.Moreover, the differences between numerical and analytical methods have been investigated.Finite element method is beneficial but expensive to implement as it is time consuming.Isvandzibaei et al. (2014) study on vibration behavior of functionally graded material cylindrical shell.The FGM shell equations have been established based on strain-displacement relationship using Love-Kirchhoff shell theory.The governing equations of motion were solved by energy functional and applying Ritz method.
On the other hand, meshless approaches have found efficient for FG materials since they require nodal connectivity.Betwixt current numerical approaches, differential quadrature method (DQM) has been considerably attended by researchers since it is straightforward to implement as well as quick to analyze.Moreover, it can overcome the most common programming difficulties, such as Latin American Journal of Solids and Structures 14 (2017) 886-905 complex algorithm and excessive use of memory and calculation time.In DQM, the term quadrature applies to the appropriate evaluation of the integral.This technique could be utilized in a simple and systematic fashion to obtain the computational solution of nonlinear differential-integral equations.The DQ method approximates a derivative at a point as a linear weighted sum of all the functional values on the domain.How to evaluate the weighting coefficients is one of the key points of this technique [Tornabene et al. (2015)].However, DQM is not always stable increasing the number of collocation points.For this reason Generalized Differential Quadrature (GDQM) [Tornabene et al. (2014)].
The applied boundary condition type directly influences on the solving method.In GDQ, solving method for conditions like SSSS can be one dimensional while for SCSC it can be 2D.But when 3D elasticity theory is used, the only appropriate solving method for FFFF, CCCC and CFCF conditions is to discretize the governing equation by 3D GDQ.Although for these cases, some are some other theories such as shear deformation theory in which 2D or 1D methods are feasible can be implement, the solving method would be more complicated.
There is not proper papers deal with thermo mechanical vibration analysis of fiber reinforced FGMs.There are several efforts in which mechanical and thermo mechanical dynamic response of FGMs have been studied; whereas in the case of vibration analysis of fiber reinforced FGM plate with different boundary condition and effect of temperature via 3D numerical solution, lack of knowledge is perceived.
In the presented work, it is desired to study on the free vibration of fiber reinforced composite plates with functionally graded volume fraction of fiber resting on different types of boundary condition under thermal-mechanical load.Fiber reinforced functionally graded plates are new type of FGMs in which fiber orientation or/and its density varies functionally through desired direction.The material properties of fiber reinforced FG orthotropic plate are estimated through micromechanical model based on simple rule of mixture and are assumed to be function of temperature.
The governing equations are based on three dimensional elasticity theory (Sadd ( 2009)) and general Von-Karman type strain-displacement equations (Reddy (2004)) which include plate foundation interaction.As various types of boundary conditions were investigated in this paper, 3D general differential quadrature method (GDQ) has been utilized to discretize the governing equation and consequently obtain natural frequencies of FG plate under different types of thermal conditions.Then, the effects of different parameters such as boundary condition, volume fraction as well as temperature are studied.

MICROMECHANICS OF ORTHOTROPIC FGM
The micro mechanic of FG fiber reinforced composites is obtained based on a micromechanical model which is able to calculate effective properties when the volume fraction is changed.Here, two functions are considered to present the variations of fibers' volume fraction through the thickness (Nejati et al. (2015)).The first model is according to Eq. ( 1) in which V depends on P and z simultaneously.For even power (n), fibers volume fraction variations are symmetric: Where V which has values that range from 0 to 1.The exponent n controls the volume fraction profile through the plate's thickness.According to symmetrical model, changes in fiber`s volume fraction occurs from 0.75 in internal surface ( 0 z = ) to 0 in mid-surface, and again reaches to 0.75 in periphery z h = of the plate.The variations incline toward 0 volume fraction of fibers in total thickness of the plate for high even power indices cause the plate to be similar to a homogenous orthotropic plate with 0% volume fraction of fiber and 100% matrix in total thickness.In the second model distribution of fibers is according to conventional (classic) form according to Eq. ( 2): Based on Classical model, fibers` volume fraction continuously reaches from zero in one surface of the plate to 0.75 in the other one.According to this relation, i V and o V which have values that range from 0 to 1, denote the volume fractions (matrix or fiber) on the inner and outer surfaces, respectively.Fibers volume fraction is zero in 0 z = and 0.75 in z h = .The alteration inclines toward 0 volume fraction of fiber for high power indices, and the plate looks like a homogenous orthotropic plate with 0% of fibers and 100% of the matrix volume fraction for the total plate thickness in this model.
In the present analysis, it is assumed that fiber angle of the fiber reinforced plate are constant with respect to the x-axis in the x y surface.The effective mechanical properties of the fiber rein- forced plate are obtained based on a micromechanical model as follows (Shen (2009), Vasiliev et al. (2001)): Where , ,

CONSTITUTIVE EQUATIONS
The equations of motion for free vibration analysis of fiber orientation FG rectangular plate can be obtained by using 3D elasticity theory in conjunction with Hamilton principle which is (Reddy (2004)).
Where KE is the kinetic energy of the rectangular plate and PE is the elastic potential energy of it.
The kinetic energy of the plate is obtained as following: The elastic potential energy of plate is comprised of strain energy due to vibration and primary stresses resulting from temperature increase respectively which both are obtained as follows: Latin American Journal of Solids and Structures 14 (2017) 886-905 Where PEP is the strain energy due to vibratory stresses and PET is the strain energy caused by the initial stresses due to temperature rise 0 ( , , , ) ij i j x y z s = is pre-stress components due to applied temperature field and ij s is stress tensor.Strain-displacement relations are considered according to the 3D elasticity theory including nonlinear terms that originally was proposed by Karman and Von (Reddy (2004)).The constitutive relation for orthotropic rectangular plate is according to general form of Hook's law.
The plate is initially stressed free at temperature to hence pre-stress components should be included by increasing the temperature.The pre-stress components due to temperature rise ΔT are defined by: It is assumed that temperature varies uniformly and non-uniformly through the thickness.The elements of Off-axis stiffness matrix [C] are given in (Reddy (2004)).

The Equations of Motion
By substituting strain-displacement and constitutive relations and Eq. ( 9) into elastic potential energy, Eqs. ( 5) and ( 6) respectively then combining with Hamilton principle the equations of motion with boundary conditions are extracted for fiber reinforced FG plate as follow: Also, the boundary conditions of the plate are described at the ends z=0 and h as follow: at x=0, a: at y=0, b: Where, x F , y F and z F are generalized forces along x, y and z directions.For normal mode analysis, the following solution is assumed for the displacement components x, y and z.Where coefficients U, V and W, are unknown displacement functions of x, y and z orientations, respectively.Combination of Eqs.(10-15) with Eq. ( 16) will result the equations of motion as a set of coupled partial differential equations in terms of displacement components.These set of equations have been solved by the general differential quadrature (GDQ) discretization procedure for spatial derivations.
Using the 3D-GDQ method for the spatial derivatives, the discretized form of the equations of motion (Eqs.10-12) at each domain grid point can be obtained as Appendix A.
For fully clamped boundary condition on four edges (CCCC_built-in edge), the discretized boundary conditions are obtained as: For SCSC type boundary condition, two opposite sides are fully clamped and two others are simply supported, we have: In the case of simply supported plate or when all edges are rested on simply supported foundation (SSSS), boundary condition is expressed on the x-constant and y-constant edges as follows: Nodes (or grids) distribute along three axes non-uniformly which is defined as Chebyshev-Gauss-Lobatto distribution as follow: Fig. 2 shows DQ meshing with uniform distribution of nodes ( Using this procedure for all of the governing equations and converting them into set of equation in terms of series form, leads to discretized equations.Then governing equations in the form of differential quadrature relations is produced in terms of series and unknown variables U, V and W. The first and second derivatives of variable U with respect to x, y and z are computed as follows: In the above relations, , , The weighting coefficients for the first-order derivatives in x-direction are thus obtained as (Bert et al. (1996)): Where x L is the length of domain along the x-direction.
1, ( ) ( ) According to DQ method, for higher-order derivatives, higher-order weighting coefficient should be employed, i.e. for second order derivatives, weighting coefficient can be expressed as below: In a similar way, the weighting coefficients for the y-direction and z-direction can be obtained.To perform the eigenvalue system of equations, the degree of freedom is separated into the domain and the boundary degree freedom.Structures 14 (2017) 886-905 In order to carry out the eigenvalue analysis, the domain and boundary degrees of freedom are separated into vector forms.They are denoted as { } d and { } b , respectively.Based on this definition, the discretized form of the motion equations and the related boundary conditions take the following forms (Bert et al. (1996)): Equations of motion (Appendix A. (A1-A3)):

Latin American Journal of Solids and
Boundary conditions (Eqs.(17-19)): By eliminating the boundary degrees of freedom in Eq. ( 25), using Eq. ( 26), this equation turn into: Where . The above eigenvalue system of equations can be solved to find the natural frequencies and mode shapes of the plates.

THERMAL ANALYSIS
Here, the influence of a temperature filed on the behavior of the FGM is investigated.The temperature varies in the thickness direction.Two types of thermal boundary conditions are applied: uniform and non-uniform.
In the case of uniform temperature rise, the temperature difference from the reference temperature o T , for both upper surface In non-uniform distribution, two different temperatures is imposed on the two surfaces while other sides are isolated.In this case, the temperature distribution along the thickness can be obtained by solving a steady-state heat transfer equation through the thickness of the plate.The equation for the temperature through the thickness is given by: Where ( ) z K z is the thermal conductivity.This equation is solved by imposing boundary condition as following: Latin American Journal of Solids and Structures 14 (2017) 886-905

RESULTS AND DISCUSSION
To validate the current three dimensional GDQ formulations, the results for an isotropic and homogenous rectangular FGM plate are compared with similar ones obtained by FEM and in the literatures by 2D GDQ.

Validation
First of all, convergence and validation study of the normalized first six natural frequencies is considered for an isotropic FGM plate in Table 1 to compare with Malekzadeh (2009) and Matsunaga (2008).
Fast rate of convergence of the method is evident for exponent index ''p=1" and it is found that only 11 DQ grid for simply supported FGM plate in the three dimensions can yield accurate results.Good agreement is observed between the results.
In order to evaluate accuracy of 3D GDQ method, non-dimensional natural frequency resulted from this method are compared with 2D GDQ formulations by Malekzadeh (2009) for as well as FEM method in Table 2 for different exponent index (p) for isotropic FGM simply supported plate and 1D GDQ formulation by Yas et al. (2010) in Table 3 in different thickness to width ratio (h/b) for orthotropic FGM simply supported plate while exponent index of 1 (p=1) and length to width ratio of 1 (a/b=1).Good agreement between all the data in Tables 2 and 3 proves the reliability of presented method.It should be noted that 13500 elements of second order 8 nodes with 250 sublayers were utilized in order to obtain slightly results.Moreover, in this example 13 GDQ grids were used to achieve adequate convergence and concordant results with Malekzadeh (2009).
The comparison shows that the present results agreed well with those in the literature.
Latin American Journal of Solids and Structures 14 (2017) 886-905 Temperature dependent properties of copper and tungsten are described in Eq. 30 by Gasik (2000).) Solids and Structures 14 (2017) 886-905 Convergence of current three dimensional GDQ approach is evaluated by obtaining first natural frequency of fiber reinforced orthotropic FGM plate with classic distribution, non-uniform temperature and temperature dependent properties.The plate is on SCSC foundation whose two opposite edges have same boundary conditions.Fig. 3 shows convergence rate of results achieved by different values of Nz when Nx and Ny has specified constant values.Reasonable convergence results and accurate natural frequency is obtained in at least a 13×13×15 or 15×15×21 or 17×17×19 nodes grid.

Effect of Temperature on the Natural Frequency
Here natural frequency results, Ω, versus ΔT are evaluated; which ΔT is divergence between maximum temperature of the plate and ambient temperature.According to part 4, maximum temperature of an orthotropic FGM plate in non-uniform temperature distribution occurs in the upper side at z=h while in uniform temperature distribution, whole the plate has a unique temperature which is considered as maximum plate temperature.Fig. 4 demonstrates effect of temperature variation on the first natural frequency in the plate for uniform and non-uniform (see Eq. 28) temperature distributions in presence of three different types of boundary conditions (clamped-clamped, simply supported-simply supported and simply supported-clamped).Besides, effect of considering material property as temperature dependent or independent is evaluated.In this figure, natural frequency would reduce when the temperature rises.Also, including dependency of material property to temperature generally causes the natural frequency to reduce more.Indeed, the results of natural frequency analysis for such these plates which would be considered without effect of temperature dependency in material property should be dealt with more carefully.

Latin
In another point of view, regardless boundary condition type and temperature dependency, assuming uniform temperature distribution results in lower natural frequency.Comparing the graphs Latin American Journal of Solids and Structures 14 (2017) 886-905 in Fig. 4 shows CC-CC boundary condition causes that stiffness of the plate increases and natural frequency of SS-SS type boundary condition is the least.

Effect of Profile of Fiber Volume Fraction on the Natural Frequency
Fig. 5 compares the natural frequencies between classic and symmetric distribution of volume fraction of fiber through the thickness.As the results show, natural frequency of the plate with a symmetric distribution of fiber is considerably higher than one with classic distribution.The other parameters such as distribution of temperature and type of boundary condition affected the plate with symmetric variation of volume fraction just as the classic one.
Latin American Journal of Solids and Structures 14 (2017) 886-905 Tables 4 and 5 investigate over fiber volume fraction distribution influence on the natural frequencies in another aspect.In these tables, the volume fraction power index of fiber varies for orthotropic FG reinforced plate subjected to different B.C. and in different situation of uniform and non-uniform temperature distributions.The natural frequencies were found for different temperature conditions.The results show that according to equation 2, when its power index increases, the natural frequencies would decrease; a result which is held for the first natural frequencies in all the considered conditions.

CONCLUSION
In this study, natural frequency of a functionally graded fiber reinforced thick plate was investigated.Response of this plate to variation of some parameters such as volume fraction profile, boundary conditions and mainly temperature distribution was examined.Because the plate is thick, a 3D GDQ method was used to solve the equations.The following conclusions are drawn:  The results proved credibility and fast convergence of 3D GDQ method for heat transfer equation and equations of motions for such these plates.Besides, replacing this method by regular solving method such as FEM and combination of high order shear deformation (instead of 3D elasticity) theory and 2D GDQ for thick FG fiber reinforced plates would not only save the time but also relived the solution from unnecessary complexity. Alteration of temperature in FG fiber reinforced plate significantly affects its natural frequency.An increase of temperature leads to a reduction in natural frequency for different types of boundary conditions.As property of the plate transforms in a wide range of temperature, evaluation of the FGFR plate would be complicated when the temperature changes.The results shows including dependency of material property to temperature generally causes the natural frequency to reduce more; in other words, FGFR plates which are more independent from temperature have higher natural frequency in high temperatures. In comparison between uniform and non-uniform temperature distribution in the plate, the uniform one results in lower natural frequency. CC-CC boundary condition causes stiffness of FGFR plate under temperature variation increases that means natural frequencies are the highest among all considered boundary conditions while SS-SS condition graphs demonstrates the lowest one. Volume fraction profile affect the natural frequency of FGFR plates subjected to thermal loads.Natural frequency of the plate with a symmetric distribution of fiber is considerably higher than one with classic distribution and for higher power indexes, the natural frequencies would be lower.

Figure 1 :
Figure 1: Schematic of fiber reinforced FG rectangular plate with fibers` volume fraction variation along thickness.
and f r are elasticity modulus, shear modulus, Poisson's ratio and density of the fiber respectively, and , , and m r are corresponding properties for the matrix.f V and m V are the fiber and matrix volume fractions in order and are related by

Figure 2 :
Figure 2: Non-uniform distribution of GDQ meshes in three directions.(a) X-Y plan, (b) X-Z plan, (c) Y-Z plan.
are weighted coefficients of first and second order derivatives along the X, Y and Z directions respectively.

Figure 3 :
Figure 3: Convergence of first frequency parameter of fiber orientation with constant value of Nx and Ny

Figure 4 :
Figure 4: Variation of first natural frequency of an orthotropic FG plate in different conditions of material property and temperature distribution by ΔT in presence of (a) CC-CC, (b) SC-SC, (c) SS-SS

Figure 5 :
Figure 5: Variation of first natural frequency of an orthotropic FG plate in different conditions of volume fraction distribution, material property and temperature distribution by ΔT in presence of a) CC-CC, b) SC-SC, c) SS-SS

Table 1 :
Convergence and accuracy of first three frequency parameters of

Table 2 :
Validation of natural frequencies for a FGM plate with SS-SS boundary

Table 3 :
Comparison of the first three non-dimensional natural frequencies of a simply supported FGM orthotropic plate with classic model (p=1, a/b = 1).

Table 4 :
First natural frequencies of an orthotropic FG plate in uniform temperature distribution and different B.C. for different volume fraction profiles in different ΔT

Table 5 :
First natural frequencies of an orthotropic FG plate in non-uniform temperature distribution and different B.C. for different volume fraction profiles in different temperatures