文 秀 秀
(四川大學水力學與山區(qū)河流開發(fā)保護國家重點實驗室,成都 610065)
自20世紀80年代起,許多學者對河道突發(fā)事故污染源反問題進行了研究。金忠青等利用脈沖譜-優(yōu)化方法解決了一維以及二維情況的對流-擴散方程源項識別以及邊界條件控制反問題[3,4];WANGER B J[5]將水流和污染物輸運模擬與非線性極大似然估計相結合,識別處置污染源及過程參數(shù);閔濤[6,7]等給出了利用遺傳算法求解對流-擴散方程源項識別反問題的一種新方法;SINGH R M[8]等將遺傳算法與仿真模型相結合,并用于污染源識別;朱嵩[9]等利用貝葉斯推理建立了二維水源對流擴散方程參數(shù)估計的數(shù)學模型,通過貝葉斯定理,獲得了模型參數(shù)的后驗分布,從而獲得反問題的解。
本文通過寬淺河道瞬時點源二維水質模擬分析的數(shù)學模型的解析解的變形,在降解系數(shù)已知的情況下,根據(jù)最小二乘原理求出待求參數(shù)的綜合值,然后根據(jù)參數(shù)間的內在關系,聯(lián)立方程得出待求參數(shù)。這種方法易于理解,所需數(shù)據(jù)簡單,能夠同時反演得到河道離散系數(shù)、河道平均流速、污染源初始位置、污染源排放時間、污染源強度等參數(shù)。
由相關資料[10]可以查閱二唯瞬時點源排放的解析解為:
(1)
式中:C(x,y,t)代表在t=t0時刻,在平面二維坐標(x0,y0)處瞬時有質量為M的污染物進入水體,則因此產(chǎn)生的不同時間不同位置的濃度,kg/m3;ux和uy分別代表水流在縱向、橫向上的流速分量,m/min;M為瞬時點源派發(fā)的污染物數(shù)量,kg;K為污染物的降解系數(shù),min;Dx、Dy分別代表污染物的縱向、橫向離散系數(shù),m2/min。
我們假設事故是發(fā)生在水域為河流寬深比足夠大的所謂寬淺河道,并且假定河道之前未被污染,即污染物的背景濃度值為0。若初始監(jiān)測時刻不為0,設為t0,則解析解變?yōu)椋?/p>
C(x,y,t)=
(2)
對式(2)兩邊取對數(shù),因為實際中uy一般取0,所以ux為水流平均速度,對上式進行變形:
如今,小學數(shù)學的主要教學目的仍是提高升學率。任課教師也把教學的主要精力放在培養(yǎng)學生的解題能力,從而在考試中取得較好的成績。很多教師在課上忽視了數(shù)學文化的滲透,注重在有限的時間內給學生灌輸更多的知識,做更多的練習題。因此,加強小學課堂數(shù)學文化的滲透應該成為小學教學課堂的重要內容之一。我們對數(shù)學的教育應該為學生數(shù)學文化素養(yǎng)的養(yǎng)成提供有效途徑,培養(yǎng)學生解決問題的思維模式,為學生解決問題提供方法。
ax2+bx+cy2+dy+e
(3)
其中:
(4)
如果我們間隔一定時間對濃度進行再次采樣,可以得到另一組a、b、c、d、e這5個表征著未知參數(shù)綜合信息的擬合參數(shù)。設第1組擬合參數(shù)的下標為1,第2組擬合參數(shù)的下標為2。
我們根據(jù)2個不同時刻的a、b、c、d、e可以得到:
(5)
經(jīng)由以上分析可以發(fā)現(xiàn),對二維寬淺河道瞬時點源的參數(shù)估計問題就轉化為了求解2個不同時刻擬合參數(shù)a、b、c、d、e的問題。
根據(jù)最小二乘法進行參數(shù)擬合的步驟,首先要選出相應于式(3)的基函數(shù),本文中為線性多項式系:{1,x,y,x2,y2}。
若假設有n個不同位置的監(jiān)測點,那么根據(jù)這些監(jiān)測點的坐標,我們可以得到一個n行5列的矩陣,即矩陣A:
(6)
式中:x,y分別為這n個監(jiān)測點各自對應的橫縱坐標。
同時,為與式(3)中所作的對數(shù)變形相對應,對某一采樣時刻下這n個監(jiān)測點的濃度數(shù)據(jù)進行對數(shù)處理,得到矩陣B:
(7)
式中:C1,C2,C3,…,Cn為各個監(jiān)測點的濃度值。
根據(jù)最小二乘法進行曲線擬合的步驟,此時問題轉化為求解線性方程組:
(ATA)S=ATB
(8)
式中:S=(S1,S2,S3,S4,S5)T,為待求的擬合參數(shù)。
通過Doolitle分解將ATA分解為一個上三角矩陣L,一個下三角矩陣U,進行求解,結果{S1、S2、S3、S4、S5}分別為式(4)中的e、b、d、a、c的值。
假定有一寬淺河道研究區(qū)域,平均水深h=1.0 m,縱向流速ux=17 m/min,橫向流速uy=0,離散系數(shù)Dx=20 m2/min,Dy=30 m2/min。在x0=-2 100 m,y0=25 m處發(fā)生瞬時點源污染事故,在t0=-120 min時刻該點源釋放污染物總量M=5 000 kg,降解系數(shù)K=0.001 min-1。
由公式(2)計算得到第i個監(jiān)測點在j時刻的污染物濃度的理論值Cobs(xi,yi,tj)(i=1,2,…,n;j=1,2,…,k),為了使結果更貼近實際情況,考慮監(jiān)測誤差,對理論值附加一個高斯數(shù)據(jù)噪聲[11],即:
Cobs(xi,yi,tj)=Cobs(xi,yi,tj)+Ewi,jCobs(xi,yi,tj)
(9)
式中:wi,j為在(-1,1)區(qū)間內的標準正態(tài)分布隨機數(shù);E為監(jiān)測誤差水平。
在此算例中,我們選取2個不同時刻:t1=10 min,t2=20 min,每個時刻選取16個監(jiān)測點,即n=16,k=2,詳細監(jiān)測點坐標數(shù)據(jù)見表1。
表1 監(jiān)測點坐標 mTab.1 Monitoring point coordinates
選擇不同的監(jiān)測誤差水平,即不同的E值, 本次算例中取E為5%、10%、20%,得到的不同監(jiān)測誤差水平下參數(shù)反演結果見表2。
表2 不同監(jiān)測誤差下污染源的參數(shù)反演結果Tab.2 Parameter inversion results for pollution sources under different monitoring errors
通過不同監(jiān)測誤差下污染物參數(shù)反演的絕對相對誤差(Absolute Relative Error)ARE來衡量反演的結果,詳見圖1。相關計算公式如下:
(10)
式中:rE為不同監(jiān)測誤差水平下的參數(shù)反演值;r0為參數(shù)的給定值。
圖1 不同E值下各參數(shù)反演結果的ARE值Fig.1 ARE for the inversion result of each parameter under different monitoring error levels
由表2可以看出:監(jiān)測誤差與參數(shù)的反演結果總的趨勢表現(xiàn)為監(jiān)測誤差越大,參數(shù)反演值ARE就越大,也就是反演結果越不精確。但是即使監(jiān)測誤差達到20%,相對誤差也在10%以內。并且由圖1可以看出,在同樣的誤差水平下,Dy、M的反演結果最差,平均流速ux的反演結果最好。這是由于參數(shù)Dy、M是在其余參數(shù)的基礎上計算的,相對于其他的參數(shù)并不獨立,所以Dy、M的反演結果受到誤差累積的影響,相對誤差較大。
為了進一步對監(jiān)測誤差帶來的影響進行分析,利用公式(2)和給定的各參數(shù)取值得到解析解濃度值,并將其在平面坐標系中畫出,得到不同位置在t1=10 min時的濃度值。分別將E=5%、10%、20%時得到的反演參數(shù)帶入公式(2),得到相應的濃度值,并與參數(shù)給定值得到的濃度值做差,將結果在平面直角坐標系中表現(xiàn)出來,見圖2、圖3。
圖2 E=0時濃度分布Fig.2 Concentration distribution under the condition of E=0
圖3 不同E值下濃度值與E=0時濃度值差值分布Fig.3 Concentration difference distribution between different value of E and E=0
可以看出,不同E值下,根據(jù)反演參數(shù)求得的濃度值均與給定參數(shù)下求得的濃度值存在差距,且表現(xiàn)為濃度值大的地方偏差大,濃度值小的地方偏差??;并且由于參數(shù)反演存在誤差,濃度值偏差存在偏態(tài)分布的特性,而不是嚴格符合正態(tài)分布特性。
從理論上講,隨著監(jiān)測點數(shù)量的增加,反演結果也會更加精確。為了探究這一點,設計了一系列算例進行研究。算例中各參數(shù)給定值與上述條件相同,只是將監(jiān)測點數(shù)量分別減少到12組以及增加到24組,監(jiān)測誤差設定為E=5%。各參數(shù)的反演結果詳見表3。
采用上文提及的ARE(絕對相對誤差)對以上的反演結果進行分析,將分析比較結果直觀地顯示在圖4中。
由圖4反演結果對比分析可得:減小測點數(shù)量,反演參數(shù)的相對誤差較之前發(fā)生了較為明顯的增大,而增加測點數(shù)量,反演參數(shù)的精度略有提高,但是提高不明顯,這主要是由于本身使用測點數(shù)量已經(jīng)足夠反演濃度分布信息,因此再增加監(jiān)測點數(shù)目對反演結果的精確性影響不大,但是如果減少測點數(shù)目,就造成了測點數(shù)量較少,在具有監(jiān)測誤差時不能真實反映濃度分布信息。同時,我們通過橫向對比,不難發(fā)現(xiàn),依然是Dy、M的反演結果最差,平均流速 的反演結果最好??梢?,該方法在流速的反演上具有較高精度,而在Dy的反演上精度較低,但是也不超過10%的相對誤差。
表3 E=5%時不同測點數(shù)量下參數(shù)反演結果Tab.3 Parameter inversion results of different number of monitoring points under the condition of E=5%
圖5 不同測點數(shù)目下各參數(shù)反演結果的ARE比較Fig.5 ARE comparison of parameters inversion results under different number of monitoring points
通過以上的算例,可以得到:對于寬淺河道中的污染源參數(shù)反演問題,本文提出的基于最小二乘法的采用同時刻不同空間點的濃度數(shù)據(jù)的反演方法是可行的,該方法能夠精確地反演出各個源項參數(shù)(釋放時間、釋放位置、釋放量)以及水流參數(shù)(流速、縱向離散系數(shù));該方法對監(jiān)測誤差具有一定的敏感性,誤差越大,反演結果越差,但是即使監(jiān)測誤差為20%時,反演結果相對誤差也不超過10%,具有實際應用的價值;監(jiān)測點數(shù)量對參數(shù)反演的結果具有一定影響,在合理范圍內,隨著監(jiān)測點數(shù)量的增加,參數(shù)反演的精度也在提高,但是監(jiān)測點數(shù)量增加所帶來的參數(shù)反演結果的精確性提升是逐漸減弱的。本研究由于缺乏實測資料,故采用假想算例進行驗證,在下一步的工作中將搜集實測資料對方法進行驗證并作出改進。