380 likes | 448 Views
GWF BCF6 PACKAGE. Computes the conductance components of the finite-difference equation which determine flow between adjacent cells. Computes the terms that determine the rate of water movement to and from storage.
E N D
GWF BCF6 PACKAGE • Computes the conductance components of the finite-difference equation which determine flow between adjacent cells. • Computes the terms that determine the rate of water movement to and from storage. • Calculates flow correction terms that are added to the difference equations when an underlying aquifer becomes partially saturated. • Requires the node be located in the center of the cell, hence the name Block-Centered Flow. BCF6
GWF BCF6 PACKAGE Basic Conduction Equations Review • Darcy’s law define one-dimensional flow • Conductance is defined as, • Darcy’s law can be written, BCF6
GWF BCF6 PACKAGE • For a set of conductances arranged in series, the inverse of the equivalent conductance equals the sum of the inverses of the individual conductances • When there are only 2 sections, the equivalent conductance reduces to, Note: the above is call the harmonic mean of the conductances C1 and C2 BCF6
GWF BCF6 PACKAGE Applying C1C2/(C1+C2) Horizontal Conductance under confined conditions • Conductances are defined between nodes of adjacent cells rather than within a cell • CR (conductance along rows) and CC (conductance along columns) are calculated between adjacent horizontal nodes. • The subscript ½ is used designate conductance between nodes (e.g. CRij+½k represents the conductance between nodes i,j,k and i,j+1,k) • Transmissivity is uniform over a cell, but may vary from cell to cell. • If the transmissivity of both cell is zero, the conductance between the nodes of the two cells is set to zero BCF6
Horizontal Conductance Under Water Table Conditions If a model layer is confined, the horizontal conductance will be constant If a layer is unconfined or potentially unconfinined, new values of the horizontal conductance must be calculated as the head fluctuates. The transmissivity is calculated (product of hydraulic conductivity and saturated thickness). When the head drops below the aquifer bottom, the cell is considered dewatered. Transmissivity within a cell in the row direction is calculated, if HNEWijk ≥ TOPijk then TRijk= (TOPijk-BOTijk)HYRijk if TOPijk > HNEWijk > BOTijk then TRijk= (HNEWijk-BOTijk)HYRijk if HNEWijk ≤ BOTijk then TRijk= 0 HYRijk is the row directed hydraulic conductivity TOPijk is the top elevation of cell i,j,k BOTijk is the bottom elevation of cell i,j,k GWF BCF6 PACKAGE BCF6
GWF BCF6 PACKAGE Vertical Conductance Formulation • The vertical interval between two nodes, i,j,k and i,j,k+1, may be considered to contain n geohydrologic layers. • Let the layers have hydraulic conductivity K,1, K2 … Kn and thickness Δz1, Δz2, … Δzn. • The map area of the cells around the nodes i,j,k and i,j,k+1 is DELRj×DELCi • The vertical conductance of an individual geohydrologic layer g is, BCF6
GWF BCF6 PACKAGE • The equivalent vertical conductance, Cij+½k, for the vertical interval between i,j,k and i,j,k+1, • Rearranging above terms gives, • The above quantity is called the vertical leakance and is designated Vcontij+½k BCF6
GWF BCF6 PACKAGE • Vcont is the input to the BCF file and is calculated externally to the model. • Vcont values are read as two dimensional arrays for each layer. • Each value of Vcontij is the vertical leakance for the interval between cell i,j,k and cell i,j,k+1 (for the interval between the layer for which the array is read, and the layer below it). • The Vcont array is not read for the lowermost layer. BCF6
GWF BCF6 PACKAGE Case 1-Single Geohydrologic Units If i,j,k and i,j,k+1 both fall with in a single hydrogeologic unit, having a vertical hydraulic conductivity Kzij that is uniform within the two cell region, then where Δzk+½ is the vertical distance between nodes and is the sum of Δvk/2 and Δvk+1/2 and Δv represents the layer thickness. BCF6
GWF BCF6 PACKAGE Case 2-Different Geohydrologic Layers If i,j,k and i,j,k+1 fall at the midpoint of two adjacent hydrogeologic units, each having its own vertical hydraulic conductivity that is uniform within its cell region, then where Δvk is the thickness of layer k, Δvk11 is the thickness of layer k+1, Kijk is the hydraulic conductivity of cell i,j,k Kijk+1 is the hydraulic conductivity of cell i,j,k+1 BCF6
GWF BCF6 PACKAGE Case 3-Quazi 3D Two aquifers separated by a confining unit, Where Δzu, Δzc, ΔzL are the thickness of the upper, confined, and lower units, respectively, and Kzu, Kzc, KzL are the their hydraulic conductivities. If Kzc is much smaller than Kzu and KzL then the above equation reduces to This formulation is called the quasi-3D approach. BCF6
GWF BCF6 PACKAGE Vertical Flow Calculations Under Dewatered Conditions • The basic flow equation for cell i,j,k is, This term gives the flow into cell i,j,k through its lower face BCF6
Cell ijk+1 Dewatered Assume cell i,j,k and the confining unit are fully saturated, then head at upper surface of confining layer is hijk. Below the confining zone is unsaturated, and the head at the base of the confining unit is atmospheric, and the head there is hijk+1 = Topijk+1 and the flow through the confining bed is given by, instead of the term from equation 1, To avoid asymmetry, we leave equation 2 in equation 1 and add, to the RHS of the flow equation. (notice the iteration parameter is set at the previous iteration, n-1) GWF BCF6 PACKAGE BCF6
Cell ijk Dewatered A correction must also be applied for the dewatered cell itself. Let the dewatered cell be i,j,k an consider flow into i,j,k from overlying cell i,j,k-1. The computed flow into i,j,k from the cell above is, whereas the actual flow into the cell is, This time we add, to the RHS of the flow equation. (notice the iteration parameter is set to present iteration because the correction term does not affect the symmetry—it is added only to a diagonal element of equation 1) GWF BCF6 PACKAGE BCF6
Storage Formulation There are two types of layers that are considered, Layer whose storage values remain constant Layer whose storage properties may convert from confined to unconfined or vice-versa. If a layer’s storage values remain constant, the rate of accumulation of water in the cell, ΔV/Δt, is given by, where SSijk(ΔrjΔciΔvk) confined SC1ijk = SY(ΔrjΔci ) unconfined The SC1ijkis called the primary storage capacity for cell i,j,k. Note: The primary storage capacity is adequate if the water level in the cell remains either above the top of the cell or below the top of the cell through out the simulation. GWF BCF6 PACKAGE BCF6
GWF BCF6 PACKAGE Storage Term Conversion • During any time step, there are four possible storage conditions for each cell • The cell is confined for the entire time step • The cell is unconfined for the entire time step • The cell converts from confined to unconfined • The cell converts from unconfined to confined • The following expression for the rate of accumulation in storage in cell i,j,k is used, where Top is the elevation of the top of the model cell, SCA is the storage capacity in effect in the cell at the start of the time step, and SBC is the “current” storage capacity (current iteration) BCF6
GWF BCF6 PACKAGE • Storage Values • If and • then, • SCA = SSijk(ΔrjΔciΔvk) and • SCB = SSijk(ΔrjΔciΔvk) • giving, i,j,k BCF6
GWF BCF6 PACKAGE Storage Values 2. If and then, SCA = SY(ΔrjΔci) and SCB = SY(ΔrjΔci) giving, i,j,k BCF6
GWF BCF6 PACKAGE 3. If and then, SCA = SSijk(ΔrjΔciΔvk) and SCB = SY(ΔrjΔci) giving, BCF6
GWF BCF6 PACKAGE 4. If and then, SCA = SY(ΔrjΔci) and SCB = SSijk(ΔrjΔciΔvk) giving, BCF6
GWF BCF6 PACKAGE • As part of the simulation of unconfined aquifers and aquifers that can convert between confined and unconfined, MODFLOW can change a variable-head cell to a no-flow cell. • If the saturated thickness becomes zero, MODFLOW converts the cell to no-flow—called “drying” the cell. • Based on heads in surrounding cells, MODFLOW will attempt to wet cells that are dry. • The user can specify the cells for which wetting is attempted. • Wetting capability useful for, • Recovery of water levels when wells are turned off, • Modeling mounds of recharge water from irrigation application, and • Situations where cells incorrectly go dry (convert to no-flow) as part of the iterative solution process. BCF6
GWF BCF6 PACKAGE Recall that MODFLOW calculates transmissivity as if HNEWijk ≥ TOPijk then TRijk= (TOPijk-BOTijk)HYRijk if TOPijk > HNEWijk > BOTijk then TRijk= (HNEWijk-BOTijk)HYRijk if HNEWijk ≤ BOTijk then TRijk= 0 • The transmissivity values are then used to calculate row and column conductances. • Vertical conductance is constant till the cell becomes dry, at which point is changed to zero. • When a cell becomes dry, IBOUND is set to zero, all conductance to the cell are set to zero, and head is set to a very large value to serve as a visual indicator. BCF6
GWF BCF6 PACKAGE • A dry cell is allowed to become wet if the head from the previous iteration in a neighboring cell is greater than or equal to a turn-on threshold, TURNON = BOT + THRESH where, BOT is the bottom elevation of a dried-out cell, THRESH is a user-specified constant called the wetting threshold BCF6
GWF BCF6 PACKAGE Option 1 • There are two options to select which neighboring are checked to see if the turn-on threshold has been reached, • Check the cell immediately below the dry cell and the four horizontally adjacent cells, or • Check only the cell immediately below the dry cell. • If the neighboring cell is either no-flow or constant-head, then that cell is not checked for TURNON. hn dry cell cell n cell i hk hi BOT + THRESH BOT cell k Option 2 dry cell hk BOT + THRESH BOT Cell k BCF6
GWF BCF6 PACKAGE • When a cell is wetted, IBOUND for the cell is set to 1, the vertical conductance for the cell are set to their origin values, and the head is set to either, h = BOT + WETFCT(hn – BOT) or h = BOT + WETFCT(THRESH) where, BOT is the bottom elevation of the dry cell, hn is the head at the neighboring cell that caused the cell to wet, WETFCT is a user-specified constant called the wetting factor. Note: The head assigned to a cell might exceed the wetting threshold of a neighboring dry cell, however a neighboring cell cannot become wet in the same iteration. BCF6
GWF BCF6 PACKAGE • There is a non-uniqueness associated with the wetting routine Threshold h starting h final h final BOT h starting Remains Dry Remains Wet The method of wetting and drying cells can cause problems with the convergence of iterative solvers used in MODFLOW. BCF6
GWF BCF6 PACKAGE IBCFB—is a flag and a unit number If IBCFCB > 0, it is a unit number to which cell-by-cell terms are written when SAVE BUDGET or a non-zero value for ICBCFL is specified in Output Control. IBCFCB = 0, cell-by-cell flow terms will not be written. IBCFCB < 0, cell-by-cell flows for constant-head cells will be written in the LIST FILE when SAVE BUDGET or a non-zero value for ICBCFL is specified in Output Control. Cell-by-cell flow to storage and between adjacent cells will not be written to any file. BCF6
GWF BCF6 PACKAGE HDRY—isthe head that is assigned to cell that are converted to dry during simulation. HDRY is similar to HNOFLO, it is an indicator. IWDFLG—is a flag that determines if the wetting capability is active. IWDFLG = 0, the wetting capability is inactive. IWDFLG ≠ 0, the wetting capability is active. WETFCT—is a factor that is included in the calculation of head that is initially establish at a cell when it is converted from dry to wet. BCF6
GWF BCF6 PACKAGE IWETIT—is the iteration interval for attempting to wet cells. Wetting is attempted every IWETIT iteration (outer iterations if PCG). If IWETIT is 0, it is set to 1. IHDWET—is a flag that determines which equation is used to define the initial head to cells that become wet: If IHDWET = 0, then h = BOT + WETFCT(hn − BOT) If IHDWET≠ 0, then h = BOT + WETFCT(WETDRY) note: the absolute value of WETDRY is the wetting threshold (THRESH) BCF6
GWF BCF6 PACKAGE Ltype—contains a combined code for each layer that specifies both the layer type (LAYCON) and the method of computing interblock conductance. Values are two digit numbers. The left digit defines method of calculating interblock transmissivity 0 or blank Harmonic mean 1 Arithmetic mean 2 Logarithmic mean 3 Arithmetic mean of saturated Thickness and Logarithmic-mean hydraulic conductivity. BCF6
GWF BCF6 PACKAGE The right digit defines layer type (LAYCON) 0 confined—Transmissivity and storativity of the layer are constant the entire simulation. 1 unconfined—Transmissivity is calculated from hydraulic conductivity and saturated thickness. Storativity is constant. Type code valid only in first layer. 2 confined/unconfined—Transmissivity is constant. Storativity may alternate between confined and unconfined values. 3 confined/unconfined—Transmissivity is calculated from hydraulic conductivity and saturated thickness. Storativity may alternate between confined and unconfined values. BCF6
GWF BCF6 PACKAGE TRPY—is a one-dimensional real variable containing a horizontal anisotropy factor for each layer. It is the ratio of transmissivity or hydraulic conductivity along the column to the comparable values along the row. Set to 1.0 for isotropic values. There is a single value for each layer and is read using U1DREL. BCF6
GWF BCF6 PACKAGE A subset of the following 2-Dimensional are used to describe each layer. • The variable needed for each layer depends on LAYCON,whether the simulation has any transient stress periods, and if the wetting capability is active (IWDFLG≠0). • If a variable is not needed , it must be omitted. • All variables are read for layer 1 first, then layer 2 and so forth. BCF6
GWF BCF6 PACKAGE Sf1—is the primary storage. Read only if there are one or more transient stress periods. For LAYCON=1, Sf1 will always be specific yield. For LAYCON=2 or 3, Sf1 will always be the storage coefficient. For LAYCON=0, Sf1 would normally be the storage coefficient. However, if it is assumed that the drawdowns in an unconfined aquifer always remain small compared to the saturated thickness (transmissivity does not vary with head), Sf1 may be specific yield. BCF6
GWF BCF6 PACKAGE Tran—is the transmissivity along the rows. Tran is multiplied by TRPY to obtain transmissivity along the columns. Read only for layers where LAYCON is 0 or 2 HY—is the hydraulic conductivity along the rows. HY is multiplied by TRPY to obtain transmissivity along the columns. Read only for layers where LAYCON is 1 or 3 BCF6
GWF BCF6 PACKAGE Vcont—is the vertical hydraulic conductivity divided by the thickness from a layer to a layer below. Vcont is not specified for the bottom layer Sf2—is the secondary storativity value. Read only for layers where LAYCON=2 or 3, and only if there are one or more transient stress periods. The secondary storativity value is always specific yield. BCF6
GWF BCF6 PACKAGE WETDRY—is a combination of the wetting threshold and a flag to indicated which neighboring cells can cause a dry cell to become wet. WETDRY<0, only the cell below the dry cell can cause the dry cell to become wet. WETDRY>0, the cell below the dry cell and the four horizontally adjacent cells can cause the dry cell to become wet. WETDRY=0, the cell can not be wetted. The absolute value of WETDRY is the wetting threshold. Read only if LAYCON=1 or 3 and IWDFLG≠0 BCF6
GWF BCF6 PACKAGE BCF Data Set BCF6