The extended timescales involved in the decay of radioactive wastes to safe levels mean that geological disposal facilities must continue to function effectively long into the future. It is therefore essential to consider long-term climate evolution in post-closure performance assessments in order to evaluate a geological disposal system's response and robustness to a variety of potential environmental changes, driven by both natural and anthropogenic forcings. In this paper, we illustrate the multiple decay components that characterize the primary driver of climate change – atmospheric CO2 – in response to fossil fuel carbon emissions. We perform a multi-exponential analysis on a series of atmospheric CO2 decay curves predicted by an Earth system model and create an empirical response function that encapsulates the long-term (>1 kyr) removal of excess CO2 from the atmosphere. We present this response function as a simple tool for rapidly projecting the future atmospheric CO2 concentration resulting from any plausible cumulative release of CO2. We discuss the implications of the long ‘tail’ to this atmospheric CO2 decay curve, both in terms of future climate evolution as well as potential impacts on radioactive waste repositories.


The significant timescales involved in the decay to safe levels of radionuclides incorporated in radioactive wastes means that geological disposal facilities containing low- and intermediate-level wastes must continue to function effectively for up to 100,000 years, as in the case of the Low Level Waste Repository (LLWR) in the UK (LLWR, 2011). For high-level wastes and spent nuclear fuel, repositories such as the proposed KBS-3 facility at Forsmark in Sweden (SKB, 2011) must remain functional for up to 1 million years. It is therefore essential to consider long-term climate evolution in post-closure performance assessments in order to evaluate a geological disposal system's response to and robustness against a variety of potential environmental changes, driven by both natural (e.g. orbital variations) and anthropogenic (e.g. fossil fuel emissions) forcings.

In order to identify what climate changes may be of relevance to different waste repositories and how they might be impacted by these changes, facilities can be categorized based on various criteria. These categories include the mode of construction (e.g. excavated from the surface or at depth), the geological context, the hydrogeological context (i.e. saturated or unsaturated), the coastal context and the potential for extreme climates. For repositories in glaciated environments for example, the environmental change that is likely to have the most significant impact is the passage of the ice sheet margin across the site (Becker et al., 2014). Whether the accompanying changes in the hydrological regime are likely to have a large impact can be inferred from the hydrogeological and geological context, while the mode of construction will give an indication of how vulnerable the facility is to surface denudation and permafrost formation. For high-level waste repositories the timing and severity of the glaciation may be of particular importance, as this will determine the remaining radionuclide inventory that can be impacted by the altered hydrogeological conditions that occur during passage of the ice sheet. It will also determine the duration of the period of permafrost conditions prior to glaciation and, therefore, be a key factor in determining whether the permafrost develops to repository depth (SKB, 2011). Thus, different repositories have different key timescales that need to be considered, based in part on the climate change processes that are expected to have the most significant impact on the facility. In the case of the UK's LLWR for example, the coming centuries will be key as sea level rise and the resultant coastal erosion are expected to affect the site on these timescales (Fish et al., 2010), whereas permafrost is projected to impact Sweden's shallow disposal facility for operational and decommissioning wastes (SFR) on timescales of tens of millennia (SKB, 2013). In contrast, the key timescale for the proposed KBS-3 repositories in Sweden and Finland is more than 100 kyr. In this context, the timing of the next glacial inception and amplitude of future glacial-interglacial cycles are key.

In general, the primary driver of global and regional climate change on timescales from decades to >1 Myr is variation in the concentration of atmospheric greenhouse gases and principally, carbon dioxide (CO2). It is in this way that anthropogenic activities are now increasingly impacting the climate system (IPCC, 1990, 1995, 2001, 2007). Between 1750 and 2011 ~555 ± 85 PgC of CO2 were emitted, raising atmospheric CO2 concentrations by 40% from 278 to 390.5 ppm in 2011 (Ciais et al., 2013). This concentration is now higher than has occurred at any time in at least the last 800 kyr as recorded in ice cores (Luthi et al., 2008) and equal to the maximum values reached about 2–3 million years ago (Martinez-Boti et al., 2015).

A sound understanding of the global carbon cycle and associated controls on the removal of anthropogenic CO2 from the atmosphere, particularly geological processes with timescales comparable to the lifetimes of some of the longer-lived major radionuclides, is hence essential. This paper focuses on the range of processes that comprise the long-term carbon cycle and how they affect atmospheric CO2, and presents a simple empirical representation of the decay of an atmospheric CO2 perturbation. We also discuss the implications of the projected atmospheric CO2 response for climate and post-closure performance assessments. This work contributes to that of Working Group 6 of the International Atomic Energy Agency (IAEA) MODARIA Programme, which has been tasked with developing a common framework for addressing climate change in post-closure radiological assessments of solid radioactive waste disposal.

The fate of fossil fuel CO2

As CO2 is emitted to the atmosphere via the combustion of fossil fuels and land-use change, its concentration is elevated above the long-term average which, for the last ten thousand years of the Holocene, lies in the range ~260–280 ppm. This excess anthropogenic CO2 is subsequently removed by a range of carbon cycle processes operating on a variety of timescales, summarized in Fig. 1. These are:

  1. On the shortest timescales of months to years (Fig. 1-I), gaseous CO2 is mixed between the two hemispheres of the atmosphere and transferred across the air–sea interface of the oceans, dissolving in surface waters and decreasing seawater pH (Turley et al., 2010). The net outcome of this process is approximated by reaction 1 in Fig. 1-I. CO2 fertilization of terrestrial primary productivity also results in a further uptake of CO2 from the atmosphere (Sitch et al., 2008). Climate feedbacks modify the rate of uptake by reducing the solubility of CO2 in seawater (warming) as well as mixing (ocean stratification).

  2. Over decades to centuries, surface ocean waters rich in dissolved CO2 are transported down into the ocean interior, primarily through the formation of intermediate and deep waters at high latitudes (Archer and Brovkin, 2008; Archer et al., 1998). In terrestrial environments, a new equilibrium will be established with enhanced primary productivity driving a larger soil carbon pool, illustrated in Fig. 1-II. Climate feedbacks again modify the strength of carbon uptake from the atmosphere for instance via soil warming, with increased metabolic activity driving an increased rate of return of CO2 back to the atmosphere.

  3. The chemical impact of dissolving excess CO2 in seawater results in a reduction in carbonate (CO32) concentration in ocean waters causing a reduction in the stability of biogenic calcium carbonate (CaCO3) minerals (principally calcite). Sedimentary CaCO3 that was originally formed at the surface and transported by gravitational settling to the seafloor may now start to dissolve. This process (reaction 2 in Fig. 1-III) is known as seafloor neutralization (Archer et al., 1997) and has the effect of partially regenerating CO32 concentrations and hence the buffering capacity of the ocean, allowing additional draw-down of atmospheric CO2 via reaction 1 (Ridgwell and Hargreaves, 2007). This all occurs on timescales of millennia and can be modified by any changes in ocean circulation induced by a warming climate.

  4. At the same time as rates of CaCO3 burial in marine sediments are reduced (by process #3 above), the supply of solutes from terrestrial rock weathering continues. This creates an imbalance in the supply vs. removal of Ca2+ and HCO3 to the ocean that on multi-millennial timescales results in further uptake of CO2 (Lenton and Britton, 2006; Ridgwell and Hargreaves, 2007). This is known as ‘terrestrial CaCO3 neutralization’ and is illustrated in Fig. 1-IV. Here climate feedbacks, by enhancing the weathering of carbonate rocks (e.g. via more acidic rainwater), will act to increase the recovery of ocean carbonate chemistry.

  5. The final process which occurs on a timescale of hundreds of thousands to millions of years – the weathering of terrestrial silicate rocks – will act to remove the remaining atmospheric CO2 perturbation and restore atmospheric CO2 back to its pre-industrial state (Lenton and Britton, 2006) via reaction 3 in Fig. 1-V.

A helpful parallel can be drawn between the decay of an atmospheric CO2 perturbation via the action of the various processes described above, and the radioactive decay of waste. Both types of decay can be characterized by a series of decay curves, each with a different characteristic timescale – each either associated with individual radionuclides or different carbon cycle processes. For CO2, the combined operating timescale and uptake capacity of each of these carbon cycle processes controls the atmospheric lifetime of the carbon perturbation, and hence the timescale and magnitude of the resultant climate changes. The radioactive isotopes of most interest in waste disposal exhibit a similar range of decay timescales to that of excess atmospheric CO2, ranging from decades (e.g. Cs-137 – 30.1 yr) to hundreds of years (e.g. Am-241 – 458 yr) to tens of thousands of years (e.g. Pu-239 – 24.1 kyr) to hundreds of thousands of years or more (e.g. U-234 – 245.5 kyr). This comparative similarity highlights the importance of consideration of long-term atmospheric CO2 in the future integrity of waste repositories.

Characterizing the timescales of atmospheric CO2 decay

The decay of an atmospheric CO2 perturbation, and hence the long-term carbon cycle processes, can be captured in the form of an empirical function comprising a series of exponential decay curves, fitted to a composite CO2 decay curve that is generated by means of a mechanistic numerical model. This approach has been undertaken in a number of past studies which produced response functions for specific sized CO2 releases (Archer et al., 1997; Colbourn, 2011), but is extended here by including a number of climate-carbon feedbacks on ocean circulation and terrestrial weathering and by using a larger range of CO2 emissions.

To characterize the timescales of CO2 decay, we utilize an Earth system model (cGENIE). This is based on the efficient climate model of Edwards and Marsh (2005) and incorporates a 2-D Energy-Moisture Balance atmosphere (EMBM), coupled to a 3-D frictional geostrophic ocean circulation model together with a dynamic-thermodynamic sea-ice component. The version of the cGENIE model used here also includes a representation of the global carbon cycle, including ocean cycling of carbon and nutrients (Ridgwell et al., 2007), geochemical interactions with marine sediments (Ridgwell and Hargreaves, 2007) and the weathering of carbonate and silicate rocks on land (Colbourn et al., 2013). The model is configured on a 36 × 36 equal horizontal grid with 8 vertical levels in the ocean. It is non-seasonally insolation forced (i.e. with an annual average latitudinal distribution of solar insolation) and because of the absence of a dynamical atmospheric circulation model component, requires that the spatial (and annual average) distribution of wind stress on the ocean surface and transport by circulation patterns, be prescribed (see Edwards and Marsh, 2005). The time step of the ocean circulation model component is approximately 7 days, with ocean and terrestrial carbon cycle processes updated a little less frequently (see Ridgwell et al., 2007) and deep-sea sediment composition only annually (Ridgwell and Hargreaves, 2007).

The Earth system model is sufficiently computationally efficient to undertake simulations of the order of 1 Myr, yet capable of representing to a first-order the basic responses of: (a) CO2 uptake by the ocean and interaction with marine sediments, (b) the global climate response to elevated CO2 and primary feedbacks of climate on ocean circulation and carbon cycling, and (c) weathering response to climate. However, there is no inter-annual variability in climate in this model and nor is the interaction of the terrestrial biosphere with CO2 and climate included. There is currently a great deal of uncertainty associated with the overall response of the terrestrial biosphere to future climate changes, such that even the sign of the change is not known for certain. There is also no representation of ice sheets in the model and associated feedbacks such as were important during the past glacial-interglacial cycles (Kohfeld and Ridgwell, 2009).

We perform a series of ten experiments, each of 1 Myr duration, in which idealized pulse emissions of CO2 are applied to the model. In these, all input of CO2 to the atmosphere occurs at a uniform rate and over a single model year. The total carbon release spans 1000 to 10,000 PgC, deemed to be a reasonable range of fossil fuel CO2 emissions that could be expected to be released in the next few hundred years. It is estimated that approximately 1000–2000 PgC remain in fossil fuel reserves that are presently economically exploitable (Ciais et al., 2013). However, fossil fuel resources (where economic extraction may be feasible in the future) are estimated to total 8543–13,649 PgC, while non-conventional resources such as methane clathrates could be as high as 20–25,000 PgC (Rogner, 1997). Pulse rather than time-dependent emissions have been used as previous studies have found that total emissions, rather than the rate of release, is the dominant control on the long-term response of atmospheric CO2, assuming emissions occur over no longer than about the next 300 years (Eby et al., 2009). This also appears to be true for the evolution of the ice sheets in response to fossil fuel emissions (Charbit et al., 2008). However, whilst a year-long pulse release facilitates the extraction of the underlying process timescales, our results will not be applicable to the evolution of atmospheric CO2 over the coming centuries as emissions continue and evolve in rate.

The long-term response of the climate to CO2 emissions of 1000 and 5000 PgC is shown in Fig. 2, as predicted by cGENIE. While the initial atmospheric CO2 peak subsides relatively quickly (Fig. 2a), it is clear that the overall perturbation is long-lived, with atmospheric pCO2 taking hundreds of thousands of years to return to pre-industrial levels, particularly for the 5000 PgC scenario. As a consequence, the wider issues associated with anthropogenic CO2 emissions, including warming (Fig. 2b), weakening of ocean circulation (Fig. 2c), ocean acidification (Fig. 2d) and increased terrestrial weathering rates expected due to an overall warmer and wetter climate (Fig. 2e) also last for hundreds of thousands of years.

In order to explore the different timescales of decay of the atmospheric CO2 perturbations, a multi-exponential function was fitted to the decay curves for the pulse experiments, taking the form:  
where CO2(t) is the atmospheric CO2 concentration (ppm) at a given time (t), E is the total emissions (PgC) which are assumed to be released at time t = 0, and n is the number of exponentials used in the fit. The factor 0.469 simply converts PgC of instantaneous CO2 emissions into the equivalent atmospheric concentration in ppm, and 278.0 is the assumed pre-industrial CO2 concentration. For each (i) of the exponentials; τi is the timescale of decay and Ai represents the fraction of total emissions removed from the atmosphere over this timescale.

An equation consisting of the sum of five exponentials was chosen as the Akaike Information Criterion (AIC) value, a measure of the relative quality of the fit which favours models with fewer parameters (a low value is best), indicated that this provided a better fit (had a lower value) than the four exponential model. Adding more exponentials does not improve the overall fit, as generally five different exponential curves are fitted, with the additional curve(s) being very similar to one of the other curves or else not having meaningful coefficient values.

The fitting coefficients (Ai and τi) were made linear functions of the total amount of carbon released (μ PgC) via a regression analysis, hence the full empirical function took the form:  
Using the parameters listed in Table 1, this empirical response function is capable of projecting the atmospheric decay of CO2 for a range of total emissions. For further details of the above methodology please see Colbourn et al. (2015).

Long-term future CO2

The dynamics of the global climate and carbon cycle system with its myriad of processes and attendant feedbacks cannot be accurately partitioned into just a few idealized decay curves. However, each of the five exponentials that were fitted to the CO2 decay curves will tend to have some correspondence with one or more of the global carbon cycle processes described earlier (and summarized in Fig. 1) as the individual processes span a wide range of timescales at which they will tend to dominate. The following section discusses the values of the fitting coefficients for each of the five exponentials, giving the median values and range (minimum to maximum) for each parameter. It is likely that the first and second exponentials, with median lifetimes of 6.5 and 50.37 yr (5.3 to 7.4 yr and 39.6 to 57.7 yr) respectively, represent invasion of CO2 into the upper water column (I). The first exponential captures processes with shorter timescales, such as dissolution of gaseous CO2 into the surface ocean, removing 5% (3 to 16%) of the total CO2 emitted. The second represents mixing of dissolved CO2 into the upper water column, resulting in the removal of 10% (7 to 20%) of CO2, meaning that the process(es) of ocean invasion removes 15% (10 to 36%) of the total CO2 emitted to the atmosphere in this model. The third exponential captures the transport of CO2 into the ocean interior (II). It has a timescale of 580 yr (214 to 734 yr), which is comparable to that of ocean mixing, and accounts for the removal of a further 41% (35 to 45%) of the CO2 perturbation. The fourth exponential has an e-folding timescale of 7 kyr (4 to 11 kyr), and represents a combination of seafloor and terrestrial CaCO3 neutralization (III and IV). 28% (13 to 42%) of total CO2 emissions are removed by these processes. Finally, enhanced silicate (and carbonate) weathering (V) results in the uptake of the remaining 7% (7 to 8%) of the excess CO2, operating on an average timescale of 256 kyr (243 to 267 kyr).

Our empirical response function (equation 2) represents a useful tool for rapidly projecting the long-term response of atmospheric CO2 to emissions of any realistic size, removing the need for long simulations using computationally expensive models. However, it should be borne in mind that the first several hundred years of the simulation cannot be considered as realistic, as the emissions are released as a pulse rather than as a time-dependent flux, which would occur in reality. Despite this, atmospheric CO2 concentrations are predicted to within ±11% by year 1000 compared to the equivalent simulation by the Earth system model. The response function is hence suitable for CO2 projections on timescales of a millennium and beyond.

The empirical model can be applied to help further illustrate the evolution of atmospheric CO2 by means of a continuous (interpolated) graphical function of both time and total emissions size, shown in Fig. 3. Readily apparent in this analysis is that while the initial large perturbations subside relatively early on, the remaining atmospheric CO2 anomaly decreases very slowly. For emissions of more than 3000 PgC, even after 1 Myr atmospheric pCO2 not yet equilibrated at natural concentrations (c. 280 ppm). The long atmospheric lifetime of perturbations is in good agreement with the findings of previous studies (Archer, 2005; Colbourn, 2011; Colbourn et al., 2015) and has important implications for future climate evolution, both over the coming centuries but also hundreds of thousands of years from now.

Potential implications for post-closure performance assessments

Glacial-interglacial cycles, largely driven by the waxing and waning of ice sheets in the Northern Hemisphere, have been a dominant feature of the Earth's climate system over the past several million years. The timing and periodicity of Northern Hemisphere glaciation is closely linked to the geometry of the Earth's orbit and its rotation through its impact on insolation (solar radiation received at the top of the atmosphere). Summertime Northern Hemisphere is known to be particularly important (Hays et al., 1976; Kawamura et al., 2007; Milankovitch, 1941; Ruddiman, 2006) as it exerts a strong control on whether winter snow cover persists into the next year and hence accumulates, or completely melts. However, on its own, an insolation-only description of the glacial-interglacial cyclicity is incomplete (Ridgwell et al., 1999) and changes in carbon dioxide (CO2) and hence the degree of greenhouse gas warming also appear intimately linked (Petit et al., 1999). This raises the potential for future anthropogenic activities to impact on the glacial-interglacial cycles by modifying atmospheric CO2 and also, because the current or proposed locations of many nuclear waste sites in the Northern Hemisphere are often at latitudes higher than or close to glacial maximum ice margins, to impact on waste repositories.

Variations in orbital parameters and insolation for the last and next 500 kyr are shown in Fig. 4. In the absence of elevated atmospheric CO2 concentrations, the relatively muted cyclical minima in 65°N June insolation (Fig. 4d) suggests that the next glacial inception may not occur for some tens of thousands of years in the future. For instance, Archer and Ganopolski (2005) found that glacial inception was initiated in the ‘CLIMBER-2’ coupled climateice sheet model when Northern Hemisphere insolation fell below a critical value (Fig. 4d). Under natural conditions of pre-industrial CO2 (~280 ppm) glacial inception occurs periodically with the first occurrence not until after 50 kyr (Fig. 4d; top panel), in good agreement with previous studies (Berger and Loutre, 2002; Paillard, 2001).

The picture changes if CO2 emissions are taken into account. Archer and Ganopolski (2005) used the same model to examine the potential timing of the next glaciation following anthropogenic CO2 emissions of 300, 1000 and 5000 PgC released over 150 years. They found that CO2 emissions had the capacity to affect the future evolution of climate for hundreds of thousands of years, increasing radiative forcing and so decreasing the critical insolation threshold from its natural value. The 300 PgC scenario glaciated after 50 kyr, similar to an experiment forced with natural CO2 evolution, although glacial termination in the emissions scenario took several millennia longer. A release of 1000 PgC prevented glaciation until 130 kyr from present, due to the long atmospheric lifetime of emissions. Anthropogenic emissions of 5000 PgC resulted in a termination of the glacial-interglacial cycles, with the current interglacial continuing for at least the next 500 kyr (Fig. 4d). A different study by Berger and Loutre (2002), forcing a different coupled climate-ice sheet model but with simplified assumptions regarding CO2, also predicted an exceptionally long current interglacial, lasting for a further 5 kyr (for constant pCO2 of 210 ppm) to 50 kyr (with an initial pCO2 of 750 ppm) after present, with the next glacial maximum occurring in 100 kyr. An atmospheric concentration of less than 220 ppm was required to force an early glacial inception.

Although influential and including a representation of coupled climate-ice sheet dynamics not available here, these previous studies did not include a fully explicit description of the long-term decline of fossil fuel CO2. For instance, Berger and Loutre (2002) assumed that atmospheric pCO2 returned to its pre-industrial state following emissions within 1000 years. Archer and Ganopolski (2005) considered different timescales of CO2 decay, but made assumptions regarding the timescale of the silicate weathering feedback. They also fixed the shorter timescales and characteristic removal fractions according to a single magnitude of CO2 release rather than allowing them to become functions of the CO2 release, as is encapsulated in equation 2.

Using our empirical response function, in Fig. 4d we illustrate the potential reduction in the critical insolation threshold due to radiative warming for the 1000, 5000 and 10,000 PgC emissions scenarios, and its decline as a function of time compared to future variations in June 65°N insolation. This was calculated based on the relationship between global radiative forcing from CO2 and the model-predicted critical insolation value found by Archer and Ganopolski (2005). The addition of anthropogenic CO2 to the atmosphere reduces the threshold such that the next glaciation may be delayed until ~130 kyr from now for the 1000 PgC scenario, and beyond 500 kyr for the 5000 and 10,000 PgC scenarios. Positive feedbacks with carbon release from methane hydrates and thawing permafrost, not included in our analysis here but included in the case of hydrates by Archer and Ganopolski (2005), may act to increase the CO2 maximum reached and/or delay the overall rate of decay, hence further displacing future glacial-interglacial cycles. It should also be noted that a significant amount of uncertainty remains in our understanding of the long-term carbon cycle. In particular, past atmospheric CO2 has varied naturally by 80–100 ppm in association with the glacial-interglacial cycles (Petit et al., 1999), but the driving mechanisms for these changes are not yet known (Kohfeld and Ridgwell, 2009).

Our analysis illustrating the importance of CO2 emissions on glacial inception has important implications for waste repositories. During glacial episodes, the advance and retreat of ice sheets can lead to surface erosion (Clayton, 1994), deposition and isostatic effects. Major changes in the hydrological and hydrogeological regime close to the ice margin may also occur (SKB, 2011). Permafrost may develop during glacial or periglacial conditions, reaching depths of tens or hundreds of metres (French, 2007). Hence, a significantly extended interglacial or a delayed glaciation with an extended prior period of periglacial conditions can have important implications, which may be either detrimental or advantageous to safety, for the long-term performance of repositories located in previously glaciated regions.

Other, both natural and anthropogenically-driven environmental changes may also be important. For example, sea-level rise is expected to have a significant impact on the safety of the LLWR (LLWR, 2011). This is affected by the atmospheric lifetime of CO2 emissions through its impact on radiative forcing, which in turn affects thermal expansion of the oceans (Williams et al., 2012) and the risk of melting and/or collapse of the Greenland and West Antarctic ice sheets (Church et al., 2013). In contrast, the proposed Yucca Mountain facility for deep geological disposal of spent fuel and high-level waste, which is located in the unsaturated zone in fractured rock, is thought to be vulnerable to high-intensity fluvial episodes in a warmer wetter climate through their impact on infiltration (Thorne, 2013). This highlights that, while regional changes in temperature and precipitation are important, modifications to atmospheric and oceanic circulation regimes are also of interest due to their influences on local and regional climates.


The comparable decay timescales of radioactivity and atmospheric CO2 emissions highlight the importance of consideration of the latter in post-closure performance assessments for radioactive waste repositories. The characteristic response times of the various processes of the global carbon cycle together with their respective carbon uptake capacities, control the overall timescale and severity of the CO2 perturbation. This can have important implications for future climate changes, including potentially affecting the timing of the future glacial-interglacial cycles, the occurrence (and extent and depth) of permafrost, and the hydrological and hydrogeological regime at different spatial scales. We provide an empirical response function as a simple tool for rapidly projecting the long-term response of atmospheric CO2 to emissions of any realistic size, removing the need for long simulations using computationally expensive models.


The work was funded by RWM via a framework contract with AMEC and Quintessa. It has benefited from discussions with and material produced by Working Group 6 of the IAEA-sponsored MODARIA Programme.

The publication of this research has been funded by the European Union's European Atomic Energy Community's (Euratom) Seventh Framework programme FP7 (2007–2013) under grant agreements n°249396, SecIGD, and n°323260, SecIGD2.
P1Freely available online through the publisher-supported open access option.