鮮少華，李清鵬，羅振先

(1.武漢市政工程設(shè)計研究院有限責(zé)任公司， 湖北 武漢 430023；2.河南大學(xué) 土木建筑學(xué)院， 河南 開封 475004)

當(dāng)季節(jié)凍土區(qū)地下水埋深較淺時，在冬季凍結(jié)過程中水分會在溫度梯度的作用下遷移至路基土體中，引起凍脹變形，繼而產(chǎn)生路面開裂破壞；融化過程中，大量的遷移水不易消散，造成土體強(qiáng)度降低，進(jìn)而引起路面沉陷等病害。溫度梯度的變化導(dǎo)致水分遷移的過程即是一種水熱耦合的過程。水熱耦合理論一直是國內(nèi)外學(xué)者的研究熱點(diǎn)，Harlan[1]通過對比分析非飽和土與部分凍結(jié)土中的水分遷移規(guī)律，最早建立了描述水熱遷移的數(shù)學(xué)模型。隨后，Taylor等[2]、Jame等[3]、Newman等[4]不同程度地改進(jìn)和發(fā)展了水熱遷移模型。Mu等[5]引入固體力學(xué)理論，在計算體應(yīng)變時考慮了水分凍脹引起的土體應(yīng)變，提出了水熱力三場耦合模型。毛雪松等[6]利用水熱力耦合模型計算了青藏公路唐南段路基溫度場、濕度場和應(yīng)力分布規(guī)律，得到的結(jié)果能夠體現(xiàn)道路凍脹的主要誘因。張玉芝等[7]聯(lián)合土體凍結(jié)過程中的熱傳導(dǎo)方程和熱彈性力學(xué)理論建立了高鐵路基的變形和內(nèi)力計算模型，預(yù)測了路基服役期內(nèi)在凍融作用下可能產(chǎn)生的凍脹量和融沉量。邰博文等[8]在凍土水熱耦合方程的基礎(chǔ)上通過土體凍脹率與含冰量的關(guān)系改進(jìn)了土體凍脹模型，并利用寒區(qū)路基變形監(jiān)測數(shù)據(jù)驗(yàn)證了模型的準(zhǔn)確性。Zheng等[9]將一維凍脹模型擴(kuò)展到了二維情形并考慮了凍結(jié)速率和約束應(yīng)力的影響。Li等[10-11]在水熱力三場耦合模型將凍融循環(huán)作用對土體的影響考慮了進(jìn)去，并利用室內(nèi)試驗(yàn)和現(xiàn)場監(jiān)測進(jìn)行了驗(yàn)證。褚志成等[12]借助溫度場-應(yīng)力場耦合模型研究了多年凍土區(qū)邊坡在凍融作用下的局部穩(wěn)定性。呂龍等[13]在Harlan水熱耦合模型的基礎(chǔ)上考慮水分遷移和冰水相變，認(rèn)為水分遷移是導(dǎo)致基坑工程冬季凍脹破壞的主要原因之一。

綜上，國內(nèi)外諸多研究者在土體凍脹機(jī)理和水熱耦合理論方面已經(jīng)進(jìn)行了大量研究工作，并且取得了一定的成果，但上述研究均未考慮有外界水分持續(xù)補(bǔ)給的情況。鑒于此，基于現(xiàn)有的水熱耦合理論，針對季節(jié)凍土區(qū)路基填土的水熱遷移問題，采用室內(nèi)試驗(yàn)和數(shù)值計算的手段，研究有外界持續(xù)補(bǔ)水條件下土體內(nèi)部的溫度場和濕度場變化規(guī)律。

1 水-熱耦合理論

1.1 基本假定

基于現(xiàn)有理論[14]，作以下基本假定：(1) 土體均勻連續(xù)，屬于各向同性體；(2) 土中無鹽分影響；(3) 凍結(jié)區(qū)和未凍區(qū)水分僅以液態(tài)形式進(jìn)行遷移；(4) 土顆粒不可壓縮；(5) 已凍土和未凍土為彈性體；(6) 無溶質(zhì)遷移；(7) 水分遷移符合Darcy定律。

1.2 溫度場方程

忽略對流傳熱和熱輻射，考慮土體中冰水相變的熱傳導(dǎo)方程表示為：

(1)

1.3 水分場基本方程

利用水動力學(xué)模型來描述土體凍結(jié)過程中水分的遷移過程。對于飽和土和非飽和土，水分遷移方程均可表示如下：

(2)

式中：θw為水分總體積含量，即θw=θu+θiρi/ρw，其中θu為未凍水的體積含量；ρw為水的密度；kx，ky為水平向和垂直向的滲透系數(shù)；Ψ為土水勢，參考文獻(xiàn)[3]將土體凍結(jié)和融化過程中的土水勢均定義為Ψ=Ψm+Y，其中Ψm為基質(zhì)勢，Y為重力勢。

(3)

1.4 水熱聯(lián)系方程

上述公式(1)和公式(3)共有T，θu和θi三個未知量，只有兩個方程，因此需要補(bǔ)充一個聯(lián)系方程才能解方程組。聯(lián)系方程為未凍水含量與溫度變化的關(guān)系，如下所示：

θu=θ0|Ttrans|B|T|-B

(4)

式中：θ0為初始含水率Ttrans為冰水相變溫度；B為與土質(zhì)相關(guān)的常數(shù)。

給定式(1)、式(3)和式(4)邊界條件和初始值就組成了路基工程水熱耦合的數(shù)學(xué)模型。

將未凍水含量θu對時間和空間求偏導(dǎo)得：

(5)

將水分遷移方程代入熱傳導(dǎo)方程并聯(lián)系上式可消去含冰率θi：

(6)

式(6)實(shí)現(xiàn)了方程組的解耦。將上式與水熱聯(lián)系方程(4)聯(lián)立即可求得溫度和水分的分布。

凍結(jié)區(qū)和未凍區(qū)之間有一個固態(tài)和液態(tài)分界界面，而且這個界面不固定。為了將該界面簡化為固定狀態(tài)，引入Heaviside階梯函數(shù)H(T)[15]利用COMSOL Multiphysics軟件內(nèi)部自帶的表達(dá)式表示如下：

H(T)=flc2hs(T-Ttrans,dT)

(7)

式中：Ttrans為相變點(diǎn)溫度；dT為轉(zhuǎn)變溫度間隙取Ttrans=-0.3℃，dT=0.3℃，建立的函數(shù)變化曲線如圖1所示。

圖1 階梯函數(shù)曲線

根據(jù)式(7)，在整個求解區(qū)域內(nèi)表示未凍水含量為：

θu=θ0|Ttrans|B|T|-B[1-H(T)]+θ0H(T)

(8)

容積熱容為：

C=Cf+(Cu-Cf)H(T)

(9)

式中：Cu為未凍土容積熱容量；Cf為已凍土容積熱容量。

滲透系數(shù)為：

k=kf+(ku-kf)H(T)

(10)

式中：ku為未凍土滲透系數(shù)；kf為已凍土滲透系數(shù)。

擴(kuò)散系數(shù)為：

D=Df+(Du-Df)H(T)

(11)

式中：Du為未凍土水分?jǐn)U散系數(shù)；Df為已凍土水分?jǐn)U散系數(shù)。

參考以往研究成果[2,10]，定義阻抗系數(shù)I:

I=1010θi

(12)

那么滲透系數(shù)和水分?jǐn)U散系數(shù)定義如下：

(13)

(14)

(15)

(16)

2 土體開放系統(tǒng)單向凍結(jié)試驗(yàn)

為驗(yàn)證土體開放系統(tǒng)中水熱耦合理論的有效性，選取路基填土制作圓柱形試樣，開展開放系統(tǒng)條件下單向凍結(jié)室內(nèi)試驗(yàn)。試樣直徑為38 mm，高度為76 mm。試樣初始含水率為20.2%，初始干密度為1.57 g/cm3。試樣制作完成后，用防水薄膜將試樣密封，僅保留底面以方便補(bǔ)水。試樣放置在圓柱形筒中，并安裝在如圖2所示的凍脹裝置中。該裝置能夠從上向下單向凍結(jié)試樣，最低溫度可達(dá)-40℃。此外，凍結(jié)過程中保持試樣底部恒溫在1℃左右，方便從外部給試樣補(bǔ)水。沿試樣高度方向每隔1.5 cm插入一個小型熱電阻，實(shí)時監(jiān)測試樣內(nèi)部溫度變化。凍結(jié)試驗(yàn)共進(jìn)行8 h，試樣取出后，將土柱沿高度方向分層切片取樣烘干測含水率。試樣不同高度處的溫度隨凍結(jié)時間的變化及溫度分布如圖3所示?？梢钥吹皆嚇釉诤銣貎鼋Y(jié)的情況下，前2 h溫度劇烈變化，之后逐漸穩(wěn)定，凍結(jié)5 h后試樣全部達(dá)到0℃以下。并且隨著凍結(jié)時間的增加溫度分布呈現(xiàn)由非線性分布向線性分布過渡的趨勢，表明溫度梯度逐漸穩(wěn)定。

圖2 土體凍融循環(huán)試驗(yàn)原理圖

3 水熱耦合理論的驗(yàn)證

3.1 計算模型及參數(shù)

本次研究借助于第二節(jié)中開放系統(tǒng)土柱單向凍結(jié)試驗(yàn)的溫度和水分變化結(jié)果驗(yàn)證水熱耦合理論在開放系統(tǒng)季節(jié)凍土區(qū)應(yīng)用的合理性。采用有限元軟件COMSOL Multiphysics根據(jù)圓柱土樣尺寸建立有限元模型如圖4所示，計算模型為二維軸對稱模型，其中高度為7.6 cm，半徑為1.9 cm。

圖3 凍結(jié)過程中不同時刻試樣的溫度分布

圖4 有限元計算模型

將試驗(yàn)開始時室內(nèi)溫度作為試樣的初始溫度，即T0=17.5℃。試樣初始含水率為20.2%，因此初始體積含水率設(shè)定為θ0=0.32。

試樣頂部即冷端(AD)溫度邊界條件按照試樣實(shí)際的溫度變化(冷卻溫度為-30℃時)進(jìn)行設(shè)置。根據(jù)圖5可以發(fā)現(xiàn)試樣頂部(h=7.6 cm)溫度變化趨勢為先線性下降然后逐漸穩(wěn)定，為簡化計算過程，將試樣頂部溫度變化擬合為分段函數(shù)即式(17)，擬合結(jié)果如圖5所示。

(17)

圖5 試樣頂部溫度變化

AB為軸對稱邊界，CD為絕熱邊界，BC為暖端恒溫邊界，溫度為試驗(yàn)過程中的實(shí)際溫度，Tb=2℃。對于水分邊界，AD和CD均為不透水邊界，AB為軸對稱邊界，BC為補(bǔ)水邊界，含水率恒定且為試樣的飽和體積含水率，即θs=0.42。參考徐學(xué)祖等[16]總結(jié)的亞黏土熱參數(shù)，綜合取值如表1所示[17]。那么可以得到土體的未凍水含量隨溫度的變化曲線如圖6所示。

表1 土體熱參數(shù)取值

圖6 未凍水含量隨溫度的變化曲線

COMSOL Multiphysics軟件自帶有系數(shù)型偏微分方程接口，可以將上述計算條件、參數(shù)及水熱耦合方程導(dǎo)入軟件中對方程重新編譯并進(jìn)行計算。

3.2 結(jié)果與分析

根據(jù)室內(nèi)試驗(yàn)結(jié)果，試樣凍結(jié)持續(xù)5 h后，內(nèi)部溫度全部達(dá)到零下并逐漸穩(wěn)定。因此以試樣凍結(jié)時間分別為1 h、2 h、4 h和5 h的溫度變化、凍結(jié)鋒面位置和水分分布對計算結(jié)果進(jìn)行驗(yàn)證。

(1)試樣溫度變化對比。試樣在凍結(jié)過程中的溫度變化云圖如圖7所示?？梢悦黠@看到隨著凍結(jié)時間的持續(xù)試樣凍結(jié)區(qū)逐步向下發(fā)展。數(shù)值計算時試樣底部保持恒溫2℃，因此圖7中凍結(jié)5 h后試樣的底部仍高于0℃。而試驗(yàn)過程中雖然底板溫度也保持在2℃左右，但是凍融循環(huán)試驗(yàn)系統(tǒng)中酒精循環(huán)液的冷卻效率明顯高于其保持恒溫的效率，因此出現(xiàn)5 h后試樣被完全凍結(jié)。

圖7 試樣凍結(jié)過程中的溫度變化云圖(單位：℃)

假定試樣凍結(jié)鋒面處土體的溫度為0℃，那么試樣在不同凍結(jié)時間計算凍結(jié)鋒面高度和實(shí)際凍結(jié)鋒面高度對比如圖8所示，其中試樣實(shí)際的凍結(jié)鋒面用白線標(biāo)記。可以看到凍結(jié)持續(xù)1 h后凍結(jié)鋒面的計算高度和實(shí)際高度均為55 mm左右；凍結(jié)2 h后高度均為30 mm左右；凍結(jié)4 h后凍結(jié)鋒面的高度出現(xiàn)一定差異，實(shí)際高度為20 mm左右，而計算高度為15 mm左右；凍結(jié)5 h后室內(nèi)試驗(yàn)的試樣已經(jīng)被完全凍結(jié)，而模擬結(jié)果的凍結(jié)鋒面的高度為5 mm左右。因此，在試樣完全凍結(jié)以前凍結(jié)鋒面的計算高度基本上和實(shí)際高度相差不大。

試樣不同高度處的計算溫度和實(shí)測溫度隨凍結(jié)時間的變化對比如圖9所示?？梢钥吹?，計算溫度和實(shí)測溫度變化趨勢一致，前2 h溫度迅速降低，凍結(jié)5 h后逐漸穩(wěn)定，且數(shù)值也比較接近。上述結(jié)果表明利用水熱耦合理論得到的開放系統(tǒng)土體凍結(jié)過程中的溫度計算結(jié)果是可信的。

圖8 試樣凍結(jié)過程中凍結(jié)鋒面的變化

圖9 試樣不同高度計算溫度和實(shí)測溫度變化對比

(2)試樣含水率變化對比。圖10為不同凍結(jié)時間試樣總體積含水率的計算值與實(shí)測值的對比。由于室內(nèi)試驗(yàn)中5 h后試樣即被完全凍結(jié)，而數(shù)值計算中溫度邊界條件的影響，土柱不會完全凍結(jié)。因此僅研究凍結(jié)時間分別為1 h、2 h和4 h時試樣內(nèi)部的總體積含水率分布狀況。所述總體積含水率包含了未凍水含水率和冰晶融化為液態(tài)水后的體積含水率。從圖10中可以看到計算含水率與實(shí)測含水率僅在試樣頂部分布差異較大其余位置數(shù)值大小接近，整體的變化趨勢基本一致，均呈現(xiàn)試樣上部含水率高于下部并在中部發(fā)生突變，并且計算得到的含水率突變位置與實(shí)際的突變位置基本一致。將圖10與圖8對比，可以發(fā)現(xiàn)含水率突變位置即凍結(jié)鋒面位置，凍結(jié)鋒面以上為凍結(jié)區(qū)，凍結(jié)鋒面以下為未凍區(qū)。由于土體中的水分是由未凍區(qū)向凍結(jié)區(qū)遷移的，因此試樣上部含水率要高于下部。上述結(jié)果表明，利用水熱耦合理論計算土體的含水率變化是符合實(shí)際的，可以應(yīng)用于季節(jié)凍土區(qū)受地下水遷移影響的巖土體的溫度和含水率計算中。

4 結(jié) 論

(1) 開放系統(tǒng)中土體的單向凍結(jié)試驗(yàn)表明，在凍結(jié)過程中隨著凍結(jié)時間的增加，沿試樣高度方向的溫度分布呈現(xiàn)由非線性分布向線性分布過渡的趨勢。

(2) 室內(nèi)試驗(yàn)結(jié)果和水熱耦合模型的計算結(jié)果均顯示，土體在凍結(jié)期間未凍區(qū)水分在溫度梯度作用下向已凍區(qū)遷移，并在凍結(jié)鋒面處水分聚積，含水率發(fā)生突變。

圖10 不同凍結(jié)時間試樣的含水率分布

(3) 建立的開放系統(tǒng)條件下的土體水-熱耦合理論模型能夠較好地預(yù)測土體凍結(jié)過程中的溫濕度變化，可以應(yīng)用于實(shí)際工程的溫度和含水率計算。