EUROPEAN JOURNAL OF PURE AND APPLIED MATHEMATICS
2025, Vol. 18, Issue 3, Article Number 6229
ISSN 1307-5543 – ejpam.com
Published by New York Business Global

New Multi-Step Collocation Methods for Non-Standard
Volterra Integral Equations

Hayder Abdullah Thabbah1, Saeed Pishbin1, Parviz Darania1,∗,
1 Department of Mathematics, Urmia University, Urmia, Iran

Abstract. Non-standard Volterra integral equations are in the category of non-linear Volterra in-
tegral equations which appear in certain mathematical modelling processes. In the current study,
new multi-step collocation techniques have been employed to numerically address non-standard
Volterra integral equations. Several established outcomes concerning convergence of the numer-
ical solution for this problem are investigated. In conclusion, two test problems are thoroughly
examined to substantiate theoretical accomplishments in practical terms. Also, we apply one-step
and multi-step collocation methods to solve Non-standard Volterra integral equations and compare
numerical results to show the efficiency and high accuracy of the proposed method.

2020 Mathematics Subject Classifications: 65R20, 45G10

Key Words and Phrases: Non-standard Volterra integral equations (NVIEs), Volterra integral
equations (VIEs), Error analysis, Multi-step collocation methods.

1. Introduction

In general, the kernels of integral equations are in the form K(t, s, x(s)). But in some
scientific and engineering applications, especially in physical sciences, a special case of
nonlinear integral equations appears, where the kernels of the equation is dependent on
time t as well as space s. This dependence of the kernels is shown in the general form as
K(t, s, x(t), x(s)). Therefore, integral equations including these types of kernels are called
non-standard integral equations.

A non-linear VIE as

y(σ) = f(σ) +

∫ σ

0
K(σ, τ, y(σ), y(τ))dτ, σ ∈ [0, T ] = I, (1)

is called a NVIEs. Different kind of physical and biological phenomena that lead to NVIEs
(1) can be found in [1–5]. Some researchers [6–12] have investigated a number of these
non-standard integral equations such as auto convolution equations and implicit VIEs.

∗Corresponding author.
DOI: https://doi.org/10.29020/nybg.ejpam.v18i3.6229

Email addresses: haider20206021@yahoo.com (H. A. Thabbah),
s.pishbin@urmia.ac.ir (S. Pishbin), p.darania@urmia.ac.ir (P. Darania)

https://www.ejpam.com 1 Copyright: © 2025 The Author(s). (CC BY-NC 4.0)



H. A. Thabbah, S. Pishbin, P. Darania / Eur. J. Pure Appl. Math, 18 (3) (2025), 6229 2 of 16

Certain established outcomes rooted in regularity, existence, uniqueness, as well as the
convergence of the solution for the general auto-convolution Volterra integral equations
studied by Zhang et al [12]. Guan et al. in [13], studied NVIEs (1) and obtained several
positive outcomes concerning regularity, existence, uniqueness, and convergence of the
one-step collocation method for this equation. In [14], the authors by adding conditions
for interpolation in the previous r step points, applied multi-step collocation scheme to
solve the equation (1) and increased the order of convergence for the m points and r steps
to m+ r compared with the one-step collocation method (Order of convergence is m [13]).

Actually, nonlinear VIEs originate in various scientific realms and this is the primary
reason why they were a subject of refined investigation from the both theoretical and
numerical perspectives so far. Specifically, they span a wide range of applied science such
as engineering, mathematical physics, economics, etc[15–19]. One of the particular kind
of non-linear VIEs is Volterra-Hammerstein IEs which have been studied in [20]. Usually,
finding an analytical solution for the non-linear VIEs is not possible, then researchers
try to provide numerical methods to avoid this complicated situation( See [16, 21, 22]).
Moreover, the great contributions of recent developments in the conceptually promising
and challenging areas of the nonlinear VIEs are all gathered in the book due to Brunner
[21].

Here, we consider new multi-step collocation scheme to NVIEs (1) and try to increase
the order of convergence compared with the one-step and multi-step collocation methods
in [13, 14].

The paper is structured in the following manner: Section 2 is dedicated to introduc-
ing a new multi-step polynomial collocation method for discretizing NVIEs using widely
recognized interpolation polynomials. In the sequel, we analyze the convergence of the
numerical solutions, considering Peano theorems for interpolation, in section 3. Section 4
then delves into the examination of two test problems to validate the theoretical results.

2. Numerical methods

The multi-step method in the context of integral equations serves a similar purpose
as in ordinary differential equations improving accuracy and computational efficiency but
is tailored to handle integral operators. Its primary aim is to solve integral equations
numerically by leveraging past computed values to approximate the solution at subsequent
steps more effectively.

The fundamental purposes of multi-step techniques in integral equations include:

1. Higher-Order Approximation;
By using information from multiple previous steps (e.g., past values of the unknown
function and its integrals), multi-step methods achieve better accuracy than single-
step quadrature-based approaches.

2. Reduced Computational Cost;
Reusing previously computed integral evaluations (e.g., in Volterra or Fredholm



H. A. Thabbah, S. Pishbin, P. Darania / Eur. J. Pure Appl. Math, 18 (3) (2025), 6229 3 of 16

equations) avoids redundant calculations, making the method more efficient than
naive quadrature rules.

3. Handling Non-Local Dependence;
Integral equations involve global dependence (due to the integral term), and multi-
step methods help approximate the history integral more efficiently by exploiting
smoothness and past data.

4. Improved Stability for Volterra Equations;
Multi-step discretizations can stabilize solutions for weakly singular or nonlinear
integral equations.

We will now introduce and examine two numerical methods: the standard multi-step
collocation method and an advanced version referred to as the new multi-step collocation
method.

2.1. Multi-step collocation scheme for NVIE

The examination conducted in [14], we consider the multi-step collocation scheme for
the equation (1). Consider {tn = nh, n = 0, 1, . . . , N (tN = T )}, as an equidistant grid
on I. Suppose that Sh be given by

Sh := {tn,i = tn + cih : 0 < c1 < · · · < cm ≤ 1, 0 ≤ n ≤ N − 1},

where {tn,i} are collocation points with {ci} as collocation parameters. Take the colloca-
tion solution on the interval ϵn = (tn, tn+1] for s ∈ [0, 1] by

y(tn + sh) =
r−1∑
k=0

Lk(s)xn−k +
m∑
j=1

L̂j(s)Yn,j , n = r, · · · , N − 1, (2)

with
Yn,j := y(tn,j), xn−k := y(tn−k),

Lk(s) =

m∏
i=1

s− ci
−k − ci

·
r−1∏

i=0,i ̸=k

s+ i

−k + i
, L̂j(s) =

r−1∏
i=0

s+ i

cj + i
·

m∏
i=1,i ̸=j

s− ci
cj − ci

.

Note that the starting values x1, . . . , xr can be obtained by a suitable one-step method.
In the sequel, we find y so that the collocation equation is hold for t ∈ Sh,

y(t) = f(t) +

∫ t

0
K
(
t, s, y(t), y(s)

)
ds.

Hence, on every subinterval ϵn, we have the system of algebraic equations for the



H. A. Thabbah, S. Pishbin, P. Darania / Eur. J. Pure Appl. Math, 18 (3) (2025), 6229 4 of 16

unknowns Yn,j , j = 1, . . . ,m, n = r, . . . , N − 1,

Yn,j = f(tn,j) + h
r−1∑
k=0

∫ 1

0
K
(
tn,j , tk + sh, Yn,j , y(tk + sh)

)
ds

+h
n−1∑
k=r

∫ 1

0
K
(
tn,j , tk + sh, Yn,j ,

r−1∑
j=0

Lj(s)xk−j +

m∑
i=1

L̂i(s)Yk,i

)
ds

+h

∫ cj

0
K
(
tn,j , tn + sh, Yn,j ,

r−1∑
k=0

Lk(s)xn−k +
m∑
i=1

L̂i(s)Yn,i

)
ds.

(3)

Also, the convergence theorem results can be presented as:

Theorem 1. [14] Consider
i) f ∈ Cm+r(I) and for D := {(t, s) : 0 ≤ s, t ≤ T}, K ∈ Cm+r(D × R× R),

ii) |K(t, s, r1, s1)−K(t, s, r2, s2)| ≤ Q(t, s)|r1−r2|+M |s1−s2| for all (t, s) ∈ D, ri, si ∈ R

with 0 < Q(s, s) ≤ M and

∫ t

0
Q(t, s)ds ≤ 1−∆ for ∆ ∈ (0, 1),

iii) The order of primary error is m + r, and the spectral radius of the matrix A is less
than 1, (ρ(A) < 1).

Then the NIVE (1) has a unique solution x ∈ Cm+r(I) and for the multi-step colloca-
tion method, we have:

∥x− y∥∞ ≤ Chm+r,

such that C is a constant not depending on N .

2.2. New multi-step collocation methods

To construct the new multi-step methods, we seek a collocation polynomial of the form

y(tn + sh) =
r−1∑
k=0

Lk(s)xn−k +

m∑
j=1

L̂j(s)Yn,j +

m∑
j=1

L̃j(s)Yn+1,j , n = r − 1, . . . , N − 2, (4)

where
Yn,j = y(tn,j), Yn+1,j = y(tn+1,j),

Lk(s) =
r−1∏
i=0
i ̸=k

s+ i

−k + i
·

m∏
i=1

s− ci
−k − ci

·
m∏
i=1

s− ci − 1

−k − ci − 1
, (5)

L̂j(s) =
r−1∏
i=0

s+ i

cj + i
·

m∏
i=1
i ̸=j

s− ci
cj − ci

·
m∏
i=1

s− ci − 1

cj − ci − 1
, (6)



H. A. Thabbah, S. Pishbin, P. Darania / Eur. J. Pure Appl. Math, 18 (3) (2025), 6229 5 of 16

L̃j(s) =
r−1∏
i=0

s+ i

1 + cj + i
·

m∏
i=1

s− ci
1 + cj − ci

·
m∏
i=1
i ̸=j

s− ci − 1

cj − ci
, (7)

also the collocation parameters are assumed to satisfy c1 ̸= 0 and ci ̸= cj for i ̸= j. In the
equations (1) by substituting t = tn+1,i and approximate solution (4), we will have

Yn+1,j = f(tn+1,j) + h
r−2∑
z=0

∫ 1

0
K
(
tn+1,j , tz + sh, Yn+1,j , y(tz + sh)

)
ds

+

n−1∑
z=r−1

h

∫ 1

0

(
K
(
tn+1,j , tz + sh, Yn+1,j ,

r−1∑
k=0

Lk(s)xz−k +

m∑
i=1

L̂i(s)Yz,i

+
m∑
i=1

L̃i(s)Yz+1,i

))
ds+ h

∫ 1

0

(
K
(
tn+1,j , tn + sh, Yn+1,j ,

r−1∑
k=0

Lk(s)xn−k

+
m∑
i=1

L̂i(s)Yn,i +
m∑
i=1

L̃i(s)Yn+1,i

))
ds+ h

∫ cj

0

(
K
(
tn+1,j , tn+1 + sh,

Yn+1,j ,
r−1∑
k=0

Lk(1 + s)xn−k +
m∑
i=1

L̂i(1 + s)Yn,i +
m∑
i=1

L̃i(1 + s)Yn+1,i

))
ds,

(8)
with

yn+1 =
r−1∑
k=0

Lk(1)xn−k +

m∑
j=1

L̂j(1)Yn,j +

m∑
j=1

L̃j(1)Yn+1,j .

In general, the integrals in the system (8), can be computed approximately by using
quadrature formulas to yield the discretized version of this non-linear system.

3. Convergence analysis

Consider the NVIE (1), and it’s corresponding collocation equation

y(t) = f(t) +

∫ t

0
K(t, s, y(t), y(s))ds. (9)

By subtracting equation (9) from (1), we will get

e(t) =

∫ t

0

(
K(t, s, x(t), x(s))−K(t, s, y(t), y(s))

)
ds, (10)

where e(t) = x(t)− y(t), is the error function.
Before investigating the convergence analysis of the new multi-step collocation solu-

tions, Taylor’s theorem to multivariate functions [16, 23, 24] is considered as:



H. A. Thabbah, S. Pishbin, P. Darania / Eur. J. Pure Appl. Math, 18 (3) (2025), 6229 6 of 16

Lemma 1. Assume that the function g : Rn → R, for some r > 0 in a closed ball
A = {y ∈ Rn : ||x0 − y|| ≤ r}, g ∈ Ck+1(A). Then we can derive the following formula
for the remainder in A, Namely,

g(x) =
∑
|α|≤k

Dαg(x0)

α!
(x− x0)

k +
∑

|β|=k+1

Rβ(x)(x− x0)
β,

where

Rβ(x) =
|β|
β!

∫ 1

0
(1− t)|β|−1Dβg(x0 + t(x− x0))dt.

Theorem 2. Let the second condition of Theorem 1 is hold and for the known functions
in (1):

A1. K ∈ Cp(D × R× R), for some m ≥ 1, and f ∈ Cp(I).

A2. Order of primary error is p and the spectral radius of the matrix à is less than 1
(ρ(Ã) < 1) where

à =

[
000r−1,1 Ir−1

Lr−1(1) Lr−2(1), . . . , L0(1),

]
. (11)

Then, there exist constant C not depending on N such that

∥x− y∥∞ ≤ Chp,

where p = 2m+ r.

Proof. Now, by Peano’s theorem [16, 25], for s ∈ [0, 1], the error e(t) := x(t)− y(t) is
written as:

e(tλ + sh) =
r−1∑
k=0

Lk(s)ελ−k +
m∑
j=1

L̂j(s)ελ,j +
m∑
j=1

L̃j(s)ελ+1,j + h2m+rRm,r,λ(s), (12)

where Rm,r,n(s) is the Peano remainder and ελ,j := Xλ,j − Yλ,j . From (8) and (10), for
λ = n+ 1, we have



H. A. Thabbah, S. Pishbin, P. Darania / Eur. J. Pure Appl. Math, 18 (3) (2025), 6229 7 of 16

εn+1,j =

∫ tn+1,j

0
K(tn+1,j , s, xn+1,j , x(s))ds−

∫ tn+1,j

0
K(tn+1,j , s, yn+1,j , y(s))ds

=

∫ tn+1,j

0

(
K(tn+1,j , s,Xn+1,j , x(s))−K(tn+1,j , s, Yn+1,j , x(s))

)
ds

+h
n∑

l=0

∫ 1

0

(
K(tn+1,j , tl + sh, Yn+1,j , x(tl + sh))

−K(tn+1,j , tl + sh, Yn+1,j , y(tl + sh))

)
ds

+h

∫ cj

0

(
K(tn+1,j , tn+1 + sh, Yn+1,j , x(tn+1 + sh))

−K(tn+1,j , tn+1 + sh, Yn+1,j , y(tn+1 + sh))

)
ds.

(13)

Now, assume that

K̂l(tη,j , s) :=

∫ 1

0

∂K

∂v
(tη,j , tl + sh, Uη,j , u(tl + sh) + ζ(y(tl + sh)− u(tl + sh)))dζ, (14)

where for l = 0, 1, ..., n − 1 the function K̂l(tη,j , s) ∈ C[0, 1]. Likewise, for l = n,
K̂n(tη,j , s) ∈ C[0, cj ].
Also,

K̂(tη,j , s) :=

∫ 1

0

∂K

∂u
(tη,j , s, Uη,j + ζ(Yη,j − Uη,j), y(s))dζ, (15)

then K̂(tη,j , s) ∈ C[0, tη,j ]. Using A1, possess

|K̂l(tη,j , s)| ≤ M,

and
|K̂(tη,j , s)| ≤ H(tη,j , s).

For further detail see [13].
By using of this assumptions, for η = n+ 1, equation (13) can be written



H. A. Thabbah, S. Pishbin, P. Darania / Eur. J. Pure Appl. Math, 18 (3) (2025), 6229 8 of 16

εn+1,j =

∫ tn+1,j

0
K̂(tn+1,j , s)dsεn+1,j

+h
n∑

z=0

∫ 1

0

(
K̂z(tn+1,j , s)(x(tz + sh)− y(tz + sh))

)
ds

+h

∫ cj

0

(
K̂n+1(tn+1,j , s)(x(tn+1 + sh)− y(tn+1 + sh))

)
ds

=

∫ tn+1,j

0
K̂(tn+1,j , s)dsεn+1,j + h

n∑
z=0

∫ 1

0
K̂z(tn+1,j , s)e(tz + sh)ds

+h

∫ cj

0

(
K̂n+1(tn+1,j , s)e(tn+1 + sh)

)
ds.

(16)

From A2, for z = 0, . . . , r − 1, and s ∈ [0, 1), we have

e(tz + sh) = hpqz(s), ||qz||∞ ≤ C1, (17)

and then using the (17), equation (16) reduced to the following form



H. A. Thabbah, S. Pishbin, P. Darania / Eur. J. Pure Appl. Math, 18 (3) (2025), 6229 9 of 16

εn+1,j −
∫ tn+1,j

0
K̂(tn+1,j , s)dsεn+1,j

−h
m∑
i=1

∫ cj

0
K̂n+1(tn+1,j , s)L̃i(s+ 1)dsεn+1,j

−h
m∑
i=1

∫ 1

0
K̂n(tn+1,j , s)L̃i(s)dsεn+1,j

−h
m∑
i=1

∫ cj

0
K̂n+1(tn+1,j , s)L̂i(s+ 1)dsεn,j

−h
m∑
i=1

∫ 1

0
K̂n−1(tn+1,j , s)L̃i(s)dsεn,j

−h

m∑
i=1

∫ 1

0
K̂n(tn+1,j , s)L̂i(s)dsεn,j

= hp+1
r−1∑
z=0

∫ 1

0
K̂z(tn+1,j , s)qz(s)ds

+hp+1
n∑

z=r

∫ 1

0
K̂z(tn+1,j , s)Rm,r,z(s)ds

+hp+1

∫ cj

0
K̂n(tn+1,j , s)Rm,r,n(1 + s)ds

+h
n∑

z=r

r−1∑
k=0

∫ 1

0
K̂z(tn+1,j , s)Lk(s)dsεz−k

+h

n−1∑
z=r

m∑
i=1

∫ 1

0
K̂z(tn+1,j , s)L̂i(s)dsεz,j

+h

n−2∑
z=r

m∑
i=1

∫ 1

0
K̂z(tn+1,j , s)L̃i(s)dsεz+1,j

+h

r−1∑
k=0

∫ cj

0
K̂n+1(tn+1,j , s)Lk(s+ 1)dsεn−k.

(18)



H. A. Thabbah, S. Pishbin, P. Darania / Eur. J. Pure Appl. Math, 18 (3) (2025), 6229 10 of 16

Now, to get the matrix form of equation (18), we define the following vectors

ājn+1 = 1−
∫ tn+1,j

0
K̂(tn+1,j , s)ds, j = 1, 2, ...,m,

E1,z = [εz−r+1, εz−r, ..., εz]
T , E2,z = [εz,1, εz,2, ..., εz,m]T ,

(19)

and matrices
Ān+1 = diag(ā1n+1, ..., ā

m
n+1), (20)

(B
(z)
n+1)j,k =



∫ 1

0
K̂z(tn+1,j , s)Lk(s)ds, z = r, ..., n,

∫ cj

0
K̂z(tn+1,j , s)Lk(1 + s)ds, z = n+ 1,

(21)

(C
(z)
n+1)j,i =



∫ 1

0
K̂z(tn+1,j , s)L̂i(s)ds, z = r, ..., n,

∫ cj

0
K̂z(tn+1,j , s)L̂i(1 + s)ds, z = n+ 1,

(22)

and

(D
(z)
n+1)j,i =



∫ 1

0
K̂z(tn+1,j , s)L̃i(s)ds, z = r, ..., n,

∫ cj

0
K̂z(tn+1,j , s)L̃i(1 + s)ds, z = n+ 1.

(23)

Therefore, by replacing these vectors and matrices in equation (18) and after some com-
putations, we will have

(
Ān+1 − h

(
D

(n+1)
n+1 +D

(n)
n+1

))
E2,n+1 − h

(
C

(n+1)
n+1 + C

(n)
n+1 +D

(n−1)
n+1

))
E2,n

= h
n−1∑
z=r

B
(z)
n+1E1,z + hB

(z)
n+1E1,n+1 + hB

(z)
n+1E1,n + h

n−1∑
z=r

C
(z)
n+1E2,z

+h
n−2∑
z=r

D
(z)
n+1E2,z +O(hp+1).

(24)

Similar to what we did before for equation (13) at point t = tn+1,j , the matrix form of
this equation at point t = tn,j can be expressed as

−hD
(n)
n E2,n+1 +

(
Ān − h

(
C

(n)
n +D

(n−1)
n

))
E2,n

= h
n∑

z=r
B

(z)
n E1,z + h

n−2∑
z=r

D
(z)
n E2,z+1 + h

n−2∑
z=r

C
(z)
n E2,z + C

(n−1)
n E2,n−1 +O(hp+1).

(25)



H. A. Thabbah, S. Pishbin, P. Darania / Eur. J. Pure Appl. Math, 18 (3) (2025), 6229 11 of 16

Now, in the equation (12), we set n = l − 1, then we get the difference equation

E1,l = ÃE1,l−1 + S̃E2,l−1 + T̃E2,l +O(hp), (26)

with exact solution

E1,l = Ãl−r+1E1,r−1 +
l−1∑

j=r−1

Ãl−j−1(S̃E2,j + T̃E2,j+1) +O(hp), (27)

where

S̃ =

 000r−1,m

L̄(1), . . . , L̄m(1)

 , T̃ =

 000r−1,m

L̂1(1), . . . , L̂m(1)

 . (28)

Inserting (27) in equations (24) and (25), we will get the general form of the system
of their equations in the form of the following system

(
Ān+1 − h

(
D

(n+1)
n+1 +D

(n)
n+1 +B

(z)
n+1T̃

))
E2,n+1

−h

(
C

(n+1)
n+1 + C

(n)
n+1 +D

(n−1)
n+1 +B

(z)
n+1S̃ +B

(z)
n+1ÃT̃ +B

(z)
n+1T̃

)
E2,n

= h
n−1∑
z=r

B
(z)
n+1E1,z + hB

(z)
n+1Ã

n−r+2E1,r−1 + h
n−1∑

j=r−1
B

(z)
n+1Ã

n−jS̃E2,j

+h
n−2∑

j=r−1
B

(z)
n+1Ã

n−j T̃E2,j+1 + hB
(z)
n+1Ã

n−r+1E1,r−1

+h
n−1∑

j=r−1
B

(z)
n+1Ã

n−j−1S̃E2,j + h
n−2∑

j=r−1
B

(z)
n+1Ã

n−j−1T̃E2,j+1

+h
n−1∑
z=r

C
(z)
n+1E2,z + h

n−2∑
z=r

D
(z)
n+1E2,z +O(hp+1),

−hD
(n)
n E2,n+1 +

(
Ān − h

(
C

(n)
n +D

(n−1)
n −B

(z)
n T̃

))
E2,n

= h
n−1∑
z=r

B
(z)
n E1,z + hB

(z)
n Ãn−r+1E1,r−1 + h

n−1∑
j=r−1

B
(z)
n Ãn−j−1S̃E2,j

+h
n−2∑

j=r−1
B

(z)
n Ãn−j−1T̃E2,j+1 + h

n−2∑
z=r

D
(z)
n E2,z+1 + h

n−2∑
z=r

C
(z)
n E2,z

+C
(n−1)
n E2,n−1 +O(hp+1).

(29)

Note that, the coefficients matrix on the left hand of (29) has the form:

AAA =

[
a11 a12
a21 a22

]
, (30)



H. A. Thabbah, S. Pishbin, P. Darania / Eur. J. Pure Appl. Math, 18 (3) (2025), 6229 12 of 16

where

a11 = Ān+1 − h
(
D

(n+1)
n+1 +D

(n)
n+1 +B

(z)
n+1T̃

)
,

a12 = −h

(
C

(n+1)
n+1 + C

(n)
n+1 +D

(n−1)
n+1 +B

(z)
n+1S̃ +B

(z)
n+1ÃT̃ +B

(z)
n+1T̃

)
,

a21 = −hD
(n)
n ,

a22 = Ān − h
(
C

(n)
n +D

(n−1)
n −B

(z)
n T̃

)
.

(31)

Now, we set

ζζζ =

[
E2,n+1

E2,n

]
, (32)

and AAAh = AAA+O(h), then we have
AAAhζζζ = bbb, (33)

where vector bbb, is clear and

AAAh =

[
Ān+1 0
0 Ān

]
. (34)

Then, matrix AAAh is invertible if for η = n, n+1, the matrices Āη, are invertible. According
to theorem 3.6 from reference [13], matrices Ān, and Ān+1, are invertible. Therefore, the
bound for ζζζ can be obtained by the same way as described in [13] (Theorem 3.6), then we
obtain

||ζζζ|| ≤ Chp,

and the proof has completes here.

4. Numerical results

Now, we aim to validate the previously established analytical findings by employing
numerical solutions to address two test problems. At the same time, we will tabulate the
attained satisfactory outcomes concerning L∞ errors and optimal convergence orders for
some values of N with m, r = 2. This is done for the purpose of numerically verifying the
underlying theorem 2. We consider two following cases:

1) c1 =
1

2
, c2 = 1 ⇒ ρ(Ã) = 0 < 1.

2) c1 =
1

3
, c2 =

7

8
⇒ ρ(Ã) = 0.505485 < 1.

To ensure clarity, it’s important to mention that all implementations have been carried out
in the Mathematica ® software. Continuing in this section, the primary values have been
gotten from well-established exact solutions. The comparison of the obtained numerical
results in Tables 1- 4, for m, r = 2, shows that new multi-step method (NMSM) is more
accurate than the one-step method (1-SM) and multi-step methods (MSMs) used in [13,
14]. Also in Fig. 1 and Fig. 2, we plot error functions for the one-step, multi-step and



H. A. Thabbah, S. Pishbin, P. Darania / Eur. J. Pure Appl. Math, 18 (3) (2025), 6229 13 of 16

(c1, c2) = (12 , 1) (c1, c2) = (12 , 1) (c1, c2) = (12 , 1)
N 1-SM [13] MSM [14] NMSM

4 5.20× 10−6 5.58× 10−10 1.30× 10−11

8 1.30× 10−6 6.52× 10−11 2.29× 10−13

16 3.25× 10−7 5.34× 10−12 4.21× 10−15

Order 2.00 3.60 5.76

Table 1: Numerical report for example 1.

(c1, c2) = (12 , 1) (c1, c2) = (12 , 1) (c1, c2) = (12 , 1)
N 1-SM [13] MSM [14] NMSM

4 1.70× 10−5 2.04× 10−9 1.77× 10−10

8 4.18× 10−6 4.21× 10−10 9.04× 10−12

16 1.03× 10−6 4.82× 10−11 1.82× 10−13

Order 2.01 3.12 5.62

Table 2: Numerical report for example 2.

new multi-step scheme with N = 16. As can be seen, the accuracy of the numerical results
obtained by the new multi-step method is higher than the other ones.

Example 1. [13] Take for the equation (1) with T = 1, as

x(t) = f(t) +

∫ t

0

s+ t+ 5

(2s+ 4t+ 10)(x2(t) + 1 + x2(s))
ds,

whose exact solution is x(t) =
√
2 + t and

f(t) = (2 + t)
1
2 − ln(

5 + 3t

5 + 2t
)
1
2 .

Example 2. [13]Take for the equation (1) for 0 ≤ t ≤ 1, as

x(t) = (1 + t)
1
2 − 1

2
exp(−(t+ 2t)) +

1

2
exp(−2(t+ 1)) +

∫ t

0

1

2
exp

(
− (y2(t) + y2(s))

)
ds,

whose exact solution is x(t) = (1 + t)
1
2 .

5. Conclusion

In this paper, we studied NMSM to solve the nonstandard Volterra integral equation,
numerically. We have displayed that this collocation approach yields an effective and very



H. A. Thabbah, S. Pishbin, P. Darania / Eur. J. Pure Appl. Math, 18 (3) (2025), 6229 14 of 16

0.2 0.4 0.6 0.8 1.0
-20

-15

-10

-5

t

L
o
g
1
0
(
E
rr
o
r)

:New multi-step

:Multi-step

:One-step

0.2 0.4 0.6 0.8 1.0

-16

-14

-12

-10

-8

-6

-4

-2

t

L
o
g
1
0
(
E
rr
o
r)

:New multi-step

:Multi-step

:One-step

Figure 1: Plot of error functions for the one-step, muti-step and new multi-step scheme with N =
16, c1 = .5, c2 = 1 in Example 1 (left). Plot of error functions related for the one-step, muti-step and
new multi-step scheme with N = 16, c1 = 1

3 , c2 = 7
8 in Example 1 (right).

(c1, c2) = (13 ,
7
8) (c1, c2) = (13 ,

7
8) (c1, c2) = (13 ,

7
8)

N 1-SM [13] MSM [14] NMSM

4 2.01× 10−4 1.23× 10−6 1.97× 10−9

8 5.37× 10−5 1.37× 10−7 4.95× 10−11

16 1.38× 10−5 1.13× 10−8 1.00× 10−12

Order 1.95 3.52 5.62

Table 3: Numerical report for example 1.

0.2 0.4 0.6 0.8 1.0
-20

-15

-10

-5

t

L
o
g
1
0
(
E
rr
o
r)

:New multi-step

:Multi-step

:One-step

0.2 0.4 0.6 0.8 1.0

-16

-14

-12

-10

-8

-6

-4

-2

t

L
o
g
1
0
(
E
rr
o
r)

:New multi-step

:Multi-step

:One-step

Figure 2: Plot of error functions for the one-step, muti-step and new multi-step scheme with N =
16, c1 = 0.5 and c2 = 1 in Example 2 (left). Plot of error functions related for the one-step, muti-step
and new multi-step scheme with N = 16, c1 = 1

3 and c2 = 7
8 in Example 2 (right).



H. A. Thabbah, S. Pishbin, P. Darania / Eur. J. Pure Appl. Math, 18 (3) (2025), 6229 15 of 16

(c1, c2) = (13 ,
7
8) (c1, c2) = (13 ,

7
8) (c1, c2) = (13 ,

7
8)

N 1-SM [13] MSM [14] NMSM

4 5.04× 10−4 6.36× 10−6 3.92× 10−8

8 1.42× 10−4 8.73× 10−7 1.39× 10−9

16 3.80× 10−5 9.21× 10−8 1.39× 10−11

Order 1.90 3.25 5.31

Table 4: Numerical report for example 2.

accurate numerical method for the approximation of solutions to NVIE. The proposed
scheme exhibited satisfactory global convergence, with the observation that the conver-
gence orders were primarily influenced by the collocation parameter and the conditions
for ρ(Ã).

References

[1] M. Bohner and A. Peterson. Advances in Dynamic Equations on Time Scales.
Birkhäuser, Boston, MA, 2003.

[2] F. Brauer and C. Castillo-Chavez. Mathematical Models in Population Biology and
Epidemiology, volume 40 of Texts in Applied Mathematics. Springer, New York, 2001.

[3] D. B. Pachpatte. Properties of solutions to nonlinear dynamic integral equations on
time scales. Electronic Journal of Differential Equations, 136:1–8, 2008.

[4] D. B. Pachpatte. On a nonstandard volterra type dynamic integral equation on time
scales. Electronic Journal of Qualitative Theory of Differential Equations, 72:1–14,
2009.

[5] C. C. Tisdell and A. Zaidi. Basic qualitative and quantitative results for solutions to
nonlinear, dynamic equations on time scales with an application to economic mod-
elling. Nonlinear Analysis, 68:3504–3524, 2008.

[6] P. Darania, S. Pishbin, and A. Ebadi. Convergence analysis of multi-step collocation
methods to solve generalized auto-convolution volterra integral equations. Numerical
Analysis and Applications, 16(2):123–134, 2023.

[7] U. Kosel and L. von Wolfersdorf. Nichtlineare integralgleichungen. In Seminar
Analysis, pages 93–128, Berlin, Akademie der Wissenschaften der DDR, 1986. Karl-
Weierstrass-Institut für Mathematik.

[8] L. Liu and J. Ma. Collocation boundary value methods for auto-convolution volterra
integral equations. Applied Numerical Mathematics, 177:1–17, 2022.

[9] L. von Wolfersdorf. Einige klassen quadratischer integralgleichungen. Sitzungs-
berichte der Sächsischen Akademie der Wissenschaften zu Leipzig, Mathematisch-
Naturwissenschaftliche Klasse, 128(2):1–34, 2000.

[10] L. von Wolfersdorf. Autoconvolution equations and special functions. Integral Trans-
forms and Special Functions, 19:677–686, 2008.



H. A. Thabbah, S. Pishbin, P. Darania / Eur. J. Pure Appl. Math, 18 (3) (2025), 6229 16 of 16

[11] L. von Wolfersdorf. Autoconvolution equations of the third kind with abel integral.
Journal of Integral Equations and Applications, 23:113–136, 2011.

[12] R. Zhang, H. Liang, and H. Brunner. Analysis of collocation methods for generalized
auto-convolution volterra integral equations. SIAM Journal on Numerical Analysis,
54(2):899–920, 2016.

[13] Q. Guan, R. Zhang, and Y. Zou. Analysis of collocation methods for nonstandard
volterra integral equations. IMA Journal of Numerical Analysis, 32:1755–1785, 2012.

[14] S. Pishbin and A. Ebadi. High-order convergence of multistep collocation methods for
nonstandard volterra integral equations. International Journal of Computer Mathe-
matics, 100(4):824–837, 2023.

[15] R. P. Agarwal, D. O’Regan, and P. J. Y. Wong. Positive Solutions of Differential,
Difference and Integral Equations. Kluwer Academic, Dordrecht, 1999.

[16] H. Brunner. Collocation Methods for Volterra Integral and Related Functional Differ-
ential Equations. Cambridge University Press, Cambridge, 2004.

[17] T. A. Burton. Volterra Integral and Differential Equations. Academic Press, New
York, 1983.

[18] A. J. Jerri. Introduction to Integral Equations with Applications. Wiley-Interscience,
New York, 1999.

[19] D. O’Regan and M. Mehhan. Existence Theory for Nonlinear Integral and Integro-
Differential Equations. Kluwer Academic, Dordrecht, 1998.

[20] G. Gripenberg, S. O. Londen, and O. Staffans. Volterra Integral and Functional
Equations. Cambridge University Press, Cambridge, 1990.

[21] H. Brunner. Volterra Integral Equations: An Introduction to Theory and Applications.
Cambridge University Press, Cambridge, UK, 2017.

[22] F. Ghoreishi and M. Hadizadeh. Numerical computation of the tau approximation
for the volterra-hammerstein integral equations. Numerical Algorithms, 52:541–559,
2009.

[23] T. Apostol. Mathematical Analysis. Addison-Wesley, Reading, MA, 1974.
[24] D. Conte, Z. Jackiewicz, and B. Paternoster. Two-step almost collocation methods

for volterra integral equations. Applied Mathematics and Computation, 204:839–853,
2008.

[25] S. Fazeli, G. Hojjati, and S. Shahmorad. Super implicit multistep collocation methods
for nonlinear volterra integral equations. Mathematical and Computer Modelling,
55:590–607, 2015.