Spotlight #41 – Geophysics & Seismology: Earth’s Interior, Seismic Waves, Plate Tectonics and the Geomagnetic Dynamo

We live on a thin silicate shell above a churning iron ocean. Seismic waves have mapped the planet’s interior to kilometre precision, revealing a phase-layered structure driven by a heat engine left over from accretion 4.5 billion years ago. The same convective engine drives the motion of tectonic plates and sustains the geomagnetic field — a magnetohydrodynamic dynamo whose pole reversals are recorded in cooling lava.

Geophysics is observational physics pushed to its extreme: we cannot drill to the mantle (the deepest borehole, Kola SG-3, reached only 12.2 km), so almost everything we know about Earth’s interior comes from remote sensing: seismology, gravity, magnetics, and heat flow. The agreement between these independent data types and the Preliminary Reference Earth Model (PREM) is one of the great achievements of 20th-century science.

1. Earth’s Interior: The PREM Model

The Preliminary Reference Earth Model (Dziewonski & Anderson 1981) is a 1D spherically symmetric model that matches a global average of seismological observations. Despite enormous lateral heterogeneity in the real Earth (imaged by seismic tomography), PREM provides accurate travel times for teleseismic phases and remains the standard reference.

PREM Layers, Densities & Seismic Discontinuities

PREM layer structure:
  Layer              Depth (km)    v_P (km/s)  v_S (km/s)  ρ (g/cm³)
  Upper crust        0–15          6.0          3.5         2.60
  Lower crust        15–35         6.8          3.9         2.90
  Mantle (upper)     35–220        8.1          4.5         3.32
  Transition zone    410–660       9.5–10.3     5.1–5.6     3.7–3.9
  Lower mantle       660–2891      10.0–13.7    5.6–7.3     4.4–5.6
  Outer core (liq)  2891–5150     8.1–10.4        0        9.9–12.2
  Inner core (sol)  5150–6371    11.3          3.6         12.8–13.1

Major discontinuities:
  Mohorovičić (Moho): crustal base at ~35 km (continental), ~8 km (oceanic)
    Marked by sharp P-wave speedup from 6.8 → 8.1 km/s
  410 km: olivine → wadsleyite phase transition (exothermic; Clapeyron slope +3 MPa/K)
  660 km: ringwoodite → perovskite + ferropericlase (endothermic; −2 MPa/K)
    Larger density jump; thermochemical boundary partially blocking mantle circulation
  Core-mantle boundary (CMB, 2891 km):
    S-waves cannot propagate in liquid outer core → S-wave shadow zone
    Ultra-low velocity zones (ULVZs): partial melt patches at CMB, δV_P ≈ −10%
  Lehmann discontinuity (5150 km): liquid → solid inner core
    Inner core boundary: Pkikp vs Pkjkp phases probe inner core anisotropy

Pressure and temperature profile:
  Pressure: 0 (surface) → 24 GPa (660 km) → 136 GPa (CMB) → 364 GPa (centre)
  Temperature: 0°C (surface) → 1600°C (base lithosphere) → 2,500°C (CMB)
               → jump to 3,600°C (top outer core) → ~5,400°C (inner core boundary)
               [comparable to photosphere of Sun]
          

2. Seismic Waves and Travel Times

Seismic body waves fall into two families: P-waves (primary, compressional, longitudinal) that travel through fluids and solids, and S-waves (secondary, shear, transverse) that cannot propagate in the liquid outer core. Surface waves — Rayleigh and Love — are guided along the free surface and carry most destructive injury from large earthquakes.

Wave Velocities, Snell’s Law & Shadow Zones

Wave velocity definitions:
  v_P = √((K + 4G/3)/ρ)    [K = bulk modulus, G = shear modulus, ρ = density]
  v_S = √(G/ρ)
  v_P/v_S = √((K/G + 4/3)) ≥ √(4/3) ≈ 1.155 (solid); → ∞ in fluid (G=0)

Typical mantle velocities (PREM, 400 km depth):
  v_P ≈ 9.1 km/s,  v_S ≈ 5.0 km/s,  v_P/v_S ≈ 1.82
Outer core (pure liquid iron):
  v_P ≈ 8.1–10.4 km/s,  v_S = 0 (no shear → liquid)

Snell's law at curved interfaces (ray parameter p = const along ray):
  p = sin(i)/v(z) = r sin(i)/v(r)    (r = radial distance)
  Rays curve back upward when v(z) increases with depth (normal Earth case)
  Ray turning depth: p = r_turn/v(r_turn) → solving for r_turn given p

Shadow zones:
  P-wave shadow zone: 103°–142° angular distance from epicentre
    Outer core refracts P-waves strongly inward → gap in direct P arrivals
    PKP phase: P-waves that traverse the core (arrives >142°)
    PKIKP: penetrates inner core; arrives 0°–180° with travel time delay
  S-wave shadow zone: 103°–180°
    S-waves cannot enter liquid outer core → absent beyond 103°
    Discovery of liquid outer core: Oldham (1906), Gutenberg (1914)

Surface waves (dispersive, frequency-dependent depth sampling):
  Rayleigh: retrograde elliptical particle motion in sagittal plane (SV + P coupling)
    Group velocity ~0.92 v_S(depth ≈ λ/3)
  Love: horizontal transverse motion; requires velocity layering for guidance
    Phase velocity c_L: 2kH tan(k_L H) = (ρ₁v_L1)/(ρ₀v_S0)...  (layer dispersion relation)
  Surface wave tomography: 20–300 s period → 30–600 km depth sensitivity
          

3. Seismology: Magnitude, Focal Mechanisms and Statistics

Quantifying earthquake size requires separating the wave amplitude recorded at a receiver (a function of distance, radiation pattern, attenuation, and instrument response) from the intrinsic strength of the source. The moment magnitude scale accomplishes this by anchoring magnitude to the physical seismic moment M0, which integrates stress drop over the rupture area.

Magnitude Scales, Seismic Moment & Gutenberg-Richter

Richter local magnitude:
  M_L = log₁₀(A/A₀)    [A = max displacement amplitude in µm at 100 km, Wood-Anderson]
  Saturates at M_L ≈ 7 (seismic moment grows faster than amplitude)

Seismic moment:
  M₀ = μ · Ā · D̄    [μ = shear modulus ~30 GPa, Ā = rupture area, D̄ = mean slip]
  Derived from long-period body-wave or static GPS deformation

Moment magnitude (Hanks & Kanamori 1979):
  M_w = (2/3)(log₁₀ M₀ − 9.1)    [M₀ in N·m]
  Does not saturate; anchored to energy release
  Examples:
    Sumatra 2004: M_w = 9.1, M₀ ≈ 4×10²² N·m, fault area ≈ 1300×200 km²
    Valdivia 1960: M_w = 9.5 (largest instrumental; Chile)
    Chicxulub impact: M_w ≈ 10.5–11 (estimated)

Focal mechanism (fault-plane solution):
  Moment tensor M_ij: 2nd-order symmetric tensor describing force couples
  Eigenvalues: T (tension), N (null), P (pressure) axes
  Beach ball: lower hemisphere stereographic projection of first-motion polarities
    Quadrant compression (down-first)/dilatation (up-first) pattern
  Fault types:
    Normal:       hanging wall down (divergent margin; extensional normal stress)
    Thrust/Reverse: hanging wall up (convergent; compressional)
    Strike-slip:  horizontal motion (transform fault)

Gutenberg-Richter relation:
  log₁₀ N(M) = a − bM
  b ≈ 1.0 (globally; range 0.6–1.3 locally)
  Means: for every M7.0, there are ~10 M6.0 and ~100 M5.0 earthquakes
  Physical interpretation: scale-free fractal fault geometry

Omori-Utsu aftershock decay:
  n(t) = K/(t + c)^p    [n = aftershock rate, K,c,p = constants; p ≈ 1.0]
  Aftershock productivity: N_total ∝ 10^{α(M_main − M_min)}, α ≈ 1
          

4. Plate Tectonics and Mantle Dynamics

The theory of plate tectonics, consolidated in the 1960s (Hess, Vine, Matthews, Wilson, McKenzie, Parker), unified the previously separate fields of continental drift, seafloor spreading, seismicity, and volcanism under a single kinematic framework. The driving engine — mantle convection — remains an active research area, as it depends on the rheology of silicate rocks under extreme pressure and temperature.

Mantle Convection, Seafloor Spreading & GPS Plate Velocities

Mantle rheology:
  Effective viscosity: η_eff ≈ 10^{21} Pa·s (upper mantle post-glacial rebound)
                        10^{22}–10^{23} Pa·s (lower mantle)
  Dislocation + diffusion creep: ε̇ = A σⁿ exp(−E_a/RT) (power-law creep, n≈3.5)
  Transition from elastic to viscous: t_Maxwell = η/G ≈ 10^{21}/70×10^9 ≈ 450 years

Mantle convection dimensionless numbers:
  Rayleigh: Ra = ρ_0 g α ΔT d³/(κη)
    ρ₀ ≈ 4000 kg/m³; g = 10 m/s²; α ≈ 2×10⁻⁵ K⁻¹; ΔT ≈ 2500 K
    d ≈ 2900 km; κ ≈ 10⁻⁶ m²/s; η ≈ 3×10²¹ Pa·s
    Ra ≈ 10⁷ → vigorous convection (onset ~Ra_c = 1100 for plane layer)
  Stagnant lid vs mobile lid: depends on yield strength of cold lithosphere

Seafloor spreading and magnetic anomalies:
  Rate: mid-Atlantic ~2.5 cm/yr (slow); East Pacific Rise ~15 cm/yr (fast)
  Vine-Matthews-Morley (1963): new oceanic crust acquires geomagnetic field direction
    → alternating normal/reversed stripes parallel to ridge axis
  Age of oceanic crust: max ~200 Ma (Jurassic); older subducted
  Seafloor age = (distance from ridge) / (spreading rate)

GPS plate velocities (ITRF2020 / NNR-MORVEL56):
  Plate pair              Rate (mm/yr)  Azimuth
  Pacific–North America   ~50           NNW
  Indo-Australian–Eurasia ~45           NNE
  Africa–Eurasia          ~7            NNW
  North America–Eurasia   ~24           ENE (Atlantic opening)
  Nubia–Somalia           ~6            ENE (East African Rift)

Wilson cycle (complete ocean opening and closing):
  1 – Continental rifting (East African Rift today: early stage)
  2 – Proto-ocean (Red Sea: narrow ocean)
  3 – Open ocean (Atlantic: mature passive margins)
  4 – Subduction initiation → island arc (Nazca–South America)
  5 – Ocean closure + continental collision (Alps, Himalaya, Appalachians)
  Duration: ~500 Myr full cycle
          

5. Earth’s Gravity Field and the Geoid

The geoid is the equipotential surface of Earth’s gravity field that coincides (on average) with mean sea level. Deviations from an ellipsoid reflect density anomalies in the crust and mantle. The GRACE satellite mission (2002–2017) and its successor GRACE-FO (2018–present) measure the gravity field to spherical harmonic degree 120 (±330 km resolution) by tracking the changing distance between two spacecraft.

Gravity Corrections, Geoid & Isostasy

Gravity field expansion:
  V(r,θ,λ) = (GM/r)∑_l ∑_m (R/r)^l P_l^m(cosθ) [C_lm cos(mλ) + S_lm sin(mλ)]
  EGM2008 (NGA): complete to l_max = 2159 (~10 km resolution)
  C₂₀ = J₂ (Earth's oblateness): −1.0826 × 10⁻³  (oblate spheroid)

Gravity corrections (reducing to a standard surface):
  Free-air correction: +0.3086 μGal/m upward (−3.086 μGal/m observation altitude)
  Bouguer correction: −2πGρh ≈ −0.1119ρh mGal/m  (mass of slab between obs. & geoid)
    Complete Bouguer: terrain correction Tc removes local topography roughness

  Free-air anomaly: Δg_FA = g_obs − g_normal + FAC         [removes altitude effect]
  Bouguer anomaly:  Δg_B  = Δg_FA − 2πGρh + Tc             [removes crust slab]

Isostasy (Airy vs Pratt models):
  Airy (George Biddell Airy, 1855): constant density ρ_c; mountain root floats
    H_root = h · ρ_c/(ρ_m − ρ_c) ≈ 7h   (h = mountain height, excess density sunk)
    Prediction: negative Bouguer anomaly under mountains (observed globally)
  Pratt (John Henry Pratt, 1855): variable density σ(x), constant base depth T_0
    σ(x) · (T_0 + h) = σ₀ · T_0   (compensating by lateral density variation)

GRACE applications:
  Greenland ice mass: −270 Gt/yr (2003–2022 trend); Antarctica: −150 Gt/yr
  Groundwater depletion: Indus Basin −24 mm/yr equivalent water height
  Post-glacial rebound: Fennoscandia +10 mm/yr gravity increase
  Hydrological cycle: seasonal 10–20 mm geoid height variations Amazon/Congo basins
          

6. The Geomagnetic Field and Dynamo Theory

Earth’s magnetic field is generated by turbulent convection in the liquid iron outer core. The field is primarily dipolar (axial dipole ≈ 80% of surface field), but changes on timescales from milliseconds (secular variation from electromagnetically coupled core dynamics) to millions of years (polarity reversals). Palaeomagnetic records in seafloor basalt preserve the history of these reversals back to the Jurassic.

MHD Dynamo, Secular Variation & Palaeomagnetic Reversals

Geomagnetic field sources:
  Core field: 95% of total; generated by MHD dynamo; varies on years–Myr timescales
  Crustal (lithospheric) field: remanent magnetism; anomalies ±2000 nT; quasi-static
  External field: magnetosphere + ionosphere currents; hours-days variation; quiet Sq ≈ 30 nT

MHD induction equation:
  ∂B/∂t = ∇×(u×B) − ∇×(η ∇×B)
  η = 1/(μ₀σ) = magnetic diffusivity [m²/s]
  Balance: magnetic Reynolds number Rm = u L/η
  Earth's outer core: u ≈ 0.3 mm/s, L ≈ 2000 km, η ≈ 1 m²/s → Rm ≈ 600 >> 1
  High Rm: diffusion negligible → field frozen into fluid (Alfvén's theorem)

Busse (1975) columnar convection model:
  Planetary rotation (Ω = 7.3×10⁻⁵ rad/s) →  Coriolis organises convection into
  cylindrical rolls aligned with rotation axis
  Taylor-Proudman: in rapidly rotating fluids, flow is 2D perpendicular to Ω
  Convective rolls generate helical flow → α-effect (mean-field dynamo):
    ⟨u×B⟩ = αB + ... (Steenbeck, Krause, Rädler 1966)
  Ω-effect: differential rotation stretches poloidal → toroidal B
  Together: αΩ dynamo → sustains both components against diffusion

IGRF (International Geomagnetic Reference Field):
  Spherical harmonic representation of internal field to degree/order 13
  Current dipole moment: 7.94×10²² A·m² (down 6% since 1840)
  South Atlantic Anomaly (SAA): weak-field region; auroral particle flux at low altitude
  Geomagnetic North Pole (2024): ~86°N, 133°W (drifting toward Siberia ~40 km/yr)

Palaeomagnetic reversals:
  Geomagnetic Polarity Time Scale (GPTS) from oceanic magnetic anomalies + biostratigraphy
  Current: Brunhes Normal Chron (since 780 ka)
  Last reversal: Matuyama-Brunhes at 781 ka (10–20 kyr transition duration)
  Average reversal rate: 4–5 per Myr (Cenozoic); Cretaceous Normal Superchron: 40 Myr no reversals
  Precursors: excursions, intensity minima, multipolar field structure during transition
  Potential consequences of reversal: 10× SAA expansion, aurora at mid-latitudes
          

Seismic tomography and mantle plumes: 3D tomographic models (LLSVP structures — Large Low-Shear-Velocity Provinces beneath Africa and the Pacific) image two large anomalous regions at the CMB, 1000–1500 km across, with δvS ≈ −3%. Whether these are chemical heterogeneities or purely thermal (or both) is actively debated. They appear to anchor hotspot tracks (Hawaii, Iceland, Yellowstone) over hundreds of millions of years.

Try These Simulations