Analytical Solution for Free Vibrations of a Moderately Thick Rectangular Plate

In the present thick plate vibration theory, governing equations of force-displacement relations and equilibrium of forces are reduced to the system of three partial differential equations of motion with total deflection, which consists of bending and shear contribution, and angles of rotation as the basic unknown functions.The system is starting one for the application of any analytical or numerical method. Most of the analytical methods deal with those three equations, some of them with two (total and bending deflection), and recently a solution based on one equation related to total deflection has been proposed. In this paper, a system of three equations is reduced to one equation with bending deflection acting as a potential function. Method of separation of variables is applied and analytical solution of differential equation is obtained in closed form. Any combination of boundary conditions can be considered. However, the exact solution of boundary value problem is achieved for a plate with two opposite simply supported edges, while for mixed boundary conditions, an approximate solution is derived. Numerical results of illustrative examples are compared with those known in the literature, and very good agreement is achieved.


Introduction
Rectangular plate is an important structural element in many modern engineering structures.Spectrum of plate natural frequencies as a continuum is very dense, and its response to an excitation can easily fall in resonance.Therefore, knowledge about plate dynamic behaviour is very important to structural designers.
Theory of thin plate vibrations is very well developed, and the achievements are collected in the fundamental Szilard book [1].Dynamics of thick plate is a more complex problem, due to the effect of shear stiffness and rotary inertia, and is still a subject of investigation.The first works on this subject are those of Reissner [2] and Mindlin [3].They assumed that plate cross-section remains in a straight line but not normal to the plate middle surface due to shear release.Accordingly, constant shear stress distribution through the plate thickness is considered, and concept of shear correction factor is introduced in order to achieve an effective value of shear stiffness.As a result of the introduced assumptions, application of the Reissner-Mindlin theory is limited to the thick plate of moderate thickness and is known as the firstorder shear deformation theory.
Natural vibrations of thick plate have been a challenging problem to many researchers during the past decades.A large number of articles have been published in the meantime and a comprehensive survey of the literature up to 1994 can be found in [4].Generally speaking, there are two main approaches for the solution of thick plate natural vibrations, that is, analytical methods for solution of differential equations of motion and numerical procedures based on the Rayleigh-Ritz energy method and FEM.In the former approach system of the governing equations of internal forces, inertia forces, displacements, and deformations is reduced to three equilibrium equations expressed with plate deflection and angles of rotations.Application of analytical methods is limited to simply supported plate at two opposite edges.Actually, total plate deflection consists of bending deflection and shear relaxation, and angles of rotation are dependent only on the bending.Analytical methods are different depending on which functions are kept as fundamental ones in further reduction of differential equations of motion.Hence, it is possible to operate with three, two, or even one function, that is, total deflection and two angles of rotation, total deflection and bending deflection, and only with total deflection, respectively.
Among a large number of articles, it is worthy to mention the Wang paper [5] in which an explicit formula for natural frequencies of simply supported thick plate in the form of that for thin plate is derived.Brunelle offers a procedure for determination of natural frequencies for a simply supported plate on two opposite edges [6].In the Xiang papers, the well-known Levy approach in thin plate theory for a simply supported plate on two opposite edges is employed to the thick plate vibration analysis [7,8].The semianalytical solution for natural vibration is proposed by Malekzadeh et al. [9].The same problem is considered in the Endo and Kimura article with total deflection and bending deflection as two fundamental functions [10].A refined plate theory with two functions, that is, bending deflection and shear deflection, is proposed by Shimpi and Patel, and it is applied only for simply supported plate at all edges [11].The rigorous analytical solution for case of simply supported edges on two opposite edges is given in the paper of Hashemi and Arsanjani [12].An interesting analytical solution for any combination of simply supported and clamped edges is presented by Xing and Liu [13].It is developed under assumption that angles of rotation in vertical longitudinal plane vary in transverse direction proportionally to the total deflection, instead of actual bending deflection, and vice versa.In spite of that approximation, very good results are achieved.Xing and Liu proposed also a new variant of their plate theory in which three classical differential equations are reformulated to arrive at two equations with total deflection and generalized curvature as the basic functions [14].
Since application of the analytical methods can be applied only to some combinations of boundary conditions, numerical methods are used for more complex problems of mixed boundary conditions, as well as for elastically supported edges.Mostly, the Rayleigh-Ritz method is applied, and accuracy of results depends on the chosen set of orthogonal functions for assumed natural modes.In the paper of Liew et al., two-dimensional polynomials are used for natural modes [15].The static Timoshenko beam functions are very often applied for plate longitudinal and transverse directions, for instance, Dawe and Roufaeil [16] and Cheung and Zhou [17].Recently, vibration problems of elastically restrained edges are in the focus [18][19][20].As it is well known, the finite element method is a very powerful tool for structural analysis, and several finite elements for the Mindlin plate have been developed and incorporated in commercial FEM software.A recent formulation of a new finite element for a thick plate is proposed by Falsone and Settineri [21].
In order to extend applicability of analytical methods for the solution of natural vibrations of a thick plate for any combination of boundary conditions, in the present paper, an advanced moderately thick plate theory is presented.The single differential equation of motion is derived with bending deflection as the basic function.The method of separation of variables is used, and a rigorous solution is obtained.Application of the method is illustrated in the case of different boundary conditions, and results are compared with those known in the relevant literature.

Basic Equations
A thick rectangular plate is considered in the Cartesian coordinate system with corresponding displacements, that is, deflection  and angles of rotation in longitudinal and transverse vertical plane, like in the Timoshenko beam theory [22], Figure 1.Deflection of thick plate can be decomposed into bending deflection and shear deflection, according to the beam theory [23] and plate theory [24]  (, , ) =   (, , ) +   (, , ) .
Angles of rotation of plate cross-sections in the lower frequency domain are predominantly caused by bending deflection: while shear deflection makes contribution only to rotation of plate generatrix.Hence, bending moments and twist moments are functions of plate bending curvature: where is plate flexural rigidity, and ℎ, , and ] are plate thickness, Young's modulus of elasticity, and Poisson's ratio, respectively.The shear forces are results of shear deformations: where  = ℎ is shear rigidity and  is shear coefficient.The vibrating plate is loaded with normal inertia load and distributed bending moments: where  = ℎ is the plate mass per unit area and  =  = (ℎ 3 /12) is the mass moment of inertia of cross-section per unit breadth.

Differential Equations of Motion
In order to derive the differential equation of motion, let us consider the equilibrium of moments and forces: Substituting expressions for forces and loads, (3a), (3b), (3c), (5), and (6), respectively, into (7) and taking (1) into account, after some rearrangement, yield Equations ( 8), (9), and (10) represent the system of three partial differential equations with two unknowns, that is,   and   .However, (8) and ( 9) can be reduced to one equation by integrating the former per  and the latter per : It is obvious that the undetermined integration functions have to be the same, that is,   (, ) =   (, ) = (), where () represents rigid body translation.Furthermore, by substituting (8), (9), and (11) or (12) with () into (10), one differential equation of motion is obtained with the unknown bending deflection: where Once   is determined, the total deflection according to ( 1) and ( 11) yields Referring to (13), bending deflection consists of the homogenous solution and a particular integral due to the disturbing function ().Putting   =   () one finds from (13), Furthermore, by substituting ( 15) into ( 14), the disturbing function () is cancelled.Since it does not appear in the expression for the total deflection (14) and derivatives of the bending deflection, by which angles of rotation and sectional forces are determined, one can get () = 0. Hence, differential equation ( 13) becomes that homogenous one.

Boundary Conditions
The following boundary conditions are considered at the transverse edge at  = .
Simply supported edge: Fixed edge: Free edge: where   is the shearing force comprising ordinary transverse shear force and the effect of torsional moment as in the thin plate theory [1].According to (3c), (5), and ( 8), the shearing force expressed with bending deflection reads Boundary conditions for longitudinal edge at  =  can be specified in an analogous way.

Differential Equation of Natural Vibrations
Natural vibrations are assumed to be harmonic, and therefore, where  is the natural frequency and   is the mode amplitude.In this case, differential equation of motion ( 13) is reduced to the frequency domain.At the same time, the second term in ( 13) can be expressed in a more appropriate way: so that the final form of ( 13) reads Amplitude of the total deflection, according to (14), yields and amplitude of rotation angles according to (2) Equation ( 22) can be reduced to the modified Timoshenko beam equation expressed with bending deflection [23].

Analytical Solution of Natural Vibrations
Differential equation of natural vibrations (22) can be solved by the method of separation of variables [25].In that case, the bending deflection is presented as a product of two functions, each depending only on one variable.By substituting ( 25) into (22), one obtains where Furthermore, according to the procedure presented in [13], each of the single variable functions can be assumed in exponential form By inserting ( 28) into (26), a characteristic equation is obtained It can be reduced to one unknown by substitution of  2 + 2 =  2 .Hence, characteristic equation ( 29) takes the following forms Two roots of (30) read where Hence, there are two solutions for  2 +  2 , that is, and each of the unknown parameters  and  has four solutions, two real and two imaginary: Parameters , , , and  are real and positive quantities.Their mutual dependence is obtained by inserting   and   in (33a) and (33b) in all combinations as follows. 1,2 and  1,2 in (33a),  3,4 and  3,4 in (33b), 1,2 and  3,4 in (33b),  3,4 and  1,2 in (33a), 1,2 and  3,4 in (33a),  3,4 and  1,2 in (33b), Substitution of  1,2 and  1,2 into (33b) and  3,4 and  3,4 into (33a) does not give real results, and therefore, that combination is not acceptable.By the determined parameters , , , and , the single variable functions (25) can be presented in the following form: =  1 sinh  +  2 cosh  +  3 sin  +  4 cos .(40)

Boundary Conditions in terms of Separated Variables
Based on known functions () and (), the total deflection (23) takes the following form; For specification of boundary conditions at the transverse and longitudinal edges, where  =  and  = , respectively, the total deflection can be presented in the following ways: Angles of rotation, (2), have simple expressions: In a similar way, the bending moments, (3a), (3b), and (3c), take the following form: The shearing forces, according to (19), yield In the previous formulae for the  and  directions, ratios   / and   / are functions which make specifying boundary conditions as algebraic equations impossible.Their scalar values are obtained only in the special case of simply supported edges [26].Since one obtains where  and  are the total number of mode half waves in  and  direction, respectively.Equation (47) can be used for arbitrary boundary conditions as the first rough approximation.Furthermore, according to (39) and (40), functions  and  consist of the hyperbolic and trigonometric part, that is, By taking into account (47) and (48), the expressions for the geometric boundary conditions read as follows.Deflection: where Mathematical Problems in Engineering Angles of rotation: In a similar way, the bending moments (44) can be written in the following form: where Finally, for the shearing forces (45), one can write where (57)

Natural Vibrations of Clamped Plate, CCCC
Each of the four rectangular plate edges can be simply supported, clamped, or free.That results in a quite large number of possible combinations of boundary conditions.However, in order to illustrate applicability of the developed theory, it is sufficient to analyse a few typical boundary value problems.
Hence, let us first consider vibrations of a clamped plate.Boundary conditions at transverse edges read According to (39), (49), and (52), (58) lead to the following expressions: From (60) it yields that Substitution of (61) into (59) gives a system of two algebraic equations which can be written in the following matrix form The nontrivial solution is obtained if determinant of ( 62) is zero, after some manipulations, which leads to the following frequency equation: In the homogenous system of algebraic equations (62), one of the integration constants  1 and  2 is arbitrary, while the other is determined by the first or the second equation.From the first equation one can write Constants  3 and  4 are obtained by employing (61).In that, way function , (39), for mode shapes in longitudinal direction is defined.Boundary conditions for the clamped longitudinal edges read Following the previous procedure, the final frequency equation can be written in an analogous way: For integration constants, one can write In that way, function , (40), is defined.Now, the bending deflection mode can be determined   =  and then the complete deflection mode by employing (41).
In frequency equations ( 63) and (66), five quantities are unknown, that is, , , , , and , and therefore, it is necessary to specify three additional conditions.For that purpose, six relations are on disposal, that is, (36a), (36b), (37a), (37b), (38a), and (38b).Frequency  is the basic variable, and it is convenient to choose parameter  of the trigonometric functions of the longitudinal plate generatrix as the second basic variable, since it remains in the case of simply supported edges.The other parameters , , and  can be expressed with .Summing up (36a) and (38b) or (36b) and (38a) yields Furthermore, one finds from (37a) and (36b) that It is interesting, that summation of square of (69a) and (69b) gives an expression similar to (68), that is, The previous relations of vibration parameters are geometrically dependent based on Pythagoras' theorem and therefore can be constructed as shown in Figure 2. Assuming the value of ,  1 and  2 are determined and presented as catheti of the right triangle ABC with hypotenuse .The circumscribed circle of rectangle ABCD is drawn.Furthermore, value of parameter  is assumed, and  is defined as cathetus of triangle AEC, according to (68).Value of  = CE is rotated so that  is obtained as the hypotenuse of triangle CDG according to (69a).Parameter  is determined as the cathetus of triangle BCF, according to (69b).Finally,  = DG and  = BF are rotated and their intersection, point , lies at the circle according to (70).In that way, the plan of the vibration parameters is closed.It changes depending on boundary conditions and vibration mode.Now, , , and  from (68) and (69a) and (69b) can be substituted into frequency equations ( 63) and (66).In the previous way, system of five nonlinear algebraic equations is reduced to two equations with unknown parameters  and .It can be solved by a Newton-Raphson iteration procedure.For that purpose, it is necessary to specify the numbers of mode half waves  and  in  and  directions, respectively.Hence, ,  is actually mode identification number.In case of clamped plate,  and  take values 2, 3, 4,. ... Furthermore, an expected value of  has to be prescribed as the initial value in the iteration procedure.Also, a realistic value for  has to be assumed.A set of few mathematical solutions for  and  is usually obtained, but only that solution which manifests the same numbers of the mode half waves  and  as the prescribed ones is physically correct.
Accuracy of the natural frequencies can be increased if calculation is repeated with ratios (47) of the condensed functions.For that purpose, the average values based on the 1-norm can be applied: where m and ñ are numbers analogous to  and .The previous integrals can be solved numerically.Calculation of natural frequencies is performed for rectangular plate of aspect ratio / = 0.4, thickness ratio ℎ/ = 0.001 and 0.1, and shear coefficient  = 5/6.The obtained values of frequency parameter Λ = ( 2 / 2 )√ℎ/ are listed in Tables 1 and 2 for thin and thick plates, respectively.PS 0 is the present solution obtained with given  and , while PS 1 and PS 2 are the first and the second iterative solution determined with average values of the functions of separated variables, (71), respectively.In that procedure,  and  change values to floating point values m and ñ.
Results for the thin plate, Table 1, are obtained directly without iteration.They are compared with Xing's solution [13]   and Liew's solution [15].The former is determined analytically in closed form, and for the latter the Rayleigh-Ritz method is employed.PS 0 and Xing's values are the same, while Liew's values are somewhat higher.
For thick plate, Table 2, stable values of vibration parameter Λ are obtained in two iteration steps.They are very close to the Xing's values, while Liew's values are somewhat higher.Values of  and  are moderately reduced in the first step of iteration and very little in the second one.Even the first approximate values of natural frequencies, PS 0 , are very close to the Liew values.Generally speaking, if natural frequency does not converge to the expected value, it is necessary to adjust values of  and .
For illustration, the first natural mode of the thick plate is presented by deflection of longitudinal generatrix in the middle of the plate, Figure 3. Similar diagrams can be drawn for transverse generatrix.Total deflection (, /2) consists of bending deflection   (, /2) and shear deflection   (, /2).Their boundary values exist as a result of ignoring integration function (), Section 3.However, the boundary values are of the close magnitudes but opposite signs, and therefore, they are almost cancelled.As a result, the boundary conditions are approximately satisfied, that is, (0, ) = (, ) = (), where the shape of discrepancy function () is similar to that of shear deflection   (  , ).The same consideration is valid for the longitudinal edges.Actually, the boundary conditions for deflection are satisfied in the mean according to the introduced 1-norms (71).Furthermore, it is obvious from the boundary slope of the bending deflection that the angles of rotation are zero as specified by the boundary conditions.The shape of the shear deflection is similar to the bending moment of clamped beam, since it depends on curvatures.Indeed, the main shear part in the formula for total deflection, (23), can be expressed with the bending moments, (3a), (3b), and (3c): Shear contribution to total deflection is increased for higher modes due to shortening of wave length.The shape of the first natural mode is shown in Figure 4, where negligible small hold-down of the plate corners can be noticed due to the previous explained reasons.(74) The previous system of algebraic equations leads to the following frequency equation: The frequency equation for transverse direction can be written in the analogous way since the boundary conditions are the same in both cases: The integration constants   and   ,  = 1, 2, 3, 4, are expressed by the same formulas as in the previous cases, (61), (64), and (67).Since one of two opposite edges is simply supported and another is clamped, numbers of mode half waves take values  = 1.5, 2.5, . . .and  = 1.5, 2.5, . ... The numerical procedure for determining natural frequencies is the same as in the case of clamped plate.
Values of vibration parameter Λ for thin and thick plates are presented in Tables 3 and 4, respectively.In the former cases results obtained by all three methods are the same.For thick plates the initial values of Λ are very close to the final values.Also, PS 2 values are very close to Xing's values.However, Liew's values are somewhat higher as in the previous example.The second natural mode is shown in Figure 5, where only slight release of the corner of clamped two adjacent edges is registered.
(77) The expanding form of (77) is the following: (78) The previous system of algebraic equations leads to the following frequency equation: The integration constants are the following: Based on the expected mode shape, the following values of modal parameters are taken into account:  = 1.5, 2.5, . . .and  = 0.5, 1.5, . ... The obtained results for thin and thick plates are listed in Tables 5 and 6, respectively.Values of the frequency parameter Λ for thin plate are determined by iteration since in this case free edge convergence is not as fast as for the other boundary conditions.Very good agreement with Liew's values can be noticed.Calculation of natural vibration is also performed by NASTRAN [27] with mesh density of the finite element model 8×16.That is sufficient to determine the first 8 natural frequencies accurately.The obtained results are quite close to those of PS 2 and Liew's values.
The fifth and sixth natural modes determined by the present theory are shown in Figures 6 and 7, where the mode identification number is ,  = 1.5, 4.5 and 2.5, 0.5, respectively.Intention is to point out the transition of natural vibrations from a dominantly large number of mode half waves in transverse direction to longitudinal direction.A slight release of the corner of clamped and free edge can be noticed.The same natural modes obtained by FEM analysis are presented in Figures 8 and 9 [27].There is no other visual difference between the mode shapes determined in these two ways.They can be distinguished mainly by the smooth and straight lines of the deformed mesh.

Conclusion
An advanced vibration theory of moderately thick plate vibrations is presented.The three traditional partial differential equations of motion with total deflection and angles of rotation are reduced to one equation with bending deflection as the basic unknown acting as a potential function.Its general solution is determined exactly by the method of separation of variables.The frequency equations are exact in the case of simply supported plate and simply supported two opposite edges.These two special problems are elaborated in detail in [26].An approximate but very reliable solution is obtained for arbitrary boundary conditions which are satisfied in the mean, except angles of rotations that are satisfied exactly in any case.For the first time, an analytical solution is obtained for the case of a free edge where Levy's approach is not applicable.
As pointed out, the new theory can be used for vibration analysis of moderately thick plate.Very thick plates, where coupling of bending shear and in-plane shear modes appears at higher natural frequencies, can be analysed by the method described in [28].Actually, very thick plates are elastic bodies, and 3D analysis is more effective.By the proposed method, results obtained by numerical methods, as the Rayleigh-Ritz method, can be verified.Also, finite elements for thick plate can be tested.The differential equation of motion is derived in the time domain and can be applied for determining both harmonic and transient force responses.
Difference between the proposed moderately thick plate theory and the known theories is obvious as can be seen from the outline of the latter presented in the Appendix.