## 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.

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

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}}$ ,

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)$ :

while at $y=0$ ,

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$
.

### 3.2. Nondimensionalisation

We nondimensionalise the model using the scalings

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$

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

and on $y=0$ ,

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.

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

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.

### 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

with boundary conditions at $y=h$ ,

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

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

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

with, on $y=h$ ,

(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

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).

### 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

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

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

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

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

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

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

We look for an expansion of the form

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

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

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

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.

The leading-order decontamination time is

so that, in terms of our original time variable,

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).

### 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

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

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

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

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}}$ .

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

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

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

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

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

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

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

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

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

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

while at $t=0$ ,

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

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

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$ .

### 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$ .

### 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

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

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

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

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

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])

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

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$ ,

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.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

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

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

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

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

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

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

and $S$ satisfies

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

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

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

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$ .

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

with solution

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

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

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.

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

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

so that, at leading order, we have

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

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

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

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

so that

with $Y_*$ the solution of

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

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

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$ )

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

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

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)$ ,

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 })$ ,

with solution

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

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