Past exam of the mathematics course of the University of Cambridge / 2018 / iii / Paper 345 (Environmental Fluid Dynamics)

Let $p$ denote the fluid pressure perturbation divided by the reference mass density, and let $b=-g{\rho }^{\prime }/{\rho }_0$ be the buoyancy perturbation. The buoyancy frequency satisfies $N^2=-g{\hat{\rho }}_z/{\rho }_0$. Linearization requires small boundary slope $k|{\eta }_0|\ll 1$, small displacements compared with the vertical scales of the disturbance and background, and advection small compared with the oscillatory acceleration. For a propagating internal gravity wave this includes $|m{\eta }_0|\ll 1$; boundary slope alone is insufficient when $|m|\gg k$. The Boussinesq approximation also requires small relative mass density differences over the region of interest. These conditions must hold for the resulting disturbance, including any amplification by resonance.

The inviscid Linearized Boussinesq equations and linear kinematic boundary condition areWrite the velocity field as $(u,w)=\mathrm{Re}\{(U(z),W(z))e^{i(kx-\omega t)}\}$. Incompressibility and horizontal momentum balance giveThe vertical momentum balance therefore yieldsFor a real vertical wavenumber $m$, this gives the dispersion relationWhen $0<\omega <N$, define $\ell =k\sqrt{N^2/{\omega }^2-1}$. The radiation condition selects $m=-\ell$, because energy must leave the boundary upward. The boundary-forced internal gravity wave isIts wavevector is $\mathbf{K}=(k,m)$, its phase velocity is ${\mathbf{c}}_p=\omega \mathbf{K}/(k^2+m^2)$, and its group velocity isThus the phase velocity points down and right, while the group velocity points up and right. They are perpendicular. Internal-wave polarization makes the particle motion an oscillation along the group velocity direction, with zero first-order mean transport. The constant-phase lines of an internal gravity wave are also parallel to the group velocity. Their inclination $\vartheta$ above the horizontal satisfies $\sin \vartheta =\omega /N$.

When $\omega >N$, set $\kappa =k\sqrt{1-N^2/{\omega }^2}$. Boundedness at infinity selects the evanescent waveThe horizontal and vertical velocity components are in quadrature: fluid particles describe small ellipses, and the response decays over ${\kappa }^{-1}$. There is no upward time-averaged energy flux, because $P$ and $W$ are in quadrature. A real vertical group velocity is not defined for this evanescent wave. The pattern travels horizontally with phase velocity $\omega /k$.

For $N=0$ the response is the evanescent wave with $\kappa =k$; the ratios $\omega /N$ should not be used. At the cutoff $\omega =N>0$, the bounded harmonic solution has $W=-i\omega {\eta }_0$, $U=0$: it neither decays nor has nonzero upward group velocity. It is the limiting cutoff response, rather than a localized radiating disturbance. A uniform nonzero buoyancy frequency in an infinitely deep Boussinesq approximation is itself a local idealization of the background mass density.

The arrows distinguish the phase velocity, group velocity, and oscillatory particle motion of the internal gravity wave. The right panel shows the decay envelope and particle ellipses of the evanescent wave.

At the change in buoyancy frequency, continuity of vertical velocity and fluid pressure givesThe background mass density is continuous, so there is no interfacial density jump or separate surface restoring force. Since $P=i\omega W^{\prime }/k^2$, these internal-wave transmission across a stratification step conditions are continuity of $W$ and $W^{\prime }$. The buoyancy perturbation need not be continuous, because the background density gradient jumps.

Let $W_0=-i\omega {\eta }_0$ and $\ell =k\sqrt{N_2^2/{\omega }^2-1}$. The lower unstratified layer obeys $W^{\prime \prime }-k^{2}W=0$ and the upper radiation condition givesPropagating these matching data downward givesThe kinematic boundary condition at the moving lower boundary therefore givesThe upper vertical-displacement complex amplitude at $H$ isIn particular, $|{\eta }_H/{\eta }_0|=[{\cosh }^2(kH)+(\ell /k{)}^2{\sinh }^2(kH){]}^{-1/2}$. The unstratified layer attenuates and phase-shifts the transmitted internal gravity wave. The matching height is $H$ in the original PDF; the TeX transcription's $H^{\prime }$ is a typo.

The central stratified layer supports upward and downward internal gravity waves, whose counterpropagating components form a standing wave. The unstratified regions support evanescent waves, and the upper disturbance decays as $e^{-k(z-2H)}$. Consequently the central layer is a stratified internal-wave guide, rather than a source of propagating energy at infinity.

The largest response occurs at its trapped normal mode frequencies. They can be specified without solving the forced problem. Set the boundary forcing to zero to find a free normal mode. The lower solution is proportional to $\sinh (kz)$, so continuity of fluid pressure and vertical velocity givesIn the central layer writeThe upper Robin boundary condition then gives the exact trapped-mode conditionThe left side increases strictly from zero to infinity, so there is one positive root for each $n$; the frequencies accumulate at zero. This is constructive phase matching after reflection at both ends. The phase shifts from the evanescent waves matter: simply imposing integer half-wavelengths across the stratified layer is generally incorrect.

In the ideal inviscid model, exact resonant forcing has no bounded steady harmonic solution. The undamped normal mode grows secularly under sustained forcing. Weak viscosity or other losses would produce large finite peaks near the displayed frequencies; the linear approximation eventually fails if the disturbance becomes too large.

The illustrated free normal mode has an oscillatory central region and evanescent wave tails. Its lower tail reaches the fixed zero-displacement boundary; its upper tail decays to infinity.

The exact mixture equation of state is $\rho =(1-\varphi ){\rho }_f+\varphi {\rho }_p$. Keeping first-order terms in the small particle volume fraction and thermal expansion givesThe omitted term is ${\rho }_0\varphi \theta$. The downward solid-volume particle deposition flux is $J_v=W_s\varphi$; the corresponding particle mass flux is $J_m={\rho }_p{W}_s\varphi$. Uniform vertical mixing and constant layer depth implyFor $W_s>0$, set $r=W_s/H$. Thenand henceThe heated particle-laden layer loses its stable density stratification when this contrast vanishes. For $\beta >0$, solving for the neutral-buoyancy time givesThe layer is neutral at $T_s$ and statically unstable for $t>T_s$. The subsequent uniformly mixed lower-layer solution cannot represent the resulting overturning.

If $W_s=0$, then $\varphi ={\varphi }_0$, $\theta =\beta {\varphi }_{0}t$, and the continuous limit is $T_s=\gamma /\beta$. If $\beta =0$, the layer stays denser than its surroundings at every finite time and becomes neutral only asymptotically when $W_s>0$; thus $T_s=\infty$. These limiting cases must replace the printed formula when its denominator vanishes.

Use the shallow-water approximation, a hydrostatic vertically mixed layer, a deep motionless ambient, and no ambient fluid entrainment. Particle loss changes reduced gravity but changes total layer volume only at the discarded dilute-particle order. Likewise thermal expansion changes the equation of state while the leading Boussinesq approximation retains volume conservation. DefineDepth integration of mass conservation, horizontal momentum balance, particle transport, and the heating law givesHere the integrated excess hydrostatic pressure is $b{h}^2/2$. There is no particle source from the bed. Bed drag and the particle contribution to inertia are neglected at this order. The equivalent variable-buoyancy shallow water equations areThe term $h{b}_x/2$ is necessary when heating or particle loss makes the reduced gravity vary horizontally.

For $v=(h,u,\varphi ,\theta {)}^T$, the quasilinear system $v_t+A(v)v_x=S(v)$ hasIts characteristic polynomial factors asThus the characteristic curves have slopesOn the two repeated contact characteristic curves, the particle and heating laws are $D\varphi =-W_s\varphi /h$ and $D\theta =\beta \varphi$. Writing $c_g=\sqrt{bh}$ and $D_{\pm }={\partial }_t+(u\pm c_g){\partial }_x$, the gravity-wave compatibility relations areFor $h>0$ and $b>0$, gravity-wave right eigenvectors can be taken as $(1,\pm c_g/h,0,0{)}^T$. The contact eigenspace has dimension two: $\delta u=0$ and $b\delta h+(h/2)\delta b=0$, with independent $\delta \varphi ,\delta \theta$. Therefore $A$ has a full set of real eigenvectors.

The system is hyperbolic throughout strict static stability, but not strictly hyperbolic because the contact speed is repeated. At $b=0$ the eigenvalues coalesce and this diagonalizability is lost; when $b<0$ the gravity-wave speeds become complex, reflecting the loss of a statically stable lower layer.

Let $\mathcal{V}=H{L}_0$ be the layer volume per unit transverse width. The gravity-current box model imposes $hL=\mathcal{V}$; it approximates the front and scalar balances, rather than an exactly uniform solution of the local momentum equation. For a specified constant front Froude number $\mathrm{Fr}$, the heated particle-laden gravity current satisfiesThe initial data are $L=L_0$, $\varphi ={\varphi }_0$, and $\theta =0$. The gravity-current front condition applies only while $\gamma \varphi -\theta \ge 0$.

For heating without settling, $\varphi ={\varphi }_0$ and $\theta =\beta {\varphi }_{0}t$. Integration of the front equation givesHence the heated gravity-current runout isIt is reached at the neutral-buoyancy time $\gamma /\beta$ within this front model. A finite-momentum current could subsequently coast or lift off: the formula is the maximum predicted by the prescribed gravity-current front condition, which contains no independent front inertia.

For settling without heating, $\theta =0$. Set $K=\mathrm{Fr}\sqrt{g\gamma \mathcal{V}}$, so $\dot{L}=K\sqrt{\varphi }{L}^{-1/2}$. Eliminating time yieldsThus the runout length of a gravity current isThis length is approached as $t\to \infty$, because positive concentration cannot disappear at a finite time under the settling law. If distance means displacement from the removed barrier, the answer is $L_{\max }-L_0$ in either case.

It is important to fix the mass flux convention. With literal mass units and the corresponding density-weighted buoyancy flux, the Boussinesq approximation and circular top-hat plume model giveReplacing $\rho$ by ${\rho }_0$ in inertial fluxes is the Boussinesq approximation; the density difference is retained in the buoyancy flux. Define the kinematic plume fluxesIf the symbols $Q,\mathbf{M},F$ instead denote the customary kinematic fluxes, simply omit the factors ${\rho }_0$ in all conversions below. Keeping the two conventions consistent is essential for the jet length and dimensional prefactors.

The Batchelor entrainment hypothesis prescribes inward edge velocityThe constant entrainment coefficient relates ambient inflow to the local axial velocity. A centreline segment of arc length $ds$ entrains volume $2\pi b{u}_{e}ds$. Hence, setting $E=2\alpha \sqrt{\pi }$,The buoyancy flux is conserved because the homogeneous entrained ambient has zero density deficit, and mixing has no buoyancy source or sink. The imposed zero source mass flux with nonzero source momentum flux is a singular point-source idealization, not a finite-radius nozzle with zero emitted fluid. Indeed $V=m/q$ and the plume reduced gravity $f/q$ diverge as $q\to 0$ with nonzero source fluxes. The Boussinesq approximation cannot remain valid arbitrarily close to that point. The top-hat plume model is applied outside the unresolved source region.

For a neutrally buoyant horizontal turbulent round jet, a segment's horizontal momentum flux changes only if an external horizontal force or ambient momentum input exists. Neither is present, soPositive buoyancy adds a vertical force, not a horizontal one. Therefore $d{M}_x/ds=0$ remains true for an inclined forced plume even when $F_0\ne 0$; in general $M_x=M_0\cos {\theta }_0$.

The net upward force on a segment is $g({\rho }_0-\rho )\pi b^{2}ds$. Since $F=g({\rho }_0-\rho )\pi b^{2}V$ and $V=|\mathbf{M}|/Q$, vertical momentum balance givesIn kinematic plume fluxes, the equations are $m_x^{\prime }=0$, $m_z^{\prime }=fq/m$, where $m=(m_x^2+m_z^2{)}^{1/2}$. The centreline geometry satisfiesFor $F_0>0$, $Q>0$ beyond the source, and an initially nonvertical turbulent round jet pointing right, $M_z$ increases and the centreline turns toward the vertical. The magnitude of momentum flux is not conserved once the tangent has a vertical component.

Put $m_0=M_0/{\rho }_0$ and $f_0=F_0/{\rho }_0>0$ under the literal mass convention. Their dimensions are $[m_0]=L^4{T}^{-2}$ and $[f_0]=L^4{T}^{-3}$. Dimensional analysis givesThe jet length compares source momentum with buoyancy. The region near the source behaves as a turbulent round jet; sufficiently far away, the forced plume behaves like a pure plume. With the constant entrainment coefficient retained explicitly, substantial turning occurs over a length of order ${\ell }_J=L_J/\sqrt{E}$, where $E=2\alpha \sqrt{\pi }$. Statements involving $s\ll L_J$ or $s\gg L_J$ normally hold this dimensionless coefficient fixed.

Near the point source, $m\sim m_0$ and the Batchelor entrainment hypothesis gives $q\sim E\sqrt{m_0}s$. ThereforeIntegrating the centreline tangent givesTo leading order the centreline is a straight turbulent round jet at its source inclination. For the horizontal source, the first curvature is the horizontal forced-plume trajectoryFor a vertical source, $\cos {\theta }_0=0$, the trajectory stays vertical rather than developing this cubic transverse displacement. If $F_0=0$, there is no finite buoyancy-induced jet length and the turbulent round jet remains straight.

For a horizontal source, $m_x=m_0$ throughout. Far downstream $m_z\gg m_0$, so the kinematic plume fluxes obey the vertical pure plume equations to leading order,Set $q=q_{\ast }{f}_0^{1/3}{\hat{s}}^{5/3}$ and $m_z=m_{\ast }{f}_0^{2/3}{\hat{s}}^{4/3}$ with $\hat{s}=s+s_0$. Matching exponents and coefficients gives $(5/3)q_{\ast }=E\sqrt{m_{\ast }}$ and $(4/3)m_{\ast }^2=q_{\ast }$. Consequentlyand the requested literal mass and density-weighted buoyancy fluxes areto leading order. The radius is $b\sim (6\alpha /5)\hat{s}$ and the axial velocity decreases as ${\hat{s}}^{-1/3}$.

The centreline becomes nearly vertical but retains finite horizontal drift:Since $dz/ds=1-O({\hat{s}}^{-8/3})$, there is a finite positive geometric deficit $D={\lim }_{s\to \infty }(s-z)$, and $z=s-D+o(1)$. The equivalent vertical plume virtual origin is therefore located atFor the horizontal point-source forced plume, the virtual source is downstream and below the actual source. The horizontal turbulent round jet first entrains substantial ambient fluid while rising only a little. Thus at the height where it turns upward it already has a finite radius and volume flux; an equivalent pure plume must have begun rising from below to acquire these. The turning displacement, entrained radius divided by its far-field spreading angle, and virtual vertical depth all scale asAt fixed entrainment coefficient this is simply $O(L_J)$. The arclength origin $s=-s_0$ should not be interpreted as a physical height: the curved near field supplies the additional shift $D$. The numerical constants depend on the chosen top-hat plume model and entrainment coefficient; scaling does not fix them.

The inclined forced plume bends from a cubic near-field trajectory toward a vertical asymptote. Dashed far-field spreading lines extrapolate to the equivalent plume virtual origin; the marked location is obtained from the integral model for this illustration, not a universal numerical prediction.

Let $\Delta \rho ={\rho }_p-{\rho }_f>0$. A spherical grain has submerged weightNeglect lift and contact torque, and use a sliding static friction model with normal reaction $R=G$. The inertial quadratic drag iswhere the last equality fixes the convention $\tau ={\rho }_f{U}^2$. Equivalently $U$ is the bed shear velocity and $C_D$ is an effective drag coefficient referred to it. A literal grain-level flow speed can differ from the shear velocity; that conversion must then be absorbed into $C_D$.

Downstream sliding begins when $F_D>{\mu }_{s}R$. Defining the Shields parameter by $\Theta =\tau /(\Delta \rho gd)$, the threshold force balance isThe grain moves downstream above this threshold in the stated sliding model. Real grain motion can instead involve lift, rolling, irregular contacts, or viscous drag; the printed constant belongs to the particular inertial-drag convention and force balance above.

The normal contact force and resisting static friction balance the drag force and submerged weight at impending motion. The right panel uses locally bed-tangent drag force.

Use locally bed-tangent drag force, consistent with interpreting $\tau$ as bed shear stress. Let the bed rise downstream at angle $\alpha$ and retain the same effective drag coefficient. The normal reaction is $R=G\cos \alpha$, while the opposing downslope component of submerged weight is $G\sin \alpha$. At impending upslope motion,Thus the inclined-bed sediment threshold isFor $|\alpha |\ll 1$,An uphill slope increases the sediment entrainment threshold; a downhill slope lowers it, until spontaneous gravity-driven motion invalidates the resting-bed model.

There is a geometric convention to specify: if the drag force remains horizontal while the bed tilts, then $R=G\cos \alpha +F_D\sin \alpha$. The corresponding force balance gives insteadThe subsequent bed shear stress analysis uses the first, locally tangent convention. The two interpretations should not be mixed.

The Exner equation is solid-volume mass conservation: a downstream increase in sediment transport flux removes material and lowers the bed. The factor ${\varphi }_b$ converts bed-height change into solid-volume change. The saturation length equation models downstream adjustment of actual sediment transport to the equilibrium sediment flux. Grain acceleration and entrainment/deposition need a finite distance, so $q$ lags $q_{\mathrm{sat}}$. The last relation is an empirical transport law above the sediment entrainment threshold, with $\chi >0$ and normally $\gamma >0$; below threshold it must be interpreted with a positive part rather than raising negative excess stress to an arbitrary power. It is a constitutive closure, not a consequence of mass conservation.

Take the physical parameters $L_{\mathrm{sat}}>0$, ${\varphi }_b>0$, $\chi >0$ and $\gamma >0$. Linearize about a horizontal, uniformly transporting bed with $\Delta ={\tau }_0-{\tau }_{\mathrm{th},0}>0$ and $q_0={\varphi }_b\chi {\Delta }^{\gamma }$. Set $\lambda =\sigma -ikc$ and use the real parts ofThe locally planar approximation uses the small instantaneous bed slope to compute the inclined-bed sediment threshold. Since $\alpha ={\eta }_x+O({\eta }_x^3)$,A sinusoidal disturbance has both uphill and downhill slopes: the printed $0<\alpha \ll 1$ describes the uphill derivation, and its first-order continuation applies to both signs. The zero-slope reference is the one implicit in the printed decomposition with constant term ${\tau }_{\mathrm{th},0}$; a finite mean inclined bed would require a shifted base threshold.

Define $C=\chi \gamma {\Delta }^{\gamma -1}$. The linearized transport equations areThe linearization requires $|{\eta }_x|\ll 1$ and $|\hat{\tau }-{\hat{\tau }}_{\mathrm{th}}|\ll \Delta$.

A fluid-dynamical shear-response closure is missing from the printed question. The three sediment equations and local slope correction do not determine $\hat{\tau }$ from $\hat{\eta }$. Write the general bed shear response as $\hat{\tau }=\mathcal{T}(k)\hat{\eta }$. A frequently used scale-invariant closure for $k>0$ isHere $A$ is the component in phase with bed height and $B$ represents an upstream phase lead in bed shear stress. They require an independent flow model and can depend on $k$. Introducing them explicitly makes the linear stability analysis complete conditional on a specified flow response; it does not turn them into data supplied by the question. An example of this hydrodynamic closure is given in Fourrière, Claudin and Andreotti's bedform-instability analysis.

Eliminate the transport amplitudes in the linear stability analysis. For a general bed shear response $\mathcal{T}(k)$, the dispersion relation isThis is the most specific result available without adding the missing hydrodynamic closure. With $\mathcal{T}(k)={\tau }_{0}k(A+iB)$, setMultiplication by $1-ix$ gives the growth rate and migration velocityThe packed-bed volume fraction cancels. For a real physical closure the negative wavenumber response is the complex conjugate; the displayed $k>0$ formulas cover the independent normal modes.

If $A>0$, $B>0$ are treated as constants, the upstream shear stress phase lead promotes growth, while the inclined-bed sediment threshold reduces it and the finite saturation length delays sediment transport. There is an unstable band precisely when $b>0$: $0<k<b/(a{L}_{\mathrm{sat}})$. Within this band $c>0$, so the growing bedforms migrate downstream. The fastest normal mode satisfiesThis last maximization assumes constant $A,B$; if the hydrodynamic coefficients vary with wavenumber, their derivatives also enter the selection condition. When $b\le 0$ and $a>0$, all nonzero normal modes decay.

Just above the sediment entrainment threshold,Thus the common regime $B<1/{\mu }_s$ is stable near onset: the slope penalty overwhelms the hydrodynamic phase lead. This conclusion is conditional on that inequality, not universal. For constant $B>0$, instability begins only when ${\tau }_0>{\tau }_{\mathrm{th},0}/({\mu }_{s}B)$ as well as ${\tau }_0>{\tau }_{\mathrm{th},0}$. The transport sensitivity $C=\chi \gamma {\Delta }^{\gamma -1}$ vanishes at onset for $\gamma >1$, is finite for $\gamma =1$, and is singular for $0<\gamma <1$. Exactly at threshold, the positive-part transport law needs separate treatment; the above linearization with perturbations small compared with $\Delta$ is unavailable.

Far above the sediment entrainment threshold, $b\simeq {\tau }_{0}B$. If $A,B>0$, the unstable cutoff tends to $B/(A{L}_{\mathrm{sat}})$, and the dominant wavelength is set by the saturation length. Both growth rate and migration velocity scale with $C{\tau }_0\sim \chi \gamma {\tau }_0^{\gamma }$, apart from their length factors. Their numerical values and any detailed dependence on wavenumber still require a flow closure.

For example, prescribing uniform shear stress independently of the bed gives $\mathcal{T}=0$ andThis equally admissible closure illustrates why instability cannot be asserted from the printed sediment equations alone.

The plots use constant $a>0$ and three illustrative ratios $b/a$. They show how the sign of the shear phase lead minus the slope correction determines stability; they are not numerical predictions for unspecified flow conditions.