戚春香，李 瑤，馬 琳，董 彬

（中國(guó)民航大學(xué)機(jī)場(chǎng)學(xué)院，天津 300300）

根據(jù)中國(guó)機(jī)場(chǎng)建設(shè)規(guī)劃，東北地區(qū)將開(kāi)展大量的機(jī)場(chǎng)建設(shè)以及改擴(kuò)建工作，其中包括大小興安嶺機(jī)場(chǎng)群的建設(shè)。大小興安嶺地區(qū)分布有大面積的高緯度凍土，而目前對(duì)于凍土地基的研究主要集中于極地地區(qū)及中國(guó)青藏高原地區(qū)。Goering 等[1]對(duì)片塊石凍土路基溫度場(chǎng)進(jìn)行了數(shù)值模擬；Andrew[2]研究了加拿大北極地區(qū)多年凍土地基溫度場(chǎng)狀況；Kudriavtcev 等[3]分析了俄羅斯凍土區(qū)鐵路路基溫度場(chǎng)的分布。由于溫度場(chǎng)上邊界條件復(fù)雜多變，建立模型時(shí)有諸多不便，朱林楠[4]基于青藏公路實(shí)測(cè)地基溫度提出了附面層理論，用穩(wěn)定的附面層代替了溫度場(chǎng)復(fù)雜的上邊界條件，在多年凍土區(qū)溫度場(chǎng)研究中得到日益廣泛的應(yīng)用；金龍等[5]根據(jù)青藏高原地基溫度資料及附面層原理，將青藏公路典型路基溫度場(chǎng)模型上邊界條件表示為溫度函數(shù)；劉偉[6]應(yīng)用青藏高原附面層理論分析了黑色寬幅路基溫度場(chǎng)規(guī)律；曹元兵等[7]應(yīng)用附面層理論建立了瑪多地區(qū)多年凍土地基溫度場(chǎng)模型；白青波[8]根據(jù)西部地區(qū)主要路段地基溫度日波動(dòng)狀況自定義了附面層厚度的選取標(biāo)準(zhǔn)，并根據(jù)該標(biāo)準(zhǔn)探究了附面層厚度及溫度增量的取值。楊簡(jiǎn)[9]應(yīng)用附面層理論建立了青藏高原多年凍土區(qū)機(jī)場(chǎng)跑道的溫度場(chǎng)模型。

大小興安嶺地區(qū)氣候環(huán)境及地基溫度狀況都與青藏高原不同，青藏高原多年凍土環(huán)境及地基溫度狀況確定的附面層理論參數(shù)是否適用于東北地區(qū)多年凍土地區(qū)還需深入研究，為此選取中國(guó)目前唯一已進(jìn)行機(jī)場(chǎng)建設(shè)的多年凍土區(qū)——漠河作為研究對(duì)象，利用ABAQUS 建立漠河地區(qū)機(jī)場(chǎng)跑道多年凍土地基溫度場(chǎng)模型，以確定漠河地區(qū)附面層理論參數(shù)取值，為東北地區(qū)多年凍土區(qū)機(jī)場(chǎng)跑道建設(shè)提供理論基礎(chǔ)。

1 有限元模型

1.1 地基溫度場(chǎng)模型結(jié)構(gòu)

采用二維對(duì)稱模型，DCAX8 計(jì)算單元，其1/2 結(jié)構(gòu)如圖1所示。整體模型道面板總寬為60 m，兩側(cè)道肩寬度為40 m，此寬度可忽略邊界效應(yīng)的影響[9-10]，土體厚度取40 m。

圖1 地基溫度場(chǎng)模型Fig.1 Foundation temperature field model

Ⅰ區(qū)為道面結(jié)構(gòu)層，根據(jù)《民用機(jī)場(chǎng)瀝青道面設(shè)計(jì)規(guī)范》[11]的相關(guān)規(guī)定，其組成如表1所示。道面結(jié)構(gòu)總厚度為1.4 m，滿足最小防凍層厚度要求[11]。Ⅱ、Ⅲ、Ⅳ區(qū)為土層結(jié)構(gòu)，Ⅱ區(qū)為1 m 厚亞粘土層（35%含水量），Ⅲ區(qū)為4 m 厚亞粘土層（15%含水量），Ⅳ區(qū)為35 m 厚碎石亞粘土層。

表1 道面結(jié)構(gòu)層Tab.1 Pavement structure layers

1.2 熱參數(shù)

由于水分遷移的影響，模型Ⅱ、Ⅲ、Ⅳ區(qū)土層熱參數(shù)應(yīng)采用等效比熱容c（T）和等效導(dǎo)熱系數(shù)[10]λ（T），分別為

其中：T為地基溫度（℃）；C（T）為體積熱容量；ρ為土體干密度（kg·m-3）；c1、c2分別為融土和凍土的比熱容（J·kg-1·℃-1）；w0為初始含水率（%）；c0、c3分別為水和冰的比熱容，c0=4 182 J·kg-1·℃-1，c3=2 090 J·kg-1·℃-1；T1為土體凍結(jié)溫度（℃）；A、B、b為經(jīng)驗(yàn)系數(shù)，A = 2.3 521，B=3.2 659，b與土質(zhì)有關(guān)，可通過(guò)試驗(yàn)測(cè)得；L為水的相變潛熱，L=334 560 J·kg-1；λx（T）、λy（T）分別為橫向和縱向等效導(dǎo)熱系數(shù)，λ1、λ2分別為融土和凍土的導(dǎo)熱系數(shù)（W·m-1·℃-1）。

各參數(shù)取值如表2、表3所示[12-14]，在ABAQUS 中利用Fortran 語(yǔ)言編寫UMATHT 子程序。

表2 土層基本參數(shù)Tab.2 Basic parameters of soil layers

表3 土層熱參數(shù)Tab.3 Thermal parameters of soil layers

模型Ⅰ區(qū)為道面結(jié)構(gòu)層，由于其保水性差可忽略水分遷移，因此無(wú)需采用等效熱參數(shù)，其各層熱參數(shù)如表4所示[10，15]。

表4 道面結(jié)構(gòu)層參數(shù)Tab.4 Parameters of pavement structure layers

1.3 邊界條件

由于附面層理論是基于上邊界條件的研究，因此只考慮模型上邊界條件的施加，利用廣泛采用的綜合考慮第二類、第三類上邊界條件，在ABAQUS 中分別定義模型上邊界的蒸發(fā)耗熱、對(duì)流換熱條件以及太陽(yáng)輻射條件。

研究[15]表明：瀝青混凝土道面的蒸發(fā)耗熱小于太陽(yáng)輻射的5%～10%，對(duì)地基影響相對(duì)很??；同時(shí)天然地表在凍結(jié)狀態(tài)幾乎沒(méi)有蒸發(fā)耗熱。因此，在計(jì)算地基溫度場(chǎng)時(shí)可忽略蒸發(fā)耗熱的影響。

1）對(duì)流換熱條件

道面與大氣之間的對(duì)流換熱主要與氣溫及兩者之間的對(duì)流換熱系數(shù)有關(guān)，通過(guò)FILM 子程序定義氣溫函數(shù)以及各月對(duì)流換熱系數(shù)即可實(shí)現(xiàn)對(duì)流換熱條件的定義。根據(jù)中國(guó)氣象數(shù)據(jù)網(wǎng)的資料，漠河地區(qū)1981—2010 年各月平均氣溫如表5所示。

表5 漠河地區(qū)各月氣溫Tab.5 Monthly atmospheric temperature in Mohe district ℃

考慮0.05 ℃/年的升溫效應(yīng)，則漠河地區(qū)月平均氣溫為

其中：t為時(shí)間（月），t=0 時(shí)為4 月。

研究[16]表明：地表與大氣間的對(duì)流換熱系數(shù)h 主要受風(fēng)速的影響，其表達(dá)式為

其中：vw為平均風(fēng)速（m·s-1）。

將中國(guó)氣象數(shù)據(jù)網(wǎng)提供的漠河地區(qū)1981—2010年各月平均風(fēng)速代入式（4）可得其各月對(duì)流換熱系數(shù)，如表6所示。

表6 對(duì)流換熱系數(shù)Tab.6 Convective heat transfer coefficient W·m-·2℃-1

2）太陽(yáng)輻射條件

太陽(yáng)凈輻射為吸收的太陽(yáng)總輻射與地面有效輻射的差值，模擬時(shí)應(yīng)用DFLUX 子程序定義各月太陽(yáng)總輻射，如表7所示。同時(shí)定義太陽(yáng)輻射吸收率和瀝青混凝土面層太陽(yáng)輻射吸收率為0.87，天然地面太陽(yáng)輻射吸收率為0.71。

表7 各月太陽(yáng)總輻射量Tab.7 Monthly total solar radiation 108 J·m-2

在ABAQUS 中通過(guò)定義氣溫幅值曲線，同時(shí)在Edit Attributes 中定義Stefan-Boltzmann 常數(shù)以及絕對(duì)零度值，可實(shí)現(xiàn)地面有效輻射的定義，即

其中：ε為發(fā)射率（黑度）；σ為Stefan-Boltzmann 常數(shù)，σ=5.669 7×10-8W·m-2·K-4；T2為地表溫度（℃）；T3為絕對(duì)零度值，T3=273.15 ℃；T4為大氣溫度（℃）。

1.4 計(jì)算分析步

采用Heat transfer 瞬態(tài)分析步，計(jì)算時(shí)長(zhǎng)以月為單位，總時(shí)長(zhǎng)設(shè)為240，步長(zhǎng)設(shè)為1，即模型共計(jì)算20年，每步計(jì)算時(shí)長(zhǎng)為1個(gè)月，共需計(jì)算240 步。

2 漠河地區(qū)附面層參數(shù)分析

2.1 附面層理論

根據(jù)相關(guān)文獻(xiàn)[5-9]，對(duì)于多年凍土區(qū)地基溫度場(chǎng)模型，上邊界可由厚度為H的附面層代替，上邊界條件可用附面層底溫度函數(shù)表示，即

其中：T為地基溫度（℃）；T0為附面層底地基溫度初始值（℃）；A0為附面層底地基溫度振幅（℃）；t為計(jì)算時(shí)間（月），t=0 時(shí)代表第1 年4 月；α為大氣溫度升溫速率，取α=0.05 ℃·a-1。

2.2 參數(shù)分析

由式（6）可知，附面層參數(shù)為附面層厚度H 及附面層底地基溫度初始值T0、振幅A0。提取漠河地區(qū)溫度場(chǎng)數(shù)值，利用OriginPro 9 軟件對(duì)地基溫度數(shù)據(jù)進(jìn)行處理分析。

1）提取數(shù)據(jù)

根據(jù)相關(guān)資料，附面層通常位于道面結(jié)構(gòu)中，模型道面結(jié)構(gòu)總厚度為1.4 m，因此分別提取瀝青混凝土道面和天然平整地面地基0～2 m 深度內(nèi)的溫度，每隔0.1 m 提取1 組數(shù)據(jù)，共40 組。為避免瀝青混凝土道面與天然平整地面地基溫度之間的相互影響，提取數(shù)據(jù)時(shí)根據(jù)溫度場(chǎng)提取等溫線近似直線處的地基溫度，分別提取了模型左、右邊界處地基溫度代表瀝青混凝土道面、天然平整地面下部地基溫度。

2）編寫函數(shù)

在OriginPro 9 軟件中打開(kāi)Fitting Function Organizer 編輯式（6）所示函數(shù)，參數(shù)為T0、A0，在Parameter Setting 中將其Value 值均設(shè)為1。

3）確定參數(shù)

將瀝青混凝土道面地基0～2 m 深度內(nèi)共21 組沿時(shí)間變化的地基溫度數(shù)據(jù)導(dǎo)入OriginPro 9 中，利用所編寫函數(shù)分別對(duì)每一組數(shù)據(jù)進(jìn)行擬合，其中時(shí)間為自變量，地基溫度為因變量。比較21 組數(shù)據(jù)的擬合度，找出擬合度最高的兩組數(shù)據(jù)，如圖2（a）所示，擬合度最高的為深度0.1 m與0.2 m 處的地基溫度。

圖2 不同深度處的地基溫度擬合度Fig.2 Fitting degree of ground temperature at different depths

以0.01 m為單位提取瀝青混凝土道面地基0.1～0.2 m 深度內(nèi)共11 組溫度隨時(shí)間變化的數(shù)據(jù)，比較其擬合度，如圖2（b）所示，深度0.12 m 處的地基溫度擬合度最高，附面層厚度即為0.12 m，擬合函數(shù)的參數(shù)T0、A0即為附面層底溫度函數(shù)的參數(shù)T0、A0，如表8所示。

同理，由天然平整地面地基0～2 m 深度內(nèi)溫度隨時(shí)間變化的數(shù)據(jù)亦可得到其對(duì)應(yīng)的附面層參數(shù)，如表8所示。

表8 漠河地區(qū)附面層參數(shù)Tab.8 Boundary layer parameters in Mohe district

3 結(jié)語(yǔ)

1）由于地域環(huán)境差異，不同凍土地區(qū)應(yīng)用附面層理論計(jì)算地基溫度場(chǎng)時(shí)其附面層參數(shù)的取值有所差異。

2）根據(jù)當(dāng)?shù)貧夂驓鉁刭Y料，應(yīng)用有限元軟件建立綜合二、三類邊界條件的溫度場(chǎng)模型，比較不同深度地基溫度與附面層理論的擬合度可確定附面層參數(shù)。

3）漠河地區(qū)跑道多年凍土地基溫度場(chǎng)瀝青混凝土道面對(duì)應(yīng)的附面層厚度為0.12 m，附面層底地基溫度初始值為1.90 ℃，振幅為25.09 ℃；天然地面對(duì)應(yīng)的附面層厚度為0.06 m，附面層底地基溫度初始值為1.03 ℃，振幅為25.20 ℃。

4）漠河古蓮機(jī)場(chǎng)正式通航之后，在多年凍土區(qū)建設(shè)機(jī)場(chǎng)成為可能，根據(jù)其附面層參數(shù)確定跑道有限元模型，利用有限元模擬分析其未來(lái)的溫度場(chǎng)變化趨勢(shì)，為及時(shí)發(fā)現(xiàn)跑道可能存在的問(wèn)題提供依據(jù)。