CN115330132B - 一种用于突发污染事故中宽浅型河流的水质分布逆时反演的方法 - Google Patents

Abstract

本发明公开了一种用于突发污染事故中宽浅型河流水质分布逆时反演的方法，基于反问题的理念，将河流水质初始浓度分布确定问题提为地表水环境二维初始条件反问题。针对初始条件反问题的不适定性难题，采用摄动量正则化方法求解二维维初始条件反问题，用一族与原问题相邻近的适定问题的解去逼近原问题的解，通过算例验证了方法的可靠性。本发明全部过程均可由计算机编程完成，人工干预少，适用于突发水污染事故中宽浅型河流初始浓度分布复演，为河流水质保护与管理提供关键技术支撑。

Description

一种用于突发污染事故中宽浅型河流的水质分布逆时反演的方法

技术领域

本发明涉及水环境规划与管理技术领域，具体涉及一种用于突发污染事故中宽浅型河流的水质分布逆时反演的方法。

背景技术

在水环境保护领域，依据河流水质的初始浓度分布(t＝0)，根据污染物在水体中的对流扩散混合输运规律，利用水环境数学模型可以预测未来某一时刻(t＝T)污染物浓度分布信息，达到水质预测的目的。水质预测是进行水环境质量影响评价、污染物排放总量控制指标制订以及水污染控制系统规划与管理中必需的工作，在环境保护规划、环境影响评价、突发污染事件风险评估和预警预报、容量总量控制、海绵城市建设、黑臭水体治理等领域得到广泛的应用。由于这一过程符合自然演变过程，属于正问题的范畴。然而，在科学研究中，我们经常遇到这样问题：知道了某个事物的现在状态，希望了解其过去状态，这往往提为逆时反问题。在水环境保护领域，已知环境系统控制方程结构、参数、边界条件以及当前(t＝T)污染物浓度分布的部分信息，推求t<T时刻浓度分布，称为初始条件反问题，即时间反演问题，或称逆时问题。包含两种问题：求t＝0时的初始浓度分布和求过去0<t<T时浓度。由于当前一种问题解决时，后一种问题的求解可以通过正问题的求解得到解决。初始条件反问题的求解对于水环境管理领域具有一定的实用价值。例如，河流突发污染事故使得大量污染物质在极短的时间内排入河流，对河流水质与生态环境造成了巨大影响，甚至威胁人群健康和生命财产安全。为了掌握事故对河流水质造成的影响，有关部门会迅速组织应急监测，获得事故发生后的浓度分布。然而这一浓度分布是河流初始浓度与事故排放影响的双重作用结果，如何区分二者的贡献对于事故的责任认定和生态补偿非常重要。因此，科学、准确、快速地识别事故发生时的水质分布对于水污染事故的处置处理工作具有重要意义与价值。

初始条件反问题常常是不适定的，从物理上看，污染物输运扩散过程是不可逆的，其逆时问题求解是高度病态的。同时，污染物的降解等作用更加增加了逆时反演的难度。由于初始条件反问题存在较强的不适定性，目前很少有人从反问题求解角度来复演历史水质状况，需要发展相关技术方法。

针对不适定问题，人们发展了一些稳定的数值求解方法。其中，最具有普适性、在理论上最完备而且行之有效的方法，就是由著名学者Tikhonov以第一类算子(特别是积分算子)方程为基本数学框架，于20世纪60年代初创造性地提出、后来得到深入发展的正则化方法。其基本思想是：用一族与原问题相邻近的适定问题的解去逼近原问题的解。

摄动量正则化就是基于正则化方法的反问题求解方法，并且可以归结为算子理论的最优化方法，其基本理论如下。

考虑如下偏微分方程的初边值问题：

其中u为一向量函数，L为微分算子，B为边值条件算子，E为初值条件算子，c(x)为待定向量函数，L依赖于g(x),Ω为区域，为Ω的边界。

该类问题的初始条件反问题为：由附加条件来确定未知向量函数c(x)。该类反问题很容易转化为如下非线性算子方程的求解问题：

利用Tikhonov正则化方法可将其化为下述非线性最优化问题的解：

其中α为正则化参数，D为L²(Ω)上的稳定泛函。利用数值方法求解非线性泛函问题的解，便可求得反问题的数值解。

摄动量正则化法是根据算子识别的摄动法、线性化技术和函数逼近论提出的一种数值迭代方法，其核心过程为：

(1)建立迭代过程：

c_n+1(x)＝c_n(x)+δc_n(x)

其中，摄动量δc_n(x)由下列非线性最优化问题来确定：

(2)对上述最优化问题进行离散化，并采用线性化方法求得δc_n(x)的数值解，即求解非线性化最优化问题的局部极小值。

目前摄动量正则化已被应用于图像处理、参数识别等领域，尚未见将摄动量正则化方法用于开展宽浅型河流水质逆时反演的相关报导，这主要需要解决两个关键问题：一是宽浅型河流不同于一般的小型河流，需要考虑污染物浓度在纵向和横向的变化，如何构造与宽浅型河流二维污染物输运扩散方程相对应的正则算子R(u,α)；二是如何选择正则参数α＝α(δ)，使之与原始数据的误差水平δ相匹配。

发明内容

发明目的：本发明目的在于针对现有技术的不足，提供一种用于突发污染事故中宽浅型河流水质分布逆时反演的方法，本发明可以为河流历史水质复演提供关键技术支撑和科学依据，可推广应用于水源地保护、水功能区划、水环境综合整治、水安全格局优化以及水资源有效保护和合理利用等规划设计与研究管理工作，克服了河流水质逆时反演中存在的高度不适定性、污染物不守恒的问题。

技术方案：本发明所述的一种用于突发污染事故中宽浅型河流水质分布逆时反演的方法，包括如下步骤：

S1、收集事故发生河流水文资料，确定宽浅型河流的研究边界为x∈[0,l]，y∈[0,b]， x为沿着河长方向的纵向坐标，y为沿着河宽方向的横向坐标，确定河流纵向流速u，横向流速v；水深h；断面面积A；确定污染物的纵向扩散系数E_x、横向扩散系数E_y和污染物的降解系数为K₁；

S2：收集突发污染事故资料，瞬时排放源强M，确定研究的时段为[0,T]，t＝T时刻的浓度数据C(x,y,T)；

S3：宽浅型河流突发污染事故水质分布逆时反演问题，可提为如下二维对流扩散系统初始条件反问题：

水质分布逆时反演问题即为已知t＝T时刻的浓度分布C(x,y,T)，推求t<T时刻的浓度分布。构造如下算例，利用摄动量正则化方法求解。T为末时刻。

S4：采用摄动量-正则化方法求解初始条件反问题，获得C(x,y,0)，即为初始时刻的浓度分布C_h(x,y)；

S5：根据正问题的解析解

即可求得t<T任意时刻的浓度分布。

进一步地，S3具体的为：

式中，δ(x)、δ(y)分别是狄拉克函数。

进一步地，S4具体的为：

4.1确定正则化系数a及求解精度EPS，对空间坐标x∈[0,l]，y∈[0,b]进行等步长离散，得到离散点坐标(x_m、y_j)，其中m＝0,1,2,…,L，j＝0,1,2,…B。

L为将步骤1中的[0,l]分成L等分，共计有L+1个节点，节点坐标为x₁，x_2～x_m；

B为将步骤1中的[0.b]分成B等分，共计有B+1个节点，节点坐标为y₁，y_2～y_m；

4.2确定基函数族为使得函数式中k_i是实系数；

4.3取有限项逼近C(x,y,0)，确定一个n维实向量K^T＝(k₁,…,k_n)∈Rⁿ,其中n的大小取决于逼近精度的要求，一般取3项即可满足精度要求；

4.4任意给定初始的K值，是优化搜索的初始点；

4.5根据当前的K值，求解正问题，计算末时刻(t＝T)各离散点的对应浓度 u(x_m,y_j,T,k_i)；所有离散点浓度组成了矩阵U；

4.6由计算导数矩阵的值A，式中τ为微小扰动量，可取0.01；

4.7计算摄动量δK_i，δK_i＝(A^TA+a)^-1A^T(V-U)，其中A^T代表A的转置矩阵，V为已知的网格点浓度分布所对应的矩阵，U为计算得到的网格点浓度分布所对应的矩阵；

4.8计算K_i+1＝K_i+δK_i，当||δK_j||>EPS时返回步骤4.3，重复执行上述步骤，直到使范数||δK_j||≤EPS为止；

4.9得到符合精度要求的K值，从而求出未知项

有益效果：与现有技术相比，本发明的优点在于：

(1)本发明基于反问题的理念，将河流水质初始浓度分布确定问题提为地表水环境初始条件反问题。针对初始条件反问题的不适定性难题，采用摄动量正则化方法求解二维初始条件反问题，用一族与原问题相邻近的适定问题的解去逼近原问题的解，解决了河流水质逆时反演的不适定性难题；

(2)本发明污染物输运过程考虑了随流作用、扩散作用和污染物降解作用等主要过程，适用于保守型污染物和以一级生物降解为主的常规有机污染物，水质指标涵盖了当前水环境管理的大多数指标，如无机盐、COD、氨氮、TP、TN等，通过算例验证了方法的可靠性；

(3)本发明不仅适用于大江大河等宽浅型二维水体，也可适用于宽深比较小的小型河流，适用范围广泛；

(4)根据方法流程，全部过程均可由计算机完成，人工干预较少，提高了计算精度；

(5)本发明不仅为河流水质逆时反演提供了关键技术支撑和科学依据，还可推广应用于水功能区划、水环境综合整治、水资源有效保护和合理利用等规划设计与研究管理工作。

附图说明

图1为本发明的流程图；

图2为实施例1的t＝60s时的浓度分布；

图3为实施例2的t＝30s时的浓度分布；

图4为当δ＝0时，实施例2的计算值与精确值对比图；

图5为当δ＝0.01时，实施例2的计算值与精确值对比图；

图6为当δ＝0.1时，实施例2的计算值与精确值对比图；

图7为当δ＝0.3时，实施例2的计算值与精确值对比图。

具体实施方式

下面通过附图对本发明技术方案进行详细说明，但是本发明的保护范围不局限于所述实施例。

按照权利要求书中提到的方法，做出以下实施例1和2。

实施例1

已知某污染物降解系数为K＝4.2d^-1，河流纵向流速为1.5m/s，横向流速为0m/s，纵向扩散系数为50m²/s，横向扩散系数为10m²/s，河宽为30m，平均河深为2.0m。如果河流突发污染事故中瞬时向河流中心排放质量为200g的污染物，已知初始浓度 C(x,y,0)＝0.1mg/l，求下游污染物浓度的时空变化。这是典型的正问题，容易得到t＝60s时的浓度分布分别如图2所示。假设已知t＝60s时的浓度分布，尝试利用本发明根据图中网格点所对应的浓度数据推求初始浓度C(x,y,0)。

由于初始浓度分布为常数，基函数族设为{1}，C(x,0)＝k₁。正则参数取0.00001，利用摄动量正则化方法反演初始分布函数系数见表1

表1初始分布函数系数反演

由计算结果可以看到，当扰动较小时，均能得到精确解，当扰动较大时，误差值往往随着正则参数的减小而减小，达到某一最优值后，随后又开始增大。本算例中正则参数取0.001较优。

实施例2

河流水文条件与事故源于实施例1相同，即某污染物降解系数为K＝4.2d^-1，河流纵向流速为1.5m/s，横向流速为0m/s，纵向扩散系数为50m²/s，横向扩散系数为10m²/s，河宽为30m，平均河深为2.0m。河流突发污染事故中瞬时向河流中心排放质量为200g 的污染物。但初始浓度分布不是常数，以指数分布为例，C(x,y,0)＝0.1*exp(-0.01x). 容易求得30s后的浓度分布分别如图3所示。试利用本发明方法，根据如图所示的浓度分布(图中网格点所对应的数据)推求初始浓度分布C(x,y,0)。

取基函数族为{1,x,x²}，C(x,0)＝k₁+k₂x+k₃x²。正则参数取0.001，利用摄动量正则化方法反演初始分布函数系数见表2。初始状态浓度分布不同噪声水平下计算值与精确值的对比见图4～7。

表2初始分布函数系数反演

由图可见，计算值与精确值吻合较好，δ＝0.3时平均相对误差最大为6.92％。δ＝0.3 时正则化参数不同取值对反演结果的影响见表3。

表3正则化参数对反演结果的影响

由图可见，正则化参数的选取对反演结果有一定影响，通常存在一个较优的正则化参数。

Claims (1)

1.一种用于突发污染事故中宽浅型河流水质分布逆时反演的方法，其特征在于包括如下步骤：

S1、收集事故发生河流水文资料，确定宽浅型河流的研究边界为x∈[0,l]，y∈[0,b]，x为沿着河长方向的纵向坐标，y为沿着河宽方向的横向坐标，确定河流纵向流速u，横向流速v；水深h；断面面积A；确定污染物的纵向扩散系数E_x、横向扩散系数E_y和污染物的降解系数为K₁；

S2：收集突发污染事故资料，瞬时排放源强M，确定研究的时段为[0,T]，t＝T时刻的浓度数据C(x,y,T)；T为末时刻；

S3：宽浅型河流突发污染事故水质分布逆时反演问题，提为二维对流扩散系统初始条件反问题：

所述S3具体的为：

式中，δ(x)、δ(y)分别是狄拉克函数；

S4：采用摄动量-正则化方法求解初始条件反问题，获得C(x,y,0)，即为初始时刻的浓度分布C_h(x,y)；

S5：根据正问题的解析解

即可求得t<T任意时刻的浓度分布；

所述S4具体的为：

4.1确定正则化系数a及求解精度EPS，对空间坐标x∈[0,l]，y∈[0,b]进行等步长离散，得到离散点坐标(x_m、y_j)，其中m＝0,1,2,…,L，j＝0,1,2,…B；

L为将步骤1中的[0,l]分成L等分，共计有L+1个节点，节点坐标为x₁，x_2～x_m；

B为将步骤1中的[0.b]分成B等分，共计有B+1个节点，节点坐标为y₁，y_2～y_m；

4.2确定基函数族为

使得函数

式中

k_i是实系数；

4.3取有限项逼近C(x,y,0)，确定一个n维实向量K^T＝(k₁,…,k_n)∈Rⁿ，其中n的大小取决于逼近精度的要求，取3项即可满足精度要求；

4.4任意给定初始的K值，

是优化搜索的初始点；

4.5根据当前的K值，求解正问题，计算末时刻(t＝T)各离散点的对应浓度u(x_m,y_j,T,k_i)，所有离散点浓度组成了矩阵U；

4.6由

计算导数矩阵的值A，式中τ为微小扰动量，取0.01；

4.7计算摄动量δK_i，δK_i＝(A^TA+a)^-1A^T(V-U)，其中A^T代表A的转置矩阵，V为已知的网格点浓度分布所对应的矩阵，U为计算得到的网格点浓度分布所对应的矩阵；

4.8计算K_i+1＝K_i+δK_i，当||δK_j||>EPS时返回步骤4.3，重复执行上述步骤，直到使范数||δK_j||≤EPS为止；

4.9得到符合精度要求的K值，从而求出未知项

Images