A NUMERICAL METHOD ON BAKHVALOV SHISHKIN MESH FOR VOLTERRA INTEGRO-DIFFERENTIAL EQUATIONS WITH A BOUNDARY LAYER

We construct a finite difference scheme for a first-order linear singularly perturbed Volterra integro-differential equation (SPVIDE) on Bakhvalov-Shishkin mesh. For the discretization of the problem, we use the integral identities and deal with the emerging integrals terms with interpolating quadrature rules which also yields remaining terms. The stability bound and the error estimates of the approximate solution are established. Further, we demonstrate that the scheme on Bakhvalov-Shishkin mesh is O(N−1) uniformly convergent, where N is the mesh parameter. The numerical results are also provided for a couple of examples.


Introduction
In this present work, we are specifically consider the following class of the singularly perturbed linear Volterra integro-differential equations (SPVIDEs) subject to where 0 < ε ≪ 1 is a small perturbation parameter. We assume a(x) ≥ α > 0, f (x)(x ∈ I) and K(x, t)((x, t) ∈ I × I) are sufficiently smooth functions such finite difference schemes are developed to solve linear first order SPVIDEs with delay in [4], [27]. In [3], using a fitted difference operator a second-order difference scheme is constructed on a piecewise uniform mesh to solve linear SPVIDEs.
In this present work, we mainly construct a uniform convergent difference scheme on a Bakhvalov-Shishkin mesh for the problem (1)- (2). Bakhvalov-Shishkin mesh is a mixed version of the Shishkin mesh and Bakhvalov mesh which are known to yield accurate results for singularly perturbed problems with boundary layers. In [30], the author demonstrated that the results from an upwind difference scheme on Bakhvalov-Shishkin mesh applied to a linear convection-diffusion equation are more accurate than the results from the upwind scheme on a Shishkin mesh. Further, a finite difference scheme on Bakhvalov-Shishkin mesh is utilized to deal with a singularly perturbed boundary value problem in [10].
The rest of the paper is organized in the following order. In Section 2, the asymptotic estimates on the exact solution to (1)- (2) are established. In Section 3, we define the Bakhvalov-Shishkin mesh points according to the boundary layer conditions of the problem (1)-(2) and derive a finite difference scheme utilizing the integral identities with exponential basis functions and then applying interpolating quadrature rules provided in [1] to the integral terms. In Section 4, we establish the stability bounds and the error estimates of the numerical solution and as a result we show that the scheme demonstrates O(N −1 ) uniform convergence with respect to the perturbation parameter. We also provide the numerical results in Section 5.

Asymptotic Behavior of the Solution
In the following lemma, we establish a priori estimates for the asymptotic behavior of the solution to the problem (1)-(2).

Lemma 1.
Let a, f ∈ C(I) and K ∈ C(I × I). The solution u to the problem where In addition, if a, f ∈ C 1 (I) and K ∈ C 1 (I × I) with then the solution u(x) satisfies Proof. To establish the first estimate given in (3) we start by rewriting (1) as where Solving the equation (6) with u(0) = A yields x ξ a(s)ds dξ, and further we calculate Here, by the definition of F (x) in (7), we get Substituting (9) into (8) yields We integrate by parts the last term with double integral here An application of the Gronwall's inequality to (10) provides which leads to the desired result in (3).
For the next estimate provided in (5), we first differentiate the equation (1) and have we have In a similar manner to the previous work above, we solve (12) v x ξ a(s)ds dξ.
Then, we have Here, by the formula of g(x) given in (11), from (3) and knowing that a, f ∈ C 1 (I), K ∈ C 1 (I × I) and from (4) we obtain which implies ||g|| ∞ ≤ C * for a C * ∈ R. Hence, utilizing this estimate on g(x) in (13) provides On the other hand, inserting x = 0 in (1) and since a, f ∈ C 1 (I) it follows that Substituting this into (15) yields which provides the desired result. □

3.1.
Notation. Before we proceed to the definition of the mesh points and discretization of the problem we provide the notation we use throughout the paper.
For any continuous mesh function v i defined on ω h we use the notation 3.2. Discretization. In this section, we construct our difference scheme based on Bakhvalov-Shishkin mesh. According to this mesh construction, we divide the domain into two subintervals [0, σ] and [σ, ℓ], where σ is the transition parameter. For a positive even discretization parameter N , we determine the transition parameter σ as We assume ε ≪ N −1 as it is used in practice. We define a set of mesh points as the following To derive the difference approximation, we use the following integral identity with the exponential basis function We remark that φ i solves the equation To obtain the difference scheme from (18), we proceed by evaluating the integrals term by term applying the interpolating quadrature rules with weight functions and obtain the remainder terms as provided in [1]. In the following, we handle the differential term on the left-hand side of (18), and Further, applying the first quadrature rules provided in [1] to the integral term in (18) twice we obtain where and T 0 (λ) = 1 for λ ≥ 0 and T 0 (λ) = 0 for λ < 0. Here, we apply the composite right-side rectangle rule to the integral term in the right-hand side of (27) and get where Then, inserting (25) in (23) provides i .
(27) On the other hand, the right-hand side of (18) gets the in the form where Inserting the relations (20), (27) and (28) in (18), we obtain the difference problem for the problem (1)-(2) as where As a result, neglecting the error term R i in (30) provides the following difference scheme where θ i defined by (21).

Stability, Error Estimates and Convergence Results
Here, we establish the stability bound and the error estimates of the approximate solution y. Further, the convergence of the difference scheme provided in (32)-(33) is analyzed.
Lemma 2. Assume that |F i | ≤ F i and F i be a non-decreasing function. The solution to the problem Proof. The proof follows from the maximum principle for difference operators. Details can be found in [27]. □ Lemma 3. Let y i be the solution of the problem (32)- (33). Then, y i satisfies Proof. The difference scheme equation given in (32) can be rewritten in the form where For F i , we have the estimate Then, applying Lemma 2 to (35) and utilizing this estimate provide Further, applying the difference analogue of the Gronwall's inequality to (36) we have which yields the result in (34). □ The error of the difference problem is given by the solution to the problem Lemma 4. Suppose that z i be the solution of (37)- (38). Then, z i holds the estimate Proof. The result follows from Lemma 3 taking A = 0 and f = R. □ Lemma 5. Let a, f ∈ C 1 (I) and K ∈ C 1 (I × I) with and Then, the truncation error R i satisfies the estimate Proof. To establish the estimate given in (43), we proceed by bounding each term in R i provided in (31). For R where s ∈ [x, x i ] comes from the Mean Value Theorem. Then, since a ∈ C 1 (I) and from (3) we get |R Further, for R (2) i we take into account of (40), (41) and |T 0 (λ)| ≤ 1, so Then, applying the Leibnitz formula to (45) yields On the other hand, by the Leibnitz formula and from (40), (42) and (5) we have Then, by the Mean Value Theorem applied to the exponential term in (47) with where h * = max 1≤j≤i h j . Lastly, for R (4) i , similarly to the work above and since f ∈ C 1 (I) we have where s ∈ [x i−1 , x i ] by the Mean Value Theorem. Further in the proof, we need to evaluate each estimate above on the sub-intervals [0, σ] and [σ, ℓ]. For this, we first establish the bounds on the step-size h i on each interval. In the first sub-interval [0, σ] with σ ≤ ℓ 2 , and hence, Then, we apply the Mean Value Theorem to h i with i * ∈ [i − 1, i] and get In the second sub-interval [σ, ℓ], we have where σ ≤ ℓ 2 and Inserting the bounds (50) and (51)  Let u be the exact solution of (1)-(2) and y be the solution of (32)- (33). If the assumptions on the functions a, f and K from Lemma 5 hold, then Proof. The proof follows from Lemma 4 and Lemma 5. □

Algorithm and Numerical Results
In this section, we present the numerical results on an example with an exact solution and an example with an unknown solution. The results include graphs of the approximate solutions, error estimates and the convergence values of the approximate solution to the exact solution. In our algorithm, we consider the following elimination method where y (0) i is the initial process.
Example 1. We study the following initial value problem The exact solution of this problem is The exact error is calculated by the formula where y N is the numerical approximation of u for different N and ε values. We compute the convergence rate by r N = ln e N /e 2N ln 2 .
In Table 1, we provide the errors e N , e 2N and the convergence rates of the approximate solution for various N and ε = 2 −i values.

Example 2. Consider the following test problem
The exact solution to this problem is not known. To compute the approximate solution and estimate the errors, we utilize the double mesh principle, that is calculating the error of the approximate solution on mesh size N with the approximate solution where y N is the approximate solution on mesh N and y 2N is the approximate solution on mesh 2N . The convergence rate is calculated as it is in Example 1.
In Table 2, the errors and the convergence rates of the approximate solution for various N and ε = 2 −i values are presented.

Conclusion
To sum up, we constructed a finite difference scheme on a Bakhvalov-Shishkin mesh to obtain the numerical solution of an initial value problem for a linear firstorder singularly perturbed Volterra integro-differential equation with a boundary layer. We proved that the method is first-order uniformly convergent with respect to the perturbation parameter. As we can see in Table 1, Table 2 and Figure 1, the numerical results of the test problems are also consistent with the analysis on the error estimates and convergence order and hence, it is confirmed that the convergence order of the scheme O(N −1 ). For future work, we suggest that this difference scheme method on Bakhvalov-Shishkin mesh can be applied to the singularly perturbed linear or non-linear problems with delay to obtain accurate numerical solutions. Further, our proposed scheme can be modified to handle integro-differential equations with fractal derivatives which are studied in [8].
Author Contribution Statements All authors contributed equally to this work, and they read and approved the final manuscript.
Declaration of Competing Interests The authors declare that they have no known competing financial interest or personal relationships that could have appeared to influence the work reported in this paper.