Past exam of the astrophysical fluid dynamics course of the University of Cambridge / 2022 / Section B

The buoyancy frequency, also called the Brunt–Väisälä frequency, is the angular frequency of a small vertical displacement restored by buoyancy in a stably stratified fluid. Take $z$ upward, gravitational acceleration $-g\hat{\mathbf{z}}$, and a background in hydrostatic equilibrium,Displace a small gas element upward by $\xi$, sufficiently slowly for its pressure to match its surroundings but without heat exchange. For this adiabatic process, $p/{\rho }^{\gamma }$ is constant along the element's motion, so its mass density changes by $\rho {\gamma }^{-1}(d\ln p/dz)\xi$. The surrounding mass density changes by $(d\rho /dz)\xi$. The element's excess mass density is consequentlyIts acceleration is $\ddot{\xi }=-g\Delta \rho /\rho =-N^2\xi$. ThusFor a gas of fixed composition obeying the ideal gas law, $\rho \propto p/T$. The adiabatic buoyancy frequency of an ideal gas is thereforewhere $c_p=\gamma \mathcal{R}/(\gamma -1)$ is the specific heat capacity at constant pressure, and $\mathcal{R}=p/(\rho T)$ is the specific gas constant. Here $N^2$ has units ${\mathrm{s}}^{-2}$; $N$ is an angular frequency, not a cyclic frequency. Positive $N^2$ gives simple harmonic motion, negative $N^2$ gives convective instability, and $N^2=0$ is marginal. The temperature form assumes fixed composition; the preceding density form is more general.

Write the background mass density as $\bar{\rho }(z)={\rho }_0+{\rho }_0^{\prime }z$ and its hydrostatic pressure as $\bar{p}(z)$, whereApply an Eulerian perturbation of a fluid variable to the background fields:Thus $\pi$ is the Eulerian pressure perturbation, in pascals; it is not the background hydrostatic pressure.

For incompressible flow, mass conservation becomes $D{\rho }_{\mathrm{tot}}/Dt=0$. Linearizing about rest givesThe Euler equations are ${\rho }_{\mathrm{tot}}D\mathbf{u}/Dt=-\nabla p_{\mathrm{tot}}-{\rho }_{\mathrm{tot}}g\hat{\mathbf{z}}$. The background pressure gradient and weight cancel. Dropping products of perturbations givesUnder $|{\rho }_0^{\prime }z|\ll {\rho }_0$, replace $\bar{\rho }$ by ${\rho }_0$ in inertia but retain the background mass density gradient in the displacement-induced buoyancy. This is the Boussinesq approximation. Its linear Boussinesq equations here areThe absolute value makes explicit that the small-density-variation condition concerns magnitude, including the usual negative ${\rho }_0^{\prime }$ of stable density stratification.

Use plane wave perturbations proportional to $e^{i(k_{x}x+k_{z}z-\omega t)}$, with amplitudes $\hat{u},\hat{w},\hat{\pi },\hat{\varrho }$. The linear equations becomeTo eliminate pressure without dividing by $k_x$, multiply the vertical momentum equation by $k_x$ and subtract $k_z$ times the horizontal momentum equation:Use incompressible flow and the equation for mass density to obtain ${\omega }^2{\rho }_0(k_x^2+k_z^2)\hat{w}=-g{\rho }_0^{\prime }{k}_x^2\hat{w}$. Hence, for a nontrivial internal gravity wave,For stable density stratification, ${\rho }_0^{\prime }<0$ and $N=\sqrt{-g{\rho }_0^{\prime }/{\rho }_0}$ is the buoyancy frequency. Then $0\le {\omega }^2\le N^2$ and the normal modes oscillate rather than grow. For ${\rho }_0^{\prime }>0$, a mode with $k_x\ne 0$ has ${\omega }^2<0$ and grows exponentially: denser fluid lies above lighter fluid. Uniform density, ${\rho }_0^{\prime }=0$, is marginal for this restoring mechanism.

The case $k_x=0$, $k_z\ne 0$, has zero frequency in the displayed relation. Incompressible flow then requires $\hat{w}=0$, so it is a degenerate non-propagating case rather than a finite-frequency vertical oscillation. These equations concern the incompressible model of part (b); its $N^2$ does not contain the adiabatic compression correction needed in part (a).

For $N^2>0$, write $K=\sqrt{k_x^2+k_z^2}$ and choose either signed branch $\omega =sN{k}_x/K$, with $s=\pm 1$. These two branches cover the roots of the squared dispersion relation. The group velocity of an internal gravity wave isIf instead the positive-frequency root is chosen for both signs of $k_x$, set $s=+1$ and multiply the displayed components by $\mathrm{sgn}(k_x)$, away from the cusp $k_x=0$.

Taking the dot product with the wave vector givesThus the group velocity, which describes the travel of a narrow wave packet and its energy, is perpendicular to the wave vector, and hence perpendicular to the phase velocity ${\mathbf{c}}_p=\omega \mathbf{k}/K^2$. This is the internal-wave phase and group velocity geometry. It also follows from phase-group orthogonality for degree-zero dispersion: the frequency depends only on the direction of $\mathbf{k}$, not its magnitude. For $k_z=0$, the group velocity is zero, consistent with the stationary extremum $|\omega |=N$.

The Rankine-Hugoniot conditions follow by integrating conservation laws across a thin moving layer, then letting its macroscopic thickness tend to zero. For example, normal momentum conservation in one dimension has the conservative formFor a shock wave moving at speed $V$, integrate over a control interval whose endpoints move with the discontinuity. The Reynolds transport theorem giveswhere $[f]=f_R-f_L$ denotes the difference of limiting states. Mass conservation similarly gives $V[\rho ]=[\rho u]$. If $v=u-V$ is the normal velocity in the shock frame, combining the two relations yieldsThe same construction applied to conservation of energy, for a non-radiating normal shock wave with nonzero mass flux, gives $[h+v^2/2]=0$, where $h$ is the specific enthalpy. An equation of state closes these relations, and the physical compression branch has increasing entropy.

Finite dynamic viscosity replaces the ideal discontinuity by a thin smooth transition. A viscous stress contributes a derivative to the momentum flux inside the layer and converts directed kinetic energy into internal energy. If its gradients vanish in the uniform far states, the same limiting Rankine-Hugoniot conditions still hold there. Balancing nonlinear transport with viscous diffusion gives a thickness of order $\nu /|u_L-u_R|$ in the viscous Burgers equation, where $\nu$ is the kinematic viscosity.

The pressure-free equation has the scalar conservation law formA smooth travelling ansatz $u=U(x-Vt)$ obeys $(U-V)U^{\prime }=0$, so it cannot continuously connect two distinct constant states. The inviscid shock wave must be interpreted as a discontinuous weak solution. Its Rankine-Hugoniot condition isFor $u_1\ne u_2$,The result is independent of which state is called upstream. For a physical compressive entropy shock of the Inviscid Burgers equation, the left and right states satisfy $u_L>u_R$ and therefore $u_L>V>u_R$: characteristic curves run into the shock wave from both sides. The ordering $u_L<u_R$ instead produces a rarefaction wave.

For $\nu >0$ and positive $\Psi$, set $f=\log \Psi$, so the Cole-Hopf transformation is $u=-2\nu f_x$. Substitution into the viscous Burgers equation givesHence the most general equation implied by this substitution isThe time normalization in the Cole-Hopf transformation accounts for the arbitrary function $C(t)$: multiplying $\Psi$ by a nonzero function of time leaves ${\partial }_x\log \Psi$ and therefore $u$ unchanged. Replacing $\Psi$ by $e^{-{\int }^{t}C(s)ds}\Psi$ sets $C=0$. With this conventional normalization,the heat equation.

For each $j=1,2$, the function ${\Psi }_j=e^{\nu {\kappa }_j^{2}t+{\kappa }_{j}x}$ obeys ${\partial }_t{\Psi }_j=\nu {\kappa }_j^2{\Psi }_j=\nu {\partial }_x^2{\Psi }_j$. By the superposition principle for this linear partial differential equation,is a solution, positive for all real $x,t,{\kappa }_1,{\kappa }_2$. Here ${\kappa }_j$ has units ${\mathrm{m}}^{-1}$ so both exponents are dimensionless. The positivity ensures that its logarithm and the resulting velocity are regular.

Let $a_j=\nu {\kappa }_j^{2}t+{\kappa }_{j}x$. Differentiating the logarithm givesFor definiteness order ${\kappa }_1<{\kappa }_2$. Then $e^{a_1}$ dominates as $x\to -\infty$, while $e^{a_2}$ dominates as $x\to +\infty$. The limiting states areWithout this ordering convention, $u(-\infty )=-2\nu \min ({\kappa }_1,{\kappa }_2)$ and $u(+\infty )=-2\nu \max ({\kappa }_1,{\kappa }_2)$. Thus the states $u_{1,2}$ correspond to $-2\nu {\kappa }_{1,2}$ after matching each label to its dominant exponential; the propagation speed does not depend on the naming of the two states.

Their exponent difference isConsequently, the Viscous Burgers travelling shock moves atSet $\Delta u=u_L-u_R=2\nu ({\kappa }_2-{\kappa }_1)>0$ and $\xi =x-Vt$. Factoring out the mean exponential gives the hyperbolic tangent profileA precise measure of the width must be stated. The natural length in the hyperbolic tangent argument isThus the layer has thickness of order $\nu /\Delta u$. A different common convention, jump divided by maximum gradient, gives $\delta =\Delta u/|u_x(0)|=8\nu /\Delta u=2\ell$; the far-tail exponential e-folding length is $\ell /2$. These differ only by their stated definitions, rather than by different solutions.

Two distinct limits require ${\kappa }_1\ne {\kappa }_2$. The source permits arbitrary real constants, but in the degenerate case ${\kappa }_1={\kappa }_2$ the solution is uniform, $u=-2\nu {\kappa }_1$, and there is no shock wave or finite shock width.

Treat the cylinder as infinitely long, with a vacuum exterior and no surface tension. Let ${\Phi }_0$ be its gravitational potential. Integrating the gravitational Poisson equation over a coaxial cylinder of radius $r<R$ and length $L$, then applying the divergence theorem, givesThus the inward gravitational acceleration is $-2\pi G\rho r\hat{\mathbf{r}}$. Hydrostatic equilibrium requiresImpose zero pressure on the free surface, $p_0(R)=0$, and integrate inward:This is the hydrostatic self-gravitating incompressible cylinder. Its central pressure is $\pi G{\rho }^2{R}^2$, and the units of the expression are $\mathrm{Pa}$. A nonzero external pressure would add that constant to $p_0(r)$.

Write total pressure and gravitational potential asHere $\pi$ is the Eulerian perturbation of a fluid variable applied to pressure, while $\psi$ is the corresponding gravitational potential perturbation. They have units $\mathrm{Pa}$ and ${\mathrm{m}}^2{\mathrm{s}}^{-2}$, respectively. The constant interior mass density is $\rho$.

The Euler equations areSubtract the hydrostatic equilibrium relation $\nabla p_0/\rho +\nabla {\Phi }_0=0$ and discard terms quadratic in the perturbation. Mass conservation gives incompressible flow, and the gravitational Poisson equation has no bulk density perturbation. Thus, strictly inside the fluid,The last equation is Laplace equation, not a statement that self-gravity has disappeared. The displaced surface carries a first-order surface mass density $\rho \eta$, where $\eta$ is its radial displacement, and this supplies the source and matching conditions for $\psi$ at the moving interface. The axisymmetric gravitational potential of a displaced cylindrical interface, expressed using modified Bessel functions, incorporates that boundary source.

It is useful to define $W=\psi +\pi /\rho$. Taking the divergence of the momentum equation gives ${\nabla }^{2}W=0$. For an axisymmetric Fourier mode with axial wavenumber $k>0$, regularity at $r=0$ selects radial dependence $I_0(kr)$, the Modified Bessel function of the first kind, since $K_0(kr)$ is singular there. In what follows $k=|k_z|$ can always be chosen nonnegative.

Use complex amplitudes and take their real parts at the end. Write the radial surface displacement as $\eta =ϵe^{i(kz-\omega t)}$ and the combined potential asThe linear kinematic boundary condition is $u_r(R)={\partial }_t\eta$. Differentiate it once more in time and use the radial momentum equation, ${\partial }_t{u}_r=-{\partial }_{r}W$, evaluated at the unperturbed radius $R$:With $\kappa =kR$, the mode amplitudes therefore satisfyThe prime differentiates the Modified Bessel function of the first kind with respect to its dimensionless argument, rather than with respect to $r$. The result has units $\mathrm{m}{\mathrm{s}}^{-2}$, as required because $ϵ{\Pi }_0$ is a potential. It applies to $k>0$; the uniform $k=0$ limit requires separate treatment.

The vanishing total pressure must be imposed at the displaced free surface, not at fixed $r=R$. The first-order Taylor expansion isSince $p_0(R)=0$ and $p_0^{\prime }(R)=-2\pi G{\rho }^{2}R$, the Lagrangian pressure perturbation condition isThe supplied gravitational perturbation and the combined potential, evaluated at $R$, then giveHere $K_0$ is the Modified Bessel function of the second kind. Equating mode amplitudes and using part (c) yieldsRearrange to obtain the dispersion relation for the axisymmetric gravitational instability of an incompressible cylinder:Every quantity on the right is dimensionless, consistent with $G\rho$ having units ${\mathrm{s}}^{-2}$. The stabilizing term $1/2$ comes from displacing the background hydrostatic pressure profile; omitting the Lagrangian pressure perturbation would lose it.

Define the dimensionless right-hand side of part (d) by $F(\kappa )$. Then ${\omega }^2=-4\pi G\rho F(\kappa )$. Because $I_0^{\prime }(\kappa )=I_1(\kappa )>0$ and $I_0(\kappa )>0$ for $\kappa >0$, its sign is the sign of $I_0(\kappa )K_0(\kappa )-1/2$. The original plot is positive for sufficiently long wavelengths and negative at short wavelengths; the nonzero crossing is approximately ${\kappa }_c\simeq 1.07$. Numerically solving the displayed equation gives ${\kappa }_c\simeq 1.0667$.

For $0<\kappa <{\kappa }_c$, ${\omega }^2<0$. Write $\omega =i\sigma$ with $\sigma =\sqrt{4\pi G\rho F}>0$. The time factor $e^{-i\omega t}$ then contains an exponentially growing solution. For $\kappa >{\kappa }_c$, ${\omega }^2>0$ and the ideal inviscid modes oscillate stably, without damping. At the crossing the frequency is zero and the mode is marginal. In terms of the axial wavelength ${\lambda }_z=2\pi /k$,The growth rate tends to zero as $kR\to 0^{+}$, so the very longest wavelengths are not the fastest growing. The source plot has its maximum near $kR\simeq 0.58$, corresponding to ${\lambda }_z\simeq 11R$ and $\sigma \simeq 0.25\sqrt{4\pi G\rho }$. Thus a characteristic fragmentation spacing is several cylinder diameters; a quantitative nonlinear fragment mass is beyond this linear stability calculation.

For unstable modes the physical displacement must be interpreted as $\mathrm{Re}[ϵe^{i(kz-\omega t)}]$ with complex amplitudes handled consistently throughout. The paper's real-part notation is harmless for the squared dispersion relation, but it should not be read as imposing real $\omega$ in a gravitational instability calculation. The exactly uniform mode $k=0$ is excluded from the finite-$k$ derivation; the limit of the curve alone does not establish a separate volume-preserving uniform-radius mode.