Design of a Bubble Reactor for Altitude Simulators Used to Humidify a Combustion Air Stream by Means of CFD Multi-Phase Models

: In this paper, a procedure for the design of a bubble reactor which allows the control of the humidity of a gas stream used as combustion air is presented. This reactor is designed to be used as a component of an altitude simulator test facility for the optimization, homologation and calibration of new hybrid engines. The design has been carried out by means of Computational Fluid Dynamics (CFD) multi-phase models and validated against the experimental data obtained from the developed prototype. A discussion about the adequate mesh topology and cell size is presented, as well as a comparison between the two available models for the air–water interphase. Lastly, a validation of the CFD results using experimental data shows that the model that should be used is the multi-regime interaction model, from which the ﬁnal design for the bubble reactor was obtained.


Introduction
The difficulty and cost of the calibration process of automotive engines has increased since the inclusion of extended testing conditions during RDE (real driving emissions) tests in the homologation process of the new Euro6 [1]: altitude up to 1300 m over sea level and −7 • C. In addition, the Worldwide Harmonized Light-Duty Test Cycle (WLTC) has been included in the homologation procedure. This new driving cycle is more demanding and closer to real driving conditions than the previous New European Driving Cycle (NEDC).
This situation become even more important when considering new control strategies of future automotive powertrains that combine different energy sources. The hybrid powertrain (electric motor + internal combustion engine (ICE)) seems to be the technology selected for the short to medium term [2]. The goal is to optimize the hybrid vehicle pollutant emissions along real driving cycles that are unknown in advance [3,4].
Atmosphere simulators test facilities are important in order to develop these methodologies to find the optimal control strategies of future hybrid powertrains in a real environment. Therefore, the development of test facilities that could be able to reproduce real driving conditions efficiently and accurately is a current need for the development of new hybrid engines.
These installations are very complex [5,6] and demand detailed studies on the control strategies in different ambient conditions. The effect of the ambient conditions on the hybrid engine is a topic that has not been extensively studied, and both the electric motor and the ICE of the hybrid power train are affected. On the one hand, the dielectric response of the electric motor insulation varies greatly with the level of humidity at which it is exposed [7]. Furthermore, the ambient humidity has an effect over maintenance procedures of the stator of electric motors during partial discharge procedures [8], whose repairing effect decreases when increasing the humidity.
On the other hand, the decrease of air density caused by the altitude conditions [9] affects greatly the performance and emissions of the engine [10], due to the change in the combustion chamber and boundary conditions of the engine, which leads to changes in the behavior of some processes. One of them is the increase in backpressure generated by the variable geometry turbine, due to a higher need in compression rate to achieve the same boost pressure than at sea level. The increase in back pressure can lead to backflows in the exhaust valves, increasing the internal EGR (exhaust gases recirculation) rate as well as the pumping losses of the engine, which reduces the performance [11].
Moreover, the lower density generates smaller fuel drop size inside the cylinder, which reduces the coalescence [12,13] thus decreasing the jet-to-wall impingement [14] as the fuel mixes better with the air and evaporates faster inside the combustion chamber. However, the penetration of the fuel spray increases, due to the decrease of the spray angle caused by the low air density [14].
In order to study the effect of the altitude on the engine performance, the use of a hypobaric chamber is the solution commonly used during the early stages of the calibration process instead of real driving tests, more used at the final stages of the engine calibration. The use of hypobaric chambers has become a standard due to the high costs and low repeatability of real driving tests. Nevertheless, the hypobaric chambers have some important drawbacks. This kind of installations are expensive and voluminous and requires extensive planning to build one. This has meant that few hypobaric chambers are available at research and development (R&D) centers, which increases greatly the development time. Moreover, the human body cannot withstand continuous pressure changes [15] so operating inside hypobaric chambers entails a serious health risk what limits the number of testing hours in chassis dyno test cells.
Therefore, to avoid the drawbacks of hypobaric chambers, devices that generate altitude only inside the intake and exhaust of the engine have been developed. These devices, called altitude simulators (AS), can generate the altitude for the engine while keeping the rest of the test cell at room pressure.
Moreover, some AS have additional conditioning modules that allow for controlling control temperature and humidity for a complete simulation of ambient conditions. One example of this family of equipment is the Multifunctional Efficient and Dynamic Altitude Simulator (MEDAS), MEDAS Temperature Module (MTM) and MEDAS Humidity Module (MHM) product from HORIBA [16][17][18][19], which allows for simulating the atmospheric conditions for almost every location on the planet in which a vehicle can be driven.
These two features of the ambient conditions are also important with respect to the effect on the performance of an ICE. The ambient temperature has a great effect on the performance of the engine, since the lower the temperature of the air the higher its density. This leads to an improvement of the volumetric efficiency of the engine, which is a measure of the effectiveness of the air intake system [20]. In addition, this temperature is also important during the combustion process. If the ambient temperature lowers, the maximum temperature during the combustion decreases, which reduces the NOx formation [21]. However, this also increases the emission of hydrocarbons, due to lower combustion efficiency and a higher quantity of fuel attached to the cylinder walls [22].
In the case of the humidity, its effect on the combustion process has not been studied in much detail, but some studies on gas turbines [23] show the negative effects of the increase of water content in the efficiency of the combustion, and thus on the efficiency and performance of the engine.
With respect to controlling the humidity for the combustion air, two technologies are mostly used: atomizing water or an air current or passing air through a bubbling reactor. The atomizing allows the amount of water injected in the air to be controlled accurately, but it needs piping length to let enough time for the water droplets to break up and dissolve into the air [24]. Contrarily, the reactor size is limited by the amount of air that has to pass through it, to control the maximum air velocity [25] and avoid an excessive increase in the water level that could lead to overflowing. The reactor solution is quite similar to fluidized bed technologies, widely used in the chemical industry. However, the procedure of design of this kind of reactor is focused on the interaction between the particles and the fluid [26]. Instead, the increase in height of the emulsion, which is the main factor that leads to the reactor overflowing is not studied in detail [27].
Moreover, these existing technologies are focused on large-scale installations, as industrial chemistry plants, with a high working mass flow rate over 1 kg/s [28], or small installations such as dedicated laboratories with mass flow rates below 0.01 kg/s [29]. Hence, the utility of the development of the installation is described in this paper. It covers the medium installation range, in which almost all the propulsion systems for land transport are located. Either fuel cells, ICE with liquid or gas fuels, or other kind of propulsion systems that require an air stream provided to the engine in order for it to work.
Therefore, a detailed analysis of a bubble reactor has been carried out. The objective of the reactor is to humidify an air stream up to a water content set point, which afterwards is going to be expanded and cooled to reach temperature and pressure set points. In Section 2 the design requirements for the installations will be explained. In Section 3 the designed experimental facility and the testing methodology will be described. Afterwards, the set-up parameters for the CFD code used will be shown in Section 4. Also, the study of the effect of the mesh size and air-water interaction models is explained, as well as the results of the parametric studies carried out to obtain the final design of the bubble reactor. Lastly, the conclusions obtained through this work will be presented in Section 5.

Equipment Requirements
Up to now there have not been many studies about the effect of the humidity in the propulsion systems of land transport, either due to lax normative in this regards or lack of proper testing facilities to carry them out. The device described in this paper has been designed to cover this gap in the literature; by developing a device that can generates and control dynamically high humidity levels even at altitude conditions with high temperature. When talking about water content in the air, we must highlight the great impact that psychometric variables pressure and temperature have over the water content in the air. The water content or specific humidity (g water/kg dry air) increases exponentially when the temperature goes up or when the pressure decreases while keeping constant the relative humidity value [30].
In Figure 1, the greenish blue area represents the specifications for the water content of the equipment. The specification has been designed given every combination of pressure, temperature, and relative humidity inside the operational range of the installation: 0.55 bar to 1.3 bar, 0 • C to 40 • C and 10% to 100%, respectively. Furthermore, the line of maximum specific humidity can be seen depicted in red, which represent 100% relative humidity at a pressure of 0.55 bar. Finally, the maximum humidity at sea level is represented in this figure by an orange line, which clearly shows how it is impossible to achieve the required humidity, in order to emulate high altitudes, just by humidifying at sea-level conditions.
The most unfavorable specifications for the humidifying process, specifically 2.5 bar and 100% RH is represented by a black line in Figure 1. These operation conditions justify the selection of a bubble reactor as a humidifier instead of atomizers, which have a maximum working pressure too restrictive for this application as well as need an extensive length of piping to fully evaporate the water injected into the air stream. In addition, due to the humidity level of the most demanding point of the installation, depicted by a black cross in Figure 1 (0.55 bar, 40 • C and 80% RH), at the most unfavorable pressure conditions (2.5 bar) it is necessary to moisten this at a temperature of around 70 • C, as shown at Figure 1. The reason why 2.5 bar could be needed is to boost combustion air and to overcome pressure losses along the rest of the elements of the system in some operative conditions [16]. Figure 2 shows the location of the bubble reactor in the installation scheme (13) of Figure 2, which has been highlighted for easier identification, with the rest of the key elements responsible for carrying out the conditioning of the air. The most unfavorable specifications for the humidifying process, specifically 2.5 bar and 100% RH is represented by a black line in Figure 1. These operation conditions justify the selection of a bubble reactor as a humidifier instead of atomizers, which have a maximum working pressure too restrictive for this application as well as need an extensive length of piping to fully evaporate the water injected into the air stream. In addition, due to the humidity level of the most demanding point of the installation, depicted by a black cross in Figure 1 (0.55 bar, 40 C and 80% RH), at the most unfavorable pressure conditions (2.5 bar) it is necessary to moisten this at a temperature of around 70 C, as shown at Figure  1. The reason why 2.5 bar could be needed is to boost combustion air and to overcome pressure losses along the rest of the elements of the system in some operative conditions [16]. Figure 2 shows the location of the bubble reactor in the installation scheme (13) of Figure 2, which has been highlighted for easier identification, with the rest of the key elements responsible for carrying out the conditioning of the air.  The airflow depicted in Figure 2 starts with a mechanical compressor (2) that is electrically driven, which sucks air form the ambient and increases the pressure of the installation. This mechanical compressor is protected from surge by means of a control valve (3). Then, a control valve (4) controls the temperature at point 'g' mixing cold air, cooled on the water cooler (5) and chiller (6), and hot air from the compressor outlet. Downstream, the valve (9) controls the humidity at point 'l', forcing part of the airflow to pass through the bubble reactor (13). Then, the flow is expanded up to the desired pressure on control valve (15) and the variable geometry turbine (17). Additionally, control valve (16) regulates the air that passes through (17) to control the temperature at point 'l', where the engine (21) sucks the amount of air needed and the rest flows through the dilution pipe. Lastly, the dilution flow mixes with the exhaust gases of the engine (22), gets cooled down on the exhaust cooler (23) and is sucked out of the installation with a mechanical compressor (24). The airflow depicted in Figure 2 starts with a mechanical compressor (2) that is electrically driven, which sucks air form the ambient and increases the pressure of the installation. This mechanical compressor is protected from surge by means of a control valve (3). Then, a control valve (4) controls the temperature at point 'g' mixing cold air, cooled on the water cooler (5) and chiller (6), and hot air from the compressor outlet. Downstream, the valve (9) controls the humidity at point 'l', forcing part of the airflow to pass through the bubble reactor (13). Then, the flow is expanded up to the desired pressure on control valve (15) and the variable geometry turbine (17). Additionally, control valve (16) regulates the air that passes through (17) to control the temperature at point 'l', where the engine (21) sucks the amount of air needed and the rest flows through the dilution pipe. Lastly, the dilution flow mixes with the exhaust gases of the engine (22), gets cooled down on the exhaust cooler (23) and is sucked out of the installation with a mechanical compressor (24).
Due to this positioning of the bubble reactor, the humidifying process is carried out before the key expansion and final tune of the temperature. This means that the air will be humidified under pressure conditions over sea level, since the air is compressed in the mechanical compressor (2), but this will allow the temperature control range of (17) to be increased in order to compensate. Since following psychometric theory, the higher the pressure the lower the amount of water that the air can withstand; whilst, with respects to the temperature, the higher the temperature the higher the possible water content. Therefore, to achieve the relative humidity levels while reproducing altitude conditions it is necessary to moisten the air either after generating the desired vacuum or before generating it but at a much higher temperature in order to compensate the increase of water capacity that the air has when decreasing the pressure.
An example of this phenomenon has been represented in the psychometric diagram of Figure 3, in which can be seen the effect of the pressure when trying to achieve the black dot conditions, 0.54 bar, 30 • C, and 70% RH (relative humidity). If the humecting process is carried out at 1.013 bar (red square of Figure 3), is not possible to achieve the water content set point while keeping the temperature at 30 • C. The set point has a higher water content than the maximum water that the air can contain at that pressure (dotted dark blue line). To achieve this objective, it would be necessary to increase the temperature up to 41.3 • C while keeping 1.013 bar (green triangle), and then expand and cool down the air. The other option is to humidify the air after lowering the pressure to 0.54 bar, which would allow the temperature to be kept constant as the maximum water at this pressure is higher than the required (dark blue line).
Appl. Sci. 2021, 11, x FOR PEER REVIEW 6 of 16 content set point while keeping the temperature at 30 °C. The set point has a higher water content than the maximum water that the air can contain at that pressure (dotted dark blue line). To achieve this objective, it would be necessary to increase the temperature up to 41.3 °C while keeping 1.013 bar (green triangle), and then expand and cool down the air. The other option is to humidify the air after lowering the pressure to 0.54 bar, which would allow the temperature to be kept constant as the maximum water at this pressure is higher than the required (dark blue line). Lastly, there are some additional problems involved in the design of the bubble reactor and the optimal water level at rest. On the one hand, if the water column height at rest is too high the risk of overflowing the reactor increases since the bubbling of the air stream creates foams in the water, but without a benefit for the humidifying process [22]. On the other hand, if the water column height at rest is too low the shorter residence time of the air inside the water also reduces the water absorption capacity of the air.
Therefore, the accurate prediction of the behavior of the air inside the water in the reactor as well as the foam level, as a function of the water column height at rest and the air stream velocity is indispensable to maximize the humidifying capacity of the installation while keeping controlled the stabilized height of the water column. Lastly, there are some additional problems involved in the design of the bubble reactor and the optimal water level at rest. On the one hand, if the water column height at rest is too high the risk of overflowing the reactor increases since the bubbling of the air stream creates foams in the water, but without a benefit for the humidifying process [22]. On the other hand, if the water column height at rest is too low the shorter residence time of the air inside the water also reduces the water absorption capacity of the air.

Development of Experimental Facility
Therefore, the accurate prediction of the behavior of the air inside the water in the reactor as well as the foam level, as a function of the water column height at rest and the air stream velocity is indispensable to maximize the humidifying capacity of the installation while keeping controlled the stabilized height of the water column.

Development of Experimental Facility
By means of some preliminary calculations run with CFD code, which will be detailed in following sections, a prototype of the bubble reactor has been designed and used to carry out the experimental validation of the CFD model used.
The tests were carried out controlling the mass flow of air that passed through the water as well as the initial height of the water volume. The air mass flow tested was between 0.16 and 0.25 kg/s, and the water level was kept at the 0 cm mark in the measurement rule when at rest. Furthermore, the density of the air at the reactor outlet was kept at a similar value for every test, with an air pressure of 2.4 bar and a temperature of 50 • C. It should be noted that the range of mass flow selected, as well as the pressure and temperature conditions, were the available controlling range of the equipment (AS HORIBA MEDAS) in which the prototype designed was installed.
The prototype built for the testing and validation is the one that appears in Figure 4, in which the external view of the bubbling reactor is depicted (Figure 4a). In this figure, the air inlet to the water can be seen (number 1 of Figure 4), as well as the air outlet (2). Furthermore, for the validation tests, the air circuit through the heat exchangers, present inside the reactor, has been closed, as shown in numbers 3 and 4 (inlet and outlet respectively). inside the reactor, has been closed, as shown in numbers 3 and 4 (inlet and outlet respectively).
(a) (b) (c) The reactor has been designed with a complex wall shape to help withstand the pressure difference between the high inside pressure and the ambient pressure, since it reduces the size of the flat areas. The bigger area inside the reactor has a footprint of 0.81 m 2 , while the smaller one has an area of 0.64 m 2 .
In Figure 4b, a cut section of the reactor is shown, and the different internal components of the design are depicted. Firstly, a holed duct in the bottom (5), which is submerged in water and works as an air distributor. The holes have a diameter of 18 mm and follow a grid distribution in the first half (point 11 of Figure 4c) with a distance between holes of 25 mm in the axial direction and 34 mm in the tangential direction. In the second half of the duct (point 12 of Figure 4c) the distribution is scattered, which allows a higher hole density, as the distance in the tangential direction decreases to 17 mm. Next, a first heat exchanger (6) which also works as a flow homogenizer to avoid all the air following a single narrow path. Then, there is a pair of reinforcing crosses (7) to help withstand the pressure difference in the walls. Four visual access points are inside the reactor (8), and a measurement rule (9), with which the water level measurements during testing were car- The reactor has been designed with a complex wall shape to help withstand the pressure difference between the high inside pressure and the ambient pressure, since it reduces the size of the flat areas. The bigger area inside the reactor has a footprint of 0.81 m 2 , while the smaller one has an area of 0.64 m 2 .
In Figure 4b, a cut section of the reactor is shown, and the different internal components of the design are depicted. Firstly, a holed duct in the bottom (5), which is submerged in water and works as an air distributor. The holes have a diameter of 18 mm and follow a grid distribution in the first half (point 11 of Figure 4c) with a distance between holes of 25 mm in the axial direction and 34 mm in the tangential direction. In the second half of the duct (point 12 of Figure 4c) the distribution is scattered, which allows a higher hole density, as the distance in the tangential direction decreases to 17 mm. Next, a first heat exchanger (6) which also works as a flow homogenizer to avoid all the air following a single narrow path. Then, there is a pair of reinforcing crosses (7) to help withstand the pressure difference in the walls. Four visual access points are inside the reactor (8), and a measurement rule (9), with which the water level measurements during testing were carried out. A second heat exchanger (10) also serves as a droplet separator to avoid liquid water to flow out of the reactor. The heat exchangers are used to heat up the water using by-passed dry air at higher temperature than the water (see Figure 2).
Finally, the testing procedure for the developed installation was as follows: • In the first place, before passing air through the water, the water level was adjusted to the 0 cm mark of the measurement rule • Then, the installation was started and the mass flow, pressure, and the temperature of the air at the outlet of the bubble reactor was controlled • Lastly, once the three controlled variables are stabilized, the water column height was measured through the sight glasses with the help of the measurement rule. In following sections, details about the CFD code used for the final design of the bubble reactor will be explained, considering the type of sub models and solver used, as well as the validation of the choices made.

Solver and Models
The segregated solver of a finite volume commercial flow program was used to solve the discretized Navier-Stokes equations (StarCCM+ ® v.15.02.007) [31]. The coupling between the momentum and continuity equations was achieved with a predictor-corrector approach. The flow was unsteady, incompressible, and turbulent. Unsteady flow calculations (URANS) were performed by means of first order implicit formulation with an adequate time-step to keep Courant-Friedrichs-Levy around 1.0 to ensure converge and stability (time step = 0.001 s). However, explicit formulations are not available for the Eulerian multiphase model used in this work. Therefore, the Shear-Stress Transport (SST) k−ω two-equations turbulence model has been used for turbulence. This approach is the only one suitable for the two-fluid modelling developed in this work. Table 1 summarizes the main aspects of the parameters describing the CFD model.  In following sections, details about the CFD code used for the final design of the bubble reactor will be explained, considering the type of sub models and solver used, as well as the validation of the choices made.

Solver and Models
The segregated solver of a finite volume commercial flow program was used to solve the discretized Navier-Stokes equations (StarCCM+ ® v.15.02.007) [31]. The coupling between the momentum and continuity equations was achieved with a predictor-corrector approach. The flow was unsteady, incompressible, and turbulent. Unsteady flow calculations (URANS) were performed by means of first order implicit formulation with an adequate time-step to keep Courant-Friedrichs-Levy around 1.0 to ensure converge and stability (time step = 0.001 s). However, explicit formulations are not available for the Eulerian multiphase model used in this work. Therefore, the Shear-Stress Transport (SST) k−ω two-equations turbulence model has been used for turbulence. This approach is the only one suitable for the two-fluid modelling developed in this work. Table 1 summarizes the main aspects of the parameters describing the CFD model.
The two-fluid model [33,34] implemented in StarCCM+ ® provides a flexible framework where dispersed, stratified, or mixed two-phase flow (air + water) can be modeled by using different closure relations. For a given pair of phases, the multiphase segregated flow model allows different types of phase modelling topology to be considered. The continuous-dispersed phase interaction [35] topology provides the means to model the interaction of a phase which is dispersed within another phase throughout the whole domain. Furthermore, the interface length scale is taken to be an effective mean diameter of the dispersed phase particles. Moreover, where particles are not spherical, this is absorbed as a correction factor into interphase transfer models such as for bubble drag and lift force.
In addition, the multiple flow regime [32] topology is an extended definition where all types of multiphase flow, that is, dispersed (bubble, droplet or particle) and stratified (freesurface) flows, can be modelled. The solver employs a volume-fraction based criterion to classify the flow regime into dispersed and separated two-phase flow. In this scenario, the multiple flow regime model treats bubbly flow in the first dispersed regime, droplet flow in the second dispersed regime, and separated water-air flow in the large-scale interface regime.

Computational Domain
The computational domain is a simplification of the bubble reactor consisting of a 2D column. The width of the system is 900 mm and includes a distributed grid located at 200 mm from the bottom side of the reactor, with channels of 18 mm that restrict the air flow at the reactor. In addition, the multiple flow regime [322] topology is an extended definition where all types of multiphase flow, that is, dispersed (bubble, droplet or particle) and stratified (free-surface) flows, can be modelled. The solver employs a volume-fraction based criterion to classify the flow regime into dispersed and separated two-phase flow. In this scenario, the multiple flow regime model treats bubbly flow in the first dispersed regime, droplet flow in the second dispersed regime, and separated water-air flow in the largescale interface regime.

Computational Domain
The computational domain is a simplification of the bubble reactor consisting of a 2D column. The width of the system is 900 mm and includes a distributed grid located at 200 mm from the bottom side of the reactor, with channels of 18 mm that restrict the air flow at the reactor.  The mesh is performed by means of a Cartesian-based grid with different embedding local refinements near the region of the grid channels where higher gradients of velocity are found. Furthermore, a mesh independence study was carried out to ensure accuracy of the CFD calculations by changing the reference cell size from 3 mm to 25 mm (minimum embedding cell size 0.75 mm and 6.25 mm respectively). Figure 7 shows the contours of volumetric fraction of water at t = 8 s obtained with the different meshes. The mesh is performed by means of a Cartesian-based grid with different embedding local refinements near the region of the grid channels where higher gradients of velocity are found. Furthermore, a mesh independence study was carried out to ensure accuracy of the CFD calculations by changing the reference cell size from 3 mm to 25 mm (minimum embedding cell size 0.75 mm and 6.25 mm respectively). Figure 7 shows the contours of volumetric fraction of water at t = 8 s obtained with the different meshes. Coarser meshes show less interaction between phases and the air appears to be diffused in the water. By contrast, finer meshes detect more heterogeneous flow in the reactor, and different scale of bubbles can be observed. Furthermore, these contours show that water column obtained by finer meshes tend to be lower.
This statement is confirmed when observing the temporal variation of the water column position (Figure 8a). In this Figure, ∆z = 0 m means initial position at t = 0 before air injection starts. After 3 s it is evident that results provided by finer meshes lead to final lower water level ∆z = 0.25 m. Nonetheless, it is observed that water level trends to decrease after 5 s which should not be expected, since the operation conditions and the height of the reactor are adequate to avoid the water arriving at the top boundary condition. Figure 8b shows the water mass imbalance in the computational domain after 20 s for every cell size mesh. The CFD solver produces a numeric mass imbalance error that increases when reducing the computational cell size. These statements are in accordance with those from other authors [36,37]. In multiphase flows the computational cell volumes have to be greater than the volume of the liquid for the Eulerian continuum approach to provide reasonable results. Therefore, an equilibrium between both features has been reached to achieve the most accurate possible solution. The cell size used in the solver mesh is 18 mm, which is considered a coarse mesh but without reaching cells too small to avoid surpassing 1% of mass imbalance. Coarser meshes show less interaction between phases and the air appears to be diffused in the water. By contrast, finer meshes detect more heterogeneous flow in the reactor, and different scale of bubbles can be observed. Furthermore, these contours show that water column obtained by finer meshes tend to be lower.
This statement is confirmed when observing the temporal variation of the water column position (Figure 8a). In this Figure, ∆z = 0 m means initial position at t = 0 before air injection starts. After 3 s it is evident that results provided by finer meshes lead to final lower water level ∆z = 0.25 m. Nonetheless, it is observed that water level trends to decrease after 5 s which should not be expected, since the operation conditions and the height of the reactor are adequate to avoid the water arriving at the top boundary condition. Figure 8b shows the water mass imbalance in the computational domain after 20 s for every cell size mesh. The CFD solver produces a numeric mass imbalance error that increases when reducing the computational cell size. These statements are in accordance with those from other authors [36,37]. In multiphase flows the computational cell volumes have to be greater than the volume of the liquid for the Eulerian continuum approach to provide reasonable results.
Therefore, an equilibrium between both features has been reached to achieve the most accurate possible solution. The cell size used in the solver mesh is 18 mm, which is considered a coarse mesh but without reaching cells too small to avoid surpassing 1% of mass imbalance. Coarser meshes show less interaction between phases and the air appears to be diffused in the water. By contrast, finer meshes detect more heterogeneous flow in the reactor, and different scale of bubbles can be observed. Furthermore, these contours show that water column obtained by finer meshes tend to be lower.
This statement is confirmed when observing the temporal variation of the water column position (Figure 8a). In this Figure, ∆z = 0 m means initial position at t = 0 before air injection starts. After 3 s it is evident that results provided by finer meshes lead to final lower water level ∆z = 0.25 m. Nonetheless, it is observed that water level trends to decrease after 5 s which should not be expected, since the operation conditions and the height of the reactor are adequate to avoid the water arriving at the top boundary condition. Figure 8b shows the water mass imbalance in the computational domain after 20 s for every cell size mesh. The CFD solver produces a numeric mass imbalance error that increases when reducing the computational cell size. These statements are in accordance with those from other authors [36,37]. In multiphase flows the computational cell volumes have to be greater than the volume of the liquid for the Eulerian continuum approach to provide reasonable results. Therefore, an equilibrium between both features has been reached to achieve the most accurate possible solution. The cell size used in the solver mesh is 18 mm, which is considered a coarse mesh but without reaching cells too small to avoid surpassing 1% of mass imbalance.

Phase Interaction Models
Both phase interaction models have been also tested in this work to determine their feasibility and accuracy for this bubbles reaction system. The reference operation condition is selected with the initial water column of 0.5 m, and air velocity through the reactor 0.1 m/s. Figure 9 shows the volume fraction contours of water at different time steps from 0 to 20 s. Both phase interaction models (presented in Section 4.1) provide significantly different results after 2 s. On the one hand, the continuous phase interaction model predicts monotonic increase of the water level, compared to the stable solution observed with the Multiple flow regime model [32]. The mixing zone is formed by big air bubbles and zones with low water volume fraction (average volume fraction through wet volume is 0.36) that leads to high water level for this point of operation, compared to results obtained with the multi-regime interaction model. On the other hand, the multi-regime interaction model shows a steady solution after 10 s, with distributed liquid zones (dark blue, average volume fraction through wet volume is 0.57).

Phase Interaction Models
Both phase interaction models have been also tested in this work to determine their feasibility and accuracy for this bubbles reaction system. The reference operation condition is selected with the initial water column of 0.5 m, and air velocity through the reactor 0.1 m/s. Figure 9 shows the volume fraction contours of water at different time steps from 0 to 20 s. Both phase interaction models (presented in Section 4.1) provide significantly different results after 2 s. On the one hand, the continuous phase interaction model predicts monotonic increase of the water level, compared to the stable solution observed with the Multiple flow regime model [32]. The mixing zone is formed by big air bubbles and zones with low water volume fraction (average volume fraction through wet volume is 0.36) that leads to high water level for this point of operation, compared to results obtained with the multi-regime interaction model. On the other hand, the multi-regime interaction model shows a steady solution after 10 s, with distributed liquid zones (dark blue, average volume fraction through wet volume is 0.57).

Phase Interaction Models
Both phase interaction models have been also tested in this work to determine their feasibility and accuracy for this bubbles reaction system. The reference operation condition is selected with the initial water column of 0.5 m, and air velocity through the reactor 0.1 m/s. Figure 9 shows the volume fraction contours of water at different time steps from 0 to 20 s. Both phase interaction models (presented in Section 4.1) provide significantly different results after 2 s. On the one hand, the continuous phase interaction model predicts monotonic increase of the water level, compared to the stable solution observed with the Multiple flow regime model [32]. The mixing zone is formed by big air bubbles and zones with low water volume fraction (average volume fraction through wet volume is 0.36) that leads to high water level for this point of operation, compared to results obtained with the multi-regime interaction model. On the other hand, the multi-regime interaction model shows a steady solution after 10 s, with distributed liquid zones (dark blue, average volume fraction through wet volume is 0.57).

Results and Validation
To check the validity of the CFD results shown above, several water-level measurements were performed on the prototype of the humidifier presented in this paper. Different tests were carried out at inlet gas velocities around the reference operation point (V_gas_inlet = 0.10 m/s) to obtain the level of the steady water column. Lower gas velocities lead to surge on the compressor that feeds the reactor, while higher gas velocities lead to choke. Figure 11 represents-in dots-the value obtained in the experiments, and the comparison with those results from the CFD calculations, both for the continuous and multi-regime interaction models, respectively. The water level predicted by the Continuous interaction models does not capture properly the behavior of the bubble reactor humidifier prototype in the range of operation. With very low gas velocity (below 0.05 m/s) only, both models provide equivalent results.

Results and Validation
To check the validity of the CFD results shown above, several water-level measurements were performed on the prototype of the humidifier presented in this paper. Different tests were carried out at inlet gas velocities around the reference operation point (V_gas_inlet = 0.10 m/s) to obtain the level of the steady water column. Lower gas velocities lead to surge on the compressor that feeds the reactor, while higher gas velocities lead to choke. Figure 11 represents-in dots-the value obtained in the experiments, and the comparison with those results from the CFD calculations, both for the continuous and multiregime interaction models, respectively. The water level predicted by the Continuous interaction models does not capture properly the behavior of the bubble reactor humidifier prototype in the range of operation. With very low gas velocity (below 0.05 m/s) only, both models provide equivalent results. As the multi-regime interaction model is the one that better fits the experimental data, it is the one that has been used to further the design of the bubble reactor. Figure 12 is showing the results of studying the effect of the air stream velocity and the initial water column height over the increase of the water column height. In Figure 12a, different air current velocity have been calculated while keeping a given initial water column height of 0.5 m. Similarly, in Figure 12b different initial water column height have been calculated keeping the same air current velocity of 0.1 m/s. As the multi-regime interaction model is the one that better fits the experimental data, it is the one that has been used to further the design of the bubble reactor. Figure 12 is showing the results of studying the effect of the air stream velocity and the initial water column height over the increase of the water column height. In Figure 12a As can be seen in Figure 12a, the air stream velocity must be controlled as the water column height increases exponentially when increasing the velocity. This parameter is the deciding factor for the minimum footprint area of the bubble reactor, to keep the air velocity in reasonable values when the maximum mass flow rate specified for the applica- As can be seen in Figure 12a, the air stream velocity must be controlled as the water column height increases exponentially when increasing the velocity. This parameter is the deciding factor for the minimum footprint area of the bubble reactor, to keep the air velocity in reasonable values when the maximum mass flow rate specified for the application is passing through the reactor (0.25 kg/s).
In the case of the initial water column height, it mainly affects the residence time of the air inside the water causing a direct effect on the final humectant capacity of the installation. Hence, it should be as high as possible, but as shown in Figure 12b, has a direct effect on the final height of the water column, so a compromise solution must be reached.
In the case of the application described in this paper, a 2D column simplification is proposed as a first step. Afterwards, the results with different phase models are studied to see the effect on the final water column height of each one. Then, with the help of an experimental facility the models are validated and the one that fits the experimental data better is selected. In the next step, the parametric studies for both the initial water height and air stream velocity are carried out, from which the optimal values for the design to avoid overflowing are 0.12 m/s of air stream velocity and 0.5 m of initial water column height.
With these two parameters the final geometry of the bubble reactor is obtained. The minimum footprint is calculated as the minimum area that assures that the air stream velocity at the outlet the bubble reactor does not increase above 0.12 m/s for the design point of the installation. In the case of the study presented in this paper, the footprint obtained was 0.81 m 2 for a design point of 2.4 bar, 50 C and 0.25 kg/s. Also, the reactor height was selected as 1.3 m to consider a security margin of 0.4 m over the nominal stabilized water column height of 0.9 m (initially 0.5 m with an increase of 0.4 m due to the 0.12 m/s of air velocity), in order to avoid overflows due to unexpected perturbances that could increase the bubbles height punctually.
Lastly, these design parameters are used on a 3D CFD model calculation of the final system to validate the design obtained, considering a wide range of operation points around the design point. In Figure 13, the results of the 3D calculations are shown for a given operation point. On the one hand, Figure 13a shows point of lower mass flow of the design point (0.2 kg/s), in which the final water column height is kept at the optimal height of 90 cm. This means that the simplification of a 2D column represents with enough accuracy the 3D behavior, validating the pre-design with the simplified geometry. On the other hand, Figure 13b shows how when increasing the mass flow way over the design (0.42 kg/s), the final water column height increases above the optimal height, but does not reach overflowing conditions, thus validating the reactor height selected previously.

Conclusions
Along the development of this paper, some conclusions came to light with respect to the calculation by means of CFD code of this kind on a bubbling case.
In the first place, it has been shown how the continuous phase interaction model overestimates the increase in height that the water can reach when an air current bubbles through, which leads to an oversized solution of the reactor obtained. Nevertheless, this

Conclusions
Along the development of this paper, some conclusions came to light with respect to the calculation by means of CFD code of this kind on a bubbling case.
In the first place, it has been shown how the continuous phase interaction model overestimates the increase in height that the water can reach when an air current bubbles through, which leads to an oversized solution of the reactor obtained. Nevertheless, this phase interaction model can be used to obtain a valid geometry for bubble reactors that avoids problems of overflowing, as shown in this paper. However, it must be noted that with this model it is not possible to design the optimal geometry for the application that is being studied.
By contrast, the multi-regime interaction model can predict accurately the behavior of the air bubbling and the water level inside the reactor for the range of air velocity tested. Therefore, it can be concluded that this is the phase interaction model that should be used to design equipment for the application described in this paper.
It should be noted that the cell size around the phase interaction is crucial for achieving a correct solution. If the cell size is too small there would be increased errors in the calculations, i.e., higher ratio of mass imbalance, and if it is too big both phases will fade which would reduce the resolution and accuracy of the calculation. Therefore, for applications like those depicted in this paper, the cell size for the interphase calculation should be around 18 mm to achieve an accurate result.
Lastly, to optimize the design it is necessary to fit the footprint to the mass flow range used in the application, in order to keep the air velocity through the air controlled. Also, a compromise should be reached on the initial water column height considering the maximum water content in the air current expected and the maximum total height available for the bubble reactor. In the case of the application described in this paper, an initial water column height of 0.5 m, footprint of 0.81 m 2 and a reactor height of 1.3 m have been adopted as a final design.