余自業(yè)，陳家豪，徐存東，3，4，李 準(zhǔn)，徐 慧，王海若，曹 駿

（1. 寧夏水利工程建設(shè)中心，寧夏 銀川 750000； 2. 華北水利水電大學(xué)水利學(xué)院，河南 鄭州 450046； 3. 浙江省農(nóng)村水利水電資源配置與調(diào)控關(guān)鍵技術(shù)重點(diǎn)實(shí)驗(yàn)室，浙江 杭州 310018； 4. 河南省水工結(jié)構(gòu)安全工程技術(shù)研究中心，河南 鄭州 450046）

0 引 言

在我國西北地區(qū)渠道表面大多采用混凝土襯砌的形式進(jìn)行防護(hù)，襯砌渠道有顯著的防滲作用，能夠減少了70%～90%的滲漏損失，是渠道防滲的重要措施之一［1］。由于寧夏固海灌區(qū)大部分為季節(jié)性凍土，襯砌渠道普遍出現(xiàn)受凍害、水體滲漏等現(xiàn)象，且灌區(qū)內(nèi)部分渠基土具有凍脹性的濕陷性黃土，遇水后在自重或受外力等情況易發(fā)生下沉變形［2］，嚴(yán)重影響灌區(qū)的安全穩(wěn)定運(yùn)行。因此迫切需要研究襯砌渠道凍脹破壞的機(jī)理，從而抓住問題的關(guān)鍵提出有效的治理方法。

面對(duì)襯砌渠道凍脹破壞的情況，相關(guān)學(xué)者已經(jīng)從凍脹機(jī)理、力學(xué)模型等方面開展了工作。STEPHEN 等［3］通過現(xiàn)場(chǎng)觀察和試驗(yàn)結(jié)果指出土壤凍脹的原因以及飽和土、黏土的水-熱-力耦合基本原理；LI 等［4］研究建立了混凝土襯砌渠道凍脹損傷的數(shù)值模型，通過等效的熱傳導(dǎo)方程精確模擬溫度場(chǎng)，同時(shí)對(duì)其采用有限元結(jié)合熱耦合的方法進(jìn)行損傷分析，模擬裂紋的發(fā)展和渠道的破壞模式，進(jìn)而揭示了其凍害機(jī)理。銀英姿等［5］通過實(shí)地試驗(yàn)，討論了在保溫條件下渠道凍脹破壞的影響因子的變化過程和基本規(guī)律；王玉寶等［6］利用渠道凍脹規(guī)律以及力學(xué)模型分析傾角與渠道抗凍脹能力之間的關(guān)系；肖雯等［7］將梯形襯砌渠道與地下水位影響因子結(jié)合，總結(jié)出計(jì)算襯砌渠道不同部位的凍脹力公式以及凍結(jié)強(qiáng)度和深度分布的方法；WANG 等［8］基于凍土的熱-水分-力耦合理論，用COMSOL 軟件對(duì)30 組不同截面的襯砌渠道進(jìn)行了研究，其研究結(jié)果揭示了襯砌渠道凍脹破壞的尺寸效應(yīng)；石姣等［9］建立了不同材料的襯砌渠道本構(gòu)模型，通過利用ABAQUS 熱力耦合模擬并比較了不同襯砌材料防止渠道凍害的效果。然而在這些襯砌渠道數(shù)值模型研究中，缺少了基于濕陷性黃土地基與渠道襯砌之間接觸作用的研究以及對(duì)高地下水位的考慮較少，因此本次研究加入濕陷性黃土地基這一特征條件，通過建立三相力學(xué)模型對(duì)渠道凍脹過程進(jìn)行模擬和分析。

鑒于此，本文基于凍土地區(qū)襯砌渠道受凍害及破壞的基本特點(diǎn)，考慮了濕陷性黃土地基這一特殊條件以及襯砌與基土的相互作用，以干渠改造前后斷面為梯形和弧底梯形結(jié)構(gòu)，建立三相相變?yōu)榛A(chǔ)的渠道凍脹力學(xué)模型，并將基土和渠道視為均勻的各向同性體，利用ABAQUS 軟件構(gòu)建襯砌渠道數(shù)值模型，旨在探明襯砌渠道各物理場(chǎng)的具體分布以及變形規(guī)律，并通過分析不同受力破壞模式進(jìn)而揭示其凍脹機(jī)理，以期在類似條件下減少灌區(qū)渠道受凍損傷和保護(hù)提供理論方法。

1 襯砌渠道凍害數(shù)學(xué)模型

1.1 渠道概況

灌區(qū)深居內(nèi)陸，多年平均氣溫為12 ℃，一般在9月至次年6月左右出現(xiàn)霜凍，一月份平均氣溫最低達(dá)到-15.5 ℃，土壤凍結(jié)期在120 d 以上。在本次研究對(duì)象為固海灌區(qū)的東一支干渠，原梯形渠道以挖、填方間斷形式布置，采用混凝土襯砌進(jìn)行防護(hù)。渠底寬200 cm，渠高250 cm，渠坡坡度為1∶1.5，混凝土襯砌厚度為6 cm，地表被風(fēng)積沙覆蓋，兩側(cè)砌護(hù)脫落、塌滑嚴(yán)重，渠道上部壤土等均為凍脹性土，場(chǎng)地為Ⅱ級(jí)自重濕陷，易發(fā)生滲透變形，需對(duì)其采取適當(dāng)防滲補(bǔ)缺措施。改造后渠道斷面為弧底梯形，兩岸護(hù)坡為8 cm 混凝土預(yù)制板，為加強(qiáng)防滲能力，在預(yù)制板下設(shè)3 cm 水泥砂漿，增設(shè)4 cm 苯板及一布一膜，干渠改造前后橫斷面如圖1所示。

圖1 渠道斷面圖（單位：mm）Fig.1 Channel cross-sectional view

梯形襯砌渠道渠底板兩端由于受到渠坡襯砌板的限制，使得底板兩端變形小而中間變形大；渠坡板由于受到兩端的約束限制，導(dǎo)致其在中間至靠近坡腳的位置變形大，兩端小；由于坡腳受到兩端襯砌板的限制以及對(duì)土體凍脹的限制，導(dǎo)致坡腳的變形量和凍脹量較小，此時(shí)渠坡板變形可視為偏心受壓情況，渠底板變形可視為壓彎情況［10］。而弧底梯形復(fù)合襯砌法向和切向凍脹力渠道渠底板分布均勻，在渠坡板上線性分布，在坡腳與弧底交接處達(dá)到最大值，渠道具體受力簡圖如圖2所示。

圖2 渠道受力簡圖Fig.2 Sketch of channel forces

1.2 數(shù)學(xué)模型構(gòu)建

1.2.1 基本假設(shè)

根據(jù)第一、第二凍脹理論［11，12］得出襯砌渠道凍脹是一個(gè)緩慢滲透的過程。主要有以下兩種方式：第一種方式為土體孔隙中的水分子被凍結(jié)的過程中發(fā)生相變并不斷釋放潛熱，進(jìn)而改變了土體的溫度場(chǎng)及凍結(jié)深度，從而影響凍脹量的變化；第二種方式為水分子在毛細(xì)力的引導(dǎo)下向結(jié)晶體處移動(dòng)，使結(jié)晶體的體積不斷擴(kuò)大，導(dǎo)致在上部土體中產(chǎn)生一層凍土層。當(dāng)溫度持續(xù)下降時(shí)，凍土層會(huì)隨著擴(kuò)大，進(jìn)而向襯砌板施力最終被破壞。

基于此，當(dāng)凍土層形成時(shí)，三相之間的耦合作用過程十分繁雜，使得精確模擬整個(gè)凍脹過程更加繁瑣，為方便在實(shí)踐中應(yīng)用，在構(gòu)建力學(xué)模型時(shí)針對(duì)凍結(jié)過程中主要特征進(jìn)行合理的簡化和假設(shè)［13，14］，其基本假設(shè)如下所示：①假設(shè)凍土是均勻、連續(xù)且各向同性體；②假設(shè)在土壤、受力條件相同的情況下相變溫度為常數(shù)，本次模型暫設(shè)為零度；③本次模型暫不考慮孔隙水壓力、土中含鹽量對(duì)渠道的影響；④假設(shè)本次模型為封閉系統(tǒng)，無水分補(bǔ)給。

1.2.2 渠基土熱傳導(dǎo)方程

由于該渠道緊鄰清水河，濕陷性壤土天然密度為1.78 g/cm3，含水量ω的范圍為15.2%～26.4%之間，場(chǎng)地處存在高地下水位的情況，受到天氣及灌溉的影響，地下水埋深最高在2～4 m。根據(jù)凍土層形成的基本理論，土體凍脹量與土體類型有關(guān)，地下水位的高低很大程度上決定了渠基土各部位的含水率，且不同土體類型不同含水率的熱膨脹系數(shù)不同，因此最終土體凍脹所造成的凍深也不同?？紤]到渠道不同部位的凍土、地下水埋深和凍深之間的關(guān)系，參考設(shè)計(jì)規(guī)范［15］給出了三者之間的關(guān)系，如圖3所示，并結(jié)合現(xiàn)場(chǎng)勘測(cè)得出該渠道改造前后凍脹情況如表1所示。

表1 渠道觀測(cè)結(jié)果Tab.1 Channel observation results

圖3 凍脹量-地下水埋深-凍深關(guān)系圖Fig.3 Frozen swelling amount-groundwater burial depth - freezing depth relationship

渠基土由于凍脹時(shí)間長且進(jìn)程緩慢，因此在進(jìn)行熱力學(xué)計(jì)算過程中，視為穩(wěn)態(tài)熱傳導(dǎo)進(jìn)行分析，則熱傳導(dǎo)方程見式（1）。

式中：T為溫度，應(yīng)滿足凍結(jié)邊界條件；λ下標(biāo)x、y分別為凍土在水平和豎直方向的導(dǎo)熱系數(shù)。

1.2.3 材料本構(gòu)方程

（1）凍土本構(gòu)方程。渠基土受凍膨脹是由于土體中的孔隙水在負(fù)溫條件下被凍結(jié)且結(jié)晶體體積不斷擴(kuò)大，加之溫度與應(yīng)力的影響導(dǎo)致土體受剪切破壞。為了方便在實(shí)際工程中進(jìn)行應(yīng)用，在ABAQUS 中將渠基土視為滿足D-P 強(qiáng)度準(zhǔn)則的彈塑性材料。因此，考慮土體應(yīng)變問題的本構(gòu)方程見式（4）：

式中：L為相變潛熱；σ為應(yīng)力張量；ε為應(yīng)變張量；γxy為剪應(yīng)變；τxy、τyz、τzx為剪應(yīng)力；α為線膨脹系數(shù)；Δt為溫度變化量；E為彈性模量；μ為泊松比。

（2）混凝土本構(gòu)方程。在ABAQUS 軟件中提供了混凝土連續(xù)、基于塑性的損傷模型（CDP），可用于模擬混凝土受壓和受拉過程，雖然CDP 模型僅適用于各向同性的彈塑性材料，但由于它的失效機(jī)理很好地響應(yīng)了混凝土關(guān)鍵的特征。CDP 模型的本構(gòu)方程為［16］：

式中：d為損傷因子；Del0、Del分別為初始、受損后的彈性剛度；εpl為塑性應(yīng)變。

2 渠道凍脹數(shù)值模型建立

2.1 數(shù)值模型

根據(jù)渠道觀測(cè)結(jié)果簡化襯砌渠道數(shù)值模型，渠頂左右兩端的邊界取100 cm，渠基土體凍深按照表1 進(jìn)行布置。對(duì)于梯形渠道，設(shè)置襯砌板與渠基土的接觸行為應(yīng)包括黏結(jié)、滑移、分離狀態(tài)。襯砌板與渠基土之間的摩擦系數(shù)取為0.554，剪應(yīng)力上限取接觸面間的凍結(jié)強(qiáng)度，允許最大彈性滑移變形為單元特征尺寸的0.5%，法向行為選取“硬接觸”關(guān)系，即當(dāng)拉應(yīng)力超過設(shè)計(jì)凍結(jié)強(qiáng)度后視為接觸面分離。對(duì)于弧底梯形渠道，將混凝土襯砌板、水泥、復(fù)合土工膜視為一體，與渠基土接觸面的設(shè)置同上述相仿。在ABAQUS 中選擇熱力耦合模塊，并將計(jì)算過程視為穩(wěn)態(tài)，網(wǎng)格劃分采用六面體、線性劃分的方式，設(shè)C3D8T 為此次模型的單元類型，模型網(wǎng)格劃分如圖4所示。

圖4 渠道網(wǎng)格劃分圖Fig.4 Channel grid division diagram

2.2 材料參數(shù)選取

本模型將凍土、襯砌板以及苯板等都視為線彈性材料且各向同性，考慮到濕陷性黃土在不同的溫度下彈性模量不同，因此需要測(cè)定和計(jì)算每個(gè)部分的導(dǎo)熱系數(shù)和熱膨脹系數(shù)，其具體的參數(shù)設(shè)置結(jié)果，見表2、表3。

表2 凍土不同溫度下的物理參數(shù)Tab.2 Physical parameters of permafrost at different temperatures

表3 其他材料計(jì)算參數(shù)Tab.3 Other material calculation parameters

2.3 邊界條件

定義好材料屬性并分配給對(duì)應(yīng)的渠道部位后，需要對(duì)其邊界條件進(jìn)行設(shè)定。由于本次模型模擬襯砌渠道熱-力耦合的過程得出3 個(gè)物理場(chǎng)，因此針對(duì)3 個(gè)物理層分別對(duì)襯砌渠道施加溫度、位移約束，具體設(shè)置如下：

（1）溫度約束：模型的上邊界溫度取實(shí)際渠道表面溫度，下邊界溫度設(shè)為0 ℃，兩側(cè)邊界設(shè)成絕熱。

（2）位移約束：對(duì)于梯形渠道，將縱向約束條件U2=0加在渠床底部土體的下邊界處，對(duì)于弧底梯形渠道的下邊界施加水平約束U1=0和縱向約束U2=0；在兩個(gè)渠道兩側(cè)邊界上設(shè)置水平方向約束，即U1=0。

3 結(jié)果分析

3.1 溫度場(chǎng)分析

圖5 為改造前后襯砌渠道模型的溫度場(chǎng)云圖，從中可以發(fā)現(xiàn)渠道表面溫度最低，土層越深，溫度越高。陰坡的溫度差和凍深要比陽坡大，分布規(guī)律一致，從而導(dǎo)致渠道土體凍脹的不均勻程度增加。且改造后的弧底梯形渠道，由于增加了防滲保溫措施，可清晰地看出等溫線分布較為集中，保溫效果較好，能在一定程度上減輕土體凍脹程度，符合工程實(shí)際情況。

圖5 溫度場(chǎng)云圖Fig.5 Temperature field cloud map

3.2 位移場(chǎng)分析

圖6（a）為改造前梯形渠道模型位移場(chǎng)輸出的法向凍脹量，其沿渠道斷面展開依次為陰坡、底板、陽坡。從位移場(chǎng)的整體分布來看，最大法向凍脹量位于陰坡板中間位置附近，最小法向凍脹量位于渠底板中間，法向凍脹量按圖中從左至右的順序分別為3.87、2.16、3.04 cm，這與數(shù)理模型得出的最大值的位置相似。這是由于渠底板兩端受到渠坡襯砌板的限制，使得底板兩端變形小而中間變形大；渠坡板由于兩端受到約束，導(dǎo)致其在中間的位置變形大，兩端小。由于坡腳受到兩端襯砌板的限制以及對(duì)土體凍脹的影響，其變形量和凍脹量較小。同時(shí)將法向凍脹量的模擬值和實(shí)際值進(jìn)行對(duì)比發(fā)現(xiàn)模型在陰坡出現(xiàn)了較為明顯的誤差，從整體來看模擬值要比實(shí)際值略低，這是由于本次模型為封閉系統(tǒng)，沒有考慮水分子遷移的情況。當(dāng)有外界水源補(bǔ)給時(shí)，隨著空氣溫度的持續(xù)下降，外界的水分子會(huì)被土體孔隙中的毛細(xì)力吸引至凍結(jié)面，導(dǎo)致該位置的凍結(jié)土厚度不斷增加。

圖6 法向凍脹量沿渠道斷面展開圖Fig.6 Normal frost heave expansion along the channel cross section

圖6（b）為改造后弧底梯形渠道模型位移場(chǎng)輸出的法向凍脹量。由圖可知法向凍脹量按圖中從左至右的順序分別為2.21、0.74、1.62 cm，在渠底與邊坡板的交接處發(fā)生凍脹的情況較小，在陰、陽坡中間附近出現(xiàn)最大值，同時(shí)由于渠底土壤含水率更高使得土體凍脹率較大，導(dǎo)致在底板中間處出現(xiàn)峰值。對(duì)比兩圖，陰坡總體法向凍脹量要大于陽坡，這是由于朝向不同導(dǎo)致渠道各部位受太陽照射和溫度的變化情況都不相同，使得陽坡凍結(jié)深度更小，加之混凝土材料在凍脹力、重力等外力條件下易受破壞，導(dǎo)致渠道中的水滲透進(jìn)土體中加劇土體的凍脹程度，使得渠道襯砌底板出現(xiàn)破壞，而改造后具有保溫防滲措施的弧底渠道在各部位的凍脹程度都明顯降低，與實(shí)際情況相符。

3.3 應(yīng)力場(chǎng)分析

（1）法向凍脹力分布。圖7為法向凍脹力具體分布圖，可以通過利用垂直于襯砌板下表面的應(yīng)力張量來獲得。根據(jù)圖7（a）分析得出渠坡上部存在少數(shù)法向凍脹力，自上而下逐漸增加，坡腳處可以看出存在應(yīng)力集中的現(xiàn)象，陰、陽兩坡的法向凍脹力分別達(dá)到2.99、1.8 MPa；渠底兩端受到諸多的限制作用，法向凍脹力較大，而中部較小，這是由于渠坡和渠底的襯砌板相互限制以及渠底板對(duì)土體凍脹的限制，導(dǎo)致在坡腳處雖然變形量以及凍脹量很小，但凍脹力卻很大。

圖7 渠道受法向凍脹力模擬值Fig.7 Simulated values of the channel subjected to normal frost swelling force

根據(jù)圖7（b）分析得，陰坡上端受壓，而后法向凍脹力由正轉(zhuǎn)負(fù)，說明襯砌板從受壓轉(zhuǎn)為受拉，襯砌板、苯板等與渠基土在相互作用力條件下有分開的趨勢(shì)，法向凍脹力最大值達(dá)到1.1 MPa。陰坡與渠底板交接處法向凍脹力出現(xiàn)峰值0.56 MPa，結(jié)合圖7（a）得出凍脹量小并不代表凍脹力也小，土壤不均勻凍脹的程度才是影響凍脹力大小的關(guān)鍵。渠底板法向凍脹力峰值為0.41 MPa，由于渠底為弧形，受到渠基土的頂托作用，因此法向凍脹力在渠底的整體分布為向內(nèi)凹。陽坡部分受拉，出現(xiàn)應(yīng)力集中的現(xiàn)象，其法向凍脹力最大值為0.87 MPa，且因?yàn)殛幤率艿降膬雒浟Υ笥陉柶?，所以渠道整體有偏向陽坡的趨勢(shì)。

（2）切向凍結(jié)力分布。圖8為切向凍結(jié)力分布圖，由圖8（a）所示梯形渠道切向凍脹力在渠坡呈現(xiàn)出上部小、下部大的分布規(guī)律，且陰坡的切向凍脹力要大于陽坡，方向均為斜向下。到了渠底呈線性分布，渠底兩側(cè)方向都朝向渠底中心處，這是由于渠床受兩端襯砌板的約束，陰坡提供拉應(yīng)力，陽坡提供壓應(yīng)力，且在渠底中心處的切向凍脹力接近于零。根據(jù)圖8（b）可知，在渠底中心出同樣出現(xiàn)切向力為0的情況，同時(shí)整個(gè)渠道受到指向陽坡方向的切向力。且在陰陽兩坡中下部均出現(xiàn)應(yīng)力峰值的情況，而在邊坡與弧底板接觸處由于切向作用不強(qiáng)，因此逐漸在降低。

圖8 渠道受切向凍脹力模擬值Fig.8 Simulated values of channel subjected to tangential frost heave

4 結(jié) 論

（1）根據(jù)襯砌渠道凍脹的不同形式以及基本特征，考慮到濕陷性黃土地基以及高地下水位的情況下，在ABAQUS 軟件中對(duì)陰坡、渠底、陽坡分別賦予不同的土質(zhì)系數(shù)參數(shù)以及對(duì)襯砌板賦予彈塑性本構(gòu)。通過進(jìn)行相應(yīng)渠道的熱力耦合數(shù)值模擬，得出的結(jié)果與實(shí)測(cè)數(shù)據(jù)相比吻合度好，能夠更好地了解預(yù)制混凝土板襯砌形式的渠道受凍脹時(shí)應(yīng)力以及變形基本規(guī)律。

（2）通過兩個(gè)渠道模型模擬的溫度場(chǎng)分布得出陽坡凍結(jié)深度和溫度梯度都要比陰坡??；從位移場(chǎng)分布得出襯砌板的不均勻變形是由于兩側(cè)渠坡的法向凍脹量較大，渠底較小。兩個(gè)渠道最大法向凍脹量在渠底板以及渠坡中心位置附近，且陽坡的凍脹程度要小于陰坡；從應(yīng)力場(chǎng)分布得出法向凍脹力在渠底兩端較大，而中部較小，特別地弧底梯形渠道的法向凍脹力在渠底處呈現(xiàn)“凹”型分布；切向凍脹力在渠底呈線性分布，且在渠底中心處接近于零，與渠道受力分析簡圖基本保持一致。

（3）通過兩個(gè)耦合模型模擬得出的3 個(gè)物理場(chǎng)結(jié)果對(duì)比發(fā)現(xiàn)改造后的剛?cè)嵋r砌渠道從保溫效果、凍脹情況以及受力大小3 個(gè)方面都得到了大幅度的改善，較好地避免了出現(xiàn)應(yīng)力集中的情況，模型整體可靠度較高，可為渠道抗凍脹設(shè)計(jì)提供參考。