# Rates of consumption of atmospheric CO2 through the weathering of loess during the next 100 yr of climate change

Y. Goddéris, S. L. Brantley, L. M. François, J. Schott, D. Pollard, M. Déqué, M. Dury

► To cite this version:

Y. Goddéris, S. L. Brantley, L. M. François, J. Schott, D. Pollard, et al.. Rates of consumption of atmospheric CO2 through the weathering of loess during the next 100 yr of climate change. Biogeo-sciences, 2013, 10, pp.135-148. 10.5194/bg-10-135-2013. insu-03620501

HAL Id: insu-03620501

https://insu.hal.science/insu-03620501

Submitted on 27 Mar 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>

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

# Rates of consumption of atmospheric CO₂ through the weathering of loess during the next 100 yr of climate change

Y. Goddéris<fnref n="1" />, S. L. Brantley<fnref n="2" />, L. M. François<fnref n="3" />, J. Schott<fnref n="1" />, D. Pollard<fnref n="2" />, M. Déqué⁴, and M. Dury<fnref n="3" />

*Correspondence to:* Y. Goddéris (yves.godderis@getobs-mip.fr)

Received: 9 July 2012 – Published in Biogeosciences Discuss.: 14 August 2012

Revised: 29 November 2012 – Accepted: 5 December 2012 – Published: 9 January 2013

**Abstract.** Quantifying how C fluxes will change in the future is a complex task for models because of the coupling between climate, hydrology, and biogeochemical reactions. Here we investigate how pedogenesis of the Peoria loess, which has been weathering for the last 13 kyr, will respond over the next 100 yr of climate change. Using a cascade of numerical models for climate (ARPEGE), vegetation (CARAIB) and weathering (WITCH), we explore the effect of an increase in CO₂ of 315 ppmv (1950) to 700 ppmv (2100 projection). The increasing CO₂ results in an increase in temperature along the entire transect. In contrast, drainage increases slightly for a focus pedon in the south but decreases strongly in the north. These two variables largely determine the behavior of weathering. In addition, although CO₂ production rate increases in the soils in response to global warming, the rate of diffusion back to the atmosphere also increases, maintaining a roughly constant or even decreasing CO₂ concentration in the soil gas phase. Our simulations predict that temperature increasing in the next 100 yr causes the weathering rates of the silicates to increase into the future. In contrast, the weathering rate of dolomite – which consumes most of the CO₂ – decreases in both end members (south and north) of the transect due to its retrograde solubility. We thus infer slower rates of advance of the dolomite reaction front into the subsurface, and faster rates of advance of the silicate reaction front. However, additional simulations for 9 pedons located along the north–south transect show that the dolomite weathering advance rate will increase in the central part of the Mississippi Valley, owing to a maximum in the response of vertical drainage to the ongoing climate change.

The carbonate reaction front can be likened to a terrestrial lysocline because it represents a depth interval over which carbonate dissolution rates increase drastically. However, in contrast to the lower pH and shallower lysocline expected in the oceans with increasing atmospheric CO₂, we predict a deeper lysocline in future soils. Furthermore, in the central Mississippi Valley, soil lysocline deepening accelerates but in the south and north the deepening rate slows. This result illustrates the complex behavior of carbonate weathering facing short term global climate change. Predicting the global response of terrestrial weathering to increased atmospheric CO₂ and temperature in the future will mostly depend upon our ability to make precise assessments of which areas of the globe increase or decrease in precipitation and soil drainage.

## 1 Introduction

Today, continental weathering consumes about 0.3 Gt of atmospheric CO₂ each year (Gaillardet et al., 1999). This value is close to the preindustrial net carbon exchange between the atmosphere and the terrestrial biosphere (0.4 Gt C yr⁻<fnref n="1" />) (Solomon et al., 2007). Continental weathering is sensitive to land use change (Raymond et al., 2008) and ongoing climate change (Gislason et al., 2009; Beaulieu et al., 2012). Despite its magnitude, apparent sensitivity, and acknowledged importance over geologic timescales, continental weathering is generally neglected in studies dealing with the evolution of the global carbon cycle and climate over century timescales.

Continental weathering is a function of the mineralogy of the exposed surfaces (Amiotte-Suchet et al., 2003), among

<footnote n="1">¹Géosciences Environnement Toulouse, CNRS-Observatoire Midi-Pyrénées, Toulouse, France</footnote>

<footnote n="2">²Earth and Environmental Systems Institute, Pennsylvania State University, University Park, PA, USA</footnote>

<footnote n="3">³Institut d'Astrophysique et de Géophysique, Université de Liège, Belgium</footnote>
many other parameters such as physical erosion (West et al., 2005), climate (White and Blum, 1995), and vegetation cover (Le Hir et al., 2011; Taylor et al., 2012). This multi-parameter dependence makes it difficult to assess the response of continental weathering to climate change on a global basis. Parametric laws can be used to predict weathering from mean annual temperature and runoff for several average lithological types (Hartmann et al., 2009). These laws were established from field observations. However, as pointed out in Goddéris et al. (2009), such a parametric method only gives a snapshot of the weathering system, and does not capture the decadal to centennial dynamics. In this contribution, we use mechanistic models to predict the response of weathering to increasing CO₂ in the future by focusing on the natural experiment of loess that was deposited approximately 13 kyr ago along the Mississippi valley transect in the USA. We have previously shown that mechanistic weathering models could successfully reproduce mineralogical patterns observed in the loess along the climosequence (Williams et al., 2010; Goddéris et al., 2010). This latter modeling effort lends credibility to our attempt here to predict weathering forward in time.

In addition, loess and loess-derived sediments cover large areas of the continents, especially in the Northern Hemisphere. Ten percent of the continents are presently covered by loessic formations derived from glacial eolian deposits (Rousseau, 2001; Haase et al., 2007). Given their abundance at the Earth's surface, they should contribute significantly to the global weathering flux in the ocean for at least two reasons: (1) loess minerals often have only slightly weathered surfaces, and (2) the specific surface area is high owing to the fine-grained nature of loess. Furthermore, since loess contains highly reactive carbonate phases, loess weathering could potentially respond rapidly to environmental change. Indeed, researchers have suggested that the high surface area of loess may cause fast weathering that can be an important cause for the long-term drawdown of CO₂ (Anderson et al., 1997). Understanding how these loessic formations are weathered might thus be key toward quantifying the future evolution of the global carbon cycle.

In North America, loess is one of the most extensive superficial deposits, especially in the central Great Plains region where it is known as the Peoria loess (Rousseau, 2001; Muhs et al., 2008). Minerals include primary silicates (potassium feldspar, plagioclase feldspar, quartz), secondary clay minerals (dominated by smectites), and carbonate phases (calcite and dolomite) (Muhs et al., 2008). Given this mixed mineralogy, a quantitative description of the evolution of loess weathering as a function of climate change is not an easy task.

Using the WITCH numerical model forced by the GENE-SIS climate model, we previously calculated the weathering rates of the Peoria loess over the last 10 kyr for two sites, one located at the northern end of the Mississippi Valley, and the other located at the southern end (Goddéris et al., 2010). The model starts with the initial mineralogical composition

of the aeolian deposits 10 kyr ago. It reproduced the time evolution of the mineralogical composition of the two pedons, simulating mineralogical concentrations close to measured values (Goddéris et al., 2010). Features that were simulated included the deepening of the reaction fronts of dolomite, albite, and potassium feldspar. These results suggest that the calculated weathering rates of the pedons were correctly simulated as a function of climate by the WITCH model. They also show that the present mineralogical composition of the loessic formation is the result of the weathering of the loess under variable climatic conditions. In other words, climate has left its imprint on the loessic profiles. The question we ask in this contribution is: what will be the effect of the climate during the next century on the weathering processes at play in the Peoria loess? We thus start the numerical modeling from the present-day state of the weathering system and project it forward to 2100 for the same two pedons to explore how weathering of Mississippi Valley loess will respond to future climate changes.

## 2 Methods

### 2.1 Model settings

To project weathering into the future, we ran an atmospheric global circulation model (ARPEGE, Salas y Mélia et al., 2005; Gibelin and Déqué, 2003) under the A1B IPCC carbon emission scenario (Solomon et al., 2007). The output of this simulation (mean air temperature and diurnal temperature range, rainfall, relative air humidity, cloud cover and surface horizontal wind speed) was then used to force a global dynamic vegetation model (CARAIB, Dury et al., 2011). CARAIB predicts the vegetation cover, the biospheric productivity, the below-ground hydrology and CO₂ production along the N-S transect of the Mississippi Valley. These physical variables were then used to force WITCH, a process-based model of weathering reactions (Goddéris et al., 2006; Goddéris et al., 2010; see Appendix A for more details).

The three numerical models have been previously validated. The climate model ARPEGE has been part of the international stretched-grid model intercomparison project (SG-MIP). ARPEGE was able to reproduce the North American climate evolution over 12 yr (1987–1998 period) (Fox-Rabinovitz et al., 2006). The CARAIB biospheric model is a dynamic vegetation model initially developed for the global scale (Warnant et al., 1994). As such it has been used and validated over all continents (Nemry et al., 1996, 1999). Model vegetation distribution maps have also been validated globally (Otto et al., 2002). Here, a globally validated classification of plant functional types (PFT) is used (updated from François et al., 2011). Finally, the WITCH weathering model has been used to simulate water-mineral interactions in various environments, from polar to tropical settings (Goddéris et al., 2006; Roelandt et al., 2010; Violette et al., 2010; Beaulieu
et al., 2011). WITCH was able to reproduce the concentration in base cations and silica in soil solutions and at the outlet of the simulated catchments.

The available ARPEGE climate simulation covers the 1950–2100 period. Instead of running the model cascade only into the future (from 2012 to 2100), we decided to perform the simulation over the whole 1950–2100 period for two reasons: (1) since the weathering response is a priori assumed to be low (specifically for silicate minerals), a longer run will maximize the weathering changes, and (2) the simulation must be started from the less anthropogenically perturbed state, for the purpose of relaxation of the initial condition.

## 2.2 Mineralogical composition of the pedons

In our previous study (Goddéris et al., 2010), we relied on published chemical analyses and X-ray diffraction data for the Peoria loess sampled from pedons along the Mississippi transect (Muhs et al., 2001). Williams et al. (2010) also reported new quantitative X-ray diffraction (XRD) analyses of Muhs' samples from pedons 1 (Clayton County, IL), 13 (Lauderdale County, TN), and 22 (St. Francisville, LA). Those samples were interpreted to contain montmorillonite, illite, plagioclase, potassium feldspar, quartz, and minor kaolinite. In addition, dolomite was observed with XRD in the deepest sample of the northern-most pedon but not in pedons 13 or 22 (i.e. in those latter pedons, dolomite < 2 or 3 wt %).

Williams (2008) developed a mineral model to calculate mineral abundances from bulk elemental analyses for the samples from pedons 1, 13, and 22 that matched the mineralogical data derived from quantitative X-ray diffraction. The mineral model included montmorillonite, illite, plagioclase, potassium feldspar, quartz, dolomite, and kaolinite plus minor calcite, chlorite, ferrihydrite. Using this mineral model, we calculated the mineralogy of the bulk elemental data from Muhs et al. (2001) for all the pedons along the north–south transect. Based on these calculations, the extent of depletion of plagioclase feldspar along the transect was observed to increase from north to south (Williams et al., 2010). In addition, evidence for a dolomite reaction front at depths between 200 and 400 cm was strong for pedons 1 through 6 (Mammoth, IL) and 8 (Greenbay Hollow, IL). The model was also consistent with a small amount of dolomite at depth in pedons 9, 17, 20–22. Consistent with this, other researchers have observed dolomite in deep loess samples near pedons 21 and 22 (Ruhe, 1984; Pye and Johnson, 1988; Muhs et al., 2001). To ascertain definitively if dolomite was present in any of these more southern profiles, we attempted to get samples in collaboration with D. Muhs (US Geological Survey). Sample volume limitations precluded analysis of most of the deepest samples from the southern part of the transect. In agreement with the mineral model, however, the one sample

(from 350 cm depth in pedon 20) that was large enough for micro-XRD revealed the presence of dolomite.

Given these observations, here we chose to follow Goddéris et al. (2010) by modeling pedons 1 in the north and 20 in the south. In both cases, we used mineralogy based on the mineral model of Williams (2008) (Table 1). This model predicts a small concentration of calcite; however, Goddéris et al. (2010) documented that calcite removal in the WITCH model was rapid. Given that such small amounts of calcite might simply reflect error in the model attribution of mineralogy, we removed all calcite from the starting pedon mineralogy. For both pedons, the initial porosity was set equal to 43 % (Bettis et al., 2003). Furthermore, because the samples reported by Muhs et al. (2001) only allowed us to estimate that the dolomite reaction front was located between about 200 cm and 350 cm depth in pedon 20, whereas it was located at approximately 350 cm in pedon 1, we decided to fix the dolomite weathering front at 280 cm depth for both pedons. By locating the two fronts at the same depth, we can compare weathering in the north and south directly. Furthermore, moving this front by a few tens of cm up or down is not critical for this study. The front was also assumed to be a sharp front, meaning that we assume the total absence of dolomite above this depth, and an abundance of 11 % vol below this level (i.e. dolomite comprises a completely developed reaction front; Brantley and White, 2009).

# 3 Results

## 3.1 Complex behavior of the driving parameters of weathering

It is generally assumed that weathering increases with runoff (runoff being the difference between rainfall and evapotranspiration) and with air temperature, at least for silicate rocks (White and Blum, 1995; Oliva et al., 2003; Dessert et al., 2003). Our simulation shows that the link between climate, vegetation and weathering is complex, especially when considered over a short timescale (here centuries).

For example, the overall temperature rise from 4 to 6 °C calculated for the 1950–2100 period along the north–south transect, as predicted by ARPEGE, should promote silicate dissolution and slow carbonate dissolution because carbonates are generally at equilibrium with soil solutions and, unlike silicates, their solubility decreases with temperature (Schott et al., 2009).

However, increases in temperature are accompanied by a change in rainfall which is likely to also play a role. For the simulations, rainfall was predicted to rise moderately from 120 to 140 mm month⁻¹ in the southern part of the Mississippi valley (31°N) over the same time period, and from 75 to 85 mm month⁻¹ for the northern area (42°N). The vertical drainage through the loess formation will also vary as rainfall increases. However, drainage does not vary in a simple way
<table><caption>**Table 1.** Pedon mineralogy (%) assumed for pedon north (42°59′ N) and pedon south (30°47′ N).</caption><thead><tr><th colspan="8">Pedon north</th></tr><tr><th>Depth (cm)</th><th>Quartz</th><th>Albite</th><th>K-spar</th><th>Montmorillonite</th><th>Dolomite</th><th>Kaolinite</th></tr></thead><tbody><tr><td>5.0</td><td>54.5</td><td>11.4</td><td>14.6</td><td>6.6</td><td>0.5</td><td>8.6</td></tr><tr><td>19.0</td><td>55.6</td><td>11.3</td><td>14.5</td><td>6.2</td><td>0.5</td><td>8.8</td></tr><tr><td>39.0</td><td>53.3</td><td>10.6</td><td>14.1</td><td>6.8</td><td>0.8</td><td>11.2</td></tr><tr><td>57.0</td><td>47.9</td><td>9.4</td><td>13.0</td><td>7.6</td><td>0.7</td><td>16.5</td></tr><tr><td>84.0</td><td>43.3</td><td>9.5</td><td>13.0</td><td>7.9</td><td>1.1</td><td>19.0</td></tr><tr><td>127.0</td><td>43.4</td><td>10.8</td><td>13.2</td><td>8.1</td><td>1.4</td><td>17.8</td></tr><tr><td>181.0</td><td>43.9</td><td>10.4</td><td>13.0</td><td>7.9</td><td>1.5</td><td>17.8</td></tr><tr><td>237.0</td><td>46.4</td><td>12.1</td><td>12.6</td><td>7.9</td><td>1.4</td><td>14.6</td></tr><tr><td>280.0</td><td>45.4</td><td>12.8</td><td>14.5</td><td>8.1</td><td>1.3</td><td>12.9</td></tr><tr><td>330.0</td><td>46.4</td><td>12.5</td><td>13.5</td><td>8.9</td><td>1.6</td><td>12.5</td></tr><tr><td>377.0</td><td>41.2</td><td>11.8</td><td>12.6</td><td>7.8</td><td>11.5</td><td>10.7</td></tr></tbody><tbody><tr><th colspan="8">Pedon south</th></tr><tr><td>1.0</td><td>63.1</td><td>6.8</td><td>10.0</td><td>4.7</td><td>0.1</td><td>11.0</td></tr><tr><td>4.0</td><td>63.6</td><td>7.0</td><td>10.7</td><td>4.6</td><td>0.1</td><td>10.8</td></tr><tr><td>19.0</td><td>62.9</td><td>6.6</td><td>11.5</td><td>4.5</td><td>0.1</td><td>11.2</td></tr><tr><td>42.0</td><td>50.9</td><td>5.5</td><td>10.9</td><td>6.9</td><td>0.0</td><td>19.7</td></tr><tr><td>66.0</td><td>48.6</td><td>6.2</td><td>11.2</td><td>7.6</td><td>0.0</td><td>20.1</td></tr><tr><td>86.0</td><td>51.0</td><td>7.5</td><td>10.4</td><td>7.7</td><td>0.0</td><td>17.8</td></tr><tr><td>103.0</td><td>50.1</td><td>8.4</td><td>11.9</td><td>7.4</td><td>0.0</td><td>16.8</td></tr><tr><td>119.0</td><td>50.1</td><td>8.7</td><td>11.6</td><td>7.5</td><td>0.0</td><td>16.6</td></tr><tr><td>144.0</td><td>49.8</td><td>9.2</td><td>11.7</td><td>8.1</td><td>0.0</td><td>16.2</td></tr><tr><td>178.0</td><td>51.1</td><td>9.8</td><td>11.1</td><td>8.1</td><td>0.3</td><td>15.0</td></tr><tr><td>350.0</td><td>43.7</td><td>9.6</td><td>10.4</td><td>7.7</td><td>11.3</td><td>10.2</td></tr></tbody></table>

since it also depends on changes in the actual evapotranspiration, which responds in a complex fashion owing to the rise in temperature and atmospheric CO₂ (including fertilization effects, stomatal closure, etc.). In the south, the yearly averaged vertical drainage calculated at the base of the pedon increases by 14% (from about 350 mm yr⁻¹) to a value slightly above 400 mm yr⁻¹ – over the simulation period (1950–2100). In contrast, the deep drainage decreases in the north by more than 40% (from 170 to < 100 mm yr⁻¹) despite an increase in rainfall. Drainage decreases in the north because evapotranspiration (AET) increases significantly as temperature increases by 6 °C (Fig. 1).

In addition to changes in temperature and precipitation, the soil CO₂ concentrations also change in a complex fashion as global temperature and atmospheric CO₂ increase. As calculated by the CARAIB model, the temperature rise promotes soil respiration, but the higher temperature also accelerates CO₂ diffusive loss to the atmosphere (Van Bavel, 1951). Therefore, even though the model-calculated soil CO₂ production rate rises from a mean annual value of 95 g C m⁻² yr⁻¹ to > 130 g C m⁻² yr⁻¹ in the south, the soil CO₂ concentration increases only from ~20 × to 22 × atmospheric levels. This lack of change in the south is the outcome of faster CO₂ diffusion to the atmosphere that competes with the CO₂ fertilization effect, which promotes

net primary productivity and hence both autotrophic and heterotrophic respiration. In contrast, this compensation does not occur in the northern site where the temperature effect on diffusion outweighs the fertilization effect causing soil CO₂ concentrations to decline from ~30 × to 20 × atmospheric levels over the 1950–2100 period. This CO₂ decline happens despite a calculated increase in the soil CO₂ production rate from 70 to 110 g C m⁻² yr⁻¹ (Fig. 2).

## 3.2 Weathering results

### 3.2.1 CO₂ consumption by dolomite weathering

The instantaneous atmospheric CO₂ consumption by weathering is calculated as the flux of HCO₃⁻ required to precisely balance the cation release by the dissolution of silicate minerals, integrated over the whole weathering profile. In addition, the occasional precipitation of smectite contributes negatively to this budget. Finally, the bicarbonate flux equal to the Ca²⁺ + Mg²⁺ release by dolomite dissolution must also be added to account for atmospheric CO₂ consumed by dolomite weathering.

Despite the increase in temperature, the total CO₂ consumption due to weathering decreases in both the north and the south from 1950 to 2100. In the south, the model predicts a slight decrease from 1 to 0.9 mol m⁻² yr⁻¹ (trend calcu-
<figure><img src="image_4.png" /></figure>

**Fig. 1.** The calculated water budget for the southern ((a) and (b)) and northern ((c) and (d)) sites. The increase in evapotranspiration driven by the temperature rise offsets the drainage response to the global warming. The 1-yr running means of the monthly outputs are plotted. The dotted lines are the mean square linear fits of the vertical drainage evolution at the basis of the pedons.

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

**Fig. 2.** Calculated yearly averaged CO$_2$ concentrations below the root zone (plotted in units of PAL, i.e., Pre-industrial Atmospheric Level) for the southern and northern sites (solid line). The dotted lines are the 10-yr running average of the annual output.

lated from a root mean square linear fit of the yearly averaged modeled CO$_2$ consumption). In the north, CO$_2$ consumption goes down by a factor of 2.3, from 0.7 to 0.3 mol m$^{-2}$ yr$^{-1}$ (Fig. 3).

Dolomite dissolution is the main weathering reaction that largely controls these net CO$_2$ consumption trends: the CO$_2$ consumption flux by dolomite dissolution is about 10 × greater than consumption by any silicate phase in the pedons. In all the pedons, dolomite is absent at shallow depths and present only beneath 2.80 m depth. This 2.80 m depth therefore represents the sharp reaction front where downward-flowing fluids that are undersaturated with respect

to dolomite begin to dissolve this mineral. The sharp reaction front and the close-to-equilibrium porefluids with dolomite are consistent with transport-controlled dissolution. Changes in the integrated dolomite weathering rate (in moles of mineral dissolved per time unit, summed over the whole weathering profile) are strongly correlated to the vertical drainage at both locations. In fact, as drainage increases over any time interval, dolomite dissolution responds linearly. This behavior is explained by the fact that waters draining the pedons are equilibrated with respect to dolomite below the dolomite front. Thus, the dissolution flux of dolomite is essentially equal to the product of the solubility of dolomite and the
<figure><img src="image_6.png" /><figcaption>**Fig. 3.** Calculated yearly averaged CO₂ consumption by weathering reactions for the southern and northern sites. In black: reference runs. In green: simulations where the air temperature is held at the averaged 1950–1960 constant value. The straight lines are the mean square linear fits of the model outputs.</figcaption></figure>

drainage flux (Brantley and White, 2009). If drainage increases, the flux increases, i.e., more water implies more dissolution. On the other hand, yearly averaged dolomite dissolution is negatively correlated to the temperature as dolomite solubility decreases with temperature. Since temperature displays a clear increasing trend from 1950 to 2100 period (4 °C warming) while drainage only increases slightly in the south or decreases drastically in the north (Fig. 1), the overall effect is a decrease in the CO₂ consumed by dolomite dissolution at both locations.

### 3.2.2 CO₂ consumption by primary silicate phases

In contrast to dolomite, the CO₂ consumption by dissolution of the primary mineral phases albite and K-feldspar either increase or stay at roughly constant values from 1950 to 2100. In the south, equilibrium with respect to albite in the soil solution is episodically reached below 1 m depth, and permanently reached below the dolomite front at ~2.80 m depth. In the plots of soil solution saturation state with respect to albite versus time, pulses of high drainage are observed to be interspersed with periods of low drainage. High-drainage spikes lead to dilution below equilibrium in pore fluids down to the dolomite front, but the soil solution remains saturated with respect to albite below 2.80 m depth (Fig. 4a). The rate of occurrence of the dry events (low drainage periods when albite may effectively stop dissolving high in the profiles) increases after 2000 (Fig. 4a). But at the same time, the annually averaged drainage slightly increases in the south over the 1950–2100 period (14% increase), promoting albite dissolution. The temperature rise (4 °C) forces albite dissolution to further increase. This results in an increasing trend in the calculated Na⁺ export for the pedon from about 380 mol ha⁻¹ yr⁻¹ (1950) to 470 mol ha⁻¹ yr⁻¹ (2100).

For the northern pedon, the soil solution gets saturated with respect to albite mainly at the dolomite front even before 2000, i.e., at 2.80 m depth (Fig. 4b). This saturation level

<figure><img src="image_7.png" /><figcaption>**Fig. 4.** Calculated time evolution of the saturation state of the soil solutions with respect to albite, plotted as a function of depth: (a) southern pedon, (b) northern pedon.</figcaption></figure>
is reached at shallower depth after 2000. Several drought
events strongly inhibit the vertical drainage, and these events
become more and more common after 2000, owing to the
overall decrease in drainage. As a result, the models predict
an almost permanent rise of the saturation depth by about 2 m
after 2030. However, this line cannot rise above the rooting
zone because the through-flux of water is higher in that zone,
as water that enters at the land surface is taken up by roots in-
creasingly with depth. In contrast, beneath the rooting zone,
the water flux is smaller and the water quickly equilibrates
with respect to the low-solubility feldspar minerals (Fig. 4b).

Like the dolomite behavior described above, the albite
dissolution flux varies directly with vertical drainage in the
north and in the south. In other words, just as dolomite dis-
solution is driven by the drainage flux, because dolomite is
maintained at chemical equilibrium at the base of the pedons,
the deep-soil solutions are at chemical equilibrium or even
oversaturated with respect to albite. Any increase in drainage
thus enhances the observed dissolution flux. In the north, the
calculated decrease in drainage over the next century com-
pensates for the temperature rise and the yearly averaged
albite dissolution fluxes stay constant (or even slightly de-
crease) from 1950 to 2100. In contrast, in the south, where
the drainage flux increases over the simulation period, the
albite dissolution flux rises as the global temperature rises.

Similar arguments can be made for the behavior of K-
feldspar dissolution. For both pedons, soil solutions are satu-
rated with respect to K-feldspar below the root zone at about
1 m. In the northern profile, dissolution increases from 1950
to 2100 and CO2 consumption rises by 15%, whereas in the
south it rises by 25%.

Not surprisingly, the flux from K-feldspar dissolution is
positively correlated with vertical drainage at both sites. In-
deed, like many two-feldspar systems (White et al., 2001),
the saturation states of the draining waters with respect to
K-feldspar are much higher than the saturation state with re-
spect to albite. As a consequence, the dissolution flux from
K-feldspar depends directly on the availability of water be-
cause everywhere in the profile dissolution occurs close to
equilibrium with respect to the mineral phase: K-feldspar
shows transport-controlled behavior.

3.2.3 CO₂ consumption by secondary silicate phases

Draining waters of both northern and southern sites are gen-
erally undersaturated with respect to smectite. Equilibrium
with this phase is only attained during dry events in the
northern site. For this site, the saturation state of the drain-
ing waters generally increases from 1950 to 2100 because
of the occurrence of dry events related to the decreasing
rainfall (Fig. 5). These dry events have a strong impact on
the overall CO2 consumption. As smectite precipitates, base
cations are removed from the solution, limiting the atmo-
spheric CO2 consumption. Furthermore, these dry events re-
duce dolomite dissolution because drainage decreases down

<figure><img src="image_8.png" /><figcaption>**Fig. 5.** Calculated smectite dissolution rates plotted as a function of depth and time for the northern pedon.</figcaption></figure>

to 0 mm yr⁻¹ below the root zone (Fig. 1b). As a consequence, CO₂ consumption falls to 0 mol m⁻² yr⁻¹ during extremely dry events such as the one indicated around 2025 (Fig. 3b).

The aridity developing in the north of the Mississippi val-
ley promotes a general increase in the pH of the soil solu-
tions (Fig. 6) (from 6.5 to 6.8 over 150 yr). There is a strong
positive correlation between temperature and mean pH of
the weathering profile. The correlation between mean annual
drainage and mean annual pH is conversely negative, i.e.,
low-draining periods are characterized by higher pH. Both
the temperature rise and the drainage decrease from 1950
to 2100 thus promote a pH rise in the soil solutions in the
north. On the contrary, pH oscillates around an average value
of 6.65 over the whole simulation period for the south pe-
don. It is also interesting to note that porewaters available
to plants change in composition and pH forward in time. In
other words, the chemistry at rooting depth is affected by
climate change and plants have access to fluids of different
chemistry over time. Of course, there is no feedback of the
soil chemistry on the biospheric model in the present study,
and such feedback should be considered in the future.

3.2.4 Relative contribution of the dissolution of each mineral to CO₂ consumption

Overall, the relative contribution of each mineral phase
to total CO2 consumption depends largely on the vertical
drainage at both sites. The dolomite contribution to total
CO2 consumption is generally above 90% but can decrease
to less than 40% during dry events in the north (Fig. 7).
Conversely, the relative silicate contribution increases during
low drainage events. Indeed, the associated CO2 consump-
tion flux goes down when drainage decreases, but at a slower
rate than the flux due to dolomite alone.
<figure><img src="image_9.png" /><figcaption>Fig. 6. Calculated soil solution pH plotted as a function of depth and time for the northern pedon.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Fig. 7. Relative contribution of calculated dolomite dissolution to the CO₂ consumption in the north plotted as a function of vertical drainage. Each point stands for a yearly averaged value.</figcaption></figure>

### 3.2.5 Calculated CO₂ consumption compared to previous estimates

Our simulations predict a CO₂ consumption of about 1 mol CO₂ m⁻² yr⁻¹ for the southern pedon and 0.3–0.7 mol CO₂ m⁻² yr⁻¹ for the northern pedon. Since most of the CO₂ consumption is the result of dolomite dissolution, these numbers can be compared to the CO₂ consumption by outcropping carbonate globally. Carbonate outcrops (including marls, dolomites, and limestones predominantly) cover about 13.4% of the continental surface (20.1 × 10⁶ km²) (Amiotte-Suchet et al., 2003). They consume 12.3 × 10¹² mol of atmospheric CO₂ yr⁻¹ (Gaillardet et al., 1999). This gives a mean area-averaged CO₂ consumption of 0.62 mol m⁻² yr⁻¹, close to our estimated values for the loess.

Results can be also compared with the modeled values (using parametric laws calibrated on field data) from Moosdorf et al. (2011). For unconsolidated sediments of North America, including loess, CO₂ consumption is estimated at 0.2 mol CO₂ m⁻² yr⁻¹ for the southern pedon using the Moosdorf et al. (2011) correlation between CO₂ consumption and runoff (taken as the difference between rainfall and evapotranspiration), and about 0.05 mol CO₂ m⁻² yr⁻¹ for the northern pedon. For carbonate rocks of North America, these numbers rise to 0.8 and 0.25 mol CO₂ m⁻² yr⁻¹, respectively, in Moosdorf et al. (2011). These numbers are of the same order of magnitude as our results.

## 4 Sensitivity tests

The key role of temperature on dolomite dissolution and the associated trends in CO₂ consumption can be illustrated by a simulation where monthly air temperature is fixed at its monthly averaged 1950–1960 value. In the south, the overall CO₂ consumption now increases at a secular rate of 1 mmol CO₂ m⁻² yr⁻², while it was decreasing at a rate of −1 mmol CO₂ m⁻² yr⁻² in the reference run. In the north, the decreasing trend of −2.5 mmol CO₂ m⁻² yr⁻² (reference run) is converted into an increasing rate of 0.25 mmol CO₂ m⁻² yr⁻² (sensitivity test; Fig. 3).

This does not mean that temperature is a stronger driver of the CO₂ consumption by weathering than drainage for these soil profiles. Indeed, for both sites, the base cation export fluxes at the base of each pedon are strongly correlated to the vertical drainage. But the vertical drainage only slightly increases in the south while it strongly decreases in the north, owing to enhanced evapotranspiration in a warmer world. The response in terms of weathering is strongly dependent on the climate and biospheric model behavior: temperature appears as a key factor because the ARPEGE GCM predicts a rather strong global warming in the Mississippi valley (4 to 6 °C), which impacts the biospheric model in such a way that the evapotranspiration increase compensates for the rainfall increase, limiting the drainage response.

## 5 Carbonate fronts along the Mississippi transect

The detailed analysis presented above was focused on the most southern and northern pedons because those are the two sites where a precise mineralogical composition has been calculated as a function of depth. Nevertheless, these two mineralogical profiles are not widely different (see Table 1). We thus expand our calculations to the whole Mississippi latitudinal transect, assuming the mineralogical composition of the southern pedon for all sites. As discussed by Williams et al. (2010), the Peoria loess is fairly invariant in composition. Our goal is to explore the response to climate change of the most reactive mineral phase with respect to the CO₂ budget
(here dolomite, as demonstrated above), along a latitudinal
transect.

The dolomite front corresponds to the depth below which soil solutions are supersaturated with respect to dolomite, and above which they are undersaturated. The dolomite re-action front can be viewed as a “terrestrial lysocline”, i.e. a depth in the soil where the rate of dissolution of carbonate changes rapidly, like the oceanic lysocline. But contrary to the downward-flowing soil solutions, seawater is oversaturated with respect to carbonate phases (either calcite or aragonite) above the oceanic lysocline and undersaturated at depth. With the ongoing oceanic acidification process driven by the rising CO₂, the oceanic lysocline is moving upward, reducing the thickness of the ocean layer saturated with respect to carbonate minerals (Archer et al., 1997; Boudreau et al., 2010).

However, the two lysoclines are coupled. When the ocean
lysocline shallows, the terrestrial lysocline is deepening such
that after some time lag, a new steady-state ocean pH and at-
mospheric CO₂ concentration must be achieved. In effect,
the enhanced continental carbonate dissolution produces a
flux of alkalinity to the oceans which eventually can counter-
act acidification (Archer, 2005). In other words, under some
conditions, the depth of carbonate depletion on the continents
would increase to compensate for the thinning of the carbon-
ate stability zone in the oceans. Such a feedback is common
in studies of the geological evolution of the climate, where
the surficial ocean stays oversaturated with respect to car-
bonate despite atmospheric CO₂ levels reaching 10 times or
more the present-day level. The question is whether this feed-
back can work at the secular or millennial timescale. The an-
swer to this question is not obvious, since carbonate dissolu-
tion (here dolomite) is heavily dependent on the temperature
and drainage evolution into the future, and this latest param-
eter displays contrasting behavior along the Mississippi tran-
sect.

Model calculations for pedons over the entire transect document that the dolomite front deepens 1.9 cm over the 1950–2100 period in the northern part of the Mississippi (40.9° N), and 6.1 cm over the same time period in the central part of the transect (34.2° N) (the greatest value observed for the transect). These rates correspond to a retreat rate of 13 cm kyr⁻¹ in the north, and 41 cm kyr⁻¹ in the central part. The amplitude of the retreat is approximately constant (≈ 5 to 6 cm retreat over 150 yr) from 30° N to 37° N. Above 37° N, the dolomite front retreat rate suddenly decreases below 2 cm over the 150 yr of the study. This sudden drop correlates well with a marked decrease in rainfall around 37° N. The averaged rainfall over the whole period of simulation falls from about 105–130 mm month⁻¹ below 37° N to 80 mm month⁻¹ above this latitude.

This rate of retreat of the front is correlated to the local
climatic conditions, but does not carry any information about
the response of the terrestrial lysocline to climate change. But
the question remains: does the rate of retreat accelerate or not

over the next century? We first calculated the time derivative
of the depth of the terrestrial lysocline for 9 pedons. These
rates of retreat display ample fluctuations as a function of
time for each location. We isolate the long-term trends by
fitting these rates (in cm yr⁻¹) with a linear fit. The slope of
this linear fit (cm yr⁻²) is a measurement of the long-term ac-
celeration or deceleration of the retreat rate of the terrestrial
lysocline (Fig. 8). Over the 1950–2100 period, the acceler-
ation of the retreat rate is maximum in the central part of
the transect. In contrast, it is at a minimum, reaching a neg-
ative value, in the north, suggesting that the dolomite front
is still retreating downward, but at a rate that decreases over
the 1950–2100 period. In the south, the downward shift of
the dolomite front is decelerating, but at a moderate rate.

The causes for this contrasting behavior are multiple. In the north, temperature is increasing at the fastest rate for any of the sites along the transect, soil CO₂ is decreasing rapidly, and drainage is decreasing. All these factors tend to slow down the dolomite dissolution and hence the downward rate of retreat of the terrestrial lysocline. In the central part of the watershed, in contrast, the drainage increase is at a maximum, accelerating the rate of retreat of the dolomite front. In the south, processes are more intertwined, with the rate of increase in soil CO₂ playing an important role for this environment where the changes in drainage and temperature are similar. It is important to note that the drainage rise in the extreme south is not large enough to compensate for the temperature rise. As a consequence, the dolomite front is retreating there at a decreasing rate for the 1950–2100 period.

We note again that the initial position of the dolomite front was fixed at 2.8 m depth for all pedons. Moving this initial position upward or downward does not alter significantly our calculated retreat rates. Indeed, dolomite front retreat rates are mostly affected by changes in pH (Fig. 6) and drainage, and these changes largely occur within the root zone. pH and drainage remain fairly constant with depth below the root zone (until the dolomite front where pH rapidly rises). As long as the dolomite front is located below the root zone, which is the case for all pedons according to Williams et al. (2010), the geographical pattern of the terrestrial lysocline retreat rate is not dependent on the initial dolomite front depth.

These factors (drainage, temperature and soil CO₂) thus
produce complex behavior that cannot easily be integrated
or averaged over the transect to answer the question we have
asked. In the central part of the watershed, dolomite weather-
ing is accelerating, while it is slowing rapidly in the north and
moderately in the extreme south. The response of the terres-
trial lysocline to the ongoing climate changes is not straight-
forward. Indeed, Fig. 8 illustrates the high sensitivity of the
terrestrial lysocline deepening to parameters which are them-
selves sensitive to the modality of climate change (drainage,
temperature and soil CO₂ levels).
<figure><img src="image_11.png" /></figure>

**Fig. 8.** Increasing rates of vertical drainage (below the root zone), ground temperature, and soil CO$_2$ level over the Mississippi transect, simulated for the 1950–2100 period. In red, acceleration of the dolomite retreat rate along the transect simulated over the same period. The dolomite front is moving downwards at each location, but this retreat rate either accelerates (positive value) or decelerates (negative value). Maximum acceleration of the retreat rate is reached in the central section of the Mississippi Valley, owing to the rate of change in climatic and geochemical conditions.

# 6 Limitations and perspectives

Contrary to the temperature and soil hydrology, the soil CO$_2$ level does not emerge amongst the critical driving factors of weathering in the future in our simulations. However, our model for soil CO$_2$ is rather basic. Diffusion only depends on air temperature (Van Bavel, 1951; Goddéris et al., 2010). The impact of soil hydrology on diffusion is not included, e.g. the expected decrease in air-filled porosity with increasing water content, which would limit diffusion back to the atmosphere when the profiles get more humid. A better model should be implemented in the future. However, a striking output of the simulation is the decrease in soil CO$_2$ in the future for the northern pedon, as diffusion rises with temperature. As the simulated profile also dries out, upward diffusion would be further promoted. Thus the soil CO$_2$ trend will still display a decrease into the future even if the role of hydric state on diffusion is included.

Another limitation arises when considering the terrestrial lysocline behavior. In our simulation, we do not allow dolomite to precipitate when solutions become supersaturated. Indeed, the kinetics of dolomite precipitation is very low at ambient temperature. Globally, dolomite precipitation is usually only observed where promoted by bacterial activity in specific environments (Vasconcelos et al., 1995; Sánchez-Román et al., 2011), but calcite is not allowed to precipitate either. Supersaturation might occur in the future, especially since the weathering profiles get warmer and possibly dry out, as shown in the simulations. If calcite were to precipitate, the terrestrial lysocline might cease deepening and actually move upward, consuming alkalinity and releasing CO$_2$ back to the atmosphere. Although atmospheric carbon would be trapped in mineral form if calcite precipitates, there would be less capture of carbon in the dissolved HCO$_3^-$ form. This latter capture stores carbon in the world oceans at the millen-

nial timescale (Beaulieu et al., 2012). An upward movement of the terrestrial lysocline would negatively impact the short term atmospheric CO$_2$ consumption by weathering.

Here, below the dolomite front in both north and south sites, the saturation state of the soil solution with respect to calcite is predicted to increase in the future (Fig. 9). In the south, the mean annual saturation state of the soil solution with respect to calcite below the dolomite front is seldom above 1 before about 2015. Then the occurrence of saturation increases rapidly after 2015. In contrast in the north, saturation of the soil solution stays generally below 1, because temperature is lower, although it increases with time. Although not critical here, future work should account for calcite precipitation.

Finally, our work shows that the weathering of carbonate lithologies may respond quickly to climate change in the near future. Because of the very high reactivity of carbonates, their dissolution is mainly controlled by transport, stressing the need for an accurate modeling of the water fluxes in these environments. Modeling of carbonate precipitation, which may occur more often in a drier and warmer future world, would require accurate knowledge of the carbonate precipitation kinetics close to equilibrium.

# 7 Conclusions

We use a cascade of numerical models (simulating climate, the continental biosphere, and the weathering processes) to extrapolate the weathering of the Mississippi Valley loess forward in time to 2100. Our interest in the weathering of those loessic formations arises because they are globally important due to their high surface area and due to the fact that they contain both carbonate and silicate mineralogies, thus exemplifying the contrasting behaviors of weathering
<figure><img src="image_12.png" /><figcaption>Fig. 9. Calculated saturation state of the soil solutions with respect to calcite plotted as a function of depth and time: (a) south pedon, (b) north pedon.</figcaption></figure>

for these two dominant rock-forming mineral types. However, loess is also of importance in that it represents a fast-weathering component of the terrestrial surface that is produced in response to glaciations worldwide (Anderson et al., 1997). The simulations we present here are an outgrowth of our previous investigation of the Peoria loess deposited along the Mississippi valley of the USA interpreted as a natural experiment for pedogenesis over the last 13 kyr (Williams et al., 2010; Goddéris et al., 2010).

We found that yearly averaged drainage is expected to slightly increase in the southern end of the Mississippi Valley (14% increase), while it will decrease drastically in the northern end (40% decrease). The decrease in the north is triggered by the rise in evapotranspiration driven by enhanced temperatures. Coupled to the rise in temperature linked to the rise in CO₂ (700 ppmv in 2100), the dolomite dissolution rate slows down in the northern and southern ends of the Mississippi Valley. Dolomite dissolution is further reduced by the accelerated diffusion of soil CO₂ to the atmosphere under a warmer climate. Conversely, primary feldspar dissolution rates are enhanced by the temperature rise.

Dry events occur more often after about 2025, particularly in the north, forcing pH to rise by 0.3 units in the weathering profile. Smectite precipitates during very low drainage events.

When simulating the whole Mississippi transect (9 sites, from 30.7° to 42.2° N), we found that the dolomite front retreat decelerates at both ends of the transect, but accelerates in the middle part of the watershed. This area experiences the strongest enhancement in vertical drainage as dictated by the delicate balance between rainfall and evapotranspiration along the transect. In the middle part of the watershed, this balance maximizes the drainage increases. North of 37.4° N, the drainage rapidly decreases, forcing the dolomite retreat rate to slow down into the future.

The dolomite front defines a terrestrial lysocline for the Mississippi Valley loess. The deepening of this lysocline defines the way carbonate weathering responds to climate change. The deepening accelerates only in the central part of the watershed, and decelerates elsewhere. Acceleration of the retreat of the terrestrial lysocline is equivalent to increasing rates of delivery of alkalinity to the oceans. Defining how this lysocline will respond at continental or global scales is thus a key issue regarding predictions of impacts of processes such as ocean acidification that are occurring as a result of increasing atmospheric CO₂.

# Appendix A

## Model details

### A1 WITCH general description

WITCH is a vertical box model that calculates the time evolution of the soil solutions as a function of time. Soil solution chemical composition is calculated for 40 vertical levels extending from the surface down to 4 m depth. The vertical resolution is 0.2 m in the upper horizons, and decreases below the dolomite weathering front to match the fluctuations in the advance rate of this front correctly. A differential equation describes the budget of the main conservative species (not involved in fast speciation reactions) and is solved at each time step for each box:

$$ \frac{d(z \cdot \theta \cdot C_j)}{dt} = F_{\text{top}} - F_{\text{bot}} + \sum_{i=1}^{n_m} F_{\text{weath.}}^i - \sum_{i=1}^{n_m} F_{\text{prec.}}^i, \quad (A1) $$

where $C_j$ is the concentration of the species j. WITCH includes a budget equation for Ca²⁺, Mg²⁺, K⁺, Na⁺, SO₄²⁻, total alkalinity, total aluminum, total silica, and total phosphorus. z is the thickness of the layer under consideration,
and $\theta$ is the water volumetric content. $F_{top}$ is the input of species from the layer above the given layer in moles per unit time through drainage, and $F_{bot}$ is the removal of water through downward drainage. $n_m$ is the total number of minerals considered in the simulation. $R_j$ is the exchange term between the soil solution and the clay–humic complex. $F_{weath, j}^i$ is the release of species j to solution through dissolution of the mineral i in the box, and $F_{prec, j}^i$ is the removal of j through precipitation when the solution becomes supersaturated with respect to mineral i. The dissolution $F_{weath, j}^i$ and precipitation $F_{prec, j}^i$ terms for silicate minerals are calculated using kinetic laws and parameters derived from transition state theory (TST; Eyring, 1935) and laboratory experiments, respectively (Schott et al., 2009; Goddéris et al., 2010):

$$F_g = A_g' \left[ \sum_l k_l \cdot g \cdot \exp\left(\frac{-E_{a, g}^l}{R \cdot T}\right) \cdot a_l^{n_l \cdot g} \cdot f_{inh} \right] \cdot (1 - \Omega_g^{1/s}), \quad (A2)$$

where $F_g$ is the overall dissolution rate of mineral g inside a given layer. The sum accounts for the four parallel rate-controlling elementary reactions that are assumed to describe dissolution promoted by H⁺, OH⁻, H₂O and organic ligands. In these relations, $a_l$ and $n_l \cdot g$ stand for the activity and the reaction order with respect to the i-th species, respectively. For organic ligands, $a_l$ equals the activity of a generic organic species RCOO⁻ (conjugate of RCOOH). $k_l \cdot g$ is the rate constant of mineral g dissolution reaction promoted by species l, and $E_{a, g}^l$ is the activation energy of this reaction. $f_{inh}$ stands for inhibitory effects (i.e. by aqueous Al). The last factor in Eq. (A2) describes the effect on the rate of the departure from equilibrium, where $\Omega_g$ is the solution saturation state with respect to mineral g ($\Omega_g = Q/K_g$ where Q is the activity quotient and $K_g$ is the equilibrium constant for mineral g dissociation reaction), and s, the Temkin's coefficient number, is the stoichiometric number of moles of activated complex formed from one mole of the mineral. Precipitation terms are calculated in the same way, but $\Omega_g$ is then > 1 (only clay secondary phases are allowed to precipitate). $A_g'$ is the reactive surface. Dolomite dissolution has been added to WITCH and modeled within the framework of TST and surface coordination chemistry concepts (Pokrovsky and Schott, 1999, 2001; Pokrovsky et al., 1999). Dolomite dissolution rate $F_{dol}$ is modeled as:

$$F_{dol} = \left[ k_H^{\text{dol}} \cdot a_H^{0.75} + k_M^{\text{dol}} \cdot \left( \frac{1.575 \cdot 10^{-9}}{1.575 \cdot 10^{-9} + 3.5 \cdot 10^{-5} \cdot a_{\text{CO}_3} + a_{\text{CO}_3} \cdot a_{\text{Ca}}} \right) \right] \cdot (1 - \Omega_{\text{dol}}^{1.9}), \quad (A3)$$

where $a_H$, $a_{\text{CO}_3}$, $a_{\text{Ca}}$ stands respectively for the H⁺, CO₃²⁻, and Ca²⁺ activities. $k_H^{\text{dol}}$ equals 10⁻³ mol m⁻² s⁻¹, and $k_M^{\text{dol}}$ 10⁻⁸.² mol m⁻² s⁻¹ at 25 °C. The activation energy for $k_H^{\text{dol}}$ and $k_M^{\text{dol}}$ is set to 40 and 60 kJ mol⁻¹, respectively (Pokrovsky and Schott, 1999, 2001; Pokrovsky et al., 1999). Dolomite is not allowed to precipitate.

Once the mass balance is calculated, the speciation of the soil solutions is calculated at each time step together with pH. In addition, WITCH includes a mass balance for each mineral phase, accounting for the amount of mineral being dissolved and precipitated within each time step.

## A2 Soil hydrology

The water volumetric content $\theta$ at each weathering profile level and the vertical drainage as a function of depth are taken for each location from the simulation of the global dynamic vegetation model CARAIB, as mentioned in the main text (Dury et al., 2011). CARAIB accounts for only one level below ground, spanning the whole root zone, while the vertical resolution of WITCH is much finer. Regarding $\theta$, we fix the water volumetric content of the first layer of WITCH below ground at the water volumetric content calculated by CARAIB. At the base of the root zone (a depth which is calculated at each time step by CARAIB), we fix $\theta$ at the field capacity. Inside the root zone, the water volumetric content is linearly interpolated.

Regarding the vertical drainage, the water flux at the ground level and entering the weathering profile is fixed at the rainfall calculated by the ARPEGE GCM minus the direct evaporation. At the base of the root depth, the drainage in WITCH is fixed at the drainage calculated by CARAIB, accounting for the water uptake by the roots. Below that level, drainage is held constant. Inside the root zone, drainage is linearly interpolated at each time step. Apart from the superficial runoff at the top of the weathering profile, the model does not account for horizontal transfer inside the profile. Once entering the profile, the water is either taken up by roots or transferred downward.

## A3 Calculating soil CO₂ levels

The CARAIB model calculates autotrophic and heterotrophic respiration. We assume that all the heterotrophic respiration and 1/3 of the autotrophic respiration occurs below ground. The sum of these two terms defines the CO₂ production inside the root zone, the thickness of which is calculated by CARAIB. This production allows calculation of soil CO₂ at the base of the root zone, accounting for the temperature effect on CO₂ upward diffusion, as detailed in Goddéris et al. (2010). We then assume that the CO₂ pressure at ground level (0 m depth) equals the atmospheric pressure. CO₂ is assumed to increase linearly inside the root zone until it reaches its maximum at the base of the root zone. Below this level, it is held constant.

Edited by: J. Middelburg
<figure><img src="image_13.png" /></figure>

The publication of this article is financed by CNRS-INSU.

## References

Amiotte-Suchet, P., Probst, J.-L., and Ludwig, W.: Worldwide distribution of continental rock lithology: implications for the atmospheric/soil CO$_2$ uptake by continental weathering and alkalinity river transport to the oceans, Global Biogeochem. Cy., 17, 1038, doi:10.1029/2002GB001891, 2003.

Anderson, S. P., Drever, J. I., and Humphrey, N. F.: Chemical weathering in glacial environments, Geology, 25, 399–402, 1997.

Archer, D.: Fate of fossil fuel CO$_2$ in geologic time, J. Geophys. Res., 110, C09S05, doi:10.1029/2004JC002625, 2005.

Archer, D., Kheshgi, H., and Maier-Reimer, E.: Multiple timescales for neutralization of fossil fuel CO$_2$, Geophys. Res. Lett., 24, 405–408, 1997.

Beaulieu, E., Goddéris, Y., Labat, D., Roelandt, C., Calmels, D., and Gaillardet, J.: Modeling of water-rock interaction in the Mackenzie basin: competition between sulfuric and carbonic acids, Chem. Geol., 289, 114–123, 2011.

Beaulieu, E., Goddéris, Y., Donnadieu, Y., Labat, D., and Roelandt C.: High sensitivity of the continental-weathering carbon dioxide sink to future climate change, Nature Climate Change, 2, 346–349, 2012.

Bettis, E. A., Muhs, D. R., Roberts, H. M., and Wintle, A. G.: Last glacial loess in the conterminous USA, Quaternary Sci. Rev., 22, 1907–1946, 2003.

Boudreau, B. P., Middelburg, J. J., Hofmann, A. F., and Meysman, F. J. R.: Ongoing transients in carbonate compensation, Glob. Biogeochem. Cy., 24, GB4010, doi:10.1029/2009GB003654, 2010.

Brantley, S. L. and White, A. F.: Approaches to modeling weathered regolith, in: Thermodynamics and Kinetics of Water-Rock Interaction, edited by: Oelkers E. H. and Schott, J., Rev. Min. Geochem., 70, 435–484, 2009.

Dessert, C., Dupré, B., Gaillardet, J., François, L. M., and Allègre, C. J.: Basalt weathering laws and the impact of basalt weathering on the global carbon cycle, Chem. Geol., 202, 257–273, 2003.

Dury, M., Hambuckers, A., Warnant, P., Henrot, A., Favre, E., Ouberdous, M., and François L.: Responses of European forest ecosystems to 21st century climate: assessing changes in interannual variability and fire intensity, iForest, 4, 82–99, doi:10.3832/ifor0572-004, 2011.

Eyring, H.: The activated complex in chemical reactions, J. Chem. Phys., 3, 107–115, 1935.

Fox-Rabinovitz, M., Côté, J., Dugas, B., Déqué, M., and McGregor, J. I.: Variable resolution general circulation models: stretched-grid model intercomparison project (SGMIP), J. Geophys. Res., 111, D16104, doi:10.1029/2005JD006520, 2006.

François, L. M., Utescher, T., Favre, E., Henrot, A.-J., Warnant, P., Micheels, A., Erdei, B., Suc, J. P., Cheddadi, R., and Mosbrugger, V.: Modelling Late Miocene vegetation in Europe: Results of

the CARAIB model and comparison with palaeovegetation data, Paleogeogr. Paleoclim., 304, 359–378, 2011.

Gaillardet, J., Dupré, B., Louvat, P., and Allègre, C. J.: Global silicate weathering and CO$_2$ consumption rates deduced from the chemistry of large rivers, Chem. Geol., 159, 3–30, 1999.

Gibelin, A. L. and Déqué, M.: Anthropogenic climate change over the Mediterranean region simulated by a global variable resolution model, Clim. Dynam., 20, 327–339, 2003.

Gislason, S. R., Oelkers, E. H., Eiriksdottir, E. S., Kardjilov, M. I., Gisladottir, G., Sigfusson, B., Snorrason, A., Elefsen, S., Hardardottir, J., Torssander, P., and Oskarsson, N.: Direct evidence of the feedback between climate and weathering, Earth Planet. Sci. Lett., 277, 213–222, 2008.

Goddéris, Y., François, L. M., Probst, A., Schott, J., Moncoulon, D., Labat, D., and Viville, D.: Modelling weathering processes at the catchment scale: the WITCH numerical model, Geochim. Cosmochim. Ac., 70, 1128–1147, 2006.

Goddéris, Y., Roelandt, C., Schott, J., Pierret, M.-C., François, and L. M.: Towards an integrated model of weathering, climate, and biospheric processes, in: Thermodynamics and Kinetics of Water-Rock Interaction, edited by: Oelkers E.H. and Schott J., Rev. Mineral. Geochem., 70, 411–434, 2009.

Goddéris, Y., Williams, J. Z., Schott, J., Pollard, D., and Brantley, S. L.: Time evolution of the mineralogical composition of the Mississippi Valley loess over the last 10 kyr: climate and geochemical modeling, Geochim. Cosmochim. Ac., 74, 6357–6374, 2010.

Haase, D., Finke, J., Haase, G., Ruske, R., Pecsi, M., Richter, H., Altermann, M., and Jäger, K. D.: Loess in Europe – its spatial distribution based on a European Loess Map, scale 1:2,500,000, Quaternary Sci. Rev., 26, 1301–1312, 2007.

Hartmann, J., Jansen, N., Dürr, H. H., Kempe, S., and Köhler, P.: Global CO$_2$-consumption by chemical weathering: what is the contribution of highly active regions?, Global Planet. Change, 69, 185–194, 2009.

Le Hir, G., Donnadieu, Y., Goddéris, Y., Meyer-Berthaud, B., Ramestein, G., and Blakey, R. C.: The climate change caused by the land plant invasion in the Devonian, Earth Planet. Sci. Lett., 310, 203–212, 2011.

Moosdorf, N., Hartmann, J., Lauerwald, R., Hagedorn, B., and Kempe, S.: Atmospheric CO$_2$ consumption by chemical weathering in North America, Geochim. Cosmochim. Ac., 75, 7829–7854, 2011.

Muhs, D. R., Bettis, E. A., Been, J., and McGeehin J. P.: Impact of climate and parent material on chemical weathering in loess-derived soils of the Mississippi River Valley, Soil Sci. Soc. Am. J., 65, 1761–1777, 2001.

Muhs, D. R., Bettis, E. A., Aleinikoff, J. N., McGeehin, J. P., Beann, J., Skipp, G., Marshall, B. D., Roberts, H. M., Johnson, W. C., and Benton, R.: Origin and paleoclimatic significance of late Quaternary loess in Nebraska: evidence from stratigraphy, chronology, sedimentology, and geochemistry, Geol. Soc. Am. Bull., 120, 1378–1407, 2008.

Nemry, B., François, L. M., Warnant, P., Robinet, F., and Gérard, J.-C.: The seasonality of the CO$_2$ exchange between the atmosphere and the land biosphere: a study with a global mechanistic vegetation model, J. Geophys. Res., 101, 7111–7125, 1996.

Nemry, B., François, L. M., Gérard J.-C., Bondeau A., and Heimann, M.: Comparing global models of terrestrial net pri-
mary productivity (NPP): analysis of the seasonal atmospheric CO$_2$ signal., Glob. Change Biol., 5, Supplement 1, 65–76, 2003.

Oliva, P., Dupré, B., and Viers, J.: Chemical weathering in granitic environments, Chem. Geol., 202, 225–256, 2003.

Otto, D., Rasse, D., Kaplan, J., Warnant, P., and François, L. M.: Biospheric carbon stocks reconstructed at the Last Glacial Maximum: comparison between general circulation models using prescribed and computed sea surface temperatures, Global Planet. Change, 33, 117–138, 2002.

Pokrovsky, O. S. and Schott, J.: Processes at the magnesium-bearing carbonates/solution interface. II. Kinetics and mechanism of magnesite dissolution, Geochim. Cosmochim. Ac., 63, 881–897, 1999.

Pokrovsky, O. S. and Schott, J.: Kinetics and mechanism of dolomite dissolution in neutral to alkaline solutions revisited, Am. J. Sci., 301, 597–626, 2001.

Pokrovsky, O. S., Schott, J., and Thomas, F.: Dolomite surface spe-
ciation and reactivity in aquatic systems, Geochim. Cosmochim.
Ac., 63, 3133–3143, 1999.

Pye, K. and Johnson, R.: Stratigraphy, geochemistry, and thermolu-
minescences ages of lower Mississippi Valley loess, Earth Surf.
Proc. Land., 13, 103–124, 1988.

Raymond, P. A., Oh, N.-H., Turner, R. E., and Broussard, W.: An-
thropogenically enhanced fluxes of water and carbon from the
Mississippi River, Nature, 451, 449–452, 2008.

Roelandt, C., Goddéris, Y., Bonnet, M.-P., and Sondag, F.: Coupled modeling of biospheric and chemical processes at the continental scale, Global Biogeochem. Cy., 24, GB2004, doi:10.1029/2008GB003420, 2010.

Rousseau, D.-D.: Loess biostratigraphy: new advances and ap-
proaches in mollusk studies, Earth Sci. Rev., 54, 157–171, 2001.

Ruhe, R. V.: Loess-derived soils, Mississippi Valley region: II, Soil-
climate system, Soil Sci. Soc. Am. J., 48, 864–867, 1984.

Salas y Mélia, D., Chauvin, F., Déqué, M., Douville, H., Guérémy, J.-F., Marquet, P., Planton, S., Royer, J.-F., and Tyteca, S.: Description and validation of CNRM-CM3 global coupled climate model, Note de centre GMCEC, CNRM, 103, 2005.

Sánchez-Román, M., Romanek, C. S., Fernández-Remolar, D. C., Sánchez-Navas, A., McKenzie, J. A., Amils Pibernat, R., and Vasconcelos, C.: Aerobic biomineralization of Mg-rich carbonates: implication for natural environments, Chem. Geol., 281, 143–150, 2011.

Schott, J., Pokrovsky, O. S., and Oelkers, E. H.: The link between mineral dissolution/precipitation kinetics and solution chemistry, in: Thermodynamics and kinetics of water-rock interaction, edited by: Oelkers, E. H. and Schott, J., Rev. Mineral Geochem. 70, 207–258, 2009.

Solomon, S., Qin, D., Manning, M., Chen, Z., Marquis, M., Averyt, K. B., Tignor, M., and Miller, H. L.: Contribution of working group I to the Fourth Assessment Report of the Intergovernmen- tal Panel on Climate Change, Cambridge University Press, Cam- bridge, United Kingdom and New York, NY, USA, 2007.

Taylor, L. L., Banwart, S. A., Valdes, P. J., Leake, J. R., and Beer-
ling, D. J.: Evaluating the effects of terrestrial ecosystems, cli-
mate and carbon dioxide on weathering over geological time: a
global-scale process-based approach, Philos. T. Roy. Soc. B, 367,
565–582, 2012.

Van Bavel, C. H. M.: A soil aeration theory based on diffusion, Soil Sci., 72, 33–46, 1951.

Vasconcelos, C., McKenzie, J. A., Bernasconi, S., Grujic, D., and Tien, A. J.: Microbial mediation as a possible mechanism for natural dolomite formation at low temperatures, Nature, 377, 220–222, 1995.

Violette, A., Goddéris, Y., Maréchal, J.-C., Riotte, J., Oliva, P., Mo-
han Kumar, M. S., Sekhar, M., and Braun, J.-J.: Modelling the
chemical weathering fluxes at the watershed scale in the Trop-
ics (Mule Hole, South India): relative contribution of the smec-
tite/kaolinite assemblage versus primary minerals, Chem. Geol.,
277, 42–60, 2010.

Warnant, P., François, L. M., Strivay, D., and Gérard, J.-C.: CARAIB: a global model of terrestrial biological productivity, Global Biogeochem. Cy., 8, 255–270, 1994.

White, A. F. and Blum, A. E.: Effects of climate on chemical weath-
ering in watersheds, Geochim. Cosmochim. Ac., 599, 1729–
1747, 1995.

White, A. F., Bullen, T. D., Schultz, M. S., Blum, A. E., Huntington, T. G., and Peeters N. E.: Differential rates of feldspar weathering in granitic regoliths, Geochim. Cosmochim. Ac., 65, 847–869, 2001.

Williams, J. Z.: The effect of temperature and precipitation on sodium depletion fronts in soils developed on Peoria loess, Pennsylvania State University, University Park, 113 pp., 2008.

Williams, J. Z., Bandstra, J. Z., Pollard, D., and Brantley S. L.: The temperature dependence of feldspar dissolution determined using a coupled weathering-climate model for Holocene-aged loess soils, Geoderma, 156, 11–19, 2010.

West, A. J., Galy, A., and Bickle, M.: Tectonic and climatic controls on silicate weathering, Earth Planet. Sci. Lett., 235, 211–228, 2005.