1 INTRODUCTION

The intensity of stress concentrations and stress effects due to the large mismatch in material properties can be substantially reduced if the microstructure transition behaviour is graded. Advances in material synthesis technologies have spurred the development of functionally graded materials (FGM) with promising applications in aerospace, transportation, energy, cutting tools, electronics, and biomedical engineering (^{Chakraborty et al., 2003}). An FGM comprises a multi-phase material with volume fractions of the constituents varying gradually in a predetermined profile, thus yielding a non-uniform microstructure in the material with continuously graded properties (^{Jin et al., 2003}).

Elastic and elastic-plastic analyses of thick-walled pressure vessels have always attracted a lot of research interest because of their importance in engineering applications. One of the first studies in analytical elastic solutions for functionally graded cylindrical vessels was presented by (^{Horgan and Chan, 1999}). (^{Figueiredo et al., 2008}) proposed a numerical methodology in order to predict the elastic-plastic stress behaviour of functionally graded cylindrical vessels subjected to internal pressure. It was assumed that the structures undergo small strain and that the material properties of the graded layer are modeled by the modified rule of mixtures approximation. Furthermore, the plastic domain for ductile phases was defined through the von-Mises yield criterion. They proposed an iterative method for solving the nonlinear system, combining a finite element approximation and an incremental-iterative scheme. (^{Haghpanah et al., 2009}, ^{2010}) extended the Variable Material Property (VMP) method developed by (^{Jahed and Dubey, 1997}) for materials with varying elastic and plastic properties. In the VMP method, the linear elastic solution of a boundary value problem is used as a basis to generate the inelastic solution. Through iterative analysis, the VMP method was used to obtain the distribution of material parameters which were considered as field variables. The application of the VMP method, generally applied to homogeneous elastic-plastic materials (^{Jahed et al., 2001}, ^{2005}, ^{2006}), was extended to materials with varying elastic-plastic properties in order to calculate the residual stresses in an autofrettaged FGM cylindrical vessel.

Although there are several papers in the elastic analysis of FGM spherical pressure vessels in the literature (^{You et al., 2005}; ^{Dai et al., 2006}; ^{Chen and Lin, 2008}), elastic-plastic stress analysis of FGM spherical pressure vessel is not such a customary study. (^{Sadeghian and Ekhteraei, 2011}) studied the thermal stress analysis for an FGM spherical pressure vessel made of elastic-perfectly plastic and power law material model.

Similar FGM cylindrical and spherical vessels, much of the studies on FGM rotating disk has been carried out numerically (^{Durodola and Attia, 2000}; ^{Bayat et al., 2008}). (^{Haghpanah et al., 2012}) applied VMP method to solve for estimating the elasto-plastic stresses in a rotating disk with varying elastic and plastic properties in the radial direction.

In this paper a new analytical method is proposed for predicting stress components of a strain-hardening cylinder based on the von-Mises yield criterion under plane-stress conditions by assuming an isotropic material model. Results obtained from finite element analysis using the commercial software, ABAQUS (v 6.10), were used to validate the proposed analytical method. The method was further extended to obtain solutions for FGM spherical vessels and rotating disks.

Elastic-plastic governing equations of functionally graded cylindrical and spherical pressure vessels and rotating disks are presented in the next section. Section 3 describes the material properties of the graded layer modelled by the modified rule of mixtures approximation. Solution procedures and results of the elastic-plastic analyses are presented in section 4. Finally in section 5 key conclusions are pointed out.

2 GOVERNING EQUATION FOR CASE STUDIES

2.1 Cylindrical Pressure Vessel

Consider the axisymmetric problem of a cylinder subjected to uniform pressure on the inner surface. The problem can be studied in polar coordinates (*r*,θ). Due to axisymmetric deformations, the only displacement component is *u _{r}
* =

*u*(

_{r}*r*), i.e. the radial displacement. The equilibrium equations with zero body force reduce to the following single equation:

Where σ*
_{r}
* and σ

*are the radial and hoop stresses respectively.*

_{o}Considering the geometric relations:

Where *E _{o}
* and

*E*are the hoop and radial strains in terms of the radial displacement,

_{r}*u*must satisfy the compatibility equation:

_{r}

Considering the quasi-static plasticity and small deformations hypothesis, the total strains may be split into elastic and plastic components as follows:

where *E*(*r*) and *v* are the Young modulus and the Poisson's ratio respectively. It is assumed that *v* is a material constant while *E*(*r*) varies with the position across the wall thickness of the cylindrical vessel.

The governing equation for the functionally graded cylindrical vessel is obtained by substituting Eq. ^{(4)} into Eq. ^{(3)}. This equation shows the relationship between stress and plastic strain (Eq. ^{(5)}).

Stress and plastic strain in the cylindrical vessel may be obtained by solving the differential equations ^{(1)} and ^{(5)} simultaneously. However solving this system of differential equation is not possible because the number of variables (σ*
_{r}
* , σ

_{θ}, ε

*, ε*

_{r}^{p}_{θ}

*) is more than the number of equations. Substituting plastic strains in terms of stresses, this system of differential equations can be solved. For a von Mises material with associated flow rule, plastic strain increment is defined as (*

^{r}^{Dunne, 2006}):

where dε*
_{e}
^{p}
* is the equivalent plastic strain increment. σ

*, the equivalent stress and*

_{e}*S*, the deviatoric stress for plane stress are defined as:

and

Assuming small strains, the plastic stress-strain relation proposed by Hencky may be written as Eq. ^{(9)}, better known as 'total strain theory' (Chakrabarty, 2011). (^{Mendelson, 1968}) explained that the total strain theory is applicable to cylindrical problems (disks), and provides the results that match well with the actual material response.

For a linear strain hardening material (Fig. 1), the equivalent stress σ*
_{e}
* is determined by:

and the equivalent plastic strain is obtained from:

where *h _{p}
* (

*r*) is plasticity modulus (i.e. the gradient of the stress-plastic strain curve) and σ

*(*

_{yo}*r*) is the initial yield stress of the FGM material. Both

*h*(

_{p}*r*) and σ

*(*

_{yo}*r*) are functions dependent on the radial position,

*r*.

The plastic strain components can be determined by substituting Eqs. ^{(7)}, ^{(8)}, and ^{(11)} into Eq. ^{(9)}:

Furthermore, differentiating the hoop plastic strain with respect to the radius, *r*, gives:

Finally, the system of differential equations ^{(14)} is obtained by substituting Eqs. ^{(12)} and ^{(13)} into Eq. ^{(5)}:

By applying the appropriate boundary conditions, the distributions of radial and tangential stresses in the plastic region may be obtained from the above equations. The distribution of stress components in the elastic zone may also be obtained from Eqs. ^{(15)}.

2.2 Spherical Pressure Vessel

In order to obtain stress distributions for a spherical thick-walled functionally graded pressure vessel, the equilibrium equations with zero body force reduce to the following single equation in the spherical coordinate system (*r*, θ, Φ)):

where σ_{θ} and σ_{φ} are the circumferential stresses and σ*
_{r}
* is the radial stress component Considering the geometric relations:

where ε_{θ} and ε_{φ} are the circumferential strains and ε*
_{r}
* is the radial strain in terms of the radial displacement,

*u*. The compatibility equation must be satisfied:

_{r}

Following similar procedures to those explained in the previous section, i.e. considering the quasi-static plasticity and small deformations hypothesis, the total strains may be split into its elastic and plastic components as follows:

By substituting Eq. ^{(19)} into Eq. ^{(18)}, the governing equation for the functionally graded spherical vessel is obtained (Eq. ^{(20)}).

σ*
_{e}
* , the equivalent stress and

*S*, the deviatoric stress for spherical vessel are defined as

and

The plastic strain components can be determined by substituting Eqs. ^{(21)}, ^{(22)}, and ^{(11)} into Eq. ^{(9)}:

Furthermore, differentiating the hoop plastic strain with respect to radius, *r*, gives:

Finally, the system of differential equations ^{(25)} is obtained by substituting Eqs. ^{(23)}, and ^{(24)} into Eq. ^{(20)}:

The distributions of radial and circumferential stresses in the plastic region are obtained by applying the appropriate boundary conditions. Also the distribution of stress in the elastic zone may be obtained from Eq. ^{(26)}.

2.3 Rotating Disk

The equilibrium equation for a thin rotating disk as an axi-symmetric problem is:

where σ_{θ}, σ*
_{r}
* , ω and ρ(

*r*)are the tangential and radial stress components, the angular velocity of the disc and its mass density respectively.

Due to axisymmetric deformations, the only displacement component is the radial displacement, *u _{r}
* =

*u*(

_{r}*r*) and the strain field is given by:

where ε_{θ} and ε*
_{r}
* are the strains in tangential and radial directions. As before the radial displacement,

*u*must satisfy the compatibility equation:

_{r}

Also the total strain may be split into its elastic and plastic components as follows:

Substituting Eq. ^{(30)} into Eq. ^{(29)}, the governing equation for the functionally graded rotating disk is obtained.

The equivalent stress, σ*
_{e}
* , and the deviatoric stress,

*S*, for plane stress case in a rotating disk are the same equations as for the cylindrical vessel can be used. Therefore, the governing system of differential equations for the rotating disk can be obtained by substituting Eqs.

^{(11)},

^{(12)}and

^{(13)}into Eq.

^{(31)}:

One can also obtain the stress distribution in the elastic zone from Eq. ^{(33)}.

3 MECHANICAL BEHAVIOUR OF THE FUNCTIONALLY GRADED MATERIAL

It was assumed that the functionally graded metal-ceramic composite was locally isotropic and followed the von-Mises yield criterion. Three important material properties for elastic-plastic analysis are the elastic modulus *E*(*r*), the initial yield stress σ*
_{y0}
* (

*r*) and the tangent modulus

*H*(

*r*). These properties can be calculated using the modified rule of mixtures for composites (

^{Suresh, 1998}).

where the subscripts 'c' and 'm' indicate ceramic and metal material respectively. The volume fraction of ceramic particles is denoted by *f _{c}
* , and

*q*is the ratio of stress to strain transfer as follows, where σ

*, ε*

_{c}*and σ*

_{c}*, ε*

_{m}*are the average stress and strain in ceramic and metal respectively (Figure 2).*

_{m}

The empirical parameter *q* depends on many factors including material composition, microstructural arrangements, and the internal constraints. For example *q* → ∞ if the constituent elements deform identically in the loading direction, while *q =* 0 if the constituent elements experience the same stress level. In present analysis the ceramic particle reinforcement was assumed to have a volume fraction that varies from 0 at the inner radius, *r _{i}
* , to

*f*at the outer radius,

_{c0}*r*according to the following relationship:

_{o}

where *n* is the reinforcement distribution exponent (*n =* 0 denotes uniformly-reinforced metal-ceramic)

The material properties for each constituent phase are listed in Table 1. The parameter *q* may be approximated by calibration of experimental data from tensile tests performed on monolithic composite specimens. For example a value of *q = 4.5 Gpa* was used for a TiB/Ti FGM (^{Carpenter, 1999}) whereas the Poisson ratio was taken constant and equal to 0.3.

4 ANALYSIS OF RESULTS

By substituting Eqs. ^{(34)} - ^{(36)} into the system of differential equations ^{(14)}, ^{(25)} and ^{(32)}, a set of complicated equations are obtained that are difficult to solve parametrically. Therefore, numerical methods may be used as an alternative. ODE45 solver in (MATLAB, 2008) was used for this purpose. ODE45 is a function for the numerical solution of ordinary differential equation and employs variable step size Runge-Kutta integration methods. ODE45 varies the size of the step of the independent variable in order to meet the accuracy specified at any particular point along the solution. If ODE45 can take big steps and still meet this accuracy, it will do so and will therefore move quickly through regions where the solution does not change greatly. In regions where the solution changes more rapidly, ODE45 will take smaller steps. This solver uses a 4th and 5th order pair for higher accuracy (^{Jamshidi, 2001}).

To verify the accuracy of the analytical solution, finite element analyses were also performed using the commercial finite element code ABAQUS (^{Karlsson, 2008}). The conventional method of modelling FGM vessels or rotating disks in the software is to subdivide the thick wall into several thin layers with equal thicknesses. This method of modelling leads to a discontinuity in the mechanical properties of FGM materials and is both difficult and time-consuming. (^{Setoodeh et al., 2008}) proposed a new approach for analysing the FGM material in the elastic zone without the need for dividing the thickness of the cylinder into thin strips. They applied a virtual temperature distribution in the cylinder wall using the facility available within the software to assign continuously variable properties across the wall thickness and then create a one to one relationship between temperature and mechanical properties. Correlation of the distribution of temperature and mechanical properties of FGM material, allowed for the variation of the FGM properties in the cylinder to be modelled. Note that conductivity factor and other thermal parameters are zero, and temperature does not change in during the analysis time. Indeed, the analysis is mechanical (no thermo mechanical analysis). This method allowed for the analysis of the elastic-plastic FGM cylindrical and spherical vessels as well as the rotating disk. The three-dimensional 8-noded linear coupled temperature displacement family of finite elements in ABAQUS was used to model the cylinder. Mesh sensitivity analysis was also performed to ensure the results remained insensitive to the element size.

In order to evaluate the analytical method, a set of results from the finite element calculation and the VMP method (^{Haghpanah, 2009}) obtained for the plane strain conditions were compared with results obtained from the analytical method for an FGM cylindrical vessel subjected to an autofrettage internal pressure 100 MPa (Fig. 3). The results indicated that the proposed analytical method was capable of calculating stress components in a cylindrical vessel with great accuracy (in ODE45, the size of the step is variable and Number of solution points is 142 points. In the finite element calculation, Number of elements in the thickness of cylinder is 400 elements). The developed method was also used to study the stress components in a cylindrical FGM vessel with various ceramic particle reinforcement distributions. All the results are presented for the plane stress condition and the material properties listed in Table 1. Fig. 4 illustrates the distribution of von-Mises stresses across the thickness in the cylindrical vessel subjected to different levels of internal pressure where *n* = 2 and *f _{c0}
* = 0.8. The magnitude of von Mises stresses, indicate that with increasing the internal pressure, the plastically deformed region extends across the thickness of the cylinder from both the inner and the outer surfaces. In cylindrical vessels made of homogeneous material the plastic region grows only from the inner surface. Fig. 5 shows the von-Mises stress distribution in an FGM vessel subjected to an internal pressure of 300 MPa with

*f*= 1 and different reinforcement distribution exponents,

_{c0}*n*. By increasing

*n*, the plastic region gradually spreads from the inner surface of the cylinder.

For the spherical vessel subjected to an internal pressure of 500 MPa results obtained from the analytical method were compared with the finite element analysis results. Excellent agreement was observed as shown in Fig. 6 The distribution of von Mises stresses across the thickness of the spherical vessel subjected to different levels of internal pressure is shown in Fig. 7 Similar to the case of the cylindrical vessel, with increasing internal pressure, the plastically deformed region extends across the thickness of the spherical vessel from both the inner and the outer surfaces.

Similar to the case of the cylindrical and spherical vessels, Fig. 8 presents the excellent agreement of the analytical and finite element predictions of elastic-plastic stresses for a rotating disc. The distribution of von-Mises stresses across the thickness in a disk rotating at different angular velocities where *n* = 2 and *f _{c0}
* = 1 is presented in Fig. 9 Results clearly indicate that the growth of the plastic zone across the thickness is initiated from both sides of the disc.

5 CONCLUSIONS

Using a new analytical method, the elastic-plastic stress distributions in a cylindrical and spherical pressurized vessels and rotating disks made of an FGM material were determined. Solutions previously presented for this problem were based on numerical methods such as finite element analysis (^{Figueiredo et al., 2008}; ^{Durodola, 2000}) and the Variable Material Property (VPM) method (^{Haghpanah et al., 2009}, ^{2010}, ^{2012}) where as the technique used in this paper is based on the application of basic plasticity equations only. Despite some the limitations in the application of this method, the proposed analytical approach provided very efficient, yet accurate results for the determination of the elastic-plastic stress strain distribution in the cylindrical pressure vessels (similarly spherical vessels and rotating disks) of FGM materials. Also in this method, modelling and analysis time costs are greatly reduced compared to FEM because geometric modelling, meshing and assembling steps are eliminated. Other strength of this method includes the ability to modelling the isotropic strain-hardening rule based on the von-Mises yield criterion.

Finite element analysis of the problem using ABAQUS commercial code was used for verification of the proposed analytical solution technique. An alternative modelling approach was developed and used for this purpose. The numerical analysis within the software was performed by the application of a "virtual thermal load" that enabled the continuous variation of material behaviour through the wall thickness of the FGM cylinder in the finite element model. An accurate solution was therefore obtained from the FE analysis. This was achieved without the need to consider a multi-layered cylinder with a stepwise solution as used in most previous solutions in the existing literature (^{Haghpanah et al., 2009}, ^{2010}, ^{2012}; ^{Figueiredo et al., 2008}).

The analysis results obtained in this work also indicate the possibility of formation and growth of the plastic region within the wall thickness from the external surface of the FGM vessels or rotating disks whereas in a cylindrical (spherical) vessels and rotating disks made of homogeneous materials, plasticity starts essentially from the inner surface.