1、Your First Groundwater Model with PMWIN实例:一个由两种不同地层构成的含水系统,北部和南部为无流量边界,东部和西部为河流,与地下水有充分的水力联系,可视为定水头边界,西边界与东边界的水头分别为9m,8m。Fig. 8-1 Configuration of the sample problem含水系统为非承压各向同性,第一层、第二层的水力传导系数分别为 0.0001m/s 和0.0005m/s,两层的纵向传导系数假设为横向的 10%,有效孔隙率为 25%,地面标高为10m。第一层、第二层的厚度分别为 4m、6m。含水层接受面状补给量为 8109 mm。在第一
2、层中靠近西部边界有一污染源。现在我们的任务是,如何用一个位于东部边界附近的完整井来隔离该污染源。那么,在该区域内,需要建立一个数学模型来逐步计算抽水井所需的抽水量。只有当抽水量足够大时,被污染的地下水才能被抽水井全部截获。我们将使用 PMWIN 来建立数学模型,使用 PMPATH 来计算抽水井的截获带。基于所计算的地下水流场,使用 MT3D 和MOC3D 来模拟污染源的运移情况。并演示如何使用 PEST 和 UCODE 来校正水流模型。最后创建一个动态演示方式,演示地下水被污染的发展变化过程。为了示范运移模型,我们假想污染物溶解于地下水中的速度为 110-4g/s/m2,含水层纵向弥散和橫向弥
3、散系数分别为 10m 和 1m,延迟系数为 2,初始浓度、分子扩散系数及衰减速率均为 0。我们将计算三年之后污染浓度分布情况,显示在两个含水层中点x,y=290,310,390,310的浓度时间突变曲线。稳定流模拟稳定流模拟中的六个主要步骤:1. 建立模拟模型2. 模型赋值3. 执行水流模拟4. 模拟结果检验5. 子区水均衡计算6. 输出结果。步骤一:建立模型1、打开 File 菜单,选择 New Model。New Model 对话框弹出。选择一个保存模型数据的文件夹,如:C:PM5DATASAMPLE,键入文件名:SAMPLE 作为示例模型。模型的文件扩展名必须是“.PM5”。在 Wind
4、ow95/98/NT 下,文件名的有效使用字符数为 120。最好在一个单独的文件夹中保存每一个模型及输出数据。可以同时运行多个模型(多任务处理) 。2、点击 OK。PMWIN 只需几秒钟即可建立一个新模型,模型的文件名显示在标题栏中。步骤二:模型赋值包括:创建模型网格、定义边界条件、模型单元赋值。在整个模型中,PMWIN 要求单位使用统一。例如,若以 m 表示长度单位,s 作为时间单位,水力传导系数必须以 m/s 表示,抽水量必须以 m3/s 表示,弥散系数必须以 m 表示。在 MODFLOW 中,一个含水系统可以由一个包含一系列节点的离散区域或一些关联的有限差分块取代。图 8-2,显示了一个
5、包含网格和节点的含水系统的空间离散化情况,在每一节点均可计算出水头值。节点网格形成了数值模型的框架。一个地质单元可由一个或多个模型层表示。每一计算单元格的厚度、长度、宽度均可改变。计算单元所处的位置可以用列、行、层来表示,PMWIN 使用索引符号J,I,K来指定计算单元的位置。如:位于第 2列、第 6 行、第 1 层的计算单元记为2,6,1。Fig. 8-2 Spatial disretization of an aquifer system and the cell indices一、 创建模型网格1、打开 Grid 菜单,选择 Mesh Size。Model Dimension 对话框弹出
6、。2、输入含水层的层数 3,行数 30、列数 30,行高、列宽 20。第一含水层和第二含水层分别由一个模型层和两个模型层来代替。3、点击 OK。PMWIN 的界面发生变化,此时显示了模型的整体网格。PMWIN 允许用户改变和旋转模型网格。用户可以随意改变模型网格行、列宽度,也可以增加/删除行和列。在该实例中,不需更改模型的网格。4、从 File 菜单中选择 Leave Editor 或点击 Leave Editor 按钮退出。二、定义含水层类型1、打开 Grid 菜单,选择 Layer Type。Layer Options 对话框弹出。2、点击标签为 Type 列中的某一单元格,单元格中将出现
7、一个带下箭头的按钮,点击该按钮,出现一个可供选择含水层类型的列表。3、为第一层含水层选取 1:Unconfined,其它含水层选取 0:Confined,然后点击 OK,关闭对话框。三、定义边界条件在模型中,边界上的各个计算单元均用不同的指定代码表示不同的边界类型:(1) “1”表示变水头边界, (2) “1”表示定水头边界, (3) “0”表示无流量边界。1、从 Grid 菜单中选择 Boundry ConditionIBOUND(Modflow)。PMWIN 的 Data Editor 界面弹出,界面显示一个模型网格平面图,单元格指针位于1,1,1 ,即第一层的左上角。当前格的数值显示在状
8、态栏底部。IBOUND 默认值为 1。单元格指针可以使用键盘的箭头键移动或使用鼠标点击你所需要的位置。使用 PgUp、PgDn,或在工具栏的编辑框中输入含水层的层号后按 Enter 键,可以从一层移到另一层。2、击鼠标右键,PMWIN 显示一个 Cell Value 对话框。3、在对话框中输入1,点击 OK。计算单元1,1,1被定义为定水头计算单元。4、点击 Duplication 按钮,打开复制开关。如果 Duplication 键下陷,表示复制开关处于开启状态。移动网格指针,可将当前计算单元的数值复制到单元格指针经过的所有计算单元,再次点击 Duplication 键,关闭复制开关。5、从
9、左上角1,1,1移动单元格指针至左下角1,30,1,-1 值被复制到模型西部边界的所有计算单元中。6、移动单元格指针到右上角30,1,1。7、再将单元格指针从右上角30,1,1移至右下角30,30,1,-1 值被复制到模型东部边界的所有计算单元中。8、点击 Layer Copy 按钮,打开层复制开关。该按钮下陷,表示层复制开关处于开启状态。此时,从一层移至另一层,当前层的所有值将被复制到另一层中。再次点击该按钮,将关闭层复制开关。9、按 PgDn 键,从第一层移至第二层、第三层。第一层的值被复制到第二层、第三层。10、点击 Leave Editor 按钮或从 File 菜单中选择 Leave
10、Editor 退出定义边界条件界面。四、模型几何参数赋值(一)模型层的顶板高度赋值1、打开 Grid 菜单,选择 Top of Layer(Top)。PMWIN 显示模型网格。2、打开 Value 菜单选择 Reset Matrix 或按 Ctrl+R 键。Reset Matrix 对话框出现。3、在对话框中输入 10,点击 OK。第一层的顶板高度定义为 10。4、按 PgDn 键移至第二层。5、重复步骤 2、3,输入第二层、第三层顶板高度 6、3。6、从 File 菜单中选择 Leave Editor 或点击 Leave Editor 按钮退出。(二)模型层的底板高度赋值1、从 Grid 菜
11、单中选择 Bottom of Layer(BOT)。2、重复上述步骤,输入第一层、第二层、第三层的底板标高值 6、3、0。3、点击 Leave Editor 按钮退出。五、时间及空间参数赋值在本实例中,时间参数包括时间单位、应力期数、时段、水质运移时段。空间参数包括初始水头、水平及垂直方向的水力传导系数、有效孔隙率。(一)时间参数赋值1、打开 Parameters 菜单,选择 Time。Time Parameters 对话框出现。在 MODFLOW 中,模拟时间被划分为几个应力期,也就是说,在外部应力作用下,时间间隔为常数,应力期又划分为时段。在大多数水质运移模型中,每一水量模型的时段又被分成
12、为更小的运移时段,时段的长度并不与稳定流模型对应。然而,当我们使用 MT3D 和 MOC3D 来运行污染运移模拟模型时,实际的时段长度需在表中给定。2、输入第一应力期的长度 9.46728E+07(秒) 。3、点击 OK 接受其它默认值。(二)初始水头赋值1、打开 Parameters 菜单,选择 Initial Hydrolic Head。PMWIN 显示模型网格。2、从 Value 菜单中选择 Reset Matrix(或按 Ctrl+R 键) ,在对话框中输入 8,点击 OK。3、将网格指针移至左上角1,1,1。4、点鼠标右键,在对话框 Cell Value 中输入 9。5、点击 Dup
13、lication 按钮,打开复制开关。6、从左上角向下拖动鼠标,当前计算单元中的值 9 将被复制到鼠标所经过的所有计算单元中。7、点击 Layer Copy 按钮,打开层复制开关。8、按 PgDn 键,移至第二层、第三层。9、点击 Leave Editor 按钮,退出。(三)导水系数赋值1、从 Parameters 菜单中选择 Vertical Hydraulic Conductivity。2、打开 Value 菜单,选择 Reset Matrix(或按 Ctrl+R) ,在对话框中输入 0.00001,点击OK。3、按 PgDn 键,移至第二层。4、重复步骤 2、3,分别在第二层、第三层中输
14、入 0.00005。5、点 Leave Editor 按钮退出。(四)有效孔隙率赋值1、从 Parameters 菜单中选 Effective Porocity。输入有效孔隙率值,并保存。模型的默认值为 0.25。2、退出。(五)面状补给量的赋值1、打开 Models 菜单,选择 MODFLOWRecharge。2、从 Value 菜单中选择 Reset Matrix,在对话框的 Recharge FluxL/T编辑栏中输入 8E-9。3、退出。(六)定义抽水井位置及抽水量赋值在 MODFLOW 中,一个注水井或抽水井由一个节点(或计算单元)取代,用户可用一个节点定义一个注水井或一个抽水井,并
15、且假设井穿透整个含水层。MODFLOW 可以模拟穿透多个含水层的抽水井、或抽水井在每一层中有不同的抽水量。多层井的总抽水量等于各单层抽水量之和。每一层的抽水量可以根据下式计算: TtoalkkQTk:含水层的导水系数;T:被井穿透所有含水层的导水系数之和。但是,当第一层含水层为非承压时,如果不假设一个饱水厚度,我们无法确切知道井所在位置含水层的饱水厚度和导水系数,上述公式无法应用。模拟一个多层井的另一种可能是:对井所在的每一个计算单元设一个很大的垂向导水系数(或垂向渗漏) ,如 1m/s,总抽水量赋予井的最底部计算单元。为了演示,在井中的其它计算单元中赋一个非常小的抽水量(110 -10m3/
16、s) ,这样,来自每个被穿透含水层的水量将被 MODFLOW 计算出,并且该值可通过 Water Budget Calculator 来获得。由于我们不能确定截获所有污染物所需的抽水量,我们假设一抽水量 0.0012m3/s 来进行调试。1. 打开 Models 菜单,选择 MODFLOWWell。2. 将网格指针移到25,15,1。3. 点鼠标右键,输入-1E-10,点击 OK。负值表示抽水井。4. 按 PgDn 键,移到第二层。5. 点鼠标右键,输入-1E-10,点击 OK。6. 按 PgDn 键,移到第三层。7. 点鼠标右键,输入-0.0012 ,点击 OK。8. 退出。步骤三:执行水流
17、模拟1、打开 Models 菜单,选择 MODFLOWRun.。Run Modflow 对话框弹出。2、点击 OK,开始水流模拟。在运行 MODFLOW 之前,PMWIN 将利用用户自定义的数据为 MODFLOW(及可选MODPATH)产生一个输入文件,该文件列在 Run Modflow 对话框的列表中。如果Generate 标记为 ,表示只产生了一个输入文件。可以点击该按钮锁定 Generate 标记。一般情况下,不需更改该标记。步骤四:模型结果检验在水流模型模拟中,MODFLOW 在列表文件 pathOUTPUT.DAT 中有一个详细运行记录。如果成功地完成了一个水流模拟,MODFLOW
18、将以多种非格式化文件(二进制)保存模拟结果。如下表。在运行 MODFLOW 之前,用户可以控制这些非格式化(二进制)文件的输出。具体操作:从 Models 菜单中选择 ModflowOutput Control。若 Interbed Storage Package 被激活,则产生 pathINTERBED.DAT 输出文件。为了检验模拟结果的质量,MODFLOW 在每个时段结尾为整个模型作了水量均衡计算,并保存在 OUTPUT.DAT 文件中。一个水均衡提供了一个数学解的所有可解性迹象。MODFLOW 的输出文件文件 内容PathOUTPUT.DAT 详细运行记录和模拟报告PathHEADS.
19、DAT 水头PathDDOWN.DAT 降深、初始水头与计算水头之差PathBUGET.DAT 各计算单元的流量项PathINTERBET.DAT 每一层中整个含水层的沉降、压密及压密之前的水位PathMT3D.FLO 与 MT3D/MT3DMS 的界面文件,这个文件由MT3D/MT3DMS 提供的软件包 LKMT 创建。Path 为保存模型数据的文件夹。步骤五:子区水均衡计算为了计算方便,每个计算单元的流量项被保存在文件 pathBUDGET.DAT 中,这些单一计算单元的流量项实际上是一个单元流向一单元的流量,有四种类型:(1)一个单元向一个单元的应力流量;或受外部应力影响,流入或流出一个
20、单元的流量,如抽水井或面状补给;(2)计算单元间的储存项:在各计算单元中,储存量的增加和减少的速率。 (3)计算单元间的定水头流量项:流向或来自各个定水头计算单元的净流量。 (4)计算单元间内部的流量项:通过各计算单元截面的流量,也即相邻计算单元间的流量。Water Budget Calculator利用逐个计算单元的流量来计算整个模型的水量均衡、子区及子区间的水量均衡。子区水均衡计算1打开 Tool 菜单,选择 Water Budget。Water Budget 对话框弹出。在稳定流计算中,不需改变 Time 的设置。2点击 Zones。一个 Zone 代表一个水均衡计算子区,由 050 间
21、的一个代号(区号)表示,每个区的区号必须赋予每一个计算单元。区号 0 表示该计算单元与其它区无关。3 Value 菜单打开 Reset Matrix对话框,输入 1,点击 OK。4 PgDn 键,移到第二层。5 从 Value 菜单中,选择 Reset Matrix,输入 2,点击 OK。6 点击退出按钮,退出。7 点击 OK,退出 Water Budget 对话框。PMWIN 在文件 pathWATERBDG.DAT 中计算并保存流量项。流量单位为L 3/T。对每一时段、每一层的每一个区分别计算。流入区内的流量记为 IN,子区之间的流量记在 Flow Matrix 中,通过子区边界的水平方向
22、流量由 HORIZ EXCHANGE 给出。EXCHANGE( UPPER)给出了来自( IN)或流向(OUT)上部相邻层的流量,EXCHANGE( LOWER)给出了来自(IN)或流向(OUT)下部相邻层的流量,如LAYER1,ZONE 1 的 EXCHANGE(LOWER) (表示第一层、第一区与下部相邻层的流量交替项) ,从第一层到第二层为流量为 2.587256E-03m3/s,计算结果中的 PERCENT DISCREPANCY 可由下式计算出:2/)(0OUTIN在该实例中,整个模型中每层、每区的流入、流出量的百分比误差很小,这意味着模型方程有正确的解。为了计算流向抽水井的确切流量
23、,重复上述步骤计算水量均衡。这次,只需对计算单元25,15,1定义为 1 区,计算单元25,15,2定义为 2 区,计算单元25,15,3定义为 3 区,所有其它计算单元均定义为 0 区,水量计算结果列在表 2.4 中。抽水井在第一层中的抽水量为7.7992809E-05m3/s,第二层为 5.603538E-04m3/s,第三层为 5.5766129E-04m3/s,几乎所有抽水量均来自第二层,与实际含水层结构吻合。步骤六:结果输出除水量均衡计算之外,PMWIN 提供了多种模拟结果检验、创建图形输出的可能。流线及速度矢量可以通过 PMPATH 显示。利用 Result Extractor,可
24、将任一时段任一层的模拟结果从非格式化(二进制)结果文件中读取,或保存在 ASCII Matrix 文件中。一个 ASCII Matrix 文件含有一个模型层的计算值(每计算单元一个值) 。PMWIN 能将 ASCII Matrix 文件装入模型网格中。以下给出操作步骤:1 用 Results Extractor 读取或保存计算的水头。2 于第一层的计算水头生成等水位线图。3使用 PMPATH 计算水流运行轨迹及抽水井的截获带。(一)读取及保存计算的水头1、从 Tool 菜单中选取 Results Extractor。Results Extractor 对话框弹出。对话框中含有四个选项标签,分别
25、为MODFLOW、MOC3D 、MT3D、MT3DMS。在 MODFLOW 工作表中,可从 Result Type下拉列表中选择结果类型,也可以指定要读取数据结果的层、应力期、时间步长。该电子表格显示了一系列的行和列,一行与一列的交汇处为一个单元格。其中的每一单元格与模型层中计算单元一一对应。下列步骤 2、6 以三个 ASCII Matrix 文件为每一层保存了水头值。2、点击 Result Type 的下拉箭头按钮,选择 Hydraulic Head。3、在 Layer 编辑框中输入 1。该实例中(稳定流模拟只有一个应力期一个时段)应力期数及时段数均为 1。4、点击 Read 按钮。第一层、
26、第一应力期,时段 1 的水头值将被读入电子表中。点击滚动条,可游览整体结果。5、点击 Save 按钮。Save Matrix 对话框出现。在 Save as Type 选项中,可根据用户需要将数据结果保存为ASCII 格式或 SUFFER 格式。给定文件名 H1.DAT,保存文件,点击 OK。6、重复步骤 3、4、5,分别保存第二层、第三层的水头计算结果 H2.DAT,H3.DAT 。7、点击 Close,关闭对话框。(二)生成计算水位等值线图1、从 Tools 菜单中选取 Presentation。在 Presentation 中定义的数据并不能被 PMWIN 使用,但可以用 Present
27、ation 来临时保存数据或显示模拟结果的图形。2、从 Value 菜单中选择 Matrix(Ctrl+B) 。Brows Matrix 对话框弹出。点击 Load 按钮,可将 ASCII Mtrix 文件装入表中,或点击Save 按钮,将表中数据保存为一个 ASCII Matrix 文件。还可通过菜单 Value 中的 Result Extractor 读取水位数据,或使用 Apply 按钮在 Presentation 数据表中添加数据。3、点击 Load 按钮。Load Matrix 对话框出现。4、点击 Open 按钮,选取由 Result Extractor 保存下的文件 H1.DAT
28、,点击 OK,H1.DAT 被装入。5、点击 OK,关闭 Browse Matrix 对话框。6、打开 Option 菜单,选择 Environment(Ctrl+E) 。Environment Option 对话框出现。对话框中包括三个选项列表。Appearance 和 Coordinate System 用来修改模型网格的外观及位置,Contours 用来生成等值线图。7、点击 Contour 标签,选中 Visible 复称选框,点击 Restore Default 按钮。点 Restore Default 按钮,PMWIN 默认等值线数为 11,并使用当前层中的最大和最小值定义等值线的最大、最小值。如果 Fill Contours 被选取,点击列标签 Fill,可自选颜色填充等值线图。用 Label Format,可定义合适的数据格式。注意:当退出编辑时,PMWIN 将清除 Visible 复选框。8、点击 OK,退出 Environment Options 对话框。9、从 File 菜单中选择 Save Plot As或 Print Plot,保存或打印图形文件。10、按 PaDn 移至第二层,重复步骤 29,装入文件 H2.DAT,显示和保存图形文件。11、点击 Leave Editor 按钮,然后点击 Yes,退出并保存 Presentation 结果。