Application of Fuzzy Fractional Kinetic Equations to Modelling of the Acid Hydrolysis Reaction

,


Introduction
A rigorous kinetic study of the acid hydrolysis of lignocellulosic materials is a complicated issue due to several factors: (1) the structure of the whole cells should be protected against the attacks of chemicals, (2) lignin hydrophobicity creates an arduous access of protons to the raw material, (3) interaction with other components, (4) presence of strong bonds in the raw material, and (5) variable exposition of hemicelluloses surface to the chemical attack along the reaction [1].Various researchers adapted Saeman's kinetics to describe the hydrolysis of hemicellulose and the formation of byproducts at high temperatures [2,3].These models begin with the work of Saeman for the hydrolysis of wood using sulfuric acid [4].
In the present work, our focus is to find the approximate solution of the fuzzy fractional model of the kinetics of the hemicelluloses hydrolysis reaction catalysed by sulphuric acid.We apply a family of orthogonal functions to derive a numerical method for solving this type of kinetic equation of the fuzzy Caputo fractional order.
1.1.Acid Hydrolyzing Reaction.First of all, in 1945, Faith [5] investigated the dilute acid hydrolysis process by using wood material in a 0.5 sulphuric acid solution.Neureiter et al. [6] found that the temperature has an important role in the formation of sugar degradation products.For hydrolyzing the hemicellulose and cellulose, the temperature ranges are from 100 ∘ C to 140 ∘ C. The acidity (concentration of acid) of the system is another parameter that affects the hydrolysis process.For the hydrolyzing of acid, both penetration into the biomass and dispersion in the reactor can influence not only the general reaction but also the reactor performance.The sulphuric acid diffusivity depends on the nature of lignocellulosic substances.It is revealed that the diffusivity of sulfuric acid in agricultural residues is extraordinary higher than in hard wood [7].It is with this motivation to investigate the effect of sulfuric acid concentration on both of hydrolysis of oil palm empty fruit bunch (OPEFB) fiber hydrolyzing and analyzing the appropriate kinetic model.The OPEFB biomass includes cellulose, hemicellulose, and lignin.In this study, to obtain an in depth understanding to the hydrolysis reaction and an optimum reaction condition for the process, the reaction kinetics at a complete set of reaction conditions were investigated.The respective rate constants were determined as functions of temperature and sulfuric acid concentration by a nonlinear regression analysis.

Fractional Kinetic Equation Model and Numerical Methods.
It is now well established that fractional kinetic equations represent an appropriate model to describe physical phenomena such as diffusion in porous media with fractal geometry, kinematics in viscoelastic media, relaxation processes in complex systems (including viscoelastic materials, glassy materials, synthetic polymers, and biopolymers), propagation of seismic waves, and the rate of change of chemical composition of a star [8][9][10][11][12][13][14][15][16][17][18].
It is somewhat extraordinary that recently the orthogonal functions received remarkable utilizing for the fractionalorder differential equations [49,50].Much efforts have been made to develop accurate algorithms using tau and collocation based on operational matrices of some orthogonal polynomials such as block pulse functions [51], Legendre polynomials [52,53], Chebyshev polynomials [54][55][56], Jacobi polynomials [57,58], and Laguerre polynomials [59].The main characteristic behind the approach using this technique is that it reduces FDEs to those of solving a system of algebraic equations thus notably simplifying the problem.
Furthermore, spectral methods have been found very robust tools for solving many types of fractional differential equations which have inspired many authors to apply them for these kinds of equations.The shifted Legendre operational matrix for fractional derivatives was introduced by [53] and applied with spectral methods for numerical solution of multiterm linear and nonlinear fractional differential equations subject to initial conditions which were developed by Kazem et al. [52] who applied fractional Legendre orthogonal functions for solving these types of equations.Subsequently, Doha et al. [55] introduced a new efficient Chebyshev spectral algorithms for solving linear and nonlinear multiterm fractional orders differential equations.Thereafter, Bhrawy et al. [60] propose a method to approximate multiterm fractional differential equations with variable coefficients using a quadrature shifted Legendre tau approach.Consequently, this way has been followed by several authors [54,57,59,61].
On the other hand, the study of fuzzy differential equations is rapidly developing as a new field of fuzzy mathematics.The fuzzy differential equations have been studied by several authors [62][63][64][65][66][67][68][69][70].Lately, Agarwal et al. [71] proposed the concept of solutions for fractional differential equations with uncertainty.They have considered Riemann-Liouville's differentiability with a fuzzy initial condition to solve FFDEs.Afterward, several authors have studied the existence and uniqueness of the solution of the fuzzy fractional differential equations (FFDEs) under different types of fuzzy fractional differentiability [72][73][74][75]; nevertheless, it is a little bit surprising that few papers reported numerical methods for solving FFDEs [75][76][77][78].
The scope of this paper is to derive an explicit formula for fuzzy fractional-order derivative of shifted fractional Legendre polynomials of any degree in terms of shifted fractional Legendre polynomials themselves, in the fuzzy Caputo sense with order (0 < V ≤ 1).Also, we are concerned with the direct solution technique for solving the fuzzy fractional kinetic equation (FFKE) which is extracted by using the kinetics data of the acid hydrolysis reaction subject to nonhomogeneous initial conditions.
In this paper, the fuzzy fractional derivative of the proposed FFKE is approximated based on the shifted fractional Legendre polynomials presented in [52], and then the operational matrix of the fractional Caputo derivative of order (0 < V ≤ 1) is specified, and we apply the shifted fractional Legendre spectral tau (SCT) method to construct the spectral solution for such problem.To the best of the authors knowledge, such approach has not been employed for solving linear fractional kinetic equations under uncertainty.
The paper is organized as follows.In Section 2, we introduce some necessary definitions of fuzzy sets, fractional calculus theory, kinetic equations, relevant properties of Legendre polynomials, and some of the main properties of fuzzy fractional derivatives.Section 3 is devoted to presentation of the governing fractional kinetic equation.In Section 4, the proposed method is explained for numerical solution of the FFDEs.In Section 5, the derived FFKE is solved based on the different values of constant coefficients in the equation by applying the presented technique, and the error of the approximate solution is depicted to demonstrate the effectiveness of the method.Finally, some conclusions are drawn.

Preliminaries
In this section, some basic definitions and notations related to the fuzzy sets and functions, Legendre polynomials, and fractional kinetic equations are presented which will be used throughout the paper.

Basic Definitions of Fuzzy Sets.
The basic definitions presented in this section are given in [79][80][81][82].
(iii)  is upper semicontinuous: for any  0 ∈ R, it holds that In this paper, the set of all fuzzy numbers is denoted by that contains all elements with membership degree in  greater than or equal to , that is, For a fuzzy number , its -cuts are closed intervals in R, and we denote them by According to Zadeh's extension principle, the operation of addition on R F is defined as follows: and scalar multiplication of a fuzzy number is given by where 0 ∈ R F .
Definition 3 (see [80]).The distance (, V) between two fuzzy numbers  and V is defined as where is the Hausdorff distance between []  and [V]  .
It is easy to see that  is a metric in R F and has the following properties (see [80,83]): Definition 4 (see [84]).Let  and  be the two fuzzy-numbervalued functions on the interval [, ], that is, ,  : [, ] → R F .The uniform distance between fuzzy-number-valued functions is defined by ( () ,  ()) . ( Remark 5 (see [83]).Let  : [, ] → R F be fuzzy continuous.Then, from property (iv) of Hausdorff distance, we can define Definition 6 (see [85]).Let ,  ∈ R F .If there exists  ∈ R F such that  =  ⊕ , then  is called the H-difference of  and , and it is denoted by  ⊖ .
In this paper, the sign "⊖" always stands for H-difference, and note that  ⊕  ̸ =  + (−).Also, throughout the paper, the Hukuhara-diference and generalized Hukuharadifferentiability are assumed to be existed.
The matrix form of the earlier equations is where the coefficient matrix  = (  ), 1 ≤ ,  ≤  is a crisp  ×  matrix and   ∈ R F , 1 ≤  ≤ .This system is called a fuzzy linear system (FLS).
In this part, a brief approach to the fractional calculus based on the crisp concept is provided.For more details, see [19,25,30].
Definition 10.The Caputo type fractional derivative of order V > 0 of a function  : (0, ∞) → R is defined by One of the most important advantages of using a Caputo type fractional derivative is that the Caputo derivative of a constant is zero, which means that this kind of derivative can be used to model the rate of change.
For the Caputo derivative, we have The ceiling function ⌈V⌉ is used to denote the smallest integer greater than or equal to V and floor function ⌊V⌋ to denote the largest integer less than or equal to V. Also, N = {1, 2, . ..} and N 0 = {0, 1, 2, . ..}.
Definition 11 (see [30]).Similar to the differential equation of integer order, the Caputo's fractional differentiation is a linear operation, that is, where  and  are constants.

Fractional Legendre Polynomials.
The shifted Legendre polynomials are generated from the three-term recurrence relation: The analytic form of the shifted Legendre polynomial   () of degree  is given by in which where The fractional-order Legendre functions (FLFs) can be defined by introducing the change of variable  =   and  > 0 on shifted Legendre polynomials.Let the FLFs   (  ) be denoted by FL   () [52].The fractional-order Legendre functions are a particular solution of the normalized eigenfunctions of the singular Sturm-Liouville problem Then, FL   () by using ( 20) can be obtained as follows The analytic form of FL   () of degree  given by where Note that FL   (0) = (−1)  and FL   (1) = 1.The FLFs are orthogonal with respect to the weight function () =  −1 in the interval (0, 1] with the orthogonality property A function (), square integrable in (0, 1], may be expressed in terms of FLFs as where the coefficients   are given by In practice, only the first ()-terms FLFs are considered.Then, we have with Theorem 12 (see [52]).Let the function   () ∈ (0, 1] for  = 0, 1, . . ., , (2 + 1) ≤  and is the best approximation to  from P  , then the error bound is presented as follows: where The later theorem proved that the approximate function based on the fractional Legendre polynomials converges to the function .

Fractional Kinetic Equations.
Chemical kinetics as a science began in the middle of the 19th century, when Wilhelmy [88] was apparently the first to recognize that the rate at which a chemical reaction proceeds follows definite laws, and although his work paved the way for the law of mass action of Waage and Guldberg [89], it attracted little attention until it was taken up by Ostwald towards the end of the century, as discussed by Laidler [90].Wilhelmy realized that chemical rates depended on the concentrations of the reactants [91].
Here, we provide a brief definition to illustrate the chemical kinetic process.One can find more details in [91,92].Definition 13.Chemical kinetics is the study of the rate at which a chemical process occurs.Besides information about the speed at which reactions occur, kinetics also sheds light on the reaction mechanism (exactly how the reaction occurs).
The production and destruction of species is described by kinetic equations governing the change of the number density   of species  over time [32][33][34], that is, where  is a constant.The destruction rate of the particles of type  is given by where  is a constant.Then, the residual effect can be taken as where  =  − .The solution of ( 36) can be seen as a simple function A fractional production-destruction equation can be obtained from a standard production-destruction equation by considering a fractional integral in place of a classical integral.A fractional production-destruction equation can be obtained, after dropping , as where 0  −V  represents the Riemann-Liouville fractional integral.
It is worth to note here the alternative approach of Mainardi [93], for solving the fractional relaxation-oscillation equation, that is, the fractional production-destruction equation with the Caputo fractional-order derivative.
As it can be seen, the fractional kinetic equation is obtained by replacing the standard integral with the fractional Riemann-Liouville integral.In the present research, firstly, we replace the Caputo fractional-order derivative with the classical derivative in the kinetic equation, and then the derived fractional kinetic equation in the sense of the fuzzy setting is solved by the proposed technique.

Fuzzy Fractional Differentiability.
In the fractional literature, Mittag-Leffler function plays an important role in the theory of fractional calculus/fuzzy fractional calculus and fractional differential equations/fuzzy fractional differential equations which is defined by [19,30,94]: Also, it is useful to review some basic definitions and theorems related to the FFDEs.For more details see [72,74,75].Firstly, we present some notations which are used later in the paper.(iii)  R F  [, ] indicates the set of all fuzzy-valued functions which are continuous up to order .
(iv)  R F [, ] denotes the set of all fuzzy-valued functions which are absolutely continuous.
Note that one can find easily these notations in the crisp context in [30] and references therein.Definition 14 (see [74]).
The Riemann-Liouville integral of fuzzy-valued function  is defined as Definition 15 (see [74]).
Then  is said to be Caputo's H-differentiable at  when Definition 16 (see [74]).
We say that () is fuzzy Caputo fractional differentiable of order 0 < V ≤ 1 at  0 , if there exists an element or or or For the sake of simplicity, we say that the fuzzy-valued function   [(1) − V] is differentiable if it is differentiable as in Definition 16 case (i), and  is  [(2) − V] differentiable if it is differentiable as in Definition 16 case (ii), and so on for the other cases.

The Governing Fractional Kinetic Equation
In this section, the kinetic data experimented based on the hydrolysis of OPEFB fiber at 120 ∘ C under low acid concentration conditions [95] is modelled by means of the fractional calculus.
In what follows are the detailed equation for the classic kinetic model that will be used to analyze the experimental data presented by [95].

Model
In 1945, Saeman [4] worked on dilute sulfuric acid that was catalyzed hemicellulose hydrolysis.The reaction of hemicellulose hydrolysis was modelled as a consecutive firstorder reaction with two steps.Firstly, xylan is hydrolyzed to xylose ( 1 ), and secondly the furfural was produced by degradation of xylose ( 2 ) when it is released from the solid matrix and entered the acidic solution.This classic kinetic where  = xylan;  = xylose;  = Decomposition products;  1 = sugar release rate;  2 = sugar decomposition rate.Material balance for components "" and "" for the firstorder kinetics gives in which the initial concentration at  = 0 is presented by   =   0 .Also, we have the same way for material : in which the initial concentration at  = 0 is presented by   =   0 .Equation ( 51) can be integrated, and using the given boundary condition results in Substituting ( 51) into (52) gives In this part, we apply the fractional Caputo-type derivative in (53) instead of the classical ones.So, we have

Determination of Reaction Rate Constants.
Determinations of the rate constants ( 1 ,  2 ) depend on the order of reaction.By using experimental data, the correct order would be calculated by which function of rate equation best fit the linear demand.When the order is established, the rate constants are determined from the slope of the linear plot (Table 1) by using Figures 1, 2, and 3 and (55).For a first-order reaction, the rate can be determined from the rate of formation of xylose which can be explained by equation

Shifted Fractional Legendre Spectral Tau Method
The main objective of this section is to apply the shifted fractional Legendre spectral tau (SFLT) method to construct the spectral solution for the fractional linear kinetic equation presented in Section 3, in the fuzzy concept.We now derive the operational matrix for the fuzzy fractional calculus and formulate the problem in terms of the Caputo fractional derivative.
Please notice that one can easily see why we choose such orthogonal function for this contribution in [52,53] and references therein.Now, definitions of fuzzy approximation function are established based on the application of fractional Legendre polynomials.For more details about the fuzzy approximation functions, see [96][97][98][99].

Approximation of Fuzzy Function
Definition 18.For  ∈  R F  (0, 1) ∩  R F (0, 1) and fractional Legendre polynomial FL   () a real-valued function over (0, 1), the fuzzy function is approximated by where the fuzzy coefficients   are obtained by in which FL   () is as the same in (26), and ∑ * means addition with respect to ⊕ in R F .Remark 19.Practically, only the first ()-terms shifted fractional Legendre polynomials are considered.So, we have hence were the fuzzy shifted fractional Legendre coefficient vector    and shifted fractional Legendre vector Φ  () are defined as Definition 20 (see [97]).A fuzzy-valued polynomial p * ∈ ∏ is the best approximation to fuzzy function  on  = { 0 , in which ∏ is the set of all fuzzy-valued polynomials.
The problem is addressed to the best shifted fractional Legendre approximation, as we use fractional Legendre's nodes.

Theorem 21. The best approximation of a fuzzy function based on the fractional Legendre nodes exists and is unique.
Proof.The proof is an immediate result of Theorem 4.2.1 in [97].
The following theorem defines an upper error bound for the fuzzy approximation function based on the fractional Legendre polynomials.Concerning this theorem, one can find that the fuzzy approximate function based on the FLFs is convergent to the desired fuzzy function.
Proof.Taking into account the analytic form of the FLFs explained in Section 2.2 and (26) and the fuzzy fractional Caputo definition, we have Now, by utilizing the properties of the Caputo fractional derivative explained in Section 2.1, the lemma can be proved.
The fuzzy Caputo operational matrix based on the FLFs is stated the same as the crisp concept.For more details, see [52,76].So, we have where   (V) is the × operational matrix of fuzzy fractional Caputo's derivative of FLFs and   V Φ() ∈  R F [, ].So, using (67) and Lemma 23, we extend the operational matrix of Caputo fractional-order derivative of FLFs in the sense of the fuzzy setting as follows.
Theorem 24.Let one assume that Φ() is the FLF vector.
is the -square operational matrix of fuzzy fractional Caputo's derivative of order 0 < V ≤ 1.Then, the elements of  (V) , are achieved as Proof.If we consider the fractional derivative in the crisp context, then the proof of this theorem is obtained from Theorem 2 in [52].In the same way, if we consider the fuzzy Caputo fractional derivative, then analogously to the demonstration of Theorem 1 in [52] and by using Lemma 23, we can prove the relation (69).
Also, we can exhibit the operational matrix  (V) in the alternative form as Remark 25.If V = 1, then Theorem 24 gives the operational matrix for the shifted Legendre polynomials function as described in Section 2.2.

The Application of the FLFs for Solving Linear FFKE.
In this section, we are concerned with providing a numerical solution to linear FFKE which was described in Section 3. We approximate the fuzzy function () by means of the FLFs and then substitute in the FFKE to derive the approximate solution.The method reduces the fuzzy fractional initial problem to a system of fuzzy algebraic equations.The large systems of algebraic equations may lead to greater computational complexity and large storage requirements.However, the operational matrix for the FLFs is structurally spare.This reduces the computational complexity of the resulting algebraic system.Fractional kinetic differential equation was derived in the deterministic case in Section 3. Now, in order to investigate the mentioned problem in a real case, we use the fuzzy initial value   0 , the fuzzy-valued function (), and the concept of Caputo's H-differentiability for fractional derivative of   (), (   V 0+   )() and generalized H-differentiability [67] for firstorder derivative of   ().
Let us consider the fuzzy version of the linear fractional kinetic differential equation in the general form: where   () :  R F ∩  R F is a continuous fuzzy-valued function,   V 0 + () denotes the fuzzy Caputo fractional derivative of order V, and () : [0, 1] → R F .
Firstly, we state the unknown fuzzy functions   (), (   V 0 +   )() and known fuzzy function () in terms of the FLFs as follows: The following theorem provides a suitable way to reach the fuzzy unknown coefficient of the fuzzy approximate function     (), by means of the fuzzy residual function (  ()) of the problem stated.

Numerical Results
In order to illustrate the effectiveness and accuracy of the proposed method, we carry it out for solving the fractional kinetic equation which was derived in Section 3 in the sense of the fuzzy concept.Comparison of the results obtained with the exact solution reveals that our method is efficient and convenient.Now, we consider (54) with the fuzzy initial condition.It means that the initial concentration of xylose is a fuzzy number.So, we have Example 27.Consider the following FFDE: in which ] is a continuous fuzzy function,   V 0 + denotes the fuzzy Caputo fractional derivative of order V, and  1 ,  2 are constant coefficients specified from Table 1.
Regarding the definition of the Caputo fuzzy fractional differentiability for  [1 − V] differentiability and Theorem 17, we can state (81) in a parametric form as follows: where   (0; ) = [0.5 + 0.5, 1.5 − 0.5].The exact solution of (81) can be gained using (82) as Employing the proposed method described in Section 4, we can derive the operational matrix of the FLFs based on the fuzzy fractional derivative of the Caputo type and the approximate fuzzy function of the solution by means of the FLFs.So, we can express (81) as   Using Theorem 26, we can derive the fuzzy residual function for (81).Afterward, employing the inner product explained in Section 4.3 and (81), we can generate fuzzy linear algebraic system with (−1) equations which is in the matrix form as follows: Additionally, we approximate fuzzy initial conditions (81) by using (80) to produce the last equation which is needed to complete our fuzzy algebraic linear equations system.So, we have From ( 85) and ( 86), ()-fuzzy linear equations are produced which on can solve it easily using Definitions 8 and 9 and the method presented in [86] to find the unknown fuzzy coefficients, {  } −1 =0 , of the fuzzy approximate solution.We solved the problem by applying the technique described in Section 4. With  = 3, V =  = 0.75 and from Table 1, we assume that  1 = 0.012,  2 = 0.001.We approximate the solution in terms of the lower and upper approximate fuzzy functions as where our aim is to find the unknown fuzzy coefficients {  } 2 =0 .The operational matrix of the FLFs is as The FLFs for  = 2 are as follows: Here, we assume that -cut = 1, so the right-hand side coefficients vector in terms of the lower and upper parametric form of the fuzzy number is as follows: Putting ( 88) and ( 90) in ( 85), we reach a fuzzy linear algebraic equations system.Unknown coefficients vectors  12 and  22 are obtained by solving this system the values for -cut = 1 are as and from the earlier vector, we can attain the approximate fuzzy function as Now, in the following results, we will reach the approximate fuzzy solutions for (81) regarding different values of  1 and  2 from Table 1.Actually, we achieve the concentration of xylose after a specific time from the stated fuzzy fractional kinetic equation (81).
In Table 2, we present the approximate solution and the absolute error with  1 = 0.012,  = 7, and  2 = 0.001 at  = 1.As expected, the absolute error demonstrates the spectral accuracy of the proposed method.The absolute error for V =  = 0.5, 075, 0.85, 0.95 is plotted in Figure 4.As we can see, the absolute error increases gradually but remains in the interval [10 −10 , 10 −12 ].This error can be explained by the computer round-off errors that prevent any further accuracy improvements.This error can be explained by the computer round-off errors that prevent any further accuracy improvements.Also, Figure 5 shows the absolute error for different number of the FLFs which is clear that decreasing the the absolute error decreases with the increasing number of functions occurring.Moreover, in Figure 6, the absolute  error obtained by our method for V = 0.75 and  = 7.As one can see, the increase in the amount of  cannot reduce the absolute error noticeably.The approximated solutions are evaluated for V =  = 0.5, 0.75, 0.85, 0.95 and  = 8.The results of the numerical simulations are plotted in Figure 7.
The approximate solutions for V =  = 0.85,  = 9,  1 = 0.014, and  2 = 0.001 are obtained in Table 3.The obtained numerical results are in very good agreement with the exact solutions.In Figure 8, the absolute errors for the problem of V =  = 0.5, 0.75, 0.85, 0.95 by the proposed method with  = 8 are given.It is obvious that the improvement of the value of the V and  could not affect the absolute error considerably, but the method could achieve a suitable error in the interval [10 −9 , 10 −14 ].The absolute error for V = 0.85 and the various choices of  and  by the presented method are shown in Figures 9 and 10, respectively, to make it easier to compare with the analytic solution.We show the approximate solution in Figure 11 for V =  = 0.5, 0.75, 0.85, 0.95 with  = 12.
Table 4 exhibits the error of analytical and numerical solutions.Numerical results demonstrate the good accuracy of the proposed method.The absolute error is evaluated for various choices of  in Figure 12.We see that as V,  approaches 1, the solution of the fuzzy fractional kinetic equation approaches that of the integer order fuzzy kinetic equations.In Figure 13, again we see that the method could reach a high accuracy for the approximate solution, but the increase in the number of the FLFS could not affect the absolute error considerably.Also, in Figure 14, the method could hit the better accuracy when the values of the fractional derivative and fractional order of the Legendre functions were the same.In Figure 15, we show the fuzzy approximate solutions which were evaluated for various values of V = .

Conclusion
In this paper, a spectral method has been adapted for numerical solution of fuzzy fractional kinetic equation of order 0 < V ≤ 1.A specific family of the orthogonal polynomials which is so-called FLFs was used as the approximation basis.The proposed method is characterized by its simplicity, efficiency, and high accuracy.The accuracy and validity of the presented method were demonstrated through numerical simulation for the derived FFKE.Using the proposed method, we could reach the suitable approximation of the amount of the concentration value of xylose after a determined time that is important to analyze the kinetic data in the chemical process.
A direction for further research would be to extend the presented method to the other types of the fuzzy fractional differential equations like fuzzy fractional oscillation differential equations of the distributed order.

Figure 3 :
Figure 3: Comparison between calculated and experimental value of xylose in hydrolysis of OPEFB with 6% H 2 SO 4 at 120 ∘ C [95].
Development.By the addition of water to the xylan, the hydrolysis occurs and xylose is produced