MeshFree|Reissner-Mindlin板问题
ABSTRACT
Implicit boundary method enables the use of background mesh to perform finite element analysis while using solid models to represent the geometry. This approach has been used in the past to model 2D and 3D structures. Thin plate or shell-like structures are more challenging to model. In this paper, the implicit boundary method is shown to be effective for plate elements modeled using Reissner-Mindlin plate theory. This plate element uses a mixed formulation and discrete collocation of shear stress field to avoid shear locking.
The trial and test functions are constructed by utilizing approximate step functions such that the boundary conditions are guaranteed to be satisfied. The incompatibility of discrete collocation with implicit boundary approach is overcome by using irreducible weak form for computing the stiffness associated with essential boundary conditions. A family of Reissner-Mindlin plate elements is presented and evaluated in this paper using several benchmark problems to test their validity and robustness.
INTRODUCTION
The desire to avoid mesh generation difficulties has motivated much research in meshless, meshfree and mesh independent approaches for analysis. In this paper, we are interested in extending a mesh independent approach to model thin plate-like structures. A mesh independent approach uses standard solid models created in CAD software to represent the geometry of the analysis domain and uses a background mesh to approximate the solution. The key advantage is that the geometry does not have to be reconstructed as a mesh and the background mesh is easy to generate because it does not have to conform to the geometry. While this method has been demonstrated for 2D and 3D elements [1]-[3] it is more
challenging to extend a mesh independent approach to thin plate and shell-like structures. One approach to modeling shells is to use 3D background mesh with geometry being a surface passing through these elements. This approach was shown to work with cubic B-Spline elements [3] even for very thin shells. However, these elements are computationally far more expensive than the traditional lower order shell elements.This provides the motivation to seek mesh independent plate/shell elements that are based on the traditional plate and shell theories that yield acceptable solutions even with lower order elements. Another reason to formulate elements based on plate theory is to model composite shell-like structures where plate theory allows one to compute the effective properties when the number of layers and ply orientations are known.
In Reissner-Mindlin (RM) plate theory or first-order shear deformation theory of plates, the transverse shear strains are taken into consideration so that moderately thick plates can be accurately modeled. However, RM plate elements exhibit shear locking when the thickness of the plate tends to zero.
Shear locking is the term used to describe unrealistically high stiffness and lower deflections predicted for thin shells due to error in shear strain. Ideally, shear strain should tend to zero as the thickness of the plate is reduced. But the interpolations of displacement and slope are approximate due to which the computed shear strains do not reduce to zero as the plate thickness reduces resulting in significant shear strain energy and higher stiffness. Over the last three decades, extensive efforts to devise an effective approach to overcome shear locking have resulted in several approaches including reduced integration or selective reduced integration methods [19], [29],mixed/hybrid method [8], [21], assumed strain and discrete
collocation methods [10], [30], enhanced assumed strain method [17], [26], twist-Kirchhoff theory plate elements [15],[25], and a few other approaches. Bathe et al [11], [27], [28]developed a family of MITC elements for plates and shells,which were theoretically and numerically shown to have good performance.
Recent endeavors have focused on developing plate/shell elements using various new techniques, such as isogeometric method [14], [20], smoothed finite element method [23]-[24],discontinuous Galerkin finite element method [6][7], meshless or meshfree finite element method [13], [22] and implicit boundary finite element method [3]. Of particular interest are
techniques that are able to avoid the use of traditional mesh to model the geometry. The isogeometric method uses NURBS or T-Spline to represent both the geometry and the displacement fields. Meshless or MeshFree methods avoid the use of a mesh altogether by using interpolation or approximation schemes that do not need elements or connectivity between nodes. Yet another approach to avoid using a traditional mesh is the mesh independent approach which represents the geometry using equations of the boundary while using a background mesh to approximate or interpolate the displacement fields. The term “mesh independent” refers to the idea that the geometry is defined independent of the mesh and it was first used in the context of X-FEM for modeling a crack propagating through a continuum [12].
In this paper we investigate ways to extend the implicit boundary method to develop a mesh independent Reissner-Mindlin plate element, using mixed formulation and discrete collocation constraints [18], [30].Discrete collocation has been used to develop plate elements that satisfy the necessary conditions (or count conditions) for locking free elements. In this approach, the shear constraint is imposed only at the shear nodes and the shear stress field within an element is interpolated using the values at these nodes. In this paper, it is shown that even though the discrete collocation method appears incompatible with the implicit boundary method, it is possible to use these two methods in conjunction, if the irreducible weak form is used for computing the stiffness associated with the essential boundary conditions.
BASIC CONCEPTS AND DISCRETE COLLOCATION CONSTRAINTS
Plate element formulations using first order shear deformation plate theory and mixed formulation is now well established [9], [30]. In this section, we briefly describe the underlying concepts and establish a nomenclature that closely follows Zienkiewicz and Taylor [30]. Figure 1 illustrates the main assumption made in first order shear deformation plate theory which is that sections normal to the mid plane remain planar during deformation but do not necessary remain normal to the mid plane.
The weak form for Reissner-Mindlin plates can be derived from the governing equations either by using Galerkin’s weighted residual form or the principle of minimization of potential energy and stated as follows in a mixed form involving transverse deflection, rotations and shear stress fields:
and S contains the rotation and shear stress components respectively. 
is the bending rigidity matrix and 
where k is the shear correction factor, G is the shear modulus and h is the plate thickness. The applied loads include, q , the transverse distributed load, 
the shear force and 
the moment applied along the boundary. The rotation component 
is the rotation about the y-axis 
and 
is negative of rotation about x-axis 
. The moments and shear forces can be expressed as:
Where,
and the bending rigidity matrix for plates is
In the mixed formulation, the shear stress S is treated as one of the field variables to be computed but its components are interpolated within each element using different shape functions than the displacement and rotation components. For shear stress components, the shape functions must be selected such that the necessary conditions for stability and solvability are satisfied [30]. In the discrete collocation approach, the shear stress at the nodes are computed using equation (3), so that the nodal values of shear stresses is computed using the displacement and rotation fields and can be expressed in terms of the nodal values of deflection and rotations as:
The matrices
and
are constant matrices that depend on the shape functions used for interpolating deflection and rotations. The discrete collocation method was studied by Hinton and Huang [18] using patch test. Zienkiewicz et al [30] have provided count conditions as necessary conditions for the stability and solvability of the elements, and Bathe et al [9]-[10] have utilized inf-sup test to examine the stability and consistency of these types of elements.
IMPLICIT BOUNDARY METHOD FOR MINDLIN PLATES
Trial and test functions
For mesh independent analysis, a background mesh, as shown in Figure 2, is used to interpolate the trial and test functions. The geometry of the analysis domain is defined using equations of its boundaries that may be imported from CAD software. The background mesh may not have nodes along these boundaries therefore essential boundary conditions cannot be imposed by assigning nodal values as in traditional finite element mesh. In the implicit boundary approach,equations of the boundary are used to construct the trial and test function such that the essential boundary conditions(EBC) are guaranteed to be imposed.
Let (
) , where,
and
, be the implicit equation of a boundary (denoted as essential boundary in Figure 2(b)) on which the deflection w is prescribed. If only parametric equations of the boundaries of the analysis domain are available, the characteristic function
may be constructed as a signed distance function.Using this implicit equation, the trial and test functions for the deflection of a plate element may be constructed as:
is an approximate step function constructed such that 
and
The transition from 0 to 1 happens over a narrow band of width
referred to as the transition band. At the boundary,where
, the displacement field is equal to the prescribed boundary condition 
.Approximate step functions can be constructed to be polynomial or trigonometric in the transition region. A quadratic approximate step function H can be constructed as
In the step function,
is the implicit equation of the boundary or 
could be the signed distance function constructed using available parametric equations of the boundary. For elements that do not have an essential boundary passing through it,
. These fields are interpolated within elements using shape functions.
Here 
is a row matrix containing the shape function used for interpolating 
and 
.
is a column matrix containing the unknown nodal values of the field variable 
.
is a similar matrix containing the nodal values of the function a w which are prescribed such that at the boundary is equal to the applied displacement boundary condition
Where n is the number of nodes per element. Similarly,for the rotation components the trial and test functions for an element that contains an essential boundary can be defined as
In the above equations,
is a 
matrix containing the shape function used for interpolating 
and 
. The nodal values of and are stored in the column matrices 
and
.
In a mixed formulation, since the shear stress S is treated as nodal variables, its components are interpolated within the elements. A different set of shape functions and nodes are used for each component of the shear stress.
Where,
is a column matrix that contains the nodal values of the shear stress components and 
is a matrix of shape functions.
m is the number of nodes used for interpolating the shear components. As stated in equation (3), the relation between shear stress and displacements is 
. The shear stress should ideally go to zero as the plate thickness tends to zero. As a result at each integration point, the equilibrium equations are trying to enforce the constraint
as the plate thickness tends to zero. If this constraint is not possible to impose due to the approximation used for w and 
, the error due to non-zero shear leads to a shear locking.Often the number of integration points needed for accurate quadrature is too many and therefore it is not possible to accurately impose the shear stress constraints at all these points. In the discrete collocation approach, the shear stresses are computed at a set of points which are the nodal values for the shear stress interpolation. The number of these nodes has to be selected such that it is possible to make the shear stress zero at these points. The count conditions suggested in past work [30] are the necessary conditions for solvability and they serve as a useful guide in selecting the nodal points for shear interpolation. The shear stress at each collocation point is calculated by using equation (3) at these points.
Where,
is the parametric coordinate of the shear stress nodes. Substituting the trial and test functions for displacements, equation (7), (11) and (12), into equation (17),the shear stresses can be calculated as
where,
and
.Note that for all elements that do not have an essential boundary passing through it 
=1.
Therefore, for all elements which do not have an essential boundary we can set 
Similarly, we can write the variation of stresses as follows
More details on the selecting the collocation points or the location of shear nodes for various interpolation schemes can be found in Zienkiewicz et al [30]. For boundary elements that have essential boundaries, we cannot assume that 
everywhere and we need an alternate approach that is discussed in a later section.
Modified weak form
Using the trial and test functions in equations (7), (11)and (12) as well as the collocation relations from equations(19) and (20) into the weak form (equation (1)), a modified set of weak forms can be derived as:
Equations (21)and(22) can be converted to a matrix form as
Where,
Numerical integration
To evaluate the stiffness matrix and force vector, the integrals in equations (24)-(26) must be evaluated numerically. The stiffness matrix for internal elements is easy to compute by evaluating the following integral using Gauss quadrature.
For elements that are on the boundary, the same integral can be evaluated by first partitioning the region that is within the boundary into triangles and then integrating over each triangle. However, if the element has an essential boundary passing through it then an alternate approach is needed. The trial functions in equation (7) and (11) include a step function
that has very high gradients near the boundary. In the discrete collocation approach we compute the shear stress values using the displacement fields only at the nodes used for shear interpolation. Therefore, the sharp gradients associated with the step functions will be missed. To avoid this, we use the irreducible weak form [30] for computing the stiffness associated with the essential boundary conditions. In the irreducible weak form we replace shear stress with its definition 
to get
Using the equation (7) and (11), the left hand side of the above weak form can be converted into the standard form as
, where
Where,
Note that in 
we have set
and
because the approximate step functions are equal to one inside most of the domain and the narrow transition band, where the step function values change from 0 to 1 near the boundary, can be ignored when computing area integrals.
contains the derivative of the step function which have large values within the transition band but is zero outside this band. By making this band very narrow,
, it is possible to assume that the entire band is entirely within the element through which the boundary passes. The stiffness matrix computation can therefore be simplified as
Computation of 
is done by triangulating the boundary element and integrating only over triangles (whose area are 
,i=1,..,
) inside the geometry.
and
are computed by setting 
because
=0 in most the domain 
except in the transition band
along the boundary 
. Here we are assuming that 
is the signed distance function for the boundary. Using these definitions, the weak form can be rewritten in the standard discrete form as
NUMERICAL TEST
Several numerical benchmark problems are employed for the purpose of testing the implicit boundary method for the following elements: bilinear 4-node (Q4), serendipity 8-node(Q8), biquadratic 9-node (Q9) and bicubic 16-node (Q16). The results are compared with shell elements, S4R (4-node element) and S8R5 (8-node element), in ABAQUS [5] that use reduced integration technique to overcome shear locking. For simplicity, we utilize the following properties for all the examples: Young’s modulus,
, Poisson’s ratio 
=0.3 , a uniform pressure load p=100 . The thickness ofplate is denoted as ‘t’ and it is assigned values of 0.1 and 1 to test the element for different length/thickness (a/t) ratios.
Square Plate
Square plate has often been used in finite element analysis to assess the performance of plate elements because its analytical solution is available in literature. Both fullyclamped and simply-supported cases are studied here. The edge length of the square plate is 10.
In order to better verify the fact that the analysis results do not depend on the constructed background mesh, we utilize two kinds of background mesh pattern, as shown in Figure 3 and Figure 4. Both meshes are non-conforming but in Figure 4 the plate is rotated by 45 degrees so that the edges of the plate are not parallel to the element edges. Using identical boundary conditions for both, the results obtained from these two mesh using mesh independent elements converge to the analytical solution. The central deflections for various elements and mesh densities are studied in this example and the convergence of strain energy is illustrated.
Results obtained for the fully clamped and simplysupported square plate for both the thin plate (length to thickness ratio of 100) and the thick plate (length to thickness ratio of 10) are provided. The deflection distribution of clamped thin case for the two mesh patterns are shown in Figure 5. The distributions of stress (
) are shown in Figure 6. The convergence plots for the total strain energy are shown in Figure 7.
As the results indicate all the elements converge to the same answer. The four node elements Q4 and S4R are slow in converging to the correct answer compared to the higher order elements.
Circular Plate
A circular plate under uniformly distributed load with fully clamped and simply supported boundary conditions is considered here. For this geometry, the background mesh with uniform regular shaped elements has to be nonconforming and is therefore a good example to test any mesh independent method. The radius of the plate is 5. A typical background mesh is shown in Figure 8.
The deflection distribution of the clamped thick case is shown in Figure 9. The convergence plots for total strain energy are shown in Figure 10.
The performance of the higher order elements are similar to the ABAQUS elements while the Q4 element is very slow to converge in this case.
Flange plate
The earlier examples are benchmark problems used for testing plate elements. In this example a flange plate which has many curved boundaries is modeled and a non-zero essential boundary condition is applied on the central circular hole. The prescribed deflection boundary condition at the center hole is 
,
and
, where
for thick plate (thickness=1) and 
for thin plate (thickness =0.1). The plate is clamped at the four patterned holes such that
,and
.
The geometry of this plate is shown in Figure 11 and a typical background mesh is presented in Figure 12. The deflection of the plate is shown in Figure 13.The convergence plots for the total strain energy are shown in Figure 14.
For this example also all elements converge to the same result despite many curve boundaries along which the mesh independent approach uses partial elements. The boundary conditions are accurately enforced and the convergence rates are similar except for the Q4 element
CONCLUDING REMARKS
REFERENCES
[1] A.V. Kumar, S. Padmanabhan, R. Burla. Implicit boundary method for finite element analysis using nonconforming mesh or mesh, Int. J. Numer. Methods Engrg. 74(9) (2008) 1421-1447.
[2] Kumar, A. V., Burla R., Padmanabhan S. and Gu L.,2008, “Finite element analysis using non-conforming mesh”,Journal of Computing and Information Science in Engineering,Transactions of ASME, Vol. 8., No. 3, 031005(11 pages) doi:10.1115/1.2956990.
[3] R. Burla, A.V. Kumar, Implicit boundary method for analysis using uniform B-spline basis and structured mesh, Int.J. Numer. Methods Engrg. 76 (13) (2008) 1993-2028.
[4] A.V. Kumar, P.S. Periyasamy, Mesh independent analysis of shell-like structures, Int.J.Numer. Methods Engrg. (2012).In press.DOI: 10.1002/nme.4278.
[5] ABAQUS, 6.10 version.
[6] D.N. Arnold, F. Brezzi, L.D. Marini, A family of discontinuous Galerkin finite elements for the Reissner-Mindlin plate, J. Sci. Comput. 22 (2005) 25-45.
[7] D.N. Arnold, F. Brezzi, R.S. Falk, L.D Marini, Lockingfree Reissner-Mindlin elements without reduced integration,Comput. Methods Appl. Mech. Engrg. 196 (2007) 3660-3671.
[8] F. Auricchio, R.L. Taylor, A shear deformable plate element with an exact thin limit,Comput. Meth. Appl. Mech.Engng. 118 (1994) 393-412.
[9] K.J. Bathe, Finite Element Procedures in Engineering Analysis, Prentice-Hall, Englewood Cliffs, NJ. 1982.
[10] K.J. Bathe and E.N. Dvorkin, A four-node plate bending element based on Mindlin/Reissner plate theory and a mixed interpolation, Int. J. Numer. Methods.Engrg. 21 (2) (1985)367-383.
[11] K.J. Bathe, P.S. Lee, Measuring the convergence behavior of shell analysis schemes,Comput. Struct.89 (2011) 285-301.
[12] T. Belytschko and C. Parimi and N. Moës and N.Sukumar and S. Usui. Structured extended finite element methods for solids defined by implicit surfaces. Int. J. Numer.Methods Engrg. 56(4)(2003) 609-635.
[13] T. Belytschko, Y. Krongauz, D. Organ, M. Fleming, P.Krysl. Meshless methods: an overview and recent developments. Comput. Methods Appl. Mech. Engrg. 1996;139 : 3-47.
[14] D.J. Benson, Y. Bazilevs, M.C. Hsu, T.J.R. Hughes,Isogeometric shell analysis: The Reissner-Mindlin shell,Comput. Methods Appl. Mech. Engrg. 199 (2010) 276-289.
[15] F. Brezzi, J.A. Evans, T.J.R. Hughes, L.D. Marini, New rectangular plate elements based on twist-Kirchhoff theory,Comput. Methods Appl. Mech. Engrg. 200 (2011) 2547-2561.
[16]C. Cartensen, K. Weinberg, Adaptive mixed finite element method for Reissner-Mindlin plate, Comput. Methods Appl. Mech. Engrg, 190 (2001) 6895-6908.
[17] S. Glaser, F. Armero, On the formulation of enhanced strain finite elements in finite deformations, Engrg. Comput,14 (7) (1997) 759-791.
[18] E. Hinton, H.C. Huang, A family of quadrilateral Mindlin plate elements with substitute shear strain fields, Comput.Struct.23 (3) (1986) 409-431.
[19] T.J.R. Hughes, M. Cohen, M. Haroun, Reduced and selective integration techniques in finite element method of plates, Nuclear Engrg. Des. 46 (1978) 203-222.
[20] T.J.R. Hughes, J.A. Cottrell, Y. Bazilevs, Isogeometric analysis: CAD, finite elements,NURBS, exact geometry and mesh refinement, Comput. Methods Appl. Mech. Engrg. 194(2005) 4135-4195.
[21]D. Karamanlidis, S.N. Atluri, Mixed finite element models for plate bending analysis-Theory, Comput. Struct. 19(3) (1984) 431-445.
[22]P. Krysl, T. Belytschko, Analysis of thin plate by element-free Galerkin method, Comput. Mech. 17 (1995) 26-35.
[23] H. Nguyen-Xuan, T. Rabczuk, S. Bordas, J.F. Debongnie,A smoothed finite element method for plate analysis, Comput.Methods Appl. Mech. Engrg. 197 (2008) 1184-1203.
[24] H. Nguyen-Xuan, T. Rabczuk, N. Nguyen-Thanh, T.Nguyen-Thoi, S. Bordas, A node-based smoothed finite element method (NS-FEM) for analysis of Reissner-Mindlin plates, Comput. Mech. 46 (2010) 679-701.
[25] H.A.F.A. Santos, J.A. Evans, T.J.R. Hughes,Generalization of the twist-Kirchhoff theory of plate elements to arbitrary quadrilaterals and assessment of convergence,Comput. Methods Appl. Mech. Engrg. 209-212 (2012) 101-114.
[26] J.C. Simo, M.S. Rifai, A class of mixed assumed strain methods and the method of incompatible modes,Int. J. Numer.Methods Engrg. 29 (1990) 1595-1638.
[27] P.S. Lee, K.J. Bathe, Development of MITC isotropic triangular shell finite elements, Comput. Struct.82 (2004) 945-962.
[28] P.S. Lee, K.J. Bathe, The quadratic MITC plate and MITC shell elements in plate bending, Advance in Engrg.Software.41 (2010) 712-728.
[29]O.C. Zienkiewicz, R.L. Taylor, J.M. Too, Reduced integration technique in general analysis of plates and shells,Int. J. Numer. Methods Engrg. 3 (1971) 275-290.
[30] O.C. Zienkiewicz and R.L. Taylor, The Finite Element Method for Solid and Structural Mechanics, sixth ed.,Butterworth-Heinemann. 2005, NY.
Author
From
https://www.researchgate.net/publication/267490668
