Page 102 - 2023年第54卷第10期

(
ｍ；Ｈ ｘ，ｙ，ｚ，ｔ ) 为第一类边界条件已知水头函数，ｍ；Ｑ ｘ，ｙ，ｚ，ｔ ) 为第二类边界条件已知单位
１ (
  →
、
、
ｎ
面积流量函数，ｍ?ｄ；Ｋ为边界法线方向ｎ的渗透系数，ｍ?ｄ；α和 β为已知函数；Γ １ Γ ２ Γ ３ 分别为
研究区第一类、第二类和第三类边界。
根据监测井地下水水位的观测值，记为 ｄ ＝ （ｄ，ｄ，…，ｄ），以模型模拟值 Ｆ（ θ ）与观测值 ｄ之
１ ２ ｎ
间误差（ ε ＝ ｄ － Ｆ（ θ ））为目标函数，推求 ε 最小时的模型参数 θ （如 Ｋ等），即为参数反演问题。
２．２ 随机参数场的 Ｋａｒｈｕｎｅｎ － Ｌｏèｖｅ展开 对于某一非均质孔隙介质含水层，如渗透系数 Ｋ值空间变
化，形成参数反演的高维问题。相比于裂隙和岩溶介质显著的非均质性，孔隙介质在一定尺度上可以
认为相对连续并具有自相关性。对此，可利用 Ｋａｒｈｕｎｅｎ － Ｌｏèｖｅ（Ｋ － Ｌ）展开 ［２０］ ，将 Ｋ值的空间变化表
述为确定性函数和随机场的线性叠加，以达到降低维度，即：
!

＾ 珔 ＋
Ｋ＝Ｋ ∑ 槡 ｆ（Ｓ） （２）
ξ ｎ λ ｎ ｎ
ｎ ＝１
式中： 珔 为服从标准高斯分布的随机数，彼此相互独立，ｎ为随机数项数；
Ｋ为参数空间分布均值；ξ ｎ
和 ｆ（ｘ）分别协方差函数的特征向量和特征方程。
λ ｎ ｎ
在地下水数值模拟中，Ｋ的自然对数（ｌｎＫ）的协方差函数 Ｃ（Ｓ，Ｓ）服从指数分布，即任意两点
１ ２
（ Ｓ 和 Ｓ）满足：
２
１
η )
Ｃ（Ｓ，Ｓ） ＝Ｃ（ｘ，ｙ；ｘ，ｙ） ＝ σｅ ( ｘ－ ｘ ＋ ｙ－ ｙ （３）
１ ２
１ ２
２ －
η
２
２
２
１
１
２
１
１
２
分别为 ｘ和 ｙ方向上的自相关长度；σ 为方差。
式中：η １ 和 η ２
２．３ 基于 Ｋｒｉｇｉｎｇ插值的替代模型 为了提高参数反演的效率，使用替代模型模仿地下水数值模型输
入－ 输出关系。Ｋｒｉｇｉｎｇ替代模型是全局模型加上局部的偏差的组合。任意一点 Ｋｒｉｇｉｎｇ空间插值 ＾ｙ 与自
变量 θ 之间的关系式为：
ｙ
＾ ＝ Ｇ（ θ ） ＋ Ｚ（ θ ） （４）
式中：Ｇ（ θ ）为 θ 的已知全局函数；Ｚ（ θ ）为高斯随机函数，为系统局部偏差，其满足 Ｅ（Ｚ） ＝０，Ｖａｒ（Ｚ） ＝
２
σ，且任一两点 Ｚ（ θ ）的相关函数可表示为：
（ｊ）
２
（ｉ）
（ｉ）
（ｊ）
ｃｏｖ［Ｚ（ θ ），Ｚ（ θ ）］ ＝ σＲ（ θ ，θ ） （５）
（ｊ）
（ｉ）
式中：Ｒ（ θ ，θ ）为任意两个样本之间的协方差函数（常用有线性函数、高斯函数、指数函数等）。
Ｋｒｉｇｉｎｇ替代模型实现步骤如下：
（１）建立模型输入（参数 θ ） － 输出（水头 ｙ）集合。在各参数先验分布上抽取 Ｎ组样本作为输入集
合，调用地下水数值模型，模拟每个监测井不同时段水头值 ｙ作为样本输出集合。
（ ２）Ｋｒｉｇｉｎｇ替代模型训练。将步骤（１）样本的输入和输出集合代入替代模型进行训练，建立每个
监测井的水头替代模型
（ ３）替代模型验证。重新随机抽取 Ｍ组参数样本作为检验样本输入集合，分别调用地下水数值模
型、替代模型模拟得到第 ｉ个监测井的水头（分别为 ｙ 和 ＾ｙ ）。如果 ｙ 与 ＾ｙ 十分接近，则替代
ｉ，ｏｕｔ ｉ，ｏｕｔ ｉ，ｏｕｔ ｉ，ｏｕｔ
模型可行。
２．４ 参数反演方法
２．４．１ 贝叶斯方法 从观测数据中推断出未知参数的可能取值，即后验概率分布 ｐ（ θ"ｄ）表述为：
ｐ（ｄ "θ ）ｐ（ θ ）
ｐ（ θ"ｄ） ＝ （６）
ｐ（ｄ）
式中：θ 为待反演的参数；ｄ为地下水水头观测值；ｐ（ θ ）为 θ 的先验概率密度函数；ｐ（ｄ "θ ）为似然函
数，表示模拟值与观测值 ｄ的接近程度；ｐ（ｄ）为归一化积分常数。
通常根据野外调查、试验等信息得出参数某一取值范围，并假设其服从某一分布，即参数的先验
， 的
分布 ｐ（ θ ）。假设模型中未知参数共有 ｎ个，即 θ ＝ （ θ １ θ ２ ，…，θ ｎ ），服从均匀分布，则参数 θ ｉ
ｐ（ θ ）为：
３
— １ ２ ８ —