Groundwater:Boundary Conditions

Each cell in the active domain is assigned a boundary condition with a map of integer values specified with the GW_BOUNDFILE project card. Boundary conditions around the border of the watershed may be constant head (type 2), constant flux (type 0) or a combination of both. Flux boundaries are represented by setting the transmissivity of flow boundary cells to zero and entering the flux in the source/sink term, W. For all non-head boundary cells (type 1), the unsaturated zone provides a flux to the source/sink term, W. The lower boundary is a flux term (cm/hr) that represents a leaky aquifer, specified with the GW_LEAKAGE_RATE project card. The default is zero. This single user-specified value is used in all cells.

Possible internal boundary conditions include constant head (type 2), constant flux (type 6), variable flux (type 3), river head (type 5) and river flux (type 4), and reservoir cell (type 7). The constant flux and variable flux boundaries are used to represent pumping wells of constant and variable rate, respectively. Water added/subtracted by pumping is added to the source/sink term, W. The possible boundary types are:

0 – no flow

1 – regular infiltration cell (no special boundary condition)

2 – specified head, taken from WATER_TABLE file

3 – dynamic flux (well)

4 – stream cell with calculated flux between stream node and groundwater cell

5 – stream cell with a specified head, value from stream routing solution

6 – static flux (well) value taken from the GW_FLUX_BOUNDTABLE

7 - reservoir cell

8.5.1 Stream Boundary Cells – Type 4 & 5

When 1-D channel flow is simulated, grid cells containing stream channels may be represented as either head (type 5) or flux boundaries (type 4). For a head boundary, the elevation of the water surface in the stream node is used as a specified head in the groundwater solution. For a river flux boundary condition the flux between the groundwater cell and stream node is calculated during the channel routing update. Fluxes accumulated over multiple channel routing updates are accumulated and then added to the groundwater cell in the source/sink term, W. Additional inputs are required to simulate groundwater/stream interactions:

thickness of the riverbed material, M_rb, (cm), and
vertical hydraulic conductivity of riverbed material, K_rb (cm/hr).

These may be specified as uniform values using the M_RIVER and K_RIVER project cards, respectively. These may also be distributed along the stream network in the channel input file (cip file), specified with the CHAN_INPUT card, as described in Section 5.1 Channel Routing.

The flux throught the river bed is calculated based on the difference in elevation between the groundwater surface and the water surface elevation in the stream node. McDonald and Harbaugh (1988) developed a simple method to compute flows between the channel network and saturated groundwater based on the Darcy equation. If the groundwater surface elevation is above the riverbed elevation but below the river water surface elevation, the river discharges to the groundwater:

(47)

where:

f = per unit area discharge (m/s),

K_rb = hydraulic conductivity of the river bed material (cm/hr),

M_rb = depth of the river bed material (cm),

E_r = elevation of the river water surface (m), and

E_ws = elevation of the groundwater water surface (m).

The negative flux indicates that the discharge is into the groundwater. If the groundwater surface elevation is below the riverbed elevation the river leaks to the groundwater at the rate

ƒ = -K_rb

(48)

When the STREAM_LOSS card is used in the absence of a known or calculated water table location, Equation 48 is used to compute the stream loss.

If the groundwater surface elevation is above the river water surface elevation, the groundwater discharges to the river according to Equation 46. In this case, the flux will be positive, indicating flow is from the groundwater to the river. The unit flux, ƒ, is multiplied by the top width and length of the channel segment. In calculating the top width an effective depth is defined. If the water surface elevation in the channel is higher than the water surface elevation of the groundwater, the channel depth is used in the calculation. When the groundwater level is higher than the water level in the stream, the effective depth used is the groundwater elevation minus the bed elevation.

The channel side slope properties are assumed to be the same as those of the channel bed. If properties vary significantly between the bed and the side slope an average or weighted average values should be input. This approach is highly empirical and the riverbed property values will require calibration to the system being modeled.

8.5.2 Static Pumping Wells – Type 6

Static pumping rate wells may be placed in any cell in the watershed domain. Pumping wells are located by placing a value of 6 in the GW_FLUX_BOUNDTABLE at each desired well location. The GW_FLUX_BOUNDTABLE specifies the name of a file with the pumping rates. The GW_FLUX_BOUNDTABLE has the following format.

I location well 1          J location well 1          Pumping rate (m³/d) well 1
I location well 2          J location well 2          Pumping rate (m³/d) well 2
I location well 3          J location well 3          Pumping rate (m³/d) well 3
…
…
…
I location well N-1    J location well N-1     Pumping rate (m³/d) well N-1 I location well N        J location well N         Pumping rate (m³/d) well N

Values are separated by spaces or tabs, not commas. Groundwater extractions have a positive pumping rate; injections have a negative pumping rate. The number and I, J locations of the wells in table must match the number and location of wells in the GW_BOUNFILE file.

8.5.3 Reservoir Cells - Type 7

Any cell that is within an active reservoir is assigned a boundary condition of 7. As the reservoirs change size within the GSSHA model, type 7 boundary conditions are assigned internally to GSSHA. If a cell has a type 7 boundary a flux between the reservoir and the groundwater will be calculated as the Darcy flux based on the head between the reservoir and the groundwater, the depth of sediments in the bottom of the reservoir, and the hydraulic conductivity of the sediments. The equations are the same as those described for groundwater channel interaction. Values of hydraulic conductivity for each cell are taken as the hydraulic conductivities prescribed for the infiltration model. The thickness of the bed sediments is specified with the "M_LAKE" card in the CHAN_INPUT file. A uniform value, for all lakes, can be specified in the project file with the M_LAKE card. The value is specified in cm. There is no default value. A different value can be entered for any or all of the lakes. Groundwater interaction will not be calculated for reservoir links without the "M_LAKE" card in the CHAN_INPUT file or the M_LAKE card in the project file.

8.5.4 Static and Dynamic Wells using the Mapping Table

In addition to specifying static wells above, both static and dynamic wells can be created by using the WELL_TABLE table of the mapping table file. An index map of wells with unique pumping characteristics is created (0's elsewhere.) This index map and the associated table is then used to modify the groundwater boundary map. This table has two parameters. The first is a flag (1 for yes and 0 for no) to indicate whether the well is dynamic (1) or static (0). If the well is static, the second parameter is the pumping rate (positive values) or injection rate (negative values). If the well is dynamic then the second value refers to a time series index; you will need to also set up the time series index table.

To create a dynamic well, there are four steps that need to be followed.
1) set up a time series of the pumping rates (positive values) or injection rates (negative values) (see Time series file format)
2) reference the time series in the project file (TIME_SERIES_FILE "filename.ts", as indicated in the above link to time series) preferably before the mapping table.
3) set up the time series index mapping table
4) set up the WELL_TABLE table, with the is_dynamic value set to 1 and the value set to the time series index for that particular time series.