Legendre-Gauss-Lobatto collocation method for solving multi-dimensional systems of mixed Volterra-Fredholm integral equations

1.   Introduction

Integral equations ^[1,2,3,4,5] can accurately describe many different phenomena in engineering and science. Although there are several powerful techniques to numerically solve integral equations, most reveal limitations with multi-dimensional problems. The systems of mixed Volterra-Fredholm integral equations (SMVF-IEs) ^{[6,7,8,9,10,11,12]} appear in the scope of parabolic boundary value models, including in physics, mathematics, biology, and other subjects. The solution of SMVF-IEs has been a matter of substantial interest. However, solving SMVF-IEs is a challenging issue, for which effective methods are still lacking. Population dynamics, parabolic boundary value problems, spatio-temporal evolution of epidemics, and other phenomena lead to SMVF-IE-based models.

In reference [13], a variational iteration technique was proposed for solving systems of integro-differential equations. In paper ^[14], a homotopy perturbation technique was utilized to numerically solve nonlinear SMVF-IEs, while in ^[15], it was adopted to solve nonlinear two-dimensional SMVF-IEs. In ^[16], hybrid mixtures were used to solve the three-dimensional $L$-shaped channel problem and to analyze the properties of heat generation via forced convection. In ^[17], the Keller Box approach was utilized to simulate nonlinear problems encountered in developing liquid and supplementary algebraic dynamics domains. In ^[18,19], numerical and analytical methods were developed to address SMVF-IEs. In ^[20], the Galerkin finite element method was used to find a closed-form solution for a nonlinear coupled partial differential equation, whereas the author in ^[21] used a quasi-digital technology called the differential transformation method to address the control of complex PDE devices. In ^[22], the Caputo–Fabrizio and Atangana–Baleanu techniques were used to solve fractional dimensionless systems, and performed a theoretical analysis via the Chebyshev spectral approach ^[23].

Spectral methods are effective for solving several types of differential and integral problems ^{[24,25,26,27]}. Particular kinds of spectral approaches include the Galerkin ^[28,29], collocation ^{[30,31,32,33,34,35,36,37]} and tau techniques ^{[39,40,41,42]}. Contrasting with other approaches, namely the finite difference and finite element methods, spectral techniques achieve superior accuracy, even with few nodes, thus involving a smaller computational burden. Indeed, they are characterized by exponential convergence rates. The main idea is to express the solution to the original equation by a finite sum of a certain basis function, and then to choose the functions' coefficients such that the error between the exact and the numerical solutions is minimized. In the spectral collocation variant ^{[43,44,45,46,47,48,49,50]}, the numerical solution is compelled to closely satisfy the original problem, thus the residuals may approach zero at specific collocation points.

This paper addresses the solution of SMVF-IEs by means of a Legendre-Gauss-Lobatto collocation method. The one-dimensional case is firstly treated. Afterwards, the method is extended to two-dimensional linear and nonlinear SMVF-IEs. Several numerical examples are presented to demonstrate the effectiveness of the approach, and to illustrate its superiority in comparison with alternative techniques for solving SMVF-IEs.

The paper is organized into sections. In section 2, some mathematical preliminaries are outlined. In section 3, one-dimensional SMVF-IEs are solved. In section 4, the novel algorithm is extended for solving two-dimensional SMVF-IEs. In section 5, error analysis of the method is addressed. In section 6, some numerical examples assess and compare the proposed approach with other techniques. In section 7, the main conclusions are summarized.

2.   Mathematical preliminaries

The Legendre polynomials $\mathcal{L}_{\jmath}(\varrho)$ $(\jmath = 0, 1\ldots)$ comply with the Rodrigues' expression ^[51]:

where $\jmath$ stands for degree.

Accordingly, the ${n_1}$th derivative of $\mathcal{L}_\jmath(\varrho)$, denoted by $\mathcal{L}^{{n_1}}_\jmath (\varrho)$, is given by the following:

where

Let us denote the norm and inner product by $\left\Vert{\Lambda}\right\Vert$ and $(\Lambda, \gamma)$, respectively, of space $L^2[-1, 1]$. The collection of $\mathcal{L}(\varrho)$ is a whole orthogonal system in $L^2[-1, 1]$ ^[52]

where $h_{i} = \frac{2}{2i+1}$, and $\delta_{{j_1}\jmath}$ is the Dirac function. Hence, for each $\gamma\in L^2[-1, 1]$,

Assume that $S_{N_1}[-1, 1]$ is a collection of all polynomials of degree at the utmost $N_1$ ($N_1\geq0$). Hence, for each $\varphi\in S_{2{N_1}-1}[-1, 1]$, we have

where $\varrho_{{N_1}, \jmath}$ ($0\leq \jmath\leq {N_1}$) and $\varpi_{{N_1}, \jmath}$ ($0\leq \jmath\leq N_1$) are the nodes and Christoffel numbers of the Legendre-Gauss-Lobatto interpolation on the classical interval $[-1, 1]$, respectively. The discrete norm and inner product correspond to

Denote the shifted Legendre polynomials specified on the interval $[0, L]$ by $\mathcal{L}_{L, \jmath}{(\varrho)}$. These polynomials are obtained by the recurrence ^[51]

Then, we can write

The integral $I^{\nu}{\mathcal{L}}_{L, {j_1}}(\varrho)$ may be obtained from

where $\mathcal{L}_{L, {j_1}}{(0)} = {(-1)}^{j_1}$. The equation of the orthogonality condition is

where $w_L(\varrho) = 1\ {\text{and}}\ h_\jmath^{L} = \dfrac{L}{2\jmath+1}.$

If function $\Lambda(\sigma)\in L^2[0, L]$, then it can be expressed by $\mathcal{L}_{L, i}(\sigma)$ as

with $c_i$ given by

In the approximation, $\Lambda(\varrho)$ can be written as

3.   One-dimensional SMVF-IEs

Let us consider the one-dimensional SMVF-IEs

where $\Delta_{1}(\varrho)$, $\Delta_{2}(\varrho)$, $\jmath_{1}(\varrho, \sigma)$, $\jmath_{2}(\varrho, \sigma)$, $\jmath_{3}(\varrho, \sigma)$ and $\jmath_{4}(\varrho, \sigma)$ are given as real valued functions, while $\Lambda(\varrho)$ and $\gamma(\varrho)$ are unknown functions.

Using $\sigma = \frac{\varrho}{L}\eta$ to write the integrals $\int\limits_{0}^{\varrho}{\jmath_{1}(\varrho, \sigma)\Lambda(\sigma)}d\sigma$, $\int\limits_{0}^{\varrho}{\jmath_{3}(\varrho, \sigma)\Lambda(\sigma)}d\sigma$ in the interval $[0, L]$, for the variable $\eta$, we perform the shifted Legendre-Gauss-Lobatto integration

We rely on the shifted Legendre-Gauss-Lobatto collocation technique to turn the previous SMVF-IEs into a system of algebraic equations. Thus, we collocate independent variables at $\varrho_{L, {N_1}, {j_1}}$ points, yielding the approximate solution

Assume $S_{N_1}(0, L)$ is the collection of all polynomials of degree at utmost $N_1$ for any positive integer $N_1$. It follows for each $\phi\in S_{2{N_1}-1}(0, L)$, based on the shifted Legendre-Gauss-Lobatto quadrature,

where $\varpi_{L, {N_1}, i}$ are known on the interval [0, L], meaning the Christoffel numbers of the shifted Legendre-Gauss-Lobatto interpolation.

Based on Eqs (3.3) and (3.4), one can write Eq (3.2) as

In the shifted Legendre-Gauss-Lobatto collocation technique presented herein, the residual of (3.5) is made zero at the ${N_1}+1$ shifted Legendre-Gauss points

After the coefficients $a_{j}, \, b_{j}$ are specified, the approximate solution $\Lambda_{N_1}(\varrho), \, \gamma_{N_1}(\varrho)$ at any value of $\varrho\in[0, L]$ in the specific domain can be easily computed from the equations

4.   Two-dimensional SMVF-IEs

The preceding numerical algorithm is extended to solve linear and nonlinear two-dimensional SMVF-IEs. The collocation points are chosen at the shifted Legendre-Gauss-Lobatto interpolation nodes. The idea is to discretize the SMVF-IEs and to construct a system of algebraic equations.

4.1.   Linear SMVF-IEs

Let us consider the two dimensional SMVF-IEs

where $\Lambda(\varrho, \sigma)$ and $\gamma(\varrho, \sigma)$ are unknown functions, whilst $\Delta_{1}(\varrho, \sigma), \, \Delta_{2}(\varrho, \sigma), \, \jmath_{1}(\varrho, \sigma)$ and $\jmath_{2}(\varrho, \sigma, \chi, \psi),$ are given as real valued functions.

Using the change of variable $\psi = \frac{\sigma}{\tau}\eta$, we can transform the integrals $\int\limits_{0}^{\sigma}\int\limits_{0}^{L}{\jmath_{1}(\varrho, \sigma, \chi, \psi) \Lambda(\chi, \psi) d\chi d\psi}$, $\int\limits_{0}^{\sigma}\int\limits_{0}^{L}{\jmath_{2}(\varrho, \sigma, \chi, \psi) \gamma(\chi, \psi) d\chi d\psi, }$ into the interval, $[0, \tau]$, for the variable $\eta$, to immediately execute the shifted Legendre-Gauss-Lobatto integration,

We extend the dependent variable by the model

In virtue of the Eqs (4.3) and (3.4), we can rewrite Eq (4.2) as

where

The residual of (4.4) is set to be zero in the suggested shifted Legendre-Gauss-Lobatto collocation technique at $({N_1}+1)\times({N_2}+1)$ of shifted Legendre-Gauss-Lobatto points

where ${n_1} = 0, \ldots, {N_1}$ and ${n_2} = 0, \ldots, {N_2}.$

Finally, Eq (4.5) is enforced to exactly satisfy (4.1) at the shifted Legendre-Gauss-Lobatto interpolation nodes $\varrho_{L, {N_1}, {n_1}}, \, \sigma_{L, {N_1}, {n_2}}$. This provides $2({N_1}+1)({N_2}+1)$ equations for $a_{i, {j_1}}, \, b_{i{j_1}}; \, i = 0, \cdots, {N_1}, \ {j_1} = 0, \cdots, {N_2}$. Consequently, the approximate solution (4.3) can be evaluated as

4.2.   Nonlinear SMVF-IEs

We expand the technique for numerically handling nonlinear SMVF-IEs

where $\Lambda(\varrho, \sigma)$ and $\gamma(\varrho, \sigma)$ are unknown functions, whilst $f(\varrho, \sigma)$, $k(\varrho, \sigma, \chi, \psi, (\Lambda(\chi, \psi))), \jmath_1(\varrho, \sigma, \chi, \psi, \Lambda(\chi, \psi), (\gamma(\chi, \psi)))$ and $\jmath_2(\varrho, \sigma, \chi, \psi, \Lambda(\chi, \psi), \gamma(\chi, \psi))$ are given as functions.

Using the change of variable $\psi = \frac{\sigma}{\tau}\eta$, we can transform the integrals $\int\limits_{0}^{\sigma}\int\limits_{0}^{L}{\jmath_{1}(\varrho, \sigma, \chi, \psi) \Lambda(\chi, \psi) d\chi d\psi}$, $\int\limits_{0}^{\sigma}\int\limits_{0}^{L}{\jmath_{2}(\varrho, \sigma, \chi, \psi) \gamma(\chi, \psi) d\chi d\psi, }$ into the interval, $[0, \tau]$, for the variable $\eta$, and apply the shifted Legendre-Gauss-Lobatto integration

We select the approximate solution from the model

Proceeding as in the previous subsection, we can rewrite the problem in the form

where

The residual of (4.10) is set to be zero at $({N_1}+1)\times({N_2}+1)$ for the shifted Legendre-Gauss-Lobatto points

where

where ${n_1} = 0, \ldots, {N_1}$ and ${n_2} = 0, \ldots, {N_2}.$

By utilizing Newton's iterative technique, we can solve the preceding nonlinear system of algebraic equations. As a result, the coefficients $a_{ij}, \, b_{ij}$ are specified, and the approximate solution can be accounted for $\Lambda_{{N_1}, {N_2}}(\varrho, \sigma), \, v_{{N_1}, {N_2}}(\varrho, \sigma)$ at any value of $(\varrho, \sigma)$ in the known domain, by means of the next equation

5.   Error analysis

In this section, we investigate an error analysis of the present method.

Definition 5.1. For a nonnegative integer $\rho$, we have ^[51,53]

whereas $\partial^{\imath}_{z}\Lambda(z) = \frac{\partial^{\imath}\Lambda(z)}{\partial{z}^{\imath}}$, and

Lemma 5.2. For $\Lambda\in {B^{q}}(\omega^{d})$ with $d\leq q\leq N+1$ ^[54]

Theorem 5.3. Let ${\boldsymbol{I}}_{N}{\boldsymbol{\Lambda(\varrho)}}$ and ${\boldsymbol{\Lambda(\varrho)}}$ be the spectral approximation and the exact solution of the Volterra-Fredholm system. Thus, we have

Proof. We can write the system in Eq (3.1) as a multivariate system

When utilizing the approximate solution we have

Subtracting (5.4) from (5.3) yields

where

We can write (5.5) by using Gronwall inequality

Then, the term $\parallel B_{1}\parallel$ is estimated as the following:

By using the Cauchy inequality, we can get,

Hence,

Now we estimate the term $\|B_{2}\|$. We use the Legendre-Gauss integration formula (2.3) to obtain

We obtain it by using the Cauchy-Schwarz inequality

By using the Lipschitz condition, we can write

By using the triangle inequality, we derive that

Furthermore, from Lemma 5.2, we can deduce that

Then, the term $\parallel B_{3}\parallel$ is estimated as the following:

By using the Cauchy inequality, we can get,

Hence,

Now, we estimate $\|B_{4}\|$. We use Eq (2.3) to obtain

We obtain it by using the Cauchy-Schwarz inequality

By using the Lipschitz condition, we can write

By using the triangle inequality, we have that

Additionally, from Lemma 5.2, we can deduce that

6.   Numerical results

To illustrate the performance of the proposed scheme and the thoroughness of the results, we present some numerical examples. The results with the new method are compared with those yielded by others ^{[11,15,18,19]}. For assessing accuracy, the difference between the exact, $\Lambda(\varrho)$, and the approximate, $\Lambda_{N, M}(\varrho)$, solutions at the point $\varrho$, meaning the absolute error $E$ is adopted, that is,

Additionally, the maximum absolute error ($M_E$) is specified by

Example 1. Let us consider the one-dimensional linear SMVF-IEs ^[18]:

where $\Delta_{1}(\varrho)$ and $\Delta_{2}(\varrho)$ are given functions, founded by the exact solution $\Lambda_{1}(\varrho) = 1-\varrho$, $\gamma_{2}(\varrho) = \varrho$.

In order to illustrate the convergence rate of our technique, we list the $M_E$ for several options of $N_1$ in Table 1, while comparing our method with the fixed point collocation approach (FPCA) ^[18]. We verify that the new technique leads to a better numerical solution with far fewer nodes, and that our numerical solutions are very close to the exact ones. In Figure 1, we illustrate the $M_E$ (i.e., $\log_{10} M_E$) obtained with the new technique for diverse values of $N$, {where $E_{1}$ and $E_{2}$ correspond to $\Lambda_{1}(\varrho)$ and $\gamma_{2}(\varrho)$, respectively, in logarithmic graphs being computed as in Eq (6.1)}.

Example 2. We solve the two-dimensional linear SMVF-IEs ^[19]:

that has the exact solution $\Lambda_{1}(\varrho, \sigma) = \sigma\sin(\varrho)$, $\gamma_{2}(\varrho, \sigma) = \sigma\cos(\varrho)$.

For several choices of $N$, we verify that the new technique is more accurate than the homotopy method (HAM) ^[19]. Table 2 summarizes the values of the absolute errors $E$ of problem 2.

Figures 2 and 3, illustrate the evolution of the absolute errors along each dimension, $E_1$ and $E_2$, for $N = 12$. To emphasize the high thoroughness and convergence rate of the new scheme, we depict the maximum absolute errors in log scale in Figure 4. Based on the results, we may conclude that the proposed technique yields excellent approximations and exponential convergence rates.

Example 3. We now consider the nonlinear two-dimensional SMVF-IEs ^[15]:

where $(\varrho, \sigma)\in[0, 1]\times[0, 1]$, and with the exact solution $\Lambda_{1}(\varrho, \sigma) = \varrho+\sigma$, $\gamma_{2}(\varrho, \sigma) = \varrho^{2}-\sigma^{2}$.

Table 3 shows the absolute errors by utilizing the new algorithm and those with the homotopy perturbation method (HPM) ^[15]. It is observed that the proposed technique is more precise than the HPM.

Example 4. We consider the nonlinear two-dimensional SMVF-IEs ^[19]:

where $\Delta_{1}(\varrho, \sigma)$ and $\Delta_{2}(\varrho, \sigma)$ are the given real valued functions, and the exact solution is $\Lambda_{1}(\varrho, \sigma) = -2\varrho+2\varrho\sigma$, $\gamma_{2}(\varrho, \sigma) = 1+2\varrho\sin(\sigma)$.

Table 4 compares the absolute error $E$ resulting from the application of the new proposed method with that by the approach in ^[19] for several values of $N\ {\text{and}}\ M$. The numerical results show that the solutions are extremely precise, even with small values of $N\ {\text{and}}\ M$.

7.   Conclusions

In this paper, a shifted Legendre-Gauss-Lobatto collocation scheme was proposed for numerically solving SMVF-IEs. First, the one-dimensional case was solved, and then the method was extended to address two-dimensional linear and nonlinear SMVF-IEs. Different numerical examples revealed the superiority of the proposed approach when compared with alternative techniques for treating SMVF-IEs. The SMVF-IEs play a crucial role in many scientific and engineering problems and, thus, methods for solving then accurately are crucial.

Use of AI tools declaration

The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

Acknowledgements

The authors would like to thank the Deanship of Scientific Research at Umm Al-Qura University for supporting this work by Grant Code: (22UQU4380165DSR01).

Conflict of interest

The authors declare no conflicts of interest.