The Scientific World Journal

Volume 2014 (2014), Article ID 520958, 11 pages

http://dx.doi.org/10.1155/2014/520958

Research Article

Numerical Manifold Method for the Forced Vibration of Thin Plates during Bending

Ding Jun,¹ Chen Song,¹ Wen Wei-Bin,² Luo Shao-Ming,² and Huang Xia¹

¹College of Mechanical Engineering, Chongqing University of Technology, Chongqing 400054, China

²Zhongkai University of Agriculture and Engineering, Guangzhou 510225, China

Received 8 February 2014; Revised 24 March 2014; Accepted 8 April 2014; Published 30 April 2014

Academic Editor: Yumin Cheng

Copyright © 2014 Ding Jun et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

A novel numerical manifold method was derived from the cubic B-spline basis function. The new interpolation function is characterized by high-order coordination at the boundary of a manifold element. The linear elastic-dynamic equation used to solve the bending vibration of thin plates was derived according to the principle of minimum instantaneous potential energy. The method for the initialization of the dynamic equation and its solution process were provided. Moreover, the analysis showed that the calculated stiffness matrix exhibited favorable performance. Numerical results showed that the generalized degrees of freedom were significantly fewer and that the calculation accuracy was higher for the manifold method than for the conventional finite element method.

1. Introduction

According to the concepts of topological manifold and differentiable manifold and based on the discontinuous dynamics of DDA block systems, the numerical manifold method (NMM) adopts the finite cover technique to set up a unified calculation form for the finite element method, the DDA, and the analytical method [1]. NMM combines the advantages of the finite element method and DDA and can deal with various types of complex geometric and material boundaries. Thus, NMM has an obvious advantage in the simulation of large geometric deformations and joint deformations, as well as dynamic and static crossings [2]. The application research on NMM mainly focuses on certain fields, such as geotechnical engineering [3, 4] and crack propagation [5–8]. However, NMM has recently been successfully applied to fatigue failure [9], fluid calculations [10], and seepage [11]. In theoretical research, most NMMs adopt shape functions in the finite element method as the weight functions for the manifold method and use the local cover functions of high-order polynomials to improve the continuity and interpolation accuracy of elements [12–14]. Notably, a polynomial form based on the global coordinate system will increase the calculated condition number of the stiffness matrix or even cause a solution failure. For this reason, Lin et al. [15] improved the form of local cover functions to enhance the computational efficiency and interpolation accuracy of elements. However, the number of generalized degrees of freedom remains large. The mesh-less numerical manifold method is more flexible in the selection of mathematical meshes and can avoid problems resulting from mesh generation [16]. However, plenty of matching points have to be used in the calculation, and matrix inversions have to be intensively performed. Thus, the number of calculations remains large, but the integral accuracy is low, which restricts the application [17, 18]. Recently, Wen and Luo [19] studied the numerical manifold method for polygonal elements, which has a stronger adaptability to physical meshes. The newly developed numerical manifold method based on quadratic B-spline interpolation can further improve the continuity of elements and prevent the linear correlation problem in the conventional numerical manifold method [20, 21].

The discrete finite element solution is the most widely adopted solution for structural dynamics. The interpolation function required by the calculation is based on the static node structure, which is inconsistent with actual problems. The interpolation shape function is related to the element shape and the coordinate system, such that the mesh requirements are relatively high, and the interpolation accuracy is obviously reduced as the structure changes significantly. For the finite elements with mesh structures, the shape functions usually have a low order, which indicates that they cannot meet the requirements of the continuity of displacement by structural deformation. NMM based on finite cover can solve the above problems effectively. The local cover function corresponding to the node displacement of the finite element is no longer relevant to the nodes. If the order of the local cover function is increased, the continuity and interpolation accuracy of the elements can easily be improved. In particular, the manifold method can adopt fixed mathematical meshes that are independent of the physical meshes (material boundary) in setting up the unified interpolation scheme. This method has the advantages of the Lagrangian and Eulerian methods. The manifold method can also avoid problems such as mesh distortion, moving boundary, and migration terms that are difficult to deal with in the finite element and conventional numerical manifold methods [22]. Based on the concept of fixed meshes and considering the deficiency of the existing manifold method in terms of interpolation accuracy, continuity, and solution efficiency, a manifold element with a new cover form and interpolation is constructed. To apply the new element to the solution of structural dynamics, the dynamic theory based on the numerical manifold method is investigated. The numerical manifold scheme for the bending vibration of thin plates is derived, and the final calculation results are comparatively analyzed.

2. Basic Principles of NMM

2.1. Basic Concept of NMM

NMM consists of block dynamics, simplex integration, and a finite cover system. Figure 1 shows the concrete structure [23]. NMM adopts two sets of independent meshes: a mathematical mesh (cover) and a physical mesh (cover). The mathematical mesh is used to define the interpolation accuracy of the solution based on user-defined forms, and the selected mathematical covers can overlap with one another. Consistency with the physical cover is unnecessary, but the mathematical cover has to include the whole material under study. The physical cover is used to describe the geometric properties of the material and the domain of integration, such as the material boundaries, cracks, blocks, and the interfaces of different materials.

Figure 1: Basic structure of the NMM.

In the finite cover system, the intersection of the mathematical cover and the physical mesh is defined as a manifold element; that is, the mathematical cover is subdivided by a physical mesh. The weight function in the manifold method is equivalent to a shape function in the finite element method. If the proper weight function is selected, the weighted summation of the local cover functions can generate the overall cover function of the whole solution domain.

2.2. Generation of Overall Cover Function

For simplicity, the manifold method adopts a mathematical cover with a regular structure or a finite element mesh to generate the desired weight function. The generated weight function should meet the property of element decomposition. Meanwhile, the local cover function corresponding to each mathematical cover is not directly relevant to the nodes. This function can be a constant, linear, and high-order polynomial or take other local series forms. The two-dimensional problem can be taken as an example. The local cover function of any mathematical cover is defined as follows:

The local cover function can be expressed as the linear combination of the basis function with a given -order and the corresponding unknown coefficients and . The matrix form is as follows:

Any manifold element is formed by the intersection of the cover ) subdivided by a physical cover. The weight function corresponding to each cover is used to solve the overall cover function in the domain, which can be represented as follows:

3. Two-Dimensional Manifold Method for Cubic B-Spline Interpolation

3.1. Definition and Property of the Cubic B-Spline Basis Function

The B-spline basis function can be defined in many ways to facilitate computer programming. In this paper, the reverse derivation method is used to define it [24]. is assumed to be a monotonically nondecreasing sequence; that is, and , where is the B-spline node. is the th B-spline curve with the following definition:

where is the node interval of th B-spline. During the process of reverse derivation, is required.

B-spline has a local support feature. In other words, when , , when , . All B-spline curves in the B-spline node interval meet the property of the element decomposition [19]; that is,

According to the definition of the B-spline basis function, the mathematical expression of the cubic B-spline curve expressed in piecewise can be derived as

Figure 2 shows the cubic B-spline curves when the B-spline nodes are evenly distributed. Formula (6) can be used to solve all the B-spline curves in the figure. Figure 3 shows the B-spline curves within the B-spline interval or within the element. The mathematical expression is as follows:

Figure 2: Cubic B-spline curve in piecewise form.

Figure 3: Cubic B-spline curve within the element .

3.2. Two-Dimensional Manifold Element of the Cubic B-Spline

A method similar to bilinear interpolation is adopted to construct the B-spline overall interpolation function (as shown in Formula (8)) within the rectangular element domain of and :

In the formula, is the generalized displacement vector of the manifold element , such that

Each term of is , such that

In the same way, cubic B-spline curves of , , , and within can be obtained in piecewise form. The calculated B-spline functions are substituted into Formula (8) to obtain the final overall cover function of the manifold element.

We suppose that the numbers of manifold elements in the directions of , are, respectively, , , and the corresponding element nodes are . Through Formulas (7), a certain number of virtual nodes have to be added to construct the B-spline basis functions. We suppose that the added virtual nodes are , , , ; , , , . The selected nodes only need to have the monotonically nondecreasing property.

Based on the continuity of the B-spline surface, the spline surface formed by uniform B-spline basis functions has a -order continuity [25] of , in the element boundary. The above manifold element adopts the model of the generalized node distribution, which is the same as that on the B-spline surface, such that second-order continuity exists in the element boundary. Figure 3 shows that the B-spline has a desirable property of local support, and its values are small within the interval to ensure the reasonability of the calculated condition number of the stiffness matrix. The number of generalized degrees of freedom required by the calculation is , which is significantly reduced compared with that of the conventional manifold method. Figure 2 shows that the B-spline basis functions have no property function in the nodes. Thus, these functions can reflect the dynamic property of the structure.

4. Manifold Scheme of the Bending Vibration of Simply Supported Thin Plates

4.1. Dynamic Equation of the Manifold Element Based on the Principle of Minimum Instantaneous Potential Energy

Suppose that a thin plate has a thickness of , a density of , an elasticity modulus of , a bending deflection (displacement) of , an actuating force in the unit area on the plate surface of , a normal direction for the simply supported boundary of , and corresponding fixed deflection and bending moments of and . According to the Kirchhoff plate theory [26] and based on the processing method of the boundary condition of the manifold element [1], the principle of the minimum instantaneous potential energy of the numerical manifold of the bending vibration of thin plates can be derived as

where is the area of the thin plate, is the simply supported boundary, is the damping coefficient, and and are the stiffness coefficients of the penalty function. is the bending stiffness matrix, and is the curvature and torsion of the elastic surface represented by vectors in the directions of , that meets the following equations:

in Formula (11) can be expressed as the following equation by the displacement function, where is the angle from the direction to the direction:

The manifold element of Formula (8) is substituted into Formula (11) along with Formulas (12) and (14) to obtain the discretization form of the minimum instantaneous potential energy:

The first variation is taken from Formula (15) that is the instantaneous variation that satisfies and . Therefore, the dynamic equation for any manifold element is given by

where

To facilitate the integration of the stiffness matrix, the transformation matrix is introduced and can be defined by the following formula, where is the total freedom of degree vector:

The transformation matrix corresponding to Formula (8) is

where is the identity matrix.

Formula (18) is used to integrate Formula (16) into the overall dynamic equation

where

It can be easily observed that (20) has an analogous form as the overall dynamic equation generated by FEM. Meanwhile, considering that and are symmetric positive matrices, we can solve the natural frequencies of the system by use of the techniques employed in FEM.

4.2. Initialization of the Dynamic Equation

Unlike the finite element method, , , and in the manifold method are no longer relevant to the nodes and thus need to be initialized. If the manifold element to be calculated is and the number of degrees of freedom is , the noncoincident points need to be regenerated in the directions of , to obtain the final nodes , , and . If the manifold element corresponding to node is and the real initial displacement is initialized with the velocity of and , we obtain

where

Formulas (23) can be used to obtain the generalized displacement and the generalized velocity , which are substituted into Formula (20) to obtain the generalized initial acceleration:

4.3. Integration Scheme of the Dynamic Equation

When the domain of area integration is a regular rectangular element, Gaussian integration can be adopted. When the domain is an irregular geometric configuration, it can be subdivided into several triangular elements. Then, the simplex integration or Hammer integration is used to obtain the calculation results with sufficient accuracy [27]. When the boundary integration is a regular curve, Gaussian integration can be directly performed. If the boundary integration is a complex curve, it should be subdivided into several, straight-line segments before Gaussian integration is performed.

For time integration, the Wilson- method can be selected in this paper. Formula (20) is then converted into the linear equation T at time . Consider

where

After solving , the generalized total displacement, velocity, and acceleration are solved:

Formulas (18) and (8) are used to obtain the function of the displacement, velocity, and acceleration of any manifold element at time :

5. Calculation Examples and Numerical Analysis

5.1. Calculation Example

As shown in Figure 4, the thin plate that is simply supported on the four edges has a length of m, a width of  m, and a thickness of  m. The transverse dynamic load is  Pa,  Pa, and . In addition, the elasticity modulus is  GPa, the Poisson’s ratio is , the stiffness coefficient of the penalty function is , the material density is  kg/m³, the damping coefficient is , and the coefficient of the Wilson- method is . Different time steps are selected for the calculation. As shown in Figure 5, the fixed rectangular mesh is used as the mathematical cover. The shaded area in the figure denotes any calculated manifold element.

Figure 4: Scheme and configuration of the thin plate.

Figure 5: Cover system and manifold element of the thin plate.

With the theory illustrated in Section 4.1, the generalized instantaneous potential energy for this example can be expressed as

where the latter four terms are the generalized energy form of boundary conditions.

After a large amount of discretization calculation, the matrices , , , and can be calculated by

where and are given in (12) and (13), respectively.

The frequency and the displacement for this case are as follows:

5.2. Numerical Analysis

Table 1 shows the frequency of the manifold method, the accuracy of which will improve along with the increasing density of meshes. The relative error of frequency in the table is defined as . When the mesh is adopted, the first-order frequency error is only 0.0041%, whereas the fifth-order frequency error is 0.24%.

Table 1: First five natural frequencies and their relative errors of manifold calculation.

Element
First frequency
Second frequency
Third frequency
Fourth frequency
Fifth frequency

/Hz

/Hz

/Hz

/Hz

/Hz

17.2564
0.084
22.5501
0.60
65.7306
111.79
70.9070
64.50

17.2427
0.0046
22.4204
0.026
31.2165
0.58
44.9846
4.36
63.9337
9.06

17.2413
0.0035
22.4141
0.0022
31.0619
0.085
43.3865
0.65
60.3066
2.87

17.2411
0.0046
22.4130
0.0071
31.0406
0.016
43.1758
0.16
59.0368
0.71

17.2412
0.0041
22.4126
0.0089
31.0348
0.0026
43.1275
0.052
58.7638
0.24

Theoretical solution
17.2419
0.00
22.4146
0.00
31.0356
0.00
43.1049
0.00
58.6227
0.00

, , and are defined as the numerical displacement, velocity, and acceleration, respectively. In the same way, , , and are the partial derivatives of numerical displacement. The corresponding relative errors have the same definitions as in Table 1. Table 2 shows the velocities and acceleration of the given points at certain time point. The accuracy of the numerical solution is high. Considering that the given interpolation functions has the same continuity and coordination for and , Table 3 only shows part of the partial derivatives of displacement. The numerical solution shows that both stress and strain in the manifold method have good accuracy.

Table 2: Velocity and acceleration of the manifold calculation and their relative errors.

Elements
, ,  s,

(×10⁻¹ m/s)

(m/s²)

6.2981
0.57
−33.705
1.81

6.2088
1.98
−33.693
1.85

6.2127
1.92
−34.199
0.38

6.2228
1.76
−34.437
0.32

Theoretical solution
6.3340
0.00
−34.328
0.00

Table 3: Manifold results of the displacement partial derivatives and their relative errors.

Element
, ,  s,

()

()

4.2774
7.85

10.17

3.9085
1.45

0.66

3.9712
0.13

0.71

3.9721
0.15

0.98

Theoretical solution
3.9660
0.00

0.00

Figures 6, 7, and 8 are, respectively, the displacements, velocities, and acceleration at the points of and for each moment. With the reduction in time step , the accuracy of the numerical solution improves. In particular, with the increase in time steps, smaller step lengths of can significantly improve the accuracy of the numerical solution. Figure 9 shows some partial derivations at the points of and . The figure presents the calculated values of the first several periods. As the meshes become denser, the calculation accuracy significantly improves. Therefore, the numerical solution is consistent with the theoretical solution.

Figure 6: Numerical result of displacement   with manifold elements.

Figure 7: Numerical result of velocity with manifold elements.

Figure 8: Numerical result of acceleration with manifold elements.

Figure 9: Numerical result of with  s.

Since the manifold method in this paper adopts an overall interpolation approximate, the integrated error norm has to be introduced to reflect the interpolation accuracy fully. is defined as the integrated error of displacement, where is the number of discrete points equally distributed in the and directions within the thin plate. During calculation, is taken as discrete points within the thin plate. The integrated errors of the other quantities in Tables 4 and 5 are defined by a similar method. Table 4 shows the bilinear Hermite finite element. The mathematical software MATLAB is strongly recommended for the calculation, similar to the manifold method. By comparing Tables 4 and 5, the overall degrees of freedom for the manifold method involving the same number of elements are significantly smaller than those for the finite element method, but its calculation accuracy is higher than that of the latter. Furthermore, Table 6 gives the global errors from the cubic B-spline finite element. By comparing Tables 6 and 5, we can safely conclude the cubic B-spline manifold element produces more accurate results than cubic B-spline finite element. All numerical results confirm the validity of the manifold method in solving structural dynamics problems.

Table 4: Global errors of the finite method with Hermite element.

Time
, , DOFs = 256,  s

20
1.22
0.66
3.21
23.26
9.59
1.52

40
0.34
0.83
2.54
12.50
3.27
0.74

60
0.47
0.47
4.11
6.67
1.41
0.88

80
0.42
0.53
2.49
10.97
2.83
0.78

100
0.50
1.02
2.26
10.25
2.50
0.77

Table 5: Global errors of the manifold method.

Time
, , ,  s

20
0.55
0.51
3.18
7.51
7.02
0.99

40
0.20
0.83
2.14
4.09
2.47
0.70

60
0.46
0.41
4.12
4.23
1.14
0.84

80
0.36
0.51
1.98
4.28
2.22
0.75

100
0.45
0.94
2.23
2.66
1.87
0.74

Table 6: Global errors of the finite method with Cubic B-spline element.

Time
, , DOFs = 256,  s

20
1.41
0.77
3.18
12.46
9.59
1.32

40
0.37
0.88
2.14
8.37
3.27
0.72

60
0.58
0.55
4.12
4.69
1.41
0.87

80
0.49
0.56
1.98
6.85
2.83
0.75

100
0.61
1.24
2.23
7.62
2.50
0.76

6. Conclusions

In this paper, based on the nonuniform cubic B-spline basis function, a new numerical manifold method was derived. By using the same structural model as the B-spline surface, the second-order coordinate can be obtained at the boundary of the manifold element. The local support property of the B-spline basis function enables the good performance of the stiffness matrix, which facilitates the numerical calculation. The new method was used to solve the bending problems of thin plates. The linear elasticity dynamic equation was derived based on the principle of the minimum instantaneous potential energy. The method for the initialization of the dynamic equation and the time integration method were also derived. The fixed mathematical mesh was used to implement the manifold method. The numerical results agreed well with the theoretical solution. Compared with the finite element method, the overall degrees of freedom of the cover calculated with the new manifold method were significantly reduced. The calculation accuracy was obviously improved, which indicates the superiority of the manifold method in solving structural dynamics problems. The theoretical and application studies in this paper are preliminary, confined to linear elasticity problems, and limited to situations in which only the time integration method is adopted. Further research on nonlinear problems should be conducted in the future to perfect the theory.

Conflict of Interests

The authors declare that they have no conflict of interests.

Acknowledgments

This work is financially supported by the Natural Science Foundation of China (11302272), the Natural Science Foundation of China (11102169), and the Science and Technology Project Affiliated with the Education Department of Chongqing Municipality (KJ120801).

References

J. M. Pei, “Numerical manifold method and discontinuous deformation analysis,” Chinese Journal of Rock Mechanics and Engineering, vol. 16, no. 3, pp. 279–292, 1997.

G. Ma, X. An, and L. He, “The numerical manifold method: a review,” International Journal of Computational Methods, vol. 7, no. 1, pp. 1–32, 2010.

G. Zhang, Y. Zhao, and X. Peng, “Simulation of toppling failure of rock slope by numerical manifold method,” International Journal of Computational Methods, vol. 7, no. 1, pp. 167–189, 2010.

Y. J. Ning, X. M. An, and G. W. Ma, “Footwall slope stability analysis with the numerical manifold method,” International Journal of Rock Mechanics and Mining Sciences, vol. 48, no. 6, pp. 964–975, 2011.

H. H. Zhang, L. X. Li, X. M. An, and G. W. Ma, “Numerical analysis of 2-D crack propagation problems using the numerical manifold method,” Engineering Analysis with Boundary Elements, vol. 34, no. 1, pp. 41–50, 2010.

S. C. Li, S. C. Li, and Y. M. Cheng, “Enriched meshless manifold method for two-dimensional crack modeling,” Theoretical and Applied Fracture Mechanics, vol. 44, no. 3, pp. 234–248, 2005.

S. C. Li, S. C. Li, J. W. Zhang, and Z. Q. Wang, “Numerical manifold method from mathematical theory and its application,” Engineering Mechanics, vol. 24, no. 6, pp. 36–42, 2007.

H. Gao and Y. Cheng, “A complex variable meshless manifold method for fracture problems,” International Journal of Computational Methods, vol. 7, no. 1, pp. 55–81, 2010.

G. Wei, K. Li, and H. Jiang, “Incompatible numerical manifold method for fracture problems,” Acta Mechanica Sinica, vol. 26, no. 2, pp. 247–255, 2010.

Z. Zhang, X. Zhang, and J. Yan, “Manifold method coupled velocity and pressure for Navier-Stokes equations and direct numerical solution of unsteady incompressible viscous flow,” Computers and Fluids, vol. 39, no. 8, pp. 1353–1365, 2010.

Q. H. Jiang, S. S. Deng, C. B. Zhou, and W. B. Lu, “Modeling unconfined seepage flow using three-dimensional numerical manifold method,” Journal of Hydrodynamics, vol. 22, no. 4, pp. 554–561, 2010.

J. Song and Y. Ohnishi, “High order rectangular element of manifold method,” Chinese Journal of Rock Mechanics and Engineering, vol. 22, no. 6, pp. 932–936, 2003.

S. Li, Y. Cheng, and Y.-F. Wu, “Numerical manifold method based on the method of weighted residuals,” Computational Mechanics, vol. 35, no. 6, pp. 470–480, 2005.

S. C. Li, Y. M. Cheng, and S. C. Li, “Meshless manifold method for dynamic fracture mechanics,” Acta Physica Sinica, vol. 55, no. 9, pp. 4760–4766, 2006.

S. Z. Lin, Y. F. Qi, and H. D. Su, “Improved local function of numerical manifold method and its application,” Journal of Yangtze River Scientific Research Institute, vol. 23, no. 6, pp. 55–58, 2006.

H. F. Gao and Y. M. Cheng, “Complex variable numerical manifold method for elasticity,” Acta Mechanica Sinica, vol. 41, no. 4, pp. 480–488, 2009.

S. C. Li and Y. M. Cheng, “Meshless numerical manifold method based on unity partition,” Acta Mechanica Sinica, vol. 36, no. 4, pp. 496–500, 2004.

H. Gao, Y. Cheng, and H. Jiang, “A meshless manifold method with complex variables for elasticity,” Chinese Journal of Applied Mechanics, vol. 27, no. 1, pp. 15–19, 2010.

W. B. Wen and S. M. Luo, “Polygonal manifold element for thin plate-bending analysis,” Engineering Mechanics, vol. 29, no. 10, pp. 249–256, 2012.

W. B. Wen, K. L. Jian, and S. M. Luo, “2D numerical manifold method based on quartic uniform B-spline interpolation and its application in thin plate bending,” Applied Mathematics and Mechanics, vol. 34, no. 8, pp. 1017–1030, 2013.

C. X. Guo and H. Zeng, “Study on linear dependence problem in high-order numerical manifold method,” Engineering Mechanics, vol. 29, no. 12, pp. 228–232, 2012.

H. D. Su, “Study on numerical manifold method with fixed meshes,” Chinese Journal of Theoretical and Applied Mechanics, vol. 43, no. 1, pp. 169–178, 2011.

S. C. Li and Y. M. Cheng, “Numerical manifold method and its applications in rock mechanics,” Advances in Mechanics, vol. 34, no. 4, pp. 446–454, 2004.

R. Sevilla, S. Fernández-Méndez, and A. Huerta, “NURBS-Enhanced Finite Element Method (NEFEM),” Archives of Computational Methods in Engineering, vol. 18, no. 4, pp. 441–484, 2011.

L. Piegl and W. Tiller, The NURBS Book, Springer, London, UK, 1996.

J. R. Barber, Elasticity Book, Kluwer Academic Publishers, New York, NY, USA, 2004.

S. M. Luo, X. W. Zhang, W. G. Lu, and D. R. Jiang, “Theoretical study of three-dimensional numerical manifold method,” Applied Mathematics and Mechanics, vol. 26, no. 9, pp. 1127–1132, 2005.