The dependency of the d 18 O discrepancy between ice cores and model simulations on the spatial model resolution

This study examines outputs of a regional climate model (RCM) enabled to compute fractionation of water isotopes over the Greenland ice sheet. The COSMO_iso RCM is forced at the lateral boundaries with isotope enabled GCM simulations with atmospheric nudging. Outputs of COSMO_iso simulations for the present day and the mid-holocene (at a 50 km spatial resolution)are compared against ice core isotopic measurements. For the present-day simulations the RCM simulations generally improved the agreement with observations compared to the GCM results, with the improvements generally occurring in regions with coarser GCM resolution. Higher-resolution RCM simulations at 7 km did not further improve the agreement, producing a worse agreement in some instances. For the mid-Holocene simulations, there was not a large improvement resulting from the RCM simulations (although data were available only from four ice cores). The authors note that the higher-resolution simulations provide a range of spatial variability for the coarse resolution grid that can be used to generate a distribution for comparison against ice core measurements. They also examine gradients of isotope ratio relative to temperature, finding higher variability in temperature and isotope ratios along the ice sheet margins.

here in the usual δ notation, δD and δ 18 O with respect to the Vienna Standard Mean Ocean Water V-SMOW) reflect the 30 atmospheric conditions under which the fractionating process took place (Dansgaard, 1964;Merlivat and Jouzel, 1979;Gat, 1996). This process is generally utilized to reconstruct paleo-climate conditions like past changes of temperature, out of isotopic ratios stored in climate archives (Dansgaard et al., 1969;Masson-Delmotte et al., 2005;Jouzel, 2013).
In Arctic regions like Greenland, ice cores constitute an exceptional climate archive. Over thousands of years, accumulated snow was solidified to ice, preserving at some locations the water isotopic ratios since the last interglacial period. Climate 35 reconstructions based on these ice cores show that the climate conditions changed considerably in Greenland during the Holocene (here defined as the period between present-day and 12 ka; Marcott et al., 2013). Temperature was steadily rising during the early Holocene, until the Holocene Thermal Maximum was reached in the mid-Holocene (6 ka). Since then, temperatures were steadily decreasing until the late Holocene (Marcott et al., 2013;Moossen et al., 2015). In this context, the mid-Holocene is a period of particular interest, as by that time an Arctic warming took place due to orbital forcing variations 40 and their related feedbacks on large-scale climate variations, which exhibits similarities to the strong recent Arctic warming.
For Greenland, the mid-Holocene provides the opportunity to investigate the processes, leading to this warming, in more detail and to potentially obtain new insights about the future development of the Arctic region (Yoshimori and Suzuki, 2019).
While General Circulation Models (GCMs) are generally able to reflect the direction and large-scale patterns of past climate changes (e.g. Timm and Timmermann, 2007;Smith and Gregory, 2012), they often fail to reproduce the magnitude of regional 45 changes (Braconnot et al., 2012;Harrison et al., 2014), which are documented in various local climate archives. Thus, a scale gap might exist between the measured point information and the large-scale climate information generated by GCMs. The comparison of observational and GCM data can therefore be subject to considerable uncertainties (Felzer and Thompson, 2001).
Especially for structured landscapes, the spatial resolution in GCMs is often too coarse to resolve relevant local factors (Jost 50 et al., 2005;Fischer and Jungclaus, 2011). Important properties like topography and surface conditions are consequently only represented in a generalized and imprecise form in climate simulations. In most cases, this does not meet the complex characteristics of the land surface and its associated interactions with the atmosphere. For stable water isotopes, key physical processes of isotope fractionation are therefore also often not entirely resolved in coarsely resolved GCMs, leading to differences between simulated and observational isotope data, especially in complex terrains (Sturm et al., 2005;Werner et 55 al., 2011). Consequently, in many cases GCM results exhibit larger deviations to observational data than the results of corresponding Regional Climate Model (RCM) simulations (Sjolte et al., 2011;Russo and Cubasch, 2016;Ludwig et al., 2018).
In the presented study, isotope-enabled GCM simulation results for the Arctic region are dynamically downscaled with an isotope-enabled RCM to a higher temporal and spatial resolution. By means of such regional simulations, the spatial variability 60 of the isotopic ratios in the Arctic is potentially increased, accounting for the heterogeneity of local conditions at the different ice core locations. In this way, the impact of highly resolved local conditions on the spatial variability of isotopic ratios is https://doi.org/10.5194/cp-2019-156 Preprint. Discussion started: 20 January 2020 c Author(s) 2020. CC BY 4.0 License. investigated, and it is analysed, if such small-scale spatial variability can potentially explain the discrepancy between simulated and observed paleo-climate conditions in the Arctic region.
To explore this, the isotope-enabled version of the RCM COSMO-CLM (Rockel et al., 2008), COSMO_iso (Pfahl et al., 2012;65 Christner et al., 2018), is used. The results of the standard variables from COSMO in the Arctic region under present-day conditions are evaluated in a separated study (Karremann et al., under preparation). In the scope of ERA-Interim reanalysis (Dee et al., 2011) driven simulations, the model showed that it is generally able to simulate reasonable near-surface temperatures and precipitation rates for Greenland and can therefore be used for isotope applications in this region. The capability of COSMO_iso to simulate realistic water isotopic ratios for Greenland is tested by downscaling a global present-70 day simulation with an isotope-enabled GCM for the Arctic region. The GCM and RCM results are subsequently compared to measured water isotope ratios from snow pit samples. Afterwards, the tested isotope-enabled COSMO_iso model system is used to downscale an isotope-enabled GCM simulation for a mid-Holocene time-slice. The simulated isotopic ratios are evaluated against Greenland ice core data. Such a dynamical downscaling of global isotope simulations for Greenland under mid-Holocene conditions, is performed for the first time in the framework of this study. 75

Model Description
In this study, simulated stable water isotope concentrations of HD 16 O and H2 18 O with isotope-enabled GCMs (section 2.1.2), 80 are regionally downscaled with COSMO_iso (Pfahl et al., 2012), an isotope-enabled version of the numerical weather prediction model COSMO (Baldauf et al., 2011) (version 4.18). For the purpose of long-term climate simulations, isotoperoutines of COSMO_iso were implemented in COSMO-CLM (Rockel et al., 2008), the climate version of COSMO. In this

Snow layer thickness
In the standard configuration of COSMO, the Greenland ice sheet is treated as a constant mass of ice, which is covered by a single snow layer. But in this model structure, dynamical processes within the ice sheet (flow, basal melt) are not included. As a result, the depth of the snow layer is constantly increasing and thus also its heat capacity. The model is therefore not able to reproduce the observed diurnal cycles of the 2m temperature at CIRES stations in Central Greenland. To avoid this spurious 100 model behavior, the snow layer depth is limited to 5 cm in this study. Using this value, realistic diurnal cycles of the 2m temperature could be simulated.

Marine regions with sea ice cover
To be able to simulate reasonable fractionation processes for marine regions with sea ice cover, a snow layer is also 105 implemented on top of the sea ice (e.g., as suggested in Bonne et al., 2019). The isotopic composition of this surface snow layer is in this case set to the isotopic composition of the most recent precipitation.

Fractionation at snow covered surfaces
Isotope fractionation during sublimation from a surface snow layer is poorly understood. Several different processes are 110 suggested to be involved, which are not yet taken into account in the models, such as non-fractionating layer-by-layer sublimation, kinetic fractionation during sublimation into sub-saturated air, a diurnal cycle of sublimation combined with fractionating vapor deposition on the snow, and fractionating melt water sublimation combined with recrystallization of residual melt water have been suggested (see e.g. discussion in Christner et al., 2017). In this study, best agreement of simulation results with data of ground-based in situ measurements of near-surface water vapour at a site near Karlsruhe, 115 Germany (measured with a Picarro water isotopologue analyser L2120-i), was found for the assumption of equilibrium fractionation during sublimation from surface layer snow and sea ice. However, the authors are aware that this description may only realistically approximate a much more complex interplay between different surface processes.

Model Simulation Setup 120
In a first step, the capability of COSMO_iso to realistically reflect the fractionating processes of stable water isotopes in Greenland is evaluated. For this, the nudged simulation outputs (standard and isotopic) from an isotope-enabled atmospheric model ECHAM5-wiso (Werner et al., 2011) simulation are dynamically downscaled with COSMO_iso for the whole Arctic region. The data from the same ECHAM5-wiso simulation have been already used as boundary conditions for COSMO_iso simulations over Europe by Christner et al. (2018). The simulation outputs from ECHAM5-wiso are at a T106 horizontal 125 spatial resolution (1.1° x 1.1°) and on 31 atmospheric vertical levels. The dynamical fields were nudged every 6 hours towards ERA-Interim reanalysis data (Dee et al., 2011). The simulation period is 2008-2014. The spatial resolution of COSMO_iso is https://doi.org/10.5194/cp-2019-156 Preprint. Discussion started: 20 January 2020 c Author(s) 2020. CC BY 4.0 License. 0.44° x 0.44°, corresponding to 50 km x 50 km in rotated coordinates (COSMO_iso_50km). Afterwards, the COSMO_iso_50km simulation is further downscaled to a spatial resolution of 0.0625° x 0.0625° (corresponding to about 7 km × 7 km) for Greenland (COSMO_iso_7km). This high-resolution simulation covers only the year 2011, which was chosen 130 for technical reasons. In the COSMO_iso runs, the horizontal wind fields above the 850 hPa level are spectrally nudged (von Storch et al., 2000) towards the reanalysis-based dynamical fields of ECHAM5-wiso. This method ensures that consistent atmospheric boundary conditions build the framework for the fractionating processes simulated in COSMO_iso.
The same model chain is applied to the mid-Holocene period. As no isotope-enabled ECHAM5-wiso simulation has been available for this period, atmospheric fields have been retrieved from a mid-Holocene simulation of the fully coupled model 135 MPI-ESM-wiso (Cauquoin et al., 2019), whose atmospheric component is ECHAM6-wiso. The major ECHAM6 model changes compared to ECHAM5 include an improved representation of radiative transfer in the solar part of the spectrum, an improved representation of surface albedo, a new aerosol climatology and an improved representation of the middle atmosphere (Stevens et al., 2013). With COSMO_iso, a representative time slice of 30 years is simulated for this climate period, only, since the regional COSMO_iso simulations are computationally very expensive. The greenhouse gas 140 concentrations and the orbital parameters are adapted, according to the Paleoclimate Modelling Intercomparison Project 4 experiment design (PMIP4, Kageyama et al., 2018). The model domain of the COSMO_iso simulations is identical to the present-day simulations.

Observations 145
The capability of the isotope-enabled regional climate model COSMO_iso to reproduce measured isotopic ratios in Greenland is evaluated by comparing the simulation results to different observational data. In this context, simulated δ 18 O ratios are compared to snow pit samples collected during the North Greenland Traverse (Fischer et al., 1998;Weißbach et al., 2016) and top core samples from four ice core locations (Renland (Vinther et al., 2008), Neem (Masson-Delmotte et al., 2015), GISP2 (Grootes and Stuiver, 1997), Summit (Fischer, 2003) Ice core samples are also used to evaluate the simulated isotopic ratios for the mid-Holocene. Beside the already mentioned Renland and GISP2 samples, two more ice core samples, namely GRIP and NGRIP (Vinther et al., 2006), are used for the 155 model evaluation. The mid-Holocene δ 18 O values (blue numbers in Table 1) are calculated as an average of the measured δ 18 O values in ice cores over the period 5.5 ka -6.5 ka.  In general, COSMO_iso in a 50 km x 50 km spatial resolution is able to reflect the observed isotopic ratios at the 16 snow pit samples and improves the simulation results of ECHAM5-wiso. In both simulations, the δ 18 O ratios are high near the coastline and low in Central Greenland. But in COSMO_iso_50km, the δ 18 O ratios decline stronger from the coastline to the inland  (Table 1). plateau than in ECHAM5-wiso. The spatial δ 18 O differences are consequently stronger pronounced and the general overestimation of δ 18 O ratios, which occurs in ECHAM5-wiso, is reduced in COSMO_iso_50km. Thus, the regional simulation 180 reaches a better agreement with the observations. Especially for the snow pit samples at which ECHAM5-wiso exhibits strong deviations to the observed δ 18 O values (1,3,4,6,7,8,9,10,16; see Table 1), a regional downscaling with COSMO_iso_50km reduces the bias considerably ( Figure 2). But for snow pit samples at which ECHAM5-wiso has already a high agreement with the observations (2,5,11,13,14), COSMO_iso_50km tends to increase the bias. A further downscaling with COSMO_iso to a spatial resolution of 7 km x 7 km does not improve the simulation results anymore. The only exception constitutes the snow 185 pit sample from Renland (1). Here, a considerable model bias in ECHAM5-wiso and COSMO_iso_50km is strongly reduced in COSMO_iso_7km.

Present-Day
The coastal area of Renland is characterized by a complex terrain and constitutes a special case for isotope-enabled modeling in Greenland. The snow pit sample is located in a transition zone from the homogeneous inland glaciation to the rugged coastline, where the glaciers calve into the sea. Thus, within short distances large differences in altitude and land surface 190 characteristics occur in this region. The isotopic ratios in the snow pit sample are therefore strongly affected by these heterogeneous local conditions, which are insufficiently represented in the coarse model resolution of ECHAM5-wiso. By increasing the spatial resolution with regional climate modeling, also the representation of the associated small-scale processes is improved. This leads generally to a higher agreement of the simulation results with observations, as seen for the COSMO_iso_7km run for Renland (Figure 2). 195 However, an increase in spatial resolution is, especially in complex terrains, also associated with an enlarged heterogeneity of the surface characteristics and the related small-scale processes. This is because the GCM grid boxes are further divided in smaller RCM grid boxes and consequently higher as well as lower values (for e.g. altitude) are included within the same area. In the case that a snow pit sample is located right at the border of two RCM grid boxes, the higher resolved model information does not necessarily have to match the local conditions in complex terrains. The characteristics of the adjacent RCM grid box 200 might fit much better to the conditions of the snow pit sample. The differences between the model data and the observations can consequently get larger than for a coarser model resolution. Thus, in the following, snow pit samples are not anymore solely compared to the model grid boxes covering the samples location. Instead, it is investigated whether the δ 18 O range of all adjacent RCM grid boxes to a snow pit location is consistent with the observed δ 18 O value of the same site. For this, all RCM grid boxes located within the corresponding GCM grid box are included in the comparison with the observations. 205

210
The spatial isotopic ratio variability of the COSMO_iso_50km grid boxes surrounding the 16 snow pit samples is shown as Box-Whiskers-plot in Figure 3a. The spatial isotopic ratio variability of the COSMO_iso_7km is shown in Figure 3b. In this spatial isotopic ratio variability, the δ 18 O values of all COSMO_iso (50 km) and COSMO_iso (7 km) grid boxes within the ECHAM5_wiso grid box closest to the snow pit sample, are included. For 11 of the 16 snow pit samples (1,2,3,4,5,6,7,8,14,15,16) the range of the COSMO_iso_50km grid box variability is consistent with the observed δ 18 O values. These stations are all located in the north of Greenland. In these latitudes, longitudes are converging and the differences in spatial resolution between ECHAM5-wiso and COSMO_iso_50km with its rotated pole become smaller. Thus, in the high 225 latitudes less COSMO_iso_50km grid boxes are covered by one ECHAM5-wiso grid box, than in the south of Greenland. The spread of the COSMO_iso_50km results is consequently reduced and the COSMO_iso_50km δ 18 O range at these stations is not consistent anymore with the observed δ 18 O values.
A downscaling to 7 km does slightly increase the spread of the COSMO_iso results. But still, the observed δ 18 O values from 5 of 16 snow pit samples are not covered within the modelled COSMO_iso_7km grid box variability (Figure 3b). Thus, a 230 further downscaling to a spatial resolution of 7 km does not increase the accuracy of the simulated isotopic ratio spread within a ECHAM5-wiso grid box. In accordance with the missing benefits of the COSMO_iso_7km simulation and its increased computing time costs, only a COSMO_iso_50km simulation is performed for the mid-Holocene.

Mid-Holocene 235
In Figure 4a the differences for the simulated MPI-ESM-wiso (blue) and COSMO_iso_50km (green) grid box results to the observed δ 18 O ratios at the corresponding ice cores are presented for the mid-Holocene. As in Figure 3, the spatial isotopic ratio variability of the COSMO_iso_50km grid boxes surrounding four Greenland ice core samples is shown as Box-Whiskersplot. MPI-ESM-wiso properly reflects the isotopic ratios of the mid-Holocene from ice core data. For the inland ice cores looking at the absolute biases, the downscaling does not seem to bring an added value to the MPI-ESM-wiso results for mid-255 Holocene conditions. But taking also into account the isotopic ratio variability in the COSMO_iso_50km simulation, the model results are again in agreement with the isotopic ratios of the ice core samples. The fact that, in contrast to the present-day simulations, only four observational data sets are available for the mid-Holocene, makes the assessment of the simulation results difficult. Moreover, the GRIP and GISP2 ice cores being located very close to each other (Figure 5), only three local isotope distributions clearly different from each other are available. Therefore in Figure   5a, the spatial isotopic ratio variability of the COSMO_iso_50km simulation is illustrated for whole Greenland. This spatial isotopic ratio variability is calculated as the standard deviation of all COSMO_iso_50km grid boxes within the respective MPI-270 ESM-wiso grid boxes. Figure 5a shows that at the coastline, the isotopic ratio variability of COSMO_iso_50km is considerably increased within the MPI-ESM-wiso grid boxes. In Central Greenland the increase in the isotopic ratio variability is lower.
Thus, the simulated variability is high in regions where large orographic differences occur within short distances, like the coastal areas of Greenland, and lower for homogeneous terrain like the inland plateau. The Renland ice core is therefore located in an area of a high isotopic variability, the GRIP and GISP2 ice cores in an area of low isotopic variability. Nevertheless, 275 regions with higher isotopic variability occur also in the inland plateau of Greenland. The NGRIP ice core, for instance, is located in an area of a moderate isotopic ratio variability. The four ice core drill sites are therefore located the three regions of Greenland with substantially different isotopic ratio variabilities.
The spatial surface temperature variability in the COSMO_iso_50km simulation, which is calculated in the same way as the isotopic ratio variability, is shown in Figure 5b. Near the coastline, the surface temperature variability is high, highlighting 280 how important the land surface characteristics are for the regional temperature variability. But in Central Greenland almost no variability occurs. This is in contrast to the water isotope ratio, exhibiting also regional variations over the inland plateau ( Figure 5a). The spatial distribution of δ 18 O does consequently not only depend on land surface processes, but also on dynamic atmospheric processes. In this way, isotopic ratios based on atmospheric fractionation processes along the trajectory of an air mass, are transported to Central Greenland and increase there the isotopic variability. 285 The same pattern occurs for the spatial isotope-temperature gradient (Figure 6a), representing the relation between the surface temperatures and the δ 18 O values. The spatial gradient is thereby calculated between the simulated δ 18 O ratios and surface temperatures at all COSMO_iso_50km grid boxes within the respective MPI-ESM-wiso grid boxes. At the coastline, the δ 18 Otemperature gradient is low, reflecting the high surface temperature and δ 18 O variability in this region. In Central Greenland the gradient is high, which is the result of the low surface temperature variability and the simultaneously increased δ 18 O 290 variability. This pattern does not only apply for the mid-Holocene, but also for the present-day COSMO_iso_50km simulation (Figure 6b). In the present-day simulation the contrast between the coastal regions and the inland plateau is even more clearly pronounced. The spatial variabilities of the surface temperature and the δ 18 O ratio in the mid-Holocene follow consequently the same mechanisms as for present-day.

Discussion and Conclusions 310
The results of several global paleo-climate simulations exhibit considerable deviations to the observed regional climate patterns during the Holocene (Braconnot et al., 2012). In the presented study, for the first time, regional climate simulations with an isotope-enabled RCM are performed for Greenland to potentially improve the agreement to climate observations in this region for the mid-Holocene. In a first step, the capability of the isotope-enabled RCM COSMO_iso to reproduce observed isotopic ratios for Greenland is demonstrated. 315 The COSMO_iso simulation results show that a spatial resolution of 50 km already leads to reasonable δ 18 O values. Especially in regions where the global ECHAM5-wiso model, which has been used to derive necessary forcing fields for the COSMO_iso simulations, deviates strongly from the observed δ 18 O values, the bias is considerably reduced by the regional climate simulation with COSMO_iso. In complex terrain like the coastal areas of Greenland, the results can be further improved with an additional downscaling to a spatial resolution of 7 km. In such simulations with high spatial resolution, small-scale processes 320 are described in more detail (e.g. Torma et al., 2015;Coppola et al., 2018) and thus the local characteristics at ice core sites are better taken into account (Sturm et al., 2005;Werner et al., 2011). But for northern Greenland, regional climate simulations with COSMO_iso even increase the deviations to the observations and the local COSMO_iso δ 18 O variability at the ice core drill sites is not consistent with the observed δ 18 O values. In this region, the gain in resolution is small for a downscaling and consequently the simulation results are not improved. But all in all, the results show that COSMO_iso is able to provide 325 reasonable isotopic ratios for Greenland and thus, can be applied for paleo-climate simulations.
For the mid-Holocene, MPI-ESM-wiso is in good agreement with observed ice core data in Greenland, as already described by Cauquoin et al. (2019). The model bias is, in this context, not further reduced by a downscaling with COSMO_iso. But an increase in the spatial model resolution leads also to an increase in the models degrees of freedom. This in turn can lead, in the case of an already good performance of the driving model, to additional noise and thus, a deviating RCM behavior with even 330 an increase in the absolute model bias, as seen for the mid-Holocene and the present-day simulations in Central Greenland.
Another consequence of these increased degrees of freedom in the COSMO_iso simulation is that the spatial variability of the simulated δ 18 O ratios is enhanced. Therefore, by analyzing the spatial variability of the COSMO_iso results, it can be demonstrated that most of the observed δ 18 O values lie within a range of spatially varying δ 18 O values, which can be derived in a physically consistent way by a regional downscaling of the spatially coarse MPI-ESM-wiso model results. This applies 335 for both, the present-day runs and the regional paleo-climate simulations for the mid-Holocene in Greenland. The deviation between the coarser resolved MPI-ESM-wiso results and the finer resolved observations is therefore potentially caused by the missing representation of important small-scale processes, which are induced by e.g. the surface conditions or orographic effects over Greenland.
As δ 18 O ratios are used as an indicator for temperatures in past climates (Dansgaard et al., 1969;Masson-Delmotte et al., 2005;340 Jouzel, 2013), it is important to understand how the presented COSMO_iso simulations might be able to improve these isotopebased temperatures reconstructions. In general, the regional surface temperature variability and the regional δ 18 O variability https://doi.org/10.5194/cp-2019-156 Preprint. Discussion started: 20 January 2020 c Author(s) 2020. CC BY 4.0 License.
show similar patterns for Greenland. In both cases the variability is high at the coast and low on the inland plateau. The same patterns as in the mid-Holocene can also be seen for the present-day simulations. In this context, the simulated spatial δ 18 Otemperature gradient reflects the same patterns as simulated by Sjolte et al. (2011) under present-day conditions for Greenland. 345 This spatial variability patterns indicate that the regional surface temperature variability highly depends on the surface characteristics. However, for the regional isotopic ratio variability, this dependence appears to be less pronounced. At the coastline, a clear relationship between surface temperatures and measured δ 18 O ratios in ice cores can be derived, while in Central Greenland this relation is less pronounced. These spatial differences might be explained by the fact that isotope changes are an integrated signal of the meso-scale variability of atmospheric processes (Dansgaard, 1964;Merlivat and Jouzel, 1979;350 Gat, 1996), which might partially be decoupled from surface temperature changes in homogeneous terrain.
The presented study demonstrates that the isotope-enabled MPI-ESM-wiso -COSMO_iso model chain with realistically implemented stable water isotope fractionation processes constitutes a useful supplement to reconstruct regional paleo-climate conditions during the mid-Holocene in Greenland. This approach might also be very helpful for other isotope-enabled GCMs and their deviations to observed isotope ratios in different paleo-time periods and regions. Particularly in regions, in which 355 large differences occur between simulated and observed δ 18 O ratios, like Europe and North America (Cauquoin et al., 2019;Comas-Bru et al., 2019), an improved representation of small-scale processes can potentially reduce these biases, and consequently, the reconstruction of regional paleo-climate patterns can become more reliable. To prove this hypothesis, in follow-up studies, more time slices will be simulated with the presented MPI-ESM-wiso -COSMO_iso model chain for different periods and different regions. 360