What is the Priestley–Taylor wet-surface evaporation parameter? Testing four hypotheses
This study compares four different hypotheses regarding the nature of the Priestley–Taylor parameter α. They are as follows:
α is a universal constant.
The Bowen ratio (, where H is the sensible heat flux, and LE is the latent heat flux) for equilibrium (i.e., saturated air column near the surface) evaporation is a constant times the Bowen ratio at minimal advection (Andreas et al., 2013).
Minimal advection over a wet surface corresponds to a particular relative humidity value.
α is a constant fraction of the difference from the minimum value of 1 to the maximum value of α proposed by Priestley and Taylor (1972).
Formulas for α are developed for the last three hypotheses. Weather, radiation, and surface energy flux data from 171 FLUXNET eddy covariance stations were used. The condition 0.90 was taken as the criterion for nearly saturated conditions (where LEref is the reference, and LEp is the apparent potential evaporation rate from the equation by Penman, 1948). Daily and monthly average data from the sites were obtained. All formulations for α include one model parameter which is optimized such that the root mean square error of the target variable was minimized. For each model, separate optimizations were done for predictions of the target variables α, wet-surface evaporation (α multiplied by equilibrium evaporation rate) and actual evaporation (the latter using a highly successful version of the complementary relationship of evaporation). Overall, the second and fourth hypotheses received the best support from the data.
On a globe dominated by ocean surfaces, wet-surface evaporation has obvious global importance (e.g., Brutsaert, 2023, p. 142; Andreas, 2013; McMahon et al., 2013; Szilagyi et al., 2014; Yang and Roderick, 2019; Tu et al., 2022). But estimates of wet-surface evaporation can be valuable over land surfaces as well. For example, the Global Land Evaporation Amsterdam Model (GLEAM) evaporation product (Miralles et al., 2011; Martens et al., 2017) uses the Priestley–Taylor (1972) wet-surface evaporation equation as the starting point for land surface evaporation. While climatic influences on wet-surface evaporation rates can differ from those on transpiration (e.g., Schymanski and Or, 2015), in the GLEAM model, adjustments for water stress are made using a “multiplicative stress factor” (Martens et al., 2017). Many other models and data products use some form of the Penman (1948) or the related Penman–Monteith equation (Monteith, 1965; Allen et al., 1998; McMahon et al., 2013). The advection–aridity version (Brutsaert and Stricker, 1979) of Bouchet's (1963) complementary relationship (CR) between actual and apparent potential evaporation (Brutsaert, 2005, p. 136) makes use of both the Priestley–Taylor and the Penman equation to estimate actual land surface evaporation from non-saturated surfaces.
Both the Penman (1948) and the Priestley–Taylor (1972) equation are estimates of potential evaporation, the hypothetical evaporation rate that one would get from a land surface if the surface were saturated (Brutsaert; 2005, p. 136). The Penman equation consists of a radiative term and an advective term. Slatyer and McIlroy (1961, p. 3–73) noted that the advective term would be zero if the surface and the lower atmosphere were fully saturated (see further discussion in Sect. 2.1). Evaporation under this condition is known as equilibrium evaporation. The Priestley–Taylor equation multiplies the equilibrium evaporation rate by a factor α (where α>1) to account for the presence of some vapor pressure deficit, even under conditions of minimal advection – what Priestley and Taylor (1972) termed “the absence of advection”.
Some work regarding α for wet surfaces treated it as a global constant to be found through field experiments (e.g., Priestley and Taylor, 1972). Other work made use of mixed-layer models of the atmospheric boundary layer (ABL), linked to a surface layer model in order to assess the role of ABL development in the value of α (de Bruin, 1983; McNaughton, 1976; McNaughton and Spriggs, 1989; Lhomme, 1997a, b; Raupach, 2000). As a whole, this work suggests that ABL processes result in variability of the value of α. Since Priestley and Taylor (1972) found a single central value for α, such variability could cast doubt on the concept of a minimal-advection, wet-surface evaporation rate that can be reliably estimated. Eichinger et al. (1996) derived an explicit equation for α. Szilagyi et al. (2014) used global ocean data products to investigate α over the oceans of the world. They found a discernible relationship between α and temperature (see also Yang and Roderick, 2019). Han et al. (2021) began with the sigmoid generalized complementary equation developed by Han and Tian (2018) and modified it to apply to wet-surface evaporation.
The context of this work is the estimation of wet-surface evaporation, whether under actually wet conditions or the hypothetical evaporation rate if an unsaturated surface were actually saturated (Thornthwaite, 1948). Such hypothetical wet-surface evaporation estimates are commonly used, particularly in models based on the CR (e.g., Brutsaert and Stricker, 1979; Szilagyi and Joszsa, 2008; Crago et al., 2016; Han and Tian, 2018, 2020). In this context, it is not immediately apparent whether formulae derived for wet surfaces will also provide good results (for hypothetical wet-surface evaporation) from unsaturated surfaces. A different context is sometimes seen in the literature, in which α is essentially a moisture-availability factor (e.g., De Bruin, 1983). This latter use of α will not be further considered here.
The objective here is to gain a better conceptual understanding of α (cf., Crago and Qualls, 2013). Similar to the limitations applied by Priestley and Taylor (1972; cf. Andreas et al., 2013), only cases where both sensible and latent heat fluxes are positive will be considered.
Four hypotheses regarding α will be examined:
Hypothesis 1. The ratio (α) between wet-surface evaporation under minimally advective conditions and under equilibrium conditions (i.e., a saturated atmospheric column near the wet surface) is a global constant.
Hypothesis 2. There is a globally constant ratio (Andreas et al., 2013; Yang and Roderick, 2019) between (1) the Bowen ratio that occurs under minimal advection with a saturated surface and (2) the Bowen ratio that would occur under equilibrium evaporation conditions, and α can be derived using this constant ratio.
Hypothesis 3. There is a globally constant relative humidity value that can be used to derive an estimate of α that corresponds to minimally advective conditions
Hypothesis 4. The parameter α is a globally constant fraction of the gap between the minimum value of α and the maximum-allowable value of α proposed by Priestley and Taylor (1972).
All four hypotheses will be examined under actually saturated conditions, but they will also be evaluated under unsaturated conditions. Because there is no measured or reference value for the hypothetical wet-surface evaporation rate, instead the hypothetical rate will be included in a well-tested CR model for actual evaporation. That is, the CR model accuracy will be taken as an indirect measure of the method's ability to estimate the hypothetical wet-surface evaporation rate.
2.1 Wet-surface evaporation equations
Latent heat flux LE (W m−2) is related to evaporation rate E (kg m−2 s−1) as LE =lvE, where lv is the latent heat of evaporation (J kg−1). Penman's (1948) equation for apparent potential evaporation (LEp) from a wet surface can be written as follows:
where (Pa K−1) is evaluated at the air temperature Ta (K) at height zT (m), e (Pa) is vapor pressure, and e* (Pa) is saturated vapor pressure, where both e and e* are calculated using the formulations given by Andreas et al. (2013) which are valid for temperatures both above and below freezing. The net radiation is Rn (W m−2); G (W m−2) is the ground heat flux; the latent heat of evaporation, lv, is also calculated with a formulation given by Andreas (2013); ) (Pa K−1) is the psychrometric constant; p (Pa) is atmospheric pressure; and cp is the specific heat of air at constant pressure (J kg−1 K−1). The formulae adapted from Andreas et al. (2013) have been included in the Supplement. The drying power of the air EA (kg m−2 s−1) is defined by
where ea is the vapor pressure at height zT and f(u) (s m−1) is a function of wind speed. The wind function can be calculated (Brutsaert, 2015) using Monin–Obukhov similarity theory (MOS theory):
where k=0.4 (dimensionless) is von Karman's constant; Rd (J kg−1 K−1) is the ideal gas constant of dry air; u (m s−1) is wind speed measured at height zT; d0 (m) is the displacement height; and z0 (m) and z0v (m) are the roughness lengths for momentum and sensible heat, respectively. Equation (3) is based on MOS theory, the standard formulation of flux–gradient relationships in the lower atmosphere (Stull, 1988, p. 376; Brutsaert, 2005, p. 128), but Penman (1948) recommended the form f(u)=c1(1+c2u2), where c1 and c2 are empirical constants, and u2 (m s−1) is wind speed at 2 m (m s−1). This latter formulation (not used here) is preferred by some authors (e.g., Szilagyi et al., 2019) because information about the roughness of the surface (needed for z0, z0v and d0) is not needed.
Note that other versions of LEp are available, including one (e.g., Qualls and Crago, 2020; Crago and Qualls, 2021) which is based on the surface energy budget with mass and energy transport functions for the latent and sensible heat fluxes, respectively. While Eqs. (1)–(3) are based on the same principles, Penman's (1948) derivation involved his well-known approximation that Δ for a wet surface is approximately equal to the ratio of the difference in vapor pressure between the surface and measurement height to the difference in temperature between the same two levels, which allowed the simple two-term Eq. (1). Only Eq. (1) will be used for LEp in this project.
As described by Brutsaert (2005, p. 129), air in the lowest layers blowing for a long distance over a wet surface would likely become increasingly humid. If it should approach saturation, the second term of (1) would go to zero, leaving the first term of (1) as an effective “lower limit” or “equilibrium value” for wet-surface evaporation (Slatyer and McIlroy, 1961, p. 3–73). This evaporation rate is often termed equilibrium evaporation (e.g., Brutsaert, 2005, p. 129). Here, this lower limit LEe (W m−2) is calculated as
where is Δ evaluated at the wet-surface temperature T0 (to be defined shortly). While Δ is commonly estimated at Ta, (e.g., Brutsaert, 2005, p. 126), Eq. (4) corresponds to the definition of equilibrium evaporation suggested by Andreas et al. (2014) and Qualls and Crago (2020). Namely, it is the lowest wet-surface evaporation rate possible for a given available energy value (Rn−G) with a surface temperature of T0 (K). It is a minimum because lower evaporation rates would require the vapor pressure to exceed the saturation value (e.g., Philip, 1987; Andreas et al., 2013; Qualls and Crago, 2020). The fact that super-saturation cannot occur during evaporation explains why wet-surface evaporation is limited by Eq. (4) rather than simply by (Rn−G) (see Qualls and Crago, 2020). The Bowen ratio (Bo ; dimensionless) corresponding to Eq. (4) is Bo (.
Priestley and Taylor (1972) introduced the parameter (dimensionless) so that
where LEPT (W m−2) estimates minimum-advection wet-surface latent heat flux. Because of some dry advection even over extensive saturated surfaces, they found α>1. Their data suggested a typical value of α≈ 1.26. Because LEPT≥ LEe and H≥ 0, the limits on α are (Priestley and Taylor, 1972)
Hypothesis 1 suggests that Eq. (5), with α a global constant, defines minimal-advection wet-surface evaporation.
Andreas et al. (2013) examined thousands of measurements taken over extensive water and ice surfaces for which H>0 and LE >0 and suggested that Bo is related to Bo* by
where aA (dimensionless) was found to be a global constant of about 0.4. This is equivalent to a Priestley–Taylor α of
In Eq. (8), aA is a constant, and Δ is a function of the skin temperature T0 so that a discernible relationship between α and T0 is implied by Eq. (8) (cf., Szilagyi et al., 2014). Hypothesis 2 suggests that Eq. (8) captures the foundational concept of α.
Yang and Roderick (2019) made a similar proposal to Eq. (7), resulting in aA=0.24 based on global ocean data products. However, they noted that, in practice, LE and Rn cannot be known independently of each other over oceans, since increased LE reduces the ocean surface skin temperature, which reduces outgoing longwave radiant fluxes, thereby increasing Rn. Their value of aA accounts for adjustments in the available energy resulting from this linkage. The present study assumes that Rn−G is known via measurements at each site.
Eichinger et al. (1996) had already proposed a dimensionless variable ] for use in an explicit method (their Eq. 7) to estimate α for wet surfaces. Plans to include an additional hypothesis based on their Eq. (7) in this study were abandoned when it became apparent that their C (taken as a constant model parameter rather than calculated with the definition given in the previous sentence) is mathematically equivalent to (1−aA). While we will refer to Eq. (8) as the Andreas et al. (2013) formula, we acknowledge the prescient contribution of Eichinger et al. (1996).
As an alternative to Eq. (8), if there is minimum advection over a wet surface, both Eqs. (1) and (5) should give the correct evaporation rate. By setting them equal to each other, one arrives at
where RH (dimensionless) is the relative humidity of the air, and αRH is dimensionless. The values of lv, , and e* could all be evaluated at the wet-surface skin temperature T0. Equation (9) gives the correct value of α within the accuracy of Penman's (1948) assumption regarding Δ, provided RH is the measured relative humidity. However Eq. (9) is proposed here as a parameterization of α for both actually and hypothetically saturated surfaces, where RH is the model parameter representing the relative humidity under saturated surface and minimal advection conditions. Small values of (Rn−G) could result in unreasonably large values of α. Therefore, the limits given by Eq. (6) are applied to estimates of αRH. That is, if , it is set to , and if αRH<1, then it is set to 1.
The limits on α given by Eq. (6) suggest that perhaps α takes a constant intermediate position in between the limits. Thus, the parameter m (dimensionless) is
Hypothesis 4 suggests that Eq. (11) is the best explanation of α.
Szilagyi and Jozsa (2008; see also Szilagyi and Schepers, 2014; Szilagyi et al., 2017) suggested T0 could be found by setting two expressions for the Bowen ratio (here, given by ) equal to each other:
where the equation used for e*(T0) (from Andreas et al., 2014) is given in the Supplement. The wet-surface temperature in Eq. (12) is T0, which can be easily found from Eq. (12) with a numerical root finder. Equation (12), thus solved, provides the wet-surface temperature T0 from data taken from either saturated or unsaturated surfaces (Szilagyi and Schepers, 2014).
2.2 The complementary relationship (CR) of evaporation
In the complementary relationship (CR) between actual and potential evaporation (Bouchet, 1963), regional evaporation from a saturated surface, the apparent-potential evaporation rate, and the actual evaporation rate are all identical (Brutsaert, 2015, p. 136). According to the advection–aridity approach (Brutsaert and Stricker, 1979), apparent potential evaporation corresponds to Penman's equation (Eq. 1) or to the evaporation from a small wet patch, and the wet regional surface rate corresponds to the Priestley –Taylor equation (1972), shown in Eq. (5). As the surface dries, less water is available to evaporate, so actual evaporation decreases. This results in a drier and warmer lower atmosphere, which increases apparent potential (wet patch) evaporation. Conversely, if the lower atmosphere becomes dry and warm (in the absence of significant dry advection), this implies that regional evaporation rates are low. Thus, evaporation and apparent potential evaporation change in opposite directions – they complement each other. An estimate of the Priestley–Taylor α is an integral part of most CR models, and the performance of CR models making use of the four different hypotheses regarding α can serve as a further test of the hypotheses. Note that Han et al. (2021) took a different approach, by adapting the CR model of Han and Tian (2018) to estimate evaporation from wet surfaces; this results in a non-linear dependence of wet-surface evaporation on equilibrium evaporation.
As formulated by Brutsaert (2015, p. 136), the CR can be formulated in terms of x= and , where LEw is given by Eq. (5) and LEp by Eqs. (1)–(3). Both x and y are dimensionless. Values of α to be used in Eq. (5) will be discussed in Sect. 3. Brutsaert (2015) used physical reasoning to suggest that at x= 0, the boundary conditions are y= 0 and d, while at x= 1, they are y= 1 and d 1. Crago et al. (2016), however, noted that y can approach zero when there is no water available to evaporate, but x cannot approach zero unless LEw goes to zero. The smallest x can get is
where (W m−2) is given by
In Eq. (15), the subscript “d” means the variable is evaluated at Td, the “dry air temperature”. A straight line with slope d (where Tg is a generic temperature variable) represents an isenthalp (line of constant available energy) through (Ta, ea) on a graph of temperature (x axis) versus vapor pressure (e on the y axis). The temperature at which this isenthalp reaches e= 0 is Td (Szilagyi et al., 2017; Crago and Qualls, 2021). That is,
Crago et al. (2016) suggested x could be “rescaled” through the following transform:
A simple formulation suggested by Crago et al. (2016) is
Crago et al. (2022) considered data from seven FLUXNET sites in Australia as well as global, gridded ERA5 data (Hersbach, 2020) produced by ECMWF (European Centre for Medium-Range Weather Forecasts; https://www.ecmwf.int/, last access: 23 August 2023). With the FLUXNET data, Eq. (17) consistently performed best at predicting reference (eddy covariance) latent heat fluxes. Since FLUXNET data are used here as well, Eq. (17) will be the assumed CR formula for this study. From Eq. (17), latent heat flux estimates can be found with LEest=X(LEp).
Equations (1)–(5) and CR methods are generally considered applicable at timescales ranging from daily to monthly, with monthly being most common (McMahon et al., 2013). Equations (1)–(5) require homogeneous surfaces corresponding to the spatial extent of the flux footprint (e.g., Schuepp, 1990), typically corresponding to several hundred meters, while CR formulations are best suited for homogeneous conditions at the “regional” scale (Brutsaert, 2005, p. 136) of perhaps tens of kilometers.
3.1 Data sources and processing
Monthly- and daily-average data pre-processed by FLUXNET were downloaded as CSV files from the fluxnet.org website for 171 eddy covariance stations (listed in Table S1 in the Supplement). At least minimally adequate fetches are assumed at all sites included in FLUXNET. Measurement heights, latitudes, longitudes, IGBP (International Geosphere-Biosphere Programme) land classes (https://www.igbp.net, last access: 23 August 2023), and canopy heights were provided for these sites by Wang et al. (2020; see their supporting information). Wang et al. (2020) assumed z0=0.123hc, d0=0.67hc, and z0v=0.1z0, where hc (m) is the reported canopy height. These values were all adopted herein. Separate wind speed, temperature, and humidity measurement heights were not included by Wang et al. (2020), so it is assumed here that they are all measured at the single given height. Net radiation, ground heat flux, sensible and latent heat fluxes, air pressure, air temperature, vapor pressure deficit, and wind speed were included in the FLUXNET downloads. All variables employed some gap filling using the MDS (marginal distribution sampling; Reichstein, 2005) method as described by Pastorello et al. (2020). Data flagging and quality assurance and control for all the variables also followed the procedures outlined by Pastorello et al. (2020). The data and the software used in this project are available through Zenodo and GitHub (see Crago et al., 2023a, b).
Following the procedures outlined by Pastorello et al. (2020), the half-hourly or hourly energy fluxes are also gap-filled using the MDS (Reichstein, 2005) method, and these are used by FLUXNET to derive the daily or monthly reference values used here. The FLUXNET dataset includes the variables “H_CORR” and “LE_CORR”, which indicate corrected values, that is, values that correspond to energy budget closure. However, the corresponding uncorrected variables (“LE_F_MDS” and “H_F_MDS”) are available for more sites and times. These latter surface fluxes were used for this study for the reference values LEref and Href, respectively. Issues regarding energy budget closure with eddy covariance fluxes are complicated, as discussed by Mauder et al. (2020). In this study, the downloaded values of sensible heat flux are taken to be the final reference values Href, while downloaded latent heat flux values are adjusted so that monthly (daily) energy budget closure is obtained: . This was the procedure recommended by Wang et al. (2020; see also Tu et al., 2022, 2023).
Months (or days) with eddy covariance values of H and LE less than zero or were screened out of the dataset; this eliminates periods of strong dry advection that result in negative Href. The ground heat flux G was not measured at all for some of the sites, and missing values of G also occurred. When measurements of G were not available, a value of zero was assumed. Over a 24 h period, G “is often near zero” (Stull, 1988), so this assumption is not unreasonable. Over the daily to monthly timescales at which Eqs. (1) and (5) are commonly used (McMahon, 2014), the assumption likely improves as the averaging time increases.
The CR has been used at timescales from hours to years (Brutsaert, 2023, p. 147), but the CR, the Penman equation (Eq. 1), and the Priestley–Taylor equation (Eq. 5) most typically use daily- to monthly-average values (McMahon, 2013). It is true that use of time averages of variables as inputs to non-linear equations can lead to “significant errors” (Slatyer and McIlroy, 1961, p. 3–58). However, CR and the Priestley–Taylor wet-surface equation both assume that the land surface conditions and the temperature and humidity in the lower atmosphere are well adjusted to each other (Brutsaert, 2023, p. 147). The diurnal cycle makes this adjustment unlikely over periods less than 24 h (McMahon, 2013). Therefore, the approach here is to use daily-average (monthly-average) input values to produce daily (monthly) energy fluxes (e.g., Penman, 1948; McMahon, 2013; Brutsaert, 2023). That is, daily to monthly timescales are suited to these equations, as spatial scales corresponding to small watersheds are suited to saturation-excess runoff (e.g., Chow et al., 1988).
With this dataset, the monthly (daily) mean of reference latent heat flux is 61 W m−2 (62 W m−2), the median is 58 W m−2 (56 W m−2), and the standard deviation is 41 W m−2 (48 W m−2). Thus, the central tendencies for monthly and daily values are similar, but the daily values have more spread about the mean.
The decisions described above regarding data inclusion, screening, and correction reflect a desire to obtain a broad range of climates, land covers, and seasons, so as to test the four hypotheses under as wide a range of conditions as possible. While these decisions do entail some risk of including lower quality data, we think they are defendable as outlined above. However, much of the analysis was repeated after removing data for which the FLUXNET quality control index (which ranges from 0 for very poor to 1 for excellent) was less than 0.9, with little difference in numerical results and no difference in qualitative results (such as which methods performed better than other methods) compared to the results without filtering for data quality.
A total of 11 different IGBP land surface classes (e.g., Loveland et al., 1999) are in the original data. They include Wooded Savanna (WSA), Grassland (GRA), Evergreen Broadleaf Forest (EBF), Cropland (CRO), Evergreen Needleleaf Forest (ENF), Savanna (SAV), Deciduous Broadleaf Forest (DBF), Closed Shrubland (CSH), Mixed Forest (MF), Open Shrubland (OSH), and Permanent Wetland (WET). Classes for each site were provided by Wang et al. (2020).
3.2 Estimates of wet-surface α
Determination of wet-surface values of α requires determination of data representing saturated surface conditions. Saturated or nearly saturated land surface conditions were assumed when LEref>T(LEp), where T (dimensionless) is a threshold value which had to be determined. A value of T was sought such that the linear regression between LEref and LEp (for data for which the condition above is met) falls nearly on the 1 : 1 line, and at the same time root mean square (rms) errors between LEref and LEp are very small. After filtering for wet-surface conditions, αref was calculated as .
Trial values of the parameters αc, aA, RH, and m, where αc is a constant (global) value of α, were selected randomly from a range of reasonable values (that is, α ranged from 1 to 1.6, and aA, RH, and m all varied from 0 to 1). These values were used at all sites and times, satisfying the wet-surface condition. A total of 3000 values drawn randomly from these ranges were evaluated to determine the optimal parameter values. These optimal parameter values were used to estimate different versions of αest, namely αc and aA from Eq. (8), αRH from Eq. (9), and αm from Eq. (11). These αest values were compared to αref= . The trial parameter values that minimized the root mean square difference (RMSD) between αest and αref were taken to the be tuned parameter values; the values of αest will be called the best-fit values.
3.3 Estimates of wet-surface evaporation
Randomly selected parameter values were chosen 3000 times (over the same range as in Sect. 3.2) to estimate αc, αA, αRH, and αm for use in Eq. (5) to estimate LEest=αestLEe, where αest is estimated from the parameter values using Eqs. (8), (9), and (11). This was done for all wet-surface measurements. Tuned parameter values were those giving minimum RMSD between the resulting LEest and LEref, and the resulting LEest values are considered best-fit values. The same tuned parameter values were used for all stations and all times. Note that the tuned parameter values found in this analysis may differ from those found in Sect. 3.2.
3.4 Complementary relationship for actual evaporation
Next, actual evaporation was estimated with the CR using all available data. That is, the analysis was not limited to only wet-surface data. Specifically, 3000 new samples of the parameter values were chosen from the same range as above. Those parameter values were substituted into Eqs. (8), (9), or (11) and were used to calculate estimates of α. Those values in turn were used in Eqs. (1)–(5) along with Eqs. (14)–(16) in Eq. (17) and finally in . The RMSD between LEest and LEref was found for each of the four methods of estimating LEPT. The tuned parameter values were those minimizing the RMSD, and the corresponding LEref values are the best-fit values.
For monthly averaging times, the wet-surface threshold was established to be T=0.90, which resulted in a regression equation with slope of about 1 and intercept very near 0, while rms errors were small (Fig. 1d). The process by which this value of T was established is described using Fig. 1, and the corresponding statistics for both monthly and daily data are found in Table 1. The second row of text in each panel identifies the range of values incorporated into the graph. If the lower limit of accepted values was T, then in the left column (Fig. 1a, c, e) an upper limit of 2−T was imposed. In the right column (Fig. 1b, d, f) no upper limit was imposed (as indicated by an upper limit denoted “−”). There seems to be no compelling reason to impose an upper limit on wet-surface , even though the upper limit improved many of the statistics (i.e., comparing Fig. 1a to b, c to d, and e to f). Figure 1b has a slope somewhat below 1a, and Fig. 1f has a large RMSD. Figure 1d, where all points with were included, was taken as a reasonable compromise. As shown in Fig. 1d, wet-surface evaporation defined in this way occurred in 430 months and from 50 of the sites. The sites included IGBP classes CRO, ENF, GRA, DBF, WET, OSH, and EBF. For daily averaging times, a similar process was followed (statistics included in Table 1). Wet-surface evaporation (defined again by T=0.90) occurred on 22 998 d, involving 158 of the 171 stations, and included IGBP classes WSA, GRA, EBF, CRO, ENF, SAV, DBF, CSH, OSH, WET, and MF. Figure 2 shows the location of sites having at least 1 month of wet-surface measurements (top panel) and those having none (bottom panel).
* R: correlation coefficient. Slope and intercept: S and I in . NSE: Nash–Sutcliffe (1970) efficiency. RMSD: root mean square difference. Range values ending with “–” indicate no upper limit to the range.
As described in Sect. 3.2, using only wet-surface measurements, tuned values of αc, aA, RH, and m were found that minimized RMSD between reference and estimated values of α. Different tuned values of these parameters were also found for use in estimating wet-surface evaporation using the Priestley–Taylor equation (Eq. 5). Results for estimating α under wet-surface conditions are found in Fig. 3a, and results for estimating wet-surface evaporation itself are shown in Fig. 3b; both sets of results are also provided in Table 2. Finally, still-different tuned values of αc, aA, RH, and m were those which produced the minimum values of RMSD between LEest and LEref using the CR formulation in Eq. (17). That is, LEest is found by taking y found with Eq. (17), where LEw is given by Eq. (5), LEp by Eqs. (1)–(3), and y by Eqs. (13)–(17). This y was multiplied by LEp from Eq. (1) to get LEest. The tuned parameters that result when the goal is to obtain the best fit between αest and αref are different than those when the goal is to obtain the best fit between wet-surface LEPT and LEref, and in turn these are different than when the goal is to obtain the best fit between LEest and LEref using Eq. (17). Reasons and implications for these differences will be discussed in Sect. 5.
a RMSD: root mean square difference. R: correlation coefficient. Slope and intercept: S and I in . NSE: Nash–Sutcliffe (1970)
b RMSD(αest, αref) is the RMSD between α estimates found with the prescribed method and the reference values found from αref= for times with wet-surface evaporation conditions.
c RMSD(LEw_est, LEref) is the RMSD between LEest estimates found with the prescribed method and the reference values LEref for times with wet-surface evaporation conditions.
d RMSD(LEest, LEref) is the RMSD of all wetness conditions between estimates LEest and LEref.
5.1 General trends
For convenience, the use of α with a single global value will be called the “αc method” (corresponding to Hypothesis 1), with (8) it will be the “aA method” (Hypothesis 2), with Eq. (9) it will be the “RH method” (Hypothesis 3), and with Eq. (11) it will be the “m method” (Hypothesis 4). After discussing trends found in the results, the four hypotheses will be evaluated based on the results.
Figure 1 shows that LEref and LEp (Eq. 1) are very similar when the threshold for wet surfaces at a monthly timescale is taken to be T= 0.90. For daily data the same threshold provides good results, with low rms error and linear regression close to the 1 : 1 line (Table 1), so was chosen as the indicator of wet-surface evaporation, with no upper threshold to . Figure 2 shows the geographical location of sites that had some wet-surface evaporation months (upper panel) and sites that did not (bottom panel).
Figure 3 gives the results from the four methods in terms of prediction of αref itself (Fig. 3a) and of wet-surface evaporation estimates from Eq. (5) (Fig. 3b). Results from the use of the four methods when used in the CR model y=X are shown in Fig. 4. Even when estimates of α differ considerably from the reference values (Fig. 3a), the methods still provide good estimates of wet-surface evaporation (LEw – Fig. 3b) and actual evaporation (Fig. 4).
Table 2 provides much of the same data as Figs. 3 and 4 (for monthly averaging), while Table 3 provides the same information for daily averaging. A large number of FLUXNET sites spanning a wide range of climates and land cover classes were included in this study; such a diverse and large number of sites provides some confidence that the trends discussed here would apply to other sites and regions.
Different tuned parameter values were found to calculate αest (Fig. 3a), LEPT (Fig. 3b), and LEest (Fig. 4). Ideally, the tuned parameter values would remain nearly identical in the three cases. A likely explanation for this difference is as follows: in the top panel of Fig. 3, all αref values count equally in determination of tuned parameter values that produce the best-fit αest for each of the methods. But in the bottom panel of Fig. 3, αest values that correspond to small values of (Rn−G) have far less influence on the RMSD of LEPT than those corresponding to larger Rn−G. If 5 W m−2, an increase in α from 1.1 to 1.3 only increases α(Rn−G) from 5.5 to 6.5 W m−2, whereas if Rn−G is 200 W m−2, it increases from 220 to 260 W m−2; therefore the larger Rn−G would influence RMSD of LEPT more. Similarly, when moving to the CR estimate LEest in Fig. 4, the CR estimates again apply different weight to the various estimates of α; therefore different tuned parameter values result here as well. Actually, Brutsaert (2023, p. 149) treats the parameter α in Eq. (5) as a completely different parameter from the α embedded in Eq. (17). While the present authors consider both to be the same parameter, we recognize that its tuned value could vary depending on the context.
For reference, when the parameter values used in Fig. 3b are used in the CR (namely α= 1.29 for the αc method, αA= 0.31 for the αA method, RH = 0.76 for the RH method, and m= 0.58 for the m method), rms errors increased (from 19.37 to 20.12, from 18.61 to 19.25, from 21.13 to 26.84, and from 18.65 to 19.29, respectively). Note that the αA and m methods still provide the lowest RMSD values.
Figure 3 and Tables 2 and 3 provide evidence that all four methods provide acceptable estimates of actual wet-surface evaporation rates. But what about estimation of hypothetical wet-surface evaporation rates? Szilagyi and Jozsa (2008) and Szilagyi and Schepers (2014) have provided good evidence that a small wet patch within a drying region with wind speed and available energy held constant should maintain a constant surface temperature. The use of Eq. (12) to get the temperature of an actually saturated region is straightforward. Crago and Qualls (2021), Qualls and Crago (2020), and Szilagyi (2021) show graphically, using e–T graphs, how air and ground-surface isenthalps (lines of constant available energy on a (T, e) graph) can be determined, and they show that T0 is simply the intersection of the ground-surface isenthalp with the saturation vapor pressure curve.
a RMSD: root mean square difference. R: correlation coefficient. Slope and intercept: S and I in LEref+I. NSE: Nash–Sutcliffe (1970)
b RMSD(αest, αref) is the RMSD for between α estimates found with the prescribed method and the reference values found from αref= for times with wet-surface evaporation conditions.
c RMSD(LEw_est, LEref) is the RMSD for between LEest estimates found with the prescribed method and the reference values LEref for times with wet-surface evaporation conditions.
d RMSD(LEest, LEref) is the RMSD for all wetness conditions between estimates LEest and LEref.
There are three distinct explanations for why a best-fit estimate LEest using the CR (17) might differ from the LEref values: first, α may not be estimated correctly (see the preceding paragraphs); second, the LEPT estimate may not adequately represent the hypothetical wet-surface evaporation rate; and third, the CR formulation may be inadequate. No method to distinguish the effects of these is apparent; therefore the results in Fig. 4 and Tables 2 and 3 only provide an indirect test of the adequacy of the four methods to estimate α for hypothetical wet surfaces.
Under drying surface conditions, since the wet-surface temperature remains constant during drying (assuming Rn−G and wind speed are constant – see Szilagyi and Schepers, 2014), T0 found with Eq. (12) under drying conditions should still be the correct wet-surface temperature, that is, the temperature at which Δ in Eq. (5) should be evaluated to estimate the hypothetical wet-surface evaporation rate (but note that Δ in Eq. (1) is always taken at air temperature). During the regional drying process, Crago and Qualls (2021) showed that ea slides down the air isenthalp as drying progresses, while T0 is found just as it would be for a saturated surface, namely using Eq. (12). So, use of Eq. (12) to determine T0 for either saturated or unsaturated surfaces seems to have good support, and this T0 value can be used to predict wet-surface evaporation rates from Eq. (5) with Eq. (4). The process is described above as a temporal drying of the region, but the analysis of Crago and Qualls (2021) is only concerned with the current status of the land and lower atmosphere and not with the drying or wetting pathway to that status.
With respect to the best formulation for the CR, there is unfortunately no consensus (e.g., Crago et al., 2022; Han and Tian, 2018). However, with FLUXNET data from Australia, Crago et al. (2022) found that the y=X formulation was the best overall for predicting latent heat fluxes under many conditions. Given the wide range of methods represented in Figs. 3 and 4, it is actually striking how little variation there is among the various methods with respect to RMSD (Fig. 4; Tables 2 and 3).
Comparison of the four methods associated with hypotheses 1–4 suggests that all the hypotheses can give good estimates in many cases. In Sect. 5.2 we will compare the methods in the context of an examination of these hypotheses.
5.2 Examination of the hypotheses
As discussed in Sect. 1, the objective of this study is to evaluate different hypotheses or conceptualizations regarding α, by using them to estimate α itself, actual wet-surface evaporation (Eq. 5), and hypothetical wet-surface evaporation as a part of a CR model that predicts actual regional evaporation rates (Eq. 17). As discussed in Sect. 5.1, outcomes from these conceptualizations are used to evaluate the hypotheses stated in the Introduction. Including a range of hypotheses in this process makes it more likely that the correct conceptualization will be included and identified as the best.
Hypothesis 1, based on the αc method, has been the default hypothesis in the majority of work with Eq. (5) and within CR formulations (e.g., Brutsaert, 2015, 2023, p. 148; Crago et al., 2016, 2022). Growing evidence that wet-surface, minimal-advection α actually has a fairly wide range of values (e.g., McNaughton and Spriggs, 1989; Lhomme, 1997a, b; Raupach, 2000) might raise doubt regarding our ability to accurately estimate wet-surface evaporation. Clearly, unexplained variability is a real challenge, but the αc estimate performs quite well in predicting actual wet-surface evaporation and in the CR model (Figs. 3 and 4 and Tables 2 and 3).
A possible explanation for this surprisingly good performance begins with work by Szilagyi et al. (2014) and Andreas et al. (2013), who showed that much of the variability of α is due to temperature, with α increasing with decreasing temperature. Several formulations for α in terms of temperature are given in Fig. 5, including αA, αm, αmax (defined in the Fig. 5 caption), and a multi-term polynomial developed by Szilagyi et al. (2014) for α over saturated land surfaces. Because α is a function of temperature, it is likely that many of the large values of αref in Fig. 3a correspond to cold temperatures, which typically imply low available energy. Because available energy is small, relatively large errors in α result in only small absolute errors in wet-surface evaporation using Eq. (5). Thus, the fact that the global constant value of α is too small for these low-temperature sites does not result in large absolute errors in wet-surface evaporation rates. Nevertheless, this is clearly not the best supported of the four hypotheses.
Hypothesis 2, based on the αA method, assumes a constant ratio (aA) of the Bowen ratio between equilibrium and minimal-advection conditions. The resulting equation for α (Eq. 8) is able to account for much of the systematic variability of αref due to temperature variability because of the variable in Eq. (8). Figure 3a shows that αA estimates of α do increase as αref increases, but not as quickly as the reference values. While the trend is not matched perfectly, the αA method is clearly an improvement over the αc method in terms of predicting αref. The method also performs well at predicting actual wet-surface evaporation and actual evaporation (Figs. 3b and 4), and it provides the smallest RMSD for estimating actual evaporation from Eq. (17) for both monthly (Table 2) and daily (Table 3) data. With a clear definition and consistently good performance, Hypothesis 2 has considerable support. However, it is not obvious (based on physical principles) why aA in Eq. (7) ought to be a constant. Overall, Hypothesis 2 gains considerable support from the data presented here.
Hypothesis 3, based on Eq. (9), also captures much of the variability of αref. Equation (9) is correct to within the accuracy of Penman's (1948) well-known assumption regarding Δ, provided that RH is the actual measured relative humidity. When measured RH is replaced with the parameter RH, Eq. (9) provides an estimate for αref. Based on Eq. (1) combined with Eq. (5), Eq. (9) suggests that the optimal value of RH should ideally represent the relative humidity that characterizes wet-surface evaporation with minimal advection. Note that in Eq. (9), α depends on f(u), Rn−G, and temperature.
As seen in Tables 2 and 3, the RH method does not rank highly for prediction of α, wet-surface evaporation, or actual evaporation. Hypothesis 3 makes a very intuitive claim regarding wet-surface minimal-advection evaporation, namely, that it is associated with a particular value of relative humidity. While this method is conceptually appealing, and it performs relatively well with some subsets of the data (not shown), its performance in this study is not as good as that of hypotheses 2 and 4. Thus, this study does not provide much support for Hypothesis 3.
Hypothesis 4 assumes that minimal advection has been achieved when α is a specified fraction (m) of the distance from α= 1 to the maximum physically realistic value of (Priestley and Taylor, 1972). The idea of m being this fraction is clear and understandable, but it is not immediately obvious that it must be true on physical grounds. Overall, this method gives the lowest rms error for estimating α and LEw, it performs nearly as well as the αA method in estimating LE using the CR (Eq. 17), and it does this at both monthly and daily timescales. Furthermore, the other statistics included in Tables 2 and 3 are consistently favourable for this method. The data examined here seem to provide support for this method comparable to the αA method.
Note that another hypothesis was considered for inclusion, based on the αS curve, developed by Szilagyi et al. (2014) for saturated land surfaces and included in Fig. 5. The fact that α is a strong function of temperature is an important insight. However, the valid temperature range of their curve is more limited (from 0 to 28∘) than the temperatures in the dataset, and variability of α with temperature is already included in hypotheses 2, 3, and 4. Also, these hypotheses can be stated in terms of the parameters aA, RH, and m, respectively, which have well defined and physically meaningful definitions. Therefore, no fifth hypothesis was evaluated.
Different data points might be assumed to represent wet-surface conditions depending on the threshold value of T as illustrated in Fig. 1. But different points could also result for a given T value for different values of z0 used in the wind function (Eq. 3). The data presented here have used the Wang et al. (2020) z0 and d0 formulations as described above. But eddy covariance measurements of friction velocity u* (m s−1) are available for most of the sites and measurement periods. This means the logarithmic wind profile )log[(] (e.g., Brutsaert, 2023, p. 41) can be solved for z0 for each measurement period for which u* is available. This value of z0 is specific to a particular site and a particular month or day, so it accounts for roughness variations with season and wind direction. With these data, the values of z0 calculated in this way are somewhat smaller than those found with the Wang et al. (2020) formulation, which causes LEp values to be smaller and more data points to be identified as wet-surface values. Nevertheless, a figure similar to Fig. 1 but using these new z0 values (not shown) suggests that T= 0.9 is still appropriate. This method results in root mean square differences (see the Supplement, Tables S2 and S3) comparable to those in Tables 2 and 3. The number of data points differ because not all time periods had u* measurements and because different z0 values resulted in different data points qualifying as wet-surface values. However, the key points remain unchanged. Results from the Wang et al. (2020) version of z0 and d0 are shown herein because they represent the way the roughness of land surfaces is usually estimated.
Four hypotheses regarding the Priestley–Taylor (1972) parameter α were considered. Each of them has a different assumption regarding the nature and variability of α. In the first hypothesis α is constant; in the second it represents a ratio of two Bowen ratios; in the third, it represents conditions at a given relative humidity value; and in the last, it can be seen as a midpoint between theoretical maximum and minimum values. Using FLUXNET data from a total of 171 stations, α, LEPT, and actual evaporation values are compared to reference values in an attempt to determine which hypotheses best explain the data.
The second and fourth hypotheses generally produce the best results. In both of these, α is dependent on temperature, although the functional forms of the relationship are different. The third hypothesis has a very intuitive physical interpretation, but it tends not to work as well as the αA and m methods. But overall, the data in this study provide the most support for hypotheses 2, the αA method, and 4, the m method. According to Hypothesis 2, the actual Bowen ratio under wet-surface conditions with minimal advection is a constant fraction of the Bowen ratio under equilibrium conditions (Eq. 4). According to Hypothesis 4, α for wet surfaces remains at a constant fraction (m) of the distance between the minimum of value of 1 and the maximum value of . Since is a function of the wet-surface temperature T0, so is αm.
Without a need for any additional data, the temperature dependence of α can be included in evaporation equations. It seems appropriate to include this dependence in applications of the Priestley–Taylor (1972) equation and in particular in the use of CR models to estimate actual evaporation over drying surfaces. It is striking that four distinct hypotheses for how to understand the physical meaning of α can be stated clearly and that they have real implications regarding the nature and the numerical value of α.
All data were downloaded from the FLUXNET website (https://fluxnet.org/login/?redirect_to=/data/download-data/ (last access: 1 September 2023). The data and code used to process them are freely available (Crago et al., 2023a) through Zenodo at https://doi.org/10.5281/zenodo.8172604. They are also available (Crago et al., 2023b) from GitHub at https://github.com/r-crago/FLUXNET_Wet_sfc_evap (last access: 1 September 2023).
The supplement related to this article is available online at: https://doi.org/10.5194/hess-27-3205-2023-supplement.
RDC: conception of study, acquisition and analysis of data, and rough draft and edits. RJQ: further development of concepts and edits. JS: further development of concepts and edits.
The contact author has declared that none of the authors has any competing interests.
Publisher’s note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
The authors appreciate the helpful suggestions from the external reviewers and the editors.
Support was provided by the Ministry of Innovation and Technology of Hungary from the National Research, Development and Innovation Fund, financed under the TKP2021 funding scheme (project no. BME-NVA-02).
This paper was edited by Stan Schymanski and reviewed by two anonymous referees.
Allen, R. G., Pereira, L. S., Raes, D., and Smith, M.: Crop Evapotranspiration, FAO Irrigation and Drainage Paper No. 56, Food and Agriculture Organization of the United Nations, Rome, ISBN 92-5-104219-5, 1998.
Andreas, E. L., Jordan, R. E., Mahrt, L., and Vickers, D.: Estimating the Bowen ratio over the open and ice-covered ocean, J. Geophys. Res.-Oceans, 118, 4334–4345, 2013.
Bouchet, R. J.: Evapotranspiration reelle, evapotranspiration potentielle, et production agricole, Annal. Agronom., 14, 743–824, 1963.
Brutsaert, W.: Hydrology: An Introduction, Cambridge University Press, Cambridge, ISBN 13 978-0-521-82479-8, 2005.
Brutsaert, W.: A generalized complementary principle with physical constraints for land-surface evaporation, Water Resour. Res., 51, 8087–8093, 2015.
Brutsaert, W.: Hydrology: An Introduction, 2nd edn., Cambridge University Press, Cambridge, ISBN 978-1-107-13527-7, 2023.
Brutsaert, W. and Stricker, H.: An advection-aridity approach to estimate actual regional evapotranspiration, Water Resour. Res., 15, 443–450, 1979.
Chow, V. T., Maidment, D. R., and Mays, L. W.: Applied Hydrology, McGraw_hill, New York, NY, ISBN 0-07-010810-2, 1988.
Crago, R., Szilagyi, J., and Qualls, R. J.: Collection of FLUXNET Data, Zenodo [data set and code], https://doi.org/10.5281/zenodo.8172604, 2023a.
Crago, R., Szilagyi, J., and Qualls, R. J.: Collection of FLUXNET Data, GitHub [data set and code], https://github.com/r-crago/FLUXNET_Wet_sfc_evap (last access: 1 September 2023), 2023b.
Crago, R., Szilagyi, J., Qualls, R., and Huntington, J.: Rescaling the complementary relationship for land surface evaporation, Water Resour. Res., 52, 8461–8471, 2016.
Crago, R., Qualls, R. J., and Szilagyi, J.: Complementary relationship for evaporation performance at different spatial and temporal scales, J. Hydrol., 608, 2022WR127575, https://doi.org/10.1016/j.jhydrol.2022.127575, 2022.
Crago, R. D. and Qualls, R. J.: The value of intuitive concepts in evaporation research, Water Resour. Res., 49, 6100–6104, 2013.
Crago, R. D. and Qualls, R. J.: A graphical interpretation of the rescaled complementary relationship for evapotranspiration, Water Resour. Res., 57, 2021WR026766, https://doi.org/10.1029/2020WR028299, 2021.
deBruin, H. A. R.: A model for the Priestley-Taylor parameter α, J. Clim. Appl. Meteorol., 22, 572–578, 1982.
Eichinger, W. E., Parlange, M. B., and Stricker, H.: On the concept of equilibrium evaporation and the value of the Priestley-Taylor coefficient, Water Resour. Res., 32, 161–164, 1996.
Han, S. and Tian, F.: Derivation of a sigmoid generalized complementary function for evaporation with physical constraints, Water Resour. Res., 54, 5050–5068, 2018.
Han, S. and Tian, F.: A review of the complementary principle of evaporation: from the original linear relationship to generalized nonlinear functions, Hydrol. Earth Syst. Sci., 24, 2269–2285, https://doi.org/10.5194/hess-24-2269-2020, 2020.
Han, S., Tian, F., Wang. W., and Wang, L.: Sigmoid generalized complementary equation for evaporation over wet surfaces: A nonlinear modification of the Priestley-Taylor equation, Water Resour. Res., 57, 2020WR028737, https://doi.org/10.1029/2020WR028737, 2020.
Hersbach, H.: The ERA5 global reanalysis, Q. J. Roy. Meteor. Soc., 146, 1999–2049, https://doi.org/10.1002/qj.3803, 2020.
Lhomme, J.-P.: An examination of the Priestley-Taylor equation using a convective boundary layer model, Water Resour. Res., 33, 2571–2578, 1997a.
Lhomme, J.-P.: A theoretical basis for the Priestley-Taylor coefficient, Bound.-Lay. Meteorol., 82.2, 179–191, 1997b.
Loveland, T. R., Zhu, Z., Ohlen, D. O., Brown, J. F., Reed, B. C., and Yang, L.: An analysis of the IGBP global land-cover characterization process, Photogramm. Eng. Rem. S., 65, 1021–1032, 1999.
Martens, B., Miralles, D. G., Lievens, H., van der Schalie, R., de Jeu, R. A. M., Fernández-Prieto, D., Beck, H. E., Dorigo, W. A., and Verhoest, N. E. C.: GLEAM v3: satellite-based land evaporation and root-zone soil moisture, Geosci. Model Dev., 10, 1903–1925, https://doi.org/10.5194/gmd-10-1903-2017, 2017.
Mauder, M. Foken, T., and Cuxart, J.: Surface-energy budget closure over land: A review, Bound.-Lay. Meteorol., 177, 395–426, 2020.
McMahon, T. A., Peel, M. C., Lowe, L., Srikanthan, R., and McVicar, T. R.: Estimating actual, potential, reference crop and pan evaporation using standard meteorological data: a pragmatic synthesis, Hydrol. Earth Syst. Sci., 17, 1331–1363, https://doi.org/10.5194/hess-17-1331-2013, 2013.
McNaughton, K. G.: Evaporation and advection, Q. J. Roy. Meteor. Soc., 102, 181–191, 1976.
McNaughton, K. G. and Spriggs, T. W.: An evaluation of the Priestley and Taylor equation and the complementary relationship using results from a mixed-layer model of the convective boundary layer, Estimation of Areal Evapotranspiration, IAHS P., 177, 89–104, 1989.
Miralles, D. G., Holmes, T. R. H., De Jeu, R. A. M., Gash, J. H., Meesters, A. G. C. A., and Dolman, A. J.: Global land-surface evaporation estimated from satellite-based observations, Hydrol. Earth Syst. Sci., 15, 453–469, https://doi.org/10.5194/hess-15-453-2011, 2011.
Monteith, J. L.: Evaporation and environment, Sym. Soc. Exp. Biol., 19, 205–234, 1965.
Nash, J. E. and Sutcliffe, J. V.: River flow forecasting through conceptual models part I – A discussion of principles, J. Hydrol., 10, 282–290, 1970.
Pastorello, G., Trotta, C., Canfora, E., Chu, H., Christianson, D., Cheah, Y.W., Poindexter, C., Chen, J., Elbashandy, A., Humphrey, M., and Isaac, P.: The FLUXNET2015 dataset and the ONEFlux processing pipeline for eddy covariance data, Sci. Data, 7, 1–27, 2020.
Penman, H. L.: Natural evaporation from open water, bare soil, and grass, P. R. Soc. London, 193, 120–145, 1948.
Philip, J. R.: A physical bound on the Bowen ratio, J. Clim. Appl. Meteorol., 26, 1043–1045, 1987.
Priestley, C. H. B. and Taylor, R. J.: On the assessment of surface heat flux and evaporation, Mon. Weather Rev., 106, 81–92, 1972.
Qualls, R. J. and Crago, R. D.: Graphical interpretation of wet surface evaporation equations, Water Resour. Res., 56, e2019WR026766, https://doi.org/10.1029/2019WR026766, 2020.
Raupach, M. R.: Equilibrium evaporation and the convective boundary layer, Bound.-Lay. Meteorol., 96.1, 107–142, 2000.
Reichstein, M., Falge, E., Baldocchi, D., Papale, D., Aubinet, M., Berbigier, P., Bernhofer, C., Buchmann, N., Gilmanov, T., Granier, A.. and Grünwald, T.: On the separation of net ecosystem exchange into assimilation and ecosystem respiration: review and improved algorithm. Glob. Change Biol., 11, 1424–1439, 2005.
Schuepp, P. H., Leclerc, M. Y., MacPherson, J. I., and Desjardins, R. L.: Footprint prediction of scalar fluxes from analytical solutions of the diffusion equation, Bound.-Lay. Meteorol., 50, 355–373, 1990.
Schymanski, S. J. and Or, D.: Wind effects on leaf transpiration challenge the concept of “potential evaporation”, Proc. IAHS, 371, 99–107, https://doi.org/10.5194/piahs-371-99-2015, 2015.
Slatyer, R. O. and McIlroy, I. C.: Practical Microclimatology, CSIRO, Melbourne, Australia, Q. J. Roy. Meteor Soc., 80, 559–560, https://doi.org/10.1002/qj.49708837822, 1961.
Stull, R. B.: An Introduction to Boundary Layer Meteorology, Kluwer Academic Publishers, Boston, ISBN 90-277-2768-6, 1988.
Szilagyi, J. and Jozsa, J.: New findings about the complementary relationship-based estimation methods, J. Hydrol., 354, 171–186, 2008.
Szilagyi, J., and Schepers, A.: Coupled heat and vapor transport: The thermostat effect of a freely evaporating land surface, Geophys. Res. Lett., 41, 435–41, 2014.
Szilagyi J., Parlange, M. B., and Katul, G. G.: Assessment of the Priestley-Taylor parameter value from ERA-Interim global reanalysis data, J. Hydro-Environ. Res., 2, 1–7, 2014.
Szilagyi, J., Crago, R., and Qualls, R.: A calibration-free formulation of the complementary relationship of evaporation for continental-scale hydrology, J. Geophys. Res.-Atmos., 122, 264–278, 2017.
Thornthwaite, C. W.: An approach toward a rational classification of climate, Geogr. Rev., 38, 55–94, https://doi.org/10.2307/210739, 1948.
Tu, Z., Yang, Y., and Roderick, M. L.: Testing a maximum evaporation theory over saturated land: implications for potential evaporation estimation, Hydrol. Earth Syst. Sci., 26, 1745–1754, https://doi.org/10.5194/hess-26-1745-2022, 2022.
Tu, Z., Yang, Y., Roderick, M. L., and McVicar, T. R.: Potential evaporation and the complementary relationship, Water Resour. Res., 59, 2022WR033763, https://doi.org/10.1029/2022WR033763, 2023.
Wang, L., Tian, F. Han, S., and Wei, Z.: Determinants of the asymmetric parameter in the generalized complementary principle of evaporation, Water Resour. Res., 56, 2020WR026570, https://doi.org/10.1029/2019WR026570, 2020.
Yang, Y. and Roderick, M. L.: Radiation, surface temperature and evaporation over wet surfaces, Q. J. Roy. Meteor. Soc., 145, 1118–1129, https://doi.org/10.1002/qj.3481, 2019.