Constructing Extreme Heatwave Storylines with Differentiable Climate Models

TL;DR

This work tackles the challenge of identifying physically plausible upper bounds of extreme heatwaves without resorting to prohibitively large ensembles. It introduces a differentiable hybrid climate model, NeuralGCM, and uses gradient-based optimization to perturb initial conditions and generate extreme heatwave trajectories. Applied to the 2021 Pacific Northwest event, the method yields trajectories that exceed the most extreme member of a 75-member ensemble by up to about 3.7 degrees Celsius, with signatures of intensified atmospheric blocking and amplified Rossby waves. The approach offers a computationally efficient, process-based tool for rapid risk assessment and scenario storylines, with potential applicability to other extremes and differentiable climate models.

Abstract

Understanding the plausible upper bounds of extreme weather events is essential for risk assessment in a warming climate. Existing methods, based on large ensembles of physics-based models, are often computationally expensive or lack the fidelity needed to simulate rare, high-impact extremes. Here, we present a novel framework that leverages a differentiable hybrid climate model, NeuralGCM, to optimize initial conditions and generate physically consistent worst-case heatwave trajectories. Applied to the 2021 Pacific Northwest heatwave, our method produces heatwave intensity up to 3.7 $^\circ$C above the most extreme member of a 75-member ensemble. These trajectories feature intensified atmospheric blocking and amplified Rossby wave patterns-hallmarks of severe heat events. Our results demonstrate that differentiable climate models can efficiently explore the upper tails of event likelihoods, providing a powerful new approach for constructing targeted storylines of extreme weather under climate change.

Figures (10)

• Figure 1: (a) Histograms of 6-hourly temperature values from a 40-year NeuralGCM simulation (green) and from ERA5 (orange) reanalysis over the study domain outlined in Figure \ref{['fig:timeseries_main_var']}. Dashed line indicates the 95th and 99th percentiles. (b) Time series of temperature forecasts at 1000-hPa from NeuralGCM with 10-, 6-, 4-, and 2-day lead times (green colored lines) compared with ERA5 reanalysis data (orange line) for the PN2021 heatwave. Grey area highlights the targeted time range for the optimization process.
• Figure 2: Top row: Evolution of the difference in 500-hPa geopotential height ($\Delta Z$, in km) between the optimized simulation and the control run for EXP75. Black contours (optimized run) outline the amplified Rossby wave pattern, with deeper troughs and higher ridges compared to the control. Middle row: The difference in 1000-hPa temperature (${}^\circ$C) between the optimized and control simulations. Bottom row: Th difference between the 500-hPa geopotential height and 1000-hPa temperature averaged over the target domain.
• Figure 3: Top row: cross section of 500-hPa geopotential along a fixed latitude for EXP75. Bottom row: the amplitude difference of the Fourier spectrum, including wavenumbers 1–5 highlighted in grey. Red dots highlight wavenumbers 2-5.
• Figure 4: (a) Time series of 1000-hPa temperature for two optimized trajectories (50 and 75 steps) and a 75-member ensemble with its mean. (b) Spatial map of average temperature anomalies from the 75-step run. (c) Time series of 500-hPa geopotential height for the same set of simulations. (d) Spatial map of 500-hPa geopotential height anomalies during the event period. (e, f) Time series of heatwave intensity (defined by Eq. \ref{['eq:hw_def']}) and duration for the ensemble and optimized cases.
• Figure 5: Time series of (a) surface pressure, (b) near-surface specific humidity, (c) U-component of wind, (d) V-component of wind, (e) and (f) temperature advection for each component at 1000-hPa . Data from the optimized trajectory are shown alongside the individual ensemble members (in grey) and the ensemble mean (thick black line).