不同時空格式在求解污染物對流擴散方程中的應用

劉忠波，房克照，孫昭晨

（大連理工大學海岸和近海工程國家重點試驗室，遼寧 大連 116023）

不同時空格式在求解污染物對流擴散方程中的應用

劉忠波，房克照，孫昭晨

（大連理工大學海岸和近海工程國家重點試驗室，遼寧 大連 116023）

為了研究污染物對流擴散方程中不同時空格式的適用性，針對對流擴散方程的一維﹑二維和三維3種情況，分別建立了預報-校正的有限差分數值模型。在時間步進格式上分別采用了Crank-Nicolson格式或混合4階Adams-Bashforth-Moulton格式，對對流項分別采用2階精度或4階精度，對擴散項采用了2階精度。利用建立的數值模型求解了經典的污染物濃度場對流擴散，通過數值解與解析解的比較討論了不同時空格式對數值模型計算結果的影響。結果表明：對空間一次導數采用4階精度可以避免采用2階精度帶來的誤差。采用混合4階Adams-Bashforth-Moulton格式或Crank-Nicolson格式數值計算結果均與解析解吻合程度較好，但對于數組為[40,40,40]的三維對流擴散問題，前者比后者省時20.7%。

對流擴散方程；Crank-Nicolson格式；Adams-Bashforth-Moulton格式；4階精度

進入21世紀后，伴隨現代工業的發展,越來越多的污染物被排放到大氣、河流、湖泊和海洋環境中,從而使所涉及區域的環境質量遭受不同程度惡化,這已嚴重影響了工農業的可持續發展，為此有必要展開污染物的運動規律及其歸宿的研究。

描述污染物在海洋環境中運動過程的數學模型是對流擴散方程，對簡單問題有解析解，但對于復雜問題只能借助于數值方法獲得相應的數值解，國內外許多學者為解決該問題做出了大量的工作[1-5]。鄭永紅等利用SOWMAC格式求解了二維對流擴散方程[4]，汪守東和沈永明將SOWMAC格式、Crank-Nicloson有限元格式等應用到了三維對流擴散方程中[5]，并指出迎風等多種格式在求解三維對流擴散問題時存在不足。為研究不同時空格式在求解對流擴散方程中的適用性，本文采用了Crank-Nicolson格式和混合4階Adams-Bashforth-Moulton格式兩種時間步進格式和兩種空間精度格式（2點2階精度和4點4階精度）處理對流項，并利用所建立的數值模型模擬了經典污染物濃度場的對流擴散問題。

1 數學模型

1.1 三維對流擴散方程

式中：c為污染物濃度；u，v，w 為 x，y，z方向的水流速度；Dx，Dy，Dz為擴散系數；Sm為源、匯項；t為時間。

1.2 模型求解

式（1）可寫成如下形式：

F中的對流項采用4點中心差分格式或2點中心格式（對于中間點）和采用2階精度（對應于邊界點），擴散項采用3點2階中心格式，這里僅給出關于x求導的表達式：

針對式（2）的求解，可選擇的方法很多，這里僅選用劉忠波等在計算高階Boussinesq水波方程時采用的Crank-Nicolson 格式或混合 4 階 Adams-Bashforth-Moulton 格式［6］，這兩種格式的具體表達形式如下：

（1）Crank-Nicolson格式（C-N格式）

該格式預報和校正格式的構造如下：

式中：n+1表示當前要計算的時間層；n表示上一個時間層。

（2）混合4階Adams-Bashforth-Moulton(ABM)格式

預報時采用3階時間格式，校正時采用4階格式，其表達如下：

無論是采用C-N格式或ABM格式，均需多次迭代過程。當兩次計算濃度誤差小于給定的數（本文取10-8），迭代過程結束，進入下一步的計算，若不滿足以上條件，繼續回到校正步進行重新校正，直到滿足條件。

2 算例與討論

以單位高斯脈沖在一維﹑二維和三維空間內的對流擴散為研究對象，該問題在三維空間下的初始條件和對應解析解的具體表達式如下：

式中：

對于一維或二維情況，可根據上式簡化得到，這里就不再贅述。

2.1 對流項的空間精度對計算結果的影響

考慮一維情況，初始投放點x0=5 m，水流速度u=1 m/s，Dx=0.01 m2/s。數值計算時，計算域0 m≤x≤45 m，空間步長△x=0.05 m，時間步長△t=0.005 s，在計算中時間步進上則采用ABM格式，對對流項采用2階精度或4階精度，討論了精度對數值計算結果的影響。圖1給出了t=5 s時刻的數值解與解析解結果的比較。

考慮二維情況，初始投放點(x0,y0)=(0.5，0.5)m，水流速度u=0.8 m/s，v=0.8 m/s，Dx=0.01 m2/s，Dy=0.01 m2/s。數值計算時，計算域0 m≤x≤2 m，0 m≤y≤2 m，空間步長△x=△y=0.05 m，時間步長△t=0.005 s，圖2分別給出了t=1 s時刻的x=y斷面上數值解與解析解結果的比較。

考慮三維情況，初始投放點(x0，y0，z0)=（0.5，0.5，0.5）m，水 流 速 度 u=0.8 m/s，v=0.8 m/s，w=0.1 m/s，Dx=0.01 m2/s，Dy=0.01 m2/s，Dz=0.01 m2/s。數值計算時，計算域 0 m≤x≤2 m，0 m≤y≤2 m，0 m≤z≤2 m，空間步長△x=△y=△z=0.05 m，時間步長△t=0.005 s。圖3分別給出了t=1.25 s時刻z=0.2 m水深處x=y上數值解與解析解結果的比較。圖4給出了t=1.25 s時刻z=0.2 m水深處全場濃度分布圖。

圖1 一維對流擴散方程計算結果與解析解的比較

圖2 二維對流擴散方程計算結果與解析解的比較

圖3 三維對流擴散方程計算結果與解析解的比較（z=0.2 m水深x=y線上）

由圖1～圖3可見，采用4階精度計算結果與解析解結果吻合程度很好，而采用2階精度則明顯與解析解結果存在差異，這種差異主要表現在兩個方面：（1）最大濃度峰值均比解析解最大濃度峰值有所減少；（2）最大濃度峰值的位置落后。以上兩點集中反映出采用2階精度求解對流項是不足的，而應采用4階精度。當進一步采用7點6階精度來求解空間一次導數時，我們發現數值計算結果與采用4階精度的計算結果差別不大，這說明提高差分精度格式達到6階精度已無必要。

圖4 z=0.2 m處三維污染物濃度計算結果與解析解的比較（右圖為解析解）

2.2 時間步進格式對計算結果的影響

利用C-N格式和ABM格式對時間導數項差分處理，針對以上算例進行計算，這里僅給出一維情況下的對比，具體對比情況見圖5。由圖可見，兩種格式計算結果幾乎重合。為此，我們考察了兩種格式的計算效率，以上面給出的三維情況為例，在現有電腦配置為3.4雙核CPU奔4系列，內存1G條件下，采用C-N格式計算 [40,40,40]空間下計算250步（1.25/0.005=250）需要時間為29 s，而采用ABM格式所需時間為23 s，省時20.7%。因此對于污染物對流擴散方程的求解我們推薦采用ABM格式。

3 結論

本文將C-N格式和混合4階ABM格式的預報-校正的有限差分法應用到求解污染物對流擴散方程，通過研究得到以下結論：

（1）C-N格式和混合4階ABM格式均能給出較好的預報結果，并與解析解吻合程度良好，但ABM格式在計算時間上存在較大優勢，因此在此類問題求解中推薦混合4階ABM格式。

（2）對于對流項，建議采用4階中心格式，從而可以避免2階中心格式帶來的較大誤差。

本文給出的兩種格式編程簡便，數值精度也較高，可期用于模擬實際復雜多變的海洋環境中污染物運動，這將在下一步進行深入研究。

[1]Komatsu T,Ohgushi K,Asai K.Refined numerical scheme for advective transport in diffusion simulation[J].Journal of Hydraulic Engineering,1997,123(1):41-50.

[2]Sankaranayanan S,Shankar N J,CheongH F.Three-dimensional finite difference model for transport of conservative pollutants[J].Ocean Engineering,1998,25(6):425-442.

[3]梁書秀.潮汐水域中污染物輸移擴散的數值模擬研究[D].大連:大連理工大學，2000.

[4]鄭永紅,沈永明,王利生,等.污染物移流擴散方程的高精度分裂格式[J].水利學報,2002(2):41－46.

[5]汪守東,沈永明.三維對流擴散方程的三種高精度分裂格式[J].應用數學和力學，2005，26(8):921-928.

[6]劉忠波，孫昭晨，鄒志利.不同時空格式在求解Boussinesq水波模型中的應用[J].港工技術，2009,46(1):4-8.

The Application of Different Time&Space Schemes in Pollutant Convectivediffusive Equation

LIU Zhong-bo,FANG Ke-zhao,SUN Zhao-chen
（State Key Laboratory of Coastal and Offshore Engineering,Dalian University of Technology,Dalian Liaoning 116023,China）

In order to understand the application of the different time and space derivative schemes in pollutant convective-diffusive equation,one-dimension,two-dimension and three-dimension models were established based on predict-correct finite difference method.In these models,time marching schemes included Crank-Nicolson scheme and a composite fourth order Adams-Bashforth-Moulton scheme.Second order accuracy or fourth order accuracy in space derivatives for convective terms and second order accuracy for diffusive terms were considered.Numerical simulations were carried out upon a classical pollutant concentration problem with these models,through the comparisons among the numerical results and the analytical solution,the effects of different schemes of time and space derivative were investigated.The results show that:Fourth order accuracy in space derivatives for convective terms can better simulate pollutant convective-diffusive problem,while second order accuracy for convective terms can not well simulate this problem.The computed results by a composite fourth order Adams-Bashforth-Moulton scheme or Crank-Nicloson scheme can better agree with the analytical solution,but for a three-dimension[40,40,40]space array,the first scheme can save 20.7%time.

convective-diffusive equation;Crank-Nicolson scheme;Adams-Bashforth-Moulton scheme;fourth-order accuracy

X143

A

1003-2029(2012)01-0096-04

2011-07-22

國家自然科學基金資助項目（51009018）；大連理工大學海岸和近海工程國家重點實驗室開放基金資助項目；河海大學海岸災害及防護教育部重點實驗室開放基金資助項目

劉忠波(1976-)，博士后，從事海岸波浪理論與數值模擬研究。