第四章 运行选项 c;l!i-
本章将详细给出: NR4+&d
¨ 以工程或批处理方式运行模型; jK|n^5\
¨ 编辑层宽变量; 0SQ!lr
¨ 编辑运行时间选项。 (MwRe?Ih
4.1运行模块菜单 Ji4p6$ .j-
从主菜单的下拉菜单中选择Run进入运行屏幕后,会出现下列对话框: m,.Y:2?*V
用户必须选定当前进行的模拟是瞬变流还是稳定流。选定后,用户可以看到流动运行模块的下拉式菜单的选项(见图)。 +
G@N
4.1.1流动的运行模块菜单 Y;i=c6
File Basic Solver Recharge Layers BCF OC Pathlines Run/Translate MT3D Tools Help Kk3+ ]W<
File 保存或返回主菜单。 2+1ybOwb
Basic 设立模型时间步参数(如果模拟瞬变流),指定初始水头(或者是每层保持不变,或者从Surfer,ASCII文件,或前一次Visual MODFLOW 运行结果中输入)。 $8tk|uh
Solver 选择MODFLOW 求解器包。 '&IGdB I
Recharge 给模型顶层或最上面的活动单元指定初给量。 K~W(ZmB
Layers 按MODFLOW 手册所述指定模型中各层的类型。 ?Ga8.0Z~KT
BCF 定义块中心流动选项:单元二次湿润和层的各向异性。 9YzV48su#
OC 指定保存结果(输出控制)的时间和格式。 s9- qR_
Pathlines 设置用于USGS MODPATH的质点选项。 C6!F6Stn]g
Run/Translate 将Visual MODFLOW 文件转换成MODFLOW、MODPATH、ZONE BUDGET和MT3D格式的文件以及按要求运行模型。
RD$:.
MT3D 转至Visual MODFLOW 的MT3D运行模块。 zakhJ
Tools 激活单元检查器。 dlu*s(O"
Help 运行菜单的通用帮助屏幕。 wJp<ZL
4.1.2 MT3D的运行模块菜单 xS*UY.>
从流动运行模块的菜单上选择MT3D菜单项,设置进行MT3D模拟的选取项。用户将会见到如下的MT3D运行模块的下拉式菜单(见图)。 HsY5wC
File Advection Output/Time Steps Initial Concs. Run/Translate Flow Tools Help !_fDL6a-
File 保存或返回主菜单。 6o't3Peh
Advection 设置MT3D的求解方法的选项。 U4D7@KY +m
Output/Time Steps设置将MT3D结果保存到文件的时间。 l;-Ml{}|0
Initial Concs. 告诉MT3D使用哪种初始浓度(在MT3D模块中指定、从SURFER、ASCII文件或前一次MT3D解中输入)。 t7=D$ua
Run/Translate 将Visual MODFLOW 文件转换成MODFLOW、MODPATH、ZONE BUDGET和MT3D格式的文件以及按要求运行模型。 2Tp2{"sB>A
Flow 转至Visual MODFLOW 的流动运行模块。 S?~0)EXj(
Tools 激活单元检查器。 /%@;t@BK4
Help 运行菜单的通用帮助屏幕。 >eJ<-3L;
4.2文件(File)菜单 gZ @+62
在运行菜单中选择File,会出现一下拉菜单: J8ni}\f
Save 以原数据集名保存选项。 I G1];vX
Export... 将Visual MODFLOW 图形窗口输出到AutoCAD的DXF文件(世界坐标或局部坐标)或WMF文件。 V_7\VKR
Main Menu 返回主菜单。 P9v(5Z00|d
4.3基本资料(Basic)菜单 mLCDN1UO{
在运行菜单中选择Basic,会出现一下拉菜单: 0ho;L 0Nr'
Time 定义每个应力期的时间步数和时间步乘数。本项只在瞬变流模式下才有效。 (0-Ol9[
Initial Heads 定义起始水头。其选项有:每层水头保持不变、从SURFER的.grd文件中输入水头、以及使用前一次的Visual MODFLOW 运行轮次、SURFER的.grd文件或ASCII文件中的水头。 .j}]J:{%
4.3.1时间(Time) (x=$b(I
如果选择的是瞬变流模式,那么你的Basic下拉菜单下就会有Time选项。你会看到如下的应力期(Stress Period)窗口,你可以输入每个应力期的时间步数和时间步长乘数。 F*72g)hVh
MODFLOW 中的模拟时间被Visual MODFLOW 分为应力期。应力期的定义为系统上所有的应力(边界条件、抽水速率等)保持不变的时间期。Visual MODFLOW 利用用户加入的边界条件来确定瞬变流模拟的每个应力期的长度。应力期限的长度不能由用户直接设置。接着每个应力期又可以分为时间步。在每个应力期内,时间步组成几何级数。用户要定义时间步数的时间步长乘数,或者说是每个时间步与前一时间步的比例。Visual MODFLOW 利用这些定义来计算应力期内每个时间步的长度。 RQVu~7d[
MODFLOW 在各应力期、时间步和打印时间上是有区别的。每个时间步都会计算水头场。但水头场只在每个打印时间(printout time)结束时才保存到某个文件中去。在默认情况下,打印时间就是应力期的结束,而用户可以通过OC菜单的选项对此进行修改。但在输入模块中指定的观测点的水头和降深值却是在每个时间步都计算并加以保存。 3j7FG%\
一般情况下,时间步长越小,其对流场的瞬时变化的刻画就越好。因此增加时间步数能使水头或降深的曲线更加光滑。 e@D_0OZ
4.3.2初始水头(Initial Heads) '|8dt "C
对稳定流模式,MODFLOW 需要一个水头分布的初始猜测值,而对瞬变流模拟,MODFLOW 需要起始水头分布。降深值也从初始水头中算得。对稳定流模拟,一个好的猜测值可以极大地减少所需的运行时间。从Basic菜单中选择Initial Heads后,会出现如下弹出窗口。 EPm~@8@"j?
它允许用户指定初始水头分布是否是按层保持不变,或者从SURFER的.grd文件、ASCII文件还是前一次的Visual MODFLOW 运行结果中输入。 : auR0FE
4.3.2.1按层保持不变 4XkI? l
按层保持不变是默认选择。如果没有选择Initial Heads,程序就会用它来模拟。如果选择了“Constant by Layer ”,就会出现如下的窗口。 7|,L{~
该窗口允许用户按层设置初始水头的猜测值。 : |'(T[~L
4.3.2.2从SURFER的.grd文件或ASCII文件中输入 (r|m&/
Visual MODFLOW 允许初始水头从已有的数据中插值求得。如果水位资料可用,它就能从SURFER的.grd文件或ASCII文件中插值出来,将有一个与如下窗口类似的窗口出现。 sJ6.3=
c
如果选择了ASCII文件,其类型可以是ASC,TXT或XYZ。然而,文件必须以空格为分隔符,每组X,Y,Z坐标占一行。对SURFER,文件必须是GRD格式,其扩展名必须是GRD。 8>KUx]AN
ASCII和SURFER选项都是用距离倒数插值方法确定每个单元的初始水头的。 1lw%RM
4.3.2.3前一次Visual MODFLOW 运行结果 ~\":o:qyc
“Previous Visual MODFLOW Run”选项允许用户将上一次的HDS文件作为初始水头的猜测值。选择该项后,将出现如下窗口。 DDE-$)lf>
HDS文件必须与当前的模型具有相同的行、列和层数。这在模型校正和敏感性分析时是十分有用的,因为它为模型提供了接近于解的猜测值。如果网格作了重新剖分,旧的HDS文件就不能再用了。 4 m:h&^`N
4.4求解器(Solver)菜单 X[B P0:`t
选择Solver后就会出现求解器的选项表列。对每个求解器都会出现一个不同的“Solver Parameter”对话框。你可根据所需的计算水头分布的精度在其中修改默认的收敛标准。一般情况下,默认值就可以得到较好的结果。然而,我们推荐用户去试验这些数,以确信迭代的解确实是收敛的。这也是用户的责任所在。 R)NSJ-A!2
选择Solver后就会出现如下对话框。 @MSmg3&
4.4.1 前承条件共轭斜量程序包(PCG2) C- .;m
PCG2用前承条件共轭斜量法(Preconditioned Conjugate-Gradient)来求解模型的联立方程。可以模拟线性和非线性的流动条件。PCG2有两个前承条件选项:改进的不完全Cholesky前承条件,它在标量计算机上效率很高;多项式前承条件,它对计算机的存储要求不高,并且通过针对计算机的优化,在矢量计算机上效率极高。求解器的收敛条件采用水头变化和残差判别标准。非线性问题采用Picard迭代算法。对PCG2的说明见USGS的Water-Resource Investigations Report 90-4048,它也已包含在用户参考手册中。 s.J4&2Q
对求解器参数的解释 K%q5:9m
PCG2求解器在每个时间步上用双层迭代法求解。外层迭代用于向解的方向改变前承条件参数矩阵。每进行一次外层迭代,前承条件参数矩阵中的流动系统的水文地质参数(即水力传导系数、饱水厚度、储水系数)就作了一次更新。 `/O`%6,f1!
最大外层迭代次数(Maximum Number of Outer (non-linear) Iterations)默认值=50 n!)$e;l
本参数指定外层迭代次数的上限。这个值只有在迭代结束前不能达到收敛解的时候才用到。对大多数问题而言,25次就足够了。然而,如果达到了最大迭代次数,而不能得到较好的质量平衡差异,就应加大该值。 R%UTYRLUn
最大内层迭代次数(Maximum Number of Inner Iterations)默认值=10 Gwd38
本参数指定内层迭代次数的上限。这个值只有在迭代结束前不能达到当前外层迭代的矩阵条件下的收敛解的时候才用到。对大多数问题而言,10次就足够了。超过10次的迭代往往对解并没有什么用处,因为在返回外层迭代时,解要进行更新。 .}IW!$
dq
收敛的水头变化判别标准(Head Change Criterion for Convergence):默认值=0.01 O}M-6!%<,
在所有的外层迭代完成后,求解器就会检查每个单元的解,找出最大的变化量。如果这个值小于收敛标准(以英尺或米为单位)。就认为解已收敛,求解器就停止运行,否则开始新的外层迭代过程。 W[2]$TwT
通常对大多数问题来说,解的精度达到0.01ft或m就足够了,除非模拟区域的水头的最大差值小于1(ft或m)。如果内外层迭代的次数没有达到上面所设置的最大值,却没有得到较好的质量平衡结果,可以将这个值减少一个数量级,如0.001。 aODh5
收敛的残差判别标准(Residual Criterion for Convergence):默认值=0.01 pz%s_g'
水头变化的判别标准用于判断整个求解器的收敛,而残差判别标准则用于判断求解器内层迭代的收敛情况。如果所有节点上的残差的最大绝对值小于这里指定的值(单位:长度3/时间),那么就进行下一个外层迭代。 7l *
&Fh9;
通常,水头变化和残差判别标准是等价的。一般0.01ft3/s或m3/s就足够了。然而,如果你注意到所有的外层迭代都只运行了一些内层迭代,而却没有得到较好的质量平衡结果,则可以将该值减少一个或几个数量级。 TgiZ
% G
注意:残差判别标准与所用的单位是有关的。如果你所用的长度单位是英尺或米,时间单位是秒,则0.01就足够了。如果不是,你就应该将默认值乘上你所用的时间单位的秒数(例如,如果时间单位是天,则残差判别标准应采用864)。 2<D| {
打印间隔(Printout Interval):默认值=10 X^\D"fmE.
打印间隔的值说明了将水头资料写入表列文件(.LST)之间要进行的迭代次数。 :ld~9
4.4.2 强隐式求解程序包(SIP) { 'b;lA]0
强隐式求解法(SIP)是一种通过迭代来求解大型联立方程组的方法。其优点是非常稳定且常常可以收敛于某个解,但通常很慢。它不如PCG法快速,但占用的计算机内存要少。由于每个方程最多可涉及7个未知水头以及在整个网格中某个方程与其下一个方程的未知量是不同的,所以整个网格上的方程必须在每个时间步上都同时求解。本程序包的说明见MODFLOW 手册第十二章。 UtQj<18<
4.4.3 分层连续超松弛程序包(SOR) <)7aNW.
分层连续超松弛法(Slice-Successive Overrelation)是一种用迭代求解大型线性方程组的方法。在SOR中通过将有限差分网格分成垂向的层,并将节点方程重新组成分开的方程组,每组对应于一层。在每次迭代中将这些方程按顺序求解,就得到各层的一组新的水头估计值。当处理各层方程时,它们首先写成两次相连迭代间的计算水头的变化的形式。然后,通过将相邻层的项当作已知量,用高斯消去法直接求解相应层的方程组。接着,计算出来的水头变化值乘上一个加速因子T,T一般在1-2之间;这个结果就作为该层在此次迭代中的水头变化值。按顺序对每个层重复这个过程,直到三维阵列中的所有层都计算完成,就完成了一次整个模拟区域的迭代。然后再重复这个求解序列,直到相邻两次迭代之间的水头计算值的差值小于网格中所有的节点的判别标准为止。SOR程序包详见MODFLOW 手册第十三章。 qTJhYxm
参数说明 (&}[2pb!
最大迭代次数(Maximum Number of Iterations):默认值=50 C4+DZ<pE
本参数指定迭代次数的上限。这个值只有在迭代结束前不能达到收敛解的时候才用到。对大多数问题而言,50次就足够了。然而,如果达到了最大迭代次数,而不能得到较好的质量平衡差异,就应加大该值。 o/pw=R/):
加速因子(Acceleration Factor):默认值=1 PR8nJts W5
加速因子用于控制两次迭代间的水头变化的大小。它必须是正数。该值如大于1将会加大迭代间的水头的变化;可能可以较快得到解但也易使解失稳。小于1的值会使水头变化减少,因而增加了迭代次数。 Xf
u0d1b
收敛的水头变化判别标准(Head Change Criterion for Convergence):默认值=0.01 <KMCNCU\+
在每次迭代完成后,求解器就会检查每个单元的解,找出最大的变化量。如果这个值小于收敛标准(以英尺或米为单位)。就认为解已收敛,求解器就停止运行,否则开始新的迭代过程。通常对大多数问题来说,解的精度达到0.01ft或m就足够了,除非模拟区域的水头的最大差值小于1(ft或m)。如果迭代的次数还没有达到上面所设置的最大值,却没有得到较好的质量平衡结果,可以将这个值减少一个数量级。 *b{IWOSe^
打印间隔(Printout Interval):默认值=10 ] Q5:JV
打印间隔的值说明了将水头资料写入表列文件(.LST)之间要进行的迭代次数。 .psb#4
4.4.4 WHS的Visual MODFLOW求解器(WHS) ,`geOJn'
WHS使用了一种双共轭斜量稳定(Bi-CGSTAB)加速程序,它由用于预处理地下水流动偏微分方程的Stone不完全分解作为补充。这个求解器如同其它所有的迭代算法的求解器一样通过对近似解的迭代来求得大型偏微分方程组的解。由于地下水流动的矩阵方程最初是“病态的”,因而要求得一个有效解,就必须对这些矩阵进行有效的预处理。WHS中有两种等级的因式分解。程序在一级因式分解时收敛所需的迭代次数要少一些。而随着因式分解等级升高,运行程序所需内存变大。而且,在一级因式分解时,每次迭代的工作量也加大了,因而总计算时间未必比采用0级因式分解少。 s%)f<3=a
参数说明 U,g8:M
xHK
WHS求解器在每个时间步上用双层迭代法求解。 H4g8
1V=
外层迭代用于向解的方向改变因式分解参数矩阵。每进行一次外层迭代,因式分解参数矩阵中的流动系统的水文地质参数(即水力传导系数、饱水厚度、储水系数)就作了一次更新。不同等级的因式分解允许这些矩阵进行不同的初始化,以增加求解效率和模型稳定性。内层迭代用于迭代求解外层迭代建立的矩阵。 1 Pk+zBJ$
最大外层迭代次数(Maximum Number of Outer (non-linear) Iterations)默认值=100 mnm
ZO}
本参数指定外层迭代次数的上限。这个值只有在迭代结束前不能达到收敛解的时候才用到。对大多数问题而言,50次就足够了。然而,如果达到了最大迭代次数,而不能得到较好的质量平衡差异,就应加大该值。 ]Lv3XMa
最大内层迭代次数(Maximum Number of Inner Iterations)默认值=25 o[Ffa#sE
本参数指定内层迭代次数的上限。这个值只有在迭代结束前不能达到当前外层迭代的矩阵条件下的收敛解的时候才用到。对大多数问题而言,500次就足够了。然而,如果所有的内层迭代都达到了最大迭代次数,而不能得到较好的质量平衡差异,就应加大该值。 |A&;m}(Mt
收敛的水头变化判别标准(Head Change Criterion for Convergence):默认值=0.05 Oe5rRQ$O
在所有的外层迭代完成后,求解器就会检查每个单元的解,找出最大的变化量。如果这个值小于收敛标准(以英尺或米为单位)。就认为解已收敛,求解器就停止运行,否则开始新的外层迭代过程。通常对大多数问题来说,解的精度达到0.01ft或m就足够了,除非模拟区域的水头的最大差值小于1(ft或m)。如果内外层迭代的次数没有达到上面所设置的最大值,却没有得到较好的质量平衡结果,可以将这个值减少一个数量级。 $d<NN2
收敛的残差判别标准(Residual Criterion for Convergence):默认值=0.005 }3
xkA
水头变化的判别标准用于判断整个求解器的收敛,而残差判别标准则用于判断求解器内层迭代的收敛情况。如果相邻两次迭代的变化小于这里指定的值(单位:长度3/时间),那么就进行下一个外层迭代。 'f( CN3.!
残差标准对大多数问题来说0.001是比较合适的。然而,如果你注意到所有的外层迭代都只运行了一些内层迭代,而却没有得到较好的质量平衡结果,则可以将该值减少一个或几个数量级。 X1#Ar)
外层迭代阻尼因子(Dampening Factor for the Outer Iterations):默认值=0.5 <>HtXn/
这个因子允许用户在每两次相邻的迭代中减少(阻减)计算的水头变化。对大多数“井型”或物理意义上的地下水流动问题,阻尼因子选择1就可以了。这个参数可以用于非收敛(摆动或发散)的问题中,使其稳定下来,可以得到一个解。这是将该因子减少到0-1之间(很少小于0.6)来运行的。该因子类似于其它求解器中的加速因子。 9c@M(U@Yh
相对残差判别标准(Relative Residual Criterion):默认值=0 ng}C$d . I
本参数是另一种检查内层迭代的方法。该方法是将最近一次内层迭代的残差与最开始时的内层迭代的残差相比较。一旦最近一次内层迭代的残差小于最开始时的内层迭代的残差乘上相对残差判别标准,当前外层迭代就完成了,开始进行新的一次外层迭代。例如: +?J N_aR
如果:最近一次内层迭代的残差 )Zq'r L<
< 最开始时的内层迭代的残差 ´ 相对残差判别标准 A@V$~&JCL5
因式分解等级(Factorization Level):默认值=1 M}8P _<,
0级(Level 0)需要 1级(Level 1)需要 |] 7c&`
较多的外层迭代次数 较少的外层迭代次数 M9iX_4
较少的内存 较多的内存 #,#`<h!
4.5补给(Recharge)菜单 w6BBu0,KC
在Visual MODFLOW 中补给仅限于最上一层。然而,补给允许通过某列的无效单元(默认条件)。常水头边界总是截断补给,防止更深的渗透。因此,当你从运行屏幕上选择了Recharge时,就会有以下的对话框让你选择要用哪种补给。 Rqe.=+Qs
4.6层(Layers)菜单 nxZz{&
当你在运行屏幕中选择Layers时,就会有以下的对话框让你选择每个层的类型。 C19N0=
Type 0 Confined 在整个模拟期间,该层的导水系数和储水系数保持不变。 A8-[EBkK
Type 1 Unconfined 层的导水系数是变化的,它从饱水层厚度和水力传导系数中计算出来。储水系数保持不变;本类型只适用于第1层。 6KddHyFz
Type 2 Confined/Uncofined 层的导水系数和储水系数保持不变。储水系数可以有承压和非承压两种值。 Ci`o;KVj
Type 1 Unconfined/Confined 层的导水系数是变化的,它从饱水层厚度和水力传导系数中计算出来。储水系数可以有承压和非承压两种值。如果是含水层稀释,则限制从上面来的垂向渗透。 f@i#Znkf*?
4.7块中流(BCF)菜单 Ark]>4x>
在块中流(Block-Centered-Flow BCF)选项下,你可以指定单元二次湿润参数和层各向异性比率。 8T1`9ITl:
4.7.1二次湿润 &%2^B[{
最初的USGS的MODFLOW 是不允许在一次模拟中变干的单元再次饱水;它们简单地处理为无效单元。然而,USGS后来扩展了BCF包(BCF2),使得用户可以让干的单元二次湿润。其缺点是数值解的稳定性会下降。关于单元二次湿润和BCF包更详细的说明可见MODFLOW 程序包参考手册。 |Y3w6 !$
BCF菜单下的Rewetting选项可以让用户激活BCF2程序包中的湿润选项并指定二次湿润参数。通过单击“Activiate cell wetting”框可激活单元二次湿润选项。 |=0vgwd"S
单元湿润控制(Cell Wetting Control) 78l);/E{v
单元的湿润的控制有两种:通过单元底部的水头直接控制或通过四个相邻的水平单元加上其底下的单元的水头控制。第一种方法更稳定一些,并且在相邻单元难以指示是否一个单元应该湿润的时候就特别有用(例如在很陡的垂向水力梯度的地区或底板高程变化的地区)。在某些情况下,会导致干的单元处在湿的单元下。 yCQvo(V[F
湿/干阈值(Wet/Dry threshold) HV a9b;
对一个要湿润的单元来说,相邻单元的水头必须大于当前单元底部高程与湿/干阈值之和。例如,有一个底部高程为10的干单元,如果湿/干阈值设置为2,单元湿润由下面的单元控制,那么,该单元只有在下面的单元的水头达到12后才能被二次湿润。 Yq ]sPE92
总是在使用,以确定单元在何时湿润。 D;en!.[Z
单元湿润间隔(Cell Wetting Interval) m.D8@[y
单元湿润间隔表明了MODFLOW 多久就测试是否要湿润单元。例如,如果单元湿润间隔为2,则每隔一次迭代就看看是否要进行单元湿润检查。 x?S86,RW
单元湿润计算(Cell Wetting Calculations) 5*44QV
单元湿润时,新水头可用以下两个方程之一来计算: vl:~&I&y;R
水头= Z底部+湿润因子(相邻水头—Z底部) 9]eG|LFD
水头= Z底部+湿润因子(湿/干阈值) m)A:w.o
其中Z底部是当前单元的底部高程。 ?lC>E[
通常,第一个方程认为更合理一些,因为单元的新水头是随相邻水头而变化的,而正是后者使其湿润的。然而在MODFLOW 在迭代时高估水头变化的情况下,该方程有可能会不收敛。此时可用方程2来获得一个较稳定的解。 gTj,I=3$?e
解的收敛 =@U5/J
单元二次湿润常产生不收敛或不稳定的解,这可从单元不断地湿润排干看出。如果发生这种情况,我们推荐你试试以下方法: OBWb0t5H?
a. 将你已知肯定不会湿润的单元设置为无效单元。 D!.c??
b. 加大湿/干阈值。这会使单元的湿润更加困难,因而让MODFLOW 从不停地重复开关单元的过程中停止下来。然而,这可能会降低解的精度,因为应该湿润的单元会一直处在干的状态。 Y(UK:LZ'
c. 改变湿润因子。这会增加或减少湿润单元的水头。 ?t'V5$k\
d. 改变控制单元湿润的相邻单元。 \c2x
udU
e. 改变计算湿润单元新水头的方程。 cZVx4y%kz
f. 试用SIP或PCG求解器,并修改求解器参数。 \,13mB6
g. 对稳定态解,用符合最终条件的较好的初始猜测值开始求解。这将很好地指示哪个单元会湿润和排干,因而将单元在干湿之间的变换降至最小。 i^DMnvV.
4.7.2各向异性 ,C,nNaW
水平各向异性是指列向的导水系数或水力传导系数与行向的分量的比值。从BCF选择该项时,会出现如下对话框。 U'=8:&
各向异性系数可逐层设定或保留水力传导系数数据库(Hydraulic Conductivity Database)中的值。如果你选择逐层设置,则必须为每层输入一个Ty/Tx系数。选择逐层设定并不会冲去Ky的原始值,但将用Ty/Tx系数重新对它们进行计算。 h$8h@2%
4.8输出控制(OC)菜单 3t-STk?
当你从运行屏幕上选择OC(输出控制)菜单项时,会出现一个包括Output Control和List File Opts选项的下拉菜单。 JCcYFtW
4.8.1输出控制(Output Control) "^&H9