基于梯度正則化聯合重構河流水質模型多項參數

邢利英, 孫三祥，張國珍

(1.蘭州交通大學環境與市政工程學院，甘肅 蘭州 730070；2.南陽師范學院土木與建筑工程學院，河南 南陽 473061)

基于梯度正則化聯合重構河流水質模型多項參數

邢利英1,2, 孫三祥1，張國珍1

(1.蘭州交通大學環境與市政工程學院，甘肅 蘭州 730070；2.南陽師范學院土木與建筑工程學院，河南 南陽 473061)

為識別河流水質模型參數，首先采用Crank-Nicolson有限差分格式離散控制方程，其次通過附加終端觀測值，設計梯度正則化算法重構河流水質模型的多項參數。結果表明，不管是順直河流的常系數模型，還是天然的線性相關或者線性無關的變系數模型，梯度正則化方法均能快速有效識別模型參數；初值和測量誤差對參數重構影響較小，充分證明梯度正則化算法是一種有效穩定地反演河流水質模型參數的方法。

水質模型；模型參數；參數重構；Crank-Nicolson；梯度正則化

隨著社會經濟的快速發展，河流污染日趨嚴重。及時查找河流污染的具體原因，建立河流的水質模型，掌握河流污染物濃度的時空變化，對于水環境預測、水污染控制治理和生態修復等諸多方面有重要的意義。其中水質模型參數估計是建立水質模型的一個非常重要環節，而水質模型應用的成敗很大程度上取決于參數估計是否正確[1-5]。通常模型參數可以通過模型試驗或經驗值法確定，然而模型試驗往往工作量較大、效率較低，而經驗值法很難保證其普適性。目前較為常用的方法是參數調試，實際上是根據實測的數據來調試參數，使得預測結果與實測數據吻合，即是基于初始邊界條件識別控制方程的參數[6-7]。

參數識別對于實際測量中污染源的位置、污染物排放時間、污染物擴散系數的測量等具有一定的理論參考和實際應用價值。目前在河流污染水質模型中，參數反演已經取得一定的成功。參數反演方法主要有隨機性優化算法和確定性算法[8]，隨機性優化算法主要包括遺傳算法、基于馬爾科夫鏈的貝葉斯法(Bayesian-Markov ChainMonte Carlo，B-MCMC)、模擬退火法、單純形法等；確定性算法有共軛梯度法、時空全域函數配點法(Globalspace-time Multi-quadric，MQ)、Landweber迭代法、Tikhonov正則化法等。朱嵩等[9]提出了有限體積法結合混合遺傳算法的參數識別算法。薛紅琴等[10]采用有限差分法結合單純形法的參數識別算法來反演天然河流的水質參數。江思珉等[11]分別利用單純形模擬退火法和Hooke-Jeeves吸引擴散粒子群混合算法求解地下水污染源反演的模擬—優化模型。邢利英等[12]應用改進的共軛梯度法識別地下水污染源。Li等[13-14]提出基于時空全域MQ函數配點法構造反演模型，求解地下水污染源識別初值反問題。李功勝等[15]成功地應用一種最佳迭代正則化算法反演污染源項，數值結果表明在測量數據存在干擾時，這種算法仍可以精確穩定地反演源項。范小平等[16]結合最佳攝動思想的遺傳算法確定了淄博灃水南部地區的地下水硫酸污染源強度。

上述研究中，參數識別方法均能有效地識別河流水質模型的某項參數，而對于多個參數的聯合識別研究較少。劉進慶等[17]應用梯度正則化方法重構地下水污染強度，該算法精度高，收斂速度快，穩定性好。本文應用具有二階精度的Crank-Nicolson格式離散控制方程，采用梯度正則化方法進行河流污染水質模型的參數識別，分別聯合重構常系數和變系數河流水質模型中的多項參數，并且考慮參數之間的耦合和非耦合情況。

1 Crank-Nicolson格式對控制方程的離散化

1.1 河流污染的控制方程及反問題

假定河流中污染物按一級動力學衰減，河流一維水質模型控制方程為

(1)

式中：Q為研究區域，Q=[0,L]×[0,T]；ρ(x,t)為測點(x,t)處的污染物濃度；u為河流的平均流速；Ex為x方向擴散系數；K為綜合降解系數;f(t)為邊值條件。

如果上述方程中的各參數u、Ex、K、f(t)均已知，方程(1)則構成了一個正問題，可以通過模型試驗和數值模擬方法預測河流下游污染物濃度ρ(x,t)，從而可以推斷河流的污染情況。然而大多數情況下，u、Ex、K往往不能準確地測得，因此上述問題也就轉化為參數識別反問題，即如何確定一組〈ρ(x,t),u,K,Ex〉滿足方程(1)。求解上述反問題，需要附加終端觀測值

ρT(x)=φ(x)

(2)

式中：ρT(x)為T時刻的終端觀測值；φ(x)為終端觀測函數。

由此，式(1)和(2)構成一個關于河流水質模型參數識別問題，可以應用迭代的方法求解。

1.2 Crank-Nicolson離散控制方程

xi=iΔxi=0,1,…,M

(3)

tn=nΔtn=0,1,…,N

(4)

Crank-Nicolson格式是一種加權隱式格式，其截斷誤差為O(Δx2+Δt2)，即具有二階精度，被廣泛地應用在偏微分方程的求解中[1,19]。上述河流污染的控制方程建立Crank-Nicolson差分格式為

(5)

經過簡化，式(5)可轉化為

(6)

方程(6)的矩陣表達式為

Gρn+1+Z=Bρn+Qn=1,2,…,N

(7)

其中

(8)

(9)

(10)

(11)

G和B是典型的三對角矩陣

(12)

(13)

對于所有的內部節點都可列出方程(7), 利用初邊界條件可得到一個主對角元占優的三對角方程組，可以應用追趕法求解上述方程組。

1.3 Crank-Nicolson差分方程的穩定性

穩定性是判斷反問題是否適定的條件之一，本文應用極值原理證明Crank-Nicolson差分格式的穩定性[20-21]。方程(6)可變形為

(14)

(2-2s+2KΔt)‖ρn‖∞+(r+s)‖ρn‖∞

0≤i≤M-1,0≤n≤N-1

(15)

(16)

(17)

當0≤i≤M-1則有

‖ρn+1‖∞≤(1+KΔt)‖ρn‖∞

(18)

式(16)對于任意的0≤i≤M得

‖ρn‖∞≤(1+KΔt)‖f‖∞

(19)

(20)

(21)

基于上述穩定性的證明，同理有

‖εn‖∞≤(1+KΔt)‖Ψ‖∞0≤i≤M

(22)

2 梯度正則化法

a. 建立迭代過程：

Dn+1(x)=Dn(x)+δDn(x)

(23)

其中，攝動量δDn(x)是由下列非線性最優化問題來確定：

αE[δDn]

(24)

式中：α為正則化參數；E[δDn]為[δDn]的穩定化泛函；A為導數矩陣；Jα為非線性泛函。

b. 對上述最優化問題進行離散，采用線性化方法求得δDn(x)的數值解，即求解非線性化最優化問題的局部極小值。

根據以上的算法，梯度正則化法的具體迭代步驟[3,15]如下：

步驟1：確定正則化參數，步長及初始值，(x,t)的離散點以及求解精度EPS；

步驟2：根據已知參數求解方程(1)，得到ρ(xi,T,Wj)；

步驟3：設在計算區間有M個離散點xm(m=1,2,…,M)，由

步驟4：計算δWi=(ATA+a)-1AT(V-D), 其中V為附加數據；

3 常系數模型的多參數重構

選取3個典型的數值算例驗證梯度正則化算法的有效性和普適性，采用孿生實驗，即設定準確參數，利用同一模型生成觀測數據，設參數未知進行識別試驗。

3.1 正問題的求解

設一均勻河段長30 km，計算時間取1 h，離散的空間步長取1.5 km，時間步長取1/15 h，污染物的質量濃度邊界值ρ0取1 mg/L，河流的平均流速u取5 km/h，彌散系數Ex取2 km2/h，污染物的一級降解速率K取0.015 h-1，算例數據引自文獻[9]。此種離散情況下，Peclet數為3.75，屬于非對流占優問題，對流項采用中心加權差分格式是比較精確的，數值彌散與振蕩不明顯。采用Crank-Nicolson的有限差分格式求解的計算結果見圖1，可以看出，計算結果較為光滑，這是由于Crank-Nicolson是二階精度，截斷誤差為O(Δx2+Δt2)。

3.2 多參數重構

常系數模型參數估計結果見表1，其中求解精度EPS取10-3，初值相同，真值為u=5 km/h，Ex=2 km2/h，K=0.015 h-1，與相關文獻的共軛梯度法[12]和混合遺傳算法[9]的反演結果相比，發現梯度正則化方法比其他2種方法更快、更有效地重構河流污染的模型參數，因此下面的模型參數識別中，均采用梯度正則化方法重構模型參數。

表1 幾種方法常系數模型參數估計結果

3.3 初值對參數重構的影響

選取不同初值，應用梯度正則化方法重構參數的結果見表2，其中x*為精確值，從表中數據可以看出，初值對于模型參數的反演結果影響不大。當選取初值偏差達到50%時，污染物的降解系數K與真值相當接近，相對誤差僅為0.7%；而當初值偏差達到100%時，K的反演誤差為8.6%，在允許的誤差范圍內，充分說明初值對于梯度正則化反演模型參數影響較小。

3.4 測量誤差對參數重構的影響

即使是精密儀器測量的觀測數據也會有誤差，為了更真實地重構模型參數，考慮不同隨機噪聲對反演結果的影響。采用下述形式表示含有測量誤差的觀測數據[21]：

表2 不同初值下常系數模型參數重構結果

ρδ(x,T)=ρ(x,T)[1+δrandom(x)]

(25)

式中：ρδ(x,T)為帶有測量誤差的觀測數據；random(x)為一隨機函數。初值取u=10 km/h，Ex=6 km2/h，K=0.02 h-1，試驗30次的平均結果見表3。

不同噪聲水平的常系數模型參數重構結果表明，當測量誤差δ為5%時，降解系數K均值的相對誤差為0.4%，標準差為0.002 6；當噪聲水平δ為20%時，降解系數K反演結果均值的相對誤差為20%，標準差為0.005 5，反演結果在合理的范圍內，并且迭代次數對于反演結果影響不大，故梯度正則化對于一維常系數河流模型多參數重構是有效的。

4 變系數模型的多參數重構

4.1 線性相關的變系數模型多參數重構

對于順直河道，水質模型中的參數往往都為常數，然而對于天然河道，水質模型參數都是變化的，并且參數之間可能存在較強的耦合性。假設河道的平均流速為空間位置的線性函數(如明渠漸變流)u=ax+b，而根據費希爾的估算公式，彌散系數為河道平均流速的二次函數,即Ex=βu2。如果平均流速函數已知，則需要重構的模型參數有系數β和污染物的降解速率K，河流污染模型其他參數參考常系數模型的參數值。采用孿生實驗，假設u=0.2x+5，β=0.08，采用Crank-Nicolson的有限差分格式求解的結果(圖2)作為終端觀測值基礎。

圖2 Crank-Nicolson格式求解的變系數模型不同時刻濃度

4.1.1 初值對參數重構的影響

選擇不同的初值，求解精度EPS取10-4，變系數模型多參數重構的結果見表4，其中x*為精確值。從表4中可以得出，不同的初值和迭代次數對反演結果的影響較小。究其原因可能是反演的參數本身數量級較小，對于重構結果影響不大。

4.1.2 測量誤差對參數重構的影響

由于觀測數據可能存在一定的測量誤差，應用梯度正則化迭代結果也會產生相應的誤差。為確定測量誤差對參數重構的影響，選擇初值β=0.12,K=0.03 h-1，表5為不同測量誤差下變系數模型參數重構結果。當測量誤差為20%時，重構參數的均值為β=0.082 0，K=0.1 587 h-1；均值的相對誤差分別為2.5%、5.8%；標準差也較小，分別為0.006 2、0.001 2，充分說明梯度正則化方法有較大的容錯性，適用于天然河道污染水質線性相關模型參數重構。

表3 不同測量誤差下常系數模型參數重構結果

表4 不同初值下線性相關變系數模型參數重構結果

表5 不同測量誤差下線性相關的變系數模型參數重構結果

表6 不同初值下線性無關的變系數模型參數重構結果

4.2 線性無關的變系數模型多參數重構

如果變系數模型參數之間不存在耦合性，假設河流的平均流速為u=ax+b，河流污染物的擴散系數為Ex=cx2+dx+e，污染物的降解速率K=0.015,求解精度EPS取10-4，應用梯度正則化方法研究線性無關的變系數模型多參數重構。

4.2.1 初值對參數重構的影響

表6列出了不同初值下線性無關的變系數模型參數的重構結果，從表6中可以看出，初值對反演結果的影響不大。

4.2.2 測量誤差對參數重構的影響

類似地，考慮測量誤差的影響。初值選擇u=0.5x+5，Ex=0.02x2+0.3x+1.0,不同測量誤差下線性無關變系數模型參數重構結果見表7。從表中可以得出，當測量誤差為5%時，河流平均流速u=0.399x+2.98，污染物彌散系數Ex=0.0167x2+0.243x+0.905；當測量誤差為20%時，平均流速u=0.398x+3.04，彌散系數Ex=0.017 1x2+0.23x+0.915，均與真實解u=0.4x+3，Ex=0.016x2+0.24x+0.9很接近。可見，測量誤差的影響很小，進一步地證明梯度正則化法穩定性好。

表7 不同測量誤差下線性無關的變系數模型參數重構結果

5 結 論

應用Crank-Nicolson格式離散一維河流水質方程，比較梯度正則化法、共軛梯度法以及混合遺傳算法重構順直河流水質模型參數的結果，得出梯度正則化法是一種有效地參數重構方法，進而應用梯度正則化方法重構順直和天然河流水質模型的多參數。結果表明，不管是順直常系數河道還是天然變系數(線性相關或線性無關)河道，梯度正則化方法均能有效地聯合重構多參數。此外，模擬結果表明，初值和測量誤差對模擬結果影響較小，進一步證明梯度正則化方法是一種快速有效的重構河流水質模型參數的方法。

梯度正則化方法不但適用于一維河流水質模型參數的重構，而且可以推廣到二維情形，具有較強的通用性。具體問題應用成功的關鍵在于河流水質模型正問題的求解和優化參數的選取。

[1] 孫培德.環境系統模型及數值模擬[M].北京：中國環境科學出版社,2005.

[2] 王紅旗,秦成,陳美陽. 地下水水源地污染防治優先性研究[J]. 中國環境科學, 2011, 31(5): 876-880. (WANG Hongqi, QIN Cheng, CHEN Meiyang.Study on priorities of prevention and control of groundwater source pollution [J]. China Environmental Science, 2011,31(5):876-880.(in Chinese))

[3] 李功勝, 姚德. 擴散模型的源項反演及其應用[M]. 北京：科學出版社, 2014.

[4] YEH W G. Review: optimization methods for groundwater modeling and management [J]. Hydrogeology Journal, 2015, 23(6): 1-15.

[5] 軒曉博, 逄勇, 李一平, 等．金屬礦區重金屬遷移對水體影響的數值模擬[J]. 水資源保護，2015, 31(2): 30-35. (XUAN Xiaobo, PANG Yong, LI Yiping, et al. Numerical simulation of influence of heavy metal migration on water in metallic mining areas [J]. Water Resources Protection, 2015, 31(2):30-35. (in Chinese))

[6] 江思珉,張亞力,蔡奕,等. 單純形模擬退火算法反演地下水污染源強度[J]. 同濟大學學報(自然科學版), 2013, 41(2): 253-257. (JIANG Simin, ZHANG Yali, CAI Yi, et al. Groundwater contamiant identification by hybrid simplex method of simulated annealing[J]. Journal of Tongji University (Nature Science), 2013, 41(2): 253-257. (in Chinese))

[7] HAZART A, DUBOST S, GIOVANNELLI J F, et al. Erratum to “inverse transport problem of estimating point-like source using a Bayesian parametric method with MCMC” [J]. Signal Processing, 2014, 96(5):346-361.

[8] 肖傳寧,盧文喜,安永凱,等. 基于兩種耦合方法的模擬-優化模型在地下水污染源識別中的對比[J]. 中國環境科學, 2015, 35(8): 2393-2399. (XIAO Chuanning, LU Wenxi, AN Yongkai, et al. Simulation-optimization model of identification of groundwater pollution sources based on two coupling method[J]. China Environmental Science. 2015,35( 8): 2393-2399. (in Chinese))

[9] 朱嵩,毛根海,劉國華.基于FVM-HGA的河流水質模型多參數識別[J].水力發電學報, 2007, 26(6):91-95. (ZHU Song, MAO Genhai, LIU Guohua. Parameters identification of river water quality model based on finite volume method-hybrid genetic algorithm[J]. Journal of Hydroelectric Engineering, 2007, 26(6): 91-95. (in Chinese))

[10] 薛紅琴,趙塵,劉曉東,等.確定天然河流縱向離散系數的有限差分-單純形法[J].解放軍理工大學學報(自然科學版), 2012, 13(2): 214-218. (XUE Hongqin, ZHAO Chen, LIU Xiaodong, et al. Finite difference method: simplex method for the determination of longitudinal dispersion coefficient in natural rive[J]. Journal of PLA University of Science and Technology (Natural Science Edition), 2012, 13(2): 214-218. (in Chinese))

[11] 江思珉,王佩,施小清,等.地下水污染源反演的Hooke-Jeeves吸引擴散粒子群混合算法[J]. 吉林大學學報(地球科學版),2012,42(6):1866-1872.(JIANG Simin, WANG Pei, SHI Xiaoqing, et al. Groundwater contaminant source identification by hybrid Hooke-Jeeves and attractive repulsive particle swarm optimization method [J].Journal of Jilin University (Earth Science Edition), 2012,42(6):1866-1872.(in Chinese))

[12] 邢利英,張國珍.基于改進的共軛梯度法重構地下水污染源項 [J]. 水資源保護,2017, 33(3):42-46.(XING Liying, ZHANG Guozhen. Reconstruction of the groundwater pollution source with improved conjugate gradient algorithm [J]. Water Resources Protection, 2017, 33(3):42-46. (in Chinese))

[13] LI Z, MAO X Z. Global multi-quadric collocation method for groundwater contaminant source identification [J]. Environmental Modelling & Software, 2011, 26(12): 1611-1621.

[14] 李子, 毛獻忠, 周康. 污染物非恒定輸運逆過程反演模型研究[J]. 水力發電學報, 2013, 32(6): 115-121. (LI Zi, MAO Xianzhong, ZHOU Kang. Inversion model for the inverse process of unsteady pollutant transport [J]. Journal of Hydroelectric Engineering, 2013, 32(6): 115-121. (in Chinese))

[15] 李功勝, 姚德, 馬昱, 等. 一維溶質運移源 (匯) 項系數反演的迭代正則化算法[J]. 地球物理學報, 2008, 51(2): 582-588. (LI Gongsheng, YAO De, MA Yu, et al. An iterative regularization algorithm for determining source (sink) coefficientin 1-D solute transportation [J]. Chinese Journal of Geophysics, 2008, 51(2): 582-588. (in Chinese))

[16] 范小平, 李功勝. 確定地下水污染強度的一種改進的遺傳算法[J]. 計算物理, 2007, 24(2): 187-191. (FAN Xiaoping, LI Gongsheng. An improved genetic algorithm for groundwater pollution [J]. Chinese Journal of Computational Physics, 2007, 24(2): 187-191.(in Chinese))

[17] 劉進慶,李功勝,馬昱. 確定地下水污染強度的梯度正則化方法[J]. 山東理工大學學報(自然科學版), 2007, 21(2): 17-20. (LIU Jinqing, LI Gongsheng, MA Yu. Gradient-regulation method for determining a pollution source term in groundwater [J].Journal of Shandong University of Technology (Natural Science Edition), 2007, 21(2): 17-20. (in Chinese))

[18] 王福軍. 計算流體動力學分析：CFD軟件原理與應用[M].北京：清華大學出版社，2004.

[19] 劉繼軍. 不適定問題的正則化方法及其應用[M]. 北京：科學出版社，2005.

[20] 孫志忠.偏微分方程數值解法[M]. 北京：科學出版社，2005.

[21] DENG Zuicha, QIAN Kun, RAO Xiaobo, et.al. An inverse problem of identifying the source coefficient in a degenerate heat equation [J]. Inverse Problems in Science and Engineering, 2014, 23(3): 498-517.

Reconstruction of river water quality model parameters with gradient regularization approach

XING Liying1, 2, SUN Sanxiang1, ZHANG Guozhen1

(1.SchoolofEnvironmentalandMunicipalEngineering,LanzhouJiaotongUniversity,Lanzhou730070,China; 2.AcademyofCivilEngineeringandArchitecture,NanyangNormalUniversity,Nanyang473061,China)

In order to identify the parameters of the river water quality model, we adopted the Crank-Nicolson finite difference scheme to discrete the control equations first, and then appended several additional terminal observations, and designed an effective gradient regularization algorithm to reconstruct the model parameters. Numerical results show that for both the constant coefficient model of straight channels and the natural coupled or non-coupled variable coefficient model, the proposed method can quickly and effectively identify the model parameters. In addition, the influences of the initial value and measuring error are insignificant, which further demonstrates that the gradient regularization algorithm is effective to inverse the parameters of the river water quality model.

water quality model; model parameter；parameter reconstruction; Crank-Nicolson；gradient regularization

10.3880/j.issn.1004-6933.2017.04.009

國家自然科學基金(51468033)

邢利英(1980—)，女，講師，博士研究生，主要從事環境水力學的反問題研究。E-mail: lyxing2011@163.com

張國珍，教授。E-mail: guozhenzhang163@163.com

X523

A

1004-6933(2017)04-0055-07

黑臭水體評價方法比較(

2016-11-03 編輯：王 芳)