Skip to content
40 changes: 36 additions & 4 deletions src/bmi_cfe.c
Original file line number Diff line number Diff line change
@@ -1,3 +1,4 @@
#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <stdbool.h>
Expand Down Expand Up @@ -975,13 +976,45 @@ int read_init_config_cfe(const char* config_file, cfe_state_struct* model)
return BMI_FAILURE;
}
if (is_gw_storage_set == FALSE) {
Log(SEVERE, "Config param 'gw_storage' not found in config file. Aborting...n");
Log(SEVERE, "Config param 'gw_storage' not found in config file. Aborting...\n");
return BMI_FAILURE;
}

// compute gw storage in meters
// compute gw storage in meters
if ((is_gw_storage_set == TRUE) && (is_gw_max_set == TRUE)) {
model->gw_reservoir.storage_m = model->gw_reservoir.gw_storage * model->gw_reservoir.storage_max_m;
if (!isfinite(model->gw_reservoir.gw_storage) ||
!isfinite(model->gw_reservoir.storage_max_m) ||
!isfinite(model->gw_reservoir.coeff_primary) ||
!isfinite(model->gw_reservoir.exponent_primary) ||
model->gw_reservoir.gw_storage < 0.0 ||
model->gw_reservoir.storage_max_m <= 0.0) {

Comment thread
cmaynard-ngwpc marked this conversation as resolved.
static int gw_bad_param_warned = 0;

if (gw_bad_param_warned == 0) {
Log(WARNING,
"Invalid CFE groundwater reservoir parameters from config. "
"gw_storage=%lf, max_gw_storage=%lf, Cgw=%lf, expon=%lf. "
"This indicates an upstream parameter-generation issue that should be corrected.\n",
model->gw_reservoir.gw_storage,
model->gw_reservoir.storage_max_m,
model->gw_reservoir.coeff_primary,
model->gw_reservoir.exponent_primary);

gw_bad_param_warned = 1;
}
else if (gw_bad_param_warned == 1) {
Log(WARNING,
"Invalid CFE groundwater reservoir parameters occurred again; "
"subsequent occurrences of this message will be suppressed.\n");

gw_bad_param_warned = 2;
}
}

model->gw_reservoir.storage_m =
model->gw_reservoir.gw_storage *
model->gw_reservoir.storage_max_m;
}

if (is_num_timesteps_set == FALSE && strcmp(model->forcing_file, "BMI")) {
Expand Down Expand Up @@ -2297,7 +2330,6 @@ static int Set_value (Bmi *self, const char *name, void *src)
return BMI_SUCCESS;
}


static int Get_component_name (Bmi *self, char * name)
{
strncpy (name, "The CFE Model", BMI_MAX_COMPONENT_NAME);
Expand Down
22 changes: 22 additions & 0 deletions src/cfe.c
Original file line number Diff line number Diff line change
Expand Up @@ -72,7 +72,29 @@ extern void cfe(
//##################################################
// partition rainfall using Schaake function
//##################################################
static int potentialEtIssueWarned = 0;

if (!isfinite(evap_struct->potential_et_m_per_s) ||
evap_struct->potential_et_m_per_s < 0.0) {

if (potentialEtIssueWarned == 0) {
Log(WARNING,
"Invalid CFE potential ET input detected. "
"water_potential_evaporation_flux=%lf, time_step_size=%lf. "
"This indicates an upstream forcing/BMI provider issue.\n",
evap_struct->potential_et_m_per_s,
time_step_size);

potentialEtIssueWarned = 1;
}
else if (potentialEtIssueWarned == 1) {
Log(WARNING,
"Invalid CFE potential ET input occurred again; "
"subsequent occurrences of this message will be suppressed.\n");

potentialEtIssueWarned = 2;
}
}
evap_struct->potential_et_m_per_timestep = evap_struct->potential_et_m_per_s * time_step_size;
evap_struct->reduced_potential_et_m_per_timestep = evap_struct->potential_et_m_per_s * time_step_size;

Expand Down
71 changes: 66 additions & 5 deletions src/conceptual_reservoir.c
Original file line number Diff line number Diff line change
@@ -1,6 +1,7 @@
#ifndef _CONCEPTUAL_RESERVOIR_C
#define _CONCEPTUAL_RESERVOIR_C

#include "logger.h"
#include "conceptual_reservoir.h"


Expand Down Expand Up @@ -40,20 +41,80 @@ extern void conceptual_reservoir_flux_calc(struct conceptual_reservoir *da_reser
// *****************************************************************************
// ------------------ Conceptual Ground Water Reservoir -------------------------
// single outlet reservoir like the NWM V1.2 exponential conceptual gw reservoir
if (da_reservoir->is_exponential == TRUE) {
if (da_reservoir->is_exponential == TRUE) {
static int gw_bad_state_warned = 0;

if (!isfinite(da_reservoir->storage_m) ||
!isfinite(da_reservoir->storage_max_m) ||
!isfinite(da_reservoir->coeff_primary) ||
!isfinite(da_reservoir->exponent_primary) ||
da_reservoir->storage_m <= 0.0 ||
da_reservoir->storage_max_m <= 0.0) {

Comment thread
cmaynard-ngwpc marked this conversation as resolved.
if (gw_bad_state_warned == 0) {
Log(WARNING,
"Invalid CFE groundwater reservoir state or parameters during flux calculation; "
"returning zero groundwater flux for this timestep. "
"storage_m=%lf, storage_max_m=%lf, Cgw=%lf, expon=%lf.",
da_reservoir->storage_m,
da_reservoir->storage_max_m,
da_reservoir->coeff_primary,
da_reservoir->exponent_primary);

gw_bad_state_warned = 1;
}
else if (gw_bad_state_warned == 1) {
Log(WARNING,
"Invalid CFE groundwater reservoir state or parameters occurred again; "
"subsequent occurrences of this message will be suppressed.");

gw_bad_state_warned = 2;
}

return;
}

// calculate the one flux and return.
double exp_term = exp(da_reservoir->exponent_primary * da_reservoir->storage_m / da_reservoir->storage_max_m);

*primary_flux_m = da_reservoir->coeff_primary * (exp_term - 1.0);
*secondary_flux_m = 0.0;
double calculated_flux_m = da_reservoir->coeff_primary * (exp_term - 1.0);

if (!isfinite(calculated_flux_m) || calculated_flux_m < 0.0) {
static int gw_bad_flux_warned = 0;

if (gw_bad_flux_warned == 0) {
Log(WARNING,
"Invalid CFE groundwater flux calculation; "
"returning zero groundwater flux for this timestep. "
"storage_m=%lf, storage_max_m=%lf, Cgw=%lf, expon=%lf, raw_flux=%lf.",
da_reservoir->storage_m,
da_reservoir->storage_max_m,
da_reservoir->coeff_primary,
da_reservoir->exponent_primary,
calculated_flux_m);

gw_bad_flux_warned = 1;
}
else if (gw_bad_flux_warned == 1) {
Log(WARNING,
"Invalid CFE groundwater flux calculation occurred again; "
"subsequent occurrences of this message will be suppressed.");

gw_bad_flux_warned = 2;
}

return;
}

*primary_flux_m = calculated_flux_m;

// cap so flux never exceeds storage, but commented out here per OWP request (9/11/2025)
//if (*primary_flux_m > da_reservoir->storage_m) {
// *primary_flux_m = da_reservoir->storage_m;
//}

return;
}

// *****************************************************************************

// code goes past here iff it is not a single outlet exponential deep groundwater reservoir of the NWM variety
Expand Down