Climate of the Past Exploring the climatic impact of the continental vegetation on the Mezosoic atmospheric CO 2 and climate history

In this contribution, we continue our exploration of the factors defining the Mesozoic climatic history. We improve the Earth system model GEOCLIM designed for long term climate and geochemical reconstructions by adding the explicit calculation of the biome dynamics using the LPJ model. The coupled GEOCLIM-LPJ model thus allows the simultaneous calculation of the climate with a 2-D spatial resolution, the coeval atmospheric CO 2, and the continental biome distribution. We found that accounting for the climatic role of the continental vegetation dynamics (albedo change, water cycle and surface roughness modulations) strongly affects the reconstructed geological climate. Indeed the calculated partial pressure of atmospheric CO 2 over the Mesozoic is twice the value calculated when assuming a uniform constant vegetation. This increase in CO 2 is triggered by a global cooling of the continents, itself triggered by a general increase in continental albedo owing to the development of desertic surfaces. This cooling reduces the CO 2 consumption through silicate weathering, and hence results in a compensating increase in the atmospheric CO 2 pressure. This study demonstrates that the impact of land plants on climate and hence on atmospheric CO 2 is as important as their geochemical effect through the enhancement of chemical weathering of the continental surface. Our GEOCLIM-LPJ simulations also define a climatic baseline for the Mesozoic, around which exceptionally cool and warm events can be identified.


Introduction
Since the early work of Walker et al. (1981), global climatic evolution of the Earth is described as the result of the long term balance between continental silicate weathering (Fsw) Correspondence to: Y. Donnadieu (yannick.donnadieu@lsce.ipsl.fr)and solid Earth CO 2 degassing (Fvol).This balance defines the multi-million years atmospheric CO 2 base level, around which fluctuations are allowed at shorter timescales (Berner and Kothaval, 2001;Donnadieu et al., 2006a;Goddéris and Franc ¸ois, 1995).The term Fsw is expressed as a function of numerous parameters, such as temperature (Brady, 1991;Dessert et al., 2001), continental runoff (Oliva et al., 2003), physical erosion (Millot et al., 2002) and continental vegetation cover (Berner, 1997), these parameters being all intertwined.Nevertheless, a comprehensive formula, physically accounting for each parameter, is far from being reached.Indeed, the most standard mathematical expression for silicate weathering is: where runoff and T are respectively the mean annual continental runoff and temperature.k is a constant.Then, in order to introduce the others factors influencing Fsw, one can add several dimensionless parameters or feedback functions representing the erosion changes driven by the uplift occurring in Earth's history, the evolution of the continental area modulated by the sea-level changes, the impact of continental vegetation, and so on (Berner, 2004).An interesting example of this method is the simulation of the effect of land plants on weathering.Land plants can activate the silicate weathering rate.With their roots, they make the soil more acidic to enhance the uptake of nutrients.They also increase the pCO 2 of the soil since there is more organic matter, and CO 2 is released by the root and microbial respiration.They finally increase the residence time of water in the weathering profile, thus promoting further mineral dissolution.Studies were made to quantify the impact of plants by measuring the silicate weathering rate over vegetated area versus area sporadically covered by mosses and lichens (Drever and Zobrist, 1992;Moulton and Berner, 1998).Those authors found a ratio of weathering rates between vegetated and non-vegetated area varying between 2 and 8 for young weathering systems Published by Copernicus Publications on behalf of the European Geosciences Union.
Y. Donnadieu et al.: Biophysical effect of vegetation on Mesozoic atmospheric CO 2 (with intense physical erosion such as in mountainous area, or on recent volcanic outcrop as in Iceland).
Based on these results, Berner (1990) introduced a multiplying parameter in Eq. (1) (f E(t), where t is geological time) that expresses the effects of the rise and evolution of vascular land plants on the rates of silicate weathering.The lowest value of f E(t) is 0.25 and is used for a non vegetated world (before 380 Ma).Then, a value of f E(t)=0.875 for weathering by pteridophytes and gymnosperms (350-130 Ma) is used with a transition to angiosperm controlled weathering taking place between 130 and 80 Ma (Berner, 2004).After 80 Ma, the value of f E(t) is set equal to one.Based on such parameterization, Berner (2004) has demonstrated that the evolution of land plants is a first-order control parameter of the long term climate evolution through its impact on silicate weathering rates.In this contribution, we will explore an effect of the land plants on the Earth climate evolution neglected up to now: its physical impact on the global water cycle, on the albedo of the Earth surface, and on the continental roughness.All those parameters should directly influence the continental temperature and runoff, and hence the continental weathering, in addition to the direct effect of land plants on weathering described above.These climatic impacts have been neglected up to now because the climate model used in long-term geochemical and climatic evolution were too basic to allow a physical quantification.Indeed, previous climatic-geochemical coupled models did not explicitly and physically include the climate as they either use a parameterization linking the mean global temperature to atmospheric CO 2 , and global continental runoff to this calculated air temperature (Berner, 1990), or 1-dimensional energy balance models with no physical water cycle representation (Godderis and Joachimski, 2004;Simon et al., 2007).However spectacular progresses have been made in the climate modelling area.In particular, owing to the ongoing improvement of supercomputers, we are now in a position to use General Circulation Models for performing such a study with an interactive vegetation model.Such types of models would allow a physical quantification of the effect of land plants on the climate.
We recently developed an Earth System Model, called GEOCLIM, capable of multi-million year integrations of long-term variables (i.e.oceans, ice sheets, geochemical cycles) accounting for relatively low frequency forcing (e.g.paleogeography, carbonate factory) and internal feedbacks contributing to non-linear behaviour (Donnadieu et al., 2004).This model coupled an atmospheric general circulation model to a global biogeochemical cycle model (including a description of the carbon and alkalinity cycles).This model has been used to explicitly quantify the impact of the continental configurations on the Mesozoic climate history (251-65 Ma).We showed that the break-up of the Pangaea supercontinent has triggered an increase in continental runoff, resulting in enhanced atmospheric CO 2 consumption through silicate weathering (Donnadieu et al., 2006a).As a result, at-mospheric CO 2 fell from values above 3000 ppmv during the Triassic to rather low levels during the Cretaceous (around 400 ppmv), resulting in a decrease in global mean annual continental temperatures from about 20 • C to 10 • C over the course of the Mesozoic.Most of the CO 2 drop is concentrated at the end of the Triassic within a 20 Myr time window in response to a moderate northward drift of the Pangea supercontinent (Godderis et al., 2008), exposing large continental surface to intense weathering in the warm and humid equatorial climate.
All those GEOCLIM simulations were performed assuming uniform physical characteristics of the continental surface, typical of a globally distributed deciduous forest.Because of this rough assumption, the climatic feedback of the biome distribution was not explored.After having unravelled the tectonic forcing over the Mesozoic, we focus on the impact of vegetation on climate and hence on weathering rates.This objective requires the estimate of the biome distribution in response to climate change, and the feedback of this distribution on the climate itself.For that reason, we have coupled the Dynamic Vegetation Model LPJ (Sitch et al., 2003) to the Fast Ocean Atmosphere Model FOAM (Donnadieu et al., 2006b).We thus quantify the feedback of changing vegetation cover on the climate, which in turn influences the silicate weathering through the runoff and the temperature variables used to calculate weathering rates (Eq.1).New predicted atmospheric CO 2 and climate are discussed regarding the most recent 18 O constraints.

Generality
We use the GEOCLIM model allowing to calculate simultaneously the global carbon cycle (including the atmospheric CO 2 partial pressure), the climate and the silicate rock weathering, the two last output with a 2-dimensional spatial resolution above continental surfaces.A full description of the model can be found in Donnadieu et al. (2006) (see also references therein).In summary, a full coupling between FOAM and COMBINE (the geochemical module) cannot be achieved owing to excessive computational times.Hence we adopt an indirect coupling that employs look-up tables from a catalogue of simulations.Climatic parameters, temperature and runoff, calculated on the 48×40 grid of FOAM are used to evaluate the weathering rates on the same grid according to Eq. (1).Then, fixing the CO 2 degassing rate, a numerical feedback loop between FOAM and COMBINE is run until a steady-state PCO 2 is reached.Each continental configuration is thus finally characterized by a steady-state atmospheric CO 2 level and a distribution of the climatic variables (mean annual temperature and runoff).For the climatic simulations, the GCM FOAM is used.The atmospheric component of FOAM is a parallelized version of NCAR's Community Climate Model 2 (CCM2) with the upgraded radiative and hydrologic physics incorporated in CCM3 v. 3.2.The atmosphere runs at R15 spectral resolution (4.5 • ×7.5 • ) with 18 levels.The ocean component (OM3) is a z-coordinate ocean model that has been optimized for performance and scalability on parallel computers.OM3 contains 24 vertical layers, a 128×128 grid (1.4 • ×2.8 • ) and uses simple second order differencing and a fully explicit time step scheme for the barotropic and baroclinic modes.FOAM successfully simulates many aspects of the present-day climate and compares well with other contemporary medium-resolution climate models; it has also been used previously to investigate Cretaceous and Neoproterozoic climates (Donnadieu et al., 2006b;Poulsen et al., 2003Poulsen et al., , 2002Poulsen et al., , 2001)).For this study, we use FOAM in mixed-layer mode, i.e., the atmospheric model is linked to a 50-m mixed-layer ocean, which parameterizes heat transport through diffusion, mainly for computation time considerations (each GEOCLIM simulation requires up to 12 GCM simulations, as described below).Still for climatic simulations with FOAM, the Earth's orbit around the Sun is circular (eccentricity = 0) and the Earth's obliquity is 23.5 • (this setting leads to an equal annual insolation for both hemispheres).Solar luminosity is assumed to evolve through time according to the stellar evolution models (from 2.5%reduction in the Early Triassic to 0.6% in the Maastrichtian (Gough, 1981)).

Vegetation
The main improvement of GEOCLIM is the coupling with the vegetation model LPJ (Lund-Potsdam-Jena) so that vegetation distribution can be coevally calculated with the climatic reconstruction.LPJ is a Dynamic Global Vegetation Model of vegetation biogeography and vegetation/soil biogeochemistry.Starting from the climate, soil and atmospheric information as input, it dynamically computes spatially explicit transient vegetation composition in terms of plant functional groups, and their associated carbon and water budgets.To couple FOAM with LPJ, we have written a python program that changes the climatic parameters in the FOAM boundary files according to the vegetation given by LPJ (python is a dynamic object-oriented programming language).The influence of vegetation on climate occurs through the fluxes of energy, water vapour and momentum between land and atmosphere.Vegetation absorbs part of the heat which modifies the sensible heat flux.Plants also transpire and change the latent heat flux.Furthermore, they hinder air movements which transfer momentum.The calculated vegetation distribution directly impacts on the following physical parameters in FOAM: (1) the albedo of each continental grid element, (2) its thermal properties, (3) its roughness, (4) its potential evaporation.

Approximations for deep time utility
Even if LPJ has been designed to be used for present and past simulations, some of the plant functional types (PFTs) are not relevant for the Mesozoic simulations.Indeed, C4 grass were only well established by the Miocene period, and C3 grass probably appeared during the Late Cretaceous (Chapman, 1993).Instead of grass, ferns would probably have grown.Thus, the temperate herbaceous (C3 grass) and the tropical herbaceous (C4 grass) must be understood as Pteridophytes.Nevertheless, the creation of a new PFT corresponding to ferns, which remains to be done, would allow us to have a more precise idea of where ferns would have developed.This study should really be understood as a first step toward a better understanding of the climatic role of land plants on weathering.Indeed, in addition to the uncertainties tied to the use of modern PFTs for deep time problems, representation of the explicit effect of plants on weathering in the soil will not be investigated here as mentioned in the introduction.Our goal here is to isolate the climatic role of land plants.Future advances on coupling LPJ with an explicit soil weathering model such as WITCH (Godderis et al., 2006) will be undertaken in the future, so that the impact of continental vegetation on chemical weathering could be explicitly modelled.
A single soil type (medium) is prescribed in LPJ because adequate empirical constraints on the types of soil during the Mesozoic are lacking.A sensitivity analysis on the type of soil prescribed in LPJ for the present has been done; comparing the vegetation distribution obtained with a single soil type to the results obtained with the real soil types.It appears that it only changes the vegetation distribution by 12% in the worst case (with the fine, vertisol soil type).As the vegetation changes are less than 2% with the medium soil type, the latter is chosen as the single soil type used for the Mesozoic simulations.

Boundary conditions
The two only forcing functions are the slowly increasing solar constant and the continental configuration.The carbon cycle, the climate and the continental vegetation evolve freely until a steady-state is reached.Since no firm constraint exists on the solid Earth degassing rate, we kept it at its present-day constant value in all simulations (see the brief discussion in Donnadieu et al., 2006).Once steadystate is reached, the total CO 2 consumption through silicate weathering is thus the same in all simulations to maintain the multi-million year balance of the carbon cycle.But the spatial distribution will differ from one geological period to the other.Seven time slices are simulated, corresponding to seven geographical configurations along the course of the Pangea breakup.These configurations are exactly the same as those used in Donnadieu et al. (2006) and in Godderis et al. (2008), from the end Permian towards the latest  Cretaceous.The only important change is the coupling of FOAM with the terrestrial biosphere model, LPJ.In our coupling of FOAM and LPJ, FOAM is run for one year, during which monthly average temperature, precipitation, cloud cover and insolation are saved.These fields are then used to force LPJ, after which the vegetation is updated in FOAM.This process is repeated throughout the experiment integration defined by geography and an atmospheric CO 2 level.Each experiment is integrated for 30 model years; results used in the GEOCLIM interpolation procedure have been averaged over the final ten years.Duration of the simulations is long enough to ensure that both models are at steadystate.The relative proportions of silicate (and among them of granitic and basaltic lithologies) and carbonate outcrops is assumed to be the same in each grid cell because of the current lack of precise lithological control.These relative proportions are not expressed in terms of relative area, but rather in terms of contribution of each lithology to the CO 2 consumption flux, so that total silicate weathering reaches 13.6×10 12 moles/yr (with 30% due to basaltic weathering), and carbonate weathering 23×10 12 moles/yr.This requires that under the same climatic conditions (in terms of mean annual air temperature and continental runoff), the contribution to the CO 2 sink per m 2 through granitic or basaltic and carbonate weathering is identical for each continental grid cell changed for all past simulations.In summary, while spatial resolution of the weathering has been introduced, we chose to keep the fractions of silicate, basalt and carbonate outcrops constant for all time slices and for each grid element.Indeed, (Gibbs et al., 1999) have demonstrated that the changes in lithology since 250 Ma weakly influence the calculation of the weathering fluxes.This result, although surprising, is mainly the consequence of the rather constant zonal relative abundance of each lithological type.

Results
Calculated atmospheric CO 2 for the 7 time slices of the geological past show the same general trend as the one simulated in Donnadieu et al. (2006) (noted D06 thereafter) (Fig. 1), a more or less monotonic decrease punctuated by a large fall at the end of the Triassic.Nevertheless, implementing vegetation as a function of climate induces changes in the absolute values of PCO 2 .Atmospheric CO 2 is always higher when the terrestrial vegetation model LPJ is used.
The ratio between atmospheric CO 2 calculated for each set of simulations (with or without interactive vegetation) remains more or less constant around 2, with a maximum value of 2.5 for the Toarcian geography and a minimum value of 1.7 for the Cenomanian geography.In  A substantial portion of continental surface becomes covered by deserts (Fig. 2).The maximum desert surface is found in the middle of the Triassic (58 millions of km 2 ) corroborating the picture of a warm and arid Triassic climate.As a consequence, FOAM-LPJ runs are systematically colder and drier than FOAM-Constant runs.Indeed, the appearance of large desert surfaces in our Mesozoic runs results in an increase of the mean continental albedo (Fig. 3) mainly located over the tropical latitudes.Other vegetation types than the mixed deciduous evergreen forest take place in the FOAM-LPJ runs (Fig. 2).However, their effect is harder to quantify because biophysical values of the various vegetation types in FOAM are similar (Table 2).In order to understand what happens with the free evolving vegetation runs, we have plotted the changes in the mean zonal continental temperature occurring between FOAM-Constant and FOAM-LPJ runs (Fig. 4).
The largest cooling in FOAM-LPJ runs are clearly located over the latitudes where the barren surfaces have replaced the previously artificial imposed mixed deciduous vegetation type.Nevertheless, the cooling occurs at all latitudes.This is due to the mixing of the atmosphere through Hadley cells and mid-latitude eddies.The large cooling occurring over the high latitudes seems to be linked to the snow albedo feedback.Indeed, the cooling induced by the appearance of desert leads to a greater snow cover (Fig. 5).In summary, using the previously modelled atmospheric CO 2 of D06, FOAM-LPJ systematically predicts lower temperature and runoff mainly as a result of the appearance of desert.Hence, total silicate weathering is reduced, and CO 2 rises until a new balance with the volcanic degassing flux is reached again.We would like to emphasize here that the primary effect of the vegetation is to warm and make the Earth wetter.Hence, if our FOAM-LPJ runs had been compared to FOAM runs in which barren surface was imposed on every grid land point, we would have found a warmer and wetter world in FOAM-LPJ runs.Unfortunately, previous GEO-CLIM runs were performed with mixed deciduous evergreen forest on every grid land point.
Figure 6 shows the time evolution of the mean annual global and continental temperature for both sets of simulations (with and without vegetation dynamics).From a general point of view, the long term trend remains the same with a decrease in mean global temperature mainly located Fig. 4. Difference in mean annual continental temperature as simulated between FOAM and FOAM-LPJ runs using the same orbital parameters for all simulations.Atmospheric CO 2 are those used in the calculation of Table 1.
at the Triassic-Jurassic boundary though the cooling from the Carnian to Rhaetian is enhanced in GEOCLIM-LPJ.In previous runs, the cooling reached 4 • C contrasting with the almost 6.5 • C cooling found here.Origin of this massive atmospheric CO 2 decrease can still be found in the northward migration of Pangea as explained in detail in Godderis et al. (2008).However, the largest cooling simulated here is due to the enhanced radiative cooling directly related to the CO 2 (the ratio of CO 2 between both time periods is 4.8 here compared to 3.5 in D06).In order to figure out the relationship between atmospheric CO 2 and global temperature in GEOCLIM-LPJ, we have plotted the atmospheric CO 2 in Fig. 6a using a logarithmic scale.We then move arbitrarily "the CO 2 curve" until a match with the simulated temperature is obtained.Indeed, mean global temperature are linked to the atmospheric CO 2 via relationship such as T =k×ln (pCO 2 ) and this should theoretically yield a perfect match between the terms "ln(pCO 2 )" and the mean global temperature.Nevertheless, other processes also influence the radiative budget of the Earth such as the albedo change through tectonic drifting for example.If we match Toarcian to lower Cretaceous CO 2 and their respective global temperature (Fig. 6), one can note that 1) the middle and upper Cretaceous are warmer than they should be in a only CO 2 forcing case and 2) the lower Triassic is colder than expected from CO 2 forcing.Several feedbacks exist in the climate system that may affect the relation between pCO 2 and the climate history.For example, Donnadieu et al. (2006) found that changes in geography from the early to mid-to-late Cretaceous, through the response of water cycle and the changes in thermal gradient, induce a 3.8 • C warming at constant atmospheric CO 2 level.Processes described in this paper match very well "the warm bias" of the mid-to-late Cretaceous seen here.In summary, the relationship between the atmospheric CO 2 and the mean temperature is conserved over the course of our Mesozoic runs though internal feedbacks of the climatic system modulate the relationship.Still on Fig. 6, we have also represented the mean oceanic and the mean continental temperature time evolution in GEOCLIM-LPJ runs.These temperatures diverge after the Toarcian time period with continental temperatures colder than the global ones and oceanic temperatures remaining relatively stable.This divergence of continental and marine curves exists with or without vegetation dynamics.To understand the underlying mechanism, let's take a look at the mean annual temperature distribution for the three time periods during which a diverging trend between global and continental temperatures takes place (Fig. 6).We can first note that this period corresponds to the acquisition of two major characters of the modern geography i.e. the existence of pole-centred continents (see also Fig. 2 to have a general overview of geographies) and the existence of a Northern Polar Ocean surrounded by land masses.The accentuated cooling over land is mainly due to the drifting of the Antarctica continent over the South Pole.
In addition, general cooling induced by decreasing CO 2 increases the surface of continents favourable to winter snow cover resulting in a positive feedback.Figure 7 clearly shows the polar thermal amplification of the cooling first generated by the CO 2 decrease.In detail, the cooling signal from Toarcian to Berriasian follows the drifting of Antarctic and Australia toward the South Pole.The continuous cooling from Berriasian to Cenomanian is more due to the snow-albedo feedback, with a larger surface of the Northern Hemisphere being snow covered.

Discussions
Implementation of LPJ in GEOCLIM allows us to account for the dynamic effects of vegetation on climate.Main conclusions regarding the previous results obtained with a unique vegetation type ("mixed deciduous evergreen forest") on every land point are: 1.For the same atmospheric CO 2 level, free evolving vegetation induces a cooler and a drier world.Indeed, large barren surfaces replace the previously imposed PFT mixed deciduous evergreen forest and lead to a cooling and to a modification of the hydrological cycle.
2. However, and surprisingly, changes in atmospheric CO 2 level following the inclusion of vegetation dynamics remains broadly constant with a doubling of CO 2 and appears only weakly dependant on the various continental positions occurring during the Mesozoic.As a result the general trend of the CO 2 history is not sensitive to the implementation of dynamic biomes, but the absolute level is.
Fig. 5. Difference in mean annual snow depth as simulated between FOAM and FOAM-LPJ runs using the same orbital parameters for all simulations.Atmospheric CO 2 are those used in the calculation of Table 1.This result demonstrates that the biome distribution impacts heavily on the global climatic evolution of the Earth.Even without accounting for the direct effect of land plants on weathering, the biome distribution can be responsible for atmospheric CO 2 changes by a factor of 2 compared to a simulation where no biome reactivity to climate changes is modelled.This emphasizes the need for explicit modelling of continental vegetation in the modelling of past CO 2 evolution, including the impact on continental surface properties.Because CO 2 globally increases by a factor of 2 when biome dynamic is considered, compared to the D06 set of simula- tions where a global temperate forest was prescribed, it might be expected that without any land plant cover, CO 2 would further increase (mainly because of an albedo increase when vegetation is replaced by barren surface, resulting in a global cooling of the continental surface and a collapse of the CO 2 consumption through continental weathering).This effect would further contribute to the global CO 2 drawdown at the end of the Devonian.Indeed, Berner and Kothavala (2001) have calculated that the rise of land plants would decrease atmospheric CO 2 from about 3300 ppmv down to 1200 ppmv (a factor of 2.75) between 380 and 350 Ma through silicate weathering enhancement.This is the maximum effect since organic carbon burial is increasing at the same time in the GEOCARB III model, a process further reducing atmospheric CO 2 .We demonstrate here that the climatic effect of land plants (albedo, roughness, water cycling changes) might be potentially of the same order of magnitude, further enhancing the cooling linked to the rise of land plants.Modelling of this key period of the Earth history (end Devonianearly Carboniferous) should be undertaken accounting for the global effect of the radiation of plants on land.While interesting, the absolute atmospheric CO 2 levels calculated here remain very theoretical and need to be discussed regarding data constraints.As discussed in Godderis et al. (2008), paleo-proxies of atmospheric CO 2 still suffer from major uncertainties.The stomatal index method fully relies on present-day correlation between stomatal indexes and the CO 2 partial pressure for living species (such as Gingko; Retallack, 2001).Such relationships are then applied to the geological past and applied to CO 2 levels very far beyond any observed values.The method based on measurement of carbon isotopic composition of pedogenic car- bonates mainly suffers from the estimate of the link existing between the soil respired CO 2 level and the atmospheric pCO 2 (Royer et al., 2001).Soil CO 2 levels are not well constrained for present-day environments, and they largely depend on the soil hydrology, soil temperature and on the vegetation type and production, parameters that cannot easily be defined for past environments.As a result of these various uncertainties, for the period spanning the first 60 millions of years following the Permo-Triassic boundary, CO 2 paleosol reconstructions show a completely different pattern than the CO 2 stomatal reconstructions, the latter reconstruction being in better agreement with our CO 2 simulation (Godderis et al., 2008).Models also suffer from major uncertainties due to the complexity of the Earth climatic system.We have seen here that the inclusion of a dynamic vegetation model results in higher atmospheric CO 2 concentrations in the coupled carbon-climate model, compared to simulations were a global vegetation type is prescribed.Nevertheless, various unaccounted processes remain to be implemented in order to come out with an integrated picture of the CO 2 evolution over the Mesozoic.For example, a work, still in progress, has been undertaken to explicitly simulate the chemical weathering within soil rather than to continue to use parametric law providing us with a broad estimate of CO 2 consumption through weathering.In addition, another difficulty arises from the short response time of carbon in the exospheric system (200 kyr, see in Francois and Godderis, 1998).Within a few million years, atmospheric CO 2 levels may rise or fall by one order of magnitude, apart from any www.clim-past.net/5/85/2009/Clim.Past, 5, 85-96, 2009 long-term trends that are explored in this contribution, and the proxies may record these abrupt changes.Indeed, the GEOCLIM-LPJ experiments designed here are only helpful to understand what forces the long-term trend of the climate.However, the paleo-climate evolution integrates short, mid, and long-term processes which can destabilize the climate for a given time period.In order to illustrate our purpose, let us now consider the Mesozoic climate evolution using the last updated compilation of carbonate marine oxygen isotope data from published literature (Prokoph et al., 2008).The database has been structured according to the paleolatitudinal band where fossils were found.Here, we use the midlatitude δ 18 O database because the mean number of samples per Ma is the most regular and the most abundant for the Mesozoic time period.Quality and default of the carbonate δ 18 O isotopic signal as a proxy for marine temperatures have been discussed intensively in the literature (Veizer et al., 1999)  (2) local evaporation/precipitation ratio and continental runoff that influence both the salinity and the local δ 18 O of seawater (Roche et al., 2006); (3) interactions between the hydrosphere and the silicate minerals on continents and on the seafloor.The Mesozoic period has often been considered to be "ice-free" because of the absence of glacial deposits during its major part.This would correspond to a δ 18 O for ocean water close to −1% SMOW, as suggested for an icefree Earth by Kennett and Shackelton (Kennett and Shackelton, 1976).However, ice-caps (though probably smaller than today) may have been episodically present during the Mesozoic (Frakes et al., 1992).Here, in order to calculate marine surface paleotemperatures, we use the palaeothermometer of Kim and O'Neil (Kim and O'Neil, 1997) and assume a water δ 18 O value of -1 %.Hence, temperatures plotted in Fig. 6 should be seen as the maximum values for the mid-latitudes as accounting for possible small ice-sheets will increase the sea water δ 18 O.The traditional view of a uniformly warm Mesozoic appears untenable regarding the isotopic oxygen record.Brief icy intervals lie throughout the 200 Ma record and corroborate previous studies having put forward evidences for cold snaps within the alleged warm mode of the Mesozoic (Frakes et al., 1992;Price, 1999).The largest cooling event occurring in the Prokoph database is located at the Callovian-Oxfordian boundary.Sea surface temperatures from the northern mid-latitudes record a 7 • C cooling.Using δ 18 O measured on apatite from fish teeth, Dromart et al. (2003) suggested a cooling of 3 to 5 • C depending on the value of the δ 18 O of sea-water.They also noted that the seawater cooling indicated by the O-isotope paleothermometry coincided exactly with the massive southward movement of boreal ammonites.The documentation of this cold snap with three independent climate proxies shows that it is a real and firm event but also that the Prokoph database produces a trustable isotopic signal.The second largest cooling event of the Mesozoic database is the one occurring at the Berriasian -Valanginian boundary.Accounting for these events and considering the Prokoph et al. (2008) database, the alleged warm Cretaceous is considerably reduced in length and appears to have taking place within the Late Cretaceous (100-70 Ma).In detail, the amplitude of the variability becomes important within this period, including a short cold snapshot in the Campanian, and makes it hard to assimilate the Late Cretaceous to a mean "warm" period.Nevertheless, it seems that the 20 Myr spanning the Cenomanian to the Santonian stages were really warm and only interrupted by brief intervals of minor cooling.We use "minor" because when looking at the 200 Ma paleotemperature record, minimum values of the sea surface temperature remain well above those found during the Callovian-Oxfordian and the Valanginian.Then, a short cooling event occurs at the beginning of the Campanian (83.5-70.6Ma) followed by a warm period.The Maastrichtian records a cooling that settles back the temperature to the mean base level of the Early Cretaceous.Fish teeth δ 18 O data gathered by Puceat et al. (2007) have already demonstrated that marine temperatures of the climatic optimum of the middle Cretaceous and of the cooler latest Cretaceous were separated by 7 • C on average.Amplitude of the cooling estimated here is of the same order.In summary, the Prokoph database show a climate highly unstable at the million years time step punctuated by warm and cool events.This first conclusion leads us to question the notion of a mean 10 Myr time step climate used in our longterm climate model.Indeed, it seems that the climate is never equilibrated and shifts from one instability to the other.Nevertheless, a way to validate our model is to see where it sets the mean climate state (see Fig. 8).For that, we calculate the mid-latitude oceanic surface temperatures in the Tethys realm for each paleogeography and steady-state CO 2 simulated.Green crosses are the temperatures as simulated by GEOCLIM-LPJ using a constant Earth degassing flux.The calculated temperatures define a multi-million year baseline, standing for the climate that should be reached if the Earth system was close to steady-state and not affected by the various processes occurring at shorter timescales than 10 Myr.Despite the simplicity of our modelling system, our midlatitudes temperature estimates fall on the mean trend of the dataset and are in good agreement in terms of absolute values.These results support the conclusion that the continental configuration is the first order process controlling the longterm climatic evolution.It is interesting to note that our mean climate trend allows to qualify the Callovian-Oxfordian and the Berriasian-Valanginian events as cool events according to our modelled baseline.In the same way, the early Jurassic and the late Cretaceous temperature estimates are well above our simulated temperatures and as such can be qualified as exceptionally warm periods (exceptionally here meaning a warming exceeding the baseline value).As already envisioned before, processes we have focused on here can not account for millions years time step instabilities, continental configurations through its impact on silicate weathering only fix the long-term mean atmospheric CO 2 level trend.Then, in the real world, several processes can modify the long-term trend.Among them, let us cite the perturbations due to the LIPs, to the OAE or to changes in global ocean dynamics (Courtillot and Renne, 2003;Jenkyns, 1988;Knoll et al., 1996).The next steps toward a better understanding of the forcing factors controlling the Mesozoic climate evolution will thus imply the consideration of these processes.This will consist in developing specific modules accounting for the sulphur cycle if we want to be able to accurately simulate the impact of LIPs or better accounting for the watersediment interface processes for OAE.In summary, our study provides the first explicit quantification of the effect of plants on weathering through climatic changes.Data/model comparisons allow to validate our simulated climate long-term trend as a plausible baseline and demonstrate the importance of the tectonic forcing factor.More explicit modelling will help to reveal the importance of the tectonic forcing over the other factors cited above.
Edited by: R. Zeebe The publication of this article is financed by CNRS-INSU.

Fig. 3 .
Fig. 3. Difference in percentage in the surface albedo of the continental surface between FOAM-LPJ and FOAM for the Carnian.The albedo of continental surface in runs from D06 was 0.06.The albedo of continental surface in this study is now provided by LPJ.

Y.
Donnadieu et al.:  Biophysical effect of vegetation on Mesozoic atmospheric CO 2 snow cover depth (cm) Changes in annual snow cover depth (cm) Changes in annual snow cover depth (cm)

Fig. 6 .
Fig. 6.Mean annual global temperature evolution over the Mesozoic as simulated by the GEOCLIM-LPJ version and by the GEOCLIM-only version from Donnadieu et al. (2006a).Continental, marine and global trends are shown.Red crosses are the atmospheric CO 2 estimates simulated with GEOCLIM-LPJ (see Fig. 1) and have been arbitrarily positioned in Fig. 6a in order to fit the early Jurassic and the early Cretaceous temperature estimates from the GEOCLIM-LPJ runs.

Fig. 7 .
Fig. 7. Geographical distribution of the mean annual temperature for three time-periods (a) Toarcian (b) Berriasian and (c) Cenomanian.Atmospheric CO 2 levels are those estimated by GEOCLIM-LPJ.

Fig. 8 .
Fig.8.Northern mid-latitudes sea surface temperature estimates from the δ 18 O database published byProkoph et al. (2008) using the palaeothermometer ofKim and O'Neil (1997) and a sea water δ 18 O value of −1% (see text for details).The brown dots are the original data, the black curve represents a 2 Myr time-window running average and the red curve represents a 5 Myr time-window running average.Green crosses show the temperature of the northern mid-latitudes as simulated by LPJ-GEOCLIM.

Table 1 .
Mean Annual Runoff and Continental Temperature as calculated by FOAM when using either the same vegetation type (here mixed deciduous evergreen forest) on all continental grid points as done in D06 or the dynamic global vegetation model LPJ as done in this study.Note that solar constant and orbital parameters are those used in D06 ; atmospheric CO 2 levels used are shown in the first column.

Table 1
, we have summarized the mean continental temperature and the mean runoff as estimated by FOAM-LPJ and FOAM-Constant (i.e.mixed deciduous evergreen forest imposed on every continental points) using atmospheric CO 2 levels closest as possible to those calculated in D06.Continental temperature and

Table 2 .
Parameters changed in the FOAM surface scheme as a function of LPJ results.
and have to be kept in mind when interpreting the isotopic signal in terms of seawater temperatures.The 18 O/ 16 O ratio recorded in carbonaceous fossils depends on both temperature and the isotopic composition of the water.Knowledge of the δ 18 O of seawater is therefore crucial for the temperature calculation.Variations in the oxygen isotope composition of seawater are governed by several factors including (1) the evolution of the mass of continental ice that modifies the δ 18 O of seawater by storing preferentially 16 O in ice;