Fluid instabilities and convection

Can a warmer lower layer stay at rest?

Imagine a horizontal fluid layer only 5 mm deep. Its bottom is held at 293 K and its top at 292.5 K. For a material whose density decreases when warmed, the warmer fluid lies below cooler, denser fluid. If a tiny warm patch is pushed upward, will the motion die away, remain small, or grow into convection?

Heating from below supplies a possible buoyancy feedback, not an automatic promise of visible motion. Viscosity resists velocity differences, and heat diffusion can erase the temperature contrast that drives the patch. Geometry and the wall conditions determine which disturbances fit. We will compare the feedback with those opposing effects before calculating an onset.

First follow the diagrams, signed tables and physical explanations: a displaced parcel, a disturbed interface and a heated layer. The final mode calculation is an optional return path. It explains one exactly solvable set of boundaries, not every laboratory container.

Useful published preparation is Pressure, buoyancy, and fluid statics and Geometry, trigonometry, and coordinate systems. Derivatives, Taylor series, and local approximation explains the gradients used below. For the optional mode discussion, return to Fixed points and stability and Complex numbers and complex exponentials. Temperature transport and the new symbols are defined locally; later thermodynamics lessons are not required for a first reading.

Use a fixed frame with x horizontal and z vertically upward. Gravity is −g e_z, with g positive. A prime marks a small departure from a stated base state; it does not label a different kind of fluid. Our example parameters are prescribed ideal-model inputs, not measured properties of a named substance.

Ask what a perturbation does

A base state is the flow, density, temperature and pressure about which we ask the question. It might be a resting hydrostatic layer, two flat uniform streams or a motionless conducting layer. A perturbation is a small displacement or change in these fields. Stability concerns the evolution of that departure under specified equations and conditions.

Bounded motion is not the same as decay. A frictionless displaced parcel can oscillate without returning permanently to rest. A decaying disturbance becomes smaller. A growing disturbance becomes larger. In linear theory, a mode with amplitude A(t) = A₀ exp(st) grows exponentially if the real part of its rate s is positive; its imaginary part describes oscillation.

Linearisation retains terms proportional to the small disturbance and omits products of disturbances. It can predict early growth rates and an onset condition. Once a disturbance is no longer small, the terms we omitted can change the result: motion may saturate, reorganise, overturn or develop more complicated dynamics. Extrapolating an exponential forever is not a nonlinear prediction.

A neutral mode has no exponential growth or decay in the chosen linear calculation. It need not be motionless, and a marginal double root can require special care. Showing that selected modes do not grow also does not prove nonlinear stability against every finite disturbance. State both the model and the disturbances tested.

Displace a parcel in a density gradient

Stratification means that a fluid property varies with height. Consider a resting, locally linear density profile ρ̄(z), and move a small parcel upward by η. In the present incompressible model the parcel retains its original density while pressure adjusts to its surroundings. If density decreases upward, the parcel is denser than its new surroundings, so the net buoyancy acceleration points downward.

Write G = dρ̄/dz, the background density gradient, with units kg/m⁴. For a reference density ρ₀, the linear net upward acceleration is (gG/ρ₀)η. Define N² = −(g/ρ₀)G, in s⁻². The parcel equation becomes d²η/dt² = −N²η. For N² > 0, N = √N² is the buoyancy angular frequency in rad/s, not N² itself.

Here ρ₀ = 1000 kg/m³ and g = 9.81 m/s². With G = −2 kg/m⁴, N² = +0.019620 s⁻². A 1 mm upward displacement initially gives −0.000019620 m/s²: a restoring acceleration. With G = +2 kg/m⁴, the same displacement gives the opposite acceleration and the local inviscid parcel model contains a growing solution.

For the same initial η = 1 mm and zero vertical speed, the restoring case is η(t) = η₀ cos(Nt), the uniform-density case stays at η₀, and the adverse-gradient case is η(t) = η₀ cosh(√(−N²)t). The last expression is a combination of growing and decaying exponentials chosen to satisfy the zero initial speed; it is not one pure growing mode.

These curves describe the assigned small-displacement parcel diagnostic without entrainment, viscosity or scalar diffusion. The profile is assigned on −0.1 ≤ z ≤ 0.1 m; the shown displacements remain inside it. Full internal-wave frequencies also depend on spatial structure and pressure coupling, and compressible air parcels can change density as they move. Do not apply this density-retaining parcel formula indiscriminately to the atmosphere.

A statically stable density arrangement resists this buoyancy displacement. It does not rule out shear-driven instability, and an adverse density gradient alone does not settle the onset of a confined viscous, thermally diffusive layer. Those additional dynamics matter below.

Compare gravity and capillarity at an interface

A sharp interface is a different model from a smooth density gradient. Put two immiscible fluids on either side of a flat horizontal interface: density ρ_l below and ρ_u above. Disturb the interface slightly with a wave of positive wavenumber k; its wavelength is ℓ = 2π/k. Gravity tends to restore a bump when the lower fluid is denser, but drives overturning when the upper fluid is denser. The latter mechanism is Rayleigh–Taylor instability in this gravity-loaded setting.

Surface tension σ penalises extra curved interfacial area and supplies a restoring pressure jump. Its effect strengthens more rapidly with k than the gravity term. Thus an adverse density arrangement can have growing long waves while sufficiently short waves remain oscillatory in the ideal capillary model. Saying “heavy above light” identifies a mechanism, not the outcome for every wavelength and boundary condition.

For two deep, inviscid, incompressible, irrotational fluids at rest, the linear result is ω² = S = [gk(ρ_l − ρ_u) + σk³]/(ρ_l + ρ_u). Here ω is angular frequency and S has units s⁻². If S > 0, the disturbance supports inviscid oscillatory modes with no decay in this model. If S < 0, one mode grows at rate √(−S). At S = 0 the roots coalesce; it is a marginal case, not a claim of robust damping.

Take g = 9.81 m/s², σ = 0.030 N/m and the two assigned densities 1000 and 900 kg/m³. Compare the denser fluid below with the same fluids exchanged. The table uses identical material inputs and wavelengths in both cases. It does not specify viscous stresses, finite depths, walls, phase change or nonlinear mixing.

For the exchanged, heavy-above case, the zero-shear marginal wavenumber follows from σk² = g(ρ_u − ρ_l): k = √[g(ρ_u − ρ_l)/σ]. Our inputs give approximately 180.831413 m⁻¹. The 300 m⁻¹ sample lies on the short-wave oscillatory side despite the adverse ordering. The 20 and 100 m⁻¹ samples have growing modes.

Add a velocity jump

Now keep the denser fluid below, but let the two uniform streams have different horizontal velocities. A small interface disturbance deflects these streams and changes their pressure response. Energy can be drawn from the relative flow: this is the shear-driven Kelvin–Helmholtz mechanism. A stable density ordering supplies buoyant restoration, but it does not remove the shear energy source.

Define ΔU = U_l − U_u and the density-weighted mean Ū = (ρ_lU_l + ρ_uU_u)/(ρ_l + ρ_u). For our deep inviscid interface, (ω − kŪ)² = S, where S = [gk(ρ_l − ρ_u) + σk³]/(ρ_l + ρ_u) − [ρ_lρ_u/(ρ_l + ρ_u)²](ΔU)²k². The first two terms restore this stable ordering; the last is destabilising. A common extra velocity shifts the real propagation frequency, not the growth rate.

At k = 100 m⁻¹, keep 1000 kg/m³ below 900 and σ = 0.030 N/m. Set U_l = +ΔU/2 and U_u = −ΔU/2. Their arithmetic mean is zero, but their density-weighted mean is not zero unless the densities are equal or ΔU vanishes. In the table, even the 0.500 m/s jump overwhelms the two restoring terms and gives a growing mode.

The velocity jump is an ideal vortex-sheet limit. A smooth shear layer, viscosity, finite depth, three-dimensional disturbances or a diffuse interface changes the stability problem. There is no one universal velocity-jump threshold for all flows. Neither a specified velocity difference nor a Reynolds number alone predicts the nonlinear roll-up or a turbulent mixing rate.

With positive surface tension and stable density ordering, sufficiently small k is protected by the gravity term and sufficiently large k by the capillary term for any fixed finite jump. An intermediate unstable band can exist. Our three jump samples at one k do not trace the complete band; the optional dispersion analysis tells which k values belong to it.

Identify the thermal feedback

Return to the heated layer. Conduction transfers heat through a temperature gradient without requiring bulk fluid motion. Convection transports temperature with moving fluid. Natural convection is driven by buoyancy; forced convection is driven by an imposed flow. Here there is no mean imposed flow: the question is whether the motionless conducting state admits growing disturbances.

Let h be layer height and ΔT = T_b − T_t the held bottom-to-top temperature difference. The motionless base temperature is T̄(z) = T_b − ΔT z/h. This linear profile solves the steady diffusion equation with the assigned temperatures. Pressure is hydrostatic; there is no base circulation to draw.

A rising parcel from a warmer lower position can be warmer than its new surroundings. Write its temperature departure as θ = T − T̄. For positive thermal expansion coefficient β, its buoyancy acceleration is approximately gβθ upward. Upward velocity w transports the warmer background into that position, helping sustain θ. This closes the positive thermal feedback when ΔT > 0.

Heat diffusion, with diffusivity κ in m²/s, reduces the temperature contrast; viscosity, with kinematic viscosity ν in m²/s, damps the motion. A horizontal heated layer can therefore have an adverse base density gradient while its full confined, diffusive linear problem is below onset. This is not a contradiction with the earlier undiffused parcel diagnostic.

For this small-temperature-difference Boussinesq model, use reference density ρ₀ except in buoyancy, where retain ρ′ ≈ −ρ₀βθ. The flow remains divergence-free. This approximation does not mean that all densities are physically identical or that thermal buoyancy disappears. Large density changes, compressibility or strongly varying material properties require another model.

Compare driving with two diffusions

The Rayleigh number Ra = gβΔT h³/(νκ) compares thermal buoyancy driving with viscous and thermal diffusion in this layer problem. Every factor must carry its stated units; Ra is dimensionless. The Prandtl number Pr = ν/κ compares the two diffusivities. These numbers organise the specified model, rather than supplying universal onset thresholds on their own.

Our assigned layer has h = 0.005 m, g = 9.81 m/s², β = 0.0002 K⁻¹, ν = 0.000001 m²/s and κ = 0.0000001 m²/s. Thus Pr = 10. The diffusion times across h are h²/ν = 25 s for momentum and h²/κ = 250 s for temperature. They are characteristic scales, not guaranteed waiting times for observable convection.

For the thermal cases, the bottom remains 293 K and the three top temperatures are 292.98, 292.5 and 292 K. Their Ra values are shown below. Even the largest βΔT is only 0.0002 in the prescribed linear density law, consistent with using a small density-contrast approximation for this example. This consistency check does not turn the assigned coefficients into measured data.

We next impose an exactly solvable boundary model: both boundaries are impermeable and stress-free and hold their temperatures fixed; the horizontal extent is unbounded. Periodic Fourier disturbances in that extended layer can have continuously admissible wavenumbers. A particular finite periodic box would instead restrict them. Stress-free means zero tangential viscous stress, not zero viscosity in the bulk. These are ideal boundaries, not ordinary no-slip rigid plates.

Within that model, the lowest stationary neutral threshold is Ra_c = 27π⁴/4 ≈ 657.511364. The first case is below onset; the other two admit a growing linear mode. This threshold comes from the equations and boundaries derived next. It must not be transferred to rigid walls, a heat-flux boundary, a finite box, a free-surface stress-driven problem or another fluid approximation.

Return path: derive a conditioned neutral curve

This section is optional on a first reading. Scale height by h, time by h²/κ, velocity by κ/h and temperature departure by ΔT. Write Z = z/h. For the conducting base, the linear disturbance equations are Pr⁻¹ ∂u′/∂t = −∇π′ + ∇²u′ + Ra θ e_z, ∇·u′ = 0 and ∂θ/∂t = w + ∇²θ. The symbols and derivatives in these equations are dimensionless; π′ is the correspondingly scaled pressure departure.

Use a horizontal Fourier mode proportional to exp(iαX + λt), where α = kh > 0 is dimensionless horizontal wavenumber and λ is dimensionless growth rate. The vertical velocity amplitude is W(Z), the temperature amplitude is Θ(Z), and D means d/dZ. Eliminating pressure gives (D² − α²)(D² − α² − λ/Pr)W = Ra α²Θ, while (D² − α² − λ)Θ = −W.

At each boundary, impermeability gives W = 0, stress-free motion gives D²W = 0, and fixed temperature gives Θ = 0. For stationary neutral motion λ = 0, eliminating Θ yields (D² − α²)³W = −Ra α²W. The eliminated boundary conditions may be written W = D²W = D⁴W = 0 at Z = 0 and 1.

The vertical eigenfunctions W = sin(nπZ), with positive integer n, satisfy those conditions. Substitution gives Ra_n(α) = (α² + n²π²)³/α². This is a neutral curve: at a given α, it is the stationary Ra where that mode changes sign of growth. The curve is not a plot of a parcel trajectory or a nonlinear roll amplitude.

For n = 1, differentiating the neutral curve gives its minimum at α_c = π/√2 and Ra_c = 27π⁴/4. Higher n have minima 27n⁴π⁴/4, so the first branch supplies the smallest threshold when the horizontal wavenumber is continuous. The associated physical wavelength is ℓ_c = 2πh/α_c = 2√2 h = 0.014142 m for our 5 mm layer. It is a wavelength, not the height or the width of one rising half-roll.

A rigid no-slip wall instead imposes zero tangential velocity, leading here to W = DW = 0 together with Θ = 0. The sine modes above do not satisfy that extra velocity condition. Its neutral problem must be solved with the different boundaries; we do not import a generic rigid-wall critical value into this stress-free example.

For finite containers, only particular horizontal structures may fit, so minimising over every α can be inappropriate. Rotation, composition-dependent density, surface-tension gradients or unequal scalar diffusivities introduce other mechanisms and equations. These modifications are not small licence to reuse the same curve unchanged.

Return path: distinguish growth from propagation

For the stress-free thermal mode, put Q = α² + n²π². The mode equations yield (λ + Pr Q)(λ + Q) = Pr Ra α²/Q. With Pr > 0 and Ra ≥ 0, the two roots are real. The larger root crosses zero exactly at Ra_n(α); above it there is a growing mode. Its dimensional rate is λκ/h². The neutral value is independent of Pr in this stationary boundary problem, but the dimensional growth rate is not.

At α = α_c and n = 1, our three thermal cases have dimensional larger-root rates approximately −0.054354, +0.043642 and +0.123537 s⁻¹. These describe one admitted mode, not a claim that the critical wavenumber remains the fastest-growing mode above onset. Its eventual amplitude and heat transport are not determined by this linear eigenvalue.

For the interface, write a disturbance as the real part of A exp(ikx − iωt). Under this convention a positive imaginary part of ω means growth. Decaying potential disturbances away from the interface, its material-displacement condition and the gravity/capillary pressure jump give ρ_l(ω − kU_l)² + ρ_u(ω − kU_u)² = gk(ρ_l − ρ_u) + σk³. Completing the square gives the S formula used earlier.

When S > 0, ω = kŪ ± √S describes two oscillatory branches. When S < 0, ω = kŪ ± i√(−S) contains a growing and a decaying branch. The advection frequency kŪ and the growth rate are distinct. At S = 0, the double root has no exponential growth, but a general initial-value solution can contain a term linear in time; do not relabel it guaranteed decay.

For stable density ordering and σ > 0, S < 0 is equivalent to (ΔU)² > [(ρ_l + ρ_u)/(ρ_lρ_u)][g(ρ_l − ρ_u)/k + σk]. Minimising the right-hand side over positive k gives its least jump at k = √[g(ρ_l − ρ_u)/σ]. This is a threshold for some mode of this deep vortex-sheet model, not a criterion for an arbitrary smooth shear profile.

Separate onset from nonlinear outcomes

At the onset stage we ask whether a small disturbance can extract enough energy from an adverse density arrangement, relative stream motion or maintained thermal contrast to overcome the restoring or dissipative effects in the model. These are related stability questions, but their energy sources and conditions differ.

If thermal convection develops, rising and falling fluid also modify the temperature distribution. That alters the feedback which started the motion. A growing interface can change its shape and the shear distribution; viscosity and mixing may then matter. A linear growth rate alone cannot choose the eventual rolls, plumes, interface roll-up or transport rate.

Convection is not a synonym for turbulence. Organised or steady convective motion is possible; more complex regimes require the actual nonlinear dynamics and boundaries. Conversely, forced fluid motion can transport heat without a thermal instability. The following permanent subject, Turbulence, addresses that further regime; it is not a preparation requirement here.

Layered oceans and atmosphere, heated equipment and experiments offer applications, but each must be modelled with its own compressibility, material properties, boundaries, forcing and relevant additional physics. Our thin ideal layer and deep two-fluid examples are mechanisms to reason with, not forecasts for an ocean, a storm or a cooling device.

Check the most common misreadings

“There is motion, so the state is unstable” confuses a disturbance with its evolution. A bounded wave can move indefinitely without exponential growth. “Stable means every disturbance decays” excludes the undamped stable modes present in our inviscid examples.

“Warm fluid is below cool fluid, so convection must start immediately” ignores diffusion, viscosity and admitted spatial modes. “Ra exceeds 657.5, so any container convects” transfers a stress-free, fixed-temperature threshold beyond its boundary model.

“Denser below means every flow is stable” ignores the energy in relative streams. “Denser above means every wavelength grows” ignores capillary restoration. “A big ΔU predicts mixing” turns a conditioned linear onset into a nonlinear transport forecast.

A reliable conclusion names the base state, the perturbation, the model/domain, the boundary conditions, the signed restoring or driving terms, and the time/amplitude range in which the calculation is valid.

Worked examples

Try it yourself

1. A displaced inviscid parcel oscillates at constant amplitude. Is it decaying? Does its motion alone demonstrate exponential instability?

   Show the answer and reasoning

   No to both. It is bounded oscillatory motion in the stated model.

   A disturbance can move without becoming larger or smaller. Damping would be additional physics; exponential growth concerns amplitude evolution rather than the mere presence of velocity.

2. Reverse the parcel gradient from −2 to +2 kg/m⁴ while retaining an upward 0.001 m displacement, ρ₀ = 1000 kg/m³ and g = 9.81 m/s². Find the initial acceleration and its direction.

   Show the answer and reasoning

   +0.000019620 m/s², upward and reinforcing the displacement.

   N² becomes −0.019620 s⁻², and −N²η is positive. This local undiffused parcel instability is not by itself a thermal-layer onset calculation.

3. For the heated-layer Ra definition, double h while holding g, β, ΔT, ν and κ fixed. By what factor does Ra change? Does this establish a universal container threshold?

   Show the answer and reasoning

   Ra becomes eight times larger; no universal threshold follows.

   Ra is proportional to h³, so 2³ = 8. Stability still depends on boundaries, geometry, admissible modes and the governing approximation.

4. Can Ra_c = 27π⁴/4 be used unchanged when the two stress-free boundaries are replaced by ordinary rigid no-slip plates? Explain.

   Show the answer and reasoning

   No. The velocity boundary conditions and the mode problem change.

   The sine solution uses W = D²W = 0. Rigid no-slip plates require W = DW = 0 with the fixed-temperature condition. Reusing a threshold without re-solving those conditions is not justified.

5. An assigned mode with dimensionless amplitude begins at 0.001 and follows A = A₀ exp(0.2t), with t in seconds. What is A at 5 s? If the linear estimate is admitted only up to amplitude 0.01, when does it reach that limit?

   Show the answer and reasoning

   A(5 s) ≈ 0.002718; the estimate reaches 0.01 at approximately 11.512925 s.

   Substitute 0.2 × 5 = 1. For the limit solve exp(0.2t) = 10, giving t = ln(10)/0.2. The limit is a stipulated validity cutoff, not a predicted nonlinear saturation amplitude or turbulent transition.

6. Add the same horizontal velocity to both interface streams without changing their densities, surface tension or k. What happens to the velocity jump and linear growth rate?

   Show the answer and reasoning

   The jump and growth rate are unchanged; the real advection frequency shifts.

   ΔU is a difference. The S expression depends on ΔU, not a common added velocity, whereas ω contains kŪ. This is a frame/propagation distinction within the given uniform-stream model.

What to carry forward

• Stability is about a perturbation evolving from a named base state; bounded oscillation, decay and growth are different outcomes.
• A density-decreasing-upward parcel profile supplies buoyant restoration in the local undiffused Boussinesq diagnostic, not a guarantee against every shear or finite-amplitude disturbance.
• At an ideal interface, gravity and surface tension compete with relative-flow energy; density order and wavelength both matter.
• Thermal convection needs buoyancy feedback to compete with viscosity and temperature diffusion under the actual boundaries.
• The stress-free fixed-temperature neutral curve and its 657.511364 threshold are conditioned derivations, not universal convection rules.
• A linear growth rate predicts early amplification, not final mixing, heat transport or turbulence.