Deep Magma Oceans: Dynamics & Evolution

• Deep magma oceans are extensive molten silicate layers that regulate planetary thermal histories and drive dynamic, episodic tidal dissipation.
• They are modeled using shallow-water formalisms and eigenmode analyses to capture complex fluid dynamics and multimodal responses.
• Observationally, these oceans manifest in variable light curves and thermal emissions, offering insights into interior-exterior coupling.

Deep magma oceans are extensive regions of subterranean or surface molten silicate, characterized by high melt fractions and dynamic fluid behavior, driven by tidal and radiative heating. These oceans are prevalent in the early evolution of rocky planets and moons, particularly in the aftermath of giant impacts and among short-period exoplanets subjected to extreme insolation and sustained orbital eccentricity. Deep magma oceans regulate planetary heat transport, control short- and long-term thermal histories, and fundamentally alter spin, orbital, and geophysical evolution, with observational signatures in planetary light curves, thermal emission, and volatile loss (Farhat et al., 2024, Farhat et al., 11 Jan 2026).

1. Governing Physical Framework

The dynamics of deep magma oceans are governed by the interplay of gravitational and tidal potentials, rotational effects, viscous damping, and thermodynamic gradients. For tidally forced magma oceans on synchronously rotating exoplanets or early Earth, the star’s eccentricity tide is expressed as: $$U_T(\theta, \lambda, t) = \mathrm{Re} \left\{ \sum_{p=2}^\infty \sum_{q=0}^p \sum_{k=-\infty}^{\infty} U_{p,q}^{T;k} P_p^q(\cos \theta) e^{i(q \lambda + \sigma_q^k t)} \right\}$$ where $U_{p,q}^{T;k}$ incorporates planetary and stellar parameters, and $\sigma_q^k$ specifies forcing frequencies (Farhat et al., 11 Jan 2026).

Within a shallow-water (Laplace-Tidal) formalism, the behavior of the ocean is modeled by: $$\partial_t u + f \hat{z} \times u + \sigma_R u = g \nabla (\zeta_\text{eq} - \zeta)$$

$\partial_t \zeta + \nabla \cdot (H u) = 0$

where $u$ is horizontal velocity, $f$ the Coriolis parameter, $\sigma_R$ the Rayleigh drag rate, $H$ the magma ocean depth, $g$ gravity, and $U_{p,q}^{T;k}$0 the actual/modeled surface displacement (Farhat et al., 11 Jan 2026).

Imposing a coastline at the day-night terminator yields mixed Dirichlet-Neumann boundary conditions:

• $U_{p,q}^{T;k}$1 at $U_{p,q}^{T;k}$2
• $U_{p,q}^{T;k}$3 for the scalar potential $U_{p,q}^{T;k}$4
• $U_{p,q}^{T;k}$5 for the toroidal potential $U_{p,q}^{T;k}$6

These boundary conditions are critical for reproducing observable hemispheric variability.

2. Normal Modes and Multimodal Response

The dayside magma ocean supports a discrete spectrum of horizontal eigenmodes calculated by solving: $U_{p,q}^{T;k}$7

$U_{p,q}^{T;k}$8

for eigenvalues $U_{p,q}^{T;k}$9 parameterized for mode index $\sigma_q^k$0. The velocity field and elevation are then expanded as: $\sigma_q^k$1 The surface elevation follows from $\sigma_q^k$2 (Farhat et al., 11 Jan 2026).

Damping (large $\sigma_q^k$3) and rapid planetary rotation suppress sharp modal resonances, resulting instead in broad, overdamped, nonperiodic responses.

3. Tidal Dissipation, Wave Interference, and Spatial–Temporal Variability

Tidal dissipation in the magma ocean is quantified locally as: $\sigma_q^k$4 with $\sigma_q^k$5 the magma density. The total velocity field, $\sigma_q^k$6 from various $\sigma_q^k$7 tidal components, yields interference terms $\sigma_q^k$8 introducing beat and sum frequencies $\sigma_q^k$9. For eccentricities $$\partial_t u + f \hat{z} \times u + \sigma_R u = g \nabla (\zeta_\text{eq} - \zeta)$$0, multiple tidal harmonics drive aperiodic, irregular patterns of tidal dissipation, temporally and longitudinally shifting hotspots, and stochastic, nontrivial thermal maps (Farhat et al., 11 Jan 2026).

Exemplar parameters for a 55 Cnc e analogue ($$\partial_t u + f \hat{z} \times u + \sigma_R u = g \nabla (\zeta_\text{eq} - \zeta)$$1, $$\partial_t u + f \hat{z} \times u + \sigma_R u = g \nabla (\zeta_\text{eq} - \zeta)$$2, $$\partial_t u + f \hat{z} \times u + \sigma_R u = g \nabla (\zeta_\text{eq} - \zeta)$$3, $$\partial_t u + f \hat{z} \times u + \sigma_R u = g \nabla (\zeta_\text{eq} - \zeta)$$4, $$\partial_t u + f \hat{z} \times u + \sigma_R u = g \nabla (\zeta_\text{eq} - \zeta)$$5, $$\partial_t u + f \hat{z} \times u + \sigma_R u = g \nabla (\zeta_\text{eq} - \zeta)$$6):

• $$\partial_t u + f \hat{z} \times u + \sigma_R u = g \nabla (\zeta_\text{eq} - \zeta)$$7 (tidal heating dominates over insolation)
• $$\partial_t u + f \hat{z} \times u + \sigma_R u = g \nabla (\zeta_\text{eq} - \zeta)$$8
• $$\partial_t u + f \hat{z} \times u + \sigma_R u = g \nabla (\zeta_\text{eq} - \zeta)$$9
• Hotspot longitude varies $\partial_t \zeta + \nabla \cdot (H u) = 0$0 per orbit (Farhat et al., 11 Jan 2026).

4. Light-Curve Predictions and Observational Diagnostics

The surface temperature at each point, $\partial_t \zeta + \nabla \cdot (H u) = 0$1, is controlled by both local absorption of stellar flux and tidal heating: $\partial_t \zeta + \nabla \cdot (H u) = 0$2 where $\partial_t \zeta + \nabla \cdot (H u) = 0$3 is the Stefan-Boltzmann constant and $\partial_t \zeta + \nabla \cdot (H u) = 0$4 is the emissivity. The disk-integrated observable thermal flux is

$\partial_t \zeta + \nabla \cdot (H u) = 0$5

Temporal and spatial modulation of $\partial_t \zeta + \nabla \cdot (H u) = 0$6 driven by mode interference yields orbit-to-orbit light-curve variability, with peak-to-peak temperature fluctuations up to $\partial_t \zeta + \nabla \cdot (H u) = 0$7 and disk-integrated flux spikes by factors of several within one orbital period for $\partial_t \zeta + \nabla \cdot (H u) = 0$8 (Farhat et al., 11 Jan 2026).

For $\partial_t \zeta + \nabla \cdot (H u) = 0$9, $u$0, while at $u$1, detectable $u$2 pulses of $u$3 are predicted if phase-curve precision reaches $u$4.

A plausible implication is that light-curve analysis can reveal both the presence and dynamical behavior of deep magma oceans and constrain global mantle properties.

5. Coupling to Mantle Convection and Thermal Equilibria

The vertical stratification beneath a magma ocean comprises three principal layers, determined by comparing the mantle adiabat, $u$5, to the solidus, $u$6, and the critical-melt curve $u$7:

• Fully molten ocean: $u$8, thickness $u$9
• Mushy two-phase layer: $f$0
• Solid viscoelastic mantle: $f$1

Typical Earth-mass planets with mantle potential temperature $f$2 present $f$3 (Farhat et al., 11 Jan 2026).

The total convective heat flux through fluid magma follows

$f$4

where Nu is the Nusselt number, $f$5, $f$6. For solid-state stagnant lids, $f$7.

Thermal equilibria are established when global tidal input $f$8 balances integrated convective output. Two branches emerge:

• Low-$f$9 (≤1600 K): Solid-state, weakly tidal, no fluid ocean.
• High-$\sigma_R$0 (≥2500 K): Fluid ocean, tidal power peaks for $\sigma_R$1–$\sigma_R$2 km, $\sigma_R$3–$\sigma_R$4 K for $\sigma_R$5 d and $\sigma_R$6. Fluid branch equilibria predict stable fluid oceans of $\sigma_R$7–$\sigma_R$8 km for Gyr timescales if eccentricity is replenished (Farhat et al., 11 Jan 2026).

6. Dynamical Consequences and Evolutionary Context

In early post-impact Earth-Moon history, the molten surface of the Earth enabled efficient angular momentum transfer, leading to rapid lunar recession up to $\sigma_R$9 Earth radii within $H$0–$H$1 years, in contrast to previous estimates based on solid-Earth models (Farhat et al., 2024). For exoplanets, the presence of a fluid magma ocean modifies tidal dissipation, removing planets from stable spin-orbit resonances and dramatically accelerating tidal synchronization—the path to synchronous rotation is shortened from Gyr to Myr timescales when compared to weakly viscous solid models (Farhat et al., 2024).

A dominant regime change occurs when tidal heating in the fluid layer overtakes stellar insolation in controlling surface thermodynamics. For certain short-period, eccentric systems, the surface temperature, day–night flux ratio, and thermal phase curve become direct probes of interior-exterior dynamical coupling (Farhat et al., 11 Jan 2026).

7. Longevity, Feedback, and Observability

Thermal settling times for the fluid-ocean branch are rapid ($H$2–$H$3 yr), while solid-state evolution is much slower ($H$4–$H$5 yr). Orbital inspiral due to tidal dissipation occurs on timescales $H$6–$H$7 yr in the fluid regime, confirming molten mantles can persist up to Gyr durations as long as eccentricity is continually supplied by planetary companions or secular effects (Farhat et al., 11 Jan 2026).

Time-variable thermal emission and aperiodic light-curve "spikes" are predicted for lava worlds with deep oceans ($H$8–$H$9 km) and eccentricities $g$0, such as analogues to 55 Cnc e, making these features compelling diagnostics for exoplanet characterization.

Table: Characteristic Parameters and Outcomes for Deep Magma Oceans on a 55 Cnc e Analogue

Quantity	Typical Value	Context/Implication
Eccentricity ($g$1)	0.03 – 0.05	Maintained by companions
Ocean depth ($g$2)	50 – 200 km	Fluid branch equilibrium
$g$3	$g$46 m s$g$5	Multimodal flow velocities
$g$6	$g$7 5	Tidal dominates insolation
$g$8 variation	$g$9800 K	Orbit-to-orbit, suborbital timescales

The combination of multimodal wave interference, episodic and dominant tidal dissipation, and convective–radiative coupling defines the unique phenomenology and evolutionary trajectories of planets and moons hosting deep magma oceans (Farhat et al., 2024, Farhat et al., 11 Jan 2026).