Do increased flow rates in displacement ventilation always lead to better results?

Abstract Indoor ventilation is essential for a healthy and comfortable living environment. A key issue is to discharge anthropogenic air contamination such as CO$_2$ gas or, of potentially more direct consequence, airborne respiratory droplets. Here, by employing direct numerical simulations, we study mechanical displacement ventilation with a wide range of ventilation rates $Q$ from 0.01 to 0.1 m$^3$ s$^{-1}$ person$^{-1}$. For this ventilation scheme, a cool lower zone is established beneath a warm upper zone with interface height $h$, which depends on $Q$. For weak ventilation, we find the scaling relation $h\sim Q^{3/5}$, as suggested by Hunt & Linden (Build. Environ., vol. 34, 1999, pp. 707–720). Also, the CO$_{2}$ concentration decreases with $Q$ within this regime. However, for too strong ventilation, the interface height $h$ becomes insensitive to $Q$, and the ambient averaged CO$_2$ concentration decreases towards the ambient value. At these values of $Q$, the concentrations of pollutants are very low and so further dilution has little effect. We suggest that such scenarios arise when the vertical kinetic energy associated with the ventilation flow is significant compared with the potential energy of the thermal stratification.


Introduction
Ventilation is important to optimise the air quality in indoor spaces (Linden 1999;Awbi 2003;Morawska et al. 2020). There are many reasons for achieving effective ventilation, such as thermal comfort, energy savings and to minimise pollutants such as carbon dioxide (CO 2 ) or aerosols contributing to the spread of infectious diseases (Fisk 2000;Li et al. 2007;WHO 2020;CDC 2021). Recently, indoor ventilation has become increasingly important considering that the majority of COVID-19 infections are caused by asymptomatic or presymptomatic infected individuals, and that, by speaking and breathing, tiny saliva droplets are released into the environment, accumulate (Abkarian et al. 2020;Pöhlker et al. 2021) and persist in the air for long periods (Bourouiba 2020(Bourouiba , 2021Chong et al. 2021). To measure the amount of respiratory contaminants in a room, one often relies on the fact that tiny droplets, aerosols or droplet nuclei (i.e. dried-out droplets) behave like tracers and are simply carried along by the air flow. Therefore, they evolve similarly to another tracer produced through exhalation, namely CO 2 . Since CO 2 is easily measured, it is often regarded as a surrogate for indoor airborne contaminants (von Pohle, Anholm & McMillan 1992;Rudnick & Milton 2003).
There are two limiting ventilation types (Linden 1999;Awbi 2003;Chenvidyakarn 2013;. First, for mixing ventilation, air is circulated throughout the room, leading to a uniform distribution of temperature and contaminants (Linden 1999), and the purpose of this ventilation type is usually thermal comfort. The second type is displacement ventilation, which has the benefit of displacing contaminants to the upper level of the room where they are extracted, so that occupants breathe fresher air . In this case, the air inlet is at a low level and the air outlet is at a high level. This arrangement results in a flow through the room interior from bottom to top, with the exhaled droplets or aerosols being discharged upwards. The flow in the interior is stratified: each occupant produces an upward thermal plume (the typical heat output of a person at rest is 80 W) due to body heat (Craven & Settles 2006;, which leads to a density stratification, with the warm, lighter air layer above the cooler, heavier counterpart (see e.g. figure 1a). Respiratory contaminants exhaled together with a warm vapour puff that are not entrained by the main body plume tend to accumulate between the stratified air layers (Bolster & Linden 2007; to form the so-called lock-up layer (Qian et al. 2006;Zhou et al. 2017).
Early studies of displacement ventilation were mainly focused on experiments and theories. Sandberg & Lindstrom (1990) were the first to discuss displacement ventilation with a source of buoyancy, which is referred to as natural displacement ventilation since there is no mechanical extraction. Linden, Lane-Serff & Smeed (1990) further investigated this phenomenon in an enclosure with inflow at the bottom and outflow at the top. They demonstrated the presence of thermal stratification with two layers of uniform temperatures separated by a horizontal interface (see e.g. figure 1). They also derived how the layer height h of the lower cold zone depends on the effective area of the vent. The stratification is also affected by the number of occupants. With separately distributed buoyancy sources, multiple layers can form, leading to a more gradual change of temperature in the domain (Cooper & Linden 1996;Livermore & Woods 2007).
When the displacement ventilation is driven by mechanical extraction or an input wind (Hunt & Linden 1999, 2001, 2005, the interface height h is obtained from matching the total extraction rateQ and the buoyancy sources, giving the hallmark formula: Equation (1.1) is obtained by scaling arguments of a self-similar buoyant plume in an unstratified background (Schmidt 1941;Morton, Taylor & Turner 1956;Linden et al. 1990). Here, B is the buoyancy flux produced by each of n occupants and h v is the 'virtual origin' of the thermal plume, i.e. the height at which the plume would start if it were a point buoyant plume. The empirical constant c has an approximate value of c ≈ 0.105 (Morton et al. 1956;Linden 1999). Given that the formula is based on some simple scaling arguments and is derived using some assumptions, one wonders how well the clean zone height h prediction works for practical extraction rates and in what range it holds. For example, when the inflow is too strong, the interface can break down and lead to a short-circuit flow between the inlet and outlet (Partridge & Linden 2017).
In this study, we use direct numerical simulation (DNS) to study the mechanical displacement ventilation with a wide and practically relevant range of ventilation rate Q. To reproduce indoor air flows, the flow is fully coupled to the temperature, and to the CO 2 and water vapour concentration fields. Previous work on indoor ventilation has mainly focused on the design for energy efficiency, whereas, due to the COVID-19 pandemic, the recent objective has been on how ventilation can effectively remove the respiratory contaminant.
Previous numerical studies mainly employed large-eddy simulation (LES) or turbulence models (Davidson 1989;Durrani, Cook & McGuirk 2015;van Hooff, Blocken & Tominaga 2017), while numerical results obtained by DNS are scarce. However, LES and turbulence models perform differently in different air flow cases (Zhang et al. 2007). Simulating flows in the near-wall region is more challenging than those in the free flow region (Piomelli & Balaras 2002) due to the dominance of smaller vortices. LES also becomes more questionable once temperature and CO 2 are coupled, since with LES the turbulent heat flux can be significantly underestimated (Taylor & Sarkar 2008). Thus the dataset from our DNS can be seen as a benchmark for comparison with results by LES or turbulence models, which are unavoidable for large parameter space studies.
The paper is organised as follows. In § 2, the flow set-up and the governing equations are introduced. In § 3, we discuss how the lock-up effects change with increasing Q ( § 3.1), explain the change in the layer formation by an energy balance ( § 3.2), analyse globally and locally averaged concentration values ( § 3.3), and discuss the influence of different ventilation set-ups on the lock-up effect ( § 3.4). Finally, in § 4, we summarise our findings and provide an outlook for future work.

Flow set-up
Owing to computational constraints, we consider a room of dimensions 3 m × 3 m × 3 m, with a single occupant, having an idealised lower-level inlet flow and higher-level outlet flow on opposite walls, as sketched in figure 1(a). For simplicity, both inlet and outlet have the same heights (0.3 m) and widths (3 m). The ambient temperature, relative humidity (RH) and CO 2 volume fraction of the room are considered fresh and homogeneous initially and set to values acceptable for indoor comfort, i.e. T amb = 22 • C, 40 % and 0.04 %, respectively. Correspondingly, the inlet flow is a uniform flow and set to be the same as the initial ambient conditions, whereas an outflow boundary condition is imposed at the outlet. For the walls, no-slip, impermeable and adiabatic boundary conditions are imposed.
A seated occupant within the room produces thermal buoyancy (from body temperature), high vapour concentration (from exhalation) and respiratory contaminants (represented by CO 2 ). The fluid is assumed to be incompressible, ∂ i u i = 0. To capture complete interactions of these anthropogenic sources, the flow is fully coupled with the heat, water vapour concentration and CO 2 concentration by applying the Boussinesq approximation. The governing equations read The scalar variable φ represents either temperature T, vapour mass fraction γ , or CO 2 volume fraction c, all of which are simultaneously taken into consideration. We define ρ ref as the reference density, β φ as the isobaric thermal expansion coefficient of the fluid for each scalar, ν is the kinematic viscosity and κ φ is the diffusivity for each scalar, for example, κ T is the thermal diffusivity. All of these quantities are assumed to be independent of temperature. Accordingly, the gas Prandtl number is Pr ≡ ν/κ T = 0.71. The thermal Rayleigh number of the occupant, Ra ocpt ≡ gβ T Δh 3 ocpt /(νκ T ) = 8.8 × 10 8 , where Δ ≡ T ocpt − T amb and h ocpt = 1.2 m. The grid resolutions are 768 3 , with a clipped Chebyshev clustering in the z-direction to resolve gradients at the upper and lower walls. Our resolutions have been validated by performing grid refinement studies and by comparing with the body thermal plume velocity profile from experiments (Craven & Settles 2006) (see the supplementary material and movies, which are available at https://doi.org/10.1017/jfm.2021.949). Equations (2.1a,b) are solved by DNS using the multi-scalar second-order finite-difference method with a fractional third-order Runge-Kutta scheme (Verzicco & Orlandi 1996;Ostilla-Mónico et al. 2015). The body in the room is modelled by the immersed boundary method with no-slip boundary, fixed body temperature and zero CO 2 and vapour concentrations.
The interactions are expected to be non-trivial because of the competing effects of hot humid air, which is lighter than the ambient, with the heavier CO 2 . The occupant body temperature is set as T ocpt = 27 • C (Houdas & Ring 2013), which is the mean of the periphery and ambient temperature. The exhaled breath is set to 27 • C (given that the temperature drops from the mouth at 37 • C), with RH = 100% and CO 2 volume fraction of 4 %. The tidal period and volume of the breath is set to 4 s and 2 l, respectively, matching experimental measurements (Gupta, Lin & Chen 2010). Since our focus is on the large-scale flows in the room, and also to make our simulations tractable, the breath is modelled as a Gaussian spatial source at a distance of ≈20 cm from the mouth and angled 60 • below horizontal. The simulations were performed for at least 500 free-fall time units (equivalent to 30 min in physical time), and we sample the last 200 free-fall units for our statistical analysis. To ensure that the steady state is reached, we compared the last 100 free-fall units to the full sampling duration, and the variation is ∼1 % difference. A visualisation of the flow field is available in the supplementary movie. The main control parameter of the simulations is the ventilation rateQ per person (in the following, p −1 represents per person), which we denote as Q (having units m 3 s −1 p −1 ). We chose to use Q as the ventilation unit instead of air change rate (commonly with units of air changes per hour, ACH), since Q includes the number of occupants as the key factor. For reference, 1 ACH = 0.0075 m 3 s −1 p −1 or, equivalently, an average inlet speed is 0.0083 m s −1 in our set-up. The range of Q is simulated from 0.01 to 0.1 m 3 s −1 p −1 corresponding to ACH from 1.3 to 13, which adopts the ventilation rate in .

Influence of Q on the lock-up effect
We begin by analysing the clean-zone height, h, defined as the height below which the air maintains the same properties as at the inflow (see figure 1). Typically, clean zones are designed to be taller than occupants to ensure a contaminant-free occupied space.
In the presence of thermal sources, a stably stratified environment can exist within a room (Bolster & Linden 2007) where lower, cooler air is separated by upper, warmer air by an interface at height h. Contaminants such as aerosols and fomites that are entrained by the thermal sources can accumulate at this interface (Bolster & Linden 2007), and this is commonly referred to as the lock-up effect (Qian et al. 2006;Zhou et al. 2017). To illustrate this point, in figure 2, height h, the mean temperature profiles and CO 2 concentrations are plotted for various Q.
In figure 2(a), we find that the interface height h increases from Q = 0.01 to 0.05 m 3 s −1 p −1 . This increasing trend is in remarkably good agreement with (1.1), shown as the dashed line in figure 2(a), indicating that the thermal body plume can be well approximated by assuming a simplified buoyant plume source in accordance with Morton et al. (1956). However, when Q 0.05 m 3 s −1 p −1 , relation (1.1) no longer holds and h remains roughly constant instead. This implies that the simplified assumptions cannot be extended to the largest ventilation rates in our set-up. We will see later that this is because of the dominance of the inflow resulting in vastly different flow structures in the room. We call this large-Q regime the inflow-dominant regime. The weak influence of h on Q for high Q values is also reflected in the temperature and CO 2 profiles (figure 2b,c), which remain roughly unchanged for Q 0.05 m 3 s −1 p −1 . To further test the robustness of these results, we conducted three additional simulations with two occupants in the room (shown as crosses in figure 2a), which are also in agreement with the general trend. It should be noted that the transition point of ventilation rate may change as the number of occupants increases, since the potential energy depends on the number of occupants. It also depends, of course, on the vent size and position, the volume of the room and the geometry. These dependences are not analysed in the present work and deserve further studies.
Using the CO 2 concentration as an indicator of respiratory contaminants (von Pohle et al. 1992), in figure 2(c), we find that the lock-up effect is evident for Q 0.05 m 3 s −1 p −1 , as shown by the peaks in the CO 2 concentration. At Q 0.05 m 3 s −1 p −1 , the CO 2 concentration is substantially reduced as compared to Q = 0.01 m 3 s −1 p −1 , suggesting that the contaminant layer depletes at large Q and as expected for mechanical displacement ventilation . However, although the lock-up layer diminishes at large Q, there is a modest increase in CO 2 concentration close to head height at z ≈ 1.2 m (see inset in figure 2c). This implies that, at Q 0.05 m 3 s −1 p −1 , contaminants become dispersed at lower heights, which contributes to poorer air quality. The findings from this simplified set-up highlight that careful design of indoor ventilation configuration is necessary to ensure efficient removal of contaminants. As stated above, the transition point of the flow rate Q also depends, of course, on other factors, such as the room geometry and layout, the vent size and position.

Regime transition explained by balance of potential and kinetic energy
We next explain why there is transition to the inflow-dominant regime by examining the relative strength of the energies in the system. As depicted in figure 3(a), three types of energies (per unit mass) can be identified: (i) kinetic energy, E k = 1 2 u 2 ; (ii) potential energy due to the stable stratification, E p ≡ N 2 h 2 layer ; and (iii) energy losses from friction and blockage, E loss . Here, u is the inflow velocity, N ≡ √ β T g|dT/dz| max is the buoyancy frequency, with β T the thermal expansion coefficient and g the gravitational acceleration, and h layer is the effective thickness of the upper layer (with height of the outlet vent w being subtracted) expressed as h layer ≡ H − w − h with h the height of the clean zone. Thus, the potential energy can be rewritten as  The energy loss in (iii) consists of two parts, which are the friction from the ground, E friction , and the drag due to the blockage of the occupant, E block . We estimate E friction by 1 2 C f u 2 with the friction coefficient for turbulent flow Pope 2000), where Re y ( y) ≡ uy/ν. The other part of the energy loss is E block , which is expressed as E block = (α/2)C d u 2 . Here, the geometrical factor α takes into account the partial blockage by the occupant (estimated by the product of the cross-sectional area of the occupant and the height of the occupant over the domain volume) and its value is approximately 0.06. The drag coefficient C d ≈ 0.7 (Wang, Zhou & Mi 2012). The dissipated energy is two orders of magnitude smaller than the kinetic energy itself and can thus be neglected.
In figure 3(b), we compare the residual kinetic energy E k to the potential energy E p for various Q. When Q < 0.05 m 3 s −1 p −1 the inflow is still relatively weak compared to the strong thermal stratification established in the room, and therefore E p is dominant. Upon increasing Q, the influence of the inflow becomes more dominant. At around Q = 0.06 m 3 s −1 p −1 , the residual kinetic energy becomes larger and overtakes E p . This energy argument explains the transition to the inflow-dominant regimes in figure 2. As can be seen in figure 2(c), the layer height remains the same for large enough Q. The physical explanation is that the strong inflow breaks the stratified layer near the outlet and directly leaves the room (see figure 3a), which has been shown experimentally by Partridge & Linden (2017). Therefore, in our set-up the excessive kinetic energy no longer contributes to making the upper stratified layers thinner. Crucially, the particular value of the ventilation rate at which this 'short-circuiting' occurs is related to the kinetic energy of the inflow via the vent size and further depends on, at least, the room layout and geometry, and the position and orientation of the vents.
3.3. Globally and locally averaged temperature and CO 2 concentration While spatial profiles are useful for quantifying the interface height h and the lock-up effect, it is also instructive to compute integral-averaged quantities, since they represent measurable metrics for real-world cases. Therefore, in figure 4, we plot the spatially averaged temperature and CO 2 concentrations versus Q. The plane averaging is performed on the longitudinal plane and the region near the occupant is omitted. We consider four different ways of averaging; see the inset sketches in figure 4(a). From figure 4, in all these cases the averaged temperature and CO 2 concentration eventually decrease towards the ambient values with increasing Q. In particular, local and pointwise measurements are reduced to ambient values when Q 0.02 m 3 s −1 p −1 . In contrast, the global values are larger than the ambient values and reach a plateau at Q 0.05 m 3 s −1 p −1 . These trends highlight two crucial points. Firstly, measurements at head height can provide an inadequate view of the temperature and CO 2 concentration values (Mahyuddin & Awbi 2012). For our case, the inaccuracy is mainly caused by the non-trivial flow organisation, where a stably stratified background flow is separated by a temperature interface. Secondly, increasing Q above 0.05 m 3 s −1 p −1 does not alter the global temperature value and CO 2 concentration since they are almost at their ambient values, which further supports the fact that a transition point of Q exists for mechanical displacement ventilation if short-circuiting occurs.

Influence of different flow configurations
To test the sensitivity to different inflow set-ups, we simulated additional cases with air inflow and outflow on the same side (see figure 5f -j). As seen from these, the stratified temperature field for this case and the original case (figure 5a-e) are quite similar (figure 5a,f ). Just as before, the region of high CO 2 concentration is closely related to the inflection point in the mean temperature profile T(z) (indicated, for example, by the isocontour of 500 ppm in figure 5b,g). On closer inspection, some key differences can be identified. For example, h is much lower in the one-sided case (≈1.9 m, see circled marker in figure 5j) and as a result the CO 2 concentration is also higher at that z-level. Some differences between the profiles are indeed expected since the large-scale flow structures are clearly different for the two-sided and one-sided ventilation cases. In particular, when we inspect the vertical velocity in the periphery above the body (see figure 5i), we find that there is a significant upward velocity contribution from the inflow. This upflow  contribution influences the effective buoyant flux from the body through mixed convection, causing a stronger temperature mixing and thus lowering the interface height h.

Conclusions and outlook
In summary, we have shown and discussed a regime transition of mechanical displacement ventilation by studying a wide range of Q. When the ventilation rate is low, the observed flow structure agrees with the classical stratified flow structure of displacement ventilation, where there is a clean and cooler zone established beneath the warm contaminated zone. In this case, the height h for the lower clean zone scales with Q as h ∼ Q 3/5 , cf. (1.1). However, when the ventilation rate becomes high , we observe the inflow-dominant regime in which h becomes insensitive to Q. There is a major change in flow structure where the stratified layer partially breaks down due to the kinetic energy of the flow exceeding the potential energy of the stratification. By measuring the ambient mean CO 2 concentration, we find that the local and global concentrations remain unchanged with increasing Q in the inflow-dominant regime, implying that too large ventilation strengths (i.e. too large flow rates Q) do not help to remove more contaminants. Our work gives clear evidence that, for the cross-flow displacement ventilation studied, the benefits of the displacement ventilation strategy break down for too strongly forced inflows. In such cases part of the ventilating flow passes to the outlet without being entrained into the upper layer via thermal convection. We call this mechanism 'short-circuiting'. Once the ventilation exhibits such short-circuiting, increasing the interface height and the further removal of pollutants from the upper layer both become more challenging.
In the presence of internal heat loads, one expects that the potential energy, or stack pressure, associated with the stratification is decreased as the ventilation flow is increased. For vents of fixed size, this is reflected in an increase of the kinetic energy. We have found no prior evidence of the energy balance being expressed in this manner in the literature, but from our point of view this is a useful approach. For the cases that we examined, once the flow short-circuited, further increasing the kinetic energy of the driving does not result in an enhanced potential energy of the stratification. Moreover, the point at which short-circuiting occurred approximately coincided with the crossing point between the kinetic energy of the inflow and the potential energy of the stratification. We hope that the energy balance analysis as done in this paper will also be applied to other ventilation flows, for example to examine when the transition between a displacement flow and a 'well-mixed' state will occur. Clearly, the essence of our findings will also hold for other ventilation scenarios and geometries.
Motivated by mitigating the COVID-19 pandemic, many studies on indoor ventilation have emerged recently Morawska & Milton 2020;Morawska et al. 2020). Given that the parameter space is vast, it is computationally demanding to obtain a generalised view on indoor ventilation. Nevertheless, by examining the energetic balance in the room with simple formulae, one can gain insight into the transition point of Q for different scenarios, such as cases with different-sized rooms, number of occupants, etc. Note that, in our set-up, we chose idealised inlet and outlet vent shapes and scheme to simplify the energy balance argument. The value of the transition point from the low-Q regime to the high-Q regime may depend on the geometry and other factors. Further work using LES and other turbulent models are more suitable for exploring the wide parameter space, and are crucial to further generalise the methodology proposed in our work. Our DNS results of ventilation flows can be used for validations of these simulation methods. Lastly, other factors could also influence our results and merit further investigation, such as when occupants are moving between rooms (Mingotti et al. 2020), or when occupants are performing different respiratory activities (Pöhlker et al. 2021).