A novel meshfree approach for free vibration and buckling analysis of thin laminated composite plates

A novel meshfree radial point interpolation approach which employs a new numerical integration scheme is introduced. The new integration scheme, namely Cartesian Transformation Method, transforms a domain integral into a double integral including a boundary integral and a one-dimensional integral, and thus allowing integration without discretizing domain into sub-domains usually called background mesh in traditional meshfree analysis. A new type of radial basis function that is little sensitive to user-defined parameters is also employed in the proposed approach. The present approach is applied to free vibration and buckling analysis of thin laminated composite plates using the classical Kirchhoff’s plate theory. Various numerical examples with different geometric shapes are considered to demonstrate the applicability and accuracy of the proposed method.


INTRODUCTION
inite element method (FEM) [1] is well-known in the engineering communities due to its advantages in solving partial differential equations. The method has many (advantages?) advatages, such as simplicity and high accuracy with not-sohigh computational cost. However, it is not The main idea of FEM is discretizing the problem domains into non-overlapping subdomains called elements. Each element usually has a common geometric shape such as triangle, quadrilateral (for two-dimensional domains), tetrahedron and hexahedron (for three-dimensional domains). A "good quality" element usually has to satisfy certain requirements such as size and convexity. In cases of large deformation, elements could be distorted and become a source of considerable error. Furthermore, in problems where the mesh has to be updated such as in moving-boundary problems and crack-propagation problems, re-meshing is always a challenging task.
The class of meshfree methods, on the other hand, does not require elements. The problem domain is represented only by nodes, including nodes on boundaries and nodes inside the domain [2]. Hence, the difficulties related to elements are avoided. Most of the meshfree methods are developed upon basis functions that do not possess Kronecker-delta property, requiring extra techniques such as Lagrangian multipliers and penalty method to enforce boundary conditions. In contrast, the Radial Point Interpolation method (RPIM) [3] satisfies the Kronecker-delta property, allowing direct imposition boundary conditions. Since the introduction, the method has been intensively investigated and applied to various engineering problems, such as structural dynamics [4], plate analysis [5], heat transfer [6], fracture mechanics [7] and unsaturated flow [8].
One drawback of the RPIM is the influence of user-defined parameters on numerical results, and it seems that each problem requires a distinct "optimum" set of parameters. Recently a new quartic radial basis function was introduced by [9], F in which the user-defined parameters are less sensitive to the numerical results. A second aspect is that during numerical integration, the problem domain still has to be divided into sub-domains called background cells. The situation could be improved by the novel integration scheme, namely Cartesian Transformation Method (CTM), proposed by [10]. In CTM scheme, the domain integral is transformed into a double integral including a boundary integral and a 1D integral. The double integral can be numerically evaluated without the creation of background cells as usual (usually?) seen in traditional meshfree methods.

BRIEF ON RPIM FORMULATION.
Consider a 2D elastic body Ω bounded by the boundary Γ = ∂Ω. A function u(x) defined in Ω can be approximated by The vector of polynomial basis p(x) is usually chosen as a complete second order polynomials  (5) R(xi, xj) is the radial basis function (RBF) and can be defined in many forms [2]. For example, the multiquadric form [3] (6) and the quartic form recently proposed by [9] In Eq. (6) and (7), rij is the distance between node i and node j; αc, q and θ are the user-chosen parameters. The parameter rs in Eq. (7) is the maximum distance between a pair of nodes in the support domain.

CARTESIAN TRANSFORMATION METHOD FOR EVALUATION OF TWO-DIMENSIONAL DOMAIN INTEGRALS.
The method was originally reported by [10] as an alternative numerical integration scheme to the popular Gaussian quadrature. The main idea of the scheme is to transform a domain integral into a double one-dimensional integral, hence it is named as Cartesian Transformation Method (CTM). Consider a domain integral defined over a domain Ω as follows , I f x y d (8) where f is an arbitrary regular function. Next, an auxiliary domain ΩR that contains the integration domain Ω is defined. The domain integral in Eq. , The double integral in Eq. (10) can then be easily evaluated by Gaussian composite scheme, as illustrated in Fig 2. The one-dimensional integral along the y-direction is first evaluated by dividing the vertical direction into k intervals. Within each interval, a certain number of Gauss points are selected. From each Gauss point on vertical direction, a horizontal ray is created. Again, each horizontal ray is divided into a certain number of intervals, and Gauss points are selected within each interval, so that the line integral in Eq. (11) can be evaluated. For an illustration of applying CTM into a specific problem, please refer to Fig 2. It is noted here that the number of intervals in both ydirection and x-direction directly relates to the number of integration points. In any numerical integration scheme, increasing the number of integration points will increase the accuracy of the evaluation, but computational time also increases. In the case of standard Gauss quadrature for integrands in form of polynomials, an optimum number of integration points can be determined, see [1]. Determination of "optimum" distribution of integration points for CTM scheme is an interesting topic but it is not in the scope of this paper and thus is scheduled for future research.

FREE VIBRATION AND BUCKLING
ANALYSIS OF THIN LAMINATED COMPOSITE PLATES.
Let us consider a thin laminated composite plate, as depicted in Fig 3, showing the fiber orientation of a layer denoted by φ. The displacements of the plate in the x-, y-and zdirection are denoted as u, v and w, respectively. Following the Kirchhoff theory for thin plates, the displacement fields can be defined as  (13) with D being the material stiffness matrix. Details on determination of the matrix D for thin laminated composite plate can be found in [10].

Free vibration analysis
Based on Kirchhoff theory for thin plate, the deflection of the plates can be approximated from the nodal deflection, wI, in the following form The constrained Galerkin weak formulation for undamped elasto-dynamic problems of thin plates without consideration of external force can be written as: where μ1, μ2 are defined as ratio between the loads μ1= Nyy/Nxx and μ2= Nxy/Nxx

NUMERICAL EXAMPLES.
To investigate the applicability and accuracy of the proposed method on free vibration and buckling analyses of thin laminated composite plates, three numerical examples with different geometrical shape are considered. In order to demonstrate the efficiency of the novel techniques, i.e. the quartic radial basis function and the CTM integration scheme, only symmetric configuration of laminated composite plates is considered for simplicity. The essential boundary conditions are restricted to simply supported on all external boundaries, as rotations are not included in the variables. Constraints related to rotations, such as an edge being clamped, may be treated by the suggestion in [12], but it is not within the scope of the present work. In all examples, the term "standard RPIM" denotes the RPIM that employ multiquadric basis function with parameters q = 1.03 and αc = 1, and the standard Gaussian quadrature for numerical integration.

Free vibration analysis of a laminated composite elliptical plate
A laminated composite plate in elliptical shape is considered in this example. The major radius and minor radius of the elliptical plate are a = 5 m and b = 2.5 m, respectively. Other geometrical and material parameters are given by: thickness t = 0.06 m, mass density ρ = 8000 kg/m 3 , ratio of elastic constants E1/ E2 = 2.45 and G12/E2 = 0.48, Poisson's ratios ν12 = 0.23 and ν21 = ν12.E2/ E1. The natural frequencies are normalized by in which  Table 1. The value of userdefined parameter θ, regardless as small as 1 or as big as 10000, seems do not affect the numerical results. This is indeed an advantage of the present method, compared with the Moving Kriging interpolation [5], where the correct results depend heavily on the "right choice" of user-defined parameter. Further observation reveals that the computational time for the present method is close to that for the standard RPIM, as shown in Table 2. As the time needed to compute the quartic basis function is not more than the multiquadrics basis function, it could be inferred from Table 2 that the CTM integration scheme is equivalent to or even faster than standard Gaussian quadrature. Accuracy could also be assumed as equivalent due to good agreement between the different methods, as shown in Table 1. However, it should be noted that in standard Gaussian quadrature, a system of background cells has to be created beforehand, which can be considered as a kind of "mesh" and thus is not favored in application of meshfree analysis. On the other hand, the CTM scheme requires no background cell.  Table 3 presents a comparison results obtained by present method with other meshfree methods, where good agreement can be observed. Further investigation on computational time again shows that CTM scheme is potentially faster than Gaussian scheme in evaluation of numerical integration, see Table 4.

Buckling analysis of a plate with a hole of complicated shape
In the last example, buckling analysis buckling for a plate with a hole of a complicated shape is investigated, see Fig 6. Plate thickness is t = 0.06. The material properties are the same as that mentioned in Example 5.1. The plate is simply supported on all four edges and loading conditions are similar to Example 5.2. Application of the CTM integration scheme for this problem was illustrated in Fig 7. The procedure of CTM is presented in Section 3 and will not be repeated. Given approximately equivalent number of integration points (2040 points for CTM and 1888 points for Gauss quadrature), computational time required by CTM is less than 10s, while that by standard Gauss quadrature is more than 11s. The dimensionless critical buckling load factor for various configurations of fiber orientation is presented in Table 5. Analytical solution for this problem is not available, therefore results calculated by finite element methods with a fine mesh of 3148 elements (3406 nodes) are taken as reference.  The quartic radial basis function has been used in the meshfree Radial Point Interpolation Method to develop an approach for free vibration and buckling analysis of thin laminated composite plate. Obtained results does not heavily depend on the user-chosen parameter θ. Hence, it is reasonable to select a default value θ = 1. Insensitivity to user-defined parameters would be a desirable property that broadens the applicability of meshfree method in practical problems. For both types of analyses considered in this paper, good agreement between results obtained by the proposed method and other results reported in literature. The second remark is the employment of CTM integration scheme. The scheme has been shown to be equivalent to the well-known Gaussian quadrature in accuracy. Collected data implies that the computational time in case of CTM scheme is potentially less than the Gaussian scheme. However, this is only preliminary observation o case of thin plate analysis. The hypothesis that CTM scheme is more efficient than Gauss quadrature in term of computational time shall be further investigated in future works. Nevertheless, it is worth noting that the CTM scheme is more practical than the Gauss quadrature, as it requires no background cells during numerical integration and thus is closer to the definition of "mesh free" methods.  Tóm tắt-Bài báo giới thiệu phương pháp không lưới mới sử dụng một kỹ thuật tích phân mới. Kỹ thuật này, với tên gọi Cartesian Transformation Method, biến đổi phép tích phân miền thành một phép tích phân biên và một phép tích phân một chiều, từ đó cho phép tính tích phân số mà không cần chia miền bài toán thành các ô tích phân, thường gọi là các ô nền trong phương pháp không lưới truyền thống. Cùng với đó, phương pháp đề xuất được tích hợp một dạng hàm nội suy hướng kính mới với đặc tính ít phụ thuộc vào các tham số tùy chọn. Phương pháp mới phát triển được ứng dụng vào phân tích dạng dao động riêng và bất ổn định tấm mỏng composite lớp theo lý thuyết tấm cổ điển Kirchhoff. Các ví dụ tính toán được phân tích và so sánh để làm rõ tính chính xác và hiệu quả của phương pháp.