在稳定流边界元方程（4.4.35）和非稳定流边界元方程（4.4.51）中，都含有源汇项部分
???Gp?(x,y)Tdxdy；
???Gp?(x,y)Tdxdy （4.4.52）
这部分积分的计算仍然需要采用离散化网格，如三角形网络，进行数值积分运算。因此，含有源汇项的地下水流问题，虽然渗透矩阵的建立只要进行边界离散化即可，但面源补给项的计算却难以避免使用模型区域内部的离散化网格。
井孔也是地下水流重要的源汇项，需要在有限元模型中加以处理。对于已知流量的井孔，设其坐标位置为为
2??-Qw/?rw,?w(x,y)????0,(xw,yw)，半径为rw，抽水流量为Qw，则这一井孔形成的补给强度可表示
d?rwd?rw （4.4.53）
其中：
d?(x?xw)2?(y?yw)2 （4.4.54）
当井孔趋于无穷小时，井孔附近源汇项的积分可计算为
???Gp?w(x,y)Tdxdy??lim[?rw?0rw0?2?0Gp(rp)Qwrdrd?] （4.4.55） 2?rwT式中：r的起点为井点；?的起算方向是从点p指向井点坐标的方向。函数GP重新表示为
Gp(rp)?ln(rp,w?rcos?) （4.4.56）
式中：rp,w为点p到井点的距离。运用极限分析可以得到
???Gp?w(x,y)Tdxdy??QwGp(rp,w) （4.4.57） T因此，在边界元模型中，已知流量井孔形成的源汇项可以直接以式（4.4.57）的形式加入到边界元方程中，不会增加边界元方程的数目。由于式（4.4.57）中并不包含井孔半径，因此在计算流场水头分布时，井孔半径并不是需要的参数。在这种情况下，如果直接调用式（4.4.14）计算井孔水头，由于井孔为奇点而无法进行。计算实际的进孔水位必须引入井孔半径，即计算距离井心半径为rw处的水头，有
Hw?Qw1lnrw?2?T2???[HB?G?H?(x,y)?G]ds???Gdxdy （4.4.58）
??n?n2?T式中：?(x,y)为除井孔之外的源汇项强度。
如果是已知水头的井孔，则需要把井孔流量Qw作为待求变量，而式（4.4.58）中的Hw为已知变量。式（4.4.58）因此加入边界节点的积分方程组，一起求解。
第五节 地下水流数值模型的反演
当含水层的时空离散模型已经建立，边界条件和初始条件也得到处理之后，输入一定
的含水层参数即能得到模型运行的结果，这种模拟过程被称为模型的正演。含水层的参数一般通过抽水试验等物理测试方法获取，这种方式得到的参数值在局部范围内可能是准确可靠的，但未必适合大范围的条件。通过模型正演人们往往发现数值模拟的结果与野外观测的地下水位或径流量不一致，而改变了模型的参数值之后，结果有所改善。于是，模型参数的确定问题就称为数值模拟可靠性的关键。这种校正型参数使模拟结果与观测结果保持一致的过程，就是模型反演，也被称为模型参数识别。
一、地下水流的逆问题
地下水流场的动态演变过程除了边界条件、初始条件和各种水文关系之外，参数也是重要的控制因素。当其他条件相同时，参数的差异也可能导致不同的流场特征。如果我们知道了流场的水头分布和水头变化特征，是否可以反过来推求边界条件、参数的取值和分布等控制因素呢？这就是地下水流的逆问题。
以一个河间地块模型为例（图4.16），设渗透系数为K的潜水含水层被两条无限长的河流切割，河流水位分别为H1和H2，含水层顶部有均匀的大于零的入渗补给强度?，在距离左侧河流x?a处有一个观测孔，其观测地下水位为ha。首先考虑稳定流的情况，假设河流为含水层的定水头边界，则描述地下水流的数学模型可表示为
K??h(h)???0 （4.5.1） ?x?xh(x?0)?H1;h(x?L)?H2 （4.5.2）
模型的正演结果是上述定解问题的解析解
2?L2x2?L2x2h?H?(H?H?)?() （4.5.3）
KLKL2212211在已知观测孔水位ha的情况下，可以求出含水层的渗透系数K，这就是上述定解问题的一个逆问题
h观测孔φH1H2ox=aLx 图4.16 两侧为已知水头边界的河间地块模型
K?L?h?H2a2?aL?L?a?21??a?H22?H21? （4.5.4）
在已知?、H1、H2、ha的情况下，可以求出一个唯一的渗透系数。
然而，地下水流的逆向问题往往是非唯一的，不同的含水层参数在某些条件下将形成完
全相同的流场。例如，在上述河间地块模型中，如果含水层没有入渗补给，即??0，则根据式（4.5.4）得到K?0。这是错误的，因为模型的正演结果实际上变为
2h2?H12??H2?H12?x （4.5.5） L不论含水层的渗透系数有多大，水头的分布都相同。在这种情况下，求含水层参数的逆向问题就有无穷个解，无法通过观测孔德水位来识别。当入渗补给强度?很小时，逆向问题也将接近这种非唯一状态。
地下水的逆向问题还有可能是不稳定的，即对观测误差具有敏感性。仍然以河间地块模型为例。如图4.17所示，有一渠道其水位高于地下水位，向两侧对称渗漏，右侧渗漏水量流向观测孔A和观测孔B，单宽流量为q，属未知。含水层的渗透系数为K，无入渗补给。该地下水稳定流的正向问题可描述为
K?h?x上述定解问题的解析式为
???h??h??0 （4.5.6） ?x??x?x?0?q；h?x?b??hB （4.5.7）
2q?b?x? （4.5.8） K2h2?hB?渗漏hABhAhBqOxx=ax=b
图4.17 渠道对称渗漏模型
由于单宽流量q属于未知变量，其逆向问题是已知hA和hB计算单宽流量，根据式（4.5.8）有
22hA?hB （4.5.9） q?K2?b?a?如果出于种种原因，观测孔A和B分别有观测误差?A和?B，则含有误差?q的单宽流量为
q??q?略去误差的二次项，有
K??hA??A?2??hB??B?2? （4.5.10）
?2?b?a??hA?A?hB?B （4.5.11）
b?a?q?K则相对误差约为
?r??qq?2?hA?A?hB?B? （4.5.12） 22hA?hB当含水层水力梯度很小，而饱和带厚度又很大时，根据式（4.5.12）可能会产生巨大的相对误差。例如，水力梯度为1%，A、B观测孔相距200m，hA=30m、hB=28m，观测孔A正好有0.5%的误差，而观测孔B正好有0.5%的负误差，则观测孔反求的单宽流量具有以下的相对误差
302?282?r?2?2?0.5%?14.5% （4.5.13）
30?282这意味着观测孔即使是只有0.5%的相对误差，也可能导致反演模型高达14.5%的相对误差，误差被放大了29倍。如果用数值模型来求解地下水流问题，则除了观测误差之外，还会引入模型本身的截断误差等计算误差，这种现象将更加严重。
因此，地下水流的逆问题可能是有唯一解的，但可能反演结果对观测误差和模型计算误差特别敏感。在某些情况下，逆向问题还可能具有非唯一性。这些因素增加了地下水流逆问题求解的复杂性，几十年来，国内外虽然做了不少研究工作，但在理论和方法上取得的进展并不大。
二、模型识别的依据和原则
既然地下水流的正演问题是根据含水层的边界条件、初始条件和非参数获得地下水的水头分布以及地下水的径流量，则一般来说地下水流的逆问题，是根据已知的水头分布和地下水径流量确定含水层的参数、边界条件和初始条件。而含水层的水头分布和地下水径流量，在绝大多数情况下只能通过调查钻孔、泉等方式观测到。因此，钻孔水位、泉流量等式模型识别的主要依据。其中又以钻孔水位作为最主要的依据。
当以钻孔的观测水位作为地下水流数值模型的反演依据时，人们经常会遇到的问题是：若干个观测孔能反求出多少个参数？或者，修要多少个观测孔来反演若干个参数？这涉及数值模型反演问题的唯一性。只有当反演问题具有唯一性时，反求参数才是有意义的。为了保证反演模型的唯一性和稳定性，模型识别需要注意以下的原则。
1）不用稳定流模型反演给水度、储水系数等只有非稳定流模型才使用的参数。
2）对于稳定流模型，尽量避免同时反求补给参数和渗透性参数。平面二维稳定流模型的偏微分方程一般可表示为
???H????H?Km?Km?y????0 （4.5.14） ?x??x??x??y??y?式中：m为含水层的饱和厚度；?为补给强度。在数值模型中，参数往往分区处理，在一个参数区内含水层假定为均质，因此上述水流方程可在同一参数内改写为
???H???Ky?H??m?0 （4.5.15） ???m????x??x??y?Kx?y?Kx可见参数Kx总是和补给强度参数联系在一起，难以分开进行模型识别。在上节与图4.17对应的稳定流模型中，侧向补给量q和渗透系数K总是联系在一起的，也不能同时对两者进行识别。在这种情况下，只有渗透系数和补给强度参数的综合参数才能够作为反演的目标。 3）在稳定流模型中，应尽量做到每个参数分区至少有一个作为识别依据的观测孔，不妨再用一维模型说明这个问题，如图4.18所示，分区模型中含水层在x?L处渗透系数由K1转变为K2，观测孔4水位已知，则正演模型的解析式可表示为
2h2?h4?2q?x4?x?,L?x?x4 （4.5.16） K22h2?h4?2q2q?x4?L???L?x?,0?x?L （4.5.17） K2K1在已知q、h4的情况下，根据式（4.5.16）和式（4.5.17）可知，必须依靠观测孔1、3或2、3才能同时求出K1、K2。如果用观测孔3，仅能反演K2。只用观测孔1或2、或同时用观测孔1和2也只能反演K1，而无法同时反演两个分区的渗透系数。
渗漏h1234qOK1K2xx=L
图4.18 分段潜水含水层模型
4）对于非稳定流模型，补给强度参数、给水度（或储水系数）和渗透系数3个参数中，只能同时对其中的2个参数进行识别。平面二维非稳定流模型的偏微分方程一般可表示为
???H????H??H （4.5.18） Km?Km???S?y??x??x??x??y??y??t在一个均质参数分区中，方程可以改写为
???H???Ky?H??S?H （4.5.19） m?m???????x??x??y?Kx?y?KxKx?t