CC BY
17Magazine of Civil Engineering. 2019. 85(1). Pp. 107-122 Инженерно-строительный журнал. 2019. № 1(85). С. 107-122
Magazine of Civil Engineering ISSN
° 2071-0305
journal homepage: http://engstroy.spbstu.ru/
DOI: 10.18720/MCE.85.9
Calculation method of bending plates with assuming shear deformations
Yu. Ya. Tyukalov*,
Vyatka State University, Kirov, Russia * E-mail: [email protected]
Keywords: bending plates, shear deformations, stresses, finite elements.
Abstract. The problem of calculating bending plates by the finite element method with considering of shear deformations is considered. The bending plates are widely used as common structures of various objects of civil and industrial construction. The solution was obtained on the basis of the principles of the minimum of additional energy and possible displacements. For approximation of moment fields, piecewise constant functions are used. Shear forces can be approximated by constant or piecewise constant functions. The necessary relations for rectangular and triangular finite elements are obtained. It is shown that the proposed method can be used in combination with traditional finite elements for thin plates obtained by the finite element method in displacements. A comparison of the solutions, obtained by the proposed method, with other known solutions for bending plates, with assuming of shear deformations, is given. A numerical estimate of the accuracy and convergence of the proposed method, when crushing the finite element mesh, is given.
1. Introduction
The bending plates are widely used as common structures of various objects of civil and industrial construction. Often, in modern constructions thick and multi-layered slabs are used. In calculating, such thick slabs, we should consider, besides the bending deformations, the shear strains, which can significantly affect the values of the plate displacements. The classical theory of Kirchhoff plate bending is based on assumption of the direct normals assumption and, therefore, does not allow for the shear deformations. The finite elements, which are developed based on the Kirchhoff theory, are can used only for the calculation of thin plates [1-2]. Therefore, Timoshenko-Mindlin theory of bending plates [3-4] are widely application for calculating thick plates. According to this theory, angles of rotation of the normals and the vertical displacements are considered as independent variables. Such an approach lowers the maximum order of derivatives in the strains energy functional and makes it possible to use the first-order function-forms for approximating the displacements. Studies have shown that the direct use of the Timoshenko-Mindlin theory for constructing finite elements in displacements leads to the effect of «locking», or to the impossibility of using these finite elements to calculate thin plates.
To overcome the «locking» effect, various procedures are used, such as putting the assumption of direct normals at discrete points or applying high order shift theories [5-6]. The finite elements based on the putting the assumption of direct normals in the middle points of the finite element sides are widely used in program complexes. The new methods for considering for shear deformations, based on equations in displacements, are also offered in [7-8]. In [9-10], nonlinear solutions for rod systems with considering for shear deformations of cross sections are considered. The high order shear theories, when deformations along a cross section change by to a law other than linear, are used in constructing analytical solutions to bending problems of rectangular plates [11-12]. To construct finite elements, that considering shear deformations, the theory of the third order is successfully applied [13-14]. In [14], the quadrangular finite element is presented, that has seven degrees of freedom at each node: three displacements along the axes of coordinates, two shear angles and
Tyukalov, Yu.Ya. Calculation method of bending plates with assuming shear deformations. Magazine of Civil Engineering. 2019. 85(1). Pp. 107-122. DOI: 10.18720/MCE.85.9.
Тюкалов Ю.Я. Метод расчета изгибаемых плит с учетом деформаций сдвига // Инженерно-строительный журнал.
2019. № 1(85). С. 107-122. DOI: 10.18720/МСЕ.85.9
("ОЕ^^^И This open access article is licensed under CC BY 4.0 (lUtps://crcativcconmrons.org/1 iccnscs/bv/4.0/)
two angles of rotation of the normals. This approach allows us to more accurately consider shear deformations, when the properties of the material change in various directions.
The procedure for introduction shear deformations in existing finite elements designed for the calculation of thin plates is proposed in [15]. The deformations of the transverse shear can be accounting since the direct use of the equations of the three-dimensional theory of elasticity. This method and the Galerkin's method in a weak form is used to construct a quadrilateral finite element [16]. Also, the Galerkin's method is used to construct triangular and quadrangular finite elements according to the Timoshenko-Mindlin theory in [17].
Another way to construct finite elements with accounting shear deformations is the use of mixed and hybrid variational formulations [18-20]. This approach, on the one hand, simplifies the consideration of shear deformations due to the use of transverse forces and moments as unknowns, together with displacements. On the other hand, in order to ensure convergence of solution, it is necessary to agree on approximations of displacements and forces. Note the work [21-23], in which the solution of the problem of plate bending with consider the shear deformations, the modified Mindlin's theory is used. The modification of the Mindlin's theory consists in the introduction of an additional unknown parameter in the form of an angle of rotation in the plane of the plate.
Thus, construction the models with considering shear deformations, which are alternatives by the finite element method in displacements, is actual for the bending plates. The purpose of this work is to develop the method for calculating the plates with accounting shear deformations based on the functional of additional energy and the principle of possible displacements [24-27], as well as comparing the solutions obtained for plates with different support conditions with solutions of the other methods.
2. Methods
Solving the problems of plate bending with considering the shear deformations due to transverse forces, we will obtain based on the functional of additional energy for an isotropic plate (for simplicity, we assume that there are no specified displacements) [1]:
nc =
^+ Ml -2MMy + 2(1 + y)Ml)dQ +1^2^)^2 + Q2^
min. (1)
E is the modulus of elasticity of the material; t is the plate thickness; v is Poisson's ratio; k is coefficient, which considering the parabolic law of change of the tangential stresses across the thickness of the plate. The functional (1), also called the Castigliano's functional, is also considered in [2]. In [1] it was
shown that in the linear theory of elasticity the value nc for the equilibrium state is minimal.
We write the functional (1) in matrix form that is more convenient for solving by the finite element method:
nc = 2f [E ]-1 {M }dQ + 2\{Q}T [Esh ]-1 (Q}dQ ^ min.
2j v , L , v , -+ 2.
In expression (2) the following notation is entered:
(2)
(M } =
M M M
xy
(q}=
\qv
[E]-1 =
12
E • t3
1 -v
- V 1
0 0
0 0
2(1 + v)
[Eh ]- =
5E • t
10 01
. (3)
In the functional (2), the first member is associated with the bending deformations of the plate, the second - with shear deformations by transverse forces. Mx and Qx are the bending moment directed along
the Xaxis and the corresponding shear force; My and Qy are the bending moment directed along the Yaxis
y ^y
and the corresponding shear force; My is torque. The bending moments are positive if the lower fibers of the plate are stretched.
In accordance with the principle of minimum of additional energy, the functions of moments and shear forces must satisfy the corresponding differential equations of equilibrium and static boundary conditions. Since, in the general case, it is almost impossible to select such functions, we will operate as follows. Divide the plate into rectangular or triangular finite elements. On the region of the finite element, we approximate the moment and shear force fields by piecewise constant functions (Figure 1a). Below we show that transverse forces can also be approximated by the functions, which are constant over the finite element region (Figure 1b).
Figure 1. Approximation of moments and shear forces in the region of finite elements: a) piecewise constant moments and shear forces; b) constant transverse forces.
Then the functional (2) can be written in the following form:
Пc = 2{M}T [D]M}+ 2{Q}T [Dshmin. (4)
{M} is vector of unknown nodal moments for the whole system;
[d] is flexibility matrix for the whole system under bending;
[Dsh ] is flexibility matrix for the whole system under shear.
Then, using the principle of possible displacements, we are construct algebraic equilibrium equations of the nodes of the grid of finite elements. In this case, independently, possible displacements causing only a bending state and possible displacements causing only a shift are considered. In Figure 2, it is showing such possible displacements by the example of beam elements. Under the possible bend (Figure 2a), the cross sections rotate and remain perpendicular to the neutral axis, and the nodal moments do the work. Under shear (Figure 2b), the cross sections remain vertical and only transverse forces perform a work. For triangular and rectangular finite elements of a plate, possible states are like to those shown in Figure 2.
Figure 2. Possible displacements of node 1: a) bending state; b) shear state.
Under bending, the equilibrium equations are expressed through the nodal moments and can be represented in the following matrix form:
C }t M}+p=o
Static boundary conditions we write in the following form:
г e„
(5)
Mnj = MXJ cos2 ъ + MVJ sin2 ъ -2M г sin^ cos^ -Mnj = 0,
(- Mx,i + mJ
y,г
sin^i cos^i
+MyA
22 cos ъ - sin ъ
)-Mn
= 0,
г e
Under shear, the equilibrium equations are expressed only through the transverse forces and are represented in follow form:
{CSJ {Q }+ P = 0, ieS z.
(6)
Static boundary conditions for transverse forces are as follows:
Qm = Qx,, cos^ - Qni = 0, /eHr,
cpt is the angle between the tangent to the border at node i and the axis X; Mn,, Qn, are the values of setting moments and shear forces, which are normal to the boundary; Mnsi is the values of the setting moments, which are tangent to the boundary; Mxi,Myi,Myi are unknown nodal moments; Qxi, Qyi are unknown nodal shear forces; {M,}, {Q,} are vectors of unknown nodal moments and shear forces of all finite elements adjoining to node i; S z are set of nodes that have free displacement along the vertical axis Z; Sr
z
are set of nodes lying on the border; Pi is the generalized force corresponding to the potential of external loads
under single possible displacement of the node i along the Z axis; {Ci}, {Csh i} are vectors, which contain
coefficients at unknown nodal moments and at shear forces in the equilibrium equations of node i along the vertical axis Z. Algebraic equilibrium equations provide the equilibrium of moments and shear forces in a discrete sense.
Thus, we have obtained the problem of minimizing quadratic function of several variables (4) with constraints in the form of system of linear algebraic equations. Unknown parameters are nodal moments and shear forces. To solve this problem, we use the well-known Lagrange multipliers method for account the equilibrium equations and static boundary conditions. Then, we get the following expression of the extended functional:
Пc = 2 M }T [d]m}+1 Q}T D ]Q}+
X ^ ({Q }T {M}+Pi)+ X w,, i (С, i } Q}+Pi)+
ie— 2 ie—z
XX (mxJ cos2 срг + Myi sin2 срг - 2Mxy,i sin Pi cosp, -Mn.)+ (7a)
ie=r
XX ((Myi - mx. )sinp cos Pi + Mxyi (cos2 p. - sin2 P )-M^)+
ie=r
XX Qx,i cosP +Qyi sinP - Qn min-
ie=r
Static boundary conditions can also be considered using the penalty function method, then we get
Пc = 1M }T [D]{M}+2 Q}T D ]Q}+ X Wi ({Ci }T {Mi}+ Pi)+ X w^ ({Csh,i }T Q.}+ Pi)+
ie^z ie—z
Xa(Mx,i cos2 P, + My.i sin2 P, - 2 Mxy,i sin Pi cos p. - Mn. )2 + (7b)
'e=r
Xa((My,i - mx. )sin P.- cos P. + Mxy,i (cos2 p1 - sin2 P.)-MnSJ )2
ie=r
XaQx,i cos Pi +Qy ,i sin Pi - Qn )2 ^ min-
2
+
*X,i
— Г
wt is vertical displacement of the node i, associated with the bend of the plate; wshi is vertical
displacement of the node i, associated with the shear of the plate cross sections; a is penalty parameter (large number). Obviously, the total displacement will be equal to the sum of these two values.
The use of penalty functions to account for static boundary conditions eliminates the introduction of additional unknowns, as compared to the Lagrange multipliers method. The calculation of the derivatives of the penalty functions along the unknown nodal forces leads to the appearance of additional addends to the
elements of the flexible matrix [D] and to the elements of the load vectors \fm } and {fq }. If the static boundary conditions are zero, then the elements of the vectors {fm } and {fq } are equal to zero. With considering the static boundary conditions, the matrix of flexibility will be denoted with the index r - [d r ]. Then, the expression of the functional (7) can be represented in a more compact matrix form:
n - = 2 M }t [d r]M 1+ 2 Q}T D ]Q}+M i fm }+fci fq }+ „ 2 2 (8)
{w}T (F }- [L]M })+ wh }T (F1- [Lh ]Q1) ^ m'n.
{w}, {wh} are vectors of global displacements of nodes, associated with bending and shearing, respectively; {F} is the loads vector, whose elements are equal to the works of external forces on the single
Magazine of Civil Engineering, 85(1), 2019
vertical displacements of nodes; [l], [lh ] are matrices of equilibrium equations, the rows of which are formed from the vectors {C,} and Csh t}, respectively.
Equating to zero the derivatives of the functional nc along the vectors {M} and {w}, we obtain the system of equations, consisting of the equations of the compatibility of deformations and the equilibrium equations for bending:
D] -[L]' liWll i-(FM ¡1 ...
-[L] [0] JlMn-{FH <9>
If piecewise constant approximations are used for the moment fields, then the we easily get inverse matrix [dr]-1 analytically. Therefore, we can easily express the vector {M} from the first matrix equation:
{M }= [d r]-[L]r {w}. + [d f]-|fm }. (10)
Substituting expression (10) into the second matrix equation of (9), we obtain the system of linear algebraic equations for determining the vector {w}:
MM={F}+[L]M-1 {fm} [k]=[L][Dr]-1 [L]T. <11>
The expressions for the elements of the vector {F}, matrices [d] , [l] , the algorithm of their formation and examples of the calculation of bent plates corresponding to the Kirchhoff theory, are given in [27].
To determine the displacement vector {w h} and the vector of transverse forces {Q} associated with
shear deformations, we equate to zero the derivatives of the functional nc along the vectors {q} and {w h} . Then we get the following system of equations, like the system of equations (9):
"D] -Lh]TUfeH J-F
-[Lh ] [0]
K }J l-{F }
(12)
The solution of the system of equations (12) can be performed in the following sequence:
[Kh ] = [Lsh D ]-1 [Lh ], [Ksh ]K }= {F }+[Lh D ]-1 F } {Q}= [or ]-1 [Lh ]T {wsh }+D ]-1.{FQ }
(13)
Thus, solving the minimization problem of the functional (7) is reduced to solving two independent systems of linear algebraic equations (9) and (12). This allows, if necessary, to analyze the state of bending and state of shear of plates separately from each other. For example, one can solve the problem of plate bending, without considering the shear, by the finite element method in displacements, and the additional displacements, associated with the shear of cross sections, can be determined from the solution of the system of equations (12).
We will obtain the necessary expressions for the elements of the matrix [Dsh ], [Lh ] and vector {F},
when using the rectangular and triangular finite elements for the discretization of the subject area (Figure 3). To approximate the fields of transverse forces on a finite element region, we consider piecewise-constant (Figure 1a) and constant (Figure 1b) functions.
Consider the variant of piecewise constant approximations of transverse forces. We introduce the following notation: {Qt }= <j Qx,t 1 is vector of shear forces at node t in the local coordinate system X1OY1.
[Qy,t J
{Qt}= |Qx,t| is vector of shear forces at node t in the global coordinate system XOY.
a)
o
Xt X3 X2
Figure 3. Triangular and rectangular finite elements.
Figure 4. The division of triangular finite element into areas with constant transverse forces: a) arbitrary triangle; b) orthogonal triangle.
b)
If the shear forces are piecewise constant, the expression of the additional strain energy in the global coordinate system can be written as simple sum:
m is total number of nodes;
n R is the number of rectangular elements adjoined to node i;
nt T is the number of triangular elements adjoined to node i;
As is area of the s-th finite element;
AS is the area of the part s-th of the triangular finite element with constant transverse forces, adjoined to the node i (Figure 4a).
For rectangular finite element, the division of the finite element into regions with constant moments is uniquely - into four equal regions. For triangular element, each side must be divided equally, but there must also be the point, inside the element, at which the three areas connected to the nodes intersect. This point can be defined as the intersection point of perpendiculars drawn from the middle of the sides (Figure 4a). If the greatest angle of the triangle is more than 90 degrees, then such the point will lie outside the triangle. In this case, the triangle is divided into zones by lines passing through the midpoints of the sides. These lines will be parallel to the sides of the triangle. The proposed division of the triangular element into regions with constant moments, allows to obtain more accurate results, than when simply dividing into three equal
In Figure 4a point O is the center of the circle, described around triangle. OA, OB, OC are perpendiculars, which drawn from the midpoints of the corresponding sides of the triangle. Denote the lengths of the sides of the triangle as ll2, l23, l3l. From geometric scheme in Figure 4a we get:
(14)
(15)
parts - As.
(16)
- As, A3 = 1 4 2
corners of the triangle is greater than 90°, then, as well as for orthogonal triangle, we get:
Obviously, for orthogonal triangle (Figure 4b) we get - A1s = A2 = — As, A3 = — As. If one of the
1
1
AS = AS = — As, AS = — As. In addition, it is obvious, that it is better not to use too long triangles.
i 2 4 ' 3 2
We introduce the notation for the flexible matrix \Dshi] of "neighborhoods" of node i and for the global flexible matrix for the whole system [Dsh ], which consists of matrices for all nodes of the system:
riAh,i ] "
K ]= A, [Eh]-1, [Dsh]= '■■ . (17)
KJ
The matrix [Dsh ] has the block-diagonal form and is easily invertible analytically.
Di I-1
D I-1 =
D h1
(18)
If the shear forces are constant on the finite element region (Figure 1b), then
U * = 2 ZI= A, Q } K № }
where Ak is the area of the finite element; n is the total number of finite elements;
{Qk }=<j Qx k i is vector of unknown transverse forces of the finite element in the global coordinate
\Qy,k J
system. In this case, the flexible matrix for the whole system is determined by formulas like formulas (17)—(18), with the replacement of the index m by the index n.
We obtain the equilibrium equation for the possible displacement of a node of the rectangular finite element in the local coordinate system (Figure 3). For the rectangular finite element, we also introduce the local coordinate system, associated with its center, and the functions, which are expressed in normalized local coordinates in the following form:
2 x
a
2y b '
n = ^, i = 1,2,3,4.
(19)
The index i denotes the local node of the finite element; x, y are coordinates of the node along the
axes X1 and Y1, respectively; are local normalized coordinates of node i, taking values of 1 or -1.
Nodes are numbered counterclockwise, starting with the lower left node.
Possible displacements of the points of the finite element, for the possible shearing displacement of a node, are expressed in the following form:
=
(1 + + nn)
4
(20)
As a result of the possible displacement of the node, such shear deformations will arise in sections:
8v = d^shj = 6(1 + nn) Y = ^hj =nÎ1±M)
dx 2a ' 5y 2b '
(21)
Then, for the case of piecewise constant approximations, the work of internal transverse forces for finite element k, on possible displacements of node , can be expressed as follows:
U = ab r1 r1
M 4
J- J- (^YxzQx + 5YyzQy Wn = T,=1
V 8 v
1 +
\
2
Qx,f +
an
8
1+
2
Qy, i ) )
(22)
Substituting in (22) i equal to from 1 to 4, we obtain, for the considered finite element, expressions for the work of internal transverse forces for possible displacements of nodes, from the first to the fourth.
Unite the nodal transverse forces, expressed in the local coordinate system, for the finite element k into the vector Q k }.
Qk F = (q 1 qy,! q*,2 qy,2 qx,3 qy,3 Qx,4 Qy,4) (23)
Also, we introduce vector combining the values of the work of internal transverse forces for possible displacements of all nodes of the finite element:
^Uk }T = (sukz SUkz SUl SUl). Then, we can write the following expression:
SUk Wsh Q }
Using (22) we get the expressions of the elements of matrix \_Lk ].
(24)
(25)
[Li ]
16
- 3b - 3a - 3b -a -b -a -b - 3a
3b - a 3b - 3a b - 3a b -a
b a b 3a 3b 3a 3b a
- b 3a -b a - 3b a - 3b 3a
(26)
Nodal forces {Q/}, expressed in the local coordinate system, and Qf}, expressed in the global coordinate system, are connected by the matrix of direction cosines
[i 1=
cos a sin a - sin a cos a
(27)
[Sf ]-
(28)
a - the angle between the Y\ axis and the Yaxis (Figure 3). Using (27) we obtain the matrix of direction cosines for the finite element
[l ]
[l ]
_ [l] The work of the internal forces (25) can be represented as follows:
ssuk }=[l]h ]qk } [lh ]=[lksh ] (29)
The matrix [l^], conditionally, can be called as local «equilibrium» matrix of finite element for shearing in the global coordinate system. From matrices for finite elements [l^ ], in accordance with the numbering of the nodes and elements, the global matrix [Lsh] for the whole system is formed (see (12)).
Consider the case, when the transverse forces are approximated, in finite element region, by constant functions (Figure 1b). Then, the work of internal transverse forces for finite element k on the possible displacement of node i can be expressed as follows:
SUl =
ab r1 r1
4
J1J1 (Sr*Q +S/yZQy Wn = b-Y QXlt
+
+ 2 Qy,f •
(30)
1
Qxk, Qyk are transverse forces for finite element, expressed in the local coordinate system, are
combined into the vector of unknowns for finite element {q k }T = Qx k Qyk ). The matrix \_Lk ], in this case, will have the following form:
- a - b
[Lkh ] = 2
(31)
a - b a b - a b
The matrix of directional cosines will coincide with the matrix [/ ] (see (27)).
The potential of the external concentrated and uniformly distributed loads, for possible displacements of the node i along the global coordinate axis, is determined by (32).
5Vl = P +4qkab = R.. (32)
Pi is force concentrated in node; qk is uniformly distributed load.
The generalized forces Ri, in accordance with the numbering of nodes, are placed in the vector {F} (see (12)).
We obtain the equilibrium equations for triangular finite elements. The equilibrium equation for possible displacement of node of the triangular finite element can be obtained directly in the global coordinate system. For this, the possible displacements of the points of the finite element k on shearing will express using the triangular coordinates:
w, i (x, y) = Ti, i = 1 , 2 , 3 T = (a + bx+cty). (33)
a = x+iy+2- x+2y+1, b = y+1- y+2, c = x+2- x+1. (34)
Ak is area of the triangular element;
xi, yi are coordinates of the node i (Figure 3). Triangular coordinates are natural coordinates of triangular area.
The function Tt takes the value 1 at node i and the value zero at other two nodes. With the possible
displacement of node, constant shear deformations arise in cross sections:
= = d(dw^)
dx 2Ak' yz dy 2Ak
8y = ' sKl' = , 8y = sKl' ^^V. (35)
' xz ^ ~ .k yz ^ 0 4 k v '
Consider the case of piecewise constant approximations of transverse forces in finite element region. The vector of nodal forces for triangular finite element, expressed in the global coordinate system, will have the following form:
{qk }t =(qx,1 qy,1 qx,2 qy,2 qx,3 qy,3) (36)
The work of the internal transverse forces of the k-th finite element on the possible displacement of the node i is expressed as an integral:
5Uk = MSv Q +Sy Q )dA = XV3 Q Ak +Y3 Q A. (36)
i,z JV fxz^-x lyz^-y^ 2Ak ¿—'J x,} } 2A j = y,1 j Ak
Ak is area of part of k-th triangular finite element is adjoined to the node j (Figure 4a) and is determined
by (16). We introduce the vector, that combines the values of the works of transverse forces on possible displacements of finite element nodes:
Then we get
{suk }T =(sukz sul su 3 z).
№ }=[Lkh ]Qk ¡i
[Lh ]
2 Ak
bi Ak b2 Ak
ci Ak C2 A1
b3 Ai c^ Ai
b Ak
b2 Ak
h Ak
Ci Ak C2 A2 C3 A2
bi A3
b2 Ak
b3A3
Ci Ak
C2 A3 C3 A3
(37)
(38)
(39)
In the case of using constant approximations of transverse forces, the vector of unknowns for the triangular element will have the following form:
{Qk}T =(Qxk Qyk)
(40)
Calculating the integral (36), we get:
[Lkh ]
(41)
From matrices [Lksh] for triangular finite elements, in accordance with the numbering of the nodes and the elements, the global matrix [Lh] for the whole system is formed (see (12)).
The potential of the external concentrated and uniformly distributed loads for possible displacements of the node i, along the global axis of coordinates, is calculated by (42).
5Vl = p + ^ qkAk = R.
(42)
The global «equilibrium» matrix [Lh ] for the whole system will have tape structure of nonzero elements. The number of rows of the matrix [Lh] is equal to of the number of loose nodes of the system. Numbering of unknown are assigned according to the numbering of nodes and finite elements. Therefore, the width of the tape of nonzero elements of the matrix [Lh ] will be determined by the maximum difference in the numbers of nodes of all finite elements adjoined to the node. After calculating the width of the tape for each of the rows, its maximum value lmax is determined. Then, for the elements of the matrix [Lh ], you can use a rectangular
array consisting of columns lmax and m rows. The tape structure of nonzero elements is used in constructing the matrix multiplication algorithm, when we calculate the elements of matrix [Kh ]. Note, using constant approximations of transverse forces in finite element region, the value of lmax and, thus, the width of the tape of matrix will be significantly smaller.
3. Results and discussion
To assess the accuracy of the proposed method, rectangular plates were calculated with different conditions for supporting the sides (Figure 5) on the action of uniformly distributed load.
In Figure 5, the dashed line and the letter S denote the hinged supported side along the X axis, the skew hatch and the letter C denote the clamped side, the letter F denote the free side. Table 1 presents the results of calculations of the SS plate, given in [8] for various theories, and the results obtained by the proposed method - SFEM. For approximation of shear forces, piecewise constant functions were used (Figure 1a). For crushing the plates, square finite elements with a side size of 0.05 m. were used. And the size b = 6 m. Poisson's ratio is v = 0.3 . The results of calculations in Table 1 are presented in the dimensionless form:
w =
Etk
4
qa
w
a b 2 i 2
Tyz =
qa
a
T = —T
xz x
qa
i o,b
I 2
T yz =-
3Qy
2t
t =
3Qxz 2t
(43)
i
b
C
b
C
2
2
b
C
3
3
t
t
0
2
Figure 5. Conditions of supports and the sizes of sides of the Levi's plates. Table 1. Displacements and stresses for the SS plate under the action of uniformly distributed
load.
b / a
a /1
Methods w
ABACUS 0.0536
CRT 0.0444 0.4909 0.4909
FSDT 0.0536 0.4909 0.4909
TSDT 0.0535 0.3703 0.3703
Analytical [8] 0.0536 0.4909 0.4909
SFEM 0.0536 0.5064 0.5064
ABACUS 0.0467
CRT 0.0444 0.4909 0.4909
FSDT 0.0467 0.4909 0.4909
TSDT 0.0467 0.4543 0.4543
Analytical [8] 0.0467 0.4909 0.4909
SFEM 0.0467 0.5064 0.5064
ABACUS 0.0444
CRT 0.0444 0.4909 0.4909
FSDT 0.0444 0.4909 0.4909
TSDT 0.0444 0.4909 0.4905
Analytical [8] 0.0444 0.4909 0.4909
SFEM 0.0444 0.5064 0.5064
ABACUS 0.1248
CRT 0.1106 0.5240 0.6813
FSDT 0.1248 0.5240 0.6813
TSDT 0.1248 0.4569 0.5615
Analytical [8] 0.1248 0.5240 0.6813
SFEM 0.1249 0.5541 0.6969
ABACUS 0.1141
CRT 0.1106 0.5240 0.6813
FSDT 0.1142 0.5240 0.6813
TSDT 0.1142 0.5051 0.6448
Analytical [8] 0.1142 0.5240 0.6813
SFEM 0.1146 0.5541 0.6969
ABACUS 0.1106
CRT 0.1106 0.5240 0.6813
FSDT 0.1106 0.5240 0.6813
TSDT 0.1106 0.5238 0.6809
Analytical [8] 0.1106 0.5240 0.6813
SFEM 0.1107 0.5541 0.6969
10
100
10
100
1
5
2
5
Table 1 shows the results of calculations performed: by the ABACUS program; according to the classical Kirchhoff plate theory (CRT); according to the theory of plate bending with first order shear theory (FSDT); according to the theory of plate bend with third order shear (TSDT); by an analytical method, based on the third-order shear theory [8].
Comparison of calculation results shows, that the proposed calculation method (SFEM) in stresses, allows us to obtain solutions of the same accuracy for both thin and thick plates. This indicates the absence of a «locking» effect, which does not allow to obtain correct results for thin plates, when we use some types of finite elements. The values of displacements, obtained by the proposed method in stresses, practically coincide with the values obtained by other methods, that considering the shear deformations.
It should also be noted that the values of tangential stresses (in dimensionless form), obtained by the proposed method, do not depend on the thickness of the plate. For square plates, the maximum tangential stresses, obtained by the proposed method, are larger than the corresponding values, obtained by other methods, by about 3 %. For rectangular plates: the maximum tangential stresses, obtained along the short side of the slab, are 5.7 % larger, than the values obtained by other methods, and values along the long side are 2.2 % larger, than the values, obtained by other methods. Accounting for shear deformations is most important, when calculating flexural plates on stability, and, when determining the frequencies of free vibrations of plates. The values of the critical forces and frequencies of free oscillations will be affected by a decrease in the rigidity of the plates, due to additional shear deformations, therefore the accuracy and correctness of considering shear deformations is important. The proposed method for considering shear deformations is based on the fundamental principles of the minimum of additional energy and possible displacements. No additional techniques, such as the satisfaction of the Kirchhoff hypothesis at individual points (DKT -elements), or double approximation of displacements (MITS - elements), are not used in the proposed method.
Table 2. Displacement the center of the plate w = w x 100D /(qa 4 ) for different conditions of supporting of sides, under the action of uniformly distributed load (b = 2 a = 6 m ).
a /1
Bounders conditions
Methods CC CS SS CF SF FF
ABACUS 1.0000 1.0703 1.1429 1.2089 1.2842 1.4280
FCDT 1.0000 1.0704 1.1430 1.2090 1.2844 1.4283
Analytical [8] 0.9357 1.0373 1.1430 1.1757 1.2849 1.4293
SFEM 0.9792 1.0617 1.1474 1.1992 1.2876 1.4325
ABACUS 0.8850 0.9637 1.0453 1.0980 1.1827 1.3225
FCDT 0.8850 0.9637 1.0454 1.0981 1.1829 1.3228
Analytical [8] 0.8673 0.9546 1.0454 1.0893 1.1834 1.3239
SFEM 0.8813 0.9637 1.0495 1.0964 1.1848 1.3248
ABACUS 0.8511 0.9329 1.0180 1.0663 1.1545 1.2935
FCDT 0.8511 0.9330 1.0181 1.0664 1.1547 1.2938
Analytical [8] 0.8481 0.9314 1.0181 1.0651 1.1550 1.2944
SFEM 0.8539 0.9363 1.0221 1.0676 1.1561 1.2947
ABACUS 0.8445 0.9270 1.0128 1.0604 1.1494 1.2884
FCDT 0.8445 0.9270 1.0129 1.0605 1.1496 1.2887
Analytical [8] 0.8445 0.9270 1.0129 1.0605 1.1496 1.2887
SFEM 0.8485 0.9311 1.0169 1.0613 1.1496 1.2877
10
25
1000
5
Table 2 presents the results of determining the displacement of the center of a rectangular slab, for different variants of supporting the sides, and, for different ratios of the slab thickness to the size of short side of the slab. For crushing the slab, square finite elements with a side size of 0.05m were also used. For ratio a /1 = 5, when the effect of shear deformations is greatest, the displacements, obtained by the proposed method (SFEM), differ from the displacements, obtained by the ABACUS program and the FCDT method by, no more than 2 %. At the same time, for some boundary conditions, the obtained values of displacements, are larger, and for others, less, than the compared ones. For thinner plates, the values of displacements for all the methods listed in Table 2, differ a little.
Table 3. Displacement of the center of the hinge plate under the action of uniformly distributed load q = 10 kN / m2.
b/a na nb Constant Q Piecewise-constant Q
Rectangular elements Triangle elements Rectangular elements Triangle elements
1 5 5 0.28875 0.28610 0.28842 0.28534
10 10 0.28214 0.28184 0.28205 0.28155
20 20 0.28049 0.28050 0.28047 0.28042
30 30 0.28020 0.28021 0.28019 0.28017
2 5 10 0.69531 0.69810 0.69515 0.69724
10 20 0.68753 0.68869 0.68749 0.68840
20 40 0.68559 0.68539 0.68558 0.68531
30 60 0.68523 0.68544 0.68523 0.68536
Table 3, to assess the convergence of the proposed solution, for different numbers of finite elements along the sides (na and nb) presents the results for plate 0.6 m thick. The modulus of elasticity
E = 10000 kN /m2 and Poisson's ratio v = 0.3 . Piecewise-constant and constant approximations of shear forces, as well as rectangular and triangular finite elements were considered. The grid of triangular finite elements was formed by dividing rectangular finite elements into two triangular elements. For the square hinged plate (b / a = 1) in [28], an analytical solution was obtained for displacement of the center of plate with considering for shear deformations:
4
a q
f f ^ 1 + 4.6|
v a,
w = 0.00406-^-^-= 0.2785 m. (44)
Et3
Comparison with the values from Table 3 shows, that the values, getting by the proposed method, differs from the analytical value by less than 1 % for fine grid, and by about 4 %, for the coarsest grid. Note, that the results differ very slightly, both when using rectangular and triangular finite elements. The influence of the choice of the type of approximation of transverse forces is also insignificant. Thus, for solving using the proposed method, one can use constant approximations of transverse forces, which are more convenient for calculating branched systems.
4. Conclusion
1. The method for considering shear deformations, when bending plates are calculating by the finite element method, is proposed. The method is based on the fundamental principles of the minimum of additional energy and possible displacements and is applicable for the calculation of both thin and thick plates.
2. For approximation of shear forces, piecewise constant and constant approximations in finite element region can be used. To calculate branched systems, one can use constant approximations of shear forces without loss of accuracy.
3. Displacements from shear deformations are determined independently of bending-related displacements, therefore, the proposed method can be used in combination with traditional finite elements for thin plates, which was obtained by the finite element method in displacements.
4. Comparison of the solutions, obtained by the proposed method, with other known solutions for bending plates, considering shear cross sections, shows its good accuracy and convergence, when crushing the finite element mesh.
•12(1 -v2 )
References
1. Gallager, R. Metod konechnykh elementov. Osnovy [Finite element method. The basics]. Moskow: Mir, 1984. 428 p. (rus)
2. Slivker, V.I. Stroitelnaya mekhanika. Variatsionnyye osnovy [Structural mechanics. Variational bases]. Moskow: Izdatelstvo Assotsiatsii stroitelnykh vuzov, 2005. 736 p. (rus)
3. Timoshenko, S.P., Voynovskiy-Kriger, S. Plastiny i obolochki [Plates and shells]. Moskow: Nauka, 1966. 636 p. (rus)
4. Fallah, N. On the use of shape functions in the cell centered finite volume formulation for plate bending analysis based on Mindlin— Reissner plate theory. Computers & Structures. 2006. Vol. 84. Pp. 1664-1672.
5. Katili, I., Batoz, J.L., Maknun, I.J., Hamdouni, A., Millet, O. The development of DKMQ plate bending element for thick to thin shell analysis based on the Naghdi/Reissner/Mindlin shell theory. Finite Elements in Analysis and Design. 2015. Vol. 100. Pp. 12-27.
6. Khezri, M., Gharib, M., Rasmussen, K.J.R. A unified approach to meshless analysis of thin to moderately thick plates based on a shear-locking-free Mindlin theory formulation. Thin-Walled Structures. 2018. Pp. 161-179.
7. Yaghoubshahi, M., Alinia, M.M. Developing an element free method for higher order shear deformation analysis of plates. Thin-Walled Structures. 2015, Vol. 94. Pp. 225-233.
8. Thai, H.-T., Nguyen, T.-K., Vo, T.P., Ngo, T. A new simple shear deformation plate theory. Composite Structures. 2017. Vol. 171. Pp. 277-285.
9. Lalin, V.V., Zdanchuk, Ye.V., Kushova, D.A., Rozin, L.A. Variatsionnyye postanovki nelineynykh zadach s nezavisimymi vrashchatelnymi stepenyami svobody [Variational statements of nonlinear problems with independent rotational degrees of freedom]. Magazine of Civil Engineering. 2015. No. 4(56). Pp. 54-65. (rus)
10. Lalin, V.V., Belyayev, M.O. Izgib geometricheski nelineynogo konsolnogo sterzhnya. Resheniye po teoriyam Kirkhgofa i Kossera -Timoshenko [Bending geometrically nonlinear cantilever rod. Decision on theories of Kirchhoff and Cossera - Tymoshenko]. Inzhenerno-stroitelnyy zhurnal. 2015. No. 1(53). Pp. 39-55. (rus)
11. Park, M., Choi, D.-H. A two-variable first-order shear deformation theory considering in-plane rotation for bending, buckling and free vibration analysis of isotropic plates. Applied Mathematical Modelling. 2018. Vol. 61. Pp. 49-71.
12. Kumara, R., Lala, A., Singhb, B.N., Singhc, J. New transverse shear deformation theory for bending analysis of FGM plate under patch load. Composite Structures. 2019. Vol. 208. Pp. 91-100.
13. Thai, H.-T., Choi, D.-H. Finite element formulation of various four unknown shear deformation theories for functionally graded plates. Finite Elements in Analysis and Design. 2013. Vol. 75. Pp. 50-61.
14. Doa, T.V., Nguyenb, D.K., Ducc, N.D., Doanc, D.H., Buie, T.Q. Analysis of bi-directional functionally graded plates by FEM and a new third-order shear deformation plate theory. Thin-Walled Structures. 2017. Vol. 119. Pp. 687-699.
15. Mulcahy, N.L., McGuckin, D.G. The addition of transverse shear flexibility to triangular thin plate elements. Finite Elements in Analysis and Design. 2012. Vol. 52. Pp. 23-30.
16. Karttunen, A.T., Hertzen, R., Reddy, J.N., Romanoff, J. Shear deformable plate elements based on exact elasticity solution. Computers & Structures. 2018. Vol. 200. Pp. 21-31.
17. Duan, H., Ma, J. Continuous finite element methods for Reissner-Mindlin plate problem. Acta Mathematica Scientia. 2018. Vol. 38. Pp. 450-470.
18. Zenkour, AM. Exact mixed-classical solutions for the bending analysis of shear deformable rectangular plates. App. Math. Model. 2003. Vol. 27. Pp. 515-534.
19. Erkmen, R.E. Shear deformable hybrid finite-element formulation for buckling analysis of thin-walled members. Finite Elements in Analysis and Design. 2014. Vol. 82. Pp. 32-45.
20. Ozutok, A., Madenci, E. Static analysis of laminated composite beams based on higher-order shear deformation theory by using mixed-type finite element method. International Journal of Mechanical Sciences. 2017. Vol. 130. Pp. 234-243
21. Senjanovic, I., Vladimir, N., Hadzic, N. Modified Mindlin plate theory and shear locking free finite element formulation. Mechanics Research Communications. 2014. Vol. 55. Pp. 95-104.
22. Sukhoterin, M.V., Baryshnikov, S.O., Potekhina, Ye.V. O raschetakh plastin po sdvigovym teoriyam [About plate calculations based on shear theories]. Vestnik gosudarstvennogo universiteta morskogo i rechnogo flota im. admirala S.O. Makarova. 2015. No. 2 (30). Pp. 81-89. (rus)
23. Sukhoterin M.V., Baryshnikov, S.O., Knysh, T.P. Stress-strain state of clamped rectangular Reissner plates. Magazine of Civil Engineering. 2017. 76(8). Pp. 225-240. DOI: 10.18720/MCE.76.20
24. Tyukalov, Yu.Ya. Finite element models in stresses for plane elasticity problems. Magazine of Civil Engineering. 2018. 77(1). Pp. 23-37. DOI: 10.18720/MCE.77.3
25. Tyukalov, Yu.Ya. The functional of additional energy for stability analysis of spatial rod systems. Magazine of Civil Engineering. 2017. 70(2). Pp. 18-32. DOI: 10.18720/MCE.70.3
26. Tyukalov, Yu.Ya. Stress finite element models for determining the frequencies of free oscillations. Magazine of Civil Engineering. 2016. 67(7). Pp. 39-54. DOI: 10.5862/MCE.67.5
27. Tyukalov, Yu.Ya. Finite element models in stresses for bending plates. Magazine of Civil Engineering. 2018. 82(6). Pp. 170-190. DOI: 10.18720/MCE.82.16
28. Donnell, L.G. Balki, plastiny i obolochki [Beams, plates and shells]. Moskow: Nauka, 1982. 568 p. (rus)
Contacts:
Yury Tyukalov, +7(912)8218977; [email protected]
© Tyukalov, Yu.Ya., 2019
Инженерно-строительный журнал
сайт журнала: http://engstroy.spbstu.ru/
ISSN
2071-0305
DOI: 10.18720/MCE.85.9
Метод расчета изгибаемых плит с учетом деформаций сдвига
Ю.Я. Тюкалов*,
Вятский государственный университет, г. Киров, Россия * E-mail: [email protected]
Ключевые слова: изгибаемые плиты, деформации сдвига, напряжения, конечные элементы.
Аннотация. Рассмотрена задача расчета изгибаемых плит методом конечных элементов с учетом деформаций сдвига. Изгибаемые плиты широко применяются для различных объектов гражданского и промышленного строительства. Решение задачи получено на основе принципов минимума дополнительной энергии и возможных перемещений. Для аппроксимации полей моментов используются кусочно-постоянные функции. Поперечные силы могут быть аппроксимированы постоянными или кусочно-постоянными функциями. Получены необходимые соотношения для прямоугольных и треугольных конечных элементов. Показано, что предлагаемый метод может использоваться в сочетании с традиционными конечными элементами для тонких пластин, полученными методом конечных элементов в перемещениях. Приведено сравнение решений, полученных по предлагаемому методу, с другими известными решениями для изгибаемых плит с учетом сдвига. Дана численная оценка точности и сходимости предлагаемого метода при измельчении сетки конечных элементов.
Литература
1. Галлагер Р. Метод конечных элементов. Основы. М.: Мир, 1984. 428 с.
2. Сливкер В. И. Строительная механика. Вариационные основы. М.: Изд-во Ассоциации строительных вузов. 2005. 736 с.
3. Тимошенко С.П., Войновский-Кригер С. Пластины и оболочки. М.: Наука, 1966. 636 с.
4. Fallah N. On the use of shape functions in the cell centered finite volume formulation for plate bending analysis based on Mindlin— Reissner plate theory // Computers & Structures. 2006. Vol. 84. Pp. 1664-1672.
5. Katili I., Batoz J.L., Maknun I.J., Hamdouni A., Millet O. The development of DKMQ plate bending element for thick to thin shell analysis based on the Naghdi/Reissner/Mindlin shell theory // Finite Elements in Analysis and Design. 2015. Vol. 100. Pp. 12-27.
6. Khezri M., Gharib M., Rasmussen K. J. R. A unified approach to meshless analysis of thin to moderately thick plates based on a shear-locking-free Mindlin theory formulation. Thin-Walled Structures. 2018. Pp. 161-179.
7. Yaghoubshahi M., Alinia M.M. Developing an element free method for higher order shear deformation analysis of plates // Thin-Walled Structures. 2015, Vol. 94. Pp. 225-233.
8. Thai H.-T., Nguyen T.-K., Vo T.P., Ngo T. A new simple shear deformation plate theory // Composite Structures. 2017. Vol. 171. Pp. 277-285.
9. Лалин В.В., Зданчук Е.В., Кушова Д.А., Розин Л.А. Вариационные постановки нелинейных задач с независимыми вращательными степенями свободы // Инженерно-строительный журнал. 2015. № 4(56). С. 54-65.
10. Лалин В.В., Беляев М.О. Изгиб геометрически нелинейного консольного стержня. Решение по теориям Кирхгофа и Коссера-Тимошенко // Инженерно-строительный журнал. 2015. №1(53). С. 39-55.
11. Park M., Choi D.-H. A two-variable first-order shear deformation theory considering in-plane rotation for bending, buckling and free vibration analysis of isotropic plates // Applied Mathematical Modelling. 2018. Vol. 61. Pp. 49-71.
12. Kumara R., Lala A., Singhb B.N., Singhc J. New transverse shear deformation theory for bending analysis of FGM plate under patch load // Composite Structures. 2019. Vol. 208. Pp. 91-100.
13. Thai H.-T., Choi D.-H. Finite element formulation of various four unknown shear deformation theories for functionally graded plates // Finite Elements in Analysis and Design. 2013. Vol. 75. Pp. 50-61.
14. Doa T.V., Nguyenb D.K., Ducc N.D., Doanc D.H., Buie T.Q. Analysis of bi-directional functionally graded plates by FEM and a new third-order shear deformation plate theory // Thin-Walled Structures. 2017. Vol. 119. Pp. 687-699.
15. Mulcahy N.L., McGuckin D.G. The addition of transverse shear flexibility to triangular thin plate elements // Finite Elements in Analysis and Design. 2012. Vol. 52. Pp. 23-30.
16. Karttunen A.T., Hertzen R., Reddy J.N., Romanoff J. Shear deformable plate elements based on exact elasticity solution // Computers & Structures. 2018. Vol. 200. Pp. 21-31.
17. Duan H., Ma J. Continuous finite element methods for Reissner-Mindlin plate problem // Acta Mathematica Scientia. 2018. Vol. 38. Pp. 450-470.
18. Zenkour AM. Exact mixed-classical solutions for the bending analysis of shear deformable rectangular plates // App. Math. Model. 2003. Vol. 27. Pp. 515-534.
19. Erkmen R.E. Shear deformable hybrid finite-element formulation for buckling analysis of thin-walled members // Finite Elements in Analysis and Design. 2014. Vol. 82. Pp. 32-45.
20. Ozutok A., Madenci E. Static analysis of laminated composite beams based on higher-order shear deformation theory by using mixed-type finite element method // International Journal of Mechanical Sciences. 2017. Vol. 130. Pp. 234-243
21. Senjanovic I., Vladimir N., Hadzic N. Modified Mindlin plate theory and shear locking free finite element formulation // Mechanics Research Communications. 2014. Vol. 55. Pp. 95-104.
22. Сухотерин М.В., Барышников С.О., Потехина Е.В. О расчетах пластин по сдвиговым теориям // Вестник государственного университета морского и речного флота им. адмирала С.О. Макарова. 2015. № 2 (30). С. 81-89.
23. Сухотерин М.В., Барышников С.О., Кныш Т.П. Напряженно-деформированное состояние защемленной прямоугольной пластины Рейсснера // Инженерно-строительный журнал. 2017. № 8 (76). С. 225-240. DOI: 10.18720/MCE.76.20
24. Тюкалов Ю.Я. Конечно-элементные модели в напряжениях для задач плоской теории упругости // Инженерно-строительный журнал. 2018. № 1(77). С. 23-37. DOI: 10.18720/MCE.77.3
25. Тюкалов Ю.Я. Функционал дополнительной энергии для анализа устойчивости пространственных стержневых систем // Инженерно-строительный журнал. 2017. № 2(70). С. 18-32. DOI: 10.18720/MCE.70.3
26. Тюкалов Ю.Я. Определение частот свободных колебаний методом конечных элементов в напряжениях // Инженерно-строительный журнал. 2016. № 7(67). С. 39-54. DOI: 10.5862/MCE.67.5
27. Тюкалов Ю.Я. Конечно элементные модели в напряжениях для изгибаемых пластин // Инженерно-строительный журнал. 2018. № 6(82). С. 170-190. DOI: 10.18720/MCE.82.16
28. Доннелл Л.Г. Балки, пластины и оболочки. М.: Наука, 1982. 568 с.
Контактные данные:
Юрий Яковлевич Тюкалов, +7(912)8218977; эл. почта: [email protected]
© Тюкалов Ю.Я., 2019
