the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
All aboard! Earth system investigations with the CH2O-CHOO TRAIN v1.0
Tyler Kukla
Daniel Ibarra
Kimberly V. Lau
Jeremy K. C. Rugenstein
Abstract. Models of the carbon cycle and climate on geologic (>104 year) timescales have improved tremendously in the last 50 years due to parallel advances in our understanding of the Earth system and the increase in computing power to simulate its key processes. Despite these advances, balancing the Earth System's vast complexity with a model's computational expense is a primary challenge in model development. Running longer simulations spanning hundreds of thousands of years or more generally requires reducing the complexity of the modeled climate system. However, simpler model frameworks often leave out certain features of the climate system, such as radiative feedbacks, shifts in atmospheric circulation, and the expansion and decay of ice sheets, which can have profound effects on the long-term carbon cycle. Here, we present a model for climate and the long-term carbon cycle that captures many fundamental features of global climate while retaining the computational efficiency needed to simulate millions of years of time. The Carbon-H2O Coupled HydrOlOgical model with Terrestrial Runoff And INsolation, or CH2O-CHOO TRAIN, couples a one-dimensional (latitudinal) moist static energy balance model of climate with a model for rock weathering and the long-term carbon cycle. The key advantages of this framework are (1) it simulates fundamental climate forcings and feedbacks; (2) it accounts for geographic configuration; and (3) it is highly customizable, equipped to easily add features, change the strength of feedbacks, and prescribe conditions that are often hard-coded or emergent properties of more complex models, such as climate sensitivity and the strength of meridional heat transport. The CH2O-CHOO TRAIN is capable of running million-year-long simulations in about thirty minutes on a laptop PC. This paper outlines the model equations, presents a sensitivity analysis of the climate responses to varied climatic and carbon cycle perturbations, and discusses potential applications and next stops for the CH2O-CHOO TRAIN.
- Preprint
(3250 KB) -
Supplement
(1195 KB) - BibTeX
- EndNote
Tyler Kukla et al.
Status: final response (author comments only)
-
RC1: 'Referee comment on egusphere-2022-1000', Pierre Maffre, 13 Jan 2023
Kukla et al. present here a model of global geochemical cycles (though the manuscript focus on carbon cycle) and climate at multi-million years timescale. They justify their approach among others by the integration of a recent improvements in the representation of hydrological cycle in energy balance models (Siler, 2018, 2019), which is a novelty in Earth system models of that level of complexity. They provide an adequate discussion of the advantages and limitations of the model framework (though some further limitations can be touched upon). One important advantage being the compromise between representing some key features of the interactions between climate dynamics and continental fluxes, while being able to simulate long timescales needed to investigate geochemical cycles. This is particularly relevant if one consider further developments (or "stops") concerning elements with long residence time, such as oxygen and sulfur. The authors also provide an appreciable discussion of the model’s potential applications, and applications for which it is not optimized. The presented experiments give a good illustration these applications, and of the model’s behavior.
The manuscript suffers from several incorrectnesses and missing information that the authors need to address. Apart from these, I found the model well described and its results well explained.
Therefore, I recommend minor revision.Â
## Main comments
* There is one key missing element in the description: "rmax" (Eq. 12) is scaled to "keff" (line 280 of both Initialization/Code/CH2O-CHO_p_func_bistable-trackBC.R and Initialization/Code/CH2O-CHO_p_func_bistable-trackBC-INSOL.R), and not held constant, as stated by the authors (lines 207–208). This scaling was, indeed, already put forward by Maher & Chamberlain (2014).
It is critical because, without this scaling, the only effect of temperature is to reduce the weathering zone reactivity (fw in Maher & Chamberlain, Eq. 12). So weathering rate would be decreasing with temperature, whatever the value of pCO2 and runoff (if held constant). This would contradict all the discussion about weathering sensitivity to temperature and runoff.The authors need to update Eqs. 12 and 13 (and Table S2), also to avoid confusion between keff and reff.
* The ability and limitations of the model to simulate the hydrological cycle in the various experiments can be more discussed. Siler (2018) showed this current MEBM reproduces the link between hemispheric asymmetries in Qnet and the position of the ITCZ and the associated cross-equatorial heat transport. It will be meaningful to show how well (or not) the Northland experiments reproduce the southward shift of the ITCZ observed by of Laguë et al. (2021) in NorthlandBright configuration.
The authors indicate that the E-P formulation "captures the general trends on land". While this statement is supported with modern continental configuration, they omit to mention that the Budyko framework assumes a limitation of evaporation by available water on land (precipitation), but no limitation of precipitation by the distance to the ocean (continentality). In that regard, the Northland configuration, with almost an entire hemisphere continental, pushes the model towards the edge of its domain of validity. The hydrological cycle of Northland here probably resembles more the one of ThreePatchLand configuration from Laguë et al. (2021), where water bodies help carrying moisture up to the North pole.
It could also be reminded in the discussion that the hydrological cycle is, to some extent, constrained by the imposed profiles of w, u and R_G, meaning that some features of climate changes under different continental configurations, or pCO2, may not be seen. Is there a reference for the idealized profiles of R_G and u? It seems to be a determinant choice regarding the hydrological cycle.
Finally, a substantial element of this article is the suggestion that silicate weathering turns into a positive feedback if land is concentrated between 10° and 30° of latitude, because of the "arid" behaviour of drier conditions when global temperature increase. On modern continental configuration, most of India and South-East Asia, including the Himalaya, is between 10°N and 30°N. Yet, they don’t fall into an arid climate zone, and the reason for this is monsoon. One should expect that the concentration of land between 10° and 30° of latitude will generate a monsoonal circulation, that cannot be captured with this MEBM. The authors need to touch upon this potential challenge to their weathering feedback breakdown hypothesis.* The model description is missing one equation to describe Qnet and connect it to LWout (Eq. 9). Siler (2018) computed Qnet from modern climate reanalysis (or climate simulations), and used it directly as a forcing. However, it is not what the author are doing (in Initialization/Code/MEBM_ODEfun.R), they compute Qnet as "(1-albedo)I – LWout", "I" being the insolation, which is the forcing term.
Following that, I recommend to change the description "difference between top of atmosphere (TOA) and surface net downward energy fluxes" (line 115). Given the construction of the model, there is no net surface flux, so Qnet is simply the TOA net downward energy flux.
Finally, the author only mention insolation line 421, for the "reference" ("Cat-eye") simulation. I suppose that the insolation forcing the same for all the other experiments?* The authors never explain what is "net C emission". I suppose that it is the sum of all fluxes in Eq. 19, but it should be specified.
Following that remark, it is not easy to follow what the authors mean by "initial increase/decrease/drop of weathering" throughout the whole section 4. Lines 522–523 gives a good example: "Northland world cools due to an initial weathering increase". Yet, looking at Fig. 7H, what one can see is, at the time of the perturbation (increased D), a peak of positive net carbon emission (so towards the "less weathering" direction) for Northland. This contradiction is found in nearly all the experiments.
Assuming that "net C emission" is indeed the sum of all fluxes in Eq. 19, the "very initial" positive peak in net C emission is likely due to increased carbonate weathering bringing more C to the system (whereas silicate weathering has no instantaneous C effect strictly speaking). Only after the ocean chemistry respond to the perturbation can the increased alkalinity flux translate into negative C emission through carbonate burial. The author should indicate that the "very initial" instantaneous response is due to uncompensated changes in carbonate weathering, and is always neglected (justifiably so) in the interpretation. What the reader need to focus on is the "second initial" increase or decrease following the instantaneous response.* How is weathering of organic carbon computed? Lines 354–355 suggest that it is constant, but it could be more explicitly said. If it is indeed constant, while organic carbon burial is scaled to carbonate burial, I suspect that there is a net imbalance of organic C cycle at the end of all simulations: since carbonate weathering does not have the same sensitivity to climate than silicate weathering, Fw,sil getting back to its initial value (Fvolc) – as climate evolves towards its new equilibrium – does not necessarily mean that Fw,carb goes back to its initial value. Consequently, neither would Fb,carb and Fb,org .
Is that imbalance small enough to be neglected?* The model resolution (discretization) is not given.
* The units of land area in Fig. 4, discharge in Figs. 4B, S4C and G, and silicate weathering in Figs. 3C, 4C and S4D and H, are confusing. It is unclear whether the authors plotted the extensive properties (m2, m3/yr, mol/yr) for each latitudinal element (or grid cell) of the model, or whether they plotted their density functions per units of x (which happened to be unitless). In the first case, the values shown in the y axis would depend on the model’s discretization, and if the discretization was not regular in x, latitudinal variations on those figures could be simply due to grid cells having different area.
I recommend to specify on the figure "density of water land area/water discharge/silicate weathering per unit of x" (provided that it is indeed what they plotted), or to simply show the land fraction (continental fluxes are trickier because specific fluxes, like mol/m2/yr, hide the modulation by the land fraction).
Indicate also on those figures the latitudinal scale, which seems to be equal-area (x).
With those 2 conditions, the area under the curve will truly represent the Earth-integrated flux, regardless of the chosen discretization.
Finally, the units "% of global discharge" (or "% of global weathering" on Fig. S4D and G) suggests that the flux is normalized by the current Earth-integrated flux, and the curves should always sum to 100%. Yet, discharge on Fig 4B is virtually everywhere higher at 2360ppm than at 320ppm. Similarly, weathering seems everywhere lower (higher) on Fig S4D (H) if we compare the light green to the dark blue curve.Â
## Mathematical mistakes
* Eq. 7 is wrong, the rightmost part should be "-E0 /P + [1 + (E0 /P)^ω]^(1/ω)".
This mistake is not in the code (line 86 of Initialization/Code/MEBM_solve-bistable-glwx.R): the final "-1" of Eq. 7 multiplies the preceding part, instead of being added to.* There are some inconsistencies with Eq. 6 (Supplement Eq. 4).
It should be "cP /(Lv q)" (The code is OK: "cp/Lv/q", line 69 of Initialization/Code/MEBM_hydrofun.R)
As it is, evaporation is expressed in W/m2, which could be indicated, given that E, P and runoff are after discussed in m/yr.
Finally, Siler et al. (2019) derived this equation using q* (saturation specific humidity), which would be q/rh here. Is there a reason why the authors used q instead?* I don’t understand the justification for the "last term" (pCO2 - pCO2,0) of Eq. 15.
Eq. 15 without the last term (which is the original equation from Volk, 1987), can be rewritten in:
WZCO2 = pCO2 + RGPP*(WZ_CO2,0 – pCO2,0)
This already ensures that WZCO2 > pCO2 , unless the reference one (WZ_CO2,0) is not.
I don't expect, however, this modification to bring any significant change, since the pCO2_soil-pCO2 relationship looks pretty similar with or without the extra term.* Eq. 20: There's a sign mistake here, it should be "(δ13Cflx - δ13C)" in all terms ("flx" being "volc", "w,org", …). Same thing line 285, at least given how it is written in Eq. 20.
This mistake is not in the code (lines 596 and 802 of Initialization/Code/CH2O-CHO_p_func_bistable-trackBC.R and lines 604 and 804 of Initialization/Code/CH2O-CHO_p_func_bistable-trackBC-INSOL.R)Â
## Details
I recommend adding a table synthesizing the World configurations, instead of the information being dispersed in the whole section 3.
Moreover, there's no information about the position of "tropicslice" and "polarslice" (though the figures and text suggest that they are, respectively, centered on the equator and reaching the North pole).
There is no mention either of the latitude boundaries of "Patchland", the value of the proportion of land of "Patchland" and "Cat-eye" and whether this proportion of land is constant in "Patchland" (as it appears to be on Fig. 2).
Finally, the author never present the configuration "Polar[slice|hat] XL", from Supplementary Fig. S4."Polarslice" and "Tropicslice" have become "Polarhat" and "Tropicbelt" on Fig. 6, Fig. S3 and its caption
- Line 101: typo: "used TO calculate"
- Line 114: "F" is the northward flux, not the divergent flux. "dF/dx" is the divergence.
Also, it would be more accurate to describe F as the column-integrated AND zonally-integrated northward energy flux (which is consistent with its units being W).- Line 120: units of cp is "J/kg/K", not "J/kg"
- Line 121: It should be less confusing to introduce "q" in kg/kg, since it is multiplied by "Lv" in J/kg.
- Line 158: Don’t you mean "E is limited by P"?
- Line 178 (Eq. 9) and Supplementary Table S1: the meaning of "B" seems to be the Planck's negative feedback, that is, Earth emitting more longwave radiation as its temperature increase. It takes into account the water vapor positive feedback, as the value of B is lower than the pure blackbody emission, but it cannot be described as "Water vapor feedback coefficient".
- Lines 189–193: The authors could refer to the Supplement (section 2.3) here.
- Lines 212–213: This is confusing. "kreac" doesn’t seem to encapsulate the effects of mineral surface area and molar mass, since they are explicitly separated from this parameter. Besides, its units is mol/m2/y, not y^-1.
- Line 407: There's a misstatement here, Kump (2018) does not argue about a positive weathering feedback, they suggest an inefficient (but still negative) weathering feedback because of weathering being limited by the amount of exhumed material.
- Lines 441–448: Elevation is another key contributor too colder temperature in the South pole, because of lapse-rate and because it favors ice-sheet inception.
- Section 4.5: A word of caution should be added: changing D imply changing the Hadley cell and eddies-driven transport of moist static energy without changing their relative contribution (w). The choice of this sensitivity test therefore have consequences on the tropical/subtropical runoff response.
- Line 565: There is a simplification here: solute concentration [C] increases with temperature, but also decreases with runoff.
- Lines 710–720: Another element that could be mentioned here is that designing experiments with constant F_volc and adjusted W coefficient – rather than the opposite – is equivalent to keeping constant the residence time of C, and therefore, the equilibration time of C cycle.
- Fig. 1: Does "S" refer to sulfur, that is also simulated in the full geochemical model?
- Fig. 2: the "Northland" cartoon is confusing because, compared to "Tropicslice" and "Business belt" just next to it, it seems that lands start around 30°N, instead of 12°N stated in the text. I suggest to redraw that cartoon.
- Fig. 5: Though it is mentioned in the text, I suggest to specify the geography (continental configuration) on the figure, or its caption.
- Eq. 13: "keff" from previous equation and text has turned into "reff". Also, these equations (12–13) need to be update to describe the scaling of r_max (see my first comment).
- Eq. 16: While being exactly how Volk (1987) introduced the equation and its parameters, I found it oddly formulated, since if "pCO2" is replaced by "pCO2,half", GPP is not GPPmax/2 (though close to).
- Eq. 19: To be fully consistent, please add commas (",") in fluxes variables: Fw,org , Fw,carb , Fb,org and Fb,carb.
Â
#### Supplement
- Line 13: "ψ" is the southward transport (i.e., positive southward), not equatorward.
- Line 16: "1.5x10^4" seems to be a typo ("x" instead of "\times"?). Please add also its units (J/kg). Besides, based on Siler et al. (2018), g should here be 0.06*heq, to match Eq. 8 of Siler et al. (2018).
- Line 72: Missing reference to main text figure.
Â
#### Supplement tables (S1–S3)
- Title of Table S3 is incorrect. It should be something like "carbon cycle and other geochemical parameters", not "weathering model".
- "kreac" and "rmax" are a little confusing since they're both described as "reaction rate" but have different units, and do not represent the same thing (dissolution rate per unit of mineral surface VS dissolution rate per unit of water volume).
- What is "ε" in Table S3? How does it relate to "εb,org" and "εb,carb"? It seems that one value is missing to get the value of both "εb,org" and "εb,carb".
* Missing parameters
- "A" (specific mineral surface area, from Eq. 12)
- "m" (molar mass, from Eq. 12)
- "tz" (soil, or weathering zone, age) in equation 12 is named "Ts" in Table S2
- "T0" (reference temperature in Eq. 13). It is not mentioned in the main text either.
- "δ13Cw,org"
* wrong or missing units
- "cP" should be in J/kg/K, not J/kg
- "Rv" should be in J/kg/K, not J/kg
- "α" should be in K^-1, not unitless
- "CLW" should be in W/m2. not unitless
- "B" should be in W/m2/K, not unitless
- "M" should be in W/m2, not unitless
- "Lφ" should be in m, not unitless
- the units of "D" is spelled "sec" instead of "s"
Citation: https://doi.org/10.5194/egusphere-2022-1000-RC1 - RC2: 'Comment on egusphere-2022-1000', Anonymous Referee #2, 06 Mar 2023
Tyler Kukla et al.
Model code and software
CH2O-CHOO TRAIN v1.0 | Model code, instructions for running the model, and accessory scripts for plotting and generating new model input files Kukla, Ibarra, Lau, Rugenstein https://zenodo.org/record/7072803#.YzLgV6TMK3C
Tyler Kukla et al.
Viewed
HTML | XML | Total | Supplement | BibTeX | EndNote | |
---|---|---|---|---|---|---|
534 | 82 | 10 | 626 | 32 | 3 | 2 |
- HTML: 534
- PDF: 82
- XML: 10
- Total: 626
- Supplement: 32
- BibTeX: 3
- EndNote: 2
Viewed (geographical distribution)
Country | # | Views | % |
---|
Total: | 0 |
HTML: | 0 |
PDF: | 0 |
XML: | 0 |
- 1