基于一維水質(zhì)模型沙潁河沈丘—界首河段水質(zhì)模擬分析

陳曉青　孫毅松　唐靜娟

摘 要：為深入分析沙潁河沈丘-界首斷面的水質(zhì)情況，根據(jù)河段2013年4月11日～2013年5月10日的水動(dòng)力數(shù)據(jù)和水質(zhì)監(jiān)測(cè)資料，建立了斯崔特-菲爾普斯BOD-DO模型，采用隱式差分法來(lái)模擬計(jì)算研究河段內(nèi)幾個(gè)斷面的水質(zhì)變化情況，并根據(jù)界首斷面的實(shí)測(cè)數(shù)據(jù)對(duì)模型的參數(shù)進(jìn)行調(diào)試率定，率定標(biāo)準(zhǔn)為界首斷面各計(jì)算時(shí)刻高錳酸鉀指數(shù)和溶解氧計(jì)算值與實(shí)測(cè)值的相對(duì)誤差之和最小，最后率定結(jié)果為K1=0.06d-1，K2=0.01d-1，并以此推求計(jì)算時(shí)間內(nèi)各斷面的水質(zhì)資料，對(duì)產(chǎn)生的誤差進(jìn)行了一定的分析。

關(guān)鍵詞：沙潁河水質(zhì)；一維水質(zhì)模型；托馬斯算法

1.1 一維水質(zhì)模型的求解

河流一維水質(zhì)模型一般需要用數(shù)值方法求解。數(shù)值法是隨著計(jì)算機(jī)的出現(xiàn)而迅速發(fā)展起來(lái)的一種近似計(jì)算方法，是用離散方法對(duì)數(shù)學(xué)模型進(jìn)行離散，進(jìn)而求出其數(shù)值解的方法。常用的數(shù)值法有有限差分法、有限單元法等。本文采用有限差分法中的隱式差分體系求解。

1.2 K1、K2率定過(guò)程

K1、K2的率定思路如下：首先，根據(jù)初始時(shí)刻（2013年4月11日）沈丘斷面和界首斷面的溶解氧和高錳酸鉀指數(shù)，線性插值分別求出初始時(shí)刻各計(jì)算斷面的溶解氧和高錳酸鉀指數(shù)；其次，根據(jù)提供的沈丘-界首水動(dòng)力資料，計(jì)算各時(shí)刻各斷面的離散系數(shù)E；最后，假設(shè)不同的K1、K2，由托馬斯法求解一維水質(zhì)模型，得到各時(shí)刻各斷面的計(jì)算溶解氧和高錳酸鉀指數(shù)，將之于界首斷面的實(shí)測(cè)值進(jìn)行誤差比較，選出誤差最小又符合實(shí)際意義的K1、K2的值，完成率定。

1.2.1 離散系數(shù)E

E可通過(guò)三種途徑計(jì)算：①由斷面流速分布資料推求。②由現(xiàn)場(chǎng)示蹤劑試驗(yàn)推求。③由經(jīng)驗(yàn)公式估算。在缺乏斷面流速分布資料和示蹤劑試驗(yàn)時(shí)，可用經(jīng)驗(yàn)公式估算。本文運(yùn)用費(fèi)希爾于1975年提出的公式：

式中：H為平均水深，m；J為水力坡降；u為斷面平均流速，m/s；B為河段平均水面寬，m；g重力加速度，m/s2 。

根據(jù)提供的沈丘-界首水動(dòng)力資料，計(jì)算各時(shí)刻各斷面的離散系數(shù)E。

1.2.2 K1 、K2率定結(jié)果

污染物降解系數(shù)K1是計(jì)算污染負(fù)荷的一項(xiàng)重要參數(shù)。與污水特性、PH值、水溫、水力特征等諸多因素有關(guān)。有機(jī)污染物質(zhì)在降解過(guò)程中消耗掉水中的溶解氧，同時(shí)水體又不停地從大氣中獲得氧。雖然水生生物的光合作用也是某些水體溶解氧的來(lái)源之一，但是對(duì)大多數(shù)水體而言，溶解氧的主要來(lái)源是通過(guò)水與空氣接觸界面上發(fā)生的氧傳遞，即大氣復(fù)氧。大氣降解系數(shù)K2即是表示這一物理過(guò)程的速度系數(shù)。

通過(guò)查閱文獻(xiàn)，得出K1、K2的取值范圍均為0-0.5d-1，對(duì)給定K1的每一個(gè)值求解在該范圍內(nèi)的K2的所有值（步長(zhǎng)均為0.01），如此循環(huán)直到取遍K1、K2在取值范圍內(nèi)的所有值。根據(jù)以上計(jì)算出的E和初始時(shí)刻各斷面的水質(zhì)資料，結(jié)合各時(shí)刻各斷面的水動(dòng)力資料以及沈丘斷面的各時(shí)刻水質(zhì)資料，即可求得假設(shè)K1、K2下各時(shí)刻各斷面的計(jì)算溶解氧和高錳酸鉀指數(shù)。將界首斷面的計(jì)算溶解氧和高錳酸鉀指數(shù)與界首斷面的各個(gè)時(shí)刻實(shí)測(cè)值對(duì)比，當(dāng)溶解氧和高錳酸鉀指數(shù)相對(duì)誤差之和最小時(shí)，此時(shí)的K1和K2即為最優(yōu)的率定結(jié)果。

2.1.溶解氧模擬結(jié)果分析

為分析誤差趨勢(shì)，以界首斷面計(jì)算時(shí)段內(nèi)溶解氧實(shí)測(cè)值與計(jì)算值的相對(duì)誤差為縱坐標(biāo)，計(jì)算時(shí)段為橫坐標(biāo)，繪制相對(duì)誤差曲線，如下圖1所示：

根據(jù)所繪制的圖像1可知，4月17日和4月24日溶解氧相對(duì)誤差異常高，分別達(dá)到80.8%、171.3%，所以在誤差分析小節(jié)內(nèi)需針對(duì)這些節(jié)點(diǎn)進(jìn)行誤差分析。但除此之外，其余溶解氧的計(jì)算濃度值與實(shí)測(cè)濃度值的相對(duì)誤差大部分都在40%的范圍內(nèi)，并且預(yù)測(cè)曲線的變化趨勢(shì)與實(shí)測(cè)的變化趨勢(shì)大致一致。說(shuō)明降解系數(shù)K1和復(fù)氧系數(shù)K2的率定基本滿足精度和曲線變化趨勢(shì)的要求。

2.2 誤差分析

由于本文采用K1、K2耦合率定，所以單獨(dú)分析高錳酸鉀指數(shù)或溶解氧的誤差結(jié)果都不客觀，應(yīng)該綜合考慮。從高錳酸鉀和溶解氧的模擬結(jié)果來(lái)看，凡是誤差較大的節(jié)點(diǎn)，均是實(shí)測(cè)值突變的節(jié)點(diǎn)，這是由于我們使用的差分方法不會(huì)突變，模型的靈敏度不足導(dǎo)致。如高錳酸鉀誤差較大的4月24日和5月2日，溶解氧誤差較大的4月17日和4月24日。其中4月24日數(shù)據(jù)不是很正常，高錳酸鉀指數(shù)7.95，溶解氧2.07，分別為計(jì)算時(shí)段內(nèi)的最大值和最小值，并且與4月24日前后的實(shí)測(cè)值也相差較大。而且2013年4月24日是周三，正好在一周五天工作日中間，又無(wú)特殊節(jié)日或發(fā)生特殊事件，所以這樣的突變存在觀測(cè)誤差的可能。

中可能會(huì)引起數(shù)值彌散和偽震蕩。如果在求解的過(guò)程中能夠考慮△t、△x、μ三者之間的關(guān)系，并且盡量保證 =1，那么計(jì)算結(jié)果與實(shí)測(cè)結(jié)果的擬合程度會(huì)更好。

參考文獻(xiàn)：

[1] 雒文生.水環(huán)境保護(hù)[M].中國(guó)水利水電出版社，2009

作者簡(jiǎn)介：

陳曉青（1996.03.29-）女，漢族，山東省樂(lè)陵市，身份證號(hào)：37148119960329004X，本科生，研究方向：水利水電工程

孫毅松（1996.07.26-）女，漢族，山東省鄒城市，身份證號(hào)： 370883199607260466，本科生，研究方向：水利水電工程

唐靜娟（1995.07.05-）女，漢族，河南省平頂山市，身份證號(hào)： 410481199507053521，本科生，研究方向：水利水電工程