Frequency Domain Spectral Element Model for the Vibration Analysis of a Thin Plate with Arbitrary Boundary Conditions

We propose a new spectral elementmodel for finite rectangular plate elements with arbitrary boundary conditions.Thenew spectral elementmodel is developed bymodifying the boundary splittingmethod used in our previous study so that the four corner nodes of a finite rectangular plate element become active.Thus, the new spectral element model can be applied to any finite rectangular plate element with arbitrary boundary conditions, while the spectral element model introduced in the our previous study is valid only for finite rectangular plate elements with four fixed corner nodes. The new spectral element model can be used as a generic finite element model because it can be assembled in any plate direction. The accuracy and computational efficiency of the new spectral element model are validated by a comparison with exact solutions, solutions obtained by the standard finite element method, and solutions from the commercial finite element analysis package ANSYS.


Introduction
The plate is a representative structural element that is widely used in many engineering fields such as mechanical, civil, aerospace, shipbuilding, and structural engineering.Severe or unwanted vibration of a plate is a very important engineering problem.Thus, it is required to accurately predict the vibration characteristics of a plate during the design phase.Exact solutions are available only for Levy-type plates [1,2].Thus, numerous computational methods have been developed for the vibrations of plates during the last two centuries.
The finite element method (FEM) is one of the most widely used computational methods that can be applied to various complex structures including the plates.The FEM in general provides reliable solutions in the low frequency range, but poor solutions in the high frequency range.Thus, to improve the solution accuracy in the high frequency range, a finite element must be divided into many smaller finite elements so that their sizes are smaller than the wavelengths of the highest vibration mode of interest.However, this will result in a significant increase in computation cost.Thus, as an alternative to FEM, we can consider the spectral element method (SEM) for the vibration analysis of plates.
The SEM considered in this study is the fast Fourier transform-(FFT-) based frequency domain analysis method [3,4].The spectral element matrix (or exact dynamic stiffness matrix) used in the SEM is formulated from free wave solutions that satisfy governing differential equations of motion in the frequency domain.Thus, compared with FEM, the SEM can provide exact solutions by representing a uniform structure as a single finite element, regardless of the size of the uniform structure.Accordingly the SEM is known as an exact solution method that has the flexibility of FEM and the exactness of continuum elements [3].
Despite the outstanding features of the SEM, it is mostly used in one-dimensional (1D) structures [3,4].The SEM application to two-dimensional (2D) structures such as plates has been limited to very specific geometries and boundary conditions, for example, Levy-type plates [6][7][8][9][10] and infinite or semi-infinite plates [11][12][13][14][15].Some researchers [16,17] have introduced the spectral super element method (SSEM) for rectangular plates with prespecified boundary conditions on two parallel edges in one direction (say, the -direction).However, as their spectral element models can be assembled only in another direction (the -direction), their applications must be limited to very specific boundary conditions.Recently, Park et al. [5,18,19] developed spectral element models for rectangular membrane, isotropic plate, and orthotropic laminated composite plate elements by using two methods in combination: the boundary splitting method [20] and the spectral super element method (SSEM) [16].They derived frequency-dependent shape functions by applying a Kantorovich method-based finite strip element method in one direction and the SEM in another direction in combination, and vice versa.Accordingly, their spectral element models can be assembled in both the and directions.However, the spectral element models by Park et al. [5,18,19] still have some limitations because they are valid only for finite rectangular membranes or plate elements whose four corner nodes are fixed.To the authors' best knowledge, there have been no reports on a generic type of spectral element model that can be assembled in any direction of a plate subjected to arbitrary boundary conditions.
Thus, the purpose of this study is to develop a new spectral element model for finite rectangular plate elements that can be applied to any plate subjected to arbitrary boundary conditions.The new spectral element model is developed by modifying the boundary splitting technique used in our previous study [5] so that the four corner nodes of a finite rectangular plate element become active.The performance of the new spectral element model is evaluated by comparison with exact solutions, FEM solutions, and solutions using the commercial finite element analysis package ANSYS [21].

Spectral Element Model for a Finite Plate Element
2.1.Governing Equations in the Frequency Domain.The time domain equation of motion and the boundary conditions of plate structures with transverse vibrations are described in [2].The time domain equation of motion of a plate can be transformed into a frequency domain equation of motion of the plate by using the FFT as follows [5]: where (, ) is the transverse displacement in spectral form, (, ) is the external force in spectral form,  is the mass per unit area of the plate, and  = ℎ 3 /[12(1 − ] 2 )] is the flexural bending rigidity of the plate where  is the modulus of elasticity, ] is Poisson's ratio, and h is the plate thickness.Similarly, the time domain boundary conditions can be transformed into frequency domain boundary conditions as follows: where   and   are the dimensions of a finite plate in the and -directions, respectively;   and   are the resultant moments; and   and   are the resultant transverse shear forces defined by And   and   are the slopes defined by We need to obtain frequency domain free wave solutions for a homogeneous equation of motion in order to formulate : Boundary splitting method used in the previous study [5] to derive (, ) for a rectangular plate element subjected to arbitrary boundary conditions (e: active nodes; O: fixed nodes).the spectral element model for a finite plate element.To realize this, the homogeneous equation of motion is considered by removing external force (, ) in (1) as follows: The weak form of (5) can be obtained in the following form: A free vibration solution satisfying the weak form given in (6) can be obtained approximately by using two combined methods: the boundary splitting method [20] and the spectral super element method (SSEM) [16].The SSEM uses a combination of the Kantorovich method (based on the finite strip element method) and the frequency domain 1D spectral element method.
The concept of the boundary splitting method is illustrated in Figures 1 and 2. Figure 1 indicates the concept used in our previous study [5].Figure 2 indicates the concept used in the present study.The original problems, shown in Figures 1(a) and 2(a), are represented by the sum of two partial problems, Problem  and Problem .In Figures 1 and 2, the geometric boundary conditions of the original problems are presented in simple forms by using the following definitions: In our previous study [5], Problem , shown in Figure 1(b), has fixed (null) boundary conditions on the parallel edges at  = −  /2 and   /2.Problem , shown in Figure 1(c), has fixed (null) boundary conditions on the parallel edges at  = −  /2 and   /2.As a result, the spectral element model developed in [5] is valid only for finite rectangular plate elements whose four corner nodes are fixed, as shown in Figure 1(d).Accordingly, an application of this approach should be limited to very specific problems as considered in [5].
We propose a new boundary splitting method by modifying the boundary splitting method used in [5] such that the four corner nodes of a finite plate element become active.Problem , shown in Figure 2 Accordingly, compared to the spectral element model developed in our previous study [5] based on the boundary splitting shown in Figure 1, the present spectral element model that was developed based on the boundary splitting shown in Figure 2 has four active corner nodes.Thus, it can be used as a generic finite element model that can be assembled in both the and -directions of a plate with arbitrary boundary conditions.

Derivation of 𝑤 𝐴 (𝑥, 𝑦).
To obtain the solution   (, ) for Problem  by using the SSEM, a rectangular finite plate element is divided into   finite strip elements in the direction, as shown in Figure 3(a).The th finite strip element, which has a width of  ()  =   −  −1 in the -direction, is shown in Figure 3(b).
The displacement field in the th finite strip element can be represented by (9) where Y yA (x) y j−1 y j w (2)  yA (x) w (1)  yA (x) where By using ( 9), the displacement field   (, ) over the entire domain of the finite plate element can be represented as where w  () = {w (0)  () , w (1)  () , . . ., w with In ( 16), ℎ () are functions defined by where () is the Heaviside unit step function.Substituting ( 13) into (5) yields where with the following definitions: The constant matrices Λ 1 , Λ 2 , Λ 3 , and Λ 4 are provided in Appendix A.
Next, we assume solutions of ( 18) to be in the following form: where   is a constant,   is the wavenumber in the direction, and Substituting ( 21) into (18) gives the following eigenvalue problem: or The dispersion relation (i.e., the frequency-wavenumber relationship) can be obtained from (24) as follows: From ( 25), the wavenumbers can be computed as By using the wavenumbers  () given by (26), we can write the general solution of (18) in the following form: where with In (29), r () is the th eigenvector, which can be readily computed from (23) using   =  () .
The nodal DOFs at the nodes defined on the edges at  = −  /2 and   /2 can be written in vector form as where Here, the superscripts L and R denote the nodal values on the left edge (i.e., at  = −  /2) and the right edge (i.e., at  =   /2) of the plate, respectively.Superscript A denotes the quantities related to or contributed by   (, ).By substituting (27) into (31), the nodal DOF vector d  can be written in terms of the constant vector a  as follows: ] where The constant vector a  can be removed from (27) by using (32) to obtain the following expression: where From (12), by using the nodal DOFs defined in Figure 4, we obtain the following expressions: (3) By applying (36) to ( 14), we obtain the following expressions: We rearrange the order of nodal DOFs in (38) to define a new nodal DOF vector as follows: (39) The nodal DOF vector d  can be related to the new nodal DOF vector d  as follows: where T 1 is the transformation matrix defined by where The matrices T 1 and T 2 in (42) are (  +1)-by-(  +1) block diagonal matrices defined by where Applying (40) to (34) gives where Finally, substituting (45) into (13) gives where By following the solution procedure used for Problem , the displacement field in the ith finite strip element, which has a width of  ()   =   −  −1 , can be written in the following form: where X ()  () is a one-by-four interpolation function matrix and w ()   () are the nodal line DOF functions defined by where By using (49), the displacement field   (, ) in the whole domain of the finite plate element can be represented by where w  () = {w (1)  () , w (2)  () , . . ., w ()  () , . . ., In (54), the following definition is used: where and () is the Heaviside unit step function.Note that the null boundary conditions at  = −  /2 and  =   /2 have been applied to obtain (52).
x y (1) ( 2) Substituting (52) into (5) gives where with the following definitions: The constant matrices Λ 1 , Λ 2 , Λ 3 , and Λ 4 are provided in Appendix B. Now we assume solutions to (57) in the following form: where   is a constant,   is the wavenumber in the direction, and Substituting (60) into (57) gives the following eigenvalue problem: or The dispersion relation can be obtained from (63) as follows: From (64), the wavenumbers can be readily computed in the following forms: ( = 1, 2, . . ., 4 (  − 1)) . (65) By using the 8 × (  − 1) wavenumbers computed from (64), we can write the general solution of (57) in the following form: where with In (61), r () is the th eigenvector that can be readily computed from (63) using   =  () .
By using (47), the nodal values contributed by   (, ) at the th nodes on the bottom edge at  = −  /2 and the upper edge at  =   /2 can be related to the nodal DOF vector d  as follows: ] ] or where the primes (  ) denote the derivatives with respect to  or , and The superscripts B and U denote the quantities at the bottom edge (i.e., at  = −  /2) and the upper edge (i.e., at  =   /2) of the plate, respectively.The superscripts A and B denote the quantities related to or contributed by   (, ) and   (, ), respectively.
By using (70), the nodal values contributed by   (, ) at all nodes on the bottom and upper edges, except for four corner nodes, can be written in vector form as where By using (70), (72) can be written in terms of d  as follows: where with Similarly, the nodal values contributed by   (, ) at all nodes on the bottom and upper edges can be computed from (66), and they can be written in the vector form as where By substituting (66) into (78), the nodal values d can be written in terms of the constant vector a  as follows: ] where

Mathematical Problems in Engineering
We rearrange the order of nodal values in (77) to define a new vector as where (2) , . . .,   () ,   () ,   () ,   () , . . .,   (  −2) , By introducing a proper transformation matrix, the nodal DOF vector d can be related to the new vector d in the following form: where T 2 is the transformation matrix defined by where and T 1 and T 2 are (  − 1)-by-(  − 1) block diagonal matrices defined by where By substituting (80) into (84), we obtain The nodal DOFs defined at the nodes on the bottom and upper edges of the finite plate element (see Figure 4) can be written in vector form as where The nodal DOF vector d  must be identical to the sum of the nodal values contributed by   (, ) and   (, ).Thus, Substituting ( 74) and ( 89) into (92) gives From (93), we obtain the constant vector a  as The constant vector a  can be removed from (66) by using (94) to obtain where Finally, substituting (95) into (52) yields where where d() is the 8(  +   )-by-one spectral nodal DOF vector defined by where  1 (),  2 (),  1 (), and  2 () are the resultant transverse shear forces acting on the four boundary edges.
By substituting (101) into (102), we obtain the spectral element equation in the following form: where In the preceding equations, the following definitions are used: where I  is the 8(  + 1)-by-8(  + 1) identity matrix and I  is the 8(  − 1)-by-8(  − 1) identity matrix.The matrix S() in ( 103) is the 8(  +   )-by-8(  +   ) symmetric dynamic stiffness matrix (often called the "spectral element matrix" in the literature), and matrix D() is defined by where where ] , ] , with the use of the following definitions: In (109), the symbol vect( ) denotes a vector-form representation of the components of a diagonal matrix, and the symbol (. * ) denotes the elementwise matrix multiplication defined in MATLAB5 [23] as follows: where the components of matrix C are defined by

Numerical Results and Discussion
To evaluate the performance of the present spectral element model, a square plate is considered as a numerical example.The plate is made of aluminum, and its material properties are as follows: Young's modulus E = 69 GPa, Poisson's ratio ] = 0.33, and mass density  = 2700 kg/m 3 .The size of the square plate is L = 1 m, and its thickness is h = 0.001 m.For numerical studies, three types of boundary conditions are considered: (1) Example 1: a square plate with simplesimple-simple-simple (S-S-S-S) edge support, as shown in Figure 6; (2) Example 2: a square plate with free-clampedfree-clamped (F-C-F-C) edge support, as shown in Figure 7; and (3) Example 3: a square plate with free-free-free-clamped (F-F-F-C) edge support, as shown in Figure 8.An evaluation of the present spectral element model (denoted by "SEM") was conducted by comparing the natural frequencies obtained from the SEM with those obtained by using exact theory, the standard finite element method  (FEM), and the commercial finite element analysis package ANSYS [21].For the FEM results, a four-node 12-DOF conforming rectangular finite element model [22] is used.It is well known that exact solutions are available in analytical forms only for Levy-type plates.This is why we considered S-S-S-S square plate (Example 1) as the first example problem to evaluate the accuracy of the present spectral element model.The exact natural frequencies of S-S-S-S square plate are given by [24] where  is the dimension of the square plate and (, ) indicates the mode number.Table 1 compares the natural frequencies of the S-S-S-S square plate (Example 1) obtained by exact theory, FEM [22], and the present SEM.The SEM results obtained by using a single finite element are found to be identical to the exact solutions given by (112).The FEM results are found to converge to the exact solutions or the SEM results as the number of finite elements is increased to more than 100×100.
From Table 1, we can compare the range of natural frequencies that can be calculated with reasonable accuracy by the present SEM and FEM by comparison with exact solutions.The mesh size of the FEM when 50 × 50 elements are used (see the 4th where   are the natural frequencies provided by (112).
Figure 10 shows that the analytical solution by the MAM is very close to that by the SEM and that the FEM results approach the result by the SEM as the number of finite elements is increased to more than 50 × 50. Figure 11 compares the transient dynamic responses at  = 0.1 m on the -axis obtained by the MAM [24], the FEM [22], and the present SEM.Similar to Figure 10, the transient dynamic response by the MAM is very close to the SEM results, and the FEM results approach the SEM results as the number of finite elements is increased to more than 100×100.
The computation times (CPU times) for the SEM and FEM are compared in Table 4.To measure the CPU times, we used a standard desktop PC equipped with two sockets of Intel Xeon E5-2630 v3 processors and 320 GB of DDR4 RAM memory clocked at 2133 MHz.Table 4 shows that the CPU time for the present SEM is much smaller than that required to obtain the same accuracy level using the FEM.
Based on aforementioned observations, we conclude that the present spectral element model has the capability to provide extremely accurate natural frequencies and dynamic responses very efficiently.

Conclusion
We propose a new frequency domain spectral element model for finite rectangular plate elements with arbitrary boundary conditions.The new spectral element model is developed by using the boundary splitting method and the spectral super element method in combination.The high solution accuracy and computational efficiency of the proposed new spectral element model are validated by comparison with exact theory, the standard FEM, and the ANSYS commercial finite element analysis package.The conclusions drawn from our results are as follows: (1) The proposed new spectral element model can be applied to any finite rectangular plate element with arbitrary boundary conditions, whereas the spectral element model introduced in the authors' previous work [5] can be applied only to finite rectangular plate elements whose four corner nodes are fixed.Furthermore, the proposed new spectral element model can be assembled in the and -directions to represent a plate structure, whereas most existing spectral element models are valid only for plates with very specific boundary conditions such as Levy-type plates, infinite or semi-infinite plates, and plates in which assembly is allowed only in one direction.(2) Through numerical studies, we show that the proposed spectral element model provides highly accurate solutions by using a relatively small number of finite elements.
(3) The proposed spectral element model is found to be computationally efficient compared to the standard FEM and ANSYS, because the proposed spectral element model has nodal DOFs only on four edges of a rectangular finite element to significantly reduce the total number of DOFs used in the spectral element analysis.

Figure 1
Figure1: Boundary splitting method used in the previous study[5] to derive (, ) for a rectangular plate element subjected to arbitrary boundary conditions (e: active nodes; O: fixed nodes).

Figure 2 :
Figure 2: Boundary splitting method used in this study to derive (, ) for a rectangular plate element subjected to arbitrary boundary conditions (e: active nodes).
) is a one-by-four interpolation function matrix and w ()  () are the nodal line degree of freedom (DOF) functions defined by

Figure 3 :
Figure 3: Finite strip element representation of a rectangular finite plate element subjected to arbitrary boundary conditions at  = −  /2 and   /2 (e: nodes; grey circles: nodal values).

Figure 4 :
Figure 4: Spectral nodal degrees of freedom (DOFs) of a rectangular finite plate element (e: nodes).

Figure 5 :
Figure 5: Finite strip element representation of a rectangular finite plate element subjected to null boundary conditions at  = −  /2 and   /2 (e: nodes).

Figure 11 :
Figure 11: Comparison of the impulse-induced transient dynamic responses at  = 0.1 m on the -axis of a simply supported square plate obtained by the present SEM, the modal analysis method (MAM), and the FEM(n), where n is the number of finite elements used in the analysis.
(b), has arbitrary boundary conditions rather than fixed boundary conditions on the parallel edges at  = −  /2 and   /2, and its solution is represented by   (, ).Problem , shown in Figure 2(c), has fixed (null) boundary conditions on the parallel edges at  = −  /2 and   /2.However, the boundary conditions at  = −  /2 and   /2 in Problem  must be specified such that the sum of the boundary conditions at  = −  /2 and   /2 in Problem  and those in Problem  is identical to the boundary conditions at  = −  /2 and   /2 in the original problem.The solution of Problem  is represented by   (, ).Then, the solution (, ) of the original problem can be obtained by summing the solutions to Problem  and Problem  as follows:  (, ) =   (, ) +   (, ) .

)
2.3.Derivation of   (, ).We can find the solution   (, ) to Problem  by using a procedure similar to that used to obtain the solution for   (, ) in Problem .Problem  can be obtained from Problem  by rotating the coordinate system (, ) 90 ∘ clockwise.However, the differences for Problem  are as follows: (1) the fixed (null) boundary conditions are placed at  = − /2 and   /2; (2) the finite plate element is divided into   finite strip elements in the -direction, as shown in Figure5; and (3)   (, ) should be determined to satisfy the boundary conditions at  = −  /2 and   /2 in combination with the boundary values contributed by   (, ).

Table 3 :
[22]arison of the natural frequencies (Hz) of a square plate with free-free-free-clamped (F-F-F-C) edge support obtained by ANSYS[21], FEM[22], and SEM.  : total number of finite elements;   : total number of degrees of freedom.Figure 9: Lowest six mode shapes of a square plate with free-free-free-clamped (F-F-F-C) edge support obtained by SEM, Example 3. Note:

Table 4 :
Comparison of the computation times (CPU times) required to obtain the dynamic responses shown in Figure11.  : total number of finite elements used in the analysis. Note: