Development of the Dynamic Stiffness Method for the Out-of- Plane Natural Vibration of an Orthotropic Plate

In this present paper, the dynamic stiffness method (DSM) was formulated to investigate the out-of-plane natural vibration of a thin orthotropic plate using the classical plate theory (CPT). Hamilton’s principle was implemented to derive the governing differential equation of motion for free vibration of the orthotropic plate for Levy-type boundary conditions. The Wittrick–Williams (W–W) algorithm was used as a solution technique to compute the natural frequencies of a thin orthotropic plate for different boundary conditions, aspect ratios, thickness ratios, and modulus ratios. The obtained results are compared with the results by the finite element method using commercial software (ANSYS and those available) in the published literature. The presented results by the dynamic stiffness method can be used as a benchmark solution to compare the natural frequencies of orthotropic plates.


Introduction
A rectangular orthotropic plate has many applications in designing different components in the engineering field, such as aerospace, mechanical, and civil. The orthotropic response of a given material is due to the existence of its constitutive relations. Various composite plates have been modeled analytically as orthotropic plates in recent years. So, one should have knowledge of the free vibration of such structures for efficient dynamic structural analyses. Starting from the earlier works on natural vibration of the plate by Rayleigh [1] and Ritz [2], the past few decades have witnessed different numerical and analytical methods, such as the Kantorovich method, the superposition method, the Rayleigh-Ritz method, and the iterative reduction method, for the investigation of natural vibrations of rectangular orthotropic plates [3][4][5][6][7][8][9][10][11][12][13][14][15][16][17][18][19][20]. The usually adopted finite element method (FEM) has also proved its popularity in the dynamic analysis of orthotropic rectangular plates [21][22][23]. However, the main drawback of traditional FEM and other approximate methods is the discretization technique of the given structure, which depends on the number of elements taken. The dynamic stiffness method (DSM) provides accurate results independent of the number of elements implemented in the investigation. The 1. The DSM was formulated to investigate the natural vibration response of thin orthotropic plates. 2. The W-W algorithm was applied to compute the natural frequency of the orthotropic plate. 3. The DSM results were compared with the published literature and the finite element method. 4. A new set of DSM results is reported for different aspect ratios, thickness ratios, and modulus ratios, which may be used as benchmark solutions for comparison. Figure 1 shows the coordinate system of a thin orthotropic plate where the plate length is , with , and thickness is ℎ. The Levy-type solution is applied where two opposite edges are simply supported along the y-axis (i.e., along the edges = 0, and = ), while the other two edges may be free ( ), simply supported ( ), or clamped ( ), as represented in Figure 2. The material orthotropic axes of the plate are parallel in the directions of x and y.

Description of Geometrical Property
Applying the CPT assumption, the transverse displacements of an arbitrary point of the plate are expressed by Equation (1) where the displacement components of the plate are represented by , , in , , directions, respectively, and = ∅ , = ∅ represents the rotational displacements of the and axes at the plate middle surface, respectively. represents the thickness (transverse) displacement in the z-direction.
In the case of orthotropic plates, the orientation of the constitutive material is such that orthotropic axes 1 and 2 are equal to axes x and y, respectively, and the material constants, are given by Equation (2).
where , are Young's moduli, along with the orthotropic directions 1 and 2, respectively; , are major and minor Poisson's ratios, and is the shear modulus [43].

Equations of Motion
The well-known Hamilton's principle is used to derive the governing differential equation for the natural vibration of the orthotropic plate based upon CPT as given by Equation (3).
where the mass density and thickness of the orthotropic plate are represented by and ℎ, respectively. , , and are flexural rigidities, which are given by Equation (4).

Boundary Conditions
The natural boundary conditions (BCs) are obtained by applying Hamilton's principle, given by Equation (5).
where , indicate the shear force and bending moments of the plate.

Formulation of Dynamic Stiffness (DS) Matrix with Levy Solution
The generalized differential Equation (3) is solved by applying force and displacement boundary conditions to develop the DS matrix. A Levy solution is implemented to solve Equation (3), and to satisfy the BCs given by Equation (5), it is sought in the given form as presented in Equation (6) [33].
where represents the unknown frequency. A generalized fourth-order ordinary differential equation is determined by substituting Equation (6) into Equation (3); it can be expressed as The developed Equation (7) produces the standard four roots; based on its nature, there are only two feasible solutions obtained and given by cases 1 and 2. Case 1.
In the above case, all roots are real ( , − , , − ) and can be expressed as The solution is given by Equation (8).
In the above case are two real roots and two imaginary roots ( , − , , − ), and can be expressed as The solution is given by Equation (9).
For case 1, the formulation of the DS matrix is explained below. A similar pattern is implemented for case 2, but is not explained here for brevity.
( , ) = ( )sin ( ) where = / / The BCs for the Levy-type plate are represented in Figure 3. The BCs for displacements are: The BCs for the forces are:  The displacement BCs in Equation (13) are substituting into Equations (8) and (10); the following equations can be obtained as: This expression can be rewritten in the matrix form and expressed by Equation (15). i.e., Similarly, the force BCs are applied; substituting Equation (14) into Equations (11) and (12), the following matrix relationship can be developed and expressed by Equation (17). i.e., By excluding the constant vector value of , the following relationship can be formed as where Thus, a square 4×4 symmetric DS matrix from Equation (20) is developed, including independent terms( , , , , , ). Therefore, the generated DS matrix of the single plate element can be expressed as The mathematical expressions of Equation (21) are explained in Appendix A.

Dynamic Stiffness (DS) Matrix Assembly Procedure with Boundary Conditions
The dynamic stiffness matrix given by Equation (21) was developed for a plate element. By obtaining the natural frequencies of a given orthotropic plate assembly, we considered four elements of the given plate geometry. Each element of the plate was connected through nodal lines. Since there were five nodal lines and the degree of freedom was two per element, a 10×10 global master stiffness matrix was formulated. The assembly procedure is similar to the finite element method, schematically shown in Figure 4. Boundary conditions can be applied the same way as we applied in the finite element method. The penalty method was applied as boundary conditions to suppress a particular degree of freedom. In this method, a considerable value of stiffness is added to the appropriate term on the leading diagonal of the dynamic stiffness matrix. The procedure for applying the boundary conditions is summarized as follows: 1. Displacement ( ) is penalized for simply supported (S) boundary conditions. 2. Displacement ( ) and rotation ( ) are penalized for clamped (C) boundary conditions. 3. No penalty is implemented for the free (F) boundary condition.
where represents the suppressed node.

Application of the Wittrick and Williams (W-W) Algorithm
One way of determining the natural frequencies is by using the zeros of the global dynamic stiffness matrix of the structure under study. However, this method has its limitations due to the transcendental behavior of the DS elements, which makes the plot of the frequency determinant tedious. Moreover, sometimes this may lead to missing the coincident frequencies. So, to avoid these difficulties, the Wittrick-Williams algorithm [33,50] was used, which ensured that no natural frequencies were missed. The procedure to follow the W-W algorithm can be represented in a flow chart, shown in Figure 5.

Results and Discussions
This section presents the non-dimensional fundamental frequencies of an orthotropic rectangular plate with Levy-type BCs, i.e., two simply supported opposite edges and the other two edges having arbitrary boundary conditions. To simplify the problem, a twoletter notation was used to describe the boundary conditions of the remaining edges, as shown in Figure 2. For example, SC signifies that one edge is simply supported (S) and the other is clamped (C). The results obtained by DSM were compared with published The bisection method is implemented to the bracket any natural frequency between its upper and lower bounds to any required

Stop
To avoid the calculation of ( ) * , a sufficiently fine mess is applied so that the clamped-clamped natural frequencies are not exceeded.
results and the FEM using ANSYS. A detailed study and discussion were conducted to analyze the effects of the boundary conditions and the variations of the modulus ratio, aspect ratio, and thickness ratio on the fundamental natural frequency of an orthotropic plate. The following material properties were used for analysis [51]: / is varied from 3 to 50, The non-dimensional fundamental natural frequency is calculated as [52]:

Comparative Study
A program in MATLAB was constructed to compute the fundamental natural frequencies of the plate. The natural frequencies obtained by DSM were compared against the results reported by Reddy and Phan [52] using CPT and with those reported by Thai and Kim [15] using the variable refined plate theory (RPT), as shown in Table 1; errors incurred by DSM when compared with Ref. [15] are reported in the brackets. The natural frequency changes for different boundary conditions (due to a change in stiffness) were noticed in Table 1. The maximum value of natural frequencies is reported for the CC boundary condition; the minimum value was obtained for the FF boundary condition. In the CC boundary condition, more constraints were introduced at the edges of the plate, which increased the stiffness of the plate, resulting in a higher natural frequency. On the other hand, in the FF boundary condition, no constraint was applied at the edges of the plate, which decreased the stiffness of the plate, and lowered the natural frequency, as shown in Table 1.  Table 1 that the maximum error is −1.7436% for the CC boundary condition at b/L = 0.5 and E1/E2 = 40, whereas the minimum error encountered is 0.00% at b/L = 0.5 and E1/E2 = 25 for the FF boundary condition. It can be seen from Table 1 that most of the errors lie within 2%. It is observed that the DSM results of a thin orthotropic plate in Table 1 are in excellent agreement with the published literature and the FEM results.

Parameter Studies
To study the effects of different boundary conditions, variations of thickness ratio, modulus ratio, and aspect ratio on the non-dimensional fundamental frequencies of thin orthotropic plates, parametric studies were conducted. The DSM results are also compared with the FEM obtained by ANSYS. For the modeling of the plate, shell element 181 was used in ANSYS. To test the accuracy of FEM, the non-dimensional fundamental frequencies were obtained by DSM based on CPT for a simply supported Levy-type orthotropic plate (b/L = 1, b/h = 20) at different modulus ratios, compared to the values obtained by FEM with different mesh sizes, as shown in Table 2. A very good convergence between the FEM and DSM results can be observed. The final choice is a mesh with 20 × 20 plate elements, which is a good compromise between the need for accuracy and limiting the computing time. In Tables 3-8, the DSM results are compared against the published results reported by Mukhtar [53] using the differential transform method (DTM) and Taylor collocation method (TCM), and with those reported by authors [15] for different Young's modulus ratios (E1/E2), aspect ratios (b/L), and thickness ratios (b/h).. It could be observed from Tables 3-8 that the natural frequencies, in general, decrease with an increase in Young's modulus, aspect ratio, and thickness ratio for the orthotropic plate. It is observed in Tables 3-7 that with boundary conditions changing from SC to SF, the natural frequencies increase, with the increase in the modulus ratio (E1/E2) from 3 to 50 because of an increase in the stiffness of the plate. Compared with Ref. [15], the relative errors of DSM are reported in brackets in Tables 3-8.    Table 6. The next highest percentage of errors is −1.7539 (for the SC plate at b/L = 1, b/h = 20 and E1/E2 = 20) in Table 8. Thus, for all cases reported in Tables 3-8, the error is less than 2%.
The effect of Young's modulus ratio (E1/E2) and the thickness ratio (a/h) on natural frequencies for all boundary conditions are presented in Figure 6 and Figure 7, respectively. The following observations are obtained from Figure 6 and Figure 7, respectively.
(a) The natural frequency increases with increases in the modulus ratio for a given value of the side-thickness ratio (b/L = 0.5), as shown in Figure 6 for all boundary conditions, except FF and FS boundary conditions. The natural frequency values are nearly constant for FF and FS boundary conditions. (b) The fundamental natural frequency values increase with increases in the side to thickness ratio (b/h) for a given Young's modulus ratio , as shown in Figure 7 for all boundary conditions except FF and FS boundaries.

Conclusions
In the present paper, the dynamic stiffness method, as a novel method, was implemented to analyze the out-of-plane free vibration of the rectangular orthotropic plate, where two opposite edges are simply supported. The classical plate theory was used to formulate the dynamic stiffness matrix for the rectangular orthotropic plate. The Wittrick-Williams algorithm was applied to solve the transcendental nature of the global dynamic stiffness matrix to extract the natural frequencies of the overall plate. The complete procedure, starting from developing the dynamic stiffness matrix and calculating natural frequencies, was implemented in a computer program using MATLAB. This enabled computation of any number of exact natural frequencies of the orthotropic plate for the Levytype boundary conditions. The computed natural frequencies were compared against published results obtained by the finite element method using commercial software AN-SYS. A new dataset on the natural frequencies for different aspect ratios, modulus ratios, and thickness ratios was computed and compared with the published literature.
The new results on natural frequency obtained by the DSM method can be implemented as a benchmark solution for future comparison purposes.