From 8bb459daec560d4748fc46676d9482d1d5243dac Mon Sep 17 00:00:00 2001
From: KyleKlenk <kyle.c.klenk@gmail.com>
Date: Wed, 9 Aug 2023 14:41:19 -0600
Subject: [PATCH] Added HRU rel and abs tolerances to be set by file

---
 .../file_access_actor/file_access_actor.hpp   |  2 +
 build/includes/global/global.hpp              |  7 ++-
 build/includes/global/settings_functions.hpp  |  8 ++-
 .../hru_actor_subroutine_wrappers.hpp         |  1 +
 build/includes/job_actor/GRU.hpp              | 12 +++-
 .../cpp_code/file_access_actor.cpp            | 18 +++---
 .../fortran_code/write_to_netcdf.f90          |  7 ++-
 .../actors/global/settings_functions.cpp      |  6 ++
 .../actors/hru_actor/cpp_code/hru_actor.cpp   | 62 +++++++++++--------
 .../hru_actor/fortran_code/hru_actor.f90      | 23 +++++++
 build/source/actors/job_actor/GRU.cpp         | 20 +++++-
 build/source/actors/job_actor/job_actor.cpp   | 16 +++--
 build/source/dshare/data_types.f90            |  4 ++
 build/source/netcdf/def_output.f90            |  8 +--
 14 files changed, 145 insertions(+), 49 deletions(-)

diff --git a/build/includes/file_access_actor/file_access_actor.hpp b/build/includes/file_access_actor/file_access_actor.hpp
index 7bffc19..500ef4b 100644
--- a/build/includes/file_access_actor/file_access_actor.hpp
+++ b/build/includes/file_access_actor/file_access_actor.hpp
@@ -21,6 +21,8 @@ struct netcdf_gru_actor_info {
     
     int state_var_id; // The success of the GRU 1 = pass, 0 = fail
     int num_attempts_var_id;
+    int rel_tol_var_id;
+    int abs_tol_var_id;
 };
 
 
diff --git a/build/includes/global/global.hpp b/build/includes/global/global.hpp
index 4730894..c5ac1a0 100644
--- a/build/includes/global/global.hpp
+++ b/build/includes/global/global.hpp
@@ -17,6 +17,9 @@ struct serializable_netcdf_gru_actor_info {
 
     int successful; // 0 = false, 1 = true
     int num_attempts;
+
+    double rel_tol;
+    double abs_tol;
 };
 
 template<class Inspector>
@@ -27,7 +30,9 @@ bool inspect(Inspector& f, serializable_netcdf_gru_actor_info& x) {
                               f.field("run_physics_duration", x.run_physics_duration),
                               f.field("write_output_duration", x.write_output_duration),
                               f.field("successful", x.successful),
-                              f.field("num_attempts", x.num_attempts));
+                              f.field("num_attempts", x.num_attempts),
+                              f.field("rel_tol", x.rel_tol),
+                              f.field("abs_tol", x.abs_tol));
 }
 
 
diff --git a/build/includes/global/settings_functions.hpp b/build/includes/global/settings_functions.hpp
index 706b7cf..f50a628 100644
--- a/build/includes/global/settings_functions.hpp
+++ b/build/includes/global/settings_functions.hpp
@@ -96,13 +96,15 @@ bool inspect(Inspector& inspector, Job_Actor_Settings& job_actor_settings) {
 Job_Actor_Settings readJobActorSettings(std::string json_settings_file);
 
 // ####################################################################
-//                          SUMMA Actor Settings
+//                          HRU Actor Settings
 // ####################################################################
 
 struct HRU_Actor_Settings {
     bool print_output;
     int output_frequency;
     int dt_init_factor; // factor to multiply the initial timestep by
+    double rel_tol;
+    double abs_tol;
 };
 
 template<class Inspector>
@@ -110,7 +112,9 @@ bool inspect(Inspector& inspector, HRU_Actor_Settings& hru_actor_settings) {
     return inspector.object(hru_actor_settings).fields(
                 inspector.field("print_output",     hru_actor_settings.print_output),
                 inspector.field("output_frequency", hru_actor_settings.output_frequency),
-                inspector.field("dt_init_factor",   hru_actor_settings.dt_init_factor));
+                inspector.field("dt_init_factor",   hru_actor_settings.dt_init_factor),
+                inspector.field("rel_tol",          hru_actor_settings.rel_tol),
+                inspector.field("abs_tol",          hru_actor_settings.abs_tol));
 }
 
 HRU_Actor_Settings readHRUActorSettings(std::string json_settings_file);
diff --git a/build/includes/hru_actor/hru_actor_subroutine_wrappers.hpp b/build/includes/hru_actor/hru_actor_subroutine_wrappers.hpp
index df1dab9..c9fff76 100644
--- a/build/includes/hru_actor/hru_actor_subroutine_wrappers.hpp
+++ b/build/includes/hru_actor/hru_actor_subroutine_wrappers.hpp
@@ -88,5 +88,6 @@ extern "C" {
       void computeTimeForcingHRU(void* handle_timeStruct, void* handle_forcStruct, double* fracJulDay, 
                               int* yearLength, int* err);
 
+      void setIDATolerances(void* handle_mparStruct, double* relTol, double* absTol);
 
 }
\ No newline at end of file
diff --git a/build/includes/job_actor/GRU.hpp b/build/includes/job_actor/GRU.hpp
index 8123f21..4d57edf 100644
--- a/build/includes/job_actor/GRU.hpp
+++ b/build/includes/job_actor/GRU.hpp
@@ -29,6 +29,8 @@ class GRU {
 
     // Modifyable Parameters
     int dt_init_factor; // The initial dt for the GRU
+    double rel_tol; // The relative tolerance for the GRU
+    double abs_tol; // The absolute tolerance for the GRU
 
     // Status Information
     int attempts_left; // The number of attempts left for the GRU to succeed
@@ -44,7 +46,8 @@ class GRU {
     
   public:
     // Constructor
-    GRU(int global_gru_index, int local_gru_index, caf::actor gru_actor, int dt_init_factor, int max_attempts);
+    GRU(int global_gru_index, int local_gru_index, caf::actor gru_actor, int dt_init_factor, 
+        double rel_tol, double abs_tol, int max_attempts);
 
     // Deconstructor
     ~GRU();
@@ -60,12 +63,16 @@ class GRU {
     double getRunPhysicsDuration();
     double getWriteOutputDuration();
 
+    double getRelTol();
+    double getAbsTol();
+
     double getAttemptsLeft();
     gru_state getStatus();
 
     bool isFailed();
 
 
+
     // Setters
     void setRunTime(double run_time);
     void setInitDuration(double init_duration);
@@ -73,6 +80,9 @@ class GRU {
     void setRunPhysicsDuration(double run_physics_duration);
     void setWriteOutputDuration(double write_output_duration);
 
+    void setRelTol(double rel_tol);
+    void setAbsTol(double abs_tol);
+
     void setSuccess();
     void setFailed();
     void setRunning();
diff --git a/build/source/actors/file_access_actor/cpp_code/file_access_actor.cpp b/build/source/actors/file_access_actor/cpp_code/file_access_actor.cpp
index 32afa0c..02d2a32 100644
--- a/build/source/actors/file_access_actor/cpp_code/file_access_actor.cpp
+++ b/build/source/actors/file_access_actor/cpp_code/file_access_actor.cpp
@@ -57,10 +57,11 @@ behavior file_access_actor(stateful_actor<file_access_state>* self, int start_gr
             set_var_d(bpar_struct, params->handle_bpar_struct);
             // write the populated data to netCDF
             writeParamToNetCDF(self->state.handle_ncid, &index_gru, &index_hru, 
-                params->handle_attr_struct, 
-                params->handle_type_struct, 
-                params->handle_mpar_struct, 
-                params->handle_bpar_struct, &err);
+                               params->handle_attr_struct, 
+                               params->handle_type_struct, 
+                               params->handle_mpar_struct, 
+                               params->handle_bpar_struct, 
+                               &err);
         
 
             self->state.file_access_timing.updateEndPoint("write_duration");
@@ -224,10 +225,13 @@ behavior file_access_actor(stateful_actor<file_access_state>* self, int start_gr
         },
 
 
-        [=](deallocate_structures, std::vector<serializable_netcdf_gru_actor_info> &netcdf_gru_info) {
+        [=](finalize, std::vector<serializable_netcdf_gru_actor_info> &netcdf_gru_info) {
             int num_gru = netcdf_gru_info.size();
-            WriteGRUStatistics(self->state.handle_ncid, &self->state.gru_actor_stats, 
-                    netcdf_gru_info.data(), &num_gru, &self->state.err);
+            WriteGRUStatistics(self->state.handle_ncid, 
+                               &self->state.gru_actor_stats, 
+                               netcdf_gru_info.data(), 
+                               &num_gru, 
+                               &self->state.err);
 
             
             // call output_container deconstructor
diff --git a/build/source/actors/file_access_actor/fortran_code/write_to_netcdf.f90 b/build/source/actors/file_access_actor/fortran_code/write_to_netcdf.f90
index dc51f1f..91bbaaf 100644
--- a/build/source/actors/file_access_actor/fortran_code/write_to_netcdf.f90
+++ b/build/source/actors/file_access_actor/fortran_code/write_to_netcdf.f90
@@ -305,7 +305,8 @@ subroutine writeGRUStatistics(handle_ncid,      &
   real(c_double), dimension(num_gru)                  :: forcing_time_array
   real(c_double), dimension(num_gru)                  :: run_physics_time_array
   real(c_double), dimension(num_gru)                  :: write_output_time_array
-
+  real(c_double), dimension(num_gru)                  :: rel_tol_array
+  real(c_double), dimension(num_gru)                  :: abs_tol_array
   integer(c_int), dimension(num_gru)                  :: successful_array
   integer(c_int), dimension(num_gru)                  :: num_attempts_array
 
@@ -322,6 +323,8 @@ subroutine writeGRUStatistics(handle_ncid,      &
     forcing_time_array(i) = gru_stats_vector(i)%forcing_duration
     run_physics_time_array(i) = gru_stats_vector(i)%run_physics_duration
     write_output_time_array(i) = gru_stats_vector(i)%write_output_duration
+    rel_tol_array(i) = gru_stats_vector(i)%rel_tol
+    abs_tol_array(i) = gru_stats_vector(i)%abs_tol
     successful_array(i) = gru_stats_vector(i)%successful
     num_attempts_array(i) = gru_stats_vector(i)%num_attempts
   end do
@@ -335,6 +338,8 @@ subroutine writeGRUStatistics(handle_ncid,      &
     err = nf90_put_var(ncid%var(iFreq), gru_var_ids%write_output_duration_var_id, write_output_time_array)
     err = nf90_put_var(ncid%var(iFreq), gru_var_ids%state_var_id, successful_array)
     err = nf90_put_var(ncid%var(iFreq), gru_var_ids%num_attempts_var_id, num_attempts_array)
+    err = nf90_put_var(ncid%var(iFreq), gru_var_ids%rel_tol_var_id, rel_tol_array)
+    err = nf90_put_var(ncid%var(iFreq), gru_var_ids%abs_tol_var_id, abs_tol_array)
   end do
 
 end subroutine writeGRUStatistics
diff --git a/build/source/actors/global/settings_functions.cpp b/build/source/actors/global/settings_functions.cpp
index ad69d90..4ec4c02 100644
--- a/build/source/actors/global/settings_functions.cpp
+++ b/build/source/actors/global/settings_functions.cpp
@@ -116,6 +116,12 @@ HRU_Actor_Settings readHRUActorSettings(std::string json_settings_file) {
     hru_actor_settings.dt_init_factor = getSettings(json_settings_file, parent_key,
         "dt_init_factor", hru_actor_settings.dt_init_factor).value_or(1);
 
+    hru_actor_settings.rel_tol = getSettings(json_settings_file, parent_key,
+        "rel_tol", hru_actor_settings.rel_tol).value_or(1e-6);
+
+    hru_actor_settings.abs_tol = getSettings(json_settings_file, parent_key,
+        "abs_tol", hru_actor_settings.abs_tol).value_or(1e-6);
+
     return hru_actor_settings;
 }
 
diff --git a/build/source/actors/hru_actor/cpp_code/hru_actor.cpp b/build/source/actors/hru_actor/cpp_code/hru_actor.cpp
index 6d6f40a..383cf4c 100644
--- a/build/source/actors/hru_actor/cpp_code/hru_actor.cpp
+++ b/build/source/actors/hru_actor/cpp_code/hru_actor.cpp
@@ -10,7 +10,9 @@
 namespace caf {
 
 behavior hru_actor(stateful_actor<hru_state>* self, int refGRU, int indxGRU,
-    HRU_Actor_Settings hru_actor_settings, caf::actor file_access_actor, caf::actor parent) {
+                   HRU_Actor_Settings hru_actor_settings, 
+                   caf::actor file_access_actor, 
+                   caf::actor parent) {
     
     // Actor References
     self->state.file_access_actor = file_access_actor;
@@ -241,40 +243,50 @@ behavior hru_actor(stateful_actor<hru_state>* self, int refGRU, int indxGRU,
 void Initialize_HRU(stateful_actor<hru_state>* self) {
 
     setupHRUParam(&self->state.indxHRU, 
-            &self->state.indxGRU,
-            self->state.handle_attrStruct, 
-            self->state.handle_typeStruct, 
-            self->state.handle_idStruct,
-            self->state.handle_mparStruct, 
-            self->state.handle_bparStruct, 
-            self->state.handle_bvarStruct,
-            self->state.handle_dparStruct, 
-            self->state.handle_lookupStruct,
-            self->state.handle_startTime, 
-            self->state.handle_oldTime,
-            &self->state.upArea, &self->state.err);
+                  &self->state.indxGRU,
+                  self->state.handle_attrStruct, 
+                  self->state.handle_typeStruct, 
+                  self->state.handle_idStruct,
+                  self->state.handle_mparStruct, 
+                  self->state.handle_bparStruct, 
+                  self->state.handle_bvarStruct,
+                  self->state.handle_dparStruct, 
+                  self->state.handle_lookupStruct,
+                  self->state.handle_startTime, 
+                  self->state.handle_oldTime,
+                  &self->state.upArea, 
+                  &self->state.err);
     if (self->state.err != 0) {
-        aout(self) << "Error: HRU_Actor - SetupHRUParam - HRU = " << self->state.indxHRU <<
-        " - indxGRU = " << self->state.indxGRU << " - refGRU = " << self->state.refGRU << std::endl;
+        aout(self) << "Error: HRU_Actor - SetupHRUParam - HRU = " << self->state.indxHRU
+                   << " - indxGRU = " << self->state.indxGRU 
+                   << " - refGRU = " << self->state.refGRU << "\n";
         self->quit();
         return;
     }
             
     summa_readRestart(&self->state.indxGRU, 
-            &self->state.indxHRU, 
-            self->state.handle_indxStruct, 
-            self->state.handle_mparStruct, 
-            self->state.handle_progStruct,
-            self->state.handle_diagStruct, 
-            self->state.handle_fluxStruct, 
-            self->state.handle_bvarStruct, 
-            &self->state.dt_init, &self->state.err);
+                      &self->state.indxHRU, 
+                      self->state.handle_indxStruct, 
+                      self->state.handle_mparStruct, 
+                      self->state.handle_progStruct,
+                      self->state.handle_diagStruct, 
+                      self->state.handle_fluxStruct, 
+                      self->state.handle_bvarStruct, 
+                      &self->state.dt_init, 
+                      &self->state.err);
     if (self->state.err != 0) {
-        aout(self) << "Error: HRU_Actor - summa_readRestart - HRU = " << self->state.indxHRU <<
-        " - indxGRU = " << self->state.indxGRU << " - refGRU = " << self->state.refGRU << std::endl;
+        aout(self) << "Error: HRU_Actor - summa_readRestart - HRU = " << self->state.indxHRU
+                   << " - indxGRU = " << self->state.indxGRU 
+                   << " - refGRU = " << self->state.refGRU << "\n";
         self->quit();
         return;
     }
+
+
+    // Set HRU Tolerances
+    setIDATolerances(self->state.handle_mparStruct, 
+                     &self->state.hru_actor_settings.rel_tol, 
+                     &self->state.hru_actor_settings.abs_tol);
             
 }
 
diff --git a/build/source/actors/hru_actor/fortran_code/hru_actor.f90 b/build/source/actors/hru_actor/fortran_code/hru_actor.f90
index 3b815dd..4714c63 100644
--- a/build/source/actors/hru_actor/fortran_code/hru_actor.f90
+++ b/build/source/actors/hru_actor/fortran_code/hru_actor.f90
@@ -14,6 +14,7 @@ public::getFirstTimestep
 public::setTimeZoneOffset
 public::prepareOutput
 public::updateCounters
+public::setIDATolerances
 
 real(dp),parameter  :: verySmall=1e-3_rkind      ! tiny number
 real(dp),parameter  :: smallOffset=1.e-8_rkind   ! small offset (units=days) to force ih=0 at the start of the day
@@ -549,4 +550,26 @@ subroutine updateCounters(handle_timeStruct, handle_statCounter, handle_outputTi
  elapsedWrite = elapsedWrite + elapsedSec(startWrite, endWrite)
 end subroutine updateCounters
 
+! Set the HRU's relative and absolute tolerances
+subroutine setIDATolerances(handle_mparStruct, rtol, atol) bind(C, name="setIDATolerances")
+  USE data_types,only:var_dlength
+  USE var_lookup,only:iLookPARAM
+
+  implicit none
+
+  type(c_ptr), intent(in), value           :: handle_mparStruct !  model parameters
+  real(c_double),intent(in)               :: rtol              ! relative tolerance
+  real(c_double),intent(in)               :: atol              ! absolute tolerance
+  ! local variables
+  type(var_dlength),pointer                :: mparStruct        ! model parameters
+
+  call c_f_pointer(handle_mparStruct, mparStruct)
+
+  mparStruct%var(iLookPARAM%relErrTol_ida)%dat(1) = rtol
+  mparStruct%var(iLookPARAM%absErrTol_ida)%dat(1) = atol
+
+
+
+end subroutine setIDATolerances
+
 end module hru_actor
\ No newline at end of file
diff --git a/build/source/actors/job_actor/GRU.cpp b/build/source/actors/job_actor/GRU.cpp
index 520c36d..9b43067 100644
--- a/build/source/actors/job_actor/GRU.cpp
+++ b/build/source/actors/job_actor/GRU.cpp
@@ -5,11 +5,14 @@
 
 
 
-GRU::GRU(int global_gru_index, int local_gru_index, caf::actor gru_actor, int dt_init_factor, int max_attempt) {
+GRU::GRU(int global_gru_index, int local_gru_index, caf::actor gru_actor, 
+         int dt_init_factor, double rel_tol, double abs_tol, int max_attempt) {
   this->global_gru_index = global_gru_index;
   this->local_gru_index = local_gru_index;
   this->gru_actor = gru_actor;
   this->dt_init_factor = dt_init_factor;
+  this->rel_tol = rel_tol;
+  this->abs_tol = abs_tol;
   this->attempts_left = max_attempt;
   this->state = gru_state::running;
 }
@@ -48,6 +51,14 @@ double GRU::getWriteOutputDuration() {
   return this->write_output_duration;
 }
 
+double GRU::getRelTol() {
+  return this->rel_tol;
+}
+
+double GRU::getAbsTol() {
+  return this->abs_tol;
+}
+
 double GRU::getAttemptsLeft() {
   return this->attempts_left;
 }
@@ -77,6 +88,13 @@ void GRU::setWriteOutputDuration(double write_output_duration) {
   this->write_output_duration = write_output_duration;
 }
 
+void GRU::setRelTol(double rel_tol) {
+  this->rel_tol = rel_tol;
+}
+void GRU::setAbsTol(double abs_tol) {
+  this->abs_tol = abs_tol;
+}
+
 void GRU::setSuccess() {
   this->state = gru_state::succeeded;
 }
diff --git a/build/source/actors/job_actor/job_actor.cpp b/build/source/actors/job_actor/job_actor.cpp
index fc86f75..0726291 100644
--- a/build/source/actors/job_actor/job_actor.cpp
+++ b/build/source/actors/job_actor/job_actor.cpp
@@ -112,6 +112,8 @@ behavior job_actor(stateful_actor<job_state>* self,
                                                      local_gru_index, 
                                                      gru, 
                                                      self->state.dt_init_start_factor, 
+                                                     self->state.hru_actor_settings.rel_tol,
+                                                      self->state.hru_actor_settings.abs_tol,
                                                      self->state.max_run_attempts));    
           }
         }, // end init_gru
@@ -169,11 +171,11 @@ behavior job_actor(stateful_actor<job_state>* self,
               auto global_gru_index = GRU->getGlobalGRUIndex();
               auto local_gru_index = GRU->getLocalGRUIndex();
               auto gru_actor = self->spawn(hru_actor, 
-                                            global_gru_index, 
-                                            local_gru_index, 
-                                            self->state.hru_actor_settings,
-                                            self->state.file_access_actor, 
-                                            self);
+                                           global_gru_index, 
+                                           local_gru_index, 
+                                           self->state.hru_actor_settings,
+                                           self->state.file_access_actor, 
+                                           self);
               self->state.gru_container.gru_list[local_gru_index-1]->setGRUActor(gru_actor);
             }
           }
@@ -186,7 +188,7 @@ behavior job_actor(stateful_actor<job_state>* self,
 
             self->request(self->state.file_access_actor, 
                           infinite,
-                          deallocate_structures_v, netcdf_gru_info)
+                          finalize_v, netcdf_gru_info)
               .await(
                 [=](std::tuple<double, double> read_write_duration) {
                 
@@ -258,6 +260,8 @@ std::vector<serializable_netcdf_gru_actor_info> getGruNetcdfInfo(int max_run_att
         
         gru_info.num_attempts = max_run_attempts - gru->getAttemptsLeft() + 1;
         gru_info.successful = success(gru->getStatus());
+        gru_info.rel_tol = gru->getRelTol();
+        gru_info.abs_tol = gru->getAbsTol();
 
         gru_netcdf_info.push_back(gru_info);
 
diff --git a/build/source/dshare/data_types.f90 b/build/source/dshare/data_types.f90
index fa0b480..a21593b 100755
--- a/build/source/dshare/data_types.f90
+++ b/build/source/dshare/data_types.f90
@@ -87,6 +87,8 @@ MODULE data_types
     integer(C_INT) :: write_output_duration_var_id
     integer(C_INT) :: state_var_id
     integer(C_INT) :: num_attempts_var_id
+    integer(C_INT) :: rel_tol_var_id
+    integer(C_INT) :: abs_tol_var_id
   end type netcdf_gru_actor_info
 
   type,public,bind(C) :: serializable_netcdf_gru_actor_info
@@ -97,6 +99,8 @@ MODULE data_types
     real(C_DOUBLE) :: write_output_duration
     integer(C_INT) :: successful
     integer(C_INT) :: num_attempts
+    real(C_DOUBLE) :: rel_tol
+    real(C_DOUBLE) :: abs_tol
   end type serializable_netcdf_gru_actor_info
 
   ! ***********************************************************************************************************
diff --git a/build/source/netcdf/def_output.f90 b/build/source/netcdf/def_output.f90
index 97f9328..0f3e4f7 100755
--- a/build/source/netcdf/def_output.f90
+++ b/build/source/netcdf/def_output.f90
@@ -233,11 +233,9 @@ subroutine def_output(handle_ncid,startGRU,nGRU,nHRU,actor_info,err) bind(C, nam
     err = nf90_def_var(ncid%var(iFreq),"write_output_duration",outputPrecision,(/gru_DimID/),actor_info%write_output_duration_var_id)
     err = nf90_def_var(ncid%var(iFreq),"successful",nf90_int,(/gru_DimID/),actor_info%state_var_id)
     err = nf90_def_var(ncid%var(iFreq),"num_attempts",nf90_int,(/gru_DimID/),actor_info%num_attempts_var_id)
-    if(err/=0) then 
-      message=trim(message)//trim(cmessage)
-      print*, message
-      return
-    end if
+    err = nf90_def_var(ncid%var(iFreq),"rel_tol",outputPrecision,(/gru_DimID/),actor_info%rel_tol_var_id)
+    err = nf90_def_var(ncid%var(iFreq),"abs_tol",outputPrecision,(/gru_DimID/),actor_info%abs_tol_var_id)
+    if(err/=0) then; message=trim(message)//trim(cmessage); print*, message; return; end if
 
   end do
 end subroutine def_output
-- 
GitLab