1 Introduction

Bending of solid circular plates is one of the classical problems in elasticity^[1]. Meanwhile, as a basic structural configuration, circular plates have been commonly used in infrastructures and other engineering applications. As a result, problems related to solid circular plates (including bending) have attracted a lot of research attention from both scientists and engineers. Moreover, functionally graded materials (FGMs) are a novel type of inhomogeneous materials whose composition and mechanical properties are graded along one or more directions. FGMs can be designed to fulfil various requirements for structural components in numerous engineering fields. Since the birth of this type of composite materials, lots of research works have been done to reveal the elastic behaviors of FGM plates under various loads and boundary conditions^[2-5].

In order to analyze the bending of FGM plates, simplified plate theories have usually been employed. For instance, using the relationships between the solutions of the classical plate theory (CPT) and the first-order shear deformation plate theory (FSDT), Reddy et al.^[6] studied the bending of FGM circular and annular plates subject to a uniform load. Analytical solutions were presented by Nosier and Fallah^[7] for both axisymmetric and asymmetric behaviors of FGM circular plates based on the FSDT. Using the third-order shear deformation plate theory (TSDT), Saidi et al.^[8] investigated the axisymmetric bending and buckling problems of FGM circular plates. Heydari et al.^[9] explored two new approaches for buckling analysis of FGM circular plates having linear and quadratic thickness variations resting on Pasternak elastic foundation with simply-supported and clamped boundary conditions; CPT and higher-order shear deformation theory were employed and compared.

Besides the above-mentioned methods based on various plate theories, exact analytical solutions directly based on three-dimensional (3D) elasticity theory have also been sought since they can be used as benchmarks for assessing the solutions obtained from plate theories or numerical methods. However, only a few relevant problems have been solved. Li et al.^[10] derived the elasticity solution of axisymmetric bending of FGM circular plates subject to even-order polynomial loads. Using the direct displacement method, Wang et al.^[11-12] studied bending responses of transversely isotropic FGM circular elastic and magneto-electro-elastic plates under arbitrary axisymmetric transverse loads. Lu et al.^[13] obtained an elasticity solution for transversely isotropic FGM circular plates subject to axisymmetric loads by the extended Plevako's general solution. Numerical or semi-analytical methods based on the elasticity theory have also been proposed or adopted in the literature. For instance, Wu and Liu^[14] recently presented a state-space meshless method to analyze the axisymmetric bending of FGM circular plates with simply-supported and clamped edges. The main idea behind the state-space meshless method is similar to that of the state-space-based differential quadrature^[15-16]. The latter has been employed by Alibeigloo^[17] to investigate axisymmetric thermoelastic responses of sandwich circular plates with an FGM core imbedded in metal and ceramic facing sheets. To the best of the authors' knowledge, however, the general asymmetric behavior of FGM solid circular plates based on 3D elasticity theory has not been undertaken yet.

In the authors' previous work^[18], England's method^[19] was generalized from isotropic materials to transversely isotropic ones, which was then used to solve the bending problem of FGM rectangular plates. In the present study, the generalized England's method is further utilized to study the bending of transversely isotropic FGM solid circular plates subject to biharmonic loads, which are expanded in Fourier series. For each item in the series, the four analytic functions α (ζ), β (ζ), ϕ (ζ), and ψ (ζ) that are the key in the generalized England's method are constructed analytically. Consequently, 3D analytical solutions are presented for a transversely isotropic FGM solid circular plate subject to axisymmetric and asymmetric loads, as will be illustrated in the paper.

2 Formulation

Employing a Cartesian coordinate system, the differential equation of equilibrium in the absence of body forces is

(1)

in which the subscript comma denotes differentiation with respect to the coordinate variable that follows.

The relations between components of stress and displacement for transversely isotropic materials are given as^[20]

(2)

where u, v, and w are the displacement components. c_ij (2c₆₆=c₁₁-c₁₂) are the elastic constants. As for FGMs, they are functions of z, namely, c_ij=c_ij (z). If c₁₁=c₃₃, c₁₂=c₁₃, and c₄₄=c₆₆, the material will degenerate to the isotropic ones. The z-axis coincides with the material axis of symmetry, and is perpendicular to the mid-plane (i.e., the xy-plane) of the plate.

The following forms of solutions to Eqs. (1) and (2) are adopted^[19]:

(3)

where R_j (j=0, 1, ..., 4) and T_k (k=1, 2, ..., 4) are functions of z; , , and are defined as the mid-plane displacements,

(4)

Consider the following stress conditions on the surface of the plate: σ_zx=σ_zy=0 at z=±h/2, σ_z=0 at z=-h/2, and σ_z=-p(x, y) at z=h/2, where p is a biharmonic load (i.e., ∇²∇²p=0). Substituting Eq. (3) into Eq. (2), the result in turn into Eq. (1), and making use of the above mentioned stress boundary conditions, we can determine the expressions of functions R₀, T₁, T₂, R₁, R₂, T₃, R₃, T₄, R₄. The readers can refer to Appendix A of Ref. [18]. Consequently, the complete solution of the mid-plane displacements is^[18]

(5)

(6)

where a prime denotes derivative with respect to the complex variable ζ=x+iy, , α (ζ), β (ζ), ϕ (ζ), and ψ (ζ) are four analytic functions which are the key in the generalized England's method^[18], κ₁, κ₂, κ₃, and κ₄ are constants, and

(7)

where S₂₁ is a constant; Q(ζ) and P(ζ) are two analytic functions related to the load p(x, y) through the following function:

(8)

The corresponding 3D displacement and stress components expressed by four analytic functions α (ζ), β (ζ), ϕ (ζ), and ψ (ζ) can then be given by substituting Eqs. (5) and (6) into Eqs. (3) and (2),

(9)

(10)

(11)

(12)

(13)

(14)

3 FGM circular plates subject to transverse loads

Figure 1 shows an FGM solid circular plate subject to a transverse load with the radius a and thickness h. In the cylindrical coordinate system (r, θ, z), take u_r, u_θ, and w as the displacement components and σ_r, σ_θ, σ_z, σ_rθ, σ_rz, and σ_θz as the corresponding stress components.

According to the Fourier series expansion law, the load p(r, θ) can be expanded in the circumferential direction as follows:

(15)

where p₀=T₀ (r), and .

The following equation can be obtained by letting Eq. (15) satisfy the biharmonic equation:

(16)

Let T_k (r)=a_kr^ℓ and substitute it into the above equation. We can get

(17)

Consequently, the solution of T_k (r) is

(18)

Notice that for a solid circular plate, the general solution to (18) corresponding to k=0 and k=±1 takes the following forms:

(19)

where a₀₁ and a₀₃ are real constants, and a_kj (j=1, 2, 3, 4; k=±1, ±2, ...) are complex constants.

As a result, each kind of loads p_k (r, θ) can be derived which constitutes a class of periodic solutions in the θ direction satisfying the biharmonic equation (see Ref. [1]). After that, using the principle of superposition, bending of transverse loads p(r, θ) can be studied by solving the response corresponding to p_k (r, θ). In what follows, the expressions of the functions α (ζ), β (ζ), ϕ (ζ), ψ (ζ) and the mid-plane displacements will be derived for different k.

Case 1   k=0

This case corresponds to the axisymmetric loading, and we have

(20)

Comparing Eq. (20) with Eq. (8) leads to

(21)

Integrating Eq. (21) gives

(22)

where

(23)

Substituting Eq. (22) into Eq. (7) leads to

(24)

where

(25)

Let

(26)

where α₀, β₁, and ϕ₁ are the unknown real constants which can be determined from the cylindrical boundary conditions at r=a for the circular plate. Substitution of Eqs. (26) and (24) into Eqs. (5) and (6) yields

(27)

(28)

where the expression of the function D₀ (r) is given in Appendix A. It is found from Eqs. (27) and (28) that and are real functions of r and .

Case 2   k=1

In this case, we have

(29)

Substituting the second term of Eq. (19) into Eq. (29) leads to

(30)

Comparing Eq. (30) with Eq. (8) leads to

(31)

which can be integrated to

(32)

where

(33)

Substituting Eq. (32) into Eq. (7) leads to

(34)

where

(35)

Let

(36)

where α₁, β₂, ϕ₂, and ψ₀ are unknown complex constants which can also be determined from the cylindrical boundary conditions at r=a for the circular plate. The corresponding expressions of mid-plane displacements are

(37)

(38)

in which the expressions of the functions B₁ (r), C₁ (r), and D₁ (r) can be seen in Appendix A.

Case 3   k=2

In this case, we have the following expression from Eqs. (15) and (18):

(39)

Comparison of Eq. (39) with Eq. (8) leads to

(40)

The integration of the above equation gives rise to

(41)

where

(42)

Substituting Eq. (41) into Eq. (7) leads to

(43)

where

(44)

Let

(45)

where α₂, β₃, ϕ₃, and ψ₁ are unknown complex constants. The corresponding expressions of mid-plane displacements are

(46)

(47)

in which the expressions of the functions B₂ (r), C₂ (r), and D₂ (r) are listed in Appendix A.

Case 4   k=3, 4, 5, ...

Based on Eqs. (15) and (18), the following expression can be given:

(48)

for k (k=3, 4, 5, ...). Comparing Eq. (48) with Eq. (8) leads to

(49)

Integrating Eq. (49) gives

(50)

where

(51)

Substitution of Eq. (50) into Eq. (7) yields

(52)

where

(53)

Let

(54)

where α_k, β_k+1, ϕ_k+1, and ψ_k-1 are unknown complex constants. The corresponding expressions of mid-plane displacements are

(55)

(56)

where the expressions of functions B_k (r), C_k (r), and D_k(r) can be found in Appendix A.

4 Determination of unknown constants

The unknown real or complex constants involved in the four analytical functions corresponding to every kind of biharmonic loads can be determined by virtue of the boundary conditions at r=a for a solid circular plate.

4.1 Resultant forces in cylindrical coordinates

Integrating Eqs. (11)-(13) gives the expressions of the resultant forces and moments,

(57)

(58)

(59)

where a_k (k=1, 2, ..., 9), b_j (j=0, 1, ..., 9), Q_z1, Q_z2, Q_z3, and Q_z4 are real constants (see Appendix B of Ref. [18].

4.2 Boundary conditions

For an FGM solid circular plate, there are two typical boundary conditions at r=a, namely, simply supported (S) and clamped (C):

(60)

(61)

In the following, the case of k=1 is taken as an example to illustrate the procedure of determining these constants in the four analytical functions.

It can be found from Eq. (37) that the following equations are equal to the conditions of and :

(62)

As for the conditions of , we can obtain the following equivalent equations from Eq. (38):

(63)

The equivalent equation of the condition of is

(64)

The transform relations between Cartesian coordinates and cylindrical coordinates for the resultant forces will be used as follows:

(65)

(66)

(67)

(68)

(69)

where the expressions of N₁ (r), N_-1 (r), N_-3 (r), M₁(r), M_-1 (r), M_-3 (r), Q₀ (r), and Q_-2 (r) are supplemented in Appendix B.

It can be found from Eqs. (65) and (66) that the following equation is equal to the condition of N_r(a, θ)=0:

(70)

The equivalent equation of the condition of M_r(a, θ)=0 is obtained from Eqs. (67) and (68) as follows:

(91)

Based on the above equations obtained for different boundary conditions, the involved unknown constants α₁, β₂, ϕ₂, and ψ₀ can be fixed. Thus, the 3D components of displacement and stress at any position in the solid circular plate can be determined for Case 2.

5 Numerical examples Example 1 An isotropic FGM solid circular plate subject to uniformly distributed load

Consider an isotropic FGM solid circular plate subject to uniform load q in this example to validate the present method. This particular case corresponds to case k=0 and the corresponding elasticity solution can be obtained by letting a₀₁=q and a₀₃=0 in Eq. (20). The same problem was studied by Reddy et al.^[6], Li et al.^[10], and Wang et al.^[11], for which the Poisson's ratio ν=0.288 remains unchanged, and the elastic modulus varies along the thickness in the following pattern:

(72)

where E_m=110.25 GPa and E_c=278.41 GPa. The parameter λ is the graded factor. When λ=0, the material becomes the homogeneous one.

The following deflection factor w₀ is introduced for numerical illustration:

(73)

Unless stated otherwise, a=1 m, h=0.15 m, and q=10⁶ N/m² are assumed in this and the next numerical examples. Numerical results at point r=0 and z=0 are presented and compared in Tables 1 and 2 for the simply supported and clamped plates, respectively. It is found that the present solution agrees well with those available in the literature. Compared with those obtained based on the FSDT^[6], the present solution is closer to the 3D elasticity solution^[10-11], thus validating the present method. In Table 2, the deflection factor w₀ of the plate with this kind of clamped boundary is independent of the thickness to radius ratio h/a which has been revealed in Ref. [10]. It is emphasized that letting a₀₁=0 and a₀₃=q in Eq. (20), we can further solve the axisymmetric problem with transverse load qr².

Example 2 A transversely isotropic FGM solid circular plate subject to load qr cos θ

A transversely isotropic FGM circular plate subject to load qr cos θ is investigated in this example to show the asymmetric bending behavior of the plate. The corresponding elasticity solution can be obtained by letting a₁₁=q/2 and a₁₄=0 in Eq. (29). The following FGM model is considered:

(74)

where C_ij^0(A) are those of Al₂O₃ at the bottom surface (z=-h/2), and C_ij^0(T) are those of Titanium at the top surface (z=h/2) whose values are listed in Table 3. It can be observed that the bending rigidity of the plate decreases with λ^[21].

Table 4 lists the dimensionless deflection w₁=w(a/2, π/4, z)/h and normal stress σ_r^*=σ_r (a/2, π/4, z)/q of the simply-supported circular plate subject to load qr cos θ. It can be observed from the results that there is little difference in the deflection along the thickness direction of the plate, regardless of the value of λ. It is also confirmed that the deflection increases with λ, just as expected. Compared with the homogeneous plate, the compressive zone of the FGM plate slightly extends to the lower half of the plate and the tensile stress on the bottom surface (z=-h/2) of the plate is far greater than the compressive stress at the top surface (z=h/2) of the plate.

Figure 2 plots the distribution of the dimensionless deflection (w(a/2, 0, z)/h)×10⁴ along the thickness direction of the plate subject to load qr cos θ. Only the simply-supported boundary is considered for brevity hereinafter. It is noted that the real through-thickness deflection distributions of the homogeneous and FGM plate can be well captured by the 3D elasticity solution, which are almost invariant along the thickness. This is just the starting point of the classical plate theories.

Figure 3 draws the distribution of the dimensionless radial normal stress σ_r (a/2, 0, z)/q along the thickness of the plate subject to load qr cos θ. It is observed that the distribution is linear for a homogeneous plate and nonlinear for an FGM plate. Among four values of λ, the radial normal stress in the homogeneous plate is the smallest on the bottom surface (z=-h/2) of the plate and that in the FGM plate with λ=2 is the smallest on the top surface (z=h/2) of the plate.

Figure 4 presents the distribution of the dimensionless shear stress σ_rz (a/2, 0, z)/q along the thickness direction of the plate subject to load qr cos θ. It is seen that the shear stress exhibits a parabolic distribution feature. It is symmetric and its maximum value occurs at the mid-plane for a homogeneous plate and its peak moves toward the lower half of the plate for an FGM plate. It is also observed that the maximum value of the shear stress decreases with λ.

Figure 5 depicts the distribution of the dimensionless deflection (w(r, 0/π, 0)/h)×10⁴ along the radial direction of the plate subject to load qr cos θ. Standard sine function distributions are observed for all four cases, among which there are three zero points and the peak and valley values both occur at r=a/2. Once again, the maximum peak and valley values are obtained for the FGM plate with λ=10 and the minimum ones are for the homogeneous plate.

Due to the components of displacement and stress in the plate in an asymmetric state are all dependent on the circumferential coordinate, which is quite different from an axisymmetric state, the following three figures focus on variations of displacements and stresses with the circumferential coordinate. Figures 6-8 depict respectively the distribution of the dimensionless deflection (w(a/2, θ, h/2)/h)×10⁴, radial normal stress σ_r (a/2, θ, h/2)/q, and shear stress σ_rz (a/2, θ, 0)/q along the circumferential direction of the plate subject to load qr cos θ. It is seen that affected by this kind of load, all field variables vary with the circumferential coordinate in the form of cosine function, with each being zero when θ=π/2 and θ=3π/2. Compared with the plate deflection and radial normal stress, the impact of λ on the distribution of the shear stress is very small.

Figures 9 and 10 show respectively the effect of thickness to radius ratio on the dimensionless deflection w(a/2, 0, z)/h×10⁴ and radial normal stress σ_r (a/2, 0, z)/q of the plate subject to load qr cos θ. Here a=1 m is fixed and h changes with the thickness to radius ratio h/a. It is revealed that the dimensionless deflection decreases with h/a, which is due to the increasing of the bending rigidity of the plate; as for radial normal stress, the maximum tensile and compressive stresses also decrease with h/a and the maximum tensile stress is much greater than the maximum compressive stress. The neutral plane occurring at z=-0.1h is also observed.

6 Conclusions

This paper studied the asymmetric bending of transversely isotropic functionally graded solid circular plates by adopting a generalization of the England's method which is based on the 3D elasticity theory. The material coefficients can vary arbitrarily and continuously along the thickness of the plate. 3D analytical solutions corresponding to every kind of biharmonic loads are thus constructed analytically. The unknown constants are determined from the cylindrical boundary conditions of the plate. The validity and accuracy of the present method are verified by comparing the numerical results with those for the axisymmetric bending of circular plates reported in the literature. It is found from the numerical examples that the gradient index, load pattern, and thickness to radius ratio have important effects on the asymmetric behavior of the FGM circular plate.

It should be emphasized that the present elasticity solutions exactly satisfy the 3D equilibrium equations and the stress boundary conditions on the top and bottom surfaces of the plate. The cylindrical boundary conditions in the plate are adopted in the Saint-Venant sense. Therefore, the proposed 3D analytical solutions can serve as a benchmark for the asymmetric bending solutions of FGM solid circular plates based on various simplified plate theories or numerical methods.

Appendix A Expressions of functions B_k (r), C_k (r), and D_k (r) (k=0, 1, 2, ...) Appendix B Expressions of functions contained in the resultant forces for Case 2