复杂函数边界控制下的潜水非稳定流模型及解的应用

吴　丹，陶月赞，林　飞

（合肥工业大学 土木与水利工程学院，安徽 合肥　230009）

1　研究背景

将完全切割潜水含水层的河渠水位作为第一类边界，边界一侧潜水非稳定渗流问题，是地下水渗流力学中的最经典问题之一［1-3］。该问题也属带自由面渗流问题，所以一直是地下水渗流力学领域中的研究难点与热点之一［4］。上述问题的解，可为灌溉渠系设计提供重要的计算依据［1］：通过调控渠系水位将影响区地下水水位控制在合理区间，是现代许多灌区的基本要求［5-6］。

就边界条件而言，河渠水位变化过程不仅受河渠自身水力特征和区域水文气象条件控制，而且受到明显的人工调控影响，具强烈的随机性［7-8］；对这种复杂多变的时变函数所构成的边界条件，随着对其概化与处理的方法不同，模型的求解方法以及所获得的解也有所差别。

这类问题的经典模型——J.G.Ferris模型（1950年），假定河渠水位变化瞬时完成之后长期保持稳定［1，3，9-10］，即下文模型（Ⅰ）中的边界条件H（t）为常数H0；对河渠水位复杂多变条件下的模型，利用常用积分变换方法求解时，按对边界处理方式与求解方法不同，现行方法主要可分为三大类：

方法一。H（t）是如时间线性函数等相对简单的函数，可通过常用的积分变换方法对模型直接求解［1，3，9-10］；该方法对一些特定情形（如经典模型所假设的情形）有较好的适用性，但不具普遍应用意义。

方法二。H（t）采用水位或水位变速“分段等值”方法离散化处理，依据经典模型的解，利用叠加原理构建问题的解［1］。这类方法就求解过程而言，方法比较简便；但当H（t）比较复杂时，划分的时段数较多，计算偏繁琐。

方法三。H（t）形式复杂，包括一些初始条件较复杂问题，在求解过程中或解的表达式中，需借助不常用的特殊函数或是构建专用的特殊函数［4-6，11-14］，这为应用带来不便。综合上述方法优点，针对这类问题的数学模型特点，在不依赖边界函数的变换过程的条件下，利用Fourier变换并依据卷积定义和卷积的微分性质，给出模型的理论解；对实际河渠水位过程，进行Lagrange线性插值离散，将插值函数代入理论解，给出关于问题的实际解；这一求解过程比较简单，所获得的问题实际解，不仅形式相对简单，而且是由常用函数表达。

2　经典模型及基本解

图1　河渠附近潜水渗流场

河渠及附近地段的水文地质条件（图1），可概括如下：

①具水平隔水底板的潜水含水层均质各向同性、平面无限延展；

②河渠水位变动过程，为时变函数H（t）；

③初始潜水位h（x，0）与河渠水位水平；

④河渠顺直且在剖面上基本完整切割含水层，潜水水流可视为一维流。该水文地质概念模型所对应的数学模型，可写成模型（Ⅰ）：

式中：μ为含水层给水度；h为潜水水位，m；K为渗透系数，m/d；H（t）为第一类边界（河渠水位）随时间变化的函数，m。

当h（x，t）-h（x，0）≤0.1hm（hm为计算期间的潜水流平均厚度，实际中的潜水问题多能符合该条件）时，可利用 Boussinesq方程的第一线性化方法；同时，令u（x，t）=h（x，t）-h（x，0），a=khm/μ，a为潜水含水层的导压系数，m2/d；可将模型（Ⅰ）转化为模型（Ⅱ）。

当H（t）=H0，H0为已知常数时，即为J.G.Ferris模型；此时模型的解为［1-3］：

式中，erfc（z）是余误差函数，z=x/2（at）1/2。式（9）也是这类问题的基本解。

3　模型理论解

对于模型（Ⅱ），在不依赖H（t）具体变换过程的条件下，建立这类模型的理论解。

由边界条件与变量的变化范围，根据Fourier变换的性质与特点，对模型（Ⅰ）求关于x的Fourier正弦变换。为u关于x的Fourier变换，ω、F和F-1分别为变换算子、变换与逆变换算符。有：

所以，由式（5），有：

上式通解为：

对式（10），求正弦逆变换，注意正弦与余弦逆变换关系，并适时交换积分次序：

式（12）是在H（t）还未明确条件下模型（Ⅱ）的解；也即对任意H（t），式（12）是模型的理论解。具体到实际问题应用，还需要将已知的H（t）代入，再进一步展开，才能获得实际问题的解。

为获得更为便捷的求解方法，并结合本领域的应用习惯，以下依据积分变换性质及卷积定理，给出以概率密度函数表达的解。

据卷积定义，由式（12），有：

式中*为卷积算符。

由卷积的微分性质：

式（15）实质上也是模型（Ⅱ）的另外一种形式的理论解。

上述求解过程中，边界条件H（t）一直未进行形式上的变换，也即在H（t）属未知的情况下进行的；这一解法可为H（t）复杂条件下问题求解提供支持。

需要说明的是，上述方法求解过程中，虽然H（t）未进行形式上的变换，但H（t）实质上参与了变换与逆变换过程；所以H（t）必须满足Fourier变换要求，即：函数在任意区间满足Dirichlet条件，在无限区间上绝对可积；就河渠水位实际变化过程而言，通常都可满足上述要求的。

4　河渠水位过程离散及其解

4.1Lagrange线性插值实际中，复杂的水位变化过程H（t）显然难以给出统一的数学表达式；此条件下，不依赖H（t）的具体函数形式，只根据H（t）的实测过程，采用Lagrange线性插值方法对其进行离散化处理。

根据Lagrange线性插值原理，在时段Ti=ti-ti-1，有：

则有：

ε（t-ti-1）系Heaviside 函数［12］， 具有如下性质：当t＜ti-1、ε（t-ti-1）=0，当t≥ti-1、ε（t-ti-1）=1。

4.2线性插值函数对应的解将式（17）带入式（15），有：

注意ε（t）函数的性质及u（x，t）=h（x，t）-h（x，0），整理式（18），可得：

式（19）是在河渠边界控制下的半无限潜水含水层中，对水位变化过程H（t）采用Lagrange线性插值离散化的条件下，潜水非稳定渗流过程的解析解。

当然，也可将式（17）带入式（13），再通过分部积分等步骤进一步展开，以获得最终形式的解；但由式（15）而获得式（19），求解过程显然要简明得多。

5　解的应用

河渠潜水渗流模型研究的重要目的之一是以其为工具，利用潜水位变动数据计算模型参数a。

5.1模型参数计算潜水位变动速度φ（x，t）=∂h()x，t ∂t，由式（16）：

以下仅讨论H0、λ1的情形。

（1）当H0=0∩λ1≠0的配线法。由式（20）：

对距离边界为x的观测孔（x为确定值），z=x/2（at）1/2，首先建立不同a值的erf（z）～t理论曲线图族；由实测地下水水位，制出φ（x，t）～t的曲线；由式（19），当φ（x，t）～t所依赖的含水层参数a值，与erf（z）～t理论曲线图族中某条曲线a值相等时，则，两条曲线形态完全相同，仅相差一个常数λ1倍，也即两条曲线应该完全重合。

因此，根据地下水水位实测数据，建φ（x，t）～t曲线，将之与erf（z）～t理论曲线族进行适线，就可确定含水层的a值（图2）。

图2　适线法求a

（2）当H0≠0∩λ1=0。式（20）转为J.G.Ferris模型及其解，解的利用有许多文献，在此不作赘述。

5.2河渠与潜水之间交换量计算

5.2.1算式河渠与潜水之间的补排关系，由河渠水位与潜水水位关系所决定，这也表明，河渠附近潜水非稳定渗流模型，可以用来计算河渠与潜水之间的水量交换；仅讨论H0、λ1的情形。

河渠与潜水之间的交换量，用单位河渠长度上的交换量表达，即单宽流量Q（t）（m3·d-1·km-1）；t时刻的交换强度q（t），由Darcy定理：

式中：hm是潜水流厚度，a=khm/μ。

式（22）是河渠与潜水之间在t时刻的交换强度算式；q＞0，表示潜水接受河渠补给；q＜0，表示潜水排泄河渠。

5.2.2λ段的影响规律由式（19），去掉等式右端第2项，转化为J.G.Ferris模型的交换强度算式，对应的是河渠水位瞬间变动H0所获得的水量交换强度，记为qH；式（19）去掉等式右端第1项，对应的是H0=0条件下、仅由河渠水位变化段所获得的交换强度，记为qλ。

λ段的影响，可用qλ/qH来反映：

假设水位在变动段获得的水位升幅也为H0（即λt等于H0），由式（23），对交换量的贡献，λt是H0的2倍；也即水位变动段的水位升幅虽然也为H0，但其对交换量的贡献，是不考虑变化过程影响而计算出效果的2倍；这是水位在t时段内缓慢上升至λt的过程中，被不断抬升的水位持续增强补给作用的过程累积效应所至。

6　实例研究

安徽淮北平原中部的蒙城县境内茨淮新河灌区，以粉细砂为主的潜水含水层分布广泛，厚度8 m左右，底部一般发育有不完全连续黏性土层；由于潜水位埋深浅，为2.5～3.0m；区内农田灌溉渠系比较完善，干渠基本深切至隔水底板、干渠的渠间距为2 km，干渠渠首多有节制闸控制；一眼国家级地下水位自记观测井，距离干渠直线距离为60 m处，观测井附近地面标高31.02 m。

2014年7月中下旬，干渠长时间处于非运行状态，水位基本保持不变；25日至30日，一直处于无雨期，也无灌溉回渗影响，计算中不考虑潜水垂向交换作用；例中，时间段较长，数据量也较大，观测孔水位按12h摘录，如表1。

表1　潜水水位动态数据(2014.7.25—2014.7.30)与“配线法”计算过程

该时段，H0=0，不符合用拐点法求参数的条件，但时间段较长，数据量足够大，可依据这段时间潜水变动速度随时间变化曲线，采用配线法求算；由上表数据，计算出的φ（x，t）～t与erf（z）～t理论曲线图中的曲线族进行适线（如图2），适线结果a值为890 m2/d；这与该地区相关文献［10］的计算结果（854 m2/d）基本一致。

7　结语

对水位变化过程复杂的河渠边界附近潜水非稳定渗流问题，在对水位过程离散化的基础上，利用Fourier变换给出问题的解，并对解的应用进行了相应的探讨，所得结论如下：

（1）河渠水位变化过程往往比较复杂且难以给出具体函数表达式，利用Fourier变换，在不依赖边界函数变换的条件下，获得理论模型通用解；再采用Lagrange线性插值方法，对实际水位过程离散化处理，将插值函数代入理论模型通用解，可获得由比较常用函数组成的解，且形式较简单。

（2）在利用Fourier变换求解过程中，可充分利用如卷积定理和卷积的微分性质等变换性质，由此获得如文中式（11）的理论通用解，可为实际问题求解提供更简便的途径。

（3）河渠水位在其变动过程中，被不断抬升的河渠水位（以上升为例），使河渠对地下水的补给作用呈持续增大并伴生有过程累积效应；在计算河渠与潜水之间水量交换时，采用水位稳定或分段稳定的计算方法，将因忽略累积效应而导致计算偏差，偏差幅度为2λ1t/H0～λ1t/（H0+λ1t）。

（4）依据模型解，利用地下水水位动态监测数据，可计算含水层参数a；在H0对地下水水位作用不明显时段，利用地下水位变动速度随时间变化过程与理论曲线进行配线，实例表明，方法可行。

值得指出的是，文中的模型及其解，仅适用模型（Ⅰ）所对应的水文地质条件，即一条直线河渠控制的半无限域潜水渗流问题；随着水文地质条件变化而应采用相应的渗流问题模型，如有多条平行河渠（如灌区内的支渠）相互干扰时，应采用河间地块潜水渗流问题模型。