Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Update to code for Esoil in sparse vegetation. #315

Merged
merged 5 commits into from
Dec 9, 2015
Merged

Update to code for Esoil in sparse vegetation. #315

Show file tree
Hide file tree
Changes from 1 commit
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
649243c
Update to code for Esoil in sparse vegetation. I had intended to add…
tbohn Nov 17, 2015
d958240
Added comments.
tbohn Nov 30, 2015
c3b772b
Added comments.
tbohn Nov 30, 2015
19ce733
Replaced hard-coded numbers with constants defined in vicNl_def.h.
tbohn Dec 6, 2015
c48f178
Replaced variable rarc_soil with constant RARC_SOIL.
tbohn Dec 7, 2015
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
3 changes: 2 additions & 1 deletion src/arno_evap.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -85,6 +85,7 @@ double arno_evap(layer_data_struct *layer,
double max_infil;
double Evap;
double tmpsum;
double rarc_soil = 100.0;
Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@jhamman commented: why 100.0?
@tbohn's reply: 100 gives better results than 0 or 200, from tests I conducted at eddy covariance towers. I didn't do a fancy Bayesian parameter estimation, not enough time for that. I apologize for making this a hard-wired value - I considered making it another input soil parameter, but we have so many of those already. And this likely won't noticeably affect anyone who has been using the traditional VIC set-up, since bare soil rarely occurs in the traditional parameters.


Evap = 0;

Expand All @@ -97,7 +98,7 @@ double arno_evap(layer_data_struct *layer,

/* Calculate the potential bare soil evaporation (mm/time step) */

Epot = penman(air_temp, elevation, rad, vpd, ra, 0.0, 0.0) * delta_t / SEC_PER_DAY;
Epot = penman(air_temp, elevation, rad, vpd, ra, 0.0, rarc_soil) * delta_t / SEC_PER_DAY;

/**********************************************************************/
/* Compute temporary infiltration rate based on given soil_moist. */
Expand Down
40 changes: 30 additions & 10 deletions src/func_surf_energy_bal.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -300,6 +300,9 @@ double func_surf_energy_bal(double Ts, va_list ap)
double tmp_displacement[3];
double tmp_roughness[3];
double tmp_ref_height[3];
double ga_veg;
double ga_bare;
double ga_average;

/************************************
Read variables from variable list
Expand Down Expand Up @@ -688,6 +691,32 @@ double func_surf_energy_bal(double Ts, va_list ap)
else
Ra_used[0] = HUGE_RESIST;

if (veg_var->vegcover < 1) {
if (Ra_used[0] > 0) {
ga_veg = 1/Ra_used[0];
tmp_wind[0] = wind[0];
tmp_wind[1] = -999.;
tmp_wind[2] = -999.;
tmp_height = soil_con->rough/0.123;
tmp_displacement[0] = calc_veg_displacement(tmp_height);
tmp_roughness[0] = soil_con->rough;
tmp_ref_height[0] = 10;
Error = CalcAerodynamic(0,0,0,soil_con->snow_rough,soil_con->rough,0,Ra_bare,tmp_wind,tmp_displacement,tmp_ref_height,tmp_roughness);
Ra_bare[0] /= StabilityCorrection(tmp_ref_height[0], tmp_displacement[0], TMean, Tair, tmp_wind[0], tmp_roughness[0]);
if (Ra_bare[0] > 0) {
ga_bare = 1/Ra_bare[0];
ga_average = veg_var->vegcover*ga_veg + (1-veg_var->vegcover)*ga_bare;
Ra_used[0] = 1/ga_average;
}
else {
Ra_used[0] = 0;
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Should this be zero or infinite (or does it not matter because it means that the value is not used?)

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It should be zero. Ra of 0 means infinite conductance. Perhaps some of the confusion here is due to the naming of the variables. Before I added this feature, Ra_used was the Ra of the entire veg tile. Since vegcover was assumed == 1.0 in the past, Ra_used was also the Ra of the vegetated portion (which was 100% of the tile). When I added this block, I intended to adjust Ra_used for the presence of exposed soil. In this case, the old Ra_used now initially only pertains to the vegetated portion (which is < 100% of the tile), and gets reassigned to the value pertaining to the entire tile later on.

So, what we have is:
If vegetated Ra is 0, then vegetated conductance (ga) is infinite, so total tile ga is also infinite, so Ra_used = 1/ga is 0.
If exposed soil Ra is 0, Ra_used also = 0.

What if I add another intermediate variable to make this clearer: Ra_veg. Ra_veg will be assigned = Ra_used at the beginning of this block. Then, ga_veg = 1/Ra_veg. ga_avg = weighted average as currently stated. Then, Ra_used = 1/ga_avg. This will hopefully make it clearer that Ra_used = 0 if either Ra_veg or Ra_bare are 0.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes - that would be much easier to understand - please go ahead.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

OK, done.

}
}
else {
Ra_used[0] = 0;
}
}

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@jhamman commented: Can you add some comments here. There are a few hard coded literals and a lot of moving pieces without any explanation.
@tbohn's reply: OK, will do.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Still needs those comments. The hardcoded values need to be moved to header files.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

But I added several lines of comments... It needs more? Which lines in particular need comments?

/*************************************************
Compute Evapotranspiration if not snow covered

Expand All @@ -709,20 +738,11 @@ double func_surf_energy_bal(double Ts, va_list ap)
transp[i] = layer[i].evap;
layer[i].evap = 0;
}
tmp_wind[0] = wind[0];
tmp_wind[1] = -999.;
tmp_wind[2] = -999.;
tmp_height = soil_con->rough/0.123;
tmp_displacement[0] = calc_veg_displacement(tmp_height);
tmp_roughness[0] = soil_con->rough;
tmp_ref_height[0] = 10;
Error = CalcAerodynamic(0,0,0,soil_con->snow_rough,soil_con->rough,0,Ra_bare,tmp_wind,tmp_displacement,tmp_ref_height,tmp_roughness);
Ra_bare[0] /= StabilityCorrection(tmp_ref_height[0], tmp_displacement[0], TMean, Tair, tmp_wind[0], tmp_roughness[0]);
Evap *= veg_var->vegcover;
Evap += (1-veg_var->vegcover)
* arno_evap(layer, surf_atten*NetBareRad, Tair, vpd,
depth[0], max_moist * depth[0] * 1000.,
elevation, b_infilt, Ra_bare[0], delta_t,
elevation, b_infilt, Ra_used[0], delta_t,
resid_moist[0], frost_fract);
for (i=0; i<options.Nlayer; i++) {
layer[i].evap = veg_var->vegcover*transp[i] + (1-veg_var->vegcover)*layer[i].evap;
Expand Down