Performance of functionally graded plates under localised transverse loading

This paper presents a study of the bending of an isotropic functionally graded plate under localised trans- verse load through a combination of analytical and computational means. The analytical modelling is based on the recently developed three-dimensional elasticity solution, expanded to cover different load- ing types, whilst the Finite Element model uses graded isoparametric elements. The plate under consid-eration is assumed to be simply supported, with Young’s and shear moduli varying exponentially through the thickness and the Poisson’s ratio constant. Comparative analysis of stress and displacement ﬁelds in functionally graded and homogeneous plates subjected to uniformly distributed and patch loadings is carried out.


Introduction
Functionally Graded Materials (FGMs), described in detail by Suresh and Mortensen [1], are a type of heterogeneous composite materials exhibiting gradual variation in volume fraction of their constituents from one surface of the material to the other, resulting in properties which vary continuously across the material. The idea of a Functionally Graded Material is not a new one, there are in fact many natural materials which exhibit this property. Study of bone, shell, balsawood and bamboo shows that they are all graded with their greatest strength on the outside, in areas where the greatest protection is required. However it was not until the 1980s in Japan [2] that the idea of a Functionally Graded Material was actively researched in order to gain advances in heat resistant materials for use in aerospace and nuclear fission reactors.
There are currently two main methods of processing Functionally Graded Materials. The first, known as constructive processes, involve layering up the material in layers of varying volume fraction (i.e. the ratio of volume present of each of these two constituents) to create the variation in properties. Examples of constructive processes include solid state powder consolidation and liquidphase sintering. Transportation processes are the second method used to process Functionally Graded Materials. These processes according to Suresh and Mortensen [1] ''rely on natural transport phenomena; the flow of fluid, the diffusion of atomic species or the conduction of heat to create gradients in local microstructures''. It is this gradient in microstructure which provides the required graduation of properties. Examples of transport processes are carburisation, nitriding and ion implantation. All of the examples outlined above are covered in detail in Suresh and Mortensen [1].
Functionally Graded Materials are usually designed specifically for their application, current applications include dental implants, heat exchanger tubes, engine components and to eliminate mismatch of thermal properties in metal and ceramic bonding. Unlike laminated plates, functionally graded plates contain no interfaces, eliminating stress concentrations at boundaries and problems such as delamination.
In order to understand the performance of Functionally Graded Materials and structures a number of different approaches, both analytical and numerical, have been undertaken by researchers and a review of the main developments up to 2006 was published by Birman and Byrd [3]. Stress and deformation of statically loaded functionally graded rectangular plates were analysed by a number of researchers, most recently in [4][5][6][7][8][9][10][11][12][13][14][15]. Based on classical plate theory, Chi and Chung [4] developed a solution for functionally graded plates with power, sinusoidal or exponential variation in stiffness properties. An exponentially graded plate was also considered by Zenkour [5] who developed both a higher order trigonometric and 3-D elasticity solution and presented results for stresses and displacements for a variety of plate geometries. Nguyen et al. [6] presented a first-order shear deformation model of a functionally graded plate, comparing results with those available in the literature. Gilhooley et al. [7] and Matsunaga [8] developed two-dimensional higher-order deformation theories and show that results compare favourably with Finite Element and 3-D solutions for moderately thick plates. However, Matsunaga [8] concluded that higher orders of expanded 2D theories may be necessary to obtain reasonably accurate solutions for very thick plate. A two dimensional theory, utilising a mixed variational approach was presented by Fares et al. [9] and results were compared with those obtained through use of three-dimensional approaches. Zhong and Shang [10], Liu and Zhong [11] and Chun and Zheng [12] presented three dimensional bending analyses of simply supported functionally graded plates and studied the effect of varying the distribution of Young's modulus through the plate thickness. Through a number of examples it is shown that changes in stiffness gradient have a significant effect on the total plate stiffness and its deflection. A three dimensional elasticity solution for functionally graded transversely isotropic rectangular plates under uniform loading was presented by Yang et al. [13] and good agreement is found with results from the literature. Through use of displacement functions for inhomogeneous transversely isotropic media Woodward and Kashtalyan [14] presented a three-dimensional elasticity solution for transversely isotropic graded plate. Xu and Zhou [15] derived expressions for stresses and displacements of an exponentially graded plate with variable thickness.
In addition to analytical methods of modelling Functionally Graded Materials, a number of numerical methods, such as the 3-D Finite Element method have also been employed. The main issue encountered in application of the Finite Element method is how to model a material with continuously varying properties. The simplest and crudest method involves the use of homogeneous elements each with different properties, giving a stepwise change in properties in the direction of the material gradient. This method has already been used by several researchers and can give reasonably accurate results. Etemadi et al. [16] modelled a sandwich panel with functionally graded core subject to low velocity impact, Zhang et al. [17] modelled the contact response of a functionally graded coating and Tilbrook et al. [18] and Wang and Nakamura [19] modelled the propagation of cracks in Functionally Graded Materials. There are however several problems in using this approach, which are discussed by Buttlar et al. [20]. Due to the model approximating the continually varying properties of the material with stepwise changes, a computational error is always introduced, particularly under cases of high stiffness gradient. In order to minimise this error, a very fine mesh in the direction of the property gradient is often used. This can however lead to extremely long computation times. A more advanced method of including property variation into a Finite Element model is to use elements which themselves contain gradient in properties. The difference between graded and stepwise elements is shown in Fig. 1. The formulation and implementation of graded elements is described in detail by Kim and Paulino [21] who also compare the behaviour of graded and homogenous elements under various loading conditions and compare with analytical solutions in literature. They conclude that graded elements give far greater accuracy when modelling FGMs.
In the current paper, the elastic deformation of a functionally graded plate subject to localised transverse loading is analysed through a combination of analytical and computational means. The analytical model employs a 3-D elasticity solution for stress and displacement fields in functionally graded plates subjected to a one-term sinusoidal loading recently developed by Kashtalyan [22]. As any loading can be expressed as a Fourier series involving terms of this kind, the solution [22] provides a basis for further study of the behaviour of graded plates under various types of loading. The authors expand solution [22] to examine three different types of loading: uniformly distributed loading (UDL), large centralised patch and small centralised patch. The computational model involves a Finite Element analysis of the same problem, using ABAQUS software and user implemented graded elements. Stresses and displacements are calculated through the thickness of the graded plate and compared with those in a homogenous plate.

Problem formulation
A functionally graded plate (Fig. 2) of length a, width b and thickness h is referred to a Cartesian co-ordinate system The plate is assumed to be simply supported at the edges such that where r ij are components of the stress tensor and u i are components of the displacement vector. The boundary conditions, Eqs.
(1), are representative of roller supports and analogous to simply supported edges in the plate theories (see Kashtalyan [22]).
The loading is applied transversely to the upper surface of the plate and provides the final six boundary conditions, such that: The loading Q(x 1 , x 2 ) is assumed to allow expansion into a Fourier series where m and n are wave numbers and q mn is the loading coefficient in the double trigonometric series expansion for a given load type and is defined in Eq. (2c).
The bottom surface of the plate is assumed to be load-free, i.e.
The plate material is a FGM with constant Poisson's ratio m = const, whilst the shear modulus is assumed to vary exponentially through the thickness from G 1 , the value at the lower surface to G 2 , the value at the upper surface, according to: where c is the inhomogeneity parameter, which can be expressed in terms of the ratio between the shear modulus at the top and bottom surfaces of the plate. The problem stated above will be solved using two methods, a 3-D elasticity solution and Finite Element solution, their results compared and conclusions drawn.

3-D elasticity solution
Using the displacement function method employed by Kashtalyan [22], the following representation of displacement and stress fields in a functionally graded plate has been obtained: For any pair of m and n, A j,mn are sets of six arbitrary constants to be determined from Eqs. (3)-(5), U 1,jmn , U 2,jmn , U 3,jmn , P 33,jmn , P 13,jmn , P 23,jmn , P 11,jmn , P 22,jmn and P 12,jmn are functions given in Appendix A. Boundary conditions at the edges of the plate, Eqs. (1), are satisfied exactly.
As previously mentioned, a number of different loads can be investigated, as long as they can be expressed as a sum of Fourier series. For a UDL as seen in Fig. 3 Q ðx 1 ; Substituting Eq. (7a) into Eq. (2c) and integrating yields For a patch of load of area 2c Â 2d (see Fig. 4) a similar method employed for the UDL can be carried out just over a different range of integration Substituting Eq. (8a) into Eq. (2c) and integrating gives

Finite Element analysis
The Finite Element method used in the present paper uses elements with stiffness gradient in the thickness direction. The method of element formulation is discussed in detail in a number of papers [20,21,23,24] and an overview is given here. Beginning with an assumed set of shape functions it can be written that the displacements within an element are interpolated as where N i (x) is a matrix of shape functions corresponding to each of the n nodes of the element and U i are nodal displacements corresponding to each of the n nodes.
To obtain the element strains differentiation of the displacements is carried out, hence where B i (x) is a matrix of derivatives of the shape functions. For linear elastic behaviour standard stress strain relations can now be used where C(x) is the material property matrix. Traditionally this matrix contains constant material properties but for a Functionally Graded Material can be set to spatially variable functions, for example exponential variation in Young's modulus. The element stiffness matrix k e , mapping the nodal displacements U i , to the nodal forces f i , can be written in the same way as for a standard Finite Element, that is Now applying the theorem of virtual work which states that the work done by nodal forces must be equal to the work of deformation within the element, allows the following to be written where V e is the volume of the element. Analytical integration of this expression is very difficult or indeed impossible, so instead a numerical integration scheme is applied. Through application of Gauss quadrature (see [25]), Eq. (12) can be calculated in the following manner where W i , W j , W k are Gauss weights (specified in [26]), i, j and k are the Gaussian integration points and J ijk is the determinant of the Jacobian matrix. In carrying out this integration, the constitutive matrix is evaluated at each Gaussian integration point and since the applied loading is known then displacements can be found through application of Commercial Finite Element software such as ABAQUS does not support graded elements directly, they have to be defined separately in a user subroutine. In ABAQUS this is carried out through a UMAT subroutine (see [27]). This subroutine is called at all material calculation points of elements for which the material definition includes a user-defined material and is used to define the mechanical constitutive behaviour of the material and provide the material Jacobian matrix, @Dr @De . Most importantly the stresses and solution-dependent state variables are calculated and updated to their values at the end of the increment for which the subroutine is called. Due to the symmetry of the problem and of the applied loads only one quarter of the plate is modelled, with 50 graded 8 node quadrilateral elements in the thickness direction required for convergence to the analytical solution.

Validation
Validation of the analytical solution is given through comparison with elasticity solutions available in literature. Reddy [28] gives numerical results for the transverse deflections and stresses in an isotropic homogeneous plate subject to uniform distributed loading. In order to reproduce the results of Reddy, the inhomogenity ratio c is set sufficiently close to zero and the geometry is fixed as a/h = 10. Table 1 shows normalised maximum transverse deflections and stresses in square and rectangular plates under sinusoidal loading Q(x 1 , x 2 ) = Àq 0 sin(px 1 /a)sin(px 2 /b) and UDL loading of intensity q 0 . Following Reddy, stresses and displacement are normalised as follows: There is good agreement between results based on the present analytical solution and those obtained by Reddy [28]. To achieve convergence, up to 29 terms in m and n are required for UDL.
The developed Finite Element model for graded materials is validated through comparison with the available three-dimensional solutions for isotropic functionally graded plates. Table 2 shows normalised mid-plane displacements u 3 ¼ G 2 u 3 q 0 h at the centre of a thick square a/h = b/h = 3 isotropic graded plate with exponential variation of the shear modulus through the thickness based on the present Finite Element solution and the 3-D elasticity solution. The plate is simply supported on its edges and loaded by transverse loading Q(x 1 , x 2 ) = Àq 0 sin(px 1 /a)sin(px 2 /b) at the top surface and results are given for a range of inhomogeneity parameters. The Finite Element solution is found to be in excellent agreement with the 3-D elasticity solution, the difference in transverse displacements being at its largest 0.009%. In Table 3, the number of elements used through the thickness in the Finite Element analysis of the functionally graded plate subject to sinusoidal loading is shown and the normalised transverse displacement is considered. It was found that using 50 elements through the thickness of the plate provides excellent agreement with the 3-D elasticity solution throughout the plate.

Comparative study of loading types
In this study the stress and displacement fields under the different loading types outlined above are examined. The stiffness gradient is chosen as G 2    In order to study the effect of load localisation, the intensity of the load for each of the different load types must be modified such that the overall load applied is the same in all cases. The factors by which the loading intensity must be increased for each of the different loading types can be calculated using double integrals and are specified in Table 4. The stress and displacements fields in a functionally graded plate (FGM) and a homogenous plate (H) plate are compared under each of the loadings outlined above, using both the 3-D elasticity analysis and Finite Element model with graded elements. Both plates are thick, with a/h = b/h = 3. It can be seen from Fig. 5, that the through thickness variation of normalised out-of-plane normal stress is similar for both the FGM and H plates, showing that grading has almost no effect on this stress component. However when load concentration ( Fig. 5B and C) is increased, corresponding increases in this stress component are observed. Excellent agreement is seen between the 3-D elasticity solution and the Finite Element model. Plots of through thickness variation of transverse shear stress r 13 (Fig. 6A-C) show that under all loading types, the stress in the upper half of the FGM plate is higher than in the H plate, whilst the stress in the lower half of the FGM is lower. There  is again excellent agreement between the 3-D elasticity solution and the Finite Element analysis.
Plots of normalised in-plane normal stress r 11 (Fig. 7A-C) and normalised in-plane shear stress r 12 (Fig. 8A-C) show that stresses at the upper surface of plate are higher for FGM plate, whilst lower at the bottom surface. As the load concentration is increased the magnitude of these stress components in the upper half of the plate substantially increase, with a much smaller increase observed for the lower half of the plate. Analysis of the through-thickness variation of the in-plane displacements u 1 (Fig. 9A-C) and transverse displacement u 3 ( Fig. 10A-C) show that as the load becomes more concentrated, displacements through the plate increase. When comparing FGM and H plates, it can be seen that the increased stiffness of the graded plate provides a reduction in both displacements. The plot   of in-plane displacement is highly non-linear, once again emphasising the importance of the 3-D stress analysis.

Conclusions
In the current paper, a simply supported functionally graded plate subjected to uniformly distributed and patch loading has been analysed using 3-D elasticity and a Finite Element model with graded elements. Results were then compared with an equivalent homogenous plate. Agreement between the 3-D elasticity solution and Finite Element model was excellent, with both methods showing that the increased stiffness provided by the Functionally Graded Material reduces displacements throughout the plate. Many of the plots produced were highly non-linear through the thickness, showing the importance of 3-D stress analysis.

Acknowledgement
Financial support of this research by EPSRC DTA is gratefully acknowledged.