For Each Simulation:

Data Set 1a

An optional character variable used for activating the option to specify text files containing relations among lake stage, surface area, and volume. A separate text file is specified for each lake; each file contains exactly 151 lines; and each line consists of one value for lake stage, volume, and surface area, respectively. This file must contain exactly 151 lines of data because this number of values is consistent with the number of values that are calculated internally by the Lake Package if the “TABLEINPUT” option is not used, and ensures that results will be consistent with results for simulations that calculate lake bathymetry information internally. For example, Lake Package test problems can be run without using the “TABLEINPUT” option and the Lake Package will output the lake bathymetry tables as 151 values to the List file. Bathymetry values printed to the List file can then be used directly for specifying a separate bathymetry input file using the keyword option “TABLEINPUT” and the solution will be the same. However, for difficult problems convergence can be improved by smoothing the lake bathymetry data (for example, by fitting polynomial regression equations to the original bathymetry information). If the keyword “TABLEINPUT” is entered on the first line (record) of the data set, then lake bathymetry data will be read from separate input files. These tables are used for calculating water balances within lakes; however, calculation of seepage between lakes and groundwater is unchanged from the original Lake Package on the basis of the lake and groundwater discretization (Merritt and Konikow, 2000). See IUNITTAB


Data Set 1b

NLAKES ILKCB
If data are being read using the fixed format mode, then each field should be entered using I10 format.
NLAKES Number of separate lakes.
Sublakes of multiplelake systems are considered separate lakes for input purposes. The variable NLAKES is used, with certain internal assumptions and approximations, to dimension arrays for the simulation.

ILKCB Whether or not to write cellbycell flows (yes if ILKCB> 0, no otherwise). If ILKCB< 0 and "Save Budget" is specified in the Output Control or ICBCFL is not equal to 0, the cellbycell flows will be printed in the standard output file.
ICBCFL is specified in the input to the Output Control Option of MODFLOW.



Data Set 2

THETA [NSSITR SSCNCR] [SURFDEPTH]
If data are being read using the fixed format mode, then the data should be entered using format (F10.4,I10,F10.4).
THETA Explicit (THETA = 0.0), semiimplicit (0.0 < THETA < 1.0), or implicit (THETA = 1.0) solution for lake stages.
SURFDEPTH is read only if THETA is assigned a negative value (the negative value of THETA is then changed to a positive value internally by the code).
•A new method of solving for lake stage uses only the timeweighting factor THETA (Merritt and Konikow, 2000, p. 52) for transient simulations. THETA is automatically set to a value of 1.0 for all steadystate stress periods. For transient stress periods, Explicit (THETA = 0.0), semiimplicit (0.0 < THETA < 1.0), or implicit (THETA = 1.0) solutions can be used to calculate lake stages. The option to specify negative values for THETA is supported to allow specification of additional variables (NSSITER, SSCNCR, SURFDEP) for simulations that only include transient stress periods. If THETA is specified as a negative value, then it is converted to a positive value for calculations of lake stage.
•In MODFLOW2000 and later, ISS is not part of the input. Instead NSSITR or SSCNCR should be included if one or more stress periods is a steady state stress period as defined in Ss/tr in the Discretization file.
•SSCNCR and NSSITR can be read for a transient only simulation by placing a negative sign immeditately in front of THETA. A negative THETA sets a flag which assumes input values for NSSITR and SSCNCR will follow THETA in the format as described by Merritt and Konikow (p. 52). A negative THETA is automatically reset to a positive value after values of NSSITR and SSCNCR are read.
•If the absolute value of THETA is less than 0.5, THETA will be set to 0 after data set 2 is read. 
Maximum number of iterations for Newton’s method of solution for equilibrium lake stages in each MODFLOW iteration for steadystate aquifer head solution. Only read if ISS (option flag input to DIS Package of MODFLOW indicating steadystate solution) is not zero or if THETA is specified as a negative value.
•NSSITR and SSCNCR may be omitted for transient solutions (ISS = 0).
•In MODFLOW2000 and later, ISS is not part of the input. Instead NSSITR or SSCNCR should be included if one or more stress periods is a steady state stress period as defined in Ss/tr in the Discretization file.
•SSCNCR and NSSITR can be read for a transient only simulation by placing a negative sign immeditately in front of THETA. A negative THETA sets a flag which assumes input values for NSSITR and SSCNCR will follow THETA in the format as described by Merritt and Konikow (p. 52). A negative THETA is automatically reset to a positive value after values of NSSITR and SSCNCR are read.
•If NSSITR = 0, a value of 100 will be used instead. 
Convergence criterion for equilibrium lake stage solution by Newton’s method. Only read if ISS is not zero or if THETA is specified as a negative value.
•NSSITR and SSCNCR may be omitted for transient solutions (ISS = 0).
•In MODFLOW2000 and later, ISS is not part of the input. Instead NSSITR or SSCNCR should be included if one or more stress periods is a steady state stress period as defined in Ss/tr in the Discretization file.
•SSCNCR and NSSITR can be read for a transient only simulation by placing a negative sign immeditately in front of THETA. A negative THETA sets a flag which assumes input values for NSSITR and SSCNCR will follow THETA in the format as described by Merritt and Konikow (p. 52). A negative THETA is automatically reset to a positive value after values of NSSITR and SSCNCR are read.
•If SSCNCR = 0, a value of 0.0001 will be used instead. 
The height of small topological variations (undulations) in lakebottom elevations that can affect groundwater discharge to lakes.
SURFDEPTH decreases the lakebed conductance for vertical flow across a horizontal lakebed caused both by a groundwater head that is between the lakebed and the lakebed plus SURFDEPTH and a lake stage that is also between the lakebed and the lakebed plus SURFDEPTH. This method provides a smooth transition from a condition of no groundwater discharge to a lake, when groundwater head is below the lakebed, to a condition of increasing groundwater discharge to a lake as groundwater head becomes greater than the elevation of the dry lakebed. The method also allows for the transition of seepage from a lake to groundwater when the lake stage decreases to the lakebed elevation. Values of SURFDEPTH ranging from 0.01 to 0.5 have been used successfully in test simulations. SURFDEP is read only if THETA is specified as a negative value.

Notes:
1. If data are being read using the fixed format mode, then the data should be entered using F10.4 format for floating point numbers and I10 format for integer numbers.
2. If the model simulation includes a steadystate stress period, then the number of iterations (NSSITER) and the closure tolerance (SSCNCR) defined by Merritt and Konikow (2000, p. 52) for ending the Newton iteration method is used for all subsequent transient or steadystate stress periods. If the model simulation only includes transient stress periods, a default value of 0.0001 is assigned to SSCNCR and a default value of 100 is assigned to NSSITR. An option was added in which values of SSCNCR and NSSITR can be read for a transientonly simulation by specifying a negative value for THETA. If THETA is specified as a negative value then it is converted to a positive value for calculations of lake stage.


For the First Stress Period Only:

Data Set 3




STAGES [SSMN SSMX] [IUNITTAB] [CLAKE(1)..........CLAKE(NSOL)]

This data set should consist of one line for each lake, where line 1 includes data for lake 1, and line n includes data for lake n. There must be exactly NLAKES lines of data.
If data are being read using the fixed format mode, then each field should be entered using F10.4 format.
STAGES The initial stage of each lake at the beginning of the run.

SSMN Minimum stage allowed for each lake in steadystate solution.
•SSMN and SSMX are not needed for a transient run and must be omitted when the solution is transient.
•When the first stress period is a steadystate stress period, SSMN is defined in record 3. For subsequent steadystate stress periods, SSMN is defined in record 9a. 
SSMX Maximum stage allowed for each lake in steadystate solution.
•SSMN and SSMX are not needed for a transient run and must be omitted when the solution is transient.
•When the first stress period is a steadystate stress period, SSMN is defined in record 3. For subsequent steadystate stress periods, SSMN is defined in record 9a. 
The unit number of the text files containing data describing relations among lake stage, surface area, and volume. Only read if the keyword “TABLEINPUT” is entered on the first line (record) of the input file.

CLAKE The initial concentrations in each lake at the beginning of the model run. Values are entered for NSOL constituents. The value of NSOL is passed from MOC3D or GWT. CLAKE values are ignored if entered in MODFLOW runs.

Notes:
1. This data set should consist of one line for each lake, where line 1 includes data for lake 1, and line n includes data for lake n. There must be exactly NLAKES lines of data.
2. SSMN and SSMX must be omitted when the solution is transient. For simulations that do not include a steadystate stress period, the minimum and maximum values of lake stage for calculated based on the lowest elevation of the lakebed and the and highest celltop elevation for all lake cells for a particular lake, respectively.
3. If data are being read using the fixed format mode, then each field should be entered using F10.4 format.
4. IUNITAB is used to indicate that the table of values of lake stage, surface area, and volume will be read from a separate input file. These files and associated file unit numbers must be included within the MODFLOW Name file as type ‘Data.’


For Each Stress Period:

Data Set 4

ITMP ITMP1 LWRT
If data are being read using the fixed format mode, then each field should be entered using I10 format.
Lake definition data are restricted to cells for which IBOUND and WETDRY values have been set to zero.
ITMP > 0, read lake definition data (records 57, and, optionally, records 8 and 9);
ITMP = 0, no lake calculations this stress period;
ITMP < 0, use lake definition data from last stress period.

ITMP1 ≥ 0, read new recharge, evaporation, runoff, and withdrawal data for each lake, and associated concentrations if needed for MOC3D runs;
ITMP1 < 0, use recharge, evaporation, runoff, and withdrawal data, and concentrations, if needed, from last stress period.

LWRT > 0, suppresses printout from the lake package.
ICBCFL ≤ 0 or not specifying "Save Budget" also suppresses printout from the lake package. ICBCFL or "Save Budget" are specified in the input to the Output Control option of MODFLOW.



If ITMP > 0 read lake definition data (records 57, and, optionally, records 8 and 9).

Data Set 5

LKARR(NCOL,NROW)
A NCOL by NROW array is read for each layer in the grid by MODFLOW module U2DINT.
LKARR A value is read in for every grid cell.
If LKARR(I,J,K) = 0, the grid cell is not a lake volume cell.
If LKARR(I,J,K) > 0, its value is the identification number of the lake occupying the grid cell. LKARR(I,J,K) must not exceed the value NLAKES. If it does, or if LKARR(I,J,K) < 0, LKARR(I,J,K) is set to zero.
Lake cells cannot be overlain by nonlake cells in a higher layer.
Lake cells must be inactive cells (IBOUND = 0) and should not be convertible to active cells (WETDRY = 0).
The Lake package can be used when all or some of the model layers containing the lake are confined. The authors recommend using the LayerProperty Flow Package (LPF) for this case, although the BCF and HUF Packages will work too. However, when using the BCF6 package to define aquifer properties, lake/aquifer conductances in the lateral direction are based solely on the lakebed leakance (and not on the lateral transmissivity of the aquifer layer). As before, when the BCF6 package is used, vertical lake/aquifer conductances are based on lakebed conductance and on the vertical hydraulic conductivity of the aquifer layer underlying the lake when the wet/dry option is implemented, and only on the lakebed leakance when the wet/dry option is not implemented.


Data Set 6

BDLKNC(NCOL,NROW)
A NCOL by NROW array is read for each layer in the grid by MODFLOW module U2DREL.
BDLKNC A value is read in for every grid cell. The value is the lakebed leakance that will be assigned to lake/aquifer interfaces that occur in the corresponding grid cell.
If the wetdry option flag (IWDFLG) is not active (cells cannot rewet if they become dry), then the BDLKNC values are assumed to represent the combined leakances of the lakebed material and the aquifer material between the lake and the centers of the underlying grid cells, i. e., the vertical conductance values (CV) will not be used in the computation of conductances across lake/aquifer boundary faces in the vertical direction.
IBOUND and WETDRY should be set to zero for every cell for which LKARR is not equal to zero. IBOUND is defined in the input to the Basic Package of MODFLOW). WETDRY is defined in the input to the BCF or other flow package of MODFLOW if the IWDFLG option is active.
When used with the HUF package, the Lake Package has been modified to compute effective lakeaquifer conductance solely on the basis of the userspecified value of lakebed leakance; aquifer hydraulic conductivities are not used in this calculation. An appropriate informational message is now printed after the lakebed conductances are written to the main output file.


Data Set 7

NSLMS
If data are being read using the fixed format mode, then NSLMS should be entered using format I5.
NSLMS The number of sublake systems if coalescing/dividing lakes are to be simulated (only in transient runs). Enter 0 if no sublake systems are to be simulated.


If ITMP > 0 and NSLMS > 0:
A pair of records (records 8a and 8b) is read for each multiplelake system, i.e., NSLMS pairs of records. However, IC = 0 will terminate the input.

Data Set 8a

IC ISUB(1) ISUB(2) ............ ISUB(IC)
If data are being read using the fixed format mode, then each field of Record 8a should be entered using I5 format
IC The number of sublakes, including the center lake, in the sublake system being described in this record.

ISUB The identification numbers of the sublakes in the sublake system being described in this record. The center lake number is listed first.



Data Set 8b

SILLVT(2) ............. SILLVT(IC)
If data are being read using the fixed format mode, then each field of Record 8b should be entered using F10.4 format
SILLVT Sill elevation that determines whether the center lake is connected with a given sublake. One value is entered in this record for each sublake in the order the sublakes are listed in the previous record.


If ITMP1 ≥ 0
At least one of the above records will be read for each lake; i.e., NLAKES records, or sets of records, will be read. If MODFLOW is being run, only the first record is read. If MOC3D or GWT is being run, a set of two or more records will be read for each lake (see next note).
If record 9b is included because solute transport is being simulated, then 9b should consist of one record (line) for each solute; each record must contain two or three values; and there must be as many records as the number of solutes being simulated (NSOL). The order of records must be that all necessary lines for 9b are listed for a given lake before line 9a for the next lake. For example, if the Lake Package is representing three lakes and the solute transport package is representing two solutes, then the order of data for Record 9 would be 9a, 9b, 9b, 9a, 9b, 9b, 9a, 9b, and 9b.

Data Set 9a

PRCPLK EVAPLK RNF WTHDRW [SSMN] [SSMX]
PRCPLK The rate of precipitation per unit area at the surface of a lake (L/T).

EVAPLK The rate of evaporation per unit area from the surface of a lake (L/T).

RNF Overland runoff from an adjacent watershed entering the lake. If RNF > 0, it is specified directly as a volumetric rate, or flux (L3 /T). If RNF < 0, its absolute value is used as a dimensionless multiplier applied to the product of the lake precipitation rate per unit area (PRCPLK) and the surface area of the lake at its full stage (occupying all layer 1 lake cells).
When RNF is entered as a dimensionless multiplier (RNF < 0), it is considered to be the product of two proportionality factors. The first is the ratio of the area of the basin contributing runoff to the surface area of the lake when it is at full stage. The second is the fraction of the current rainfall rate that becomes runoff to the lake. This procedure provides a means for the automated computation of runoff rate from a watershed to a lake as a function of varying rainfall rate. For example, if the basin area is 10 times greater than the surface area of the lake, and 20 percent of the precipitation on the basin becomes overland runoff directly into the lake, then set RNF = 2.0.

WTHDRW The volumetric rate, or flux (L3 /T), of water removal from a lake by means other than rainfall, evaporation, surface outflow, or groundwater seepage. A negative value indicates augmentation. Normally, this would be used to specify the rate of artificial withdrawal from a lake for human water use, or if negative, artificial augmentation of a lake volume for esthetic or recreational purposes.

SSMN Minimum stage allowed for each lake in steadystate solution.
•SSMN and SSMX are not needed for a transient run and must be omitted when the solution is transient.
•When the first stress period is a steadystate stress period, SSMN is defined in record 3. For subsequent steadystate stress periods, SSMN is defined in record 9a. 
SSMX Maximum stage allowed for each lake in steadystate solution.
•SSMN and SSMX are not needed for a transient run and must be omitted when the solution is transient.
•When the first stress period is a steadystate stress period, SSMN is defined in record 3. For subsequent steadystate stress periods, SSMN is defined in record 9a. 


Data Set 9b

If solute transport is also being simulated (Ftype “CONC” exists), then for each solute the following data are read:
CPPT(NSOL) CRNF(NSOL) [CAUG(NSOL)]
If data are being read using the fixed format mode, then each field should be entered using F10.4 format.
CPPT The concentration of solute in precipitation onto the lake surface.

CRNF The concentration of solute in overland runoff directly into the lake.

CAUG The concentration of solute in water used to augment the lake volume.
It is implicitly assumed that no solute is present in water evaporated from the lake surface.


