Solving the Linear Integral Equations Based on Radial Basis Function Interpolation

The radial basis function (RBF) method, especially the multiquadric (MQ) function, was introduced in solving linear integral equations. The procedure of MQ method includes that the unknown function was firstly expressed in linear combination forms of RBFs, then the integral equation was transformed into collocation matrix of RBFs, and finally, solving the matrix equation and an approximation solution was obtained. Because of the superior interpolation performance of MQ, the method can acquire higher precision with fewer nodes and low computations which takes obvious advantages over thin plate splines (TPS) method. In implementation, two types of integration schemes as the Gauss quadrature formula and regional split technique were put forward. Numerical results showed that the MQ solution can achieve accuracy of 1E − 5. So, the MQ method is suitable and promising for integral equations.


Introduction
In mechanics, electromagnetic theory, thermal, radiation science, and other fields, the mathematical models of science or engineering problems can be attributed to the integral equations.Comparing to the differential equation, the integral equation embodies the initial value or boundary information just in a brief form.Moreover, the numerical integral brings smaller relative error than the numerical differentiation.So it can obtain higher accuracy, and the research on numerical solution of integral equations is of great importance.The general form for linear integral equation is  () =  ∫ Ω  (, )  () d +  () ,  ∈ Ω. (1) In the above, () is free function, (, ) is kernel function, () is the unknown function, Ω is integral region, and  is the parameter where the (), (, ), and  are known.The numerical methods for solving integral equations mainly include the undetermined coefficient approximation method, the direct numerical integration method, the successive approximation method, and the approximation method for specific kernel function.
Currently, a variety of new numerical methods were proposed, for instance, the combination method of Galerkin and collocation with Haar wavelet function for the first kind of linear Fredholm equation [1], the Galerkin method with Daubechies wavelet for the second kind of linear Volterra equations [2], however, the calculation error of wavelet methods is relatively large.The Sinc-collocation method [3] and the modified block pulse function [4] are for one dimension problems.The Chebyshev collocation method [5], block pulse method [6], and polynomial approximation with smooth kernel function [7] are for solving two-dimensional Fredholm or Volterra equations.In the above unknown function approximation methods, the polynomial, Chebyshev polynomials, or wavelet function, were adopted, respectively.However, the polynomial interpolation is unstable due to the strong rigidity of the polynomial function space, the interpolation error is relatively large for wavelets, the suitable interpolation basis functions for two or multidimensional problems are hard to construct, and so on.The interpolation basis function with high accuracy and convenience for multidimensions is particularly necessary.
And, fortunately, the radial basis function (RBF) not only is insensitive to the space dimensions but also has high interpolation accuracy.The radial basis functions were first developed by Hardy [8] in 1971 as a multidimensional scattered interpolation method in modeling of the earth's gravitational field.It was not recognized by most of the academic researchers until Franke [9] published a review paper in the evaluation of two-dimensional interpolation methods.The RBF method in solving integral equations was initially proposed in 2006 [10][11][12][13].The related research attracted a lot of attention recently.Therefore, this paper attempted to adopt the RBF method for solving linear integral equations.The RBF and interpolation principle were described in Section 2. And, in Section 3, the RBF method for solving linear integral equation and the numerical integration scheme for coefficients matrix were proposed.And, finally, the RBF method was applied in one-and two-dimensional integral equations; its performance was analyzed and compared.

Radial Basis Function and Interpolation Principle
In the above, c means the center of basis function and  is the shape parameter which is generally associated with the distance between adjacent centers.So  = ||c  − c  ||, and  is also called the shape parameter.The principle of RBF interpolation regards the unknown function as linear combination of radial functions.So the approximation function can be obtained after calculating the coefficients.For a series of known data points (  , (  )),  = 1, 2, . . ., , the approximation function can be constructed by RBF as Then, substituting the interpolation data points into the above equation, we obtain With its abstract matrix form, where [f  ] = [( 1 ), ( 1 ), . . ., (  )]  is the known data, [] = [ 1 ,  2 , . . .,   ]  is the coefficients, and [Φ  ] is the coefficient matrix whose element calculation formula is The MQ function is globally defined which results in a full resultant coefficient matrix.And, consequently, the coefficient matrix in ( 5) is nonsingular and usually ill conditioned.The accuracy of the MQ solution depends heavily on the shape parameter.In general, for a fixed number of centers , smaller shape parameters produce more accurate approximations, but they also are associated with a poorly conditioned interpolation matrix.Also, the condition number grows with , for fixed values of the shape parameter .The choice of this optimal value is still under intensive investigation.Therefore, solving (5) can acquire The approximation value f() at any point  is

Solving Linear Integral Equations Based on the RBF Interpolation
The basic idea of RBF method for solving integral equation is to employ a linear combination of RBF to approximate the unknown function; thus, the integral equation is transformed into the combination form of RBFs and their coefficients.Then, the weight coefficients of RBF can be calculated by using the collocation type weighted residual method.And, finally, an approximate representation of the unknown function can be obtained.Specifically, the RBF method for solving linear integral equations is proposed in the following.

RBF Interpolation to Approximate One-Dimensional Linear Integral Equation.
The general form of one-dimensional linear integral equation is Here,  |  represents the upper limit of integral.The constant  or variable  is, respectively, corresponding to Fredholm or Volterra equation.In the following discussion, we choose the upper limit  for demonstration.In RBF approximation, the unknown function () is expressed as a combination of RBF.Thus, the approximate formula ( 9) is obtained: Furthermore, for any collocation point   in [, ], by using the linear property of integral, we obtain the specific expression of the above equation as when choosing  collocation points, hence forming the collocation equations in matrix form as where The notations   and   represent correspondingly the collocation point and RBF center, while the element   is the definite integral of function ℎ() which can be calculated through the following Gauss quadrature formula: The integral range can be transformed with formula  = ( + )/2 + ( − )/2 = (), so the final calculation formula is However, the traditional regional split technique was adopted in this paper.Firstly, the irregular region was split geometrically.Then, the entire domain definite integral can be obtained by adding up all the integral of subregions which is indicated as the following: where Ω  represents the th subregion.Specifically, the triangular element is chosen because of its strong adaptability in conforming two-dimensional irregular domain.The triangulation can be implemented automatically by the PDE toolbox in Matlab.Then, let the vertex coordinates denote by (  ,   ),  = 1, 2, 3, so the element area   and the three midpoints (  ,   ) of element boundary were easily obtained.And, finally, the quadratic interpolation formula was applied to calculate the definite integral as

The Implementation Procedure of RBF for Solving Linear
Integral Equation.In summary, the implementation procedure of RBF for solving linear integral equation was listed as the following.
(1) Setting the centers c of RBFs and shape parameters , then setting the collocation points x  , the numbers of RBFs and collocation points are, respectively,  and  in general.Note that, to simplify and speed up the computations, we can identify the set of centers with the set of collocation points.
(2) Substituting approximate solution based on RBF into the linear integral equation, then the collocation conditions are imposed to form the coefficient matrix.
(3) Matrix elements calculation: the Gauss quadrature formula is adopted for one-dimensional and two-dimensional regular domain problems, while for irregular domain case, the numerical integration scheme based on regional split technique was employed.
( (5) RBF solution for unknown function then can be achieved as

Error Analysis.
Utilizing the same way in error analysis as [13], let the operator  :  1 →  1 for one-dimensional case be defined as Now, we can rewrite (8) in abstract form as Let  be a bounded operator with |||| < 1, and then the operator  −  is a contraction operator; by the Banach contraction mapping principle, (8) has a unique solution.In such a case, the geometric series theorem implies that the operator  −  has a bounded inverse with So, the error of the presented RBF method is mainly based upon the RBF interpolation error which emerges in (3) and the error caused by performing the numerical integration scheme over domain which appears in (18).For sufficiently large nodes , the error of the RBF interpolation is dominated over the integration error and, thence, increasing the number of nodes in the numerical integration method has no significant effect on the error.

Numerical Experiments and Analysis
In this section, the RBF method is applied to some numerical examples involving one-dimensional and two-dimensional integral equations on the regular and nonrectangular regions.Two types of RBF as multiquadric (MQ) function and thin plate splines (TPS) function were adopted simultaneously in the following simulations.Their parameters were set as the same except that the additional shape parameter was provided for MQ.In order to measure the accuracy of the method, the maximum absolute error (MAE), the maximum relative error (MRE), and the root mean squared error (RMSE) have been used, respectively.The RMSE has the following definition as where   (  ), (  ) denote the exact and approximation value at point   and symbol   is the total testing points.
Example 1.Consider the following Fredholm integral equation in [1]: The exact solution for this equation is () =  2 .
Setting RBF centers with interval ℎ = 0.2, the collocation points are the same as the center points  = 6.For MQ, the shape parameter  = 12.Then, the 20-point Gauss quadrature formula is adopted in matrix element calculation.And, finally, the numerical result was shown in Table 1.
Then, the shape parameter was analyzed with discrete interval of 0.1; the RMSE and the condition number of matrix [Φ − Κ] were calculated in case of ℎ = 0.2.It is obvious to see in Figure 1 that the optimal shape parameter lies in the interval of [10,20] showed in red dashed box.And, in the whole range of [0.1,100], the MQ method can maintain the RMSE at the level of 1 − 3. The curve of condition number was shown in Figure 2, and it fluctuates in the range of [1 + 16, 1 + 18] when  is greater than 20.So the MQ method is stable when the shape parameter varies in a relative wide range.
And, finally, the RBF's center spacing was analyzed in case of  = 12.The RMSE curve was shown in Figure 3.It can   be seen clearly that the optimal center spacing was ℎ = 0.2 and 0.0833.But for mostly center spacing, the MQ method can maintain RMSE at the level of 1 − 3.As for condition number, whose curve also fluctuates in the range of [1 + 13, 1 + 20], it was shown in Figure 4. Example 2. Consider the one-dimensional Volterra integral equation in [4]: The exact solution is () = 2 sin().Firstly, setting the RBF center spacing ℎ = 0.1, the collocation points and centers are the same as  = 11.The shape parameter  = 10 for MQ and the 20-point Gauss quadrature formula is adopted.The total test points are 101 with equidistant ℎ  = 0.01.We focus on the center spacing variation; the approximation solutions were shown in Figure 5 for TPS method.Generally, the TPS can acquire better approximation as the center spacing decreases which means that the solution of ℎ = 0.05 was better than that of ℎ = 0.1.However, numerical simulation also indicated that keeping on decreasing ℎ may result in unstable solution.
The comparison between MQ and TPS was shown in Table 2.
Therefore, we can conclude from Table 2 that the MQ method takes obvious advantages over TPS whose accuracy is improved about 3 or 4 orders of magnitude.And the MQ method is also much stable for wide range of center spacing ℎ.The MQ especially can acquire high precision even with larger center spacing as ℎ = 0.25.Compared to the method in [4], the MQ also showed its superiority; the highest precision of modified block pulse functions method is about the level of 1 − 2. ( The exact solution is (, ) =  cos().
Setting the RBF center spacing ℎ = 0.1, the collocation points and centers are  = 121.The shape parameter is  = 35 for MQ and the 20-point Gauss quadrature formula was also adopted in both  and  directions.The total testing points are 441 with equidistant ℎ  = 0.05.The RMSE, MAE, MRE, and condition number for MQ are 2.6583 − 6, 5.5309 − 6, 1.6077 − 4, and 8.0241 + 18, respectively.However, the above values for TPS are 1.2273 − 3, 9.9312 − 3, 3.4146−1, and 1.9795+3.And comparing to the rational Haar function method in [14] whose accuracy was about level 1 − 2 and 1 − 3, the MQ acquires much higher accuracy than others.The error surfaces of MQ and TPS were shown in Figure 6.Then, analyzing the shape parameter of MQ with discrete step of 1.0, the RMSE and the condition number were calculated in case of ℎ = 0.1.The local optimal shape parameters in this instance were  = 7, 35, 58, and 64 revealed in Figure 7.The overall RMSE of MQ method was about the level of 1 − 4 and 1 − 5. And, correspondingly, the curve of condition number was similar to that in Figure 2 except that the curve fluctuates in the range of [1 + 18, 1 + 20] when  is greater than 10.
The exact solution is (, ) =  cos().Firstly, setting the RBF centers and collocation points as Figure 8, the total number was  = 35.The shape parameter was  = ℎ = 2.3 for MQ and the regional splitting with triangular.
Then, the PDE toolbox in Matlab was utilized in triangulation.The above irregular region Ω was split into 2272 triangular elements with 1225 nodes.In simplicity, the testing points were chosen at the element nodes.And the quadratic form of interpolation was employed in calculating the definite integral of subregions in (20).The results showed that the RMSE, MAE, MRE, and condition number are 2.6260 − 5, 3.0345−5, 1.5783−3, and 3.3527+14.However, the above values for TPS method are 3.9563−3, 1.9925−2, 8.6492− 1, and 4.3509 + 2. The error curves of MQ and TPS are shown in Figure 9. Comparing to [13] where 33 collocation points were selected, the RMSE is 1.96−4, and for GA weight function the RMSE is 1.65−4.So the RBF method proposed in this paper can obtain higher accuracy.

Conclusion
In this paper, the RBF method, specifically, the MQ, was introduced in solving linear integral equations.The following conclusions can be drawn.
(1) The procedure of MQ method includes that the unknown function was firstly expressed in linear combination forms of RBFs, then the integral equation was transformed into collocation matrix of RBFs, and finally, solving the matrix equation and an approximation solution was obtained.
(2) Due to the excellent interpolation performance of MQ, the method can acquire higher precision with fewer nodes and low computations.Numerical results showed that the MQ method can achieve solution accuracy of 1 − 5. Numerical results showed that the MQ method takes obvious advantages over TPS.In addition, the method was very convenient for solving higher dimensional integral equations because the only modification in implementation is the distance formula.
(3) One of the key issues in MQ method is the definite integral calculation.Two types of integration schemes as the Gauss quadrature formula and regional split technique were put forward.Specifically, the former is adopted in one-dimensional and two-dimensional regular domain problems and the latter is employed in irregular domain case.

Figure 2 :
Figure 2: The condition number curve versus the shape parameter .

Figure 5 :
Figure 5: The TPS solutions for different RBF center spacing ℎ.

Figure 7 :Figure 8 :
Figure 7: The RMSE curve versus the shape parameter .

Table 1 :
Comparison between RBF methods and Haar wavelet methods.

Table 2 :
Error comparison between MQ and TPS.