2.1.1 Boundary conditions

On the channel walls we impose no-flux and no-slip conditions on the exterior flow and a no-flux condition on the Darcy flow, so that

(2.3a,b ) $$\begin{eqnarray}\boldsymbol{u}=\mathbf{0}\quad \text{on }z=0,1\text{ for }|x|>1;\quad W=0\quad \text{on }z=0,1\text{ for }|x|<1.\end{eqnarray}$$

On the inflow ( $x=-1$ ) and outflow ( $x=1$ ) interfaces, continuity of flux, continuity of pressure and a no tangential slip condition (Levy & Sanchez-Palencia Reference Levy and Sanchez-Palencia1975; Carraro et al. Reference Carraro, Goll, Marciniak-Czochra and Mikelić2015) are given by

(2.4a-c ) $$\begin{eqnarray}u=U,\quad p=P,\quad w=0\quad \text{on }x=\pm 1\text{ for }0<z<1.\end{eqnarray}$$

The flow is driven by imposing a constant, upstream pressure gradient in the far field, with

(2.5) $$\begin{eqnarray}\displaystyle \frac{\text{d}p}{\text{d}x}\rightarrow -12\quad \text{as }x\rightarrow -\infty , & & \displaystyle\end{eqnarray}$$

ensuring that the cross-sectionally averaged horizontal component of velocity has magnitude $1$ in the far field.

2.3.1 Outer regions

In the outer regions I, IV and VII, we begin by posing asymptotic expansions of the form

(2.6) $$\begin{eqnarray}\displaystyle (\boldsymbol{u},\boldsymbol{Q},p,P)=(\boldsymbol{u}_{0},\boldsymbol{Q}_{0},p_{0},P_{0})+{\it\epsilon}(\boldsymbol{u}_{1},\boldsymbol{Q}_{1},p_{1},P_{1})+O({\it\epsilon}^{2})\quad \text{as }{\it\epsilon}\rightarrow 0. & & \displaystyle\end{eqnarray}$$

Substituting (2.6) into (2.1), we find that the leading-order governing equations for the exterior flow are simply the lubrication equations

(2.7a-c ) $$\begin{eqnarray}0=-p_{0x}+u_{0zz},\quad 0=-p_{0z},\quad 0=u_{0x}+w_{0z}\quad \text{for }|x|>1,0<z<1.\end{eqnarray}$$

Substituting (2.6) into (2.2), the leading-order governing equations for the interior flow are

(2.8a-c ) $$\begin{eqnarray}U_{0}=-KP_{0x},\quad 0=-KP_{0z},\quad 0=U_{0x}+W_{0z}\quad \text{for }|x|<1,0<z<1.\end{eqnarray}$$

At leading-order, the channel wall boundary conditions (2.3) become

(2.9a ) $$\begin{eqnarray}\displaystyle & \displaystyle \boldsymbol{u}_{0}=\mathbf{0}\quad \text{on }z=0,1\text{ for }|x|>1, & \displaystyle\end{eqnarray}$$

(2.9b ) $$\begin{eqnarray}\displaystyle & \displaystyle W_{0}=0\quad \text{on }z=0,1\text{ for }|x|<1, & \displaystyle\end{eqnarray}$$

(2.10) $$\begin{eqnarray}\displaystyle \frac{\text{d}p_{0}}{\text{d}x}\rightarrow -12\quad \text{as }x\rightarrow -\infty . & & \displaystyle\end{eqnarray}$$

We are unable to impose at leading-order the interfacial boundary conditions (2.4) in the outer regions; these conditions are satisfied within the inner regions described in § 2.3.2. For the outer problems, we must instead impose matching conditions (using the matching law given in Van Dyke (Reference Van Dyke1975)). These are determined in § 2.3.2, and we state them here to close the leading-order outer problem. As the outer pressure is independent of $z$ and the fluid velocity is proportional to the pressure gradient, we find that we may impose continuity of total flux and pressure at each interface (corresponding to Neumann and Dirichlet conditions on the pressure, respectively). That is,

(2.11a,b ) $$\begin{eqnarray}\int _{0}^{1}u_{0}\,\text{d}z=\int _{0}^{1}U_{0}\,\text{d}z,\quad p_{0}=P_{0}\quad \text{on }x=\pm 1.\end{eqnarray}$$

We note that whilst the conditions (2.11) close the outer problem, they give no information about, for example, the leading-order stress acting on the interfaces.

The leading-order solutions are readily obtained as follows: in region I,

(2.12a,b ) $$\begin{eqnarray}\boldsymbol{u}_{0}=6z(1-z)\boldsymbol{e}_{x},\quad p_{0}=-12x;\end{eqnarray}$$

in region IV,

(2.13a,b ) $$\begin{eqnarray}\boldsymbol{U}_{0}=\boldsymbol{e}_{x},\quad P_{0}=-\frac{1+x}{K}+12;\end{eqnarray}$$

and in region VII,

(2.14a,b ) $$\begin{eqnarray}\boldsymbol{u}_{0}=6z(1-z)\boldsymbol{e}_{x},\quad p_{0}=-12x+\left(24-\frac{2}{K}\right).\end{eqnarray}$$

The pressure is determined up to an arbitrary constant which (at this order) we fix by supposing $p_{0}+12x\rightarrow 0$ as $x\rightarrow -\infty$ , without loss of generality. We remark that the presence of the porous plug introduces a jump in pressure of $24-2/K\sim -2/K$ (as $K$ is small). Equivalently, an $O(1)$ pressure drop creates an $O(K)$ flow through the porous plug.

We note that the horizontal components of velocity (2.12a ) and (2.13a ) are incompatible with the continuity of flux condition (2.4a ) on the interfacial inflow boundary for the full problem. Similarly for outflow, we find that (2.13a ) and (2.14a ) are incompatible with (2.4a ). In the next section we introduce inner regions to resolve this issue, but for the remainder of this section we consider higher-order terms in the outer regions as the higher-order pressure is coupled to the leading-order flow in the inner regions.

Using the leading-order solutions (2.12) and (2.14), the $O({\it\epsilon})$ terms in the governing equations for the exterior flow (2.1) become

(2.15a-c ) $$\begin{eqnarray}0=-p_{1x}+u_{1zz},\quad 0=-p_{1z},\quad 0=u_{1x}+w_{1z}\quad \text{for }|x|>1,0<z<1.\end{eqnarray}$$

Similarly, using the leading-order solutions (2.13), the $O({\it\epsilon})$ terms in the governing equations for the interior flow (2.2) become

(2.16a-c ) $$\begin{eqnarray}U_{1}=-KP_{1x},\quad 0=-KP_{1z},\quad 0=U_{1x}+W_{1z}\quad \text{for }|x|<1,0<z<1.\end{eqnarray}$$

The channel wall boundary conditions (2.3) at this order are as follows

(2.17a ) $$\begin{eqnarray}\displaystyle & \displaystyle \boldsymbol{u}_{1}=\mathbf{0}\quad \text{on }z=0,1\text{ for }|x|>1, & \displaystyle\end{eqnarray}$$

(2.17b ) $$\begin{eqnarray}\displaystyle & \displaystyle W_{1}=0\quad \text{on }z=0,1\text{ for }|x|<1, & \displaystyle\end{eqnarray}$$

(2.18) $$\begin{eqnarray}\displaystyle \frac{\text{d}p_{1}}{\text{d}x}\rightarrow 0\quad \text{as }x\rightarrow -\infty . & & \displaystyle\end{eqnarray}$$

The matching conditions across each interface are determined in § 2.3.2, and we give them here to close the problem. For inflow, we have

(2.19a,b ) $$\begin{eqnarray}\int _{0}^{1}u_{1}\,\text{d}z=\int _{0}^{1}U_{1}\,\text{d}z,\quad p_{1}-P_{1}={\it\Pi}^{-}(Re,K)\quad \text{on }x=-1,\end{eqnarray}$$

and for outflow we have

(2.20a,b ) $$\begin{eqnarray}\int _{0}^{1}u_{1}\,\text{d}z=\int _{0}^{1}U_{1}\,\text{d}z,\quad p_{1}-P_{1}={\it\Pi}^{+}(Re,K)\quad \text{on }x=1,\end{eqnarray}$$

where the pressure jumps ${\it\Pi}^{-}$ and ${\it\Pi}^{+}$ are functions of the parameters $Re$ and $K$ . In § 2.3.2, we specify the problems that need to be solved to determine these functions, while in §§ 2.4 and 2.5 we determine their asymptotic behaviour for small and large $Re$ .

Given ${\it\Pi}^{\pm }$ , the solution to (2.15)–(2.20) is readily determined: in region I,

(2.21a,b ) $$\begin{eqnarray}\boldsymbol{u}_{1}=\mathbf{0},\quad p_{1}={\it\Pi}^{-};\end{eqnarray}$$

in region IV,

(2.22a,b ) $$\begin{eqnarray}\boldsymbol{Q}_{1}=\mathbf{0},\quad P_{1}=0;\end{eqnarray}$$

and in region VII,

(2.23a,b ) $$\begin{eqnarray}\boldsymbol{u}_{1}=\mathbf{0},\quad p_{1}={\it\Pi}^{+}.\end{eqnarray}$$

As with the leading-order solutions (2.12)–(2.14), the pressure field is unique up to an arbitrary constant. For convenience, we choose the constant such that $P_{1}(0,0.5)=0$ , without loss of generality. We proceed by considering the inner regions, highlighting the problems that would need to be solved numerically to determine the functions ${\it\Pi}^{-}$ and ${\it\Pi}^{+}$ .

2.3.2 Inner regions

The appropriate scalings in the inflow and outflow inner regions are $x=\mp 1+{\it\epsilon}X^{\mp }$ , respectively, with $(w,W)=(\overline{w},\overline{W})/{\it\epsilon}$ . To make it clear that we are working in the inner regions we introduce overbars on the other dimensionless variables, writing $(u,U)=(\overline{u},\overline{U})$ , $(p,P)=(\overline{p},\overline{P})$ and $(\overline{\boldsymbol{u}},\overline{\boldsymbol{Q}})=(\overline{u},\overline{U})\boldsymbol{e}_{x}+(\overline{w},\overline{W})\boldsymbol{e}_{z}$ . Using vector operators in terms of $X^{\mp }$ and $z$ as appropriate, the exterior fluid governing equations (2.1) become

(2.24a,b ) $$\begin{eqnarray}Re(\overline{\boldsymbol{u}}\boldsymbol{\cdot }\boldsymbol{{\rm\nabla}})\overline{\boldsymbol{u}}=-{\it\epsilon}^{-1}\boldsymbol{{\rm\nabla}}\overline{p}+{\rm\nabla}^{2}\overline{\boldsymbol{u}},\quad 0=\boldsymbol{{\rm\nabla}}\boldsymbol{\cdot }\overline{\boldsymbol{u}}\quad \text{in regions II and VI},\end{eqnarray}$$

which are valid when $X^{-}<0$ for inflow (region II) and $X^{+}>0$ for outflow (region VI), with $0<z<1$ in both cases. The interior fluid governing equations (2.2) become

(2.25a,b ) $$\begin{eqnarray}\overline{\boldsymbol{Q}}=-{\it\epsilon}^{-1}K\boldsymbol{{\rm\nabla}}\overline{P},\quad 0=\boldsymbol{{\rm\nabla}}\boldsymbol{\cdot }\overline{\boldsymbol{Q}}\quad \text{in regions III and V},\end{eqnarray}$$

which are valid for $X^{-}>0$ for inflow (region III) and $X^{+}<0$ for outflow (region V), with $0<z<1$ in both cases.

Similarly to (2.6), we form inner expansions in powers of ${\it\epsilon}$ as follows:

(2.26) $$\begin{eqnarray}\displaystyle (\overline{\boldsymbol{u}},\overline{\boldsymbol{Q}},\overline{p},\overline{P})=(\overline{\boldsymbol{u}}_{0},\overline{\boldsymbol{Q}}_{0},\overline{p}_{0},\overline{P}_{0})+{\it\epsilon}(\overline{\boldsymbol{u}}_{1},\overline{\boldsymbol{Q}}_{1},\overline{p}_{1},\overline{P}_{1})+O({\it\epsilon}^{2})\quad \text{as }{\it\epsilon}\rightarrow 0. & & \displaystyle\end{eqnarray}$$

At leading order, it follows from (2.24a ) that

(2.27) $$\begin{eqnarray}\displaystyle \mathbf{0}=-\boldsymbol{{\rm\nabla}}\overline{p}_{0}\quad \text{in regions II and VI}, & & \displaystyle\end{eqnarray}$$

while (2.25a ) gives, similarly,

(2.28) $$\begin{eqnarray}\displaystyle \mathbf{0}=-K\boldsymbol{{\rm\nabla}}\overline{P}_{0}\quad \text{in regions III and V}. & & \displaystyle\end{eqnarray}$$

The appropriate interfacial conditions on $X^{\mp }=0$ for the leading-order pressure are given by the leading-order terms in (2.4b ), which become

(2.29) $$\begin{eqnarray}\displaystyle \overline{p}_{0}=\overline{P}_{0}\quad \text{on }X^{\mp }=0. & & \displaystyle\end{eqnarray}$$

Thus $p_{0}$ and $P_{0}$ are constant and equal in regions II and III and in regions V and VI, leading to the matching condition (2.11b ). Using the leading-order outer solutions (2.12)–(2.14) for the pressure, we find

(2.30a,b ) $$\begin{eqnarray}\overline{p}_{0}=12\quad \text{in region II},\quad \overline{P}_{0}=12\quad \text{in region III},\end{eqnarray}$$

for inflow, and

(2.31a,b ) $$\begin{eqnarray}\overline{P}_{0}=12-\frac{2}{K}\quad \text{in region V},\quad \overline{p}_{0}=12-\frac{2}{K}\quad \text{in region VI}\end{eqnarray}$$

for outflow. We note that these leading-order pressures are constant on the interface, as predicted by Levy & Sanchez-Palencia (Reference Levy and Sanchez-Palencia1975).

As the leading-order fluid pressures are constant within each inner region, the leading-order inner flow is driven by the first correction to pressure, and these equations are given by the $O(1)$ terms in (2.24) and (2.25). We illustrate the leading-order flow problems in figure 4 for inflow, and in figure 5 for outflow. In both cases, the inner flow transition problem is governed by the two-dimensional steady Navier–Stokes equations coupled to the Darcy equations. To fully determine the inflow and outflow pressure jumps ${\it\Pi}^{-}$ and ${\it\Pi}^{+}$ , we would have to solve the problems indicated in figures 4 and 5 numerically for two parameters: $Re$ and $K$ .

Figure 4. The inner inflow problem for $Re=O(1)$ . The flow is from left to right.

Figure 5. The inner outflow problem for $Re=O(1)$ . The flow is from left to right.

However, to couple the outer regions, we require another condition on the pressure. We note that we can analytically determine suitable coupling conditions for the horizontal components of velocity averaged over the channel height by appealing to global conservation of mass, rather than solving the full problem. By integrating the continuity equations across the channel height, using the conditions for the flow in the normal direction at each boundary and matching with the outer solutions, we deduce the following coupling conditions for the outer horizontal components of velocity

(2.32a,b ) $$\begin{eqnarray}\int _{0}^{1}u_{0}\,\text{d}z=\int _{0}^{1}U_{0}\,\text{d}z,\quad \int _{0}^{1}u_{1}\,\text{d}z=\int _{0}^{1}U_{1}\,\text{d}z\quad \text{on }x=\pm 1.\end{eqnarray}$$

We proceed by determining the inner flow in the asymptotic sublimits $Re\ll 1$ and $1\ll 1/K\ll Re\ll 1/{\it\epsilon}$ , starting with the former. This will allow us to calculate quantities that are unobtainable from sole consideration of the outer problem, such as the stress acting on the interface and the limiting behaviour of the functions ${\it\Pi}^{-}$ , ${\it\Pi}^{+}$ .

2.5.1 Inflow

The inflow inner problem (within regions II and III) is shown in figure 4. The limit as ${\it\delta}\rightarrow 0$ yields leading-order inviscid flow in region II away from the boundaries, with five additional boundary layers being required to satisfy all of the boundary conditions. Two of these viscous boundary layers are on the channel walls, one is on the porous interface, whilst the remaining boundary layers occur in the corners between the channel wall and porous interface. The boundary layer structure for the inflow problem is shown in figure 7.

Figure 7. Schematic diagram of the inner regions IIa–c within the inflow boundary region II. The light region is the exterior fluid and the dark region is the porous obstacle. The size of the boundary layers has been exaggerated for illustrative purposes.

We start by considering the coupled regions II and III, and show that the transition to Poiseuille flow occurs solely within region III at leading order in ${\it\delta}$ .

2.5.2 Regions II and III

Within regions II and III, we make the following asymptotic expansions

(2.36a ) $$\begin{eqnarray}\displaystyle & \displaystyle \overline{\boldsymbol{u}}_{0}=6z(1-z)\boldsymbol{e}_{x}+{\it\delta}\overline{\boldsymbol{u}}_{01}+{\it\delta}^{4/3}\overline{\boldsymbol{u}}_{02}+O({\it\delta}^{5/3}), & \displaystyle\end{eqnarray}$$

(2.36b ) $$\begin{eqnarray}\displaystyle & \displaystyle \overline{\boldsymbol{Q}}_{0}=\overline{\boldsymbol{Q}}_{00}+{\it\delta}\overline{\boldsymbol{Q}}_{01}+{\it\delta}^{4/3}\overline{\boldsymbol{Q}}_{02}+O({\it\delta}^{5/3}), & \displaystyle\end{eqnarray}$$

(2.36c ) $$\begin{eqnarray}\displaystyle & \displaystyle (\overline{p}_{1},\overline{P}_{1})=(\overline{p}_{11},\overline{P}_{11})+{\it\delta}^{1/3}(\overline{p}_{12},\overline{P}_{12})+O({\it\delta}^{2/3}) & \displaystyle\end{eqnarray}$$

${\it\delta}\rightarrow 0$

${\it\delta}$

Substituting the asymptotic expansions (2.36) into the governing equations presented in figure 4, the $O(1/{\it\delta})$ terms in region II are trivially satisfied, but so far $\overline{\boldsymbol{Q}}_{00}$ is unknown. Proceeding with the expansions, we find at $O(1)$ that $\overline{\boldsymbol{u}}_{01}=(\overline{u}_{01},\overline{w}_{01})$ is governed by

(2.37a-c ) $$\begin{eqnarray}\overline{u}_{00}\overline{u}_{01X}+\overline{u}_{00z}\overline{w}_{01}=-\overline{p}_{11X}-12,\quad \overline{u}_{00}\overline{w}_{01X}=-\overline{p}_{11z},\quad \boldsymbol{{\rm\nabla}}\boldsymbol{\cdot }\overline{\boldsymbol{u}}_{01}=0\quad \text{in region II},\end{eqnarray}$$

where $\overline{u}_{00}=6z(1-z)$ , while $\overline{\boldsymbol{Q}}_{00}=(\overline{U}_{00},\overline{W}_{00})$ is governed by

(2.38a,b ) $$\begin{eqnarray}\overline{\boldsymbol{Q}}_{00}=-K\boldsymbol{{\rm\nabla}}\overline{P}_{11},\quad 0=\boldsymbol{{\rm\nabla}}\boldsymbol{\cdot }\overline{\boldsymbol{Q}}_{00}\quad \text{in region III}.\end{eqnarray}$$

The boundary conditions we impose at this order in region II are again obtained via matching. On the channel walls (via matching with region IIb in appendix A), we have no flux, so that

(2.39) $$\begin{eqnarray}\displaystyle \overline{w}_{01}=0\quad \text{on }z=0,1\text{ for }X^{-}<0, & & \displaystyle\end{eqnarray}$$

while in the far field we match with region I, to obtain

(2.40a,b ) $$\begin{eqnarray}\overline{\boldsymbol{u}}_{01}\rightarrow \mathbf{0},\quad \overline{p}_{11}\sim -12X^{-}+{\it\Pi}_{1}^{-}\quad \text{as }X^{-}\rightarrow -\infty ~\text{for }0<z<1.\end{eqnarray}$$

Similarly, the boundary conditions in region III are

(2.41) $$\begin{eqnarray}\displaystyle \overline{W}_{00}=0\quad \text{on }z=0,1\text{ for }X^{-}>0, & & \displaystyle\end{eqnarray}$$

while matching with region IV gives

(2.42) $$\begin{eqnarray}\displaystyle \overline{P}_{11}\sim -\frac{X^{-}}{K}\quad \text{as }X^{-}\rightarrow \infty ~\text{for }0<z<1. & & \displaystyle\end{eqnarray}$$

We note that the far-field matching condition $\overline{\boldsymbol{Q}}_{00}\sim \boldsymbol{e}_{x}$ as $X^{-}\rightarrow \infty$ follows consistently from (2.42) using (2.38a ).

The interfacial conditions at $X^{-}=0$ for (2.37)–(2.38) are formally derived from consideration of the inner–inner boundary layer (region IIa) in § 2.5.3, within which we are able to satisfy the no tangential slip condition, but additionally find that the leading-order flow and the first correction to the pressure $\overline{p}_{11}$ are unchanged. From § 2.5.3, the correct conditions to couple regions II and III are continuity of flux and pressure, with

(2.43a,b ) $$\begin{eqnarray}6z(1-z)=\overline{U}_{00},\quad \overline{p}_{11}=\overline{P}_{11}\quad \text{on }X^{-}=0\text{ for }0<z<1.\end{eqnarray}$$

The system to solve in regions II and III (up to the $O({\it\delta})$ terms) is given by (2.37)–(2.43). We see that regions II and III have decoupled: the boundary condition (2.43a ) allows us to solve for $\overline{\boldsymbol{Q}}_{00}$ and $\overline{P}_{11}$ using (2.38) with (2.41)–(2.42), and then we can use the boundary condition (2.43b ) to determine $\overline{\boldsymbol{u}}_{01}$ and $\overline{p}_{11}$ using (2.37) with (2.39)–(2.40). The governing equations (2.38), with boundary conditions (2.41)–(2.42) and (2.43a ), are solved by

(2.44) $$\begin{eqnarray}\displaystyle -K\overline{P}_{11}=X^{-}+\mathop{\sum }_{k=1}^{\infty }A_{k}\exp (-2k{\rm\pi}X^{-})\cos (2k{\rm\pi}z),\quad A_{k}=\frac{3}{({\rm\pi}k)^{3}}. & & \displaystyle\end{eqnarray}$$

Since $\overline{\boldsymbol{u}}_{00}$ is independent of $X^{-}$ within region II, the leading-order transition from Poiseuille to plug flow occurs entirely within the porous medium in region III, and is governed by (2.44).

We now complete the solution in region II by determining $\overline{u}_{01}$ and $\overline{p}_{11}$ . As well as allowing us to determine ${\it\Pi}_{1}^{-}$ , knowledge of $\overline{p}_{11}$ gives the normal stress acting on the interface up to $O({\it\epsilon})$ , the first order at which the stress varies in the $z$ -direction. A rearrangement of the exterior flow (2.37) allows us to decouple the region II system (where $-\infty <X^{-}<0$ and $0<z<1$ ) into

(2.45a ) $$\begin{eqnarray}\displaystyle & \displaystyle z(1-z){\rm\nabla}^{2}\overline{w}_{01}+2\overline{w}_{01}=0, & \displaystyle\end{eqnarray}$$

(2.45b ) $$\begin{eqnarray}\displaystyle & \displaystyle {\rm\nabla}^{2}\overline{p}_{11}=-2\overline{u}_{00z}\overline{w}_{01X^{-}}. & \displaystyle\end{eqnarray}$$

We solve (2.45a ) for $\overline{w}_{01}$ first, using boundary conditions (2.39) on $z=0$ , $1$ and (2.40a ) as $X^{-}\rightarrow -\infty$ . We derive a consistent boundary condition for $\overline{w}_{01}$ on $X^{-}=0$ using the boundary condition (2.43b ) for the pressure on $X^{-}=0$ , in combination with the governing equation (2.37b ), giving

(2.46) $$\begin{eqnarray}\displaystyle \overline{w}_{01X^{-}}=-\frac{1}{K\overline{u}_{00}}\mathop{\sum }_{k=1}^{\infty }2k{\rm\pi}A_{k}\sin (2k{\rm\pi}z)\quad \text{on }X^{-}=0\text{ for }0<z<1. & & \displaystyle\end{eqnarray}$$

The boundary condition (2.46) ensures that $\overline{w}_{01}$ does not vanish on $X^{-}=0$ , inducing a slip velocity on this interface within region II. This issue is taken care of in an inner–inner boundary layer closer to the interface (region IIa in figure 7) in § 2.5.3. The linearity of the governing equation (2.45a ), with boundary conditions (2.39), (2.40a ) and (2.46), also allows us to deduce that $\overline{w}_{01}$ scales with $K^{-1}$ . Further, from the continuity equation (2.37c ) it follows that $\overline{u}_{01}$ also scales with $K^{-1}$ .

We solve for $\overline{w}_{01}$ using a standard central finite-difference scheme. By numerically integrating $\overline{w}_{01}$ with respect to $X^{-}$ , we are able to plot the streamlines in figure 8. In figure 8(a), we show the streamlines for the first correction to the velocity, and in figure 8(b), we show the streamlines for the full velocity up to this first correction. We see that the fluid moves from the centre of the channel towards the channel walls, as expected. Therefore, adjustment from Poiseuille to plug flow occurs in region II, but at higher order. This movement induces a horizontal slip velocity on the channel walls, which is resolved in appendix A by introducing a boundary layer near the channel walls in which $z=O({\it\delta}^{1/3})$ on the bottom wall, with a similar layer near the top wall. These regions are labelled as IIb in figure 7. In appendix A, we show that the flow near $z=0$ has the form

(2.47) $$\begin{eqnarray}\displaystyle \overline{u}_{0}\sim 6z(1-z)+{\it\delta}\mathop{\sum }_{k=1}^{\infty }C_{k}\exp ({\it\mu}_{k}X^{-})g_{k}^{\prime }(0)\frac{\displaystyle \int _{0}^{z/{\it\delta}^{1/3}}\text{Ai}(a_{k}v)\,\text{d}v}{\displaystyle \int _{0}^{\infty }\text{Ai}(a_{k}s)\,\text{d}s}\quad \text{as }z\rightarrow 0^{+}, & & \displaystyle\end{eqnarray}$$

where $C_{k}$ , ${\it\mu}_{k}$ , and $g_{k}$ are described in (A1)–(A 2), $\text{Ai}$ is the Airy function of the first kind and $a_{k}$ is defined below (A 10). From (2.47), we see that the leading-order solution is preserved near the channel wall. In fact, the slip is resolved within higher-order terms in the boundary layer, and does not affect the leading-order solution in region II.

Figure 8. The streamlines for inflow near the interface. (a) The streamlines for the scaled first-correction velocity $K\overline{\boldsymbol{u}}_{01}$ , which is independent of $K$ . The fluid moves from the centre of the channel to the channel walls. (b) The streamlines for the full velocity up to first correction $\overline{\boldsymbol{u}}_{00}+{\it\delta}\overline{\boldsymbol{u}}_{01}$ (with a large value of ${\it\delta}/K=1$ taken to emphasise the effect of the correction). The fluid moves from left to right.

The problem for $\overline{p}_{11}$ is given by (2.45b ) with boundary conditions (2.40b ) and (2.43b ). We obtain consistent boundary conditions for $\overline{p}_{11}$ on the channel walls by substituting the leading-order flow solution into the governing equation (2.37b ) to obtain

(2.48) $$\begin{eqnarray}\displaystyle \overline{p}_{11z}=0\quad \text{on }z=0,1\text{ for }X^{-}<0. & & \displaystyle\end{eqnarray}$$

We solve the full system for $\overline{p}_{11}$ , given by (2.45b ) with boundary conditions (2.40b ), (2.43b ), and (2.48) numerically using a central finite-difference scheme. Using a similar scaling argument as for $\overline{w}_{01}$ , we find that ${\it\Pi}_{1}^{-}$ scales with $1/K$ . We are able to calculate ${\it\Pi}_{1}^{-}$ by plotting $(\overline{p}_{11}+12X^{-})K$ , and using the matching condition (2.40b ). As ${\it\Pi}_{1}^{-}$ is independent of $z$ , and to more easily show the pressure drop on a graph, we plot $(\int _{0}^{1}\overline{p}_{11}\,\text{d}z+12X^{-})K$ for varying $X^{-}$ in figure 9. We find thereby that ${\it\Pi}_{1}^{-}={\it\alpha}/K$ , where ${\it\alpha}=-0.117$ to three significant figures.

Figure 9. The function $(\int _{0}^{1}\overline{p}_{11}\,\text{d}z+12X^{-})K$ for varying $X^{-}$ in region II, which allows us to determine the outer pressure drop for inflow.

2.5.3 Inner–inner region IIa

We now investigate the inner–inner boundary layer closer to the interface, region IIa in figure 7. This boundary layer is required to satisfy the no tangential slip condition on the interface, and captures a balance between the relevant inertial terms and the viscous terms induced by the presence of the porous interface.

The relevant scalings are $X^{-}={\it\delta}\widehat{X}$ and $\overline{w}_{0}={\it\delta}\widehat{w}_{0}$ , and we introduce the new variables $\overline{u}_{0}=\widehat{u}_{0}$ and $\overline{p}_{1}=\widehat{p}_{1}$ to signify that we are working in region IIa. Under these scalings, in region IIa (where $-\infty <\widehat{X}<0$ and $0<z<1$ ) the exterior flow equations in figure 4 become

(2.49a ) $$\begin{eqnarray}\displaystyle & \displaystyle {\it\delta}^{-1}\widehat{u}_{0}\widehat{u}_{0\widehat{X}}+{\it\delta}\widehat{w}_{0}\widehat{u}_{0z}=-\widehat{p}_{1\widehat{X}}+{\it\delta}^{-1}\widehat{u}_{0\widehat{X}\widehat{X}}+{\it\delta}\widehat{u}_{0zz}, & \displaystyle\end{eqnarray}$$

(2.49b ) $$\begin{eqnarray}\displaystyle & \displaystyle \widehat{u}_{0}\widehat{w}_{0\widehat{X}}+{\it\delta}^{2}\widehat{w}_{0}\widehat{w}_{0z}=-{\it\delta}\widehat{p}_{1z}+\widehat{w}_{0\widehat{X}\widehat{X}}+{\it\delta}^{2}\widehat{w}_{0zz}, & \displaystyle\end{eqnarray}$$

(2.49c ) $$\begin{eqnarray}\displaystyle & \displaystyle 0=\widehat{u}_{0\widehat{X}}+{\it\delta}^{2}\widehat{w}_{0z}, & \displaystyle\end{eqnarray}$$

(2.50a-c ) $$\begin{eqnarray}\widehat{u}_{0}=\overline{U}_{0},\quad \widehat{p}_{1}=\overline{P}_{1},\quad \widehat{w}_{0}=0\quad \text{on }\widehat{X}=0\text{ for }0<z<1.\end{eqnarray}$$

The purpose of this boundary layer is to determine the tangential velocity $\widehat{w}_{0}$ , and to ensure that it is possible to satisfy the no tangential slip boundary condition on $X^{-}=0$ . Therefore, we form asymptotic expansions as follows:

(2.51a ) $$\begin{eqnarray}\displaystyle & \displaystyle \widehat{u}_{0}=6z(1-z)+O({\it\delta}), & \displaystyle\end{eqnarray}$$

(2.51b ) $$\begin{eqnarray}\displaystyle & \displaystyle \widehat{w}_{0}=\widehat{w}_{01}+O({\it\delta}^{4/3}), & \displaystyle\end{eqnarray}$$

(2.51c ) $$\begin{eqnarray}\displaystyle & \displaystyle \widehat{p}_{1}=\overline{p}_{11}(0,z)+O({\it\delta}^{1/3}) & \displaystyle\end{eqnarray}$$

${\it\delta}\rightarrow 0$

$\overline{p}_{11}$

${\it\delta}$

(2.52) $$\begin{eqnarray}\displaystyle 6z(1-z)\widehat{w}_{01\widehat{X}}=\widehat{w}_{01\widehat{X}\widehat{X}}. & & \displaystyle\end{eqnarray}$$

Using (2.44), the leading-order interfacial boundary conditions (2.50) become

(2.53a-c ) $$\begin{eqnarray}6z(1-z)=-K\overline{P}_{11\widehat{X}},\quad \overline{p}_{11}=\overline{P}_{11},\quad \widehat{w}_{01}=0\quad \text{on }\widehat{X}=0\text{ for }0<z<1.\end{eqnarray}$$

The relevant far-field condition arises via matching with region II, giving

(2.54) $$\begin{eqnarray}\displaystyle \widehat{w}_{01}\rightarrow \overline{w}_{01}(0,z)\quad \text{as }\widehat{X}\rightarrow -\infty ~\text{for }0<z<1, & & \displaystyle\end{eqnarray}$$

where $\overline{w}_{01}$ is determined in § 2.5.2.

The interfacial conditions (2.53a,b ) are used to couple regions II and III in § 2.5.2. The remaining terms in the system (2.52)–(2.54) yield $\widehat{w}_{01}$ , and are solved by

(2.55) $$\begin{eqnarray}\displaystyle \widehat{w}_{01}=\overline{w}_{01}(0,z)(1-\exp (6z(1-z)\widehat{X})). & & \displaystyle\end{eqnarray}$$

We see that (2.55) can only be a solution for inflow, as we require $6z(1-z)\widehat{X}<0$ as we move back into region II. For outflow, this inner–inner region is not a valid boundary layer, and we must tackle the problem in a different manner.

We note that this boundary layer persists for alternate tangential slip conditions. We can showcase the robust nature of the boundary layer via two examples. The first example uses the general tangential slip condition derived in Carraro et al. (Reference Carraro, Goll, Marciniak-Czochra and Mikelić2015). That is, instead of (2.53c ), we consider a condition of the form

(2.56) $$\begin{eqnarray}\displaystyle \overline{w}={\it\lambda}\overline{u}\quad \text{on }X^{-}=0,\text{ for }0<z<1, & & \displaystyle\end{eqnarray}$$

which reduces to (2.53c ) when ${\it\lambda}=0$ . Here, ${\it\lambda}$ is a slip coefficient, obtained by solving the cell problem given in Carraro et al. (Reference Carraro, Goll, Marciniak-Czochra and Mikelić2015). Using this condition, the transverse velocity is

(2.57) $$\begin{eqnarray}\displaystyle \overline{w}\sim 6{\it\lambda}z(1-z)\exp (6z(1-z)\widehat{X})+{\it\delta}\overline{w}_{01}(0,z)(1-\exp (6z(1-z)\widehat{X}))+O({\it\delta}^{4/3}),\quad & & \displaystyle\end{eqnarray}$$

in region IIa. In the second example, we use a Beavers and Joseph slip condition (Beavers & Joseph Reference Beavers and Joseph1967) of the form

(2.58) $$\begin{eqnarray}\displaystyle \overline{w}_{\widehat{X}}=\frac{\hat{{\it\lambda}}}{\sqrt{K}}(\overline{W}-\overline{w})\quad \text{on }X^{-}=0,\text{ for }0<z<1, & & \displaystyle\end{eqnarray}$$

where $\hat{{\it\lambda}}$ is the experimentally determined slip coefficient. In this case, the transverse velocity would be

(2.59) $$\begin{eqnarray}\displaystyle \overline{w}\sim {\it\delta}\left(\overline{w}_{01}(0,z)+\frac{\hat{{\it\lambda}}\overline{W}_{0}}{6\sqrt{K}z(1-z)}\right)+O({\it\delta}^{4/3}). & & \displaystyle\end{eqnarray}$$

Hence the boundary layer structure persists if the no-slip condition is modified in either of these two ways. Physically, we can conclude that the thickness of each boundary layer involved in inflow is unchanged for any of the three boundary conditions applied at the interface. Moreover, the flow induced by each of the three boundary conditions is only different in region IIa, the boundary layer closest to the interface.

Returning to the original no-slip boundary condition (2.50c ), it is straightforward to calculate that the components of stress on the interface at $x=-1$ are given by

(2.60a ) $$\begin{eqnarray}\displaystyle & \displaystyle {\it\sigma}_{xx}=-(12+{\it\epsilon}\overline{p}_{11}(0,z))+O({\it\delta}^{1/3}{\it\epsilon}), & \displaystyle\end{eqnarray}$$

(2.60b ) $$\begin{eqnarray}\displaystyle & \displaystyle {\it\sigma}_{xz}=-{\it\epsilon}(6z(1-z)\overline{w}_{01}(0,z))+O({\it\delta}^{1/3}{\it\epsilon}) & \displaystyle\end{eqnarray}$$

${\it\epsilon},{\it\delta}\rightarrow 0$

${\it\epsilon}\ll {\it\delta}\ll 1$

$O(1)$

$P(0,0.5)=0$

$O({\it\epsilon})$

$z$

$O({\it\epsilon})$

$z$

$z$

$\overline{p}_{11}$

$\overline{w}_{01}$

$1/K$

$1/K$

$O({\it\epsilon})$

Figure 10. The interfacial stress predicted by the inner flows (from (2.60)) at $O({\it\epsilon})$ . (a) The normal stress (we plot $K\overline{p}_{11}(0,z)$ ) and (b) the shear stress (we plot $6z(1-z)K\overline{w}_{01}(0,z)$ ). As both the normal and shear stress scale with $1/K$ at this order, both plotted terms are independent of $K$ .

Additionally, since we have imposed $P(0,z)=0$ , we have no $O({\it\epsilon}/{\it\delta})$ constant term from the pressure in the normal stress (2.60a ). In three dimensions, the tangential (that is, in the $(\boldsymbol{e}_{x}\times \boldsymbol{e}_{z})$ -direction) variation in pressure ensures that there will be an $O({\it\epsilon}/{\it\delta})$ term in the corresponding normal stress component no matter how we define the pressure (as we show in § 3.6.1).

We have now solved for the flow variables within the inflow boundary layers for the regions that affect the stress and the pressure drop. As shown in figure 7, there are further boundary layers which are important for inflow, and we discuss these in appendix A.

2.5.4 Summary

The main part of the leading-order inflow may be summarized as follows. The leading-order axial flow and pressure are unchanged until they reach the interface, after which they transition to plug flow inside the obstacle in region III, whose length is comparable to the channel height. This has an effect on the correction to normal flow and pressure in region II. Meanwhile, the tangential velocity satisfies the no tangential slip condition on the obstacle surface in region IIa – a boundary layer of even smaller width, which is comparable to the length of the channel height multiplied by the reciprocal of the Reynolds number.

Combining the original asymptotic expansion (2.26) with the asymptotic expansion in region II, namely (2.36c ), we deduce that the pressure jump between outer inflow regions I and IV is $O({\it\epsilon})$ . Specifically, we have determined that

(2.61) $$\begin{eqnarray}\displaystyle {\it\Pi}^{-}\sim \frac{{\it\alpha}}{K}\quad \text{as }{\it\epsilon},{\it\delta}\rightarrow 0~\text{with }{\it\epsilon}\ll {\it\delta}\ll 1, & & \displaystyle\end{eqnarray}$$

where ${\it\alpha}\approx -0.117$ .

The argument presented in this section will not hold if the fluid is leaving the obstacle rather than entering it. This is because the boundary layer region IIa is only valid for inflow (as can be seen from (2.55)). We reconsider our approach for outflow in the next section.

2.5.5 Outflow

The outflow inner problem (within regions V and VI) is shown in figure 5. We again seek a solution as ${\it\delta}\rightarrow 0$ with ${\it\epsilon}\ll {\it\delta}\ll 1$ using the method of matched asymptotic expansions.

The asymptotic structure for outflow is shown schematically in figure 11. The flow in region V is unchanged at leading order and reaches the interface as uniform plug flow with magnitude $1$ . There is an inviscid core in the centre of region VI, flanked by Prandtl boundary layers near $z=0,1$ whose thicknesses are of $O(({\it\delta}X^{+})^{1/2})$ . These boundary layers (region VIb) grow downstream, and the transition to Poiseuille flow occurs when their thicknesses become of $O(1)$ . This occurs on a horizontal length scale of $O(1/{\it\delta})$ in a transition region we denote region VIa. There is a third boundary layer, denoted region VIc, which occurs in the corner between the channel wall and the porous interface. This region deals with a singularity at the corner that arises in region VIb, as in the classic Prandtl problem.

Figure 11. Schematic diagram of the inner regions for outflow – VIa, VIb, and VIc. The size of the boundary region VIb has been exaggerated for illustrative purposes, is of height $O({\it\epsilon}{\it\delta}^{1/2})$ , and grows proportional to the square root of the distance from the porous obstacle. Region VIc lies within region VIb and has height and length $O({\it\epsilon}{\it\delta})$ .

In this section, we investigate the regions VI and VIa to understand the transition from plug to Poiseuille flow and to calculate the pressure jump between the outer regions IV and VII. Regions VIb and VIc are discussed in appendix B.1 since they are not important to the transition. Recall that, from (2.61), we calculated an $O({\it\epsilon})$ pressure jump between outer regions I and IV for inflow. We show that the pressure jump between outer regions IV and VII is of $O({\it\epsilon}/{\it\delta})$ for outflow.

We start by investigating how the fluid moves through region VI. The leading-order behaviour is plug flow, and we determine the effect of the interfacial no tangential slip condition on the correction to the plug flow. In particular, we determine the far-field behaviour of the flow, which allows us to impose the correct matching conditions for the problem of transition to Poiseuille flow in region VIa.

2.5.6 Inner regions V and VI

To determine the outflow pressure jump (2.35b ), we must consider first-order flow terms. Therefore, we pose the asymptotic expansions

(2.62a-d ) $$\begin{eqnarray}\overline{\boldsymbol{u}}_{0}\sim \boldsymbol{e}_{x}+{\it\delta}^{1/2}\overline{\boldsymbol{u}}_{0a},\quad \overline{\boldsymbol{Q}}_{0}\sim \boldsymbol{e}_{x}+{\it\delta}^{1/2}\overline{\boldsymbol{Q}}_{0a},\quad \overline{p}_{1}\sim {\it\delta}^{-1/2}\overline{p}_{1a},\quad \overline{P}_{1}\sim -\frac{X^{+}}{K}+{\it\delta}^{1/2}\overline{P}_{1b},\end{eqnarray}$$

as ${\it\delta}\rightarrow 0$ , where the non-integer powers of ${\it\delta}$ arise due to the correction terms from the boundary layers in region VIb near the channel walls. We determine the leading-order behaviour within the coupled regions V and VI, but do not fully solve for the higher-order terms due to their complexity. Instead, we analyse the higher-order problems to obtain the information required for the matching condition with the transition region VIa.

Substituting the asymptotic expansions (2.62) into the system outlined in figure 5, the leading-order equations are trivially satisfied. The $O({\it\delta}^{1/2})$ exterior flow equations are

(2.63a-c ) $$\begin{eqnarray}\overline{u}_{0aX^{+}}=-\overline{p}_{1aX^{+}},\quad \overline{w}_{0aX^{+}}=-\overline{p}_{1az},\quad 0=\overline{u}_{0aX^{+}}+\overline{w}_{0az}\quad \text{in region VI},\end{eqnarray}$$

and the $O({\it\delta}^{1/2})$ interior flow equations are

(2.64a,b ) $$\begin{eqnarray}\overline{\boldsymbol{Q}}_{0a}=-K\boldsymbol{{\rm\nabla}}\overline{P}_{1b},\quad 0=\boldsymbol{{\rm\nabla}}\boldsymbol{\cdot }\overline{\boldsymbol{Q}}_{0a}\quad \text{in region V}.\end{eqnarray}$$

In region V, the first correction to the boundary condition on the channel walls is given by

(2.65) $$\begin{eqnarray}\displaystyle \overline{W}_{0a}=0\quad \text{on }z=0,1\text{ for }X^{+}<0, & & \displaystyle\end{eqnarray}$$

while the far-field conditions are obtained by matching with region IV, to yield

(2.66a,b ) $$\begin{eqnarray}\overline{P}_{1b}\rightarrow 0,\quad \overline{\boldsymbol{Q}}_{0a}\rightarrow \mathbf{0}\quad \text{as }X^{+}\rightarrow -\infty ~\text{for }0<z<1.\end{eqnarray}$$

In region VI, the boundary conditions on the channel walls are obtained by matching with region VIb. We show in appendix B.1 that, at this order, the flow problem within region VIb is equivalent to the standard Prandtl problem of uniform flow past a flat plate (Prandtl Reference Prandtl and Krazer1905). The matching conditions are therefore given by

(2.67a ) $$\begin{eqnarray}\displaystyle & \displaystyle \overline{w}_{0a}=\tilde{{\it\beta}}(4X^{+})^{-1/2}\quad \text{on }z=0\text{ for }X^{+}>0, & \displaystyle\end{eqnarray}$$

(2.67b ) $$\begin{eqnarray}\displaystyle & \displaystyle \overline{w}_{0a}=-\tilde{{\it\beta}}(4X^{+})^{-1/2}\quad \text{on }z=1\text{ for }X^{+}>0, & \displaystyle\end{eqnarray}$$

$\tilde{{\it\beta}}$

${\approx}\,$

${\it\beta}$

$\tilde{{\it\beta}}$

$\tilde{{\it\beta}}=2^{1/2}{\it\beta}$

$X^{+}=0$

(2.68a-c ) $$\begin{eqnarray}\overline{u}_{0a}=\overline{U}_{0a},\quad \overline{p}_{1a}=0,\quad \overline{w}_{0a}=0\quad \text{on }X^{+}=0\text{ for }0<z<1.\end{eqnarray}$$

We do not solve the coupled system of equations (2.64)–(2.68) here, though we note that it is possible to obtain an implicit representation for $\overline{w}_{0a}$ via a Fourier transformation. Instead, we acquire the necessary information for the matching with region VIa, which allows us to determine ${\it\Pi}_{0}^{+}$ in the next section. We note that (2.63) can be rearranged to deduce the following three results: firstly, that

(2.69) $$\begin{eqnarray}\displaystyle {\rm\nabla}^{2}\overline{w}_{0a}=0\quad \text{in region VI}; & & \displaystyle\end{eqnarray}$$

secondly, that $\overline{p}_{1a}$ is the harmonic conjugate of $\overline{w}_{0a}$ ; and, thirdly, that $\overline{p}_{1a}+\overline{u}_{1a}$ is independent of $X^{+}$ . We deduce the leading-order far-field behaviour of $\overline{w}_{0a}$ by considering the governing equation (2.69) together with the channel wall boundary conditions (2.67) as $X^{+}\rightarrow \infty$ , to obtain the expansion

(2.70) $$\begin{eqnarray}\displaystyle \overline{w}_{0a}\sim \frac{\tilde{{\it\beta}}(1-2z)}{(4X^{+})^{1/2}}\quad \text{as }X^{+}\rightarrow \infty ~\text{for }0<z<1. & & \displaystyle\end{eqnarray}$$

Using the far-field condition (2.70) in the governing equations (2.63a,c ), we further deduce that

(2.71a,b ) $$\begin{eqnarray}\overline{u}_{0a}\sim \tilde{{\it\beta}}(4X^{+})^{1/2},\quad \overline{p}_{1a}\sim -\tilde{{\it\beta}}(4X^{+})^{1/2}\quad \text{as }X^{+}\rightarrow \infty ~\text{for }0<z<1.\end{eqnarray}$$

We note that the far-field expansions (2.70)–(2.71) have no dependence on the permeability $K$ . As these far-field expansions are used to form matching conditions, which transfer information into the transition region VIa (which is where we calculate the pressure drop ${\it\Pi}_{0}^{+}$ ), we note that ${\it\Pi}_{0}^{+}$ will be independent of $K$ . We carry out this matching in appendix B.2.

The components of stress on the interface at $x=1$ are given by

(2.72a ) $$\begin{eqnarray}\displaystyle & \displaystyle {\it\sigma}_{xx}=-\hspace{-3.0pt}\left(12-\frac{2}{K}+{\it\epsilon}{\it\delta}^{-1/2}\overline{p}_{1a}\right)+\mathit{o}({\it\epsilon}{\it\delta}^{-1/2}), & \displaystyle\end{eqnarray}$$

(2.72b ) $$\begin{eqnarray}\displaystyle & \displaystyle {\it\sigma}_{xz}={\it\epsilon}{\it\delta}^{1/2}\left(\overline{u}_{0az}+\overline{w}_{0aX^{+}}\right)+\mathit{o}({\it\epsilon}{\it\delta}^{1/2}) & \displaystyle\end{eqnarray}$$

${\it\epsilon},{\it\delta}\rightarrow 0$

${\it\epsilon}\ll {\it\delta}\ll 1$

$\overline{p}_{1a}$

$\overline{u}_{0a}$

$\overline{w}_{0a}$

$z$

$O({\it\delta}^{-1/2})$

$O({\it\delta}^{1/2})$

2.5.7 Transition region VIa

We now consider the region in which the flow transitions from plug to Poiseuille flow. This occurs when the growing Prandtl boundary layers with thickness of $O(({\it\delta}X^{+})^{1/2})$ become of $O(1)$ , that is, when $X^{+}=O(1/{\it\delta})$ , bringing about a balance in the exterior flow equations between fluid momentum and viscosity across the full channel width. Recall that we have called this region VIa, as illustrated in figure 11.

The relevant scalings are $X^{+}={\it\delta}^{-1}\widetilde{X}$ and $\overline{w}_{0}={\it\delta}\widetilde{w}_{0}$ , and we introduce the new variables $\overline{u}_{0}=\widetilde{u}_{0}$ and $\overline{p}_{1}=\widetilde{p}_{1}$ to signify that we are working in region VIa. Then, in region VIa (where $0<\widetilde{X}<\infty$ and $0<z<1$ ) the exterior flow equations given in figure 5 become

(2.73a ) $$\begin{eqnarray}\displaystyle & \displaystyle \widetilde{u}_{0}\widetilde{u}_{0\widetilde{X}}+\widetilde{w}_{0}\widetilde{u}_{0z}=-{\it\delta}\widetilde{p}_{1\widetilde{X}}+{\it\delta}^{2}\widetilde{u}_{0\widetilde{X}\widetilde{X}}+\widetilde{u}_{0zz}, & \displaystyle\end{eqnarray}$$

(2.73b ) $$\begin{eqnarray}\displaystyle & \displaystyle {\it\delta}^{2}(\widetilde{u}_{0}\widetilde{u}_{0\widetilde{X}}+\widetilde{w}_{0}\widetilde{w}_{0z})=-{\it\delta}\widetilde{p}_{1z}+{\it\delta}^{4}\widetilde{w}_{0\widetilde{X}\widetilde{X}}+{\it\delta}^{2}\widetilde{w}_{0zz}, & \displaystyle\end{eqnarray}$$

(2.73c ) $$\begin{eqnarray}\displaystyle & \displaystyle \widetilde{u}_{0\widetilde{X}}+\widetilde{w}_{0z}=0. & \displaystyle\end{eqnarray}$$

We pose the asymptotic expansions

(2.74a-c ) $$\begin{eqnarray}\widetilde{u}_{0}=\widetilde{u}_{00}+O({\it\delta}^{1/2}),\quad \widetilde{w}_{0}=\widetilde{w}_{01}+O({\it\delta}^{1/2}),\quad \widetilde{p}_{1}={\it\delta}^{-1}\widetilde{p}_{10}+O({\it\delta}^{-1/2}),\end{eqnarray}$$

as ${\it\delta}\rightarrow 0$ , and substitute them into the governing equations (2.73) to yield the following leading-order governing equations in region VIa:

(2.75a-c ) $$\begin{eqnarray}\widetilde{u}_{00}\widetilde{u}_{00\widetilde{X}}+\widetilde{w}_{01}\widetilde{u}_{00z}=-\widetilde{p}_{10\widetilde{X}}+\widetilde{u}_{00zz},\quad 0=-\widetilde{p}_{10z},\quad \widetilde{u}_{00\widetilde{X}}+\widetilde{w}_{01z}=0.\end{eqnarray}$$

The boundary conditions on the channel walls become

(2.76) $$\begin{eqnarray}\displaystyle \overline{\boldsymbol{u}}_{00}=\mathbf{0}\quad \text{on }z=0,1\text{ for}~\widetilde{X}>0. & & \displaystyle\end{eqnarray}$$

The leading-order matching conditions are derived in appendix B.2. They are obtained by forming composite expansions between regions VI and VIb, then using Van Dyke’s matching rule (Van Dyke Reference Van Dyke1975). We use the multiplicative composite expansions (Van Dyke Reference Van Dyke1975) to ensure that the matching condition for the velocity satisfies the no-slip condition on the channel walls. This allows a greater accuracy when we solve numerically the resulting system. The matching conditions are given by

(2.77a ) $$\begin{eqnarray}\displaystyle & \displaystyle \widetilde{u}_{00}\sim (1+2\tilde{{\it\beta}}\widetilde{X}^{1/2})f^{\prime }(\widetilde{{\it\eta}}_{(0)})f^{\prime }(\widetilde{{\it\eta}}_{(1)}), & \displaystyle\end{eqnarray}$$

(2.77b ) $$\begin{eqnarray}\displaystyle & \displaystyle \widetilde{w}_{01}\sim \frac{(1-2z)}{2\tilde{{\it\beta}}\widetilde{X}^{1/2}}(\widetilde{{\it\eta}}_{(0)}f^{\prime }(\widetilde{{\it\eta}}_{(0)})-f(\widetilde{{\it\eta}}_{(0)}))(f(\widetilde{{\it\eta}}_{(1)})-\widetilde{{\it\eta}}_{(1)}f^{\prime }(\widetilde{{\it\eta}}_{(1)})), & \displaystyle\end{eqnarray}$$

(2.77c ) $$\begin{eqnarray}\displaystyle & \displaystyle \widetilde{p}_{10}\sim -2\tilde{{\it\beta}}\widetilde{X}^{1/2}, & \displaystyle\end{eqnarray}$$

$\widetilde{X}\downarrow 0$

$f$

$\widetilde{{\it\eta}}_{(0)}=z/\widetilde{X}^{1/2}$

$\widetilde{{\it\eta}}_{(1)}=(1-z)/\widetilde{X}^{1/2}$

$\widetilde{X}$

$z$

(2.78a-c ) $$\begin{eqnarray}\widetilde{u}_{00}\sim 6z(1-z),\quad \widetilde{w}_{01}\rightarrow 0,\quad \widetilde{p}_{10}\sim -12\widetilde{X}+{\it\Pi}_{0}^{+}\quad \text{as}~\widetilde{X}\rightarrow \infty ~\text{for }0<z<1,\end{eqnarray}$$

where ${\it\Pi}_{0}^{+}$ is determined from the solution in region VIa. The full system for the flow in region VIa is given by (2.75)–(2.78).

In Bodoia & Osterle (Reference Bodoia and Osterle1961), this same transition region was considered, though was not derived formally through asymptotic methods. The authors used the boundary conditions

(2.79a-c ) $$\begin{eqnarray}\widetilde{u}_{00}=1,\quad \widetilde{w}_{01}=0,\quad \widetilde{p}_{10}=0\quad \text{on}~\widetilde{X}=0\text{ for }0<z<1,\end{eqnarray}$$

instead of the matching conditions (2.77). We note that, unlike (2.77b ), the condition (2.79b ) is not consistent with the small $\widetilde{X}$ approximation to the governing equations (2.75), from which it can be deduced that $\widetilde{w}_{01}$ is singular at the corners at $\widetilde{X}=0$ , $z=0,1$ . With the boundary conditions (2.79), it was found in Bodoia & Osterle (Reference Bodoia and Osterle1961) that the pressure jump ${\it\Pi}_{0}^{+}\approx -0.338$ .

We now calculate ${\it\Pi}_{0}^{+}$ numerically using the formal matching conditions (2.77). We start the numerical simulation at $\widetilde{X}=10^{-4}$ to avoid the inverse square root singularities in $\widetilde{w}_{01}$ and in $\widetilde{p}_{10\widetilde{X}}$ , and use the finite-difference scheme outlined in Bodoia & Osterle (Reference Bodoia and Osterle1961). In figure 12, we show the numerical solutions for $\widetilde{p}_{10}+12\widetilde{X}$ as a function of $\widetilde{X}$ , and we emphasise that $\widetilde{p}_{10}$ is independent of $z$ . From the matching condition (2.78c ), we determine that the value of the leading-order outflow pressure jump between outer regions IV and VII is given by ${\it\Pi}_{0}^{+}=-0.331$ to 3 significant figures. We note that the value of ${\it\Pi}^{+}$ found by Bodoia & Osterle (Reference Bodoia and Osterle1961) using their inconsistent boundary conditions (2.79) only differs by approximately 2 % from the value of ${\it\Pi}^{+}$ obtained using formal matching conditions.

Figure 12. Calculating the outflow pressure drop along the channel by plotting $\widetilde{p}_{10}+12\widetilde{X}$ . See text for further details.

Combining the original asymptotic expansion (2.26) with the asymptotic expansion in region VIa (2.74c ), we deduce that the pressure jump between outer outflow regions IV and VII is $O({\it\epsilon}/{\it\delta})$ . Specifically, we have determined that

(2.80) $$\begin{eqnarray}\displaystyle {\it\Pi}^{+}\sim \frac{{\it\Pi}_{0}^{+}}{{\it\delta}}\quad \text{as}~{\it\epsilon},{\it\delta}\rightarrow 0~\text{with}~{\it\epsilon}\ll {\it\delta}\ll 1, & & \displaystyle\end{eqnarray}$$

where ${\it\Pi}_{0}^{+}\approx -0.331$ .

3.6.1 Inflow

The inflow behaviour is similar to the flow considered in § 2.5.1, for which the asymptotic structure is given in figure 7. We show that the leading-order tangential velocity is independent of $X^{-}$ through region II and region III, but the next-order flow depends on $X^{-}$ in region II. In this section we also show that the no tangential slip condition (3.29) is satisfied in region IIa, and we calculate the asymptotic behaviour of ${\it\Lambda}^{s-}$ as ${\it\delta}_{s}\rightarrow 0$ .

In region II, the tangential flow is unchanged at leading order and satisfies the asymptotic expansion

(3.38) $$\begin{eqnarray}\displaystyle \overline{v}_{0}(X^{-},z;s,t)=6z(1-z)+{\it\delta}_{s}\overline{v}_{01}+O({\it\delta}_{s}^{4/3})\quad \text{as}~{\it\delta}_{s}\rightarrow 0. & & \displaystyle\end{eqnarray}$$

Substituting (3.38) into the governing equation (3.27) and equating terms of $O(1)$ , gives

(3.39) $$\begin{eqnarray}\displaystyle z(1-z)\overline{v}_{01X^{-}}+(1-2z)\overline{w}_{01}=0, & & \displaystyle\end{eqnarray}$$

with the matching condition

(3.40) $$\begin{eqnarray}\displaystyle \overline{v}_{01}\rightarrow 0\quad \text{as }X^{-}\rightarrow -\infty . & & \displaystyle\end{eqnarray}$$

Equation (3.39), with far-field condition (3.40), can be integrated numerically with respect to $X^{-}$ to determine $\overline{v}_{01}$ . The resulting correction to the tangential velocity induces a slip on the cell walls, which is remedied by a boundary layer in region IIb (as illustrated in figure 7), the details of which are omitted from this paper for the sake of brevity.

The coordinate scaling in region IIa is $X^{-}={\it\delta}_{s}\widehat{X}$ , and we use hats to denote all variables in this region. We use the flow variables calculated in § 2.5.3, and pose an asymptotic expansion for the tangential velocity of the form

(3.41) $$\begin{eqnarray}\displaystyle \widehat{v}_{0}(\widehat{X},z;s,t)=\widehat{v}_{00}(\widehat{X},z;s,t)+O({\it\delta}_{s}), & & \displaystyle\end{eqnarray}$$

where the leading-order term satisfies the governing equation

(3.42) $$\begin{eqnarray}\displaystyle 6z(1-z)\widehat{v}_{00\widehat{X}}=\widehat{v}_{00\widehat{X}\widehat{X}}\quad \text{for}~\widehat{X}<0,\,0<z<1, & & \displaystyle\end{eqnarray}$$

with the no tangential slip condition (3.29) yielding

(3.43) $$\begin{eqnarray}\displaystyle \widehat{v}_{00}=0\quad \text{on}~\widehat{X}=0,\;0<z<1, & & \displaystyle\end{eqnarray}$$

and a matching condition with region II given by

(3.44) $$\begin{eqnarray}\displaystyle \widehat{v}_{00}(\widehat{X},z;s,t)\sim 6z(1-z)\quad \text{as}~\widehat{X}\rightarrow -\infty ,\text{ for }0<z<1. & & \displaystyle\end{eqnarray}$$

The system (3.42)–(3.44) is solved by

(3.45) $$\begin{eqnarray}\displaystyle \widehat{v}_{00}=6z(1-z)(1-\exp (6z(1-z)\widehat{X})). & & \displaystyle\end{eqnarray}$$

The leading-order contribution to ${\it\Lambda}^{s-}$ in (3.32a ) is of $O({\it\delta}_{s})$ , and comes from both the $O({\it\delta}_{s})$ term in region II and the $O(1)$ term in region IIa, as region IIa is $O({\it\delta}_{s})$ smaller than region II. Therefore, substituting the asymptotic expansion (3.38) into (3.32a ), and using the governing equation (3.39) to relate the velocity $\overline{v}_{01}$ to $\overline{w}_{01}$ (which has already been calculated numerically in § 2.5.2), along with the leading-order solution for the tangential velocity in region IIa (3.45), the leading-order contribution to (3.32a ) becomes

(3.46) $$\begin{eqnarray}\displaystyle {\it\Lambda}^{s-}\sim \frac{1}{\overline{Re}}\frac{\partial }{\partial s}\!\left(\!\frac{p_{0s}^{s}(0,s,t)}{{\it\nu}(s,t)}\left(\!1-\int _{-\infty }^{0}\int _{0}^{1}\left(\!\frac{2z-1}{z(1-z)}\int _{-\infty }^{X^{-}}\overline{w}_{01}({\it\xi},z;s,t)\,\text{d}{\it\xi}\!\right)\!\,\text{d}z\,\text{d}X^{-}\!\right)\!\right)\!, & & \displaystyle \nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

as ${\it\epsilon},1/\overline{Re}\rightarrow 0$ with ${\it\epsilon}\ll 1/\overline{Re}\ll K\ll 1$ . We note that the integral is finite because $\overline{w}_{01}\rightarrow 0$ as $z\rightarrow 0,1$ .

In contrast to the two-dimensional case, the pressure now varies in the $s$ -direction. Therefore, the component of the interfacial stress in the normal direction is different from the two-dimensional case (2.60a ), and in three dimensions is given by

(3.47) $$\begin{eqnarray}\displaystyle {\it\sigma}_{nn}=-(P_{0}^{s}(0,s,t)+{\it\epsilon}(\overline{Re}P_{10}^{s}(0,s,t)+\overline{p}_{11}^{s}(0,s,z,t)))+O({\it\epsilon}/\overline{Re}^{1/3}) & & \displaystyle\end{eqnarray}$$

as ${\it\epsilon},1/\overline{Re}\rightarrow 0$ with ${\it\epsilon}\ll 1/\overline{Re}\ll K\ll 1$ . We note that $P_{0}^{s}$ and $P_{10}^{s}$ have no $z$ -dependence. Therefore, as in two dimensions, the first order at which the normal stress varies in $z$ is at $O({\it\epsilon})$ . Finally, we determine the leading-order term in the component of the shear stress on the obstacle boundary in the $s$ -direction for inflow, ${\it\sigma}_{ns}$ . From the tangential flow solution in region IIa given by (3.45), and recalling the scaling (3.26), we deduce that

(3.48) $$\begin{eqnarray}\displaystyle {\it\sigma}_{ns}\sim \frac{{\it\epsilon}}{{\it\delta}_{s}}\frac{p_{0s}^{s}(0,s,t)}{12}\widehat{v}_{0\widehat{X}}(0,z;s,t)=-{\it\epsilon}\overline{Re}(3{\it\nu}(s,t)p_{0s}^{s}(0,s,t)(z(1-z))^{2}), & & \displaystyle\end{eqnarray}$$

as ${\it\epsilon},1/\overline{Re}\rightarrow 0$ with ${\it\epsilon}\ll 1/\overline{Re}\ll K\ll 1$ .

3.6.2 Outflow

Outflow occurs in a similar manner to § 2.5.5, for which the asymptotic structure is illustrated in figure 11. The tangential velocity within region V is plug flow (as given by (3.25)) and the tangential flow in region VI is of $O({\it\delta}_{s})$ . The tangential flow then transitions to Poiseuille flow in region VIa, where it becomes of $O(1)$ . In this section we solve for the leading-order tangential flow for outflow and calculate the asymptotic behaviour of ${\it\Lambda}^{s+}$ as ${\it\delta}_{s}\rightarrow 0$ .

In region VI, we use the flow variables calculated in § 2.5.6 as well as the asymptotic expansion

(3.49) $$\begin{eqnarray}\displaystyle \overline{v}_{0}(X^{+},z;s,t)={\it\delta}_{s}(12X^{+})+O({\it\delta}_{s}^{3/2})\quad \text{as}~{\it\delta}_{s}\rightarrow 0, & & \displaystyle\end{eqnarray}$$

which satisfies the no tangential slip condition (3.29). The no-slip condition on the cell walls is satisfied in region VIb (as illustrated in figure 11) via a similarity solution, as described in appendix B.1. As with the three-dimensional inflow case considered in the previous section, the component of the interfacial stress in the normal direction is now

(3.50) $$\begin{eqnarray}\displaystyle {\it\sigma}_{nn}=-(P_{0}^{s}(0,s,t)+{\it\epsilon}\overline{Re}P_{10}^{s}(0,s,t))+O({\it\epsilon}/{\it\delta}_{s}^{1/2}), & & \displaystyle\end{eqnarray}$$

as ${\it\epsilon},1/\overline{Re}\rightarrow 0$ with ${\it\epsilon}\ll 1/\overline{Re}\ll 1$ , and where $P_{0}^{s}$ and $P_{10}^{s}$ have no $z$ -dependence. Therefore, as in two dimensions, the first order at which the normal stress varies in $z$ is at $O({\it\epsilon}/{\it\delta}_{s}^{1/2})$ . From (3.49) and the scaling (3.26), we deduce that the stress on the interface in the tangential direction is

(3.51) $$\begin{eqnarray}\displaystyle {\it\sigma}_{ns}\sim -{\it\epsilon}\frac{p_{0s}^{s}(0,s,t)}{12}\widehat{v}_{0X^{+}}(0,z;s,t)=\frac{{\it\epsilon}}{{\it\nu}(s,t)\overline{Re}}p_{0s}^{s}(0,s,t), & & \displaystyle\end{eqnarray}$$

as ${\it\epsilon},1/\overline{Re}\rightarrow 0$ with ${\it\epsilon}\ll 1/\overline{Re}\ll K\ll 1$ .

The transition to Poiseuille flow occurs in region VIa, within which the relevant coordinate scaling is $X^{+}={\it\delta}_{s}^{-1}\widetilde{X}$ . We further use $\overline{v}_{0}=\widetilde{v}_{0}$ to indicate that we are in this intermediate region, and use the variables introduced in § 2.5.7. We substitute the asymptotic expansion

(3.52) $$\begin{eqnarray}\displaystyle \widetilde{v}_{0}=\widetilde{v}_{00}+O({\it\delta}_{s}^{1/2})\quad \text{as}~{\it\delta}_{s}\rightarrow 0, & & \displaystyle\end{eqnarray}$$

into the governing equation (3.27) to obtain the leading-order governing equation

(3.53) $$\begin{eqnarray}\displaystyle \widetilde{u}_{00}\widetilde{v}_{00\widetilde{X}}+\widetilde{w}_{01}\widetilde{v}_{00z}=12+\widetilde{v}_{00zz}. & & \displaystyle\end{eqnarray}$$

The cell wall boundary conditions (3.28) become

(3.54) $$\begin{eqnarray}\displaystyle \widetilde{v}_{00}=0\quad \text{on }z=0,1\text{ for}~\widetilde{X}>0, & & \displaystyle\end{eqnarray}$$

the matching condition with region VI is

(3.55) $$\begin{eqnarray}\displaystyle \widetilde{v}_{00}\rightarrow 0\quad \text{as}~\widetilde{X}\downarrow 0\text{ for }0<z<1, & & \displaystyle\end{eqnarray}$$

and the far-field behaviour is determined via the matching condition (3.30), and is given by

(3.56) $$\begin{eqnarray}\displaystyle \widetilde{v}_{00}\sim 6z(1-z)\quad \text{as}~\widetilde{X}\rightarrow \infty ~\text{for }0<z<1. & & \displaystyle\end{eqnarray}$$

The system (3.53)–(3.56) details the transition from plug to Poiseuille tangential flow. To determine the flux jump ${\it\Lambda}^{s+}$ at leading order, the system must be solved numerically. We note that the system (3.53)–(3.56) is independent of $K$ , and hence only needs to be solved once to determine ${\it\Lambda}^{s+}$ at leading order.

Accounting for the change of variable $X^{+}={\it\delta}_{s}^{-1}\widetilde{X}$ being dependent on $s$ , we find

(3.57) $$\begin{eqnarray}\displaystyle \frac{\partial }{\partial s}\mapsto \frac{\partial }{\partial s}-\frac{\widetilde{X}}{{\it\nu}}\frac{\text{d}{\it\nu}}{\text{d}s}\frac{\partial }{\partial \widetilde{X}}, & & \displaystyle\end{eqnarray}$$

and it follows that the leading-order behaviour of ${\it\Lambda}^{s+}$ is given by

(3.58a ) $$\begin{eqnarray}\displaystyle & \displaystyle {\it\Lambda}^{s+}\sim -\overline{Re}{\it\Lambda}_{a}\frac{\partial }{\partial s}({\it\nu}p_{0s}^{s}(0,s,t)), & \displaystyle\end{eqnarray}$$

where

(3.58b ) $$\begin{eqnarray}\displaystyle & \displaystyle {\it\Lambda}_{a}=\int _{0}^{\infty }\int _{0}^{1}(1-\widetilde{v}_{00}(\widetilde{X},z))\,\text{d}z\,\text{d}\widetilde{X}. & \displaystyle\end{eqnarray}$$

We calculate $\widetilde{v}_{00}$ numerically using a second-order central finite-difference scheme, then determine ${\it\Lambda}_{a}$ using the trapezium rule. We find that ${\it\Lambda}_{a}\approx 0.125$ to three significant figures.

Combining the outflow flux jump (3.58) with its inflow equivalent (3.46), and recalling that ${\it\nu}=KP_{0n}^{s}(0,s,t)$ , we find that the leading-order (in ${\it\delta}$ ) coupling condition (3.16b ) in the large Reynolds number asymptotic sublimit is given by

(3.59a ) $$\begin{eqnarray}\displaystyle \frac{\partial }{\partial n}\left(p_{10}^{s}-\frac{1}{10}p_{0t}^{s}+\frac{3}{560}|\boldsymbol{{\rm\nabla}}p_{0}^{s}|^{2}\right)-12K\frac{\partial P_{10}^{s}}{\partial n}=-K{\it\Lambda}_{a}\frac{\partial }{\partial s}\left(P_{0n}^{s}p_{0s}^{s}\right)H(\boldsymbol{u}_{0}\boldsymbol{\cdot }\boldsymbol{n})\quad \text{on}~\partial D. & & \displaystyle \nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

The leading-order pressure jump condition (3.16a ) in the large Reynolds number sublimit is given by

(3.59b ) $$\begin{eqnarray}\displaystyle & \displaystyle p_{10}^{s}-P_{10}^{s}=(KP_{0n}^{s})^{2}{\it\Pi}_{0}^{+}H(\boldsymbol{u}_{0}\boldsymbol{\cdot }\boldsymbol{n})\quad \text{on}~\partial D. & \displaystyle\end{eqnarray}$$

We see that the coefficient premultiplying the Heaviside function in (3.59a ) does not necessarily vanish when $\boldsymbol{u}_{0}\boldsymbol{\cdot }\boldsymbol{n}=0$ , whereas the coefficient premultiplying the Heaviside function in (3.59b ) does vanish. Thus, at these points the pressure is continuous but the pressure derivative along the boundary will have a $\log$ singularity.