The European summer heatwave of 2019 – a regional storyline perspective

The number and intensity of heatwaves have increased in the recent past, along with anthropogenic climate change. This poses challenges to many communities and emphasises the need to develop adaptation measures based on more accurate information regarding regional-to-local changes in temperature extremes and their impacts. While the general increase in global mean temperature is well established, current global climate projections show a large model spread concerning possible future circulation changes. To isolate the more certain thermodynamic response from the less certain dynamical response to anthropogenic climate change, we employ an event-based storyline approach and focus the present study on the 2019 summer heatwaves that occurred over central Europe. Our approach comprises three steps. Firstly, the large-scale circulation in the free troposphere was spectrally nudged to the ERA5 reanalyses within the global coupled climate model AWI-CM-1.1-MR for the recent period of 2017–2022, corresponding to global warming of +1.4 K. This process was then repeated under pre-industrial conditions, as well as under +2, +3, and +4 K global-warming climates. Secondly, the global storylines were dynamically downscaled using the regional ICOsahedral Nonhydrostatic model in Climate Limited-area Mode (ICON-CLM) to a EURO-CORDEX domain with a horizontal resolution of 12 km, and, thirdly, they were downscaled to a central European (German) domain with a resolution of 3 km. We provide evidence that the downscaling of global storyline integrations significantly improved the representation of present-day temperature patterns and reduced errors in daily 2 m temperatures relative to observations from central Europe. The magnitude of the heatwave temperature response significantly exceeds the globally modelled background warming, with distinct spatial and temporal variations in the regional increments. Our simulations indicate a general linear dependence of the 2 m temperature response on global-warming levels: during the July 2019 heatwave, warming rates ranged between factors of 2 and 3 in central Europe, resulting in anthropogenic warming of 8 to 12 °C in the +4 K climate. The spatial extent and duration of the heatwave are also amplified in the warmer climates. With this three-step downscaling approach, we gain new insights into possible future changes in heat extremes in central Europe, which apparently surpass global-warming trends. Along with its scientific value, our method provides ways to facilitate the communication of regional climate change information to users.

1 Introduction

Heatwaves are a major natural hazard worldwide, with the heatwaves of 2003, 2010, and 2018 serving as prominent examples for the European continent (e.g. Fink et al., 2004; Barriopedro et al., 2011; Miralles et al., 2014; Spensberger et al., 2020). In the last 2 decades, Europe has witnessed an increase in the frequency, duration, and intensity of extreme heat events, which, in turn, has led to an increase in mortality rates, food and water insecurity, and long-term economic and cultural stress (Robine et al., 2008; García-Herrera et al., 2010; Perkins-Kirkpatrick and Lewis, 2020; Becker et al., 2022; Calvin et al., 2023; Knutzen et al., 2025). The general effects of climate change on heatwave characteristics have been shown to be robust in recent studies, which are often aimed at estimating the impact of human-induced greenhouse gas (GHG) forcing on recent extreme events using observational records and at projecting future changes by means of climate modelling (e.g. Barriopedro et al., 2011, 2023). The traditional approaches, which commonly involve multi-model averaging and probabilistic event attribution, can provide estimates of trends in the frequency, intensity, and persistence of extreme events, but they lack clarity regarding the physical processes causing the changes (Shepherd, 2014).

One current scientific challenge involves disentangling the relative roles of dynamical and thermodynamic contributions to future heatwave characteristics in attribution and projection studies (Shepherd, 2014, 2021; Sousa et al., 2020; Sánchez-Benítez et al., 2018). Thermodynamic effects, such as near-surface warming, moistening of the atmosphere, and the partitioning of radiative and turbulent fluxes, show a relatively robust response to anthropogenic GHG forcing in models and generally tend to exhibit less internal variability (Deser et al., 2014; Wehrli et al., 2018). On the other hand, dynamic effects, which include changes in the position, strength, and meandering of the jet stream, as well as changes in the occurrence of weather regimes and more localised circulation patterns, are subject to larger uncertainties (Deser et al., 2014; Shepherd et al., 2018; Zappa, 2019). First, this is associated with inherent model uncertainties and differences in the parameterisation of unresolved processes in the models (Shepherd, 2014). Second, the internal variability in the dynamical component of the atmosphere is responsible for a low signal-to-noise ratio in studies aimed at quantifying the regional response of extreme events to global warming (Deser et al., 2014; Shepherd, 2014, 2021; Wehrli et al., 2018; Barriopedro et al., 2023).

An alternative to circumventing the uncertainties associated with different atmospheric-circulation changes under enhanced GHG forcing is to use an event-based storyline approach, in which the dynamical conditions are constrained to the present-day state in a specified way (Caviedes-Voullième and Shepherd, 2023). Along with the uncertainty in changing dynamics, the internal variability is considerably reduced in event-based storylines (Sánchez-Benítez et al., 2022), which improves the signal-to-noise ratio, providing us with the opportunity to better quantify the actual event-specific thermodynamic response. Assuming the robustness of this quantification, we improve our understanding of the potential impacts of future extreme events and communicate them to the public and authorities in a more understandable manner.

Constraining the dynamical conditions on a regional scale can be achieved with the pseudo-global-warming (PGW) approach, i.e. by perturbing boundary conditions obtained from reanalyses using the average climate change signal (deltas) from global climate models (GCMs) (Schär et al., 1996; Aalbers et al., 2023; Ludwig et al., 2023; Vries et al., 2024). As only smoothed, multi-year-averaged GCM fields are used to modify the boundary conditions, this method is computationally effective, representing a significant advantage for multi-model and multi-ensemble studies of regional thermodynamic responses to global warming (Brogli et al., 2023).

While the inter-annual variability in delta fields obtained from GCMs is intentionally suppressed in the PGW approach, it can be consistently taken into account in spectrally nudged storylines. In this storyline method, a GCM is run by nudging the upper-tropospheric winds to force the large-scale circulation to follow the reanalysis state, while the background climate corresponds to a specific warming level (Sánchez-Benítez et al., 2022; Athanase et al., 2024; van Garderen et al., 2021; Wehrli et al., 2020). If the GCM is a coupled model, as in Sánchez-Benítez et al. (2022), no assumptions about the deltas of sea surface temperature and sea ice content need to be made (van Garderen et al., 2021).

On the other hand, one advantage of the PGW approach over the nudged storyline approach is its potential for avoiding GCM-specific biases by repeating the experiment with deltas derived from multi-model ensemble means or different single-model multi-member ensemble means (see e.g. Aalbers et al., 2023; Vries et al., 2024). In our work, we follow the path prescribed by a single GCM. In the context of the storyline approach, this unfolding of events is physically self-consistent and plausible, adhering to the storyline definition introduced by Shepherd et al. (2018) and allowing for a process-oriented evaluation of the obtained responses.

Nevertheless, GCMs have difficulty representing regional climate means and variability due to unresolved orography and shortcomings in model parameterisations associated with their coarse horizontal resolution (Giorgi and Gutowski, 2015). For example, using global nudged storylines, Sánchez-Benítez et al. (2022) found a strong amplification of the July 2019 European heatwave under global warming; however, their global simulations underestimated the high-temperature extremes reached during the heatwave. An effective method for addressing the lack of precision in GCMs involves performing dynamical downscaling with a regional climate model (RCM) (Feser et al., 2011; Giorgi, 2019; Vautard et al., 2021). Commonly, the horizontal resolution of RCMs applied to the European domain is on the order of 10–12 km (e.g. Jacob et al., 2014; Giorgi and Gutowski, 2015). Being computationally efficient, this resolution allows for the production of large ensembles of simulations, significantly improving the representation of relevant climatological variables compared to GCMs (Vautard et al., 2021). Nevertheless, convective processes can only be resolved by regional convection-permitting models (CPMs) operating at resolutions finer than 4 km (Prein et al., 2015; Giorgi, 2019; Hundhausen et al., 2023). Associated with explicitly resolved deep convection and better-captured processes in regions with complex topography, CPMs have been shown to add further value to RCMs in their representation of precipitation and near-surface temperature, especially in regions with complex topography (Prein et al., 2015; Giorgi, 2019).

In this study, we go a step beyond the global storyline approach by providing a regional perspective on the heatwave that occurred in Europe during July 2019 and its unfolding in colder (pre-industrial) and warmer climates (see e.g. Sánchez-Benítez et al., 2022; Sousa et al., 2020). With this aim, we dynamically downscale the global spectrally nudged storylines for the summer 2019 heatwave across five different background climates, ranging from pre-industrial conditions to +4 K global warming, for (central) Europe. The approach utilises a global-to-regional (GCM–RCM–CPM) model chain comprising global spectrally nudged storyline simulations obtained from the Alfred Wegener Institute's global coupled climate model AWI-CM-1.1-MR (hereafter referred to as AWI-CM1; Semmler et al., 2020), with large-scale horizontal winds spectrally nudged to ERA5 (Hersbach et al., 2020), and the ICOsahedral Nonhydrostatic model in Climate Limited-area Mode (ICON-CLM) (Pham et al., 2021) for dynamical downscaling – first to a 12 km horizontal grid spacing over Europe and subsequently to a 3 km horizontal grid spacing over central Europe. This approach permits the derivation of climate data at a high resolution, thereby providing detailed information for attribution and impact studies.

We address the following research questions:

1. How accurately can a regional event-based storyline simulation represent a recent event, and what improvement does it offer compared to the global spectrally nudged storyline simulation?

2. What was the effect of climate change on the 2019 European heatwave based on regional and convection-permitting ICON-CLM simulations?

3. What local-to-regional extreme temperature scaling occurs in response to global warming for an event like the 2019 heatwave?

The paper is structured as follows. Section 2 describes the global and regional model setups, as well as the datasets used for model evaluation. The main results are presented in Sect. 3, with the first research question addressed in Sect. 3.1 and the regional storylines analysed in Sect. 3.2 and 3.3. Section 4 summarises and discusses the results, formulates the main conclusions for each research question, and provides future research ideas.

2 Data and methods

2.1 Global spectrally nudged storylines

The global spectrally nudged simulations are based on the global coupled climate model AWI-CM1 (Semmler et al., 2020). This model has contributed to the sixth phase of the Coupled Model Intercomparison Project (CMIP6; Eyring et al., 2016). It consists of the atmospheric model ECHAM6.3.04p1 from the Max Planck Institute for Meteorology (MPI-M) (Stevens et al., 2013), coupled with version 1.4 of the Finite-volumE Sea ice–Ocean Model (FESOM) for the ocean component (Wang et al., 2014). The atmospheric component is run at a T127L95 spectral resolution, i.e. with truncation at wavenumber 127, which corresponds to a horizontal resolution of about 100 km in the tropics and 95 vertical levels extending up to ∼0.01 hPa. The ocean model FESOM uses an unstructured mesh that allows for fine resolutions in energetically active areas such as the Gulf Stream (Sidorenko et al., 2015; Sein et al., 2017). Consequently, the horizontal resolution of the ocean ranges from 80 km in the subtropical Pacific to 8–10 km in the North Sea and 8–20 km in the Arctic (see Fig. 1 in Semmler et al., 2020).

In the storyline experiments, the evolution of the AWI-CM1 large-scale atmospheric circulation is constrained by spectrally nudging the model's vorticity and divergence (representing large-scale horizontal winds) to ERA5 reanalysis data (Hersbach et al., 2020), with an e-folding time of 24 h and a spectral truncation of 20 for zonal wavenumbers. Nudging is applied only to vertical levels between 700 and 100 hPa. This configuration has been shown to optimally constrain large-scale events, such as heatwaves (Sánchez-Benítez et al., 2022), warm and moist intrusions in the Arctic (Pithan et al., 2023), and marine heatwaves (Athanase et al., 2024), while preserving some freedom in the boundary layer and at small spatiotemporal timescales.

A series of nudged storyline simulations were conducted using ESM-Tools (Barbi et al., 2021) for a range of climate states based on the configuration described above. Specifically, nudging experiments were branched off from the historical CMIP6 runs (Semmler et al., 2020) for 1 January 1851 to produce pre-industrial climate conditions. Meanwhile, climate scenarios for the present day (+1.4 K) and for +2, +3, and +4 K (relative to pre-industrial) conditions were made by spawning the CMIP6 experiments corresponding to the Shared Socioeconomic Pathway scenario SSP370 (Semmler et al., 2019) on 1 January of the years 2017, 2038, 2065, and 2093, respectively. These years were selected based on when the global-warming levels were reached. Each of the storylines comprises five ensemble members, spawned from the five respective CMIP6 simulations and thus started using slightly different initial conditions. Each storyline is simulated continuously from 1 January 2017 to 30 September 2022. Throughout the text, for all storylines, we refer to the years corresponding to the present-day circulation inferred from ERA5 as the “dynamical years” of 2017–2022.

2.2 Dynamical downscaling with ICON-CLM

Next, data from the global spectrally nudged storylines described in Sect. 2.1 were used as initial and boundary conditions to drive version 2.6.5.1 of the ICON (ICOsahedral Nonhydrostatic) model (Zängl et al., 2015), implemented using a regional climate configuration known as ICON-CLM (Pham et al., 2021). We used the runtime environment “Starter Package for ICON-CLM Experiments” (SPICE v2.0) (Rockel and Geyer, 2022) to conduct all the simulations. For each storyline, we simulated the full 2017–2022 period, as described in Sect. 2.1. The downscaling was performed on a EURO-CORDEX domain (Jacob et al., 2014) with an R12B5 resolution, meaning that each edge of the base icosahedron was initially divided into 12 parts, followed by five bisections of each resulting triangle's edges. This corresponds to a horizontal grid spacing of 0.11° (or 12 km), hereafter referred to as EUR-12 (Fig. 1a; see e.g. Prill et al., 2023, for a description of icosahedral grid-spacing conventions). Subsequently, we ran a nested ICON-CLM simulation at an R13B7 resolution, corresponding to a horizontal grid spacing of 0.0275° (or 3 km), on ensemble member 1 for the extended German domain (hereafter referred to as GER-3), including the peripheral hydrological catchment areas (see Fig. 1b). In this study, we only focus on ensemble member 1 for the GER-3 domain. A detailed analysis of the full GER-3 ensemble is thus reserved for another study.

Figure 1 (a) The EUR-12 and GER-3 domains used for downscaling the global AWI-CM1 storylines with ICON-CLM. (b) The GER-3 domain, showing the locations of the DWD (German Weather Service) observation stations and the area used for spatial averaging (48–51° N, 6–10° E). Colours correspond to the orography used in the simulations. The shown domains include the lateral boundary zone.

At the upper boundary, grid-point nudging was applied in the ICON-CLM simulations to maintain proximity to the present-time circulation represented by the global spectrally nudged AWI-CM1 runs. In ICON, this nudging was implemented as an additional forcing term, added to the prognostic equations at each fast-physics time step (Prill et al., 2023):

$$\begin{array}{}\text{(1)}& \mathit{\psi }\left(t\right)={\mathit{\psi }}^{\ast }\left(t\right)+{\mathit{\alpha }}_{\mathrm{nudge}}{N}_{\mathrm{ds}}\left[{\mathit{\psi }}_{\mathrm{bc}}\right(t)-{\mathit{\psi }}^{\ast }(t\left)\right],\end{array}$$

where ψ_bc(t) is the value of the prognostic variable ψ at time t, taken from the driving model; ψ^∗ is the value of the variable ψ before the nudging; α_nudge refers to the nudging strength; and N_ds is the number of dynamics substeps per fast-physics step (Prill et al., 2023). Upper-boundary nudging is applied to a sponge layer of a chosen thickness, where the nudging strength increases quadratically with height, starting with zero at the nudging start height (z_start) and reaching the maximum nudging coefficient (B₀) at the model's top height:

$$\begin{array}{}\text{(2)}& {\mathit{\alpha }}_{\mathrm{nudge}}=B_{\mathrm{0}}{\left({\displaystyle \frac{z-z_{\mathrm{start}}}{z_{\mathrm{top}}-z_{\mathrm{start}}}}\right)}^{\mathrm{2}}.\end{array}$$

The nudging coefficient for the thermodynamic prognostic variables θ_v, ρ, and q_v was set to zero to prevent overfitting and, even more importantly in this context, to allow for the free development of thermodynamics. Thus, upper-boundary nudging was applied only to the horizontal velocity. We kept the maximum coefficient (B₀) at its default value of 0.04. The nudging start height (z_start) in the EUR-12 domain was set to 5000 m, while for the GER-3 domain, it was left at 10 500 m to prevent any interaction with deep convection. Further information can be found in Prill et al. (2023).

Soil temperature and soil moisture data from ERA5 were used for soil initialisation in all EUR-12 simulations due to the partial unavailability of soil temperature data in AWI-CM1 outputs. In ICON, the initialising soil data are pre-processed and remapped onto the eight-layer mesh by the built-in algorithm (Prill et al., 2023; Pham et al., 2021). To account for the possible discrepancy in soil types between ICON and ERA5, the volumetric soil moisture was transformed to align with the universal soil moisture index (SMI), making it independent of the soil type (Prill et al., 2023). For the present-time experiment, the initial simulation year (2017) was considered for spin-up. In the storyline simulations, an additional year was required for the soil to adapt to the warmer climate. Therefore, we ran the dynamical year 2017 twice. Additionally, the temperature of the bottom soil layer, which is not prognostic in the TERRA land module of ICON but is set to the climatological annual mean near-surface temperature (T_CL) based on Climate Research Unit data (Schulz et al., 2016; Mitchell and Jones, 2005), was adjusted to reflect the respective global-warming level for each storyline. The lower-boundary condition for soil moisture was provided by a free-drainage formulation and thus did not require additional adjustments (Prill et al., 2023; Chen et al., 2018; Zeng and Decker, 2009).

2.3 Model evaluation approach

We compare the obtained daily 2 m temperature fields to the ERA5 reanalyses (Hersbach et al., 2020; Hersbach et al., 2023) and to E-OBS (v28), a daily gridded, land-only observational dataset for Europe (see e.g. Cornes et al., 2018; Copernicus Climate Change Service, 2023). To evaluate the simulations at the station level, we use 20 stations from the German Weather Service (DWD) (DWD, 2023; Kaspar et al., 2013; see the full list in Table S1 in the Supplement). The root mean square difference (RMSD) corresponding to the observational datasets (DWD and E-OBS) and its change between simulations of different resolutions (ΔRMSD) are chosen as metrics to compare the representations of near-surface temperatures provided by the models in our GCM–RCM–CPM chain for the present-day storyline. The significance of ΔRMSD between ICON EUR-12 and GER-3 is computed with the paired-difference test (Rubin, 1973).

3 Results

3.1 Evaluation of the present-day storylines

The results of the ICON regional model simulations are evaluated by comparing the model output with observational E-OBS and ERA5 reanalysis data. A comparison of the 2 m temperature fields and 500 hPa geopotential-height fields for the heatwave peaks of June and July 2019 shows good agreement between the ICON EUR-12 simulations (Fig. 2c and f) and ERA5 reanalyses (Fig. 2a and d). Furthermore, the ICON simulations show a clear improvement in the representation of daily 2 m temperature compared to the global AWI-CM1 simulations (Fig. 2b and e). A comparison of the daily temperature fields with E-OBS yielded similar results (see Fig. S1 in the Supplement).

Figure 2 Mean 2 m temperature (colours) and geopotential height at 500 hPa (contours) on 26 June 2019 (the first peak of the June heatwave) and 25 July 2019 (the peak of the July heatwave) for ERA5 (a, d), AWI-CM1 (b, e), and ICON EUR-12 (c, f).

In Fig. 3, time series of daily maximum, mean, and minimum 2 m temperatures are compared with E-OBS over a given latitude–longitude area (48–51° N, 6–10° E) to investigate the underestimation of maximum temperatures during the exceptionally hot periods of June and July 2019, as mentioned in Sánchez-Benítez et al. (2022). The ICON simulations show an improvement in both the daily maximum and minimum temperatures for July in this area. This indicates a more accurate representation of the diurnal temperature range in the regional ICON simulations.

Figure 3 Comparison with E-OBS of daily (a) maximum, (b) mean, and (c) minimum 2 m temperatures obtained with AWI-CM1, ICON EUR-12, and ICON GER-3 (ensemble member 1 – EM1), averaged across a latitude–longitude box with the following boundaries: 48–51° N, 6–10° E (dashed box in Fig. 1b). Shading indicates the minimum–maximum range of values obtained from the five-member ensembles. The numbers shown in brackets in the legends display the RMSD in ensemble member 1 relative to E-OBS for June, July, and August.

We compared the performance of the simulations within the model chain by calculating the root mean square difference (RMSD) between the simulated 2 m temperature and DWD observations (DWD, 2023) for June–August, using data from 20 selected stations (locations are shown in Fig. 1b, with the full list of stations given in Table S1). The RMSD is significantly reduced by the dynamical downscaling of AWI-CM1 data to the EUR-12 domain (see Table 1). In the case of the GER-3 simulation, a further reduction in the RMSD could only be achieved for the daily maximum 2 m temperatures, whereas the daily minimum and mean temperatures slightly deteriorated but remained clearly improved compared to those from AWI-CM1. A similar result can be obtained when the RMSD is calculated between the simulations and the E-OBS time series shown in Fig. 3 (see the RMSD values in the legend).

Table 1 The root mean square differences (RMSDs), expressed in degrees Celsius, of daily maximum, mean, and minimum 2 m temperatures (TMAX_2M, TMEAN_2M, and TMIN_2M, respectively) relative to DWD station observations recorded during summer (June–August) 2019, averaged over 20 stations (for locations, see Fig. 1).

Figure 4a–c shows the spatial distribution of the RMSD between EUR-12 and E-OBS for the daily maximum, mean, and minimum temperatures corresponding to the summer of 2019. The RMSD varies between 1 and 2 °C in central Europe for daily mean temperatures and between 2 and 3.5 °C for daily minimum and maximum temperatures, which aligns with the values shown in Table 1. We interpolated the ICON EUR-12 data to the grid of AWI-CM1 and compared the RMSDs for both models to E-OBS. The green colours in Fig. 4d–f indicate an improvement in the performance of the regional simulations compared to the global AWI-CM1 simulation for the respective temperatures. The overall improvement in central and southern Europe is robust, with the exception of the maximum temperatures found along the coast of the Iberian Peninsula and north of the Black Sea.

Figure 4 Performance assessment of the ICON EUR-12 simulations for June, July, and August 2019 for daily (a, d) maximum, (b, e) mean, and (c, f) minimum 2 m temperatures. (a–c) Root mean square difference (RMSD) in the simulated daily 2 m temperatures from ICON EUR-12 relative to E-OBS. (d–f) Change in RMSD achieved by dynamical downscaling; the green colours correspond to the reduced squared error in daily temperatures.

To identify potential systematic biases, the observed E-OBS temperature fields were compared to the EUR-12 ones. From Figs. S2 and S3, it can be inferred that the negative bias in daily maximum temperatures in central Europe, introduced in the AWI-CM1 simulations, is reduced (see Fig. S3a and d), along with the positive bias in daily minimum temperatures (Fig. S3c and f). According to Fig. S3a and g, the increased error in daily maximum temperatures in western Iberia and eastern Europe, mentioned above, can be attributed to the positive bias introduced by ICON-CLM.

A similar comparison was conducted for all simulated summer seasons between 2018 and 2022 (Figs. S4 and S5). This comparison demonstrates that the bias patterns persist consistently across all the simulated summers, suggesting that the errors are not flow-dependent.

The nested convection-permitting GER-3 simulation was assessed by comparing the RMSD in the 2 m temperature to the RMSD corresponding to the driving EUR-12 simulation (see Fig. 5). To this end, both temperature datasets were interpolated to the E-OBS grid (0.1° horizontal resolution). We also evaluated the average bias between the nested simulations by subtracting the EUR-12 seasonal mean temperature fields from the GER-3 fields, which were coarsened to the EUR-12 domain (see Fig. S6a–c). The bias of the GER-3 fields relative to E-OBS is shown in Fig. S6d–f.

Figure 5 Evaluation of the nested convection-permitting GER-3 simulations for June, July, and August 2019 for daily (a, d) maximum, (b, e) mean, and (c, f) minimum 2 m temperatures. (a–c) Root mean square difference (RMSD) in daily 2 m temperature relative to E-OBS for the GER-3 simulation. (d–f) Change in RMSD with respect to (w.r.t.) the EUR-12 simulation, with significant differences indicated by hatching (p<0.05).

The nested GER-3 simulation is between 0.5 and 2 K warmer than the EUR-12 simulation, which further reduces the negative bias in the daily maximum temperature over Germany. However, the positive bias in the daily minimum temperature increases, and the representation of the daily mean temperature shows no robust improvement over Germany. The most substantial bias in daily maximum temperature occurred in the Po Valley and in the westernmost parts of the domain over land. However, these areas are located within the lateral boundary zone and are not used in the analyses (see Fig. S6). Due to refined topography, the nesting significantly reduced the RMSD in maximum and mean 2 m temperatures over the Alps (see Fig. 5d–f).

Analogous to the assessment of the EUR-12 simulation, the evaluation of the GER-3 simulation over the full simulation period demonstrates that the patterns shown in Fig. 5 remain consistent over the years for all modelled summer seasons between 2018 and 2022 (see Figs. S7 and S8), with a more prominent improvement in the RMSD observed for the maximum temperatures (Figs. 5d and S7d).

3.2 Storylines for the summer 2019 heatwaves

Given the improved performance of dynamical downscaling with ICON-CLM under present-day conditions, we now analyse the regionalised past and future analogues of the July 2019 European heatwave.

We first consider the period corresponding to the peak of the July 2019 heatwave. According to Fig. 6, the maximum temperature on 25 July exceeds 40 °C over a considerable area of western Europe in a +4 K climate. The area affected by temperatures exceeding 40 °C is projected to increase significantly in the EUR-12 simulations, from 21 000 ± 14 000 km² in the pre-industrial climate to 290 000 ± 40 000 km² in the present-day climate and 1 110 000 ± 70 000 km² in the +4 K climate (see Fig. S9). Moreover, the 45 °C threshold is projected to be surpassed across a large area of western France in the +3 K climate and across Benelux and the Rhine Valley in the +4 K warmer world (Fig. 6i and j).

Figure 6 Daily maximum 2 m temperatures on 25 July for the ICON EUR-12 simulations (a, c, e, g, i) and ICON GER-3 simulations (b, d, f, h, j) with respect to pre-industrial (a, b), present-time (c, d), +2 K (e, f), +3 K (g, h), and +4 K (i, j) climates (based on ensemble member 1). For the ensemble spread of the EUR-12 simulations, see Fig. S11.

The increase in 2 m temperature is not spatially homogeneous, with regions located outside of the heatwave's core experiencing stronger warming, thereby contributing to the increasing spatial extent of the heatwave. This is exemplified in Fig. S12, where the maximum temperature differences between the +4 K and pre-industrial climates reach 12 °C in Luxembourg, southern Belgium, western Germany, and the easternmost parts of France. The temperature increase in the core of the heatwave is close to 8 °C, which corresponds to a doubling of the global-warming level (+4 K).

Figure 7 displays the time series of daily maximum, mean, and minimum temperatures, averaged over the given latitude–longitude area (48–51° N and 6–10° E; depicted in Fig. 1b), for the extended summer season (May to September) of 2019. During July and August, an increase in spread is found between the temperature curves corresponding to the different warming levels. Conversely, the spread between the time series of temperatures is much smaller in May and early June. This finding was previously confirmed in the analyses of global storylines conducted by Sánchez-Benítez et al. (2022) and will be further explored in the following section. We also observe a larger spread in the daily maximum curves for July and August compared to the mean and minimum temperatures.

Figure 7 Daily (a) maximum, (b) mean, and (c) minimum temperatures averaged over a given latitude–longitude box (48–51° N, 6–10° E; see Fig. 1b) for the May–September (MJJAS) period of 2019, based on the EUR-12 storyline simulations. Shading indicates the minimum–maximum range of values obtained from the five-member ensembles. The three highlighted periods (orange) are discussed in detail in Sect. 3.3. For the GER-3 simulations, see Fig. S13. PI: pre-industrial climate. PR: present-day climate. K2: +2 K climate. K3: +3 K climate. K4: +4 K climate.

3.3 Temperature scaling in response to global warming

In order to gain a deeper insight into the spread of temperature curves depicted in Fig. 7 and to address the question of temperature scaling in response to global warming, three 5 d periods were selected for detailed analysis (highlighted in orange in Fig. 7). The first period corresponds to mid-June, when no heatwave was observed; the second period corresponds to late June during the first heatwave; and the third period occurs around the peak of the July heatwave. The daily maximum, mean, and minimum 2 m temperatures – averaged spatially over the area 48–51° N, 6–10° E, and temporally over the aforementioned periods – show a clear linear dependence on the global-warming level (see Fig. 8a–c). Therefore, we express anthropogenic change in 2 m temperature per 1 °C of global warming as the slope of this line. This slope will be referred to in the following text as the “warming rate”.

Figure 8 (a–c) Daily maximum (red), mean (orange), and minimum (blue) 2 m temperatures over a given latitude–longitude box (48–51° N, 6–10° E), averaged over three 5 d periods and plotted against the global-warming level. The numbers in the legends represent the ensemble means of the slopes of the respective lines. (d) Warming rates for the rolling averages (5 d windows) of daily maximum, mean, and minimum temperatures over the same box. Shading indicates the minimum–maximum range of values obtained from the five-member ensembles. The three highlighted periods are discussed in detail in Sect. 3.3. This figure is based on the EUR-12 simulations. For the GER-3 simulations, see Fig. S14. WR: warming rate.

As shown in Fig. 8a, the average warming rate over the studied area is close to a factor of 1 in mid-June for all three curves. This indicates that the warming rate is comparable to the global-warming level in the absence of an extreme event in early summer. The warming rates increase during the first heatwave and approach or even exceed a factor of 2 for maximum temperatures during the July heatwave (Fig. 8b and c). This suggests that the expected differences in maximum temperature for a similar event in a future climate would be exacerbated and twice as high as the corresponding global-warming level.

Using this method, we computed the warming rates for 5 d running-mean maximum, mean, and minimum temperatures (see Fig. 8d) for the extended summer season of 2019. Between May and early June, the warming rates fluctuate around 1. The values increase with the onset of the first heatwave in late June, with the maximum temperature responding more strongly than the minimum temperature. Following a relative minimum in mid-July, the warming rates strongly increase before 25 July, reaching a factor of 2.2 for maximum temperatures ahead of and during the late July heatwave. After a short decrease immediately following the temperature maximum, the warming rates increase again to values over 2.0 at the beginning of August. This is followed by a gradual decrease before the last peak of the warming rates, which has an even higher magnitude and occurs in late August, showing a further temperature increase during the summer. We observe a general tendency for the warming rates to rise towards late summer, possibly associated with the gradual decrease in soil moisture. The increase in global-warming amplification in central Europe was also identified by Sánchez-Benítez et al. (2022) in both free and nudged AWI-CM1 simulations.

To comment on the broader warming-rate peak preceding the July heatwave seen in Fig. 8d, we estimate the duration of this heatwave across different climates based on the exceedance of the 90th percentile of the modelled maximum 2 m temperature, which is 30 °C when computed for all days in July over the simulated period of 2018–2022. According to our estimation, the duration of this heatwave would increase non-linearly from 4 d in the present-day climate to 5 d in the +3 K climate and 9 d in the +4 K warmer world (not shown).

While the response of the maximum temperature in May, June, and the second half of September is smaller or comparable to the response of the mean and minimum temperatures, the warming rate for the maximum temperature is very large during July and August. Thus, with global warming, the diurnal temperature range tends to increase during the mid- and late-summer periods of the dynamical year 2019. In contrast, the temperature response appears to be distributed more uniformly during the day in spring, early summer, and early autumn. We obtained similar behaviour in the other summer seasons (May–September) of the simulated period.

The warming rates are now computed for each grid point in the study area. Figure 9a displays the warming rates for the 5 d period around the peak of the heatwave in late July 2019 (the third shaded area in Fig. 7) at each grid point in both the EUR-12 and GER-3 domains. According to the goodness-of-fit maps in Fig. S15, the assumption of linear growth in the areas affected by the heatwave is valid in all cases. For comparability with the GER-3 simulations, Fig. 9 shows the warming rates obtained from ensemble member 1 for both domains. The ensemble mean of the temperature response in the EUR-12 simulations is shown in Fig. S16, which depicts similar patterns with small differences in detail.

Figure 9 (a) Warming rates for the period from 23 to 27 July 2019. (b) Warming rates for 25 July. The contours represent the geopotential heights from the EUR-12 simulations. This figure corresponds to ensemble member 1.

The warming rates for the daily maximum temperatures exceed a factor of 2 over large areas in central and southern Europe (see Fig. 9a). In contrast, the minimum temperatures increase at a comparatively lower rate. This also indicates the presence of an enhanced diurnal temperature range during extreme heat events, as discussed above. Consistent with the values shown in Fig. 8a, the warming rates for a “neutral” period in early summer are much lower and below a factor of 1 across central Europe (Fig. S17a).

A closer look at the response of the 2 m temperature on 25 July reveals that the warming rates during the peak of the event reached a factor of 3.0 east of the heatwave's core, underscoring the fact that these areas became up to 12 °C warmer in the +4 K climate compared to the pre-industrial conditions in our simulations (see Fig. 9b). The mean warming rates of the EUR-12 ensemble show a similar pattern: the black contour in Fig. S16, encompassing the core of the heatwave in the present-day climate, does not coincide with the 2.5 K K⁻¹ contour of the temperature response. Accordingly, the spatial extent of the heatwave is expected to be subject to a rapid increase in future climates. As shown in Fig. S9, the area affected by maximum temperatures over 40 °C corresponding to 25 July grows with the global-warming level at an approximate rate of 270 000 ± 1300 km² K⁻¹, slightly accelerating towards the +4 K climate. The lack of warming over the British Isles may be associated with the fact that this area is surrounded by sea and located at the edge of the subtropical ridge triggering this event, making it less affected.

To compare the scaling of temperatures during the extreme event to the mean scaling, we estimated the average monthly response of the 2 m temperatures over Europe for five simulated summers from 2018 to 2022. According to Fig. 10, the warming rates are close to 1 in central Europe in June, increase to 1.5 in July, and approach 2.0 in August. Such an intensification of the temperature response towards late summer indicates a higher risk of heatwave development in a warmer world. This aspect was investigated in Hundhausen et al. (2023), where, based on an ensemble of regional high-resolution climate simulations, the probability of large heatwaves was found to gradually increase during the summer.

Comparing the mean scaling in Fig. 10 to the warming rates in Fig. 9b, we find that the scaling of the maximum temperature at the peak of the July 2019 heatwave is nearly twice as high in western France and eastern Germany. Nevertheless, it is worth mentioning that all the years of the simulated period fall within the span of the European multi-year drought and heat event of 2018–2022 (Knutzen et al., 2025). Thus, the average warming rates in Fig. 10 may differ from those for the years unaffected by severe drought conditions.

Figure 10 Mean warming rates for maximum 2 m temperature with respect to June, July, and August 2018–2022, based on the EUR-12 storylines.

4 Summary and discussion

In this study, we follow an event-based storyline approach using a GCM–RCM–CPM model chain to analyse the thermodynamic response of the European summer heatwaves of 2019 to global warming. We obtained our storylines using global spectrally nudged AWI-CM1 simulations and subsequent dynamical downscaling with the regional model ICON-CLM to resolutions of 12 km (EUR-12) and 3 km (GER-3).

The outcomes of the simulations allow us to answer the three key research questions as follows:

1. How accurately can a regional event-based storyline simulation represent a recent event, and what improvement does it offer compared to the global spectrally nudged storyline simulation?

   The daily spatial and temporal patterns of the 2 m temperature fields obtained with the AWI-CM1–ICON model chain show good agreement with ERA5, E-OBS, and DWD station observations for the summer of 2019. Compared to the driving AWI-CM1 simulations, the dynamical downscaling significantly reduced the RMSD in 2 m temperature over most of Europe (by about 1.5 °C in central Europe). The GER-3 simulation provides more spatial detail than the EUR-12 simulation for daily maximum temperatures; however, it does not necessarily improve the representation of daily minimum temperatures.

2. What was the effect of climate change on the 2019 European heatwave based on regional and convection-permitting ICON-CLM simulations?

   Based on the simulations for our 2019 case study with different thresholds for global warming, the peak temperatures of the July heatwave increased significantly beyond the underlying global-warming level, with the magnitude of the temperature response depending on the location. The increments in daily maximum temperature in the +4 K climate, relative to the pre-industrial climate, vary between 8 °C (in the centre of the heatwave) and 12 °C (to the east of the heatwave), leading to an increased spatial extent of the heatwave in a warmer world. In the context of the present-day climate, where we took the 40 °C isoline as a benchmark, the affected area was encompassed by the 45 °C isoline in the +4 K world. Consequently, the areas affected by temperatures exceeding 40 °C experienced a significant expansion, exceeding 1×10⁶ km² in the +4 K storyline. When considering the time series of 2 m temperature over an area ranging from 48–51° N, 6–10° E, for all five storylines, the spread between the curves appeared to be larger in the mid- and late-summer periods compared to early summer. This aligns with previous findings suggesting the intra-seasonal dependence of anthropogenic warming (intensification towards late summer; see e.g. Hundhausen et al., 2023), based not only on the nudged storylines but also on the free CMIP6 runs (Sánchez-Benítez et al., 2022).

3. What local-to-regional extreme temperature scaling occurs in response to global warming for an event like the 2019 heatwave?

   Our findings reveal a linear dependence of the 2 m temperature response on the global-warming level, with the derived warming rate determined by the slope of the linear regression, demonstrating spatial and intra-seasonal variations. Quantifying the smaller spacing between the temperature curves mentioned above, the warming rates over the studied area in early summer fluctuate around 1, indicating that regional warming aligns with global warming during this period. However, in July and August, the warming rates for daily maximum temperatures vary between 1.5 and 2.5, reaching higher values with each successive heatwave during the study period (see Fig. 8d).

Furthermore, the broadening of the warming-rate peak during the July heatwave indicates an extension of the heatwave's duration in a warmer world. We also observe the broadening of the diurnal temperature range in future climates, which is indicated by much lower warming rates for the minimum temperature compared to those for the daily maximum. This difference does not occur in early summer and disappears by late September. Considering that the first heatwave in late June modified the soil moisture for the rest of the summer (Sousa et al., 2020; Sánchez-Benítez et al., 2022), this case exemplifies the dependence of global-warming amplification on the event-specific regional evolution of thermodynamic conditions.

On 25 July, the response of the maximum 2 m temperature (warming rate) reached a factor of 3 in some areas (Fig. 9b). However, the highest warming rates were not located over the heatwave centre but instead shifted eastwards (see Figs. S12 and S16). The areas less affected by unprecedented temperatures but still located within the influence of the event triggered by the subtropical ridge tend to experience stronger warming in a warmer world. Thus, along with the increasing duration, the area affected by the heatwave expands. While the GER-3 domain effectively captures the regions with the strongest temperature responses, it does not fully encompass the core area affected by the July heatwave. We admit that for a more comprehensive kilometre-scale investigation of this heatwave, this domain should have covered a larger fraction of France. Still, as this study is part of the InnoPool SCENIC (Storyline Scenarios of Extreme Weather, Climate, and Environmental Events along with their Impacts in a Warmer World) project (Helmholtz Changing Earth, 2024), which focuses on the impacts of climate change within Germany, our CPM simulations provide the necessary insights for this context.

Several factors may have influenced the spatially variable magnitude of the 2 m temperature response, as seen both for a single member (Fig. 9) and for all ensemble members (Fig. S16). One possible explanation for this behaviour is the heterogeneous response of the coupling between soil moisture and temperature within and outside the heatwave's core (see e.g. Gevaert et al., 2018; Miralles et al., 2014). We hypothesise that the overall amplification of the warming rates during heatwave events and the extension of the diurnal temperature range may have been exacerbated due to soil–atmosphere feedback. However, a dedicated and detailed analysis would be necessary to demonstrate or dismiss this hypothesis as other factors, such as small changes in atmospheric dynamics, may have also played an important role.

The same event was analysed by Vries et al. (2024) using a PGW approach, yielding similar temperature responses of 1.5 to 2.5 K K⁻¹ during heatwaves. Both studies show an increased response of maximum temperatures compared to minimum temperatures. However, in contrast to our findings, Vries et al. (2024) observed no significant response dampening within the heatwave core. Nevertheless, both studies indicate that areas in France impacted by extreme temperatures on 25 July do not show higher scaling than the surrounding areas. While Vries et al. (2024) focus on the southeastern Netherlands, where higher temperatures yield stronger responses (see Fig. S10 in Vries et al., 2024), this relationship may not be directly applicable to central France for events as extreme as the July 2019 heatwave, given the region's distinct climate, and further investigation is needed to determine whether such scaling applies. The complementary results from both studies enrich our understanding of heatwave dynamics and provide a broader context for future investigations.

5 Conclusions

The aim of this study was to provide a regional perspective on global spectrally nudged storylines for the summer heatwaves of 2019 in central Europe. We addressed the unfolding of these heatwaves on regional-to-local spatial scales and followed the evolution of near-surface temperatures throughout the entire summer season in five dynamical analogues of the summer of 2019 using a dynamical downscaling approach. We observed that the late June heatwave triggered higher warming rates and an extension of the diurnal temperature range in a warmer world, continuing for the rest of the summer. Additionally, we identified higher warming rates over regions east of the centre of the July heatwave, as well as a broadening of the warming-rate peaks associated with both 2019 heatwaves. This demonstrates that our approach not only allows for the estimation of the possible impacts of extreme heat events in a warmer world but also facilitates the investigation of the mechanisms and conditions that lead to different rates of response to background warming.

Our regional storylines can be used to drive hydrology and land surface models, as well as other impact models, providing relevant information for developing adaptation measures. Finally, the insights gained from storyline-based regional impact studies are more tangible than probabilistic estimates and, thus, have the potential to raise public awareness about the significance of the effects of climate change at the community level.