BROOK90 - EVP - INTERCEPTION AND TRANSPIRATION

INTER - INTER24 - PLNTRES - TBYLAYER

September 30, 2015

BROOK90 routines in this module relate to interception and actual transpiration. Subroutines INTER and INTER24, which handle the interception of rain or snow by the plant canopy are equivalent routines. INTER24 is used when parameter NPINT is 1 and precipitation is input once a day; it assumes that the daily precipitation all occurs within DURATN hours in the middle of the day. INTER is used when NPINT > 1 and precipitation is input more than once a day; it assumes that precipitation rate is constant through the precipitation time step. INTER and INTER24 are used both for rain and snow, with different calling parameters and variables. PLNTRES calculates parameters related to rhizosphere, root, and xylem resistance; it is called once at the beginning of each day. These parameters affect only soil water supply rate, not potential transpiration. TBYLAYER calculates the daily transpiration from each soil layer from the potential transpiration (PTRAN) and the total soil water potential in each layer (PSITI).

Subroutine INTER - interception for precipitation interval

Older studies of rain and snow interception regressed throughfall on precipitation, but such interpretation ignored the fact that energy supply rather than water supply may limit interception and also ignores storm duration/intensity and interstorm interval. Detailed simulation models of rain interception over time through a storm have been developed (Rutter et al. 1972) but these are too complex to include in a hydrologic model like BROOK90. Much less is known about the complicated process of evaporation of intercepted snow.

Subroutine INTER accounts in the simplest way for the concepts of catch rate, evaporation rate, and canopy capacity. The same algorithm is applied to both rain and snow, which are considered to behave independently with respect to their interception. INTER is used when precipitation data are input more than once a day in a precip. interval file (PINT > 1). If only daily precipitation is input, then the modified procedure of INTER24 is used.

The conservation of mass equation for rain interception can be written as

dS/dt = C - I - D

where S is the amount of water stored on the canopy (mm), C is the catch rate, or rate of water input to the canopy, I is the rate of evaporation of intercepted water, and D is the drip rate, or rate of transfer of liquid water to the ground. The same equation applies to snow or mixed snow and rain, when any solid-liquid phase change is ignored and D includes all rain or snow blowing or falling from canopy to ground.

BROOK90 ignores D by defining C as a net catch rate (C - D), or only the portion of the catch that will sooner or later evaporate, so, from the Flow Chart,

d INTR / dt = RINT - IRVP

for rain, and

d INTS / dt = SINT - ISVP

for snow, where INTR and INTS are the canopy storages, RINT and SINT are the net catch rates, and IRVP and ISVP are the evaporation rates, for rain and snow respectively.

BROOK90 assumes that interception catch rates, RINT and SINT, are a constant fraction of rainfall or snowfall until the canopy reaches a storage capacity. Until the capacity is reached, RINT and SINT are assumed to be linear functions of LAI and SAI, so that

RINT = (FRINTL * LAI + FRINTS * SAI) * RFAL

and

SINT = (FSINTL * LAI + FSINTS * SAI) * SFAL

where RFAL and SFAL are rainfall rate and snowfall rate as determined from subroutine SNOFRAC, FRINTL and FSINTL are the catch fraction per unit LAI for rain and snow, respectively, and FRINTS and FSINTS are the catch fraction per unit SAI for rain and snow, respectively.

The canopy has capacities or maximum values of INTR and INTS that depend on LAI and SAI. In BROOK90 these dependencies are assumed linear. The parameters CINTRL and CINTRS are the capacities for intercepted rain per unit LAI and SAI respectively, so that INTRMX, the capacity for rain, is

INTRMX = CINTRL * LAI + CINTRS * SAI.

For snow,

INTSMX = CINTSL * LAI + CINTSS * SAI,

and the capacity parameters are generally larger than for rain. The eight interception parameters, FRINTL, FRINTS, FSINTL, FSINTS, CINTRL, CINTRS, CINTSL, and CINTSS, only control interception loss in small storms; interception loss in large storms is controlled by the evaporation rate of intercepted water (PINT) and the storm intensity and duration.

The rate at which intercepted water evaporates (PINT) is calculated from the Shuttleworth-Wallace equations by calling subroutine SWPE (Section PET) with the canopy resistance r_c = 0. The soil surface resistance (RSS) is not reduced for the PINT calculation. The MSBDAYNIGHT routine does this separately for daytime and for nighttime weather variables and the results are weighted by daylength (DAYLEN) to produce PINT. PINT is considered to be constant throughout the daily time step; its actual diurnal variation is ignored.

The canopy is considered to be either completely wetted or completely dry. Partial canopy wetting and drying is not treated in BROOK90, though it is a key component of specific models of the interception process (Rutter et al. 1972). Subroutine INTER determines the actual catch rate (RINT or SINT) and the actual evaporation rate (IRVP or ISVP) for the precipitation time step in the three cases that the canopy dries during the timestep, the canopy wets but does not reach capacity, and the canopy reaches capacity. The routine appropriately handles the case of a wet canopy with decreasing capacity because of decreasing LAI or SAI by allowing RINT or SINT to be negative.

Subroutine INTER24 - daily interception

Proper representation and integration of the interception process is a problem for hydrologic models that use a daily interval for precipitation input (NPINT = 1), because the storm duration is not known. For a brief, intense storm, the canopy wets once and the interception loss is limited primarily by canopy capacity. For a low intensity, all day storm, the canopy stays wet and the interception loss is limited primarily by the potential interception, PINT. This problem is worst when only daily precipitation is known, and decreases as precipitation is given at shorter intervals.

Subroutine INTER24 was developed because the use of subroutine INTER for daily precipitation consistently produced too much interception. INTER24 is a modification of INTER that loops through the procedure every hour,using the PINT rate for each hour. DURATN is a parameter that specifies the average hourly duration of precipitation for each month of the year. INTER24 truncates DURATN to the next lower even integer, and then centers the "storm" on noon. Thus if DURATN is input as 7.5, the daily precipitation is assumed to occur at a constant rate from time 0900 to 1500. Centering on noon is only used to see how much interception carries over into the next day. The algorithm for each hourly loop is the same as for INTER, except that rates are in mm/hr and amounts are summed over the day. The interception catch rate (RINT or SINT), and the evaporation rate (IRVP or ISVP) are returned to MSBPREINT as average rates over the day.

To determine appropriate values of DURATN I examined hourly precipitation data for 4 years (one year at Hubbard Brook) from the SAMSON data set. Averaging the number of hours per day of precipitation of 0.02 inch (0.5 mm) or greater over days with such precipitation gave the following results after a little smoothing

                    J   F   M   A   M   J   J   A   S   O   N   D
San Juan PR         3   2   2   2   2   2   2   3   3   3   3   3 
Atlanta GA          5   5   5   5   4   4   3   3   4   4   5   6
Caribou ME          4   4   5   5   4   4   4   4   4   6   6   5
Madison WI          4   4   5   3   3   2   3   3   4   4   5   5
Lake Charles LA     5   4   3   3   3   3   2   2   3   3   4   5
Phoenix AZ          4   4   4   4   4   2   2   2   2   2   4   4
Rapid City SD       3   3   3   4   4   3   2   2   2   2   4   4
Tacoma WA           6   6   5   4   4   4   4   4   4   4   6   6
Fairbanks AK        3   3   4   4   4   4   3   3   4   4   4   3
Hubbard Brook NH    5   5   5   4   4   4   4   4   4   5   5   5

Apparently a default DURATN of 4 hours is appropriate for all months anywhere in the U.S.

Subroutine PLNTRES

Subroutine PLNTRES is called at the beginning of each day to obtain resistivities to liquid water flow: rhizosphere resistivity for each soil layer, root resistivity in each soil layer, and xylem resistivity. These parameters, together with soil water potential in each layer (PSITI) and critical plant water potential (PSICR) control the supply of water to transpiring leaves and thus the reduction of actual transpiration below potential transpiration. As defined by Hunt et al. (1991) the resistances used here are "potential difference resistivities", because the transpiration flux rate is in units of mm/d and the potential gradient is in MPa. The resistivities have units of MPa d mm^-1.

In most of this subsection and the following subsection (TBYLAYER), algebraic notation is used instead of variable names. The correspondence is:

The following additional parameters or constants are input to the routines

Several other algebraic variables occur in the derivations below, but are not needed in the program.

If soil water potential is uniform through the root zone and rhizosphere resistivity is negligible, a bulk plant resistivity, r_p, can be defined by

r_p = ( ψ_t - ψ - ρ_w gd ) / T'

where T' is the transpiration rate, ψ_t is the soil water potential, and ψ is the leaf water potential. The ρ_wgd term is the gravity potential difference above the ground surface, where ρ_wg is the density of water times the acceleration of gravity, and d is the effective height of canopy evaporation, taken here as the closed canopy zero-plane displacement (DISPC). Change in water storage within the plant is ignored in BROOK90, as Hunt et al. (1991) have shown that it does not matter to total daily transpiration.

In BROOK90, r_p is determined primarily by the maximum bulk plant conductivity (MXKPL), which is an input parameter. MXKPL is the water uptake rate for a closed canopy per unit of soil to leaf potential difference. When the soil is wet so ψ_t is close to 0, many plants have a ψ of around -1.5 MPa when transpiration rate is about 0.5 mm/hr, a normal midday rate for a sunny day in temperate regions. This gives a MXKPL of 8 mm d^-1 MPa^-1. Abdul-Jabbar et al. (1984) found literature values to range only from 7 to 30 mm d^-1 MPa^-1, which seems a surprisingly narrow range for plant canopies of any species. Actual plant resistivity, r_p, is determined in subroutine CANOPY.

Fig. EVP-1. Resistances and potentials in the liquid flow path for transpiration, for 3 layers with roots. ψ is the leaf potential, ψ_x is potential at gorund level, and ψ_ti are total soil water potentials. r_s is xylem resistance, r_ri are root resistances, and r_si are rhizosphere resistances.

Figure EVP-1 shows the resistance network used in BROOK90. Each soil layer is considered to have a rhizosphere resistivity, r_si, and a root resistivity, r_ri, in series. Each layer is considered in parallel with the others. An additional resistivity in series accounts for resistance to flow through the xylem above ground level, r_x. The total soil-water potential, ψ_ti, differs among layers. The leaf water potential, ψ, is assumed constant through the canopy. The potential at the ground surface is ψ_x. This system and the parameterization of it was developed by Federer (1979) and Hunt and others(1991).

The xylem resistance, r_x, is

r_x = f_x r_p ;

where f_x (FXYLEM) is an input parameter, which is probably close to zero for short canopies, but may be 0.5 or higher for forests (Hunt et al. 1991). The use of f_x here specifies the amount of r_p (RPLANT) that should not be divided among soil layers, so it is effectively the fraction of the plant resistance that is above ground level.

The plant resistance that is below ground level, r_p - r_x, is the parallel combination of the individual layer resistances, r_ri (Fig. EVP-1). BROOK90 assumes that the root resistivity per unit length of root is constant, so r_ri (RROOTI) is inversely proportional to the fraction of total root length that is in each layer, f_i

r_ri = ( r_p - r_x ) / f_{i .}

BROOK90 parameterizes root distribution by the relative density of roots in each soil layer d_i (RELDEN_i). RELDEN_i is obtained from ROOTDEN in subroutine RTDEN . Then f_i is

f_i = d_i D_i / S ( d_i D_i )

where D_i is the stone-free layer thickness, THICK(i) * (1 - STONEF(i)).

When the soil is at uniform water potential and rhizosphere resistances, r_si , are negligible and f_x = 0, then the relative transpiration withdrawal from each layer is proportional to f_i (see subroutine TBYLAYER). Increasing f_x makes the withdrawal less dependent on f_i .

Usually the rhizosphere resistance only becomes significant when the soil is dry. Following Federer (1979) and Cowan (1965), the rhizosphere resistance, r_si , calculated in subroutine TBYLAYER is r_si = α_i / K_i , where K_i is the hydraulic conductivity of the rhizosphere. The value of α_i (ALPHA_i) for a layer is

α_i = A_i / ρ_wg D_i

where A_i from Cowan (1965) is

A_i = ( 1 / 8 π L_i ) [ δ_i - 3 - 2 ( ln δ_i ) / ( 1 - δ_i ) ]

where L_i is the root density in the layer (mm/mm³) and δ_i is the root volume fraction in the layer, obtained as

δ_i = π R₁² L_i

where R₁ is the average radius of the absorbing roots (RTRAD) , which is an input parameter. (Note that δ_i has no relation to the daylength, δ ). The root density is

L_i = f_i L_r / D_i

where L_r is the total length of absorbing roots per unit ground area in mm/mm². The value of L_r in m/m² is found in subroutine CANOPY as the parameter MXRTLN reduced by DENSEF and RELHT. The dependence on the seasonal RELHT assumes that root length increases proportionally with height growth. MXRTLN is an input parameter expressing the length of absorbing roots per unit ground area in m/m². L_r and R₁ only affect rhizosphere resistance and thus are only important for dry soil or when L_r is small.

Subroutine TBYLAYER

TBYLAYER determines the rate at which liquid water can be supplied to transpiring leaves, compares this rate with the potential transpiration rate, sets actual transpiration equal to the lesser of the two, and then allocates the transpiration among soil layers. This routine is based on the model of Federer (1979), which has been widely used, e.g. Wetzel and Chang (1987), Levine and Salvucci (1999).

The routine requires iteration when outflow from roots is prevented for layers in which xylem potential ψ_x is greater than soil water potential ψ_ti(Fig. EVP-1).

The resistance to uptake from each layer, r_i, is

r_i = r_ri + r_si = r_ri + α_i / K_i

(Fig. EVP-1) where r_ri (RROOTI) and α_i (ALPHAi) are variables from subroutine PLNTRES and K_i is the rhizosphere conductivity of the layer (Cowan 1965, Federer 1979). Estimating rhizosphere conductivity requires an iterative solution as described in Federer (1979), so BROOK90 uses the K_i (KK_i) of the bulk soil instead. This tends to underestimate the rhizosphere resistance, but the error is probably unimportant unless r_p is quite small. The total root resistance to uptake is

r_t = 1 / Σ ( 1 / r_i ) .

The transpiration rate for each layer, T_i , is

T_i = ( ψ_ti - ψ_x ) / r_i

where ψ_ti is the total water potential for the layer (PSITI) , and ψ_x is the xylem potential. The total transpiration from all layers is

T = ΣT_i = Σ( ψ_ti / r_i ) - ψ_x Σ ( 1 / r_i ) = ( ψ_t - ψ_x ) / r_t

where ψ_t is a weighted mean soil water potential defined as

ψ_t = ψ_t Σ (ψ_ti / r_i ) .

Eliminating ψ_x from the T_i and T equations gives

T_i = ( ψ_ti - ψ_t + r_t T ) / r_i

which is the equation used to distribute T among layers.

In BROOK90 actual transpiration rate, T, is the lesser of potential transpiration rate, P, and the maximum possible rate of water uptake from the soil, S. Following Federer (1979) and Lynn and Carlson (1990), BROOK90 assumes that the S obtains when the leaf potential, ψ, is the critical potential at which stomates close, ψ_c , which is an input parameter (PSICR). The assumption that the T is the lesser of the P and S is equivalent to an assumption that stomatal opening is not limited by leaf water potential, ψ, until the critical potential, ψ_c, is reached. The stomata are then assumed to close as much as necessary to maintain ψ = ψ_c. This abrupt switch from a demand-limited to supply-limited regime is oversimplified, but is convenient for modeling. This type of behavior induces a flat-topped diurnal transpiration curve (Lynn and Carlson 1990). The critical potential, ψ_c, represents the maximum suction that the plant can exert to get water from the soil, the minimum soil-water potential that the plant can induce, and the lower limit of soil water availability.

The water supply rate, S, is the potential difference between leaf water potential, ψ, and the xylem potential, ψ_x, divided by the xylem resistance, r_x , when ψ is equal to the critical potential, ψ_c (see Fig EVP-1). With allowance for the gravity potential between the ground and the canopy zero-plane displacement height, d (DISP), this is

S = ( ψ_x - ψ_c - ρ_w g d ) / r_x

S below ground must be T = ( ψ_t - ψ_x ) / r_t from above, and equating the two S equations to eliminating ψ_x, gives

S = ( ψ_t - ψ_c - ρ_w g d) / ( r_t + r_x )

which is assumed to be constant throughout a day.

Fig. EVP-2. Transpiration for the day, T' (shaded area), as the lesser of a constant water supply rate, S, and a half-sine potential transpiration, P'.

Following Federer (1982), BROOK90 assumes that the daytime potential transpiration rate, P', varies as a half sine wave. The actual transpiration rate, T' , is the lesser of S and P' (Fig. EVP-2). The average value of T' over the daytime, T, is

T = P ( 1 + R cos^-1 R - sin (cos^-1 R) )  R < 1

T = P       R >= 1

where

R = 2 S / π P

where P is the average of P' over the daytime. (Notation here differs slightly from Federer (1982) who used S' for S, D/d for P, and T/d for T where d is daylength.) At night, P is assumed constant and T is the lesser of S and P.

Normally, all layers with roots are included in the above calculations. However, when some layers are wet and others are dry, it is possible for ψ_ti (PSITI) in one or more layers to be lower than ψ_x so that the roots in those layers are releasing water to the soil. The question of whether such outflow occurs or is somehow prevented by the roots is controversial (Hunt et al. 1991). More recent work shows that outflow does occur (Dawson 1993). Richards and Caldwell (1987) name the process "hydraulic lift" because it moves water upward from wetter deeper soil layers through the plant to shallow drier layers. In BROOK90 outflow from roots is prevented when the parameter NOOUTF is set to 1, and is allowed if NOOUTF is 0. In general, BROOK90 usually will work better when NOOUTF = 1; this is the default. When NOOUTF = 1 and any T_i is negative, the layer with the most negative T_i is eliminated and new values of r_t , ψ_t , T, and T_i are obtained. If any T_i are still negative, the elimination process is repeated. This elimination procedure causes transpiration from a layer to cease when its potential is still greater than PSICR.

COMPASS BROOK