HomeMy WebLinkAboutGround Water Contamination Model for Diagnosing 1985,SEP 3 0
A PRAGMATIC MODEL for DIAGNOSING and FORECASTING
GROUND WATER CONTAMINATION
T~MMO S. STEENHUtS*, MARJOLEIN VAN DER MAREL*, and STEVEN PACENKA**
Cornell University
Department of Agricultural Engineering* and
Center for Environmental Research**
Ithaca, NY 14853
Abstract
A mathematical management model for tracking the movement and fate of
a soluble chemical in the unsaturated and saturated zones is described and
tested against field and'laboratory data. The model is named MOUSE which is
the acronym for "Method Of Underground Solute Evaluation". MOUSE's soft-
ware runs on an IBM PC microcomputer.
MOUSE divides the surface and subsurface into three parts consisting
of the atmosphere, the unsaturated zone and the aquifer. The model simu-
lates these three parts using four linked submodels each of which can also
be used in a stand alone mode. They perform the following functions:
generation of synthetic climate patterns; calculations of moisture content
and fluxes in the unsaturated zone; simulation of degradation and movement
of the solute In the unsaturated zone and finally the movement of water and
the solute and attent~ation of pesticides in the aquifer.
MOUSE is a management model with intended use as a training tool for
proiessionals and sit,dents. It will.give them a better understanding of
ground water contaminattun. Another use envisioned is screening pesticides
under different environmental conditions. The predictions of the model
are compared against recharge measurements and pesticide concentrations in
the ,msaturated and saturated zones. In all cases model predictions are in
reasonable agreement with the observed values.
Proceedings of Practical Application of Ground Water Models
August 1984. National Water Well Association
Introduction
Ground water provides half of the U.S. population with their drinking
water. About 95 percent of the nation's rural population are similarly
dependent (U.S. Water Resources Council, 1978). Despite this substantial
reliance, we collectively know little about the relationships between our
act£vtties on the land surface and the quantity and quality of the water
which recharges this vital resource.
Relatively recent national concern about the negative impacts some
land uae activities have on ground water has spurred a great deal of
governmental activity relative to this resource. Numerous communities are
seeking ways to safeguard (or restore, in too amny cases) their local
ground water supplies. At the Federal level, the U.S. EPA has proposed to
reallocate existing funds to give ground water a much mere significant
place in their operations (U.S. EPA, 1984). What motivates this activity
is a recognition of the number of unsolved contamination cases coupled with
an acknowledgment of our historical state of ignorance about ground water
in general.
One major contribution to this new effort to understand aquifers for
management has been research into mathematical models which represent the
behavior of contaminants introduced on and into the land. Such models are
essential in understanding such normally invisible, expensive-to-observe,
and significantly complex systems.
This paper describes one such mathematical model, termed MOUSE, which
is the acronym for "Method Of Underground Solute Evaluation."
MOUSE packages mathematical algorithms describing the soil-water-
chemical system into a highly interactive computer program. The model
represents the complete subsurface domain, tracking the movement and
behavior of a chemical from its introduction at the land surface until its
possible eventual reappearance in a surface water body. The computer
program is designed to teach its user about key factors affecting the fate
of the introduced chemical, by allowing rapid changes of parameter values
and by portraying the simulated behavior of the system dynamically using a
graphical display.
The two main intended uses of this program are:
to train planners, public officials and students (as future decision
makers) about general aspects of the environmental processes which
affect the transport and fate of a chemical introduced at the land
surface; and
to help perform preliminary evaluations ("screening") of the potential
for ground water contamination by specific chemicals used under
different soil-climate-management regimes.
Mathematical Models for Water and Contaminant Accounting
MOUSE is part of a rich context. There are many published models
which trace soluble chemicals through the vadose zone to ground water.
There are also numerous models available which route chemicals through the
saturated zone of an aquifer. MOUSE does not break new ground in modeling
theory; rather, it packages selected, proven techniques into a comprehen-
sive and usable framework.
Types of Models
Addiscott and Wagenet's (1984) classification scheme divides the
continuum of models in several ways:
- deterministic versus stochastic;
- mechanistic versus functional; and
- rate versus capacity.
Deterministic models presume that the changes of state in the system
are uniquely governed by one given set of input parameters. Stochastic
models presuppose that the path followed will be uncertain and incorporate
this uncertainty explicitly in their structures. Much of the impetus for
development of stochastic models arises from the extreme soil variability
measured in the field (e.g., Bigger and Nielsen, 1976). Most deterministic
models cannot account for small-scale spatial variability. In theory,
stochastic approaches can represent this variation, but they remain
severely limited by the lack of quantitative knowledge about spatial
variability in specific field soils. Stochastic models and patterns of
field variability are the foci of a major wave of current research.
Mechanistic models consist of a network of equations which represent
fundamental physical, chemical, and biological processes in rich detail.
Conceptually elegant, these models can be difficult to apply to an
arbitrary field situation and also may require rather substantial computer
expertise and resources to operate (Wagenet and Rao, 1984). Functional
modeling approaches sacrifice this detail for sufficiency of results.
Simplifying assumptions in functional models limit their resolution and
precision to gain ease and economy of use. These models are diverse
because they embrace a variety of simplifying assumptions.
Rate models treat fluxes and changes of storage of water and chemicals
using differential equations solved over short time intervals. In some
cases they can employ analytical solutions. Capacity models use difference
equations to define changes in amounts of solute and water content over
finite intervals. These typically use time intervals of a day or longer,
whereas some numerical rate models must resort to time steps on the order
of minutes to follow the transient flow of infiltrating water.
By this classification system, MOUSE is a deterministic, functional,
rate model. Like all deterministic models, it allows the consequences of
variability in parameters to be assessed by way of repeated conditional
simulations. MOUSE fits into the functional category since its governing
differential equations of water and solute flux are simplified to allow
closed-form solutions. As a result the model can use large time steps,
thereby reducing the execution time significantly. The simplifications
introduced probably do not degrade computational accuracy beyond the limits
imposed by uncertainties in the field data with which it will be applied
and tested.
The choice of model characteristics depends heavily upon the intended
uses. Two broad types of usage include research and agricultural and
envlronn~ntal management. Research models aim at advancing the state of
knowledge about fundamental environmental processes, and thus tend to fit
into the mechanistic and rate categories. Management models sacrifice
detail for practicality of regular application and use. Simpler and more
functional models can be used in a ~rlder variety of situations, especially
where the less detailed results are sufficient to satisfy the manager's
needs. MOUSE's aims and simplifications make it a menagement model.
Design Criteria
Management models used by environmental planners ideally have the
following characteristics (modified after Haith et al., 1984):
1. Input data are generally available from standard published
2. The amount of computer expertise required is minimal.
3. Computational needs are sufficiently modest to permit long-term
simulations.
4. The algorlthn~ are based on sound physical principles and can be
used with a minimum of locally-~easured field data.
These criteria also apply in training settings where the aim is to
transfer appropriate conceptual models in users' minds. MOUSE's use of
dynamic graphical displays of simulated solute and water processes and Its
interactive control structure both assist in transferring such informa-
tion. These features work together to permit repeated trial runs with
different asau~ptions and quick recognition of their consequences.
Screening is facilitated by a general structure which permits widely
varying chemical, soil, and climate characteristics to be represented using
different parameter values. Most ~arameters needed may easily be estimated
from standard sources, such as SCS soil survey reports and National Weather
Service climatic records. By changing parameter values, conditions in most
of the continental U.S. can be represented.
The Overall Structure of MOUSE
For simplicity of construction and operation, the total model is
broken up into four loosely-connected submedels, each of which can also be
used by itself. A submodel obtains its parameters from one disk file,
obtains [ts boundary inflows from a different file created by its predeces-
sor submodel (if any), and saves its own results into a third file for use
by the next module in the chain. Each submodel also contains a screen-
oriented parameter editor and displays its dynamic results in graphical
form. The modular structure permits eubstitution of submodels using
different algorithms into any slot, and also allows mmasured data for
boundary inflows to be inserted when testing a submodel.
The submodels are:
- Synthetic Climate Generator
Vadose Water Balancer
Vadose Solute Transporter
Aquifer Water and Solute Transporter
The Synthetic Climate Generator generates daily precipitation, air
temperature, and soil temperature at a 10 cm depth. The Vadose Water
Balancer simulates water flows and storages within the unsaturated zone,
calculating evapotranspiration, moisture movement between vertical soil
strata, and recharge to the ground water. The Vadose Solute Transporter
superimposes solute movement on the calculated water fluxes, adjusted for
adsorption and desorption of the solute on the porous medium. Other
processes it represents include chemical diffusion and dispersion and
possible degradation due to biochemical action. The end reeult of this
submodel is a time series of solute addition to the aquifer. The Aquifer
Water and Solute Transporter simulates the movement and the degradation of
the solute in a vertical, two-dimenelonal cross-section of an aquifer.
Water movement is assumed to obey steady state Dupuit equations. Chemical
transport is based on water movement, linear adsorption/dasorption, and
first order degradation.
The following sections of this paper provide overviews of the four
submodels. In the final section the realism of each euhmodel is assessed
relative to the degree needed for the intended uses.
Synthetic Climate Generator
The synthetic climate generator calculates daily precipitation, and
air temperature. The cumbersome task of entering daily climate data can be
avoided with this submodel. The model uses monthly climate statistics from
historical records to produce daily weather patter'ne with statistical
properties similar to the historical pattern.
The input data needed are the average total precipitation amounts by
month; the average number of days per month with more than 0.25 mm of
precipitation; the average monthly temperatures for January and July; and
the number of years to be simulated. The climate data are eaeily obtainable
in summarized form from records compiled by the National Weather Service.
The output consists of a plot and file written to disk consisting of
daily precipitation, daily air temperature and soil temperature. The plot
is generated on the monitor during the simulation (Figure 1). The disk file
is used by the Vadose Water Balancer.
SIMULATED CLIMATE: Bridgehampton, N. Y.
oc
10-
0-
-10-
C
M
5-
O
J F M A M J J A S O N D
Air,Soil Temp. ,~:~
! ..~.. ..:"": ' '"': :~.- ..;.
Precipitation
Figure 1: Screen Display of the Synthetic Climate Generator.
Methods Used for Simulation
The simulation of the daily precipitation and temperature is based on
the work of Pickertng (Picketing, 1982; Haith et al. [984). This model was
chosen because of its simplicity and availability.
Daily precipitation is simulated using a Markov chain method which is
driven by conditional probabilities of precipitation given that a previous
day had rain or not based on historical records. This determines whether
any precipitation at all falls on a day. The equations for the probability
of precipitation are:
pdd - 0.1718 + 0.8642*WD/MD
(t)
pdw = WD*(I - pdd) (2)
where pdd is the probability of precipitation after a dry day, pdw is the
probability of precipitation after a previous wet day, W~ is the average
number of wet days in the month, and ~ is the total number of days in the
month. (WD can be taken from historical statistics). The coefficients in
the first equation are from Pickering's work. The occurrence of
precipitation uses these conditional probabilities and a random number from
a uniform distribution between 0 and 1.
The amount of precipitation falling on a "wet" day is also random,
based on an exponential distribution:
P = 0.0127 + I (3)
(PM - 0.0127)*1n(1 - U)
where P is the amount of precipitation (cm), PM is the mean amount of
precipitation on a wet day in the current month (from historical statis-
tics), and U is a random number from a uniform distribution between 0 and
1. The 0.0127 is half of the minimal amount of precipitation on a wet day.
Daily air temperatures are simulated by assuming that the average
temperature follows an annual sinusoidal pattern with random perturbations:
t[d] = u[d] + cl*(c[d-l] - u[d]) + /' ( I - cl*cl )*s[d]*V (4)
cl is a constant estimated by Picketing (1982) as 0.65; d is the Julian
day; t[d] is the day's air temperature; and V is a random number between 0
and I from a normal distribution; u[d] is the average temperature for
this day, computed by
u[d] Ujuly - UJan * { 1 + sin( (d - 100_,) ~) } + u]an (5)
2 180
which fits a sine wave between the average ~nthly temperatures of July and
January (Ujul¥ and Ujan respectively), s[d] is the standard
deviation of ~he day's temperature, based on
s[d] = 5.72 - O.122*u[d]
(6)
Soil temperature at 10 cm depth is computed as the average of the air
temperature of the previous 20 days.
The Vadose Water Balancer
The Vadose Water Balancer simulates the movement and retention of
water in the unsaturated zone above the gronnd water table for a given
crop, soil and climate. The submodel accounts for snowfall, snowmalt.
rainfall, runoff, infiltration, evaporation, transpiration, downward
movement (percolation and recharge) and changes in snow and liquid water
storage.
The input consists of daily precipitation and temperature, which is
obtained from a disk file created by the Synthetic Climate Generator. Other
climatic input consists of potential evaporation. The user may supply
observed data. The evapotranspiration may also be simulated by supplying
the program with the average monthly potential evapotranspiration for
January and July.
Soil hydrologic parameters are also required. For each of the three
zones the saturated hydraulic conductivity and the moisture contents at
saturation, wilting point and when air dry, should be entered. To compute
the overland flow on soils without a hardpan the curve numbers for the SCS
runoff equation need to be known. Finally the depth of each of the three
zones (including the depth to the ground water) has to be specified. The
soil hydrologic data may be found in the more recent S.C.S. soil surveys.
The curve numbers are tabulated in most hydrology textbooks. Pan evapora-
tion is measured at several locations in each state.
The output consists of a disk file containing the amounts of infiltra-
tion, the percolation and the recharge to the ground water. These amounts
are summarized over a user-defined time interval. The file is used by the
Vadose Solute Transporter. During the execution of the program the
recharge, infiltration and evapotranspiration are shown graphically on the
monitor. As an example Figure 2 gives a printed copy of what appears on
the screen.
cm Water Balance Summaries
J F M A M J J A S O N D
incremental
2
1
°
cumulative
150
~0 -,
0 .~-
Inflow [] Evap ~ Recharge ~
Figure 2. Display of Screen of the Vadose Water Balancer.
Methods Used for Simulation
The theory on which the Vadose Water Balancer is based is treated in
som~ detail because of the importance of mass flow in the transport of
chemicals in the soil profile.
To meet the design criteria for n~nagement models as set forth earlier
-- permitting long-term simulation and requiring that input data is
generally available from published sources~ the description of the soil is
kept simple. The soil profile is broken up in four zones as shown in
Figure 3. The first zone is the active evaporation zone (ZONE l) and has a
depth which changes with time in response to the seasonal plant growth,
Below the evaporation zone and at a constant depth is the boundary between
the upper and lower transmission zones (ZONES 2 and 3. respectively). The
distinction is made for computational reasons as well as to provide sensi-
tivity to properties of different soil horizons. Below the lower (unsatu-
rated) transmission zone is the saturated zone or aquifer. The unsaturated
movement in zones 1 through 3 is vertical while the flow of water and
solutes in the saturated is mainly in the horizontal direction.
Zone I
] Active Evopo-
'~lraznSo~r°ti°n J
Zone 2
Upper tronsmission zone
Zone ~3
Lower Ironsmis$ion zon~ ·
Soturoted zone or impermeable Ioyer
Figure 3. Division of Soil Profile for Vadose Water Balancer.
To facilitate long-term simulation a variable time step is used. The
time step is short (less than one day) when the profile is wet and the
fluxes are high. The time step is large when the profile is dry and fluxes
are small. The meximum time step is set in the program and depends on the
interval over which the fluxes are summarized and is always smaller than
two consecutive rainfall events.
The main components of the Vadose Water Balance are shown in Figure
4. The overall structure is in essence very similar to other functional as
well as mechanistic models such as ARM (Donigan and Crawford, 1976),
DRAINMOD (Skaggs, 1979), CREAMS (Kuisel, 1980), ANSWERS (Beasley, 1977) and
CPM (Haith, et al., 1984) but the methematical basis of most of the
components is very different from these earlier models.
The water mass balance equations for each zones are:
Zone 1: dlt Am! = Rt + It + Mt - Qt - Et -Plt (7)
Zone 2: d2t Am2 = Pit- P2t (8)
Zone 3: d3 Am3 = P2t - P3t (9)
where dlt, d2t and d3 are the thicknesses of each zone (em). dlt
and d2t vary with season, but d3 is constant; see figure 3. Ami is
SHALLOW
PECOLATION
MOISTU~[ · INT[RF~OW
ZONE 2
RECHARGE.O WHEN THE SOIL HAS
AN ~MPERMEAgLE LAYER
INTERFLOW.O WHEN THE SOIL IS
HOMOGENEOUS AND
DEEP
Figure 6.
I DEEP
PERCOLATION
I SOl L
MOISTURE · INTERFLOW
ZONE
IRECHARGE
GROUNDWATER OR
IMPERMEABLE LAYER
Hydrologic Model Structure.
the change in moisture content per unit volume in layer J. Rt is the
rainfall (tm), It is irrigation (tm), Mt is enowmelt (tm) and Qt is
runoff (tm). Et is evapotranspiratton (tm), which occurs only in Zone
and Pjt is the amount of water conducted downward through Zone J.
Snowmelt is computed by a degree-day factor (U.S. Army Corps of
Engineers, 1960). Irrigation is applied whenever soil moisture falls below
some specified limit provided that the elapsed time since irrigation is
greater than or equal to the length of the irrigation cycle.
In applying equations 7-9, water inputs to each soil zone, which
include rainfall, irrigation and snowmelt for Zone 1, are assumed to infil-
trate at a ~uch higher rate than that of evapotransptration and percola-
tion. In this model water is added at the beginning of the timestep and
then allowed To redistribute slowly during the timeetep as a function of
moisture content.
Unlike precipitation, irrigation, snowmelt and runoff, evapotranspira-
tion and percolation are modeled as continuous processes occurring
10
throughout the time step. Because simultaneous evaporation and percolation
occur in Zone I while only percolation occurs in Zones 2 and 3, a separate
mathematical treatment is required for the active evaporation and transmis-
sion zones. Thus, water movement in Zone 1 is discussed separately from
Zones 2 and 3.
Evaporation and percolation are slow continuous processes. 8otb rates
depend on the volumetric moisture content of the soil. When the moisture
content is greater than or equal to mrsd (a meisture content below which
no do~raward flow occurs in the root zone), evapotranspiration occurs at its
potential rate, Pg. For the model, potential evaporation was estimated by
a simplified Penman Equation (Merva and Fernandez, 1982; Norman, 1984) or
if no daily evaporation data are available by
PE: ( PEJul - PEJan) { I + sin (<d - 100) x ~ )} + PE%an <10)
2 180
where PEjan and PEju1 are the average potential evapotranspiration during
January and July respectively.
At moisture contents above mred percolation is significant; below
Mred moisture content of redistribution, hydraulic conductivity is so
small that downward flow is assumed to be zero. Stated mathematically, the
relationships for redistribution following rain, malt, and irrigation
inputs to the root zone are as follows:
d1 dml = -PE - k(ml), for m1 > mred (11)
dt
dI dml ~ -E(m1) for m1 ~ mred (12)
dt
where m1 is the updated moisture content in Zone I taking into account the
water added by precipitation, snowmalt and/or irrigation, k(m1) is the
unsaturated hydraulic conductivity at moisture content mi, PE is the poten-
tial evapotranspiration and E is the actual evapotranspiration rate.
Equation Il is based on the assumption proposed in Baver, et al.
(1972) that gravity dominates over matrlc forces making the hydraulic
potential gradient close to unity. For the unsaturated conductivity in
equation 11, an exponential relationship was used because of the necessity
to obtain closed form solutions; e.g.
k(m1) - ks exp[~(&-l)]
(13)
9 is the reduced moisture content, given by:
e = (14)
ms - md
In the above equations k(m) is the hydraulic conductivity at moisture
content m (cm/day), ks is the saturated hydraulic conductivity (cm/day),
11
md is the air dry moisture content, and K is a constant which is based on
the work of Bresler, et al. (1978) has an universal value of 13.
Assuming that the moisture content decreases uniformly with depth when
the soil is drying out, equations 11 and 13 can be combined to solve for
t', the time it would take for a given moisture content, ml, to be reached
starting at an initial moisture content for Zone I of mil.
t' = -d1 [ml - md - ms - md In [P~ + ks (exp(-m))(exp(~01))]] + t*(15)
PE PE ~
where t* is an integration constant determined by the initial conditions
expressed by the equation below:
t*~ dx [~il - md ms - md
PE PE ~
In + ks (16)
When moisture content was smaller than mrsd the downward flux was
assumed zero amd the equatlon for the actual rate of evapotranspiration
assumed a linear relationship between the actual evapotransptration rate
and the moisture content between the moisture where soil water first limits
evapotransptratton, mred, and the moisture content where the soil is so
dry that evapotransptration ceases, mwp. Note, mwp is the permanent
wilting point and not necessarily equal to md, the air dry moisture
content. Stated mathematically,
w!
E - ~ PE (17)
W
where w1 = dl(mI - mwp), or moisture in layer at any time (tm),
w* ~ d1 (mrsd - mwp).
Rearranging equation 17, substituting it in equation 12 and solving
the differential equation, the evapotranspiration under water limiting
conditions is obtained.
Et = dl(mil - m~p) (1 - exp ( - PE At)) (18)
dl(mred - mwp)
The corresponding moisture content is given by
dlmtl - Et
mi = (19)
dl
Beneath the active evaporation zone there is no plant uptake of water,
and the flux of water is equal to the hydraulic conductivity (assuming unit
hydraulic gradient).
12
dj dm3 = -k(mj) (20)
dt
where the subscript J refers to Zone 2 or 3. The relationship for k(m)
(equation 13) is substituted into equation (20) and the resulting equation
is solved by separation of variables to yield an expression for the water
content as a function of time and initial water content.
msj - mdj [ksj exp(- )At + exp(- eij)]
mj = mdj ~ dj(msJ - mdj)
(21)
Moisture content at the end of the time interval can be found by
substituting the length of the time interval for At in equation (21). In
the unsaturated zone below the root zone, the downward flow of water Pit
is simply equal to the change in moisture content over the interval ti~es
the thickness of the soil zone.
rjr ~ dj{mj - mlj)
(22)
where mij is the initial moisture content in Zone J.
The overland flow is calculated with the modified SCS curve number
method. In addition runoff was assumed to occur on the shallow soil when
the profile was saturated (i.e., when no water can infiltrate into the
soil). Thus, when the soil was saturated the curve number Jumps to almost
100 (i.e., all precipitation becomes runoff).
The Vadose Water Balance also has a component that calculates inter-
flow in soils with a hardpan.
13
The Vadose Solute Transporter
The Vadose Solute Transporter simulates the movement and the attenua-
tion of chemicals within the unsaturated zone above the water table. The
submodel includes the processes of first order degradation, dispersion and
diffusion as well as mass movement.
The input of the model consists of a file generated by the Vadose
Water Balancer. This file contains not only the fluxes of water such as
infiltration, percolation and recharge but also common hydrologic soil
parameters that are shared between the two models such as the equilibrium
moisture content for the three zones as well as the depth of each zone.
Additional data required for this submodel are the application dates,
amounts and incorporation depths of the applied contaminants. Adsorption
parameters for the root zone and lower transmission zones are also needed.
Finally, the degradation rates as a function of depth from the surface are
required to run the submedel.
The output consists of a dynamic display of the solute concentration
pattern over the depth of the vadose zone (Figure 5). The contributions of
the chemical to the ground water are written to a file on the disk for use
by the Aquifer Water and Solute Transporter.
CONTAMINANT MOVEMENT FROM
THE LAND SURFACE TO THE AQUIFER
Year
2
Day
293
Conc.
at
aquifer
362
d
e
6
Land surface
0 (solute, g/cm ha) 300.0
14
Figure 5. Screen Display for the Vadose Solute Transporter.
Methods Used for Simulation
The simulation of the movement of the solutes through the soil profile
is complicated by the adsorption-desorption and decay of pesticides as
well as the simultaneous convection and dispersion. Mechanical modeling
approaches solve the exact differential equation describing these processes
with a finite difference or finite element technique. This technique is
time consuming considering the intended uses. In this submodel the simula-
tion is simplified by decoupling the transport from the degradation. In the
discussion in the following section this two step approach is maintained.
The movement and dispersion are discussed first followed by the attenuation
of the chemical.
The solute movement is superimposed on the water movement and is
partly based on simplifying assumptions proposed by Gardner (1965) and Rao
et al. (1976, 1981). The depth of the solute peak is calculated first and
then the dispersion of the s61ute around the peak location is obtained
with an analytical solution to the convective-dispersive transport equation
as original proposed by Scheidegger (1960) and used by Walter (1974).
The program keeps track of the peak location for each of the pesticide
messes that are applied. The depth of the solute peak may be found as
(Gardner, 1965; Rao et al.,1976);
zt = zt-1 + AP (23)
meq + Kd P
where zt_1 is the location of the pesticide at time "t-l", zt
is the location at time t, meq is the moisture content at equilibrium and
may be found from equating the average yearly flux to the unsaturated
conductivity in equation 13, AP is the quantity of water flowing past the
peak and is equal to the average of the in-out going fluxes of the zone in
which the peak is located. Kd is the adsorption partition coefficient
(based on a linear adsorption isotherm) and p is the density of the soil.
The solute dispersion around the peak is found as
M(x) = 0.5 Mo erf x
2/ 8 zt
(2~)
where M(x) is the mass of the chemical between the peak and a distance x
from the peak. Mo is the total amount of pesticide. ~ is the
dispersivity divided by the average downward velocity of the peak.
Walter (1974) found that a value of three for ~ fits most conditions, zt is
the distance of the peak from the surface. Thus z ia defined in the
Euclidian coordinate system. The x coordinate is used in a moving coordi-
nate system with an origin at the peak location. This coordinate system
moves at the velocity of the peak as calculated with equation 23.
The first step in the computational procedure is to obtain the flow of
water through each zone simulated by the Vadose Water Balancer. The flux is
then substituted in equation 23 to find the solute peak location. On the
day of pesticide application, the pesticide peak is at a depth equal to the
incorporation depth, where it remains until the first precipitation or
irrigation event. With the first infiltration event, the peak location is
determined by equation 23. When the pesticide concentration with depth is
displayed, the position of the band is substituted into equation 24 t.o find
the distribution of the pesticide around the peak. After distribution
several pesticide applications may come to overlap. Contributions of ail
pesticides at each depth interval are summed and displayed. The process is
shown schematically in figure 6. When the pesticide band is near the
surface, it is possible that equation 24 will give a positive mass of
pesticide above the surface. When this happens the portion above the
surface is repositioned equally between the peak and the surface.
The preceding equations do not account for degradation of the
pesticide. As the pesticide is in the soil, however, the concentration
decreases due to break down of the molecule. The rate of break down
decreases with depth. In the root zone the degradation is higher than in
the lower transmission zone. In the root zone, breakdow~ is microbial in
character while in the subsoil hydrolysis takes over in importance. Usually
microbial degradation is much more effective than hydrolysis.
To account for varying rates of degradation with depth the pesticide
mass for each application is divided up in portions with the same percent-
age of equal mass in each. Each pesticide mass is divided up in a maximum
of 10 portions. The division occurs shortly after application. During the
solute Journey the pesticide mass remaining in each of these original
partitions is being tracked. This is accomplished by calculating the
location of the partition and then finding the degradation rate for the
15
Figure 6.
PESTICIDE CONCENTRATION
Schematic of Summing up of Separate Pesticide Bands.
pesticide at that location. Assuming first order decay, the mass of
pesticide in any partition at any time is given as
Mpt = Mpo exp[-( kxt! + k2t2 + .......... + kntn)] (25)
where Mpt is the mass of the pesticide in partition p after elapsed time
t, MpolS the original amount of pesticide present in partition p at time
of application, tj is the length of time the pesticide in the partition
has a decay rate of kj. The degradation rate kj is found by averaging
degradation rates over which the portion extend~ itself.
In the beginning of the program the degradation rates for each of the
partitions associated with various peak locations are calculated and stored
in memory as a lookup table. The advantage is that during the simulation,
there is no need to distribute the pesticide around the peak for each time
step to find the degradation rate, k, for each partition; but, by finding
the nearest peak in the lookup table, the degradation rate for each of the
partitions is known too. This process significantly reduces the execution
time.
As stated before, the pesticide input to the ground water is is one of
the output variables written to a disk file. To calculate this an imaginary
line is placed at the depth of ground water table. The soil profile is
thought to extend below this imaginary line. By setting the pesticide
degradation rate to zero below the imaginary line, and keeping track of the
total pesticide amounts below the line, the quantity of pesticide passing
into the ground water is equal to difference in total pesticide amounts
below the line.
16
Aquifer Water and Solute Transporter
The Aquifer Water and Solute transporter simulates the movement of
water and the movement and degradation of a solute in a vertical, two-
dimensional cross-section of an unconfined aquifer. The water flow is
based on steady-state, Dupuit fl~w assumptions. The chemical movement is
superimposed on the water flow and includes first-order degradation and
linear adsorption/desorption.
The input consists of recharge and solute inputs aa a function of time
and is obtained from a file on the disk, created I~ the Vadose Solute
Transporter. Other input data required are the characteristics of the
aquifer consisting of the width of the cross-section, elevation of water at
outflow point, locations of flow path starting points, and the saturated,
hydraulic conductivity in the horizontal and vertical direction. The
chemical degradation and adsorption characteristics are also needed.
The output consists of solute concentration over time and along one
flow path and is displayed on the screen. An example is given in Figure 7.
t7
SOLUTE TRANSPORT IN THE SATURATED ZONE
Solute concentration (ug/l}
Flowpath 1 2' 0--['--~~ ~~ ~--400
I sochrone
Clock: ~
9.9yr. th
Cell:
Conc.: (rr -----:--
0.0
Max age is 22.2 years
0 distance (m) 1000
Figure 7. Screen Display for Aquifer Water and Solute Transporter.
Methods Used for Simulation
The elevation of the water table in a two-dimensional, steady-state
ground water flow pattern is calculated using the hydraulic approach of
Dupuit and Forchheimer as expressed by Gelhar and Wilson (1974). The
boundary conditions (see figure 7) include a hollow boundary at the left
and an impervious layer at the bottom, uniform recharge along the top, and
a constant-head outlet at the right. The equation below gives the height
of the water table at any location:
18
-- /{ h°2 + q y (2 L - y)
head(y)
)
(26)
ks
where head(y) is the elevation of the water table at position y measured
from the impervious left edge, ho is the constant head at the right edge,
q is the recharge rate averaged over the period of simulation, Ks
conductivity, and L is the width of the cross-section from the left to the
right edges.
Flow directions and speeds are represented by plotting streamlines and
isochrones. Flow paths are plotted using the equations:
y(t) = L + (Yo - L) ex~ [~otI
(27)
z(t) = head(yo) exp [~ot]
(28)
where q, ho, and L are as defined above, t is a "travel time" measured from
the entry of a "tracer particle" at the water table, Yo is the location
where the path intersects the water table, y(t) is the horizontal location
of the particle at time t, and z(t) is the vertical location of the tracer
particle at the same time. A whole path may be plotted by varying the
parameter t up to the point where y(t) is over the right boundary. An
isochrone Joins the points on adjacent flow paths that correspond to the
same t value. Isochrones are plotted such that they reflect the solute's
lagged movement (due to adsorption) rather than the water's movement.
Chemical movement and degradation are calculated for successive time
intervals. It is assumed that there is no chemical in the system at the
start and that there is no dispersion of the chemical as it moves.
The equations used to represent the effects of chemical adsorption and
degradation are equivalent to those used in the Vadose Solute Transporter.
They are far simpler due to the omission of dispersion. A "band" of
chemical which enters in one time interval thus stays together until it
exits from the saturated system at its downstream end. Degradation reduces
the m~ss in the band depending only on the amount of time since the band
entered the saturated zone.
Realism
A mathematical model must be realistic to be used in practical appli-
cations. The degree to which it should simulate the reality depends on the
intended use. If a model mast reproduce day to day variations in measured
field observations for research purposes, the simulation should represent
the reality for the particular area in detail. However, a model used for
for training, diagnosis of contamination problems in the soil and/or
screening out combinations of soil chemical and climate posing an
environmental problem ~st be qualitatively realistic. More important the
latter kind of the model (such as MOUSE) must perform reasonably over a
wide range values, rather than be specialized within tight limits.
In all cases, comparisons of a model's predictions with field measure-
ments provides an indication of consistency and accuracy. Recognizing the
variability in field measurements and the limited capacity of management
models to represent this variability, no attempt has been made to compare
statistically the observed versus the real data. Instead, for various
experiments, the predicted and observed data are shown and it is up to the
user to judge how mach confidence she/he is willing to invest in this tool
MOUSE for a given application.
Great variability in solute and water flow has been observed in the
vadose zone. For a summary of these observations see for example Rao and
Wagenet (1984). However, in predicting groundwater concentration the
spatial distribution of solute movement in the vadose zone is of less
concern. Water resides long enough under the field that the variation in
arrival time of a pesticide mass (as affected by the variability in the
vadose zone), are negated by the the dispersion and diffusion in the
saturated zone, thus resulting in averaging of the many variations intro-
duced in the unsaturated zone. Temporal trends on the order of the travel
time in the saturated zone under the field rm~st be preserved in the simula-
tion results. Consequently, management mudele, which are concerned with
the eventual fate of pollutants in the aquifer, predict areal average
concentrations and fluxes based on their intended use.
The submodels are tested separately. For each model physical measured
data are used. If they do not exist data are estimated based on data avail-
able for other similar areas. In all field cases historical weather records
Climate Generator
Picketing (1982) tested the algorithm for precipitation and air
temperature for 25-year periods at weather stations in New York (Aurora),
Iowa (Ames), and Georgia (Athens). The validation is described in more
detail by Haith and others (1984). The soil temperature equations have
been validated by Steenhuis (1979). Test results for the Climate Generator
submodel are presented in these references.
The Vadose Water Balancer
This submodel was tested against precipitation and recharge or deep
percolation (measured at the 90 depth) data for a site on the grounds of
Cornell University's Horticultural Research Farm located 4 km north of
Riverhead, N.Y. and 1 km south of the Long Island Sound. The farm lies on
a glacial outwash plain with soils typical of the region. Soil at the site
is a Haven sandy loam. The soil profile has a sharp distinction between
the top 30 cm and the rest of the profile. The plow layer has up to 40%
clay and silt and has been compacted by traffic. The soil below 30 cm is a
yellowish-brown, gravelly sand with less than 10% silt and clay. The sand
has a high saturated conductivity of approximately 10 m/day.
19
The direct and the closure methods were used to measure recharge. The
direct method consisted, first, of measuring simultaneously the matric
potential and moisture content at 90 and 120 cm depths. By relating
moisture content to hydraulic conductivity the vertical flux of water could
be calculated directly with Darcy's Law. The closure method calculated
recharge from the hydrologic budget equation with evaporation computed from
the mtcrometeorologic data. A detailed description can be found in
Steenhuis, et al. (1984).
Recharge at the site was measured by both methods in 1980 from January
through June and during November and December. Precipitation and minimum
and maximum temperature were measured at the site. Average monthly cloudi-
ness index and wind speed were taken from the airports in New York City.
Soil parameters used in the simulation model are given in Table I and were
measured at the site itself (Steenhuis, et al., 1984).
The simulated recharge as well as the observed recharge with the
direct method are compared in figures 8 and 9. The results of the direct
method were assumed to be more accurate because the recharge estimate does
not contain the errors in the measurement of individual components of the
Figure 8. Precipitation and Predicted and Observed Recharge on a Haven
Sandy Loam from January through June 1980.
hydrologic budget as with the closure method. The closure method measured
a 10 to 20% higher recharge than the direct method except for January and
December when the closure method was much higher than the direct method.
As can be seen from Fig. 8 and 9 the simulated recharge was generally
slightly higher than the directly measured recharge. The difference
between the direct method and the simulated recharge on a month by month
basis was always less than the difference between the direct and closure
methods. The model therefore gives good results of recharge based on
physical measured data.
2O
21
Table 1. Soil Parameter Values for Water Balance.
Zones
1 and 2
Zone
3
Saturated hydraulic conductivity (m/day)
Saturated moisture content (cmS/tm3)
Moisture content at wilting point (cmS/tm~)
Moisture content at which actual evapo-
transpiration is smal}er ~han potential
evapotranspiraton (cm /tm~)
Dry moisture content (cm~/cm~)
Depth of zone (m)
0.6
0.50
0.08
0.24
0. O05
0.26
10.0
0.45
0.005
0.64
*Not applicable.
Figure 9.
"1
..... ..=.-
Precipitation and predicted and observed recharge on a Haven
Sandy Loam during November and December 1980.
The model was also tested for watersheds with soils underlain by a
hardpan at shallow depth in upstate New York. Predicted and observed flows
from watersheds up to 100 sq km agreed well (Norman, 1984; Steenhuis et
al., 1983).
Vadose Solute Transporter
The Solute Transporter was tested against a laboratory soil column
experiments, a mechanistic model and soil cores taken under two potato
fields on Long Island.
[n the laboratory columas the distribution of radioactivity with depth
was measured from radio-labelled butylate, alachlor, and metolachlor after
they were leached with 15 cm of water. The application rates ~ere equiva-
lent to 4.5. 2.5 and 2.5 kg/ha, respectively, and the adsorption partition
coefficients were 3.14, 1.8 and 1.98 cc/g respectively. The soil used was a
Keeton sandy loam. More information may be found in Spillner et al.
([983). To test the Vadose Solute Transporter the 15 cm of recharge as
well as pesticide and soil characteristics were supplied as input param-
eters to the model. The comparison of predicted and observed distribu-
tions for the three pesticides are shown in figure 10. A satisfactory
agreement was obtained for these pesticides. Although metolachlor had
approximately the same characteristics as alachlor, its observed
distribution with depth was different. The model underestimated slightly
the movement of matolachlor.
P~STICID~ CONC£NTRATION (MG/KG)
O,S 0,9 1.5 2.1 2.7
0
0
0
BUTYLATE
22
j 0.3 0.6
~. 20
ALACHLOR
0.9 1.2 1.5
10
20
O 0.~ 0.6
0 METOLACHLOR
0.9 1.2 1,5
Figure 10.
SIMULATED PESTICIDE CONCENTRATION
0 OBSERVED PESTIC[DE CONCENTRATION
A comparison of Predicted and Observed Distribution of
Butylate Alachlor and Metolachlor.
The second test of ~OUSE compared the model against the predictions of
a mechanistic finite difference model developed by Intern (1980). The model
was used to predict aldicarb leaching under a potato field on Long Island.
The field was located in Cutchogue, NY, which is within a intensive potato
growing area on the North Fork of Long Island. Biannual aldicarb applica-
tions were made from 1977 through 1979 at a total rate of 3200, 4000 and
5000 grams of active ingredient per hectare per year for the three years.
Five cores were taken for sampling in December, 1979. The soil was compa-
rable to the soil on the Horticulture Experimental Farm as described
before. The depth to ground water was 2.5 m. Intera calibrated their model
and found that a half life of 77 days in the whole profile and an adsorp-
tion partition coefficient of 0.16 cc/g best fitted the observations iu two
of the five cores. In MOUSE we substituted the same parameters (including
the same weather information from the Greenport Weather Station). As can be
seen from figure 11 we obtained almost a complete match between the
aldicarb contents simulated b~ the t~o models. This is not surprising as
both models are based on the same original mathe~atical formulations. It
also implies that one can substitute simpler functional models for mecha-
nistic models when the system's parameters are not very well known, The
Mouse simulation run was done on a micro computer in less than 15 minutes
including the time it took to enter the data, The finite difference
approach required a larger computer and considerable expense,
23
0
,,
ALDICARB CONCENTRATION ( ~G/L)
20 ~0 60 80 100 120 1~0
~L2. O~ o I
2,aJ
ol
0 I
O SIMULATED ALDICARB CONCENTRATION BY INTERA (1980)
SIMULATED ALDICARB CONCENTRATION BY ~0US[
Figure 11. Comparison Between the Intera (1980) and MOUSE Simulation of
Aldicarb Leaching.
The third and most rigorous test was for soil cores under potato
fields taken in December 1983 at t~o locations on Long Island. The cores
were analyzed for aldtcarb content b~ the Union Carbide Corporation (Hughes
and Porter, 1984). The first location was near the site of the Intera
study. The depth of the groundwater was approximately 10 meters. The
second site was near Wading River with a depth to groundwater of 31-33
meters. In both cases the soil was very similar to the Haven soil at the
experimental farm. The aldicarb applications were made from 1975 through
1979 after which the pesticide was withdrawn from the Long Island market.
During each year an application of 3000 g AI/ha was made at planting. In
addition during the last three years a later application of aldicarb was
made during June at a rate of 1200 or 1800 g AI/ha.
The greatest uncertainty in modeling aldtcarb is its rate of degrada-
tion in the vadose zone and saturated zone. Unfortunately, laboratory and
field experiments giving the rate of degradation are few, and the results
are inconclusive. A previous review summarizing studies on aldicarb degra-
dation in a variety of soils showed ranges in half lives from 1 to 231 days
(Kain and Steenhuis, 1984). In our simulation a 45 day half live was used
for the root zone. This was increased to 90 days for the next 50 cm to
reflect the decrease of microbial activity with depth (Pacenka and Porter,
1981). For the remainder of the vadose zone, a half life of I0 years was
used. This was obtained from work of Lemley (Personal communication, 1983)
and analysis of aldicarb content data from a deep soil core taken at the
Long Island Horticulture Research Laboratory (unpublished data).
The adsorption portion coefficient was set at 0.10 in the root zone
and zero in the subsoil. The numbers were calculated from a linear regres-
sion relating adsorption coefficients to soil organic matter derived from
the results of four field and laboratory experiments. (Kain and Steenhuis,
1984).
Based on the above pesticide properties and soil characteristics
similar to those used for the recharge study, the pesticide distributien
with depth was simulated. For the field on the North Fork predictions
agreed with that what was observed: i.e. all aldicarb had moved out of the
unsaturated zone (into the ground water). For the field near Wading River
the model and the observed concentration started near the same depth.
(Figure 12). The variation in concentration with depth was larger for the
predicted than the observed results. The dispersion parameter which deter-
mines the spread of the band was based on experiments done with laboratory
columns. It has been noted that in ground water the dispersivity tends to
increase with distance. Thus the simulation was repeated with a ~ value
(equation 24) increasing linearly fro~ 3 at 2 m to 10 at 30 m. The result-
~ng prediction fitted slightly better.
Overall, as long as no better information about the subsoil parameters
is available, the simulation results have to be considered acceptable.
The Aquifer Water and Solute Transporter
The aquifer n~del was tested for two pesticides for a transect on the
North Fork of Long Island, near Depot Lane. In 1981 the Suffolk County
24
ALDICARB CONCENTRATION (~,G/I(G)
8 12 .. 16 20
2~
25
0
0
--- DISPERSION PARAMETER INCREASING NITH DEPTH
DISPERSION PARAMETER CONSTANT WITH DEPTH
0 OBSERVED ALDICARB CONCENTRATION ~N T~O CORES
Figure 12. Comparison of Simulated and Measured Aldtcarb Content of a
Potato Field near Wading River.
Department of Health Services drilled a series of geological test holes and
observation wells along Depot Lane in orde~ to obtain information about how
aldicarb and other contaminants move through the North's Fork aquifer. The
transect runs from the Long Island Sound to the Cutchogue harbor and runs
parallel to the ground water flow. Along the transect the water table
height, aldicarb and dichloropropane concentrations were measured in a time
span from spring 1981 until spring 1982. Additional information about the
sampling procedure and results can be found in Baler and Robbins (1982).
The prediction of the water table height is first compared with the
observed water table height. The hydraulic conductivity was taken as
90 m/day horizontal with a horizontal-to-vertical anisotropy ratio of 5.0
(Baler and Robbins, 1982). The depth of the aquifer was taken as 33 m near
the Long Island Sound and 40 m near the Cutchogue Harbor. The predicted
and observed water table height are compared in figure 13. The predicted
water table height is a "steady state" height based on the average recharge
during a nine year simulation period from 1975 - 1983. The variation in
observed ground water was the difference between the readings in the summer
of 1981 and following wLnter. As expected the ground water levels near the
salt water interface were slightly off. The simplified flow equation
employed with MOUSE does not account for the interface between salt and
fresh water and the upward curving of streamlines near the interface. This
causes a higher loss of energy than predicted with the Dupuit assumptions
and results in a higher water table than predicted.
I I
1.0
0,5
0
LONG
Figure 13.
1000 2000 3000 4000 COTCHOQuE
HARBOR
DISTANCE (M)
Observed and Simulated Water Table Height for Depot Lane
Transect.
For the aldicarb simulation pesticide application data along the
transect by Trautmann et al. (1983) were used. The degradation and adso~p-
tion constanta for aldicarb were the same as used in the Vadose Solute
Transporter test near Wading River. Soil parameters were equal to those
used in the Vadose Water Balancer validation. The model simulations in
figure 14 revealed that the depth over which aldicarb was found in the
ground water resembles closely the calculated distribution. The model
predicted a penetration of approximately 6 meter below the ground water
table. Aldicarb had only been found below this depth in a few locations by
1981 (Figure 14). The comparison between predicted and observed concentra-
tions was hampered because the exact locations and depths of aampled wells
were not known. In general the model predicted concentrations ranging from
10 ppb to 50 ppb with a very few exceptions when the concentratione
exceeded 100 ppb. In general this was at the same locations where high
concentrations were observed in the ground water.
The final test involved a comparison of calculated times of ground
water travel on the North Fork since the early 1950's. Dichloropropane, a
soil fumigant which has been detected in the ground water in the area was
first used in 1951 in the Depot Lane transect and can serve as tracer.
Assuming that as in case of aldicarb the adsorption partition coefficient
is zero, the observed distribution of dichloropropane (figure 15) is
roughly consistent with the ground water travel times predicted with the
Aquifer Water and Solute Transporter. This lends additional confidence to
the applicability of MOUSE to this kind of setting.
26
27
( ,,T,A } ~0 ~ ,T,¥/~'~~r
28
(,~,J) i(o z &¥A.1"l
Conclusions
Several tests of simulation results against observed field and
laboratory data all indicate that MOUSE can reproduce the general measured
patterns of subsurface water and selected pesticide concentrations quite
well within the limits of knowledge of field conditions. This degree of
realism is adequate for planning and management applications and ensures
that MOUSE users will obtain simulation results that are quantitatively
reasonable for the chemicals examined.
As with most models of the functional, management type, MOUSE's
realism is coarse. Only with fine- tuning of parameters can it reproduce
short-term small-scale details of a particular field or laboratory
experiment. This suggest that our intended restriction8 -- to training and
management uses -- are appropriate.
29
References
Addtscott, T.M. and R.J. Wagenet. 1984. Concepts of Solute Leaching in
Soils: A Review of Modelling Approaches. Journal of Soil Science (in
press).
Baler, J.H. and S.F. Robbins. 1982. Report on the Occurrence and Movement
of Agricultural Chemicals in Ground Water: North Fork of Suffolk
County. Suffolk County Department of Health Services, Bureau of Water
Resources.
Baver, L.P., W.H. Gardner and W.R. Gardner. 1972. Soil Physics. 4th Ed.,
John Wiley, NY, NY.
Beasley, D.B. 1977. ANSWERS: A Mathematical Model for Simulating the
Effects of Land Use and Management on Water Quality. Unpublished Ph.D.
Thesis, Purdue University, Lafayette, IN.
Biggar, J.W. and D.R. Nielsen. 1976. Spatial Variability of the Leaching
Characteristics of a Field Soil. Water Resour. Res. 12:78-84.
Bresler, E., D. Russo and R.D. Miller. 1968. Rapid Estimate of Unsatu- rated HYdraulic Conductivity Function. Soil Sci. Soc. Am. J. 42:170-172.
Donigan, A.S. and N.H. Crawford. 1976. Modeling Pesticides and Nutrients
on Agricultural Land. EPA 600/2-76-043. U.S. Environmental Protection
Agency, Athens, GA.
Gardner, W.R. 1965. Movement of Nitrogen in Soil. IN: Bartholomew,
W.V. and F.E. Clark, editors: Soil Nitrogen. Agronomy, A Series of
Monographs, V. 10. American Society of Agronomy. Madison, WI.
Gelhar, L.W. and J.L. Wilson. 1974. Ground Water Quality Modeling.
Ground Water, 12:399-408.
Haith, D.A., L.J. Tubbs and N.B. Picketing. 1984. Simulation of Pollution
by Soil Erosion and Soil Nutrient Loss. PUDOC, Wagntngen, The
Netherlands.
Hughes, H.B.F. and K.S. Porter. 1984. Tracking Aldtcarb Levels in Long
Island Ground Water. Unpublished. Water Resources Program. Center for
Environmental Research. Cornell University, Ithaca, NY.
Intera. 1980. Mathematical Simulation of Aldicarb Behavior on Long
Island. Unsaturated Flow and Ground Water Transport. Inters Environ-
mental Consultants, Inc., Houston, TX.
Kain, D.P. and T.S. Steenhuis. 198&. Adsorption Partition Coefficients
and Degradation Rate Constants and Half Life of Selected Pesticides
compiled from literature. Department of .Agricultural Engineering
Professional Paper 84-3. Cornell University, Ithaca, NY.
Knisel, W.G. 1980. CREAMS: A Field-Scale Model for Chemicals Runoff and
Erosion from Agricultural Management Systems. USDA Cons. Res. Rep. 26.
Merva, G. and A. Fernandez. 1982. Simplified Application of Penman's
Equation. ASAE Paper No. 82-2013, American Society of Agricultural
Engineers, St. Joseph, MI.
Norman, W.R. 1984. Drought-flow Analysis and Prediction in Small Water-
sheds. M.S. Thesis, Cornell University, Department of Agricultural
Engineering, Ithaca, NY.
Pickering, N.B. 1982. Operational Stochastic Meteorologic Models for Non-
point Source Pollution Modeling. M.S. Thesis, Cornell University,
Department of Agricultural Engineering, Ithaca, NY.
Rao, P.S.C. and R.J. Wagenet. 1984. Spatial Variability of Pesticides in
Field Soils: Methods for Data Analysis and Consequences. In review.
Weed Sci. Soc. Amer. J.
Rao, P.S.C., J.M. Davidson and R.E. Jessup. 1981. Simulation of Nitrogen
Behavior in Cropped Land Areas Receiving Organic Wastes:33-37. IN:
M.J. Frissel and J.A. vanVeen (eds.). Nitrogen Behavior of Soil-Plant
Systems. PUDOC, Wageningen, The Netherlands.
Rao, P.S.C., J.M. Davidson and L.C. Hammond. 1976. Estimation of Non-
reactive and Reactive Solute Front Locations in Soils. In: Fuller, W.H.
(ed): Residual Management by Land Dispos~l. Proc. of Hazardous Wastes
Research Symp., Tucson, AZ. EPA 600/9-76-015.
Scheidegger, A.E. 1960. The Physics of Flow Through Porous Media.
MacMillan Co., NY. pp. 256-274.
Skaggs, R.W. 1979. A Water Management Model for Shallow Water Table
Soils. Report No. 134, Water Resources Research Inst., North Carolina
Stats University.
Spillner, C.J., V.M. Thomas, D.C. Takahashi and H.B. Scher. 1983. A
Comparative Study of the Relationships Between the Mobility of Aiachlor,
Butylate and Metolachlor in Soil and Their Physico Chemical Properties.
IN: Fate of Chemicals in the Environment, American Chemical Society.
Steenhuts, T.S., C. Jackson, S. Kung and W.G. Brutsaert. 1984. Measure-
ment of Ground Water Recharge on Eastern Long Island. Journal of
Hydrology (in press).
Steenhuis, T.S., R. E. Muck and M.F. Walter. 1983. Predictions of Water
Budgets for Soils with or Without a Hardpan. IN: Advances in Infiltra-
tion. ASAE Monograph.
Steenhuis, T.S. 1979. Simulation of Soil and Water Conservation
Practices. Action in Controlling Pesticides. IN: D.A. Haith and R.C.
Loehr (eds.), Effectiveness of Soil and Water Conservation Practices for
Pollution Control. EPA 600/3-79-106. U.S. Environmental Protection
Agency.
Trsutman, N.M., K.S. Porter and H.B.F. Hughes. 1983. Southhold Demonstra-
tion Site: New York State Fertilizer and Pesticide Demonstration
Project. Unpublished. Center for Environmental Research, Cornell
University, Ithaca, NY.
U.S. Army Corps of Engineers. 1960. Runoff from Snowmalt. Manual
1110-2-1406. Washington, DC.
3O
U.S. Environmental Protection Agency. 1984, A Ground Water Protection
Strategy for the Environmental Protection Agency, Draft. Washington,
DC.
U,S, Water Resources Council, 1978, The Nation's Water Resources,
1975-2000. Vol. 1: Summary. Second National Water Assessment.
Wagenet, R.J. and P.S.C. Rao. 19874. Basic Concepts of Modeling Pesticide
Fate in the Crop Root Zone (in review), Weed Sci. Soc. Am. J.
Walter, H,F. 1974. Nitrate Hove~nt in Soil Under Early Spring Condi-
tions. Unpublished Ph.D, Thesis, University of WEsconein, Madison, WI.
3!