Published online 24 January 2007
Published in Vadose Zone J 6:77-92 (2007)
DOI: 10.2136/vzj2006.0001
© 2007 Soil Science Society of America
677 S. Segoe Rd., Madison, WI 53711 USA
ORIGINAL RESEARCH
Climatic Forecasting of Net Infiltration at Yucca Mountain Using Analogue Meteorological Data
Boris Faybishenko*
Earth Sciences Division, Lawrence Berkeley National Lab., Berkeley, CA 94720
* Corresponding author (bfayb{at}lbl.gov)
Received 6 January 2006.
 |
ABSTRACT
|
|---|
At Yucca Mountain, NV, future changes in climatic conditions will probably alter net infiltration, drainage below the bottom of the evapotranspiration zone within the soil profile, or flow across the interface between soil and the densely welded part of the Tiva Canyon Tuff. The objectives of this study were to: (i) develop a semiempirical model and forecast average net infiltration rates, using the limited meteorological data from analog meteorological stations, for interglacial (present day), and future monsoon, glacial transition, and glacial climates over the Yucca Mountain region; and (ii) corroborate the computed net infiltration rates by comparing them with the empirically and numerically determined groundwater recharge and percolation rates through the unsaturated zone from published data. This study approached calculations of net infiltration, aridity, and precipitation-effectiveness indices using a modified Budyko's water-balance model, with reference-surface potential evapotranspiration determined from the radiation-based Penman formula. Results of calculations show that net infiltration rates are expected to generally increase from the present-day climate to monsoon climate, to glacial transition climate, and then to the glacial climate, following a power law relationship between net infiltration and precipitation. The forecasting results indicate the overlap between the ranges of net infiltration for different climates. Forecasting of net infiltration for different climate states is subject to numerous uncertainties associated with selecting climate analog sites, using relatively short analog meteorological records, neglecting the effects of vegetation and surface runoff and run-on on a local scale, as well as possible anthropogenically induced climate changes.
Abbreviations: DRI, Desert Research Institute G, glacial IG, interglacial IM, intermediate or glacial transition M, monsoon
 |
INTRODUCTION
|
|---|
PRESENT-DAY and potential future net infiltration is a hydrologic parameter that controls the rate of deep percolation, groundwater recharge, radionuclide transport, and seepage into tunnelswhich are all, in turn, parameters for the total system performance assessment of the nuclear waste repository at Yucca Mountain, Nevada. Net infiltration is defined as water drainage below the bottom of the evapotranspiration zone within the soil profile or flow across the interface between soil and the densely welded part of the Tiva Canyon tuff at Yucca Mountain. Because net infiltration is largely dependent on climatic conditions, future changes in climatic conditions will potentially alter net infiltration into the deep unsaturated zone (at Yucca Mountain, the depth of the unsaturated zone is on the order of 600 m [Bodvarsson et al., 2003a]).
Although a variety of sophisticated numerical models are being used for predictions of soil infiltration, a key point in selecting an adequate prediction model is to start with the simplest linear or nonlinear functions to describe the structure in the data. Then, if required, more complex models could be used, but they should not be used unnecessarily to preclude increasing the uncertainty of predictions performed to answer engineering or scientific questions. The reasonable accuracy of estimates using simple functions is demonstrated here by corroboration of predicted net infiltration rates with the results of other field and modeling studies as obtained from published sources.
Because of the limited amount of meteorological information (such as precipitation, temperature, dew point, and wind velocity records) from meteorologically analogous sites, it is reasonable to apply a relatively simple soil-water-budget approach, which has been broadly used for watershed- and regional-scale hydrological and climatological predictions (e.g., Thornthwaite, 1948; Thornthwaite and Mather, 1955; Budyko, 1948, 1951, 1974; Rasmusson, 1971; Mather, 1978; Alley, 1984; Willmott et al., 1985; Mintz and Walker, 1993; Mintz and Serafini, 1992; Milly and Dunne, 2002). Such an approach has been used successfully for annual (Mather, 1978) and long-term predictions (Brutsaert, 1982).
Conventional models for forecasting changes in the waterenergy balance usually require using such meteorological parameters as precipitation, solar radiation flux, diurnal and seasonal temperature cycles, evapotranspiration, and relative humidity. Because these parameters are not known for future climates, changes in future climatic conditions at Yucca Mountain could be forecast using meteorological records from analog meteorological stations (Sharpe, 2003; Bechtel SAIC Company, 2004a, 2004b). In particular, precipitation and temperature can generally be considered as proxy parameters affecting other processes involved in water and energy transfer in an atmosphericshallow subsurface system.
The objectives of this study were to: (i) develop a semiempirical model and forecast average net infiltration rates, using the limited meteorological records from analog meteorological stations, for interglacial (present-day), monsoon, intermediate (glacial transition), and glacial climates across the Yucca Mountain region expected for the next 500 000 yr; and (ii) corroborate the forecast net infiltration rates by comparing them with empirically and numerically determined groundwater recharge and percolation fluxes through the unsaturated zone at different field sites as gathered from published data.
First, the data characterizing present-day and future climates are described, reconciling the Desert Research Institute (DRI) (Sharpe, 2003) and USGS (Thompson et al., 1999; USGS, 2001; Bechtel SAIC Company, 2004a) reports and records from analog meteorological stations. The conceptual model and main assumptions of the semiempirical approach used for net-infiltration forecasting for Yucca Mountain's analog meteorological stations are discussed. The results of calculations of net infiltration and the aridity and precipitation-effectiveness indices for these meteorological stations are presented. The types of uncertainties involved in climatic forecasting of net infiltration are summarized and the results of corroboration studies are presented in comparison with published data.
Forecasting of net infiltration for different climate states is subject to numerous uncertainties: selection of climate analog sites, the use of relatively short meteorological records from the analog meteorological stations, neglecting the effects of vegetation and surface runoff and run-on on a local scale, as well as possible anthropogenic climate changes. A detailed analysis of how these factors would affect net infiltration is beyond the scope of this study, however.
 |
CHARACTERIZATION OF PRESENT-DAY AND FORECASTING FUTURE CLIMATES
|
|---|
Types of Climatic Data and Climate Timing
Characterization of climatic conditions at Yucca Mountain is based mainly on the results of the USGS (USGS, 2001; Thompson et al., 1999) and DRI (Sharpe, 2003) paleogeographic and paleoclimatic investigations of the fossil records, specifically the ostracode and diatom assemblages recovered from Owens Lake, California (Sharpe, 2003), and Devils Hole, Nevada (Winograd et al., 1992), as well as Vostok Station, Antarctica (Petit et al., 1999) and orbital cycle periods (Milankovitch theory). Sharpe (2003) identified the sequence and duration of past climate states during a period of 500 000 yr, including: (i) interglacial climate (IG, present-day); (ii) monsoon (M); (iii) intermediate (IM, glacial transition); (iv) glacial 4/2 (G 4/2, which corresponds to two equivalent oxygen isotope stages [OIS] 4 and 2), (v) glacial 10/8 (G 10/8, which corresponds to two equivalent OIS 10 and 8), and (vi) glacial 16/6 (G 16/6, which corresponds to two equivalent OIS 16 and 6).
Table 1 presents the duration of past climate states, indicating that the total duration of glacial climate states was 18.3%, with the longest total duration (63.6%) for the IM climate. The common approach to forecasting future climate states is based on the assumption that the sequence and duration of past climate states will recur in the future (Knox, 1991). For each climate, Sharpe (2003) identified two types of climatic conditions: the lower-bound climate, causing lower net infiltration; and the upper-bound climate, causing higher net infiltration.
Present-Day Climate
Both USGS (2001) and DRI reports (Sharpe, 2003) indicate the existence of a long-term, present-day interglacial climate state for at least the last 9000 yr before the present. The present-day climate is estimated to last
600 more years. The present-day meteorological conditions of the Yucca Mountain region feature a mean annual precipitation of 125 mm and a mean annual temperature of 13.4°C (Thompson et al., 1999, Table 4, Fig. 16 and 17). The special distribution of meteorological parameters across the Yucca Mountain region has been characterized using the data collected from a network of nine automated weather stations (Bechtel SAIC Company, 2004b); however, the meteorological conditions are changing with elevation and time. For example, evidence has recently accumulated that one of the most important features of the present-day climate is that the world climate has begun to warm since the early 1900s. Temperature increased nearly 1°C during the 20th century. Although the causes of this warming are not fully understood, one of the possible reasons for warming is the release of CO2 and other greenhouse gases into the atmosphere (Muller and MacDonald, 2000). The pattern of increasing temperature and precipitation during the past century indicate that the mean temperature and precipitation calculated from the last 30 to 60 yr of observations at analog meteorological stations may not be statistically representative of the future interglacial climate, if temperature and precipitation continue to increase with time.
Future Climates
The future interglacial climate states are assumed, in general, to be comparable to the relatively warm present-day climate state. (For the last 420 000 yr, brief periods of interglacial peaks lasted typically from a few thousand to perhaps 20000 yr [Muller and MacDonald, 2000].)The monsoon climate state is characterized by hot summers with increased summer rainfall relative to the present-day climate. This monsoon climate is somewhat similar to the climate in the equatorial region, because of a similar abundant precipitation (rainfall is distributed seasonally as in tropical climates) and temperature regime, even though annual excursion is higher by about 7 to 8°C. Monsoon climate conditions can presently be found in the southwestern USA (Wright et al., 2001; Cavazos et al., 2002; Douglas et al., 1993).
The glacial-transition climate state is expected to have cooler and wetter summers and winters relative to the present-day climate. The future glacial climate is expected to be wetter (pluvial) and cooler than the present-day climate. According to analog-based precipitation estimates, the mean annual precipitation for the last glacial maximum was from 266 to 321 mm/yr, which is within the range of the upper-bound present-day precipitation; and the mean annual temperature was 7.9 to 8.5°C, which is near the lower bound of the present-day temperature range for Nevada District 3 (Thompson et al., 1999).
Analog Meteorological Station Data
The locations of the analog meteorological stations (identified by Sharpe, 2003) for the Yucca Mountain future climates are shown in Fig. 1
. Individual meteorological stations provide meteorological records, which are obtained at a point scale and for a limited duration of monitoring, only in a few instances exceeding 100 yr (Table 2). The relationships between the mean annual precipitation and temperature for present-day, monsoon, intermediate, and glacial climates, using data from analog meteorological stations, are summarized in Fig. 2
. The monthly meteorological data for analog meteorological stations were taken from the database of the Water Regional Climate Center (WRCC) of the DRI, Reno, NV, at http://www.wrcc.dri.edu/ (verified 23 Oct. 2006).

View larger version (35K):
[in this window]
[in a new window]
|
Fig. 2. Changes in mean annual temperature (MAT) and mean annual precipitation (MAP), using data from analog meteorological stations (see Tables 2 and 3). Open circles are the MAT and MAP data for interglacial (present-day) climate from the Yucca Mountain Meteorological Sites 2 and 5.
|
|
View this table:
[in this window]
[in a new window]
|
Table 3. Results of calculations of reference potential evapotranspiration (Eo), net infiltration, and net infiltration, precipitation-effectiveness (P-E), and aridity indices for analogue meteorological stations.
|
|
The precipitation and temperature data from the analog meteorological stations are assumed to be constant for each climate state. In other words, these data do not take into account the dynamic pattern of temperature changes with time, as determined from the Devils Hole (Winograd et al., 1992) and Vostok ice core (Muller and McDonald, 2000) data analysis.
 |
SOIL-WATER-BALANCE MODEL FOR CALCULATIONS OF NET INFILTRATION USING CLIMATIC DATA
|
|---|
Soil-Water Balance and Main Assumptions
The general form of the water-balance equation for the evaluation of net infiltration can be given by:
 | [1] |
where In is the net infiltration, P is the total precipitation, including the snowmelt, ET is the evapotranspiration, S is the change in soil water storage, Roff is the runoff, and Ron is the run-on. Semiempirical formulae are generally used to determine different components of the water-balance equation for large-scale characterization of soil moisture balance (Rasmusson, 1971; Milly and Dunne, 2002).
Depth and time intervals of the soilrock profile, for which the components of Eq. [1] are calculated, are generally dependent on the investigation objectives. Despite large values of net radiation (largely affecting potential evapotranspiration) at Yucca Mountain, episodic infiltration (of precipitated and snowmelt water) into the subsurface may cause preferential and transient flow through the upper portion of a deep unsaturated zone (Scanlon et al., 1997). Walvoord et al. (2002b) incorporated into their vapor transport model observations of temporally invariant matric potentials at 3- to 5-m depths during
5-yr monitoring periods, and simulated the presence of net upward water movement from 3- to
10- or 20-m depths. Yet, the conventional Cl mass balance approach indicated an overall downward advective liquid flux into a deep unsaturated zone.
In general, all terms of Eq. [1] are likely to vary with time, as affected by changes in climatic conditions. The time step may vary from 1 d to tens of years or longer, and the depth may vary from the topsoil depth to the depth of seasonal fluctuations of moisture content or the depth of evapotranspiration. Using the water-balance approach, which was developed for large-scale investigations (Dooge, 1988), we assumed a steady-state (time-averaged) net-infiltration regime for each climate. The errors that could be caused by this assumption should be further evaluated, because modeling of the coupled liquidgasheat movement through a deep unsaturated zone in arid environments indicates the presence of unsteady water flow even after 10000 to 15000 yr of continuous drying (Walvoord et al., 2002a). For the first-order estimation of long-term average net infiltration for future climates, we also assumed that (i) soil water storage does not change, (ii) lateral water motion within the shallow subsurface is negligible, and (iii) the terms of the surface water runoff and run-on in a regional-scale water-balance model simply cancel each other out and need not be included in the large-scale, regional water-balance model for the net-infiltration estimation. The latter is based on the results of field monitoring within the arid and semiarid areas of the southwestern USA, indicating that stream runoff at the mountain front is generally ephemeral and almost always disappears within the mountain front zone. Consequently, downstream runoff beyond the mountain front could be considered negligible, leading to a simplification of the water-balance model (Wilson and Guan, 2004). The surface runoff and run-on are likely to affect net infiltration at the local scale, such as the crest of Yucca Mountain, and could change with changes in climatic conditions. The estimates of surface runoff and run-on under the influence of climate, however, are beyond the scope of this study. Therefore, in our study, we assumed that the surface runoff and run-on within the watershed cancel each other out, so that all surplus water presents a source of net infiltration.
Semiempirical Budyko's Hydrological Model
For long-term estimates, at least for 1 yr, assuming that the change in moisture storage in the soil and the net ground heat flux are small, and that a sensible heat flux is positive, the evapotranspiration, E, can be expressed as a function of the aridity index,
= Eo/P, where Eo is the potential evapotranspiration (Arora, 2002):
 | [2] |
Budyko (1974) used net radiation as a surrogate for potential evapotranspiration Eo, and stated that if Eo = R/L (where R is the net radiation and L is the latent heat of evaporation) then the following conditions should satisfy:
 | [19] |
These conditions would determine the form of the function f(
). Several formulae were developed to describe the empirical relationship between precipitation and the aridity index. Schreiber (1904 [in Arora, 2002]) was probably the first to propose an exponential relationship to express the relation between E, P, and the aridity index,
, given by
 | [3] |
Then Ol'dekop (1911 [in Arora, 2002]) developed a hyperbolic tangent relationship, given by
 | [4] |
Using the water-balance data from a number of catchments around the world, Budyko (1974) found that empirical data were scattered between the curves described by the exponential relationship (Eq. [3]) of Schreiber (1904 [in Arora, 2002]) and the hyperbolic tangent relationship (Eq. [4]) of Ol'dekop (1911 [in Arora, 2002]). To describe experimental data, Budyko (1974) used the geometric mean of the right-hand sides of Eq. [3] and [4] given by
 | [5] |
or, in a simpler form,
 | [6] |
Equation [6] was initially tested for 29 European river basins (Budyko, 1951) and then for 1200 regions with known precipitation and runoff data (Budyko and Zubenok, 1961). Although the original Budyko's model was developed for the determination of surface runoff, the Budyko-like approach was also used to assess an infiltrationrunoff component of the water balance and the catchment-scale soil moisture capacity (Potter et al., 2005). Several studies have been published in which the researchers described experimental data obtained on the watershed scale using various relationships analogous to that of Budyko. For example, Milly and Dunne (2002) conducted their studies for large river basins (10000 km2 and greater) and Sankarasubramanian and Vogel (2003) incorporated the soil moisture storage capacity into their Budyko-like model, based on the results of observations at 1337 watersheds throughout the USA with at least 10 yr of records.
Several other Budyko-like models have been used for hydrological calculations. For example, the generalized TurcPike equation is given by
 | [7] |
and was tested using data from 250 catchments from different climatic zones (Pike, 1964). (In the original Turc [1954 (in Brutsaert, 1982)] equation, the first coefficient is 0.9.) Zhang et al. (2001) implemented the "plant-available water coefficient" (introduced by Milly, 1994) to represent soil moisture transpiration by plants. The rational function equation developed by Zhang et al. (2001) is given by
 | [8] |
where w is the plant water-availability coefficient, which is proportional to the root-zone depth. To take into account Budyko's idea of using net radiation to represent the value of potential evaporation, Zhang et al. (2001) used the Priestley and Taylor (1972) formula for calculating Eo.
Figure 3
shows close agreement between various curves relating the evaporation ratio (E/P) and the aridity index,
= Eo/P, using the Budyko (1974), Turc (1954 [in Brutsaert, 1982]), and Zhang et al. (2001) formulae. This figure shows two curves calculated using the Zhang et al. (2001) formula, given by Eq. [8]: for w = 0.5 (for pasture) and w = 2 (for forests). The statistical analysis of curves shown in Fig. 3 indicates that the mean relative error when using Budyko's curve is only 0.7%, in comparison with the average of all other curves shown in Fig. 3. An example of the comparison of experimental data and calculated curves from the Zhang et al. (2001) study is shown in Fig. 4
.

View larger version (31K):
[in this window]
[in a new window]
|
Fig. 3. Plots of the relationship between the evaporation (E)/precipitation (P) ratio and the aridity index (Eo/P) calculated from different semiempirical formulae, illustrating that Budyko's curve (Eq. [6]) is in the middle of curves from other formulae.
|
|
Figures 3 and 4 show that the E/P vs.
curves approach unity asymptotically as the aridity index increases. The straight segments A and B reflect the physical constraints of a water-balance model: the straight line A presents an asymptote for energy-limited evapotranspiration, and the straight line B presents an asymptote for water-limited evapotranspiration. The annual and seasonal cycling of climate may cause the transition between segments A and B (Budyko and Zubenok, 1961; Milly, 1994; Milly and Dunne, 2002).
Budyko and Zubenok (1961) showed that the mean discrepancy between the evapotranspiration calculated from Eq. [6] and that derived by the water balance was about 10%. Budyko (1974) also stated that this relationship could be applied to most mountainous basins (but not for the highest mountain basins) and to watersheds with runoff that does not vary appreciably across the area. The departure from the classical Budyko curve could be caused by biases in estimations of precipitation, discharge, net radiation, or potential evaporation, and human disturbance of natural water fluxes in arid basins (Milly and Dunne, 2002).
Although Budyko (1974) hypothesized that radiative energy supply is equivalent to the upper bound of the latent heat flux, Milly and Dunne (2002) showed that actual evaporation could exceed that determined from net radiative energy supply. Milly and Shmakin (2002, p. 302) indicated that "[O]verall, no model performed substantially better than Budyko's equation, and most models performed much worse. The superior performance of Budyko's equation was found despite the fact that most or all of the models had the advantage of using information on the global distribution of surface characteristics."
Calculations of Net Infiltration and the Precipitation-Effectiveness Index
Based on the assumptions introduced above, for large spatial and long-term temporal scales, all surplus water calculated from the water-balance equation will leave the system as net infiltration, which can be determined from
 | [9] |
or
 | [10] |
where the ratio In/P can be called a net infiltration index (dimensionless value or a percentage of the total precipitation, i.e., the sum of precipitation and snowmelt). Using E/P = f(
) calculated from Eq. [6] as an example, Fig. 5
demonstrates the variations of net infiltration for different values of Eo. The approach to calculating the value of Eo for the evaluation of the aridity index is described below.

View larger version (32K):
[in this window]
[in a new window]
|
Fig. 5. Net infiltration index (net infiltration as a percentage of precipitation) for different reference potential evapotranspiration (Eo) values (given on the curves in millimeters per year), calculated from Eq. [10] using the Budyko model (Eq. [6]).
|
|
Using precipitation and temperature as proxies representing climatic processes, the moisture conditions can be characterized using the Thornthwaite precipitation-effectiveness (P-E) index (NRCS, 2002). The P-E Index is calculated using monthly precipitation and temperature values (Thornthwaite, 1931):
 | [11] |
where the monthly P-E ratio is 11.5P/(T 10)]10/9, P is the average monthly precipitation (in inches, with 0.5 being the minimum value), T is the average monthly temperature (°F, minimum temperature of 28.4°F is used in calculations), and summation is provided for 12 mo of the year. (The results of calculations of the relationship between the P-E and net-infiltration indices are given below).
Evaluation of Reference Potential Evapotranspiration
Rationale for Selecting a Method for the Evaluation of Potential Evapotranspiration
Evapotranspiration, which combines bare-soil evaporation and transpiration by plants, is a dominant water-balance component in arid and semiarid areas. The potential evapotranspiration is often determined using various experimental methods and mathematical formulae, which, however, may often produce inconsistent results (Lu et al., 2005), especially for interannual predictions (Sankarasubramanian and Vogel, 2003). The determination of evapotranspiration is particularly difficult for mountain areas with varying elevation, vegetation, and runoff (Wilson and Guan, 2004). Furthermore, significant uncertainty and ambiguity in estimating potential evapotranspiration are caused by limited meteorological data (Brutsaert, 1982).
Semiempirical methods used for the evaluation of potential evapotranspiration can be grouped into two categories: (i) reference-surface potential evapotranspiration (for example, temperature-based HargreavesSamani, Thornthwaite, Hamon, JensenHaise, and Turc models, and radiation-based PriestleyTaylor and Penman methods), and (ii) surface-dependent potential evaporation (for example, radiation-based PenmanMonteith and Shuttleworth and Wallace [1985] methods). The reference-surface potential evapotranspiration is defined as evapotranspiration that would occur from a land surface with a "reference crop," which is usually a short, uniform, green plant cover (such as alfalfa [Medicago sativa L.] or grass) under designated weather conditions and well-moist soil (Federer et al., 1996). Although empirical reference-surface Eo relationships take into account the effect of meteorological factors, they do not explicitly include the effect of vegetation. The surface-dependent Eo depends on the surface and aerodynamic resistances, which are used to account separately for transpiration and soil evaporation. Because the reference-surface Eo is a climatic parameter and is computed from meteorological data, it expresses the evaporation rate generated by the atmosphere at a specific location and time, with no effects of crop characteristics and soil factors (Allen et al., 1998).
To calculate reference-surface potential evapotranspiration to represent the effect of net radiation in the Budyko model, this study used the Penman (1948) model, which is known to produce accurate results (Thom et al., 1981; Jensen et al., 1990). Another reason for using this formula is the fact that the WRCC database contains practically all meteorological parameters from observations at analog meteorological stations, which are needed for calculations using the Penman model. The meteorological records in the WRCC database contain the following types of average-monthly data, which we used in our calculations: total precipitation (precipitation plus snow melt); minimum, maximum, and mean air temperature; dew point temperature; wind speed; solar radiation; and pan evaporation (determined using Class A evaporation pans). The types of meteorological data used in our calculations are summarized in Table 2.
Estimates of Reference-Surface Potential Evapotranspiration
Penman Model
Penman's equation (Penman, 1948) combines the two main processes affecting the evaporation rate, or evapotranspiration rate from a well-watered surface: (i) the energy input, and (ii) the aerodynamic exchange between the surface and atmosphere. Accordingly, the common two-term form of the Penman (1948) equation for the evaluation of Eo is given by
 | [12] |
where
is the slope of the saturation vapor pressuretemperature curve,
is the psychometric constant, Rn is the net radiation expressed in water-depth units (equivalents of energy), G is the soil heat flux, which can be assumed zero for annual (or longer) predictions, and Ea is the aerodynamic transport term, which is commonly given by
 | [13] |
where f(u) is the wind speed (u) function, es is the saturation vapor pressure, and ed is the saturation vapor pressure corresponding to the dew point temperature. Various forms of the wind function f(u) (depending on crop types, the height of measurements, and other factors) were described by Hatfield and Allen (1996). In this study, we used the function
 | [14] |
with coefficients a = 1 and b = 0.56, originally proposed by Penman. The Penman formula estimates reference-surface evapotranspiration from unvegetated (or sparsely vegetated) areas.
Assuming that under abundant water-supply conditions evapotranspiration would eventually attain an equilibrium rate, the actual evapotranspiration rate would be equal to the Penman potential evapotranspiration. To express the aerodynamic effect on the potential evapotranspiration, Priestley and Taylor (1972) multiplied the first term of the Penman equation by a factor
= 1.26; however, this factor could vary depending, for instance, on the surface roughness and soil moisture content, and may underestimate both peak and seasonal evapotranspiration because of neglecting the advection term in the heat balance equation. For example, Fisher et al. (2005), who compared several evapotranspiration models for a Sierra Nevada forest ecosystem, concluded that the PriestleyTaylor model performed well to estimate the actual evapotranspiration using a modified
factor as a function of the actual soil-moisture content.
In areas with no or a small water deficit, approximately 95% of the annual evaporative demand is supplied by radiation (Stagnitti et al., 1989). Shuttleworth and Calder (1979) reported that the difference in estimates of Eo produced using Penman and PriestleyTaylor equations is within
5% of each other. Although the Penman equation may produce accurate results (Jensen et al., 1990, p. 249), uncertainties of meteorological data for future climates may create commensurate uncertainty in predicting potential evaporation for future climates. (Note that Penman formula estimates of Eo closely match those from Class A evaporation pans with corrections involving the pan-adjusted coefficient for dry areassee below).
Conversion of Pan Evaporation to Reference Evapotranspiration
Direct measurements of the evaporation rate from shallow water pans at meteorological stations are commonly used for estimating potential evaporation. Evaporation-pan rates depend on the pan's geometry, latitude, elevation, solar declination, and the cloud coverage, and usually overestimate the potential evapotranspiration under arid climate conditions (Linacre, 1994; Allen et al., 1998). To obtain realistic estimates of potential evapotranspiration in arid climates, the results of pan evaporation measurements should be adjusted by taking into account the pan's geometry, environmental setting, and operation conditions (Rosenberg et al., 1983; Allen et al., 1998). Pan coefficients also depend on the size and state of the upwind buffer zone (fetch): the larger the upwind buffer zone, the more the air moving across the pan will be in equilibrium with the buffer zone. The equation for the evaporation-pan adjustment coefficient for dry fetch (which is more likely to represent the unvegetated or sparsely vegetated Yucca Mountain area) is given by (Allen et al., 1998, Chapter 4):
 | [15] |
where RHmean is the mean relative humidity, u2 is the wind speed at the 2-m elevation, and FET is the fetch distance, which varies from 50 to 2000 m. In our calculations, FET was 1000 m. The Kp values vary typically from 0.5 to 1.0. It will be illustrated below that calculations of Eo using Penman's formula for Yucca Mountain analog meteorological stations show a good agreement with the corrected values of Eo determined using Class A evaporation pans, as well as PriestleyTaylor's formula.
 |
RESULTS OF CALCULATIONS
|
|---|
Net infiltration for Analog Meteorological Stations
Table 3 presents the results from calculating the potential evapotranspiration and net infiltration for different climates. Using the calculated net infiltration rates, Fig. 6
illustrates the relationship between net infiltration, In, and the mean annual precipitation, Pm (both in mm/yr), given by
 | [16] |
with R2 = 0.93.

View larger version (22K):
[in this window]
[in a new window]
|
Fig. 6. Relationship between calculated net infiltration and precipitation, showing the names of analog meteorological stations and climates. Also shown are the forecast data (black dots) and the power-law regression line (Eq. [16]); L indicates the lower bound climate state, and U, the upper bound climate state.
|
|
Figure 7
present the plots of climatic ranking of the net-infiltration index (% of precipitation) and net infiltration rates (mm/yr). These plots demonstrate a general trend of increasing net infiltration from the present-day climate to monsoon, glacial transition, and then to glacial climate. For the glacial climate, net infiltration during the G 16/6 climate (its duration is only 2.5% of the total duration of future climates, see Table 1) ranges from 39.9 to 213 mm/yr, which is the maximum net infiltration for the glacial climate. Net infiltration for the G 4/2 climate (its duration is 7.3% of the total duration of future climates) is from 5.5 to 71.1 mm/yr, and it overlays the lower bound of the G 16/6 net infiltration. At the same time, it roughly corresponds to the net infiltration rate for the glacial transition climate. The G 10/8 (its duration is 8.5% of the future climates) net infiltration rate generally exceeds that of the G 4/2 climate; its lower bound is within that for the glacial transition climate, and its upper bound exceeds that for the glacial transition.

View larger version (46K):
[in this window]
[in a new window]
|
Fig. 7. Climatic ranking of ranges of forecast net infiltration index (upper panel), and net infiltration (lower panel) for different climates. On the lower panel, for the present-day (interglacial) climate, red dashed lines show the ranges of the percolation flux from calculations using a Cl mass balance model (solid diamonds), calcite mass model (open diamonds), temperature data (closed circles), and experts' evaluation (solid triangles).
|
|
Aridity and Precipitation-Effectiveness Indices
The aridity index can be used to classify climatic regimes (Ponce et al., 2000): arid (12 >
5), semiarid (5 >
2), subhumid (2 >
0.75), and humid (0.75 >
0.375). Figure 8a
depicts the ranking of the annually averaged aridity indices, which is generally consistent with that from the net-infiltration ranking shown in Fig. 7. Figure 8a shows that for the present-day climate, the aridity index ranges from that typical for arid (lower bound arid climate net infiltration) and semiarid climates (upper bound arid climate net infiltration); the monsoon climate is characterized by the aridity index spanning from the arid climate (lower bound monsoon infiltration) to the border between the semiarid and subhumid climates (upper bound monsoon net infiltration). For the intermediate (glacial transition) climate, the aridity index spans the range from the middle of the semiarid climate to the low aridity subhumid indices. Finally, for the glacial climate, the aridity index is mostly within the range typical for a subhumid climate, and it even decreases to that for a humid climate for the G 16/6 climate.

View larger version (27K):
[in this window]
[in a new window]
|
Fig. 8. Climatic ranking of (a) the annual average aridity index, and (b) precipitation-effectiveness (P-E) index. Relationships of (c) the aridity index vs. net infiltration index, and (d) the P-E index vs. net infiltration index.
|
|
Climatic ranking of the P-E indices, shown in Fig. 8b, has essentially the same trend as that of the net-infiltration indices, because there is virtually a linear relationship between the P-E and net-infiltration indices. Figures 8c and 8d show the fitting curves for the net-infiltration vs. the P-E and aridity indices, which can be used for forecasting net infiltration if these indices are known.
 |
CORROBORATION OF THE FORECASTING RESULTS
|
|---|
Sources of Uncertainties and Approach to Corroboration
An often-encountered difficulty in the evaluation of model-predicted components of the water balance, including evapotranspiration and net infiltration, is the lack of widespread field observations that can be used to compare model predictions at the spatial and temporal scales. It is apparent that a significant error (or uncertainty) in evaluating net infiltration from the regional water-balance model could result from net infiltration being the smallest component of the water-balance equation. In other words, net infiltration is computed as the difference between other, much greater values of the water-balance equation (e.g., precipitation, evapotranspiration, and runoff/run-on). Moreover, the difficulty in validating computed values of net infiltration for future climates at Yucca Mountain arises from there being no reliable direct (field) measurements of net infiltration representing different climatic conditions.
As part of establishing confidence in the results of this study, the approach developed here was corroborated by comparing the results of evapotranspiration and net-infiltration calculations with other independently determined estimates. Below, the estimates of Eo from the Penman formula are compared with measurements conducted using Class A evaporation pans and calculations using the PriestleyTaylor formula for different meteorological stations. Then the estimates of net infiltration are compared with local and area-averaged groundwater recharge and percolation flux through the unsaturated zone from different sites, using published data. This approach is based on the assumption of steady-state water flow through the unsaturated zone, in spite of the results of modeling that show that deep flow and transport processes are still responding slowly to large shifts in PleistoceneHolocene climatic and vegetation changes that occurred about 10000 to 15000 yr ago (Walvoord et al., 2002b).
Comparison of Computed and Experimentally Determined Evapotranspiration Rates
To establish confidence in the results of the evaluation of the reference-surface potential evapotranspiration, we compare the estimates of potential evaporation using the Penman (1948) and Priestley and Taylor (1972) formulae with field observations conducted using Class A evaporation pans at different meteorological stations. The measured Class A evaporation rates were corrected using the correction coefficient suggested in FAO56 recommendations for dry surfaces (Allen et al., 1998, Ch. 4)see above. Figure 9
illustrates a good agreement between the estimates using the semiempirical Penman and PriestleyTaylor formulae and corrected evaporation pan measurements from analog meteorological stations. Our results agree with the conclusions of comprehensive experimental and theoretical studies by Thom et al. (1981), who showed a good comparison of the results of corrected evaporation pan measurements with those computed using the Penman formula.
Comparison of Net Infiltration with Groundwater Recharge
One of the widely used methods for estimating recharge is the MaxeyEakin method (Maxey and Eakin, 1950). This method was used in several previous water-balance studies of the Death Valley region to estimate groundwater basins' recharge. According to Maxey and Eakin (1950): (i) for precipitation <203 mm/yr, no groundwater recharge occurs; (ii) for precipitation from 203 to 304 mm/yr, groundwater recharge is 3% (this estimate corresponds to the results of the water-balance calculations of discharge measurements from springs south of Yucca Mountain near the NevadaCalifornia border by Winograd and Thordarson [1975]); (iii) for precipitation from 305 to 380 mm/yr, groundwater recharge is 7%; (iv) for precipitation from 381 to 507 mm/yr, groundwater recharge is 15%; and (v) for precipitation of 508 mm/yr and greater, groundwater recharge is 25%. The MaxeyEakin recharge rates were determined from groundwater balance estimates of the recharge and discharge, depending on the depth to the water table, for 13 valleys in east-central Nevada. By comparing the MaxeyEakin estimates with 40 estimates of recharge obtained from the southern Great Basin, using a basinwide water-budget analysis, and 27 estimates of recharge obtained using geochemical and numerical modeling approaches, Avon and Durbin (1994) and Harrill and Prudic (1998) concluded that the MaxeyEakin method provides reasonable estimates of recharge for basins in Nevada. Several studies have presented modified and updated versions of the MaxeyEakin method, based on recent precipitation data, geochemical data, and basinwide water-balance data (D'Agnese et al., 1997; Donovan and Katzer, 2000).
In the MaxeyEakin method, the areas with annual precipitation of <200 mm are not considered to recharge the groundwater. At Yucca Mountain, however, recharge is known to occur within areas where annual precipitation is <200 mm. Therefore, the comparison of the calculated net infiltration with that from the MaxeyEakin coefficients for the annual precipitation of <200 mm is invalid. Moreover, estimates of net infiltration for the Yucca Mountain area may not correspond directly to recharge because of the time lag between the net infiltration and groundwater recharge in the thick unsaturated zone.
Figure 10
summarizes the results of comparing forecast net infiltration for analog meteorological stations with estimation of groundwater recharge determined using various independent field methods and modeling, including:
- The MaxeyEakin recharge rates
- Groundwater recharge estimates, using a Cl-balance method, for two small, upland watersheds in central and south-central Nevada310 mm/yr, or about 50% of the estimated average annual precipitation of 639 mm, and 33 mm/yr, or 9.8% of the average precipitation of 336 mm/yr (Lichty and McKinley, 1995, Table 15)
- Groundwater recharge rates for Fenner Basin of the eastern Mojave Desert, California (Davisson and Rose, 2000)
- Assessments of mountain front recharge for various locations, from Table 2 of Wilson and Guan (2004)
- Groundwater recharge rates for Huntington Valley in northern Nevada (Czarnecki, 1985)
- Groundwater recharge rates for northeastern Arizona determined from 14C and Cl data (Zhu, 2000)
- An empirical power-law relationship given by Wilson and Guan (2004):
 | [17] |
where Rg is the groundwater recharge, Pm is the mean annual precipitation (both in mm/yr). Figure 10 shows that this equation deviates from Maxey-Eakin estimates for Pm >600 mm/yr.
- An empirical power-law relationship for subsurface flow and surface runoff in mountain areas, which potentially become the groundwater recharge, at Carson Basin, Nevada, given by Maurer and Berger (1997):
 | [18] |

View larger version (20K):
[in this window]
[in a new window]
|
Fig. 10. Comparison of climatic forecasting of net infiltration vs. precipitation with groundwater recharge from published data.
|
|
In Eq. [18], Rg and Pm are also in millimeters per year. Figure 10 shows that calculations using this equation exceed the results of the MaxeyEaking estimates for Pm <350 mm yr1.
To provide confidence in the results of calculations of net infiltration, Fig. 7 (lower panel) also includes the estimates of percolation rates through the Yucca Mountain unsaturated zone from several independent corroborative studies: Cl mass balance, from 0.73 to 10.6 mm/yr (Liu et al., 2003); calcite data, from 2 to 6 mm/yr (Xu et al. (2003); temperature measurements in boreholes at the crest of Yucca Mountain, 5 to 10 mm/yr (Bodvarsson et al., 2003b); and the results of the experts' evaluation of net infiltration, from 3.9 to 12.7 mm/yr (Civilian Radioactive Waste Management System Management & Operating Contractor, 1997b).
Thus, Fig. 7 and 10 demonstrate that computed net infiltration rates vs. precipitation for analog meteorological stations correspond relatively well to independently determined empirical and numerical estimates of groundwater recharge and percolation rates from published data.
 |
SUMMARY AND CONCLUSIONS
|
|---|
It is essential to forecast the range of (or to bound) net infiltration across the Yucca Mountain areafor both the present-day climate state and future climatic conditions representing the monsoon, glacial transition, and glacial climatesto assess long-term repository performance. These climate conditions are represented using temporally limited meteorological records of monthly averaged total precipitation, temperature, solar radiation, dew point temperature, and evapotranspiration from analog meteorological stations at Yucca Mountain.
The developed semianalytical model is based on computing net infiltration from Budyko's empirical water-balance model, using the estimates of reference-surface potential evapotranspiration from the Penman (1948) formula (for the analog meteorological stations, the estimates of potential evapotranspiration from the Penman formula are in close agreement with PriestleyTaylor and adjusted Class A pan evaporation measurements).
The results of calculations were used for ranking net infiltration, along with aridity and precipitation-effectiveness indices, for future climatic scenarios. We determined a general power law trend of increasing net infiltration from the present-day climate to monsoon, to intermediate (glacial transition), and then to glacial climates. The ranking of the aridity and P-E indices is practically the same as that of net infiltration. The calculated net infiltration rates for the Yucca Mountain analog meteorological stations have yielded a good match with other field and modeling study results pertaining to groundwater recharge and percolation flux through the unsaturated zone. This comparison indicates the robustness of the simple water-balance approach used here.
Future research should include the evaluation of uncertainties related to selecting analog meteorological sites spanning the anticipated range of meteorological conditions within each climatic state, calculations using relatively short meteorological records (for example, only precipitation and temperature) from the analog stations, and accounting for possible anthropogenic climate changes. Future research should also include the evaluation of uncertainties and deviations from the regional-scale Budyko curve (Potter et al., 2005) as affected by the soil plant-available water-holding capacity, various seasonality parameters (Milly, 1994), vegetation and plant-available water coefficient (Zhang et al., 2001), soil-moisture storage capacity (Rasmusson, 1971; Sankarasubrumanian and Vogel, 2002), surface runoff (Rasmusson, 1971; Sharif and Miller, 2006), and anthropogenic climate effects. Since infiltration rates affect the percolation flux through the unsaturated zone and groundwater recharge, it would be desirable to perform an uncertainty analysis to address how sensitive unsaturated and saturated zone contaminant transport is to the variability of infiltration.
 |
ACKNOWLEDGMENTS
|
|---|
This work was supported by the Director, Office of Civilian Radioactive Waste Management, U.S. Department of Energy, through Memorandum Purchase Order EA9013MC5X between Bechtel SAIC, LLC, and the Ernest Orlando Lawrence Berkeley National Laboratory (Berkeley Lab). The support is provided to Berkeley Lab through the U.S. Department of Energy Contract no. DE-AC03-76SF00098. Thanks to Stefan Finsterle, H.H. Liu, Norman Miller, and Dan Hawkes of LBNL for review of this paper, and Diana Swantek for help in drawing Fig. 1.
 |
REFERENCES
|
|---|
- Allen, R.G., L.S. Pereira, D. Raes, and M. Smith. 1998. Crop evapotranspiration: Guidelines for computing crop water requirements. Irrig. Drain. Pap. 56. FAO, Rome.
- Alley, W.M. 1984. On the treatment of evapotranspiration, soil moisture accounting, and aquifer recharge in monthly water balance models. Water Resour. Res. 20:11371149.
- Arora, V.K. 2002. The use of the aridity index to assess climate change effect on annual runoff. J. Hydrol. 265:164177.[CrossRef]
- Avon, L., and T.J. Durbin. 1994. Evaluation of the MaxeyEakin method for estimating recharge to ground-water basins in Nevada. Water Resour. Bull. 30(1):99111.
- Ball, R.A., L.C. Purcell, and S.K. Carey. 2004. Evaluation of solar radiation prediction models in North America. Agron. J. 96:391397.[Abstract/Free Full Text]
- Bechtel SAIC Company. 2004a. Future climate analysis. ANL-NBS-GS-000008, Rev. 01. Bechtel SAIC Co., Las Vegas, NV.
- Bechtel SAIC Company. 2004b. Yucca Mountain site description. TDR-CRW-GS-000001, Rev. 02 ICN 01. Bechtel SAIC Co., Las Vegas, NV.
- Bodvarsson, G.S., C.K. Ho, and B.A. Robinson (ed.). 2003a. Preface. Special issue: Yucca Mountain Project. J. Contam. Hydrol. 6263:12.
- Bodvarsson, G.S., E. Kwicklis, C. Shan, and Y.S. Wu. 2003b. Estimation of percolation flux from borehole temperature data at Yucca Mountain, Nevada. J. Contam. Hydrol. 6263:322.
- Brutsaert, W. 1982. Evaporation into the atmosphere: Theory, history, and applications. D. Reidel Publ. Co., Dordrecht, the Netherlands.
- Budyko, M.I. 1948. Evaporation under natural conditions. (IPST translation, Jerusalem, 1963.) GIMIZ, Leningrad.
- Budyko, M.I. 1951. On climatic factors of runoff. Probl. Fiz. Geogr. 16. (In Russian.)
- Budyko, M.I. 1974. Climate and life. Academic Press, San Diego.
- Budyko, M.I., and L.I. Zubenok. 1961. The determination of evaporation from the land surface. (In Russian.) Izv. Akad. Nauk SSSR Ser. Geogr. 6:317.[Medline]
- Cavazos, T., A.C. Comrie, and D.M. Liverman. 2002. Intraseasonal variability associated with wet monsoons in southeast Arizona. J. Clim. 15:24772490.[CrossRef]
- Civilian Radioactive Waste Management System Management & Operating Contractor. 1997a. Engineering design climatology and regional meteorological conditions report. B0000000001717570700066 Rev. 00. TRW Environ. Safety Syst. Inc., Las Vegas, NV.
- Civilian Radioactive Waste Management System Management & Operating Contractor. 1997b. Unsaturated zone flow model expert elicitation project. TRW Environ. Safety Syst. Inc., Las Vegas, NV.
- Czarnecki, J.B. 1985. Simulated effects of increased recharge on the ground-water flow system of Yucca Mountain and vicinity, NevadaCalifornia. USGS-WRI-844344. Water-Resour. Invest. Rep. 84-4344. U.S. Geol. Surv., Denver, CO.
- D'Agnese, F.A., C.C. Faunt, A.K. Turner, and M.C. Hill. 1997. Hydrogeologic evaluation and numerical simulation of the Death Valley regional ground-water flow system, Nevada and California. Water-Resour. Invest. Rep. 96-4300. U.S. Geol. Surv., Denver, CO.
- Davisson, M.L., and T.P. Rose. 2000. A calibrated MaxeyEakin curve for the Fenner Basin of the eastern Mojave Desert, California. Rep. UCRL-ID-139030. Lawrence Livermore Natl. Lab., Livermore, CA.
- Donovan, D.J., and T. Katzer. 2000. Hydrologic implications of greater ground-water recharge to Las Vegas Valley, Nevada. J. Am. Water Resour. Assoc. 36:11331148.
- Dooge, J.C.I. 1988. Modelling the behaviour of water. Ch. 7. In E.R.C. Reynolds and F.B.Thompson (ed.) Forests, climate, and hydrology: Regional impacts. United Nations Univ. Press, Tokyo.
- Douglas, M.W., R.A. Maddox, K. Howard, and S. Reyes. 1993. The Mexican monsoon. J. Clim. 6:16651677.[CrossRef]
- Federer, C.A., C. Vörösmarty, and B. Fekete. 1996. Intercomparison of methods for calculating potential evaporation in regional and global water balance models. Water Resour. Res. 32:23152321.[CrossRef]
- Fisher, J.B., T.A. DeBiase, Y. Qi, M. Xu, and A.H. Goldstein. 2005. Evapotranspiration models compared on a Sierra Nevada forest ecosystem. Environ. Modell. Software 20:783796.[CrossRef]
- Hargreaves, G.H., and Z.A. Samani. 1982. Estimating potential evapotranspiration. J. Irrig. Drain Eng. Am. Soc. Civ. Eng. 108:223230.
- Harrill, J.R., and D.E. Prudic. 1998. Aquifer systems in the Great Basin region of Nevada, Utah, and adjacent statesSummary report. Prof. Pap. 1409-A. U.S. Geol. Surv., Denver, CO.
- Hatfield, J.L., and R.G. Allen. 1996. Evapotranspiration estimates under deficient water supplies. J. Irrig. Drain. Eng. 122:301308.[CrossRef]
- Jensen, M.E., R.D. Burman, and R.G. Allen (ed.). 1990. Evapotranspiration and irrigation water requirements. Manuals and Rep. on Eng. Pract. 70. Am. Soc. Civ. Eng., Reston, VA.
- Knox, J.B. 1991. Global climate change: Impacts on California. An introduction and overview. p. 125. In J.B. Knox and A.F. Scheuring (ed.) Global climate change and California: Potential impacts and responses. Univ. of Calif. Press, Berkeley.
- Lichty, R.W., and P.W. McKinley. 1995. Estimates of ground-water recharge rates for two small basins in central Nevada. Water-Resour. Invest. Rep. 944104. U.S. Geol. Surv., Denver, CO.
- Linacre, E.T. 1994. Estimating U.S. Class A pan evaporation from few climate data. Water Int. 19:514.
- Liu, J., E.L. Sonnenthal, and G.S. Bodvarsson. 2003. Calibration of Yucca Mountain unsaturated zone flow and transport model using porewater chloride data. J. Contam. Hydrol. 6263:213235.
- Lu, J., G. Sun, S.G. McNulty, and D.M. Amatya. 2005. A comparison of six potential evapotranspiration methods for regional use in the southeastern United States. J. Am. Water Resour. Assoc. 41:621633.
- Mather, J.R. 1978. The climatic water budget in environmental analysis. Lexington Books, Lexington, MA.
- Maurer, D.K., and D.L. Berger. 1997. Subsurface flow and water yield from watersheds tributary to Eagle Valley hydrographic area, west-central Nevada. Water-Resour. Invest. Rep. 974191. U.S. Geol. Surv., Denver, CO.
- Maxey, G.B., and T.E. Eakin. 1950. Ground water in White River Valley, White Pine, Nye, and Lincoln counties, Nevada. Water Resour. Bull. 8. Office of the State Engineer, Carson City, NV.
- Milly, P.C.D. 1994. Climate, soil water storage, and the average annual water balance. Water Resour. Res. 30:21432156.[CrossRef]
- Milly, P.C.D., and K.A. Dunne. 2002. Macroscale water fluxes: 2. Water and energy supply control of their interannual variability. Water Resour. Res. 38(10):1206, doi:10.1029/2001WR000760.[CrossRef]
- Milly, P.C.D., and A.B. Shmakin. 2002. Global modeling of land water and energy balances. Part II: Land-characteristic contributions to spatial variability. J. Hydrometeorol. 3:301310.[CrossRef]
- Mintz, Y., and Y.V. Serafini. 1992. A global monthly climatology of soil moisture and water balance. Clim. Dyn. 8:1327.
- Mintz, Y., and G.K. Walker. 1993. Global fields of soil moisture and land surface evapotranspiration derived from observed precipitation and surface air temperature. J. Appl. Meteorol. 32:13051334.[CrossRef]
- Muller, R.A., and G.J. MacDonald. 2000. Ice ages and astronomical causes. Springer Praxis, Chichester, UK.
- NRCS. 2002. Wind erosion: Climatic factor C. p. 502-12502-13. In National agronomy manual. 3rd ed. USDA-NRCS, Washington, DC.
- Ol'dekop, E.M. 1911. On evaporation from the surface of river basins. Trans. Met. Obs. lur-evskogo, Univ. Tartu 4 (in Russian).
- Penman, H.L. 1948. Natural evaporation from open water, bare soil and grass. Proc. R. Soc. London, Ser. A 193:120145.
- Petit, J.R., J. Jouzel, D. Raynaud, N.I. Barkov, J.-M. Barnola, I. Basile et al. 1999. Climate and atmospheric history of the past 420000 years from the Vostok ice core, Antarctica. Nature 399:429436.[CrossRef]
- Pike, J.G. 1964. The estimation of annual runoff from meteorological data in a tropical climate. J. Hydrol. 2:116