環(huán)境科學(xué)研究第13卷第1期esearch of environmental sciencesVol.13No.12000大氣化學(xué)模式的數(shù)值計(jì)算方法劉峻峰,李金龍,白郁華北京大學(xué)技術(shù)物理系北京100871)摘要∶大氣化學(xué)模式計(jì)算的一個(gè)重要研究熱點(diǎn)是如何快速、精確地對(duì)相互關(guān)聯(lián)的非線性常微分方程組(OEs進(jìn)行數(shù)值積分求解。由于大氣中督物種壽命長(zhǎng)短相差很大,導(dǎo)致ωDEs具有了很強(qiáng)的剛性因而在數(shù)值求算過程中存在著計(jì)算精度和計(jì)算效率的平衡問題。筆者介紹了當(dāng)前幾種常見的大氣化學(xué)模式的計(jì)算方法包括 GEAR QSSA( MQSSA)Y&BMY&B),RADM, LSODE SMVGEAR小參數(shù)法 Gong and Cho及 TWOSTEP等并結(jié)合當(dāng)前本領(lǐng)域內(nèi)的部分研究成果對(duì)一些重要的方法進(jìn)行了比較和分析。關(guān)鍵詞:大氣化學(xué)模式;剛性常微分方程組;數(shù)值積分中圖分類號(hào):X11文獻(xiàn)標(biāo)識(shí)碼:A文章編號(hào):1001-692×200001-0044-06The Numerical Solvers for Atmospheric Chemistry ModelsLIU Jun-feng, LI Jin-long, BAI Yu-huaDepartment of Technical Physics Peking University Beijing 100871)Abstract: The numerical integration of nonlinear systems of ordinary differential equations( ODEs )is a major problem to atmosphechemistry models. Because of the great disparities in the rates of the atmospheric chemistry reactions it is importantd a properway to accurately and efficiently solve the stiff systems of ODEs. In this paper several solvers including GEAR QSSA( MQSSA)Y&E MY&B) RADM SODE SMVGEAR SPM ong and Cho and TWOSTEP aare presented and compared the based on current research achievements in this fieldKey words: atmospheric chemistry model stiff ODE numerical integration自70年代以來,模式研究已成為大氣環(huán)境化學(xué)領(lǐng)域中或(1)代表了一組數(shù)量較大、非線性且互相關(guān)聯(lián)的常微不可缺少的組成部分。一—個(gè)性能優(yōu)異的大氣化學(xué)模式應(yīng)具分方程細(xì)ODEs)故無法求出其解析解,只能在給定每種化備2個(gè)必要條件⑩合理的模式機(jī)理。此機(jī)理能較準(zhǔn)確地反合物初值后用各物種濃度對(duì)時(shí)間數(shù)值積分的方法,求岀各映實(shí)際大氣中物種數(shù)量及彼此間的化學(xué)反應(yīng)同時(shí)給出較準(zhǔn)物種在特定時(shí)刻的濃度。確的動(dòng)力學(xué)數(shù)據(jù)。②高效精確的計(jì)算方法。該方法在滿足最早的方法是對(duì)具有初值的式(1)在固定時(shí)間步長(zhǎng)的條定誤差范圍的前提下應(yīng)盡可能快速地計(jì)算出最終結(jié)果件下進(jìn)行 Euler積分。但由于各反應(yīng)物種在大氣中的壽命長(zhǎng)根據(jù)研究的對(duì)象及目的不同,人們?cè)O(shè)計(jì)出了不同的機(jī)短不一,物種間的反應(yīng)速率相差高達(dá)幾個(gè)數(shù)量級(jí),造成了理。常見的幾種大氣化學(xué)模式機(jī)理及計(jì)算方法見表1。ODEs的剛性。為了保證求解方程組的準(zhǔn)確性和穩(wěn)定性由表1可以看出,般機(jī)理包含的化合物種類由幾十種須將時(shí)間步長(zhǎng)定得很小造成計(jì)算效率極低。因此在數(shù)值積到上百種反應(yīng)方程式由幾百個(gè)到上千個(gè)。由這些化合物及分過程中存在著積分時(shí)間步長(zhǎng)同各物種的壽命是否相匹配反應(yīng)方程式根據(jù)質(zhì)量作用定律可以寫出每種化合物的反的問題。若積分時(shí)間步長(zhǎng)取的較大,對(duì)于壽命短的化合物會(huì)應(yīng)速率方程產(chǎn)生較大誤差時(shí)間步長(zhǎng)取得過小,又會(huì)使計(jì)算過慢造成dh=F(C)=P(Ct)-L(CC(1)對(duì)計(jì)算資源的浪費(fèi)中國(guó)煤化工中新數(shù)值求解OEs的方法式中:P(Cd)第i種化合物的產(chǎn)生速率HCNMHG應(yīng)用,但是周看人們研究物L(fēng)(C,)第i種化合物的去除速率。質(zhì)量模式中所取的網(wǎng)絡(luò)數(shù)成干上萬,一般又要模擬較長(zhǎng)收稿日期:1999-09-14段時(shí)間的反應(yīng)因而對(duì)計(jì)算速度的要求也在提高,即需要作者簡(jiǎn)介劉峻峰(1974-)男北京市人碩士研究生種快速且較為精確的計(jì)算方法。GEAR法逐漸成為人們對(duì)新研究的快速方法精確度的檢驗(yàn)手段。1977年Y&B法161第1期劉峻峰等大氣化學(xué)模式的數(shù)值計(jì)算方法45和1978年QSsA"擬穩(wěn)態(tài)假設(shè))法的提出標(biāo)志著快速計(jì)算較18-201,并提出了一批新的或是改進(jìn)的方法,如法的誕生。此后直到現(xiàn)在又有許多科學(xué)家對(duì)大氣化學(xué)模式 LSODEL211,RADM21, Gong andTWOSTEP 241中剛性O(shè)DEs數(shù)值解法的效率和精度進(jìn)行了研究和比 SMVGEAR25] MQSSA和MY&B201表1幾種常見大氣化學(xué)模式及其計(jì)算方法Table 1 Atmospheric chemistry models and their numerical solvers機(jī)理縮寫名參考物種數(shù)反應(yīng)方程數(shù)計(jì)算法RADM2-IFUStockwell et al.[(1990)RADM2-FZKStockwell et al. [( 1990)RADMRADM2-KFAStockwell et al. [I(1990)Implicit EulaEuro-RADMStockwell and Kley?( 1994)CBM-ⅣLOTOGery et al. [3( 1989 )Builties+( 1992)63CBM-ⅣCBM-N-TNOGery, et al. 3(1989 Sterr'( 1994)GearC4.1Gery et al. [3(1989EMEPEMEPSimpson et al. [7( 1993), Simpson&( 1995140MCMPilling et al. [(1998)>2500>7000ADOM-ⅡlADOMLurmann et al. 10( 1986)47114Strand and HoyI( 1994)Andersson-Skold 12(1995RuhnkeRuhnke 13( 1995)145342資料主要來源 Kuhn et al.41998筆者在前面簡(jiǎn)要介紹了GEAR, LSODE, SMVGEARhF(Cm)Cm-1為前一步的C值,CQSSA, MQSSA,Y&B,MY&B,小參數(shù)法, RADM Gong and為本步的C值Cho及 TWOSTEP等計(jì)算方法。在下面結(jié)合各方面的研究對(duì)這些方法進(jìn)行了分類并對(duì)重要方法的精確性和計(jì)算效率n——計(jì)算步數(shù)做了比較和分析每步校正次數(shù)1方法介紹積分開始時(shí),Z0=(CC′o),C0為初值,C'o=dC/d1.1GEAR法和改進(jìn)GEAR法(C)第一步步長(zhǎng)取值較小P取1以后根據(jù)計(jì)算的誤GEAR法及與GEAR法類似的 LSODE, LSODES,差大小可按一定的方法自動(dòng)調(diào)整h和P。SMVGEAR等代表了一類較精確的數(shù)值計(jì)算方1.1.2 LSODE與 LSODES1.1.1GEAR法LSODE( Livermore solver for ordinary differential equaGEAR數(shù)值計(jì)算法151具有適用范圍廣、可以自動(dòng)調(diào)節(jié)ons521是建立在GEAR法基礎(chǔ)上,被廣為應(yīng)用計(jì)算ODEs時(shí)間步長(zhǎng)與誤差,以及計(jì)算結(jié)果準(zhǔn)確度高等特點(diǎn),是一種的較精確方法。此方法允許使用者根據(jù)ODEs的剛性程度為理想的計(jì)算方法。但其最大的缺點(diǎn)是計(jì)算步驟繁瑣要頻選擇Adam多步計(jì)算公式(對(duì)于非剛性體系)或BDF公式繁計(jì)算 Jacobian行列式,當(dāng)計(jì)算的物種比較多時(shí)會(huì)耗費(fèi)大 backward difference multistep formulas對(duì)于剛性體系)大量的機(jī)時(shí)和內(nèi)存。GEAR法的數(shù)值計(jì)算公式如下氣化學(xué)模式的ODEs體系剛性較大故多用BDF法。BD)F的a.預(yù)測(cè)公式計(jì)算公式為(2)C)+△tbf(Cb.校正公式式中:C"+1—當(dāng)前t"1的濃度Zm+1=Zm+1-[(l1-lnh)(xZm+1)](3)中國(guó)煤化工C.最終公式(4)CNMHG邇以刊酊炬陣結(jié)構(gòu)稀少的不同來源(是式中zn—(k)×(P+1矩陣來自于用戶直接輸入,還是來自于調(diào)用 Jacobian子程序或是A-(P+1)×(P+1)矩陣l≤P≤6多次調(diào)用方程評(píng)估子程序)分別采用不同方法簡(jiǎn)化計(jì)算步Ⅰ——單位矩陣驟提高計(jì)算效率26。(xm})k×1矩陣,G(Zm+1)=AF(Cm1.1. 3 SMVGEAR46環(huán)境科學(xué)研究第13卷SMVGEAR( Sparse-Matrix Vectorized Gear Code25是另(△tm≤Mt≤△tmx)(9)種GEAR式的解ODEs的方法。它與 LSODE有一定的相式中:最大允許△t的改變量取值為1%~100%。似性但也有其特點(diǎn)1.3Y&B和MY&B其最大的特點(diǎn)是利用了矢量式計(jì)算機(jī)在計(jì)算內(nèi)循環(huán)時(shí)1.3.1Y&B可以顯著提高浮點(diǎn)運(yùn)算效率這一現(xiàn)象設(shè)計(jì)計(jì)算代碼的。矢Y&B數(shù)值計(jì)算法又稱 Hybrid法l6,由 Young和Bori量化是指計(jì)算機(jī)在計(jì)算內(nèi)部DO…LOOP循環(huán)時(shí),可以提速于1977年提出其目的是建立一種先預(yù)測(cè)再迭代校正的積到其最大速度的90%以上。當(dāng)然若循環(huán)中包含一些特殊分方法。此法后被 McRae odman等采用并改進(jìn)。 Yamani作如遞歸、矩陣運(yùn)算時(shí))矢量化也會(huì)減慢因而此法也減1o等將其應(yīng)用于 CALGRID大尺度模式計(jì)算中271少了不必要的矩陣操作以提高速度Y&B法開始時(shí)按反應(yīng)物的特征時(shí)間與積分時(shí)間步長(zhǎng)的對(duì)于化學(xué)模式計(jì)算中的ODEs可以認(rèn)為是二維的即關(guān)系將反應(yīng)物分成2類:①剛性物質(zhì)(τ<△t),其在大氣中維是化學(xué)物種數(shù),一維是網(wǎng)格數(shù)。 Jacobson等研究認(rèn)為油于的壽命較短非剛性物質(zhì)(r;>△t)在大氣中的壽命較長(zhǎng)般情況下網(wǎng)格數(shù)會(huì)遠(yuǎn)遠(yuǎn)大于物種數(shù)故以網(wǎng)格數(shù)矢量化的公式如下速度比以化學(xué)物種數(shù)矢量化的速度快。計(jì)算每步的階和時(shí)a.剛性(x;≥△tn)間步長(zhǎng)的方法同于 LSODE。減少矩陣運(yùn)算的措施是引入預(yù)測(cè)值JSPARSE法。此法先對(duì)物種進(jìn)行排序,將出現(xiàn)頻率較少、和其他物種關(guān)系簡(jiǎn)單的物種排在矩陣的頂部而出現(xiàn)頻率C:= Ci+(P(10)多、關(guān)系復(fù)雜的物種排在矩陣的底部這樣處理可以最大限修正值度地提高矩陣和零相乘的幾率洏而后 JSPARSE將所有零和矩陣相乘的項(xiàng)刪除。這是其相對(duì)于 LSODE法的主要差別。(11)當(dāng)然有一點(diǎn)SMVGEAR運(yùn)行在非矢量化的計(jì)算b.非剛性(r<△tn)機(jī)上時(shí)和 LSODE法在性能上是相差不大的預(yù)測(cè)值1.2QSsA和 MQSSAC(2r"-4tn)+2P?r"△1.2.1 QSSAQSSA數(shù)值計(jì)算法17又稱擬穩(wěn)態(tài)近似法是 Hesstwedu修正值等于1978年提出的。其主要思想是按照積分時(shí)間步長(zhǎng)△t△t(P"+P;r"+r;y2與壽命τ的比將計(jì)算公式分成3種r+r;-△tna.Δt/τ>10,此時(shí)化合物的壽命很短按穩(wěn)態(tài)假設(shè)處(13)理即將P;和L;看成是常數(shù)1.3.2MY&BPMathur20等人指出倘若對(duì)n的物種應(yīng)用穩(wěn)態(tài)假(6)設(shè)則會(huì)大大提高計(jì)算效率。他們提出了將反應(yīng)物分成3類b.0.01<△t/r<10,此時(shí)化合物的壽命與所取步長(zhǎng)的方法即快速物質(zhì)τ;<0.2△n中速物質(zhì)0.2△tn≤r;接近5△tn和慢速物質(zhì)r2>5△tn。對(duì)于中速和慢速物質(zhì)應(yīng)用公式(7)(10~13)快速物質(zhì)則應(yīng)用公式(6) Mathur等同樣將自動(dòng)c,△t/x<0.01此時(shí)化合物的壽命較長(zhǎng)時(shí)間步長(zhǎng)控制機(jī)制引入到了Y&B法見式(9+1=CLn+△(P-LC)(8)14小參數(shù)法(SPM)公式(7是建立在P和L;是常數(shù)的基礎(chǔ)上的公式6廂(8)小參數(shù)法由我國(guó)韓天敏提出281適用于中等剛性的是為了增加計(jì)算速度。當(dāng)△t/r>10時(shí)用公式(6算所引ODEs。計(jì)算公式為起的誤差小于公式(7)所引起誤差的4.5×10-3%而公式M= p[ef(tn+I Ch+1)+Ch+I-C] (14)(8X簡(jiǎn)單 Euler積分公式)雖引入了誤差項(xiàng)(△t/r;)/2但引+1=Cn+(1-a)xn+ax+1(15)起的誤差小于公式(7)所引起誤差的5.0×103%。式中1.2. 2 MQSSAp=hh+∈)Mh2等引入矯正因子P和L:P=-2,中國(guó)煤化工L;=L;。每次計(jì)算出的結(jié)果用上式校正后將P和L;代入1.5CNMHG公式(6~8)中繼續(xù)迭代RADM法22是針對(duì)RADM( Regional acid depositon mod此外 Mathur等還引入了時(shí)間步長(zhǎng)控制機(jī)制:el)模式提出。其對(duì)18種自由基作了穩(wěn)態(tài)假設(shè),與QSSA法類似對(duì)時(shí)間步長(zhǎng)也規(guī)定了上下界。計(jì)算公式如下a.OH自由基濃度由解下式方程得到第1期劉峻峰等大氣化學(xué)模式的數(shù)值計(jì)算方法47al oh]+ b[ OH]+c=16)各不相同但最終目的都是為了提高計(jì)算速度和精確程度對(duì)于其他自由基應(yīng)用穩(wěn)態(tài)假設(shè)見公式13然而就象測(cè)不準(zhǔn)原理中描述的時(shí)間和空間不可能同時(shí)測(cè)準(zhǔn)C,H,O,與HO計(jì)算公式一樣剛性O(shè)DEs的數(shù)值解法也存在著計(jì)算速度和計(jì)算精度(L:)(17)間難以同時(shí)滿足的矛盾。計(jì)算速度的提高是建立在犧牲計(jì)算精度的基礎(chǔ)上的而高精度的計(jì)算方法必然伴隨著占有大d.其他物種量的計(jì)算機(jī)資源。在模式的實(shí)際應(yīng)用過程中有時(shí)對(duì)計(jì)算速C(1-Ltn2)+P△t度要求很高而計(jì)算精度只要滿足一定的誤差范圍即可故L"△tn/2(18)采用一定的方法簡(jiǎn)化計(jì)算過程可以大大提高計(jì)算速度1.6 Gong and Cho筆者介紹的方法基本上可以分成3類:第1類是GEARng and Cho法23將時(shí)間步長(zhǎng)定為0.5h根據(jù)物質(zhì)的類包括 GEAR LSODE LSODES) SMVGEAR等方法其基壽命將其分為2類:①r≥104s稱為慢物質(zhì)熜τ<104s稱本思路來源于GEAR法;第2類是穩(wěn)態(tài)假設(shè)類,此類方法的為快物質(zhì)。計(jì)算公式為基本思路是將物種按其壽命大間步長(zhǎng)之比分類,對(duì)不C=(C1C2)同類的物種分別采用不同的計(jì)算公式,包括QSSA(MQSSA),Y&BMY&B) RADM Gong and Cho等方法第3類是(19)d迭代法包括小參數(shù)法和 TWOSTEP法慢物質(zhì)在策(n+1)s濃度用式(20)算在各種數(shù)值計(jì)算方法中GEAR法是一種較為精確的計(jì)Cr= Ci +f(C )At(20)算方法也是研究其他快速方法的基礎(chǔ)。在研究模式機(jī)理階對(duì)于快物質(zhì)則用牛頓迭代公式計(jì)算出n-(n+1)s的濃度段對(duì)計(jì)箅結(jié)果的精確程度要求較高,般使用該方法這點(diǎn)增量δC可從表1看出。但是該法的計(jì)算速度較慢不適于大尺度的C2+-1+C2+f(Cm1+1,C2+-1)氣質(zhì)模式計(jì)算。 LSODE和 SMVGEAR兩者為了提高計(jì)算速(21)度都減少了計(jì)算 Jacobian矩陣的量,即只有在收斂失敗的情況下才計(jì)算。 SMVGEAR則又利用了矢量計(jì)算機(jī)的特點(diǎn)來其中C2=C21+8C2而送代計(jì)算出的C2又可以利用加快計(jì)算速度?？傮w而言第1類的3種方法基本上可以滿式(22增加C1的精度足較高的計(jì)算精度,其時(shí)間步長(zhǎng)可根據(jù)計(jì)算過程動(dòng)態(tài)改變。C+=C+(C1)(22)穩(wěn)態(tài)假設(shè)法代表了一類高速計(jì)算方法。由于大氣體系中物1.7 TWOSTEP種的壽命相差很大是造成ODEs剛性的來源故其核心內(nèi)容WOSTEP法2428是專為解大氣化學(xué)中的ODEs設(shè)計(jì)是對(duì)壽命很小和很大的物種采用較為簡(jiǎn)單的公式而對(duì)壽命的。其用 Gauss-Seidel迭代代替了牛頓迭代公式如下適中的物種進(jìn)行復(fù)雜計(jì)算這樣可以提高計(jì)算速度。但是此C+l= Y+ yf( tn+I Cn+I)(23)類方法的精度依賴于時(shí)間步長(zhǎng)取值大小,當(dāng)Δt的取值較大式中+5=時(shí)計(jì)算精度較差而Δt取值較小雖然提高了精度但計(jì)算c+2速度較低。QSSA法和Y&B法從形式上可以看到是比較相似的,尤其經(jīng)過改進(jìn)后的 MQSSA和MY&B法可近似認(rèn)為是一種方法。( ong and Chou法介于 LSODE法和Y&B法之間其對(duì)壽命較長(zhǎng)的物種采用一般的積分方法,而對(duì)壽命較利用 Gauss-Seidel迭代式(23)可寫成短物質(zhì)的計(jì)算方法類似于 LSODE法。RADM法建立在C=C)Y+mxtt11、(24)RADM模式的基礎(chǔ)上,它雖然針對(duì)特定模式提出但也有其普遍的適用意義特別是它提出了18種自由基的計(jì)算順序TWOSTEP的迭代計(jì)算公式為眾所周知在數(shù)值計(jì)算過程中不同的計(jì)算順序會(huì)產(chǎn)生不同的結(jié)果,該算法給出的合理計(jì)算順序?qū)笕说难芯繒?huì)帶來一定的幫助。第3類方法主要是用迭代法來求解,其中(k=12(25) TWOSTEP法假定了本步所用濃度是通過前2次的濃度計(jì)式中:m—物種數(shù)。算得汁汗田工內(nèi)竿剛性的體系當(dāng)然式(25)中物種的計(jì)算順序?qū)k的計(jì)算結(jié)果有影中國(guó)煤化工因而有必要對(duì)其進(jìn)行誤響但試驗(yàn)證明采用式25廂(26)對(duì)結(jié)果的影響是不大的。差和CNMHG交方法是計(jì)算各物種的均C=F(C1C2…Ck=12,Ck…,,Cm)(26)方根(RMS)對(duì)于第k種物質(zhì)其RMS定義為2各算法間的比較RISk(27)以上所列方法雖然在基本原理、公式結(jié)構(gòu)及條件參數(shù)上式中C",C"—分別代表精確值和計(jì)算值環(huán)境科學(xué)研究第13卷總步數(shù),m=12,,M[5] Stern R M. Entwichlung und Anwendung eines ells mit ver在滿足了一定的誤差范圍前提下看其對(duì)相同計(jì)算機(jī)所schiedenen chemischen Mechanismen[ A ] In: Meteorologische用時(shí)間來比較速度的高低。Say19等對(duì) LSODE Gong&Abhandlungen Vol. 8. Institut Fur Meteorologie der Freien UrCho,Y& B QSSA4種方法作了比較。其模式選用了CBMersitat Berlin[ C ] Verlag von Dietrich Reimer, Berlin, ISSNⅣ機(jī)理(包含33種化合物和81個(gè)反應(yīng)),對(duì)10幾種物質(zhì)采0342-4324,1994用不同方法時(shí)的RMS與CPU耗時(shí)情況進(jìn)行了比較?？蓺w[6] Milford J B Gao D Russell A G et al. Use of sensitivity analysis納如下:①隨著計(jì)算誤差的降低,各種方法的CPU耗時(shí)也在to compare chemical mechanisms for air-quality modeling J ]. Environ Sci Technol 1992 26: 1179-1 189成倍增加。② LSODE法的RMS隨著CPU耗時(shí)的增加降低[ 7] Simpson D Andersson-Skold Y Jenkin M E. Updating the chemi-最快而QSSA法降低最慢。③當(dāng)對(duì)CPU耗時(shí)要求較少時(shí)4cal scheme for the EMEP MSC-w oxidant model: current status種算法的RMS近似QSA法的性能較高。 ODMAN等研[R]. EMEP MSC-W Note 2/93, Norwegian Meteorological Insti究認(rèn)為使用Y&B和改進(jìn)后的QSA法可以使計(jì)算精度得到tute Oslo ,Norway 1993保證的前提下大幅提高計(jì)算速度。他們同樣采用CBM-Ⅳ[8] Simpson D. Biogenic emissions in Europe2: implication for ozone機(jī)理(一些條件不同于上者),以Y&B和改進(jìn)的QSSA法分control strategie J ] Journal of Geophysical Rescarch 1995 .100別與 LSODE法作了比較結(jié)果發(fā)現(xiàn)在與 LSODE的誤差控22891-22906在1%之內(nèi)時(shí),Y&B和QSSA法分別是 LSODE法的40與20[9] Pilling M Saunders s,Jenkin M,etal. Tropospheric Chemistry低計(jì)算過程中對(duì)QSSA法使用了質(zhì)量守恒的控制技術(shù),此Modelling[Eb/oLI.Downloadedfromhttp:/www.chem技術(shù)占時(shí)不大此外 MATHUR等與 JACOBSON等分別leeds. ac, uk/Atmospheric/MCM/ main. html. 1998對(duì) MQSSA MY& B RADM及 SMVGEAR LSODE等方法作[10] Lurmann F W, Lloyd A C ,Atkinson R, A chemical mechanism了比較結(jié)論較為相似。for use in long-range transport/ acid deposition computer modeling結(jié)論[J]. J Geophys Res 1986 91: 10905--1093611] Strand A Hov O. A two-dimensional global study of tropospheric通過對(duì)以上各種方法的介紹和討論可以看出大氣化學(xué)ozone production[ J]. j Geophys Res 1994 99 22877-2289模式中剛性O(shè)DFs的計(jì)算問題的解決方法是多種多樣的但[12] Andersson-Skold Y, Updationg the chemical scheme for the IVL大多都圍繞在計(jì)算精度和計(jì)算時(shí)間的取舍上建立的可分為tochemical trajectory model[ R]. IVL B-1151, IVL,Gote-2類,一類在可以接受的誤差范圍內(nèi)盡可能地提高計(jì)算速borg sweden 1995另一類則追求在可以忍受的計(jì)算速度下提高計(jì)算精度。顯[13] Ruhnke r. Ein Verfabren zur Analyse von chemischen Reak然隨著三維模式和大尺度模式的廣泛應(yīng)用,研究前者的人占tionssystemen und zur Erstellung von Chemienodulen fur atmo-多數(shù)spharische Modelld Ph-D]. Thesis University Essen germany在各種方法中 GEAR LSODE QSSA和Y&B法一直是人們研究的熱點(diǎn)。其中GEAR和 LSODE等方法代表了一類[14] Kuhn M Builtjes P J H Poppe D et al. Intercomparison for the精度較高的算法,是人們研究其他方法時(shí)精確度的參照法gas-phase chemistry in several chemistry and transport models而 LSODE因其在計(jì)算速度上的優(yōu)勢(shì)已被廣泛應(yīng)用。QSSA[J]. Atmos Environ 1998 32 593--709和Y&B法作為快速計(jì)算方法其計(jì)箅精度與系統(tǒng)穩(wěn)定性能15 Gear C W. Numerical Initial Value Problems in Ordinary Differ-也在不斷提高,些大型模式計(jì)算就采用了該方法。其他方ential Equation R ] Prentiedhill Englewood Cliffs , J, 197116] Young T R Boris J P. A numerical technique for solving stiff or法的建立原則也很明顯,一種是尋找新的原理,另一種是對(duì)dinary differential equations associated with the chemical kineticEAR與QSSA等已有方法進(jìn)行綜合改進(jìn)其中有許多地方blems[J]. J Phys Chem, 1977, 81: 2424值得后人借鑒。2427參考文獻(xiàn)[17] Hesstvedt E Hov O,Isaksen I S A. Quasi-steady-state approxi[ 1] Stockwell WR, Middleton P Chang J S et al. The second generamations in air pollution modeling: Comparison of two numericaltion regional acid deposition model chemical mechanism for regionalschemes for oxidant predictio[ J ] Int J Chem Kinet, 1978, 10air quality modeling J ]. J Geophys Res 1990 95: 16343-16367[2 Stockwell W R, Kley D. The Euro-RADM Mechanism: A Gas. 18 J Odman M T, Kumar N Russell A G. A comparison of fast chemPhase Chemical Mechanism for European Air Quality Studie[ C]eling[j]. Atmos Environ 1992Forschungszentrum Julich gmbH KFA) Julich Germany 1994中國(guó)煤化工[ 3] Gery M w Whitten G Z Killus J P et al. A photochemical kinet- [19]CN MH Gmparison of numerical methodsics mechanism for urban and regional scale computer modelinfor the integration of kinetic equations in atmospheric chemistry[J]. I Geophys Res19899412925-12956and transport models[ J ]. Atmos Environ 1995,29: 2585[ 4] Builtjes P J. The LOTOS- Long tem ozone simulation project summary report[ R ] IMW-TNO Report 92/240, TNO, Delft, The 20 Mathur R Young JO Schere K L , et al. A comparison of niNetherlands 1992merical techniques for solution of atmospheric equation[ J ].At第1期劉峻峰等大氣化學(xué)模式的數(shù)值計(jì)算方法49mos Environ 1998 32: 1535--1553torized gear code for atmospheric models[ J ]. Atmos Environ[21] Hindmarsh A C. LSODE and LSODI, Two new initial value or199428273-284dinary differential equation solver J]. ACM Newsl 1980, 15: 10 [26] Hindmarsh A C. ODEPACK: a systematized collection of ode22] Chang I S Brosr R A Isaksen I S A, et al. A three-dimensionalNorth-Holland Amsterdam 1983 55-64acid deposition model: physical concepts anion[J]. [27 Yamartino RJ Scire J S Carmichael GR et al. The CALGRIDGeophys Res98792:14681-14700rid model: I. Model formulation[J[ 23] Gong W Cho HR. A numerical scheme for the integration of the中國(guó)煤化工1512gas-phase chemical rate equations in three-dimensional atmo- [28 1CNMHG境化學(xué)M]北京高等教育spheric model![ J ] Atmos Environ 1993 27A 2147--2160出版社1990.156~158[24] VerwerJ(. Gauss-Seidel iteration for stiff ODEs from chemical[29]關(guān)治陸金甫.數(shù)值分析基砒[M].北京高等教育出版社kinetic J ]. SIAM J Sci Comput 1994, 15: 1243-12501998.445-45125] Jacobson M Z, Turco R P. SMVGEAR: A sparse-matrix ve