李 平，馬 放

(哈爾濱工業(yè)大學城市水資源與水環(huán)境國家重點實驗室，150090哈爾濱)

污染物進入河流，在河流水動力的作用下發(fā)生稀釋、擴散、沉降、降解，這一過程的模擬和預測，需要根據(jù)河流水文特點、污染物的性質(zhì)建立相應水質(zhì)響應模型，來描述水體中污染物隨時間和空間遷移轉(zhuǎn)化規(guī)律［1－3］.水質(zhì)模擬過程是利用水體歷史數(shù)據(jù)、污染物的物理化學演變規(guī)律，建立相應的數(shù)學、物理、化學及實驗模型，進而描述水體中污染物質(zhì)與河流水質(zhì)之間的響應關(guān)系，對河流水環(huán)境系統(tǒng)進行模擬、預測和評價，水質(zhì)模擬是水環(huán)境監(jiān)督管理的重要工具［4-12］.水質(zhì)模型研究是從1925年 H.Streeter和E.Phelps提出的BOD和DO之間耦合關(guān)系一維數(shù)學方程開始，1970年，EPA開始推出QUAL水質(zhì)綜合模型，先后開發(fā)出了QUAL－Ⅱ、QUAL2E、QUAL2E－UNCAS、QUAL2K 等一維平移－彌散質(zhì)量遷移水質(zhì)模型［13］.1983年，US EPA 又開發(fā)了應用于河流、水庫、河口、海岸水質(zhì)分析模擬的通用WASP模型體系，可進行一維、二維、三維水質(zhì)模擬［14］.1993年，丹麥開發(fā)了 MIKE模型體系，包括多用于模擬一維河流水體的一維動態(tài)模型MIKE11、多用于具有垂向變化湖庫的二維動態(tài)模型 MIKE21 和三維空間模型 MIKE3［15-17］.水質(zhì)模型規(guī)格化、模塊化，集成化程度不斷加深拓展，但水質(zhì)模擬參數(shù)的研究是準確分析仿真的前提和必要基礎(chǔ)，是進一步降低水質(zhì)模型失真度、提高水質(zhì)模型預報精度的必要工作.

本研究利用水團追蹤及羅丹明B示蹤實驗獲得的數(shù)據(jù)，重點探討了佳木斯江段河流水質(zhì)擴散特性，研究了佳木斯江段平水期和冰封期的COD模型，識別了相應特征參數(shù)，驗證了模型的有效性，確定了松花江佳木斯段橫向、縱向擴散系數(shù)，為進一步精確可靠的流域水質(zhì)模擬研究提供參考.

1 實驗

1.1 水質(zhì)模型建立基本方法

建立水質(zhì)模型應首先根據(jù)研究對象的建模目的、研究范圍、河流特點，分析河段水質(zhì)模擬需求和水質(zhì)的時空分布規(guī)律，按照模擬需求進行數(shù)據(jù)收集和資料整理工作.結(jié)合河流水質(zhì)的特點，確定COD指標作為控制因子，結(jié)合河流水文水質(zhì)特點，確定相應水質(zhì)模型描述污染物輸入和水質(zhì)輸出的關(guān)系.模型建立后，利用相關(guān)水文和實驗數(shù)據(jù)進行計算機求解模型，確定水質(zhì)模型參數(shù)后對模擬的可靠性進行檢驗.

1.2 水團追蹤及羅丹明B示蹤實驗

淡水、河水、海水中的彌散僅從理論上分析較困難，采用現(xiàn)場示蹤試驗，選用易溶于水、稀溶液呈強熒光的工業(yè)染料堿性羅丹明B為示蹤劑，采用熒光分光光度法進行測定，確定縱向彌散系數(shù)，研究干流擴散規(guī)律.示蹤劑羅丹明B(分析純)瞬時投放，在河流的下游布置斷面采樣，樣品送回化驗室檢測，波長578 nm定量，水質(zhì)指標采用國家標準方法［15］測定.

利用人工投放的方式在河道中造成一個含有機物的污水團進行水團跟蹤監(jiān)測，利用污水團沿程濃度衰變情況推定河流對污染物質(zhì)的降解規(guī)律.BOD20為現(xiàn)場測定，枯水期、平水期在佳木斯江段前后兩次共進行水團追蹤作業(yè)，1月水溫為0℃，各斷面采樣時間為佳木斯下1:00，樺川6:00，新城14:30，富錦18:00.6月水溫為14℃，各斷面采樣時間為佳木斯下8:00，樺川19:00，新城8:00，富錦4:00.

2 結(jié)果與分析

2.1 河流水質(zhì)擴散研究

天然河流的寬深比一般較大，污染物在垂向一般比較均勻，可用二維擴散模型分析其擴散過程.若河道均勻平直，橫向流速可以忽略，并認為縱向流速和縱、橫方向的擴散系數(shù)皆為常數(shù)，則擴散模型基本方程為

在無限寬邊界中瞬時點源條件下的解析解為

其中:x、y為縱、橫距離(m);ρ為垂向平均質(zhì)量濃度(mg/m3);μ為流速(m/s);Dx、Dy為縱、橫向擴散系數(shù);t為時間;h為水深(m).

2.1.1 擴散系數(shù)估值

河流擴散系數(shù)取決于河流深度和摩擦速度，摩擦系數(shù)又與河深、河面坡度有關(guān)，另外河流寬度和河道河曲彎度對河流橫向擴散有影響.擴散系數(shù)可由下面公式求出:

式中:D為河流擴散系數(shù)(m2/s);h為河深(m);u*為摩擦速度(m/s);g為重力加速度(m/s2);s為河面坡降;a為經(jīng)驗常數(shù)，求橫向、縱向擴散系數(shù)時，其經(jīng)驗常數(shù)分別為0.23和5.96.

研究江段中深度變化較小，最深處約為5 m，淺處為2 m，平均水深為3.36 m.整個研究江段河底較為平坦，江寬平均為800～1 200 m，江段底質(zhì)為砂礫，江面坡度為0.015%.

1)信號重疊是植物提取物中復雜光譜的一個問題，妨礙了提取物的識別和量化，因此需要更好的信號分辨率。多維核磁共振波譜法的共振信號分別在多個獨立的頻率軸上傳播，從而減少了信號重疊，但多維核磁共振波譜法的局限是采集時間較長。近年來報道了一種稱作Hadamard譜的快速多維譜方法[49]，在相同的分辨率和信噪比的前提下，基于Hadamard的信號采集時間比傳統(tǒng)方法減少約20倍，因此Hadamard譜或許會成為新的研究趨勢。

根據(jù)式(2)求出整個研究河段摩擦速度平均值u*=0.070 3(m/s)，粗略估算研究江段的橫向擴散系數(shù)為 0.054 m2/s，縱向擴散系數(shù)為1.41 m2/s.由于估算未計入河寬、河曲以及流速變化的影響，一般較試驗求得數(shù)值要小1～2個數(shù)量級.

2.1.2 試驗江段水文特征

松花江佳木斯江段全長201 km，流經(jīng)佳木斯、樺川、綏濱、富錦、同江等市縣，并與黑龍江匯合于同江市.松花江下游地處三江平原，河道寬，分叉較多.歷史上徑流量年變化幅度較大，存在明顯的豐、枯周期.枯水年與豐水年徑流量比為1∶1～6(佳木斯斷面).試驗江段為松花江下游佳木斯－富錦江段，江段長167 km，河道彎曲系數(shù)為1.17.

根據(jù)26年水文資料統(tǒng)計，研究江段多年平均流量為2 150 m3/s，多年最大流量為18 400 m3/s是多年最小流量的100余倍，表明該江段流量變化很大.在染料釋放期間，即6月份流量為3 840 m3/s，研究江段平均流速為1.05 m/s.通過平水期在松花江佳木斯江段染料(羅丹明B)示蹤試驗，采用瞬時點源形式釋放，以人工一次性傾倒式投藥，距離排放點(高程76.2 m)約4 km處設監(jiān)測斷面(高程75.6 m)，對其濃度進行判斷，以多元回歸方法求出擴散基本方程解析式，進而確定相應擴散系數(shù).

2.1.3 擴散特征研究

染料自釋放點向下游縱向遷移的同時亦進行橫向擴散.據(jù)資料報道，染料混合后，不論橫向或縱向擴散，開始為偏態(tài)分布，以后接近正態(tài)分布.由圖1可以看出，位于接近河流中泓線的點染料質(zhì)量濃度分布也基本趨于正態(tài)，染料在斷面橫向的質(zhì)量濃度分布基本接近正態(tài)(圖2).

圖1 河流中泓線斷面質(zhì)量濃度分布

圖2 橫斷面質(zhì)量濃度分布

示蹤試驗數(shù)據(jù)顯示，示蹤劑質(zhì)量濃度分布基本呈正態(tài)分布，以擴散方程解析式來推求擴散系數(shù).將式(1)兩端同乘以t，然后分別取對數(shù)，可變化成下列形式:

2.1.4 擴散系數(shù)計算結(jié)果

研究江段的擴散系數(shù)按染料實驗求得，表1列出松花江佳木斯江段擴散系數(shù)計算結(jié)果.最終確定松花江佳木斯江段橫向擴散系數(shù)為0.41m2/s，縱向擴散系數(shù)為43.51m2/s，為進一步精確可靠的流域水質(zhì)模擬研究奠定基礎(chǔ).因橫向擴散系數(shù)較縱向擴散系數(shù)小兩個數(shù)量級，影響較小，本文只考慮水體縱向擴散作用.

表1 松花江佳木斯江段擴散系數(shù) m2·s-1

2.2 化學耗氧量模型

2.2.1 模型的選擇

河流水質(zhì)模擬計算中，綜合衰減系數(shù)起著非常重要的作用［18］.通過對松花江水系部分江段COD監(jiān)測，發(fā)現(xiàn)COD在最初階段質(zhì)量濃度下降較快，但隨著水團流動距離的增加，COD質(zhì)量濃度下降速度逐漸減緩，最終趨于穩(wěn)定.質(zhì)量濃度負荷變化規(guī)律符合負指數(shù)方程，選用一級反應速率方程來表征河段COD水質(zhì)變化特征.

微分方程形式為

此微分方程解為

式(4)和(5)為穩(wěn)態(tài)條件下基本模型，并未考慮河流彌散作用［19］.綜合河流的穩(wěn)態(tài)彌散作用，以對松花江佳木斯江段進行更精確的模擬，將式(4)右側(cè)增加相應的彌散作用表達式，即

式中:L為水質(zhì)質(zhì)量濃度(mg·L-1);L0為初始COD質(zhì)量濃度(mg·L-1);t為流程時間(s);K1為綜合衰減速率(d-1);D為擴散系數(shù)(m2/s)，初值條件L(x=0)=L0.

2.2.2 模型參數(shù)

1)參數(shù)識別方法.流速與縱向彌散系數(shù)是反映河流特征的重要參數(shù)，河流縱向彌散系數(shù)與流速之間的關(guān)系可以采用以下公式表達

式中:u*為摩擦速度(m/s);u為平均速度(m/s);b為水面寬(m);h為水深(m).

表2列出了佳木斯江段1月及6月兩次化學耗氧量測試結(jié)果.共計8組實測數(shù)據(jù)，5組用于參數(shù)識別，3組用于參數(shù)驗證.其中1、2號數(shù)據(jù)用于冰封期參數(shù)識別;5～7號數(shù)據(jù)用于平水期參數(shù)識別，3、4、8號數(shù)據(jù)用于參數(shù)驗證.

2)參數(shù)識別結(jié)果.將表3的冰封枯水期(水溫0℃)1～4號數(shù)據(jù)和平水期(水溫14℃)8～10號水團追蹤數(shù)據(jù)輸入式(7)，識別流速和擴散系數(shù)結(jié)果見表4.

表2 松花江佳木斯段COD實測結(jié)果

表3 模型參數(shù)識別輸入

表4 參數(shù)結(jié)果

將表4結(jié)果和冰封枯水期(水溫0℃)1～2號數(shù)據(jù)和平水期(水溫14℃)水團追蹤數(shù)據(jù)輸入，利用已知的x(距離)、質(zhì)量濃度L以及初始質(zhì)量濃度L0、彌散系數(shù)D、平均流速u，擬合式(6)計算識別化學耗氧量水質(zhì)模型參數(shù)K1，冰封期和平水期K1，結(jié)果見表5.

表5 參數(shù)識別結(jié)果

利用模型參數(shù)修正方程基本型式K(T)=K(20)·θT-20，對計算所得衰減系數(shù)進行溫度修正.式中θ為溫度修正系數(shù)，K(T)為任意溫度(T)下的參數(shù);K(20)為20℃ 下的參數(shù).根據(jù)表5的參數(shù)識別結(jié)果進行溫度修正，得到水溫20℃ 下衰減系數(shù)K為0.188 d，溫度系數(shù)θ為1.109 8，任意溫度下衰減系數(shù)經(jīng)驗式為

3)參數(shù)驗證結(jié)果.由表4結(jié)果可知，冰封期流速不足平水期流速1/2，冰封期擴散系數(shù)更遠較平水期擴散系數(shù)小，水溫直接導致冰封期COD衰減系數(shù)遠小于平水期COD衰減系數(shù).分別利用冰封期3、4號和平水期8號數(shù)據(jù)進行參數(shù)驗證，相應參數(shù)輸入見表6，驗證結(jié)果如表7所示.

冰封期流速0.41m/s、擴散系數(shù)6.65m2/s、平水期流速1.05m/s、擴散系數(shù)43.51m2/s條件下，通過模型輸入計算，最小相對誤差0.31%，最大相對誤差8.55%，實測COD和模擬值的相對誤差控制在10%以內(nèi)，滿足《水污染物排放許可證技術(shù)報告編寫指南》對模型誤差應控制在±40%精度要求［20］.所獲取的水質(zhì)模型參數(shù)值模擬精度較高，應用在松花江佳木斯段水質(zhì)模擬中較合適.文中所獲得的擴散系數(shù)和降解系數(shù)將為進一步提高河段水質(zhì)模擬精度和COD總量控制提供數(shù)據(jù)和理論支撐.

表6 模型參數(shù)驗證輸入

表7 松花江化學耗氧量模型驗證結(jié)果

3 結(jié)論

1)通過羅丹明B示蹤實驗分析，確定了松花江佳木斯段橫向擴散系數(shù)為0.41 m2/s，縱向擴散系數(shù)為43.51 m2/s，比估算值高一個數(shù)量級，為進一步精確可靠的流域水質(zhì)模擬研究奠定基礎(chǔ).

2)確定了松花江佳木斯段COD衰減系數(shù)K(20)=0.187，計算模擬結(jié)果和實測COD相對誤差都在10%以內(nèi)，表明該模型應用在松花江佳木斯段水質(zhì)模擬中較合適.

3)對松花江佳木斯段COD模型參數(shù)進行了溫度參數(shù)修正，得到任意溫度下衰減系數(shù)經(jīng)驗式K(T)=0.188 ×1.109 8T-20.

［1］李蘭.水質(zhì)多參數(shù)辨識與反演算法［J］.水利學報，1998(6):46－52.

［2］MACAVOY D C，MASSCHELEYN P，PENG C，et al.Risk assessment approach for untreated wastewater using the QUAL2E water quality model［J］.Chemosphere，2003，52:55－66.

［3］郭德偉，韓光明，田秀麗，等.復雜邊界河流水質(zhì)預測計算［J］.山東農(nóng)業(yè)大學學報，2008，39(4):553－556.

［4］劉偉，劉洪超，徐海巖.基于MIKE11模型計算河流水功能區(qū)納污能力方法［J］.東北水利水電，2009，8:69－72.

［5］金春久，王超，范曉娜，等.松花江干流水質(zhì)模型在流域水資源保護管理中的應用［J］.水利學報，2010，41(1):86－92.

［6］LEE J H W，JAYAWARDENA A W，WANG Z Y.Proceeding of international conference on environmental hydraulics［M］.Northland:A A Balkman，1999:387 －392.

［7］張波，王橋，李順，等.基于系統(tǒng)動力學模型的松花江水污染事故水質(zhì)模擬［J］.中國環(huán)境科學，2007，27(6):811－815.

［8］張荔，王曉昌.水文水質(zhì)綜合模型與城市小流域河流污染分析［J］.西安建筑科技大學學報，2007，39(4):569－574.

［9］PORTER A L，F(xiàn)ITTIPALDI J J.Environmental methods review:retooling impact assessment for the new century［M］.Fargo，North Dakota，USA:The Press Club，1998:100-106.

［10］MARSILI-LIBELLI S，GIUSTIA E.Water quality modelling for smallriver basins［J］.Environmental Modeling＆ Software，2008，23(4):451－463.

［11］MARTINEZ C J，WISE W R.Analysis of constructed treatment wetland hydraulics with the transient storage model OTIS［J］.Ecological Engineering，2003，20:211 －222.

［12］付彩菊.河流水污染定量研究進展［J］.環(huán)境科學與管理，2007，32(9):59－61.

［13］曹永中，周孝德，吳秋平，等.河流水質(zhì)模型研究概述［J］.水利科技與經(jīng)濟，2008，14(3):197 －200.

［14］TORO di D M，SIFITZPATRICK J J.Documentation for water quality analysis simulation program(WASP)and model verification program(MVP)［M］.Duluth，MN:US Environmental Protection Agency，1983.

［15］Danish Hydraulics Institute.MIKE11:user guide and reference manual［R］.Horsholm，Denmark:［s.n.］，1993.

［16］Danish Hydraulics Institute.MIKE21:user guide and reference manual［R］.Horsholm，Denmark:［s.n.］，1996.

［17］Danish Hydraulics Institute.MIKE3:eutrophication module，user guide and reference manual［R］.Horsholm，Denmark:［s.n.］，1996.

［18］馬常仁，徐得潛，周慧，等.基于最小二乘法與最速下降法的綜合衰減系數(shù)的率定［J］.水利科技與經(jīng)濟，2012，18(1):20 －24.

［19］張解放，史忠科.河流水質(zhì)分析的狀態(tài)方法［J］.重慶環(huán)境科學，2002，24(2):40 －42.

［20］操群，柳炳祥，林洋.BOD－DO耦合模型在鄱陽湖入湖河流水質(zhì)預測中的應用［J］.科技資訊，2011，7:132-133.