Transfer of mass and momentum at rough and porous surfaces

TL;DR

This paper develops a framework to characterize mass and momentum transfer at rough and porous surfaces, introducing measures like transpiration length and resistance coefficients, validated through canonical flow problems.

Contribution

It introduces a set of effective boundary conditions for non-smooth surfaces, linking slip, transpiration, and pressure jump to measurable surface properties.

Findings

Transpiration velocity arises from slip variations and mass conservation.

Normal momentum transfer causes a pressure jump characterized by resistance coefficients.

Effective boundary conditions are validated on canonical flow problems.

Abstract

The surface texture of materials plays a critical role in wettability, turbulence and transport phenomena. In order to design surfaces for these applications, it is desirable to characterise non-smooth and porous materials by their ability to exchange mass and momentum with flowing fluids. While the underlying physics of the tangential (slip) velocity at a fluid-solid interface is well understood, the importance and treatment of normal (transpiration) velocity and normal stress is unclear. We show that, when slip velocity varies at an interface above the texture, a non-zero transpiration velocity arises from mass conservation. The ability of a given surface texture to accommodate for a normal velocity of this kind is quantified by a transpiration length. We further demonstrate that normal momentum transfer gives rise to a pressure jump. For a porous material, the pressure jump can be characterised by so called resistance coefficients. By solving five Stokes problems, the introduced measures of slip, transpiration and resistance can be determined for any anisotropic non-smooth surface consisting of regularly repeating geometric patterns. The proposed conditions are a subset of effective boundary conditions derived from formal multi-scale expansion. We validate and demonstrate the physical significance of the effective conditions on two canonical problems -- a lid-driven cavity and a turbulent channel flow, both with non-smooth bottom surfaces.

Tables

Table 1. Table 1: The slip length L x ​ x subscript 𝐿 𝑥 𝑥 L_{xx} and the transpiration length M x ​ x subscript 𝑀 𝑥 𝑥 M_{xx} for a range of interface locations z i subscript 𝑧 𝑖 z_{i} above the textured surface. The effective tangential velocity u x subscript 𝑢 𝑥 u_{x} is sampled at ( 0.5 ​ H , z i ) 0.5 𝐻 subscript 𝑧 𝑖 (0.5\,H,z_{i}) . The effective transpiration velocity u z subscript 𝑢 𝑧 u_{z} is sampled at ( 0.25 ​ H , z i ) 0.25 𝐻 subscript 𝑧 𝑖 (0.25\,H,z_{i}) . The model predictions are finally normalised using the ensemble averaged results from DNS u ¯ x subscript ¯ 𝑢 𝑥 \bar{u}_{x} and u ¯ z subscript ¯ 𝑢 𝑧 \bar{u}_{z} .

$z_i/l$	$L_{xx}/l$	$M_{xx}/l$	$u_x/U_0$	$u_z/U_0$	$u_x/{\overline{u}}_x$	$u_z/{\overline{u}}_z$
$0.0$	$0.018$	$0.025$	$-1.24\cdot {10}^{-3}$	$1.09\cdot {10}^{-5}$	$0.989\pm 0.001$	$0.820\pm 0.087$
$0.1$	$0.118$	$0.061$	$-7.81\cdot {10}^{-3}$	$1.67\cdot {10}^{-4}$	$0.975\pm 0.000$	$0.933\pm 0.001$
$0.2$	$0.218$	$0.110$	$-1.39\cdot {10}^{-2}$	$5.22\cdot {10}^{-4}$	$0.961\pm 0.000$	$0.914\pm 0.001$
$0.3$	$0.318$	$0.160$	$-1.96\cdot {10}^{-2}$	$1.05\cdot {10}^{-3}$	$0.951\pm 0.000$	$0.895\pm 0.001$
$0.4$	$0.418$	$0.210$	$-2.50\cdot {10}^{-2}$	$1.74\cdot {10}^{-3}$	$0.943\pm 0.000$	$0.880\pm 0.001$
$0.5$	$0.518$	$0.259$	$-3.02\cdot {10}^{-2}$	$2.58\cdot {10}^{-3}$	$0.937\pm 0.000$	$0.869\pm 0.001$

Table 2. Table 2: Geometrical properties of the porous media considered in this work (see graphical representation in figure 9 b). The subscript i 𝑖 i corresponds to the interface and the subscript b 𝑏 b corresponds to the bulk. The porosity ϕ italic-ϕ \phi is defined as the ratio between the solid volume and the fluid volume. The last column shows the interior permeability tensor 𝐊 𝐊 \bm{\mathrm{{K}}} .

Configuration	Porosity ($\varphi$)	Geometry details	Permeability tensor ($𝐊$)
(i)	0.75	$r=0.28l$
(ii)	0.78	$\begin{array}{c}a=0.36l,b=0.19l,\\ \alpha ={45}^o\end{array}$
(iii)	$\begin{array}{c}{\varphi }_b=0.95\\ {\varphi }_i=0.80\end{array}$	$\begin{array}{c}{r}_b=0.13l\\ r_i=0.25l\end{array}$

Table 3. Table 3: The slip length, the transpiration length, and the resistance coefficients 𝐟 ( 1 ) superscript 𝐟 1 \bm{\mathrm{f}}^{(1)} and 𝐟 ( 2 ) superscript 𝐟 2 \bm{\mathrm{f}}^{(2)} for the porous medium geometries shown graphically in figure 9 (b). The last three columns show the ratio between the model and the DNS results for the slip and the transpiration velocities.

Config.	$L_{xx}/l$	$M_{xx}/l$	$f_x^{\left(1\right)}/\left(\mu /l\right)$	$f_z^{\left(1\right)}/\left(\mu /l\right)$	$f_x^{\left(2\right)}/\left(\mu /l\right)$	$u_x/{\overline{u}}_x$	$u_z^{-}/{\overline{u}}_z$	$u_z/{\overline{u}}_z$
(i)	$0.1516$	$0.0856$	$0.000$	$-10.43$	$\mathrm{\hspace{0.17em}0.000}$	$0.958$	$0.629$	$1.127$
(ii)	$0.1563$	$0.0885$	$2.125$	$-7.948$	$-1.541$	$0.958$	$0.661$	$1.137$
(iii)	$0.1538$	$0.0866$	$0.000$	$-38.23$	$\mathrm{\hspace{0.17em}0.000}$	$0.986$	$0.828$	$1.153$

Table 4. Table 4: Friction Reynolds numbers ( R ​ e τ 𝑅 subscript 𝑒 𝜏 Re_{\tau} ) at the bottom wall of the turbulent channel flow with smooth, rough and slip boundary conditions at the bottom wall. Here k 𝑘 k and δ 𝛿 \delta represent roughness elements’ height and channel half height, respectively. L x ​ x subscript 𝐿 𝑥 𝑥 L_{xx} and M x ​ x subscript 𝑀 𝑥 𝑥 M_{xx} are coefficients used in equations ( 5 , 12 ). For the considered roughness geometry L x ​ x = L y ​ y subscript 𝐿 𝑥 𝑥 subscript 𝐿 𝑦 𝑦 L_{xx}=L_{yy} , M x ​ x = M y ​ y subscript 𝑀 𝑥 𝑥 subscript 𝑀 𝑦 𝑦 M_{xx}=M_{yy} and L x ​ y = L y ​ x = M x ​ y = 0 subscript 𝐿 𝑥 𝑦 subscript 𝐿 𝑦 𝑥 subscript 𝑀 𝑥 𝑦 0 L_{xy}=L_{yx}=M_{xy}=0 .

	$k/\delta$	$L_{xx}/\delta$	$M_{xx}/\delta$	$R{e}_{\tau }$
No-slip	-	-	-	178.76
Slip only	0.04	0.01146	0	172.32
TR model	0.04	0.01146	0.01602	188.02
Geometry-resolved	0.04	-	-	184.70

Table 5. Table 5: Summary of the TR model boundary condition terms at the orders, at which corresponding terms emerge from the formal multi-scale expansion.

	Tangential velocity	Normal velocity	Pressure jump
	${𝐮}_t$	$u_n$	$p^{-}-p$
Leading order	${\mu }^{-1}𝐋\cdot 𝝉$	$0$	$-2\mu {\partial }_z{u}_z+{𝐟}^{\left(2\right)}\cdot {𝐮}_t$
Next order	$0$	$u_n^{-}-𝐌:{\nabla }_2{𝐮}_t$	${𝐟}^{\left(1\right)}\cdot {𝐮}^{-}$

Equations

$$u_x=L{\partial }_z{u}_x\text{mbox}onz=z_i.$$

$$p^{-}-p=-\mu C_{\pi }{\partial }_z{u}_x-2\mu {\partial }_z{u}_z.$$

$$u_z=u_z^{-}-M{\partial }_x{u}_x.$$

$$-p+2\mu {\partial }_z{u}_z=-p^{-}+f{u}_z^{-}.$$

$$\left(u_x,u_y\right)={\mathrm{u}}_t=\frac{\mathrm{L}}{\mu }\cdot \tau \text{mbox}onz=z_i,$$

$$-\nabla \hat{p}+\mu \Delta \hat{\mathrm{u}}$$

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

$$\hat{\mathrm{u}}={\hat{\mathrm{R}}}_{\tau }\cdot \tau .$$

$${\hat{\mathrm{u}}}^{(\tau x)}=\mu {\hat{\mathrm{R}}}_{\tau }\cdot {\mathrm{e}}_x\text{and}{\hat{\mathrm{u}}}^{(\tau y)}=\mu {\hat{\mathrm{R}}}_{\tau }\cdot {\mathrm{e}}_y.$$

$$-\nabla {\hat{p}}^{(\tau x)}+\mu \Delta {\hat{\mathrm{u}}}^{(\tau x)}=-\mu \delta \left(\hat{z}-{\hat{z}}_i\right){\mathrm{e}}_x,$$

$$-\nabla {\hat{p}}^{(\tau y)}+\mu \Delta {\hat{\mathrm{u}}}^{(\tau y)}=-\mu \delta \left(\hat{z}-{\hat{z}}_i\right){\mathrm{e}}_y,$$

$$\mathrm{u}=⟨\hat{\mathrm{u}}{⟩}_i=⟨{\hat{\mathrm{u}}}^{(\tau x)}{⟩}_i\frac{{\tau }_x}{\mu }+⟨{\hat{\mathrm{u}}}^{(\tau y)}{⟩}_i\frac{{\tau }_y}{\mu }\text{mbox}onz=z_i.$$

$$⟨\hat{a}{⟩}_i=⟨\hat{a}{⟩}_p\left({\hat{z}}_i\right),\text{mbox}where⟨\hat{a}{⟩}_p\left(\hat{z}\right)=\frac{1}{l^2}\underset{0}{\overset{l}{\int }}\underset{0}{\overset{l}{\int }}\hat{a}\left(\hat{x},\hat{y},\hat{z}\right)d\hat{x}d\hat{y}.$$

$$L_{xx}=⟨{\hat{u}}_x^{(\tau x)}{⟩}_i,L_{yx}=⟨{\hat{u}}_y^{(\tau x)}{⟩}_i,L_{xy}=⟨{\hat{u}}_x^{(\tau y)}{⟩}_i,L_{yy}=⟨{\hat{u}}_y^{(\tau y)}{⟩}_i.$$

$$u_n=u_z=-\mathrm{M}:{\nabla }_2{\mathrm{u}}_t\text{mbox}onz=z_i,$$

$$Q_1+Q_2+Q_3+Q_4+Q_5=0,$$

$$Q_i={\int }_{S_i}\mathrm{u}\cdot \mathrm{n}dS=\underset{{\hat{z}}_b}{\overset{{\hat{z}}_i}{\int }}\underset{s_1}{\overset{s_2}{\int }}\mathrm{u}\cdot \mathrm{n}d\hat{z}ds,$$

$$Q_i$$

$$=\underset{{\hat{z}}_b}{\overset{{\hat{z}}_i}{\int }}\underset{s_1}{\overset{s_2}{\int }}\left[⟨{\hat{\mathrm{R}}}_{\tau }{⟩}_p\left(\hat{z}\right)\cdot \tau \left(x,y\right)\right]\cdot \mathrm{n}d\hat{z}ds=\left({\mathrm{R}}_{\tau }\cdot \underset{s_1}{\overset{s_2}{\int }}\tau ds\right)\cdot \mathrm{n},$$

$${\mathrm{R}}_{\tau }=\underset{{\hat{z}}_b}{\overset{{\hat{z}}_i}{\int }}⟨{\hat{\mathrm{R}}}_{\tau }{⟩}_{p}d\hat{z}=\frac{1}{l^2}\underset{0}{\overset{l}{\int }}\underset{0}{\overset{l}{\int }}\underset{z_b}{\overset{z_i}{\int }}{\hat{\mathrm{R}}}_{\tau }\left(\hat{x},\hat{y},\hat{z}\right)d\hat{x}d\hat{y}d\hat{z}.$$

$$Q_5=\underset{x_1}{\overset{x_2}{\int }}\underset{y_1}{\overset{y_2}{\int }}{u}_{z}dxdy.$$

$$\underset{x_1}{\overset{x_2}{\int }}\underset{y_1}{\overset{y_2}{\int }}{u}_{z}dxdy=$$

$$-\left({\mathrm{R}}_{\tau }\cdot \underset{x_1}{\overset{x_2}{\int }}\left[\tau \left(x,y_2\right)-\tau \left(x,y_1\right)\right]dx\right)\cdot {\mathrm{e}}_y.$$

$$u_z\Delta x\Delta y=$$

$$-\left[{\mathrm{R}}_{\tau }\cdot \left(\tau \left(x,y+\Delta y\right)-\tau \left(x,y\right)\right)\Delta x\right]\cdot {\mathrm{e}}_y,$$

$$u_z=-\left({\mathrm{R}}_{\tau }\cdot {\partial }_x\tau \right)\cdot {\mathrm{e}}_x-\left({\mathrm{R}}_{\tau }\cdot {\partial }_y\tau \right)\cdot {\mathrm{e}}_y.$$

$$u_z=-{\mathrm{R}}_{\tau }:{\nabla }_2\tau .$$

$$\mathrm{M}=\mu {\mathrm{R}}_{\tau }\cdot {\mathrm{L}}^{-1}.$$

$$\left(\begin{array}{c}{u}_x\\ u_y\\ u_z\end{array}\right)=\left(\begin{array}{ccc}{L}_{xx}& 0& 0\\ 0& L_{xx}& 0\\ 0& 0& M_{xx}\end{array}\right)\cdot \left(\begin{array}{c}{\partial }_z{u}_x\\ {\partial }_z{u}_y\\ {\partial }_z{u}_z\end{array}\right),$$

$$u_n=u_n^{-}-\mathrm{M}:{\nabla }_2{\mathrm{u}}_t\text{mbox}onz=z_i.$$

$$u_n^{-}=\left(-\frac{\mathrm{K}}{\mu }\cdot \nabla p^{-}\right)\cdot \mathrm{n}.$$

$$-p+2\mu {\partial }_z{u}_z=-p^{-}+{\mathrm{f}}^{(1)}\cdot {\mathrm{u}}^{-}+{\mathrm{f}}^{(2)}\cdot {\mathrm{u}}_t\text{mbox}onz=z_i.$$

$$-\nabla \hat{p}+\mu \Delta \hat{\mathrm{u}}$$

Full text

1 Introduction

The physical behaviour of a number of fluid systems is dramatically modified by the presence of a small-scale surface roughness. For example, in wetting (figure 1a) – that is when a liquid in contact with a solid reaches the balance of surface tensions – the resulting apparent contact angle $\theta$ is very sensitive to the details of the surface texture (Wenzel, 1936; Quéré, 2008). At high Reynolds numbers (of order $1000$ and above), the pressure loss in turbulent pipes is a function of the wall roughness (figure 1b) (Nikuradse, 1950; Jiménez, 2004). Yet another example is the transport phenomena involving porous media, where the exchange of mass, momentum, energy, and other passive scalars between a free flowing fluid and a porous medium depends very much on the roughness at the interface between the two domains (figure 1c).

Engineers take an advantage of the sensitivity to the surface texture to modify large-scale flow features and to enhance the transport phenomena. Efficiency of heat exchangers (Mehendale et al., 2000; Agyenim et al., 2010) is highly dependent on the surface texture. In scaffold design for a bone regeneration, the cell growth on the implant (a porous biomaterial such as a calcium phosphate cement) depends on the interaction between the surrounding liquid and the surface texture of the implant (Dalby et al., 2007; Perez & Mestres, 2016). The performance of fuel cells depends on the ability of gas flow to efficiently transport the water vapour away from the cathode, a thin porous medium (Prat, 2002; Haghighi & Kirchner, 2017). Turbulent skin friction on wings or turbine blades can be reduced by using riblets, which are able to push quasi-streamwise vortices away from the wall (Walsh & Lindemann, 1984).

The design of the surface texture in the examples mentioned above is based on a trial and error procedure, that may require tremendous amount of effort, time and expensive surface manufacturing equipment. The formulation presented in this paper provides a framework for modelling the interaction between free flows and various textured and porous surfaces. Our modelling approach provides a direct relationship between the microscopic geometrical details of a complex surface and the associated macroscopic transport of mass and momentum. Thus, it has the potential to replace the trial and error procedure in the design phase.

Due to the multiscale nature of the problems described above, fully resolved numerical investigations – of both the complex surface and the free flow above it – are practically impossible to perform in applied settings. Effective approaches are actively pursued to circumvent this difficulty. In this way, one can capture the averaged effect of the microscale features on the macroscopic processes, and hence avoid resolving microscopic geometric details. Some recent examples of effective modelling applied to drying, cell growth and heat exchange can be found in works by Mosthaf et al. (2014); Vaca-González et al. (2018); Laloui et al. (2006); Wang et al. (2018). The main challenge for effective models describing fluid-surface interaction is the specification of a boundary condition at an artificially created interface between the free-fluid region and the complex surface. Despite the recent advancements, we still lack interface conditions that capture the dominant physical features associated with complex anisotropic surfaces.

Before highlighting the main ingredients of our model, we make a brief account of the current state-of-the-art of effective boundary conditions of textured and porous surfaces. A two-dimensional configuration is sufficient for this purpose. The streamwise and wall-normal coordinates are denoted by $x$ and $z$, where the effective boundary conditions are imposed at a planar interface at coordinate $z=z_{i}$. For rigid textured surfaces with a characteristic size $l$, one may impose the slip velocity condition (Navier, 1823) as an effective boundary condition,

[TABLE]

Here, $u_{x}$ is the tangential velocity component at the interface and $L\sim l$ is the slip length. Geometrically – as shown in figure 2(a) – the slip length is the distance that the velocity profile has to be linearly extrapolated to reach zero value. There has been extensive development of the slip boundary condition for textured and porous surfaces (Saffman, 1971; Sahraoui & Kaviany, 1992; Miksis & Davis, 1994; Sarkar & Prosperetti, 1996; Gupte & Advani, 1997; Jäger & Mikelić, 2001; Stroock et al., 2002; Bolanos & Vernescu, 2017). In current approaches, the interface normal (transpiration) velocity is typically set either to zero $u_{z}=0$ (for textured surfaces) or to the interior flow $u_{z}=u_{z}^{-}$ (for porous surfaces) due to mass conservation arguments or as the leading-order boundary condition (Mohammadi & Floryan, 2013; Lācis & Bagheri, 2016; Jiménez Bolaños & Vernescu, 2017).

Configurations with porous surfaces require additional boundary conditions. If the bulk of the surface is governed by the Darcy-Brinkmann equation, which is typical for works considering method of volume averaging (Whitaker, 1998), stress jump conditions are often derived (Ochoa-Tapia & Whitaker, 1995; Valdés-Parada et al., 2009, 2013; Angot et al., 2017). In the current work, however, we consider only Darcy’s law within the bulk of the porous surface. Consequently, a condition for the Darcy pressure or the pore pressure $p^{-}$ is needed. The pressure continuity $p=p^{-}$, where $p$ is the free fluid pressure, has been a common choice in the past (Ene & Sanchez-Palencia, 1975; Levy & Sanchez-Palencia, 1975; Hou et al., 1989; Lācis & Bagheri, 2016). The most notable recent theoretical and numerical developments (Marciniak-Czochra & Mikelić, 2012; Carraro et al., 2013, 2018) have resulted in the pressure jump condition

[TABLE]

Here, $\mu$ is the fluid dynamic viscosity and $C_{\pi}$ is a stabilisation parameter derived from matching boundary layer solutions with exterior solutions. The coefficient $C_{\pi}$ is non-zero only for anisotropic porous surfaces. The pressure interface condition – as well as the velocity interface condition – for porous media have been a subject of many investigations and is still debated (Beavers & Joseph, 1967; Han et al., 2005; Le Bars & Grae Worster, 2006; Rosti et al., 2015; Zampogna & Bottaro, 2016; Mikelić & Jäger, 2000; Jäger & Mikelić, 2009; Carraro et al., 2015, 2018; Lācis & Bagheri, 2016; Zampogna et al., 2019).

In this work, we extend the above conditions with new terms for the wall-normal velocity condition and the pressure condition. Our proposed set of boundary conditions is called transpiration-resistance (TR) model, and it is applicable for any textured or porous surface consisting of regular repeating geometric entities. The TR model captures the transport of interface tangential momentum as well as the transport of mass and interface normal momentum. It is a homogenised boundary condition, valid for configurations with a scale separation $\epsilon=l/H\ll 1$, where $H$ is the characteristic length scale of the free fluid. It consists of the slip boundary condition (1) for the interface tangential velocity. The wall normal velocity in the TR model is

[TABLE]

The first term is the seepage Darcy velocity, given by $u^{-}_{z}=\left(K/\mu\right)\,\partial_{z}p^{-}$, where $K$ is the interior permeability. The second term quantifies how much a surface texture allows exchange of mass with the surrounding fluid due to a streamwise variation of the slip velocity. Using continuity, the above condition for textured surfaces (where $u^{-}_{z}=0$) can be written as $u_{z}=M\,\partial_{z}u_{z}$. Geometrically (figure 2b), the transpiration length $M$ is thus the distance below the interface for which a non-zero transpiration velocity $u_{z}$ can exist. This depth is obtained under an assumption of a linear decay of the velocity with a slope $\partial_{z}u_{z}$. For a porous surface, the TR model provides the pressure condition,

[TABLE]

Here, the left hand side is the normal stress of the outside free flow on the interface plane, and the right hand side is the normal stress from the porous material. The resistance coefficient $f$ quantifies the friction force that the Darcy seepage velocity generates while passing through the interface (figure 2c).

Three assumptions underlies the the proposed TR model.

A1

Creeping flow assumption Re $<1$, which allows to solve a given flow problem near the interface with help of linear decomposition.

A2

Scale separation assumption $\epsilon\ll 1$, which leads to constant macroscopic flow field variables over the characteristic length $l$ of the surface.

A3

The surface is homogeneous i.e., it consists of repeating geometric entities or elements, which allows to consider only single structure to determine surface properties.

Under these assumptions, the slip length, the transpiration length and the resistance coefficients are properties of the surface texture only, and can be computed by solving five fundamental Stokes problems. For a given texture, the knowledge of these effective coefficients provides important information of the diffusive/advective transport into the material as well as the ability of the solid skeleton to resist externally imposed shear stress. The TR model is based on conditions derived from a formal multi-scale expansion in the small parameter $\epsilon$. By including higher-order terms for the transpiration velocity and for the pressure, we will show using numerical simulations that the error of the TR model is close to $\mathcal{O}\left(\epsilon^{2}\right)$.

This paper is organized as follows. In sections 2 and 3, we describe and validate the TR model for textured surfaces and porous surfaces, respectively. In section 4, we show using the turbulent channel flow that the transpiration velocity in the TR model – despite being a higher-order term from an asymptotic viewpoint – is essential from a physical viewpoint. In section 5, the TR model is discussed in the context of formal multi-scale expansion and, finally, we provide conclusions in section 6.

2 A model for textured surfaces

In this section, we present the transpiration-resistance (TR) model for 3D textured surfaces in contact with a free flowing fluid, under the assumptions (A1–A3). First, we explain the boundary condition for the interface tangential velocity (the slip condition) and show how to obtain the associated slip length tensor. Then, we introduce the transpiration velocity condition and demonstrate how to determine the transpiration length tensor by making use of mass conservation. Finally, we compute the slip and transpirations tensors and validate the model by using fully resolved numerical simulations. The effect of the interface location on the accuracy of the TR model is discussed in the last subsection.

2.1 Tangential interface velocity and slip length

The tangential velocity condition in the TR model is provided by the standard slip condition, which for 3D textured surfaces reads

[TABLE]

where $\bm{\mathrm{\tau}}=\mu\left(\partial_{z}u_{x},\partial_{z}u_{y}\right)$ and $\bm{\mathrm{{L}}}=\left(L_{xx},L_{xy};L_{xy},L_{yy}\right)$ is the symmetric positive definite (Kamrin & Stone, 2011) surface slip length tensor. Here, $\bm{\mathrm{u}}_{t}$ is the tangential velocity vector. The tangential $t$ subscript is used interchangeably with the $x$ and $y$ components.

Let us consider a patterned wall and a vortical flow over it, as illustrated in figure 3, left. The scale separation assumption (A2) allows us to introduce two different spatial coordinates: $x_{i}$ and $\hat{x}_{i}$. The former is used for describing spatial variations over large length scales ($x_{i}\sim H$). The latter is used to describe microscopic variations over much smaller roughness scale ($\hat{x}_{i}\sim l$). The effective boundary condition (5) is a macroscopic condition; the microscopic features of the texture are embedded in an averaged sense in the slip length tensor $\bm{\mathrm{{L}}}$.

To determine $\bm{\mathrm{{L}}}$, we consider a small volume near the surface of the texture with a cross section $l\times l$. This volume contains one representative surface structure, see figure 3. Within this interface cell, the scale separation assumption (A2) allows us to treat the shear stress from the free fluid $\bm{\mathrm{\tau}}$ as spatially constant external parameter. Due to the creeping flow assumption (A1), the equations governing the flow response to the free fluid shear stress are the Stokes equations,

[TABLE]

This set of equations is equivalent to a two-domain description employing velocity continuity and stress jump at the interface (appendix B). Additionally, equations (6–7) are the same as previously used and derived by Luchini et al. (1991); Kamrin et al. (2010); Luchini (2013). The imposed boundary conditions are no-slip and no-penetration at the surface of the solid structure ($\hat{z}=\hat{z}_{b}$). We impose periodic conditions at the vertical faces of the interface cell (due to the assumption A3). At the top surface of the cell, we impose zero-stress condition to keep the shear stress at the interface as the only driving force of the problem.

The linearity assumption (A1) allows us to write the solution as a product between a response operator $\bm{\mathrm{{\hat{R}}}}_{\tau}$ and the free fluid shear stress,

[TABLE]

This can be expanded as {subeqnarray} ^u &= ( ^R_τ⋅e_x )  τx + ( ^Rτ⋅e_y )  τ_y

= ^u^(τx)   τxμ + ^u^(τy)   τyμ,

where we have defined

[TABLE]

Thus, the velocity fields $\hat{\bm{\mathrm{u}}}^{(\tau x)}$ and $\hat{\bm{\mathrm{u}}}^{(\tau y)}$ are solutions of the following two fundamental problems,

[TABLE]

The fundamental problems are forced in $x$ and $y$ directions, respectively, with unit shear at the plane $z_{i}$.

Taking the surface average of expression (8$b$) at the interface, we obtain

[TABLE]

No average is carried out for the free fluid shear stress, because it is constant within the interface cell (A2). The surface average of an arbitrary quantity $\hat{a}$ is defined as

[TABLE]

By comparing the surface averaged velocity in the interface cell (9) with the slip boundary condition (5), we observe that the components of the slip length tensor can be obtained as

[TABLE]

In terms of the response operator, the slip tensor becomes $\bm{\mathrm{{L}}}=\mu\,\langle\bm{\mathrm{{\hat{R}}}}_{\tau}\rangle_{i}$. Note that dimension of the vector fields $\bm{\mathrm{\hat{u}}}^{(\tau x)}$ and $\bm{\mathrm{\hat{u}}}^{(\tau y)}$ is velocity per shear, which gives unit of length.

2.2 Interface normal velocity and transpiration length

We begin with a simple motivation for the transpiration velocity based on the mass conservation. We consider a two dimensional rough surface and define a control volume (CV) below the interface $z=z_{i}$ as shown in figure 4. We assume that there is a slip velocity variation from $u_{x,1}$ at the left side of the CV to $u_{x,2}$ at the right side of the CV. The mass fluxes at the left and the right boundaries of the CV are proportional to the slip velocities at the interface (direct consequence of A1). Consequently, mass conservation requires a non-zero transpiration velocity $u_{z}$ at the interface. If the slip velocity is increasing $u_{x,2}>u_{x,1}$, the generated transpiration velocity is, therefore, negative.

More generally, the interface normal velocity condition in the TR model for a 3D textured surface is provided by a linear law relating the normal velocity with the tangential variation of slip velocity,

[TABLE]

where $\bm{\mathrm{{M}}}=\left(M_{xx},M_{xy};M_{xy},M_{yy}\right)$ is the transpiration length tensor – exhibiting the same symmetry properties as the slip length tensor – and $\nabla_{2}=\left(\partial_{x},\partial_{y}\right)$ is gradient operator containing the two tangential directions. We use the normal $n$ subscript interchangeably with the $z$ component. The proposed expression, motivated from the mass conservation, also emerges from a formal multi-scale expansion (see section 5 and appendix A). Bottaro (2019) has recently used the multi-scale expansion to confirm the transpiration velocity condition proposed here.

To determine $\bm{\mathrm{{M}}}$, we make use of mass conservation in a 3D setting. We define a CV with size $\left(x_{2}-x_{1}\right)\times\left(y_{2}-y_{1}\right)\times\left(\hat{z}_{i}-\hat{z}_{b}\right)$ over a number of texture elements as shown in figure 5. Taking into account that there can be no flux through the impermeable bottom surface, mass conservation requires

[TABLE]

where $Q_{i}$ are the volumetric flux through faces of the CV (figure 5). The flux through the vertical faces ($i=1,\dots,4$) of the CV can be evaluated as

[TABLE]

where $\bm{\mathrm{u}}$ is the effective velocity field at the CV face, $\bm{\mathrm{n}}$ is the unit normal vector of the surface and $s$ is either $x$ or $y$, depending on which surface the integral is carried over. Note that the integral in the wall-normal direction is carried out over the microscale $\hat{z}$, because macroscopically the textured surface is infinitesimal and variations in depth does not exist.

Next, we use equation (10) in conjunction with the solution of fundamental problems (FP1,FP2) to rewrite (14),

[TABLE]

where we have defined the response tensor $\bm{\mathrm{{R}}}_{\tau}$ as

[TABLE]

The flux through the top wall can be expressed as

[TABLE]

Inserting the expressions for the fluxes through the CV faces into the mass conservation identity (13) we obtain

[TABLE]

To proceed towards the effective boundary condition (12), we take an infinitesimal CV limit, which gives us

[TABLE]

were we have $\Delta y=y_{2}-y_{1}$, $\Delta x=x_{2}-x_{1}$ and $x_{1}=x$ and $y_{1}=y$. Dividing both sides by $\mathit{\Delta x}\mathit{\Delta y}$ and using the definition of a derivative, we obtain

[TABLE]

This expression can be rewritten using double contraction as

[TABLE]

To obtain the transpiration length tensor, we express the tangential shear stress from equation (5) and insert the result into (21). Comparing the final result with equation (12) yields

[TABLE]

Recall that the tensor $\bm{\mathrm{{R}}}_{\tau}$ can be obtained as a post-processing step from the fundamental problems (FP1, FP2) using the volume integral (16).

It is interesting to note that the velocity conditions (5,12) can be written in a more compact form,

[TABLE]

valid for an incompressible flow over isotropic geometries or for incompressible two-dimensional flows. The upper left $2\times 2$ block corresponds to the slip length tensor $\bm{\mathrm{{L}}}$ introduced before, while the lower right element $M_{xx}$ is the first term of the transpiration length tensor $\bm{\mathrm{{M}}}$. The equivalence with the previous formulation can be seen through the application of continuity, i.e., $\partial_{z}u_{z}=-\partial_{x}u_{x}-\partial_{y}u_{y}$. The form (23) can be useful in practice, for example, if boundary conditions are imposed weakly in a finite element method. Such a set of boundary conditions was numerically investigated by Gómez de Segura et al. (2018). In their work, the focus was on elucidating the turbulent flow response to the boundary condition (23) where all the coefficients for the slip and the transpiration lengths could take different values.

2.3 Numerical validation of velocity conditions

We consider a lid-driven cavity whose bottom surface is made of a texture with the characteristic length scale $l$ (figure 6a). The macroscopic length scale $H$ corresponds to the cavity length and the cavity height. The scale separation parameter is set to $\epsilon=l/H=0.1$. A no-slip condition is applied on all surfaces except the top wall, which moves with a prescribed velocity $(U_{0},0)$. Details about the numerical solver can be found in appendix C.1.

The moving upper wall generates a clock-wise rotating vortex. This vortex imposes a negative shear on the rough surface. It also induces a downward mass flux at the right half of the cavity and an upward mass flux at the left half of the cavity. Near the surface texture, one can observe velocity fluctuations with a wavelength corresponding to the texture size $l$. To obtain macroscopic flow fields from DNS, we average out the microscale oscillations by creating an ensemble of $50$ DNS simulations. The ensemble consists of configurations in which the textured surface at the bottom of the cavity is incrementally shifted in the $x$ direction. The tangential and the transpiration velocities from the ensemble averaged DNS along $z=0.2\,l$ are shown with black lines in figure 7(a,b), respectively.

For comparison, we set up an effective simulation of the problem with the domain and boundary conditions shown in figure 6(b). We position the interface at the previously selected coordinate $z_{i}=0.2\,l$. Coefficients for the boundary conditions (5,12) – the slip and the transpiration lengths – are obtained as described in sections 2.1 and 2.2, using a FreeFEM++ open-source code (Lācis & Bagheri, 2016–2019). For the chosen configuration and interface location, we have $L_{xx}=0.218\,l$ and $M_{xx}=0.110\,l$. The velocities at $z=0.2l$ from the effective simulation are compared to the ensemble averaged DNS in figure 7(a,b). It is clear that the employed boundary conditions accurately predict the ensemble average (or the macroscopic variation) of both velocity components.

The results obtained using the TR model are not sensitive to the exact interface location. To show this, we repeat the previous effective computation for a range of interface locations $z_{i}=\left(0.0,0.1,0.3,0.4,0.5\right)\,l$. For each interface coordinate, the effective coefficients $L_{xx}$ and $M_{xx}$ are recomputed using the fundamental problems (FP1–FP2), see table 1. To quantitatively present the TR model predictions, we select two streamwise positions at the interface plane, shown with vertical dashed lines in figure 7. In table 1, we show the model predictions of $u_{x}$ sampled at point $\left(0.5\,H,z_{i}\right)$ for all interface locations. As the interface moves upwards – further away from the solid structures –, the value of the predicted slip velocity increases due to a larger distance over which the viscous friction can bring the velocity to the no-slip value at the wall. This effect is correctly captured by the model through the linear increase of the slip length $L_{xx}$ (table 1). In other words, the information about the interface location is provided to the effective model through the adjustment of the coefficients. Similar behaviour can be observed also for the interface normal velocity component $u_{z}$ and transpiration length $M_{xx}$.

In the last two columns of table 1 we present the ratio between model predictions and the ensemble averaged DNS result. The error estimates indicate uncertainty in the averaged DNS due to the presence to some unfiltered small-scale fluctuations. As one can see, this uncertainty is significant only for transpiration velocity at coordinate $z=0$. From table 1, we observe that for all interface locations the relative error ($1-u/\bar{u}$) is below $7\%$ for the slip velocity. The relative error for the transpiration velocity is below $14\%$. There is a trend of an increasing error as the interface is moved upward, with an exception for interface location $z_{i}=0$. This exception arises due to the large uncertainty in the reference result. Despite the trend of increasing error with interface location, the TR model has a remarkably good accuracy taking into account that the transpiration velocity is varied over two orders of magnitude (see 5th column of table 1).

We have carried out similar numerical computations on equilateral triangular surface texture and obtained the same behaviour as reported above. This investigation shows that it is possible to adjust the interface height over distances $\mathcal{O}\left(l\right)$ without a significant loss of accuracy. Such invariance of interface location has already been demonstrated numerically by Lācis & Bagheri (2016) and theoretically by Marciniak-Czochra & Mikelić (2012) for the slip velocity alone. However, as a “rule-of-thumb”, we suggest to place the interface as close to the solid structure as possible without intersecting the solids.

3 The TR model for porous surfaces

In this section, we extend the TR model to 3D porous surfaces by augmenting the set of boundary conditions from the previous section with a pressure condition. This is achieved by considering the transfer of normal momentum between the free flow region and the porous surface.

To determine the coefficients appearing in the interface boundary conditions, we will adopt a similar interface-cell approach as for the textured surface (see figure 8, leftmost frame). The bottom coordinate of the interface cell, $\hat{z}_{b}$, is chosen such that all flow variations have decayed inside the porous medium, where only the interior (Darcy) flow remains. As a rule-of-thumb, the interface cell should contain around four solid skeleton entities $\hat{z}_{b}\approx\hat{z}_{i}-4\,l$. From scale separation assumption (A2) it follows that, in the interface cell, shear stress from the free fluid and the pore pressure gradient from the porous material are both constant.

3.1 Velocity boundary conditions

For a porous surface, the tangential velocity boundary condition is identical to the textured surface, i.e. the slip condition (5). The interface normal velocity condition, on the other hand, is

[TABLE]

Here, $u_{n}^{-}$ is the interface normal Darcy velocity, satisfying

[TABLE]

This term is induced by mass conservation between free fluid and pore flow.

The tensors $\bm{\mathrm{{M}}}$ and $\bm{\mathrm{{L}}}$ are determined by solving the fundamental problems (FP1–FP2). The only difference from the textured surface is that the solid structures within the interface cell represent the porous material. Consequently, all the elements in $\bm{\mathrm{{L}}}$ and $\bm{\mathrm{{M}}}$ can be obtained through expressions (11) and (22), respectively. The interior permeability tensor ($\bm{\mathrm{{K}}}$) of the porous medium is computed through a set of Stokes equations in a bulk unit cell (Whitaker, 1998; Mei & Vernescu, 2010).

3.2 Pressure boundary condition

The pressure boundary condition for a general 3D porous surface, obeying assumptions (A1–A3), is obtained through a balance between the normal free-fluid stress and the stress from the porous material, i.e.

[TABLE]

The normal stress from the porous material consists of the pore pressure $p^{-}$ and two friction coefficients, $\bm{\mathrm{f}}^{(1)}$ and $\bm{\mathrm{f}}^{(2)}$. The coefficient $\bm{\mathrm{f}}^{(1)}$ describes the interface normal resistance that the Darcy flow $\bm{\mathrm{u}}^{-}$ must overcome to transport mass and momentum across and along the interface. The coefficient $\bm{\mathrm{f}}^{(2)}$ provides the interface normal force due to the slip velocity near the interface. It exists only for anisotropic surface geometries, similarly to the stabilisation parameter (see equation 2) derived by Marciniak-Czochra & Mikelić (2012); Carraro et al. (2018).

The friction coefficients are again determined by considering the interface cell (figure 8). The difference from the textured surface is the existence of the pore pressure resulting in an additional forcing term in the governing equations of the interface cell, yielding:

[TABLE]

Darcy’s law is valid only in the porous material, therefore the pressure gradient forcing is considered only below the interface. A 1D Heaviside step function $H\left(\hat{z}\right)$ is used to distinguish between regions above and below the interface. Boundary conditions for the interface cell are the same as for the equations (6–7) except at the bottom of the domain, where we have to impose the interior solution corresponding to the Darcy flow due to the same pressure gradient $\nabla p^{-}$ (Whitaker, 1998; Mei & Vernescu, 2010; Lācis & Bagheri, 2016).

To continue, we make use of the linearity assumption (A1) and write the pressure as

[TABLE]

Here, $\bm{\mathrm{\hat{r}}}_{\tau}$ and $\bm{\mathrm{\hat{r}}}_{p}$ are the response operators related to the shear stress $\bm{\mathrm{\tau}}$ and the pressure gradient $\nabla p^{-}$, respectively. This expression is expanded as

[TABLE]

where we have defined

[TABLE]

and

[TABLE]

Note that $\hat{p}^{(\tau x)}$ and $\hat{p}^{(\tau y)}$ are the pressure fields appearing in the fundamental problems (FP1,FP2); they are the pressure responses to the interface shear forcing in the $x$ and $y$ directions, respectively. Furthermore, $\hat{p}^{(px)}$, $\hat{p}^{(py)}$ and $\hat{p}^{(pz)}$ are the pressure fields associated with the following three fundamental problems

[TABLE]

These problems describe the response to the pressure gradient forcing along the three coordinates and have been previously derived by Lācis & Bagheri (2016) using formal multi-scale expansion. Keep in mind that in equations (FP3–FP5) the fields $\hat{\bm{\mathrm{u}}}^{(px)}$, $\hat{\bm{\mathrm{u}}}^{(py)}$ and $\hat{\bm{\mathrm{u}}}^{(pz)}$ have the dimension of length squared, similar as the permeability of a porous medium.

From the five fundamental problems (FP1–FP5), we can determine the resistance vectors $\bm{\mathrm{f}}^{(1)}$ and $\bm{\mathrm{f}}^{(2)}$ for the pressure condition (26). We begin with $\bm{\mathrm{f}}^{(2)}$. It generates a pressure jump

[TABLE]

where $\tilde{\bm{\mathrm{f}}}^{(2)}=\bm{\mathrm{f}}^{(2)}\cdot\bm{\mathrm{{L}}}$. The pressure field response (30) due to the shear is

[TABLE]

Next, we need to relate the effective pressures in porous and free fluid regions with the linear pressure responses $\hat{p}^{(\tau x)}$ and $\hat{p}^{(\tau y)}$ in the interface cell. For the velocity, a simple plane average at the interface (9) was sufficient. However, for the pressure condition a single pressure value will not provide the necessary information about the pressure jump. Therefore, we define the effective pressure in the interior and free fluid as

[TABLE]

Here, $V_{f}$ corresponds to fluid volume in the integration region. To neglect any transition effects of the pressure field near the interface, these volume averages are taken at the bottom and at the top of the interface cell. In this way, the averaging operation is sufficiently far away from the interface to obtain a representative pressure value for the interior and the free fluid.

Now we insert the pressure field decomposition (34) into equation (35) and we take the difference between the interior pressure and the free fluid pressure,

[TABLE]

By comparing the above to equation (33), we obtain

[TABLE]

We emphasise that $\hat{p}^{(\tau x)}$ and $\hat{p}^{(\tau y)}$ are pressure fields in the interface cell generated due to shear stress forcing (figure 8, middle frame) and can be computed from the fundamental problems (FP1,FP2). Finally, the resistance vector $\bm{\mathrm{f}}^{(2)}$, appearing in the front of the slip velocity in equation (26), is obtained from

[TABLE]

The procedure to get this friction coefficient is similar to the one reported by Marciniak-Czochra & Mikelić (2012); Carraro et al. (2013).

We turn our attention to the resistance coefficient $\bm{\mathrm{f}}^{(1)}$. The pressure jump condition (26) due to the Darcy velocity is

[TABLE]

where $\tilde{\bm{\mathrm{f}}}^{(1)}=-\bm{\mathrm{f}}^{(1)}\cdot\bm{\mathrm{{K}}}$. The pressure field response (30), corresponding to the pore pressure gradient forcing, is

[TABLE]

Using (35), we can express the pressure jump as

[TABLE]

By comparing equations (38) and (40), we identify the friction vector components as

[TABLE]

We recall that the pressure fields in the interface cell are generated by the pore pressure gradient forcing below the interface (figure 8), and they are computed from fundamental problems (FP3–FP5). The final form of the friction coefficient is,

[TABLE]

This friction coefficient term, which to the best of authors’ knowledge is reported for the first time, is particularly important for capturing the correct pressure jump across the interface for layered problems, as we demonstrate in the next section.

3.3 Validation of the TR model for porous surfaces

We consider the same flow configuration as in section 2.3, but we replace the bottom textured wall with a porous surface, see figure 9(a). The porous medium consists of a periodic distribution of solid inclusions with a characteristic length scale $l$ (see $l\times l$ square in figure 9a). The width and the height of the cavity is $H$, while the depth of the porous material is $H/2$. The scale separation parameter is again set to $\epsilon=l/H=0.1$. The flow reaches the interior seepage velocity quickly (Lācis & Bagheri, 2016; Lācis et al., 2017); therefore, a porous material containing only five repeating structures in depth is sufficient for using Darcy equation in the interior.

To demonstrate the generality of the TR model, we consider three kinds of porous geometries, shown in figure 9(b). Configuration (i) has circular solid inclusions, which results in an isotropic porous medium. The anisotropic elliptic inclusions considered in configuration (ii) are the same as investigated by Carraro et al. (2013). The last geometry (iii) has isotropic circular inclusions with the interface layer different from the interior. The porosity and geometrical details of the three configurations are listed in table 2.

We carry out fully resolved DNS. For each configuration, the mean over an ensemble of $50$ shifted porous beds is computed. The free flow is similar to the flow in the lid driven cavity with the textured bottom (section 2.3) except that there is a mass flux in and out of the porous material.

The averaged DNS will be compared to effective representations of the porous bed (figure 9c). Within the porous domain, we employ the Darcy’s law, where the only unknown quantity is the pore pressure $p^{-}$. The interior permeability tensors ($\bm{\mathrm{{K}}}$) for all geometries are listed in the last column of table 2. They were obtained by solving a set of Stokes equations in a periodic unit cell in the bulk (Whitaker, 1998; Mei & Vernescu, 2010). A Neumann condition on pore pressure $\nabla p^{-}\cdot\bm{\mathrm{n}}=0$ is enforced at solid boundaries of the cavity. This condition corresponds to zero fluid flux through the wall. Boundary conditions for the free fluid remain the same as in section 2.3. We place the interface at a distance $z_{i}=0.1\,l$ above the solid structures; the interface location $z_{i}=0.0$ was not accessible due to meshing issues. We compute the effective parameters appearing in boundary conditions (5,24,26) using the procedure explained in sections 2.1,2.2 and 3.2. The slip and the transpiration lengths, reported in table 3, have nearly the same value for the three configurations because the porous materials have similar porosity near the interface (table 2).

To validate the velocity conditions (5,24), we sample the ensemble-averaged DNS and the effective model at coordinates $(0.5\,H,z_{i})$ and $(0.25\,H,z_{i})$ for slip and transpiration velocities, respectively. The ratio between the model predictions and the DNS are given in last columns of table 3. It is clear that the slip velocity is predicted as accurately as for the textured surfaces. In the second to last column of table 3 we list the ratio between the Darcy transpiration velocity $u_{z}^{-}$ – sampled just below point $(0.25\,H,z_{i})$ – and the DNS result. The Darcy velocity alone is a rather inaccurate predictor and has a relative error up to $37\%$. The agreement between the TR model (24) – that augments the Darcy contribution with the term containing the transpiration length – and the DNS is better as the relative error is smaller than $15\%$.

To validate the pressure condition (26), we analyze the ensemble-averaged DNS of configuration (iii). The pressure $\hat{p}$ along the vertical dashed lined in figure 9(a) is shown in figure 10(a) with a solid black line. The region corresponding to porous domain is shadowed. We observe that the pressure field undergoes a sharp variation when transitioning from the free fluid region to the porous medium, as shown by an inset in figure 10(a). The sharp variation indicates that there could be a pressure jump in the effective representation. Note that there are microscale oscillations in the pressure field; the employed ensemble average filters out only the microscale variations in the $x$ direction. For comparison, the pressure obtained from the effective model is shown using dotted blue symbols in figure 10(a). We can observe that the agreement between model predictions and DNS is good both in free fluid and interior, while the sharp variation in the near vicinity of the interface is not modelled. This is a direct consequence of having an infinitely thin interface in the effective representation, which condenses all variations near the boundary to a single line. The same quantitative agreement is observed for configurations (i) and (ii), but are not shown here.

To show the importance of the resistance coefficients $\bm{\mathrm{f}}^{(1)}$ and $\bm{\mathrm{f}}^{(2)}$, we carried out two more effective simulations. In the first one – called “leading order” – we set $\bm{\mathrm{f}}^{(1)}=0$ in the pressure condition (26). This essentially corresponds to pressure condition proposed by Marciniak-Czochra & Mikelić (2012); Carraro et al. (2013, 2018). In the second one – called “pressure continuity” – we impose $p=p^{-}$, which has been a common approach in the past (Ene & Sanchez-Palencia, 1975; Levy & Sanchez-Palencia, 1975; Hou et al., 1989; Lācis & Bagheri, 2016). Results from the leading order model and the pressure continuity model are reported in figure 10(a) using crosses and a dashed curve, respectively. We observe that both conditions result in a poor agreement between the model and the DNS if compared to the TR model.

An inaccurate pressure condition can also influence the flow field. To illustrate this, we provide model results of interface normal velocity in figure 10(b), using the same symbols. We observe that the error in pressure condition can lead to significantly different – and inaccurate – vertical velocity predictions. The coefficient $\bm{\mathrm{f}}^{(1)}$ imposes larger resistance for the wall-normal velocity, and, thus, it decreases the transpiration by precisely the correct amount.

4 Role of the transpiration in a turbulent channel flow

In this section, we demonstrate with a specific example – fully-developed turbulent flow over a textured surface – that a small transpiration velocity can be crucial to capture the correct physics of the problem. Here, the domain consists of a periodic channel whose bottom surface is covered with ordered cuboid roughness elements (see inset in figure 11a).

We define a flat interface on the crest plane of the cuboids. The region below this interface is discarded in the effective representation of the textured wall. We impose three different boundary conditions on the interface; (i) no-slip condition corresponding to a smooth wall, (ii) slip condition (5) and (iii) the TR model, including also the transpiration velocity (12). Since the cuboids have the same geometry in both $x$ and $y$ directions and they are aligned with the chosen coordinate system, we have $L_{xx}=L_{yy}$, $M_{xx}=M_{yy}$ and $L_{xy}=M_{xy}=0$. The values of $L_{xx}$ and $M_{xx}$ are provided in table 4. They were computed a priori by solving the fundamental problems (FP1–FP2) for a cuboid roughness element in the interface unit cell using the procedure described in section 2.

The simulations are carried out under conditions which lead to $Re_{\tau}=u_{\tau}\delta/\nu\approx 180$, where $u_{\tau}$ is the friction velocity. For all simulations, we impose a constant mass flux; the driving pressure gradient is continuously adjusted. For more simulation details, see appendix C.2. Figure 11(a) shows the time and space averaged mean velocity profiles for the three effective simulations in plus (or wall) units. The averaged velocity on the interface $U_{s}^{+}$ is subtracted from the mean flow $U^{+}$, such that all profiles have zero value at the crest plane of the roughness. We note a downward shift of the logarithmic part of mean velocity profiles that increases from (i) to (iii). The logarithmic part of the mean flow can be represented by

[TABLE]

where $\kappa=0.392$ and $B=4.48$ (Millikan, 1939; Luchini, 2017). Moreover, $\Delta U^{+}$ is the roughness function that quantifies the shift in the mean velocity profile.

From figure 11(a), we observe that the slip velocity boundary condition produces a relatively small shift ($\Delta U^{+}_{s}$) compared to the simulations where both slip and transpiration velocity are imposed. The transpiration induce an additional shift $\Delta U^{+}_{2}$ which is in fact larger than $\Delta U^{+}_{s}$, despite that the transpiration velocity is formally a higher-order boundary condition. This illustrates the sensitivity of the turbulent channel flow to the transpiration velocity as recognised in earlier studies (Jiménez et al., 2001; Orlandi & Leonardi, 2006; García-Mayoral & Jiménez, 2011).

We also carried out DNS using an immersed boundary method to resolve the flow around cuboids (Breugem et al., 2006).

Figure 11(a) shows the time and space averaged mean flow of the geometrically resolved DNS. In order to make a comparison with effective simulations, we subtracted the mean slip velocity at the crest plane for the DNS as well. We observe that there is a relatively good agreement with the mean velocity profile of the TR model. A good agreement is also observed for the rms velocity fluctuations (figure 11b).

The friction Reynolds number $Re_{\tau}$ for all simulations is given in table 4. For channel flow with rectangular cuboid roughness elements, the friction velocity (and thus the skin-friction drag) at the rough wall is larger compared to that of a smooth wall (Orlandi & Leonardi, 2006). However, effective simulations without the transpiration predict a reduction in the skin-friction drag: this can be observed by a smaller $Re_{\tau}$ for the simulations denoted as “Slip only”. In contrast, the TR model is able to modify the near wall turbulence in a correct way, predicting the roughness-induced drag increase.

Additional insight into the role of the transpiration can be gained from the total wall stress, given by

[TABLE]

where $\rho$ is the fluid density. Here, the over-bar denotes time and space averaged quantities, and $(\bm{\cdot})^{\prime}$ represent turbulent fluctuating quantities. The first term ($\tau_{V}$) is the viscous stress, and the second term ($\tau_{R}$) is the Reynolds stress. For channel flow with slippage only, the wall normal fluctuations are zero ($u_{z}^{\prime}=0$) at the wall. Hence, equation (42) simplifies to $\tau_{w}=\tau_{V}$ and only the viscous stress is modified by the boundary condition. In contrast, when also the transpiration condition is imposed, $u_{z}^{\prime}\neq 0$, which allows for a direct modification of the Reynolds stress at the wall. More thorough analysis of transpiration velocity effect on the turbulence can be found in papers by Gómez de Segura et al. (2018); García-Mayoral et al. (2019).

The TR model is based on the creeping flow or linearity assumption (A1). Naturally the model is expected to work, if the roughness size is below the size of the viscous sub-layer. The roughness considered here is so-called transitional roughness (Jiménez, 2004), which is slightly larger than the viscous sub-layer ($k^{+}\approx 7$). The model, however, still provides a reasonable approximation of the DNS, despite the fact that there are inertial effects present in the flow between textured elements. For even larger roughness elements, inertial effects inside the textured surface will become more important, rendering the fundamental problems (FP1–FP5) inaccurate. From our experience, the TR model will fail for roughness heights around $k^{+}\approx 20$.

An empirical model for describing substantially larger surface textures would be required to capture non-linear effects, such as sweeps and ejections in the boundary layer (Breugem et al., 2006) and transpiration velocity due to the pressure fluctuations (García-Mayoral et al., 2019). This is however out in scope of this work.

5 Comparison to the multi-scale expansion

We compare the effective conditions in the TR model (5,12,24,26) with a set of conditions obtained from a multi-scale expansion (MSE). Appendix A provides the essential components of the derivation. A more in depth analysis and derivations of the MSE model will be presented in a separate paper. Here, we focus our attention on two aspects; first, one-to-one comparison between the terms of the TR model and the MSE and, second, the accuracy of the TR model compared to the MSE.

5.1 One-to-one comparison between the TR model and the MSE

An expansion up to the order $\mathcal{O}\left(\epsilon\right)$ results in the following boundary conditions (75) for the free fluid velocity at the interface,

[TABLE]

where $\bm{\mathrm{\tau}}_{e}=\mu\left(\partial_{z}u_{x}+\partial_{x}u_{z},\partial_{z}u_{y}+\partial_{y}u_{z},2\partial_{z}u_{z}\right)$ is the fluid shear stress (containing the symmetric part) of the free fluid and $U^{s}$ is a characteristic magnitude of the slip velocity. The tensor $\bm{\mathrm{{M}}}_{e}$ is a 3rd-rank tensor with $81$ elements, while the tensors $\bm{\mathrm{{L}}}_{e}$ and $\bm{\mathrm{{K}}}_{e}$ are both 2nd-rank tensors with $9$ elements. The double dot operation between a 3rd-rank tensor $\bm{\mathrm{{A}}}$ and a 2nd-rank tensor $\bm{\mathrm{{B}}}$ is defined as $\bm{\mathrm{{A}}}:\bm{\mathrm{{B}}}=A_{ijk}B_{jk}$, where summation over repeating indices is implied.

The leading order term ${\mathcal{O}\left(U^{s}\right)}$ of the velocity condition (43) is the slip term. The slip tensor $\bm{\mathrm{{L}}}$ in the TR model (5) corresponds to the upper left $2\times 2$ block of $\bm{\mathrm{{L}}}_{e}$. In other words, the TR model contains the leading order MSE term with a reduced expression of the shear stress tensor. From mass conservation arguments, it can be shown that the last row and column of tensor $\bm{\mathrm{{L}}}_{e}$ appearing in (43) are zero. Consequently there is no transpiration velocity at ${\mathcal{O}\left(U^{s}\right)}$.

There are two higher order ${\mathcal{O}\left(\epsilon U^{s}\right)}$ terms in the velocity condition (43); a Darcian term related to the pore pressure gradient and a term related to the variation of the shear stress. In the TR model, the Darcy contribution to the tangential velocity components at the interface is neglected. In other words, the TR model has the Darcy contribution only for the wall-normal transpiration component (24). One may again show from mass conservation, that the last row of $\bm{\mathrm{{K}}}_{e}$ in equation (43) is equal to the last row of the interior permeability tensor $\bm{\mathrm{{K}}}$. Consequently, the Darcy term $u_{z}^{-}$ in the TR model (25) corresponds to $(\bm{\mathrm{{K}}}_{e}\cdot\nabla p^{-}/\mu)\cdot\hat{n}$.

Finally, the term related to the variation of the shear stress in equation (43) is compared with the transpiration boundary conditions in the TR model (12). We observe that the TR model contains some of the next order terms arising from the variation of the shear stress; the TR model contains only $2$nd-rank tensor corresponding to tangential shear stress variations in the tangential directions, which is in contrast to the full $3$rd-rank tensor in the MSE corresponding to all shear stress variations in all directions.

The boundary condition for the pressure derived using MSE (76) is

[TABLE]

where $\Delta P$ is a characteristic magnitude of a pressure drop in the system. Here, $\bm{\mathrm{b}}$ and $\bm{\mathrm{a}}$ are vectors and $\bm{\mathrm{{C}}}$ is a 2nd-rank tensor.

The leading term ${\mathcal{O}\left(\Delta P\right)}$ of expression (44) induces a pressure jump proportional to the shear stress. It can be shown through mass conservation and force balance that the last element in $\bm{\mathrm{b}}$ (corresponding to shear stress $2\mu\,\partial_{z}u_{z}$) is always equal to $-1$. Therefore this term can be transferred to the left hand side and grouped together with “$-p$” to yield the total free fluid stress, as appearing in the TR model (26). Further, we assert that the first two elements of the vector $\bm{\mathrm{b}}$ in expression (44) corresponds to the friction factor $\bm{\mathrm{f}}^{(2)}$ in equation (26). This can be confirmed by replacing the shear stress in expression (44) with $\bm{\mathrm{{L}}}^{-1}\bm{\mathrm{u}}_{t}$. Consequently there is a full overlap of the leading order pressure jump terms between the TR and the MSE. The two higher order ${\mathcal{O}\left(\epsilon\Delta P\right)}$ terms in the pressure condition (44) are considered next. The vector $\bm{\mathrm{a}}$ corresponds to the friction factor $\bm{\mathrm{f}}^{(1)}$, which is confirmed by replacing the pore pressure in equation (44) with $\bm{\mathrm{{K}}}^{-1}\bm{\mathrm{u}}^{-}$. Finally, the term corresponding to variations of the shear stress is completely neglected for the pressure condition in the TR model.

In table 5 we group the different terms of the TR model according to the order at which the corresponding terms emerge in the multi-scale expansion. The slip velocity in the TR model contains only the leading order term with a reduced shear stress vector, while the transpiration condition and the pressure condition contains all the leading order contributions and some of next order corrections.

5.2 Accuracy of the TR model compared to the MSE

We compare the flow and the pressure in the lid-driven cavity with a textured and a porous surface computed from fully resolved ensemble averaged DNS with three different effective models; (i) the zeroth-order model, containing only leading order terms from the MSE condition (43–44), (ii) the first-order model, containing all terms from the MSE condition (43–44) and (iii) the TR model.

The transpiration velocity along the interface – at the same textured wall as discussed in section 2.3 – is shown in figure 12(a). The pressure distribution along the vertical slice – of the same layered geometry from section 3.3 – is shown in figure 12(b). The interface in both cases is located at $z_{i}=0.1\,l$. As expected, the TR model has a clear improvement over the zeroth-order MSE model. It is observed that the first-order MSE model and the TR model nearly overlap for both the transpiration and the pressure fields. We can thus conclude that the TR model provides nearly as good approximation as the MSE model, but with significantly reduced complexity.

As a final remark, we note that the TR model is not mathematically (or asymptotically) fully consistent. If mathematical rigour is sought, all the next order terms for the slip velocity, transpiration velocity and pressure condition should be taken into account. However, as we have demonstrated using the turbulent channel flow, a small transpiration velocity gives a similar magnitude of shift in the mean velocity profile as the large slip velocity (figure 11a). Consequently, small corrections of slip velocity would not change results significantly, while the introduction of the small transpiration velocity leads to notable modification of the mean flow. In other words, it can be argued that – for certain physical problems – some higher-order terms are more important than others.

6 Conclusions

The TR model provides a set of accurate effective boundary conditions suitable for modelling free-fluid interacting with rough and porous surfaces. These boundary conditions can be incorporated into computational fluid dynamics codes and thus enable investigations of how fully anisotropic textures and porous materials interact with external fluids. We have validated the TR model for creeping flows over textured and porous surfaces. Moreover, based on our investigations, we suggest to place the interface as close to the solid structures as possible without intersecting the solid structures.

The values of the coefficients within the TR model – the slip length tensor $\bm{\mathrm{{L}}}$, the transpiration length tensor $\bm{\mathrm{{M}}}$, the resistance vectors $\bm{\mathrm{f}}^{(1)}$ and $\bm{\mathrm{f}}^{(2)}$ – provide direct information of transfer of mass and momentum that can be expected when the surface interacts with a flow. These coefficients can be computed for any surface topology using five fundamental Stokes problems and are properties only of the surface itself. To understand the significance of this, we make an analogy with a bulk porous material. In this case, the permeability is an established measure that characterises the ability of the porous material to transmit fluids. This measure is invaluable in understanding and designing porous materials for applications. In a similar way, we believe that slip, transpiration and resistance coefficients have a physical meaning on their own.

The TR model is derived by making three assumptions; (i) a creeping flow near the surface texture; (ii) scale separation between the texture size and the flow length scale and; (iii) a repeating surface geometry. This means that the proposed model is a type of homogenised interface condition. By using a formal multi-scale expansion (MSE), we have identified the theoretical orders of different terms present in the TR model. We have also shown how the TR model and MSE model predictions compare using the lid-driven cavity as a test bed. Based on this comparison, we have justified that certain higher-order terms of the MSE model can be neglected. This results in a set of effective boundary conditions that are much simpler to implement and use compared to the full MSE conditions.

For configurations where there exists an intrinsic hydrodynamic sensitivity to wall-normal velocity, the transpiration velocity may become as important as the slip velocity, although the former is, from an asymptotic viewpoint, a higher order correction. We have shown one such example here, namely the turbulent channel flow, for which the friction at the rough wall has a direct contribution from wall-normal velocity fluctuations. Similarly, the resistance terms in the normal stress balance condition for porous media can be shown to have different sizes via scaling analysis, but the relevance of the terms can only be determined when the targeted application is taken into consideration. Here, we have shown that a so-called layered porous material need a higher-order resistance coefficient in order to physically capture the “layering effect”. In nature, there is an abundance of porous materials with inhomogeneous layers; one example is the otolith structure inside human ear, which is part of our vestibular apparatus. Otoliths (calcium carbonate crystals) are located on top of a gel membrane, in which hairy sensory structure is located. The transfer of external fluid into these types of complex materials thus requires the higher order description based on transpiration length and resistance coefficients.

The generalisation of the TR model to elastic and poroelastic surfaces is relatively straightforward by applying the model locally, fixed to the displacing solid (Lācis, 2019). Furthermore, the TR model can also be used for curved interfaces using a coordinate transformation, provided that the curvature of the interface is larger than characteristic surface length $l$. Another interesting direction is the extension of the TR model by considering the slip length, the transpiration length and the resistance coefficients as spatially varying, time dependent as well as flow dependent (for example, shear or Reynolds number dependent). Indeed, extensions of this kind may open up exciting modelling opportunities in applied problems including turbulent flows, heat transfer, nutrition transport, etc.

Acknowledgements

U.L. and S.B. acknowledge funding from Swedish Research Council (INTERFACE center and grant nr. VR-2014-5680). S.Y. acknowledges funding from the European Union’s Horizon 2020 research and innovation programme under the Marie Sklodowska-Curie grant agreement number 708281. S.P. acknowledges funding from the Swiss National Science foundation (project nr. P2ELP2_181788). Numerical simulations were performed on resources provided by the Swedish National Infrastructure for Computing (SNIC) at PDC and HPC2N. Authors acknowledge Prof. Luca Brandt and Dr. Marco E. Rosti for sharing their turbulent code and assisting with implementation of boundary conditions. Authors thank Prof. Alessandro Bottaro for his valuable comments about computing the transpiration length.

Appendix A Multi-scale analysis for boundary conditions

In this appendix, we provide a derivation of the boundary condition terms at different orders using multi-scale expansion (MSE). The derivation follows the approach previously used by Lācis & Bagheri (2016).

A.1 Dimensionless Navier-Stokes equations

The starting point is incompressible dimensional Navier-Stokes (NS) equations, defined in all space filled by fluid. These equations are rendered dimensionless using relationships

[TABLE]

where primed variables are dimensionless. Here, $U^{s}$ is characteristic slip velocity near the surface, $\Delta P$ is characteristic pressure drop in the system and $t_{s}$ is characteristic time scale. Using expression (45), the dimensionless Navier-Stokes equations read

[TABLE]

To obtain these equations, we have used an estimate $l\Delta P/\mu\sim U^{s}$. Reynolds and Strouhal numbers are defined as

[TABLE]

respectively. Here, $\rho$ is fluid density and $\mu$ is fluid dynamic viscosity.

A.2 Fast macroscale flow decomposition

We employ so called fast macroscale flow decomposition in the region, where there is only free fluid. The decomposition for velocity and pressure is

[TABLE]

respectively. Here, $\hat{\bm{\mathrm{u}}}^{\prime+}$ and $\hat{p}^{\prime+}$ are velocity and pressure perturbation caused by textured or porous surface, while $\bm{\mathrm{U}}^{\prime}$ and $P^{\prime}$ is flow field assumed to obey NS equations exposed to homogeneous no-slip condition $\bm{\mathrm{U}}^{\prime}=0$ at the artificial interface with the surface ($z=z_{i}$). In the surface, we assume that there is only slow flow denoted by $\hat{\bm{\mathrm{u}}}^{\prime-}$ and $\hat{p}^{\prime-}$ satisfying same NS equations. Gathering resulting governing equations for flow variables $\hat{\bm{\mathrm{u}}}^{\prime\pm}$ and $\hat{p}^{\prime\pm}$ we obtain

[TABLE]

where we have set continuity of velocity and stress at the artificial interface. The forcing from macroscopic fast flow is defined as

[TABLE]

while the stress tensors are

[TABLE]

Taking sum of equations (50–54) with governing equations for $\bm{\mathrm{U}}^{\prime}$ and $P^{\prime}$, one recovers the single set of equations for whole domain.

A.3 Multi-scale expansion

To carry out the MSE, we introduce two different dimensionless coordinates, so called macroscale and microscale, as

[TABLE]

respectively. Note that the second coordinate is identical as the one introduced in the section A.1. In these coordinates, there are two derivatives appearing due to chain rule

[TABLE]

where $\nabla_{1}$ and $\nabla_{0}$ corresponds to derivatives with respect to $\bm{\mathrm{\hat{x}}}^{\prime}$ and $\bm{\mathrm{x}}^{\prime}$, respectively. Here, $\epsilon=l/H\ll 1$ is scale separation parameter. In addition, the standard amplitude expansion is employed for perturbation velocity and pressure fields as

[TABLE]

We insert the chain rule (58) and amplitude expansions (59–60) in the governing equations for perturbation velocity (50–54). We assume that we are working with small Reynolds number $Re\leq\mathcal{O}\left(\epsilon^{2}\right)$ and small Strouhal number $St\leq\mathcal{O}\left(1\right)$ and group the terms appearing at different orders.

A.3.1 $\mathcal{O}\left(1\right)$ problem and solution Ansatz

The problem at order $\mathcal{O}\left(1\right)$ reads

[TABLE]

One can observe, that this problem is forced by stress at the interface from the fast macroscopic flow $\bm{\mathrm{U}}^{\prime}$. Therefore we anticipate that the solution for velocity perturbations will take form

[TABLE]

where $\bm{\mathrm{{L^{\prime}}}}_{e}$ is unknown $3\times 3$ tensor field, and $\bm{\mathrm{\tau}}^{\prime}_{e}=\epsilon\left(\nabla_{0}\bm{\mathrm{U}}^{\prime}+\left(\nabla_{0}\bm{\mathrm{U}}^{\prime}\right)^{T}\right)\cdot\bm{\mathrm{n}}$ is viscous stress vector at the surface. The pressure, on the other hand, should take form

[TABLE]

where $\bm{\mathrm{b^{\prime}}}^{\pm}$ is unknown vector field and $p^{\prime-(0)}$ is zeroth order pressure field below the interface containing only macroscale variations. This term is needed to match the free flow pressure $P^{\prime}$.

A.3.2 $\mathcal{O}\left(\epsilon\right)$ problem and solution Ansatz

The problem at order $\mathcal{O}\left(\epsilon\right)$ is obtained by collecting all the terms at appropriate order and inserting solution of fields from the order $\mathcal{O}\left(1\right)$ problem. Resulting equations are

[TABLE]

where we the double contraction between $3\times 3\times 3$ tensor $\bm{\mathrm{{A}}}$ and $3\times 3$ tensor $\bm{\mathrm{{B}}}$ is defined as $\bm{\mathrm{{A}}}:\bm{\mathrm{{B}}}=A_{ijk}B_{jk}$. Here one observes that the problem is forced using volume forcing and mass sources which are proportional to macroscopic gradient of pressure within the surface as well as macroscopic gradient of free flow shear stress. The volume forcing factor standing in front of shear stress variations is defined as $g\left(\bm{\mathrm{{L^{\prime}}}}_{e}^{\pm},\bm{\mathrm{b^{\prime}}}^{\pm}\right)=\bm{\mathrm{b^{\prime}}}^{\pm}\bm{\mathrm{{\delta}}}+2\nabla_{1}\bm{\mathrm{{L^{\prime}}}}_{e}^{\pm}$. Therefore we assume that the solution for next order velocity perturbations is

[TABLE]

where $\bm{\mathrm{{K^{\prime}}}}_{e}^{\pm}$ is unknown $3\times 3$ tensor field, and $\bm{\mathrm{{M^{\prime}}}}_{e}^{\pm}$ is unknown $3\times 3\times 3$ tensor field. The pressure, on the other hand, should take the following form

[TABLE]

where $\bm{\mathrm{a^{\prime}}}^{\pm}$ is an unknown vector field and $\bm{\mathrm{{C^{\prime}}}}^{\pm}$ is unknown $3\times 3$ tensor field.

A.4 Boundary conditions up to order $\mathcal{O}\left(\epsilon\right)$

In this section, we present resulting boundary conditions. The result and derivation holds both for textured and porous surface.

A.4.1 Order $\mathcal{O}\left(1\right)$ boundary conditions

To determine the boundary condition with an error of $\mathcal{O}\left(\epsilon\right)$, we use Ansatzes for solution of $\mathcal{O}\left(1\right)$-problem (64–65) and insert them back into the amplitude expansion (59–60).

For velocity condition, we take the surface average at the interface, neglect all higher order terms, relate the no-slip solution with the corrected flow field and go back to the dimensional quantities to obtain

[TABLE]

where the dimensionless $3\times 3$ tensor $\bm{\mathrm{{\mathcal{L}}}}=\langle\bm{\mathrm{{L}}}^{\prime+}_{e}\rangle_{i}$ is a surface average of the microscale tensor field $\bm{\mathrm{{L}}}^{\prime+}_{e}$.

For pressure jump condition, we follow the same approach as for velocity, but instead of single surface average we take instead two volume averages (one in the free fluid region, and second in the interior region) and take a difference for estimating the pressure jump condition. For the estimation of the pressure jump, we have to remember also that the pressure in the free fluid contain both the perturbation and the no-slip solution (49). Taking all this into account, we obtain

[TABLE]

where vector with $3$ components $\bm{\mathrm{\bar{b}}}=\langle\bm{\mathrm{b}}^{\prime-}\rangle^{-}-\langle\bm{\mathrm{b}}^{\prime+}\rangle^{+}$ is the difference in volume averaged fields $\bm{\mathrm{b}}^{\pm}$.

A.4.2 Order $\mathcal{O}\left(\epsilon\right)$ boundary conditions

For boundary condition with an error of $\mathcal{O}\left(\epsilon^{2}\right)$, we repeat the same procedure as in section A.4.1 by taking additionally into account the Ansatzes for solution of $\mathcal{O}\left(\epsilon\right)$-problem (71–72).

Velocity boundary condition then becomes

[TABLE]

where the dimensionless $3\times 3$ tensor $\bm{\mathrm{{\mathcal{K}}}}=\langle\bm{\mathrm{{K}}}^{\prime+}_{e}\rangle_{i}$ is a surface average of the microscale tensor field $\bm{\mathrm{{K}}}^{\prime+}_{e}$ and the dimensionless $3\times 3\times 3$ tensor $\bm{\mathrm{{\mathcal{M}}}}=\langle\bm{\mathrm{{M}}}^{\prime+}_{e}\rangle_{i}$ is the surface average of the microscale tensor field $\bm{\mathrm{{M}}}^{\prime+}_{e}$.

The pressure boundary condition with next order corrections becomes

[TABLE]

where dimensionless vector with $3$ components $\bm{\mathrm{\bar{a}}}=\langle\bm{\mathrm{a}}^{\prime-}\rangle^{-}-\langle\bm{\mathrm{a}}^{\prime+}\rangle^{+}$ is the difference in volume averaged fields $\bm{\mathrm{a}}^{\pm}$ and dimensionless $3\times 3$ tensor $\bm{\mathrm{{\mathcal{C}}}}=\langle\bm{\mathrm{{C}}}^{\prime-}\rangle^{-}-\langle\bm{\mathrm{{C}}}^{\prime+}\rangle^{+}$ is the difference in volume averaged fields $\bm{\mathrm{{C}}}^{\prime\pm}$. The boundary conditions (75–76) are presented in the main paper (43–44) and discussed in the context of the TR model.

Appendix B Equivalence to a two-domain description

In this appendix, we elaborate on how the Dirac delta function is used for surface forcing in equations (6–7). In essence, this notation is equivalent to having a two-domain description and enforcing continuity of velocities and jump in stress, as appearing in multi-scale expansion (62) and also as reported in work by Lācis & Bagheri (2016).

Let us consider the equations (6–7) in three different regions; (i) above the interface, (ii) below the interface and (iii) in a close vicinity of the interface. Introducing plus notation for variables above the interface and minus notation for variables below the interface, we rewrite (6–7) as

[TABLE]

where we have used the property of the Dirac delta that it is zero everywhere except at $\hat{z}=\hat{z}_{i}$. Now, however, there is additional surface that requires new boundary conditions. From the continuity of solution at the single domain, we state that the first condition at $\hat{z}=\hat{z}_{i}$ must be continuity of velocities, i.e., $\hat{\bm{\mathrm{u}}}^{+}=\hat{\bm{\mathrm{u}}}^{-}$. This is, however, not sufficient and a stress condition is also required. Consider equations (6–7) in the near vicinity of the interface $\hat{z}=\hat{z}_{i}$ before introduction of plus and minus notation. We rewrite the momentum equation (6) by making use of the Newtonian fluid stress tensor $\bm{\mathrm{{\Sigma}}}$ as

[TABLE]

where $\bm{\mathrm{{I}}}$ is identity tensor. We integrate the equation in the interface normal direction from $\hat{z}_{i}-\delta\hat{z}$ to $\hat{z}_{i}+\delta\hat{z}$ and get

[TABLE]

where we have explicitly integrated out the divergence part in $z$-direction as well as Dirac delta function, which gives one as long as the integration interval is encapsulating the interface coordinate from both sides. Assuming that the fluid stress tensor varies smoothly in $x$ and $y$ directions, we have

[TABLE]

as the integration interval shrinks to zero $\delta\hat{z}\rightarrow 0$. Consequently the moment equation integral (80) can be rewritten as

[TABLE]

where fluid stress tensors are evaluated at the interface from the negative side $\left.\bm{\mathrm{{\Sigma}}}\right|_{\hat{z}_{i}^{-}}$ (function of $\hat{p}^{-}$ and $\bm{\mathrm{\hat{u}}}^{-}$) and from the positive side $\left.\bm{\mathrm{{\Sigma}}}\right|_{\hat{z}_{i}^{+}}$ (function of $\hat{p}^{+}$ and $\bm{\mathrm{\hat{u}}}^{+}$). By taking into account the orientation of unit normal of the interface $\bm{\mathrm{e}}_{z}=\bm{\mathrm{n}}$, we have

[TABLE]

which is the final boundary condition needed for the two domain formulation of the interface cell problem.

Appendix C Description of numerical methods

In this appendix, we describe more details of the numerical methods we have used through this work.

C.1 Laminar flow

To discretise the domain for the rough configuration (figure 6), we use node spacing $\Delta s_{t}=0.05\,l$ at the top wall and $\Delta s_{b}=0.005\,l$ at the surface texture. The very fine mesh for the rough configuration was chosen to make sure that the variations in the final averaged data are not due to resolution issues. The node spacing for effective textured simulations is $\Delta s_{t}=\Delta s_{b}=0.05\,l$.

To discretise the domain for the porous configurations (figure 9), we use node spacing $\Delta s_{t}=0.125\,l$ at the top wall and $\Delta s_{b}=0.05\,l$ at the porous structures. For effective simulations of porous configurations, we use node spacing $\Delta s_{t}=\Delta s_{b}=0.083\,l$ at all walls.

We solve the incompressible Stokes equations with finite element solver FreeFEM++ (Hecht, 2012). We choose a monolithic approach, i.e., the momentum and continuity equations are treated at the same time, which leads to natural treatment of boundary conditions, which mix velocities and pressures.

C.2 Turbulent flow

For turbulent simulations, we use a periodic uniform grid in the spanwise ($x$) and streamwise ($y$) directions, and a stretched grid in wall normal ($z$) direction. The solver is based on a staggered grid with a third order Runge-Kutta time scheme combined with a splitting technique. A semi-explicit sub-iteration scheme is used to implement the model velocity boundary conditions. The results are made dimensionless using the viscous length and the friction velocity defined as $u_{\tau}/\nu$ and $u_{\tau}$, respectively. These quantities are computed at the cuboids crest plane to be compared to the TR model. The mesh spacings for the geometry resolving simulation in viscous units are $\Delta x^{+}=\Delta y^{+}=1.794$, $\Delta z_{w}^{+}=0.358$, and $\Delta z_{c}^{+}=3.096$; subscripts ’w’ and ’c’ denote near-wall and channel center line respectively. Within the textured layer, we use a constant near-wall mesh spacing $\Delta z_{w}^{+}$. We also impose 5 layers of uniform $\Delta z_{w}^{+}$ at both channel walls. For the effective (and smooth wall) simulations, the mesh spacing is $\Delta x^{+}=7.8$, $\Delta y^{+}=5.1$, $\Delta z_{w}^{+}=0.3$, and $\Delta z_{c}^{+}=4.8$.

The Reynolds number for the geometry-resolved case is defined as $Re_{rough}=U(\delta+k/2)/\nu=2856$ based on the bulk velocity $U=1$, the half-channel height and the kinematic viscosity $\nu$. As the height of the domain is truncated for the TR model, we consider the Reynolds number computed on the reduced domain, cutting off the roughness part. It becomes $Re=U\delta/\nu=2839.2$ with $U_{c}=1.014$ the bulk velocity computed on this truncated domain. This Reynolds number is kept constant for all DNSs, which leads to $Re_{\tau}=u_{\tau}\delta/\nu\approx 180$ for smooth channel, where $u_{\tau}$ is the friction velocity.