A.1 Slip-modified friction law and inertial correction

Firstly, this section gives the derivation of the slip-modified friction law (2.10), together with its correction in the presence of near-wall inertial effects (3.17). The baseline relation follows the angular-momentum defect theory for turbulent TC flow with wall slip (Panton Reference Panton2007; Srinivasan et al. Reference Srinivasan, Kleingartner, Gilbert, Cohen, Milne and McKinley2015). The correction is then introduced by noting that, at larger $w^+$ , the near-wall convective transport effectively displaces the inner overlap layer. In the reduced-order model, this effect is represented by a near-wall averaged turbulent-viscosity correction $\widetilde {\nu _t}$ .

Following Srinivasan et al. (Reference Srinivasan, Kleingartner, Gilbert, Cohen, Milne and McKinley2015), the mean angular momentum in the bulk is written in defect form as

(A1) \begin{align} \langle \gamma \rangle = \langle \gamma _0 \rangle + \left (\frac {u_\tau }{V_i}\right )\langle \gamma _1 \rangle , \end{align}

where $\gamma =r u_\theta$ , $V_i=r_i\omega _i$ is the inner-cylinder velocity, and $\langle \gamma _0 \rangle$ and $\langle \gamma _1 \rangle$ denote the leading-order bulk angular momentum and its defect correction, respectively. Near the inner wall, the overlap-layer variable is defined as

(A2) \begin{align} \gamma ^+ = \frac {r_i(V_i-U_s)-\langle \gamma \rangle }{r_i u_\tau }, \end{align}

with $U_s$ the effective slip velocity. Matching the wall-layer and outer representations, following the argument of Panton (Reference Panton2007), leads to a logarithmic overlap relation at the inner wall,

(A3) \begin{align} (1-\varGamma _0)\frac {V_i}{u_\tau } = \frac {1}{\alpha _g}\ln \left (\frac {u_\tau d}{\nu }\right )+C_1+C_2+U_s^+, \end{align}

where $\varGamma _0=\langle \gamma _0\rangle /(r_iV_i)$ and $U_s^+=U_s/u_\tau$ . At the outer wall, where a no-slip condition is imposed, the corresponding relation is

(A4) \begin{align} \varGamma _0\frac {V_i}{u_\tau } = \frac {1}{\beta _g}\ln \left (\frac {u_\tau d}{\nu }\right )+C_3+C_4. \end{align}

Summing (A3) and (A4) with ${V_i}/{u_\tau }=\sqrt {{2}/{C_f}}$ gives

(A5) \begin{align} \sqrt {\frac {2}{C_f}} = C_{g,I}\ln \left (Re_i\sqrt {\frac {C_f}{2}}\right )+C_{g,\textit{II}}+U_s^+, \end{align}

where the geometry-dependent constants are $C_{g,I}={1}/{\alpha _g}+{1}/{\beta _g}$ and $ C_{g,\textit{II}}=C_1+C_2+C_3+C_4$ . Equation (A5) is the baseline slip-modified friction law quoted in (2.10) of the main text.

The baseline relation (A5) assumes that the inner-wall overlap layer is controlled by diffusion and the overall wall slip alone. As shown in §§ 3.1 and 3.2, this assumption breaks down at larger $w^+$ , where the near-wall convective transport becomes non-negligible. In the present reduced-order model, this missing transport is represented by a near-wall averaged turbulent-viscosity correction $\widetilde {\nu _t}$ . The correction is not introduced as a local pointwise modification of the wall law; rather, it is treated as an integral effect associated with the overall displacement of the inner overlap region. Accordingly, $\widetilde {\nu _t}$ denotes a near-wall averaged quantity over the modulated layer, rather than a local viscosity profile.

With this interpretation, the effective inner-wall overlap scale becomes

(A6) \begin{align} \delta _e=\frac {\nu +\widetilde {\nu _t}}{u_\tau }, \end{align}

while the outer-wall scaling remains unchanged, since the patterned slip and its associated convective correction are confined to the inner-wall side. Accordingly, the inner-wall overlap variable is shifted from $y^+=yu_\tau /\nu$ to $y/\delta _e$ , which corresponds to a downward shift of the logarithmic overlap region, enhanced momentum transport, and thus, an increase in drag. Replacing $\nu$ by $\nu +\widetilde {\nu _t}$ only in the inner-wall matching relation (A3), while retaining the outer-wall relation (A4), gives

(A7) \begin{align} \sqrt {\frac {2}{C_f}} = C_{g,I}\ln \left (Re_i\sqrt {\frac {C_f}{2}}\right )+C_{g,\textit{II}}+U_s^+ - \frac {C_{g,I}}{2}\ln \left (1+\frac {\widetilde {\nu _t}}{\nu }\right ). \end{align}

At sufficiently high Reynolds numbers, or in the limit of weak curvature, the inner- and outer-wall logarithmic coefficients both approach values related to the von Kármán constant $\kappa$ . It is therefore natural to approximate the inner-wall contribution as one half of the total logarithmic coefficient, i.e. $1/\alpha _g \approx C_{g,I}/2$ . Here, $\widetilde {\nu _t}$ denotes the representative near-wall turbulent viscosity, obtained by averaging $\nu _{t,\textit{eff}}$ at $y=0$ and $y=y_{adv}$ , as the overall shift of the velocity profile results from the cumulative effect of momentum mixing across the near-wall layer. Moreover, results in table 2 show that $\widetilde {\nu _t}$ is insensitive to the precise definition of the layer-averaging process. The endpoint average $\widetilde {\nu _t}^{(e)}$ , spatial average $\widetilde {\nu _t}^{(s)}$ and convective-flux-weighted average $\widetilde {\nu _t}^{(w)}$ yield comparable values, whereas the local maximum $\nu _{t,\textit{eff}}(y_{adv})$ gives a noticeably larger estimate. In this table, differentiating the relation (A7) with respect to $\chi _\nu =\widetilde {\nu _t}/\nu$ gives the local sensitivity

(A8) \begin{align} \frac {\partial \ln C_f}{\partial \chi _\nu } = \frac {C_{g,I}} {(1+\chi _\nu ) \left [\sqrt {2/C_f}+C_{g,I}\right ]}, \end{align}

which provides a first-order estimate of how uncertainty in $\widetilde {\nu _t}$ propagates to the predicted skin-friction coefficient.

Table 2.

Sensitivity of the representative near-wall turbulent viscosity and the propagated effect on the predicted skin-friction coefficient for three representative cases. The endpoint average is used in the main model.

Table 2. Long description

A table comparing parameters and angular-momentum transport for different cases. The table has 3 rows and 8 columns. Column headers are: Case, v_t ⁺ (e)/v_t ⁺, v_t ⁺ (s)/v_t ⁺, v_t ⁺ (w)/v_t ⁺, v_t ⁺ (max)/v_t ⁺, (ΔC_f/C_f)_s (%), (ΔC_f/C_f)_w (%), (ΔC_f/C_f)_max (%). Row 1: Case W16O30, 0.204, 0.201, 0.208, 0.287, -0.05, +0.06, +1.20. Row 2: Case W28O40, 0.548, 0.510, 0.538, 0.814, -0.38, -0.10, +2.65. Row 3: Case W40O50, 0.838, 0.764, 0.798, 1.254, -0.57, -0.31, +5.31.

Equation (A7) reduces to the baseline law (A5) when $\widetilde {\nu _t}\rightarrow 0$ . For finite $\widetilde {\nu _t}$ , the additional logarithmic term lowers the effective drag-reducing contribution of slip, thereby capturing the counteracting role of near-wall inertial transport against overall wall slip.

A.2 Critical condition and crossover range

Having established the physical basis of the two scaling branches (2.9) and (3.16), the remaining task is to quantify the finite range over which the slip-dominated $3/4$ scaling gives way to the inertia-affected $1/4$ scaling, which is depicted as the first grey-shaded region in figure 10.

The transition between the two slip-velocity scalings occurs when the wall-induced convective transport is no longer a small correction to molecular viscous transport. The relative magnitude of this contribution is measured quantitatively by $\nu _{t,\textit{eff}}/{\nu }$ . The limit $\nu _{t,\textit{eff}}/\nu \ll 1$ corresponds to a slip-dominated near-wall response, whereas $\nu _{t,\textit{eff}}/\nu =O(1)$ indicates that the inertial contribution must be retained. This order-of-magnitude condition does not imply an exact or sharp threshold at $\nu _{t,\textit{eff}}/\nu =1$ , since the inertial correction enters the friction relation logarithmically. Therefore, $\nu _{t,\textit{eff}}/\nu =O(1)$ should be regarded as a qualitative order-of-magnitude condition, and a more quantitative operational criterion is required to estimate the finite transition range.

From the corrected friction law (A7), the displacement associated with near-wall inertial transport is

(A9) \begin{align} \Delta U_I^+ = \frac {C_{g,I}}{2} \ln \left ( 1+\frac {\widetilde {\nu _t}}{\nu } \right ). \end{align}

For $\Delta U_I^+\ll 1$ , the inertial term remains a small perturbation to the slip-only friction relation. By contrast, when $\Delta U_I^+=O(1)$ , the inertial contribution produces an order-unity velocity profile displacement in wall units. It is therefore no longer asymptotically small and must be retained in the reduced description.

To obtain a quantitative estimate of the transition range, we take

(A10) \begin{align} \Delta U_I^+=1, \qquad \frac {\widetilde {\nu _t}}{\nu } = \exp \left (\frac {2}{C_{g,I}}\right )-1 \simeq 0.445, \end{align}

where $C_{g,I}=5.43$ . This value provides a practical quantitative realisation of the qualitative condition $\widetilde {\nu _t}/\nu =O(1)$ .

The corresponding texture scales are estimated using the mixing-length relation (3.11). Applying this condition (A10) to the fitted $3/4$ and $1/4$ branches gives

(A11) \begin{align} w_{3/4}^+\simeq 94.6, \qquad w_{1/4}^+\simeq 68.9. \end{align}

Here, $w_{3/4}^+$ estimates where the upward extrapolation of the slip-only branch loses self-consistency because the inertial correction has reached an order-unity level, whereas $w_{1/4}^+$ estimates where the downward extrapolation of the inertia-affected branch first becomes self-consistent. These two values therefore provide the upper and lower estimates of the finite crossover range, respectively.

Generally, the two estimates bracket a finite crossover interval,

(A12) \begin{align} w_{\textit{trans}}^+ \approx 70{-}100 . \end{align}

This interval is further supported by the intersection of the two scaling branches in figure 10, $w_{\textit{cross}}^+ \simeq 85$ , which lies within the same range. Hence, the transition should be understood as a gradual overlap between two scaling descriptions, rather than as a sharp critical point. Below this range, the slip-only $(w^+)^{3/4}$ scaling is appropriate. Above this range, the wall-induced secondary vortices provide an appreciable inertial correction, and the $(w^+)^{1/4}$ scaling becomes the appropriate reduced description. Data within $70\lt w^+\lt 100$ are therefore treated as transitional and are not used to define either asymptotic scaling branch.

A.3 Upper applicability bounds

The upper applicability range of the reduced model is set by two underlying assumptions. First, the modulated near-wall layer is assumed to remain governed primarily by the wall-induced secondary vortices and separated from the canonical near-wall shear-production region, namely (2.5). In terms of the production components, this requires

(A13) \begin{align} \left \langle \mathcal{P}_{r\theta }\right \rangle _V^\lt \ll \left \langle \mathcal{P}_{\theta z}\right \rangle _V^\lt , \end{align}

where $\langle \boldsymbol{\cdot }\rangle _V^\lt$ denotes a volume average over the modulated layer. Second, the gas–liquid interface is assumed to remain approximately rigid, pinned and shear free, without persistent deformation, depinning or wetting. The corresponding dynamical and interfacial limits are estimated below.

A.3.1 Near-wall dynamical limit

The former dynamical upper bound is estimated by comparing the modulated-layer thickness $y_{adv}^+$ with the location $y_p^+$ of the canonical near-wall production peak. The scaling derivation assumes that the dominant production within the modulated layer is $\mathcal{P}_{\theta z}$ , generated by the slip/no-slip heterogeneity and the associated secondary vortices. This assumption remains appropriate while $y_{adv}^+$ lies below the conventional shear-production region.

As $w^+$ increases, however, the modulated layer thickens and progressively approaches this region, as shown in figure 4(b–d). Using the vortex scale in (3.6), with prefactor $A_y=0.126$ , together with the canonical estimate $y_p^+=O(10)$ (Kim et al. Reference Kim, Moin and Moser1987; Huisman et al. Reference Huisman, Scharnowski, Cierpka, Kähler, Lohse and Sun2013), the dynamical limit is defined by

(A14) \begin{align} y_{adv}^+ \simeq y_p^+. \end{align}

For the present smooth-wall reference case, figure 4(a) gives $y_p^+\simeq 10.8$ , yielding

(A15) \begin{align} w_P^+ \simeq \left ( \frac {y_p^+}{A_y} \right )^{4/3} \simeq 380. \end{align}

For $w^+\gtrsim w_P^+$ , the conventional radial–azimuthal production $\mathcal{P}_{r\theta }$ becomes appreciable within the modulated layer. The near-wall response is then no longer governed solely by the secondary-vortex scales used in the $1/4$ derivation. Additional shear-induced mixing and form-drag-like losses are therefore expected to weaken the effective slip response and lead to an overprediction of $U_s^+$ by the reduced scaling.

A.3.2 Interfacial stability limit

A second upper bound arises from the stability of the pinned gas–liquid interface. Even if the additional contribution of $\mathcal{P}_{r\theta }$ were modelled, the prescribed shear-free interface would become questionable once the turbulent pressure load became comparable to the available capillary restoring pressure.

Following the Young–Laplace stability framework for superhydrophobic interfaces (Rastegari & Akhavan Reference Rastegari and Akhavan2019), the relevant condition is

(A16) \begin{align} \Delta P_{\textit{load}} \lesssim \Delta P_{\sigma ,max } \simeq C_\sigma \frac {\sigma }{g}, \end{align}

or, equivalently,

(A17) \begin{align} \Delta P_{\textit{load}}^{+}g^{+}We_\tau \lesssim C_\sigma , \end{align}

where $g=(1-\phi _s)w$ is the gas-interface width and $We_\tau =\rho u_\tau \nu /\sigma$ is the friction Weber number.

The dimensionless coefficient $C_\sigma$ represents the maximum capillary pressure supported by the pinned meniscus relative to the scale $\sigma /g$ . It therefore depends on the limiting contact angle, meniscus morphology and local pinning geometry. For a trench with sidewalls, the Young–Laplace pinning condition gives $C_\sigma =-2\cos \theta _a$ , where $\theta _a$ is the advancing contact angle measured through the liquid phase (Yu et al. Reference Yu, Li, McClelland, del Campo Melchor, Lee, Lee and Kim2023). The value $\theta _a=120^\circ$ adopted by Rastegari & Akhavan (Reference Rastegari and Akhavan2019) gives $C_\sigma =1$ , whereas re-entrant edges can increase the limiting capillary pressure to $2\sigma /g$ , corresponding to the upper value $C_\sigma =2$ (Yu et al. Reference Yu, Li, McClelland, del Campo Melchor, Lee, Lee and Kim2023). Surface treatments and chemically patterned pinning sites, including the hydrophilic ridge tips employed here, can also modify the accessible range of pinning angles and the associated capillary resistance (Yu et al. Reference Yu, Li, McClelland, del Campo Melchor, Lee, Lee and Kim2023; Cao et al. Reference Cao, Liu, Zhang, Wu, You, Li, Duan and Huang2024; Liu et al. Reference Liu2024). Accordingly, we adopt $1\leqslant C_\sigma \leqslant 2$ as a literature-based sensitivity range, rather than assigning a unique value to the present surface.

Using $w^{+}=wu_\tau /\nu$ , (A16) gives the implicit stability estimate

(A18) \begin{align} \left (w_{\textit{stab}}^{+}\right )^2 \simeq \frac {C_\sigma \sigma w} {\rho \nu ^2(1-\phi _s) \Delta P_{\textit{load}}^{+} \left (w_{\textit{stab}}^{+},Re_\tau \right )}. \end{align}

Notably, the interfacial upper bound is therefore not universal, but depends on the physical texture width, gas fraction, surface tension, pressure statistics and pinning geometry.

In the pressure-PDF-based estimate of Rastegari & Akhavan (Reference Rastegari and Akhavan2019), the root-mean-square wall-pressure fluctuation is

(A19) \begin{align} p_{\textit{rms}}^{\prime +} = \left \{ 2.32\ln \left [(g^+)^{3/4}\right ] + 2.31\ln (Re_\tau ) - 14 \right \}^{1/2}. \end{align}

For a conservative instantaneous plastron-failure criterion, the pressure load in (A17) is represented by the upper-tail value

(A20) \begin{align} \Delta P_{\textit{load,inst}}^{+} \simeq 4p_{\textit{rms}}^{\prime +}, \end{align}

which accounts for rare, intense pressure events capable of causing local depinning or wetting. In addition, as discussed by Yu et al. (Reference Yu, Li, McClelland, del Campo Melchor, Lee, Lee and Kim2023), this estimate corresponds to the limiting case of sufficiently deep and long gas-filled grooves, for which the gas phase provides negligible compensating pressure fluctuations.

By contrast, the present objective is to estimate the applicability of a time-averaged shear-free model, rather than an instantaneous plastron-failure threshold. The forcing is therefore characterised by the root-mean-square level fluctuation

(A21) \begin{align} \Delta P_{\textit{load}}^{+} = p_{\textit{rms}}^{\prime +}, \end{align}

which assumes that the mean friction response is altered primarily by sufficiently strong and persistent interfacial deformation or depinning. It should therefore be interpreted as an indicative time-averaged applicability limit, rather than as a strict instantaneous-failure condition.

For $1\leqslant C_\sigma \leqslant 2$ , substituting (A19) into (A18) gives the case-specific ranges

(A22) \begin{align} w_{\textit{stab}}^{+}\simeq 213-293,\quad 283-389,\quad 339-466,\quad 387-532 \end{align}

for $w=1.6,\ 2.8,\ 4.0$ and $5.2 \ \textrm {mm}$ , respectively. Over the parameter range considered here, the unweighted mean of the four interval midpoints is approximately $363$ . For compact comparison with the other applicability estimates, we therefore adopt

(A23) \begin{align} w_{\textit{stab}}^{+}\simeq 360 \end{align}

as a representative central estimate here. This value is not a universal threshold; the case-specific intervals should be retained because of the strong dependence of the stability limit on the physical texture width and pinning geometry.

Above the corresponding case-specific interfacial limits, persistent interface deformation, depinning or fluctuations may invalidate the prescribed shear-free boundary condition. The resulting interfacial motion can enhance near-wall momentum transport and reduce the effective slip velocity $U_s^{+}$ , or equivalently, the log-law shift $\Delta U^{+}$ , consistent with the observations of Rastegari & Akhavan (Reference Rastegari and Akhavan2019).

A.3.3 Combined upper limit

The dynamical estimate $w_P^+\simeq 380$ is close to the interfacial estimate $w_{\textit{stab}}^+\simeq 360$ . This proximity indicates that enhanced cross-plane mixing and persistent meniscus fluctuations may become important over a similar large- $w^+$ range. Both mechanisms increase near-wall momentum transport and weaken the effective slip response.

The dominant limiting mechanism may nevertheless depend on the physical texture size. For larger grooves, shear-induced mixing and form-drag-like effects may arise first, whereas for smaller grooves, the capillary-stability limit may be reached before a comparable dynamical penalty develops.

Combining the dynamical and interfacial estimates, a representative upper applicability limit here is taken as

(A24) \begin{align} w_{\textit{upper}}^+ \simeq 360, \end{align}

with an indicative upper-limit transition interval $w_{\textit{upper}}^+ \in [340,400]$ . Within and above this range, indicated by the second grey-shaded region in figure 10, the interfaces may undergo strong and persistent deformation or fluctuations, making the assumption of a stable, shear-free interface questionable. The resulting interfacial motion can enhance near-wall momentum transport and reduce the effective slip velocity $U_s^+$ , or equivalently, the log-law shift $\Delta U^+$ , as also reported by Rastegari & Akhavan (Reference Rastegari and Akhavan2019).

Taken together, the scaling behaviour is divided into the following physically distinct ranges:

(A25) \begin{align} \begin{aligned} w^+ &\lesssim 70, && \text{viscous-dominated $3/4$ scaling};\\ 70 \lesssim w^+ &\lesssim 100, && \text{finite crossover range};\\ 100 \lesssim w^+ &\lesssim 360, && \text{inertia-affected $1/4$ scaling};\\ w^+ &\gtrsim 360, && \text{breakdown-prone regime}. \end{aligned} \end{align}

Beyond the upper-limit range, persistent interfacial fluctuations and additional form-drag-like losses are expected to weaken the effective slip response. The $1/4$ scaling should therefore not be extrapolated into this regime.

B.1 Energy-budget framework

The instantaneous field is decomposed into a base and a residual component as

(B1) \begin{align} [\boldsymbol{u},p]=[\boldsymbol{U},P]+[\boldsymbol{u}_d,p_d], \end{align}

where the choice of base flow $\boldsymbol{U}$ depends on the analysis. In the physical-space analysis, $\boldsymbol{U}=\overline {\boldsymbol{u}}$ is the time-averaged field. In the spectral analysis, $\boldsymbol{U}=\langle \boldsymbol{u}\rangle _{A,t}$ is the mean field averaged over a cylindrical surface and time; quantities defined using this base flow are denoted by the superscript $(\boldsymbol{\cdot })^\dagger$ . Subtracting the averaged Navier–Stokes equations from the full equations gives the residual equations

(B2a) \begin{gather} \frac {\partial \boldsymbol{u}_d}{\partial t} +\boldsymbol{U}\boldsymbol{\cdot }\boldsymbol{\nabla }\boldsymbol{u}_d +\boldsymbol{u}_d\boldsymbol{\cdot }\boldsymbol{\nabla }\boldsymbol{U} +\boldsymbol{u}_d\boldsymbol{\cdot }\boldsymbol{\nabla }\boldsymbol{u}_d -\langle \boldsymbol{u}_d\boldsymbol{\cdot }\boldsymbol{\nabla }\boldsymbol{u}_d\rangle = -\frac {1}{\rho }\boldsymbol{\nabla }p_d+\nu \boldsymbol{\nabla} ^2\boldsymbol{u}_d, \\[-28pt] \nonumber \end{gather}

(B2b) \begin{gather} \boldsymbol{\nabla }\boldsymbol{\cdot }\boldsymbol{u}_d=0, \\[0pt] \nonumber \end{gather}

where $\langle \boldsymbol{\cdot }\rangle$ denotes the averaging operator associated with the chosen base flow.

Taking the inner product of (B2) with $\boldsymbol{u}_d$ and averaging yields the Reynolds–Orr equation (Schmid & Henningson Reference Schmid and Henningson2001) quoted in the main text,

(B3) \begin{align} \underbrace {\partial E_d/\partial t + \boldsymbol{U} \boldsymbol{\cdot }\boldsymbol{\nabla }E_d}_{\bar {D} E_d/\bar {D} t} & = \underbrace {-\langle \boldsymbol{u}_d \otimes \boldsymbol{u}_d \rangle \colon \boldsymbol{\nabla }\boldsymbol{U}}_{\mathcal{P}} \underbrace {-\langle \boldsymbol{u}_d \otimes \boldsymbol{u}_d \colon \boldsymbol{\nabla }\boldsymbol{u}_d \rangle }_{\mathcal{N}} \nonumber \\& \quad - \underbrace {\nu \langle \boldsymbol{\nabla }\boldsymbol{u}_d \colon \boldsymbol{\nabla }\boldsymbol{u}_d \rangle }_{\mathcal{D}} - \boldsymbol{\nabla }\boldsymbol{\cdot }\langle p_d \boldsymbol{u}_d \rangle + \nu \boldsymbol{\nabla} ^2 E_d , \end{align}

where $E_d=\langle \boldsymbol{u}_d\boldsymbol{\cdot }\boldsymbol{u}_d\rangle /2$ . Here $\mathcal{P}$ denotes production by the base-flow shear, $\mathcal{N}$ denotes nonlinear redistribution of residual kinetic energy and $\mathcal{D}$ denotes viscous dissipation. Equation (B3) describes the local kinetic-energy balance associated with the wall-induced structures and turbulent fluctuations.

Under statistically stationary conditions, the volume-averaged form of (B3) gives

(B4) \begin{align} \langle \mathcal{P}\rangle _V - \langle \mathcal{D}\rangle _V =0, \end{align}

because the temporal term vanishes, while the pressure-transport and viscous-diffusion terms reduce to boundary contributions that do not enter the net volume balance. The nonlinear term is conservative under volume averaging and redistributes energy without contributing to the net gain or loss. Equation (B4) shows that, in the integrated sense, production and dissipation act as the source and sink of residual kinetic energy, respectively.

To resolve scale-specific contributions, a generic variable $f(\theta ,z,t;r)$ is Fourier-transformed in the azimuthal, axial and temporal directions as $\widehat {f}(\boldsymbol{k},\varOmega ;r)$ , where $\boldsymbol{k}=[k_\theta ,k_z]^T$ . In the present groove-aligned configuration, we focus on the steady axial family $[k_\theta =0,k_z,\varOmega =0]$ , which captures the wall mode and its harmonics. Applying the Fourier transform to (B2), multiplying by $\widehat {\boldsymbol{u}_d}^{\,*}$ and averaging yields the spectral kinetic-energy budget (Pope Reference Pope2000)

(B5) \begin{align} \boldsymbol{U} \boldsymbol{\cdot }\boldsymbol{\nabla }\widehat {E_d} & = \underbrace {- \langle \widehat {\boldsymbol{u}_d} \otimes \widehat {\boldsymbol{u}_d}^{\,*}\rangle \colon \boldsymbol{\nabla }\boldsymbol{U}}_{\widehat {\mathcal{P}}} \underbrace {- \langle \widehat {\boldsymbol{u}_d}^{\,*}\boldsymbol{\cdot }\widehat {\boldsymbol{u}_d \boldsymbol{\cdot }\boldsymbol{\nabla }\boldsymbol{u}_d} \rangle }_{\widehat {\mathcal{N}}} \nonumber \\& \quad - \underbrace {\nu \langle \boldsymbol{\nabla }\widehat {\boldsymbol{u}_d} \colon \boldsymbol{\nabla }\widehat {\boldsymbol{u}_d}^{\,*} \rangle }_{\widehat {\mathcal{D}}} - \langle \widehat {\boldsymbol{u}_d}^{\,*} \boldsymbol{\cdot }\boldsymbol{\nabla }\widehat {p_d} \rangle + \nu \boldsymbol{\nabla} ^2 \widehat {E_d}. \end{align}

Here, $\widehat {\mathcal{P}}$ represents modal production from the mean shear, $\widehat {\mathcal{N}}$ represents nonlinear redistribution among modes and $\widehat {\mathcal{D}}$ represents viscous dissipation of modal kinetic energy.

After volume averaging, and neglecting the pressure-transport and viscous-diffusion terms that are subdominant to the leading production, transfer and dissipation terms (Aubry et al. Reference Aubry, Holmes, Lumley and Stone1988; Symon et al. Reference Symon, Illingworth and Marusic2021), the spectral balance reduces to

(B6) \begin{align} \langle \widehat {\mathcal{P}}\rangle _V + \langle \widehat {\mathcal{N}}\rangle _V - \langle \widehat {\mathcal{D}}\rangle _V =0. \end{align}

Thus, $\langle \widehat {\mathcal{P}}\rangle _V$ and $\langle \widehat {\mathcal{D}}\rangle _V$ act as the source and sink of modal kinetic energy, whereas $\langle \widehat {\mathcal{N}}\rangle _V$ redistributes energy across modes without altering the total volume-averaged amount.

To reveal the detailed pathways embedded in $\widehat {\mathcal{N}}$ , the volume-averaged nonlinear transfer is further decomposed as

(B7) \begin{align} \langle \widehat {\mathcal{N}}^\dagger (k_z) \rangle _V = \sum _s \langle \widehat {\mathcal{M}}^\dagger (k_z,s) \rangle _V , \end{align}

where the first argument $k_z$ denotes the recipient mode and the second argument $s$ denotes the source mode. The mode-to-mode transfer tensor is defined as

(B8) \begin{align} \widehat {\mathcal{M}}^\dagger (k_z,s) = -\widehat {\boldsymbol{u}_d}^{\,*}(k_z)\boldsymbol{\cdot }\big [ \widehat {\boldsymbol{u}_d}(k_z-s) \boldsymbol{\cdot }\widehat {\boldsymbol{\nabla }\boldsymbol{u}_d}(s) \big ]. \end{align}

Equivalently, since the physical velocity field is real, $\widehat {\boldsymbol{u}_d}^{\,*}(k_z)=\widehat {\boldsymbol{u}_d}(-k_z)$ . This quantity measures the signed transfer associated with the triad $(k_z,s,k_z-s)$ rooted in the nonlinear term $\boldsymbol{u}\boldsymbol{\cdot }\boldsymbol{\nabla }\boldsymbol{u}$ (Domaradzki & Rogallo Reference Domaradzki and Rogallo1990; Alexakis et al. Reference Alexakis, Mininni and Pouquet2005; Duvvuri & McKeon Reference Duvvuri and McKeon2015). A positive value of $\widehat {\mathcal{M}}^\dagger (k_z,s)$ means that the recipient mode $k_z$ gains energy from the source mode $s$ , whereas a negative value means that the recipient mode $k_z$ loses energy to the source mode $s$ .

B.2 Energy transfer properties

The volume-averaged mode-to-mode transfer satisfies two basic conservative properties (Ding et al. Reference Ding, Chung and Illingworth2025). First, the transfer between any two modes is antisymmetric,

(B9) \begin{align} \langle \widehat {\mathcal{M}}^\dagger (k_z,s) \rangle _V = - \langle \widehat {\mathcal{M}}^\dagger (s,k_z) \rangle _V , \end{align}

so that energy gained by one mode is lost by another. Second, there is no net self-transfer,

(B10) \begin{align} \langle \widehat {\mathcal{M}}^\dagger (k_z,k_z) \rangle _V =0 . \end{align}

These properties ensure that the nonlinear term redistributes energy among modes without changing the total volume-averaged energy, consistent with (B6). For the regional analysis, the volume-averaged mode-to-mode transfer is separated into the modulated layer and the outer region:

(B11) \begin{align} & \langle \widehat {\mathcal{M}}^\dagger (k_z,s) \rangle _V^{\lt } = \frac {1}{V_\lt } \int _{y\leqslant y_{adv}} \widehat {\mathcal{M}}^\dagger (k_z,s)\,\text{d}V , \nonumber \\ & \langle \widehat {\mathcal{M}}^\dagger (k_z,s) \rangle _V^{\gt } = \frac {1}{V_\gt } \int _{y\gt y_{adv}} \widehat {\mathcal{M}}^\dagger (k_z,s)\,\text{d}V . \end{align}

Here, $\langle \widehat {\mathcal{M}}^\dagger (k_z,s)\rangle _V^{\chi }$ with $ \chi \in \{\lt ,\gt \}$ denotes either of the two regional averages.

The classification of a pairwise modal transfer as forward/downscale or inverse/upscale depends on both the sign of $\langle \widehat {\mathcal{M}}^\dagger (k_z,s)\rangle _V^{\chi }$ and the relative magnitude of $k_z$ and $s$ . Forward, or downscale, transfer denotes energy redistribution from a lower-wavenumber, larger-scale mode to a higher-wavenumber, smaller-scale mode. Inverse, or upscale, transfer denotes the opposite direction. Therefore,

(B12) \begin{align} (k_z-s) \langle \widehat {\mathcal{M}}^\dagger (k_z,s) \rangle _V^{\chi } \gt 0 \quad \text{for forward/downscale transfer} \end{align}

and

(B13) \begin{align} (k_z-s) \langle \widehat {\mathcal{M}}^\dagger (k_z,s) \rangle _V^{\chi } \lt 0 \quad \text{for inverse/upscale transfer}. \end{align}

The transfer matrix above identifies direct pairwise modal exchange. However, an individual matrix element does not by itself represent a complete cascade.

As a result, to quantify organised transfer pathways beyond the direct pairwise transfer, we introduce two additional diagnostics: the transfer-distance distribution $\mathcal{H}_q^\chi (\Delta k)$ and the cutoff-integrated modal-transfer flux $\varPi _q^\chi (K)$ . Here $ \chi \in \{\lt ,\gt \}$ denotes the radial region and $q=F,I$ denotes forward/downscale and inverse/upscale transfer, respectively (Symon et al. Reference Symon, Illingworth and Marusic2021; Ding et al. Reference Ding, Chung and Illingworth2025). In general, the transfer-distance distribution identifies the modal separation and locality of the transfer, while the cutoff-integrated flux quantifies the cumulative scale-to-scale consequence of an organised transfer pattern.

First, given the modal separation $\Delta k=|k_z-s|$ with the total modal-exchange activity in region $\chi$ written as

(B14) \begin{align} A^{\chi } = \sum _{k_z\gt s} \big | \langle \widehat {\mathcal{M}}^\dagger (k_z,s) \rangle _V^{\chi } \big |, \end{align}

the forward/downscale and inverse/upscale-transfer-distance distributions are then defined as

(B15) \begin{align} \mathcal{H}_F^{\chi }(\Delta k) = \frac {1}{A^{\chi }} \sum _{\substack {k_z\gt s\\ k_z-s=\Delta k}} \big | \langle \widehat {\mathcal{M}}^\dagger (k_z,s) \rangle _V^{\chi } \big | \mathbf{1}\big [ (k_z-s) \langle \widehat {\mathcal{M}}^\dagger (k_z,s) \rangle _V^{\chi } \gt 0 \big]\end{align}

and

(B16) \begin{align} \mathcal{H}_I^{\chi }(\Delta k) = \frac {1}{A^{\chi }} \sum _{\substack {k_z\gt s\\ k_z-s=\Delta k}} \big | \langle \widehat {\mathcal{M}}^\dagger (k_z,s) \rangle _V^{\chi } \big | \mathbf{1}\big [ (k_z-s) \langle \widehat {\mathcal{M}}^\dagger (k_z,s) \rangle _V^{\chi } \lt 0 \big ], \end{align}

where $\mathbf{1}[\boldsymbol{\cdot }]$ is the indicator function. Thus, $\mathcal{H}_F^\chi$ and $\mathcal{H}_I^\chi$ measure the fractions of the total modal-exchange activity carried by forward/downscale and inverse/upscale transfers at a given modal separation $\Delta k$ . Small $\Delta k$ , especially $\Delta k=1$ , corresponds to local modal transfer, whereas a peak at a larger prescribed separation, such as $\Delta k=k_w$ , indicates non-local transfer locked to a specific periodicity.

Figure 11 shows $\mathcal{H}_q^\chi (\Delta k)$ for SMO50 and W28O50 in the modulated and outer regions. Inside the modulated layer, the two cases exhibit distinct transfer characteristics. In SMO50 the forward component is dominated by local transfer at small separation $\Delta k=1$ (figure 11 a). By contrast, in W28O50 the structured wall produces a pronounced non-local peak near $\Delta k_*\simeq k_w$ in the forward component, confirming a wall-harmonic transfer pathway locked to the wall periodicity (figure 11 c). However, outside the modulated layer, the two cases share a common organised feature: the inverse/upscale component $\mathcal{H}_I^\gt (\Delta k)$ exhibits a finite-separation peak at $\Delta k=\Delta k_*$ associated with an energetic scale. In SMO50, this peak occurs near $\Delta k_*\simeq k_s/2$ , indicating a subharmonic separation related to streaks (figure 11 b). In W28O50 the corresponding peak occurs near $\Delta k_*\simeq k_w$ , indicating that the outer-region upscale feedback is associated with the wall-induced secondary vortices (figure 11 d). Notably, the energy-containing streak signature appears at a subharmonic scale rather than exactly at $k_s$ , which is consistent with the fact that streaks are self-organised flow structures that may naturally merge or modulate over longer wavelengths, unlike the wall-imposed secondary vortices that are externally locked to $k_w$ ; the volume-averaged transfer may further reflect their growth and merging stages. Generally, although the near-wall transfer changes from local streak-related transfer in SMO50 to non-local wall-harmonic transfer in W28O50, the outer region in both cases exhibits finite-separation inverse/upscale transfer that links energetic structures to larger-scale motions.

Moreover, the dominant peak of $\mathcal{H}_I^\gt (\Delta k)$ at $\Delta k=\Delta k_*$ (figure 11 b,d) provides a quantitative link to the organised inverse-transfer bands marked by the dash-dotted lines in figures 6(d,h). These bands have an approximately constant offset from the diagonal, corresponding to source–recipient pairs with $|s-k_z|=\Delta k_*$ . They are therefore consistent with the definition of $\mathcal{H}_I^\gt (\Delta k_*)$ and can also be interpreted as triadic-compatible transfer patterns,

(B17) \begin{align} \left \{s,\ k_z,\ |s-k_z|=\Delta k_*\right \} , \end{align}

through which energetic modes are connected to larger-scale recipient modes as part of the outer-region cascade-like redistribution. As a result, the two outer-region transfer maps are not pointwise identical, but they share a common organised feature: finite-separation inverse/upscale transfer links energetic structures to larger scales, with the relevant separation set by the streak scale in SMO50 and by the wall mode in W28O50.

Figure 11.

Forward and inverse-transfer-distance distributions, i.e. $\mathcal{H}_{F}^{\chi }$ and $\mathcal{H}_{I}^{\chi }$ , in terms of the modal separation $\Delta k=|k_z-s|$ , for (a,b) SMO50 and (c,d) W28O50, in the modulated and outer regions. The dotted line marks the local-transfer threshold $\Delta _c=2$ , while the dashed lines mark streak-related modes $k_s$ or $k_s/2$ for SMO50, and wall modes $n k_w (n=1,2,\ldots )$ for W28O50.

Figure 12.

Cutoff-integrated modal-transfer fluxes in terms of cutoff wavenumber $K$ for (a,b) SMO50 and (c,d) W28O50, in the modulated and outer regions, including the forward flux $\varPi _F^\chi (K)$ , inverse flux $\varPi _I^\chi (K)$ and net flux $\varPi _{\textit{net}}^\chi (K)=\varPi _F^\chi (K)-\varPi _I^\chi (K)$ . The dashed lines mark streak-related modes $k_s$ or $k_s/2$ for SMO50, and wall modes $n k_w (n=1,2,\ldots )$ for W28O50.

Figure 12. Long description

Panel A: A line graph shows the near-wall flux for SMO50 in the modulated layer. The x-axis represents the cutoff wavenumber K, ranging from 0 to 20. The y-axis represents the near-wall flux Π(K), ranging from -0.05 to 0.20. Three lines are plotted: forward flux (blue circles), inverse flux (red squares), and net flux (yellow triangles). The dashed lines mark streak-related modes k_s/2 and k_s. Panel B: A line graph shows the outer flux for SMO50 outside the modulated layer. The x-axis represents the cutoff wavenumber K, ranging from 0 to 20. The y-axis represents the outer flux Π(K), ranging from 0 to 0.10. Three lines are plotted: forward flux (blue circles), inverse flux (red squares), and net flux (yellow triangles). The dashed lines mark streak-related modes k_s/2 and k_s. Panel C: A line graph shows the near-wall flux for W28O50 in the modulated layer. The x-axis represents the cutoff wavenumber K, ranging from 0 to 20. The y-axis represents the near-wall flux Π(K), ranging from 0 to 0.6. Three lines are plotted: forward flux (blue circles), inverse flux (red squares), and net flux (yellow triangles). The dashed lines mark wall modes kw and 2kw. Panel D: A line graph shows the outer flux for W28O50 outside the modulated layer. The x-axis represents the cutoff wavenumber K, ranging from 0 to 20. The y-axis represents the outer flux Π(K), ranging from 0 to 0.20. Three lines are plotted: forward flux (blue circles), inverse flux (red squares), and net flux (yellow triangles). The dashed lines mark wall modes kw and 2kw.

Second, the cutoff-integrated modal-transfer flux is defined using the same signed pairwise transfer. For a cutoff $K$ , modes with wavenumbers not larger than $K$ form the lower-wavenumber group, whereas modes with wavenumbers larger than $K$ form the higher-wavenumber group. The normalised forward and inverse fluxes are

(B18) \begin{align} \varPi _F^{\chi }(K) = \frac {1}{A^{\chi }} \sum _{\substack {k_z\gt s\\ s\leqslant K\lt k_z}} \big | \langle \widehat {\mathcal{M}}^\dagger (k_z,s) \rangle _V^{\chi } \big | \mathbf{1}\big [ (k_z-s) \langle \widehat {\mathcal{M}}^\dagger (k_z,s) \rangle _V^{\chi } \gt 0 \big ]\end{align}

and

(B19) \begin{align} \varPi _I^{\chi }(K) = \frac {1}{A^{\chi }} \sum _{\substack {k_z\gt s\\ s\leqslant K\lt k_z}} \big | \langle \widehat {\mathcal{M}}^\dagger (k_z,s) \rangle _V^{\chi } \big | \mathbf{1}\big [ (k_z-s) \langle \widehat {\mathcal{M}}^\dagger (k_z,s) \rangle _V^{\chi } \lt 0 \big ]. \end{align}

Here, $\varPi _F^{\chi }(K)$ measures the cumulative transfer from lower to higher wavenumbers across the cutoff, whereas $\varPi _I^{\chi }(K)$ measures the cumulative transfer from higher to lower wavenumbers. The net cutoff-integrated flux is $\varPi _{\textit{net}}^{\chi }(K) =\varPi _F^{\chi }(K)-\varPi _I^{\chi }(K)$ .

Figure 12 shows the corresponding cutoff-integrated modal-transfer fluxes. In the modulated layer of SMO50, the flux varies smoothly with $K$ , indicating a broadly distributed streak-related transfer. By contrast, W28O50 exhibits a clear wall-harmonic modulation: $\varPi _F^\lt (K)$ and $\varPi _{\textit{net}}^\lt (K)$ increase sharply near $K=k_w$ and then approach a plateau towards $2k_w$ . This behaviour indicates that the structured wall injects a wall-harmonic-locked pathway into the cumulative scale-to-scale redistribution. The inverse component is finite, but less than the forward component; the near-wall process is therefore interpreted as harmonic-locked transfer, rather than as a broadband net inverse cascade. Outside the modulated layer, the detailed cutoff dependence remains different between SMO50 and W28O50, but both cases exhibit the same organised tendency: a forward-dominated net flux coexists with finite inverse/upscale feedback. This supports the interpretation that the two outer-region transfer maps are not pointwise identical, but they are similar in their collective cascade-like redistribution towards larger scales through finite upscale feedback.

In summary, the physical-space budget (B3) clarifies the near-wall production and maintenance of coherent structures induced by the patterned boundary. The spectral-space budget (B5) and the mode-to-mode transfer tensor (B8) identify the wall mode, its harmonics and their interaction with large-scale roll motions. The additional diagnostics $\mathcal{H}_q^\chi (\Delta k)$ and $\varPi _q^\chi (K)$ quantify the locality and cumulative direction of the modal transfer. Together they show that the structured wall reorganises the near-wall pathway from a primarily streak-related local transfer in SMO50 to a wall-harmonic-locked non-local transfer in W28O50, while the outer region retains a forward-dominated cascade-like redistribution with finite upscale feedback. These physical interpretations of the modal-transfer diagnostics are summarised in table 3.

Table 3.

Summary of modal organisation and transfer pathways in SMO50 and W28O50. The near-wall modulated layer and the outer region are listed separately to distinguish wall-harmonic transfer from outer-region cascade-like redistribution. The entries connect each physical interpretation with the corresponding evidence from nonlinear transfer $\mathcal{N}$ and $\widehat {\mathcal{N}}^{\dagger }$ , mode-to-mode transfer $\langle \widehat {\mathcal{M}}^\dagger (k_z,s)\rangle _V^{\chi }$ , transfer-distance distributions $\mathcal{H}_q^\chi (\Delta k)$ and cutoff-integrated fluxes $\varPi _q^\chi (K)$ .

Table 3. Long description

A table comparing the near-wall modulated layer and outer region in smooth (SMO50) and patterned (W28O50) cases. The table has four rows and three columns. Column headers are Physics, Feature, and Evidence. Row labels are SMO50 (smooth) Near-wall modulated layer, SMO50 (smooth) Outer region, W28O50 (patterned) Near-wall modulated layer, and W28O50 (patterned) Outer region. Each row describes the physical interpretation, feature, and evidence for nonlinear transfer, mode-to-mode transfer, transfer-distance distributions, and cutoff-integrated fluxes. The table distinguishes wall-harmonic transfer from outer-region cascade-like redistribution.

C.1 Numerical setting

The DNS are performed using the finite-volume framework OpenFOAM v2412. The incompressible Newtonian flow between two coaxial cylinders is solved in cylindrical geometry. The inner cylinder rotates at angular velocity $\omega _i$ , while the outer cylinder is kept stationary. The radius ratio is fixed at $\eta ={r_i}/{r_o}=0.714$ , where $r_i$ and $r_o$ are the inner- and outer-cylinder radii, respectively, with the gap width of $d=r_o-r_i$ . The computational domain uses a six-fold rotational symmetry, corresponding to an azimuthal sector $\Delta \theta ={\pi }/{3}$ and an axially periodic domain with aspect ratio $\varGamma ={L_z}/{d}={2\pi }/{3}$ . This reduced-domain setting follows previous TC DNS practice, where azimuthal subdomains are used to reduce computational cost (Brauckmann & Eckhardt Reference Brauckmann and Eckhardt2013; Spandan et al. Reference Spandan, Verzicco and Lohse2018; Xu et al. Reference Xu, Su, Lan, Zhao, He, Sun and Wang2023; Su et al. Reference Su, Yi, Zhao, Wang, Xu, Wang and Sun2024).

The governing equations are advanced using the PIMPLE algorithm for pressure–velocity coupling. The temporal term is discretised using an implicit Euler scheme, and the spatial terms are discretised using Gauss-type finite-volume schemes with linear interpolation. The pressure equation is solved using a geometric-algebraic multigrid (GAMG) solver with a DIC smoother, while the velocity equation is solved using an iterative solver with a symmetric Gauss–Seidel smoother. The residual tolerance is set to $10^{-6}$ for the flow variables. The time step is adjusted to keep the maximum Courant number below the prescribed limit, so that the temporal resolution remains consistent with the spatial resolution. The mesh is refined near the wall to resolve the viscous sublayer and the near-wall modulated structures induced by the patterned surface. In all cases, the near-wall grid spacing is checked in wall units based on the inner-wall friction velocity, ensuring DNS-level resolution in the viscous region. After the statistically stationary state is reached, the statistics are collected over at least $50$ turnover times $t_c={d}/({r_i\omega _i})$ , which is sufficient for the converged mean transport and energy-budget quantities used in the present analysis (Xu et al. Reference Xu, Su, Lan, Zhao, He, Sun and Wang2023; Su et al. Reference Su, Yi, Zhao, Wang, Xu, Wang and Sun2024).

Table 4.

Parameters, numerical resolutions and dimensionless angular-momentum transport of the smooth- and patterned-wall cases. Here $ \textit{Nu}^\omega _0$ denotes the corresponding reference value (first row from Ostilla et al. (Reference Ostilla, Stevens, Grossmann, Verzicco and Lohse2013) and others from experiments). The cases W28O50z and W28O50 $\theta$ are validation cases with doubled axial and azimuthal computational domains, respectively.

Table 4. Long description

The table presents data for different cases of fluid dynamics experiments. It has 15 rows and 9 columns. The columns are labeled as follows: Case, w/d, w+, Re_i, Ta, N_r x N_theta x N_z, Nu^ω, Nu₀ ^ω, and Err. (percent). Each row provides specific values for these parameters. For example, Row 1: Case SMO10, w/d Smooth, w+ not applicable, Re_i 2.50 times 10 to the power of 3, Ta 9.52 times 10 to the power of 6, N_r x N_theta x N_z 200 x 120 x 120, Nu^ω 6.389, Nu₀ ^ω 6.422, Err. (percent) 0.51. The table captures the sensitivity of the representative near-wall turbulent viscosity and its effect on the predicted skin-friction coefficient for different cases.A table comparing parameters and angular-momentum transport for different cases.

The boundary conditions are imposed as follows. Periodic conditions are applied in the axial and azimuthal directions. The outer cylinder is no-slip. For smooth-wall cases, the inner cylinder is also no-slip. For patterned-wall cases, no-slip conditions are imposed on the solid ridges of the inner cylinder, whereas the gas–liquid interfaces are modelled as shear free and impermeable. This idealised shear-free-interface condition is used to isolate the inertial effect of the patterned slip boundary without resolving interfacial deformation or multiphase dynamics. Such a boundary modelling strategy has been widely used in previous studies of gas–liquid or slip interfaces (Park et al. Reference Park, Park and Kim2013; Crowdy Reference Crowdy2021; Jeganathan et al. Reference Jeganathan, Shannak, Alba and Ostilla-Mónico2023), and is further supported here by the comparison with reference or experimental values in table 4.

The numerical accuracy is assessed using three checks. First, the computed $ \textit{Nu}^\omega (r)$ is required to be nearly independent of the radial location, consistent with the conservation of angular-momentum transport in statistically stationary TC flow. Second, smooth-wall and reference patterned-wall cases are compared with available benchmark or experimental values, as summarised in table 4. For all cases where a reference value is available, the relative difference remains below $3\,\%$ . Third, because the present work focuses on secondary-flow organisation and modal energy transfer, direct axial and azimuthal domain-sensitivity checks are performed for the representative case W28O50, as discussed in Appendix C.3.

C.2 Parameter range and cases of simulations

The simulations are designed to vary the main control parameters while monitoring both the global transport response and the local flow structures. The primary control parameters are the texture size and the driving strength. The texture size is characterised by the groove width $w/d$ and its inner-unit value $w^+$ normalised by viscous scale $\delta _\nu$ . The driving strength is characterised by the inner-cylinder Reynolds number $ \textit{Re}_i$ , or equivalently, by the Taylor number $ \textit{Ta}$ . The main response parameter is the dimensionless angular-momentum transport $ \textit{Nu}^\omega$ , which measures the global torque relative to the laminar value. In addition to this global response, the simulations also provide detailed local fields for the angular-momentum decomposition, TKE budget and modal energy-transfer analysis.

For the patterned-wall cases, the solid fraction is fixed at $\phi _s=1/4$ , while the normalised groove width is varied over $w/d\in \{0.16,\ 0.28,\ 0.40,\ 0.52\}$ . The driving strength in the simulations covers $ \textit{Re}_i\in [2.50\times 10^3,\ 1.50\times 10^4]$ , or equivalently, $ \textit{Ta}\in [9.52\times 10^6,\ 3.43\times 10^8]$ . Accordingly, the wall-unit texture size of the DNS patterned-wall cases lies in

(C1) \begin{align} w^+\in [47.01,\ 217.34], \end{align}

as listed in table 4. This range covers the small-texture, transitional and inertial-affected regimes discussed in the main text. We note that the experimental data and the literature data (Park et al. Reference Park, Park and Kim2013; Rastegari & Akhavan Reference Rastegari and Akhavan2019) used for comparison in figures 2 and 10 may extend beyond this DNS range.

Table 4 summarises the control parameters, numerical resolutions and global angular-momentum transport for all smooth- and patterned-wall simulations. The case names encode the texture size and driving condition. For example, W28O50 denotes a patterned-wall case with $w/d=0.28$ and $\omega _i=50\ \textrm {rad s}^{-1}$ , while SM denotes a smooth-wall benchmark. The additional cases W28O50z and W28O50 $\theta$ are domain-sensitivity checks, in which the axial and azimuthal extents are doubled, respectively. For cases with available reference data, $ \textit{Nu}^\omega _0$ denotes the corresponding reference value from either literature (Ostilla et al. Reference Ostilla, Stevens, Grossmann, Verzicco and Lohse2013; Su et al. Reference Su, Yi, Zhao, Wang, Xu, Wang and Sun2024) or experiments, and the relative error is defined as

(C2) \begin{align} \textrm {err.} =\left | \frac {Nu^\omega -Nu^\omega _0}{Nu^\omega _0} \right |. \end{align}

For the smooth-wall and validation cases, the relative difference remains below $3\,\%$ , confirming the accuracy of the present numerical set-up.

C.3 Domain-sensitivity checks

The baseline simulations are performed in a reduced TC domain with a six-fold rotational symmetry, $\Delta \theta =\pi /3$ , and an axially periodic length $\varGamma =L_z/d=2\pi /3$ . This setting is used to reduce computational cost while retaining the dynamically relevant near-wall and bulk structures. The validity of this reduced domain is assessed from three complementary aspects: comparison with previous TC DNS using similar reduced domains, the texture-induced origin of the secondary vortices and direct domain-sensitivity checks in the present simulations.

Reduced azimuthal domains have been widely used in turbulent TC DNS when the relevant azimuthal correlations are shorter than the full circumference. Previous studies showed that suitable reduced domains can reproduce the essential turbulent dynamics and global torque, provided that the relevant length scales are resolved (Brauckmann & Eckhardt Reference Brauckmann and Eckhardt2013; Su et al. Reference Su, Yi, Zhao, Wang, Xu, Wang and Sun2024). For example, Brauckmann & Eckhardt (Reference Brauckmann and Eckhardt2013) used axially periodic domains with a default aspect ratio $\varGamma =L_z/d\simeq 2$ , which supports one Taylor-vortex pair. Because the full azimuthal circumference is large for $\eta =0.71$ , they used periodically repeated azimuthal subdomains of size $2\pi /n_{\textit{sym}}$ for high- $ \textit{Re}$ simulations, with $n_{\textit{sym}}$ up to 9. Their local angular-momentum-current analysis at $ \textit{Re}_i\simeq 1.69\times 10^4$ was shown in a one-sixth annulus, corresponding to $n_{\textit{sym}}=6$ and $\Delta \theta =2\pi /6=\pi /3$ . A reference calculation for the stationary-outer-cylinder case showed that the torque difference caused by the azimuthal reduction was below $1\,\%$ . More recent interface-resolved TC simulations also used a six-fold rotational sector and $\varGamma =2\pi /3$ , with this choice validated against previous single- and multiphase TC simulations (Su et al. Reference Su, Yi, Zhao, Wang, Xu, Wang and Sun2024). These results support the use of the present reduced sector, provided that the transport and energy-transfer diagnostics are checked.

The secondary vortices analysed here are also not expected to be box-scale artefacts. Their origin is local and texture induced: the axial alternation of slip and no-slip regions produces an axial gradient of the mean azimuthal velocity and, hence, the surface-modulation production $\mathcal{P}_{\theta z}$ . The dominant near-wall mode is therefore locked to the texture wavenumber $k_w$ and its harmonics, rather than to the lowest axial or azimuthal wavenumber permitted by the computational domain. This scale separation is evident in the spectral diagnostics, where the near-wall transfer and cutoff flux exhibit clear features near $K=nk_w$ $(n=1,2,\ldots )$ , while the large-scale roll modes remain at much lower wavenumbers.

Direct domain-sensitivity checks were performed for the representative case W28O50. The baseline calculation uses $\Delta \theta ={\pi }/{3}$ , $\varGamma ={2\pi }/{3}$ . To test the axial extent, the W28O50z case keeps the same azimuthal sector but doubles the axial length to $\varGamma ={4\pi }/{3}$ . To test the azimuthal extent, the W28O50 $\theta$ case keeps the same axial length but doubles the azimuthal sector to $\Delta \theta ={2\pi }/{3}$ , corresponding to $n_{\textit{sym}}=3$ . The resulting global transport values, selected from table 4, are summarised in table 5. The change in $ \textit{Nu}^\omega$ is only $1.26\,\%$ for the doubled axial domain and $0.28\,\%$ for the doubled azimuthal domain, indicating that the global angular-momentum transport is insensitive to the tested domain extensions.

Table 5.

Domain-sensitivity check for W28O50. Here W28O50z doubles the axial extent, whereas W28O50 $\theta$ doubles the azimuthal extent. The difference is computed relative to the baseline W28O50 case.

Table 5. Long description

The table presents data for three cases: W28O50, W28O50z, and W28O50θ. It includes columns for case names, width to diameter ratio (w/d), positive width (w+), Reynolds number (Re_i), grid resolution (Nr x Nθ x Nz), Nusselt number (Nu^ω), and error percentage (Err.%). Row 1: W28O50, 0.28, 109.65, 12500, 500 x 240 x 240, 7.512, no error percentage. Row 2: W28O50z, 0.28, 110.34, 12500, 500 x 240 x 480, 7.607, 1.26 percent. Row 3: W28O50, 0.28, 109.80, 12500, 500 x 480 x 240, 7.533, 0.28 percent.

Figure 13.

Domain-sensitivity check of the TKE budget for W28O50, W28O50z and W28O50 $\theta$ . The three cases show nearly identical profiles of production, dissipation and nonlinear transfer, confirming that the near-wall modulated layer and the relative roles of $\mathcal{P}_{\theta z}$ and $\mathcal{P}_{r\theta }$ are not artefacts of the reduced computational domain.

The local energy-budget diagnostics give the same conclusion. As shown in figure 13, W28O50, W28O50z and W28O50 $\theta$ exhibit nearly identical distributions of $\mathcal{P}^+$ , $\mathcal{D}^+$ , $\mathcal{N}^+$ , $\mathcal{P}_{r\theta }^+$ and $\mathcal{P}_{\theta z}^+$ . In all three cases, the surface-modulation production $\mathcal{P}_{\theta z}$ remains localised in the same near-wall region, while the conventional radial–azimuthal production $\mathcal{P}_{r\theta }$ becomes important farther from the wall. The location of $y_{adv}$ and the relative balance among production, dissipation and nonlinear transfer are therefore essentially unchanged by the domain extension.

Figure 14.

Domain-sensitivity check of the cutoff-integrated modal-transfer flux. The comparison between W28O50 and W28O50z shows that the wall-harmonic-locked near-wall flux structure and the forward-dominated outer-region redistribution with finite upscale feedback are preserved when the axial domain is extended.

The modal-transfer diagnostics are also robust. Since W28O50 $\theta$ differs only weakly from W28O50 in the global transport and local budget statistics, the modal-flux comparison is focused on W28O50z, which tests the sensitivity to the axial extent. Figure 14 compares the cutoff-integrated fluxes $\varPi (K)$ for W28O50 and W28O50z inside and outside the modulated layer. In the modulated layer, both cases show the same wall-harmonic-locked flux structure, with a sharp increase near $K=k_w$ and a plateau extending towards $2k_w$ . Outside the modulated layer, both cases recover a forward-dominated cutoff flux with finite inverse/upscale feedback. The amplitudes vary slightly because the larger domain contains more realisations of large-scale motions, but the identified forward-dominated redistribution with finite inverse feedback is preserved.

Overall, previous reduced-domain TC DNS, the texture-induced origin of the near-wall secondary vortices and the present sensitivity checks consistently show that the global angular-momentum transport, secondary-flow organisation, energy-budget features and modal-transfer pathways are not artefacts of the imposed six-fold rotational symmetry or the finite axial domain.