An Improved Bernoulli Collocation Method for Solving Volterra Integral Equations

1. Introduction

Most of the problems of science and technology can be treated with the aid of theories of ordinary and partial differential equations (PDEs). However, there are better methods called the theory of integral equations to solve these problems. Years ago, these types of equations were a hot topic in the minds of mathematicians. If the unknown function is given inside the integration symbol in an equation, it is named an integral equation, which is regarded as a common type of functional equation. The theory of integral equations is one of the most operative mathematical tools for solving problems in pure and applied mathematics, mathematical physics, mechanical vibrations, and fields related to science and technology. However, the real development of the theory of integral equations started with the attempts of the Italian mathematician Vito Volterra and the Swedish mathematician Fredholm. Volterra was the first person who realized the significance of the theory of integral equations and systematically paid attention to it. In 1896, he presented the first general scheme for solving a category of linear integral equations characterized by a variable that appears at the upper bound of the integral. Interest in VIEs has been growing in recent years. These equations arise in many physical applications, for example see [1,2,3,4,5,6,7,8,9,10,11,12]. Here, a numerical scheme based on Bernoulli polynomials is developed to solve VIEs of the second kind.

Bernoulli polynomials have a significant role in many areas of mathematical analysis, like the theory of modular forms [14], the theory of distributions in p-adic analysis [13], the polynomial expansions of analytic functions [15], and so on. Some applications of these polynomials in mathematical physics are related with the theory of the KdV equation [16], solving Lamé equation [17], and in the field of vertex algebra [18].

The collocation method as a numerical method for solving all kinds of functional equations and real world problems has always been the interest and attention of mathematicians. Many authors have presented different kinds of collocation method during past decades. For example, Doha et al. [19] proposed a Jacobi–Gauss–Lobatto collocation scheme, used in relationship with the fourth order implicit Runge–Kutta method as a numerical algorithm for approximating solutions of nonlinear Schrödinger equations (NLSE). Nemati [20] solved Volterra–Fredholm integral equations applying Legendre collocation method. His method is based upon shifted Legendre polynomials approximation. Then, utilizing the shifted Gauss–Legendre as collocation nodes reduces the Volterra–Fredholm integral equations to the solution of a matrix equation. Mirzaee and Hoseini [21] proposed a new matrix method on the basis of Fibonacci polynomials and collocation points for numerically solving the Volterra–Fredholm integral equations. Ren and Tian [22] proposed a scheme to solve a boundary value problem for Kirchhoff type of nonlinear integro-differential equation numerically. Gouyandeh et al. [23] by using the Tau-Collocation method approximated solution of the nonlinear Volterra-Fredholm-Hammerstein integral equations. Aziz et al. [24], based on Haar wavelet, provided a new collocation method for numerical solution of 3D elliptic PDEs with Dirichlet boundary conditions. Çelik [25], by utilizing Chebyshev wavelet collocation method, studied free vibration problems of non-uniform Euler–Bernoulli beam under different supporting conditions. Samadyar and Mirzaee [26] presented orthonormal BCM to approximate the linear singular stochastic Itô-VIEs. Biçer and Yalçinbaş [27] an approximate solution of the telegraph equation applying BCM. Alijani et al. [28] investigated systems of fuzzy fractional differential equations numerically with a lateral type of the Hukuhara derivative and its generalization. Wang et al. [29] presented a new collocation method for evaluating the 2D elliptic PDEs. Singh [30] used Jacobi collocation method to solve the fractional advection-dispersion equation arising in porous media. Kumbinarasaiah et al. [31] presented an integration operational matrix applying the Bernoulli wavelet and suggested a new scheme named the Bernoulli wavelet collocation method (BWCM). In [32], the authors investigated a space-time Sinc-collocation method for treating the fourth-order nonlocal heat model appearing in viscoelasticity. Laib et al. [33], based on the using Taylor polynomials, suggested an algorithm to construct a collocation solution for approximating the solution of 2D-VIEs. Wang et al. [34], by utilizing the zeros of Chebyshev polynomial as collocation points, proposed a new collocation method to solve the second kind VIE.

In this work, based on Bernoulli polynomials, a new collocation method is suggested to approximate numerically VIEs of the second kind. The main goal of the suggested scheme is to improve the results of the classic BCM by dividing the interval into some sub-intervals and considering the collocation points on each of them. Here, the zeros of the SCPs are considered as collocation points. Then, BCM is applied step by step from the first sub-interval to the last one. By this process, a system of algebraic equations is attained for each sub-interval that could be easily solved using computing software. At last, the approximate solution is obtained as a piece-wise function. This idea is very effective. Although, we have tested this idea on BCM in this work but we guess that all collocation methods mentioned above can be improved by this idea.

The rest of the paper is organized as follows. In section 2, Bernoulli polynomials are introduced and their features are stated. In section 3, the proposed scheme is explained. The convergence analysis is discussed in in section 4. In section 5, numerical results are presented. At last, conclusions are given in section 5.

2. Bernoulli Polynomials

The traditional Bernoulli polynomials $B_n\left(t\right)$ is often characterized by the following exponential generating functions [35]:

$$\frac{s{e}^{sx}}{e^s-1}={\sum }_{z=0}^{\mathrm{\infty }} B_z\left(x\right)\frac{s^k}{z!}$$

(1)

The Bernoulli polynomials of $n$th degree are determined in the interval $\left[\mathrm{0,1}\right]$ as follows [36]

$$B_N\left(x\right)=\sum _{z=0}^N \left(\begin{array}{l}N\\ z\end{array}\right)B_z{x}^{N-\mathrm{z}},$$

(2)

where $B_z=B_z\left(0\right)$ is the Bernoulli number for each $k=\mathrm{0,1},\dots N$. We write some polynomials like

$$\begin{array}{ll}& B_0\left(x\right)=1,\\ & B_1\left(x\right)=x-\frac{1}{2},\\ & B_2\left(x\right)=x^2-x+\frac{1}{6},\\ & B_3\left(x\right)=x^3-\frac{3}{2}{x}^2+\frac{1}{2}x,\\ & B_4\left(x\right)=x^4-2x^3+x^2-\frac{1}{30},\\ & B_5\left(x\right)=x^5-\frac{5}{2}{x}^4+\frac{5}{3}{x}^3-\frac{1}{6}x,\\ & B_6\left(x\right)=x^6-3x^5+\frac{5}{2}{x}^4-\frac{1}{2}{x}^2+\frac{1}{42},\\ & B_7\left(x\right)=x^7-\frac{7}{2}{x}^6+\frac{7}{2}{x}^5-\frac{7}{6}{x}^3+\frac{1}{6}x.\end{array}$$

Bernoulli polynomials have following important properties [37].

• $\frac{d}{dx}{B}_N\left(x\right)=N{B}_{N-1}\left(x\right),N\ge 1$.
• $B_N(x+1)-B_N\left(x\right)=N{x}^{N-1}$.
• ${\int }_0^1 B_N\left(x\right)dx=0,N\ge 1$.
• ${\int }_0^1 B_N\left(x\right)B_M\left(x\right)dx=(-1{)}^{N-1}\frac{M!N!}{(M+N)!}{B}_{N+M}$.
• ${\int }_a^x B_N\left(t\right)dt=\frac{B_{N+1}\left(x\right)-B_{N+1}\left(a\right)}{(N+1)}$.
• $B_N(1-x)=(-1{)}^N{B}_N(x)$.

2.1. Function Approximation

A square integrable function $u\left(t\right)$ can be expressed according to of Bernoulli polynomials as

$$u\left(t\right)=\sum _{i=0}^{\mathrm{\infty }} U_i{B}_i\left(t\right),$$

and the truncated series is

$$\tilde{u}\left(x\right)\simeq \sum _{i=0}^M U_i{B}_i\left(x\right)=U^{T}B\left(x\right),$$

(3)

where

$$U={\left[\begin{array}{lllll}{U}_0& U_1& U_2& \dots & U_M\end{array}\right]}^T,$$

(4)

is the vector of unknown coefficients and

$$B\left(x\right)={\left[\begin{array}{lllll}{B}_0\left(x\right)& B_1\left(x\right)& B_2\left(x\right)& \dots & B_M\left(x\right)\end{array}\right]}^T,$$

(5)

is the vector of Bernoulli polynomials.

3. Description of the Proposed Scheme

In this section, the proposed scheme is described to deal with the VIEs of the second kind. The main purpose of the suggested method is to improve the results of the classic collocation method by dividing the interval into some sub-intervals and applying the collocation method in each of them.

3.1. Solving VIEs of the Second Kind by BCM

Regard the following VIE of the second kind.

$$u\left(x\right)=f\left(x\right)+{\int }_0^x k\left(x,t\right)N\left(u\left(t\right)\right)dt, x\in \left[\mathrm{0,1}\right],$$

(6)

where $u$ is the unknown function, while $f$ and the kernel $k$ are known functions, and $N$ is a given continuous function which is nonlinear with respect to $u$.

Substituting Eq. (3) in Eq. (6) leads to

$$\sum _{i=1}^M{U_{i}B}_i\left(x\right)=f\left(x\right)+{\int }_0^{x}k\left(x,t\right)N\left(\sum _{i=1}^M{U_{i}B}_i\left(t\right)\right)dt.$$

(7)

Now, we use the collocation method to determine unknowns $U_i$, $i=1,\dots ,M$. Let $c_k, k=1,\dots ,M$ be the collocation points. Here, we apply the zeros of SCPs of degree $M$ in the interval $\left[\mathrm{0,1}\right]$ as collocation points. For example for $M=3$, the collocation points are $c_1=0.0670$, $c_2=0.5$, $c_3=$

0.9330. Then, we will a system of nonlinear algebraic equation that could be easily solved computer software.

In the case of $N\left(u\left(t\right)\right)=u\left(t\right)$, that is the equation is linear, by rearranging this equation in terms of $U_i$, we have

$$\sum _{i=1}^M\left(B_i\left(x\right)-{\int }_0^{x}k\left(x,t\right)B_i\left(t\right)dt\right)U_i=f\left(x\right).$$

(8)

By substituting collocation points $c_k, k=1,\dots ,M$ in Eq. (8) we get

$$\sum _{i=1}^M\left(B_i\left(c_k\right)-{\int }_0^{c_k}k\left(c_k,t\right)B_i\left(t\right)dt\right)U_i=f\left(c_k\right), k=1,\dots ,M.$$

(9)

The last equation is a system of $M$ algebraic equations with unknown coefficients $U_i$, $i=1,\dots ,M$ that could be stated in the following matrix form.

$AU=F,$

In which

$A=\left[\begin{array}{ccc}{B}_1\left(c_1\right)-{\int }_0^{c_1}k\left(c_1,t\right)B_1\left(t\right)dt& \cdots & B_M\left(c_1\right)-{\int }_0^{c_1}k\left(c_1,t\right)B_M\left(t\right)dt\\ ⋮& \ddots & ⋮\\ B_1\left(c_M\right)-{\int }_0^{c_M}k\left(c_M,t\right)B_1\left(t\right)dt& \cdots & B_M\left(c_M\right)-{\int }_0^{c_M}k\left(c_M,t\right)B_M\left(t\right)dt\end{array}\right],$vector $U$ is defined as Eq. (4) and

$F_1={\left[f\left(c_1\right),f\left(c_2\right),\dots ,f\left(c_M\right)\right]}^T.$

3.2. Solving VIEs of the Second Kind by IOBCM

In order to apply the idea of the suggested scheme, first we divide the interval $\left[\mathrm{0,1}\right]$ into $N$ sub-intervals as $I_j=\left[\left(j-1\right)h,jh\right]$ where $h=\frac{1}{N}$. Then, we consider the approximation of $u$ by Bernoulli polynomials of degree $M$ in each sub-interval as follows.

$$u\left(x\right)\simeq \sum _{i=1}^M{U}_{i,j}{B}_i\left(x\right)=B\left(x\right)U_j, x\in I_j, j=1,\dots ,N,$$

(10)

where $U_{i,j}$, $i=1,\dots ,M, j=1,\dots ,N$, are unknown coefficients to be determined and $B\left(x\right)$ is define as Eq. (4) and

$$U_j={\left[U_{1,j},U_{2,j},\dots ,U_{M,j}\right]}^T, j=1,\dots ,N.$$

(11)

According to Eq. (10), the approximate solution is considered as a piece-wise function in the proposed method.

In general, there are $MN$ unknowns, $U_{\mathrm{1,1}},U_{\mathrm{2,1}},\dots ,U_{M,1},\dots , U_{1,N},U_{\mathrm{2,1}},\dots ,U_{M,N}$, to be determined. To find these unknowns, we do as follows.

For finding the unknowns $U_{\mathrm{1,1}},U_{\mathrm{2,1}},\dots ,U_{M,1}$, suppose that $x\in I_1$. Then, according to Eq. (10) the approximate solution in the interval $I_1=\left[0,h\right]$ is

$$u\left(x\right)\simeq \sum _{i=1}^M{U}_{i,1}{B}_i\left(x\right).$$

(12)

Substituting Eq. (12) in Eq. (6) leads to

$$\sum _{i=1}^M{U_{i,1}B}_i\left(x\right)=f\left(x\right)+{\int }_0^{x}k\left(x,t\right)N\left(\sum _{i=1}^M{U_{i,1}B}_i\left(t\right)\right)dt.$$

(13)

We use the collocation method to determine unknowns $U_{i,1},$ $i=1,\dots ,M$. Let $c_{k,1}, k=1,\dots ,M$ be the zeros of SCPs of degree $M$ in the interval $\left[0,h\right]$ as collocation points. For example for $M=3$, and $N=4$ the collocation points are $c_{\mathrm{1,1}}=0167$, $c_{\mathrm{2,1}}=0.1250$, $c_{\mathrm{3,1}}=0.2333$, $c_{\mathrm{1,2}}=0.2667$, $c_{\mathrm{2,2}}=0.3750$, $c_{\mathrm{3,2}}=0.4833$, $c_{\mathrm{1,3}}=0.5167$, $c_{\mathrm{2,3}}=0.6250$, $c_{\mathrm{3,3}}=0.7333$, $c_{\mathrm{1,4}}=0.7667$, $c_{\mathrm{2,4}}=0.8750$, $c_{\mathrm{3,4}}=0.9833$. Then, a system of nonlinear algebraic equation is produced that could be easily solved.

In the case of $N\left(u\left(t\right)\right)=u\left(t\right)$, that is the equation is linear, equation (13) could be stated as follows.

$$\sum _{i=1}^M\left(B_i\left(x\right)-{\int }_0^{x}k\left(x,t\right)B_i\left(t\right)dt\right)U_{i,1}=f\left(x\right).$$

(14)

By substituting collocation points $c_{k,1}, k=1,\dots ,M$ in Eq. (13) we will have

$$\sum _{i=1}^M\left(B_i\left(c_{k,1}\right)-{\int }_0^{c_{k,1}}k\left(c_{k,1},t\right)B_i\left(t\right)dt\right)U_{i,1}=f\left(c_{k,1}\right).$$

(15)

This equation is a system of $M$ algebraic equations with unknown coefficients $U_{i,1}$, $i=1,\dots ,M$ which can be written in the following matrix form.

$A_1{U}_1=F_1,$

In which

$A_1=\left[\begin{array}{ccc}{B}_1\left(c_{\mathrm{1,1}}\right)-{\int }_0^{c_{\mathrm{1,1}}}k\left(c_{\mathrm{1,1}},t\right)B_1\left(t\right)dt& \cdots & B_M\left(c_{\mathrm{1,1}}\right)-{\int }_0^{c_{\mathrm{1,1}}}k\left(c_{\mathrm{1,1}},t\right)B_M\left(t\right)dt\\ ⋮& \ddots & ⋮\\ B_1\left(c_{M,1}\right)-{\int }_0^{c_{M,1}}k\left(c_{M,1},t\right)B_1\left(t\right)dt& \cdots & B_M\left(c_{M,1}\right)-{\int }_0^{c_{M,1}}k\left(c_{M,1},t\right)B_M\left(t\right)dt\end{array}\right],$vector $U_1$ is defined as Eq. (11) for $j=1$ and

$F_1={\left[f\left(c_{\mathrm{1,1}}\right),f\left(c_{\mathrm{2,1}}\right),\dots ,f\left(c_{M,1}\right)\right]}^T.$

By using the coefficients $U_{i,1}$, $i=1,\dots ,M$ which have been determined in the previous stage, we can find the unknowns $U_{i,2}$, $i=1,\dots ,M$ as follows.

Suppose that $x\in I_2$. Then, Substituting Eq. (10) for $j=2$ in Eq. (6) gives

$$\begin{array}{cc}\hfill {\displaystyle \sum _{i=1}^M}{B}_i\left(x\right)U_{i,2}& =f\left(x\right)+{\int }_0^{h}k\left(x,t\right)\left({\displaystyle \sum _{i=1}^M}{B}_i\left(t\right)U_{i,1}\right)dt\hfill \\ & +{\int }_h^{x}k\left(x,t\right)\left({\displaystyle \sum _{i=1}^N}{B}_i\left(t\right)U_{i,2}\right)dt.\hfill \end{array}$$

(16)

This equation can be written as follows.

$$\begin{gathered} \sum _{i=1}^M\left(B_i\left(x\right)-{\int }_h^{x}k\left(x,t\right)B_i\left(t\right)dt\right)U_{i,2}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em} \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}=f\left(x\right)+{\int }_0^{h}k\left(x,t\right)\left(\sum _{i=1}^M{U_{i,1}B}_i\left(t\right)\right)dt. \end{gathered}$$

(17)

Let $c_{k,2}, k=1,\dots ,M$ be the zeros of SCPs of degree $M$ in the interval $\left[h,2h\right]$ as collocation points. Substituting the collocation points in Eq. (17) leads to the following relation.

$$\sum _{i=1}^M\left(B_i\left(c_{k,2}\right)-{\int }_h^{c_{k,2}}k\left(c_{k,2},t\right)B_i\left(t\right)dt\right)U_{i,2}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}$$

(18)

This equation is a system of $M$ algebraic equations with unknown coefficients $U_{i,2}$, $i=1,\dots ,M$. This system could be stated in the following matrix form.

$A_2{U}_2=F_2,$

where

$A_2=\left[\begin{array}{ccc}{B}_1\left(c_{\mathrm{1,2}}\right)-{\int }_h^{c_{\mathrm{1,2}}}k\left(c_{\mathrm{1,2}},t\right)B_1\left(t\right)dt& \cdots & B_M\left(c_{\mathrm{1,2}}\right)-{\int }_h^{c_{\mathrm{1,2}}}k\left(c_{\mathrm{1,2}},t\right)B_M\left(t\right)dt\\ ⋮& \ddots & ⋮\\ B_1\left(c_{M,2}\right)-{\int }_h^{c_{M,2}}k\left(c_{M,2},t\right)B_1\left(t\right)dt& \cdots & B_M\left(c_{M,2}\right)-{\int }_h^{c_{M,2}}k\left(c_{M,2},t\right)B_M\left(t\right)dt\end{array}\right],$vector $U_2$ is defined as Eq. (11) for $j=2$ and

$F_2=\left(\begin{array}{c}f\left(c_{\mathrm{1,2}}\right)+{\int }_0^{h}k\left(c_{\mathrm{1,2}},t\right)\left({\displaystyle \sum _{i=1}^M}{U_{i,2}B}_i\left(t\right)\right)dt\\ ⋮\\ f\left(c_{M,2}\right)+{\int }_0^{h}k\left(c_{M,2},t\right)\left({\displaystyle \sum _{i=1}^M}{U_{i,2}B}_i\left(t\right)\right)dt\end{array}\right).$

Now, we can present a general formula to find the unknowns in in the interval $I_j=\left[\left(j-1\right)h,jh\right]$.

Suppose that $x\in I_j$. Then, according to Eq. (12) the approximate solution in the interval $I_j$ is

$$u\left(x\right)\simeq \sum _{i=1}^M{U}_{i,j}{B}_i\left(x\right).$$

(19)

By substituting Eq. (19) in Eq. (6) we will have

$$\begin{array}{cc}\hfill {\displaystyle \sum _{i=1}^M}{B}_i\left(x\right)U_{i,j}=& f\left(x\right)+{\int }_0^{h}k\left(x,t\right)\left({\displaystyle \sum _{i=1}^M}{U_{i,1}B}_i\left(t\right)\right)dt+\dots \hfill \\ & +{\int }_{(j-2)h}^{(j-1)h}k\left(x,t\right)\left({\displaystyle \sum _{i=1}^M}{U_{i,j-1}B}_i\left(t\right)\right)dt+{\int }_{(j-1)h}^{x}k\left(x,t\right)\left({\displaystyle \sum _{i=1}^M}{U_{i,j}B}_i\left(t\right)\right)dt,\hfill \end{array}$$

This equation would be stated as follows.

$$\begin{array}{c}{\displaystyle \sum _{i=1}^M}\left(B_i\left(x\right)-{\int }_{(j-1)h}^{x}k\left(x,t\right)B_i\left(t\right)dt\right)U_{i,j}\hfill \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}=f\left(x\right)+{\int }_0^{h}k\left(x,t\right)\left({\displaystyle \sum _{i=1}^M}{U_{i,1}B}_i\left(t\right)\right)dt+\dots \hfill \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}+{\int }_{(j-2)h}^{(j-1)h}k\left(x,t\right)\left({\displaystyle \sum _{i=1}^M}{U_{i,j-1}B}_i\left(t\right)\right)dt.\hfill \end{array}$$

(20)

Let $c_{k,j}, k=1,\dots ,M$ be the zeros of SCPs of degree $M$ in the interval $I_j$ as collocation points. Substituting the collocation points in Eq. (20) yields

$$\begin{array}{c}{\displaystyle \sum _{i=1}^M}\left(B_i\left(c_{k,j}\right)-{\int }_{(j-1)h}^{c_{k,j}}k\left(c_{k,j},t\right)B_i\left(t\right)dt\right)U_{i,j}\hfill \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}=f\left(c_{k,j}\right)+{\int }_0^{h}k\left(c_{k,j},t\right)\left({\displaystyle \sum _{i=1}^M}{U_{i,1}B}_i\left(t\right)\right)dt+\dots \hfill \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}+{\int }_{(j-2)h}^{(j-1)h}k\left(c_{k,j},t\right)\left({\displaystyle \sum _{i=1}^M}{U_{i,j-1}B}_i\left(t\right)\right)dt.\hfill \end{array}$$

This equation is a system of $M$ algebraic equations with unknown coefficients $U_{i,j}$, $i=1,\dots ,M$. It has the following matrix form.

$A_j{U}_j=F_j,$

where

$A_j=\left[\begin{array}{ccc}{B}_1\left(c_{1,j}\right)-{\int }_{(j-1)h}^{c_{1,j}}k\left(c_{1,j},t\right)B_1\left(t\right)dt& \cdots & B_M\left(c_{1,j}\right)-{\int }_{(j-1)h}^{c_{1,j}}k\left(c_{1,j},t\right)B_M\left(t\right)dt\\ ⋮& \ddots & ⋮\\ B_1\left(c_{M,j}\right)-{\int }_{(j-1)h}^{c_{M,j}}k\left(c_{M,j},t\right)B_1\left(t\right)dt& \cdots & B_M\left(c_{M,j}\right)-{\int }_{(j-1)h}^{c_{M,j}}k\left(c_{M,j},t\right)B_M\left(t\right)dt\end{array}\right],$vector $U_j$ is defined as Eq. (11) and

$F_j=\left(\begin{array}{c}f\left(c_{1,j}\right)+{\displaystyle \sum _{r=1}^{j-1}}{\int }_{(r-1)h}^{rh}k\left(c_{1,j},t\right)\left({\displaystyle \sum _{i=1}^M}{U_{i,r}B}_i\left(t\right)\right)dt\\ ⋮\\ f\left(c_{M,j}\right)+{\displaystyle \sum _{r=1}^{j-1}}{\int }_{(r-1)h}^{rh}k\left(c_{M,j},t\right)\left({\displaystyle \sum _{i=1}^M}{U_{i,r}B}_i\left(t\right)\right)dt\end{array}\right).$Therefore, if we continue this process step by step and find the unknowns in the next subintervals using the coefficients which have been determined so far, then all unknowns will be obtained finally.

For finding the unknowns in the last interval, suppose that $x\in I_N$. Then according to Eq. (10) the approximate solution in this interval is

$$u\left(x\right)\simeq \sum _{i=1}^M{U}_{i,N}{B}_i\left(x\right).$$

(21)

Substituting Eq. (21) in Eq. (6) gives

$\begin{array}{cc}\hfill {\displaystyle \sum _{i=1}^M}{B}_i\left(x\right)U_{i,N}=& f\left(x\right)+{\int }_0^{h}k\left(x,t\right)\left({\displaystyle \sum _{i=1}^M}{U_{i,1}B}_i\left(t\right)\right)dt+\dots \hfill \\ & +{\int }_{(N-2)h}^{(N-1)h}k\left(x,t\right)\left({\displaystyle \sum _{i=1}^M}{U_{i,N-1}B}_i\left(t\right)\right)dt\hfill \\ & +{\int }_{(N-1)h}^{x}k\left(x,t\right)\left({\displaystyle \sum _{i=1}^M}{U_{i,N}B}_i\left(t\right)\right)dt.\hfill \end{array}$

This equation can be written in terms of $U_{i,N}$ as follows.

$$\begin{array}{c}{\displaystyle \sum _{i=1}^M}\left(B_i\left(x\right)-{\int }_{(N-1)h}^{x}k\left(x,t\right)B_i\left(t\right)dt\right)U_{i,N}\hfill \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}=f\left(x\right)+{\int }_0^{h}k\left(x,t\right)\left({\displaystyle \sum _{i=1}^M}{U_{i,1}B}_i\left(t\right)\right)dt+\dots \hfill \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}+{\int }_{(N-2)h}^{(N-1)h}k\left(x,t\right)\left({\displaystyle \sum _{i=1}^M}{U_{i,N-1}B}_i\left(t\right)\right)dt.\hfill \end{array}$$

(22)

Let $c_{k,N}, k=1,\dots ,M$ be the zeros of SCPs of degree $M$ in the interval $I_N$ as collocation points. By substituting these points in Eq. (22) we will have

$\begin{array}{c}{\displaystyle \sum _{i=1}^M}\left(B_i\left(c_{k,N}\right)-{\int }_{(N-1)h}^{c_{k,N}}k\left(c_{k,N},t\right)B_i\left(t\right)dt\right)U_{i,N}\hfill \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}=f\left(c_{k,N}\right)+{\int }_0^{h}k\left(c_{k,N},t\right)\left({\displaystyle \sum _{i=1}^M}{U_{i,1}B}_i\left(t\right)\right)dt+\dots \hfill \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}+{\int }_{(N-2)h}^{(N-1)h}k\left(c_{k,N},t\right)\left({\displaystyle \sum _{i=1}^M}{U_{i,N-1}B}_i\left(t\right)\right)dt.\hfill \end{array}$

This equation is a system of $M$ algebraic equations with unknown coefficients $U_{i,N}$, $i=1,\dots ,M$ which can be stated in a matrix form as follows.

$A_N{U}_N=F_N,$

In which

$A_N=\left[\begin{array}{ccc}{B}_1\left(c_{1,N}\right)-{\int }_{(N-1)h}^{c_{1,N}}k\left(c_{1,N},t\right)B_1\left(t\right)dt& \cdots & B_M\left(c_{M,N}\right)-{\int }_{(N-1)h}^{c_{1,N}}k\left(c_{1,N},t\right)B_M\left(t\right)dt\\ ⋮& \ddots & ⋮\\ B_1\left(c_{M,N}\right)-{\int }_{(N-1)h}^{c_{M,N}}k\left(c_{M,N},t\right)B_1\left(t\right)dt& \cdots & B_M\left(c_{M,N}\right)-{\int }_{(N-1)h}^{c_{M,N}}k\left(c_{M,N},t\right)B_M\left(t\right)dt\end{array}\right],$vector $U_N$ is defined as Eq. (11) for $j=N$ and

$F_N=\left(\begin{array}{c}f\left(c_{1,N}\right)+{\displaystyle \sum _{r=1}^{N-1}}{\int }_{(r-1)h}^{rh}k\left(c_{1,N},t\right)\left({\displaystyle \sum _{i=1}^M}{U_{i,r}B}_i\left(t\right)\right)dt\\ ⋮\\ f\left(c_{M,N}\right)+{\displaystyle \sum _{r=1}^{N-1}}{\int }_{(r-1)h}^{rh}k\left(c_{M,N},t\right)\left({\displaystyle \sum _{i=1}^M}{U_{i,r}B}_i\left(t\right)\right)dt\end{array}\right).$Eventually, we can calculate the solution by the following piece wise function.

$$\tilde{u}=\left\{\begin{array}{c}{\displaystyle \sum _{i=1}^M}{U}_{i,1}{B}_i\left(x\right), x\in I_1,\\ {\displaystyle \sum _{i=1}^M}{U}_{i,2}{B}_i\left(x\right), x\in I_2,\\ ⋮\\ {\displaystyle \sum _{i=1}^M}{U}_{i,N}{B}_i\left(x\right), x\in I_N,\end{array}\right.$$

(23)

where $\tilde{u}$ is the approximation of the exact solution $u$.

4. Convergence Analysis

Here, the convergence of the suggested scheme is discussed. For this purpose, we first recall one of the important theorems related to the residual interpolation error by Chebyshev nodes.

Theorem 1. Let $u$ be a sufficiently smooth function on $I=[0, 1]$ and ${\mathsf{\Pi }}_M$ be the space of polynomials of order $M.$ Also, let $P_M\in {\mathsf{\Pi }}_M$ be the interpolating polynomials of $u$ at points $c_1,\dots ,c_{M+1}$ which are the zeros of the SCP of degree $M+1$ on $I$. Then, the following relation is established.

$$u\left(t\right)-P_M\left(t\right)=\frac{{\partial }^{M+1}u\left(\xi \right)}{\partial x^{M+1}(M+1)!}\prod _{i=0}^M \left(t-c_i\right),$$

(24)

where $\xi \in I$.

Proof. [38].

According to the last theorem, we can write

$$\left|u\left(t\right)-P_M\left(t\right)\right|\le \underset{x\in I}{\mathrm{m}\mathrm{a}\mathrm{x}} \left|\frac{{\partial }^{M+1}u\left(t\right)}{\partial x^{M+1}}\right|\frac{\prod _{i=0}^M \left|t-c_i\right|}{(M+1)!}.$$

(25)

Now, Assume that

$$\underset{x\in I}{\mathrm{m}\mathrm{a}\mathrm{x}} \left|\frac{{\partial }^{M+1}u\left(t\right)}{\partial x^{M+1}}\right|\le \eta .$$

(26)

Applying this upper bound to Eq. (25) and considering the approximations for Chebyshev interpolation nodes [39] leads to

$$\left|u\left(t\right)-P_M\left(t\right)\right|\le \frac{\eta }{(M+1)!2^{2M+1}}.$$

(27)

Theorem 2. Suppose that $\tilde{u}$ defined in Eq. (3), be the best approximation of real sufficiently smooth function $u$ by Bernoulli polynomials. Then a real constant $\eta$ exists such that

$${∥u\left(t\right)-\tilde{u}\left(t\right)∥}_2\le \frac{\eta }{(M+1)!2^{2M+1}}.$$

(28)

Proof. According to the definition, $\tilde{u}$ is the best approximation of $u$ provided that

$$\forall v\left(t\right)\in {\mathsf{\Pi }}_N;{∥u\left(t\right)-\tilde{u}\left(t\right)∥}_2\le \parallel u\left(t\right)-v\left(t\right){\parallel }_2.$$

(29)

Particularly, if $v\left(t\right)=P_M\left(t\right)$ then according to Eq. (27), we get

$$\begin{array}{rr}{∥u\left(t\right)-\tilde{u}\left(t\right)∥}_2^2& \le {∥u\left(t\right)-P_M\left(t\right)∥}_2^2={\int }_0^1 {\left|u\left(t\right)-P_M\left(t\right)\right|}^{2}dt\\ & \le {\int }_0^1 {\left[\frac{\eta }{(M+1)!2^{2M+1}}\right]}^{2}dt={\left[\frac{\eta }{(M+1)!2^{2M+1}}\right]}^2.\end{array}$$

(30)

Hence, Eq. (28) is proved.

According to Eq. (28), it can be written

$${∥u\left(t\right)-\tilde{u}\left(t\right)∥}_2=\mathcal{O}\left(\frac{1}{(M+1)!2^{2M+1}}\right).$$

(31)

So, $\frac{1}{(M+1)!2^{2M+1}}\to 0$ when $M\to \mathrm{\infty }$ which implies that $\tilde{u}\to u$. Therefore, the collocation method based on the Bernoulli polynomials is convergent.

• Theorem 3. Assume that $u_{M,j}$ be the approximate solution of Eq. (6) in the interval $I_j=[\left(j-1\right)h,jh]$, and

  $$\left(1-L_1{\lambda }_1\right)\left(1-L_2{\lambda }_2\right)\dots \left(1-L_j{\lambda }_{\mathrm{j}}\right)>0,$$

  where $j=1,\dots .,N$. Also, the nonlinear term satisfies the Lipschitz condition as follows:

  $$‖N\left(u\left(x\right)\right)-N\left(u_{M,j}\left(x\right)\right)‖\le L_j‖u\left(x\right)-u_{M,j}‖.$$

  (32)

Then, there is an upper error bound as follows.

$$∥u\left(t\right)-u_{M,j}\left(t\right)∥\le \frac{{\lambda }_1}{\left(1-{\lambda }_1\right)\left(1-{\lambda }_2\right)\dots \left(1-{\lambda }_{\mathrm{j}}\right)}$$

(33)

where

$$\mathrm{m}\mathrm{a}\mathrm{x}\left|k\left(x,t\right)\right|={\lambda }_{\mathrm{j}}, x\in I_j$$

(34)

Proof. The approximate solution of Eq. (6) in the interval $I_1$ could be stated as

$$u_{M,1}\left(x\right)=f\left(x\right)+{\int }_0^x k\left(x,t\right)N\left(u_{M,1}\left(t\right)\right)dt.$$

(35)

From Eqs. (6) and (35), we get

$u\left(x\right)-u_{M,1}\left(x\right)={\int }_0^x k\left(x,t\right)\left(N\left(u\left(t\right)\right)-N\left(u_{M,1}\left(t\right)\right)\right)dy.$Then, we have

$∥u\left(x\right)-u_{M,1}\left(x\right)∥\le L_1∥k\left(x,t\right)∥∥u\left(t\right)-u_{M,1}\left(t\right)∥$By using Eq. (34), for $i=1$, we have

$∥u\left(x\right)-u_{M,1}\left(x\right)∥\le L_1{\lambda }_1∥u\left(x\right)-u_{M,1}\left(x\right)∥.$So,

$$∥u\left(x\right)-u_{M,1}\left(x\right)∥\le \frac{1}{1-L_1{\lambda }_1}.$$

(36)

Now, regard the approximate solution of Eq. (6) in the interval $I_2$as follows.

$$u_{M,2}\left(x\right)=f\left(x\right)+{\int }_0^h k\left(x,t\right)N\left(u_{M,1}\left(t\right)\right)dt+{\int }_h^x k\left(x,t\right)N\left(u_{M,2}\left(t\right)\right)dt.$$

(37)

From Eqs. (6) and (36), we get

$\begin{array}{cc}\hfill u\left(x\right)-u_{M,2}\left(x\right)& \\ & ={\int }_0^h k\left(x,t\right) \left(N\left(u\left(t\right)\right)-N\left(u_{M,1}\left(t\right)\right)\right)dt\hfill \\ & +{\int }_h^{2h} k\left(x,t\right) \left(N\left(u\left(t\right)\right)-N\left(u_{M,2}\left(t\right)\right)\right)dt.\hfill \end{array}$Then, we have

$∥u\left(x\right)-u_{M,2}\left(x\right)∥\le ∥k\left(x,t\right)∥∥u\left(x\right)-u_{M,1}\left(x\right)∥+∥k\left(x,t\right)∥∥u\left(x\right)-u_{M,2}\left(x\right)∥$Using Eq. (34), for $i=2$ and also Eq. (36) leads to

$∥u\left(x\right)-u_{M,2}\left(x\right)∥\le {L_1\lambda }_1\left(\frac{1}{1-{L_1\lambda }_1}\right)+L_2{\lambda }_2∥u\left(x\right)-u_{M,2}\left(x\right)∥.$Therefore, we have

$$∥u\left(x\right)-u_{M,2}\left(x\right)∥\le \frac{L_1{\lambda }_1}{\left(1-L_1{\lambda }_1\right)\left(1-{L_1\lambda }_2\right)}.$$

(38)

For the approximate solution of Eq. (6) in the interval $I_3$ we can write

$$\begin{array}{cc}\hfill u_{M,3}\left(x\right)=f\left(x\right)& +{\int }_0^h k\left(x,t\right)N\left(u_{M,1}\left(t\right)\right)dt+{\int }_h^{2h} k\left(x,t\right)N\left(u_{M,2}\left(t\right)\right)dt+\hfill \\ & {\int }_{2h}^x k\left(x,t\right)N\left(u_{M,3}\left(t\right)\right)dt.\hfill \end{array}$$

(39)

From Eqs. (6) and (36), we can write

$\begin{array}{cc}\hfill u\left(x\right)-u_{M,3}\left(x\right)& \\ & ={\int }_0^h k\left(x,t\right)\left(N\left(u\left(t\right)\right)-N\left(u_{M,1}\left(t\right)\right)\right)dt\hfill \\ & +{\int }_h^{2h} k\left(x,t\right)\left(N\left(u\left(t\right)\right)-N\left(u_{M,2}\left(t\right)\right)\right)dt\hfill \\ & +{\int }_{2h}^x k\left(x,t\right)\left(N\left(u\left(t\right)\right)-N\left(u_{M,3}\left(t\right)\right)\right)dt.\hfill \end{array}$

Then, we have

$\begin{array}{cc}\hfill ∥u\left(x\right)-u_{M,3}\left(x\right)∥& \\ \hfill \le & L_1∥k\left(x,t\right)∥∥u\left(t\right)-u_{M,1}\left(t\right)∥+L_2∥k\left(x,t\right)∥∥u\left(t\right)-u_{M,2}\left(t\right)∥\hfill \\ \hfill +& L_3∥k\left(x,t\right)∥∥u\left(t\right)-u_{M,3}\left(t\right)∥\hfill \end{array}$Using Eq. (34), for $i=3$ and also Eq. (38) leads to

$\begin{array}{c}∥u\left(x\right)-u_{M,3}\left(x\right)∥\le \frac{L_1{\lambda }_1}{1-L_1{\lambda }_1}+\frac{L_1{\lambda }_2{\lambda }_1}{\left(1-{\lambda }_1\right)\left(1-{\lambda }_2\right)}+{\lambda }_3∥u\left(t\right)-u_{M,3}\left(t\right)∥\\ =\frac{L_1{\lambda }_1}{\left(1-{\lambda }_1\right)\left(1-{\lambda }_2\right)}+{\lambda }_3∥u\left(x\right)-u_{M,3}\left(x\right)∥.\end{array}$Finally, we have

$$∥u\left(x\right)-u_{M,3}\left(x\right)∥\le \frac{L_1{\lambda }_1}{\left(1-{L_1\lambda }_1\right)\left(1-{L_2\lambda }_2\right)\left(1-L_3{\lambda }_3\right)}.$$

(40)

By comparing the upper error bounds obtained in previous steps, it can be concluded that an upper error bound for $u_{M,j}$ is as follows:

$∥u\left(x\right)-u_{M,j}\left(x\right)∥\le \frac{L_1{\lambda }_1}{\left(1-{L_1\lambda }_1\right)\left(1-L_2{\lambda }_2\right)\dots \left(1-L_j{\lambda }_{\mathrm{j}}\right)}.$

Therefore, Eq. (33) is established.

5. Numerical Results

In this section, several examples are presented to demonstrate the validity, applicability, and efficiency of the suggested scheme. Whole the numerical calculations were carried out utilizing Matlab software (R2018b).

• Example 1. Consider the following linear the second kind VIE [40].

  $$u\left(x\right)=\cos x-e^x\sin x+{\int }_0^x{e}^{x}u\left(t\right)dt,$$

  (41)

  with exact solution $u\left(x\right)=\cos x.$.

The numerical results for example one are presented in Table 1, Table 2 and Table 3. In Table 1, the absolute errors of the IBCM with $M=5$ are compared with those of the BCM for two different values of $N$ ($N=3$ and $=5$) at five points in the interval $\left[\mathrm{0,1}\right]$. Table 2 and Table 3 present comparisons similar to what has been presented in Table 1. Although, in these tables, the number of subintervals in the IBCM is increased. In Table 2, $M$ is 10 while in Table 3, it is doubled. Investigating Table 1, Table 2 and Table 3 reveals that in both methods, the higher the polynomial degree ($N$), the higher the accuracy. On the other hand, the IBCM is more accurate than the BCM. For example, for $N=5$, the order of error in the IBCM with $M=20$ (Table 3) is e-12 while the order of error in the BCM is e-3 (Table 1). In Table 3, for $N=7$, the order of error in the IBCM is e-16 while the order of error in the BCM is e-8.

The precision of the IBCM can be ameliorated by adding the number of subintervals, $M$ while $N$ is fixed. For example, for $N=5$, the order of error in the IBCM with $M=5$ (Table 1) is e-9 while it is e-10 and e-12 for $M=10$ (Table 2) and $M=20$ (Table 3), respectively. Thee error of IBCM is plotted in Figure 1 for $M=20$ and $N=7$.

This equation was solved in [12] by Bernstein’s approximation. The authors computed errors for $n=2,\dots ,9$ where $n$ is the polynomial degree. The order of error was e-10 for $n=9$ while in our propose method, the order of error is e-16 for $n=7$.

Example 2. Regard the following nonlinear the second kind VIE [23].

$$u\left(x\right)=e^x-\frac{1}{3}\left(e^{3x}+1\right)+{\int }_0^x{u\left(t\right)}^{3}dt,$$

(42)

with exact solution $u\left(x\right)=e^x.$

The results for the above example are presented in Table 4 and Table 5. In Table 4, the absolute errors of the IBCM with $M=10$ are compared with those of the BCM for two different values of $N$ ($N=4$ and $=6$ ). Investigating this table reveals that in both methods, the higher the polynomial degree ($N$), the higher the accuracy. On the other hand, the IBCM is more precise than the BCM. In Table 5, numerical results for both IBCM and BCM are presented for $N=8$ and they are compared with the results reported in [23]. According to this table the accuracy of the suggested method is more than the method of [23]. Thee error of IBCM is plotted in Figure 2 for $M=5$ and $N=8$.

Example 3. Regard the following nonlinear the second kind VIE [23].

$$u\left(x\right)=-\frac{1}{10}{x}^4+\frac{5}{6}{x}^2+\frac{3}{8}+{\int }_0^x\frac{1}{2x}{u\left(t\right)}^{2}dt,$$

(43)

with exact solution $u\left(x\right)=x^2+\frac{1}{2}.$The numerical results for above example are presented in Table 6. In this table, numerical results for both IBCM and BCM are presented for $N=11$ and they are compared with the results reported in [23]. According to this table the accuracy of the suggested method is more than the method of [23]. Thee error of IBCM is plotted in Figure 3 for $M=5$ and $N=11$.

Example 4. Regard the following nonlinear the second kind VIE [41].

$$u\left(x\right)=\sin \left(\pi x\right)+{\int }_0^x\sin \left(\pi t\right)\cos \left(\pi x\right){u\left(t\right)}^{3}dt,$$

(43)

with exact solution $u\left(x\right)=\sin \left(\pi x\right)+\frac{20-\sqrt{391}}{3}\cos \left(\pi x\right).$The numerical results for this example are presented in Table 7 and Table 8. In Table 7, exact and approximate solutions of Eq. (43) by IBCM and Modification of hat functions [41] are presented. Comparison the results shows that IBCM is more accurate. In Table 8, the absolute errors of the IBCM with $M=10$ are compared with those of the BCM for $N=5$. Investigating this table reveals that in both methods, the higher the polynomial degree ($N$), the higher the accuracy. On the other hand, the IBCM is more accurate than the BCM. Thee error of IBCM is plotted in Figure 1 for $M=10$ and $N=10$.

Figure 3. The absolute error of the IBCM with $M=10$ and $N=10$ for example 3.

Figure 3. The absolute error of the IBCM with $M=10$ and $N=10$ for example 3.

6. Conclusions

In this work, an improved collocation method based on the Bernoulli polynomials was presented to solve the VIE of the second kind. In classic collocation methods, regardless of the type of polynomial used, collocation points are scattered throughout the whole interval and numerical computations are performed at once on the given interval. The main goal of the suggested method is to improve the results of the classic BCM by dividing the interval into some sub-intervals and considering the collocation points on each of them. Here, the zeros of the SCPs are considered as collocation points. Then, BCM is applied step by step from the first sub-interval to the last one. By this process, a system of algebraic equations is attained for each sub-interval which can be lightly solved using computing software. At last, the approximate solution is obtained as a piece-wise function. Convergence of the scheme was also analyzed. Several numerical examples are presented in order to illustrate the validity, applicability, and efficiency of the suggested method. Numerical results show that in both methods, the higher the accuracy. On the other hand, IBCM is more accurate than the BCM. The precision of the IBCM can be ameliorated by adding the number of subintervals while the degree of the polynomial is fixed. We suggest that to test idea on other kinds of collocation methods.