feat(ngen): add mass balance check for all bmi modules during runtime

Author: hellkite500

CMakeLists.txt

@@ -321,6 +321,7 @@ add_subdirectory("src/forcing")
 add_subdirectory("src/utilities/mdarray")
 add_subdirectory("src/utilities/mdframe")
 add_subdirectory("src/utilities/logging")
+add_subdirectory("src/utilities/bmi")
 
 target_link_libraries(ngen
     PUBLIC
@@ -333,6 +334,7 @@ target_link_libraries(ngen
         NGen::forcing
         NGen::core_mediator
         NGen::logging
+        NGen::bmi_protocols
 )
 
 if(NGEN_WITH_SQLITE)

include/core/Layer.hpp

@@ -119,6 +119,8 @@ namespace ngen
                 double response(0.0);
                 try{
                     response = r_c->get_response(output_time_index, simulation_time.get_output_interval_seconds());
+                    // Check mass balance if able
+                    r_c->check_mass_balance(output_time_index, simulation_time.get_total_output_times(), current_timestamp);
                 }
                 catch(models::external::State_Exception& e){
                     std::string msg = e.what();
@@ -127,6 +129,13 @@ namespace ngen
                              +" at feature id "+id;
                     throw models::external::State_Exception(msg);
                 }
+                catch(std::exception& e){
+                    std::string msg = e.what();
+                    msg = msg+" at timestep "+std::to_string(output_time_index)
+                             +" ("+current_timestamp+")"
+                             +" at feature id "+id;
+                    throw std::runtime_error(msg);
+                }
                 std::string output = std::to_string(output_time_index)+","+current_timestamp+","+
                                     r_c->get_output_line_for_timestep(output_time_index)+"\n";
                 r_c->write_output(output);

include/realizations/catchment/Bmi_Module_Formulation.hpp

@@ -7,6 +7,7 @@
 #include "Bmi_Adapter.hpp"
 #include <DataProvider.hpp>
 #include "bmi_utilities.hpp"
+#include "bmi/protocols.hpp"
 
 using data_access::MEAN;
 using data_access::SUM;
@@ -250,6 +251,12 @@ namespace realization {
         const std::vector<std::string> get_bmi_input_variables() const override;
         const std::vector<std::string> get_bmi_output_variables() const override;
 
+        virtual void check_mass_balance(const int& iteration, const int& total_steps, const std::string& timestamp) const override {
+            //Create the protocol context, each member is const, and cannot change during the check
+            models::bmi::protocols::Context ctx{iteration, total_steps, timestamp, id};
+            bmi_protocols.mass_balance.run(ctx);
+        }
+
     protected:
 
         /**
@@ -419,6 +426,7 @@ namespace realization {
         int next_time_step_index = 0;
 
     private:
+        models::bmi::protocols::NgenBmiProtocols bmi_protocols;
         /**
          * Whether model ``Update`` calls are allowed and handled in some way by the backing model for time steps after
          * the model's ``end_time``.

include/realizations/catchment/Bmi_Multi_Formulation.hpp

@@ -42,6 +42,15 @@ namespace realization {
 
         virtual ~Bmi_Multi_Formulation() {};
 
+        virtual void check_mass_balance(const int& iteration, const int& total_steps, const std::string& timestamp) const final {
+            for( const auto &module : modules ) {
+                // TODO may need to check on outputs form each module indepdently???
+                // Right now, the assumption is that if each component is mass balanced
+                // then the entire formulation is mass balanced
+                module->check_mass_balance(iteration, total_steps, timestamp);
+            }
+        };
+
         /**
          * Convert a time value from the model to an epoch time in seconds.
          *

include/realizations/catchment/Formulation.hpp

@@ -44,6 +44,7 @@ namespace realization {
             virtual void create_formulation(boost::property_tree::ptree &config, geojson::PropertyMap *global = nullptr) = 0;
             virtual void create_formulation(geojson::PropertyMap properties) = 0;
 
+            virtual void check_mass_balance(const int& iteration, const int& total_steps, const std::string& timestamp) const = 0;
         protected:
 
             virtual const std::vector<std::string>& get_required_parameters() const = 0;

src/realizations/catchment/Bmi_Module_Formulation.cpp

@@ -382,6 +382,9 @@ namespace realization {
                     if (bmi_var_names_map.find(output_var_name) != bmi_var_names_map.end())
                         available_forcings.push_back(bmi_var_names_map[output_var_name]);
                 }
+                //Initialize all NgenBmiProtocols with the valid adapter pointer and any properties
+                //provided in the read configuration.
+                bmi_protocols = models::bmi::protocols::NgenBmiProtocols(get_bmi_model(), properties);
             }
         }
         /**

src/realizations/catchment/CMakeLists.txt

@@ -21,5 +21,6 @@ target_link_libraries(realizations_catchment PUBLIC
         NGen::geojson
         NGen::logging
         NGen::ngen_bmi
+        NGen::bmi_protocols
         )