Mathematical Modelling of Face Coverings for Virus Protection

TL;DR

This paper develops a mathematical model to analyze airflow and aerosol transport through different face coverings, aiming to optimize mask permeability for better virus protection.

Contribution

It introduces a novel mathematical model for porous face masks, comparing traditional masks and neck gaiters, and identifies optimal permeability for protection.

Findings

Neck gaiters have lower leakage than traditional masks.

An optimal permeability exists that maximizes protection.

The model provides insights into airflow and aerosol dynamics through masks.

Abstract

Traditional face masks used extensively, for example, during the COVID-19 pandemic, are relatively impermeable, forcing aerosol-containing breath to flow into the surrounding atmosphere by ``leaking" out from the lateral edges of the mask. This protects close-proximity persons but results in a high viral load being delivered into the air. An alternative face covering, the Virustatic Shield, developed by Virustatic, involves a much more permeable material coated in an active agent, with the goal of increasing the transmission through the mask and thus reducing the viral load released. In this paper, we build and solve a model that describes the air flow and aerosol transport through a porous face mask. Our model for air flow is based on a low-Reynolds-number approximation to the Navier--Stokes equations, and reduces to a two-dimensional scalar partial differential equation for the local flux through the mask. We use the model to explore the interplay between leakage from the mask sides and flow directly through. We compare and contrast the leakage for traditional face masks and for neck gaiters, and find that neck gaiters have a lower leakage. We utilise a simple model for aerosol transport that balances advection with capture by the mask, and find that there is an ideal permeability, which optimises real-world protection offered by the mask.

Tables

Table 1. Table 1: Typical values of dimensional parameter values in the problem. We will use the reference values for all parameters except gap width for the remainder of the paper. Separate values are given for the permeabilities of different face covering materials.

Parameter		Range	Reference Value	Source
Mask width	$\hat{a}$		$1.8\text{×}{10}^{-1}m$
Neck Gaiter width	$\hat{a}$		$5.4\text{×}{10}^{-1}m$	111Whilst worn, due to stretching
Mask height	$\hat{b}$		$9\text{×}{10}^{-2}m$
Mask thickness	$\hat{T}$		$4\text{×}{10}^{-4}m$	[11]
Mouth area	$\hat{A}$	$1.2\text{×}{10}^{-4}{m}^2\text{ to }1.8\text{×}{10}^{-4}{m}^2$	$1.2\text{×}{10}^{-4}{m}^2$	[36]
Mean gap thickness	${\hat{h}}_0$	$3.6\text{×}{10}^{-4}m\text{ to }1.1\text{×}{10}^{-3}m$	$4.9\text{×}{10}^{-4}m$	[24]222Mean perimeter gap thickness for face masks with different fits predicted by simulations.
Density of air	$\hat{\rho }$		$1.2\mathrm{kg}{m}^{-3}$	[24]333At $20\mathrm{°C}$.
Viscosity of air	$\hat{\mu }$		$1.8\text{×}{10}^{-5}\mathrm{Pa}s$	[24]333At $20\mathrm{°C}$.
Exhalation velocity	${\hat{w}}^{\ast }$	$0.3m{s}^{-1}\text{ to }0.7m{s}^{-1}$	$0.5m{s}^{-1}$	[7]
Mask permeability (surgical mask)	$\hat{K}$	$9.9\text{×}{10}^{-12}{m}^2\text{ to }1.4\text{×}{10}^{-11}{m}^2$	$1.0\text{×}{10}^{-11}{m}^2$	[24]
Mask permeability (neck gaiter)	$\hat{K}$	$2.4\text{×}{10}^{-11}{m}^2\text{ to }5.2\text{×}{10}^{-11}{m}^2$	$5.0\text{×}{10}^{-11}{m}^2$	Supplied by Virustatic
Virus denature constant	$\hat{\gamma }$		Unknown

Table 2. Table 2: Values of the dimensionless parameter values in the model based on reference values for the dimensional parameters in Table 1 and the given dimensional range of h ^ \hat{h} .

Dimensionless parameter	Formula	Value
$ϵ$	$\frac{{\hat{h}}_0}{\hat{b}}$	$4\text{×}{10}^{-3}\text{ to }1.2\text{×}{10}^{-2}$
$l$ (surgical mask)	$\frac{\hat{a}}{\hat{b}}$	$2.0$
$l$ (neck gaiter)	$\frac{\hat{a}}{\hat{b}}$	$6.0$
$T$	$\frac{\hat{T}}{{\hat{h}}_0}$	$0.36\text{ to }1.1$
${ϵ}^2\mathrm{Re}$	$\frac{\hat{\rho }\hat{A}{\hat{w}}^{\ast }{\hat{h}}_0}{{\hat{b}}^2\hat{\mu }}$	$0.18\text{ to }0.54$
$\Gamma$	$\frac{\hat{\gamma }\hat{b}}{\varphi }$	Unknown
$M$	$\frac{\hat{A}}{{\hat{b}}^2}$	$1.4\text{×}{10}^{-2}$
$K$ (surgical mask)	$\frac{\hat{K}{\hat{b}}^2}{{\hat{h}}^4}$	$5.5\text{×}{10}^{-2}\text{ to }4.8$
$K$ (neck gaiter)	$\frac{\hat{K}{\hat{b}}^2}{{\hat{h}}^4}$	$2.8\text{×}{10}^{-1}\text{ to }2.4\text{×}{10}^1$

Table 3. Table 3: Parameter values so that ( 71 ) approximates the shape of face f f and a face mask g g .

	${\lambda }_1$	${\lambda }_2$	${\lambda }_3$	${\lambda }_4$	${\lambda }_5$	${\lambda }_6$	${\lambda }_7$	${\lambda }_8$
${𝝀}^f$	$1.0$	$0.07$	$0.0$	$50.0$	$0.6$	$25.0$	$0.55$	$0.0$
${𝝀}^g$	$1.0$	$0.02$	$0.1$	$51.0$	$0.605$	$15.0$	$0.48$	$0.01$

Equations

$$\hat{\nabla }\cdot \hat{\mathbf{u}}$$

$$\hat{\rho }\left(\frac{\partial \hat{\mathbf{u}}}{\partial \hat{t}}+\hat{\mathbf{u}}\cdot \hat{\nabla }\hat{\mathbf{u}}\right)$$

$$\hat{\nabla }\cdot \left(\varphi \hat{\mathbf{U}}\right)$$

$$\varphi \hat{\mathbf{U}}$$

$$\hat{\mathbf{u}}\times {\mathbf{n}}_f=0,\hat{\mathbf{u}}\cdot {\mathbf{n}}_f={\hat{w}}_{\mathrm{in}}\left(\hat{x},\hat{y}\right)\text{ at }\hat{z}=\hat{f}\left(\hat{x},\hat{y}\right),$$

$${\hat{w}}_{\mathrm{in}}\left(\hat{x},\hat{y}\right)=\{\begin{array}{ll}{\hat{w}}^{\ast }& \text{ if }\left(\hat{x},\hat{y}\right)\in {\Omega }_m,\\ 0& \text{ if }\left(\hat{x},\hat{y}\right)\notin {\Omega }_m\end{array}$$

$$\hat{\mathbf{u}}\times {\mathbf{n}}_g=0,\hat{\mathbf{u}}\cdot {\mathbf{n}}_g=\varphi \hat{\mathbf{U}}\cdot {\mathbf{n}}_g,\hat{p}=\hat{P}\text{ at }\hat{z}=\hat{g}\left(\hat{x},\hat{y}\right).$$

$$\hat{P}={\hat{p}}_{\mathrm{atm}}\text{ at }\hat{z}=\hat{g}\left(\hat{x},\hat{y}\right)+\hat{T}\sqrt{1+{\left(\frac{\partial \hat{g}}{\partial \hat{x}}\right)}^2+{\left(\frac{\partial \hat{g}}{\partial \hat{y}}\right)}^2}$$

$$\hat{p}={\hat{p}}_{\mathrm{atm}},\hat{P}={\hat{p}}_{\mathrm{atm}}\text{ if }\hat{x}=0,\hat{a},\text{ or }\hat{y}=0,\hat{b}.$$

$$\hat{p}={\hat{p}}_{\mathrm{atm}},\hat{P}={\hat{p}}_{\mathrm{atm}}$$

$$\hat{\mathbf{u}},\hat{p},\hat{P}\text{ periodic }$$

$$\hat{\mathbf{u}}=0\text{ at }\hat{t}=0.$$

$$\frac{\partial \hat{c}}{\partial \hat{t}}+\hat{\mathbf{u}}\cdot \hat{\nabla }\hat{c}=0.$$

$$\frac{\partial \hat{C}}{\partial \hat{t}}+\varphi \hat{\mathbf{U}}\cdot \hat{\nabla }\hat{C}=-\frac{\hat{\gamma }\left|\hat{\mathbf{U}}\right|\hat{C}}{\varphi },$$

$$\hat{c}=\{\begin{array}{ll}{\hat{c}}_{\mathrm{in}}& \text{ if }\left(\hat{x},\hat{y}\right)\in {\Omega }_m\\ 0& \text{ if }\left(\hat{x},\hat{y}\right)\notin {\Omega }_m\end{array},\text{ at }\hat{z}=\hat{f}\left(\hat{x},\hat{y}\right),$$

$$\hat{c}\hat{\mathbf{u}}\cdot {\mathbf{n}}_g=\varphi \hat{C}\hat{\mathbf{U}}\cdot {\mathbf{n}}_g,\text{ at }\hat{z}=\hat{g}\left(\hat{x},\hat{y}\right).$$

$$\hat{c}=0,\hat{C}=0\text{ at }\hat{t}=0.$$

$$\left(\hat{x},\hat{y},\hat{z}\right)$$

$$\left(\hat{u},\hat{v},\hat{w}\right)$$

$$\nabla \cdot \mathbf{u}$$

$${ϵ}^2\mathrm{Re}\left(\frac{\partial u}{\partial t}+\mathbf{u}\cdot \nabla u\right)$$

$${ϵ}^2\mathrm{Re}\left(\frac{\partial v}{\partial t}+\mathbf{u}\cdot \nabla v\right)$$

$${ϵ}^2\mathrm{Re}\left(\frac{\partial w}{\partial t}+\mathbf{u}\cdot \nabla w\right)$$

$$\frac{\partial c}{\partial t}+\mathbf{u}\cdot \nabla c$$

$${ϵ}^2\left(\frac{\partial }{\partial x}\left(\varphi U\right)+\frac{\partial }{\partial y}\left(\varphi V\right)\right)+\frac{\partial }{\partial z}\left(\varphi W\right)$$

$$\varphi U$$

$$\varphi V$$

$$\varphi W$$

$$\frac{\partial C}{\partial t}+{ϵ}^2\left(\varphi U\frac{\partial C}{\partial x}+\varphi V\frac{\partial C}{\partial y}\right)+\varphi W\frac{\partial C}{\partial z}$$

$$K=\frac{\hat{K}}{{ϵ}^4{\hat{b}}^2},\Gamma =\frac{\hat{\gamma }\hat{b}}{\varphi }$$

$$\mathbf{u}\times {\mathbf{n}}_f=0,\mathbf{u}\cdot {\mathbf{n}}_f=w_{\mathrm{in}}(x,y),c=\{\begin{array}{ll}1& \text{ if }(x,y)\in {\Omega }_m\\ 0& \text{ if }(x,y)\notin {\Omega }_m\end{array}\text{ at }z=f(x,y),$$

$$\mathbf{u}\times {\mathbf{n}}_g=0,-\left(u\frac{\partial g}{\partial x}+v\frac{\partial g}{\partial y}\right)+w=-{ϵ}^2\left(\varphi U\frac{\partial g}{\partial x}+\varphi V\frac{\partial g}{\partial y}\right)+\varphi W,p=P,c=C\text{ at }z=g(x,y),$$

$$P=0\text{ at }z=g(x,y)+T\sqrt{1+{ϵ}^2{\left(\frac{\partial g}{\partial x}\right)}^2+{ϵ}^2{\left(\frac{\partial g}{\partial y}\right)}^2},$$

$$w_{\mathrm{in}}\left(x,y\right)=\{\begin{array}{ll}\frac{1}{M}& \text{ if }\left(x,y\right)\in {\Omega }_m\\ 0& \text{ if }\left(x,y\right)\notin {\Omega }_m\end{array},$$

$$p=0,P=0\text{ if }x=0,l,\text{ or }y=0,1,$$

$$p=0,P=0$$

Full text

I Introduction

During the COVID-19 pandemic, the use of face coverings greatly expanded [1] and evidence shows that they reduced the spread of COVID-19 and other respiratory diseases [2, 3, 4]. Face coverings may be broadly divided into three categories [5]. The highest protection is offered by respiratory protection devices (respirators), designed to provide a perfect seal between the covering and the face. Such devices are intended to protect the wearer from infection in high-risk environments such as COVID wards in hospitals, and require fit tests by each worker to ensure a particular design fits their face. Respirators are governed by regulations specifying a minimum filtration efficiency and a maximum allowed total inward leakage of air [5, 1]. Surgical (or medical) face masks offer the next highest level of protection and are also subject to regulated standards. These masks were originally intended for source control of large droplets, and so no specification of the amount of either inward or outward leakage is given [5, 1]. Finally, community face coverings (cloth masks, neck gaiters) grew in popularity in the wider population during the pandemic [1, 5, 6]. Regulations for community face coverings were introduced in 2020 due to the pandemic and are again primarily focused on source control, and to a lesser extent on the wearer’s protection [5]. However, evidence exists demonstrating that these masks leak a significant proportion of the exhaled air and particles [7, 8].

There has been considerable interest in understanding and characterising the key factors in the design and use of face coverings to optimise the protection they offer to the wearer and others. To quantify the performance of face masks, studies have focused on four key metrics, two of which relate to the material and two to the real-world performance [9]. The properties of the mask material are measured by the permeability (often reported as the pressure drop), which quantifies the ability of air to flow through the mask material, and the filtration efficiency given by the proportion of particles entering the mask material that are captured. The real-world performance of face masks is quantified by the leakage ratio, $\eta$, the proportion of air exhaled that does not pass through the mask material, and the face-fitted filtration efficiency, $\mathrm{FFE}$, which is the proportion of exhaled particles captured by the face covering. Note that these definitions can be naturally altered to describe the amount of inward protection offered to the wearer, but we will focus on exhalation in this paper.

The growth in popularity of community face coverings led to a great deal of innovation, for example, by Virustatic developing the Virustatic Shield${}^{\textrm{\textregistered}}$. Virustatic’s aim was to design a face covering, based on a neck gaiter, which was more permeable than a traditional cloth mask, and made of a face conforming elastic material, with the aim of reducing the leakage. This ensures a greater proportion of exhaled air is filtered through the face covering, which is treated with an active coating to maintain a high filtration efficiency [10]. However, as the regulations intended for surgical masks focus only on material properties, and do not account for leakage, there is no way for their design to progress up the regulatory pathway, despite aiming to offer better real-world protection than a traditional face mask.

There have been a number of experimental studies on the effect of face coverings worn by human volunteers, including [7, 11] who measured the effect of face masks on the distance jets of air travel from the wearer during talking, breathing, and coughing. The measurements presented in [11] indicate that the velocity of air just as it leaks out the top of the mask is $0.1\text{\,}\mathrm{m}\,\mathrm{s}^{-1}$ during normal respiration. In [12], direct measurements of the filtration efficiency and inward face-fitted filtration efficiency of coverings are presented for reusable fabric masks and disposable surgical and N95 masks. They find a variation in the filtration efficiency between the different materials, and that ensuring an effective seal to minimise leakage has a much greater effect on protecting the wearer than the material choice. Other papers [13, 14, 6, 15, 16] have focused on using headforms to improve reproducibility. Schlieren imaging is used in [13] to visualise the jets from headforms wearing face coverings to identify where leakage occurs and see the main leakage jets out of the top and sides of the mask. The filtration efficiency and inward and outward face-fitted filtration efficiencies of eleven different face coverings worn on a mannequin were measured in [14]. They found the difference between the material and face-fitted value was greatest for loose-fitting coverings. In [6], the outwards face-fitted efficiency during a simulated coughing cycle was measured for a variety of face masks, a neck gaiter, and a face shield worn by a mannequin. They found no significant statistical difference between a cloth mask, a surgical mask, and a neck gaiter, but all performed significantly better than the face shield. In [15], a headform was also used to measure the aerosol concentration for a surgical mask, an N-95 mask, and a homemade cloth mask, and found that all face coverings reduced the concentration of aerosols in front of the wearer, but increased it above and to the sides, indicating leakage of particles was occurring. They also found that the differences in concentration between different masks decreased with distance from the wearer. In [17], the leakage flow was estimated to be 25% for cloth and surgical masks, whereas [16] estimate a leakage between 43% and 95% at low flow rates, with a decrease of about 10% if the flow rate increases. These were calculated in both papers by (i) measuring the pressure drop across the face covering in the case where there is leakage, and (ii) when there is a perfect seal, then used an empirical model to estimate the leakage. These results show the wide range of values for leakage available and illustrate how it appears to be sensitive to a particular wearer’s facial geometry.

Several studies have measured the permeability and filtration efficiency of a wide range of materials used in the manufacture of cloth and surgical masks [18, 12, 19, 20]. Analysis of these results has demonstrated that no simple relationship exists between permeability and filtration efficiency due to the many factors affecting them, including porosity, weave pattern, thread size, and electrostatic properties [21]. These measurements have been used to parametrise a variety of studies modelling the airflow and aerosol transport of the respiration of a face covering wearer, which then allow the more relevant leakage ratio and face-fitted filtration efficiency to be computed. The models can broadly be divided into two categories: empirical models and computational fluid dynamics (CFD) models, which we discuss separately.

The simplest empirical model is an equivalent circuit model with a resistor representing the face mask in parallel with a resistor representing flow through a gap [22, 23]. The resistance through the gap is governed by Darcy’s law, while the resistance to flow through the gap is given by Bernoulli’s equation [23], or a resistor governed by Bernoulli’s equation in series with another given by the Hagen–Poiseuille equation [22]. Both models are used to predict the amount of leakage out of the sides of the mask and show that this may be a significant fraction of the total flux. The solution presented in [22] is validated in two dimensions against the numerical solution of the Navier–Stokes equations and found good agreement in the proportion of leakage predicted. The model in [22] has been extended in [24] to contain additional spatial information by dividing the boundary into segments, represented by parallel resistors, with the resistance in each segment varying according to simulations of the gap width on the perimeter of a face mask. Leakage was again found to be significant. The basic two-resistor model has also been coupled to an equivalent circuit model for the human respiratory system [25]. Here, it was found that the material of a surgical or N95 mask provided a sufficiently high resistance to breathing that it would become difficult to do so without the effect of leakage lowering the pressure required. Two-resistor equivalent circuit models have been used to quantify the amount of leakage during experiments on mannequin heads in [17] and [16]. The common feature of all these models is that they admit analytic or fast numerical solutions, allowing for rapid parameter sweeps. However, the equivalent circuit models assume a uniform pressure within the gap behind the mask, making it impossible to predict spatial variation in the flux through the mask.

The other existing modelling approach involves computational fluid dynamics (CFD) simulations, where the full Navier–Stokes equations are solved in realistic three-dimensional geometries representing the face covering, the wearer’s head, and the surrounding air. These models have been used to simulate the behaviour during normal breathing [26, 23, 27, 28], exercise [26], and coughing [21, 29, 28]. While [26, 23] consider laminar flow models, the remainder all include a turbulence model to account for the high Reynolds number of the flow. The different models are reviewed in [28], along with a discussion of numerical issues that arise when implementing these models. The porous face covering is either modelled as a separate region [26, 23] or as a spatially varying body force in the fluid momentum equations [29, 28, 30, 27, 21]. Resistance to flow through the mask is commonly assumed to be linear, that is, to obey Darcy’s law, for example, in [26, 23, 27, 31]. However, in [29, 28, 30], they assume a nonlinear resistance, motivated by experimental evidence that the pressure drop–flux relationship may become quadratic at high flow rates [11]. Given a solution for the airflow, different methods are adopted to estimate the resulting face-fitted filtration efficiency. The simplest model, adopted by [21, 23], assumes perfect adherence of the particles to the streamlines, i.e., an advective model, an assumption which is only true for aerosol-sized droplets (diameter less than $5\text{\,}\mathrm{\SIUnitSymbolMicro}\mathrm{m}$) [21]. To simulate larger particles, their inertia must be accounted for, which [26, 29] do using Lagrangian particle dynamics. In particular, [29] show that a large proportion of larger particles will impact the surface of a face mask even when the leakage ratio is high. A convection–diffusion model is considered in [30]. These simulations all provide a large amount of detail about features in the flow, such as the strength of secondary jets and the distance particles may be transmitted from the wearer. However, the computational resources required to compute the fluid flow limit the results reported to a small number of simulations, with [30] reporting that five high-fidelity simulations of a single coughing cycle required, on average, four months of CPU time on a 1000 CPU cluster, while [29] reported using three days on a 32 CPU cluster to simulate 10 coughing cycles. The parameter sweep of different mask permeabilities in [21] required two million CPU hours of computational resources from the National Science Foundation’s XSEDE [32].

In this paper, we build and solve a simple mathematical model for the exhalation airflow in the gap between a face covering and the wearer’s face, and through the face covering. We will use asymptotic analysis to systematically reduce the Navier–Stokes equations and Darcy’s law to obtain a reduced model. Our reduced model has the advantage of satisfying physical laws, giving interpretability, and allowing fine spatial resolution, but may be solved in a similar computational time as an empirical model. Using our reduced model, we explore how the amount of air leaking out of the mask depends on the permeability of the face covering, the geometry of the mouth, and the gap between the face covering and the wearer’s face. We will also compare two designs of face covering: the standard face mask, representing surgical and cotton masks, and a neck gaiter. The former is open on all sides and attached to the wearer’s face by elastic straps around the ears. The latter comprises a closed loop of fabric that wraps around the entire head, and so is only open to the air at the top and bottom. The neck gaiter is held to the wearer’s face by the elasticity of the material itself. We couple the flow model to a simple advection-capture model for aerosolised viral particles, which we will use to predict face-fitted filtration efficiency. By proposing a simple model linking the particle capture rate to the permeability of the mask material, we will explore the trade-offs in selecting a material to minimise leakage or maximise filtration efficiency.

The remainder of the paper is structured as follows. In Section II we present and nondimensionalise a mathematical model for the airflow and aerosol transport through a face covering and in the gap between the face covering and the wearer’s face. In Section III we will exploit the small size of (i) the aspect ratio of the gap between the face covering and the mouth, and (ii) the reduced Reynolds number, to obtain a reduced problem consisting of a single scalar partial differential equation (PDE) for the flux through the mask at each point. In Section IV we solve this reduced model explicitly for a constant gap width, and we simulate the model in the case where the gap has variable size. Using these solutions we will explore how the leakage ratio varies with the properties of the material, the size of the gap between the face and the mask, and between a face mask and a neck gaiter. In Section V we will solve the aerosol transport problem and calculate the face-fitted filtration efficiency of a face covering and explore how the parameters in the problem can be varied to optimise this. We present our conclusions in Section VI.

II Mathematical model

II.1 Airflow modelling

We begin by constructing a model for the airflow through the gap between the face covering and the wearer’s face, and through the face covering. For simplicity, we pose the model only for exhalation, but discuss its modification to inhalation in Section VI. We first consider the geometry shown in Fig. 1, where we assume the curvature of the face-covering-wearer’s head is negligible, so that the features of their face can be described by a function of the form $\hat{z}=\hat{f}\left(\hat{x},\hat{y}\right)$, where $\left(\hat{x},\hat{y},\hat{z}\right)$ form a three-dimensional Cartesian coordinate system. We assume that the rear of the face covering is positioned at $\hat{z}=\hat{g}\left(\hat{x},\hat{y}\right)$, so that the width of the gap between the face and the face covering is given by $\hat{h}=\hat{g}-\hat{f}$. We assume that the face covering is inextensible and has constant thickness $\hat{T}$, so that the front of the face covering is located at $\hat{z}=\hat{g}\left(\hat{x},\hat{y}\right)+\hat{T}/({\bf n}_{g}\cdot{\bf k})$, where ${\bf n}_{g}$ is the unit normal to the surface at $z=g$ pointing into the mask and $\bf k$ is a vector pointing in the $z$-direction. We further assume that the mask is rectangular in plane view, so that it occupies the region $[0,\hat{a}]\times[0,\hat{b}]$ in the $\hat{x}\hat{y}$-plane, as shown in Fig. 1.

We assume that the air is a Newtonian fluid. Thus, in the gap between the face and the face covering, fluid velocity, $\hat{\bm{u}}=\left(\hat{u},\hat{v},\hat{w}\right)$, and pressure, $\hat{p}$, satisfy the Navier–Stokes equations

[TABLE]

where $\hat{\rho}$ is the density, and $\hat{\mu}$ the dynamic viscosity of the air, both assumed constant.

We model the flow through the face covering using Darcy’s law, which has been shown experimentally to be valid for face coverings at low breathing rates [11]. For simplicity, we will assume that the material of the face covering is homogeneous and isotropic. Thus, the governing equations that describe the (local average) velocity of fluid within the pore-space, $\hat{\bm{U}}=\left(\hat{U},\hat{V},\hat{W}\right)$, and the pressure, $\hat{P}$, are

[TABLE]

where $\phi$ is the porosity of the mask, and $\hat{K}$ is the permeability of the mask, both assumed constant. We will assume that there is a constant unidirectional flux out of the mouth, in the normal direction to the surface defining the face, and no-flux and no-slip conditions hold on the remainder of the face. These conditions may be written compactly as

[TABLE]

where $\bm{n}_{f}$ is the unit normal to $\hat{f}$ pointing out of the wearer’s face, and

[TABLE]

is a piecewise constant function, where $\hat{w}^{*}$ is the exhalation velocity, and $\Omega_{m}\subset\left[0,\hat{a}\right]\times\left[0,\hat{b}\right]$ is the location of the wearer’s mouth on the surface $z=f(x,y)$. On the rear of the mask, we impose continuity of flux, pressure, and following the discussion in [33], a no-slip condition,

[TABLE]

We assume that the air surrounding the mask and wearer is at atmospheric pressure. Consequently, we set the pressure in front of the mask

[TABLE]

We also set the pressure to be atmospheric at any open boundaries of mask material. A traditional face mask will be open on all sides, so we set

[TABLE]

A neck gaiter is open only on the top and bottom, so we set

[TABLE]

To complete the specification of the fluid problem, we assume that the fluid velocity is initially zero, i.e. we write

[TABLE]

II.2 Aerosol transport modelling

We now consider a simple model for the transport of viral-containing aerosols by the airflow. We assume that the aerosols have negligible volume fraction and do not affect the airflow. We further assume that all the aerosols are identical, have constant physical properties, and are transported by advection alone. Thus, in the gap between the face covering and the wearer’s face, their distribution is described by the transport equation

[TABLE]

where $\hat{c}$ is the concentration of aerosols, measured in \unitmol.m^-3. In the mask, we will assume that all the pores are identical, and the number of aerosol particles is small compared with the number of pores, so by [34], the appropriate transport equation is

[TABLE]

where $\hat{C}$ is the (local average) concentration of aerosols and $\hat{\gamma}\left|\hat{\bm{U}}\right|$ is the rate of particle capture by the filter. We further assume for simplicity that the capture of aerosols does not result in a change in the porosity, permeability, or capture rate of the mask. Finally, we must impose boundary conditions for the aerosol. Since we are considering the airflow during exhalation, we will focus on the outwards protection offered by the face covering, and so assume that the viral aerosols are being generated by the wearer. Consequently, we impose the concentration at the mouth,

[TABLE]

and continuity in the flux of aerosols across the rear surface of the mask,

[TABLE]

Conditions on the remaining boundaries are not necessary due to the hyperbolic nature of (12)–(13). Finally, we assume that there are no viral aerosols in the gap or face covering initially, and so we write

[TABLE]

II.3 Nondimensionalisation

Our model comprises the six equations (1)–(4) and (12)–(13), subject to the boundary conditions (5)–(8) and (14)–(15), along with (9) for a face mask or (10) for a neck gaiter, and initial conditions (11) and (16) for the six variables $\hat{\bm{u}}$, $\hat{p}$, $\hat{\bm{U}}$, $\hat{P}$, $\hat{c}$, and $\hat{C}$. We nondimensionalise the model using the following scalings

[TABLE]

Here, $\epsilon=\hat{h}_{0}/\hat{b}$ is the aspect ratio of the gap between the face covering and the face, where $\hat{h}_{0}$ is a typical thickness (equal to the mean of $\hat{h}(x,y)$ over the face). The velocity scale has been chosen so that there is a unit total dimensionless flux through the mouth, and we have chosen the pressure scaling to balance the pressure gradient along the gap between the face covering and the wearer’s face with viscous forces in the $z$-direction. Applying the scalings (18) to the space between the mask and the wearer’s face, (1)–(2) and (12), gives

[TABLE]

where $\mathrm{Re}=\hat{\rho}\hat{A}\hat{w}^{*}/\left(\hat{\mu}\hat{h}_{0}\right)$ is the Reynolds number. The dimensionless governing equations within the mask, (3)–(4) and (13), are

[TABLE]

where

[TABLE]

are the dimensionless permeability and the dimensionless capture number, respectively. The boundary conditions, (5)–(8) and (14)–(15), become

[TABLE]

[TABLE]

[TABLE]

where

[TABLE]

where $M=\hat{A}/\hat{b}^{2}$ is the dimensionless area of the mouth and $T=\hat{T}/\hat{h}_{0}$ is the ratio of the mask width to the reference gap width. (As mentioned previously, our scalings result in $\iint_{\Omega_{m}}w_{\mathrm{in}}\,\rm{d}\Omega_{m}=1$, since $\Omega_{m}$ has dimensionless area $M$). The conditions specific to a face mask, (8), become

[TABLE]

where $l=\hat{a}/\hat{b}$ is the lateral aspect ratio and those for a neck gaiter, (10), become

[TABLE]

The model (24)–(35) depends on seven dimensionless parameters: the aspect ratio $\epsilon$, the reduced Reynolds number $\epsilon^{2}\mathrm{Re}$, the dimensionless permeability $K$, the dimensionless capture number $\Gamma$, the ratio of the mask area to the mouth area $M$, the dimensionless thickness of the mask $T$, and the mask’s lateral aspect ratio $l$. In addition, the definitions of $f(x,y)$, $g(x,y)$, and $\Omega_{m}$ will contain dimensionless geometric factors. In Table 1, we present the range of typical values of the dimensional parameters, representing a wearer at rest. With the exception of mean gap width, we will use fixed reference values to estimate the value of the dimensionless parameters, given in Table 2. Since variation in mean gap width can have a significant impact on the behaviour of our model, we will consider how the behaviour changes over the range of gap widths from an excellently fitting to a badly fitting face mask, based on mean gap sizes reported from simulations in [24]. The reference value for $\hat{h}_{0}$ corresponds to the mean gap width reported for a typical fit with a gap adjacent to the nose in [24]. Since no corresponding simulations exist for a neck gaiter, we will use these values for both designs of face covering. No data is currently available on the virus denature constant, $\hat{\gamma}$.

From Table 2, we see that $\epsilon$ and $\epsilon^{2}\mathrm{Re}$ are small and so we will work henceforth in the limit where $\epsilon\ll 1$ and $\epsilon^{2}\mathrm{Re}\ll 1$, the standard approximations that lead to lubrication flows [35]. For the other parameters (including the dimensionless virus denature constant), we assume that they are $\textit{O}(1)$, although we note that $M$ is small (but not as small as $O(\epsilon^{2})$). There are geometric parameters contained in the descriptions of the face and face covering. We will assume that variations in $f$ and $g$ are long wavelength, so that gradients are at most $O(1)$.

III Leading-order flow model

We now seek a steady solution of (19)–(35) in the limit of small $\epsilon$ and $\epsilon^{2}\mathrm{Re}$, and we neglect terms of these sizes henceforth. As we only consider the leading-order behaviour, we will not introduce extra notation to distinguish this. We will first solve for the flow through the mask and then for the flow in the gap between the face and face covering.

III.1 Flow through the mask

The leading-order versions of (24)–(27) for flow in the mask are

[TABLE]

for $g(x,y)<z<g(x,y)+T$, and are subject to the leading order of the boundary condition at the front of the mask, (32),

[TABLE]

We solve the first and last equation in (36), using (37), to find

[TABLE]

where $Q(x,y)$ is the unknown (Darcy) flux of air through the mask at each point. We note that the remaining equations in (36), provide expressions for $U$ and $V$ in terms of $Q$.

III.2 Flow in the gap between face and face covering

The leading-order versions of the equations for the flow in the gap between the face and the face covering (19)–(22) are

[TABLE]

for $f(x,y)<z<g(x,y)$, and, using the solution for the mask, the leading-order dimensionless boundary conditions (30)–(31) are

[TABLE]

The leading-order boundary conditions on the perimeter of the mask are

[TABLE]

for a face mask, using (34), and

[TABLE]

Equation (42) tells us that the pressure in the gap does not vary with $z$ and so, using (44), we find that

[TABLE]

Consequently, we integrate (40) and (41) and apply (43) and (44) to find

[TABLE]

We thus calculate the total flux parallel to the face in the gap between the face and the face covering, $\bm{q}=\left(q_{x},q_{y}\right)$, to be given by, in component form,

[TABLE]

where $h(x,y)=g(x,y)-f(x,y)$ is the (dimensionless) gap width. Integrating the conservation of mass equation (39) and applying the flux boundary conditions (43) and (44), yields the relationship

[TABLE]

Substituting (47) into (49) and combining with (50), we arrive at a single partial differential equation for the flux $Q$, namely

[TABLE]

where $k^{2}=12K/T$. Equation (51) is an inhomogeneous modified Helmholtz equation [37]. Using (45) and (47), we see that appropriate boundary conditions for $Q$ for a face mask are

[TABLE]

while (46) and (47) give the appropriate conditions for a neck gaiter,

[TABLE]

We therefore have a closed problem (51)–(53) for $Q$, and using its solution, may subsequently determine $W$, $P$, $p$, $w$, $u$, and $v$ by using, (38), (39), (47), and (49). Given a solution to (51)–(53), we calculate the total flux through each boundary in the problem. We therefore define

[TABLE]

as the magnitude of the total flux through the mask, out of the mouth, and through the left, right, bottom and top boundaries, respectively. By our choice of scaling we know that $\eta_{\mathrm{in}}=1$, implying that $\eta_{\textrm{mask}}$, $\eta_{\textrm{left}}$, $\eta_{\textrm{right}}$, $\eta_{\textrm{top}}$, and $\eta_{\textrm{bottom}}$, are just the proportion of the total flux out their boundaries, respectively, as shown in Fig. 2.

It follows that the leakage ratio, $\eta$, is given by

[TABLE]

Often, we will consider configurations with horizontal or vertical symmetry, and therefore define

[TABLE]

as the total proportion of the flux leaking out the horizontal and vertical boundaries, respectively (and note that $\eta=\eta_{\mathrm{hoz}}+\eta_{\mathrm{vert}}$). Integrating (51) over the mask implies $1=\eta_{\textrm{mask}}+\eta_{\textrm{left}}+\eta_{\textrm{right}}+\eta_{\textrm{bottom}}+\eta_{\textrm{top}}$, which expresses the total conservation of mass, and provides the alternate expression $\eta=1-\eta_{\textrm{mask}}$.

IV Results for leakage

IV.1 One-dimensional results

We briefly consider the solution to our reduced model for a face mask in one dimension to gain intuition about its behaviour. For simplicity, we consider only the case $l=2$, with the mouth, of width $2\mathcal{M}$, centred on $x=1$, allowing us to use symmetry to only consider half of the domain. We also assume the gap width $h$ is constant ($h=1$) so that the problem is linear, and (51)–(52) reduce to

[TABLE]

where

[TABLE]

subject to

[TABLE]

By requiring continuity of $Q$ and its derivative at $x=1-\mathcal{M}$, we solve (57)–(59) to find

[TABLE]

From this, we use conservation of mass to calculate a closed-form expression for the leakage ratio,

[TABLE]

In Fig. 3 we plot $\eta$ as a function of $k$ for different values of $\mathcal{M}$. We see that varying $\mathcal{M}$ only has a small effect on the amount of leakage, and, once $\mathcal{M}<1/3$ (which corresponds to a mask more than three times the width of the mouth), the leakage is less than 2% from the value as $\mathcal{M}\to 0$, given by $\eta\sim 1/\cosh(k)$, and shown in the figure by a dashed line.

Examining the variation with $k$, we see that $\eta\to 1$ as $k\to 0$ and $\eta\to 0$ as $k\to\infty$. From (61) it is clear that the transition between these two limits happens exponentially when $k\approx 1$. In terms of the dimensional parameters, $k^{2}=12\hat{b}^{2}\hat{K}/\left(\hat{h}_{0}^{3}\hat{T}\right)$, which is the ratio between the resistance to flow through the porous material $\hat{K}\hat{b}^{2}/\hat{T}\hat{\mu}$, and the resistance to fully developed viscous flow through a two-dimensional channel $\hat{h}_{0}^{3}/12\hat{\mu}$. Therefore, small $k$ corresponds to a low permeability and large gap width, while large $k$ corresponds to high permeability and low gap width. A key observation is that $k^{2}$ scales linearly with the permeability and the mask thickness, but with the inverse cube of the gap width. This means that the leakage ratio will be more sensitive to changes in fit than to changes in material. Although this one-dimensional solution provides an explicit approximation for the leakage, we cannot use it to investigate the differences between a mask and a neck gaiter, as we cannot simultaneously have two ends joined to form a closed loop, and still have the possibility of leakage. We therefore now investigate the full two-dimensional version of (51)–(53).

IV.2 Two-dimensional results: constant gap width

We again consider the solution under the restriction that $h$ is constant, i.e. $h\equiv 1$. We begin by seeking a solution to (51) for a traditional face covering, so that (52) holds on all four boundaries. It is then clear from the boundary conditions that the appropriate series expansion is

[TABLE]

with the coefficients given by

[TABLE]

We also represent the velocity through the mouth in the same basis, and so write

[TABLE]

where

[TABLE]

which depends on the given shape of the mouth. By multiplying (51) by the $mn$-th term in (62), integrating by parts and using orthogonality, it is straightforward to show that $\alpha_{mn}$ satisfies

[TABLE]

We note that $\beta_{mn}$ is independent of $k$, so only needs computing once per mouth. For the simplest case of a rectangular mouth with centre $\left(E_{x},E_{y}\right)$, and dimensions $2L_{x}$ and $2L_{y}$ in the $x$- and $y$-directions respectively, equation (65) yields

[TABLE]

Given (62), we calculate the depth-integrated flux of the air behind the mask via (47) and (49), to find

[TABLE]

which are straightforwardly shown to be convergent. By evaluating (68) along the boundary we calculate the flux of air leaking out of the sides at each point. We also substitute (62) and (68) into (54) to find expressions for the total fluxes through each boundary, to find that

[TABLE]

allowing us to evaluate the leakage ratio.

In Fig. 4, we show an example solution for the flux of exhaled air from a face mask with constant gap width and $k=1$, found using (62) and (68). For numerical purposes, the infinite sums are truncated at $m=200$ and $n=100$, to account for the size of the domain in different directions; beyond these, the solutions were indistinguishable as $m$ and $n$ increase. The central panel shows the flux through the mask at each point $Q(x,y)$, while the side panels show the flux out of the sides, given by the non-zero component of $\bm{q}$. The position of the mouth is shown by the red-dashed box. As the mouth is positioned at the centre of the mask, the flux out of the top and bottom, and the left and right of the mask are symmetric. We observe that the flux through the mask varies spatially, with a larger flux closer to the mouth. This shows that treating a mask as a single resistor with uniform spatial flux, as in an equivalent circuit model, neglects important behaviour that is required for accurately predicting the flux. The figure also shows that the amount of leakage at each point on the boundary decreases as the distance from the mouth increases. This is because the gap width is uniform, and so the viscous resistance to flow in the gap is simply proportional to the distance the air flows through the gap. We also see that the leakage flux from the corners of the mask is zero. This arises from our assumption that the air everywhere outside the region covered by the mask is at atmospheric pressure and hence $\bm{q}$ becomes zero at the sharp corner at leading order (since $q_{x}=0$ along the horizontal boundaries and $q_{y}=0$ along the vertical boundaries in our approximation). We explore the effect of rounding at these corners in Appendix A.

We note that, when $h=1$, so that $g(x,y)=f(x,y)+1$, it is possible to find $w$ explicitly, by substituting (47), (48), and (51) into (39), then integrating with respect to $z$, and applying the boundary condition (43). We find that

[TABLE]

IV.3 Neck gaiters

In the neck gaiter case, we solve (51) subject to the boundary conditions (53), so that periodic conditions hold at $x=0$, $x=l$, yielding similar expressions to (62)–(68) as given by (83)–(89) in Appendix B. In Fig. 5, we plot the flux of exhaled air from a neck gaiter, predicted by (83)–(89), when $k=1$ in the case where the mouth is at the centre of the mask. As the material must wrap around the wearer’s head, we assume the domain is larger and take $l=6$. From this solution, we see that the majority of flux out through the neck gaiter material and the boundaries is occurring in the central portion of the domain, and so not all the material is being used.

To compare the two mask designs, in Fig. 6 we show the leakage ratio $\eta$ as a function of $k$ for both a face mask with $l=2$ (as in Fig. 4) and for a neck gaiter with $l=6$ (as in Fig. 5), keeping all other geometric parameters identical. (We plot $\eta$ rather than $1-\eta$ for consistency with previous literature, e.g. [24].) We also show the horizontal leakage ratio $\eta_{\mathrm{hoz}}$ and the vertical leakage ratio $\eta_{\mathrm{vert}}$, defined in (56) for the face mask. This allows us to see that the majority of the leakage occurs through the top and bottom boundaries. This is to be expected as these boundaries are closer to a mouth positioned in the centre of the mask, and as the gap width is constant, resistance to flow behind the mask is equal at every point.

Comparing the total leakage ratio $\eta$ for a face mask and a neck gaiter as a function of $k$, we see that the two functions are almost indistinguishable. However, based on the parameter ranges for $K$ in Table 2, we estimate the value of $k$ for a face mask and a neck gaiter to be in $1.47.2$ and $3.116.0$, respectively, depending on the fit. These ranges are shown in Fig. 6 as grey and green regions, respectively, and the value of $k$ for a nominal fit are shown by grey and green vertical lines. We see that our model predicts that a face mask can leak as little as 5% to as much as 81% of the airflow, depending on the fit, with 21% leakage with a typical fit, whereas the neck gaiter will have at most 40% leakage, and only 1% for a typical fit. This suggests that a well-fitted face mask may have less leakage than a badly fitted neck gaiter, but on average a neck gaiter will leak less than a face mask. However, our estimate of the gap width is based on simulations of the fit of a face mask, whose material is much less elastic than that of a face mask, so the worst-case leakage of 40% arising from the gap width of an oversized mask is likely to be an overestimate. To see why $\eta$ is very close for a mask and neck gaiter when $k$ is equal, we can compare Fig. 4 and Fig. 5, which are both shown for $k=1$. Although the neck gaiter removes the possibility of air escaping from the sides, this comes at the expense of larger upper and lower boundaries than the face mask from which air can leak out. As a result, approximately the same total leakage can occur in both the face mask and the neck gaiter.

Until now, we have only considered the amount of leakage when the mouth is positioned in the centre of the mask. In Fig. 7, we show how the total leakage ratio varies as the centre of the mouth, $\left(E_{x},E_{y}\right)$, varies across the $xy$-plane. The results show that the lowest leakage occurs when the mouth is positioned at the centre of the mask. This is expected as it maximises the distance from the mouth to the sides of the mask. We also see that the increase in leakage ratio as the mouth moves towards the edge is nonlinear, with small perturbations from the centre having a small effect on the leakage, but a large increase in the leakage when the mouth approaches the edge of the mask. We could also consider varying the shape of the mouth, but we find this has minimal effect on the amount of leakage (see Appendix C for details).

IV.4 Two-dimensional results: variable gap width

In reality, the width of the gap between a face covering and the wearer’s face is not constant, so to predict realistic behaviour we must solve (51)–(52) for variable $h$. We do this by employing the finite element method (also used in Section IV.2 to solve the reduced model with constant $h$ for a mask with rounded corners). We will solve for a given functional form for $h$. We note that an analytic expression is not necessary and the problem could be solved for a given set of measurements [26] or using the output from of simulations [38, 39, 40, 41, 31]. We focus on using simple functions to obtain qualitatively correct results. We define $f(x,y)=F\left(x,y,\bm{\lambda}^{f}\right)$, and $g(x,y)=F\left(x,y,\bm{\lambda}^{g}\right)$ where $\bm{\lambda}^{f}$ and $\bm{\lambda}^{g}$ are different sets of parameters and the function $F\left(x,y,\bm{\lambda}\right)$ for $\bm{\lambda}=\bm{\lambda}^{f}$ or $\bm{\lambda}^{g}$ is defined by

[TABLE]

where $\bm{\lambda}=\left(\lambda_{i}\right)_{i=1}^{8}$. The function $F$ consists of a Gaussian distribution in the $x$-direction multiplied by a function representing the silhouette of a nose in side profile. The parameter $\lambda_{1}$ controls the horizontal position of the nose, and $\lambda_{2}$ and $\lambda_{3}$ control the horizontal spread, with $\lambda_{3}$ controlling the variation of the spread in the $y$ direction. The parameter $\lambda_{4}$ controls the maximum size of the nose, and $\lambda_{5}$ controls the rate of decrease from the tip to the bridge of the nose. The parameters $\lambda_{6}$ and $\lambda_{7}$ control the position of the bottom of the nose and the rate at which it grows to the tip, and $\lambda_{8}$ allows for an offset at the edge of the mask. The values of $\bm{\lambda}^{f}$ and $\bm{\lambda}^{g}$ used throughout this paper are given in Table 3 and chosen such that $h(x,y)=g(x,y)-f(x,y)>0$ for all $(x,y)\in\Omega$, and the mean of $h$ is one. Using the scalings (18), we find that the choice of parameters for $f$ gives a nose extending $3\text{\,}\mathrm{c}\mathrm{m}$ from the face, while the choice of parameters for $g$ gives a gap size either side of the nose peaking at $4\text{\,}\mathrm{m}\mathrm{m}$, which corresponds to the maximum gap size around the nose in the profiles given in [24]. The resulting dimensional profile is shown in Fig. 8.

In Fig. 9a, we show a heat map of $h(x,y)$ using the parameter values given in Table 3. We see that there is a small ‘u’-shaped region of gap width $h\approx{2.0}$ due to the nose (which joins the mouth to the top boundary) surrounding a region of very small gap width representing close fitting over the nose. We also show in Fig. 9a the lateral flux in the gap behind the mask, $\bm{q}$, (as black arrows), in the case when $k=6.0$, obtained by solving (51)–(52) and using (49), along with the position of the mouth as a red dashed line. We keep our attention focused on the flow from the mouth rather than adding in complications from the additional nasal flow. We see that the majority of the airflow is upwards in the gap along the ‘u’-shaped channel, with the speed decreasing away from the mouth, and as $h$ decreases. In Fig. 9b, we show the corresponding flux through the mask $Q$, and the leakage out of the sides. The position of the mouth is again shown as a red dashed box. We see that the flow through the ‘u’- shape channel creates a double-peaked distribution in the flow out of the top of the mask, and the presence of the nose creates an asymmetry in the flow through the mask, as expected. The flux through the mask is approximately uniform above the large-gap-width region located below the nose. This is because the gap is so large there is little resistance to the flow at that location, so it is easier for the air to spread out. We also see that the peak (i.e. the bright yellow region) of the flow rate through the mask is located beneath the mouth, rather than centred on it as seen in the constant-width solution shown in Fig.4. Physically, we can interpret this as air travelling along the side of the nose finding less resistance to leaking out of the sides compared with flowing through the mask, whereas the flow over the chin finds less resistance when flowing through the mask compared with leaking out.

In Fig. 10, we present the leakage ratio for the variable gap width case as $k$ varies. As our choice of mouth position and gap width is symmetric horizontally, we also show $\eta_{\textrm{top}}$, $\eta_{\textrm{bottom}}$, and $\eta_{\textrm{hoz}}$. The results show the same qualitative behaviour as the constant-$h$ case, with $\eta\to 1$ as $k\to 0$, and $\eta\to 0$ as $k\to\infty$, with a transition as $k$ goes between approximately $0.1$ and $10$. The different components of the leakage ratio show, as expected from Fig. 9, that all of the leakage occurs out of the top boundary for all $k$. In the figure, we also show the leakage ratio predicted by the constant-gap-width model, when $h\equiv 1$, the mean of the variable-gap-width model. We see that the constant-width model gives a reasonable approximation for the shape of the variable-width curve, but at lower values of $k$. This discrepancy gives a difference in predicted values up to 30%. These results suggest that, although the constant-width model cannot accurately predict where leakage will occur around the perimeter of a mask, it gives a reasonable approximation of the more complicated model provided the gap width is correctly parametrised. The sweep presented in Fig. 10 was computed for 250 data points and took only $50\text{\,}\mathrm{s}$ to run on a standard desktop, so each simulation takes only $0.2\text{\,}\mathrm{s}$ on average. Individual simulations take $2.5\text{\,}\mathrm{s}$ on average, suggesting mesh creating is much more costly than the problem solution. Therefore, even with variation in $h$, our reduced model is six orders of magnitude faster to run than any of the CFD simulations described in the Introduction.

V Aerosol transport

We now solve for the leading-order steady aerosol concentration. The leading-order transport equation in the gap (23) is

[TABLE]

which we solve subject to the boundary condition

[TABLE]

on $z=f(x,y)$. Equation (72) implies that $c$ is constant along streamlines and, since all streamlines originate from the mouth, it follows that $c=1$ everywhere in the gap and, in particular, on $z=g(x,y)$. The leading-order form of the governing equation of concentration inside the face covering (28) is

[TABLE]

subject to $C=1$ at $z=g(x,y)$ from (31). The solution is

[TABLE]

We straightforwardly calculate the filtration efficiency, to find that

[TABLE]

which has equivalent form to the formula given in [18]. Since we are considering an advection-only model, it is straightforward to see that the face-fitted filtration efficiency

[TABLE]

So far, we have assumed that the permeability and capture rate of a material are independent properties. However, we would intuitively expect a less permeable material to have a higher capture rate, as it will contain more microscale obstacles that impede the aerosol particles. To obtain an expression that captures this relationship, we assume that the permeability is related to the porosity, $\phi$, through the Kozeny–Carmen equation [42], a commonly used empirical model, in which

[TABLE]

where $\hat{D}_{p}$ is the effective average diameter of the microscale fibers (\unitm). Assuming the capture rate $\hat{\gamma}$ to be proportional to the available surface area on the microstructure, it is straightforward to deduce that

[TABLE]

where $\hat{\gamma}_{0}$ (\unitm^-1) is the surface reaction rate constant, which depends on the size and shape of the fibres and their packing and any coatings applied [10]. The $\phi$-dependent versions of $k$ and $\Gamma$ are

[TABLE]

$k_{0}=\hat{D}_{p}/\left(\sqrt{15}\epsilon^{2}\hat{b}\right)$ is a shape parameter that captures the effect of the microscale geometry and $\Gamma_{0}=\hat{\gamma}_{0}\hat{b}$ is the dimensionless surface reaction rate constant, capturing the effect of the surface chemistry on the macroscale capture rate. Combining these expressions allows us to link the capture number to the dimensionless permeability, which, for convenience, we write in terms of $k$ as

[TABLE]

Relationship (81) has the property that, as the permeability becomes infinite, the capture number tends to zero, whereas if the permeability tends to zero then the capture rate becomes infinite. The derivation of these expressions could be made more formal by specifying a particular microstructure model and geometry, then applying the method of multiple scales [43], but our expressions exhibit the correct qualitative behaviour, so we do not pursue this further here.

In Fig. 11, we use (81) to show how $\eta$, $\mathrm{FE}$ and $\mathrm{FFE}$ vary with $k$, with each panel showing a different value of $\Gamma_{0}$. For simplicity, we use the values of $\eta$ obtained by solving the reduced model for a face mask, (51)–(52), assuming constant gap width $h\equiv 1$, and a rectangular mouth at the centre of the mask. We see that $\eta$ is the same in each of the panels, since $\eta$ is independent of $\Gamma_{0}$. We do not show an equivalent curve for a neck gaiter, since, as already discussed, $\eta$ as a function of $k$ is indistinguishable for a face mask and a neck gaiter. The ranges of $k$ for a face mask and a neck gaiter, as the fit of the face covering changes, are shown as grey and green regions, respectively, with the value of a typical fit shown by vertical lines. All the panels show the existence of an optimal value, $k^{*}$ that gives the maximum $\mathrm{FFE}^{*}$. It is clear that, as $\Gamma_{0}$ increases, the values of $k^{*}$ and $\mathrm{FFE}^{*}$ both increase, so that the optimal point provides greater protection and is more permeable and/or tighter fitting as the surface reaction rate constant increases. We also see that, at lower values of $\Gamma_{0}$, $k^{*}$ lies in the range of possible values for both a face mask and a neck gaiter (with different mean gap widths). However, as $\Gamma_{0}$ increases, the optimum moves out of the range of $k$s possible for a face mask. Comparing the typical values of $k$ shown by the vertical lines in Fig. 11, we see that, for low $\Gamma_{0}$, the face mask is less permeable than the optimum and the neck gaiter is more permeable than the optimum, and thus both under-filter the aerosols compared with the optimum. However, as $\Gamma_{0}$ increases, the value of $k$ that gives the optimal efficiency moves away from the mean value associated with face masks and gets much closer to the mean value for neck gaiters. We illustrate how $k^{*}$ and $\mathrm{FFE}^{*}$, along with the corresponding leakage ratio $\eta^{*}$ and filtration efficiency $\mathrm{FE}^{*}$, vary with $\Gamma_{0}$ in in Fig. 12. We see that, to achieve a face-fitted filtration efficiency above 90%, the leakage ratio must be less than 2%.

The range of possible values of $k$ for a face mask and a neck gaiter are also shown in Fig. 12, and suggest that, by adjusting the fit, the optimal value $k^{*}$ can be achieved by a face mask for $\Gamma_{0}<1.2$, whereas a neck gaiter can achieve the optimal value for $\Gamma_{0}<7.2$. This means that, with these parameter values, our model predicts that a face mask can have a face-fitted filtration efficiency of at most 80%, whereas a neck gaiter can have a face-fitted filtration efficiency up to 99%. It is reasonable to also ask about the dependence of $\mathrm{FFE}$ on the only other remaining parameter in the problem, $k_{0}$. However, the results presented depend on $k_{0}$ only through (81) and hence variation in $k_{0}$ is equivalent to variation in $\Gamma_{0}$.

VI Conclusions, discussion, and extensions

In this paper, we have built and solved a model for the airflow and aerosol transport in and behind a face mask and a neck gaiter during exhalation. By working in the asymptotic limit of small aspect ratio of air gap to mask size, and small reduced Reynolds number, we were able to obtain a reduced model (51)–(53) for the airflow. This model was considerably easier to solve than a full CFD simulation, whilst retaining the effect of spatial variation that is not captured in previous reduced models. We used this model to understand how the various parameters in the problem affect two key real-world performance metrics: the leakage ratio $\eta$, and the face-fitted filtration efficiency $\mathrm{FFE}$, previously identified in the literature [9]. The behaviour of this reduced model was determined principally by a key dimensionless parameter, $k$, given by $k^{2}=12\hat{b}^{2}\hat{K}/\left(\hat{h}_{0}^{3}\hat{T}\right)$, where $\hat{K}$ is the permeability of the material of the face covering, and $\hat{h}_{0}$ is a reference value for the size of the gap between the face covering and the wearer’s face. Thus, $k^{2}$ is the ratio between the resistance to flow through the mask and the resistance to flow through the gap between the face and the face covering.

For the simplest case of a constant gap between face covering and the wearer’s face, it was possible to solve (51)–(53) explicitly, and find expressions for the leakage ratio for both a mask and a neck gaiter. Our results show that, as $k\to 0$ (i.e. $\hat{K}\to 0$, or $\hat{h}\to\infty$), $\eta\to 1$, whereas as $k\to\infty$ (i.e. $\hat{K}\to\infty$, or $\hat{h}\to 0$), $\eta\to 0$. In practice, this transition from all the flow leaking out of the sides to all the flow passing through the mask occurs as $k$ ranges from approximately $0.1$ to $10$. Besides $k$, the other parameters in the reduced model relate to the geometry of the mouth. We found that the mouth shape and position only had a small effect on the amount of leakage for any plausible choice of geometry. Therefore, the leakage ratio is determined almost entirely by the value of $k$. Using appropriate values for the permeability, we were able to estimate the possible range of $k$ for both a surgical face mask and a neck gaiter, which suggests that a neck gaiter will have less leakage for the same gap width. However, as neck gaiters are made of a more elastic material, they may conform better to the face and create a smaller gap. Both of these factors will increase $k$, and hence lead to a lower leakage ratio, highlighting the potential superiority of neck gaiters over face masks in minimising viral spread.

We also solved (51)–(53) numerically (using finite element methods) in the more general case of a variable gap size. Each solution is computed in less than $2.5\text{\,}\mathrm{s}$, representing a very significant saving in time and computational resources compared with a CFD simulation. By choosing an appropriate functional form for the gap size, we were able to reproduce realistic leakage profiles with a channel of larger gap size connecting the mouth to either side of the nose. We saw that the leakage rates depended on the dimensionless permeability in a qualitatively similar manner to that for a constant gap width, and was reasonably well approximated by the constant gap width model using the mean of the variable gap width. This represents a step improvement over previous reduced models in the literature, which could not capture the variation in the gap size away from the perimeter of the mask.

Our results show that increasing the permeability of the material will increase the proportion of air that flows through the face covering. However, increasing the permeability is also likely to decrease the filtration efficiency of the material, as there will be less material to impede the particles. To investigate this, we developed a simple model for the advection of aerosol particles by the airflow, with a capture rate that depends linearly on the local concentration of aerosol particles, the cross-sectional area of the material microstructure, and a surface reaction rate constant. By using the empirical Kozeny-Carman equation, we related the capture rate to the permeability of the material. Using this relationship, for a face covering of constant gap width we predicted the existence of an optimal material permeability which maximises the face-fitted filtration efficiency. By varying the surface rate constant of the material, we were able to show that the optimal value of face-fitted filtration efficiency always requires a face covering with only a small amount of leakage. These results also show that increasing the surface rate constant, for example by applying a coating [10], can substantially increase the maximum value of $\mathrm{FFE}$ and increase the corresponding value of $k$ at which the maximum occurs, to the point where they can be achieved by only a neck gaiter not a face mask.

Throughout this paper, we have considered the performance of the face covering during exhalation. However, by applying the transformation $Q\to-Q$ to (51)–(52), it is easy to see that leakage ratios predicted by our models for constant exhalation also describe the proportion of the flux through the mask during constant inhalation (since we are in a Stokes’ flow regime). Further, by assuming that there is a constant background concentration of aerosol particles outside the mask, it is again straightforward to show that the appropriately re-defined $\mathrm{FE}$ and $\mathrm{FFE}$ for inward protection [9] correspond exactly to the expressions we predict for exhalation. However, experimental results in [13] may be used to suggest that a small difference exists between $\eta$ for inhalation and exhalation; either the effect of the flow on the gap size, or inertial effects must distinguish one process from another.

In future, it would be interesting to consider a more detailed particle transport model accounting for additional physical processes, such as inertia, molecular diffusion, and electrostatic forces, as these can all be important for different-sized aerosol particles [9]. It would also be interesting to employ systematic upscaling techniques, such as in [43], to derive effective macroscale equations for the behaviour of the mask, where the macroscale parameters depend on microscale geometry and properties, such as weave pattern, thread size, and chemistry. The model could also be formulated in curvilinear coordinates to better capture the underpinning curvature of the head. However, we anticipate that, as the size of the gap is small compared with the curvature of the head, we would end up with the same model (51)–(53), since the same result is found for the standard lubrication flow equations [37]. In our analysis, we focused on exhalation from the mouth and neglected any nasal flow. However, none of our calculations require the region where air exits the face to be a connected region, and so the model is also valid for exhalation from disjoint nostrils, albeit with a different $\Omega_{m}$. In particular, for a face mask with constant gap width $h$, the expression for $\beta_{mn}$ (65) may be split into integrals over each discrete region, allowing, for example, the solution for two rectangular nostrils to be calculated using (67) for each nostril. An identical result also holds for neck gaiters. In our model, assume that the face covering material is inextensible. The model may be naturally extended to include an elastic face covering with permeability that may vary with tension. Spatial variations in permeability might affect the leakage calculation in interesting ways. We note that, inevitably, the constant-thickness case will be less relevant when we consider extensible materials. Another key extension will be to consider dynamic breathing (i.e. $w_{\mathrm{in}}$ is a function of time) and mouth opening and closing ($\Omega_{m}$ is a function of time). This has the potential to alter $h$ dynamically, especially if combined with the elastic model. Finally, this model was derived for velocities arising for light breathing, so it would be interesting to compare it with numerical simulations incorporating inertia in the airflow, which will become important during talking and coughing.

Our model provides a systematic framework for scientists and innovators to study permeable face masks and to determine how their operation depends on their physical parameters, aerosol capture rate, and fit. This will enable the expensive and time-consuming experiments to be focused on the relevant parameter regimes. It is a stepping stone towards re-evaluating how face masks protect the population from viruses released by infected wearers, which goes beyond current legislation, providing the opportunity for minimising transmission in future pandemics.

VII Acknowledgements

The main part of this work was funded by an Innovate UK Analysis for Innovators Award to Virustatic and the Isaac Newton Institute, grant reference 10074419.

J. Houghton, L. Hope, P. Hope, I. M. Griffiths, and C. J. W. Breward contributed to the conceptualization of the problem. J.Houghton and P. Hope provided background technical expertise, while L. Hope provided legislative background. The mathematical model was formulated, solved, and interpreted by M. D. Shirley, I. M. Griffiths, and C. J. W. Breward. The write-up was completed by M. D. Shirley, I. M. Griffiths, and C. J. W. Breward.

VIII Copyright

For the purpose of open access, the authors have applied a CC BY public copyright licence to any author-accepted manuscript (AAM) arising from this submission.

Appendix A Rounding the corners of a mask

In Fig. 4, we saw that the flux tends to zero at the corners of the mask. To investigate whether this effect is also present for masks with rounded corners, we solve the reduced model for $Q$, (51)–(52), numerically using the finite element package Firedrake [44], imposing a mask shape with rounded corners of radius $r$. In Fig. 13, we vary the amount of rounding of the mask boundary, and see that the fluxes are unaffected away from the corners, and only change appreciably in a small region close to the corners of the mask, where the flux becomes small but non-zero when the rounding is applied. This demonstrates that the flux approaching zero in the corners, while on the face of it unphysical, is only a local effect, with minimal impact on the overall behaviour of the system.

Appendix B Solution for a neck gaiter

We now derive the solution to (51) subject to the Dirichlet condition (52) at $y=0$ and $y=1$, and the periodic condition (53) at $x=0$ and $x=l$. We will assume the solution is sufficiently smooth that the periodic condition becomes

[TABLE]

where throughout we use $\widetilde{\cdot}$ notation to distinguish expressions for a neck gaiter from those for a face mask. Subject to these boundary conditions, the series expansion for $\widetilde{Q}$ is

[TABLE]

where the coefficients $\widetilde{\alpha}_{mn}$ are given by

[TABLE]

for $m=0,1,\dotsc$ and $n=1,2,\dotsc$, where, as in the case of a face mask, $\widetilde{\beta}_{mn}$ are the coefficients of the series representation of $w_{\mathrm{in}}$ in this basis, i.e.

[TABLE]

and are given by

[TABLE]

for a general mouth. For a rectangular mouth with centre $\left(E_{x},E_{y}\right)$ and dimensions $2L_{x}$ and $2L_{y}$, in the $x$- and $y$-directions respectively, then

[TABLE]

for $n=1,2,\dotsc$ and

[TABLE]

if $m,n=1,2,\dotsc$. The depth-integrated flux of the air behind the neck gaiter is found by (49), giving

[TABLE]

Finally, we consider the total fluxes, defined by (54). As we assume the solution is periodic in the $x$-direction, $\widetilde{\eta}_{\textrm{left}}=0$, $\widetilde{\eta}_{\textrm{right}}=0$. For the remaining fluxes, we find that

[TABLE]

Note that these global quantities only depend on the zeroth mode in the $x$-direction (since all other modes integrate to zero). Using these expressions, the leakage ratio may be calculated using (55).

Appendix C Effect of mouth shape

In this appendix, we consider the effect of varying the shape of the mouth, $\Omega_{m}$. For simplicity, we only present results for a face mask with constant gap width, but the findings generalise to a neck gaiter, and to a face covering with variable gap width.

If $\Omega_{m}$ is elliptical with centre $\left(E_{x},E_{y}\right)$ and width and height $\left(2R_{x},2R_{y}\right)$, respectively, we find that the expression for the coefficients of the series representation of $w_{\mathrm{in}}$, (65), are

[TABLE]

The integral must be evaluated numerically, which we perform using Simpson’s rule. Using this expression, we plot in Fig. 14a, how the total, horizontal, and vertical leakage ratios for a face mask vary with $k$ for both a rectangular and an elliptical mouth, positioned at the centre of the mask, and with equal areas. We also keep equal the aspect ratios of the mouths, defined as their width divided by their height, i.e. $L_{x}/L_{y}$ and $R_{x}/R_{y}$ for a rectangle and an ellipse, respectively. We conclude that, at this aspect ratio, the leakage ratios are indistinguishable.

In Fig. 14b we again show the leakage curves for rectangular and elliptical mouths, but varying their aspect ratio and keeping $k$ fixed. This plot shows that the leakage for the two mouth shapes are indistinguishable for moderate aspect ratios that encompass all possible realistic mouth shapes, and the leakage ratio remains approximately constant. For extremely small or large aspect ratios, the total leakage increases as the edge of the mouth becomes close to a boundary, greatly increasing the vertical or horizontal leakage, as appropriate. Since $\eta$ is increasing as the aspect ratio tends either to zero or infinity, there is also an aspect ratio that minimises the leakage. We find this to be 3.31 for a rectangular mouth and 4.07 for an elliptical mouth. Interestingly, this is slightly higher than 2, the equivalent aspect ratio of the standard face mask considered in this paper.