Alexandre Pohl
email: pohl@cerege.fr

David A T Harper
email: david.harper@durham.ac.uk

Yannick Donnadieu
email: donnadieu@cerege.fr

Guillaume Le Hir
email: guillaumelehir[lehir@ipgp.fr

Elise Nardin

Thomas Servais
email: thomas.servais@univ-lille1.fr

Possible patterns of marine primary productivity during the Great Ordovician Biodiversification Event

Following the appearance of numerous animal phyla during the 'Cambrian Explosion', the 'Great Ordovician Biodiversification Event' (GOBE) records their rapid diversification at the lower taxonomic levels, constituting the most significant rise in biodiversity in Earth's history. Recent studies suggest that the rapid rise in phytoplankton diversity observed at the Cambrian-Ordovician boundary may have profoundly restructured marine trophic chains, paving the way for the subsequent flourishing of plankton-feeding groups during the Ordovician. Unfortunately, the fossil record of plankton is incomplete. Its smaller members represent the bulk of the modern marine biomass, but they are usually not documented in Palaeozoic sediments, preventing any definitive assumption with regard to an eventual correlation between biodiversity and biomass at that time. Here we use an up-to-date ocean general circulation model with biogeochemical capabilities (MITgcm) to simulate the spatial patterns of marine primary productivity throughout the Ordovician, and we compare the model output with available palaeontological and sedimentological data.

Introduction

The 'Great Ordovician Biodiversification Event' (GOBE) was arguably the most important and sustained increase of marine biodiversity in Earth's history (e.g. [START_REF] Sepkoski | The Ordovician Radiations: Diversification and extinction shown by global genus level taxonomic data[END_REF][START_REF] Webby | Introduction[END_REF][START_REF] Harper | The Ordovician biodiversification: Setting an agenda for marine life[END_REF]. During the so-called 'Cambrian explosion' most, if not all, animal phyla first appeared in the fossil record. Subsequently, during the Ordovician Period, an 'explosion' of diversity at the order, family, genus, and species level occurred.

The search for the triggers of this biodiversification is ongoing, but most probably there was no unique cause, but the combined effects of several geological and biological processes that helped generate the GOBE.

During the Ordovician Period, a unique palaeogeographical scenario existed, with the greatest continental dispersal of the Palaeozoic. High sea levels that were the highest during the Palaeozoic, if not of the entire Phanerozoic, allowed marine waters to cover large epicontinental areas and flooded large tropical shelf areas, enabling diversification (e.g. [START_REF] Servais | Understanding the Great Ordovician Biodiversification Event (GOBE): Influences of paleogeography, paleoclimate, or paleoecology?[END_REF][START_REF] Servais | The Great Ordovician Biodiversification Event (GOBE): the palaeoecological dimension[END_REF]. The climate was warm, although recent studies indicate a significant, long cooling trend throughout the Ordovician, followed up by an abrupt cooling at the end of the period, triggering the Late Ordovician glacial peak [START_REF] Trotter | Did cooling oceans trigger Ordovician biodiversification? Evidence from conodont thermometry[END_REF][START_REF] Nardin | Modeling the early Paleozoic long-term climatic trend[END_REF], and Harper et al. 2014).

During the Ordovician, important ecological evolutionary changes occurred, beginning with the 'explosion' of the phyto-and zooplankton, leading to the 'Ordovician plankton revolution' [START_REF] Servais | The Ordovician Biodiversification: revolution in the oceanic trophic chain[END_REF]. The onset of the GOBE is actually, at least partly for the planktonic groups, rooted in the late Cambrian, when most of the planktonic organisms started to rapidly diversify [START_REF] Servais | The onset of the 'Ordovician Plankton Revolution' in the late Cambrian[END_REF]. In this context, [START_REF] Saltzman | Pulse of atmospheric oxygen during the late Cambrian[END_REF] already considered that during the latest Cambrian a major increase of atmospheric oxygen concentration (pO 2 ) had already taken place. The temporal correlation between the Cambrian oxygenation event and the concomitant rise in plankton diversity led [START_REF] Saltzman | Pulse of atmospheric oxygen during the late Cambrian[END_REF] and [START_REF] Servais | The onset of the 'Ordovician Plankton Revolution' in the late Cambrian[END_REF] to hypothesize that the increase of pO 2 might have triggered the increase of plankton diversity, which may be related to changes in macro-and micronutrient abundances in increasingly oxic marine environments. The higher amount of available nutrients in the oceans possibly triggered the development of pico-and phytoplankton, i.e. the basis of modern marine trophic chains [START_REF] Saltzman | Pulse of atmospheric oxygen during the late Cambrian[END_REF][START_REF] Servais | The onset of the 'Ordovician Plankton Revolution' in the late Cambrian[END_REF]. Logically, the trophic chain is then assembled with additional tiers provided by the suspension-feeding benthos (e.g., brachiopods, bryozoans and corals) and nekton, themselves prey to a range of predators, including the trilobites, with orthoconic cephalopods and fishes at and near the top of the food chain.

Great advances have been made during the last decades concerning palaeogeographical reconstructions for the Early Palaeozoic, including the Ordovician (e.g. Torsvik & Cocks 2013). These more reliable global palaeogeographical reconstructions have allowed the formulation of simple conceptual and more complex numerical models for ancient atmospheric and ocean circulation patterns, or to hypothetically locate upwelling zones (e.g. [START_REF] Wilde | Oceanography in the Ordovician[END_REF][START_REF] Christiansen | Oceanic circulation as an element in palaeogeographical reconstructions: the Arenig (early Ordovician) as an example[END_REF]Hermann et al. 2004;[START_REF] Pohl | Effect of the Ordovician paleogeography on the (in)stability of the climate[END_REF][START_REF] Servais | Possible oceanic circulation patterns, surface water currents and upwelling zones in the Early Palaeozoic[END_REF]. More recently, new constraints on the palaeobiogeography of marine living communities were provided by the publication of maps showing much more precisely the ocean surface circulation modelled at various atmospheric CO 2 levels during the Early, Middle and Late Ordovician (Pohl et al. 2016b).

The aim of the present paper is to model possible patterns of biomass production in the Ordovician seas, in order to attempt to understand where and how the diversification originated. We use here an ocean-atmosphere general circulation model with biogeochemical capabilities (MITgcm) in order to simulate the changing spatial patterns of marine primary productivity in response to the palaeogeographical evolution throughout the Ordovician. We subsequently attempt to compare the model output with available palaeontological and sedimentological data.

Ordovician biomass and the sedimentary record

Primary productivity in ancient oceans

Many palaeontologists have focussed on the analyses of palaeobiodiversity during the Phanerozoic (e.g., [START_REF] Sepkoski | Phanerozoic marine diversity and the fossil record[END_REF][START_REF] Alroy | The shifting balance of diversity among major marine animal groups[END_REF][START_REF] Harper | The Great Ordovician Biodiversification Event: Reviewing two decades of research on diversity's big bang illustrated by mainly brachiopod data[END_REF]. Such studies lead to the understanding of the major trends in biodiversification during Earth History, and to the discovery of the major extinction phases, including the Big-Five mass extinctions. Specialists in macroevolution and macroecology usually apply various statistical methods in order to better understand and interpret the changing palaeobiodiversities (e.g. [START_REF] Sepkoski | The Ordovician Radiations: Diversification and extinction shown by global genus level taxonomic data[END_REF][START_REF] Bambach | Phanerozoic Biodiversity Mass Extinctions[END_REF][START_REF] Alroy | The shifting balance of diversity among major marine animal groups[END_REF][START_REF] Stanley | Estimates of the magnitudes of major marine mass extinctions in earth history[END_REF]. The GOBE has also been recognized based on these studies.

Biodiversity is a measure of the number of biological organisms present at a given moment, but usually provides no information about the abundance of these organisms. Little attention has been paid to the evolution of the biomass in the oceans (e.g. [START_REF] Franck | Causes and timing of future biosphere extinctions[END_REF][START_REF] Kallmeyer | Global distribution of microbial abundance and biomass in subseafloor sediment[END_REF].

Similarly, the evolution of the abundance of nutrients available in the oceans during the history of the Earth is only poorly known (e.g. [START_REF] Allmon | Seafood through time revisited: the Phanerozoic increas in marine trophic resources and its macroevolutionary consequences[END_REF]). In addition, biodiversity is not (at least not directly) linked to the biomass produced [START_REF] Irigoien | Global biodiversity patterns of marine phytoplankton and zooplankton[END_REF][START_REF] Finnegan | Body size, energetics, and the Ordovician restructuring of marine ecosystems[END_REF]. As a result, little information is available about biomass, or the primary productivity in ancient oceans, including that during the GOBE.

A few authors attempted to understand the evolution of the marine biomass during geological time. [START_REF] Martin | The fossil record of biodiversity: nutrients, productivity, habitat area and differential preservation[END_REF], for instance, analysed the fossil record of biodiversity in relation to nutrients, productivity and habitat area, whereas [START_REF] Martin | Marine biodiversification in response to evolving phytoplankton stoichiometry[END_REF] investigated the evolution of ocean stoichiometry (nutrient content) in order to understand the biodiversification of the Phanerozoic marine biosphere. Concerning the Ordovician, [START_REF] Payne | Controls on marine animal biomass through geological time[END_REF] considered that during the GOBE the increase in the complexity of the marine trophic chains and in the efficiency of marine organisms in removing available food, from both the water column and the sediment, appears to account for a secular increase in animal biomass. Based on [START_REF] Martin | Marine biodiversification in response to evolving phytoplankton stoichiometry[END_REF], [START_REF] Servais | The onset of the 'Ordovician Plankton Revolution' in the late Cambrian[END_REF] further considered that the increasing presence of planktonic organisms in the late Cambrian -Early Ordovician must coincide with increasing nutrient supply, increased primary productivity and expanded biomass production, that resulted during the initiation of the GOBE with a higher diversity and increased abundance of plankton-feeding groups during the Ordovician. [START_REF] Nowak | Phytoplankton dynamics from the Cambrian Explosion to the onset of the Great Ordovician Biodiversification Event: a review of Cambrian acritarch diversity[END_REF] effectively observed a dramatic increase in diversity of the acritarchs, i.e. the organic-walled fraction of the phytoplankton, in the late Cambrian -Early Ordovician interval. But do the higher diversities of phytoplanktonic organisms also indicate an increased abundance of phytoplankton, and an increased biomass?

Ocean productivity largely refers to the production of organic matter by 'phytoplankton' (e.g. [START_REF] Sigman | The Biological Productivity of the Ocean[END_REF]. In summary, most of the single-celled phytoplankton are 'photoautotrophs' that use nutrients and light to convert inorganic to organic carbon. They are subsequently consumed by the 'heterotrophs,' that include the 'zooplankton', the 'benthos,' and the 'nekton.' The most important nutrients necessary for the phytoplankton are nitrogen (N), phosphorus (P), iron (Fe), and silicon (Si), while sunlight is the basic energy source needed. Until recently, it was assumed that the larger parts of the phytoplankton (between 5 and 100 µm in diameter) account for most phytoplankton biomass and productivity. This larger phytoplankton is partly preserved in the fossil record, and corresponds in the Palaeozoic to the informal group of the acritarchs. However, recent studies indicate that more than half of the biomass in modern oceans is actually produced by the much smaller fraction (< 2 µm in diameter), referred to the picoplankton (e.g. [START_REF] Buitenhuis | Picoheterotroph (Bacteria and Archaea) biomass distribution in the global ocean[END_REF]. Recent estimates indicate that 30 % of the modern oceanic biomass is constituted by picoheterotrophs and 25 % of picophytoplankton (e.g., [START_REF] Buitenhuis | MAREDAT: towards a world atlas of MARine Ecosystem DATa[END_REF]Moriarty & O'Brien 2013). This small fraction is almost entirely unknown from the fossil record, the picoplankton not being observed by palaeopalynologists because it falls out of the range of classical observational methods. However, most recent studies indicate that it is this fraction that is the most diverse in modern oceans [START_REF] De Vargas | Eukaryotic plankton diversity in the sunlit ocean[END_REF].

Thus, our understanding of the fossil record of the base of the trophic chain is highly incomplete in the Ordovician, as only a minor fraction of the phytoplankton is recorded under the informal grouping of, for example, the acritarchs. A key approach, therefore, to tentatively understand the biomass production in the Early Palaeozoic is numerical models. In the present study we apply such a model for three different time intervals in the Ordovician to tentatively understand marine productivity during the GOBE interval. Our analysis specifically focuses on large-scale upwelling systems, which generate high rates of biomass production and leave a distinctive fingerprint in the sedimentary record.

Secondary production in upwelling systems

Modern coastal upwelling zones have been a focus of investigation for a number of years as loci of bioproductivity. Here, benthic and nektonic diversity is influenced by a wide range of environmental factors together with low-oxygen concentrations and biotic interactions such as competition and predation. Upwelling zones promote bioproductivity through the delivery of nutrient-rich, deep water onto the shelf, igniting the growth and abundance of phytoplankton. During the GOBE and early stages of the establishment of the Palaeozoic Evolutionary Fauna, communities were dominated by suspension feeders such as the brachiopods, bryozoans and corals and potentially could benefit directly as primary consumers. There is a wealth of biodiversity data available for all the major benthic and nektonic groups (e.g. [START_REF] Webby | Introduction[END_REF]) but there are relatively few regional studies related to specific geographic areas.

Recent studies of some of the key coastal upwelling zones, e.g. along the Namibian Coast (Eisenbrath & Zettler 2016), associated with the Benguela Current Large Marine Ecosystem (BCLME), together with those on coastal upwelling along the Peru coast [START_REF] Rosenberg | Benthos biomass and oxygen deficiency in the upwelling system off Peru[END_REF], the effects of El Nino on the benthos of the Benguela, California and Humboldt upwelling ecosystems [START_REF] Arntz | El Niño and similar perturbation effects on the benthos of the Humboldt, California, and Benguela Current upwelling ecosystems[END_REF] and upwelling along the NW Africa coast [START_REF] Thiel | Zoobenthos of the CINECA area and other upwelling regions[END_REF]) have established some key properties for these zones. Primary and secondary production is substantial, however, it also generates low-oxygen conditions commonly moving the Oxygen Minimum Zone (OMZ) into shallower-water environments.

Biotas are dominated by soft-bodied taxa, there is a reduced diversity and evenness and fewer calcified forms [START_REF] Levin | Oxygen minimum zone benthos: adaptation and community response to hypoxia[END_REF]. In many cases assemblages are dominated by pioneer communities populated by opportunistic species forming dense accumulations. Moreover, the poor oxygen conditions encourage the migration of taxa into shallower water, extending too their geographic ranges along the shelf. Organisms are abundant, generating substantial biomass but not necessarily high diversities.

Methods: model description

Ocean, atmosphere and sea ice

We used a coupled ocean-atmosphere-sea ice setup of the Massachusetts Institute of Technology general circulation model (MITgcm). An isomorphism between ocean and atmosphere dynamics is exploited to allow a single hydrodynamical core to simulate both fluids [START_REF] Marshall | Atmosphere-ocean modeling exploiting fluid isomorphisms[END_REF]). The oceanic and the atmospheric components also share the same cubed-sphere grid with 32 x 32 points per face (cs32), yielding a mean equatorial resolution of 2.8° x 2.8°. The cubed-sphere grid avoids polar singularities resulting from the convergence of the meridians at the poles, thus ensuring that the model dynamics there is treated with as much fidelity as elsewhere [START_REF] Adcroft | Implementation of an atmosphere-ocean general circulation model on the expanded spherical cube[END_REF]).

The oceanic component is an up-to-date, hydrostatic, implicit free-surface, partial step topography ocean general circulation model (Marshall et al. 1997a, b). Twenty-eight levels are defined vertically, the thickness of which gradually increases from 10 m at the surface to 1300 m at the bottom. Effects of mesoscale eddies are parameterised as an advective process [START_REF] Gent | Isopycnal mixing in ocean circulation models[END_REF]) and an isopycnal diffusion [START_REF] Redi | Oceanic isopycnal mixing by coordinate rotation[END_REF]. The nonlocal K-Profile Parameterisation (KPP) scheme of [START_REF] Large | Oceanic vertical mixing: a review and a model with a nonlocal boundary layer parameterization[END_REF] accounts for vertical mixing in the ocean's surface boundary layer, and the interior. The atmospheric physics is based on the Simplified Parameterisations, Primitive-Equation Dynamics (SPEEDY) scheme [START_REF] Molteni | Atmospheric simulations using a GCM with simplified physical parametrizations. I: Model climatology and variability in multi-decadal experiments[END_REF]. The latter comprises a four-band longwave radiation scheme, a parameterisation of moist convection, diagnostic clouds, and a boundary layer scheme. A low vertical resolution is used. Five levels are defined: one level represents the planetary boundary layer, three layers are placed in the troposphere and the fifth layer is placed in the stratosphere. The pressure coordinate p is employed. Sea ice is simulated using a thermodynamic sea-ice model based on the [START_REF] Winton | A reformulated three-layer sea ice model[END_REF] two and a half layer enthalpy-conserving scheme. Sea-ice growth occurs when the ocean temperature falls below the salinity dependent freezing point.

Fluxes of momentum, freshwater, heat, and salt are exchanged every 20 minutes in the model (i.e., the ocean time step). The resulting coupled model can be integrated for ca. 100 years in one day of dedicated computer time. Relatively similar model configurations were used in the past [START_REF] Enderton | Explorations of atmosphere-ocean-ice climates on an aquaplanet and their meridional energy transports[END_REF][START_REF] Ferreira | Localization of deep water formation: role of atmospheric moisture transport and geometrical constraints on ocean circulation[END_REF]2011), including those for palaeoceanographical purposes [START_REF] Brunetti | Modeling the Middle Jurassic ocean circulation[END_REF][START_REF] Pohl | The climatic significance of Late Ordovicianearly Silurian black shales[END_REF].

Primary productivity

The MITgcm includes a biogeochemistry model that simulates the net primary productivity (NPP) in the ocean. NPP is computed based on Michaelis-Menten equations as a function of available photosynthetically active radiation (PAR) and phosphate concentration (PO 4 ), where α = 2 × 10 -3 mol m -3 yr -1 is the maximum community productivity, K PAR = 30 W m -2 the half saturation light constant, and K PO4 = 5 × 10 -4 mol m -3 the half saturation phosphate constant. In this configuration, phosphate is the single limiting nutrient. It is consumed in the photic zone to fuel the marine primary productivity, regenerated by remineralisation throughout the water column based on the empirical law of [START_REF] Martin | VERTEX: carbon cycling in the northeast Pacific[END_REF], redistributed within the ocean using the velocity and diffusivity fields provided by the general circulation model and ultimately returned back to the ocean surface in upwelling zones. Because phosphate is assumed to have an oceanic residence time much longer than the oceanic turnover time scale (i.e., 10-40 kyr; [START_REF] Ruttenberg | Reassessment of the oceanic residence time of phosphorus[END_REF][START_REF] Wallmann | Feedbacks between oceanic redox states and marine productivity: A model perspective focused on benthic phosphorus cycling[END_REF], its global oceanic concentration is fixed in the model [START_REF] Dutkiewicz | Interactions of the iron and phosphorus cycles: A three-dimensional model study[END_REF]. Iron is known as another major factor in productivity [START_REF] Falkowski | Ocean science: the power of plankton[END_REF]). Because it is mainly delivered to the ocean surface as dust from deserts, the emissions of which are difficult to quantify today [START_REF] Bryant | Recent advances in our understanding of dust source emission processes[END_REF]); providing the model with seasonal maps of iron input in the Ordovician is challenging. Although our biogeochemical model has the provision to account for cycling of iron, we prefer not to consider iron fertilisation here. The PAR is computed at the ocean surface as a fraction of the incident shortwave radiation provided by the atmospheric component of the MITgcm. It is then attenuated throughout the water column assuming a uniform extinction coefficient. Similar configurations of this biogeochemistry model have been used in the past [START_REF] Friis | Possible overestimation of shallow-depth calcium carbonate dissolution in the ocean[END_REF]2007), including those for the Ordovician [START_REF] Pohl | The climatic significance of Late Ordovicianearly Silurian black shales[END_REF].

Boundary and initial conditions

We ran our model on three palaeogeographical reconstructions representative of the Early Ordovician (480 Ma, Tremadocian), the Middle Ordovician (460 Ma, Darriwilian) and the early Silurian (440 Ma, Aeronian). The location of the continental masses is taken from the reconstructions by [START_REF] Torsvik | BugPlates: linking biogeography and palaeogeography, software manual[END_REF]. The topography and the bathymetry are reconstructed based on published global reconstructions, with additional information from regional studies for Gondwana Siberia (e.g. [START_REF] Cocks | Siberia, the wandering northern terrane, and its changing geography through the Palaeozoic[END_REF] and Asia (e.g. Cocks & Torsvik 2013). Because the location and depth of Ordovician ocean ridges is not well constrained, they are not included in the model. We use a flat-bottom ocean, the depth of which is set to present-day mean seafloor depth, i.e., -4000 m [START_REF] Pohl | Effect of the Ordovician paleogeography on the (in)stability of the climate[END_REF]. The flat bottom is not expected to constitute a major bias. Several studies on late Palaeozoic oceans suggest that it does not critically impact the large-scale patterns of simulated ocean circulation [START_REF] Montenegro | Climate simulations of the Permian-Triassic boundary: Ocean acidification and the extinction event[END_REF][START_REF] Osen | Sensitivity of Late Permian climate to bathymetric features and implications for the mass extinction[END_REF]. The resulting palaeomaps are very similar to those used by Pohl et al. (Pohl et al. 2016b).

A gradual greening of the continents occurred throughout the Ordovician. The first, non-vascular land plants are documented from the Middle Ordovician Dapingian [START_REF] Rubinstein | Early Middle Ordovician evidence for land plants in Argentina (eastern Gondwana)[END_REF]. There is no evidence of plants on land before that date, including during the Tremadocian, i.e., the first time slice used in the present study. In addition, the spatial cover of this primitive vegetation is difficult to estimate for the remainder of the period [START_REF] Edwards | Could land-based early photosynthesizing ecosystems have bioengineered the planet in mid-Palaeozoic times?[END_REF][START_REF] Porada | High potential for weathering and climate effects of non-vascular vegetation in the Late Ordovician[END_REF]. As a consequence we here follow previous studies [START_REF] Nardin | Modeling the early Paleozoic long-term climatic trend[END_REF][START_REF] Pohl | Effect of the Ordovician paleogeography on the (in)stability of the climate[END_REF]) and impose a rocky desert landscape on the continents (ground albedo of 0.24, which is potentially modified by snow).

During the Ordovician, the atmospheric partial pressure of CO 2 (pCO 2 ) was significantly higher than today (e.g. [START_REF] Berner | GEOCARBSULF: A combined model for Phanerozoic atmospheric O 2 and CO 2[END_REF]). However, the MITgcm does not account for varying pCO 2 levels. It is tuned to the present-day pCO 2 . We therefore increased the solar forcing in the model to simulate various climatic states [START_REF] Ferreira | Climate determinism revisited: multiple equilibria in a complex climate model[END_REF]) and subsequently compared the simulated temperatures with Ordovician estimates. A solar constant of 350 W m -2 (instead of 342 W m -2 today) induces an increase in tropical sea-surface temperatures (SSTs) up to 32.5 °C to 33.7 °C (depending on the time slice considered). These values compare well with the SSTs reconstructed by [START_REF] Trotter | Did cooling oceans trigger Ordovician biodiversification? Evidence from conodont thermometry[END_REF] for the Middle Ordovician based on ∂ 18 O measurements (see their Fig. 3 in particular). Alternative values would better fit the Early and the Late Ordovician SSTs, but we here aim at quantifying the impact of the palaeogeographical changes on the spatial patterns of Ordovician primary productivity, all other things kept equal. We therefore conduct our three simulations using a solar constant of 350 W m -2 , which are here selected to provide the best match with Middle Ordovician SST estimates. The orbital configuration is defined with an obliquity of 23.45° and an eccentricity of 0°.

We use identical initial conditions in all simulations, including an homogeneous salinity of 35 psu (practical salinity units) and a theoretical latitudinal gradient of ocean temperature characterised by equatorial and polar SSTs of respectively 35 °C and 6 °C and an ocean bottom potential temperature of 3 °C. These values ensure that the ocean is sea ice free at the beginning of each model run.

Phosphate is initialised with its present-day depth profile. For each simulation, the physical oceanatmosphere-sea ice model is first run until deep-ocean equilibrium is reached (≥ 2000 model years).

It is subsequently restarted with marine biogeochemistry for 550 additional years and climatic fields used for analysis are averaged over the last 50 years of the simulation.

Results: simulated Ordovician marine productivity

In space: wind belts and Ekman pumping

The main patterns of simulated Ordovician surface primary productivity (Fig. 1) reflect in simple terms the ocean phosphate concentration, while the PAR only imposes a hemispheric-scale decrease of NPP with latitude. The phosphate concentration in surface seawaters results, in turn, from the large-scale patterns of Ekman pumping. The concentration in PO 4 is lower in shallow waters because nutrients are consumed during photosynthesis and it increases with depth due to remineralisation of sinking particles. Upwelling systems allow nutrient-rich deep waters to be transported back to the surface and they are therefore associated with locally high phosphate content. On the contrary, down welling areas are poor in nutrients. The spatial patterns of upwelling and down welling are essentially driven by the direction of the wind blowing over the ocean surface through Ekman pumping and suction. The Trade Winds induce large-scale upwelling systems on the western margin of tropical landmasses, and the Westerlies cause high phosphate concentrations at the mid-latitudes (40° -60°). Between the Trade Winds and the Westerlies, the down welling of surface waters along the tropics (30°) leads to low phosphate levels. A local minimum in phosphate concentration and thus primary productivity occurs in the Palaeo-Tethys, between the eastern coasts of Baltica and the western margin of tropical Gondwana (Fig. 1A-C). Here, the down welling of surface waters combines with a strong freshwater input from the continent. The latter results from the intense orographic precipitation that occurs when the moisture-laden Westerlies intercept the coastal topography of Gondwana. [START_REF] Pohl | The climatic significance of Late Ordovicianearly Silurian black shales[END_REF] demonstrated that this strong runoff to the ocean is a robust model result that it not overly model-dependent. While the polar latitudes (60° -90°) are dynamically isolated from the global ocean and thus depleted in phosphate in the Northern Hemisphere [START_REF] Pohl | The climatic significance of Late Ordovicianearly Silurian black shales[END_REF], deep-water convection along the coast of Gondwana (e.g. [START_REF] Poussart | Late Ordovician glaciation under high atmospheric CO 2 : A coupled model analysis[END_REF][START_REF] Herrmann | Response of Late Ordovician paleoceanography to changes in sea level, continental drift, and atmospheric pCO 2 : potential causes for long-term cooling and glaciation[END_REF]) drives the local enrichment of surface waters in nutrients, ensuring relatively high productivity levels there.

Although the small-scale spatial patterns of simulated NPP may be model-dependent to some extent, the fact that the first-order signal results from a fundamental characteristic of Earth's climate (the zonal wind belts) suggests that these results are relatively robust [START_REF] Pohl | The climatic significance of Late Ordovicianearly Silurian black shales[END_REF].

In time: throughout the Ordovician

The three studied time slices share a certain number of common features (Fig. 1). The western margin of Laurentia, first, is associated with high levels of simulated NPP. In both hemispheres, the mid-latitudes are further characterised by zonal currents inducing Ekman pumping and the upwelling of nutrients fuelling a strong productivity. On the contrary, low marine productivity levels typify the tropics (30°), the margin of Gondwana situated over the South Pole and the Northern high-latitudes.

A local minimum persists throughout the period at 30° S between the western coast of Gondwana and the tropical landmasses.

Nevertheless, major changes can be observed from 480 Ma to 440 Ma (Fig. 1). The most considerable alteration of the NPP patterns between the Early and the Late Ordovician resulted from the gradual drift of Gondwana to the North. Confined in the Southern Hemisphere at 480 Ma, the supercontinent reached 30° N at 440 Ma. The direct consequence of this continental drift was the appearance of a major upwelling system at tropical latitudes in the Late Ordovician along the coasts of Gondwana (Australia and India) and South China-Annamia. The contrast between the eastern and western tropical coasts of Gondwana increased at the same time, with simulated NPP significantly decreasing along the eastern coast of Gondwana at 440 Ma. Elsewhere, the spatial patterns of simulated NPP remained relatively stable as the tropical continental masses slowly migrated to the North. This continental drift prompted Siberia to slowly shift, from the Early to the Late Ordovician, from the zone of minimum NPP reported previously at 30° S to the highly productive Panthalassic circumpolar current. On the contrary, Baltica underwent a migration from nutrient-rich mid-latitude waters in the Southern Hemisphere at 480 Ma to the depleted water masses centred on 30° S at 440 Ma.

Discussion

Here we compare our modelling results with the Ordovician geological record in order to validate our simulations. We also discuss some enigmatic biotic events in the light of our model runs, such as the sudden and widespread bioherm development that punctuated the Late Ordovician of Baltoscandia.

Geological evidence for upwelling systems

During the Ordovician two key coastal areas are typified by high levels of primary productivity in the model: the west coasts of Laurentia and tropical Gondwana (Fig. 1), in particular Australia and South China. These modelling results are supported by geological data. Indeed, a wide range of sedimentary indicators has been confirmed for a prolonged phase of upwelling across much of Laurentia during the Middle and Late Ordovician related to glaciation [START_REF] Pope | Widespread, prolonged late Middle to Late Ordovician upwelling in North America: A proxy record of glaciation?[END_REF]. A careful reassessment of the ages of many of the units [START_REF] Leslie | Widespread, prolonged late Middle to Late Ordovician upwelling in North America: A proxy record of glaciation?: Comment and Reply[END_REF] indicates that the formation of cherts may be much more widespread, suggesting the possibility that this is related to more general phases of upwelling across the Laurentian continent. Associated biotic indicators are sparse and under-developed. Graptolites, however, were most common in upwelling zones along continental margins. In the Vinini Formation, Roberts Mountains, Nevada, for example, changes within the graptolite zooplankton have been associated with fluctuating oceanographic conditions and the upwelling of anoxic waters. During high-stands graptolite faunas diversified within the OMZ but the ecosystem collapsed with the substantial fall of sea level associated with the end Ordovician glaciation and the retreat of the OMZ [START_REF] Finney | The influence of denitrifying seawater on graptolite extinction and diversification during the Hirnantian (latest Ordovician) mass extinction event[END_REF]. Elsewhere in the Great Basin the macroshelly fauna is abundant throughout much of the Ordovician, forming locally shell beds and concentrations [START_REF] Finnegan | Relative and absolute abundance of trilobites and rhynchonelliform brachiopods across the Lower/Middle Ordovician Boundary, Eastern Basin and Range[END_REF], with individuals increasing in size [START_REF] Payne | Controls on marine animal biomass through geological time[END_REF] and shell thickness [START_REF] Pruss | Carbonates in skeleton-poor seas: new insights from Cambrian and Ordovician strata of Laurentia[END_REF]; the faunas are never highly diverse but they are abundant.

Bryozoan-rich deposits, formed by another suspension feeder, have also been associated with the upwelling of nutrients in parts of Laurentia during the later Ordovician [START_REF] Taylor | Latitudinal distribution of bryozoan-rich sediments in the Ordovician[END_REF].

In detail the location of the large-scale upwelling systems simulated in the model is in relatively good agreement with most of the evidence of upwelling (cherts and phosphate deposits) documented by [START_REF] Pope | Widespread, prolonged late Middle to Late Ordovician upwelling in North America: A proxy record of glaciation?[END_REF] (Fig. 1A). Best match is observed on the western margin of the tropical landmasses, i.e., Laurentia and Gondwana. Some data points are more difficult to explain (e.g., in southeastern Laurentia), and the most outstanding model-data mismatch is observed in Baltica.

Clearly, our simulation provides no explanation for the preservation of cherts in that precise location.

In order to explain the discrepancy, we emphasize that the spatial resolution of our model (ca. 300 km) does not allow us to capture small-scale processes. In addition, the land-sea mask interpolated on the model's grid constitutes a crude approximation of Ordovician coastlines. Finally, current palaeogeographical reconstructions only provide first-order indications of Early Palaeozoic bathymetry and topography and they generate large uncertainties in the position of the continental masses (up to 15° in latitude) [START_REF] Lees | Quantifying paleogeography using biogeography: a test case for the Ordovician and Silurian of Avalonia based on brachiopods and trilobites[END_REF], which makes any straightforward model-data comparison challenging.

To explain the appearance of these unusual carbonates with abundant chert and phosphate, [START_REF] Pope | Widespread, prolonged late Middle to Late Ordovician upwelling in North America: A proxy record of glaciation?[END_REF] required the strengthening of the meridional overturning circulation and associated increased nutrient supply to the ocean surface in response to glacial onset during the late Middle Ordovician. Because similarly high levels of primary productivity are simulated on the western margin of Laurentia during each of the three studied time slices (Fig. 1A-C), with no ice sheet over the South Pole, our model suggests that such coastal upwelling may instead have been a persistent characteristic of Ordovician oceans. These results are supported by the reappraisal of the age of many of the deposits originally reported by [START_REF] Pope | Widespread, prolonged late Middle to Late Ordovician upwelling in North America: A proxy record of glaciation?[END_REF], a number of them being potentially Early Ordovician in age [START_REF] Leslie | Widespread, prolonged late Middle to Late Ordovician upwelling in North America: A proxy record of glaciation?: Comment and Reply[END_REF]. Together this raises serious doubts about the climatic implications of the phosphatic rocks proposed by [START_REF] Pope | Widespread, prolonged late Middle to Late Ordovician upwelling in North America: A proxy record of glaciation?[END_REF]. Glacial onset does not seem to be a necessary condition for strong upwelling to occur along the margin of Laurentia.

Another possible indicator of high levels of primary productivity is the preservation of sediments enriched in organic matter. [START_REF] Melchin | Environmental changes in the Late Ordovician-early Silurian: Review and new insights from black shales and nitrogen isotopes[END_REF] published a compilation of Late Ordovician -early Silurian black shale occurrence. They illustrated three time slices immediately before, during and right after the Hirnantian glacial peak. Using simulations similar to ours, but focusing on the Late Ordovician, [START_REF] Pohl | The climatic significance of Late Ordovicianearly Silurian black shales[END_REF] recently demonstrated a striking correlation between the regions of high (low) primary productivity simulated in their model and the preservation of organic-rich (organic-poor) sediments in the geological record of the late Katian (i.e., the first time slice of [START_REF] Melchin | Environmental changes in the Late Ordovician-early Silurian: Review and new insights from black shales and nitrogen isotopes[END_REF]. Similar to the cherts and phosphate deposits of [START_REF] Pope | Widespread, prolonged late Middle to Late Ordovician upwelling in North America: A proxy record of glaciation?[END_REF], the black shales are documented on the western margins of equatorial Laurentia and equatorial Gondwana, thus matching simulated upwelling systems (Fig. 1A). More interestingly, the deposits depleted in organic matter are found around Baltica and along the coast of Gondwana over the South Pole [START_REF] Melchin | Environmental changes in the Late Ordovician-early Silurian: Review and new insights from black shales and nitrogen isotopes[END_REF], precisely where the model simulates local NPP minima (Fig. 1A). This modeldata agreement supports the spatial patterns of NPP simulated in the present study. It also suggests that the preservation of black shales in the late Katian may have been driven by the patterns of primary production at the ocean surface [START_REF] Pohl | The climatic significance of Late Ordovicianearly Silurian black shales[END_REF].

Baltoscandian reefs and mounds

In their comprehensive review, [START_REF] Kröger | Ordovician reef and mound evolution: the Baltoscandian picture[END_REF] They postulated climatic and eustatic drivers and further identified the northward drift of Baltica as the main determining factor for the timing of the start of the Baltic reef and mound development.

More specifically, they suggested that 'the entry of Baltica in a geographical zone that allowed for a widespread bioherm formation was a major factor for the Sandbian radiation of reefal and reefrelated organisms'. Although the most direct effect of the latitudinal shift of Baltica was probably the establishment of climatic conditions more favourable to the shallow-water carbonate factory, our model runs also indicate that this northward migration made Baltica enter the nutrient-depleted tropical region during the Late Ordovician (Fig. 1A). Such unparalleled oligotrophic conditions may have provided the appropriate environmental background for the explosion of the bioherms on Baltoscandia.

Limitations

Our simulations provide an overview of possible spatial patterns of Ordovician NPP in time and space. However, it does not tell us anything about the total amount of NPP. The latter is a relatively direct function of nutrient availability in the ocean, which depends in turn on the intensity of continental weathering. Estimating continental weathering requires major assumptions on both the lithology and vegetation cover of emerged continental masses and global climate (e.g. [START_REF] Goddéris | The role of palaeogeography in the Phanerozoic history of atmospheric CO 2 and climate[END_REF]. In particular, weathering may have significantly varied throughout the Ordovician as a result of changes in global climate [START_REF] Trotter | Did cooling oceans trigger Ordovician biodiversification? Evidence from conodont thermometry[END_REF], continental configuration [START_REF] Nardin | Modeling the early Paleozoic long-term climatic trend[END_REF], volcanic activity [START_REF] Lefebvre | Did a Katian large igneous province trigger the Late Ordovician glaciation? A hypothesis tested with a carbon cycle model[END_REF]), land-ice cover (Pogge von Strandmann et al. 2017) or in response to the advent of the first land plants [START_REF] Lenton | First plants cooled the Ordovician[END_REF][START_REF] Stanley | Estimates of the magnitudes of major marine mass extinctions in earth history[END_REF][START_REF] Porada | High potential for weathering and climate effects of non-vascular vegetation in the Late Ordovician[END_REF]. Such mechanisms do lie beyond what we are able to resolve using our ocean-atmosphere model. As a result, we are unable to predict the trend towards an increase or a decrease in NPP throughout the period.

Conclusion

Numerical simulations conducted with an ocean-atmosphere general circulation model with biogeochemical capabilities (MITgcm) predict the position and longevity of upwelling zones in the Ordovician Earth system. Upwelling zones host specific types of ecosystems, characterised by high organic productivity, abundant organisms (commonly opportunist species, often soft-bodied) but not necessarily high diversities. Nevertheless, high productivity levels may have provided the primary resources needed by superior consumers and thus paved the way for the development of more complex trophic chains featuring more diverse taxonomic assemblages. Fossil data have, to date, not been aligned with these predictions although accumulations of organic matter, bone and skeletal concentrations are identifiable in the fossil record and can be correlated with characteristic sediments, such as cherts, that signpost upwelling zones in deep time. The present study therefore targets in-depth analysis of the smaller members of the fossil record as the logical next step towards the integrated understanding of the diversification patterns throughout the early Palaeozoic.

  (e.g. Torsvik & Cocks 2013), Laurentia (e.g. Cocks & Torsvik 2011), Baltica (e.g. Cocks & Torsvik 2005),

  demonstrated a critical change in the mode of carbonate production across the Baltic Basin during the latest Sandbian -earliest Katian (Late Ordovician, ca. 453 Ma). This boundary marked the beginning of a protracted period of widespread development of reefs in shallow-water areas and mud mounds in deeper epicontinental settings. The authors showed that biotic factors do not explain the initiation of Late Ordovician bioherm growth.

Figure 1 :

 1 Figure 1: Surface primary productivity simulated at 440 Ma (A), 460 Ma (B) and 480 Ma (C). The black mask indicates continental masses. In subfigure A, white points represent

Acknowledgments

The findings in this study are based on climatic fields simulated by the ocean-atmosphere general circulation model MITgcm. Code for the climate model MITgcm can be accessed at http://mitgcm.org.

Requests for the climate model output can be sent to A.P (pohl@cerege.fr). The authors thank Peter Doyle and Alan Owen for editorial handling. We also thank Tom Challands and Axel Munnecke for helpful and constructive reviews, and the guest-editors of this volume for the invitation to submit. The authors acknowledge the financial support from the CNRS (INSU, action SYSTER). A.P. and Y.D. thank the CEA/CCRT for providing access to the HPC resources of TGCC under the allocation 2014-012212 made by GENCI. A.P. thanks Maura Brunetti from the University of Geneva and David Ferreira from the University of Reading, who provided expertise that greatly assisted the research. This research was funded through a CEA PhD grant CFR. This is a contribution to the IGCP Project-653, 'The onset of the Great Ordovician Biodiversification Event'. D.A.T.H. acknowledges financial support from the Leverhulme Trust and the Wenner Gren Foundation.

evidence of upwelling. The latter are plotted on the Late Ordovician time slice following [START_REF] Pope | Widespread, prolonged late Middle to Late Ordovician upwelling in North America: A proxy record of glaciation?[END_REF], but see discussion in the main text. L: Laurentia; B: Baltica; S: Siberia; G: Gondwana; Sc: South China; An: Annamia; Aus: Australia.