δ 18 O water isotope in the i LOVECLIM model ( version 1 . 0 ) – Part 1 : Implementation and verification

Introduction Conclusions References


Introduction
Water isotopes are widely used tracers of the hydrological cycle.With fractionation occurring at phase changes (evaporation, condensation, freezing, e.g.Craig and Gordon, 1965) and through diffusive processes at smaller scale, water isotopes are faithful recorders of the complex processes at work within the hydrological cycle.They have been used for decades in the field of palaeoclimate research to infer climatic conditions from the ice-cores (Dansgaard, 1964;Dansgaard et al., 1993;EPICA community members, 2004;North Greenland Ice Core Project members, 2004) or from oceanic sediment cores (Emiliani, 1955, for example) but may also be used at much smaller timescale and spatial scale to link climate variability and water isotope compositions (Kurita et al., 2011) or even to infer the mixing properties within rain events (Risi et al., 2010a).
From a climatic modelling perspective, the inclusion of water isotopes enable a thorough evaluation of the hydrological cycle in climate models, not only against precipitation amount or evaporation amount observations, but also on the actual transport of water through the atmospheric model.Applied to palaeoclimate simulations, it enables an accurate comparison of the model results with palaeoproxies, avoiding intermediate steps through derivation of temperature or salinities.Finally, it is an important procedure to assess the closure and adequate transport of water within models since any non-conserving process will likely lead to unrealistic δ 18 O of water.Since the pioneering work of Joussaume et al. (1984), much progress has been achieved in atmospheric general circulation models (Hoffmann et al., 1998;Noone and Simmonds, 2002;Mathieu et al., 2002;Tindall et al., 2009, for example) that can simulate quite accurately the δ 18 O of precipitation, even at fine resolution (Werner et al., 2011).Some secondary parameters like the deuterium excess have proven to be more challenging (Risi et al., 2010b).The development of water isotope-enabled isotopic models (Schmidt, 1998;Delaygue et al., 2000, for example) have further enabled the use of coupled isotope-enabled climate models that are then applied to palaeoclimate science questions (Legrande and Schmidt, 2008, for example).
Published by Copernicus Publications on behalf of the European Geosciences Union.

D. M. Roche: Water isotopes in the iLOVECLIMmodel
Given the computing resources needed to run coupled climate models, applying less complex coupled climate models with water isotopes to long-term palaeoclimate perspectives is still promising.The requirement, given the nature of isotope fractionation and distillation processes within the atmosphere (Craig and Gordon, 1965), is to explicitly compute the transport of water isotopes within the atmosphere.Earlier attempts (Roche et al., 2004, for example) have shown that such a perspective, though clearly not applicable for kilometre scale issues, could help towards a better understanding of palaeoproxy records.More recent developments (Brennan et al., 2012) showed that for present-day conditions, simplified modelling approaches could yield large-scale fields in accordance with data.The approach taken in the present study is somewhat similar, albeit with a more comprehensive treatment of water advection, precipitation and evaporation within the atmospheric component.Our study comprises three parts: (1) development and verification, (2) evaluation against water isotopes observations, and (3) evaluation against carbonate isotopes proxy data.
In the present manuscript, I present the design and verification of a δ 18 O water isotopes module in the iLOVECLIM climate model.I start with the equations needed to simulate the water isotopes in our simplified coupled climate model and proceed to the verification against well-known relationships of the climate system.In the two companion manuscripts, we present the model validation and evaluation, at first from the perspective of δ 18 O in water from present-day observations (Roche and Caley, 2013) and second with a palaeoperspective against late Holocene carbonate proxy data (Caley and Roche, 2013).

Technical description of the water isotopic scheme used in iLOVECLIM
The iLOVECLIM model is a code fork of the LOVECLIM-1.2 climate model extensively described in Goosse et al. (2010).From the original model, we retain the atmospheric (ECBilt), oceanic (CLIO), vegetation (VECODE) and landsurface (LBM) components and developed a complete, conservative, water isotope cycle through all cited components.With regards to water isotopes, the main development lies in the atmospheric component in which evaporation, condensation and existence of different phases (liquid and solid) all affect the isotopic conditions of the water isotopes.Hence, in the following I first describe extensively the method used to trace the water isotopes in ECBilt and only briefly their treatment in other components afterwards.

ECBilt-wiso: water isotopes tracking in the atmosphere
ECBilt is the simplified component of the iLOVECLIM earth system model.It is a quasi-geostrophic atmospheric model with some additional correction terms, described in details in Opsteegh et al. (1998) and Goosse et al. (2010).The atmosphere runs at a T21L3 resolution, that is approximately 5.6 • resolution in latitude and longitude.Of main interest for the purpose of developing a water isotopic module is the water cycle dynamics.ECBilt contains a full description of the water cycle from the evaporation to precipitation through condensation.The vertical structure is on three levels with only one humid layer (troposphere) and two dry layers (stratosphere).A schematic representation of the water cycle in ECBilt is given in Fig. 1.Evaporative water fluxes are added to the humid layer.Then vertical advection is computed.Since the two upper layers are dry, water fluxes crossing the boundary between the troposphere and the stratosphere are rained out as convective rain.If the water specific humidity of the humid layer is greater than a specific q sat value (set in ECBilt as 80 % of the saturation humidity at given temperature), the excess water is removed as large scale precipitation.Finally, if large scale precipitation occurs with negative temperatures, excess precipitation is removed as large scale snowfall.

Prognostic variable for water isotopes
For water isotopes, I follow roughly the same procedure.The prognostic variable for humidity in ECBilt is the quantity of precipitable water for the whole atmospheric column (in metres), q.It may be written as where m H 2 O is the mass of water in the given cell A the surface area of the cell, ρ the water density, n the number of moles of water and M H 2 O the molar mass.
The water isotope variable to be used is, by analogy, written here for 18 O: However, since the interaction of the different water isotopic species that form the water will not be dealt with, the previous formulation may be simplified to have the water isotopes as a simple tracer of water as in Merlivat and Jouzel (1979).Thus, the water isotopic quantity is expressed as where R 18 refers to the classical simplified form R 18 = n 18 n .
Let us now describe the isotopic changes throughout the water cycle, from evaporation to precipitation.

Isotopic evaporation
At the evaporation stage, the formulation for water is E = (q s − q v ), where q s is the surface specific humidity immediately above the ocean, q v the humidity of the first The abbreviation are as follow: q sat stands for saturation humidity and q for humidity; T i stands for limit temperature for snowfall and T 1 is the temperature of the first layer of the atmosphere; "Snowf."stands for Snowfall and "Precip" for precipitation.
atmospheric layer above the ocean and the drag coefficient depending at least on wind speed.For the water isotopes, the equivalent formulation is where the 18 denote the oxygen-18 related terms.To simulate the water isotopes in the evaporation, we need to determine q 18 s .However, in ECbilt, the terms q s and q v are computed from a climatological discretization on the vertical to take into account the effect of the planetary boundary layer.Since there is no equivalent vertical discretization for water isotopes, I cannot use the same procedure and need to rely on an approximate solution, computing first the water isotopic ratio in the evaporation.With that obtained and using the property, with the given definitions, we can write the computation of E 18 follows logically.Computing R 18 E requires some assumption on the processes occurring between the ocean and the atmosphere.I chose here to use the method introduced by Cappa et al. (2003) with a slight modification to account for our context.Cappa et al. (2003) assume that at the interface of the ocean, there is a thin layer in equilibrium with the ocean, overlaid by the planetary boundary layer that exchanges moisture with the free atmosphere and the previous thin equilibrated layer.I will not repeat the equations developed by these authors here since they apply to our case and only repeat the resulting final formulation for R 18 E : where h a is the relative humidity and R 18 a the isotopic ratio in the free atmosphere, R 18 eq is the isotopic ratio at equilibrium with the ocean and α * diff is the kinetic fractionation factor for the isotope considered.To obtain the previous formulation, Cappa et al. (2003) assume implicitly that the saturated humidity is the same at the ocean surface and in the free atmosphere.This is not the case for the ECBilt model.I thus need to slightly modify the previous equation, introducing the saturated specific humidity in both the free atmosphere and in the thin equilibrated layer above the ocean.An "apparent" relative humidity for the free atmosphere may be defined as where h a is the model free atmosphere relative humidity, q eq a the specific free atmosphere humidity and q eq s the specific humidity of the equilibrated thin layer above the ocean.Using the apparent relative humidity, Eq. ( 6) can finally be used to compute the isotopic ratio of evaporating moisture as where M. Roche: Water isotopes in the iLOVECLIMmodel D denoting the molecular diffusivity of water, D 18 the one of the respective isotope and n a coefficient that varies with turbulence and evaporative surface (Brutsaert, 1975;Mathieu and Bariac, 1996).

The
D 18 D coefficient is determined experimentally.In iLOVECLIM, the values determined by Merlivat (1978) for 18 O are used, i.e.
which fully determine the isotopic evaporation in our model.

Water isotopes in precipitation
As described previously, we have only one moist level and only two types of precipitation.A full description of the fractionation for precipitation would need to describe the processes from the distillation between the liquid droplets formation and atmospheric moisture to the re-equilibration of falling precipitation with its surrounding moist air.Since the ECBilt model does not allow such an implementation, I rely on a very simple approach.I assume that precipitation always forms in isotopic equilibrium with the surrounding moisture with instantaneous rainout to the surface.For convective liquid precipitation, I use the equilibrium isotopic ratios at the altitude of the tropopause whereas for large scale liquid precipitation, I use the mid-troposphere conditions.In the case of solid precipitation (snow), I consider it to be always in equilibrium with isotopic moisture at the tropopause.However, to account for enhanced kinetic fractionation that was reported in high latitude regions, the formula of Merlivat and Jouzel (1979) is used.Thus our fractionation scheme for large-scale or convective precipitation, and snow, may be summarised as where the fractionation coefficient α 18 l-v is taken from Majoube (1971b) as and kin α 18 s-v derives from α 18 s-v as in Merlivat and Jouzel (1979): where S is a function depending on temperature as where T is in • C and λ is a tunable parameter generally taken between 2 × 10 −3 and 4 × 10 −3 , taken here to be 4 × 10 −3 .
The equilibrium fractionation coefficient between water vapor and solid water is taken from Majoube (1971a) as which entirely determines our system of equations.The rationale behind this formulation is that the liquid precipitation, while formed at high altitude, does re-equilibrate partially with the surrounding water vapor during its fall (as shown in paired vapor/precipitation measurements, Araguas-Araguas et al., 2000).Hence, the apparent fractionation equilibrium is not the altitude of precipitation formation but somewhat lower in the atmospheric column.Since we do have only three layers, the mid-troposphere is the most appropriate choice, as confirmed a posteriori by our results.In the case of solid precipitation (snow), I consider it to be always in equilibrium with isotopic moisture at the tropopause.The case of solid precipitation is different to that of the liquid precipitation since the re-equilibration time is expected to be longer.

Water isotopes in other components
In the other components of the earth system, I assume to a first order approximation that the water isotopes act as passive tracers in the ocean and under equilibrium fractionation for the other components.

Land surface model
As precipitation falls on land surface, the water and water isotopes are added to the bucket water model.If reevaporation occurs, it is assumed to be formed under isotopic equilibrium with environmental conditions.Thus the ratio of oxygen-18 isotopes in re-evaporation is where T is the local surface temperature.If the amount of (isotopic) water in the soil bucket exceeds a threshold then water is routed instantaneously to the ocean following a simple routing scheme with pre-defined river basins.There is no fractionation associated with that process since there is no phase changes.Similarly, evapotranspiration occurs with equilibrium fractionation from the soil bucket water.Departing from this hypothesis would require to model what is occurring for leaf water: from the bucket water uptake in roots to the transpiration in leafs.This is clearly beyond the modelling scale we are attempting here.Further evolution of the model will test the necessity to use a simplified parametrization going beyond the presented simplistic assumption.
Finally, the snow layer is also represented as a bucket type: snow piles up until a threshold is reached.When additional snow is added, the snow is routed directly in the ocean following the same routing as for liquid water.As I do not deal with the accumulation of snow in different layers through the course of winter, the snow layer is assumed to be one well-mixed layer: additional snow precipitation modifies the δ 18 O content of all the layer.In turn, snow sublimation produces moisture with δ 18 O at the snow δ 18 O, that is the evaporated snow is assumed to be isolated from the rest of the snow layer, as is expected.

Ocean model
Water isotopes are passive tracers in the ocean.Since CLIO is a free surface Oceanic General Circulation Model, I took care to implement the isotopes so as to be mass conserving, following exactly what is done for salinity.This is especially important if one wants to conserve the water isotopes to salinity relationship.In the present initial version, it is assumed that there is no isotopic effect in relationship with sea ice.The actual measured fractionation is relatively small (on the order of 2 per mil (Craig and Gordon, 1965;Melling and Moore, 1995) with regard to the surface ocean) in comparison with the salinity effect.Since there might be a local effect in regions where sea ice is formed we plan to implement it in a later version.However, the current state of the sea-ice model does not allow to easily trace the water isotopic content and overcoming this limitation would need relatively extensive model development.

Simulation set-up
In the following, I present results of a 5000 yr equilibrium run under fixed pre-industrial boundary conditions.The atmospheric pCO 2 is chosen to be 280 ppm, methane concentration is 760 ppb and nitrous oxide concentration is 270 ppb.The orbital configuration is calculated from Berger (1978) with constant year 1950.The run is performed using presentday land-sea mask, freshwater routing and interactive vegetation.
With regards to the water isotopes, the atmospheric moisture is initialised at −12 per mil and the ocean at 0 per mil for δ 18 O.

Verification: atmospheric component
In order to assess implementation of all above fractionation factors, we now conduct a verification step for the atmosphere by checking if the model is able to reproduce the expected relationships between simulated δ 18 O and selected simulated climatological variables.

Rayleigh distillation
One of the simplest transformations that may occur is the socalled Rayleigh distillation.It is described as follows: starting from a moist air mass with a certain composition 0 R 18 v , the air mass progressively loses its water by equilibrium precipitation and immediate removal from the air mass.The air mass is assumed to be isolated through this process.Such a process allows to simply relate the humidity of the air mass at a certain point in the drying process to its δ 18 O compositionor the equilibrium δ 18 O composition of the next precipitation to be formed as where f is the fraction of remaining water (see Appendix A for the derivation of the equations).
In the following I compare the results for δ 18 O in precipitation to the results obtained in a theoretical Rayleigh distillation assumed to start from a 0 per mil as first condensate.

Comparison with model results
Figure 2 shows such a Rayleigh distillation for a constant temperature of 15 • C. Model results (plotted in a colour code representing their latitude) are, as expected, largely above the theoretical Rayleigh distillation line.This is expected since the moist air in the model is not isolated at first and thus is recharged over its course from the evaporative regions to the δ 18 O depleted regions with higher δ 18 O content from the surface (oceanic mainly).Overall, the evolution of δ 18 O in precipitation is following the Rayleigh distillation, which shows that the δ 18 O module computes the ratios in precipitation reasonably.
The points present below the theoretical line are not problematic since there are dry regions with mean temperatures for condensation lower than 15 • C that could be approximated by a colder theoretical Rayleigh distillation line.Moreover, the implementation presented is more complex than a simple distillation and thus does not have to fit exactly one particular theoretical line.
A more problematic part that clearly shows up in the modelling results are the points below −25 per mil from Antarctica.Contrary to the Rayleigh distillation theoretical curve that show a steady decrease in humidity with decreasing δ 18 O, our modelling results are showing an increase in humidity with decreasing δ 18 O.Since those points are coming from Antarctica, it is hard to imagine a likely moisture source with decreasing δ 18 O content over the continent.Antarctica is probably indeed what is the closest at large scale on earth from a Rayleigh distillation.Hence, this points to an inconsistency in the modelling set-up that is discussed further below.

Dansgaard relationship
Another well-known feature observed initially by Dansgaard (1964) is the local relationship between δ 18 O in precipitation and the mean annual temperature at the site.Using essentially high latitude sites (low mean annual temperatures) he found that the relationship was well approximated by the following linear approximation for mean annual values: where T surf. is in • C.He also noted that the relationship was not linear anymore for annual mean temperature above 15 • C. The Dansgaard equation was used extensively for palaeotemperature evaluation from δ 18 O measurements in palaeoclimate proxy data.However, it is unlikely to always be stable through time (Werner et al., 2000).
Using available data for δ 18 O in precipitation from the Global Network for Isotopes in Precipitation (IAEA, 2006), a linear fit on this larger dataset (not shown) can be computed, limiting also the temperature up to 15 • C. It results in a slightly lower slope and intercept obtained with an R 2 value of 0.96.The updated equation is δ 18 O precip.= 0.61 T surf.− 15.6. (21) Hence, the result is very close to the traditional Dansgaard equation.Using a second-order polynomial fit on all data, the fit can even be better and not limited to low annual mean temperatures.The obtained equation is then: with an R 2 value of 0.977.Figure 3 presents our modelling results within this framework.The two previous equations are represented together with a second-order polynomial fit on the modelled δ 18 O from our simulation with equation: with an R 2 value of 0.911.At first glance, the two secondorder polynomials are quite similar in shape over most of the range of the data, although our results are biased towards heavier δ 18 O values at identical mean annual temperature.Also, the model results yield flatter relationship for high annual mean temperature.Overall, this indicates that our implementation of δ 18 O in atmospheric moisture yields too low fractionation from oceanic source moisture towards drier regions, as shown from the lower slope of the fitted polynomial.Additionally, it is quite obvious that the modelled δ 18 O values below −35 • C are at odd with all fitted lines.All these datapoints are in Antarctica and highlight a clear issue with the fractionation and advection of water isotopes along the path from the lower latitude.Assuming that a realistic estimation of the modelled "Dansgaard relationship" is given by the fitted second-order polynomial on the modelled values at temperatures higher than −35 • C, it seems that the atmospheric moisture -and hence the δ 18 O in precipitation -  Dansgaard (1964).Coloured dots are from the model simulation, coloured after their latitude, given on the scale to the right ; grey dots is the fitted second-order polynomial to the model data for temperature above −30 • C; magenta dots are constructed from the temperature data in iLOVECLIM, using the regression from Dansgaard (1964); yellow dots are constructed from the temperature data in iLOVECLIM using a second order polynomial fitted to the GNIP dataset.R 2 values for the fit to the iLOVECLIM model is 0.911 (grey dots); R 2 values for the fit to the GNIP data is 0.977 (yellow dots).
is modified by a source with higher δ 18 O content.Assuming that these anomalous datapoints are on a mixing line between the fitted polynomial and a source of moisture implies that the contaminating source has an isotopic signal of +5 to +10 per mil.Since there is no such δ 18 O rich source over Antarctica, it is necessary to conclude that the mixing is only apparent and that the cause is to be found in the numerical advection scheme of ECBilt.Such analysis is reinforced by the already noted bias in humidity over Antarctica at very low δ 18 O content in precipitation (cf.Fig. 2).
Additional checks performed (not reported here) show that indeed in the case of very low humidity content, the numerical advection scheme is not fully conservative, in isotopes and in water moisture and results in absurd δ 18 O values in precipitation.Correction of such bias will need relatively thorough analysis of the numerical scheme of ECBilt that is beyond the scope of the present study.It is however noteworthy that the presented implementation of δ 18 O yields a very positive result enabling detection of some defects in extremely dry climatic regions, a fact that was ignored so far.
To summarise, we can state that -apart from Antarcticathe relationship modelled between δ 18 O in precipitation and surface temperature bears a strong resemblance with what is expected from data inferences, albeit with a lower slope.This lower slope is probably a consequence of the simplicity of the model with the discretization of atmospheric moisture in a single layer.The fit is also better with a second-order polynomial in the model world than with a linear relationship, as is the case for observations.

Annual δ 18 O amplitude
So far, only annual mean values were reported.Since many of the proxy records may have a seasonal bias or record seasonal changes, it is important to have a look at basic features of the seasonal cycle, to verify whether the model is capable of reproducing some aspects of the yearly variations.I choose to present the maximum δ 18 O amplitude in precipitation from monthly data to evaluate the geographical variations of the amplitude of the seasonal cycle.
The first evident feature from Fig. 4 is the contrast between the continental and the oceanic regions.Over the oceans, the amplitude remains low.Tropical areas are characterised by 1 per mil amplitude, since the very active evaporation brings moisture in all year round, buffered by the ocean δ 18 O around 0 per mil.The amplitude increases towards the higher latitudes, but also in the equatorial regions.The maximum amplitude is reached over the Arctic ocean where cold and dry winters, with sea-ice covered ocean yield very low δ 18 O values and the retreat of sea ice and influence of warmer sea-surface conditions increases the δ 18 O during the rest of the year.
Looking at the longitudinal evolution, we see that the amplitude over the continent is several per mil higher than in the adjacent ocean, an expected result since over the continent there is no buffering effect from a large water body as the ocean.Over Africa for example, the amplitude is comprised between 4 and 9 per mil, while over the neighbouring ocean is between 2 and 4 per mil.Such a result is encouraging and shows that our model capture correctly some expected largescale patterns of climate and δ 18 O variations.Again, it shows that the implementation of the water isotopes seems to function correctly in the atmospheric component.
Comparing the modelled results with the amplitude derived from the GNIP database (IAEA, 2006), there is an overall good agreement between the two.Both observational data and simulation show an enhanced seasonal cycle in δ 18 O of precipitation over the continents and a dampened seasonal cycle over the oceans.There is also a continentality gradient observed both in model and observations, with a tendency to higher seasonal amplitude for higher continentality.The model even obtain a good representation of the minimum values over the ocean (around 0.9 per mil over the tropical oceans) and a high amplitude for ice covered regions like Antarctica.

Verification: oceanic component
As the focus is on verification of the implementation within the oceanic model and not on validation of the oceanic model itself, I will concentrate on a well-known relationship of the observed ocean: the δ 18 O-salinity relationship.Since both salinity and δ 18 O are affected by the balance between evaporation (that extract freshwater from the ocean) and input terms (precipitation, runoff, etc. providing freshwater), it is logical that a certain relationship exists between the two terms (Craig and Gordon, 1965).However, the water extracted by evaporation does not include salt at all (it is pure freshwater) whereas it contains 18 O.Thus the relationship is likely to break in hydrologically very active regions that is where the local hydrological balance is dominant over the surface advection and mixing terms.Similarly in the small regions where the sea-ice formation is dominant, the rejection of large amounts of salt into the surface oceanic waters with no similar 18 O counterpart may alter this relationship.
Figure 5 present results from the simulation for the surface Atlantic data.The distribution within the δ 18 O-salinity space shows a good correlation between the two variables within the range simulated.Using a linear regression, we obtain a slope of 0.43 to be compared to 0.52 when the same linear regression is performed on modern observations (Schmidt et al., 1999) of the GISS database (not shown).The agreement between the two slopes is excellent, showing our ability to simulate correctly both salinities and δ 18 O as described in the implementation.The modelled intercept is −14.9 per mil while the one calculated on data is −17.9.The end member of the modelling results are thus three per mil too high.Since the δ 18 O of precipitation is also positively biased around values of −20 per mil by about 5 per mil (see Fig. 3), the positive bias toward depleted δ 18 O values is expected to show also in the oceanic values.The fact that our slope is also underestimated by 0.1 per mil (δ 18 O)/per mil (salinity) is also related to the same phenomenon.Linear fit to the model data is given in the equation in the bottom right box, with a R 2 coefficient of 0.877.Performing the same analysis on Atlantic ocean data from the GISS database, we obtain a linear regression of 0.52 × −17.9 with a R 2 of 0.81.The colour scale applied on modelled points is the given latitude of each point.

Conclusions
In the present study, I have presented the design, implementation and verification of a δ 18 O water isotopic module in a simplified, single moisture layer, atmospheric model.The verification step showed that the implementation was successful and that well-known relationships of δ 18 O with other climatic variables (temperature, humidity) are well represented.I have shown as well that the moisture content was not fully conserved over Antarctica causing unrealistic results for δ 18 O in precipitation over that region.Since the replacement of the advection module of ECBilt on its Gaussian grid is beyond the scope of the present study, we have analysed the biases that such an non-conservation caused.Though problematic for future use of the isotopic model over Antarctica where some of the most well-known data is recorded in ice cores (e.g.EPICA community members, 2004), this example is interesting to point out the benefits of water isotope modelling already for the validation of climate models.
The iLOVECLIM climate model was developed further with the introduction of an additional scheme to simulate the same δ 18 O water isotope in the ocean, land-surface and vegetation part of our fully coupled climate model.The analysis of the results from an oceanic perspective showed a good accordance with observation-derived relationships, though still presenting similar biases as the ones detected in the atmospheric part.
Overall, the model seems to perform adequatedly when its simplicity is taken into account.The availability of water isotopes in such a fast model open wide prospects for longterm palaeoclimate simulations.
Fig. 2. δ 18 O -humidity relationship in iLOVECLIM compared to a theoretical Rayleigh distillation model.The colour scale applied to model points indicates the latitude.The Rayleigh distillation curve is given for a temperature of 15 • C starting from an oceanic-like value of 0 per mil, as an example.The inset details the anomalous humidity values discussed in the text.

Fig. 3 .
Fig. 3. Annual δ 18 O -Temperature relationship in iLOVECLIM compared to the original regression fromDansgaard (1964).Coloured dots are from the model simulation, coloured after their latitude, given on the scale to the right ; grey dots is the fitted second-order polynomial to the model data for temperature above −30 • C; magenta dots are constructed from the temperature data in iLOVECLIM, using the regression fromDansgaard (1964); yellow dots are constructed from the temperature data in iLOVECLIM using a second order polynomial fitted to the GNIP dataset.R 2 values for the fit to the iLOVECLIM model is 0.911 (grey dots); R 2 values for the fit to the GNIP data is 0.977 (yellow dots).

Fig. 5 .
Fig.5.Near surface ocean annual δ 18 O-Salinity relationship in iLOVECLIM for the Atlantic ocean.Linear fit to the model data is given in the equation in the bottom right box, with a R 2 coefficient of 0.877.Performing the same analysis on Atlantic ocean data from the GISS database, we obtain a linear regression of 0.52 × −17.9 with a R 2 of 0.81.The colour scale applied on modelled points is the given latitude of each point.

One humid layer Two dry layers C o n v e c t i v e P r e c i p i t a t i o n Large scale Precip. Evaporation. Large scale Snowf. T 1 >T i T 1 <T i q>q sat q>q sat Ocean surface layer Fig. 1. Schematic
representation of the water cycle in ECBilt.