# A Series Solution for the Vibration of Mindlin Rectangular Plates with Elastic Point Supports around the Edges

Authors: Fuzhen Pang; Haichao Li; Yuan Du; Shuo Li; Hailong Chen; Ning Liu
Journal: Shock and Vibration
 (2018)

---

## Abstract

A series solution for the transverse vibration of Mindlin rectangular plates with elastic point supports around the edges is studied. The series solution for the problem is obtained using improved Fourier series method, in which the vibration displacements and the cross-sectional rotations of the midplane are represented by a double Fourier cosine series and four supplementary functions. The supplementary functions are expressed as the combination of trigonometric functions and a single cosine series expansion and are introduced to remove the potential discontinuities associated with the original admissible functions along the edges when they are viewed as periodic functions defined over the entirex-y plane. This series solution is approximately accurate in the sense that it explicitly satisfies, to any specified accuracy, both the governing equations and the boundary conditions. The convergence, accuracy, stability, and efficiency of the proposed method have been examined through a series of numerical examples. Some numerical examples about the nondimensional frequency and mode shapes of Mindlin rectangular plates with different point-supported edge conditions are given.


---

## Body

## 1. Introduction

Rectangular plate is of great importance in various engineering branches, such as aerospace, electronics, and mechanical, nuclear, and marine engineering. A better understanding of its dynamic characteristics is meaningful when dealing with the design of a plate-type structure. In the early stage of research about rectangular plate, the majority of research was about thin plate based on classical Kirchhoff hypothesis. However, for the shortage of consideration about transverse shear deformation and rotary inertia, vibration frequencies of thick plates are overestimated. More and more investigations about Mindlin plates began to attract attention after the Mindlin first-order plate theory was proposed by Mindlin [1]. The vibration characteristics of Mindlin plates have been well investigated by researches with classical and elastic edge support [2–5]. Few literatures focused on the vibration behavior of Mindlin plates with elastic point supports around the edges. However, the boundary conditions of the plates are not always classical and elastic around the edges in practical engineering applications. And the boundary conditions of elastic point supports around the edges do exist. So it is much of great significance to study the vibration behavior of Mindlin plates with elastic point supports around the edges.In the field of vibration analysis of structures subject to elastic edge/point supports, Takashi and Jin [6] used the collocation method to investigate Mindlin plates with constant thickness and two opposite edges simply supported. Jiarang and Jianqiao [7] established the three-dimensional state equation for the laminated thick orthotropic plate with simply supported edges and obtained the numerical results. Ohya et al. [3] investigated the free vibration characteristics of rectangular Mindlin plates which have simultaneous elastic edge and internal supports via the superposition method. Hosseini-Hashemi et al. [8] proposed an exact closed-form procedure to solve free vibration of moderately thick rectangular plates with two opposite edges simply supported. Dozio et al. [9] developed a Ritz method using a set of trigonometric functions to obtain accurate modal properties of rectangular plates with arbitrary thickness. Maretic [10] analyzed the transverse vibration of a circular plate loaded by uniform pressure along its edge. Based on the shear deformable plate theory, Bahmyari [11] used the Element-Free Galerkin Method to analyze the free vibration of inhomogeneous moderately thick plates with point supports resting on a two-parameter elastic foundation. The results show that the vibration results obtained are in a very good agreement with the available literatures in spite of using low numbers of nodes. Bashmal et al. [12] investigated the in-plane free vibrations of an elastic and isotropic annular disk with elastic constraints at the inner and outer boundaries. And the inner and outer boundaries are applied either along the entire periphery of the disk or at a point. Esendemir [13] studied polymer-matrix composite beam of arbitrary orientation supported from two ends acted upon with a force at the midpoint by an analytical elastoplastic stress analysis. Foyouzat et al. [14] presented an analytical-numerical approach to determine the dynamic response of thin plates resting on multiple elastic point supports with time-varying stiffness. Gan et al. [15] studied the effect of intermediate elastic support on the vibration of functionally graded Euler-Bernoulli beams excited by a moving point load. Kucuk et al. [16] carried out the analytical elastic-plastic stress analysis on an aluminum metal-matrix composite beam reinforced by unidirectional steel fibers supported at the ends acted upon with a force at the midpoint. Hosseini-Hashemi et al. [17–19] used the generalized differential quadrature method to study the buckling analysis and dynamic transverse vibration characteristics. Setoodeh and Karami [20] employed a three-dimensional elasticity based layerwise finite element method (FEM) to study the static, free vibration, and buckling responses of general laminated thick composite plates. In this paper, the elastic line and point supports are successfully incorporated for thick plates. On the basis of three-dimensional elasticity theory, Xu and Zhou [21] studied the bending of simple-supported rectangular plate on point supports, line supports, and elastic foundation. Wu [22] analyzed the free vibration of arbitrary quadrilateral thick plates with internal columns and elastic edge supports by using the powerful pb-2 Ritz method and Reddy’s third-order shear deformation plate theory. Wu and Lu [23] studied the free vibration of rectangular plates with internal columns and elastic edge supports using the powerful pb-2 Ritz method. Other literatures related to this field are shown in [24–28].From the literatures review, we can know that most researches focused on the classical plates. Only few literatures studied the vibration characteristics of Mindlin plates with point supports around the edges. And these literatures mainly used the collocation method, analytical-numerical method, and so on. To the best of the authors’ knowledge, there are no published literatures focused on vibration characteristics of Mindlin plates with point supports around the edges by the method of modified Fourier. Thus, a unified, efficient, and accurate formulation to deal with the free vibration of Mindlin plates subjected to arbitrary point-supported boundary condition is necessary and much of great significance.In this paper, a modified Fourier solution for the free vibration of Mindlin rectangular plates with elastic point supports around the edges is proposed. The vibration displacements and the cross-sectional rotations of the midplane are represented by a double Fourier cosine series and four supplementary functions. The supplementary functions are in the form of the combination between trigonometric functions and a single cosine series expansion. By importing supplementary functions, the discontinuity of the boundary condition is overcome. The change of the boundary conditions can be easily achieved by only varying the stiffness of the boundary springs around all edges of the plates without involving any change to the solution procedure. The natural frequencies of Mindlin rectangular plates are obtained by using the Rayleigh–Ritz method. Effectiveness and stability of this method will be verified by comparing with the results of FE modeling. Some numerical examples of free vibration for Mindlin rectangular plates with different aspect ratio and thickness are conducted under different point-supported conditions.


## 2. Theoretical Formulations

### 2.1. Point-Supported Edge Conditions

The rectangular Mindlin plate with arbitrary elastic point edge supports is shown in Figure1. The boundary conditions are presented by three kinds of restraining springs [29–34], namely, translational, rotational, and torsional springs. Springs are evenly arranged on each edge of Mindlin plate. Through changing the stiffness of springs, different boundary conditions can be achieved [35, 36]. The governing differential equations of a Mindlin plate about its free vibration are as follows:(1)kGh∂2w∂x2+∂2w∂y2+∂ψx∂x+∂ψy∂y-ρh∂2w∂t2=0(2)D∂2ψx∂x2+1-μ2∂2ψx∂y2+1+μ2∂2ψy∂x∂y-κGh∂w∂x+ψx-ρh∂2ψx∂t2=0(3)D∂2ψy∂y2+1-μ2∂2ψy∂x2+1+μ2∂2ψx∂x∂y-κGh∂w∂y+ψy-ρh∂2ψy∂t2=0Figure 1
A Mindlin plate with arbitrary elastic point edge supports.In the formulas above, transverse displacement is presented byw, ψx is the slope along x direction, and ψy is the slope along y direction. ρ means the mass density, μ is Poisson’s ratio, κ is the shear correction factor, and the thickness of the plate is h. G=E/(2(1+μ)) is the shear modulus, μ is Poisson’s ratio, and D=Eh3/(12(1-μ2)) is the flexural rigidity.Besides, the bending moment is expressed by (4) and (5), twisting moment is calculated by (6), and the transverse shearing forces in plates are expressed by (7) and (8).(4)Mxx=D∂ψx∂x+μ∂ψy∂y(5)Myy=D∂ψy∂y+μ∂ψx∂x(6)Mxy=D1-μ2∂ψx∂y+∂ψy∂x(7)Qx=κGhψx+∂w∂x(8)Qy=κGhψy+∂w∂yThere are three kinds of forces along every edge, namely, the bending moment, the twisting moment, and the shearing forces. Rotational, torsional, and translational springs along every edge are the counterparts to these three forces in this paper. The boundary conditions for an elastically restrained rectangular plate are as follows.Atx=0,(9)kx0yw=-Qx(10)Kx0yψx=-Mxx(11)Kyx0yψy=-Mxyk x 0 ( y ), Kx0(y), and Kyx0(y) represent the stiffness of the translation spring, rotational spring, and torsional spring, respectively, on the boundary edge x=0; the corresponding expression results in (12)kx0y=∑i=1Nkx0iδy-yx0i(13)Kx0y=∑i=1NKx0iδy-yx0i(14)Kyx0y=∑i=1NKyx0iδy-yx0iwhere N is the elastic point supports numbers on edge x=0, δ(x) is the Dirac delta function, and kx0i, Kx0i, and Kyx0i are, respectively, the stiffness of the translational, rotational, and torsional springs at the ith support located at yx0i. Similarly, the remaining boundary conditions on the other three edges can be expressed as follows.Atx=a,(15)kxayw=Qx(16)Kxayψx=Mxx(17)Kyxayψy=Mxy(18)kxay=∑i=1Nkxaiδy-yxai(19)Kxay=∑i=1NKxaiδy-yxai(20)Kyxay=∑i=1NKyxaiδy-yxaiAt y=0,(21)ky0xw=-Qy(22)Ky0xψy=-Myy(23)Kxy0xψx=-Mxy(24)ky0x=∑i=1Nky0iδx-xy0i(25)Ky0x=∑i=1NKy0iδy-yy0i(26)Kxy0x=∑i=1NKxy0iδx-xy0iAt y=b,(27)kybxw=Qy(28)Kybxψy=Myy(29)Kxybxψx=Mxy(30)kybx=∑i=1Nkybiδx-xybi(31)Kybx=∑i=1NKybiδx-xybi(32)Kyx0x=∑i=1NKyx0iδy-yybiEquations (9)–(32) represent a set of general boundary conditions by setting the spring stiffness to appropriate values. Based on (4)–(8), the boundary condition can be finally written as follows.Atx=0,(33)∑i=1Nkx0iδy-yx0iw=-κGhψx+∂w∂x(34)∑i=1NKx0iδy-yx0iψx=-D∂ψx∂x+μ∂ψy∂y(35)∑i=1NKyx0iδy-yx0iψy=-D1-μ2∂ψx∂y+∂ψy∂xAt x=a,(36)∑i=1Nkxaiδy-yxaiw=κGhψx+∂w∂x(37)∑i=1NKxaiδy-yxaiψx=D∂ψx∂x+μ∂ψy∂y(38)∑i=1NKyxaiδy-yxaiψy=D1-μ2∂ψx∂y+∂ψy∂xAt y=0,(39)∑i=1Nky0iδx-xy0iw=-κGhψy+∂w∂y(40)∑i=1NKy0iδy-yy0iψy=-D∂ψy∂y+μ∂ψx∂x(41)∑i=1NKxy0iδx-xy0iψx=-D1-μ2∂ψx∂y+∂ψy∂xAt y=b,(42)∑i=1Nkybiδx-xybiw=κGhψy+∂w∂y(43)∑i=1NKybiδx-xybiψy=D∂ψy∂y+μ∂ψx∂x(44)∑i=1NKyx0iδy-yybiψx=D1-μ2∂ψx∂y+∂ψy∂xAccording to the Mindlin plate theory, the transverse displacement of the plate middle surface and the rotations of the cross section, respectively, along thex direction and the y direction are utilized. Based on the traditional solution method, the admissible functions usually express a Fourier series expansion, because the Fourier functions constitute a complete set and exhibit an excellent numerical stability in the previous study [37–43]. We found the conventional Fourier series expression to have some defects which contain the convergence problem along the boundary edges except for a few simple boundary conditions, and the derivatives of a Fourier series cannot be obtained simply through term-by-term differentiation.In this study, in order to overcome the shortcoming with the conventional Fourier series expression, the admissible functions will be expressed as a more wholesome form of Fourier series expansion:(45)wx,y=∑m=0∞∑n=0∞Amncos⁡λamxcos⁡λbny+∑l=12ζbly∑m=0∞amlcos⁡λamx+∑l=12ζalx∑n=0∞bnlcos⁡λbny(46)ψxx,y=∑m=0∞∑n=0∞Bmncos⁡λamxcos⁡λbny+∑l=12ζbly∑m=0∞cmlcos⁡λamx+∑l=12ζalx∑n=0∞dnlcos⁡λbny(47)ψyx,y=∑m=0∞∑n=0∞Cmncos⁡λamxcos⁡λbny+∑l=12ζbly∑m=0∞emlcos⁡λamx+∑l=12ζalx∑n=0∞fnlcos⁡λbnywhere Amn, aml, bnl, Bmn, cml, dnl, Cmn, eml, and fnl represent the unknown coefficients, λam=mπ/a, λbn=nπ/b. The specific expressions of the auxiliary functions ζal and ζbl are defined as(48)ζa1x=a2πsin⁡πx2a+a2πsin⁡3πx2a(49)ζa2x=-a2πcos⁡πx2a+a2πcos⁡3πx2a(50)ζb1y=b2πsin⁡πy2b+b2πsin⁡3πy2b(51)ζb2y=-b2πcos⁡πy2b+b2πcos⁡3πy2bAs shown in (45)–(47), the supplementary functions ζa1, ζa2, ζb1, and ζb2 are used for the x and y direction displacement expressions. It is easy to find that(52)ζa10=ζa1a=ζa1′a=0,ζa1′0=1(53)ζa20=ζa2a=ζa2′0=0,ζa2′a=1Similar conditions are present in the supplementary function iny-direction, namely, ζb1y and ζb2(y). Whereas these conditions are not necessary, the existence of these conditions is purely for the sake of simplifying follow-up mathematical expressions and derivation process.Examining the admissible functions in (45)–(47), one will notice that, besides the standard 2D Fourier series defined over the domain R2:0,a⊗0,b, four additional 1D Fourier expansions are also included in the admissible functions expressions. In light of (52), it is not difficult to see that the fourth term (or the third single Fourier series) on the right side of (45) is equal to the normal derivative of the displacement function w(x,y) at edge x=0. Thus, this term will actually inherit any potential discontinuity associated with the normal derivative at the boundary x=0. Similarly, the other three terms are used to take care of the possible discontinuities at the remaining edges. Thus, the 2D Fourier expansion will now represent a periodic Π2 (residual displacement) function defined over the entire x-y plane. The convergence characteristic for the Fourier series representation of such a function has been well established in mathematics. Namely, the Fourier series will converge with a speed of mπ3 at least.By substituting (45)–(47) into boundary condition (33) and expanding all the y-related terms, except for those containing kx0(y), into Fourier cosine series, one will have(54)∑i=1Nk-x0iδy-yx0i∑m=0∞∑n=0∞Amncos⁡λbny+ζb1y∑m=0∞am1+ζb2y∑m=0∞am2=-∑m=0∞∑n=0∞Bmncos⁡λbny+∑n=0∞β1ncos⁡λbny∑m=0∞cm1+∑n=0∞β1ncos⁡λbny∑m=0∞cm2+∑n=0∞bn1cos⁡λbnyThe expressions for Fourier cosine expansion coefficient of related functions can be found in AppendixA. To establish the relationship between the Fourier coefficients in (54), the terms on the left side shall also be expressed in the form of cosine series, resulting in(55)∑m=0∞∑t=0∞Amt∑i=1Nk-x0icos⁡λbtyx0icos⁡λbnyx0i+∑m=0∞am1∑i=1Nk-x0iζb1yx0icos⁡λbnyx0i+∑m=0∞am2∑i=1Nk-x0iζb2yx0icos⁡λbnyx0i=-δ-bn∑m=0∞Bmn-δ-bnβ1n∑m=0∞cm1-δ-bnβ2n∑m=0∞cm2-δ-bnbn1(56)δ-bn=b21+δ0n(57)δ-am=a21+δ0mwhere δ0n and δ0m are Kronecker delta. Using similar procedures, the other seven boundary condition equations can be obtained. Thus, a total of twelve constraint equations can be derived as follows.Atx=0,(58)∑m=0∞∑t=0∞Amt∑i=1Nk-x0icos⁡λbtyx0icos⁡λbnyx0i+∑m=0∞am1∑i=1Nk-x0iζb1yx0icos⁡λbnyx0i+∑m=0∞am2∑i=1Nk-x0iζb2yx0icos⁡λbnyx0i=-δ-bn∑m=0∞Bmn-δ-bnβ1n∑m=0∞cm1-δ-bnβ2n∑m=0∞cm2-δ-bnbn1(59)∑m=0∞∑t=0∞Bmt∑i=1NK-x0icos⁡λbtyx0icos⁡λbnyx0i+∑m=0∞cm1∑i=1NK-x0iζb1yx0icos⁡λbnyx0i+∑m=0∞cm2∑i=1NK-x0iζb2yx0icos⁡λbnyx0i=-δ-bndn1-μδ-bn∑m=0∞∑q=0∞Cmtλbqκbnq+δ-bnμη1n∑m=0∞em1+δ-bnμη2n∑m=0∞em2(60)∑m=0∞∑t=0∞Cmt∑i=1NK-yx0icos⁡λbtyx0icos⁡λbnyx0i+∑m=0∞em1∑i=1NK-yx0iζb1yx0icos⁡λbnyx0i+∑m=0∞em2∑i=1NK-yx0iζb2yx0icos⁡λbnyx0i=--δ-bn∑m=0∞∑q=0∞Bmqλbqκbnq+δ-bnη1n∑m=0∞cm1-δ-bnη2n∑m=0∞cm2+δ-bnfn1At x=a,(61)∑m=0∞∑t=0∞-1mAmt∑i=1Nk-x0icos⁡λbtyx0icos⁡λbnyx0i+∑m=0∞-1mam1∑i=1Nk-x0iζb1yx0icos⁡λbnyx0i+∑m=0∞-1mam2∑i=1Nk-x0iζb2yx0icos⁡λbnyx0i=δ-bn∑m=0∞-1mBmn+δ-bnβ1n∑m=0∞-1mcm1+δ-bnβ2n∑m=0∞-1mcm2+δ-bnbn1(62)∑m=0∞∑q=0∞-1mBmq∑i=1NK-x0icos⁡λbnyx0icos⁡λbqyx0i+∑m=0∞-1mcm1∑i=1NK-x0iζb1yx0icos⁡λbnyx0i+∑m=0∞-1mcm2∑i=1NK-x0iζb2yx0icos⁡λbnyx0i=δ-bndn2-μδ-bn∑m=0∞∑q=0∞-1mCmqλbqκbnq-δ-bnμη1n∑m=0∞-1mem1-δ-bnμη2n∑m=0∞-1mem2(63)∑m=0∞∑q=0∞-1mCmq∑i=1NK-yx0icos⁡λbqyx0icos⁡λbnyx0i+∑m=0∞-1mem1∑i=1NK-yx0iζb1yx0icos⁡λbnyx0i+∑m=0∞-1mem2∑i=1NK-yx0iζb2yx0icos⁡λbnyx0i=-δ-bn∑m=0∞∑q=0∞-1mBmqλbqκbnq+δ-bnη1n∑m=0∞-1mcm1+δ-bnη2n∑m=0∞-1mcm2+δ-bnfn2At y=0,(64)∑s=0∞∑n=0∞Asn∑i=1Nk-y0icos⁡λasxy0icos⁡λamxy0i+∑n=0∞bn1∑i=1Nk-y0iζa1xy0icos⁡λamxy0i+∑n=0∞bn2∑i=1Nk-y0iζa2xy0icos⁡λamxy0i=-δ-am∑n=0∞Cmn+δ-amα1m∑n=0∞fn1+δ-amα2m∑n=0∞fn2+δ-amam1(65)∑s=0∞∑n=0∞Cmn∑i=1NK-y0icos⁡λasxy0icos⁡λamxy0i+∑n=0∞fn1∑i=1NK-y0iζa1xy0icos⁡λamxy0i+∑n=0∞fn2∑i=1NK-y0iζa2xy0icos⁡λamxy0i=-δ-amem1-μδ-am∑s=0∞∑n=0∞Bmnλasτams-μδ-amγ1m∑n=0∞dn1-μδ-amγ2m∑n=0∞dn2(66)∑s=0∞∑n=0∞Bsn∑i=1NK-xy0icos⁡λasxy0icos⁡λamxy0i+∑n=0∞dn1∑i=1NK-xy0iζa1xy0icos⁡λamxy0i+∑n=0∞dn2∑i=1NK-xy0iζa2xy0icos⁡λamxy0i=--δ-am∑s=0∞∑n=0∞Csnλasτams+δ-amγ1m∑n=0∞fn1+δ-amγ2m∑n=0∞fn2+δ-amcm1At y=b,(67)∑s=0∞∑n=0∞-1nAsn∑i=1Nk-y0icos⁡λasxy0icos⁡λamxy0i+∑n=0∞-1nbn1∑i=1Nk-y0iζa1xy0icos⁡λamxy0i+∑n=0∞-1nbn2∑i=1Nk-y0iζa2xy0icos⁡λamxy0i=δ-am∑n=0∞-1nCmn+δ-amα1m∑n=0∞-1nfn1+δ-amα2m∑n=0∞-1nfn2+δ-amam2(68)∑s=0∞∑n=0∞-1nCsn∑i=1NK-y0icos⁡λasxy0icos⁡λamxy0i+∑n=0∞-1nfn1∑i=1NK-y0iζa1xy0icos⁡λamxy0i+∑n=0∞-1nfn2∑i=1NK-y0iζa2xy0icos⁡λamxy0i=δ-amem2-μδ-am∑s=0∞∑n=0∞-1nBmnλasτams-μδ-amγ1m∑n=0∞-1ndn1-μδ-amγ2m∑n=0∞-1ndn2(69)∑s=0∞∑n=0∞-1nBsn∑i=1NK-xy0icos⁡λasxy0icos⁡λamxy0i+∑n=0∞-1ndn1∑i=1NK-xy0iζa1xy0icos⁡λamxy0i+∑n=0∞-1ndn2∑i=1NK-xy0iζa2xy0icos⁡λamxy0i=-δ-am∑s=0∞∑n=0∞-1nCsnλasτams+δ-amγ1m∑n=0∞-1nfn1+δ-amγ2m∑n=0∞-1nfn2+δ-amcm2In this paper, all the series expansions have been uniformly truncated tom=R and n=S. Equations (58)–(69) can be rewritten in matrix form as(70)HB=QGwhere (71)H=H1-1H1-2H1-3H1-4⋯H1-12H2-1H2-2H2-3H2-4⋯H2-12H3-1H3-2H3-3H3-4⋯H3-12H4-1H4-2H4-3H4-4⋯H4-12⋮⋮⋮⋮⋱⋮H12-1H12-2H12-3H12-4⋯H12-12Q=Q1-1Q1-2Q1-3Q2-1Q2-2Q2-3Q3-1Q3-2Q3-3Q4-1Q4-2Q4-3⋮⋮⋮Q12-1Q12-2Q12-3(72)B=a01,a11,⋯aM1,a02,a12⋯aM2,b01,b11,⋯bN1,b02,b12⋯bN2,c01,c11,⋯cM1,c02,c12⋯cM2,d01,d11,⋯dN1,d02,d12⋯dN2,e01,e11,⋯eM1,e02,e12⋯eM2,f01,f11,⋯fN1,f02,f12⋯fN2Tand (73)G=A00,A01,⋯Am′0,Am′1⋯Am′n′,⋯AMN,B00,B01,⋯Bm′0,Bm′1⋯Bm′n′,⋯BMN,B00,B01,⋯Bm′0,Bm′1⋯Bm′n′,⋯BMNAll the elements in matricesH and Q can be directly derived from (58)–(69). For example, Appendix B shows the explicit expressions that correspond to the (n+1)th equation in (58). Equations (58)–(69) represent 6(R+S+2) linear algebraic equations against a total of 3×(R+1)×(S+1)+6(R+S+2) unknown Fourier coefficients. To solve for all these unknowns, additional 3×(R+1)×(S+1) equations will have to be derived from the governing equations as described below.


### 2.2. Solving Governing Differential Equations

By substituting (45)–(47) into governing differential equations (1)–(3) and then expanding all the resulting equations and sine functions into cosine series and comparing the like terms, the following equations in matrix form will be deduced:(74)CG+DB+ρω2κGEG+FB=0where (75)C=C1-1C1-2C1-3C2-1C2-2C2-3C3-1C3-2C3-3E=E1-1E1-2E1-3E2-1E2-2E2-3E3-1E3-2E3-3(76)D=D1-1D1-2D1-3D1-4⋯D1-12D2-1D2-2D2-3D2-4⋯D2-12D3-1D3-2D3-3D3-4⋯D3-12and (77)F=F1-1F1-2F1-3F1-4⋯F1-12F2-1F2-2F2-3F2-4⋯F2-12F3-1F3-2F3-3F3-4⋯F3-12The elements in the matricesC, D, E, and F can be directly obtained from the results of substitution process. The coefficients in B cannot be treated as independent variables and need to be eliminated from (74) by making use of (70). By doing so, the final characteristic equation can be written as (78)K+ρω2κGMA=0where K=C+DH-1Q and M=E+FH-1Q, and the modal frequencies and the corresponding eigenvectors can now be easily determined by solving a standard matrix eigenproblem. The elements in each eigenvector are actually the Fourier coefficients for the corresponding mode whose mode shape in the physical space can be simply calculated using (45), (46), (47), and (70).


## 2.1. Point-Supported Edge Conditions

The rectangular Mindlin plate with arbitrary elastic point edge supports is shown in Figure1. The boundary conditions are presented by three kinds of restraining springs [29–34], namely, translational, rotational, and torsional springs. Springs are evenly arranged on each edge of Mindlin plate. Through changing the stiffness of springs, different boundary conditions can be achieved [35, 36]. The governing differential equations of a Mindlin plate about its free vibration are as follows:(1)kGh∂2w∂x2+∂2w∂y2+∂ψx∂x+∂ψy∂y-ρh∂2w∂t2=0(2)D∂2ψx∂x2+1-μ2∂2ψx∂y2+1+μ2∂2ψy∂x∂y-κGh∂w∂x+ψx-ρh∂2ψx∂t2=0(3)D∂2ψy∂y2+1-μ2∂2ψy∂x2+1+μ2∂2ψx∂x∂y-κGh∂w∂y+ψy-ρh∂2ψy∂t2=0Figure 1
A Mindlin plate with arbitrary elastic point edge supports.In the formulas above, transverse displacement is presented byw, ψx is the slope along x direction, and ψy is the slope along y direction. ρ means the mass density, μ is Poisson’s ratio, κ is the shear correction factor, and the thickness of the plate is h. G=E/(2(1+μ)) is the shear modulus, μ is Poisson’s ratio, and D=Eh3/(12(1-μ2)) is the flexural rigidity.Besides, the bending moment is expressed by (4) and (5), twisting moment is calculated by (6), and the transverse shearing forces in plates are expressed by (7) and (8).(4)Mxx=D∂ψx∂x+μ∂ψy∂y(5)Myy=D∂ψy∂y+μ∂ψx∂x(6)Mxy=D1-μ2∂ψx∂y+∂ψy∂x(7)Qx=κGhψx+∂w∂x(8)Qy=κGhψy+∂w∂yThere are three kinds of forces along every edge, namely, the bending moment, the twisting moment, and the shearing forces. Rotational, torsional, and translational springs along every edge are the counterparts to these three forces in this paper. The boundary conditions for an elastically restrained rectangular plate are as follows.Atx=0,(9)kx0yw=-Qx(10)Kx0yψx=-Mxx(11)Kyx0yψy=-Mxyk x 0 ( y ), Kx0(y), and Kyx0(y) represent the stiffness of the translation spring, rotational spring, and torsional spring, respectively, on the boundary edge x=0; the corresponding expression results in (12)kx0y=∑i=1Nkx0iδy-yx0i(13)Kx0y=∑i=1NKx0iδy-yx0i(14)Kyx0y=∑i=1NKyx0iδy-yx0iwhere N is the elastic point supports numbers on edge x=0, δ(x) is the Dirac delta function, and kx0i, Kx0i, and Kyx0i are, respectively, the stiffness of the translational, rotational, and torsional springs at the ith support located at yx0i. Similarly, the remaining boundary conditions on the other three edges can be expressed as follows.Atx=a,(15)kxayw=Qx(16)Kxayψx=Mxx(17)Kyxayψy=Mxy(18)kxay=∑i=1Nkxaiδy-yxai(19)Kxay=∑i=1NKxaiδy-yxai(20)Kyxay=∑i=1NKyxaiδy-yxaiAt y=0,(21)ky0xw=-Qy(22)Ky0xψy=-Myy(23)Kxy0xψx=-Mxy(24)ky0x=∑i=1Nky0iδx-xy0i(25)Ky0x=∑i=1NKy0iδy-yy0i(26)Kxy0x=∑i=1NKxy0iδx-xy0iAt y=b,(27)kybxw=Qy(28)Kybxψy=Myy(29)Kxybxψx=Mxy(30)kybx=∑i=1Nkybiδx-xybi(31)Kybx=∑i=1NKybiδx-xybi(32)Kyx0x=∑i=1NKyx0iδy-yybiEquations (9)–(32) represent a set of general boundary conditions by setting the spring stiffness to appropriate values. Based on (4)–(8), the boundary condition can be finally written as follows.Atx=0,(33)∑i=1Nkx0iδy-yx0iw=-κGhψx+∂w∂x(34)∑i=1NKx0iδy-yx0iψx=-D∂ψx∂x+μ∂ψy∂y(35)∑i=1NKyx0iδy-yx0iψy=-D1-μ2∂ψx∂y+∂ψy∂xAt x=a,(36)∑i=1Nkxaiδy-yxaiw=κGhψx+∂w∂x(37)∑i=1NKxaiδy-yxaiψx=D∂ψx∂x+μ∂ψy∂y(38)∑i=1NKyxaiδy-yxaiψy=D1-μ2∂ψx∂y+∂ψy∂xAt y=0,(39)∑i=1Nky0iδx-xy0iw=-κGhψy+∂w∂y(40)∑i=1NKy0iδy-yy0iψy=-D∂ψy∂y+μ∂ψx∂x(41)∑i=1NKxy0iδx-xy0iψx=-D1-μ2∂ψx∂y+∂ψy∂xAt y=b,(42)∑i=1Nkybiδx-xybiw=κGhψy+∂w∂y(43)∑i=1NKybiδx-xybiψy=D∂ψy∂y+μ∂ψx∂x(44)∑i=1NKyx0iδy-yybiψx=D1-μ2∂ψx∂y+∂ψy∂xAccording to the Mindlin plate theory, the transverse displacement of the plate middle surface and the rotations of the cross section, respectively, along thex direction and the y direction are utilized. Based on the traditional solution method, the admissible functions usually express a Fourier series expansion, because the Fourier functions constitute a complete set and exhibit an excellent numerical stability in the previous study [37–43]. We found the conventional Fourier series expression to have some defects which contain the convergence problem along the boundary edges except for a few simple boundary conditions, and the derivatives of a Fourier series cannot be obtained simply through term-by-term differentiation.In this study, in order to overcome the shortcoming with the conventional Fourier series expression, the admissible functions will be expressed as a more wholesome form of Fourier series expansion:(45)wx,y=∑m=0∞∑n=0∞Amncos⁡λamxcos⁡λbny+∑l=12ζbly∑m=0∞amlcos⁡λamx+∑l=12ζalx∑n=0∞bnlcos⁡λbny(46)ψxx,y=∑m=0∞∑n=0∞Bmncos⁡λamxcos⁡λbny+∑l=12ζbly∑m=0∞cmlcos⁡λamx+∑l=12ζalx∑n=0∞dnlcos⁡λbny(47)ψyx,y=∑m=0∞∑n=0∞Cmncos⁡λamxcos⁡λbny+∑l=12ζbly∑m=0∞emlcos⁡λamx+∑l=12ζalx∑n=0∞fnlcos⁡λbnywhere Amn, aml, bnl, Bmn, cml, dnl, Cmn, eml, and fnl represent the unknown coefficients, λam=mπ/a, λbn=nπ/b. The specific expressions of the auxiliary functions ζal and ζbl are defined as(48)ζa1x=a2πsin⁡πx2a+a2πsin⁡3πx2a(49)ζa2x=-a2πcos⁡πx2a+a2πcos⁡3πx2a(50)ζb1y=b2πsin⁡πy2b+b2πsin⁡3πy2b(51)ζb2y=-b2πcos⁡πy2b+b2πcos⁡3πy2bAs shown in (45)–(47), the supplementary functions ζa1, ζa2, ζb1, and ζb2 are used for the x and y direction displacement expressions. It is easy to find that(52)ζa10=ζa1a=ζa1′a=0,ζa1′0=1(53)ζa20=ζa2a=ζa2′0=0,ζa2′a=1Similar conditions are present in the supplementary function iny-direction, namely, ζb1y and ζb2(y). Whereas these conditions are not necessary, the existence of these conditions is purely for the sake of simplifying follow-up mathematical expressions and derivation process.Examining the admissible functions in (45)–(47), one will notice that, besides the standard 2D Fourier series defined over the domain R2:0,a⊗0,b, four additional 1D Fourier expansions are also included in the admissible functions expressions. In light of (52), it is not difficult to see that the fourth term (or the third single Fourier series) on the right side of (45) is equal to the normal derivative of the displacement function w(x,y) at edge x=0. Thus, this term will actually inherit any potential discontinuity associated with the normal derivative at the boundary x=0. Similarly, the other three terms are used to take care of the possible discontinuities at the remaining edges. Thus, the 2D Fourier expansion will now represent a periodic Π2 (residual displacement) function defined over the entire x-y plane. The convergence characteristic for the Fourier series representation of such a function has been well established in mathematics. Namely, the Fourier series will converge with a speed of mπ3 at least.By substituting (45)–(47) into boundary condition (33) and expanding all the y-related terms, except for those containing kx0(y), into Fourier cosine series, one will have(54)∑i=1Nk-x0iδy-yx0i∑m=0∞∑n=0∞Amncos⁡λbny+ζb1y∑m=0∞am1+ζb2y∑m=0∞am2=-∑m=0∞∑n=0∞Bmncos⁡λbny+∑n=0∞β1ncos⁡λbny∑m=0∞cm1+∑n=0∞β1ncos⁡λbny∑m=0∞cm2+∑n=0∞bn1cos⁡λbnyThe expressions for Fourier cosine expansion coefficient of related functions can be found in AppendixA. To establish the relationship between the Fourier coefficients in (54), the terms on the left side shall also be expressed in the form of cosine series, resulting in(55)∑m=0∞∑t=0∞Amt∑i=1Nk-x0icos⁡λbtyx0icos⁡λbnyx0i+∑m=0∞am1∑i=1Nk-x0iζb1yx0icos⁡λbnyx0i+∑m=0∞am2∑i=1Nk-x0iζb2yx0icos⁡λbnyx0i=-δ-bn∑m=0∞Bmn-δ-bnβ1n∑m=0∞cm1-δ-bnβ2n∑m=0∞cm2-δ-bnbn1(56)δ-bn=b21+δ0n(57)δ-am=a21+δ0mwhere δ0n and δ0m are Kronecker delta. Using similar procedures, the other seven boundary condition equations can be obtained. Thus, a total of twelve constraint equations can be derived as follows.Atx=0,(58)∑m=0∞∑t=0∞Amt∑i=1Nk-x0icos⁡λbtyx0icos⁡λbnyx0i+∑m=0∞am1∑i=1Nk-x0iζb1yx0icos⁡λbnyx0i+∑m=0∞am2∑i=1Nk-x0iζb2yx0icos⁡λbnyx0i=-δ-bn∑m=0∞Bmn-δ-bnβ1n∑m=0∞cm1-δ-bnβ2n∑m=0∞cm2-δ-bnbn1(59)∑m=0∞∑t=0∞Bmt∑i=1NK-x0icos⁡λbtyx0icos⁡λbnyx0i+∑m=0∞cm1∑i=1NK-x0iζb1yx0icos⁡λbnyx0i+∑m=0∞cm2∑i=1NK-x0iζb2yx0icos⁡λbnyx0i=-δ-bndn1-μδ-bn∑m=0∞∑q=0∞Cmtλbqκbnq+δ-bnμη1n∑m=0∞em1+δ-bnμη2n∑m=0∞em2(60)∑m=0∞∑t=0∞Cmt∑i=1NK-yx0icos⁡λbtyx0icos⁡λbnyx0i+∑m=0∞em1∑i=1NK-yx0iζb1yx0icos⁡λbnyx0i+∑m=0∞em2∑i=1NK-yx0iζb2yx0icos⁡λbnyx0i=--δ-bn∑m=0∞∑q=0∞Bmqλbqκbnq+δ-bnη1n∑m=0∞cm1-δ-bnη2n∑m=0∞cm2+δ-bnfn1At x=a,(61)∑m=0∞∑t=0∞-1mAmt∑i=1Nk-x0icos⁡λbtyx0icos⁡λbnyx0i+∑m=0∞-1mam1∑i=1Nk-x0iζb1yx0icos⁡λbnyx0i+∑m=0∞-1mam2∑i=1Nk-x0iζb2yx0icos⁡λbnyx0i=δ-bn∑m=0∞-1mBmn+δ-bnβ1n∑m=0∞-1mcm1+δ-bnβ2n∑m=0∞-1mcm2+δ-bnbn1(62)∑m=0∞∑q=0∞-1mBmq∑i=1NK-x0icos⁡λbnyx0icos⁡λbqyx0i+∑m=0∞-1mcm1∑i=1NK-x0iζb1yx0icos⁡λbnyx0i+∑m=0∞-1mcm2∑i=1NK-x0iζb2yx0icos⁡λbnyx0i=δ-bndn2-μδ-bn∑m=0∞∑q=0∞-1mCmqλbqκbnq-δ-bnμη1n∑m=0∞-1mem1-δ-bnμη2n∑m=0∞-1mem2(63)∑m=0∞∑q=0∞-1mCmq∑i=1NK-yx0icos⁡λbqyx0icos⁡λbnyx0i+∑m=0∞-1mem1∑i=1NK-yx0iζb1yx0icos⁡λbnyx0i+∑m=0∞-1mem2∑i=1NK-yx0iζb2yx0icos⁡λbnyx0i=-δ-bn∑m=0∞∑q=0∞-1mBmqλbqκbnq+δ-bnη1n∑m=0∞-1mcm1+δ-bnη2n∑m=0∞-1mcm2+δ-bnfn2At y=0,(64)∑s=0∞∑n=0∞Asn∑i=1Nk-y0icos⁡λasxy0icos⁡λamxy0i+∑n=0∞bn1∑i=1Nk-y0iζa1xy0icos⁡λamxy0i+∑n=0∞bn2∑i=1Nk-y0iζa2xy0icos⁡λamxy0i=-δ-am∑n=0∞Cmn+δ-amα1m∑n=0∞fn1+δ-amα2m∑n=0∞fn2+δ-amam1(65)∑s=0∞∑n=0∞Cmn∑i=1NK-y0icos⁡λasxy0icos⁡λamxy0i+∑n=0∞fn1∑i=1NK-y0iζa1xy0icos⁡λamxy0i+∑n=0∞fn2∑i=1NK-y0iζa2xy0icos⁡λamxy0i=-δ-amem1-μδ-am∑s=0∞∑n=0∞Bmnλasτams-μδ-amγ1m∑n=0∞dn1-μδ-amγ2m∑n=0∞dn2(66)∑s=0∞∑n=0∞Bsn∑i=1NK-xy0icos⁡λasxy0icos⁡λamxy0i+∑n=0∞dn1∑i=1NK-xy0iζa1xy0icos⁡λamxy0i+∑n=0∞dn2∑i=1NK-xy0iζa2xy0icos⁡λamxy0i=--δ-am∑s=0∞∑n=0∞Csnλasτams+δ-amγ1m∑n=0∞fn1+δ-amγ2m∑n=0∞fn2+δ-amcm1At y=b,(67)∑s=0∞∑n=0∞-1nAsn∑i=1Nk-y0icos⁡λasxy0icos⁡λamxy0i+∑n=0∞-1nbn1∑i=1Nk-y0iζa1xy0icos⁡λamxy0i+∑n=0∞-1nbn2∑i=1Nk-y0iζa2xy0icos⁡λamxy0i=δ-am∑n=0∞-1nCmn+δ-amα1m∑n=0∞-1nfn1+δ-amα2m∑n=0∞-1nfn2+δ-amam2(68)∑s=0∞∑n=0∞-1nCsn∑i=1NK-y0icos⁡λasxy0icos⁡λamxy0i+∑n=0∞-1nfn1∑i=1NK-y0iζa1xy0icos⁡λamxy0i+∑n=0∞-1nfn2∑i=1NK-y0iζa2xy0icos⁡λamxy0i=δ-amem2-μδ-am∑s=0∞∑n=0∞-1nBmnλasτams-μδ-amγ1m∑n=0∞-1ndn1-μδ-amγ2m∑n=0∞-1ndn2(69)∑s=0∞∑n=0∞-1nBsn∑i=1NK-xy0icos⁡λasxy0icos⁡λamxy0i+∑n=0∞-1ndn1∑i=1NK-xy0iζa1xy0icos⁡λamxy0i+∑n=0∞-1ndn2∑i=1NK-xy0iζa2xy0icos⁡λamxy0i=-δ-am∑s=0∞∑n=0∞-1nCsnλasτams+δ-amγ1m∑n=0∞-1nfn1+δ-amγ2m∑n=0∞-1nfn2+δ-amcm2In this paper, all the series expansions have been uniformly truncated tom=R and n=S. Equations (58)–(69) can be rewritten in matrix form as(70)HB=QGwhere (71)H=H1-1H1-2H1-3H1-4⋯H1-12H2-1H2-2H2-3H2-4⋯H2-12H3-1H3-2H3-3H3-4⋯H3-12H4-1H4-2H4-3H4-4⋯H4-12⋮⋮⋮⋮⋱⋮H12-1H12-2H12-3H12-4⋯H12-12Q=Q1-1Q1-2Q1-3Q2-1Q2-2Q2-3Q3-1Q3-2Q3-3Q4-1Q4-2Q4-3⋮⋮⋮Q12-1Q12-2Q12-3(72)B=a01,a11,⋯aM1,a02,a12⋯aM2,b01,b11,⋯bN1,b02,b12⋯bN2,c01,c11,⋯cM1,c02,c12⋯cM2,d01,d11,⋯dN1,d02,d12⋯dN2,e01,e11,⋯eM1,e02,e12⋯eM2,f01,f11,⋯fN1,f02,f12⋯fN2Tand (73)G=A00,A01,⋯Am′0,Am′1⋯Am′n′,⋯AMN,B00,B01,⋯Bm′0,Bm′1⋯Bm′n′,⋯BMN,B00,B01,⋯Bm′0,Bm′1⋯Bm′n′,⋯BMNAll the elements in matricesH and Q can be directly derived from (58)–(69). For example, Appendix B shows the explicit expressions that correspond to the (n+1)th equation in (58). Equations (58)–(69) represent 6(R+S+2) linear algebraic equations against a total of 3×(R+1)×(S+1)+6(R+S+2) unknown Fourier coefficients. To solve for all these unknowns, additional 3×(R+1)×(S+1) equations will have to be derived from the governing equations as described below.


## 2.2. Solving Governing Differential Equations

By substituting (45)–(47) into governing differential equations (1)–(3) and then expanding all the resulting equations and sine functions into cosine series and comparing the like terms, the following equations in matrix form will be deduced:(74)CG+DB+ρω2κGEG+FB=0where (75)C=C1-1C1-2C1-3C2-1C2-2C2-3C3-1C3-2C3-3E=E1-1E1-2E1-3E2-1E2-2E2-3E3-1E3-2E3-3(76)D=D1-1D1-2D1-3D1-4⋯D1-12D2-1D2-2D2-3D2-4⋯D2-12D3-1D3-2D3-3D3-4⋯D3-12and (77)F=F1-1F1-2F1-3F1-4⋯F1-12F2-1F2-2F2-3F2-4⋯F2-12F3-1F3-2F3-3F3-4⋯F3-12The elements in the matricesC, D, E, and F can be directly obtained from the results of substitution process. The coefficients in B cannot be treated as independent variables and need to be eliminated from (74) by making use of (70). By doing so, the final characteristic equation can be written as (78)K+ρω2κGMA=0where K=C+DH-1Q and M=E+FH-1Q, and the modal frequencies and the corresponding eigenvectors can now be easily determined by solving a standard matrix eigenproblem. The elements in each eigenvector are actually the Fourier coefficients for the corresponding mode whose mode shape in the physical space can be simply calculated using (45), (46), (47), and (70).


## 3. Numerical Examples and Discussion

In this section, some results and discussions about the free vibration of Mindlin rectangular plates are presented to verify the accuracy and flexibility of the proposed method. Based on that, some new results are obtained for Mindlin rectangular plates with various aspect ratios and thicknesses subjected to different point-supported conditions. The discussion is arranged as follows. Firstly, the convergence of the present method is checked when the plate is square rectangular and point-supported. Then, the accuracy of this method is compared with those obtained by FEM commercial program ABAQUS (S4R model). And the following results obtained by FEM are all based on the commercial program ABAQUS (S4R model) unless otherwise stated. In addition, point-supported plates with different geometry parameters and various thicknesses are calculated by the present method. Then, the method is extended to multipoint-supported plate. After the convergence of the multipoint-supported plate is verified, more multipoint-supported plates with different geometry parameters and number of clamped points are calculated through this method.If not mentioned, the material constants are chosen as follows in this study:E=211 Gpa, ρ=7800 kg/m3, and μ=0.3. To avoid round-off results, nondimensional frequency is used and stated as follows: Ω=ωb2(ρh/D)1/2/π2. The shear correction factor is selected as κ=5/6. In addition, the nondimensional translational, rotational, and torsional springs parameters are defined as ξk=ki/ϗ and ξK=Ki/ϰ, ξK=Kij/ϰ, where ϗ=104 N/m and ϰ=104 N·m/rad.


### 3.1. Convergence Studies

A Mindlin plate with arbitrary elastic point edge supports is shown in Figure1. The first example is about point-supported square rectangular plate, namely, case 1 in Figure 2. In Table 1, the first eight nondimensional frequency parameters Ω=ωb2(ρh/D)1/2/π2 for the square plate of case 1 are given with different restraining stiffness (ξki=ξKi=ξ, ξKij=0) for each point support, and the results are compared with the FEA results. FEA results of Tables 1 and 2 are both obtained by FEM commercial program ABAQUS (S4R model, 18720 elements).Table 1
Frequency parametersΩ=ωb2(ρh/D)1/2/π2 for the square plate of case 1 with restraining stiffness, ξki=ξKi=ξ, ξKij=0, for each point support.
ξ Method Model sequence 1 2 3 4 5 6 7 8 0 Present 1.2897 1.9201 2.3639 3.2356 3.2356 5.6176 5.6176 5.6459 FEA 1.2898 1.9212 2.3662 3.2379 3.2379 5.6231 5.6231 5.6522 1 Present 1.2893 1.9204 2.3640 3.2351 3.2351 5.6176 5.6176 5.6447 FEA 1.2898 1.9216 2.3663 3.2379 3.2379 5.6231 5.6231 5.6522 101 Present 1.2919 1.9235 2.3650 3.2363 3.2363 5.6198 5.6198 5.6471 FEA 1.2898 1.9247 2.3673 3.2385 3.2385 5.6244 5.6244 5.6547 102 Present 1.2893 1.9522 2.3741 3.2415 3.2415 5.6322 5.6322 5.6697 FEA 1.2898 1.9508 2.3756 3.2439 3.2439 5.6364 5.6364 5.6751 103 Present 1.2893 2.1249 2.4258 3.2771 3.2771 5.7010 5.7010 5.8077 FEA 1.2898 2.1093 2.4052 3.2671 3.2671 5.6817 5.6817 5.7570 104 Present 1.2898 2.4948 2.4988 3.3767 3.3767 5.8271 5.8271 6.0942 FEA 1.2898 2.4551 2.4870 3.3610 3.3410 5.8060 5.8060 6.0210 106 Present 1.2893 2.4807 2.4889 3.3696 3.3696 5.8123 5.8123 6.0496 FEA 1.2898 2.4834 2.5011 3.3787 3.3787 5.7769 5.7769 6.0209 107 Present 1.2887 1.7945 1.7945 2.1151 2.5149 4.9471 4.9471 6.1398 108 Present 1.2887 1.8018 1.8018 2.1325 2.5156 4.9694 4.9694 6.1398 109 Present 1.2887 1.8026 1.8026 2.1342 2.5156 4.9716 4.9716 6.1398Table 2
Frequency parametersΩ=ωb2(ρh/D)1/2/π2 for the square plate of case 1 with restraining stiffness, ξki=ξKi=104, ξKij=0, for each point support.
R, S Model sequence 1 2 3 4 5 6 7 8 R = S = 5 1.2972 2.5625 2.6799 3.4185 3.4185 5.9387 5.9387 6.1675  R = S = 6 1.2965 2.5340 2.6052 3.3997 3.3997 5.8960 5.8960 6.1669  R = S = 7 1.2919 2.5334 2.6040 3.3969 3.3969 5.8861 5.8861 6.1509  R = S = 8 1.2917 2.5164 2.5574 3.3860 3.3860 5.8612 5.8612 6.1305  R = S = 9 1.2899 2.5011 2.5562 3.3845 3.3845 5.8561 5.8561 6.1291  R = S = 10 1.2898 2.4948 2.4988 3.3767 3.3767 5.8271 5.8271 6.0942  R = S = 11 1.2891 2.4796 2.4978 3.3720 3.3720 5.8098 5.8098 6.0657  R = S = 12 1.2887 2.4795 2.4877 3.3680 3.3680 5.8096 5.8096 6.0471  R = S = 13 1.2887 2.4794 2.4877 3.3679 3.3679 5.8094 5.8094 6.0466  FEA 1.2898 2.4551 2.4870 3.3610 3.3410 5.8060 5.8060 6.0210Figure 2
Two-case elastic point support arrangement considered in the calculation.The geometric dimensions of case 1 are as follows:a=b=1 m, h=0.1 m. Besides, series expansion truncated number R=S=15. It is noted that the current results compare well with the FEA results. As shown in Table 1, when nondimensional springs parameters of translational and rotational springs are in the range of ξki=ξKi=107~109, nondimensional frequency parameter is of convergence. In other words, nondimensional springs parameters at this range can guarantee rigid boundary conditions. When nondimensional springs parameters are less than 10, the boundary conditions can be seen free. To understand the modes of square rectangular plate in free boundary conditions better, the first six mode shapes of square rectangular plates with completely free boundary conditions calculated by the present method and FEA are presented in Figures 3 and 4. And when the nondimensional springs parameters are in the range of 102 and 106, the boundary conditions are elastic. As shown in Table 2, the first eight nondimensional frequency parameters become quickly converged at R=S=10 for the given 5-digit precision. For simplicity, the displacement expansion will be truncated to R=S=10 in all the subsequent calculations.Figure 3
The first six mode shapes of square rectangular plates with completely free boundary conditions calculated by the present method.
(a)
The 1st mode (b)
The 2nd mode (c)
The 3rd mode (d)
The 4th mode (e)
The 5th mode (f)
The 6th modeFigure 4
The first six mode shapes of square rectangular plates with completely free boundary conditions calculated by FEA method.
(a)
The 1st mode (b)
The 2nd mode (c)
The 3rd mode (d)
The 4th mode (e)
The 5th mode (f)
The 6th modeTo eliminate the particularity of square rectangular plate, the second example is about point-supported rectangular plate 1b in Figure2. The geometric dimensions of case 2 are as follows: a=2 m, b=1 m, and h=0.1 m. Besides, series expansion truncated number R=S=10. In Table 3, the first eight nondimensional frequency parameters about 1b calculated by the new method and FEA are given. Meanwhile, the comparison in Table 3 shows the accuracy of the new method to ordinary rectangular Mindlin plate. FEA results of plate 1b in Table 3 and Figure 8 are obtained by FEM commercial program ABAQUS (S4R model, 37500 elements).Table 3
Frequency parametersΩ=ωb2(ρh/D)1/2/π2 for the rectangular plate of case 2 with restraining stiffness, ξki=ξKi=ξ, ξKij=0, for each point support.
ξ Method Model sequence 1 2 3 4 5 6 7 8 0 Present 0.5383 0.6459 1.4067 1.4669 2.1484 2.4072 2.4924 2.8472 FEA 0.5396 0.6458 1.4062 1.4685 2.1514 2.4047 2.4948 2.8528 1 Present 0.5387 0.6462 1.4074 1.4677 2.1496 2.4085 2.4937 2.8488 FEA 0.5388 0.6458 1.4062 1.4687 2.1512 2.4057 2.4948 2.8529 101 Present 0.5387 0.6462 1.4074 1.4677 2.1496 2.4085 2.4937 2.8488 FEA 0.5401 0.6458 1.4067 1.4698 2.1521 2.4057 2.4949 2.8543 102 Present 0.5525 0.6462 1.4127 1.4794 2.1593 2.4085 2.4949 2.8633 FEA 0.5515 0.6458 1.4111 1.4793 2.1599 2.4057 2.4960 2.8661 103 Present 0.6283 0.6462 1.4473 1.5438 2.2165 2.4085 2.5031 2.9408 FEA 0.6045 0.6458 1.4320 1.5171 2.1908 2.4057 2.5019 2.9113 104 Present 0.6462 0.8315 1.5615 1.6646 2.3366 2.4085 2.5410 3.0559 FEA 0.6458 0.8335 1.5597 1.6775 2.3375 2.4057 2.5351 3.0695 106 Present 0.6462 0.7899 0.9192 1.1024 1.6561 2.0563 2.4085 2.6004 108 Present 0.6462 0.8028 0.9396 1.1561 1.6957 2.0763 2.4085 2.7058 1010 Present 0.6462 0.8043 0.9419 1.1624 1.7001 2.0786 2.4085 2.7180


### 3.2. Point-Supported Plate

The accuracy and nice convergence characteristic of this method have been proved when the plate is point-supported. In Table4, the first eight nondimensional frequency parameters are given with different stiffness values, aspect ratios, and thickness ratios. The first six mode shapes of square rectangular plates with completely elastic boundary conditions (ξki=ξKi=102, ξKij=0) calculated by the present method and FEA are presented in Figures 5 and 6. FEA results of plate 1a in Figure 6 are obtained by FEM commercial program ABAQUS (S4R model, 18720 elements). The first six mode shapes matched well between current method and FEA. It can be seen from Table 4 that nondimensional frequency parameters under elastic boundary conditions tend to decrease with the aspect ratio, whereas the first three modes in each case basically represent the rigid-body motions. The change of these three modes is smaller than higher modes. With the increase of nondimensional spring parameters, the boundary conditions change from elastic to clamped. When the spring stiffness is set as ξki=ξKi=108, ξKij=0, the first six mode shapes of 1a and 1b which are displayed in Figures 7 and 8 can be seen as point clamped.Table 4
Frequency parametersΩ=ωb2(ρh/D)1/2/π2 for the restraining stiffness point support plate with different geometry parameters, ξki=ξKi=ξ, ξKij=0, for each point support.
ξ a / b h / b Model sequence 1 2 3 4 5 6 7 8 102 a / b = 1 h / b = 0.1 1.2890 1.9514 2.3730 3.2405 3.2405 5.6196 5.6196 5.6683 h / b = 0.2 1.1847 1.7660 2.1468 2.7913 2.7913 4.5648 4.5766 4.5766 a / b = 2 h / b = 0.1 0.5522 0.6459 1.4120 1.4786 2.1582 2.4072 2.4936 2.8618 h / b = 0.2 0.5257 0.6078 1.2908 1.3669 1.9617 2.1320 2.2271 2.5177 a / b = 5 h / b = 0.1 0.4743 0.5256 0.7779 0.8091 1.1345 1.1556 1.4961 1.6006 h / b = 0.2 0.4572 0.4914 0.7426 0.7581 1.0495 1.0856 1.3696 1.4769 104 a / b = 1 h / b = 0.1 1.2890 2.4948 2.4988 3.3720 3.3720 5.8218 5.8218 6.0670 h / b = 0.2 1.1847 1.9651 2.1982 2.8358 2.8358 4.6511 4.6511 4.7226 a / b = 2 h / b = 0.1 0.6459 0.8311 1.5607 1.6637 2.3354 2.4072 2.5397 3.0544 h / b = 0.2 0.6078 0.6183 1.3357 1.4389 2.0232 2.1320 2.2380 2.6037 a / b = 5 h / b = 0.1 0.3967 0.6293 0.6891 0.8091 0.8794 1.2422 1.2678 1.4961 h / b = 0.2 0.2805 0.5012 0.5478 0.7581 0.7766 1.0786 1.1222 1.3696 106 a / b = 1 h / b = 0.1 1.2890 1.7336 1.7336 1.9936 2.5169 4.7623 4.7623 5.6910 h / b = 0.2 1.1585 1.1763 1.1763 1.1847 2.2452 3.2096 3.3393 3.3393 a / b = 2 h / b = 0.1 0.6459 0.7895 0.9188 1.1018 1.6552 2.0553 2.4072 2.5990 h / b = 0.2 0.6078 0.6419 0.7129 0.7410 1.2653 1.7296 1.8630 2.1320 a / b = 5 h / b = 0.1 0.1945 0.2536 0.2787 0.2843 0.6492 0.7673 0.8091 0.8821 h / b = 0.2 0.1807 0.2407 0.2524 0.2540 0.5723 0.6156 0.7581 0.7802Figure 5
The first six mode shapes of square rectangular plates with completely elastic boundary conditions,ξki=ξKi=102, ξKij=0, calculated by the present method.
(a)
The 1st mode (b)
The 2nd mode (c)
The 3rd mode (d)
The 4th mode (e)
The 5th mode (f)
The 6th modeFigure 6
The first six mode shapes of square rectangular plates with completely elastic boundary conditions,ξki=ξKi=102, ξKij=0, calculated by FEA method.
(a)
The 1st mode (b)
The 2nd mode (c)
The 3rd mode (d)
The 4th mode (e)
The 5th mode (f)
The 6th modeFigure 7
The first six mode shapes of square rectangular plates with completely elastic boundary conditions,ξki=ξKi=108, ξKij=0, calculated by the present method.
(a)
The 1st mode (b)
The 2nd mode (c)
The 3rd mode (d)
The 4th mode (e)
The 5th mode (f)
The 6th modeFigure 8
The first six mode shapes of rectangular plates with completely elastic boundary conditions and the aspect ratio asa/b=2, h/b=0.1, ξki=ξKi=108, and ξKij=0, calculated by the present method.
(a)
The 1st mode (b)
The 2nd mode (c)
The 3rd mode (d)
The 4th mode (e)
The 5th mode (f)
The 6th mode


### 3.3. Multipoint-Supported Plate

It is proved that the present method is convergent and accurate when the plate is point-supported from the previous study. In this subsection, nondimensional frequency parametersΩ=ωb2(ρh/D)1/2/π2 of the multipoint-supported square plate 1a (a=1 m, b=1 m, and h=0.1 m) with different numbers of terms in the series expansion (R and S) are shown in Table 5. The nondimensional spring parameters translational and rotational springs are set as 108 N/m and 108 N·m/rad, respectively; besides, there are 26 clamped points on every edge. This shows that the results are very accurate when R and S are small numbers. When truncated numbers R and S are larger than 10, results are almost invariant. The displacement expansion will be truncated to R=S=10 in the subsequent calculations of multipoint-supported plate. For the lack of relevant literature, the FEM data is given as a comparison. The FEM data is obtained through ABAQUS (S4R); each edge is divided into 100 pieces, which is considered adequately fine to capture the spatial variations of these lower order modes. Besides, the mode shapes calculated through FEM are displayed in Figure 10. The first six mode shapes of square rectangular plates with 26 clamped points on every edge which are calculated by the present method are presented in Figure 9. Both nondimensional frequency and mode shapes of square rectangular plate matched well between the new method and FEA. The accuracy and nice convergence characteristic of this method have been proved when the plate is multipoint-supported.Table 5
Frequency parametersΩ=ωb2(ρh/D)1/2/π2 for the multipoint-supported square plate; the numbers of clamped points is 26 on every edge.
R, S Model sequence 1 2 3 4 5 6 7 8 R = S = 6 3.2837 6.2606 6.2606 8.7450 10.3330 10.4515 12.4545 12.4545 R = S = 7 3.2835 6.2569 6.2569 8.7372 10.3311 10.4495 12.4501 12.4501 R = S = 8 3.2826 6.2562 6.2562 8.7369 10.3221 10.4411 12.4411 12.4411 R = S = 9 3.2823 6.2552 6.2552 8.7349 10.3218 10.4408 12.4400 12.4400 R = S = 10 3.2823 6.2550 6.2550 8.7348 10.3201 10.4386 12.4376 12.4376 R = S = 12 3.2823 6.2550 6.2550 8.7348 10.3194 10.4386 12.4375 12.4375 R = S = 14 3.2823 6.2550 6.2550 8.7348 10.3194 10.4386 12.4375 12.4375 CC1 3.2870 6.2578 6.2578 8.7352 10.3247 10.4418 12.4412 12.4412 ∗CC1 refers to the uniformly clamped plate with a/b=1, h/b=0.1 and calculated by ABAQUS.Figure 9
The first six mode shapes of square rectangular plates with multipoint support, calculated by the present method.
(a)
The 1st mode (b)
The 2nd mode (c)
The 3rd mode (d)
The 4th mode (e)
The 5th mode (f)
The 6th modeFigure 10
The first six mode shapes of square rectangular plates with uniform completely clamped boundary conditions, calculated by FEA method.
(a)
The 1st mode (b)
The 2nd mode (c)
The 3rd mode (d)
The 4th mode (e)
The 5th mode (f)
The 6th modeIn Table6, the first eight nondimensional frequency parameters are given with different number of clamped points, aspect ratios, and thickness ratios. Similarly, nondimensional frequency parameters tend to decrease with the aspect ratio when the plate is multipoint-supported.Table 6
Frequency parametersΩ=ωb2(ρh/D)1/2/π2 for the multipoint-supported plate with different geometry parameters and numbers of clamped points.
geometry parameters No. of clamped points Model sequence 1 2 3 4 5 6 7 8 a / b = 3, h/b=0.1 4 1.5833 1.6832 2.0553 2.0597 2.4456 2.7332 3.0820 3.1456 8 2.1798 2.4123 2.6437 3.2877 4.1133 5.1177 5.4172 5.5959 12 2.1833 2.4178 2.8495 3.4943 4.3430 5.3744 5.4539 5.6759 16 2.1834 2.4181 2.8503 3.4958 4.3451 5.3773 5.4542 5.6769 20 2.1834 2.4181 2.8503 3.4958 4.3451 5.3773 5.4542 5.6769 22 2.1834 2.4181 2.8503 3.4958 4.3451 5.3773 5.4542 5.6769 26 2.1834 2.4181 2.8503 3.4958 4.3451 5.3773 5.4542 5.6769 28 2.1834 2.4181 2.8503 3.4958 4.3451 5.3773 5.4542 5.6769 30 2.1834 2.4181 2.8503 3.4958 4.3451 5.3773 5.4542 5.6769 CC1 2.1810 2.4181 2.8536 3.4987 4.3445 5.3795 5.4602 5.6776 a / b = 3, h/b=0.2 4 1.3105 1.4309 1.6865 1.7451 2.1095 2.1927 2.3602 2.4646 8 1.8389 2.0139 2.1817 2.6793 3.2958 4.0280 4.0619 4.1648 12 1.8424 2.0205 2.3551 2.8496 3.4791 4.1086 4.2146 4.2583 16 1.8426 2.0212 2.3569 2.8521 3.4824 4.1089 4.2180 4.2592 20 1.8426 2.0212 2.3569 2.8521 3.4824 4.1089 4.2180 4.2592 22 1.8426 2.0212 2.3569 2.8521 3.4824 4.1089 4.2180 4.2592 26 1.8426 2.0212 2.3569 2.8521 3.4824 4.1089 4.2180 4.2592 28 1.8426 2.0212 2.3569 2.8521 3.4824 4.1089 4.2180 4.2592 30 1.8426 2.0212 2.3569 2.8521 3.4824 4.1089 4.2180 4.2592 CC2 1.8449 2.0268 2.3604 2.8596 3.4840 4.1151 4.2190 4.2672 a / b = 5, h/b=0.2 4 0.6581 0.8020 0.9433 0.9906 1.1117 1.1914 1.3123 1.4223 8 1.6839 1.7695 1.8092 1.8715 1.9115 2.0716 2.5176 2.5228 12 1.8057 1.8566 1.9511 2.0984 2.3047 2.5709 2.8942 3.2685 16 1.8058 1.8571 1.9521 2.1004 2.3074 2.5745 2.8979 3.2729 20 1.8058 1.8571 1.9521 2.1004 2.3074 2.5745 2.8979 3.2729 22 1.8058 1.8571 1.9521 2.1004 2.3074 2.5745 2.8979 3.2729 26 1.8058 1.8571 1.9521 2.1004 2.3074 2.5745 2.8979 3.2729 28 1.8058 1.8571 1.9521 2.1004 2.3074 2.5745 2.8979 3.2729 30 1.8058 1.8571 1.9521 2.1004 2.3074 2.5745 2.8979 3.2729 CC3 1.8097 1.8629 1.9553 2.1043 2.3078 2.5626 2.9013 3.2779 ∗CC1: refers to the uniformly clamped plate with a/b=3, h/b=0.1 and calculated by ABAQUS. CC2: refers to the uniformly clamped plate with a/b=3, h/b=0.2 and calculated by ABAQUS. CC3: refers to the uniformly clamped plate with a/b=5, h/b=0.2 and calculated by ABAQUS.Figure11 shows the difference of nondimensional frequency parameters between different numbers of clamped points and clamped boundary condition, in which three kinds of aspect ratios and thickness ratios are considered. The first three frequency parameters of multipoints plate and clamped plate are, respectively, calculated by the current method and ABAQUS. With the increase of the number of clamped points on each edge, the percentage error between multipoints plate and clamped plate narrows quickly. When the number of clamped points is bigger than 16 at each edge, the discrepancy between nondimensional frequency parameters of multipoints plate and clamped plate is tiny and invariant. The maximum percentage error is 0.31%, which means multipoints boundary condition is equal to clamped boundary condition.Figure 11
The discrepancy of nondimensional frequency parameters between different numbers of clamped points and clamped boundary condition (b=1 m).
(a)
a / b = 3 , h / b = 0.1 (b)
a / b = 3 , h / b = 0.2 (c)
a / b = 5 , h / b = 0.2Based on different plate theories, the first frequency and error of multipoint clamped square plates with different thickness ratios and number of clamped points are listed in Figure12. The error is defined as Error(%)=Ωα-Ω3D/Ω3D×100%, in which α represents classical plate theory, Mindlin plate theory, and three-dimensional (3D) elastic theory, respectively. It is easy to find that classical thin plate theory overestimates the plate frequencies for the neglect of transverse shear and rotary inertia. The discrepancy of thin plates using Mindlin and classical plate theory is tiny to neglect when the thickness ratio is less than 0.05. However, the error of classical plate theory increases rapidly when the thickness ratio is more than 0.1 and is generally more than 10%. Meanwhile, the error of Mindlin plate theory is only 1.1% when the thickness ratio is 0.2. The introduction of Mindlin theory is necessary when the plate is not a thin plate (thickness ratio is less than 0.1).Figure 12
The frequency and error (%) of multipoint-supported plates with different thickness ratiosh/b and theories.
(a)
The frequency of 4 clamped points’ plate with different thickness ratiosh/b and theories (b)
The frequency of 16 clamped points’ plate with different thickness ratiosh/b and theories (c)
Error (%) of 4 clamped points’ plate with different thickness ratiosh/b and theories (d)
Error (%) of 16 clamped points’ plate with different thickness ratiosh/b and theories


## 3.1. Convergence Studies

A Mindlin plate with arbitrary elastic point edge supports is shown in Figure1. The first example is about point-supported square rectangular plate, namely, case 1 in Figure 2. In Table 1, the first eight nondimensional frequency parameters Ω=ωb2(ρh/D)1/2/π2 for the square plate of case 1 are given with different restraining stiffness (ξki=ξKi=ξ, ξKij=0) for each point support, and the results are compared with the FEA results. FEA results of Tables 1 and 2 are both obtained by FEM commercial program ABAQUS (S4R model, 18720 elements).Table 1
Frequency parametersΩ=ωb2(ρh/D)1/2/π2 for the square plate of case 1 with restraining stiffness, ξki=ξKi=ξ, ξKij=0, for each point support.
ξ Method Model sequence 1 2 3 4 5 6 7 8 0 Present 1.2897 1.9201 2.3639 3.2356 3.2356 5.6176 5.6176 5.6459 FEA 1.2898 1.9212 2.3662 3.2379 3.2379 5.6231 5.6231 5.6522 1 Present 1.2893 1.9204 2.3640 3.2351 3.2351 5.6176 5.6176 5.6447 FEA 1.2898 1.9216 2.3663 3.2379 3.2379 5.6231 5.6231 5.6522 101 Present 1.2919 1.9235 2.3650 3.2363 3.2363 5.6198 5.6198 5.6471 FEA 1.2898 1.9247 2.3673 3.2385 3.2385 5.6244 5.6244 5.6547 102 Present 1.2893 1.9522 2.3741 3.2415 3.2415 5.6322 5.6322 5.6697 FEA 1.2898 1.9508 2.3756 3.2439 3.2439 5.6364 5.6364 5.6751 103 Present 1.2893 2.1249 2.4258 3.2771 3.2771 5.7010 5.7010 5.8077 FEA 1.2898 2.1093 2.4052 3.2671 3.2671 5.6817 5.6817 5.7570 104 Present 1.2898 2.4948 2.4988 3.3767 3.3767 5.8271 5.8271 6.0942 FEA 1.2898 2.4551 2.4870 3.3610 3.3410 5.8060 5.8060 6.0210 106 Present 1.2893 2.4807 2.4889 3.3696 3.3696 5.8123 5.8123 6.0496 FEA 1.2898 2.4834 2.5011 3.3787 3.3787 5.7769 5.7769 6.0209 107 Present 1.2887 1.7945 1.7945 2.1151 2.5149 4.9471 4.9471 6.1398 108 Present 1.2887 1.8018 1.8018 2.1325 2.5156 4.9694 4.9694 6.1398 109 Present 1.2887 1.8026 1.8026 2.1342 2.5156 4.9716 4.9716 6.1398Table 2
Frequency parametersΩ=ωb2(ρh/D)1/2/π2 for the square plate of case 1 with restraining stiffness, ξki=ξKi=104, ξKij=0, for each point support.
R, S Model sequence 1 2 3 4 5 6 7 8 R = S = 5 1.2972 2.5625 2.6799 3.4185 3.4185 5.9387 5.9387 6.1675  R = S = 6 1.2965 2.5340 2.6052 3.3997 3.3997 5.8960 5.8960 6.1669  R = S = 7 1.2919 2.5334 2.6040 3.3969 3.3969 5.8861 5.8861 6.1509  R = S = 8 1.2917 2.5164 2.5574 3.3860 3.3860 5.8612 5.8612 6.1305  R = S = 9 1.2899 2.5011 2.5562 3.3845 3.3845 5.8561 5.8561 6.1291  R = S = 10 1.2898 2.4948 2.4988 3.3767 3.3767 5.8271 5.8271 6.0942  R = S = 11 1.2891 2.4796 2.4978 3.3720 3.3720 5.8098 5.8098 6.0657  R = S = 12 1.2887 2.4795 2.4877 3.3680 3.3680 5.8096 5.8096 6.0471  R = S = 13 1.2887 2.4794 2.4877 3.3679 3.3679 5.8094 5.8094 6.0466  FEA 1.2898 2.4551 2.4870 3.3610 3.3410 5.8060 5.8060 6.0210Figure 2
Two-case elastic point support arrangement considered in the calculation.The geometric dimensions of case 1 are as follows:a=b=1 m, h=0.1 m. Besides, series expansion truncated number R=S=15. It is noted that the current results compare well with the FEA results. As shown in Table 1, when nondimensional springs parameters of translational and rotational springs are in the range of ξki=ξKi=107~109, nondimensional frequency parameter is of convergence. In other words, nondimensional springs parameters at this range can guarantee rigid boundary conditions. When nondimensional springs parameters are less than 10, the boundary conditions can be seen free. To understand the modes of square rectangular plate in free boundary conditions better, the first six mode shapes of square rectangular plates with completely free boundary conditions calculated by the present method and FEA are presented in Figures 3 and 4. And when the nondimensional springs parameters are in the range of 102 and 106, the boundary conditions are elastic. As shown in Table 2, the first eight nondimensional frequency parameters become quickly converged at R=S=10 for the given 5-digit precision. For simplicity, the displacement expansion will be truncated to R=S=10 in all the subsequent calculations.Figure 3
The first six mode shapes of square rectangular plates with completely free boundary conditions calculated by the present method.
(a)
The 1st mode (b)
The 2nd mode (c)
The 3rd mode (d)
The 4th mode (e)
The 5th mode (f)
The 6th modeFigure 4
The first six mode shapes of square rectangular plates with completely free boundary conditions calculated by FEA method.
(a)
The 1st mode (b)
The 2nd mode (c)
The 3rd mode (d)
The 4th mode (e)
The 5th mode (f)
The 6th modeTo eliminate the particularity of square rectangular plate, the second example is about point-supported rectangular plate 1b in Figure2. The geometric dimensions of case 2 are as follows: a=2 m, b=1 m, and h=0.1 m. Besides, series expansion truncated number R=S=10. In Table 3, the first eight nondimensional frequency parameters about 1b calculated by the new method and FEA are given. Meanwhile, the comparison in Table 3 shows the accuracy of the new method to ordinary rectangular Mindlin plate. FEA results of plate 1b in Table 3 and Figure 8 are obtained by FEM commercial program ABAQUS (S4R model, 37500 elements).Table 3
Frequency parametersΩ=ωb2(ρh/D)1/2/π2 for the rectangular plate of case 2 with restraining stiffness, ξki=ξKi=ξ, ξKij=0, for each point support.
ξ Method Model sequence 1 2 3 4 5 6 7 8 0 Present 0.5383 0.6459 1.4067 1.4669 2.1484 2.4072 2.4924 2.8472 FEA 0.5396 0.6458 1.4062 1.4685 2.1514 2.4047 2.4948 2.8528 1 Present 0.5387 0.6462 1.4074 1.4677 2.1496 2.4085 2.4937 2.8488 FEA 0.5388 0.6458 1.4062 1.4687 2.1512 2.4057 2.4948 2.8529 101 Present 0.5387 0.6462 1.4074 1.4677 2.1496 2.4085 2.4937 2.8488 FEA 0.5401 0.6458 1.4067 1.4698 2.1521 2.4057 2.4949 2.8543 102 Present 0.5525 0.6462 1.4127 1.4794 2.1593 2.4085 2.4949 2.8633 FEA 0.5515 0.6458 1.4111 1.4793 2.1599 2.4057 2.4960 2.8661 103 Present 0.6283 0.6462 1.4473 1.5438 2.2165 2.4085 2.5031 2.9408 FEA 0.6045 0.6458 1.4320 1.5171 2.1908 2.4057 2.5019 2.9113 104 Present 0.6462 0.8315 1.5615 1.6646 2.3366 2.4085 2.5410 3.0559 FEA 0.6458 0.8335 1.5597 1.6775 2.3375 2.4057 2.5351 3.0695 106 Present 0.6462 0.7899 0.9192 1.1024 1.6561 2.0563 2.4085 2.6004 108 Present 0.6462 0.8028 0.9396 1.1561 1.6957 2.0763 2.4085 2.7058 1010 Present 0.6462 0.8043 0.9419 1.1624 1.7001 2.0786 2.4085 2.7180


## 3.2. Point-Supported Plate

The accuracy and nice convergence characteristic of this method have been proved when the plate is point-supported. In Table4, the first eight nondimensional frequency parameters are given with different stiffness values, aspect ratios, and thickness ratios. The first six mode shapes of square rectangular plates with completely elastic boundary conditions (ξki=ξKi=102, ξKij=0) calculated by the present method and FEA are presented in Figures 5 and 6. FEA results of plate 1a in Figure 6 are obtained by FEM commercial program ABAQUS (S4R model, 18720 elements). The first six mode shapes matched well between current method and FEA. It can be seen from Table 4 that nondimensional frequency parameters under elastic boundary conditions tend to decrease with the aspect ratio, whereas the first three modes in each case basically represent the rigid-body motions. The change of these three modes is smaller than higher modes. With the increase of nondimensional spring parameters, the boundary conditions change from elastic to clamped. When the spring stiffness is set as ξki=ξKi=108, ξKij=0, the first six mode shapes of 1a and 1b which are displayed in Figures 7 and 8 can be seen as point clamped.Table 4
Frequency parametersΩ=ωb2(ρh/D)1/2/π2 for the restraining stiffness point support plate with different geometry parameters, ξki=ξKi=ξ, ξKij=0, for each point support.
ξ a / b h / b Model sequence 1 2 3 4 5 6 7 8 102 a / b = 1 h / b = 0.1 1.2890 1.9514 2.3730 3.2405 3.2405 5.6196 5.6196 5.6683 h / b = 0.2 1.1847 1.7660 2.1468 2.7913 2.7913 4.5648 4.5766 4.5766 a / b = 2 h / b = 0.1 0.5522 0.6459 1.4120 1.4786 2.1582 2.4072 2.4936 2.8618 h / b = 0.2 0.5257 0.6078 1.2908 1.3669 1.9617 2.1320 2.2271 2.5177 a / b = 5 h / b = 0.1 0.4743 0.5256 0.7779 0.8091 1.1345 1.1556 1.4961 1.6006 h / b = 0.2 0.4572 0.4914 0.7426 0.7581 1.0495 1.0856 1.3696 1.4769 104 a / b = 1 h / b = 0.1 1.2890 2.4948 2.4988 3.3720 3.3720 5.8218 5.8218 6.0670 h / b = 0.2 1.1847 1.9651 2.1982 2.8358 2.8358 4.6511 4.6511 4.7226 a / b = 2 h / b = 0.1 0.6459 0.8311 1.5607 1.6637 2.3354 2.4072 2.5397 3.0544 h / b = 0.2 0.6078 0.6183 1.3357 1.4389 2.0232 2.1320 2.2380 2.6037 a / b = 5 h / b = 0.1 0.3967 0.6293 0.6891 0.8091 0.8794 1.2422 1.2678 1.4961 h / b = 0.2 0.2805 0.5012 0.5478 0.7581 0.7766 1.0786 1.1222 1.3696 106 a / b = 1 h / b = 0.1 1.2890 1.7336 1.7336 1.9936 2.5169 4.7623 4.7623 5.6910 h / b = 0.2 1.1585 1.1763 1.1763 1.1847 2.2452 3.2096 3.3393 3.3393 a / b = 2 h / b = 0.1 0.6459 0.7895 0.9188 1.1018 1.6552 2.0553 2.4072 2.5990 h / b = 0.2 0.6078 0.6419 0.7129 0.7410 1.2653 1.7296 1.8630 2.1320 a / b = 5 h / b = 0.1 0.1945 0.2536 0.2787 0.2843 0.6492 0.7673 0.8091 0.8821 h / b = 0.2 0.1807 0.2407 0.2524 0.2540 0.5723 0.6156 0.7581 0.7802Figure 5
The first six mode shapes of square rectangular plates with completely elastic boundary conditions,ξki=ξKi=102, ξKij=0, calculated by the present method.
(a)
The 1st mode (b)
The 2nd mode (c)
The 3rd mode (d)
The 4th mode (e)
The 5th mode (f)
The 6th modeFigure 6
The first six mode shapes of square rectangular plates with completely elastic boundary conditions,ξki=ξKi=102, ξKij=0, calculated by FEA method.
(a)
The 1st mode (b)
The 2nd mode (c)
The 3rd mode (d)
The 4th mode (e)
The 5th mode (f)
The 6th modeFigure 7
The first six mode shapes of square rectangular plates with completely elastic boundary conditions,ξki=ξKi=108, ξKij=0, calculated by the present method.
(a)
The 1st mode (b)
The 2nd mode (c)
The 3rd mode (d)
The 4th mode (e)
The 5th mode (f)
The 6th modeFigure 8
The first six mode shapes of rectangular plates with completely elastic boundary conditions and the aspect ratio asa/b=2, h/b=0.1, ξki=ξKi=108, and ξKij=0, calculated by the present method.
(a)
The 1st mode (b)
The 2nd mode (c)
The 3rd mode (d)
The 4th mode (e)
The 5th mode (f)
The 6th mode


## 3.3. Multipoint-Supported Plate

It is proved that the present method is convergent and accurate when the plate is point-supported from the previous study. In this subsection, nondimensional frequency parametersΩ=ωb2(ρh/D)1/2/π2 of the multipoint-supported square plate 1a (a=1 m, b=1 m, and h=0.1 m) with different numbers of terms in the series expansion (R and S) are shown in Table 5. The nondimensional spring parameters translational and rotational springs are set as 108 N/m and 108 N·m/rad, respectively; besides, there are 26 clamped points on every edge. This shows that the results are very accurate when R and S are small numbers. When truncated numbers R and S are larger than 10, results are almost invariant. The displacement expansion will be truncated to R=S=10 in the subsequent calculations of multipoint-supported plate. For the lack of relevant literature, the FEM data is given as a comparison. The FEM data is obtained through ABAQUS (S4R); each edge is divided into 100 pieces, which is considered adequately fine to capture the spatial variations of these lower order modes. Besides, the mode shapes calculated through FEM are displayed in Figure 10. The first six mode shapes of square rectangular plates with 26 clamped points on every edge which are calculated by the present method are presented in Figure 9. Both nondimensional frequency and mode shapes of square rectangular plate matched well between the new method and FEA. The accuracy and nice convergence characteristic of this method have been proved when the plate is multipoint-supported.Table 5
Frequency parametersΩ=ωb2(ρh/D)1/2/π2 for the multipoint-supported square plate; the numbers of clamped points is 26 on every edge.
R, S Model sequence 1 2 3 4 5 6 7 8 R = S = 6 3.2837 6.2606 6.2606 8.7450 10.3330 10.4515 12.4545 12.4545 R = S = 7 3.2835 6.2569 6.2569 8.7372 10.3311 10.4495 12.4501 12.4501 R = S = 8 3.2826 6.2562 6.2562 8.7369 10.3221 10.4411 12.4411 12.4411 R = S = 9 3.2823 6.2552 6.2552 8.7349 10.3218 10.4408 12.4400 12.4400 R = S = 10 3.2823 6.2550 6.2550 8.7348 10.3201 10.4386 12.4376 12.4376 R = S = 12 3.2823 6.2550 6.2550 8.7348 10.3194 10.4386 12.4375 12.4375 R = S = 14 3.2823 6.2550 6.2550 8.7348 10.3194 10.4386 12.4375 12.4375 CC1 3.2870 6.2578 6.2578 8.7352 10.3247 10.4418 12.4412 12.4412 ∗CC1 refers to the uniformly clamped plate with a/b=1, h/b=0.1 and calculated by ABAQUS.Figure 9
The first six mode shapes of square rectangular plates with multipoint support, calculated by the present method.
(a)
The 1st mode (b)
The 2nd mode (c)
The 3rd mode (d)
The 4th mode (e)
The 5th mode (f)
The 6th modeFigure 10
The first six mode shapes of square rectangular plates with uniform completely clamped boundary conditions, calculated by FEA method.
(a)
The 1st mode (b)
The 2nd mode (c)
The 3rd mode (d)
The 4th mode (e)
The 5th mode (f)
The 6th modeIn Table6, the first eight nondimensional frequency parameters are given with different number of clamped points, aspect ratios, and thickness ratios. Similarly, nondimensional frequency parameters tend to decrease with the aspect ratio when the plate is multipoint-supported.Table 6
Frequency parametersΩ=ωb2(ρh/D)1/2/π2 for the multipoint-supported plate with different geometry parameters and numbers of clamped points.
geometry parameters No. of clamped points Model sequence 1 2 3 4 5 6 7 8 a / b = 3, h/b=0.1 4 1.5833 1.6832 2.0553 2.0597 2.4456 2.7332 3.0820 3.1456 8 2.1798 2.4123 2.6437 3.2877 4.1133 5.1177 5.4172 5.5959 12 2.1833 2.4178 2.8495 3.4943 4.3430 5.3744 5.4539 5.6759 16 2.1834 2.4181 2.8503 3.4958 4.3451 5.3773 5.4542 5.6769 20 2.1834 2.4181 2.8503 3.4958 4.3451 5.3773 5.4542 5.6769 22 2.1834 2.4181 2.8503 3.4958 4.3451 5.3773 5.4542 5.6769 26 2.1834 2.4181 2.8503 3.4958 4.3451 5.3773 5.4542 5.6769 28 2.1834 2.4181 2.8503 3.4958 4.3451 5.3773 5.4542 5.6769 30 2.1834 2.4181 2.8503 3.4958 4.3451 5.3773 5.4542 5.6769 CC1 2.1810 2.4181 2.8536 3.4987 4.3445 5.3795 5.4602 5.6776 a / b = 3, h/b=0.2 4 1.3105 1.4309 1.6865 1.7451 2.1095 2.1927 2.3602 2.4646 8 1.8389 2.0139 2.1817 2.6793 3.2958 4.0280 4.0619 4.1648 12 1.8424 2.0205 2.3551 2.8496 3.4791 4.1086 4.2146 4.2583 16 1.8426 2.0212 2.3569 2.8521 3.4824 4.1089 4.2180 4.2592 20 1.8426 2.0212 2.3569 2.8521 3.4824 4.1089 4.2180 4.2592 22 1.8426 2.0212 2.3569 2.8521 3.4824 4.1089 4.2180 4.2592 26 1.8426 2.0212 2.3569 2.8521 3.4824 4.1089 4.2180 4.2592 28 1.8426 2.0212 2.3569 2.8521 3.4824 4.1089 4.2180 4.2592 30 1.8426 2.0212 2.3569 2.8521 3.4824 4.1089 4.2180 4.2592 CC2 1.8449 2.0268 2.3604 2.8596 3.4840 4.1151 4.2190 4.2672 a / b = 5, h/b=0.2 4 0.6581 0.8020 0.9433 0.9906 1.1117 1.1914 1.3123 1.4223 8 1.6839 1.7695 1.8092 1.8715 1.9115 2.0716 2.5176 2.5228 12 1.8057 1.8566 1.9511 2.0984 2.3047 2.5709 2.8942 3.2685 16 1.8058 1.8571 1.9521 2.1004 2.3074 2.5745 2.8979 3.2729 20 1.8058 1.8571 1.9521 2.1004 2.3074 2.5745 2.8979 3.2729 22 1.8058 1.8571 1.9521 2.1004 2.3074 2.5745 2.8979 3.2729 26 1.8058 1.8571 1.9521 2.1004 2.3074 2.5745 2.8979 3.2729 28 1.8058 1.8571 1.9521 2.1004 2.3074 2.5745 2.8979 3.2729 30 1.8058 1.8571 1.9521 2.1004 2.3074 2.5745 2.8979 3.2729 CC3 1.8097 1.8629 1.9553 2.1043 2.3078 2.5626 2.9013 3.2779 ∗CC1: refers to the uniformly clamped plate with a/b=3, h/b=0.1 and calculated by ABAQUS. CC2: refers to the uniformly clamped plate with a/b=3, h/b=0.2 and calculated by ABAQUS. CC3: refers to the uniformly clamped plate with a/b=5, h/b=0.2 and calculated by ABAQUS.Figure11 shows the difference of nondimensional frequency parameters between different numbers of clamped points and clamped boundary condition, in which three kinds of aspect ratios and thickness ratios are considered. The first three frequency parameters of multipoints plate and clamped plate are, respectively, calculated by the current method and ABAQUS. With the increase of the number of clamped points on each edge, the percentage error between multipoints plate and clamped plate narrows quickly. When the number of clamped points is bigger than 16 at each edge, the discrepancy between nondimensional frequency parameters of multipoints plate and clamped plate is tiny and invariant. The maximum percentage error is 0.31%, which means multipoints boundary condition is equal to clamped boundary condition.Figure 11
The discrepancy of nondimensional frequency parameters between different numbers of clamped points and clamped boundary condition (b=1 m).
(a)
a / b = 3 , h / b = 0.1 (b)
a / b = 3 , h / b = 0.2 (c)
a / b = 5 , h / b = 0.2Based on different plate theories, the first frequency and error of multipoint clamped square plates with different thickness ratios and number of clamped points are listed in Figure12. The error is defined as Error(%)=Ωα-Ω3D/Ω3D×100%, in which α represents classical plate theory, Mindlin plate theory, and three-dimensional (3D) elastic theory, respectively. It is easy to find that classical thin plate theory overestimates the plate frequencies for the neglect of transverse shear and rotary inertia. The discrepancy of thin plates using Mindlin and classical plate theory is tiny to neglect when the thickness ratio is less than 0.05. However, the error of classical plate theory increases rapidly when the thickness ratio is more than 0.1 and is generally more than 10%. Meanwhile, the error of Mindlin plate theory is only 1.1% when the thickness ratio is 0.2. The introduction of Mindlin theory is necessary when the plate is not a thin plate (thickness ratio is less than 0.1).Figure 12
The frequency and error (%) of multipoint-supported plates with different thickness ratiosh/b and theories.
(a)
The frequency of 4 clamped points’ plate with different thickness ratiosh/b and theories (b)
The frequency of 16 clamped points’ plate with different thickness ratiosh/b and theories (c)
Error (%) of 4 clamped points’ plate with different thickness ratiosh/b and theories (d)
Error (%) of 16 clamped points’ plate with different thickness ratiosh/b and theories


## 4. Conclusions

In this paper, a modified Fourier method has been presented to study the free vibration behaviors of moderately thick rectangular plates with different point-supported conditions. The first-order shear deformation plate theory is adopted to formulate the theoretical model. The displacement and rotation fields of plates, regardless of boundary conditions, are generally sought as a new form of trigonometric series expansions in which several supplementary closed functions are introduced to ensure and accelerate the convergence of the series expansion. Not only is the series expansion representation of solution applicable to any boundary conditions, but also the convergence of the series expansion can be substantially improved. Rayleigh–Ritz method is employed to obtain solution by the energy description of the plates. The convergence of the present solution is examined and the excellent accuracy is validated by comparison with FEM data. Excellent agreements are obtained from these comparisons. The proposed method provides a unified means for extracting the modal parameters and predicting the vibration behaviors of moderately thick plates with arbitrary point-supported edge restraints. A variety of free vibration results for moderately thick rectangular plates with different aspect ratios, thickness ratios, and boundary conditions are presented. From the results in this paper, we can find that nondimensional frequency tends to decrease with the aspect ratio, whereas the change of these three modes is smaller than higher modes under elastic boundary conditions. In addition, with the increase of the number of clamped points on edge of rectangular plate, the boundary condition converges to fully clamped edge condition. Finally, we have verified that classical thin plate theory overestimates the plate frequencies for the neglect of transverse shear and rotary inertia. It is necessary to introduce Mindlin theory when considering the vibration of moderately thick plate.



---

*Source: 8562079-2018-06-28.xml*