H2 effects on CH4 lifetime

DESCRIPTION

@@ -1,6 +1,6 @@
 Package: hector
 Title: The Hector Simple Climate Model
-Version: 3.2.0
+Version: 3.3.0.9999999
 Authors@R: c(person("Kalyn", "Dorheim", 
 	        email = "kalyn.dorheim@pnnl.gov",
 	        role = c("aut", "cre"),

NEWS.md

@@ -1,3 +1,7 @@
+# hector 3.3.0.9999999
+
+* under development for the H2 integration, increasing the version index though will be ensure that the right version of hector is being run with gcam
+
 # hector 3.2.0 
 [![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.10698028.svg)](https://doi.org/10.5281/zenodo.10698028)
 * Correct aerosol forcing coefficients based on Zelinka et al. (2023)

R/RcppExports.R

@@ -4,6 +4,9 @@
 #' @describeIn emissions Emissions hydrogen gas (\code{"Tg H2"})
 NULL
 
+#' @describeIn oh coefficent for the h2 emissions
+NULL
+
 #' @describeIn msgtype Message type for retrieving data from a component
 #' @keywords internal
 GETDATA <- function() {
@@ -1571,6 +1574,10 @@ HEAT_FLUX <- function() {
     .Call('_hector_HEAT_FLUX', PACKAGE = 'hector')
 }
 
+COEFF_H2 <- function() {
+    .Call('_hector_COEFF_H2', PACKAGE = 'hector')
+}
+
 newcore_impl <- function(inifile, loglevel, suppresslogging, name) {
     .Call('_hector_newcore_impl', PACKAGE = 'hector', inifile, loglevel, suppresslogging, name)
 }

inst/include/component_data.hpp

@@ -281,6 +281,7 @@
 #define D_COEFFICENT_CH4 "CCH4"
 #define D_COEFFICENT_NMVOC "CNMVOC"
 #define D_COEFFICENT_CO "CCO"
+#define D_COEFFICENT_H2 "CH2"
 #define D_EMISSIONS_H2 "H2_emissions"
 
 

inst/include/h_util.hpp

@@ -26,7 +26,7 @@
  * \brief The model version number to be included in logs and outputs.
  * \note  Manually update the git tag to match this.
  */
-#define MODEL_VERSION "3.2.0"
+#define MODEL_VERSION "3.3.0.9999999"
 
 #define OUTPUT_DIRECTORY "output/"
 

inst/include/oh_component.hpp

@@ -71,6 +71,8 @@ class OHComponent : public IModelComponent {
   double CNMVOC; // coefficent for NMVOC
   double CNOX;   // coefficent for NOX
   double CCH4;   // coefficent for CH4
+  unitval CH2;   // coefficent for CH4
+
 
   // logger
   Logger logger;

inst/input/hector_ssp119.ini

@@ -138,6 +138,7 @@ CNOX=0.0042			; coefficent for NOX
 CCO=-0.000105		; coefficent for CO
 CNMVOC=-0.000315	; coefficent for NMVOC (non methane VOC)
 CCH4=-0.32			; coefficent for CH4
+CH2=-0.00044625 ; coefficent for H2
 
 ;------------------------------------------------------------------------
 [ozone]

inst/input/hector_ssp126.ini

@@ -137,6 +137,8 @@ CNOX=0.0042			; coefficent for NOX
 CCO=-0.000105		; coefficent for CO
 CNMVOC=-0.000315	; coefficent for NMVOC (non methane VOC)
 CCH4=-0.32			; coefficent for CH4
+CH2=-0.00044625 ; coefficent for H2
+
 
 ;------------------------------------------------------------------------
 [ozone]

inst/input/hector_ssp245.ini

@@ -138,6 +138,8 @@ CNOX=0.0042			; coefficent for NOX
 CCO=-0.000105		; coefficent for CO
 CNMVOC=-0.000315	; coefficent for NMVOC (non methane VOC)
 CCH4=-0.32			; coefficent for CH4
+CH2=-0.00044625 ; coefficent for H2
+
 
 ;------------------------------------------------------------------------
 [ozone]

inst/input/hector_ssp370.ini

@@ -138,6 +138,8 @@ CNOX=0.0042			; coefficent for NOX
 CCO=-0.000105		; coefficent for CO
 CNMVOC=-0.000315	; coefficent for NMVOC (non methane VOC)
 CCH4=-0.32			; coefficent for CH4
+CH2=-0.00044625 ; coefficent for H2
+
 
 ;------------------------------------------------------------------------
 [ozone]

inst/input/hector_ssp434.ini

@@ -138,6 +138,8 @@ CNOX=0.0042			; coefficent for NOX
 CCO=-0.000105		; coefficent for CO
 CNMVOC=-0.000315	; coefficent for NMVOC (non methane VOC)
 CCH4=-0.32			; coefficent for CH4
+CH2=-0.00044625 ; coefficent for H2
+
 
 ;------------------------------------------------------------------------
 [ozone]

inst/input/hector_ssp460.ini

@@ -138,6 +138,8 @@ CNOX=0.0042			; coefficent for NOX
 CCO=-0.000105		; coefficent for CO
 CNMVOC=-0.000315	; coefficent for NMVOC (non methane VOC)
 CCH4=-0.32			; coefficent for CH4
+CH2=-0.00044625 ; coefficent for H2
+
 
 ;------------------------------------------------------------------------
 [ozone]

inst/input/hector_ssp534-over.ini

@@ -138,6 +138,8 @@ CNOX=0.0042			; coefficent for NOX
 CCO=-0.000105		; coefficent for CO
 CNMVOC=-0.000315	; coefficent for NMVOC (non methane VOC)
 CCH4=-0.32			; coefficent for CH4
+CH2=-0.00044625 ; coefficent for H2
+
 
 ;------------------------------------------------------------------------
 [ozone]

inst/input/hector_ssp585.ini

@@ -137,6 +137,8 @@ CNOX=0.0042			; coefficent for NOX
 CCO=-0.000105		; coefficent for CO
 CNMVOC=-0.000315	; coefficent for NMVOC (non methane VOC)
 CCH4=-0.32			; coefficent for CH4
+CH2=-0.00044625 ; coefficent for H2
+
 
 ;------------------------------------------------------------------------
 [ozone]

src/RcppExports.cpp

@@ -2627,6 +2627,16 @@ BEGIN_RCPP
     return rcpp_result_gen;
 END_RCPP
 }
+// COEFF_H2
+String COEFF_H2();
+RcppExport SEXP _hector_COEFF_H2() {
+BEGIN_RCPP
+    Rcpp::RObject rcpp_result_gen;
+    Rcpp::RNGScope rcpp_rngScope_gen;
+    rcpp_result_gen = Rcpp::wrap(COEFF_H2());
+    return rcpp_result_gen;
+END_RCPP
+}
 // newcore_impl
 Environment newcore_impl(String inifile, int loglevel, bool suppresslogging, String name);
 RcppExport SEXP _hector_newcore_impl(SEXP inifileSEXP, SEXP loglevelSEXP, SEXP suppressloggingSEXP, SEXP nameSEXP) {
@@ -3035,6 +3045,7 @@ static const R_CallMethodDef CallEntries[] = {
     {"_hector_FLUX_MIXED", (DL_FUNC) &_hector_FLUX_MIXED, 0},
     {"_hector_FLUX_INTERIOR", (DL_FUNC) &_hector_FLUX_INTERIOR, 0},
     {"_hector_HEAT_FLUX", (DL_FUNC) &_hector_HEAT_FLUX, 0},
+    {"_hector_COEFF_H2", (DL_FUNC) &_hector_COEFF_H2, 0},
     {"_hector_newcore_impl", (DL_FUNC) &_hector_newcore_impl, 4},
     {"_hector_shutdown", (DL_FUNC) &_hector_shutdown, 1},
     {"_hector_reset", (DL_FUNC) &_hector_reset, 2},

src/ch4_component.cpp

@@ -165,7 +165,7 @@ void CH4Component::run(const double runToDate) {
         << std::endl;
 
 // Permafrost thaw produces CH4 emissions
-#define PG_C_TO_TG_CH4 (1000.0 * 16.04 / 12.01)
+    #define PG_C_TO_TG_CH4 (1000.0 * 16.04 / 12.01);
     const double rh_ch4 =
         core->sendMessage(M_GETDATA, D_RH_CH4).value(U_PGC_YR) * PG_C_TO_TG_CH4;
 

src/oh_component.cpp

@@ -57,12 +57,15 @@ void OHComponent::init(Core *coreptr) {
 
   // Inform core what data we can provide
   core->registerCapability(D_LIFETIME_OH, getComponentName());
+  core->registerCapability(D_COEFFICENT_H2, getComponentName());
   // Register inputs accepted.  Note that more than one component
   // can accept an input
   core->registerInput(D_EMISSIONS_CO, getComponentName());
   core->registerInput(D_EMISSIONS_NMVOC, getComponentName());
   core->registerInput(D_EMISSIONS_NOX, getComponentName());
   core->registerInput(D_EMISSIONS_H2, getComponentName());
+  core->registerInput(D_COEFFICENT_H2, getComponentName());
+
 
 }
 
@@ -100,8 +103,8 @@ void OHComponent::setData(const string &varName, const message_data &data) {
       H_ASSERT(data.date != Core::undefinedIndex(), "date required");
       NMVOC_emissions.set(data.date, data.getUnitval(U_TG_NMVOC));
     } else if (varName == D_EMISSIONS_H2) {
-        H_ASSERT(data.date != Core::undefinedIndex(), "date required");
-        H2_emissions.set(data.date, data.getUnitval(U_TG_H2));
+      H_ASSERT(data.date != Core::undefinedIndex(), "date required");
+      H2_emissions.set(data.date, data.getUnitval(U_TG_H2));
     } else if (varName == D_INITIAL_LIFETIME_OH) {
       H_ASSERT(data.date == Core::undefinedIndex(), "date not allowed");
       TOH0 = data.getUnitval(U_YRS);
@@ -117,6 +120,9 @@ void OHComponent::setData(const string &varName, const message_data &data) {
     } else if (varName == D_COEFFICENT_NOX) {
       H_ASSERT(data.date == Core::undefinedIndex(), "date not allowed");
       CNOX = data.getUnitval(U_UNDEFINED);
+    } else if (varName == D_COEFFICENT_H2) {
+      H_ASSERT(data.date == Core::undefinedIndex(), "date not allowed");
+      CH2 = data.getUnitval(U_UNDEFINED);
     } else {
       H_THROW("Unknown variable name while parsing " + getComponentName() +
               ": " + varName);
@@ -149,6 +155,8 @@ void OHComponent::run(const double runToDate) {
   unitval current_nox = NOX_emissions.get(runToDate);
   unitval current_co = CO_emissions.get(runToDate);
   unitval current_nmvoc = NMVOC_emissions.get(runToDate);
+  unitval current_h2 = H2_emissions.get(runToDate);
+
 
   // get this from CH4 component, this is last year's value
   const double previous_ch4 =
@@ -169,7 +177,13 @@ void OHComponent::run(const double runToDate) {
         CNMVOC *
         ((1.0 * +current_nmvoc) -
          NMVOC_emissions.get(NMVOC_emissions.firstdate()).value(U_TG_NMVOC));
-    toh = a + b + c + d;
+    const double f =
+         CH2 *
+          ((1.0 * +current_h2) -
+           H2_emissions.get(H2_emissions.firstdate()).value(U_TG_H2));
+
+
+    toh = a + b + c + d + f;
     H_LOG(logger, Logger::DEBUG)
         << "Year " << runToDate << " toh = " << toh << std::endl;
   }
@@ -201,9 +215,12 @@ unitval OHComponent::getData(const std::string &varName, const double date) {
     H_ASSERT(date != Core::undefinedIndex(),
              "Date required for NMVOC emissions");
     returnval = NMVOC_emissions.get(date);
+  } else if (varName == D_COEFFICENT_H2) {
+    H_ASSERT(date == Core::undefinedIndex(), "Date not allowed for H2 coefficent");
+    returnval = CH2 ;
   } else if (varName == D_EMISSIONS_H2) {
-      H_ASSERT(date != Core::undefinedIndex(), "Date required for H2 emissions");
-      returnval = H2_emissions.get(date);
+    H_ASSERT(date != Core::undefinedIndex(), "Date required for H2 emissions");
+    returnval = H2_emissions.get(date);
   } else {
     H_THROW("Caller is requesting unknown variable: " + varName);
   }

src/rcpp_constants.cpp

@@ -1503,3 +1503,9 @@ String FLUX_INTERIOR() { return D_FLUX_INTERIOR; }
 //' @export
 // [[Rcpp::export]]
 String HEAT_FLUX() { return D_HEAT_FLUX; }
+
+
+//' @describeIn oh coefficent for the h2 emissions
+ //' @export
+ // [[Rcpp::export]]
+ String COEFF_H2() { return D_COEFFICENT_H2; }

tests/testthat/test-h2.R

@@ -25,9 +25,9 @@ test_that("tau ", {
 
 test_that("H2 emissions ", {
 
-    # Check to make sure that can fetch and set the H2 emissions although
-    # at this point changing the emissions will have no impact on
-    # [ch4] dynamics but in the future it should...
+    # Check to make sure that can fetch and set the H2 emissions and
+    # that the effects of adding some H2 emissions change thigns
+    # now that we have the CH4 indrect effects implement.
 
     hc <- newcore(ini)
     run(hc)
@@ -45,12 +45,13 @@ test_that("H2 emissions ", {
 
     diff <- abs(out$value - out2$value)
 
-    # As of now the [CH4] and tau oh should not change
-    expect_equal(mean(diff[out$variable == CONCENTRATIONS_CH4()]), 0)
-    expect_equal(mean(diff[out$variable == LIFETIME_OH()]), 0)
+    # Now we do expect to see changes in CH4 concentrations in response
+    # to H2 emissions
+    expect_gt(mean(diff[out$variable == CONCENTRATIONS_CH4()]), 0)
+    expect_gt(mean(diff[out$variable == LIFETIME_OH()]), 0)
 
-    # But if we can change the H2 emissions we should see a difference in
-    # H2 emissions between the two new runs!
+    # Since we are working off of some default H2 emissions = 0, then
+    # changing the H2 emissions we should see a difference.
     expect_equal(mean(diff[out$variable == EMISSIONS_H2()]), new_val)
 
     })