Moved SPATIAL_SNOW from compile-time to run-time options.

Author: tbohn

samples/global.param.sample

@@ -16,14 +16,20 @@ STARTHOUR	00	# hour model simulation starts
 ENDYEAR 	2000	# year model simulation ends
 ENDMONTH	12	# month model simulation ends
 ENDDAY		31	# day model simulation ends
+
+#######################################################################
+# Energy Balance Parameters
+#######################################################################
 FULL_ENERGY 	TRUE	# TRUE = calculate full energy balance; FALSE = compute water balance only
-FROZEN_SOIL	TRUE	# TRUE = calculate frozen soils
+#CLOSE_ENERGY	FALSE	# TRUE = all energy balance calculations (canopy air, canopy snow, ground snow,
+                        # and ground surface) are iterated to minimize the total column error
 
 #######################################################################
 # Soil Temperature Parameters
 # VIC will choose appropriate values for QUICK_FLUX and IMPLICIT depending on values of FULL_ENERGY and FROZEN_SOIL; the user should only need to override VIC's choices in special cases.
 # The other options in this section are only applicable when FROZEN_SOIL is TRUE and their values depend on the application.
 #######################################################################
+FROZEN_SOIL	TRUE	# TRUE = calculate frozen soils
 #QUICK_FLUX	FALSE	# TRUE = use simplified ground heat flux method of Liang et al (1999); FALSE = use finite element method of Cherkauer et al (1999)
 #IMPLICIT	TRUE	# TRUE = use implicit solution for soil heat flux equation of Cherkauer et al (1999), otherwise uses original explicit solution.
 #QUICK_SOLVE	FALSE	# TRUE = Use Liang et al., 1999 formulation for iteration, but explicit finite difference method for final step.
@@ -60,8 +66,6 @@ FROZEN_SOIL	TRUE	# TRUE = calculate frozen soils
 #			# AR_406_FULL 	= multiply by 10 for latent heat AND sensible heat and do NOT apply stability correction; additionally, always use overstory aero_resist for ET (as in 4.1.0).
 #			# AR_410 	= apply stability correction but do NOT multiply by 10 (as in VIC 4.1.0); additionally, always use overstory aero_resist for ET (as in 4.1.0).
 #			# Default 	= AR_406_FULL
-#CLOSE_ENERGY	FALSE	# TRUE = all energy balance calculations (canopy air, canopy snow, ground snow,
-                        # and ground surface) are iterated to minimize the total column error
 
 #######################################################################
 # Meteorological Forcing Disaggregation Parameters
@@ -94,6 +98,11 @@ FROZEN_SOIL	TRUE	# TRUE = calculate frozen soils
 #LOG_MATRIC	FALSE	# TRUE = use the linear interpolation of the logarithm of the matric potential
                         # from the two surrounding layers to estimate the soil moisture drainage from each
                         # layer (Boone and Wetzel, 1996).
+#SPATIAL_SNOW	FALSE	# TRUE = use a uniform distribution to simulate the partial coverage of the
+                        # surface by a thin snowpack.  Coverage is assumed to be uniform after snowfall
+                        # until the pack begins to melt.  If TRUE, VIC will expect an additional column
+                        # in the soil paramter file containing the snow distibution slope parameter
+                        # (= 2 * snow depth below which coverage < 1).
 
 #######################################################################
 # State Files and Parameters

src/README.txt

@@ -619,11 +619,13 @@ Cleanup of compile-time options
 	calc_atmos_energy_bal.c
 	calc_forcing_stats.c (removed)
 	calc_longwave.c
+	calc_surf_energy_bal.c
 	display_current_settings.c
 	frozen_soil.c
 	full_energy.c
 	func_surf_energy_bal.c
 	get_global_param.c
+	ice_melt.c
 	initialize_atmos.c
 	initialize_global.c
 	initialize_model_state.c
@@ -638,11 +640,13 @@ Cleanup of compile-time options
 	runoff.c
 	soil_conduction.c
 	soil_thermal_eqn.c
+	SnowPackEnergyBalance.c
+	solve_snow.c
 	surface_fluxes.c
 	user_def.h
 	vicNl_def.h
 	vicNl.h
-	write_model_state.h
+	write_model_state.c
 
 	Description:
 
@@ -662,6 +666,8 @@ Cleanup of compile-time options
 
 	* Moved CLOSE_ENERGY from compile-time to run-time options.
 
+	* Moved SPATIAL_SNOW from compile-time to run-time options.
+
 	* Moved max array dimensions such as MAX_VEG, MAX_LAYERS, etc. to
 	vicNl_def.h.
 

src/SnowPackEnergyBalance.c

@@ -82,6 +82,7 @@ static char vcid[] = "$Id$";
 	      two elements (surface and overstory); added
 	      options.AERO_RESIST_CANSNOW.				TJB
   2009-Sep-19 Added Added ground flux computation consistent with 4.0.6.	TJB
+  2013-Dec-27 Moved SPATIAL_SNOW from compile-time to run-time options.	TJB
 
 *****************************************************************************/
 double SnowPackEnergyBalance(double TSurf, va_list ap)
@@ -253,18 +254,19 @@ double SnowPackEnergyBalance(double TSurf, va_list ap)
 
   *SensibleHeat = AirDens * Cp * (Tair - TMean) / Ra_used[0];
 
-#if SPATIAL_SNOW  
-  /* Add in Sensible heat flux turbulent exchange from surrounding 
-     snow free patches - if present */
-  if ( SnowCoverFract > 0 ) {
-    *(AdvectedSensibleHeat) = advected_sensible_heat(SnowCoverFract, 
-						     AirDens, Tair, TGrnd, 
-						     Ra_used[0]);
+  if (options.SPATIAL_SNOW) {
+    /* Add in Sensible heat flux turbulent exchange from surrounding 
+       snow free patches - if present */
+    if ( SnowCoverFract > 0 ) {
+      *(AdvectedSensibleHeat) = advected_sensible_heat(SnowCoverFract, 
+						       AirDens, Tair, TGrnd, 
+						       Ra_used[0]);
+    }
+    else (*AdvectedSensibleHeat) = 0;
+  }
+  else {
+    (*AdvectedSensibleHeat) = 0;
   }
-  else (*AdvectedSensibleHeat) = 0;
-#else
-  (*AdvectedSensibleHeat) = 0;
-#endif // SPATIAL_SNOW
 
   /* Convert sublimation terms from m/timestep to kg/m2s */
   VaporMassFlux = *vapor_flux * Density / Dt;

src/calc_surf_energy_bal.c

@@ -172,6 +172,7 @@ double calc_surf_energy_bal(double             Le,
 	      and clarified the descriptions of the SPATIAL_SNOW
 	      option.							TJB
   2013-Jul-25 Added photosynthesis.					TJB
+  2013-Dec-27 Moved SPATIAL_SNOW to options_struct.			TJB
 ***************************************************************/
 {
   extern veg_lib_struct *veg_lib;
@@ -745,21 +746,21 @@ double calc_surf_energy_bal(double             Le,
 
       /** Check for Thin Snowpack which only Partially Covers Grid Cell
 	  exists only if not snowing and snowpack has started to melt **/
-#if SPATIAL_SNOW
-      snow->coverage = calc_snow_coverage(&snow->store_snow, 
-					  soil_con->max_snow_distrib_slope, 
-					  snow_coverage, snow->swq,
-					  old_swq, snow->depth, old_depth, 
-					  (*melt) - snow->vapor_flux, 
-					  &snow->max_snow_depth, snowfall, 
-					  &snow->store_swq, 
-					  &snow->snow_distrib_slope,
-					  &snow->store_coverage);
-      
-#else
-      if ( snow->swq > 0 ) snow->coverage = 1.;
-      else snow->coverage = 0.;
-#endif // SPATIAL_SNOW
+      if (options.SPATIAL_SNOW) {
+        snow->coverage = calc_snow_coverage(&snow->store_snow, 
+					    soil_con->max_snow_distrib_slope, 
+					    snow_coverage, snow->swq,
+					    old_swq, snow->depth, old_depth, 
+					    (*melt) - snow->vapor_flux, 
+					    &snow->max_snow_depth, snowfall, 
+					    &snow->store_swq, 
+					    &snow->snow_distrib_slope,
+					    &snow->store_coverage);
+      }
+      else {
+        if ( snow->swq > 0 ) snow->coverage = 1.;
+        else snow->coverage = 0.;
+      }
 
       if ( snow->surf_temp > 0 ) 
 	energy->snow_flux = ( energy->grnd_flux + energy->deltaH 
@@ -775,9 +776,8 @@ double calc_surf_energy_bal(double             Le,
       snow->surf_temp  = 0;
       snow->pack_temp  = 0;
       snow->coverage   = 0;
-#if SPATIAL_SNOW
-      snow->store_swq = 0.;
-#endif // SPATIAL_SNOW
+      if (options.SPATIAL_SNOW)
+        snow->store_swq = 0.;
     }
     snow->vapor_flux *= -1;
   }

src/display_current_settings.c

@@ -57,6 +57,7 @@ void display_current_settings(int                 mode,
 	      run-time option.						TJB
   2013-Dec-26 Moved CLOSE_ENERGY from compile-time to run-time options.	TJB
   2013-Dec-26 Removed EXCESS_ICE option.				TJB
+  2013-Dec-27 Moved SPATIAL_SNOW from compile-time to run-time options.	TJB
 **********************************************************************/
 {
 
@@ -116,11 +117,6 @@ void display_current_settings(int                 mode,
 #else
   fprintf(stderr,"SPATIAL_FROST\t\tFALSE\n");
 #endif
-#if SPATIAL_SNOW
-  fprintf(stderr,"SPATIAL_SNOW\t\tTRUE\n");
-#else
-  fprintf(stderr,"SPATIAL_SNOW\t\tFALSE\n");
-#endif
 
   fprintf(stderr,"\n");
   fprintf(stderr,"Maximum Array Sizes:\n");
@@ -268,6 +264,10 @@ void display_current_settings(int                 mode,
     fprintf(stderr,"QUICK_SOLVE\t\tTRUE\n");
   else
     fprintf(stderr,"QUICK_SOLVE\t\tFALSE\n");
+  if (options.SPATIAL_SNOW == TRUE)
+    fprintf(stderr,"SPATIAL_SNOW\t\tTRUE\n");
+  else
+    fprintf(stderr,"SPATIAL_SNOW\t\tFALSE\n");
   if (options.SNOW_ALBEDO == USACE)
     fprintf(stderr,"SNOW_ALBEDO\t\tUSACE\n");
   else if (options.SNOW_ALBEDO == SUN1999)

src/get_global_param.c

@@ -107,6 +107,7 @@ global_param_struct get_global_param(filenames_struct *names,
   2013-Dec-26 Added LOG_MATRIC option.						TJB
   2013-Dec-26 Moved CLOSE_ENERGY from compile-time to run-time options.	TJB
   2013-Dec-26 Removed EXCESS_ICE option.				TJB
+  2013-Dec-27 Moved SPATIAL_SNOW from compile-time to run-time options.	TJB
 **********************************************************************/
 {
   extern option_struct    options;
@@ -380,6 +381,11 @@ global_param_struct get_global_param(filenames_struct *names,
 	  options.PLAPSE = TRUE;
 	}
       }
+      else if(strcasecmp("SPATIAL_SNOW",optstr)==0) {
+        sscanf(cmdstr,"%*s %s",flgstr);
+        if(strcasecmp("TRUE",flgstr)==0) options.SPATIAL_SNOW=TRUE;
+        else options.SPATIAL_SNOW = FALSE;
+      }
       else if(strcasecmp("SW_PREC_THRESH",optstr)==0) {
         sscanf(cmdstr,"%*s %s",flgstr);
         options.SW_PREC_THRESH = atof(flgstr);

src/ice_melt.c

@@ -115,6 +115,7 @@ static char vcid[] = "$Id$";
   2012-Feb-08 Renamed depth_full_snow_cover to max_snow_distrib_slope
 	      and clarified the descriptions of the SPATIAL_SNOW
 	      option.								TJB
+  2013-Dec-27 Moved SPATIAL_SNOW from compile-time to run-time options.	TJB
 *****************************************************************************/
 int ice_melt(double            z2,
 	      double            aero_resist,
@@ -623,21 +624,22 @@ int ice_melt(double            z2,
   if (lake->ice_water_eq <= 0.0) {
     lake->ice_water_eq = 0.0;
   }
-#if SPATIAL_SNOW
+  if (options.SPATIAL_SNOW) {
   /*  snow->coverage = calc_snow_coverage(&snow->store_snow,
-      soil_con->max_snow_distrib_slope,
-      old_coverage, snow->swq,
-      old_swq, snow->depth, old_depth,
-      melt + snow->vapor_flux,
-      &snow->max_snow_depth, snowfall,
-      &snow->store_swq,
-      &snow->snow_distrib_slope,
-      &snow->store_coverage);
-  */
-#else
-  if ( snow->swq > 0 ) snow->coverage = 1.;
-  else snow->coverage = 0.;
-#endif
+    soil_con->max_snow_distrib_slope,
+    old_coverage, snow->swq,
+    old_swq, snow->depth, old_depth,
+    melt + snow->vapor_flux,
+    &snow->max_snow_depth, snowfall,
+    &snow->store_swq,
+    &snow->snow_distrib_slope,
+    &snow->store_coverage);
+    */
+  }
+  else {
+    if ( snow->swq > 0 ) snow->coverage = 1.;
+    else snow->coverage = 0.;
+  }
 
   /* Mass balance test */
   MassBalanceError = (InitialSwq - snow->swq) + (InitialIce - LakeIce)

src/initialize_global.c

@@ -125,6 +125,7 @@ void initialize_global() {
   2013-Jul-25 Added CARBON, SHARE_LAYER_MOIST, and VEGLIB_PHOTO.		TJB
   2013-Dec-26 Added LOG_MATRIC option.						TJB
   2013-Dec-26 Moved CLOSE_ENERGY from compile-time to run-time options.	TJB
+  2013-Dec-27 Moved SPATIAL_SNOW from compile-time to run-time options.	TJB
 *********************************************************************/
 
   extern option_struct options;
@@ -172,6 +173,7 @@ void initialize_global() {
   options.SNOW_BAND             = 1;
   options.SNOW_DENSITY          = DENS_BRAS;
   options.SNOW_STEP             = 1;
+  options.SPATIAL_SNOW          = FALSE;
   options.SW_PREC_THRESH        = 0;
   options.TFALLBACK             = TRUE;
   options.VP_INTERP             = TRUE;

src/read_soilparam.c

@@ -122,6 +122,7 @@ soil_con_struct read_soilparam(FILE *soilparam,
 	      option.							TJB
   2013-Jul-25 Added calculation of soil albedo in PAR range.		TJB
   2013-Dec-26 Removed EXCESS_ICE option.				TJB
+  2013-Dec-27 Moved SPATIAL_SNOW from compile-time to run-time options.	TJB
 **********************************************************************/
 {
   void ttrim( char *string );
@@ -598,16 +599,19 @@ soil_con_struct read_soilparam(FILE *soilparam,
       temp.FS_ACTIVE = (char)tempint;
 
       /* read minimum snow depth for full coverage */
-#if SPATIAL_SNOW
-      token = strtok (NULL, delimiters);
-      while (token != NULL && (length=strlen(token))==0) token = strtok (NULL, delimiters);
-      if( token == NULL ) {
-        sprintf(ErrStr,"ERROR: Can't find values for SPATIAL SNOW in soil file\n");
-        nrerror(ErrStr);
+      if (options.SPATIAL_SNOW) {
+        token = strtok (NULL, delimiters);
+        while (token != NULL && (length=strlen(token))==0) token = strtok (NULL, delimiters);
+        if( token == NULL ) {
+          sprintf(ErrStr,"ERROR: Can't find values for SPATIAL SNOW in soil file\n");
+          nrerror(ErrStr);
+        }
+        sscanf(token, "%lf", &tempdbl);
+        temp.max_snow_distrib_slope = tempdbl;
+      }
+      else {
+        temp.max_snow_distrib_slope = 0;
       }
-      sscanf(token, "%lf", &tempdbl);
-      temp.max_snow_distrib_slope = tempdbl;
-#endif // SPATIAL_SNOW
 
       /* read slope of frozen soil distribution */
 #if SPATIAL_FROST
@@ -716,12 +720,12 @@ soil_con_struct read_soilparam(FILE *soilparam,
       /**********************************************
         Validate Spatial Snow/Frost Params
       **********************************************/
-#if SPATIAL_SNOW
-      if (temp.max_snow_distrib_slope < 0.0) {
-        sprintf(ErrStr,"max_snow_distrib_slope (%f) must be positive.\n", temp.max_snow_distrib_slope);
-        nrerror(ErrStr);
+      if (options.SPATIAL_SNOW) {
+        if (temp.max_snow_distrib_slope < 0.0) {
+          sprintf(ErrStr,"max_snow_distrib_slope (%f) must be positive.\n", temp.max_snow_distrib_slope);
+          nrerror(ErrStr);
+        }
       }
-#endif // SPATIAL_SNOW
 
 #if SPATIAL_FROST
       if (temp.frost_slope < 0.0) {

src/read_soilparam_arc.c

@@ -131,6 +131,7 @@ soil_con_struct read_soilparam_arc(FILE *soilparam,
 	      and clarified the descriptions of the SPATIAL_SNOW
 	      option.								TJB
   2013-Dec-26 Removed EXCESS_ICE option.				TJB
+  2013-Dec-27 Moved SPATIAL_SNOW from compile-time to run-time options.	TJB
 **********************************************************************/
 {
   extern option_struct options;
@@ -406,12 +407,15 @@ soil_con_struct read_soilparam_arc(FILE *soilparam,
 
     /** Minimum Snow Depth of Full Coverage **/
     fscanf(soilparam,"%s",tmpstr);
-#if SPATIAL_SNOW
-    strcpy(namestr,soilparamdir);
-    strcat(namestr,"/");
-    strcat(namestr,tmpstr);
-    temp.max_snow_distrib_slope = read_arcinfo_value(namestr,temp.lat,temp.lng);
-#endif // SPATIAL_SNOW
+    if (options.SPATIAL_SNOW) {
+      strcpy(namestr,soilparamdir);
+      strcat(namestr,"/");
+      strcat(namestr,tmpstr);
+      temp.max_snow_distrib_slope = read_arcinfo_value(namestr,temp.lat,temp.lng);
+    }
+    else {
+      temp.max_snow_distrib_slope = 0;
+    }
 
     /** Slope of Frozen Soil Distribution **/
     fscanf(soilparam,"%s",tmpstr);