第5章空间离散化和参数赋值

___Y1

已于 2024-04-09 12:29:41 修改

阅读量1.6k

点赞数 29

分类专栏：地下水数值模拟文章标签：算法

于 2023-12-08 10:20:47 首次发布

本文链接：https://blog.csdn.net/weixin_44259522/article/details/134871720

版权

我可以被限制在一个核桃壳里，却算自己是无限空间的国王。

——《哈姆雷特》，第二幕第二场

在前面的章节中，我们介绍了将问题域离散化为在有限差分（FD）网格中定义单元/块或在有限元（FE）网格中定义元素的节点的概念。本章重点讨论了空间离散化和将初始参数值分配给网格/网格的过程。时域离散化在第7.6节中讨论。

5.1 空间离散化

设计有限差分（FD）网格或有限元（FE）网格是模型设计中最重要且最困难的步骤之一。网格/网格中的节点数量影响解决方案的准确性、求解模型所需的计算时间以及将生成的输出量。通常使用图形用户界面（GUI）（第3.6节）或FE网格生成器来生成网格/网格，但建模者必须输入节点间距并指定FE网格的元素类型。一些GUI允许在设计网格/网格之后添加或减少节点；然而，重新设计节点网络需要重新分配参数值，并通常需要调整边界条件。因此，花些时间思考如何最好地设计网格/网格以满足建模目标是明智的。此外，相同的网格/网格将用于建模项目的所有阶段，可能包括校准、粒子跟踪和预测模拟。因此，提前计划确保网格/网格将满足建模项目中进行的所有类型模拟的目标至关重要。

在数值模型中，使用节点和FD单元/块或有限元（图3.4，3.5，3.10和4.11）离散化问题域。节点网络形成数值模型的计算框架，节点间距决定输出的分辨率，因为在节点计算头，并且必须在节点之间插值。虽然有时将术语网格和网格互换使用，但我们在指FD节点网络时使用网格，而在指FE网络时使用网格。概念模型和边界的位置确定了水平网格/网格的总体尺寸。面向二维（2D）的模型只有一个由水平节点间距定义的层。三维（3D）模型有两个或更多层，但通常模型的所有层中水平节点网络是相同的。然而，一些代码允许水平节点网络按层变化。

总体而言，FE网格在离散化问题域方面提供了更大的灵活性，因为一个网格可以包含三角形和/或四边形（即四边形，不限于矩形）元素，而标准FD网格则要求矩形网格。

网格和网格都可以是结构化的或非结构化的。在结构化网格/网格中，3D空间中的一个节点最多与六个相邻节点相连，但在非结构化网格/网格中，连接可能少于或多于六个（图5.1(b)），并且连接的数量可能因节点而异。FD代码使用控制体积有限差分（CVFD）方法来适应非结构化网格（Panday等，2013年）。早期使用此方法形成非结构化网格的FD方法被称为集成有限差分（IFD）。数值模型首次应用于地下水问题是使用非结构化FD网格的（Tyson和Weber，1964；也见Narasimhan和Witherspoon，1976）。MODFLOW-USG，其中USG代表非结构化网格（Panday等，2013），是一个允许使用非结构化网格的CVFD代码的最新示例。为了明晰起见，我们区分使用CVFD的FD方法，该方法通过非结构化网格实现更灵活的网格设计，以及传统（标准）FD方法，该方法需要结构化矩形网格。

图5.1 节点之间的连接。（a）在3D结构化网格或网格中，一个节点最多与六个相邻节点连接；在2D中，最多有四个连接（http://doi.

ieeecomputersociety.org/cms/Computer.org/dl/mags/cs/2012/03/figures/mcs20120300483.gif）。（b）水平2D非结构化FD网格，其中中央节点连接到其他六个节点；在3D中，此网格将具有八个连接（修改自Tyson和Weber，1964年。此材料仅可用于个人用途。任何其他用途均需获得美国土木工程师协会的事先许可）。

垂直方向的离散化取决于模型层的厚度，其中每个层通常代表一个水文地质单元（第5.3节）。因此，在大多数数值模型中，垂直方向的离散化尺度与水平节点网络中的节点间距不同。我们推迟对垂直离散化的讨论到第5.3节。以下是关于2D水平网格/网格的一些一般性评论。水平节点间距的选择在第5.2节中讨论。

5.1.1 网格/网格的定向

理想情况下，网格/网格的方向应使x或y坐标轴与地下水流的主导方向平行。此外，坐标轴应与水力传导率（K）张量的主导方向共线（Box 3.1）。

主导流向可能由结构性特征控制，如断层和节理组（图5.2），水文地质单元的走向，或层理（图5.3）。当没有一个显著的流向时，网格/网格通常与GIS或GUI中使用的投影对齐，这些软件用于组装代码的输入文件。

在罕见的情况下，各向异性的方向必须被模拟为非均匀的，即K张量的主导方向在空间中改变方向，因此无法将模型的全局坐标轴与K张量的主导方向共线（Box 3.1中的图B3.1.1(b)）。例如，在模拟倾斜的水文地质单元（第5.3节）时，模型会出现表示非均匀各向异性的问题。在这种情况下，所提出的方程必须包括K张量的非对角分量（Box 3.1中的Eqn（B3.1.4）），并且最好使用FE代码解决该方程的形式。FD代码理论上可以修改以处理非均匀各向异性，对于K的水平分量来说这更容易实现（Anderman等，2002）。理论上可以校正非均匀垂直各向异性，但在实践中编程难度较大（Hoaglund和Pollard，2003；Yager等，2009）。

与FD方法（包括CVFD）相比，有限元方法更容易适应各向异性方向的变化，因为在组装FE解的系数矩阵时，每个元素都被单独考虑（第3.5节），从而元素内的局部坐标轴可以与水力传导率张量的局部方向对齐（参见Box 3.1中的图B3.1.1(b)）。许多FE代码允许各向异性方向的非均匀性，但并非所有FE代码都支持。在FE代码FEFLOW（Diersch，2014）中，各向异性可以与全局坐标轴平行，与倾斜的地层和地质结构平行（例如图5.3中的情况），或与元素内的局部坐标系平行，这是最一般的情况。FE代码SUTRA（Voss和Provost，2002）也允许非均匀各向异性的最一般情况。

5.1.2 有限差分网格

标准的有限差分（FD）网格是矩形和结构化的（图5.4(a)和(b)）。非结构化的FD网格在网格设计上提供了更大的灵活性，并提供了FE网格的许多灵活性。

图5.3 FE网格中的坐标轴被定位以与K张量的主要分量一致，如左侧的地层显示。未显示与挡土墙附近的详细网格（修改自Townley和Wilson，1980）。

图5.3 FE网格中的坐标轴与K张量的主要分量一致，如左侧的地层所示。未显示靠近挡土墙的详细网格（修改自Townley和Wilson，1980）。

5.1.2.1 结构化网格

结构化FD网格具有矩形单元/块（图4.11(b)和3.4(c)），可以是块居中的（图5.5(a)和4.11(b)），也可以是点居中的（图5.5(b)和4.11(c)）。在块居中的网格中，网格线形成单元的边缘，节点位于单元的中心；而在点居中的网格中，网格线穿过节点。在任一情况下，节点周围的矩形区域是FD单元。结构化FD网格的整体形状是矩形的，但定义问题域的地下系统通常不是矩形的（图5.2）；问题域的边界在表示模型域范围的较大矩形区域内表示（图5.2）。在代码计算头值的模型域内，节点是活动节点，而落在问题域边界之外的节点是非活动节点。非活动节点不参与解决方案，代码不计算这些节点的头值（但它们仍然占用代码使用的数组中的存储空间）。

在结构化FD网格中，节点的位置使用两维的行和列索引约定进行标识，并使用三维的行、列和层约定（图3.4）。通过在水平平面上指定Δx和Δy的值数组，生成2D结构化FD网格。用户应注意设置网格的代码特定约定。例如，在MODFLOW中，变量DELR（Δr）表示列之间的水平节点间距（等效于Dx），而DELC（Δc）表示行之间的水平间距（等效于Dy）。节点使用i索引表示行，使用j表示列（Section 3.5.1），因此DELR数组对应网格中的每一列具有一个值（Drj），DELC对每一行（Δc_i）具有一个值（图5.5）。在3D模型中，Δz值的数组可以直接输入，或者可以从顶部和底部层的高程计算Dz。

5.1.2.1 结构化网格

Figure 5.4 结构化和非结构化网格设计（Panday等人，2013）。

在规则间隔的网格中（图5.2a），水平平面上的节点间距是均匀的（Δx=Δy = 常数），并且节点位于单元的中心。在不规则网格中，节点在一个或两个水平方向上的间距是变化的；在块居中的网格中，节点仍位于单元的中心，但在点居中的网格中节点不在中心（图5.5）。

规则间隔的水平网格在计算上更精确和操作上更高效；许多后处理例程设计用于规则间隔的水平节点网络。然而，在仅在局部需要细分以表示地表水特征、抽水井、参数变化、陡峭水力梯度和其他水文学上重要的特征时，规则网格可能不切实际。当只在模型域的一部分需要细分时，设计不规则网格，或者在较大尺度网格中嵌入一个局部尺度网格（图5.4(h)、(i); Section 4.4, Fig. 4.20）。从概念上讲，嵌入的局部尺度网格通常更受欢迎，但可能需要更多的时间来设置和运行。这里我们讨论标准FD结构化网格中使用的不规则网格。

Figure 5.5 使用MODFLOW约定显示的不规则网格，其中i = 行，j = 列。列之间的间距（行间的间距）为Drj，行之间的间距（列间的间距）为Dci：(a) 块居中的网格；(b) 点居中的网格（修改自McDonald和Harbaugh，1988）。

不规则网格在整个网格中必须具有相同数量的行和列（图5.6），这几乎总会导致在焦点区域之外进行细化。在构建不规则网格时，两个水平维度上的节点间距应该逐渐增加或减小。一个常见的准则是通过增加节点间距来扩展网格，使其不超过先前节点间距的1.5倍。

例如，如果两个节点之间的距离为1米，则下两个节点之间的距离不应超过1.5米，下一个节点间距不应超过2.25米，依此类推。对于少数几行或列，可以使用2的倍数，但最好不要使用超过1.5倍的因子扩展整个网格。不规则网格的准则是因为对于不规则节点间距，二阶导数的FD表达式的误差较大（Eqn (3.25)）。不规则网格的截断误差在一阶上是正确的，而对于规则网格，FD表达式在二阶上是正确的。通过泰勒级数展开可以展示不规则网格截断误差的差异（Remson等人，1971，pp. 65e67）。对于不同的网格，数值误差较大的一种启发式解释见Box 5.1。当x-和y-方向都存在可变节点间距时，数值误差会叠加，但只要在两个水平方向都遵循1.5规则，就是可以接受的。

尽管规则网格更受欢迎，因为它对数值精度更高，并且输入数据和处理输出更为方便，但如果在整个问题域中使用细分节点间距，则需要更多的节点。更多的节点需要更大的输入和输出数组，更多的计算机存储和数据管理。然而，随着今天计算机的计算能力和存储能力的增加，对于许多地下水问题，使用均匀节点间距的全局细化网格是可能的。

Figure 5.6 不规则FD网格的设计，旨在在新泽西州（美国）一个超级基金现场附近提供细致的节点间距（Lewis-Brown等人，2005）。

**框5.1** 不规则FD网格固有的数值误差

在本框中，我们提供了一个启发性解释，说明相对于规则FD网格（第5.1节）而言，基于不规则FD网格的FD逼近存在更大误差。考虑图B5.1.1中的一维块中心网格。在节点i周围，二阶导数的FD表达式是在点 $i+1/2$ 和 $i-1/2$ 处计算的一阶导数的差异：

当 $\Delta x$ 是常数时（图B5.1.1(a)），点 $i+1/2$ 和 $i-1/2$ 与块的边缘重合，使得节点点(i)位于 $i+1/2$ 和 $i-1/2$ 之间。然而，在不规则网格中（图B5.1.1(b)）， $i+1/2$ 不与包含节点i的块的边缘重合，因此节点点(i)不再位于 $i-1/2$ 和 $i+1/2$ 之间。

FD解决方案计算了在 $i+1/2$ 和 $i-1/2$ 之间的位置的水头。

当节点不位于 $i+1/2$ 和 $i-1/2$ 之间时，计算水头的位置与节点不重合，从而引入了一些误差到解决方案中。

**图B5.1.1** （a）规则FD网格；（b）不规则FD网格。两者均以一维形式展示。节点用填充的圆表示；位于节点中间的位置由红色的X标出。位置 $i + 1/2$ 位于节点 $i$ 和 $i + 1$ 之间的中点；位置 $i -1/2$ 位于节点 $i -1$ 和 $i$ 之间的中点； $\Delta x_{i+1/2}$ 是节点i和i + 1之间的距离； $\Delta x_{i-1/2}$ 是节点 $i -1$ 和 $i$ 之间的距离； $\Delta x_i$ 是围绕节点i的单元的宽度。在规则网格中， $\Delta x_{i-1/2} = \Delta x_{i+1/2 }= \Delta x$ 。出于说明目的，不规则网格扩展了四倍，而不是推荐的1.5倍（第5.1节）。

5.1.2.2 无结构网格

在CVFD模型中，非结构化FD网格（第3.5节）通过使用嵌套网格（图5.4（h），（i）和5.7）提供了与FE网格相似的节点设计灵活性，这是结构化不规则网格的理想替代方案。嵌套网格的优势在于在问题域焦点区域之外的区域不需要细节点间距。这个概念类似于望远镜网格细化（TMR）和局部网格细化（LGR）（第4.4节），其中在焦点区域的详细离散化在局部尺度模型中实现，该模型的边界条件从更大尺度的模型中提取。然而，TMR/LGR方法需要链接多个数值模型和多组边界条件来表示嵌套问题域。相反，在非结构化嵌套网格中，小的节点间距嵌入到具有较粗节点间距的较大尺度网格中，只有一个问题域，一个边界条件集和一个解。

MODFLOW-USG（Panday等，2013）提供了嵌套网格（图5.4（h），（i）和5.7）以及围绕地质或工程特征使用四叉树细化（图5.4（k）和（l））的详细网格分辨率。通过四叉树网格化（在三维中称为八叉树网格化），相对较粗的节点间距的区域网格通过将每个矩形单元细分为最多四个较小单元而进行多次本地细化，以形成围绕诸如抽水井、河流和渠道等特征的嵌套网格（图5.8）。

非结构化网格的灵活性带来了构建网格所需输入的额外复杂性。与FE网格一样，节点是按顺序编号的，而不是使用标准FD的高效i、j、k索引约定定位。每个单元也被编号（图5.7）。相邻节点被连接，用户必须为每个单元连接输入共享单元面积。鉴于组装非结构化网格的数据的复杂性，建议使用预处理器代码（Panday等，2013）。

基于笛卡尔坐标设计非结构化FD网格的一个重要准则是相邻单元中心连线应该垂直穿过共享面（图5.9），并且交点应该大致位于面的中央。对于圆柱网格，交点应与半径的对数平均值重合。遵循这一指导方针，称为控制体积有限差分（CVFD）要求，确保FD方程正确到二阶。严格来说，对于嵌套网格，该指南被违反（有关详细信息，请参见Panday等，2013），但可以容忍小的违规行为。图5.4（a）至（e）和（j）中的网格符合CVFD要求，而图中的其他网格则不符合。MODFLOW-USG包括一个名为Ghost Node Correction（GNC）Package的模块，旨在最小化CVFD要求违规引起的错误。GNC Package会自动放置一个幽灵节点，使连接幽灵节点与相邻单元中的活动节点的线在两个单元之间的面上呈直角（图5.10）。跨越面的通量基于幽灵节点处的水头计算，该水头通过插值来近似。幽灵节点还用于计算FD近似中地下水流方程（方程（3.27））的单元间导水率项的连接长度。幽灵节点本身不作为解中的活动节点包括在内；这些节点的水头是通过插值计算的。