Numerical modeling has been the primary tool for analyzing the interaction of groundwater and surface water since the mid-1960s. Most models continue to use a simple Darcy calculation (where discharge = hydraulic conductivity * hydraulic gradient * cross-sectional area) as well as highly idealized stream geometry, to transfer water (seepage) through the stream sediments based on head differences between the surface water and groundwater. In recent years additional surface-water modules have been added to the commonly used modular groundwater flow model MODFLOW [ McDonald and Harbaugh, 1988]. The original release of MODFLOW had a `river module' that had the simple seepage calculation that was just described, with the further assumption that the stream stage remained constant over one model stress period. Prudic [1989] modified the river module so streamflow routing could be accomodated by the MODFLOW code, but its use is limited to steady flow in rectangular, prismatic channels. More recently, Swain and Wexler [1992] developed another module for MODFLOW, termed MODBRANCH. The code links BRANCH, which is a 1-dimensional numerical model designed to simulate unsteady flow in open-channel networks, to MODFLOW. In these recently developed modules, streams are simulated more realistically than in the original river module, but the transfer of water between surface water and groundwater is still based on the Darcy calculation, relying on head differences between surface water and groundwater and hydraulic conductivity of the streambed. Investigators outside the U.S. Geological Survey, such as Schenk et al. [1990] also have developed modifications to MODFLOW for simulating the interaction of groundwater and surface water.
Other studies related to basic development of numerical models for simulating the interaction of groundwater and surface water include the following. Jorgensen et al. [1989] developed a procedure to account for intracell flow in models that consider recharge to a water-table aquifer that includes a stream. Their study indicates that intercepted intracell flow, where the cell includes recharge and a stream sink, can result in significant model error if the intracell flow is not accounted for. To reduce this error, for model cells that contain both recharge and a sink, recharge needs to be reduced by the ratio of the area of influence of the sink within the cell to the area of the cell.
Modica [1993] used numerical modeling to investigate the configuration of the boundary separating the stream subsystem from regional groundwater for the Upper Rancocas watershed, New Jersey. The study also included determining the flow patterns within the stream subsystem. He found that stream subsystem geometry and its flow patterns are controlled largely by the quantity and distribution of stream discharge. The study also indicated that streams are line sinks that induce complex flow patterns whereby flow from far and near sources are drawn into a common discharge area.
Computer codes have also been developed for simulating the interaction of groundwater and surface water for variably saturated groundwater systems. Cooley [1983] developed a finite-element code that was used by Winter [1983] to evaluate the dynamic groundwater flow conditions caused by nonuniform areal distribution of recharge directly adjacent to surface-water bodies. Using the model for 2-dimensional vertical sections, the study indicated that recharge is initially focused directly adjacent to the surface water. These near-shore flow paths cause increased flow where groundwater seeps into the surface water, and may cause flow to reverse where surface water seeps out to groundwater prior to recharge. Tracy and Marino [1987] developed a Galerkin finite-element code for simulating the interaction of surface water with a variably saturated groundwater system. The model solves for the transient position of surface-water stage and the distribution of the phreatic surface (water table) in the porous medium, as well as providing estimates of seepage rates.
Numerical models have been used as routine analytical tools for numerous field studies in virtually every conceivable hydrogeologic setting and for virtually all scales of systems. Recent studies concerned with understanding the interaction in natural systems have been done for a wide variety of hydrogeologic and climatic environments; such as, the Chalk aquifer of the Berkshire Downs, England [ Rushton et al., 1989], Alsace, France [ Ackerer et al., 1990], a small stream on Long Island, New York [ Prince et al., 1989], a large alluvial aquifer in Kansas [ Wolf and Helgesen, 1993], and high altitude basins in Wyoming [ Hasfurther et al., 1991]. Many other recent studies have used numerical models to evaluate the effect of groundwater pumping on streamflow, such as that by Eberts and Blair [1990].
Numerical models are increasingly being used for studies of solute transport between groundwater and surface water. Duffy and Lee [1992] conducted numerical experiments to evaluate the response of baseflow chemistry to nonpoint source contamination in an idealized stream-aquifer setting. Jakeman et al. [1989] used such a model to study salinity of the Murray River, Australia, as related to a contaminated aquifer. Sorek et al. [1992] modeled flow and transport in a shallow aquifer that was affected by surface reservoirs.