Comparison of Selected Numerical Methods for Solving Integro-Differential Equations with the Cauchy Kernel

Abstract

The integro-differential equation with the Cauchy kernel is used in many different technical problems, such as in circuit analysis or gas infrared radiation studies. Therefore, it is important to be able to solve this type of equation, even in an approximate way. This article compares two approaches for solving this type of equation. One of the considered methods is based on the application of the differential Taylor series, while the second approach uses selected heuristic algorithms inspired by the behavior of animals. Due to the problem domain, which is symmetric, and taking into account the form of the function appearing in this equation, we can use this symmetry in some cases. The paper also presents numerical examples illustrating how each method works and comparing the discussed approaches.

1. Introduction

The integro-differential equation with Cauchy kernel appears in many different technical problems, such as circuit analysis, gas infrared radiation studies, molecular conduction, modeling of chemical processes or a model of the interaction of neurons [1,2,3,4,5].

The Singular integro-differential Cauchy equation discussed in this article has the following form

$$2{\phi }^{\prime }\left(x\right)+\lambda {\int }_{-1}^1\frac{\phi \left(t\right)}{t-x}dt=f\left(x\right),-1<x<1,$$

(1)

in which $\phi$ is the unknown function, $\lambda \in \mathbb{R}$ and f is known function. Due to the problem domain, which is symmetric around zero, and taking into account the form of the f function appearing in this equation, we can use this symmetry in some cases. The function $\phi$ is assumed to be many times differentiable in the interval $(-1,1)$.

Equation (1) with appropriate assumptions has a solution; however, apart from the null solution ($\phi \left(x\right)=f\left(x\right)=0$), the analytical form of such a solution is not known in general. There are some analytical methods for solving the Equation (1) (see, e.g., [6,7,8]), but in many cases they are general because they can only be used in certain specific cases of the form of the Equation (1). There are also known numerical methods that can be used to solve the Equation (1) or equations with forms similar to it. These include the following: a projection method or the Galerkin method using the orthogonal basis of Legendre polynomials or using the Taylor series, Bernstein polynomial (see, e.g., [9,10,11,12,13,14,15]). However, these methods also have their drawbacks and limitations.

Equations such as integro-delay differential equations (IDDEs) are considered in [16]. The authors focus on certain perturbed and un-perturbed nonlinear systems of continuous and discrete integro-delay differential equations. The work discussed various types of solution stability. It also includes relevant theorems as well as numerical applications. In the paper [17], Tunç et al. consider a class of scalar nonlinear integro-differential equations with fading memory. The article is devoted to examining concepts such as convergence, asymptotic stability, uniform stability, boundedness and square integrability according to considered equations. The authors present a numerical example illustrating the theorem included in the paper.

Taylor series have a wide range of applications. As mentioned above, they can also be used for issues similar to the problem discussed in this work. Of course, apart from this, there are many problems that can be solved using Taylor series. It might seem that this method has become forgotten, but many researchers still use it in their research. Interesting applications include the following: research related to determining the spectral shift [18], solving fractional equations [19], image analysis [20], market analysis [21], biological and chemical research [22], research on bivalent functions [23], engineering design [24], controller design [25], constructing alloy models [26] and many others.

This article attempts to solve the problem (1) using meta-heuristic algorithms. For this purpose, a functional describing the error of the approximate solution is constructed, and selected meta-heuristic algorithms are responsible for finding the minimum of the error function. Algorithms of this type are very popular due to their accessible description, low requirements for the objective function and effectiveness in many optimization problems. For example, in [27], one of such method was compared to an iterative method to solve the inverse problem. The problem considered in the article concerned the anomalous diffusion model with the Caputo fractional derivative, while the inverse problem involved the identification of the thermal conductivity coefficient. Dai et al. [28] use the Whale Optimization Algorithm (WOA) to improve the efficiency of indoor mobile robots. Mate-heuristic was used in the field of path planning of mobile robots. The results presented in the article showed that the application of WOA turned out to be effective.

In this article two meta-heuristic algorithms were used, namely Artificial Bee Colony (ABC) and Whale Optimization Algorithm (WOA). The applications of these algorithms are wide, and some of them can be found in the papers [29,30,31,32,33,34]. However, the articles [35,36,37,38] present the applications of ABC and WOA in various types of symmetry problems.

2. Proposed Solution Methods

The Taylor series method discussed in this section is similar to the DTM method, which has a wide range of applications, e.g., for solving the ordinary differential equations and their systems, for solving integro-differential equations with delayed argument. Also in this case, use of DTM would be beneficial, but the specificity of the task does not allow it, or more precisely—the initial condition (the $\phi \left(0\right)$ value) is missing in the Equation (1). In some special cases the value of $\phi \left(0\right)$ or $\phi (x_0\ne 0)$ (which makes it difficult but not impossible to use the DTM method) may be known (see e.g., [11,39]) and then the DTM method can be used. However, in this article, in order to expand the class of Equation (1), we discuss the use of the Taylor series in the next subsection.

2.1. Method Using Taylor Series

In this subsection, description of the method for solving the considered equation is presented. For this purpose, the function is expanded into a Taylor series. First, expand the function $\phi \left(t\right)$ into a Taylor series around the point x. Then we get

$$\phi \left(t\right)=\sum _{i=0}^{\infty }\frac{{\phi }^{\left(i\right)}\left(x\right)}{i!}{(t-x)}^i.$$

(2)

Taking into account the Formula (2), left side of Equation (1) takes the following form

$$2{\phi }^{\prime }\left(x\right)+\lambda {\int }_{-1}^1\frac{\phi \left(t\right)}{t-x}dt=2{\phi }^{\prime }\left(x\right)+\lambda {\int }_{-1}^1\left(\frac{\phi \left(x\right)}{t-x}+\sum _{i=1}^{\infty }\frac{{\phi }^{\left(i\right)}\left(x\right)}{i!}{(t-x)}^{i-1}\right)dt.$$

After integration and applying appropriate transformations, we obtain

$$\begin{array}{cc}\hfill 2{\phi }^{\prime }\left(x\right)& +\lambda {\int }_{-1}^1\frac{\phi \left(t\right)}{t-x}dt=\lambda \phi \left(x\right)ln\left|\frac{x-1}{x+1}\right|+2(\lambda +1){\phi }^{\prime }\left(x\right)\hfill \\ \hfill & +\lambda \sum _{i=2}^{\infty }\frac{{\phi }^{\left(i\right)}\left(x\right)}{i·i!}\left({(1-x)}^i+{(-1)}^{i+1}{(x+1)}^i\right).\hfill \end{array}$$

For case $-1<x<1$, we have

$$ln\left|\frac{x-1}{x+1}\right|=-2\sum _{i=1}^{\infty }\frac{x^{2i-1}}{2i-1},$$

and expanding the functions $\phi \left(x\right)$ and ${\phi }^{\prime }\left(x\right)$ into Maclaurin series, we obtain

$$\phi \left(x\right)=\sum _{i=0}^{\infty }\frac{{\phi }^{\left(i\right)}\left(0\right)}{i!}{x}^i,{\phi }^{\prime }\left(x\right)=\sum _{i=0}^{\infty }\frac{{\phi }^{(i+1)}\left(0\right)}{i!}{x}^i.$$

Similarly, expanding the subsequent components of the last sum, using the products of appropriate series, we present the left side of the Equation (1) as a polynomial of the variable x with unknown coefficients depending on the value at zero of the function $\phi$ and its subsequent derivatives. Expanding the function $f\left(x\right)$ (the right side of the discussed Cauchy equation) into Maclaurin series

$$f\left(x\right)=\sum _{i=0}^{\infty }\frac{f^{\left(i\right)}\left(0\right)}{i!}{x}^i$$

and comparing the corresponding coefficients at subsequent powers of the variable x, we construct a system of equations with unknowns $\phi \left(0\right)$, ${\phi }^{\prime }\left(0\right)$, ${\phi }^{″}\left(0\right)$, ….

Of course, since we cannot solve an infinite system of equations, instead of the corresponding series, we take their n-th partial sums. For example, for $n=5$ and $\lambda =1$, we would have the following system of equations

$$\left\{\begin{array}{c}4{\phi }^{\prime }\left(0\right)+\frac{1}{9}{\phi }^{\left(3\right)}\left(0\right)+\frac{1}{300}{\phi }^{\left(5\right)}\left(0\right)=f\left(0\right),\hfill \\ -2\phi \left(0\right)+3{\phi }^{″}\left(0\right)+\frac{1}{36}{\phi }^{\left(4\right)}\left(0\right)=f^{\prime }\left(0\right),\hfill \\ -2{\phi }^{\prime }\left(0\right)+\frac{4}{3}{\phi }^{\left(3\right)}\left(0\right)+\frac{1}{18}{\phi }^{\left(5\right)}\left(0\right)=\frac{f^{″}\left(0\right)}{2},\hfill \\ -\frac{2}{3}\phi \left(0\right)-{\phi }^{″}\left(0\right)+\frac{1}{6}{\phi }^{\left(4\right)}\left(0\right)=\frac{f^{‴}\left(0\right)}{6},\hfill \\ -\frac{2}{3}{\phi }^{\prime }\left(0\right)-\frac{1}{3}{\phi }^{\left(3\right)}\left(0\right)+\frac{1}{6}{\phi }^{\left(5\right)}\left(0\right)=\frac{f^{IV}\left(0\right)}{24},\hfill \\ -\frac{2}{5}\phi \left(0\right)-\frac{1}{3}{\phi }^{″}\left(0\right)-\frac{1}{12}{\phi }^{\left(4\right)}\left(0\right)=\frac{f^V\left(0\right)}{5!}.\hfill \end{array}\right.$$

Solving the above system of equations using the nth partial sums, we find the subsequent values $\phi \left(0\right)$, ${\phi }^{\prime }\left(0\right)$, …, ${\phi }^{\left(n\right)}\left(0\right)$, thanks to which we obtain corresponding coefficients of the expansion of the sought function $\phi \left(x\right)$ into the Maclaurin series. The approximate solution is then assumed to be of the following form

$$\phi \left(x\right)\approx \sum _{i=0}^n\frac{{\phi }^{\left(i\right)}\left(0\right)}{i!}{x}^i.$$

2.2. Approach with Metaheuristic Optimization Algorithms

The previously discussed method, as well as many other numerical methods, uses polynomials to solve the Equation (1). Based on this assumption, we use a heuristic algorithm to find the appropriate polynomial. The sought polynomial has a degree of n and it approximates the solution function $\phi$

$$\phi \left(x\right)\approx {\phi }_n\left(x\right)=\sum _{i=0}^n{a}_i{x}^i,$$

where $a_i\in \mathbb{R}$, $i=0,1,\dots ,n$ are unknown coefficients of the sought polynomial and, at the same time, arguments of the functional $\mathcal{F}$ are minimized by heuristic algorithms. This functional measures the fit of the polynomial (approximate solution) to the exact solution of the Equation (1). Due to this definition, the functional $\mathcal{F}$ takes the following form:

$$\mathcal{F}=\frac{1}{2}{\int }_{-1}^1\left|2{\phi }_n^{\prime }\left(x\right)+\lambda {\int }_{-1}^1\frac{{\phi }_n\left(t\right)}{t-x}dt-f\left(x\right)\right|dx.$$

Recently, heuristic algorithms have attracted great interest from researchers and are used for many different types of tasks. For example, the ABC algorithm was used in an inverse continuous casting problem, image and signal processing [40]. In this article, we apply and compare two metaheuristic swarm algorithms for searching for the global minimum of the objective function. These algorithms are Artificial Bee Colony (ABC) and Whale Optimization Algorithm (WOA). In Section 3, these two metaheuristic algorithms are described in detail.

3. Metaheuristic Algorithms for Optimization Problems

This section is devoted to the description of two algorithms that are popular in the scientific literature heuristic. Both algorithms have their background in nature and are inspired by the communication of a group of animals—bees and whales. These algorithms have been widely used in optimization problems, hence their selection in this article.

3.1. Whale Optimization Algorithm

This section describes the Whale Optimization Algorithm (WOA), which is used to solve the considered optimization problem. The inspiration for this algorithm is the behavior of whales when searching for prey. The mechanism these mammals use when hunting for prey seems interesting and served as a background for the algorithm. In the literature, the hunting method of these animals is called the bubble-net feeding method. To put it simply, when hunting, whales first dive into the water and then move along a spiral path, circling the prey and creating bubbles. Then they use one of several techniques: coral loop, lobtail, and capture loop. This behavior has only been observed in humpback whales. Rules and principles governing the algorithm are presented next in this section.

• In each of iteration (t is number of iteration), the best individual in the population ${\mathbf{x}}_{best}^t$ is determined. This individual’s position is closest to the prey, and other individuals move towards it according to the following formulas

  $$\mathbf{d}=|\mathbf{c}\ast {\mathbf{x}}_{best}^t-{\mathbf{x}}^t|,$$

  (3)

  $${\mathbf{x}}^{t+1}={\mathbf{x}}_{best}^t-\mathbf{a}\ast \mathbf{d},$$

  (4)

  where $\mathbf{a},\mathbf{c}$ are vectors of random coefficients, ${\mathbf{x}}^t$ is an individual in population t, and $\mathbf{d}$ is a scaled vector of distance between an individual ${\mathbf{x}}^t$ and the best individual ${\mathbf{x}}_{best}^t$ in the population. The operation ∗ means multiplying vectors element by element. Vectors $\mathbf{a}$, $\mathbf{c}$ are calculated in following way

  $$\mathbf{a}=2\mathbf{w}\ast \mathbf{r}-\mathbf{w},$$

  (5)

  $$\mathbf{c}=2\mathbf{r},$$

  (6)

  where $\mathbf{r}$ is a vector of random numbers in the range $[0,1]$, and $\mathbf{w}$ is a vector of numbers decreasing from iteration to iteration from 2 to 0. The values of the vector $\mathbf{w}$ in a fixed iteration t are assumed to be $2\frac{T-t}{T}$, where T is the maximum number of iterations. Decreasing the $\mathbf{w}$ parameter in each iteration simulates the process of narrowing the space to search for a victim. This technique is called shrinking encircling mechanism.
• The next important stage is the mechanism of spiral-shaped movement of whales. Mathematically, we describe this process with the following formula

  $${\mathbf{x}}^{t+1}={\mathbf{d}}^{\prime }{e}^{br}cos\left(2\pi r\right)+{\mathbf{x}}_{best}^t,$$

  (7)

  where ${\mathbf{d}}^{\prime }=|{\mathbf{x}}_{best}^t-{\mathbf{x}}^t|$ is the distance vector between the current individual and the best individual in the population, b is a constant (parameter of the algorithm), and r is a random number in the range $[-1,1]$.
• Whales move using both a shrinking encircling mechanism and a spiral-shaped movement. In the algorithm, this behavior is simulated by the formula

  $${\mathbf{x}}^{t+1}=\left\{\begin{array}{cc}{\mathbf{x}}_{best}^t-\mathbf{a}\ast \mathbf{d}\mathrm{Equation}\left(4\right),\hfill & \mathrm{for}p\le 0.5,\hfill \\ {\mathbf{d}}^{\prime }{e}^{br}cos\left(2\pi r\right)+{\mathbf{x}}_{best}^t\mathrm{Equation}\left(7\right),\hfill & \mathrm{for}p>0.5.\hfill \end{array}\right.$$

  (8)
  In the above formula $p\in [0,1]$, therefore, there is a $50\%$ chance of making each move. It is possible to control behavior of the algorithm by the level of probability p. Standard approach assumes a level of $0.5$.
• During the exploration phase, whales behave similarly to in Equation (4), with the difference that for the vector $\mathbf{a}$, we assume $\left|\mathbf{a}\right|>1$, which simulates the exploration phase and the whales move along with a random individual in the population, not the best individual. This stage of the algorithm is described mathematically by the formulas

  $$\mathbf{d}=|\mathbf{c}\ast {\mathbf{x}}_{rand}^t-{\mathbf{x}}^t|,$$

  (9)

  $${\mathbf{x}}^{t+1}={\mathbf{x}}_{rand}^t-\mathbf{a}\ast \mathbf{d},$$

  (10)

  where ${\mathbf{x}}_{rand}^t$ denotes the position of a random whale in the population t.
• The vector $\mathbf{w}$, whose values change from iteration to iteration in the range from 2 to 0, is responsible for the transition from the exploration phase to the exploitation phase. For values of the vector $\mathbf{w}$ in the range $(1,2]$, there is an exploration phase, while for the range $[0,1]$, there is an exploitation phase.
• The way of decreasing the value of the vector $\mathbf{w}$ in the Formula (5) and the level of probability appearing in the Formula (8) are tuning parameters of the WOA and can take different values.

Algorithm 1 presents the following steps of the WOA algorithm. Details about the described algorithm, as well as its application to selected problems, can be found in the articles [41,42,43].

Algorithm 1: Pseudocode of WOA

Figure 1 presents the block diagram of WOA.

3.2. Artificial Bee Colony

The Artificial Bee Colony (ABC) optimization algorithm draws its inspiration from the behavior of bees. This subsection is devoted to the mathematical description of the ABC algorithm.

Bees looking for a food source can be divided into two groups

• working bees—these are bees whose job is to look for a food source. Important information for these bees consists of the following: the distance between the hive and the food source, the direction the bee should follow to reach the food source and the amount of nectar in the source.
• bees unclassified—these are bees that search for new food sources. We can divide them into two groups: scouts and onlookers. Scouts, after leaving a food source, look for another one in randomly way, while viewers look for visited sources based on the information provided.

Bees communicate with each other to exchange information about their food source through dancing. Unclassified bees belonging to the onlookers group can choose the most relevant food source by observing the dance. Once a bee locates a food, it remembers its location and immediately starts exploring it. After obtaining a sufficient amount of nectar, it returns to the hive and unloads the collected food. Then it is faced with choosing one of three options:

• abandons the source, becomes an onlooker and watches the bees conveying information,
• transmits information through dance and recruits other bees,
• continues to explore on its own, without hiring other bees.

Based on the above assumptions, the bee algorithm was created, the main assumptions of which are as follows

• the locations of food sources correspond to potential solutions of the optimized problem.
• the quantity of nectar in the source corresponds to the quality of the solution.
• the number of working bees is equal to the number of onlookers, which is denoted by SN.
• Food sources are modified according to the formula

  $${\mathbf{v}}_i^t={\mathbf{x}}_i^t+\mathbf{r}\ast ({\mathbf{x}}_i^t-{\mathbf{x}}_k^t),$$

  (11)

  where $k\in \{1,2,\dots ,SN\}\wedge k\ne i$ is the random index, $\mathbf{r}$ is a vector of randomly generated numbers in the range $[-1,1]$, the operation ∗ means multiplying the vector “element by element”, ${\mathbf{x}}_i^t$ denotes i-th solution from the population in the t-th iteration and ${\mathbf{v}}_i^t$ denotes i-th the modified solution from the population in the t-th iteration.
  Compare positions ${\mathbf{v}}_i^t$ with  ${\mathbf{x}}_i^t$. If $\mathcal{F}\left({\mathbf{v}}_i^t\right)\le \mathcal{F}\left({\mathbf{x}}_i^t\right)$. The position of ${\mathbf{x}}_i^t$ in the population t is then replaced by ${\mathbf{v}}_i^t$. Otherwise, the position ${\mathbf{x}}_i^t$ remains in the population.
• Each item in population is assigned a probability according to the formula

  $$p_i=\frac{fit\left({\mathbf{x}}^i\right)}{{\displaystyle \sum _{j=1}^{SN}}fit\left({\mathbf{x}}^j\right)},i=1,2,\dots ,SN,$$

  (12)

  where

  $$fit\left({\mathbf{x}}_i\right)=\left\{\begin{array}{cc}\frac{1}{1+\mathcal{F}\left({\mathbf{x}}_i\right)},\hfill & \mathrm{if}\mathcal{F}\left({\mathbf{x}}_i\right)\ge 0,\hfill \\ 1+\left|\mathcal{F}\right({\mathbf{x}}_i\left)\right|,\hfill & \mathrm{if}\mathcal{F}\left({\mathbf{x}}_i\right)<0.\hfill \end{array}\right.$$

  (13)
• Each onlooker bee selects one source according to the probability $p_i$ and starts searching near it according to the Formula (11). Then the bee compares two locations—the new and previous one.
• If, after performing the previous step of the algorithm, any of the food sources have not changed their position, then they are omitted and replaced with a new random source

  $${\mathbf{x}}_i={\mathbf{x}}_{min}+\mathbf{r}\ast ({\mathbf{x}}_{max}-{\mathbf{x}}_{min}),i=1,2,\dots ,SN,$$

  (14)

  where $\mathbf{r}$ is a vector of random numbers in the range $[0,1]$, and ${\mathbf{x}}_{min}$, ${\mathbf{x}}_{max}$ are vectors of constraints, respectively lower and upper limits of the search space domain.

Algorithm 2 shows the following steps of the ABC method. More about this algorithm and its applications can be found in [44].

Algorithm 2: Pseudocode of ABC

Figure 2 presents the block diagram of the ABC algorithm.

4. Numerical Examples

This section presents computational examples illustrating the effectiveness of the methods described in this paper. The considered algorithms were compared with each other.

4.1. Example 1

Let us consider the equation

$$2{\phi }^{\prime }\left(x\right)+{\int }_{-1}^1\frac{\phi \left(t\right)}{t-x}dt=-\frac{x}{2},$$

which is considered important in the scientific literature (see, e.g., [39]), due to its occurrence in works related to heat conduction and radiation research.

Since we do not know the analytical (exact) solution to this problem, we plot absolute error $\Delta \left(x\right)$ based on the obtained approximate solution

$$\Delta \left(x\right)=\left|2{\phi }_n^{\prime }\left(x\right)+\lambda {\int }_{-1}^1\frac{{\phi }_n\left(t\right)}{t-x}dt-f\left(x\right)\right|,$$

(15)

which, in this case, comes down to the following form

$$\Delta \left(x\right)=\left|2{\phi }_n^{\prime }\left(x\right)+\frac{x}{2}+{\int }_{-1}^1\frac{{\phi }_n\left(t\right)}{t-x}dt\right|.$$

Figure 3 presents the errors of the solution obtained from the method using the Taylor series and for different values of the number of terms n: for $n=8$, left side and for $n=10$, right side.

The error curves in both cases are similar in nature, i.e., the highest error values are achieved at the edges of the interval. Taking $n=10$ reduces the size of the error by one order compared to $n=8$.

Now, we present a solution to this example by using the WOA and ABC algorithms. Since we do not know the exact solution, an approximate solution was adopted in the form of polynomials of degree 3, 4 and 5. Figure 4, Figure 5 and Figure 6 show the distribution of errors in the approximate solution for WOA and ABC. For example Figure 5 shows the distribution of errors in the case of a fourth-degree polynomial. As we can see in Figure 5, smaller errors were obtained in the case of the WOA solution. The highest errors are obtained in the case of ABC, but they are also at a satisfactory level. In the case of polynomials of degree 3 and 5, analogous error graphs are included in Figure 4 and Figure 6.

4.2. Example 2

Let us consider equation

$$2{\phi }^{\prime }\left(x\right)+{\int }_{-1}^1\frac{\phi \left(t\right)}{t-x}dt=f\left(x\right),$$

where $\lambda =1$,

$$f\left(x\right)=-\frac{2}{3}(17+38x-12x^2+30x^3)+(2-3x-4x^2+x^3-2x^4)ln\left|\frac{1-x}{1+x}\right|.$$

Exact solution for this problem has following form

$$\phi \left(x\right)=2-3x-4x^2+x^3-2x^4.$$

Since we know the exact solution, we do not need to check the accuracy of the approximate solution using (15). In this case, the absolute error ${\Delta }_n^{\prime }$ is defined more simply

$${\Delta }_n^{\prime }\left(x\right)=\left|\phi \left(x\right)-{\phi }_n\left(x\right)\right|.$$

(16)

First, the Taylor series method was used to solve the equation. The plots of the exact solution $\phi$ and the approximate solution ${\phi }_n$ and the errors of these approximations are provided in Figure 7 (for $n=4$) and Figure 8 (for $n=6$). As we can see, only for $n=6$, the errors of the obtained solution are minimal, approximately ${10}^{-15}$. In the case of $n=4$, the approximate solution maintains properties of the exact solution, but the differences between the two functions are significant.

Now, we discuss the obtained solutions when using ABC and WOA algorithms. It is assumed that the form of the sought function would be a polynomial of the fourth degree, i.e., $\overline{\phi }\left(x\right)=a_0+a_{1}x+a_2{x}^2+a_3{x}^3+a_4{x}^4$. Therefore, five coefficients are sought: $a_0,a_1,\dots ,a_4$. As a result of using the metaheuristic algorithms, the following results were obtained (see

Table 1. Sought coefficient values obtained using WOA and ABC methods. $a_i(i=0,1,2,3,4)$—obtained coefficient value, $\mathcal{F}$—objective function value.

	${\mathit{a}}_0$	${\mathit{a}}_1$	${\mathit{a}}_2$	${\mathit{a}}_3$	${\mathit{a}}_4$	$\mathcal{F}$
WOA	$2.6771$	$-3.0077$	$-3.9007$	$1.0142$	$-1.7346$	$0.1631$
ABC	$2.0051$	$-2.9996$	$-3.9976$	$0.9996$	$-1.9981$	$0.0034$

). As we can see, the value of the fitness function $\mathcal{F}$ in the case of ABC is much smaller than in the case of WOA. This proves that the ABC algorithm works better for this example, which is especially visible when comparing results of finding parameters $a_0$ and $a_4$. The exact values of these coefficients are 2 and $-2$, respectively. WOA returned values of $a_0\approx 2.67$, $a_4\approx -1.73$, while ABC returned $a_0\approx 2.00$, $a_4\approx -1.99$. Figure 9 and Figure 10 present the solutions found by WOA and ABC, respectively, compared with the exact solution, along with the distribution of errors in the domain. It is clear that ABC solves the problem very well, while the result obtained from WOA is not so good.

4.3. Example 3

In this example, we again consider the following equation

$$2{\phi }^{\prime }\left(x\right)+{\int }_{-1}^1\frac{\phi \left(t\right)}{t-x}dt=f\left(x\right),$$

where function f has the form

$$f\left(x\right)=-\frac{2\left(4x^4-5x^2+1\right)}{\sqrt{1-x^2}}-\frac{1}{8}\pi \left(8x^4-12x^2+3\right).$$

The exact solution to this problem is a function

$$\phi \left(x\right)=x(x^2-1)\sqrt{1-x^2}.$$

As in the case of the previous examples, the solution is sought in the form of a polynomial. To estimate the absolute errors Formula (16) was used.

Plots of the exact solution $\phi$ and the approximate solution ${\phi }_n$ and the errors of these approximations (for the method with Taylor series expansion) are presented in Figure 11 (for $n=6$), and Figure 12 (for $n=16$). We can notice that inside the considered area the approximation, errors are minimal and when closing to the ends of interval $[-1,1]$, these errors increase slightly.

We now present solutions obtained using metaheuristic optimization algorithms. The solution is sought in the form

$$\overline{\varphi }\left(x\right)=\sum _{i=0}^n{a}_i{x}^i.$$

In the case of metaheristic algorithms, it is assumed $n=3,4,5$. As can be seen in the Figure 13, Figure 14, Figure 15, Figure 16, Figure 17 and Figure 18, in the case of metaheuristic algorithms, the best results are obtained for the ABC algorithm and the $n=3$ case. Definitely, worse results are obtained in the case of WOA.

5. Conclusions

In the article, two different approaches for solving the integro-differential equation with the Cauchy kernel were compared. The first method used Taylor series expansion of functions, and the second group consists of two selected metaheuristic optimization algorithms—the classic ABC algorithm and the younger WOA algorithm. Both are used in recent papers in many problems.

Numerical examples show that both methods are effective in solving this important and difficult task. For this purpose, three numerical examples of varying complexity were presented, both for equations in which we know the exact solution and for those where such a solution is not known. The research has shown that the first of the discussed methods works better with solving the considered type of equation and generates noticeably smaller errors. Heuristic methods also solve this type of task.

To sum up, we can say that the main advantages of the method based on the Taylor series is higher in accuracy than in case of the metaheuristic approach. A slight disadvantage of this method is the specific form of the system of equations being solved, which may be indeterminate for certain degrees of the sought polynomial. However, the main advantages of metaheuristic methods are their universality. Once the method has been implemented, the form of the equation being solved does not affect the difficulty of adopting it to the algorithm. It only requires specifying the class of functions (e.g., polynomials) in which the solution will be sought. Then, an error function is constructed that depends on the coefficients of the function being sought. By finding the minimum of the error function, we obtain an approximate solution. The disadvantages of these methods include the total running time. In order to effectively search for solutions, it is necessary to spend time in advance to ’tune’ the algorithm parameters. We also obtained errors in the approximate solution that were higher than in the case of the method with the Taylor series. But in many engineering problems, the obtained solution is sufficiently satisfactory.

It should be mentioned here that thanks to the symmetric domain and the appropriate form of the f function, certain conclusions can be drawn in some cases. For example, in example 1, the given function f is even (its graph is symmetrical about the Y axis), so we can expect that the solution will be an odd function (its graph is symmetrical about the origin), and in example 3, we have the opposite situation—the graph of the f function is symmetrical with respect to the origin of the coordinate system, so it is probable that the graph of the sought function $\phi$ will be symmetrical with respect to the Y axis. In both cases, this assumption turned out to be true, and the correct assumption of the truth of the hypothesis about the form of the function $\phi$ can significantly shorten the time of searching for this function and improve the quality of the solution.

Our future research plans focus on exploring a broader class of meta-heuristic algorithms. Another idea is to take into account hybrid methods, i.e., combining a heuristic method with a certain deterministic method (e.g., Nelder–Mead or Hooke–Jeeves). We also plan to generalize the class of the considered integro-differential equation.