Articles | Volume 16, issue 6
Clim. Past, 16, 2573–2597, 2020
Clim. Past, 16, 2573–2597, 2020

Research article 23 Dec 2020

Research article | 23 Dec 2020

The middle to late Eocene greenhouse climate modelled using the CESM 1.0.5

The middle to late Eocene greenhouse climate modelled using the CESM 1.0.5
Michiel Baatsen1, Anna S. von der Heydt1, Matthew Huber2, Michael A. Kliphuis1, Peter K. Bijl3, Appy Sluijs3, and Henk A. Dijkstra1 Michiel Baatsen et al.
  • 1Department of Physics, Institute for Marine and Atmospheric Research, Utrecht University, Princetonplein 5, 3584CC Utrecht, the Netherlands
  • 2Department of Earth, Atmospheric, and Planetary Sciences, Purdue University, 610 Purdue Mall, West Lafayette, IN 47906, USA
  • 3Department of Earth Sciences, Marine Palynology and Paleoceanography, Laboratory of Palaeobotany and Palynology, Utrecht University, Princetonlaan 8a, 3584CB Utrecht, the Netherlands

Correspondence: Michiel Baatsen (


The early and late Eocene have both been the subject of many modelling studies, but few have focused on the middle Eocene. The latter still holds many challenges for climate modellers but is also key to understanding the events leading towards the conditions needed for Antarctic glaciation at the Eocene–Oligocene transition. Here, we present the results of CMIP5-like coupled climate simulations using the Community Earth System Model (CESM) version 1. Using a new detailed 38 Ma geography reconstruction and higher model resolution compared to most previous modelling studies and sufficiently long equilibration times, these simulations will help to further understand the middle to late Eocene climate. At realistic levels of atmospheric greenhouse gases, the model is able to show overall good agreement with proxy records and capture the important aspects of a warm greenhouse climate during the Eocene.

With a quadrupling of pre-industrial concentrations of both CO2 and CH4 (i.e. 1120 ppm and ∼2700 ppb, respectively, or 4 × PIC; pre-industrial carbon), sea surface temperatures correspond well to the available late middle Eocene (42–38 Ma;  Bartonian) proxies. Being generally cooler, the simulated climate under 2 × PIC forcing is a good analogue for that of the late Eocene (38–34 Ma;  Priabonian). Terrestrial temperature proxies, although their geographical coverage is sparse, also indicate that the results presented here are in agreement with the available information.

Our simulated middle to late Eocene climate has a reduced Equator-to-pole temperature gradient and a more symmetric meridional heat distribution compared to the pre-industrial reference. The collective effects of geography, vegetation, and ice account for a global average 5–7 C difference between pre-industrial and 38 Ma Eocene boundary conditions, with important contributions from cloud and water vapour feedbacks. This helps to explain Eocene warmth in general, without the need for greenhouse gas levels much higher than indicated by proxy estimates (i.e. ∼500–1200 ppm CO2) or low-latitude regions becoming unreasonably warm. High-latitude warmth supports the idea of mostly ice-free polar regions, even at 2 × PIC, with Antarctica experiencing particularly warm summers. An overall wet climate is seen in the simulated Eocene climate, which has a strongly monsoonal character.

Equilibrium climate sensitivity is reduced (0.62 C W−1 m2; 3.21 C warming between 38 Ma 2 × PIC and 4 × PIC) compared to that of the present-day climate (0.80 C W−1 m2; 3.17 C per CO2 doubling). While the actual warming is similar, we see mainly a higher radiative forcing from the second PIC doubling. A more detailed analysis of energy fluxes shows that the regional radiative balance is mainly responsible for sustaining a low meridional temperature gradient in the Eocene climate, as well as the polar amplification seen towards even warmer conditions. These model results may be useful to reconsider the drivers of Eocene warmth and the Eocene–Oligocene transition (EOT) but can also be a base for more detailed comparisons to future proxy estimates.

1 Introduction

The Eocene–Oligocene transition (EOT) is one of the most dramatic climate transitions of the Cenozoic, thought to be associated with the formation of a continental-scale ice sheet on Antarctica (Zachos et al.1994; Coxall et al.2005; Lear et al.2008). A possible cause for the inception of ice is a long-term decline of greenhouse gas concentrations through the middle Eocene, eventually crossing a threshold for glaciation (DeConto and Pollard2003; DeConto et al.2008; Gasson et al.2014). Following the early Eocene (∼50 Ma), a gradual cooling levelled off during the middle Eocene (43–42 Ma) and eventually reversed into a warming (Zachos et al.2001, 2008; Bijl et al.2009; Cramwinckel et al.2018) towards the Middle Eocene Climatic Optimum at ∼40 Ma (MECO; Bohaty and Zachos2003; Bijl et al.2010; Sluijs et al.2013). The cooling trend continued during the late Eocene (∼38–34 Ma), with a cold interval at ∼37.3 Ma characterised by the Priabonian oxygen isotope maximum (PrOM; Scher et al.2014). While these temperature changes may have caused some ice growth as early as the middle Eocene, they did not allow the formation of a continental-scale Antarctic ice sheet to occur until after 34 Ma (Scher et al.2014; Passchier et al.2017; Carter et al.2017). It remains a question to what extent continental geometry (e.g. opening of Southern Ocean gateways), next to gradual shifts in both the atmospheric and oceanic circulation, was a driver for both regional and global climate change during the Eocene (Bijl et al.2013; Bosboom et al.2014; Goldner et al.2014; Sijp et al.2014, 2016; Toumoulin et al.2020).

The climate throughout most of the Eocene was characterised by a reduced Equator-to-pole temperature gradient compared to present day (Bijl et al.2009; Hollis et al.2012; Douglas et al.2014; Evans et al.2018). This aspect of the Eocene greenhouse climate has proven challenging to adequately simulate with global climate models (Huber and Sloan2001; Huber and Caballero2011; Cramwinckel et al.2018). Very high greenhouse gas concentrations were often needed to reproduce high-latitude warmth at the expense of equatorial temperatures being significantly higher than indicated by proxy data (Huber and Caballero2011; Lunt et al.2012). Estimates of the meridional temperature gradient during the Eocene have come up with respect to some earlier studies, mainly due to warmer equatorial temperatures (Pearson et al.2007; Schouten et al.2013; Inglis et al.2015; Evans et al.2018; Cramwinckel et al.2018). Meanwhile, global climate models have been under continuous development by including more processes (especially cloud properties; e.g. Abbot et al.2009; Kiehl and Shields2013) and using a higher spatial resolution with better resolved palaeogeographies (Baatsen et al.2016; Lunt et al.2016; Hutchinson et al.2018) to improve the overall model–proxy comparison.

Many modelling studies have focused on the early Eocene (Huber and Caballero2011; Lunt et al.2012; Herold et al.2014; Zhu et al.2019) or looked at the latest Eocene–early Oligocene (Hill et al.2013; Ladant et al.2014; Kennedy et al.2015; Elsworth et al.2017; Hutchinson et al.2018; Kennedy-Asser et al.2019, 2020; Tardif et al.2020; see also the overview of Gasson et al.2014). However, there has been limited attention to the middle Eocene. Licht et al. (2014) performed simulations with both a 40 and 34 Ma geography reconstruction using a reduced-complexity global climate model (FOAM). A set of comprehensive model studies with several time slices covering the Eocene (Inglis et al.2015; Lunt et al.2016; Farnsworth et al.2019) use a lower-resolution version of the HadCM3 model. Cramwinckel et al. (2018) provide an overview of different climatic states throughout the middle Eocene cooling using the CCSM3 simulations of Goldner et al. (2014).

Here, we present the results of a set of coupled atmosphere–ocean simulations with the Community Earth System Model (CESM) version 1.0.5 using 38 Ma boundary conditions from Baatsen et al. (2016). The 38 Ma boundary conditions differ from those used to study the EOT, primarily in the treatment of Southern Ocean gateways (i.e. narrow, shallow passages) and an open Turgai Strait. With time-specific geographic boundary conditions and adequate equilibration, the aim is to show a more detailed and representative model-based overview of the middle to late Eocene (i.e. Bartonian–Priabonian; ∼42–34 Ma) climate. The considered period is suitable to investigate both the warm greenhouse climate of the middle Eocene and the conditions leading up to the EOT. The focus will therefore be on the general features of the modelled climate, a comparison to proxies and other model results, and the similarities and differences with respect to the present-day climate (specifically regarding climate sensitivity).

The model set-up and experimental design of the CESM simulations are first explained in Sect. 2. The main results are then presented in Sect. 3 including the equilibrium climate of each simulation (Sect. 3.1), a model–proxy comparison (Sect. 3.2), and a model–model comparison (Sect. 3.3). This is followed by an analysis of equilibrium climate sensitivity (Sect. 3.4) with a focus on the main changes involved in the radiative balance. Finally, the main findings are summarised and discussed in Sect. 4.

2 Methods

2.1 The CESM 1.0.5

The Community Earth System Model (CESM; Hurrell et al.2013) version 1 is a fully coupled atmosphere–land–ice–ocean general circulation model (GCM) that was developed at the National Center for Atmospheric Research (NCAR) in Boulder, Colorado. For use in palaeoclimate modelling, version 1.0.5 of the CESM is a suitable choice motivated by a trade-off between increasing model complexity and computational cost. This version of the CESM is equivalent to the latest version of the Community Climate System Model (CCSM4; Blackmon et al.2001; Gent et al.2011). In the configuration used here, the horizontal resolution of both the atmosphere and ocean is doubled compared to that of Herold et al. (2014), while that of the atmosphere is about 50 % higher than in Hutchinson et al. (2018), with comparable ocean grids.

The atmospheric component of the CESM is the Community Atmosphere Model (CAM4; Neale et al.2013), which uses a finite-volume grid at a nominal resolution of 2 (2.5× 1.9) and 26 vertical levels with a hybrid sigma vertical coordinate extending upward to 2 hPa. In this configuration, the model has a reported a warming response of 3.13 C to a doubling of CO2 starting from pre-industrial conditions (compared to ∼2.5C in CCSM3; Bitz et al.2012; Kiehl et al.2006).

The physical, chemical, and biological processes taking place on land are represented in the Community Land Model (CLM4; Oleson et al.2010; Lawrence et al.2011). A static rather than dynamic vegetation model is used here to avoid runaway feedback effects, which can become an issue, especially in warm greenhouse climates (e.g. dieback of vegetation at low latitudes; Loptson et al.2014; Herold et al.2014). The considered biomes are translated into fractions of the corresponding CLM4 plant functional types (PFTs), from which a set of monthly forcing files is finally used in the model. Freshwater runoff is treated by a simple river-routing scheme, in which all runoff is transported to one of the surrounding eight model grid cells until the ocean is reached. The direction is determined by the local topography gradient and manually adjusted where the runoff scheme would otherwise form closed loops.

The sea ice component consists of the Los Alamos National Laboratory (LANL) Community Ice Code version 4 (CICE4; Hunke and Lipscomb2004). For simplicity, sea ice only forms when the sea surface cools down to −1.8C, after which its dynamical behaviour (e.g. melt and advection) is specifically treated by the model.

The CESM1 uses the LANL Parallel Ocean Program version 2 (POP2; Smith et al.2010) for the ocean model component. The standard configuration is applied here, with a nominal 1 (1.25× 0.9) horizontal resolution on a curvilinear grid placing the northern pole over Greenland. The POP2 model is set up with 60 layers of varying thickness between 10 m near the surface and 250 m at greater depth. Horizontal viscosity is considered anisotropic (Smith and McWilliams2003) and horizontal tracer diffusion follows the parameterisation of Gent and Mcwilliams (1990). The model further uses the K-profile parameterisation (KPP) scheme to determine vertical mixing coefficients (Large et al.1994). More information and discussion on the ocean model physics and parameterisations can be found in Danabasoglu et al. (2008, 2012).

2.2 Model experiments

To allow comparisons of various climatic features between Eocene and pre-industrial conditions within a similar framework, a pre-industrial reference run with the same version of the CESM is performed (using the geography shown in Fig. 1b, e, f). Since the atmospheric component used here consists of the CAM4 at a nominal 2 horizontal resolution, this reference is similar to the 2 1850 control of Gent et al. (2011). The solar constant in the simulation is 1361.27 W m−2, and the atmospheric levels of CO2 and CH4 are 280 ppm and 671 ppb, respectively (i.e. pre-industrial carbon; PIC). Astronomical orbital parameters are set to their present-day configurations: eccentricity 0.0167, obliquity 23.44, and precession placing the aphelion in Northern Hemisphere summer. Both vegetation and atmospheric aerosols are kept fixed at their respective pre-industrial distributions (Fig. S1), meaning that all anthropogenic influences apart from land cover changes are disregarded. Overflow and tidal mixing parameterisations are switched off for a better comparison to the Eocene cases. Background vertical diffusivity is a constant 0.17 cm2 s−1 nearly everywhere, increasing to 0.3 cm2 s−1 around 30 N or S and 1 cm2 s−1 in the Banda Sea, following Jochum (2009).

Figure 1Topography (atmosphere) and bathymetry (ocean) grid used in the model simulations for the (a) 38 Ma Eocene cases and (b) pre-industrial reference. Black lines are drawn every 20 model grid boxes for the ocean, and the thick grey line shows the coastlines defined as the 0.5 land fraction. Polar stereographic projections of the Eocene model geography are shown for (c) the Arctic and (d) the Antarctic region, similarly for the pre-industrial one in (e) and (f), respectively. Dashed black lines indicate the coastlines of the alternative geography, and light blue contours show the edges of ice sheets.

For deep-time climate simulations there are a number of model parameters and settings that need to be reconsidered, while others are left unchanged with respect to the pre-industrial reference. Apart from a newer model version and overall increases in resolution, our Eocene CESM configuration is similar to that of the CCSM3 simulations from Herold et al. (2014) and based on the NCAR-suggested deep-time standards (found at:, last access: 29 November 2020). This means that overflow and tidal mixing parameterisations in the ocean are switched off and the background vertical diffusivity κw is horizontally homogeneous but dependent on model level depth z as follows: κw=vdc1+vdc2tan-1((|z|-dpth)linv), where vdc1= 0.524 cm2 s−1, vdc2= 0.313 cm2 s−1, dpth= 1000 m, and linv= 4.5×10−3 m−1. The solar constant is also slightly reduced from its pre-industrial value to 1360.89 W m−2 to reflect Eocene conditions. Since our Eocene model simulations are designed to reconstruct the mean climate over longer periods in time (4–8 Myr), there is not a constant set of orbital parameters that is representative for this entire time interval. A set of parameters generally conducive for Antarctic ice growth is chosen: minimum eccentricity (i.e. 0; cancelling the effect of precession) and present-day obliquity (23.44). The low eccentricity choice is motivated by the conjunction of two such minima occurring around the EOT, as shown by Coxall et al. (2005) and DeConto et al. (2008).

Our Eocene simulations use the 38 Ma geography reconstruction from (Baatsen et al.2016; Fig. 1a, c, d), which in contrast to most previously used geographies is based on a palaeomagnetic (“PaleoMag”) reference frame. This method prioritises the reconstruction of exact palaeolatitudes, crucial for palaeoclimate simulations (van Hinsbergen et al.2015). As for any other geography reconstruction, it comes with its own limitations and uncertainties (Baatsen et al.2016). The model vegetation used here (Fig. S1) is largely based on reconstructed biomes of Sewall et al. (2000) and comparable to the early Eocene vegetation used by Herold et al. (2014). The mostly zonal bands in Sewall et al. (2000) are adjusted to the new geography reconstruction and expected maritime influences before being translated into the CLM4 plant functional types (PFTs; see also Table S1). This translates into most of the Earth being covered by various types of forests and shrubs (tropical forest or savannah in low and mixed forests at middle to high latitudes). Neither desert regions nor any land ice coverage are incorporated, while elevated surfaces are covered with a separate biome. These choices are made specifically for the middle to late Eocene period, for which proxy data suggest the absence of any persistent large-scale deserts or ice sheets.

Despite being an important contribution to the radiative forcing in present and future climate simulations, atmospheric aerosols are tricky to include and are therefore often either omitted or assumed to be similar to pre-industrial in palaeoclimate model simulations. A considerable improvement can be made by running a bulk aerosol model (BAM; Heavens et al.2012) version of the CAM4 to determine a more realistic distribution of aerosols. The sources of naturally formed aerosols (mainly dust, sea salt, and organic carbon, excluding volcanic emissions) are adjusted from pre-industrial levels and redistributed based on the new land surface properties. Using these sources the (stand-alone) CAM4 is run for 50 years, at the end of which a monthly climatology of aerosol distributions is derived (see Fig. S1 for the corresponding annual mean aerosol optical depth).

Two 38 Ma simulations are carried out with 2 × PIC and 4 × PIC, respectively, to cover the most likely range of atmospheric greenhouse gases during the middle to late Eocene (∼500–1200 ppm; Beerling and Royer2011; Anagnostou et al.2016;, last access: 29 November 2020) and to estimate climate sensitivity. CO2 and CH4 concentrations are increased simultaneously, in line with middle to late Eocene estimates (Beerling et al.2009, 2011; Goldner et al.2014). According to Etminan et al. (2016) the radiative forcing of 2 × PIC and 4 × PIC is equivalent to that of 2.15 × CO2 and 4.69 × CO2, respectively. An overview of the different model cases analysed in this paper is given in Table 1.

(Baatsen et al.2016)(Sewall et al.2000)

Table 1Overview of characteristics for all CESM 1.0.5 simulations that were performed (PFT: plant functional type, BAM: bulk aerosol model).

Download Print Version | Download XLSX

2.3 Spin-up procedure

Our 38 Ma Eocene simulations are both initialised using a stagnant ocean with a horizontally homogeneous temperature distribution, decreasing with depth from 15 C at the surface to 9 C at the bottom. The pre-industrial reference is initialised using present-day temperature and salinity fields from the PHC2 dataset (Steele et al.2001). A long spin-up (see Table 1) is performed to allow the deep ocean to equilibrate sufficiently. An overview of absolute (ΔT, ΔS) and normalised (ΔTT, ΔSS) drifts over the last 200 model years is given in Table 2 for each spin-up. Drifts are generally 10-5–10−4 K per year for the global average volume-weighted ocean temperature (values of <10−4 are often regarded as well equilibrated; see e.g. Goldner et al.2014; Hutchinson et al.2018). Globally averaged salinity drifts at the end of all three spin-ups are 10-7 psu per year.

Table 2Overview of drifts in globally averaged (full-depth, volume-weighted) ocean temperature and salinity over the last 200 model years for all three CESM spin-up simulations. Normalised drifts ΔTT and ΔSS, with T temperature in Kelvin and S salinity in practical salinity units (psu), are also shown for each case using the same interval.

Download Print Version | Download XLSX

Time series of globally averaged upper ocean (<1 km), deep-ocean (>2 km), and full-depth ocean temperatures are shown for all three spin-up runs in Fig. 2a. Starting from present-day initial conditions, the pre-industrial simulation cools down by about 0.5 C globally. The 38 Ma 2 × PIC simulation is seen to equilibrate much faster than the 4 × PIC one, probably because the deep-ocean equilibrium temperature is close to that of the initial state (∼9C). With higher greenhouse gas concentrations, the ocean experiences additional heating at the surface, causing it to become more stably stratified which consequently reduces vertical mixing into the deep ocean. This reduced mixing in the 4 × PIC case causes the global average deep-ocean temperature to still increase about 0.1 C over the last 500 years. As expected, changes in the upper 1 km temperatures settle more quickly and are small (<0.1C) during the last 2000 model years of each simulation.

Figure 2(a) Time series of global upper 1 km (red), below 2 km (black), and full-depth (purple) volume-weighted average temperature for the complete pre-industrial reference (dotted), 38 Ma 2 × PIC (dashed), and 4 × PIC (solid) spin-up simulations. (b) Global average near-surface air temperature for the 38 Ma cases (markers: annual, grey crosses: 2 × PIC, black dots: 4 × PIC, lines: smoothed–dashed dark blue: 2 × PIC, solid light blue: 4 × PIC), with mean values over the last 200 years given. (c) Globally integrated top-of-model net fluxes for both 38 Ma Eocene cases using the same marker and line conventions as in (b).


While generally adjusting more quickly than the ocean, near-surface air temperatures can also be seen to equilibrate faster in the 38 Ma 2 × PIC case compared to the 4 × PIC one. A small but abrupt temperature drop is seen at model years 3000 and 4000 for the 2 × PIC and 4 × PIC case, respectively. This is caused by a regional shift in vegetation cover that was necessary to adjust regions of erroneously switched PFTs in the Eocene boundary conditions (Fig. S1a shows the correct vegetation cover used at the end of the 38 Ma simulations). Since most of the land surface is covered with different types of forest and the complexity of interactions with the modelled climate system is limited, the overall effect is minor. Both 38 Ma simulations were extended by another 600 years, showing little change after a quick initial adjustment. Correcting the vegetation cover mainly affects the albedo (a 1.8 × 10−3 increase globally), lowering the net fluxes at the top of the model atmosphere (with a net radiative forcing of about −0.1 W m−2 and a cooling of 0.2–0.3 C; Fig. 2c). Time series of these fluxes also indicate that the atmosphere is close to radiative equilibrium towards the end of both simulations, with a <0.1 W m−2 imbalance between globally integrated shortwave and longwave fluxes. The slightly positive radiative balance agrees with the deep ocean continuing to warm up slowly, while both the upper ocean and atmosphere show little change after model year ∼1500 for the 38 Ma 4 × PIC simulation.

In addition, global patterns of ideal age tracers (a measure for oceanic ventilation timescales; Thiele and Sarmiento1990; England1995) are observed to equilibrate over the last 1000 years of each simulation (Supplement Fig. S2) as are the meridional overturning strength and oceanic gateway transports (Supplement Fig. S3). Average values of ideal age still have considerable trends because of further ageing in stagnant deep-ocean regions (Fig. S2a) and upwelling of older water masses into the upper ocean (Fig. S2b). Such changes will likely continue for many thousands of years, but the associated circulation pattern has equilibrated by the end of each simulation, which is used for further analysis. The difference in spin-up between the two 38 Ma cases is once again highlighted in the evolution of maximum overturning strength. Whereas the 2 × PIC simulation has a stable southern overturning cell within 500 model years, its equivalent in the 4 × PIC case only materialises after 2750 years (Fig. S3b).

2.4 Proxy compilations

The 38 Ma CESM simulations presented here are validated using both qualitative and quantitative measures based on proxy records of the middle to late Eocene. After a general assessment of the model results we look more specifically at sea surface temperatures (SSTs); those from the 38 Ma 4 × PIC simulation are assumed to best represent the 42–38 Ma proxy records and from the 2 × PIC case for the 38–34 Ma interval. An overview of the considered SST proxies is given in Table S3 for 42–38 Ma and Table S4 for 38–34 Ma using data from Pearson et al. (2001, 2007), Tripati et al. (2003), Kobashi et al. (2004), Bijl et al. (2009), Liu et al. (2009), Okafor et al. (2009), Douglas et al. (2014), Hines et al. (2017), Cramwinckel et al. (2018), and Evans et al. (2018). A point-by-point comparison with proxy records is made using their estimated 38 Ma positions in accordance with the PaleoMag geography from Baatsen et al. (2016). Proxy-derived values for annual mean temperatures are considered with their calibration uncertainty for different methods (UK37, TEX86H, Mg  Ca, Δ47, and δ18O). Modelled SSTs at the reconstructed proxy locations are given with error bars representing the variance in a surrounding 4× 4 box, covering the uncertainty associated with the palaeogeographic reconstruction (, last access: 29 November 2020). Estimates using UK37 are disregarded at low latitudes as the calibration saturates at ∼28C (Conte et al.2006). There is an ongoing discussion on which calibration to use for TEX86 (O'Brien et al.2017; Cramwinckel et al.2018; Hollis et al.2019). In addition to SST estimates using the exponential TEX86H calibration from Kim et al. (2010), we show those using the linear relation of Kim et al. (2008) (light blue markers in Fig. 6). There is also uncertainty in the calibration of Mg  Ca, for which we use the constraints on seawater chemistry from Evans et al. (2018).

In addition to SSTs, we compare modelled near-surface air temperatures to terrestrial proxies. Similarly, the distinction is made between the late middle Eocene (∼42–38 Ma) and late Eocene (∼38–34 Ma) as a reference for the 38 Ma 4 × PIC and 2 × PIC simulations, respectively. The considered records are listed in Tables S5 and S6, containing data from Greenwood and Wing (1995), Gregory-Wodzicki (1997), Smith et al. (1998), Wolfe et al. (1998), Greenwood et al. (2004), Retallack et al. (2004), Hinojosa and Villagrán (2005), Uhl et al. (2007), Boyle et al. (2008), Schouten et al. (2008), Prothero (2008), Eldrett et al. (2009), Quan et al. (2012), and Passchier et al. (2013). Terrestrial proxies used here consist of vegetation-based indicators using pollen, nearest living relative (NLR), leaf physiognomy (LMA, CLAMP, ELPA, and LMA; Wing and Greenwood1993; Yang et al.2011; Traiser et al.2005; Kowalski and Dilcher2003), and MBT–CBT (Peterse et al.2012).

2.5 Previous model results

The results of our 38 Ma simulations are compared to those of Goldner et al. (2014) (hereafter: GH14) and Hutchinson et al. (2018) (H18). The latter provides an opportunity to test the robustness of the overall temperature distribution in the modelled middle to late Eocene cases using very similar boundary conditions but a different model. Here, we are mainly looking at zonally averaged sea surface and near-surface air temperatures and how they correspond to the previously presented proxies. A side-by-side overview of oceanic and atmospheric fields from this study and GH14 can be found in the Supplement.

The simulations of GH14 were carried out with the CCSM3, a predecessor of the CESM1 (the former having CAM3 rather than CAM4 with a spectral versus finite-volume core, respectively). GH14 used a 45 Ma hot-spot-referenced reconstruction for their geographical boundary conditions, in which Antarctica is shifted by ∼6 latitude compared to the 38 Ma PaleoMag reconstruction used here. They applied four subsequent doublings of pre-industrial CO2 levels, which correspond to the EO1–EO4 cases of Cramwinckel et al. (2018). In terms of horizontal resolution GH14 implemented ∼2.5 and 3.7 for the model's ocean and atmosphere grids, respectively, versus ∼1 and ∼2 here.

A model set-up similar to ours was implemented by H18 using the 38 Ma geography reconstruction from Baatsen et al. (2016), a similar oceanic resolution (∼1 – atmosphere: ∼3), and a vegetation cover based on the same biome distribution used here (Fig. S1a). They studied mostly the late Eocene using the GFDL CM2.1 and subsequent doublings of present-day rather than pre-industrial CO2 concentrations (i.e. ∼400 ppm).

A more quantitative assessment of the model's performances is made following a procedure similar to the one presented by Lunt et al. (2012) using the SST proxy compilations considered here (Tables S3 and S4). The model–proxy discrepancy σ is defined as the average difference between model-predicted and proxy-induced SSTs at reconstructed locations (using 38 Ma PaleoMag for this study and H18 versus 45 Ma hot spot for GH14). A set of different comparisons is considered using either all proxies or a subset (i.e. global, equatorial, or extra-tropical) and modelled annual mean versus summertime temperatures. An absolute error |σ| is also introduced (rather than e.g. RMS) to rule out better scores by compensating for errors without exaggerating the impact of a single (possibly unrealistic or not representative) proxy value.

2.6 Climate sensitivity

A number of additional simulations are carried out to determine the model's response to an altered radiative forcing when only the concentration of atmospheric greenhouse gases is changed. These include a 2000-year 4 × CO2 (Fig. S16) and four 20-year perturbation experiments listed below. Using the extrapolation method introduced by Gregory et al. (2004) on the results of those shorter simulations (see Fig. S17 for 4 × CO2), the model-derived radiative forcing from consecutive doublings of either CO2 or both CO2 and CH4 is determined as (bracketed values: theoretical estimates from Etminan et al.2016; E16)

  • 1 × CO2 2 × CO2: ΔRF2×CO2= 3.49 W m−2
    (E16: 3.80 W m−2),

  • 1 × CO2 4 × CO2: ΔRF4×CO2= 7.93 W m−2
    (E16: 7.96 W m−2),

  • 1 × PIC  2 × PIC: ΔRF2×= 4.18 W m−2
    (E16: 4.23 W m−2),

  • 1 × PIC  4 × PIC: ΔRF4×= 9.33 W m−2
    (E16: 8.96 W m−2).

The value of ΔRF2×CO2 determined here is slightly lower than the 3.8 W m−2 estimated by E16 but very close to the 3.5 W m−2 reported by Kay et al. (2012) using offline radiative transfer calculations for CAM4 specifically. Note that the model-derived radiative forcing is generally close to the corresponding one from E16 but featuring a stronger nonlinear behaviour towards higher perturbations, especially when including the effect of CH4. Using the values shown above, the radiative forcing that results only from a second PIC doubling can thus be estimated as ΔRF2×4×= 5.15 W m−2 (RF4×-ΔRF). Through a (log-scaled) quadratic fit, these estimates can also be used to show that our 2 × PIC is equivalent to 2.25 × CO2, while 4 × PIC corresponds to 4.85 × CO2 (rather than 2.15 and 4.69 based on E16, respectively).

Starting from the equilibrated pre-industrial reference (at model year ∼3100), an instant quadrupling of atmospheric CO2 was applied and continued for another 2000 model years. Using the same procedure as Gregory et al. (2004), the equilibrium climate sensitivity (ECS) of the model was assessed at 3.17 C per CO2 doubling. The first 100 model years are disregarded in the extrapolation to better capture the effect of slow feedbacks (see Fig. S16). Using the same atmospheric configuration and a slab ocean model, Bitz et al. (2012) reported a comparable ECS of 3.1 C in the CCSM4. With the model-derived radiative forcing of 7.93 W m−2 and an extrapolated equilibrium response of 6.34 C, the normalised equilibrium climate sensitivity of the pre-industrial reference is therefore SPI= 0.80 C W−1 m2.

Since we carry out a long spin-up simulation for each of the 38 Ma cases, we can determine the equilibrium climate response as the difference between the actual equilibrium states rather than making an estimate through extrapolation. The mean over the last 200 (rather than 50 for the other results) model years of each simulation is taken to exclude as much of the internal variability as possible. When comparing global average temperature differences between the pre-industrial reference and both Eocene simulations, not all of the warming is a result of higher greenhouse gas concentrations. To separate the internal effect of palaeogeography from the externally (i.e. greenhouse gas) driven warming, the 38 Ma Eocene equilibrium climate sensitivity SEO can also be calculated using the combined radiative forcing due to greenhouse gas and palaeogeography changes by assuming that

(1) S EO = Δ T Δ RF + G ,

where ΔT is the temperature difference between two climatic states, ΔRF the net radiative forcing from greenhouse gases, and G the radiative forcing due to the (integral) effect of palaeogeography. Note that the latter also includes changes in the distribution of e.g. land ice and vegetation, as those are applied boundary conditions in the model set-up used here. In order to be compatible with SPI, this formulation returns S in normalised units of [C W−1 m2] rather than [C per CO2 doubling] as used by Royer et al. (2012).

Since both 38 Ma Eocene simulations use the same boundary conditions (except CO2 and CH4 concentrations), it is reasonable to assume that S and G estimated from the comparison with the pre-industrial case should be similar (provided that their nonlinear contribution is small). By comparing the temperatures of both Eocene runs to those of the pre-industrial reference, G can then be estimated from

(2) G = Δ T 2 × Δ RF 4 × - Δ T 4 × Δ RF 2 × Δ T 4 × - Δ T 2 × ,

where ΔT and ΔT denote the temperature difference with respect to the pre-industrial climate for the 38 Ma 4 × PIC and 2 × PIC case, respectively. Within this formulation, globally averaged near-surface air temperatures can then be substituted by equatorial (<23.5 N or S) SSTs (with a 3∕2 ratio between global and equatorial change suggested by Royer et al.2012) or deep-sea temperatures, which should both be more compatible with ECS estimates from proxy data.

3 Results

3.1 The simulated middle to late Eocene equilibrium climate

Taking the mean over the last 50 years of the 38 Ma 4 × PIC simulation, a number of annual mean oceanic fields are shown in Fig. 3 (see also Fig. S4 for the pre-industrial reference and Fig. S5 for the 38 Ma 2 × PIC case; seasonal fields are shown in Fig. S8). Sea surface temperatures are quite warm at low latitudes (34 C on average, regionally >36C) but also mild across high latitudes (Fig. 3a). The Southern Ocean, especially its Pacific part, is mainly characterised by annual mean temperatures of 10–20 C, while those of the Arctic Ocean are mostly 6–10 C. The former Neo-Tethys (Pacific-Indian-Mediterranean) Ocean still mostly acts as one basin in the equatorial region, with an expansive warm pool and cold tongue across its western and eastern part, respectively.

Figure 3Annual mean for the 38 Ma 4 × PIC simulation with (a) sea surface temperature and (b) salinity, (c) barotropic stream function (positive for clockwise flow) and zonal wind stress (contours every 4 × 10−2 Pa, solid positive and dashed negative, thick line at 0 Pa), and (d) upper 200 m average potential density and mixed layer depth (contours at 100 and 250 m, thick white line at 500 m).

Mild temperatures in the sub-polar South Pacific are accompanied by relatively high salinities of around 35 psu (Fig. 3b). Higher salinities are present across low-latitude evaporative regions around the world, in contrast to the much fresher Arctic Ocean (∼20 psu). The latter is geographically almost isolated and only connected to the shallow Para-Tethys by the Turgai Strait (<100 m depth) and through the Nordic Seas (<1 km depth) to the North Atlantic Ocean where consequently low surface salinities are seen as well.

The surface temperature and salinity patterns are reflected by upper 200 m potential density (Fig. 3d), with high values throughout the Southern Ocean and much lower densities and thus more stably stratified waters across the North Atlantic and Arctic Ocean. Deepwater formation occurs only in the Pacific sector of the Southern Ocean in winter (at ∼11.4C and 34.8 psu), but high upper-level densities suggest that this may take place virtually anywhere around the Antarctic continental slope. The only part of the Northern Hemisphere where deepwater formation could occur is the high-latitude Pacific Ocean, but it is never seen in either of the 38 Ma simulations.

The global ocean circulation consists mostly of expansive sub-tropical gyres, a geographically restricted Antarctic Circumpolar Current (ACC), and some rather weak sub-polar gyres (Fig. 3c). This proto-ACC is associated with a sharp frontal zone separating warm sub-tropical from cooler sub-polar waters. Even with a shallow Drake Passage and Tasmanian Gateway (200–500 m; see Fig. 1a, d), the integrated zonal flow is about 25 %–30 % of its pre-industrial equivalent (45–50 Sv vs. 180 Sv; see also Fig. S3a). Strikingly, the temperature front is located at 55–60 S, which is almost 10 poleward of where it is found in the pre-industrial reference (Fig. S4). The zonal variability in the location of the front also has profound implications for regional meridional temperature contrasts. This location is strongly bound by the zonal wind stress maximum, indicating that it is fixed by the latitudes where both atmospheric and oceanic flow are the least obstructed by continents.

Despite differences in spin-up time, patterns of the equilibrium ocean circulation state are generally similar for our 38 Ma 2 × PIC and 4 × PIC simulations (see also Fig. S5). The annual mean global average sea surface temperature (SST) is 28.4 C in the 4 × PIC case versus 25.7 C in the 2 × PIC one and 18.4 C in the pre-industrial reference (Table ). Similar temperature differences of 2.5–3 C are seen in both the upper and deep ocean between the 38 Ma cases (Fig. 2a). Globally averaged deep-sea temperatures (below 2 km) reach ∼11.5C in the 38 Ma 4 × PIC case and ∼8.6C in the 2 × PIC one, which is much warmer than for the pre-industrial ocean (∼0.7C). The difference in equatorial (<23.5 N or S) SSTs between the Eocene cases and the pre-industrial reference is smaller (i.e. 2–3 C less than the other measures listed) but still considerable.

Table 3Mean equilibrium temperatures over the last 200 model years of each simulation, showing MATglob: global average air temperature (at 2 m reference height), SSTglob: global average sea surface temperature, SSTeq: equatorial (<23.5 N or S) average SST, and Tdeep: global average deep-ocean temperature (below 2 km).

Download Print Version | Download XLSX

Also for the atmospheric fields, the mean over the last 50 years of the 38 Ma 4 × PIC simulation is taken and shown in Fig. 4 (see also Fig. S6 for the pre-industrial reference and Fig. S7 for the 38 Ma 2 × PIC case; seasonal fields of 4 × PIC are shown in Fig. S9). As for the ocean, the modelled 38 Ma 2 × PIC atmosphere (Fig. S7) is quite similar to the 4 × PIC one. At 23.34 C the globally averaged annual mean near-surface air temperature at 2 × PIC is ∼3.2C cooler than at 4 × PIC (26.55 C). This is, however, still ∼9.5C warmer than in the pre-industrial reference (13.82 C).

Figure 4Annual mean for the 38 Ma 4 × PIC simulation with (a) near-surface (2 m) air temperature (shading) and min–max temperature (contours; blue for Tmin<0C and white for Tmax>40C), (b) precipitation (shading) and mean sea level pressure (MSLP; contours every 5 hPa, thick black lines every 20 hPa, ≥1000 hPa in black and <1000 hPa in white), (c) 700 hPa temperature (shading) and 200 hPa wind speed (contours every 5 m s−1 starting at 10 m s−1; white for <20 m s−1 and black for ≥20 m s−1, thick black lines every 20 m s−1), and (d) 500 hPa geopotential height (shading) and tropopause height (contours every 1 km; white for <15 km and black for ≥15 km, thick black line at 15 km).

The simulated Eocene climate is overall warmer and wetter than that of the pre-industrial reference, with reduced contrasts between low and high latitudes. Near-surface air temperatures (at 2 m; Fig. 4a) under 4 × PIC are above zero and up to 40 C over land at low latitudes in the annual mean (with summer maxima regionally >50C; see Fig. S9a, b). Annual mean (daily) minimum temperatures are only below freezing in East Antarctica and northeast Siberia (the latter only at 2 × PIC). Extreme seasonality is seen over the continental interior of Antarctica (Fig. S9a, b), with 35–45 C differences between daily mean summer and winter temperatures.

Simulated precipitation and mean sea level pressure (MSLP) patterns (Fig. 4b) indicate the presence of a prominent tropical trough and expansive sub-tropical ridges in the Eocene. The inter-tropical convergence zone (ITCZ) consists of two precipitation maxima extending over the Indo-Pacific basin on both sides of the Equator. A pronounced double ITCZ over the Pacific Ocean is also seen in the pre-industrial reference (Fig. S6b) and a known model-related issue (Song and Zhang2009; Bellucci et al.2010). Effects of orographic lift are evident on westward-facing coastlines and mountain ranges across middle and high latitudes (e.g. southern Andes and northern Rocky Mountains). Regardless of an overall reduced Equator-to-pole temperature gradient in the 38 Ma 4 × PIC case, middle-latitude storm tracks are prominent with increased precipitation and a poleward expansion with respect to the pre-industrial reference. Strong seasonality in the precipitation patterns also indicates the importance of monsoons in this warm Eocene climate (Fig. S9c, d).

The modelled 700 hPa temperature (Fig. 4c) of the 38 Ma 4 × PIC case highlights warm mid-level air masses in persistent continental high-pressure regions in the sub-tropics. The air over Antarctica is substantially warmer than over the Arctic because of its elevation and continental climate. Like the Tibetan Plateau today, the Antarctic continent therefore acts as an elevated heat island in summer (Hoskins and Karoly1981; Ye and Wu1998; Fig. S9f). Wind speeds at 200 hPa highlight the positions of both sub-tropical and polar jet streams. Several wind maxima, related to (mostly topographically induced) preferred Rossby wave activity, are evident but most pronounced at ∼90 E. A comparable jet stream pattern with regional wind maxima is seen in the pre-industrial reference, but it is more zonally uniform and shifted equatorward.

Reduced latitudinal and hemispheric differences in the modelled 38 Ma 4 × PIC compared to the pre-industrial reference are also seen in 500 hPa geopotential and tropopause heights (Fig. 4d). The latter is defined using the WMO definition of a −2C km−1 lapse rate. Both fields show a nearly meridionally symmetric pattern for the Eocene, with the sharpest gradient across middle-latitude regions. The overall warmer air column in the Eocene case accounts for an increase in the 500 hPa surface by 100–200 m (up to 400 m over Antarctica) and the tropopause by 1–2 km.

We assess the global heat budget in our model results by looking at both the oceanic and atmospheric heat transports, along with their link to global circulation patterns, presented in Fig. 5. The globally integrated oceanic meridional overturning stream function is similar in pattern and extent (Fig. 5a) for both 38 Ma Eocene simulations but different to that of the pre-industrial reference (Fig. 5b). Little change in shallow wind-driven cells around the Equator is seen among all three cases. A deep meridional overturning cell is present in the Southern Hemisphere rather than the Northern Hemisphere in the 38 Ma configuration, which is slightly stronger and deeper in the 2 × PIC case. The effect of upwelling due to (Ekman pumping by) the zonal flow in the (proto-)ACC can be seen at around 55 S in both 38 Ma cases as well as in the pre-industrial reference. The related overturning cell in the latter is much deeper and stronger as a result of a more developed ACC but mostly the unrestricted zonal flow in the present geographical configuration. Note that the average latitudinal shift of these overturning cells between the modelled Eocene and pre-industrial circulation is less pronounced than the more regional variations of the surface polar front in Figs. 3 and S4.

Figure 5Global oceanic meridional overturning stream function averaged over the last 50 model years of the (a) 38 Ma 4 × PIC (shading) and 2 × PIC (contours) simulations as well as the (b) pre-industrial reference (shading; 38 Ma 2 × PIC again in contours). Contours are drawn using the same colour scale at levels indicated by white lines in the colour bar to the right. (c) Integrated meridional heat transport in the ocean, globally (greyscale) and Indo-Pacific only (coloured) for the three model cases. Note that the horizontal (latitude) scale covers the range [75 S, 75 N], with latitude increasing from left to right, for a better comparison with oceanic fields in (a) and (b). The top-of-model (TOM) required total meridional heat transport is shown in (c) for 38 Ma 4 × PIC (solid black) and the pre-industrial reference (dotted grey) with the corresponding atmospheric component in red (4 × PIC) and dotted blue (pre-industrial). Differences with respect to 38 Ma 4 × PIC for both total (dashed dark grey) and atmospheric (dashed orange) heat transport are given for the 2 × PIC case, magnified 10-fold.


An otherwise predominantly wind-driven gyre circulation in the 38 Ma simulations is reflected by symmetric oceanic meridional heat transport into both hemispheres (Fig. 5c). This is in contrast to the pre-industrial situation, in which a difference of up to 1 PW is seen between hemispheres, making the Northern Hemisphere relatively warm (Trenberth and Caron2001). A major part of the 38 Ma oceanic heat transport occurs in the Neo-Tethys–Indo-Pacific Ocean, indicating the dominant role of this basin within the Eocene circulation. The presence of a deep overturning cell in the North Atlantic Ocean in the pre-industrial reference substantially reduces the relative Indo-Pacific contribution. The largest relative differences between hemispheres are seen at high latitudes, where about 0.5 PW is transported southward at 45 S, while the heat transport is close to zero at 45 N in both 38 Ma cases. This difference can be explained by the presence of a deep southern overturning cell, pulling warm waters into the southern high latitudes.

To determine the total required meridional heat transport, the top-of-model net radiative flux is integrated zonally (Fig. 5d). Overall differences are small between the 38 Ma cases and pre-industrial reference, with a slight shift towards the Southern Hemisphere, making the Eocene pattern more symmetric. Changes seen in oceanic meridional heat transport are less pronounced in the total transport and are thus compensated for by the atmospheric component, albeit only partially. This suggests that the reduced latitudinal temperature gradient seen in the modelled 38 Ma cases is thus not sustained by increased meridional heat transport, but rather induces a reduced transport with respect to the pre-industrial reference.

3.2 Model–proxy comparisons

In terms of the general ocean circulation, proxy-based information is limited but mostly agrees well with what is shown by the model. Persistent deepwater formation in the South Pacific Ocean throughout most of the Eocene and the possibility of a North Pacific or South Atlantic source are consistent with the findings of Cramer et al. (2009) and Hague et al. (2012). A strongly stratified North Atlantic Ocean with low-salinity waters at middle to high latitudes was also suggested by Coxall et al. (2018) prior to ∼36 Ma.

The all but complete absence of Arctic sea ice in the 38 Ma 4 × PIC model results differs from indications of seasonal sea ice during the middle Eocene by Stickley et al. (2009). Although limited in extent, the model does show sea ice during the winter months under 2 × PIC (contours in Fig. S5). In addition to the sea ice indications being highly localised, age restrictions on the considered section of the ACEX core make it difficult to rule out whether the sea ice actually occurred in the late Eocene. Still, the model's limitations to represent sea ice need to be considered here, as the applied scheme is quite simple. Low salinities across the Arctic Ocean would have made it easier for sea ice to form, which is not reflected by the implemented fixed SST threshold (−1.8C).

A southward migration of the temperature front in the Southern Ocean can help to explain some of the high SSTs induced from proxies, especially across the southwest Pacific region (Bijl et al.2009; Hines et al.2017). Additionally, zonal variations in the front induce a thermal heterogeneity in the Southern Ocean. The 6–8 C (Fig. 3a) difference between Tasmania and the tip of the Antarctic Peninsula (Seymour Island) agrees well with proxy indications from Douglas et al. (2014). A 10–15 latitude difference between Tasmania and the tip of the Antarctic Peninsula in the (PaleoMag) geography reconstruction used (Baatsen et al.2016) here may help to explain the temperature differences seen in proxies. On the other hand, eastward flow through the Tasmanian Gateway limits the strength and extent of the Antarctic countercurrent around most of Antarctica (Bijl et al.2013; Houben et al.2019). The model results disagree with indications of a strong sub-polar gyre in the Ross Sea (limited in the model by the northern extent of East Antarctica) into the late Eocene, with predominantly westward flow through the Tasmanian Gateway and Antarctic-derived surface waters found in the southwest Pacific Ocean (Stickley et al.2004; Huber et al.2004; Bijl et al.2013; Cramwinckel et al.2020). The absence of a clear temperature discrepancy between the two sides of the gateway still suggests a strong connection as seen in the model (Bijl et al.2013), possibly explained by atmospheric influences at shallow and near-coastal locations.

In Fig. 6a point-by-point comparison with available SST proxies (see also Table S3) is complemented by zonally averaged SSTs from the model along with their longitudinal spread for both the annual and summer mean (for 38 Ma 4 × PIC; see also Fig. S10 for 2 × PIC).

Figure 6Zonal average annual mean (black) and summer mean (red) sea surface temperature (SST), with shaded regions showing minimum and maximum values for each latitude, from the 38 Ma 4 × PIC simulation. Blue markers indicate estimates from 42–38 Ma SST proxies (see also Table S3; light blue: linear, dark: log calibration from Kim et al.2008, for TEX86H), whereas the white (annual) and red (summer) circles depict model values at the corresponding 38 Ma locations. Error bars are obtained using proxy calibration errors and the spatial variation within a 4× 4 box surrounding the corresponding location in the model. The inset shows a scatter plot comparing proxy and model SST, with dashed lines indicating a linear fit using annual mean (black) or summer mean (red; only poleward of 23 N or S) model temperatures.


For low-latitude regions, model results and proxies are generally in good agreement (i.e. within uncertainty). Some near-equatorial sites (e.g. Atlantic ODP Sites 925, 929, and 959) show cooler temperatures than suggested by the model of around 30 C (when considering TEX86H). Both the model's limitations representing the sharp upwelling region and uncertainty in the palaeolatitude can explain this offset. Cooler near-equatorial temperatures (down to <30C) are in fact indicated by the increased zonal spread, showing that these are possible in the simulated Eocene climate. TEX86H estimates from Tanzania and recent clumped isotope measurements from Java suggest that SSTs of ∼35C or more are possible outside upwelling regions, which is well in line with the 38 Ma 4 × PIC results.

Mixed agreement is seen at other latitudes, with model results seemingly too warm at northern middle latitudes, too cold at southern middle latitudes, and good at high latitudes. The large discrepancy with proxy estimates from the Gulf of Mexico and Neo-Tethys (Mediterranean) stands out. While it seems unlikely that these regions would be colder in the Eocene than they are today, the model probably underestimates strong seasonal cooling in shallow coastal waters and other localised effects near the coastline due to the complex palaeogeography at both locations. SST estimates from other locations at a similar latitude, such as ODP Site 1052 on Blake Nose (east of Florida), are much higher, bringing them back into good agreement with the model.

An improved match is generally found at most southern middle and high southern latitudes when modelled summertime temperatures are considered. Since SST proxies are based on past living organisms, they possibly have a bias towards the warm season as their activity and sedimentation are directly or indirectly affected by the available amount of sunlight (Sluijs et al.2006, 2008; Bijl et al.2009; Hollis et al.2012; Schouten et al.2013). Confidence in the modelled annual mean SSTs is boosted by nearly perfect agreement in the Arctic and clumped isotope indicators elsewhere (UK and Seymour Island).

The modelled 38 Ma 2 × PIC SSTs similarly show overall good agreement with 38–34 Ma proxy estimates, with a higher spread in the latter (Fig. S10 and Table S4). The best overlap is again seen in low-latitude regions, as is a better match with modelled summer temperatures at high latitudes (except for clumped isotopes). Considering the model's limitations and uncertainties or possible biases in proxy-derived SST estimates, these 38 Ma simulations can reconstruct the middle to late Eocene (42–34 Ma) temperature distribution well.

Despite the limited latitudinal coverage of terrestrial proxies, overall good agreement is also seen between the modelled temperatures on land (Fig. S11 for 38 Ma 4 × PIC; see also Fig. S12 for 2 × PIC). Proxy estimates from China, however, indicate little change in temperature, while the model shows considerably warmer conditions at lower latitudes. Quan et al. (2012) discuss the limitations of capturing high summer temperatures, as it is unlikely to see similar conditions over a >20 latitude range. There are currently no proxy indications to assess the possibility of very high model-based temperatures (exceeding 50 C; see Fig. S9a, b) in low-latitude continental interiors. An important limitation of the model, especially on land, is its resolution, which smoothens the topography and therefore underestimates local temperature effects. A correction using the original 0.1 geography reconstruction from Baatsen et al. (2016) and a free tropospheric lapse rate of −6.5 K km−1 significantly improves the agreement between modelled and proxy-based temperatures, especially in North America (small squares versus filled markers in the inset of Fig. S11).

3.3 Model–model comparisons

Zonally averaged SSTs of our 38 Ma 4 × PIC and 2 × PIC simulations are compared to those of the 4 × CO2 (i.e. 1120 ppm) case from GH14 and the 800 ppm CO2 (i.e. 2× PD) case from H18 in Fig. 7, along with the pre-industrial reference and the available middle to late Eocene proxies. A similar overview considering near-surface air temperatures and terrestrial proxies can be found in Fig. 8. All of the Eocene cases feature an overall warming compared to the pre-industrial reference, which is the strongest over southern high latitudes, along with a reduction of the Equator-to-pole temperature gradient. Both of our 38 Ma cases and the 800 ppm simulation of H18 are considerably warmer than the 45 Ma 4 × CO2 of GH14, indicating a higher sensitivity to the applied Eocene boundary conditions. It is important to mention that the simulations of H18 show deepwater formation in the North Pacific Ocean, in addition to the Southern Ocean as found here and in GH14. This results in a zonal variation of northern high-latitude SSTs, but little change is seen (outside the Arctic Ocean) in zonally averaged profiles in Fig. 7.

Figure 7Annual mean zonal average sea surface temperature for the 38 Ma 4 × PIC (black; this study), 2 × PIC (green; this study), 45 Ma 4 × CO2 (dashed red; GH14 – Goldner et al.2014), 38 Ma 2 × CO2 (i.e. 800 ppm, dash-dotted magenta; H18 – Hutchinson et al.2018), and pre-industrial reference (dotted blue; this study) simulations. Estimates from proxy data are represented by blue squares for 42–38 Ma and circles for the 38–34 Ma period. Shading indicates the zonal range in temperatures for both the 38 Ma 4 × PIC (grey) and 45 Ma 4 × CO2 (yellow) case; thin lines show the highest summer temperatures at each latitude (light grey shading added for 38 Ma 4 × PIC).


Figure 8Annual mean zonal average near-surface (2 m) air temperature for the 38 Ma 4 × PIC (black; this study), 2 × PIC (green; this study), 45 Ma 4 × CO2 (dashed red; GH14 – Goldner et al.2014), 38 Ma 2 × CO2 (i.e. 800 ppm, dash-dotted magenta; H18 – Hutchinson et al.2018), and pre-industrial reference (dotted blue; this study) simulations. Estimates from proxy data are represented by crosses for 42–38 Ma and circles for the 38–34 Ma period. Shading indicates the seasonal (DJF, JJA) range in temperatures for both the 38 Ma 4 × PIC (grey) and 45 Ma 4 × CO2 (yellow) case; thin lines show the highest summer temperatures at each latitude (light grey shading added for 38 Ma 4 × PIC).


The most prominent differences in geography between the 45 Ma GH14 and our 38 Ma simulations are in the representation of Antarctica, the Tasmanian Gateway opening, and the position of India. The Southern Ocean in particular shows large differences in circulation and the resulting temperatures, mainly related to the different model geographies. The formation of a proto-ACC and southern extent of the East Australia current act to shift the polar front in the South Pacific southward for the 38 Ma cases, while the opposite happens for 45 Ma (Fig. S13). The southward expansion of sub-tropical gyres and migration of the associated temperature front seen here are consistent with the findings of Viebahn et al. (2016) in response to a Drake Passage closure under present-day conditions. These changes can thus be directly linked to the continental configuration and associated shifts in zonal wind stress (maximum at 55 S versus 45 S in GH14). Generally, western boundary currents (e.g. Kuroshio, Agulhas, East Australia Current) and the effects of ocean bathymetry are more pronounced in the 38 Ma results. An issue in the 45 Ma results with very low (negative) salinities in the Arctic Ocean, although having seemingly little impact on the general circulation, is mostly resolved (lowest salinities down to ∼10 psu) in the 38 Ma case by having several shallow passages.

With several indications of near-equatorial temperatures as high as 34–36 C (Tanzania, Java, and Saint Stephens Quarry), the 38 Ma 4 × PIC case is able to match those proxies while still allowing equatorial upwelling zones to be <30C. The 800 ppm CO2 case of H18 shows similar tropical warmth but a steeper Equator-to-pole temperature gradient, especially in the Northern Hemisphere. Cooler low-latitude proxies of ∼30C are better matched by both the 45 Ma 4 × CO2 (GH14) and 38 Ma 2 × PIC results. Southern Hemisphere high-latitude proxies are difficult to meet by any model because of their large spread, with some of the higher estimates only matched by the warmest summertime temperatures in both the 38 Ma 4 × PIC and 800 ppm CO2 (H18) case. Still, most of the lower estimates are best reconstructed by the annual mean SSTs in the 38 Ma 4 × PIC case while also meeting the higher estimates when considering summer conditions.

Similar to the ocean results, all of the considered Eocene simulations show a comparable zonal average near-surface air temperature distribution that differs greatly from the pre-industrial reference (Fig. 8). Extreme summertime heat over low- to middle-latitude continental regions is seen in our 38 Ma 4 × PIC cases, but it is also present in the 38 Ma 2 × CO2 from H18 and to a lesser extent in the 45 Ma 4 × CO2 from GH14. Part of this warmth can be explained by the use of fixed vegetation types in all of the models, creating an efficient evaporative surface, while most plants would otherwise perish. The large seasonality over Antarctica is enhanced in our 38 Ma simulations compared to both GH14 and H18, with summer temperatures >30C in the 4 × PIC results (versus 20 C in GH14). This Antarctic summertime warmth, together with cloud–albedo feedbacks over the Arctic, keeps annual mean zonally averaged temperatures above 5 C in our 38 Ma 4 × PIC simulation across all latitudes.

Global average temperatures under similar radiative forcing are strongly model-dependent for the Eocene conditions studied here. While the 45 Ma 4 × CO2 (GH14) case is cooler than our 38 Ma 2 × PIC one, the H18 study simulates temperatures similar to those of our 4 × PIC case at 800 ppm CO2 (albeit with a steeper Equator-to-pole gradient). Although our 38 Ma simulations are the only ones with elevated CH4 concentrations, the difference in global temperature with GH14 is much larger than what would be expected solely from the resulting radiative forcing. The grid and resolution used here for the atmospheric component (CAM4) are both shown to slightly increase the sensitivity to a doubling of CO2 (Bitz et al.2012). Using the CAM4 finite-volume dynamical core instead of a spectral one greatly influences cloud radiative forcing and causes a warming compared to GH14 (there is a globally averaged net cloud forcing between the simulations of 13.8  6.4 = 7.4 W m−2; see also Table S2). Indeed, similar global temperatures are found under comparable radiative forcing by H18, who also use a finite-volume dynamical core.

Part of the temperature increase between pre-industrial and Eocene conditions can be explained by the global land fraction being reduced from 29.2 % to 26.2 % in our CESM simulations. A lower land fraction generally reduces the Earth's albedo and therefore induces a warming response as shown by Farnsworth et al. (2019). This cannot, however, explain the lower temperatures found by GH14, who implement a further reduced land fraction of 25.1 % in their model boundary conditions. A warm initialisation of the deep ocean, general circulation changes, and the applied aerosol distribution further add to temperature differences between models. It is worth noting here that recent model results from Zhu et al. (2019) also found a ∼5C warming in the CESM1.2 compared to their pre-industrial reference under the same external forcing, albeit for early Eocene (55 Ma) conditions.

A more detailed comparison of the different components in the radiative balance between the considered simulations can be found in Table S2. The Eocene boundary conditions implemented here mainly induce a warming through the radiative forcing from greenhouse gases and a reduced albedo. The latter effect is mostly compensated for by a negative cloud–albedo feedback in GH14, but not in our 38 Ma Eocene simulations. The warming is further enhanced by a strong positive water vapour feedback, explaining the much larger clear-sky longwave flux differences (i.e. 22.4 W m−2 globally) compared to what would be expected based on greenhouse gas concentrations alone.

An overview of model skill scores is given in Table 4, comparing modelled SSTs from the 45 Ma 4 × CO2 (GH14), 38 Ma 800 ppm (2×) CO2 (H18), and our 38 Ma 4 × PIC cases to 42–38 Ma proxies and those of our 38 Ma 2 × PIC to 38–34 Ma proxies. Although all of the considered model cases simulate on average cooler SSTs than those suggested by the proxies (σglob<0), considerable improvement is seen for the new 38 Ma cases (including H18). The Equator-to-pole temperature gradient is still being overestimated by the models, but the discrepancy (σeqσex) is smallest for the 38 Ma 4 × PIC case presented here. As expected, better model skill is seen when considering summertime temperatures but only if those are highly underestimated by the annual mean. This is thus not the case for our 38 Ma 4 × PIC and H18 800 ppm CO2 results, as the improvement at some sites is compensated for by others at which the proxies better agree with annual mean SST (see also Figs. 6 and 7). Those two model cases even underestimate the meridional gradient, taking summertime temperatures at extra-tropical locations (dashed red line in the inset of Fig. 6). This indicates that we cannot simply assume either an overall summer bias in those proxies or a cold bias in the models at extra-tropical sites. Indeed, further improvement of the model–proxy agreement is seen when considering modelled summer temperatures only at those sites where they provide a better match to the corresponding proxy estimate. For nearly all of the respective skill measures the 38 Ma 4 × PIC results show the best agreement with the available 42–38 Ma SST proxies, followed closely by those of H18 at 800 ppm CO2.

Table 4Model skill scores following a procedure comparable to the one presented by Lunt et al. (2012) for the model results of GH14, H18, and those presented here using the 42–38 Ma SST proxy compilation from Fig. 6 (see also Table S3). The last column considers the 38 Ma 2 × PIC case and 38–34 Ma SST proxies (see Fig. S10 and Table S4). Different measures are used (top to bottom): global, equatorial (<23.5 N or S), extra-tropical (>23.5 N or S), summer temperatures (only for >23.5 N or S), meridional temperature gradient (i.e. equatorial and extra-tropical; annual or summer), global absolute error, and summer absolute error. Bracketed values indicate model skill using modelled summer temperatures only when they improve agreement with the corresponding proxy estimate. An average skill score for each model–proxy comparison is also determined by taking the root mean square of all the above measures and highlighted in bold.

Download Print Version | Download XLSX

3.4 Climate sensitivity of the middle to late Eocene

The difference in the global average equilibrium temperature between our 38 Ma 2 × PIC and 4 × PIC simulations (3.21 C; see Table ) is similar to the model's response to a CO2 doubling in the pre-industrial reference (3.14 C). Yet, the estimated radiative forcing from a second PIC doubling (ΔRF2×4×=5.15 W m−2) results in a considerably lower equilibrium climate sensitivity for our Eocene simulations: SEO=0.62C W−1 m2. Even in this virtually ice-free world significant polar amplification is found in the warming response to a PIC doubling, which ranges between 2 C in low-latitude regions to as much as 8 C at high latitudes (∼0.4–1.6 C W−1 m2; see Fig. 9a–b). Zonal variation in the warming signal is highest at middle latitudes because of land–ocean contrasts, while differences in seasonal response are largest in polar regions. The polar amplification is reflected by differences in outgoing longwave radiation, being increasingly positive towards higher latitudes as a result of the Planck feedback (Fig. 9c–d). Overall net negative longwave fluxes over equatorial regions can be explained by tropopause responses and stratospheric cooling, which reduce outgoing radiation. On a more regional scale, alternating signals of differences in shortwave and longwave radiation occur with a pronounced minimum at the Equator and maxima over the ITCZ locations. This suggests the importance of cloud cover, which is seen to differ substantially (Fig. 9e–f) between the 38 Ma 2 × PIC and 4 × PIC cases. An increase in deep convective clouds over the equatorial Pacific corresponds to lower fluxes (more reflection of incoming shortwave and less outgoing longwave), while the opposite is seen off the Equator. Differences in cloud cover are seen to have a warming effect over most of the high-latitude regions. Increased wintertime low cloud cover over Antarctica and reduced summertime low cloud cover over the Arctic, along with increased high cloud cover, tend to warm the surface. The overall slight net heat loss (shortwave < longwave difference) at low latitudes versus a gain at middle to high latitudes agrees with a reduced meridional heat transport seen in the warmer 4 × PIC climate (Fig. 5c, d).

Figure 9(a) Annual mean temperature response, normalised by the global average TOM forcing (in W m−2) of the 38 Ma 4 × PIC compared to the 2 × PIC equilibrium climate with contours showing the winter season only (using the same colour scale). (b) Zonal average normalised temperature response; annual (black, shading for zonal variation), December–January–February (blue) and June–July–August (red). The (area-weighted) global average response (i.e. climate sensitivity) is 0.62 C W−1 m2 as indicated by the black dashed line. (c) Clear-sky component of the net longwave flux change at the TOM for 38 Ma 4 × PIC vs. 2 × PIC, similar for the clear-sky shortwave flux in contours (black to white; 1, 2, 5, and 10 W m−2). (d) Zonal average TOM flux change for longwave (black) and shortwave (red) fluxes, corresponding clear-sky components are shown using dashed lines. Note that longwave fluxes are defined positive upward while shortwave fluxes are positive downward. (e) Change in the total cloud forcing; longwave component (shading) and shortwave (contours every 5 W m−2; solid positive and dashed negative, thick black line at 0 W m−2). (f) Zonal average changes in cloud cover (in %) at low (solid red), medium (dashed green), high (dotted blue) and all (thick solid black) atmospheric levels.

As in the present-day climate the polar amplification is due in part to albedo effects from clouds, vegetation, snow cover, and some seasonal sea ice, but it is also the result of a strong water vapour feedback and meridionally dependent lapse-rate feedback (less negative towards higher latitudes). Low-latitude regions generally warm less than the global average in response to a positive radiative forcing, as they are strongly governed by moist processes (i.e. cloud cover and negative lapse-rate feedback) and more tied to SSTs (low land fraction). Additionally, the total column water vapour induces a negative shortwave feedback that is strongest in the tropics and helps to explain the smaller warming between our 38 Ma 2 × PIC and 4 × PIC than would be expected based on radiative forcing alone. An overview of zonally averaged atmospheric temperature changes between some of the different model cases considered here can be found in Fig. S15.

A list of averaged differences in our 38 Ma 4 × PIC and 2 × PIC simulations with respect to the pre-industrial reference is presented in Table 5, considering global near-surface air temperature (MATglob), equatorial SST (SSTeq), and deep-ocean temperature (Tdeep). Using Eq. (2), we can determine both the climate sensitivity S and integral geography effect G for the different measures. Not surprisingly, using MATglob yields the same value of SEO=0.62C W−1 m2, with G=11.07W m−2 (so G∼6.89C). The estimate of G decreases when using only oceanic temperatures as the direct effect of cooler temperatures over pre-industrial land ice is removed. While considering global MAT differences will likely overestimate the effect of global geography changes (mainly through lapse-rate effects on ice sheets and mountains), the opposite is true for equatorial SSTs. The corresponding results for SEO show less variation between different methods as they are tied to the respective temperature changes between the 38 Ma 2 × PIC and 4 × PIC cases. A possible range of G=6.08–11.10 W m−2 is found, equivalent to ∼1.6–2.6 CO2 doublings or a 4.2–6.9 C warming globally.

Table 5Temperature differences comparing the 38 Ma 4 × PIC (ΔT) and 2 × PIC (ΔT) Eocene climate to the pre-industrial reference as well as derived values for equilibrium climate sensitivity (S) and forcing from (integral) geography changes (G) using Eq. (2). Results are shown using global average 2 m temperature, equatorial SST (with a 3∕2 ratio, as discussed in Royer et al.2012), and deep-sea temperature (below 2 km).

Download Print Version | Download XLSX

A more direct estimate of G can be obtained by comparing a pre-industrial climate under 2 × PIC and 4 × PIC with the results from the modelled 38 Ma 2 × PIC case. Using the model's pre-industrial ECS (SPI=0.80C W−1 m2) and radiative forcing of two consecutive PIC doublings (ΔRF and ΔRF) yields an expected warming of 3.34 C and 7.46 C, respectively. As the modelled change in global average temperature with respect to the pre-industrial reference is either 9.50 C (2 × PIC) or 12.71 C (4 × PIC), an additional warming of 5.25–6.16 C takes place owing to integrated global geography changes (the effect of which is thus not entirely independent of the climatic state).

A similar exercise can be done for all of the different cases shown in Fig. 10 using SEO on one hand and SPI in combination with the previously defined values of RF and the available model results on the other hand. The projected temperatures suggest a contribution of ∼10 % at 2×CO2 PIC and ∼20 % at 4 × CO2 PIC from slow feedbacks (e.g. deep ocean, water vapour) to the global average warming response. The respective equilibrium temperatures for either pre-industrial or Eocene climates then yield a consistent estimate of G=8.6 W m−2, corresponding to a 5.3–6.9 C global average temperature difference.

Figure 10Overview of globally averaged near-surface air temperatures in the different pre-industrial reference (blue) and 38 Ma Eocene (red) simulations. Blue squares show the extrapolated temperature response of fast feedbacks from 20-year perturbations (4 × CO2, 2×CO2, 2 × PIC, and 4 × PIC), while grey circles also include that of slow feedbacks. Filled markers indicate results from equilibrated simulations, and open markers show values that were estimated from the model-derived climate sensitivity and radiative forcing (red axis to the right; note the change in scale between pre-industrial and Eocene using SPI and SEO, respectively). Offsets in both temperature and radiative forcing (i.e. G) between pre-industrial and 38 Ma Eocene cases are specified for all of the considered atmospheric compositions.


The possible range of G obtained here (∼6–11 W m−2; 4–7 C) is considerably higher than either the estimated 1.8 C warming related to geography (using HadCM3L) found by Lunt et al. (2012) or the ∼2–4 W m−2 (i.e. 0.6–1.1 CO2 doublings) suggested by Royer et al. (2012). A comparably high warming of ∼5C from integral geography effects in the Eocene was proposed by Caballero and Huber (2013) and was found more recently for the early Eocene by Zhu et al. (2019). Higher estimates of 2–6 C can also be deduced using the different models considered by Lunt et al. (2012).

4 Summary and conclusions

Using version 1.0.5 of the Community Earth System Model (CESM), we have presented the results of a simulated 38 Ma Eocene climate under both high (4 × PIC) and moderate (2 × PIC) concentrations of CO2 and CH4. These are among the first simulations with a fully coupled and detailed (CMIP5-like) climate model to specifically study the middle to late Eocene climate using a recent 38 Ma palaeogeography reconstruction.

The 38 Ma 4 × PIC case shows a warm climate with a global average near-surface air temperature of 26.5 C (pre-industrial: 13.8 C) and a low Equator-to-pole temperature gradient. The global heat budget is approximately meridionally symmetric, which is reflected by the zonal average temperature pattern. Deepwater formation occurs in the South Pacific Ocean, while the North Atlantic is stably stratified and stagnant due to the outflow of brackish Arctic waters. A shallow and rather weak precursor of an Antarctic Circumpolar Current is found at a more southerly latitude than its present equivalent, mainly driven by the position of zonal wind stresses. Continental low- to middle-latitude regions are characterised by high seasonality in both hemispheres and strong summer monsoons. Middle and high latitudes mostly have mild winters, warm summers, and pronounced storm tracks. The Arctic is rather cool due to its geographic isolation, and the Antarctic continental interior shows strong seasonality with especially hot summers.

Comparing our 38 Ma 4 × PIC model results to the available 42–38 Ma sea surface and terrestrial temperature proxy records shows an overall good match at middle to high latitudes without low latitudes being too hot. This indicates that the CESM is able to simulate the warm greenhouse climate of the late middle Eocene ( Bartonian) without the need for extremely high (>1200 ppm CO2) greenhouse gas forcing. In the simulated 38 Ma 2 × PIC climate, patterns of the oceanic and atmospheric circulation are qualitatively very similar to those of the 4 × PIC one. Based on a similar comparison between model results and 38–34 Ma proxy temperature estimates, our 38 Ma 2 × PIC case is a good analogue for the late Eocene climate ( Priabonian). While the model results are able to explain Southern Ocean heterogeneity and southwest Pacific warmth suggested by proxy records, they disagree with proxy-based indications of predominant westward flow through the Tasmanian Gateway during most of the Eocene.

Previous Eocene simulations (at 4× pre-industrial CO2) with a similar model but a different (45 Ma) continental configuration and lower resolution resulted in overall similar sea surface temperature distributions (Goldner et al.2014). However, our 38 Ma 4 × PIC case is about 4–5 C warmer globally in both SST and land temperature. Comparable results were found by Hutchinson et al. (2018) for the late Eocene using the same 38 Ma palaeogeography reconstruction of Baatsen et al. (2016) and the GFDL CM2.1. Higher resolution and more time-specific geographic boundary conditions allow for a better representation of regional features, including equatorial upwelling, zonal heterogeneity in the Southern Ocean, and Antarctic summer warmth. Warm Eocene temperatures in our simulations compared to most previous Eocene modelling studies are thus the combined result of a quadrupling of atmospheric methane concentrations (radiative forcing 4.85× pre-industrial CO2), the atmospheric dynamical core (finite-volume versus spectral), different cloud parameterisations, a strong water vapour feedback, and higher spatial resolution in combination with a newer palaeogeography reconstruction.

An equilibrium climate sensitivity of SEO=0.62C W−1 m2 is found between the 38 Ma 2 × PIC and 4 × PIC cases, which is lower than the same model's pre-industrial value (SPI=0.80C W−1 m2). While some of the altered model boundary conditions for the Eocene cases add to increased temperatures with respect to the pre-industrial reference, they also act to reduce climate sensitivity (e.g. land–sea distribution and cloud cover). Water vapour and lapse-rate feedbacks play crucial roles in both the reference state and sensitivity of the Eocene greenhouse climate, reducing the Equator-to-pole gradient but also the overall warming response to increased atmospheric greenhouse gas concentrations. The model-derived radiative forcing from a second PIC doubling (5.15 W m−2) is much higher than that of a first CO2 doubling (3.49 W m−2), starting from the pre-industrial reference. Even after 4600 model years the 38 Ma 4 × PIC simulation continues to warm up slightly, which suggests a possible underestimation of SEO.

When also taking the pre-industrial reference simulation into consideration, a fixed forcing G≈9 W m−2 from (integral) geography changes in the 38 Ma Eocene cases is estimated, corresponding to a ∼6C warming globally. Previous studies have noted this effect in terms of an offset in global average temperature between pre-industrial and (pre-EOT) palaeoclimate simulations (Lunt et al.2012; Farnsworth et al.2019). Similar to what was found by Caballero and Huber (2013) the direct effect of ice sheet coverage is limited, leaving a considerable warming due to other geography-related changes. When using oceanic instead of atmospheric model temperatures, the influence of topography and land surface changes between the Eocene and pre-industrial cases is indeed reduced, mostly by excluding the direct effect of the ice sheets and vegetation changes. Although smaller, the estimate of G obtained through oceanic temperatures is still larger (∼6 W m−2; 4 C) than suggested in most previous studies. This indicates a major contribution to G from changes in continental geometry, land surface properties, and the related circulation patterns, which our 38 Ma simulations resolve in more detail.

Several other notable (extreme) phenomena are found in our 38 Ma CESM simulations. Extremely high (∼50C) summer temperatures occur in the sub-tropics under 4 × PIC and are possibly related to fixed vegetation types, but this cannot be disproven by the available proxies. Strong seasonality is seen on Antarctica, where summer temperatures reach up to 35 C in the 4 × PIC case. The absence of an ice sheet, together with warm waters surrounding the continent and summertime insolation, causes the Antarctic continent to become a heat island. Sea ice coverage is nearly nonexistent and only occurs sporadically during the winter months in the 38 Ma 2 × PIC case.

As the simulated middle to late Eocene (38 Ma 2 × PIC and 4 × PIC) climate is in good agreement with estimates from currently available proxy records, the results presented here can be used to interpret (also novel) proxy records in more detail using the modelled circulation patterns. The results mainly show that a warm and equable Eocene climate can be simulated under realistic levels of atmospheric greenhouse gases. In line with another recent study from Hutchinson et al. (2018), our results can provide a basis to reconsider the conditions leading up to the EOT.

Data availability

The model data used to generate the main figures in this paper are publicly available on the Utrecht University Yoda platform. 38 Ma 4×PIC: (Baatsen2020a); 38 Ma 2×PIC: (Baatsen2020b); pre-industrial reference: (Baatsen2020c).


The supplement related to this article is available online at:

Author contributions

The model experiments were designed by ASvdH, HAD, and MB. Simulations were set up by MAK and MB, and data management and post-processing were performed by MAK. Analysis of the model results and generation of the figures and tables were performed by MB. The Introduction results from a joint effort of all authors, who also made various contributions to the discussion of results. SST proxy data were provided by AS and PKB, who joined MB in the processing and assessment of the model–data comparison. MH provided SAT proxy data, provided results of previously published CCSM3 Eocene simulations, and joined in the overall assessment of the model results.

Competing interests

The authors declare that they have no conflict of interest.


The authors would like to thank Nicholas Herold for his extensive help setting up the CESM and assisting to fix numerous bugs and difficulties emerging with deep-time palaeoclimate simulations. David Hutchinson is thanked for providing his late Eocene model results and his overall help through many discussions. Both anonymous reviewers and the editor are also thanked for their valuable input to improve the paper. The computations using CESM1.0.5 were done on the Cartesius at SURFsara in Amsterdam.

Financial support

This work was carried out under the programme of the Netherlands Earth System Science Centre (NESSC), financially supported by the Dutch Ministry of Education, Culture and Science. The use of the SURFsara computing facilities was sponsored by NWO-EW (Netherlands Organisation for Scientific Research, Exact Sciences) under project 15508. Anna S. von der Heydt acknowledges travel support to network partners from the EPSRC-funded Past Earth Network (grant no. EP/M008363/1), and Appy Sluijs thanks the European Research Council for a consolidator grant (grant no. 771497).

Review statement

This paper was edited by Yves Godderis and reviewed by two anonymous referees.


Abbot, D. S., Huber, M., Bousquet, G., and Walker, C. C.: High-CO2 cloud radiative forcing feedback over both land and ocean in a global climate model, Geophys. Res. Lett., 36, L05702,, 2009. a

Anagnostou, E., John, E. H., Edgar, K. M., Foster, G. L., Ridgwell, A., Inglis, G. N., Pancost, R. D., Lunt, D. J., and Pearson, P. N.: Changing atmospheric CO2 concentration was the primary driver of early Cenozoic climate, Nature, 533, 380–384,, 2016. a

Baatsen, M.: CESM data for Baatsen et al. 2020, 38 Ma 4×PIC, Utrecht University,, 2020a. a

Baatsen, M.: CESM data for Baatsen et al. 2020, 38 Ma 2×PIC, Utrecht University,, 2020b. a

Baatsen, M.: CESM data for Baatsen et al. 2020, Pre-Industrial Reference, Utrecht University,, 2020c. a

Baatsen, M., van Hinsbergen, D. J. J., von der Heydt, A. S., Dijkstra, H. A., Sluijs, A., Abels, H. A., and Bijl, P. K.: Reconstructing geographical boundary conditions for palaeoclimate modelling during the Cenozoic, Clim. Past, 12, 1635–1644,, 2016. a, b, c, d, e, f, g, h, i, j

Beerling, D., Berner, R. A., Mackenzie, F. T., Harfoot, M. B., and Pyle, J. A.: Methane and the CH4-related greenhouse effect over the past 400 million years, Am. J. Sci., 309, 97–113,, 2009. a

Beerling, D. J. and Royer, D. L.: Convergent Cenozoic CO2 history, Nat. Geosci., 4, 418–420,, 2011. a


Bellucci, A., Gualdi, S., and Navarra, A.: The double-ITCZ syndrome in coupled general circulation models: The role of large-scale vertical circulation regimes, J. Climate, 23, 1127–1145,, 2010. a

Bijl, P. K., Schouten, S., Sluijs, A., Reichart, G.-J., Zachos, J. C., and Brinkhuis, H.: Early Palaeogene temperature evolution of the southwest Pacific Ocean, Nature, 461, 776–779,, 2009. a, b, c, d, e

Bijl, P. K., Houben, A. J. P., Schouten, S., Bohaty, S. M., Sluijs, A., Reichart, G. J., Sinninghe Damsté, J. S., and Brinkhuis, H.: Transient Middle Eocene atmospheric CO2 and temperature variations, Science, 330, 819–821,, 2010. a


Bitz, C. M., Shell, K. M., Gent, P. R., Bailey, D. A., Danabasoglu, G., Armour, K. C., Holland, M. M., and Kiehl, J. T.: Climate sensitivity of the community climate system model, version 4, J. Climate, 25, 3053–3070,, 2012. a, b, c

Blackmon, M., Boville, B., Bryan, F., Dickinson, R., Gent, P., Kiehl, J., Moritz, R., Randall, D., Shukla, J., Solomon, S., Bonan, G., Doney, S., Fung, I., Hack, J., Hunke, E., Hurrell, J., Kutzbach, J., Meehl, J., Otto-Bliesner, B., Saravanan, R., Schneider, E. K., Sloan, L., Spall, M., Taylor, K., Tribbia, J., and Washington, W.: The Community Climate System Model, B. Am. Meteorol. Soc., 82, 2357–2376,<2357:TCCSM>2.3.CO;2, 2001. a

Bohaty, S. M. and Zachos, J. C.: Significant Southern Ocean warming event in the late middle Eocene, Geology, 31, 1017–1020,, 2003. a

Bosboom, R. E., Abels, H. A., Hoorn, C., van den Berg, B. C. J., Guo, Z., and Dupont-Nivet, G.: Aridification in continental Asia after the Middle Eocene Climatic Optimum (MECO), Earth Planet. Sc. Lett., 389, 34–42,, 2014. a

Boyle, B., Meyer, H., Enquist, B., and Salas, S.: Higher taxa as paleoecological and paleoclimatic indicators: A search for the modern analog of the Florissant fossil flora, in: Paleontology of the Upper Eocene Florissant Formation, Colorado, edited by: Meyer, H. W. and Smith, D. M., The Geological Society of America, Boulder, Colorado, 33–51, 2008. a

Caballero, R. and Huber, M.: State-dependent climate sensitivity in past warm climates and its implications for future climate projections, P. Natl. Acad. Sci. USA, 110, 14162–14167,, 2013. a, b

Carter, A., Riley, T. R., Hillenbrand, C. D., and Rittner, M.: Widespread Antarctic glaciation during the Late Eocene, Earth Planet. Sc. Lett., 458, 49–57,, 2017. a

Conte, M. H., Sicre, M. A., Rühlemann, C., Weber, J. C., Schulte, S., Schulz-Bull, D., and Blanz, T.: Global temperature calibration of the alkenone unsaturation index (UK′37) in surface waters and comparison with surface sediments, Geochem. Geophy. Geosy., 7, Q02005,, 2006. a

Coxall, H. K., Wilson, P. A., Pälike, H., Lear, C. H., and Backman, J.: Rapid stepwise onset of Antarctic glaciation and deeper calcite compensation in the Pacific Ocean, Nature, 433, 53–57,, 2005. a, b

Coxall, H. K., Huck, C. E., Huber, M., Lear, C. H., Legarda-Lisarri, A., O'Regan, M., Sliwinska, K. K., van de Flierdt, T., de Boer, A. M., Zachos, J. C., and Backman, J.: Export of nutrient rich Northern Component Water preceded early Oligocene Antarctic glaciation, Nat. Geosci., 11, 190–196,, 2018. a


Cramwinckel, M. J., Huber, M., Kocken, I. J., Agnini, C., Bijl, P. K., Bohaty, S. M., Frieling, J., Goldner, A., Hilgen, F. J., Kip, E. L., Peterse, F., Van Der Ploeg, R., Röhl, U., Schouten, S., and Sluijs, A.: Synchronous tropical and polar temperature evolution in the Eocene letter, Nature, 559, 382–386,, 2018. a, b, c, d, e, f, g

Cramwinckel, M. J., Woelders, L., Huurdeman, E. P., Peterse, F., Gallagher, S. J., Pross, J., Burgess, C. E., Reichart, G.-J., Sluijs, A., and Bijl, P. K.: Surface-circulation change in the southwest Pacific Ocean across the Middle Eocene Climatic Optimum: inferences from dinoflagellate cysts and biomarker paleothermometry, Clim. Past, 16, 1667–1689,, 2020. a

Danabasoglu, G., Ferrari, R., and McWilliams, J. C.: Sensitivity of an ocean general circulation model to a parameterization of near-surface eddy fluxes, J. Climate, 21, 1192–1208,, 2008. a

Danabasoglu, G., Bates, S. C., Briegleb, B. P., Jayne, S. R., Jochum, M., Large, W. G., Peacock, S., and Yeager, S. G.: The CCSM4 ocean component, J. Climate, 25, 1361–1389,, 2012. a

DeConto, R. M. and Pollard, D.: Rapid Cenozoic glaciation of Antarctica induced by declining atmospheric CO2, Nature, 421, 245–249,, 2003. a



Eldrett, J. S., Greenwood, D. R., Harding, I. C., and Huber, M.: Increased seasonality through the Eocene to Oligocene transition in northern high latitudes, Nature, 459, 969–973,, 2009. a

Elsworth, G., Galbraith, E., Halverson, G., and Yang, S.: Enhanced weathering and CO2 drawdown caused by latest Eocene strengthening of the Atlantic meridional overturning circulation, Nat. Geosci., 10, 213–216,, 2017. a

England, M. H.: The Age of Water and Ventilation Timescales in a Global Ocean Model, J. Phys. Oceanogr., 25, 2756–2777,<2756:TAOWAV>2.0.CO;2, 1995. a


Evans, D., Sagoo, N., Renema, W., Cotton, L. J., Müller, W., Todd, J. A., Saraswati, P. K., Stassen, P., Ziegler, M., Pearson, P. N., Valdes, P. J., and Affek, H. P.: Eocene greenhouse climate revealed by coupled clumped isotope-Mg  Ca thermometry, P. Natl. Acad. Sci. USA, 115, 1174–1179,, 2018. a, b, c, d

Farnsworth, A., Lunt, D. J., Brien, C. L. O., and Foster, G. L.: Climate Sensitivity on Geological Timescales Controlled by Nonlinear Feedbacks and Ocean Circulation, Geophys. Res. Lett., 46, 9880–9889,, 2019. a, b, c

Gasson, E., Lunt, D. J., DeConto, R., Goldner, A., Heinemann, M., Huber, M., LeGrande, A. N., Pollard, D., Sagoo, N., Siddall, M., Winguth, A., and Valdes, P. J.: Uncertainties in the modelled CO2 threshold for Antarctic glaciation, Clim. Past, 10, 451–466,, 2014. a, b




Greenwood, D. R. and Wing, S. L.: Eocene continental climates and latitudinal temperature gradients, Geology, 23, 1044–1048,<1044:ECCALT>2.3.CO;2, 1995. a

Greenwood, D. R., Wilf, P., Wing, S. L., and Christophel, D. C.: Paleotemperature Estimation Using Leaf-Margin Analysis: Is Australia Different?, PALAIOS, 19, 129–142,<0129:PEULAI>2.0.CO;2, 2004. a


Gregory-Wodzicki, K. M.: The Late Eocene House Range Flora, Sevier Desert, Utah: Paleoclimate and Paleoelevation, PALAIOS, 12, 552–567,, 1997. a

Hague, A. M., Thomas, D. J., Huber, M., Korty, R., Woodard, S. C., and Jones, L. B.: Convection of North Pacific deep water during the early Cenozoic, Geology, 40, 527–530,, 2012. a

Heavens, N. G., Shields, C. A., and Mahowald, N. M.: A paleogeographic approach to aerosol prescription in simulations of deep time climate, J. Adv. Model. Earth Sy., 4, M11002,, 2012. a



Hines, B. R., Hollis, C. J., Atkins, C. B., Baker, J. A., Morgans, H. E., and Strong, P. C.: Reduction of oceanic temperature gradients in the early Eocene Southwest Pacific Ocean, Palaeogeogr. Palaeocl., 475, 41–54,, 2017. a, b

Hinojosa, L. F. and Villagrán, C.: Did South American Mixed Paleofloras evolve under thermal equability or in the absence of an effective Andean barrier during the Cenozoic?, Palaeogeogr. Palaeocl., 217, 1–23,, 2005. a

Hollis, C. J., Taylor, K. W., Handley, L., Pancost, R. D., Huber, M., Creech, J. B., Hines, B. R., Crouch, E. M., Morgans, H. E., Crampton, J. S., Gibbs, S., Pearson, P. N., and Zachos, J. C.: Early Paleogene temperature history of the Southwest Pacific Ocean: Reconciling proxies and models, Earth Planet. Sc. Lett., 349–350, 53–66,, 2012. a, b

Hollis, C. J., Dunkley Jones, T., Anagnostou, E., Bijl, P. K., Cramwinckel, M. J., Cui, Y., Dickens, G. R., Edgar, K. M., Eley, Y., Evans, D., Foster, G. L., Frieling, J., Inglis, G. N., Kennedy, E. M., Kozdon, R., Lauretano, V., Lear, C. H., Littler, K., Lourens, L., Meckler, A. N., Naafs, B. D. A., Pälike, H., Pancost, R. D., Pearson, P. N., Röhl, U., Royer, D. L., Salzmann, U., Schubert, B. A., Seebeck, H., Sluijs, A., Speijer, R. P., Stassen, P., Tierney, J., Tripati, A., Wade, B., Westerhold, T., Witkowski, C., Zachos, J. C., Zhang, Y. G., Huber, M., and Lunt, D. J.: The DeepMIP contribution to PMIP4: methodologies for selection, compilation and analysis of latest Paleocene and early Eocene climate proxy data, incorporating version 0.1 of the DeepMIP database, Geosci. Model Dev., 12, 3149–3206,, 2019. a

Hoskins, B. J. and Karoly, D. J.: The Steady Linear Response of a Spherical Atmosphere to Thermal and Orographic Forcing, J. Atmos. Sci., 38, 1179–1196,<1179:TSLROA>2.0.CO;2, 1981. a

Houben, A. J., Bijl, P. K., Sluijs, A., Schouten, S., and Brinkhuis, H.: Late Eocene Southern Ocean Cooling and Invigoration of Circulation Preconditioned Antarctica for Full-Scale Glaciation, Geochem. Geophy. Geosy., 20, 2214–2234,, 2019. a

Huber, M. and Caballero, R.: The early Eocene equable climate problem revisited, Clim. Past, 7, 603–633,, 2011. a, b, c

Huber, M. and Sloan, L. C.: Heat transport, deep waters, and thermal gradients: Coupled simulation of an Eocene greenhouse climate, Geophys. Res. Lett., 28, 3481–3484,, 2001. a

Huber, M., Brinkhuis, H., Stickley, C. E., Döös, K., Sluijs, A., Warnaar, J., Schellenberg, S. A., and Williams, G. L.: Eocene circulation of the Southern Ocean: Was Antarctica kept warm by subtropical waters?, Paleoceanography, 19, PA4026,, 2004. a

Hunke, E. C. and Lipscomb, W. H.: The Los Alamos sea ice model, documentation and software, Version 3.1., Tech. Rep. Los Alamos National Laboratory, 2004. a


Hutchinson, D. K., de Boer, A. M., Coxall, H. K., Caballero, R., Nilsson, J., and Baatsen, M.: Climate sensitivity and meridional overturning circulation in the late Eocene using GFDL CM2.1, Clim. Past, 14, 789–810,, 2018. a, b, c, d, e, f, g, h, i

Inglis, G. N., Farnsworth, A., Lunt, D. J., Foster, G. L., Hollis, C. J., Pagani, M., Jardine, P. E., Pearson, P. N., Markwick, P., and Galsworthy, A. M. J.: Descent toward the Icehouse: Eocene sea surface cooling inferred from GDGT distributions, Paleoceanography, 30, 1000–1020,, 2015. a, b


Kay, J. E., Holland, M. M., Bitz, C. M., Blanchard-Wrigglesworth, E., Gettelman, A., Conley, A., and Bailey, D.: The influence of local feedbacks and northward heat transport on the equilibrium arctic climate response to increased greenhouse gas forcing, J. Climate, 25, 5433–5450,, 2012. a

Kennedy, A. T., Farnsworth, A., Lunt, D. J., Lear, C. H., and Markwick, P. J.: Atmospheric and oceanic impacts of Antarctic glaciation across the Eocene-Oligocene transition, Philos. T. R. Soc. A., 373, 2054,, 2015. a

Kennedy-Asser, A. T., Lunt, D. J., Farnsworth, A., and Valdes, P. J.: Assessing Mechanisms and Uncertainty in Modeled Climatic Change at the Eocene-Oligocene Transition, Paleoceanography and Paleoclimatology, 34, 16–34,, 2019. a

Kennedy-Asser, A. T., Lunt, D. J., Valdes, P. J., Ladant, J.-B., Frieling, J., and Lauretano, V.: Changes in the high-latitude Southern Hemisphere through the Eocene–Oligocene transition: a model–data comparison, Clim. Past, 16, 555–573,, 2020. a

Kiehl, J. T. and Shields, C. A.: Sensitivity of the Palaeocene-Eocene Thermal Maximum climate to cloud properties, Philos. T. R. Soc. A., 371, 2001,, 2013. a

Kiehl, J. T., Shields, C. A., Hack, J. J., and Collins, W. D.: The climate sensitivity of the Community Climate System Model version 3 (CCSM3), J. Climate, 19, 2584–2596,, 2006. a


Kim, J. H., van der Meer, J., Schouten, S., Helmke, P., Willmott, V., Sangiorgi, F., Koç, N., Hopmans, E. C., and Damsté, J. S.: New indices and calibrations derived from the distribution of crenarchaeal isoprenoid tetraether lipids: Implications for past sea surface temperature reconstructions, Geochim. Cosmochim. Ac., 74, 4639–4654,, 2010. a

Kobashi, T., Grossman, E. L., Dockery, D. T., and Ivany, L. C.: Water mass stability reconstructions from greenhouse (Eocene) to icehouse (Oligocene) for the northern Gulf Coast continental shelf (USA), Paleoceanography, 19, PA1022,, 2004. a

Kowalski, E. A. and Dilcher, D. L.: Warmer paleotemperatures for terrestrial ecosystems, P. Natl. Acad. Sci. USA, 100, 167–170,, 2003. a



Lawrence, D. M., Oleson, K. W., Flanner, M. G., Thornton, P. E., Swenson, S. C., Lawrence, P. J., Zeng, X., Yang, Z.-L., Levis, S., Sakaguchi, K., Bonan, G. B., and Slater, A. G.: Parameterization improvements and functional and structural advances in Version 4 of the Community Land Model, J. Adv. Model. Earth Sy., 3, M03001,, 2011. a

Lear, C. H., Bailey, T. R., Pearson, P. N., Coxall, H. K., and Rosenthal, Y.: Cooling and ice growth across the Eocene-Oligocene transition, Geology, 36, 251–254,, 2008. a

Licht, A., van Cappelle, M., Abels, H. A., Ladant, J.-B., Trabucho-Alexandre, J., France-Lanord, C., Donnadieu, Y., Vandenberghe, J., Rigaudier, T., Lécuyer, C., Terry, D., Adriaens, R., Boura, A., Guo, Z., Soe, A. N., Quade, J., Dupont-Nivet, G., and Jaeger, J.-J.: Asian monsoons in a late Eocene greenhouse world, Nature, 513, 501–506,, 2014. a

Liu, Z., Pagani, M., Zinniker, D., Deconto, R., Huber, M., Brinkhuis, H., Shah, S. R., Leckie, R. M., and Pearson, A.: Global Cooling During the Eocene-Oligocene Climate Transition, Science, 323, 1187–1190, 2009. a

Loptson, C. A., Lunt, D. J., and Francis, J. E.: Investigating vegetation–climate feedbacks during the early Eocene, Clim. Past, 10, 419–436,, 2014. a


Lunt, D. J., Farnsworth, A., Loptson, C., Foster, G. L., Markwick, P., O'Brien, C. L., Pancost, R. D., Robinson, S. A., and Wrobel, N.: Palaeogeographic controls on climate and proxy interpretation, Clim. Past, 12, 1181–1198,, 2016. a, b


O'Brien, C. L., Robinson, S. A., Pancost, R. D., Sinninghe Damsté, J. S., Schouten, S., Lunt, D. J., Alsenz, H., Bornemann, A., Bottini, C., Brassell, S. C., Farnsworth, A., Forster, A., Huber, B. T., Inglis, G. N., Jenkyns, H. C., Linnert, C., Littler, K., Markwick, P., McAnena, A., Mutterlose, J., Naafs, B. D. A., Püttmann, W., Sluijs, A., van Helmond, N. A., Vellekoop, J., Wagner, T., and Wrobel, N. E.: Cretaceous sea-surface temperature evolution: Constraints from TEX86 and planktonic foraminiferal oxygen isotopes, Earth-Sci. Rev., 172, 224–247,, 2017. a

Okafor, C. U., Thomas, D. J., Wade, B. S., and Firth, J.: Environmental change in the subtropics during the late middle Eocene greenhouse and global implications, Geochem. Geophy. Geosy., 10, Q07003,, 2009. a

Oleson, K. W., Lawrence, D. M., Bonan, G. B., Flanner, M. G., Kluzek, E., Lawrence, P. J., Levis, S., Swenson, S. C., Thornton, P. E., Dai, A., Decker, M., Dickinson, R., Feddema, J., Heald, C. L., Hoffman, F., Lamarque, J.-F., Mahowald, N., Niu, G.-Y., Qian, T., Randerson, J., Running, S., Sakaguchi, K., Slater, A., Stöckli, R., Wang, A., Yang, Z.-L., Zeng, X., and Zeng, X.: Technical Description of version 4.0 of the Community Land Model (CLM), NCAR Tech. Note, NCAR/TN-47, 2010. a

Passchier, S., Bohaty, S. M., Jiménez-Espejo, F., Pross, J., Röhl, U., Van De Flierdt, T., Escutia, C., and Brinkhuis, H.: Early Eocene to middle Miocene cooling and aridification of East Antarctica, Geochem. Geophy. Geosy., 14, 1399–1410,, 2013. a

Passchier, S., Ciarletta, D. J., Miriagos, T. E., Bijl, P. K., and Bohaty, S. M.: An antarctic stratigraphic record of stepwise ice growth through the eocene-oligocene transition, Bull. Geol. Soc. Am., 129, 318–330,, 2017. a

Pearson, P. N., Ditchfield, P. W., Singano, J., Harcourt-Brown, K. G., Nicholas, C. J., Olsson, R. K., Shackleton, N. J., and Hall, M. A.: Warm tropical sea surface temperatures in the Late Cretaceous and Eocene epochs, Nature, 413, 481–487,, 2001. a

Pearson, P. N., van Dongen, B. E., Nicholas, C. J., Pancost, R. D., Schouten, S., Singano, J. M., and Wade, B. S.: Stable warm tropical climate through the Eocene Epoch, Geology, 35, 211–214,, 2007. a, b

Peterse, F., van der Meer, J., Schouten, S., Weijers, J. W. H., Fierer, N., Jackson, R. B., Kim, J. H., and Sinninghe Damsté, J. S.: Revised calibration of the MBT-CBT paleotemperature proxy based on branched tetraether membrane lipids in surface soils, Geochim. Cosmochim. Ac., 96, 215–229,, 2012. a

Prothero, D. R.: Magnetic stratigraphy of the Eocene-Oligocene floral transition in western North America, in: Paleontology of the Upper Eocene Florissant Formation, Colorado, edited by: Meyer, H. W. and Smith, D. M., The Geological Society of America, Boulder, CO, 71–87, 2008. a

Quan, C., Liu, Y. S. C., and Utescher, T.: Eocene monsoon prevalence over China: A paleobotanical perspective, Palaeogeogr. Palaeocl., 365–366, 302–311,, 2012. a, b

Retallack, G. J., Orr, W. N., Prothero, D. R., Duncan, R. A., Kester, P. R., and Ambers, C. P.: Eocene-Oligocene extinction and paleoclimatic change near Eugene, Oregon, Bull. Geol. Soc. Am., 116, 817–839,, 2004. a

Royer, D. L., Pagani, M., and Beerling, D. J.: Geobiological constraints on Earth system sensitivity to CO2 during the Cretaceous and Cenozoic, Geobiology, 10, 298–310,, 2012. a, b, c, d

Scher, H. D., Bohaty, S. M., Smith, B. W., and Munn, G. H.: Eocene glaciation, Paleoceanography, 29, 628–644,, 2014. a, b

Schouten, S., Eldrett, J., Greenwood, D. R., Harding, I., Baas, M., and Damsté, J. S.: Onset of long-term cooling of Greenland near the Eocene-Oligocene boundary as revealed by branched tetraether lipids, Geology, 36, 147–150,, 2008. a


Sewall, J., Sloan, L., Huber, M., and Wing, S.: Climate sensitivity to changes in land surface characteristics, Global Planet. Change, 26, 445–465,, 2000. a, b, c

Sijp, W. P., von der Heydt, A. S., Dijkstra, H. A., Flögel, S., Douglas, P. M. J., and Bijl, P. K.: The role of ocean gateways on cooling climate on long time scales, Global Planet. Change, 119, 1–22,, 2014. a

Sijp, W. P., von der Heydt, A. S., and Bijl, P. K.: Model simulations of early westward flow across the Tasman Gateway during the early Eocene, Clim. Past, 12, 807–817,, 2016. a

Sluijs, A., Schouten, S., Pagani, M., Woltering, M., Brinkhuis, H., Damsté, J. S., Dickens, G. R., Huber, M., Reichart, G. J., Stein, R., Matthiessen, J., Lourens, L. J., Pedentchouk, N., Backman, J., Moran, K., Clemens, S., Cronin, T., Eynaud, F., Gattacceca, J., Jakobsson, M., Jordan, R., Kaminski, M., King, J., Koc, N., Martinez, N. C., McInroy, D., Moore, T. C., O'Regan, M., Onodera, J., Pälike, H., Rea, B., Rio, D., Sakamoto, T., Smith, D. C., St John, K. E., Suto, I., Suzuki, N., Takahashi, K., Watanabe, M., and Yamamoto, M.: Subtropical Arctic Ocean temperatures during the Palaeocene/Eocene thermal maximum, Nature, 441, 610–613,, 2006. a

Sluijs, A., Röhl, U., Schouten, S., Brumsack, H. J., Sangiorgi, F., Sinninghe Damsté, J. S., and Brinkhuis, H.: Article late Paleocene–Early Eocene paleoenvironments with special emphasis on the Paleocene-Eocene thermal maximum (Lomonosov Ridge, Integrated Ocean Drilling Program Expedition 302), Paleoceanography, 23, PA1S11,, 2008. a

Sluijs, A., Zeebe, R. E., Bijl, P. K., and Bohaty, S. M.: A middle Eocene carbon cycle conundrum, Nat. Geosci., 6, 429–434,, 2013. a

Smith, G. A., Manchester, S. R., Ashwill, M., McIntosh, W. C., and Conrey, R. M.: Late Eocene-early Oligocene tectonism, volcanism, and floristic change near Gray Butte, central Oregon, Bull. Geol. Soc. Am., 110, 759–778,<0759:LEEOTV>2.3.CO;2, 1998. a

Smith, R. D. and McWilliams, J. C.: Anisotropic horizontal viscosity for ocean models, Ocean Model., 5, 129–156,, 2003. a

Smith, R. D., Jones, P., Briegleb, B., Bryan, F., Danabasoglu, G., Dennis, J., Dukowicz, J., Eden, C., Fox-Kemper, B., Gent, P., Hecht, M., Jayne, S., Jochum, M., Large, W., Lindsay, K., Maltrud, M., Norton, N., Peacock, S., Vertenstein, M., and Yeager, S.: The Parallel Ocean Program (POP) reference manual: Ocean component of the Community Climate System Model (CCSM), available at: (last access: 1 December 2020), Los Alamos National Laboratory Tech. Rep., LAUR-10-01853, 141, 1–141, 2010. a

Song, X. and Zhang, G. J.: Convection parameterization, tropical pacific double ITCZ, and upper-ocean biases in the NCAR CCSM3. Part II: Coupled feedback and the role of ocean heat transport, J. Climate, 23, 800–812,, 2009. a

Steele, M., Morley, R., and Ermold, W.: PHC: A global ocean hydrography with a high-quality Arctic Ocean, J. Climate, 14, 2079–2087,<2079:PAGOHW>2.0.CO;2, 2001. a


Stickley, C. E., St John, K., Koç, N., Jordan, R. W., Passchier, S., Pearce, R. B., and Kearns, L. E.: Evidence for middle Eocene Arctic sea ice from diatoms and ice-rafted debris, Nature, 460, 376–379,, 2009. a

Tardif, D., Fluteau, F., Donnadieu, Y., Le Hir, G., Ladant, J.-B., Sepulchre, P., Licht, A., Poblete, F., and Dupont-Nivet, G.: The origin of Asian monsoons: a modelling perspective, Clim. Past, 16, 847–865,, 2020. a

Thiele, G. and Sarmiento, J.: Tracer Dating and Ocean Ventilation, J. Geophys. Res., 95, 9377–9391,, 1990. a

Toumoulin, A., Donnadieu, Y., Ladant, J. B., Batenburg, S. J., Poblete, F., and Dupont-Nivet, G.: Quantifying the Effect of the Drake Passage Opening on the Eocene Ocean, Paleoceanography and Paleoclimatology, 35, e2020PA003889,, 2020. a

Traiser, C., Klotz, S., Uhl, D., and Mosbrugger, V.: Environmental signals from leaves–a physiognomic analysis of European vegetation, New Phytol., 166, 465–484,, 2005. a


Tripati, A. K., Delaney, M. L., Zachos, J. C., Anderson, L. D., Kelly, D. C., and Elderfield, H.: Tropical sea-surface temperature reconstruction for the early Paleogene using Mg  Ca ratios of planktonic foraminifera, Paleoceanography, 18, 1101,, 2003. a

Uhl, D., Klotz, S., Traiser, C., Thiel, C., Utescher, T., Kowalski, E., and Dilcher, D. L.: Cenozoic paleotemperatures and leaf physiognomy – A European perspective, Palaeogeogr. Palaeocl., 248, 24–31,, 2007. a


Viebahn, J. P., von der Heydt, A. S., Le Bars, D., and Dijkstra, H. A.: Effects of Drake Passage on a strongly eddying global ocean, Paleoceanography, 31, 564–581,, 2016.  a

Wing, S. and Greenwood, D.: Fossils and fossil climate: The case for equable continental interiors in the Eocene, Philos. T. R. Soc. B, 341, 243–252,, 1993. a

Wolfe, J. A., Forest, C. E., and Molnar, P.: Paleobotanical evidence of Eocene and Oligocene paleoaltitudes in midlatitude western North America, Bull. Geol. Soc. Am., 110, 664–678,<0664:PEOEAO>2.3.CO;2, 1998. a

Yang, J., Spicer, R. A., Spicer, T. E., and Li, C. S.: 'CLAMP Online': A new web-based palaeoclimate tool and its application to the terrestrial Paleogene and Neogene of North America, Palaeobio. Palaeoenv., 91, 163–183,, 2011. a


Zachos, J. C., Stott, D., and Lohmann, K. C.: Evolution of early Cenozoic marine temperatures at Equator, Paleoceanography, 9, 353–387,, 1994. a

Zachos, J. C., Pagani, M., Sloan, L., Thomas, E., and Billups, K.: Trends, rhythms, and aberrations in global climate 65 Ma to present, Science, 292, 686–693,, 2001. a


Zhu, J., Poulsen, C. J., and Tierney, J. E.: Simulation of Eocene extreme warmth and high climate sensitivity through cloud feedbacks, Science Advances, 5, eaax1874,, 2019. a, b, c

Short summary
Warm climates of the deep past have proven to be challenging to reconstruct with the same numerical models used for future predictions. We present results of CESM simulations for the middle to late Eocene (∼ 38 Ma), in which we managed to match the available indications of temperature well. With these results we can now look into regional features and the response to external changes to ultimately better understand the climate when it is in such a warm state.