A numerical scheme for the solution of neutral integro-differential equations including variable delay

In this study, an effective numerical technique has been introduced for finding the solutions of the first-order integro-differential equations including neutral terms with variable delays. The problem has been defined by using the neutral integro-differential equations with initial value. Then, an alternative numerical method has been introduced for solving these type of problems. The method is expressed by fundamental matrices, Laguerre polynomials with their matrix forms. Besides, the solution has been obtained by using the collocation points with regard to the reduced system of algebraic equations and Laguerre series.


Introduction
Delay differential equations in neutral type and integro-differential equations are of an attractive interest in many applications in science and engineering. In applied mathematics, they have an increasing enthusiasm by their implementation in dynamical systems, electrodynamics and mechanics as well. In recent years, numerical treatment of the neutral type integro-differential equations has been arised [1]. These type of equations occur in mechanics, physics, technical problems such as progress for showing cutting and infeed grinding. Besides, they describe some procedures in chemistry and physics for reactors. Moreover, some well-known biological processes growth, death and birth are determined by neutral type equations. In ecological phenomena, some models with respect to the neutral type integro-differential equations are used for evolution equations of single species [2]. Furthermore, many applications of these types of equations exist in medicine. For instance, sugar quantity in blood is modeled by using them; immunology, epidemiology, cancer chemotherapy may explain different aspects of human body interaction with diseases. As another example, in Fig. 1 model of pressure regulation of model of arterial blood is represented by the formulation of functional differential equations. In this figure, arterial vessels are shown by A and B. These arterial vessels have connections to each other, blood flow from A to B has a rate which is shown by Q, R represents peripheral resistance, and the heart productivity is denoted by Q h . The incoming liquid rate is Q in while the outcoming rate is Q out [3]. Furthermore, epidemic of the human immunodeficiency virus (HIV) is modeled by the system of functional integro-differential equations. Besides, fishing process, river pollution control can be described and these examples are ecological applications in the ecology field.
Neutral type delay integro-differential equations under the initial condition are of numerical solutions which have been reached by many authors. Such problems have difficulties in motivation but also often appear at numerical investigations. Functional integro-differential equations including Volterra type integrals have been searched with regard to numerical analysis aspect by Brunner [4]. Collocation methods have been applied for solving Volterra functional integral equations including non-vanishing delays [5]. Besides, in the literature, continuous spline collocation methods [6], Lagrange interpolation and Chebyshev interpolation [7], Adams-Moulton method [8], backward substitution method [9], continuous Runge-Kutta method [10][11][12], Spectral method [13] have been implemented in order to find the solutions of neutral type integro-differential equations including variable delays.
In this work, the following first-order integro-differential equation including neutral terms and variable delays is considered as under the initial condition where P 0 (t), P 1 (t), u(t), v(t), g(t) and the delay term (t) are defined as continuous functions for ∞ > b > t ≥ 0 . Here the aim is about finding a numerical solution y N (t) by using the truncated Laguerre series of given problem (1)- (2): where a n , n = 0, 1, ..., N are unknown coefficients; L n (t) are the Laguerre polynomials for n = 0, 1, ..., N and defined as

Numerical method
In this section, the numerical method based on Laguerre polynomials is introduced. The main advantage of the method lies in its straightforwardness since it has no aim for discretization as a reliable tool. On the other hand, Laguerre polynomials give powerful solutions, especially on the positive interval since its applications in several fields with these types of properties.

Fundamental relations
In here, the matrix forms of Eq. (1) are composed. First of all, it is organized as where Then, the matrix form of the numerical solution (3) is considered as where here (t) is defined as in the matrix form: where Thus, the connection between (t) and � (t) is defined as in which Therefore, the matrix relations (7) and (8) are used and � (t) is defined as follows So that, from (7) and (9) is obtained in which is described as Thus, from (6), (7), and (10) [14]: By replacing t → t − (t) into (11) and by using (9), we obtain where Diversely, the kernel function K(t, s) is obtained by using the Taylor series as Then, its matrix form is written as By means of the relations (9) and (13), the integral part of Eq. (1) has a matrix form as where and Similarly, from (11), the initial condition in Eq. (2) has the matrix form as Therefore, we obtain the matrix forms of D 1 , D 2 , D 3 and I in Eq. (5), from (6), (12), (11), and (14), respectively, as So that, Eq. (1) can be represented by the following matrix equation as

Method of solution
In this section, the collocation points are defined by Thus, the collocation points (17) are substituted into Eq. (16) and the fundamental matrix is obtained asor where Briefly, Here we have the augmented matrix form of Eq. (18). Besides, we consider the initial condition which is given in Eq. (2). Its matrix form is defined in Eq. (15). By putting the row from (15) in place of the last row of (19), Then, we construct the new augmented matrix as In Eq. (20), if rank ̃ = rank ̃ ;̃ = N + 1 , then the coefficients matrix is uniquely determined with the help of Gauss Elimination procedure [15][16][17][18][19][20][21][22]. Then, by using Eq. (3), the problem (1)-(2) is solved numerically and its numerical solution is obtained as in the form

Analysis of the method
In this section, some properties of the method are introduced. The stability of the collocation methods has been investigated previously by Brunner et al. [23,24]. Besides, the stability of the related problem has been also considered by some authors [25,26]. However, existence and uniqueness theorems and convergence of the method are given as follows.

Definition 1 Let y N (t) be the approximate solution of the problem (1)-(2) which has an exact solution as y(t) . Then, the collocation method is said to be convergent if an only if
where h is the step size and h → 0, N → ∞ [27][28][29][30][31].

Existence and uniqueness
Theorem 3 Consider the first-order integro-differential equation including neutral terms and variable delays in Eq. (1) under the initial condition (2). Assume that g(t) and K(t, s) are continuous functions which satisfy the Lipschitz conditions. Then, for every | | t − t 0 | | ≤ a and | | s − t 0 | | ≤ a for any initial value t 0 and the positive constant a > 0 , ‖ ‖ y 1 ‖ ‖ < ∞ and ‖ ‖ y 2 ‖ ‖ < ∞ . Then, the IVP has a unique solution.

Accuracy
In this section, the accuracy of the approximate solution is investigated. As an important factor for the numerical methods in the literature, the approximate solutions are corrected with regard to the residual error correction procedure. On the other hand, a brief error analysis is given in order to reach the approximation for the problem (1)-(2).

Residual correction
In here, the residual correction is given in order to improve the solutions and for the comprehensive error analysis for the approximate solution of the problem (1)-(2) [17,18,[32][33][34][35][36][37]. Now, let us consider the er ror function as e N (t) = y(t) − y N (t) . Then, we construct an error problem in the form Subsequently, the numerical method is applied on the error problem in Eq. (26) and we have the approximate solution for the error function in the form as follows Thus, the corrected approximate solution is obtained as

Error analysis
In this section, we investigate the absolute error function E N (t) for t = t p ∈ [a, b], p = 0, 1, 2, .... and the accuracy of the numerical solutions is checked. Specifically, if E N (t q ) ≤ 10 −k q ( k q any positive integer) is small enough, then the approximation has its reliability.

Algorithm
In here, the present method is shown by its algorithm. The steps are explained clearly in order to see the implementation of the computer programming part of the work.
Step 0 Input initial data: Construct the matrices such as (t), , Step 3 Replace in the fundamental equation Step 4 Put the collocation points, Compute the matrix for the initial condition Step 7 Substitute the outcome from S6 into the matrix in S5 and get ̃ ;̃ Step 8 Determine the system in S7. Output: y N (t) Step 9 Accuracy check:E N (t) Else back S1

Numerical experiments
Numerical experiments section gives us an idea about the method and its applicability, validity and reliability. The applications of this method have been implemented by using some numerical illustrations. Here Maple and MATLAB computer programs are used for the calculation algorithm and plotting.

Example 1
Firstly, the first-order Volterra integro-differential equation including neutral term and variable delay.
under the initial condition Hereby approximate solution of (30)-(31) is found by using the algorithm of the present method for N = 4.
The collocation points (17) are set for N = 4 on the interval [0,1] which are as follows.
Then, Eq. (1) has its fundamental matrix equation as   So, the augumented matrix is found as Now, the initial condition is defined in (31) has the matrix form as So that the new augmented matrix is constructed as By following the procedure, the system is solved, and then, a n Laguerre coefficients are found as follows Consequently, these calculated coefficients are established and they are replaced in Eq. (32). Then, the exact solution of (30)-(31) is acquired as (33) y(t) = t + 1.  Example 2 As a second example, the first-order Volterra integro-differential equation including neutral term and delay is considered [37] and the initial condition corresponding to the exact solution y exact (t) = e t . Correspondingly, the problem (34)- (35) is solved by the similar approach which is followed in Example 1. In Table 1, absolute errors for different N truncation values are demonstrated for their comparison. Herein this error comparison with different numerical methods: Laguerre collocation method (LCM) for M = 10 and N = 9 and mixed spline/spectral method (MSSM) with h = 0.2 in Table 2 and in Fig. 1. From these comparisons, we can see apparently that more suitable and efficient results are obtained for the smaller when we have increasing N value.

Example 3
Consider the first-order Volterra integro-differential equation including neutral term and delay from Eq. (1) with P 1 (t) = v(t) = g(t) = 0 , K(t, s) = u(t) = 1 and including a delay term y(t − 1) together with P 0 (t) = 1 . Moreover, the initial condition is given as y(0) = 1 and the exact solution is y exact (t) = e t . Then, the absolute errors for different N truncation values are demonstrated for their comparison in Table 3 and Fig. 2. Besides, we can see the numerical solutions with regard to Laguerre approach and exact solution of the problem in Fig. 3. From these comparisons, we can see apparently that more suitable and efficient results are obtained for the smaller when we have increasing N value. Subsequently, there are improved results with related to the (34) Sect. 4.1. Residual correction helps for a better approximation and give smaller values for the error function which is shown in Table 4.

Conclusion
In this study, a powerful numerical technique to determine the numerical solutions of first-order integro-differential equations including neutral terms and variable delays is proposed. The technique affords approximate solutions of the problem which are mainly close enough to the exact solutions with respect to N . Accuracy and applicability of the method have been proved by the visible results at the tables and the figures. The main advantages of this technique are including but not limited to straightforward coding, its apparent algorithm and accessible matrix calculations. Moreover, the error analysis including the residual correction supports the results with the additional error problem solution which is explained and applied in an example.
As a future outlook, this numerical study and the technique can be extended to other models with related to Volterra integro-differential equations including retarded term. However, some modifications are required [38].
Funding Open Access funding enabled and organized by Projekt DEAL.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.