Multivariable Tracking Control of a Bioethanol Process under Uncertainties

Bioethanol is one of the most studied alternative fuels nowadays. Due to its production process complexity and the low quality of the mathematical models that describe it, a reliable controller is needed to maximize the fuel production and minimize its environmental impact, even in the presence of uncertainty. Here, a controller for tracking optimal profiles considering model errors and external perturbations is proposed. +is work is an improvement of a previously presented technique. To reduce the earlier mentioned uncertainties’ effect during the fermentation, some tracking error integrators are added in the control action calculation. +is simple modification ensures the tracking error convergence to zero, even in the presence of uncertainties (demonstration available). Different tests are carried out and a performance comparison with the original controller is shown to highlight improvements in the tracking error of up to 98% when integrators are incorporated. Furthermore, a classical PI controller is contrasted with the proposed techniques.


Introduction
Excessive consumption of fossil fuels has led to a depletion of this natural resource and an environmental contamination caused by the greenhouse gases generated in their combustion [1,2]. Some of the most notorious consequences are the global warming, acid rain, and urban smog [3]. For this reason, renewable and nonpolluting alternative fuels have become an investigation focus. One of the most important research lines in this topic is the bioethanol production, which is less toxic, readily biodegradable, and produces less air-borne pollutants than petroleum fuel [4]. Its production process involves a fermentation that transforms carbohydrates into alcohol. Moreover, it can be obtained from different raw materials, which classify it into three categories [4]: (i) first generation, involves feedstock like sugar cane, sugar beet, wheat, sorghum, fruits, corn, potato, rice, sweet potato, or barley [5,6]; (ii) second generation, uses lignocellulosic biomass from wood, straw, grasses, olive mill [7], or pineapple waste [8], among others [9]; and (iii) third generation comes from macro-and microalgae biomass [10].
Furthermore, there are new investigations oriented to obtain bioethanol from fecal waste [1,11]. However, the second generation is the most studied nowadays because it does not generate a competition with the food industry as the first one and is an alternative of final disposal for many types of waste. e bioethanol production process consists of a series of stages; some of them depend on the raw material used: (i) Pretreatment: it allows altering the lignin and hemicellulose structure of lignocellulosic materials to leave cellulose exposed; this stage facilitates the next one and increases its efficiency.
(ii) Hydrolysis: with this treatment, long chains of carbohydrates are reduced into monomeric sugars.
(iii) Fermentation: biological process ensures the transformation of the monomeric sugars into alcohol molecules by some microorganism (yeast or bacteria).
(iv) Distillation: the alcohol is separated from the culture medium. is method requires a significant amount of energy, which makes the process more expensive.
Many researchers have devoted their efforts to improve each one of the stages mentioned. Herrero, et al. [7] presented the optimization results of the sulfuric acid pretreatment variables applied to improve sugars bioavailability. Aimaretti et al. [12], studied two different enzymatic hydrolysis strategies to increase fermentable sugar concentration in the must. Kumar et al. [13] demonstrated how to increase the efficiency of sugar conversion and the ethanol yield of a fermentation applying gas stripping at high temperature. In Tgarguifa et al.'s study [14], the model and a distillation column optimization to produce ethanol from wine are proposed. e fermentation step can be improved with many strategies, such as the specific microorganism selection or genetic modification of those microorganisms [3,4]. Fedbatch bioprocesses are intensively studied nowadays because of their main advantages [15]. For example, the medium substrate concentration can be regulated by an appropriate feed rate profile [16], obtaining better production yields and minimizing the production costs [15,17]. However, bioprocesses control is required to follow a certain feed flow rate and get stability and the best productivity [18]. Moreover, a real-time bioprocess detailed control is a complex target but is necessary to ensure raw materials optimal use, water and energy-saving, final product consistent quality, a reduction in wastes and process cycle time, replacement of costly and slow laboratory testing, and continuous learning, which opens up the possibility of bioprocess innovation. Furthermore, the mathematical representation of the process is the key to achieve good results.
A mathematical model provides a map from inputs to responses. A model quality depends on how closely its responses match those of the true plant. However, a model set which includes the true physical plant can never be constructed [19]. Generally, a bioprocess modelling presents particular difficulty in their parameters determination caused by the poorly understood microorganism's dynamics (multivariable and highly nonlinear dynamics), the strongly coupled variables, and the presence of numerous external disturbances, which leads to many modelling uncertainties [20]. Furthermore, sometimes those parameters are determined without a previous model analysis, or their values are not informed with their respective confidence intervals [21,22]. Moreover, the time-varying parameters are usually assumed constant [23]. Also, uncertainties associated with processing technologies' parameters are rarely considered [24]. Dismissal of aspects leads to a reality misrepresentation and, as a result, a deficient performance with serious and imminent risks [25,26]. Hence, the main task to assure the bioprocess quality implies finding a way to control these flaws [27][28][29]. For this reason, the development of new control schemes that reduce the effect of uncertainties in the tracking error has become an attention focus for the scientific community [21]. e problem of model parameters determination has motivated the development of many strategies. Donoso-Bravo et al. [30] discussed existing methodologies for parameter estimation and model validation in the field of anaerobic digestion processes. Balsa-Canto and Banga [31] proposed an iterative procedure to improve practical identifiability and to help in the design of informative experiments. Vilas, et al. [32] presented the theoretical background of a parameter identification protocol intended to deal with those mentioned challenges. To reduce the complexity of tuning a complex anaerobic digestion model, a particle swarm optimizationbased smart algorithm was developed to estimate all parameters [33]. Elenchezhiyan and Prakash [34] formulated state estimation schemes for nonlinear hybrid dynamic systems subjected to stochastic state disturbances and random errors in measurements using interacting multiple-model algorithms. In Pantano et al.'s study [35], the problem of optimal profiles tracking control under uncertainties for a fed-batch bioprocess with two control actions is addressed. In Lara-Cisneros et al.'s study [36], an extreme-seeking scheme based on an approach to variable structure control for fed-batch bioreactors, which deals with uncertainty on the specific growth rate, is proposed. For some other examples, see [37] and [38].
On the other hand, many scientists have developed several feedback control strategies to deal with bioprocess uncertainties, such as optimal control, nonlinear model predictive control, fuzzy control, hybrid control, adaptive control, tracking control, and neural network [37,[39][40][41][42][43]. Due to the online implementation difficulty, the high computational cost, the imprecise mathematical models, and online solutions, their use for bioprocesses is limited [17] and is still a research topic [44].
Generally, for disturbances treatment, strategies to reject or limit their effect are used. Among the most frequent approaches can be found prefeeding, disturbance observer, disturbance and uncertainty estimation, or control for active rejection of disturbances [45][46][47][48][49]. Disturbance estimation needs a filter in their structure design (to ensure infinite norm less than one), which increases the mathematical complexity.
For the specific case of ethanol production defined by Hunag, et al. [50], Fernández, et al. [51] presented a controller design focused on looking for control actions, to track predefined concentration profiles. As the controller structure comes from the process mathematical model, it can be implemented in many systems.
is procedure is simple, versatile, and precise. Besides, only basic knowledge of numerical methods and linear algebra are needed to implement it. e technique was tested against different disturbances, showing an excellent performance and the tracking error convergence to zero. However, its weakness is handling model uncertainties and external perturbations.
In this paper, the control strategy presented in Fernández, et al. [51] is improved.
is new alternative allows tracking predefined optimal profiles considering additive uncertainty in an underactuated system with only one control action. To achieve this goal, a new term is added to the original controller design equation. is term represents the model uncertainties, which are avoided using tracking error integrators. It is important to highlight that as the number of integrators increases, the approximation improves and the tracking error decreases. Furthermore, this adjustment ensures signals uniformity and the tracking error convergence to zero. In this way, a positive answer to the tracking control challenging problem in multivariable nonlinear systems with additive uncertainty is provided. e paper is organized in four sections. Section 2 summarizes the ethanol bioprocess and the previously presented controller technique. Section 3 explains the new controller design, which takes into account the uncertainties. Section 4 shows simulations testing the new controller performance and a comparison with other control techniques. Finally, conclusions are exposed.

Ethanol Bioprocess Description.
e following mathematical model describes a fed-batch ethanol fermentation. It was originally proposed by Hunag, et al. [50]. Here, the yeast Saccharomyces diastaticus is used. e system input, feed rate (U), is a 50-50 glucose and fructose mixture, while the states are cells (X), ethanol (P 1 ), glycerol (P 2 ), glucose (S 1 ), and fructose (S 2 ) concentration inside the reactor, whose profiles variation over time is wanted to be tracked: where V is the culture volume, µ 1 and µ 2 are the specific yeast cells growth rates, q S1/P1 and q S2/P1 are the specific ethanol production rates, q S1/P2 and q S2/P2 are the specific glycerol production rates, in all cases from glucose and fructose, respectively. U is the control variable. Hunag, et al. [50] applied a run-to-run optimization procedure in order to estimate kinetic model parameters, fed glucose concentration, feed strategy, and fermentation time to maximize the ethanol production rate. ose results are used as references in this paper. A reactor with a total working volume of 5 L was used. Temperature was maintained at 35.8°C by controlling the circulation of the cooling water. e airflow was fixed at 1.5 vvm, and the pH was kept at 5.0 by adding 1N NaOH. e biomass concentration was determined spectrophotometrically at 540 nm, and the corresponding dry weight of cells was obtained from a calibration curve. e concentrations of glucose, fructose, ethanol, and glycerol were analyzed with a high-performance liquid chromatograph [50]. Table 1 shows system's initial conditions, while parameters' nomenclature, description, and values are exposed in Table 2.

Controller Design.
e methodology described in [51] aims to find U in order to make the system follow preestablished state variables profiles (references) with minimum error. For the controller design is assumed that reference profiles and states at each sampling instant are known. Figure 1 shows the reference concentration profiles and the feed flow rate. Next, the mentioned methodology is summarized.
First, equation (1) is discretized using numerical methods. Here, Euler is applied due to its simplicity and good results [51].
where σ symbolizes each state variable, σ n is the current σ value measured from the reactor (online), σ n+1 is σ value in the next measurement instant, and T S is the sampling time (0.1 h) [52]. e process lasts 15.7 h (T f ). e state variables in n + 1 are approximated as follows: where σ ref is the reference state variables and k σ represents the controller parameter for the variable σ. For this system, the controller parameters are k X , k P1 , k P2 , k S1 , and k S2 . As expected is the error reduction in each sampling time, it must be fulfilled that 0 ≤ k σ < 1. en, replacing equation (4) in (3),

Mathematical Problems in Engineering
Applying equation (5) in (1), To find U, the equation system (7) must have an exact solution. So, b has to be a linear combination of A columns [53], in other words, A and b must be parallel. One way to satisfy this condition is as follows: where the operation between < > and ||.|| represents the inner product and the vectors norm in R n space, respectively.
At this point, a "sacrificed variable" selection is necessary. is variable, defined as S 1ez , ensures that (7) has an exact solution. For more details on its selection and calculation, see [51]. Finally, U n is obtained using least squares [53].
Yield coefficient for glycerol from glucose 0.5331 Yield coefficient for glycerol from fructose 0.4462

K S1
Saturation coefficient for cell growth on glucose (g/L) 159.75

K S1I
Inhibition coefficient for cell growth on glucose (g/L) 94.233

K P1
Saturation coefficient for cell growth on ethanol (g/L) 238.39

K P1I
Inhibition coefficient for cell growth on ethanol (g/l) 2.7378

K S2
Saturation coefficient for cell growth on fructose (g/L) 0.0726

K S2I
Inhibition coefficient for cell growth on fructose (g/L) 9.0048

K P2
Saturation coefficient for cell growth on glycerol (g/L) 35.958

K P2I
Inhibition coefficient for cell growth on glycerol (g/L) 9.9722 K S1P1 Saturation coefficient for ethanol production on glucose (g/L) 1.3409 k S1P1 Inhibition coefficient for ethanol production on glucose (g/L) 18.612

K S2P1
Saturation coefficient for ethanol production on fructose (g/L) 0.9129

k S2P1
Inhibition coefficient for ethanol production on fructose (g/L) 1000 K S1P2 Saturation coefficient for glycerol production on glucose (g/L) 6.7116 k S1P2 Inhibition coefficient for glycerol production on glucose (g/L) 0.5863

K S2P2
Saturation coefficient for glycerol production on fructose (g/L) 0.4310

k S2P2
Inhibition coefficient for glycerol production on fructose (g/L) Sugar total feed concentration (g/L) 300 Table 1: Initial conditions for ethanol fermentation.

Variable
Initial value X (g/L) 1.5 Mathematical Problems in Engineering

Controller Tuning.
In [51], the Monte Carlo algorithm was used to find the k σ combination that minimizes the accumulated error. is methodology consists in simulating the process N times using random k σ [54]: where δ is the confidence and ε the accuracy. en, two new terms are introduced, "tracking error" and "total error", equations (11) and (12), respectively: where p represents the simulation in progress, p � 1, 2, . . . , N and nT S the measurement instant, nT S � 1, 2, . . ., J. Equation (12) is the function cost to be minimized with the Monte Carlo Algorithm. en, k σ that minimize the total error are selected.
Sequence of the Monte Carlo algorithm is as follows [55]: (1) Define the controller's parameters to be optimized.
(2) Determine the number of simulations to be performed (N). e values of δ and ε are chosen depending on the desired accuracy, for example, δ � 0.01 and ε � 0.005, resulting in N � 1000. Theorem. If the discrete system is given by equations (1), the control action is calculated with equations (9), and k σ takes values between zero and one (0 < k σ < 1), then, the tracking error convergence to zero when n tends to infinity is achieved.

Controller Implementation.
e control system diagram is shown in Figure 2, and the mathematical procedure is summarized in the following steps: Step 1: Define T s and σ ref . Read σ n .
Step 4: Define and calculate the sacrificed variable.

Controller Design under Uncertainties
As a main contribution, the following procedure suggests an important improvement for Fernández et al.'s [51] controller performance.

Controller Design under Uncertainty.
In this section, the uncertainties effect in the tracking error is evaluated. With the aim of representing the additive uncertainty effect in the process mathematical model, the terms E σ,n are added in equation (7): e additive uncertainty (E n ) models several kinds of mismatches as well as perturbed systems. E n might depend on the states variables and the system input. Moreover, considering a real plant like z n+1 � f (z n ,u n ), the additive uncertainty can be expressed as E n � f (z n ,u n ) -f (z n ,u n ), where f (z n ,u n ) is the discrete-time nonlinear system model.
If z and u are assumed to be bounded, and f is Lipschitz [57,58], then E σ,n can be modelled as a bounded uncertainty [59,60].
Taking into account this uncertainty and following the same procedure of Section 2.3, the next final expression is obtained: Comparing (24) with (22), the error nonconvergence to zero is noticed due to E n presence. erefore, the next step is E n approximation to reduce their effect in the tracking error, achieving the convergence to zero.

Integral Action.
Here, a novel approach is presented to reduce the effect of E n in the tracking error by adding some integrators in the control action calculus equation. e advantage is that new terms do not add extra mathematical complexity and maintain the simplicity of the scheme presented by Fernández et al. [51].
As it was explained, E n represents the additive uncertainty, typical of mathematical model formulation errors and the influence of external perturbation. ose uncertainties can not be measured, that is to say, they are unknown, but they can be assumed as a polynomial so that all kinds of possible uncertainties can be considered, which can be constant, linear, quadratic, or other. In this way, by increasing the order of the polynomial, the uncertainty approximation is improved and there is a tendency to reduce the tracking error more and more.
is statement is mathematically demonstrated below. E n differences are defined as

Constant Uncertainty.
If E n is assumed as a constant, then δE n � 0. Denoting generally the error as e σ,n , An integrator could force the tracking error convergence to zero. is integrator is defined as follows: nT S e σ,n (t)dt � U σ,n + e σ,n T S . (27) Consequently, equation (4) can be rewritten as where k σ and l σ are the proportional and integral controller parameters, respectively. Demonstration of error convergence to zero when uncertainties are considered is given below.

References
Control action U n Figure 2: Control system diagram.

Linear Uncertainty.
If E n is supposed as a linear function, then δ 2 E n � 0. Following the same reasoning as in 3.2.1, now two integrators are needed to compensate the added uncertainty. ose integrators are Accordingly, equation (4) can be rewritten as where k σ and l σ1 and l σ2 are the proportional and both integral controller parameters, respectively.
After the same steps of 3.2.1 demonstration, the error convergence to zero is achieved.

Polynomial Uncertainty.
In this case, E n is supposed as an m-order polynomial function. So, if q > m, δ q E n � 0. Using the same steps before mentioned, it should be added m integrators in order to compensate the uncertainties considered. Figure 3 shows the block diagram of the new control system.

Controller Parameters Selection.
For the error to decrease progressively and remain at low values, the linear term of equation (24), L, must tend to zero. In Fernández et al.'s study [51], this state was possible choosing k σ values between zero and one. However, when integrators are considered, two types of parameters are involved, proportional and integral, what changes the parameters selection condition to achieve the error convergence to zero.
In this way, starting from the characteristic equation of the linear term of equation (33), e σ,n � − NL n + NL n+1 , If the characteristic equation roots (r 1 and r 2 ) are between zero and one, it is guaranteed that the linear term tends to zero when n tends to infinity. Hence, the root values must be determined to calculate parameters.
As many authors recommend, the Monte Carlo algorithm is used to find the best roots values [61][62][63]. e steps to follow are similar to those described in 2.3. Considering that only one integrator is added, the process must be simulated N times (equation (10)), using random values of r 1 and r 2 (between zero and one). en, k σ and l σ are obtained with For each simulation, the cost function (equation (12)) is calculated. Finally, the parameters corresponding to the simulation with minimal total error are selected.
When more than one integrator are added, the process is repeated, but now the characteristic equation order increases, and consequently, the expressions for parameters calculation change. For example, using two integrators, the parameters are calculated as follows: Table 3 shows the results of the three controllers tuning, the original controller [51] (C1), a controller with one integrator (C2), and another with two integrators (C3).

Results and Discussion
e next section demonstrates the effectiveness of the improvement proposed through four important tests. In the first one, the operation of the three controllers is evaluated with a step perturbation in the control action. In the second one, a step and a ramp disturbance are added in the control action. In the third, the controllers are evaluated under parametric uncertainty using the Monte Carlo Algorithm. Finally, the controllers are tested against all earlier mentioned perturbations together. In each situation, the better performance of controllers C2 and C3 is evidenced. Furthermore, a comparison with a typical industry controller is shown.

Simulation Adding a Step Perturbation in the Control
Action. In this test, a hypothetical situation that may produce an unexpected variation in the production is simulated. According to this, a − 30% step perturbation in the bioreactor feed rate is added to evaluate the response of the controllers. Table 4 summarizes simulations results. Figure 4 shows the control action variation compared to the reference and the percentage error, considering C1 error as 100%. Note how the Total Error improves when the new algorithm is applied, it is reduced by 89.91% just adding one integrator and by 94.55% with two integrators.

Simulation Adding
Step and Ramp Perturbations in the Control Action. In this test, a ramp disturbance is added to the step perturbation previously presented. Under those conditions, controllers C1, C2, and C3 are evaluated. In Table 5, the Total Errors obtained in the simulations are exposed. Figure 5 shows how the feed rate varies in each case. Besides, the decrease of the Total Error is represented in the bar graph, demonstrating the improvement of the controller when one (91.4%) and two (97.1%) integrators are incorporated.

Parametric Uncertainties.
In all bioprocess, model parameters may change unexpectedly [64]. If this variation is not foreseen, there could be significant instabilities in the process [26]. erefore, to achieve good performance and efficiency, bioreactors require advanced regulation procedures [29].
As Monte Carlo Randomized Algorithm has been used for uncertainty quantification in many applications [65,66], it is used for this test, following the procedure described in 2.3. e number of simulations, N � 1000, is obtained with equation (10), adopting a confidence (δ) of 0.01 and an accuracy (ε) of 0.005. erefore, the following test demonstrates the technique success from a statistical point of view [67][68][69].
In simulation, a way to quantify uncertainties and perturbations is to specify the parameters' real range of variation instead of using a constant value with greater error [70]. Hunag et al. [50] specified, in Table 1 of their work, the parameters' confidence interval for the ethanol process under study. For this test, the information of the mentioned table is used for control purposes. Using C1, 1000 simulations are carried out. In each simulation, all the system parameters are randomly changed by ±10% of their original range value (Table 1 of Hunag et al.'s study [50]). en, the Total Error is calculated with equation (12), and the system parameters for the worst situation are selected. en, C2 and C3 are tested with those system parameters. e results are then compared considering the Total Error (Table 6). Note how, even considering the most problematic situation of the system, the Total Error decreases for each controller, and note the performance improvement for the controllers C2 (51.27%) and C3 (56.56%) over C1. Moreover, Figure 6 illustrates the accumulated error for each simulation and evidences the improvement with a bar graph.
Considering eorem 1 of Tempo and Ishii [54], as the tracking error (equation (11)) for the 1000 simulations remain bounded, the operation of the controllers C2 and C3 will be satisfactory with 99% of probability while the parameters vary within a range of ±10%.

Simulation Adding a Step and a Ramp Perturbation in the Control Action and Parametric Uncertainty.
is last test aims to demonstrate how the controllers can fix an important perturbation. In this way, the perturbations described in 4.1, 4.2, and 4.3 are simultaneously considered. Table 7 compares the resulting Total Errors. Also, Figure 7 shows the feed rate evolution in the three simulations, the accumulated error obtained with each controller, and a visible comparison of the performance improvement by a bar graph, when controllers C2 and C3 are used. In this case,

Linear algebra controller System
Integral action

References
Control action U n X Figure 3: New control system block diagram.         the performance error can be reduced by 87.87% just adding one integrator and by 98.07% with two integrators. To summarize, the addition of integrators in the controller design has advantages, such as performance improvement, no significant mathematical complexity increase, and additive uncertainty consideration, although it has a disadvantage as integrators are added, the equations order increases, making mathematical development less simpler.

Comparison with a Typical
Controller. For their simple structure and easy parameter adjustment, PID or PI controllers are still the most widely used in factories [71], although their weakness with tracking variable set points is well known [72][73][74]. Figure 8 shows a performance comparison between the proposed controllers and a traditional PI in the same ethanol production system. Here, a step and a ramp perturbation in the control action are added. e best PI parameters were selected with a Monte Carlo algorithm. PI accumulated error is much higher than that of C1, due to the increase in error over time.
It is important to note that a rigorous comparison between C1 and PI is detailed in [51].

Conclusions
A previously presented control technique for tracking optimal profiles in nonlinear and multivariable systems was improved. e proposed methodology allows tracking predefined concentration profiles, even in the presence of typical bioprocess uncertainties like unknown nonlinearities or parameters uncertainties. To take into account the mentioned uncertainties, one or more integrators of the tracking error are added in the controller design. e main contribution of this paper is to reduce additive uncertainties influence on the tracking error without increasing the controller complexity. Note that C2 and C3 preserve the original controller advantages (simplicity, fast design, etc.) but manage to reduce uncertainties effect.
e Monte Carlo Randomized Algorithm is used to tune the controllers and carry out the test with parametric uncertainty.
e different tests carried out show the effectiveness of the proposed procedure, reaching an improvement of 98% in the worst case considered.
Comparing with other previous strategies that deal with analogous uncertain control problems, the proposed controller presents the advantage of avoiding the stochastic modelling needed to deal with parameters under perturbation of white noise. Besides, this nonlinear control does not require a great mathematical effort and does not add significant complexity to the original controller.

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

Conflicts of Interest
e authors declare that they have no conflicts of interest.