Hostname: page-component-cd9895bd7-gvvz8 Total loading time: 0 Render date: 2024-12-29T12:11:25.862Z Has data issue: false hasContentIssue false

The effect of pore-scale contaminant distribution on the reactive decontamination of porous media

Published online by Cambridge University Press:  08 August 2023

Ellen K. Luckins*
Affiliation:
Mathematical Institute, University of Oxford, Oxford, UK
Christopher J. W. Breward
Affiliation:
Mathematical Institute, University of Oxford, Oxford, UK
Ian M. Griffiths
Affiliation:
Mathematical Institute, University of Oxford, Oxford, UK
Colin P. Please
Affiliation:
Mathematical Institute, University of Oxford, Oxford, UK
*
Corresponding author: Ellen Luckins; Email: [email protected]
Rights & Permissions [Opens in a new window]

Abstract

A porous material that has been contaminated with a hazardous chemical agent is typically decontaminated by applying a cleanser solution to the surface and allowing the cleanser to react into the porous material, neutralising the agent. The agent and cleanser are often immiscible fluids and so, if the porous material is initially saturated with agent, a reaction front develops with the decontamination reaction occurring at this interface between the fluids. We investigate the effect of different initial agent configurations within the pore space on the decontamination process. Specifically, we compare the decontamination of a material initially saturated by the agent with the situation when, initially, the agent only coats the walls of the pores (referred to as the ‘agent-on-walls’ case). In previous work (Luckins et al., European Journal of Applied Mathematics, 31(5):782–805, 2020), we derived homogenised models for both of these decontamination scenarios, and in this paper we explore the solutions of these two models. We find that, for an identical initial volume of agent, the decontamination time is generally much faster for the agent-on-walls case compared with the initially saturated case, since the surface area on which the reaction can occur is greater. However for sufficiently deep spills of contaminant, or sufficiently slow reaction rates, decontamination in the agent-on-walls scenario can be slower. We also show that, in the limit of a dilute cleanser with a deep initial agent spill, the agent-on-walls model exhibits behaviour akin to a Stefan problem of the same form as that arising in the initially saturated model. The decontamination time is shown to decrease with both the applied cleanser concentration and the rate of the chemical reaction. However, increasing the cleanser concentration is also shown to result in lower decontamination efficiency, with an increase in the amount of cleanser chemical that is wasted.

Type
Papers
Creative Commons
Creative Common License - CCCreative Common License - BY
This is an Open Access article, distributed under the terms of the Creative Commons Attribution licence (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted re-use, distribution and reproduction, provided the original article is properly cited.
Copyright
© The Author(s), 2023. Published by Cambridge University Press

1. Introduction

Following the spill or malicious release of a hazardous chemical agent, it is crucial for public health that any materials contaminated with the agent are thoroughly cleaned. When the agent has seeped into a porous material, such as a brick or concrete block, the decontamination process typically consists of applying a cleanser chemical to the surface of the porous material and allowing this to react into the agent-contaminated material, thus neutralising the agent in a chemical reaction. The two fluids, the cleanser and the agent, are often immiscible. For example, the agent is often an oily substance, while the cleanser is often in aqueous solution. Thus, the decontamination reaction often occurs at the interfaces between these two fluids, within the pores of the porous material.

If the porous material is initially saturated with agent, a reaction front develops between the two fluids. We will refer to such a decontamination scenario as a ‘sharp-interface’ case. This scenario has been studied in [Reference Dalwadi, Dubrovina and Eisenträger2] and [Reference Dalwadi, O’Kiely, Thomson, Khaleque and Hall4], under the assumption that the decontamination occurred at a flat macroscale reaction front that moved into the porous medium as the decontamination process evolved. In the limit of a deep agent spill, the decontamination rate is found [Reference Dalwadi, O’Kiely, Thomson, Khaleque and Hall4] to be strongly affected by the partition coefficient of the reaction product between the two phases, as well as by the concentration of the applied cleanser and the rate of the chemical reaction. Two parameter regimes were identified, corresponding to whether (i) the decontamination was limited by the supply of cleanser or (ii) by the removal of the reaction product.

The effect of the pore-scale processes on the motion of such a reaction front was explored in [Reference Luckins, Breward, Griffiths and Wilmott10] (and a minor correction to the homogenised model is given in [Reference Luckins, Breward, Griffiths and Please11]), in the case that the reaction product is only soluble in the cleanser phase. A homogenisation analysis was used to derive an effective model, which has the same form as the model used in [Reference Dalwadi, O’Kiely, Thomson, Khaleque and Hall4], but with a modified reaction rate taking into account the effect of the microstructure on the motion of the reaction front. This sharp-interface model is equivalent to a Stefan problem with kinetic undercooling, as arises in models for many other physical systems. A mathematical analysis of such Stefan problems is given in [Reference Evans and King6].

The sharp-interface decontamination models are valid in scenarios when agent initially saturates the pores of the porous material, but it is not clear whether or not this is the case in reality. It is plausible that, initially, the agent only partially saturates the pore space. A second decontamination scenario was proposed in [Reference Luckins, Breward, Griffiths and Wilmott10], in which the agent initially coats the walls of the pores and does not saturate the pores. This scenario is an example of an unsaturated porous material, and we refer to it as the ‘agent-on-walls’ model. A homogenised model for this unsaturated microscale geometry was also derived in [Reference Luckins, Breward, Griffiths and Wilmott10] and takes the form of a reaction–diffusion equation for the concentration of cleanser, coupled with a reaction equation describing the evolution of the agent-layer thickness. This agent-on-walls model is similar to homogenised models describing the growth of biofilms in porous media [Reference Dalwadi, Wang, King and Minton5, Reference Schulz and Knabner12, Reference Schulz and Knabner13, Reference Wood, Golfier and Quintard15], reactive filtration [Reference Kiradjiev, Breward, Griffiths and Schwendeman9] and the accumulation of particles within a filter [Reference Dalwadi, Griffiths and Bruna3]. We note that, in both the sharp-interface and agent-on-walls models derived in [Reference Luckins, Breward, Griffiths and Wilmott10], the neat agent is stationary, and that the transport of cleanser chemical through solution is purely diffusive.

In this paper, we will examine solutions of the two homogenised models derived in [Reference Luckins, Breward, Griffiths and Wilmott10] describing the decontamination of porous media under two different initial pore-scale distributions of agent, namely the sharp-interface and agent-on-walls scenarios. In order to directly compare between the two scenarios, we suppose the total initial volume of agent within the porous material is the same. In particular, we examine how the overall decontamination time and the amount of cleanser required differ between the two scenarios. In Section 2, we describe the two decontamination scenarios and the common assumptions made about geometry and the decontamination procedure. In Section 3, we state and solve the sharp-interface model. Although the sharp-interface model has been analysed in detail before [Reference Dalwadi, O’Kiely, Thomson, Khaleque and Hall4, Reference Evans and King6], it is helpful to state some key results in different asymptotic limits, in order to compare directly between the two decontamination scenarios. In Section 4, we then state and analyse the agent-on-walls model, drawing analogy to the sharp-interface case throughout. Discussion and concluding remarks are given in Section 5.

2. Common geometry and model setup

We consider the reactive decontamination of a porous material, for two canonical initial distributions of agent within the pores, one in which the agent initially saturates a region of pore space (the sharp-interface case), and the other in which the agent coats the walls of the pores in a layer (the agent-on-walls case). As in [Reference Dalwadi, O’Kiely, Thomson, Khaleque and Hall4], the agent fluid is assumed to be stationary and does not flow or move through the porous material during the decontamination, in either the sharp-interface or agent-on-walls cases.

To directly compare between these two cases, we will assume the same pore-scale geometry of the porous material, consisting of circles, radius $r_0$ , in a square lattice, with centres a distance $d$ apart, in keeping with the pore-scale geometry used in [Reference Luckins, Breward, Griffiths and Wilmott10] for the agent-on-walls model. We consider a one-dimensional porous material, with spatial variable $y$ pointing into the porous material, with the external surface at $y=0$ . We suppose that agent initially resides uniformly within the region $y\in [0,L]$ in the porous medium. In order to directly compare the two models, we will also assume that there is the same initial total volume of agent (per unit cross-sectional area), $v\,$ [m $^3\,$ m $^{-2}$ ], within the porous medium in each of the two cases, and so the spill depth, $L$ , will be different for the two cases (as discussed below). A schematic illustrating the geometry is given in Figure 1.

Figure 1. Schematic illustrating the two canonical agent distributions within the pore space: the sharp-interface case (left) with the agent initially saturating the pore space and the agent-on-walls case (right).

In both cases, as in [Reference Luckins, Breward, Griffiths and Wilmott10] we assume the agent is neat, with (constant) molar volume $\chi \,$ [m $^{3}\,$ mol $^{-1}$ ], while the cleanser is in solution, with local concentration $c(y,t)$ , where $t$ is the time variable. For the agent-on-walls case, the solvent (likely water) of the cleanser solution is assumed to initially fill the remaining pore space around the agent layers. The agent and cleanser fluids are assumed to have the same density so that there is no flow generated by the reaction process, and the transport of cleanser chemical is therefore purely by diffusion within the solution. For simplicity, one mole of cleanser is assumed to react with one mole of agent, to produce the reaction product, although this may be easily generalised. This product of the chemical reaction is taken to be only soluble in the cleanser phase, so that the agent remains neat at all times. The reaction rate is modelled using the law of mass action, as in [Reference Dalwadi, O’Kiely, Thomson, Khaleque and Hall4, Reference Luckins, Breward, Griffiths and Wilmott10], so that the rate of reaction depends on the local concentration of cleanser at the reacting interface(s).

The same cleaning protocol is assumed to be applied in both scenarios: a cleanser is applied to the surface at $y=0$ , at a known concentration $c_0$ , and the cleanser is continuously replenished so that this surface concentration remains constant for all time.

3. Sharp-interface case

In this case, the agent initially saturates the pore space from the top surface of the porous material at $y=0$ to a fixed depth $L_{\text{SI}}$ . This depth of the spill is determined by the amount of agent, $v$ , that has been spilled, and the available pore space, and is thus given by

(3.1) \begin{align} L_{\text{SI}}=\frac{v}{\phi }, \end{align}

where $\phi =1-\pi r_0^2/d^2$ is the porosity of the medium. As in [Reference Dalwadi, O’Kiely, Thomson, Khaleque and Hall4, Reference Luckins, Breward, Griffiths and Wilmott10], we assume that the decontamination reaction occurs at a flat macroscale interface at $y=h(t)$ , which we call the ‘sharp interface’, which moves down from the top surface into the porous medium as the decontamination occurs, so that $h\in [0,L_{\text{SI}}]$ is an increasing function of $t$ . The cleanser chemical diffuses through the region $y\in [0,h)$ . We are only interested in the behaviour of the system until the time at which $h(t)=L_{\text{SI}}$ , when all the agent has been removed. Beyond that time, our model no longer holds.

3.1. Model statement

As derived in [Reference Luckins, Breward, Griffiths and Wilmott10] (and including the minor correction given in [Reference Luckins, Breward, Griffiths and Please11]), the dimensional model for the cleanser concentration $c$ in the sharp-interface case is, for $0\lt y\lt h\lt L_{\text{SI}}$ ,

(3.2a) \begin{align} c_t=D_0c_{yy}, \end{align}

where subscripts $t$ and $y$ denote partial derivatives with respect to time and space, respectively, along with boundary conditions at the reacting interface at $y=h(t)$ :

(3.2b) \begin{align} ch_t+D_0 c_y = -\beta kc, \end{align}
(3.2c) \begin{align} h_t =\chi \beta kc, \end{align}

while at $y=0$ ,

(3.2d) \begin{align} c=c_0. \end{align}

Here, the diffusivity $D_0$ is the effective diffusivity of the cleanser through the (constant porosity) porous material, given by solving a ‘cell problem’ (see [Reference Luckins, Breward, Griffiths and Wilmott10]), $k$ [ms $^{-1}$ ] is the surface reaction rate, $\chi$ [m $^3\,$ mol $^{-1}$ ] is the molar volume of the neat agent and $\beta$ is a dimensionless geometric parameter which is the ratio of the average length of the microscale reacting interface to the porosity of the porous medium that arises from averaging the microscale problem (see [Reference Luckins, Breward, Griffiths and Please11]). The dimensional parameters in our models are summarised in Table 1. We assume that, initially, the agent fills the porous medium, so that $h=0$ at $t=0$ . An initial condition for $c$ is not required since there is no cleanser-occupied region at $t=0$ .

Table 1. Dimensional parameters in the decontamination models (SI is sharp-interface model, AW is agent-on-walls model)

3.2. Nondimensionalisation

We nondimensionalise the model using the scalings

(3.3) \begin{align} c=c_0c^{\prime}, \qquad y=L_{\text{SI}}y^{\prime}, \qquad h=L_{\text{SI}}h^{\prime}, \qquad t=\frac{L_{\text{SI}}}{\chi c_0 \beta k}t^{\prime}, \end{align}

where we have fixed the cleanser concentration scaling, $c_0$ , by the boundary condition at $y=0$ , the length scaling to be the initial depth of the agent spill, and we have chosen the timescale associated with the motion of the interface due to the chemical reaction by balancing (3.2c). This is the timescale of a reaction-limited decontamination, which we will find is the shortest timescale the overall decontamination can take.

Substituting these scalings (3.3) into the model (3.2) and dropping the prime notation, we obtain the dimensionless problem, for $0\lt y\lt h(t)\leq 1$

(3.4a) \begin{align} \alpha \gamma c_t=c_{yy}, \end{align}

with boundary conditions at $y=h(t)$ ,

(3.4b) \begin{align} \alpha \gamma ch_t+c_y=-\gamma c, \end{align}
(3.4c) \begin{align} h_t=c, \end{align}

and on $y=0$ ,

(3.4d) \begin{align} c=1. \end{align}

We have introduced two dimensionless parameters, $\alpha$ and $\gamma$ , defined in Table 2. Since $\chi$ is the molar volume of the neat agent, $\alpha$ characterises the ratio of cleanser concentration in solution to the agent concentration (in its neat form). Since the cleanser is dilute while the agent is neat, it is likely that $\alpha$ is less than one. This justifies our choice of the reaction timescale over the diffusion timescale, since for relatively dilute cleanser ( $\alpha \lt 1$ ) the reaction will often be the controlling process. We note that $\alpha ^{-1}$ may be interpreted as the Stefan number, as defined in [Reference Evans and King6]. Meanwhile, $\gamma$ is a ratio of the lengthscale of the spill to be decontaminated, $L_{\text{SI}}$ , to the lengthscale of diffusion over the reaction timescale. We may therefore consider $\gamma$ to be a proxy for the depth of the spill: $\gamma$ is small for shallow spills and large for deep spills.

Table 2. Dimensionless parameters in the sharp-interface model (3.4)

This problem (3.4) is a simplified version of the decontamination model studied in [Reference Dalwadi, O’Kiely, Thomson, Khaleque and Hall4], in the case that the reaction product is only soluble in the cleanser fluid (or $\mathcal{K}\rightarrow 0$ in the notation of [Reference Dalwadi, O’Kiely, Thomson, Khaleque and Hall4]). As noted in [Reference Dalwadi, O’Kiely, Thomson, Khaleque and Hall4], the problem (3.4) is a one-phase Stefan problem with kinetic undercooling and has been studied in detail, for instance in [Reference Evans and King6]. In the remainder of this section, we solve (3.4) numerically and discuss some relevant limits, which will be helpful when we compare with the agent-on-walls model in Section 4.

3.3. Numerical method and early-time behaviour

To solve (3.4) numerically, we first transform to a fixed domain by making the change of variables $y=hY$ , for $Y\in [0,1]$ , and we then use the method of lines on the transformed problem with 100 spatial meshpoints, which we find to be sufficient for convergence of the solution. For the spatial discretisation, we use a central-difference scheme for the diffusive term and a first-order upwind scheme for the artificial advection due to change of variables to a fixed domain, so that the overall scheme is first order. We then use the inbuilt ODE solver ode15s in Matlab for the timestepping, noting that the system is stiff for $\alpha \gamma \ll 1$ . This is a multistep solver, using numerical differentiation formulas of order 1–5, specifically designed for stiff systems [Reference Shampine and Reichelt14].

The numerical method is not valid at $t=0$ since initially we have $h=0$ . In order to initialise the numerical simulations, we therefore use asymptotic expansions for early time. Specifically, we consider $t\ll 1$ to be small. As in, for instance, [Reference Dalwadi, O’Kiely, Thomson, Khaleque and Hall4, Reference Evans and King6], we see that in order to balance (3.4c) for small $t$ we must have $h$ and $y$ of order $t$ . In (3.4a) we see that the transport of cleanser is fast at early times, with $c_{yy}=0$ . Further, the boundary condition (3.4b) requires that $c_y=0$ at $y=h$ , and thus $c$ is uniform at early times, and equal to 1 using the boundary condition at $y=0$ . At next order in $t$ ( $\ll 1$ ), we find linear equations for the first-order correction terms, which we solve to find

(3.5) \begin{align} c= 1-\gamma (1+\alpha )y+O(t^2), \qquad h=t-\frac{\gamma (1+\alpha )}{2}t^2+O(t^3). \end{align}

These early-time expansions (3.5) are valid so long as $t\ll \max (1,\gamma ^{-1},(\alpha \gamma )^{-1})$ . We use these early-time expansions with appropriately chosen $t$ to initialise numerical simulations of (3.4).

In the remainder of this section, we solve the model both numerically and asymptotically in relevant limits of small and order-one $\alpha$ . For ease of later comparison with the agent-on-walls case, we summarise these asymptotic limits in Figure 2 and refer to them as case i (defined by $\alpha \ll 1$ and with the sublimits cases ia and ib for which $\gamma \ll 1$ and $\gamma \gg 1$ , respectively), case ii ( $\alpha =O(1)$ and $\gamma \ll 1$ ) and case iii ( $\alpha =O(1)$ and $\gamma \gg 1$ ), respectively. We will compare the analytic solutions in these asymptotic limits with full numerical solutions throughout.

Figure 2. Schematic of the asymptotic limits considered for the sharp-interface model. Case i is $\alpha \ll 1$ , with sublimits ia for $\gamma \ll 1$ and ib for $\gamma \gg 1$ . Case ii is $\alpha =O(1)$ and $\gamma \ll 1$ , and case iii is $\alpha =O(1)$ and $\gamma \gg 1$ .

3.4. Case i: dilute cleanser, $\alpha \ll 1$

As discussed in Section 3.2, we expect that a physically relevant limit is $\alpha \ll 1$ , in which the cleanser chemical is dilute in relation to the neat agent. In this section, we therefore study asymptotic solutions of (3.4) in the limit $\alpha \ll 1$ , which we refer to as case i in the schematic of Figure 2. We assume first that $\gamma =O(1)$ , before summarising the behaviour in the sublimits $\gamma \ll 1$ (case ia) and $\gamma \gg 1$ (case ib).

We assume first that $\gamma =O(1)$ . To leading order in $\alpha \ll 1$ , the system (3.4) is simply

(3.6a) \begin{align} c_{yy}=0, \end{align}

with boundary conditions at $y=h$ ,

(3.6b) \begin{align} c_y=-\gamma c, \qquad h_t=c, \end{align}

along with $c=1$ on $y=0$ , at $h=0$ at $t=0$ . This reduced problem has solution

(3.7) \begin{align} c(y,t)=1-\frac{\gamma y}{\sqrt{1+2\gamma t}}, \qquad h=\frac{1}{\gamma }\left ({-}1+\sqrt{1+2\gamma t}\right ), \end{align}

for $y\in [0,h]$ and $h\in [0,1]$ . In particular, the leading-order time at which the decontamination is completed, defined to be the time at which $h=1$ and the reacting interface has moved through the entire spill, is given by

(3.8) \begin{align} t_D=1+\frac{\gamma }{2}. \end{align}

We note that the diffusion of cleanser is quasi-steady in this limit.

This approximate solution (3.7) is valid for the sublimit $\gamma \ll 1$ , case ia, corresponding to a thin spill depth, and we note that in this case $c\approx 1$ at all times, as the cleanser consumed in the reaction is quickly replenished from the surface of the porous medium by diffusion. Additionally, $h\approx t$ in this sublimit, that is, the front moves linearly with time through the material, rather than like $\sqrt{t}$ . The decontamination is limited by the rate of reaction in this case.

We now consider large spill depths, $\gamma \gg 1$ , referred to as case ib in Figure 2. It is not immediately clear that the case-i solution is valid in this limit, since for $\gamma \gg 1/\alpha$ , we might expect the timescale of diffusion to re-enter the problem. We will, however, show that the case-i solution (3.7) is indeed valid for all $\gamma \gg 1$ .

For large $\gamma$ , the case-i solution is valid for early $O(1)$ time. However, the timescale of decontamination must be longer. Making the change of variables $t=\gamma T$ , we see that the model (3.4) becomes

(3.9a) \begin{align} \alpha c_T=c_{yy}, \end{align}

with, on $y=h$ ,

(3.9b) \begin{align} h_T=-\frac{c_y}{1+\alpha c}, \qquad \frac{1}{\gamma }h_T=c, \end{align}

(where we have substituted (3.4c) into (3.4b)), and as usual $c=1$ on $y=0$ . For $\alpha \ll 1$ and $\gamma \gg 1$ (for the leading-order behaviour it is not relevant whether $\alpha \gamma$ is large or small), we see that, at leading order in $\gamma ^{-1}\ll 1$ and $\alpha \ll 1$ , $c$ is linear and equal to zero at $y=h$ . The leading-order interface speed $h_T=-c_y|_{y=h}$ is determined by the flux of cleanser into the interface. The leading-order solution is therefore

(3.10) \begin{align} c=1-\frac{y}{\sqrt{2T}}, \qquad h=\sqrt{2T}, \end{align}

which we note is the sublimit $\gamma \gg 1$ of (3.7). The fact that $\alpha \ll 1$ means that the cleanser concentration profile evolves instantaneously with the interface motion. Indeed both diffusion and reaction are quasi-steady processes in the $\gamma \gg 1$ sublimit; the decontamination is limited by the flux of cleanser into the reacting interface.

The position of the reaction front $h$ against time and the concentration $c$ of cleanser against $y$ for various times through the decontamination process, in case i, given by (3.7), are plotted in Figure 3, along with numerical solutions of (3.4) for various $\alpha$ and $\gamma$ . We see that the cleanser diffuses into the porous medium, with a monotonic decreasing concentration front the surface $y=0$ to the reaction front $y=h(t)$ . The reaction front $h(t)$ progresses into the porous material, with a speed that is non-increasing in time. We observe excellent agreement between the numerical and analytic solutions for the relatively small value of $\alpha =0.1$ . For larger $\alpha =0.5$ , corresponding to more-concentrated cleanser solutions, we see that the dimensionless decontamination time increases, since the diffusive transport of cleanser to the reaction front is no longer quasi-steady. (We note that our chosen timescale, given in (3.3), is inversely proportional to $\alpha$ ; the dimensional decontamination time is considered in Section 3.6 below). For non-zero $\alpha \gamma$ , we also have a nonlinear interface condition (3.4b), which slows the interface motion. We observe in Figure 3b, d and f that $c$ is nearly linear in $y$ , although for larger $\alpha$ and $\gamma$ the profile becomes slightly convex. We note that for larger $\gamma$ there is $O(1)$ variation in $c$ over the domain (for instance in Figure 3b), while for smaller $\gamma$ (smaller spill depths) $c$ remains close to 1 throughout the domain (as in Figure 3f).

Figure 3. Numerical and analytic solutions of the sharp-interface model (3.4), for small $\alpha$ , with $\gamma =5$ (top), $\gamma =1$ (middle) and $\gamma =0.2$ (bottom). On the left we show the position, $h(t)$ , of the reacting interface against time $t$ for various $\gamma$ and for various small values of $\alpha$ , along with the leading-order asymptotic solution (3.7) in the limit of small $\alpha$ . On the right, we show the cleanser concentration profiles $c(y)$ at times $t_D/4$ , $t_D/2$ , $3t_D/4$ and $t_D$ (the decontamination time, at which $h=1$ ), with the arrows showing the direction of increasing time, for the corresponding $\gamma$ and $\alpha$ , and again showing the asymptotic small- $\alpha$ solution (3.7).

3.5. Strong cleanser, $\alpha =O(1)$

With $\alpha =O(1)$ the cleanser applied at the surface of the porous medium has a similar concentration to the (neat) agent. We can make analytical progress in both the limit of shallow ( $\gamma \ll 1$ ) and deep ( $\gamma \gg 1$ ) agent spills, corresponding to cases ii and iii, respectively, in Figure 2. The case $\alpha,\gamma =O(1)$ is intractable analytically, but we will also show numerical solutions in this regime in Figure 5.

3.5.1. Case ii: shallow spills, $\gamma \ll 1$

For shallow spills, with $\gamma \ll 1$ , we look for a solution of (3.4) of the form

(3.11) \begin{align} c=c_0+\gamma c_1+\ldots, \qquad h=h_0+\gamma h_1+\ldots . \end{align}

Substituting these expansions (3.11) into (3.4), we find that, at leading order, $c_0$ and $h_0$ are given by

(3.12) \begin{align} c_0=1, \qquad h_0=t. \end{align}

Continuing to next order in $\gamma$ , we find that

(3.13) \begin{align} c= 1-\gamma (1+\alpha )y+O\left(\gamma ^2\right), \qquad h=t-\frac{\gamma (1+\alpha )}{2}t^2+O\left(\gamma ^2\right), \end{align}

which are in fact the early-time expansions (3.5). Thus, for sufficiently shallow spills, the linear early-time approximation is valid for the entirety of the decontamination process.

We therefore see that, in this shallow-spill case, the reaction front moves with constant velocity to leading order (in agreement with case ia), and that the concentration profile is independent of time at both leading and first order in $\gamma$ . The leading-order solution is shown in Figure 5a and b and discussed further in Section 3.5.2 below.

The decontamination time, at which $h=1$ , is given by

(3.14) \begin{align} t_D=1+\frac{\gamma (1+\alpha )}{2}+O\left(\gamma ^2\right). \end{align}

We note that we regain the dilute-cleanser solution (3.8) in the limit $\alpha \rightarrow 0$ .

3.5.2. Case iii: deep spills $\gamma \gg 1$

In the case of deep agent spills, so that $\gamma \gg 1$ , we note from (3.4a) that the diffusive flux of cleanser through the domain is a limiting process over the long lengthscale of the spill. The decontamination therefore occurs over a longer timescale, and so we make the change of variables

(3.15) \begin{align} t=\gamma T, \end{align}

for $T=O(1)$ , and we write $h(t)=H(T)$ . On this new timescale, the model (3.4) becomes

(3.16a) \begin{align} \alpha c_T&=c_{yy} & &\text{ for }0\lt y\lt H(T), \end{align}
(3.16b) \begin{align} c&=1 & &\text{ at }y=0, \end{align}
(3.16c) \begin{align} H_T&=-\frac{c_y}{1+\alpha c} \text{ and } H_{T}=\gamma c & &\text{ at }y=H(T), \end{align}

where we have combined the boundary conditions at $H(T)$ .

We look for an expansion of the form

(3.17) \begin{align} c=c_0+\frac{1}{\gamma }c_1+O\left (\frac{1}{\gamma ^2}\right ), \qquad H=H_0+\frac{1}{\gamma }H_1+O\left (\frac{1}{\gamma ^2}\right ). \end{align}

Substituting (3.17) into (3.16), we find that, at leading order in $\gamma ^{-1}$ ,

(3.18a) \begin{align} \alpha c_{0T}&=c_{0yy} & &\text{ for }0\lt y\lt H_0(T), \end{align}
(3.18b) \begin{align} c_0&=1 & &\text{ at }y=0, \end{align}
(3.18c) \begin{align} c_{0y}&=-H_{0T} \text{ and } c_0=0 & &\text{ at }y=H_0(T). \end{align}

This leading-order model admits a similarity solution, of the form

(3.19) \begin{align} c_0=1-\frac{\textrm{erf}\left (\frac{\sqrt{\alpha }y}{2\sqrt{T}}\right )}{\textrm{erf}\left (\frac{\sqrt{\alpha }\eta _*}{2}\right )}, \qquad H_0=\eta _*\sqrt{T}, \end{align}

where the constant $\eta _*$ satisfies the algebraic equation

(3.20) \begin{align} 2\sqrt{\frac{\alpha }{\pi }}=\eta _*\exp \left (\frac{\alpha \eta _*^2}{4}\right )\textrm{erf}\left (\frac{\sqrt{\alpha }\eta _*}{2}\right ). \end{align}

Equation (3.20) may be solved numerically for $\eta _*$ for any given $\alpha$ . The solution for a range of $\alpha$ is shown in Figure 4. In particular, we note that as $\alpha \rightarrow 0$ , $\eta _*\rightarrow \sqrt{2}$ , (indeed for $\alpha \lt 1$ we find $\eta _*$ very close to $\sqrt{2}$ ) and so we recover the solution (3.10) for the dilute-cleanser, deep-spill case, in this limit.

Figure 4. The numerically computed solution $\eta _*$ of (3.20), as a function of $\alpha$ .

The leading-order decontamination time is

(3.21) \begin{align} T_D=\frac{1}{\eta _*^2}+O\left (\frac{1}{\gamma }\right ), \end{align}

so that, in terms of our original time variable,

(3.22) \begin{align} t_D=\frac{\gamma }{\eta _*^2}+O(1). \end{align}

In Figure 5, we show solutions of the sharp-interface model (3.4) for $\alpha =1$ , and various values of $\gamma$ . We see excellent agreement between the numerical solutions and both the analytic solution (3.13) in the small- $\gamma$ limit (in Figure 5a and b) and the analytic solution (3.19) in the large- $\gamma$ limit (in Figure 5c and d).

Figure 5. Numerical and analytic solutions of the sharp-interface model (3.4), examining the effect of $\gamma$ , with $\alpha =1$ throughout. In (a) and (b), we show the small- $\gamma$ numerical solutions of (3.4), with the leading-order small- $\gamma$ asymptotic solution (3.12). In (c) and (d), we show large- $\gamma$ solutions with the leading-order large- $\gamma$ asymptotic solution (3.19). In (a) and (c), we show the motion of the reaction front $h(t)$ . In (b) and (d), we show the cleanser concentration profiles at times $t_D/4$ , $t_D/2$ , $3t_D/4$ and $t_D$ (the decontamination time, at which $h=1$ ), with arrows indicating increasing time. In (b), we mark the end of each line with a circle for the cases $\gamma =0.1$ and $\gamma =0.01$ , as otherwise the curves are difficult to distinguish.

3.6. Decontamination times and decontamination efficiency

There are two practical metrics we can evaluate to assess the quality of the decontamination process. These are:

(i) the decontamination time, which is the time $t_D$ at which all the agent has been reacted away, that is, when $h(t_D)=1$ . In the various asymptotic limits above, we found the analytical expressions (3.8), (3.14) and (3.22) for the leading-order decontamination time in cases i, ii and iii, respectively.

(ii) the efficiency of the decontamination process, in terms of the amount of cleanser chemical used. We would ideally like to minimise the excess cleanser introduced to the system, since (although less toxic than the agent) the cleanser can also be harmful to people or to the environment. Since we have assumed the reaction has stoichiometry 1 (one mole of agent reacts with one mole of cleanser), the minimum number of moles of cleanser that are required to neutralise the entire agent spill (per unit surface area of the spill) is

(3.23) \begin{align} N=\frac{v}{\chi }\quad \text{mol}\,\text{m}^{-2}. \end{align}

The amount of cleanser per unit cross-sectional area that enters the porous medium over the entire decontamination process, however, is

(3.24) \begin{align} \text{total flux in}= - \frac{D_0\phi }{\chi \beta k}\int _{t=0}^{t_D}c_y|_{y=0}\,\mathrm{d}t \quad \text{mol}\,\text{m}^{-2}\,. \end{align}

We define the decontamination efficiency to be the ratio of the minimum amount of cleanser required to the flux of cleanser into the porous medium, namely

(3.25) \begin{align} \text{efficiency}=\frac{N}{\text{total flux in}}=\frac{\gamma }{-\int _{t=0}^{t_D}c_y|_{y=0}\,\mathrm{d}t}\,, \end{align}

where we have made use of (3.1) and (3.23). The efficiency defined by (3.25) is a dimensionless value between 0 and 1, with 1 being perfectly efficient decontamination. We note that we may compute the leading-order decontamination efficiencies in the various asymptotic limits studied above and find

(3.26) \begin{align} \text{efficiency}=\begin{cases}1+O(\alpha ) &\text{ for case i, }\alpha \ll 1,\\[5pt] \dfrac{1}{1+\alpha }+O(\gamma ) & \text{ for case ii, }\alpha =O(1),\gamma \ll 1,\\[9pt] \exp \left ({-}\dfrac{\alpha \eta _*^2}{4}\right )+O\left (\gamma ^{-1}\right ) &\text{ for case iii, }\alpha =O(1),\gamma \gg 1, \end{cases} \end{align}

where $\eta _*(\alpha )$ is the solution of (3.20).

In Figure 6, we show the variation of the decontamination time (left) and decontamination efficiency (right) with the system parameters and compare numerical solutions of (3.4) to these limiting analytic expressions. In Figure 6a, we show the variation of the numerically computed decontamination times as a function of $\gamma$ . In particular, since $\gamma$ is proportional to $L_{\text{SI}}$ , we may interpret this plot as the effect of the spill depth $L_{SI}$ on the decontamination time (since for a particular agent, cleanser and porous material, $\beta k$ and $D_0$ are constant). For small $\alpha$ , we observe that $t_D$ behaves like the case-i solution, $1+\gamma/2$ , as given by (3.8). For larger $\alpha =O(1)$ the limiting behaviours (3.14) and (3.22) capture the behaviour of the numerical solution well. In particular, we see that, for both $\alpha \ll 1$ and $\alpha =O(1)$ , the dimensionless decontamination time appears to increase linearly with spill depth. Since the time scaling used in (3.3) is also linear in $L_{\text{SI}}$ , we therefore conclude that the dimensional decontamination time increases quadratically with spill depth $L_{\text{SI}}$ .

Figure 6. Sharp-interface model decontamination time $t_D$ or scaled decontamination time $t_D/\alpha \gamma$ (left) and decontamination efficiency (right) as the system parameters $\alpha$ and $\gamma$ are varied. Numerical solutions of (3.4) are shown in colour, and the analytical expressions (3.8), (3.14) and (3.22) for the decontamination time, and (3.26) for the efficiency, in the various asymptotic limits, are shown in black (dashed, dotted or solid).

We note that the time scaling used in (3.3) also depends on both the reaction rate $k$ and the applied cleanser concentration $c_0$ . In order to draw physical conclusions about the impact of these parameters on the decontamination time, in Figure 6c and e we therefore plot the scaled decontamination time

(3.27) \begin{align} \frac{t_D}{\alpha \gamma }=\frac{D_0}{\chi c_0\beta k L_{\text{SI}}}t_D, \end{align}

which has the same dependence on $c_0$ and $ k$ as the time scaling used in (3.3). Thus, the dependence of $t_D/\alpha \gamma$ on $ k$ is the same as the dependence of the dimensional decontamination time on $k$ . In Figure 6c, we see that this scaled decontamination time decreases with $\gamma$ so that (as we might expect) faster reaction rates $k$ result in faster decontamination. However, we observe limiting returns, with the scaled decontamination time approaching an ( $\alpha$ -dependent) constant for larger values of $\gamma$ . This is because as $\gamma$ increases the reaction becomes almost instantaneous and the problem is no longer reaction-limited, but instead becomes limited by the supply of cleanser. In Figure 6e, we show the variation of the scaled decontamination time $t_D/\alpha \gamma$ with $\alpha$ (which is proportional to the applied cleanser concentration $c_0$ , and so we may interpret increasing $\alpha$ with increasing the applied cleanser concentration). For larger $\alpha$ , and so stronger cleanser concentrations, we observe faster decontamination. As $\alpha$ becomes small, the scaled decontamination time grows like $\alpha ^{-1}$ : the small- $\alpha$ expression, which from (3.8) is given by $\alpha ^{-1}(\gamma ^{-1}+1/2)$ , appears to be a good approximation for $\alpha \lt 1$ .

In Figure 6b and d, we show the variation of the decontamination efficiency (3.25) with $\gamma$ and $\alpha$ , respectively. For this sharp-interface case, the efficiency can also be thought of as a measure of how much cleanser is left inside the porous medium, in $y\in [0,1]$ , at time $t_D$ . We see in Figure 6b that the efficiency increases with $\gamma$ , since we recall (from Figure 5d, say) that, for larger $\gamma$ , we expect lower values of $c$ , and thus less cleanser is contained in $y\in [0,1]$ . The variation with $\gamma$ is most noticeable for larger $\alpha$ . For $\alpha \ll 1$ , the efficiency approaches 1. This is because for small applied cleanser concentrations $c_0$ , there are fewer moles of cleanser chemical left within the porous material at the decontamination time. From Figure 6d, we see that the small- $\gamma$ and large- $\gamma$ expressions (3.26) capture the decrease of the efficiency with $\alpha$ well. Moreover, we observe that the numerically simulated efficiencies shown in Figure 6b and d appear to be bounded below and above by the small- and large- $\gamma$ asymptotic formulae, respectively.

In summary, we have seen that increasing the reaction rate $\beta k$ decreases the decontamination time (albeit with limited returns) and also improves the efficiency of the system. Increasing the concentration of the applied cleanser $c_0$ results in faster decontamination times but reduced decontamination efficiency. In practice, therefore, there is a pay-off between cleanser wastage and speed of decontamination. To illustrate this, we show the variation of the efficiency with the decontamination time (parameterised by $\alpha$ or $c_0$ ) in Figure 6f. The applied cleanser concentration should be chosen to strike a reasonable balance between fast decontamination times and limited waste of the cleanser.

4. Agent-on-walls case

For the second of our two decontamination models, the agent is stuck to the solid porous structure in a layer of thickness $R$ , as shown in Figure 1 (right). The agent-coated solid circles are surrounded by cleanser solution, with cleanser concentration again given by $c$ .

4.1. Model statement

In [Reference Luckins, Breward, Griffiths and Wilmott10], we showed that in regions with $R=0$ the cleanser concentration $c$ satisfies

(4.1a) \begin{align} c_t=D_0c_{yy}, \end{align}

while if $R\gt 0$ we have the reaction–diffusion problem

(4.1b) \begin{align} c_t&=\frac{1}{\mathcal{V}(R)}\left (\mathcal{D}(R)\mathcal{V}(R)c_y\right )_y-\mathcal{F}(R)k\left (1+\chi c\right )c, \end{align}
(4.1c) \begin{align} R_t&=-\chi k c, \end{align}

where $\mathcal{D}(R)$ is the diffusivity of the cleanser through the available pore space when the agent thickness is $R$ (this is found by the numerical solution of a cell problem which is stated in Appendix A and solved in, for instance, [Reference Bruna and Chapman1]), and $D_0=\mathcal{D}|_{R=0}$ , $\mathcal{V}(R)=1-\pi (r_0+R)^2/d^2$ , is the fraction of the cell volume occupied by the cleanser solution, and $\mathcal{F}(R)=2\pi (r_0+R)/d^2\mathcal{V}(R)$ is the ratio of the length of the agent–cleanser interface on the microscale to the volume of the cell occupied by cleanser.

Initially, we suppose that $c=0$ is uniform throughout $y\gt 0$ , and that the agent is uniformly distributed with initial layer thickness $R=R_0$ , through the region $0\lt y\lt L_{\text{AW}}$ , where

(4.2) \begin{align} L_{\text{AW}}=\frac{vd^2}{\pi \left(2R_0r_0+R_0^2\right)}, \end{align}

is the spill lengthscale. We note that, since we want the cleanser-occupied region to be initially connected, we require $R_0+r_0\lt d/2$ . Since the agent layer is being removed, the cleanser region is then connected for all time. In higher dimensions or for different microscale geometries, we would require a similar constraint on the initial agent-layer thickness to ensure this connectedness. For all $t\gt 0$ , as in the sharp-interface case, we impose a constant concentration $c=c_0$ of cleanser at $y=0$ , but as $y\rightarrow \infty$ we require $c\rightarrow 0$ .

4.2. Nondimensionalisation

In addition to the spill lengthscale $L_{AW}$ , there are various additional lengthscales in this model, compared with the sharp-interface model, namely the initial agent-layer thickness $R_0$ , the radius of the circular solid inclusions, $r_0$ , and the separation distance between these circular solid inclusions, $d$ . We will assume that $\rho _0\,:\!=\,r_0/d=O(1)$ (a property of the porous material only) is not particularly small. However, we will allow that ratio

(4.3) \begin{align} \delta =\frac{R_0}{d}, \end{align}

to be $O(1)$ or smaller, so that the layer of agent may be thin relative to the pore space. We also define the ratio of the pore lengthscale to the spill lengthscale as

(4.4) \begin{align} \epsilon =\frac{d}{L_{\text{AW}}}. \end{align}

We require $\epsilon \ll 1$ in order that our homogenised model (4.1) is valid, and so we will assume this to be true throughout this paper.

We use the scalings

(4.5a) \begin{align} c=c_0\hat{c}, \qquad R=R_0\hat{R}, \qquad y=L_{\text{AW}}\hat{y}, \qquad t=\frac{R_0}{\chi k c_0}\hat{t}, \end{align}

where we have chosen the timescale that balances (4.1c). The functions $\mathcal{D},\mathcal{V}$ and $\mathcal{F}$ are simply functions of geometry, although $\mathcal{D}$ must be calculated numerically by solving the cell problem given in Appendix A (as in, for instance, [Reference Bruna and Chapman1, Reference Dalwadi, Griffiths and Bruna3]). We make the scalings

(4.5b) \begin{align} \mathcal{D}&=D_0\hat{D}\left(\rho _0+\delta \hat{R}\right), \end{align}
(4.5c) \begin{align} \mathcal{V}&=V_0\hat{V}\left(\hat{R};\delta,\rho _0\right)=V_0\left (1-\delta \frac{\pi \left(2\rho _0\hat{R}+\delta \hat{R}^2\right)}{V_0}\right ), \end{align}
(4.5d) \begin{align} \mathcal{F}&=F_0\hat{F}\left(\hat{R};\delta,\rho _0\right)=F_0\left (1+\frac{\delta \hat{R}}{\rho _0}\right )\left (\frac{V_0}{1-\pi \left(\rho _0+\delta \hat{R}\right)^2}\right ), \end{align}

where $D_0$ , $V_0=1-\pi \rho _0^2$ and $F_0=2\pi \rho _0/dV_0$ are the values of $\mathcal{D}$ , $\mathcal{V}$ and $\mathcal{F}$ , respectively, when the agent thickness $\hat{R}=0$ .

Substituting (4.5) into (4.1) and dropping the hat notation, our model becomes

(4.6a) \begin{align} A\Gamma c_t&=\frac{1}{V(R)}\left (D(R)V(R)c_y\right )_y-\Gamma H(R)F(R)\left (1+\frac{\delta A}{f_1} c\right )c, \end{align}
(4.6b) \begin{align} R_t&=-H(R)c, \end{align}

in $y\gt 0$ , where $H(R)$ is the Heaviside function. The boundary conditions are

(4.6c) \begin{align} c=1 \qquad &\text{ at } \qquad y=0, \end{align}
(4.6d) \begin{align} c\rightarrow 0 \qquad &\text{ as } \qquad y\rightarrow \infty, \end{align}

while at $t=0$ ,

(4.6e) \begin{align} c=0 \quad \text{for }y\gt 0, \qquad \text{ and } \qquad R=\begin{cases}1 & \text{if }y\in [0,1],\\ 0 & \text{otherwise}. \end{cases} \end{align}

The dimensionless parameters in the model are the lengthscale ratios $\delta$ and $\epsilon$ , given by (4.3) and (4.4), and $\Gamma$ and $A$ , given by

(4.7) \begin{align} \Gamma =\frac{\gamma f_2}{\epsilon \delta }=\frac{kL_{\text{AW}}}{\epsilon f_1D_0}, \qquad A=\frac{\alpha f_1}{\delta }=\frac{c_0\chi f_1}{\delta }, \end{align}

where the sharp-interface parameters $\gamma$ and $\alpha$ are as defined in Table 2. The model also depends on the $O(1)$ geometric factors

(4.8) \begin{align} f_1=\frac{1-\pi \rho _0^2}{2\pi \rho _0}, \qquad f_2=\frac{2 \rho _0}{\beta (2\rho _0+\delta )}. \end{align}

The dimensionless parameters in (4.6) are summarised in Table 3. We plot $D$ , $V$ and $F$ as functions of $R$ for given $\rho _0$ and $\delta$ in Figure 7, using $D(R)$ as computed in [Reference Dalwadi, Griffiths and Bruna3]. We note that $V$ and $D$ decrease with $R$ , while $F$ increases with $R$ .

Table 3. Dimensionless parameters in the agent-on-walls model (4.6)

Figure 7. The variation of (a) $V$ and $D$ and (b) $F$ with agent thickness $R$ . We take $\rho _0=0.2$ throughout, so that the model is valid for $\delta \in (0,0.3)$ . Three values $\delta =0.01$ (dash-dotted curves), $\delta =0.1$ (dashed curves) and $\delta =0.3$ (solid curves) are shown.

4.3. Comparison with the sharp-interface model

Before analysing solutions of (4.6), we first compare the model structure, scalings and dimensionless parameters with that of the sharp-interface model (3.4).

Firstly, because we are comparing the models where the initial volume of agent is the same, we note that the depth of porous medium to be decontaminated is smaller in the sharp-interface case than the agent-on-walls case, and $L_{\text{SI}}/L_{\text{AW}}=O(\delta )$ . This is because the agent occupies only a fraction of the pore space in the agent-on-walls case, rather than saturating it. The timescales of the decontamination are likewise different with $[t_{\text{SI}}]/[t_{\text{AW}}]=O(\delta/\epsilon )$ , (where $[t_{\text{SI}}]$ and $[t_{\text{AW}}]$ are the time scalings used in (3.3) and (4.5), respectively). We recall that we require $\epsilon \ll 1$ for our agent-on-pores model to be valid, while $\delta$ may be small or of order one. If $\delta/\epsilon \gg 1$ , then the decontamination timescale is much shorter for the agent-on-walls case than for the sharp-interface case. This is despite the fact that the spill is deeper in the agent-on-walls case than for the sharp-interface case. The reason the agent-on-walls decontamination is faster is because the surface area of the interfaces between the fluids is much greater in the agent-on-pores case, and so the overall decontamination can occur much more quickly. We note, however, that with our choice of timescales we have assumed that the decontamination is limited by the reaction rate, rather than by diffusion of the cleanser. In each case, this choice is valid for sufficiently shallow agent spills but for deeper spills, where the transport of cleanser to the reacting interface(s) is also important, the timescale of decontamination is longer. Indeed for deep spills we will see that the decontamination timescale may actually be longer for the agent-on-walls model than for the sharp-interface case.

The parameters $A$ and $\Gamma$ in the agent-on-walls model are closely linked to $\alpha$ and $\gamma$ (defined in Table 2), respectively, in the sharp-interface model; the differences are purely geometric scaling factors. Specifically, $A$ and $\alpha$ in the agent-on-walls and sharp-interface models, respectively, characterise the amount of cleanser in the system compared with the local amount of agent to be removed. In the sharp-interface case, this is simply $\alpha =c_0\chi$ (essentially a ratio of the chemical concentrations), but in the agent-on-walls case we have $A=\alpha f_1/\delta$ , since the volume fraction of agent within the porespace is reduced by a factor of $\delta/f_1$ . Since we expect that $c_0\lesssim \chi$ , we only consider $\alpha$ to be $O(1)$ or smaller in the sharp-interface case. However, in the agent-on-walls case, $A$ may be small, order one or large, depending on the relative sizes of $\alpha$ and $\delta$ : if $\delta \ll \alpha$ so that $A\gg 1$ , we will find in Section 4.7 that diffusion is the limiting process of the decontamination. Both $\gamma$ and $\Gamma$ are a measure of the depth of the spill or speed of the reaction: they are the ratio of the spill depth to the diffusion lengthscale of the cleanser, over the reaction timescale.

As for the sharp-interface model, we can make asymptotic progress in certain parameter regimes. These are illustrated schematically in Figure 8, and referred to as cases I–III, which, as we will show, bear some resemblance to the limits i–iii in the sharp-interface case. We note that, for all of our asymptotic analyses of the agent-on-walls model, we will use $\delta \ll 1$ .

Figure 8. Schematic of the asymptotic limits considered for the agent-on-walls model. Case I is $A\ll 1$ , with sublimits $\Gamma \ll 1$ and $\Gamma \gg 1$ referred to as cases Ia and Ib, respectively. Case II is $A\gg 1$ and $\Gamma \ll 1$ , while case III is $A\gg 1$ and $\Gamma \gg 1$ . For all these asymptotic limits, we additionally assume $\delta \ll 1$ .

4.4. Numerical method

As for the sharp-interface model (3.4), we solve the model (4.6) using the method of lines, with 100 spatial meshpoints (beyond which the results are indistinguishable), implemented in Matlab with the inbuilt Ordinary Differential Equation (ODE) solver ode15s for the timestepping. To discretise the spatial derivatives of $c$ , we use a central finite difference scheme.

We note that the lengthscale of the diffusion of cleanser in $y\gt 1$ (where $R=0$ ) may vary drastically from the lengthscale of the agent-occupied region which has been nondimensionalised to reside between $y=0$ and $y=1$ , depending on the size of the parameter group $A\Gamma$ and also on the time $t$ . Specifically, the distance that cleanser has travelled beyond $y=1$ is on the order of $\sqrt{t/A \Gamma }$ , which may be large or small depending on the parameter regime of interest. In order that our numerical scheme is accurate and efficient in both regions $y\in [0,1]$ and $y\gt 1$ , we therefore discretise (4.6) in the two regions $y\lt 1$ and $y\gt 1$ separately and impose conditions of continuity of cleanser concentration and cleanser flux between the two (at $y=1$ ). For $y\lt 1$ , we use a uniform mesh and impose the boundary condition (4.6c) straightforwardly. For $y\gt 1$ , we make the change of variables

(4.9) \begin{align} y=1+\sqrt{\frac{t}{A\Gamma }}\eta, \qquad c(y,t)=\bar{c}(\eta,t). \end{align}

The cleanser concentration $\bar{c}$ in $y\gt 1$ (or $\eta \gt 0$ ) then satisfies

(4.10) \begin{align} t\bar{c}_t-\frac{\eta }{2}\bar{c}_\eta = \bar{c}_{\eta \eta }. \end{align}

The two conditions imposed at $y=1$ ( $\eta =0$ ) to couple the two systems together are continuity of cleanser concentration and flux, namely

(4.11) \begin{align} c(1,t)=\bar{c}(0,t), \qquad D(R)V(R)c_y(1,t)= \sqrt{A\Gamma } \frac{\bar{c}_\eta (0,t)}{\sqrt{t}}. \end{align}

We solve this problem (4.10) and (4.11) on the domain $\eta =[0,5]$ , imposing the Dirichlet condition (4.6d) at $\eta =5$ (which we find to be sufficiently large for numerical accuracy). To do this, we discretise the equation (4.10) with central differences for the diffusion term, and a first-order upwind (forward) scheme for the artificial advection term due to the growing coordinate system, to ensure stability. (For the practical implementation, we rescale $\eta$ onto $[0,1]$ and so use the same mesh for $\eta$ as for $y$ .) We note that the apparent singularity at $t=0$ in (4.10) and (4.11) does not cause numerical difficulty, since $\bar{c}=\bar{c}_\eta =0$ when $t=0$ .

4.5. Preliminary remarks

4.5.1. The first time at which $R=0$ , and the emergence of a decontamination front

We first consider the early-time behaviour of (4.6) at $y=0$ , where we have $c=1$ . Initially, $R=1$ for $y\geq 0$ and, while $R\gt 0$ everywhere, from (4.6b) we have

(4.12) \begin{align} R_t=-c. \end{align}

Thus $R_t|_{y=0}=-1,$ and so $R=1-t$ while $t\lt 1$ at $y=0$ . The first occurrence of $R=0$ is therefore at time $t=1$ , regardless of the size of the dimensionless parameters. This justifies our choice of timescale used in (4.5), as this is the minimum time that the decontamination can take.

We define the position, $y=y_*(t)$ to be the largest value of $y$ for which $R=0$ , so that $R=0$ for $y\lt y_*$ and $R\gt 0$ for $y\gt y_*$ . This position is analogous to the position of $y=h(t)$ of the reacting interface in the sharp-interface model, and we view it as the decontamination front. In particular, the decontamination time, that is, the time when the decontamination process is complete, occurs at the time $t_D$ at which $y_*=1$ . We note that $y_*=0$ when $t=1$ , unlike the reaction front $y=h(t)$ for the sharp-interface case, which satisfies $h=0$ at $t=0$ .

4.5.2. Behaviour for $y\gt 1$

For $y\gt 1$ where $R=0$ for all time, we may make analytical progress. Here, $c$ satisfies the diffusion equation

(4.13) \begin{align} A\Gamma c_t=c_{yy}, \end{align}

with the requirement that $c=0$ at $t=0$ and $c\rightarrow 0$ as $y\rightarrow \infty$ . This problem has solution (using the Green’s function approach, [Reference Keener8])

(4.14) \begin{align} c(y,t)=\sqrt{\frac{A\Gamma }{4 \pi }}\int _{\tau =0}^{t}\frac{c(1,\tau )(y-1)}{(t-\tau )^{3/2}}\exp \left({-}\frac{A\Gamma (y-1)^2}{4(t-\tau )}\right)\,\mathrm{d}\tau, \end{align}

in terms of the value of $c$ at $y=1$ at all previous times. As well as continuity at $y=1$ , we require that the solution has a continuous diffusive flux of cleanser at $y=1$ . Since $R$ is discontinuous at $y=1$ (while $R\gt 0$ for $y\lt 1$ ), the diffusivity and pore-volume occupied by cleanser are discontinuous. Continuity of diffusive flux requires

(4.15) \begin{align} D(R)&V(R)c_y|_{y=1^-}=c_y|_{y=1^+}\nonumber \\ &=\sqrt{\frac{A\Gamma }{4 \pi }}\,\lim _{y\rightarrow 1}\left (\int _{\tau =0}^{t}\frac{c(1,\tau )}{(t-\tau )^{3/2}}\exp \left({-}\frac{A\Gamma (y-1)^2}{4(t-\tau )}\right)\left (1-\frac{A\Gamma (y-1)^2}{2(t-\tau )}\right )\,\mathrm{d}\tau \right ),\end{align}

This expression (4.15) will be helpful in our analytic analyses to come, as it quantifies the size of the flux of cleanser leaving the agent-occupied region through $y=1$ , namely it is on the order of $\sqrt{A\Gamma/ t}$ .

4.5.3. The limit $\delta \ll 1$

In order to make analytic progress in the following sections, we regularly consider the limit $\delta \ll 1$ , in which the initial agent-layer thickness $R_0$ is small relative to the separation between the solid circles making up the microstructure. With $\delta \ll 1$ ,

(4.16) \begin{align} D(R)=1+O(\delta ), \qquad V(R)=1+O(\delta ), \qquad F(R)=1+O(\delta ), \end{align}

are all constant at leading order in $\delta$ , and – so long as $A\delta \ll 1$ – the model equations (4.6a)–(4.6b) reduce at leading order to

(4.17a) \begin{align} A\Gamma c_t&=c_{yy}-\Gamma H(R)c, \end{align}
(4.17b) \begin{align} R_t&=-H(R)c. \end{align}

4.6. Case $A\ll 1$ , reaction-limited decontamination

In this section, we suppose $\alpha \ll \delta$ so that $A\ll 1$ . This is the dilute-cleanser limit, since $\delta$ is order 1 or smaller, and so we have $\alpha \ll 1$ . As in the sharp-interface case $\alpha \ll 1$ , we will see that for $A\ll 1$ equation (4.6a) is dominated by diffusion (regardless of the size of $\Gamma$ ), so that diffusion is a quasi-steady process.

4.6.1. Case I: analytical solutions in the limit $\delta \ll 1$

In order to make analytical progress, we further assume that $\delta \ll 1$ . In this limit, where $A\ll 1$ and $\delta \ll 1$ (case I as shown in the schematic of Figure 8), the leading-order model is (4.17) in $y\gt 0$ , with boundary and initial conditions (4.6c)–(4.6e). We note that in the analysis below we will replace the far-field boundary condition (4.6d) with the condition (4.15), and solve (4.17) on $y\in [0,1]$ only. In particular, we note that the problem (4.6c) for $c$ decouples from the problem for $R$ (so long as we know the region for which $R\gt 0$ ).

We now explore solutions of (4.17) for $A\ll 1$ , for various sizes of $\Gamma$ . We begin with case I ( $\Gamma =O(1)$ ) before then looking at the sublimit $\Gamma \ll 1$ which is case Ia. We then find additional structures in case Ib for which $\Gamma \gg 1$ .

With $\Gamma =O(1)$ , to leading order in $A\ll 1$ and $\delta \ll 1$ we have

(4.18a) \begin{align} c_{yy}&=\Gamma H(R)c,\\[-10pt]\nonumber \end{align}
(4.18b) \begin{align} R_t&=-H(R)c. \end{align}

The flux of cleanser through $y=1$ is given by (4.15) and is of order $\sqrt{A\Gamma }$ . Thus, to leading order in $\sqrt{A}\ll 1$ , we replace (4.6d) with the boundary condition

(4.18c) \begin{align} c_y=0 \qquad \text{ on }y=1. \end{align}

This, along with the boundary condition (4.6c), and the initial conditions (4.6e), closes the leading-order problem, and we need not consider $y\gt 1$ explicitly.

For $t\lt 1$ when $R\gt 0$ throughout $y\in (0,1)$ , we find that (4.18) has solution

(4.19) \begin{align} c=\frac{\cosh \left(\sqrt{\Gamma }(1-y)\right)}{\cosh \left(\sqrt{\Gamma }\right)}, \qquad R=1-\frac{\cosh \left(\sqrt{\Gamma }(1-y)\right)}{\cosh \left(\sqrt{\Gamma }\right)}t. \end{align}

We note, in particular, that $c$ is independent of $t$ , since the reaction–diffusion process is quasi-steady, and for $\delta \ll 1$ the depletion of agent has no effect on the transport of $c$ at leading order. (We also note that there is an additional early-time problem when $t=O(A)$ that describes how $c$ evolves to its quasi-steady profile, but which is not otherwise of interest.)

For $t\gt 1$ , we also have a decontaminated region $y\lt y_*(t)$ for which $R=0$ , in addition to a region $y\in [y_*(t),1]$ for which $R\gt 0$ . Solving (4.18a) for $c$ in $y\lt y_*(t)$ and $y\gt y_*(t)$ separately, and imposing the boundary conditions at $y=0$ and $y=1$ , as well as continuity of both $c$ and $c_y$ at $y=y_*(t)$ , we find the solution

(4.20) \begin{align} c=\begin{cases} 1-\dfrac{\sqrt{\Gamma }y}{\coth \left(\sqrt{\Gamma }(1-y_*)\right)+\sqrt{\Gamma }y_*} & \text{ if }y\in [0,y_*(t)],\\[18pt] \dfrac{\cosh \left(\sqrt{\Gamma }(1-y)\right)}{\cosh \left(\sqrt{\Gamma }(1-y_*)\right)+\sqrt{\Gamma }y_*\sinh \left(\sqrt{\Gamma }(1-y_*)\right)} & \text{ if }y\in [y_*(t),1], \end{cases} \end{align}

given in terms of $y_*(t)$ . The leading-order equation (4.18b) for $R$ is therefore

(4.21) \begin{align} R_t=-\frac{\cosh \left(\sqrt{\Gamma }(1-y)\right)}{\cosh \left(\sqrt{\Gamma }(1-y_*(t))\right)+\sqrt{\Gamma }y_*(t)\sinh \left(\sqrt{\Gamma }(1-y_*(t))\right)}, \end{align}

for $y\gt y_*(t)$ , which must be solved along with the initial condition $R=R_1(y)\,:\!=\,R(1,y)$ at $t=1$ from (4.19), and the fact that $R=0$ at $y=y_*(t)$ . We look for a separable solution to (4.21) of the form

(4.22) \begin{align} R(t,y)=R_1(y)-S(t)U(y). \end{align}

Substituting (4.22) into (4.21), we see that

(4.23) \begin{align} U(y)=\cosh \left(\sqrt{\Gamma }(1-y)\right), \end{align}

and $S$ satisfies

(4.24) \begin{align} S^{\prime}(t)=\frac{1}{\cosh \left(\sqrt{\Gamma }(1-y_*)\right)+\sqrt{\Gamma }y_*\sinh \left(\sqrt{\Gamma }(1-y_*)\right)} . \end{align}

Furthermore, in order that $R=0$ at $y=y_*$ , we must have that

(4.25) \begin{align} S(t)=\frac{R_1(y_*(t))}{U(y_*(t))}=\textrm{sech}\left(\sqrt{\Gamma }(1-y_*)\right)-\textrm{sech}\left(\sqrt{\Gamma }\right) . \end{align}

We note that $S=0$ when $y_*=0$ , as required for the initial condition. Substituting (4.25) into (4.24) gives an ODE for $y_*(t)$ , namely

(4.26) \begin{align} y^{\prime}_*(t)=\frac{\coth \left(\sqrt{\Gamma }(1-y_*)\right)}{\sqrt{\Gamma }\left (1+\sqrt{\Gamma }y_*\tanh \left(\sqrt{\Gamma }(1-y_*)\right)\right )} . \end{align}

This ODE can be solved numerically, with the initial condition $y_*=0$ at $t=1$ . The leading-order solutions for $c$ and $R$ are then given in terms of $y_*$ by (4.20) and

(4.27) \begin{align} R=1-\frac{\cosh \left(\sqrt{\Gamma }(1-y)\right)}{\cosh \left(\sqrt{\Gamma }(1-y_*)\right)} \qquad \text{ for }y\in [y_*(t),1] . \end{align}

We present solutions to (4.26) for various $\Gamma$ in Figures 9 and 11. Solutions with $\Gamma =1$ are shown in Figure 9a, c and e, along with numerical solutions of (4.17) with $A=10^{-3}$ and $\delta \ll 1$ , with good agreement. For $t\gt 1$ , we see in Figure 9a that the speed of the decontamination front $y^{\prime}_*$ is increasing. We note from Figure 9c that initial transients in $c$ are not captured by the quasi-steady expression (4.20), but that this is a good approximation thereafter, both for $t\lt 1$ and $t\gt 1$ . The agent-layer thickness, $R$ , shown in Figure 9e decreases from its initial value of 1 and decreases most quickly in regions with higher cleanser concentration (close to $y=0$ ). The analytic solution (4.27) for $R$ is in good agreement with the numerical solution, although the error appears to increase with time. We note that the largest source of error in the small- $A$ approximation is the $O(\sqrt{A})$ error in the flux boundary condition applied at $y=1$ (since the other error from the time-derivative term in (4.17) is of $O(A)$ ); we accordingly see much better agreement at $y=0$ than at $y=1$ , and better agreement at earlier times. We note that the solution (4.27) for $R$ bears some resemblance to a travelling wave. However, it is not a travelling wave in the traditional sense: as well as moving with the decontamination front $y_*(t)$ , the profile $R$ changes with time, decreasing in amplitude. Furthermore, the wave-speed $y^{\prime}_*(t)$ is not constant in time: the wave-speed may initially decrease for $t\gt 1$ , with $y^{\prime\prime}_*(t)\lt 0$ (this is true for the large- $\Gamma$ case considered below), but $y^{\prime}_*$ is always seen to be increasing as $y_*$ approaches the end of the contaminated region at $y=1$ .

Figure 9. Numerical and asymptotic solution of (4.6) in the case $A\ll 1$ , $\delta \ll 1$ and for $\Gamma =1$ (left) and $\Gamma =0.1$ (right). We show the decontamination front $y_*(t)$ (top), cleanser profiles $c$ (middle) and agent-layer thickness $R$ (bottom). The numerical solutions (coloured curves) are computed in the limit $\delta \ll 1$ , making the assumption that $D=V=F=1$ , and taking $A=10^{-3}$ . In (b)–(f), solutions are shown at times $0.05t_D$ , $t_D/4$ , $t_D/2$ , $3t_D/4$ and $0.95t_D$ (where $t_D$ is the decontamination time, at which $y_*=1$ ), and arrows show the direction of increasing time. Analytical solutions are shown in black. In (c) and (d), thick dotted black curves are the analytical solution (4.20) and in (e) and (f) thick dashed curves are (4.27), where $y_*(t)$ is the (numerically computed) solution of (4.26). The thin dashed lines in (f) are the approximation (4.31) to $R$ , and the thick black line $c=1$ in (d) is the approximation (4.31) to $c$ in the small- $\Gamma$ sublimit.

We may also make further analytical progress in case Ia for which $\Gamma \ll 1$ . This is a sublimit of case I ( $\Gamma =O(1)$ ) considered above. For $\Gamma \ll 1$ we see that (4.26) reduces to

(4.28) \begin{align} y^{\prime}_*(t)\sim \frac{1}{\Gamma (1-y_*)}, \end{align}

with solution

(4.29) \begin{align} y_*(t)\sim 1-\sqrt{1+\frac{2}{\Gamma }(1-t)}. \end{align}

This expression is valid only for $1\lt t\lt t_D$ where

(4.30) \begin{align} t_D=1+\frac{\Gamma }{2}, \end{align}

is the (approximate) decontamination time in this limit. We further note that, when $\Gamma \ll 1$ , the analytic solution given by (4.20) and (4.27) reduces to simply

(4.31) \begin{align} c=1, \qquad R=1-t, \end{align}

at leading order, for both $t\lt 1$ and $t\in [1,t_D]$ .

In Figure 9 (right), we compare this solution (4.31) in the $\Gamma \ll 1$ sublimit to the full asymptotic solution for $\Gamma =O(1)$ , given by (4.20) and (4.27), and to numerical solutions of (4.17) in the limit $\delta \rightarrow 0$ , for small $A=10^{-3}$ , taking the value $\Gamma =0.1$ . In Figure 9b, we see that the small- $\Gamma$ approximation (4.29) for the decontamination front $y_*(t)$ agrees well with the numerical solution and the full solution of (4.26). The small- $\Gamma$ solutions (4.31) also approximate the $c$ and $R$ profiles well for $\Gamma =0.1$ , since $c$ is close to 1 and at each time the $R$ -profile is nearly flat.

We note that the limit $A\ll 1$ is analogous to the limit $\alpha \ll 1$ in the sharp-interface model. In particular, we see that, for small $\Gamma$ (analogous to small $\gamma$ in the sharp-interface model), the form of the decontamination time $t_D$ given by (4.30) is closely related to (3.14), in the limit $\alpha \ll 1$ of the sharp-interface model, despite the very different motion of the decontamination front in each case.

Case Ib: the limit $\Gamma \gg 1$

The case $\Gamma \gg 1$ corresponds to the deep-spill limit and is referred to as case Ib in Figure 8. We might expect that the time derivative in (4.17a) would become important, but – as in Section 3.4 for the sharp-interface case – we will see that this is not the case. For $\Gamma \gg 1$ . we will show that the timescale of the overall decontamination is longer than the reaction timescale, by a factor of $\Gamma$ . We explore the decontamination process in this limit in detail by examining the behaviour separately at ‘early’ times (which is the reaction timescale previously considered), at $O(\Gamma )$ times, and near the ‘end’ time through a temporal boundary layer analysis at the end of the decontamination. A schematic of the solution structure in this case Ib is given in Figure 10.

Figure 10. An illustration of the asymptotic solution structure for case Ib. Left: the decontamination front $y_*(t)$ against time, showing the three regions “early” ( $t=O(1)$ ), “late” ( $t=O(\Gamma )$ ) and “end” ( $1-y_*(t)=O(1/\sqrt{\Gamma })$ ) times. Right: the cleanser concentration $c$ and agent-layer thickness $R$ during the late-time stage.

To understand the early-time behaviour, we first consider the sublimit of (4.17) for which $\Gamma \gg 1$ , on an $O(1)$ timescale. Initially, $R\gt 0$ throughout $y\in [0,1]$ , and so, at leading order, we see from (4.17a) and (4.6e) that we must have $c=0$ and therefore, from (4.17b), that $R=1$ . This is incompatible with the boundary condition $c=1$ at $y=0$ , and so there must be a boundary layer near $y=0$ in which $c$ and $R$ vary. Making the change of variables

(4.32) \begin{align} y=\frac{Y}{\sqrt{\Gamma }}, \end{align}

equation (4.17a) becomes (at leading order in $\delta/\Gamma$ )

(4.33) \begin{align} A c_t=c_{YY}- H(R)c, \end{align}

so that, at leading order, we have

(4.34) \begin{align} c_{YY}=H(R)c, \end{align}

along with (4.17b) describing the evolution of $R$ . As $Y\rightarrow \infty$ , we must have $c\rightarrow 0$ and $R\rightarrow 1$ to match with the far-field behaviour. For $t\lt 1$ (when $R\gt 0$ everywhere in $Y\in [0,\infty )$ ), we find the solution

(4.35) \begin{align} c=\exp ({-}Y), \qquad R=1-t\exp ({-}Y). \end{align}

For $t\gt 1$ , we also have a decontaminated region $Y\lt Y_*(t)$ for which $R=0$ . Solving for $c$ in $Y\lt Y_*(t)$ and $Y\gt Y_*(t)$ separately, and imposing continuity of both $c$ and $c_Y$ at $Y=Y_*(t)$ , we find the solution

(4.36) \begin{align} c(Y,t)&=\begin{cases} 1-\dfrac{Y}{1+Y_*(t)} & \text{if }Y\lt Y_*(t),\\[15pt] \dfrac{\exp ({-}(Y-Y_*(t)))}{1+Y_*(t)} & \text{if }Y\gt Y_*(t). \end{cases} \end{align}

In $Y\gt Y_*(t)$ , the agent thickness $R$ then satisfies

(4.37) \begin{align} R_t=-\frac{\exp ({-}(Y-Y_*(t)))}{1+Y_*}. \end{align}

As in the $\Gamma =O(1)$ case, we look for a separable solution $R=1-U(t)S(Y)$ . Following the same method, we find

(4.38) \begin{align} S(Y)=\exp ({-}Y), \qquad U(t)=\frac{1}{S(Y_*)}=\exp \left (Y_*(t)\right ), \end{align}

so that

(4.39) \begin{align} R(Y,t)=1-\exp \left ({-}(Y-Y_*(t))\right ) \qquad \text{ for }Y\in [Y_*(t),\infty ). \end{align}

with $Y_*$ the solution of

(4.40) \begin{align} Y^{\prime}_*(t)=\frac{1}{1+Y_*(t)}. \end{align}

We note that this ODE for $Y_*$ is much simpler than the equivalent equation (4.26) in the case $\Gamma =O(1)$ , where effects from the end of the agent-occupied domain impacted on the cleanser concentration profile. The explicit solution to (4.40) with initial condition $Y_*=0$ at $t=1$ is

(4.41) \begin{align} Y_*(t)=-1+\sqrt{2t-1} \qquad \text{ for }t\gt 1. \end{align}

For $t\gt 1$ , the solution is therefore

(4.42a) \begin{align} c&=\begin{cases}1-\dfrac{Y}{Y_*(t)+1} &\text{if }Y\lt Y_*(t),\\[10pt] \dfrac{1}{Y_*(t)+1}\exp \left ({-}(Y-Y_*(t))\right ) &\text{if }Y\gt Y_*(t),\end{cases}\end{align}
(4.42b) \begin{align} R&=\begin{cases} 0 &\hspace{0.6cm}\text{if }Y\lt Y_*(t),\\[5pt] 1-\exp \left ({-}(Y-Y_*(t))\right ) &\hspace{0.6cm}\text{if }Y\gt Y_*(t).\end{cases}\end{align}

Unlike the $\Gamma =O(1)$ case given by (4.20)–(4.27), the wave-like solution $R$ does not change shape as it moves into the domain, although the speed of the wave motion, $Y^{\prime}_*(t)$ , is still non-uniform.

We now consider $O(\Gamma )$ time. We note that, for large $t$ , the early-time decontamination front given by (4.41) behaves like $Y_*\sim \sqrt{2t}$ , so that the position $y_*=Y_*/\sqrt{\Gamma }$ of the decontamination front becomes $O(1)$ after an $O(\Gamma )$ time. To examine the decontamination of the entire spill, we therefore change to an $O(\Gamma )$ timescale.

Making the change of variables $t=\Gamma T$ , equations (4.17a)–(4.17b) become (at leading order in $\delta$ )

(4.43a) \begin{align} A c_T&=c_{yy}-\Gamma H(R)c,\\[-5pt]\nonumber \end{align}
(4.43b) \begin{align} R_T&=-\Gamma H(R)c. \end{align}

The position, $y=y_*(T)$ , of the decontamination front is now of $O(1)$ . For $y\lt y_*$ we see that, to leading order (in $A$ and $\delta$ ), $c$ is linear, and given by

(4.44) \begin{align} c=1-B(T)y \end{align}

for some $B(T)$ to be determined.

For $y\gt y_*(T)$ , where $R\gt 0$ , we see that (4.43a) requires $c=O(1/\Gamma )$ to be small. We therefore again look for an $O(1/\sqrt{\Gamma })$ boundary layer in which $c=O(1)$ , this time centred on the moving interface, by making the change of variables

(4.45) \begin{align} y=y_*(T)+\frac{Y}{\sqrt{\Gamma }}, \qquad \text{for }Y\gt 0. \end{align}

This boundary layer structure is illustrated in Figure 10 (right). Substituting (4.45) into (4.43), we see that, in this boundary layer at $y=y_*(T)$ ,

(4.46a) \begin{align} A \left (c_T-\sqrt{\Gamma }y^{\prime}_*(T)c_Y\right )&=\Gamma \left (c_{YY}- H(R)c\right )\\[-5pt]\nonumber \end{align}
(4.46b) \begin{align} R_T-\sqrt{\Gamma }y^{\prime}_*(T)R_Y&=-\Gamma H(R)c. \end{align}

In particular, we see from (4.46b) that $c$ must be small within this boundary layer and of $O(1/\sqrt{\Gamma })$ . Setting $c=C(Y,T)/\sqrt{\Gamma }$ in (4.46), we find that, at leading order in $O(1/\sqrt{\Gamma })$ ,

(4.47) \begin{align} C_{YY}=C \qquad \text{ for }Y\gt 0, \end{align}

with solution

(4.48) \begin{align} C=D(T)e^{-Y}, \end{align}

in order that $C\rightarrow 0$ as $Y\rightarrow \infty$ . The two integration constants $B(T)$ and $D(T)$ are fixed by imposing continuity of both the cleanser concentration and its derivative at $y=y_*$ , corresponding to $Y=0$ . Specifically, we require $c=0$ and $c_y=C_Y$ at $y=y_*$ ( $Y=0$ ), which fix

(4.49) \begin{align} B(T)=D(T)=\frac{1}{y_*(T)}. \end{align}

In summary, and returning to the original spatial variable $y$ , the leading-order solution, $c$ , is therefore

(4.50) \begin{align} c(y,T)=\begin{cases} 1-\dfrac{y}{y_*(T)} &\text{ if }y\lt y_*(T),\\[10pt] \dfrac{1}{\sqrt{\Gamma }y_*(T)}\exp ({-}\sqrt{\Gamma }(y-y_*(T))) &\text{ if }y\gt y_*(T). \end{cases} \end{align}

The leading-order equation for $R$ is therefore

(4.51) \begin{align} y^{\prime}_*(T)R_Y=\frac{1}{y_*(T)}\exp ({-}Y), \end{align}

which has the wave-like solution

(4.52) \begin{align} R=1-\exp ({-}Y)=1-\exp \left (\sqrt{\Gamma }(y-y_*(T))\right ) \qquad \text{ with } \qquad y_*(T)=\sqrt{2T}. \end{align}

We note this solution is precisely the large- $t$ limit of the early-time solution (4.42); in particular the shape of $R$ is exactly the same as (4.42b). Indeed, the early-time solution (4.42) is a higher-order approximation of the late-time behaviour than (4.50) and (4.52), since there is no new physics entering the problem for $O(\Gamma )$ times. This late-time solution is illustrated schematically on the right of Figure 10.

As an aside, we notice that the leading-order solution given by (4.50) and (4.52) satisfies a reaction-limited, deep-spill, sharp-interface problem, namely

(4.53a) \begin{align} c_{yy}&=0 & \text{ for }y\lt y_*(T), \end{align}
(4.53b) \begin{align} c&=0 &\text{ at }y=y_*(T), \end{align}
(4.53c) \begin{align} y^{\prime}_*(T)&=-c_y, &\text{ at }y=y_*(T), \end{align}

with $c=1$ at $y=0$ , which is precisely the form of the sharp-interface model in case ib. Thus, for $A\ll 1$ , $\delta \ll 1$ , and $\Gamma \gg 1$ (case Ib), the agent-on-pores model behaves precisely like a sharp-interface model (in case ib).

The leading-order decontamination time predicted by (4.52), in the large- $\Gamma$ limit, is $t_D=\Gamma/2$ . Since, as we have already noted, the early-time solution (4.42) is, in fact, valid so long as $y_*$ is an $O(1)$ distance from 1, we may use (4.42) to find the first few correction terms to the decontamination time. Setting $Y_*=\sqrt{\Gamma }$ in (4.41), we find that the decontamination time is

(4.54) \begin{align} t_D^*=\frac{\Gamma }{2}+\sqrt{\Gamma }+1. \end{align}

However, we have assumed in our analysis that the decontamination front is an $O(1)$ distance away from the end of the agent spill at $y=1$ . When $1-y_*(T)$ becomes the same size as the width of the boundary layer, $O(1/\sqrt{\Gamma })$ , we will see end effects similar to those seen in the $\Gamma =O(1)$ case studied above, and we expect that this will provide an $O(\sqrt{\Gamma})$ correction to the decontamination time predicted by (4.54). Indeed, we study the decontamination process in this end-time temporal boundary layer near to $t_D$ , in Appendix B, and find that the $O(\sqrt{\Gamma})$ correction to the decontamination time is $-\sqrt{\Gamma }$ . The improved decontamination time is therefore

(4.55) \begin{align} t_D&=t_D^*+\sqrt{\Gamma }\tau _D^*=\left (\frac{\Gamma }{2}+\sqrt{\Gamma }+1\right )-\sqrt{\Gamma }\nonumber \\ &=\frac{\Gamma }{2}+1.\end{align}

Surprisingly, this is the same as in the small- $\Gamma$ sublimit (4.30), as well as analogous to the sharp-interface model (3.8) (case i).

In summary, we have found that, in the sublimit $\Gamma \gg 1$ , the leading-order solution of (4.17) is, for $t\lt 1$ ,

(4.56a) \begin{align} c=\exp ({-}\sqrt{\Gamma }y), \qquad R=1-t\exp ({-}\sqrt{\Gamma }y), \end{align}

while for $1\lt t\lesssim \Gamma/2-\sqrt{\Gamma }$ ,

(4.56b) \begin{align} c&=\begin{cases} 1-\dfrac{\sqrt{\Gamma }y}{\sqrt{2t-1}} &\hspace{1.8cm}\text{ if }y\lt y_*=\dfrac{1}{\sqrt{\Gamma }}\left ({-}1+\sqrt{2t-1}\right ),\\[20pt] \dfrac{\exp ({-}\sqrt{\Gamma }y-1+\sqrt{2t-1})}{\sqrt{2t-1}} &\hspace{1.8cm}\text{ if }y\gt y_*=\dfrac{1}{\sqrt{\Gamma }}\left ({-}1+\sqrt{2t-1}\right ), \end{cases}\end{align}
(4.56c) \begin{align} R&=\begin{cases} 0 &\text{ if }y\lt y_*=\dfrac{1}{\sqrt{\Gamma }}\left ({-}1+\sqrt{2t-1}\right ),\\[20pt] 1-\exp ({-}\sqrt{\Gamma }y-1+\sqrt{2t-1}) &\text{ if }y\gt y_*=\dfrac{1}{\sqrt{\Gamma }}\left ({-}1+\sqrt{2t-1}\right ). \end{cases}\end{align}

For completeness, we also note that, from the analysis in Appendix B, the leading-order solution near the end-time, for $\Gamma/2-\sqrt{\Gamma }\lesssim t\leq t_D=\Gamma/2+1$ , as

(4.56d) \begin{align} c&=\begin{cases} 1-y+\dfrac{y}{\sqrt{\Gamma }}\left (\dfrac{\coth (Z_*)-Z_*}{1-Z_*/\sqrt{\Gamma }}\right ) &\text{ if }y\lt y_*=1-\dfrac{Z_*}{\sqrt{\Gamma }},\\[15pt] \dfrac{1}{\sqrt{\Gamma }}\dfrac{\left (\exp \left(\sqrt{\Gamma }(1-y)\right)+\exp ({-}\sqrt{\Gamma }(1-y))\right )}{e^{Z_*}-e^{-Z_*}} &\text{ if }y\gt y_*=1-\dfrac{Z_*}{\sqrt{\Gamma }}, \end{cases}\end{align}
(4.56e) \begin{align} R&=\begin{cases} 0 &\text{ if }y\lt y_*=1-\dfrac{Z_*}{\sqrt{\Gamma }},\\[15pt] 1-\dfrac{\exp \left(\sqrt{\Gamma }(1-y)\right)+\exp ({-}\sqrt{\Gamma }(1-y))}{e^{Z_*}+e^{-Z_*}} &\text{ if }y\gt y_*=1-\dfrac{Z_*}{\sqrt{\Gamma }}, \end{cases}\end{align}

where $Z_*(\tau )=Z_*((t-t_D^*)/\sqrt{\Gamma })$ is the solution of (B.9).

In Figure 11, we show solutions of (4.17) in the limit $\delta \rightarrow 0$ , for small $A=10^{-3}$ , and fairly large $\Gamma$ . We compare the numerical solution of (4.17) (with $\delta =0$ ) to the full asymptotic solution, (4.20) and (4.27), as well as to the reduced asymptotic solutions in the large- $\Gamma$ sublimit, given by (4.56). For brevity, here and throughout the remainder of the paper, we show $c$ and $R$ on the same axes. We note that, as in the small- $\Gamma$ case presented in Figure 9, the full asymptotic solution is in very good agreement with the numerical solution. The asymptotic solutions for the large- $\Gamma$ sublimit also capture the profiles $R$ and $c$ well but, in particular, are good approximations for the decontamination time, $t_D$ , even for $\Gamma =20$ , which is not particularly large. We note that we have chosen a very small value of $A$ ( $=10^{-3}$ ) here; this is not necessary to see good agreement between the asymptotic and numerical solutions, but by doing so we can more clearly see the errors in the asymptotic approximations that are due to the $O(\Gamma ^{-1/2})$ or $O(\Gamma ^{-1})$ correction terms.

Figure 11. Numerical and asymptotic solution of (4.6) in the case $A\ll 1$ , $\delta \ll 1$ and for $\Gamma =20$ . The numerical solutions (solid curves) are computed in the limit $\delta \ll 1$ , making the assumption that $D=V=F=1$ , and taking $A=10^{-3}$ . In (a), we show the decontamination front $y_*(t)$ , with the full, early-time and end-time asymptotic approximations of the solution of (4.26), (4.41) and (B.9), respectively. In (b), the solid lines are the numerical solution $c$ (red) and $R$ (blue) at times $0.01t_D$ , $0.1t_D$ , $t_D/4$ , $t_D/2$ , $3t_D/4$ and $0.95t_D$ , where $t_D$ is the time at which $y_*(t_D)=1$ , and the arrows indicate the direction of increasing time. Black dashed ( $R$ ) and dotted ( $c$ ) curves are the full analytical solutions (4.20) and (4.27), where $y_*(t)$ is the (numerically computed) solution of (4.26). The magenta dashed ( $R$ ) and dotted ( $c$ ) curves are the approximations (4.56) in the large- $\Gamma$ sublimit.

4.6.2. Thick agent layers: case $\delta =O(1)$

The asymptotic analysis presented in Section 4.6.1 is in the limit $\delta \ll 1$ , so that the layer of agent $R_0$ initially coating the walls of the pore is thin relative to the pore lengthscale $d$ . In that limit $\delta \ll 1$ , we saw that the equations (4.6a)–(4.6b) for $c$ and $R$ partially decouple, since variations in the thickness of the agent layer do not affect the effective cleanser diffusivity $D$ , or the geometric factors $F$ and $V$ at leading order.

Figure 12. Variation of solutions with $\delta$ in the limit $A\ll 1$ , for $\Gamma =0.01$ (top), $\Gamma =1$ (middle) and $\Gamma =100$ (bottom). On the left, we show the decontamination front $y_*(t)$ and on the right the cleanser $c$ (red) and agent-layer thickness $R$ (blue) profiles. Numerical solutions of (4.6) are shown as solid lines with $\delta =0.25$ , and in the limit $\delta \ll 1$ as dotted lines, for comparison. Throughout we take $A=0.01$ , and $\rho _0=0.2$ . Solutions in (b), (d) and (f) are shown at the times $0.1t_D$ , $0.2t_D$ , $0.4t_D$ , $0.6t_D$ , $0.8t_D$ and $t_D$ , where $t_D$ is the decontamination time for the $\delta =0.25$ solution.

We now investigate the effect of non-negligible $\delta$ in the limit of dilute cleanser $A\ll 1$ . In Figure 12, we show solutions of (4.6) for $\delta =0.25$ (solid lines) compared with solutions in the limit $\delta \ll 1$ (dotted lines), for three different values of $\Gamma$ . We note that due to the geometrical constraints we require $\delta \in (0,0.5-\rho _0)=(0,0.3)$ . We see qualitatively similar behaviours of both the decontamination front $y_*(t)$ and the profiles of $c$ and $R$ to those observed in the small- $\delta$ limit, in all of the cases $\Gamma =0.01$ , $\Gamma =1$ and $\Gamma =100$ . This is unsurprising, as, although the local values of $D$ , $F$ and $V$ vary with $\delta$ , the asymptotic balances are not expected to change until $D$ becomes very small, which occurs when $\delta \rightarrow 0.5-\rho _0$ . However, we notice that, since $D$ now decreases with $R$ , cleanser cannot diffuse so easily into regions with large amounts of agent, and so we observe lower $c$ values in regions with large $R$ for the $\delta =0.25$ solutions, compared with the $\delta \ll 1$ solutions at the same time-points. Since $c$ is reduced, the rate of reaction is slower in these regions, and so the agent-depletion rate is reduced, leading to steeper spatial gradients in the $R$ -profiles. These effects are particularly evident when $\Gamma$ is large, since $R$ exhibits greater spatial variation than for $\Gamma \ll 1$ .

We note that we have kept $\Gamma$ constant in each of the cases shown in Figure 12, in order to observe the effect of the spatial variation in $D$ , $V$ and $F$ . However, since

(4.57) \begin{align} \Gamma =\frac{vk}{\epsilon f_1D_0\pi }\cdot \frac{1}{\delta (2\rho _0+\delta )}, \end{align}

in order to hold the total initial amount of agent, $v$ , constant, we should vary $\Gamma$ simultaneously with $\delta$ . In this case, the qualitative changes as $\delta$ is varied are mainly due to the corresponding change in $\Gamma$ and are similar to those observed in the $\delta \ll 1$ limit previously studied. We therefore do not show these solutions here.

4.6.3. Comparison with the sharp-interface model

When $A\ll 1$ (case I), we have seen that the decontamination process is limited by the reaction, with quasi-steady diffusion of cleanser on the reaction timescale. This is identical to the case-i sharp-interface model when $\alpha \ll 1$ . For both models, this is because the reaction rate is proportional to the cleanser concentration, so that for small $A$ or $\alpha$ (corresponding to a weak applied cleanser concentration compared with the local amount of agent), the reaction is slow relative to the diffusive transport of cleanser.

When $\Gamma \ll 1$ and $A\ll 1$ (case Ia) the motion of the decontamination front $y_*(t)$ is seen to be very different from the $\sqrt{t}$ behaviour of the reaction front in the sharp-interface model (case ia), since the decontamination reaction occurs almost uniformly throughout $y\in [0,1]$ . Nevertheless, we found an analogous leading-order expression for the decontamination time, $t_D=1+\Gamma/2$ , compared with $t_D=1+\gamma/2$ for the sharp-interface case. It is surprising that these expressions are so directly analogous: it is not clear to us precisely what the underlying reason for this is. Of course, we have used different time scalings for the two models, and so these do not correspond to equal dimensional decontamination times.

When $\Gamma \gg 1$ and $A\ll 1$ (case Ib), we recovered the same $\sqrt{t}$ behaviour of the leading-order decontamination front, $y_*(t)$ as seen for the motion of the reaction front $h(t)$ in the sharp-interface model in case i. The approximate decontamination time in case Ib was also found to be $\Gamma/2+1$ , again analogous to the sharp-interface case ib. Indeed, the agent-on-walls model in this limit exhibited behaviour like a sharp interface, with only a narrow $O(\Gamma ^{-1/2})$ reaction region over which $R$ varied, and we showed that the leading-order solution satisfies a sharp-interface model (4.53).

4.7. The case $A\gg 1$ : diffusion-limited decontamination

As discussed previously, we expect the dilution of the cleanser to be $\alpha = O(1)$ or smaller for physically relevant decontamination scenarios, but since $A=\alpha f_1/\delta$ we may have $A\gg 1$ if $\delta \ll \alpha$ , that is, if the layer of agent coating the walls of the porous medium is sufficiently thin. In this section, we suppose $A\gg 1$ and, since $\delta$ is therefore small, the model (4.6) reduces to

(4.58a) \begin{align} A\Gamma c_t&=c_{yy}-\Gamma H(R)\left (1+\frac{\delta A}{f_1}c\right )c,\\[-10pt]\nonumber \end{align}
(4.58b) \begin{align} R_t&=-H(R)c, \end{align}

with the initial and boundary conditions given by (4.6c)–(4.6e).

Since $A\gg 1$ , we immediately notice that the time derivative of $c$ in (4.58a) dominates over the reaction term. Thus, we expect that the system is diffusion-limited when $A\gg 1$ (compared to reaction-limited as we found for the case $A\ll 1$ ). We consider the behaviour in subcases depending on the size of $\Gamma$ , from small $\Gamma =O( A^{-1})\ll 1$ to large $\Gamma \gg 1$ .

Figure 13. Asymptotic behaviour of the decontamination front $y_*$ and decontamination time $t_D$ in case II. In (a), we show the solution $\sqrt{A\Gamma }y_*$ of (4.60c) against $t$ (black) and the large- $t$ (large- $A\Gamma$ ) limit (4.62) (red dashed). In (b), we show the solution $t_D$ of (4.61) (blue) and the large $A\Gamma$ limit (4.63) (red dotted).

4.7.1. Case II: shallow spills, $\Gamma = O(A^{-1})$

We first suppose the spill depth is shallow, taking the distinguished limit $\Gamma =O (A^{-1})\ll 1$ .

With $A\Gamma =O(1)$ , we see that, at leading order in $A^{-1}$ , (4.58a) becomes

(4.59) \begin{align} A\Gamma c_t=c_{yy} \end{align}

regardless of whether $R$ is zero or positive. The cleanser and agent equations therefore fully decouple. Solving (4.59) subject to (4.6c)–(4.6e), the leading-order cleanser concentration $c$ is given by the similarity solution

(4.60a) \begin{align} c=\text{erfc}\left (\frac{\sqrt{A\Gamma }y}{2\sqrt{t}}\right ). \end{align}

Substituting (4.60a) into (4.58b) and solving the resulting equation, we find that $R$ takes the form

(4.60b) \begin{align} R=1-\text{erfc}\left (\frac{\sqrt{A\Gamma }y}{2\sqrt{t}}\right )\left (\frac{A\Gamma y^2}{2}+t\right )+\frac{y\sqrt{A\Gamma t}}{\sqrt{\pi }}\exp \left ({-}\frac{A\Gamma y^2}{4t}\right ), \end{align}

for $y\in [y_*(t),1]$ (and $R=0$ otherwise), with $y_*(t)$ defined as the point at which $R$ is first zero, namely the solution of

(4.60c) \begin{align} 1+\frac{y_*\sqrt{A\Gamma t}}{\sqrt{\pi }}\exp \left ({-}\frac{A\Gamma y_*^2}{4t}\right )=\text{erfc}\left (\frac{\sqrt{A\Gamma }y_*}{2\sqrt{t}}\right )\left (\frac{A\Gamma y_*^2}{2}+t\right ). \end{align}

We show how $\sqrt{A\Gamma } y_*$ , given by the solution of (4.60c), varies with $t$ in Figure 13a, and note that (4.60c) has a solution $y_*\gt 0$ only for $t\gt 1$ . In particular, the decontamination time $t_D$ (at which $y_*=1$ ) is the solution of

(4.61) \begin{align} 1+\sqrt{\frac{A\Gamma t_D}{\pi }}\exp \left ({-}\frac{A\Gamma }{4t_D}\right )=\text{erfc}\left (\frac{\sqrt{A\Gamma }}{2\sqrt{t_D}}\right )\left (\frac{A\Gamma }{2}+t_D\right ). \end{align}

We show how $t_D$ , given by the solution of (4.61), varies with $A\Gamma$ in Figure 13b. This has error of order $\max (A^{-1},\delta )$ (since it is the reaction term of (4.58a), and the $O(\delta )$ corrections to $D,F,$ and $V$ that we have neglected. This solution in the distinguished limit $A\Gamma =O(1)$ is in fact valid for all $\Gamma \ll 1$ , since the loss of cleanser to the chemical reaction only affects the leading-order solution of (4.58a) when $\Gamma =O(1)$ or larger.

For very shallow spills $\Gamma \ll A^{-1}$ , the solutions $c$ and $R$ are approximately uniform in $y$ , and indeed from (4.61) we see that $t_D\rightarrow 1$ as $A\Gamma \rightarrow 0$ . For deeper spills and so larger $\Gamma$ (with $A^{-1}\ll \Gamma \ll 1$ ), we see from (4.58a) that the decontamination process is slower, since the diffusive timescale is longer. Indeed, for $A\Gamma \gg 1$ (with $A^{-1}\ll \Gamma \ll 1$ ) we find that the solution $y_*$ of (4.60c) behaves like

(4.62) \begin{align} y_*\sim \sqrt{\frac{6t}{A\Gamma }W_0\left (\frac{2t^{2/3}}{3\pi ^{1/3}}\right )}, \end{align}

where $W_0$ is the principal branch of the Lambert-W function, and from (4.61) we see that

(4.63) \begin{align} t_D\sim \frac{A\Gamma }{10 W_0\left (\frac{2}{5}\left (\frac{A\Gamma }{4\sqrt{\pi }}\right )^{2/5}\right )}. \end{align}

In Figure 14 (left), we show numerical solutions of (4.58) and the asymptotic approximations (4.60), for $\Gamma =0.1\in (A^{-1},1)$ , and we see very good agreement even for this moderate value of $\Gamma$ . Unlike in the case $A\ll 1$ studied in Section 4.6, the cleanser profile is independent of the agent thickness, with no change in behaviour between $y\lt 1$ and $y\gt 1$ .

Figure 14. Numerical and asymptotic solution of (4.6) in the case $A\gg 1$ (and $\delta =O(A^{-1})\ll 1$ ). We take $A=100$ throughout with $\Gamma =0.1$ (left) and $\Gamma =100$ (right). We show the $c$ and $R$ profiles against the spatial coordinate $y$ in (a)–(b) and against the similarity variable $\sqrt{A\Gamma } y/2\sqrt{t}$ in (c)–(d). The decontamination front $y_*(t)$ is shown in (e)–(f). The numerical solutions (solid curves) are computed in the limit $\delta \ll 1$ , making the assumption that $D=V=F=1$ . In (a)–(d), the solid curves are the numerical solution of (4.58) ( $c$ red, $R$ blue) and the black curves are the analytical solution (4.60a) for $c$ (dotted) and (4.60b) for $R$ (dashed), where $y_*(t)$ is the (numerically computed) solution of (4.60c), shown in (e)–(f), along with the numerical decontamination front.

4.7.2. Case III: deep spills, $\Gamma =O(1)$ and larger

Alternatively if the agent depth is so deep that $\Gamma =O(1)$ or larger, we see from (4.58a) that the consumption of cleanser in the reaction will slow down the decontamination process, and so the solution $t_D$ given by (4.61) will underestimate the decontamination time. At early, $O(1)$ times, the same solution (4.60a)–(4.60c) holds, with $c$ and $R$ only varying within an $O(1/\sqrt{A\Gamma })$ boundary layer at $y=0$ .

In Figure 14 (right), we show solutions for $\Gamma =100$ . We see that (4.60) still provides a good approximation at early time, and the $c$ profile is still close to the similarity solution (4.60a) for all time (although for $y\gt y_*$ we observe $c$ decaying more quickly with $y$ , since cleanser is consumed appreciably in the reaction). We see that the cleanser concentration is very low for regions where $R\gt 0$ , and $R$ varies rapidly for $y\gt y_*$ , because the system is diffusion-limited. The non-negligible depletion of $c$ has a major impact on the decontamination behaviour, with the solution $R$ diverging significantly from the analytic expression (4.60b), and the shape of the decontamination front $y_*(t)$ is very different to that given by (4.60c).

4.8. Moderate cleanser concentration, $A=O(1)$

Having studied the model (4.6) in the limits of $A\ll 1$ and $A\gg 1$ we now investigate the case $A=O(1)$ . Since $A=\alpha f_1/\delta$ , the case $A=O(1)$ requires that $\alpha/\delta =O(1)$ , so that the local amounts of cleanser and agent are approximately balanced. Specifically, this is true for both $\alpha,\delta$ small (so that the cleanser is dilute, but also there is a small volume fraction of agent within the pore space), or equally for both $\delta$ and $\alpha$ of order 1. We show numerical solutions of (4.6) with $A=1$ in Figure 15, for various sizes of $\Gamma$ , and also for both $\delta \ll 1$ and $\delta =0.25$ . With $A=O(1)$ , we see that both the diffusion of cleanser and the reaction rate play a role, and the system is neither diffusion-limited nor reaction-limited. As in the case $A\ll 1$ , we see qualitative agreement between the cases $\delta \ll 1$ and $\delta =O(1)$ , although the non-uniform cleanser diffusivity $D(R)$ for $\delta =O(1)$ results in slower transport of $c$ through regions where $R$ is close to 1, and correspondingly steeper spatial profiles of $R$ , as well as a clear discontinuity in $c_y$ at $y=1$ .

Figure 15. Variation of solutions with $\delta$ , holding $A=1$ , for $\Gamma =0.01$ (top), $\Gamma =1$ (middle) and $\Gamma =100$ (bottom). On the left, we show the decontamination front $y_*(t)$ and on the right the cleanser $c$ (red) and agent-layer thickness $R$ (blue) profiles. Numerical solutions of (4.6) are shown as solid lines with $\delta =0.25$ , and in the limit $\delta \ll 1$ (making the assumption $D=V=F=1$ ) as dotted lines, for comparison. Throughout we take $\rho _0=0.2$ . Solutions in (b), (d) and (f) are shown at the times $0.1t_D$ , $0.2t_D$ , $0.4t_D$ , $0.6t_D$ , $0.8t_D$ and $t_D$ , where $t_D$ is the decontamination time for the $\delta =0.25$ solution. The solutions for $\delta =0.25$ show clear discontinuities in $c_y$ at $y=1$ due to the discontinuous diffusivity $D(R)$ .

4.9. Decontamination times and decontamination efficiency

As in the sharp-interface case, the decontamination time and the decontamination efficiency are of practical interest. The decontamination time $t_D$ is the time when all the agent has been removed and occurs when $y_*(t_D)=1$ . We have found the analytic expressions (4.30) and (4.54) for $t_D$ in the limits $\Gamma \ll 1$ and $\Gamma \gg 1$ , respectively, in both cases assuming that $A\ll 1$ , $\delta \ll 1$ . The efficiency is defined (analogously to the sharp-interface case) as the ratio of the minimum amount of cleanser required to react away all the agent, to the total amount of cleanser that has entered the porous medium, namely

(4.64) \begin{align} \text{efficiency}=\frac{N}{\text{total flux in}}=\left (\frac{2\rho _0+\delta }{2\rho _0}\right )\frac{\Gamma }{-\int _{t=0}^{t_D}\left (DVc_y\right )|_{y=0}\,\mathrm{d}t}\,, \end{align}

using the expression (3.23) for $N$ , and noting that

(4.65) \begin{align} \text{total flux in}= - \frac{R_0D_0V_0}{\chi k L_{\text{AW}}}\int _{t=0}^{t_D}\left (DVc_y\right )|_{y=0}\,\mathrm{d}t \quad \text{mol}\,\text{m}^{-2}\,. \end{align}

In the limits in which we have approximate analytical solutions, we compute the approximate efficiency, finding that

(4.66) \begin{align} \text{efficiency}\sim \begin{cases}1 & \text{ if }A\ll 1, \delta \ll 1, \Gamma \ll 1,\text{ case Ia},\\[4pt] 1 & \text{ if }A\ll 1, \delta \ll 1, \Gamma \gg 1,\text{ case Ib},\\[11pt] \dfrac{2\sqrt{\pi \Gamma }}{\sqrt{A}} & \text{ if }A\gg 1, \delta \ll 1, \Gamma \ll 1,\text{case II}. \end{cases} \end{align}

In Figure 16, we explore the effect of the model parameters on the decontamination time. Broadly speaking, we observe similar dependence of the decontamination time on $A$ and $\Gamma$ as we saw for $\alpha$ and $\gamma$ , respectively, in the sharp-interface case discussed in Section 3.6. In Figure 16a, we show the variation of the numerically computed decontamination times as a function of $\Gamma$ , which we may interpret as the effect of the spill depth $L_{AW}$ on the decontamination time. We recall that for $A\ll 1$ (case I), for both $\Gamma \ll 1$ (case Ia) and $\Gamma \gg 1$ (case Ib) we have found the same expression for the decontamination time, namely $t_D=1+\Gamma/2$ . We see that, for small values of $A$ , this approximation $\Gamma/2+1$ is a good estimate for the decontamination time, even when $\Gamma =O(1)$ . For moderate values of $A=O(1)$ , we see that, while the approximation remains reasonable for $\Gamma \ll 1$ , the decontamination time is greater than $\Gamma/2+1$ for intermediate and large values of $\Gamma$ . For large values of $A$ , the small- $\Gamma$ approximation (case II) to $t_D$ (the solution of (4.61)) shows good agreement while $\Gamma$ is small. For $A\gt 1$ , the growth of $t_D$ with the spill depth appears to be roughly linear in $\Gamma$ , although we note that the scaling $A\Gamma$ suggested in Section 4.7.2 for large $A$ and $\Gamma$ gives an upper bound.

Figure 16. The behaviour of the decontamination time $t_D$ (in (a)), or the scaled decontamination time $t_D/A\Gamma$ (in (b)–(c)) or (4.70) (in (d)) for the agent-on-walls model, computed from the numerical solution of (4.6). In (a)–(c), we take the small- $\delta$ limit and investigate the effect of spill depth and reaction rate on the decontamination time via $\Gamma$ , and of the applied cleanser concentration via $A$ . In (d), we show the effect of the initial agent-layer thickness via $\delta$ .

As in the sharp-interface case, we note that the time scaling used in (4.5) depends on both the reaction rate $k$ and the applied cleanser concentration $c_0$ . In order to draw physical conclusions about the impact of these parameters on the decontamination time, in Figure 16b and c we therefore plot the scaled decontamination time

(4.67) \begin{align} \frac{t_D}{A\Gamma }=\frac{R_0D_0}{\chi k c_0f_2L_{\text{AW}}^2}t_D, \end{align}

which has the same dependence on $c_0$ and $k$ as the time scaling used in (4.5). In Figure 16b, we see that the decontamination time decreases with reaction rate $k$ , so that faster reaction rates result in faster decontamination. As in the sharp-interface case, we observe limiting returns. With $A\leq O(1)$ we see that for very large $k$ (and so large $\Gamma$ ), the decontamination time approaches a constant value, $1/2A$ . With $A\gt 1$ we also see much slower variation of $t_D$ with $\Gamma$ as $\Gamma$ becomes large. In Figure 16c, we show the variation of $t_D$ with $A$ , which is proportional to the applied cleanser concentration $c_0$ , and so we may interpret increasing $A$ with increasing the applied cleanser concentration. For stronger cleanser concentrations, we observe faster decontamination. The large- $A$ asymptotic solution given by (4.61) appears to be a good approximation at large $A$ values, even for large values of $\Gamma$ .

We now consider the effect of $\delta$ on the decontamination time. We note that for a fixed total volume of agent $v$ , as $\delta$ increases the depth $L_{\text{AW}}$ of the spill is decreased, since, from the definition (4.2) of $L_{\text{AW}}$ ,

(4.68) \begin{align} L_{\text{AW}}=\frac{v}{\pi \delta (2\rho _0+\delta )}. \end{align}

Accordingly, $\Gamma$ decreases as $\delta$ increases, since

(4.69) \begin{align} \Gamma =\frac{vk}{\pi \epsilon f_1D_0}\cdot \frac{1}{\delta (2\rho _0+\delta )}. \end{align}

In order to understand how the agent thickness impacts on the physical decontamination time while $v$ is held constant, we therefore compute solutions holding $\Gamma \delta (2\rho _0+\delta )$ constant. In order to understand the effect of $\delta$ on the dimensional decontamination time, we define another scaled decontamination time, namely

(4.70) \begin{align} \frac{t_D}{A\Gamma \delta ^2(2\rho _0+\delta )^3}=\frac{D_0 R_0\pi ^2 g}{2\rho _0\chi k c_0}t_D, \end{align}

which has the same dependence on $R_0$ , and so $\delta$ , as the time scaling used in (4.5). We plot this second scaled decontamination time (4.70) against $\delta$ in Figure 16d. The dominant effect of increasing $\delta$ is to significantly reduce the spill depth $L_{\text{AW}}$ , and thus we see that for larger $\delta$ and thicker agent layers on the pore walls, the scaled decontamination time decreases. We therefore conclude that the dimensional decontamination time decreases with the agent layer thickness $\delta$ . This behaviour is, of course, limited by the requirement that the cleanser-occupied region of the porespace remains connected, which is true so long as $\delta \lt 0.5-\rho _0$ .

In Figure 17, we explore the effect of the system parameters on the efficiency (4.64). Like the decontamination time, the efficiency also behaves mostly similarly with $\Gamma$ and $A$ as the sharp-interface efficiency behaved with $\gamma$ and $\alpha$ . In Figure 17a, we show the variation of the efficiency (4.64) with $\Gamma$ , observing that, for all values of $A$ , the efficiency increases with $\Gamma$ . Thus, for faster-reacting cleansers (larger values of $k$ and so $\Gamma$ ), less cleanser is wasted. Unlike in the sharp-interface case, for small $\Gamma$ we see the efficiency decreasing to zero, rather than approaching a constant as we saw in Figure 6d. This is because, in this agent-on-walls case, cleanser can diffuse through $y=1$ and be lost from the system, unlike in the sharp-interface case where cleanser was confined to $y\in [0,h(t)]$ . From (4.15), we find that the flux of cleanser through $y=1$ scales with $\sqrt{A\Gamma/t}$ . For small $\Gamma$ (when $t_D\rightarrow 1$ ), we therefore expect the cleanser lost through $y=1$ to be of order $\sqrt{A\Gamma }$ .

Figure 17. The behaviour of the decontamination efficiency (4.64) for the agent-on-walls model, computed from the numerical solution of (4.6). In (a)–(c), we take the small- $\delta$ limit. In (a), (b) and (d), we show the effect of $\Gamma$ , $A$ and $\delta$ , respectively, on the efficiency. In (c), we show the variation of the efficiency with the decontamination time, parameterised by $A$ .

From Figure 17b, we observe that, as in the sharp-interface case, the efficiency decreases with $A$ and so, for stronger applied cleanser concentrations (larger $c_0$ and so larger $A$ ), the efficiency worsens. Thus, as in the sharp-interface case, there is a pay-off between improving the decontamination time (by increasing the applied cleanser concentration $c_0$ ) and improving the efficiency (by decreasing $c_0$ ). In Figure 17c, we illustrate this pay-off by plotting the decontamination time against the efficiency, and we observe very similar qualitative behaviour to the sharp-interface case (in Figure 6f). In Figure 17d, we explore the effect of $\delta$ on the efficiency. As for Figure 16d, we hold the total initial volume of agent fixed, by holding $\Gamma \delta (2\rho _0+\delta )$ constant. We see that, as $\delta$ increases and $\Gamma$ is reduced, the efficiency drops.

5. Discussion and conclusions

In this paper, we investigated the effect of the initial local distribution of agent within the pores of a porous material on the reactive decontamination of the porous material. Our macroscale models for two canonical agent distributions (sharp-interface and agent-on-walls) have very different structures, with the reaction in the sharp-interface case occurring on the macroscale domain boundary, while the average effect of the reaction for the agent-on-walls case appears as a bulk macroscale reaction term. We solved both models numerically and investigated relevant parameter regimes asymptotically, as summarised in Figures 2 and 8 for the sharp-interface and agent-on-walls cases, respectively. Key dimensionless parameters were found to be the dilution of the cleanser compared with the local amount of agent ( $\alpha$ for the sharp-interface and $A$ agent-on-walls cases), and the depth of the spill or speed of the chemical reaction, captured by $\gamma$ and $\Gamma$ in the sharp-interface and agent-on-walls models, respectively.

Despite the fundamental differences between the two models, some aspects of the decontamination were qualitatively similar. In particular, in the specific case (case Ib, shown in Figure 8) of a dilute applied cleanser we saw that for sufficiently deep agent spills or sufficiently fast reaction rates, the reaction zone in the agent-on-walls model becomes very narrow. The agent-on-walls model reduces to a Stefan problem in this limit, and thus behaves precisely like a sharp-interface model (case ib, as in Figure 2).

Assuming the same initial total volume of agent for both models, we directly compared the decontamination between the two scenarios. We found that, for the same total initial amount of agent, the agent distribution in the pores has a significant effect on the timescale of decontamination. If the agent spill is sufficiently shallow or the reaction rate is sufficiently slow so that $\gamma,\Gamma \ll 1$ , we saw that the decontamination of both models is reaction-limited, and the sharp-interface decontamination is slower than the agent-on-walls case ( $O(\delta/\epsilon )=O(R_0L_{\text{AW}}/d^2)$ slower), because the surface area on which the reaction can occur is so much greater in the agent-on-walls case. Conversely, for sufficiently deep spills or fast reaction rates so that $\gamma,\Gamma \gg 1$ , the fact that the lengthscale of the agent spill for the agent-on-walls model is so much deeper than in the sharp-interface case becomes important. For instance, in the dilute limit $\alpha, A\ll 1$ for which we have analytical expressions for the decontamination time, in the case $1\ll \gamma \ll \Gamma$ , the sharp-interface decontamination is significantly faster than the agent-on-walls decontamination ( $O(\epsilon \delta )$ faster, equal to the ratio of the agent-layer thickness on the pore walls to the depth of the overall spill). Although the dimensional decontamination times vary significantly between the two models, we found analogous dimensionless expressions for the sharp-interface and agent-on-walls model in the dilute-cleanser limit, namely $t_D=1+\gamma/2$ and $t_D=1+\Gamma/2$ , respectively. This hints at some deeper similarity between these problems, which would be interesting to explore further.

By investigating the decontamination time and decontamination efficiency, we saw that, for both the agent-on-walls and sharp-interface scenarios, by choosing a cleanser with a faster reaction rate we can both reduce the decontamination time and the amount of cleanser that is wasted. We also found that increasing the strength of the applied cleanser solution speeds up the decontamination but also resulted in a less efficient decontamination, with more cleanser chemical left unreacted in the environment. In practice, an appropriate cleanser concentration should be chosen to balance these effects as best as possible.

Throughout this paper, we have assumed a two-dimensional, circular pore geometry. Different geometries, and in particular a three-dimensional pore structure, would alter the geometric factors $f_1$ , $f_2$ and $F$ , as well as the effective diffusivity $D$ , volume fraction $V$ and porosity $\phi$ . Although the solutions would change quantitatively, the qualitative behaviour of the models is not expected to change with different pore-scale geometry.

In our decontamination models, we have assumed that the concentration of cleanser applied at the surface of the porous material is constant. An interesting avenue for future work might be to prescribe a time-varying cleanser concentration at the surface of the porous medium and investigate how to both minimise the decontamination time and maximise the decontamination efficiency. Our analysis suggests that applying a strong cleanser initially, and then a progressively weaker concentration, might be beneficial. The transport of cleanser in both cases is assumed to be purely diffusive. If the agent and cleanser solution had different densities then, as the agent material is reacted away to form the reaction product, a flow would be generated and we would expect both the diffusion and advection of chemicals with such a flow to be important. For the sharp-interface case, such a flow in the cleanser phase might be incorporated into the homogenised model using a similar analysis to that developed for an evaporation front in [Reference Luckins, Breward, Griffiths and Please11]. The effect of such a flow is investigated in [Reference Geng, Kamilova and Luckins7] for the sharp-interface case, but it is not clear what the effect would be on solutions of the agent-on-walls model. Finally, we have assumed two highly idealised agent distributions in the pore space. It is not clear how the agent might be distributed in practice; this may depend on the porous medium and the fluid properties of the agent, and we might expect spatial variation in the initial agent distribution. Investigation of the contamination problem, or how the agent infiltrates the porous medium initially, is an important area for further work.

Data availability statement

In compliance with the University of Oxford’s open access initiative, the data in this paper are available from http://dx.doi.org/10.5287/ora-nyvkkkxdn.

Competing interest

None.

Acknowledgements

EKL is grateful for funding from the Industrial Fund of the Centre For Doctoral Training in Industrially Focused Mathematical Modelling. IMG is grateful to the Royal Society for funding through a University Research Fellowship. We also thank U. Beuscher and V. Venkateshwaran (W.L. Gore & Associates, Inc.), L. Hetherington (Defra), G. Anderson (Beko) and S. Tarkuç (Arçelik) for their advice and guidance.

Appendix A: The cell problem for the effective diffusivity $\mathcal{D}$

As shown in [Reference Luckins, Breward, Griffiths and Wilmott10], the effective diffusivity for the agent-on-walls problem $\mathcal{D}(R)$ may be found by the solution of a cell problem on the microscale of the porous medium. Specifically, $\mathcal{D}$ is given by

(A.1) \begin{align} \mathcal{D}=1+\frac{1}{\mathcal{V}(R)}\iint _{\omega (R)}W_X\,\mathrm{d}X\mathrm{d}Y, \end{align}

where $W$ satisfies the cell problem

(A.2) \begin{align} \nabla ^2 W=0 \qquad &\text{ for }X,Y\in \omega (R),\\[-10pt]\nonumber \end{align}
(A.3) \begin{align} \left (\nabla W+\boldsymbol{e}_1\right )\cdot \boldsymbol{n}=0 \qquad &\text{ on }|\boldsymbol{X}|=r_0+R,\\[-10pt]\nonumber \end{align}
(A.4) \begin{align} W\text{ periodic in both }X\text{ and }Y \qquad &\text{ over }\omega (R), \end{align}

where $\boldsymbol{n}$ is the unit normal to the agent–cleanser interface at $|\boldsymbol{X}|=r_0+R$ , and $\omega (R)$ is the cleanser-occupied region, namely $\omega (R)=[-d/2,d/2]^2\setminus \{|\boldsymbol{X}|\lt r_0+R\})$ . This definition (A.1) of the effective diffusivity $\mathcal{D}$ is common to many diffusion problems in multi-scale geometries. The solution $\mathcal{D}$ has been computed for our circular geometry in, for instance, [Reference Bruna and Chapman1, Reference Dalwadi, Griffiths and Bruna3], and we will use their solution in this paper.

Appendix B: The end-time behaviour of the agent-on-walls model in case Ib

In this appendix, we study the behaviour of the agent-on-walls model in the case where $\delta,A\ll 1$ , and $\Gamma \gg 1$ (that is case Ib), at times for which $y_*(t)$ is within $O(1/\sqrt{\Gamma })$ of the end of the contaminated region located at $y=1$ . This allows us to compute the $O(\sqrt{\Gamma})$ correction to the decontamination time in this limit.

In the region $y\lt y_*$ in which $R=0$ , we have the same linear profile for $c$ as found in the late-time studied in Section 4.6.1, so that, to leading order in $1/\sqrt{\Gamma }$ , in $y\lt y_*$ ,

(B.1) \begin{align} c=1-y. \end{align}

We note that we have fixed $c=0$ at $y=1$ for this outer solution, since, as in the late-time behaviour, we expect $c=O(1/\sqrt{\Gamma })$ for $y\gt y_*$ , and we have $y_*=1+O(1/\sqrt{\Gamma })$ . For $y\gt y_*$ , we make the change of variables $y=1-Z/\sqrt{\Gamma }$ and define the decontamination front to be located at $Z=Z_*$ , where $y_*=1-Z_*/\sqrt{\Gamma }$ . We also use the same change of variables $c=C/\sqrt{\Gamma }$ as in the boundary layer located at $y=y_*$ in the late-time case. By balancing the terms in (4.17b), we see that the timescale of the motion of the reaction front is $O(\sqrt{\Gamma})$ in this end-zone. We therefore make the change of time-variables $t=t_D^*+\sqrt{\Gamma }\tau$ , where $\tau =O(1)$ and $t_D^*$ is our prediction (4.54) for the decontamination time before taking end effects into account. (We note that this timescale $O(\sqrt{\Gamma})$ is different from the $O(\Gamma )$ timescale in the late-time analysis presented in Section 4.6.1 because we are interested in the motion of the decontamination front over a shorter, $O(1/\sqrt{\Gamma })$ , distance.) With these changes of variable, equations (4.17) become

(B.2a) \begin{align} \frac{A}{\sqrt{\Gamma }}C_\tau &=C_{ZZ}-H(R)C+O(\delta/\Gamma ),\end{align}
(B.2b) \begin{align} R_\tau &=-H(R)C.\end{align}

At $Z=Z_*$ we impose continuity of the cleanser flux (matching with the outer solution (B.1)) to find that

(B.2c) \begin{align} C_Z=1 \qquad \text{ at }\qquad Z=Z_*(\tau ). \end{align}

At the end of the agent-occupied domain, $Z=0$ (equivalently $y=1$ ), we see that the flux through this point, given by (4.15), is

(B.2d) \begin{align} C_Z|_{Z=0}=O\left (A^{1/2}\Gamma ^{-1/4}\right )\ll 1. \end{align}

Initial conditions are given by matching with the late-time solution as $\tau \rightarrow \infty$ , and we will give details of this later as needed.

From (B.2), the leading-order cleanser concentration $C$ therefore satisfies, for $0\lt Z\lt Z_*(\tau )$ ,

(B.3a) \begin{align} C_{ZZ}=C, \end{align}

with

(B.3b) \begin{align} C_Z=0 \qquad \text{ at }Z=0, \qquad \text{ and } \qquad C_Z=1 \qquad \text{ at }Z=Z_*, \end{align}

with solution

(B.4) \begin{align} C=\frac{\cosh (Z)}{\sinh (Z_*)} . \end{align}

(As an aside, we note that we may easily compute a correction term, linear in $y$ to the outer solution (B.1) to ensure that the solution is continuous at $y=y_*$ to $O(1/\sqrt{\Gamma })$ .) From (B.4), we see that the leading-order equation for $R$ , (B.2b), is

(B.5) \begin{align} R_\tau =-\frac{\cosh (Z)}{\sinh (Z_*)} . \end{align}

This must be solved for along with the matching requirements that $R\rightarrow 1$ and $Z_*\rightarrow \infty$ as $\tau \rightarrow -\infty$ , and the requirement that $R=0$ at $Z=Z_*(\tau )$ . We look for a solution of the form

(B.6) \begin{align} R=1-r(\tau )\cosh (Z). \end{align}

This has the correct behaviour as $\tau \rightarrow -\infty$ , and the constraint at $Z=Z_*$ requires that

(B.7) \begin{align} r(\tau )=\frac{1}{\cosh (Z_*)}. \end{align}

Substituting into (B.5), we obtain the ODE

(B.8) \begin{align} Z^{\prime}_*(\tau )=-\frac{\cosh ^2(Z_*)}{\sinh ^2(Z_*)}, \end{align}

with implicit solution

(B.9) \begin{align} Z_*-\tanh (Z_*)+a=-\tau, \end{align}

where $a$ is an integration constant, which must be fixed by more careful matching with the late-time solution as $\tau \rightarrow -\infty$ . We note from (B.9) that as $\tau \rightarrow -\infty$ , $Z_*\sim -(\tau +a-1)$ , and thus the decontamination front is at

(B.10) \begin{align} y_*\sim 1+\frac{\tau +a-1}{\sqrt{\Gamma }}+\text{exponentially small terms} \end{align}

as $\tau \rightarrow -\infty$ . Conversely, the late-time solution

(B.11) \begin{align} y_*(t)=\frac{1}{\sqrt{\Gamma }}\left ({-}1+\sqrt{2t-1}\right ), \end{align}

in terms of the $\tau$ variable, behaves like

(B.12) \begin{align} y_*\sim 1+\frac{\tau }{\sqrt{\Gamma }}+O\left (\frac{1}{\Gamma }\right ) \end{align}

as $\tau \rightarrow -\infty$ given by (4.54). Comparing the two we conclude that we must set $a=1$ for the expansions to match.

The solution for $Z_*(\tau )$ may be computed by solving (B.9) numerically with $a=1$ . The leading-order solutions $C$ and $R$ are then given by (B.4) and (B.6), respectively. In particular, we note from (B.9) (with $a=1$ ) that $Z_*=0$ when $\tau =\tau _D=-1$ . Thus, we have found the correction to the predicted decontamination time $t_D^*$ , due to the end effects, namely $\sqrt{\Gamma }\tau _D=-\sqrt{\Gamma }$ .

We note that we have implicitly assumed that $A\ll \Gamma ^{-1/2}$ in order to match between the late and end times. Were instead $A=O(\Gamma ^{-1/2})$ then we expect an additional contribution at $O(\Gamma ^{-1/2})$ in the late-time behaviour, which would affect the matching and thus the estimate of the decontamination time.

References

Bruna, M. & Chapman, S. J. (2015) Diffusion in spatially varying porous media. SIAM J. Appl. Math. 75(4), 16481674.CrossRefGoogle Scholar
Dalwadi, M., Dubrovina, E., Eisenträger, A., et al. (2021) Toxic chemicals and their neutralising agents in porous media. ESGI100. http://www.maths-in-industry.org/miis/671/ CrossRefGoogle Scholar
Dalwadi, M. P., Griffiths, I. M. & Bruna, M. (2015) Understanding how porosity gradients can make a better filter using homogenization theory. Proc. R. Soc. A Math. Phys. Eng. Sci. 471(2182), 20150464.Google Scholar
Dalwadi, M. P., O’Kiely, D., Thomson, S. J., Khaleque, T. S. & Hall, C. L. (2017) Mathematical modeling of chemical agent removal by reaction with an immiscible cleanser. SIAM J. Appl. Math. 77(6), 19371961.CrossRefGoogle Scholar
Dalwadi, M. P., Wang, Y., King, J. R. & Minton, N. P. (2018) Upscaling diffusion through first-order volumetric sinks: A homogenization of bacterial nutrient uptake. SIAM J. Appl. Math. 78(3), 13001329.CrossRefGoogle Scholar
Evans, J. D. & King, J. R. (2000) Asymptotic results for the Stefan problem with kinetic undercooling. Q. J. Mech. Appl. Math. 53(3), 449473.CrossRefGoogle Scholar
Geng, Y., Kamilova, A. A. & Luckins, E. K. (2023) Fluid-flow effects in the reactive decontamination of porous materials driven by chemical swelling or contraction. J. Eng. Math. Accepted.CrossRefGoogle Scholar
Keener, J. P. (2018). Principles of Applied Mathematics: Transformation and Approximation, CRC Press, Boca Raton.Google Scholar
Kiradjiev, K. B., Breward, C. J. W., Griffiths, I. M. & Schwendeman, D. W. (2021) A homogenized model for a reactive filter. SIAM J. Appl. Math. 81(2), 591619.CrossRefGoogle Scholar
Luckins, E., Breward, C. J. W., Griffiths, I. M. & Wilmott, Z. (2020) Homogenisation problems in reactive decontamination. Eur. J. Appl. Math. 31(5), 782805.CrossRefGoogle Scholar
Luckins, E. K., Breward, C. J. W., Griffiths, I. M. & Please, C. P. (2023) A homogenised model for the motion of evaporating fronts in porous media. Eur. J. Appl. Math. 34(4), 806837.CrossRefGoogle Scholar
Schulz, R. & Knabner, P. (2017) Derivation and analysis of an effective model for biofilm growth in evolving porous media. Math. Methods Appl. Sci. 40(8), 29302948.CrossRefGoogle Scholar
Schulz, R. & Knabner, P. (2017) An effective model for biofilm growth made by chemotactical bacteria in evolving porous media. SIAM J. Appl. Math. 77(5), 16531677.CrossRefGoogle Scholar
Shampine, L. F. & Reichelt, M. W. (1997) The MATLAB ODE suite. SIAM J. Sci. Comput. 18(1), 122.CrossRefGoogle Scholar
Wood, B. D., Golfier, F. & Quintard, M. (2011) Dispersive transport in porous media with biofilms: local mass equilibrium in simple unit cells. Int. J. Environ. Waste Manag. 1(2), 2448.CrossRefGoogle Scholar
Figure 0

Figure 1. Schematic illustrating the two canonical agent distributions within the pore space: the sharp-interface case (left) with the agent initially saturating the pore space and the agent-on-walls case (right).

Figure 1

Table 1. Dimensional parameters in the decontamination models (SI is sharp-interface model, AW is agent-on-walls model)

Figure 2

Table 2. Dimensionless parameters in the sharp-interface model (3.4)

Figure 3

Figure 2. Schematic of the asymptotic limits considered for the sharp-interface model. Case i is $\alpha \ll 1$, with sublimits ia for $\gamma \ll 1$ and ib for $\gamma \gg 1$. Case ii is $\alpha =O(1)$ and $\gamma \ll 1$, and case iii is $\alpha =O(1)$ and $\gamma \gg 1$.

Figure 4

Figure 3. Numerical and analytic solutions of the sharp-interface model (3.4), for small $\alpha$, with $\gamma =5$ (top), $\gamma =1$ (middle) and $\gamma =0.2$ (bottom). On the left we show the position, $h(t)$, of the reacting interface against time $t$ for various $\gamma$ and for various small values of $\alpha$, along with the leading-order asymptotic solution (3.7) in the limit of small $\alpha$. On the right, we show the cleanser concentration profiles $c(y)$ at times $t_D/4$, $t_D/2$, $3t_D/4$ and $t_D$ (the decontamination time, at which $h=1$), with the arrows showing the direction of increasing time, for the corresponding $\gamma$ and $\alpha$, and again showing the asymptotic small-$\alpha$ solution (3.7).

Figure 5

Figure 4. The numerically computed solution $\eta _*$ of (3.20), as a function of $\alpha$.

Figure 6

Figure 5. Numerical and analytic solutions of the sharp-interface model (3.4), examining the effect of $\gamma$, with $\alpha =1$ throughout. In (a) and (b), we show the small-$\gamma$ numerical solutions of (3.4), with the leading-order small-$\gamma$ asymptotic solution (3.12). In (c) and (d), we show large-$\gamma$ solutions with the leading-order large-$\gamma$ asymptotic solution (3.19). In (a) and (c), we show the motion of the reaction front $h(t)$. In (b) and (d), we show the cleanser concentration profiles at times $t_D/4$, $t_D/2$, $3t_D/4$ and $t_D$ (the decontamination time, at which $h=1$), with arrows indicating increasing time. In (b), we mark the end of each line with a circle for the cases $\gamma =0.1$ and $\gamma =0.01$, as otherwise the curves are difficult to distinguish.

Figure 7

Figure 6. Sharp-interface model decontamination time $t_D$ or scaled decontamination time $t_D/\alpha \gamma$ (left) and decontamination efficiency (right) as the system parameters $\alpha$ and $\gamma$ are varied. Numerical solutions of (3.4) are shown in colour, and the analytical expressions (3.8), (3.14) and (3.22) for the decontamination time, and (3.26) for the efficiency, in the various asymptotic limits, are shown in black (dashed, dotted or solid).

Figure 8

Table 3. Dimensionless parameters in the agent-on-walls model (4.6)

Figure 9

Figure 7. The variation of (a) $V$ and $D$ and (b) $F$ with agent thickness $R$. We take $\rho _0=0.2$ throughout, so that the model is valid for $\delta \in (0,0.3)$. Three values $\delta =0.01$ (dash-dotted curves), $\delta =0.1$ (dashed curves) and $\delta =0.3$ (solid curves) are shown.

Figure 10

Figure 8. Schematic of the asymptotic limits considered for the agent-on-walls model. Case I is $A\ll 1$, with sublimits $\Gamma \ll 1$ and $\Gamma \gg 1$ referred to as cases Ia and Ib, respectively. Case II is $A\gg 1$ and $\Gamma \ll 1$, while case III is $A\gg 1$ and $\Gamma \gg 1$. For all these asymptotic limits, we additionally assume $\delta \ll 1$.

Figure 11

Figure 9. Numerical and asymptotic solution of (4.6) in the case $A\ll 1$, $\delta \ll 1$ and for $\Gamma =1$ (left) and $\Gamma =0.1$ (right). We show the decontamination front $y_*(t)$ (top), cleanser profiles $c$ (middle) and agent-layer thickness $R$ (bottom). The numerical solutions (coloured curves) are computed in the limit $\delta \ll 1$, making the assumption that $D=V=F=1$, and taking $A=10^{-3}$. In (b)–(f), solutions are shown at times $0.05t_D$, $t_D/4$, $t_D/2$, $3t_D/4$ and $0.95t_D$ (where $t_D$ is the decontamination time, at which $y_*=1$), and arrows show the direction of increasing time. Analytical solutions are shown in black. In (c) and (d), thick dotted black curves are the analytical solution (4.20) and in (e) and (f) thick dashed curves are (4.27), where $y_*(t)$ is the (numerically computed) solution of (4.26). The thin dashed lines in (f) are the approximation (4.31) to $R$, and the thick black line $c=1$ in (d) is the approximation (4.31) to $c$ in the small-$\Gamma$ sublimit.

Figure 12

Figure 10. An illustration of the asymptotic solution structure for case Ib. Left: the decontamination front $y_*(t)$ against time, showing the three regions “early” ($t=O(1)$), “late” ($t=O(\Gamma )$) and “end” ($1-y_*(t)=O(1/\sqrt{\Gamma })$) times. Right: the cleanser concentration $c$ and agent-layer thickness $R$ during the late-time stage.

Figure 13

Figure 11. Numerical and asymptotic solution of (4.6) in the case $A\ll 1$, $\delta \ll 1$ and for $\Gamma =20$. The numerical solutions (solid curves) are computed in the limit $\delta \ll 1$, making the assumption that $D=V=F=1$, and taking $A=10^{-3}$. In (a), we show the decontamination front $y_*(t)$, with the full, early-time and end-time asymptotic approximations of the solution of (4.26), (4.41) and (B.9), respectively. In (b), the solid lines are the numerical solution $c$ (red) and $R$ (blue) at times $0.01t_D$, $0.1t_D$, $t_D/4$, $t_D/2$, $3t_D/4$ and $0.95t_D$, where $t_D$ is the time at which $y_*(t_D)=1$, and the arrows indicate the direction of increasing time. Black dashed ($R$) and dotted ($c$) curves are the full analytical solutions (4.20) and (4.27), where $y_*(t)$ is the (numerically computed) solution of (4.26). The magenta dashed ($R$) and dotted ($c$) curves are the approximations (4.56) in the large-$\Gamma$ sublimit.

Figure 14

Figure 12. Variation of solutions with $\delta$ in the limit $A\ll 1$, for $\Gamma =0.01$ (top), $\Gamma =1$ (middle) and $\Gamma =100$ (bottom). On the left, we show the decontamination front $y_*(t)$ and on the right the cleanser $c$ (red) and agent-layer thickness $R$ (blue) profiles. Numerical solutions of (4.6) are shown as solid lines with $\delta =0.25$, and in the limit $\delta \ll 1$ as dotted lines, for comparison. Throughout we take $A=0.01$, and $\rho _0=0.2$. Solutions in (b), (d) and (f) are shown at the times $0.1t_D$, $0.2t_D$, $0.4t_D$, $0.6t_D$, $0.8t_D$ and $t_D$, where $t_D$ is the decontamination time for the $\delta =0.25$ solution.

Figure 15

Figure 13. Asymptotic behaviour of the decontamination front $y_*$ and decontamination time $t_D$ in case II. In (a), we show the solution $\sqrt{A\Gamma }y_*$ of (4.60c) against $t$ (black) and the large-$t$ (large-$A\Gamma$) limit (4.62) (red dashed). In (b), we show the solution $t_D$ of (4.61) (blue) and the large $A\Gamma$ limit (4.63) (red dotted).

Figure 16

Figure 14. Numerical and asymptotic solution of (4.6) in the case $A\gg 1$ (and $\delta =O(A^{-1})\ll 1$). We take $A=100$ throughout with $\Gamma =0.1$ (left) and $\Gamma =100$ (right). We show the $c$ and $R$ profiles against the spatial coordinate $y$ in (a)–(b) and against the similarity variable $\sqrt{A\Gamma } y/2\sqrt{t}$ in (c)–(d). The decontamination front $y_*(t)$ is shown in (e)–(f). The numerical solutions (solid curves) are computed in the limit $\delta \ll 1$, making the assumption that $D=V=F=1$. In (a)–(d), the solid curves are the numerical solution of (4.58) ($c$ red, $R$ blue) and the black curves are the analytical solution (4.60a) for $c$ (dotted) and (4.60b) for $R$ (dashed), where $y_*(t)$ is the (numerically computed) solution of (4.60c), shown in (e)–(f), along with the numerical decontamination front.

Figure 17

Figure 15. Variation of solutions with $\delta$, holding $A=1$, for $\Gamma =0.01$ (top), $\Gamma =1$ (middle) and $\Gamma =100$ (bottom). On the left, we show the decontamination front $y_*(t)$ and on the right the cleanser $c$ (red) and agent-layer thickness $R$ (blue) profiles. Numerical solutions of (4.6) are shown as solid lines with $\delta =0.25$, and in the limit $\delta \ll 1$ (making the assumption $D=V=F=1$) as dotted lines, for comparison. Throughout we take $\rho _0=0.2$. Solutions in (b), (d) and (f) are shown at the times $0.1t_D$, $0.2t_D$, $0.4t_D$, $0.6t_D$, $0.8t_D$ and $t_D$, where $t_D$ is the decontamination time for the $\delta =0.25$ solution. The solutions for $\delta =0.25$ show clear discontinuities in $c_y$ at $y=1$ due to the discontinuous diffusivity $D(R)$.

Figure 18

Figure 16. The behaviour of the decontamination time $t_D$ (in (a)), or the scaled decontamination time $t_D/A\Gamma$ (in (b)–(c)) or (4.70) (in (d)) for the agent-on-walls model, computed from the numerical solution of (4.6). In (a)–(c), we take the small-$\delta$ limit and investigate the effect of spill depth and reaction rate on the decontamination time via $\Gamma$, and of the applied cleanser concentration via $A$. In (d), we show the effect of the initial agent-layer thickness via $\delta$.

Figure 19

Figure 17. The behaviour of the decontamination efficiency (4.64) for the agent-on-walls model, computed from the numerical solution of (4.6). In (a)–(c), we take the small-$\delta$ limit. In (a), (b) and (d), we show the effect of $\Gamma$, $A$ and $\delta$, respectively, on the efficiency. In (c), we show the variation of the efficiency with the decontamination time, parameterised by $A$.