A Novel RBF Collocation Method Using Fictitious Centre Nodes for Elasticity Problems

: The traditional radial basis function collocation method (RBFCM) has poor stability when solving two-dimensional elastic problems, and the numerical results are very sensitive to shape parameters, especially in solving elastic problems. In this paper, a novel radial basis function collocation method (RBFCM) using ﬁctitious centre nodes is applied to the elastic problem. The proposed RBFCM employs ﬁctitious centre nodes to interpolate the unknown coefﬁcients, and is much less sensitive to the shape parameter compared with the traditional RBFCM. The details of the shape parameters are discussed for the novel RBFCM in elastic problems. Elastic problems with and without analytical solutions are given to show the effectiveness of the improved RBFCM. 35Q68


Introduction
Over the past few decades, one type of the popular numerical methods that have been used for solving scientific and engineering problems has been the mesh-based methods, such as the finite element method (FEM) [1,2], finite difference method [3], finite volume method [4], and so on [5]. However, mesh-based methods face huge challenges and difficulties when dealing with large deformations. In recent years, meshless methods have been developed for their simplicity and effectiveness. Without grid division, the meshless methods have obvious advantages in solving complex domains, moving boundaries, and high-dimensional problems.
The Kansa method or RBFCM is one type of meshless method based on the collocation technique [6]. The uniqueness and convergence of the RBFCM have been investigated [7,8]. However, in the RBFCM simulation process, an asymmetric and fully populated matrix of a system of the linear equations is generated, which may cause a high condition number and affect the stability of the method. In order to deal with the asymmetric matrix, some methods have been proposed and developed to avoid the difficulties associated with the asymmetric matrix, such as the RBF Hermite collocation method [9] and the modified Kansa's method (MKM) [10]. In the traditional RBFCM, the centre and collocation nodes are the same. The matrix formulated in the RBFCM is ill-conditioned, and the traditional RBFCM is not stable, especially in dealing with Neumann boundary conditions [11]. The boundary conditions should be treated properly by using the weighted RBFCM and other numerical techniques [12,13].
The shape parameter of the traditional RBFCM should really be taken care of, and the accuracy of the traditional RBFCM is related to the shape parameters. Many numerical algorithms have been put forward to pick out the optimal shape parameters, such as the leave-one-out cross-validation (LOOCV) algorithm [14], the genetic algorithm [15], and where µ is Poisson's ratio and E is the elasticity modulus, u i (x, y), i = 1, 2 are the displacements, and f i (x, y), i = 1, 2 are the forcing terms. Ω is a domain bounded by a segmented smooth surface. The traction boundary conditions can be expressed as follows where n i (x, y), i = 1, 2 are the unit normal vectors at the boundary node. ∂Ω N is the boundary that satisfies traction boundary conditions, f i (x, y), i = 1, 2 are the tractions. The displacement boundary conditions can be given as g i (x, y), i = 1, 2 represent the displacement at boundary nodes. ∂Ω D is the boundary that satisfies Dirichlet boundary conditions. For the plane strain problem, the governing equation and boundary conditions can be obtained only by changing E to E/ 1 − µ 2 and µ to µ/(1 − µ) in the Equations (1)- (6).

Three-Dimensional Cases
The governing equations of the three-dimensional elasticity problem [28] can be expressed as where θ = ∂u 1 ∂x + ∂u 2 ∂y + ∂u 3 ∂z . µ is Poisson's ratio and E is the elasticity modulus. u i , i = 1, 2, 3 are the displacements. σ ij , i, j = 1, 2, 3 are the stresses. The governing equations of the 3D elasticity problem can be expressed as where f i , i = 1, 2, 3 are the forcing terms, Ω is a 3D computational domain. The traction boundary conditions can be expressed as follows where n i , i = 1, 2, 3 are the unit normal vectors at the boundary node. ∂Ω N is the boundary that satisfies traction boundary conditions, f i (x, y), i = 1, 2, 3 are the tractions. The displacement boundary conditions can be given as g i , i = 1, 2, 3 represent the displacement at boundary nodes. ∂Ω D is the boundary that satisfies Dirichlet boundary conditions.

The RBFCM
Here we take the 2D elastic problem as an example. In the RBFCM [6], the displacements of the 2D elastic problem can be approximated as Mathematics 2022, 10, 3711 4 of 15 where N is the number of centre nodes, ξ 1,j and ξ 2,j are the unknown coefficients related to the numerical solutions u n 1 and u n 2 , respectively; φ j is the RBF, in this work the multi-quadric (MQ) RBF is given as follows The RBF is associated with centre nodes (x j , y j ), which are the same as the collocation nodes (x k , y k ) in the traditional RBFCM, as shown in Figure 1. c is the shape parameter.
on ∂Ω D , and K = K i + K b1 + K b2 is the number of total nodes. In the traditional RBFCM K = N. By substituting Equations (7) and (8) back to the governing equations in Equations (1) and (2), the following discretized form can be obtained (1) (2) 1, , , is the index of the number of the inner nodes. Equations (13) and (14) can be obtained by substituting Equations (7) and (8) back to traction boundary conditions in Equations (3) and (4). By substituting Equations (7) and (8) back to the governing equations in Equations (1) and (2), the following discretized form can be obtained where where k = 1, . . . , K i , is the index of the number of the inner nodes. Equations (13) and (14) can be obtained by substituting Equations (7) and (8) back to traction boundary conditions in Equations (3) and (4).
Mathematics 2022, 10, 3711 j,k = n 1 µ ∂φ j (x k ,y k ) ∂y j,k = n 2 µ ∂φ j (x k ,y k ) ∂x j,k = n 2 ∂φ j (x k ,y k ) ∂y is the index of the number of the boundary nodes that satisfy the traction boundary conditions. Equations (16) and (17) can be obtained by substituting Equations (7) and (8) back to the displacement boundary conditions in Equations (5) and (6).
is the index of the number of the boundary nodes that satisfy the Dirichlet boundary conditions. Generally, we choose K ≥ N. The unknown coefficients ξ 1,n N n=1 and ξ 2,n N n=1 can be calculated from the linear system obtained from Equations (10)- (17). Then, the displacements of Equations (1)-(6) can be approximated by Equations (7) and (8). The collocation and the centre nodes are the same in the RBFCM, and the size of the system matrix is 2K × 2N, where K = N here.

The Improved RBFCM with Fictitious Centre Nodes
In the improved RBFCM, the collocation nodes and centre nodes are different. A set of ghost nodes are taken as fictitious centres, as the red nodes shown in Figure 2. The number of the fictitious nodes is defined as K. The fictitious centre nodes can be placed inside and outside the domain, as shown in Figure 2. The size of the fictitious centre nodes can be controlled by a radius of R. The relationship of the distance in the RBFCM is changed with such a simple adjustment. To study the performance of the improved RBFCM, the shape parameters and the values of radius R require further study. •" and fictitious centre nodes "○").

Shape Parameters
Over the past two decades, various algorithms have been proposed to predict the optimal shape parameter in the RBFCM. The shape parameter is vital for the performance of the RBFCM in different problems. In this work, three different approaches are considered for choosing shape parameters in the improved RBFCM for the elastic problem.
Approach 1: Brute force (BF) A typical method to find an optimal c is the "brute force" method, in which the value Figure 2. The node distribution with fictitious centre nodes (boundary nodes "•", interior nodes "•" and fictitious centre nodes " ").

Shape Parameters
Over the past two decades, various algorithms have been proposed to predict the optimal shape parameter in the RBFCM. The shape parameter is vital for the performance of the RBFCM in different problems. In this work, three different approaches are considered for choosing shape parameters in the improved RBFCM for the elastic problem.

Approach 1: Brute force (BF)
A typical method to find an optimal c is the "brute force" method, in which the value of the shape parameter starts from c = 0.01, increasing by 0.01 each time until c arrives at a certain value. Then an optimal shape parameter can be decided by using the trial errors obtained from different shape parameters. The BF costs a lot of time in the repeating process, and it is not an efficient way to find the optimal shape parameter. The time costs of the "brute force" method can be dramatically reduced by using the modified Franke formula.

Approach 2: The modified Franke formula (MFF)
In 1982, an experimental formula based on the density of the interpolation nodes was proposed by Franke to estimate the optimal shape parameter [21]. Recently, owing to the popularity of double precision arithmetic, the Franke formula has been revised to where D is the diameter of the smallest circle containing all fictitious centre nodes and N is the number of the collocation nodes. The modified Franke formula has been demonstrated to be a satisfactory prediction for a reasonably good shape parameter of the multi-quadric (MQ) radial basis function.
Approach 3: The sample solution approach (SSA) When the analytical solutions are unavailable, the sample solution approach (SSA) is employed to validate the numerical results [23]. In the procedure of the SSA, an exact solution to the pseudo-problem is set up as the sample solution. The pseudo-problem has the same geometry and the same number of degrees of freedom as the current problem. The optimal shape parameters of the pseudo-problem can be found with the help of the exact solutions. Then the optimal shape parameters obtained from the pseudo-problem are used to solve the problem where the exact solution is not available.

Numerical Results
In this section, an example with an exact solution is presented to illustrate the effectiveness of an improved RBFCM, and then two other numerical examples of two-dimensional elastic problems without exact solutions are presented by comparing them with the numerical results of the finite element method. The fictitious centre nodes are distributed by using the Halton quasi-random number generator, which is available via the MATLAB command haltonset. The relative error is defined as follows where u n and u are the numerical and analytical solutions of the displacement, respectively. Example 1. In this example, the plan strain problem described in Equations (1) The modulus of elasticity is taken to be E = 2.1 × 10 5 Mpa, and µ = 0.2. The computational domain with node distributions is shown in Figure 3, where a 2 m × 2 m square elastic plate with 100 circular holes (radius = 0.05 m) is presented. The traction boundary conditions are applied to the four boundaries of the square domain, and the displacement boundary conditions are imposed on the inner boundaries of 100 circular holes. The analytical solutions are given as tions is shown in Figure 3, where a 2 m × 2 m square elastic plate with 100 circular holes (radius = 0.05 m) is presented. The traction boundary conditions are applied to the four boundaries of the square domain, and the displacement boundary conditions are imposed on the inner boundaries of 100 circular holes. The analytical solutions are given as The node distribution of the improved RBFCM is shown in Figure 4, where 1296 fictitious centre nodes (red circle "○") are given to evaluate the numerical results. Unlike the traditional RBFCM, the fictitious centre nodes are placed in a larger circular area with R = 3. The numerical results of the improved RBFCM are compared with the results obtained with the traditional RBFCM in Figure 5. The node distribution of the improved RBFCM is shown in Figure 4, where 1296 fictitious centre nodes (red circle " ") are given to evaluate the numerical results. Unlike the traditional RBFCM, the fictitious centre nodes are placed in a larger circular area with R = 3. The numerical results of the improved RBFCM are compared with the results obtained with the traditional RBFCM in Figure 5.   In Figure 5, the relative errors of the displacements computed by the traditional and the improved RBFCM are shown with the red and blue colors, respectively. Different shape parameters are used to evaluate the relative errors of both traditional and improved RBFCM. The numerical results show that the numerical results of the improved RBFCM are always much better than the traditional RBFCM. Moreover, when the shape parameter changes, the relative errors of the improved RBFCM are much smaller than 10 −5 , which indicates that the improved RBFCM is not sensitive to the variations in the shape parameters.
In Figure 5, the three different approaches to choosing shape parameters are tested for the improved RBFCM. The c BF , c SSA , and c MFF are the optimal shape parameters calculated by the BF, SSA, and MFF, respectively. The relative errors of c BF , c SSA , and c MFF in Figure 5 are less than 10 −9 . The BF, SSA, and MFF can also be applied to obtain the optimal shape parameter in the improved RBFCM.   In Figure 5, the relative errors of the displacements computed by the traditional and the improved RBFCM are shown with the red and blue colors, respectively. Different shape parameters are used to evaluate the relative errors of both traditional and improved RBFCM. The numerical results show that the numerical results of the improved RBFCM are always much better than the traditional RBFCM. Moreover, when the shape parameter changes, the relative errors of the improved RBFCM are much smaller than 10 −5 , which indicates that the improved RBFCM is not sensitive to the variations in the shape parameters.
In Figure 5, the three different approaches to choosing shape parameters are tested for the improved RBFCM. The cBF, cSSA, and cMFF are the optimal shape parameters calculated by the BF, SSA, and MFF, respectively. The relative errors of cBF, cSSA, and cMFF in Figure 5 are less than 10 −9 . The BF, SSA, and MFF can also be applied to obtain the optimal shape parameter in the improved RBFCM.   As the improved RBFCM is influenced by the fictitious centre nodes, the influence of R is studied in Table 1.  As the improved RBFCM is influenced by the fictitious centre nodes, the influence of R is studied in Table 1. In Table 1, 700 interior nodes and 600 boundary nodes are considered, thus 1300 ghost nodes are used to guarantee the size of the discretized matrix to be square. R changes from 2 to 8, and c BF , c SSA , and c MFF are used to evaluate the optimal relative errors. The relative errors in Table 1 show that, as the value of R is larger than 2, the improved RBFCM can always obtain accurate results, which further indicate the stability of the improved RBFCM. Compared with the relative errors of c BF , c SSA , and c MFF , the relative errors of c BF are the best results; however, the relative errors of c BF and c SSA are also less than 10 −9 . The improved RBFCM is not sensitive to the shape parameter when R changes in a certain range. Figure 7, a problem with the partition wall is considered. The density of the partition wall is ρ = 7800 kg/m 3 , the height is H = 1 m, and the thickness L = 0.5 m. The density of water is ρ 1 = 1000 kg/m 3 . The plane strain problem described in Equations (1)-(6) is considered, where E and µ are replaced with E/(1 − µ 2 ) and µ/(1 − µ), respectively. Then E = 2.1 × 10 5 Mpa, and µ = 0.3. In the numerical calculation, 361 interior nodes and 80 boundary nodes are considered, thus, 441 fictitious centre nodes are used to guarantee the discretized matrix is square.     The numerical error is defined as

MFF SSA
where u n = u n 1,1 , . . . , u n 1,N t , u n 2,1 , . . . , u n 2,N t and u f = u Firstly, the displacements obtained from the improved RBFCM by considering different R are given in Table 2. As the trial errors in the BF are not easy to decide without exact solutions, only the MFF and SSA are used to choose the optimal shape parameters. The numerical results are compared with the results of FEM, where 50,530 degrees of freedom (DOF) are used.  Table 2 show that, as R changes from 4 to 8, the relative errors are always less than 5%. The shape parameter obtained from both SSA and MFF always leads to good numerical results for the improved RBFCM. To show the similarity of the results, the displacements on the x-axis and y-axis are given in Figures 8 and 9, respectively. The numerical results are obtained from the improved RBFCM with R = 8.    Table 2, the color maps of the displacements still show high similarity.

Results in
Example 3. In this case, a rectangular cross-section of a vertical column with a fixed bottom is given in Figure 10. The density ρ = 7800 kg/m 3 , and a uniform shear q = 1 × 10 7 N is imposed on the right side of the column. The length of the rectangular cross-section is L = 0.5 m and the height    Table 2, the color maps of the displacements still show high similarity.
Example 3. In this case, a rectangular cross-section of a vertical column with a fixed bottom is given in Figure 10. The density ρ = 7800 kg/m 3 , and a uniform shear q = 1 × 10 7 N is imposed on the right side of the column. The length of the rectangular cross-section is L = 0.5 m and the height   Figures 8c and 9c show the results of SSA. Figures 8a and 9a are the results of the FEM. Although the relative errors obtained by SSA are close to 5% in Table 2, the color maps of the displacements still show high similarity.

Example 3.
In this case, a rectangular cross-section of a vertical column with a fixed bottom is given in Figure 10. The density ρ = 7800 kg/m 3 , and a uniform shear q = 1 × 10 7 N is imposed on the right side of the column. The length of the rectangular cross-section is L = 0.5 m and the height is H = 1 m. E = 2.1 × 10 5 Mpa, and µ = 0.3. In the numerical calculations, 361 interior nodes and 80 boundary nodes are considered, and the number of ghost nodes is 441. The results in Table 3 show that, as R changes from 4 to 7, the relative errors are always less than 8%. The relative errors obtained from MFF lead to worse numerical results for the improved RBFCM compared with SSA. When R = 8, the relative errors of MFF are close to 10%. The relative errors are larger compared with example 2; this is because the stress concentrations in this case are larger. To show the similarity of the results, the displacements on the x-axis and y-axis are given in Figures 11 and 12, respectively. The numerical results are obtained from the improved RBFCM with R = 8.  The results in Table 3 show that, as R changes from 4 to 7, the relative errors are always less than 8%. The relative errors obtained from MFF lead to worse numerical results for the improved RBFCM compared with SSA. When R = 8, the relative errors of MFF are close to 10%. The relative errors are larger compared with example 2; this is because the stress concentrations in this case are larger. To show the similarity of the results, the displacements on the x-axis and y-axis are given in Figures 11 and 12, respectively. The numerical results are obtained from the improved RBFCM with R = 8. Figures 11b and 12b show the displacements obtained with the MFF. Figures 11c and 12c show the results of SSA. Figures 11a and 12a are the results of the FEM. Although the relative error of the MFF in Table 3 is close to 10%, the color map of the displacements still shows a high similarity, which further validates the effectiveness of the improved RBFCM.

Example 4.
In this case, a 3D elastic problem is given as shown in Figure 13, with the length L = 1 m. Displacement boundary conditions are given to the lower bottom surface, and traction boundary conditions are given to the other surfaces. For more details of the problem please refer to [29]. The ghost nodes are evenly distributed inside the sphere, as shown in Figure 14. To show the similarity of the results, the displacements on the x-axis and y-axis are given in Figures 11 and 12, respectively. The numerical results are obtained from the improved RBFCM with R = 8.   Table 3 is close to 10%, the color map of the displacements still shows a high similarity, which further validates the effectiveness of the improved RBFCM. Figure 13, with the length L = 1 m. Displacement boundary conditions are given to the lower bottom surface, and traction boundary conditions are given to the other surfaces. For more details of the problem please refer to [29]. The ghost nodes are evenly distributed inside the sphere, as shown in Figure 14.   Table 3 is close to 10%, the color map of the displacements still shows a high similarity, which further validates the effectiveness of the improved RBFCM.

Example 4.
In this case, a 3D elastic problem is given as shown in Figure 13, with the length L = 1 m. Displacement boundary conditions are given to the lower bottom surface, and traction boundary conditions are given to the other surfaces. For more details of the problem please refer to [29]. The ghost nodes are evenly distributed inside the sphere, as shown in Figure 14.  The relative errors of the improved RBFCM are given and compared with the traditional RBFCM in Figure 15; the accuracy of the improved RBFCM is much higher than that of the traditional RBFCM, and as the number of nodes increases, the convergence rate of the improved RBFCM is much higher. In Table 4, the errors of different R using the optimal shape parameters obtained with the MFF and the SSA are given. In total, 729 interior nodes, 485 boundary nodes, and 1214 ghost nodes are considered, as shown in Figure 14. As shown in Table 4, with the MFF and the SSA used in the improved RBFCM, we can obtain an accurate solution.  The relative errors of the improved RBFCM are given and compared with the traditional RBFCM in Figure 15; the accuracy of the improved RBFCM is much higher than that of the traditional RBFCM, and as the number of nodes increases, the convergence rate of the improved RBFCM is much higher. The relative errors of the improved RBFCM are given and compared with the traditional RBFCM in Figure 15; the accuracy of the improved RBFCM is much higher than that of the traditional RBFCM, and as the number of nodes increases, the convergence rate of the improved RBFCM is much higher. In Table 4, the errors of different R using the optimal shape parameters obtained with the MFF and the SSA are given. In total, 729 interior nodes, 485 boundary nodes, and 1214 ghost nodes are considered, as shown in Figure 14. As shown in Table 4, with the MFF and the SSA used in the improved RBFCM, we can obtain an accurate solution.  In Table 4, the errors of different R using the optimal shape parameters obtained with the MFF and the SSA are given. In total, 729 interior nodes, 485 boundary nodes, and 1214 ghost nodes are considered, as shown in Figure 14. As shown in Table 4, with the MFF and the SSA used in the improved RBFCM, we can obtain an accurate solution. The displacement contours are shown in Figure 16; the displacement contours of the improved RBFCM are similar to the exact solution. This further validates that the improved RBFCM can be extended to the 3D cases. The displacement contours are shown in Figure 16; the displacement contours of the improved RBFCM are similar to the exact solution. This further validates that the improved RBFCM can be extended to the 3D cases.

Conclusions
In this work, a novel RBFCM is proposed for elastic problems with high stability and accuracy. Elastic problems with and without analytical solutions are given to validate the improved RBFCM. Three different approaches to choosing shape parameters are tested in our examples. The numerical results show that the improved RBFCM has a high convergence rate and is less sensitive to the shape parameters. Moreover, the modified Franke formula and the sample solution approach adopted in this paper effectively settle the difficulty of choosing ideal shape parameters for the complex problems without analytical solutions, and further improve the accuracy, efficiency, and stability of the proposed method. One of our future research projects is to extend the proposed ghost node method to a local method for solving more practical real-life problems.

Conflicts of Interest:
The authors declare no conflict of interest.

Conclusions
In this work, a novel RBFCM is proposed for elastic problems with high stability and accuracy. Elastic problems with and without analytical solutions are given to validate the improved RBFCM. Three different approaches to choosing shape parameters are tested in our examples. The numerical results show that the improved RBFCM has a high convergence rate and is less sensitive to the shape parameters. Moreover, the modified Franke formula and the sample solution approach adopted in this paper effectively settle the difficulty of choosing ideal shape parameters for the complex problems without analytical solutions, and further improve the accuracy, efficiency, and stability of the proposed method. One of our future research projects is to extend the proposed ghost node method to a local method for solving more practical real-life problems.