FRACTIONAL HERMITE DEGENERATE KERNEL METHOD FOR LINEAR FREDHOLM INTEGRAL EQUATIONS INVOLVING ENDPOINT WEAK SINGULARITIES∗

In this article, the Fredholm integral equation of the second kind with endpoint weakly singular kernel is considered and suppose that the kernel possesses fractional Taylor’s expansions about the endpoints of the interval. For this type kernel, the fractional order interpolation is adopted in a small interval involving the singularity and piecewise cubic Hermite interpolation is used in the remaining part of the interval, which leads to a kind of fractional degenerate kernel method. We discuss the condition that the method can converge and give the convergence order. Furthermore, we design an adaptive mesh adjusting algorithm to improve the computational accuracy of the degenerate kernel method. Numerical examples confirm that the fractional order hybrid interpolation method has good computational results for the kernels involving endpoint weak singularities.


Introduction
Fredholm integral equations of the second kind arise in many scientific and engineering applications. Some kinds of boundary value problems can convert to these integral equations. The general form of linear Fredholm integral equation is applications, see [11,13]. There are many numerical methods developed to solve some kinds of Fredholm integral equations, see for example [3,4]. In this paper, we focus on the efficient computation of linear Fredholm integral equations of the second kind with endpoint weakly singular kernels. For approximating (1.1) with singular kernels, it is difficult to achieve a satisfactory result via traditional numerical methods. Vainikko et al. [21] discussed the properties of solutions of weakly singular integral equations of the second kind in detail. There have been many attempts to overcome the difficulties caused by the weak singularities of the kernel function [2, 6-10, 12, 16-18, 22, 25]. Recently, a modified Galerkin method was studied to show the superconvergence property for nonlinear Fredholm integral equations with weakly singular kernels [2,10]. In [6], Cao et al. successfully constructed Galerkin methods with high convergence order for the second kind Fredholm integral equations involving algebraic or logarithmic weak singularity by using singularity preserving projection methods. Collocation methods [18,25] were also used to solve linear Fredholm integral equations. Guebbai et al. [12] presented a modified degenerate kernel method by transforming the kernel function into the regular form. As for the kernel function containing weak singularities of the type (y − a) α (b − y) β , where α, β > −1, Fermo et al. [8,9] presented a modified Nyström's method by "moving" the singularities from the kernel and then regularizing the equation by means of a suitable polynomial type transformation.
It is imperative to note that one important reason for traditional numerical methods having lower convergence order for weakly singular equations is that the kernels can not be expanded as standard Taylor's series about the singularities. In order to approximate a function near its singularity, we need to introduce fractional Taylor's series [19,20], which are a generalization of Taylor's series. Liu et al. [15] discussed a general fractional Taylor's formula and its computation for insufficiently smooth functions. Like traditional Taylor's series, fractional Taylor's series usually get highly accurate approximation only in the region near the singular point, and when the argument is far away from the expanding point the series gradually become less accurate. Hence, if we want to obtain a uniform approximation for the function with a singularity, the best choice is using the fractional Taylor's expansion in a small interval involving the singularity and adopting interpolation in the remaining part of the interval. Following this idea, Wang et al. [22] successfully constructed fractional order degenerate kernel methods based on the fractional Taylor's expansion and piecewise quadratic interpolation for Fredholm integral equations of the second kind with endpoint weak singularities.
In the present paper, we aim to construct a hybrid fractional Hermite interpolation to approximate the kernel and then design a fractional degenerate kernel method to solve Fredholm integral equations of the second kind. As mentioned above, the fractional Taylor's expansion is only applied to a small interval involving the singularity and piecewise cubic Hermite interpolation is adopted in the remaining part of the interval. In order to guarantee that the fractional Taylor's expansion and the piecewise Hermite interpolation are continuously differentiable at their adjacent nodes, we further modify the coefficients of the last two terms of the fractional Taylor's expansion. In Section 2, we develop this hybrid fractional Hermite interpolation and provide the expression of the remainder for the interpolation. In Section 3, we construct a kind of fractional degenerate kernel method based on this hybrid Hermite interpolation. We discuss the condition that the piecewise hybrid interpo-lation method can converge and further estimate the convergence order. In order to improve the computational accuracy, we also present a method to generate an adaptive mesh based on the error analysis in this section. In Section 4, two typical numerical examples are given to show that the piecewise hybrid Hermite interpolation method has good computational results for the kernel functions with endpoint weak singularities. Numerical examples also confirm that the numerical convergence order of the piecewise hybrid Hermite interpolation method is consistent with the theoretical analysis. Finally, a conclusion is presented in Section 5.

Hybrid fractional Hermite interpolation
In this section, we develop a piecewise hybrid interpolation method by combining fractional Hermite interpolation with traditional cubic Hermite interpolation to approximate the kernel function. Suppose that the kernel function k (x, y) in (1.1) is not sufficiently smooth at the endpoints of the interval [a, b], and it can further be expanded as local fractional Taylor's series about y = a and y = b It is well known that fractional Taylor's series have local approximation, which means that the series are only accurate at the neighbourhood of the expanding point. Hence, it is necessary to introduce interpolation on the whole interval. Divide the interval [a, b] into a grid with nodes a = a 0 < a 1 < · · · < a n = b and step length In the subintervals [a 0 , a 1 ] and [a n−1 , a n ], k (x, y) has fractional Taylor's expansions (2.1) and (2.2), respectively. By truncating them into finite terms, we havẽ In other subintervals, we want to construct piecewise cubic Hermite interpolation. In order to connectk a ,k b with the cubic Hermite interpolation smoothly, we need to modify them first. Let where κ i (x) , i = 1, 2, 3, 4 are unknown functions to be determined by using the following interpolation conditions k n (x, a n−1 ) = k (x, a n−1 ) , ∂k n (x, y) ∂y y=an−1 = k y (x, a n−1 ).
It is straightforward to show that As for the remainders of the interpolants (2.4), (2.5), we have the following theorem.
which satisfy Proof. We only prove the formula (2.6). Subtracting (2.4) from (2.1), we obtain A straightforward computation shows Hence, we have Noting that α 1 < α 2 < · · · < α ma+1 < · · · , it is easy to show that |τ 1 (y)| ≤ 2 when y ∈ [a, a + h 1 ]. The proof of (2.7) follows in a similar way, and so is omitted. For simplicity of presentation, we still denote is sufficiently smooth with respect to y over the other subintervals [a i−1 , a i ], i = 2, 3, . . . , n − 1. Hence, we can construct piecewise cubic Hermite interpolation For the error estimate of this interpolation, we introduce the following lemma. y) is the cubic Hermite interpolating polynomial defined by (2.8). Assume that k(x, y) is four times continuously differentiable with respect to y ∈ [a 1 , a n−1 ], then the interpolation error is given by

The degenerate kernel method via hybrid Hermite type interpolation
In this section, we develop a degenerate kernel method using the above hybrid Hermite interpolation. By discussing the convergence of the method, we present an algorithm to adjust the grids to improve the computational accuracy.

The construction of the degenerate kernel method
The integral equation (1.1) can be rewritten as Approximating the kernel function k(x, y) by its hybrid Hermite interpolation and denoting the corresponding approximate solution by u n (x) yield We further rewrite (3.2) as Obviously, B i , i = 1, 2, . . . , n are the unknown vectors to be determined. Multiplying both sides of (3.2) by (x − a) αν , ν = 1, 2, . . . , m a and integrating on the interval [a 0 , a 1 ], we can obtain Here Analogously, for the last subinterval, multiplying both sides of (3.2) by (b−x) βγ , γ = 1, 2, . . . , m b and integrating on the interval [a n−1 , a n ], we acquire λ an an−1 Likewise, we rewrite (3.6) in vector form where Multiplying both sides of (3.2) by l j,γ , γ = 1, 2, 3, 4, and integrating on the subinterval [a j−1 , a j ], we obtain (3.8) The vector form of (3.8) reads where Further by letting the integral equation (1.1) can be discretized as a system of linear equations in matrix form where B is the unknown vector to be solved. Obviously, if λ is not an eigenvalue of the matrix A, the linear system (3.10) has only one unique solution. It should be pointed out that some entries of the matrix A and vector F are weakly singular integrals, which can be effectively evaluated by the modified Gauss-Legendre rule [23,24].

Convergence analysis
In this subsection, we analyze the convergence of the degenerate kernel method by using the framework of theoretical analysis introduced in [4,5]. We transform (1.1) to an equivalent operator equation by introducing an operator in a suitable Banach space  Then the operators (λ − K n ) −1 exist from X onto X for all sufficiently large n, say n ≥ N , and

Lemma 3.3 ( [4]). For the operator equations
It can be seen from (3.13) that if ∥ K − K n ∥ converges to zero rapidly, then the same is true for ∥ u − u n ∥, independent of the differentiability of u, as well as the smoothness of the force term f . This is an advantage of degenerate kernel methods. We now define a hybrid interpolation operator K n corresponding to (3.1). For u ∈ C[a, b] where k 1 (x, y), k n (x, y) are the fractional Hermite interpolations (2.4), (2.5), respectively, and k i (x, y), i = 2, 3, . . . , n − 1 are the cubic Hermite interpolations defining by (2.8). The interpolating points a 1 < a 2 < · · · < a n−1 can be arbitrary.
In order to analyze the convergence of the degenerate kernel method, we need to evaluate ∥K − K n ∥, denoted by ρ n . By using Theorem 2.1 and Lemma 2.1, we obtain are bounded, the following theorem holds.
It should be noted that even though we remove the weak singularities from the kernel k(x, y), the higher order derivatives of k(x, y) may become large when y is closed to the endpoint a or b. In order to obtain uniform precision, we should choose the subinterval lengths h 1 and h n suitably large and partition [a + h 1 , b − h n ] into non-uniform grids, which will be discussed in next subsection.

Adaptive mesh based on error analysis
In this subsection, we present a method to generate interpolating points based on the error analysis to get an optimal order of global convergence.
To this end, we describe an adaptive partition of the interval in terms of the exponents α i , β i defined by (2.1), (2.2), and a reference step length h. As we mentioned above, the length of the subinterval containing the singularity cannot be too short. We fix the length of [a 0 , a 1 ] being h 1 and [a n−1 , a n ] being h n , respectively. Furthermore, we define the distance between the node and the singularity by d i = a i − a, ord i = b − a i and the threshold value d = min{1, b−a 2 }. When d i is greater than 1, the function is smooth enough at the node a i , hence uniform step size can be used for the remaining nodes. In the following, we adjust the step lengths h i , i = 2, 3, . . . , n − 1. We . .}, which reflects the strength of the singular behavior, andξ(x) is the coefficient corresponding to (y − a)ᾱ. Sinceᾱ < 4 in most cases, the error defined by (3.16) will tend to infinity when y → a. Although we have removed the singularity from (3.16), (y − a)ᾱ −4 is still large when y is near to a. Hence, we need to adjust h i to balance the interpolation error.
Give a reference step length h. First, let d 1 = h 1 . Second, on the interpolation interval [a 1 , a 2 ], let which means that the interpolation precision on [a 1 , a 2 ] is approximately the same as the one on subinterval [a + d − h, a + d]. Hence, we can obtain Recursively, for the interval [a i−1 , a i ], i = 2, 3, . . ., we have from which we know The above procedure is repeated until d l ≥ d, and let d l = d, h l = d l − d l−1 , a l = a l−1 + h l . Until now, we successfully adjust the step lengths of the interval [a, a + d]. For the interval [b − d, b], we can also adjust the step lengths using the same procedure. For the remaining part [a + d, b − d], the function is sufficiently smooth and we choose uniform partition as follows a + d, a + d + h, · · · , a + d + ih, · · · , b − d.
Through the above argument, we obtain the interpolating nodes which have uniform interpolation precision on the whole interval. In next section we will give some numerical examples to show the superiority of the method just described.

Numerical examples
In this section, two typical examples are provided to illustrate the high accuracy of the algorithm in Section 3. Since Mathematica can easily formulate the fractional Taylor's expansions of a function about some special points and achieve arbitrary precision in numerical computation, we choose Mathematica as a platform to implement the fractional degenerate kernel method in this paper.

Example 4.1. Compute the following second kind Fredholm integral equation using the fractional Hermite interpolation method
The exact solution of this example is e √ x .
In this example, the partial derivative of the kernel function k(x, y) is weakly singular at the left endpoint x = 0 and y = 0. Referring to (2.1), the parameters of the local fractional Taylor's expansion about y = 0 are as follows We can derive from (3.15) that It follows that ∥ K − K n ∥→ 0 as α ma → ∞ (h 1 < 1) and h i → 0, i = 2, 3, . . . , n, since the series in the above formula is convergent. Hence Theorem 3.1 implies that the approximate solution u n (x) based on the fractional hybrid Hermite interpolation uniformly converges to the exact solution u(x).
In this example, we choose α ma = 7 and fix h 1 = 0.2 and hence m a = 15. First, give a reference step length h = 0.2 and generate the adaptive mesh using the algorithm in Section 3.3. We plot the absolute error of the hybrid fractional Hermite interpolation at x = 2 with logarithmic scale, shown in Fig.1. In this figure, we also plot the errors generated by two uniform meshes h = 0.2 and h = 0.1. We point out that we only plot the errors on the subinterval [0.2, 1] since uniform grid is generated on the remaining interval [1,4]. It is clear from Fig.1 that the interpolation errors are distributed evenly throughout the interval by using adaptive mesh, whereas the errors generated by uniform mesh vary from large to small as y gradually becomes large. Just as the non-uniform distribution of interpolation errors for uniform meshes, we shall see the accuracy of the degenerate kernel method is lower on the whole internal since the degenerate method always formulates a global approximate solution.
Second, we solve the Fredholm integral equation using the degenerate kernel method in Section 3. Since the expression of the approximate solution is very complicated, we only plot the logarithmic error between the accurate solution and the approximate solution in Fig.2 when the reference step length is h = 0.2. In this figure, we also plot the logarithmic error with uniform mesh h = 0.2. It can be seen  Fig.2 that the errors of the two methods have at least 10 −6 order accuracy, which means that our method is successful to solve this integral equation, and the method using adaptive mesh achieves better result on the whole interval.
Finally, we check the convergence order of our method. Define maximal absolute error E n = ∥u (x) − u n (x)∥ ∞ , numerical convergence order O n = log 2 (E n /E 2n ), where n is the number of nodes corresponding to the reference step length h. We note that a non-uniform mesh is generated in computation and hence the nodal number n is definitely increased, denoted byn. In addition, the spectral condition number of the matrix A formulated by (3.10) is discussed, denoted by Cond n . In comparison, we list the results obtained by using fractional Hermit type (FHT) degenerate kernel method, as well as the results computed by Nyström's method in Table 1. We note that Nyström's method [4,14,26] is a numerical integration method to discretize the integral equation and the convergence order of this method is determined by the accuracy of the numerical quadrature formula. In this example, we choose composite Simpson's rule. For sufficiently smooth kernel function, the convergence order of composite Simpson's rule is 4. It can be seen from Table 1 that the standard Nyström method has lower accuracy for this non-smooth kernel function.
The results of this example show that our method gets high accuracy and the convergence order is nearby 4, which is well matched with the theoretical one. The adjusted nodal numbern is not much bigger than n, which means that the method in this paper achieves high accuracy at less additional cost. Table 1 also shows that the growth of the condition number of matrix A with respect to spectral norm is relatively slow, which illustrates that the method is well-conditioned. Finally, we point out that the selection α ma = 7 is enough for all the computations, since the error of the first subinterval is 1.00553 × 10 −14 .
In next example, we will consider an integral equation with two-endpoint weak singularities and the standard Nyström's method is invalid in that case.
The exact solution of this equation is (x(1 − x)) 1 5 . In this example, the kernel function k(x, y) is weakly singular at the two endpoints y = 0, 1. Referring to (2.1) and (2.2), the parameters of the local fractional Taylor's expansions about y = 0 and y = 1 are as follows where Γ(x) is the gamma function. We can derive from (3.15) that It shows that the series ∞ p=ma+1 α p M p h 1 αp−αm a +1 is convergent when h 1 < 1.
Analogously, the series ∞ p=m b +1 β pMp h n βp−βm b +1 is convergent when h n < 1. According to Theorem 3.1, we can conclude that the approximate solution u n (x) based on fractional hybrid Hermite interpolation converges uniformly to the exact solution u(x). In this example, we take h 1 = h n = 0.1, α ma = β m b = 10, from which we know the error of the fractional Hermite interpolation is smaller than 9.80 × 10 −13 and it's accurate enough for all computations. In Fig.3, we plot the absolute errors with logarithmic scale for the hybrid Hermite interpolation based on adaptive mesh and uniform mesh using h = 0.05 when x = 0.5. It shows that the error is distributed evenly throughout the interpolation interval for the adaptive mesh. In addition, it can be easily seen that the fractional Taylor's approximation to the kernel function near its singularities is very successful.  In order to illustrate the approximate accuracy of u n (x) to u(x), we plot the absolute error function with logarithmic scale in Fig.4 when the reference step length is h = 0.0125. We can see that the method we developed has high accuracy in the whole interval and the method using adaptive mesh gets better results. Comparing Fig.3 with Fig.4, it is necessary to point out that although uniform mesh has large error only near singularities, its effect on the accuracy of degenerate kernel method is global, which leads to the overall accuracy reduction of degenerate kernel method. It also shows that the adaptive mesh for this example is necessary. In addition, the singularities of the kernel function make the standard Nyström's method invalid for this weakly singular Fredholm integral equation because k(x, y) is infinite at the two endpoints about y. Table 2 shows that the fractional Hermite interpolation method still has good computational results and the numerical convergence order is nearly 4 in this example.

Conclusion
In this paper, we design an efficient algorithm to solve the second kind Fredholm integral equation involving endpoint weak singularities for the kernel. We successfully separate the singularities of the kernel function by using fractional Taylor's series and develop a hybrid fractional Hermite interpolation method to approximate the kernel function. Based on this interpolation, we design an efficient degenerate kernel method for solving second kind Fredholm integral equation with two-endpoint weak singularities. We discuss the condition that the method can converge and give the convergence order estimation. In addition, we present an adaptive mesh generating algorithm to improve the computational accuracy. The computation shows that the numerical convergence order is well matched with the theoretical one. For some types of Fredholm integral equations with non-smooth kernels or weakly singular force terms at the endpoints, our method can get better accuracy than the standard numerical methods such as Nyström's method.