Solution of non-linear boundary value problems in immobilized glucoamylase kinetics

A mathematical model to describe the enzyme reaction, mass transfer and heat effects in the calorimetric system is discussed. The model is based on non-stationary diffusion Equation containing a non-linear term related to immobilize liver esterase by flow calorimetry. This paper presents the complex numerical methods (Adomian decomposition method, Homotopy analysis and perturbation method) to solve the nonlinear differential Equations that depict the diffusion coupled with a non-linear reaction terms. Approximate analytical expressions for substrate concentration have been derived for all values of parameters α, β and γE. These analytical results are compared with the available numerical results and are found to be in good agreement.


INTRODUCTION
Flow reaction calorimetry has several advantages over a batch calorimetry method.The operation at a calorimetric experiment can be made exceedingly simple and equilibration time prior to the experiment can be omitted.Mixing of reactants can be achieved without the presence of a gaseous phase which is of great importance when experiments are performed with volatile liquids and in micro-calorimetric expriments where very small conden-sation-evaporation effects may affect the result.Surface adsorption effects which may cause seious systematic errors in micro calorimetry can be neglected if a steady liquid flow is allowed to continue until possible wall reactions have occurred [1].Immobilized biocatalysts (IMB)enzymes or whole cells are used in various areas of analytical, medical, and industrial applications.Basically, enzyme kinetic parameters cannot be determined experimental data.For this purpose many experimental techniques can be used, that are more or less laborious and time consuming.Reaction kinetics of carboxyl esterase's depends strongly on the nature of substrate.The hydrolysis of different substrate activation [2] or it can follow the simple Michaelis-Menten kinetics [3].
Stefuca et al. [4] have described the principles and applications of flow calorimetry (FC) in the investigation of the IMB properties.One of the last improvements of this technique was the introduction of an "auto calibration" principle based on reaction solution recirculation enabling to determine true reaction rate of biocatalyst reaction without any requirement of an additional analytical technique [5].Vladimir Stefuca et al. [6] have derived the experimental data were treated by mathematical modelling based on material and heat balances of the reaction system.Recently, Fedor Malik [7] has developed the mathematical model describing the enzyme reaction, mass transfer and heat effects in the calorimetric system.
To my knowledge no rigorous analytical solutions of the substrate of phenyl acetate hydrolysis with steadystate conditions for all values of parameters ,   and E  have been reported.The purpose of this communication is to derive approximate analytical expressions for the steady-state concentration of substrate using Adomian decomposition method, Homotopy analysis and perturbation method.

MATHEMATICAL FORMULATION OF THE PROBLEM
The experimental set-up used for the capacity is depicted in (Figure 1).The main part of the system thermostatic cell through immobilized enzyme column.The column was operated packed bed reactor The temperature difference between the column input and output, , is measured by thermistors and registered by a personal computer.The system was kept at temperature of 303.15 K, while the buffer solution was continuously pumped through the column with constant flow rate of 1 ml/min.The experiment was started by replacing the buffer solution by the substrate solution containing 1 -200 g/l of MDX in 0.1 M acetate buffer (pH 4.7).Two techniques of measurement were applied: single flow mode and total recirculation mode.The single flow mode was performed with the switching valve 2 opened to the waste [7].The substrate concentration gradient on the particle surface was calculated by the equation of substrate balance in the particle [7]: The equation must be solved subject to the following initial and boundary conditions: where SP is the substrate concentration in the particle, S is the phenyl acetate concentration, e is the diffusion coefficient, are the kinetic parameters and r is the particle radial co-ordinate, is the particle radius.We can write the steady state equation as [7]: The system governs the substrate concentration SP when there is no competitive inhibition in the reaction.The non-linear ODE (Eq.5) is made dimensionless by defining the following parameters: where E  denote the reaction diffusion parameter, x is the dimensionless distance and

 
U x is the dimensionless concentration.Here  and  denotes the saturation parameters.The above Eq.5 reduces to the following dimensionless form The corresponding boundary conditions are

SOLUTION OF BOUNDARY VALUE PROBLEM USING ADOMIAN DECOMPOSITION METHOD
Adomian's decomposition method has been successfully applied to linear and nonlinear problems.One of its advantages is that it provides a rapid convergent solution series [8].However, the method applied to nonlinear Equations does not seem to be fast enough to be a efficient method to solve these kind of equations and one can find in the open literature some modifications proposed by several authors [9][10][11][12][13].By applying the Adomian's decomposition method, a new iterative method to compute nonlinear equations is developed and is presented in this work.The Adomian decomposition method is an extremely simple method [9][10][11][12][13] to solve non-linear differential equations.First iteration is enough.Furthermore, the obtained result is of high accuracy.Using this Adomian decomposition method (see appendix A and B), the solution of Eq.7 becomes:

SOLUTION OF BOUNDARY VALUE PROBLEM USING HOMOTOPY ANALYSIS METHOD
Perturbation methods are the most famous analytic techniques for nonlinear problems, which are widely applied in science and engineering.In 1992, the homotopy, a traditional concept in topology, was used by Liao [14] to propose an approximation technique for nonlinear problems, namely the homotopy analysis method (HAM).Using the concept of the homotopy, a nonlinear problem is transformed into a sequence of linear sub-problems that are easy to solve by means of the symbolic computation software.In 1997 Liao [14] further generalized the HAM by introducing an auxiliary nonzero parameter (called today the convergence-control parameter).Different from perturbation techniques, the HAM does not depend upon any small physical parameters, and besides provides great freedom to choose different base functions to approximate nonlinear problems.Especially, different from all other analytic approximation methods, the socalled convergence-control parameter of the HAM provides us a convenient way to ensure the convergence of series solution.Thus, the HAM overcomes the restrictions of the perturbation methods and therefore is more general.With these advantages and having the aid of high-performance computer and symbolic computation software, the HAM has been widely applied to solve many types of nonlinear differential Equations in science, engineering and finance [15].Using this HAM (see Appendix C and D) we obtain, the concentration of substate as follows 3 120 The analytical solution represented by Eq.11 contains the auxiliary parameter h, which gives the convergence region and rate of approximation for the homotopy ana-lysis method.The analytical solution should converge.It should be noted that the auxiliary parameter h controls the convergence and accuracy of the solution series.In order to define region such that the solution series is independent of h, a multiple of h curves are plotted in Figure 2.

SOLUTION OF BOUNDARY VALUE PROBLEM USING HOMOTOPY PERTURBATION METHOD
The Homotopy perturbation method which doesn't need small parameter is implemented for solving the differential Equations and it is predicted that HPM can be founded widely applicable in engineering and in cases that don't have exact solution this method can be used as semi-exact solution.Homotopy perturbation method yields a very rapid convergence and usually, one iteration leads to high accuracy of solution [17][18][19][20][21][22][23][24][25].The Homotopy perturbation method is a high accuracy method.Using this method (see Appendix E and F) we obtain

NUMERICAL SIMULATION
The non-linear differentials Eqs.7-9 are also solved by numerical methods.The function bvp4c in Matlab software which is a function of solving two-point boundary value problems (BVPs) for ordinary differential equations is used to solve this Equation.The Matlab program is also given in appendix G. Numerical values of the parameters used in Fedor Malik et al. [7] in this work is given in Table 1.Its numerical solution is compared with Adomian decomposition method, Homotopy analysis and perturbation method in Table 2-5 and it gives satisfactory result when α ≤ 1 and β ≤ 1.

RESULTS AND DISCUSSION
The primary result of this work is the first approximate and simple expression of concentrations of substrate (Eqs.10(ADM), 11 (HAM) and 12 (HPM)).From this figures, it is apparent that the value of the concentration of substrate increases for various values of  increases.

CONCLUSION
A non-linear ordinary differential equation in the investigation of kinetics of immobilized liver esterase by flow calorimetry has been solved using Adomian decomposition method, Homotopy analysis and perturbation method.The simple approximate expression of concentration of substrate for all values of parameters ,   and E  are reported.These methods can be easily extended to find the solution of all other non-linear reaction diffusion equations for immobilized enzymes with reversible Michaelis-Menten kinetics for various complex boundary conditions.These analytical results are useful for design and optimization of immobilized liver esterase by flow calorimetry.

APPENDIX A Basic Concept of the Adomian Decomposition Method (ADM)
Adomian decomposition method [9][10][11][12][13] depends on the non-linear differential equation into the two components where L and N are the linear and non-linear parts of F respectively.The operator L is assumed to be an invertible operator.Solving for leads to Applying the inverse operator L on both sides of Eq.A.3 yields where   x  y is the constant of integration which satisfies the condition Now assuming that the solution can be represented as infinite series of the form Furthermore, suppose that the non-linear term   N y can be written as infinite series in terms of the Adomian polynomials n A of the form where the Adomian polynomials n   N y are evaluated using the formula: Then substituting Eqs.A. 5 and A.6 in Eq.A. 4 gives Then equating the terms in the linear system of Eq.A.8 gives the recurrent relation However, in practice all the terms of series in Eq.A.7 cannot be determined, and the solution is approximated by the truncated series .This method has been

APPENDIX B Analytical Solutions of Concentrations of Substrate Using ADM
In this appendix, we derive the general solution of nonlinear Eq.7 by using Adomian decomposition method.We write the Eq.7 in the operator form, where . Applying the inverse operator 1 L  on both sides of Eq.B.1 yields where A and B are the constants of integration.We let, We identify the zeroth component as and the remaining components as the recurrence relation where A n are the Adomian polynomials of .We can find the first few A n as follows: The remaining polynomials can be generated easily, and so, Adding (B.11) to (B.13) we get Eq.11 in the text.

APPENDIX C Basic Concept of Homotopy Analysis Method
Consider the following differential equation [15]: where Ν is a nonlinear operator, t denotes an independent variable; u(t) is an unknown function.For simplicity, we ignore all boundary or initial conditions, which can be treated in the similar way.By means of generalizing the conventional homotopy method, Liao constructed the so-called zero-order deformation equation where is the embedding parameter, 0 h  is a nonzero auxiliary parameter, is an auxiliary function, L is an auxiliary linear operator, , is an unknown function.It is an important, that one has great freedom to choose auxiliary unknowns in HAM.Obviously, when p = 0 and p = 1, it holds: where If the auxiliary linear operator, the initial guess, the auxiliary parameter h, and the auxiliary function are so properly chosen, the series (C.4) converges at p = 1 then we have: Define the vector Differentiating Eq.C.2 for m times with respect to the embedding parameter p, and then setting p = 0 and finally dividing them by m! we will have the so-called -order deformation Equation as: where L  on both side of Eq.A8, we get In this way, it is easily to obtain for at when , we get an accurate approximation of the original Eq.C.1.For the convergence of the above method we refer the reader to Liao [15].If Eq.C.1 admits unique solution, then this method will produce the unique solution.If Eq.C.1 does not possess unique solution, the HAM will give a solution among many other (possible) solutions.

APPENDIX D Approximate Analytical Solutions of the System of Eqs.7-9 Using HAM
In this appendix, we indicate how Eq.11 in this paper is derived.The Homotopy analysis method was constructed to determine the solution of Eqs.7-9.
In order to solve Eq.D.1 by means of the HAM, we first construct the zeroth-order deformation Equation by taking The approximate solutions of Eq.D.2 are as follows Substituting the series (D.3) in Eq.D.2 and equating the like powers of p we get The boundary conditions becomes Now applying the boundary conditions Eq.D.7 in Eq.D. 4 we get Substituting the values of 0 in Eq.D.5 and solving the Equation using the boundary conditions Eq.D.8 we obtain the following result: Substituting the values of U 0 and U 1 in Eq.D.6 and solving the Equation using the boundary conditions Eq. D.8 we obtain the following result: To find few iteration we get, the solution of

APPENDIX E Basic Concept of the Homotopy Perturbation Method (HPM)
We outline the basic idea of Homotopy perturbation method.This method has eliminated the limitations of the traditional perturbation methods.On the other hand it can take full advantage of the traditional perturbation techniques, so there has been a considerable deal of research in applying homotopy technique for solving various strongly non-linear equations.To explain this method, let us consider the following function with the boundary conditions of where A , , B   f r and  denote a general differential operator, a boundary operator, a known analytical function and the boundary of the domain  , respectively.Generally speaking, the operator A can be divided into a linear part L and a non-linear part N. Eq.2.1 can therefore be written as By the homotopy perturbation technique, we construct where is an embedding parameter, and 0 is an initial approximation of Eq.E.1, which satisfies the boundary conditions.Obviously from Eqs.E. 4 and E.5, we will have when p = 0 Eq.E.4 or E.5 becomes a linear Equation; when p = 1 it becomes a non-linear Equation.So the changing process of p from zero to unity is just that of We can first use the embedding parameter p as a "small parameter", and assume that the solutions of Eqs.E. 4 and E.5 can be written as a power series in p Setting p = 1, results in the approximate solution of Eq.E.1: The combination of the perturbation method and the Homotopy method is called the Homotopy perturbation method.

APPENDIX F
Approximate Analytical Solutions of the system of Eqs.7-9Using HPM Solution of the Eqs.7-9 using Homotopy perturbation method.In this appendix, we indicate how Eq.12 in this paper is derived.Furthermore, a Homotopy was constructed to determine the solution of Eqs.7.
The approximate solutions of Eq.F.1 is Substituting Eq.F.2 into Eq.F.1, and comparing the coefficients of like powers of p The boundary conditions are Solving the Eqs.F.3 to F.5 and using the boundary conditions (F.6) and (F.7), we can find the following results

 
6 120 According to the HPM, we can conclude that Using Eqs.F.8-F.10 in Eq.F.11, we obtain the final results are described in Eq.12.

APPENDIX G
In this appendix, we derive the solution of Eq.F.4 by using reduction of order.To illustrate the basic concepts of reduction of order, we consider the Equation where P, Q, R are function of x.Eq.F.4 can be simplified to Using reduction of order, we have   2 2 1 2 ; 0; 6 Substituting the value of P in the above Eq.G.5 becomes The given Eq.G.2 reduces to 6 (G.9) ql = 1; pr = ur-1; qr = 0; Integrating Eq.G.9 twice, we obtain

APPENDIX I Determining the Region of h for Validity
Substituting (G.6) and (G.10) in (G.4) we have, The analytical solution should converge.It should be noted that the auxiliary parameter h controls the convergence and accuracy of the solution series.The analytical solution represented by Eq.11 contains the auxiliary parameter h, which gives the convergence region and rate of approximation for the Homotopy analysis method.In order to define region such that the solution series is independent of h, a multiple of h-curves are plotted.The region where the distribution of and Using the boundary conditions Eqs.F.6 and F.7, we can obtain the value of the constants A and B. Substituting the value of the constants A and B in the Eq.G.11 we obtain the Eq.F.10.Similarly we can solve the other differential Eqs.B.11, D.4-D.6, F.3 and F.5 using the reduction of order method.

Figures 3- 6 Figures 3 - 6 , 7 - 9 .
show the analytical expression of concentration of substrate for various values of dimensionless reac-it is inferred that the value of the concentration of substrate increases when dimensionless reaction diffusion parameter in these Figures3 to 6, it is known that the value of the concentration of substrate increases gradually and attains the maximum at the boundary   1 x  .The normalized numerical simulation of three dimensional substrate concentration is shown in Figures The time independent concentration slowly decreasing when E  is increasing.Then the concentration of u x   1  when 1 x  and also for all values of E  ,  and  .In these figure, it should be noted that the value of the concentration of substrate decreases for all values of E  .
be very efficient in solving various types of non-linear boundary and initial value problems.
versus h is a horizontal line is known as the convergence region for the corresponding function.The common region among   U x and its derivatives are known as the over all convergence region.To study the influence of h on the convergence of solution, the h-curves of clearly indicate that the valid region of h is about (−2 to −0.5).Similarly we can find the value of the convergence control parameter h for different values of constant parameters.