From 3a870caa33276cdbf503eddb8518871c3b7fdbb9 Mon Sep 17 00:00:00 2001
From: KyleKlenk <kyle.c.klenk@gmail.com>
Date: Tue, 15 Aug 2023 10:06:04 -0600
Subject: [PATCH] Add ability to pass in state variable tolerances:

Can use json to pass state variable tolerances.
Can also just pass in atol and rtol to have all state variables set the same.
---
 build/includes/global/settings_functions.hpp  | 14 +++
 .../hru_actor_subroutine_wrappers.hpp         | 16 +++-
 .../actors/global/settings_functions.cpp      | 87 ++++++++++++++++++-
 .../actors/hru_actor/cpp_code/hru_actor.cpp   | 17 +++-
 .../hru_actor/fortran_code/hru_actor.f90      | 53 +++++++++--
 .../hru_actor/fortran_code/hru_setup.f90      |  1 -
 build/source/netcdf/read_icond.f90            |  2 +-
 7 files changed, 173 insertions(+), 17 deletions(-)

diff --git a/build/includes/global/settings_functions.hpp b/build/includes/global/settings_functions.hpp
index f50a628..1b9f05a 100644
--- a/build/includes/global/settings_functions.hpp
+++ b/build/includes/global/settings_functions.hpp
@@ -105,6 +105,20 @@ struct HRU_Actor_Settings {
     int dt_init_factor; // factor to multiply the initial timestep by
     double rel_tol;
     double abs_tol;
+    double relTolTempCas;
+    double absTolTempCas;
+    double relTolTempVeg;
+    double absTolTempVeg;
+    double relTolWatVeg;
+    double absTolWatVeg;
+    double relTolTempSoilSnow;
+    double absTolTempSoilSnow;
+    double relTolWatSnow;
+    double absTolWatSnow;
+    double relTolMatric;
+    double absTolMatric;
+    double relTolAquifr;
+    double absTolAquifr;
 };
 
 template<class Inspector>
diff --git a/build/includes/hru_actor/hru_actor_subroutine_wrappers.hpp b/build/includes/hru_actor/hru_actor_subroutine_wrappers.hpp
index c9fff76..cee5ba4 100644
--- a/build/includes/hru_actor/hru_actor_subroutine_wrappers.hpp
+++ b/build/includes/hru_actor/hru_actor_subroutine_wrappers.hpp
@@ -88,6 +88,20 @@ 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);
+      void setIDATolerances(void* handle_mparStruct, 
+                            double* relTolTempCas,
+                            double* absTolTempCas,
+                            double* relTolTempVeg,
+                            double* absTolTempVeg,
+                            double* relTolWatVeg,
+                            double* absTolWatVeg,
+                            double* relTolTempSoilSnow,
+                            double* absTolTempSoilSnow,
+                            double* relTolWatSnow,
+                            double* absTolWatSnow,
+                            double* relTolMatric,
+                            double* absTolMatric,
+                            double* relTolAquifr,
+                            double* absTolAquifr);
 
 }
\ No newline at end of file
diff --git a/build/source/actors/global/settings_functions.cpp b/build/source/actors/global/settings_functions.cpp
index 1cf22f4..e09fbb4 100644
--- a/build/source/actors/global/settings_functions.cpp
+++ b/build/source/actors/global/settings_functions.cpp
@@ -116,11 +116,76 @@ 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);
 
+
+    /*
+    Set Tolerances
+    ---------------
+    We can use rel_tol and abs_tol to set the tolerances for all the state variables.
+    If we set rel_tol and abs_tol in the config.json file then we just don't include 
+    the other tolerance values and they will be set to the value of rtol and atol.
+    */
     hru_actor_settings.rel_tol = getSettings(json_settings_file, parent_key,
-        "rel_tol", hru_actor_settings.rel_tol).value_or(1e-6);
+        "rel_tol", hru_actor_settings.rel_tol).value_or(-9999);
 
     hru_actor_settings.abs_tol = getSettings(json_settings_file, parent_key,
-        "abs_tol", hru_actor_settings.abs_tol).value_or(1e-6);
+        "abs_tol", hru_actor_settings.abs_tol).value_or(-9999);
+
+    double local_rtol;
+    double local_atol;
+
+    if (hru_actor_settings.rel_tol > 0) {
+        local_rtol = hru_actor_settings.rel_tol;
+    } else {
+        local_rtol = 1e-6;
+    }
+
+    if (hru_actor_settings.abs_tol > 0) {
+        local_atol = hru_actor_settings.abs_tol;
+    } else {
+        local_atol = 1e-6;
+    }
+
+    hru_actor_settings.relTolTempCas = getSettings(json_settings_file, parent_key,
+        "relTolTempCas", hru_actor_settings.relTolTempCas).value_or(local_rtol);
+
+    hru_actor_settings.absTolTempCas = getSettings(json_settings_file, parent_key,
+        "absTolTempCas", hru_actor_settings.absTolTempCas).value_or(local_atol);
+
+    hru_actor_settings.relTolTempVeg = getSettings(json_settings_file, parent_key,
+        "relTolTempVeg", hru_actor_settings.relTolTempVeg).value_or(local_rtol);
+
+    hru_actor_settings.absTolTempVeg = getSettings(json_settings_file, parent_key,
+        "absTolTempVeg", hru_actor_settings.absTolTempVeg).value_or(local_atol);
+
+    hru_actor_settings.relTolWatVeg = getSettings(json_settings_file, parent_key,
+        "relTolWatVeg", hru_actor_settings.relTolWatVeg).value_or(local_rtol);
+
+    hru_actor_settings.absTolWatVeg = getSettings(json_settings_file, parent_key,
+        "absTolWatVeg", hru_actor_settings.absTolWatVeg).value_or(local_atol);
+
+    hru_actor_settings.relTolTempSoilSnow = getSettings(json_settings_file, parent_key,
+        "relTolTempSoilSnow", hru_actor_settings.relTolTempSoilSnow).value_or(local_rtol);
+
+    hru_actor_settings.absTolTempSoilSnow = getSettings(json_settings_file, parent_key,
+        "absTolTempSoilSnow", hru_actor_settings.absTolTempSoilSnow).value_or(local_atol);
+
+    hru_actor_settings.relTolWatSnow = getSettings(json_settings_file, parent_key,
+        "relTolWatSnow", hru_actor_settings.relTolWatSnow).value_or(local_rtol);
+
+    hru_actor_settings.absTolWatSnow = getSettings(json_settings_file, parent_key,
+        "absTolWatSnow", hru_actor_settings.absTolWatSnow).value_or(local_atol);
+
+    hru_actor_settings.relTolMatric = getSettings(json_settings_file, parent_key,
+        "relTolMatric", hru_actor_settings.relTolMatric).value_or(local_rtol);
+
+    hru_actor_settings.absTolMatric = getSettings(json_settings_file, parent_key,
+        "absTolMatric", hru_actor_settings.absTolMatric).value_or(local_atol);
+
+    hru_actor_settings.relTolAquifr = getSettings(json_settings_file, parent_key,
+        "relTolAquifr", hru_actor_settings.relTolAquifr).value_or(local_rtol);
+
+    hru_actor_settings.absTolAquifr = getSettings(json_settings_file, parent_key,
+        "absTolAquifr", hru_actor_settings.absTolAquifr).value_or(local_atol);
 
     return hru_actor_settings;
 }
@@ -149,7 +214,21 @@ void check_settings_from_json(Distributed_Settings &distributed_settings,
               << "************ HRU_ACTOR_SETTINGS ************\n"
               << hru_actor_settings.print_output << "\n"
               << hru_actor_settings.output_frequency << "\n"
-              << "rel_tol: " << hru_actor_settings.rel_tol << "\n"
-              << "abs_tol: " << hru_actor_settings.abs_tol << "\n\n\n"; 
+              << "rel_tol: "            << hru_actor_settings.rel_tol << "\n"
+              << "abs_tol: "            << hru_actor_settings.abs_tol << "\n"
+              << "relTolTempCas: "      << hru_actor_settings.relTolTempCas << "\n"
+              << "absTolTempCas: "      << hru_actor_settings.absTolTempCas << "\n"
+              << "relTolTempVeg: "      << hru_actor_settings.relTolTempVeg << "\n"
+              << "absTolTempVeg: "      << hru_actor_settings.absTolTempVeg << "\n"
+              << "relTolWatVeg: "       << hru_actor_settings.relTolWatVeg << "\n"
+              << "absTolWatVeg: "       << hru_actor_settings.absTolWatVeg << "\n"
+              << "relTolTempSoilSnow: " << hru_actor_settings.relTolTempSoilSnow << "\n"
+              << "absTolTempSoilSnow: " << hru_actor_settings.absTolTempSoilSnow << "\n"
+              << "relTolWatSnow: "      << hru_actor_settings.relTolWatSnow << "\n"
+              << "absTolWatSnow: "      << hru_actor_settings.absTolWatSnow << "\n"
+              << "relTolMatric: "       << hru_actor_settings.relTolMatric << "\n"
+              << "absTolMatric: "       << hru_actor_settings.absTolMatric << "\n"
+              << "relTolAquifr: "       << hru_actor_settings.relTolAquifr << "\n"
+              << "absTolAquifr: "       << hru_actor_settings.absTolAquifr << "\n\n\n";
 
 }
\ No newline at end of file
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 82bcb4c..171f994 100644
--- a/build/source/actors/hru_actor/cpp_code/hru_actor.cpp
+++ b/build/source/actors/hru_actor/cpp_code/hru_actor.cpp
@@ -282,11 +282,22 @@ void Initialize_HRU(stateful_actor<hru_state>* self) {
         return;
     }
 
-
     // Set HRU Tolerances
     setIDATolerances(self->state.handle_mparStruct, 
-                     &self->state.hru_actor_settings.rel_tol, 
-                     &self->state.hru_actor_settings.abs_tol);
+                     &self->state.hru_actor_settings.relTolTempCas,
+                     &self->state.hru_actor_settings.absTolTempCas,
+                     &self->state.hru_actor_settings.relTolTempVeg,
+                     &self->state.hru_actor_settings.absTolTempVeg,
+                     &self->state.hru_actor_settings.relTolWatVeg,
+                     &self->state.hru_actor_settings.absTolWatVeg,
+                     &self->state.hru_actor_settings.relTolTempSoilSnow,
+                     &self->state.hru_actor_settings.absTolTempSoilSnow,
+                     &self->state.hru_actor_settings.relTolWatSnow,
+                     &self->state.hru_actor_settings.absTolWatSnow,
+                     &self->state.hru_actor_settings.relTolMatric,
+                     &self->state.hru_actor_settings.absTolMatric,
+                     &self->state.hru_actor_settings.relTolAquifr,
+                     &self->state.hru_actor_settings.absTolAquifr);
             
 }
 
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 4714c63..1ad1316 100644
--- a/build/source/actors/hru_actor/fortran_code/hru_actor.f90
+++ b/build/source/actors/hru_actor/fortran_code/hru_actor.f90
@@ -551,22 +551,61 @@ subroutine updateCounters(handle_timeStruct, handle_statCounter, handle_outputTi
 end subroutine updateCounters
 
 ! Set the HRU's relative and absolute tolerances
-subroutine setIDATolerances(handle_mparStruct, rtol, atol) bind(C, name="setIDATolerances")
+subroutine setIDATolerances(handle_mparStruct,  &
+                            relTolTempCas,      &
+                            absTolTempCas,      &
+                            relTolTempVeg,      &
+                            absTolTempVeg,      &
+                            relTolWatVeg,       &
+                            absTolWatVeg,       &
+                            relTolTempSoilSnow, &
+                            absTolTempSoilSnow, &
+                            relTolWatSnow,      &
+                            absTolWatSnow,      &
+                            relTolMatric,       &
+                            absTolMatric,       &
+                            relTolAquifr,       &
+                            absTolAquifr) 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
+  type(c_ptr), intent(in), value          :: handle_mparStruct !  model parameters
+  real(c_double),intent(in)               :: relTolTempCas
+  real(c_double),intent(in)               :: absTolTempCas
+  real(c_double),intent(in)               :: relTolTempVeg
+  real(c_double),intent(in)               :: absTolTempVeg
+  real(c_double),intent(in)               :: relTolWatVeg
+  real(c_double),intent(in)               :: absTolWatVeg
+  real(c_double),intent(in)               :: relTolTempSoilSnow
+  real(c_double),intent(in)               :: absTolTempSoilSnow
+  real(c_double),intent(in)               :: relTolWatSnow
+  real(c_double),intent(in)               :: absTolWatSnow
+  real(c_double),intent(in)               :: relTolMatric
+  real(c_double),intent(in)               :: absTolMatric
+  real(c_double),intent(in)               :: relTolAquifr
+  real(c_double),intent(in)               :: absTolAquifr
   ! local variables
-  type(var_dlength),pointer                :: mparStruct        ! model parameters
+  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
+  mparStruct%var(iLookPARAM%relTolTempCas)%dat(1)       = relTolTempCas 
+  mparStruct%var(iLookPARAM%absTolTempCas)%dat(1)       = absTolTempCas
+  mparStruct%var(iLookPARAM%relTolTempVeg)%dat(1)       = relTolTempVeg
+  mparStruct%var(iLookPARAM%absTolTempVeg)%dat(1)       = absTolTempVeg
+  mparStruct%var(iLookPARAM%relTolWatVeg)%dat(1)        = relTolWatVeg
+  mparStruct%var(iLookPARAM%absTolWatVeg)%dat(1)        = absTolWatVeg
+  mparStruct%var(iLookPARAM%relTolTempSoilSnow)%dat(1)  = relTolTempSoilSnow
+  mparStruct%var(iLookPARAM%absTolTempSoilSnow)%dat(1)  = absTolTempSoilSnow
+  mparStruct%var(iLookPARAM%relTolWatSnow)%dat(1)       = relTolWatSnow
+  mparStruct%var(iLookPARAM%absTolWatSnow)%dat(1)       = absTolWatSnow
+  mparStruct%var(iLookPARAM%relTolMatric)%dat(1)        = relTolMatric
+  mparStruct%var(iLookPARAM%absTolMatric)%dat(1)        = absTolMatric
+  mparStruct%var(iLookPARAM%relTolAquifr)%dat(1)        = relTolAquifr
+  mparStruct%var(iLookPARAM%absTolAquifr)%dat(1)        = absTolAquifr
+
 
 
 
diff --git a/build/source/actors/hru_actor/fortran_code/hru_setup.f90 b/build/source/actors/hru_actor/fortran_code/hru_setup.f90
index 60e5fc7..323c499 100644
--- a/build/source/actors/hru_actor/fortran_code/hru_setup.f90
+++ b/build/source/actors/hru_actor/fortran_code/hru_setup.f90
@@ -92,7 +92,6 @@ subroutine setupHRUParam(&
   ! subroutines and functions
   use time_utils_module,only:elapsedSec                       ! calculate the elapsed time
   USE mDecisions_module,only:mDecisions                       ! module to read model decisions
-  USE ffile_info_module,only:ffile_info                       ! module to read information on forcing datafile
   ! USE read_attrb_module,only:read_attrb               ! module to read local attributes
   USE paramCheck_module,only:paramCheck                       ! module to check consistency of model parameters
   USE pOverwrite_module,only:pOverwrite                       ! module to overwrite default parameter values with info from the Noah tables
diff --git a/build/source/netcdf/read_icond.f90 b/build/source/netcdf/read_icond.f90
index 676f20f..c95af1c 100644
--- a/build/source/netcdf/read_icond.f90
+++ b/build/source/netcdf/read_icond.f90
@@ -264,7 +264,7 @@ subroutine read_icond(&
     endif
 
     ! make sure canopy water is positive
-    if(progData%var(iLookPROG%scalarCanopyliq)%dat(1) < 0.0001_rkind)then
+    if(progData%var(iLookPROG%scalarCanopyWat)%dat(1) < 0.0001_rkind)then
       progData%var(iLookPROG%scalarCanopyliq)%dat(1) = 0.0001_rkind
     endif
 
-- 
GitLab