Downscaling Industrial-Scale Syngas Fermentation to Simulate Frequent and Irregular Dissolved Gas Concentration Shocks

In large-scale syngas fermentation, strong gradients in dissolved gas (CO, H2) concentrations are very likely to occur due to locally varying mass transfer and convection rates. Using Euler-Lagrangian CFD simulations, we analyzed these gradients in an industrial-scale external-loop gas-lift reactor (EL-GLR) for a wide range of biomass concentrations, considering CO inhibition for both CO and H2 uptake. Lifeline analyses showed that micro-organisms are likely to experience frequent (5 to 30 s) oscillations in dissolved gas concentrations with one order of magnitude. From the lifeline analyses, we developed a conceptual scale-down simulator (stirred-tank reactor with varying stirrer speed) to replicate industrial-scale environmental fluctuations at bench scale. The configuration of the scale-down simulator can be adjusted to match a broad range of environmental fluctuations. Our results suggest a preference for industrial operation at high biomass concentrations, as this would strongly reduce inhibitory effects, provide operational flexibility and enhance the product yield. The peaks in dissolved gas concentration were hypothesized to increase the syngas-to-ethanol yield due to the fast uptake mechanisms in C. autoethanogenum. The proposed scale-down simulator can be used to validate such results and to obtain data for parametrizing lumped kinetic metabolic models that describe such short-term responses.


Introduction
Syngas fermentation is nowadays an established process for the conversion of waste gases into chemicals [1,2]. The company LanzaTech successfully commercializes the fermentation of synthesis gas (containing CO, H 2 and CO 2 ) into ethanol, and is currently exploring other products, such as acetone and isopropanol [3]. Although mass transfer limitations have often been accounted as a limiting factor for scale-up of syngas fermentation, such limitations could highly relieved greatly by making products which are bubble coalescence-suppressing, such as ethanol [4].
High product specificity towards ethanol (>90%) is required for successful commercialization [1]. In a process called solventogenesis, Clostridium autoethanogenum, the workhorse of industrial-scale syngas fermentation, produces ethanol from syngas (e.g., using Reaction (1) and (2)), while during acetogenesis, syngas is converted into acetate [5,6]. Solventogenesis can be triggered by low extracellular pH [7], by high extracellular concentrations of acetate [8], or by H 2 supplementation [9]. 6CO + 3H 2 O → C 2 H 6 O + 4CO 2 (1) Since industrial-scale reactors are of considerable size (e.g., 5 m diameter by 25 m height, or~500 m 3 , is not exceptional), the occurrence of spatial gradients is more of a rule than an exception. Usually substrate gradients occur when the characteristic time of reaction τ rxn is significantly lower than the characteristic time of transport, which is related to mixing for substrates in the liquid phase via the circulation time (t c ), and to mass transfer for gaseous-phase substrates (τ MT ) [10,11]. For large-scale syngas fermentation, τ rxn is expected to be much lower (~0.3 s; [12]) than τ MT (around 10-20 s; [4]) and t c (~40 s; Figure S1). Furthermore, both hydrostatic pressure differences (3.5 bar at the bottom vs. 1 bar near the headspace) and gas mole fraction differences (e.g., y CO decreases from 50% to 5% from bottom to top due to consumption) cause a gradient of around factor 35 in saturation concentration, while the volumetric mass transfer coefficient k L a might vary locally due to turbulent fluctuations, differences in bubble size and gas hold-up [4]. We hypothesize that all of this leads to sizeable dissolved CO and H 2 concentration gradients, which might have implications for the syngas fermentation performance.
The impact of such concentration differences on C. autoethanogenum can be studied with Euler-Lagrangian CFD modelling. In this way, environmental changes, for example in substrate concentration, temperature, and shear stress, can be recorded from the perspective of the microbe (the so-called "lifeline") [13][14][15]. The cells are simulated as Lagrangian flowfollowers (particles) and, when they do not interact with the flow or concentration field (one-way coupling), are used for analyzing the environmental fluctuations occurring in the bioreactor [16][17][18]. Such analyses could be used for the development of scale-down simulators [19,20] or to study cell population heterogeneity [21]. Two-way coupling has to be realized when studying the influence of biomass on the flow or concentration fields [22,23]. This method requires the use of a structured metabolic-kinetic model that could be coupled with the CFD model in a computationally viable fashion [15]. Although very detailed genome-scale metabolic models and kinetic ensemble models are currently available for C. autoethanogenum and other acetogens [9,24,25], two-way coupling of these models is currently too computationally intensive for practical application. Development of lessdetailed, yet structured, kinetic models by metabolite lumping [26][27][28] is key in studying the influence of C. autoethanogenum on the flow and concentration fields more accurately.
Previously was studied how C. ljungdahlii would respond to CO gradients in a severely mass transfer limited bubble column reactor [17]. It was hypothesized, based upon Euler-Lagrangian results and in analogy with Escherichia coli, that in such a case transcriptional changes were very likely (>84%) to occur, because long-lasting CO limitations would lead to a maintenance-dominated metabolism. Redox-controlled oscillations in biomass-specific uptake and production rates were observed in C. autoethanogenum [29], within the timescale of hours, while substrate fluctuations in the order of seconds (~t c ) or minutes are expected at the large-scale. With scale-down simulators (e.g., based on a single vessel, multiple vessels such as stirred tank reactors coupled with plug flow reactors, or microfluidics) the physiological cell response on such short term fluctuations could be studied, so that the large-scale environment as experienced by the cell is reproduced at bench scale [10,11,13,15,30]. The obtained metabolite fluctuations can be used for parametrization of the lumped metabolic models [26].
Several scale-down simulators have been developed and used in recent decades, but the requirements of a scale-down simulator for syngas fermentation have not yet been identified. Since there are many unknowns in the scientific literature regarding kinetics and short-term cell responses, the execution of scale-down experiments that are representative of the large-scale behavior is crucial for advancing the syngas fermentation field. In this study, we propose a scale-down simulator to study industrial-scale syngas fermentation at lab-scale. To exemplify the distinctive applicability of the proposed scale-down simulator, a wide range of industrial biomass concentrations were studied, since this is a major determinant for the dissolved gas concentration. The impact of gas (CO 2 ) production on the dissolved gas concentration gradients, and thus possible fluctuations for the microbe, was studied. We used our previous CFD model of an industrial-scale external-loop gas-lift reactor (EL-GLR) and our lab observations to develop and analyses lifelines representative for large-scale syngas fermentation [4,31].

Geometry and Flow Field
As a starting point for the simulations, the 3D reactor geometry and computed flow field of the EL-GLR were used, for which the modelling approach was validated on pilotscale data, as described in [4]. The only change was in the syngas composition, from a 50% CO, 50% N 2 mixture to a 50% CO, 20% H 2 , 30% CO 2 (v/v) mixture [29]. Since the average molar mass of these compositions is similar and the ideal gas law applies, the mass-flow inlet boundary condition of 2.11 kg s −1 was kept the same, as well as the absolute headspace pressure of 101 kPa.
Next to the equations for gas and liquid flow, volume fraction and turbulence, the species equations were solved transiently for both phases to obtain the concentration fields, by implementing user-defined functions for both mass transfer and biological reaction ( Figure 1) in ANSYS Fluent 2021R1. The 3D geometry, based on publicly available pictures of the Shougang-LanzaTech plant, the hydrodynamic model CFD of the EL-GLR and its validation on pilot-scale data, are extensively described in our previous work [4]. The scale-down simulator is based upon a 3 L bioreactor (2 L liquid volume), with varying durations of high and low stirrer speed after the start-up phase. See Table S1 for details of the geometry of this reactor [32]. Created with BioRender.com (accessed on 30 August 2022).

Mass Transfer Model
The mass transfer coefficient k L was computed by taking the maximum value derived from either the Higbie [33] (Equation (3)) or the Lamont-Scott relation [34], the result of the latter corrected for the underestimation of the energy dissipation rate ε by the k-ε model [20] using the pneumatic power input derived from standard correlations [35] and the liquid volume integral of ε (Equation (4)). The maximum k L of species i was used to account for both surface layer renewal mechanisms, since high k L might be obtained in zones with high energy dissipation [36], but transfer in low-turbulent conditions is better approximated by the Higbie relation [37].
In order to obtain a realistic mass transfer rate for industrial-scale syngas fermentation, spherical bubbles with constant diameter (3 mm) were assumed, based upon our previous work [4] (Equation (5)). Since coalescence could be suppressed by the presence of surfaceactive compounds (e.g., ethanol, salts) in syngas-to-ethanol fermentations, small bubbles can be obtained, leading to high mass transfer rates [18,31,38,39]. While our multiphase model accounted for mass loss through interphase mass transfer and gas expansion using the ideal gas law in Fluent's volume fraction equation [40], we acknowledge that bubble coalescence, break-up, shrinkage by consumption, and pressurebased bubble expansion were not considered by assuming a constant bubble size. Although these factors could have potentially improved the accuracy of the gas phase description, we chose to prioritize realistic gas mass transfer rates and to focus on the biological aspects of our study. Therefore, we opted for a simplified set of equations, similar to those in [17], that were within the scope of our work.
The saturation concentration was calculated considering the local gas phase mole fraction. The pH equilibrium of CO 2 with carbonate species could increase the gas-toliquid mass transfer rate in neutral (pH 6-8) and basic conditions (pH > 8); however, this effect can be neglected, as the syngas fermentation process is operated at pH 5. To ensure complete saturation of CO, H 2 , and CO 2 , and achieve steady-state conditions (i.e., statistically stationary) in the average flow and concentration fields, the mass transfer model was run for 1000 s. Although short-time fluctuations occurred, there were no longterm dynamics, as evidenced by the constant rolling average of the global gas hold-up and dissolved gas concentrations.

Biological Reaction Modelling
The biomass-specific CO and H 2 uptake rates (q i ) were modelled using a recently derived kinetic model [12]. The kinetic models for both CO and H 2 uptake account for CO inhibition (Equations (6) and (7)), and are based on the models derived by [41] and [42], respectively.
The overall reaction rate r i is the product of q i and the biomass concentration c X , the latter assumed to be spatio-temporally constant, as a continuous process is considered in steady state. The reaction rates were enabled once CO and H 2 concentrations reached a steady saturation value. Once statistically stationary flow and concentration fields were obtained (after 600 s, the rolling averages of the global dissolved gas concentrations and hold-up remained constant), time-averaged fields were collected over an averaging period of 200 s. Parameters used for computing the Eulerian concentration fields are provided in Table 1.  [12] Inhibition coefficients K I , K I,CO 0.246 0.025 -mol 2 m −6 , mol m −3 [12] The influence of microbial CO 2 production was examined by modelling two extreme cases at 25 g L −1 biomass: (1) only CO 2 consumption by H 2 catabolism, q CO 2 = − 1 3 q H 2 , and (2) also including production by CO catabolism: q CO 2 = 4 6 q CO − 1 3 q H 2 . The case with the most extensive dissolved gas concentration gradient was subsequently used to study a wide range of industrially relevant conditions, by running simulations with varying biomass concentrations (2, 5, 7.5, 10 and 25 g L −1 ).
The obtained dissolved gas concentrations from the Eulerian simulations for the different biomass concentrations were compared to approximations obtained by a simple ideal-mixing model (c L,i from Equation (8)), wherein it was assumed that all transferred gas is directly consumed. The ideal-mixing model used the same parameters (Table 1) and uptake kinetics (Equations (6) and (7)), an average pressure (274 kPa), and the volumeaverage k L a obtained from the CFD simulations.

Lifeline Analysis
Microbial lifelines were obtained for three cases with different biomass concentrations (5, 10 and 25 g L −1 ). Massless Lagrangian particles were injected at the sparger and tracked for a certain number of circulation times N t c and particles N p . To account for turbulence effects, Fluent's discrete random walk model was enabled. While tracking the particles, the current time, position, concentrations and biomass-specific uptake rates were stored in text format every 0.1 s. Data obtained for the first 90 s (approximately one 95% mixing time t m , Figure S1) were discarded to ensure that the particles were evenly dispersed over the whole reactor volume during the entire analysis.
The lifelines revealed distinct periods of maxima and minima (i.e., peaks and valleys) in both dissolved CO and H 2 concentrations. Peak and valley periods were defined from the lifelines by comparing the transient concentrations in the lifeline with the Eulerian average dissolved gas concentration: in case the transient concentration was 2 times (for the 5 g L −1 case) or 1.5 times (for the other cases) higher or lower than the Eulerian average concentration for at least 1 s, then a peak or valley was assigned, for which the residence time and the average gas concentration were stored. Probability-normalized histograms were calculated subsequently using 100 linearly distributed bins over the whole parameter space (e.g., residence time or average concentration in peak), except for time in the valleys, where the maximum value was capped at 150 s. The circulation time t c was calculated as the average time between two peaks (Equation (9)): For the case with 5 g L −1 biomass, lifelines were obtained during t lifeline = 1000 s (around 23 circulation times) and for N p = 160,000 Lagrangian trajectories. This resulted in extensive simulation time and data usage, so that the analysis of the full dataset was computationally unwieldy. We determined how many Lagrangian trajectories (N p ) and circulation times (N t c ) were needed to ensure statistical independence using the Kullback-Leibler divergence (see Supplementary Information, Figures S2 and S3).

Design of a Scale-Down Simulator
The scale-down simulator was designed based on the results of the lifeline analysis (i.e., the probability density functions of concentrations and residence times in peaks and valleys). The goal of this scaled-down system is to reproduce to the best possible degree the residence times and concentrations experienced by microbes in the full-scale system. The starting point was a continuously operated bench-scale stirred tank reactor (CSTR) (see Figure 1 and Table S1), for which operational conditions were varied to mimic the large-scale environment at several biomass concentrations.
Mass transfer, dilution and consumption rates were modelled for CO, H 2 , CO 2 and biomass while assuming ideal mixing in the liquid phase (Equation (10)). The evolution of the gas composition in the dispersed phase y D,i and in the headspace y H,i were also considered (Equations (11) and (12)) since these could be highly variable during operation at low gas flow rates. The dispersed gas volume V G,D was determined by approximating the gas hold-up using the method proposed by [45], while volume balancing was done to calculate the headspace gas volume V G,H .
The volumetric mass transfer coefficient k L a of compound i is dependent on the superficial gas velocity and the stirrer speed [46] and was estimated by considering mass transfer enhancement by the broth composition ( f broth = 1.5; [31]), the temperature and the compound-specific diffusion coefficient in water (Equation (13)). The power input was estimated for a Rushton impeller with P 0 = N Po n 3 d 5 i and the geometry used [32,45].
The overall gas consumption rate was determined using the local concentrations in the liquid phase via Equations (6) and (7) and the biomass concentration. The biomass growth rate µ · c X was determined using the model parameters derived in [47] for solventogenic conditions (Equation (14)), while neglecting the maintenance requirements of the biomass. Biomass retention in the system was assumed (e.g., [48]) and varied by adjusting the biomass recycling rate R rec .
The modelled bench-scale reactor was operated with a constant dilution rate of 0.021 h −1 , inlet gas flow of 0.05 vvm, temperature of 37 • C, pressure of 101 kPa and a stirrer speed during start-up of 75 rpm. Initial concentrations of CO, H 2 , CO 2 and biomass in the liquid (Table 2) were assumed to solve the system with the ode15s function in MAT-LAB. After the start-up period, the concentration oscillations were repeatedly imposed by varying the stirrer speed. The obtained scale-down lifelines were analyzed using the same routine as for the industrial-scale reactor but considering no threshold factor to discriminate between the peaks and valleys, since these were manually imposed. Inlet gas fraction

Influence of Gas Production
The results of Eulerian simulations with 25 g L −1 biomass in two CO 2 production cases were compared in terms of the dissolved CO concentration c L,CO distribution in the reactor (Figure 2a,b). Although the dissolved gas concentrations in both simulations were in the same range (as expected, since the biomass concentration was kept constant), the spatial distribution of c L,CO within the riser was completely different. In both cases, the highest CO concentrations appeared at the base of the riser, where the mass transfer rates are high due to the hydrostatic pressure and high CO and H 2 gas fractions. As the gas rises, the pressure and gas fractions decrease, leading to lower mass transfer rates. More mass transfer was observed in the top separator due to the locally increased gas hold-up ( Figure S4), leading to increased c L,CO . In the downcomer, the long biomass residence time and poor gas renewal caused low CO concentrations. Effect of CO 2 production on mass transfer and dissolved gas distribution. Time-averaged (200 s) dissolved CO concentrations c L,CO in the zy-plane (x = 0) of the EL-GLR for 25 g L −1 biomass (a) without CO 2 production and (b) by considering CO 2 production. Time-dependence of the dispersion volume-averaged (c) ε G and (d) k L,CO a during the CFD-simulations. Until 2200 s, only gas-liquid mass transfer was included (black line). At 2200 s, gas consumption was switched on in the model in three cases: 2 g L −1 biomass without CO 2 production (blue line), 25 g L −1 biomass without CO 2 production (red line), 25 g L −1 biomass including CO 2 production (green line).
The gas plume is pushed towards the left side by the liquid exiting the downcomer, causing high dissolved gas concentrations at the left side in the case without CO 2 production, which is due to reduced oscillations in the gas plume. When CO 2 is considered, the gas, and thus the dissolved syngas, concentrate towards the middle (Figure 2a,b), in a similar manner to the case without gas consumption [4]. Additional gas is generated halfway along the riser by microbial reaction and is transferred back to the gas phase due to CO 2 oversaturation at decreased hydrostatic pressure. The additional gas in the riser (cf. Figure S4a,b) leads to transport dissolved CO towards the right side of the riser (cf. Figure 2a,b) and homogenizes the dissolved gas distribution (i.e., the variation of c L,CO and c L,H2 ) within the whole reactor volume ( Figure S5).
Next to c L,CO , the evolution of gas hold-up ε G and consequently k L a in the EL-GLR are highly affected by gas consumption (Figure 2c,d). As the mass transfer simulation starts without dissolved gas at t = 1200 s, there is an initial drop in ε G due to gas dissolution. The lower ε G causes a drop in k L a because of their linear dependence (Equation (5)). After about 400 s, the liquid saturates with dissolved gas and ε G and k L a stabilize. When the reaction is switched on at t = 2200 s, both ε G and k L a suffer a significant drop in cases with high biomass concentration, since high amounts of gas are being consumed (Figure 2c), even when CO 2 production is included. Similar decreases in ε G were also visible in the model in [17]. Interestingly, with little biomass in the reactor (2 g L −1 ), the gas conversion decreases significantly (from 0.67 kg s −1 at 25 g L −1 to 0.16 kg s −1 ) due to inhibiting CO concentrations (see Section 3.1.2 and Figure S6), and increasing ε G and k L a, compared to the cases with more biomass.
Although ε G could be well predicted with empirical relations in cases without gas consumption [4], the 33% decrease in ε G by microbial gas consumption makes the prediction of ε G , and thus k L a, even more challenging in operational EL-GLRs. This observation is especially relevant for gases rich in carbon source or electron donors, like the used syngas, in contrast to air, where the dilution with inert N 2 , and typically near equimolar conversion of O 2 into CO 2 results in negligible volume changes due to mass transfer.
The reduced gradients when considering CO 2 production make the c L,CO variations less impactful for the micro-organisms. Due to uncertainties regarding the metabolism, e.g., the possibility of simultaneous CO and H 2 consumption [49], the modelled cases would either under-or overestimate the CO 2 production rate. Since the case without CO 2 production appears to generate larger fluctuations and thus complicate the design of the scale-down simulator (and this is, dependent on the syngas composition, the ideal gas fermentation process from an environmental point of view), we chose to further examine this scenario.

Influence of Biomass Concentration
Dissolved gas concentration fields in the large-scale reactor were computed with 2, 5, 7.5, 10, and 25 g L −1 biomass, as shown in Figures S5 and S6. The variability that the microbes experience in dissolved CO and H 2 concentrations, as well as the corresponding biomass-specific uptake rates q CO and q H2 , are displayed in Figure 3.
The mean values for the Eulerian fields follow the same trend as the ideal-mixing model (Equation (8)), indicating that the concentration range is predominantly c L,i < K S,i . For c X below 5 g L −1 , the high potential mass transfer capacity compared to the reaction rate leads to strong CO inhibition, while in the 5-10 g L −1 range, the mass transfer rate is in equilibrium with microbial syngas consumption at lower dissolved gas concentrations. At high c X (25 g L −1 ), gas uptake is fast, leading to low c L,i and thus also to low uptake rates. This decrease in q i is compensated by the greater c X , causing the volumetric reaction rate and gas conversion to remain similar to cases with less biomass and higher biomass-specific uptake rates ( Figure S6). The ideal-mixing model suggests that an optimum q i could be obtained at a certain biomass concentration, but the exact biomass concentration remains difficult to be determined using the CFD models, considering the wide concentration distribution, the non-linear kinetics, and that iteratively running these models is very time-consuming. As there is less inhibition at higher c X , there could be a possibility of increasing gas conversion by supplying more gas, providing that coalescence remains suppressed by the broth components [31]. There is a large volumetric spread in the dissolved gas concentrations obtained by the CFD models. The highest quartile of concentrations is often a factor of 10 higher than the concentrations in the second quartile (e.g., at 5 g L −1 Q2 of c L,CO starts at around 10 −2 mol m −3 , while Q4 starts at 10 −1 mol m −3 ). This would imply that the micro-organism could experience regular concentration fluctuations of around one order of magnitude. However, due to the non-linear nature of the CO and H 2 uptake kinetics, such fluctuations only lead to minor oscillations in biomass-specific uptake rates. Here, the observed concentration gradients are significantly smaller than those in sugar fermentations with similar τ rxn [16], due to the continuous gaseous substrate supply. However, the spread in the concentration fields may cause an overestimation of uptake rates by the ideal mixing models.
Overall, the ideal mixing model was able to describe the concentration range reasonably well, especially in the limitation regime, and could still be used for quick estimations of dissolved gas concentrations at varying conditions (e.g., increased mass transfer, pressure, or with different kinetics). Then, the spatio-temporal variations c L,i , which can only be obtained by CFD modelling, are to be estimated as ± half-an-order of magnitude around the derived concentrations from the ideal mixing model.

Lifeline Analysis
From the lifelines obtained in cases with 5, 10 and 25 g L −1 biomass, it appears that the micro-organisms could experience frequent fluctuations in solute concentrations (in 5 to 30 s), as visible from Video S1. To quantify the microbial experience, the residence times in the peaks and valleys of substrate concentrations were determined, as well as the average dissolved CO and H 2 concentration during a peak or valley. From the resulting probability density functions, we determined the joint probability of observing a specific residence time and concentration in a peak or valley (Figure 4). The dynamic behavior of the EL-GLR causes a large spread in the observed concentrations and residence times. This makes it impossible to standardize a concentration profile of a lifeline. The differences in concentration between the peaks and valleys are around a factor of 5, but within these there are significant deviations (up to 50%) around their specific mean values.
Although the average gas concentrations during the oscillations are very different for the three biomass concentrations, the microbial residence time distributions are quite similar. This is caused by the similar hydrodynamic behavior in the three cases, resulting from similar superficial gas flow velocity and gas conversion rates, while neglecting the influence of biomass concentration on fluid properties (density and viscosity). Interestingly, in the cases with 5 and 10 g L −1 the dips in concentration lasted sometimes for a very long time (>100 s). This could be due to a recirculation pattern in the downcomer. Since, in the 25 g L −1 case, the concentration difference is small between peaks and valleys, and some gas pockets with relative high CO concentration still exist, such moments were not observed.
During the CO peaks at 5 g L −1 , cells can spend short moments (between 5 and 15 s) at inhibitory concentrations (since K I = 0.25 mol 2 m −6 ). In this case c L,CO still remains at around the values needed for an optimum q CO , so that strong inhibition would not be expected, since t peak < t c ≈ 40 s, based upon the used kinetic model. The precise microbial response to such short moments of potential inhibition is unclear.
Similar results were obtained for H 2 ( Figure S9). The residence times in peaks and valleys were in the same ranges as for CO, with high concentration fluctuations of about a factor of 5 noticed. Because H 2 uptake was inhibited at a relatively low c L,CO (K I,CO = 0.025 mol m −3 ), this resulted in high levels of CO inhibition during peaks. When interested in H 2 (and thus CO 2 ) conversion, CO levels should be kept well below the inhibitory values, which could be achieved by adjusting the inlet gas composition (e.g., by green hydrogen supplementation) and/or by increasing the biomass concentration. Due to the strong fluctuations in c L,CO and the inhibiting effect of CO on H 2 uptake, q H2 was significantly influenced, highlighting the need to study the mutual effect of CO and H 2 fluctuations on q H2 .
The ratio between the average dissolved concentrations of CO and H 2 (c L,CO /c L,H2 ) increases with an increased biomass concentration: c L,CO /c L,H2 ≈ 2 at 5 g L −1 , 3 at 10 g L −1 and 6 at 25 g L −1 . This is caused by the faster decrease of q CO compared to q H2 with biomass concentration (cf. Figure 3c,d), due to decreased CO inhibition. To inspect the level of inhibition by CO in the determined ranges for peaks and valleys, the biomass-specific CO and H 2 uptake rates were calculated for each case using its respective c L,CO /c L,H2 ratio ( Figure 5).
From determining the specific gas uptake rates, it became clear that the reactor should be operated in the limitation regime, when increasing c L,i would result in a greater q i (e.g., at 25 g L −1 ), while inhibitory concentrations are avoided. At low biomass concentration (5 g L −1 ), CO inhibition is already problematic, leading to decreased H 2 uptake rates in the peaks. With 10 g L −1 , a significant increase in q CO is observed when transitioning from a valley to a peak (from 0.3 to 0.8 mol mol x −1 h −1 ), but small increases in c L,CO during the peaks could worsen overall performance, since q CO is close to optimum. From the oscillatory dataset in [29], it was derived that fluctuations in q CO and q H2 in the timescale of hours lead to large increases in q EtOH, and thus the ethanol yield [12]. Scale-down experiments with imposed concentration fluctuations could be informative as to whether this observation also holds for the circulation timescale.
Too low dissolved gas concentrations would cause a thermodynamically infeasible catabolism and thus no syngas uptake at all. Such concentrations were estimated to be around 4 × 10 −4 and 3 × 10 −3 mol m −3 for CO and H 2 , respectively, assuming independent consumption of CO and H 2 /CO 2 for solventogenesis [47]. Since such low c L,CO was not obtained in our analysis, we do not expect such problems for CO consumption. However, for H 2 , values below the thermodynamic limit were attained in the valleys for 10 and 25 g biomass L −1 , so that a coupling with CO consumption is potentially required to supply enough electrons for H 2 uptake. As this may come at the expense of the product yield, further scale-down studies are required to determine how C. autoethanogenum may react to such short-term fluctuations in H 2 concentration. Figure 5. Overall biomass-specific (a) CO and (b) H 2 uptake rates computed for biomass concentrations of 5 (blue), 10 (red), and 25 g L −1 (green). c L,CO in the H 2 uptake kinetics was calculated by using a c L,CO /c L,H2 of 2, 3 and 6 for each respective case. The full arrows indicate the concentration ranges of the peaks, while the dashed arrows indicate the ranges for the valleys. CO-uptake is independent of the biomass concentration, hence the single line.
It has been estimated that a starvation regime could occur when c L,CO < 3 × 10 −3 mol m −3 [17]. Since, in this regime, a major portion of the energy might be spend on maintenance catabolism, lower growth rates can be expected, leading to higher product yield. In the configuration they studied, such CO shortages were highly likely to occur, causing a probable shift towards a starvation regime. In our simulations, this situation may only occur in the valleys when operating with high biomass concentrations. Due to higher c L,CO in our other cases, reaching the maintenance catabolism was very unlikely to occur, since k L,CO a of the EL-GLR was a factor 5 higher than in the BCR studied in [17] (with k L a~0.01 s −1 and c X = 10 g L −1 ).
Our results suggest that even higher biomass concentrations may be advantageous, considering the current operation in the limitation regime and that high mass transfer could be obtained due to bubble coalescence suppression in the fermentation broths. Operation at very low c L,CO would enable operational flexibility and a high product yield, without sacrificing gas conversion. Caution is needed to prevent the dissolved gas concentrations becoming so low that the reaction becomes thermodynamically unfeasible, or that the high biomass concentrations hampers mass transfer and mixing by increasing broth viscosity [50].
In the LanzaTech process, k L,CO a could well be around 3-4 times higher [4] than the final one obtained in our model (650 h −1 vs. 180 h −1 ). This could be due to mass transfer intensification (e.g., by introduction of perforated plates [51]) or by achieving smaller bubbles (~1 mm). Although bubbles would become more rigid in such a case, mass transfer might still be enhanced by cell monolayers around the bubbles, especially in case of operating at high biomass concentrations [45]. In such high mass transfer cases, substantial CO inhibition might be expected, stressing the need to operate at high (>25 g L −1 ) biomass concentrations.

Development of Scale-Down Simulator
Based upon the analysis of the CFD data, we used numerical simulations to propose a conceptual design of a scale-down simulator, in order to experimentally replicate the dissolved gas concentrations which were estimated to be experienced by micro-organisms in the industrial-scale syngas fermentation. A single-vessel system, with a 2 L working volume CSTR, was chosen as a basis for the scale-down simulator. Since rapid and irregular fluctuations in peaks and valleys should be obtained, we did not consider multi-vessel systems with forced circulation, because rapid consumption of dissolved gas in the tubes connecting the vessels would be detrimental to performance. In addition to the advantage of no substrate depletion in the tubes, clogging or high shear stress by pump action could also be avoided. In a well-mixed stirred-tank, there are little spatial variations in dissolved gas concentration (unlike plug-flow systems) so that the dissolved gas concentrations can be better controlled. Another potential scale-down system would be a two-stage STR, with a perforated plate separating two well-mixed zones [52]. In this way, the dynamic interchange between two concentrations could be reproduced, although the step transitions might be unrealistic for the observed large-scale behavior. A significant disadvantage of using a single-vessel is the lack of population heterogeneity in cellular experience, which would definitely be present in multi-vessel systems, as well as the poor incorporation of dead (or low concentration) zones, such as the downcomer [53,54]. Furthermore, slow mixing at low power input might possibly lead to local concentration gradients [55]. Operation at smaller scales with better mixing (e.g., 200 mL) might, however, lead to practical problems regarding sampling.
To mimic the large-scale successfully, we should make sure that the microbes experience similar peak/valley durations and the same concentration differences as in the large-scale bioreactor. Although we could argue from Figure 4 that most peaks last between 5 and 15 s, and valleys 10 to 30 s, this argument does not account for frequently occurring irregularities. There are peaks and valleys that largely exceed these times, e.g., peaks of 25 s and valleys of 60 s are not exceptional, and the scale-down simulator should replicate such outliers in terms of time and concentration. With the probability distributions of the residence times derived from the CFD lifelines, variations in the stirrer speed were imposed to obtain corresponding peaks and valleys in the scale-down simulator. It was determined that around 2000 oscillations, lasting~15 h in total, should be applied in the scale-down simulator to make sure that enough variation in peak or valley residence time is imposed ( Figure S10). To account for the varying biomass concentrations from the three CFD cases (5, 10 and 25 g L −1 ), the biomass recycling rate R rec was altered between the different cases to adjust the reaction rate.
With this computational set-up and the iteratively derived operational conditions (Table 3), lifelines were simulated in the conceptual bench-scale reactor for the three different biomass concentrations that roughly correspond to the large-scale lifelines ( Figure 6). The pulses in stirring speed are well captured and provide the same peak-valley frequencies as is expected at the large-scale. The concentrations that the microbes would experience are similar as in the large-scale bioreactor within the peaks and valleys. For example, for CO and H 2 in the 5 g L −1 cases, the upper concentrations are always in the same order of magnitude and the experienced valleys are very similar to those in the large-scale reactor (Figure 6a,b). For the cases with high c X , it is more challenging to represent the deep concentration valleys well (c L,i << 10 −3 ), since the increased c X generally requires more mass transfer in the valleys. Since the impact of such concentrations on the biomass-specific uptake rates is small (Figure 5), negligible influence was expected. Table 3. Operational conditions of the scale-down simulator to obtain an acceptable fit of the lifelines obtained by the scale-down simulator with CFD-derived lifelines.

Peak
Valley Recycle  The rate of increase in dissolved gas concentration during the transition from a valley to a peak in the scale-down simulator is very similar to that in the large-scale bioreactor: instantaneously, the microbes experience concentration increases of up to around 1-2 orders of magnitude in a matter of seconds (1-5 s). The decrease in the slope at the beginning of the peaks in dissolved CO concentrations can equally be identified in some of the peaks of the CFD lifelines. In the large-scale, this rapid increase is due to the micro-organism travelling instantly into a zone with high mass transfer and thus dissolved gas concentrations, while, in the scale-down simulator, the mass transfer increase responds to the step increase in stirring speed.
The transition from the peak to the valley was found to be more problematic to reproduce in the bench-scale reactor for cases with 10 and 25 g L −1 . Although the dissolved gas concentration decay is more plug flow-like at the large-scale, immediate decreases back to a representative "valley-baseline" were observed in the scale-down simulator. Simulating a ramped decrease in stirring speed could be helpful in obtaining a more realistic decay in concentrations.
A major factor varying between the two scales is the frequency and magnitude of cellular exposure to shear forces. This was quantified with the energy dissipation/circulation function (EDCF = P total V eff 1 t c ) [56,57] while approximating t c as 1 4 t m [55]. In the scale-down simulator, there is highly varying exposure to shear between peaks and valleys when the cells are close to the impellers (EDCF varies from 50 kW m −3 s −1 and 1 × 10 −4 kW m −3 s −1 , respectively). In the EL-GLR. there is only a high shear region around the gas plume (EDCF 0.06 kW m −3 s −1 ), without considering bubble burst. The significantly varying EDCF between two scales could be impactful when C. autoethanogenum is shear-sensitive, but this was not expected, due to its small size (around 3 µm) compared to the Kolmogorov scale (~10 µm) [58,59].
The correspondence of the results of the scale-down simulator with the large-scale reactor was determined by performing a lifeline analysis. In this way, the probability distributions for the residence times and the concentrations in the peaks and valleys could be compared quantitatively (Figures 7 and S11-S15). Generally, a very good correspondence of the residence time distributions was obtained. To some extent this is logical, as CFD results of these are the inputs of the scale-down simulator, although the limitations of a bench-scale reactor do not guarantee sufficiently good correspondence to be feasible; this indicates the feasibility of imposing rapid stirring speed fluctuations in a well-mixed bench-scale system.
Corresponding concentration distributions were more challenging to obtain, since the ranges of the large-scale CO peak concentrations are very large (0.1-0.25 mol m −3 , Figure 7a). In the scale-down simulator, the CO peak concentration could not become that high (maximum 0.2 mol m −3 ) so that a narrower range was obtained, which was more skewed towards the lower concentrations. The assumptions and parameters used in the mass transfer and kinetic models make it challenging, however, to rely purely on quantitative results for the concentration fluctuations in the scale-down setup. Using more accurate mass transfer and kinetic models would increase the reliability of our quantitative predictions and thus our conceptual scale-down simulator.
Despite all these limitations, we showed that. with a conceptually relatively simple scale-down simulator, the large-scale dissolved gas concentration gradients for a wide range of biomass concentrations could be reproduced at lab-scale. Model-based tuning of the operational conditions (e.g., stirrer speed, gas flow rate, gas composition) of the scale-down set-up on the probability distribution functions of the large-scale reactor, is a possible strategy to maximize correspondence between the two scales and thereby provide a fruitful basis for representative scale-down of syngas fermentation.

Outlook
Further improvement of the scale-down simulator could lead to even better representations of the large-scale behavior. An optimization routine could help in obtaining the best-fit parameters with the CFD-derived data. The ideal scale-down simulator has as few as possible variable parameters and represents the large-scale behavior for a wide range of conditions. The effect of parameters (e.g., stirrer speed during the start-up phase) could be derived using tools such as principal component analysis during the optimization and in helping to decide whether or not to use the parameter in further analyses. The quantitative results of the CFD study and the proposed scale-down study are strongly dependent on the process conditions (e.g., headspace pressure, gas fraction), as well as on the kinetic model for CO and H 2 uptake. For example, the current model for COuptake significantly underestimates the maximum specific growth rate (µ max model = 0.03 h −1 ) compared to experimental values (µ max exp = 0.12 h −1 ) [60]. Since the used kinetic models for gas uptake are parametrized using insufficient data, the accuracy of our simulations is decreased and is therefore a major drawback of the study. Development of accurate kinetic models is crucial for reliably modelling bioreactors, and we hope that our work motivates further research in this area. The MATLAB scripts describing the conceptual scale-down simulator are openly available and can be used for further development with updated sub-models.se.
Despite the limitations, the proposed set-up and method are still applicable to a wide range of conditions. Even without using CFD, but with an ideal mixing model (Equation (8)), one can estimate the effect of process variables (e.g., increased mass transfer rates) on the average concentration in the large-scale reactor. In case the reactor is operated in the gas limitation regime (c L,i < K S,i ), spatial concentration differences of a factor of 5 around the mean could be expected based upon our CFD simulations. Since, in the scenarios with varying biomass concentrations the residence time distributions in the peaks and valleys do not greatly differ, neither are such differences expected in other operational cases in the limitation regime. A scale-down simulator could then be tuned based upon the estimated concentration differences and the residence time distributions.
Although the developed scale-down simulator is conceptually easy to understand, practical installation and operation might be challenging. The repeated variations in stirrer speed at high frequencies should be controlled, along with possible ramp phases when increasing or decreasing the stirring speed. Ideally, rapid-sampling and/or online measurement for CO and H 2 [61,62] should be applied to make sure that peaks and valleys are obtained in the intended manner. As the peaks and valleys are applied in the secondtimescale, probe lag should be taken into account when analyzing the experimental data. Furthermore, the influence of broth components on k L a in real fermentations should be considered for better predictions [31]. If the in-situ concentrations cannot be measured, they could be predicted using a precisely determined k L a [63].
Our conceptual scale-down simulator makes it possible to simulate a statistically representative lifeline of the EL-GLR within a fraction of the time that it would take for the CFD model to be run and analyzed. Such a lifeline could then be used to study interactions between the extra-and intracellular environment by coupling with a metabolic kinetic model [12]. The obtained response should then be similar to the large-scale response, due to the correspondence of both lifelines. By making variations in such lifelines, the peak and valley residence time and concentration distributions can be obtained, which could lead to a desired large-scale response (e.g., high ethanol specificity).
Ramp and feast-famine studies in the scale-down set-up could be used to parametrize kinetic models that describe the short-term response of C. autoethanogenum [12,28,64], by rapid sampling of metabolite and enzyme concentrations. Ramp studies would be helpful in determining whether the instantaneous electron supply in the peaks would indeed lead to increased ethanol production, as was expected from long-term oscillations [29]. If the gas uptake rates are product-independent, then such scale-down simulators could be used for engineered strains to produce higher-value products [3], proteins [65], or for coupled reactions with other micro-organisms, such as chain-elongators [66] or PHA production [67]. With the scale-down simulator, the microbe could be adapted to large-scale conditions, so that fewer scale-up problems might be expected [14].
The analyses of the industrial syngas fermentation process in this and our previous study [4] are all model-based and only slightly tuned, based upon the scarcely available literature data of the full-scale LanzaTech operation [1,2,68]. To advance the syngas fermentation process, for model validation and the execution of highly representative scale-down simulators, the publication of real industrial data would be required, such as large-scale circulation times, operational k L a values, a range of dissolved gas concentrations and their gradients. All of this could, for example, enable the utilization of a broader range of gas compositions, the development of processes towards higher-value products and intensified fermentation equipment.
In our analyses we showed that high biomass concentrations (e.g., c X > 10 g L −1 ) might be advantageous for both product yield and gas conversion. Since the highest reported biomass concentration in syngas fermentation reactors is around 9 g L −1 [48], experiments should target the influence of increased biomass concentrations and its viability on gas uptake, broth viscosity and mass transfer. The precise operating interval in terms of dissolved gas (CO and H 2 ) concentrations should be retrieved experimentally, so that the thermodynamically infeasible range can be avoided while operating in the maintenance-dominated regime, which would enable high product-to-substrate yields. Our results show that, with the currently used syngas composition (20% H 2 ) and high biomass concentration (25 g L −1 ), H 2 catabolism may be thermodynamically infeasible, although co-consumption with CO might occur [9]. Thus, in a future with intermittent green hydrogen supply from renewable resources [69], supplementation of hydrogen might be a good option to valorize excess electricity and increase the CO-to-product yield.

Conclusions
The effect of biomass concentration and dissolved gas concentration fluctuations in large-scale syngas fermentation was studied with Euler-Lagrangian simulations. Based upon these numerical simulations, we recommend industrial operation at relatively high biomass concentrations, as this would reduce the effects of CO inhibition, could increase the product yield and would provide high operational flexibility. Simulations indicate that, in large-scale syngas fermentation, C. autoethanogenum will experience frequent oscillations (peaks and valleys) in a dissolved gas (CO, H 2 ) concentration of about one order of magnitude, in a timescale of seconds (5 to 30 s). Such concentration fluctuations may occur irrespective of the biomass concentration and were hypothesized to favor the ethanol yield.
The large-scale concentration fluctuations should be simulated during small-scale experiments to study how C. autoethanogenum adapts to industrial-scale conditions. We proposed a single-vessel scale-down simulator that theoretically replicates the fluctuations in dissolved gas concentrations by varying the stirrer speed based on the large-scale oscillations. Numerical analysis shows that the durations of the oscillations could be well replicated, but the settings might be adjusted to achieve higher similarities for the variations in concentration. The obtained lifelines in the proposed scale-down simulator well represent the large-scale reactor for a wide range of biomass concentrations and operational conditions. Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/bioengineering10050518/s1, Figure S1: Mixing study in the EL-GLR; Table S1: Characteristics of the scale-down simulator; Figure S2: Kullback-Leibler divergence of the residence time distribution in a CO peak; Figure S3: Kullback-Leibler divergence of the residence time distribution in a CO valley, H 2 peak, and valley; Figure S4: Surface plots of the time-averaged gas hold-up in the EL-GLR; Figure S5: Variation of dissolved CO and dissolved H 2 concentrations within the EL-GLR reactor volume; Figure S6: CO and H 2 conversion and consumption rates for varying biomass concentrations; Figure S7: Surface plots of the dissolved CO concentration in the EL-GLR; Figure S8: Surface plots of the dissolved H 2 concentration in the EL-GLR; Figure S9: The probability of a microbe to experience a specific H 2 concentration peak or valley. Figure S9: Kullback-Leibler divergence for the number of peaks used as input for the scale-down simulator, and the operational time of the scale-down simulator; Figure S10: Probability density functions comparison for H 2 with 5 g L −1 biomass; Figure S11: Probability density functions comparison for CO with 10 g L −1 biomass; Figure S12: Probability density functions comparison for H 2 with 10 g L −1 biomass; Figure S13: Probability density functions comparison for CO with 25 g L −1 biomass; Figure S14: Probability density functions comparison for H 2 with 25 g L −1 biomass. Video S1: Video of Euler-Lagrangian simulation. (References [70][71][72][73] are cited in Supplementary Materials).  Acknowledgments: The authors would like to thank colleagues at the TU Delft Bioprocess Engineering group for their critical input in discussions, with special mention to Rob van der Lans.

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