Past exam of the mathematics course of the University of Cambridge / 2018 / iii / Paper 321 (Dynamics of Astrophysical Discs)

The specific orbital energy released in a black-hole accretion flow reaching an inner radius $r_{\mathrm{in}}$ is of order $GM/r_{\mathrm{in}}$. For a circular Keplerian orbit, the binding energy is more precisely $GM/(2r_{\mathrm{in}})$, with the other half of the potential-energy decrease retained as orbital kinetic energy. Multiplying by the mass accretion rate gives the accretion luminosity scaleFor an inner radius of order the Schwarzschild radius, $r_{\mathrm{in}}\sim 2GM/c^2$, this becomesThis is an energy-release scaling, not unit radiative efficiency of black-hole accretion. The actual inner orbit, stresses and advected energy determine the coefficient and the emitted luminosity; the innermost stable circular orbit of a Schwarzschild black hole lies outside its horizon.

Take $v>0$ to mean the magnitude of the inward radial velocity. In the Gaussian units used for the magnetic pressure, approximate the surface density of a disk by $\Sigma \simeq 2H\rho$. The mass accretion rate and free-fall speed then giveAt the magnetospheric truncation radius, balance the magnetic pressure of the magnetic dipole field with the specified ram pressure:ThusDiscarding the convention-dependent numerical factor givesThe magnetic pressure grows inward as $r^{-6}$, faster than the ram pressure $\propto r^{-5/2}$, so the field disrupts the flow inside this radius. The numerical coefficient depends on the vertical averaging and on treating the inflow as radial free fall, as stipulated here; a detailed thin disk boundary model need not have the same coefficient.

The corotation radius is where the stellar spin matches Keplerian rotation:A parcel forced into corotation at $r_m$ has outward centrifugal acceleration ${\Omega }_m^2{r}_m$. This exceeds the inward gravitational acceleration whenThis is the onset condition for a centrifugal propeller effect. It does not by itself guarantee escape to infinity. For example, a parcel released with only azimuthal speed ${\Omega }_m{r}_m$ has nonnegative specific orbital energy only ifAdditional magnetic work, radial motion or pressure can alter that escape condition.

The magnetic torque transfers angular momentum from the faster rotator to the slower one. When ${\Omega }_m>{\Omega }_K(r_m)$, the neutron star supplies angular momentum to the inner accretion disk and possible propeller flow, so it tends to spin down. When the inner accretion disk rotates faster than the star, the magnetic connection and accreted angular momentum tend to spin it up.

For slowly varying mass accretion rate and magnetic dipole moment, these tendencies favour a magnetic spin equilibrium with the corotation radius near the magnetospheric truncation radius:The precise balance depends on accreted and ejected angular momentum and on the extent of magnetic coupling across the accretion disk; corotation is a scaling estimate, not an exact zero-torque theorem.

Use $\dot{M}>0$ for inward mass accretion rate, and $T_m>0$ for magnetic torque adding angular momentum to the disk. Combining the continuity equation with conservation of angular momentum identifies the inward mass flux as $(dh/dr{)}^{-1}d\mathcal{G}/dr$. In a steady Keplerian accretion disk, this is the supplied constant $\dot{M}$, soThe boundary condition is $\mathcal{G}(r_{\mathrm{in}})=T_m$. Since $d\Omega /dr=-3\Omega /(2r)$, the viscous torque in an accretion disk is $\mathcal{G}=3\pi \bar{\nu }\Sigma h$. Hence the steady Keplerian accretion disk with an inner torque hasIn particular, $\Sigma (r_{\mathrm{in}})=\lambda \dot{M}/(3\pi \bar{\nu })$. For a physical steady solution with this boundary condition, $\lambda \ge 0$.

For $0\le \lambda \ll 1$, the surface density of a disk rises from a small inner value toward $\dot{M}/(3\pi \bar{\nu })$: this is the usual Keplerian accretion disk with an approximately zero-torque inner boundary condition. For $\lambda \gg 1$, the inner surface density of a disk is large and decreases approximately as $r^{-1/2}$, giving a torque-dominated accretion disk. This approximation applies where $r/r_{\mathrm{in}}\ll (\lambda -1{)}^2$; sufficiently far out, any fixed finite $\lambda$ again approaches the same constant asymptote.

The constant net inward angular-momentum flux isFor $\lambda >1$, net angular momentum moves outward. The strong-torque profile is therefore decretion disk-like in its stress-dominated structure, but the specified mass flux remains inward: a genuine net decretion disk requires an outward mass flux, not merely large $\lambda$.

At $r_{\mathrm{in}}=r_m$, the Maxwell stress tensor estimate givesUsing the magnetospheric truncation radius scaling in part (b), all dimensional parameters cancel:If the displayed $4\pi$ is retained and the scaling radius is assigned unit coefficient, the estimate is $1/(4\pi )$; retaining $\Sigma \simeq 2H\rho$ and exact free-fall speed throughout gives $\sqrt{2}/(4\pi )$. Those prefactors are not controlled by the scaling argument. In particular, it does not generate a parametrically large $\lambda \gg 1$.

Distinguish the cylindrical radius $R$ from the spherical radius $r=\sqrt{R^2+z^2}$. A Taylor expansion of the central gravitational potential for $|z|/R\ll 1$ givesThus the usual thin-disk expansion of the stellar gravitational potential isThe vertical acceleration is $-{\partial }_z\Phi =-{\Omega }_K^{2}z$. The first term can be omitted when considering only vertical forces at fixed $R$, but must be retained for radial dynamics.

Take $\gamma >1$ for the ideal gas. At fixed cylindrical radius, put $\mathcal{R}=k/(\mu m_p)$, where $k$ is the Boltzmann constant, $\mu$ the mean molecular weight and $m_p$ the proton mass. Also write $a=9\alpha \Omega /4$. Vertical hydrostatic equilibrium, thermal equilibrium and the ideal gas law requireConsequently, with $K=A/a$,Integrating from the midplane yieldsThe polytropic vertical structure in stellar gravity therefore has, for $|z|\le H$,whereThis is a polytropic equation of state $P\propto {\rho }^{1+1/\beta }$ with polytropic index $\beta$. For $|z|>H$, set $\rho =P=0$; temperature has no physical gas value in that vacuum. Here $H$ is the finite semi-thickness, rather than an exponential disk scale height.

The vertical acceleration in the thin disk is of order ${\Omega }^{2}H$. Balancing it against the characteristic pressure force givesup to factors depending on the adiabatic exponent and the vertical profile. A displacement of size $\delta z$ has a restoring acceleration of order ${\Omega }^2\delta z$, so $\delta \ddot{z}\sim -{\Omega }^2\delta z$. Equivalently the sound-crossing time is $H/c_s$. Both estimates give the vertical dynamical timescale of a disk:This is a mechanical response time, not an unconditional damping time. An inviscid stable layer can exhibit a vertical breathing mode of an astrophysical disk without settling; its adiabatic frequency satisfies ${\omega }_H^2=(\gamma +1){\Omega }^2$. Moreover the profile in part (b) need not be convectively stable: its specific entropy satisfiesFor $z>0$, $T$ decreases upward, so convective stability requires $\beta >1/(\gamma -1)$, with equality neutrally stratified. Thus interpreting the timescale as re-establishment of a stable equilibrium presupposes suitable stability and damping; the stated equations alone do not guarantee either.

For slow cooling of an astrophysical disk, let $t_{\mathrm{cool}}\gg {\Omega }^{-1}$ and assume a vertically symmetric laminar column with no significant radial mass exchange. Its vertical velocity is of order $H/t_{\mathrm{cool}}$, so the vertical material derivative of that velocity is of order $H/t_{\mathrm{cool}}^2$. Relative to gravity or pressure acceleration $\sim {\Omega }^{2}H$, the inertial correction is $O[(\Omega t_{\mathrm{cool}}{)}^{-2}]$.

The hydrostatic approximation therefore holds to leading order while the column slowly cools and contracts. Non-turbulent evolution removes the prescribed turbulent heating, but not the compressional work in the ideal gas energy equation:The continuity equation must still determine the slow vertical flow; $D_t$ is not generally ${\partial }_t$ at fixed height. This is a quasi-static approximation for a well-prepared, mechanically stable column, with fast free oscillations or growing convection excluded from the assumed slow solution.

At fixed cylindrical radius, let $w=u_z$ and use the continuity equation ${\partial }_t\rho +{\partial }_z(\rho w)=0$. For the similarity solution $z=\xi \eta (t)$ and $\rho =\tilde{\rho }(\xi )/\eta$, the continuity equation reduces toHere $\eta$ is dimensionless and $\xi$ labels the initial height. Midplane symmetry sets the integration constant to zero, so the flow is a homologous vertical motion of an astrophysical disk:The hydrostatic approximation and ideal gas law becomeThe material derivatives are $D_{t}P=\tilde{P}\dot{\eta }$ and $D_t\rho =-\tilde{\rho }\dot{\eta }/{\eta }^2$, so the energy equation givesSeparation of variables requires a positive constant $K$ such thatIntegrating with $\eta (0)=1$ gives the slow cooling of an astrophysical disk similarity factorA closed set of profile equations is thereforeThese have the same polytropic vertical structure in stellar gravity as part (b). Writingwith $P_0=\mathcal{R}{\rho }_0{T}_0$, fixesReturning to physical height, the explicit profiles and evolving semi-thickness for $|z|\le H(t)$ areandOutside this moving surface $\rho =P=0$. The surface density of a disk is conserved because $\rho dz=\tilde{\rho }d\xi$. The column cools and contracts without changing its scaled profile: its midplane temperature falls as ${\eta }^2$, its mass density rises as ${\eta }^{-1}$ and its pressure falls as $\eta$.

The value $T_0$ determines $H_0$, but the mass density normalization also enters $C$. If the initial column is specifically the heated equilibrium of part (b), switching off its heating givesOtherwise $P_0$, or equivalently the conserved column mass, is additional initial data. The contraction timescale grows with $1+Ct$, so sufficiently slow initial evolution remains slow relative to the fixed vertical dynamical timescale of a disk. This is an exact solution of the reduced hydrostatic approximation, rather than the full momentum equation: its omitted vertical acceleration is $\ddot{\eta }\xi$, with $\ddot{\eta }/(\eta {\Omega }^2)=2(1+\beta )C^2/[(1+2\beta {)}^2{\Omega }^2(1+Ct{)}^2]$, which remains small and decreases if $C\ll \Omega$.

Let $\mathbf{\zeta }=\nabla \times \mathbf{u}$ and $\mathbf{\omega }=\mathbf{\zeta }+2\Omega {\mathbf{e}}_z$ be the relative vorticity and absolute vorticity. The nonlinear acceleration identity rewrites the momentum equation asTaking the curl removes the pressure, kinetic energy gradient and shearing-sheet tidal potential. ThusBoth $\nabla \cdot \mathbf{u}$ and $\nabla \cdot \mathbf{\omega }$ vanish. The vorticity equation is consequentlyFor the absolute-vorticity flux tensor $T_{ij}=u_i{\omega }_j-{\omega }_i{u}_j$, use the flux convention $(\nabla \cdot \mathbf{T}{)}_j={\partial }_i{T}_{ij}$, with the first index denoting transport direction. This givesThe scalar $\theta$ is a buoyancy displacement variable, so dimensional consistency gives it units of length; the printed potential-temperature terminology presupposes this normalization. Its curl source directly forces only the $y,z$ components of absolute vorticity. Hence radial vorticity conservation in a shearing sheet takes the local flux formIts volume integral is constant when the net boundary flux vanishes. However $D_t{\omega }_x=\mathbf{\omega }\cdot \nabla u_x$ still permits vortex stretching; the radial component is not generally a materially conserved scalar in three dimensions. Using divergence on the second index of the stated tensor would reverse the transport sign, so the tensor convention matters.

For the Keplerian shearing sheet background $\mathbf{U}=-3\Omega x{\mathbf{e}}_y/2$, the only spatial derivative is ${\partial }_x{U}_y=-3\Omega /2$. Since $U_x=0$ and all fields are independent of $y$, $\mathbf{U}\cdot \nabla \mathbf{U}=0$ and $\nabla \cdot \mathbf{U}=0$.

The Coriolis acceleration exactly balances the tidal term:With constant pressure and $\theta =0$, the buoyancy displacement variable equation also has zero advection, zero source $U_x$ and zero diffusion. Thereforeis a steady solution. Its absolute vorticity is ${\mathbf{\omega }}_0=(\Omega /2){\mathbf{e}}_z$, since the background relative vorticity is $-3\Omega {\mathbf{e}}_z/2$.

For an axisymmetric vertical mode of a shearing sheet in this Boussinesq approximation, the perturbations have no horizontal gradients, so background-shear advection vanishes. The linear Coriolis acceleration and perturbation advection of the background shear givewhile thermal conduction and incompressibility giveHere $\beta$ is the thermal diffusion rate of a Fourier mode, not the cooling exponent in Question 2. Since $k\ne 0$, incompressibility forces $u_z^{\prime }=0$, and the vertical momentum equation then forces $P^{\prime }=0$ for this nonzero-wavenumber amplitude.

The three remaining amplitudes satisfyA nontrivial solution requires the determinant to vanish, giving the convective overstability dispersion relationNo division by $s$ or $s+\beta$ was needed, so the thermal energy mode of a shearing sheet is retained. The notation $N^2$ denotes a signed buoyancy frequency squared and can be negative for an adverse stratification; it must not be restricted to the square of a real positive frequency.

Treat $n^2=N^2/{\Omega }^2$ as a signed small parameter, with $|n^2|\ll 1$, as printed in the original PDF. At $n^2=0$, the dispersion relation factors asThus the two epicyclic motion modes have $s_0=\pm i\Omega$, while the thermal energy mode of a shearing sheet hasIt decays on the thermal diffusion time ${\beta }^{-1}$. For completeness, its next correction is $s_E=-\beta +\beta N^2/({\beta }^2+{\Omega }^2)+O(\Omega n^4)$ at fixed $\beta /\Omega$.

For either epicyclic motion root, the implicit function theorem applied to the dispersion relation givesHenceSince $\beta >0$, the oscillations undergo overstability precisely whenThe convective overstability growth rate tends to zero both for very slow diffusion and for very fast diffusion. Differentiating $\beta /({\beta }^2+{\Omega }^2)$ shows that its maximum occurs at $\beta =\Omega$, corresponding to $k^2=\Omega /\xi$. ThereforeThis is the growth rate of the perturbation amplitude; a quadratic perturbation energy grows at twice that rate. The convective overstability arises when an adverse radial specific entropy gradient couples to epicyclic motion with a finite thermal conduction lag.