updated config files for cfe2.0, check if nash storage is zero then s…

…kip the iteration.

configs/README.md

@@ -25,7 +25,7 @@ Example configuration files are provided in this directory. To build and run the
 | nash_storage | *double* |   |   | parameter_adjustable |  | Nash Config param - secondary reservoir   |
 | giuh_ordinates   | *double* |   |   | parameter_adjustable |  | Giuh ordinates in dt time steps   |
 | num_timesteps  | *int* |   |  | time_info |  | set to `1` if `forcing_file=BMI`   |
-| verbosity | *int* | `0`-`3`  |   | option |   |  prints various debug and bmi info  |
+| verbosity | *int* | `0`-`3`  |   | optional |   |  prints various debug and bmi info (defaults to 0) |
 | surface_water_partitioning_scheme | *char* | `Xinanjiang` or `Schaake`  |  | parameter_adjustable | infiltraton exces |    |
 | surface_runoff_scheme | *char* | GIUH or NASH_CASCADE | | parameter_adjustable | surface runoff | also supports 1 for GIUH and 2 for NASH_CASCADE; default is GIUH |
 | N_nash_surface | *int* |   |   | parameter_adjustable | surface runoff | number of Nash reservoirs for surface runoff   |

configs/cfe_config_cat_87.txt

@@ -17,16 +17,15 @@ soil_storage=0.585626[m/m]
 K_nash=0.03[]
 K_lf=0.01[]
 nash_storage=0.0,0.0
-giuh_ordinates=0.06,0.51,0.28,0.12,0.03
 num_timesteps=1
-verbosity=2
-surface_runoff_scheme=GIUH
+surface_runoff_scheme=NASH_CASCADE
+N_nash_surface=2
+K_nash_surface=0.83089
+nsubsteps_nash_surface=10
+nash_storage_surface=0.0,0.0
+#surface_partitioning_scheme=Schaake
 surface_water_partitioning_scheme=Xinanjiang
 a_Xinanjiang_inflection_point_parameter=1
 b_Xinanjiang_shape_parameter=1
 x_Xinanjiang_shape_parameter=1
 urban_decimal_fraction=0.0
-DEBUG=0
-#surface_water_partitioning_scheme=Schaake
-#ice_fraction=0
-#ice_content_threshold=0.15

configs/cfe_config_cat_87_pass.txt

@@ -17,11 +17,12 @@ soil_storage=0.585626[m/m]
 K_nash=0.03[]
 K_lf=0.01[]
 nash_storage=0.0,0.0
-giuh_ordinates=0.06,0.51,0.28,0.12,0.03
-num_timesteps=1
-verbosity=1
-surface_runoff_scheme=GIUH
-surface_water_partitioning_scheme=Schaake
+surface_runoff_scheme=NASH_CASCADE
+N_nash_surface=2
+K_nash_surface=0.83089
+nsubsteps_nash_surface=10
+nash_storage_surface=0.0,0.0
+surface_partitioning_scheme=Schaake
 #surface_water_partitioning_scheme=Xinanjiang
 #a_Xinanjiang_inflection_point_parameter=1
 #b_Xinanjiang_shape_parameter=1

configs/cfe_config_laramie_pass.txt

@@ -17,10 +17,14 @@ soil_storage=0.585626[m/m]
 K_nash=0.03[]
 K_lf=0.01[]
 nash_storage=0.0,0.0
-giuh_ordinates=0.06,0.51,0.28,0.12,0.03
-num_timesteps=1
-verbosity=1
+surface_runoff_scheme=NASH_CASCADE
+N_nash_surface=2
+K_nash_surface=0.83089
+nsubsteps_nash_surface=10
+nash_storage_surface=0.0,0.0
 surface_partitioning_scheme=Schaake
+
+
 #surface_partitioning_scheme=Xinanjiang
 #a_Xinanjiang_inflection_point_parameter=1
 #b_Xinanjiang_shape_parameter=1

configs/cfe_config_laramie_pass_aet_rz.txt

@@ -17,16 +17,12 @@ soil_storage=0.585626[m/m]
 K_nash=0.03[]
 K_lf=0.01[]
 nash_storage=0.0,0.0
-giuh_ordinates=0.06,0.51,0.28,0.12,0.03
-num_timesteps=1
-verbosity=1
-surface_runoff_scheme=GIUH
+surface_runoff_scheme=NASH_CASCADE
+N_nash_surface=2
+K_nash_surface=0.83089
+nsubsteps_nash_surface=10
+nash_storage_surface=0.0,0.0
 surface_water_partitioning_scheme=Schaake
-#surface_water_partitioning_scheme=Xinanjiang
-#a_Xinanjiang_inflection_point_parameter=1
-#b_Xinanjiang_shape_parameter=1
-#x_Xinanjiang_shape_parameter=1
-#urban_decimal_fraction=0.0
 aet_rootzone=true
 soil_layer_depths=0.1,0.4,1.0,2.0
 max_rootzone_layer=2

src/bmi_cfe.c

@@ -17,7 +17,7 @@
   infiltration by Ahmad Jan Khattak, May 2024.
 *******************************************************************************************************/
 
-#define CFE_DEBUG 1
+#define CFE_DEBUG 0
 
 #define INPUT_VAR_NAME_COUNT 5
 #define OUTPUT_VAR_NAME_COUNT 14
@@ -460,7 +460,6 @@ int read_init_config_cfe(const char* config_file, cfe_state_struct* model)
     //    - number of Nash lf reservoirs (optional, defaults to 2, ignored if storage values present)
     //    - K_nash
     //    - initial Nash storage values (optional, defaults to 0.0 for all reservoirs according to number)
-    //    - GIUH ordinates
 
     char config_line[max_config_line_length + 1];
 
@@ -1004,10 +1003,11 @@ int read_init_config_cfe(const char* config_file, cfe_state_struct* model)
       return BMI_FAILURE;
     }
     if (is_verbosity_set == FALSE) {
+#if CFE_DEBUG >= 1
       printf("Config param 'verbosity' not found in config file\n");
-      printf("setting verbosity to a high value\n");
-      model->verbosity = 10;
-      return BMI_FAILURE;
+      printf("setting verbosity to a low value\n");
+#endif
+      model->verbosity = 0;
     }
     if (is_infiltration_excess_method_set == FALSE) {
 #if CFE_DEBUG >= 1
@@ -3301,21 +3301,22 @@ extern void mass_balance_check(cfe_state_struct* cfe_ptr){
     direct_residual = cfe_ptr->vol_struct.volin - cfe_ptr->vol_struct.vol_runoff - cfe_ptr->vol_struct.vol_infilt -
                       cfe_ptr->vol_struct.vol_et_from_rain;
 
-    printf("********************* DIRECT RUNOFF MASS BALANCE ***************\n");
-    printf(" Surface runoff   = %8.4lf m\n",cfe_ptr->vol_struct.vol_runoff);
-    printf(" Infiltration     = %8.4lf m\n",cfe_ptr->vol_struct.vol_infilt);
-    printf(" Vol_et_from_rain = %8.4lf m\n",cfe_ptr->vol_struct.vol_et_from_rain);
-    printf(" Direct residual  = %6.4e m\n",direct_residual);  // should equal 0.0
+    printf("******************* PRECIPITATION MASS BALANCE *****************\n");
+    printf(" Volume input        = %8.4lf m\n",cfe_ptr->vol_struct.volin);
+    printf(" Infiltration Excess = %8.4lf m\n",cfe_ptr->vol_struct.vol_runoff);
+    printf(" Infiltration        = %8.4lf m\n",cfe_ptr->vol_struct.vol_infilt);
+    printf(" Vol_et_from_rain    = %8.4lf m\n",cfe_ptr->vol_struct.vol_et_from_rain);
+    printf(" Precip residual     = %6.4e m\n",direct_residual);  // should equal 0.0
     if(!is_fabs_less_than_epsilon(direct_residual,1.0e-12))
       printf("WARNING: DIRECT RUNOFF PARTITIONING MASS BALANCE CHECK FAILED\n");
 
     giuh_residual = cfe_ptr->vol_struct.vol_runoff - cfe_ptr->vol_struct.vol_out_surface - cfe_ptr->vol_struct.vol_end_surface -
                     cfe_ptr->vol_struct.vol_runon_infilt;
     printf("********************* SURFACE MASS BALANCE *********************\n");
     printf(" Volume into surface  = %8.4lf m\n",cfe_ptr->vol_struct.vol_runoff);
-    fprintf(stderr," Volume out surface   = %8.4lf m\n",cfe_ptr->vol_struct.vol_out_surface);
-    fprintf(stderr," Volume end surface   = %8.4lf m\n",cfe_ptr->vol_struct.vol_end_surface);
-    fprintf(stderr," Runon infiltration   = %8.4lf m\n",cfe_ptr->vol_struct.vol_runon_infilt);
+    printf(" Volume out surface   = %8.4lf m\n",cfe_ptr->vol_struct.vol_out_surface);
+    printf(" Volume end surface   = %8.4lf m\n",cfe_ptr->vol_struct.vol_end_surface);
+    printf(" Runon infiltration   = %8.4lf m\n",cfe_ptr->vol_struct.vol_runon_infilt);
     printf(" Surface residual     = %6.4e m\n",giuh_residual);  // should equal zero
     if(!is_fabs_less_than_epsilon(giuh_residual,1.0e-12))
       printf("WARNING: GIUH MASS BALANCE CHECK FAILED\n");
@@ -3330,14 +3331,14 @@ extern void mass_balance_check(cfe_state_struct* cfe_ptr){
     printf(" Volume into soil     = %8.4lf m\n",cfe_ptr->vol_struct.vol_infilt);
     printf(" Volume soil2latflow  = %8.4lf m\n",cfe_ptr->vol_struct.vol_soil_to_lat_flow);
     printf(" Volume soil to GW    = %8.4lf m\n",cfe_ptr->vol_struct.vol_soil_to_gw);
-    printf(" Final volume soil    =  %8.4lf m\n",vol_soil_end);
-    printf(" Volume et from soil  =  %8.4lf m\n",cfe_ptr->vol_struct.vol_et_from_soil);
+    printf(" Final volume soil    = %8.4lf m\n",vol_soil_end);
+    printf(" Volume et from soil  = %8.4lf m\n",cfe_ptr->vol_struct.vol_et_from_soil);
     printf(" Volume soil residual = %6.4e m\n",soil_residual);
     if(!is_fabs_less_than_epsilon(soil_residual,1.0e-12))
       printf("WARNING: SOIL CONCEPTUAL RESERVOIR MASS BALANCE CHECK FAILED\n");
 
     nash_residual = cfe_ptr->vol_struct.vol_in_nash - cfe_ptr->vol_struct.vol_out_nash - vol_in_nash_end;
-    printf("********* NASH CASCADE CONCEPTUAL RESERVOIR MASS BALANCE *******\n");
+    printf("************* NASH CASCADE (SUBSURFACE) MASS BALANCE ***********\n");
     printf(" Volume to Nash   = %8.4lf m\n",cfe_ptr->vol_struct.vol_in_nash);
     printf(" Volume from Nash = %8.4lf m\n",cfe_ptr->vol_struct.vol_out_nash);
     printf(" Final vol. Nash  = %8.4lf m\n",vol_in_nash_end);

src/nash_cascade.c

@@ -60,12 +60,12 @@ double nash_cascade_surface(double runoff_m, double soil_storage_deficit_m,
   double dS    = 0.0;            // change in reservoir storage
   double dS_infil = 0.0;         // change in reservoir storage due to infiltration
   double Q_r;                    // discharge from reservoir
-  double Q_out;                  // discharge at the outlet (the last reservoir) per subtimestep
+  double Q_out = 0.0;            // discharge at the outlet (the last reservoir) per subtimestep
   double Q_infil;                // discharge from reservoirs to soil
-  double Q_to_channel_m = 0.0;  // total outflow to channel per timestep
-  double Q_to_soil_m    = 0.0;  // runon infiltration (losses from surface runoff to soil)
+  double Q_to_channel_m = 0.0;   // total outflow to channel per timestep
+  double Q_to_soil_m    = 0.0;   // runon infiltration (losses from surface runoff to soil)
   double soil_deficit_m = soil_storage_deficit_m; // local variable to track the soil storage deficit
-  double infil_m;
+  double infil_m = 0.0;
 
   nash_params->nash_storage[0] += runoff_m;
 
@@ -75,6 +75,10 @@ double nash_cascade_surface(double runoff_m, double soil_storage_deficit_m,
     //Loop through reservoirs
     for(int i = 0; i < N_nash; i++) {
 
+      // if storage of ith reservoir is zero, move to the next reservoir
+      if (nash_params->nash_storage[i] == 0.0)
+	continue;
+
       // First: infiltration capacity should be satisfied before routing water through Nash reservoirs
 
       // compute runon infiltration