陳銘儒，胡軍國(guó)，*，崔武峰，俞 平

(1.浙江農(nóng)林大學(xué) 信息工程學(xué)院，浙江 杭州 311300； 2.中共杭州市臨安區(qū)委統(tǒng)一戰(zhàn)線(xiàn)工作部，浙江 杭州 311300)

土壤是生態(tài)環(huán)境中最大的碳庫(kù)，其總儲(chǔ)量約為1 394 Pg，每年向大氣中排放的二氧化碳量在68～98 Pg[1]；因此，土壤中二氧化碳濃度的細(xì)微變化都會(huì)對(duì)生態(tài)環(huán)境產(chǎn)生巨大的影響。研究土壤呼吸作用排放二氧化碳的變化規(guī)律對(duì)于緩解全球氣候變暖具有關(guān)鍵作用[2]。在這一背景下，準(zhǔn)確測(cè)量土壤呼吸排放的二氧化碳通量成為了生態(tài)學(xué)研究的一項(xiàng)重要課題。

目前，土壤呼吸碳通量的監(jiān)測(cè)方法主要包括微氣象法[3]和氣室法[4]2大類(lèi)：微氣象法不能夠測(cè)量小范圍的碳通量，且儀器成本高昂；氣室法成本較低、便于監(jiān)測(cè)，廣泛應(yīng)用于土壤呼吸研究。氣室法又可細(xì)分為動(dòng)態(tài)氣室法和靜態(tài)氣室法。其中，靜態(tài)氣室法利用堿液吸收二氧化碳而產(chǎn)生碳酸氫鈉的原理，將密閉氣室置于土壤上方，并在其內(nèi)部放入氫氧化鈉溶液，通過(guò)測(cè)定單位時(shí)間氣室內(nèi)部碳酸氫鈉的產(chǎn)生量而計(jì)算出二氧化碳通量，此方法為最原始的二氧化碳通量監(jiān)測(cè)方法[5-6]，測(cè)量過(guò)程復(fù)雜緩慢，而且計(jì)算誤差較大[7]，現(xiàn)已不再使用。動(dòng)態(tài)氣室法的原理是利用高精度氣體傳感器來(lái)采集氣室內(nèi)部的二氧化碳濃度，通過(guò)得到單位時(shí)間內(nèi)的積累量或氣室內(nèi)的流量差，從而計(jì)算出二氧化碳的通量值，是目前碳匯研究中使用最多的方法[8]。但動(dòng)態(tài)氣室法的缺點(diǎn)也較為明顯：其一，氣室在密閉期間會(huì)使二氧化碳在其內(nèi)部堆積，從而抑制了土壤向大氣中排放二氧化碳[9]：其二，弧形氣室會(huì)干擾二氧化碳向上擴(kuò)散時(shí)的湍流作用[10]，因此會(huì)對(duì)擴(kuò)散規(guī)律產(chǎn)生影響。這就使得動(dòng)態(tài)氣室法無(wú)法對(duì)二氧化碳擴(kuò)散通量進(jìn)行精確的測(cè)量。針對(duì)上述問(wèn)題，本研究提出一種使用Maxwell-Stefan擴(kuò)散模型，通過(guò)開(kāi)放型氣室計(jì)算土壤呼吸時(shí)排放二氧化碳通量的方式，旨在實(shí)現(xiàn)對(duì)土壤碳通量的長(zhǎng)時(shí)間精確監(jiān)測(cè)。

Maxwell-Stefan擴(kuò)散模型由英國(guó)物理學(xué)家、數(shù)學(xué)家James Clerk Maxwell和奧地利物理學(xué)家Josef Stefan共同提出。在此模型提出之前，F(xiàn)ick擴(kuò)散模型是計(jì)算流體通量最主要的方法，并已得到大量應(yīng)用[11-13]。Maxwell-Stefan擴(kuò)散模型在Fick擴(kuò)散模型的基礎(chǔ)上對(duì)擴(kuò)散系數(shù)進(jìn)行了改進(jìn)和擴(kuò)展，使之成為適用于計(jì)算多組分混合物擴(kuò)散通量的數(shù)學(xué)模型。1962年，研究者Duncan和Toor在進(jìn)行氮?dú)?、氫氣和二氧化碳?個(gè)儲(chǔ)器介質(zhì)之間擴(kuò)散的實(shí)驗(yàn)時(shí)發(fā)現(xiàn)了氮?dú)饽嬗跐舛忍荻冗M(jìn)行擴(kuò)散的現(xiàn)象[14]，從而證明了多組分系統(tǒng)和二元擴(kuò)散系統(tǒng)的擴(kuò)散規(guī)律是不同的。在此之后，多項(xiàng)研究證實(shí)Maxwell-Stefan擴(kuò)散模型在多組分?jǐn)U散領(lǐng)域可以取得很好的計(jì)算效果[15-17]。使用該模型計(jì)算土壤呼吸碳通量主要有以下幾點(diǎn)優(yōu)勢(shì)：第一，Maxwell-Stefan擴(kuò)散模型的擴(kuò)散系數(shù)使用矩陣來(lái)表示，矩陣中的元素反映了分子間作用力對(duì)其中2種組分?jǐn)U散系數(shù)的影響，該模型能夠很好地描述多組分?jǐn)U散系統(tǒng)的狀態(tài)；第二，該擴(kuò)散模型非常適用于氣相中的擴(kuò)散過(guò)程，并在此領(lǐng)域已得到了廣泛的使用，因此利用該模型來(lái)計(jì)算土壤呼吸碳通量更為適合；第三，Maxwell-Stefan擴(kuò)散模型考慮到了分子的熱力學(xué)運(yùn)動(dòng)，因此在對(duì)高溫高壓下的物質(zhì)擴(kuò)散進(jìn)行計(jì)算時(shí)同樣能夠取得較好的計(jì)算結(jié)果。Moldrup等[18]運(yùn)用Maxwell-Stefan擴(kuò)散模型對(duì)土壤氣體的擴(kuò)散活動(dòng)進(jìn)行了仿真研究，并且取得了很好的效果，更加說(shuō)明該方法適用于土壤碳通量研究。

本研究使用Maxwell-Stefan擴(kuò)散模型計(jì)算土壤呼吸排放二氧化碳的通量，主要目的如下：其一，改進(jìn)氣室法存在的抑制土壤呼吸的問(wèn)題，從而實(shí)現(xiàn)長(zhǎng)時(shí)間的監(jiān)測(cè)；其二，改進(jìn)Fick擴(kuò)散模型不能精確計(jì)算多組分?jǐn)U散系統(tǒng)通量的問(wèn)題，以此取得更高的二氧化碳擴(kuò)散通量精度；其三，為基于開(kāi)放氣室的土壤呼吸監(jiān)測(cè)裝置研發(fā)提供理論基礎(chǔ)。

1 材料與方法

1.1 半經(jīng)驗(yàn)擴(kuò)散系數(shù)

擴(kuò)散系數(shù)是描述流體擴(kuò)散能力的物理量，它表示某組分沿其擴(kuò)散方向，在單位時(shí)間和單位濃度條件下，垂直通過(guò)單位面積的質(zhì)量或物質(zhì)的量。在Maxwell-Stefan擴(kuò)散模型中，擴(kuò)散矩陣由半經(jīng)驗(yàn)擴(kuò)散系數(shù)計(jì)算得到；因此，擴(kuò)散系數(shù)是用于計(jì)算的關(guān)鍵因子。本文使用半經(jīng)驗(yàn)法計(jì)算各組分?jǐn)U散系數(shù)(D)。

(1)

式(1)中：T為熱力學(xué)溫度；P為大氣壓強(qiáng)；MA、MB為氣體A、B的相對(duì)分子質(zhì)量，本研究中二氧化碳?xì)怏w的擴(kuò)散介質(zhì)為氮?dú)?，因此MA、MB分別取值MA=44、MB=28；VA、VB為氣體A、B在正常沸點(diǎn)時(shí)的摩爾容積，VA、VB分別取值34.0和29.9?？梢钥闯觯珼為關(guān)于溫度和壓強(qiáng)的函數(shù)。利用半經(jīng)驗(yàn)擴(kuò)散系數(shù)能夠很好地描述二元組分之間的擴(kuò)散規(guī)律，這也是流體力學(xué)領(lǐng)域中最為常用的擴(kuò)散系數(shù)公式[18]。擴(kuò)散系數(shù)與濃度梯度的乘積即為擴(kuò)散通量。

1.2 Maxwell-Stefan擴(kuò)散模型

Maxwell-Stefan擴(kuò)散模型是適用于計(jì)算多組分混合物擴(kuò)散通量的數(shù)學(xué)模型，是一種基于濃度梯度進(jìn)行通量計(jì)算的方法。它在Fick擴(kuò)散模型的基礎(chǔ)上對(duì)擴(kuò)散系數(shù)進(jìn)行了改進(jìn)和擴(kuò)展，因此Fick模型是該模型在二元組分下的一種特殊形式。Maxwell-Stefan擴(kuò)散模型將半經(jīng)驗(yàn)擴(kuò)散系數(shù)拓展為擴(kuò)散系數(shù)矩陣，其中主對(duì)角線(xiàn)上的元素分別反映了物質(zhì)的自擴(kuò)散系數(shù)，主對(duì)角線(xiàn)以外的元素反映了多種物質(zhì)間的相互影響[19]。Maxwell-Stefan擴(kuò)散模型從分子之間的相互作用力出發(fā)進(jìn)行推算。促使分子發(fā)生擴(kuò)散的力被稱(chēng)為分子驅(qū)動(dòng)力[20]，由d表示，對(duì)于組分i：

(2)

式(2)中：μi為組分i的分子勢(shì)能；xi為組分i的摩爾分?jǐn)?shù)；R為氣體常數(shù)，8.314 J·mol-1·K-1；grad表示梯度算子。

(3)

式(3)中：x為氣體分子間距，μi可表示為

(4)

對(duì)于二元混合物而言，促使組分i和組分j分子發(fā)生相對(duì)運(yùn)動(dòng)的驅(qū)動(dòng)力與2種分子的速度差和組分j的物質(zhì)的量分?jǐn)?shù)成正相關(guān)，可表示為

-dμi/dz=(RT/D)·[xj(ui-uj)]。

(5)

式(5)中：z為分子移動(dòng)方向，ui、uj分別表示組分i、j的分子移動(dòng)速度，D為二元組分系統(tǒng)下的半經(jīng)驗(yàn)擴(kuò)散系數(shù)。對(duì)于多元混合物，式(5)可拓展為

-dμi/dz=RT[xij(ui-uj)/Dij+xik(ui-uk)/Dik+…+xin(ui-un)/Din]。

(6)

更普遍的形式為

(7)

又有

Ni=ctxiui=Ji+xiNt，

(8)

式中，Ni為組分i的擴(kuò)散通量，ct為混合物總濃度，Ji為組分i受分子作用力影響下的通量，Nt為混合物的總體擴(kuò)散通量。整理可得：

(J)=ct[B]-1[Γ](x)。

(9)

式(9)即為多元混合物擴(kuò)散的普遍化Maxwell-Stefan模型。其中：(J)為一維列矩陣，矩陣中每個(gè)元素為對(duì)應(yīng)組分的擴(kuò)散通量；[B]為擴(kuò)散系數(shù)矩陣；Γ為熱力學(xué)系數(shù)；x為各個(gè)組分的濃度梯度。式中：

(10)

(11)

式(11)中，δij為單位沖激函數(shù)，n為參與擴(kuò)散的組分?jǐn)?shù)量。

在本研究中，參與擴(kuò)散的介質(zhì)為氧氣和二氧化碳，擴(kuò)散介質(zhì)為氮?dú)夂推渌粎⑴c擴(kuò)散的氣體組成的混合物，由于其他氣體物質(zhì)含量較低，對(duì)實(shí)驗(yàn)結(jié)果產(chǎn)生的影響可忽略不計(jì)，因此介質(zhì)可以近似成氮?dú)?。根?jù)上述方程，土壤呼吸過(guò)程的通量模型可表示為

(12)

式(12)中：D12表示二氧化碳和氧氣之間的擴(kuò)散系數(shù)，D13表示二氧化碳和氮?dú)庵g的擴(kuò)散系數(shù)，D23表示氧氣和氮?dú)庵g的擴(kuò)散系數(shù)，x1為二氧化碳濃度，x2為氧氣濃度。根據(jù)式(12)，即可同時(shí)計(jì)算出每一時(shí)刻二氧化碳和氧氣的通量值。

1.3 數(shù)據(jù)來(lái)源

本研究所使用的各種氣體濃度數(shù)據(jù)由Ansys仿真系統(tǒng)進(jìn)行仿真獲得。Ansys為美國(guó)ANSYS公司開(kāi)發(fā)的有限元分析軟件，是國(guó)際上最受認(rèn)可的計(jì)算機(jī)輔助設(shè)計(jì)軟件，是一個(gè)非常適用于進(jìn)行流體力學(xué)仿真的軟件。本研究利用湍流模型，采用SIMPLE(Semi-Implicit Method for Pressure-Linked Equations)算法進(jìn)行仿真數(shù)據(jù)的生成。其中，湍流模型用于根據(jù)設(shè)定的初始值生成不同高度的二氧化碳濃度和氧氣濃度的時(shí)間序列，SIMPLE算法用于對(duì)產(chǎn)生的濃度進(jìn)行壓力修正。

二氧化碳和氧氣的時(shí)間序列使用Spalart-Allmaras(S-A)湍流模型仿真生成。本研究擬對(duì)圓桶狀氣室進(jìn)行仿真，S-A湍流模型適用于墻壁束縛流體運(yùn)動(dòng)，因此用其生成土壤呼吸的二氧化碳和氧氣的時(shí)間序列是可行的。在A(yíng)nsys-Fluent仿真系統(tǒng)中，需要為氣體流速提供邊界條件。氣室壁面的邊界條件由Fluent系統(tǒng)自動(dòng)生成，不需要手動(dòng)設(shè)定，因此只需要設(shè)定氣室入口處，即氣室底部的邊界條件，包括入口速度、入口壓力、入口濃度等邊界條件參數(shù)。由于入口邊界條件會(huì)在很大程度上影響到之后仿真時(shí)的擴(kuò)散過(guò)程和最終生成的仿真數(shù)據(jù)；因此，在很多情況下，設(shè)定正確的入口邊界條件都是十分重要的。本研究的邊界條件均由實(shí)驗(yàn)采集而來(lái)。

仿真實(shí)驗(yàn)與真實(shí)的土壤環(huán)境存在的差異主要包括：第一，在實(shí)測(cè)條件下，不同氣體的密度差異會(huì)導(dǎo)致擴(kuò)散規(guī)律發(fā)生一定的改變。對(duì)于這一影響因素，Ansys仿真系統(tǒng)在進(jìn)行初始化時(shí)會(huì)對(duì)各個(gè)組分的密度進(jìn)行設(shè)定，生成的仿真數(shù)據(jù)都是將重力因素計(jì)算在內(nèi)的。第二，受水平氣流的影響，氣流通過(guò)氣室出口時(shí)可能會(huì)流入氣室內(nèi)部，從而干擾二氧化碳的擴(kuò)散。為了盡可能地與真實(shí)的監(jiān)測(cè)條件接近，對(duì)于這一問(wèn)題，本文使用了S-A湍流模型。當(dāng)氣流流入氣室內(nèi)，會(huì)對(duì)內(nèi)部氣體造成低雷諾湍流運(yùn)動(dòng)，而S-A擴(kuò)散模型最適用于低雷諾指數(shù)的湍流分析。實(shí)驗(yàn)將湍流強(qiáng)度設(shè)定為1%，接近于氣流對(duì)于氣室內(nèi)部氣體的影響。

本研究分別進(jìn)行3組仿真實(shí)驗(yàn)，設(shè)定的二氧化碳通量分別是0.5、1.0、2.0 μmol·m-2·s-1，具體仿真參數(shù)見(jiàn)表1。

1.4 實(shí)現(xiàn)步驟

微生物和植物根部的呼吸作用會(huì)產(chǎn)生二氧化碳，土壤內(nèi)部的二氧化碳濃度為0.15%～0.65%，遠(yuǎn)高于大氣中0.039%的濃度，二氧化碳會(huì)隨土壤孔隙到達(dá)地表，使得地表的二氧化碳濃度最高可達(dá)0.16%。土壤表面和大氣之間形成的二氧化碳濃度差會(huì)導(dǎo)致二氧化碳由土壤表面向上方擴(kuò)散，因此二氧化碳濃度會(huì)隨高度的增加而降低，到達(dá)離地面約350 mm位置時(shí)與大氣濃度持平。當(dāng)土壤呼吸達(dá)到穩(wěn)態(tài)時(shí)，二氧化碳會(huì)產(chǎn)生有規(guī)律的濃度梯度，采用線(xiàn)性擬合法就可以得到濃度梯度，此時(shí)便于使用Maxwell-Stefan擴(kuò)散模型進(jìn)行二氧化碳通量的計(jì)算。

仿真實(shí)驗(yàn)的建立步驟如下。

(1)建立氣室模型。本研究所建立的圓筒狀氣室，直徑200 mm，高度350 mm，氣室底部與頂部均開(kāi)口，其中底部緊貼土壤，二氧化碳由氣室底部向頂部進(jìn)行擴(kuò)散。氣室筒壁厚度為5 mm，氣室頂部有3圈均勻分布的氣孔用于與外界的大氣進(jìn)行氣體交換，氣孔的直徑為5 mm，3圈通氣孔距離氣室橫截面中心的距離分別為65、85、105 mm。建立好的氣室模型如圖1所示。之后將模型導(dǎo)入到Fluent程序中。

表1 仿真參數(shù)

Table 1 Simulation parameter

通量Flux/(μmol·m-2·s-1)速率Rate/(kg·m-2)T/KP/kPa入口CO2體積分?jǐn)?shù)Inlet CO2 volumnfraction入口O2體積分?jǐn)?shù)Inlet O2 volumnfraction出口CO2體積分?jǐn)?shù)Outlet CO2 volumnfraction出口O2體積分?jǐn)?shù)Outlet O2 volumnfraction 0.52.2×10-8298.51 0100.000 810.193 50.000 5960.209 81.04.4×10-8298.51 0100.001 070.193 50.000 5960.209 82.08.8×10-8298.51 0100.002 290.193 50.000 5960.209 8

(2)劃分網(wǎng)格和建立控制方程。網(wǎng)格劃分的精密程度決定了仿真實(shí)驗(yàn)的效果。網(wǎng)格的作用在于將氣室這個(gè)整體分解成不同的區(qū)域，從而對(duì)每個(gè)區(qū)域分別進(jìn)行濃度動(dòng)態(tài)變化的推演。網(wǎng)格劃分程度越高，所占用的計(jì)算機(jī)性能也就越多。Ansys仿真系統(tǒng)支持多種類(lèi)型的網(wǎng)絡(luò)。本研究采用三角形非結(jié)構(gòu)化網(wǎng)格，其優(yōu)點(diǎn)是生成網(wǎng)格的時(shí)間較少，適用于較為復(fù)雜的實(shí)體，并且對(duì)于未知方向的氣流具有良好的計(jì)算性能，雖然計(jì)算數(shù)據(jù)的時(shí)間較長(zhǎng)，但數(shù)據(jù)精確度高。實(shí)驗(yàn)中將氣室劃分成大約60萬(wàn)個(gè)網(wǎng)格。圖2為氣室網(wǎng)格劃分示意圖。

氣室內(nèi)部的氣體屬于典型的牛頓流體，通過(guò)牛頓的本構(gòu)方程可以獲得較為全面的流體力學(xué)方程組。在流場(chǎng)求解時(shí)所用到的動(dòng)量方程可表示為

a，氣室模型；b，氣室半剖視圖。a， Chamber model; b， Half section view of chamber.圖1 氣室結(jié)構(gòu)圖Fig.1 Structural diagram of gas chamber

a，氣室壁面；b，氣室底部。a， Chamber wall; b， Chamber bottom.圖2 氣室網(wǎng)格劃分示意圖Fig.2 Schematic chart of gas chamber mesh generation

(13)

式(13)中：φ為通用因變量，Sφ為廣義源項(xiàng)，u為速度矢量，Гφ為廣義擴(kuò)散系數(shù)，ρ為密度，t為時(shí)間。

(3)邊界條件設(shè)定。邊界條件直接影響到生成數(shù)據(jù)的結(jié)果，因此應(yīng)盡量接近大氣環(huán)境值。組分氣體采用二階迎風(fēng)格式進(jìn)行離散，其他采用一階迎風(fēng)格式離散，氣室避免采用無(wú)滑移條件。采用非定常計(jì)算，時(shí)間步長(zhǎng)設(shè)定為0.1 s，滿(mǎn)足Courant-Friedrich-Levy穩(wěn)定準(zhǔn)則。氣室內(nèi)部氣體為近地表空氣，主要成分為氧氣、二氧化碳、氮?dú)夂推渌麣怏w，初始溫度為298.5 K，考慮重力加速度與能量方程，土壤呼吸界面采用質(zhì)量通量，湍流強(qiáng)度為1%，氣室出口界面設(shè)定為自由流出。不同土壤所排出的組分氣體含量不同，本研究采用的地表空氣與土壤空氣的各項(xiàng)組分氣體含量如表1所示。由于其他氣體組分復(fù)雜且含量較少，對(duì)本研究的影響可以忽略不計(jì)，因此其他氣體均由水蒸氣替代。

(4)由Fluent程序設(shè)定仿真初始化參數(shù)，并開(kāi)始生成仿真數(shù)據(jù)，產(chǎn)生0～0.35 m處不同高度氣體濃度的2 000 s時(shí)間序列。此時(shí)間段足夠使二氧化碳的擴(kuò)散在氣室內(nèi)形成穩(wěn)態(tài)。

(5)根據(jù)設(shè)定的溫度和壓強(qiáng)計(jì)算Maxwell-Stefan擴(kuò)散系數(shù)矩陣，并求得其逆矩陣。

(6)根據(jù)仿真系統(tǒng)生成的二氧化碳和氧氣濃度時(shí)間序列，每隔10 s取一組0、0.10、0.20、0.30、0.35 m高度的濃度數(shù)據(jù)，從而獲得200組二氧化碳和氧氣的濃度數(shù)據(jù)。仿照獲取的濃度數(shù)據(jù)，根據(jù)穩(wěn)態(tài)情況選取相應(yīng)的擬合方法計(jì)算出200組二氧化碳和氧氣濃度關(guān)于高度的梯度值，并將每一組的二氧化碳和氧氣濃度數(shù)據(jù)組成一個(gè)二階矩陣。

(7)將擴(kuò)散系數(shù)矩陣與每一組濃度梯度矩陣相乘，獲得每一時(shí)刻二氧化碳的通量值。

2 結(jié)果與分析

2.1 二氧化碳濃度時(shí)間分析

在碳通量的監(jiān)測(cè)過(guò)程中，二氧化碳濃度的變化規(guī)律是監(jiān)測(cè)的一項(xiàng)重要指標(biāo)，只有產(chǎn)生較為穩(wěn)定的濃度梯度，才能夠使用Maxwell-Stefan擴(kuò)散模型進(jìn)行碳通量的計(jì)算。因此，研究二氧化碳的動(dòng)態(tài)變化，揭示兩相混合氣體在氣室內(nèi)的自由擴(kuò)散機(jī)理，對(duì)提升土壤呼吸監(jiān)測(cè)采樣精度具有重大意義。圖3為氣室內(nèi)部二氧化碳濃度變化。

從圖3可以明顯看出，整個(gè)氣室內(nèi)的二氧化碳濃度隨時(shí)間延長(zhǎng)而增加。在擴(kuò)散過(guò)程中，大氣表層與土壤氣體由于缺少?gòu)?qiáng)對(duì)流作用，二氧化碳濃度的梯度變化較為平緩。經(jīng)過(guò)一段時(shí)間的擴(kuò)散，氣室內(nèi)的二氧化碳?xì)怏w逐層上升，并在600 s之后趨于穩(wěn)定。由此不難得出二氧化碳的擴(kuò)散規(guī)律，首先是二氧化碳?xì)怏w在氣室內(nèi)充分自由擴(kuò)散的階段，之后是二氧化碳達(dá)到飽和，進(jìn)入平衡階段。

仿真程序根據(jù)設(shè)定的通量值(0.5、1.0、2.0 μmol·m-2·s-1)分別進(jìn)行二氧化碳的擴(kuò)散仿真，并生成2 000 s內(nèi)每一時(shí)刻不同高度的氣體濃度值，其中二氧化碳的濃度值如圖4所示。

a～c中設(shè)定的通量值分別為0.5、1.0、2.0 μmol·m-2·s-1。下同。The set flux in a-c was 0.5， 1.0， 2.0 μmol·m-2·s-1， respectively. The same as below.圖4 不同高度下CO2濃度的時(shí)間序列Fig.4 Carbon dioxide concentration of different altitude

由圖4可知，氣室內(nèi)部二氧化碳的濃度隨時(shí)間延長(zhǎng)而增加，最后形成穩(wěn)態(tài)，即氣室內(nèi)部各個(gè)高度的二氧化碳濃度不再隨時(shí)間而顯著變化。造成這一現(xiàn)象的原因是，開(kāi)始時(shí)氣室內(nèi)部的二氧化碳濃度相比氣室入口處少，二氧化碳由高濃度向低濃度擴(kuò)散，而隨著氣室上方的二氧化碳濃度增加，下方濃度與上方濃度無(wú)法再形成能夠產(chǎn)生擴(kuò)散的濃度差，此時(shí)擴(kuò)散運(yùn)動(dòng)就趨于穩(wěn)定。當(dāng)設(shè)定的二氧化碳通量值較低時(shí)，氣室上方的二氧化碳濃度會(huì)呈現(xiàn)先下降后上升的趨勢(shì)。這是因?yàn)槎趸細(xì)怏w的密度相對(duì)空氣較大，而此時(shí)底部的二氧化碳濃度還未升高，因此二氧化碳不足以擴(kuò)散到足夠的高度。當(dāng)擴(kuò)散達(dá)到穩(wěn)態(tài)時(shí)，氣室內(nèi)部的二氧化碳分布較為均衡，因此這一現(xiàn)象并不會(huì)影響到最終的計(jì)算結(jié)果。開(kāi)放型擴(kuò)散氣室的設(shè)計(jì)能夠使高濃度處的二氧化碳以垂直方向向低濃度處擴(kuò)散而不受水平氣流的影響，從而使二氧化碳的濃度變化規(guī)律更加明顯，這樣更易于使用Maxwell-Stefan擴(kuò)散模型進(jìn)行計(jì)算。另外，這樣的氣室還能夠有效地減少氣室對(duì)土壤呼吸的抑制作用。若使用封閉型氣室，二氧化碳會(huì)在氣室內(nèi)部積累，一段時(shí)間過(guò)后氣室內(nèi)部各個(gè)高度的二氧化碳濃度將會(huì)相等，并且會(huì)隨時(shí)間延長(zhǎng)而不斷增加，直至與地表的二氧化碳濃度一致從而不再發(fā)生擴(kuò)散現(xiàn)象。從圖4的數(shù)據(jù)可以看出，各個(gè)高度的二氧化碳濃度都有不同的變化曲線(xiàn)，并沒(méi)有趨于一致，這也證明了開(kāi)放型氣室設(shè)計(jì)的優(yōu)勢(shì)。其中，0.5 μmol·m-2·s-1的擴(kuò)散通量下，二氧化碳濃度在600 s之后達(dá)到穩(wěn)態(tài)；在1.0 μmol·m-2·s-1的擴(kuò)散通量下，二氧化碳濃度在800 s之后達(dá)到穩(wěn)態(tài)；在2.0 μmol·m-2·s-1的擴(kuò)散通量下，二氧化碳濃度在1 620 s之后達(dá)到穩(wěn)態(tài)。氣室內(nèi)達(dá)到穩(wěn)態(tài)的時(shí)間隨通量的增加而延長(zhǎng)，而在濃度達(dá)到穩(wěn)態(tài)之前，濃度梯度逐漸減小，達(dá)到穩(wěn)態(tài)之后濃度與高度成線(xiàn)性關(guān)系，因此穩(wěn)態(tài)時(shí)濃度梯度的計(jì)算可以采用線(xiàn)性擬合法。這種線(xiàn)性關(guān)系也證明了在土壤呼吸系統(tǒng)中，使用Maxwell-Stefan擴(kuò)散模型計(jì)算二氧化碳擴(kuò)散通量是可行的。

2.2 Maxwell-Stefan模型計(jì)算結(jié)果分析

根據(jù)獲得的200組間隔為10 s的數(shù)據(jù)，取每一組中0.10、0.20、0.30、0.35 m高度的濃度數(shù)據(jù)進(jìn)行擬合，得出濃度梯度，并用Maxwell-Stefan擴(kuò)散模型進(jìn)行計(jì)算，得出2 000 s時(shí)間段內(nèi)的通量值，計(jì)算所得通量值的變化曲線(xiàn)見(jiàn)圖5。

圖5顯示了由Maxwell-Stefan擴(kuò)散模型計(jì)算得到的二氧化碳通量值時(shí)間序列。初始時(shí)氣室內(nèi)部的二氧化碳濃度值與大氣相等，而氣室入口處的二氧化碳濃度較高，導(dǎo)致濃度梯度大，此時(shí)計(jì)算的通量值要高于實(shí)際設(shè)定的通量值，而隨時(shí)間延長(zhǎng)，通量值逐漸減小，最終趨于穩(wěn)定，至穩(wěn)態(tài)時(shí)，設(shè)定的0.5、1.0、2.0 μmol·m-2·s-1的通量下，計(jì)算結(jié)果分別為0.547、0.969、2.122 μmol·m-2·s-1，實(shí)驗(yàn)結(jié)果與設(shè)定值的誤差分別為9.4%、6.9%、6.1%?？梢钥闯?，計(jì)算的精確度隨設(shè)定通量值的提高而升高，產(chǎn)生這一現(xiàn)象的原因主要有2點(diǎn)：第一，二氧化碳擴(kuò)散通量較低時(shí)，重力對(duì)擴(kuò)散的影響較大，向上擴(kuò)散的速度無(wú)法抵消重力的影響，導(dǎo)致一部分二氧化碳無(wú)法向上排放從而在氣室底部積累，使氣室下方的濃度梯度增大，擴(kuò)散模型所計(jì)算出的通量結(jié)果也因此增加。第二，較小的擴(kuò)散通量使得組分之間的相互影響力減小，而Mawell-Stefan模型的計(jì)算效果在組分的影響因素更多時(shí)才會(huì)更加明顯。但要說(shuō)明的是，以上這2點(diǎn)并不會(huì)對(duì)擴(kuò)散模型的計(jì)算結(jié)果產(chǎn)生較大的影響，實(shí)驗(yàn)結(jié)果也證明，利用Maxwell-Stefan擴(kuò)散模型計(jì)算的土壤呼吸時(shí)排放二氧化碳的通量與設(shè)定通量值非常相近。

圖5 Maxwell-Stefan通量值與設(shè)定通量值對(duì)比Fig.5 Comparison of calculated flux by Maxwell-Stefan with set flux

2.3 Maxwell-Stefan與Fick擴(kuò)散模型計(jì)算結(jié)果對(duì)比

Maxwell-Stefan擴(kuò)散模型是Fick模型的普遍情況，將Maxwell-Stefan模型計(jì)算得到的碳通量結(jié)果與Fick模型的計(jì)算結(jié)果進(jìn)行對(duì)比，以便說(shuō)明本研究中算法的合理性和優(yōu)勢(shì)所在。每隔10 s取一組二氧化碳濃度值，構(gòu)成一個(gè)由200組濃度構(gòu)成的時(shí)間序列，并對(duì)二氧化碳濃度進(jìn)行擬合，得到二氧化碳梯度的時(shí)間序列，之后運(yùn)用Fick模型進(jìn)行二氧化碳通量的計(jì)算，其計(jì)算結(jié)果與Maxwell-Stefan模型結(jié)果的對(duì)比情況如圖6所示。在穩(wěn)態(tài)條件下，F(xiàn)ick擴(kuò)散模型所計(jì)算的二氧化碳通量分別為0.589、1.044、2.286 μmol·m-2·s-1?？梢钥闯?，在多組實(shí)驗(yàn)下，本文算法所得的計(jì)算結(jié)果都小于Fick模型所得結(jié)果，且更接近設(shè)定值。當(dāng)二氧化碳向氣室上方擴(kuò)散時(shí)，分子間的作用力會(huì)減緩這一過(guò)程，這會(huì)使擴(kuò)散的通量值有所降低。Maxwell-Stefan擴(kuò)散系數(shù)矩陣對(duì)分子間的這一影響因子進(jìn)行了校正，在計(jì)算時(shí)減去了被減緩的通量。因此，利用這一模型所計(jì)算出的通量值為消除分子間作用之后的通量值。而Fick擴(kuò)散模型將除二氧化碳之外的氣體歸于一種組分，不能夠準(zhǔn)確表述多組分之間的相互影響和擴(kuò)散過(guò)程，這就導(dǎo)致了在計(jì)算多組分?jǐn)U散系統(tǒng)的通量時(shí)Fick模型的效果不如Maxwell-Stefan模型，并且計(jì)算結(jié)果偏大。本文3組實(shí)驗(yàn)的結(jié)果很好地證明了這一點(diǎn)。顯然，Maxwell-Stefan擴(kuò)散模型從理論上能夠更精確地計(jì)算土壤呼吸時(shí)排放二氧化碳的通量。

圖6 Maxwell-Stefan模型與Fick模型的結(jié)果對(duì)比Fig.6 Comparison of flux result calculated by Maxwell-Stefan model and Fick model

3 討論

Maxwell-Stefan擴(kuò)散模型很好地解決了計(jì)算多組分?jǐn)U散系統(tǒng)中擴(kuò)散通量的問(wèn)題，彌補(bǔ)了Fick擴(kuò)散模型的不足之處。開(kāi)放型的氣室很好地解決了密閉氣室抑制土壤呼吸時(shí)排出二氧化碳的問(wèn)題。本文將這兩者相結(jié)合，進(jìn)行仿真研究，并將計(jì)算結(jié)果與仿真設(shè)定的通量值進(jìn)行對(duì)比，根據(jù)實(shí)驗(yàn)分析與計(jì)算的結(jié)果得出以下結(jié)論。

(1)土壤呼吸作為多組分?jǐn)U散系統(tǒng)，利用Maxwell-Stefan擴(kuò)散模型計(jì)算碳通量是可行的。

(2)開(kāi)放型氣室的設(shè)計(jì)不僅解決了抑制土壤呼吸的問(wèn)題，還使氣室內(nèi)的二氧化碳濃度梯度呈線(xiàn)性趨勢(shì)，能夠?qū)崿F(xiàn)對(duì)碳通量的長(zhǎng)時(shí)間持續(xù)監(jiān)測(cè)。

(3)在利用Maxwell-Stefan擴(kuò)散模型計(jì)算碳通量時(shí)可以得到與實(shí)際相近的計(jì)算結(jié)果。

土壤呼吸作為一種由多組分參與的擴(kuò)散過(guò)程，理論上能夠運(yùn)用Maxwell-Stefan模型進(jìn)行計(jì)算。該模型基于濃度梯度計(jì)算碳通量，能夠解決傳統(tǒng)監(jiān)測(cè)方法對(duì)二氧化碳產(chǎn)生抑制的問(wèn)題，從而實(shí)現(xiàn)長(zhǎng)時(shí)間的連續(xù)監(jiān)測(cè)。為了驗(yàn)證該方法有效，本研究進(jìn)行了2個(gè)階段的實(shí)驗(yàn)；首先，分析氣室內(nèi)部二氧化碳濃度的動(dòng)態(tài)變化。從仿真實(shí)驗(yàn)的結(jié)果可以看出，氣室內(nèi)部的二氧化碳在穩(wěn)態(tài)之后呈現(xiàn)出穩(wěn)定的濃度梯度，這使得基于Maxwell-Stefan模型計(jì)算碳通量具備可行性。之后，通過(guò)程序根據(jù)S-A湍流模型生成不同高度處的濃度時(shí)間序列進(jìn)行碳通量計(jì)算，得出每一時(shí)刻下的碳通量值。通量值與設(shè)定值對(duì)比的結(jié)果表明，該方法能夠得出較為準(zhǔn)確的通量計(jì)算結(jié)果。另外，還將前述計(jì)算結(jié)果與Fick模型進(jìn)行了對(duì)比。Fick模型作為Maxwell-Stefan模型的一種二元組分特殊形式，已經(jīng)取得了較好的二氧化碳通量監(jiān)測(cè)效果[21]。Krishna等[22]利用氣井法并使用Fick模型進(jìn)行土壤呼吸碳通量的計(jì)算，計(jì)算得出的最大誤差在16%以?xún)?nèi)，這一精度在碳通量監(jiān)測(cè)領(lǐng)域是非常領(lǐng)先的。Maxwell-Stefan因其適用于多組分?jǐn)U散這一優(yōu)點(diǎn)，應(yīng)該能夠取得更加精確的計(jì)算結(jié)果[23]，本實(shí)驗(yàn)也很好地證明了這一猜想。在多組實(shí)驗(yàn)中，Maxwell-Stefan模型所計(jì)算出的二氧化碳通量值都比Fick模型更加接近設(shè)定值。本研究從理論和實(shí)驗(yàn)上同時(shí)證明了Maxwell-Stefan更加適用于計(jì)算土壤碳通量并且能夠得到很好的計(jì)算結(jié)果。

本研究證明了Maxwell-Stefan模型與開(kāi)放氣室相結(jié)合的可行性，為后續(xù)研究土壤生態(tài)環(huán)境提供了一種新的方法。但在野外實(shí)地測(cè)量時(shí)，可能還會(huì)受到許多其他因素的干擾，解決這些問(wèn)題將是今后拓展該方法實(shí)用性的重點(diǎn)。此外，今后應(yīng)著重探討的地方還有2點(diǎn)：第一，不同地區(qū)和類(lèi)型的土壤中所排放的氣體有所不同，利用本文模型計(jì)算時(shí)應(yīng)考慮到更多的擴(kuò)散組分；第二，在真實(shí)環(huán)境下測(cè)量時(shí)，水平氣流的速率并不是恒定的，因此在關(guān)于野外情況下如何對(duì)濃度梯度進(jìn)行擬合，還應(yīng)做進(jìn)一步研究。