Approximate analytical solutions for the blood ethanol concentration system and predator-prey equations by using variational iteration method

1.   Introduction

Over the past three decades, numerous authors have continued to be interested in fractional calculus ^{[1,2,3,4,5,6]}. To fulfill the demand to model real-world problems in various domains like mechanics, biology, finance and engineering ^{[7,8,9,10,11,12]}, several academics have discovered that developing new fractional derivatives with various single or non-singular kernels is crucial. There were very important studies done with the help of the Caputo operator, for example, periodic motion for piecewise Caputo derivatives of impulsive fractional functional differential equations was studied in ^[13] and the Mittag-Leffler stabilization and exponentially stable periodic oscillation for factional-order impulsive control neural networks using piecewise Caputo derivatives also was studied in ^[14]. At the same time, there are many efforts to develop the theory of fractional calculus, for example on the foundation of the conventional Caputo fractional derivatives, Caputo and Fabrizio constructed the Caputo-Fabrizio operator (CFO). This new definition has a nonsingular kernel, which sets it apart in the Caputo sense (i.e., exponential function) ^[15,16]. See for further information on the definitions of fractional derivatives and their characteristics ^[17].

The majority of the fractional differential equations (FDEs) lack an exact solution, so numerical and approximate techniques must be used ^{[18,19,20,21,22,23,24,25]}. One of the most popular of these methods is, the variational iteration method (VIM), which is a potent analytical technique, and was initially described in ^[26]. Numerous circumstances have successfully used this method: for instance, see ^[27]. In addition, the first instance of the VIM being used with FDEs was in ^[7]. The VIM was used by Odibat and Momani ^[29] to solve partial FDEs. The advantages of the VIM lie in it treating the shortcoming in some other approximate and numerical methods, such as the Adomian decomposition method ^[30], where it does not need to compute the Adomian polynomials of the nonlinear terms in the differential equation under study and does not need to divide the domain of the problem to compute the numerical solution only at the resulting nodes as in the finite difference method.

In this work, we will use the VIM to investigate the approximate solutions for two important models. The first one is the fractional BECS, whereas the second model is the fractional PPEs. The following is how the paper is set up:

● Some of the preliminaries and the studied models are presented in Section 2.

● The solution procedure is shown in Section 3.

● Some numerical examples for the two studied models are presented in Section 4.

● The paper's conclusion is presented in Section 5.

2.   Preliminaries, notations and formulations of the models

2.1.   Some definitions of fractional calculus

Definition 2.1. The fractional integration in the Riemann-Liouville sense is defined by ^[31]:

where $\Gamma(.)$ is the Gamma function.

Definition 2.2. The fractional derivative $D^{\nu}$ of order $0 < \nu\leq1$ in the Caputo sense for $\varphi(t)\in H^{1}(0, b)$ is given by:

2.2.   The proposed fractional models

2.2.1.   Model 1: blood ethanol concentration system

This part focuses on figuring out the levels of alcohol in a human body's stomach $\Phi(t)$ and blood $\Psi(t)$. An experimental study was carried out in ^[32], the major source of the actual data for the current work. The suggested model is based on the fractional derivative kinetic reaction (for $0 < \nu\leq1$) and provided by:

where we have the following descriptions for the included functions and parameters ^[32]:

$\Phi(t)$: The alcohol's concentration in the stomach at time $t$ (mg/l).

$\Psi(t)$: The alcohol's concentration in the blood at time $t$ (mg/l).

$\beta$: The rate law constant 1 (min$^{-1}$).

$\mu$: The rate law constant 2 (min$^{-1}$).

The exact solution of the above system is ^[33]:

2.2.2.   Model 2: predator-prey system

The predator and prey equations are a pair of nonlinear first-order differential equations that are frequently used to analyze the dynamics of biological systems involving interactions between two species, one of which is a predator and the other a prey. Samardzija ^[34] made an extension of this model and proposed the concept of predators and single prey for the Lotka-Volterra system. This system is what we are going to study here, and it is formulated in its fractional form as follows:

where $c_{1}$–$c_{3}$ are constants. Here, the number of predators are $\Phi(t)$ and $\Psi(t)$, and the number of its prey is $\Upsilon(t)$, through the time $t$, $\sigma_{1}$–$\sigma_{7}$ are parameters that elucidate the interaction between the three species ^[34].

3.   Implementation of variational iteration method

The fundamentals of the VIM and how it applies to different types of differential equations are described in ^[35]. It was demonstrated in ^[7] that the FDEs may also be solved using the VIM. In this section, we expand the application of the VIM to solve the two proposed fractional system models.

3.1.   Concepts of VIM via BECS

We can create the correction functionals form of the systems (2.1) and (2.2) according to the VIM as follows:

where $\lambda_{k}, \, (k = 1, 2)$ are the general Lagrange multipliers, which can be identified optimally via variational theory ^[36,37]. Some approximations must be made to determine the estimated Lagrange multipliers. It is possible to roughly express the rectification of functional Eqs (3.1) and (3.2) as follows:

where $\tilde{\Phi}_{m}$ and $\tilde{\Psi}_{m}$ are considered as restricted variations, in which $\delta\, \tilde{\Phi}_{m} = \delta\, \tilde{\Psi}_{m} = 0$. To find the optimal $\lambda_{1}$ and $\lambda_{2}$, we proceed as follows:

The following can be done to acquire the stationary conditions:

The solutions of the Eq (3.7) is given by:

We substitute from (3.8) into the functional Eqs (3.1) and (3.2) to obtain the following iteration formula:

The beginning estimations $\Phi_{0}(t)$ and $\Psi_{0}(t)$ can be freely chosen if they satisfy the initial conditions of the problem. Finally, the solutions $\Phi(t)$ and $\Psi(t)$ can be approximated by the $m$-th terms $\Phi_{m}(t)$ and $\Psi_{m}(t)$, respectively as follows:

3.2.   Concepts of VIM via PPS

By following the same procedure in subsection 3.1 but on the systems (2.5)–(2.7), we can obtain the following iteration formula:

The initial approximations $\Phi_{0}(t), \, \Psi_{0}(t)$ and $\Upsilon_{0}(t)$ can be freely chosen if they satisfy the initial conditions of the problem. Finally, we approximate the solutions $\Phi(t), \, \Psi(t)$ and $\Upsilon(t)$ by the $m$-th terms $\Phi_{m}(t), \, \Psi_{m}(t)$ and $\Upsilon_{m}(t)$, respectively as follows:

4.   Numerical applications

4.1.   Approximate solution of the BECS

Now, we are ready to get the approximate solution of the studied model BECS by considering the variational iteration formulas (3.9) and (3.10) with distinct values of $\nu, \ m$ with $\beta = 0.02873$, $\mu = \; 0.08442$ and initials conditions $\Phi_{0} = 4$ and $\Psi_{0} = 0$ through Figures 1–6.

If we start with $\Phi_{0}(t) = 4$ and $\Psi_{0}(t) = 0$ in the iteration formulas (3.9) and (3.10) we can obtain directly some of the other components of the solution as follows:

and so on. The remaining parts of the iteration formula can be produced in a similar manner. Here, in our computation, we approximated the solution to $\Phi(t)$ and $\Psi(t)$ by $\Phi(t)\cong\Phi_{m}(t)$ and $\Psi(t)\cong\Psi_{m}(t)$, respectively.

In Figure 1, both the exact and numerical solutions with $\nu = 0.97$ and $m = 4$ are compared, whereas in Figure 2, the absolute error with $\, \nu = 0.93$ and $m = 4$ is presented. In Figure 3, again the exact and numerical solutions are compared but for $\nu = 0.87$ and $m = 4$, whereas in Figure 4, the absolute error is presented with $\nu = 0.83$ and $m = 4$.

Figure 5 is given to see both the behavior of numerical (a and c) and exact solutions (b and d) with distinct values of $\nu = (0.92, \, 0.82, \, 0.72, \, 0.62)$ at $m = 5$. Finally, Figure 6 is given to show both the behavior of numerical (a and c) and exact solutions (b and d) with distinct values of $\beta = 0.05, \, 0.10, \, 0.15, \, 0.20$ and $\mu = 0.15, \, 0.20, \, 0.25, \, 0.30$ with $\nu = 0.95$ and $m = 5$.

As a result, the behavior of the numerical solution is strongly dependent on $\nu$, $\beta$ and $\mu$ which confirms that, in the case of fractional derivatives, the suggested approximation method is effectively used to solve the stated problem.

To validate the numerical solutions at ($\nu = 0.9$, $\Phi_{0} = 6$ and $\Psi_{0} = 0$), $m = 5$ and the same values of the parameters as in Figure 1, a comparison of the absolute error (AE) is presented in Table 1 for the proposed method and the Chebyshev spectral collocation method for the same model but with the Atangana-Baleanu-Caputo fractional derivative with non-singular kernel ^[38]. This comparison demonstrates how the approach suggested in this article is suitable for solving the proposed model with its fractional form.

4.2.   Approximate solution of the PPS

Here, we verify the accuracy of the presented scheme by presenting a numerical simulation on a test example, where we address the systems (2.5)–(2.8) with different values of $\nu$ and $m$ with $\sigma_{1} = \sigma_{2} = \sigma_{3} = \sigma_{4} = 1, \, \sigma_{5} = 2, \, \sigma_{6} = 3$ and $\sigma_{7} = 2.7$ and initial conditions $\Phi_{0} = 1, \, \Psi_{0} = 1.4$ and $\Upsilon_{0} = 1.$ The obtained numerical results for the studied model by applying the proposed technique are shown in Figures 7–9.

If we start with $\Phi_{0}(t) = 1, \, \Psi_{0}(t) = 1.4$ and $\Upsilon_{0}(t) = 1$ in the iteration formulas (3.12)–(3.14), we can obtain directly some of the other components of the solution as follows:

and so on. The remaining parts of the iteration formula can be produced in a similar manner. Here, in our computation, we approximated the solution to $\Phi(t), \, \Psi(t)$ and $\Upsilon(t)$ by $\Phi(t)\cong\Phi_{7}(t), \, \Psi(t)\cong\Psi_{7}(t)$ and $\Upsilon(t)\cong\Upsilon_{7}(t)$, respectively.

The behavior of the approximate solution for $\nu = 1.0, \, 0.97, \, 0.94$ with $m = 7$ is given in Figure 7. In Figure 8, we present the approximate solution for $m = 6, \, 7, \, 8$ with $\nu = 0.95$. Figure 9 is plotted to represent the residual error function (REF) ^[39] of the approximate solution at $\nu = 0.96$ with $m = 6$. From these results, it can be noted that the behavior of the numerical solution resulting from the application of the proposed method depends on $\nu$ and $m$, and this demonstrates that the suggested approach is appropriate for solving the suggested model in its fractional version.

As opposed to their counterparts, the integer order derivative operators, fractional models have been found to be more accurate at estimating the true data. The fractional operator has the fewest errors, according to the REF. These instances involve the exact same mistakes. On the other hand, under the optimum values of the computed parameters, the errors for the integer order and the Caputo operator are the same.

5.   Conclusions

The powerful variational iteration method has been used to approximate the solutions of two important nonlinear systems of the FDEs, namely the predator-prey equations and the blood ethanol concentration system. The findings demonstrate that the used technique is efficient and cost-effective for obtaining rough solutions of the two suggested models. By including new components of the solution drawn from the solution sequence, we can control and reduce the absolute approximate error. Also, to demonstrate the viability of the suggested method, approximate solutions with various values of the fractional-order $\nu$, the order of approximation $m$ and the residual error function were computed. It has been shown that the Caputo fractional operator is substantially more accurate than the integer order version of the model for estimating the amount of alcohol in a human's blood. Also, it has been noted that the performance of this fractional operator is on par with that of its integer order model. We intend to deal with these models in the future, but on a larger scale, by generalizing this research to include a modified proposed method or additional types of fractional derivatives. The Mathematica software program was used to perform numerical simulation operations.

Use of AI tools declaration

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

Acknowledgments

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

Conflict of interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.