Buckling of circular plates with functional grading in two directions

This short note considers thin circular plates, functionally graded in both axial and transverse directions and loaded in compression on the middle plane by a uniform axisymmetric load. The functional grading is based on recent literature on the subject and we deal with a direct problem for buckling, i.e., given the geometry of the plate and its constitutive properties, the critical load multiplier and the buckling mode are determined by a usual non-triviality condition. Original expressions are found in a nearly closed form and show that suitable functional grading, actually achievable in practice, may lead to gentle solutions. Numerical examples, physical interpretations and comments are provided.


Introduction
Functionally graded materials (FGM) have attracted interest in the scientific community in the last decades due to their promising possible applications in several branches of engineering and applied sciences. Among the literature we may quote [1][2][3][4], where actual realizations of functional grading are widely illustrated, with details on the technological processes involved. Thus, many investigations have been performed also on the field of structural mechanics. Without pretending to be exhaustive, and only to give a hint, we may quote the monographs [1,5,6], where a lot of problems are solved with a view towards structural applications.
One benchmark study in structural mechanics is on the possible buckling of an element, since the pioneering works by Euler on the bifurcated solutions of a compressed elastic column. About buckling of FGM structural elements, we may quote the papers [7], dealing with the so-called Timoshenko beams, and [8], dealing with circular plates. There actually is a vast literature on this subject, considering various beam models (purely flexible, fully deformable, sheardeformable at higher orders with respect to the standard theories), various plate models (according to Kirchhoff and Love, Reissner and Mindlin, Von Kármán), plus interactions with elastic foundations, porosity and other physical aspects. However, this is out of our scope, since in this short note we wish to Indeed, it is not necessary to remark how closedform solutions can be of interest in engineering; one may quote the monograph [9] and the paper [10], devoted to the search of such solutions for problems of technical interest. In a recent paper [11], one of the authors investigated the possibility of finding closedform solutions for the buckling of FGM circular thin plates under uniform radial compression. The present short technical note investigates the possibility to extend the study in [11] to a more general functional grading, following the suggestions of literature, for instance in [11][12][13][14]. The novel point of this note with respect to the quoted literature is the finding of results in a nearly closed form, which makes them quite suitable for possible applications. Since this is only a short note, paving the way for a more vast and profound investigation, we do not pretend to be exhaustive. We are satisfied to show that by recent technological suggestions we are able to find closedform solutions for buckling loads; we will not then investigate here about possible imperfections, postbuckling behaviour, sensitivity to initial imperfections, limit loads for actual equilibria.

Formulation of the problem
Let us consider a circular plate of radius R and uniform thickness h ( R. Under these hypotheses, the plate is considered thin and can be described by the so-called Kirchhoff-Love model: material fibres orthogonal to the middle plane in the reference shape remain orthogonal to the deflected middle surface.
Let us fix a cylindrical coordinate system with origin at the centre of the middle plane of the plate and let q, /, z be the radial, angular and transverse coordinates of a point of the plate, respectively. The middle plane of the plate is characterized by z = 0; the kinematic and resultant action fields are defined on it.
In order to model FGM, Young's modulus E is assumed to vary both along the radius and the thickness of the plate, i.e., E ¼ E q; z ð Þ. By contrast, Poisson's ratio v is assumed uniform over the plate, which means that the shear modulus G is made to vary accordingly; however, the law of functional grading for G will not enter the field equations and we will henceforth not give any further detail on this.
Let w be the deflection of the middle plane, describing a surface of revolution for axial symmetry; let h be the angle between the axis of revolution of the deflected surface and any normal to it. For axial symmetry, these fields depend on the radial variable q only. According to the Kirchhoff-Love model, h coincides with the rotation of a material fibre originally orthogonal to the middle plane in a meridian plane, i.e. ( [14][15][16][17]) Let v q , v / be the local radial and hoop curvatures of the middle plane, respectively. These are representative of the behaviour of the entire plate, since shear deformations are neglected, and are given by ( [14][15][16][17] where, for the sake of simplicity, the dependence of the indicated fields on the radial coordinate q is omitted. Let the material of the plate be linear elastic, with chemical composition varying smoothly across the plate. The elongation of the radial fibres e q will depend linearly on z through the local curvature v q ; analogously, that of the hoop fibres e / will depend linearly on z through the local curvature v / ( [14][15][16][17]): By Hooke's law, the normal stresses along the radial and hoop directions, r q and r / , are ([14-17]): By integrating the product of the expressions in Eq. (4) for r q and r / by the lever arm z over the thickness, we obtain the radial and hoop bending moments M q ; M / per unit of length ( [14][15][16][17]): Remark that no shearing stress, shearing strain, and mixed resultant moment exist, due to axial symmetry ( [14][15][16][17]). Moreover, in Eq. (5) the stress resultants depend on both radial and transverse coordinates due to the functional grading, which, inspired by some literature ( [8,[11][12][13][14]), is assumed to be described by: The uniform value m 0 of Poisson's ratio is admissible for metals and their alloys, and E 1 is uniform across the plate. The real parameters n; m are form factors of the distribution of Young's modulus and can be fixed, e.g., using conditions for the local composition of the plate at special places.
The functional grading in Eq. (6), suitable for our scopes of easier calculation and the seek for closedform expressions, has solid bases in the literature, as remarked above, and is application-oriented. Indeed, it foresees a stiffer, heavier material towards the upper and lower faces of the plate, to increase its global strength, and a softer, lighter material moving from the centre to the outer edge, to keep global strength while sparing weight.
Replacing the distributions of Eq. (6) into Eq. (5) yields where D q ð Þ denotes the local flexural rigidity, expressed as while in a homogeneous plate it is uniformly given by ( [14][15][16][17]) Thus, by the functional grading (6) the variation of Young's modulus along the radius and the thickness makes the flexural rigidity depend only on q, modulo the amplification factor a in Eq. . This is due to the symmetric distribution of Young's modulus with respect to z and has interesting consequences.
Indeed, symmetric mechanical properties of the plate with respect to its middle plane imply that the latter necessarily coincides with the neutral surface. Further, this makes all the standard equations for plates that we find in ( [14][15][16][17]) valid, since all the mechanical fields of interest may be reduced to the plane z ¼ 0: consequently, bulk equations and boundary conditions are well defined on the middle plane and on its outer circumference. Hence, the usual considerations for possible buckling hold, and no special theories, nor adjusted boundary conditions, are required: we will then proceed with usual means.
Inserting Eq. (2) into Eq. (7), and denoting qderivatives by primes, yields The local balance of moment along the radial direction is ( [14][15][16][17]): where QðqÞ is the transverse shearing force per unit length, which is constitutively unprescribed and can thus be determined only via the balance of transverse force.
Let the external compressive load N q be applied uniformly along the outer edge of the plate ðq ¼ RÞ, with no eccentricity with respect to the middle plane; due to the assumed symmetry of material properties, this is reasonable and represents the ''perfectly loaded'' structural element, similar to the centred compression of a Bernoulli-Euler beam. Moreover, we suppose that the plate remains plane until buckling, and there is no bending or deflection during the prebuckling phase (that is, the fundamental path is trivial): this is also a standard assumption, and is reasonable since, as already remarked, by the functional grading (6) the global strength of the plate is maintained while sparing weight.
In any possible adjacent configuration in which the presence of a deflection of the middle surface is admitted, the balance of transverse force yields the value of the shearing force ( [10][11][12][13][14]) Substituting Eq. (8) and (10) into Eq. (11) and (12) yields the novel field equation for the buckled shape In order to abstract from the actual values of the geometrical and physical parameters involved, let us perform a non-dimensional analysis. If we then let q ¼ rR, Eq. (13) assumes the dimensionless form Equation (14) is a second-order ordinary differential equation with variable coefficients, which must be completed by two boundary conditions to look for solutions. The first derives from axial symmetry: the transverse fibre at the centre of the plate must remain vertical, i.e., h r ¼ 0 ð Þ¼0, which results in where 2 F 1 c; d; e; z ½ is Gauss' hypergeometric function ( [18]) and C is an integration constant, since just one boundary condition was used; C must be imaginary in order for the rotation to assume a real value. By imposing the second boundary condition at the outer edge of the plate q ¼ R ) r ¼ 1 ð Þ , the latter may attain a configuration different from the fundamental, trivial, one, provided by C ¼ 0. This nontriviality condition, expressed in terms of the nondimensional compressive load ratio p, provides the value of the buckling load. As usual, the buckling mode may be found, but C remains undetermined [19].
Equation (14) may be turned in terms of the transverse displacement w r ð Þ, yielding Equation (16) is a third-order ordinary differential equation, and maintains the same expression in nondimensional form; again, the transverse fibre at the centre of the plate cannot rotate, i.e., w 0 r ¼ 0 ð Þ¼0, which results in In Eq. (17) j F k ½c j ; d k ; z is the generalized hypergeometric function ( [20]) and A; B are integration constants, in general complex. By imposing the two remaining boundary conditions at the outer edge of the plate, and using the non-triviality condition, the buckling load and mode can be found.

Clamped plate
If the outer edge of the plate is clamped, we use Eq. (14) and (15); the relevant boundary condition is Once again, remark that, due to the symmetry of material properties with respect to the middle plane, the usual expressions for the boundary conditions may be provided -a clamp prevents the corresponding rotation of the transverse material fibres, whence Eq. (18).
Equation (18) admits non-trivial solutions for C, hence buckled shapes for the plate, if the numerator fðm; pÞ of its left-hand side vanishes. This condition can be graphically visualized by intersecting the twodimensional surface represented by the numerator of the left-hand side of Eq. (18) depending on the parameters m and p for a given value m ¼ 1=3; with the plane f ¼ 0, as shown in Fig. 1; the curve common to the two surfaces collects all the solutions of Eq. (18).
To provide an approximate (yet meaningful) analytical formulation of the critical values p cr of p as a function of the shape factor m, at least in the range 1:1 m 5, we tabulate the values p cr corresponding to equally spaced nodes, then approximate the data with a second-order polynomial by the least square method, yielding 4 Simply supported plate If the outer edge of the plate is simply supported, the additional boundary conditions for Eq. (16) and (17) are where the usual comments on the standard expressions for the conditions provided by the support.
The first Eq. (21) provides a relationship between A and B The second Eq. (21) accounts for the constraint, ineffective towards rotation about /, thus exerting null bending couple in the same direction; due to Eq. (10), this condition equals Inserting Eq. (17) and (22) into Eq. (23), in order to have non-trivial solutions one must satisfy yielding the critical multipliers of the applied load. As before, the solution of Eq. (24) is shown graphically, intersecting the surface given by its left-hand side (evaluated for m ¼ 1=3) with the plane f ¼ 0, see Fig. 2.
Proceeding as in the case of the clamped plate, the least square method leads to p cr m ð Þ ¼ À0:0988 m 2 þ 0:8311 m þ 2:2313 ð25Þ

An example
In order to provide a numerical example, a composite Aluminium-Alumina plate is studied; these materials have Young's Modulus 70 GPa and 380 GPa, respectively. It is further assumed h ¼ 1mm and R ¼ 25mm.
The following functional material grading conditions are considered: These conditions are verified only if the parameters n, m and E 1 assume the values Once these values, plus m ¼ 1=3, are fixed, the buckling load for the clamped plate is found by Eq. (15), (20) and that for the simply supported plate is found by Eqs. (15), (25) In order to compare these results with those provided by literature for a homogeneous and isotropic material, the parameter n must be near zero, while m must be indefinitely large. Let us then pose n ! 0 in Eq. (6), so that the material is uniform along the thickness of the plate and has Young's modulus E 0 q ð Þ that varies only along the radius. However, since the least-square approximation for the buckling load was obtained for 1:1 m 5, it is improper to evaluate the buckling load for m ! 1, in order to have an actual uniform plate. Thus, the buckling load is compared for m providing the vertex of the parabolas in Eqs. (20) and (25) (thus, Young's modulus E 0 q ð Þ has a small variation with respect to E 1 ). In the branch preceding this point, the polynomial function well fits the data: after that point the parabolas in Eqs. (20) and (25)

Final remarks
We investigated the possibility to find closed-form expressions for the buckling loads of uniformly compressed circular plates, the material of which is functionally graded following a quadratic law along both the radius and thickness, inspired by some recent literature; metals and some alloys may comply with it.
The functional grading aims at grossly maintaining the strength of the plate while sparing weight; the symmetric distribution of the material properties with respect to the middle plane of the plate makes the standard bulk equations and boundary conditions valid, yet it highlights the differences with respect to the homogeneous case. Indeed, the bending stiffness varies along the radius of the middle plane and the bulk equations are without constant coefficients. The non-triviality conditions providing buckling are expressed in closed form, which is not usual in the literature. The dependence of the critical value of the load multiplier on a certain range of the characteristic parameter of the functional grading is also found in closed form by a least square method, and this is also not common in the literature. A numerical example is provided; a comparison with known results of the literature is provided, together with a mechanical interpretation. This short note paves the way for the search of other closed-form solutions of technical interest.
Funding Open Access funding provided by Università degli Studi di Roma La Sapienza. This study was partially funded by the institutional Grant Number RM11916B7ECCFCBF of the university ''La Sapienza'' of Roma, Italy.

Compliance with ethical standards
Conflict of interest The authors declare that they have no conflict of interest.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.