Multi-model evaluation of the sensitivity of the global energy budget and hydrological cycle to resolution

This study undertakes a multi-model comparison with the aim to describe and quantify systematic changes of the global energy and water budgets when the horizontal resolution of atmospheric models is increased and to identify common factors of these changes among models. To do so, we analyse an ensemble of twelve atmosphere-only and six coupled GCMs, with different model formulations and with resolutions spanning those of state-of-the-art coupled GCMs, i.e. from resolutions coarser than 100 km to resolutions finer than 25 km. The main changes in the global energy budget with resolution are a systematic increase in outgoing longwave radiation and decrease in outgoing shortwave radiation due to changes in cloud properties, and a systematic increase in surface latent heat flux; when resolution is increased from 100 to 25 km, the magnitude of the change of those fluxes can be as large as 5 W m−2. Moreover, all but one atmosphere-only model simulate a decrease of the poleward energy transport at higher resolution, mainly explained by a reduction of the equator-to-pole tropospheric temperature gradient. Regarding hydrological processes, our results are the following: (1) there is an increase of global precipitation with increasing resolution in all models (up to 40 × 103 km3 year−1) but the partitioning between land and ocean varies among models; (2) the fraction of total precipitation that falls on land is on average 10% larger at higher resolution in grid point models, but it is smaller at higher resolution in spectral models; (3) grid points models simulate an increase of the fraction of land precipitation due to moisture convergence twice as large as in spectral models; (4) grid point models, which have a better resolved orography, show an increase of orographic precipitation of up to 13 × 103 km3 year−1 which explains most of the change in land precipitation; (5) at the regional scale, precipitation pattern and amplitude are improved with increased resolution due to a better simulated seasonal mean circulation. We discuss our results against several observational estimates of the Earth's energy budget and hydrological cycle and show that they support recent high estimates of global precipitation.


Introduction
The state of the climate results from the complex balance of fluxes of energy in and out of the Earth system, as well as fluxes of energy and water between its components. The Electronic supplementary material The online version of this article (https ://doi.org/10.1007/s0038 2-018-4547-y) contains supplementary material, which is available to authorized users. 1 3 steady progress of satellite imagery in the recent decades and the development of networks of in situ measurements have undoubtedly led to more precise estimates of those budgets (Trenberth et al. 2009(Trenberth et al. , 2011Stephens et al. 2012;L'Ecuyer et al. 2015;Rodell et al. 2015). However, large uncertainties remain: surface radiative fluxes, which are difficult to measure from space, are only known within ± 10 W m −2 ; the Global Precipitation Climatology Project (GPCP), which is the best available dataset for global precipitation, appears to be underestimating precipitation by at least 10% in the light of recent CloudSat observations (Stephens et al. 2012). Climate models, despite their inherent biases, have the potential to ensure consistency between the intertwined energy budget and hydrological cycle. As such, they are used to infer more accurate estimates of energy and water budgets when observations are missing. Wild et al. (2013Wild et al. ( , 2015, for instance, used CMIP5 models to remedy the lack of coverage of the Baseline Surface Radiation Network (BSRN) stations and derive estimates of the radiative balance at the surface. In addition, all the other estimates use, to some extent, model derived data when measurement cannot be obtained from satellite data as is the case for moisture fluxes (Trenberth et al. 2011;Rodell et al. 2015).
As a new generation of advanced high-resolution climate models emerges with horizontal resolution of 25 km or finer (Haarsma et al. 2016;Kinter et al. 2013;Mizielinski et al. 2014;Roberts et al. 2018a, b), the question thus arises as to how the modelled energy and water budgets could benefit from increasing general circulation models (GCMs) resolution. As proposed by several studies, such as Pope and Stratton (2002) and Haarsma et al. (2016), studies looking at the sensitivity to resolution should try to answer the following questions: (1) are there processes emerging at higher resolution that are not accounted for by model parameterization?
(2) does higher resolution alleviate known errors associated with parameterizations? (3) do small-scale processes feedback on the large scale (4) making high resolution desirable for simulating the large-scale climate with improved fidelity? (5) is there a resolution threshold at which climate processes emerging with resolution converge?
Those remain open questions. Few studies (e.g. Pope and Stratton 2002;Bacmeister et al. 2014;Demory et al. 2014) have assessed the resolution sensitivity of the energy budget and hydrological cycle in global climate models. Those which did often focused on changes at the regional scale or on a single model.
In evaluating the resolution effect on the energy budget, only few studies applied no (or little) resolution-specific retuning, and these noted that different resolutions exhibit very similar energy budget properties on global annual time scales (Hack et al. 2006;Demory et al. 2014). In Demory et al. (2014) the main change in the energy budget related to the shortwave radiation (SW) reflected by clouds, which tends to decrease with higher resolution, leading to an erroneous increase of surface net SW radiation of 2 W m −2 . This error is compensated by a decrease of downward longwave radiation (LW) at surface (a common compensation in climate models, e.g. Wild et al. 2015).
The increase of atmospheric models resolution has not been found to produce drastic changes to the simulation of the time mean global hydrological cycle (Boyle andKlein 2010, Bacmeister et al. 2014). Most changes occur at the regional scale and are due to improved regional and mesoscale circulation, attributed to better resolved orography and land-sea contrasts. For example, recent publications have demonstrated better represented orographic effects in the Himalayas during the Indian Monsoon (Shaffrey et al. 2009;Lau and Ploshay 2009;Bacmeister et al. 2014;Bush et al. 2015), better representation of orographic jets and their associated moisture transport (de Souza Custodio et al. 2017) or diurnal forced circulation such as sea-breeze effect (Boyle and Klein 2010) and a general improvement of the simulated maritime continent hydroclimate because of better resolved coastline and orography Johnson et al. 2016). High resolution also improves the representation of synoptic systems: it produces more tropical cyclones of strong intensity (Strachan et al. 2013;Shaevitz et al. 2014;Bacmeister et al. 2014;Roberts et al. 2015) and better represented extratropical cyclones (Jung et al. 2006) which might impact the hydrological cycle through their moisture transport Terai et al. 2018). Moreover, increasing model resolution favours the organisation of convection. Vellinga et al. (2016) showed that convection was more organised within African easterly waves when model resolution was increased, giving rise to a positive feedback between moisture convergence and precipitation, involving large-scale circulation and improving regional rainfall processes over the Sahel. Several studies have noticed a role of resolution on modulating the diurnal cycle of precipitation (Bacmeister et al. 2014;Birch et al. 2015). In fact, high resolution enables the generation of more explicit precipitation, which peaks at different time from parameterized convective precipitation.
At the global scale, Demory et al. (2014) and Terai et al. (2018) found that while ocean precipitation decreases with higher resolution, land precipitation increases due to higher moisture convergence over land. The contribution of moisture transport to land precipitation also increases, whereas moisture recycling, a quantity that is known to be overestimated by state-of-the-art GCMs, tends to decrease. However, the mechanisms leading to those changes are still unclear. Pope and Stratton (2002) could not find asymptotic convergence of hydrological processes over four different resolutions ranging from N48 (2.5° × 3.75°) to N144 (0.833° × 1.25°) in HadAM3 and attributed this result to non-linearities both in the hydrological cycle and the dynamics. Demory et al. (2014) results further suggested that the global water cycle is sensitive to horizontal resolution, up to about 60 km, from where their results asymptote to a local plateau. This finding is based on two different versions of the same model and therefore needs to be tested in other GCMs for the potential to generalise.
Despite the beneficial impact of resolution on several aspects of the hydrological cycle at the global and regional scales, the global pattern of precipitation biases remains largely the same when resolution is increased Bacmeister et al. 2014;Terai et al. 2018). Boyle and Klein (2010) noted that some inherent model errors associated with physical parameterizations (e.g., tropical diabatic drying rates in their model) cannot be ameliorated with resolution only.
One question of interest is how the sensitivity to resolution changes with model formulation [e.g. dynamical core, physical parameterizations and coupling strategy (coupled/ uncoupled)]. Dynamical cores have different abilities to conserve physical properties such as mass, energy and moisture: finite volume schemes are conservative by design, finite difference schemes can conserve moisture although they are not designed specifically to do so, but spectral models do not conserve moisture. Not only do dynamical cores have a crucial impact on the closure of global budgets but they also affect the representation of processes. Landu et al. (2014), for instance, investigated the role of resolution in the double ITCZ bias and showed that the sensitivity to resolution is strongly dependent on the dynamical core used in their model.
The assessment of the energy budget and hydrological cycle changes with resolution is often made difficult by the fact that simulations with the same model at different resolutions use different model configurations (e.g. model tuning, radiative and SST forcings). In particular, the different resolutions of the same model might be retuned independently to maintain the top of atmosphere (TOA) energy balance and hydrological cycle as close to observations as possible (Pope and Stratton 2002;Roeckner et al. 2006), which introduces uncertainties as to what is actually causing the changes: tuning or resolution itself. Parameters that are typically retuned include time-step and horizontal diffusion (to keep the model stable), critical relative humidity thresholds for clouds and cloud water auto-conversion threshold (to balance energy at TOA), drag coefficient (to compensate the increase of surface latent heat flux, Duffy et al. 2003), gravity wave drag (to account for the resolution of the orography, Pope and Stratton 2002), adjustment time scale related to moist convection. Duffy et al. (2003) investigated the effect of retuning their high-resolution model and showed that, despite a reduction of cloudiness and precipitation biases, it had only moderate effects on precipitation spatial patterns. When comparing their high and low-resolution models, Pope and Stratton (2002) and Terai et al. (2018) explicitly distinguished the role of the horizontal grid resolution itself from tuning parameters, by running a simulation on the low-resolution grid using the high-resolution parameters. Terai et al. (2018) found that grid resolution contributed to two-thirds of the global precipitation increase in their model while tuning parameters were responsible for one-third. A multi-model assessment seems necessary to further identify systematic changes of the energy and water budget with resolution, to evaluate those changes quantitatively and to determine the impact of model formulation. The present study takes the opportunity of the European project PRI-MAVERA to undertake such a model intercomparison. PRI-MAVERA is a European funded project whose aim is to develop, under a coordinated protocol, a new generation of advanced and well-evaluated high-resolution global climate models, capable of simulating and predicting regional climate with unprecedented fidelity. Data and methods are presented in Sect. 2, the model intercomparison is presented in Sect. 3 and results and possible implications for our understanding of global energy budget and hydrological cycle estimates and projections are discussed in Sect. 4.

Model ensemble
In this study, we have collected 12 atmosphere-only models (AMIP type) and 6 coupled GCMs, each of them having the capability to span a large range of resolutions between 200 and 15 km. The vertical resolution is kept constant in all the configurations of a same model. The model tuning is performed on the lowest resolution version of each model, except for MRI3.2 for which the highest resolution has been tuned. No resolution-specific retuning is applied unless otherwise stated. The simulations have been performed as part of different projects: UPSCALE (HadGEM3-GA3), Febbraio (HadGEM3-GC2), EU-Horizon 2020 PRIMAVERA (CNRM-CM6-1, EC-EARTH3.1, ECMWF-IFS, HadGEM3-GC31, CMCC-CM2, MPIESM1.2), CASCADE (CAM5.1), CMIP5 (GFDL-HIRAM). The models are described below in further detail and Table 1 summarises the main information.

CMCC-CM2
CMCC-CM2 atmospheric component is CAM4.5. The model has 30 vertical levels and a model top at ~ 2 hPa. It has been run at two atmospheric resolutions (0.9° latitude × 1.25° longitude; 0.23° latitude × 0.31° longitude, respectively 1° and 0.25° hereafter). In coupled mode, the two atmospheric resolutions are coupled to NEMO3.6 at a resolution of 0.25°. The model uses MACv2-SP (produced for the EASY aerosols framework, Stevens et al. 2018), which provides a prescription of anthropogenic aerosol optical properties and an associated Twomey effect.

CNRM-CM6-1
CNRM-CM6-1 consists in the atmospheric model ARPEGE-Climat version 6.3. The ARPEGE-Climat dynamical core is derived from IFS cycle 37t1. Compared to CNRM-CM5, the atmospheric physics has been largely revisited. In particular, convection scheme, microphysics scheme and turbulent scheme have been updated. The convection scheme (Guérémy 2011;Piriou et al. 2007) provides a consistent, continuous, and prognostic treatment of convection from dry thermals to deep precipitating events. The microphysics scheme is derived from Lopez (2002) and takes into account autoconversion, sedimentation, ice-melting, precipitation evaporation and collection. The turbulence scheme represents the TKE with a 1.5-order scheme prognostic equation according to Cuxart et al. (2000). The model is operated with a T127 and a T359 truncation. In both versions there are 91 vertical levels in the atmosphere (with top level at 78.4 km). The model is forced by HadISST2 daily SST at a resolution of 0.25° (Rayner et al. 2006). Ozone concentration is a model's prognostic variable (Michou et al. 2011) and anthropogenic aerosols are interpolated from CNRM TACTIC v3 (Michou et al. 2015).

HadGEM3-GA3
This model is the atmospheric component of HadGEM3 (Hewitt et al. 2011) in the GA3.0 configuration with 85 vertical levels extending to 85 km in height (Walters et al. 2011). HadGEM3-GA3 describes a non-hydrostatic atmosphere using a semi-Lagrangian, semi-implicit formulation, which allows an increase in horizontal resolution, while keeping a relatively long time step necessary for climate integrations (Davies et al. 2005). It uses a regular latitude/ longitude grid and is developed at three horizontal resolutions (N96, N216, and N512). Simulations cover the period 1986-2011. Daily SST and sea ice forcings are derived from the Operational Sea Surface Temperature and Sea IceAnalysis (OSTIA) product (Donlon et al. 2012), which has a native resolution of 1/20°. The ensemble has 5 members for N96, 3 for N216 and 5 for N512. See Mizielinski et al. (2014) for more detail on these simulations.
A 100-year present-day control simulation with forcings from year 2000 (equivalent to experiment 2 from CMIP3) has been completed at three atmospheric resolutions (N96, N216 and N512) coupled to a 0.25° ocean and where relevant (e.g., for aerosols) emissions vary through the annual cycle.

ECMWF-IFS
ECMWF-IFS is a climate model configuration based on cycle 43r1 of the European Centre for Medium-Range Weather Forecasts Integrated Forecast System. This configuration follows forcing recommendations from CMIP6/High-ResMIP and is described in detail by Roberts et al. (2018a, b). The atmospheric dynamical core of the IFS is hydrostatic, semi-Lagrangian, and semi-implicit (SISL) with computations alternated between spectral and reduced Gaussian grid-point representations at each time step. ECMWF-IFS has 91 vertical levels on a hybrid sigma-pressure coordinate. The grid-point representation of ECMWF-IFS uses a cubic octahedral reduced Gaussian grid such that the shortest wavelengths in spectral space are represented by 4 model grid points to give resolutions of 25 km (Tco399) and 50 km (Tco199) in ECMWF-IFS-HR and ECMWF-IFS-LR configurations, respectively. There was one resolution-specific tuning of an atmospheric parameterization. In prototype atmosphere-only experiments, it was found that net surface energy balance in ECMWF-IFS-HR was about 1 W m −2 too low compared to observations. To reduce this bias, the autoconversion threshold for liquid precipitation over the ocean was modified in ECMWF-IFS-HR to give the same global surface energy balance as ECMWF-IFS-LR. All configurations include representation of the land-surface (Balsamo et al. 2009) and ocean-wave field (Janssen 2004). The ocean and sea-ice components of the coupled configurations are based on NEMO v3.4 (Madec 2008) and LIM2 (Bouillon et al. 2009;Fichefet and Maqueda 1997), respectively. In ECMWF-IFS-HR, NEMO and LIM2 use the eddy-permitting ORCA025 configuration which as a nominal resolution of ∼ 0.25°. In ECMWF-IFS-LR, NEMO and LIM2 use non-eddying ORCA1 configuration, which has a nominal horizontal resolution of ~ 1°. The ECMWF-IFS-LR version of NEMO is configured to be as close as possible to the ECMWF-IFS-HR, with differences limited to resolutiondependent parameterizations such as the Gent and Mcwilliams (1990) scheme for the effect of mesoscale eddies. Anthropogenic aerosols are prescribed from MACv2-SP (Stevens et al. 2018). It is a successor of version 2.3, which has been used for CMIP5. EC-Earth3.0.1 includes updates of its atmospheric and oceanic model components as well as higher vertical and horizontal resolutions in the atmosphere. It has been used by e.g. Hartung et al. (2017) and Koenigk and Brodeau (2017). EC-EARTH3.0.1 is based on ECM-WF's atmospheric circulation model IFS, cycle 36r4. EC-EARTH3.0.1 is provided by the Swedish Meteorological and Hydrological Institute (SMHI). The model has been run at six different spectral truncations (T159, T255, T319, T511, T799, T1279).

EC-EARTH3.1
• EC-EARTH3.1a is provided by the Consiglio Nazionale delle Ricerche (CNR) and described in Davini et al. (2017). Several changes have been brought to the model since the version 3.0.1 described above, in particular to solve a P-E imbalance issue. Changes made to the model include: the introduction of a proportional advection mass fixer for water species in the atmosphere, a retuning of oceanic diffusive albedo (increased from 0.06 to 0.05), a change of runoff flux correction and a reduc- The name and version of each model is given in column 1. Column 2 indicates whether or not the models are coupled. Column 3 indicates the model formulation (SL stands for semi-Lagrangian, SI for semi-implicit). Column 4 lists the topology of the grid used by each model. In column 5, we indicate the number of vertical levels. In column 6, we list the different model resolutions used in this article. We indicate for each model the configurations chosen at the low resolution (in italic) and high resolution (underlined) to build the low and high resolution multi-model ensemble means, respectively. In column 7, we provide the number of years used for each simulation and in column 8, the number of ensemble members used in this study (more members have been produced in the various projects but only a subset was made available to our study)  (Titchner and Rayner 2014). • EC-EARTH3.1b is provided by various members of the EC-Earth Consortium as part of the PRIMAVERA project. The atmospheric only version of the model is forced by HadISST2 daily SST at 0.25°. Anthropogenic aerosols are prescribed from MACv2-SP (Stevens et al. 2018). • In EC-EARTH3-CPL, IFS is coupled to ocean model NEMOv3.6 and sea ice model LIM3. The low resolution coupled model combines the spectral truncation T255 and the ocean grid ORCA1 at a resolution of ~ 1° and the high resolution coupled model combines the spectral truncation T511 and the ocean grid ORCA025 at a resolution of ~ 0.25°.

GFDL-HIRAM
The GFDL global High Resolution Atmospheric Model is based on AM2 (GAMDT 2004) with increased horizontal and vertical resolutions as well as simplified parameterisations for moist convection and large-scale (stratiform) cloudiness. A finite-volume core using a cubed-sphere grid topology with a quasi-uniform horizontal grid spacing with 32 vertical levels. The model is forced by HadISST 1.1 monthly SSTs (Rayner et al. 2003 Demory et al. (2014) recommends to assess the global energy and water budgets simulated by climate models against two types of budgets, either derived from observational datasets or simulated by reanalyses. Studies based on (globally consistent) observational data are constructed to satisfy closure of the energy and water budgets in principle. However, some issues remain in terms of the availability (e.g. data might be absent or scarce in some regions) and quality of observational products (e.g. error in instrument calibration, bias in the retrieval procedure). On the other hand, reanalyses provide a comprehensive and self-consistent view of the energy and water budgets, but they are not specifically designed to conserve energy and water. The energy budget simulated by models is compared to four observation-based estimates derived independently. Trenberth et al. (2011, hereafter T11)  , derive an energy budget using a variational approach that optimizes the match between several satellite products collected by NASA's Energy and Water cycle Study (NEWS) and a set of equations describing energy and water budgets closure. Their approach aims to improve the closure of the Earth's energy budget and hydrological cycle, as well as the consistency between those two.

Observations and reanalysis products
The hydrological cycle is compared against two estimates. Trenberth et al. (2011)  Those observational estimates are complemented with the values simulated by the two reanalyses which have the best ability to balance the global energy and water budgets according to Trenberth et al. (2011): the European Centre for Medium-Range Weather Forecasts (ECMWF)'s reanalysis, ERA-Interim (Dee et al. 2011) and the NASA Goddard Center's reanalysis, MERRA (Rienecker et al. 2011).
Radiative fluxes are further assessed using CERES-EBAF, covering the period 2000-2015. CERES-EBAF provides net balanced Top of Atmosphere (TOA) fluxes derived from CERES satellite, where the global net is constrained to the ocean heat storage term ) and calculates mean surface fluxes using a radiative transfer model .
Precipitation is further assessed against the global precipitation analysis GPCP (Adler et al. 2003) merging data from rain gauge stations, satellites, and sounding observations and the multi-satellite product TRMM 3B42 (Huffman et al. 2007).

Estimation of orographic precipitation
Orographic precipitation is known to be particularly sensitive to resolution and we define a quantitative metric to assess it across the simulations analysed in this paper. An obvious metric would be to choose an altitude threshold and define precipitation as orographic for all regions above this threshold. There are at least two disadvantages to this method: first, the area defined as orographic would differ a lot among the different resolutions of a same model as orography is generally lower at lower resolution; second, regions of high altitude but relatively dry, such as the Tibetan plateau or Antarctica, would be included in the composite while producing little rain. Another possible metric is to define a mask based on an estimation of the orographic enhancement in a reanalysis, which has the advantage of addressing the two issues mentioned above. In fact, applying the same mask to all simulations has the advantage of identifying precipitation occurring over the same area, thus avoiding "sampling issues". Moreover, a mask identifying areas where orographic enhancement occurs will not retain dry region at high altitudes. Finally, a reanalysis offers realistic circulation patterns which is essential in the simple model described hereafter, and we choose to use ERA-Interim.
The estimation of orographic enhancement is based on a model developed by Sinclair (1994). In this model, the vertical orographic uplift at the surface is determined as a function of the dot product of horizontal surface wind and orographic slope. The vertical wind profile is a function of the orographic uplift and set to decrease vertically up to 200 hPa following a power law of atmospheric pressure. The condensation rate, which is function of the vertical wind, specific humidity and relative humidity at each level, is integrated vertically from the Lifting Condensation Level to the top of the atmosphere, to give an estimate of orographic precipitation [please refer to Sinclair (1994) for additional detail about their model]. This model can identify regions of potential orographic precipitation only, since several processes have not been accounted for, such as the background vertical wind, hydrometeor advection by the background flow, drying downwind from mountains, re-evaporation by the entrainment of dry air subsidence. We apply this model to daily ERA-Interim data. Our calculation gives a maximum orographic enhancement of 2 mm day −1 which occurs over mountainous terrain (see Fig. S1a). This is less than the total precipitation simulated by ERA-Interim in these areas and this can be explained because of all the assumptions listed above.
A monthly climatology of the orographic precipitation estimate is then computed and a mask is created with 1 when the orographic precipitation estimate is above 0.5 mm day −1 and 0 otherwise. Regions of orographic precipitation are extended by a radius of 100 km to account for the transport of hydrometeor away from their region of formation. In the Andes, orographic precipitation is detected all year round (see Fig. S1b supplementary) whereas over the western part of the Himalayas orographic precipitation is detected only during summer as would be expected due to the strong seasonality of the monsoonal flow.

Mean state sensitivity to resolution
In this section, as a prerequisite to the evaluation of the energy budget and hydrological cycle, we assess the zonally averaged climatologies of temperature, specific humidity and zonal wind in each model of the ensemble (Fig. 1).
The most striking change associated with the increase of horizontal resolution is a stratospheric cooling, which occurs at all latitudes in most models, with the exception of HadGEM3-GC31 and CAM5.1, as well as MRI3.2 in the tropics. Pope and Stratton (2002) showed that a similar cooling occurred in dynamical core experiments which they attributed to the non-linear interactions of stratospheric gravity waves at higher resolution. The stratospheric cooling at higher resolution corresponds to an increase of the multi-model mean temperature bias with regard to ERA-Interim (Fig. 2).
There is a tropospheric warming in the midlatitudes when resolution is increased in all models except EC-EARTH3.1a and b and ECMWF-IFS, but with different vertical structures depending on models: it is largest in the upper troposphere in HadGEM3_GA3, HadGEM3-GC31 and MRI3.2 in contrast to CAM5.1, GFDL-HIRAM, CNRM-CM6-1, which exhibit a vertically uniform warming. In their model, Pope and Stratton (2002) attributed the midlatitudes warming to several different causes. In the mid-troposphere, increasing resolution resulted in enhanced transient vertical winds, Fig. 1 Difference in the zonal mean annual mean temperature (low resolution represented in black contours every 10 K; difference in colour), specific humidity (low resolution represented in black contours every 2 × 10 −3 kg kg −1 ; difference in colour) and zonal wind (low resolution represented in black contours every 5 m s −1 ; difference in colour) between the high and low resolutions of each model of the ensemble leading to more condensation and latent warming, whereas in the upper troposphere the warming was found to be due to a better representation of the dynamics. The tropospheric warming at higher resolution corresponds to a decrease of the temperature bias in the troposphere (Fig. 2). Consistently with the tropospheric warming, there is a decrease of relative humidity in the troposphere in all models of the ensemble (not shown).
There is no agreement between models on the sensitivity of specific humidity to resolution, especially in the tropics. However, most models tend to simulate a drying of the low-level troposphere in both hemispheres around 50°. This is consistent with the drying described by Pope and Stratton (2002), which they attributed to the moisture depletion of low-levels by enhanced transient vertical winds. It is unclear if this drying corresponds to an improvement or a worsening when compared with ERA-Interim (Fig. 2).
The main changes in the zonal wind consist of a reduction of the subtropical jet in most models at 200 hPa and around 30° in both hemispheres, except in ECMWF-IFS. Several models show a poleward shift of the southern eddydriven jet with resolution (HadGEM3-GA3, HadGEM3-GC31, MRI3.2 and GFDL-HIRAM). Similar changes have been described in Pope and Stratton (2002), Roeckner et al. (2006), Demory et al. (2014) and Lu et al. (2015). However, other models do simulate an equatorward shift (CNRM-CM6-1), an intensification (CAM5.1) or a weakening (EC-EARTH3.1a and EC-EARTH3.1b) of the eddydriven jet. In the northern hemisphere there is no clear and systematic response of the eddy-driven jet. However, there is a marked reduction of the multi-model mean bias of surface easterlies in the Tropics when resolution is increased (Fig. 2).
Finally, we assess the change in Hadley cell circulation in Fig. 3. Following Byrne and Schneider (2016), we recalculate vertical velocity in pressure coordinates from the meridional gradient of the meridional mass stream function. The differences in vertical velocity are mainly in the Tropics. They suggest a narrowing of the ascending region of the Inter-Tropical Convergence Zone (ITCZ) at higher resolution and a less intense subsidence in subsiding regions (i.e. the opposite response to that predicted in a warmer climate, Lu et al. 2007).

Global energy budget
In Fig. 4, we analyse the different components of the energy budget following Demory et al. (2014) against the observational estimates by T11, S12, W13, R15. We have also added the observations used in R15 to Fig. 4 (open gray diamond) as this gives an indication of how their optimization procedure adjusts the budget components to satisfy closure.
TOA net radiation in the model ensemble is contained within a range of ± 2 W m −2 , which is larger than the uncertainty of present-day TOA imbalance deduced from changes of ocean heat content 0.6 ± 0.4 W m −2 (Loeb et al. 2009;Llovel et al. 2014;Dieng et al. 2015) and used in the various observational estimates (Fig. 4a). Variations of TOA net radiation with resolution are within 2 W m − 2 for each individual model, except CAM5.1, which has a decrease of 3.5 W m −2 between 2° and 0.25° (Fig. 4a). The variations of TOA net radiation with resolution reflect overall those of the net energy absorbed at the surface (Fig. 4l). Note that there is a 2 W m −2 loss of energy at TOA in EC-EARTH3.0.1, which is not compensated by an upward net energy flux at the surface: hence the larger net energy loss of the atmosphere in this model (Fig. 4n represents the net energy input in the atmosphere, NEI, see Sect. 3.2.2). However, EC-EARTH3.0.1 NEI bias is only moderately sensitive to resolution. In contrast, in EC-EARTH3.1a, there is an increase of TOA net radiation with resolution and a decrease of the net surface energy balance so that the atmosphere is receiving more energy with increasing resolution.
All models show a decrease of outgoing SW radiation with resolution (Fig. 4b). Models tend to have values close to observations and reanalyses, i.e. 100 W m −2 , for the Fig. 2 Multi-model mean of the RMSE change between the high and low resolutions for AMIP models (top) and coupled models (bottom) and for temperature, specific humidity and zonal wind (from left to right). ERA-Interim is used as the reference. Brown indicates a decrease of the bias at higher resolution and green an increase of the bias resolution at which the tuning was performed (low-resolution for all models except high resolution for MRI3.2). There is a difference of up to 6 W m −2 in MRI3.2 between low and high resolution outgoing SW radiation. The decrease of outgoing SW radiation with resolution is almost entirely explained by SW reflected by clouds and atmosphere (Fig. 4d). The various estimates evaluate SW reflected by clouds and atmosphere to be in the range 76-80 W m −2 , but most high resolution models simulate values lower than 76 W m −2 , suggesting a lack of cloud reflectivity. Moreover, the decrease of SW reflected by clouds results in an increase of net SW radiation at the surface with higher resolution (Fig. 4f). Shortwave reflected by clouds remains in most models in the range of uncertainty of S12 but is larger than the upper estimate of W13. Note that the bias in net SW at the surface in HadGEM3-GA3 and HadGEM3-GC2, irrespective of resolution dependence, is consistent with a lack of SW absorbed by the atmosphere (Fig. 4e), which possibly reflects a bias introduced by the radiative transfer code.
At TOA, the decrease in outgoing SW radiation is largely compensated by an increase of outgoing longwave radiation (OLR) (Fig. 4c). OLR in high resolution models Zonal section of vertical velocity in pressure coordinates deduced from the meridional mass-stream function (see text for further detail). a Low (resp. high) resolution fields are in green (resp. red), individual models (resp. model ensemble mean) are represented by thin (resp. thick) lines and ERA-interim is represented by the black line, b difference between high and low resolutions for atmosphere-only models, c difference between high and low resolutions for coupled models. In b and c, the thick black line stands for the AMIP and coupled models ensemble mean, respectively tend to exceed the value of the different estimates even though the range of uncertainty of S12 encompasses most of the models' values. Longwave radiation (LW) emitted at the surface is not very sensitive to resolution in atmosphere-only models, due to prescribed SST, so that the OLR increase can be mostly attributed to changes in the atmosphere.
The incoming LW radiation lies below the range of uncertainty of the various estimates in half of the models and is underestimated by up to 10 W m −2 , which is a bias commonly found in climate models as previously reported by Stephens et al. 2012 andWild et al. 2013 (Fig. 4g). However, it has very different resolution sensitivity among models. At the surface, the increase of net SW is compensated by an increase of latent heat flux (LH) in all models (Fig. 4k). In atmosphere-only runs, the mechanism of this compensation is not obvious: as will be shown later, the largest LH increase is over the ocean where SST is fixed and insensitive to incoming SW changes. All high resolution models simulate a LH compatible with S12 but exceed the upperbound of the other estimates. The relatively high LH value of 88 W m −2 proposed by S12 is plausible, as it is obtained by correcting a systematic underestimation of tropical precipitation in GPCP. The optimization procedure used in LE15 brings LH from 75 (open diamond in Fig. 4k) to 81 W m −2 (plain diamond), underlining the call of their optimization procedure for higher LH values.
The surface sensible heat flux (SH) presents almost no change with resolution. SH is close to T11 in most models but lower than the other estimates. A slight decrease in SH in MRI3.2, GFDL-HIRAM, HadGEM3-GC2 is consistent with the tropospheric warming simulated in those models, which reduces the temperature contrast between the sea surface and low-level atmosphere.
In all models, there is little sensitivity to resolution of the SW radiation absorbed by the atmosphere (Fig. 4e), which has a large uncertainty in S12 and W13.
Overall, coupled runs (HadGEM3-GC2, HadGEM3-GC31-CPL, ECMWF-IFS-CPL and EC-EARTH3-CPL) behave similarly to their atmosphere-only counterparts. In particular, they do simulate an increase of OLR and LH and a decrease of outgoing shortwave radiation when resolution increases. Those changes have the same amplitude as in atmosphere-only runs.
Most of the changes described so far are consistent either with higher tropospheric temperature (increase of OLR, LH) or lower cloud reflectivity (decrease of SW reflected by clouds, increase in OLR and increase in net surface SW). To further distinguish these effects, we decompose radiative fluxes into cloud forcing and clear-sky conditions (Fig. 5). SW cloud forcing explains the largest change of SW reflected by clouds and atmosphere (Fig. 5a) in all models except EC-Earth3.1b. However, both LW cloud forcing and clear-sky emission contribute to the change in OLR (Fig. 5b,  d). The OLR increase is mostly due to an increase in clearsky emission in HadGEM3-GC3.1 but to a reduction in LW cloud forcing in CAM5.1.
Additional analysis showed that the compensation between SW and LW cloud forcing largely occurs within tropical regions characterized by mean ascending motion (not shown) and are suggestive of a change in how convection is organized (more compact convective systems at higher resolution) and/or a change in convective cloud properties.
In the next section, we further investigate how the changes described so far depend on latitude.

Zonal sections of the energy budget
We now consider the net energy input in the atmosphere (NEI) based on net TOA and surface energy balances. If we ignore kinetic energy and surface friction, NEI is equal to the tendency of moist static energy: where S is the moist static energy, SW, shortwave radiative flux, LW, longwave radiative flux, SH, surface sensible heat flux and LH, surface latent heat flux, SFC stands for surface. Arrows indicate whether the radiation is upward or downward positive definite.
Almost all models, with the exception of EC-EARTH3.1a, consistently show a decrease of the energy gained in the tropics at higher resolution and a decrease of the energy lost in the extra-tropics (Fig. 6f). At steady state, we expect (1)  . 6 Zonal sections of (a, f, k) energy absorbed by the atmosphere, b, g, l SW absorbed by the atmosphere, c, h, m surface latent heat flux, d, i, n surface sensible heat flux, e, j, o LW cooling. a-e Low (resp. high) resolution fields are in green (resp. red), individual models (resp. model ensemble mean) are represented by thin (resp. thick) lines, ERA-interim is represented by the black line and CERES-EBAF by the dashed black line; f-j difference between high and low resolutions for atmosphere-only models; k-o difference between high and low resolutions for coupled models. In f-j and k-o, the thick black line stands for the AMIP and coupled models ensemble mean. All fluxes are expressed in W m −2 the atmosphere to lose as much energy in the extra-tropics as it gains in the tropics. This is obviously not the case for EC-EARTH3.1a, whose atmosphere gains energy globally as explained in Sect. 3.2.1. When NEI is split into the individual components of the energy budget (SW absorbed by the atmosphere, surface latent and sensible heat fluxes and LW cooling), the agreement between models is not as good as for NEI (Fig. 6g-j) suggesting that several processes might be playing a role or that another common mechanism drives this response.
LH increases in all models in the Southern Ocean (Fig. 6h). In atmosphere-only models, this can be explained either by a decrease in specific humidity at low-level or by an increase of the module of surface wind. Terai et al. (2018) differentiated the bulk formulae of LH and found that the change with resolution was mainly driven by an increase of surface wind speed. Note that we would expect an increase of surface wind speed to enhance SH too, which is not the case (Fig. 6i), possibly because the temperature increase in the troposphere has an opposite effect on the bulk formula. It is also possible that the increase of LH is the result of a more complex change in non-linear air-sea interactions within extratropical cyclones.
The LW cooling (i.e. the net energy loss by the atmosphere due to LW) is increasing in the Tropics with increasing resolution in most models (Fig. 6j) but there is a significant spread of the resolution sensitivity among models. Despite the tropospheric warming simulated by most atmospheric models in the extratropics, there is only a moderate increase of LW cooling at those latitudes. The inter model spread (i.e. around the thick black line) in LH seems to balance in the inter model spread in LW cooling, which can be expected from radiative convective equilibrium.
As already shown, there is only little variation in the SW absorbed by the atmosphere with resolution (see Fig. 3e) and the difference between the model ensemble mean high and low resolutions is almost negligible at all latitudes. Thus, we conclude that SW absorbed by the atmosphere is not a likely contributor to the systematic latitudinal changes of NEI with resolution.
We find noticeable changes between the coupled and the atmosphere-only models. There is an increase of the energy absorbed in the tropics at higher resolution, in particular in the region of the ITCZ (Fig. 6k). This is due to a large increase of LH in the subtropics and even more so in the ITCZ with resolution (Fig. 6m), which is only partly compensated by LW cooling (Fig. 6o).
In the previous section, we demonstrated that cloudradiative forcing plays an important role in the sensitivity of the radiative fluxes with resolution, therefore we investigate in Fig. 7 the change in cloud fraction. The sensitivity of the cloud cover with resolution varies among models. For instance, HadGEM3-GA3 and HadGEM3-GC31 simulate an increase of the tropical cloud fraction of 5% between high and low resolution but CAM5.1 simulates a decrease overall. Those changes could be consistent with changes in LW cloud forcing (Fig. 5b) and must relate to high convective clouds which exert a larger LW cloud forcing. However, the change in cloud fraction does not relate in an easy way with SW cloud forcing (Fig. 5a) which is for instance decreasing in HadGEM-GC31. The change in cloud fraction only cannot fully explain the change in cloud forcing and more work is needed using high frequency data to better determine which of the different cloud properties (optical thickness, height, fraction) have the largest effect on radiation.
The decrease of NEI in the tropics must affect the amount of energy that is transported to the extratropics. This is investigated in the next section.

Meridional energy transport
At steady state, NEI has been integrated meridionally to find the meridional flux of moist static energy (VS): where R = 6380 km is the radius of the Earth and all the other quantities have been defined above.
As expected from the reduction of energy absorbed in the tropics and lost in the midlatitudes, there is a reduction of the poleward moist static energy transport in most models when resolution is increased (Fig. 8). The amplitude of this reduction is about 0.2 PW (Fig. 8b) and is maximum around 30ºN/S (which corresponds approximately to the poleward edge of the Hadley cell). In the tropics the decrease of the poleward energy transport is consistent with the weakening of the Hadley cells (Fig. 3b).
The moist static energy transport in coupled models does share some similarity with AMIP models but is shifted toward more negative values in the tropics. Thus, there is more southward transport of energy at the Equator at higher resolution (Fig. 8c) and this is consistent with the gain of energy in the ITCZ region (Fig. 6k).
In the midlatitudes, the decrease of the meridional energy transport is consistent with the reduction of the meridional temperature gradient. This is suggested by two diagnostics in the additional material: 1. the contrast of longwave cooling between the midlatitudes and the Equator is the main contributor to the contrast of NEI (Fig. S2), 2. the eddy sensible heat flux in the component of the horizontal flux that explains the best the reduction of the total energy transport and its decrease is the largest at the tropopause in the few models in which we could calculate it (Fig. S3). 1

Hydrological cycle sensitivity to resolution
In this section, we assess the sensitivity of the global water cycle to atmospheric models' resolution. In particular, we test whether the increase of moisture transport to land at higher resolution demonstrated by Demory et al. (2014) is a systematic behaviour of GCMs with a special focus on the role played by model formulation (dynamical core and coupling).

Global hydrological cycle
We first analyse if models have the ability to conserve moisture. At steady state, a model which conserves water in the atmosphere should have E-P equal to 0 globally. Zonally averaged cloud fraction (see text for further detail). a Low (resp. high) resolution fields are in green (resp. red), individual models (resp. model ensemble mean) are represented by thin (resp. thick) lines and ERA-Interim is represented by the black line, b dif-ference between high and low resolutions for atmosphere only models, c difference between high and low resolutions for coupled models. In b and c, the thick black line stands for the AMIP and coupled models ensemble mean, respectively Any deviation from this corresponds to a virtual and spurious source or sink of moisture. Global E-P is presented in Fig. 9h. CAM5.1, GFDL-HIRAM, CNRM-CM6-1 and all the HadGEM3 models maintain E-P less than 0.5 × 10 3 km 3 year −1 , at all resolutions. By contrast, MRI3.2, EC-EARTH3, EC-EARTH3.1a, ECMWF-IFS are those for which E-P deviates most significantly from 0. In EC-EARTH3.0.1, global P-E can be as large as 15 × 10 3 km 3 year −1 , which is roughly equivalent to a third of the moisture advection to land (Fig. 9g). Several factors might explain the different ability of models to conserve moisture, in particular the type of dynamical core and the use of a semi-Lagrangian scheme possibly associated to a mass fixer. Regarding the dynamical core, finite volume models, such as CAM5.1 and GFDL-HIRAM, are designed to conserve moisture. Conversely, spectral models, such as MRI3.2, EC-EARTH3.0.1, EC-EARTH3.1a and b, ECMWF-IFS, have poor performance at conserving moisture. In EC-EARTH3.0.1 and 3.1a, physical parameters are tuned on the lowest resolution, whereas in MRI3.2, they are tuned on the highest resolution. Therefore, high resolution in the different EC-EARTH3 simulations and low resolution in MRI3.2 deviates the most from balance. Despite being . a Low (resp. high) resolution fields are in green (resp. red), individual models (resp. model ensemble mean) are represented by thin (resp. thick) lines, b difference between high and low resolutions for atmosphere-only models, c difference between high and low resolutions for coupled models. In b and c, the thick black line stands for the AMIP and coupled model ensemble mean respectively. All fluxes are expressed in PW a spectral model, CNRM-CM6-1 simulates a water budget near closure and performs as well as finite volume models.
The use of a semi-Lagrangian scheme requires an interpolation to estimate the value of the field being advected at departure point: this is a source of errors in the water budget. This problem can be partly solved by the introduction of a mass fixer algorithm, applied on the tracer field, which allows to find a mass conserving solution. In ECMWF-IFS, no mass fixer was used. However, sensitivity experiments have shown that implementing the mass fixer for semi-Lagrangian advection dramatically reduces the magnitude of the non-conservation (not shown), but gives a residual of the opposite sign indicating that there are remaining processes that are also contributing to the non-conservation. A similar reduction and change in sign of E-P is simulated between EC-EARTH3.0.1 and EC-EARTH3.1a, after the introduction of a mass fixer in EC-EARTH3.1a (Davini et al. 2014). The introduction of a mass fixer in a previous version of CNRM model also contributed significantly to the reduction of the non-conservation of moisture (Voldoire et al. 2013).
Most models simulate a global precipitation increase at higher resolution (Fig. 9a), which despite the inability of some models to conserve moisture, is consistent with the overall increase of evaporation (Fig. 9d), as already mentioned in Sect. 3.1 (Fig. 4k). The increase of precipitation with resolution occurs either over land (CAM5.1 and HadGEM3 family) or over ocean (EC-EARTH family, MRI3.2). There is a fairly good agreement between estimates for land precipitation to be in the range 113 to 120 × 10 3 km 3 year −1 , but models which show an increase in land precipitation simulate values at high resolution ranging . 9 Resolution dependence of hydrological processes: a global precipitation, b land precipitation, c ocean precipitation, d global evaporation, e land evaporation, f ocean evaporation, g land E-P, h global E-P, i precipitable water, j strength of hydrological cycle. Units are a-h 10 3 km 3 year −1 , i mm and j cycles year −1 . The x-axis represents the equivalent resolution at 50°N (in km). In addition to Rodell et al. (2015) optimized estimates (plain gray diamond), we have added Rodell et al. (2015) OBS corresponding to the raw observations (open gray diamond) in order to show the effect of the optimization from 130 to 145 × 10 3 km 3 year −1 . Contrary to precipitation, the increase of evaporation with resolution occurs mostly over the ocean in all models (Fig. 9e, f). The inter-model spread in land evaporation is as large as the spread of observational estimates going from 70 to 85 × 10 3 km 3 year −1 . The spread of moisture convergence to land in the low resolution models is similar to that of the various observational estimates, i.e. 30-45 × 10 3 km 3 year −1 , however high resolution models span a much larger range with values up to 60 × 10 3 km 3 year −1 (Fig. 9g).
Models do not agree on the change of precipitable water with increasing resolution (Fig. 9i): it increases in all HadGEM3 models and ECMWF-IFS but decreases in CAM5.1. In this model, the strength of the hydrological cycle (defined by Trenberth 1998, as the ratio between global precipitation and precipitable water) increases with resolution, possibly leading to a drying of the atmosphere (Fig. 9j).

Regional precipitation changes due to resolution
In this section, we evaluate how the sensitivity of precipitation to resolution varies at different latitudes and in various regions.
From Fig. 10a, it is clear that the difference between the model ensemble mean precipitation and the observations is of the same order of magnitude as the difference between the two different observational datasets, TRMM and GPCP. Figure 10b shows the annual mean and zonally averaged precipitation difference between high and low resolution for each model. All models with the exception of GFDL-HIRAM show a decrease of precipitation at the Equator at higher resolution and an increase in the subtropics. CAM5.1 is the only model to simulate a narrowing of the northern branch of the ITCZ. This model also has a marked increase of precipitation in the Southern Ocean. For most models, the off-equatorial increase of precipitation corresponds to a deterioration of the double ITCZ bias, as already noted by Bacmeister et al. (2014).
The multimodel zonal mean precipitation around the Equator shows a decrease in the summer hemisphere and an increase in the winter hemisphere of about 0.5 mm day −1 (in Fig. 10e, h). Most models agree with this change but there is a large intermodel spread. In coupled models, there is an increase of precipitation in the ITCZ region and a decrease in the south hemisphere subtropics with increasing resolution both in DJF and JJA.

Fig. 10
Zonal sections of precipitation. a, d, g Low resolution models in green and high resolution models in red, GPCP in plain black and TRMM in dashed black. b, e, h difference between high and low resolution atmosphere-only models and multi-model mean difference in black. c, f, i Difference between high and low resolution coupled models and multi-model mean difference in black. a-c Annual mean, d-f JJA and g-i DJF. Precipitation is given in mm day −1 comparable resolution changes. Interestingly, models agree on the sign of the sensitivity to resolution in several regions. In particular, there is an increase of precipitation in regions of high orography, storm tracks, South American river basins (Amazon and La Plata) and Sahel (Fig. 11a). In DJF, there is a shift of precipitation over the Southern Ocean toward the pole in most models (Fig. 11b). There is a decrease of precipitation in the western tropical Atlantic and an increase in the eastern tropical Atlantic (in the year mean and in JJA particular, Fig. 11c). Increasing resolution has also a drying effect over the Maritime Continent (mostly ocean, Fig. 11a).
In Fig. 12, we evaluate whether the precipitation changes presented in Fig. 11a, d correspond to an improvement or a deterioration of the bias. RMSE has been calculated against two different observational datasets, GPCP and TRMM. Of all the regions with a change of precipitation with resolution, only a few are associated with a clear reduction of the precipitation bias. One of those is the La Plata basin, which is supplied in moisture from the Amazonian basin by a lowlevel orographic jet that is better represented at higher resolution (as also shown in de Souza Custodio et al. 2017). The drying of the western Atlantic basin corresponds also to an improvement, consistently with Siongco et al. (2015) who suggested that higher (lower) resolution models tend to have an eastern (western) precipitation bias. There is an overestimation of precipitation in the ITCZ in the eastern Pacific with increasing resolution.

d-f
The six coupled models included in the ensemble are HadGEM3-GA3, EC-EARTH3, EC-EARTH3.1a, GFDL-HIRAM, CAM5.1 and GFDL-HIRAM. Units are mm day −1 . Note the color bar has a log scale. Large and small dots indicates that all models and all but one model, respectively, agree on the sign of the change Several precipitation changes with resolution are the same as in AMIP. However, in coupled models, precipitation increases north of the Equator and decreases south of the Equator when resolution increases, which corresponds to a stronger asymmetry of precipitation than in AMIP models. The six models agree on the reduction of the precipitation bias throughout the tropical Atlantic and in the southeastern Pacific. There is a systematic reduction of precipitation biases in South Africa, over the Agulhas current and over Sahel in coupled models that is not observed in AMIP models. More work is needed to determine if the strong reduction of precipitation biases at higher resolution in coupled models comes from a reduction of SST biases or emerging coupled processes.

Precipitation intensity
In Fig. 13 we have represented the exceedance frequency of a given precipitation rate in the tropical band (± 30ºN). To do so, 6 hourly data have been regridded on a 1º × 1º grid and classified by bins of 0.75 mm (6 h) −1 . The figure shows a systematic divergence of precipitation rate between the high and low resolution models above 20 mm (6 h) −1 , which is consistent with Bacmeister et al. (2014), Vellinga et al. (2016) and Terai et al. (2018). Hence, models produce more intense high intensity precipitation events at higher resolution. Two models, CAM5.1 and HadGEM3-GC3.1, are able to simulate the same frequency of high intensity rain rate as TRMM. Models tend to overestimate light rain [< 1.5 mm (6 h) −1 ] in comparison to TRMM as already  Exceedance frequency of precipitation rate as a function of precipitation rate in atmospheric only models at low resolution (plain curve) and high resolution (dashed curve) and in TRMM (gray) between 30ºS and 30ºN. The inset shows the frequency of no-rain conditions. Frequencies are given for bins of width 0.75 mm (6 h) −1 described by Terai et al. (2018). We further evaluated the frequency of no rain conditions in the inset of Fig. 13. TRMM has more no rain conditions than any model but TRMM is reported to underdetect low to moderate rain rates in the Tropics (Huffman et al. 2007). We do not find a consistent change of no rain conditions with resolution among this set of models.

Land precipitation and moisture advection to land
A direct comparison of land precipitation among models is of limited relevance, as the intensity of the hydrological cycle varies among models (Fig. 9j). Instead, Demory et al. (2014) evaluated the ratio of land precipitation to global precipitation and found that it increases in two versions of their model. We reproduce their diagnostic in Fig. 14b. We find that all observational estimates have relatively close values of this ratio, at around 0.22. Interestingly, two behaviours seem to emerge among models. The ratio increases by approximately 10% from low to high resolution in HadGEM3 models, CAM5.1 and CMCC-CM2. The high resolution version of these models also get closer to the various observational estimates. On the other hand, the ratio decreases in MRI3.2 and EC-EARTH3.0.1 and 3.1. It is interesting to note that models in which the ratio increases are grid-point models, whereas those in which it decreases are spectral models. Moreover, coupling does not seem to change these findings (compare HadGEM3-GA3 and GC2 and EC-EARTH3.0.1 and EC-EARTH3.0.1 CPL).
We now test another finding of Demory et al. (2014), who showed that the fraction of land precipitation caused by moisture convergence increased with resolution. At steady state, P-E averaged over land should be equal to the moisture advected from ocean to land. Thus, the ratio (P-E) land /P land corresponds to the fraction of land precipitation explained by moisture convergence. The values of this ratio are presented in Fig. 14a. The ratio is lower in reanalysis MERRA (0.25) and ERA-Interim (0.30), indicative of a strong moisture recycling over land, and higher in T11 (0.35) and R15 (0.39). In grid-point models CAM5.1, HadGEM3-GA3, GC2 and GC31, there is an increase of about 5% or more of the ratio between the low and high resolutions. Concurrently, the fraction of precipitation explained by moisture convergence increases by only 2% or less in EC-EARTH3.0.1 and 3.1a and there is no net change between low and high resolutions in MRI3.2. High resolution models have values closer to the range of T11 and R15.
Overall, we find that models behaving as described in Demory et al. (2014), i.e. showing both an increase of P land / P total and (P-E) land /P land with resolution, are finite difference and finite volume models. On the other hand, the spectral models of the ensemble show a decrease of P land /P total and a lesser increase of (P-E) land /P land .
In the remainder of the paper, we test two hypotheses to explain why the fraction of total precipitation that falls on land and the moisture convergence to land increase more in grid-point models: first, it is due to a better representation of atmospheric eddies and their associated transport of moisture; second, it is triggered by an increase of orographic precipitation.
We tested the first hypothesis by recalculating mean and eddy moisture transport to land in HadGEM3-GA3 and EC-Earth3.0.1 and found that the component that changes the most with resolution is the mean moisture advection in the tropics (not shown). In the extratropics, moisture advection (a) (b) Fig. 14 Resolution dependence of a (P-E) land /P land (i.e., fraction of land precipitation explained by moisture convergence) and b P land /P total to land by eddies is more sensitive to resolution but explains a fraction of the global change that is one order of magnitude less than in the Tropics. We further tested whether the mean moisture advection change was due to a change in mean circulation or mean distribution of moisture by replacing subsequently low resolution wind and moisture in the calculation of high resolution fluxes (not shown) and found that it is caused by mean circulation. Hence, we conclude that it is the mean moisture advection to land in the tropics, due to a change in mean circulation, that explains the largest variation of moisture convergence to land at higher resolution.

Role of orographic precipitation
We now test the hypothesis that the orography triggers more precipitation at higher resolution and that it is balanced by an increase in moisture convergence. The impact of orography on the hydrological cycle was evaluated in two different ways: first, we compare a control simulation HadGEM3-GA6 at high resolution (N480) with a sensitivity experiment in which the resolution of the orography was decreased while that of the atmospheric model was held constant (Schiemann et al. 2018; we refer to this experiment as N480_LRO hereafter), which has the potential to evaluate changes due to the resolution of the orography only; second, for each model and resolution of the ensemble, we partition orographic and non-orographic precipitation. The partitioning is achieved using a pre-defined mask for regions of high orographic precipitation rate in ERA-Interim (see Sect. 2 for detailed information on how this mask was build), which was then applied to each individual model. The careful comparison of Figs. 11a and 15 suggests that a significant portion of regional precipitation changes with increasing resolution can be attributed to the resolution of orography itself. Among those, we can cite the increase of precipitation in regions of high orography and the drying of the Maritime Continent, which presents complex coastlines and orography. On the other hand, some changes such as the intensification of precipitation in the eastern Atlantic at high resolution are not reproduced when the resolution of the orography only is modified. This suggests that the better resolved orography over Africa and its possible positive effect on tropical cyclogenesis does not play a major role in this change. It is more likely that the increase of precipitation over the Sahel at higher resolution has to do with the better representation of convective processes (Vellinga et al. 2016). A thorough description of these biases and understanding their origin will be the object of another study.
In Fig. 16 we present the partitioning of precipitation occurring in regions of high orographic precipitation and elsewhere. Before all, note that the difference of precipitation between N480 and N480_LRO is entirely captured by orographic precipitation. This suggests that precipitation differences in non-orographic regions between those two runs (which were described in the previous paragraph) tend to cancel out and/or that the precipitation difference occurring in orographic region is at least one order of magnitude larger than in non-orographic regions (note that the color bars of Figs. 11, 15). Global orographic precipitation ranges between 43 and 55 × 10 3 km 3 year −1 in all climate models and between 43 and 49 × 10 3 km 3 year −1 in GPCP and the two reanalysis products. Non-orographic precipitation ranges between 63 and 83 × 10 3 km 3 year −1 in models and between 70 and 77 × 10 3 km 3 year −1 in the observational estimates. Thus, regions identified as producing orographic precipitation in the mask represent around 20% of land surfaces, but correspond to approximately 35-45% of total precipitation in models and in the various observational estimates. More interestingly, despite representing less than half of total precipitation, orographic precipitation captures most of the change of precipitation due to resolution: the change of orographic precipitation between low and high resolution is larger than 10 × 10 3 km 3 year −1 in HadGEM3-GA3 and CAM5.1, whereas the change is respectively 3 and 5 × 10 3 km 3 year −1 for non-orographic precipitation in those two models. Finally, the orographic precipitation change with resolution is two times larger in the grid-point models of the ensemble than in the spectral models. Finally, we should mention that the difference of global precipitation between N480 and N480_LRO (respectively 576.7 and 574.7 × 10 3 km 3 year −1 ) is smaller than the difference of land precipitation (respectively 130.5 and 122.9 × 10 3 km 3 year −1 ). Hence, the resolution of the orography affects mostly the split of precipitation between ocean and atmosphere in the HadGEM3 models, but generates only a small increase of global precipitation.

Discussion
The present study investigates the sensitivity of the energy budget and hydrological cycle to atmospheric models' horizontal resolution in an ensemble of 18 CGCMs which covers a range of resolutions going from 200 km, as typically used by state-of-the-art climate models, to 15 km. Our primary goal is to evaluate how processes emerging with resolution modify the simulated climate rather than to build better models. Thus no or little specific retuning has been performed at each resolution that would have ensured that models remain as close as possible to the observations. In other words, we do not necessarily expect all aspects of the simulated climate to be more accurate or in better agreement with the observations at higher resolution. Moreover, we take the opportunity that the models of the ensemble use different formulations and different degrees of coupling to test the robustness of changes due to increasing resolution.
Our main findings can be summarized as follows: • The radiative balance at TOA is modified, with an OLR increase and an outgoing SW decrease. At surface, the increase of net SW is compensated by an increase of LH, in particular in the midlatitudes. Overall, there is less energy absorbed by the atmosphere in the tropics at higher resolution. As a result, the poleward atmospheric energy transport decreases. The reduction of the poleward energy transport is consistent with a weakening of the Hadley circulation in the tropics and the reduction of the Equator to Pole temperature contrast. • The increase of OLR and the decrease of outgoing SW in the tropics are primarily due to a change of cloud radiative forcings in regions of mean ascending motion. A possible explanation is that, at higher resolution, high intensity precipitation events are generated by more compact and more intense convective systems, thus reducing the mean cloud fraction. A more detail analysis of cloud radiative properties per dynamical regime such as in Bony and Dufresne (2005) and using high frequency data should help understand which cloud processes are responsible for those changes. This will be the object of a future study. • There is a large increase of surface evaporation with resolution in all models but the fate of moisture was found to depend on model formulation. In grid point models, there is a large increase of moisture convergence to land and of the fraction of total precipitation that falls on land as found by Demory et al. (2014). In spectral models however, the moisture advection to land increases less with resolution and the fraction of precipitation that falls on land decreases. • We find that it is the increase of orographic precipitation that explains the large increase of land precipitation in grid-point models. Despite representing only about 40% of total precipitation, orographic precipitation was shown to capture most of the change of land precipita-(a) (b) Fig. 16 Partitioning of precipitation into a orographic and b non-orographic precipitation, using the mask presented in Sect. 2 tion due to resolution. Spectral models, which apply a smoothing to orography to avoid the spurious generation of gravity waves, show less sensitivity to resolution. Nonorographic precipitation exhibits little change with resolution but is characterized by a large intermodel spread. • We find overall little asymptotic behaviour of the fraction of land precipitation explained by moisture convergence at high resolution as described by Demory et al. (2014). Thus, moisture convergence to land might continue increasing if resolution is increased even more. • Overall, coupling plays only a small role in the sensitivity of the global moisture budget to resolution. However, at the regional scale, several systematic improvements were found in coupled models when resolution increases which do not occur in their AMIP counterpart (in particular the large increase of precipitation in the northern branch of the ITCZ, the reduction of the double ITCZ bias in the southeastern tropical Pacific and the reduction of precipitation biases in the tropical Atlantic). This suggests that identical coupled processes emerge with increasing resolution in the various models (such as tropical instability waves, Shaffrey et al. 2009).
The global increases of precipitation and LH with resolution must be understood together as the two should be exactly balanced. Demory et al. (2014) proposed that the increase of LH was caused by increasing SW radiation at surface while Terai et al. (2018) suggested that it was the result of stronger surface wind speed. Moreover, Terai et al. (2018) attributed part of the change in global precipitation to the tuning parameters (time steps). Our results show that the increase of surface evaporation occurs mostly over the ocean and is as strong in AMIP as in coupled simulations. As SST is fixed in AMIP simulations, it cannot mediate the increase of incoming shortwave to the surface latent heat flux. More dedicated simulations would be needed to conclude on the mechanism causing the increase of precipitation/LH with resolution across models.

Possible implications for the observational estimate of the energy budget and hydrological cycle
There is a large increase of global precipitation from 20 to 40 × 10 3 km 3 year −1 across the range of resolutions investigated, corresponding to an increase of 5 to 8% for each model. At the highest resolution, models support the hypothesis proposed by Stephens et al. (2012) that low-resolution global observations of precipitation (i.e. GPCP) may not be accurate enough and may lack precipitation over high heterogeneous mountains. Their correction of GPCP precipitation results in a global latent heat value larger than other observational estimates, and closer to high-resolution GCM values (see also Fig. 1 in Roberts et al. 2018a, b). The optimization of the variational approach of L'Ecuyer et al. (2015) and Rodell et al. (2015) ensures a closure of the water and energy budgets, but tends to adjust global precipitation by smaller increments than Stephens et al. (2012). Interestingly, L'Ecuyer et al. (2015) note that "while initial energy imbalances suggest that latent heating should be increased significantly, water cycle constraints limit the magnitude of the adjustments since precipitation already exceeds the sum of evaporation and runoff". Following Stephens et al. (2012) hypothesis, our results confirm that land precipitation in L'Ecuyer et al. (2015) (based on GPCP) and moisture advection to land (partly based on MERRA) might be underestimated by about 15 × 10 3 km 3 year −1 . We believe that higher values of land precipitation and moisture advection to land without modifying land evaporation would lead to a larger increase of global latent heat during the optimization. In addition, this increase would occur at the expense of sensible heat which, after the optimization, has a larger value than in any other estimate. We have shown here that the GCMs used in this study are within the large range of observational uncertainty, and we argue that making a better use of such model configurations, together with up-to-date observations, would improve our understanding of the processes, and thus the way we should treat both GCMs and observations to reduce the uncertainty range.

Caveats on model configuration
We acknowledge that, in order to assess the sensitivity to resolution of physical processes simulated in climate models, a relevant definition of what is low and high resolution should be physically motivated. In this study, we used grid properties to define resolution, i.e. the actual grid points spacing at 50ºN in grid-point models, or half the wavelength of the shortest zonal wave in spectral models. It has been suggested that a more relevant definition of the models' effective resolution is to take the wavelength at which the eddy kinetic energy departs from the expected − 5/3 slope in the inertial range of a log-log spectrum (Skamarock 2004;Evans et al. 2013). A new PRIMAVERA study (Klamer et al. 2018), comparing the same GCMs presented here, is working in this direction. Although the resolution specific retuning is kept as little as possible, some parameters need to be modified to satisfy models' stability, in particular the model's timestep or the interval at which the parameterization of convection is called. The effect of changing those parameters has not been addressed in this study but was shown to lead to a substantial increase of global precipitation according to Terai et al. (2018).
Despite these caveats, some robust changes of the energy and water budgets are shared by models across the range of resolutions investigated. This is remarkable as different tuning strategies were adopted among models by adjusting different sets of parameters and performing the retuning either on the lowest or highest resolution. For instance, despite displaying opposite sensitivities of TOA net radiation to resolution, the sensitivity to hydrological processes is similar in HadGEM3-GA3 and CAM5.1 (strong increase of the fraction of precipitation that falls on land and the fraction of land precipitation explained by moisture convergence). Moreover, following the PRIMAVERA protocol, HadGEM3-GC31, CNRM-CM6-1 and EC-EARTH3.1b use the same aerosols forcing (EASY aerosols) and SST forcing (HadISST2 daily). On the other hand, EC-EARTH3.1a and EC-EARTH3.1b are close versions of the same model and mostly differ by their forcing. Despite this, EC-EARTH3.1a and b energy and water budgets behave the same way. From these two considerations, we conclude that forcing does not affect the sensitivity of the energy and water budget to resolution.

Some recommendations to further identify the drivers of the sensitivity to resolution
We can suggest the following recommendations for future intercomparison of high-resolution GCMs, such as that performed within the framework of HighResMIP (High Resolution Model Intercomparison Project; Haarsma et al. 2016): • It is not always easy to distinguish the role of tuning parameters from that of horizontal grid resolution. To have a companion low-resolution simulation using the tuning parameters of the high-resolution simulation, as done in Pope and Stratton (2002) and Terai et al. (2018), would help distinguish their respective roles. • We demonstrated in the present study the crucial role of the resolution of orography on the hydrological cycle in most GCMs. Performing simulations at high resolution with the low-resolution orography has the potential to attribute changes due to the resolution of orography only. • Comparison with simulations using idealized configurations (e.g., aquaplanet, idealized land or orography) might be useful for hypothesis testing.