基于相關系數優化法的河流突發污染源項識別

陳媛華,王 鵬,姜繼平,郭 亮 (哈爾濱工業大學市政環境工程學院,黑龍江 哈爾濱 150090)

基于相關系數優化法的河流突發污染源項識別

陳媛華,王 鵬*,姜繼平,郭 亮 (哈爾濱工業大學市政環境工程學院,黑龍江 哈爾濱 150090)

基于相關系數優化法,結合地表水環境特征和污染物水質過程特征,推導出一維河道中單點源瞬時排放的源項反演算法,得到了污染源排放特征與河流環境特征參數的反演公式.采用假想算例進行數值試驗,綜合分析了流速信息、污染物衰減、監測距離、監測數據誤差及中間參數?T選取等因素對反演結果的影響,確定了該方法的適用條件和最優條件的尋找方式.間隔10min進行兩次監測采樣,若監測誤差小于5%,反演結果的相關系數達到-0.97,污染源位置和排放量反演結果的相對誤差均小于4%,綜合相對誤差在2%以內.并且方法具有監測布點簡單高效,數據需求低,編程簡單等優點,值得在環境應急管理中進行實際應用.

點源污染源；河流突發污染；源項反演；參數識別；相關系數優化

世界各地都受到化學品泄漏事故的危害[1],其中地表水環境污染所占比重及對社會影響和損失最為顯著.一些突發性水污染事故很有可能在事故初期不知道污染源的位置以及排放特征,甚至有些污染事故在污染團流出管轄范圍后管理者才知道污染的發生,這對流域水環境管理十分不利. 因此需要找到有效的手段減少事故危害.上述問題可以描述為,已知污染物時空分布求解模型參數、結構和邊界條件等信息,在數學上可以歸納為反問題[2-5].反問題廣泛存在于環境水力學的各個研究領域,如廢水排放污染源控制問題,水環境容量計算問題,污染物總量控制及分配問題等[6-8].

目前,地表水污染源反演研究越來越受到關注.如 Boano等[9]應用地學統計的方法識別河流中的污染源,張玉超等[10]采用支持向量機方法反演太湖葉綠素a的濃度分布,Cheng等[11]使用反向位置概率密度函數法和CCHE2D模型程序對河流污染源進行了識別,朱嵩等[12-13]采用貝葉斯方法進行了污染源識別和模型參數反演研究,牟行洋[14]利用微分進化算法研究了污染物源項識別反問題.

相關系數優化算法在地下水污染源識別問題研究中曾有應用[15],但地下水系統與地表水存在較大差異,此算法用于識別河流污染源的可行性有待研究.2007年郭建青等[16]嘗試將其引入地表水參數識別體系,但尚未考慮污染物質在水體中的遷移轉化規律(水質過程),并且算例的濃度數據直接來源于解析解的計算結果.因此極大地限制了算法的有效性和普適性.本研究在此基礎上,結合河流的環境特征和污染物水質過程,考慮水力-水質耦合過程改進算法,得到了適用于可降解污染物的一維概化的地表水污染源識別方法.此算法可以反演出污染源特征(污染物排放量、污染源位置及初始排放時刻),同時計算出擴散系數和水流速,并考察濃度監測數據存在誤差條件下算法的準確性和有效性.但本算法只適用于穩態流動條件下地表水污染物的遷移轉化過程.

1 相關系數優化反演算法及修正

一維均一地表水瞬時點源排放流場的解析解[17]在許多教科書中可以找到:

式中:C(x,t)為瞬時脈沖排污引起的污染物濃度在水體的時空分布,mg/L;M為瞬時投放的污染物質量,kg;A為河流斷面積,m2;Dx為縱向彌散系數,m2/min; k為污染物一級反應動力學衰減速率系數,min-1; u為河流平均速度, m/min; x0為污染源位置,m;t0為污染源排放初始時刻,min.

同一時刻監測水體中不同斷面 xi污染物的濃度Ci,為簡化計算將首次監測時刻設定為0時刻(t=0),污染物初始排放時刻先于監測時刻,故為負值.根據已知的一系列濃度數據求得5個參數, M、Dx、u、x0和t0(假定水的濃度ρ為1).

基于上述概化的水質模型和解析解,對以監測點與污染源位置的距離為自變量,污染物濃度對數值為變量的線性函數的相關系數取極值.求得污染源位置,然后選用另一組監測數據,求得其他參數.具體過程如下:

首先,對式(1)取對數,結果如下:

其中:

根據式(2)在直角坐標系中繪制濃度對數值關于(x–xm)2的曲線,可得到一條以mx為斜率、以bx為截距的直線.但由于污染源位置xm未知,無法確定橫軸上各點的坐標, 因此尚不能確定此方程的斜率與截距.

對于任意的 xm值,可對式(2)進行線性回歸,并且計算回歸的相關系數.目標就是選取一個適當的xm值使相關系數達到最小,即最接近-1.式(2)的相關系數定義為:

n

式中:nx為監測點數;為平均距離,其表達式為:

關于xm對相關系數最小化:

求導,可得:

其中,

可見,根據式(10)～(15)可以直接估算出污染源的位置(xm).實現了污染源識別的第一步.

確定xm后,通過式(4)可以計算得到Xi,然后,可對式(2)進行最小二乘法線性擬合.

可得,

從式(10)、式(16) 和式(17)可得xm,mx和bx.通過更進一步的檢驗發現只有 2個可以獨立的用于求解的未知參數.因此,其余參數需要通過附加的條件確定.由于未知參數不能通過同一時刻的污染物濃度數據確定,其與監測點數目及濃度的時空分布無關.這些信息一定通過不同時間監測的污染物濃度提供.

為得到污染源的信息,需要估算全部的未知參數.首先,在取樣的時刻(t=0)改寫式(1):

根據式(18)可以計算得到3個不同組合的參數:M,Dxt0和x0-ut0,其估計取值如下:

然后為了獨立地求解出各個參數,污染物濃度數據必須取自一定時間后(?t),同樣按照式(10)和式(16)求得xm′和mx′,它們的定義是:

最后,根據上式能夠計算出其余4個參數:

此外,根據式(21)可以計算 M,得到這些參數后,就可以實現反演計算過程.算法實際應用過程見圖1.

圖1 算法實際應用過程Fig.1 Process of applying the inversion algorithm

2 案例應用及反演計算影響因素分析

2.1 算例應用

2.1.1 算例來源 某點源污染源非常規瞬時排放進河流中一定數量的化學物質 α.事發一段時間后,預警系統進行響應,沿水流方向設置一系列監測點反演其排放歷史,其位置見表 1.將解析解(式1)計算得到的結果Ci作為相應斷面濃度數據的期望,加入高斯分布的白噪聲[11],此處為5%.詳細“觀測”數據見表1.其他參數如下: M=3000kg; u=18m/min;Dx=60m2/min;k=0.001min-1;t0=-85min (第一次取樣前 85min);x0=-1000m; A=20m2.為確定所有參數,在10min后進行第二次采樣.

表1 瞬時點源算例的“觀測”數據Table 1 Observed concentration data for the instantaneous source case

表2 算例污染源與參數反演結果Table 2 Estimated and real parameters for the conservative and decay instantaneous source case

圖2 t=0時刻瞬時排放點源的預測濃度數據和“觀測”濃度數據Fig.2 Predicted and observed concentration for a conservative and a decay solute at t=0

2.1.2 反演結果 基于修正的相關系數優化法的反演計算公式,估算出未知參數見表 2.同時表2列出了將上述算例在不考慮污染物降解條件下進行反演的計算結果.估計參數與真實值接近,反演效果十分理想.圖1為第一次采樣時刻(t=0),實測數據、反演數據與理論數據的曲線.

由表2和圖2可見,計算結果較接近假想算例中的真實值,相關系數絕對值達到了0.9747.

2.2 反演計算影響因素分析

2.2.1 河流流速的影響 實際研究中,河流的流速可方便獲得.因此,考察水流速度條件已知與否對反演結果的影響.

由表3可知,速度條件對相關系數R沒有影響.但其他參數的計算結果有所差別.本算例中流速已知的反演結果比未知時略差.

表3 流速對瞬時點源算例計算結果的影響Table 3 The impact of known velocity and unknown velocity for the instantaneous source case

2.2.2 污染物降解過程的影響 BOD是反映水體中污染物濃度的綜合指標之一,其在河流中的降解系數一般為 0.039～20d-1[18].參照 BOD的降解系數范圍,考察污染物在河流中遷移轉化過程對反演結果的影響,因此選取化學物質 α的降解系數為0.001min-1,即1.44d-1.由圖3可見,對于可降解的污染物質(k=0.001min-1),其在水體中的總質量隨距離的增大而減少,因而在同一時刻同一位置其濃度小于非降解污染物;對于保守物質(k=0),則濃度計算結果偏高.但從反演結果(表2)可以看出,兩者的計算誤差均在5%之內,其中保守物質的計算結果較優于非保守物質,而考慮可降解物質的算法更適合應用于實際情況.

2.2.3 采樣位置和時間的影響 構建一個評價指標,式(28)用向量的方法計算相對誤差,用于考察采樣位置和時間同污染源間的時空距離對算例中反演結果的綜合影響.相對誤差式中:為反演求得的參數向量;I為向量第i個參數;為真實值的參數向量.

由圖4可見,在反演過程中,監測位置與時間是一對耦合的參數,共同影響反演結果.在二維的時空坐標下,源排放后 65～85min內在距污染源800～1000m的范圍內進行第一次采樣, 反演效果較好.

圖4 監測位置和時間對反演結果的綜合影響Fig.4 The impact of monitoring location and time on the inversion results

2.2.4 監測數據誤差及兩次采樣時間間隔對反演結果的影響 利用相對誤差考察采樣噪聲和時間間隔(?t)對算例中反演計算結果的綜合影響.

由圖5可見,當采樣的濃度噪聲不超過5%,取樣時間間隔對反演結果的影響較小;當濃度誤差超過 5%時,隨著誤差增大,計算結果與真實值間的偏差越大.

圖5 采樣噪聲和時間間隔對反演結果的綜合影響Fig.5 The comprehensive impact of sampling noises and interval on the inversion results

3 結語

污染物在實際河流的遷移轉化規律十分復雜,這對水質預測帶了很大的困難.本文考慮了污染物水質特征,給出適于可降解污染物質的修正的相關系數優化反演算法,使本身運用于地下水的方法向地表水的特定環境下開展實際運用推進了一大步.本方法更適合實際情景的應用.

假想算例計算表明,本方法的計算結果接近真實情況,并探討和分析論證其在實際運用時可能遇到的問題,具有一定的應用價值.

[1] 楊 潔,畢 軍,張海燕,等.中國環境污染事故發生與經濟發展的動態關系 [J]. 中國環境科學, 2010,30(4):571-576.

[2] 劉家琦.數學物理中的反問題及其數值解法 [J]. 科學探索1983,(3):75-83.

[3] 王彥飛.反演問題的計算方法及其應用 [M]. 北京:高等教育出版社, 2005:3-13.

[4] 劉曉東,姚 琪,薛紅琴.環境水力學反問題研究進展 [J]. 水科學進展, 2009,20(6):885-893.

[5] Atmadja J, Bagtzoglou. A C. State of the art report on mathematical method for groundwater pollution source identification [J]. Environmental Forensics, 2001,(2):205-214.

[6] 閔 濤,周孝德.污染物一維非恒定擴散逆過程反問題的數值求解 [J]. 西安理工大學學報, 2003,19(1):1-5.

[7] Kitanidis P K, Vomvoris E G. A geostatistical approach to inverse modeling in ground water modeling and one-dimensional simulations [J]. Water Resource Research, 1983,19(3):677-690.

[8] 姜繼平.水污染危險區域判定及其原型系統開發 [D]. 哈爾濱:哈爾濱工業大學, 2008.

[9] Boano F, Revelli R, Ridolfi L. Source identification in river pollution problems: A geostatistical approach [J]. Water Resource Research, 2005,41(7):W07023.doi:1029/2004WR003754.

[10] 張玉超,錢 新,錢 瑜,等.支持向量機在太湖葉綠素 a非線性反演中的應用 [J]. 中國環境科學, 2009,29(1):78-83.

[11] Cheng W P, Jia Y. Identification of contaminant point source in surface waters based on backward location probability density function method [J]. Advances in Water Resources, 2010,33(4): 397-410.

[12] 朱 嵩,劉國華,王立忠,等.水動力-水質耦合模型污染源識別的貝葉斯方法 [J]. 四川大學學報, 2009,41(5):30-35.

[13] 朱 嵩,劉國華,毛欣煒,等.基于貝葉斯推理的標準 k-ε湍流模型參數識別 [J]. 四川大學學報, 2010,42(4):78-82.

[14] 牟行洋.基于微分進化算法的污染物源項識別反問題研究 [J].水動力學研究與進展, 2011,26(1):24-30.

[15] Sidauruk P, Cheng A H D, Ouazar D. Groundwater contaminant source and transport parameter identification by correlation coefficient optimization [J]. Ground Water, 1998:36(2):208-214.

[16] 郭建青,李 彥,王洪勝,等.突然排污時計算河流水質參數的直線解析法 [J]. 水力發電學報, 2007,26(4):61-65.

[17] 鄭 彤.環境系統數學模型 [M]. 北京:化學工業出版社, 2003: 43-46.

[18] 謝永明.環境水質模型概論 [M]. 北京:中國科學技術出版社,1996:124-130.

Contaminant point source identification of rivers chemical spills based on correlation coefficients optimization method.

CHEN Yuan-hua, WANG Peng*, JIANG Ji-ping, GUO Liang (School of Municipal and Environmental Engineering, Harbin Institute of Technology, Harbin 150090, China). China Environmental Science, 2011,31(11)：1802～1807

A novel inversion algorithm based on an optimization approach for river point pollution sources was developed. Mass transport and kinetics processes of the contaminants in surface waters were combined along with the discharge history. And other relative parameters were deduced under the scenario that singular source instantly discharges degradable and soluble chemicals into one-dimensional rivers. A series of numerical experiments were carried out based on the hypothetic cases to analyze inversion effects associated with ambient river flow rates, contaminant decay rates, monitoring sites setting, sampling data errors and time intervals between two groups of sampling. When the monitoring time interval was less than 10 minutes and sampling data errors were controlled fewer than 5% approximately, the relative errors of pollution source location, total released mass and synthetical relative error are under 4%, 4% and 2%, respectively. Results show that parameters calculated fit well with the real values. In addition, the algorithms had the advantages such as efficient sampling process, minimum data requirement as well as easy programming. It was worthwhile to utilize this method for emergency environmental management practices.

point pollution source；river chemical spill；pollution source inversion；parameter identification；correlation coefficient optimization

X703.1

A

1000-6923(2011)11-1802-06

2011-02-01

國家自然科學基金資助項目(50821002)

* 責任作者, 教授, pwang73@hit.edu.cn

致謝：文章的英文部分由美國弗吉尼亞大學土木與環境工程系龍梧生教授修改與指導,在此表示感謝.

陳媛華(1986-),女,黑龍江哈爾濱人,哈爾濱工業大學市政環境工程學院碩士研究生,主要從事環境模擬與預警研究.發表論文1篇.