Page Web du cours

    Introduction aux instabilit�s hydrodynamiques

    D. Fabre, novembre 2018-janvier 2019.

    Presentation of the course

    Objectives of the course

    The concept of Instability generally means the spontaneous transition of a system from a “regular state” (steady, symmetric, …) to a “less regular state” (unsteady, asymetric, …).

    The identification and description of instabilities is a general situation encountered in many situation in fluid dynamics… and more generally in other fields of physics.

    The mathematical study of instabilities originates from the second half of the 19th century and is associated with the name of classical scientists (Rayleigh, Kelvin, Helmholtz, …). However this discipline is still in flourishing development nowadays.

    This course is an introduction to:

    • The various kinds of instabilities encountered in fluid dynamics, with emphasis on their physical origin,
    • The mathematical formalism allowing to describe an instability,
    • The specific numerical methods used in the study of instability phenomena.

    Material provided for this course.

    This webpage gathers the online material for the course. The material consists of various documents, including lecture notes, exercices, and commented programs.

    The commented programs aim at providing simple, self-contained programs incorporating their own documentation, written in Matlab/Octave language. The proposed programs allow to study the main classes of instabilities considered in this course, and to produce most of the numerical results illustrating the course. T

    he programs are largely built using the support and the philosophy of the easystab project launched by J. Hoeffner and his colleagues from Universit� Paris VI.

    The website of the present course is actually hosted by the website of the easystab project, which is itself hosted by the website of the basilisk numerical simulation code.

    How to use the programs provided in these pages :

    • Get the basic bricks of the project (functions to build matricial operators, etc…) by downloading and unpacking
    • To download a program, open its web page, click “raw page source” in the left column, cut-paste into the matlab/octave editor window, and save it with the corresponding name and “.m” extension.
    • Once it is downloaded, just play with the programs and have fun !
    • If you have improved a program and wish to share it, you are welcome ! the collection of programs is handled under the form of a wiki, so you can easily edit, discuss, and create new programs.

    Personal work

    The interactive teaching part of this course consists of only 10 lectures of 1h30. To complement this limited number of teaching hours, a signicant personal work is required from you, to prepare and digest in the best

    Before the lecture :

    • Check the suggested preparation work,
    • Read the material for the lecture (possibly print it so that you can write on it during the lecture),
    • Download the commented programs and check that they correctly run on your computer.

    The programs will be used and modified during the lectures, so it is advised that you bring your personal laptops with you.

    After the lecture :

    • Reconstruct the results produced during the course with the programs.
    • It is recommended that you bring your personal laptops with you.
    • Do the suggested exercices.

    Bibliography for the course

    • F. Charru, Instabilit�s hydrodynamiques, CNRS editions (main source for this course)
    • Drazin & Reid, Hydrodynamic instabilities
    • Schmid & Henningson, Instabilities and transition in shear flows
    • Glendinning, Instabilities, chaos & Dynamical systems
    • Godr�che & Manneville
    • Bender & Orzag, advanced mathematical methods for scientists and engineers.

    Organisation of the course

    Lecture 1 Dynamical systems I : linear systems

    Preparation work :

    • Review your knowledge about linear systems of differential equations, linear algebra, eigenvalues of real and complex matrices.
    • Advised lectures : F. Charru, sections 1.1 and 1.2; Glendining.

    Material for the lecture :

    Personal work :

    Exercices for lecture 1

    Lecture 2 Dynamical systems II : nonlinear systems

    Preparation work:

    • Charru, section 1.3 (see in particular section 11.6)

    Material for the lecture :

    Personal work :

    • Charru, exercices 1.6.4 and 1.6.5 (warning : there is an error in equations 1.62 and 8.4 ! the sign of V(A) must be changed)

    • (Charru, exercice 11.7.11) Build the bifurcation diagram of the following amplitude equation, where \mu is the bifurcation parameter: \displaystyle \frac{d x}{d t} = x (x^2+\mu)(x^2+\mu^2-1)

    • (Charru, exercice 11.7.12) Build the bifurcation diagram of the following amplitude equation, where \mu is the bifurcation parameter: \displaystyle \frac{d x}{d t} = -4 x ( (x-1)^2-\mu-1 )

    Verify that the bifurcation diagram is the one given in figure 11.20 of Charru.

    (warning : there is a sign error in the book !)

    • Charru, exercice 11.7.13. You may use the program PhasePortrait_NonLinear.m to draw phase portraits of the system for several values of the bifurcation parameter \mu.

    • Consider the dynamical system defined as follows : \displaystyle \frac{d}{d t} \left[ \begin{array}{c} x_1 \\ x_1 \end{array} \right] = \left[ \begin{array}{c} r x_1 - x_1^3 -3 x_2^2 x_1 \\ (r-1) x_2 + x_1^2 x_2 - x_2^3 \end{array} \right]

    where r is a control parameter.

    1. Using the program PhasePortrait_NonLinear.m, draw phase portraits for various values of r.

    2. Study the number of fixed points and their stability as function of r.

    3. Regroup the results under the form of a bifurcation diagram showing the amplitude parameter A= |x_1|+|x_2| of all solutions as function of r.

    Going further:

    • Charru, chapter 11
    • Glendinning

    Lecture 3 Mathematical and numerical study of a model problem displaying instabilities : the Ginsburg-Landau equation.

    Preparation work:

    Material for the lecture:

    Commented program GinsburgLandau_Linear.m

    Commented program GinsburgLandau_NonLinear_OneMode.m

    Personal work:

    • Check the theoretical results cited in the program.
    • Use the program to study the Inhomogeneous Ginsburg-Landau Equation with \sigma(x) = \sigma_0 - \sigma_2 x^2. you may start with the parameters x0 = -5, L = 10, \sigma_2 = 1 and vary the parameter \sigma_0.

    Lecture 4 Gravity-capillary waves and Rayleigh-Taylor instabilities

    Personnal work

    This chapter is left as personal work. The situation investigated corresponds to two fluids separated by a horizontal interface. The region y>0 is occupied by a fluid of density \rho_1 and the region y<0 is occupied by a fluid of density \rho_2. We note \gamma the surface tension and we neglect viscous effects.

    1. Show that modal perturbations proportional to e^{i k x - i \omega t} are governed by the following relation between k and \omega (dispersion relation) :

    \displaystyle \omega^2 = \frac{(\rho_2-\rho_1)}{\rho_2+\rho_1} g k + \frac{\gamma}{\rho_2+\rho_1} k^3.

    1. In the case where the lower fluid is heavier (\rho_2> \rho_1), show that the dispersion relation predicts surface waves with real frequencies. Recall the definition and physical interpretation of the group velocity c_g. Justify that the group velocity has a minimum value c_m corresponding to a particular wavelength \Lambda_m = 2 \pi / k_m.

    2. In the case where the lower fluid is lighter (\rho_1> \rho_2), show that the dispersion relation predicts either surface waves with real frequencies, or unstable modes with purely imaginary frequencies \omega (hence purely real eigenvalues \lambda). Show that instabilities occur for wavelength \Lambda = 2 \pi / k > 2 \pi l_c where l_c is the capillary length defined by \displaystyle l_c = \sqrt{\frac{\gamma}{|\rho_1-\rho_2| g}}

    Suggested lectures :

    • Charru, sections 2.3 and 2.4

    • Wikipedia :


    • Rayleigh-Taylor instability with confinement effects (Charru, exercice 2.8.1)
    • Rayleigh-Taylor instability of a suspended film (Charru, exercice 2.8.2)

    • Derive the dispersion relation and instability criteria for a double layer of fluid, defined as follows: \displaystyle \rho(y) = \rho_1 (y<-L/2) ; \quad \rho(y) = \rho_2 (L/22) ; \quad \rho(y) = \rho_3 (y>L/2)

    Lecture 5 Thermal instabilities (Rayleigh-Benard)

    Preparation work:

    • Review your knowledge about thermal convection (Bousinesq approximation)
    • Charru, section 2.5

    Material for the course:

    Personal work:

    • Mathematical resolution of the Rayleigh-B�nard problem for “free-free” boundaries (sources: Drazin & Reid, sections 2.8 and 2.10 ; Rieutort ; Chandrasekar)
    1. Starting from the stability equations in primitive form for [\hat u, \hat v, \hat, p, \hat T], show that the problem can be reduced to a single differential equation of order 6 with the following form (Drazin & Reid, Eq. 8.39, p. 43) : \displaystyle (\partial_y^2 - k^2) (\partial_y^2 - k^2 - \lambda) (\partial_y^2 - k^2 - \lambda/Pr)\hat v = 0
    2. precise the set of boundary conditions to be used with this equation, in the respective cases of (a) rigid boundaries, and (b) “free-slip” boundaries.

    3. In the case of “free-slip” boundaries, show that the least-stable eigenvalue is (D&R, eq. 10.5): \displaystyle \lambda = -\frac{1}{2}(1+Pr)(\pi^2 + k^2) \pm \frac{1}{2} \sqrt{ (Pr-1)^2(\pi^2+k^2)^2+k^2 \, Ra \, Pr/(\pi^2+k^2)}

    4. Show that for a given k, the flow is unstable for Ra >Ra_j(k), with

    \displaystyle Ra_j(k) = (\pi^2+k^2)^3/k^2

    (D&R, Eq. 10.6)

    Draw this curve in the (Ra,k) plane.

    1. Show that the instability threshold is Ra_c = min [Ra_j(k)] = 27 \pi^4/4 \approx 657.5 and give the corresponding value of k.
    • Study of the Lorenz system.
    1. Find all fixed-point solutions of the Lorenz system. Show that the system undergoes a supercritical Pitchfork bifurcation for r>1.

    2. Show that if if P>b+1, the stable solution existing for r>1 becomes unstable through a Hopf bifurcation for r>r_c with \displaystyle r_c = \frac{ P ^2 + b P + 3P}{P - b - 1}. Give the numerical value for b=1 (Charru) and for b=8/3 (initial choice of Lorenz).

    Tip : Study the stability of the fixed-point solution. Write the characteristic polynomial governing its stability under the form P(\lambda) = \lambda^3 + A \lambda^2 +B \lambda + C. Remarking that if a Hopf bifurcation occurs, this polynomial has two purely imaginary roots and one real root, write a condition on the coefficients A,B,C of the polynomial.

    1. Using the program, check that the bifurcation is subcritical.

    2. Using the program, verify the existence of a “strange attractor” for r>r_s with r_s\approx 24 (try to find the best approximation of r_s !)

    3. Verify that for r_c<r<r_s both stationnary and “strange” solutions are possible depending on initial conditions.

    4. Oberve the behaviour of the Lorenz system for much larger values or r. Observe in particular the behaviour in the range $r .

    Lecture 6 Shear flow instabilities I (inviscid instabilities)

    Preparation work:

    Material for the lecture:

    External material�_de_Kelvin-Helmholtz

    Personal work:

    Exercice 1: Kelvin-Helmholtz instability of a piecewise-linear shear layer (Charru, section 4.3.2; Drazin & Reid, p. 146; Schmid & Henningson).

    Consider a base flow defined as follows :

    \displaystyle U(y) = \left\{ \begin{array}{ll} U_1 & \quad (y<-\delta) \\ U_m + \Delta U y/\delta & \quad (-\delta <y<\delta) \\ U_2 & (y>\delta) \\ \end{array} \right. with U_m = (U_1+U_2)/2 and \Delta U = (U_2-U_1)/2.

    1. Starting from the solution of the Rayleigh Equation in the three zones and using suitable matching conditions, show that modal perturbations are governed by the dispersion relation with the following form : \displaystyle 4 k^2 \delta^2 (c-U_m)^2 - \Delta U^2 \left[ (2 k \delta -1)^2- e^{-4 k \delta} \right]= 0.

    2. Justify that the flow is unstable for k\delta <1 and that for k \delta \ll 1 the dispersion relation reduces to the classical one for a shear layer of zero thickness.

    3. Plot the amplification rate \omega_i = k c_i as function of k. Compare with the case of a contious profile (tanh profile).

    Exercice 2: Combined Rayleigh-Taylor-Kelvin-Helmoltz instability (Charru, exercice 4.5.1)

    Consider two superposed fluid layers with different densities and velocities :

    \displaystyle y<0 : \quad \rho = \rho_1 ; \, u = - U

    \displaystyle y>0 : \quad \rho = \rho_2 ; \, u = U

    Show that the dispersion relation governing modal perturbations is given by : \displaystyle \rho_1 (U+c)^2 + \rho_2 (U-c)^2 - \frac{(\rho_1-\rho_2)g}{k} - \gamma k = 0

    Show that this dispersion relation generalizes the cases previously considered.

    Discuss the stability conditions as function of the parameters (see Charru, p. 132)

    Exercice 3: Demonstrate the Fjortoft theorem

    Lecture 7 Shear-flow instabilities II (viscous instabilities)

    Preparation work:

    Charru, sections 5.2 and 5.3.

    Material for the lecture:

    • Viscous stability analysis of the tanh shear layer

    Commented program KH_temporal_viscous.m

    • Viscous stability analysis of the plane Poiseuille flow

    Commented program /sandbox/easystab/poiseuille_uvp.m for explaining the numerical method and validation.

    Program /sandbox/easystab/david/TS_PlanePoiseuille.m for parametric study.

    • Viscous stability analysis of the Blasius boundary layer

    (Program /sandbox/easystab/blasius.m to compute the base-flow solution for a boundary layer (Blasius equation))

    (Program /sandbox/easystab/blasiusStability.m to compute the stability properties of the Blasius boundary layer (maybe next year…)

    Personal work:

    Lecture 8 Shear-flow instabilities IV (three-dimensional perturbations and transient growth mechanisms)

    Preparation work:

    • Review the theory of dynamical systems (lectures 1 and 2)

    • Study the “green” sections in the document

    Material for the lecture:

    Personal work:

    Exercice 1 Transient growth in a nonnormal linear system.

    Sources : Rieutort, p. 228 ; Schmid & Henningson, p. 522.

    Warning : this exercice is close to the case considered by Charru (p. 31-32) but there is an error in Eq. 1.43. Will you find it ?

    We consider the following two-dimensional system :

    \displaystyle \frac{dx}{dt} = -\epsilon x

    \displaystyle \frac{dy}{dt} = - 2 \epsilon y + x.

    Or in matricial form :

    \displaystyle \frac{dX}{dt} = \left[\begin{array}{cc} -\epsilon & 0 \\ 1 & - 2 \epsilon \end{array}\right] X.

    1. Determine the eigenvalues (\lambda_1, \lambda_2) and normalised eigenvectors (\hat X_1, \hat X_2) of the system.

    2. Represent the eigenvectors (\hat X_1, \hat X_2) in the (x,y) plane. What is the angle \alpha between them ?

    3. Determine the adjoint eigenvectors (\hat X_1^\dag, \hat X_2^\dag), and represent them graphically. What is the norm of these vectors ?

    4. Give the solution x(t);y(t) with initial conditions (x_0,y_0).

    5. Among the initial conditions of ||X_0|| = \sqrt{x_0^2+y_0^2}, which one maximizes the projection upon the least unstable mode ? (this initial condition is called the optimal perturbation). Give a geometrical interpretation.

    6. Show that the energy growth associated to the optimal perturbation is given by

    \displaystyle G = \frac{x(t)^2+y(t)^2}{x_0^2+y_0^2} = e^{-2 \epsilon t} + \left(\frac{e^{-\epsilon t} - e^{-2 \epsilon t}}{\epsilon} \right)^2

    1. Show that when \epsilon \ll 1 the energy growth reaches a maximum for t_{opt} \approx ln(2)/\epsilon, with G_{max} = G(t_{opt}) \approx 1/(16 \epsilon^2).

    2. Show that for t \ll t_{opt} the energy growth is algebraic. Interpret this by comparing with the solution of the linear system with \epsilon = 0.

    Exercice 2 Subcritical transition in a stable, nonnormal system.

    Sources : Schmid & Henningson, p. 525; Trefethen et al. , Science, 1993.

    We consider the following two-dimensional system :

    \displaystyle \frac{dX}{dt} = \left[\begin{array}{cc} -\epsilon & 0 \\ 1 & - 2 \epsilon \end{array}\right] X + ||X|| \left[\begin{array}{cc} 0 & 1 \\ -1 & 0 \end{array}\right] X

    Using the program PhasePortrait_NonLinear.m, study the phase portrait of this problem.

    Show that the combined effect of transient growth and nonlinearities can lead to transition to a non-zero solution, even for very small initial conditions.

    Lecture 9 Shear-flow instabilities III (spatial an spatio-temporal theory)

    Preparation work:

    Charru, Chapter 3.

    Material for the lecture:

    Personal work:

    Lecture 10 Bluff-body Wake instabilities

    Preparation work:

    Material for the lecture:

    This final lecture is based on a paper available here:

    The programs used for this lecture are part of the StabFem project, available here.

    To get this project, you will have to do the following steps : - Install the FreeFem++ software (free, available here ) - open a terminal (linux/unix/mac systems) or a console (windows) and type the following command :

    git clone

    Personal work