# Understanding the early Paleozoic carbon cycle balance and climate change from modelling

Chloé M. Marcilly, Pierre Maffre, Guillaume Le Hir, Alexandre Pohl, Frédéric Fluteau, Yves Goddéris, Yannick Donnadieu, Thea H. Heimdal, Trond Torsvik

► To cite this version:

Chloé M. Marcilly, Pierre Maffre, Guillaume Le Hir, Alexandre Pohl, Frédéric Fluteau, et al. Understanding the early Paleozoic carbon cycle balance and climate change from modelling. Earth and Planetary Science Letters, 2022, 594, pp.117717. 10.1016/j.epsl.2022.117717. hal-03870225

HAL Id: hal-03870225

https://cnrs.hal.science/hal-03870225

Submitted on 26 Nov 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

ELSEVIER

Contents lists available at [ScienceDirect](https://www.sciencedirect.com)

Earth and Planetary Science Letters

[www.elsevier.com/locate/epsl](https://www.elsevier.com/locate/epsl)

<figure><img src="image_3.png" /></figure>

# Understanding the early Paleozoic carbon cycle balance and climate change from modelling

<figure><img src="image_4.png" /></figure>

Check for
updates

Chloé M. Marcilly<sup>a,*</sup>, Pierre Maffre<sup>b</sup>, Guillaume Le Hir<sup>c</sup>, Alexandre Pohl<sup>d,e</sup>, Frédéric Fluteau<sup>c</sup>, Yves Goddéris<sup>f</sup>, Yannick Donnadieu<sup>g</sup>, Thea H. Heimdal<sup>a</sup>, Trond H. Torsvik<sup>a,h</sup>

## ARTICLE INFO

**Article history:**

Received 11 February 2022

Received in revised form 21 June 2022

Accepted 4 July 2022

Available online 20 July 2022

Editor: B. Wing

**Keywords:**

Ordovician

paleoclimate

carbon modelling

paleogeography

solid Earth degassing

## ABSTRACT

The Ordovician global cooling trend observed by several temperature proxies, which coincides with one of the most significant evolutionary diversifications on Earth, is yet to be fully understood. This study presents new simulations of pCO₂ and surface temperatures using a spatially resolved climate-carbon cycle Earth system model fed with refined continental reconstructions and new estimates of solid Earth degassing. First, in order to quantify the respective roles of paleogeography and degassing, the impact of continental drift alone is investigated. This is done by calculating its imprint on continental weathering rates throughout the Ordovician under a constant pCO₂, and by testing different topography inputs. Secondly, the sensitivity of the Earth's climate to paleogeography and carbon degassing changes is investigated by coupling the climate and carbon cycles in the GEOCLIM model. Based on our experiments we show that, although Early Ordovician high temperatures can be replicated within error margins, our new constraints cannot explain the intense cooling over the Mid to Late Ordovician, even if a progressive enhancement in Earth surface weatherability is taken into account. By using GEOCLIM in an inverse modelling approach, we calculate that the theoretical degassing necessary to reach proxy-derived temperatures for the Early Ordovician is three-to-five times higher compared to modern values. Further, in order to simulate the following Ordovician cooling trend, the solid Earth degassing must be reduced to modern-day values in only 30 Myrs. We conclude that, if accepting the veracity of the high Early Ordovician temperatures, alternative sources but also sinks of carbon must be considered to explain the climatic shift over the period.

© 2022 The Author(s). Published by Elsevier B.V. This is an open access article under the CC BY license ([http://creativecommons.org/licenses/by/4.0/](http://creativecommons.org/licenses/by/4.0/)).

## 1. Introduction

The Ordovician long-term global cooling trend is a key feature of the early Paleozoic climate (Trotter et al., 2008), and is also observed in recent studies (e.g. Edwards et al., 2021; Goldberg et al., 2021; Song et al., 2019). Global cooling initiated in the late Cambrian, and proceeded rapidly during the Early Ordovician, before tropical sea-surface temperatures (SSTs) reached

a plateau during the Middle to Late Ordovician and subsequently stabilized at modern-like levels (Fig. 1a) (i.e., ~30°C; Goldberg et al., 2021; Song et al., 2019; Trotter et al., 2008; Scotese et al., 2021). The striking feature of this climatic transition is its magnitude, with tropical temperatures plummeting by around 10°C, from ~40°C to 30°C, excluding the short-lived end-Ordovician Hirnantian cooling (445-444 Ma) (Finnegan et al., 2011; Trotter et al., 2008) (Fig. 1a). During the Ordovician cooling phase, the global climate transitioned from a greenhouse state to the first major continental glaciation of the Phanerozoic, with continental glacial onset argued to have occurred sometime around the Darriwilian (Fig. 1a,b) (~467-458 Ma; e.g. Pohl et al., 2016; Rasmussen et al.,

<footnote><sup>a</sup> Centre for Earth Evolution and Dynamics (CEED), University of Oslo, 0315 Oslo, Norway</footnote>

<footnote><sup>b</sup> Department of Earth and Planetary Science, University of California, Berkeley, CA 94720, USA</footnote>

<footnote><sup>c</sup> Université Paris, Institut de Physique du Globe de Paris, CNRS, 1 rue Jussieu, 75005 Paris, France</footnote>

<footnote><sup>d</sup> Biogéosciences, UMR 6282, UBFC/CNRS, Université Bourgogne Franche-Comté, 6 boulevard Gabriel, F-21000 Dijon, France</footnote>

<footnote><sup>e</sup> Department of Earth and Planetary Sciences, University of California, Riverside, CA, USA</footnote>

<footnote><sup>f</sup> Géosciences-Environnement Toulouse, CNRS-Université Paul Sabatier, Toulouse, France</footnote>

<footnote><sup>g</sup> Aix-Marseille Univ, CNRS, IRD, INRA, Coll. France, CEREGE, Aix-en-Provence, France</footnote>

<footnote><sup>h</sup> School of Geosciences, University of Witwatersrand, Johannesburg 2050, South Africa</footnote>

<footnote>* Corresponding author.
E-mail address: c.f.m.marcilly@geo.uio.no (C. M. Marcilly).</footnote>
2016) or the Katian (453-445 Ma; Cocks and Torsvik, 2020; Loi et al., 2010), but peak glacial activity occurred during the short-lived Hirnantian cooling (Finnegan et al., 2011; Ghienne et al., 2014).

Recent studies suggest an intricate coupling between climate and marine biological evolution throughout the Ordovician. The Great Ordovician Biodiversification Event (GOBE) marked the radiation of the marine animal phyla that had emerged during the Cambrian Explosion (Cocks and Torsvik, 2020). Biodiversification may have been driven by cooling to more habitable temperatures alone (Trotter et al., 2008; Rasmussen et al., 2016), or in conjunction with the ensuing ocean oxygenation (Edwards et al., 2017), both being known to act synergistically on marine ectotherm habitability (Deutsch et al., 2015; Stockey et al., 2021). Ordovician cooling was not only associated with a significant increase in the number of different complex animal forms during the Mid Ordovician (~470-455 Ma; GOBE), but the later Hirnantian cooling (with SSTs averaging $22 \pm 6^{\circ}C$) coincided with the only known mass extinction to have occurred during icehouse conditions (Torsvik et al., 2021).

The cause(s) of the long-term Ordovician cooling trend remain(s) poorly understood due to the sparsity of data and the large spread of temperature estimates (Fig. 1a). First, temporal trends from oxygen isotope data ($\delta^{18}\text{O}$) may be distorted by diagenesis, although meticulous screening for pristine samples and trace element analyses and the use of material supposedly resistant to diagenesis (such as biogenic apatite) permit gaining confidence in reconstructed temperatures (Hearing et al., 2018). Second, the heterogeneous spatial sampling of proxy data can create an underestimation of climatic perturbations, which has been shown for the Ordovician (Jones and Eichenseer, 2021). Absolute temperature values are even more challenging to constrain compared to temporal trends, due to large error margins associated with functions used to convert isotopic values to temperature (Lécuyer et al., 2013), as well as uncertain constraints on early Paleozoic seawater isotopic compositions. Some studies suggest that the seawater isotopic composition has remained unchanged (e.g. Henkes et al., 2018) while others follow the secular trend idea of the Veizer hypothesis (Veizer and Prokoph, 2015) of a gradual enrichment of the $\delta^{18}\text{O}$ isotope (e.g. Galili et al., 2019; Hodel et al., 2018; De fliese, 2020).

Several hypotheses have been proposed to explain the Ordovician long-term cooling, either by increasing the weatherability of the continents or by decreasing the amplitude of the solid Earth degassing. On the one hand, Earth's surface weatherability increases when continental drift brings continents into the warm and humid equatorial area (Nardin et al., 2011, Fig. 1b) in response to early land-plant colonization (Lenton et al., 2012; Porada et al., 2016), or, for example, through increased weathering of mafic silicate rock linked to physical erosion through changes in topography (Kump et al., 1999; Young et al., 2009; Swanson-Hysell and Macdonald, 2017; Macdonald et al., 2019). On the other hand, volcanic outgassing rates may have decreased through the early Paleozoic (McKenzie et al., 2016). The search for a unique causative mechanism for the long-term Ordovician cooling trend is elusive and these various factors may all have played a role. Unfortunately, these different hypotheses have been proposed in isolation using different methods and tools, and the importance and contribution of each of these factors remain to be quantified in an integrated numerical framework.

Here, we use the climate-carbon model GEOCLIM (Donnadieu et al., 2006-updated) to simulate the long-term Ordovician climatic trend in response to changes in paleogeography and associated silicate weathering activity as well as changes in volcanic degassing rates.

The model is fed by climatic fields that have been simulated using the atmospheric general circulation model LMDz6 (see Supplementary material) based on improved continental reconstructions (Marcilly et al., 2021). The GEOCLIM model accounts for the competing impacts of continental weathering and volcanic outgassing on the atmospheric partial pressure of $CO_2$ ($pCO_2$) and climate. This allows us to quantify the magnitude and temporal trend in global cooling, resulting from these competing mechanisms, and evaluate their relative contributions. In order to quantify the respective roles of paleogeography and degassing, we first explore the impact of continental drift alone, by calculating its imprint on continental weathering rates throughout the Ordovician under a constant $pCO_2$. Then, the sensitivity of the Earth's climate to paleogeography and carbon degassing variations is investigated by coupling the climate and carbon cycles in the GEOCLIM model. This study highlights the current limits of our understanding of factors driving the early Paleozoic climate and presents conceivable directions for future studies.

<figure><img src="image_5.png" /><figcaption>**Fig. 1.** (a) Compilation of temperature reconstructions covering the Late Cambrian to early Silurian. These datasets are from Song et al. (2019) in red, filtered to reflect only fossils found within $\pm 30^{\circ}$ latitude with error margins shaded. The stippled purple curve is from Goldberg et al. (2021) and represents the average of the least diagenetically altered samples. The purple shaded area represents the spread of data of the mid-to-least diagenetically altered data presented in Fig. 4 of Goldberg et al. (2021). Finally, the green dots are the temperatures reconstructed by Trotter et al. (2008) with their error band in shaded green. (b) Model simulations of atmospheric $CO_2$ levels. The green line is simulated using the forward model COPSE (Mills et al., 2019) whereas the stippled blue curve is obtained using GEOCARBSULF with updated paleogeography and degassing parameters (Marcilly et al., 2021). The dotted black curve is from Nardin et al. (2011) using GEOCLIM with the GCM FOAM. The numbers refer to the theoretical glacial inception levels with 3 and $6^{\circ}C$ climate sensitivity assuming 500 ppm inception level for modern times and a fainter sun in deep time (Marcilly et al., 2021). (For interpretation of the colours in the figure(s), the reader is referred to the web version of this article.)</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>**Fig. 2.** (a) Modelled Global Mean Surface Temperatures (GMST), continental mean, tropical mean and tropical (21°N-21°S) sea-surface temperature (SST) though time. (b) Mean continental runoff (solid line, left axis), which corresponds to precipitation minus evaporation and globally cumulated continental discharge (dashed line, right axis) corresponding to the area-integrated runoff through time. (c) Total silicate weathering fluxes predicted by GEOCLIM (relative to pre-industrial control weathering flux). Climatic variables from LMDz simulations conducted at 12 PAL $pCO_2$ (3360 ppmv) were used to force the weathering model.</figcaption></figure>

## 2. The paleogeographical imprint on Ordovician climate

### 2.1. Response of weathering to the continental configuration

To explore the interplays between Ordovician paleogeography and chemical weathering, we ran the LMDz6 climate model (see Supplementary Information) with late Cambrian to early Silurian paleogeographic maps ranging from 490 Ma to 430 Ma at 10 Myr intervals (7 time slices from Marcilly et al., 2021, see Supplementary Fig. S1). Several models for the tectonic and paleogeographic evolution exist for the Ordovician (e.g. Scotese, 2016; Blakey, 2009; Swanson-Hysell and Macdonald, 2017), but we here chose to work with the global model developed by Torsvik and collaborators (2012, 2016), which forms the basis for the exposed land maps produced in Marcilly et al. (2021; see details in supplementary section G).

In order to isolate the paleogeographical imprint, the $pCO_2$ in our simulations is fixed at 12 PAL (3360 ppmv) (i.e., Pre-industrial Atmospheric Level; 1 PAL = 280 ppmv) to avoid strongly nonlinear climate responses associated with sea-ice feedbacks encountered at lower $pCO_2$ levels (Pohl et al., 2014, 2016). We vary the model solar luminosity with geological age (from ~1313 W m$^{-2}$ at 490 Ma to ~1320 W m$^{-2}$ at 430 Ma) according to the model of stellar evolution (Gough, 1981).

At the global scale, the global mean surface temperature (GMST) remains relatively constant at ~22.5 °C (Fig. 2a) in our simulations. However, the continental climate (i.e., temperature, runoff, and discharge) is more sensitive to paleogeographical changes (Fig. 2a, b; Fig. 3). The most important observations include significant continental warming between 490 Ma and 480 Ma (Fig. 2a), as well as local maxima of global dryness at 480 Ma and 450 Ma (Fig. 2b, Fig. S2 in SI). From 490 to 480 Ma, the

global mean temperature over the continents increases by 4°C (from 14.5°C to 18.6°C), followed by fluctuations between 17°C and 19.3°C during the rest of the investigated time period (Fig. 2a). The warmest episodes occur at 480 Ma, 460 Ma, and 430 Ma when the Gondwana continent extends less over the South Pole region.

In the tropical band (i.e., from 21°S to 21°N), the annual mean temperature increases by 2°C (from 30.1°C to 32.1°C) from 480 to 430 Ma (Fig. 2a), mainly due to the solar constant increase. Despite the large fluctuations in emerged land surface, reaching a peak of ~36 × 10⁶ km² in the tropics by 450 Ma (40% of the total land exposed, Table S1 in SI), the paleogeography has apparently little influence on the modelled terrestrial temperature in the tropics (Fig. 2a).

In the subtropics, arid zones are larger at 480 Ma and 450 Ma compared to other periods (Fig. 3, Supplementary Fig. S2), which is when the surface of the Gondwana continent located below the descending branch of the Hadley circulation is at its maximum. The increase in aridity (i.e., areas where mean precipitation minus evaporation is inferior to 0) is also favoured by the drift of smaller blocks (Southern part of South-China and Baltica) entering into the subtropics between 460 and 440 Ma (SI Figs. S1, S2).

Precipitation and runoff reach maximum values over scattered islands in the equatorial belt and equatorial Gondwana (Fig. 3). A second belt with significant runoff is shown by all model runs at mid to high latitudes (i.e., between 40°S and 80°S), reaching maximum values over the western coast of Gondwana and the western edge of mountain ranges (Fig. 3; see also Pohl et al., 2017). The slight decoupling between the mean continental runoff and globally cumulated continental discharge at 450–440 Ma is due to an increase in the total area of emerged (exposed) lands between 490
<figure><img src="image_7.png" /><figcaption>**Fig. 3.** Runoff patterns for each time slice at the $pCO_2$ of 12 PAL (3360 ppmv). Mollweide projection (equal area), parallels plotted every 20° between 60°S and 40°N. B, Baltica, L, Laurentia, S, Siberia, La, Laurussia.</figcaption></figure>

and 450 Ma. This land area increase counteracts the drop in mean continental runoff. Collectively, these results explain why the Late Ordovician and early Silurian configurations (440 Ma, and mostly

430 Ma), and to a lesser extent, the 470 Ma configuration, between the two "dry intervals" (480 and 450 Ma), are associated with the highest weathering rates (Fig. 2c & Fig. 4, SI Fig. S4).
<figure><img src="image_8.png" /><figcaption>**Fig. 4.** Weathering patterns for each time slice with a constant $pCO_2$ of 12 PAL (3360 ppmv). Simulated weathering patterns are only driven by continental temperature and runoff, both depending on tectonic settings. B, Baltica, L, Laurentia, S, Siberia, La, Laurussia. Same projection and graticules as in Fig. 3.</figcaption></figure>

Overall, these results demonstrate a long-term increase in continental weathering in response to continental rearrangement during the Ordovician rather than a sudden rise at the end of the

Ordovician (Fig. 2c & 4). This progressive enhancement of Earth's surface weatherability (see also Marcilly et al., 2021) should lead to a $pCO_2$ reduction, a point explored in the following section.
<figure><img src="image_9.png" /><figcaption>**Fig. 5.** a) Equilibrium *pCO₂* with *constant* (black), and *variable* degassing from subduction fluxes (SBF; pink) and age-zircon distribution (Z; turquoise). The purple arrow and band display the glaciation threshold range estimated from climate models coupled with land-ice models (Herrmann et al., 2004; Lowry et al., 2014; Pohl et al., 2016). The numbers refer to the theoretical glacial inception levels with 2, 3 and 6 °C climate sensitivity assuming 500 ppm inception level from modern times and fainter sun activity in deeper time (Marcilly et al., 2021). b) Global mean surface temperature (GMST) of GEOCLIM, run with constant degassing (black) or variable degassing from SBF and Z (pink and turquoise). The temperature evolution at fixed CO₂ (12 PAL or 3360 ppmv) is added (orange) for comparison. Climatic trends shown with background colours in panel b are extracted from proxy compilations of different studies (Trotter et al., 2008; Song et al., 2019; Scotese et al., 2021; Goldberg et al., 2021) and simplified as follows: cooling until 480 Ma followed by rapid cooling until 450-440 Ma followed by a plateau.</figcaption></figure>

## 2.2. Climate-carbon response to continental drift

We explore the climatic evolution over the entire Ordovician assuming a constant degassing rate of the Earth. We couple the LMDz6 climate model with a long-term carbon cycle model to compute the steady-state climate, based on the climate-carbon negative feedback (Walker et al., 1981). The resulting climate-carbon model, GEOCLIM (Donnadieu et al., 2006-updated), is described in the Supplementary Information and available online. The calculated CO₂ consumption by silicate rock weathering is calibrated on the present-day state: this means that under present-day conditions, the model predicts the silicate weathering flux of the main world rivers (Park et al., 2020). The corresponding global silicate weathering (3.3 Tmol of CO₂ consumed each year) is close to the estimate by Gaillardet et al. (1999) (2.5 Tmol/yr). This implies that the solid Earth degassing must be close to the CO₂ consumption by silicate weathering. To reach a steady state at present day, we assume that the magmatic degassing of the present day Earth is equal to 3.3 Tmol/yr of CO₂. This value is close to the CO₂ degassing estimates by Gerlach (2011) (3.4 to 5.9 Tmol CO₂/yr). Therefore, to highlight the effect of paleogeography alone, the solid Earth degassing rate for all experiments was set to 7.1 Tmol C/yr and kept constant. This value is equivalent to about 2 times that of the present-day, which corresponds to the average degassing level over the studied period from proxies for seafloor spreading (the proxies and corresponding values are presented in section 3.1). As shown in Fig. 5a (black line), the *pCO₂* drops from 30 to 18 PAL (i.e., 8400-5040 ppm) between 490 and 430 Ma (Table S3 in SI), in response to the increasing surface weatherability (Fig. 2c).

The substantial drop in *pCO₂* observed between 480 Ma and 470 Ma induces an important global cooling event, with GMST decreasing by 2.5 °C (Fig. 5b; black line). Between 470 Ma and 430 Ma, although *pCO₂* slightly decreases (Fig. 5a), *pCO₂*-driven cooling is counteracted by the paleogeographical evolution accompanied by an increase in solar luminosity which, everything else unchanged, impose a warming from 470 Ma to 430 Ma (Fig. 5b, orange curve). As a result, the modelled global temperature slightly

increases from 470 Ma to 430 Ma (Fig. 5b; black curve). Therefore, in our model, the paleogeographical forcing alone cannot explain the long-term Ordovician cooling trend reported based on geological data (Goldberg et al., 2021; Song et al., 2019; Trotter et al., 2008).

## 2.3. Climate-carbon response to enhanced topography

In order to test the effect(s) of changing topography, we ran three sets of simulations with varying topographic inputs. The first simulation (i.e., run #1, "revised topo" in Fig. 6) aimed to test the effect of higher altitudes, and so the altitudes were doubled compared to our previous simulations. Secondly, the influence of more substantial island arcs and obducted arcs were tested (Macdonald et al., 2019). To do so, new maps with increased relief areas and same doubled altitudes as in run #1 were created (i.e., run #2, "increased topo" in Fig. 6; see SI Fig. S7 for maps). The third set of simulations (i.e., run #3, "increased topo + slope" in Fig. 6) are similar to the second, but with increased slopes to accentuate the effect of increased topography on physical erosion, and hence on chemical weathering (see SI Fig. S8). All tests were first performed with constant degassing in order to isolate the topography effect (Fig. 6a).

The higher altitude simulations (run #1) led to similar temperatures compared to the previous simulations with low topography, even slightly higher between 460-440 Ma, and no intensification of the cooling can be observed (Fig. 6a). As the influence of topography is increased by taking into account larger accretion zones and island arcs (run #2), the temperature at 440 Ma is slightly lowered (Fig. 6a), but almost negligible as showing only a -0.18 °C difference. Even with steepened slopes (run #3), the difference only reached -0.51 °C (Fig. 6a). Contrarily to our first simulations with constant degassing and lower topography (Section 2.2), which suggest a warming from 450 Ma, all the simulations with enhanced topography suggest cooling from 460 Ma and lasting until 440 Ma.

These results show that changing topography does not have a large influence on the simulated Ordovician climatic trend. Taking
<figure><img src="image_10.png" /><figcaption>**Fig. 6.** Global mean average surface temperature (GMST) simulated in GEOCLIM. Simulations with a) constant degassing, b) variable degassing from subduction fluxes, and c) variable degassing from zircon proxy. Each panel includes 4 simulations with varying topography. The solid line is the topography used in our primary set of maps, the dashed line with doubled relief altitude, the dotted line with doubled relief altitude and increased topography from Macdonald et al. (2019) and finally the dash-dotted line is the same as the dotted line but with increased slopes.</figcaption></figure>

into account larger topography and steeper slope does make the Earth slightly cooler (by ~1°C compared to the revised but not increased topography), but quasi uniformly on the whole period, since there is no significant temporal variation of global topography (growth nor decay) over the studied period. Therefore, the climatic trend is virtually unchanged, except for a slight prolonged cooling at 450-430 Ma in the increased topography and slope simulation. However, it is important to note that only modest changes in topography have been tested, as this parameter is poorly constrained through time.

# 3. Transition from a warm climate state to the ice age

## 3.1. Testing a variable degassing throughout the Ordovician

The Ordovician long-term cooling trend could potentially be attributed to fluctuations in the major long-term carbon sources. However, as direct measurements for solid Earth degassing in the

Ordovician do not exist, we must rely on proxies. The long-term solid Earth degassing mainly occurs at spreading ridges, volcanic arcs and continental rifts (Kelemen and Manning, 2015). Reconstructions of degassing fluctuations through time have traditionally been based on the past activity of seafloor spreading, and such estimates are commonly used in long-term carbon cycle models (Berner, 2006; Goddéris and Donnadieu, 2019; Mills et al., 2019; Pohl et al., 2020).

Here we use two proxies for seafloor production rate, which in turn act as a proxy for solid Earth degassing (Fig. 7). The first is based on estimates of subduction fluxes from the full-plate model of Domeier and Torsvik (2014) extended back to the early Paleozoic (Marcilly et al., 2022) (see SI for more information). Note that the early Paleozoic estimates have large uncertainties since all the oceanic lithosphere prior to the Jurassic has been subducted and absolute longitudes in the Paleozoic are strongly debated. We here estimate uncertainties associated with the subduction flux based
<figure><img src="image_11.png" /><figcaption>**Fig. 7.** Degassing relative to present-day flux derived from subduction fluxes ("subduction flux deg.") from a full-plate model (pink solid line) and age zircon distribution ("Zircon deg.") (turquoise dashed line). The proxies are plotted with their estimated error range. For comparison purposes, other estimations are plotted: Mills et al. (2019) based on hybrid subduction and rift lengths, Berner (2004) using sea-level inversion and calculated subduction flux from the Merdith et al. (2021) full plate model.</figcaption></figure>

on previous versions of our subduction flux estimates (Fig. 7). The second proxy is the age distribution of zircon presented in Marcilly et al. (2021). Considering the increasing uncertainty of subduction flux estimates with deeper time, the age-frequency relationship of detrital zircons, mimicking variations in the subduction flux, both regionally and globally has been suggested as a more objective proxy for seafloor production (Domeier et al., 2018; McKenzie et al., 2016; Marcilly et al., 2021). Uncertainties for this proxy are based on its difference with subduction flux estimates in more recent time (Fig. 7). This zircon-based proxy suggests solid Earth degassing values in the same order of magnitude compared to other proxies commonly used in long-term carbon cycle modelling (e.g., Berner, 2006; Mills et al., 2019; Goddéris and Donnadieu, 2019; Merdith et al., 2021; Fig. 7), although it suggests lower estimates towards the end-Ordovician and a slight decrease between 480 Ma and 430 Ma. In contrast, the subduction-flux-based proxy we use here suggests higher values, especially between 490 and 470 Ma.

By introducing an up-to-date variable degassing reconstruction from subduction fluxes in the GEOCLIM simulations, the GMST temporal trend is amplified (Fig. 5b; purple curve) compared to simulations performed with a solid degassing held constant (Fig. 5b; black curve). While the Early Ordovician becomes warmer compared to the constant degassing scenario ($T_{max}$ of 28.5 °C vs. 26 °C at 480 Ma; Fig. 5b), the model now predicts a cooler Late Ordovician climate (22 °C vs. 24.5 °C at 440 Ma; Fig. 5b). However, even by combining paleogeographical evolution and solid Earth degassing variations, our model fails to reproduce the full extent of the observed tropical SST drop derived from proxy data (Song et al., 2019; Goldberg et al., 2021; Fig. 8). The lower degassing scenario derived from the zircon-based proxy produces simulations generally more in line with proxy temperatures (Fig. 8c). However, the high temperatures of the Early Ordovician are well out of reach, with a maximum SST modelled at 32.5 °C at 480 Ma against

~40 °C from proxies (Fig. 8c). In summary, the simulation based on the subduction flux degassing matches more closely the proxy data between 490-480 Ma while the zircon proxy degassing curve is a better match from 470 Ma and onwards (Fig. 8b, & 8c), even if the total cooling intensity is only ~2.5 °C (Fig. 8c).

Increasing the topography with varying degassing does not lead to a more pronounced global cooling (Fig. 6b,c), and contrarily to simulations with constant degassing, this trend is not changed nor intensified between 460 and 440 Ma. Here as well, the change in values at 440 Ma is negligible with a maximum change of only -0.42 °C with subduction flux degassing (Fig. 6b) and -0.53 °C with zircon degassing (Fig. 6c).

Overall, our simulations match the Early Ordovician cooling trend, but cannot explain why the cooling continues after the Middle Ordovician (~460 Ma) and lasted for more than 40 Myrs (Fig. 8). It should be noted that, although various models and temperature reconstructions point towards a Mid-Ordovician cooling (Trotter et al., 2008; Rasmussen et al., 2016) there are no known Mid-Ordovician glacial deposits that can confirm this (Cocks and Torsvik, 2020). Furthermore, while the proxies suggest a total cooling of more than 10 °C (Fig. 8), our results show that the tectonic settings (from subduction fluxes) can drive a cooling close to 4 °C. This is in agreement with regional $\delta^{18}$O data from Baltica, suggesting a drop of 4-5 °C during the Middle Ordovician (Rasmussen et al., 2016).

## 3.2. Linking the $\delta^{18}$O record, paleogeography, and degassing

Regarding the conundrum of the Mid-Late Ordovician cooling, we complete our study by permitting the model parameters to reconcile modelled and proxy-derived SSTs. For that purpose, we adopt an inverse modelling approach by prescribing a tropical SST from the proxy reconstructions as a target in our model and calculating the degassing rate required to reproduce this Ordovician SST trend. We consider two alternative SST targets: after Goldberg et al. (2021) or after Song et al. (2019), both here averaged over 10 Myrs to fit our reconstruction intervals (see Fig. 1a & 8). The latter estimates have been checked to reflect only tropical temperatures (Marcilly et al., 2021). Further, Late Cambrian to early Silurian SSTs are almost exclusively derived from (apatite) conodonts that once lived from 12°N to 23°S (mean paleolatitude = 13±5°N/S). These SSTs should therefore reflect a reasonably well constrained cooling trend over the Ordovician with respect to the underestimation from biased geographical sampling described by Jones and Eichenseer (2021).

The two sets of temperature estimates differ mostly during the Late Cambrian-Early Ordovician, but following this time period, their trend and absolute values are similar and within error margins of each other (Figs. 1a, 8). While the estimates by Song et al. (2019) are based on phosphatic fossils (conodonts), the ones from Goldberg et al. (2021) are from bulk rock isotopes. The scarcity of estimates in both studies may represent a limitation mostly for 490 Ma where only a few points are analyzed. In particular, the spread of the estimates is the largest for the Ordovician between 490-480 Ma in Goldberg et al. (2021), and therefore might represent the least constrained interval of their study as all data for this age is extracted from the same location (i.e., Newfoundland at the equator at the time). Similarly, the data reported in Song et al. (2019) become depleted of data coverage for 490 Ma after filtering and therefore reflect only an extrapolation of the continuity between the late Cambrian and 480 Ma. Therefore, we cannot rule out a local anomaly such as a temperature low as reconstructed in Goldberg et al. (2021).

The degassing rate scenarios obtained using the inverse approach differ significantly compared to those previously used in this study (Fig. 9a). For example, the inverse approach suggests
<figure><img src="image_12.png" /><figcaption>**Fig. 8.** Equilibrium tropical sea-surface temperature (SST) simulated in GEOCLIM. Simulations with **a)** constant degassing, **b)** variable degassing from subduction fluxes, and **c)** variable degassing from zircon proxy. Shaded pink areas represent standard deviation. Purple and red solid lines represent mean tropical SST estimates from Goldberg et al. (2021) and Song et al. (2019), respectively.</figcaption></figure>

that 10.7 to 11.8 Tmol C/yr would be required in order to match the observed mean tropical SST at 480 Ma (i.e., 38-41 °C; Song et al., 2019; Goldberg et al., 2021), which is significantly higher compared to the solid Earth degassing rate set in our previous experiments (averaging to 7.1 Tmol C/yr; see section 2.2). On the other hand, to simulate temperatures cooler than 30 °C at 450 Ma, the solid Earth degassing must be reduced to modern-day values (Fig. 9a).

**4. Discussion**

The existence of extremely high temperatures in the Late Cambrian-Early Ordovician remains uncertain as interpretations of observed δ¹⁸O data are debated. Two models for the early Paleozoic temperature exist; the “hot” and “cool” models (Scotese et al., 2021). Studies supporting the “cool” model suggest that the long-term trend in seawater δ¹⁸O is due to an enrichment of the heavier

oxygen isotope (¹⁸O) through geological time and does therefore not reflect change of the water temperature. This model, often referred as the “Veizer hypothesis” (Veizer and Prokoph, 2015), argues that early Paleozoic seawater was enriched in the light oxygen isotope (¹⁶O) because of this evolution and that the high temperatures of the early Paleozoic are not accurate (Fig. 10). The main obstacle to this hypothesis is the lack of a primary mechanism that could trigger such a change, but also that other records that do not show the same feature. Those records support the “hot” model, which builds on the fact that the isotopic composition of seawater has remained constant so that the lighter δ¹⁸O values accurately reflect warmer paleo-temperatures during the early Paleozoic (e.g. Henkes et al., 2018; Bergmann et al., 2018) (Fig. 10a,b). Moreover, the peak in SST during this interval followed by a sharp cooling is well in agreement with the Cambrian radiation and mostly the GOBE demonstrating a strong (negative) correlation between biodiversity expansion and temperature (Cocks and Torsvik, 2020).
<figure><img src="image_13.png" /><figcaption>**Fig. 9.** a) CO₂ degassing rate required by GEOCLIM to match the observed tropical sea-surface temperature (SST) of Goldberg et al. (2021) (purple), Song et al. (2019) (red line) and for comparison degassing from subduction fluxes proxy (pink line) and zircon proxy (turquoise line). b) pCO₂ of the corresponding GEOCLIM runs as in a.</figcaption></figure>

If the “hot” model is correct, the state of the carbon mass balance during the Early Ordovician must have been largely different from the modern one. Indeed, the SST reconstructed from sedimentary $\delta^{18}$O records considering the “hot” hypothesis are critically high for the Early Ordovician and difficult to reproduce with a 3D-climate model. Tropical temperatures above 40 °C require very high atmospheric CO$_2$. Those levels will necessarily depend on the sensitivity of the climate model to atmospheric CO$_2$ (SI Fig. S5). But regardless of what this sensitivity was in the deep past, which has been estimated to be ~3 °C per CO$_2$ doubling during ice-free periods (greenhouse climates) and ~6 °C per CO$_2$ doubling in the presence of land ice (Royer, 2016), or fluctuating from 2 to 5 °C (IPCC working group 1: Solomon et al., 2007), simulating temperatures above 40 °C in the tropical ocean will necessarily require very high levels of CO$_2$. To generate a very warm tropical ocean, CO$_2$ levels as high as 300 times the present-day partial pressure are required in our models. By comparison, such levels are close to the post-snowball glaciation levels (Le Hir et al., 2008). We consider such CO$_2$ levels to be unrealistically high for the early Ordovician, however, they cannot be ruled out, because no reliable CO$_2$ proxies are extended further in the past than the Late Ordovician.

By targeting the proxy-based temperatures in our model, we have shown that the high CO$_2$ values prevailing in the early Ordovician would require a release flux of carbon into the ocean-atmosphere system of about 3 to 5 times that of the present-day volcanic outgassing (Fig. 9a). A factor of 3 to 5 decreasing to only 1 over 30 Myrs (Fig. 9a) would require drastic, and perhaps unrealistic, changes in the geodynamic settings since the largest reduction in degassing reconstructed from proxies is by a factor of 2.2 (observed over the past 120 Myrs; Marcilly et al., 2021).

If we consider the proxy record as being biased (such as in the Veizer hypothesis; Fig. 10), then the calculated temperatures become significantly lower for the Ordovician. If our models with both degassing scenarios fit this “cool” model until 470 Ma, the simulated temperatures plot well above the data of Veizer following this time period (Fig. 10b). In order to simulate such temperatures, an even sharper decrease in degassing needs to be invoked, which seems unrealistic considering the tectonic activity of the time period and that the zircon degassing proxy already suggests the lowest estimates among different reconstruction methods.

Reconstructing actual CO₂ outgassing levels of the early Paleozoic is challenging. Considering the tectonic settings of the early

Ordovician (490-480 Ma), with the initial rifting of the Rheic Ocean between Gondwana and Avalonia, a high rate of seafloor spreading is plausible. A change in volcanic outgassing has previously been argued to represent the principal driver for Ordovician climate change (McKenzie et al., 2016). However, no full-plate models, nor other methods to assess seafloor spreading (e.g., age-zircon distribution, sea level inversion, subduction zone lengths, continental rift lengths, or combinations of these) have ever suggested such high rates during the Phanerozoic or even the Neoproterozoic (Goddéris and Donnadieu, 2019; Marcilly et al., 2021; Mills et al., 2019).

Adding fluxes from continental rifts (Brune et al., 2017) and explicitly considering the intersection of continental magmatic arcs with carbonate-filled basins (e.g., Lee and Lackey, 2015) could certainly increase solid Earth degassing during the early Ordovician, but an additional release of carbon required to sustain high temperatures may come from other sources than solid Earth degassing. Imbalances in the organic carbon cycle may theoretically release CO$_2$ into the ocean and atmosphere if the oxidation of reduced sediments exposed on the continents exceeds the burial of organic carbon. The overall increase in the carbonate $\delta^{13}$C from the base of the Ordovician to the Upper Ordovician sections suggests a progressive predominance of burial versus oxidation (Veizer et al., 1999; Edwards and Saltzman, 2016; Hu et al., 2021) which is compatible with a decreasing carbon release from the Early to the Late Ordovician. However, the amplitude of the decrease needed to explain the Ordovician cooling is extremely large. If we assume that the “missing” Early Ordovician carbon source originates from the organic carbon sub-cycle instead of the outgassing flux, then the early Paleozoic must have been characterized by a very low organic carbon burial rate, or by significant oxidation of reduced sediments on the continents. Simple mass balance calculations show that if this were the case, the $\delta^{13}$C of seawater would be as low as -20‰ PDB in the Early Ordovician, far below the observed value of about -2‰ PDB. Finally, another process able to release carbon into the ocean-atmosphere system is the oxidative weathering of sedimentary pyrite exposed in the uplifted regions. This process generates sulfuric acid, which dissolves carbonate rocks, releasing carbon in the rivers and ultimately in the ocean (Torres et al., 2014). However, it has been shown that this carbon source is only transitory and cannot drive the long-term evolution of the carbon cycle (Maffre et al., 2020).
<figure><img src="image_14.png" /><figcaption>**Fig. 10.** (a) Oxygen isotope data from Veizer et al. (1999) in black. Changes following the secular trend from Veizer and Prokoph (2015) are indicated in grey. The dataset of the "warm" models from Bergmann et al. (2018) and Henkes et al. (2018) are represented with their error margins in shaded corresponding colours. (b) Reconstructed paleotemperatures from phosphatic fossils from Song et al. (2019), Bulk-rock δ¹⁸O from Goldberg et al. (2021) and carbonate fossil database of Veizer et al. (1999). The grey data points represent the calculated temperature of the Veizer dataset following changes due to secular trend as suggested by the Veizer hypothesis.</figcaption></figure>

Similar to the high temperatures of the Early Ordovician, our simulations do not replicate the cooling trend reconstructed based on proxy data (Figs. 1 & 8b,c). However, it should be noted that only the mean tropical SST has been considered, while proxy data have been collected at various geographical locations. Accounting for the spread of simulated SSTs at various longitudes within the tropical band reduces the model-data discrepancy (Fig. 8b,c). However, the corresponding ±1 σ SST envelope, more representative of common SST conditions found in the simulated climatic states, remains far from capturing the proxy-derived temperatures (Fig. 8b,c). Moreover, when invoking changing degassing, the correlation between our models and the latest temperature reconstructions (i.e., Song et al., 2019; Scotese et al., 2021 & Goldberg et al., 2021) is reasonable with an average Pearson correlation coefficient of 0.85 for both scenarios (Table 1).

Our simulations reveal that the paleogeographical evolution on its own cannot explain observed Ordovician cooling, and if no change in solid Earth degassing is accounted for, the continen-

tal drift and masses distribution imposes a warming from 470 Ma onward (Fig. 5b). This is highlighted by the very low Pearson coefficient correlation between the latest temperature reconstructions and our model (0.2; Table 1). This observation differs slightly from the findings of Nardin et al. (2011), even if, globally, both studies show a similar flat trend in temperature (i.e., maximum variation of ~3°C). However, these authors used the GEOCLIM model with a fixed degassing, as well as a simpler representation of continental weathering. Their climatic fields were derived from a previous generation of the model, based on a lower spatial resolution (FOAM) of different continental reconstructions (Nardin et al., 2011). Further, the temporal resolution of this study was also lower compared to our study, with one GEOCLIM simulation every 20 Myrs (vs. 10 Myrs). Discrepancies between Nardin et al. (2011) and the results presented here reflect therefore differences in the continental reconstructions, the climate model, and the representation of continental weathering (Park et al., 2020). The major difference between these two GEOCLIM studies is that the model
<table><caption>**Table 1**

Pearson correlation coefficients for our GEOCLIM simulations with different degassing scenarios compared to two proxy reconstructions of sea surface temperatures Song et al. (2019) and Goldberg et al. (2021).</caption><thead><tr><th>PEARSON</th><th>GEOCLIM no degassing</th><th>GEOCLIM zircon degassing</th><th>GEOCLIM subduction flux degassing</th></tr></thead><tbody><tr><td>Song et al. (2019)</td><td>0.24</td><td>0.89</td><td>0.83</td></tr><tr><td>Scotese et al. (2021)</td><td>0.19</td><td>0.84</td><td>0.87</td></tr><tr><td>Goldberg et al. (2021)</td><td>0.17</td><td>0.82</td><td>0.86</td></tr><tr><td>Average</td><td>0.20</td><td>0.85</td><td>0.85</td></tr></tbody></table>

of Nardin et al. (2011) is able to create a cooler climate by the end of the Ordovician, whereas the update provided here, benefiting from new constraints on paleogeography, does not support global cooling triggered by changes in land distribution.

Constraining $pCO_2$ thresholds for the glacial onset during the Ordovician has been the topic of a range of studies, going back to Crowley and Baum (1991). Some of these studies coupled climate models to land-ice models, which effectively simulate ice-sheet nucleation (Herrmann et al., 2004; Lowry et al., 2014; Pohl et al., 2016). Due to differences in the models and boundary conditions, values reported by these studies are widespread (1120 ppm to 3360 ppm; Herrmann et al., 2004; Lowry et al., 2014; Pohl et al., 2016). Pohl et al. (2016) used a similar climate modelling setup as our study (LMDz6 forced with FOAM-derived SSTs, see Supplementary Information), so these results might be the most comparable. The simulations presented here yield $CO_2$ levels below the 3360 ppm threshold as suggested by Pohl et al. (2016) by 440 Ma, at the end of a 40-Myr cooling (Fig. 5a). Therefore, we demonstrate that $CO_2$ levels for the Ordovician are compatible with the onset of a glaciation such as the short-lived Hirnantian glaciation (445-444 Ma), even though the intense global cooling as shown by proxies cannot be replicated. For comparison, other carbon cycle models, such as GEOCARBSULFvolc, estimate $CO_2$ levels almost similar to the simulations presented here when applying the same degassing parametrization (variable age-zircon distribution degassing; Marcilly et al., 2021). This study estimated theoretical glacial inception levels at 1035 ppm and 2142 ppm, based on climate sensitivities of 6 °C and 3 °C. Assuming an even lower climatic sensitivity (2 °C) yields a threshold as high as 4434 ppm (Fig. 5a).

Model-data disagreement in our simulations may come, at least in part, from the limitations of our models. Regarding the climatic component, results are necessarily model-dependent. More particularly, the double intertropical convergence zone, which represents a very general bias of ocean-atmosphere general circulation models, and is amplified in our simulations by the continental configuration featuring vast ocean basins at low latitudes, reduces precipitation rates over tropical emerged lands (Figs. 2, 3). This bias alters the climate-weathering relationship. In addition, due to the large number of simulations required to force the carbon cycle component, we did not run the ocean-atmosphere version of our CMIP6-class model. Instead, the atmospheric component (LMDz6) was forced with SSTs derived from an ocean-atmosphere model of a previous generation run under identical boundary conditions (FOAM see Supplementary Information). Such setup, used in previous work (Ladant et al., 2014; Licht et al., 2014; Pohl et al., 2016), may also lead to bias in our results, although this impact remains difficult to quantify. Importantly, we emphasize that uncertainty in the model climatic sensitivity would only alter the simulated $pCO_2$ and not affect simulated temperatures, which constitute our focus. Another source of uncertainty could come from the computation of weathering rates. The recent developments embedded in the version of GEOCLIM used here (Park et al., 2020) largely rely on continental topography, which is used to derive local slope and physical erosion rates (see Supplementary Information). Therefore, uncertainties in topography on land, inherent in every deep-time study due to the lack of robust and widespread proxies for paleo-elevation, constitute a limitation in our work.

The inability of our climate-carbon cycle model to simulate the long-term Ordovician cooling trend may also reflect missing mechanisms. In particular, including the impact of land plant colonization on weathering could help reconcile models and proxies (Porada et al., 2016; Lenton et al., 2012). However, the impact of non-vascular land plants on weathering is debated (Edwards et al., 2015; Porada et al., 2016). In addition, recent studies pinpoint the impact of the colonization scenario on weathering rates (Maffre et al., 2022), and this colonization remains poorly documented in the Ordovician. Consequently, the impact of non-vascular plants on Ordovician climate evolution remains elusive.

Another possible mechanism to explain the intensification of climate change during the Ordovician could reside in the weathering of obducted volcanic arcs as recently demonstrated by Conwell et al. (2022). Indeed, this study argues, using $^{87}$Sr/$^{86}$Sr seawater records, that a shift in the weathering source lithology is synchronous with the Ordovician cooling, and this can be linked to enhanced weathering of uplifted mafic oceanic crust associated with the Taconic orogeny. This hypothesis has also been suggested by several authors (i.e., Swanson-Hysell and Macdonald, 2017; Macdonald et al., 2019). Macdonald et al. (2019) placed Laurentia about 10° more northward between 480 and 444 Ma (Fig. 11a) compared with our reconstructions (see SI Fig S6), motivated by exhuming mafic and ultramafic lithologies during arc-continent collisions within the warm and wet tropical latitudes (±10°). Swanson-Hysell and Macdonald (2017) argued that a Taconic-related increase in silicate weatherability contributed to long-term Ordovician cooling, but this was not quantified in their study. One way of estimating changes in silicate weatherability is to examine the amount of exposed land where chemical weathering is assumed to be the most intense. Exposed land within ±10° was steadily increasing during the Ordovician (SI Table S1), peaking at around 450 Ma with ~22 × 10⁶ km² or ~24% of total land exposed between 10°S and 10°N. In the Late Ordovician, we estimate about 3000 km globally of ophiolite-bearing sutures in our paleogeographic model or ~5500 km in the Swanson-Hysell and Macdonald (2017) model. These estimates are based on Macdonald et al. (2019) but here revised since no Ordovician ophiolites are known in Greenland and all Ordovician ophiolites in Baltica (Fig. 11a) were obducted in the Silurian (e.g., Andersen and Andresen, 1994; Jakob et al., 2022). Fig. 1 in Macdonald et al. (2019) at 445 Ma, therefore, requires a major revision in suture lengths. In order to compare suture lengths with exposed land we (generously) assume an obduction width of 100 km, and as a simple theoretical exercise we multiply the resulting area by a factor of 10 to account for potentially higher weatherability of mafic/ultramafic rocks and elevated topography (Dessert et al., 2003; West et al., 2005). Note, however, that hydrothermal basaltic alteration and ultramafic serpentinization would have lowered their weathering efficiency compared to direct weathering of juvenile volcanics. This exercise shows that the potentially added weathering effect of obducted rocks in the tropics (compared to the amount of exposed land, mostly ‘granitic’) amounts to 27% at most, and as low as 3% considering the Scotese (2016) reconstructions (Fig. 11).

Adding the potential effect of obduction and increased weathering using the GEOCARBSULF model (e.g., increasing the $f_{AW}/f_A$ parameter, i.e. "fraction of land area undergoing chemical weathering"
<figure><img src="image_15.png" /><figcaption>**Fig. 11.** a) Late Ordovician reconstruction of the lapetus bordering continents. We place the southern margin of Laurentia at ~20°S (Model-I, Torsvik et al., 2012) whilst Swanson-Hysell and Macdonald (2017; see their figure DR1 at 445 Ma) place Laurentia ~10° more northward (Model-II, see supplementary section G). Thick dark brown lines in Model-II reconstruction were named 'ophiolite-bearing sutures' in Macdonald et al. (2019; see their Fig. 1 at 445 Ma) but there are no known ophiolites in Greenland and ophiolites in western Baltica were all obducted in the Silurian. This has been accounted for in Model-I reconstruction where no ophiolites in the lapetus realm were obducted within 10° from the equator (yellow band). Globally, reconstruction in Swanson-Hysell and Macdonald (2017) and Macdonald et al. (2019) are identical to those of Torsvik & Cocks (2016), except for locating Laurentia 10° more northward during the Ordovician (Supplementary Fig. G). This result in a much wider lapetus Ocean at 445 Ma, the Caroline Terrane has yet not collided with Laurentia and the Model-II requires fast latitudinal plate velocities (18 cm yr⁻¹) to close the lapetus at around 430 Ma. (b) Exposed land globally within 10° from the equator based on our Model I, Model-II (as Model-I except modified Laurentia) and Scotese (2016; Model III). All models are based on exposed land mapped out in 10 Myrs intervals in Marcilly et al. (2021) and peaks in exposed land in all models occur at 450 Ma. Swanson-Hysell and Macdonald (2017) argued for a peak in obducted arc lengths in the Late Ordovician but did not attempt to quantify the effect in terms of atmospheric CO₂ drawdown. As a simple back-of-the-envelope calculation, we estimate the effect of increased weathering by considering an obduction width of 100 km, which have been multiplied by 10 in order to account for a higher theoretical weathering efficiency (mafic versus felsic rocks) and elevated topography. Compared with the total exposed land (10°S-10°N) the added weathering potential can be (generously) estimated to 27% (Model-I), 15% (Model-II) and 3% (Model II). In terms of CO₂ drawdown that amounts to 275 ppm (Model-I), 150 ppm (Model-II) and 30 ppm (Model-III) in the GEOCARBsulf model. The red stippled line is based on the original estimates in Macdonald et al. (2019).</figcaption></figure>

ering”) results in a theoretical CO₂ reduction of 150 ppm in our preferred paleogeographic model (Fig. 11b). This exercise however, assumes that the climate gradients during the Ordovician were similar to modern times, but the occurrences of evaporites in Laurentia at low latitudes suggest a dry climate, rather than warm and wet tropics (e.g., Boucot et al., 2013; Torsvik & Cocks 2016; Torsvik et al., 2021).

When integrating island arcs and increased topography linked to obducted arcs as suggested by Macdonald et al. (2019) (Fig. 6, SI Fig. S6 and SI Fig. S7), our simulations did not yield any large effect linked to topography. Moreover, considering the trend in cooling enhancement seen in the simulations with a shift from low to moderate topography, only a considerable change towards very high relief can trigger large fluctuations.

## 5. Conclusion

By using a spatially resolved climate-carbon cycle Earth system model fed with updated continental reconstructions and new estimates of volcanic degassing, we show that the Early Ordovician high temperatures can be modelled within error margins of the proxy data. However, our new constraints cannot explain the intense cooling over the Mid to Late Ordovician even if revealing a progressive enhancement in Earth's surface weatherability. The calculated, theoretical degassing scenario necessary to reproduce the proxy-derived Ordovician long-term cooling trend features degassing rates reaching 3 to 5 times that of the current value, which

is not supported by geodynamical models nor geological proxies. Therefore, we interpret the inability of our up-to-date climate-carbon cycle model to reproduce Ordovician cooling as the existence of a missing parameter or a faulty influence of a present one. Accepting the veracity of the high Early Ordovician temperatures, alternative sources but also sinks of carbon must be considered. The contribution of the organic carbon sub-cycle and oxidative weathering of pyrite is currently not supported by proxy data. The impact of the colonization of land surfaces by non-vascular plants may help bring models and data into a closer agreement, but the overall impact of vegetation on climate remains difficult to quantify. Our simulations with increased topography did not intensify the Ordovician cooling, but note that the model used here holds large limitations considering the integration of relief and slopes. A more detailed analysis on this subject is necessary to draw any firm conclusions.

## CRediT authorship contribution statement

**Chloé M. Marcilly:** Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Software, Visualization, Writing – original draft. **Pierre Maffre:** Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Resources, Software, Visualization, Writing – review & editing. **Guillaume Le Hir:** Conceptualization, Formal analysis, Investigation, Resources, Visualization, Writing – review & editing. **Alexandre Pohl:** Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Visualization, Writing – review & editing.
ology, Resources, Software, Writing – review & editing. **Frédéric Fluteau:** Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Resources, Software, Visualization, Writing – review & editing. **Yves Goddéris:** Conceptualization, Data curation, Formal analysis, Investigation, Resources, Software, Supervision, Writing – review & editing. **Yannick Donnadieu:** Conceptualization, Data curation, Formal analysis, Investigation, Resources, Software, Writing – review & editing. **Thea H. Heimdal:** Supervision, Writing – review & editing. **Trond H. Torsvik:** Methodology, Resources, Supervision, Writing – review & editing.

## Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

## Acknowledgements

Calculations were performed using HPC resources from DNUM CCUB (Centre de Calcul de l'Université de Bourgogne). This work was granted access to the HPC resources of TGCC under the allocations 2019-A0070107601 and 2020-A0090107601 made by GENCI. This project has received funding from the European Union's Horizon 2020 research and innovation program under the Marie Skłodowska-Curie grant agreement No. 838373, and the Research Council of Norway (RCN), through its Centres of The Excellence funding scheme, project number 223272 (CEED). This is a contribution to UNESCO project IGCP 735 "Rocks and the Rise of Ordovician Life" (Rocks n' ROL).

## Appendix A. Supplementary material

Supplementary material related to this article can be found online at https://doi.org/10.1016/j.epsl.2022.117717.

## References

Andersen, T.B., Andresen, A., 1994. Stratigraphy, tectonostratigraphy and the accretion of outboard terranes in the Caledonides of Sunnhordland, W. Norway. Tectonophysics 231, 71–84.

Bergmann, K.D., Finnegan, S., Creel, R., Eiler, J.M., Hughes, N.C., Popov, L.E., Fischer, W.W., 2018. A paired apatite and calcite clumped isotope thermometry approach to estimating Cambro-Ordovician seawater temperatures and isotopic composition. Geochim. Cosmochim. Acta 224, 18–41. https://doi.org/10.1016/j.gca.2017.11.015.

Berner, R.A., 2004. The Phanerozoic carbon cycle: CO₂ and O₂. Oxford University Press on Demand.

Berner, R.A., 2006. GEOCARBSULF: a combined model for Phanerozoic atmospheric O₂ and CO₂. Geochim. Cosmochim. Acta 70 (23), 5653–5664. https://doi.org/10.1016/j.gca.2005.11.032.

Blakey, R.C., 2009. Gondwana paleogeography from assembly to break-up—a 500 my odyssey. Spec. Pap., Geol. Soc. Am. 441, 1–28.

Boucot, A.J., Xu, C., Scotese, C.R., 2013. Phanerozoic paleoclimate: an atlas of lithologic indicators of climate. In: SEPM Concepts in Sedimentology and Paleontology, p. 11.

Brune, S., Williams, S.E., Müller, R.D., 2017. Potential links between continental rifting, CO₂ degassing and climate change through time. Nat. Geosci. 10, 941–946. https://doi.org/10.1038/s41561-017-0003-6.

Cocks, L.R.M., Torsvik, T.H., 2020. Ordovician palaeogeography and climate change. Gondwana Res. 100, 53–72. https://doi.org/10.1016/j.gr.2020.09.008.

Conwell, C.T., Saltzman, M.R., Edwards, C.T., Griffith, E.M., Adiatma, Y.D., 2022. Nd isotopic evidence for enhanced mafic weathering leading to Ordovician cooling. Geology.

Crowley, T.J., Baum, S.K., 1991. Toward reconciliation of Late Ordovician (440 Ma) glaciation with very high CO₂ levels. J. Geophys. Res., Atmos. 96, 22597–22610. https://doi.org/10.1029/91JD02449.

Defliese, W.F., 2020. The impact of Snowball Earth glaciation on ocean water δ18O values. Earth Planet. Sci. Lett. 554, 116661. https://doi.org/10.1016/j.epsl.2020.116661.

Dessert, C., Dupré, B., Gaillardet, J., François, L.M., Allegre, C.J., 2003. Basalt weathering laws and the impact of basalt weathering on the global carbon cycle. Chem. Geol. 202 (3–4), 257–273.

Deutsch, C., Ferrel, A., Seibel, B., Pörtner, H.O., Huey, R.B., 2015. Climate change tightens a metabolic constraint on marine habitats. Science 348 (6239), 1132–1135. https://doi.org/10.1126/science.aaa1605.

Donnadieu, Y., Goddéris, Y., Pierrehumbert, R., Dromart, G., Fluteau, F., Jacob, R., 2006. A GEOCLIM simulation of climatic and biogeochemical consequences of Pangea breakup. Geochem. Geophys. Geosyst. 7 (11), Q11019. http://doi.wiley.com/10.1029/2006GC001278.

Domeier, M., Torsvik, T.H., 2014. Focus review paper: plate kinematics of the Late Paleozoic. Geosci. Front. 5, 303–350. https://doi.org/10.1016/j.gsf.2014.01.002.

Domeier, M., Magni, V., Hounslow, M.W., Torsvik, T.H., 2018. Episodic zircon age spectra mimic fluctuations in subduction. Sci. Rep. 8, 1747. https://doi.org/10.1038/s41598-018-35040-z.

Edwards, D., Cherns, L., Raven, J.A., 2015. Could land-based early photosynthesizing ecosystems have bioengineered the planet in mid-Paleozoic times? Palaeontology 58, 803–837. https://doi.org/10.1111/pala.12187.

Edwards, C.T., Saltzman, M.R., 2016. Paired carbon isotopic analysis of Ordovician bulk carbonate (δ13Ccarb) and organic matter (δ13Corg) spanning the Great Ordovician Biodiversification Event. Palaeogeogr. Palaeoclimatol. Palaeoecol. 458, 102–117. https://doi.org/10.1016/j.palaeo.2015.08.005.

Edwards, C.T., Saltzman, M.R., Royer, D.L., Fike, D.A., 2017. Oxygenation as a driver of the Great Ordovician Biodiversification Event. Nat. Geosci. 10, 925–929. https://doi.org/10.1038/s41561-017-0006-3.

Edwards, C.T., Jones, C.M., Quinton, P.C., Fike, D.A., Sciences, P., Louis, S., 2021. Oxygen isotope (18 O) trends measured from Ordovician conodont apatite using secondary ion mass spectrometry (SIMS): implications for paleo-thermometry studies. GSA Bull. 1, 1–14. https://doi.org/10.1130/B35891.1/5294555/b35891.pdf.

Finnegan, S., Bergmann, K., Eiler, J.M., Jones, D.S., Fike, D.A., Eisenman, I., Hughes, N.C., Tripati, A.K., Fischer, W.W., 2011. The magnitude and duration of late Ordovician-early Silurian glaciation. Science 331 (6019), 903–906. https://doi.org/10.1126/science.1200803.

Gaillardet, J., Dupré, B., Louvat, P., Allegre, C.J., 1999. Global silicate weathering and CO₂ consumption rates deduced from the chemistry of large rivers. Chem. Geol. 159 (1–4), 3–30. https://doi.org/10.1016/S0009-2541(99)00031-5.

Galili, N., Shemesh, A., Yam, R., Brailovsky, I., Sela-Adler, M., Schuster, E.M., Halevy, I., 2019. The geologic history of seawater oxygen isotopes from marine iron oxides. Science 365 (6452), 469–473. https://doi.org/10.1126/science.aaw9247.

Gerlach, T., 2011. Volcanic versus anthropogenic carbon dioxide. Eos Trans. AGU 92 (24), 201–202. https://doi.org/10.1029/2011EO240001.

Ghienne, J.F., Desrochers, A., Vandenbroucke, T., et al., 2014. A Cenozoic-style scenario for the end-Ordovician glaciation. Nat. Commun. 5, 4485. https://doi.org/10.1038/ncomms5485.

Goddéris, Y., Donnadieu, Y., 2019. A sink-or a source-driven carbon cycle at the geological timescale? Relative importance of palaeogeography versus solid Earth degassing rate in the Phanerozoic climatic evolution. Geol. Mag. 156 (2), 355–365. https://doi.org/10.1017/S0016756817001054.

Goldberg, S.L., Present, T.M., Finnegan, S., Bergmann, K.D., 2021. A high-resolution record of early Paleozoic climate. Proc. Natl. Acad. Sci. U.S.A. 118 (6), e2013083118. https://doi.org/10.1073/pnas.2013083118.

Gough, D.O., 1981. Solar interior structure and luminosity variations. Sol. Phys. 74 (13), 21–34. https://doi.org/10.1007/BF00151270.

Hearing, T.W., Harvey, T.H.P., Williams, M., Leng, M.J., Lamb, A.L., Wilby, P.R., Gabbott, S.E., Pohl, A., Donnadieu, Y., 2018. An early Cambrian greenhouse climate. Sci. Adv. 4 (5), eaar5690. https://doi.org/10.1126/sciadv.aar5690.

Henkes, G.A., Passey, B.H., Grossman, E.L., Shenton, B.J., Yancey, T.E., Pérez-Huerta, A., 2018. Temperature evolution and the oxygen isotope composition of Phanerozoic oceans from carbonate clumped isotope thermometry. Earth Planet. Sci. Lett. 490, 40–50. https://doi.org/10.1016/j.epsl.2018.02.001.

Herrmann, A.D., Patzkowsky, M.E., Pollard, D., 2004. The impact of paleogeography, pCO₂, poleward ocean heat transport and sea level change on global cooling during the Late Ordovician. Palaeogeogr. Palaeoclimatol. Palaeoecol. 206 (1), 59–74. https://doi.org/10.1016/j.palaeo.2003.12.019.

Hodel, F., Macouin, M., Trindade, R.I.F., Triantafyllou, A., Ganne, J., Chavagnac, V., Agrinier, P., 2018. Fossil black smoker yields oxygen isotopic composition of Neoproterozoic seawater. Nat. Commun. 9 (1), 1–7. https://doi.org/10.1038/s41467-018-03890-w.

Hu, D., Zhang, X., Li, M., Xu, Y., Shen, Y., 2021. Carbon isotope (δ13Ccarb) stratigraphy of the Lower-Upper Ordovician of the Yangtze Platform, South China: Implications for global correlation and the Great Ordovician Biodiversification Event (GOBE). Glob. Planet. Change, 103546. https://doi.org/10.1016/j.gloplacha.2021.103546.

Jakob, J., Andersen, T.B., Mohn, G., Kjøll, H.J., Beyssac, O., 2022. A revised tectono-stratigraphic scheme for the Scandinavian Caledonides and its implications for our understanding of the Scandian orogeny. Spec. Pap., Geol. Soc. Am. In New Developments in the Appalachian-Caledonian-Variscan Orogen.

Jones, L.A., Eichenseer, K., 2021. Uneven spatial sampling distorts reconstructions of Phanerozoic seawater temperature. Geology. https://doi.org/10.1130/G49132.1.

Kelemen, P.B., Manning, C.E., 2015. Reevaluating carbon fluxes in subduction zones, what goes down, mostly comes up. Proc. Natl. Acad. Sci. U.S.A. 112 (30), E3997–E4006.
Kump, L.R., Arthur, M.A., Patzkowsky, M.E., Gibbs, M.T., Pinkus, D.S., Sheehan, P.M., 1999. A weathering hypothesis for glaciation at high atmospheric pCO2 during the Late Ordovician. Palaeogeogr. Palaeoclimatol. Palaeoecol. 152 (1-2), 173-187.

Ladant, J.B., Donnadieu, Y., Lefebvre, V., Dumas, C., 2014. The respective role of atmospheric carbon dioxide and orbital parameters on ice sheet evolution at the Eocene-Oligocene transition. Paleogeography 29, 810-823. https://doi.org/10.1002/2013PA002593.

Le Hir, G., Donnadieu, Y., Goddéris, Y., Pierrehumbert, R.T., Halverson, G.P., Macouin, M., Ramstein, G., 2008. The snowball Earth aftermath: exploring the limits of continental weathering processes. Earth Planet. Sci. Lett. 277 (3-4), 453-463. https://doi.org/10.1016/j.epsl.2008.11.010.

Lécuyer, C., Amiot, R., Touzeau, A., Trotter, J., 2013. Calibration of the phosphate δ18O thermometer with carbonate-water oxygen isotope fractionation equations. Chem. Geol. 347, 217-226. https://doi.org/10.1016/j.chemgeo.2013.03.008.

Lee, C.T.A., Lackey, J.S., 2015. Global continental arc flare-ups and their relation to long-term greenhouse conditions. Elements 11 (2), 125-130. https://doi.org/10.2113/gselements.11.2.125.

Lenton, T., Crouch, M., Johnson, M., et al., 2012. First plants cooled the Ordovician. Nat. Geosci. 5, 86-89. https://doi.org/10.1038/ngeo1390.

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. Jr, Adriaens, R., Boura, A., Guo, Z., Soe, A.N., Quade, J., Dupont-Nivet, G., Jaeger, J.J., 2014. Asian monsoons in a late Eocene greenhouse world. Nature 513, 501-506. https://doi.org/10.1038/nature13704.

Loi, A., Ghienne, J.-F., Dabard, M.P., Paris, F., Botquelen, A., Christ, N., Elaouad-Debbaj, Z., Gorini, A., Vidal, M., Videt, B., Destombes, J., 2010. The Late Ordovician glacio-eustatic record from a high-latitude storm-dominated shelf succession: The Bou Ingar section (Anti-Atlas, Southern Morocco). Palaeogeogr. Palaeoclimatol. Palaeoecol. 296, 332-358. https://doi.org/10.1016/j.palaeo.2010.01.018.

Lowry, D.P., Poulsen, C.J., Horton, D.E., Torsvik, T.H., Pollard, D., 2014. Thresholds for Paleozoic ice sheet initiation. Geology 42, 627-630. https://doi.org/10.1130/G35615.1.

Macdonald, F.A., Swanson-Hysell, N.L., Park, Y., Lisiecki, L., Jagoutz, O., 2019. Arc-continent collisions in the tropics set Earth's climate state. Science 364 (6436), 181-184.

Maffre, P., Swanson-Hysell, N.L., Goddéris, Y., 2020. Limited carbon cycle response to increased sulfide weathering due to oxygen feedback. Geophys. Res. Lett. 48 (19), e2021GL094589. https://doi.org/10.1029/2021GL094589.

Maffre, P., Goddéris, Y., Pohl, A., Donnadieu, Y., Carretier, S., Le Hir, G., 2022. The complex response of continental silicate rock weathering to the colonization of the continents by vascular plants in the Devonian. Am. J. Sci. 322 (3), 461-492.

Marcilly, C.M., Torsvik, T.H., Domeier, M., Royer, D.L., 2021. New paleogeographic and degassing parameters for long-term carbon cycle models. Gondwana Res. https://doi.org/10.1016/j.gr.2021.05.016.

Marcilly, C.M., Torsvik, T.H., Conrad, C.P., 2022. Global Phanerozoic sea levels from paleogeographic flooding maps. Gondwana Res.

Merdith, A.S., Williams, S.E., Collins, A.S., Tetley, M.G., Mulder, J.A., Blades, M.L., Müller, R.D., 2021. Extending full-plate tectonic models into deep time: linking the Neoproterozoic and the Phanerozoic. Earth-Sci. Rev. 214, 103477. https://doi.org/10.1016/j.earscirev.2020.103477.

McKenzie, N.R., Horton, B.K., Loomis, S.E., Stockli, D.F., Planavsky, N.J., Lee, C.-T.A., 2016. Continental arc volcanism as the principal driver of icehouse-greenhouse variability. Science 352 (6284), 444-447. https://doi.org/10.1126/science.aad5787.

Mills, B.J., Krause, A.J., Scotese, C.R., Hill, D.J., Shields, G.A., Lenton, T.M., 2019. Modelling the long-term carbon cycle, atmospheric CO2 and Earth surface temperature from late Neoproterozoic to present day. Gondwana Res. 67, 172-186. https://doi.org/10.1016/j.gr.2018.12.001.

Nardin, E., Goddéris, Y., Donnadieu, Y., Le Hir, G., Blakey, R.C., Pucéat, E., Aretz, M., 2011. Modeling the early Paleozoic long-term climatic trend. Bull. Geol. Soc. Am. 123 (5), 1181-1192. https://doi.org/10.1130/B30364.1.

Park, Y., Maffre, P., Goddéris, Y., Macdonald, F.A., Anttila, E.S.C., Swanson-Hysell, N.L., 2020. Emergence of the Southeast Asian islands as a driver for Neogene cooling. Proc. Natl. Acad. Sci. 117 (41), 25319-25326. https://doi.org/10.1073/pnas.2011033117.

Pohl, A., Donnadieu, Y., Le Hir, G., Buoncristiani, J.F., Vennin, E., 2014. Effect of the Ordovician paleogeography on the (in)stability of the climate. Clim. Past 10 (6), 2053-2066. https://doi.org/10.5194/cp-10-2053-2014.

Pohl, A., Donnadieu, Y., Le Hir, G., Ladant, J.B., Dumas, C., Alvarez-Solas, J., Vandembroucke, T.R.A., 2016. Glacial onset predated Late Ordovician climate cooling. Paleogeography 31, 800-821. https://doi.org/10.1002/2016PA002928.

Pohl, A., Donnadieu, Y., Le Hir, G., Ferreira, D., 2017. The climatic significance of Late Ordovician-early Silurian black shales. Paleogeography 32, 397-423. https://doi.org/10.1002/2016PA003064.

Pohl, A., Donnadieu, Y., Godderis, Y., Lanteaume, C., Hairabian, A., Frau, C., Michel, J., Laugié, M., Reijmer, J.J.G., Scotese, C.R., Jean, B., 2020. Carbonate platform production during the Cretaceous. GSA Bull. 132, 2606-2610. https://doi.org/10.1130/B35680.1.

Porada, P., Lenton, T.M., Pohl, A., Weber, B., Mander, L., Donnadieu, Y., Beer, C., Pöschl, U., Kleidon, A., 2016. High potential for weathering and climate effects of non-vascular vegetation in the Late Ordovician. Nat. Commun. 7, 12113. https://doi.org/10.1038/ncomms12113.

Rasmussen, C.M. Ø, Ullmann, C.V., Jakobsen, K.G., Lindskog, A., Hansen, J., Hansen, T., Eriksson, M.E., Dronov, A., Frei, R., Korte, C., Nielsen, A.T., Harper, D.A.T., 2016. Onset of mid Phanerozoic marine radiation sparked by emerging Mid Ordovician icehouse. Sci. Rep. 6, 18884. https://doi.org/10.1038/srep18884.

Royer, D.L., 2016. Climate sensitivity in the geologic past. Annu. Rev. Earth Planet. Sci. 44. https://doi.org/10.1146/annurev-earth-100815-024150.

Scotese, C.R., 2016. PALEOMAP Paleozo Atlas for GPlates and the PaleoData plotter program. Abstr. Program – Geol. Soc. Am. 48 (5). http://www.scotese.com/newpage1.htm.

Scotese, C.R., Song, H., Mills, B.J., van der Meer, D.G., 2021. Phanerozoic paleotemperatures: the Earth's changing climate during the last 540 million years. Earth-Sci. Rev. 103503. https://doi.org/10.1016/j.earscirev.2021.103503.

Solomon, S., Manning, M., Marquis, M., Qin, D., 2007. Climate change 2007-the physical science basis: Working group I contribution to the fourth assessment report of the IPCC (Vol. 4). Cambridge University Press.

Song, H., Wignall, P.B., Song, H., Dai, X., Chu, D., 2019. Seawater temperature and dissolved oxygen over the past 500 million years. J. Earth Sci. 30 (2), 236-243. https://doi.org/10.1007/s12583-018-1002-2.

Stockey, R., Pohl, A., Ridgwell, A., Finnegan, S., Sperling, E.A., 2021. Decreasing Phanerozoic extinction intensity is a predictable consequence of Earth surface oxygenation and metazoan ecophysiology. Proc. Natl. Acad. Sci. 118, e2101900118. https://doi.org/10.1073/pnas.2101900118.

Swanson-Hysell, N.L., Macdonald, F.A., 2017. Tropical weathering of the Taconic orogeny as a driver for Ordovician cooling. Geology, G38985.1. http://geology.gsapubs.org/lookup/doi/10.1130/G38985.1.

Torres, M.A., West, A.J., Li, G., 2014. Sulphide oxidation and carbonate dissolution as a source of CO 2 over geological timescales. Nature 507 (7492), 346-349. https://doi.org/10.1038/nature13030.

Torsvik, T.H., Van der Voo, R., Preeden, U., Mac Niocaill, C., Steinberger, B., Doubrovine, P.V., Cocks, L.R.M., 2012. Phanerozoic polar wander, palaeogeography and dynamics. Earth-Sci. Rev. 114 (3-4), 325-368. https://doi.org/10.1016/j.earscirev.2012.06.007.

Torsvik, T.H., Svensen, H.H., Steinberger, B., Royer, D.L., Jerram, D.A., Jones, M.T., Domeier, M., 2021. Connecting the deep Earth and the atmosphere. In: Marquardt, H., Ballmer, M., Cottaar, S., Konter, J. (Eds.), Mantle Convection and Surface Expressions. In: Geophysical Monograph, vol. 263. American Geophysical, Union, pp. 413-453.

Trotter, J.A., Williams, I.S., Barnes, C.R., Lécuyer, C., Nicoll, R.S., 2008. Did cooling oceans trigger Ordovician biodiversification? Evidence from conodont thermometry. Science 321 (5888), 550-554. https://doi.org/10.1126/science.1155814.

Veizer, J., Ala, D., Azmy, K., Bruckschen, P., Buhl, D., Bruhn, F., Strauss, H., 1999. 87Sr/86Sr, δ13C and δ18O evolution of Phanerozoic seawater. Chem. Geol. 161 (1-3), 59-88. https://doi.org/10.1016/S0009-2541(99)00081-9.

Veizer, J., Prokoph, A., 2015. Temperatures and oxygen isotopic composition of Phanerozoic oceans. Earth-Sci. Rev. 146, 92-104. https://doi.org/10.1016/j.earscirev.2015.03.008.

Walker, J.C., Hays, P.B., Kasting, J.F., 1981. A negative feedback mechanism for the long-term stabilization of Earth's surface temperature. J. Geophys. Res., Oceans 86 (C10), 9776-9782. https://doi.org/10.1029/JC086iC10p09776.

West, A.J., Galy, A., Bickle, M., 2005. Tectonic and climatic controls on silicate weathering. Earth Planet. Sci. Lett. 235 (1-2), 211-228.

Young, S.A., Saltzman, M.R., Foland, K.A., Linder, J.S., Kump, L.R., 2009. A major drop in seawater 87Sr/86Sr during the Middle Ordovician (Darriwilian): links to volcanism and climate? Geology 37 (10), 951-954.