Collocation Discrete Least Squares (CDLS) Method for Elasticity Problems and Grid Irregularity Effect Assessment

A meshless approach, collocation discrete least square (CDLS) method, is extended in this paper, for solving elasticity problems and grid irregularity effect is assessed. In the present CDLS method, the problem domain is discretized by distributed field nodes. The field nodes are used to construct the trial functions. The moving least-squares interpolant is employed to construct the trial functions. Some collocation points that can be independent of the field nodes are used to form the total residuals of the problem. The least-squares technique is used to obtain the solution of the problem by minimizing the summation of the residuals for the collocation points. The final stiffness matrix is symmetric and therefore can be solved directly via efficient solvers. The boundary conditions are easily enforced by the penalty method. The present method does not require any mesh so it is a truly meshless method. Numerical examples are studied in detail, which show that the present method is stable and possesses good accuracy, high convergence rate and high efficiency for both regular and irregular point distribution.


INTRODUCTION
The finite element method (FEM) has been the most frequently used numerical method in engineering during the three past decades. It has been used in most fields of applied sciences such as computational solid mechanics [1, 2, and 3] and so on. Mesh-based methods are not well suited to the problems associated with extremely large deformation and problems associated with frequently remeshing. To avoid these drawbacks of the FEM, a new class of numerical methods, meshless methods (also called mesh-free methods) have been developing [4,5] in the recent decade. These methods have become an important tool in computational solid mechanics, owing to their advantages over the traditional finite element method (FEM), finite-volume method (FVM), and finitedifference method (FDM). Meshless methods rely only on a group of scatter points, which means not only that the burdensome work of mesh generation is avoided, but also more accurate description of irregular complex geometries can be achieved. Furthermore, the meshless approximation has higher smoothness, and no additional post-processing is needed.
In the field of meshless methods for solving elasticity problems, Krysl and Belytschko [6] employed Element-Free Galerkin Method (EFGM) to analyze thin plates; Onate et al. [7] proposed a stabilization technique by introducing new terms in both the governing equations and the traction boundary conditions to solve elasticity problems; Kwon et al. [8] presented a leastsquares meshfree method for solving linear elastic problems; Zhang et al. [9] proposed a meshless weighted least-squares (MWLS) method, to solve problems of elastostatics; Atluri et al. proposed a MLPG mixed collocation method [10] and MLPG mixed finite difference method [11] for solid mechanics.
All the above meshless methods can be categorized into two groups according to their discretization scheme. The first group is Galerkin-based meshless methods (GBMMs), of which the EFGM proposed by Belytschko in 1994 [12] is a famous representative. In GBMMs, the highest order of derivatives is lowered by using a weak form of the original partial differential equations (PDEs). The accuracy of GBMMs is high, and good stability can always be obtained. The main shortcoming of GBMMs is that the integrals in the weak form must be evaluated properly. One way of evaluating integrals is to use a background mesh, which makes the method not truly meshless; another is to use nodal integration [13] , which results in significant errors because the divergence theorem used in the establishment of the weak form demands accurate integration [14] . In addition, because meshless shape functions are too complex to be expressed in closed form, a delicate background mesh and a large number of quadrature points are always employed, which decreases the efficiency seriously. As a consequence, GBMMs are much more computationally expensive than the FEM [15] .
The other group of meshless methods is built on collocation schemes. The SPH, FPM, least-square collocation meshless method [16] , and point weighted least-square (PWLS) method [17] all belong to this group. These methods are very efficient and easy to program, but they usually suffer from poor stability, and the accuracy often goes down near the boundary.
The universal law of least squares can be used for discretization. In fact, it has been introduced into the FEM successfully [18] . A truly meshless method based on the least-squares approach, the collocation discrete least-squares (CDLS) method, was proposed to solve Poisson's equation [19] and free surface seepage problem [20] and also was presented for error estimation and adaptive refinement in one dimensional fluid mechanics [21] . Because of using the least-squares technique and more collocation points the CDLS method is not bothered by instability as collocation-based meshless methods. In this research the CDLS method is extended for elasticity problems and effect of grid irregularity is assessed.

MATERIALS AND METHODS
In the present CDLS method, the problem domain is discretized by distributed field nodes. The field nodes are used to construct the trial functions by employing the moving least-squares interpolant. Some collocation points that can be independent of the field nodes are used to form the total residuals of the problem. The least-squares technique is used to obtain the solution of the problem by minimizing the summation of the residuals for the collocation points.
Moving least square shape functions: Among the available meshless approximation schemes, the moving least squares (MLS) method [22] is generally considered to be one of the best methods to interpolate random data with a reasonable accuracy, because of its completeness, robustness and continuity [10,23] . With the MLS interpolation, the unknown function is approximated by: where P T (X) is a polynomial basis in the space coordinates, and m is the total number of the terms in the basis. For a 2D problem we can specify T 1 for 6. is the vector of coefficients and can be obtained by minimizing a weighted discrete L 2 norm as follows: The weight function is usually built in such a way that it takes a unit value in the vicinity of the point j where the function and its derivatives are to be computed and vanishes outside a region Ω surrounding the point and N is the number of nodes in the domain. In this research the cubic spline weight function is considered as follows: Comparing equation (4) with the well known form of equation (7) yields to equation (8) where T contains the shape functions of nodes at point which are called moving least square (MLS) shape functions.

Collocation discrete least square (CDLS) method:
Consider the following partial differential equation 0 Ω 9 Subject to appropriate Dirichlet and Neumann boundaries 0 10 0 11 where L and B are partial differential operators and f represents external forces or source term on the problem domain.
Upon discretization of the problem domain and its boundaries using Equation (7) defined as the residual of partial differential equation at a typical collocation point k is: , 1~ 12 The residual of Neumann boundary condition at typical collocation k on the Neumann boundary can also be written as: , 1~ 13 and finally the residual of Dirichlet boundary condition at collocations on the Dirichlet boundary could be stated by: The stiffness matrix K in Eq. (16) is square , symmetric and positive definite. Therefore, the final system of equations can be solved directly via efficient solvers.

RESULTS AND DISCUSSION
In this section, two 2D numerical examples, which are solved by the CDLS method, are presented. The examples include: 1) a cantilever beam under end point load, and 2) an infinite plate with a circular hole under uniaxial load which is solved with both regular and irregular node distribution. Cantilever beam: A cantilever beam under a point load at the end (see Fig. 1) is solved as the first example. For this problem, the exact stress and displacement solution for plane stress is given in Timoshenko and Goodier [24] Fig. 2. This problem is simulated using the MLS with the second order polynomial basis. Fig. 3 shows the vertical displacement along the central line of the beam for the three nodal configurations. The simulation results agree with the analytical solution very well. Fig. 4 shows the stress at stations along upper surface of the cantilever beam for the three nodal configurations and a very good agreement with the analytical solution is obtained. The contours of stress displayed over the deformed shape for the nodal configuration with d = 0.25 are shown in Fig. 5.  The convergence rate is studied with three nodal configurations (d =1.0, 0.5, and 0.25). The following error norm is used for showing the convergence rate in Fig. 6: where N is the total number of nodes, and are the approximation and exact stress values at node , respectively.
The results clearly show that a stable convergence rate is obtained for the present CDLS method.  respectively. In the above equations, G is the shear modulus and 3 1 ⁄ with v the Poisson's ratio. Due to symmetry, only the upper right square quadrant of the plate is modeled [see Fig. 8]. The edge length of the square is 5a, with a being the radius of the circular hole. The exact analytical displacements solution is imposed on the left and bottom edges, and the tractions obtained from the analytical solution [Eq. 22] are applied to the top and right edges. The periphery of the circular hole is traction-free.
The problem is solved using the CDLS method, under a plane stress condition, with the following constants: t =1, E = 1000, and v=0.3. For assessment of grid irregularity effect on results of CDLS method, three regular nodal configurations with 53, 183 and 721 nodes, respectively, and also three irregular nodal distributions with 49, 179 and 657 nodes, respectively, are used. The regular and irregular grids with 183 and 179 nodes, respectively, are shown in Fig. 8. The MLS with quadratic basis is used in the simulation. The horizontal displacement along the bottom edge (y = 0), and the stress component along the left edge (x = 0) obtained with the regular and irregular grids are shown in Fig. 9 and Fig. 10, respectively. The contours of stress obtained with the regular and irregular nodal configurations of 721 and 657 nodes are plotted in Fig. 11. Compared with the analytical solutions, good agreements are obtained for both the displacements and stresses. These results also show that the CDLS method possesses good accuracy for both regular and irregular grids. Fig. 12 shows convergence rate with three regular nodal configurations (n =53, 183, and 721) and three irregular nodal configurations (n =49, 179, and 657). A stable and monotonic convergence rate is observed for the problem. It is remarkable that similar convergence rates are obtained for both regular and irregular grids.

CONCLUSION
In this article, the collocation discrete least square (CDLS) method is extended to solve elasticity problems and grid irregularity effect is assessed. The final coefficient matrix obtained by the present method is symmetric and therefore can be solved directly via efficient solvers. The boundary conditions are easily enforced by penalty method. Because of no need to any mesh, this method is a truly meshless method. Results of the numerical examples show that the present method is stable and have good accuracy, high convergence rate and high efficiency for both regular and irregular grids. All of these advantages of the CDLS ensure that it is a very potential meshless method for the applications in elasticity problems.