AA277 Instability and transition to turbulence Andr V

AA-277 Instability and transition to turbulence André V. G. Cavalieri Peter Jordan (Visiting researcher, Ciência Sem Fronteiras)

Instability and transition to instability Laminar flows Turbulent flows turbulence ❑ ❑ ❑ Lower skin friction and heat transfer Prone to flow separation Steady flows, no noise generation ❑ Increased skin friction and heat transfer ■ ❑ ❑ Turbulent mixing Less prone to separation Aerodynamic noise Challenges: 1 – Predict transition to turbulence 2 – Understand mechanisms of transition 3 – Flow control: delay, or simply avoid turbulence

Instability and transition to Are there steady (laminar) solutions to the Navier-Stokes turbulence equations? ■ ❑ Yes! ■ ■ Are there small (or sometimes not so small) disturbances in the environment? ❑ Yes! ■ ■ ■ ■ Example: zero-pressure-gradient boundary layer (Blasius): solution for arbitrary Re Surface roughness Vibration Sound waves Incoming vorticity A mosquito flying nearby Are the steady solutions stable to these disturbances? ❑ Answer: it depends. . . strongly on Re (but also on M) Unstable laminar solutions will not prevail on nature

This course ■ Theory ❑ Suggested literature: ■ ■ Criminale, Jackson & Joslin, Theory and computation of hydrodynamic stability, Cambridge University Press, 2003 Schmid & Henningson, Stability and transition in shear flows, Springer, 2001 + various articles to be provided during the course Computation ❑ Matlab programming ■ ■ Trefethen, Spectral methods in Matlab Boyd, Chebyshev and Fourier spectral methods Lectures with assignments Goals: ❑ Build a toolbox to solve some canonical stability problems ❑ Learn with examples Grades: ❑ Theoretical and computational assignments ❑ ■ ■ ■ Material can be found in ftp: //161. 24. 15. 247/Andre/AA 277

Videos! Rayleigh-Taylor instability http: //www. youtube. com/watch? v=NI 85 o. C-3 m. J 0 Pipe flow: Reynolds experiment (1883) http: //www. youtube. com/watch? v=xi. X 5 Pf. Fxm. Is Boundary-layer transition: Tollmien-Schlichting waves http: //www. youtube. com/watch? v=nb 22 g 6 ky 2 XE Kelvin-Helmholtz instability http: //www. youtube. com/watch? v=Ub. Afvca. Yr 00 Rayleigh-Bénard convection http: //www. youtube. com/watch? v=Uh. Im. CA 5 Ds. Q 0

Videos! Rayleigh-Taylor instability http: //www. youtube. com/watch? v=NI 85 o. C-3 m. J 0 Pipe flow: Reynolds experiment (1883) http: //www. youtube. com/watch? v=xi. X 5 Pf. Fxm. Is Boundary-layer transition: Tollmien-Schlichting waves http: //www. youtube. com/watch? v=nb 22 g 6 ky 2 XE WAVES!!! Kelvin-Helmholtz instability http: //www. youtube. com/watch? v=Ub. Afvca. Yr 00 Rayleigh-Bénard convection http: //www. youtube. com/watch? v=Uh. Im. CA 5 Ds. Q 0 Fourier

Basic equations: Rayleigh and Orr-Sommerfeld equations

Basic equations: Rayleigh and Orr-Sommerfeld equations Basic equations for two-dimensional incompressible flow with constant viscosity Mass conservation Navier-Stoles equation

Basic equations: Rayleigh and Orr-Sommerfeld equations

Basic equations: Rayleigh and Orr-Sommerfeld equations Mass conservation, non-dimensional Navier-Stoles equation, non-dimensional

Basic equations: Rayleigh and Orr-Sommerfeld equations

Basic equations: Rayleigh and Orr-Sommerfeld equations We consider infinitesimal disturbances to the base flow, such that << << Linearisation is formally justified (but only for initial stages of transition) Tools for linear PDEs: Superposition of solutions Fourier transforms Linear is simpler! Linear algebra

Basic equations: Rayleigh and Orr-Sommerfeld equations tedious but straightforward algebra

Basic equations: Rayleigh and Orr-Sommerfeld equations Other variables:

Basic equations: Rayleigh and Orr-Sommerfeld equations Coefficients are • constant in x and t • variable in y Normal mode Ansatz: Wave propagating in x, with: • phase velocity c • wavenumber α • frequency ω=αc All complex-valued!

Basic equations: Rayleigh and Orr-Sommerfeld equations Coefficients are • constant in x and t • variable in y Normal mode Ansatz:

Basic equations: Rayleigh and Orr-Sommerfeld equations Coefficients are • constant in x and t • variable in y Normal mode Ansatz: Neglect viscous effects:

Basic equations: Rayleigh and Orr-Sommerfeld equations Normal mode Ansatz: 4 th order, 4 BCs 2 nd order, 2 BCs Boundary conditions: Rigid wall: v = 0 (normal velocity is zero) – Rayleigh & Orr-Sommerfeld dv/dy = 0 (since u=0 and du/dx=0, no-slip condition) – only O-S Infinity: bounded v, dv/dy

Basic equations: Rayleigh and Orr-Sommerfeld equations Normal mode Ansatz: Temporal stability: • wavenumber α is real-valued • frequency ω= ωR + ωI (=αc) is complex-valued Space-periodic disturbances are added to the base-flow Will they grow or decay in time?

Basic equations: Rayleigh and Orr-Sommerfeld equations Normal mode Ansatz: Spatial stability: • frequency ω is real-valued • wavenumber α=αR+iαI is complex-valued Time-periodic disturbances are added to the base-flow Will they grow or decay in space? (see Matlab demo)

Temporal stability Real-valued wavenumber is given Base-flow U is given Reynolds number is given c and v are unknowns Rearrange O-S equation as or

Temporal stability Rearrange O-S equation as or or Generalised eigenvalue problem with linear operators L and F c is the eigenvalue Numerics: discretise operators as matrices, then eig(L, F)

A simpler eigenvalue problem or Let’s first look at a simpler eigenvalue problem with BC: v=0 for y=-1 and y=1 This problem appears in: • Linear acoustics • Vibration of membranes • Heat transfer etc

A simpler eigenvalue problem or Let’s first look at a simpler eigenvalue problem with BC: v=0 for y=-1 and y=1 Analytical solution:

A simpler eigenvalue problem Let’s first look at a simpler eigenvalue problem with BC: v=0 for y=-1 and y=1 Analytical solution: Numerical solution: needs discretization of derivative in matrix form i. e. a differentiation matrix Choices to build differentiation matrices: • Finite differences (example: EVPExample_Finite. Differences. m) • Spectral method (example: EVPExample_Chebyshev. m )

Finite differences v is sampled at N gridpoints, spacing=Δx interior points: centered differences boundary: forward or backward differences (less accurate) Derivative is estimated from values of function on the neighbourhood of xi Differentiation matrix is sparse Freedom to choose the grid

Pseudo-spectral method – Chebyshev polynomials Tn(x) Orthogonal basis in [-1, 1] Functions in [-1, 1] can be expanded in Chebyshev polynomials Derivatives:

Pseudo-spectral method v is sampled at N+1 Chebyshev gridponts (clustered near boundaries) 1. 2. 3. 4. Find the coefficients of N+1 Chebyshev polynomials from N+1 values of v Calculate derivatives of Chebyshev polynomials Calculate values of dv/dy at the gridpoints Gather everything inside a single matrix Trefethen 2001 Derivative is estimated from polynomial fit applied to N+1 samples Differentiation matrix is full Grid must be of Chebyshev points See chebyshev_introduction. m

Finite differences x pseudo-spectral method Error in the calculation of the derivative Finite differences: algebraic convergence Spectral: exponential convergence machine precision

Finite differences x pseudo-spectral method 2 nd-order finite differences Pseudo-spectral 0. 01 error • Some eigenvalues (and eigenfunctions) well converged, but some (at least half) should not be trusted! • Pseudo-spectral method better suited for accurate predictions • But: revival of finite differences due to sparse-matrix algorithms (compensate with higher N) – see Gennaro et al. 2013

Exercise – numerical solution of Rayleigh equation (temporal stability) Rayleigh equation, rearranged for temporal stability or Exercise: #1: Adapt the code to solve the temporal stability problem for the Rayleigh equation for a mixing layer (U=0. 5(1+tanh(y))), with BC: v 0 for y ±∞ (numerically: v=0 for y= ± H, H>>1) Start with α=0. 1. Compare with Michalke (1964).

Exercise – numerical solution of Rayleigh equation (temporal stability) Useful Matlab commands (help can be obtained with comands “help” and “doc”): diag(U): diagonal matrix, main diagonal = vector U II = eye(size(D)): identity matrix, same size as matrix D ZZ = zeros(size(D)): matrix of zeros, same size as matrix D H = 20; y = y*H; D = D/H; D 2 = D 2/(H^2); : stretches domain from [-1: 1] to [-H: H] A*B is the standard matrix multiplication (number of lines of A should be equal to number of columns of B) A. *B is an element-to-element multiplication (dimensions of A and B shoud be equal) A^2 is A*A (A should be a square matrix) A. ^2 is A. *A (A can have any size; could be a vector) A(1, : ) means all elements in line 1 of matrix A A(: , 5) means all elements in column 5 of matrix A

Expected results Kelvin-Helmholtz (K-H) mode Conjugate of K-H

Expected results If the c, v pair is a solution (eigenvalue-eigenfuncion) of the Rayleigh equation the c*, v* pair is also a solution. Is this true for the Orr-Sommerfeld equation?

Expected results, α=0. 1, absolute value of eigenfunction exponential decay Kelvin-Helmholtz (K-H) mode artificial boundary condition (v=0)

Expected results, α=0. 1 In the region where U is constant, Rayleigh equation reduces to If c ≠ U then Plus or minus signs are chosen such that solution be bounded at infinity.

Expected results, absolute value of eigenfunction α=0. 1 α=0. 4 Artificial boundary condition is more problematic at low α. • Domain size should be large compared to the wavelength

Expected results, absolute value of eigenfunction α=0. 1 α=0. 4 Solution: mapping from infinite (z) to finite domain (y) See appendix E of Boyd for expressions of derivatives

Expected results Kelvin-Helmholtz (K-H) mode Temporal growth rate: ωi=αci What is the fastest growing wavenumber? What is the wavelength of fastest growing disturbances?

Expected results Temporal growth rate of Kelvin-Helmholtz mode This result shows: • Wavenumber of most amplified disturbances • Growth rate of each amplified wavenumber • Range of amplified wavenumbers

Expected results What is the “preferred” wavelength?

Expected results Kelvin-Helmholtz (K-H) mode Phase speed of disturbances: cr What is the phase speed (convection velocity) of the amplified disturbances?

Expected results Phase speed of Kelvin-Helmholtz mode For temporal stability, phase speed = ½ • Average speed between two streams What happens for U(y)=tanh(y)? (check the average phase speeds and growths, observe that the nondimensionalization changes)

Expected results

Expected results Eigenfunctions allow study of shapes of growing disturbances (streamlines, vorticity, etc. . . ) Vorticity Streamlines (Michalke 1964) Initial stages of roll-up of shear layer into “vortices”

Expected results Eigenfunctions allow study of shapes of growing disturbances (streamlines, vorticity, etc. . . ) Vorticity Experiment (Winant & Browand 1974) Obs: faster stream for negative y (Michalke 1964) Initial stages of roll-up of shear layer into “vortices”

Discussion ■ In unsteady, turbulent flows, clear understanding is exception rather than rule ❑ Nonlinear ■ system of PDEs, 4 dimensions, 5 dependent variables. . . Linear stability allows identification of a well-defined flow feature: an instability wave ❑ Appropriate ■ for initial stages of transition Eigenfunctions of the linearised stability problem form a complete basis (Salwen & Grosch 1981) ❑ It is possible to project flow data (from calculations or experiments) onto the eigenfunctions to obtain amplitudes of each mode ■ Problem is much cheaper to solve numerically than a DNS or an LES ❑ Simplicity through linearisation and use of normal modes - PDE becomes an ODE

Further remarks ■ Rayleigh equation neglects viscous effects ❑ ❑ ❑ ■ No information on the critical Reynolds number Asymptotic behaviour for Re → ∞ Predicts initial stage of transition for sufficiently high Re Rayleigh equation becomes singular for c=U ❑ ❑ Critical layer Special care to treat singularity ■ Easy way out: account for viscous effects This branch should not be trusted!

Critical layer modes in Couette flow

Critical layer modes in Couette flow Discontinuity in dv/dy, discontinuity in u

Critical layer modes in Couette flow Phase opposition in u across the critical layer

Exercise – numerical solution of Orr-Sommerfeld equation (temporal stability) or Orr-Sommerfeld equation, rearranged for temporal stability Exercise: #2: Write a code to solve the temporal stability problem for the Orr. Sommerfeld equation for plane Poiseuille flow (U=1 - y 2). Start with α=1. 0 and Re=10000. Compare with Orszag (1971).

Some further numerical tricks If v 1=0 and v. N+1 = 0 (homogeneous Dirichlet boundary conditions), the easiest way to impose BCs is L=L(2: N); F=F(2: N); Trefethen 2001 Implicitly, we take only Chebyshev polynomials satisfying the boundary conditions

Some further numerical tricks Orr-Sommerfeld: we need to impose v 1=0 and v. N+1 = 0 (homogeneous Dirichlet boundary conditions) dv 1/dy=0 and dv. N+1/dy=0 (homogeneous Neumann boundary conditions) Here is the trick: If q(-1)=q(1)=0, then: v(-1)=v(1)=0 AND dv/dy(-1)=dv/dy(1)=0 Write a problem with q as the unknown, and impose q(-1)=q(1)=0.

Exercise #2, expected results Exercise: #2: Write a code to solve the temporal stability problem for the Orr. Sommerfeld equation for plane Poiseuille flow (U=1 - y 2). Start with α=1. 0 and Re=10000. Compare with Orszag (1971). Eigenspectrum ? ? ?

Exercise #2, expected results Recall that not all eigenvalues are converged Orr-Sommerfeld: no analytical solutions How to decide which eigenvalues are accurate?

Exercise #2, expected results Recall that not all eigenvalues are converged Compare results for different discretisations

Exercise #2, expected results Boyd, section 6. 14: Checking

Exercise #2, expected results Recall that not all eigenvalues are converged Compare results for different discretisations Some errors still present; region of high sensitivity of the O-S operator (see Reddy, Schmid & Henningson 1993)

Exercise #2, expected results The Orr-Sommerfeld operator is non-normal • Adjoint Orr-Sommerfeld ≠ Direct Orr-Sommerfeld • Matrix L is not Hermitian (i. e. L*T ≠ L) • High sensitivity of some eigenvalues to numerical errors More on that later!

Exercise #2, expected results Convergence of most unstable mode, α=1, Re=10000 Orszag 1971:

Exercise #2, expected results Growth rate of most unstable mode: Effect of Reynolds number

Exercise #2, expected results Growth rate of most unstable mode: Effect of Reynolds number Low Re is stabilising (effect of viscosity) Below a certain Re, all modes are stable Above a certain Re, one mode becomes unstable (and that is enough for transition!)

Exercise #2, expected results Temporal growth rate of most unstable mode of Poiseuille flow Green, yellow, red: unstable Cyan, blue, white: stable

Exercise #2, expected results Stable Critical Reynolds number (Recr): Unstable Stable Above these value a range of wavenumbers will grow exponentially with time Initial stage of transition to turbulence

Exercise #2, expected results Instability waves are dispersive: different wavenumbers travel at different speeds Most unstable mode in plane Poiseuille fow Each mode has a dispersion relationship

Experimental results Nishioka, Iida & Ichikawa 1975 – channel flow, background turbulence=0. 05%

Experimental results Nishioka, Iida & Ichikawa 1975 – channel flow, background turbulence=0. 05% Note: other experiments report transition at lower Re Patel & Head 1969: Recr = 2500 Karnitz, Potter & Smith 1974: Recr = 5000 (background turbulence: 0. 3%) Subcritical transition More on that later!

Further exercises ■ Exercise #3: Evaluate the Reynolds number effect on the temporal growth rate for a mixing layer (tanh profile, Michalke 1964). Is it possible to determine the critical Reynolds number? ❑ Note: to normalise the Reynolds number with the momentum thickness of the mixing layer, it is preferrable to use U=0. 5(1+tanh(y/2)) – unit momentum thickness ■ Exercise #4: Study the temporal stability of the Blasius boundary layer. ■ Exercise #5: Study the effect of favourable and adverse pressure gradients in the stability of boundary layers using the Falkner-Skan family of velocity profiles ❑ The code in “OS_Falkner_Skan. m” may help in the last two exercises.

Expected results ■ Exercise #3: Evaluate the Reynolds number effect on the temporal growth rate for a mixing layer (tanh profile, Michalke 1964). Is it possible to determine the critical Reynolds number? Recrit ≈ 2!!! Does this make sense?

Expected results Recrit ≈ 2!!! Does this make sense? Recall: from boundary layer theory but we assume U = U(y) Parallel flow hypothesis only applicable for Re>>1

Plane channel flow x Mixing layer Plane Poiseuille flow High Recrit Small range of unstable α Stable for Re → ∞ (check numerically) Low Recrit Greater range of unstable α Unstable for Re → ∞ (Michalke 1964) WHY?

Expected results Blasius boundary layer β=0 ■ Exercise #4: Study the temporal stability of the Blasius boundary layer. ■ Exercise #5: Study the effect of favourable and adverse pressure gradients in the stability of boundary layers using the Falkner-Skan family of velocity profiles Tollmien 1931, 1936 Schlichting 1932, 1933, 1935 Tollmien-Schlichting (T-S)waves

Expected results Visualisation of Tollmien-Schlichting waves

Historical note: Schubauer & Skramstad 1943 Note: Reynolds number based on displacement thickness

Historical note: Schubauer & Skramstad 1943 Eigenfunctions Full lines: experiment Dashed lines: theory (Schlichting) All this with only vibrating ribbons and hot wires. . .

Historical note: Schubauer & Skramstad 1943 Hans W. Liepmann memories, 1997 “Sometime in 1941 with the war in Europe in its second year, I came up the narrow substandard stairs in the Guggenheim laboratory. At the top of the stairs I met Clark Millikan handing a sheaf of papers to von Kármán with the words, ‘It’s a complete German victory!’ I was stunned. Fortunately, however, the victory was not another one won by Hitler but referred to the experimental verification of the [German] theory. . . by Schubauer and Skramstad. . . I still vividly remember the impact these experiments had on me; they forcefully demonstrated. . . the beginning of a new area in transition research. ”

The effect of pressure gradient Blasius boundary layer β=0 Favourable pressure gradient β = 0. 1

Expected results Blasius boundary layer β=0 Adverse pressure gradient β = -0. 1 Reynolds number defined using δ (boundary layer thickness).

Expected results Blasius boundary layer β=0 Adverse pressure gradient β = -0. 15 WHY?

A hint: base flows

Three-dimensionality: Squire transformation and Squire equation Consider a parallel base-flow U(y), with small disturbances u, v, w(x, y, z, t), p(x, y, z, t) Which are the homogeneous directions? What is the appropriate Ansatz?

Three-dimensionality: Squire transformation and Squire equation Which are the homogeneous directions? What is the appropriate Ansatz?

Three-dimensionality: Squire transformation and Squire equation Orr-Sommerfeld equation, three-dimensional disturbances Eigenvalue problem: given Re, α, β, find non-trivial eigenfunction v(y) and eigenvalue ω, as previously.

The Squire equation Besides Orr-Sommerfeld, from linearised Navier-Stokes we can derive an equation for the normal vorticity called Squire equation

The Orr-Sommerfeld-Squire system O-S Squire B. C. : η=0 on a wall Note that O-S is decoupled from Squire and can be solved separately. Two sets of modes: • Orr-Sommerfeld modes: solve eigenvalue problem for v (O-S), then solve ODE for η • Squire modes: take v=0 and solve Exercise: obtain the Squire eigenvalues for plane Poiseuille flow. Can you find an unstable Squire mode?

The Orr-Sommerfeld-Squire system All Squire modes are stable (see demonstration) • this is why Orr. Sommerfeld alone is sufficient to look for unstable modes The Squire equation will become important later in this course.