HTML
--> --> -->As an alternative, the complementary relationship (CR) (Bouchet, 1963) of evaporation inherently accounts for the dynamic feedback mechanisms found in the soil?vegetation?atmosphere interface, and thus provides a viable, robust alternative for land ET estimation relying solely on standard atmospheric forcing without the need for any soil or vegetation data. The description in the next two paragraphs of the applied CR method parallels that of Ma and Szilagyi (2019).
The generalized nonlinear version of the complementary relationship (GCR) by Brutsaert (2015) relates two scaled variables, x = Ew Ep?1 and y = E Ep?1 as
Here, E (mm d?1) is the actual ET rate, while Ep (mm d?1) is the potential ET rate, i.e., the ET rate of a plot-sized wet area in a drying (i.e., not completely wet) environment, typically specified by the Penman (1948) equation as
where

where u2 (m s?1) is the 2-m horizontal wind speed.
The so-called wet-environment ET rate, Ew (mm d?1), of a well-watered land surface having a regionally significant areal extent, is specified by the Priestley and Taylor (1972) equation:
The dimensionless Priestley?Taylor (PT) coefficient, α, in Eq. (4), normally attains values in the range of [1.1?1.32] (Morton, 1983). For large-scale model applications of gridded data, Szilagyi et al. (2017) proposed a method of finding a value for α, thus avoiding the need for any calibration.
Very soon after the publication of the GCR, Crago et al. (2016) and Szilagyi et al. (2017) introduced a necessary scaling into Eq. (1) by means of a time-varying wetness index, w = (Ep,max ? Ep)(Ep,max ? Ew)?1, to define the dimensionless variable, X, as X = wx, by which Eq. (1) keeps its original nonlinear form, i.e.,
Note that Eq. (4) in Priestley and Taylor (1972) was designed with measurements under wet environmental conditions; therefore, Δ should be evaluated at the wet-environment air temperature, Tw (°C), instead of the typical drying-environment air temperature, T (Szilagyi and Jozsa, 2008; Szilagyi, 2014). By making use of a mild vertical air temperature gradient (de Vries, 1959; Szilagyi and Jozsa, 2009; Szilagyi, 2014) observable in wet environments (as Rn is consumed predominantly by the latent heat flux at the expense of the sensible one, and water representing an unusually high latent heat of the vaporization value found in nature), Tw can be approximated by the wet surface temperature, Tws (°C). Note that Tws may still be larger than the drying-environment air temperature, T, when the air is close to saturation, but the same is not true for Tw, due to the cooling effect of evaporation. In such cases, the estimated value of Tw should be capped by the actual air temperature, T (Szilagyi, 2014; Szilagyi and Jozsa, 2018). Szilagyi and Schepers (2014) demonstrated that Tws is independent of the size of the wet area. Thus, Tws can be obtained through iterations from the Bowen ratio (β) of the sensible and latent heat fluxes (Bowen, 1926) when applied over a small, plot-sized, wet patch (by the necessary assumption that the available energy for the wet patch is close to that for the drying surface) the Penman equation is valid for, i.e.,
Ep,max in the definition of X within Eq. (5) is the maximum value that Ep can achieve (under unchanging available energy for the surface) during a complete dry-out (i.e., when e becomes close to zero) of the environment, i.e.,
in which Tdry (°C) is the so-achieved dry-environment air temperature. The latter can be estimated from the (isoenthalp) adiabat of an air layer in contact with the drying surface (Szilagyi, 2018a), i.e.,
where Twb (°C) is the wet-bulb temperature. Twb can be obtained with the help of another iteration of writing out the Bowen ratio for adiabatic changes (e.g., Szilagyi, 2014), such as
For a graphical illustration of the saturation vapor pressure curve, the different temperatures and the related ET rates defined, please refer to Ma and Szilagyi (2019). The same source also includes a brief description of how the CR evolved into Eq. (5) over the past 40 years. Additionally, it plots selected historical CR functions over sample data, and explains how assigning a value of α is performed without resorting to any calibration. A sensitivity analysis of the ET rates in Eq. (5) to their atmospheric forcing is found in Ma et al. (2019).
While Brutsaert et al. (2020) realized the necessity of scaling x with the help of a static aridity index, Crago et al. (2016), Szilagyi et al. (2017), Szilagyi (2018a, b), Szilagyi and Jozsa (2018), Ma and Szilagyi (2019), and Ma et al. (2019) performed one (and the same one) via a dynamic wetness index. Whereas the wetness index assigns increasing values to wetter environmental conditions, the aridity index does the same to drier ones. Brutsaert et al. (2020) did not include this dynamic wetness index method in their study, and therefore the present work was initiated to fill this gap.
An alternative, static scaling of x by Brutsaert et al. (2020) takes place via an adjustable parameter, αc, that acts as the PT-α value for wet environments. Since Eq. (4) can also be written as Ew = αEe, where Ee is the equilibrium evaporation rate of Slatyer and Mcllroy (1961), thus the scaled variable, X, becomes X = αc Ee Ep?1 = αcxα?1. The spatially variable (but constant through time at a given location) value of αc was then related to a long-term aridity index by Brutsaert et al. (2020), with the latter defined as the ratio of the mean annual Ep and rain depth, and globally calibrated with the help of additional water-balance data, requiring altogether seven parameters in highly nonlinear equations.
Note that the X = wx scaling by Crago et al. (2016) and Szilagyi et al. (2017) requires only the forcing variables (Rn, G, T, u2, and VPD), without the need for external precipitation/rain data, which is significant as precipitation is possibly the most uncertain meteorological variable to predict in climate models. It is important to mention that w changes with each value of x, unlike αc. As Szilagyi et al. (2017) demonstrated, a (temporally and spatially) constant value of the PT α, necessary for x, can be set by the sole use of the forcing variables, without resorting to additional water-balance data of precipitation and stream discharge, thus making the approach calibration-free on a large scale (Szilagyi, 2018b; Ma et al., 2019; Ma and Szilagyi, 2019) where wet-environmental conditions, necessary for setting the value of α, can likely be found. Note that setting a constant value of α is also necessary for Brutsaert et al. (2020) in order to force their spatially variable but temporally constant αc values to reach a predetermined value of about 1.3 under wet conditions. Despite almost half a century of research following publication of the Priestley and Taylor (1972) equation, there is still no consensus about what environmental variables (atmospheric, radiative, and/or surficial properties), and exactly how their spatial and temporal averaging, influence the value of the PT α. Until compelling information is available on these variables, a spatially and temporally constant α value may suffice for modeling purposes.
As was found by Szilagyi (2018b), the value of the PT α depends slightly on the temporal averaging of the forcing data, i.e., whether or not the monthly values are long-term averages [yielding α = 1.13 (Szilagyi et al., 2017) and 1.15 (Szilagyi, 2018b), respectively]. Therefore, here, it is tested if such is the case for the globally calibrated model of Brutsaert et al. (2020). Namely, if its performance is affected by similar changes (from long-term mean monthly values to monthly values) in the input/forcing variables, then some caution must be exercised during its routine future application, and recalibration of its seven parameters may be necessary. Note that besides the different scaling of x, everything is the same (including data requirements) in the two GCR model versions applied here, except that Δ in Ee is evaluated at the measured air temperature in Brutsaert et al. (2020) while the same in Ew (= αEe) is evaluated at an estimated wet-environment air temperature (Szilagyi et al., 2017) explained above.
Both model versions (denoted for brevity by GCR-stat and GCR-dyn, respectively) were tested over the conterminous United States, first with long-term averages (1981?2010) of monthly 32-km resolution North American Regional Reanalysis (NARR) (Mesinger et al., 2006) radiation and 10-m wind (u10) data [reduced to 2-m values via u2 = u10 (2/10)1/7 (Brutsaert, 1982)], as well as with 4.2-km PRISM air, and dewpoint temperature values (Daly et al., 1994) followed by a continuous 37-year simulation of monthly values over the 1979?2015 period. The NARR data were resampled to the PRISM grid by the nearest neighbor method. Monthly soil heat fluxes were considered negligible. Evaluation of the model estimates were performed by water-balance estimates of basin-representative evaporation rates (Ewb) with the help of United States Geological Survey two- and six-digit Hydrologic Unit Code (HUC2 and HUC6) basin (Fig. 1) discharge data (Q) together with basin-averaged PRISM precipitation (P) values as Ewb = P ? Q, either on an annual (for trend analysis) or long-term mean annual basis. The application of a simplified water balance is justifiable as soil-moisture and groundwater-storage changes are typically negligible over an annual (or longer) scale (Senay et al., 2011) for catchments with no significant trend in the groundwater-table elevation values.


However, the picture changes when switching from long-term mean monthly forcing values to monthly values in a continuous simulation (Fig. 3). GCR-dyn, with a slightly changed PT-α value [from 1.13 to 1.15, using the procedure of Szilagyi et al. (2017)] continues to produce unbiased estimates of basin-averaged mean annual evaporation values. However, the globally calibrated GCR-stat model underestimates the water-balance-derived values by about 10% [i.e., relative bias (RB) of ?0.09 for both basin scales] and produces reduced interannual variability (see the horizontally elongated “crosses” for the HUC2 basins in Fig. 3) in comparison with GCR-dyn. Reduced model performance of GCR-stat is also apparent in the long-term linear tendencies (obtained as least-squares fitted linear trends) of the basin-averaged annual evaporation values (Fig. 4) by being less effective than GCR-dyn in reproducing the observed linear trends in the water-balance data.


As on a mean-annual basis GCR-stat performs better than (with mean monthly values) or about equal to (in a continuous simulation) GCR-dyn by exploiting precipitation data (which on the watershed scale naturally serves as an upper bound for land ET), its weakened performance in trends can only be explained by the same reliance on the long-term means of the precipitation (and Ep) rates in the (therefore) static αc values that will hinder its response to slow (decadal) changes in aridity. The same problem cannot occur in GCR-dyn, since its wetness index (w) is updated in each step of calculations.
The current GCR-dyn model has already been shown to (a) yield correlation coefficient values in excess of 0.9 with local measurements of latent heat fluxes across diverse climates in China (Ma et al., 2019) in spite of large differences in spatial representativeness (i.e., grid resolution vs flux measurement footprint), and (b) outperform several popular large-scale ET products over the conterminous United States (Ma and Szilagyi, 2019). These products include three LSMs—namely, Noah (Chen and Dudhia, 2001), VIC (Liang et al., 1994), and Mosaic (Koster and Suarez, 1996); two reanalysis products—namely, NCEP-II (Kanamitsu et al., 2002), and ERA-Interim (Dee et al., 2011); another two remote-sensing based approaches—namely, GLEAM (Miralles et al., 2011; Martens et al., 2017) and PML (Zhang et al., 2017; Leuning et al., 2008); and a spatially upscaled eddy-covariance measurement product, FLUXNET-MTE (Jung et al., 2011). In a comparison with water-balance data, GCR-dyn turns out to produce even better statistics on the HUC2 scale than the spatially interpolated eddy-covariance measured ET fluxes (Fig. 5), which is remarkable from a calibration-free approach. GCR-dyn especially excels in the long-term linear tendency estimates of the HUC2 ET rates, demonstrated by Figs. 6 and 7. As FLUXNET-MTE employs several temporally static variables for its spatial interpolation method, its ability to capture long-term trends is somewhat limited (Jung et al., 2011). On the contrary, the dynamic scaling inherent in GCR-dyn automatically adapts to such trends and identifies them rather accurately.



Among the different popular large-scale ET models, GCR-dyn produces multi-year mean annual ET rates closest in its spatial distribution to those of FLUXNET-MTE (Fig. 8), with a spatially averaged ET value almost identical (both in its spatial average and standard deviation) to that of GLEAM (Fig. 8), which is a remote-sensing based approach. Note that all models of the comparison (except GCR-dyn) rely on precipitation data as input, which greatly aids ET estimations, as on a regional scale and long-term basis precipitation forms an upper bound for terrestrial ET rates; plus, it may drive any required soil-moisture water-balance calculations.

The spatial distribution of the modeled multi-year linear ET trends is displayed in Fig. 9. Again, GCR-dyn produces results closest in spatial distribution to FLUXNET-MTE in terms of the statistically significant trends and to GLEAM in general. For a more detailed discussion of model comparisons (including additional model descriptions), please refer to Ma and Szilagyi (2019).

In conclusion, it can be stated that the GCR of evaporation (Brutsaert, 2015) is a very effective tool in land ET modeling, as it requires only a few, largely meteorological forcing variables, and avoids the need for detailed soil-moisture and land-cover information. Although attractive, as its (altogether seven) parameters have already been globally precalibrated, the GCR model version (GCR-stat) of Brutsaert et al. (2020) may, however, not perform optimally in estimating long-term tendencies in basin-wide ET rates. This is particularly the case in comparison to an earlier, calibration-free version (GCR-dyn), having no precalibrated parameter values but requiring that its sole, temporally?and spatially?constant parameter (i.e., the PT coefficient) be set using the actual forcing data through a largely automated method, described in Szilagyi et al. (2017). Since in a continuous monthly simulation both models performed about the same, while the GCR-dyn produced better long-term tendencies, a dynamic scaling of Ew Ep?1 is recommended over a static one in future applications of the GCR of evaporation.
As has been recommended before (Szilagyi, 2018b; Szilagyi and Jozsa, 2018; Ma and Szilagyi, 2019), it is reiterated here that GCR-dyn, due to its minimal data requirement, calibration-free nature and dynamic scaling, may continue to serve as a diagnostic and benchmarking tool for more complex and data-intensive models of terrestrial ET rates, including calibration/verification (for past values) and reality checking (for future scenario values) of the LSM-predicted ET rates of existing regional and general circulation models.
Acknowledgements. All data used in this study can be accessed from the following websites. NARR data:
Open Access This article is distributed under the terms of the creative commons attribution 4.0 international license (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the creative commons license, and indicate if changes were made.