Exploring the uncertainties in the aviation soot–cirrus effect

A global aerosol–climate model, including a two-moment cloud microphysical scheme and a parametrization for aerosol-induced ice formation in cirrus clouds, is applied in order to quantify the impact of aviation soot on natural cirrus clouds. Several sensitivity experiments are performed to assess the uncertainties in this effect related to (i) the assumptions on the ice nucleation abilities of aviation soot, (ii) the representation of vertical updrafts in the model, and (iii) the use of reanalysis data to relax the model dynamics (the so-called nudging technique). Based on the results of the model simulations, a radiative forcing from the aviation soot–cirrus effect in the range of −35 to 13 mW m⁻² is quantified, depending on the assumed critical saturation ratio for ice nucleation and active fraction of aviation soot but with a confidence level below 95 % in several cases. Simple idealized experiments with prescribed vertical velocities further show that the uncertainties on this aspect of the model dynamics are critical for the investigated effect and could potentially add a factor of about 2 of further uncertainty to the model estimates of the resulting radiative forcing. The use of the nudging technique to relax model dynamics is proved essential in order to identify a statistically significant signal from the model internal variability, while simulations performed in free-running mode and with prescribed sea-surface temperatures and sea-ice concentrations are shown to be unable to provide robust estimates of the investigated effect. A comparison with analogous model studies on the aviation soot–cirrus effect show a very large model diversity, with a conspicuous lack of consensus across the various estimates, which points to the need for more in-depth analyses on the roots of such discrepancies.

1 Introduction

The aviation sector contributes about 2.4 % of the global anthropogenic CO₂ and is one of the fastest growing anthropogenic sectors, which makes it one of the key targets for mitigating the anthropogenic impact on climate (Lee et al., 2010; Grewe et al., 2017; Lee et al., 2021). In the last decades, civil aviation has experienced a steady growth in activity, resulting in increasing CO₂ emissions at an average rate of 2 % yr⁻¹ between 1970 and 2012, further accelerating in recent years (2013–2018) to a rate of 5 % yr⁻¹ (Lee et al., 2021). Most of the Shared Socioeconomic Pathways (SSPs; Riahi et al., 2017), developed in the context of the assessments by the Intergovernmental Panel for Climate Change (IPCC), project this growth in aviation emissions to continue until at least 2040 (Gidden et al., 2019), although the unexpected worldwide outbreak of the SARS-CoV-2 pandemic in the early months of 2020 temporarily reduced such an increase (e.g. Forster et al., 2020; Gettelman et al., 2021).

In addition to the well-understood impact of CO₂ and the related mitigation measures (Fuglestvedt et al., 2008; Dahlmann et al., 2016), aircraft also emit a number of non-CO₂ components, whose climate impact is still uncertain (Grewe et al., 2017). This concerns, for instance, the role of nitrogen oxides (${\mathrm{NO}}_x=\mathrm{NO}+{\mathrm{NO}}_{\mathrm{2}}$), which control ozone formation and affect methane lifetime (Grewe et al., 2019), aerosol particles, and their interactions with clouds (e.g. Gettelman and Chen, 2013; Righi et al., 2013; Penner et al., 2018), as well as the formation and growth of contrails and contrail cirrus (Burkhardt and Kärcher, 2011; Chen and Gettelman, 2013; Bock and Burkhardt, 2016).

Among these various aviation effects, the impact of aviation soot on natural cirrus clouds has gained attention in recent years due to its potentially large climate impact, possibly exceeding the contribution of most of the aforementioned components, including CO₂ (Hendricks et al., 2005; Penner et al., 2018). Cirrus clouds cover about 30 % of the globe and have an overall warming impact on the Earth radiative balance, as on average their longwave warming effect dominates over shortwave cooling (Hartmann et al., 1992; Hong et al., 2016; Chen et al., 2000; Gasparini and Lohmann, 2016; Heymsfield et al., 2017). The aerosol-induced formation of ice crystals (ICs) in cirrus clouds can occur either via homogeneous freezing of supercooled liquid solution aerosols or via heterogeneous freezing on the surface of a so-called ice-nucleating particle (INP; Vali et al., 2015). During the formation of cirrus clouds, both mechanisms can occur, and their competition for the available supersaturated water vapour has crucial effects on the microphysical structure of these clouds, since it controls the number concentration and the size of ICs and hence their optical and radiative properties (Kärcher, 2017). Several aerosol types have been shown to act as INPs at cirrus temperatures, including mineral dust and soot (Möhler et al., 2005; Hoose and Möhler, 2012; Cziczo et al., 2013; Kanji et al., 2017). Since aircraft soot emissions are directly released at cirrus altitudes (i.e. above ∼400 hPa), several studies focused on their impact on cirrus properties and the resulting climate effect (Hendricks et al., 2005, 2011; Liu et al., 2009; Penner et al., 2009; Gettelman and Chen, 2013; Zhou and Penner, 2014; Kärcher, 2017; Penner et al., 2018; Zhu and Penner, 2020; McGraw et al., 2020).

Nevertheless, the above-mentioned studies have indicated very large uncertainties on the magnitude, and even on the sign, of the impact of aviation soot on the radiative forcing (RF) exerted by natural cirrus. Part of these uncertainties derives from the assumptions on the still poorly understood ice-nucleating properties of aviation soot. This concerns in particular the critical saturation ratio (S_crit) at which soot can initiate freezing and the active fraction (f_act) of the soot particle population that can act as INPs. Applying a global aerosol model coupled with a two-moment ice microphysical scheme (CAM3-IMPACT), Liu et al. (2009) estimated an aviation soot effect of −110 mW m⁻², assuming a critical ice saturation ratio for soot nucleation S_crit=1.2–1.3 but a positive RF of 260 mW m⁻² if soot is considered to be a worse INP with S_crit=1.4. In both cases, a very high active fraction f_act=100 % was assumed for aviation soot. Using the same model, Penner et al. (2009) applied two different ice nucleation parametrizations (Kärcher and Lohmann, 2002; Liu and Penner, 2005) and quantified an aviation soot–cirrus RF of −161 and −124 mW m⁻². In the first case, they assumed S_crit=1.3 for soot, while in the latter the occurrence of heterogeneous freezing on soot particles was determined as a function of temperature. In both cases, the whole soot population was assumed to act as INPs (i.e. f_act=100 %). Using the ECHAM4 model and the Kärcher et al. (2006) parametrization for ice nucleation which is also applied in the present study, Hendricks et al. (2011) found no statistically significant climate effect when assuming that 10 % of aviation soot acts as INPs, with S_crit=1.2. No statistically significant effects were also reported by Gettelman and Chen (2013), who applied the CAM5 model with the assumption that soot has similar ice-nucleating properties as mineral dust (i.e. S_crit=1.2–1.3) but with lower (and presumably more realistic) active fractions. Zhou and Penner (2014) again used the CAM5 model with the Liu and Penner (2005) ice nucleation parametrization but provided an explicit calculation of the soot active fraction by considering its preprocessing in contrail cirrus. This resulted in 0.6 % of aviation soot being an efficient INP, leading to a RF quantified in the range of −350 to 90 mW m⁻², depending on the assumptions on the background sulfate and dust concentrations, both competing with soot for the formation of ICs in cirrus clouds. Using the same model with the Kärcher et al. (2006) parametrization for ice nucleation and an improved scheme for subgrid-scale vertical updrafts, Penner et al. (2018) simulated an aircraft effect on cirrus clouds of −200 mW m⁻², considering S_crit=1.35 and again f_act=0.6 % as a result of soot preprocessing in contrails. A subsequent study, using the Community Earth System Model (CESM) featuring a hybrid ice nucleation scheme combining the best features of the Liu and Penner (2005) and Kärcher et al. (2006) parametrizations, updated this estimate to −140 mW m⁻², also assuming S_crit=1.35 (Zhu and Penner, 2020). In the most recent work on this subject, McGraw et al. (2020) used the CESM2 with the cirrus nucleation scheme by Barahona and Nenes (2009) and found that aircraft soot can perturb cirrus clouds in the Northern Hemisphere, with effects both in the shortwave and in the longwave RF, but they were not able to extract a statistically significant signal from their model simulations.

This brief literature review shows that the existing quantifications of the aviation soot–cirrus effect obtained with different model approaches lead to very different results, ranging from statistically non-significant effects to potentially large RF values. The reason for this lack of consensus lies in the high complexity of the physical processes controlling this effect, which are hard to constrain by measurements and difficult to represent in global models (Kärcher, 2017). A few studies (Koehler et al., 2009; Mahrt et al., 2020) have shown that soot particles released by the combustion process in aircraft turbines could act as INPs at various ranges of atmospheric temperatures and ice supersaturations, but measurements are still limited. To further complicate this picture, ice formation in the cirrus regime is also crucially driven by atmospheric dynamics and depends on the occurrence of vertical updrafts and on their strength (Kärcher and Podglajen, 2019; Kärcher et al., 2019). In spite of recent progress on the modelling side (Podglajen et al., 2016), the representation of vertical updrafts in global models is still subject to considerable approximations and does not capture small-scale variability in sufficient detail.

In this work, we apply the EMAC model (ECHAM/MESSy Atmospheric Chemistry; Jöckel et al., 2010) with the aerosol submodel MADE3 (Modal Aerosol Dynamics for Europe, adapted for global applications, third generation; Kaiser et al., 2014, 2019) in a recently developed configuration (Righi et al., 2020) which includes a detailed parametrization for aerosol-induced ice formation in cirrus clouds (Kärcher et al., 2006; Kuebbeler et al., 2014). Rather than attempting to provide a single estimate of the aviation soot–cirrus effect, the goal of this study is to explore the uncertainties related to the microphysical and dynamic aspects of this effect, in order to provide a realistic, albeit broad, range of possible values for the resulting climate impact. The microphysical analysis focuses on the ice-nucleating properties of aviation soot, based on the results of the laboratory measurements reported in the literature. We show that the variation in both the critical saturation ratio S_crit for ice nucleation and the active fraction f_act of aviation soot can have a significant impact not only on the magnitude, but also on the sign of the resulting RF. We attempt to relate these changes to the underlying physical processes as represented in the model. Furthermore, we analyse the role of the model's representation of vertical updraft by means of mechanistic studies in which a very simple representation of such updraft is implemented in order to explore the sensitivity of the relevant microphysics to the model dynamics. This parametric approach has proven useful in a previous climate impact study (Righi et al., 2013), where we analysed the uncertainties related to the assumptions on the size distribution of aerosol particles from different transport emission sources, including aviation. On a more general level, such parametric approaches were successfully used, for instance, to constrain the uncertainties on the microphysical properties of warm clouds (Lee et al., 2013) and on the aerosol indirect effect (Carslaw et al., 2013; Regayre et al., 2020).

This paper is organized as follows: in Sect. 2 the processes controlling the aviation soot–cirrus effect and its large uncertainties are discussed, and the results of available laboratory studies are summarized. The model setup, which is largely based on Righi et al. (2020) but has been further improved here, and the performed numerical experiments are described in Sect. 3. The results are presented and discussed in Sect. 4, and we summarize the main outcomes of this study in Sect. 5.

2 Uncertainties in the soot–cirrus effect

As mentioned in the introduction, aerosol-induced ice formation in the cirrus regime ($T\mathit{\lesssim }-\mathrm{37}{}^{\circ }$C) can occur either via homogeneous freezing or by heterogeneous freezing on the surface of an INP. The latter process usually requires a lower critical supersaturation over ice and can therefore occur prior to the onset of homogeneous freezing, attenuating or even inhibiting the direct freezing of supercooled liquid solution aerosol. This has of course important consequences on the microphysical properties of cirrus, since it affects the number concentrations and size of ICs and, in turn, the lifetime and the radiative properties of the clouds. Cirrus clouds generally exert large RFs, both in the shortwave and in the longwave spectrum, with the latter being slightly larger, which results in an overall warming effect (Hartmann et al., 1992; Hong et al., 2016; Chen et al., 2000; Gasparini and Lohmann, 2016). If a cirrus cloud is dominated by homogeneous freezing, adding more INPs typically results in a decrease of IC number concentration (ICNC) and a corresponding increase in their size (the so-called negative Twomey effect; Kärcher and Lohmann, 2003). Adding INPs to a cirrus cloud where heterogeneous freezing already dominates, on the other hand, could result in a further increase of ICNC and a decrease in crystal sizes. As shown by previous studies (e.g. Zhang et al., 1999), an increasing IC size reduces the longwave cloud forcing (i.e. less warming) but also the shortwave cloud forcing in absolute terms (i.e. less cooling), and their combination can be either a net warming or a net cooling depending on the cloud ice water content. Other effects, such as a more efficient sedimentation of less abundant, larger ICs or an increased deposition of water vapour in the presence of more efficient INPs (dehydration; Jensen et al., 2001, 2013), add even more complexity to this picture and to the interpretation of the model results. For a reliable quantification of the effect of aviation soot on cirrus clouds and its climate impact, it is therefore essential not only to have a reliable estimate of the ice-nucleating properties of aviation soot, but also of the dynamic processes that control the background state of natural cirrus clouds in the model. Cirrus formed in slow (fast) updrafts are usually dominated by heterogeneous (homogeneous) nucleation, resulting in lower (higher) concentrations of ICs and larger (smaller) IC sizes (Kärcher and Lohmann, 2002; Krämer et al., 2016, 2020). Hence, the way aviation soot INPs can impact these clouds depends not only on their ice nucleation ability, but also on the dynamical background conditions. In the present study we therefore focus on both the microphysical and the dynamic aspects, by analysing the uncertainties related to the assumptions on aviation soot INP characteristics and exploring the impact of different (simplified) representations of the dynamic forcing in the model. We stress again that the main goal of this study is not to provide an updated estimate on the aviation soot–cirrus effect but to explore its sensitivity to aviation soot microphysics and, to some extent, to the underlying model dynamics.

Figure 1 summarizes the ice-nucleating properties of different soot types retrieved from the literature, including cloud-processed soot. These properties are given in terms of ice saturation ratio S_ice and active fraction f_act and are compared with the assumptions made in the modelling studies on the aviation soot–cirrus effect reviewed in the Introduction. Note that models usually assume a critical ice saturation ratio S_crit at which ice nucleation takes place, while in the laboratory measurement a wide spectrum of S_ice values is explored. Furthermore, the experimental results show a clear temperature dependence of the soot ice-nucleating properties, which most of the model parametrizations do not take into account (an exception here is the Liu and Penner (2005) parametrization). The measurements mostly report medium to low ice nucleation efficiency by different soot types (S_ice≳1.3) and activation fractions hardly above 10 %. The modelling studies, on the other hand, tend to represent soot as a better INP, especially in terms of critical supersaturation ratio. Recent laboratories studies (Marcolli, 2017; Mahrt et al., 2020) support indeed higher ice nucleation efficiency for soot particles which experienced cloud processing, for example in contrails, but this applies only to particles with diameters of about 400 nm, which implies that only a very low fraction of aviation soot can effectively be active as INPs in the upper troposphere. This was also the conclusion by Kärcher et al. (2021): in a recent process-oriented analysis using a high-resolution column model based on the parametrization by Marcolli et al. (2021), they suggested that less than 1 % of aviation soot particles can lead to the formation of ice crystals in competition with homogeneous freezing. To explore this parameter space in sufficient detail, also reproducing the assumptions of previous model studies, we therefore perform nine sensitivity simulations in this study, varying S_crit from 1.2 to 1.4 and f_act from 0.1 % to 10 %. These assumptions are marked with the black star symbols in Fig. 1 and result in nine combinations of these two parameters. Values of S_crit>1.4 would mostly exceed the homogeneous freezing threshold at relevant cirrus temperatures and are therefore not worthy to be explored for the scope of the present study.

Figure 1 Summary of ice-nucleating properties of soot measured in different laboratory studies (Möhler et al., 2005; Koehler et al., 2009; Crawford et al., 2011; Kanji et al., 2011; Chou et al., 2013; Kulkarni et al., 2016; Mahrt et al., 2018; Nichman et al., 2019; Mahrt et al., 2020) compared with the values assumed in model studies on the impact of aviation soot on cirrus (Hendricks et al., 2005; Penner et al., 2009; Liu et al., 2009; Hendricks et al., 2011; Gettelman and Chen, 2013; Zhou and Penner, 2014; Penner et al., 2018; Zhu and Penner, 2020). The parameters explored in this study are symbolized with stars. The vertical dashed lines show the homogeneous ice nucleation threshold at 210, 220, and 230 K.

All previous model-based investigations on the aviation soot–cirrus effect considered approximated representations of the vertical updraft and its subscale variability: a common approach uses the square root of the turbulent kinetic energy (TKE) as a proxy for such variability (Lohmann and Kärcher, 2002). Later studies (Kuebbeler et al., 2014) also included the contribution of orographic gravity waves generated on the lee of mountain ranges (Joos et al., 2008), while recent models (Penner et al., 2018) used methodologies based on measurements (Podglajen et al., 2016) to consider the contribution of gravity waves to the vertical velocity. In this study, the aviation soot–cirrus effect is quantified with the EMAC-MADE3 model (see Sect. 3), which follows the TKE approach for the subscale vertical velocity, also considering the impact of orographic waves in relevant regions. To further investigate how different cirrus regimes may react to the perturbation represented by aviation soot INPs, we also consider an idealized representation of the vertical velocity. We prescribe constant values of the vertical velocity in the range from 2 to 50 cm s⁻¹, to explore the full range of possible updraft regimes, including both slow and fast updrafts (Kärcher et al., 2006; Krämer et al., 2020). While this approach is of course idealized, it offers the possibility to separate the microphysical from the dynamic effect by artificially introducing a spatially uniform dynamic regime, thus allowing the competition among the different INPs considered by the model to be interpreted purely in terms of their microphysical properties. Furthermore, it enables the investigation of INP effects under possible regimes, not covered by the TKE and orographic gravity wave approaches mentioned above.

3 Model description and simulations

The global aerosol model EMAC-MADE3 used in this study has been extensively evaluated in Kaiser et al. (2019), focusing in particular on the representation of aerosol particles and their global distribution. In this work, we apply an advanced version of the model, including couplings of the aerosol submodel MADE3 with the radiation and cloud schemes. The model includes a two-moment cloud microphysical scheme based on Kuebbeler et al. (2014), which features a parametrization for aerosol-driven ice formation in cirrus clouds following Kärcher et al. (2006). This parametrization considers the competition among the different ice formation processes in increasing order of critical saturation ratio, i.e. from the most efficient heterogeneous freezing process to homogeneous freezing. The increase of supersaturation driven by vertical updrafts proceeds until it is compensated by the loss of water vapour due to the growth of ice crystals (either freshly formed or pre-existing from the previous model time step). For large enough supersaturations, both heterogeneous and homogeneous freezing can take place, and their competition is accounted for in the model.

This configuration has been tuned and evaluated with respect to the main cloud and radiation variables in Righi et al. (2020). The work by Righi et al. (2020) demonstrated the overall ability of the model to reproduce the most important cirrus properties, like ice water content (IWC) and ICNC, when compared with in situ measurements (Krämer et al., 2016; Voigt et al., 2017; Krämer et al., 2020), a necessary requirement in the context of the present work.

With respect to the EMAC-MADE3 version documented in Righi et al. (2020), some technical changes and a few improvements specific to the application for this study have been introduced. In particular, the time integration of the cloud prognostic variables, including cloud droplet number concentration (CDNC) and ICNC, has been improved, to make it consistent with the approach used by the other submodels which are part of the EMAC model system. This change required a retuning of the model configuration, which now uses slightly different values for the tuning parameters with respect to the Righi et al. (2020) setup. More specifically, the enhancement factor γ_r for the rate of rain formation by autoconversion has been reduced from γ_r=8 to γ_r=4, and the minimum CDNC (a threshold introduced in the cloud scheme to avoid unrealistically low concentrations in pristine regions) has been increased from 20 to 50 cm⁻³. This value is close to the 40 cm⁻³ value used by Neubauer et al. (2019) in the ECHAM6 model and is supported by typical CDNC values measured by satellites in pristine regions (Bennartz and Rausch, 2017; Grosvenor and Wood, 2018). The overall model performance after this correction is similar to Righi et al. (2020), and the main conclusions of that study remain valid: the model provides a good representation of cloud cover, CDNC, precipitation, and cloud radiative effects, while it is still affected by a relatively high liquid water path (LWP; at the high end of the observed range). The cirrus specific quantities, IWC and ICNC, show a slight improvement. Furthermore, the model now simulates an effective RF due to the total anthropogenic aerosol effect of −1.16 W m⁻², in better agreement with the recent assessment by Bellouin et al. (2020), who reported a range −1.6 to −0.6 W m⁻² (68 % confidence interval). The tuned model is characterized by a radiative balance of 3.4 mW m⁻² when run in nudged mode, i.e. by relaxing the meteorology towards reanalysis data. The same configuration in free-running mode has a radiative balance of 0.9 mW m⁻² (see Righi et al., 2020, for a detailed discussion on the impact of nudging on the model radiative balance).

Another relevant improvement to the model configuration applied here is the introduction of an additional tracer BCtag to which the soot emissions from the aviation sector are assigned. The BCtag tracer is distributed into the same six modes as the standard BC tracer of MADE3, namely Aitken, accumulation, and coarse mode, each with insoluble and mixed states. The BC and BCtag tracers have the same physical properties and undergo exactly the same processes in the model but allow for different ice-nucleating properties between background and aviation soot in the cirrus parametrization. The ice-nucleating properties of mineral dust and background soot are the same as in Righi et al. (2020), namely $S_{\mathrm{crit}}=\mathrm{1.1}-\mathrm{1.2}$ with a temperature-dependent active fraction for mineral dust in the deposition mode, S_crit=1.3 and f_act=5 % for mineral dust in the immersion mode (Kuebbeler et al., 2014), and S_crit=1.4 and f_act=0.25 % for background soot (Hendricks et al., 2011). The S_crit and f_act parameters for aviation (and in part also background) soot are explored in more detail in the dedicated sensitivity studies, as discussed in Sect. 2. To avoid confusion, we note here that the BC and BCtag tracers of MADE3 actually refer to black carbon, i.e. an aerosol type composed only of carbon, but we are using the term soot in this paper for consistency with most of the literature on aviation effects, although these definitions are not fully consistent (see Petzold et al., 2013, for a more detailed discussion on this terminology).

To calculate the number concentration of INPs for the different types, we use the same approach as Righi et al. (2020), while for the newly introduced BCtag tracer we derive the number concentration from the tracer mass, by assuming aviation soot to follow the bimodal size distribution measured by Petzold et al. (1999) in the plume of a B737-300 aircraft. This distribution is characterized by median diameters of 25 and 150 nm and geometric standard deviations of 1.55 and 1.65, for the Aitken and accumulation modes, respectively. The same size distribution parameters were used in Righi et al. (2013) to characterize particle number emissions from aviation. Introducing the BCtag tracer has the advantage that a lower number of simulations needs to be performed to isolate the impact of aviation soot on cirrus clouds, since only two experiments are required for that, i.e. with and without the effect of the BCtag INPs in the cirrus parametrization. The difference between these two experiments hence provides an estimate of the climate impact of aviation soot on natural cirrus clouds while excluding the effects of the interactions with clouds at lower levels (e.g. sulfate impact on liquid clouds). The statistical significance of this estimate is also improved with respect to an approach where no tagging of aviation soot is included, since in that case four experiments would be required (with and without aviation, with and without soot impact on cirrus) to isolate the effect. Another advantage of this tagging approach is that different ice nucleation abilities can be assumed for aviation soot and background soot, i.e. soot from background sources. The aviation soot–cirrus effect is estimated by calculating the difference between a given simulation and a baseline experiment (BASE), where aviation soot (i.e. the BCtag tracer) is not considered as INPs in the cirrus parametrization. A paired sample t test is applied to verify the null hypothesis that the annual mean values of a given quantity (e.g. RF) are identical in the two simulations (with and without aviation soot impact on cirrus). We express the response of the test in terms of confidence level, i.e. 100(1−p), where p is the p value. Unless otherwise specified, we regard the results as statistically significant when the null hypothesis can be rejected at a confidence level larger than 95 % (p<0.05).

Finally, to provide estimates of the aviation effects closer to the present-day conditions, we updated the model emissions setup for anthropogenic and biomass burning (or open burning) sources of aerosol and aerosol precursor species to the recent CMIP6 inventory for the year 2014 (van Marle et al., 2017; Hoesly et al., 2018; Feng et al., 2020). Aviation BC emissions in this inventory amount to about 10.5 Gg in 2014. Natural emissions (e.g. mineral dust and volcanoes) are considered as in Righi et al. (2020), where mineral dust and sea spray emissions were calculated as a function of the local wind speed. For consistency, the prescribed mixing ratios of greenhouse gases used as input to the radiation scheme of EMAC (RAD; Dietmüller et al., 2016) are also updated to the 2014 values provided in Meinshausen et al. (2017) for the CMIP6 project, while for ozone we use the geographically resolved mixing ratios generated in support of CMIP6 by Hegglin et al. (in prep.; see also Hegglin et al., 2016) for the same year. The simulations performed in this study are summarized in Table 1. They are all characterized by a T42L41 configuration, corresponding to a horizontal resolution of about ${\mathrm{2.8}}^{\circ }\times {\mathrm{2.8}}^{\circ }$ and 41 non-equidistant vertical layers. The simulated time period covers the years from 2000 to 2015, with the year 2000 as spin-up and not included in the analysis. Unless otherwise specified, the model meteorology (temperature, winds, and logarithm of the surface pressure) is nudged towards the ERA-Interim reanalysis data (Dee et al., 2011) of the European Centre for Medium-Range Weather Forecasts (ECMWF) for the same time period. The relaxation times for nudging are 6 h for wind vorticity, 24 h for surface pressure and temperature, and 48 h for wind divergency. Further sensitivity experiments are performed in free-running mode to analyse the effect of nudging on the results and cover a period of 30 years, using prescribed climatological mean (2000–2009) sea-surface temperature (SST) and sea-ice concentration (SIC) from the Met Office Hadley Centre dataset (HadISST; Rayner et al., 2003).

Table 1 Summary of the EMAC-MADE3 simulations performed in this study. S_crit and f_act refer to the ice-nucleating properties of aviation soot, while $S_{\mathrm{crit}}^{\mathrm{bg}}$ and $f_{\mathrm{act}}^{\mathrm{bg}}$ indicate the ones of background soot (i.e. from other anthropogenic and biomass burning sources). Every simulation includes an extra spin-up year which is not considered for the analysis. The nudged simulations use meteorological reanalysis data for the period 2001–2015.

The model's representation of vertical velocity follows the approach of Kuebbeler et al. (2014) as described in Righi et al. (2020), but given its importance for the present study, we recall the main aspects here. The vertical velocity in the model is described as the sum of a large-scale and a subscale component. The latter is approximated according to Lohmann and Kärcher (2002) as $\mathrm{0.7}\sqrt{\mathrm{TKE}}$, while in the vicinity of mountain ranges the subscale contribution of orographic waves to the vertical velocity is considered instead, following the parametrization of Joos et al. (2008). The vertical velocities from these three components, as simulated by EMAC-MADE3, are depicted in Fig. 2 for the BASE experiment, considering the 95th percentile of the vertical velocity distribution above 400 hPa. This allows the largest fluctuations in the vertical velocity distribution to be analysed for each component. The large-scale component (Fig. 2a) is usually characterized by small updrafts of the order of a few centimetres per second and with little spatial variation. The additional subscale contribution from the TKE (Fig. 2b) emerges mostly in the tropics, but it is generally below about 1 cm s⁻¹, with the exception of the Himalayan and Andean mountain ranges, where spikes of a few tens of centimetres per second can be expected. In these and other mountain regions, however, the model accounts for the contribution of orographic waves (Fig. 2c), which can induce large vertical velocity fluctuations of about 50 cm s⁻¹, mostly impacting the western side of North America, Europe, and large parts of Asia around the Himalayas, i.e. near the most prominent mountain ranges of the world. Such strong updrafts can induce large supersaturations, thus efficiently driving homogeneous freezing and enhancing cirrus formation.

Figure 2 The vertical velocity components considered by the model: (a) large-scale, (b) subscale using the turbulent kinetic energy (TKE) as proxy, and (c) subscale contribution of the orographic waves in the vicinity of mountain ranges. Each panel shows the 95th percentile of the vertical velocity distribution with respect to the time (over the 2001–2015 simulation period) and vertical (above 400 hPa) coordinates.

We further analyse the model mean state and, in particular, the impact of the vertical velocity on the ice formation processes in Fig. 3, depicting IWC and ICNC distributions in the BASE experiment, as well as the frequency of occurrence of homogeneous freezing events. The latter is estimated by considering the number concentration of homogeneously and heterogeneously formed ice crystals as diagnosed by the cirrus parametrization in the model. For each model grid box and output time step (11 h frequency; see Righi et al., 2020), we then assign a value 1 (0), where ICNC from homogeneous freezing is larger (smaller) than ICNC from heterogeneous freezing. By averaging this over the whole simulation time period (2001–2015), we obtain an estimate of the frequency of occurrence of homogeneous freezing events. IWC has its largest values (2–5 mg kg⁻¹) along the tropopause at all latitudes (Fig. 3a), while at 250 hPa the maxima are mostly in the tropics (Fig. 3b). A similar pattern is followed by ICNC (Fig. 3c, d), with slightly higher concentrations (200–500 L⁻¹) in the northern extra-tropics, possibly driven by anthropogenic and dust emissions contributing large concentration of INPs at these latitudes. The zonal mean profile of the homogeneous freezing fraction (Fig. 3e) shows a large contribution (≳50 %) of homogeneous freezing at lower levels in the polar regions and close to the tropopause at all latitudes, with larger values in the extra-tropics than in the tropics. As expected, the geographical distribution at 250 hPa (Fig. 3f) closely matches the pattern of vertical velocities; in particular the strong updrafts of several tens of centimetres per second induced by orographic waves (Fig. 2c) lead to enhanced homogeneous ice formation, with larger fractions found in the vicinity of the highest mountain ranges over the globe. Note that, due to the parametrization design, the analysis of the homogeneous freezing fraction only counts pristine crystals; hence it only applies to the formation process (disregarding for instance the different residence times of ice crystals of different sizes). The patterns depicted in Figs. 2 and 3 are important for the aviation effects investigated here and need to be taken into account for the interpretation of the results discussed in the next section.

Figure 3 Zonally averaged (a, c, e) and geographical distribution at 250 hPa (b, d, f) of IWC (a, b), ICNC (c, d), and the fraction of occurrence of homogeneous freezing events (e, f). All panels show results from the BASE experiment (see Table 1) and consider the multi-year average over the simulated period (2001–2015). IWC and ICNC include both cloudy and cloud-free model grid boxes.

4 Results and discussion

4.1 Geographical distribution of aviation soot and INPs

Before analysing the aviation soot radiative effects, we present in Fig. 4 the geographical distribution of aviation soot emissions from the CMIP6 dataset used here, together with the simulated mass and number concentration of aviation soot. Note that only a fraction of this number concentration is actually effective as INPs in the model, depending on the f_act fraction assumed in the different sensitivity experiments discussed above. We stress again that the aviation emissions are assigned to a tagged soot tracer in EMAC-MADE3, thus making it possible to track aviation soot in the model without the need for an extra sensitivity experiment with aviation emissions switched off, as in the standard perturbation approach applied, for example, in Righi et al. (2013). Not surprisingly, the aviation soot emissions are largest in the northern mid-latitudes (Fig. 4a), with maxima in the Northern Hemisphere at typical flight altitudes (200–250 hPa) and close to the surface, due to the climb and descent phases but also to the paths of short-range flights, which are mostly common over the continents. The geographical pattern at flight altitude (Fig. 4b) shows the major routes connecting the most populated areas of the world, in particular the North Atlantic flight corridor between Europe and the eastern United States, while the connections to eastern Asia are less marked in this inventory. The aviation soot mass concentration (Fig. 4c, d) follows a similar distribution, again with maxima at flight altitudes and close to the surface, with concentrations of the order of 0.2–0.5 ng m⁻³. Not surprisingly, much lower mass concentrations are found in the Southern Hemisphere, again with a maximum at flight altitudes of 0.05–0.1 ng m⁻³, at slightly higher level than in the Northern Hemisphere. This pattern is consistent with the results of a previous model study with the predecessor version of EMAC-MADE3 (Righi et al., 2013), where, however, lower concentrations of aviation soot were simulated. This could be due to the use of a different emission inventory (CMIP5 instead of CMIP6) and for earlier conditions (year 2000 vs. 2014) in that study but also to the introduction of the tagging method for aviation soot in the present work, which provides more accurate estimates. The difference in mass concentration between the Northern and Southern Hemisphere is evident in the geographical map at ∼250 hPa (Fig. 4d), where also a clear gradient from the tropics towards the polar regions is present. The number concentration of aviation soot particles (Fig. 4e, f) shows a similar pattern to the mass concentration, although with clearer features near the emission sources, such as the North Atlantic flight corridor, and the regions above central Europe and the eastern United States. Some flight corridors can also be identified in the Southern Hemisphere. The reason for such sharply structured patterns is that particles in the Aitken mode, which dominate total particle number, are characterized by a shorter lifetime due to particle–particle interactions, which effectively reduce their number concentration away from sources, while their mass is of course conserved. For comparison, Fig. S1 in the Supplement shows the number concentration of the other INPs considered by the model, namely mineral dust in the immersion and deposition mode, as well as background BC resulting from the other (non-aviation) emission sources.

Figure 4 Aviation soot emissions from the CMIP6 dataset used for the model simulations in this work (a, b). The other panels show the mass concentration (c, d) and the number concentration (e, f) of aviation soot simulated by EMAC-MADE3 (BASE experiment; see Table 1). Panels (a), (c), and (e) show zonal means, whereas (b), (d), and (f) show the ∼250 hPa level. All panels consider the multi-year average over the simulated period (2001–2015). Note that the number concentration of aviation soot is further multiplied by the active fraction f_act=0.1 %–10 % in the cirrus parametrization to obtain the number concentration of aviation soot INPs.

4.2 The aviation soot–cirrus effect

The RFs from the aviation soot–cirrus effect under the nine different assumptions for the ice-nucleating properties of aviation soot are presented in form of a matrix in Fig. 5. The EMAC-MADE3 simulations estimate this effect to be in the range of −35 to 13 mW m⁻² but lacking sufficient statistical significance for S_crit=1.3 and S_crit=1.4. The effect is negative (cooling) for low critical saturation ratios (higher nucleation efficiencies) and tends to increase towards a positive (warming) effect for medium to high critical saturation ratios (medium to low nucleation efficiencies). The statistical insignificance of the RF for some combinations of the parameters, however, makes it difficult to draw general conclusions from these overall numbers alone.

Figure 5 Multi-year average (2001–2015) top-of-the-atmosphere all-sky RF due to the effect of aviation soot on natural cirrus clouds, under different assumptions for the ice nucleation efficiency of aviation soot INPs (S_crit and f_act). The values in brackets within each box indicate the confidence level.

To facilitate the interpretation, we separate in Fig. 6a–d the different components of the RF, namely shortwave and longwave, all-sky and clear-sky. As expected, these graphs demonstrate how the aviation soot–cirrus effect actually results from the opposite changes in the shortwave (Fig. 6a) and longwave (Fig. 6b) all-sky RFs, which also have a much higher statistical significance than the total RF. The effect in the shortwave corresponds to a warming (i.e. a reduced cooling), meaning that aviation soot reduces the (background) cooling impact of natural cirrus clouds in the shortwave. This could be a manifestation of the negative Twomey effect (Kärcher and Lohmann, 2003): the additional INPs from aviation compete with homogeneous freezing for the available supersaturation and lead to the formation of fewer and larger ICs in these clouds, thereby decreasing their reflectivity and their shortwave cooling. This interpretation is consistent with the very small changes in the clear-sky part of the shortwave spectrum (Fig. 6c), and it is further supported by the general decrease in the homogeneous freezing fraction (Fig. 6e), although the changes in overall ICNC are very small and not statistically significant (Fig. 6f). The decrease in cloud frequency (Fig. 6h) could also contribute to these changes in the shortwave RF, especially at higher S_crit, where the competition between homogeneous and heterogeneous freezing which drives the negative Twomey effect can be expected to be less important (see again Fig. 6e). The warming effect in the shortwave spectrum is counteracted by a cooling (i.e. a reduced warming) in the longwave (Fig. 6b), due to the fact that larger ICs usually sediment more efficiently and reduce the lifetime of the cirrus clouds, resulting therefore in a decrease of their (background) warming effect and hence in a cooling. This is supported by the overall decrease in cloud frequency shown in Fig. 6h. The cooling (reduced warming) effect in the longwave is further enhanced when the nucleation efficiency of aviation soot is high (S_crit=1.2): in this case also the clear-sky RF (Fig. 6d) significantly contributes to the cooling. A possible aviation-induced (or aviation-enhanced) dehydration of the affected air masses is suggested here, resulting from the increased deposition of supersaturated water vapour on the very efficient INPs from aviation, which are then rapidly removed via sedimentation (Jensen et al., 2001, 2013). This reasoning is supported by the marked decrease in total water (sum of water vapour and ice water; Fig. 6g) and in cloud frequency (Fig. 6h): here small (≲1%) but statistically significant decreases in the simulations with S_crit=1.2 are evident. For lower ice nucleation efficiencies (S_crit≥1.3), the changes are about 1 order of magnitude smaller and mostly not statistically significant. We can therefore conclude that the aviation soot INPs can effectively enhance the dehydration of the upper troposphere and induce a statistically significant cooling effect when high nucleation efficiencies are assumed. For lower efficiencies (higher critical saturation ratios), the warming (reduced cooling) in the shortwave appears to be more important (possibly due to a reduction in cloud frequency; Fig. 6h) and leads to a slightly positive overall RF effect.

Figure 6 As in Fig. 5 but for (a) all-sky shortwave, (b) all-sky longwave, (c) clear-sky shortwave, and (d) clear-sky longwave top-of-the-atmosphere RFs. Panels (e–h) depict the aviation-soot-induced relative changes in ICNC from homogeneous freezing, ICNC, total water (as the sum of water vapour and ice water), and cloud frequency, respectively, all spatially averaged above 400 hPa and over cloudy and cloud-free model grid boxes.

The above analysis is based on globally averaged values, but large regional variations can be expected. This is because of the uneven geographical distribution of aviation soot and other INPs (Figs. 4 and S1) and of cirrus clouds but also because of the large spatial variations in the vertical velocity simulated by the model, as shown in Fig. 2. As discussed above, the latter effect is particularly important, as it controls the prevalence of homogeneous over heterogeneous ice formation regimes (Fig. 3e, f) and therefore the properties of natural cirrus clouds. The regional variations of the aviation soot RF are depicted in Fig. 7, where the RF effect is separated into three different latitude bands, namely southern extra-tropics (SH-Ext, 30–90^∘ S), tropics (Trop, 30^∘ S–30^∘ N), and northern extra-tropics (NH-Ext, 30–90^∘ N), and compared to the global values as presented in Fig. 5. For high nucleation efficiencies (S_crit=1.2; Fig. 7a, d, g), the net cooling effect is common to all regions and particularly strong in the Northern Hemisphere, where the aviation traffic is most intense, and aviation soot shows the largest values of mass and number concentrations (Fig. 4). This cooling effect in the northern extra-tropics decreases gradually with decreasing ice nucleation efficiency S_crit. For low efficiencies (S_crit=1.4; Fig. 7c, f, i), the warming effect from the southern extra-tropics and tropics dominates and results in a global warming effect. A reason for this pattern could be the smaller mean vertical velocities and the relatively clean background compared to the Northern Hemisphere. The additional aviation-emitted soot in this region could lead to enhanced heterogeneous nucleation and, due to the smaller cooling rates, to less homogeneous freezing, thus strengthening the negative Twomey effect.

Figure 7 As in Fig. 5 but averaged over different regions. Each of the panels (a–i) corresponds to a given assumption for ice nucleation efficiency of aviation soot particles (S_crit and f_act) as shown on the top. The values above or below each bar indicate the confidence level. Significant and non-significant results are further marked using filled and hatched bars, respectively. To facilitate the comparison, the bar in the left box in each panel shows the global RF values as in Fig. 5.

As above, decomposing the RF effect into its different parts helps to disentangle the physical reasons for these aviation-induced effects on the RF. This is shown in Figs. S2–S6 in the Supplement. The warming effect in the shortwave (Fig. S2) is characterized by a very noisy pattern, with a prevalence of strong local warming effects especially in the northern extra-tropics, although the generally low confidence level (<90 %) of the results hampers the identification of coherent patterns. The low values of the clear-sky RF (Fig. S4) confirm the dominance of cloud effects in the shortwave. The longwave effect (Fig. S3) shows a distinct and strong cooling over the continental regions of the northern extra-tropics, especially at high f_act, with a slight dependency on S_crit. This pattern very closely matches the one of the homogeneous freezing fraction (Fig. 3f): in these regions, homogeneous ice formation is therefore expected to dominate the total ICNC, while aviation soot can effectively compete against this process for the available water vapour. Its impact appears to be very effective regardless of the critical supersaturation, as long as the latter remains below the homogeneous freezing threshold, and a sufficient fraction (i.e. 10 %) of aviation soot particles can be active as INPs. The consequence of this is a marked decrease in the cloud frequency (Fig. S6), which then results in the reduced longwave warming. Of course, the shortwave radiation component (Fig. S2) is also affected and indeed shows a warming over the continents, but the signal is very noisy and not as evident as in the longwave. A possible explanation for this could be that other shortwave forcers like low-level clouds are contributing to the model variability, hence enhancing the noise, while in the longwave, the cloud radiative effects are mostly due to cirrus clouds only. At high nucleation efficiencies (S_crit=1.2), a significant cooling is also evident in the clear-sky longwave RF (Fig. S5), with a pretty uniform distribution in the extra-tropics, possibly due to a dehydration effect leading to a reduction in water vapour concentration and a resulting decrease in water-vapour-induced warming. The clear-sky longwave effect rapidly disappears at higher S_crit, and in the northern extra-tropics it turns to a warming, thus contributing to the decrease of the RF with S_crit in this region.

The results of our simulations show therefore that the largest (in absolute terms) and most statistically significant effect is simulated for large efficiencies of aviation soot, i.e. S_crit=1.2 and f_act=10 %, resulting in a cooling effect of −35 mW m⁻². Good nucleation abilities for aircraft soot (S_crit≃1.2) could be considered realistic for soot particles undergoing cloud preprocessing, for example, in contrails. Laboratory measurements (Mahrt et al., 2020) show, however, that only large soot particles with a diameter above ∼400 nm may gain improved ice nucleation abilities resulting from such a preprocessing effect, but particles of this size are rare in the upper troposphere. Considering the measured size distribution by Petzold et al. (1999) as in Sect. 3, one can estimate the fraction of particles with diameter larger than 400 nm to be of the order of 0.001 %, i.e. 2 orders of magnitude below the lowest f_act considered in this study (0.1 %). The parameters measured by Petzold et al. (1999) refer to a young plume, so ageing processes and soot aggregation within contrails might contribute to the increase of the fraction of larger particles in the population, but it is unlikely that a significant fraction of aviation soot will end up in the size range where preprocessing is effective, as also confirmed by aircraft measurements of ice residuals in cirrus and contrail cirrus (Voigt et al., 2017).

Our results generally point towards a relatively small aviation soot–cirrus effect, of the order of tens of milliwatts per square metre (mW m⁻²; in absolute terms), with statistically non-significant figures in several cases. This is in contrast with the estimates by Zhou and Penner (2014), Penner et al. (2018), and Zhu and Penner (2020), who reported larger effects, of the order of hundreds of milliwatts per square metre (in absolute terms), also testing various assumptions on the ice nucleation abilities of aviation soot and experimenting with different parametrizations for ice nucleation. Due to the high complexity of the involved models and the coupling between their different components (aerosol, clouds, radiation, and dynamics), it is challenging to track down the reasons for this disagreement, which could be due not only to different models' schemes and parametrizations, but also to the diverse configurations and tuning approaches. For example, Gettelman and Chen (2013) used a similar model (CAM5) to the aforementioned studies but found no statistically significant effect of aviation soot on natural cirrus clouds, thus being more consistent with the results presented here. McGraw et al. (2020) used a model version from the same family (CESM2, which is based on CAM6 for the atmospheric component), also concluding that the impact of aviation soot is not statistically significant. This was also the case in Hendricks et al. (2011), who used an ECHAM-based global climate model (GCM) as the present one and the same cirrus cloud parametrization, albeit with different aerosol and cloud microphysical schemes. Kärcher et al. (2021) also argued against large RF effects from aviation soot–cirrus interactions, pointing to the limitations of global models in representing key processes as a possible reason for overestimated effects. In conclusion, a consensus among modelling groups on the aviation soot–cirrus effect is still lacking, and future research should consider working towards a concept for a model intercomparison study with common assumptions and detailed analyses of the differences among model configurations and tuning approaches.

We finally recall that the sensitivity experiments conducted in this section are focusing on the ice nucleation abilities of aviation soot, while the properties of background soot (i.e. soot originating from other anthropogenic and biomass burning sources) are not varied and are assumed to be S_crit=1.4 and f_act=0.25 % as in Hendricks et al. (2011) (see also Sect. 3). As a further sensitivity study, we perform two simulations assuming S_crit=1.2 for background soot (simulations BASE-BG12 and S12F10-BG12; see Table 1). This experiment pair results in an aviation soot–cirrus effect of −25.7 mW m⁻² (99.2 % confidence level), which is lower (in absolute term) than the corresponding case calculated above with S_crit=1.4 for background soot (and S_crit=1.2 for aviation soot in both cases). This was to be expected, since increasing the ice nucleation abilities of background soot enhances the competition with aviation soot and the other INPs for available water vapour, thus reducing the impact of aviation soot on natural cirrus clouds. The properties of background INPs could therefore be a further source of uncertainties for the aviation soot–cirrus effect and will be the subject of a companion study.

4.3 Dependency of the aviation soot–cirrus effect on the model representation of the vertical velocity

Besides the ice-nucleating properties of soot, another major source of uncertainties in model studies attempting to quantify the climate impact of aviation on cirrus clouds is the representation of vertical velocities. In the cirrus parametrization adopted here (Kärcher et al., 2006), the vertical velocity is a key parameter, as it controls the critical supersaturation, the competition between homogeneous and heterogeneous freezing (Fig. 3e, f), and the nucleation rate in cirrus clouds. Ice formation in cirrus clouds is strongly influenced by small-scale updrafts, of the order of 1–10 cm s⁻¹ (Barahona et al., 2017), but due to their coarse spatial resolution, global models are not able to represent such phenomena in sufficient detail, and rough approximations are usually introduced to account for them. As explained in Sect. 3, in the EMAC-MADE3 configuration adopted here, the subscale vertical velocity is accounted for by adding an extra-term proportional to the square root of the turbulent kinetic energy to the large-scale, grid box mean, vertical velocity. In the vicinity of mountain ranges, this term is replaced by the contribution of orographic waves to small-scale fluctuations in the vertical velocity, based on the parametrization by Joos et al. (2008), which could lead to stronger updraft of the order of several tens of centimetres per second.

To further explore the limitations behind this approach and their possible impact on the estimated aviation-cirrus effects presented in Sect. 4.2, we perform an additional set of sensitivity experiments, prescribing a geographically uniform constant vertical velocity in the range from 2 to 50 cm s⁻¹. Such an assumption is of course not realistic, but the goal here is not to provide a refined estimate on the aviation soot–cirrus effect but rather to understand the role of dynamic forcing in the results and hence to estimate the uncertainties associated with the RF effects quantified in the previous section. A prescribed constant uniform field allows, for example, regions of the world to be explored which could be important for the cirrus effect but where the model does not simulate a sufficiently strong updraft for ice formation to occur. Due to the large amount of computational resources required, we restrict this sensitivity analysis to a single pair of assumptions for the ice nucleation abilities of aviation soot. To facilitate the analysis, we select the case S_crit=1.2 and f_act=10 %, which is more likely to return statistically significant results, as demonstrated in Fig. 5.

The simulated aviation soot–cirrus effect as a function of different prescribed vertical velocities is depicted in Fig. 8 (global averages) and Fig. 9 (zonal profiles) for RF (and its different components) and other relevant cloud properties, whereas Fig. S7 depicts the geographical distributions of the same quantities. The all-sky RF (Fig. 8a) remains negative for all explored values of the prescribed vertical velocity. The aviation soot–cirrus RF for a prescribed vertical velocity of 2 cm s⁻¹ (simulation S12F10-V2) is −35.9 mW m⁻², very close to the value of −35.4 mW m⁻² obtained with the parametrized vertical velocity (simulation S12F10). This suggests that the aviation effect discussed in this study is mostly controlled by small updrafts of a few centimetres per second. The RF is significantly larger, around −60 mW m⁻², for larger prescribed vertical velocities, except at 10 cm s⁻¹, where a smaller value is simulated. This behaviour results from the combination of the shortwave (Fig. 8b) and longwave (Fig. 8c) all-sky RF, both increasing in absolute terms to a maximum value at 20 cm s⁻¹ and significantly dropping above that value. This trend also characterizes the longwave clear-sky RF (Fig. 8e), which grows towards larger values with increasing vertical velocity and switches from cooling to warming at 20 cm s⁻¹, while the shortwave clear-sky RF (Fig. 8d) is generally small across the explored range of vertical velocities. The zonal profiles (Fig. 9) reveal that the aviation effect is particularly large in the northern extra-tropics, consistent with the distribution of aviation soot particles (Fig. 4) and with the largest values above the North Atlantic, Europe, and the western United States, especially for larger vertical velocities (Fig. S7). The zonal profiles (in particular Fig. 9h) show that the change of sign in the longwave clear-sky RF for v>20 cm s⁻¹ is driven by an increase in total water in the northern extra-tropics. Aviation-induced changes in ICNC (Fig. 9g) and cloud frequency (Fig. 9i) show a similar behaviour, with the strongest effects simulated at v=20 cm s⁻¹. At larger vertical velocities, homogeneous freezing becomes very effective and rapidly consumes the available supersaturated water vapour, so that heterogeneously formed ice crystals have less time to grow to larger sizes, and the sedimentation process becomes less effective. As shown in Figs. 8f and  9f, the aviation-induced reduction of the homogeneous freezing fraction is significantly smaller at v≳20 cm s⁻¹. As a consequence, in this somewhat extreme regime, the impact of aviation soot on the main cloud variables (Figs. 8f–i and 9f–i) loses importance and results in smaller aviation-induced changes in the radiative fluxes (Figs. 8b–e and 9b–e), although their combined effect remains large (Fig. 8a). For the same reason, also the dehydration effect, which we found to induce a pretty strong clear-sky cooling in the longwave (Sect. 4.2), becomes less effective at higher vertical velocities and explains the change of sign at v=20 cm s⁻¹ in Fig. 8e.

Figure 8 Multi-year average (2005–2015) top-of-the-atmosphere (a) all-sky, (b) all-sky shortwave, (c) all-sky longwave, (d) clear-sky shortwave, and (e) clear-sky shortwave RF due to the effect of aviation soot on natural cirrus clouds, for different values of the vertical velocity. Panels (f–i) depict the aviation-soot-induced relative changes in ICNC from homogeneous freezing, ICNC, total water (as the sum of water vapour and ice water), and cloud frequency, respectively, all spatially averaged above 400 hPa and over cloudy and cloud-free model grid boxes. The values above or below each bar indicate the confidence level. Significant and non-significant results are further marked using filled and hatched bars, respectively.

Figure 9 As Fig. 8 but showing zonal averages. Values with a statistical significance below 95 % are drawn in grey.

The sensitivity simulations analysed in this subsection confirm that the model dynamics and the representation of the vertical velocity play an essential role in the microphysics of cirrus clouds. Even relatively small updrafts of a few centimetres per second can induce large changes in the properties of cirrus clouds perturbed by aviation emissions, and the resulting estimates on the climate impact are equally sensitive to the subscale fluctuations in the vertical velocity as to the ice-nucleating properties of aviation soot particles. Our sensitivity experiments further show that the vertical velocity mostly controls the magnitude of the aviation soot–cirrus effect, while the ice-nucleating properties of aviation soot also act on the sign. Both components still represent the largest source of uncertainties for all modelling studies attempting to estimate the climate impact of aviation soot on natural cirrus clouds.

4.4 Impact of nudging

All the results discussed so far refer to model simulations performed in nudged mode, i.e. relaxing meteorological variables (temperature, divergence, vorticity and surface pressure) towards reanalysis data. This approach has been chosen in order to maximize the chances to obtain statistically significant results for the small climate effects and RFs investigated here, and it is a common practice in this kind of study. This technique ensures that the simulations to be compared (in this case with and without aviation effect on cirrus) are characterized by a similar synoptic situation. Due to the chaotic nature of the climate system, this is not the case when running the model in free (climate) mode, as each experiment will develop its own meteorology, which hinders the isolation of the effect of a small perturbation such as the one represented by the impact of aviation soot on cirrus clouds. However, nudging is known to have an impact on simulated temperature profiles (Schultz et al., 2018), in turn affecting the heating rates and all kinds of radiative adjustments in the atmosphere, which implies a potential influence on the effective RF of the climate perturbation under consideration (e.g. Forster et al., 2016; Johnson et al., 2019). In an attempt to characterize the impact of nudging on our results, we repeat the BASE, S12F10, and S14F10 simulations in free-running mode (BASE-FREE, S12F10-FREE, and S14F10-FREE, respectively) and in nudged mode but without relaxing temperature (BASE-FREE-T, S12F10-FREE-T, and S14F10-FREE-T, respectively; see Table 1). In these simulations, short-term feedbacks on temperature (and in the FREE case on other dynamical quantities) can freely evolve. To establish climate conditions comparable to the years modelled in the nudged simulations, SST and SIC are prescribed in both sets: the FREE-T experiments use SST and SIC from the ERA-Interim reanalysis (i.e, the same dataset used for nudging winds and surface pressure), while in the FREE ones climatological means of the 2000–2009 period from the Met Office Hadley Centre dataset (HadISST; Rayner et al., 2003) are used. To increase the chance of obtaining a statistically significant result, the FREE simulations are performed on a longer time period (30 years), while the FREE-T ones cover the same 15-year period as the nudged case.

The resulting all-sky RFs from the aviation soot–cirrus effect are compared in Fig. 10 for the three configurations (nudged, nudged without temperature and free) and the two values of the S_crit parameter. As expected, the free-running simulations are characterized by a much larger statistical noise, which prevents any robust and statistically significant conclusion being drawn on the investigated effect and supports our choice for a nudged configuration. The results of the nudged experiments are, however, consistent with the free-running ones, since they lie within the uncertainty ranges. The two nudging methods (with and without temperature) are highly consistent for S_crit=1.2, which further supports the robustness of the results discussed in this work, although a feedback of temperature nudging on the dynamics cannot be excluded. For S_crit=1.4, the two nudging methods show RFs with opposite sign, although both are statistically non-significant at the 95 % confidence level; therefore no conclusions can be drawn in this case. Nevertheless, this sensitivity study confirms the main conclusions of this work, that large aviation soot–cirrus effects can be simulated only under the assumptions of good ice nucleation ability of aviation soot particles (S_crit=1.2).

Figure 10 Multi-year average top-of-the-atmosphere all-sky RF due to the effect of aviation soot on natural cirrus clouds, under different model configurations: nudged, nudged without relating temperature, and free-running. The latter configuration covers a period of 30 years, while the nudged ones cover 15 years. The horizontal line at the centre of each box is the average value, and each box spans the 95 % confidence interval. Values are shown for S_crit=1.2 and S_crit=1.4, in both cases assuming f_act=10 %. Note that for the free-running experiments, the confidence levels are calculated using Welch's t test instead of a paired-sample t test as in the nudged ones.

5 Conclusions and outlook

In this work, we applied the EMAC-MADE3 global aerosol climate model coupled with a cloud microphysical scheme to quantify the impact of aviation soot on natural cirrus clouds and to estimate the uncertainties in this effect due to the assumptions on the ice-nucleating properties of aviation and background soot particles and on the representation of vertical updrafts in the model. According to our model simulations, the aviation soot–cirrus effect is in the range of −35 to 13 mW m⁻² (with a confidence level >99.9 % and 87 %, respectively) depending on the assumptions on the critical ice saturation ratio and on the fraction of active INPs for aviation soot. Further sensitivity experiments performed with a prescribed, geographically uniform, vertical velocity suggest that the uncertainties in the model representation of vertical updrafts and the resulting cooling rates can add a factor of ≲2 to the estimated effect, although this has been tested only for a specific case corresponding to high nucleation efficiency of aviation soot (S_crit=1.2 and f_act=10 %). For comparison, the estimated effective RFs of the main climate forcers from aviation reported in the recent assessment by Lee et al. (2021) are 34.3 mW m⁻² (CO₂), 57.3 mW m⁻² (contrail cirrus), and 17.5 mW m⁻² (NO_x), with the total aviation effect amounting to 110.9 mW m⁻². Therefore, the aviation soot–cirrus effect quantified here potentially represents an important contribution to the aviation climate impact.

Our results partly support the findings of Liu et al. (2009), who also reported a change of sign in the aviation soot–cirrus effect from negative to positive when increasing the critical saturation ratio from 1.2–1.3 to 1.4, although our RF values are considerably smaller, possibly due to the much lower active fraction assumed in this study: 0.1 %–10 % compared to 100 %. Systematically larger negative RFs are also found by the studies of Penner et al. (2009), Zhou and Penner (2014), Penner et al. (2018), and Zhu and Penner (2020), all assuming generally good ice nucleation abilities for aviation soot. Our study therefore agrees with these in terms of the sign of aviation soot effect, but there is a clear disagreement in the magnitude. This could be due to the details of the cloud microphysical parametrizations, to differences in the model setups (e.g. in the use of nudging), or to the representation of the vertical velocity, which can have a large impact on the resulting effects, as shown by our idealized sensitivity experiments. The model setup adopted here is very close to the one used in Hendricks et al. (2011), who used a similar base model (ECHAM4) and the same cirrus parametrization (Kärcher et al., 2006). They nevertheless found a non-significant impact of aviation soot on natural cirrus clouds, which can possibly be ascribed to the use of a free-running model setup in that study, while the nudging mode used here likely helped to extract a statistically significant information from the model. Non-significant results were also found by Gettelman and Chen (2013) and, more recently, by McGraw et al. (2020). We note that for RFs of the considered magnitude, a statistically significant signal quantification is only possible with simulations using the nudging technique, although it is not completely clear how strong this procedure may impact the resulting effective RF values.

This work helps to disentangle and quantify the main uncertainties in the aviation soot–cirrus effect, but its actual magnitude (and to a lesser extent also its sign) remain uncertain. Our model simulations show that a more precise characterization of the ice-nucleating properties of aviation soot could help to constrain at least the sign of the resulting RF. Further laboratory measurements are therefore needed, in particular concerning the role of cloud processing. In situ measurements are also essential, in order to characterize the microphysical properties of the population of aviation soot, like number and size. Recent advancements in the representation of soot-induced ice formation in cirrus clouds (Kärcher et al., 2021; Marcolli et al., 2021) should also be considered for implementation in global model parametrizations. Further on the modelling side, increasing the models intercomparability, for example, by performing commonly designed experiments in the context of an intercomparison exercise, would provide valuable insights on the results and could help interpret the model diversity in better detail. In addition to the CMIP activities (Eyring et al., 2016; Collins et al., 2017), notable examples of the advantages of such exercises are provided by the AeroCom community (e.g. Schulz et al., 2006; Quaas et al., 2009; Mann et al., 2014; Samset et al., 2014, see also https://aerocom.met.no/), last access: 26 November 2021.

The EMAC-MADE3 model adopted here was shown to provide a satisfactory representation of aerosol, aerosol-induced ice formation in cirrus, and overall cloud and radiation properties, but biases still exist. The model tends to overestimate the ice crystal number concentration in cirrus clouds at higher cirrus temperatures (≳225 K). Hence, the contribution of cirrus clouds formed in situ by heterogeneous and/or homogeneous ice nucleation and cirrus clouds resulting from transport of ice from the mixed phase (e.g. via detrainment from convective clouds) should be further explored. Moreover, the cirrus parametrization by Kärcher et al. (2006) implemented in EMAC-MADE3 has some known limitations which should be addressed in future versions of the model, for example, replacing critical supersaturation with an activation spectrum following the change in the active INP population with increasing supersaturation and considering the effects of particle size and coating on the soot nucleation process. The representation of background aerosol like soot and mineral dust was already demonstrated to be quite good by Kaiser et al. (2019), Righi et al. (2020), and Beer et al. (2020). Nevertheless, improvements are still possible, especially with the support of in situ measurements. The inclusion of additional background INP types, such as ammonium sulfate and glassy organics, in the model should also be considered, as this might change the outcome of the competition for available supersaturated water vapour among the different INPs, with possible impacts on the aviation soot effect estimated in this work.