Nonlinear Elimination of Drugs in One-Compartment Pharmacokinetic Models : Nonstandard Finite Difference Approach for Various Routes of Administration

The motivation for this study is to introduce and motivate the use of nonstandard finite difference (NSFD) schemes, capable of solving one-compartment pharmacokinetic models. These models are modeled by both linear and nonlinear ordinary differential equations. “Exact” finite difference schemes, which are a special NSFD, are provided for the linear models while we apply the NSFD rules, based on Mickens’ idea of transferring nonlinear models into discrete schemes. The method used was compared with other established methods to verify its efficiency and accuracy. One-compartment pharmacokinetic models are considered for different routes of administration: I.V. bolus injection, I.V. bolus infusion and extravascular administration.


Introduction
Pharmacokinetics modeling is the mathematical representation of the behaviour of a drug in the body or an area of the body, created to describe the pharmacologic or physiologic kinetics characteristics.Pharmacokinetics is the study of the basic processes that determine the duration and intensity of a drug effect within an organism.These models can assist in simulating the biological processes involved in the kinetic behaviour of a drug after it has been introduced into the body, leading to a better understanding of its dynamic effects.Mathematical modeling is currently a common tool used in the study of physiological and biochemical systems.It can be developed from non-compartmental representations to large scale multi-compartment models.One of the early uses of compartment models was reported by Widmark [1].He used a compartment model to describe the distribution of alcohol in the body.Such compartment models have proven to be a great advantage when screening drugs used by humans at any instant in time.
In the case of compartment models, mass-balance equations are used to represent each compartment.The number of compartments in the model depends on the rate of drug distribution to different parts of the body.Most studies use one-or two-compartment models.When a drug is eliminated, the drug concentration in the systemic circulation and in all tissues decline at the same rate because of the rapid distribution equilibrium.Drugs that follow this behaviour follow the one-compartment pharmacokinetic model, while in two-compartment models, the movement of the administered drug is distributed instantaneously to some tissues and slowly to other tissues.However, if the distribution of the drug happened at three different rates, a three-compartment model would be applicable.Our focus is on a one-compartment pharmacokinetic model, specifically aimed at different models of drug elimination.
After a drug is released from its dosage form, the drug is absorbed into the surrounding tissue and/or the body.As commented on by Shargel et al. [2], the distribution and elimination of the drug in the body varies for each patient but can be characterized using mathematical models and statistics.Being able to characterise drug distribution and elimination is an important prerequisite to be able to determine or modify the dosing regimens of individuals and groups of patients.Among the three main types of pharmacokinetic (PK) models: compartment, physiologic and non-compartmental models, compartmentally based models are known to be a very simple and useful tool in pharmacokinetic.In essence, a compartment model provides a simple way of grouping all the tissues into one or more compartments where drugs move to and from the central or plasma compartment.In this manner, we are able to model the transport processes between interconnected volumes, such as the movement of drugs and hormones in the human body.Compartment models assume that there is rapid and perfect mixing, so that the drug concentration remains the same in each compartment.The complex transport processes are approximated by assuming that the flow rates between the compartments are proportional to the concentration difference in the compartments.
Compartment models play a significant role in understanding the dynamics of drug concentration in the body.In practice, PK models seldom consider all the rate processes ongoing in the body, Shargel et al. [3].Due to the complexity of the models which incorporate such information, simplifying assumptions are often made so that solutions may be obtained.Traditional PK models, being simplified mathematical expressions, are based on the assumption of a linear relationship between the dose of a drug and its concentration-see Be ňová et al. [4].In a linear model, these rate coefficients called k are assumed to be constant.However, such assumptions regarding the linearity of the model do not necessarily describe the actual physical processes as accurately as a non-linear relationship may.In fact, the non-linearities seen in such models are related to drug absorption, distribution, metabolism and excretion and the pharmacokinetic of drug action.Since in most cases, these compartment models are described by autonomous linear or nonlinear ordinary differential equations, we choose to consider the latter in this work, focusing on three regimes of excretion.In the case when the nonlinearity or kinetics in the system are complex or when the number of compartments in the model becomes large, such as in the work of Sharma [5], exact solutions are not obtainable and hence we turn to numerical methods.Some of the well known standard numerical schemes produce unnecessary oscillations, introduce extraneous or spurious solutions, and converge to fixed-point solutions different from the corresponding derivative [6].Hence, one observes the occurrence of numerical instabilities.The nonstandard finite difference (NSFD) scheme was developed by Mickens as an alternative method providing an approximate solution to a wide range of differential equations and catering for the numerical instabilities that occur when using standard methods.NSFD methods have been well reported in recent years, mainly because they are efficient and and preserve qualitative properties, see for example Villatoro [7], Roeger [8], Ibijola & Obayomi [9], Manning & Margrave [10], Mickens [11][12][13][14][15][16] and Sunday [17] which give the relevant background materials on this topic.Some of these authors deal with the exact finite difference scheme which is a special NSFD method.
A standard finite difference (SFD) scheme is said to be "exact" for a particular differential equation if its local truncated error is exactly zero for its general solution (Gander & Meyer-Spasche [18] and Mickens [19]).In other words, when the analytical solution of a differential equation can be matched exactly with its corresponding SFD equation, then the "exact" finite difference solution exists.If exact general solutions of a differential equation are explicitly known, then an "exact" finite difference scheme exists.The idea of an "exact" finite difference scheme was first conceived by Mickens [19,20] who has shown that an exact explicit scheme is easily obtained from the knowledge of its analytical solution.Therefore, "exact" finite difference schemes are designed in such a fashion that the difference equation has the same general solution as the corresponding differential equation.In the situation where exact solutions exist, the solution can be re-structured in such a way to obtain the "exact" finite difference scheme.However, in the case where exact solutions are not possible, the rules proposed by Mickens [19] will be deployed.In such situations, the method is referred to as the NSFD method as discussed by Anguelov & Lubuma [21].The consequence of these rules is that while the scheme may not be "exact", qualitative properties of the corresponding differential equations for all step-sizes are preserved and thus elementary numerical instabilities that can arise are eliminated.
Our research is aimed at the well-known one-compartment model.We aim to investigate three forms of drug elimination from the body.In a one-compartment model, the body is assumed to be a single compartment and the drug absorbed achieves instantaneous distribution throughout the body, metabolizing between tissues.The drug output is characterized by an elimination rate.In this study, three dose regimen models are considered:
We consider these cases to illustrate the effectiveness of the NSFD method for the solution of nonlinear differential equations of this nature.The work conducted here is done with the aim of introducing a numerical method which may act as an effective tool to be employed in future research for the solution of models which are non-linear and/or describe multiple compartments.As such, we propose and illustrate the use of a numerical method of solution, namely the NSFD method, capable of efficiently obtaining solutions which are not only accurate but maintain the underlying dynamics of the system of equations.This choice of method impacts on whether we are able to consider non-compartment models; the NSFD method is not amenable to the simulation of non-compartment models as it provides a meta-analysis of the inter-compartment dynamics, whereas non-compartment models are unable to describe these meta-dynamics and instead conduct parameter estimation of the entire system as a whole through the use of experimental data.The advantage of the NSFD method is the ability to predict the concentration-time profile of a drug when there are alterations in the dosing regimen-this would not be possible were one to consider non-compartment analysis.Another advantage of the NSFD method is that it preserves significant properties of the analogous models and consequently gives reliable numerical results even when analytical solutions are not possible.The standard approaches to multi-compartment models assume linear dynamics over the duration of each time step, whereas the NSFD method assumes exponential dynamics.Hence, in the case of a linear model, the NSFD method recovers the model dynamics exactly.This paper illustrates the ability of the NSFD method to solve a one-compartment PK model with various modes of elimination, in a stable and robust fashion, with the ability to be extended to non-linear and/or multi-compartment models.
The variables of importance and their meaning are give below: The absorption rate constant for oral administration.

•
k el : Elimination rate of the drug leaving the central compartment.
The apparent volume of distribution.

Methods
While the implementation of the NSFD method is the focus of this research, we employ the Runge-Kutta as a means of comparison.This section provides an overview of NSFD and Runge-Kutta.

NSFD Modeling Fundamental Principles
NSFD methods provide numerical solutions to differential equations by constructing discrete models.They preserve the significant properties of their continuous analogues and consequently give reliable numerical results.The following rules were given by Mickens in [19] for constructing an NSFD scheme: Rule 1 The orders of the discrete representation of the derivative must be equal to the orders of the corresponding derivatives appearing in the differential equations.Rule 2 Denominator functions for the discrete representations for derivatives must, in general, be expressed in terms of more complicated functions of the step-sizes than those conventionally used.Rule 3 Nonlinear terms must, in general, be modeled by nonlocal discrete representations.Rule 4 All the special conditions that correspond to either the differential equation and/or its solutions should also correspond to the difference equation and/or its solutions.Rule 5 The discrete scheme should not introduce extraneous or spurious solutions.
Remark.Exact finite difference is a special NSFD.

Runge-Kutta Method
In a similar fashion with the finite difference scheme, we introduce the concept of the Runge-Kutta method from Taylor's theorem, where h is the step size between the values of the independent variable x.Consider Then, the Taylor's series expansion of Equation ( 1) is given by Differentiating Equation ( 1), we have x (t) is given in Equation ( 1), therefore Equation (3) becomes Substituting Equations ( 1) and (4) into Equation (2), we have With some manipulations, we have From Equation ( 6), if then classical second order Runge-Kutta method is given as The approximation given by Equation ( 9) has a local truncation error O(h 3 ).This second order Runge-Kutta method is also known as Heun's method.
The most widely used method is the fourth-order Runge-Kutta method which can be developed in a similar fashion to the second order Runge-Kutta.The local truncated error of the fourth-order Runge-Kutta method is O(h 5 ).Equation ( 1) can be solved using the classical fourth-order Runge-Kutta as follows: where The approximation given by Equation ( 10) subjected to k i , i = 1, 2, 3, 4 has an error of O(h 5 ) and thus is deemed the most accurate of all the approximations provided.

I.V. Bolus Injection
When drugs are administered by an I.V. bolus injection, the entire dose administered enters the bloodstream directly and is able to produce pharmacological effects.This is followed by the distribution of the drug through the circulatory system to all tissues in the body.Hence, we assume that a drug given by an I.V. bolus injection is rapidly mixed.Naturally, if you inject it directly into the bloodstream, the drug is immediately found in the bloodstream and does not have to be absorbed.The concentration at time t 0 = 0 corresponds to the dose given in this manner and hence, it describes the I.V. bolus injection route of administration.We consider this mode of administration alongside two different elimination processes: (1) drugs eliminated by linear pharmacokinetic and (2) drugs eliminated by nonlinear processes.We solve the differential equation arising from these elimination processes, employing the NSFD method as a means of comparison with the SFD method.
The case where the drugs are given via an I.V. bolus injection, distributed as a two-compartment model and then eliminated only by linear pharmacokinetic, is presented in Egbelowo et al. [22].In the manuscript, we did not provide the results in the case of a one-compartment model, which is of interest here.When considering the SFD method for the one-compartment model that describes the distribution and elimination after an IV bolus dose, i.e., we notice the following interesting dynamics: In turn, the corresponding NSFD scheme constructed gives accurate results without the requirement (16) as needed by the SFD scheme.

I.V. Bolus Injection: Nonlinear Pharmacokinetic Elimination
The equation that describes the elimination of a drug that is distributed in the body as a one-compartment and is eliminated by nonlinear pharmacokinetic after an I.V. bolus injection is given as per Shargel et al. [23].From the compartment diagram in Figure 1, we have the differential equation subject to where V max is the maximum elimination rate and K m is the Michaelis constant.Michaelis-Menten kinetics are also referred to as the capacity-limited metabolism, saturable metabolism, or mixed-order kinetics as discussed in Be ňová et al. [4].
V max , K m C The solution given as C(t) to the differential Equation ( 17) has no closed form expression.The explicit closed-form solution of the one-compartment I.V. bolus injection model that follows Michaelis-Menten kinetics given in Equation ( 17) is, however, expressible in terms of the Lambert W-function as Definition 1.The Lambert W-function is defined to be a multivalued (single valued in the case of PK applications) inverse of the function x → xe x satisfying [24] such that x = ye y .( 21) Equation ( 21) can then be written as where W is Lambert's function.Equation (19) satisfies the transcendental equation (22).The exact finite difference scheme of ( 19) is given as A SFD scheme for the Michaelis-Menten Equation ( 17) might take the form We provide NSFD schemes for Equation (17) for the case where Equation ( 17) is discretized using three different discretization methods, as done in the works of Mickens [25] and Chapwanya et al. [26,27] namely, semi-implicit forward-Euler, implicit forward-Euler and explicit forward-Euler.

Case 1: Semi-Implicit Forward-Euler
The NSFD scheme here is obtained from the semi-implicit forward-Euler discretization as This can be rewritten as We then use the fact that allowing us to make the replacements which implies Therefore, the denominator function for the semi-implicit forward-Euler discretization is given as The NSFD semi-implicit discretization of Equation ( 17) is The NSFD scheme (31) is compared with the corresponding SFD scheme (25).The scheme given by Equation (31) has the same qualitative behaviour as the original differential equation for all step sizes.

Case 2: Implicit Forward-Euler
When Equation ( 17) is discretized using an implicit discretization, we obtain and following the same process as described by Equations ( 25)-(29), we obtain the same denominator function given in Equation (30).Therefore, the NSFD scheme for Equation (17), when discretized using an implicit forward-Euler approximation, is

Case 3: Explicit Forward-Euler
The NSFD scheme, in this case, is obtained from an explicit forward-Euler discretization as The denominator function is given by Similarly, the NSFD discretization of Equation ( 17) is From the three cases, we have that the denominator function obtained depends on the discretization method used.Upon employing the Lambert function, the resultant scheme for Equation (17) as exhibited by Equation ( 23), is This scheme will serve as a means of comparison with the results obtained from the three cases above.

I.V. Bolus Injection: Mixed Drug Elimination
Another possible route of drug elimination is mixed drug elimination.In this elimination process, as depicted in Figure 2, drugs are eliminated by nonlinear processes.Therefore, the equation that best describes a drug that is eliminated by Michaelis-Menten kinetics after an I.V. bolus injection is given by where k el is the first-order rate constant representing the sum of all first-order elimination processes.
The second term of the Equation (38) represents the saturable process.The SFD scheme of Equation ( 38) is given by Implementing the NSFD scheme as before, we obtain a denominator function given by which provides the following NSFD scheme for Equation (38)

I.V. Bolus Infusion
I.V. bolus infusion is the process of infusing a drug at a constant rate.The drug input is constant and equal to the rate of infusion of the drug.On starting the infusion, there is no drug in the body and therefore no elimination.The concentration of the drug in the body then rises, but as the drug concentration increases, so does the rate of elimination.Thus, the rate of elimination will keep rising until it matches the rate of infusion.The concentration of the drug in the body is then constant and is said to have reached a steady state.A similar approach is used for the I.V bolus injection process, and we consider this mode of administration alongside two different elimination processes as before: (1) drugs eliminated by linear pharmacokinetic and (2) drugs eliminated by nonlinear processes.
In the case when the drug is given via I.V. bolus infusion and the drug is excreted in a linear way, we observe similar dynamics as we did for the case where an I.V. bolus injection is the means of administration (see Egbelowo et al. [22]) where both cases was considered for a two-compartment model.The one-compartment model that describes the distribution and elimination after an IV infusion dose is given by The solution obtained via the SFD scheme gives the following results: with an oscillating amplitude via an alternating sign at each step, (iv) if hk el = 2, C k oscillates with a constant amplitude 2R k el , and (v) if hk el > 2, C k oscillates with an increasing amplitude.
h is the step-size, k el is the elimination rate of the drug, C k represents the concentration of the drug at time is the flow rate of the drug, and R 1 is the infusion rate per unit time.From these results, we conclude that this model will have numerical instabilities for all cases except for cases (i) and (ii).Maintaining the requirements given by cases (i) and (ii), the same qualitative behaviour is observed as the original differential equation.The NSFD scheme constructed in turn gave accurate results for all the cases given above.

I.V. Bolus Infusion: Nonlinear Pharmacokinetic Elimination
In Figure 3, the drug is administered by constant infusion and is eliminated by nonlinear pharmacokinetic processes.The equation that describes the rate of change of the plasma concentration, as depicted in Figure 3, is given by All the parameters in Equation ( 43) are defined as for model (17).Solving Equation (43) using the concept of the W-Lambert function, we have where Equation ( 45) can be written in the form where L 1 = L(h) and M 1 = M(h).The steady-state concentration of Equation ( 43) is determined by the following equation The NSFD scheme for Equation ( 43) is structured as where φ 1 (h) is defined as before.Comparing Equation (49) with the SFD scheme given by shows that the NSFD scheme is dynamically consistent with the original differential equation for any step size.

C R 1
V max , K m

I.V. Bolus Infusion: Mixed Drug Elimination
Figure 4 describes the rate of change in the plasma drug concentration for a drug that is given by I.V. infusion and eliminated by nonlinear pharmacokinetic.This is an extension of Figure 3, which leads to Equation ( 51) All the parameters are defined as done for the model given by Equation (17).The steady-state concentration of Equation ( 51) can be determined by The NSFD scheme of ( 51) is given by where φ 3 (h) is given in Equation (40).The SFD scheme of the relevant equation is

Extravasular Administration
A drug administered via the extravascular route of administration undergoes the process of absorption before it gets to the systemic circulation.This type of drug delivery is complicated by the variable at the site of absorption.The level of absorption of a drug from the gastrointestinal tract (GIT) depends on the anatomy and physiology of the absorption site, physiochemical properties of the drug, and physiochemical properties of the dosage form.Initially, the entire drug is in the site of absorption and none has yet reached the systemic circulation [28].Most drugs administered extravascularly act systemically.In such cases, systemic absorption is a prerequisite for efficacy.This section describes the extravascular route of administration.We consider this mode of administration via two different elimination processes: (1) drugs eliminated by linear pharmacokinetic (2) drugs eliminated by nonlinear processes.

Extravasular Administration: Linear Pharmacokinetic Elimination
When a drug is administered through extravascular administration and eliminated by a linear process as shown in Figure 5, we apply Equations ( 55) and ( 56) in order to model the process.Equation (55) describes the drug at the site of absorption before it reaches the systemic circulation, while Equation (56) describes the concentration of the drug at the systemic circulation.D is the amount of drug in the GIT at any time t, k a is the first-order absorption rate constant, C is the plasma concentration of the drug in the body and k el is the elimination rate.Thus, the disappearance rate of the drug from the GIT is given by (also termed the equation for drug in GIT), The rate of change of the amount of drug in the body is given by where C is the amount of the drug available at the absorption site.The time course of the the amount of drug that follows the oral route of administration is given by where D 0 is the dose of the administered drug, k a is the constant of the absorption, and k el is the rate of elimination.At t = ∞ (at later time intervals) the above equation reduces to (i.e., when e −k a t ≈ 0) The "exact" finite difference scheme of the model is derived from the analytical solution.Since Equations ( 55) and ( 56) can be solved simultaneously, we proceed as The particular solution of Equation ( 59) is Considering the system of Equations (59) in matrix form, the corresponding matrix is The matrix has eigenvalues λ if det(M − λI) = 0 or λ 2 − tra(M)λ + det(M) = 0. Therefore, the eigenvalue equation to be which provide the eigenvalues represents the eigenvectors corresponding to the eigenvalues, then (M − λI)v = 0. Hence, we have that giving the general solution of the system as To calculate A and B, we use the initial values given in Equation (60) so that Equation (66) simplifies to From Equation (65), we have that Substituting Equation (68) into Equation (69) and after some algebraic manipulations, we obtain Substituting Equation ( 70) into (68), we obtain Substituting Equations ( 70) and (71) into Equation (66) with some manipulations we obtain where p = k a λ 2 −λ 1 .The 'exact' finite difference scheme of Equation ( 59) is obtained by making the following transformations in Equation (72) Thus, giving the "exact" finite difference scheme where The "exact" finite difference result obtained in system (76) will be compared with the SFD scheme k el C k a D

Extravascular Administration: Mixed Drug Elimination
In the situation when the drug is administered by the extravascular mode of administration and eliminated by parallel pathways, Equation ( 79) is applied to describe Figure 6.Consider with initial conditions D 0 = D(t 0 ) and C 0 = C(t 0 ).The NSFD scheme of Equation ( 79) is and may be compared to the SFD scheme of the form

Numerical Simulations and Discussion
One-compartment models with different routes of administration (I.V. bolus injection, I.V. bolus infusion and extravasular) are considered for simulations.In order to perform a useful comparison, these methods were tested under similar conditions corresponding to the intended practical application.The purpose of the tests was to compare the accuracy and stability of the various numerical schemes employed.This is done for varying step-sizes and the results are examined in the figures and numerical results presented in the next few sections.The numerical calculations are carried out in MATHEMATICA, and the results are then processed in MATLAB to generate visual representations.Thi section presents the results of the case when the drug is administered by I.V. bolus injection and eliminated by Michaelis-Menten elimination.NSFD schemes (31), ( 33) and (36), and the SFD schemes ( 25), ( 32) and (34) respectively, are compared with the analytical solution (19) in Figure 7.The analytical solution was obtained through the use of the W-Lambert function.19) and the corresponding SFD scheme (32) and (c) NSFD scheme (36) in case 3 plotted against the analytical solution ( 19) and the corresponding SFD scheme (35).

Results Describing Mixed Drug Elimination
In this section, simulations of the equation that describe a drug that is eliminated by mixed drug elimination after an I.V. bolus injection are provided.Figure 8 shows a comparison between the SFD scheme in Equation (39), the NSFD scheme in Equation ( 41) and MATLAB built-in function ODE45.

Results Describing Nonlinear Pharmacokinetic Elimination
The results for the case when the drug is administered by I.V. bolus infusion and eliminated by Michaelis-Menten elimination is presented.The NSFD scheme (49) and SFD scheme (50) are compared with MATLAB built-in, ODE45.From Figure 9 we see that regardless of the step-size, the NSFD scheme (49) converges to the steady state.Table 1 gives the simulation results of the I.V. bolus infusion case where nonlinear pharmacokinetic elimination is present.

Results Describing Mixed Drug Elimination
Figure 10 show the simulation results of the NSFD scheme (53) and SFD scheme (54) in comparison to the results obtained via the in-built function ODE45.

Results Describing Linear Pharmacokinetic Elimination
Results of the one compartment pharmacokinetic model administered by an extravascular mode of administration and following linear elimination are presented here.The NSFD scheme (76) is compared to the SFD scheme (78).The schemes obtained from the model using different methods are tested under similar conditions.Simulations are provided for h = 0.5 and h = 1 in Figure 11.Table 2 shows numerical results for the "exact" finite difference scheme for extravascular administration in comparison with standard methods (Euler, Heun and Runge-Kutta) and the analytical solution of the model.The numerical results for the 'exact' finite difference scheme are the same as the analytical solution for any value of t. Figure 12 shows the simulation results for the model describing extravascular administration along with mixed drug elimination processes.Comparison of methods for one-compartment extravascular administration that follows mixed drug elimination i.e., NSFD scheme (80) and SFD scheme (81) is compared with ODE45.

Conclusions
In this work, we presented one-compartment pharmacokinetic models with different routes of administration.We presented numerical results via a variety of schemes for each of the developed models, paying attention particularly to the efficiency of the NSFD method in comparison to standard methods.From the results obtained, we observe that the stability of the NSFD scheme is independent of the chosen step-size for the linear cases.This is not the case with standard methods such as Euler and Heun methods.With the later methods, the step-size must be chosen in a reasonable domain, otherwise numerical instabilities will occur.The numerical simulations conducted verify that NSFD schemes are efficient and accurate for the solution of problems modelling pharmacokinetic processes.Importantly, as pointed out through test cases in this work, the NSFD method is able to generate numerical schemes that are dynamically consistent with the original differential equations.

Figure 2 .
Figure 2. Schematic representation of I.V. bolus injection with mixed drug elimination.

Figure 4 .
Figure 4. Schematic representation of I.V. infusion with mixed drug elimination.

Figure 5 .
Figure 5. One-compartment pharmacokinetic model for first-order drug absorption and first-order elimination.

Figure 6 .
Figure 6.Schematic presentation of extravascular administration with both linear and Michealis-Menten elimination.

Figure 10 .
Figure 10.NSFD scheme (53) and SFD scheme (54) is comparison with ODE45.C ss2 is the steady state Equation (52).The concentration of the drug when administered via I.V. bolus infusion and eliminated by mixed drug processes for (a) h = 0.75 and (b) h = 3.3333.

Figure 12 .
Figure 12.Comparison of methods for one-compartment extravascular administration that follows mixed drug elimination i.e., NSFD scheme (80) and SFD scheme (81) is compared with ODE45.
C k tends to 0 with an oscillating amplitude via an alternating sign at each step, (iv) if hk el = 2, C k oscillates with a constant amplitude C 0 , and (v) if hk el > 2, C k oscillates with an increasing amplitude, where h is the step-size, k el is the elimination rate of the drug, and C k represents the concentration of the drug at time t k .Results indicate that the SFD scheme developed has the same qualitative behaviour as the analytical solution of the pharmacokinetic model if 0 < hk el 1.

Table 1 .
The absolute error results of Equations (43) for C with parameters values R = 0.5, K m = 4, and V max = 2.

Table 2 .
The numerical results for the extravascular administration model.