An unconfined groundwater model of the Death Valley Regional Flow Systemand a comparison to its con

The MODFLOW version of the United States Geological Survey (USGS) Death Valley Regional Flow System (DVRFS) in California and Nevada is conceptually inaccurate in that it models an unconfined aquifer as a confined system and does not accurately simulate unconfined drawdown in transient pumping simulations. The transfer of geologic and hydrologic information from the confined MODFLOW DVRFS model to an unconfined MODFLOW–SURFACT (SURFACT) version was accomplished by maintaining cell structure between models and computing effective cell properties to translate the HUF2 package used in MODFLOW to the BCF4 package used by SURFACT. The confined version of the DVRFS was compared to the unconfined SURFACT version by examining head contour maps and the ability of the SURFACT model to match the 4900 observations of hydraulic head/drawdown, 49 observations of groundwater discharge, and 15 estimates of groundwater fluxes into/out of the model domain. Resultant weighted root mean squared error (xRMSE) for the unconfined SURFACT model was lower than the USGS confined model. Despite a lower xRMSE, unconfined conditions simulated with SURFACT did produce greater heads in mountainous regions compared to the confined MODFLOW with differences most pronounced in regions where cell thickness is large, horizontal conductivity small and recharge large. Difference in computed heads reflects computation schemes employed by both models to estimate interblock conductance. Specifically, interblock conductance for the unconfined SURFACT model is dependent on the relative saturation of a modeled cell while MODFLOW’s confined system is not. Despite head differences, SURFACT simulates comparable flux estimates to MODFLOW (e.g. observed ET, groundwater spring flow, and groundwater flux across model boundaries), while significantly improving transient well drawdown estimates. SURFACT is also capable of producing more realistic estimates of water availability from proposed groundwater development and resultant potential impacts to the region.

美国地质调查局（USGS）在加利福尼亚和内华达州的死亡谷地区流系统（DVRFS）的MODFLOW版本在概念上存在不准确之处，

它将潜水含水层建模为封闭系统，

并且在瞬时抽水模拟中不能准确模拟潜水的降水。

​​​​​​​

通过保持模型之间的单元结构并计算有效单元特性，

将受限MODFLOW DVRFS模型中的地质和水文信息转移到无围限的MODFLOW-SURFACT（SURFACT）版本是可行的。

此转换涉及将MODFLOW中使用的HUF2包翻译为SURFACT使用的BCF4包。

通过检查头等值线图以及SURFACT模型匹配4900个水头/降水观测、

49个地下水排泄观测和15个模型域内地下水通量估计的能力，

将受限版本的DVRFS与无围限SURFACT版本进行了比较。

无围限SURFACT模型的结果加权均方根误差（xRMSE）低于USGS受限模型。

尽管xRMSE较低，使用SURFACT模拟的无围限条件在山区产生的水头要比受限的MODFLOW大，差异在单元厚度较大、水平导水性较小和补给量较大的地区最为显著。

计算头差异反映了两个模型用于估计单元间导水性的计算方案。

具体而言，无围限SURFACT模型的单元间导水性取决于模型单元的相对饱和度，而MODFLOW的受限系统则不是如此。

尽管存在头差异，SURFACT模拟了与MODFLOW相当的通量估算（例如观察到的蒸发蒸腾、地下水泉流和模型边界处的地下水通量），同时显著改善了瞬时井降估算。

SURFACT还能够更实际地估计提议的地下水开发的水资源可用性及其对该地区的潜在影响。

Introduction 

There have been numerous water rights applications submitted throughout southern Nevada to offset the needs of growth in Las Vegas. Several applications include groundwater withdrawals adjacent to the Nevada Test Site (NTS) where underground nuclear testing was conducted. If large quantities of groundwater are pumped adjacent to the NTS, the groundwater system could change dramatically. Potential impacts from groundwater pumping include decreasing water levels, reduction in groundwater resources on the NTS, reduction in spring flows adjacent to proposed pumping centers, and the alteration of groundwater flowpaths.

In 1998, the US Department of Energy (DOE) Nevada Site Office funded the United States Geological Survey (USGS) to improve upon two previous groundwater flow models of the Death Valley region with the initial intention of understanding groundwater flowpaths and travel times associated with potential movement of radioactive material from the NTS as well as to characterize the groundwater system in the vicinity of Yucca Mountain, and address effects on users down-gradient from the NTS and Yucca Mountain (Belcher, 2004). The first of these earlier models was developed by DOE for the National Nuclear Security Administration/Nevada Site Office (NNSA/NSO) Underground Test Area (UGTA) project (IT Corporation, 1996). The second was developed by the Office of Civilian Radioactive Waste Management’s (OCRWM) Yucca Mountain Project (YMP) and the NNSA/NSO Hydrologic Resources Management Program (HRMP).

The resultant USGS Death Valley Regional Groundwater Flow System Model (DVRFSM) has improved upon the two previous models by using newly acquired data and modeling tools (Belcher, 2004). Studies compiled and used in the DVRFSM construction include reassessing groundwater discharge via evapotranspiration (ET) (Laczniak et al, 1999, 2001; Reiner et al, 2002; DeMeo et al, 2003) and spring flow (refer to Table C-2, Belcher, 2004), cataloging historical groundwater pumping from 1913 through 1998 (Moreo et al, 2003), reinterpreting groundwater recharge as net infiltration (Hevesi et al, 2002, 2003), assessing model boundary inflows and outflows from regional hydraulic gradients, developing a water budget (Belcher, 2004), and finally incorporating hydraulic conductivity relationships as a function of depth (Belcher et al, 2001, 2002).

Carroll et al (2006) used the DVRFSM to test the impacts of proposed groundwater development south of the NTS and found that, despite increased level of geologic detail and improved water budget accounting, the DVRFSM produced drawdowns on the order of thousands of meters after only a few decades of pumping. Such large drawdowns are not realistic and serve to highlight the most significant limitation of the DVRFSM: the model was built using the confined layer assumption to improve numeric stability. The use of a confined, rather than a ‘‘convertible” (i.e., unconfined), layer type within MODFLOW assumes that the saturated thickness remains constant throughout the entire simulation. Therefore, cells are not allowed to dry or become inactive as the water level decreases below the bottom of a cell. The confined approach was adopted by the USGS because the model was computationally unstable when cells were allowed to convert between confined and unconfined. Unfortunately the confined approach, while numerically stable, allows unrealistic estimates of drawdown to occur in transient simulations.

The primary objective of this study was to produce a numerically stable, unconfined version of the Death Valley Regional Flow System (DVRFS) for future use as a tool to more accurately address potential withdrawals in the NTS region and to determine the potential differences in the solution of the groundwater flow equation given confined and unconfined assumptions.

引言

在整个内华达南部，

已经提交了大量的水权申请，

以满足拉斯维加斯地区的增长需求。

其中一些申请包括在内华达试验场（NTS）附近进行地下核试验的地方进行地下水提取。

如果在NTS附近大量抽取地下水，

地下水系统可能发生巨大变化。

地下水抽取可能导致水位下降、

NTS地下水资源减少、

在拟议的抽水中心附近泉水流量减少

以及地下水流路径的改变等潜在影响。

1998年，

美国能源部（DOE）内华达站办公室

资助了美国地质调查局（USGS）

改进死亡谷地区两个先前地下水流模型的工作，

最初的目的是理解潜在核材料从NTS流动的地下水流路径和时间，

以及表征Yucca Mountain附近的地下水系统，

并解决NTS和Yucca Mountain下游用户的影响（Belcher，2004）。

这两个早期模型中的第一个是由DOE为国家核安全局/内华达站办公室（NNSA/NSO）的地下试验区（UGTA）项目开发的（IT Corporation，1996）。

第二个是由民用放射性废物管理办公室的（OCRWM）Yucca Mountain项目（YMP）和NNSA/NSO水资源管理计划（HRMP）开发的。

由此产生的USGS死亡谷地区地下水流系统模型（DVRFSM）

通过使用新获得的数据和建模工具改进了这两个先前的模型（Belcher，2004）。

用于DVRFSM建模的研究包括重新评估通过蒸散蒸发（ET）排放的地下水（Laczniak等，1999，2001；Reiner等，2002；DeMeo等，2003）

和泉水流量（参见表C-2，Belcher，2004），

编制了1913年至1998年的历史地下水抽取情况（Moreo等，2003），

重新解释地下水充沛为净渗透（Hevesi等，2002，2003），

评估区域水力梯度的模型边界流入和流出，

制定了水量预算（Belcher，2004），

最后将水力传导率关系作为深度的函数纳入模型（Belcher等，2001，2002）。

Carroll等人（2006年）使用DVRFSM测试了NTS南部提出的地下水开发对地下水系统的影响，

并发现尽管具有更高级别的地质细节和改进的水量预算核算，

DVRFSM在仅抽水几十年后就产生了数千米的降水。

这样大的降水是不现实的，

并突显了DVRFSM最重要的局限性：

该模型采用了受限层的假设以提高数值稳定性。

在MODFLOW中使用受限层（而不是“可转换”的，即无围限的层类型）的假设，

假设整个模拟中饱和层厚度保持不变。

因此，当水位降到单元底部以下时，单元不允许变干或变为非活动状态。

USGS采用了受限层的方法，

因为当允许单元在受限和无围限之间转换时，

模型在计算上是不稳定的。

不幸的是，受限的方法虽然在数值上稳定，

但在瞬时模拟中允许出现不现实的降水估算。

本研究的主要目标是

为未来将死亡谷地区地区流系统（DVRFS）作为工具更准确地解决NTS地区的可能抽水和确定在受限和无围限假设下地下水流方程解的潜在差异，

制作一个数值稳定的无围限版本。

Numerical models 

Death Valley Regional Flow System model: MODFLOW 

The DVRFS and the associated USGS model domain (DVRFSM), with the NTS superimposed, are shown in Fig. 1. The DVRFS is approximately 100,000 km2 in Nevada and California and is bounded by latitudes 35000 N and 38150 N and by longitudes 115000 W and 118000 W. This system encompasses flow between recharge areas in the mountains of central and southern Nevada and discharge areas of wet playas and springs south and west of the NTS and in Death Valley, California. The flow is strongly influenced by a complex geologic framework and the USGS DVRFSM incorporates the distribution of the flow system’s principal aquifers and confining units as well as the principal geologic structures that might affect subsurface flow (Belcher, 2004).

The DVRFSM was built with MODFLOW-2000 (McDonald and Harbaugh, 1988; Harbaugh et al, 2000), a public domain, finite-difference code developed by the USGS that is capable of simulating groundwater flow in transient, three-dimensional, anisotropic and heterogeneous systems. MODFLOW’s governing three-dimensional flow equation for a confined aquifer combines Darcy’s Law and the principle of conservation of mass via

where x, y, and z are Cartesian coordinates (L), Kxx, Kyy and Kzz are the principal components of saturated hydraulic conductivity along the x, y, and z axes, respectively (L/T), h is hydraulic head (L), G is a volumetric source/sink term (1/T), S is stortivity (dimensionless), b is saturated thickness (L) which is equivalent to cell thickness, and t is time (T). Conversely, Eq. (2) defines groundwater flow given unconfined conditions, where h (L) is the phreatic surface elevation and Sy is specific yield (dimensionless).

\frac{\partial }{\partial x}(Kxxh \frac{\partial h}{\partial x})+\frac{\partial }{\partial y}(Kyy h\frac{\partial h}{\partial y})+\frac{\partial }{\partial z}(Kzz h\frac{\partial h}{\partial z})-G=Sy\frac{\partial h}{\partial t}\qquad(2)

Solving for h in Eq. (2) is non-linear problem because h depends on hydraulic conductivity and saturated thickness, with saturated thickness dependent on the water table’s elevation. To address non-linearity inherent within unconfined aquifer systems, MODFLOW’s updated HUF2 package maintains a confined system, with the top elevation equal to the approximate elevation of the water table, and substitutes Sy for S into Eq. (1). The HUF2 approach to simplifying the non-linearity of an unconfined solution is only applicable when there is little change in saturated thickness.

Belcher (2004) describes the DVRFSM in detail, however a consolidated review of the DVRFSM is provided here for clarity. The DVRFSM uses a constant-grid spacing in the horizontal plane of 1500 m by 1500 m. The finite-difference grid is oriented north– south, with 160 columns, 194 rows, and 16 layers for a total of 314,784 active cells (Fig. 2). Cell thicknesses range between 1 m and more than 3000 m, with layer 16 (the bottom layer) reaching 4000 m below sea level. The top elevation of layer one represents the water table (and modeled as a potentiometric surface) and not the actual land surface elevation with cell thickness in layer one increasing substantially in mountain ranges (Fig. 3). Fig. 4 shows the model layer configuration along cross section A–A0 (refer to Fig. 1 for cross section location) in comparison to land surface elevation and simulated hydrogeologic units. Twenty-seven hydrogeologic units were modeled in the DVRFSM via the MODFLOW2000 HUF2 package. A brief description of these units is provided in Table 1. In general, the model’s upper layers portray relatively shallow basin-fill sediments and volcanics, while deeper layers define the regional carbonate aquifer and confining units. Model cells containing more than one hydrogeologic unit are simulated with vertically averaged hydraulic properties.

Model calibration (Belcher, 2004) for steady state conditions prior to year 1912 and for transient conditions for years 1913– 1998 was accomplished by best matching 4899 observations of hydraulic head and drawdown as well as 49 estimates of groundwater discharge and 15 estimates of regional flux into/out of the model domain using MODFLOW-2000 parameter estimation processes (Hill et al, 2000). One hundred parameters were estimated in the calibration process including hydraulic conductivity, vertical anisotropy, drain conductance, fault-associated hydraulic conductivity, recharge zone multipliers, depth decay coefficients for hydraulic conductivity, and storage parameters (S and Sy). Fig. 5 provides an example of where boundary conditions, applied stresses (e.g. wells, drains) and faults (i.e., hydrologic flow barriers) are designated in the DVRFSM’s top layer during the final modeled timestep (year 1998). Fig. 6 shows recharge rates (m/d) applied throughout the modeled domain. Recharge is highest in mountainous regions, in particular the Spring Mountains contribute the greatest volume. Recharge was not modeled in most lowland areas, but those with designated recharge represent irrigation recharge or the re-infiltration of groundwater springs (e.g. Death Valley). Lateral flows into and out of the model boundary were specified as constant head cells that were interpolated from the regional potentiometric surface. Only the region bounding the Spring Mountains was defined as a no flow boundary. Drains were used to define both groundwater ET and spring flows. Discharge was defined by observation and drain conductances were determined by parameter estimation while drain elevations were 

set equal to 10 m below the lowest land-surface altitudes for each group of drain cells modeled. These elevations were assumed to adequately represent extinction depths for ET as well as account for springs being located in land-surface depressions (Belcher, 2004). Faults acting as barriers to flow were defined in MODFLOW-2000 with the Horizontal-Flow Barrier (HFB) package.

Faults were assumed to extend the entire thickness of the model and the width of each fault was assumed 1 m, while the hydraulic conductivity of each fault was determined by parameter estimation. Nine faults (i.e., horizontal flow barriers) were found to exert some influence on hydraulic heads and are shown in Fig. 5. The DVRFSM simulated hydraulic heads in layer one for the final time step (year 1998) are shown in Fig. 7.

数值模型

死亡谷地区流系统模型：MODFLOW

图1显示了DVRFS及其相关的USGS模型域（DVRFSM），

上面叠加了NTS。

DVRFS在内华达和加利福尼亚州约为100,000 km2，

其边界为北纬35°00′至38°15′，

西经115°00′至118°00′。

该系统包括中南内华达山区的补给区域

与NTS以南和西的湿盐滩和温泉的出流区域

以及加利福尼亚死亡谷之间的流动。

地下水流受到复杂的地质框架的影响，

USGS DVRFSM包括流系统主要含水层和封闭单元的分布，

以及可能影响地下流的主要地质结构（Belcher，2004）。

DVRFSM是使用MODFLOW-2000（McDonald和Harbaugh，1988；Harbaugh等，2000）建立的，

这是由USGS开发的一种能够模拟瞬态、三维、各向异性和异质系统中地下水流动的有限差分代码。

MODFLOW的受限含水层的三维流动方程结合了达西定律和质量守恒原理，

其表达式为：

\[ \]

其中，\(x\)、\(y\)和\(z\)是笛卡尔坐标（长度单位为L）

，\(K_{xx}\)、\(K_{yy}\)和\(K_{zz}\)分别是饱和水力传导率在\(x\)、\(y\)和\(z\)轴上的主要分量（L/T），

\(h\)是水头（L），

\(G\)是体积源/汇项（1/T），

\(S\)是储存率（无量纲），

\(b\)是饱和厚度（L，等同于单元厚度），

\(t\)是时间（T）。

相反，方程（2）定义了无受限条件下的地下水流动，

其中\(h\)（L）是井管地下水位高程，\(S_y\)是比流量（无量纲）。

\[ \]

在方程（2）中求解\(h\)是一个非线性问题，

因为\(h\)取决于水力传导率和饱和厚度，

而饱和厚度又取决于水位面的高程。

为了解决无受限含水层系统内在的非线性性，

MODFLOW的更新HUF2包维持了一个受限系统，

其顶部高程等于水位面的近似高程，

并将\(S\)替换为方程（1）中的\(S\)。

HUF2方法简化了无受限解的非线性，

仅在饱和厚度变化很小的情况下适用。

Belcher（2004）对DVRFSM进行了详细描述，

但为了清晰起见，

在此提供了对DVRFSM的综合审查。

DVRFSM在水平平面上采用了1500米×1500米的常数网格间距。

有限差分网格朝南北方向，共有160列、194行和16层，共计314,784个活动单元（图2）。

单元厚度在1米到3000多米之间变化，第16层（底层）达到海平面以下4000米。

第一层的顶部高程代表水位面（并被建模为势流面），

而在山脉中第一层的单元厚度显著增加（图3）。

图4显示了横截面A-A0（参见图1以获取横截面位置）沿线的模型层配置，

与地表高程和模拟的水文地质单元进行比较。

通过MODFLOW2000 HUF2包在DVRFSM中建模了27个水文地质单元。

表1提供了这些单元的简要描述。

总体而言，模型的上层描绘了相对浅的盆地沉积物和火山岩，

而深层则定义了区域性的碳酸盐含水层和封闭单元。

包含多个水文地质单元的模型单元使用了垂直平均的水力特性进行模拟。

在1912年之前的稳态条件和1913年到1998年的瞬态条件下，

通过使用MODFLOW-2000参数估计过程（Hill等，2000），

对DVRFSM进行了校准，

以最佳匹配4899个水头和降水观测，

以及49个地下水排泄和15个模型域内外区域通量的估计。

在校准过程中，估计了100个参数，包括水力传导率、垂直各向异性、排水导纳、与断层相关的水力传导率、补给区域倍增器、水力

传导率深度衰减系数和存储参数（S和\(S_y\)）。

图5提供了DVRFSM最终模拟时间步（1998年）中的边界条件、应用的应力（例如井、排水）和断层（即水文流动屏障）的示例。

图6显示了在整个模拟域中应用的补给率（m/d）

。补给在山区最高，特别是Spring Mountains贡献了最大的体积。

在大多数低地区域，未建模补给，但指定了补给的地区代表灌溉补给或地下水泉的再渗透（例如死亡谷）。

指定为常量头单元的模型边界处的横向流动进出由区域势流面插值得到。

只有Spring Mountains周边的区域被定义为无流边界。

排水用于定义地下水蒸散和泉水流。

观测定义了排水流量，

而排水导纳是通过参数估计确定的，

排水高程设置为模拟排水单元组的最低地表高程以下10米。

假定这些高程足以充分表示蒸散的消失深度以及泉水位于地表凹陷中（Belcher，2004）。

在MODFLOW-2000中，作为流动屏障的断层是使用水平流屏障（HFB）包定义的。

假定断层延伸到模型的整个厚度，

每个断层的宽度假定为1米，

而每个断层的水力传导率由参数估计确定。

在图5中显示了对九个断层（即水文流动屏障）产生一定影响的水力头的DVRFSM模拟结果。

图7显示了DVRFSM在最终时间步（1998年）的第一层中模拟的水力头。

SURFACT 

SURFACT (HGL, 2007) is a fully integrated flow and transport numeric code based on the USGS finite-difference groundwater model MODFLOW. SURFACT potentially improves upon MODFLOW with additional modules. Specifically, SURFACT’s modified BCF4 package performs a complex saturated-unsaturated subsurface flow analysis using a pseudo-soil function for a rigorous treatment of unconfined flow. This is done using robust numerical and matrix solution techniques for greater stability (PCG4) and an adaptive time-stepping scheme (ATO) to improve model efficiency. The three-dimensional movement of water in variably saturated media is given by (Huyakorn et al, 1986),

where additional terms relative to Eq. (1) include krw equal to relative permeability (a function of water saturation, dimensionless), u equal to a drainable porosity (i.e., specific yield Sy), and Sw representing the degree of water saturation (a function of pressure head, dimensionless). When fully saturated conditions exist, then Sw = 1.0, krw = 1.0, and Eq. (3) reduces to Eq. (1), the confined groundwater equation used by MODFLOW. Consequently, the conventional flow equation is used by SURFACT below the water table and in confined systems.

While MODFLOW’s other flow packages (BCF and LPF) support several averaging schemes (e.g. harmonic mean, arithmetic mean, arithmetic mean thickness-logarithmic mean hydraulic conductivity) for computing effective interblock transmissivity, MODFLOW’s HUF2 package is limited to only the harmonic mean (i.e., harmonic mean of transmissivity between cell i and cell j is defined as T_{i,j} = 2-T_iT_j/(T_i + T_j)). The effective transmissivity is equal to the saturated hydraulic conductivity multiplied by the saturated block thickness.

For confined systems, the saturated thickness is equal to cell thickness (b). The harmonic mean, by definition, removes the impacts of large outliers, aggregates smaller values and is best used to average  

fluxes. However, its estimation automatically underestimates the equivalent interblock conductivity (Romeu and Noetinger, 1995) by biasing toward the lower block value, and in the extreme case never allows a desaturated block (i.e saturated thickness of zero) to resaturated. In contrast, SURFACT’s BCF4 package computes interblock conductance as a product of the weighted harmonic mean of block saturated hydraulic conductivities (Kxx, Kyy, Kzz), relative permeability (krw), and mean flow area. The SURFACT approach to interblock conductance is more stable compared to MODFLOW’s unconfined solution, does not bias toward the lower block value and does not inactivate dry cells during desaturation (HGL, 2007).

SURFACT provides a fracture well (FWL) technique for more realistic simulation of groundwater pumping withdrawals compared to the original MODFLOW well (WEL) package. The original WEL package required the user to prescribe the contribution of water from each node (i.e., layer) in a multi-node well and did not automatically adjust pumping rates when water levels dropped below screened intervals. Instead, pumping produced withdrawals into the nonphysical realm. On the other hand, FWL automatically computes (and updates) relative nodal contributions of total well flux rates. For an unconfined system, the FWL package adjusts water withdrawal to feasible volumes from a potentially over-pumped well such that heads do not drop below a well’s screened interval.

SURFACT

SURFACT（HGL，2007）是一个基于美国地质调查局有限差分地下水模型MODFLOW的完全集成的流动和输运数值代码。

SURFACT通过额外的模块可能对MODFLOW进行改进。

具体而言，SURFACT的改进BCF4包使用伪土壤函数进行饱和-非饱和地下流分析，

用于对无围限流的严格处理。

这是通过使用强大的数值和矩阵解技术（PCG4）以及自适应时间步进方案（ATO）来提高模型稳定性和效率的。

在不同饱和介质中的三维水运动由以下方程给出（Huyakorn等，1986），

与方程（1）相比，

与方程（3）相关的附加项包括krw等于相对渗透率（水饱和度的函数，无量纲），

u等于可排水孔隙度（即比流量Sy），

Sw表示水饱和度的程度（压力水头的函数，无量纲）。

当存在完全饱和的条件时，Sw = 1.0，krw = 1.0，方程（3）就简化为方程（1），即MODFLOW使用的受限含水层方程。

因此，在水位下方和受限系统中，SURFACT使用传统的流动方程。

尽管MODFLOW的其他流动包（BCF和LPF）支持多种平均方案（例如，调和平均、算术平均、算术平均厚度-对数平均水力传导率）来计算有效的区块间渗透率，

但MODFLOW的HUF2包仅限于调和平均（即，单元i和单元j之间的渗透率的调和平均被定义为\[ \]）。

有效渗透率等于饱和水力传导率乘以饱和块厚度。

对于受限系统，饱和厚度等于单元厚度（b）。

调和平均，根据定义，消除了大的离群值的影响，汇总了较小的值，最适合平均流量。

然而，调和平均估计会自动低估等效的区块间导电性（Romeu和Noetinger，1995），

因为它倾向于偏向较低的块值，

并在极端情况下永远不允许非饱和块（即饱和厚度为零）重新饱和。

相比之下，SURFACT的BCF4包将区块间导电性计算为块饱和水力传导率（Kxx、Kyy、Kzz）、相对渗透率（krw）和均匀流动面积的加权调和平均的乘积。

SURFACT对区块间导电性的处理方式比MODFLOW的无围限解更稳定，

不偏向较低的块值，并在非饱和时不使干燥单元失活（HGL，2007）。

SURFACT提供了一种断裂井（FWL）技术，

相对于原始的MODFLOW井（WEL）包，更实际地模拟地下水提取。

原始的WEL包要求用户为多节点井中的每个节点（即层）规定水的贡献，

并且当水位降到筛选间隔以下时，

不会自动调整抽水率。

相反，抽水产生了进入非物理领域的提取。

另一方面，FWL自动计算（并更新）总井流速的相对节点贡献。

对于无围限系统，FWL包调整了从可能过度抽水的井中提取的水，以使水头不会降到井的筛选间隔以下。

Methods Translation of the DVRFS from MODFLOW to SURFACT Translation of MODFLOW’s geological model packages HUF2 and KDEP (Anderman and Hill, 2000; 2003), which accounts for hydraulic conductivity as a function of depth, into SURFACT’s modified BCF4 package was accomplished by extracting effective cell properties from the original DVRFSM (e.g. horizontal hydraulic conductivity, vertical hydraulic conductivity, and specific storage) using a Fortran code developed for this purpose and compiled as a subroutine of MODFLOW. Specific yield in the DVRFSM was originally calibrated between 0.19 and 0.20 depending on the hydrogeologic unit considered. A value of 0.19 was assumed reasonable and assigned across the entire SURFACT model domain.

Fig. 8 shows hydraulic conductivity (K) in the model’s top layer as an example of effective values used by SURFACT.

Comparison between MODFLOW and SURFACT The SURFACT DVRFS model was run using USGS-calibrated parameters for the MODFLOW DVRFSM. Comparison between the unconfined and confined modeling approaches included: (1) areal head distributions at the end of the transient model (year 1998), (2) ability to predict observed values of heads, drawdown, groundwater fluxes to springs (i.e., modeled drains) and across modeled boundaries (i.e., constant heads), and (3) estimated drawdowns related to 50 years additional pumping and 50 years recovery beginning in year 1998 from a hypothetical well south of the NTS, a region often considered for groundwater development to support Las Vegas growing need for water.

Results and discussion

Predicted heads

Computed heads at the end of the transient simulation (year 1998) found both MODFLOW and SURFACT consistent in their estimates of flow magnitude and direction. However, differences in predicted heads did occur. In general, the unconfined modeled heads were higher than those predicted by the confined model in mountainous regions with greatest differences occurring in mountain ranges defined with thick model cells (b), low hydraulic conductivity (K) and high recharge flux (w). Specifically, SURFACT predicted heads in excess of 100 m compared to MODFLOW in the Spring Mountains, the Sheep Range, Pintwater Range, Timpahute Range, and Last Chance Range. SURFACT predicted heads on the order of 10–100 m above MODFLOW predicted heads occurred within the Amargosa Range, Black Mountains, Greenwater Range, Belted Range, Black Mountain and the Yucca Mountain region.

These regions are defined with thinner cell thicknesses, higher hydraulic conductivities and lower recharge values than other mountain regions. SURFACT and MODFLOW predicted heads were relatively similar for most lowland regions in the model, with SURFACT predicted heads less than MODFLOW in Cactus Flat, Penoyer Valley, Tikaboo Valley, the southern end of the Amargosa Desert (i.e., Amargosa Flat) and the northern end of Death Valley. Differences were on the order of 1.0–10.0 m.

To understand discrepancies between modeling approaches it is necessary to revisit the groundwater flow equation. The analytical solution for a confined aquifer defined by Eq. (1) and assuming 1dimensional, steady state flow given homogenous, isotropic conditions is

h=h_0-\frac{wx^2}{2kb}+(\frac{h_1-h_0}{L}+\frac{2L}{Kb})x

where h0 and h1 equal heads (L) at x_0 = 0 and x_1 = L, respectively and w is recharge (L/T). Assuming h_0 = h_1 = 0 and that the groundwater divide (q'_ x=-Kb(dh/dx)= 0 )occurs at x = L/2 will produce a maximum head (hmax) equal to,

DVRFS从MODFLOW到SURFACT的方法

将MODFLOW的地质模型包HUF2和KDEP（Anderman和Hill，2000; 2003）

翻译成SURFACT的改进BCF4包，

该包考虑了水力导度随深度变化，

通过提取原始DVRFSM中的有效单元属性

（例如水平水力导度，垂直水力导度和比存储）完成，

使用专门为此目的开发的Fortran代码，

并编译为MODFLOW的子程序。

DVRFSM中的比流量最初在0.19到0.20之间进行了校准，

具体取决于考虑的水文地质单元。

假定0.19是合理的，并在整个SURFACT模型域中分配。

图8显示了作为SURFACT使用的有效值的模型顶层中的水力导度（K）的示例。

MODFLOW和SURFACT之间的比较

使用为MODFLOW DVRFSM校准的参数运行了SURFACT DVRFS模型。

无围限和有围限建模方法之间的比较包括：

（1）瞬态模型末端（1998年）的面积头分布，

（2）能够预测观测值的头，降水，对泉水的地下水通量（即，模拟排水）和模拟边界（即，常量头）以及

（3）与1998年开始的50年额外抽水和50年恢复相关的估计降水，来自NTS南部的一个假设井，这个区域通常被考虑为支持拉斯维加斯不断增长的水需求的地下水开发的地方。

结果和讨论

预测的水头

在瞬态模拟结束时（1998年），

计算出的水头发现MODFLOW和SURFACT在其流量大小和方向的估计上是一致的。

然而，确实存在预测水头的差异。

总体而言，在山区，未围限的建模水头比围限模型预测的水头高，最大的差异出现在定义有厚模型单元（b），低水力导度（K）和高补给通量（w）的山脉。

具体而言，在Spring Mountains，Sheep Range，Pintwater Range，Timpahute Range和Last Chance Range，SURFACT预测的水头超过MODFLOW的100 m。

在Amargosa Range，Black Mountains，Greenwater Range，Belted Range，Black Mountain和Yucca Mountain地区，SURFACT预测的水头在10-100 m的数量级以上的MODFLOW预测的水头。

这些地区的单元厚度较薄，水力导度较高，补给值较低。

在模型中的大多数低地区，

SURFACT和MODFLOW预测的水头相对相似，SURFACT在Cactus Flat，Penoyer Valley，Tikaboo Valley，Amargosa Desert的南端（即Amargosa Flat）和Death Valley的北端比MODFLOW低。

差异在1.0-10.0 m的数量级。

要理解建模方法之间的差异，有必要重新审视地下水流方程。

对于由方程（1）定义的受限含水层的解析解，假设均匀、各向同性条件下给定一维、稳态流是

\[\]

其中h0和h1分别等于和处的水头（L），

w是补给（L/T）。假设，

并且地下水分水岭（\[\]）

发生在\[\]，将产生一个最大水头（\[\]）等于

​​​​​​​

\[ \]