FGMs are a mixture of ceramic and metal or a combination of different metals made by gradually varying the volume fraction of the constituent materials. The FGM can be continually produced by varying the constituent multi-phase materials in a predetermined profile. The most distinct features of an FGM are its non-uniform microstructures within its continuously graded properties. FGMs can be described by variations in volume fractions, exponential function, or sigmoid function. Most researchers use the power-law function in which the material properties gradually varied from a bottom surface of pure metal to a top surface of pure ceramic, as shown in the section of the plate in Figure 2. To acquire the effective thermal and mechanical properties of an FGM plate, a simple mixture rule based on a power-law distribution for which the material properties vary throughout the thickness of the plate was assumed as

  $P(z)=P_{t}V+P_b(1-V)$                                                                                                                (1-a)

where subscripts m and c referred to the metal and ceramic constituents, respectively. P denoted a property of the material, and z was measured over the thickness of the plate. In Eq. (1-a), V was the volume fraction of the ceramic phase, and it was represented by a power-law distribution such that

   $V={(\frac{z}{h}+\frac{1}{2})}^n\begin{array}{cc}& \end{array},-\frac{h}{2}\le z\le \frac{h}{2}$                                                                                                      (1-b)

In which h was the thickness of the plate, n was the power-law index, which is always positive, and z was measured from the geometric center of the plate, in which -(h/2) < z < +(h/2). Figure 3 shows the variation of the ceramic volume fraction throughout the thickness of the plate for different volume fraction potencies. The plate is fully ceramic when n = 0, while the composition of ceramic and metal takes a linear shape when n = 1.

Figure 3 Variation of volume fraction over the thickness for different values of power index n.

Mathematical formulation the floor and ceiling FGM Plates

The final differential equation for the FGM plate, considering the loads on the plate, (q), was obtained as follows:

                                                                                $\overline{D}{\nabla }^4{w}_0=q$                                                            (2)

where D ̅ was the flexural rigidity of the FGM plate, and w was the displacement of the plate along the z-axis.

Figure 3 Variation of volume fraction over the thickness for different values of power index n.

Displacement depends on bending rigidity and the amount of uniform loading parameters perpendicular to the plate. Assuming a constant uniform load value on the floor and ceiling tank plates, the bending rigidity parameter of plates was representative of the floor and ceiling under static analysis. This section presents the flexural rigidity of the FGM floor and ceiling plates within a cubic tank.

                                                                   $\overline{D}=\frac{I_1{I}_3-I_2^2}{I_1(1-{\upsilon }^2)}$                                                                   (3)

The flexural rigidity of FGM floor plates, where the top surface of the plate faced the load and was made of ceramic, was calculated using Eq. (4) through (6):

$I_1=E_{m}h+(E_C-E_m)(\frac{h}{n+1})$                                                                                                                    (4)

$I_2=(E_C-E_m)h^2(\frac{1}{n+2}-\frac{1}{2n+2})$ (5)

$I_3=E_m\frac{h^3}{12}+(E_C-E_m)h^3(\frac{1}{n+3}-\frac{1}{n+2}-\frac{1}{4n+4})$ (6)

However, the flexural rigidity of FGM ceiling plates, where the top surface of the plate faced the load and was made of metal, can be calculated using Eq. (7) through (9):

$I_1=E_{c}h+(E_m-E_c)(\frac{h}{n+1})$ (7)

$I_2=(E_m-E_c)h^2(\frac{1}{n+2}-\frac{1}{2n+2})$ (8)

$I_3=E_c\frac{h^3}{12}+(E_m-E_c)h^3(\frac{1}{n+3}-\frac{1}{n+2}-\frac{1}{4n+4})$ (9)

where Ec and Em represented the modulus of elasticity of ceramic and metal, respectively.  

 As illustrated in Figure 4, the bending rigidity of FGM plates appeared as follows:  the metal part of the plate was Aluminum $E_m=70\begin{array}{c}\end{array}GPa$  and the ceramic part of the plate was Alumina $E_C=380\begin{array}{c}\end{array}GPa$ . The Poisson's ratio of the plate was equal to 0.3 $(\nu =0.3)$ . The thickness of the plate was 0.01 of the length of the plate. The n value varied between zero and one hundred. Figure 4 shows that by increasing the value of parameter n, the flexural rigidity of the FGM plate on the floor of the cubic tank was more than that of the ceiling plate.

Figure 4 Variation of flexural rigidity of the FGM plate used in ceiling and floor of tanks for different values of ‘n’.

Interpretation from $(\nu =0.3)$  Figure 4 shows that until n = 1, the flexural rigidity of the FGM plate in the ceiling plate of the cubic tank was lower than that of the floor plate. However, if n was more than 1, the flexural rigidity of the ceiling FGM plate was more than the bending rigidity of the floor plate.

Galerkin method for loading perpendicular to the mid - Plane

In this study, evaluation was based on an FGM rectangular plate, in the thickness as shown in , which was under constant uniform load on the surface of the plate perpendicular to the mid-plane with variation in the flexural rigidity parameter, and with fully simply supported boundary conditions for both the longitudinal and transverse edges.

Plates’ behavior analysis can be applied to trigonometric functions in the longitudinal direction of the plate that satisfies boundary conditions. To satisfy the boundary conditions in the transverse direction of the plate, a polynomial function was employed due to the difficulty in meeting the requirements of trigonometric functions. Simply supported (denoted as S) and clamped supported (denoted as C) plates are illustrated in the following figures.

Figure 5 Model rectangular plate under uniformly load perpendicular to the mid-Plane.

To solve the problem and maintain clarity of presentation, Eq. (2) was non-dimensionalized by substituting

   $\zeta =\frac{x}{a}$ $\eta =\frac{2y-b}{b},$ $\phi =\frac{a}{b},$ $W=\frac{w}{h}$                                                                                      (10)

where a and b were the length and width of the plate, respectively; W was the non-dimensional form of the displacement of the middle surface of the plate in a direction normal to the undeformed middle surface.

Thus, Eq. (2) was re-written as follows:

                   $(\frac{1}{{\phi }^2}\frac{{\partial }^{4}w}{\partial {\xi }^4}+8\frac{{\partial }^{4}w}{\partial {\xi }^2\partial {\eta }^2}+16{\phi }^2\frac{{\partial }^{4}w}{\partial {\eta }^4})-\frac{q}{\overline{D}}{a}^2{b}^2=0$                                                    (11)

The solution of Eq. (11) is taken in the form                    

$w={\displaystyle \sum _m{\displaystyle \sum _n{q}_{mn}{f}_m(\xi }})g{}_n(\eta )$ (12)

Eq. (12) was simplified to yield a fourth-order partial differential equation with variable coefficients for which a closed form solution would be difficult, if not impossible, to obtain. Therefore, Galerkin’s series method was employed to solve the problem.

The function precisely described the deflection of the transverse direction of the plate.

                   $f_m(\xi )=\{\begin{array}{l}\sin (m\pi \xi )\\ \sin [(m+1)\pi \xi ]+(\frac{m+1}{m})\sin (m\pi \xi )\\ \sin (m\pi \xi )\sin (\pi \xi )\\ \sin [(m-\frac{1}{2})\pi \xi ]\sin (\frac{\pi \xi }{2})\end{array}\begin{array}{c}\\ \\ \end{array}\begin{array}{c}\begin{array}{l}S-S\\ \end{array}\\ \begin{array}{l}S-C\\ \end{array}\\ \begin{array}{l}C-C\\ \end{array}\\ C-G\end{array}$                                                    (13)

The function was taken as polynomial with four constants (four boundary conditions), as shown in Eq. (14).

                                                               $g_n(\eta )={\eta }^{n+4}+A_n{\eta }^{n+3}+B_n{\eta }^{n+2}+C_n{\eta }^{n+1}+D_n{\eta }^n$            (14)

These constants were evaluated through applying the equations in Table 1.

To get the deflection function corresponding to the number of series, the constants were substituted into Eq. (12)

$w={{\displaystyle \sum _{n=1,2,3}{q}_n{f}_m(\xi )[\eta }}^{n+4}+A_n{\eta }^{n+3}+B_n{\eta }^{n+2}+C{}_n{}^{}\eta {}^{n+1}+D_n{\eta }^n]$ (15)

For simplicity

$w=Y.f_m(\xi )$ (16)

Where   

$Y={\displaystyle \sum _n{q}_n[g_n(\eta )}]$ (17)

The number of half waves, m in Eq. (13), depended on the aspect ratio of the plate. Indeed, m indicated the mode shape in which the plate buckles. Substituting Eq. (16) into Eq. (11) leads to an expression simplified to a fourth order differential equation. The final form of the governing differential equation is

$\overline{S}=16\frac{d^{4}Y}{d{\eta }^4}f(\xi )+\frac{8}{{\phi }^2}\frac{d^{2}Y}{d{\eta }^2}\frac{d^{2}f(\xi )}{d{\xi }^2}+\frac{1}{{\phi }^4}Y\frac{d^{4}f(\xi )}{d{\xi }^4}-a^2{b}^2\frac{q}{\overline{D}}Y=0$ (18)

${\displaystyle \underset{0}{\overset{1}{\int }}{\displaystyle \underset{-1}{\overset{1}{\int }}[\overline{S}}}]\frac{\partial [Y\times f(\xi )]}{\partial q_n}d\eta d\xi =0$ (19)

By applying the Galerkin method and using weight functions, differential Eq. (19) can be solved

$$\begin{array}{l}{\displaystyle \sum _{n=0,1,\mathrm{..}}^j{q}_n}{\displaystyle \underset{0}{\overset{1}{\int }}}{\displaystyle \underset{-1}{\overset{1}{\int }}[16(f_m}(\xi ){)}^2\frac{d^4{Y}_n}{d{\eta }^4}+\frac{8}{{\phi }^2}{f}_m(\xi )\frac{d^2{f}_m(\xi )}{d{\xi }^2}\frac{d^2{Y}_n}{d{\eta }^2}+\frac{1}{{\phi }^4}{f}_m(\xi )\frac{d^4{f}_m(\xi )}{d{\xi }^4}{Y}_n\\ -\frac{q}{\overline{D}}{f}_m(\xi )a^2{b}^2{Y}_n]Y_{j}d\xi d\eta =0\end{array}$$                                                                                                                                             (20)

Equation (20) leads to a square matrix as follows:

$\left[\begin{array}{cccccc}(L_{00}-K{M}_{00})& (L_{01}-K{M}_{01})& .& .& .& (L_{0J}-K{M}_{01})\\ (L_{10}-K{M}_{10})& .& .& .& .& .\\ .& .& .& .& .& .\\ .& .& .& .& .& .\\ .& .& .& .& .& .\\ (L_{J0}-K{M}_{J0})& .& .& .& .& (L_{JJ}-K{M}_{JJ})\end{array}\right]\left[\begin{array}{c}{q}_0\\ .\\ .\\ .\\ .\\ q_J\end{array}\right]=0$                                          (21)

In which

$L_{ij}={\displaystyle \underset{0}{\overset{1}{\int }}}{\displaystyle \underset{-1}{\overset{1}{\int }}[16}{(f_m(\xi ))}^2\frac{d^4{Y}_i}{d{\eta }^4}+\frac{8}{{\phi }^2}{f}_m(\xi )\frac{d^2{f}_m(\xi )}{d{\xi }^2}\frac{d^2{Y}_i}{d{\eta }^2}+\frac{1}{{\phi }^4}{f}_m(\xi )\frac{d^4{f}_m(\xi )}{d{\xi }^4}{Y}_i]Y_{j}d\xi d\eta$ (22)

In this relationship, matrix [L] represented the stiffness matrix for the plates. The second part of Eq. (20) or the force was

           $f={\displaystyle \underset{0}{\overset{1}{\int }}}{\displaystyle \underset{-1}{\overset{1}{\int }}\frac{q}{\overline{D}}}{f}_m(\xi )a^2{b}^2{Y}_n{Y}_{j}d\xi d\eta$                                                                                     (23)

Also:

$f=k\Delta =[L]\times \Delta$ (24)

Displacement of the plate was calculated with Eq. (25)  

                                                  $\Delta ={[L]}^{-1}\times \left[f\right]$      (25)

Finally, displacement of the plate was obtained using Eq. (26)

     $\Delta =w={[L]}^{-1}\times {\displaystyle \underset{0}{\overset{1}{\int }}}{\displaystyle \underset{-1}{\overset{1}{\int }}\frac{q}{\overline{D}}}{f}_m(\xi )a^2{b}^2{Y}_n{Y}_{j}d\xi d\eta$                                                                             (26)