## THE NONLINEAR BENDING OF SIMPLY SUPPORTED ELASTIC PLATE

**Authors:**Mathieu G., Tyekolo D., Belay S.**Issue:**Vol 18, No 1 (2017)**Pages:**58-69**Section:**Innovative software engineering research**URL:**http://journals.rudn.ru/engineering-researches/article/view/16000**DOI:**http://dx.doi.org/10.22363/2312-8143-2017-18-1-58-69

#### Abstract

In this article, assumptions in Classical Plate Theory (CPT) are explained; followed by concepts involved in Finite Element discretization for elastic Plate bending in CPT. Computer implementation aspects and Numerical Results of CPT elements are also included for analyzing nonlinear bending of simply supported elastic plates.

IntroductionA plate is a three-dimensional structural element, which is characterized by two key properties. Firstly, its geometrically three-dimensional solid whose thickness is very small (thin when compared with other dimensions of the faces (length, width, diameter, etc.). Secondly, the static or dynamic loads carried by plates are predominantly perpendicular to the plate faces.By “thin,” it is meant that the plate’s transverse dimension, or thickness h, is small compared to the length and width dimensions. A mathematical expression of “thin” in the aforementioned paragraph is: t/h < 1 thickness, and L represents a representative length or width.Prominently a plate has two special geometric features; viz,Thinness: One of the plate dimensions, called its thickness, is much smaller than the other two.Figure 1. Plate and associated (x, y, z) coordinate systemFlatness: The mid surface of the plate, which is the locus of the points located halfway between the two plate surfaces, is a plane.In this article attention is focused on the Kirchhoff model for bending of thin (but not too thin) plates. The term “thin” is to be interpreted in the engineering sense and not in the mathematical sense. For example, h/Lc is typically 1/5 to 1/100 for most plate structures.Consider first a flat surface, called the plate reference surface or simply its mid surface or midplane (see Figure 2, b). We place the axes x and y on that surface to locate its points. The third axis, z is taken normal to the reference surface forming a right-handed Cartesian system. Axis x and y are placed in the mid plane, forming a right-handed Rectangular Cartesian Coordinate (RCC) system.If the plate is shown with a horizontal mid surface, as in Figure 2 below, we shall orientz upwards.Figure 2. Idealization of plate as two-dimensional mathematical problemNext, imagine material normals, also called material filaments, directed along the normal to the reference surface (that is, in the z direction) and extending h/2 above and h/2 below it. The magnitude h is the plate thickness. We will generally allow the thickness to be a function of x, y, that is h = h(x, y), although most plates used in practice are of uniform thickness because of fabrication considerations.The end points of these filaments describe two bounding surfaces, called the plate surfaces. The one in the +z direction is by convention called the top surface whereas the one in the -z direction is the bottom surface.Assumptions in Classical Theory of PlatesThe classical plate theory (CPT) is based on the Kirchhoff hypothesis. Three assumptions involved in this hypothesis are:A cross-section perpendicular to the middle surface prior to deformation remains plane and perpendicular to the deformed middle surface (Fig. 3).The transverse normals do not experience elongation (i.e. they are inextensible).The transverse normals rotate such that they remain perpendicular to the mid- surface after deformation.Figure 3. Deformation of the cross section in the xz plane according to the Kirchhoff assumptionsFinite Elements for elastic Plate bending in CPTHistorically the first model of thin plate bending was developed by Lagrange, Poisson and Kirchhoff. It is known as the Kirchhoff plate model. In the finite element literature Kirchhoff plate elements are often called C1 plate elements because that is the continuity order nominally required for the transverse displacement shape functions.The application of the finite element method (FEM) to the analysis of Kirchhoff plate bending demands the continuity in the first derivative of the expansion of the deflectionThe reader is referred to Zienckiewicz’s excellent book for details.Conforming vs Non-Conforming Elements. Displacement should be compatible between adjacent elements. There should not be any discontinuity or overlapping while deformed. The adjacent elements must deform without causing openings, overlaps or discontinuity between the elements.On each element displacements and test functions are interpolated using shape functions and the corresponding nodal values.N p N pk k ku3(x1,x2 ) = ∑ N k (x1,k =1x2 )u3 , v3(x1,x2 ) = ∑ Nk =1(x1,x2 )v3 .k kWhere, Nk is a shape functions and u3 , v3 the Nodal values.To obtain the FE equations the preceding interpolation equations are introduced intothe weak form.Similar to Euler-Bernoulli Beam the internal virtual work depends on the second order derivatives of the deflection, u3 and virtual deflection, υ.The problem domain is partitioned into a collection of pre-selected finite elements(either triangular or rectangular).Conforming Elements. Elements that satisfy all the three convergence requirements and compatibility condition are called Compatible or Conforming elements.Figure 4. Conforming triangular elementTriangular element: A conforming triangular element due to Clough and Tocher is an assemblage of three triangles as shown in fig. 4. In each sub-triangle, the transverse deflection is represented by the polynomial (i = 1, 2, 3):i 2 2 2 2 2w0 (x, y) = ai + biξ + ciη + diξη + eiξ+ fiη+ giξ+ hiξ η + kiξη+ liη.Where, (ξ, η) are the local coordinates, as shown in fig. 4 above.Rectangular element: Conforming rectangular element with w0, w0/x, w0/y and2w0/xy as the nodal variables (see fig. 5) was developed by Bogner. This element hasfewer degrees of freedom, which are less accurate in theory, but in practice the analyticalsolution of the problems frequently are not sufficiently high to achieve the accuracy of the triangular elements. This is good choice of element for thin plate analysis.Figure 5. Conforming rectangular elementNon-Conforming Elements. Elements that violate continuity conditions are known asNon-conforming elements.The Adini-Clough element is a nonconforming element. Despite this deficiency the element is known to give good results.Triangular element: The first successful nonconforming triangular plate-bending element was the original BCIZ, (Nonconforming element of Bazeley, Cheung, Irons and Zienkeiwicz), and it consists of three degrees of freedom (DOF) (w0, θx, θy) at the vertex nodes (see fig. 6). The element performs very well in bending as well as vibration problems (with a consistent mass matrix).Figure 6. Non-conforming triangular element with 3 DOF (w0, w0/x, w0/y) per nodeRectangular element: Non-conforming rectangular element has w0, θx and θy as the nodal variables (see fig. 7).Non-conformity: The element is C0 continuous, since the functions are cubic along an edge and we have four degrees of freedom to specify the function. The normal derivatives are not continuous across inter element boundaries in general.Figure 7. Non-conforming rectangular elementComputer implementation aspects and Numerical Results of CPT elementsThe conforming and non-conforming rectangular finite elements discussed in this article are implemented into a computer program using bilinear interpolation of (u0, υ0) and Hermite cubic interpolation of w0.Results of Linear Analysis. We consider the bending of rectangular plates with various edge conditions to evaluate the elements developed herein. The foundation modulus k is set to zero in all examples. The linear stiffness coefficients are evaluated using 4 × 4 Gauss rule.Example 1: Consider a simply supported (SS-1) rectangular plate under uniformly distributed load. The geometric boundary conditions of the computational domain (see the shaded quadrant in fig. 8) are:u0 = w0/x = 0 at x = 0; υ0 = w0/y = 0 at y = 0; υ0 = w0 = w0/y = 0 at x = a/2; u0 = w0 = w0/x = 0 at y = b/2; 2w0/xy = 0 at x = y = 0(For conforming element only).Figure 8. Boundary conditions for rectangular plates with biaxial symmetryTable 1 shows a comparison of non-dimensionalized finite element solutions with the analytical solutions (see Reddy [3]) of isotropic and orthotropic square plates under uniformly distributed transverse load q0. The stresses were evaluated at the center of the element. Hence, the locations of the maximum normal stresses are (a/8, b/8), (a/16, b/16), and (a/32, b/32) for uniform meshes 2 × 2, 4 × 4, and 8 × 8, respectively, while those of σxv are (3a/8, 3b/8), (7a/16, 7b/16), and (15a/32, 15b/32) for the three meshes. The analytical solutions were evaluated using m, n = 1,3…19. The exact maximum deflection occurs at x = y = 0, maximum stresses σxx and σyy occur at (0,0,h/2), and the maximum shear stress σxv occurs at (a/2, b/2, -h/2).Table 1Maximum transverse deflections and stresses* of simply supported plates under a uniformly distributed load q0 (linear analysis)VariableNon-conformingConforming2 × 24 × 48 × 82 × 24 × 48 × 8Analytical solutionIsotropic plate (v = 0,25)ω × 1024.85714.64254.58834.76194.59524.57394.5698End of Table 1VariableNon-conformingConforming2 × 24 × 48 × 82 × 24 × 48 × 8Analytical solutionσxx0.24050.26730.27400.22390.26370.27310.2762σxy0.17130.19640.20500.16880.19350.20400.2085Orthotropic plate (E1/E2 = 25, G12 = G13 = 0,5E2, v = 0,25)ω × 1020.70820.66350.65310.77100.66510.65220.6497σxx0.71480.77090.78280.55600.73880.77430.7866σyy0.02960.02530.02460.02780.02490.02450.0244σxy0.03370.04210.04440.03750.04160.04480.04633 2ω= ω0 E2h , σ = σh .q a4q a20 0Example 2: Here we consider a clamped square plate under uniformly distributed load. The boundary conditions are taken to be:u0 = w0/x = 0 at x = 0; υ0 = w0/y = 0 at y = 0; at x = a/2; u0 = υ0 = w0 = w0/x = w0/y = 0 at y = b/2; u0 = υ0 = w0 = w0/x = w0/y = 02w0/xy = 0 on clamped edges (For conforming element only)Table 2Maximum transverse deflections and stresses* of clamped (CCCC), isotropic and orthotropic, square plates (a = b) under a uniformly distributed load q0 (linear analysis)VariableNon-conformingConforming2 × 24 × 48 × 82 × 24 × 48 × 8Isotropic plate (v = 0,25)ω × 1021.57311.46531.43421.47781.43701.4249σxx0.09870.12380.13010.08610.11971.1288σxy0.04870.02220.00670.04890.02240.0068Orthotropic plate (E1/E2 = 25, G12 = G13 = 0,5E2, v = 0,25)ω × 1020.14340.13320.13140.14020.13300.1311σxx0.19620.24910.25980.15590.23580.2576σyy0.00850.00460.00420.00660.00470.0043σxy0.00760.00460.00190.00830.00480.00203 2ω= ω0 E2h , σ = σh .q a4q a20 0Table 2 contains the nondimensionalized deflections and stresses. The locations of the normal stresses reported for the three meshes are:2 × 2: (a/8, b/8); 4 × 4: (a/16, b/16); 8 × 8: (a/32, b/32)and shear stresses reported for the three meshes are:2 × 2: (3a/8, 3b/8); 4 × 4: (7a/16, 7b/16); 8 × 8: (15a/32, 15b/32).These stresses are not necessarily the maximum ones in the plate. For example for an 8 × 8 mesh, the maximum normal stress in the isotropic plate is found to be 0.2300 at (0.46875a, 0.03125b, -h/2) and the maximum shear stress u is 0.0226 at (0.28125a, 0.09375b, -h/2) for the non-conforming element. The conforming element yields slightly better solutions than the non-conforming element for deflections but not for the stresses, and both elements show good convergence.Results of Nonlinear Analysis. Here we investigate geometrically nonlinear response of plates using the conforming and non-conforming plate finite elements. The nonlinear terms are evaluated using reduced integration. Full integration (F) means 4 × 4 Gauss rule and reduced integration (R) means 1 × 1 Gauss rule.Example 3: Having established the credibility of the finite element for the linear analysis of CPT plates, we now employ the element in the nonlinear analyses. First, results are presented for single-layer isotropic square plates under uniform loading. The essential geometric boundary conditions used are (BC3 and BC5 in Table 1), simply supported (SS-3):SS-3: u = v = w = 0 on all edges.Clamped (CC-l):u = v = w = 0 on all edges,Figures 9 and 10 show the nondimensiona- lized center deflection, w = w/h, and non- dimensionalized center stress, σ = σa2/Eh2, as a function of the load parameter, P = q0a4/Eh4 for.ψy = 0 along edges parallel to x-axis,ψx = 0 along edges parallel to y-axis.simply supported (SS-3) square plate, under uniformly distributed load. Figure 11 shows similar results for clamped (CC-l) square plate under uniformly distributed loadThe results are compared with the Ritz solution of Way [5], double Fourier-series solution of Levy [6], the finite-difference solution of Wang [7], the Galerkin solution of Yamaki [8], and the displacement finite-element solution of Kawai and Yoshimura [9]. Finite-element solutions were computed for the five degrees of freedom (NDF = 5), and for three degrees of freedom (NDF = 3); in the latter case, the in-plane displacements were suppressed. Since suppressing the in-plane displacements stiffens the plate, the deflections are smaller and stresses are larger than those obtained by including the in-plane displacements. Solutions of the other investigators were read from the graphs presented in their papers. The present solutions are in good agreement with the results of other investigators.Figure 9. Load-deflection curves by various investigators for simply supported (SS-3) square plate under uniformly distributed loadFigure 10. Nondimensionalized center stress versus the load parameter for simply supported (SS-3) square plate under uniformly distributed load(a) Load-deflection curves (b) Center stress versus load parameterFigure 11. Nondimensionalized center deflection and stress versus the load parameter for clamped (CC-l) square plate under uniformly distributed loadTable 2 shows nondimensionalized center deflection w = w/h, center stressA 0 2 2B 2 2x = σx a / Eh , and edge stress = x= σx(a/2, 0)a /Eh , of aclamped (CC-2) squareplate under uniformly distributed load, qo. The boundary conditions are of type BC6,CC-2: u = v = w = ψx = ψy = 0 along all edges.The present solution is obtained using a 2 by 2 uniform mesh (in quarter plate) ofnine-node isoparametric elements (2Q9) with (R) and without (F) reduced integration. The present solution for center deflection and stresses agree very closely with the finite element solution of Pica et al. [10], and the analytical solution of Levy [11]. The edge stress, for some reason, does not agree with the other two results.Types of boundary conditions used in the present studyTable 1Notation (Type)Description of essential boundary conditionsSide 1Side 2BC1 (SS-I)v = w = ψy = 0u = w = ψx = 0BC2 (SS-2)u = w = ψy = 0v = w = ψx = 0BC3 (SS-3)u = v = w = 0u = v = w = 0BC4 (SS-4)u = v = w = ψy = 0u = v = w = ψx = 0BC5 (CC-i)u = v = w = ψx = 0u = v = w = ψy = 0BC6 (CC-2)u = v = w = ψx = ψy = 0u = v = w = ψx = ψy = 0BC7 (CC-3)u = w = ψx = 0v = w = ψy = 0BC8 (CC-4)w = ψx = 0w = ψy = 0Table 2A BNondimensionalized center deflection (w), center stress (x ), and edge stress (x )p == q0a4/Eh4ω = ω0/hA A 2 2B B 2 2PicaPicaPicapresent F(2Q9) Ret al.1 [70]Levy [71]present F(2Q9)2 Ret al. [70]Levy [71]present F(2Q9) Ret al. [70]Levy [71]17,790,19040,24550,23680,2372,2392,4592,63192,60,89040,5585,31635,5838,30,38810,47840,36990,4714,8395,1295,48165,21,6921,39411,21611,5263,40,58970,70450,69150,6957,7677,8348,32588,02,3732,67217,72618,0395,00,79090,91470,90290,91210,9710,4611,10311,12,9154,38924,96725,32134,90,98621,11891,10631,12114,3813,0913,85713,33,3166,60633,04533,5for clamped (CC-2) square plate under uniformly distributed load, (q0) ¶σx = σx a /Ehσx = σx a /Ehp == q0a4/Eh4ω = ω0/hA A 2 2B B 2 2PicaPicaPicapresent F(2Q9) Ret al.1 [70]Levy [71]present F(2Q9)2 Ret al. [70]Levy [71]present F(2Q9) Ret al. [70]Levy [71]189,01,17911,31891,30091,32317,9815,7516,49715,93,5689,32541,88542,4245,01,37811,51551,49281,52121,8918,4819,22519,23,66512,5951,71952,8318,01,56721,70201,67861,71425,8921,2121,99421,93,62716,2962,32563,9402,01,74691,87601,85551,90229,9123,8924,78025,13,46420,3073,40775,8End of Table 2σx = σx a /Ehσx = σx a /Eh4Q9 (4 by 4 mesh of 9-node elements).computed at the nearest Gauss points.Figure 12. Load-deflection curves for simply supported and clamped orthotropic square plates under uniform loadExample 4: This example is concerned with the bending of simply supported (SS-I) orthotropic square plate under a uniformly distributed transverse load q0. The geometry and material properties used are given below.a = b = 12 in., h = 0.138 in., E1 = 3 × 106 psi, E2 = 1.28 × 106 psi,G12 = G23 = G13 = 0.37 × 106 psi, υ12 = υ23 = υ13 = 0.32.A load increment of Δq0 0.2 psi and a uniform mesh of 4 × 4 in a quarter plate was used Figure 12 shows plots of the center deflection versus the intensity of the distributed load. The finite element results are in close agreement with the experimental results of Zaghloul and Kennedy [12].Example 5: In this last example of this section we consider the bending of a clamped square plate under uniformly distributed load, q0. The geometry and material parameters used are the same as those in Example 4. The boundary conditions of a clamped edge are taken to be:u0 = υ0 = w0 =∂w0 =∂z∂w0∂y= 0.Of course, for conforming element, one may also impose (2w0/xy = 0). A uniform mesh of 8 × 8 non-conforming elements in a quarter plate is used, and a load increments of:{Δq0} = {0.05, 0.05, 0.1, 0.2, 0.2 … 0.2}psi was used. The linear solution at q0 0.05 is found to be w0(0,0) = 0.00302 in. A plot of the center deflection versus the intensity of the distributed load for the clamped orthotropic plate is included in Figure 12 above.ConclusionA finite element formulation based on Kirchhoff’s plate model assumptions has been developed for the analysis of nonlinear bending of simply supported elastic plates.Different continuity conditions, according to displacement gradients, can be introduced into the formulation of geometrically non-linear non-conforming elastic plate elements to ensure convergence.We can also conclude with a note that the plate bending elements of the CPT discussed here are adequate for most engineering applications, which involve thin, and isotropic plate structures, which can be used to analyze both thin and thick plates.© Gil-oulbé Mathieu, Dau Tyekolo, Soresa Belay, 2017

### Gil-oulbé Mathieu

Peoples’ Friendship University of Russia
**Author for correspondence.**

Email: giloulbem@mail.ru

Miklukho-Maklaya str., 6, Moscow, Russia, 117198

Associate professor of the Department of Architecture and Construction

### Dau Tyekolo

Peoples’ Friendship University of Russia
Email: tiek.d@hotmail.com

Miklukho-Maklaya str., 6, Moscow, Russia, 117198

Assistant professor of the Department of Architecture and Construction

### Soresa Belay

Peoples’ Friendship University of Russia
Email: soresably@gmail.com

Miklukho-Maklaya str., 6, Moscow, Russia, 117198

Graduate student of the Department of Architecture and Construction

- Eduard Ventsel Theodor Krauthammer (2001) Thin Plates and Shells-theory, analysis and applications (The Pennsylvania State University, Pennsylvania).
- IT Kharagpur NPTEL (National Program on Technology Enhanced Learning) Web Course. #.Module 1.
- URL: http://www.colorado.edu/engineering/CAS/courses.d/AFEM.d/Home.html (21.12.2016).
- Reddy, J.N. A Penalty Plate-Bending Element for the Analysis of Laminated AnisotropicComposite Plates. Int. J. Numer. Meth. Engng., Vol. 15, pp. 1187-1206, 1980.
- Way, S. Uniformly Loaded, Clamped, Rectangular Plates with Large Deformation. Proc. 5th Int.Congr. Appl. Mech. (Cambridge, Mass., 1938), John Wiley, pp. 123-238.
- Levy, S. Bending of Rectangular Plates with Large Deflections. Report No. 737, NACA, 1942.
- Wang, C.T. Bending of Rectangular Plates with Large Deflections. Report No. 1462, NACA,1948.
- Yamaki, N. Influence of Large Amplitudes on Flexural Vibrations of Elastic Plates. ZAMM, Vol. 41, pp. 501-510, 1967.
- Kawai, T. and YOSHIMURA, N. Analysis of Large Deflection of Plates by the Finite Element Method. Int. J. Numer. Meth. Engng., Vol. 1, pp. 123-133, 1969.
- Pica, A., Wood, R.O. and Hinton, E. Finite Element Analysis of Geometrically Nonlinear Plate Behavior Using a Mindlin Formulation. Computers & Structures, Vol. 11, pp. 203-215, 1980.
- Levy, S. Square Plate with Clamped Edges Under Pressure Producing Large Deflections. NACA, Tech. Note 847, 1942.
- Zaghloul, S.A. and Kennedy, J.B. Nonlinear Analysis of Unsymmetrically Laminated Plates.J. Engng. Mech. Div., ASCE, Vol. 101 (EM3), pp. 169-185, 1975.