Internally heated convection (2-D)

A mantle heated from within has no hot lower boundary layer. Heat is produced everywhere, leaves only through the surface, and the flow is driven by cold sheets that sink from the top boundary layer, while the return flow rises broadly and slowly. This example runs that case in the spherical annulus at three Rayleigh numbers, with an insulating core so that every watt produced inside the shell has to cross the surface.

Ra = 10⁴

Ra = 10⁵

Ra = 10⁶

Temperature (300 to 900 K) at the end of each run, 9.5 Gyr (0.3 thermal diffusion times), Nr = 64 (64 × 512 cells). The cold downwellings get thinner and more numerous as the Rayleigh number rises; there are no plumes, because nothing heats the base.

Set-up

shell	$R_{\mathrm{in}}=1220$ km, $R_{\mathrm{out}}=2220$ km ($f=0.5495$, $D=1000$ km), spherical annulus
boundaries	free slip; 300 K at the surface; insulating core-mantle boundary (cfg.cmbHeatfluxMode = "insulating")
fluid	Boussinesq, constant viscosity set by eta_for_Ra with a temperature scale $\mathrm{\Delta }T=1000$ K
heating	uniform, $H^{\prime }=H{D}^2/(c_p\kappa \mathrm{\Delta }T)=3.4,6.6,14.0$ at Ra $={10}^4,{10}^5,{10}^6$, i.e. $H=4.25\times {10}^{-12},8.25\times {10}^{-12},1.75\times {10}^{-11}$ W/kg
initial state	a linear ramp from 1300 K to 300 K plus a 20 K perturbation of wavenumber $m=6$ and 1 K of random noise
run	0.3 thermal diffusion times, $0.3D^2/\kappa \approx 9.5$ Gyr

The heating rates and the shell are those of Hernlund & Tackley (2008, Table 3), chosen there so that each case carries the surface heat flow of the basally heated case at the same Rayleigh number. The Rayleigh number uses the temperature scale, not a temperature drop the model holds, so the internal-heating Rayleigh number is $\mathrm{Ra}{H}^{\prime }$.

Running it

julia --project=. -t 2 examples/annulus_internal_heating_2d.jl 1e6 64

Arguments are Ra Nr m (defaults 1e5 64 6); Nt = 8 Nr. Wall time on one desktop machine at Nr = 64: 30 min at Ra = 10⁴, 2.5 h at 10⁵, 5.7 h at 10⁶. At the end the driver prints the mean nondimensional temperature $⟨T{⟩}^{\prime }=(⟨T⟩-T_{\mathrm{s}})/\mathrm{\Delta }T$ over the $r^2$ measure and $V_{\mathrm{rms}}$ in units of $\kappa /D$, averaged over the last third of diagnostics.csv.

What to expect

Ra = 10⁶ at 262, 1812, 4575 and 9484 Myr. The initial ramp is hot at depth and cools first: the sheets that sink from the top pool over the insulating core (262 Myr), and with no heat flowing in from below the pool warms only by the heat produced inside it. Once it has, the shell settles into a time-dependent state of sinking sheets and broad, slow upwelling.

At statistical steady state the surface heat flow equals the heat produced in the shell, $q_{\mathrm{surf}}=H{\rho }_0(R_{\mathrm{out}}^3-R_{\mathrm{in}}^3)/(3R_{\mathrm{out}}^2)$, and the energy budget in diagnostics.csv closes to a few parts in ${10}^4$:

Ra	$⟨T{⟩}^{\prime }$	$⟨T⟩$	$V_{\mathrm{rms}}$ ($\kappa /D$)	$q_{\mathrm{surf}}$ (mW/m²)	heat produced (mW/m²)	state
${10}^4$	0.278	578 K	19.8	8.67	8.66	steady
${10}^5$	0.351	651 K	69.4	16.73	16.80	weakly time dependent (68 to 70)
${10}^6$	0.341	641 K	251.5	35.57	35.64	time dependent (190 to 360)

The interior temperature rises steeply from Ra = 10⁴ to 10⁵ and then levels off, because a thinner top boundary layer removes the same heat with a smaller temperature drop. The velocity grows roughly as ${\mathrm{Ra}}^{0.55}$. The same cases against the published spherical-annulus values are listed on Benchmarks and verification.