Appendix A. Solution to fluid flow problem

We describe fluid flow through the lumen and ECS by the steady Stokes equations and continuity equation for an incompressible viscous Newtonian fluid of constant dynamic viscosity, $\mu ^{*}$. In the membrane and cell layer, we model fluid flow using Darcy's law for an incompressible Newtonian fluid, where $\phi _i, k^*_i$ denote the porosity and permeability of porous region $i = m,c$. We let $k^*_m = k^*_m(z^*)$ vary across the length of the domain. We take advantage of the long, thin bioreactor geometry and non-dimensionalise the fluid flow problem using the following scalings:

(A1)\begin{equation} \left.\begin{array}{c} (r^ * , z^ * ) = L^ * (\varepsilon r, z), \quad (U^ * _l, U^ * _e, W^ * _l, W^ * _e ) = \bar{W}^ * (\varepsilon U_l, \varepsilon U_e, W_l, W_e), \\ (U^ * _m, U^ * _c, W^ * _m, W^ * _c ) = \varepsilon\bar{W}^ * (U_m, U_c, \varepsilon W_m, \varepsilon W_c), \quad P^ * _i = \dfrac{\mu^{ * }\bar{W}^ * }{\varepsilon^2L^ * }P_i, \quad Q^ * _{i,in} = 2\pi\varepsilon^2{L^ * }^2\bar{W}^ * Q_{i,in}, \end{array}\right\} \end{equation}

where $\varepsilon = R^*_{l0}/L^*$ is the aspect ratio of the lumen, $\bar {W}^*$ is the typical axial lumen fluid velocity and $i=l,m,c,e$. We define the cell layer thickness $\delta = \varepsilon (R_{c0} - R_{m0})$, and consider the physically relevant sublimit $0 < \delta \ll \varepsilon \ll 1$. Standard lubrication scalings are applied in the lumen and ECS (Reference BatchelorBatchelor 1967), and lubrication scalings for the coupled porous regions are used in the membrane and ECS (Reference Dalwadi, King, Dyson and ArkillDalwadi et al. 2020). The scaled Stokes equations and Darcy equations in the limit $\varepsilon \to 0$ are given by

(A2a)\begin{gather} \frac{1}{r}\frac{\partial}{\partial r}(r U_i) + \frac{\partial W_i}{\partial z} = 0, \end{gather}

(A2b)\begin{gather}\frac{\partial P_i}{\partial r} = 0, \quad \frac{\partial P_i}{\partial z} = \frac{1}{r}\frac{\partial}{\partial r}\left(r\frac{\partial W_i}{\partial r}\right), \end{gather}

for $i = l,e$, and

(A3a)\begin{gather} \frac{1}{r}\frac{\partial}{\partial r}(r U_i) = 0, \end{gather}

(A3b)\begin{gather}\phi_i U_i = {-}\kappa_i\frac{\partial P_i}{\partial r}, \quad \phi_i W_i = {-}\kappa_i\frac{\partial P_i}{\partial z}, \end{gather}

for $i = m,c$. Here, $\kappa _i = k^*_i/\varepsilon ^4L^{*2}$ is the non-dimensional permeability in region $i$. We take $\kappa _i$ to be a $O (1)$ quantity in our analysis to investigate the most general model set-up i.e. the distinguished asymptotic limit in the language of matched asymptotic expansions.

At the lumen–membrane, membrane–cell and cell–ECS interfaces, we prescribe interfacial conditions to couple the fluid flow in each region. Typically, flow past porous media can be modelled by imposing a Beavers–Joseph slip condition on the tangential component of fluid stress, where the slip coefficient is experimentally derived (Reference Beavers and JosephBeavers & Joseph 1967). However, Shipley et al. show that accounting for slip at a porous interface in a typical bioreactor set-up has negligible impact on the flow solution (Reference Shipley, Waters and EllisShipley et al. 2010). Motivated by this, we impose a no axial slip condition and a continuity of normal flux condition at each interface. Secondly, we specify that the normal component of fluid stress in the lumen is balanced by the fluid pore pressure in the membrane at the lumen–membrane interface similarly to Reference Chapman, Whiteley, Byrne, Waters and ShipleyChapman et al. (2017). Additionally, we specify that the normal component of fluid stress in the ECS is balanced by the fluid pore pressure in the cell layer at the cell–ECS interface, and we prescribe continuity of fluid pressure at the membrane–cell layer interface. The scaled interfacial conditions are then given by

(A4a)\begin{gather} -P_l + 2\varepsilon^2\frac{\partial U_l}{\partial r} = {-}P_m, \quad U_l = \phi_mU_m, \quad W_l = 0 \quad \mbox{at} \ r = 1, \end{gather}

(A4b)\begin{gather}P_m = P_c, \quad \phi_mU_m = \phi_cU_c, \quad W_m,W_c = 0 \quad \mbox{at} \ r = R_{m0}, \end{gather}

(A4c)\begin{gather}-P_c = {-}P_e + 2\varepsilon^2\frac{\partial U_e}{\partial r}, \quad \phi_cU_c = U_e, \quad W_e = 0 \quad \mbox{at} \ r = R_{c0}. \end{gather}

Given that we are investigating the limit in which cell layer thickness is small (i.e. $\delta = \varepsilon (R_{c0} - R_{m0}) \ll \varepsilon \ll 1$), we can further reduce our flow problem. From inspection of (A3) and (A4) for $i=c$, we see that in the limit $\delta \to 0$ (i.e. $R_{c0} \to R_{m0}$)

(A5a,b)\begin{equation} P_c|_{R_{c0}} \to P_c|_{R_{m0}}, \quad U_c|_{R_{c0}} \to U_c|_{R_{m0}}. \end{equation}

Utilising (A5a,b) and taking the limit $\varepsilon \to 0$, we can replace (A4) with the following boundary conditions:

(A6a)\begin{gather} P_l = P_m, \quad U_l = \phi_mU_m, \quad W_l = 0 \quad \mbox{at} \ r = 1, \end{gather}

(A6b)\begin{gather}P_m = P_e, \quad \phi_mU_m = U_e, \quad W_e = 0 \quad \mbox{at} \ r = R_{m0}. \end{gather}

Therefore, our flow problem reduces to flow through the lumen, membrane and ECS only, and the cell layer is replaced by an interfacial condition at the membrane–ECS interface (A6b).

The membrane is sealed at each axial end, therefore we prescribe no flux through these ends in the membrane

(A7)\begin{equation} W_m = 0 \quad \mbox{at} \ z = 0,1. \end{equation}

To drive flow along the lumen and ECS, we prescribe a pressure drop along each channel

(A8a)\begin{gather} P_i = P_{i,in} \quad \mbox{at} \ z = 0, \end{gather}

(A8b)\begin{gather}P_i = P_{i,out} \quad \mbox{at} \ z = 1, \end{gather}

for $i = l,e$. We note that, while we prescribe fluid pressures at the inlet and outlet of the lumen and ECS in our model, we determine the appropriate inlet pressures to prescribe in (A8) from a set of volumetric flow rates from literature summarised in table 1. We relate the inlet flux and inlet pressures via the relationship

(A9)\begin{equation} Q_{i,in} = \int_{R_1}^{R_2} r W_i \, \mbox{d} r \quad \mbox{at} \ z = 0, \end{equation}

for $i = l,e$, where $R_1, R_2$ define the inner and outer radial boundaries of the respective region.

Lastly, we prescribe a symmetry condition at the lumen centreline

(A10)\begin{equation} U_l = \frac{\partial W_l}{\partial r} = 0 \quad \mbox{at} \ r = 0, \end{equation}

and at the outer wall of the bioreactor, we prescribe no slip and no flux

(A11)\begin{equation} \boldsymbol{U}_e = \boldsymbol{0} \quad \mbox{at} \ r = R_{e0}. \end{equation}

Integrating the lumen and ECS flow equations (A2) with respect to $r$ and applying the no axial slip boundary conditions in (A6), (A11), we find that $P_l = P_l(z), P_e = P_e(z)$ and

(A12a)\begin{gather} W_l = \frac{1}{4}\frac{\mbox{d}P_l}{\mbox{d}z}(r^{2} - 1), \quad U_l = \frac{r}{16}\frac{\mbox{d}^{2}P_l}{\mbox{d}z^{2}}(2 - r^{2}), \end{gather}

(A12b)\begin{gather}W_e = \frac{1}{4}\frac{\mbox{d}P_e}{\mbox{d}z}\left(r^{2}-R_{e0}^{2} + \frac{(R_{e0}^{2}-R_{m0}^{2})\ln (r/R_{e0})}{\ln (R_{m0}/R_{e0})}\right), \end{gather}

(A12c)\begin{align} U_e & = \frac{1}{16r\ln (R_{m0}/R_{e0})}\frac{\mbox{d}^{2}P_e}{\mbox{d}z^{2}} \Big(2r^{2}(R_{m0}^{2} - R_{e0}^{2})\ln (r/R_{e0})\nonumber\\ &\quad -(r^{2}-R_{e0}^{2})(R_{m0}^{2} - R_{e0}^{2} + (r^{2} - R_{e0}^{2})\ln (R_{m0}/R_{e0}))). \end{align}

Substituting (A3b) into (A3a) for $i=m$ and solving subject to continuity of pressure in (A6) we find that

(A13)\begin{equation} P_m = \frac{(P_e-P_l)}{\ln R_{m0}}\ln r + P_l. \end{equation}

We substitute (A13) into (A3b) to find

(A14a,b)\begin{equation} \phi_mU_m = \frac{\kappa(P_l-P_e)}{r\ln R_{m0}}, \quad \phi_mW_m = \kappa\left(\frac{\ln r}{\ln R_{m0}}\left(\frac{\mbox{d}P_l}{\mbox{d}z} - \frac{\mbox{d}P_e}{\mbox{d}z}\right) - \frac{\mbox{d}P_l}{\mbox{d}z}\right). \end{equation}

Substituting into the continuity of flux conditions in (A6), we end up with a system of coupled ODEs for $P_l$ and $P_e$

(A15a,b)\begin{equation} \frac{A}{\kappa}\frac{\mbox{d}^{2}P_l}{\mbox{d}z^{2}} = P_l - P_e, \quad \frac{B}{\kappa}\frac{\mbox{d}^{2}P_e}{\mbox{d}z^{2}} = P_l - P_e, \end{equation}

where

(A16a,b)\begin{equation} A = \frac{\ln R_{m0}}{16}, \quad B = {-}\frac{\ln R_{m0}}{16\ln (R_{m0}/R_{e0})}[(R_{e0}^{2}-R_{m0}^{2})^{2} + (R_{e0}^4-R_{m0}^4)\ln (R_{m0}/R_{e0})]. \end{equation}

A.1 Case I: homogeneous permeability

Substituting the ODE for $P_l$ in to the ODE for $P_e$ in (A14a,b), we obtain a fourth-order ODE in $P_l$ with constant coefficients, the characteristic equation of which is solved by $\lambda _{1,2} = \pm \sqrt {\kappa (B-A)/BA}, \lambda _{3,4}=0$. Noting that $\kappa (B-A)/BA > 0$ for all possible values of $R_{m0}$, $R_{e0}$ and $\kappa$, we can substitute our general solution for $P_l$ back in to (A14a,b) and solve subject to (A8) to find

(A17a)\begin{align} P_l & = \frac{B(P_{l,in} - P_{e,in})}{(B-A)({\rm e}^{\lambda_1} - {\rm e}^{\lambda_2})}[(z-1)({\rm e}^{\lambda_1} - {\rm e}^{\lambda_2}) - ({\rm e}^{\lambda_1(z-1)} - {\rm e}^{\lambda_2(z-1)})] \nonumber\\ &\quad - \frac{B(P_{l,out} - P_{e,out})}{(B-A)({\rm e}^{\lambda_1} - {\rm e}^{\lambda_2})}[z({\rm e}^{\lambda_1} - {\rm e}^{\lambda_2}) - ({\rm e}^{\lambda_1z} - {\rm e}^{\lambda_2z})] + (P_{l,out} - P_{l,in})z + P_{l,in},\end{align}

(A17b)\begin{align} P_e & = \frac{(P_{l,in} - P_{e,in})}{({\rm e}^{\lambda_1} - {\rm e}^{\lambda_2})}\left[\frac{B}{(B-A)}(z-1)({\rm e}^{\lambda_1} - {\rm e}^{\lambda_2}) - \frac{A}{(B-A)}({\rm e}^{\lambda_1(z-1)} - {\rm e}^{\lambda_2(z-1)})\right] \nonumber\\ &\quad - \frac{(P_{l,out} - P_{e,out})}{({\rm e}^{\lambda_1} - {\rm e}^{\lambda_2})}\left[\frac{B}{(B-A)}z({\rm e}^{\lambda_1} - {\rm e}^{\lambda_2}) - \frac{A}{(B-A)}({\rm e}^{\lambda_1z} - {\rm e}^{\lambda_2z})\right]\nonumber\\ &\quad + (P_{l,out} - P_{l,in})z + P_{l,in}, \end{align}

giving us a fully analytic solution for $P_l$ and $P_e$.

A.2 Case II: uniform flux

A common objective of bioreactor systems is to deliver nutrient uniformly to cells. For a homogeneous cell population with uniform seeding density at the outer wall of our fibre membrane, this corresponds to an axially uniform nutrient flux across the membrane wall. For nutrient delivered via advection, this can be achieved with a radial membrane flux that is independent of axial position at the outer wall.

From (A14a,b), we see that the axial dependency in membrane flux comes from the difference in pressures in the lumen and ECS. Axial pressure gradients drive fluid through our bioreactor system and thus cannot be eliminated. Instead, we seek a solution to our membrane permeability such that

(A18)\begin{equation} \kappa(z) = \frac{1}{P_l(z) - P_e(z)}, \end{equation}

where $P_l, P_e$ are yet to be solved for. Therefore (A14a,b) becomes

(A19)\begin{equation} \phi_mU_m = \frac{1}{r\ln R_{m0}}. \end{equation}

The flux delivered to the outer membrane wall where cells are seeded can then be determined by evaluating (A19) at $r = R_{m0}$.

With (A18), the system (A14a,b) decouples, simplifying to

(A20a,b)\begin{equation} \frac{\mbox{d} ^2P_l}{\mbox{d} z^2} = \frac{1}{A}, \quad \frac{\mbox{d} ^2P_e}{\mbox{d} z^2} = \frac{1}{B}. \end{equation}

Solving subject to inlet and outlet pressures (A8), (A20a,b) becomes

(A21a)\begin{gather} P_l = \frac{z(z-1)}{2A} + (P_{l,out} - P_{l,in})z + P_{l,in}, \end{gather}

(A21b)\begin{gather}P_e = \frac{z(z-1)}{2B} + (P_{e,out} - P_{e,in})z + P_{e,in}. \end{gather}

From (A18) it is clear that we encounter a singularity when $P_l = P_e$. That is, we cannot have uniform flow if the pressure drop across the membrane vanishes at any axial position. Additionally we cannot have negative permeability, therefore this analysis only holds for $P_l > P_e$ at each value of $z$. For this to hold, we must also ensure $P_{l,in} > P_{e,in}$ and $P_{l,out} > P_{e,out}$.

Appendix B. Cell layer reduction

As justified by the parameter values in table 1, we investigate the physically relevant asymptotic limit $\delta \ll \varepsilon \ll 1$. We transform the radial coordinate $r$ such that $r = n + \varepsilon (R_{m0} + R_{c0})/2$. In doing so, we take the axial centreline of the cell layer to be $n = 0$, where $n \in (-\delta /2, \delta /2)$ in the cell layer. A schematic is shown in figure 8, where we let $\boldsymbol {Q}_i = \mathcal{P}_l c _i\boldsymbol {U}_i - \phi _i\boldsymbol {\nabla } c_i$ for $i =m,c$, and $\boldsymbol {Q}_e = \mathcal{P}_l c _e\boldsymbol {U}_e - \boldsymbol {\nabla } c_e$. Given the cell layer is $O (\delta )$, we introduce an inner scaling $n = \delta N$ such that $N = {O}(1)$. For a distinguished asymptotic limit where uptake appears at leading order, we also take $\mathcal {R} = {O}(\delta ^{-1})$, and define $\bar {\mathcal {R}} = {O}(1)$ such that $\mathcal {R} = \bar {\mathcal {R}}/\delta$. We then pose an asymptotic expansion for $c_i$ in powers of $\delta$ i.e.

(B1)\begin{equation} c_i \sim {c_i}_0 + \delta{c_i}_1 + \cdots \quad \mbox{as} \ \delta \rightarrow 0. \end{equation}

Substituting (B1) into (2.2) and (2.3) and collecting $O (\delta ^{-2})$ terms, we obtain

(B2)\begin{equation} \frac{\partial^2 {c_m}_0}{\partial N^2} = \frac{\partial^2 {c_c}_0}{\partial N^2} = \frac{\partial^2 {c_e}_0}{\partial N^2} = 0. \end{equation}

Equation (B2) has general solution

(B3a,b)\begin{equation} {c_c}_0 = A_c(z)N + B_c(z), \quad {c_i}_0 = A_i(z)N + B_i(z), \quad i = m,e. \end{equation}

We note that in the cell layer, $N = {O}(1)$ and bounded. In the membrane and ECS, however, $N$ can tend to $\pm \infty$. Since concentrations in the outer regions will be bounded as $\delta \to 0$, matching these inner solutions for ${c_m}_0, {c_e}_0$ as $N \to \pm \infty$ with the outer solutions as $n \to 0$, we obtain the matching condition $A_i \equiv 0$ for $i = m,e$, so that

(B4a,b)\begin{equation} {c_m}_0 = B_m(z), \quad {c_e}_0 = B_e(z). \end{equation}

From the continuity of concentration and flux boundary conditions (2.5), we then find $A_c(z) = 0$, and $B_m(z) = B_c(z) = B_e(z)$. Therefore, concentration $c_c$ is constant at leading order in the boundary layer, and therefore

(B5)\begin{equation} {c_c}_0 = {c_m}_0|_{n = 0} = {c_e}_0|_{n = 0}. \end{equation}

Collecting terms at $O (\delta ^{-1})$, (2.2) and (2.3) become

(B6a)\begin{gather} \frac{\partial}{\partial N}\left(\mathcal{P}_c(\boldsymbol{U}_c\boldsymbol{\cdot} \hat{\boldsymbol{n}}){c_c}_0 - \frac{\partial {c_c}_1}{\partial N}\right) = {-}\bar{\mathcal{R}}({c_c}_0), \end{gather}

(B6b)\begin{gather}\frac{\partial}{\partial N}\left(\mathcal{P}_i(\boldsymbol{U}_i\boldsymbol{\cdot} \hat{\boldsymbol{n}}){c_i}_0 - \frac{\partial {c_i}_1}{\partial N}\right) = 0, \end{gather}

for $i = m,e$, where $\hat {\boldsymbol {n}}$ is the unit normal vector in the positive $N$-direction. From (B6) we see that the flux is conserved in the inner regions of the membrane and cell layer domains. Multiplying (B6a) by the square of the porosity of the cell layer, $\phi _c^2$, and noting that $\mathcal {P}_i = \mathcal {P}_l/\phi _i$, we integrate over $N \in (-1/2,1/2)$ to find

(B7)\begin{equation} \phi_c\left.\left(\mathcal{P}_l(\boldsymbol{U}_c\boldsymbol{\cdot} \hat{\boldsymbol{n}}) {c_c}_0 - \phi_c\frac{\partial {c_c}_1}{\partial N}\right) \right|_{-{1}/{2}}^{{1}/{2}} = {-}\phi_c^2\bar{\mathcal{R}}({c_c}_0). \end{equation}

Noting continuity of fluid flux in (2.5), and that $\mathcal {P}_i = \mathcal {P}_l/\phi _i$, we see

(B8a)\begin{gather} \left.\phi_c\left(\mathcal{P}_l(\boldsymbol{U}_c\boldsymbol{\cdot} \hat{\boldsymbol{n}}){c_c}_0 - \phi_c\frac{\partial {c_c}_1}{\partial N}\right)\right|_{-{1}/{2}} = \left.\phi_m\left(\mathcal{P}_l(\boldsymbol{U}_m\boldsymbol{\cdot} \hat{\boldsymbol{n}}){c_m}_0 - \phi_m\frac{\partial {c_m}_1}{\partial N}\right)\right|_{-{1}/{2}}, \end{gather}

(B8b)\begin{gather}\left.\phi_c\left(\mathcal{P}_l(\boldsymbol{U}_c\boldsymbol{\cdot} \hat{\boldsymbol{n}}){c_c}_0 - \phi_c\frac{\partial {c_c}_1}{\partial N}\right)\right|_{{1}/{2}} = \left.\left(\mathcal{P}_l(\boldsymbol{U}_e\boldsymbol{\cdot} \hat{\boldsymbol{n}}){c_e}_0 - \frac{\partial {c_e}_1}{\partial N}\right)\right|_{{1}/{2}}. \end{gather}

We also note that, by continuity of fluid flux and continuity of concentration, the advective terms in (B8) cancel out. Matching into the outer problem in the membrane and ECS, we find that

(B9a)\begin{gather} \left.\frac{\partial {c_m}}{\partial n}\right|_{n\uparrow 0} \sim \frac{1}{\delta}\frac{\partial}{\partial N}({c_m}_0 + \delta {c_m}_1)|_{N = {-}{1}/{2}} = \left.\frac{\partial {c_m}_1}{\partial N}\right|_{N = {-}{1}/{2}}, \end{gather}

(B9b)\begin{gather}\left.\frac{\partial {c_e}}{\partial n}\right|_{n\downarrow 0} \sim \frac{1}{\delta}\frac{\partial}{\partial N}({c_e}_0 + \delta {c_e}_1)|_{N = {1}/{2}} = \left.\frac{\partial {c_e}_1}{\partial N}\right|_{N = {1}/{2}}. \end{gather}

Therefore, using the results of (B6) and (B9), we may derive the following effective boundary condition:

(B10)\begin{equation} \phi_m^2\frac{\partial {c_m}}{\partial n} - \frac{\partial {c_e}}{\partial n} = {-}\phi_c^2\bar{\mathcal{R}}({c_c}_0). \end{equation}

Figure 8.

Cell layer with thickness $\delta$. The central dotted line denotes the centreline along which we define $N=0$, and the external dotted lines above and below denote the asymptotic extent of the inner region in the outside membrane and ECS regions which we scale into.

Transforming back into the radial coordinate $r$, (B5), (B10) become

(B11)\begin{equation} c_m = c_e, \quad \phi_m^2\frac{\partial c_m}{\partial r} - \frac{\partial c_e}{\partial r} = {-}\phi_c^2\bar{\mathcal{R}}({c_m}) \quad \mbox{at} \ r = \varepsilon R_{m0}, \end{equation}

thus reducing our problem to fluid flow and nutrient transport in the lumen, membrane and ECS only. Formally, (B11) replaces (2.5b)–(2.5c) as $R_{c0} \rightarrow R_{m0}$. Physically, (B11) states that the effective jump in diffusive flux across the membrane wall is due to nutrient uptake by the cells. Similarly to the derivation of (B11), the effective boundary condition at the membrane-ECS interface for the waste transport problem becomes

(B12)\begin{equation} w_m = w_e, \quad \phi_m^2\frac{\partial w_m}{\partial r} = \frac{\partial w_e}{\partial r} + \nu\phi_c^2\bar{\mathcal{R}}(c_e) \quad \mbox{at} \ r = \varepsilon R_{m0}. \end{equation}

Appendix C. Streamfunction–arclength formulation

To mathematically reduce the dimensionality of the two-dimensional advection term in our nutrient transport model, we consider the transformation of our system from a three-dimensional Cartesian coordinate system to a streamfunction–arclength formulation defined in axisymmetric cylindrical coordinates, $(x, y, z) \mapsto (\psi, \theta, s)$. However, given that our solutions to the fluid flow problem are known in an axisymmetric cylindrical coordinate system, it is prudent to express our Cartesian coordinate system in terms of cylindrical coordinates, such that

(C1a–c)\begin{equation} x = r(\psi, s)\cos{\theta}, \quad y = r(\psi, s)\sin{\theta}, \quad z = z(\psi, s). \end{equation}

The metrics of the transformation can be expressed as

(C2a)\begin{gather} H_\psi \equiv \left|\frac{\partial\boldsymbol{x}}{\partial\psi}\right| = \sqrt{x_\psi^2 + y_\psi^2 + z_\psi^2} = \sqrt{r_\psi^2(\cos^2\theta + \sin^2\theta) + z_\psi^2} = \left|\frac{\partial\boldsymbol{r}}{\partial\psi}\right|, \end{gather}

(C2b)\begin{gather}H_\theta \equiv \left|\frac{\partial\boldsymbol{x}}{\partial\theta}\right| = r\sqrt{\cos^2\theta + \sin^2\theta} = r, \end{gather}

(C2c)\begin{gather}H_s \equiv \left|\frac{\partial\boldsymbol{x}}{\partial s}\right| = \sqrt{x_s^2 + y_s^2 + z_s^2} = \sqrt{r_s^2(\cos^2\theta + \sin^2\theta) + z_s^2} = \left|\frac{\partial\boldsymbol{r}}{\partial s}\right| , \end{gather}

where $\boldsymbol {x} = (x,y,z)$ and $\boldsymbol {r} = (r,\theta, z)$.

The streamfunction ($\psi$) relates to the fluid velocity in domain $i$ by

(C3)\begin{equation} \boldsymbol{U}_i = \boldsymbol{\nabla} \times \boldsymbol{\psi}, \end{equation}

where $\boldsymbol {\psi } = (0,0,\psi )$. Having previously solved for components of the fluid velocity in cylindrical coordinates, (C3) can be written as

(C4a,b)\begin{equation} U_{i} = {-}\frac{1}{r}\frac{\partial\psi}{\partial z}, \quad W_{i} = \frac{1}{r}\frac{\partial\psi}{\partial r}, \end{equation}

where $U_i, W_i$ are the radial and axial components of fluid velocity, respectively. Therefore

(C5)\begin{equation} \left|\frac{\partial\psi}{\partial \boldsymbol{r}}\right| = \sqrt{\left(\frac{\partial\psi}{\partial r}\right)^2 + \left(\frac{\partial\psi}{\partial z}\right)^2} = r\sqrt{U_i^2 + W_i^2} = r\tilde{U}_i. \end{equation}

We can relate (C5) to the metric of transformation $H_\psi$ by

(C6)\begin{equation} \frac{\partial\psi}{\partial \boldsymbol{r}}\boldsymbol{\cdot} \frac{\partial \boldsymbol{r}}{\partial\psi} = \left|\frac{\partial\psi}{\partial \boldsymbol{r}}\right|\left|\frac{\partial \boldsymbol{r}}{\partial\psi}\right|\cos\theta. \end{equation}

Given that both vectors are tangential and align in the direction of constant $s$, this means that $\theta =0$, and thus $\cos \theta = 1$. The left hand side of (C6) can be simplified if we consider the Jacobian of the transformation $(r,z) \mapsto (\psi, s)$

(C7)\begin{equation} \left|\frac{\partial(r,z)}{\partial(\psi,s)}\right| = \begin{pmatrix} r_\psi & r_s \\ z_\psi & z_s \end{pmatrix} = \begin{pmatrix} \psi_r & \psi_z \\ s_r & s_z \end{pmatrix}^{{-}1} = \frac{1}{\psi_r s_z - \psi_z s_r} \begin{pmatrix} s_z & -\psi_z \\ -s_r & \psi_r \end{pmatrix}. \end{equation}

Therefore, the left-hand side of (C6) becomes

(C8)\begin{equation} \frac{\partial\psi}{\partial \boldsymbol{r}}\boldsymbol{\cdot} \frac{\partial \boldsymbol{r}}{\partial\psi} = \psi_r r_\psi + \psi_z z_\psi = \frac{\psi_r s_z - \psi_z s_r}{\psi_r s_z - \psi_z s_r} = 1, \end{equation}

as a direct result of the inverse function theorem. Therefore,

(C9a)\begin{gather} \left|\frac{\partial\psi}{\partial \boldsymbol{r}}\right|\left|\frac{\partial \boldsymbol{r}}{\partial\psi}\right| = 1, \end{gather}

(C9b)\begin{gather}\Rightarrow \left|\frac{\partial \boldsymbol{r}}{\partial\psi}\right| = \left|\frac{\partial\psi}{\partial \boldsymbol{r}}\right|^{{-}1} = \frac{1}{r\tilde{U}}. \end{gather}

The normalised arclength ($s$) parameterises the distance along streamlines. Without loss of generality, we impose $0 \leq s \leq 1$ in the lumen, $1 \leq s \leq 2$ in the membrane and $2 \leq s \leq 3$ in the ECS. We let $\mathcal {L}_i$ denote the actual arclength of a streamline in domain $i$, which varies between 0 and the total length of the streamline, $l_{\psi,i}$. The arclength $\mathcal {L}_{i}$ and normalised arclength $s$ are therefore related by

(C10)\begin{equation} \mathcal{L}_{i} = l_{\psi,i}(s - n_i), \end{equation}

where $n_i = 0,1,2$ for the lumen, membrane and ECS, respectively. We also know that, along a streamline of constant $\psi$, we can relate the arclength of a streamline to axisymmetric coordinates such that

(C11)\begin{equation} {\mbox{d}\mathcal{L}_i}^{2} = {\mbox{d} r}^2 + {\mbox{d} z}^2. \end{equation}

With $\psi$ fixed, $\mathcal {L}_i$ varies only with $s$. Combining (C10) and (C11), this means that

(C12a)\begin{gather} \left(\frac{\partial\mathcal{L}_i}{\partial s}\right)^2 = \left(\frac{\partial r}{\partial s}\right)^2 + \left(\frac{\partial z}{\partial s}\right)^2 = l_{\psi,i}^2, \end{gather}

(C12b)\begin{gather}\Rightarrow \left|\frac{\partial\boldsymbol{r}}{\partial s}\right| = \sqrt{\left(\frac{\partial r}{\partial s}\right)^2 + \left(\frac{\partial z}{\partial s}\right)^2} = l_{\psi,i}. \end{gather}

Substituting (C9), (C12) in to (C2), we can determine the metrics of the transformation to be

(C13a–c)\begin{equation} H_\psi = \frac{1}{r\tilde{U}_i} , \quad H_\theta = r, \quad H_s = l_{\psi,i}. \end{equation}

In these transformed coordinates, (2.2) becomes

(C14)\begin{equation} \mathcal{P}_i\frac{\tilde{U}_i}{H_s}\frac{\partial c_i}{\partial s} = \frac{1}{H_\psi H_\theta H_s}\left[\frac{\partial}{\partial\psi}\left(\frac{H_\theta H_s}{H_\psi}\frac{\partial c_i}{\partial \psi}\right) + \frac{\partial}{\partial s}\left(\frac{H_\psi H_\theta}{H_s}\frac{\partial c_i}{\partial s}\right)\right], \quad i = l,e.\end{equation}

Finally, substituting in our metrics of transformation (C13a–c) into (C14), we recover (3.3).

To transform the boundary conditions from axisymmetric to streamfunction–arclength formulation, we make use of the fact that

(C15a)\begin{gather} \frac{\partial}{\partial r} = \frac{\partial\psi}{\partial r}\frac{\partial}{\partial\psi} + \frac{\partial s}{\partial r}\frac{\partial}{\partial s} = rW_i\frac{\partial}{\partial\psi} + \frac{U_i}{l_{\psi,i}\tilde{U}_i}\frac{\partial}{\partial s}, \end{gather}

(C15b)\begin{gather}\frac{\partial}{\partial z} = \frac{\partial\psi}{\partial z}\frac{\partial}{\partial\psi} + \frac{\partial s}{\partial z}\frac{\partial}{\partial s} = {-}rU_i\frac{\partial}{\partial\psi} + \frac{W_i}{l_{\psi,i}\tilde{U}_i}\frac{\partial}{\partial s}. \end{gather}

With this, continuity of concentration and solute flux (2.5a) becomes

(C16a)\begin{equation} c_l = c_m, \end{equation}

(C16b)\begin{align} &rW_l\frac{\partial c_l}{\partial\psi} + \frac{U_l}{l_{\psi,l}\tilde{U}_l}\frac{\partial c_l}{\partial s} = \phi_m^2\left(rW_m\frac{\partial c_m}{\partial\psi} + \frac{U_m}{l_{\psi,m}\tilde{U}_m}\frac{\partial c_m}{\partial s}\right)\nonumber\\ &\quad \mbox{at} \ s = 1, \psi \in [\psi_{l,crit},\psi_{l,max}), \end{align}

and the effective cell layer condition (B11) becomes

(C17a)\begin{equation} c_m = c_e, \end{equation}

(C17b)\begin{align} &\phi_m^2\left(rW_m\frac{\partial c_m}{\partial\psi} + \frac{U_m}{l_{\psi,m}\tilde{U}_m}\frac{\partial c_m}{\partial s}\right) = \left(rW_e\frac{\partial c_e}{\partial\psi} + \frac{U_e}{l_{\psi,e}\tilde{U}_e}\frac{\partial c_e}{\partial s}\right) - \phi_c^2\bar{\mathcal{R}}(c_e) \nonumber\\ &\quad \mbox{at}\ s = 2, \psi \in (\psi_{e,min},\psi_{e,crit}]. \end{align}

The symmetry condition (2.6) and the no-flux condition at the outer wall (2.7) become

(C18)\begin{equation} rW_i\frac{\partial c_i}{\partial\psi} + \frac{U_i}{l_{\psi,i}\tilde{U}_i}\frac{\partial c_i}{\partial s} = 0 \quad \mbox{at} \ \psi = 0, \psi_{e,max}, \end{equation}

and the no-diffusive-flux outlet conditions (2.10) become

(C19)\begin{equation} \frac{W_i}{l_{\psi,i}\tilde{U}_i}\frac{\partial c_i}{\partial s} -rU_i\frac{\partial c_i}{\partial\psi} = 0 \quad \mbox{at} \ s = 1, \psi \in (0, \psi_{l,crit}] \quad \mbox{and} \quad s = 3. \end{equation}

The inlet concentration conditions (2.9) remain the same,

(C20a)\begin{gather} c_l = 1 \quad \mbox{at} \ s = 0, \end{gather}

(C20b)\begin{gather}c_e = 0 \quad \mbox{at} \ s = 2, \psi \in (\psi_{e,crit},\psi_{e,max}]. \end{gather}

For waste metabolite transport, the governing equations (3.3) remain the same, swapping $c_i$ for $w_i$. The jump condition (B12) becomes

(C21a)\begin{equation} w_m = w_e, \end{equation}

(C21b)\begin{align} &\phi_m^2\left(rW_m\frac{\partial w_m}{\partial\psi} + \frac{U_m}{l_{\psi,m}\tilde{U}_m}\frac{\partial w_m}{\partial s}\right) = \left(rW_e\frac{\partial w_e}{\partial\psi} + \frac{U_e}{l_{\psi,e}\tilde{U}_e}\frac{\partial w_e}{\partial s}\right) + \nu\phi_c^2\bar{\mathcal{R}}(w_e) \nonumber\\ &\quad \mbox{at}\ s = 2, \psi \in (\psi_{e,min},\psi_{e,crit}], \end{align}

and the inlet concentration (2.12) condition becomes

(C22a)\begin{gather} w_l = 0 \quad \mbox{at} \ s = 0, \end{gather}

(C22b)\begin{gather}w_e = 0 \quad \mbox{at} \ s = 2, \psi \in (\psi_{e,crit},\psi_{e,max}]. \end{gather}

The remaining boundary conditions (C16), (C18) and (C19) remain the same as for nutrient transport, exchanging $c_i$, for $w_i$.

Appendix D. Model reduction of metabolite transport problem

Motivated by the scaling arguments outlined in § 3.3.1 and dropping the superscript $^\star$, (3.3) in the lumen and ECS becomes

(D1)\begin{equation} \varepsilon^2\mathcal{P}_i\frac{\partial c_i}{\partial s} = \frac{\partial}{\partial\psi}\left({r}^2 l_{\psi,i}\tilde{U}_i\frac{\partial c_i}{\partial \psi}\right) + \varepsilon^2\frac{\partial}{\partial s}\left(\frac{1}{l_{\psi,i}\tilde{U}_i}\frac{\partial c_i}{\partial s}\right), \end{equation}

for $i=l,e$ where

(D2)\begin{equation} \tilde{U}_i = \sqrt{\varepsilon^2U_i^2 + W_i^2}, \end{equation}

and the reduced Péclet number $\varepsilon ^{2}\mathcal {P}_i = {O}(1)$. Utilising the scalings outlined in § 3.3.2, (3.3) in the membrane becomes

(D3)\begin{equation} \varepsilon^2\mathcal{P}_m\frac{\partial c_m}{\partial s} = \varepsilon^2\frac{\partial}{\partial\psi}\left({r}^2 l_{\psi,m}\tilde{U}_m\frac{\partial c_m}{\partial \psi}\right) + \frac{\partial}{\partial s}\left(\frac{1}{l_{\psi,m}\tilde{U}_m}\frac{\partial c_m}{\partial s}\right), \end{equation}

where the reduced Péclet number in the membrane $\varepsilon ^{2}\mathcal {P}_m = \varepsilon ^{2}\mathcal {P}_l/\phi _m = {O}(1)$. Utilising the appropriate scalings from §§ 3.3.1 and 3.3.2, continuity of concentration and solute flux becomes

(D4a)\begin{equation} c_l = c_m, \end{equation}

(D4b)\begin{align} &r W_l\frac{\partial c_l}{\partial\psi} + \varepsilon^2\frac{U_l}{l_{\psi,l}\tilde{U}_l}\frac{\partial c_l}{\partial s} = \phi_m^2\left(\varepsilon^2r W_m\frac{\partial c_m}{\partial\psi} + \frac{U_m}{l_{\psi,m}\tilde{U}_m}\frac{\partial c_m}{\partial s}\right) \nonumber\\ &\quad \mbox{at} \ s = 1, \psi \in [\psi_{l,crit},\psi_{l,max}), \end{align}

and the effective cell layer condition becomes

(D5a)\begin{equation} c_m = c_e, \end{equation}

(D5b)\begin{align} &\phi_m\left(\varepsilon^2 r W_m\frac{\partial c_m}{\partial\psi} + \frac{U_m}{l_{\psi,m}\tilde{U}_m}\frac{\partial c_m}{\partial s}\right) = \left(r W_e\frac{\partial c_e}{\partial\psi} + \varepsilon^2\frac{U_e}{l_{\psi,e}\tilde{U}_e}\frac{\partial c_e}{\partial s}\right) - \varepsilon \phi_c^2\bar{\mathcal{R}}(c_e) \nonumber\\ &\quad \mbox{at}\ s = 2, \psi \in (\psi_{e,min},\psi_{e,crit}]. \end{align}

The symmetry condition and the no-flux condition at the outer wall (C18) become

(D6)\begin{equation} r W_i\frac{\partial c_i}{\partial\psi} + \varepsilon^2\frac{U_i}{l_{\psi,i}\tilde{U}_i}\frac{\partial c_i}{\partial s} = 0 \quad \mbox{at} \ \psi = 0, \psi_{e,max}, \end{equation}

and the no-diffusive-flux outlet conditions become

(D7)\begin{equation} \frac{W_i}{l_{\psi,i}\tilde{U}_i}\frac{\partial c_i}{\partial s} -r U_i \frac{\partial c_i}{\partial\psi} = 0 \quad \mbox{at} \ s = 1, \psi \in (0, \psi_{l,crit}] \quad \mbox{and} \quad s = 3. \end{equation}

The inlet concentration conditions remain the same

(D8a)\begin{gather} c_l = 1 \quad \mbox{at} \ s = 0, \end{gather}

(D8b)\begin{gather}c_e = 0 \quad \mbox{at} \ s = 2, \psi \in (\psi_{e,crit},\psi_{e,max}]. \end{gather}

Again, (D1)–(D4), (D6) and (D7) remain the same for waste metabolite transport exchanging $c_i$ for $w_i$. The effective boundary condition (C21) becomes

(D9a)\begin{equation} w_m = w_e, \end{equation}

(D9b)\begin{align} &\phi_m^2\left(\varepsilon^2 r W_m\frac{\partial w_m}{\partial\psi} + \frac{U_m}{l_{\psi,m}\tilde{U}_m}\frac{\partial w_m}{\partial s}\right) = \left(r W_e\frac{\partial w_e}{\partial\psi} + \varepsilon^2\frac{U_e}{l_{\psi,e}\tilde{U}_e}\frac{\partial w_e}{\partial s}\right) - \varepsilon \phi_c^2\bar{\mathcal{R}}(w_e) \nonumber\\ &\quad \mbox{at}\ s = 2, \psi \in (\psi_{e,min},\psi_{e,crit}], \end{align}

and the inlet concentration condition remains the same as (C22).

Appendix E. Solution to metabolite transport problem

E.1 Membrane solution

From (3.5), we see that $\tilde {U}_m = U_m$ at leading order in $\varepsilon$. Therefore at leading order and dropping superscript $\star$, (D3) becomes (3.6). Equation (3.6) is then solved subject to continuity of concentration (D4a) and the leading-order effective cell layer boundary condition (D5b), which becomes

(E1)\begin{equation} \frac{\phi_m^2}{l_{\psi,m}}\frac{\partial c_m}{\partial s} = {-} \varepsilon \phi_c^2\bar{\mathcal{R}}(c_m) \quad \mbox{at}\ s = 2, \psi \in (\psi_{e,min},\psi_{e,crit}], \end{equation}

noting that $W_e$ vanishes at the membrane–ECS interface as per (A4), and continuity of concentration (D5a) means we can consider $c_m = c_e$ at this interface too. For generality, we retain $\varepsilon \bar {\mathcal {R}}$ as an $O (1)$ quantity. Having previously solved the fluid transport problem analytically in Appendix A to obtain $U_m$, we know that streamlines in the membrane align with the radial axis as demonstrated in figure 3. Therefore, we can relate the radial coordinate to the normalised arclength in the membrane via the relationship

(E2)\begin{equation} r = (s-1)(R_{m0} - 1) + 1, \end{equation}

and streamline lengths throughout the membrane become,

(E3)\begin{equation} l_{\psi,m} = R_{m0} - 1. \end{equation}

On this basis, we can express the fluid velocity along a streamline in the membrane wall by the radial velocity in (A3b), having transformed the radial coordinate by (E2) to give

(E4)\begin{equation} \tilde{U}_m = U_m = \frac{\kappa(P_l - P_e)}{\phi_m[(s-1)(R_{m0}-1) + 1]\ln R_{m0}}. \end{equation}

Substituting (E2)–(E4) into (3.6), we can integrate subject to continuity of concentration (D4a) and the effective membrane-ECS interfacial condition (E1) to determine that

(E5)\begin{equation} c_m = c_l|_{[l,m]} - \mathcal{A}([(s-1)(R_{m0}-1) + 1]^{-\xi} - 1), \quad \xi = \frac{\varepsilon^2\mathcal{P}_l\kappa(P_e - P_l)}{\phi_m^2\ln R_{m0}}, \end{equation}

where $c_l|_{[l,m]}$ denotes the unknown lumen concentration at the lumen–membrane interface and $\mathcal {A}$ depends on $c_l|_{[l,m]}$. These unknowns depend on the solution to nutrient concentration in the lumen which in turn is coupled to the solution to the membrane concentration itself. Therefore, care must be taken when handling interfacial conditions coupling each region together. To understand this coupling, we must now consider nutrient concentration in the lumen and ECS.

E.2 Lumen/ECS solution

E.2.1 Boundary layer analysis

Equation (D1) results in a singular perturbation problem as $\varepsilon \to 0$, and implies the existence of a boundary layer at the end of the streamlines in the lumen and ECS where diffusion along streamlines dominates. Therefore, care must be taken in deriving an appropriate matching condition for the outer region.

We first investigate the boundary layer at the lumen and ECS exits. To maintain a dominant balance of advection and diffusion along streamlines, $s$ must scale with $\varepsilon ^2$. We introduce an $O (1)$ variable $\hat {s}$ and scale such that $s = n + \varepsilon ^2 \hat {s}$ where $n=1,3$ in the lumen and ECS, respectively. Equation (D1) then becomes

(E6)\begin{equation} \varepsilon^2\mathcal{P}_i\frac{\partial \hat{c}_i}{\partial \hat{s}} = \varepsilon^2\frac{\partial}{\partial\psi}\left({r}^2 l_{\psi,i}\tilde{U}_i\frac{\partial \hat{c}_i}{\partial \psi}\right) + \frac{\partial}{\partial \hat{s}}\left(\frac{1}{l_{\psi,i}\tilde{U}_i}\frac{\partial \hat{c}_i}{\partial \hat{s}}\right), \end{equation}

subject to

(E7)\begin{equation} \frac{W_i}{l_{\psi,i}\tilde{U}}_i\frac{\partial \hat{c}_i}{\partial \hat{s}} -\varepsilon^2 r U_i \frac{\partial \hat{c}_i}{\partial\psi} = 0 \quad \mbox{at} \ \hat{s} = 0, \quad \psi \in [0,\psi_{l,crit})\ \mbox{and} \ \psi \in (\psi_{e,min},\psi_{e,max}), \end{equation}

at the bioreactor outlet. Here, $\psi _{i,min},\psi _{i,max}$ denote the minimum and maximum streamline values in domain $i=l,e$, respectively. At leading order in $\varepsilon$, (E6) and (E7) become

(E8)\begin{equation} \varepsilon^2\mathcal{P}_i\frac{\partial \hat{c}_i}{\partial \hat{s}} = \frac{\partial}{\partial \hat{s}}\left(\frac{1}{l_{\psi,i}\tilde{U}}_i\frac{\partial \hat{c}_i}{\partial \hat{s}}\right), \end{equation}

subject to

(E9)\begin{equation} \frac{\partial \hat{c}_i}{\partial \hat{s}} = 0 \quad \mbox{at} \ \hat{s} = 0, \quad \psi \in [0,\psi_{l,crit})\ \mbox{and} \ \psi \in (\psi_{e,min},\psi_{e,max}). \end{equation}

Equation (E8) is satisfied by the general solution

(E10)\begin{equation} \hat{c}_i = A_1(\psi) - A_2(\psi)\exp{\left(\varepsilon^2\mathcal{P}_i l_{\psi,i}\int_0^{\hat{s}}\tilde{U}_i\,\mbox{d}\hat{s}\right)}, \end{equation}

and solving subject to (E9) we see that $A_2 = 0$. $A_1$ can be solved for by matching into the concentration leaving the outer region, which shows that the concentration remains constant along streamlines in the boundary layer and equal to the concentration entering the boundary layer (i.e. there is no boundary layer at leading order for outlet flow).

We now investigate the lumen boundary layer at the lumen–membrane interface. This boundary layer is of $O (\varepsilon )$ thickness, given that in this region it is consistent to scale $r$ and $W_l$ as well as $s$. Therefore, we scale $s = 1 + \varepsilon \hat {s}$, $r = 1 + \varepsilon \hat {r}$ and, due to continuity of flux, scale $W_l = \varepsilon \hat {W}_l$ in the lumen. Using these scalings, (D1) becomes

(E11)\begin{equation} \varepsilon^2\mathcal{P}_l\frac{\partial \hat{c}_l}{\partial \hat{s}} = \varepsilon^2\frac{\partial}{\partial\psi}\left((1 + \varepsilon{\hat{r}})^2 l_{\psi,i}\hat{\tilde{U}}_l\frac{\partial \hat{c}_l}{\partial \psi}\right) + \frac{\partial}{\partial \hat{s}}\left(\frac{1}{l_{\psi,i}\hat{\tilde{U}}_l}\frac{\partial \hat{c}_l}{\partial \hat{s}}\right), \end{equation}

where

(E12)\begin{equation} \tilde{U}_l = \sqrt{\varepsilon^2 {U_l}^2 + W_l^2} = \sqrt{\varepsilon^2{U_l}^2 + \varepsilon^2\hat{W}_l^2} = \varepsilon\sqrt{{U_l}^2 + \hat{W_l}^2} = \varepsilon\hat{\tilde{U}}_l. \end{equation}

The lumen–membrane interfacial condition (D4) becomes

(E13)\begin{align} &\varepsilon(1 + \varepsilon\hat{r}) W_l\frac{\partial \hat{c}_l}{\partial\psi} + \frac{U_l}{l_{\psi,l}\hat{\tilde{U}}_l}\frac{\partial \hat{c}_l}{\partial \hat{s}} = \phi_m^2\left(\varepsilon^2r W_m\frac{\partial c_m}{\partial\psi} + \frac{U_m}{l_{\psi,m}\tilde{U}_m}\frac{\partial c_m}{\partial s}\right) \nonumber\\ &\quad \mbox{at} \ \hat{s} = 0, \psi \in [\psi_{l,crit},\psi_{l,max}). \end{align}

At $O (1)$, (E11) reduces to (E8) and thus gives the general solution (E10). Equation (E13) becomes

(E14)\begin{equation} \frac{1}{l_{\psi,l}}\frac{\partial \hat{c}_l}{\partial \hat{s}} = \frac{\phi_m^2}{l_{\psi,m}}\frac{\partial c_m}{\partial s} \quad \mbox{at} \ \hat{s} = 0, \psi \in [\psi_{l,crit},\psi_{l,max}), \end{equation}

noting that $\tilde {U}_m = U_m$ at leading order, and $\hat {\tilde {U}}_l = U_l$, $\hat {W}_l = 0$ given by the no-slip condition at the lumen–membrane interface (A4). We know the solution to $c_m$ in terms of the concentration at the lumen–membrane interface $c_l|_{[l,m]}$ given by (E5). Substituting (E10) and (E5) in to (E14), we find that $A_2 = \mathcal {A}$, and from matching into the outer solution in the lumen, we know that $A_1 = c_l|_{s \uparrow 1}$, giving

(E15)\begin{equation} \hat{c}_l = c_l|_{s \uparrow 1} - \mathcal{A}\exp\left(\varepsilon^2\mathcal{P}_l l_{\psi,l}\int_0^{\hat{s}} \tilde{U}_l\,\mbox{d}\hat{s}\right). \end{equation}

From the continuity of normal fluid flux condition (A4), we can write

(E16)\begin{equation} \left.\int_0^{\hat{s}} \tilde{U}_l \, \mbox{d} \hat{s} \right|_{\hat{s} = 0} = \left.\int_0^{\hat{s}} U_l \, \mbox{d} \hat{s} \right|_{\hat{s} = 0} = \left.\phi_m\int_1^{s} U_m\, \mbox{d} s \right|_{s = 1} = 0, \end{equation}

and thus using (E15) can express the lumen–membrane interfacial concentration as

(E17)\begin{equation} c_l|_{[l,m]} = c_l|_{s \uparrow 1} - \mathcal{A}. \end{equation}

Therefore substituting (E17) into (E5) and dropping superscript notation gives (3.11).

E.3 Lumen/ECS outer solution

We now solve for the nutrient concentration in the outer regions of the lumen and ECS. From (D2), $\tilde {U}_i = W_i$ at leading order in $\varepsilon$. Dropping the superscript $\star$, (D1) at leading order becomes (3.4). Owing to our asymptotic reduction, we see that (3.4) is first order with respect to $s$. Therefore, we require initial conditions in $s$ to solve for transport in the lumen and ECS, which come from inlet conditions (D5a), (D8). Additionally, the symmetry and no-flux condition at leading order becomes

(E18)\begin{equation} c_m = 1 - \mathcal{A}((s-1)(R_{m0}-1) + 1)^{-\xi}, \quad \xi = \frac{\varepsilon^2\mathcal{P}_l\kappa(P_e - P_l)}{\phi_m^2\ln R_{m0}}, \quad \mathcal{A} = \frac{b + \sqrt{b^2 + 4ac}}{2a}, \end{equation}

which is also given as (3.8) in the main body.

We solve (3.4) using the method of lines, whereby we discretise along the $\psi$ coordinate, resulting in a linear system of ODEs with respect to $s$ which can be solved readily using MATLAB's in-built ODE solver, ode15s.

E.4 Waste transport

We solve for the waste transport problem similarly to the nutrient transport problem in this appendix. As in the derivation of nutrient transport concentration in the membrane (3.11), we solve (3.6) swapping $c_i$ for $w_i$ subject to continuity of concentration

(E19)\begin{equation} w_l = w_m \quad \mbox{at} \ s = 1, \psi \in [\psi_{l,crit},\psi_{l,max}), \end{equation}

and the leading-order effective boundary condition

(E20)\begin{equation} \frac{\phi_m^2}{l_{\psi,m}}\frac{\partial w_m}{\partial s} = \varepsilon \phi_c^2\bar{\mathcal{R}}(c_m) \quad \mbox{at}\ s = 2, \psi \in (\psi_{e,min},\psi_{e,crit}], \end{equation}

giving

(E21)\begin{equation} w_m = w_l|_{s \uparrow 1} + \mathcal{B}((s-1)(R_{m0}-1) + 1)^{-\xi}, \end{equation}

where

(E22)\begin{equation} \mathcal{B} = {-}\frac{\varepsilon\phi_c^2\bar{\mathcal{R}}(c_m)\log{R_{m0}}}{\varepsilon^2 \mathcal{P}_{w,l}\kappa(P_e - P_l)R_{m0}^{-(\xi + 1)}}. \end{equation}