University of Southern California Dissertations and Theses

Model based design of porous and patterned surfaces for passive turbulence control

(USC Thesis Other)

Model based design of porous and patterned surfaces for passive turbulence control

PDF

Download a page range

Download transcript

Copy asset link

Request this asset

Transcript (if available)

Content Model based design of por ous and p a tterned surf a ces f or p assive tur bulence contr ol b y Andrew Cha v arin A Dissertation Presen ted to the F A CUL TY OF THE GRADUA TE SCHOOL UNIVERSITY OF SOUTHERN CALIF ORNIA In P artial F ulfillmen t of the Requiremen ts for the Degree DOCTOR OF PHILOSOPHY (Mec hanical Engineering) Ma y 2021 Cop yrigh t 2021 Andrew Cha v arin A c kno wledgmen ts My journey and success would not have been possible without the support and help of my colleagues at USC, my friends, and most importantly, my loving family. I would like to personally acknowledge the following people: My close colleagues and friends at USC: Dr. Christoph Efstathiou and Dr. Yohanna Hanna, Mark Hermes, Shilpa Vijay, Vamsi Krishna, Chris Ohh, and Saakar Byahut for their help and advice over the years. I enjoyed the company, the warm sentiments, and the cups of coffee shared over the years. My close and personal friends: Raymond Earl Jackson, Daniel Acevedo, and Tony Estrada. Thank you for believing in me, being there in my time of need, and for the many memories we’ve created over the years. And my mentor Dr. Mitul Luhar for the opportunity to work with you and for your guidance and knowledge you have imparted upon me over the years. Thank you for the patience you have had with me and for believing in me. ii T able of Con ten ts A c kno wledgmen ts ii List Of T ables v List Of Figures vi Nomenclature x Abstract xiv Chapter 1: In tro duction 1 1.1 Near Wall Streaks . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 1.2 Drag Reduction via a Slip Surface . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 1.3 Riblets . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13 1.4 Anistropic Permeable Surfaces . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18 1.5 Drag Breakdown Mechanism . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23 1.6 Contributions and Outline . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26 Chapter 2: Resolv en t F ramew ork 28 Chapter 3: Resolv en t Analysis for T urbulen t Channel Flo w With Riblets 34 3.1 Overview . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 3.2 Modeling Approach . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36 3.2.1 Volume Penalization for Riblets . . . . . . . . . . . . . . . . . . . . . . . . . 36 3.2.2 Modified Resolvent Operator . . . . . . . . . . . . . . . . . . . . . . . . . . 38 3.2.3 Numerical Implementation . . . . . . . . . . . . . . . . . . . . . . . . . . . 44 3.3 Model Predictions and Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47 3.3.1 Near-Wall Cycle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48 3.3.2 Kelvin-Helmholtz Vortices . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54 3.3.3 Riblet Shape Optimization . . . . . . . . . . . . . . . . . . . . . . . . . . . 59 3.4 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61 Chapter 4: Resolv en t-Based Optimization of T w o-Dimensional Riblet Profiles 63 4.1 Overview . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 63 4.2 Geometry Construction via Parametric Bézier Curves . . . . . . . . . . . . . . . . 64 4.3 Optimization Problem . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66 4.4 Nonlinear conjugate gradient method . . . . . . . . . . . . . . . . . . . . . . . . . . 68 4.4.1 Line Search Method . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 71 4.4.2 Feasible Direction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72 4.4.3 Conjugate Direction β, and Restart Procedure . . . . . . . . . . . . . . . . 74 iii 4.4.4 Approximating Derivatives and Stopping Conditions . . . . . . . . . . . . . 76 4.5 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 76 4.6 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 83 Chapter 5: Resolv en t-Based Predictions for T urbulen t Flo w Ov er Anisotropic P ermeable Substrates 85 5.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 85 5.1.1 Contribution and Outline . . . . . . . . . . . . . . . . . . . . . . . . . . . . 86 5.2 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 88 5.2.1 Accounting for Permeable Substrates . . . . . . . . . . . . . . . . . . . . . . 88 5.2.2 Modified Resolvent Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . 91 5.2.3 Boundary and Interface Conditions . . . . . . . . . . . . . . . . . . . . . . . 95 5.2.4 Mean Velocity Profile . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 96 5.2.5 Numerical Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98 5.3 Results and Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99 5.3.1 Mean Velocity Profiles . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99 5.3.2 Near-Wall Resolvent Mode . . . . . . . . . . . . . . . . . . . . . . . . . . . 101 5.3.3 Phase-Corrected Near-Wall Resolvent Mode . . . . . . . . . . . . . . . . . . 109 5.3.4 Spanwise-Coherent Resolvent Modes . . . . . . . . . . . . . . . . . . . . . . 112 5.4 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 117 Chapter 6: Concluding Remarks 120 References 123 iv List Of T ables 4.1 Summary of the various optimization cases executed. Optimal spacing (s + ), base width (s + b ) and riblet length (l + g ) for h s ≈ 0.5 and h s ≈ 0.75. . . . . . . . . . . . . . 83 v List Of Figures 1.1 Scanning electron micrograph (SEM), left figure, of the dermal denticles covering the dorsal fin of the Shortfin mako shark. SEM shows quasi-periodic riblet like structure displaying troughs and ridges, these denticles are oriented in the stream- wise direction [63]. Manufactured riblet, right figure, emulating the properties of shark skin used in the experiments of Bechert, Bruse, and Hage [8]. . . . . . . . . . 2 1.2 Hydrogen bubble visualization of the near wall coherent streaks at y + = 4.5 by Kline et al. [58]. Photograph shows the spanwise organize of the low-speed streaks and the streamwise extend of the streaks. The streamwise and spanwise directions are denoted by x + and z + respectively. . . . . . . . . . . . . . . . . . . . . . . . . . 5 1.3 Planar snapshot aty + ≈ 20 of the near-wall streaks from the numerical simulations of Schoppa and Hussain [104]. The dark contour show the elongated low-speed streaks u ′ < 0. The dimension of the snapshot is (L + x ,L + z ) = (450,1500) and spanwise spacing of λ + z ≈ 100 can be seen here. . . . . . . . . . . . . . . . . . . . . 6 1.4 Simplified cartoon depicting the structure of the coherent near wall streaks and their relationship to the quasi-streamwise vortices (QSV’s). In this carton the low- speed (−u ′ ) and high-speed (+u ′ ) streaks are color coded red and blue respectively. Horizontal gradient indicates the transition between regions of high and low wall shear stress produced by the quasi-streamwise vortices. . . . . . . . . . . . . . . . . 8 1.5 Cartoon depiction of the displacement of the perceived origin by the mean stream- wise flow (center panel) and the perceived origin by quasi-streamwise vortices and thecross-flow(rightpanel). Panelontheleftshowsadepictionofspanwiseperiodic riblets with a spacing and height of s + and h + respectively. . . . . . . . . . . . . . 10 1.6 Drag measurements collected over triangular riblets, panel a), by Bechert et al. [9]. The optimal riblet spacing s + opt are visible and the breakdown in the linear regime can be seen over the various geomtries tested. Panel b) shows the change in drag for a general riblet geometry plotted verses the riblet spacing, displaying the linear regime, the optimal spacing and the breakdown of this regime. Figures adapted from Bechert et al. [9]. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 vi 1.7 Scaled drag reduction curves collected from Bechert et al. [9] shown in panel a). The change in drag is scaled by the slope of the linear drag regime m ℓ and is plotted against the riblet groove cross-sectionℓ + g . Panels b)-d) display the optimal riblet spacing, height and groove cross-sectional area scale. Figures adapted from García-Mayoral and Jiménez [35]. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17 1.8 Measured drag reduction over anisotropic permeable mediums with streamwise anisotropy ratio ϕ xy ≈ 3.6–11.1. Figures adapted from Gómez-de-Segura and García-Mayoral [40] . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 3.1 Meanvelocitypredictions. (a)Predictedslipvelocityasafunctionofribletspacing. Open symbols show predictions for rectangular and trapezoidal riblets. Closed symbols show DNS results from García-Mayoral and Jiménez [35] for rectangular riblets. (b) Mean streamwise velocity profile for the case identified in (a) in red. Solid contours show velocity values in increments of 2u τ . . . . . . . . . . . . . . . . 46 3.2 Gain for resolvent modes resembling the NW cycle as a function of riblet spacing. . 49 3.3 Velocity structure for resolvent modes resembling the NW cycle for riblet spacing: a) s + = 8, b) s + = 21, and c) s + = 33. Figures on the left show the wall normal velocity in the streamwise-spanwise plane at y + ≈ 6. Figures on the right depict the flow structure in the spanwise-wall normal plane. The vectors show in-plane velocities (v,w) while the red and blue shading shows the out-of-plane streamwise velocity (u). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51 3.4 Amplification of spanwise-constant structures as a function of streamwise wave- length and wave speed for riblet spacing: a) s + = 8, b) s + = 21, and c) s + = 33. Red and blue shading represent contours of log 10 (σ k+ ). The black circles show a local maximum. The vertical bar shows where λ + x = 150. . . . . . . . . . . . . . . . 55 3.5 Maximumgainforspanwise-constantresolventmodesresemblingKelvin-Helmholtz vortices as a function of riblet spacing. The inset text shows the streamwise wave- length for the most amplified mode. . . . . . . . . . . . . . . . . . . . . . . . . . . 56 3.6 Predicted flow structure for the most amplified spanwise-constant resolvent mode over riblets with spacing s + = 21. a) Distribution of wall-normal velocity in a streamwise-spanwise plane at y + ≈ 5. b) Predicted flow field above the riblet tips. c) Predicted flow field in the middle of the riblet grooves. d) Predicted flow field at an intermediate spanwise location. For b)-d) the vectors show the in-plane velocities (u,v), while the red and blue shading shows the out-of-plane spanwise velocity (w). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 58 3.7 Search for optimal trapezoidal riblet geometries with height to spacing ratioh/s = 0.5. a) Predictions for the interfacial slip velocity, U + s , as a function of inner- normalized spacing and ridge angle. b) Predicted gain for resolvent modes resem- bling the NW cycle. Solid black contours show values for the inner-normalized groove length scale, ℓ + g . The assumed geometry is shown on the right. . . . . . . . 59 vii 3.8 Distribution of the optimal riblet size ℓ + g , determined from the square root of the groove cross-sectional area for rectangular, triangular, and trapezoidal riblets. The optimum size is identified based on minimum singular values for resolvent modes resembling the NW cycle. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 60 4.1 Symmetric two-dimensional riblet geometries assembled using Bézier curve(s). On each panel the control points are labeled in order and are denoted using ( ) andP i . Panels a-c) show a quadratic, cubic and quartic Bézier curve respectively. Panel d) shows geometry assembled using a composite Bézier curves. . . . . . . . . . . . 65 4.2 Schematic of the two riblet configurations used with our optimization scheme. Panel a) shows a configuration in where base thickness is equal to the riblet spacing and panel b) shows a configuration where the base thickness is independent of the riblet spacing. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68 4.3 Optimized geometry determined using the algorithm outlined in Section 4.4. Panel (a) shows the gain normalized by the smooth wall cases,σ k+ /σ k,s plotted verses it- erations of the optimizer. Panel (b) plots the residualR k at each iteration, panel(c) shows the initial geometry fed into the optimizer. Panels (d) and (e) and shows the geometry produced after 15 and 38 iterations. In panels (c) -(e) the control points are drawn alongside their Bézier curves. . . . . . . . . . . . . . . . . . . . . . . . . 78 4.4 See figure 4.3 for caption. Panels (d) and (e) show the geometry after 8 and 32 iterations. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 79 4.5 Optimal geometries determined for wavelengths corresponding to, a) λ + z = 80, b) λ + z = 90 and c) λ + z = 100. The control points P 0 to P 3 are drawn alongside the Bézier curve they correspond to. The riblet spacing and reduction in the gain are draw above each configuration. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 81 4.6 Optimized geometry determined using the algorithm outlined in Section 4.4, in where the base width is independent of the riblet spacing. In both panels the normalized mode gain, the optimal riblet spacing (s + ) and the base width (s + b ) are shown along the optimal geometry. In panel(a) the maximum height was limited to h = 15 corresponding to an aspect ratio of h s ≈ 0.5 and in panel (b) the maximum was limited to h = 10 corresponding to h s ≈ 0.5. . . . . . . . . . . . . . . . . . . . . 82 5.1 Schematic showing the symmetric channel flow configuration considered in this chapter. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 90 5.2 Predicted slip velocity at the porous interface as a function of streamwise per- meability p K + x . Closed symbols show DNS results from Gómez-de-Segura and García-Mayoral [40]. Open symbols show predictions made using (5.13). Red sym- bolsindicatethecasesforwhichwall-normalpermeabilityfirstexceeds p K + y & 0.4. This is roughly the threshold above which drag reduction performance deteriorates in DNS. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 100 viii 5.3 Predictions for the NW resolvent mode. (a) Comparison between normalized mode gain (σ k,p /σ k,s , black symbols) and drag reduction observed in DNS (∆U + , light gray symbols) by Gómez-de-Segura and García-Mayoral [40], plotted as a function of p K + y . The resolvent-based predictions shown in this panel are obtained using thevelocityprofilesfromDNS.(b)Comparisonbetweennormalizedgainspredicted using the DNS mean profiles (filled symbols) and the synthetic mean profiles (open symbols) plotted as a function of p K + y . (c) Comparison between normalized gains predicted using the DNS mean profiles (filled symbols) and the synthetic mean profiles (open symbols) plotted as a function of p K + x − p K + y . The black and red dashed lines show linear fits to the initial decrease in normalized gain for the DNS- and synthetic-mean predictions, respectively. For all panels, the symbols represent substrates with ϕ xy = 3.6, symbols represent substrates with ϕ xy = 5.5, and symbols represent substrates with ϕ xy = 11.4. . . . . . . . . . . . . . . 102 5.4 NW mode structure predictions made using the DNS velocity profiles for (a,b) the smooth-wall case, (c,d) the substrate with ϕ xy = 5.5 and p K + y = 0.18, and (e,f) the substrate with ϕ xy = 5.5 and p K + y = 0.45, and (g,h) the substrate with ϕ xy = 5.5 and p K + y = 1.0. In (a,c,e,f), the shading shows normalized contours of positive (red) and negative (blue) streamwise velocity normalized by the maximum value. The vectors show the wall-normal and spanwise velocity fluctuations. In (b,d,f,g), the solid lines show the magnitude of the streamwise velocity for this resolvent mode,|u k | and the dashed lines show the wall-normal velocity magnitude multiplied by a factor of ten, 10|v k |. The black lines represent the smooth-wall case while the gray lines represent the permeable substrates. These plots make use of the shifted coordinate ˆ y = 1+y, such that ˆ y = 0 represents the location of the smooth wall or the porous interface. . . . . . . . . . . . . . . . . . . . . . . . . . . 106 5.5 Predictions for the NW resolvent mode with a modified phase speed of c + ≈ 10 + p K + x − p K + z . (a) Comparison between normalized mode gain (σ k,p /σ k,s , black symbols) and drag reduction observed in DNS (∆U + , light gray symbols) by Gómez-de-Segura and García-Mayoral [40], plotted as a function of p K + y . The resolvent-based predictions shown in this panel are obtained using the velocity profiles from DNS. (b) Comparison between normalized gains predicted using the DNS mean profiles (filled symbols) and the synthetic mean profiles (open sym- bols) plotted as a function of p K + y . (c) Comparison between normalized gains predicted using the DNS mean profiles (filled symbols) and the synthetic mean profiles (open symbols) plotted as a function of p K + x − p K + y . The black and red dashed lines show linear fits to the initial decrease in normalized gain for the DNS- and synthetic-mean predictions, respectively. For all panels, the symbols rep- resent substrates with ϕ xy = 3.6, symbols represent substrates with ϕ xy = 5.5, and symbols represent substrates with ϕ xy = 11.4. . . . . . . . . . . . . . . . . . 111 ix 5.6 Normalized gain for spanwise-coherent structures plotted as a function of stream- wise length λ + x and wave speed c + for substrates with an anisotropy ϕ xy = 5.5. These predictions make use of the synthetic mean profile. Panels (a)-(d) repre- sent structures with spanwise wavenumber κ z ≈ 2.3 (λ + z = 500). Panels (e)-(h) represents structures with spanwise wavenumber κ z = 0. Red and blue shading respectively represent mode amplification and suppression relative to the smooth- wall case. Wall-normal permeability increases from left to right. Structures with the largest amplification are labeled with a ( •) marker. . . . . . . . . . . . . . . . . 112 5.7 Flow structure associated with the most amplified spanwise-constant resolvent mode over the permeable substrate with ϕ xy = 5.5 and p K + y = 0.39; see• in figure 5.6(f). The red and blue shaded contours show regions of positive and neg- ative pressure, normalized by the maximum value. . . . . . . . . . . . . . . . . . . 115 5.8 Comparison of resolvent-based gain predictions for spanwise-coherent structures with anisotropy ratios ϕ xy = 3.6 ( ), 5.5 ( ), and 11.4 ( ). The maximum normalized gain obtained for resolvent modes with λ + x ∈ [50,500] and c + ∈ [4,10] is shown as a function of wall-normal permeability for (a) κ z ≈ 2.3 (λ + z = 500) and for (b) κ z = 0 (λ + z =∞). All of these predictions were obtained using the synthetic mean profiles. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 116 x Nomenclature Op erators Sym b ol Description ( ) + = ( )/δ v , length normalized by viscous length scale ( ), time-averaged quantity Sym b ols Sym b ol Description x, y, z, streamwise, wall normal, spanwise coordinate ν, kinematic viscosity ρ, fluid density h † , channel half height t, time p, pressure u, v, w, streamwise, wall normal, spanwise velocity u = (u,v,w), velocity vector u τ , friction velocity δ v =ν/u τ , viscous length xi Sym b ol Description Re τ , friction Reynolds number h, riblet height s, riblet spacing A g , cross-sectional area of riblet groove α, ridge anlge for trapazodial and triangular riblets l g = p A g , riblet length scale based on cross-sectional area ℓ x , streamwise virtual origin ℓ z , transerve virtual origin λ x , λ z , streamwise and spanwise wavelength k = (κ x ,κ z ,ω), individual wavenumber frequency combination κ s = 2π/s, wavenumber corresponding to riblet spacing κ x , κ z , streamwise and spanwise wavenumber ω, frequency c =ω/κ x , phase speed K † , permeability matrix K, dimensionless permeability matrix K x , dimensionless streamwise permeability K y , dimnesinoless wall-normal permeability K z , dimensionless spanwise permeability ϕ xy = √ K x /K y , streamwise to wall-normal anisotropy ratio ϕ xz = √ K x / √ K z , streamwise to spanwise anisotropy ratio xii Sym b ol Description α k , line search step length at step k β k , conjugate parameter at step k A,b, linear inequality and equality constraints A c ,b c , set of active linear constraints g k , gradient of the objective function at step k P i , i th Bézier curve control point p k , search direction at step k T, fesiable direction projection matrix η, parameter space b i,n (ξ), Bernstein basis polynomial f(η), objective function xiii Abstract Many turbulent flows of engineering interest involve rough, porous, or patterned surfaces that can alter the dynamics of the near-wall turbulence. In particular, surfaces such as sharkskin- inspired riblets and anisotropic permeable coatings have been shown to be effective tools for passive turbulence control. Streamwise aligned riblets have demonstrated the ability to reduce turbulent skin friction by as much as 10% in simulations and laboratory experiments. Recent idealized numerical simulations show that streamwise-preferential porous materials have the po- tential to reduce turbulent skin friction by over 20%. In addition to passive flow control for drag reduction, appropriately-designed complex surfaces could also be used to increase turbulent mixing for applications involving chemical reactors and heat exchangers. However, the immense parameter space involved with designing these surfaces makes it prohibitively expensive to op- timize via experiments or high-fidelity simulations. Testing these surfaces at Reynolds numbers relevant to engineering applications imposes additional hardships. In this thesis, we develop a reduced-complexity modeling framework grounded in linear input-output (or forcing-response) analyses of the governing Navier-Stokes equations to consider the effect of porous and patterned surfaces on wall-bounded turbulent flows. Using this framework, we investigate the effects of spanwise periodic riblets and anisotropic permeable substrates for drag reduction in turbulent flows. xiv The reduced-complexity modeling framework builds on the so-called resolvent analysis frame- work. Under the resolvent formulation, the Navier–Stokes equations are interpreted as a forcing- response system: the nonlinear convective term is interpreted as a feedback forcing on the remain- ing linear terms, which generates a velocity and pressure response. A gain-based decomposition of the linear forcing-response transfer function—the resolvent operator—yields highly amplified ve- locity and pressure modes, which can be considered key building blocks of the turbulent flow field. Previous work has shown that these high-gain modes provide substantial insight into turbulence statistics, structure, and control of smooth-walled flows. The forcing-response gain determined from this system is a measure of energy amplification and can serve as a proxy for control per- formance. The effect of spanwise periodic riblets is introduced in the governing equations using an immersed boundary method. A gain-based decomposition of this modified system reproduces important trends observed in previous simulations and experiments over riblets. The gain of a single resolvent mode, which is representative of the streaks and streamwise vortices associ- ated with the energetic near-wall cycle in turbulent flows, reproduces the drag reduction trends observed in previous simulations of flow over blade-like riblets. Moreover, this framework also predicts the appearance of energetic spanwise rollers resembling Kelvin-Helmholtz vortices, which are hypothesized to be responsible for the breakdown in drag reduction performance for larger riblets. Usingthisreduced-complexityframework,wegeneratepredictionsfortheoptimaltwo-dimensional riblet profile. Specifically, we search for the riblet geometry that leads to maximum suppression of the resolvent mode that serves as a surrogate for the near-wall cycle. To search for this optimal geometry we implement a nonlinear conjugate gradient-based optimizer. The riblet profiles are assembled using cubic Bézier splines which allow us to generate smooth and scalable shapes using a finite set of control points. The optimal profiles identified using this framework are in line with xv profiles predicted by theory and those tested in prior experiments. Predictions suggest that the optimal profile has a thin blade-like profile with an optimal spacing of roughly 20 wall units. We also make use of the resolvent-based modeling framework to study the effect of anisotropic porous substrates on wall-bounded turbulent flows. In this case, the resolvent framework is refor- mulated using the Volume-Averaged Navier-Stokes equations in which the effect of the permeable substrate is modeled using a Darcy-type permeability term. Following prior work, we focus on materials with high porosity and identical wall-normal and spanwise permeabilities. Model pre- dictions indicate that streamwise-preferential substrates with high streamwise permeability and low spanwise/wall-normal permeabilities are able to suppress the resolvent mode that serves as a surrogate for the dynamically important near-wall streaks. Further, the observed reduction in mode gain is consistent with trends predicted by theory and observed in numerical simulations. Similar to riblets the appearance of spanwise-constant rollers resembling Kelvin-Helmholtz vor- tices is hypothesized to be responsible for the breakdown of drag reduction performance at high permeabilities. The resolvent framework predicts the appearance of these spanwise rollers as the wall-normal permeability exceeds a threshold value. These efforts reemphasize the similar flow physics that give rise to drag reduction over sharkskin-inspired riblets and streamwise-preferential porous materials. Moreover, this work demonstrates the utility of the resolvent framework in de- veloping targeted reduced-complexity models that can be used to design and optimize passive control techniques for wall-bounded turbulent flows. xvi Chapter 1 In tro duction Many turbulent flows of engineering interest involve flows over complex patterned, rough, and porous surfaces that alter the dynamics of the near-wall turbulence. Depending on the characteristics of such surfaces and the flow conditions, these surfaces can modulate the near- wall turbulence thereby reducing turbulent mixing and leading to a significant decrease in drag. Examples of such drag-reducing surfaces include dermal denticles found on the skin of sharks which feature a quasi-periodic structure of troughs and valleys which are aligned in the stream- wise direction, see figure 1.1(a). The early experiments of Walsh [123] showed that sharkskin inspired surfaces riblets—featuring spanwise periodic, streamwise aligned ribs/protrusion—can reduce turbulent drag. Experiments of Bechert and Bartenwerfer [7] and Bechert, Bruse, and Hage [8] focused on cataloging the drag change over riblets which mimicked the three-dimensional pattern of shark denticles, see figure 1.1. From these experiments, a maximum drag reduction of 10% was observed over thin blade-like riblets. Other surface modifications such as superhydropho- bic coatings and streamwise preferential anisotropic permeable substrates have also shown high potential in reducing turbulent drag. Experiments of Itoh et al. [47] using seal fur—which similar to riblets displays a high streamwise and low transverse compliance—showed a drag reduction of 1 a) b) Figure 1.1: Scanning electron micrograph (SEM), left figure, of the dermal denticles covering the dorsal fin of the Shortfin mako shark. SEM shows quasi-periodic riblet like structure displaying troughs and ridges, these denticles are oriented in the streamwise direction [63]. Manufactured riblet, right figure, emulating the properties of shark skin used in the experiments of Bechert, Bruse, and Hage [8]. 12%. Numerical simulations pursed Gómez-de-Segura and García-Mayoral [40] over anisotropic streamwise preferential permeable substrate led to drag reductions as high as 20%. Although the examples listed have focused on reducing turbulent mixing and thereby reducing turbulent drag, complex porous and patterned surfaces can also lead to an increase in turbulent mixing, which can be beneficial for chemical reactors and heat exchangers [ 32]. Complex surfaces can have negative and costly effects, for example Schultz et al. [ 105] estimate that the surfaces roughness created by biofouling increases the fuel consumption of naval ships by 10% and increase the cost of operation by $1.2M per year per ship. Understanding how the characteristics of these surfaces modulate the near-wall turbulence is imperative for the design of novel passive flow control strategies for turbulent flows. Given the large parameter space, it becomes expensive, either through experiments or through computational simulations, to investigate howthegeometry of suchsurfaces affect turbulentflows. Thus far only general and sparse guidelines exist for the design of such surfaces. To overcome this problem, we present extensions to the resolvent framework of McKeon and Sharma [73] to account for patterned and permeable surfaces. Using this framework we develop low-dimensional 2 and computationally inexpensive predictive models and apply these models to optimize the char- acteristics of these surfaces for passive turbulence control. In the following sections, we describe the effects that these surfaces have in turbulent flows. Specifically, previous work indicates that the efficacy of porous and patterned surfaces in reducing drag is determined by two key features: (i) their effect on the energetic near-wall cycle [ 52, 98, 120], and (ii) the potential emergence of spanwise rollers resembling Kelvin-Helmholtz vortices over such surfaces [35]. Our reduced- complexity models aim to capture these features. The following sections describe the near-wall cycle as well as the mechanisms dictating the drag-reducing performance of riblets and anisotropic porous materials. 1.1 Near W all Streaks Our journey in understanding why passive control strategies such as riblets and anisotropic permeable wall coatings have been successful begins with understanding the nature and impor- tance of the energetic near-wall cycle in wall-bounded turbulent flows [ 48, 98, 120]. The near-wall region encompasses the region located between 0≤ y + ≤ 100 and for canonical turbulent flows, this region is the most active portion of the turbulent boundary layer[98, 113]. Here,y + represents the distance away from the wall measured in so-called wall units or viscous units. Specifically length scales presented using ( +) have been made dimensionless using the viscous length scale given by δ v = ν/u τ . Here ν is the kinematic viscosity of the fluid and u τ is the friction velocity which is determined from the shear stress τ measured at the wall as u τ = p τ/ρ. The friction Reynolds number is given as the ratio between the characteristic length scale h (e.g., channel half-height or boundary layer thickness) and the viscous length scale, i.e. Re τ =h † /δ v =u τ h † /ν. Within the near-wall region, the largest velocity amplitudes are present and the production of turbulence is at its highest level. For canonical channel, pipe and boundary layer flows, a 3 significant portion of the total turbulence production occurs in the near-wall region and the peak turbulence production occurs very close to the wall at y + ≈ 12. The early experiments Kline et al. [58] provide detailed evidence for the presence of persistent and detectable low-speed fluid streaks present in the viscous and buffer layer (i.e. y + < 30), see figure 1.2. The low and high-speed streaks refer to regions of the flow where the local streamwise velocityishigher(high-speed)orlower(low-speed)thanthelocalmeanvelocity. UsingaReynolds decomposition the flow velocity u(x,t) is decomposed into its temporal mean component denoted by U(x) and a fluctuating velocity component denoted by u ′ (x,t) =u(x,t)−U(x). Using this decomposition, a high-speed streak has streamwise velocity fluctuation u ′ > 0 and a low-speed streak has streamwise velocity fluctuation u ′ < 0. The low and high-speed fluid streaks feature a high degree of spatial and temporal organization; these streaks are quasi-periodic in the spanwise direction having an average spanwise spacing of λ + z ≈ 100 which has been shown to be invariant with the Reynolds number [58, 98, 109]. The typical reported streamwise extent of the near-wall streaks is approximately λ + x ≈ 1000 [13, 48, 50, 59]. The near-wall region is also home to quasi- streamwise vortices (QSV’s) which plays an important role in the development and sustainment of the low and high-speed fluid streaks [ 14, 48]. Advancements in high-performance computing and high fidelity numerical simulations have advancedour understanding of the autonomous, orself-sustaining natureof thenear-wallcoherent streaks. Numerical simulations of Jiménez and Pinelli [48] isolated a single streak by constraining the spanwise dimension of their computational domain to L + z ≈ 100. In these minimal domains simulations, the turbulent statistics generated by the singular autonomous streak matched the statistics generated from larger domains. For more constricted domains with L + z < 100 a tur- bulent flow could not be sustained and the flow would eventually laminarize. Early simulations 4 x + z + Figure 1.2: Hydrogen bubble visualization of the near wall coherent streaks at y + = 4.5 by Kline et al. [58]. Photograph shows the spanwise organize of the low-speed streaks and the streamwise extend of the streaks. The streamwise and spanwise directions are denoted byx + and z + respectively. of Jiménez and Moin [51] and Hamilton, Kim, and Waleffe [ 44] observed similar relaminariza- tion effects for domains in where L + z < 100. Following these results, numerical simulations of Jiménez and Pinelli [48] investigated the regeneration mechanism behind the wall streaks. Min- imal channel simulations of Jiménez and Pinelli [48] and simulations of Jiménez and Moin [51] and Hamilton, Kim, and Waleffe [ 44] suggested that the generation of wall streaks is necessary to sustain turbulence. Additionally, the simulations of Jiménez and Pinelli [48] indicate that the near-wall flow operates autonomously from the outer flow. The quasi-streamwise vortices present in the buffer layer can be best described as ‘turbulence pumps’whichliftlow-speedfluidawayfromthewallandbringinhigh-speedfluidtowardsthewall. This ejection of low-speed fluid away from the wall and sweeping—inrush—of the high-speed fluid plays a vital role in transferring momentum between the viscous sublayer and turbulent outer layer. Although we will primarily focus on describing how these ’turbulence pumps’ transfer momentum and generate wall friction, this mechanism is vital for transporting heat and mass 5 λ + z ≈ 100 Figure 1.3: Planar snapshot at y + ≈ 20 of the near-wall streaks from the numerical simulations of Schoppa and Hussain [104]. The dark contour show the elongated low-speed streaks u ′ < 0. The dimension of the snapshot is (L + x ,L + z ) = (450,1500) and spanwise spacing of λ + z ≈ 100 can be seen here. to/away from the wall as well. A significant amount of turbulence in the buffer region is produced by the ejection and bursting of the low-speed streaks [58, 98]. Kline et al. [58] initially described the sequence of steps involving the eventual breakup of the low-speed streaks, which was referred to as bursting. This process involved the gradual lifting of the low-speed streaks, an onset of violent oscillation occurring aroundy + ≈ 8−12, and the chaotic breakup—bursting—of the low- speed streaks occurring between 10≤ y + ≤ 30. The bursting phenomenon associated with the low-speed streaks was later understood to be a product of stronger ejections emanating from one or more quasi-streamwise vortices [2, 98]. The contribution by the sweeping and ejection events can be analyzed by performing quadrant decomposition of the (−u ′ v ′ ) component of the Reynolds shear stress. Using standard notation, ejection correspond to quadrant II (Q2) and represent the movement of slow speed fluid away from the wall (u ′ < 0, v ′ > 0). Sweeps correspond to quadrant IV (Q4) corresponding to high-speed fluid moving towards the wall ( u ′ > 0, v ′ < 0). The two remaining quadrants, quadrant I (Q1, u ′ > 0,v ′ > 0)andquadrantIII(Q3,u ′ < 0,v ′ < 0)arereferredtoasoutwardandinwardmotions respectively. Usingthisdecompositionitisevidentthatbothsweepingandejectionevents—which 6 haveu ′ v ′ < 0—contribute to the positive production of Reynolds shear stress (−u ′ v ′ )> 0. In the near-wall region Q2 and Q4 events make the largest contribution to the Reynolds stresses. In the vicinity of the wall,y + ≤ 15, Q4 events (sweeping) contribute more to the mean shear stress than Q2 events. Further away from the wall the contrary is true, where Q2 events contribute more to the mean shear stress than Q4 events [32, 113, 121]. Given the information extracted from a quadrant decomposition of the Reynolds stress, an interesting question to pose is if the increase in drag and wall friction observed in turbulent flows is a consequence quasi-streamwise vortices and low and high-speed streaks. Experiments and numerical simulations have shown that the skin friction near the wall is not constant but is instead separated by regions of locally low and high skin friction. Orlandi and Jiménez [85] pursed a simplified two-dimensional model to answer this question. Given that the streamwise length of the streaks is significantly larger than their spanwise spacing the variation in streamwise direction is neglected. The streamwise velocity field behaves as a passive scalar under the effects of diffusion and advection due to the wall-normal and spanwise velocities( v,w). The spanwise vorticity is given by Ω z =∂ y u and the wall shear is determined by evaluating spanwise vorticity at the wall, i.e. τ 0 = µ∂ y u 0 =µΩ z,0 . In this simplified model the streamwise momentum equations evolve like a passive scalar under the effects of viscosity and convection by the spanwise and wall-normal velocity field. Spanwise averaging and forgoing the effects of viscosity and evaluating the result at the wall yields ∂ t ⟨Ω z ⟩ =−2⟨Ω z ∂ y v⟩ 0 . Rewriting this result as ∂ t ⟨τ 0 ⟩ =−2⟨τ 0 ∂ y v⟩| 0 shows that the average shear increase in regions where the downward velocity correlates with regions of high shear and consequently in regions where the upward velocity correlates with regions of low shear [85, 113]. Additionally, this analysis indicates that Ω z,0 (t)∝ exp − R ∂ y v 0 dt . Assuming that the wall-normal velocity remains constant over some spanwise distance the spanwise vorticity and 7 y + z + QSV − v ′ + u ′ QSV + v ′ − u ′ − v ′ + u ′ y + ∼20 λ + z ≈100 Figure 1.4: Simplified cartoon depicting the structure of the coherent near wall streaks and their relationship to the quasi-streamwise vortices (QSV’s). In this carton the low-speed (−u ′ ) and high-speed (+u ′ ) streaks are color coded red and blue respectively. Horizontal gradient indicates the transition between regions of high and low wall shear stress produced by the quasi-streamwise vortices. consequently, the wall shear exponentially decays if there is an ejection (∂ y v 0 > 0) motion and exponentially grows if there is sweep (∂ y v 0 < 0) towards the wall [85]. Orlandi and Jiménez [85] argued that this unique balance of growth and decay caused by the induced velocity of the QSV’s increases average wall friction. This simple yet informative model shows that the QSV’s and the low and high-speed streaks are responsible for the higher wall friction present in turbulent flows. This model also highlights the importance of sweeping motions, which serve to increase the local wall friction and are responsible the bulk production of turbulent stress near the wall. In figure 1.4 we present a simplified cartoon that summarizes the basic structure of the quasi-streamwise vortices and the near-wall streaks. In this depiction, the wall-normal velocity direction is shown above the high-speed (+u ′ ) and low-speed (−u ′ ) streaks. This cartoon depicts how quasi-streamwise vortices resemble turbulence pumps which serve to transport high momentum flow towards the wall and low momentum fluid away from the wall. Thus, the quasi-streamwise vortices generate regions of locally high and locally low wall friction. Understandingthebasicstructureofthenear-wallsstreakswillaidusinunderstandingthesuccess of passive drag control surfaces such as riblets and anisotropic permeable surfaces. 8 1.2 Drag Reduction via a Slip Surface Analytical models of Bechert and Bartenwerfer [7] and Luchini, Manzo, and Pozzi [67] have shed insight into the mechanism responsible for the drag-reducing properties of riblets. Briefly, riblets reduce drag by preferentially damping cross-flow velocity fluctuations and weakening the quasi-streamwise vortices, while allowing the mean streamwise flow to remain relatively unim- peded. The work of Bechert and Bartenwerfer [7] attempted to connect the drag-reducing capac- ity of riblets to their protrusion height h + p . In its original definition, the protrusion height is the distance, measured from the riblet tips, in which the turbulent flow appears to perceive a no-slip wall. The concept behind this model is that for drag-reducing riblets, the riblet height is on the order of the viscous sublayer, i.e. (O(h + ) . 10), and thus the flow within and near the riblets is dominated by viscous effects. Given these conditions, the protrusion height can be calculated by solving Laplace’s equation for u ′ in the y−z plane. In other words, the protrusion height calculated by Bechert and Bartenwerfer [7] captures the virtual origin of the streamwise mean flow, the location of the virtual wall where the mean flow appears to originate from. We will refer to the virtual origin of the mean streamwise flow as ℓ + U . Luchini, Manzo, and Pozzi [67] expanded the protrusion height concept of Bechertand Barten- werfer [7] to include the effects that riblets have on damping the turbulent cross-flow. In this ex- panded model the virtual origin of the cross-flow, ℓ + z was introduced and represents the location at which the cross-flow perceives a no-slip wall. In this revised model the drag-reducing capac- ity of riblets depends on the difference in virtual origin (protrusion height) between the mean streamwise flow and the virtual origin of cross-flow, i.e. ∆ + ℓ = ℓ + U −ℓ + z . Thus, ∆ + ℓ characterizes the ability of a given riblet geometry to suppress the turbulent cross-flows (e.g., arising from the QSVs) compared to the streamwise mean flow. Unlike the protrusion height used by Bechert 9 y + x + y + x + y + z + QSV ℓ + U ℓ + z DR∝∆ + ℓ =ℓ + U −ℓ + z U + s + h + ℓ + U ℓ + z ∆ + ℓ Figure 1.5: Cartoon depiction of the displacement of the perceived origin by the mean streamwise flow (center panel) and the perceived origin by quasi-streamwise vortices and the cross-flow (right panel). Panel on the left shows a depiction of spanwise periodic riblets with a spacing and height of s + and h + respectively. and Bartenwerfer [7], ∆ + ℓ is independent of the choice of origin used. In figure 1.5 we show an illustration depicting the offset between the virtual origin of the mean streamwise flow and the virtual origin of the quasi-streamwise vortices. Figure 1.5 shows a visual depiction of how the low streamwise resistance shifts the perceived origin closer to the riblet valley. In the same fashion, the higher spanwise resistance prohibits the cross-flow virtual origin from moving towards the riblet valley. Numerical simulations of Choi, Moin, and Kim [23, 24] and Chu and Karniadakis [25] and Goldstein, Handler, and Sirovich [37] investigated the mechanism behind the drag-reducing abil- ity of riblets. Choi, Moin, and Kim [23] and Chu and Karniadakis [25] conducted laminar flow simulations over riblets and observed an increase in drag when compared to smooth wall configu- ration. Intuitively this result makes sense since riblets increase the amount of wetted area. Choi, Moin, and Kim [24] performed direct numerical simulations over the optimal triangular riblet geometries identified in the experiments of Walsh and Lindemann [ 122]. For the drag-reducing configurations results indicated a decrease in all turbulence fluctuations in the near-wall region. Ofthe three fluctuatingcomponentsthe cross-flowandwall-normalvelocitiesexhibitedthe largest 10 suppression, decreasing in magnitude by roughly 10%. The Reynolds stress reduced by 15% near the riblet vicinity. Quadrant decomposition revealed that Q2 and Q4 events, corresponding to sweep and ejection events, were reduced near the vicinity of the riblets. Similar results were obtained by Chu and Karniadakis [25] and Goldstein, Handler, and Sirovich [37]. Experiments of Suzuki and Kasagi [111] using particle tracking velocimetry (PTV) and those of Lee and Lee [65] using particle image velocimetry (PIV) also showed a decrease in the wall-normal and cross-flow velocity fluctuation for drag-reducing riblets. Goldstein, Handler, and Sirovich [37] conducted numerical experiments on primitive riblet-like elements consisting of wires placed at a location h + away from the wall. Similar to prior results, the cross-flow fluctuations were dampened and the drag was reduced. The evidence from these simulations and experi- ments supports the conclusion that riblets reduce drag by damping the cross-flow fluctuation and weakening the effects of the quasi-streamwise vortices. Given that the quasi-streamwise vortices are responsible for the increase in drag observed in turbulent flows, Jiménez [ 49] conducted numerical studying the effects of offsetting the virtual origin of mean streamwise flow and the cross-flow ℓ + U ̸= ℓ + z . The virtual origin of the mean streamwise flow was located at the no-slip wall ℓ + U = 0 and the virtual origin of the cross-flow( ℓ + z ) was prescribed by setting the transverse velocity equal to a constantγ t times the velocity at some given height h. In this setup the cross-flow virtual origin is equal to ℓ + z =−γ t h/(1−γ t ) and the offset in virtual origin is given by ∆ + ℓ =−ℓ + z . Shifting the virtual origin to lie above the wall (∆ + ℓ > 0, ℓ + z < 0) resulted in the quasi-streamwise moving away from the wall and conversely seating the virtual origin below the wall (∆ + ℓ < 0, ℓ + z > 0) resulted in the QSVs moving towards the wall. Shifting the QSVs away from the wall weakened the induced transverse flow and resulted in a decrease in skin friction. Conversely, when the quasi-streamwise vortices are brought closer to the wall the transverse flow is amplified leading to an increase in skin friction. Similar results 11 were observed in the numerical simulations of Goldstein, Handler, and Sirovich [37] using wires placed at the riblet tips which applied a negative feedback force in response to the cross-flow fluctuations. Results varying the damping strength indicated that cross-flow damping did not produce an increase in drag but instead lead to a decrease in turbulent mixing and drag. The numerical simulations of Choi, Moin, and Kim [22] analyzed the effectiveness of damping the cross-flow using active control. In this setup, wall actuators applied either blowing-suction or spanwise velocity in response to a sensor located at y + ≈ 10. From the various configurations tested pure spanwise control resulted in a 30% decrease in turbulent drag. Similar to the case of riblets a quadrant decomposition indicated that spanwise control reduces the contributions from the positive Reynolds stress-producing (Q2, Q4) events. The action of active control of Choi, Moin, and Kim [22] and primitive riblets of Goldstein, Handler, and Sirovich [37] can be equivalently viewed as moving the virtual origin of the cross-flow into the channel and shifting the quasi-streamwise vortices away from the true wall. In terms of virtual origins, this can be viewed as maintain the streamwise virtual origin at the wall (ℓ + U = 0) and moving the cross-flow virtual origin away from the wall (ℓ + z < 0) and thus producing a positive offset between the two virtual origins, ∆ + ℓ > 0. For drag-reducing riblets, the effects on the flow are limited to the near-wall region. The experiments of Lee and Lee [65] and Suzuki and Kasagi [111] over trapezoidal riblets observed that the turbulent profiles and turbulent statistics for y + & 30 resembled those of a flat plate. Choi, Moin, and Kim [24] observed a similar effect for triangular riblets with s + = 20. Two- point correlations computed by Suzuki and Kasagi [111] showed that the size of the spanwise structures is not affected by the presence of riblets. Numerical simulations of Goldstein, Handler, and Sirovich [37] showed that the spacing of the near-wall streaks remained constant (λ + z ≈ 100) for the various riblet configurations tested. Wall normal snapshots from the numerical simulation 12 of El-Samni, Chun, and Yoon [102] show the presences of the near-wall streaks even in the drag increasing cases where the riblet spacing was s + ≈ 40. For the drag increasing case, the streaks showed high spanwise coherence and were less meandering than the drag-reducing and smooth wall cases. SimilarconclusionswerereachedbyJiménezetal.[52]andJiménezandPinelli[48], suggesting that drag reduction can be achieved by disrupting the cycle by which the streaks are generated; either by weakening and damping the quasi-streamwise vortices or disrupting the coherence of the streaks. In the cases of riblets and anisotropic permeable substrates, drag reduction is achieved by shifting the quasi-streamwise vortices away from the origin of the streamwise mean flow. This offset in virtual origin ∆ + ℓ > 0 shifts the peak in Reynolds stress (−u ′ v ′ ) away from regions of high mean shear. This decreases the production of turbulent kinetic energy via the termu ′ v ′ ∂ y U and decreases the amount of energy extracted from the mean flow and effectively weakens the quasi-streamwise vortices. Note that the shift between the virtual origins (∆ + ℓ > 0) occurs due to the higher spanwise resistance than the streamwise resistance posed by riblets. In the following section, we detail the connection between the geometry of riblets and their potential to reduce drag. 1.3 Riblets Given the economic and performance advantages associated with passive turbulence control, research endeavors have focused on investigating devices such as riblets, large eddy break up de- vices and additive polymers coatings. Among these devices, riblets have been the most successful in reducing the net drag in turbulent boundary layers both experimentally and numerically [9, 24, 25, 122, 123]. In addition to this, the effectiveness of riblets has been tested over wings and on aircrafts at various flight configurations. The drag-reducing capacity of riblets has been verified 13 at low speeds and at supersonic speeds with Mach numbers as high as 2.5 [71, 119, 124]. Ex- periments in the early ’80s and ’90s carried out by Walsh and Lindemann [122] and Walsh [123] and by Bechert and Bartenwerfer [7] and Bechert et al. [9] explored the connection between the size and shape of two-dimensional (streamwise-varying) riblets and their potential to reduce drag. Experiments indicated that a drag reduction as high as 10% can be achieved using thin blade-like riblets. Although these experiments were performed on ideal canonical turbulent flow, analysis and experiments suggest that riblets can reduce 2%–7% of turbulent drag over an aircraft in flight conditions [ 71, 119]. The experiments and applications mentioned above have focused on two-dimensional riblets, i.e. streamwise constant and spanwise periodic. The design space for two-dimensional riblets is immense and in the case of three-dimensional riblets this additional degree of freedom ex- ponentially increases the size of the design space. As a consequence the literature concerning three-dimensional riblets is limited. Notable experiments include those of Bechert, Bruse, and Hage [8] over simple three-dimensional riblets mimicking the three-dimensional geometry of a shark denticles, see figure 1.1. Experiments of Sirovich and Karlsson [108] investigated the effects of streamwise alignment of simple riblets. Numerical experiments of Peet, Sagaut, and Char- ron [86] investigated triangular riblets with a streamwise sinusoidal pattern. Other numerical experiments include those of Benschop and Breugem [11] who investigated blade-like riblets in various herringbone configurations. These experiments found no performance benefit—regarding the drag-reducing capacity—for three-dimensional riblets. In many of these cases, the measured drag reduction was inferior to a similar two-dimensional riblet. Given the immense parameter space and the success of two-dimensional riblets, many exper- iments have focused on investigating the connection between riblet size and their drag-reducing potential. The early experiments of Walsh [123] consisted on optimizing the drag reduction by 14 varying the size, height, and valley curvature of two-dimensional triangular riblets. The optimal geometries identified from these experiments consist of symmetric triangular riblets with a height to spacing ratio of h + /s + = 1 and spacing of s + = 8. Other optimal geometries identified were riblets with sharp tips and large valley curvature having a h + /s + = 0.5 and spacing of s + = 16. Experiments in the late ’90s by Bechert et al. [9] systematically reevaluated the performance of various two-dimensional riblet geometries. Using an oil channel facility with an adjustable riblet surface, shear and drag was measured over triangular/v-shaped, scalloped and blade-like riblets. A maximum drag reduction of 10% was observed for thin blade-like riblets which had a height to spacing ratio of h/s = 0.5. We summarize the collected data from the various experiments conducted in figure 1.6. Not all the data are presented here but the following trends are consistent across the various experiments. First, there exists some optimal spacing, s + opt for which the largest drag reduction is achieved. The optimal spacing lies in between s + opt ≈ 10–20 and varies by geometry. The decrease in drag is initially linear with respect to the riblet size for riblets in which the spacing is smaller than the optimum (s + <s + opt ). Thirdly, we note that in the case of triangular/v-shape riblets, a lower the ridge angle α—corresponding to a more blade-like riblet, see figure 1.6(b)—generates a larger decrease in turbulent drag. From these experiments, the highest drag reduction was 10% and was measured over thin blade-like riblets. Lastly for all the configurations which reduced turbulent drag, the initial reduction in drag is linear with respect to the riblet spacings + . This trend continues up to the optimal riblet spacing. For riblet spacing larger than this s + > s + opt , the linear drag regime degrades and eventually leads to an increase in turbulent drag. As noted earlier Luchini, Manzo, and Pozzi [67] expanded the protrusion height concept of Bechert and Bartenwerfer [7] to factor in the effect of high spanwise resistance on the cross-flow. In this expanded model the drag-reducing capacity of a given riblet geometry is measured by the 15 a) b) Figure1.6: Dragmeasurementscollectedovertriangularriblets,panela),byBechertetal.[9]. The optimal riblet spacings + opt are visible and the breakdown in the linear regime can be seen over the various geomtries tested. Panel b) shows the change in drag for a general riblet geometry plotted verses the riblet spacing, displaying the linear regime, the optimal spacing and the breakdown of this regime. Figures adapted from Bechert et al. [9]. difference between the streamwise and cross-flow protrusion height, i.e. ∆ + ℓ =ℓ + U −ℓ + z . Bechert et al. [9] developed a simple predictive model using the law of the wall and Prandtl’s universal law of wall friction for smooth pipe flows. A similar model was derived by García-Mayoral and Jiménez [34] using classical wall turbulence theory and Clauser’s friction law for turbulent boundary layers. The driving idea behind these models is that for riblets within the linear drag-reducing regime, the riblets increase the thickness of the viscous sublayer. This increase can be interpreted as shifting the intercept term preset in the law of the wall by a distance proportion to the difference in the protrusion height ∆B =µ 0 ∆ + ℓ . In other words the drag reduction can be approximated as DR = ∆c f c f0 ≈ −∆B (2c f0 ) −1/2 +(2κ) −1 = µ 0 ∆ + ℓ (2c f0 ) −1/2 +(2κ) −1 . (1.1) In equation (1.1), c f0 corresponds to the skin friction coefficient for a smooth wall flow, κ is the Von Kármán constant and µ 0 is a constant with a value between 0.6− 0.8. The constant B here is the integration constant appearing in the log law and ∆B models the increase in the thickness of the viscous sublayer. This simple model predicts the initial drag decrease is a function 16 of the protrusion height ∆ + ℓ . Comparisons between the experimental data of Bechert et al. [9] and predictions by the model presented in equation (1.1) reasonably reproduced the initial decrease in drag observed for various geometries. This single parameter model over-predicts the drag reduction for riblets with finite tip widths and those with large ridge angles. Unfortunately, this model cannot predict the optimal riblet spacing and cannot predict the drag reduction and drag increase beyond the optimal riblet spacing s + opt . a) b) c) d) Figure 1.7: Scaled drag reduction curves collected from Bechert et al. [9] shown in panel a). The change in drag is scaled by the slope of the linear drag regime m ℓ and is plotted against the riblet groove cross-section ℓ + g . Panels b)-d) display the optimal riblet spacing, height and groove cross-sectional area scale. Figures adapted from García-Mayoral and Jiménez [35]. García-Mayoral and Jiménez [35] pursued a non-exhaustive scaling analysis using the data collected from the experiments of experimental of Bechert and Bartenwerfer [7]. From this search, the best collapse was achieved when drag curves were plotted against the square root of the riblet groove cross-sectional area ℓ + g = √ A g + and normalized by the slope of the linear regime, m ℓ = ∂( DR) ∂ℓ + g 0 . Analysis of the shape parameters indicated that the groove cross-section is a better predictor of optimal riblet geometry, see figure 1.7. The scaling analysis indicated that optimal riblet shape has a grove cross-section ℓ + g ≈ 10.7± 1. Across all geometries the maximum drag 17 reduction predicted using equation (1.1) is approximately DR max ≈ 0.83·m ℓ ℓ + g . For b ack-of-the- envelop e calculationsthismodelbasedongroovecross-sectionℓ + g workswellforconventionalriblet designs. García-Mayoral and Jiménez [35] note that this model based on ℓ + g is purely empirical and may not work for non-conventional riblet designs. 1.4 Anistropic P ermeable Surfaces Anisotropic permeable substrates can be designed to function similarly to riblets having low streamwise and high spanwise resistance. In this configuration, permeable materials produce a favorable offset between the streamwise and spanwise slip length. Here, ℓ + U and ℓ + z refer to the streamwise and spanwise slip lengths and carry the same definitions as the protrusion height. The streamwise and spanwise slip length corresponds to the location of virtual origin perceived by the streamwise and cross-flow respectively. The permeable materials we will consider here are characterized by their thickness (H), porosity (ε), and their permeability (K). The porosity rep- resents the fraction of empty volume to total volume in permeable material, and the permeability matrixK represents the resistance to flow through the permeable material. The permeability was first used by Henry Darcy in 1856 to characterize the resistance of a porous material to the flow of fluid. In his famous column experiments Darcy postulated his empirical relationship, known as Darcy Law, to model the flow of fluid through sand filters. In this empirical relationship, the velocity of fluid is proportional to the gradient in the hydraulic head across the material and the material permeability [107]. Here we present a generalized version of Darcy’s empirical law in where the pressure gradient and the velocity through a permeable medium are related by the 18 materials permeability and the fluid dynamic viscosity, i.e. ∇p+µK −1 u = 0. The permeability tensor inR 3 has the following form K = 2 6 6 4 k xx k xy k xz k yx k yy k yz k zx k zy k zz 3 7 7 5 (1.2) Abderrahaman-Elena and García-Mayoral [1] extended the virtual origin concept of Luchini, Manzo, and Pozzi [67] to anisotropic porous substrates using the Darcy-Brinkman equations. In this analysis, they considered materials for which the permeability matrix took the following form: K = diag(K x ,K y ,K z ) and the off diagonal terms are zero. They showed that for flow conditions in where the turbulence length scales were much larger than the size of porous-scale structure the streamwise and spanwise slip lengths can be approximated as ℓ + U = p K + x and ℓ + z = p K + z . For cases where the height is not the limiting factor i.e. H + > ℓ + U , Abderrahaman-Elena and García-Mayoral [1] showed that drag reduction can be achieved with anisotropic porous materials with streamwise to spanwise anisotropy ratio, defined as ϕ xz = p K + x / p K + z is greater than one. For these materials the drag reduction is proportional to the difference in the streamwise and spanwise permeability, i.e. DR∝ (ℓ + U −ℓ + z ) = ( p K + x − p K + z ). Physically-motivated slip length models have provided useful insight into the drag-reducing effects of anisotropic porous surfaces. Numerical simulations of Busse and Sandham [ 17] paramet- rically studied the effect of anisotropic slip lengths ( ℓ + U ̸= ℓ + z ) and observed drag reduction over streamwise preferential materials i.e. (ℓ + U >ℓ + z ). These results were later re-evaluated by García- Mayoral, Gómez-de-Segura, and Fairhall [33] who showed that the drag reduction scaled linearly with (ℓ + U −ℓ + z ). Hahn, Je, and Choi [43] employed the empirical slip length model proposed by Beavers and Joseph [6] and numerically measured a 29% decrease in drag. The slip length models used for these simulations neglect the effects of wall-normal transpiration. Additionally, these 19 slip length models might not hold for more complex surfaces where the permeability length scale is larger than the viscous length scale or in cases where permeability matrix is not isotropic or diagonally dominant. More advanced models such as those by Lācis and Bagheri [61] and Lācis et al. [62] have provided a generalized tensor-based slip model. These models are derived using a homogenized approach, i.e. using volume-averaging and multi-scale expansion techniques. The effect of a per- meable surface is introduced via a slip term and a Darcy like resistance. Unlike the slip models used by Busse and Sandham [17] and Hahn, Je, and Choi [43] the framework developed Lācis and Bagheri [61] and Lācis et al. [62] includes a transpiration model. Other homogenized mod- els, such as the framework proposed by Zampogna and Bottaro [130] model the flow within the porous substrate using an effective permeability. The effects of fluid flow within the permeable medium both magnitude and and direction are accounted via this effective permeability. This effective permeability is computed via an Oseen approximation which introduces an inertial term in the leading order permeability equations developed by Mei and Vernescu [74], see Zampogna and Bottaro [130] for this extension. In the context of modeling porous lined channels, the flow within the permeable region can be modeled using Darcy’s Law with the effective permeability and this result can be coupled to the free channel flow using stress or velocity matching boundary conditions. Otherhomogenizedmodelsderivedusingvolume-averagingtechniquesincludethoseofWhitaker [127] and Whitaker [126]. In this framework, the governing equations are averaged over a micro- scopic volume in which the characteristic dimension of the volume filter is typically larger than the characteristic length of the porous element. The resulting equations are known as the volume- averagedNavierStokesequations(VANS)whichincludeavolumefiltertermthatrequiresaclosure problem to be solved. Whitaker [126] showed that the volume filter term can be represented by a 20 Darcy drag term in conjunction with a nonlinear Forchheimer correction term. The Forchheimer term models the effects of fluid inertia at high Reynolds number. This homogenized approach has been used to model turbulent flows over permeable substrates such as the simulations over isotropic permeable substrates by Breugem, Boersma, and Uittenbogaar [15] and Samanta et al. [101], and in simulations over anisotropic substrates by Rosti, Brandt, and Pinelli [99] and Rosti, Cortelezzi, and Quadrio [100]. Other models such as the one proposed by Brinkman [16] include a diffusive correction term (˜ ν∇ 2 u) to Darcy’s law. Here the viscosity in the permeable medium is given by (˜ ν) which differs from the fluid viscosity ( ν). This diffusive term is included to account for the effects of diffusion within the permeable substrates at large permeabilities. The Brinkman equations can be coupled to a free flow by matching normal and tangential stress at the interface. The Brinkman equations have been used in the analysis of Gómez-de-Segura, Sharma, and García-Mayoral [41] and in the numerical simulations of Gómez-de-Segura and García-Mayoral [40] to model anisotropic permeable substrates. a) b) Figure 1.8: Measured drag reduction over anisotropic permeable mediums with streamwise anisotropy ratio ϕ xy ≈ 3.6–11.1. Figures adapted from Gómez-de-Segura and García-Mayoral [40] 21 Simulations over anisotropic permeable substrates have provided insight into the mechanisms responsible for drag reduction. Simulations of Rosti, Brandt, and Pinelli [99] utilizing the VANS equations measured a 20% decrease in drag. These simulations considered substrates defined by K = diag(K x ,K y ,K z = K x ), where the streamwise and spanwise permeability were equal. A reduction in drag was observed in cases where the ratio between the streamwise to wall- normal permeability was greater than 1 i.e. ϕ xy = p K + x / p K + y > 1, not to be confused with ϕ xz = p K + x / p K + z which was set to unity for these simulations. An increase in drag was ob- served for isotropic configurations ( ϕ xy = 1) and configurations where the streamwise to wall- normal permeability ratio was below unity (ϕ xy < 1). The results of Rosti, Brandt, and Pinelli [99] appear to contradict the models developed by Abderrahaman-Elena and García-Mayoral [1], since K x =K z would imply no drag reduction. This can be explained by looking at the limiting effects of the wall-normal permeability. The analysis of Abderrahaman-Elena and García-Mayoral [1] considered cases where the wall-normal permeability was not the limiting factor and where the spanwise slip length could be approximated by ℓ + z = p K + z . In cases where the wall-normal per- meability is a limiting factor the spanwise slip length will be less than the spanwise permeability. In these cases the spanwise slip will depend on the wall-normal and spanwise permeability, i.e. ℓ + z =f( p K + y , p K + z ). Numerical simulations of Gómez-de-Segura and García-Mayoral [40] considered configurations where the permeability matrix was defined by K = diag(K x ,K y =K z ,K z ) andK x >K z . In this configuration the spanwise and streamwise slip can be approximated by (ℓ + U ,ℓ + z )=( p K + x , p K + z ). Numerical simulations showed that the initial reduction in drag over these substrates is propor- tional to ∆ + ℓ = ( p K + x − p K + z ), see figure 1.8(a). These results are in agreement with the predictions established by Abderrahaman-Elena and García-Mayoral [1]. Numerical simulations indicated that a larger decrease in turbulent drag could be achieved with materials featuring a 22 higher degree of anisotropy, i.e.ϕ xz > 1. Gómez-de-Segura and García-Mayoral [40] demonstrated that the drag reduction curves collapse when plotted against the wall-normal permeability p K + y , and normalized by the initial slope of the drag reduction curves (ϕ xz −1), see figure 1.8(b). Thus, the drag reduction curves for anisotropic permeable materials show characteristics similar to those over riblets: a linear regime where the drag reduction is proportional to the offset in the virtual origin (∆ + ℓ ) and a regime where this drag reduction degrades and eventually turns into a drag increase. This breakdown occurs around p K + y ≈ 0.4, and as we will detailed in the next sec- tion, this breakdown is the consequence of energetic spanwise rollers resembling Kelvin-Helmholtz vortices emerging over the porous substrates. 1.5 Drag Breakdo wn Mec hanism As mentioned in the prior section, the scaled drag reduction curves measured over various riblet geometries, and the scaled drag reduction over anisotropic porous substrates share similar characteristics. Additionally, they share one more common feature: the degradation in the linear drag-reducing regime can be attributed to the appearance of spanwise coherent rollers arising from Kelvin-Helmholtz (KH) like instability. Our discussion of these coherent rollers begins by describing their appearance over riblets. In the case of riblets, the breakdown in the linear drag regime can be attributed to either the lodging of turbulent structures within the riblet groove or the appearance of spanwise coherent rollers [24, 35, 65, 111]. Choi, Moin, and Kim [24] hypothesized that the drag increase observed over riblets was the consequence of turbulent eddies lodging within the grooves. Simulations over drag-increasing triangular riblets confirmed this hypothesis. In the case of drag-reducing riblet, the eddies remained above the riblet and their downwash was limited to the riblet tips. Similar results were observed in the experiments of Suzuki and Kasagi [111] and Lee and Lee [65] over 23 trapezoidal and scalloped riblets. Specifically, the lodging effect was observed over drag-increasing riblets for which the increased downwash associated with these eddies and the increased exposed wetted led to an increase in drag. In terms of groove lengthℓ + g , the lodging phenomena has been observed for riblets with ℓ + g > 17 well beyond the optimal value of ℓ + g, opt ≈ 11. In order to address this García-Mayoral and Jiménez [35] numerically studied the viscous breakdown near ℓ + g, opt . Simulations at Re τ ≈ 180 showed the appearance of spanwise rollers over thin blade-like riblets with an aspect ratio of h s = 0.5 occurring forℓ + g & 11. The premultiplied en- ergy spectra showed a visible peak centered around structures with a streamwise size ofλ + x ≈ 150 and with a spanwise size greater thanλ + z > 100. This peak was also present in the premultiplied cospectra of the Reynolds stress in the region y + . 25. Simulations at higher Reynolds numbers (Re τ ≈ 550) displayed an identical trend suggesting that this phenomenon scales in inner wall units [36]. Using a momentum balance García-Mayoral and Jiménez [35] showed that structures with 65≥ λ + x ≥ 290 and λ + z ≥ 130 were responsible for the increase in drag for l + g & 11. To explain the appearance of these structures García-Mayoral and Jiménez [35] pursued an inviscid two-dimensional linear stability analysis. Consistent with the DNS observations a stability anal- ysis showed a sharp transition towards a Kelvin-Helmholtz instability for riblets with ℓ + g ≈ 11. However, the most unstable mode determined via this analysis had a streamwise wavelength λ + x ≈ 80, approximately half the streamwise size of structures observed in numerical simulations. The scaling analysis of García-Mayoral and Jiménez [35] indicates that the onset of this in- stability is dependent onl + g , and so this phenomenon is dependent on the overall riblet geometry. Recent simulations of Rastegari and Akhavan [95] over blunt groove riblets did not observe the appearance of these rollers. Recent work by Endrikat et al. [28] indicates blunt triangular and trapezoidal riblets suppress this instability and that sharp blade-like riblet geometries promote this instability. For cases where instability is not present, the increase in drag can be traced back 24 to an increase in dispersive stresses caused by secondary motions within the riblet groove. Numer- ical simulations indicate that for riblets with large grooves and aspect ratios h s . 0.5 the increase in drag is predominately associated with an increase in dispersive stress which can account for up to 40% of the total stresses [79, 80]. Similarly, the breakdown in the lineardrag regime overpermeable substrates can be attributed to the appearance of KH-like rollers caused by a relaxation in the wall-normal permeability. Sim- ilar to riblets, the onset of this instability limits the maximum achievable drag reduction. Using tools from linear stability analysis Abderrahaman-Elena and García-Mayoral [1] and Gómez-de- Segura, Sharma, and García-Mayoral [41] predicted that this instability is governed by a relax- ation in the wall-normal permeability p K + y . The numerical simulations of Gómez-de-Segura and García-Mayoral [40] over streamwise preferential substrates indicated that the maximum drag reduction occurs at p K + y ≈ 0.4, see figure 1.8(b), and drag increases beyond this value for wall- normal permeability. For a drag increasing configuration, the premultiplied energy spectra show a region of energy accumulation centered around structures with λ + x ≈ 150. In the drag in- creasing configuration, flow visualizations do not show the signature streaky structure associated with smooth-walled flows. Instead, the dominant flow structures present are streamwise periodic spanwise constant rollers resembling Kelvin-Helmholtz vortices. Similar to the cases of riblets, a momentum balance showed that structures with 70 < λ + x < 320 and with λ + z > 120 were responsible for the increase in the Reynolds stress and drag observed for p K + y & 0.4. For completeness, we note that spanwise rollers have previously been observed in the numeri- cal simulations of Breugem, Boersma, and Uittenbogaar [15] over isotropic permeable substrates and in the simulations Jiménez and Pinelli [48] for wall-normal permeable substrates. Simulations carried out by Kuwata and Suga [60] over spanwise impermeable substrates observed the devel- opment of coherent spanwise rollers resembling Kelvin-Helmholtz vortices. These rollers produce 25 strong wall-normal velocities which prevent the development of the near-wall streaks [15]. Ad- ditionally, these rollers increase the Reynolds stress near the fluid-porous interface leading to a significant increase in drag. We note that for the simulations of Breugem, Boersma, and Uitten- bogaar [15] and Kuwata and Suga [60] the wall-normal permeability length scales were p K + y & 5, an order of magnitude higher than the value determined by Gómez-de-Segura and García-Mayoral [40]. Experiments carried out by Suga, Nakagawa, and Kaneda [110] also showed the formation of these rollers over manufactured anisotropic permeable substrates featuring an anisotropy ratio of ϕ xy ≤ 1 and p K + y & 30. For the drag increasing cases of Rosti, Brandt, and Pinelli [99] flow visualizations did not show evidence of vortical structures resembling a Kelvin-Helmholtz vortices. Autocorrelation function for the wall-normal velocity showed a weak signature of said rollers near the fluid-porous interface. The presence of these weak Kelvin-Helmholtz vortices can be attributed to low porosity and low permeability. 1.6 Con tributions and Outline In the above sections, we have provided insight into the mechanism responsible for drag reduction and drag increase over riblets and streamwise preferential porous materials. Drag reduction arises due to how riblets and anisotropic permeable substrates interact with the near- wall quasi-streamwise vortices. The appearance of spanwise constant rollers resembling Kelvin- Helmholtz vortices limits the maximum achievable drag reduction. For riblets the appearance of these rollers is dependent on riblet geometry. For cases where this instability is not present, the increase in drag can be traced back to an increase in turbulent activity and secondary motions within the riblet grooves. In essence, predicting the drag performance of these substrates reduces down capturing two effects; predicting the effects of these surfaces on the near-wall streaks and predicting the appearance and impact of Kelvin-Helmholtz like rollers. 26 Unfortunately, the large parameter space involved with designing patterned surfaces makes it prohibitivelyexpensivetocharacterizethedragpropertiesviaexperimentsorthroughhighfidelity simulations. Testing such surfaces at application relevant Reynolds numbers imposes additional hardships both of experiment and numerical simulations. Over the past decades, several riblet geometries have been proposed and tested ranging from scalloped to blade-like riblets. Experi- ments have also tested various streamwise arrangements mimicking the three-dimensional pattern observed of shark denticles. Similarly, simulations have shown the potential for streamwise prefer- ential permeable substrates to reduce turbulent drag. Many of these experiments and simulations have focused on simple patterned surfaces and on porous substrates characterized by diagonally dominant permeability tensors. Unfortunately, the connection between the geometric properties of these surfaces and their interaction with the near-wall turbulence is not fully understood. It remains unclear if these surfaces can optimized to mitigate the detrimental effects of the drag increasing Kelvin-Helmholtz like rollers. To address these questions we will present a reduced-order model founded on the resolvent framework of McKeon and Sharma [73] to investigate and characterized the effects of patterned and porous surfaces on turbulent flows. These models are aimed at providing a computationally efficient,inexpensive,andfullypredictiveframework. InChapter 2wewillgiveabriefintroduction totheresolventframeworkforsmoothwall-boundedturbulentflows. InChapter 3wewilldescribe the modification to the resolvent framework to account for spanwise periodic and streamwise constant riblets. In Chapter 4 we combine this reduced-order model with a nonlinear gradient- based optimizer to optimize the shape and size for a two-dimensional riblet configuration. In Chapter 5 we present an extension to the resolvent framework for anisotropic porous substrates. Finally, in Chapter 6 we present our concluding remarks and possible future directions for these reduced-order models. 27 Chapter 2 Resolv en t F ramew ork This chapter provides a brief review of the resolvent formulation proposed by McKeon and Sharma [73] for smooth-walled turbulent channel flow. Under the resolvent formulation, the Fourier transformed Navier-Stokes equations (NSE) are interpreted as a forcing-response system in which the nonlinear convective terms are treated as internal forcing (input) to the linear system which generates a turbulent velocity and pressure response (output). This framework borrows tools from the successful non-modal analysis techniques of Trefethen et al. [118] which analyze the resolvent operator and the pseudospectra of the linearized Navier-Stokes equations. We note that the linearized Navier Stokes equations have been a powerful tool in studying behaviors such as the transition to turbulence (see Butler and Farrell [18], Farrell and Ioannou [30], Henningson [45], and Reddy et al. [96] and Jovanović and Bamieh [56]), studying the amplifications of stochastic disturbances (see Bamieh and Dahleh [5] and Jovanović [53] and Jovanović and Bamieh [56]) and reproducing structural and statistical features of wall-bounded turbulent flows (see McKeon and Sharma [73] and Hwang and Cossu [46]). The input–output approach of Jovanović [53] and Jovanović and Bamieh [56] has been a powerful tool for analyzing the resonant nature of the Tollmien–Schlichting and their importance in the transition to turbulence and in understanding the resonant behavior of the near-wall streamwise streaks. Additionally, these tools have been 28 useful in quantifying the effects of control on turbulence (see Jovanović [ 55], Luhar, Sharma, and McKeon [68, 70], and Moarref and Jovanović [76, 77] and Nakashima, Fukagata, and Luhar [81]) and in predicting second-order statistics (see Zare, Georgiou, and Jovanović [131] and Zare, Jovanović, and Georgiou [132] and Towne, Lozano-Durán, and Yang [117]). For a summary and review of the advancement of tools and models based on this input-output formalism the reader is direct to Jovanović [54]. The key difference between the resolvent framework of McKeon and Sharma [ 73] and the linearized approach of Butler and Farrell [18], Jovanović [53], and Trefethen et al. [118] and others, albeit mathematically these approaches are very similar, is in the manner in which the nonlinear terms are treated. In the resolvent framework McKeon and Sharma [73], the nonlinear term are retained and treated as endogenous forcing to the transfer function formed by the linear components of the Navier-Stokes equations. To construct the system resolvent we begin writing the nondimensionalized governing equations for incompressible flow: ∂u ∂t =−∇p−(u·∇)u+ 1 Re τ ∇ 2 u, (2.1a) 0 =∇·u. (2.1b) Using a Reynolds decomposition the velocity and pressure can be decomposed into time averaged component and into a time fluctuating component about this mean. Here ( ) denotes the time- averaging operator, and we define the time-averaged (mean) turbulent velocity profile as U(x) := u(x,t) and the time-averaged (mean) pressure as P(x) := p(x,t). Up to this point the only assumption we make is that for the flow which we are analyzing there exists a well-defined mean velocity U(x) and mean pressure P(x) profile. The fluctuating velocity components can be computed as u ′ (x,t) =u(x,t)−U(x) and equivalently the fluctuating pressure is computed as 29 p ′ (x,t) =p(x,t)−P(x). Substituting this decomposition into the incompressible Navier-Stokes equations and subtracting the mean components yields the governing equations for the fluctuating components: ∂u ′ ∂t =−∇p ′ +Lu ′ +f, (2.2a) 0 =∇·u ′ . (2.2b) Here f =− u ′ ·∇u ′ −u ′ ·∇u ′ represents the nonlinear fluctuating terms and L represents the linearized terms computed about the mean turbulent velocity profile. For the remainder of this thesis we will consider a turbulent channel flow configuration that is stationary in time ( t) and homogeneous in the streamwise (x) and spanwise directions (z). For this configuration the mean turbulent velocity profile can be expressed as U(x) = [U(y),0,0]. The velocity and pressure fluctuations can be Fourier transformed and expressed as follows: " u ′ (t,x) p ′ (t,x) # = y " u k (y) p k (y) # exp(−iωt+iκ x x+iκ z z)dωdκ x dκ z . (2.3) We can compactly express the Fourier-transformed Navier-Stokes equation (2.2) as follows: " u k p k # = −iω " I 0 # − " L k − ˜ ∇ − ˜ ∇ T 0 #! −1 " I 0 # f k =H k f k . (2.4) and in where L k = 2 6 6 6 6 4 −iκ x U(y)+Re −1 τ ˜ ∇ 2 − dU/dy 0 0 −iκ x U +Re −1 τ ˜ ∇ 2 0 0 0 −iκ x U(y)+Re −1 τ ˜ ∇ 2 3 7 7 7 7 5 (2.5) 30 In the expression above, ˜ ∇ T and ˜ ∇ represent the divergence and gradient operators, respec- tively. Thedifferentialoperator ˜ ∇isdefinedas ˜ ∇ = (iκ x ,∂/∂y,iκ z ), andtheFouriertransformed Laplacian operator ˜ ∇ 2 =−κ 2 x −κ 2 z +∂ 2 ∂y 2 . The first row in ( 2.4) corresponds to the momentum equations and the second row represents the continuity constraint. A subscript k = (κ x ,κ z ,ω) denotes a specific wavenumber-frequency combination. Note that f k = −(u ′ ·∇u ′ ) k is the Fourier-transformed nonlinear term. This nonlinear term is interpreted as an internal forcing that is mapped to a velocity (u k ) and pressure (p k ) response by the resolvent operator, H k . Bear in mind that the resolvent depends onL k , which represents the linear terms in the NS equations after Reynolds-decomposition. As a result, construction of H k requires knowledge of the turbulent mean profile, U(y). At each wavenumber-frequency combination,k, the resolvent operator maps the forcing input to a velocity-pressure response. These response modes can be seen as traveling waves with a streamwise and spanwise wavelength of λ x = 2π/k x and λ z = 2π/k z traveling at a speed of c = ω/k x . A major contribution of the resolvent framework lies in the finding that the forcing- response transfer function tends to be low-rank at wavenumber-frequency combinations that are energetic in wall-bounded turbulent flows [ 73, 78]. As a result, for many k combinations, the discretized resolvent operator can be approximated well via a rank-1 truncation given by its singular value decomposition (SVD). To ensure grid independence and enforce anL 2 energy norm under the SVD, the resolvent operator is scaled such that h W u 0 i " u k p k # = h W u 0 i H k W f −1 W f f k , (2.6) or W u u k =H w k W f f k . (2.7) 31 The block diagonal matricesW u andW f contain numerical quadrature weights which ensure that the SVD (and rank-1 truncation) of the scaled resolvent operator H w k = X m ψ k,m σ k,m ϕ ∗ k,m ≈ψ k,1 σ k,1 ϕ ∗ k,1 (2.8) where σ k,1 >σ k,2 ...>σ k,m ...> 0, ϕ ∗ k,l ϕ k,m =δ lm , ψ ∗ k,l ψ k,m =δ lm , (2.9) yields forcing modesf k,m =W f −1 ϕ k,m and velocity response modesu k,m =W u −1 ψ k,m with unit energy over the channel cross-section. In other words, the scaling ensures that the orthonormality constraints in (2.9) yield Z 2 0 f ∗ k,l f k,m dy =δ lm , Z 2 0 u ∗ k,l u k,m dy =δ lm , (2.10) where a superscript (∗) denotes a complex conjugate. The wall-normaly coordinate is normalized with the channel half-height, h † ; wherein y = 0 denotes the lower wall and y = 2 represents the upper wall of the channel. Equations (2.4)-(2.9) shows that forcing in the direction of them th singular forcing mode with unit amplitude results in a response in the direction of the m th singular response mode which is amplifiedbythesingularvalue σ k,m . Thus, forcingf k =f k,1 createsaresponseu k =σ k,1 u k,1 . As noted earlier, the resolvent operator is known to be low-rank atk combinations that are energetic in natural turbulence, often with σ k,1 is significantly larger than σ k,2 . High amplification of the rank-1 response modes implies that the flow field may be reasonably approximated using u k,1 at each k. Recent studies show that this rank-1 approximation is able to reproduce key structural 32 and statistical features of smooth-walled turbulent flows [ 72, 78]. This rank-1 truncation has also yielded useful control-oriented low-complexity models to design and asses the effects of opposition control on turbulence [68–70, 81]. In particular, it has been demonstrated that the amplification of specific rank-1 response modes that resemble dynamically-important features of wall-bounded turbulent flows (e.g., the near-wallcycle)isausefulpredictorofperformance. Motivatedbythesepreviousmodelingefforts, we retain the rank-1 approximation for the remainder of this thesis. We assume that the flow field at each k corresponds to the rank-1 responseu k,1 . The term r esolvent mo des refers to these rank-1 flow fields. The terms amplification and gain are used interchangeably to refer to the corresponding singular values, σ k,1 . In addition to the work on wall-bounded turbulence discussed above, resolvent analysis has also been used to develop low-order models for other flow configurations [e.g., 39, 103, 112, 129]. For more information on the rank-1 approximation and its validity, the readers are referred to Beneddine et al. [10]. Further details on numerical implementation for a smooth channel can be foundinMcKeonandSharma[73]andLuhar, Sharma, andMcKeon[69]. Inthefollowingchapter, we will describe the modifications to the resolvent formalism to account for two-dimensional spanwise periodic riblets. 33 Chapter 3 Resolv en t Analysis for T urbulen t Channel Flo w With Riblets 3.1 Ov erview This chapter has been published in the American Institute of Aeronautics and Astronautics Journal, see Chavarin and Luhar [21]. Over the past decades, experiments and simulations have providedinsightintothemechanismresponsibleforthedrag-reducinganddragincreasingcapacity of two-dimensional spanwise periodic riblets. Unfortunately, the large parameter space involved with the design of 2D riblets makes it prohibitively expensive to quantify the change in drag via experiments or through high fidelity simulations. The limitations on the manufacturing of small riblets for experiments and the computational demands for these numerical simulations prohibits the parametric study of riblet geometries at higher Reynolds numbers. The experiments of Walsh and Lindemann [122] and Bechert et al. [9] have cataloged the performance over a limited set of two-dimensional riblets which serve as general guidelines for designing such riblets. Given these limitations, we provide an extension to the resolvent framework for smooth wall turbulence proposed by McKeon and Sharma [73] to account for the effect of streamwise-constant riblets. High-gain resolvent modes have proven to be useful models for dynamically-important 34 flow features such as the near-wall cycle [ 73, 106]. These modes can serve as the building- blocks for reduced-order models, in which the effect of control techniques can be evaluated on a limited set of high-gain resolvent modes instead of the full turbulent flow field. In recent years this framework has been successfully used to develop computationally inexpensive models for designing and assessing active turbulence control techniques [68–70, 81]. In the present extension, the effects of riblets are introduced into the Navier-Stokes equations using a linear spatially-varying body force. The linearized equations derived using this modeling approach require special attention. We focus on riblets which are spanwise periodic and rela- tively narrow-banded. In other words the NS equations are linearized around a spanwise periodic solution which gives rise to a linear system with periodic coefficients. Linear systems that fall within this category have been the subject of investigation by Fardad, Jovanović, and Bamieh [29] and Moarref and Jovanović [76]. Fardad, Jovanović, and Bamieh [29] investigated how the periodic feedback introduced by these systems can alter the stability characteristics of linear time- invariant systems. Moarref and Jovanović [76] investigated the spatially periodic linear system derived via the linearization of the Navier-Stokes equations describing the effects of upstream and downstream traveling waves. These waves were actuated via blowing and suction at the wall with a tunable streamwise wavelength, wave velocity, and amplitude. Model predictions using the spatially periodic linear system indicated that downstream traveling waves can lead to a reduction in turbulent kinetic energy and can prevent the transition to turbulence. Downstream traveling wise waves designed using the models of Moarref and Jovanović [76] were implemented in the numerical simulations of Lieu, Moarref, and Jovanović [66] which showed that the waves with the parameters selected using this model were effective in preventing the onset of turbulence. In this chapter, we will investigate the spatially distributed linear system which describes the turbulent interaction with two-dimensional spanwise periodic riblets. We will show that 35 depending on the riblet configuration—which prescribes the coefficient of our spatially periodic linear system—can lead to a suppression of the modes which act as a surrogate model for the near-wall streak. Depending on the riblet configuration this can lead to the destabilization and amplification of spanwise coherent mode which resembles Kelvin-Helmholtz vortices. We will present predictions for spanwise-periodic and streamwise-constant riblets and show that specific high-gain modes identified from the modified governing equations reproduce observations made in prior direct numerical simulations with limited computation. In Section 3.2 we describe the volume penalization technique used to account for the pres- ences of riblets, the gain-based decomposition, as well as the numerical implementation. Model predictions are presented in Section 3.3. In particular, we consider the effect of rectangular riblets of varying size on the energetic near-wall cycle (Section 3.3.1) and we test for the emergence of spanwise-constant rollers resembling Kelvin-Helmholtz vortices (Section 3.3.2). We also compare these predictions against DNS results [34, 35]. In Section 3.3.3, we pursue limited brute optimiza- tion of triangular and trapezoidal riblets, and compare our predictions against prior results and experiments. Concluding remarks are presented in Section 3.4. 3.2 Mo deling Approac h 3.2.1 V olume P enalization for Riblets To account for the effect of patterned walls of various geometries and shapes we employ a volume penalization technique: the fictitious domain methodology of Khadra et al. [57]. Solid obstructions within the fluid domain — riblets in our case — are modeled as a spatially varying permeability function K(x,y,z), which appears as an additional body force in the momentum 36 equations. This approach allows consolidation of the coupled fluid-solid system into a single heterogeneous numerical domain. The equations governing fluid flow in this domain are given by: ∂u ∂t +u·∇u =−∇p+ 1 Re τ ∇ 2 u−K −1 u, ∇·u = 0. (3.1) The volume penalization term is modeled using a Darcy-type linear resistance, with dimensionless permeability K. The friction Reynolds number is defined as Re τ = u τ h † /ν, in which u τ is the friction velocity and ν is kinematic velocity. For our problem, the permeability function takes on two values. In the fluid region that is free from all solid obstacles,K→∞. Here, the governing equations reduce down to the standard NSE. In the region of the domain that is occupied by a solid obstruction, K → 0. Here, the permeability term becomes the lead order term in Eq. (3.1), forcing the velocity in the solid domain towards zero. Angot [3] provided a rigorous analysis of the error bounds for this method. The H 1 norm of the error in velocity computed by the penalization method, when compared to the true velocity field, is expected to be O(K 1/4 ) over the entire domain. The L 2 norm of error in velocity is expected to beO(K 3/4 ) in the solid domain andO(K 1/4 ) in the fluid domain Angot [3] and Khadra et al. [57]. SinceK cannot be driven to zero numerically, a value ofK∼O(10 −8 ) guarantees an L 2 error ofO(10 −2 ) in the fluid region of the domain, and of O(10 −5 ) in the solid region. For completeness, we note that the volume penalization approach is conceptually similar to more sophisticated immersed boundary techniques used in fluid flow simulations, albeit without a feedback component to the body force that explicitly enforces the boundary conditions [37, 38, 87]. 37 This volume penalization method allows us to avoid the additional work required to form meshes that conform to the geometry of the riblets. Instead, we utilize a simpler Cartesian mesh corresponding to the standard channel flow geometry into which the riblets are embedded via the modified governing equations. From the perspective of resolvent analysis, the key benefit of this method is that the modified equations retain the linearity of the forcing-response transfer function used for the gain-based decomposition. In other words, the resolvent analysis approach described in the previous section can be extended to account for the effects of riblets as long as the permeability function can be expressed analytically and its spatial Fourier transform can be computed. 3.2.2 Mo dified Resolv en t Op erator For the remainder of this chapter we focus on streamwise-constant riblets such that the per- meability function only varies in the wall-normal and spanwise directions, K(y,z). Moreover, assuming that the riblets are periodic in the spanwise direction with spacing s, the permeability function can be expressed analytically as a Fourier series of the form K −1 (y,z) = ∞ X n=−∞ a n (y)exp{inκ s z} (3.2) with κ s = 2π/s. The Fourier coefficients a n (y) account for any wall-normal variation in riblet geometry. Physically, the function K −1 is a periodic rectangular wave (or pulse tr ain) in the spanwise direction for riblets. K −1 is zero in the unobstructed fluid domain corresponding to the grooves and K −1 →∞ in the solid domain. The spacing between these rectangular pulses in K −1 is always s, which is why the Fourier series representation in equation (3.2) only includes harmonics of κ s . However, the width of these pulses, which represents the width of the solid domain, can vary in the wall-normal direction. As an example, for rectangular riblets, the pulse 38 width is constant iny over the height of the riblets. For triangular or trapezoidal riblets, the pulse width decreases linearly over the height of the riblets. This pulse width determines the Fourier coefficients a n (y). Since the mean velocity profile will exhibit the same periodicity as the riblets, this profile can also be expressed as a Fourier series of the form: U(y,z) = ∞ X n=−∞ U n (y)exp{inκ s z}. (3.3) To construct the resolvent operffiator, we follow a similar set of steps to those outlined in Section 3.2. The full turbulent velocity field in the channel is decomposed into a time-averaged mean component that varies in the y and z directions, and fluctuations about this mean u(x,t) =U(y,z)+u ′ (x,t). (3.4) This yields the following equation for the non-zero streamwise component of the mean flow − ∂P ∂x + 1 Re τ ∇ 2 U−K −1 U− ∇·(u ′ u ′ ) x = 0. (3.5) Note that the wall-normal and spanwise components of the mean flow are assumed to be zero, i.e., U = (U,0,0). The pressure gradient driving the mean flow is −∂P/∂x = 1 for the inner- normalized equations presented here, and the subscript x for the Reynolds stress term denotes the streamwise component. The fluctuations satisfy ∂u ′ ∂t +U·∇u ′ +u ′ ·∇U +∇·(u ′ u ′ −u ′ u ′ ) =−∇p ′ + 1 Re τ ∇ 2 u ′ −K −1 u ′ , (3.6) 39 and the continuity constraint ∇·u ′ = 0. (3.7) To estimate the mean profile, which is required for the resolvent operator, we model the Reynolds stressterm in equation (3.5) using an eddy viscosity,ν e . Combiningequation (3.2), equation (3.3), and equation (3.5) with the eddy viscosity model, the governing equation for then th Fourier mode for the mean flow can be expressed as 1 Re τ ν e d 2 dy 2 + dν e dy d dy δ(n)+ d 2 dy 2 −n 2 κ 2 s U n − X p+q=n U p a q = ∂P ∂x δ(n), (3.8) in which δ(n) = 1 for n = 0 and δ(n) = 0 for n̸= 0. As shown in equation (3.24), the eddy viscosity is assumed to act only on the spanwise-constant component of the mean flow ( U 0 ). The mean pressure gradient driving the channel flow is also assumed to be spanwise-constant. Further details on the eddy viscosity formulation and the numerical procedure used to compute the coupled set of equations for each mean flow harmonic are provided in Section 3.2.3. The fluctuating component of the flow field is Fourier-transformed in the homogeneous spatial directions (x,z) and time (t). Under the Fourier transform, the governing equation for a specific wavenumber-frequency combinationk can be expressed as −iω +iκ x U 0 +a 0 − 1 Re τ ˜ ∇ 2 u k + ˜ ∇p k +v k ∂U 0 ∂y e x + ∞ X n=−∞, n̸=0 (a n +iκ x U n )u k−nκs + ∞ X n=−∞, n̸=0 v k−nκs ∂ ∂y +inκ s w k−nκs U n e x =f k , (3.9) in which the subscriptk−nκ s represents the wavenumber frequency combination(κ x ,κ z −nκ s ,ω). The top row in equation (3.9) represents linear interactions between the Fourier mode of interest 40 and the 0 th spatial harmonic of the mean flow and the penalization term (i.e., the spanwise- constantcomponents). Thesecondrowincludesinteractionswiththeremainingspatialharmonics of the mean flow and penalization term. Equation ( 3.9) can be organized into the following block operator form, grouping together the terms involving interactions with the spanwise-constant and spanwise-varying components of the mean flow and penalization term −iω " I 0 # − " L k −a 0 I − ˜ ∇ ˜ ∇ T 0 #!" u k p k # + ∞ X n=−∞, n̸=0 F n " u k−nκs p k−nκs # = " I 0 # f k , (3.10) in which the linear operatorL k is identical to that in equation (2.4) L k = 2 6 6 4 −iκ x U 0 +Re −1 τ ˜ ∇ 2 −dU 0 /dy 0 0 −iκ x U 0 +Re −1 τ ˜ ∇ 2 0 0 0 −ik x U 0 +Re −1 τ ˜ ∇ 2 3 7 7 5 , (3.11) and the block operatorF n is defined as F n = 2 6 6 6 6 4 a n +iκ x U n dU n /dy inκ s U n 0 0 a n +ik x U n 0 0 0 0 a n +ik x U n 0 0 0 0 0 3 7 7 7 7 5 . (3.12) Equation (3.10) demonstrates that the spanwise variation in the governing equations and mean velocity profile introduced by the riblets couples together the wavenumber-frequency combination ofinterest,k, andaninfinitesetofharmonicswithwavenumber-frequencycombinations k−nκ s = (κ x ,κ z −nκ s ,ω) in whichn ={±1,±2,...±∞}. If the wall is homogeneous, containing no spatial variationinthespanwisedirection, thecouplingtermsdropoutandresolventanalysiscanproceed on a mode-by-mode basis, as shown in equation (2.4) for the smooth-wall case. However, for the case with riblets, the coupled system must be analyzed together. 41 To proceed, we truncate the coupled input-output system to only consider±n h harmonics of κ s , so that n ={±1,±2,...±n h }. In other words, the infinite Fourier series for K −1 and U are truncated to only consider±n h harmonics. This effectively means that the riblet geometries included in the model are an approximation to the true geometry. The numerical effect of this truncation is discussed in the next section. The governing equations for this truncated system of interacting modes can be expressed as 2 6 6 6 6 6 6 6 6 6 6 4 u s k+n h κs u s k u s k−n h κs 3 7 7 7 7 7 7 7 7 7 7 5 = 0 B B B B B B B B B B @ −iω 2 6 6 6 6 6 6 6 6 6 6 4 ˜ I ˜ I ˜ I 3 7 7 7 7 7 7 7 7 7 7 5 + 2 6 6 6 6 6 6 6 6 6 6 4 ˜ L k+nκs F 1 F 2n h F −1 ˜ L k F 1 F −2n h F −1 ˜ L k−nκs 3 7 7 7 7 7 7 7 7 7 7 5 1 C C C C C C C C C C A −12 6 6 6 6 6 6 6 6 6 6 4 f s k+n h κs f s k f s k−n h κs 3 7 7 7 7 7 7 7 7 7 7 5 . (3.13) Here, the superscript s is used to denote the complete system of governing equations, i.e., the momentum equations as well as the continuity constraint. In other words,u s k = [u k ,p k ] T includes both the velocity and pressure fields, while the forcing term f s k = [f k ,0] T explicitly accounts for the continuity constraint. The operators ˜ I and ˜ L k are defined as ˜ I = " I 0 # and ˜ L k =− " L k −a 0 I − ˜ ∇ ˜ ∇ T 0 # , (3.14) respectively. Equation (3.13) can be expressed more compactly as u s k+ =H k+ f s k+ , (3.15) with u s k+ = h u s k+n h κs , ··· u s k , ··· u s k−n h κs i T (3.16a) 42 and f s k+ = h f s k+n h κs , ··· f s k , ··· f s k−n h κs i T (3.16b) Thenotationk+representsthecoupledsystemofFouriermodesassociatedwiththeindividual wavenumber-frequency combination k. Similar to the steps outlined in (equations (2.6)–(2.10) for the smooth-wall case, the extended resolvent operator in equation (3.15) is weighted using the extended quadrature matrices, ˆ W u and ˆ W f , ˆ W u u s k+ = ˆ W u H k+ ˆ W −1 f ˆ W f f s k+ , (3.17) prior to carrying out the SVD and rank-1 truncation ˆ W u H k+ ˆ W −1 f = X m ψ k+,m σ k+,m ϕ ∗ k+,m ≈ψ k+,1 σ k+,1 ϕ ∗ k+,m . (3.18) This weighting step ensures that the forcing modes (right singular vectors, ϕ k+,m ) and response modes (left singular vectors, ψ k+,m ) identified from the SVD have unit energy when integrated over the entire channel, i.e., the unobstructed flow as well as the region where the riblets are located. The highest-gain flow structure, which comprises a coupled set of traveling wave Fourier harmonics, is computed from the rank-1 response mode via the inverse-weighting operation u s k+,1 = ˆ W −1 u ψ k+,1 . (3.19) The singular valueσ k+,1 is used as a measure of energy amplification for this high-gain flow struc- ture. In Section 3.3, we show that changes in the gain and structure of specific resolvent modes 43 provide substantial insight into riblet performance. Further details on numerical implementation (i.e., discretization iny direction, grid resolution, number of coupled Fourier harmonics retained) are provided in the following section. 3.2.3 Numerical Implemen tation The mean flow equations and resolvent operator are discretized in the wall-normal y direction using the split-domain approach described in Aurentz and Trefethen [4]. The upper domain corresponds to the unobstructed flow and is bounded between the riblet tips ( y = 0) and the channel centerline (y = 1). The bottom domain corresponds to the region in which the riblets are present and the governing equations include the penalization term; this domain is bounded by the channel flow ( y =−h) and the riblet tips (y = 0). This split-domain representation ensures that the riblet height,h, remains fixed even as grid resolution changes. In other words, it allows us to impose a hard wall-normal cut-off for the penalization term. Across the domains, we impose the following matching conditions for each spatial harmonic of the mean velocity profile U n | y=0 +− U n | y=0 − = 0 and dU n dy y=0 + − dU n dy y=0 − = 0. (3.20) Similarly, the following matching conditions are used for each harmonic present in the extended resolvent operator (equation (3.13)) u k | y=0 +−u k | y=0 − = 0, (3.21) p k | y=0 +− p k | y=0 − = 0, (3.22) du k dy y=0 + − du k dy y=0 − ! e x + dw k dy y=0 + − dw k dy y=0 − ! e z = 0. (3.23) 44 At the lower wall of the channel, y =−h, the standard no-slip boundary conditions are applied for the mean flow as well as the Fourier-transformed fluctuations. As noted earlier, an eddy viscosity model is used to compute the mean velocity profile needed for the resolvent operator. To solve for the coupled system of mean flow harmonics governed by equation (3.8), we assume that the eddy viscosity is zero in the lower domain occupied by the riblets. In other words, flow within the riblet grooves is assumed to remain laminar and so the velocity profile can be computed using the penalized version of the Stokes equations. In the upper unobstructed domain, we employ the eddy viscosity model proposed by Reynolds and Tiederman [97]. Thus, we have ν e = 8 < : 1 2 h 1+ Reτκ 3 (2y−y 2 )(3−4y+2y 2 ) 1−exp (|y−1|−1) Reτ α 2 i 1/2 − 1 2 y≥ 0 0 y< 0 . (3.24) The eddy viscosity shown above is normalized by the kinematic viscosity, ν. The parameter α appears in van Driest’s mixing length model and κ is the von Karman constant. We recognize that the eddy viscosity model shown in equation (3.24) is not rigorously justified. The ν e profile is likely to change due to the presence of the riblets, particularly in the near-wall region. Moreover, for larger riblets, the flow may not remain laminar within the riblet grooves. To account for this effect, the eddy viscosity profile would have to include a geometry-dependent origin below y = 0. Unfortunately, there is little prior literature on eddy viscosity models for turbulent flows over riblets. Therefore, accounting for these effects would introduce further assumptions and tuning parameters (e.g., virtual origin for ν e profile) that are not rigorously justified. For the results shown below, a total of N = 81 Chebyshev nodes is used for discretization in thewall-normaldirection; 2N/3 = 54nodesareusedintheupperunobstructeddomain(y∈ [0,1]) and N/3 = 27 nodes are used in the lower domain occupied by the riblets (y∈ [−h,0]). The 45 Figure 3.1: Mean velocity predictions. (a) Predicted slip velocity as a function of riblet spacing. Open symbols show predictions for rectangular and trapezoidal riblets. Closed symbols show DNS resultsfromGarcía-MayoralandJiménez[35]forrectangularriblets. (b)Meanstreamwisevelocity profile for the case identified in (a) in red. Solid contours show velocity values in increments of 2u τ . number of Fourier harmonics is truncated to n =±n h =±6. For the mean flow, this choice of N and n h led to convergence withinO(10 −2 ) in the predicted slip velocity at the riblet tips, the centerline velocity, and the bulk-averaged velocity compared to predictions made with N = 120 andn h = 10. Singular values for the resolvent modes considered in this chapter also converged to withinO(10 −2 ). Convergence toO(10 −2 ) was deemed sufficient for the purposes of this study for two reasons. First, the L 2 norm of the error associated with the volume penalization technique is expected to be ofO(10 −2 ). Second, compared to the smooth-wall case, riblets led to changes in singular values that were significantly greater than O(10 −2 ). Note that the extended resolvent operator in equation (3.13) is of size 4N(2n h + 1)× 4N(2n h + 1). As a result, increases in N and n h carry a substantial computational penalty. For instance, the computation time required to compute the resolvent operator and SVD for a single set of coupled Fourier modes increases from≈ 60 seconds (on a single core of a desktop workstation) for N = 81 and n h = 6 to≈ 600 seconds for N = 100 and n h = 10. 46 Sample predictions for the mean velocity profile are shown in figure 3.1 for rectangular and trapezoidal riblets of varying size at Re τ = 180. For the rectangular riblets, the riblet height is fixed at h/s = 0.5 and width is fixed at t/s = 0.25, in which s is the spacing of the riblets. This geometry is identical to that considered by García-Mayoral and Jiménez [34] in DNS. For the trapezoidal riblets, the ridge angle is set at α = 30 ◦ . Despite the substantial simplifications (penalization, truncation to±n h harmonics, eddy viscosity model) the predicted slip velocity at the tips of the rectangular riblets is in good agreement with DNS results from García-Mayoral and Jiménez [35]; see figure 3.1(a). Over all the configurations tested, the error in our predictions was less than 8% compared to the DNS results from García-Mayoral and Jiménez [35]. The relative error in predicted slip velocity increased with increasing riblet size. This increase in error could be due to the breakdown of the laminar flow assumption in the region occupied by the riblets. As riblet size increases, the effect of turbulence penetration into the grooves must be considered [ 35]. Further, the predicted slip velocity for trapezoidal riblets is larger than that for rectangular riblets of similar size, which is consistent with physical intuition. The sample mean profile predictions shown in figure 3.1(b) confirm that volume penalization is able to successfully enforce the presence of riblets. 3.3 Mo del Predictions and Discussion The literature review presented in the introduction suggests that, as a starting point, riblet performance may be assessed by (i) considering the effect of riblets on the energetic near-wall (NW) cycle, and (ii) by testing for the emergence of spanwise-constant rollers resembling Kelvin- Helmholtz vortices. With this in mind, we employ the extended resolvent framework to consider the effect of riblets on resolvent modes resembling the NW cycle and modes that are spanwise- constant. To enable direct comparison with the DNS results of García-Mayoral and Jiménez 47 [35], we focus on geometrically-similar rectangular riblets with height h/s = 0.5 and thickness t/s = 0.25 at friction Reynolds number Re τ = 180. García-Mayoral and Jiménez [35] considered riblets with inner-normalized spacing between s + = 8 and s + = 33. Here, we consider riblet spacing ranging from s + = 0.5 to s + = 40. Throughout this section, we use the term gain and amplification interchangeably to refer to the first singular value: σ k+,1 in equation (3.18). Predictions for flow structure in physical space are obtained via an inverse Fourier transform of the highest-gain response modes: u s k+,1 in equation (3.19). Note that u s k+,1 represents a coupled set of Fourier harmonics. From here on, we drop the additional subscript 1 for convenience. 3.3.1 Near-W all Cycle In smooth wall flows, the dynamically-important NW cycle is known to comprise alternating streaks of high and low streamwise velocity, with streamwise lengthλ + x ≈ 10 3 and spanwise length λ + x ≈ 10 2 . In addition, these structures are found in the buffer region of the flow, near y + ≈ 15 where the local mean velocity is U + ≈ 10 [see e.g., 48, 98]. To evaluate the effect of riblets on such structures at Re τ = 180, we consider resolvent modes corresponding to k = (κ x ,κ z ,c + ) = (2πRe τ /10 3 ,2πRe τ /10 2 ,10). Previous studies have demonstrated that resolvent modes with these length and velocity scales reproduce known features of the NW cycle [73, 106]. Recent studies have also shown that these structures provide a useful starting point for the evaluation and design of control techniques [68, 69, 81]. Figure 3.2 shows how riblets of varying size affect the gain (singular value, σ k+ ) for the NW resolvent modes. For smooth wall flows at this Reynolds number, a gain-based decomposition of the unmodified resolvent operator (equation ( 2.4)) yields a singular value ofσ k ≈ 3.56. The gain predicted by the extended resolvent operator (equation (3.13)) approaches this value for riblet 48 Figure 3.2: Gain for resolvent modes resembling the NW cycle as a function of riblet spacing. spacing s + = 0.5. Keep in mind that the extended resolvent operator utilizes a different wall- normal discretization compared to the smooth wall case. Further, the volume-penalized governing equations are used in the region occupied by the riblets. Given these differences, asymptotic convergence towards the smooth-wall value as s + → 0 serves as additional verification of our modeling approach. The gain distribution shown in figure 3.2 also reproduces two important trends that are consis- tent with drag reduction curves from previous experiments and DNS. First,σ k+ initially decreases linearly with increasing riblet size. Since the singular values are a measure of energy amplification within the resolvent framework, this reduction in σ k+ can be interpreted as mode suppression due to the presence of riblets. Thus, the monotonic decrease in gain is consistent with previous observations which show that drag decreases linearly with increasing riblet size in the viscous regime [9, 35]. Second, the singular values approach a minimum value for spacing s + ≈ 23, and above this value there is an increase in gain. This critical value for riblet spacing is remarkably consistent with previous DNS results. Specifically, García-Mayoral and Jiménez [ 35] observed that drag reduction saturated nears + ≈ 21 and performance deteriorated as riblet size increased 49 beyond this critical value. The DNS results of Choi, Moin, and Kim [24] also indicate that drag- reducing riblets weaken sweeps and ejections in the near-wall region, thereby reducing turbulent momentumtransferinthewall-normaldirection. Thisobservationisconsistentwiththepredicted suppression of the NW resolvent mode. Of course, there are quantitative differences between the present low-order model predictions and DNS results. For instance, the NW resolvent modes are suppressed by approximately 35% relative to smooth wall values for s + ≈ 21. The maximum drag reduction obtained in DNS by García-Mayoral and Jiménez [35] was< 10%. Further, model predictions show that the gain does not exceed the smooth wall value until riblet size exceeds s + ≈ 38. In contrast, DNS results indicate an increase in drag relative to the smooth wall value for riblets of size s + ≈ 30. One possible reason for this discrepancy at highs + is the breakdown of the viscous flow regime within the riblet grooves. The model used to compute the mean velocity profile for the resolvent operator isnotapplicableinthisregime. Anotherpotentialexplanationforthemoredramaticdeterioration of performance observed in DNS is the emergence of energetic spanwise rollers resembling Kelvin- Helmholtz vortices for large s + . This issue is explored further in the following section. The discussion presented above suggests that the amplification of resolvent modes resembling the NW cycle serves as a useful indicator of performance, even if this measure does not yield quantitative predictions for drag reduction. Next, we evaluate the effect of riblets on the structur e of resolvent modes resembling the NW cycle. Figure 3.3 shows the predicted mode structure over riblets with spacing s + = 8, s + = 21, and s + = 33. The spacings s + = 8 and s + = 33 represent the limiting cases tested in DNS by García-Mayoral and Jiménez [35]. Riblet spacing s + = 21 corresponds to the optimal size before performance begins to deteriorate. The predictions shown in figure 3.3 superpose contributions from resolvent modes with span- wise wavenumberκ z = +2πRe τ /10 2 as well asκ z =−2πRe τ /10 2 . Flow fields for resolvent modes 50 Figure 3.3: Velocity structure for resolvent modes resembling the NW cycle for riblet spacing: a) s + = 8, b) s + = 21, and c) s + = 33. Figures on the left show the wall normal velocity in the streamwise-spanwise plane at y + ≈ 6. Figures on the right depict the flow structure in the spanwise-wall normal plane. The vectors show in-plane velocities (v,w) while the red and blue shading shows the out-of-plane streamwise velocity (u). 51 with either positive or negative κ z would align obliquely with the x axis. In addition to moving downstream, these flow fields would also move in the positive or negative spanwise direction. Su- perposition of both oblique modes ensures that the flow field is aligned in the streamwise direction and propagates downstream. Apart from this imposed spanwise symmetry, the mode structures shown in figure 3.3 emerge naturally from the gain-based decomposition of the volume-penalized governing equations. In other words, the velocity distributions shown in figure 3.3 maximize forcing-response gain for the riblet geometries specified. Keep in mind that the resolvent modes and riblet geometries considered here exhibit differing spatial periodicity. The streamwise periodicity for the resolvent modes is set by the wavenumber κ x ,whilethespanwiseperiodicityisdeterminedbythewavenumberκ z ,andthecoupledharmonics emerging due to interaction with the riblets, κ z ±nκ s , where n ={1,2,...,n h }. The riblets are streamwise-constant and exhibit spanwise periodicity determined just by harmonics ofκ s . Due to this differing periodicity, the spatial windows used to plot figure 3.3 also play an important role in determining the relative distribution of flow features (e.g., high- and low-speed streaks) and the riblet tips or grooves. Predicted mode structures reveal similar mechanisms to those observed by Suzuki and Kasagi [111] and Goldstein and Tuan [38] in numerical simulations and by Lee and Lee [65] in flow visualization experiments, albeit for triangular or scalloped riblets with sharper tips. Specifically, for s + = 8 and s + = 21, the counter-rotating vortices and alternating regions of high and low streamwise velocity associated with the NW resolvent mode are pushed above the riblet tips (figure 3.3(a,b)). There is limited flow penetration into the riblet grooves. Thus, it appears that theseribletseffectivelylimitthewall-normaltransferofmomentumbyblockingthedownwashand crossflows associated with the NW cycle. Recall from figure 3.2 that this change in flow structure for riblets withs + = 8 ands + = 21 is associated with a substantial reduction in amplification for 52 the NW resolvent modes. One potential explanation for this reduction in gain is that riblets limit turbulent energy extraction from the mean flow by pushing the NW cycle away from the wall into a region of lower mean shear. Figures 3.3(a) and (b) also show that mode structure remains relatively narrow-banded over the smaller riblets. In other words, most of the energy associated with this resolvent mode is concentrated in the specified wavenumber-frequency combination, k. Energy content in the coupled Fourier modes k±nκ s is more limited. This concentration of energy is clearly evident in the highly periodic distribution of high and low streamwise velocity observed in the spanwise direction (right-hand side of figure 3.3(a,b)). Physically, this observation is indicative of limited interaction between the NW cycle and the riblets. In contrasttotheflowfeatures discussed above, for ribletswith s + = 33 (figure 3.3(c))the NW resolvent mode is characterized by increased flow penetration into the riblet groove region. Thus, as size increases, the riblets become less effective at limiting the turbulent transfer of momentum in the wall normal direction, which is consistent with previous observations [24, 35, 65]. The flow structure is also noticeably less narrow-banded in the spanwise direction. In other words, energy is distributed more evenly across all the different components of the wavenumber-frequency set k+. Greater energy transfer into the coupled Fourier harmonics is consistent with previous simulation results which suggest that interaction between large riblets and near-wall turbulence leads to the generation of secondary vorticity near the tips [38]. Figure 3.2 shows that fors + = 33, the riblets are also less effective at reducing gain for the NW resolvent mode. Since the predicted mode structure suggests greater flow penetration into the riblet grooves and increased interaction with additional harmonics of the mean flow (equation ( 3.3)), this deterioration of performance can be attributed to greater energy extraction from the mean flow. 53 3.3.2 Kelvin-Helmholtz V ortices Next, we use the resolvent framework to test for the emergence of spanwise rollers resembling Kelvin-Helmholtz vortices over rectangular riblets. This evaluation is motivated by the DNS results of García-Mayoral and Jiménez [35], who observed the emergence of spanwise rollers over rectangular riblets with spacing larger than s + ≈ 24 (ℓ + g ≈ 15). For these cases, spectra for the Reynolds shear stress and wall-normal velocity showed a strong accumulation of energy for structureswithstreamwiselengthscaleλ + x ≈ 150andlargespanwiseextent. Amomentumbalance showed that the Reynolds stress contribution from structures with 65≤λ + x ≤ 290 and λ + z ≥ 130 was responsible for the skin friction increase observed fors + ≥ 24. Further, conditional averaging suggested that these spanwise rollers were centered around y + ≈ 10− 15, while cospectra for Reynolds stress showed a peak at y + ≈ 5. The convective speed associated with these structures was estimated to be c + < 8. To explain the emergence of these spanwise-coherent structures, García-Mayoral and Jiménez [35] pursued inviscid, two-dimensional linear stability analyses using both a piecewise-linear velocity profile and an approximation to the turbulent mean profile. This effort only considered flow in the region between the riblet tips. Modified boundary conditions, informed by viscous flow solutions for the groove region, were used to account for the presence of the riblets in a spanwise-averaged sense. Consistent with DNS observations, these analyses showed a sharp transition in stability as riblet size increased above ℓ + g ≈ 11. However, the streamwise wavelength of the most unstable modes was predicted to beλ + x ≈ 80. García-Mayoral andJiménez[35]attributedthisunder-predictionofλ + x totheinviscidflowassumption, suggesting that viscosity may suppress structures with shorter wavelengths. Figure 3.4 shows resolvent-based predictions for the gain associated with spanwise-constant structures (i.e., κ z = 0) of varying streamwise wavelength and convective speed. Predictions are shown for riblet spacings s + = 8, s + = 21, and s + = 33. For riblets with spacing s + = 8, 54 Figure 3.4: Amplification of spanwise-constant structures as a function of streamwise wavelength and wave speed for riblet spacing: a) s + = 8, b) s + = 21, and c) s + = 33. Red and blue shading represent contours of log 10 (σ k+ ). The black circles show a local maximum. The vertical bar shows where λ + x = 150. figure 3.4(a) does not show an obvious local maximum in amplification. For this drag-reducing case, no spanwise rollers were observed in DNS. A local peak in amplification is observed for s + = 21, as shown in figure 3.4(b). For this case, there is a region of high amplification localized between λ + x ≈ 90−200 and c + ≈ 4−5. The gain in this region of spectral space is an order of magnitude larger than that predicted for the smaller riblets with s + = 8. The most amplified structure has wavelength λ + x ≈ 130 and travels at speed c + ≈ 4.6. As riblet size increases, the region of high amplification spreads in spectral space. As an example, for s + = 24 the high-gain region extends from λ + x ≈ 80−300 and c + ≈ 4−6 (data not shown here). The most amplified structure for s + = 24 has streamwise wavelength λ + x ≈ 200 and convection speed c + ≈ 4.5. Figure 3.4(c) shows that as riblet spacing increases further to s + = 33, the single region of high amplification evident in figure 3.4(b) separates into two distinct high-gain zones. The first region corresponds to structures with streamwise wavelength λ + x ≈ 70, while the second region corresponds to longer structures with λ + x ≈ 450. For a direct comparison with the NW resolvent modes considered in the previous section, figure 3.5 shows how the gain associated with the most amplified spanwise-constant structure 55 Figure 3.5: Maximum gain for spanwise-constant resolvent modes resembling Kelvin-Helmholtz vortices as a function of riblet spacing. The inset text shows the streamwise wavelength for the most amplified mode. with streamwise wavelength 50≤λ + x ≤ 500 and convection speed 2≤c + ≤ 10 varies with riblet size. Consistent with the results shown in figure 3.4, a sharp increase in gain is observed as riblet size increases above s + = 21. This increase in gain also coincides with an increase in the wavelength of the most amplified spanwise roller. A comparison between figure 3.2 and figure 3.5 shows that the amplification of spanwise-constant structures exceeds that of the NW resolvent mode for s + ≥ 24. For s + = 33, gain for the spanwise-constant structure (σ k+ ≈ 33) is an order of magnitude larger than that for the NW resolvent mode (σ k+ ≈ 3). Although the DNS results of García-Mayoral and Jiménez [35] do show that spanwise rollers become more energetic over larger riblets, the dramatic increase in gain observed for riblets with s + ≥ 21 in figure 3.5 is unlikely to be representative. Instead, the observed increase in amplification may be attributed to the eddy viscosity model shown in equation ( 3.24), which assumes that the flow within the grooves remains laminar. For larger riblets, the viscous flow regime in the riblet grooves is expected to break down. This breakdown would lead to greater turbulent momentum transfer into the grooves. The ensuing reduction in mean shear could 56 substantially lower gain for the spanwise-constant structures. As noted earlier, to account for this effect, the eddy viscosity profile would have to include a geometry-dependent origin below the riblet tips at y = 0. The predictions shown in figure 3.4 are in reasonable agreement with DNS results of García- Mayoral and Jiménez [35]. Although the spanwise rollers first become energetic for s + ≈ 24 in DNS, there is evidence of their appearance at s + = 21. The wavelength and convection speed of the most amplified resolvent mode for s + = 21 is slightly lower than that for the structures observed in DNS: (λ + x ,c + ) ≈ (130,4.6) vs. (λ + x ,c + ) ≈ (150,6). However, these gain-based predictions improve upon the inviscid stability analysis results of García-Mayoral and Jiménez [35], which identified structures with λ + x ≈ 80 as being most unstable. This improvement may be attributed to two factors: the present analysis explicitly accounts for the presence of the riblets and includes viscous effects. In addition, the predicted spr e ading of the high-gain region from λ + x ≈ 90−200 for s + = 21 to λ + x ≈ 65−450 for s + = 33 is also consistent with spectra obtained in DNS. Specifically, premultiplied spectra for wall-normal velocity obtained in DNS show that the high-energy region for structures with large λ + z widens in λ + x as riblet size increases. Finally, figure 3.6 shows the predicted flow structure for the highest-gain resolvent mode at s + = 21, i.e., the structure with (λ + x ,c + ) ≈ (130,4.6) highlighted in figure 3.4(b). Predictions shown in figure 3.6(b-d) indicate that this resolvent mode resembles a spanwise-constant roller that is centered near y + ≈ 10 and does not extend beyond y + ≈ 30. Conditionally-averaged flow fields obtainedfromDNSbyGarcía-MayoralandJiménez[35]exhibitedsimilarcharacteristics. Further, figure 3.6(c) shows that the vertical velocity associated with this resolvent mode penetrates well into the riblet grooves. Penetration into the riblet grooves is also evident in the wall-normal velocity contours shown in figure 3.6(a), which exhibit spanwise periodicity that coincides with the riblet spacing. Thus, turbulent structures resembling this high-gain resolvent mode have the 57 potential to generate substantial mixing and momentum transfer in the region occupied by the riblets — as observed in DNS. Figure 3.6: Predicted flow structure for the most amplified spanwise-constant resolvent mode over riblets with spacing s + = 21. a) Distribution of wall-normal velocity in a streamwise-spanwise plane aty + ≈ 5. b) Predicted flow field above the riblet tips. c) Predicted flow field in the middle of the riblet grooves. d) Predicted flow field at an intermediate spanwise location. For b)-d) the vectors show the in-plane velocities (u,v), while the red and blue shading shows the out-of-plane spanwise velocity (w). 58 3.3.3 Riblet Shap e Optimization Figure3.7: Searchforoptimaltrapezoidalribletgeometrieswithheighttospacingratioh/s = 0.5. a) Predictions for the interfacial slip velocity, U + s , as a function of inner-normalized spacing and ridge angle. b) Predicted gain for resolvent modes resembling the NW cycle. Solid black contours show values for the inner-normalized groove length scale,ℓ + g . The assumed geometry is shown on the right. Model predictions presented in the previous two sections show that the extended resolvent framework developed here is able to reproduce trends observed in DNS with limited computation. In particular, the amplification of a single resolvent mode resembling the energetic NW cycle serves as a useful indicator of drag reduction performance. Of course, this single parameter does not capture all the salient physics dictating riblet performance, such as the emergence of spanwise rollers. However, the emergence of high-gain spanwise-constant structures appears to coincide with minimum gain for the NW mode (see figure 3.2 and figure 3.5). Therefore, minimum σ k+ for the NW mode may serve as a useful metric for the purposes of shape optimization. For proof- of-concept, in this section, we pursue limited shape optimization for triangular, trapezoidal, and rectangular riblets using this metric. Figure 3.7 shows predictions for the interfacial slip velocity (a) and NW mode gain (b) for a family of trapezoidal riblets. The height to spacing ratio for these riblets is fixed at h/s = 0.5. The ridge angle and spacing are varied systematically over the following ranges: α∈ [10 ◦ ,90 ◦ ] and s + ∈ [0.125,40]. Atotalof576differentconfigurationsareconsideredinthisfigure. Thisextensive 59 Figure 3.8: Distribution of the optimal riblet size ℓ + g , determined from the square root of the groove cross-sectional area for rectangular, triangular, and trapezoidal riblets. The optimum size is identified based on minimum singular values for resolvent modes resembling the NW cycle. exploration is enabled by the relatively low computational cost associated with evaluating a single mode. Importantly, model predictions are broadly consistent with the results compiled by Bechert et al. [9]. For instance, experimental results show that drag reduction increases as the ridge angle decreases, i.e., as the trapezoids approach a blade-like geometry. Further, for fixed h/s = 0.5, the optimal spacing remained between s + ≈ 16− 19 for all the ridge angles considered. The predictions shown in figure 3.7 follow a similar trend. The amplification of NW modes shows a clear minimum around s + ≈ 20 for all the different ridge angles considered. In addition, for spacing s + ≤ 30, the interfacial slip velocity increases and NW mode gain decreases as the ridge angle,α, decreases. An increase in slip velocity and NW mode suppression are both indicative of improved riblet performance. Note that minimum gain for the NW modes over trapezoidal riblets coincides with the contour ℓ + g ≈ 12 (solid black line in figure 3.7(b)). We have also pursued similar parametric searches for triangular and rectangular riblets. For brevity, model predictions for all of these cases are not reproduced here. Instead, figure 3.8 shows a histogram for the optimal ℓ + g values identified through this process. These ℓ + g values correspond to the spacing that minimizes amplification of the NW resolvent mode for a specific 60 family of shapes (e.g., fixed α for trapezoidal and triangular riblets, fixed t/s ratio for rectangular riblets). The histogram peaks near ℓ + g ≈ 12, which is unsurprising given the results shown in figure 3.7(b). The mean value of ℓ + g for the 54 different optimal geometries represented in this figure is ℓ + g = 13.0±1.7. This mean value is in agreement with the findings of García-Mayoral and Jiménez [35], who compiled many previous experimental and numerical datasets to show that optimal drag reduction for a variety of riblet geometries coincides with a value of ℓ + g ≈ 10.7±1. 3.4 Conclusion Riblets are arguably some of the simplest—and most effective—control solutions proposed thus far for wall-bounded turbulent flows. Simulations, experiments, and theoretical efforts over the past three decades have provided significant insight into the physical mechanisms responsible for drag reduction over riblets. The present chapter leverages these physical insight to generate a low-order modeling framework that can guide the design and optimization of riblet geometry. Specifically, the extended resolvent framework developed here allows us to consider the effect of riblets on individual Fourier modes that resemble dynamically-important features of the flow field—features known to have a direct bearing on riblet performance. Model predictions presented in Section 3.3.1 showed that the variation in gain of a single resolvent mode resembling the NW cycle over riblets reproduced drag reduction trends from DNS. Moreover, the predicted changes in flow structure for this single mode were consistent with prior observations in simulations and experiments [24, 38, 65]. The low computational cost associated with this single mode evaluation enabled the preliminary shape optimization effort presented in Section 3.3.3. Recent work by García-Mayoral and Jiménez [35] shows that the emergence of spanwise rollers resembling Kelvin-Helmholtz vortices is another pathway that leads to riblet 61 performance deterioration. Section 3.3.2 confirmed that the resolvent framework is also able to predict the emergence of such flow structures. Of course, the modeling framework developed here does have limitations. One key limitation is the requirement of a mean velocity profile. The model employed here was able to reproduce prior mean profile observations reasonably well. However, this may not be the case for all flow conditions and riblet geometries. Related to this point, evaluation of model prediction sensitiv- ity to the mean profile is also needed. The predictive mean velocity profile which we compute using the Reynolds and Tiederman [97] enforce smo oth-like conditions and may not faithfully represent the turbulent dynamics. Recent developments of [94] using a similar volume penaliza- tion approach address this problem by using the stochastically forced Navier-Stokes equations to generate corrections to the turbulent eddy viscosity model. Predictions for the second order statistics are computed by solving Lyapunov equation with white white-in-time forcing. Fur- ther, the volume penalization framework employed here does not reproduce the riblet geometry exactly. In particular, the truncated Fourier series representation shown in equation (3.2) does not yield rapid convergence for the sharp volume penalization transitions needed to accurately represent the riblets. This limitation could be remedied by extending the resolvent framework to more sophisticated immersed boundary approaches. Also, we only considered relatively simple narrow-banded surface geometries in this chapter, i.e., spanwise-periodic and streamwise constant riblets. Extension to account for multi-scale or 3D surface features would require inclusion of a greater number of coupled Fourier harmonics in the extended resolvent operator (equation (3.13)), increasing computational expense. Despite these limitations, the model verification and prelim- inary optimization effort described in this chapter shows that the extended resolvent framework can serve as a useful low-order design tool prior to testing in more expensive simulations or experiments. 62 Chapter 4 Resolv en t-Based Optimization of T w o-Dimensional Riblet Profiles 4.1 Ov erview In Chapter 3 we have shown that reduced complexity models based on the resolvent formu- lation reproduce important trends observed over blade-like riblets. Additionally, we conducted a parametric investigation on the effects of spacing ( s + ) and ridge angle (α) for triangular riblets with a fixed aspect ratio ( h/s = 0.5) and trapezoidal riblets. This simple search, based on iden- tifying the gain minimum the near-wall mode, reproduced trends observed in the experiments of Bechert et al. [9], indicating that blade-like riblet geometries are optimal for drag reduction. Using the framework developed in Chapter 3 in this chapter we will search for the optimal the riblet geometry which minimizes the near-wall mode gain. The optimal geometry will depend on the riblet spacing (s) and height (h) and on the overall riblet profile. From prior experimental and numerical work the optimal riblet spacing ranges froms + ≈ 10−20 with the optimal groove length ranging betweenl + g ≈ 6−14 [9, 35, 122]. For our optimization we will constrain the riblet spacing to lie between 10<s + < 30. Our objective function consists of the gain of the near-wall mode computed using the parameters outlined in Section 3.3 with n h = 6 harmonics but the 63 number of Chebyshev nodes has been increased to N = 91. With this node choice, the number of nodes lying within the riblet region is N = 30, see section 3.2.3 for more information. If we consider these points as the parameters over which we optimize, the overall dimension of design (search) space is 32. A brute force search consisting enumerating all parameters by 2 candidate values would result in testing 2 32 ≈ 10 9 possible variations. Furthermore, optimizing over 32 parameters could result in non-physical geometries which couldhavediscontinuitiesandlackanydegreeofsmoothness. IncreasingthenumberofChebyshev nodes used to discretize the wall-normal direction or increasing the number of parameters would exponentially increase the time to complete a brute search. Additionally, increasing the number of candidate values for each parameter would incur unfeasible computational costs. In the following section we will describe how the number of parameters can be significantly reduced by considering geometries which are constructed using Bézier curves. To search for the optimal two-dimensional riblet profile we will implement a nonlinear conjugate gradient algorithm. 4.2 Geometry Construction via P arametric Bézier Curv es The geometries we consider are those that can be described by simple Bézier curves. The use of such curves was popularized by Pierre Bézier who at the time was working as an automobile designer for Renault. In his early years at Renault, Bézier recognized the importance of efficiently and accurately digitizing the hand-drawn models created by designers and draftsman [12]. De- scribing these curves using Bézier approach offers two significant advantages: the description of a curve can be reduced down to a finite set of control points and the curves defined by these control points are locally smooth. Ann th degreeBéziercurvecanbedescribedbythesetof(n+1)controlpointswhichwedenote as [P 0 ,··· ,P n ]. For a two-dimensional curve, each control pointP n represents a vector inR 2 and 64 y z P 0 P 1 P 2 0 y z P 0 P 1 P 2 P 3 0 y z P 0 P 1 P 2 P 3 P 4 0 y z P 0 P 1 P 2 P 3 P 4 0 a) b) c) d) Figure 4.1: Symmetric two-dimensional riblet geometries assembled using Bézier curve(s). On each panel the control points are labeled in order and are denoted using ( ) andP i . Panels a-c) show a quadratic, cubic and quartic Bézier curve respectively. Panel d) shows geometry assembled using a composite Bézier curves. can be described in a Cartesian plane as P n = (z n ,y n ). Here the z-c omp onent and y-c omp onent of the n th control point are denoted using (z n ) and (y n ) respectively. We note that P n is not restricted to R 2 and readily extended to curves defined in R 3 . The points (P 0 ,P n ) correspond to the endpoints of the curve and the curve itself is tangent at these two points. The remaining control points [P 1 ···P n−1 ] generally do not lie on the curve itself. The parametric curve defined by r(ξ) can be represented using a linear combination of Bernstein basis polynomials given by b i,n (ξ). Explicitly this can be written as follows r(ξ) = n X i=0 b i,n (ξ)·P i = n X i=0 n i ξ i (1−ξ) n−i ·P i 0≤ξ≤ 1. (4.1) Here ξ is a parameter defined 0≤ξ≤ 1 and n i is the binomial coefficient which is equal to n i = n! i!(n−1)! . AnimportantadvantageofdescribingshapesusingBéziercurveisthatthesecurves can quickly and easily be transformed by applying a linear transformation matrix (T) onto the control points, i.e.P ′ i =TP i . In figure 4.1 we show a representation of simple quadratic (n = 3), cubic (n = 4) and quartic (n = 5) Bézier curve and how such curves can be used to construct riblet geometries. For our analysis we primarly focus on smooth and symmetric riblet profiles. In figure 4.1 the control points are drawn alongside each curve, and from these points, the convex 65 hull for each curve can be readily identified. Here the convex hull is defined as the smallest convex polygon encompassing all the control points. For a Bézier curve the convex hull represents the domain where the curve is confined to; in other words the curve will exist exclusively in the convex hull. A given Bézier curve can be constrained by controlling and constraining the convex hull. We note that by using this representation geometries with sharp angles and sharp transitions can be difficult to construct given the smooth nature of Bézier curves. Sharp transitions and edges can be achieved either by one of two methods. The first involves modifying the control points on highn th order curves, these transition appear sharp but ultimately are smooth, see figure 4.1(c). The second method is to use multiple lown th order Bézier curves forming a piece-wise continuous curve, see figure 4.1(d). In the following section we will provide details on the optimizer. For our analysis we will forgo the use of high order and piece-wise Bézier curves since both of these methods require additional work to ensure that the resulting Bézier curves do not self-intersect. We instead focus on simpler geometries that showcase the power of our reduced complexity model. 4.3 Optimization Problem For the remainder of this chapter will focus on symmetric riblet geometries constructed using a cubic Bézier spline. For a cubic curve, four control points are required and in two-dimensions, these four controlpointsencompass eightindependentvalues. Since aBézier curveis encapsulated by its convex hull, and since the convex hull is defined by the control points, we can constrain the size and shape of the riblet profile by constraining the control points. In this representation the height of the riblet is given by y-c omp onent of the last control point, i.e. P 3 = (z 3 ,y 3 ) = (z 3 ,h). To control the height of the convex hull and constrain the vertical translation of the curve the y-c omp onent of the first control point is fixed to 0, i.e. P 0 = (z 0 ,0). Additionally, we enforce that y-c omp onent of the remaining control points remains below the riblet height, i.e.y i ≤h. The last 66 two constraints we impose is that the control points do not intercept the axis of symmetry and remain above our reference point i.e.z i ≥ 0 andy i ≥ 0. Having these constraints we can formally state our optimization problem. min f(η) =σ k+ (η) subject to Aη≥b (4.2) where η = n s, h =y 3 , ··· , y 0 , z 3 , ··· , z 0 o . (4.3) Since all the inequality and equality constraints are linear our constraints can be represented by a matrix equation. Here z and b are vectors in R m and A is matrix in R p×m . Here m represents the size of our parameter space and p represents the number of constraints on our system. We will consider two configuration, for the first case we constrain the base width to be equal to the riblet spacing, see figure 4.2(a), and in the second case we allow the base width of the riblet to be independent of the spacing, see figure 4.2(b). In the second configuration, the size of the parameter space is m = 7. In the first case, in where the base width is equal to the riblet spacing, i.e. x 0 = 1 2 s, the size of the parameter space is m = 6. This problem can be solved by dealing with the reduced parameter space but it can also be solved using the original parameter space of m = 7. With the full parameter space, this constraint can be imposed by including it into the active constraint matrixA c . Since we are dealing with inequality constraints the active constraint matrixA c , at every step (k) tracks the rows of matrixA in which (Aη k ) i = (b) i . Here we use the notation ( ) i to denote the i th component of a vector. In the following section, we provide the details of our optimization technique. 67 y z P 0 P 1 P 2 P 3 s y z P 0 P 1 P 2 P 3 s a) b) Figure 4.2: Schematic of the two riblet configurations used with our optimization scheme. Panel a) shows a configuration in where base thickness is equal to the riblet spacing and panel b) shows a configuration where the base thickness is independent of the riblet spacing. 4.4 Nonlinear conjugate gradien t metho d To find the optimal geometry we implement an optimization scheme based on the nonlinear conjugate gradient method. Similar to the linear conjugate-gradient (CG) method the search direction is determined as a linear combination of gradient of the objective function and the prior search direction. The linear CG method can be viewed as an iterative method to solve a linear system of equations. For example consider the system of equations given Qη = b, in which Q is strictly positive definite and thus non-singular. This system can be solved by minimizing the convex quadratic function given byf(η) = 1 2 η | Qη−b | η. In this method the new search direction p k is dependent on the residual of the linear systemr k =Qη k −b. In the context of our quadratic objective function the residual is equal to the gradient of the objective function, i.e.r k =∇f(η k ). An important distinction between the conjugate gradient and the steepest descent method is in how search direction is computed. In the steepest descent method, at every iteration, the search direction is equal to the gradient of the objective function. In the CG method the search direction is updated using information from the prior step, p k =−r k +β k p k−1 =−∇f(η k )+β k p k−1 . (4.4) 68 Here, β is chosen to ensure that each new search direction is orthogonal to the previous search direction with respect to Q, i.e. p | k Qp k−1 = 0. In other words this property ensures that each new search step is conjugate—referred to as Q-orthogonal— with the prior steps. Using this definition, β takes the following value β k = r | k Qp k−1 p | k−1 Qp k−1 (4.5) Having an updated value ofp k the new search direction can be updated such that η k+1 =η k +αp k . (4.6) Here α i can be interpreted as the step length, and is determined by solving the one-dimensional minimization problem, min α k f(η k +α k p k ). If the problem is quadratic in which Q is defined, α i has an analytical solution which is given by α i =p | k r k /p | k Qp k . If the objective function is not quadraticthesteplengthcanbedeterminedviaasimplelinesearchalgorithm. Fortheseproblems Qmaynotbereadilyavailablesuchasinourcase,seesection4.3. Byusingtheconjugacyproperty ofp k ,η k andr k , β can be calculated independently ofQ and has the following form β k = r | k Qp k−1 p | k−1 Qp k−1 = r | k (r k −r k−1 ) p | k−1 (r k −r k−1 ) = r | k r k r | k−1 r k−1 . (4.7) The advantage of the conjugate gradient method is in being able to generate a new conjugate direction at every step k using the information of the prior step k−1. The nonlinear conjugate gradientmethodextendstheideologybehindthelinearCGmethodtononlinearproblems. Similar 69 to its linear counterpart, it requires the computation of the residue and for a nonlinear and non- quadratic problem the residue is computed usingthe gradientof theobjectivefunction. Compared to other methods, the computation of higher-order derivatives is not required. This is important given that our objective function (i.e the singular value) is not analytical and each derivative will need to be approximated using a finite difference approach. Compared to other optimization methods it has low storage requirements, it does not include or depends on any matrix inversion operations and at most requires minimal and rudimentary vector operations. The distinctions between the linear CG and nonlinear CG lies in howβ k is calculated and in the line search method implemented to solve for the step length. Below we provide a skeleton of our optimization scheme. In the following section we will describe the intricacies of each step in detail. 70 Algorithm 1: Nonlinear CG with Restart Evaluate: f 0 =f(η 0 ), g 0 =−∇f(η 0 ); Setp 0 ←g 0 and k = 0; while ∥g k ∦= 0 or exit c ondition do Compute α k using inexact Line Search; Updateη k+1 =η k +α k p k ; Compute active constraintsA c ; Evaluateg k+1 =−∇f(η k+1 ); Compute β PRP k+1 = (g k+1 −g k ) | g k+1 g | k g k ; β k+1 ← min(0, β PRP k+1 ); p k+1 ← g k+1 +β k+1 p k ; if A c p k ≥ 0 then Update active constraint matrixA c ; Compute feasible direction using the projection matrixT; p k ← Tp k ; g k ← Tg k ; end end 4.4.1 Line Searc h Metho d In the above algorithm, the first search direction is similar to the steepest descent method in where the first search direction is equal to negative gradient of the objective function. Having this direction we need to compute a suitable step size which is determined by performing an inexact line search. For this problem we are interested in finding the value of α k that loosely minimizes f(η k +α k p k ). The computed α k should provide a sufficient decrease in the objective function. 71 The Wolfe conditions are employed to ensure that for the given α k a sufficient decrease in the objective function is reached. Wolfe conditions impose the following constraints on α k f(η k +α k p k )≤f(η k )+µ 1 α k p | k ∇f(η k ). (4.8) p | k ∇f(η k +α k p k ) ≤−µ 2 p | k ∇f(η k ) (4.9) In the equations aboveµ 1 andµ 2 are arbitrary constants such that 0<µ 1 <µ 2 < 1. The first Wolfe condition can be interpreted as requiringf(η k +α k p k ) to lie beneath the secant line located at pointf(η k ) and having a slope equal toµ 1 α k p | k ∇f(η k ). As the value ofµ 1 approaches 1, the secant line approaches the tangent line atf(η k ), and the conditions for a sufficient decrease in the value of the objective function become more stringent. At the opposing end as µ 1 approaches 0 the secant line becomes more horizontal and the conditions for sufficient decrease in the objective function are relaxed. In order to ensure that the step lengthα k is not too small, the second Wolfe condition, see equation (4.9), is implemented. The second Wolfe condition can be interpreted as requiring a sufficient decrease in the gradient of the objective function, since if α k is too small then the gradient at stepsk andk+1 are approximately equal, i.e.∥∇f(η k+1 )∥≈∥∇f(η k )∥. In our optimization schemeµ 1 = 1e−3 andµ 2 = 0.90, for more information on the Wolfe conditions or on the choice ofµ 1 orµ 2 the readers are referred to Nocedal and Wright [82] and Griva, Nash, and Sofer [42]. 4.4.2 F easible Direction Since we are dealing with a constrained optimization problem we need to ensure that the search directionη k+1 does not violate the active set of constraints. At every iteration, the search direction is checked against our set of constraints given by Aη k+1 ≥b. If any component of the search direction reaches the boundary defined by our constraints, such that (Aη k+1 ) i = (b) i , the 72 i th row of the matrix A and i th component of b is added to the set of active constraintsA c and b c . At every subsequent iteration, the product of A c p k ≥ 0 is checked to ensure that the new search direction does not violate the imposed constraints. If a constraint is active (A c η k =b c ) it is sufficient to check the product of A c p k , since α k > 0 and A c (η k +α k p k )≥b c , A c η k =b c , → A c p k > 0 (4.10) If the component of (A c p k ) i > 0, the constraints given by those components remain active. If any component of (A c p k ) i < 0, that component and its corresponding row are removed from the active constraint matrix. Having updated the active constraints we must ensure that the search direction lies within the feasible search space. If the new search direction violate the active constraints, i.e.A c p k > 0, the search direction must be projected to lie within the feasible search space. In order for the search direction to be valid we require that A c p k = 0. In order words we require that the new search direction to lie within the null space of A c . To achieve this the search direction is multiplied by the orthogonal projection matrix defined as T =I−A | c (A c A | c ) −1 A c . (4.11) By doing this we ensure that the projected search direction given by v k = Tp k lies in the null space ofA c . In other words A c v k =A c Tp k =A c (I−A | c (A c A | c ) −1 A c )p k = 0. (4.12) 73 The projection matrix can also be computed in terms of the basis Z which spans the null space of A c , i.e. Z = null(A c ). Having this basis the projection matrix can be computed as T = Z(Z | Z) −1 Z | . If the basis is orthonormal the projection matrix can quickly be computed via T =ZZ | . 4.4.3 Conjugate Direction β , and Restart Pro cedure One of the important distinctions between the nonlinear CG and linear CG method is in the method that the search direction is updated. In the linear CG method, β is chosen such that the new search direction is conjugate to the previous search direction. The form of β detailed in section 4.4 produced a new conjugate search direction at every iteration. For a nonlinear problem, the residue in Section 4.4 is replaced with the gradient of the objective function and has the following form β k = g | k g k g | k−1 g k−1 . (4.13) In the equation above g k =−∇f(η k ). The form of β provided above was first introduced by Fletcher and Reeves [31] to tackle nonlinear objective functions that can be well approximated by quadratic matrix function in the vicinity of the optimum value η ∗ . Modifications have been introduced to improve the convergence rate of the Fletcher-Reeves (FR) method. The most important of these modifications include periodically restarting the algorithm such that all prior information on the conjugate steps are forgotten. In this restart procedure, the value of β k is set to zero, and the subsequent search direction is equal to the direction of steepest descent, p k =−∇f(η k ). The simplest restart procedure involves restarting the algorithm after m or m+1 steps—herem is the size of the vectorη—and the convergence rate improves and becomes superlinear [27]. If no restart procedure is implemented the rate of convergence becomes linear 74 [26, 91, 92]. In our optimization scheme, we have opted for the modified Polak-Ribière-Polyak (PRP) method in where β is given by β PRP k = g | k (g k −g k−1 ) g | k−1 g k−1 , and (4.14) β k = max([0, β PRP k ]). (4.15) Equation (4.14) represent the method proposed Polyak [89] in 1962 and Polak and Ribière [88] in 1969 which we will refer to here as the PRP method. The PRP method described by equation (4.14) has been augmented using equation (4.15), suggest by Powell [90]. This provides an automatic method for restarting the algorithm and additionally it ensures that p k remains a direction of descent[82]. If our nonlinear function can be well approximated by a quadratic func- tion, the matrixQ is equal to the Hessian of the objective function and the conjugacy condition can be rewritten as p | k ∇ 2 f(η k )p k−1 . The PRP method avoids computing the hessian, since it’s inconvenient and computationally expensive, and instead relaxes the conjugacy condition such that p k (∇f(η k )−∇f(η k−1 )) = 0. This condition can be derived using a Taylor expansion of ∇f(η k−1 ) or by invoking the mean value theorem. For quadratic problems which are highly convex, the gradients g k and g k−1 are orthogonal and the value of β determined by the PRP and FR methods are equal. For nonlinear problems in where an inexact line search algorithm is implemented, the PRP method has been shown to be more robust and efficient than the FR method [82, 90]. 75 4.4.4 Appro ximating Deriv ativ es and Stopping Conditions Given that our cost function does not have an analytical form, the partial derivatives of the objective must be approximated. To do this, a finite difference approximation is implemented in where ∇f(η) | e i = ∂f ∂η i = f(η +δe i )−f(η−δe i ) δ +O(δ 2 ) (4.16) Equation (4.16) provides an approximation to the partial derivative of the objective function with respect to i th component of η. In equation above δ represents the step or perturbation length ande i is the unit vector in the direction corresponding toi th component ofη. We have chosen to use a central difference approach for the smaller truncation error at the cost of two evaluations of the objective function. Here the step size is chosen to be δ =ϵ M 1/3 , in where ϵ M represents the machine error. This choice of δ minimizes the total error consisting of the sum of the machine error and truncation error. In our case, for a machine performing double-precision floating-point operations, the optimal step size corresponds to δ =ϵ M 1/3 ≈ 6.0·10 −6 . The stopping conditions for the optimizer have been chosen if either the norm of gradient or residuals are smaller than some arbitrary value γ. In our optimization cases the stopping condition is set if∥g k ∥ 2 ≤ γ = 1.0·10 −6 . Additionally, the optimizer reaches an exit condition if the residual of the objective function, defined as R k = f(η k )−f(η k−1 ) max(1,f(η k−1 )) remains and saturates below R k ≤ 1.0·10 −10 for more than one consecutive iterations. 4.5 Results Figure 4.3 and figure 4.4 show the results from our optimizer for two different starting ge- ometries. For this case, we have chosen to minimize the gain of the near-wall mode described in Section 3.3.1. Furthermore we constrain the maximum height to h + = 15 and from preliminary 76 results, at the optimal profile results in an aspect ratio of h s ≈ 0.75. The choice of h + = 15 is arbitrary and has been implemented to combat the monotonic relationship between riblet height and near wall mode gain and improve convergence. We will describe the effects and consequence of limiting the height in the discussion below. In figure 4.3(a) we display the normalized gain for the near-wall resolvent mode plotted versus the steps taken by the optimizer and in panel (b) we plot the residual R k at every step. In figure 4.3(c) we display the initial condition/geometry fed into the optimizer. After approximately 14 steps the reduction in the gain appears to saturate at approximately 0.512, the residual for k > 14 hovers around R k = 10 −3 − 10 −4 . In figure 4.3(d) we show the geometry that the optimizer produces at k = 14. After approximately 30 steps the residuals appear to level off, where R k < 10 −10 . The geometry reached after the residuals have leveled off is displayed in figure 4.3(e). The residual shows that there are earlier steps that reach R k < 10 −10 . These iterations correspond to cases for which the inexact line search algorithm chooses a small step size so that the new search directionη k+1 does not violate the active constraints. Inspecting the final geometry produced by the optimizer, see figure 4.3(e), we note that the overall geometry is visually indistinguishable from the geometry produced at step k = 15, see figure 4.3(d). Both geometries feature a moderate radius at the base which thins out near the riblet tip. The radius at the base is constrained but this is a consequence of using a cubic Bézier spline. The minimum radius at the base is never reached. The minimum base radius is achieved when both P 1 and P 2 lie at the origin. The normalized gain achieved at the end of step 15 is σ k+ /σ k,s = 0.516 and at the end of step 36 is σ k+ /σ k,s = 0.511. The improvement achieved between both steps is less than 0.60% suggesting that the decreases in the base radius and the slight increase in groove area is not significant. 77 Figure 4.3: Optimized geometry determined using the algorithm outlined in Section 4.4. Panel (a) shows the gain normalized by the smooth wall cases,σ k+ /σ k,s plotted verses iterations of the optimizer. Panel (b) plots the residual R k at each iteration, panel(c) shows the initial geometry fed into the optimizer. Panels (d) and (e) and shows the geometry produced after 15 and 38 iterations. In panels (c) -(e) the control points are drawn alongside their Bézier curves. In figure 4.4 we repeat the optimization process but change the initial starting geometry. In this case the same trends seen in figure 4.3(a,b) are observed in figure 4.4(a,b). After 8 iteration the normalized gain appears to saturate around σ k+ /σ k,s = 0.510, representing a 49% reduction in gain of the near-wall mode. For this case the residuals begin to level off after 22 iterations remaining below R k < 10 −10 . The geometries produced after 8 and 28 iterations are plotted in figure 4.3(d,e), and similar to the case detailed above, the two profiles are indistinguishable from one another. Compared to the final profile shown in figure 4.3(e), the geometry shown in 78 Figure 4.4: See figure 4.3 for caption. Panels (d) and (e) show the geometry after 8 and 32 iterations. figure 4.4(e) features a smaller base radius and slightly larger groove area. The normalized gain at the end of stepk = 28, corresponding to the end of the optimization scheme, isσ k+ /σ k,s = 0.508, an improvement of less than 0.20% over stepk = 8. Comparing the final geometries produced by the optimizer, see figure 4.3(e) and figure 4.4(e), the difference is gain between both geometries is less than 0.30% confirming our statement above whereby this small additional decrease in gain produced achieved by decreasing the base radius is insignificant. In both cases, the optimal geometry resembles a thin-blade like riblet; the use of a higher-order Bézier curve would yield a smaller base radius and more pronounced blade-like profile. 79 For the results presented above a sensitivity analysis is performed comparing the optimal riblet geometries determine for modes which have a smaller spanwise wavelengths corresponding toλ + z = 70 andλ + z = 80. Similar to the prior cases, the maximum height was limited toh + = 15 and the optimum spacing corresponds to an aspect ratio of h s ≈ 0.75. For these cases, the streamwise wavelength and wave mode speed was kept constant and matched that of the near- wall mode. In figure 4.5 we present the optimized geometries determined using the algorithm outlined in Section 4.4. Similar to the results presented above, the optimal geometry resembles a blade-like riblet, in where the optimal riblet spacing is approximatelys + ≈ 18−21. For structures withλ + z = 80 the maximum reduction in gain achieved using this geometry is approximately 42% and λ + z = 70 the maximum reduction is approximately 37%, see figure 4.5(a,b). Comparing the riblet spacing for the modes with a smaller spanwise wavelength corresponding to λ + x = 70 and λ + z = 80 the optimal spacing lies between 18.9 and 19.9 wall units. The overall geometries remain similar, see figure 4.5, but for the smaller wavelengths, the optimizer predicts a 3-9% reduction in riblet spacing. The riblet length, measured using the groove cross-sectional area, is approximately l + g ≈ 16.0−17.0 which is significantly larger than the optimal groove length determined by García- Mayoral and Jiménez [35]. This over-prediction is a result of only modeling the near-wall mode and capping the maximum height at h + = 15. For the near-wall mode the optimal geometry maximizes the streamwise compliance and max- imizes the protrusion height ∆ l . If we momentarily ignore the drag degrading effects of the KH-like spanwise rollers, maximizing the streamwise compliance and the protrusion height ∆ l can hypothetically be achieved by spacing the riblets further away from each other, up-to-the optimum spacing, or by increasing the maximum height. Across the various modes analyzed the optimal riblet spacing has remained constant at arounds + ≈ 20, and thus in order to increase the protrusion height the riblet height needs to increase. In the context of our optimization scheme 80 Figure 4.5: Optimal geometries determined for wavelengths corresponding to, a) λ + z = 80, b) λ + z = 90 and c)λ + z = 100. The control pointsP 0 toP 3 are drawn alongside the Bézier curve they correspond to. The riblet spacing and reduction in the gain are draw above each configuration. and our objective function which is given by the near-wall mode gain, the height is non-convex; a reduction in gain can be achieved by increasing the riblet height. This is because our objective function does not account for the emergence of modes resembling Kelvin-Helmholtz rollers over riblets with largel + g . To demonstrate this the optimizer was repeated for the near-wall mode but the maximum height was limited to h + = 10. From prior results the optimal spacing is s + ≈ 20 so at this given h + the aspect ratio is h s ≈ 0.5. This choice has been motivated by the results of Endrikat et al. [28] suggesting that for riblets with an aspect ratio of h s > 0.5 the primary drag breakdown mechanism is the appearance of spanwise rollers resembling Kelvin-Helmholtz vortices. Rerunning the optimizer for this scenario ( h s ≈ 0.5) the minimum gain is σ k+ /σ k,s = 0.536 and the optimal riblet length is approximately l + g ≈ 13.38. Comparing this to our original case with h s ≈ 0.75, for a 50% increase in height (10→ 15) the reduction in gains is marginally better at σ k+ /σ k,s = 0.508. For modes with a smaller spanwise wavelength of λ + z = 70, the minimum mode gain for h s ≈ 0.5 was σ k+ /σ k,s = 0.649, which is marginally better than the minimum gain achieved with an aspect ratio of h s ≈ 0.75. Similar to the gain, the optimal riblet spacing for 81 these two configurations shows a marginal difference; for h s ≈ 0.75 the optimal riblet spacing is 18.9 wall units and for h s ≈ 0.5 the optimal riblet spacing is 18.3 wall units. Figure 4.6: Optimized geometry determined using the algorithm outlined in Section 4.4, in where the base width is independent of the riblet spacing. In both panels the normalized mode gain, the optimal riblet spacing (s + ) and the base width (s + b ) are shown along the optimal geometry. In panel(a) the maximum height was limited to h = 15 corresponding to an aspect ratio of h s ≈ 0.5 and in panel (b) the maximum was limited to h = 10 corresponding to h s ≈ 0.5. For our final case, we consider riblet geometries where the base width is independent of riblet spacing, see figure 4.2(b). Given the monotonic relationship between riblet height and the near- wall mode gain, we have limited the height to h + = 15 and h + = 10. Similar to the cases above this yields an optimal aspect ratio of h s ≈ 0.75 and h s ≈ 0.5. In figure 4.6 we show the results produced by our optimizer for these cases. For both cases, the optimal geometry is similar to the geometries presented in figure 4.5(c), where the optimal profile resembles a thin blade-like riblet with a spacing of s + ≈ 20. Comparing the results shown in figure 4.6(a) and figure 4.5(c), the minimum gain reached varies by less than 10 −4 and the optimal spacing increases by less than 1%. The results for h s ≈ 0.5 display marginal differences in terms of overall geometry and optimal spacing and reaffirm our statement above in where the additional decrease in gain achieved by a geometry which is 50% taller is insignificant. For the original cases tested with h s ≈ 0.75 and h s ≈ 0.5 the optimal riblet length was l + g ≈ 17.00 and l + g ≈ 13.4 (see table 4.1 for a summary of 82 λ + z σ k+ /σ k,s s + s + b l + g Aspect Ratio h s ≈ 0.75 case 1 100 0.511 21.19 − 17.00 case 2 100 0.508 20.72 − 17.20 case 3 80 0.576 18.90 − 16.03 case 4 70 0.626 18.90 − 16.42 case 5 100 0.508 20.89 13.56 17.21 Aspect Ratio h s ≈ 0.5 case 6 100 0.536 19.48 − 13.38 case 7 70 0.649 18.30 − 13.07 case 8 100 0.534 20.73 12.44 14.08 Table 4.1: Summary of the various optimization cases executed. Optimal spacing (s + ), base width (s + b ) and riblet length (l + g ) for h s ≈ 0.5 and h s ≈ 0.75. these cases). For the cases presented in figure 4.6 the optimal riblet length is similar atl + g ≈ 17.21 and l + g ≈ 14.08. In terms of the mode gain, the increased streamwise compliance gained by the increasing the groove area does not lead to a significant improvement. The largest improvement observed was of 0.3% in the cases with aspect ratio of h s ≈ 0.5. However, as noted previously, the increase in groove area could lead to the emergence of energetic spanwise rollers that degrade drag reduction performance. 4.6 Conclusion In this chapter we have shown how the framework developed in Chapter 3 can be combined with a nonlinear conjugate based optimizer to generate predictions for optimal riblet geometries which maximally suppress the near-wall mode. We have detailed how locally smooth and easily scalable two-dimensional riblet profiles can be constructed using low-order Bézier splines. The benefit of using such splines lies in being able to describe a smooth curve using a finite set of control points. This framework including the optimization framework can easily be extended to account for riblet profile assembled using b-spline or a composite spline formed using a set of n th 83 order Bézier splines. For our optimization scheme we have opted to use a nonlinear conjugate based optimizer with a non-exact line search method. Given the nature of the model, a reduced model for the near-wall model introduces certain issues that make our objective function non-convex. One issue, which we have described in detail involves the monotonic relationship between mode gain and riblet height. However, for our analysis concerning the optimization of riblets with heights of h + = 15 and h + = 10 the differences in the gain for the near-wall mode are marginal and insignificant especially considering the 50% increase in height when going from h + = 10 to h + = 15. The optimal riblet profiles which we have identified resemble a blade-like profile that maximizes the streamwise compliance. In the context of the protrusion height concept of [67] and Bechert et al. [9], our riblet profiles maximize the protrusion height available given the constraints imposed. To optimize the riblet height, our objective function would need to include surrogate modes for the near-wall mode in conjunction with modes modeling the effect of the Kelvin-Helmholtz like rollers observed in riblets with aspect ratio h s > 0.5. Although our model has some limitations, the optimized profiles extracted using this frame- work are in line with prior experiments and predictions. Predictions based on the near-wall mode and larger spanwise modes suggest that the optimal spacing for thin-blade like riblets is s + opt ≈ 18− 20 with a maximum mode suppression of 40-50% observed. The optimal groove length based on these cases ranges from l g opt +≈ 13−17 but the values derived may not be the optimal values since the profile height has not been optimized. More advanced models can in- clude the usage of b-splines, non-symmetric riblet profiles, and the inclusion of modes resembling Kelvin-Helmholtz vortices to predict optimal heights. Our modeling and optimization framework is flexible and such changes can be implemented with ease. 84 Chapter 5 Resolv en t-Based Predictions for T urbulen t Flo w Ov er Anisotropic P ermeable Substrates 5.1 In tro duction In this study, we provide an extension of the resolvent framework to analyze the effects of anisotropic permeable substrates on turbulent flows. Predictions developed by this framework are compared with the simulations of Gómez-de-Segura and García-Mayoral [40] over streamwise preferential anisotropic permeable substrates. The effect of permeable substrates is introduced in the resolvent framework using the Volume- Averaged Navier-Stokes equations and a generalized form of Darcy’s law. Predictions for the mode gain are computed using a synthetically mean velocity profile which is generated using an eddy viscosity model. Additionally, a comparison between the mode gain determined via this method and the mode gain computed using mean velocity profiles obtained from the numeri- cal simulations of Gómez-de-Segura and García-Mayoral [40] is presented. Substrates with high streamwise permeability and low spanwise permeability are found to suppress the forcing-response gain for the resolvent mode that serves as a surrogate for the energetic near-wall cycle. For drag- reducing cases in which p K + y < 0.4 the mode gain calculated using DNS and the synthetic mean 85 velocity profiles is indistinguishable. Additionally, the reduction in mode gain is shown to be consistent with the drag reduction trends predicted by theory and observed in numerical simula- tions. Simulation results suggest that the maximum drag reduction is limited by the emergence of spanwise rollers resembling Kelvin-Helmholtz vortices appearing above a threshold value of the wall-normal permeability. The resolvent framework also predicts the conditions for which these energetic spanwise-coherent rollers emerge. The framework presented here is the first fully predictive, physics-based modeling framework that reproduces trends observed in high fidelity simulations over drag-reducing anisotropic per- meable substrates. These findings suggest that a limited set of resolvent modes can serve as the building blocks for computationally-efficient models that enable the design and optimization of permeable substrates for passive turbulence control. Preliminary results and development of this predictive framework and its use in the model-based design of porous substrates can be found in Chavarin et al. [19] and Chavarin et al. [20]. 5.1.1 Con tribution and Outline In this study, we develop a reduced-order modeling framework grounded in resolvent analysis [72, 73] that can be used to predict the effect of substrates with known permeability on the NW cycle and test for the emergence of KH rollers. Under the resolvent formulation, the Navier-Stokes equations are interpreted as a forcing-response system in which the nonlinear convective terms are treated as the forcing to the linear system that generates a velocity and pressure response. A gain-baseddecompositionoftheresolventoperator, whichisthelineartransferfunctionthatmaps the nonlinear forcing to the velocity and pressure response, is used to identify high-gain forcing and response modes across spectral space. Specific high-gain response modes ( r esolvent mo des) have been shown to serve as useful models for dynamically-important flow features such as the 86 NW cycle [78]. This means that, as a starting point, the effect of any control can be evaluated on these individual resolvent modes instead of the full turbulent flow field. Indeed, previous studies show that the gain for the NW resolvent mode is a useful predictor of drag reduction performance for both active [68, 81, 116] and passive [21, 69] control techniques in wall-bounded turbulent flows. In particular, recent work by Chavarin and Luhar [ 21] shows that the gain for the NW resolvent mode is a useful surrogate for total drag reduction in turbulent flows over riblets. Riblet geometries that lead to drag reduction in experiments and high-fidelity simulations are found to reduce the forcing-response gain for the NW resolvent mode relative to its smooth wall value. In addition, Chavarin and Luhar [21] show that the resolvent framework is also able to predict the emergence of spanwise rollers resembling KH rollers over certain riblet geometries, which is consistent with previous DNS results [35]. Motivated by these prior modeling successes, we consider the effect of anisotropic porous substrates on resolvent modes resembling the NW cycle and spanwise-coherent structures resembling KH rollers. To enable a direct comparison with the simulation results of Gómez-de-Segura and García-Mayoral [40], we consider a symmetric channel geometry at friction Reynolds number Re τ = 180 and substrates with identical wall- normal and spanwise permeabilities, i.e., substrates with ϕ yz = p K + y / p K + z = 1. We model the flow in the substrate using the volume-averaged Navier-Stokes equations, in which the effect of the permeable substrate is included via a permeability tensor. However, this modeling framework can be extended to include more sophisticated interfacial boundary conditions [61, 62], and to account for inertial effects via the so-called Forchheimer term [ 15, 128]. The remainder of this chapter is structured as follows. Section 5.2 describes the permeable substrate model used here, the extended resolvent formulation, as well as the numerical methods used for resolvent analysis. Section 5.3 presents model predictions for the effect of anisotropic permeable substrates on the NW resolvent mode as well as spanwise-coherent resolvent modes 87 resembling KH rollers. These predictions are compared against DNS results from Gómez-de- Segura and García-Mayoral [40]. We also pursue a limited sensitivity analysis of model predictions to the exact form of the mean profile used to construct the resolvent operator. Specifically, we compare predictions made using a synthetic mean profile that is computed using an eddy viscosity model against predictions generated using the mean velocity profiles obtained in DNS by Gómez- de-Segura and García-Mayoral [40]. Section 5.4 concludes the chapter. 5.2 Metho ds In this section, we present the equations used to model flow in the porous medium (Sec- tion 5.2.1), briefly describe the resolvent formulation and discuss its extension to account for per- meable substrates (Section 5.2.2), present the boundary conditions imposed at the fluid-substrate interface(Section5.2.3), discussthemeanvelocityprofilesusedtoconstructtheresolventoperator (Section 5.2.4), and describe the numerical method used to implement the analysis (Section 5.2.5). 5.2.1 A ccoun ting for P ermeable Substrates The resolvent framework is reformulated using the volume-averaged Navier-Stokes equations. Volume-averaging gives rise to two additional terms: a term representing the sub-filter scale stresses and a surface filter term that accounts for the force exerted by the solid phase of the permeable medium on the fluid phase [ 125, 126, 128]. A typical closure model for the surface filter term involves parameterizing the flow resistance using the Darcy permeability tensor and the so-called Forchheimer correction term that accounts for inertial effects [ 126]. This model has beenusedinpreviousnumericalsimulationsofflowoverporoussubstrates[ 15, 99, 100]aswellasin linear stability analyses [114, 115]. The volume-averaged Navier-Stokes equations and continuity 88 constraint for flow through a porous medium with porosity ε, dimensionless permeabilityK, and dimensionless Forchheimer resistanceF can be expressed as: ∂⟨u⟩ ∂t + 1 ε ∇· ε⟨u⟩⟨u⟩+ετ =− 1 ε ∇(ε⟨p⟩)+ 1 εRe τ ∇ 2 (ε⟨u⟩)− ε Re τ K −1 (I +F)⟨u⟩ (5.1a) and ∇·(ε⟨u⟩) = 0. (5.1b) Intheequationsabove,⟨·⟩representsavolume-averagedquantity,⟨u⟩isthedimensionlessvelocity, ⟨p⟩ is the dimensionless pressure, and t is dimensionless time. The equations presented above have been normalized using the channel half-height h † and the friction velocity u τ . The friction Reynolds number is given by Re τ = u τ h † /ν and the dimensionless permeability defined as K = K † /h †,2 , whereK † is the dimensional permeability. The quantityτ =⟨uu⟩−⟨u⟩⟨u⟩ represents the sub-filter scale stresses which arise from volume averaging the Navier-Stokes equations. For our analysis we consider the following simplifications. Consistent with prior numerical simulations [40, 99], we focus on substrates characterized by a permeability tensor of the form K = diag(K x ,K y ,K z ) with the ratio of the wall-normal and spanwise permeabilities set to unity, i.e.,K y =K z . Notethatthepermeabilitytensorissymmetric, andsoaneigenvaluedecomposition can be used to identify its principal values and directions (or axes). The assumed form of the permeability tensor implies that its principal directions align with the streamwise, wall-normal, and spanwise directions of the flow. Second, we omit the nonlinear Forchheimer correction term, F, that is used to account for inertial effects in flows through porous media [ 15, 83, 84, 128]. This assumption is again consistent with the numerical simulations we will use to test model predictions. Third, we assume that the porosity of permeable substrates is ϵ≈ 1 to maximize any potential drag reduction [1]. Finally, since we are primarily interested in structures that are 89 y =(1+H) y =1 y =−(1+H) y =−1 y x H P ermeable Substrate F ree c hannel U(y) Figure 5.1: Schematic showing the symmetric channel flow configuration considered in this chap- ter. . much larger than the characteristic length scale of the porous medium (i.e., NW cycle and KH rollers), we neglect the sub-filter scale stresses [ 15]. With these assumptions, (5.1a) and (5.1b) can be expressed as: ∂u ∂t +∇· uu =−∇(p)+ 1 Re τ ∇ 2 u− 1 Re τ K −1 u (5.2a) and ∇·u = 0, (5.2b) where the⟨·⟩ notation has been eliminated for simplicity. The unobstructed fluid domain is characterized by infinite permeability. In this region, the permeability term goes to zero and (5.2a) reduces to the standard Navier-Stokes momentum equation. Figure 5.1 shows the symmetric channel flow configuration considered in this study. The unobstructed region corresponds toy∈ [−1,1]. The regions occupied by the permeable substrates correspond to y∈ [−(1+H),−1) and y∈ (1,1+H]. The height of the permeable substrate is H. Note that all lengths are normalized by the channel half-height (h † ). 90 5.2.2 Mo dified Resolv en t Analysis The resolvent formulation for wall-bounded turbulent flows proposed by McKeon and Sharma [73]—and employed in several recent flow control studies [ 21, 68, 69, 81, 116]—is extended to accountforthepresenceofpermeablesubstratesasfollows. Foranextendeddiscussionofresolvent analysis and its applications, the reader is referred to McKeon [72]. ConstructionofthemodifiedresolventoperatorsbeginswithastandardReynolds-decomposition of the simplified volume-averaged Navier-Stokes equations in ( 5.2). The velocity and pressure fields are decomposed into a time-averaged component (denoted by (·)) and a fluctuating compo- nent about this average (denoted by(·) ′ ). Under this decomposition, the velocity field is expressed asu(t,x) =U(x)+u ′ (t,x) and the pressure field is expressed as p(t,x) =P(x)+p ′ (t,x). Note that U(x) = [U(y),0,0] T represents the turbulent mean profile. Next, the velocity and pressure fluctuation are Fourier-transformed in the homogeneous streamwise and spanwise directions as well as in time: " u ′ (t,x) p ′ (t,x) # = y " u k (y) p k (y) # exp(−iωt+iκ x x+iκ z z)dωdκ x dκ z . (5.3) In the expression above,κ x is the streamwise wavenumber,κ z is the spanwise wavenumber, andω isthefrequency. TheFouriercoefficientsforthevelocityandpressurefieldatagivenwavenumber- frequency combination,k = (κ x ,κ z ,ω), are denotedu k andp k . Under the Fourier transform, the equations in (5.2) can be expressed compactly as: " u k p k # = −iω " I 0 # − " L k − ˜ ∇ − ˜ ∇ T 0 #! −1 " I 0 # f k =H k f k . (5.4) 91 In the expression above, the first row represents the momentum equations and the second row represents the continuity constraint. The operator L k = 2 6 6 6 6 6 6 4 −iκ x U +Re −1 τ ( ˜ ∇ 2 −K −1 ) − dU dy 0 0 −iκ x U +Re −1 τ ( ˜ ∇ 2 −K −1 ) 0 0 0 −iκ x U +Re −1 τ ( ˜ ∇ 2 −K −1 ) 3 7 7 7 7 7 7 5 (5.5) represents the linear dynamics in (5.2), andf k is the Fourier coefficient for the non-linear terms. The differential operator ˜ ∇ is defined as ˜ ∇ = (iκ x ,∂/∂y,iκ z ), and so the symbols ˜ ∇ T and ˜ ∇ essentially represent the Fourier-transformed divergence and gradient operators, respectively. Similarly, the Laplacian is defined as ˜ ∇ 2 = (−κ 2 x −κ 2 z + ∂ 2 ∂y 2 ). Note that the velocity and pressure response at a given wavenumber-frequency combination constitute a traveling wave flow field with streamwise wavelength λ x = 2π/κ x and spanwise wavelength λ z = 2π/κ z that is moving downstream at speed c = ω/κ x . The transfer function that maps the nonlinear forcing f k to the velocity and pressure response [u k ,p k ] T in (5.4) is the resolvent operator,H k . The central difference between this resolvent operator for channel flow over permeable substrates and its smooth wall counterpart is the inclusion of the Darcy permeability tensor K in (5.5). Note that the resistance exerted by the Darcy term is linear with respect tou. As detailed in prior studies [68, 69, 73, 78], a singular value decomposition (SVD) of the discretized resolvent operator (see Section 5.2.5) is used to identify a set of orthonormal forcing and response modes that are ordered based on their forcing-response gain under an L 2 energy norm. To enforce this norm, the discretized resolvent operator in (5.4) is scaled as follows: h W u 0 i " u k p k # = h W u 0 i H k W −1 f W f f k (5.6a) 92 or W u u k =H w k W f f k . (5.6b) Here, W u and W u incorporate numerical quadrature weights for the entire domain spanning y∈ [−(1+H),(1+H)]. With this weighting, the SVD of the scaled resolvent: H w k = X m ψ k,m σ k,m ϕ ∗ k,m , (5.7a) where σ k,1 >σ k,2 >...> 0, ψ ∗ k,l ψ k,m =δ lm , ϕ ∗ k,l ϕ k,m =δ lm (5.7b) yields forcing modes f k,m = W −1 f ϕ k,m and velocity responses u k,m = W −1 u ψ k,m that have unit energy when integrated across the entire domain spanning y∈ [−(1 +H),(1 +H)]. In other words, this scaling ensures that Z (1+H) −(1+H) u ∗ k,l u k,m dy =δ lm , (5.8a) and Z (1+H) −(1+H) f ∗ k,l f k,m dy =δ lm . (5.8b) In (5.7)-(5.8), the superscript (·) ∗ denotes a conjugate transpose. A major contribution of the resolvent framework lies in the finding that the forcing-response transferfunctiontendstolow-rankatk combinationsthatareenergeticinwall-boundedturbulent flows. Often, the first singular value is an order of magnitude larger than subsequent singular 93 values, σ k,1 ≫σ k,2 >..., and so the resolvent operator can be well approximated using a rank-1 truncation after the SVD [73, 78]: H w k ≈ψ k,1 σ k,1 ϕ ∗ k,1 . (5.9) The expressions in (5.6)-(5.9) show that forcing in the direction of the first forcing mode f k,1 = W −1 f ϕ k,1 generates a velocity response u k,1 = W −1 u ψ k,1 that is amplified by factor σ k,1 . Put another way, a forcing of the form f k,1 to the unscaled resolvent operator in (5.4) generates a velocity and pressure responseσ k,1 [u k,1 ,p k,1 ] T . Under theL 2 scaling used here,σ 2 k,1 is a measure of energy amplification. Recent modeling efforts for active and passive flow control techniques show that specific rank- 1 modes serve as useful surrogates for the dynamically-important NW cycle. Specifically, the ability of a control technique to suppress the forcing-response gain for modes with wavenumber- frequency combinations corresponding to λ + x ≈ 10 3 , λ + z ≈ 10 2 , and c + ≈ 10 (i.e., similar to the length and velocity scales associated with NW streaks) has been shown to be a useful predictor of drag reduction performance [21, 68, 81]. Building on these prior efforts, in this study we evaluate the effect of anisotropic permeable substrates on the rank-1 resolvent mode that serves as a surrogate for the NW cycle. A reduction in gain for this mode relative to the smooth-wall value is interpreted as mode suppression, which is indicative of drag reduction. In addition, we also test whether the permeable substrates lead to an increase in principal singular values for spanwise-coherent modes (e.g., with κ z = 0) that resemble KH rollers. Since we only consider the rank-1 truncation shown in (5.9) for the remainder of this chapter, we drop the additional subscript 1 to simplify notation. 94 5.2.3 Boundary and In terface Conditions As discussed in Section 5.2.5 below, the resolvent operator is discretized using spectral dis- cretization and rectangular block matrices as described by Aurentz and Trefethen [4]. This ap- proach enables us to use two different sets of equations in the unobstructed region and porous domain (i.e., without and with the permeability term), and couple the two via appropriate inter- facial conditions. As shown in figure 5.1, the unobstructed channel corresponds to the region cor- responding toy∈ [−1,1] and the upper and lower permeable regions correspond toy∈ (1,1+H] and y∈ [−(1+H),−1), respectively. At the lower and upper substrate walls, y =±(1+H), we apply no-slip boundary conditions. At the interfaces between the porous medium and the unob- structed flow, y =±1, we impose continuity in all three components of velocity and pressure. We also impose continuity in the streamwise and spanwise shear at the interface. These boundary conditions can be summarized as follows: u = 0 at y =±(1+H), (5.10a) u| y + = u| y − and p| y + = p| y − at y =±1, (5.10b) ∂u ∂y y + = 1 ε ∂u ∂y y − and ∂w ∂y y + = 1 ε ∂w ∂y y − at y =−1, (5.10c) ∂u ∂y y − = 1 ε ∂u ∂y y + and ∂w ∂y y − = 1 ε ∂w ∂y y + at y = 1. (5.10d) In the expressions above, y + and y − denote values taken on either side of the porous substrate- unobstructed flow interface. Note that the boundary conditions shownin ( 5.10c,d) do not penalize momentum transfer into the permeable substrate. These conditions also assume that the effective viscosity in the porous medium is ˜ ν =ν/ϵ. Per the discussion presented in Abderrahaman-Elena 95 and García-Mayoral [1] and Gómez-de-Segura and García-Mayoral [40], these assumptions are reasonable for a highly-connected medium with high porosity, ϵ≈ 1. For a poorly-connected medium, previous studies show that a stress jump boundary condition may be more appropriate, and that more complex effective viscosity models may be needed [see e.g., 75, 83]. Also keep in mind that the equations above imply a sharp transition between the porous medium and the unobstructed fluid. Previous studies employing the volume-averaged Navier-Stokes equations have typically assumed the existence of a finite transition zone between the homogeneous porous and homogeneous fluid regions [ 15, 83]. Finally, since we are interested in permeable substrates that have the highest potential for drag reduction, we assume a porosity of ϵ = 1. For this value of porosity, the boundary conditions shown in (5.10) are similar to those used in previous high-fidelity simulations [ 40]. Emulation of the channel configuration and boundary conditions used by Gómez-de-Segura and García-Mayoral [40] allows for a direct comparison between model predictions and simulation results. 5.2.4 Mean V elo cit y Profile As shown in (5.5), construction of the resolvent operator requires knowledge of the mean velocity profile U(y). Here, we use two different sets of mean velocity profiles to generate model predictions: (i) mean profiles obtained in DNS by Gómez-de-Segura and García-Mayoral [ 40], and (ii) mean profiles generated using a synthetic eddy viscosity profile. The DNS mean profile datasetconsists of 22 different cases: 8 profiles for ϕ xy = 3.6; 7 profiles for ϕ xy = 5.5; and 7 profiles for ϕ xy = 11.4. For all 22 cases tested in DNS, the ratio between the wall-normal and spanwise permeability length scales was ϕ yz = p K + y / p K + z = 1. Each of the DNS profiles contained 153 points spanning the unobstructed region, y∈ [−1,1]. Within the permeable medium, the analytical solutions developed by Gómez-de-Segura and García-Mayoral [40] were used. These 96 mean profiles were interpolated onto our Chebyshev grid using the modified Akima interpolation algorithm. For each of the 22 different DNS cases, the synthetic profiles were generated as follows. The modified governing equation ( 5.2) was Reynolds-averaged to yield the following equation for the streamwise mean flow: 1 Re τ d 2 U dy 2 − U Re τ K x − d(u ′ v ′ ) dy = dp dx , (5.11) where K x is the streamwise component of the permeability tensor. The Reynolds shear stress term in (5.11) was estimated using a modified version of the eddy viscosity formulation proposed by Reynolds and Tiederman [97]: ν e = 8 > > > > > > < > > > > > > : 1 2 1+ c 2 Re τ 3 (2y−y 2 )(3−4y+2y 2 ) × 1−exp (|y−1|−1) Re τ c 1 2 # 1/2 − 1 2 |y|≤ 1, 0 |y|> 1. (5.12) Here, theeddyviscosityhasbeennormalizedbythefluidviscosity ν. Togeneratethesynthetic profiles, we used the values c 1 = 46.2 and c 2 = 0.61, which were identified by Moarref and Jovanović [77] as yielding the best fit to the mean profiles obtained in smooth wall DNS at Re τ ≈ 180. Thus, the Reynolds shear stress term was modeled using a standard smooth-wall eddy viscosity profile in the unobstructed region of the channel and assumed to be zero in the permeable substrate. In other words, this model assumes that there is no turbulence penetration into the porous medium. In reality, the penetration of turbulent cross-flows into the permeable substrate will depend on the spanwise and wall-normal components of permeability, and so the 97 model above is only valid for smallK y andK z . Finally, equations (5.11) and (5.12) were combined to yield the following equation for the mean velocity profile: 1 Re τ (1+ν e ) d 2 dy 2 + dν e dy d dy − 1 K x U = dp dx , (5.13) which was solved numerically to yield U(y). We recognize that the procedure outlined above to estimate the synthetic mean profile involves significant assumptions. However, as we show below, the interfacial slip velocities generated using this model are in good agreement with the DNS mean profiles in the initial drag reduction regime over anisotropic permeable substrates, i.e., until the point of performance degradation. Moreover, the resolvent-based predictions obtained using these profiles are also similar to those computed using the DNS mean profiles. 5.2.5 Numerical Metho ds The resolvent operator and the equations used to synthesize the mean velocity profile are discretized in the wall-normal direction using spectral discretization methods involving rectan- gular block operators, as described by Aurentz and Trefethen [4]. Each differential operator is discretized using Chebyshev polynomials and the resulting matrices are rectangular. The size of these matrices is [N×N +n] and where n represents the dimension of the operator null space. The overall block operator is made square by appending n boundary conditions. By using these block operators, we can deal with a system of boundary value problems coupled through boundary conditions. For a more thorough discussion of these methods readers are directed to Aurentz and Trefethen [4]. As shown in figure 5.1, for our problem configuration the channel is separated into 3 regions: the lower permeable substrate spanningy∈ [−H,−1), the free channel spanningy∈ [−1,1], and 98 the upper permeable substrate spanning y∈ (1,H]. Although the configuration is symmetric across the channel centerline, we consider the entire channel to retain resolvent modes that are both symmetric and anti-symmetric across the channel centerline [69, 78]. Each of the three regions is discretized usingN Chebyshev points and so the total number of Chebyshev points for the channel is 3N. A grid convergence study showed that the normalized change in singular values is ofO(10 −6 ) for grid sizes N = 80, N = 112, and N = 200. Each resolvent mode computation takes approximately 0.7s forN = 80, approximately 2.0s forN = 112 and 7.0s forN = 200. The results presented below were generated usingN = 112, which corresponds to a total of 3N = 336 grid points across the the entire channel. 5.3 Results and Discussion In this section, we compare resolvent-based predictions for the NW mode (Section 5.3.2) and the emergence of energetic spanwise rollers (Section 5.3.4) against DNS observations. Before that, we briefly compare the mean velocity profiles obtained from DNS against the synthetic profiles generated using the eddy viscosity model (Section 5.3.1). 5.3.1 Mean V elo cit y Profiles Figure5.2comparestheinterfacialslipvelocityU + s predictedusingthesyntheticeddyviscosity profile( 5.12)againstthoseobtainedinDNSbyGómez-de-SeguraandGarcía-Mayoral[40]foreach of the 22 different cases being considered here, with anisotropy ratios ϕ xy = 3.6, 5.5, and 11.4. The slip velocities for the synthetic profiles all collapse together onto a straight line corresponding to U + s ≈ p K + x . In other words, the slip velocity only depends on the streamwise component of permeability, which is consistent with the assumptions outlined in Section 5.2.4. The synthetic 99 Figure 5.2: Predicted slip velocity at the porous interface as a function of streamwise permeability p K + x . Closed symbols show DNS results from Gómez-de-Segura and García-Mayoral [40]. Open symbolsshowpredictionsmadeusing(5.13). Redsymbolsindicatethecasesforwhichwall-normal permeability first exceeds p K + y & 0.4. This is roughly the threshold above which drag reduction performance deteriorates in DNS. mean profiles do not account for the effect of wall-normal or spanwise permeability, which are likely to determine the extent to which turbulence penetrates into the permeable substrate. The synthetic slip velocities (open symbols) show close agreement with DNS results (closed symbols) for low streamwise permeabilities but begin to deviate from DNS results at higher p K + x . Moreover, the threshold values of p K + x above which the synthetic predictions begin to deviate from DNS results depend on the anisotropy ratio, ϕ xy . As an example, for ϕ xy = 3.6, synthetic U + s begin deviating from DNS results for p K + x & 2. For ϕ xy = 11.4, the predicted slip velocities deviate significantly from DNS results for p K + x & 6. Closer inspection of this trend shows that the synthetic slip velocities deviate from DNS data above a constant threshold value for the wall- normalpermeabilityforallanisotropyratios, p K + y & 0.4. Thisvalueof p K + y correspondsclosely to the conditions in which Kelvin-Helmholtz-type rollers appear over the permeable substrates in the simulations of Gómez-de-Segura and García-Mayoral [40]. The appearance of these rollers is likely to generate significant interfacial turbulence that penetrates the porous medium. The eddy 100 viscosity model used to generate the synthetic mean profiles does not account for these effects. The eddy viscosity model shown in (5.12) assumes that no turbulence penetrates into the porous medium, and that the turbulence in the unobstructed region remains similar to that over a smooth wall regardless of the value of p K + y and p K + z . For cases in which the interfacial slip velocities agree, the synthetic mean profiles are similar to the DNS profiles across the entire channel (data not shown here for brevity). For these cases, resolvent-based predictions are not very sensitive to the choice of mean profile. Moreover, as we show in Section 5.3.4, resolvent analysis is able to predict the emergence of high-gain KH rollers over permeable substrates as the wall-normal permeability increases beyond p K + y & 0.4. These observations indicate that resolvent analysis carried out using the synthetic mean profile can be a useful design tool for the design of passive flow control using anisotropic permeable substrates. 5.3.2 Near-W all Resolv en t Mo de Previous studies have shown that the resolvent mode with streamwise wavelength λ + x = 10 3 , spanwise wavelength λ + z = 10 2 , and phase speed c + = (ω/κ x ) + = 10 serves as a useful surro- gate for the dynamically-important NW cycle characterized by the presence of quasi-streamwise vortices and alternating high- and low-speed streaks [73, 106]. The forcing-response gain (i.e., singular value) for this resolvent mode has also proven to be a useful predictor of control per- formance for both active [68, 81, 116] and passive [21, 69, 70] techniques. Here, we evaluate whether it can serve as a useful reduced-complexity tool for the evaluation of anisotropic perme- able substrates for passive drag reduction. Since the resolvent-based predictions are evaluated against DNS results from Gómez-de-Segura and García-Mayoral [40] carried out at Re τ = 180, the resolvent mode with streamwise wavenumber κ x = 2πRe τ /10 3 ≈ 1.1, spanwise wavenumber κ z = 2πRe τ /10 2 ≈ 11 and frequency ω = 10κ x ≈ 11 is used as a surrogate for the NW cycle. 101 Figure 5.3: Predictions for the NW resolvent mode. (a) Comparison between normalized mode gain (σ k,p /σ k,s , black symbols) and drag reduction observed in DNS (∆U + , light gray symbols) by Gómez-de-Segura and García-Mayoral [40], plotted as a function of p K + y . The resolvent-based predictions shown in this panel are obtained using the velocity profiles from DNS. (b) Comparison betweennormalizedgainspredictedusingtheDNSmeanprofiles(filledsymbols)andthesynthetic mean profiles (open symbols) plotted as a function of p K + y . (c) Comparison between normalized gains predicted using the DNS mean profiles (filled symbols) and the synthetic mean profiles (open symbols) plotted as a function of p K + x − p K + y . The black and red dashed lines show linear fits to the initial decrease in normalized gain for the DNS- and synthetic-mean predictions, respectively. For all panels, the symbols represent substrates with ϕ xy = 3.6, symbols represent substrates with ϕ xy = 5.5, and symbols represent substrates with ϕ xy = 11.4. 102 Figure 5.3(a) shows the forcing-gain for the NW resolvent mode over the porous substrateσ k,p normalized by the smooth-wall value σ k,s at Re τ = 180 for all 22 substrates tested by Gómez- de-Segura and García-Mayoral [40]. Resolvent-based gain predictions (filled black symbols) are compared against the drag reduction measured in DNS (filled gray symbols). Following Gómez- de-Segura and García-Mayoral [40], the outward shift in the logarithmic region of the mean profile is used to quantify the change in drag: ∆U + > 0 denotes a decrease in drag and ∆U + < 0 denotes an increase in drag. For all three anisotropy ratios, there is qualitative agreement between NW modesuppressionanddragreduction. Foragivenvalueof p K + y , materialswithhigheranisotropy ratios are predicted to yield greater suppression of the NW resolvent mode, which is consistent with the drag reduction trends observed in DNS. Further, mode suppression and drag reduction both increase monotonically up to a value of p K + y ≈ 0.4. Agreement between drag reduction and mode suppression is particularly good for substrates with anisotropy ratiosϕ xy = 3.6 and 5.5. For p K + y & 0.4, the normalized gain begins to increase again for the substrates with ϕ xy = 3.6 and 5.5, which is consistent with DNS drag reduction trend. DNS results show an increase in drag relative to smooth-wall values for p K + y & 0.55 for the substrates with ϕ xy = 3.6 and 5.5. Resolvent-based predictions show that the gain of the NW mode over the permeable substrate exceeds the smooth-wall value beyond p K + y & 0.7 for these substrates. For substrates with the higher anisotropy ratio, ϕ xy = 11.4, there are visible discrepancies between the predicted change in NW mode gain and the measured drag reduction, particularly for higher values of the wall normal permeability, p K + y . For instance, the increase in normalized gain beyond p K + y & 0.4 is less pronounced for the material with ϕ xy = 11.4. For this substrate, DNS results show an increase in drag for p K + y & 0.7. However, the resolvent predictions do not show an increase in gain beyond the smooth-wall value for the conditions tested (up to p K + y ≈ 1). 103 Note that the observed disagreement between NW mode gain and drag reduction for p K + y & 0.4 is not surprising. As discussed in prior studies [1, 40, 41], the initial decrease in drag over anisotropic permeable substrates depends on the shift in the virtual origins for the streamwise mean flow and the turbulent cross-flow fluctuations produced by the quasi-streamwise vortices associated with the NW cycle. As we show below, this effect is reproduced by the NW resolvent mode. However, the increase in drag beyond p K + y ≈ 0.4 is associated with the emergence of spanwiserollersresemblingKHvortices. Thisphenomenoncannotbecapturedbyresolventmodes that serve as surrogates for the streamwise-elongated structures constituting the NW cycle. We use the resolvent framework to test for the emergence of energetic spanwise rollers in Section 5.3.4. Figure 5.3(b) shows normalized gain predictions made using the DNS mean profile (filled symbols) and the synthetic mean profile (open symbols) for the NW resolvent mode, plotted as a function of p K + y . For values of p K + y . 0.4, the normalized gain predictions are in close agree- ment. For this range of wall-normal permeabilities, the normalized gain predictions made using the synthetic mean profile are within 6% of those made using the DNS mean profile. However, for p K + y & 0.4, the synthetic profile predictions show a continued decrease in NW mode gain, while the DNS profile predictions show an increase in gain. Discrepancy between the predictions for higher values of p K + y is consistent with the slip velocities shown in figure 5.2. The synthetic profile agrees with the DNS mean profile for p K + y . 0.4, but not beyond this value. As discussed earlier, this is because higher wall-normal permeabilities give rise to spanwise rollers resembling KH vortices in the simulations. The effect of the additional interfacial turbulence generated by these rollers is not captured in the eddy viscosity model (5.12). Instead, the eddy viscosity model assumes smooth-wall like turbulence regardless of the value of p K + y and p K + z . Figure 5.3(c) shows normalized gain predictions made using the DNS mean profile (filled symbols) and the synthetic mean profile (open symbols) for the NW resolvent mode, plotted 104 as a function of p K + x − p K + y . Since the permeable substrate model used here has identical wall-normal and spanwise permeabilities, the difference between streamwise and wall-normal per- meability is also equal to the difference between streamwise and spanwise permeabilities, i.e., p K + x − p K + y = p K + x − p K + z . Slip length arguments [1, 67] suggest that the initial decrease in drag over anisotropic permeable substrates is proportional to the difference between the stream- wise and spanwise permeability length scales, which determine the virtual origins felt by the mean streamwise flow and turbulent cross-flow fluctuations, respectively. In other words, we expect ∆U + ∝ p K + x − p K + z . These arguments are supported by the simulation results of Gómez-de-Segura and García-Mayoral [40], which show a linear decrease in drag that is propor- tional to p K + x − p K + z for substrates with small permeabilities. The NW mode gains shown in figure 5.3(c) agree well with this model. Specifically, the normalized NW-mode gains initially decrease linearly with p K + x − p K + y (= p K + x − p K + z ) across all anisotropy ratios. Moreover, there is little difference between the predictions made using the synthetic and DNS mean pro- files. Linear relationships fitted to the initial decrease in normalized gain are nearly identical for the predictions based on the DNS profiles (dashed black lines) and those based on the syn- thetic mean profiles (dashed red lines). Note that the predictions based on the synthetic mean profile all collapse together. This suggests that NW mode behavior is self-similar when appropri- ately normalized. Predictions made using the DNS profiles diverge from this self similarity when p K + y & 0.4. Figure 5.4 provides physical insight into the observed changes in mode gain. The structure of the NW resolvent mode for the smooth-wall case is shown in figure 5.4(a). By construction, the mode structure shows alternating regions of positive and negative velocity fluctuations with span- wise wavelengthλ + z = 10 2 . The wall-normal and spanwise velocity fluctuations show the presence of counter-rotating streamwise vortices. Regions of downwelling (v ′ towards wall) coincide with 105 Figure 5.4: NW mode structure predictions made using the DNS velocity profiles for (a,b) the smooth-wallcase, (c,d)thesubstratewithϕ xy = 5.5and p K + y = 0.18,and(e,f)thesubstratewith ϕ xy = 5.5 and p K + y = 0.45, and (g,h) the substrate withϕ xy = 5.5 and p K + y = 1.0. In (a,c,e,f), the shading shows normalized contours of positive (red) and negative (blue) streamwise velocity normalized by the maximum value. The vectors show the wall-normal and spanwise velocity fluctuations. In (b,d,f,g), the solid lines show the magnitude of the streamwise velocity for this resolvent mode,|u k | and the dashed lines show the wall-normal velocity magnitude multiplied by a factor of ten, 10|v k |. The black lines represent the smooth-wall case while the gray lines represent the permeable substrates. These plots make use of the shifted coordinate ˆ y = 1 +y, such that ˆ y = 0 represents the location of the smooth wall or the porous interface. 106 regions of high streamwise velocity (red shading) and regions of upwelling (v ′ away from wall) coincide with regions of low streamwise velocity (blue shading). These observations are consistent with known features of the NW cycle [48, 98, 120]. When the solid wall is replaced by anisotropic permeable substrates withϕ xy = 5.5, there are some subtle but important changes in mode struc- ture. For the permeable substrates that yield mode suppression—and drag reduction—there is a downward shift in the velocity fluctuations and the mode is c ompr esse d in the wall-normal direction (see Figs. 5.4(c,d) for p K + y = 0.18 and Figs. 5.4(e,f) for p K + y = 0.45). To explain why these changes lead to mode suppression, we briefly review two mechanisms responsible for high gain in the resolvent framework [72, 73]. First, resolvent modes tend to localize around the critical layer, y c , where the mode speed matches the local mean velocity, U(y c ) = c = ω/κ x . This minimizes the term−iω+iκ x U inside the inverted operator shown on the right-hand side of (5.4). Second, energy is extracted from the mean flow via the interaction between the wall-normal velocity fluctuations and the mean shear, i.e., the so-called lift-up mechanism that leads to energy transfer from the mean flow to the fluctuations via the term −u ′ v ′ (dU/dy) in the turbulent kinetic energy equation. As the streamwise permeability of the substrate increases, there is a shift in the virtual origin of the mean profile and so the critical layer location moves closer to the permeable interface, from ˆ y + c ≈ 13 for the smooth-wall case to ˆ y + c ≈ 9.7 for the case with p K + y = 0.45. Here, ˆ y = 1+y such that ˆ y = 0 represents the porous interface at y =−1 (see figure 5.1). As noted in the introduction, this shift in the mean profile depends on p K + x . A downward shift in the critical layer location causes the NW resolvent mode to localize closer to the porous interface (see→−← in figure 5.4). However, the downward shift for the NW mode is constrained by the low wall-normal and spanwise permeabilities, p K + y = p K + z , which leads to the observed wall-normal compression in mode structure. This wall-normal compression in mode structure reduces the energy extracted from the mean flow, which depends on the wall-normal integral of 107 the product of the Reynolds shear stress and mean shear, i.e.,∝ R Re(−u ∗ k v k )(dU/dy)dy. Note that Re(·) represents the real component and so the Reynolds shear stress contribution from the NW resolvent mode is proportional to Re(−u ∗ k v k ) [68]. A reduction in energy extraction from the mean flow leads to a reduction in mode gain. Importantly, this observation is consistent with the interpretation that drag reduction is only expected if the slip length for the mean streamwise flow ℓ + x ∝ p K + x is larger than that for the turbulent cross-flows, ℓ + z ∝ p K + z , i.e., if the virtual origin for the mean flow is offset from the virtual origin for the turbulence. Additional resolvent-based predictions for p K + z > p K + x (ℓ + z >ℓ + x ) lead to an increase in mode gain (data not shown). Figures 5.4(g,h) show predictions for NW mode structure for the anisotropic permeable sub- strate with ϕ xy = 5.5 and p K + y = 1.0 that leads to an increase in drag and forcing-response gain (see figure 5.3). Compared to the cases that lead to mode suppression and drag reduction, the wall-normal footprint for the NW mode increases substantially for this case. This can be attributed to the change in the DNS mean profiles for p K + y & 0.4, for which spanwise-coherent structures resembling KH vortices emerge over the permeable substrate. The emergence of these rollers leads to a substantial increase in turbulent mixing which causes the mean profile to flatten in the interfacial region [40]. For the substrate with ϕ xy = 5.5 and p K + y = 1.0, this pushes the critical layer for the NW mode outwards to ˆ y + c ≈ 26. This means that the NW mode is again able to localize around the critical layer, and there is no vertical compression in mode structure. This allows for greater energy extraction from the mean flow and leads to the observed increase in NW mode gain over the material with ϕ xy = 5.5 and p K + y = 1.0 in figure 5.3(a). Note that the NW mode gain decreases monotonically with increasing p K + x − p K + y for predictionsmadeusingthesyntheticmeanprofiles(seefigure 5.3(c)). Thisisbecausethesynthetic mean profiles assume no change in turbulence characteristics near the interface relative to smooth- wall conditions, i.e., they do not account for the emergence of spanwise rollers and turbulence 108 penetration into the porous medium. As shown in figure 5.3(a), NW mode gain begins to increase again for p K + y & 0.4 when the DNS mean profiles are used in the resolvent operator. In the DNS, spanwise rollers resembling KH vortices emerge beyond this value of wall-normal permeability and alter the form of the mean profile. The preceding discussion assumes that the phase speed of the near-wall structures (c + ≈ 10) remains unchanged relative to smooth wall conditions. It could also be argued that the phase speed changes over the permeable medium to reflect the interfacial slip velocity, which depends on the streamwise slip length (U + s ≈ℓ + x for smallℓ + x ), and turbulence penetration into the permeable medium, which depends on the spanwise slip length (ℓ + z ). To account for these effects, the phase speed can be modified to c + = 10+ℓ + x −ℓ + z ≈ 10+ p K + x − p K + z . Figure 5.5 in Section 5.3.3 shows additional predictions for NW mode gain made using this modified phase speed. These additional predictions are broadly similar to the predictions shown in figure 5.3 for c + = 10. Specific observations are discussed further in Section 5.3.3. Together, the predictions shown in figure 5.3 and figure 5.4 confirm that suppression of the dynamically-important NW cycle requires high streamwise permeability and low spanwise/wall- normal permeabilities. In other words, materials with high anisotropy ratios (ϕ xy ≫ 1) are expected to yield the largest suppression of NW turbulence. However, the absolute value of the wall-normal permeability also plays an important role in dictating drag reduction performance. Specifically, the emergence of energetic spanwise coherent rollers is dictated primarily by p K + y . 5.3.3 Phase-Corrected Near-W all Resolv en t Mo de As noted earlier, the predictions shown in Section 5.3.2 assume that the phase speed of the resolvent mode that serves as a surrogate for the NW cycle (c + ≈ 10) remains unchanged over the permeable substrate. Figure 5.5 shows additional predictions for mode gain assuming the 109 phase speed change with substrate permeability as c + = 10+ℓ + x −ℓ + z ≈ 10+ p K + x − p K + z . In effect, this modification accounts for the interfacial slip velocity, which depends on the streamwise slip length ℓ + x ≈ p K + x , and the virtual origin for the turbulence inside the permeable substrate, which depends on the spanwise slip length ℓ + z ≈ p K + z . The modified phase speed does not lead to significant changes in mode gain compared to the predictions shown in figure 5.3 for c + ≈ 10. Figure 5.5(a) shows that mode gain predictions made using the DNS mean profile are consistent with drag reduction trends obtained in DNS. Compared to the predictions shown in figure 5.3(a), singular value ratios for the substrates with ϕ xy = 3.6 andϕ xy = 5.5 show slightly better agreement with the drag reduction trends from DNS for p K + y < 0.6. However, for p K + y & 0.6, the predictions are more scattered with the modified phase speed. For the substrate with ϕ xy = 11.4, the modified phase speed leads to reduced mode suppression. Further, the deterioration in performance—in terms of mode suppression— commences at a lower value of wall-normal permeability, p K + y ≈ 0.2, compared to the threshold observed in DNS. Figure 5.5(b) shows that the predictions made using the synthetic mean profile are consistent with those made using the DNS mean profiles for small values of p K + y , i.e., until the emergence of spanwise rollers resembling Kelvin-Helmholtz vortices in the DNS. Figure 5.5(c) shows that the initial reduction in mode gain depends primarily on p K + x − p K + z , even with the modified phase speed. These observations are again consistent with the predictions shown in figure 5.3(b,c). However, linear fits to the initial decrease in mode gain show that the modified phase speed leads to less pronounced mode suppression; see red and black dashed lines in figures 5.3(c) and 5.5(c). 110 Figure 5.5: Predictions for the NW resolvent mode with a modified phase speed of c + ≈ 10 + p K + x − p K + z . (a) Comparison between normalized mode gain (σ k,p /σ k,s , black symbols) and drag reduction observed in DNS (∆U + , light gray symbols) by Gómez-de-Segura and García- Mayoral [40], plotted as a function of p K + y . The resolvent-based predictions shown in this panel are obtained using the velocity profiles from DNS. (b) Comparison between normalized gains predicted using the DNS mean profiles (filled symbols) and the synthetic mean profiles (open symbols) plotted as a function of p K + y . (c) Comparison between normalized gains predicted using the DNS mean profiles (filled symbols) and the synthetic mean profiles (open symbols) plotted as a function of p K + x − p K + y . The black and red dashed lines show linear fits to the initial decrease in normalized gain for the DNS- and synthetic-mean predictions, respectively. For all panels, the symbols represent substrates with ϕ xy = 3.6, symbols represent substrates with ϕ xy = 5.5, and symbols represent substrates with ϕ xy = 11.4. 111 5.3.4 Span wise-Coheren t Resolv en t Mo des Figure 5.6: Normalized gain for spanwise-coherent structures plotted as a function of streamwise lengthλ + x and wave speedc + for substrates with an anisotropyϕ xy = 5.5. These predictions make use of the synthetic mean profile. Panels (a)-(d) represent structures with spanwise wavenumber κ z ≈ 2.3 (λ + z = 500). Panels (e)-(h) represents structures with spanwise wavenumber κ z = 0. Red and blue shading respectively represent mode amplification and suppression relative to the smooth-wall case. Wall-normal permeability increases from left to right. Structures with the largest amplification are labeled with a ( •) marker. In this section, we use the resolvent analysis framework to evaluate the conditions in which energetic spanwise-coherent rollers resembling Kelvin-Helmholtz vortices emerge over anisotropic permeable substrates. Snapshots of the flow field from the simulations of Gómez-de-Segura and García-Mayoral [40] show that for p K + y < 0.4, the flow is dominated by streamwise-elongated structuresassociatedwiththeNWcycle. Asthewall-normalpermeabilityincreasesto p K + y & 0.4 the flow becomes dominated by spanwise coherent rollers with an approximate streamwise spacing of λ + x ≈ 100−300. Further, using momentum balance arguments, Gómez-de-Segura and García- Mayoral [40] showed that the Reynolds stress contribution from structures with λ + x ≈ 70−320 andλ + z ≥ 120 is responsible for the increase in drag observed for substrates with p K + y & 0.4. We 112 note that spanwise-coherent rollers have also been observed in previous simulations over isotropic porous materials [15, 100]. Linear stability analyses successfully predict the emergence of Kelvin-Helmholtz vortices in channel flows over anisotropic permeable substrates [ 1, 41]. Gómez-de-Segura, Sharma, and García-Mayoral [41] also show that the wall-normal permeability becomes the driving parameter for the onset of KH rollers for streamwise-preferential substrates of large depth, i.e., substrates with p K + x ≫ p K + y and h + ≫ p K + y . However, models based on linear stability theory do not accurately predict the threshold value of p K + y for the onset of the rollers or the wavelength of the rollers. Linear stability analyses predict structures with streamwise wavelengthλ + x ≈ 50−70 to be most unstable, while DNS results suggest that structures withλ + x ≈ 150 are most energetic initially. Figure 5.6 shows resolvent-based predictions for the normalized amplification of spanwise- coherent structures over the anisotropic permeable substrate with ϕ xy = 5.5 and varying p K + y . Based on the DNS results of Gómez-de-Segura and García-Mayoral [40], we limit ourselves to the evaluation of structures with streamwise wavelength λ + x ∈ [50,500] and mode speed c + ∈ [4,10]. In addition, we consider structures that are infinitely long in the spanwise direction ( κ z = 0) as well as structures with a finite, but large, spanwise wavelength ( κ z ≈ 2.3; λ + z = 500). Figures 5.6(a,e) show that for low wall-normal permeability, p K + y = 0.18, resolvent analysis predicts minimal amplification of spanwise-coherent structures over the anisotropic permeable substrate. The maximum increase in gain relative to smooth-wall values is less than 25%, i.e., σ k,p /σ k,s < 1.25 for both κ z = 0 and κ z ≈ 2.3. Moreover, the gain for these spanwise-coherent 113 structures is also small in absolute terms. This is consistent with the simulation results of Gómez- de-Segura and García-Mayoral [40], which show that the Reynolds stress contributions from struc- tures with λ + x ≈ 70−300 and λ + z & 120 are negligible for materials with p K + y . 0.31. In other words, spanwise-coherent modes do not have a significant effect on drag for these conditions. However, the amplification of these spanwise-coherent structures increases significantly for p K + y ≥ 0.39. Figures 5.6(b,f) show a clear increase in the size and intensity of the region showing mode amplification. For λ + z = 500, the region of highly-amplified structures is localized around (λ + x ,c + )≈ (130,7) and the normalized gain for the most amplified structure is σ k,p /σ k,s ≈ 2. For spanwise constant structures, the region of high amplification is localized around (λ + x ,c + )≈ (180,6.7) and extends to structures with λ + x ≈ 300. The normalized gain for structures in this region is as high as σ k,p /σ k,s ≈ 3, indicating a 200% increase in gain relative to smooth-wall values. The extent and intensity of these high-gain regions continue to increase as wall-normal permeability increases further, as shown in Figs. 5.6(c,d) and (g,h). Indeed, for the case with p K + y = 0.66 shown in figure 5.6(d), the region of high amplification extends from λ + x ≈ 100− 500 and c + ≈ 4− 9. The highest-amplification structures in this region show normalized gains σ k,p /σ k,s ≥ 100. The substantial increase in amplification of spanwise-coherent structures for p K + y ≥ 0.39 is again consistent with the DNS results of Gómez-de-Segura and García-Mayoral [40], which show that the drag-reducing performance of anisotropic substrates begins to degrade for p K + y ≥ 0.39 due to the additional Reynolds shear stresses contributed by spanwise coherent structures with λ + x ≈ 70−300 and λ + z > 120. The DNS results also show the emergence of a new peak in the premultiplied spectra for wall-normal velocity near λ + x ≈ 150 for p K + y = 0.39. This is very close to the streamwise wavelength of the most amplified modes in Figs. 5.6(b,f). Moreover, the streamwise extent of the high-intensity region in the DNS premultiplied spectra increases with 114 Figure 5.7: Flow structure associated with the most amplified spanwise-constant resolvent mode over the permeable substrate withϕ xy = 5.5 and p K + y = 0.39; see• in figure 5.6(f). The red and blue shaded contours show regions of positive and negative pressure, normalized by the maximum value. increasing p K + y , which is similar to the expansion of the high-amplification region evident in Figs. 5.6(c,d,g,h). Figure 5.7 shows a snapshot of the flow field associated with the highest-gain resolvent mode in figure 5.6(f) with λ + x ≈ 180 and c + ≈ 6.7. As expected for Kelvin-Helmholtz type vortices, the structure shows the presence of counter-rotating rollers in the x−y plane, localized near the porous interface. Regions of prograde rotation (i.e., in the direction of the mean shear) are associated with negative pressure fluctuations and regions of retrograde rotation are associated with positive pressure fluctuations. Moreover, the velocity and pressure fields associated with these rollers span a significant portion of the buffer region of the flow, which is consistent with DNS observations. Importantly, the model predictions shown in figure 5.6 for ϕ xy = 5.5 are representative of what is observed over the permeable substrates with ϕ xy = 3.6 and ϕ xy = 11.4 as well. This is illustrated well in figure 5.8, which shows the normalized gain for the most amplifie d resolvent mode in the window spanning λ + z ∈ [50,500] and c + ∈ [4,10], for both κ z ≈ 2.3 and κ z = 0. For all three anisotropy ratios, the normalized gain show a sharp increase for p K + y & 0.4. For structures with λ + z = 500 (κ z ≈ 2.3), model predictions for all three anisotropy ratios collapse 115 Figure 5.8: Comparison of resolvent-based gain predictions for spanwise-coherent structures with anisotropy ratios ϕ xy = 3.6 ( ), 5.5 ( ), and 11.4 ( ). The maximum normalized gain obtained for resolvent modes with λ + x ∈ [50,500] and c + ∈ [4,10] is shown as a function of wall- normal permeability for (a) κ z ≈ 2.3 (λ + z = 500) and for (b) κ z = 0 (λ + z =∞). All of these predictions were obtained using the synthetic mean profiles. together nicely (figure 5.8(a)). For spanwise-constant structures (κ z = 0), there is a bit more scatter in the normalized gain values but the overall trend remains similar. The DNS results in figure 5.3(a) clearly show drag reduction performance deteriorating for p K + y & 0.4 for all three anisotropy ratios. The flow field snapshots, velocity spectra, and momentum-based analyses pursed by Gómez-de-Segura and García-Mayoral [40] attribute this deterioration in performance to the emergence of high-gain spanwise-coherent structures resem- bling Kelvin-Helmholtz vortices. The results presented in this section confirm that resolvent analysis—based on synthetic mean profiles—is able to successfully predict the emergence of such spanwise-coherent structures over anisotropic permeable substrates. Since the synthetic mean profiles are generated used an eddy viscosity model, these results also show that the resolvent framework can generate a priori predictions for whether a porous substrate with known K will give rise to KH rollers. Moreover, the resolvent-based predictions are in better quantitative agree- ment with the DNS results than linear stability predictions. Specifically, resolvent analysis yields 116 better predictions for the wall-normal permeability threshold beyond which such structures be- come important ( p K + y ≈ 0.4) as well as the streamwise wavelength of the structures that emerge first ( λ + x ≈ 150). 5.4 Conclusion Recent theoretical efforts and numerical simulations show that anisotropic permeable sub- strates have the potential to yield significant drag reduction in wall-bounded turbulent flows. Specifically, numerical simulations by Gómez-de-Segura and García-Mayoral [ 40] indicate that a drag reduction of up to 25% can be achieved by appropriately tuning the streamwise, wall-normal, and spanwise permeabilities of such substrates. The mechanism through which these substrates reducedragissimilartothatforriblets, andcanberationalizedintermsoftheoffsetinthevirtual origin felt by the mean flow, which depends on streamwise permeability, and the virtual origin felt by the turbulent fluctuations, which depends primarily on spanwise permeability. Maximum drag reduction is limited by the appearance of spanwise rollers resembling Kelvin-Helmholtz vortices. The emergence of these vortices is linked to a relaxation in the wall-normal permeability. The extended resolvent analysis framework described in this chapter reproduces all of these features with minimal computation. As shown in Section 5.3.2, the gain for a single resolvent mode that serves as a surrogate for the near-wall cycle reproduces the initial drag reduction trends observed in numerical simulations. Model predictions show thatthe reduction in gain for thismode dependson thedifference between the streamwise and spanwise permeability length scales, p K + x − p K + z , such that substrates with higher anisotropy ratios lead to greater mode suppression. This is consistent with the trends observed in numerical simulations. Moreover, the link between resolvent mode gain and drag reduction also provides a complementary explanation to the slip-length based arguments used 117 in previous studies. Specifically, resolvent-based predictions suggest that the difference between p K + x and p K + z leads to a wall-normal compression of the quasi-streamwise vortices associated with the near-wall cycle, and this compression limits energy extraction from the mean flow. Further, as we show in Section 5.3.4, resolvent analysis also predicts the emergence of energetic spanwise rollers with streamwise wavelength λ + x ≈ 150 as the wall-normal permeability increases beyond p K + y ≈ 0.4. These features are quantitatively consistent with simulation results. One weakness of the extended resolvent framework developed here is the requirement of a mean velocity profile over the permeable substrate. For turbulent flows over smooth walls, mean profile estimates can be generated using a variety of simplified models or obtained from previous simulations and experiments. However, such models or data are not readily available for flows over anisotropic permeable substrates. The form of the mean profile can have a significant ef- fect on resolvent-based predictions. Indeed, the comparison between predictions made using the mean profiles obtained in DNS and those generated using an eddy viscosity model presented in Section 5.3.2 shows important differences in the gain for the near-wall mode for p K + y & 0.4, i.e., when drag reduction performance begins to deteriorate. This is because the emergence of the spanwise rollers leads to a substantial increase in interfacial turbulence. This effect is not captured by the eddy viscosity model used to generate the synthetic profiles, which assumes that the turbulence remains smooth-wall like and does not penetrate the permeable substrate (see Section 5.2.4). Yet, resolvent analysis with the synthetic mean profile does show the emergence of these spanwise rollers beyond p K + y ≈ 0.4. In other words, resolvent analysis can generate useful a priori predictions even with a synthetic mean profile. Near-wall mode gain can be used as a measure of initial drag reduction performance, and the emergence of spanwise-coherent rollers can be used to estimate the point at which performance is likely to deteriorate. Since evaluation of these features does not require significant computation, resolvent analysis can serve as a useful 118 design tool for more complex permeable substrates. For instance, resolvent analysis could be used to pursue formal optimization efforts for the full permeability tensor. We also recognize that the volume-averaged Navier-Stokes equations used in this study in- volved significant simplifying assumptions (see Section 5.2.1). These simplifying assumptions were made primarily to allow for a direct comparison with previous simulation results. Moving forward, resolvent analysis could also be used to consider how the inclusion of inertial effects (e.g., via the Forchheimer term) or more complex interfacial boundary conditions (e.g., involving stress jumps, or more complex slip and transpiration models) is likely to affect control performance. 119 Chapter 6 Concluding Remarks In this thesis, we have explored extensions to the resolvent framework to account for riblets and anisotropic porous substrates. For wall-bounded turbulent flows, the resolvent formulation interprets the governing equations as a forcing-response system. A gain-based-decomposition of the Fourier transformed system yields a set of highly amplified velocity and pressure response modes and the forcing that gives rise to these modes. The forcing-response gain determined from this system is a measure of energy amplification and can serve as a measure of control perfor- mance. The reduced-order models developed for riblets and anisotropic permeable substrates using this gain-based metric for performance, to our knowledge, represent the first physics-based low-order models derived directly from the governing Navier-Stokes equations. Key findings and accomplishments are summarized below. In Chapter 3 we developed an extension to the resolvent framework to account for two- dimensional riblets. The effect of such surfaces is introduced using a volume penalization tech- nique. Predictions made using this extended framework for thin rectangular riblets are in agree- ment with prior experiments and numerical simulations. Specifically, the variation in gain for the resolvent mode that serves as a surrogate for the near-wall cycle reproduces the drag reduc- tion trends observed in previous studies. This extended framework also predicts the appearance 120 of highly energetic spanwise constant rollers resembling Kelvin-Helmholtz vortices, which are hy- pothesized to be responsible for the degradation in drag performance for riblets with a high aspect ratio. In Chapter 4 we combine the framework developed in Chapter 3 with a non-linear conjugate gradient-based optimization algorithm to search for the optimal two-dimensional riblet profile. To reduce the size of the search space we consider riblet profiles which are constructed using simple Bézier curves. The benefits of using such splines lie in the ability to construct locally smooth and easily scalable profiles and in their ability to describe a continuous smooth curve using a finite set of control points. Although our model has some limitations the optimal profiles extracted using this framework are in line with prior experiments and theoretical predictions. Our predictions suggest that the optimal profile resembles a thin blade-like riblet with s + opt ≈ 18−20. Given the monotonic between the riblet height and the near near-wall mode the profiles determined have not been optimized for height. Future models and improvements include using a non-symmetric riblet profile and modifying our objective function to include KH-like modes to generate predictions for the optimal height. Lastly,inChapter5wedevelopedandexploredaframeworktoanalyzetheeffectsofanisotropic porous substrates on wall bounded turbulent flows. In this extension, the resolvent framework is formulated using the Volume-Averaged Navier-Stokes equations in which the effect of permeable substrate is modeled using a Darcy-type drag term. Similar to the cases of riblets, we shown how a single resolvent mode which acts as a surrogate model for the dynamically important near-wall streaks reproduces the drag decreasing trend observed in numerical simulations. In our analy- sis, we have focused on materials that have high porosity and equal wall-normal and spanwise permeability. Model predictions show that the decrease in the near-wall mode gain is linearly proportional to the difference between the streamwise and spanwise slip length scales given by 121 p K + x − p K + z . These trends and predictions observed in our models are in line with the trends theorized and observed in numerical simulations. Resolvent predictions indicate that materials which have high anisotropy lead to greater mode suppression. The resolvent framework also pre- dictstheappearanceofenergeticspanwisecoherentrollerswithastreamwisewavelengthλ + x ≈ 150 for materials in which the wall-normal permeability length scale is greater than p K + y ≈ 0.4. The degradation ofthe linear drag regime in anisotropicporous substratesis hypothesizedtobe caused by the appearance of these spanwise coherent structures resembling Kelvin-Helmholtz vortices. The extensions presented here to the resolvent framework have some limitations, one of these limitations which we have discussed in some detail, is the need for a mean velocity profile. The resultspresentedhereprimarilyuseapredictivemeanvelocityprofileswhichhavebeensynthesized using an eddy-viscosity model. Nonetheless, these predictive mean profiles have generated useful a priori predictions. We have demonstrated how this purely predictive framework, utilizing the gain of near-wall mode can provide a measure of drag reduction. We have shown how this frameworkcanalsopredicttheemergenceofdrag-degradingspanwiseconstantrollersasafunction of the properties of riblet geometry or porous material permeability. Importantly, the reduced- complexitymodelingframeworkdevelopedhererequiresminimalcomputation. Allthepredictions presented in this thesis were generated using a standard workstation computer. Finally, though we did not pursue these issues in the present work, the modeling framework couldalsoprovideimportantinsightsintoReynoldsnumberscalingandalternateribletgeometries or porous substrate properties. Firstly, trends at higher Reynolds numbers can be quantified and compared with the predictions generated for Re τ = 180 for both riblets and anisotropic porous substrates. Given the minimal computational power required to produce predictions, scaling laws for non-conventional riblet geometries could also be expediently evaluated. In chapter 5 we focused our analysis on anisotropic porous materials defined by their streamwise, wall-normal, 122 and spanwise permeability, K = diag(K x ,K y ,K z ) in where K + y = K + z . Using the resolvent framework, a more comprehensive analysis can be performed and can provide useful insight into the scaling trends for anisotropic materials defined by a full permeability tensor. In summary, the framework presented here represents to our knowledge, the first computa- tionally inexpensive physics-based reduced-order models to predict the effect of two-dimensional riblets and anisotropic porous substrates. Future development includes extending the resolvent framework to account for three-dimensional (streamwise varying) riblets. To tackle this problem we need to develop a mean velocity profile that represents and describes the underlying dynamics. Development of more sophisticated models, both for riblets and anisotropic porous substrates, in- cludes generating predictions for second-order statistics and using these statistics to estimating the drag change over these surfaces. 123 References [1] Abderrahaman-Elena, N. and García-Mayoral, R. “Analysis of anisotropically permeable surfacesforturbulentdragreduction”.In: Physic al R eview Fluids2(11Nov.2017),p.114609. [2] Adrian, R. J. “Hairpin vortex organization in wall turbulence”. In: Physics of Fluids 19.4 (2007), p. 041301. [3] Angot, P. “Analysis of singular perturbations on the Brinkman problem for fictitious do- main models of viscous flows”. In: Mathematic al metho ds in the applie d scienc es 22.16 (1999), pp. 1395–1412. [4] Aurentz, J. and Trefethen, L. “Block Operators and Spectral Discretizations”. In: SIAM R eview 59.2 (2017), pp. 423–446. [5] Bamieh, B. and Dahleh, M. “Energy amplification in channel flows with stochastic excita- tion”. In: Physics of Fluids 13.11 (2001), pp. 3258–3269. [6] Beavers, G. S. and Joseph, D. D. “Boundary conditions at a naturally permeable wall”. In: Journal of Fluid Me chanics 30.1 (1967), pp. 197–207. [7] Bechert, D. W. and Bartenwerfer, M. “The viscous flow on surfaces with longitudinal ribs”. In: Journal of Fluid Me chanics 206 (1989), pp. 105–129. [8] Bechert, D. W., Bruse, M., and Hage, W. “Experiments with three-dimensional riblets as an idealized model of shark skin”. In: Exp eriments in Fluids 28.5 (2000), pp. 403–412. [9] Bechert, D. W., Bruse, M., Hage, W., Van Der Hoeven, J. G. T., and Hoppe, G. “Experi- ments on drag-reducing surfaces and their optimization with an adjustable geometry”. In: Journal of Fluid Me chanics 338 (1997), pp. 59–87. [10] Beneddine, S., Sipp, D., Arnault, A., Dandois, J., and Lesshafft, L. “Conditions for validity of mean flow stability analysis”. In: Journal of Fluid Me chanics 798 (2016), pp. 485–504. [11] Benschop, H. and Breugem, W.-P. “Drag reduction by herringbone riblet texture in direct numerical simulations of turbulent channel flow”. In: Journal of T urbulenc e 18.8 (Aug. 2017), pp. 717–759. [12] Bézier, P. E. “Example of an Existing System in the Motor Industry: The Unisurf System”. In: Pr o c e e dings of the R oyal So ciety of L ondon. Series A, Mathematic al and Physic al Scienc es 321.1545 (1971), pp. 207–218. [13] Blackwelder, R. F. and Eckelmann, H. “The spanwise structure of the bursting phe- nomenon”. In: Structur e and Me chanisms of T urbulenc e I. Ed. by H. Fiedler. Berlin, Hei- delberg: Springer Berlin Heidelberg, 1978, pp. 190–204. [14] Blackwelder, R. F. and Eckelmann, H. “Streamwise vortices associated with the bursting phenomenon”. In: Journal of Fluid Me chanics 94.03 (Oct. 1979), p. 577. 124 [15] Breugem, W. P., Boersma, B. J., and Uittenbogaar, R. E. “The influence of wall perme- ability on turbulent channel flow”. In: Journal of Fluid Me chanics 562 (2006), pp. 35– 72. [16] Brinkman, H. C. “A calculation of the viscous force exerted by a flowing fluid on a dense swarm of particles”. In: Flow, T urbulenc e and Combustion 1.1 (Dec. 1949), p. 27. [17] Busse, A. and Sandham, N. D. “Influence of an anisotropic slip-length boundary condition on turbulent channel flow”. In: Physics of Fluids 24.5 (May 2012), p. 055111. [18] Butler, K. M. and Farrell, B. F. “Three‐dimensional optimal perturbations in viscous shear flow”. In: Physics of Fluids A: Fluid Dynamics 4.8 (1992), pp. 1637–1650. [19] Chavarin, A., Efstathiou, C., Vijay, S., and Luhar, M. “Contribution of dispersive stress to skin friction drag in turbulent flow over riblets”. In: Pr o c e e dings of the 11th International Symp osium on T urbulenc e and She ar Flow Phenomena, TSFP 2019 (Southampton, United Kingdom). July 2019. [20] Chavarin, A., Efstathiou, C., Vijay, S., and Luhar, M. “Resolvent-based design and ex- perimental testing of porous materials for passive turbulence control”. In: International Journal of He at and Fluid Flow 86 (2020), p. 108722. [21] Chavarin,A.andLuhar,M.“ResolventAnalysisforTurbulentChannelFlowwithRiblets”. In: AIAA Journal 58.2 (July 2019), pp. 589–599. [22] Choi, H., Moin, P., and Kim, J. “Active turbulence control for drag reduction in wall- bounded flows”. In: Journal of Fluid Me chanics 262 (1994), pp. 75–110. [23] Choi, H., Moin, P., and Kim, J. “On the effect of riblets in fully developed laminar channel flows”. In: Physics of Fluids A: Fluid Dynamics 3.8 (1991), pp. 1892–1896. [24] Choi, H., Moin, P., and Kim, J. “Direct numerical simulation of turbulent flow over riblets”. In: Journal of Fluid Me chanics 255 (1993), pp. 503–539. [25] Chu, D. C. and Karniadakis, G. E. “A direct numerical simulation of laminar and turbulent flow over riblet-mounted surfaces”. In: Journal of Fluid Me chanics 250 (1993), pp. 1–42. [26] Cohen, A. I. “Rate of Convergence of Several Conjugate Gradient Algorithms”. In: SIAM Journal on Numeric al A nalysis 9.2 (1972), pp. 248–259. [27] Crowder, H. and Wolfe, P. “Linear Convergence of the Conjugate Gradient Method”. In: IBM Journal of R ese ar ch and Development 16.4 (1972), pp. 431–433. [28] Endrikat, S., Modesti, D., García-Mayoral, R., Hutchins, N., and Chung, D. “Kelvin – Helmholtz Rollers in Turbulent Flow over Riblets”. In: 21st A ustr alasian Fluid Me chanics Confer enc e(Adelaide,SouthAustralia).Ed.byT.C.W.LauandR.M.Kelso.Australasian Fluid Mechanics Society, Dec. 2018. [29] Fardad,M.,Jovanović,M.R.,andBamieh,B.“Frequencyanalysisandnormsofdistributed spatially periodic systems”. In: IEEE T r ans. A utomat. Contr ol 53 (Nov. 2008), pp. 2266– 2279. [30] Farrell, B. F. and Ioannou, P. J. “Stochastic forcing of the linearized Navier–Stokes equa- tions”. In: Physics of Fluids A: Fluid Dynamics 5.11 (Nov. 1993), pp. 2600–2609. 125 [31] Fletcher, R. and Reeves, C. M. “Function minimization by conjugate gradients”. In: The Computer Journal 7.2 (Jan. 1964), pp. 149–154. [32] Gad-el-Hak,M. Flow Contr ol: Passive, A ctive, and R e active Flow Management.Cambridge University Press, 2007. [33] García-Mayoral, R., Gómez-de-Segura, G., and Fairhall, C. T. “The control of near-wall turbulence through surface texturing”. In: Fluid Dynamics R ese ar ch 51.1 (Jan. 2019), p. 011410. [34] García-Mayoral, R. and Jiménez, J. “Drag reduction by riblets”. In: Philosophic al T r ansac- tions of the R oyal So ciety of L ondon A: Mathematic al, Physic al and Engine ering Scienc es 369.1940 (2011), pp. 1412–1427. [35] García-Mayoral, R. and Jiménez, J. “Hydrodynamic stability and breakdown of the viscous regime over riblets”. In: Journal of Fluid Me chanics 678 (2011), pp. 317–347. [36] Garcı́a-Mayoral, R. and Jiménez, J. “Scaling of turbulent structures in riblet channels up to Re τ ≈ 550”. In: Physics of Fluids 24.10 (2012), pp. 105–101. [37] Goldstein, D., Handler, R., and Sirovich, L. “Direct numerical simulation of turbulent flow over a modeled riblet covered surface”. In: Journal of Fluid Me chanics 302 (1995), pp. 333– 376. [38] Goldstein, D. and Tuan, T.-C. “Secondary flow induced by riblets”. In: Journal of Fluid Me chanics 363 (1998), pp. 115–151. [39] Gómez, F., Blackburn, H. M., Rudman, M., Sharma, A. S., and McKeon, B. J. “A reduced- ordermodelofthree-dimensionalunsteadyflowinacavitybasedontheresolventoperator”. In: Journal of Fluid Me chanics 798 (2016). [40] Gómez-de-Segura, G. and García-Mayoral, R. “Turbulent drag reduction by anisotropic permeable substrates – analysis and direct numerical simulations”. In: Journal of Fluid Me chanics 875 (Sept. 2019), pp. 124–172. [41] Gómez-de-Segura, G., Sharma, A., and García-Mayoral, R. “Turbulent Drag Reduction Using Anisotropic Permeable Substrates”. In: Flow, T urbulenc e and Combustion 100.4 (June 2018), pp. 995–1014. [42] Griva,I.,Nash,S.,andSofer,A. Line ar and Nonline ar Optimization: Se c ond Edition.Other Titles in Applied Mathematics. Society for Industrial and Applied Mathematics (SIAM, 3600 Market Street, Floor 6, Philadelphia, PA 19104), 2009. [43] Hahn, S., Je, J., and Choi, H. “Direct numerical simulation of turbulent channel flow with permeable walls”. In: Journal of Fluid Me chanics 450 (2002), pp. 259–285. [44] Hamilton, J. M., Kim, J., and Waleffe, F. “Regeneration mechanisms of near-wall turbu- lence structures”. In: Journal of Fluid Me chanics 287 (1995), pp. 317–348. [45] Henningson, D. “Bypass Transition and Linear Growth Mechanisms”. In: A dvanc es in T urbulenc e V. Ed. by R. Benzi. Springer Netherlands, 1995, pp. 190–204. [46] Hwang,Y.andCossu,C.“Linearnon-normalenergyamplificationofharmonicandstochas- ticforcingintheturbulentchannelflow”.In: Journal of Fluid Me chanics664(2010),pp.51– 73. 126 [47] Itoh, M., Tamano, S., Iguchi, R., Yokota, K., Akino, N., Hino, R., and Kubo, S. “Turbulent drag reduction by the seal fur surface”. In: Physics of Fluids 18.6 (2006), pp. 065–102. [48] Jiménez, J. and Pinelli, A. “The autonomous cycle of near-wall turbulence”. In: Journal of Fluid Me chanics 389 (1999), pp. 335–359. [49] Jiménez, J. “On the structure and control of near wall turbulence”. In: Physics of Fluids 6.2 (1994), pp. 944–953. [50] Jiménez, J. “Near-wall turbulence”. In: Physics of Fluids 25.10 (2013), p. 101302. [51] Jiménez, J. and Moin, P. “The minimal flow unit in near-wall turbulence”. In: Journal of Fluid Me chanics 225 (1991). [52] Jiménez, J., Pinelli, A., Jiménez, J., and Pinelli, A. “Wall turbulence-How it works and how to damp it”. In: 4th S he ar Flow Contr ol Confer enc e. June 1997. [53] Jovanović, M. R. “Modeling, analysis, and control of spatially distributed systems”. Ph.D. dissertation. University of California, Santa Barbara, Sept. 2004. [54] Jovanović, M. R. “From bypass transition to flow control and data-driven turbulence mod- eling:Aninput-outputviewpoint”.In: A nnu. R ev. Fluid Me ch.53.1(2021).doi:10.1146/annurev- fluid-010719-060244; also arXiv:2003.10104. [55] Jovanović, M. R. “Turbulence suppression in channel flows by small amplitude transverse wall oscillations”. In: Physics of Fluids 20.1 (2008), p. 014101. [56] Jovanović, M. R. and Bamieh, B. “Componentwise energy amplification in channel flows”. In: Journal of Fluid Me chanics 534 (2005), pp. 145–183. [57] Khadra, K., Angot, P., Parneix, S., and Caltagirone, J.-P. “Fictitious domain approach for numerical modelling of Navier–Stokes equations”. In: International Journal for Numeric al Metho ds in Fluids 34.8 (2000), pp. 651–684. [58] Kline, S. J., Reynolds, W. C., Schraub, F. A., and Runstadler, P. W. “The structure of turbulent boundary layers”. In: Journal of Fluid Me chanics 30.4 (1967), pp. 741–773. [59] Kreplin, H.-P. and Eckelmann, H. “Propagation of perturbations in the viscous sublayer and adjacent wall region”. In: Journal of Fluid Me chanics 95.2 (1979), pp. 305–322. [60] Kuwata, Y. and Suga, K. “Direct numerical simulation of turbulence over anisotropic porous media”. In: Journal of Fluid Me chanics 831 (2017), pp. 41–71. [61] Lācis, U. and Bagheri, S. “A framework for computing effective boundary conditions at the interface between free fluid and a porous medium”. In: Journal of Fluid Me chanics 812 (2017), pp. 866–889. [62] Lācis, U., Sudhakar, Y., Pasche, S., and Bagheri, S. “Transfer of mass and momentum at rough and porous surfaces”. In: Journal of Fluid Me chanics 884 (2020), A21. [63] Lang, A., Habegger, M. L., and Motta, P. “Shark skin drag reduction”. In: Encyclop e dia of Nanote chnolo gy 19 (2012), pp. 2394–2400. [64] T. C. W. Lau and R. M. Kelso, eds. 21st A ustr alasian Fluid Me chanics Confer enc e (Ade- laide, South Australia). Australasian Fluid Mechanics Society, Dec. 2018. 127 [65] Lee, S.-J. and Lee, S.-H. “Flow field analysis of a turbulent boundary layer over a riblet surface”. In: Exp eriments in Fluids 30.2 (2001), pp. 153–166. [66] Lieu, B. K., Moarref, R., and Jovanović, M. R. “Controlling the onset of turbulence by streamwise traveling waves. Part 2: Direct numerical simulations”. In: Journal of Fluid Me chanics 663 (Nov. 2010), pp. 100–119. [67] Luchini, P., Manzo, F., and Pozzi, A. “Resistance of a grooved surface to parallel flow and cross-flow”. In: Journal of Fluid Me chanics 228 (1991), pp. 87–109. [68] Luhar, M., Sharma, A. S., and McKeon, B. J. “Opposition control within the resolvent analysis framework”. In: Journal of Fluid Me chanics 749 (2014), pp. 597–626. [69] Luhar, M., Sharma, A. S., and McKeon, B. J. “A framework for studying the effect of compliant surfaces on wall turbulence”. In: Journal of Fluid Me chanics 768 (2015), pp. 415– 441. [70] Luhar, M., Sharma, A. S., and McKeon, B. J. “On the design of optimal compliant walls for turbulence control”. In: Journal of T urbulenc e 17.8 (2016), pp. 787–806. [71] McClure, P. D., Smith, B. R., Baker, W., and Yagle, P. “Design and Testing of Conven- tional Riblets and 3-D Riblets with Streamwise Variable Height”. In: 55th AIAA A er osp ac e Scienc es Me eting. 2017, p. 0048. [72] McKeon, B. J. “The engine behind (wall) turbulence: perspectives on scale interactions”. In: Journal of Fluid Me chanics 817 (2017), P1. [73] McKeon, B. J. and Sharma, A. S. “A critical-layer framework for turbulent pipe flow”. In: Journal of Fluid Me chanics 658 (2010), pp. 336–382. [74] Mei, C. C. and Vernescu, B. Homo genization Metho ds for Multisc ale Me chanics. WORLD SCIENTIFIC, Sept. 2010. [75] Minale, M. “Momentum transfer within a porous medium. II. Stress boundary condition”. In: Physics of Fluids 26.12 (Dec. 2014), p. 123102. [76] Moarref, R. and Jovanović, M. R. “Controlling the onset of turbulence by streamwise traveling waves. Part 1: Receptivity analysis”. In: Journal of Fluid Me chanics 663 (Nov. 2010), pp. 70–99. [77] Moarref, R. and Jovanović, M. R. “Model-based design of transverse wall oscillations for turbulent drag reduction”. In: Journal of Fluid Me chanics 707 (2012), pp. 205–240. [78] Moarref, R., Sharma, A. S., Tropp, J. A., and McKeon, B. J. “Model-based scaling of the streamwise energy density in high-Reynolds-number turbulent channels”. In: Journal of Fluid Me chanics 734 (2013), pp. 275–316. [79] Modesti, D., Endrikat, S., García-Mayoral, R., Hutchins, N., and Chung, D. “Contribution of dispersive stress to skin friction drag in turbulent flow over riblets”. In: Pr o c e e dings of the 11th International Symp osium on T urbulenc e and She ar Flow Phenomena, TSFP 2019 (Southampton, United Kingdom). July 2019. [80] Modesti,D.,Endrikat,S.,García-Mayoral,R.,Hutchins,N.,andChung,D.“Form-induced stress in turbulent flow over riblets”. In: 21st A ustr alasian Fluid Me chanics Confer enc e 128 (Adelaide, South Australia). Ed. by T. C. W. Lau and R. M. Kelso. Australasian Fluid Mechanics Society, Dec. 2018. [81] Nakashima, S., Fukagata, K., and Luhar, M. “Assessment of suboptimal control for tur- bulent skin friction reduction via resolvent analysis”. In: Journal of Fluid Me chanics 828 (2017), pp. 496–526. [82] Nocedal, J. and Wright, S. J. Numeric al Optimization. second. New York, NY, USA: Springer, 2006. [83] Ochoa-Tapia, J. A. and Whitaker, S. “Momentum transfer at the boundary between a porous medium and a homogeneous fluid—I. Theoretical development”. In: International Journal of He at and Mass T r ansfer 38.14 (1995), pp. 2635–2646. [84] Ochoa-Tapia, J. A. and Whitaker, S. “Momentum transfer at the boundary between a porous medium and a homogeneous fluid—II. Comparison with experiment”. In: Interna- tional Journal of He at and Mass T r ansfer 38.14 (1995), pp. 2647–2655. [85] Orlandi, P. and Jiménez, J. “On the generation of turbulent wall friction”. In: Physics of Fluids 6.2 (1994), pp. 634–641. [86] Peet, Y., Sagaut, P., and Charron, Y. “Turbulent drag reduction using sinusoidal riblets with triangular cross-section”. In: 38th Fluid Dynamics Confer enc e and Exhibit. Seattle, Washington: American Institute of Aeronautics and Astronautics, June 2008. [87] Peskin, C. S. “The immersed boundary method”. In: A cta Numeric a 11 (2002), pp. 479– 517. [88] Polak, E. and Ribière, G. “Note sur la convergence de méthodes de directions conjuguées”. In: ESAIM: Mathematic al Mo del ling and Numeric al A nalysis-Mo délisation Mathématique et A nalyse Numérique 3.R1 (1969), pp. 35–43. [89] Polyak, B. “Some methods of speeding up the convergence of iteration methods”. In: USSR Computational Mathematics and Mathematic al Physics 4.5 (1964), pp. 1–17. [90] Powell, M. J. D. “Convergence Properties of Algorithms for Nonlinear Optimization”. In: SIAM R eview 28.4 (1986), pp. 487–500. [91] Powell, M. J. D. “Restart procedures for the conjugate gradient method”. In: Mathematic al Pr o gr amming 12.1 (Dec. 1977), pp. 241–254. [92] Powell, M. J. D. “Some convergence properties of the conjugate gradient method”. In: Mathematic al Pr o gr amming 11.1 (1976), pp. 42–49. [93] Pr o c e e dings of the 11th International Symp osium on T urbulenc e and She ar Flow Phenom- ena, TSFP 2019 (Southampton, United Kingdom). July 2019. [94] Ran, W., Zare, A., and Jovanović, M. R. “Model-based design of riblets for turbulent drag reduction”. In: Journal of Fluid Me chanics (2020). [95] Rastegari, A. and Akhavan, R. “The common mechanism of turbulent skin-friction drag reduction with superhydrophobic longitudinal microgrooves and riblets”. In: Journal of Fluid Me chanics 838 (2018), pp. 68–104. 129 [96] Reddy, S. C., Schmid, P. J., Baggett, J. S., and Henningson, D. S. “On stability of stream- wise streaks and transition thresholds in plane channel flows”. In: Journal of Fluid Me- chanics 365 (June 1998), pp. 269–303. [97] Reynolds, W. C. and Tiederman, W. G. “Stability of turbulent channel flow, with appli- cation to Malkus’s theory”. In: Journal of Fluid Me chanics 27.2 (1967), pp. 253–272. [98] Robinson, S. K. “Coherent motions in the turbulent boundary layer”. In: A nnual R eview of Fluid Me chanics 23.1 (1991), pp. 601–639. [99] Rosti,M.E.,Brandt,L.,andPinelli,A.“Turbulentchannelflowoverananisotropicporous wall–drag increase and reduction”. In: Journal of Fluid Me chanics 842 (2018), pp. 381–394. [100] Rosti, M. E., Cortelezzi, L., and Quadrio, M. “Direct numerical simulation of turbulent channel flow over porous walls”. In: Journal of Fluid Me chanics 784 (2015), pp. 396–442. [101] Samanta, A., Vinuesa, R., Lashgari, I., Schlatter, P., and Brandt, L. “Enhanced secondary motion of the turbulent flow through a porous square duct”. In: Journal of Fluid Me chanics 784 (2015), pp. 681–693. [102] El-Samni, O., Chun, H., and Yoon, H. “Drag reduction of turbulent flow over thin rectan- gularriblets”.In: International Journal of Engine ering Scienc e 45.2–8(Feb.2007),pp.436– 454. [103] Schmidt, O. T., Towne, A., Rigas, G., Colonius, T., and Brès, G. A. “Spectral analysis of jet turbulence”. In: Journal of Fluid Me chanics 855 (2018), pp. 953–982. [104] Schoppa, W. and Hussain, F. “Coherent structure generation in near-wall turbulence”. In: Journal of Fluid Me chanics 453 (Feb. 2002). [105] Schultz, M. P., Bendick, J. A., Holm, E. R., and Hertel, W. M. “Economic impact of biofouling on a naval surface ship”. In: Biofouling 27.1 (2011). PMID: 21161774, pp. 87–98. [106] Sharma, A. S. and McKeon, B. J. “On coherent structure in wall turbulence”. In: Journal of Fluid Me chanics 728 (2013), pp. 196–238. [107] Sharp Jr., J. M. and Simmons, C. T. “The compleat Darcy: New lessons learned from the first English translation of les fontaines publiques de la V il le de Dijon”. In: Gr oundwater 43.3 (2005), pp. 457–460. [108] Sirovich, L. and Karlsson, S. “Turbulent drag reduction by passive mechanisms”. en. In: Natur e 388.6644 (Aug. 1997), pp. 753–755. [109] Smith, C. R. and Metzler, S. P. “The characteristics of low-speed streaks in the near-wall region of a turbulent boundary layer”. In: Journal of Fluid Me chanics 129 (1983), pp. 27– 54. [110] Suga, K., Nakagawa, Y., and Kaneda, M. “Spanwise turbulence structure over permeable walls”. In: Journal of Fluid Me chanics 822 (2017), pp. 186–201. [111] Suzuki, Y. and Kasagi, N. “Turbulent drag reduction mechanism above a riblet surface”. In: AIAA Journal 32.9 (1994), pp. 1781–1790. 130 [112] Taira, K., Brunton, S. L., Dawson, S. T., Rowley, C. W., Colonius, T., McKeon, B. J., Schmidt, O. T., Gordeyev, S., Theofilis, V., and Ukeiley, L. S. “Modal Analysis of Fluid Flows: An Overview”. In: AIAA Journal 55.12 (2017), pp. 4013–4041. [113] Tardu, S. “Near-Wall Coherent Structures: History, Identification and Detection”. In: T r ansp ort and Coher ent Structur es in W al l T urbulenc e. John Wiley & Sons, Ltd, 2014. Chap. 3, pp. 129–198. [114] Tilton, N. and Cortelezzi, L. “The destabilizing effects of wall permeability in channel flows: A linear stability analysis”. In: Physics of Fluids 18.5 (2006), p. 051702. [115] Tilton, N. and Cortelezzi, L. “Linear stability analysis of pressure-driven flows in channels with porous walls”. In: Journal of Fluid Me chanics 604 (2008), pp. 411–445. [116] Toedtli, S. S., Luhar, M., and McKeon, B. J. “Predicting the response of turbulent channel flow to varying-phase opposition control: Resolvent analysis as a tool for flow control design”. In: Phys. R ev. Fluids 4 (7 July 2019), p. 073905. [117] Towne, A., Lozano-Durán, A., and Yang, X. “Resolvent-based estimation of space–time flow statistics”. In: Journal of Fluid Me chanics 883 (2020), A17. [118] Trefethen, L. N., Trefethen, A. E., Reddy, S. C., and Driscoll, T. A. “Hydrodynamic Sta- bility Without Eigenvalues”. In: Scienc e, New Series 261.5121 (1993), pp. 578–584. [119] Viswanath, P. “Aircraft viscous drag reduction using riblets”. In: Pr o gr ess in A er osp ac e Scienc es 38.6–7 (Aug. 2002), pp. 571–600. [120] Waleffe, F. “On a self-sustaining process in shear flows”. In: Physics of Fluids 9.4 (1997), pp. 883–900. [121] Wallace, J. M. “Quadrant Analysis in Turbulence Research: History and Evolution”. In: A nnual R eview of Fluid Me chanics 48.1 (2016), pp. 131–158. [122] Walsh, M. and Lindemann, A. “Optimization and application of riblets for turbulent drag reduction”. In: 22nd A er osp ac e Scienc es Me eting. American Institute of Aeronautics and Astronautics, Jan. 1984. [123] Walsh, M. J. “Turbulent boundary layer drag reduction using riblets”. In: 20th A er osp ac e Scienc es Me eting. American Institute of Aeronautics and Astronautics, Jan. 1982, pp. 1–8. [124] Walsh, M. J., Sellers III, W. L., and Mcginley, C. B. “Riblet drag at flight conditions”. In: Journal of A ir cr aft 26.6 (1989), pp. 570–575. [125] Whitaker, S. “Advances in Theory of Fluid Motion in Porous Media”. In: Industrial & Engine ering Chemistry 61.12 (Dec. 1969), pp. 14–28. [126] Whitaker, S. “The Forchheimer equation: A theoretical development”. In: T r ansp ort in Por ous Me dia 25.1 (Oct. 1996), pp. 27–61. [127] Whitaker, S. “Volume Averaging of Transport Equations”. In: Fluid T r ansp ort in Por ous Me dia. Ed. by J. P. D. Plessis. Computational Mechanics Publications, 1997. Chap. 1, pp. 1–60. [128] Whitaker, S. The Metho d of V olume A ver aging. Vol. 13. Theory and Applications of Trans- port in Porous Media. Springer Netherlands, 1999. 131 [129] Yeh, C.-A. and Taira, K. “Resolvent-analysis-based design of airfoil separation control”. In: Journal of Fluid Me chanics 867 (2019), pp. 572–610. [130] Zampogna, G. A. and Bottaro, A. “Fluid flow over and through a regular bundle of rigid fibres”. In: Journal of Fluid Me chanics 792 (Apr. 2016), pp. 5–35. [131] Zare, A., Georgiou, T., and Jovanović, M. “Stochastic Dynamical Modeling of Turbulent Flows”. In: A nnual R eview of Contr ol, R ob otics, and A utonomous Systems 3.1 (2020), pp. 195–219. [132] Zare, A., Jovanović, M. R., and Georgiou, T. T. “Colour of turbulence”. In: Journal of Fluid Me chanics 812 (2017), pp. 636–680. 132

Abstract (if available)

Linked assets

University of Southern California Dissertations and Theses

Conceptually similar

PDF

Modeling and analysis of parallel and spatially-evolving wall-bounded shear flows

PDF

Flow and thermal transport at porous interfaces

PDF

Developing tuned anisotropic porous substrates for passive turbulence control

PDF

Reconstruction and estimation of shear flows using physics-based and data-driven models

PDF

Hydrodynamic stability of two fluid columns of different densities and viscosities subject to gravity

PDF

Numerical study of shock-wave/turbulent boundary layer interactions on flexible and rigid panels with wall-modeled large-eddy simulations

PDF

Pattern generation in stratified wakes

PDF

Large eddy simulations of turbulent flows without use of the eddy viscosity concept

PDF

Tracking and evolution of compressible turbulent flow structures

PDF

Aerodynamics at low Re: separation, reattachment, and control

PDF

Development and applications of a body-force propulsor model for high-fidelity CFD

PDF

Numerical study of flow characteristics of controlled vortex induced vibrations in cylinders

PDF

Multiscale modeling of cilia mechanics and function

PDF

The role of counter-current flow in the modeling and simulation of multi-phase flow in porous media

PDF

Chemical and mechanical deformation of porous media and materials during adsorption and fluid flow

PDF

Efficient inverse analysis with dynamic and stochastic reductions for large-scale models of multi-component systems

PDF

Machine-learning approaches for modeling of complex materials and media

PDF

A subgrid-scale model for large-eddy simulation based on the physics of interscale energy transfer in turbulence

PDF

Multiscale and multiresolution approach to characterization and modeling of porous media: From pore to field scale

PDF

Modeling and simulation of complex recovery processes

Asset Metadata

Creator Chavarin, Andrew (author)

Core Title Model based design of porous and patterned surfaces for passive turbulence control

School Viterbi School of Engineering

Degree Doctor of Philosophy

Degree Program Mechanical Engineering

Publication Date 01/15/2021

Defense Date 10/22/2020

Publisher University of Southern California (original), University of Southern California. Libraries (digital)

Tag anisotropic porous substrates,drag reduction,OAI-PMH Harvest,optimization,passive control of turbulent flows,porous media,reduced order-modeling,resolvent analysis,riblets,spatially periodic systems,turbulence,turbulence modeling

Language English

Contributor Electronically uploaded by the author (provenance)

Advisor Luhar, Mitul (committee chair), Bermejo-Moreno, Ivan (committee member), Jovanovic, Mihailo R. (committee member), Oberai, Assad A. (committee member)

Creator Email achavari@usc.edu

Permanent Link (DOI) https://doi.org/10.25549/usctheses-c89-414845

Unique identifier UC11666720

Identifier etd-ChavarinAn-9235.pdf (filename),usctheses-c89-414845 (legacy record id)

Legacy Identifier etd-ChavarinAn-9235.pdf

Dmrecord 414845

Document Type Dissertation

Rights Chavarin, Andrew

Type texts

Source University of Southern California (contributing entity), University of Southern California Dissertations and Theses (collection)

Access Conditions The author retains rights to his/her dissertation, thesis or other graduate work according to U.S. copyright law. Electronic access is being provided by the USC Libraries in agreement with the a...

Repository Name University of Southern California Digital Library

Repository Location USC Digital Library, University of Southern California, University Park Campus MC 2810, 3434 South Grand Avenue, 2nd Floor, Los Angeles, California 90089-2810, USA