Feature/debug cesm driver (#642)

Diana Gergel · Joe Hamman · commit d1872e72936e · 2017-01-06T10:18:31.000-07:00
* update command line usage docs * update release notes * Fixed image driver history file name (#635) * Fixed image driver history file name timetamp - now is timstep-beginning * Added release notes about PR #635 (fixing image driver history filename timestamp) * Committing a fix to vic_restore.c Added an if statement so that when VIC is run on multiple nodes, only the master node will validate state file dimensions and coordinate variables. * Use multiple processors for image restart tests * dummy forcings are now created while initializing the model prior to communicating with the coupler * removed decimal point for missing value for unsigned ints * initial fix for getting state file date/time info * updated fix for getting statefile time info from dmy_current (period-ending) * Fixed a bug when calculating steps of forcing to skip (#639) * Fixed a bug when calculating steps of forcing to skip - should round to integer first before converting to int type to prevent possible 1-timestep wrong offset * Added bug fix release notes to PR#639 * fixed bug that passed fillvalue as double when it was supposed to be an int in gather_put_nc_field_int fctn * Added release notes about PR #638. (#640) (Fixing a problem with image restarts when using multiple processors). * use appropriate fillvalue datatypes in gather_put_nc_field fctns * add caseid to statefile name by passing case metadata in to vic_cesm_run function * Revert "add caseid to statefile name by passing case metadata in to vic_cesm_run function because adding this logic to vic_cesm_init instead" This reverts commit 7e28f13. * assign caseid to statefile name in vic_cesm_init * fix warnings for SNOW_BAND and AGGFREQ * fixed typo in comment * Added check to ensure mask variable in the input domain file is integer type (#645) * Added check to ensure mask varible in domain file is integer type * Cleaned up some printing lines * Added get_nc_var_type.c * Minor update of comments and ReleaseNotes * Added domain file description in docs * Small fix of table in docs * Added ncdump -h results for domain file in docs * Minor docs update * fix typo in primary reference (#646) * Fix bug make lastday (#647) * Added check to ensure mask varible in domain file is integer type * Cleaned up some printing lines * Added get_nc_var_type.c * Minor update of comments and ReleaseNotes * Added domain file description in docs * Small fix of table in docs * Added ncdump -h results for domain file in docs * Minor docs update * The input arguments to function `make_lastday` are sometimes in a wrong order - fixed the bug * Added ReleaseNotes regarding PR#647 * Revert "Merge branch 'hotfix/5.0.1' of github.com:UW-Hydro/VIC into hotfix/5.0.1" This reverts commit a9a202b, reversing changes made to bbf7d9a. * pressure was being converted to kPa, should be in Pa * Revert "pressure was being converted to kPa, should be in Pa" This reverts commit 72a56a5. * pressure was being converted to kPa, should be in Pa * convert pressure to kPa when reading it in rather than Pa since this was updated in VIC 5 * Revert "convert pressure to kPa when reading it in rather than Pa since this was updated in VIC 5" This reverts commit e7192c5. * fixed bug with using the correct local domain grid cells and units issues * address PR review comments * add change log entry
diff --git a/docs/Development/ReleaseNotes.md b/docs/Development/ReleaseNotes.md
@@ -197,6 +197,16 @@ This is a major update from VIC 4. The VIC 5.0.0 release aims to have nearly ide
 
 	Fixed a bug where volumetric heat capacity of water should be used in `func_canopy_energy_bal` (previously specific heat capacity was used).
 
+6. Fixes for the CESM driver ([GH#642](https://github.com/UW-Hydro/VIC/pull/642))
+
+	1. Using correct fill value datatypes in MPI Gather steps 
+	2. Updated state file name time step to be period-ending rather than period-beginning
+	3. Set the state file name to the RASM case ID
+	4. Removed decimal point for missing values for unsigned integers
+	5. Create dummy forcings when initializing the model (so that there is forcing data for the first time step)
+	6. Changed pressure units from kPa to Pa
+	7. Fixed bug that prevented using the correct local domain grid cells in `cesm_put_data.c`
+	8. Changed reference temperature units from Celsius to Kelvin in `cesm_put_data.c`
 ------------------------------
 
 ## VIC 4.2.d [![DOI](https://zenodo.org/badge/doi/10.5281/zenodo.56058.svg)](http://dx.doi.org/10.5281/zenodo.56058)
diff --git a/vic/drivers/cesm/src/cesm_interface_c.c b/vic/drivers/cesm/src/cesm_interface_c.c
@@ -96,6 +96,9 @@ vic_cesm_init(vic_clock     *vclock,
     // populate model state, either using a cold start or from a restart file
     vic_populate_model_state(trim(cmeta->starttype));
 
+    // initialize forcings
+    vic_force();
+
     // initialize output structures
     vic_init_output(&dmy_current);
 
@@ -137,6 +140,7 @@ vic_cesm_run(vic_clock *vclock)
 
     // if save:
     if (vclock->state_flag) {
+        // write state file
         vic_store(&dmy_current, state_filename);
         write_rpointer_file(state_filename);
     }
diff --git a/vic/drivers/cesm/src/cesm_put_data.c b/vic/drivers/cesm/src/cesm_put_data.c
@@ -140,12 +140,12 @@ vic_cesm_put_data()
 
                 // 2m reference temperature
                 // CESM units: K
-                l2x_vic[i].l2x_Sl_tref += AreaFactor * force->air_temp[NR];
+                l2x_vic[i].l2x_Sl_tref += AreaFactor * (force[i].air_temp[NR] + CONST_TKFRZ);
 
                 // 2m reference specific humidity
                 // CESM units: g/g
                 l2x_vic[i].l2x_Sl_qref += AreaFactor * CONST_EPS *
-                                          force->vp[NR] / force->pressure[NR];
+                                          force[i].vp[NR] / force[i].pressure[NR];
 
                 // Albedo Note: VIC does not partition its albedo, all returned
                 // values will be the same
@@ -156,10 +156,10 @@ vic_cesm_put_data()
                 // force->NetShortAtmos net shortwave flux (+ down)
                 // SWup = force->shortwave[NR] - energy.NetShortAtmos
                 // Set the albedo to zero for the case where there is no shortwave down
-                if (force->shortwave[NR] > 0.) {
+                if (force[i].shortwave[NR] > 0.) {
                     albedo = AreaFactor *
-                             (force->shortwave[NR] - energy.NetShortAtmos) /
-                             force->shortwave[NR];
+                             (force[i].shortwave[NR] - energy.NetShortAtmos) /
+                             force[i].shortwave[NR];
                 }
                 else {
                     albedo = 0.;
@@ -184,7 +184,7 @@ vic_cesm_put_data()
 
                 // 10m wind
                 // CESM units: m/s
-                l2x_vic[i].l2x_Sl_u10 += AreaFactor * force->wind[NR];
+                l2x_vic[i].l2x_Sl_u10 += AreaFactor * force[i].wind[NR];
 
                 // dry deposition velocities (optional)
                 // CESM units: ?
@@ -259,7 +259,7 @@ vic_cesm_put_data()
                 // upward longwave heat flux
                 // CESM units: W m-2
                 l2x_vic[i].l2x_Fall_lwup += AreaFactor *
-                                            (force->longwave[NR] -
+                                            (force[i].longwave[NR] -
                                              energy.NetLongAtmos);
 
                 // evaporation water flux
@@ -278,7 +278,7 @@ vic_cesm_put_data()
 
                 // heat flux shortwave net
                 l2x_vic[i].l2x_Fall_swnet += AreaFactor *
-                                             (force->shortwave[NR] -
+                                             (force[i].shortwave[NR] -
                                               energy.NetShortAtmos);
 
                 // co2 flux **For testing set to 0
diff --git a/vic/drivers/cesm/src/get_global_param.c b/vic/drivers/cesm/src/get_global_param.c
@@ -385,6 +385,12 @@ get_global_param(FILE *gp)
                             "control file.");
                 }
             }
+            else if (strcasecmp("SNOW_BAND", optstr) == 0) {
+	        sscanf(cmdstr, "%*s %s", flgstr);
+		if (str_to_bool(flgstr)) {
+                    options.SNOW_BAND = SNOW_BAND_TRUE_BUT_UNSET;
+                }
+            }
             else if (strcasecmp("LAKES", optstr) == 0) {
                 sscanf(cmdstr, "%*s %s", flgstr);
                 if (strcasecmp("FALSE", flgstr) == 0) {
@@ -415,6 +421,9 @@ get_global_param(FILE *gp)
             else if (strcasecmp("OUTVAR", optstr) == 0) {
                 ; // do nothing
             }
+            else if (strcasecmp("AGGFREQ", optstr) == 0) {
+		; // do nothing
+            }
             else if (strcasecmp("OUTPUT_STEPS_PER_DAY", optstr) == 0) {
                 ; // do nothing
             }
diff --git a/vic/drivers/cesm/src/vic_cesm_start.c b/vic/drivers/cesm/src/vic_cesm_start.c
@@ -46,6 +46,11 @@ vic_cesm_start(vic_clock     *vclock,
     // Driver specific settings
     if (mpi_rank == VIC_MPI_ROOT) {
         strcpy(filenames.global, GLOBALPARAM);
+
+	// assign case name to state file name
+    	strncpy(filenames.statefile, trim(cmeta->caseid),
+                sizeof(filenames.statefile));
+
         // read global settings
         filep.globalparam = open_file(filenames.global, "r");
         get_global_param(filep.globalparam);
diff --git a/vic/drivers/cesm/src/vic_force.c b/vic/drivers/cesm/src/vic_force.c
@@ -127,8 +127,8 @@ vic_force(void)
     for (j = 0; j < NF; j++) {
         for (i = 0; i < local_domain.ncells_active; i++) {
             // CESM units: Pa
-            // VIC units: kPa
-            force[i].pressure[j] = x2l_vic[i].x2l_Sa_pbot / PA_PER_KPA;
+            // VIC units: Pa
+            force[i].pressure[j] = x2l_vic[i].x2l_Sa_pbot;
         }
     }
 
diff --git a/vic/drivers/shared_image/src/vic_init_output.c b/vic/drivers/shared_image/src/vic_init_output.c
@@ -224,8 +224,8 @@ initialize_history_file(nc_file_struct *nc,
     double                    *dvar;
 
 
-    // This could be further refined but for now, I've choosen a file naming
-    // Convention that goes like this:
+    // This could be further refined but for now, I've chosen a file naming
+    // convention that goes like this:
     switch (stream->agg_alarm.freq) {
     // If FREQ_NDAYS -- filename = result_dir/prefix.YYYY-MM-DD.nc
     case FREQ_NDAYS:
diff --git a/vic/drivers/shared_image/src/vic_store.c b/vic/drivers/shared_image/src/vic_store.c
@@ -48,25 +48,45 @@ vic_store(dmy_struct *dmy_current,
     size_t                     k;
     size_t                     m;
     size_t                     p;
+    double                     offset;
+    double                     time_num;
+    double                     end_time_num;
     int                       *ivar = NULL;
     double                    *dvar = NULL;
     size_t                     d2start[2];
     size_t                     d3start[3];
     size_t                     d4start[4];
     size_t                     d5start[5];
     size_t                     d6start[6];
+    dmy_struct                 end_time_date;
     nc_file_struct             nc_state_file;
     nc_var_struct             *nc_var;
 
     set_nc_state_file_info(&nc_state_file);
 
     // only open and initialize the netcdf file on the first thread
     if (mpi_rank == VIC_MPI_ROOT) {
+
+	// advance dmy_current by one timestep since dmy_current is the 
+    	// timestep-beginning timestamp, and state file date should be 
+    	// the end of the current time step 
+    	dt_seconds_to_time_units(global_param.time_units, global_param.dt,
+                                 &offset);
+    	time_num = date2num(global_param.time_origin_num, dmy_current, 0,
+                            global_param.calendar, global_param.time_units);
+    	end_time_num = time_num + offset;
+
+    	// allocate dmy struct for end of current time step 
+    	num2date(global_param.time_origin_num, end_time_num, 0.,
+                                 global_param.calendar, global_param.time_units,
+                                 &end_time_date);
+	
+
         // create netcdf file for storing model state
         sprintf(filename, "%s.%04i%02i%02i_%05u.nc",
-                filenames.statefile, global_param.stateyear,
-                global_param.statemonth, global_param.stateday,
-                global_param.statesec);
+                filenames.statefile, end_time_date.year,
+                end_time_date.month, end_time_date.day,
+                end_time_date.dayseconds);
 
         initialize_state_file(filename, &nc_state_file, dmy_current);
 
@@ -350,7 +370,7 @@ vic_store(dmy_struct *dmy_current,
             }
             gather_put_nc_field_int(nc_state_file.nc_id,
                                     nc_var->nc_varid,
-                                    nc_state_file.d_fillvalue,
+                                    nc_state_file.i_fillvalue,
                                     d4start, nc_var->nc_counts, ivar);
             for (i = 0; i < local_domain.ncells_active; i++) {
                 ivar[i] = nc_state_file.i_fillvalue;
@@ -376,7 +396,7 @@ vic_store(dmy_struct *dmy_current,
             }
             gather_put_nc_field_int(nc_state_file.nc_id,
                                     nc_var->nc_varid,
-                                    nc_state_file.d_fillvalue,
+                                    nc_state_file.i_fillvalue,
                                     d4start, nc_var->nc_counts, ivar);
             for (i = 0; i < local_domain.ncells_active; i++) {
                 ivar[i] = nc_state_file.i_fillvalue;
@@ -813,7 +833,7 @@ vic_store(dmy_struct *dmy_current,
         }
         gather_put_nc_field_int(nc_state_file.nc_id,
                                 nc_var->nc_varid,
-                                nc_state_file.d_fillvalue,
+                                nc_state_file.i_fillvalue,
                                 d2start, nc_var->nc_counts, ivar);
         for (i = 0; i < local_domain.ncells_active; i++) {
             ivar[i] = nc_state_file.i_fillvalue;
@@ -826,7 +846,7 @@ vic_store(dmy_struct *dmy_current,
         }
         gather_put_nc_field_int(nc_state_file.nc_id,
                                 nc_var->nc_varid,
-                                nc_state_file.d_fillvalue,
+                                nc_state_file.i_fillvalue,
                                 d2start, nc_var->nc_counts, ivar);
         for (i = 0; i < local_domain.ncells_active; i++) {
             ivar[i] = nc_state_file.i_fillvalue;
@@ -972,7 +992,7 @@ vic_store(dmy_struct *dmy_current,
         }
         gather_put_nc_field_int(nc_state_file.nc_id,
                                 nc_var->nc_varid,
-                                nc_state_file.d_fillvalue,
+                                nc_state_file.i_fillvalue,
                                 d2start, nc_var->nc_counts, ivar);
         for (i = 0; i < local_domain.ncells_active; i++) {
             ivar[i] = nc_state_file.i_fillvalue;
diff --git a/vic/drivers/shared_image/src/vic_write.c b/vic/drivers/shared_image/src/vic_write.c
@@ -194,7 +194,7 @@ vic_write(stream_struct  *stream,
                 }
                 gather_put_nc_field_schar(nc_hist_file->nc_id,
                                           nc_hist_file->nc_vars[k].nc_varid,
-                                          nc_hist_file->d_fillvalue,
+                                          nc_hist_file->c_fillvalue,
                                           dstart, dcount, cvar);
             }
             else {
diff --git a/vic/vic_run/include/vic_def.h b/vic/vic_run/include/vic_def.h
@@ -52,7 +52,7 @@
 /***** Model Constants *****/
 #define MAXSTRING    2048
 #define MISSING      -99999.   /**< missing value */
-#define MISSING_USI  99999.    /**< missing value for unsigned ints */
+#define MISSING_USI  99999    /**< missing value for unsigned ints */
 #define MISSING_S    "MISSING"    /**< missing value for strings */
 #define NODATA_VH    -1        /**< missing value for veg_hist inputs */
 #define NODATA_VEG   -1        /**< flag for veg types not in grid cell */

Original file line number	Diff line number	Diff line change
`@@ -385,6 +385,12 @@ get_global_param(FILE *gp)`
`385`	`385`	`"control file.");`
`386`	`386`	`}`
`387`	`387`	`}`
	`388`	`+ else if (strcasecmp("SNOW_BAND", optstr) == 0) {`
	`389`	`+ sscanf(cmdstr, "%*s %s", flgstr);`
	`390`	`+ if (str_to_bool(flgstr)) {`
	`391`	`+ options.SNOW_BAND = SNOW_BAND_TRUE_BUT_UNSET;`
	`392`	`+ }`
	`393`	`+ }`
`388`	`394`	`else if (strcasecmp("LAKES", optstr) == 0) {`
`389`	`395`	`sscanf(cmdstr, "%*s %s", flgstr);`
`390`	`396`	`if (strcasecmp("FALSE", flgstr) == 0) {`
`@@ -415,6 +421,9 @@ get_global_param(FILE *gp)`
`415`	`421`	`else if (strcasecmp("OUTVAR", optstr) == 0) {`
`416`	`422`	`; // do nothing`
`417`	`423`	`}`
	`424`	`+ else if (strcasecmp("AGGFREQ", optstr) == 0) {`
	`425`	`+ ; // do nothing`
	`426`	`+ }`
`418`	`427`	`else if (strcasecmp("OUTPUT_STEPS_PER_DAY", optstr) == 0) {`
`419`	`428`	`; // do nothing`
`420`	`429`	`}`
Original file line number	Diff line number	Diff line change
`@@ -127,8 +127,8 @@ vic_force(void)`
`127`	`127`	`for (j = 0; j < NF; j++) {`
`128`	`128`	`for (i = 0; i < local_domain.ncells_active; i++) {`
`129`	`129`	`// CESM units: Pa`
`130`		`- // VIC units: kPa`
`131`		`- force[i].pressure[j] = x2l_vic[i].x2l_Sa_pbot / PA_PER_KPA;`
	`130`	`+ // VIC units: Pa`
	`131`	`+ force[i].pressure[j] = x2l_vic[i].x2l_Sa_pbot;`
`132`	`132`	`}`
`133`	`133`	`}`
`134`	`134`
Original file line number	Diff line number	Diff line change
`@@ -194,7 +194,7 @@ vic_write(stream_struct *stream,`
`194`	`194`	`}`
`195`	`195`	`gather_put_nc_field_schar(nc_hist_file->nc_id,`
`196`	`196`	`nc_hist_file->nc_vars[k].nc_varid,`
`197`		`- nc_hist_file->d_fillvalue,`
	`197`	`+ nc_hist_file->c_fillvalue,`
`198`	`198`	`dstart, dcount, cvar);`
`199`	`199`	`}`
`200`	`200`	`else {`