Astrophysical Fluid Dynamics course of the University of Cambridge (Option lecture)

A fluid particle follows a Lagrangian trajectory, so its position satisfies $\dot{x}=u(x,t)$ rather than the velocity evaluated permanently at its initial position. Assume $\omega \ne 0$, and writeSet $\alpha =k{X}_0$ and $\theta =\alpha -\omega t$. A Taylor expansion of the velocity about $X_0$ givesThus $\dot{a}=\sin \theta$, andAt the next order,so the required second-order trajectory in a sinusoidal velocity field isBoth corrections vanish at $t=0$. All retained terms have units of length. The secular contribution gives the Stokes driftalthough the Eulerian time average of the given velocity is zero. The particle samples the wave at its displaced position, producing this second-order drift. The expansion is in the dimensionless parameter $|k{u}_0/\omega |\ll 1$, for fixed times; it ceases to be uniform when the accumulated drift changes the wave phase by order unity.

Interpret “reach” first as remaining an energetic, supersonic blast wave, rather than as the arrival of any acoustic wave. An atmospheric blast-wave range estimate must account for the thin atmosphere: the distance $D=8.5\times 10^6\mathrm{m}$ is much greater than its height $H=10^4\mathrm{m}$.

For an optimistic estimate, suppose the deposited energy $E$ remains in a horizontally expanding atmospheric layer. If its representative mass density is ${\rho }_a$, the swept mass is of order $\pi D^{2}H{\rho }_a$. Conservation of energy gives the self-similar blast wave scalingUsing ordinary atmospheric estimates ${\rho }_a\sim 1\mathrm{kg}{\mathrm{m}}^{-3}$ and speed of sound $c_a\sim 300\mathrm{m}{\mathrm{s}}^{-1}$, this givesThese are scaling estimates, not an exact shock wave solution. The strong shock wave approximation would already have failed atEquivalently, keeping the disturbance strongly supersonic out to $D$ would require substantially more than the sonic-transition energy scalewhich exceeds $10^{22}\mathrm{J}$ by about a factor of twenty. Reasonable atmospheric averages and order-unity geometrical coefficients do not reverse the comparison. Escape of energy through the top of the atmosphere makes the fully confined estimate optimistic; Earth's curvature changes its area estimate by an order-unity factor.

The strong atmospheric blast will not survive to Lansing in this model. A weakened acoustic wave can nevertheless travel there. The supplied energy and atmospheric height do not determine its detectable amplitude, detailed attenuation or damage: those require an atmospheric propagation model. Thus “no disturbance of any kind reaches Lansing” would be an unjustified stronger conclusion.

This is magnetic flux freezing in ideal magnetohydrodynamics. In addition to the supplied ideal magnetohydrodynamic induction equation, a physical magnetic field obeys Gauss's law for magnetism, $\nabla \cdot \mathbf{B}=0$.

Parameterize the material surface by two particle labels $a,b$, so $d\mathbf{S}=({\partial }_a\mathbf{x}\times {\partial }_b\mathbf{x})dadb$ and $D\mathbf{x}/Dt=\mathbf{u}$. Differentiating these tangents gives the material surface element evolution lawwhere $(\nabla \mathbf{u}{)}_{ij}={\partial }_j{u}_i$. Thus the material magnetic-flux transport theorem isExpand the curl in the ideal magnetohydrodynamic induction equation:Since $D/Dt={\partial }_t+\mathbf{u}\cdot \nabla$ is the material derivative, the first three terms cancel the corresponding terms in the surface transport formula. The remaining integrand is $\mathbf{u}\nabla \cdot \mathbf{B}=0$. ThereforeNo incompressible flow assumption is needed: the changing surface area cancels the compression term. This proof assumes smooth velocity and magnetic field and a surface carried with the fluid; the ideal magnetohydrodynamic induction equation by itself, for an arbitrary non-solenoidal vector field, would leave the final divergence term.

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.

The first argument mixes a material fluid element with a fixed spatial volume. For a material volume element following the fluid, its volume is generally time dependent:Therefore the material derivative cannot be taken outside the product with $\delta V$ as in the proposed cancellation. Instead, mass conservation statesThe momentum balance for a material fluid element, using Newton's second law, is consequentlyFor a fixed control volume, fluid crosses its boundary and carries momentum; a momentum flux term must be included. Its local conservative form is ${\partial }_t(\rho \mathbf{u})+\nabla \cdot (\rho \mathbf{u}\otimes \mathbf{u})=-\nabla p$, equivalent to the displayed Euler equations after using the continuity equation. The erroneous equation omits the term $\rho \mathbf{u}\nabla \cdot \mathbf{u}$ needed to convert a material derivative of momentum density into this conservative form.

For the second task, take a homogeneous fluid with constant mass density $\rho$, prescribed gravitational acceleration $\mathbf{g}$, incompressible flow and inviscid flow. Its equations areUnder the replacement $\rho \mapsto \tilde{\rho }=C\rho$, the same velocity solves these equations if $p\mapsto \tilde{p}=Cp$, up to a spatially constant pressure offset. This is density scaling invariance of incompressible Euler flow:The comparison assumes the same geometry and compatible velocity or gravitational free surface conditions; an independently specified nonuniform boundary pressure would supply additional information. The gravitational field is taken as prescribed. If the fluid's own self-gravity must be included, its field changes with mass density, so the invariance need not hold. A single freely moving particle's density-independent gravitational acceleration does not justify neglecting that distinction.

A highly conducting solar wind carries the Sun's magnetic field outward by magnetic flux freezing. This can open magnetic field lines that would form closed loops in a static dipole model. Beyond an opening radius $r_o$, approximate the field as radial. Constancy of the magnetic flux in a fixed solid-angle flux tube gives the radial solar-wind magnetic flux lawwhich decays more slowly than the vacuum magnetic dipole field, proportional to $r^{-3}$.

For example, if the field becomes radial already near the solar radius $R_{\odot }$, normalizing both estimates to the same field there givesThus a dipole estimate $10^{-11}\mathrm{T}$ becomes of order $2\times 10^{-9}\mathrm{T}$, comparable to the stated hundredfold enhancement. More generally, matching dipolar decay inside $r_o$ to radial decay outside givesThis is a plausible scaling explanation, not a determination of the actual opening radius. The magnetic field need not be outward everywhere: opposite-polarity sectors can preserve zero net magnetic flux through a closed surface, consistent with Gauss's law for magnetism. The argument concerns the changed geometry of the field in a flowing plasma, rather than changing the dipole exponent by assumption.

Take the outward wavenumber $k>0$; the absolute value in the source makes the two propagation directions equivalent. The phase velocity of the gravity-capillary wave obeysDifferentiate with respect to $k$:The second derivative is positive, and $c_p\to \infty$ at both ends of the positive-$k$ range, so this is the global minimum. The minimum phase speed of a gravity-capillary wave isFor $g=9.81\mathrm{m}{\mathrm{s}}^{-2}$ and the supplied surface tension and mass density,The corresponding wavenumber is $k_{\ast }\simeq 369{\mathrm{m}}^{-1}$, giving wavelength ${\lambda }_{\ast }=2\pi /k_{\ast }\simeq 1.70\mathrm{cm}$. At this minimum, $c_g=d\omega /dk=c_p+kd{c}_p/dk=c_p$. Thus the group velocity and phase velocity agree for this particular ripple, although they differ for general gravity-capillary waves.