Model-Based Analysis of Feedback Control Strategies in Aerobic Biotrickling Filters for Biogas Desulfurization

: Biotrickling ﬁlters are one of the most widely used biological technologies to perform biogas desulfurization. Their industrial application has been hampered due to the difﬁculty to achieve a robust and reliable operation of this bioreactor. Speciﬁcally, biotrickling ﬁlters process performance is affected mostly by ﬂuctuations in the hydrogen sulﬁde (H 2 S) loading rate due to changes in the gas inlet concentration or in the volumetric gas ﬂowrate. The process can be controlled by means of the regulation of the air ﬂowrate (AFR) to control the oxygen (O 2 ) gas outlet concentration ([O 2 ] out ) and the trickling liquid velocity (TLV) to control the H 2 S gas outlet concentration ([H 2 S] out ). In this work, efforts were placed towards the understanding and development of control strategies in biological H 2 S removal in a biotrickling ﬁlter under aerobic conditions. Classical proportional and proportional-integral feedback controllers were applied in a model of an aerobic biotrickling ﬁlter for biogas desulfurization. Two different control loops were studied: (i) AFR Closed-Loop based on AFR regulation to control the [O 2 ] out , and (ii) TLV Closed-Loop based on TLV regulation to control the [H 2 S] out . AFR regulation span was limited to values so that corresponds to biogas dilution factors that would give a biogas mixture with a minimum methane content in air, far from those values required to obtain an explosive mixture. A minimum TLV of 5.9 m h − 1 was applied to provide the nutrients and moisture to the packed bed and a maximum TLV of 28.3 m h − 1 was set to prevent biotrickling ﬁlter (BTF) ﬂooding. Control loops were evaluated with a stepwise increase from 2000 ppm v until 6000 ppm v and with changes in the biogas ﬂowrate using stepwise increments from 61.5 L h − 1 (EBRT = 118 s) to 184.5 L h − 1 (EBRT = 48.4 s). Controller parameters were determined based on time-integral criteria and simple criteria such as stability and oscillatory controller response. Before implementing the control strategies, two different mass transfer correlations were evaluated to study the effect of the manipulable variables. Open-loop behavior was also studied to determine the impact of control strategies on process performance variables such as removal efﬁciency, sulfate and sulfur selectivity, and oxygen consumption. AFR regulation efﬁciently controlled [O 2 ] out ; however, the impact on process performance parameters was not as great as when TLV was regulated to control [H 2 S] out . This model-based analysis provided valuable information about the controllability limits of each strategy and the impact that each strategy can have on the process performance.


Introduction
Biogas is a valuable renewable energy source produced during organic matter fermentation in anaerobic digestors in wastewater treatment plants (WWTPs). At the top of biogas utilization pathways are the use of biogas as vehicle fuel, the production of chemicals and the biogas upgrading and injection in the natural gas grids, while sole heat production and combined heat and power production (CHP) are situated lower in the value chain of biogas utilization [1,2]. Biogas utilization in WWTPs would allow supplying, partially or entirely, their electricity and heating needs, turning these facilities into "net zero" energy consumers. Selection of the best biogas utilization pathway depends strongly on the methane (CH 4 ) composition, but also on the composition of other impurities such as hydrogen sulfide (H 2 S) (0.1-2%), which must be removed to lengthen the lifespan of cogeneration engines [3,4]. According to the Biogas Utilization Handbook [4], engines manufacturers recommend limiting H 2 S below 10 ppm v in the biogas stream, while other references recommend value up to 1000 ppm v [5]. The exact value will be dictated by the specific characteristics of the engine manufacturer and operating conditions of the engine.
H 2 S removal in biogas streams can be done by using physical-chemical technologies or by using biological technologies such as biotrickling filters [6]. According to Cox and Deshusses [7], biotrickling filters are biological scrubbers in which a polluted air stream is passed through a packed bed on which a mixed culture of pollutant-degrading organisms is naturally immobilized. Biotrickling filters have shown to be environmentally friendly and an economically and technically efficient strategy to perform biogas desulfurization [8]. Several works in the literature have shown that this technology performs satisfactorily at different H 2 S loading rate (H 2 S-LR) conditions, pH conditions, gas-liquid flow patterns, different packing materials, different oxygen (O 2 ) mass transfer devices as well as different electron acceptors [3,[9][10][11]. In desulfurizing biotrickling filters, the H 2 S contained in the biogas can be converted under aerobic or anoxic conditions into either elemental sulfur (S 0 , an inorganic solid) or sulfate (SO 2− 4 ) by S-oxidizing bacteria (SOB), depending on the S/O ratio. Equations (1) and (2) show the H 2 S oxidation reaction under aerobic conditions [12]: Experimental results show that according to Equations (1) and (2), if the S/O ratio is higher than 2, S 0 is the main oxidation product, while if the S/O ratio is lower than 0.5 the reaction continues to SO 2− 4 formation [13]. Due to the packed-bed nature of biotrickling filters, bed clogging caused by S 0 accumulation in the biotrickling filter (BTF) bed is the main drawback related with the industrial application of aerobic desulfurization. This situation it is most likely to occur when biotrickling filters are operated under variable H 2 S-LR at limited dissolved oxygen (DO) concentrations [10,11].
Therefore, it is of great relevance to avoid oxygen limitations by optimizing the Gas-Liquid O 2 mass transfer in desulfurizing biotrickling filters to maximize SO 2− 4 formation and to avoid S 0 formation in the bed but limiting the biogas dilution. This can be achieved by applying control strategies to have a robust and reliable desulfurization process regardless of any existing disturbances, and to prevent environmental and technical problems when unexpected H 2 S-LR variations occur. In full-scale desulfurizing biotrickling filters, the prevailing H 2 S-LR may vary due to varying inlet H 2 S concentrations and/or biogas flowrates. The increase in H 2 S concentration may happen when the organic sulfur content of the feedstock to an anaerobic digestor increases [14,15], while the biogas flowrate fed to the trickling filter may vary, in facilities with small gasholders for biogas storage or when large amounts of biogas are produced due to a peak in the feedstock to the anaerobic digestor [16]. The most common solution in industrial facilities to avoid excessive S 0 accumulation due to the fluctuation of H 2 S-LR is to provide additional O 2 by increasing the air flowrate [3]. Nevertheless, this strategy causes biogas dilution and leads to the subsequent loss of its volumetric energy content. For this reason, strategies based on the improvement of O 2 transfer from the gas to the liquid phase are key to optimize the complete oxidation of H 2 S to sulfate, without diminishing biogas calorific power. Lopez et al. [11] found that the regulation of the trickling liquid velocity (TLV) improved sulfate selectivity [11]. TLV does not only increases the dissolved oxygen load but also improves the penetrability of the liquid through the BTF bed, reducing DO gradients along the bed, without diminishing the biogas calorific power.
Previous work in the literature have focused on the development and implementation of feedforward control and feedback control strategies for lab-scale H 2 S desulfurizing biotrickling filters under aerobic and anoxic conditions [17,18]. However, modeling work on control strategies are scarce in the literature. Coupling modeling and simulation has great importance, since the simulation of control strategies for desulfurizing biotrickling filters allows testing a great number of operating conditions and different types of control strategies in a very short time. This type of approach can be very helpful when there are a various number of possible combinations of control strategies derived from the combination of the different process disturbances, manipulated variables, type of control strategies and type of controllers applied. Recent works in literature have developed mathematical models for biotrickling filters for biogas desulfurization under aerobic conditions [18][19][20] and also under anoxic conditions [17,21]. Such models are of great importance towards the development of control strategies to optimize the performance of biotrickling filters for biogas desulfurization. Hence, this work evaluates different feedback control loops to optimize the performance of a biogas desulfurizing biotrickling filter under different H 2 S-LR conditions. This work performs a theoretical analysis of feedback controllers and control strategies through the analysis of control limits, capabilities and added value of these strategies applied in an aerobic BTF for biogas desulfurization. Feedback control was selected over feedforward control since it is the most common strategy used in industrial processes owing to its simplicity. Additionally, feedforward control was experimentally tested in a previous work for the same application and unsatisfactory results supported this decision [18]. A feedback controller was incorporated into a mathematical model for a BTF for biogas desulfurization under aerobic conditions that was calibrated and validated under steady and dynamic conditions [22]. An analytical approach using time-integral criteria was used to determine the most suitable control parameters. Control strategies tested in this work derive from common actuations in industrial biotrickling filters as well as from experimental results from the influence of the main manipulated variables in aerobic desulfurizing biotrickling filters.

Materials and Methods
The model used in this work to simulate an aerobic biotrickling filter for biogas desulfurization is based on the work done by Lopez et al. [22]. Model development, model equations, model calibration and model validation procedure can be found elsewhere [22]. In short, the model used in this work incorporates expressions for the different mechanisms occurring in the packed bed such as mass transport by advective flow in the gas and liquid phases, mass transfer at the gas-liquid interface, mass transfer by diffusion at the liquid biofilm interface, internal diffusion in the biofilm and biological reaction in the biofilm. The model includes time-dependent differential equations describing the mass balances for the main species involved in the H 2 S oxidation process, such as hydrogen sulfide, sulfide, elemental sulfur, oxygen, sulfate and thiosulfate.

Liquid-Gas Transfer Correlations
Gas-liquid mass transfer is an important step to consider when describing the mechanisms of aerobic biogas desulfurizing in biotrickling filters. Experiments and modeling works have shown that the main process variables such as the removal efficiency (RE), S 0 and SO 2− 4 production are highly sensitive to the variation of O 2 and H 2 S gas-liquid mass transfer coefficients K L,O 2 , K L,H 2 S in biotrickling filters for biogas desulfurization [22].
The BTF model used in this work [22] was calibrated and validated at constant TLV conditions, thus the gas-liquid mass transfer coefficient was defined as a constant value for O 2 and H 2 S K L,O 2 , K L,H 2 S . Therefore, to have a better representation of the mass transfer phenomena between the gas and liquid phase, the effect of the variation of the main manipulated variables, such as the gas and liquid flowrate must be included in the estimation of K L,O 2 and K L,H 2 S . Lopez et al. [22] demonstrated that the contribution of the gas-side mass transfer resistance was negligible, and that only the liquid-side mass transfer resistance was significant. Thus, the effect of air flowrate over the K L,O 2 and K L,H 2 S values can be neglected, and focus can be placed only on the effect of TLV over mass transfer coefficients.
In this study, two empirical mass transfer correlations for packed columns, the Billet and Schultes and Onda et al. correlations [23,24] respectively, were compared to select the most suitable correlation for describing gas-liquid mass transfer in our system. Billet and Schultes and Onda's correlations are described by Equations (3) and (4): where:  (3) and (4) TLV is named differently. In Equation (3) it is shown as u L , the superficial liquid velocity, while in Equation (4) as L, the superficial mass velocity of liquid, which can be related to TLV through the liquid density and the cross-sectional area of the bed (A). Operationally, TLV is regulated through the liquid recirculation flowrate in a BTF. Equation (3) is named after the work published by Billet and Schultes [24] to predict mass transfer in packed beds, although it does not appear in their original publication. Equation (3) is the result of algebraic manipulation published in [25] of the mass transfer equations of the Billet and Schultes.
The Billet and Schultes and Onda correlations were first compared by evaluating typical TLV values for BTFs ranging between (5-30 m h −1 ). The results obtained here were compared with the K L,H 2 S obtained during model calibration. Later, the effect of the mass transfer correlation on the BTF model performance was evaluated in terms of RE and S 0 and SO 2− 4 production. The effect of the mass transfer correlation was assessed by running a set of simulations applying the operational conditions described in Table 1. Conditions presented in Table 1 are the same as those presented for the constant TLV experiment in López et al. [11]. The effect of TLV over K L values and secondly for a fixed TLV value was tested, in order to assess the effect of each correlation on the BTF performance in terms of RE and S 0 and SO 2− 4 production. Finally, the best correlation was added to the model. The main biotrickling filter model used in this work includes multiple expressions for the different mechanisms occurring in the packed bed and the incorporation of this correlation aimed to describe better the gas-liquid mass transfer phenomenon in our system. Further information about the BTF model used in this work can be found in [22].

Open-Loop Study
To characterize the response of the BTF when a disturbance occurs under Open-Loop conditions, operational conditions were varied from the reference conditions as shown in Table 2. Further information about the reference conditions can be found in López et al. [11]. Mainly, focus was placed on studying the effect of two of the most common disturbances occurring in the H 2 S-LR in full-scale desulfurizing BTFs (Table 2)  For both open-loop scenarios (OL-C and OL-EBRT), simulations were run sufficiently long to ensure steady state conditions. Sulfur and O 2 mass balances were performed to complete the evaluation of the system under. Information related to model performance obtained during the Open-Loop simulations was useful to delimit and define the operational conditions for the Closed-Loop tests.

Closed-Loop Study
The design of the Closed-Loop framework analysis was based on the results obtained in the Open-Loop framework. The flowchart shown in Figure A1 (See Appendix A) illustrates the proposed methodology followed for the design of control strategies for an aerobic desulfurizing BTF.

Definition of Closed-Loop Scenarios
Two case studies were designed to evaluate feedback control strategies in an aerobic BTF for biogas desulfurization. The first strategy studied, namely air flow regulation Closed-Loop (AFR-CL), is the traditional control strategy applied in desulfurizing Biotrickling filters [3]. Figure 1a presents a detailed description of the AFR-CL strategy. This strategy controls the O 2 outlet concentration ([O 2 ] out ) through air flowrate regulation when H 2 S-LR is varied due to H 2 S inlet concentration or biogas flowrate changes. In AFR-CL, the reference H 2 S-LR (56.3 g S-H 2 S m −3 h −1 ) increases up to a H 2 S-LR of 168.9 g S-H 2 S m −3 h −1 . This strategy aims at regulating the oxygen concentration in the biotrickling filter to avoid S 0 accumulation in the packed bed, and therefore the packed bed clogging [26]. If S 0 accumulates in the packed bed the overall performance of the biotrickling filter will decrease due to a decrease in the sulfide consumption rate and subsequent reduction in the H 2 S removal. The H 2 S inlet concentration was stepwise increased from 2000 ppm v to 6000 ppm v ( Figure 2). As in the Open-Loop analysis, each concentration was evaluated long enough to reach steady-state conditions as well as to study the long-term effect of disturbances in the system. In AFR-CL, an O 2 outlet concentration set point (SP) of 5.7% (v/v) was established to calculate the error. The SP value for the AFR-CL was selected according to the O 2 outlet concentration achieved under reference conditions when a 100% H 2 S RE was achieved. The aim of this strategy is to provide the O 2 required accordingly with the H 2 S-LR changes occurring in the BTF. With regards to the air flowrate span, the regulation range for this variable was limited between the reference air flowrate of 24.2 L h −1 and to a maximum value of 120.8 L h −1 that corresponds to biogas dilution factors of 28.2 and 66.3%, respectively. Considering the methane average composition is in the range of 50 to 70% (v/v) [4], a minimum methane content in air of 16.9% would be obtained when the above-mentioned dilution factor is applied. This value is higher than the percentage of methane required to obtain an explosive mixture. It is worthwhile mentioning that such air flowrates and their corresponding biogas dilutions are nonsense from a practical application perspective. However, the analysis from a theoretical perspective was considered as a worthwhile exercise to be done to quantify the effects of such strategy on a desulfurizing BTF. The second strategy studied, namely TLV regulation Closed-Loop (TLV-CL), is an innovative strategy to control biogas desulfurization in biotrickling filters. Figure 1b presents a detailed description of the TLV-CL control strategy, which aims at controlling the H 2 S gas outlet concentration through the regulation of the liquid recirculation flowrate, even under disturbances of H 2 S-LR due to changes of H 2 S inlet concentration or biogas flowrate. In TLV-CL, the final control element is the recycling pump (element 10 in Figure 1a,b), which regulates the liquid recirculation flowrate, and it is related to TLV by dividing this flowrate by the cross-sectional area of the BTF. The H 2 S-LR (56.3 g S-H 2 S m −3 h −1 ) was increased up to 168.9 g S-H 2 S m −3 h −1 firstly due to H 2 S inlet concentration and secondly due to biogas flowrate changes.
The effect of changes in the H 2 S inlet concentration was studied with a stepwise increase from 2000 ppm v until 6000 ppm v as shown in Figure 2a. The effect of changes in the biogas flowrate was assessed performing stepwise increments from 61.5 L h −1 (EBRT = 118 s) to 184.5 L h −1 (EBRT = 48.4 s) as shown in Figure 2b. An H 2 S outlet concentration SP of 100 ppmv of H 2 S was established to calculate the error as described in Equation (5). This SP was selected in order to produce a high-quality biogas effluent to enhance the life-time of co-generation engines or the engines and machinery where biogas is finally burned [27]. Regarding TLV regulation limits, a minimum TLV of 5.9 m h −1 was applied when the H 2 S outlet concentration was below the SP to provide the nutrients and moisture to the packed bed. A maximum TLV of 28.3 m h −1 was set to prevent BTF flooding. H 2 S-LR conditions presented in this section for the two control strategies AFR-CL and TLV-CL were also tested under open-loop conditions in order to quantify process improvement when the control loop was implemented. Table 3 shows the different control loops that can be designed when different controlled variables and disturbances are combined with different types of feedback controllers. Note that in Table 3 no derivative action was included, hence only proportional (P) and proportional integral (PI) feedback controller were evaluated. During preliminary tests, it was observed that due to the characteristics of the disturbances tested, the anticipative effect was null since the response of the process was a constant nonzero error, giving no action for the derivative term (dε/dt = 0, where ε is defined as the error between the measured value and the setpoint value). Where K c and τ I in Table 3 refer to the controller gain and the integral time for a feedback controller.

Controller Tuning
In order to tune the parameters of the feedback controller, time-integral performance criteria were used to evaluate the complete Closed-Loop response [28]. Time-integral criteria are based on the entire response of the process, unlike the simple criteria that use only isolated characteristics of the dynamic response (e.g., decay ratio, settling time, overshoot) [29]. Integral of square error (ISE), integral of the absolute value of the error (IAE) and integral of the time-weighted absolute error (ITAE) were determined as described in Equations (5)-(7).
Time-integral criteria were calculated for each Closed-Loop defined in Table 3. Although ISE, IAE and ITAE were determined for all the cases here studied, only results using IAE criteria are presented. Initially, alternative tuning methodologies were evaluated such as Cohen and Coon method using the process reaction curve method and Ziegler and Nichols open-loop response. These methodologies were discarded to be used as final tuning methodology to obtain the controller parameters. However, preliminary results obtained from these methods served to obtain the initial tuning parameters. Results obtained using ISE and ITAE resulted in similar parameter selection as IAE. Since the parameters that provide the lower time-integral criteria may lead to process instability, production of oscillations in the process response was also considered to assess the controller response to certain parameters. Process instability was used to define the parameters ranges as shown in Table 3. Parameters out of this range produced unstable or oscillatory response (data not shown). If the value of the controller parameters is analyzed, it can be observed that lower gains were needed to control the O 2 outlet concentration if compared to the gains obtained when H 2 S was controlled through TLV regulation, denoting that the process was sensitive to air flowrate changes. Although the range of the controlled variables is not the same, and the magnitude of the controller gains cannot be fairly compared, this analysis denotes that the process was less sensitive to H 2 S outlet concentration regulation rather O 2 outlet concentration. Such differences are due to the gas and liquid phase dynamics of the system since, to control the O 2 outlet gas phase concentration, a gas phase variable was manipulated while, to control H 2 S outlet gas phase concentration, a liquid phase variable was manipulated.

Added Value of Feedback Control Strategies on Biogas Desulfurization in Aerobic Biotrickling Filters
The impact of control parameters over process performance was also analyzed. It helped to have an additional criterion to select the most suitable control strategy since besides evaluating if the control objective was achieved (e.g., reaching the SP value), it is also important to evaluate if such strategy also improved process performance. In this sense, it is important to evaluate if a proper H 2 S removal was achieved and if whether elemental sulfur or sulfate was the main product from H 2 S oxidation, in order to ensure process stability in the long run [26]. Sulfur species mass balances referring to the total amount of H 2 S at the inlet and the amount of O 2 consumed during each simulation were calculated to evaluate the added value of the different control strategies. Mass balances using each control strategy were compared with mass balances under Open-Loop conditions. BTF performance at each H 2 S-LR tested was compared between Open-Loop conditions and Closed-Loop conditions using only the data for each LR step. Equations to determine sulfur mass balances and the amount of O 2 consumed are presented in Appendix A.2.

Effect of Mass Transfer Correlation on BTF Model Performance
Results from the comparison of the k L,H 2 S as a function of TLV of the two correlations are showed in Figure 3A, while Figure 3B shows the model predictions using both correlations studied, Billet and Schultes and Onda's correlation, together with experimental data profiles of the BTF performance operated under Open Loop conditions (see Table 1).  Figure 3B was obtained can be found in [11]. Figure 3B shows that model predictions using the Billet and Schultes correlation for the prediction of K L,H 2 S fit the experimental data better than when using Onda correlation. The Billet and Schultes correlation was also preferred by Van der Heyden et al. [30], during the evaluation and comparison of the most frequently used correlations for randomly or structured packed columns. The goodness of the model was evaluated using the Nash and Sutcliffe efficiency (NSE) criteria [31]. NSE is defined as one minus the sum of the absolute squared differences between the predicted (y e ) and observed data (y m ), normalized by the variance of the observed data during the period under investigation according to Equation (8). Essentially, the closer the model efficiency is to 1, the more accurate the model is.
NSE values for both correlations are shown in Table 4, showing that higher efficiencies values are achieved with the Billet and Schultes correlation, confirming that this correlation fits better to the experimental data for all the variables considered. The negative NSE value for Onda's correlation means that the average value of the experimental data was a better predictor than the model. According to these results, the Billet and Schultes correlation was chosen as the one to be used for simulating control strategies with the BTF model.

Analysis of Open-Loop Behavior: Effect of H 2 S LR Changes Due to H 2 S Inlet Concentration and Biogas Flowrate Changes on BTF Performance
Steady-state results of the BTF performance operated under Open-Loop conditions are presented in Figure 4a,b for H 2 S-LR disturbances due to stepwise changes in the H 2 S inlet concentration ( Figure 4a) and biogas flowrate (Figure 4b). Results presented in Figure 4a indicate that H 2 S RE decreased from 98.5% to 42.6%, sulfate selectivity from 99.80% to 80.30%, while elemental sulfur produced increased from 0.10% to 17.80%. When H 2 S-LR is increased by means of increasing the H 2 S inlet concentration, the O 2 /H 2 S volumetric ratio decreases, favoring the accumulation of elemental sulfur, and affecting significantly the sulfide consumption rate [32]. Thus, the BTF is operating under biological limiting conditions. If the condition causing this problem persists, sulfur tends to accumulate in the liquid phase decreasing the H 2 S concentrations gradient between the gas and liquid phase and affecting the RE of H 2 S as it is observed in the data in Figure 4a  On the other hand, when the H 2 S-LR was increased due to the increase of the biogas flowrate (Figure 4b), the loss in the system performance was even more dramatic. H 2 S RE decreased from 98.52 to 18.67%, while sulfate selectivity decreased from 99.8 to 75.5% and elemental sulfur accumulation increased from 0.1 to 23.5% for the lowest and highest H 2 S inlet concentration tested (10,000 ppm v ). Such reduction is directly related to a reduction in the contact time between the gaseous compounds (O 2 and H 2 S) and the liquid phase, which concomitantly affects the mass transfer capacity of the system and, therefore, the overall performance of the BTF [6,9]. This is reflected in the steep slopes of the H 2 S RE profile in Figure 4b. Results of the BTF operated under Open-Loop conditions show that when a disturbance affects the BTF, such as an H 2 S concentration increase or an increase in the biogas flowrate, the main process performance variables are also significantly affected; therefore, process control needs are clearly evidenced.

Controller Tuning
In feedback control, the controlled variables are compared to the setpoint to calculate the error upon which the controller applies an action into the system through the manipulated variable to minimize the error value. It is helpful, before running the system with the feedback controller coupled, to optimize the controller parameters. This step is called controller tuning and it can be done using mathematical equations or algorithms [17]. According to the flowchart of the methodology followed to design the control strategies ( Figure A1), the step after the definition of the closed-loop scenarios is the controller parameter tuning. For each closed-loop scenario, AFR-CL and TLV-CL, and for each type of disturbance of the system, a different set of controller parameters were analyzed. Table 5 shows the results for the IAE criteria obtained by applying controller parameters defined in Table 3, for a proportional (P) and proportional-integral (PI) feedback controller coupled to the BTF model to control H 2 S inlet concentration disturbances through the regulation of air flowrate (AFR-CL).  For a P controller, K c values higher than 0.12 (shaded values in gray) provided an oscillatory response of the final control element (air flowrate in the case of the AFR-CL). K c of 0.1 (in bold) was selected as a suitable value for the P controller due to the lower ISE value and the lack of oscillatory response. This was the best tuning for the scenario with air flowrate regulation (AFR-CL) during variable H 2 S-LR conditions with increased H 2 S concentration. On the other hand, for the evaluation of IAE criteria for the PI controller, three different time integral values (1, 2.5 and 5) were evaluated as shown in Table 5. In order to reduce the degrees of freedom during the controller tuning, once the suitable K c value was determined for the P controller, a 26.5% of reduction for the K c value for the PI controller was applied, as is done for Cohen and Coon reaction curve methodology [28]. In addition, K c values generating an unstable response of the system for a P controller were not included in the study of the PI controller. According to Table 5, the most suitable combination of parameters was K c = 0.074 and τ I = 1. In Figure 5a, the PI controller actuation applied to control the O 2 concentration ([O 2 ] out ) through AFR-CL with the selected parameters (K c = 0.074 and τ I = 1) is presented. While Figure 5b shows the offset between the [O 2 ] out and the SP. Although the [O 2 ] out signal is not exactly at the setpoint value, the difference between the simulated value and the SP was minimal. The average [O 2 ] out was 5.707 while the SP value was 5.70, which gives a relative error of 0.12%. Hence, it was considered that there was no offset between the signal and the SP.

Closed-Loop
The AFR-CL response of the controller to the H 2 S concentration increase is to increase the air flowrate as reflected in Figure 5a to maintain the [O 2 ] out in the SP value as shown in Figure 5b. Figure 5b also shows the [O 2 ] out of the process when operated under Open-Loop conditions. The increase in air flowrate resulted in an increase of the dilution by 2.3%. Figure 6 shows the H 2 S outlet concentration when [O 2 ] out was the controlled variable in the AFR-CL control strategy. The impact of AFR-CL on H 2 S outlet concentration is minimal and it is mainly related to biogas dilution. However, if the controller performance is analyzed, the PI controller resulted useful to control the [O 2 ] out . To have a complete analysis of the impact of the control strategy over the process, in Section 3.4 it will be assessed and discussed the added value of control strategies over the system performance in terms of sulfate selectivity, oxygen consumption and reduction of elemental sulfur production.

•
Disturbance: H 2 S inlet concentration Table 6 presents the results for IAE criteria for the P and PI controllers for the TLV-CL scenario when the H 2 S inlet concentration is the disturbance. For the P controller, it can be observed how a K c higher than 300 barely reduced the error criteria. Moreover, the higher the gain the more unstable response was obtained (data not shown). For example, a gain of 1000 produced an unstable response of the actuator (shaded values in gray). Therefore, following the IAE criteria and the simple performance criteria, a K c in the range of 200 and 300 was selected as a suitable range for the P controller in the present case. On the other hand, the evaluation of a PI controller for the TLV-CL shows little differences between the different cases studied. Overall, the smaller error was obtained with the highest K c and the lowest τ I . However, considering only the IAE criteria would be erroneous since, as mentioned before parameters that provide the smaller error can lead to an unstable/oscillatory response of the controller, as in the case of the values shaded in gray in Table 6. Hence, considering both criteria, IAE error and unstable/oscillatory response, the parameters selected were K c = 147 and τ I = 1. Figure 7 presents simulation results for this tuning for the PI controller in the scenario of H 2 S outlet concentration control facing H 2 S inlet concentration changes through TLV regulation (TLV-CL).  Figure 7a shows that during the first two steps of H 2 S concentration, no change in the controller output (TLV regulation) was performed, since the [H 2 S] out was still below the SP. Only during the last concentration step did the controller increase the TLV to reduce the offset between the controlled variable and the SP. It must be noted in Figure 7a that the controller output (TLV) did not reached its maximum as it would be expected. For a digital PI controller, the manipulated variable is expected to continuously increase until reaching its maximum value or until decreasing the controlled variable. In this case that did not occurred, as observed in Figure 7a for the controller output (TLV) and in Figure 7b for the controlled variable ([H 2 S] out ). As it can be observed in Figure 7a, at the beginning of the third H 2 S inlet concentration increase (t = 800 h), corresponding to an inlet concentration of 6000 ppm v H 2 S, the first controller response is to increase the TLV as much as possible according to what is expected for a digital PI controller. However, the controller was not capable of decreasing the error through that actuation after some hours and some overshoot and oscillatory response was produced. During that time the proportional action turned into a null action, since the increase of the error was constant (the increase in the error is not significant). Although the TLV was progressively increased during the last period, the I action by itself was not capable of decreasing the [H 2 S] out . As shown in Table 6, the tuning still provided an IAE error for the PI controller, as it also can be seen in Figure 7b by an offset of 20 ppm v . In addition, Open-Loop [H 2 S] out is presented in Figure 7b to show that the impact of the application of control strategies specifically for the last concentration step is very significant, as will be discussed further in Section 3.
Finally, if the AFR-CL control strategy is compared with the TLV-CL strategy, the latter is based on a more complex process, since, in order to control the [H 2 S] out , first a gas-liquid mass transfer step is needed to reduce the [H 2 S] out , while for the AFR-CL control strategy, the manipulated variable directly affects the controlled variable [O 2 ] out .
• Disturbance: gas flowrate In the case of considering the biogas flowrate as the disturbance for the TLV-CL scenario, results for the IAE criteria are presented in Table 7. As in the TLV-CL case, a different set of parameters might be chosen for both P and PI controllers if only time-integral criteria would be considered and simple criteria such as the stability of the controller response would not be included. In this sense, a K c value for the P controller in the range of 100 to 300 would produce a suitable result for both criteria, while a K c higher than 300 would result in unstable/oscillatory response of the actuator. Table 7. IAE criteria result for the TLV-CL when biogas flowrate is the disturbance of the system. Biogas flowrate was step-wise changed corresponding to an Empty Bed Residence Time (EBRT) decrease from 118 s to 48.4 s. Values in bold represent the selected parameters and values shaded in gray provided an oscillatory response of the final control element. For the PI controller, a K c in the range of 100 to 200 in combination with τ I between 1 and 10 would provide a satisfactory performance of the controller. Figure 8 presents the simulation results for the combination of the selected parameters for the PI controller (K c = 110 and τ I = 1) for the H 2 S outlet concentration control facing biogas flowrate changes (TLV-CL).

Closed-Loop
In contrast to the previous scenario when the H 2 S inlet concentration is the disturbance, when biogas flowrate is modified, the Open-Loop response shows that the [H 2 S] out is already higher just after the first step change. In consequence, the controller increased the TLV to reduce the [H 2 S] out , indicating that the process is much more sensitive to the effect of reducing the EBRT (biogas flowrate increase) that to changes in the H 2 S concentration. Interestingly, at the highest H 2 S LR tested (169 g S-H 2 S m −3 h −1 ), the [H 2 S] out concentration was reduced from 1100 ppm v (Open-Loop response) down to 557 ppm v , which represents a 49.37% of reduction. However, such concentration was still far from the SP (100 ppm v ) indicating that changes in the H 2 S-LR due to a biogas flowrate increase at the highest H 2 S-LR tested were not controllable by means of a TLV-CL control strategy. In addition, another reason why the SP was not achieved during the last H 2 S LR tested was because the maximum TLV was constrained to 28.3 m h −1 . However, for the second H 2 S LR tested (112.6 g S-H 2 S m −3 h −1 ), was reduced to 120 ppm v which is only 20 ppm v above the SP, indicating that under those H 2 S LR condition and EBRT the process could be controlled. It is important noticing that once the last H 2 S-LR was applied, the manipulated variable (TLV) was maxed out to try to reduce the error. Since the maximum TLV was constrained to 28.3 m h −1 , the DO provided to the biotrickling filter was maintained constant although not sufficient to the system. That is reflected in the staircase response of the controlled variable after hour 800 in Figure 8b. The staircase shape is caused by a sum of causes, being the most important ones the progressive S 0 accumulation in the bioreactor, the sulfide accumulation in the liquid phase and the loss of mass transfer capacity under high biogas flowrates (low EBRTs). Under these conditions, S was accumulated in the packed bed since the biotrickling filter was operated under O 2 limiting conditions, since the oxygen supplied by the TLV was not enough to carry out H 2 S complete oxidation to sulfate. During the last biogas flowrate step, 16.10 g of S 0 were produced, compared to only 7.90 g of S 0 produced during the first two steps. The effect of S 0 accumulation as well as the effect of sulfide accumulation over H 2 S consumption by SOB are both included in the kinetic equations of the model used in this work [22].The conditions applied to the biotrickling filter showed that the controller was not able to control the process. The uncontrollability of the process during the last H 2 S-LR was due to a limitation on the gas-liquid mass transfer under the low EBRT values. Therefore, to improve the controllability of the process under such conditions, the gas-liquid mass transfer must be improved, by using intensive gas-liquid mass transport devices such as Jet Venturi [3] to produce an intensive contact between the gas and liquid phase. In this sense, the model-based study of control strategies denotes its usefulness, since results as those shown in Figure 8 warn to the designer that under such H 2 S-LR conditions a redesign of the BTF would be needed to achieve higher H 2 S RE conditions.  Table 8 presents the improvement in the main process variables for the AFR-CL and the TLV-CL strategies facing H 2 S inlet concentration changes using P and PI controllers tuned as discussed in Section 3.3. Interestingly, Table 8 shows that the [O 2 ] out control through AFR regulation reduced the [H 2 S] out and also the elemental sulfur production, which is linked to an increase in the oxygen consumption and consequently to an increase in the sulfate production at the end of the simulation. However, the reduction and improvement percentages for AFR-CL compared to those ones obtained with TLV-CL are much lower. AFR-CL was shown to be less effective than TLV-CL since controlling the [O 2 ] out does not ensure a better gas-liquid mass transfer, which is one of the main limiting steps of the desulfurization process. Through TLV regulation, oxygen consumption was also substantially increased as shown in Table 8 without any further biogas dilution effect if compared to the AFR-CL control strategy. According to Lopez et al. [11], TLV regulation improves the penetrability of the liquid through the packed bed reducing possible DO gradients and, therefore, the formation of elemental sulfur. Regarding the comparison between P and PI controllers, a PI control law performed better in AFR-CL control strategy because the PI controller was able to remove completely the offset by means of a greater AFR actuation, consequently resulting in better performance of the system. However, the P controller provided larger improvements than the PI controller for the TLV-CL control strategy. This difference was due to the different tuning results between both controllers. Independently of the type of feedback controller evaluated, the TLV-CL control strategy offered interesting results to improve the performance of aerobic biogas desulfurization in biotrickling filters. TLV-CL strategy indicated that TLV regulation can help to achieve a reliable operation of the system, without any further biogas dilution. Table 9 shows the improvement of the process variables using P and PI controllers in a TLV-CL control strategy facing biogas flowrate changes.

Selecting the Most Suitable Feedback Controller Using a TLV-CL Control Strategy to Face Variable Biogas Flowrate
As discussed in Section 3.3, the process was out of the controllable limits using TLV as the manipulated variable when the highest flowrate was tested. Therefore, such results were excluded from this analysis. Table 9 shows that despite of operating the system under lower EBRT conditions, RE and oxygen consumption increased due to an enhancement of mass transfer along the packed bed due to TLV regulation. According to empirical correlations described in Section 3.1, a TLV increase also increases the global mass transfer coefficient (K L a) for H 2 S and oxygen, which is here validated by the improvement in RE and oxygen consumption when TLV was increased. The main differences between the two types of controller tested were related to the TLV applied in each case. For the P and PI controller, a maximum TLV of 21.9 m h −1 and 15.9 m h −1 was used, respectively. Thus, greater improvements are achieved with the P controller. It must be highlighted that differences in the performance of the controllers are strongly correlated to the type and duration of the disturbance tested in this study. Nonetheless, irrespective of the type of controller studied, process performance was positively improved through TLV regulation in order to face biogas flowrate changes. However, controlling biogas flowrate changes in BTFs at EBRT values below 69 s can hardly be done using classical controllers as it was done in this work. Table 9. Improvement on the main process performance variables for TLV control strategy facing biogas flowrate changes using a using a P and PI feedback controllers.

Conclusions
The effect of different control strategies applied in an aerobic biotrickling filter for biogas desulfurization under different loading rate conditions was studied. The Billet and Schultes correlation fitted better experimental RE and sulfur species selectivity data than the Onda correlation and was therefore implemented in the model to achieve a better description of the manipulable variables (liquid and gas flowrate) during control strategies simulation. Evaluation of aerobic biogas desulfurization in biotrickling filters under Open-Loop conditions determined the strong need to mitigate disturbances such as H 2 S inlet increase or biogas flowrate changes When evaluating the controllability of the different process disturbances, changes on H 2 S inlet concentration have a higher controllability limit than biogas flowrate changes. The latter was found to be hardly controllable with classical controllers at EBRTs below 69 s.
Regarding the type of control loop studied, the Air flowrate closed loop control strategy was inefficient in terms of process improvement when added-value variables were considered because increasing the aeration flowrate does not ensure a proper gasliquid mass transfer. The Air flowrate closed loop control strategy involves a certain biogas dilution, which reduces the calorific power of the biogas produced. Therefore, Air flowrate closed loop control strategy is not recommended in biotrickling filters for biogas desulfurization. On the other hand, the trickling liquid velocity closed loop control strategy resulted as a proper strategy to regulate H 2 S outlet concentration at H 2 S LR below 169 g S-H 2 S m −3 h −1 when changes in the inlet H 2 S concentration occur. Although this strategy was only able to control H 2 S LR changes below 113 g S-H 2 S m −3 h −1 when the biogas flowrate was the disturbance.
Overall, the TLV regulation also was a feasible strategy to improve complementary BTF performance parameters such as RE, sulfate selectivity, and O 2 consumption under variable H 2 S LR conditions, independent of the type of feedback controller (P or PI). TLV showed better results than Air flowrate regulation strategy since no biogas dilution is involved and because the percentage of process improvement compared to the Open-Loop case is much higher. Acknowledgments: To the Department of Chemical, Biological and Environmental Engineering at UAB (Universitat Autònoma de Barcelona). Also special thanks to Javier Guerrero for his continuous support in control strategies and process automation along this work.

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

Appendix A
Appendix A.1. Flowchart of the Methodology Followed for Feedback Control Strategies Evaluation Figure A1. Flowchart of the methodology followed for feedback control strategies evaluation.