Static, free vibration, and buckling analysis of plates using strain-based Reissner–Mindlin elements

A quadrilateral and a triangular element based on the strain approach are developed for static, free vibration and buckling analyses of Reissner–Mindlin plates. The four-node triangular element SBTP4 has the three essential external degrees of freedom at each of the three corner nodes and at a mid-side node; whereas the quadrilateral element SBQP has the same degrees of freedom at each of the four corner nodes. Both elements use the same assumed strain functions which are in the linear variation where bending and transverse shear strains are independent and satisfy the compatibility equations. The use of the strain approach allows obtaining elements with higher-order terms for the displacements field. The formulated elements have been proposed to improve the strain-based rectangular plate element SBRP previously published. Several numerical examples demonstrate that the present elements are free of shear locking and provide high-accuracy results compared to the available published numerical and analytical solutions.


Introduction
Analyses of static, buckling and free vibration of plate structures play a large role in structural engineering applications. Considerable research works on analysis of plates are still being conducted (Mackerle 1997(Mackerle , 2002Leissa 1969Leissa , 1987Liew et al. 1995. Designers prefer low-order Reissner-Mindlin plate elements due to their simplicity and efficiency. However, for thin plates, these elements often suffer from the shear locking phenomenon when dealing with thin plates. To overcome shear locking, many research works have been undertaken where the use of the selective reduced integration was first intervened (Zienkiewicz et al. 1971;Hughes et al. 1978;Malkus and Hughes 1978). The formulation procedure used is to divide the strain energy into two parts, one of bending and the other of shear. Then, two different integration rules for these two parts are used. For low-order polynomial elements based on displacement model, such as the four-node classical bilinear element, an exact integration (two Gauss points in each direction) is taken for the bending strain energy; whereas a reduced integration (one Gauss point) is used for the shear strain energy. This selective integration can be provided with a more efficient element but often leads to numerical instability. Considerable investigations have been oriented to develop robust elements using different improved formulations and numerical techniques to avoid shear locking such as mixed formulation, enhanced assumed strain methods, assumed natural strain methods, discrete shear gap method and smoothed finite element method (Lee and Wong 1982;Ayad et al. 1998;Lovadina 1998;César de Sá and Natal Jorge 1999;César de Sá et al. 2002;Cardoso et al. 2008;MacNeal 1982;Dvorkin 1985, 1986;Zienkiewicz et al. 1990; Batoz and Katili 1992;Bletzinger et al. 2000;Nguyen-Xuan et al. 2008;Liu and Nguyen-Thoi 2010).
The strain approach has been employed as an alternative to formulate robust plate elements (Belarbi and Charif 1999;Belounar and Guenfoud 2005;Belounar and Guerraiche 2014;Guerraiche et al. 2018;Belounar et al. 2018) to increase the accuracy and stability of the numerical solutions as well as to eliminate shear locking phenomena. The use of the strain approach (Belarbi and Charif 1999;Belounar and Guenfoud 2005;Belounar and Guerraiche 2014;Guerraiche et al. 2018;Belounar et al. 2018;Djoudi and Bahai 2004a, b;Rebiai and Belounar 2013; has several advantages where it enables to obtain efficient elements with high-order polynomial terms for the displacement functions without the need of including internal nodes. The first developed strain-based Mindlin plate element SBRP (Belounar and Guenfoud 2005) has been adopted for the linear analysis of plates having only rectangular shapes. However, this element suffers from shear locking for very thin plates . Then, the formulation of a new three-node strain-based triangular Mindlin plate element SBTMP ) has been developed for static and free vibration of plate bending. The assumed curvatures and transverse shear strains for the SBRP element (Belounar and Guenfoud 2005) are coupled and contain quadratic terms.
The key idea used in this paper is to formulate new elements to overcome shear locking for very thin plates and to improve the accuracy for plates with regular and distorted shapes.
In this paper, a quadrilateral and a triangular strain-based plate element have been formulated for static, free vibration and buckling analyses of plates using Reissner-Mindlin theory. The opportunity is taken to explore the displacements field obtained from the strain-based quadrilateral plate element (SBQP) by applying it to a four-node triangular element strain-based triangular plate with four nodes (SBTP4) having the same degrees of freedom (W, β x , and β y ) at each of the three corner nodes and a mid-side node. In the process of formulation, these elements are based on linear variation for the five strain components where bending and transverse shear strains are independent and satisfying the compatibility equations. The numerical study shows that the SBQP and SBTP4 elements pass the patch test, are free of shear locking, and can be found numerically more efficient than the SBRP element (Belounar and Guenfoud 2005).

Derivation of the displacements field
For Reissner-Mindlin plate elements (Fig. 1), the strains in terms of the displacements are given as: In matrix form, it can be given as The five strains, bending (κ x , κ y and κ xy ) and transverse shear (γ xz and γ yz ), given in Eq. (1a) cannot be considered independent, for they are in terms of the displacements W, β x and β y and therefore, they must satisfy the compatibility equations (Belounar and Guenfoud 2005) given as: The field of displacements due to the three rigid body modes is obtained by having Eq. (1a) equal to zero and the following results are obtained: The proposed quadrilateral and triangular elements (SBQP and SBTP4) have three degrees of freedom (W, β x and β y ) at each of the four nodes. Therefore, the displacements field should contain twelve independent constants and having used three (α 1 , α 2 , α 3 ) for the representation of the rigid body modes, the remaining nine constants (α 4 , α 5 , …, α 12 ) are to be apportioned among the five assumed strains of the two elements.
The interpolation of the assumed strains field for the present elements (SBQP and SBTP4) is given as: Assumed bending (κ x , κ y and κ xy ) and transverse shear (γ xz and γ yz ) strains given in Eq. (4) of the proposed elements are independent and have only linear terms contrarily for the SBRP element (Belounar and Guenfoud 2005), where bending and transverse shear strains are coupled and quadratic terms are included in the assumed shear strain components.

Element matrices
The standard weak form for free vibration and buckling can, respectively, be expressed as: By substituting Eqs. (6), (7) and (13) into Eqs. (19) and (20), we obtain: Where the element stiffness, mass and geometrical stiff- , are, respectively, as: The stress-strain relationship is given by: [D] s are, respectively, rigidity, bending rigidity, shear rigidity matrices and [T] is the matrix containing the mass material density: For static analysis, we use For free vibration, we use For the buckling analysis, we use

Numerical validation
To validate the accuracy and efficiency of the formulated quadrilateral and triangular elements (SBQP and SBTP4), several numerical examples have been investigated for static, free vibration and buckling analysis of isotropic plates where the patch test of rigid body modes and the mechanic patch test are first carried out. The obtained results of the SBQP and SBTP4 elements are compared with other numerical and analytical solutions available in the literature.

Patch test of rigid body modes
To verify that both SBQP and SBTP4 elements pass the patch test of rigid body modes, the eigenvalues of the stiffness matrix for a single element are computed for various shapes and different aspect ratio. The only three zero eigenvalues obtained correspond to the three rigid displacement modes for a plate.

Mechanic patch test
In this patch test, a rectangular plate of (L = 2a = 40) length and (2b = 20) width simply supported at the three corner 1,2 and 3 (W 1 = W 2 = W 3 = 0) is considered where the plate is modeled by several elements as shown in Fig. 2 (Batoz and Dhatt 1990) for various side-thickness L/h ratio (10,100 and  (Table 1). Whereas for the case of M ns = 1 applied on all sides (Fig. 2), the obtained results at any points of the plate are M xy = 1 ( Table 1). The results given in Table 1 confirm that both SBQP and SBTP4 elements fulfill the mechanic patch test.

Square plates
A classical benchmark is first studied of square plate bending problem (Fig. 3) with different boundary conditions and various thickness-side (h/L) ratios subjected to a uniform load (q = 1), where the shear locking free test and convergence investigation of central deflection are considered in this study.
Shear locking free test is considered for a clamped square plate with several values of ratios (L/h = 10-1,000,000) using a mesh of 12 × 12. The central deflection results of the plate illustrated in Table 2 and Fig. 4, confirm that the new formulated elements (SBQP and SBTP4) are able to solve the shear locking problem when the plate thickness becomes gradually small. However, it is observed that the SBRP element (Belounar and Guenfoud 2005) exhibits from shear locking phenomena for (L/h > 100). Now, convergence tests of a square plate are investigated with three cases of boundary conditions [clamped, soft simply supported SS1 (W = 0), and hard simply supported SS2 (W = β s = 0)]. Various values of h/L ratios of 0.1, 0.01, and 0.001 are considered for thick, thin and very thin plates, respectively. The obtained results of the vertical displacement at the center of the plate are presented in Tables 3, 4 and 5 and Figs. 5, 6 and 7, which show that: • Faster convergence towards analytical solutions (Taylor and Auricchio 1993) is obtained using only a small number of elements for all cases of ratios (h/L = 0.1, 0.01, and 0.001) and boundary conditions.

Skew plates
To show the performance of the present elements to the sensitivity of mesh distortion, two examples of thin skew  plates subjected to a uniform load (q = 1) are considered which are known in the literature as severe tests and studied by many researchers (Razzaque 1973;Morley 1963). The first is concerned with Razzaque's skew plate (Razzaque 1973) (β = 60°) with simply supported on two sides and free on the other sides (Fig. 8). The results of the vertical displacement at the center of the plate using uniform meshes N = 2, 4, 8, 12 and 16 are given in Table 6 and Fig. 9 for (h/L = 0.001). The obtained results for both elements (SBQP and SBTP4) are in quite a good agreement with the reference solution given by Razzaque (1973). But it can be seen that the SBTP4 element is a little better than the SBQP and MITC4 (Nguyen-Xuan et al. 2008) elements.
The second example treated by Morley (β = 30°) (Morley 1963) is simply supported (W = 0) on all sides (Fig. 8). Using meshes of N = 4, 8, 16 and 32, the obtained vertical displacement at the center of the plate are presented in Table 7 and Fig. 10 for h/L = 0.01 and 0.001. It can be observed that for h/L = 0.01, the results of the SBTP4 and SBQP elements are in good agreement with the reference solution (Morley 1963); whereas, for h/L = 0.001, the SBTP4 element is more efficient than the SBQP and MITC4 (Chen and Cheung 2000) elements.

Free vibration of square plates
Convergence tests of the formulated quadrilateral and triangular elements are first undertaken for simply supported (W = β s = 0) and clamped plates with two thickness-side ratios (h/L = 0.005 and 0.1) (Fig. 3). The results of the first six non-dimensional frequencies (λ = (ω 2 ρL 4 h/D) 1/4 ) using the SBQP and SBTP4 elements with four regular meshes (N = 4, 8, 16 and 22) are presented in Tables 8, 9 Table 12 and the first four mode shapes of SSSF and CFCF plates are plotted in Figs. 13 and 14. For all cases of boundary condition, the following can be concluded: • The present results are very close to analytical solutions (Leissa 1969) and are more accurate than those of the MITC4, DSG3 and ES-DSG elements (Nguyen-Thoi et al. 2012). • The two elements (SBQP and SBTP4) have similar behavior, are shear locking free and their accuracy is insensitive to boundary conditions.

Free vibration of parallelogram plates
A cantilever parallelogram plate of skew angle = 60° with two h/L ratios (0.001 and 0.2) is studied (Fig. 15) using 22 × 22 mesh. The computed six non-dimensional frequencies (λ = ωL 2 /π 2 (ρh/D) 1/2 ) and the mode shapes are illustrated in Table 13 and Fig. 16, respectively. These results are compared with other numerical (DSG3, ES-DSG3, and MITC4) (Nguyen-Thoi et al. 2012) and analytical solutions (Karunasena et al. 1996). It can be seen that the SBQP and SBTP4 elements have a good accuracy compared to exact solutions (Karunasena et al. 1996)

Buckling of square plates subjected to uniaxial compression
Square plates subjected to uniaxial compression ( Fig. 17) with h/L of 0.01 is analyzed for both simply supported (SSSS) and clamped (CCCC). The buckling load factor is defined as K h = λ cr L 2 /(π 2 /D). The results of the buckling load factor for the SBQP and SBTP4 elements using 4 × 4, 8 × 8, 12 × 12, 16 × 16 and 20 × 20 meshes are presented in Table 14 and Fig. 18. For all cases of boundary condition, the two elements (SBQP and SBTP4) have similar results and converge to analytical solutions (Timoshenko and Gere 1970). In addition, these elements have excellent accuracy compared to other elements (DSG3 and ES-DSG3) (Nguyen-Xuan et al. 2010a, b).
The results of the buckling load factor (K h ) and the relative error using 20 × 20 mesh are presented in Table 15.
Numerical results of the SBQP and SBTP4 elements are in good agreement with analytical solutions (Timoshenko and Gere 1970) and other numerical solutions (Nguyen-Xuan et al. 2010a, b;Tham and Szeto 1990;Vrcelj and Bradford 2008;Liew and Chen 2004).

Buckling of square plates subjected to biaxial compression
Square plate subjected to biaxial compression (Fig. 19) with three essential boundary conditions (SSSS, CCCC, SCSC) is considered for h/L = 0.01 using a mesh of 16 × 16. The buckling load factor results (K h = λ cr L 2 /(π 2 /D)) of the proposed elements are presented in Table 16 with analytical (Timoshenko and Gere 1970) and other numerical solutions (Nguyen-Xuan et al. 2010a, b;Tham and Szeto 1990; Morley (1963) 0.408 0.408   (Abbassian et al. 1987) 4.4430 7.0250 7.0250 8.8860 9.93500 9.93500  (Abbassian et al. 1987) 4.3700 6.7400 6.7400 8.3500 9.2200 9.2200  (Abbassian et al. 1987) 5.9990 8.568 8.568 10.4070 11.4720 11.4980  (Abbassian et al. 1987) 5.7100 7.8800 7.8800 9.3300 10.1300 10.1800     plates. The four-node strain-based triangular element SBTP4 has the three engineering external degrees of freedom at each of the three corner nodes and one mid-edge point, while the quadrilateral element SBQP has the same engineering degrees of freedom at each of the four corner nodes. These developed elements passed successfully both patch and benchmark tests for plate bending problems. Numerical results show that the SBQP and SBTP4 elements are shear locking free, stable and superior to the original strain-based rectangular plate element (SBRP) (Belounar and Guenfoud 2005) which suffers from shear locking when the plate thickness becomes progressively very thin and has less rate of convergence to analytical L L Fig. 17 Square plate subjected to axial compression  Fig. 18 Convergence of uniaxial buckling load factor (K h /K exact ) of square plates with h/L = 0.01 solutions for thick and thin plates. The obtained results using both strain-based elements (SBQP and SBTP4) show that a rapid convergence to analytical solutions can be achieved with relatively coarse meshes compared with other robust elements based on different methods. In perspective, these elements can be superposed with membrane robust elements to construct shell elements for the analysis of complex shell structures.  Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.