A simple metabolic model of glacial-interglacial energy supply to the upper ocean

Introduction Conclusions References


Introduction
A sustainable natural system, from individuals to large ecosystems, must provide energy fast enough to match its power requirements at any time.Bejan (1997Bejan ( , 2000) ) has proposed that different complex systems have analogous spatial and temporal archi- patterns during state transitions is frequent in systems of very different nature (Scheffer et al., 2009).There is no reason to expect the earth system should be an exception.The (yearly and globally integrated) Sun's radiant energy over the earth has remained nearly invariant through the last several millions of years.What has changed is the rate at which this energy is incorporated, transformed and distributed to the earth system, the result being that the earth's climate has switched between different metabolic states.These changes are reflected by substantial changes in key thermodynamic earth variables, e.g. the globally-averaged surface earth temperature increases several degrees Celsius from glacial to interglacial periods (e.g.Fig. 6.5 in IPCC, 2007).At temporal scales between order 10 yr and 100 kyr, the changes in the earth's thermodynamic variables are mainly related to the flow of energy through the ocean compartment (Sarmiento and Toggweiler, 1984;Siegenthaler and Wenk, 1984;Sigman and Haug, 2003;Sundquist and Visser, 2003).The Sun's radiant energy is the only significant energy source but, when incorporated to the earth, it takes many forms that arise from physical and chemical transformations, the principal ones being internal heat, latent heat, mechanical energy, and chemical energy following photosynthesis.Mechanical energy is the combination of kinetic and potential contributions, with kinetic energy eventually being dissipated as heat at the Kolmogorov scale.The incoming radiant energy is stored as heat or converted into organic matter after primary production, but this only takes place near the ocean's surface; the efficiency of photosynthetic conversion of solar radiation is set by the supply of nutrients.The distribution of all these energy forms to the whole system relies on the ocean circulation, both as surface and thermocline (upper-ocean) and thermohaline (deep-ocean) branches, which covers many temporal and length scales.Substantial variations in the thermohaline circulation lead to very different rates of transfer of these different energy forms between the upper and deep oceans.James Lovelock and others have searched for self-regulating mechanisms of the earth system with moderate success, undoubtedly limited by the system non-linear complexity (Lovelock, 1972(Lovelock, , 1989;;Lovelock and Margulis, 1974;Margulis and Figures Back Close Full  Lovelock, 1974;Watson and Lovelock, 1983;Charlson et al., 1987).Their most outstanding achievement has been raising to scientific debate the concept that the earth's homeostasis is the result of life's modulating force.An alternative path to explore the living earth is to ignore the individual mechanisms or subsystems and, rather, to examine variables that describe the behaviour of the whole system, in a similar way as physiology studies the nature of organs or systems in living beings.Such a macroscopic view is not new as it is already the basis of a whole discipline, allometry (Peters, 1983;Calders, 1984;Schmidt-Nielsen, 1984), and of novel ideas by diverse authors (Bejan, 1997(Bejan, , 2000;;Enquist et al., 2003;Brown et al., 2004;West and Brown, 2005;Reis and Bejan, 2006;L ópez-Urrutia et al., 2007;Kirwan, 2008).Pelegrí (2008) has also proposed to consider the macroscopic properties of the ocean, as the dominant component of the earth system at glacial-interglacial timescales, using a physiological perspective.It differs from other efforts in its utilization of analogous physiological behaviours, searching for parallelisms with physiological responses in complex living organisms.The undergoing idea is that awareness on how energy supply and metabolic response operates in living beings may be helpful to better understand the response of the ocean system.Pelegrí (2008) chose mammals because they are a profoundly studied complex system and, certainly, because everyday we evidence its macroscopic behaviour, i.e. each of us experiences how our own body responds.Here we will continue to use mammals as a source of inspiration, but recognizing there likely are other beings with a much closer earth-like behaviour.
The energy distribution in mammals and the ocean takes place through circulatory systems with substantial differences but, nevertheless, remarkable functional analogies.Firstly, both mammals and the ocean have short and long circulation branches.The short circuit in mammals is the pulmonary branch while the long one is the systemic branch; for the ocean the short circuit would be the upper-ocean branch while the long one would be the deep-ocean branch (e.g.Campbell, 1990;Randall et al., 2002;Guyton and Hall, 2005;Pelegrí, 2008, and references therein).The short circuit plays an important role for gas exchange while the long one replenishes the nutrient levels, Figures

Back Close
Full ensuring large enough nutrient delivery to the whole system.Secondly, for both systems fluxes are initially driven by surface-force gradients, set during a relatively short pumping phase.These gradients are mechanically-driven pressure gradients for mammals and induced by wind stress and density gradients for the ocean.And, thirdly, both systems deliver nutrients and gases through a combination of major advective currents at long scales with diffusion processes at short scales.The major delivering paths in the ocean are nearly-closed streams (Pelegrí and Csanady, 1991), resembling the closed blood conduits in mammals.
Another thing we learn from mammals is that an increase in energy requirements is initially supplied by a proximal reservoir through anaerobic processes, as the circulatory system takes some time to accommodate the new energy necessities.Since this reservoir is limited in size its supply cannot last very long.Eventually the energy distribution rate through the circulatory system is large enough so that aerobic takes over anaerobic supply.Pelegrí (2008) proposed the earth's state depends on the total energy supply (carbon per unit time), both as readily available energy (local sources) and as potentially available energy (advection).Semi-labile or ancient dissolved organic carbon within the upper ocean is rapidly available, although its rate of transformation slows down from labile to recalcitrant forms.Deep inorganic carbon and nutrients, on the other hand, may reach the upper ocean and get transformed through photosynthesis, the limitation being on the rate of upwelling of these deep waters.
Here we apply the physiological approach to examine glacial-interglacial changes by means of a very simple metabolic model for the ocean, inspired on the way energy flows in mammals.We call it a metabolic model as it assumes the earth system switches between two different states through variations in the rate of net energy supply (carbon flow) to the upper ocean, expressed as the product of two main parameters, flux rate and carbon-concentration difference.We pursue the argument initially proposed by Pelegrí (2008), now with a much improved conceptual justification for the model, a sensitivity analysis on how the model responds to the two controlling parameters, and a careful discussion of the results.The model aims at predicting the temporal Introduction

Conclusions References
Tables Figures

Back Close
Full changes in the upper-ocean dissolved inorganic carbon (or nutrient) concentration during the last four interglacial-glacial cycles, tuned to follow the same pattern as the Vostok time series of atmospheric CO 2 (Petit et al., 1999;Siegenthaler et al., 2005;L üthi, 2008).This allows us to place constraints on past changes in the thermohaline circulation rate.

Rationale behind a metabolic model for the upper ocean
At steady state, the metabolic rate for a heterotroph is the rate of respiration.Brown et al. (2004) argued the rate of photosynthesis is the inverse reaction to the rate of respiration and defined it to be the metabolic rate for an autotroph.These authors, however, did not clarify whether this rate of conversion of solar energy should take into account the autotroph's respiration, i.e. whether metabolic rate would be equivalent to either gross or net autotrophic primary production.The answer is not trivial: the basal metabolic rate for a heterotroph is the energy necessary to sustain its vital functions so an analogous definition for an autotroph should include respiration; with this perspective the metabolic rate of an autotroph would be gross primary production.However, if we consider the metabolic rate as the net inverse respiration, or net photosynthesis, then it should be net primary production, i.e. respiration is intrinsically linked to photosynthesis and subtracts to the final rate of conversion of solar energy.In this work we will use this last definition: the metabolic rate of an autotrophic ecosystem is the net autotrophic primary production.
We may extend these ideas to the earth system, its metabolic rate being the amount of solar energy effectively transformed, i.e. the planet's net autotrophic primary production.The transformation of solar energy relies on the supply of inorganic carbon and nutrients, which comes from new (external) and regenerated (internal) supplies.The principal internal source is respiration by the heterotrophic community as it provides regenerated inorganic carbon and nutrients so that more solar energy may be photosynthesized.Metabolic rate, or net autotrophic community production, must no Figures

Back Close
Full be confused with metabolic balance, the later being a comparison between gross primary production and community respiration, or net community production (Duarte et al., 2004;L ópez-Urrutia et al., 2007).The metabolic balance may be a negative quantity but metabolic rate cannot.Temporal changes in the metabolic state of the upper-ocean (surface and permanent thermocline waters) may result either from alterations in the use of organic and inorganic reserves within the upper ocean (proximal fast sources) or from variations in the supply of deep ocean carbon and nutrients (a relatively slow advective supply).
Hereafter we will always refer to dissolved inorganic carbon (DIC) in the upper ocean although the same ideas should equally apply to dissolved inorganic nutrients (Pelegrí, 2008).We must keep in mind that carbon is always in excess at the sea surface so that inorganic nutrients, with substantial concentration differences between the upper and deep ocean, are the truly limiting factor for primary production.Hence, we should think of carbon as accompanying the nutrient supply in proportions close to the stoichiometric Redfield ratios.The reason we aim at modelling DIC is simply because we have the paleo-record of atmospheric CO 2 which will be used to calibrate a simple two-state metabolic model.

Analogous energy supply mechanisms in mammals and the earth system
Let us start with an analysis of the similarities between energy supply and utilization in mammals and the earth system (Fig. 1).Such an analogy will show to be illuminating and becomes a key motivation to build an idealized metabolic model for the ocean (Sect.3).Mammals are heterotrophic systems which rely on intake of organic matter and its oxidation into inorganic matter to produce energy.At steady state the metabolic rate ρ of a mammal is the rate of energy transformation, through oxidation of organic matter, given by Full where q is the blood flow rate, φ d is the blood-oxygen concentration arriving through the arteries to a muscle or body subsystem, and φ is the blood-oxygen concentration leaving it, with the s subindex denoting steady-state conditions; the quantity φ d − φ is known as the arterial-venous oxygen difference.The metabolic rate depends on the rate at which unused oxygen leaves the subsystem, φ s q s .We may rewrite Eq. (1) as the following balance: the term in parenthesis being net energy supply.There are two extreme situations: (i) φ s = 0, all the oxygen flux is used and the metabolic rate equals the rate of energy supply, ρ s = q s φ d , and (ii) φ s = φ d , no energy is transformed so ρ s = 0; in a mammal φ s is substantially less than φ d but never approaches zero.When a mammal switches state, from basal to enhanced or vice versa, the new energy supply has to change abruptly to sustain the new energetic requirements.During the transition, however, the blood-oxygen delivery changes progressively through both the blood flow rate and the oxygen concentration difference, so the initial difference in energy supply must come from anaerobic energy sources or sinks ρ ana .Some of the delivered oxygen may not immediately be converted and will accumulate or, contrarily, oxygen reserves may be used, we call this term d φ d V /d t.Therefore, during the transition we may write The earth, on the other hand, is an autotrophic system whose source of energy is the Sun's radiation, which results in the transformation of inorganic into organic matter through photosynthesis.The earth's metabolic rate, or amount of solar energy transformed, is a direct function of the available inorganic nutrients.The "organ" we examine is the upper ocean, composed by surface and thermocline waters, where about half of the earth's primary production takes place.The recirculation within this upper ocean Figures

Back Close
Full is relatively fast, of the order of 10 years, so we assume this compartment gets homogeneous at glacial-interglacial time scales.Let c d be the concentration of inorganic matter that returns from the deep to the upper-ocean layers.Some, but not all, of this DIC is converted into organic matter: the unused portion has concentration c and eventually escapes during deep-water formation.Therefore, the source of new DIC is where Q is the water transported by the thermohaline circulation.The autotrophic upper-ocean photosynthesizes this deep DIC supply into organic matter at a rate P , although in the process it also uses some organic matter through autotrophic respiration, R aut ; hence, the metabolic rate, or net autotrophic community production, is given by M = P − R aut .Therefore, at steady state the supply of DIC by the deep ocean sustains the metabolic rate, where the s subindex denotes steady-state conditions.As in mammals, the steadystate metabolic rate is controlled by the DIC leaving the upper ocean, c s Q s .Equation (3) may be rewritten as M ana .During the transition, some DIC may accumulate or be extracted from an upper ocean of volume V , this term being d (cV )/d t.Therefore

Advective versus local energy sources
The advective or slow energy source for a heterotrophic system is oxidation of organic matter sustained by the blood flow while for the earth (a long-term autotrophic system) it is the photosynthetic transformation of inorganic matter arriving through the thermohaline recirculation.Pursuing the physiological analogy we name this long-term energy source to be an aerobic supply (although this term is usually associated with oxidation, here it will only imply sustained or long-run maintenance of certain energy state).As a system changes from basal to enhanced states the advective supply increases notably, in the earth this has been reflected by large glacial-interglacial changes in the intensity of the thermohaline circulation (Shackleton et al., 1983;Imbrie et al., 1992;Labeyrie et al., 1992).
We turn now to examine the proximal fast-energy sources.Since these are not directly related with the delivery by the circulatory system we call them anaerobic supply.
Let us first consider mammals.When the transition is towards a more energetic state, the initial aerobic supply is insufficient to sustain the enhanced energy requirements so the system also relies on a temporal anaerobic source, ρ ana > 0. In a mammal this essentially is reserve adenosine triphosphate ATP molecules in cells (the phosphagen energy system).During a transition to a less energetic state, the decrease in aerobic supply does not instantaneously match the new energetic expenditure; as a result the aerobic energy supply temporarily exceeds the new energy requirements and an energy sink replenishes the anaerobic reserves, ρ ana < 0. Equation ( 2) may be rewritten as the following balance:

Conclusions References
Tables Figures

Back Close
Full As steady state is approached the temporal derivative goes to zero and so does anaerobic supply, and this equation becomes Eq. ( 1a) with the new φ s , q s and ρ s values.A comparison between both equations shows the instantaneous net supply to be the parenthesis term in the right-hand-side of Eq. (2a).In order to sustain the new metabolic state the net supply has to remain approximately constant, equal to the steady-state value: An analogous response applies for the upper ocean.During increased energy requirement the initial advective supply is temporarily insufficient, so the system uses energy reserves, M ana > 0. Similarly, during decreased energy requirement, the advective supply will decrease smoothly and temporarily exceed the requirements, resulting in an energy sink, M ana < 0. We propose these processes take place through remineralization of DOC by the heterotrophic community during a transition to a more energetic state, or through accumulation of DOC reserves when changing to a less energy-demanding situation, so we may also write M ana ≡ R het .
Equation ( 4) may be expressed as When reaching steady state this equation reduces to Eq. (3a) with the new c s , Q s and M s values.Physiology tells us that as soon as the system changes state, i.e. during the transition and as long as the new state lasts, the net supply remains stable: Recapitulating, the expressions that set the energy requirements for mammals (Eqs. 1 and 2) and for the upper-ocean system (Eqs.3 and 4) are completely analogous: The steady-state rate of organic-matter respiration in mammals and photosynthetic production in the earth system are attained through sole advective supply.The out-ofsteady-state metabolic rates are attained through both advective and local supply.Figure 1, a classical figure in physiology textbooks (e.g.Campbell, 1990;Randall et al., 2002;Guyton and Hall, 2005), illustrates aerobic and anaerobic supply and metabolic expenditure in either mammals or the earth system.A change in metabolic state occurs almost instantaneously but the aerobic energy supply does not, so the system recourses to anaerobic reserves.During a transition from basal to enhanced states the energetic deficit (energy expense less aerobic supply) is much less than the energetic excess (aerobic supply less energy expense) during the opposite transition; this is because the system is best adapted to heighten its metabolic rate using proximal energetic sources.This results in a marked asymmetry between these two transitions which is reflected in many physiological variables.One relevant example in mammals is the arterial-venous oxygen concentration difference; Pelegrí (2008) has argued that the concentration of atmospheric CO 2 is the analogous variable for the earth system.

Model description
Let us formalize the above ideas for a two-box ocean model, as shown in Fig. 2 (Pelegrí, 2008).The model is extremely idealized as it ignores any chemical ocean changes, such as the deep-ocean carbonate system, and does not consider the high-nutrient and low-chlorophyll Southern Ocean.Further, it does not specifically address export of particulate organic matter.The model, therefore, focuses on the Atlantic high-latitude physical pump, with a rudimentary biological pump through deep-convection of DOC.
The model will be used to obtain the evolution of DIC in the upper ocean during the last four glacial-interglacial cycles, tuned to follow the Vostok record of atmospheric CO 2 .Since it ignores some potentially important mechanisms we should not expect very good hindcasting skills.Therefore, if the hindcast turns out to be satisfying, it is Introduction

Conclusions References
Tables Figures

Back Close
Full possible that some of the neglected mechanisms are somehow embedded in the model dynamics.We will come back to these considerations in Sect. 5.

Idealized model for inorganic carbon (or inorganic nutrients)
We idealize the ocean as composed of upper and deep compartments.The upper compartment contains surface and thermocline waters, in direct or indirect contact with the atmosphere (the latter via Ekman pumping).Water from the upper compartment sinks into the deep compartment at high latitudes and eventually recirculates back into the upper ocean, a portion through the equatorial ocean.The rate of exchange between the two compartments will change in time through variations in the intensity of the thermohaline recirculation (Stommel, 1958;Broecker, 1991).Some water also recirculates within the upper compartment but much faster than the water exchange between the two compartments, so at glacial-interglacial time scales we consider the upper-ocean recirculation is solely responsible of homogenizing the carbon and nutrient concentrations in the upper compartment.We assume DIC in the deep compartment to be constant in time but we let DIC and dissolved organic carbon (DOC) in the upper compartment to change in time.After arrival of deep water to the upper ocean, some of its inorganic carbon and nutrients are photosynthesized and converted into organic, while some are left unused and return back to the deep ocean during winter deep-water formation.The equations for DIC, c, and DOC, c r , in the upper ocean are: where Q(t) is the overturning rate, V is the upper-ocean volume, c d is the DIC of the deep ocean, and R and P are the remineralization and primary production of the upper ocean, respectively; the s subindex denotes steady-state conditions.If we split respiration into autotrophic and heterotrophic fractions, R = R aut +R het , and recall M = P − R aut and M ana = R het , then P − R = P − R aut − R het = M − M ana , so Eq. ( 7) is identical to Eq. ( 6).Therefore, we may substantially simplify Eq (7) by setting an approximately constant net energy supply throughout each metabolic state, i.e. we let the energy supply to the upper ocean S ≡ R − P + Q c d be set by the metabolic requirements, indistinctly supported by new and regenerated DIC.The DIC equation becomes For simplicity we set S proportional to the overturning flow rate Assuming the upper ocean has a constant volume V , Eq. ( 9) becomes The steady-state form of Eq. (9a) shows that a constant g = c s sets the concentration level attained by the upper-ocean after sufficiently long time, Q s c s = Q s c d − M s , which is exactly the steady-state form of Eq. ( 6).The steady-state concentration level may be written as c s = c d − M s /Q s .Since long-term autotrophy requires M s > 0 then, as expected, c s < c d .Further, the greater the recirculation rate Q s the larger the steadystate concentration level will be, in consonance with the idea that a swifter thermohaline circulation allows the deep ocean waters to reach the sea surface with less dilution.
Let the non-dimensional variables (primed) be defined as c = c d c , t = τ 0 t and Q = Q 0 Q , where τ 0 and Q 0 are the overturning basal time scale and basal rate, respectively, related through τ 0 = V/Q 0 .We may further define a non-dimensional concentration level γ, so g = c d γ.In these relations we take c d , τ 0 and Q 0 as constants.

Conclusions References
Tables Figures

Back Close
Full Equation (9a) becomes non-dimensional as follows (dropping primes): The main advantage of this equation is its simplicity; its solution depends solely on two non-dimensional quantities: concentration level γ and recirculation rate Q.The steadystate dimensional concentration level g = c s should always be less or equal than the deep-sea concentration c d , in concordance with the idea that the deep-ocean is the only long-term source of inorganic carbon, so γ has to be equal or less than one.
Similarly, the dimensional recirculation rate will always be greater than the dimensional basal recirculation rate Q 0 so that the non-dimensional variable Q in Eq. ( 10) is always greater or equal than one.We propose next that the minimum recirculation rate occurs during glacial periods and the maximum one near the interglacial peak.

Model parameters
The solution of Eq. ( 10) depends on the rate of energy supply through two non-dimensional parameters: the recirculation rate or inverse recirculation time, , and the concentration level, γ(t).To set the energy supply we pursue the physiological approach: as exercise begins the energy requirements increase abruptly while as exercise ends these requirements also decrease sharply, so S (t) must resemble a step function.The energy level is set by the aerobic supply at steady-state, the dimensional relation being S s = Q s c s (Eq.9).Let us consider one basal-enhanced energy supply cycle, or rest-exercise cycle, with the system switching from a basal energy requirement before t 1 to the enhanced state between t 1 and t 2 , to return to the basal state afterwards.We specify this switch by letting the system aim towards a constant non-dimensional concentration level, i.e. we let the normalized concentration level change as a simple step function: where γ 1 and γ 2 respectively set the rest and exercise levels.
The actual pace at which the system aims towards an aerobically-sustained level is defined by the flow rate, Q(t).Therefore we set γ 1 and γ 2 as constant, to be specified through data fitting (see below), and specify the changes in aerobic supply through the non-dimensional flow rate as follows: where the non-dimensional recirculation rate anomaly ε(t) is zero for t < t 1 and positive when t 1 ≤ t < t 2 .In a two-state cycle, the steady-state non-dimensional flow rate Q s will either be its basal value, Q 0 , or an enhanced value Q 0 (1 + ε max ), where is the maximum circulation rate anomaly which occurs at some time t max during the interglacial period; during a basal to enhanced state transition the recirculation time will decrease by the same (1 + ε max ) factor.

Modeling the Antarctic CO 2 time-series
The above metabolic model may be applied to either inorganic nutrients or carbon, but here we will focus on carbon as we can tune our output by means of the Vostok atmospheric CO 2 time series.The implicit premise is that changes in atmospheric CO 2 mirror the variations of DIC in the upper ocean (Sarmiento and Toggweiler, 1984;Siegenthaler and Wenk, 1984;Toggweiler, 1999).However, rather than modelling the absolute value of the carbon concentration we look at its anomalies referred to some background state.By doing this we are effectively looking at carbon accompanying the limiting nutrient supply.Specifically, we define c = c 0 + ĉ and g = c 0 + ĝ and substitute back into Eq.(9a), so that this equation may now also be understood as an anomaly equation (after dropping hats) and Eq. ( 10) is the corresponding nondimensional anomaly relation.
Hereafter we will assume that the dissolved inorganic carbon in the upper ocean is linearly proportional to the atmospheric CO 2 value.In order to tune our DIC anomaly

Conclusions References
Tables Figures

Back Close
Full model (Eq.10) we use a non-dimensional, or normalized, atmospheric CO 2 signal anomaly, derived from the Vostok time series for the last four interglacial-glacial cycles, i.e. between years 437 975 and 20 755 before present.To produce this normalized anomaly time series we follow three steps: (a) select the smallest and largest values of the original time series, (b) calculate the anomalies with respect to the smallest value, and (c) normalize these anomalies by the difference between the largest and smallest values.This normalized time series is shown as Fig. 3, it illustrates four interglacialglacial cycles which we name cycles 1 to 4, chronologically in time, with cycle 4 the most recent one.In this figure we also show those periods with enhanced metabolic activity, which we associate with the interglacial times.The normalized energy function γ raises from a basal-glacial γ 1 = 0 value to an enhanced-interglacial level γ 2 , the latter being constant during each interglacial period but changing from one cycle to another.Figure 3 also shows a detailed view of cycle 4, with the CO 2 normalized anomaly time series together with a smoothed time series produced with a 2 kyr moving average.This figure shows those times when the system changes from glacial to interglacial, t 1 , and back to glacial, t 2 .We use very simple criteria to delimit the interglacial period.The start time t 1 corresponds to the last time when the derivative of the smoothed time series becomes zero before the original time series overcomes some threshold, here taken as 75% of the interglacial maximum value.One exception occurs if the last two minima are separated by less than 10 kyr, in this case t 1 is taken to be the average between these last two minima.Similarly, the end time t 2 is taken to be the last time when the derivative of the smoothed time series becomes zero before the original time series decreases towards glacial values, here again set as less than 75% of the interglacial maximum value.Figure 3b also illustrates the criterion used to set γ 2 .This is taken to be the mean value of the 2-kyr smoothed time series between t peak and t 2 , with t peak being when the smoothed time series takes the maximum value.Introduction

Conclusions References
Tables Figures

Back Close
Full

Model sensitivity: test cases
In last section we argued the model depends critically on two non-dimensional parameters: the energy level γ and the recirculation rate Q.Let us now examine the sensitivity of the results to simple expressions for these two parameters.First, we let the energy level γ change as a step function between 0 and γ 2 , with each individual interglacial energy level set as explained above (Eq.11).Second, we let the circulation rate be a function of the basal recirculation time τ 0 and the circulation rate anomaly ε, Q = V/ τ 0 (1 + ε(t)) (Eq.12).The circulation rate anomaly ε increases linearly from ε(t = t 1 ) = 0 to a maximum ε(t = t max ) = ε max and then back to ε(t = t 2 ) = 0.During interglacial states, such as occurring nowadays, the maximum flow rate is Q max = Q 0 (1 + ε max ).Notice t max does not need to coincide with the time when the time series takes its maximum value t peak , the latter is directly obtained from the smoothed CO 2 time series while the former is determined through fitting the model results to the CO 2 time series as discussed below.
To carry out the above procedure we estimate the upper-ocean area to be 3.3 × 10 14 m 2 , or about 65% of the earth's surface (71% of water coverage less some 6% for the continental shelves), and set a mean depth of 1400 m to get V = 4.7 × 10 17 m 3 .For the maximum rate of deep water formation we set Q max = 22 Sv (1 Sv ≡ 10 6 m 3 s −1 ), as a rough estimate obtained from today's rate of deep water formation (for a review see Pelegrí, 2008).The corresponding minimum recirculation time is τ min = V/Q max = 700 yr.We assume Q max has been the same for all interglacial periods and apply our numerical model for different τ 0 values with Choosing such a maximum constant rate is an additional approximation, necessary for the model to run.However, it may be again justified based on physiological constraints: the heart of a mammal has a maximum cardiac output which changes slowly during its life time; in the same way the maximum deep-water formation rate is likely constrained to be within some narrow range at each time throughout the earth's history.Introduction

Conclusions References
Tables Figures

Back Close
Full Our strategy is to explore the model's sensitivity to substantially different values of τ 0 (or V/Q 0 ), with different combinations of t 1 , t max and t 2 (Fig. 4).We have tried many possibilities, but the most satisfactory results correspond either to t max = (t 1 + t 2 )/2 or t max = t 2 (hereafter to be named Cases 1 and 2, which respectively correspond to panels b and d in Fig. 4).Several other combinations, such as those shown in panels a and c in Fig. 4, are significantly less satisfactory and are not shown.Once we get a feeling of the overturning basal time scale, we will refine the results by a more specific optimization method (next Section).
Figures 5 and 6 respectively illustrate the results of Cases 1 and 2 for three different τ 0 values: 10, 40 and 100 kyr.For the sake of clarity, in these and the remaining figures we recover the dimensional time scale but plot the non-dimensional values for all dependent variables.The correlation coefficient between the CO 2 original time series and the DIC modelled values serves as a useful indicator for the goodness of the solution.We find that the best fit occurs for τ 0 = 40 kyr, with correlations of 0.769 and 0.794 for Cases 1 and 2, respectively.The modelled time series cannot match the abrupt changes in the time series, which take place at time scales between several hundred and a few thousand years, but do reproduce quite well the long-scale (10 to 100 kyr) oscillations.

Selection of parameters through genetic algorithms
From the above preliminary analysis we conclude that both Cases 1 and 2, with an overturning basal time scale of τ 0 = 40 kyr, do reproduce reasonably well the long-term temporal changes in the CO 2 original time series.We next use genetic algorithms (GA) to search the optimal τ 0 and t max values.
A GA is a search technique used to find globally optimized solutions (see reviews by Goldberg, 1989;B äck, 1996;Fogel, 2000).In our model we have one single differential equation whose output is controlled by several parameters, each parameter supposed to be within a prescribed range of values.We can define a fitness function that compares the model's output with an available empirical time-series, the goodness

ESDD Figures Back Close
Full of the comparison depending on the set of parameters chosen to run the model.The GA looks for the combination of parameters which maximizes such a fitness function, i.e. minimizes the residual.The evolution usually starts from a large population of randomly generated combinations of our parameters (individuals), covering the entire range of possible solutions.In each iteration of the algorithm (generation), the fitness (goodness of the model result) of every individual in the population is evaluated.Then, individuals among those with highest fitness are stochastically selected from the current population and modified (recombined) to form a new population, which is used in the next iteration of the algorithm.The algorithm terminates correctly when the fitness function reaches a satisfactory level for the whole population, or possibly incorrectly when it exceeds some pre-established maximum number of iterations.
We allow the GA to search for τ 0 during each individual glacial-interglacial cycle, with a starting value of 40 kyr.Inspired by the relatively good results for Cases 1 and 2, we further ask the GA to search for t max in two different ways.In the first one t max remains between t 1 and t 2 , with ε(t) increasing linearly between t 1 and t max and then decreasing linearly between t max and t 2 (panel e in Fig. 4).In the second one we allow t max to happen after t 2 , with ε(t) increasing linearly between t 1 and t max but decreasing to zero immediately afterwards (panel f in Fig. 4).
The GA is executed 400 times and the correlation between the numerical output and the normalized CO 2 anomaly time-series is computed.The high correlation simply reflects good agreement between the modelled DIC and the atmospheric CO 2 for each individual glacial-interglacial cycle.The parameters for the best run (correlation 0.813) are presented in Table 1 and the corresponding numerical solution is shown in Fig. 7.Those realizations with correlation above 0.81 (a total of 38) are chosen to get the mean values and standard deviations of the controlling parameters (Table 1).A plot of the numerical output for the mean values in Table 1 is nearly indistinguishable from that in Fig. 7 and is not shown.
The above procedure tells us the size of several important quantities for each of the four glacial-interglacial cycles: the width of the enhanced (interglacial) periods,

Conclusions References
Tables Figures

Back Close
Full t 2 −t 1 ; the relative height of these enhanced states, γ 2 ; the location and intensity of the maximum recirculation rate, t max and ε max ; and the overturning basal time scale, τ 0 .It is tempting to use one or more of these parameters to build simple functions for the nondimensional energy supply γ and the circulation rate anomaly ε. Figure 8 illustrates the results with such a simple function, its design based on the best GA execution (Table 1, Fig. 7), as explained next: (a) the same t 1 and t 2 values are selected, (b) for the interglacial ε(t) we use a Gaussian function whose width (standard deviation) and height is given by the width and height of the γ = γ 2 period, (c) the peak value of this Gaussian function is located at the middle of this enhanced state, and (d) for γ (t) we have chosen exponential growing and decaying functions with decay times given by τ 0 .
We may appreciate that such a simple choice also leads to a high correlation, close to what we found with the GA.In particular, now the solution displays two different rates of DIC change during the interglacial-glacial transitions, a rapid but relatively short initial decay followed by a much longer and slower decrease, similar to the atmospheric CO 2 pattern in the Vostok's time series.
The results obtained in Fig. 8 suggest that the numerical solution could probably be improved by properly tuning the controlling parameters.However, here we do not explore this any further as it goes beyond our initial objective of having a model as simple as possible, capable of capturing the essential mechanism controlling the glacialinterglacial transitions.thanks to the shape of the circulation rate anomaly ε, which results in effectively fast circulation rates during interglacial periods and slow rates during glacial periods.The model, however, fails at simulating the intermittent abrupt changes that take place at time scales of order 100 years to 1 kyr.These abrupt changes may become very intense, particularly when the system goes towards the deep glacial state, and are the reason why the ocean-DIC -atmospheric-CO 2 correlation cannot become much larger than about 0.8.One remarkable result is that ε max gets quite high for all four cycles (Fig. 7 and Table 1); specifically, during the interglacial periods the circulation rate increases by a factor of about 60 as compared with the basal circulation rate.This number depends on the selected values for nowadays rate of deepwater formation, as well as on the upperocean volume, so its error bars are probably as large as a factor of three.If today's rate of deep-water formation was 50 Sv (instead of 22 Sv) and the upper-ocean had a mean depth of 2000 m (instead of 1400 m) ε max would still be about 20.A large ε max value is in qualitative agreement with indications that during interglacial periods the meridional overturning circulation has been more intense than during glacial periods (Shackleton et al., 1983;Imbrie et al., 1992;Labeyrie et al., 1992), although there are no clear evidences of such a very large change.This result suggests there are several oversimplifications in our model, which we briefly consider next.

Discussion
One main approximation is the absence of the relatively well-ventilated and nutrientrich high latitude waters, which play an important role in the efficiency of the physical pump.The current interglacial Southern Ocean is apparently much less stratified than during the last glacial maximum, as a result of the formation of relatively fresh Ice Shelf Water which influences on the characteristics of the Antarctic Bottom Waters (Nicholls and Osterhus, 2004).Therefore, during an interglacial maximum the Southern Ocean deep DIC-rich waters are more directly connected to the sea surface than during the glacial maximum, leading to natural CO 2 outgassing.The effect of stratification in the Southern Ocean was included in a box-model of Pleistocene climatic oscillations by Paillard and Parrenin (2004).Despite its simplicity, this model was the first one to Introduction

Conclusions References
Tables Figures

Back Close
Full accurately simulate the timing of all the glacial terminations as well as all maxima and minima of the δ 18 O oscillations observed in the past 3.5 million years.The incorporation of a third box in our model, i.e. a high-latitude compartment such as in Sarmiento and Toggweiler (1984) and Siegenthaler and Wenk (1984), should help clarify these ideas.
Another important limitation in our model is the lack of particulate organic matter export.Such export is partly responsible for maintaining the large DIC and inorganic nutrient concentrations in the deep ocean, something we waive by simply setting constant concentrations in the deep ocean.This is an idealization consistent with our approximation of one same maximum recirculation rate for all interglacial periods (a minimum recirculation time of 700 yr); it appears to work rather well but restricts us from examining any possible changes in the concentrations of dissolved inorganic matter in the deep ocean (Sigman and Haug, 2003).The incorporation of these two elements (Southern Ocean and particulate matter) go beyond the objectives of this initial study, i.e. having the simplest model inspired on physiologic analogies, and is left to future efforts.
Let us see what the best GA execution tells us about the temporal changes in aerobic and anaerobic DIC supply to the upper ocean.For this purpose we non-dimensionalise the net community production as P − R = Q 0 c d P − R , then the non-dimensional form of Eq. ( 7) is (after dropping primes): where we have used Eq. ( 12) for the non-dimensional flow rate Q.The first and second terms in the right-hand-side of Eq. ( 13) respectively are the DIC input from the deep ocean and the unused DIC fraction that is returned back to the deep ocean, the third term is the DIC supply resulting from remineralization minus production.
The upper panel in Fig. 9 illustrates the (1 + ε) input and −(1 + ε) c output terms and the lower panel in Fig. 9 shows their difference, or net input (1 + ε) (1 − c).The lower panel in Fig. 9  The zoom at the glacial-interglacial transition in Fig. 10 is clarifying.Near steady state P is always larger than R but the difference is substantially larger by the end of the interglacial period.This grand excess of production over respiration is responsible for the storage of organic matter in the upper-ocean, these are the reserves that will be necessary to sustain the next glacial-interglacial transition.Clearly this result depends on Γ(t) but, because of the limited DOC reserves in the upper ocean, we expect this parameter to be less than the non-dimensional enhanced-concentration level γ 2 only at the beginning of the interglacial period, as both the interglacial and glacial steady states have to be sustained by the deep-ocean supply.
We view this anaerobic contribution, the M ana term discussed in Sect.2.2, as arising mainly from microbial remineralization of DOC reserves.Del Giorgio and Duarte (2002) have estimated that nowadays respiration exceeds production by 6 Pg yr −1 , so the available dissolved organic carbon in the upper ocean (several hundred Pg according to Sundquist andVisser, 2003 andHansell et al., 2009) could sustain today's remineralization rate up to about one hundred years.The above results, however, suggest these reserves might have been much greater by the end of a glacial period, i.e. they would have been capable of sustaining enhanced respiration much longer, perhaps during a few thousand years at the beginning of the glacial-interglacial transition.

Concluding remarks
A comparison of the temporal changes in the distribution of energy for the earth and mammals, two complex systems, shows inspiring analogies.For this comparison, however, we must keep in mind that the energy supply for mammals, as heterotrophic systems, relies on oxidation of organic matter; hence, the faster oxygen and organic nutrients are distributed, the faster will oxidation and sustained energy supply take place.Instead, the earth as a whole is an autotrophic system, as it relies on the photosynthetic transformation of the sun's radiative energy, so it needs inorganic carbon and nutrients; hence, the faster carbon and nutrients are distributed, the faster solar Introduction

Conclusions References
Tables Figures

Back Close
Full energy will be incorporated to the system.The first pleasant surprise is the existence of relatively short and long distribution circuits in both systems, with analogous gasexchange and nutrient-supply roles (Pelegrí, 2008).But the physiologic analogy goes further beyond.As a mammal switches from one metabolic state to another, from rest to exercise or vice versa, the energy expenditure suddenly changes one order of magnitude.During a sudden low to high switch, the circulatory distribution system tries to match the increased demand but it cannot do so instantly, the distribution system requires some finite time to increase its flow rate in order to advect enough oxygen and nutrients to sustain the enhanced requirements.Hence, there is a period when it heavily relies on fast-access anaerobic supply, i.e. nearby stored energy.During a sudden high to low change, the opposite situation takes place as advective energy flow decreases much more slowly than the energy demand and there is a surplus that goes into replenishing the exhausted proximal storages.
We propose the earth, or ocean, system follows an analogous temporal pattern.As the global thermodynamic state of the system changes, the deep ocean circulation attempts to match the required "aerobic" supply (in physiological terms, in the sense that it is provided by enhanced advection, but actually transporting DIC and inorganic nutrients).This match is not attained instantly so the system must resource on "anaerobic" storage (again in physiological terms, as now we actually mean inorganic nutrients and carbon resulting from oxidation of organic matter).Pelegrí (2008) proposed that this fast source comes from increased respiration of organic matter, possibly induced by an initial increase in temperature and microbial respiration (Rivkin and Legendre, 2001;Brown et al., 2004).Such respiration would be capable of providing energy during sufficiently long times as DOC, from labile to more ancient forms, progressively get remineralized.
Another key concept that arises from the physiological comparison is what we mean by the earth's metabolic rate.For living beings the metabolic rate is the energy per unit time necessary to sustain all the system activities.For mammals at steady state the Figures

Back Close
Full metabolic rate equals the energy supply through biomass oxidation.We have shown the earth's metabolic rate at any time is net autotrophic-community production, which at steady state becomes net community production.In this sense we may think of glacial periods as times of basal (relatively low) metabolic requirements and of interglacial periods as times of enhanced (relatively high) energy demand, periods of rest and exercise of the earth system.These ideas are coherent with observations suggesting that, during the termination of a glaciation, temperature actually leads carbon dioxide by a time difference of order 1 kyr (Fischer et al., 1999;Shackleton, 2000;Monin et al., 2001;Caillon et al., 2003;Siegenthaler et al., 2005).If the change in carbon dioxide was due to a simple temperature-related decrease in its water solubility then the two time series should change in parallel.Alternatively, an increase in the thermohaline circulation would lead to cooler waters reaching the ocean surface and, hence, to the progressive incorporation of more radiative energy.This would enhance the rates of remineralisation, increase DIC and outgassing, followed by an even faster recirculation of (DIC-rich) deep waters to the upper ocean and further outgassing.
But what drives the state change?Here enters again a physiological concept.In a mammal the heart is the centre of the energy distribution system, and the heart rate is imposed by pacemaker cells which require very minor energy input to change their pace, many orders of magnitude less than the energy involved in the metabolic response.The year-integrated amount of solar energy arriving to the earth has remained nearly constant for millions of years so the way this energy is incorporated and distributed into the system, through the thermohaline circulation, is the principal responsible for the observed global changes in the earth's state.Both mammals and the earth undergo state transitions which do not imply additional external energy, the source remains constant but the amount of energy actually supplied responds to the changes in the energy distribution patterns.
We have used these ideas to set up a very simple two-box two-state metabolic model that simulates the glacial-interglacial transitions of the earth system.Box models have Introduction

Conclusions References
Tables Figures

Back Close
Full been repeatedly used with substantial success in ocean sciences but, to our knowledge, they have never been used to simulate glacial-interglacial transitions.During such a transition the external energy supply does not change, rather there is a swifter or slower distribution of the system's own resources.The system switches between low-energy supply, the basal or glacial state, and high-energy supply, the enhanced or interglacial state.The metabolic rate closely responds to changes in energy supply, from low to high levels of energy transformation during glacial and interglacial periods, respectively.
The model, despite its extreme simplicity and limitations, reproduces well the glacialinterglacial changes.It indeed misses the abrupt changes but provides information on several glacial-interglacial key parameters: the times of change between the basal and enhanced states, the relative energetic level of each interglacial period, and the location and intensity of the maximum interglacial recirculation rates (Table 1).In particular, the model predicts a basal recirculation period that ranges between about 40 and 45 kyr (Table 1), similar to the period of the earth's obliquity cycle.The model confirms that during glacial periods production marginally exceeds respiration while at the peak of the interglacial state production greatly exceeds respiration.Only at the start of an interglacial cycle respiration exceeds production, so this regenerated inorganic supply is capable of sustaining the early stages of the enhanced metabolic state.
The results depend largely on the size and sign of the thermohaline circulation, the largest the circulation rate the fastest deep-DIC reaches the upper ocean.A rapid increase in the circulation rate from glacial to interglacial states results in an effective high-energy supply, on the contrary a rapid decrease between interglacial and glacial states results in an effective low supply.The larger the rate of change the greater is the asymmetry between the glacial-interglacial and interglacial-glacial transitions, in particular with long interglacial-glacial transitions.A similar asymmetry is observed in the aerobic supply of energy to mammals: energy is efficiently extracted to sustain fast transitions from basal to enhanced states but is not as efficiently stored during the opposite transition (Fig. 1).

ESDD Introduction Conclusions References
Tables Figures

Back Close
Full Our work suggests a number of lines for future research.One principal path is to further explore the behaviour of organic matter, not only through Eq. ( 8) but also including additional equations for dissolved and particulate materials.Another important course would be the inclusion of a third ocean compartment, for nutrient-rich highlatitude waters.We trust that future work will continue to benefit from observations of the organization and response of other complex systems, hence using physiology as a source of inspiration.Introduction

Conclusions References
Tables Figures

Back Close
Full  Full tectures, capable of optimizing the flow of energy.The existence of similar temporal Figures Back Close Full Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | ) with the net energy supply being c d Q s − M s .The two opposite situations are (i) c s = 0 when all DIC flux is exhausted and the metabolic rate equals the rate of energy supply, M s = Q s c d , and (ii) c s = c d if no energy is transformed, M s = 0.For the ocean, the DIC c s is not substantially different from c d but, as explained above, its transformation is limited by the concentration of inorganic nutrients, which indeed may approach nearzero values.During a state transition the DIC net supply has to suddenly change to meet the new requirements.As in a mammal, the deep ocean DIC extraction rate (c d − c)Q changes smoothly, through both Q and (c d − c), to meet this new energetic requirement.However, this advective extraction takes some time to reach the required level, meanwhile the system either uses energy reserves or builds up reserves, we call this transfer Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper |

Figure 7
Figure7presents the model results for the best GA execution, where the DIC model output reaches a correlation of 0.813 with the atmospheric CO 2 time-series.The model is capable of reproducing very well the overall features of all four interglacial-glacial cycles.In particular, it reproduces correctly the asymmetry between the rapid glacialinterglacial transition and the much slower interglacial-glacial decay.This is possible Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | also shows the rate of change of DIC directly computed from the Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper |

Fig. 1 .Fig. 2 .Fig. 5 .Fig. 8 .
Fig. 1.Schematics of energy supply in mammals during state transitions, a similar pattern is assumed to occur in the Earth system.Energy supply is given by S= Q c d + R − P = Q c d − M which is parameterized as S = Q g.In this plot the supply changes between a basal level Q 1 g 1 and an enhanced level Q 2 g 2 .The aerobic supply is given by Q c d and the anaerobic processes R − P account for the deficit or excess supply.

Table 1 .
The top three files show the parameters obtained from the analysis of the CO 2 time series, rows 4 to 6 show those model parameters corresponding to the best data fit as obtained after 400 GA executions (correlation 0.813), and rows 7 to 9 show the statistics from those executions (38) whose correlation with the observations is above 0.81..977± 0.028 0.923 ± 0.021 0.993 ± 0.012 1.004 ± 0.002 executions τ 0 (kyr BP) 40.4 ± 0.9 43.7 ± 0.6 44.8 ± 0.2 42.0 ± 1.2 Introduction 0