Optimized Direct Padé and HPM for Solving Equation of Oxygen Diffusion in a Spherical Cell

1National Institute for Astrophysics, Optics and Electronics, Luis Enrique Erro No. 1, Sta. Maŕıa Tonantzintla, 72840, Puebla, Mexico 2Facultad de Instrumentación Electrónica, Universidad Veracruzana, Cto. Gonzalo Aguirre Beltrán S/N, 91000, Xalapa, Veracruz, Mexico 3Centro de Investigación de Micro y Nanotecnologı́a, Universidad Veracruzana, Calzada Ruı́z Cortines 455, Boca del Rı́o, 94292, Veracruz, Mexico 4Maestŕıa en Ingenieŕıa Aplicada, Facultad de Ingenieŕıa de la Construcción y el Hábitat, Universidad Veracruzana, Calzada Ruı́z Cortines 455, 94294, Boca del Rı́o, Veracruz, Mexico 5Instituto Tecnológico de Celaya, Antonio Garcı́a Cubas Pte No. 600 Esq. Av. Tecnológico, Celaya, 38010, Guanajuato, Mexico


Introduction
Michaelis-Menten kinetics describes the rate of enzymatic reactions. This model is valid when the concentration of a substrate is higher than the concentration of the enzymes, and for steady state conditions, that is, when the concentration of the complex enzyme-substrate is constant [1]. There are several works about the Michaelis-Menten oxygen uptake kinetics. For example, in [2], the relation between Michaelis-Menten direct and inverse kinetics, chemical kinetics of approximation, and second-order kinetics is presented. Moreover, in [3], the subject of transport phenomena is approached by means of the study of transport of quantity of motion (viscous flow), transport of energy (heat conduction, convection, and radiation), and matter transport (Diffusion) [4,5]. Consequently, transport phenomena are employed for solving different problems in the area of sciences such as chemistry, biochemistry, soil physics, meteorology, biology, and semiconductors, disciplines in which the use of Bessel functions, differential equations, and Laplace transform is required [6][7][8][9][10]. In [11], the behaviour of dopamine released from a small iontophoresis electrode and its voltammetric detection by a carbon fiber sensor 100pm away is presented as a basis for developing a new paradigm for measuring dopamine kinetics in intact rat neostriatum. In [11], the presented model was derived from the work of diffusion in a punctual iontophoretic source where the lineal term of absorption is replaced by a nonlinear expression that describes a Michaelis-Menten kinetic given by a constant , the Michaelis-Menten constant . In [12], the equation that models oxygen diffusion in a spherical cell, including nonlinear absorption kinetics, is solved by transforming the Lane-Emden equation into its equivalent Volterra integral form, then Adomian decomposition was employed to solve the nonlinear Fredholm-Volterra integral. Furthermore, in [13] a fully symbolic solution was proposed to solve the aforementioned equation by means of a modified Taylor series obtaining an accuracy of 1.66×10 −15 as the lowest mean square error without the use of complicated integrals.
In [37] was presented a procedure to apply Padé method to find approximate solutions for nonlinear differential equations, which consist in that the solution of a differential equation can be directly expressed as a rational power series of the independent variable as a Padé approximant. From (6) ODP employs a polynomial-like rational expression as the proposal of approximation of the nonlinear differential equation to be solved. In general terms, it works by means of substituting the rational expression in the differential equation and then regroups the powers of the independent variable. It is important to note that due to the rational expression a large amount of algebraic operations is generated. However, this work proposes a Direct Padé (DP) modification oriented to reduce this algebraic operation by means of the Taylor expansion of the rational function.
This paper is organized as follows. In Section 2, a brief description of Michaelis-Menten kinetics and its equation is given. In Section 3, a brief description of Lane-Emden equation is presented. The equation to be solved and its boundary conditions are introduced in Section 4. Section 5 describes in detail Optimized Direct Padé. In Section 6, the HPM solution is presented. Section 7 presents the analysis of solution for the case study by HPM and ODP. In Section 8, its numerical simulations, comparisons, and discussion are presented. Lastly, in Section 9 a brief conclusion of this work is presented.

Michaelis-Menten Equation
General principles of kinetics in chemical reactions are applicable to catalysed reactions by enzymes in living things. However, this shows a characteristic side which is not observed in nonenzymatic catalyst, the substrate saturation, in terms of enzyme molecules active sites occupation [1]. The study of the effect that substrate concentration has over enzyme activity is not simple task; logically thinking, substrate concentration lowers as the reaction increases. A simplification on kinetics experiments consists in measuring initial velocity . If time is short enough, substrate lowering will be minimal and this could be considered, thus, constant. This behaviour is characteristic of most enzymes and was studied by Michaelis and Menten in 1913 [38]. Figure 1 presents the three phases of enzyme kinetics: (i) For a low substrate concentration, velocity reaction is directly proportional to substrate concentration (linear relation), first-order kinetics. (ii) For a high substrate concentration, velocity reaction is practically constant and independent of substrate concentration; kinetics is considered zero order. (iii) For medium concentrations of substrate, velocity of the process becomes nonlinear and this phase is called mixed kinetics.
The curve which expresses the relationship between substrate concentration and initial velocity has the same shape for most enzymes, a rectangular hyperbola, whose algebraic expression is given by Michaelis-Menten equation (1), which describes the relationship between initial velocity and substrate concentration, [1]: , maximum velocity, and , Michaelis constant, are two characteristic kinetics parameters of every enzyme that can be experimentally determined.
Maximum velocity is obtained when velocity of reaction becomes independent of substrate concentration; when this occurs, velocity reaches a maximum value. This value depends on the available amount of enzyme. Michaelis constant indicates substrate concentration in which reaction rate is half of maximum velocity. This parameter is independent of enzyme concentration and it is characteristic of every enzyme according to the employed substrate (or substrates). Constant indicates the affinity that the enzyme has for the substrate, being this higher at a lower value of . For a lower value , the lower the substrate needed to reach half of maximum velocity, thus, the higher the affinity an enzyme has for the substrate.
Discrete Dynamics in Nature and Society 3

Lane-Emden Equation
Lane-Emden equation (2) is a basic nonlinear equation arising in the study of stellar structure theory and has where can be chosen to describe a Cartesian geometry, = 0; cylindrical, = 1; or spherical, = 2 [39]. ( , ) is a real function and ( ) is an analytical function. Initial conditions are where and are constants.
Equation (2) is used to model physical and astrophysical phenomena such as theory of stellar structure, spherical cloud of gas thermal behaviour, isothermal gas sphere [40], and transport phenomena [5] as it occurs for the case of particle diffusion in biophysical processes [4] as well in physiological processes.
A model of oxygen diffusion and nonlinear absorption in a sphere was originally proposed and resolved in [41]. Later, the same model was proposed, analysed, and solved in [42]. The equation of the model is a steady state diffusion-reaction equation that represents oxygen transport by means of linear diffusion in a sphere.

Equation of Oxygen Diffusion in a Spherical Cell
The differential equation that models the oxygen diffusion in a spherical cell with Michaelis-Menten oxygen uptake kinetics (4) [41][42][43][44][45][46] is obtained by substituting (1) into (2) and given by with boundary conditions Variables and represent the oxygen concentration and the radial distance, respectively. Parameters , , and are constants that represent the maximum reaction rate. Michaelis constant is the concentration at half saturation and cell membrane permeability, respectively.
When boundary condition is equal to 0, oxygen distribution is symmetrical to the center of the sphere. Boundary condition at = 1 determines that oxygen flow through the cell membrane is proportional to (1) = 1 − (1) differential, which is less than the cell membrane normalized concentration [44,45].

Optimized Direct Padé
The method ODP [37] considers that a nonlinear differential equation can be expressed as follows: with boundary conditions where and are the linear and nonlinear operators, respectively; ( ) is a known analytical function; is a boundary operator; Γ is the boundary of the domain Ω [47]. 0 is the initial condition for (6) that satisfies the boundary condition and is the homotopy perturbation parameter. The solution of (6) can be written as where 0 , 1 , ⋅ ⋅ ⋅ and 0 , 1 , ⋅ ⋅ ⋅ are unknowns determined by the method; , are the numerator and denominator order.
There are no systematic criteria to choose the optimal Padé order [ / ] for a given problem. However, generally a finite number of terms are required to obtain highly precise Padé approximation. The basic procedure for Direct Padé can be described as follows: (1) Boundary conditions are substituted into (8) to generate an equation for each boundary condition. It is important to note that there exists an algebraic equation for each boundary condition.
(2) ( ) from (8) is substituted into (6), then regrouping the resulting equation in terms of powers. Afterwards, the regrouping procedure includes eliminating denominator terms arising from Padé approximant (8). In this way, the resulting expression is a power series that represent the residual error from (6).
(3) In order to reduce the residual error, from the lowest order, each coefficient from powers in the resulting power series is equalized to zero to obtain an algebraic equation in terms of the unknown coefficients of (8).
The disadvantage of the method is that it becomes as difficult as the number of terms that the expression has in the numerator as well in the denominator because it is a rational expression. Even more, the substitution into (8) gets more complicated if higher order derivatives are encountered, mainly because the first derivative of a rational expression lies in terms of another rational expression with a squared denominator and if in the nonlinear equation exist several different powers, another rational expression with even more terms is generated. For this reason, a common denominator has to be found, to formulate the system of equations 4 Discrete Dynamics in Nature and Society that allows finding and . All this process may be too cumbersome and computationally expensive. In this matter, it is concluded that if the nonlinear equation possesses higher order derivatives of (the independent variable) the number of algebraic operations is notoriously increased when powers are regrouped to form the equation system to find and . Thus, in order to avoid the algebraic complexity that implies the substitution of (8) into (6), in this paper we propose a novel strategy to circumvent such issue by obtaining the Taylor expansion of (8).

HPM Solution
Multiplying Michaelis-Menten kinetic equation (4) by [ + ( )], the following is obtained: The proposed homotopy perturbation formulation is Assuming the solution for (11) is given by a power series of substituting (12) into (11), and regrouping in terms of , Substituting the first two system solutions (13) into (12), we obtain This is the solution employing HPM, where and are constants to be determined. By substituting (14) into the second boundary condition of (5), i.e., (1) + (1) = , we obtain Solving for in (15), Substituting in (14), With the end of reducing the residual error, (17) is substituted into (4), at middle point = 0.5, between boundaries in order to cancel the error at that point, to obtain for This strategy seeks to reduce the approximation error. Substituting (18) into (17), we obtain the solution in function of , , and .

ODP Solution Procedure
The following rational solution is now proposed for (4): In order to reduce the amount of algebraic steps, it is proposed in this work, and for the first time, to expand the rational expression (21) and employ the Taylor series of order four: differentiating with respect to and to satisfy the first boundary with = 0 yields With expansion (22) we have the advantage of eliminating the denominator and thus reducing the amount of algebraic terms when substituting (21) into (4). To satisfy the second boundary condition using (21) and (5) Substituting fourth-order Taylor expansion of (21) into (10), 3 (2 2 − 2 0 2 + 2 (− 1 + 0 1 ) 1 ) + 3 (2 2 − 2 0 2 Substituting the values of , , and solving the simultaneous system of equations given by (23), (24), and (25), the numerical values of 0 , 1 , 2 , 1 , and 2 are obtained to satisfy (21).

Numerical Simulations and Discussion
For all the simulations presented in this work, with the purpose of comparison, Maple 17 numerical routines were employed. For example, numerical values of = 2.5, = 2, 6 Discrete Dynamics in Nature and Society  (19) and (21) expressed with the following values: In order to verify the precision of the solutions (26) and ODP (27) obtained by HPM, their performance is compared with ADM solutions ( 2 ( ) and 3 ( )) reported in [12], as well against the numerical solution given by Maple 17 numerical routines.
6 ( + 1) and (13) given by MTSM [13] is Equations (19) and (21) are second-order compact expressions that present a good accuracy at the interval ∈ [0, 1]. Figure 2(a) shows, in a graphical manner, that the obtained solutions with ODP and HPM present a better accuracy than the obtained for (28) and (30) In Table 1 is presented and compared the absolute relative error (A.R.E.) for each approximation according to the values of , , , for HPM and ODP solutions, including (28), (29), and (30).
The absolute relative error is defined as the absolute difference of the areas under the curve for the numerical and approximate solutions, divided by the area of the exact solution. The employed integrating method was the Rule of Simpson 1/3. In all of the presented cases, the A.R.E. obtained with ODP is much lower compared to the Adomian and MTSM methods reported in [12,13] and matches with all the graphics presented in this work.
In the same manner, for the cases ( = 4.5, = 3, = 5) and ( = 8.5, = 5, = 10) for HPM it is observed that the relative absolute error is lower than the ones presented in Adomian and MTSM. Besides, the solution obtained by HPM has the advantage of being simpler because it only consists of two terms, compared against (28), (29), and (30) that possess more terms.
Thus, for Optimized Direct Padé and HPM, the proposed results in this work coincide with [12,13,39,48,49]; that is, the oxygen concentration is higher at the substrate and diminishes as it approaches the center of spherical cell, which is due to the oxygen diffusing into the cell; it reacts to give rise to the boundary condition at the center. Likewise, results are congruent with the spherical geometry proposed by Lane-Emden (2).
Finally, it is important to remark that although the precision for ODP tends to be the same as DP, it has the advantage of having a simpler algebraic procedure, avoiding the algebraic complexity of substituting a rational expression that generates more rational expressions with higher powers in the denominator o in its case, more numerator terms, making the formulation of the system of equations difficult for finding and from the rational expression of the solution.

Conclusions
In this article were proposed two compact solutions for Michaelis-Menten oxygen uptake kinetics equation employing HPM and ODP. ODP presents a quotient of two polynomials, one of them of second order divided by a first-order one. On the other hand, HPM presents only two algebraic terms, the independent term and the second degree term. Likewise, the solutions found in this work were compared to those reported by other authors, showing that ODP has a better accuracy. HPM solution shows a good approximation in comparison to (29) obtained by Adomian; despite the fact that it only has two terms, it concluded that the solution is in good agreement. Furthermore, in this paper, we propose to replace the rational expression (8) by its Taylor series in order to reduce the amount of algebraic steps during the application of DP. The graphs presented in this work show that the solution given by ODP and HPM tends to converge as we approach = 1 unlike Adomian 3 ( ) which tends to diverge by increasing the absolute error. For the case of = 8.5, = 5, and = 10 with ODP, the A.R.E. was 9.59 × 10 −4 . Also the HPM and ODP methods were compared against the numerical one obtaining very good approximations. In addition, ODP requires fewer iterations to obtain a more compact and accurate expression than the solution proposed by the Adomian method (29).

Data Availability
The data used to support the findings of this study are included within the article.

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