Atmospheric CO2: A 2-Box model accurately tracks 14C and 13C spanning 200 years and estimates global uptake of fossil fuel emissions

1. Introduction

Unlike nett values of CO₂ exchange flux between the atmosphere and the land/ocean, the absolute each-way values are not accurately known, yet they are key to understanding the carbon cycle and to estimating the effects of anthropogenic CO₂ fossil fuel emissions, (CO2ff) (IPCC (Ed.) 2013; IPCC(Ed) 2021). The exchange flow is caused by diurnal and seasonal variations which, as the earth spins on its axis and orbits the sun, alternately cause absorption and emission, arising from cyclical diurnal and seasonal changes in seawater solubility, photosynthesis, plant decay and respiration. This complicated pattern of both positive and negative fluxes implies that while absorption is occurring at one geographical location, there is simultaneous emission at another. This two-way exchange flow is estimated at approximately 210GtC in each direction per year, easily exceeding the one way CO2ff despite its rise from 2.5 GtC in 1960 to 9 GtC in 2021 (Friedlingstein et al., 2022). However, although CO2ff is smaller than exchange flow, it is still significant when compared to the imbalance in the exchange flux. In the 1970s box models were developed (e.g., Oeschger et al., 1975) to describe the carbon cycle; however, their popularity waned as General Circulation Models (GCMs) became common which, although requiring greater computer resources, could create gridded global maps of the carbon-related processes.

Undoubtedly GCMs have a significant role to play in modeling the complexity of the global carbon cycle. However, there are advantages in adopting a strategy of top-down simplification using box models, rather than adopting the open-ended system of micro-accounting in the GCMs. This is because simpler models by their very nature are "better testable" (Popper, 1934) and are more likely to shed light on errors in formulation, focusing on the most significant factors. This paper describes a 2-box model that derives ¹⁴C and ¹³C isotopic ratios given anthropogenic emissions and CO₂ mixing ratio. The major differences between this and previous box models are; a) its use of atmospheric CO₂ level as an input, not an output, and b) its use of absolute gross flow, not nett flow. This model accurately calculates isotopic ¹⁴CO2 and ¹³CO2 time-series spanning some 200 years (including the bomb pulse). Yet according to conventional wisdom the model should not work. It is widely known that the carbonate ions form a buffer solution limiting the absorption rate of CO₂, this effect being known as the Revelle factor. It is commonly held that the ¹³C and ¹⁴C isotopic forms of carbon bypass seawater carbonate chemistry, exhibiting very different absorption properties as compared to ¹²C. (Tans 1993, Joos 1994, Harvey 2000), although this was never proven. This study demonstrates this assumption is invalid, since it uses the same residence time and reservoir mixing properties for all isotopes and yet produces accurate results. For a full discussion see Section 7.1.

1.1. Radiocarbon

The global isotopic abundance of ¹³C and ¹⁴C present within atmospheric CO₂ is accurately known, having been measured in samples from tree rings (Stuiver & Quay 1981), ice cores (Rubino et al., 2013), and the atmosphere (Stuiver et al., 1998; D6). While ¹³C is radioactively stable and forms around one percent of atmospheric carbon, ¹⁴C undergoes radioactive decay with a half-life of 5700 ± 30 years (Kutschera 2013), comprising about 1 part in 10¹². There is virtually no presence of ¹⁴C in fossil fuel since it has all radioactively decayed. The ¹⁴C/¹²C ratio decreased between 1800 and 1960 because of the lack of ¹⁴C in fossil fuel emissions. Then in 1960 atmospheric nuclear weapons testing caused a sudden doubling of the absolute level, also known as ¹⁴C activity. Since 1975 it has decreased (D6), partially due to washout from the exchange flow which contains a somewhat lower level of ¹⁴C but also because of dilution arising from fossil fuels CO₂ input, which is ¹⁴C depleted (Levin, 2010). This process is known as "Suess" dilution (Suess 1955). Turning now to ¹³C, its presence in fossil fuels is reduced when compared with the atmosphere, but only by around 2% (Stuiver and Polach, 1977); this small reduction arising from fractionation, i.e., because the larger less mobile isotopic CO₂ molecules are less likely to become embedded in leaves and fossil carbon. Measurements show that the atmospheric ¹³C/¹²C ratio has decreased during the past 200-year period (D4), due to dilution by fossil fuel carbon which contains the lower level of ¹³C and washout from the exchange flow.(see Discussion)

2. CO2 Absolute Flow Finite Reservoir Model

The model describes the movement of CO₂ and its carbon isotopes between atmosphere and terrestrial / oceanic mixing reservoir. See Figure 1. The size of this mixing reservoir is an "effective" value, corresponding to the combined effect of the ocean and terrestrial storage. The idea of effective values is commonplace in many fields of engineering, science and economics. Economic inflation represents, in a single figure, the average of many different price increases on many items. Productivity measures the effective value created by each hour of work. Effective resistance in electronic engineering or hydraulics gives the overall resistance of a complex system to current or fluid flow. The assumption in each of these examples is that a single value can account for some property which is distributed throughout a complex system. The elegance of this method is its ability, unlike other methods, to determine the effective values without requiring a detailed audit of the underlying constituent components.

Figure 1. CO₂ Finite Reservoir Two-Box Model. The model describes the CO₂ flux between the reservoir and atmosphere, including its isotopic carbon contents. Isotopic mixing occurs in the two storage regions, atmosphere and reservoir. A proportion, β of anthropogenic fossil fuel CO₂, F_a, enters the atmosphere, while the remainder, 1-β enters the reservoir by direct uptake. At each iteration, the CO₂ flux outflow F_e is proportional to atmospheric CO₂, C(t); its isotopic content corresponding to that of the atmosphere. The CO₂ flux in, F_i is determined from F_e and the rate of change ΔC/ΔT, with its isotopic content corresponding to that of the reservoir, R. Atomic weapons ¹⁴C B(t) directly enters the atmosphere at intervals and amounts determined by historical records of atmospheric atomic weapons testing and bomb yield.

Figure 1. CO₂ Finite Reservoir Two-Box Model. The model describes the CO₂ flux between the reservoir and atmosphere, including its isotopic carbon contents. Isotopic mixing occurs in the two storage regions, atmosphere and reservoir. A proportion, β of anthropogenic fossil fuel CO₂, F_a, enters the atmosphere, while the remainder, 1-β enters the reservoir by direct uptake. At each iteration, the CO₂ flux outflow F_e is proportional to atmospheric CO₂, C(t); its isotopic content corresponding to that of the atmosphere. The CO₂ flux in, F_i is determined from F_e and the rate of change ΔC/ΔT, with its isotopic content corresponding to that of the reservoir, R. Atomic weapons ¹⁴C B(t) directly enters the atmosphere at intervals and amounts determined by historical records of atmospheric atomic weapons testing and bomb yield.

The two boxes in the model are considered well-mixed. At each iteration, the new isotopic concentration within the reservoir and atmosphere caused by the incoming ¹³C and ¹⁴C is calculated. We now explain the ideas behind the model's formulation. Using a similar manner and notation as Tans et al.,(1993) and Cawley (2011) we have

dC/dT = F_a + F_i - F_e

where C is the atmospheric carbon mass, F_a is the anthropogenic flux, F_i, is the flux in and F_e is the flux exiting the atmosphere. Many authors e.g., Oeschger (1975), assume a pre-industrial equilibrium but their method is not suitable for our purposes since it relies on the calculation of a net flux; it is essential that F_i and F_e are kept separate because we will need the separate fluxes to calculate the isotopic flows. We first need an expression which connects the atmospheric CO₂ mass C, to the outflow, F_e an obvious choice being

F_e = k_e C

(1)

Absolute inflow, as opposed to outflow, can then be calculated by balancing the budget; if the atmosphere is growing by an amount ΔC in time ΔT it must be receiving the outflow plus the growth amount minus the inflow from anthropogenic input, which is anthropogenic output, F_a, times the proportion which arrives in the atmosphere, β, hence we obtain

F_i = F_e + ΔC / ΔT - βF_a

(2)

Equation (1) satisfies the obvious requirement that the exit flux F_e , is zero when CO₂ mass C is zero. A further requirement, satisfied by Equation (2) is that when F_i - F_e is zero the atmospheric growth ΔC/ΔT, equals the anthropogenic input, F_a. These requirements distinguish this absolute flow description from other descriptions such as provided by Oeschger (1975) and Cawley (2011). This formulation leads to a radically different solution from conventional box models. It avoids the necessity to include allowance for seawater carbonate chemistry, because the outflow behaviour is lumped into a single constant, k_e which is empirically derived. Similarly the inflow behaviour, being calculated by balancing the budget, eliminates the need to provide a bottom-up audit of CO₂ release from the ocean or land due to outgassing, respiration, decay and combustion. Once outflow is derived, inflow can be defined from the known anthropogenic input and CO₂ measurements. We assume the constant k_e applies equally to the calculation of outflow for all three carbon species ¹⁴CO₂, ¹³CO₂ and ¹²CO₂ (see Discussion).

The model estimates the nuclear bomb yield of ¹⁴C each year, B₁₄[i], by using a bomb yield multiplication factor, Yb, this being one of the seven optimisation parameters. The bomb yield in megatons was obtained from records of atomic weapons atmospheric test detonations [D5] and was converted using Yb. A simple time delay of one year was used to allow for atmospheric mixing. The initial values i.e., δ¹³C_init, Δ¹⁴C_init, are parameters within the model and determine the initial level of the curves in Figures 3 and 4, while δ¹³C_FF determine the curve slopes in Figures 4 and 5.

2.1. Mass-Balance and Isotopic Dilution

If a number of gases are mixed but do not react, and each component gas, i, has a mass M_i, then from the law of conservation of mass (also known as Mass Balance)

M_T = Σ M_i

If within each component gas there is a specific molecule present in a ratio to its mass, R_i, then, since that specific molecule is conserved, the resulting mixture (e.g., paint mixing) having mass M_T has a ratio R_T given by

R_T = Σ R_i M_i / M_T

(3)

Equation (3) applies to the situation of mixing of isotopic gases. For practical reasons associated with historical measurements of radioactivity, studies of isotopic mixtures normally utilize the ratio of the sample activity, A_S to a standard sample activity A_abs, this ratio being known as the relative specific activity, A. Its value is then offset by one to provide the radioactivity scale, δ written (see Strenstrom 2011) in the form

δ = (A_s/A_abs) -1 = A - 1

(4)

In this scale, although δ increases linearly with A_s, its value is offset to be zero at the value of the absolute standard (often chosen to approximate background level), as for example in the cases of δ13Cand δ14C. Applying (4) and substituting A_Si for each component R_i of the mixture, leads after some algebraic manipulation, to

δ _T = Σ (δ _i M_i ) / M_T

(5)

The similarity with Equation (3) shows that although δ is defined by an offset ratio scale we can still directly apply Equation (3) to calculate the resulting value of δ _T for the mixture. The above scheme forms the basis for the system of equations A₁₄[], R₁₄[], A_FF[] and R_FF[] shown in the implementation section below.

2.2. Fractionation

The model generates values of δ ¹⁴C but comparison is required with measured values of Δ¹⁴C, which incorporates a fractionation correction to "translate the measured activity to the activity the sample would have had if it had been wood" (Strenstrom 2011) . The correction formulae, a function of δ ¹³C, may be derived from equations by Strenstrom 2011 using Equations (25, 28, 38, or Stuiver & Polach 1977 p356, 360, Table 1 , where A_SN is the normalised specific sample activity, A_s is the specific activity of the atmosphere or reservoir, A_abs is the absolute specific standard, δ ¹³C_M is the ¹³C sample measurement , and δ ¹³C_W is the ¹³C standard for wood, giving

δ ¹⁴C = A_S / A_abs - 1

Δ¹⁴C = A_SN / A_abs - 1

A_SN = A_S [ ( 1 + δ¹³C_W) / ( 1 + δ¹³C_M) ]²

Eliminating A_S , A_abs and A_SN gives

Δ¹⁴C = [1 + δ ¹⁴C ] [ (1 + δ¹³C_W) / (1 + δ ¹³C_M) ]² - 1

The fractionation correction was applied to δ ¹⁴C, with δ ¹³C for wood as -25‰ (the standard), and δ ¹³C_M being taken from the model. The correction is quite small, Δ¹⁴C is around 10‰ less than δ ¹⁴C at the peak of the bomb pulse ( δ ¹⁴C ~ 800‰) decreasing to 0.03‰ (δ ¹⁴C ~0).

2.3. Age Correction

The ¹⁴C half-life of 5700± 30 years (Kutschera, W., 2013), which translates to a radioactive decay of around 2% over the period from 1820 to 2020. However, we know the steady production of stratospheric ¹⁴C roughly balances the amount of ¹⁴C decay since the two are in approximate equilibrium. Therefore neither ¹⁴C decay nor natural stratospheric ¹⁴C production were included in the model.

2.4. Implementation

The implementation processes several calculations at each iteration. The following 3 Equations (are based upon . The amount of CO₂ in the atmosphere is C, its removal time is α_t, the atmosphere outflow is F_e , the inflow is F_i , and the amount of CO₂ initially in the reservoir is R, with the square brackets '[i]' referring to the value at iteration i, hence the iteration relationships (using Equations (1 and 2) become:-

F_e[i] = C[i] /α_t

F_i[i] = F_e [i] + C [i] - C [i-1] - β F_a [i-1]

R[i] = R [i-1]+ F_e[i] - F_i[i] + (1- β) F_a [i-1]

The ratio of ¹⁴C / C in the atmosphere is denoted by A₁₄ and for the reservoir by R₁₄. The new amount is given by the previous amount added to the input amounts with the output amount being subtracted, all divided by the total mass, see equation 5. If necessary we include the annual ¹⁴C production due to atomic weapons bomb input, B₁₄. Hence

A₁₄[i] = (A₁₄ [i-1] . C [i-1] + F_i [i-1].R₁₄[i-1] - F_e[i-1]. A₁₄[i-1] + B₁₄[i-1]) / C [i]

R₁₄[i] = (R₁₄ [i-1] . R [i-1] - F_i [i-1].R₁₄[i-1] + F_e[i-1]. A₁₄[i-1] ) / R [i]

Similarly, the relative atmospheric and reservoir fossil fuel content A_FF and R_FF, are calculated on a scale of 0 to 1 using the values at the previous iteration levels, where A refers to the atmosphere, and R the reservoir and subscripts FL means fossil fuel level, and NL means natural non-fossil level,

A_FF[i] = (A_FF [i-1].C[i-1] + F_i[i-1].R_FF [i-1] - F_e[i-1].A_FF[i-1] + β F_a [i-1]) / C [i]

R_FF[i] = (R_FF [i-1].R [i-1] - F_i[i-1].R_FF [i-1] + F_e[i-1].A_FF[i-1] + (1- β) F_a [i-1]) / R [i]

A_FL[i] = C[i] A_FF[i], A_NL[i] = C[i] - A_FL[i]

R_FL[i] = R_FF[i] R[i], R_NL[i] = R[i] - R_FL[i]

The offset isotopic ratiometric measure for ¹³C is δ¹³C and that for ¹⁴C is Δ¹⁴C. These are calculated from the iterative series by:-

δ ¹³C = δ ¹³C_FF . A_FF + δ ¹³C_N ( 1- A_FF ), Δ¹⁴C = (A₁₄ -1)

Initial Conditions: F_i [0] = F_e[0] , A₁₄ [0] = 1 , R₁₄ [0] = 1, A_FF[0] = 0 , R_FF[0] = 0

In accounting for isotopic concentration, it is not necessary to explicitly embed an attenuation factor as discussed by Stuiver & Quay (1981), or to account separately for Suess dilution (1955) because they are implicitly represented.

3. Method

The input data was selected between 1750 to 2020 from the appropriate data source listing [D1-D6] and entered into a spreadsheet. The cell formulae were set to correspond to those iteration Equations (given in the Appendix. The solution was found by minimizing the discrepancy between observed and calculated values of δ¹³C and Δ¹⁴C using the standard "Solver" function of Microsoft Excel. The solver system was set to use the 7 parameters shown in Table 1 using the default solver options. For version details see appendix.

Figure 2. Solution Optimisation. The 7 parameter model and the deviation comparison formulae are implemented as cell formulae within Microsoft Excel. The deviation comparison formulae compares the measured time-series values of atmospheric Δ14C and δ13C with the calculated values. The Add-In Solver was used to numerically minimise the standard deviation by varying the values of the 7 parameters, providing the above outputs shown below in this study.

Figure 2. Solution Optimisation. The 7 parameter model and the deviation comparison formulae are implemented as cell formulae within Microsoft Excel. The deviation comparison formulae compares the measured time-series values of atmospheric Δ14C and δ13C with the calculated values. The Add-In Solver was used to numerically minimise the standard deviation by varying the values of the 7 parameters, providing the above outputs shown below in this study.

The optimisation-solution process is shown in Figure 2 above. The standard deviations σ₁ and σ₂ between observed and predicted was determined for δ¹³C and Δ¹⁴C. The square root of this product gives the overall geometric mean standard deviation, σ. Hence

σ_j = √1/n  (obs_i - calc_i)²

(3)

σ_ = (σ₁ x σ₂ )^1/2

(4)

The Excel Solver was set to minimise the value of the target cell containing the overall geometric mean standard deviation by variation of the seven parameters shown in Table 1. Only one true constant was used, atmospheric capacity, which was taken to have a value of 2.124 ppm GtC-1 as published by Ballantyne (2012). Thus the model uses a very different technique from standard modeling techniques for deriving anthropogenic and non-anthropogenic flux. It optimises the seven parameters to fit δ¹³C and Δ¹⁴C ratios, and the flux and mixing reservoir levels emerge as by-products from this process.

4. Results

The results divide into three main types; the optimised seven parameters themselves, the graphical output, and flux tables. Flux and cumulative flux are model outputs.

Figures 3 and 4 show actual and modelled values of ∆¹⁴C. In total 340 observed data points were used giving a total overall standard deviation of 0.39‰. The quality of fit is excellent, with a standard deviation of 3‰. In Figure 3 the slight decrease from 1820 to 1940 is due to Suess dilution, the magnitude of the dilution being somewhat reduced by reservoir inflow, giving excellent agreement with observations. Figure 4 shows the bomb pulse plotted as ∆¹⁴C. Initially the pulse shape is a decaying exponential but then, as it is influenced by Suess dilution and the re-emergence of ¹⁴C from the mixing reservoir, becomes more of a linear fall. Nevertheless, despite these complicating factors, the predicted and observed values are in excellent agreement. Figures 3 and 4 also show two hypothetical scenarios "no bombs" and "no fossil". The "no bombs" scenario shows how the decrease would have continued without nuclear atmospheric weapons testing, while "no fossils" show the situation without fossil fuel emissions.

Figure 3. Atmospheric Δ14C 1820-2020. showing 130 values before 1950 [D6] and 70 values from 1968 onwards (Hua et al.. 2022). Standard deviation of observed and predicted is σ = 3.0‰, excluding values during the rising pulse from 1950 to 1968. For "no bombs" and "no fossil". See text.

Figure 3. Atmospheric Δ14C 1820-2020. showing 130 values before 1950 [D6] and 70 values from 1968 onwards (Hua et al.. 2022). Standard deviation of observed and predicted is σ = 3.0‰, excluding values during the rising pulse from 1950 to 1968. For "no bombs" and "no fossil". See text.

Figure 4. Atmospheric Δ14C between 1950 and 2020 showing the "bomb pulse" to 2020. For "no bombs" and "no fossil" see text.

Figure 4. Atmospheric Δ14C between 1950 and 2020 showing the "bomb pulse" to 2020. For "no bombs" and "no fossil" see text.

Figure 5 shows the change in δ¹³C, with a standard deviation of 0.05‰. with the reduction being caused by Suess dilution (since fossil fuel emissions contains a lower level of ¹³C) and by re-emission of CO₂ from the reservoir. This shape is also accurately tracked by the model.

Figure 5. Atmospheric δ13C 1820-2020. Model calculated values: solid green, Observed NOAA Paleo [D4]: squares, Mauna Loa [D2]: crosses, σ_1820-2020= 0.05‰.

Figure 5. Atmospheric δ13C 1820-2020. Model calculated values: solid green, Observed NOAA Paleo [D4]: squares, Mauna Loa [D2]: crosses, σ_1820-2020= 0.05‰.

Table 1. Solved 7 Parameter Values.

Parameter	Sym.	Value	err.	Notes
Removal time (years)	αt	14.87	1.7	Compares reasonably with Revelle & Suess (1957) of approximately 10 yr, and Arnold (1957) of 10-20 yr.
Relative Reservoir Size	RCO2	6.16	1.4	IPCC Ed., (2014) of 4250 GtC. 4250/588 = 7.2.
CO2ff Direct Uptake	β	0.56	0.11	See Discussion
Atomic Bomb Yield: 14C/MT	Yb	1.60	0.1	In mills. See Footnote
Pre-industrial14C ‰	Δ14Cinit	-3.1	10	Broad agreement with Intcal20[D6] for 1820 of 0.7‰. Max level between 1820-2020 is 800‰.
Pre-industrial13C ‰	δ13Cinit	-6.7	0.2	Excellent agreement with Rubino (2013) 1820 of δ13C =6.7‰.
Fossil fuel13C ‰	δ13Cff	-20.4	4	Slightly low as compared to figures reported by Stuiver & Polach (1977) for coal of δ13C = -23 ‰.

Table 1. Solved 7 Parameter Values.

Parameter	Sym.	Value	err.	Notes
Removal time (years)	αt	14.87	1.7	Compares reasonably with Revelle & Suess (1957) of approximately 10 yr, and Arnold (1957) of 10-20 yr.
Relative Reservoir Size	RCO2	6.16	1.4	IPCC Ed., (2014) of 4250 GtC. 4250/588 = 7.2.
CO2ff Direct Uptake	β	0.56	0.11	See Discussion
Atomic Bomb Yield: 14C/MT	Yb	1.60	0.1	In mills. See Footnote
Pre-industrial14C ‰	Δ14Cinit	-3.1	10	Broad agreement with Intcal20[D6] for 1820 of 0.7‰. Max level between 1820-2020 is 800‰.
Pre-industrial13C ‰	δ13Cinit	-6.7	0.2	Excellent agreement with Rubino (2013) 1820 of δ13C =6.7‰.
Fossil fuel13C ‰	δ13Cff	-20.4	4	Slightly low as compared to figures reported by Stuiver & Polach (1977) for coal of δ13C = -23 ‰.

Nuclear Bomb Yield: The value returned is 1.6 ± 0.1 mills of 14C per MT, with a total bomb yield of 440MT. For comparison, we convert the figures to estimates of the number of 14C atoms per MT bomb yield using formulae derived from those of Svetlik (2010) and Strensom (2011). The number of carbon atoms in the atmosphere in 1950 with the CO₂ mixing level of 312.8 ppm corresponding to ACO2 = 664.43 GtC, where the molecular weight of naturally abundant carbon, MC = 12.01, and Avogadro's number N_A = 6.022×10²³, was

NC = A_CO2 x N_A / MC

The ¹⁴C/C ratio R¹⁴C can be obtained from the specific activity α = 0.226 Bq per gram, ¹⁴C half life T_1/2 = 5730 years = 1.807 10¹¹ s, molar weight of ¹⁴C is M_C14 = 14, and ln(2) as

R¹⁴C = α T_1/2 M_C14 / (N_A x ln(2))

The number of ¹⁴C atoms N¹⁴C created is given by the product of the two expressions giving

N¹⁴C = NC x R¹⁴C = (A_CO2 / MC) x ( α T_1/2 MC₁₄ / ln(2)) = 4.6 x 10²⁸.

Hence the number of 14C atoms per MT of bomb yield is 0.0016 x 4.6 x 10²⁸ = 7.4 x 10²⁵, comparing favourably with reported values by Hesshaimer et al. (1994) of 1 to 2 x 10²⁶ atoms per MT. The total bomb yield of 440MT gives a figure of 3.23 x 10²⁸ which compares with 6.1 x 10²⁸ by Naegler & Levin (2006). The discrepancy may be due to ¹⁴C which is stored but does not mix in the biosphere.

The amount of ¹⁴C, also known as activity concentration, (as opposed to Δ¹⁴C), in the reservoir and atmosphere is shown in Figure 6. Their total (blue) represents the amount of ¹⁴C in the system, the sudden step increase can be seen around 1960 due to atmospheric nuclear weapons testing. The atmospheric activity concentration shows the bomb pulse decaying to a minimum around 2000 but then rising slightly again, which has been experimentally reported by Svetlik (2010). The model shows the increase can be attributed to ¹⁴CO₂ re-emerging from the mixing reservoir, at the same time causing a slight drop in reservoir level.

Figure 6. ¹⁴CO2 Budget. Model derived absolute ¹⁴CO2 (also known as activity concentration) for the reservoir (orange), the atmosphere (green), and their total (blue) in units relative to the ¹⁴CO2 content of the atmosphere in 1950. The total, being the sum of reservoir plus atmosphere, shows the near step function in ¹⁴CO2 created by nuclear atmospheric weapons testing. After the initial fall of atmospheric content, the small rise since 2000 is confirmed by the results of Svetlik [23] (see text).

Figure 6. ¹⁴CO2 Budget. Model derived absolute ¹⁴CO2 (also known as activity concentration) for the reservoir (orange), the atmosphere (green), and their total (blue) in units relative to the ¹⁴CO2 content of the atmosphere in 1950. The total, being the sum of reservoir plus atmosphere, shows the near step function in ¹⁴CO2 created by nuclear atmospheric weapons testing. After the initial fall of atmospheric content, the small rise since 2000 is confirmed by the results of Svetlik [23] (see text).

4.1. Parameter Values

The solution of all seven parameter values, is listed in Table 1 along with their estimated error, which was derived by varying each parameter up and down, while holding the other parameter values constant, until the overall geometric mean standard deviation, σ, approximately doubled. The error figure was the mean of the two variations. For all the parameters except fossil-fuel direct uptake, the variation by the error figure causes an insignificant change in the flux figures, however in the case of fossil-fuel direct uptake the change in non-fossil input is from -10GtC to 83GtC, implying an uncertainty in non-fossil input from -3.5% sink to 29% source.

Figure 7. Atmospheric CO₂ levels and cumulative flux 1750 to 2020. (GtC) The figures in red represent the increase over the period, underscore indicates 2020 values. The concentration of CO₂ from fossil fuel emissions (CO2ff) is higher in the atmosphere than in the reservoir, therefore the atmospheric flux OUT contains more CO2ff than the flux IN, removing fossil fuel emissions from the atmosphere. Atmospheric CO2ff = 253 - 295 + 164 = 121 (without rounding) corresponding to 13.9% of the atmosphere in 2020 (See Table 2, rows 3 / 9).

Figure 7. Atmospheric CO₂ levels and cumulative flux 1750 to 2020. (GtC) The figures in red represent the increase over the period, underscore indicates 2020 values. The concentration of CO₂ from fossil fuel emissions (CO2ff) is higher in the atmosphere than in the reservoir, therefore the atmospheric flux OUT contains more CO2ff than the flux IN, removing fossil fuel emissions from the atmosphere. Atmospheric CO2ff = 253 - 295 + 164 = 121 (without rounding) corresponding to 13.9% of the atmosphere in 2020 (See Table 2, rows 3 / 9).

The study reports that, over the past 270 years, 27% of CO2ff remains in the atmosphere corresponding to 13.9% of the atmospheric CO₂ content, and 88% of atmospheric growth is caused by CO2ff (see Table 2). It is found that 44% of CO2ff is directly absorbed into terrestrial biosphere and ocean with the remainder i.e., 56%, becoming mixed into the atmosphere. A small amount, around 12% of atmospheric growth originates from natural sources although this figure had the largest uncertainty. See Discussion.

4.2. Flux and Cumulative Flux Output

Table 2 lists cumulative CO² flows and storage from 1750 to 2020, and 1960 to 2020, while fig. 7 shows the flux diagrammatically.

Table 2. Cumulative CO2 Flow and Storage.

Duration	1750-2020		1960 -2020
	GtC	%	GtC	%
CO2ff Delivered
delivered to Atmos. (β)	253	56	207	46
delivered to Rsvr. (1-β)	201	44	165	36
Total	454	100	372	82
CO2ff by Destination
In Atmos	121	27	100	26
In Rsvr	332	73	272	71
Total	454	100	382	100
Atmospheric Growth
due to CO2ff	253	88	207	102
due to non-Foss*	35	12	-3	-1
Total	287	100	204	100
Reservoir Growth
due to CO2ff	332		272
due to non-Foss*	-166		-103
Total	166		164
Atmospheric Outflow
CO2ff	295		248
non-Foss	11447		2798
Total	11742		3045
Reservoir Outflow
CO2ff	164		141
non-Foss	11614		2902
Total	11777		3042
CO2ff rel. to atmos CO2 2020†	121/876	14	100/876	11.5
CO2ff induced growth/Co2ff (AF)	253/454	55	56
Average Annual Flux from Reservoir	43.7		50.8
Average Annual Flux to Reservoir	43.5		50.9

+ Relative to CO2ff supplied; * Non-Foss means "not due to anthropogenic fossil fuel CO2 emissions". † The CO2 atmospheric mixing ratio in 2020 = 412 ppm, 876 GtC, in1750 = 277ppm, 589GtC. These values take into account the backflow of fossil fuel CO2 from the reservoir to the atmosphere. Year Intervals are chosen to align with GCB from Friedlingstein (2021) Table 8. The flux / cumulative flux tables did not form part of the optimisation loop. Note that some of the figures are not model dependant, (e.g., CO2ff Total, Atmospheric Growth Total) but are based upon subtractions of input data. See [3]. Table 8.

Table 2. Cumulative CO2 Flow and Storage.

Duration	1750-2020		1960 -2020
	GtC	%	GtC	%
CO2ff Delivered
delivered to Atmos. (β)	253	56	207	46
delivered to Rsvr. (1-β)	201	44	165	36
Total	454	100	372	82
CO2ff by Destination
In Atmos	121	27	100	26
In Rsvr	332	73	272	71
Total	454	100	382	100
Atmospheric Growth
due to CO2ff	253	88	207	102
due to non-Foss*	35	12	-3	-1
Total	287	100	204	100
Reservoir Growth
due to CO2ff	332		272
due to non-Foss*	-166		-103
Total	166		164
Atmospheric Outflow
CO2ff	295		248
non-Foss	11447		2798
Total	11742		3045
Reservoir Outflow
CO2ff	164		141
non-Foss	11614		2902
Total	11777		3042
CO2ff rel. to atmos CO2 2020†	121/876	14	100/876	11.5
CO2ff induced growth/Co2ff (AF)	253/454	55	56
Average Annual Flux from Reservoir	43.7		50.8
Average Annual Flux to Reservoir	43.5		50.9

+ Relative to CO2ff supplied; * Non-Foss means "not due to anthropogenic fossil fuel CO2 emissions". † The CO2 atmospheric mixing ratio in 2020 = 412 ppm, 876 GtC, in1750 = 277ppm, 589GtC. These values take into account the backflow of fossil fuel CO2 from the reservoir to the atmosphere. Year Intervals are chosen to align with GCB from Friedlingstein (2021) Table 8. The flux / cumulative flux tables did not form part of the optimisation loop. Note that some of the figures are not model dependant, (e.g., CO2ff Total, Atmospheric Growth Total) but are based upon subtractions of input data. See [3]. Table 8.

4.3. Correlation of Inflow and Temperature

Figure 8 shows the absolute inflow (taken from the absolute model using the same 7 optimised parameters) and global temperature (D7) plotted on suitable vertical scales. The correlation coefficient is 0.94. Although correlation does not by itself imply causation, it is difficult to conceive of a mechanism whereby a global change in temperature could be caused by the rate of CO₂ inflow. The obvious temperature dependents mechanisms for CO₂ gas release in the ocean is the decrease of solubility with temperature, while on land there is a balance between productivity, respiration and decay (Melillo et al. 2011). This evidence suggests that when global temperature increases, CO₂ absolute inflow also increases although it could be conceivable that both share a common cause. The absolute global inflow should not be confused with net figures which show a net uptake of anthropogenic CO₂.

Figure 8. CO₂ inflow (green) and observed global sea-surface temperature (blue) indicating a high degree of correlation of 0.94.

Figure 8. CO₂ inflow (green) and observed global sea-surface temperature (blue) indicating a high degree of correlation of 0.94.

5. Discussion

This section discusses the differences between this work and that of other authors. Some authors (e.g., Levin et al. 2010) have found it necessary to propose additional industrial pollution figures for ¹⁴C in order to close the ¹⁴C budget. This paper has no such requirement, yet accurately reproduces the ¹⁴C curve. A recently published work (Graven et al. 2020) used a model originally produced by M. Heimann and R. Keeling and rewritten for Matlab, to give results for a "no bombs" and "no fossil" scenarios , in a similar way to those shown in Figure 5 of this paper. The results shown by Graven are almost identical to those shown here further vindicating the model described and its approach. Recent work by Svetlik (2010) showed measurements of the ¹⁴C activity concentration. This model produces an output of activity concentration which is very similar to that in Svetlik's work showing that the calculations regarding ¹⁴C are in broad agreement.

5.1. Isotopic Behaviour

It is a widely held view that the various isotopes of carbon in their molecular form of CO₂ behave differently when the exchange flux interacts with sea-water, bypassing carbonate chemistry (Revelle). This was first suggested by Siegenthaler and Oeschger (1987) and has been repeated many time since. It is claimed to lead to a factor of ten difference in isotopic behaviour of ¹³C and ¹⁴C as compared to ¹²C, thus preventing a straightforward analysis of the ¹⁴C bomb pulse (e.g., Harvey 2000). The argument is refuted here. It suggests isotope separation occurs naturally on a global scale. This is unlikely since, from a thermodynamic point of view, the isotopic mixture has a lower Gibbs free energy and lower entropy, requiring an expenditure of energy to achieve separation of its constituents (Seader & Henley 2005). We present here a further argument showing why we believe the claimed bypass effect does not occur in practice. Given the carbon concentration C and the isotopic ratio R the flux F14C can be written as

F14C = d(CR)/dt = R.dC/dt + C.dR/dt

The term, R.dC/dt, describes the bulk flow of carbon without isotopic variation and hence suffers from carbonate chemistry while C.dR/dt represents the isotopic flow which it is claimed bypasses seawater carbonate chemistry. We can better estimate the magnitude of the contribution of variations of either term, by dividing both sides by CR giving

A B

d(CR)/CR = dC/C + dR/R

A comparison of the local values of the relative sizes of dC/C and dR/R indicates which term dominates. Fortunately one can obtain approximate estimates of these ratios in a variety of situations, over sea and land from a number of papers (e.g., Bishcof 1960, Dai et al. 2009, Faassen et al., 2022, Leinweber et al., 2009, Olsen et al. 2004, Palonen et al. 2018, Takahashi et al. 2002). It is essential to use local estimates rather than global nett figures, because nett figures would average out differences in isotopic content. The local variation of dC/C was estimated from the above references as of the order of ≈6 ppm/400ppm (≈0.015) per hour, with variation in dR/R for 14C estimated as around ≈8‰ (≈0.008) per hour and for 13C ≈ 0.2 ‰ (≈0.0002) per hour. Hence |A| amply dominates |B|. Therefore, it is irrelevant whether or not B has an ability to bypass carbonate chemistry (Revelle factor) as claimed by many, since A is the dominant effect. Regarding the A term how does it respond to different isotopes of carbon? According to (Zeebe and Wolf 2001: Section 2.5 Page 119) referring to three isotopic forms ¹²C, ¹³C and ¹⁴C isotopic components “There are no direct reactions between each of these compartments... Thus there are three independent systems with respect to carbon isotopes.” Furthermore there is only one common subsidiary chemical process, i.e., “chemical coupling of both systems brought about by the chemical reactions of the carbon compounds with "H+ and OH-". This is assumed here to be a relatively weak effect. Therefore in this study we consider the isotopic behaviour to be equal, although we accept there could be a small differences due to fractionation. The accuracy of the results for Δ14C and δ13C shown in Figure 3, Figure 4 and Figure 5 results amply vindicates these assumptions, despite the widely held belief to the contrary. Furthermore, the values of the 7 parameters in Table 1 all have reasonable values, which seem extremely unlikely by chance. Therefore, regarding global absorption of atmospheric isotopic CO₂ by seawater, the argument is rejected that there is significant bypass of the Revelle factor.

5.2. Exchange Flow

As the earth daily rotates and annually orbits the sun, the terrestrial /oceanic reservoirs cyclically sucks in and blows out CO₂. These are the exchange fluxes. In the model, the effective exchange flow is found to be around 50 GtC, considerably smaller than the IPCC figure of exchange flow at around 210 GtC. Consequently the exchange washout Residence Time (RT) in 1960 was around C/Fe= 689/46 ~ 15 years. Thus the residence time RT now aligns with the exponential decay of the atomic weapons pulse, without requiring the usual, but unlikely, hypothesis of unequal behaviour of isotopes. The difference in exchange flow arises because the model determines "effective values" which reflects the global average behaviour of molecular transfer, bearing in mind the mixing properties and cycle time lag (Oeschger 1975). Consider an analogy. Two isolated islands are separated by the services of a large passenger ferry. You may think the average residence time of a person on each island is given by the number of island inhabitants divided by the ferry flow. But the true average can be much less than this. How? Because the passenger ferry may carry much the same passengers back and forth, while many of the rural residents never take the ferry. Therefore the residence time is much longer because the "effective" washout flux is smaller. Similarly if the same exchange flow arises from decaying leaves and plant litter being cycled each year, the effective exchange flow is reduced, changing the residence time from around 4 years to 15.

5.3. Airborne Fraction

The Airborne Fraction, (AF) is said to be the ratio of "CO₂ remaining in the atmosphere" divided by "the anthropogenic CO₂ supplied" (Broeker 1975). It has been suggested that the reason for its constancy is related to the solution of an ordinary differential equation as applied to atmospheric CO2 flux (Cawley 2011). According to this view, if the time constant for the exponential rate of rise of CO2ff matches the time constant for the rate of depletion of CO2ff then the Airborne Fraction would be 0.5. However this paper suggests that the origin is quite different and is related to the direct uptake of CO₂ before the remainder is mixed into the atmosphere (Figure 1).

5.4. Direct Uptake

The presence of a direct uptake mechanism (as shown in Figure 1) significantly improves the model's accuracy, reducing the overall standard deviation by approximately 30%. Recent studies of daytime releases of ¹⁴CO₂ in the vicinity of nuclear power plants have reported significant ¹⁴CO₂ uptake (Naegler & Levin 2006) while a study of urban grasses near a major highway reported plants stored up to 13% of fossil-fuel carbon (Lichtfouse et al., 2005; Ota et al. 2016). According to Kuderer et al. (2018), "The 14C signals from such point sources are well detectable in plant samples". Therefore it seems reasonable to propose the existence of a direct uptake mechanism for CO2ff before it becomes mixed in the atmosphere. The effect of the direct uptake mechanism is to more rapidly transfer part of the rising CO2ff emissions to the reservoir, hence slightly changing the shape of the ¹⁴C and ¹³C curves.

5.5. 14 C Bomb Pulse

Regarding the shape of the bomb pulse, it would be expected that, if there is a constant level of washout exchange flux, its shape would be a reducing exponential converging back to its original background value. However simple inspection reveals this not to be the case. See Figure 4. This shape departs from an exponential decay because of three main reasons; a) incoming fossil fuel (free of ¹⁴C) is diluting the atmosphere, this being called "Suess Dilution" (Suess 1955), b) some of that incoming fossil fuel is also being washed out by the exchange flow and c) some of the incoming exchange flow has a higher level of ¹⁴C, because bomb ¹⁴C having accumulated in the oceans is now being re-emitted back to the atmosphere. If (a) was the sole cause of the deviation of the bomb curve shape from a pure exponential decay, the curve would be well below Δ¹⁴C=0, bearing in mind that the atmospheric growth of 30% is attributed to fossil fuel emissions containing no ¹⁴C. The effect of washout of fossil fuel emissions (b) raises the curve; this was noticed many years ago and termed by Stuiver and Quay (1981) the "attenuation factor". However the raise is still not quite sufficient to make a perfect match. Finally effect (c) raises the tail slightly further improving the fit to the current level shown in Figure 4.

Keeping track of all of these factors in a "back of the envelope" calculation is complex. The absolute flow model described calculates at each iteration ¹³C and ¹⁴C in reservoir and atmosphere (see Implementation section). It does not need to explicitly consider Suess dilution, the Stuiver attenuation factor, flow-back or the airborne factor as these emerge implicitly from the calculation.

5.6. Amount "Remaining" in the Atmosphere

The ambiguous meaning of "remaining" has led to some confusion between different authors; the situation is clarified here by extending the above island analogy further. Suppose one of the above islands, island A has an annual supply of people arriving via its Airport, increasing its population. The ferry continues to carry exactly the same number of people back and forth between island A and island B. The population of the island A increases, yet over time the ferry has carried many of the airport arrivals to island B, and has brought inhabitants back from B to A in their place. A census would report that the number of people on island A who flew in, is less than its population growth because some are living on island B. However others would say that airport arrivals caused 100% of the growth of A since without airport arrivals their population would not have grown. This illustrates the discrepancy. IPCC state "The combustion of fossil fuels and land-use change for the period 1750–2019 resulted in the release of 700 ± 75 PgC (likely range, 1 PgC = 10¹⁵ g of carbon) to the atmosphere, of which about 41% ± 11% remains in the atmosphere today (high confidence)" (IPCC 2023 p80). However the amount "remaining" of CO2ff has been found by Skrable to be 23% and in this paper to be 27%. The apparent discrepancy arises because a) a significant part of CO2ff causing the growth has been washed out and almost exactly replaced by the exchange flow, and b) because the IPCC figure includes an allowance for Land Use Change (LUC) which cannot be accounted for by the isotopic method because of its unknown isotopic signature (although LUC has significantly reduced in recent years). The alternative measure of molecular CO₂ remaining relative to atmospheric CO₂ content was provided by Skrable (12%), Stallinga (<10%), this paper (14%) providing reasonable agreement between the three authors.

6. Conclusions

A two box CO2 absolute flow model is described which

• accurately predicts the values of Δ14C and δ13C over 200 years.
• revises the view of the airborne fraction, proposing that it reflects the relative amount of CO2ff that is absorbed directly into the terrestrial/oceanic reservoir,
• does not use or require a consideration of carbonate seawater chemistry (Revelle),
• shows there is no practical significant difference in behaviour of different isotopes apart from fractionation
• explains how exponential decay time of the bomb pulse does, after all, relate to the residence time,
• resolves the conflicting calculations of how much CO2ff "remains" in the atmosphere.

The method produces annual estimates of both nett inflow and the individual exchange flows, making it possible to calculate the molecular CO2ff remaining and also to provide the amount of atmospheric growth due to CO2ff. Thus the model leads to new estimates for the period 1750 to 2020 of the partitioning of anthropogenic CO₂ emissions, with 88% of growth being attributable to CO₂ fossil fuel emissions, and 12% being due to net "natural" emissions from the terrestrial and oceanic source, including land use change. Relative to total CO2ff of 454 GtC, 253 GtC or 55% was accountable for the growth of the atmosphere. However, at the molecular level, because of the washout caused by exchange flow, only 27% of fossil fuel emissions remain airborne, which corresponds to 14% of the entire CO₂ atmosphere. The model uniquely provides both figures thus helping to resolve any confusion regarding the meaning and effects of fossil fuel CO₂ emissions.