Pell–Lucas collocation method for numerical solutions of two population models and residual correction

Our aim in this article is to present a collocation method to solve two population models for single and interacting species. For this, logistic growth model and prey–predator model are examined. These models are solved numerically by Pell–Lucas collocation method. The method gives the approximate solutions of these models in form of truncated Pell–Lucas series. By utilizing Pell–Lucas collocation method, non-linear mathematical models are converted to a system of non-linear algebraic equations. This non-linear equation system is solved and the obtained coefficients are the coefficients of the truncated Pell–Lucas serie solution. Furthermore, the residual correction method is used to find better approximate solutions. All results are shown in tables and graphs for different values, and additionally the comparisons are made with other methods from. It is seen that the method gives effective results to the presented model problems.


Introduction
The mathematical models in population biology have a very important role in applied mathematics and science. These models, non-linear differential equations and their systems usually, are practically beneficial. In research into living systems, mathematics, physics, biology and chemistry are essential for the development of an integrative perspective. There are biotic interactions between each organism and its environment. These interactions affect the growth and destiny of a species or community. Population dynamics are related to population growth in a short time under certain conditions. In this short period, the variation in the number of individuals in the population depends on some interactions. Examples of these interactions include reproductive and food supply, death rates, climate change and predators, competition, parasitism, mutualism, disease and social context.
The number of individuals in a population is constantly changing in the continuous time approach. The most common type of modelling is the interaction between a species and the environment. Such models are very convenient for determining the destiny of a small number of interactive species or the destiny of a single population. These models were pioneered by Pierre-Fran Verhulst [1][2][3], with the introduction of the logistic model in the nineteenth century. Also, in the first quarter of the twentieth century, these models were pioneered by Vito Volterra [4,5], with the introduction of a model to describe, the cycling behaviour of communities of herbivore and carnivore fishes, qualitatively. Recently, some model problems related to HIV infection, population growth, reactiondiffusion, gas dynamic and the non-linear porous fin problem have been solved by various numerical methods in [6][7][8][9][10]. Also, the differential equations and integral equations has been solved by numerical techniques based on various polynomial bases [11][12][13][14][15][16][17][18][19][20][21][22][23]. Besides these, lately, there are some recent developments that have obtained numerical results for non-linear PDEs [24][25][26][27][28][29][30].
In this study, we will first discuss the logistic equation Equation (1) is the logistic equation introduced firstly by Verhulst and then studied by R. Pearl and L. J. Reed (1920) [31]. The population growth of a particular species is modelled by this equation over time. In the logistic Equation (1), where the exact solution is Here, Y(t), Y (t) and K represent respectively the size of the population at the time t, the rate of change of population size and the carrying capacity of the environment. While r is the positive constant, λ is the suitable constant.
In 1910s, Lotka [32] developed population models in which two animals are in a predator-prey relation. In 1926s, Volterra [5] independently discovered the same model, making a statistical analysis of fishing in the Adriatic Sea. Secondly, we will consider the Lotka-Volterra equations and in the interval 0 ≤ t ≤ b. Here, a, b, c, d, λ 1 and λ 2 are suitable constant. Here, Y 1 (t) and Y 2 (t) express the prey and predator populations at time t, respectively. See [33] for details. Hitherto, the decomposition method [34], He's homotopy perturbation method [35], the Bessel collocation method [36], the Chebyshev collocation method [37], the Galerkin-like approach [38] have been used for solving logistic Equation (1) numerically.
Besides, the Adomian decomposition method (ADM) and the power series method [39,40], the homotopy perturbation method (HPM) and the decomposition method [34,41], the variational iteration method (VIM) [42,43], the Chebyshev collocation method [37], the Galerkin-like approach [38] and the Bessel collocation method [36] have been studied for numerical solutions of the Lotka-Volterra Equation (3). Also, the Pell-Lucas collocation method [44] has been performed to solve high-order functional differential equations. In addition, a system of delay differential equations has been solved by Boubaker collocation method [52]. There are similar models in studies such as [45][46][47].
We purpose to discover the approximate solutions of model (1)-(2), expressed as truncated Pell-Lucas series form and we purpose to discover the approximate solutions of the model (3), expressed as truncated Pell-Lucas series form Here, N ≥ 1, a n and a i,n (i = 1, 2) are the Pell-Lucas coefficients and Q n (t) are the Pell-Lucas polynomials described by The aim of this study is to enrich the existing studies in the literature for numerical solutions of the equations or systems of the equations discussed in science and engineering. Thus, the results are shedded light on studies in science and engineering. Let's summarize this article briefly: In Section 2, we present some fundamental relations about the Pell-Lucas polynomial. In Section 3, we develop a collocation technique based on Pell-Lucas polynomials to obtain approximate solutions of the single species model and the Lotka-Volterra model. In Section 4, we submit an error estimation with the help of the residual functions. Moreover, we improve the approximate solutions. We examine two numerical examples related to the presented method in Sections 3 and 4 and we show the outcome of these examples with the help of tables and graphs in Sections 5. In addition, we check the numerical conclusions in this technique by different techniques in the literature in Section 6. After all, in Section 7, we sum up the conclusions of this research.

Fundamental definitions
Evaluation of polynomials can be done easily, as they are simple and smooth functions. Also, polynomials are often handled in numerical analysis for polynomial interpolation. Lucas numbers and polynomials are widely used in the study of many topics such as combinatorics, approximation theory, algebra, geometry, graph theory and number theory. Pell-Lucas polynomials are also a special type of generalized Lucas polynomials. In order to comment on whether Pell-Lucas polynomials have advantages over other generalized Lucas classes, it is necessary to apply the method to other polynomials and discuss the results. However, it can be said that similar results can be obtained. One advantage of using this polynomial is that it offers very efficient solutions. Another advantage of using this polynomial is the ease of application. Also, the generation of these polynomials is easy.
Pell-Lucas polynomials will be used in this study. Some important features of Pell-Lucas polynomials are as follows [48][49][50]: • The Pell-Lucas polynomials Q n (t) are described by the repetition correlation • The generating function of the Pell-Lucas Q n (t) is • The Pell-Lucas polynomials Q n (t) can also be given expressly as follows • Pell-Lucas numbers are explicitly given by the Binettype formulas as • The Pell-Lucas numbers are also equal where F n (t) is a Fibonacci polynomial. • For the convergence of the Pell-Lucas expansion, see [51].

Pell-Lucas collocation method
In this section, a method will be developed to find approximate solutions of the models (1)- (2) and (3). For this purpose, firstly, the matrix representation of Pell-Lucas polynomials will be found. After, the collocation points will be chosen. Consequently, the approximate solution will be achieved by giving the relations between the matrix relations and collocation points required for the model (1)-(2) in Section 3.1. Similarly, the approximate solution will be obtained by giving the relations between the matrix relations and collocation points required for the model (3) in Section 3.2. The Pell-Lucas polynomials Q n (t) can be represented in matrix form by using Equation (6) as follows or where Now, the collocation points required for the Pell-Lucas collocation method are defined as

Solution method for the model (1-2)
In this section, Pell-Lucas collocation method will be established to obtain approximate solutions of the model (1)- (2). Therefore the single species model can be shown as Suppose the model (10) has an approximate solution in the truncated Pell-Lucas series form (4). In that case, Equation (4) can be converted into a matrix form as Now, let's construct the matrix representation of the Y (t) expression in model (10). For this purpose, if the derivative of Equation (12) is taken, there is the relation Here, if the derivative of T(t) = [ 1 t t 2 · · · t N ] is taken, it is obtained and it can be written as where As for the matrix representation of the Y 2 (t), it becomes wherē . .
As the next step for the method, the collocation points (9) are written in Equations (12), (13) and (15). So, it is obtained the relationships These three expressions, respectively, can be represented in compact form and HereȲ =TMĀ. All these matrices are expressed as . . .
If the collocation points t i is written instead of t, in model (1), it is obtained or along with matrix representations it becomes Here, If Equations (16), (17) and (18) are substituted in Equation (19), then it is obtained or briefly Here, The important point here is that the dimensions of the matrices T,T,M R, K, andĀ are (N + 1) On the other hand, if 0 is written instead of t in Equation (12), the matrix representation of the initial condition depended on the matrix of Pell-Lucas coefficients is obtained as or briefly Here, As the last step for the method, the one row of the increased matrix [W; G] is discarded. The matrix [U; λ] is written instead of the line that is discarded. Thus, we have the new increased matrix Equation (23) represents a non-linear algebraic system. This system is solved with the help of Matlab. It should be noted that Matlab solves this non-linear system of equations using the trust-region-dogleg algorithm.
The trust-region-dogleg algorithm is a subspace trustregion method and is based on the interior-reflective Newton method. Thus the Pell-Lucas coefficient matrix A is obtained. This matrix is the coefficient matrix in the solution form required for the method. As a result, Pell-Lucas polynomial solution is

Solution method for the model (3)
In this section, Pell-Lucas collocation method will be established to obtain the solutions that provide the model (3). Accordingly, the predator-prey model can be shown in the form Suppose the model (25) has an approximate solution in the form of a truncated Pell-Lucas series (6). In that case, Equation (5) can be converted into a matrix form as (26) or if Equation (8) is used, then it can be written as or briefly where and let's put the collocation points (9) in Equation (28) for s = 0, 1,..N, then we have . . .
Now, let's construct the matrix representation of the Y i (t) expression in model (25). For this purpose, if the derivative of Equation (28) is taken and Equation (14) is used, we obtain or briefly where and let's put the collocation points (9) in Equation (31) for s = 0, 1,..N, then it is obtained . . .
As the next step for the method, let's construct the matrix representation of the model (25) as where and let's put the collocation points (9) in Equation (33) for s = 0, 1, . . . N, then we have or it can be expressed as an equivalent in the form of Here, Similar operations should be done for theȲ statement. For this purpose, the matrix Y 1,2 (t s ) in Equation (34) is written as All these matrices are expressed as . . .
The matrix 0 is the zero matrix in dimensional (N + 1)x(N + 1) If Equations (29), (32), (36), (37) and (38) are written in Equation (35), then it is obtained or briefly where The important point here is that the dimensions of the matrices T,T,B,M,M,Ā 1 ,H,T,L,M, A and On the other hand, by writing t → 0 in the solution form (27), the matrix representation of the initial conditions depended on the matrix of Pell-Lucas coefficients becomes (41) or briefly Here, As the last step for the method, the two row of the increased matrix [W; G] is replaced by the matrix [U; λ]. So, the new increased matrix is obtained as Equation (42) represents a non-linear algebraic system. Pell-Lucas coefficients a i (i = 0, 1, . . . , N) required for the method are also found by solving this system. This non-linear system of equations is solved by Matlab using the trust-region-dogleg algorithm. This algorithm is also based on the interior-reflective Newton method. As a result, Pell-Lucas polynomial solutions are

Residual error estimation and improvement of solutions
In this section, a method known as residual correction will be presented for Sections 3.1 and 3.2, respectively. For this purpose, these discussed models firstly are expressed in the form of L[Y(t)] = g(t). Then, the approximate solutions are replaced by the solutions in and an error problem based on the approximate solution is generated. This error problem is solved according to the method in Section 3 and estimated errors are found. As a next step, approximate solutions to estimated errors are added and so improved approximate solutions are obtained. Finally, improved errors are obtained by taking the difference of the improved approximate solutions with the exact solution of the problem.

The case of single species model
In this section, the method described in the steps in Section 4 will be applied to the model (10). For this purpose, Equation (1) is expressed as The approximate solution Y N (t) is written in L[Y(t)] = g(t), since the approximate solution obtained in Section 3.1 provides the equation. Then, the residual function becomes where Y N (t) is the approximate solution. Therefore, when Equation (44) is written in place of Equation (45) and the necessary arrangements are made, we get Then, it can be defined the error function as Now, by using Equations (1), (2), (44) and (46), we obtain the error differential equation as with the intial condition e N (0) = 0 or clearly, by using Equation (47), we obtain the error problem Then, when the problem (48) is solved as in the previous section 3.1, the estimated error function is gained in the form Here, a * n are the new Pell-Lucas coefficients that we will determine for n = 0, 1, . . . , N. Thus, it is determined the improved approximate solution in the form Consequently, the improved error function is

The case of interacting species model
In this section, the described method in the steps in Section 4 will be applied to the model (25). For this purpose, Equation (25) is shown as The approximate solutions Y i,N (t) are written in L[Y i (t)] = g i (t) since the obtained approximate solutions in Section 3.2 provide the equation. Then, the residual function can be represented in the format Here Y i,N (t) are the approximate solutions. Thus, whensoever the Equation (50) is written instead of Equation (51) and the necessary arrangements are made, we have Then, it can be written Now, by using Equations (25), (50) and (52), the error differential equation can be obtained as together with the initial conditions Here, a * n are the new Pell-Lucas coefficients that we will find for n = 0, 1, . . . , N. Thus, it is written the improved approximate solution in the format As a result, we can expresssed the improved error functions as follows

Applications
In this section, the described method in Section 3 will be applied. Next, the described method in Section 4 will be applied to indicate the reliability and sensitivity for procedure. Thus, it will be possible to comment on the accuracy and reliability of approximate solutions based on Pell-Lucas polynomials. In addition, the obtained numerical results are presented in tables and graphs, so the results can be commented more easily. All calculations obtained in this section have been calculated with the programs established in Matlab for convenience. Firstly, we consider the following logistic equation Here, r = 1, K = 1, λ = 2 and g(t) = 0 which the required for the method. The exact solution of the problem is Here, we take N = 3 and applying the method in Section 3.1. In such case, the Pell-Lucas approximate solution Y N (t) = Y 3 (t) can be written as Y 3 (t) = 3 n=0 a n Q n (t).
From here, we want to achieve the Pell-Lucas coefficients a n . So, if N is selected as 3, the set of collocation points is written as By using Equation (20) in the method in Section 3.1 and these collocation points, the fundamental matrix equation is obtained as By using Equation (22) in the method in Section 3.1 and the collocation points, the initial condition depended on the matrix of Pell-Lucas coefficients is written as U; λ = 2 0 2 0 ; 2 . As a result, it is substituted this matrix A into Equation (12) and hence when N is selected as 3, the Pell-Lucas solution becomes Since the exact solution of the equation is known, it can be obtained the actual absolute error function e 3 In addition, if the obtained approximate solution by Pell-Lucas polynomials Y 3 (t) is written in Equation (56) and subtracted from Equation (56), the error problem is obtained as follows is written explicitly, we get By applying the procedure in Section 4.1 for this problem when M is selected as 5, it is obtained the estimated error function e 3,5 (t). Thus, it is obtained the improved approximate solution Y i, 3,5 3,5 (t). Namely, better approximate solutions are obtained than the approximate solutions Y 3 (t). Finally, it is obtained the improved error function The exact solution and the approximate solutions for Example (5.1) at some points are given in Figure 1. According to these data, it can be observed that errors decrease as N values increase. Namely, the population of animal species exhibits a decline, but the rate of this decline decreases spontaneously.
We analyse the actual, estimated and improved errors of the Example 5.1 in Table 1 and Figure 2(b). Figure 3(a) compares the actual errors with the estimated errors, while Figure 3(b) compares the actual errors with the improved errors. According to these data, it can be said that the estimated errors are very close to the actual errors and the improved errors are better than the actual errors. In the case of improved errors, it can be said that as the M value rises, the errors diminish. Thus, improved approximate solutions were obtained which gave better results than the obtained approximate solutions.  represents the improved approximate solution and e i,N,M (t) represents the estimated error function.
As another example, we consider the following prey-predator models (3) Here, a = b = c = d = 1, λ 1 = 1.3, λ 2 = 0.6 and g i (t) = 0 (i = 1, 2) which the required for the method. Now, we take N = 3 and by applying the method in Section 3.2, so it can be obtained the Pell-Lucas approximate solution Y i, By using Equation (41) in the method in Section 3.       Since the exact solution of the equation is unknown, it can not be obtained the actual absolute error func- 3 (t). However, if the obtained approximate solution by Pell-Lucas polynomials Y i,3 (t) is written in Equation (62) and subtracted from Equation (62), the error problem is obtained as follows For this problem, when M is taken as 5 in the method according to Section 4.2, it is obtained the estimated error function e i, 3,5 (t). Thus, it is obtained the improved approximate solution Y i, 3,5 3,5 (t). Namely, better approximate solutions are obtained than the approximate solutions Y i,3 (t). Table 2 shows the values of the obtained estimated absolute errors for Example 5.2 at some points for different N and M values. In Tables 3 and 4 for ADM [34] forHPM [35] forBCM [36] f o rP M  [34] forHPM [35] forBCM [36] Figure 7(a,b). According to these data, it can be said that errors decrease as N and M values increase. Thus, residual correction is a extremely effective technique in the predator-prey model, as in the single species model. Finally, the results of the obtained estimated errors for different N and M values are displayed in Figure 6(a,b). According to these data, it can be said that errors decrease as N values increase.

Results and discussion
In this section, it is interpreted by making comparisons between the outcomes of the procedure in Section 5 and the outcomes of other procedures available in the literature. For this, the results are shown in the tables.
The comparisons with ADM [34], HPM [35] and BCM [36] were made in Tables 5 and 6 for the exact solution and the approximate solutions of the Example (5.1). For Table 8. Comparison of the approximate solutions with other methods at some points of Equation (62).
for ADM [34] [34] f  Table 9. Residual errors of the approximate solutions at some points of Equation (62).  [34] f Example (5.1) in Table 6, by using the exact solution and the obtained approximate solution, the actual errors are given by comparing the ADM [34], HPM [35] and BCM [36] for N = 3 and N = 5. According to all these data, it can be said that the obtained approximate solution by the present method is very close to the exact solution and gives better results than ADM [34] and HPM [35]. Similar results are obtained when the current method is compared to BCM [36]. The obtained approximate solutions for Example (5.2) in Tables 7 and 8 were compared for different N values by ADM [34], HPM [35] and BCM [36] methods.  [38] forGM [38] forGM [38] forGM [38]  In Table 9, residual errors of the obtained approximate solutions by the current procedure and residual errors of the obtained approximate solutions by ADM [34], HPM [35] and BCM [36] methods were compared for different N values. Accordingly, the obtained residual errors by the current procedure yielded better results than other procedure. In Table 10 the residual errors are compared by the GM method [38]. It can be said that although the current method and the GM method [38] give similar results at some points, the current method gives better results compared to the GM method [38] at some points. The advantage of the method is that it can be calculated in a short time with the help of Matlab. We show these CPU times (in seconds) in Tables 9  and 10.

Conclusions
In this study, numerical solutions of continuous population models for single and interactive species have been investigated. A collocation method based on Pell-Lucas polynomials for these numerical solutions has been established in Section 3. Next, a technique known as residual correction has been given in Section 4. According to this technique, better approximate solutions have been obtained using the approximate solutions given in Section 3. This technique is based on estimating  Tables 5 and 6. According to these results, the obtained approximate solutions from the discussed method are closer to the exact solution at more points than other methods. In Table 1 and Figure 2, it is examined the results obtained from actual, estimated and improved errors. From here, while the results of the estimated errors are very close to the results of the actual errors, the results of the improved errors give better results than the results of the actual errors. Here, it can be said that both the discussed method in Section 3 and the discussed method in Section 4 have been yielded good results. The results of the considered secondly sample for the prey-predator model have been compared with those of other methods available in the literature in Tables 7, 8 and 9. According to these results, the obtained approximate solutions from the discussed method are closer to the exact solution at more points than other methods. According to the obtained results from the residual errors and the residual errors for improved approximate solutions in Table 10, Figures 4 and 7, it can be said that the errors decrease as the N and M values rise.
As a result, it can be said that as the N and M values rise in the discussed method, the approximate solutions are closer to the exact solutions and give better results than the other methods available in the literature. The advantage of the procedure is that it yields well conclusions for models without an exact solution, as can be seen from Example (5.2). Another benefit of the procedure is that all results can be achieved in a short time. These obtained CPU times (in seconds) are shown in Tables 9 and 10. For the results of this article, the codes written in Matlab program have been used. Moreover, after the necessary arrangements have been made in the discussed method in this article, the method can be developed for different models.