戚春香，韓 卓，馬 琳，張獻(xiàn)民，李 瑤

（中國民航大學(xué)交通科學(xué)與工程學(xué)院，天津 300300）

中國分布有大面積的多年凍土[1]。隨著中國民航業(yè)的不斷發(fā)展，多年凍土區(qū)的機(jī)場建設(shè)將陸續(xù)開展。多年凍土條件能否滿足建設(shè)寬幅跑道的特殊要求，是目前多年凍土地區(qū)建設(shè)機(jī)場面臨的關(guān)鍵難題。由于凍土的特殊性，在溫度變化和水分遷移作用下，多年凍土區(qū)的跑道極易出現(xiàn)融沉和凍脹等問題。因此，水熱力耦合作用下跑道結(jié)構(gòu)的溫度應(yīng)力大小對跑道的影響是多年凍土區(qū)機(jī)場設(shè)計過程中需要深入研究的問題。

針對凍土多場耦合的溫度應(yīng)力研究，國外學(xué)者最早對凍土的凍脹機(jī)理進(jìn)行分析[2-6]，為后續(xù)多場耦合研究提供了理論依據(jù)。而在國內(nèi)的相關(guān)研究中，文獻(xiàn)[7]分析了考慮水熱遷移作用下的正凍土熱傳導(dǎo)問題，并建立相關(guān)理論模型；文獻(xiàn)[8-9]開展了溫度場、滲流場和應(yīng)力場耦合問題的非線性分析，提出多場耦合條件下的水熱力耦合模型；文獻(xiàn)[10-12]基于水熱力三場耦合效應(yīng)對非飽和的孔隙凍土進(jìn)行了分析，建立了凍土水熱力三場耦合的理論框架及分析平臺；文獻(xiàn)[13-14]基于凍土路基水熱耦合模型分析了凍土路基的應(yīng)力分布；文獻(xiàn)[15-16]利用有限元模擬，并結(jié)合室內(nèi)模型試驗(yàn)分析了多年凍土水熱力場的耦合效應(yīng)，分析了凍土路基在凍結(jié)過程中的應(yīng)力場；文獻(xiàn)[17]利用有限元模擬和現(xiàn)場試驗(yàn)對青藏鐵路多年凍土路基坡腳排水渠道的應(yīng)力場進(jìn)行研究。上述關(guān)于多場耦合分析的研究多是針對多年凍土區(qū)公路與鐵路工程進(jìn)行分析的，由于機(jī)場跑道本身結(jié)構(gòu)的特殊性，其凍土地基的溫度場與公路和鐵路有很大不同[18]，因此，其應(yīng)力場也會存在較大差異。

本文依據(jù)大興安嶺地區(qū)多年凍土地區(qū)月平均溫度年變化數(shù)據(jù)，開展了冷季降溫過程下的凍土土基溫度場的室內(nèi)模型試驗(yàn)，并結(jié)合室內(nèi)試驗(yàn)?zāi)Ｐ屠肁BAQUS 軟件建立了考慮水熱耦合效應(yīng)的二維凍土地基模型，通過對比監(jiān)測試驗(yàn)數(shù)據(jù)與有限元計算結(jié)果證明了水熱耦合下溫度場計算方法的有效性。結(jié)合水、熱對力學(xué)參數(shù)的影響，建立水熱力耦合作用下的跑道地基有限元模型，并研究多年凍土區(qū)機(jī)場跑道道面結(jié)構(gòu)層和地基的水平溫度應(yīng)力分布規(guī)律，為陸續(xù)開展的多年凍土環(huán)境下機(jī)場跑道建設(shè)提供理論支持。

1 凍土土基溫度場室內(nèi)模型試驗(yàn)

1.1 模型試驗(yàn)裝置及土質(zhì)參數(shù)

模型試驗(yàn)在配有控溫系統(tǒng)、溫度采集系統(tǒng)和試驗(yàn)?zāi)Ｐ筒鄣膬鐾猎囼?yàn)室內(nèi)開展，試驗(yàn)室控溫幅度為-30～50 ℃。室內(nèi)模型槽尺寸為3.8 m×1.9 m×1.1 m，低溫試驗(yàn)室及室內(nèi)模型槽如圖1 和圖2 所示。

圖1 低溫試驗(yàn)室Fig.1 Low temperature laboratory

圖2 室內(nèi)模型槽Fig.2 Indoor model slot

首先，利用液塑限聯(lián)合測定儀測得模型土塑性指數(shù)為16.6，土類型為亞黏土；其次，在進(jìn)行模型槽填土?xí)r，為控制整體含水率與壓實(shí)度一致，采用分層填筑，每層厚度為0.22 m，共分為5 層。沿模型槽內(nèi)深度、長度方向分別取不同位置土樣10 組，每組3 個試樣，測得平均含水率及干密度，進(jìn)而依據(jù)《凍土地區(qū)建筑地基基礎(chǔ)設(shè)計規(guī)范》[19]，由土質(zhì)類型、含水率與干密度可查得土體相關(guān)熱參數(shù)，各參數(shù)如表1 所示。試驗(yàn)通過控制室內(nèi)負(fù)溫溫度來模擬凍土的地溫狀況。

表1 亞黏土參數(shù)Tab.1 Subclay parameters

1.2 溫度傳感器布置及試驗(yàn)溫度控制

1.2.1 溫度傳感器布置

模型土裝入深度為1.1 m，沿深度分別在橫斷面的中部和端部（距邊界0.3 m）利用桿布置溫度傳感器，開始工作時在ENC12/14 數(shù)據(jù)采集機(jī)箱中安裝CR1000數(shù)據(jù)采集器和AM16/32B 擴(kuò)展板，并將溫度傳感器通過線纜與采集器相連，啟動電源即可自動采集溫度數(shù)據(jù)。模型槽內(nèi)埋入溫度傳感器的布置情況如圖3 所示。

圖3 試驗(yàn)?zāi)Ｐ蜋M斷面及傳感器放置位置Fig.3 Cross section of experimental model and sensor distribution

1.2.2 確定試驗(yàn)溫度及持續(xù)時間

通過國家氣象科學(xué)數(shù)據(jù)中心查得東北某地區(qū)近30 年的月平均溫度年變化如圖4 所示。由圖4 可以看出，該地區(qū)每年的冷季大約出現(xiàn)在10 月到次年的4 月，月平均最低溫度為-29.2 ℃。

圖4 30 年東北某地區(qū)月平均溫度年變化Fig.4 Annual change of monthly average temperature in a certain area of Northeast China during 30 years

為模擬土基在冷季的降溫過程，選用該地區(qū)溫度曲線-29.2～0 ℃中的部分溫度段為室內(nèi)模型試驗(yàn)控制溫度。經(jīng)前期試驗(yàn)測得同一溫度設(shè)定時長為1 d 時，影響深度約為0.3 m。由于模型土深度為1.1 m，若使最上層溫度足以影響下部土層溫度，則同一溫度的合理持續(xù)時長至少為6 d。因此，為了能夠確保試驗(yàn)地溫控制準(zhǔn)確，設(shè)定同一溫度持續(xù)時長為7～14 d。設(shè)定的試驗(yàn)控制溫度及持續(xù)時長如表2 所示。

表2 試驗(yàn)控制溫度及持續(xù)時間Tab.2 Test control temperature and duration

1.3 模型試驗(yàn)結(jié)果

經(jīng)測得，模型土在第7、14、21、28 d 時對應(yīng)的溫度沿深度變化趨勢如圖5 所示。

圖5 模型土在不同降溫時長下的溫度變化趨勢Fig.5 Temperature trend of simulated soil under different cooling duration

從圖5 可以看出，當(dāng)外界溫度下降時，隨著降溫天數(shù)增加，各深度土層的溫度均不斷降低；隨著降溫天數(shù)逐漸增加，土層上下溫差逐漸增大，上部土層相比下部土層隨環(huán)境溫度下降更明顯。這說明地基土層深度的增加能明顯減小環(huán)境溫度的影響，且土基會存在一個溫度影響深度范圍，該深度范圍與環(huán)境溫度、土質(zhì)條件以及持續(xù)時間相關(guān)，這與實(shí)際多年凍土的地溫分布規(guī)律相同。

2 模型試驗(yàn)有限元模擬

2.1 土基模型參數(shù)及邊界條件

2.1.1 土質(zhì)參數(shù)及初始條件

土基模型試驗(yàn)利用ABAQUS 有限元模擬軟件創(chuàng)建二維土基水熱耦合溫度場模型，采用DC3D20 計算單元。模型土參數(shù)按照表1 所示設(shè)置。沿深度每層土的初始溫度可以通過試驗(yàn)測得，如表3 所示。

表3 不同深度土層的初始溫度Tab.3 Initial temperature of soil layers at different depth

2.1.2 土基模型邊界條件

室內(nèi)模型上邊界為開放條件，模型槽四周及底部均設(shè)置絕熱材料，可視為絕熱條件。因模型試驗(yàn)位于室內(nèi)，故忽略太陽輻射因素。模型中通過設(shè)置對流熱交換系數(shù)以實(shí)現(xiàn)對流換熱條件的定義，其中每階段對流熱交換系數(shù)由降溫系統(tǒng)吹風(fēng)風(fēng)速與吹風(fēng)時間決定[20]，其中風(fēng)速即為室內(nèi)降溫系統(tǒng)吹風(fēng)的平均風(fēng)速，取8 m/s。

2.2 模型試驗(yàn)與有限元模擬結(jié)果分析

在與室內(nèi)試驗(yàn)相同環(huán)境溫度及持續(xù)時間內(nèi)，采用有限元模型計算土層深度范圍溫度隨時間的變化規(guī)律。28 天內(nèi)0.1、0.4、1.0 m 深度處，土層溫度隨時間變化的有限元結(jié)果與模型試驗(yàn)結(jié)果對比如圖6 所示。

圖6 土層溫度隨時間變化對比Fig.6 Comparison of soil temperature changes over time

由圖6 可得，0.1 m 和0.4 m 深度處的有限元模擬結(jié)果與試驗(yàn)所測溫度隨時間的變化規(guī)律基本一致。這說明應(yīng)用該水熱耦合模型模擬計算凍土的溫度場是可行的。但1.0 m 深度處土層溫度伴隨試驗(yàn)持續(xù)時間的延長，二者誤差逐漸增大，其原因是模型試驗(yàn)底部的隔熱效果與有限元設(shè)置的絕熱邊界有差距導(dǎo)致的。

3 水熱力耦合地基溫度應(yīng)力分析

3.1 水熱力耦合地基有限元模型概述

3.1.1 溫度場數(shù)學(xué)模型

由文獻(xiàn)[21-22]建立的受凍融影響的物質(zhì)遷移方程，并結(jié)合三維非穩(wěn)態(tài)熱傳導(dǎo)方程、三維質(zhì)量遷移方程以及未凍水含量關(guān)于溫度的函數(shù)關(guān)系式，且計算中僅考慮水分遷移與凍土相變，可以將凍土溫度場的計算轉(zhuǎn)變?yōu)榍蠼馊S非線性方程組。

3.1.2 溫度應(yīng)力數(shù)學(xué)模型

通常溫度變化速度較慢，其熱平衡狀態(tài)每一段可以看作是均勻的，因此只需要對每個平衡狀態(tài)下的溫度分布進(jìn)行求解，然后再將其代入含溫度應(yīng)變的本構(gòu)方程就可以得出溫度應(yīng)力的分布，道面結(jié)構(gòu)選用層狀彈性理論，彈塑性本構(gòu)模型選用摩爾-庫倫本構(gòu)[23]。

3.1.3 模型結(jié)構(gòu)

結(jié)合東北大小興安嶺多年凍土區(qū)的地基典型斷面[24]，建立了多年凍土區(qū)跑道二維模型，如圖7 所示，其中，跑道寬度設(shè)為40 m，且大興安嶺地區(qū)多年凍土層受外界氣溫的影響深度最大約為18 m[25]，因此模型土體總深度設(shè)為20 m。依據(jù)有關(guān)勘測，大興安嶺多年凍土區(qū)地基典型結(jié)構(gòu)層可近似為3 層，整個模型斷面可以分為4 個區(qū)域，其中，Ⅰ區(qū)為跑道道面結(jié)構(gòu)層，參照民用機(jī)場道面設(shè)計標(biāo)準(zhǔn)，道面結(jié)構(gòu)層組成材料的性質(zhì)如表4 所示，Ⅱ區(qū)與Ⅲ區(qū)土層均為亞黏土，其中Ⅱ區(qū)和Ⅲ區(qū)土含水率分別為20%、15%；Ⅳ區(qū)土層為碎石亞黏土[24]。

表4 跑道道面結(jié)構(gòu)層Tab.4 Structure layer of runway pavement

圖7 多年凍土跑道二維有限元模型Fig.7 Two-dimensional finite element model of permafrost runway

3.1.4 初始條件及邊界條件設(shè)置1）初始應(yīng)力場和初始溫度場

在水熱力耦合模型中既要計算應(yīng)力，又要考慮溫度變化，模型需要導(dǎo)入初始應(yīng)力場和溫度場，即考慮重力作用的地應(yīng)力，將該應(yīng)力場導(dǎo)入模型，并在應(yīng)力計算中通過設(shè)定分析步來導(dǎo)入計算結(jié)果中每一相應(yīng)時間所對應(yīng)的溫度場。模型初始溫度設(shè)為-1.6 ℃[24]，初始計算時間為4 月。

2）邊界條件

模型兩側(cè)與道肩距離長為40 m，故兩側(cè)邊界可設(shè)為絕熱[26]，且左右邊界水平位移固定為0。下邊界地溫梯度可取0.04 ℃/m[25]，且豎向位移限制為0。上邊界為開放邊界，并定義太陽輻射和對流換熱條件[27]，由于瀝青混凝土道面的蒸發(fā)耗熱相較太陽輻射可忽略不計，因此忽略其對土基溫度場的影響[28]。

3.2 水熱耦合地基有限元模型參數(shù)

3.2.1 熱參數(shù)

結(jié)合模型中各土層的相關(guān)熱參數(shù)，考慮水熱耦合效應(yīng)，運(yùn)用ABAQUS 中的材料生熱模擬UMATHT 子程序完成等效熱參數(shù)的定義[1,28]。其中Ⅰ區(qū)為瀝青混凝土跑道道面結(jié)構(gòu)層，該層保水性較差，故可不考慮水分遷移影響[27]，道面結(jié)構(gòu)層相應(yīng)的熱參數(shù)[19，28-29]，如表5所示。

表5 道面結(jié)構(gòu)層熱參數(shù)Tab.5 Thermal parameters of pavement structure layer

3.2.2 道面結(jié)構(gòu)層力學(xué)參數(shù)

參考民用機(jī)場瀝青道面設(shè)計標(biāo)準(zhǔn)，道面面層下各層力學(xué)指標(biāo)如表6 所示。其中，瀝青混合料膨脹系數(shù)取1×10-5[30]。各溫度影響下瀝青面層的彈性模量[31-32]如表7 所示。

表6 道面結(jié)構(gòu)層力學(xué)參數(shù)Tab.6 Mechanical parameters of pavement structure layer

表7 瀝青面層彈性模量Tab.7 Elastic modulus of asphalt surface layer

3.2.3 凍土地基力學(xué)參數(shù)

溫度影響下凍土地基的力學(xué)參數(shù)變化明顯[23]。不同溫度條件下的彈性模量、泊松比、黏聚力、摩擦角[33-34]可分別表示如下

式中：E 為土的彈性模量（MPa）；m 為彈性模量公式中的非線性指數(shù)，可取0.6；v 為土的泊松比；T 為土體溫度（℃）；c 為土的黏聚力（MPa）；φ 為土的摩擦角（°）；a1～a4和b1～b4均為與土質(zhì)相關(guān)的經(jīng)驗(yàn)系數(shù)，通常情況下ɑ1=28，ɑ2=0.40，ɑ3=0.15，ɑ4=14，b1=27，b2=-0.008，b3=0.09，b4=0.60，但當(dāng)土體溫度高于0 ℃時，b1～ b4均為0。

根據(jù)文獻(xiàn)[35]，凍土的剪脹角可以為0.5 倍的摩擦角，模型土在不同溫度條件下的力學(xué)參數(shù)如表8 所示。忽略土骨架的凍脹效果[36]，模型土膨脹系數(shù)α 可表示如下

表8 土層力學(xué)參數(shù)Tab.8 Mechanical parameters of soil layer

式中，Vs和Vw分別表示土骨架與凍結(jié)水的體積（m3）。

模型土Ⅱ～Ⅳ區(qū)在不同溫度下的膨脹系數(shù)如表9所示。

表9 土層膨脹系數(shù)Tab.9 Soil expansion coefficient

3.3 東北多年凍土區(qū)跑道地基溫度應(yīng)力

采用本文第3.1 節(jié)建立的凍土區(qū)跑道有限元模型，并考慮水熱力耦合效應(yīng)，計算東北多年凍土跑道結(jié)構(gòu)層和地基分別在第10 年、20 年和30 年內(nèi)1、4、7、10 月的水平溫度應(yīng)力。

3.3.1 水平應(yīng)力等值線圖

由于結(jié)構(gòu)對稱性，取1/2 跑道模型進(jìn)行分析。圖8分別給出了凍土跑道地基在第30 年內(nèi)1、4、7 月和10 月對應(yīng)的水平溫度應(yīng)力等值線圖，其中“+”號表示拉應(yīng)力，“-”號表示壓應(yīng)力，應(yīng)力單位為kPa。

圖8 第30 年跑道中部水平溫度應(yīng)力等值線圖Fig.8 Horizontal temperature stress contour map in the middle of the runway in the 30th year

由圖8 可以看出，由于氣溫變化和跑道材料膨脹系數(shù)差異的影響，各月橫斷面等溫線上部更復(fù)雜，且道面結(jié)構(gòu)與道肩連接部分最為突出；跑道中部的水平溫度應(yīng)力，在道面結(jié)構(gòu)層以拉應(yīng)力為主，而地基土層中大部分為壓應(yīng)力；在一年當(dāng)中，從1 月至4 月隨著地溫的增加，地基土凍脹減弱，水平溫度應(yīng)力也減小，直至7 月，地溫明顯升高，表層道面結(jié)構(gòu)由于熱脹使得水平溫度應(yīng)力增加，到10 月份地溫再次下降，熱脹減弱，水平溫度應(yīng)力減小，隨著氣溫降至0℃以下，至次年的1 月份，地基土逐漸發(fā)生凍脹，溫度應(yīng)力再次增大。

3.3.2 水平溫度應(yīng)力沿深度變化規(guī)律

為了研究水平溫度應(yīng)力沿跑道深度如何變化，選取跑道道面中心位置，在第30 年內(nèi)對應(yīng)4 個月份沿地基深度的溫度應(yīng)力，以及第10 年、20 年、30 年1 月份的地基溫度應(yīng)力分布，如圖9 所示。

圖9 水平溫度應(yīng)力隨深度的變化Fig.9 Changes of horizontal temperature stress with depth

通過圖9（a）可以看出，在跑道中心位置處，水平溫度應(yīng)力在跑道結(jié)構(gòu)上部變化十分顯著，而在下部地基深度5 m 左右趨于一致；1 月份溫度應(yīng)力變化最大，且沿跑道深度出現(xiàn)的最大拉應(yīng)力和最大壓應(yīng)力都明顯大于其他月份，其中1 月的最大拉應(yīng)力位于跑道表層結(jié)構(gòu)，可見多年凍土區(qū)跑道在冷季需要重點(diǎn)防止凍裂病害的出現(xiàn)。此外，從圖9（b）可以看出，第30 年水平應(yīng)力最大，說明隨著時間的增加，相同月份下道面結(jié)構(gòu)的水平應(yīng)力逐漸增大。

3.3.3 道面結(jié)構(gòu)水平應(yīng)力橫向分布規(guī)律

在一年當(dāng)中，1 月和7 月對應(yīng)的跑道水平溫度應(yīng)力較大（如圖9 所示），且在道面結(jié)構(gòu)層內(nèi)變化最明顯，所以將針對道面結(jié)構(gòu)內(nèi)的水平溫度應(yīng)力進(jìn)行深入分析。分別選取第30 年內(nèi)道面結(jié)構(gòu)層在1 月和7 月沿斷面橫向的水平溫度應(yīng)力分布，如圖10 所示。

圖10 第30 年道面結(jié)構(gòu)層水平應(yīng)力Fig.10 Horizontal stress of pavement structure layer in the 30th year

由圖10 可以看出，在1 月份，溫度拉應(yīng)力主要出現(xiàn)在跑道表面和瀝青上面層底，該應(yīng)力從道中沿橫向逐漸減小，瀝青下面層底以下均為壓應(yīng)力，且沿橫向變化不大，其中最大拉應(yīng)力位于跑道表面層中心位置，大小約為1.21 MPa，約是此時溫度對應(yīng)該層劈裂強(qiáng)度[34]的30%；在7 月份，溫度拉應(yīng)力主要出現(xiàn)在跑道基層以下，而壓應(yīng)力出現(xiàn)在表層和面層，跑道基層下部溫度拉應(yīng)力較大，最大溫度拉應(yīng)力出現(xiàn)在跑道基層，大小約為0.14 MPa，約占該層材料劈裂強(qiáng)度的28%。此外不難看出，不論是冬季還是夏季，跑道結(jié)構(gòu)層中大部分還是以溫度壓應(yīng)力為主，當(dāng)同時作用飛機(jī)荷載時，跑道結(jié)構(gòu)很容易發(fā)生破壞。因此，多年凍土區(qū)機(jī)場跑道不能忽略溫度引起的水分場和應(yīng)力場，尤其是在冬季和夏季，由于溫度和荷載的長期累積作用，要預(yù)防凍土地基的融沉和跑道面層網(wǎng)裂病害的出現(xiàn)。

4 結(jié)語

（1）采用模型試驗(yàn)?zāi)M了東北多年凍土區(qū)跑道地基的溫度變化，當(dāng)外界溫度降低時，上部土層隨環(huán)境溫度下降明顯，隨著凍土地基厚度的增加，其溫度受環(huán)境溫度影響顯著降低，與基于水熱耦合有限元模型計算得到的溫度規(guī)律一致。

（2）考慮水熱力耦合作用，采用有限元模擬了東北大興安嶺多年凍土區(qū)跑道地基30 年內(nèi)的溫度應(yīng)力變化規(guī)律，研究得出，位于跑道中心位置處的水平溫度應(yīng)力在跑道上部結(jié)構(gòu)層變化十分顯著，而在下部約5 m 深度趨于穩(wěn)定；一年之中溫度引起的較大拉應(yīng)力，在冬季出現(xiàn)在跑道中心面層，最大占該層材料劈裂強(qiáng)度的30%，在夏季最大拉應(yīng)力出現(xiàn)在跑道基層，占此層材料劈裂強(qiáng)度的28%。跑道結(jié)構(gòu)層中大部分以溫度壓應(yīng)力為主，由于氣候變暖趨勢，相同時期的溫度應(yīng)力也會逐年增加。因此，多年凍土區(qū)跑道的溫度應(yīng)力不能忽略。

（3）文中重點(diǎn)研究了環(huán)境溫度引起多年凍土跑道地基溫度應(yīng)力的變化規(guī)律，未涉及跑道地基變形特性的分析，由于溫度和荷載的長期累積，跑道道面結(jié)構(gòu)與地基的變形也與溫度作用息息相關(guān)，特別是夏季和冬季尤為顯著。因此，多年凍土跑道地基的溫度場、應(yīng)力場及變形場還有待于進(jìn)一步地深入研究。