Approximate Methods for Solving Linear and Nonlinear Hypersingular Integral Equations

We propose an iterative projection method for solving linear and nonlinear hypersingular integral equations with non-Riemann integrable functions on the right-hand sides. We investigate hypersingular integral equations with second order singularities. Today, hypersingular integral equations of this type are widely used in physics and technology. The convergence of the proposed method is based on the Lyapunov stability theory of solutions of ordinary differential equation systems. The advantage of the method for linear equations is in simplicity of unique solvability verification for the approximate equations system in terms of the operator logarithmic norm. This makes it possible to estimate the norm of the inverse matrix for an approximating system. The advantage of the method for nonlinear equations is that neither the existence or reversibility of the nonlinear operator derivative is required. Examples are given illustrating the effectiveness of the proposed method.


Introduction
The importance of developing analytical and numerical methods for solving hypersingular integral equations is determined by a variety of fields of mathematics and by applications that use hypersingular integral equations.
Hadamard introduced the concept of a finite part of an integral, or the hypersingular integral in modern terminology, when studying hyperbolic equations. The Riemann boundary problem leads to hypersingular integral equations in exceptional cases. The boundary integral equations method reduces the dimensions of partial differential equations; that leads to hypersingular integral equations.
Analytical methods for solving singular and hypersingular integral equations are known only for certain particular types of equations [6][7][8]. Thus, the importance of constructing numerical solutions is clear.
Developing numerical methods for solving singular integral equations began in the middle of the last century. By now, exhaustive results have been obtained for many types of equations. A detailed account of numerical methods for solving singular integral equations as well as numerous bibliography references can be found in [9][10][11][12][13][14].
Numerical methods for solving hypersingular integral equations have been developed to a much lesser extent. Mostly numerical methods to solve hypersingular integral equations of the first kind have been developed. Numerical methods for solving hypersingular integral equations of the second kind have been much less developed. Apparently, hypersingular integral equations of the first kind are more common. Naturally, the equations of the first kind are widely used in aerodynamics (one-dimensional [15] and multi-dimensional [5,16] Prandtl equation), electrodynamics, antennae theory, etc.
The following methods are used in solving hypersingular integral equations of the first kind. Collocations, mechanical quadratures and Galerkin methods were employed to solve equations with p = 2 singularity [6,[17][18][19].
Approximate methods for solving hypersingular integral equations having singularities of order p = 2, 3, . . . , and defined on closed smooth integration contours are constructed in [20].
In [21,22] spline-collocation methods for solving hypersingular and polyhypersingular integral equations of the second kind with odd and even singularities have been developed and justified. The spline-collocation methods for solving nonlinear hypersingular and polyhypersingular integral equations have been developed and justified in [23].
An iterative projection method for solving linear and nonlinear hypersingular integral equations, and polyhyperpersingular and multidimensional hypersingular equations, was proposed in [24].
In [22] the unique solvability of hypersingular integral equations with even singularities (p = 2, 4, . . . ) was proven. Meanwhile the convergence of approximate solution to the exact one was not justified. In [24] a unique solvability of the spline-collocation method was proven. In addition, for hypersingular integral equations with bounded right-hand sides the convergence of an approximate solution sequence to the exact solution was proven under certain additional conditions.
The iterative projection method proposed here overcomes these limitations. It was shown that if the exact equation has a solution for large enough N, where N is the dimension of an approximate system of equations, an approximate solution converges to the exact one.
Hypersingular integral equations with bounded right-hand sides are a small subset of the hypersingular integral equations family. Therefore, the problem arises of constructing and justifying approximate methods for solutions for hypersingular integral equations with non-Riemann integrable functions on the right-hand sides. This paper is devoted to those issues.
A large number of works are devoted to approximate methods for solving hypersingular integral equations of the first kind To solve the Equation (1), collocation and mechanical quadrature methods [17,18], the method of orthogonal polynomials [25], the method of discrete vortices [19], the method of homotopy [26] and others are used.
In the works [27][28][29] computational schemes for the approximate solution of the Equation (1) are constructed and their justification is carried out under the assumption that the solution has the forms The hypersingular integral equations are widely used in aerodinamical problems and in the theory of antennae [30,31]. In the works [30,31] the Equation (2) is investigated under the assumption that the right-hand side has the form is the delta-function. An analytical solution of the Equation (2) with the indicated right-hand sides is obtained under the assumption that it has the form A fairly detailed review of analytical and numerical methods for solving hypersingular integral equations is given in [32].
In this paper, we propose an approach to solving linear and nonlinear hypersingular integral equations, the right parts of which contain functions with power features.
In particular, the right-hand sides of the form are considered. Here g(t) is a smooth function. Below, for simplicity of notation, we put l = 1 in (3).

Remark 1.
It can be shown that if in the hypersingular integral Equation (1) of the first kind the right side f (t) ∈ H, H is a Holder class, then the solution to this equation has the form For singular right-hand sides, the classes of solutions of (1) are unknown. Below, when constructing and justifying the computational method, we assume that the Equation (1) with a given right-hand side has a unique solution.
The proposed method has the following advantages: (1) It allows us to extend collocations and mechanical quadratures methods to hypersingular integral equations with non-Riemann integrable right sides; (2) For linear hypersingular integral equations, it allows one to verify the inverse operator existence and estimate its norm quite easily; (3) The method is stable with respect to the operator and right hand side perturbations; (4) The method does not require the existence and reversibility of the nonlinear operator derivative.
The paper is organized as follows. The continuous method for linear and nonlinear operator equations is explained in Section 2. The numerical method for solving hypersingular integral equations is presented in Section 2.

Continuous Method and Its Convergence Properties
The method we employ in the next section for solving hypersingular integral equations is based on the continuous method introduced in [33].

Continuous Method for Solving Operator Equations
The continuous method for solving operator equations is based on the Lyapunov theory of stability.
Let x(t) be a solution of the differential equation in a Banach space B which is defined for all t ≥ t 0 . The solution x(t) is said to be stable if (i) for each ε > 0 there is a corresponding δ = δ(ε) > 0 such that any solutionx(t) of (4) which satisfies the inequality |x(t 0 ) − x(t 0 )| < δ exists and satisfies the inequality |x(t) − x(t)| < ε for all t ≥ t 0 .
It is said to be asymptotically stable if in addition We will use the following notation: B(a, r) = {z ∈ B : z − a ≤ r}, S(a, r) = {z ∈ B : z − a = r}, Re(K) = (K) = (K + K * )/2, Here B is a Banach space, a ∈ B, K is a linear operator on B, Λ(K) is the logarithmic norm [34] of the operator K, K * is the conjugate operator to K and I is the identity operator.
The analytical expressions for logarithmic norms are known for operators in many spaces. We restrict ourselves to a description of the three norms.
In the n-dimensional space R n of vectors x = (x 1 , . . . , x n ) the following norms are often used: Here are analytical expressions of the logarithmic norm of n × n matrix A = (a ij ), due to the above norms of the vectors: where A * is the conjugate matrix for A. Note that the logarithmic norm of the same matrix can be positive in one space and negative in another.
The logarithmic norm has the some properties which are very useful for numerical mathematics. Let A, B be n × n matrices with complex elements; and x = (x 1 , . . . , x n ), y = (y 1 , . . . , y n ), ξ = (ξ 1 , . . . , ξ n ) and η = (η 1 , . . . , η n ) are n-dimensional vectors with complex components. Let the systems of algebraic equations Ax = ξ and By = η be given. The norm of a vector and its subordinate operator norm of the matrix are agreed upon; the logarithmic norm Λ(A) corresponds to the operator norm. Theorem 1 ([35]). If Λ(A) < 0, the matrix A is non-singular and A −1 ≤ 1/|Λ(A)|.
Some properties of the logarithmic norm in a Banach space, which are useful in numerical mathematics, are given in [34].
Let us consider in a Banach space B, the Cauchy problem Let us assume that the nonlinear operator A has a Frechet derivative and A(0) = 0.
The sufficiently satisfying conditions of asymptotically stability for the solution of the Cauchy problem (5), (6) were obtained in [36,37]. Theorem 3. Let the integral t 0 Λ(A (ϕ(τ)))dτ be non-positive (respectively, be negative and satisfy lying in a ball B(0, r) of some radius r. Then the trivial solution of Equation (5) is stable (respectively, asymptotically stable).

Remark 2.
Additionally, the Theorem is valid under r = ∞.
Let us consider in a Banach space B a nonlinear operator equation where operator A acts from B into B.
We associate Equation (7) with the Cauchy problem Let x * be a the solution of Equation (7). Let us make the change of variable x = x * + v in Equation (8). This change reduces the Cauchy problem (8), (9) to the form It is easy to see that if the trivial solution of Equation (10) is globally asymptotically stable, then lim t→+∞ v(t) → 0. So, for any initial value the solution of the Cauchy problem (8), (9) tends to x * . It follows from the next assertions which were proven in [33]. (7) have a solution x * , and let inequality

Theorem 4. Let Equation
be true on each differentiable curve g(t) lying in the Banach space B. Then the solution of the Cauchy problem (8), (9) converges to the solution x * of Equation (7) for any initial value. (7) have a solution x * , and let the following conditions be satisfied on any differentiable curve g(t) lying in the ball B(x * , r).

The inequality
holds for all t(t > 0).

Remark 3.
The sufficient condition for convergence of the Cauchy problem (8), (9) solution to the solution of the operator Equation (7) is given above. It was obtained by analysing Lyapunov stability. One of the first basic results in accretive operator theory was a relation between the solution of operator equation Au = 0, where A is a locally Lipschitzian and accretive operator, and the differential equation du dt = Au was obtained in [38]. Later, accretive operator theory and its applications for finding fixed points and constructing iterative procedures were studied by many authors. Basic results and a detailed bibliography devoted to the subject may be found in [39][40][41][42].

An Solution of Hypersingular Integral Equations with the Continuous Method
Let us consider the method of mechanical quadrature for solving hypersingular integral equation of the types and It is assumed that in the Equations (13) and (14) the right-hand sides have features of the following types In what follows, without loss of generality, we set l = 1.
Let us recall the Hadamard definition of hypersingular integrals [43].
for an integer p and 0 < α < 1, is defined as if A(x) has p derivatives in the neighborhood of the point b. Here B(x) is any function that satisfies the following two conditions: (i) The above limit exists; (ii) B(x) has at least p derivatives in the neighborhood of the point x = b.
It is easy to see [43], that the conditions (i) and (ii) are sufficient for the existence of the limit. Chikin in [44] introduced the following definition.
where ξ(v) is a function constructed so that the limit exists.

An Approximate Solution of Linear Hypersingular Integral Equations with Second Order Singularity
Consider a one-dimensional hypersingular integral equation of the type where As an approximate solution of (16), we shall seek in the form of a continuous function where ϕ k (t), k = 0, 1, . . . , N 1 , ψ k (t), k = 0, 1, . . . , N 2 is a family of basis functions. For nodes t k , k = 1, . . . , N 1 − 1, the corresponding basis elements are determined by For boundary nodes t k , k = 0 and k = N 1 the corresponding basis elements are defined as and For nodes τ k , k = 0, 1, . . . , N 2 , the corresponding basis elements ψ k , k = 0, 1, . . . , N 2 , are determined in a the similar way: For nodes τ k , k = 1, . . . , N 2 − 1, the corresponding basis elements are determined by For boundary nodes τ k , k = 0 and k = N 2 the corresponding basis elements are defined as and To simplify the description of computational scheme, we introduce the following notation: Here v i = t i , i = 0, 1, . . . , N 1 , v N 1 +i = τ i , i = 1, 2, . . . , N 2 , Let us recall that the points t N 1 and τ 0 coincide. Applying the collocation method on the knotsv k , k = 0, 1, . . . , N * + 1 to the Equation (16), we obtain the following system of algebraic equations for finding unknown coefficients {γ k } of the polygon (17) Using the definition of hypersingular integrals, we receive: Here ∑ l , ∑ l , ∑ l , ∑ l indicates a summation over l = N 1 , , respectively. Detailed calculations are given in [23].
Repeating the proof presented in [24] we see that the approximation of f (t) ∈ W 1 (( f : f 1 , f 2 ; M, c)) by piecewise linear functions constructed on the basis ζ k (t), k = 0, 1, . . . , N * + 1, has the error In this paper, we denote the constants that do not depend on N by C.
, and x * (1) Repeating the arguments given in [24], we arrive at the following statement.
Theorem 6. Let the following conditions be fulfilled: (1) Equation (16) has the unique solution x * (t) Then the system of Equations (37)-(42) has a unique solution x * N (t) and the following estimate holds:
By computing the hypersingular integrals in (50), we can rewrite the system of Equation (50) as Here ∑ , ∑ , ∑ indicate summations over l = k, l = N 1 , l = N 1 + 1, respectively. The Frechet derivative on a vector (α 0 , α 1 , · · · , α N * +1 ) in the space R N * +1 is equal to Here ∑ , ∑ ∑ indicate summations over l = k, l = N 1 , l = N 1 + 1, respectively. The notation h 3 (t, τ, u) = δh(t,τ,u) δu is used here. Let the Equation (49) has the unique solution x * (t) inside the ball B(x * , δ). We shall assume that the Frechet derivative (57) in the ball R N * +1 (x * , δ) satisfies the conditions of Theorem 5. Thus, according to statements of the Theorem 5, the solution of the system of differential equations converges to the solution of the Equation (49).
Thus, we have proven the following statement.
Then the system of Equations (51)-(56) has a unique solution inside the ball B(x * , δ), and the solution of Equation (58) converges to this solution.
The effectiveness of the presented algorithms is illustrated by solving two hypersingular integral equations modeling aerodynamics problems. Example 1. Let us illustrate the effectiveness of continuous method by solving the following linear hypersingular equation where f (γ 1 , γ 2 , t) is the given right-hand side of the equation: The exact solution of the equation is x(t) = (x 1 (t), x 2 (t)); x i (t) = a i + γ i t, i = 1, 2.
To solve the Equation (59) numerically we use the continuous method for solving operator equations and arrive to the following evolution equation Nodes v k ,v k , k = 0, 1, . . . , N * + 1, have been entered above. In Figure 1 we show the trajectories of the exact solution of the Equation (59); its approximate solution, received with continuous method; and values of error. Here a 1 = 1, a 2 = 1.5, γ 1 = 0.5, γ 2 = 0.3.
where f (γ 1 , γ 2 , t) is the given right-hand side of the equation: The exact solution of the equation is It easy to see that, if x(t) is a solution of the Equation (60), then functions −x(t), |x(t)| and −|x(t)| are solutions of this equation too.
To solve the Equation (60) numerically we use the continuous method and receive the following evolution equation At first, we take α k (0) = 0.0 as an initial condition in order to demonstrate applicability of our method in cases of the Newton-Kantorovich method, the minimal residual method and other numerical methods; using in their construction the derivative of nonlinear operator is not applicable. Indeed, in this case the Frechet derivative (57) is not only degenerate-and, therefore, not invertable-but is an identical zero.
In Figure 2 we put a 1 = 1, a 2 = 1.4, γ 1 = 0.5, γ 2 = −0.4. In Figure 2 we show the trajectories of the exact solution of the Equation (60), its approximate solution, received with continuous method and values of error.
The exact solution at t = 0 has a jump discontinuity of h = 0.4. The slopes of the exact solution also change at t = 0. In Figure 2 we demonstrate that the numerical solution approximates the exact one at [−1, 0) well. At t = 0 the approximate solution has a jumph = 0.15.

Summary and Discussion
An iterative projection method for solving linear and nonlinear hypersingular integral equations has been proposed. The method is based on the use of sufficient conditions for asymptotic stability of ODE systems. Stability conditions are expressed in terms of the logarithmic norms of the corresponding matrices. In a number of spaces often used in computational mathematics, the calculation of logarithmic norms does not cause difficulties, even for large-dimensional matrices.
What are the advantages of the presented method?
(1) The method is applicable for solving linear and nonlinear hypersingular integral equations, whose right-hand sides contain non-Riemann integrable functions. (2) In Section 3.1 the continuous method is applied to linear hypersingular integral equations with the singularities of the second order. The conditions for the unique solvability of the constructed computing scheme are obtained and the convergence of the sequence of approximate solutions to the exact one is proven. It is shown that for linear hypersingular integral equations, the method converges for sufficiently large N and for b(t) = 0, t ∈ [−1, 1]. (3) In Section 3.2 the continuous method is applied to nonlinear hypersingular integral equations with the singularities of the second order. Conditions are given for the convergence of the constructed iterative spline-collocation method to the solution of a nonlinear hypersingular integral equation. It should be noted that the method is applicable to hypersingular integral equations of the first and second kinds.
The detailed bibliography of approximation methods of hypersingular integral equations of the first and the second kinds is given in [32]. The bibliography on solving hypersingular integral equations of the first kind is presented in [45].
Mostly, papers devoted to hypersingular integral equations of the first kind focused to seek solutions in the class of functions √ 1 − t 2 φ(t), where φ(t) is a smooth function. The presented method provides solutions in a general form.
The theoretical justification of the method is based on Lyapunov stability theory. It connects convergence of the method to the sign of the approximate system matrix logarithmic norm.
Said justification has advantages that allow us 1. To obtain a set of convergence conditions owing to logarithmic norm values in various spaces; 2. To determine the norm of the inverse matrix of an approximate system; 3. To determine stability boundaries for solutions with respect to variations of kernels and right-hand sides of the equations.
The major advantage of the method for nonlinear equations is as follows. The Newton-Kantorovich method requires the Frechet derivative reversibility at each iteration step. Similar conditions are required when using other iteration methods. Our method lacks such a deficiency. It does not put any restrictions on the Frechet derivative of the nonlinear operator.

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