Mathematical Analysis of Oxygen Uptake Rate in Continuous Process under Caputo Derivative

In this paper, the wastewater treatment model is investigated by means of one of the most robust fractional derivatives, namely, the Caputo fractional derivative. The growth rate is assumed to obey the Contois model, which is often used to model the growth of biomass in wastewaters. The characteristics of the model under consideration are derived and evaluated, such as equilibrium, stability analysis, and steady-state solutions. Further, important characteristics of the fractional wastewater model allow us to understand the dynamics of the model in detail. To this end, we discuss several important analyses of the fractional variant of the model under consideration. To observe the efficiency of the non-local fractional differential operator of Caputo over its counterclassical version, we perform numerical simulations.


Introduction
Mathematical models as routes of interpreting real-life situations are included in the literature. These models display a vital update in the quantification and assessment of real-life problems and preventive measures [1][2][3]. Mathematical modeling has, in many ways, proven to be a very versatile and effective way of studying the dynamics of many situations. To build and test control measures that are persuasive, statistical analysis and numerical simulations may be used. It is also widely accepted that by using mathematical models, the prediction to control and solve real-life problems can be reached. A relevant argument here is that it is important to obtain acceptable criteria for the particular situation under consideration and to use these criteria to determine the effects of potential control measures. The dilemma, then, is how, from an optimal point of view, to implement such measures. Several noteworthy attempts have recently been made to establish results through mathematical modelling [4][5][6][7][8].
One of the key areas to see the application of mathematical models is the field of biochemical engineering, where oxygen is a key substrate in aerobic bioprocesses that is used for growth, maintenance, and other metabolic routes, including product synthesis [9]. Oxygen is continuously supplied via the gas phase due to its low solubility in broths, which are typically aqueous solutions, and thus, knowledge of the oxygen transfer rate is needed for the design and scale-up of the bioreactor, which can be addressed by a mathematical modelling approach. Both the transition of oxygen mass from the phase of liquid to gas and its microorganism intake are of decisive importance, as the transfer of oxygen is the control stage for microbial growth in many processes and may influence the evolution of bioprocesses. Early in the history of biochemical engineering, this fact was established and reviewed [10][11][12][13][14]. The oxygen uptake rate is one of the essential physiological features of culture growth, and has been used to optimize the fermentation method. The real oxygen uptake rate is usually estimated from the experimentally defined oxygen uptake rate [15,16].
We consider a wastewater treatment model with the Contois model for the growth rate. As the main feeds into a tank, the aerator is named [17][18][19], and the usage theory can be established. At this point, by biological oxidation of the substrate, one bacteria population degrades the pollutant. By absorbing oxygen, the reaction produces an aerobic environment. The mixture is transmitted to a settler tank in the second stage. Here, because of gravity, the solid pieces settle and collect at the bottom. Part of the sludge inside the aerator is recycled to promote oxidation. The topic of waste water treatment has also been studied from many points of view, such as complex observation [20,21] and control [22,23], although the anaerobic component is considered in the majority of situations [24][25][26]. However, this model has never been studied using the contemporary fractional operators with the Contois model. Additionally, it has been multifariously proved that natural phenomena can be more effectively enunciated by fractional models than by differential integer-order equations. The fractional calculus has taken on the value and popularity of modeling realistic cases, in particular those with memory effects, because of this advantage [27][28][29][30]. Fractional operators have been used in so many fields where [31][32][33][34][35] is used. We were inspired to study and evaluate a new fractional version of Caputo's operator model given this value. This is the first time that the Caputo operator has been hired for the model being considered, to the best of our knowledge.

Formation of the Model
We present the model to be investigated in this research here. The model with dimensional form is given by: The oxygen accumulation in the liquid phase can be expressed as where µ(S, X) is the Contois model for the specific growth rate, which is expressed as The rate of oxygen transfer from the gas to the liquid phase is provided by -Y x/O 2 is the overall yield of cells from oxygen; and -K Ia is the oxygen mass transfer coefficient.
A detailed definition of the model can be found in Nomenclature.

The Non-Dimensional Model
Non-dimensionalization is very important in mathematical modeling. Here, we aim to non-dimensionalize (1) and (2) to study various mathematical properties of the model. Equations (1) and (2) can be written in non-dimensional form by the following equations using non-dimensional variable C * = C S 0 for oxygen concentration, S * = S S 0 for the concentration of substrate, X * = K s X S 0 as the microorganism concentration, and t * = µ m t as time In Equation (4), τ * is the main experimental control parameter. It is assumed for the investigation of the mathematical model that the substrate concentration entering into the bioreactor is positive S 0 > 0 and the rate of change of death is positive K * d > 0. In addition, it is assumed that 0 < K * d < 1 and A1: All the parameters are nonnegative.
In Equation (4), represents the dimensionless death rate, • τ * = Vµ m F represents the dimensionless residence time, • K * Ia = K I a µ m represents the dimensionless oxygen mass transfer coefficient, and In the rest of the analysis of Model 4, we drop the star symbols from Equation (4).

Analysis of the Model in Classical Sense
In this section, our aim is to analyze the properties of the model in a classical sense.

Equilibria and Stability Analysis
We consider the number of equilibrium solutions of the model (4). It is obvious that the model has a branch of the washout given by One can obtain the model's non-free steady-state (4) by setting the right side to zero. From the model (4), we have, Substituting (6) into the third equation of (4), we obtain Next, we focus on analyzing the stability of the equilibrium of (4). The stability of the two forms of the equilibrium, E 0 and E 1 , is studied. Lemma 1. The steady-state solution E 0 is locally asymptotically stable when τ < β 1−K d , and is unstable when τ > β Proof. The Jacobian matrix of system (4) at E 0 is given by The eigenvalues of (8) is expressed as Note that λ 1,3 are negative. After some algebra, we notice λ 2 < 0 when τ < β

Lemma 2.
The solution E 1 is locally asymptotically stable for a physically meaningful solution.
Proof. For the stability of the solution E 1 , it is possible to write the Jacobian of the system as where, As a simple analysis, J (E 1 ) characteristics are expressed as where We use the Rough-Hurtwiz criterion, that is, Thus, the condition of the Routh-Huewitz theorem for the stability of E 1 is satisfied.

Numerical Simulations
The model's steady-state solution is evaluated in this section with respect to parameters using the experimental values in Table 1. Figure 1 shows stable curves of the oxygen concentration, the substrate concentration, and the microorganism concentration as functions of residence time for different values of the parameter β. Note that the parameter β represents the configuration of the reactor. β = 1 corresponds to the continuous flow reactor, while β = 0 corresponds to the membrane reactor.

Parameters
Values Unit Reference β 0-1 -Range value of the parameter Figure 1a shows that the oxygen concentration decreases as a function of residence time to minimum value, called the minimum residence time. After this minimum value, the oxygen concentration increases as the value of the residence time increases. Note that the minimum residence time for β = 1 and β = 0.1 are 1.30025 and 0.201052, respectively. Thus, the value of parameter β must be taken lower to minimise the oxygen concentration.
In Figure 1b, the substrate concentration decreases as the value of the residence time increases for different values of β. In Figure 1c, the microorganism concentration increases to the maximum value, then decreases as the value of the residence time.     Figure 1a as a function of β and the influent concentration (S 0 ). We note that the minimum residence time increases as the value of these parameters increases.

Analysis of the Model under Caputo Fractional Operator
Here, we aim to extend the model (4) by incorporating the fractional order in the sense of a Caputor fractional derivative. The model (4) in the sense of a Caputo fractional operator is given by
Proof. Consider that P : E → E is defined by Now, P is the unique solution of (13) and is well-defined and depicts the fixed-point of P. Consider sup t∈J Λ(t, 0) = M 1 and κ ≥ 0 + ΩM 1 . Now, it suffices to justify PH κ ⊂ H κ , and H κ = { ∈ E : ≤ κ}, as convex and closed. Thus, each ∈ H κ , we have We justify the results. Further, for 1 , 2 ∈ E, one attains This justifies that (P 1 ) − (P 2 ) ≤ ΩL Λ 1 − 2 . Thus, due to Banach contraction, the unique solution for (13) is reached. Now, we go for the existence of (13) solutions by Krasnoselskii's fixed-point justification.
Lemma 3. Let M = ∅ be a closed, bounded, and convex subset of a Banach Space E. Let two operators that respect the given relation be P 1 , P 2 .
Then, there is u ∈ M as far as u = P 1 u + P 2 u.

Theorem 2.
Surmising Λ : J × R 5 → R is continuous and holds for the condition (A 1 ).

Iterative and Stability Analysis via Caputo Operator
Let (B, ||.||) be a Banach space and Y * be a self-map of B. Further, consider y n+1 = h(Y * , y n ) and G(Y * ) as a fixed-point set of non-empty Y * . It can be interpreted that y n has converged to q * ∈ G(Y * ). Defining To get convergence for z n , it has to possess an upper bound. Having satisfied the aforementioned conditions for y n+1 = Y * with n representing Picard's iteration as in [40], then the iteration is Y * -stable. Thus:

Theorem 3 ([40]
). Let (B, ||.||) be a Banach space and Y * be a self-map on B. Then, for all x, y ∈ B, the inequality as comes next holds with K ≥ 0, k ∈ [0, 1). If we assume that Y * is Picard Y * -stable, we can write Theorem 4. Let F be a self-map. Then, it is defined as follows: Proof. Since F is a fixed point, thus, for each (m, n) ∈ N × N, one attains Imposing norms on both sides of the first equation in (29), we attain and imposing a triangular inequality gives F (S n (t)) − F (S m (t)) = S n (t) − S m (t) Thus, (31) becomes To achieve our goal, we assume that Inserting (33) into the relation (32), one can have Since S m (t), X m (t), and C m (t) get convergent and bounded, there exists P > 0, Q > 0, and R > 0 for every t. Thereby, we get S m (t) < P, X m (t) < Q, C m (t) < R, (m, n) ∈ N × N. (35) From (34) and (35), one reaches with κ 1 , κ 2 , and κ 3 being acquired functions by the inverse Laplace transform in (34). In a similar way, we achieve (38) where the condition (28) is valid. Now, F has a fixed point. To prove that F holds for Theorem 3.3, let (37) and (38) hold, and the following relations are satisfied: Hence, we reach the intended result.

The Solution's Positiveness
Let us decide the invariant region and indicate that for all t ≥ 0, all solutions of the model's fractional system are positive. The main objective is to implement the convenience of the analyzed model solutions by observing whether they join the Υ invariant field. Employing the benefits of the fractional derivative of Caputo, we assume that is a model solution with positive initial conditions. Further, n = (S(t), X(t), C(t)) T . We must also show that on each hyperplane that is bound by the non-negative hyperoctant, the vector field points to R+ 3 . Thus, one can write Thus, we give the convenient region as follows: Thus, in the Υ area, the fractional model is biologically acceptable and mathematically well-defined when t > 0. In addition, this area is positively invariant, that is, the underlying system's solutions are positive for all t.

Numerical Dynamics
This is where, with the aid of the best-suited parameters as obtained in the table, we present different numerical results for the proposed model by using a fractional differential operator known as the Caputo operator. To serve the purpose of model simulations, a predictor-corrector numerical technique for the fractional type of dynamical systems as implemented and analyzed in [41][42][43] has been used on an effective software called Matlab 2019a. In terms of the Caputo differential operator of order ψ, the following form of the Cauchy ordinary differential equation is taken into consideration: with p = 0, 1, . . . , n − 1, n = ψ . The Volterra equation version of the Equation (43) becomes: Applying the described predictor-corrector related with the Adam-Bashforth-Moulton solver [42], we assume h = τ/N, t j = jh, and j = 0, 1, . . . , N ∈ Z + , by taking ϕ j ≈ ϕ t j , with ϕ j standing for the approximate solution and ϕ(t j ) for the true one. We can now write Here, (p) implies the predictor portion and it can be determined as: The numerical simulations of the Caputo model (13) using the parameters mentioned in Table 2 are presented in this section. In order to solve the nonlinear model shown in Equation (13), we used the predictor-corrector Adams technique and obtained graphical results based on the given parameters values. The concentration of substrate S(t), of oxygen C(t), and the microorganism concentration X(t) are investigated using different values of the fractional order ψ in order to see the dynamical behavior of the Caputo model (13) for all three compartments. It can be observed in Figure 3 that the behaviour of the curves changes as the value of fractional order increases. Moreover, the values of one of the most important and effective parameters in the model, which is the main experimental control parameter τ, have been varied in Figure 4 to practically see its influence in the model. Comparing Figures 3 and 4, the changes in the curves can be noticed as a result of the effect of the main experimental control, τ. (c) Figure 3. Numerical simulation of the oxygen saturation concentration, the substrate concentration, and the microorganism concentration for different values of fractional order. The value of the parameters is stated in Table 1

Conclusions
In the current analysis, the wastewater model has been studied from a classical viewpoint and modelled by one of the most robust non-local fractional operators, called the Caputo operator. Firstly, some important properties of the model in a classical sense were reported. The dynamical features of the model were depicted through the steady-state region in its classical sense. Since it has been proved in the literature that the fractional operators are especially suitable for studying the dynamical behaviour of real-world problems, it is of vehement importance to analyze more critically the dynamics of the model under consideration with the robust fractional operator called Caputo, and also make some comparisons with the classical variant. The fractionalized order is ψ. As a consequence, many important characteristics of the proposed fractional version of the model have been reported, such as the formation of the model, the existence and uniqueness of the solution through the fixed-point theorem, iterative solutions, and the solutions' positivity. It should be observed that the fractional version of the model under investigation more correctly understands the behavior of the model than the variant of the integer order. This has been done by numerous numerical simulations by using an efficient numerical scheme.  Data Availability Statement: All data used in this work can be found within the manuscript.

Conflicts of Interest:
The authors declare no conflict of interest.

Symbol Description F
The flow rate through the reactor.

S
The substrate concentration within the reactor. S * The dimensionless substrate concentration within the reactor [S * = S/S 0 ]. C The dissolved oxygen concentration . C * - The dimensionless dissolved oxygen concentration C S 0 . K I a Volumetric oxygen mass transfer coefficient C 0 The dissolved oxygen solubility limit C * 0 The dimensionless dissolved oxygen solubility limit K d The death coefficient K * d The dimensionless death coefficient K s The saturation constant S 0 The concentration of substrate flowing into the reactor. V The volume of the reactor.

Y xls
The yield factor for substrate α * The dimensionless yield factor for substrate Y xlo 2 The yield coefficient for oxygen γ * The dimensionless yield coefficient for oxygen. µ m The specific growth rate parameter β The reactor model for the reactor. X The biomass concentration within the reactor X * The dimensionless biomass concentration within the reactor t The time. t * The dimensionless time. µ(S, X) The specific growth rate model. τ The residence time. τ * The dimensionless residence time.