R-Function and variation method for bending problem of clamped thin plate with complex shape

Abstract

Solving ordinary thin plate bending problem in engineering, only a few analytical solutions with simple boundary shapes have been proposed. When using numerical methods (e.g. the variational method) to solve the problem, the trial functions can be found only it exhibits a simple boundary shape. The R-functions can be applied to solve the problem with complex boundary shapes. In the paper, the R-function theory is combined with the variational method to study the thin plate bending problem with the complex boundary shape. The paper employs the R-function theory to express the complex area as the implicit function, so it is easily to build the trial function of the complex shape thin plate, which satisfies with the complex boundary conditions. The variational principle and the R-function theory are introduced, and the variational equation of thin plate bending problem is derived. The feasibility and correctness of this method are verified by five numerical examples of rectangular, I-shaped, T-shaped, U-shaped, and L-shaped thin plates, and the results of this method are compared with that of other literatures and ANSYS finite element method (FEM). The results of the method show a good agreement with the calculation results of literatures and FEM.

Keywords

Complex shape thin plate bending R-function trial function variational method

Introduction

Thin plate refers to a flat plate with the thickness significantly smaller than the size of the plate surface. The problem of thin plate bending is recognized as a critical problem in the elastic theory. On the whole, the car body floor, the engine cylinder, the gearbox housing, the floor of building structure, the bridge deck, the aircraft panel, and the missile panel in aerospace engineering belong to thin plate bending problem.¹ Given Kirchhoff assumption, the governing equation of elastic thin plate bending problem is the biharmonic equation. To work out this type of higher order partial differential equation boundary value problem, a few regular plates with simple geometric shapes and boundary conditions can be obtained analytical solutions. So most engineering practical problems should be solved by using the finite element method (FEM), the boundary element method and other numerical methods.²

Two main methods are available for solving the mentioned problems, namely, the numerical method and the analytical method. Analytical methods include the superposition method, the quasi-Green function method, the trigonometric series method, the hierarchical wavelet method, the extended separation-of-variable method, the optimized separation-of-variable method, the finite integral transform method and the up-to-date symplectic superposition method. Numerical methods cover the boundary element method, the R-function method, the difference method, the finite element method, the weighted residual method, the high-order meshless method, the spline finite point method, the three-hinge-line method, and the wave-based method. In a wide range of methods, the deflection functions should be selected in advance, and the selection of the mentioned functions is arbitrary without any definite rule. In several novel analysis methods, the bending function is not required to be selected in advance (e.g. the symplectic superposition method and the finite integral transformation). Zhang et al. proposed a two-dimensional generalized finite integral transform method for new analytic bending solutions of orthotropic rectangular thin foundation plates.³ A Hamiltonian system-based variational principle via the Lagrangian multiplier method was proposed to formulate the thin plate buckling in the symplectic space.⁴ Hu et al. presented a first attempt to extend an up-to-date symplectic superposition method to linear buckling of side-cracked rectangular thin plates.⁵ An et al. developed a double finite integral transform method for analytical bending solutions of non-Levy-type cylindrical shell panels without a free edge that were not obtained by classical semi-inverse methods.⁶ Numerous engineering calculation problems (e.g. arbitrary geometry, complex boundary conditions, diverse, and uneven material properties) can be obtained directly from the mathematical model by numerical methods.⁷ The bending, vibration, and buckling of orthotropic rectangular thin plate under various boundary conditions by employing Fourier series and superposition principle with MATLAB programming were studied by Li⁸. Li et al. adopted the quasi-Green’s function method to study vibration analysis of simply supported polygonal thin plates on Winkler foundation.^9–12 Huang obtained the general solution of the differential equation of the bending problem of anisotropic rectangular thin plates by exploiting the trigonometric series solution as complex numbers and the double sinusoidal series solution, and solved the bending problems of various boundaries under arbitrary loads.¹³ Yu applied a hierarchical wavelet homotopy methodology for studying the nonlinear bending of anisotropic thin plates.¹⁴ Xing and Wang¹⁵ adopted a method to solve closed-form analytical solutions of the free vibration problems of orthotropic rectangular thin plates with arbitrary homogeneous boundary conditions. The optimized separation-of-variable method for the free vibration of orthotropic rectangular thin plates was investigated by Xing and Wang.¹⁶ The finite integral transform method was developed to conduct the bending analysis of thin plates with the combination of simply supported, clamped, and free boundary conditions.¹⁷ Hu et al.¹⁸ developed an advanced symplectic superposition method for several novel analytic free vibration solutions of side-cracked rectangular thin plates. Feng¹⁹ employed the boundary element method to solve the bending problem of thin plate with small deflection and using variational principle, and derived the boundary integral equation of the bending vibration of the base plate and thin plate on Winkler elastic foundation. Dong²⁰ used the boundary element method to study the bending problem of elastic thin plate in inhomogeneous infinite domain. The boundary element method was adopted to reduce the dimension of the problem and simplify the calculation.²¹ Chen and Xiao²² studied the numerical characteristics of the preconditioned GMRES algorithm in solving the boundary element method equations of elastic mechanics problems.²² Kurpa and Mazur²³ proposed a method of R-function theory for studying parametric vibration of complex shape plates with irregular shapes and different boundary conditions. Kurpa et al. applied the R-function theory for the free vibration of shallow shells, combined the spline approximation method, and the Ritz method to analyze the stress-strain state of the laminated shallow shell and the vibration of the laminated shallow shell with complex notch.^24–26 The variational-difference method was used to study thin plate bending problem.²⁷ The classical theory of the thin plate bending problem according to the generalized variational principle of three types of variables in 3D elastic theory was analyzed by Zhang et al.²⁸ The bending problems of Kirchhoff plate and Winkler plate under different transverse loads effectively was researched by generalized difference method.²⁹ The approximate calculation method for rectangular thin plate by using the principle of minimum potential energy and calculating the deflection of thin plate under uniform load was investigated by Dong and Lv.³⁰ The bending vibration problems of rectangular thin plates under simply supported boundary and fixed boundary conditions was studied by conducting the finite element analysis.³¹ A method for solving the bending problems of thin plates with arbitrary geometric boundaries and supporting conditions was studied by combining the weighted residual method.³² The application of small deflection theory of elastic thin plate in engineering design was conducted by Xu et al.³³ Analytical solutions for bending, buckling, and vibration analyses of thick rectangular plates with various boundary conditions were presented using two variable refined plate theory.³⁴ Based on the parametrized mixed variational principle, a high-performance finite element model has been developed for bending and vibration analysis of thick plates.³⁵ A novel mixed-field theory with relatively low number of unknown variables was introduced for bending and vibration analysis of multi-layered composite plates.³⁶ A four-node quadrilateral partial mixed plate element with low degrees of freedom was developed to conduct the static and free vibration analysis of functionally graded material plates rested on Winkler-Pasternak elastic foundations.³⁷ A spline finite point method was developed to study the nonlinear bending behavior of functionally graded material plates with different thickness in thermal environments.³⁸ Wang et al. presented a mesh free Galerkin method to study SG thin beams and plates to satisfy the continuity and convergence requirement. Besides, moving least square or reproducing kernel shape functions were employed with the cubic approximation bases.³⁹ The three-hinge-line method was adopted to study the crushing deformation of sheet metals.⁴⁰ Liu et al.⁴¹ applied finite element method and wave-based method simulations to study data on the flexural vibration of thin plates.

The Ritz method which is a common energy method has been extensively used in mechanical research because of its simple solving process.⁴² The R-function theory can be used to describe complex regions as implicit functions. The key is to express the boundary of the problem by standardized equation $ω = 0$ , and the region of the problem is expressed by inequality $ω > 0$ .²³ The introduction of the R-function theory can easily build the deflection function w of the irregular thin plate to meet the boundary conditions, and then the specific deflection function of the concave-convex irregular thin plate can be obtained by combining the variational method.

It is required to develop efficient and accurate numerical analysis methods for mechanical problems (e.g. bending, vibration, and buckling problem of thin plate) to ensure the design and safety of plate structures in practical engineering.⁴³ Besides, the governing equation of thin plate bending problem is a fourth order partial differential equation with deflection as independent variable. Developing numerical methods for the mentioned problem is of great scientific significance for solving high-order partial differential equations.⁴⁴

The paper introduces the calculation principle and gives numerical examples of rectangular, I-shaped, L-shaped, T-shaped, and U-shaped thin plates. After deducing the variational equation, the results are calculated by MATLAB program, and the feasibility of the method is verified by comparing with the results of other literatures and ANSYS finite element method (FEM) calculation results which show a good agreement.

Calculation of small deflection bending of thin plate by variational method

To the small deflection bending problem of thin plates, according to the calculation assumption, the deformation components $ε_{z}$ , $γ_{zx}$ , $γ_{zy}$ are not considered, thus the expression of deformation potential energy is simplified to the expression as follows

V_{ε} = \frac{1}{2} \int \int \int_{A} (σ_{x} ε_{x} + σ_{y} ε_{y} + τ_{xy} γ_{xy}) dxdydz

(1)

The stress component and strain component in the mentioned equation are expressed by deflection $w$ as follows

σ_{x} = - \frac{Ez}{1 - μ^{2}} (\frac{\partial^{2} w}{\partial x^{2}} + μ \frac{\partial^{2} w}{\partial y^{2}})

(2)

σ_{y} = - \frac{Ez}{1 - μ^{2}} (\frac{\partial^{2} w}{\partial y^{2}} + μ \frac{\partial^{2} w}{\partial x^{2}})

(3)

τ_{xy} = - \frac{Ez}{1 + μ} \frac{\partial^{2} w}{\partial x \partial y}

(4)

ε_{x} = - z \frac{\partial^{2} w}{\partial x^{2}}

(5)

ε_{y} = - z \frac{\partial^{2} w}{\partial y^{2}}

(6)

γ_{xy} = - 2 z \frac{\partial^{2} w}{\partial x \partial y}

(7)

Substituting equations (2)–(7) into equation (1), it yields

\begin{matrix} V_{ε} = \frac{E}{2 (1 - μ^{2})} \\ \int \int \int_{A} z^{2} {{(\nabla^{2} w)}^{2} - 2 (1 - μ) [\frac{\partial^{2} w}{\partial x^{2}} \frac{\partial^{2} w}{\partial y^{2}} - {(\frac{\partial^{2} w}{\partial x \partial y})}^{2}]} dxdydz \end{matrix}

(8)

For equation (8) each item does not change with z, integral z, from $- δ / 2$ to $δ / 2$ , and using the equation $D = E δ^{3} / [12 (1 - μ^{2})]$ , it yields

V_{ε} = \frac{1}{2} \int \int_{A} D {{(\nabla^{2} w)}^{2} - 2 (1 - μ) [\frac{\partial^{2} w}{\partial x^{2}} \frac{\partial^{2} w}{\partial y^{2}} - {(\frac{\partial^{2} w}{\partial x \partial y})}^{2}]} dxdy

(9)

The flexural stiffness D in the equation, in the case of a plate with variable thickness, will be a function of x and y, since δ will be a function of x and y.

In the thin plate of equal thickness, $D$ is a constant, and equation (8) can be rewritten as

V_{ε} = \frac{D}{2} \int \int_{A} {{(\nabla^{2} w)}^{2} - 2 (1 - μ) [\frac{\partial^{2} w}{\partial x^{2}} \frac{\partial^{2} w}{\partial y^{2}} - {(\frac{\partial^{2} w}{\partial x \partial y})}^{2}]} dxdy

(10)

where

\begin{matrix} \int \int_{A} (\frac{\partial^{2} w}{\partial x^{2}} \frac{\partial^{2} w}{\partial y^{2}} - {(\frac{\partial^{2} w}{\partial x \partial y})}^{2}) dxdy \\ = \int \int_{A} [\frac{\partial}{\partial x} (\frac{\partial w}{\partial x} \frac{\partial^{2} w}{\partial y^{2}}) - \frac{\partial}{\partial y} (\frac{\partial w}{\partial x} \frac{\partial^{2} w}{\partial x \partial y})] dxdy \end{matrix}

(11)

According to Green’s theorem, there are

\begin{matrix} \int \int_{A} (\frac{\partial}{\partial x} P (x, y) - \frac{\partial}{\partial y} Q (x, y)) dxdy \\ = \int [Q (x, y) dx + P (x, y) dy] \end{matrix}

(12)

According the equation (12), it yields

\begin{matrix} \int \int_{A} (\frac{\partial^{2} w}{\partial x^{2}} \frac{\partial^{2} w}{\partial y^{2}} - {(\frac{\partial^{2} w}{\partial x \partial y})}^{2}) dxdy \\ = \int [(\frac{\partial w}{\partial x} \frac{\partial^{2} w}{\partial y \partial x}) dx - (\frac{\partial w}{\partial x} \frac{\partial^{2} w}{\partial y^{2}}) dy] \end{matrix}

(13)

The line integral on the right is along the boundary of the thin plate.

All the boundaries of the thin plate are fixed edges. Thus, regardless of the boundary shape, there is a first-order partial derivative on the boundary equal to zero.

The boundary conditions of the fixed edge are

\begin{matrix} {(w)}_{x = \pm a} = 0, {(\frac{\partial w}{\partial x})}_{x = \pm a} = 0; \\ {(w)}_{y = \pm b} = 0, {(\frac{\partial w}{\partial y})}_{y = \pm b} = 0 . \end{matrix}

V_{ε} = \frac{D}{2} \int \int_{A} {(\nabla^{2} w)}^{2} dxdy

(14)

According to the calculation assumption and geometric equation in the small deflection bending problem of thin plate, the displacement component $u$ and $v$ can be expressed by deflection w, which need not be taken as the basic unknown function. Thus, there’s only one w fundamentally unknown function. Now set the w expression as follows

w = Σ_{m} C_{m} w_{m}

(15)

Where, $C_{m}$ are m individual undetermined coefficients independent of each other, $w_{m}$ is the setting function satisfying the displacement boundary condition of thin plate. In this way, regardless of the value of $C_{m}$ , the deflection $w$ shown in equation (15) can always satisfy the displacement boundary condition. Note that the variation of deflection w is only realized by the variation of coefficient $C_{m}$ . As for the set function $w$ , it changes only with coordinates and is completely independent of the above variation.

To determine the coefficients $C_{m}$ , the expression is adopted as follows

\frac{\partial V_{ε}}{\partial C_{m}} = \int_{V} f_{z} w_{m} dV + \int_{S} \bar{f_{z}} w_{m} dS

(16)

To the bending problem of thin plates, both the body force and the surface force are attributed to the load $q$ , and there is $dS = dxdy$ on the panel, so equation (11) is written as

\frac{\partial V_{ε}}{\partial C_{m}} = \int \int_{A} q w_{m} dxdy

(17)

The $C_{m}$ can be obtained from the m individual linear equations, thus the deflection $w$ can be obtained by equation (15), and the internal force of the thin plate can be calculated.

R-function theory

According to R-function theory,^23–26 the $ω_{0}$ describing a complex region can be described by Boolean operations $\lor_{α}$ and $\land_{α}$ of each subdomain, which are corresponding to the union set and intersection set. The function $ω_{l}$ describing each subdomain satisfies the first-order normalization equation, and $ω_{0}$ also satisfies the first-order normalization equation.

Definition of the first-order normalized equation as follows

ω_{l} (x, y) = 0, | \nabla ω_{l} | = 1, \forall (x, y) \in \partial Ω

ω_{l} (x, y) > 0, \forall (x, y) \in Ω

Let X and Y satisfy the definition of the first-order normalized equation, then the following Boolean operations $X \lor_{α} Y$ and $X \land_{α} Y$ are also the first-order normalized equations

X \lor_{α} Y = \frac{1}{1 + α} (X + Y + \sqrt{X^{2} + Y^{2} - 2 α XY})

(18)

X \land_{α} Y = \frac{1}{1 + α} (X + Y - \sqrt{X^{2} + Y^{2} - 2 α XY})

(19)

Where $- 1 \leq α \leq 1$ .

For the thin plate bending problem, the R-function is introduced into equation (15). This paper discusses the situation of the clamped edge, according to the properties of the R-function, the obtained $ω_{0}$ is squared to satisfy the boundary conditions of the clamped edge.

Derivation of variation equations

The expression of deflection is

w = {ω_{0}}^{2} Σ_{m} C_{m} w_{m}

(20)

Thus, the boundary conditions are expressed by $ω_{0}^{2}$ .

The expression of deformation potential energy of thin plate with clamped boundary conditions as follows

V_{ε} = \frac{D}{2} \int \int {(\nabla^{2} w)}^{2} dxdy

(21)

Namely

V_{ε} = \frac{D}{2} \int \int {(\frac{\partial^{2} w}{\partial x^{2}} + \frac{\partial^{2} w}{\partial y^{2}})}^{2} dxdy

(22)

Substituting equation (18) into equation (19), we can get

\begin{matrix} \frac{\partial^{2} w}{\partial x^{2}} = \sum_{m} C_{m} \\ ((2 {(\frac{\partial ω_{0}}{\partial x})}^{2} + 2 ω_{0} \frac{\partial^{2} ω_{0}}{\partial x^{2}}) w_{m} + 4 ω_{0} \frac{\partial ω_{0}}{\partial x} \frac{\partial w_{m}}{\partial x} + {ω_{0}}^{2} \frac{\partial^{2} w_{m}}{\partial x^{2}}) \end{matrix}

(23)

\begin{matrix} \frac{\partial^{2} w}{\partial y^{2}} = \sum_{m} C_{m} \\ ((2 {(\frac{\partial ω_{0}}{\partial y})}^{2} + 2 ω_{0} \frac{\partial^{2} ω_{0}}{\partial y^{2}}) w_{m} + 4 ω_{0} \frac{\partial ω_{0}}{\partial y} \frac{\partial w_{m}}{\partial y} + {ω_{0}}^{2} \frac{\partial^{2} w_{m}}{\partial^{2}}) \end{matrix}

(24)

Take the derivative of $C_{i}$ $(i = 1, 2, 3 \dots \dots m - 1, m)$ , we can obtain

\begin{matrix} \frac{\partial}{\partial C_{i}} \frac{\partial^{2} w}{\partial x^{2}} = \\ (2 {(\frac{\partial ω_{0}}{\partial x})}^{2} + 2 ω_{0} \frac{\partial^{2} ω_{0}}{\partial x^{2}}) w_{i} + 4 ω_{0} \frac{\partial ω_{0}}{\partial x} \frac{\partial w_{i}}{\partial x} + {ω_{0}}^{2} \frac{\partial^{2} w_{i}}{\partial x^{2}} \end{matrix}

(25)

\begin{matrix} \frac{\partial}{\partial C_{i}} \frac{\partial^{2} w}{\partial y^{2}} = \\ (2 {(\frac{\partial ω_{0}}{\partial y})}^{2} + 2 ω_{0} \frac{\partial^{2} ω_{0}}{\partial y^{2}}) w_{i} + 4 ω_{0} \frac{\partial ω_{0}}{\partial y} \frac{\partial w_{i}}{\partial y} + {ω_{0}}^{2} \frac{\partial^{2} w_{i}}{\partial y^{2}} \end{matrix}

(26)

Substitute the R-function expression into equation (17), the follow expression can be obtained

\frac{\partial V_{ε}}{\partial C_{m}} = \int \int q {ω_{0}}^{2} w_{m} dxdy

(27)

Take the derivative of $C_{i}$ $(i = 1, 2, 3 \dots \dots m - 1, m)$ , we can get

\begin{matrix} \frac{\partial V_{ε}}{\partial C_{i}} = \frac{\partial}{\partial C_{i}} \frac{D}{2} \int \int {(\frac{\partial^{2} w}{\partial x^{2}} + \frac{\partial^{2} w}{\partial y^{2}})}^{2} dxdy \\ = D \int \int (\frac{\partial^{2} w}{\partial x^{2}} + \frac{\partial^{2} w}{\partial y^{2}}) \frac{\partial}{\partial C_{i}} (\frac{\partial^{2} w}{\partial x^{2}} + \frac{\partial^{2} w}{\partial y^{2}}) dxdy \\ = D \int \int (\frac{\partial^{2} w}{\partial x^{2}} + \frac{\partial^{2} w}{\partial y^{2}}) (\frac{\partial}{\partial C_{i}} \frac{\partial^{2} w}{\partial x^{2}} + \frac{\partial}{\partial C_{i}} \frac{\partial^{2} w}{\partial y^{2}}) dxdy \end{matrix}

(28)

Substituting equations (22)–(25) into equation (27), we can get

\begin{matrix} \frac{\partial V_{ε}}{\partial C_{i}} = D {\begin{matrix} \sum_{m} C_{m} (\begin{matrix} 2 w_{i} ({(\frac{\partial ω_{0}}{\partial x})}^{2} + {(\frac{\partial ω_{0}}{\partial y})}^{2}) + 2 ω_{0} w_{i} (\frac{\partial^{2} ω_{0}}{\partial x^{2}} + \frac{\partial^{2} ω_{0}}{\partial y^{2}}) + \\ 4 ω_{0} (\frac{\partial ω_{0}}{\partial x} \frac{\partial w_{i}}{\partial x} + \frac{\partial ω_{0}}{\partial y} \frac{\partial w_{i}}{\partial y}) + {ω_{0}}^{2} (\frac{\partial^{2} w_{i}}{\partial x^{2}} + \frac{\partial^{2} w_{i}}{\partial y^{2}}) \end{matrix}) \times \\ (\begin{matrix} 2 w_{j} ({(\frac{\partial ω_{0}}{\partial x})}^{2} + {(\frac{\partial ω_{0}}{\partial y})}^{2}) + 2 ω_{0} w_{j} (\frac{\partial^{2} ω_{0}}{\partial x^{2}} + \frac{\partial^{2} ω_{0}}{\partial y^{2}}) + \\ 4 ω_{0} (\frac{\partial ω_{0}}{\partial x} \frac{\partial w_{j}}{\partial x} + \frac{\partial ω_{0}}{\partial y} \frac{\partial w_{j}}{\partial y}) + {ω_{0}}^{2} (\frac{\partial^{2} w_{j}}{\partial x^{2}} + \frac{\partial^{2} w_{j}}{\partial y^{2}}) \end{matrix}) \end{matrix}} dxdy \end{matrix}

(29)

The coefficient equation system can be obtained by derivation as follow

A {[\begin{matrix} C_{1} & C_{2} & C_{3} & \dots & \dots & C_{m} \end{matrix}]}^{T} = B

(30)

Where

A = {(a_{ij})}_{m \times m},

a_{ij} = D \int \int (\begin{matrix} (\begin{matrix} 2 w_{i} ({(\frac{\partial ω_{0}}{\partial x})}^{2} + {(\frac{\partial ω_{0}}{\partial y})}^{2}) + 2 ω_{0} w_{i} (\frac{\partial^{2} ω_{0}}{\partial x^{2}} + \frac{\partial^{2} ω_{0}}{\partial y^{2}}) + \\ 4 ω_{0} (\frac{\partial ω_{0}}{\partial x} \frac{\partial w_{i}}{\partial x} + \frac{\partial ω_{0}}{\partial y} \frac{\partial w_{i}}{\partial y}) + {ω_{0}}^{2} (\frac{\partial^{2} w_{i}}{\partial x^{2}} + \frac{\partial^{2} w_{i}}{\partial y^{2}}) \end{matrix}) \times \\ (\begin{matrix} 2 w_{j} ({(\frac{\partial ω_{0}}{\partial x})}^{2} + {(\frac{\partial ω_{0}}{\partial y})}^{2}) + 2 ω_{0} w_{j} (\frac{\partial^{2} ω_{0}}{\partial x^{2}} + \frac{\partial^{2} ω_{0}}{\partial y^{2}}) + \\ 4 ω_{0} (\frac{\partial ω_{0}}{\partial x} \frac{\partial w_{j}}{\partial x} + \frac{\partial ω_{0}}{\partial y} \frac{\partial w_{j}}{\partial y}) + {ω_{0}}^{2} (\frac{\partial^{2} w_{j}}{\partial x^{2}} + \frac{\partial^{2} w_{j}}{\partial y^{2}}) \end{matrix}) \end{matrix}) dxdy

(31)

B = {(b_{ij})}_{m \times 1},

b_{ij} = \int_{Ω} q {ω_{0}}^{2} w_{i} d Ω

(32)

Here $i = 1, 2, \dots \dots, m$ ; $j = 1, 2, \dots \dots m$ , $Ω$ is the area contained in the plane of thin plate (x–y) in the Cartesian coordinate system, which $C_{m}$ can be obtained from equation (17), and $C_{m}$ can be substituted into the equation (20) to obtain the deflection expression.

Numerical examples

Example 1 Rectangular thin plate

Rectangular thin plates shown in Figure 1 with side lengths of 2a and 2b, where $a = b = 0.5$ , the elastic modulus $E$ of the material is $100 GPa$ , the Poisson’s ratio $ν$ is 0.3, the thickness $δ$ of the thin plate is $0.05 m$ , and the uniform load is $q = 10^{6} Pa$ .

Figure 1.

Rectangular thin plate.

According to the R-function theory,^23–26 $ω_{0}$ can be taken as follows

ω_{0} = ω_{1} \land_{α} ω_{2},

Where $ω_{1} = (a^{2} - x^{2}) / (2 a) \geq 0$ and $ω_{2} = (b^{2} - y^{2}) / (2 b) \geq 0$ .

Setting the deflection $w$ to be an even function of $x$ and $y$ , so $w_{m}$ take $1$ , $x^{2}$ , $y^{2}$ , $x^{4}$ , $y^{4}$ , $x^{2} y^{2}$ , $x^{6}$ , $y^{6}$ , $x^{2} y^{4}$ , and $x^{4} y^{2}$ .

Combining R-function with the trial function of variational method, in the process of calculating integral, it will appear that the denominator of integral function is zero on the boundary of integral. In this case, it is impossible to calculate the integral of integral function. Thus, the rectangular area is divided to $n \times n$ numbers by gridding.

When number of the trial function is $m = 10$ , deflection of the center point of thin plate is listed in Table 1, and the number of grids are $15 \times 15$ , $20 \times 20$ , $25 \times 25$ , $30 \times 30$ , $35 \times 35$ , $40 \times 40$ . When the number of grids is $40 \times 40$ , the deflection results of center point are listed in Table 2, and $m$ takes different values. It can be seen from Tables 1 and 2 that the deflection of thin plate calculated by the method in this paper is convergent and effective, and it is a good agreement with the analytic solution⁴⁵ and the FEM solution. These prove the feasibility and correctness of this method. When the number of grids is $40 \times 40$ , $m$ take different values, the deflection curve of the center line $y = 0$ of the plate is shown in Figure 2 which show the convergence of the method.

Table 1.

Deflection of center point with different n.

	Grid density	Deflection of center point (m)	Error (%)
The method in this paper	$15 \times 15$	0.002130	93.66
	$20 \times 20$	0.001303	18.48
	$25 \times 25$	0.001154	4.95
	$30 \times 30$	0.001121	1.92
	$35 \times 35$	0.001111	1.02
	$40 \times 40$	0.001108	0.69
Exact solution³⁷		0.0011	−
FEM solution		0.001159	5.09

Table 2.

Deflection of center point with different m.

	m	Deflection of center point(m)	Error (%)
The method in this paper	1	0.000524	52.36
	3	0.00105	4.55
	5	0.001121	1.91
	6	0.001121	1.91
	8	0.001106	0.55
	10	0.0011076	0.69
Exact solution³⁷		0.0011	–

Figure 2.

The curve of w at y = 0.

In FEM modeling, SHELL181 is selected as the element type of the thin plate, and the edge length of each element of the rectangular thin plate is 0.05 m. The rectangular thin plate is divided into a total of 400 square grids. When m = 10, the results of the paper and FEM are presented in Figure 3 for comparison, and which shows a good agreement.

Figure 3.

Comparison between the method and FEM at y = 0.

Example 2 I-shaped thin plate

A I-shaped thin plate is shown in Figure 4. The elastic modulus $E$ of the material is $100 Gpa$ , the Poisson’s ratio $ν$ is 0.3, the thickness $δ$ of the thin plate is $0.05 m$ , and the uniform load is $q = 5 \times 10^{5} pa$ . We set $a = b = 1$ , $c = d = 0.5$ .

Figure 4.

I-shaped thin plate.

According to the R-function theory,^23–26 $ω_{0}$ can be taken

ω_{α} = (ω_{1} \land_{α} ω_{2}) \land_{α} (ω_{3} \lor_{α} ω_{4}),

Where $ω_{1} = \frac{a^{2} - x^{2}}{2 a} \geq 0, ω_{2} = \frac{b^{2} - y^{2}}{2 b} \geq 0, ω_{3} = \frac{c^{2} - x^{2}}{2 c} \geq 0$ and $ω_{4} = \frac{y^{2} - d^{2}}{2 d} > 0$ .

$w_{m}$ take $1$ , $x^{2}$ , $y^{2}$ , $x^{4}$ , $y^{4}$ , $x^{2} y^{2}$ , $x^{6}$ , $y^{6}$ , $x^{2} y^{4}$ , and $x^{4} y^{2}$ .

When programming with MATLAB, thin plate is divided into three small rectangles, the first rectangular part is divided into $4 n \times n$ , the second rectangle is divided into $2 n \times 2 n$ , and the third rectangle is divided into $4 n \times n$ . When $m = 10$ , and n takes different values, the results are listed in Table 3. When $n = 9$ , and m takes different values, the results of center point are listed in Table 4, and the deflection curve of the center line $y = 0$ of the thin plate is calculated which shown in Figure 5 which show the convergence of the method. In FEM model, SHELL181 is taken as the element type of the thin plate, and the respective element of the I-shaped thin plate takes the edge length of 0.05 m. The I-shaped thin plate is divided into a total of 1200 square grids. When $n = 9$ and $m = 10$ , the results are compared with the FEM results which shows a good agreement in Figure 6, and that proves the convergence and correctness of the method in this paper.

Table 3.

Deflection of center point with different n.

	Grid density	Deflection of center point (m)	Error (%)
The method in this paper	$n = 3$	0.001685	46.39
	$n = 4$	0.001169	1.56
	$n = 5$	0.001087	5.56
	$n = 6$	0.001085	5.73
	$n = 7$	0.001096	4.78
	$n = 8$	0.001107	3.82
	$n = 9$	0.001115	3.13
	FEM solution	0.001151	−

Table 4.

Deflection of center point with different m.

	m	Deflection of center point (m)	Error (%)
The method in this paper	1	0.000304	73.59
	3	0.000584	49.26
	5	0.000947	17.23
	6	0.000987	14.25
	8	0.001089	5.39
	10	0.001115	3.13
	FEM solution	0.001151	–

Figure 5.

The curve of w at y = 0.

Figure 6.

Comparison between the method and FEM at y = 0.

Example 3 T-shaped thin plate

A T-shaped thin plate is shown in Figure 7. The elastic modulus $E$ of the material is $100 Gpa$ , the Poisson’s ratio $ν$ is 0.3, the thickness $δ$ of the thin plate is $0.005 m$ , and the uniform load $q = 10^{3}$ , and set $a = b = 0.5$ , $c = d = 0.25$ .

Figure 7.

T-shaped thin plate.

According to the R-function theory.^23–26 $ω_{0}$ can be taken

ω_{0} = (ω_{1} \land_{α} ω_{2}) \land_{α} (ω_{3} \lor_{α} ω_{4}),

Where $ω_{1} = \frac{a^{2} - x^{2}}{2 a} \geq 0, ω_{2} = \frac{b^{2} - y^{2}}{2 b} \geq 0, ω_{3} = \frac{c^{2} - x^{2}}{2 c} \geq 0$ and $ω_{4} = (d - y) > 0$ .

$w_{m}$ take $1$ , $y$ , $x^{2}$ , $y^{2}$ , $y^{3}$ , $x^{2} y$ , $x^{4}$ , $y^{4}$ , $x^{2} y^{2}$ , and $y^{5}$ .

When programming with MATLAB, the T-shaped thin plate is divided into three rectangles. The first rectangular part is divided into $2 n \times n$ grids, the second rectangular part is divided into $4 n \times 2 n$ , and the third rectangular part is divided into $4 n \times n$ . When $m = 10$ , $n$ takes different values, the results are listed in Table 5. When $n = 9$ , $m$ takes different values, the results of center point are listed in Table 6, and the deflection curve of center line $y = 0$ of thin plate is calculated as shown in Figure 8 which show the convergence of the method. In FEM model, SHELL181 is taken as the element type of the thin plate, and the element of the T-shaped thin plate is the edge length of 0.05 m. The T-shaped thin plate is divided into a total of 350 rectangular grids. When $n = 9$ and $m = 10$ , the results of the paper and FEM are compared which shows a good agreement in Figure 9, and that proves the convergence and correctness of this method.

Table 5.

Deflection of center point with different n.

	Grid density	Deflection of center point (m)	Error (%)
The method in this paper	$n = 3$	0.001250	89.68
	$n = 4$	0.000870	32.02
	$n = 5$	0.000751	13.96
	$n = 6$	0.000709	7.59
	$n = 7$	0.000690	4.70
	$n = 8$	0.000680	3.19
	$n = 9$	0.000668	1.37
	FEM solution	0.000659	−

Table 6.

Deflection of center point with different m.

	m	Deflection of center point (m)	Error (%)
The method in this paper	1	0.000282	57.21
	3	0.000401	39.15
	5	0.000597	9.41
	6	0.000607	7.89
	8	0.000674	2.28
	10	0.000665	1.37
	FEM solution	0.000659	–

Figure 8.

The curve of w at y = 0.

Figure 9.

Comparison between the method and FEM at y = 0.

Example 4 U-shaped thin plate

A U-shaped thin plate is shown in Figure 10. The elastic modulus $E$ of the material is $100 Gpa$ , the Poisson’s ratio $ν$ is 0.3, the thickness $δ$ of the thin plate is $0.005 m$ , and the uniform load is $q = 10^{3}$ , and set $a = b = 1$ , $c = d = 0.5$ .

Figure 10.

U-shaped thin plate.

According to the R-function theory,^23–26 $ω_{0}$ can be taken

ω_{0} = (ω_{1} \land_{α} ω_{2}) \land_{α} (ω_{3} \lor_{α} ω_{4}),

Where $ω_{1} = \frac{a^{2} - x^{2}}{2 a} \geq 0, ω_{2} = \frac{b^{2} - y^{2}}{2 b} \geq 0, ω_{3} = (c - x) \geq 0$ and $ω_{4} = \frac{y^{2} - d^{2}}{2 d} \geq 0$ .

$w_{m}$ take $1$ , $y$ , $x^{2}$ , $y^{2}$ , $y^{3}$ , $x^{2} y$ , $x^{4}$ , $y^{4}$ , $x^{2} y^{2}$ , and $y^{5}$ .

When programming with MATLAB, the U-shaped thin plate is divided into three rectangles. The first rectangular part is divided into $4 n \times n$ grids, the second rectangular part is divided into $3 n \times 2 n$ , and the third rectangular part is divided into $4 n \times n$ . When $m = 10$ , $n$ takes different values, the results of center point are listed in Table 7. When $n = 9$ , $m$ takes different values, the results of center point are listed in Table 8, and the deflection curve of center line $y = 0$ of thin plate is calculated as shown in Figure 11 which show the convergence of the method. In FEM model, SHELL181 is taken as the element type of the thin plate, and the element of the U-shaped thin plate is the edge length of 0.05 m. The U-shaped thin plate is divided into a total of 1400 rectangular grids. When $n = 9$ and $m = 10$ , the results of the paper and FEM are compared which shows a good agreement in Figure 12, and that proves the convergence and correctness of this method.

Table 7.

Deflection of center point with different n.

	Grid density	Deflection of center point (m)	Error (%)
The method in this paper	$n = 3$	0.0018480	69.53
	$n = 4$	0.0013582	24.59
	$n = 5$	0.0012026	10.32
	$n = 6$	0.0011463	5.16
	$n = 7$	0.0011203	2.77
	$n = 8$	0.0011061	1.47
	$n = 9$	0.0010980	0.72
	FEM solution	0.0010901	−

Table 8.

Deflection of center point with different m.

	m	Deflection of center point (m)	Error (%)
The method in this paper	1	0.000319	192.63
	3	0.000646	40.74
	5	0.000851	21.93
	6	0.000943	13.49
	8	0.001064	2.39
	10	0.001101	1.00
	12	0.001098	0.72
	FEM solution	0.0010901	–

Figure 11.

The curve of w at y = 0.

Figure 12.

Comparison between the method and FEM at y = 0.

Example 5 L-shaped thin plate

A L-shaped thin plate is shown in Figure 13. The elastic modulus $E$ of the material is $100 Gpa$ , the Poisson’s ratio $ν$ is 0.3, the thickness $δ$ of the thin plate is $0.005 m$ , and the uniform load is $q = 5 \times 10^{2} pa$ , and set $a = b = 0.5$ , $c = d = 0.125$ .

Figure 13.

L-shaped thin plate.

According to the R-function theory,^23–26 $ω_{0}$ can be taken

ω_{0} = (ω_{1} \land_{α} ω_{2}) \land_{α} (ω_{3} \lor_{α} ω_{4}),

Where $ω_{1} = \frac{a^{2} - x^{2}}{2 a} \geq 0, ω_{2} = \frac{b^{2} - y^{2}}{2 b} \geq 0, ω_{3} = (c - x) \geq 0$ and $ω_{4} = (d - y) \geq 0$ .

$w_{m}$ take $1$ , $x$ , $y$ , $x^{2}$ , $xy$ , $y^{2}$ , $x^{3}$ , $x^{2} y$ , $x y^{2}$ , and $y^{3}$ .

When programming with MATLAB, the L-shaped thin plate is divided into two rectangles. The first rectangular part is divided into $7 n \times n$ grids, the second rectangular part is divided into $8 n \times 7 n$ . When $m = 10$ , $n$ takes different values, the results of center point are listed in Table 9. When $n = 9$ , $m$ takes different values, the results of center point are listed in Table 10, and the deflection curve of center line $y = 0$ of thin plate is calculated as shown in Figure 14 which show the convergence of the method. In FEM model, SHELL181 is taken as the element type of the thin plate, and the element of the L-shaped thin plate is the edge length of 0.025 m. The L-shaped thin plate is divided into a total of 1575 rectangular grids. When $n = 9$ and $m = 10$ , the results of the paper and FEM are compared which shows a good agreement in Figure 15, and that proves the convergence and correctness of this method.

Table 9.

Deflection of center point with different n.

	Grid density	Deflection of center point (m)	Error (%)
The method in this paper	$n = 1$	0.0007416	36.07
	$n = 2$	0.0005072	6.94
	$n = 3$	0.0005129	5.89
	$n = 4$	0.0005176	5.03
	$n = 5$	0.0005195	4.68
	$n = 6$	0.0005203	4.53
	$n = 7$	0.0005206	4.48
	FEM solution	0.000545	−

Table 10.

Deflection of center point with different m.

	m	Deflection of center point (m)	Error (%)
The method in this paper	1	0.0002215	59.36
	3	0.0002214	59.38
	5	0.0003413	37.38
	6	0.0005154	5.43
	8	0.0005177	5.01
	10	0.0005206	4.48
	FEM solution	0.000545	–

Figure 14.

The curve of w at y = 0.

Figure 15.

Comparison between the method and FEM at y = 0.

In the calculation of example 2, example 3, example 4, and example 5, the final results have some errors with the results calculated by the FEM, which may be caused by three factors as follows.

Firstly, the trial function selected in this paper may meet the additional non-existence boundary conditions in the non-boundary area and impact the calculation results of this research. For instance, the trial function $x^{2}$ will constantly maintain zero on the line $x = 0$ , which reveals that the trial function will maintain zero deflection on the line $x = 0$ , whereas this obviously does not comply with the conditions under the consideration of this research.

Secondly, in this paper, the number of the trial functions may be used in the process of research is insufficient, the data of example 2 can indicate that when m is small, the results will cause a large error compared with the results of FEM. Thus, to use the variational method for obtaining the correct results, a maximal number of trial functions should be taken.

Thirdly, when the R-function is integral, the denominator is equated to zero on the boundary, which reveals that the R-function cannot be solved by the traditional symbolic integral method. It must be solved by the method of adding the grid after meshing consistent with the FEM, and this also means that the density of the grid will affect the correctness of the results. Accordingly, the R-function and the variational method should be employed to solve the problem.

Conclusions

The variational method, as an energy method, is capable of effectively calculating the deflection variation of the thin plate bending problem. The disadvantage of the variational method is that when the shape of the thin plate is very complex, the advantage of the R-function should be used. The advantage of the R-function is that the implicit function can be used to represent the complex boundary shape. Taking the advantages of the R-function, the combination of the R-function and the function of the variational method can be applied easily to calculate the deflection variation of the thin plate bending problem with complex boundary conditions.

In this paper, the effect of the combination of the R-function and the variational method is illustrated by the theoretical derivation and relevant examples. The advantage of R-function is that it can simply represent the complex boundary shape. For the thin plate often used in engineering, the advantages of R-function are reflected. As one of the energy methods, the variational method can easily calculate the deflection variation of the thin plate. Combining the variational method with the R-function theory to solve the problem with complex boundary conditions, the thin plate bending problem with complex shape can be easily calculated.

In this paper, the examples are presented to illustrate that the combination of the R-function and the variational method can be adopted to calculate the bending problem of thin plates with complex shapes. It is suggested that the combination of the R-function and the variational method can feasibly calculate the bending problem of thin plates with complex shapes. Moreover, the correlation between the results and the trial function of variational method and the number of integral grids is listed, thereby it carry forward the R-function on solving the bending problem of thin plates.

Footnotes

Handling Editor: James Baldwin

Declaration of conflicting interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: The authors gratefully acknowledge the financial support provided by the Science and Technology Scheme of Guangzhou City (no.201904010141).

ORCID iD

Shanqing Li

Data availability statement

The data used to support the findings of this study are available from the corresponding author upon request.

References

Xue

Deng

. Finite element analysis of elastic thin plate bending. J Shenyang Aerosp Univ 2011; 28(04): 9–11.

Tan

Zhang

. Double reciprocal hybrid boundary point method for bending problems of elastic thin plates. J Appl Mech 2011; 28(06): 583–588.

Zhang

Zhou

Ullah

, et al. Two-dimensional generalized finite integral transform method for new analytic bending solutions of orthotropic rectangular thin foundation plates. Appl Math Lett 2019; 92: 8–14.

Wang

. Symplectic superposition method for new analytic buckling solutions of rectangular thin plates. Int J Mech Sci 2016; 119: 432–441.

Zhaoyang

Zheng

, et al. New analytic buckling solutions of side-cracked rectangular thin plates by the symplectic superposition method. Int J Mech Sci 2021; 191: 106051.

, et al. Finite integral transform method for analytical solutions of static problems of cylindrical shell panels. Eur J Mech A Solids 2020; 83: 104033.

Chen

. The finite element analysis of rectangular thin plate under uniform load based on Matlab. J Chongqing Inst Technol 2009; 23(1): 32–35.

. Bending, vibration and buckling of orthotropic rectangular thin plates. Xi’an: Xi’an University of Architecture and Technology, 2007.

S-Q

Yuan

. Quasi-Green’s function method for free vibration of clamped thin plates on Winkler foundation. Appl Math Mech 2011; 32(3): 265–276.

10.

S-Q

Yuan

. Quasi-Green’s function method for free vibration of simply- supported trapezoidal shallow spherical shell. Appl Math Mech 2010; 31(5): 635–642.

11.

Yuan

. Green quasi-function method for free vibration of simply-supported trapezoidal shallow spherical shell on Winkler foundation. Acta Mechanica Solida Sinica 2010; 23(4): 370–376.

12.

Yuan

. Green quasi-function method for free vibration of clamped thin plates. Acta Mechanica Solida Sinica 2012; 25(1): 37–45.

13.

Huang

Liao

Qiu

. Bending solutions of anisotropic rectangular thin plates under arbitrary loads. J National Univ Defense Technol 2002; 1: 30–33.

14.

. A hierarchical wavelet method for nonlinear bending of materially and geometrically anisotropic thin plate. Commun Nonlinear Sci Numer Simul 2021; 92: 105498.

15.

Xing

Wang

. An extended separation-of-variable method for the free vibration of orthotropic rectangular thin plates. Int J Mech Sci 2020; 182: 105739.

16.

Xing

Wang

. An improved separation-of-variable method for the free vibration of orthotropic rectangular thin plates. Compos Struct 2020; 252: 112664.

17.

Zhang

Liu

Ullah

, et al. Analytical bending solutions of thin plates with two adjacent edges free and the others clamped or simply supported using finite integral transform method. Comput Appl Math 2020; 39(4): 1–23.

18.

Yang

Zhou

, et al. On the symplectic superposition method for new analytic free vibration solutions of side-cracked rectangular thin plates. J Sound Vib 2020; 489: 115695.

19.

Feng

. Establishing unified boundary integral equation of thin plate bending problem by variational principle. J Shandong Univ Archit Eng 1989; 1; 15–23.

20.

Dong

. Solutions to bending problems of inhomogeneous thin plates based on boundary element method. J Comput Mech 2011; 28(0z1): 25–28.

21.

Zhang

. A new way to calculate the bending problem of thin plates by boundary element method. J Hefei Univ Tech 1991; 1: 70–77.

22.

Chen

Xiao

. Preconditioned GMRES algorithm in boundary element analysis of thin plate bending problem. J Harbin Eng Univ 2018; 39(11): 1809–1815.

23.

Kurpa

Mazur

. Method of R-function for investigation of parametric vibrations of orthotropic plates of complex shape. Compos Struct 2011; 174(3): 269–282.

24.

Kurpa

Shmatko

Timchenko

. Free vibration analysis of laminated shallow shells with complex shape using the R-functions method. Compos Struct 2010; 93(1): 225–233.

25.

Awrejcewicz

Kurpa

Osetrov

. Investigation of the stress-strain state of the laminated shallow shells by R-functions method combined with spline-approximation. J Appl Math Mech 2011; 91(6): 458–467.

26.

Kurpa

Shmatko

Awrejcewicz

. Vibration analysis of laminated functionally graded shallow shells with clamped cutout of the complex form by the Ritz method and the R-functions theory. Latin Am J Solid Struct 2019; 16(1): 1–16.

27.

Xie

. Small curved steel Variation-difference method for thin plate bending problem. Comput Phys 2009; 26(1): 135–140.

28.

Zhang

Yang

Sun

. Another proof of classical theory of thin plate bending. J Luoyang Inst Technol 1994; 2: 61–65.

29.

Tang

Fan

. Generalized finite difference method for solving Kirchhoff and Winkler thin plate bending problems. J Solid Mech 2018; 39(4): 419–428.

30.

Dong

. Application of minimum potential energy principle in small deflection thin plate. J Jilin Inst Construct Eng 2003; 1994; 3: 59–62.

31.

Liu

. Calculation of bending problems of thin plates by finite analysis. J Hunan Univ 1994; 3: 36–40.

32.

Zhou

Liu

. A numerical method for the bending problem of thin plates. J Wuhan Univ Water Resource Electr Power 1996; 6; 53–57.

33.

Guan

Zhou

. Examples of application of small deflection theory of elastic thin plate in engineering design. Struct Eng 2016; 32(2): 40–45.

34.

Thai

H-T

Choi

D-H

. Analytical solutions of refined plate theory for bending, buckling and vibration analyses of thick plates. Appl Math Model 2013; 37(18–19): 8310–8323.

35.

Lezgy-Nazargah

. A high-performance parametrized mixed finite element model for bending and vibration analyses of thick plates. Acta Mechanica 2016; 227(12): 3429–3450.

36.

Lezgy-Nazargah

Salahshuran

. A new mixed-field theory for bending and vibration analysis of multi-layered composite plate. Arch Civil Mech Eng 2018; 18(3): 818–832.

37.

Lezgy-Nazargah

Meshkani

. An efficient partial mixed finite element model for static and free vibration analyses of FGM plates rested on two-parameter elastic foundations. Struct Eng Mech 2018; 66(5): 665–676.

38.

Huang

, et al. A spline finite point method for nonlinear bending analysis of FG plates in thermal environments based on a locking-free thin/thick plate theory. Mathematical Problems in Engineering, 2020; 21: 1–22.

39.

Wang

Fan

, et al. A stable and efficient meshfree Galerkin method with consistent integration schemes for strain gradient thin beams and plates. Thin-Walled Struct 2020; 153: 106791.

40.

Yao

Tan

, et al. A three-hinge-line method for calculating the deformation process of thin plates under quasi-static axial loading. Int J Crashworthiness 2019; 25(2): 167–174.

41.

Liu

Corradi

Ripamonti

, et al. Data on the flexural vibration of thin plate with elastically restrained edges: finite element method and wave based method simulations. Data Brief 2020; 31: 105883.

42.

Chen

Feng

. Bending vibration characteristics analysis of rectangular thin plate with complex open shape based on Chebyshev-variation method. Vib Impact 2020; 39(2): 178–187.

43.

Wang

Duan

. Node integral high order meshless method for bending analysis of thin plates. J Comput Mech 2019; 36(1): 103–109.

44.

Wang

Duan

. High order efficient meshless method for bending analysis of thin plates. J Solid Mech 2018; 39(2): 152–161.

45.

Timoshenko

Woinowsky-Krieger

. Theory of plates and shells. New York: McGraw-Hill, 1959.