基于流网的非均质含水层优势流路径识别方法技术

本发明专利技术公开了一种基于流网的非均质含水层优势流路径识别方法，该方法首先建立流管面积或体积与水流移动时间之间的定量关系；将已有的地下水流模拟程序用于非均质介质中流函数计算；通过计算流管面积或体积，高效、准确地识别时间最短的移动路径。与已有的粒子追踪法和图论技术相比，本发明专利技术提供的方法效率更高、准确性更好。

Path identification method of dominant flow in heterogeneous aquifer based on flow network

【技术实现步骤摘要】
基于流网的非均质含水层优势流路径识别方法
本专利技术属于水力学
，具体指代一种基于流网的非均质含水层优势流路径识别方法。
技术介绍
目前，地下水污染问题得到越来越多的关注；地下水溶质运移难点之一在于地下水流场中优势流的存在，即：占含水层小部分体积、导水能力强的优势通道，常常可以传输大量的水流和溶质，甚至显著改变含水层行为。在自然界中普遍存在的含水层非均质性，是导致了地下水流场非均一性和优势流的主要因素。基于含水层非均质性识别地下水流场中的优势流路径，是地下水运移模拟与水质监控的重要部分。目前优势流路径识别与分析的主要手段是传统的粒子追踪法(PT)，以及近年新出现的基于图论的最小阻力路径法(MHR)。传统PT法在流场中释放大量粒子，通过追踪这些粒子在对流和弥散作用下的运动，分析水流运动不均匀性。该方法可靠性好，但主要问题在于计算量很大，需要大量时间及存储空间，而且根据粒子当前位置计算下一时刻位置容易引入不断积累的误差。MHR法以类似于导航寻路的算法寻找整体阻力最小的路径，计算效率较高，但其不考虑水动力条件和边界条件，容易导致不符合水动力规律的搜寻结果，因此该方法的可靠性还无法令人满意。
技术实现思路
针对于上述现有技术的不足，本专利技术的目的在于提供一种基于流网的非均质含水层优势流路径识别方法，通过利用已有的工具，便捷、高效且可靠地识别非均质含水层中的优势流路径，为有效刻画自然条件下污染物监控提供有力工具。为达到上述目的，本专利技术采用的技术方案如下：本专利技术的一种基于流网的非均质含水层优势流路径识别方法，包括步骤如下：(1)根据所要研究的非均质含水层，确定渗透系数K的分布；(2)对上述渗透系数K取倒数，得阻力场R分布，即R＝1/K；(3)根据含水层条件，指定流函数Ψ的边界条件；(4)计算含水层内流函数Ψ的分布；(5)根据计算得到的流函数Ψ分布，采用等高线法计算一系列等差流线(即相邻流线的流函数之差相等)，流线包裹的空间即为流管；(6)计算各流管的面积或体积A；(7)对各流管的面积或体积A进行排序，A最小的流管对应最短的移动时间，即为优势流路径所在。进一步地，所述步骤(3)中流函数Ψ与步骤(2)中阻力场R的关系满足如下方程：其中，x为空间向量，该方程对应的边界条件需要指定如下：按照区域流场状况，上下游分别为流函数零坡度边界其他平行于主流方向的边或面为恒定流函数边界Ψ(x)＝C，C为常数。进一步地，所述步骤(4)中利用地下水流模拟软件(如MODFLOW、FEWATER等)在非均质含水层中求解上述流函数方程，得到流函数的空间分布。进一步地，所述步骤(5)利用等高线计算方法(如MATLAB中的contour函数等)计算若干等差流函数值的在空间上的分布，即等差流线。进一步地，所述步骤(6)中使用polyarea函数计算流管的面积或体积。进一步地，所述步骤(6)中计算的流管面积或体积，正比于流管内的渗流运动过程的时间t，即：其中，ds为沿着流管方向的微小距离，v为达西流速的大小，w(s)为流管的宽度(三维情况下为截面积)，|Ψi+1-Ψi-1|为构成流管的流函数之差，A为该流管的总面积(三维情况下为总体积)。本专利技术的有益效果：1、本专利技术采用地下水模拟软件直接计算流函数，再通过流函数得到流线，不涉及大量计算，精度和效率均很高。2、本专利技术既考虑了含水层非均质性，又能兼顾水动力条件，得出的优势流路径是流场中的一条流线，因此比MHR法的结果更符合物理规律、更可靠。3、本专利技术各步骤涉及的细节技术可以用已有的工具解决，实现简单、操作便捷。附图说明图1为本专利技术的执行步骤示意图。图2为流函数Ψ的边界条件示意图。图3为实施例1中含有低渗透墙体的含水层示意图。图4为实施例1中用本专利技术的SF方法得到的优势流路径结果示意图。图5为实施例1中PT、MHR及SF方法得到的三种优势流路径结果与PT法中的粒子锋面及最先抵达处对比示意图。图6为实施例2中含有低渗透墙体的含水层示意图。图7为实施例2中PT、MHR及SF方法得到的优势流路径结果与PT法中的粒子锋面及最先抵达处对比示意图。图8为实施例2中PT、MHR及SF方法得到的优势流路径结果与流场流线对比示意图。图9为实施例3中对数渗透系数(lnK)空间分布示意图。图10为实施例3中PT、MHR及SF方法得到的优势流路径结果与PT法中的粒子锋面及最先抵达处对比示意图。图11为实施例3中PT、MHR及SF方法得到的优势流路径结果与流场流线对比示意图。图12a为100个算例中观察到的第一种典型情况示意图。图12b为100个算例中观察到的第二种典型情况示意图。图12c为100个算例中观察到的第三种典型情况示意图。图12d为100个算例中观察到的第四种典型情况示意图。图13为100个算例中SF和MHR相对于PT路径的平均距离d的整体比较结果图。图14为100个算例三种方法的耗时比较结果图。具体实施方式为了便于本领域技术人员的理解，下面结合实施例与附图对本专利技术作进一步的说明，实施方式提及的内容并非对本专利技术的限定。参照图1所示，本专利技术的一种基于流网的非均质含水层优势流路径识别方法，包括步骤如下：(1)根据所要研究的非均质含水层，确定渗透系数K的分布；(2)对上述渗透系数K取倒数，得阻力场R分布，即R＝1/K；(3)根据含水层条件，指定流函数Ψ的边界条件；(4)计算含水层内流函数Ψ的分布；(5)根据计算得到的流函数Ψ分布，采用等高线法计算一系列等差流线(即相邻流线的流函数之差相等)，流线包裹的空间即为流管；(6)计算各流管的面积或体积A；(7)对各流管的A进行排序，A最小的流管对应最短的移动时间，即为优势流路径所在。参照图2所示，所述步骤(3)中流函数Ψ与步骤(2)中阻力场R的关系满足如下方程：其中，x为空间向量，该方程对应的边界条件需要指定如下：按照区域流场状况，上下游分别为流函数零坡度边界其他平行于主流方向的边或面为恒定流函数边界Ψ(x)＝C，C为常数。以下对粒子追踪法简称PT，基于图论的最小水力阻力法简称MHR，本专利技术的基于流函数的算法简称SF。实施例1：含有两个低渗透墙体的均质含水层研究对象为含有两个低渗透墙体的均质承压含水层，区域平均水流由左到右，主体与低渗墙体的渗透系数比为K/Kb＝103(如图3所示)，孔隙度为常数。采用本专利技术的SF方法，计算等差流线和各流管面积A(如图4所示，图中显示了等差流线的空间分布，即相邻流线的流函数之差相等；右侧展示了流管的A值，最小A值以空心圆标示，接近含水层对称中轴)，最小的A对应的流管即为优势流路径。所得的SF路径为一平直通过墙体空隙的流线，与常识相符。为了比较本专利技术与PT、MHR方法并进行进一步分析，图5将PT、MHR及SF方法得到的三种优势流路径结果与PT法中的粒子锋面及最先抵达处进行了对比。结果表明：SF路径(图5中虚线)与PT路径(图5中实线)几乎完全一致，均指向粒子最先抵达处；MHR路径(图5中点线)也为一通过墙体的平直线。实际上所有通过空隙的平直线对MHR算法来说都是优势路径，然而真正的最快路径只有PT路径那一条(当然，算得的PT路径本身也有一定的误差)。但结合水动力条件可知，除了对称轴y＝100和上下边界以外的流线都是弯曲的，MHR方法无法识别这些路本文档来自技高网...

【技术保护点】
1.一种基于流网的非均质含水层优势流路径识别方法，其特征在于，包括步骤如下：(1)根据所要研究的非均质含水层，确定渗透系数K的分布；(2)对上述渗透系数K取倒数，获得阻力场R分布，即R＝1/K；(3)根据含水层条件，指定流函数Ψ的边界条件；(4)计算含水层内流函数Ψ的分布；(5)根据计算得到的流函数Ψ分布，采用等高线法计算一系列等差流线，流线包裹的空间即为流管；(6)计算各流管的面积或体积A；(7)对各流管的面积或体积A进行排序，A最小的流管对应最短的移动时间，即为优势流路径所在。

【技术特征摘要】
1.一种基于流网的非均质含水层优势流路径识别方法，其特征在于，包括步骤如下：(1)根据所要研究的非均质含水层，确定渗透系数K的分布；(2)对上述渗透系数K取倒数，获得阻力场R分布，即R＝1/K；(3)根据含水层条件，指定流函数Ψ的边界条件；(4)计算含水层内流函数Ψ的分布；(5)根据计算得到的流函数Ψ分布，采用等高线法计算一系列等差流线，流线包裹的空间即为流管；(6)计算各流管的面积或体积A；(7)对各流管的面积或体积A进行排序，A最小的流管对应最短的移动时间，即为优势流路径所在。2.根据权利要求1所述的基于流网的非均质含水层优势流路径识别方法，其特征在于，所述步骤(3)中流函数Ψ与步骤(2)中阻力场R的关系满足如下方程：其中，x为空间向量，该方程对应的边界条件需要指定如下：按照区域流场状况，上下游分别为流函数零坡度边界其他平行于主流方向的边或面为恒定流函数边界Ψ(x)＝C，C...

【专利技术属性】
技术研发人员：南统超，吴吉春，王远坤，
申请(专利权)人：南京大学，
类型：发明
国别省市：江苏,32