江浩源，王正中，2，劉銓鴻，王 羿，譚志翔，孟曉棟

（1.西北農(nóng)林科技大學(xué) 旱區(qū)寒區(qū)水工程安全研究中心 旱區(qū)農(nóng)業(yè)水土工程教育部重點實驗室，陜西 楊凌 712100；2.西北生態(tài)環(huán)境資源研究院 凍土工程國家重點實驗室，甘肅 蘭州 730000）

1 研究背景

我國水資源時空分布不均，南澇北旱、南暖北寒現(xiàn)象普遍，北方寒旱區(qū)的灌區(qū)及引調(diào)水工程建設(shè)是保障旱區(qū)農(nóng)業(yè)用水及安全飲水的生命線工程［1］。作為主要輸水方式的渠道工程在溫-水-土-結(jié)構(gòu)等滲-凍互饋機制作用下［2］，常發(fā)生鼓脹、斷裂、脫空、失穩(wěn)滑塌等破壞［3-5］，嚴重影響了渠道的正常輸水功能及效益發(fā)揮。

太陽輻射是地表熱量的主要來源，是影響渠道溫度場的重要因素。因渠道走向、結(jié)構(gòu)型式、所處位置等原因，使渠道陰陽坡受到的太陽輻射不同，引起了渠道橫向溫度場的非對稱分布，產(chǎn)生了不穩(wěn)定的凍融狀態(tài)［6-8］。經(jīng)調(diào)查，陜西馮家山水庫總干退水渠［9］、新疆阜康灌區(qū)［10-11］及山東打漁張灌區(qū)［12］等部分東西走向渠段均存在顯著的陰陽坡效應(yīng)，陰陽坡凍深及凍脹變形差異明顯。太陽輻射是凍融循環(huán)的前提基礎(chǔ)，是引發(fā)渠道水分遷移、冰水相變等水-熱-力耦合凍脹破壞的內(nèi)在驅(qū)動力。太陽輻射導(dǎo)致凍結(jié)特征的不對稱、不同步變化是渠道凍脹破壞的重要原因［7-8］，探明其變化規(guī)律對了解渠道凍脹破壞的機理具有重要意義。

季節(jié)性凍土區(qū)襯砌渠道的陰陽坡效應(yīng)研究較少，王正中等［7-8］、安元等［6］分別將太陽輻射簡化為等效氣溫邊界條件，或采用HOTTEL 輻射模型建立了太陽輻射下渠道的溫度場計算模型，而模型中并未考慮渠坡對太陽光的遮蔽作用或基土凍結(jié)產(chǎn)生的相變潛熱、水分遷移等復(fù)雜的凍融循環(huán)作用。雖有學(xué)者建立了渠道水-熱-力耦合模型，但并未考慮渠道在太陽輻射作用下的陰陽坡效應(yīng)［2-3，13］，或僅是參考監(jiān)測數(shù)據(jù)賦予襯砌板溫度邊界來近似考慮太陽輻射［12］，不能反映太陽輻射的時空變化；同樣，多年凍土區(qū)的青藏公路、鐵路等工程陰陽坡效應(yīng)雖研究較多［14-15］，但亦存在上述問題。學(xué)者主要從現(xiàn)場監(jiān)測數(shù)據(jù)入手，結(jié)合凍土的地溫及路基變形，擬合太陽輻射量與地表溫度的計算公式［16］，認為土體中的不均勻熱狀況導(dǎo)致了差異性的凍融過程和不均勻沉降，使路基產(chǎn)生凍害［17-19］，但與季凍區(qū)上建筑物破壞機理有所差異?，F(xiàn)場監(jiān)測雖可準確評價凍害是否發(fā)生，但因凍土區(qū)環(huán)境惡劣、監(jiān)測成本高等原因，并不能全尺度的監(jiān)測出渠道系統(tǒng)內(nèi)溫度、水分和變形的分布規(guī)律及陰陽坡太陽輻射的空間效應(yīng)和晝夜溫差的時間效應(yīng)。從提升設(shè)計理論水平出發(fā)，需要建立季凍區(qū)考慮太陽輻射和晝夜溫差的渠道水-熱-力耦合凍脹模型，從理論上分析襯砌渠道的凍害發(fā)生原因，在已知外界環(huán)境下便可預(yù)測渠道基土的水、熱、力耦合變化，對渠道凍脹破壞進行科學(xué)分析預(yù)測、防控和設(shè)計。

目前已有考慮太陽輻射的結(jié)構(gòu)溫度場計算模型，且得到驗證。如Jin 等［20］基于光線追蹤算法，計算了任一時刻拱壩的溫度場及應(yīng)力場分布；Liu 等［21］和Chen 等［22］基于太陽輻射下構(gòu)件的非均勻溫度場試驗結(jié)果，建立了考慮太陽輻射的溫度場計算模型，并對大型穹頂結(jié)構(gòu)［23］和射電望遠鏡［24］的非均勻溫度場進行了分析。上述模型僅需考慮結(jié)構(gòu)吸收的太陽輻射能，而對于太陽輻射作用下的渠道凍脹分析而言，還需考慮溫度場、水分場、變形場的耦合作用，模型更為復(fù)雜。事實上在太陽輻射作用下的渠基土存在“夜凍晝消”“冬凍春融”現(xiàn)象，凍融循環(huán)頻繁，具有短時急變、分布不均勻、作用效應(yīng)復(fù)雜等特點，因此建立考慮太陽輻射的渠道水-熱-力耦合模型仍存在不少問題有待解決。

綜上，本文基于凍土水-熱-力耦合理論，考慮冰水相變及水分遷移作用，結(jié)合太陽輻射氣象參數(shù)模型，建立了可反映陰陽坡太陽輻射空間效應(yīng)和晝夜溫差時間效應(yīng)的渠道水-熱-力三場耦合模型，并結(jié)合原型監(jiān)測驗證了數(shù)值模型的準確性?；谀Ｐ?，通過系統(tǒng)分析渠道吸收的太陽輻射能、溫度場的時空變化規(guī)律，量化了渠道斷面的水、熱、變形差異，分析了寒旱區(qū)渠道凍脹破壞的內(nèi)在原因，并給出渠道凍脹數(shù)值計算中熱邊界的選取建議。該模型從理論上分析渠道凍脹的發(fā)展過程，不依賴建成后的現(xiàn)場監(jiān)測數(shù)據(jù)，通用性較強，可為寒旱區(qū)渠道的設(shè)計、維護與研究提供參考。

2 太陽輻射下的渠道凍脹分析模型

季凍區(qū)襯砌渠道在空氣熱對流、太陽輻射、襯砌板間熱輻射和環(huán)境熱輻射等作用下產(chǎn)生復(fù)雜的水-熱-力耦合凍脹現(xiàn)象。因此在水-熱-力耦合模型和太陽輻射模型基礎(chǔ)上，提出了渠道陰影計算方法，并基于輻射度算法對上述多種熱輻射建立了渠道邊界平衡方程，以作為溫度場計算的熱通量邊界條件，從而建立了太陽輻射下的渠道凍脹分析模型。因渠道屬線性工程，可按平面問題處理。

2.1 渠基土水-熱-力耦合凍脹模型

2.1.1 水-熱耦合控制方程 因土顆粒間隙較小而以考慮熱傳導(dǎo)為主，采用冰水相變修正后的傅里葉熱傳導(dǎo)方程如下［25-26］：

式中：T 為溫度，℃；Lf為冰水相變潛熱，kJ/kg；ρi為冰密度，kg/m3；θi為體積冰含量；Cp和λ分別為土體等效定壓熱容J/（kg·K）和等效導(dǎo)熱系數(shù)W/（m·K），其值由下述半經(jīng)驗公式估算［2，9，27］：

式中：下標s、w、i、a 分別代表土顆粒、水、冰及氣相。

變飽和多孔介質(zhì)內(nèi)的水分運動可用水頭型Richards 方程描述，添加冰相的方程如下［9，13，27］：

式中：C 為比水容量；h 為基質(zhì)勢，m；k 為土體滲透系數(shù)，m/s；i 為重力項。

采用van Genuchten 模型來描述未凍水含量與基質(zhì)勢、滲透系數(shù)的關(guān)系，方程如下［9，27］：

式中：α、m 為試驗擬合參數(shù)；θs、θr分別為飽和、殘余含水量；Se 為等效飽和度；k（Se）、ks分別為非飽和土、飽和土的滲透系數(shù)，m/s。

引入冰阻抗系數(shù)I 來近似估算因冰的存在而引起的凍土滲透系數(shù)k 的降低，方程如下［25-26］：

因顆粒表面能作用，凍土中的未凍水含量始終與溫度保持動態(tài)平衡，方程如下［28］：

式中：Wu為未凍水的質(zhì)量含水量；a、b 為實驗參數(shù)；T 為含水量為Wu時對應(yīng)的凍結(jié)溫度（Tf），℃。

2.1.2 應(yīng)力-應(yīng)變控制方程 基于增量彈塑性理論，凍土的應(yīng)力-應(yīng)變方程如下［18］：

式中：λp為塑性乘子；Q 為塑性勢函數(shù)，滿足相關(guān)聯(lián)流動法則，等于屈服面函數(shù)F，取M-C 準則。

考慮凍土垂直和平行溫度梯度的橫觀各向同性凍脹特征［9］，并依據(jù)溫度梯度方向?qū)崟r修正主凍脹方向［2］，坐標轉(zhuǎn)換至整體坐標系，凍脹應(yīng)變向量方程如下：

式中：ξ為凍脹方向的分配權(quán)重，取0.9；m、n 為溫度梯度向量的方向余弦；n0為凍土初始孔隙率。

聯(lián)立式（1）（4）（9）（10）構(gòu)成渠基凍土的溫度場、水分場和變形場的耦合控制方程。由式（1）（4）（9）可計算寒區(qū)渠道在外界環(huán)境作用下的溫度場和冰、水含量分布，結(jié)合式（12）得到凍脹應(yīng)變參數(shù)，采用式（10）便可求解得到渠道的應(yīng)力變形分布，從而實現(xiàn)渠基凍土的水-熱-力耦合凍脹分析。

2.2 渠道表面輻射換熱模型

2.2.1 太陽輻射模型 太陽輻射具有明顯的時空效應(yīng)，其位置由赤緯δ（°）和時角ω（°）組成的赤道坐標系及由太陽高度角αs和方位角γs組成的地平坐標系決定，各參數(shù)計算公式如下［6，29］：

式中：n 為日序數(shù)，其中1月1日為1，12月31日為365，表征不同季節(jié)日期的參數(shù)。

式中：t 為從0 時開始的小時數(shù)，該值表征當日時間的參數(shù)，上午為負，下午為正。

式中：φ 為地理緯度；γs偏東為負，偏西為正。

對北方寒旱區(qū)而言，Hottel 晴空模型因具有較好的適用性和通用性而被優(yōu)先采用［6，21，29］，可彌補國內(nèi)站點稀疏、數(shù)據(jù)較少等不足。該模型中太陽射線到達大氣層外切平面的太陽輻射強度G0由下式計算：

式中：Gsc為太陽常數(shù)，表示單位時間單位面積上所接收的太陽輻照度，取1367 W/m2。

任一斜面的太陽入射角i、太陽輻照度Gs（太陽直射輻照度Gtb和散射輻照度Gtd）公式如下［29-30］：

式中：β 為斜面傾角；γt為斜面方位角，指渠道表面法線在水平面上的投影與正南方向的夾角，面向東時為負，面向西時為正；斜面上的日照開始和結(jié)束時間以i=90°來計算［30］；α0、α1、k 為標準晴空大氣常數(shù)，根據(jù)文獻［6，29］計算。

2.2.2 渠道陰影計算方法 渠道邊坡相互遮擋產(chǎn)生的陰影分布影響了渠道表面的太陽輻射分布。在太陽高度角較小時，陰坡的遮擋使渠道局部形成陰影，且此陰影的邊界線與渠道的走向平行?；诖耍蓪⑶狸幤驴醋饕桓鶙U，將此桿在太陽下形成的陰影邊界點做渠道走向的平行線，此為渠道陰坡產(chǎn)生的陰影邊界。以北方東西走向渠道為例，為增加圖幅立體感，將方向逆時針旋轉(zhuǎn)，計算簡圖見圖1。

圖1 不同走向渠道陰影長度計算示意

圖中OA 為細桿，高度為渠深h，l1和l4分別為渠坡頂在水平面上的投影。若此時太陽位于渠道正南方，細桿將產(chǎn)生OB 的陰影，渠道陰影區(qū)域為l2—l4，l1—l2為受光區(qū)；若將渠道順時針旋轉(zhuǎn)一定角度，此值為面向太陽的渠道表面方位角γt。為便于圖幅描述，亦可表示為渠道走向不變，太陽逆時針旋轉(zhuǎn)此角度。此時細桿OA 將產(chǎn)生OC 的陰影，渠道陰影區(qū)域為l3—l4，陰影邊界離細桿距離為OD。因太陽高度角未發(fā)生變化，細桿產(chǎn)生的陰影長度一致，即OB=OC。綜合考慮渠道走向和太陽位置變化，其產(chǎn)生的陰影長度OD 由下式計算：

進一步結(jié)合渠道斷面進行陰影判定，如圖2所示。A、E 均為太陽光線可穿過位置，為簡化圖幅，僅以A 點為例介紹。D1、D2、D3為不同太陽位置下產(chǎn)生的陰影點，結(jié)合陰影長度OD1、OD2、OD3和渠道斷面幾何關(guān)系，確定不同時刻下渠道的受光區(qū)R，計算公式如下：

圖2 渠道陰影計算示意

聯(lián)立式（13）—式（18）及（21）（22）得到不同地區(qū)、時間下的不同走向、斷面大小的渠道表面太陽入射角和受光面積，二者共同決定渠道表面的日照時間，結(jié)合式（19）（20）得到各表面接受的太陽輻照度。2.2.3 襯砌板間輻射和環(huán)境熱輻射 因襯砌板間溫度及其與環(huán)境溫度不同而產(chǎn)生熱輻射，計算簡圖見圖3，任一點P 得到的S ′面的輻照度（Gm：W/m2）及周圍環(huán)境輻照度（Gamb：W/m2）見下式［31］：

式中：J ′、Jamb分別為S′面和環(huán)境（Samb）的輻射度，W/m2；n 和n′為外法線矢量；r 為面上兩點的距離矢量；n 為折射率，不透明物體取1；σ為Stefan-Boltzmann 常數(shù)，取5.67×10-8W/（m2·K4）；Famb為環(huán)境角系數(shù)；Tamb為環(huán)境溫度，K。2.2.4 基于輻射度方法的表面輻射換熱方程 輻射度J 由自身輻射及對輻照度（Gm、Gs、Gamb）的反射輻射組成，而每一點Gm又是其他可見點J 的函數(shù)，據(jù)此得到如下輻射平衡方程［31］：

式中：ρd為漫反射系數(shù)；ε 為發(fā)射率；T 為渠道表面溫度，K。

假設(shè)渠道表面為理想漫射灰體，發(fā)射率與吸收率α 相等，可計算得到渠道表面吸收的輻射量，并以此作為凍土水-熱-力耦合中熱模塊的第二邊界條件（熱通量邊界條件），方程如下：

圖3 角系數(shù)計算示意

3 數(shù)值模型

以新疆阜康某灌區(qū)（87oE、44oN）E-W 走向“陰陽坡效應(yīng)”顯著的大型輸水渠道為例，基于現(xiàn)場實測數(shù)據(jù)，聯(lián)立上節(jié)考慮太陽輻射及陰影遮蔽的渠道水-熱-力耦合方程，并采用COMSOL 求解因太陽輻射而導(dǎo)致的不均勻水、熱、變形分布，分析渠道凍害發(fā)生原因，并給出熱邊界的選取建議。

3.1 基本概況該梯形渠道修建于季節(jié)性凍土區(qū)，屬溫帶大陸性干旱氣候，多年平均溫度6.5 ℃，極端最低溫度-25～-32.2 ℃，多年最大凍深為1.72 m，平均凍深為1.45 m。渠基土中粒徑小于0.075 mm的土粒質(zhì)量占比大于10%，為凍脹敏感性土。

渠道設(shè)計流量為21 m3/s，加大流量26 m3/s，渠深3.0 m，正常水深2.6 m，底寬4.0 m，坡比1∶1.5。采用C20 現(xiàn)澆混凝土襯砌，板厚12 cm。姜海波等［10-11］于2013年11月5日—2014年3月28日對該渠道的溫度場（WS-4 型遙感土壤溫度計，測量深度距基土表面5～120 cm）、凍深（DTM-2 型凍土器）、水分場（土鉆取樣后烘干，測量距基土表面5～120 cm）及襯砌板凍脹量（水準儀）等進行監(jiān)測，監(jiān)測點布設(shè)及渠道有限元網(wǎng)格如圖4所示。凍深區(qū)域內(nèi)網(wǎng)絡(luò)高度最小為0.5 cm，最大為1 cm，最大時間步長120 s。其中A、B、C 為陽坡測點，分別位于1、2/3、1/3 設(shè)計水位處，D 為底板中心，E、F、G 為陰坡測點，位置同陽坡對稱。

圖4 有限元網(wǎng)格及數(shù)值模型示意（單位：m）

3.2 計算參數(shù)根據(jù)試驗段土工試驗［10-11］可知，淺層基土天然體積含水率為0.256～0.368，處于潮濕-飽和狀態(tài)，液性指數(shù)為0.53～0.82，塑性指數(shù)為5.8～6.4，屬低液限粉土。為簡化計算，計算參數(shù)近似取試驗均值，其中黏聚力為10.0 kPa，內(nèi)摩擦角為20°，干密度為1.6 g/cm3；結(jié)合基土土質(zhì)，其彈性參數(shù)與溫度的關(guān)系采用式（27）（28）計算［2，13，18］；土體各相、漫反射系數(shù)和發(fā)射率、凍結(jié)曲線及土水特征曲線的參數(shù)見表1［2，9，24，28，32］，混凝土參數(shù)見表2［2，20］。

表1 土體材料計算參數(shù)

表2 混凝土材料計算參數(shù)

3.3 邊界和初始條件渠道上表面溫度邊界包括熱輻射（式（26））及熱對流邊界，其中對流換熱方程如下：

式中：n 為渠道上表面法向向量；Tamb和T 分別為外界溫度和渠道表面溫度，℃；hc為對流換熱系數(shù)，W/（m2·℃），由風速來近似確定［2］，公式如下：

圖5 外界最高、最低溫度變化曲線

根據(jù)現(xiàn)場實測數(shù)據(jù)，渠頂和渠底風速近似取5 m/s 和1 m/s，并采用二次拋物線函數(shù)過渡。

該渠道已修建多年，基土溫度基本達到吞吐平衡，因此以年周期平衡溫度場作為計算的初始溫度場（11月1日）。根據(jù)氣象站溫度監(jiān)測值，采用正弦函數(shù)擬合溫度變化曲線［9，19］，公式如下：

式中：T0為年平均溫度，取6.5 ℃；B 為年溫度振幅，取22.5 ℃；tm為時間，月；α0為相角，取π。

選取試驗段渠道岸邊百葉箱每日溫度實測極值［10］，整理如圖5所示。日溫度變化曲線結(jié)合每日極值溫度采用正弦函數(shù)擬合［6］，以此作為計算分析的外界溫度，公式如下：

渠基下表面恒溫層溫度近似取年平均溫度6.5 ℃，距地表深度為11 m［33］。

渠道因防滲措施不當，滲漏較為嚴重，停水后，渠底下4.6 m 深處土體基本處于飽和狀態(tài)，水分場基本穩(wěn)定。此后按設(shè)計水位模擬通水一個月，并在凍結(jié)試驗開始前一周排干水分，以此作為計算分析的初始水分場。

因復(fù)合土工膜較?。?3 mm），采用無厚度土工膜滲阻計算，土工膜前區(qū)域為混凝土襯砌，膜后區(qū)域為渠基土，質(zhì)量守恒方程建立如下：

式中：下標b、u、d 分別代表土工膜、膜前和膜后區(qū)域；n 為邊界外法線矢量；k 為滲透系數(shù)，混凝土取1.2×10-9m/s，土工膜取1.0×10-11m/s；h 為壓力水頭，m；▽h為滲流域在邊界處的壓力梯度；R 為滲流阻抗，s；d 為厚度，m。

在渠道周邊施加法向位移約束，見圖4，滿足初始地應(yīng)力平衡。

4 結(jié)果分析

4.1 溫度場分析

4.1.1 渠道陰影分布、太陽輻射及溫度變化 因渠道為二維斷面，采用襯砌板吸收的太陽輻射量沿渠周分布來定量描述不同時刻渠道的陰影分布，并分析由此而導(dǎo)致的溫度差異。因東西走向渠道上、下午接收的太陽輻射量基本一致，且阜康地區(qū)比北京約晚2 h，因此取10∶00—14∶00 時段分析太陽輻射量；因襯砌板溫度存在累積作用，取10∶00—18∶00 時段分析。以12月10日為例，結(jié)果如圖6—8所示。其中，坐標原點從陽坡渠頂算起，下同。

由圖6 可知，從早上太陽初升，隨時間增加，太陽直射區(qū)域由陽坡逐漸擴大至渠底，陰影區(qū)域逐漸減少；隨太陽高度角的增加，太陽輻射量逐漸增大，至14∶00 時最大達374 W/m2，照射近1/2 渠周；而陰坡始終處于陰影區(qū)域，僅吸收太陽散射值，太陽輻射遠小于陽坡，輻射量差值達324.6 W/m2。

圖6 不同時刻襯砌吸收的太陽輻射沿渠周分布曲線

由圖7 可知，因太陽輻射區(qū)域及量值的增加，襯砌表面溫度隨之升高，陽坡溫度遠高于陰坡，最大差值為7.5 ℃，“陰陽坡效應(yīng)”明顯，且由太陽輻射引起的溫度累積效應(yīng)明顯；同時襯砌板陰、陽坡溫度均隨外界氣溫的降低而降低，且二者差值逐漸減小。

圖7 不同時刻襯砌表面溫度沿渠周分布曲線

由圖8 可進一步得知，襯砌表面溫度與外界氣溫和太陽輻射存在較好的相關(guān)性，但滯后明顯，二者共同造成了陰陽坡的橫向溫度差異。太陽輻射時段內(nèi)，陽坡溫度主要受太陽輻射影響，而陰坡則主要受外界氣溫熱對流作用；非太陽輻射時段內(nèi)，陰陽坡溫度主要受外界氣溫影響。因此以往分析中不考慮陰陽坡或?qū)⑻栞椛浒谕饨鐨鉁刂谢蛉〔町愋缘那辣砻鏈囟戎?，均不能科學(xué)合理的反映出陰陽坡溫度場的時空分布差異。

圖8 陰陽坡2/3 設(shè)計水位點溫度、輻射分布曲線

4.1.2 渠道陰陽坡襯砌板表面溫度變化 取監(jiān)測期內(nèi)陰陽坡襯砌板最低和日均溫度值來分析二者的橫向溫度差異，如圖9所示。

圖9 陰陽坡襯砌板表面溫度變化曲線

由現(xiàn)場實測值可知，襯砌板溫度隨外界氣溫的降低而降低，其變化規(guī)律與外界氣溫基本一致，但存在滯后性。陰、陽坡板最低溫度分別為-24.6 和-21.5 ℃，溫差最大為4.3 ℃。

由數(shù)值結(jié)果可知，最低溫度值略高于現(xiàn)場實測值，可能是由于太陽輻射模型中未考慮空氣透明度、陰天、降雪等因素引起的太陽輻射值偏高所致。同時溫度變化規(guī)律同現(xiàn)場實測基本一致，且陰陽坡板日均溫度差異明顯，最大為3.5 ℃，進一步表明本模型可較好的分析由太陽輻射的時空效應(yīng)和陰坡遮蔽作用而引起的襯砌板橫向溫度差異。

4.1.3 渠道陰陽坡凍深及基土溫度場變化 渠道表面溫度場的不對稱導(dǎo)致基土溫度場的差異性分布。以陰、陽坡2/3 設(shè)計水位點為例，將0 ℃等溫線作為凍深判別依據(jù)，其發(fā)展如圖10所示。

圖10 陰、陽坡凍深發(fā)展過程線

因晝夜溫差及太陽輻射作用，凍結(jié)前期陰、陽坡呈現(xiàn)出“夜凍晝消”的現(xiàn)象。由圖可知，凍深隨負積溫的累積而逐漸向下發(fā)展，“夜凍晝消”現(xiàn)象逐漸消失，其中陰坡于11月25日消失，陽坡于12月10日消失，二者前后差異15 d。結(jié)合圖5，至1月10日最低溫時刻，凍深并未發(fā)展至最大，后期雖存在大幅度升溫現(xiàn)象，但并未達到凍土的融化溫度，凍深依然向下發(fā)展；至3月16日，凍深最大，為135 cm；隨著溫度的繼續(xù)升高，陽坡于3月17日開始融化，陰坡于3月18日開始融化，二者前后差異1 d。外加太陽輻射作用，使得基土的融化速率遠大于凍結(jié)速率，至監(jiān)測期的最后一天3月28日，陰、陽坡基土表面分別融化40 和49 cm，凍深分別為97 和66 cm，差異明顯。

為進一步分析太陽輻射作用下的基土溫度場分布規(guī)律，取最大凍深時刻，如圖11所示。

圖11 3月16日渠道溫度場分布曲線

由圖可知，陰坡溫度明顯低于陽坡，溫度場呈不對稱分布；基土溫度隨外界氣溫（-8 ℃）變化存在滯后性，凍深以上溫度變化劇烈，以下則較為緩慢。結(jié)合凍深監(jiān)測結(jié)果可知，數(shù)值結(jié)果略小，最大差異為22 cm，主要原因在于太陽輻射模型未考慮陰天、降水對太陽光的阻擋所致。但反映出的凍深分布規(guī)律與監(jiān)測結(jié)果分布規(guī)律基本一致，僅在數(shù)值上存在微小差異，表明本文模型可較好的模擬太陽輻射作用下的基土溫度場分布規(guī)律。

4.2 水分場分析太陽輻射和環(huán)境溫度共同作用引起渠道的不對稱溫度場，勢必引起不對稱的水分場分布。為此取11月5日（初始時刻）和1月20日（監(jiān)測時刻）陽坡、渠底及陰坡的測點A、B、D、F、G 的總含水率沿基土深度分布結(jié)果，如圖12所示。

圖12 關(guān)鍵點渠基土體積含水率分布曲線

由圖可知，相比于初始時刻，凍結(jié)區(qū)域內(nèi)因外界負積溫作用而使未凍區(qū)內(nèi)的水分逐漸向上遷移，水分集聚現(xiàn)象明顯。由數(shù)值結(jié)果可知，因太陽輻射和晝夜溫度的疊加作用，凍土內(nèi)存在多個水分集聚位置，而非只發(fā)生在凍結(jié)鋒面處；陽坡、渠底及陰坡水分的集聚個數(shù)、位置及深度均不同，呈現(xiàn)出陰坡含水量最多而陽坡最少的分布規(guī)律。對比監(jiān)測結(jié)果可知，因現(xiàn)場選取局部點來烘干，僅能反映出水分運移規(guī)律，并不能嚴格反映出含水量最高點或凍結(jié)鋒面位置，但在凍深范圍內(nèi)呈現(xiàn)陰坡、渠底總含水量高而陽坡較少的分布規(guī)律，與數(shù)值結(jié)果一致。然數(shù)值解的表層含水量偏高，80～100 cm 深度處偏低，可能與未考慮陰天、蒸發(fā)等因素及現(xiàn)場局部地質(zhì)引起的初始含水偏高所致；但其他位置點的含水量吻合度較高，滿足工程要求，即此模型可模擬太陽輻射下的水分遷移、凍結(jié)相變等特征。

為進一步分析不對稱的冰、水含量分布，取最大凍深時刻（3月16日）和最大融深時刻（3月28日）的水分場結(jié)果，如圖13—14所示。

由圖13 可知，凍結(jié)區(qū)域內(nèi)的未凍水含量隨溫度降低而逐漸減少，最低為0.06，接近于基土的殘余含水量；而相應(yīng)位置的結(jié)冰量逐漸增加，冰含量呈條帶狀分布，無明顯的空間分布規(guī)律［2-3，12］，但與基土凍深發(fā)展相呼應(yīng)（圖10）；因凍結(jié)初期基土表層常處于凍融循環(huán)狀態(tài)，在多次溫度抽吸力下表層冰含量較大；冰、水總含量最大為0.6，已遠大于基土的孔隙率，將發(fā)生嚴重的凍脹變形。

3月16日表層基土出現(xiàn)局部融化現(xiàn)象，含水量增加而冰含量減少。隨著外界溫度的逐漸升高，凍結(jié)層內(nèi)冰透鏡體開始從兩端發(fā)生雙向融化，局部區(qū)域形成過量孔隙水，隨著時間推移，孔隙水逐漸擴散，形成隨深度增加含水量逐漸增大的分布規(guī)律（圖14（a）），因凍結(jié)區(qū)域內(nèi)冰含量的阻水作用，使含水量較初始時刻明顯增大。對比圖13（b）和圖14（b）可知，凍深范圍外冰透鏡體全部融化，且凍深線附近因溫度升高至其融化溫度，使冰含量值明顯減少，而凍深范圍內(nèi)的其余區(qū)域因熱流傳輸密度少和凍深線附近的冰融化吸熱，從而使冰含量未發(fā)生較大變化。

圖13 3月16日渠道水分場分布

圖14 3月28日渠道水分場分布

由圖13、圖14 亦可發(fā)現(xiàn)，因太陽輻射和陰坡遮蔽而產(chǎn)生的熱邊界差異產(chǎn)生了明顯的不對稱水分場分布，陰坡總含水量明顯大于陽坡。首先，陰陽坡的凍深差異特征形成了陰坡冰透鏡體總體厚度較陽坡大的現(xiàn)象；其次，陰陽坡的熱邊界差異使陰坡含冰量大于陽坡含冰量；最后，因陽面接受太陽輻射多，溫升快，從而使陽坡融化速度明顯大于陰坡，冰含量快速減少。

4.3 變形場分析 由上節(jié)可知，陰、陽坡局部位置總含水量大于孔隙率，將會產(chǎn)生差異性凍脹變形，選取多個時刻下襯砌板的法向凍脹變形沿渠周分布曲線，如圖15所示。

圖15 襯砌板法向凍脹變形沿渠周分布曲線

由圖可知，梯形渠道襯砌因基土凍脹作用，呈現(xiàn)出渠底向上隆起，坡腳受擠壓明顯，兩坡板向渠內(nèi)凸起，襯砌整體上抬的變形趨勢。其中，渠底襯砌板在偏陽坡處的法向凍脹量最大，陰、陽坡板約在1/4～1/3 坡長位置處最大。在因太陽輻射而產(chǎn)生的差異凍脹力作用下，襯砌凍脹變形量由陰坡、渠底至陽坡依次遞減，陰、陽坡板凍脹變形差異明顯，整體與基土內(nèi)凍深發(fā)展、水分分布規(guī)律一致，即渠道水熱分布的不對稱、不同步特征導(dǎo)致了陰、陽坡板不對稱的變形分布。

由上圖亦可知，隨著凍結(jié)時間的發(fā)展，襯砌板法向凍脹變形均逐漸增大，至3月16日時，凍脹量最大，陽坡、渠底和陰坡最大分別為8.47、10.97 和12.48 cm，與現(xiàn)場監(jiān)測值基本吻合，且襯砌板的整體變形趨勢與監(jiān)測值分布規(guī)律基本一致，同時，在后續(xù)融化過程中，襯砌板凍脹變形亦隨基土內(nèi)冰水總含量的減少而逐漸變小，表明本文的數(shù)值模型可較好的模擬襯砌結(jié)構(gòu)在基土凍脹作用下的凍融變形特征。

5 考慮太陽輻射渠道凍脹模型的熱邊界選取及合理性評價

太陽輻射是渠道凍脹分析中能量的來源，影響著基土內(nèi)部溫度場、水分場及變形場的耦合作用，是進行寒區(qū)渠道設(shè)計的重要因素，其中熱邊界合理選取至關(guān)重要?，F(xiàn)階段，凍脹數(shù)值模擬中常用以下兩種熱邊界：（1）以日或月為單位擬合環(huán)境溫度，只采用熱對流邊界［2-3，9，13］，但無法考慮太陽輻射產(chǎn)生的陰陽坡效應(yīng)；（2）結(jié)合現(xiàn)場實測數(shù)據(jù)，采用正弦函數(shù)擬合陰、陽坡溫度變化曲線，并將其直接賦予結(jié)構(gòu)邊界［17-18］，但卻無法充分反映出太陽輻射的空間效應(yīng)和晝夜溫差的時間效應(yīng)。目前公認更為準確的是第二種以現(xiàn)場地溫監(jiān)測為依據(jù)的擬合溫度邊界，以此為參照，對比分析本文采用的太陽輻射、面-面輻射、環(huán)境輻射的熱輻射邊界和晝夜溫差的熱對流邊界組合，對熱邊界選取及其合理性進行評價。擬合溫度邊界取上章仿真結(jié)果的陰坡、渠底和陽坡的日均表面溫度數(shù)據(jù)（圖9）來正弦擬合熱邊界，可彌補現(xiàn)場監(jiān)測點少且數(shù)據(jù)不全的缺點。下面將從溫度、水分、變形來分別闡述。

取陰、陽坡2/3 設(shè)計水位點的凍深發(fā)展過程結(jié)果，如圖16所示。由圖可知，擬合的溫度邊界使襯砌板溫度緩慢變化，凍結(jié)速率發(fā)展較為均勻；而本文模型因考慮了太陽輻射空間效應(yīng)和晝夜溫度時間效應(yīng)的耦合疊加作用，使得凍深發(fā)展在局部位置出現(xiàn)了短暫的平臺狀，凍結(jié)速率忽快忽慢，但與外界的環(huán)境變化相協(xié)調(diào)，更符合現(xiàn)場實際規(guī)律。

圖16 不同熱邊界下的凍深發(fā)展過程線

陰、陽坡含水量分布規(guī)律基本類似，取最大凍深時刻下陰坡2/3 設(shè)計水位點的總含水量沿基土深度的分布曲線，如圖17所示。由圖可知，因擬合溫度邊界下凍結(jié)速率變化較為均勻，除淺層基土外，含水率由上至下逐漸增大，至凍結(jié)峰面處最大，與單向凍結(jié)結(jié)果基本一致。而本文模型因考慮太陽輻射和晝夜溫差作用，凍結(jié)初期的“夜凍晝消”往復(fù)凍融循環(huán)作用使得基土表層含水增加；水分分布與凍土的凍結(jié)速率相對應(yīng)，凍結(jié)速率慢時水分遷移多，反之亦然，更能反映出外界環(huán)境作用引起的水熱耦合作用。

圖17 不同熱邊界下的水分場分布

總含水率和凍深差異導(dǎo)致了不同的凍脹變形，取陰、陽坡2/3 設(shè)計水位點，其襯砌板法向凍脹變形隨時間的分布規(guī)律如圖18所示。由圖可知，襯砌板的起始凍脹時刻、凍脹量發(fā)展速度和量級，以及融化速度等變量均與熱邊界有關(guān)，擬合溫度邊界方法得到的凍脹量不能反映出太陽輻射和晝夜溫差作用所引起的凍結(jié)速率差異而導(dǎo)致的水熱分布規(guī)律，導(dǎo)致結(jié)構(gòu)凍融破壞計算不準確。

圖18 不同工況下襯砌板法向凍脹量隨時間變化曲線

結(jié)合溫度、水分和變形結(jié)果可知，本文數(shù)值模型的結(jié)果與目前公認的以地溫監(jiān)測數(shù)據(jù)為基礎(chǔ)的擬合溫度邊界所得規(guī)律基本一致，結(jié)合上一節(jié)結(jié)果，可進一步說明本文模型的準確性；同時本文模型綜合考慮了太陽輻射的空間效應(yīng)和晝夜溫度的時間效應(yīng)，可反映出外界環(huán)境作用下寒區(qū)渠道水、熱、力耦合的凍脹破壞本質(zhì)，特別是能減少對現(xiàn)場監(jiān)測值的依賴性，且具有理論通用性，相比于目前使用較多的擬合溫度邊界，結(jié)果更加合理、準確，且可輔助設(shè)計單位進行前期渠道設(shè)計，滿足精細化設(shè)計、復(fù)核的要求。

6 結(jié)論

考慮太陽輻射和晝夜溫差的時空變化，提出了渠道的陰影計算方法，基于凍土的橫觀各向同性凍脹特征，建立了“太陽輻射-環(huán)境溫度-凍土-工程”相互作用的水-熱-力耦合凍脹模型，并探討了不同熱邊界對渠道水-熱-力耦合計算的影響，結(jié)合工程監(jiān)測及數(shù)值模擬結(jié)果，得到如下結(jié)論：（1）考慮太陽輻射的渠道凍脹模型計算的溫度、水分和變形結(jié)果與現(xiàn)場實測值較為吻合，分布規(guī)律一致，表明本文模型可分析太陽輻射空間效應(yīng)和晝夜溫差時間效應(yīng)下的寒區(qū)渠道凍結(jié)過程，理論通用性強，且可減少對現(xiàn)場監(jiān)測的依賴，可為寒區(qū)渠道的凍脹模擬、預(yù)測及災(zāi)害防控設(shè)計提供理論依據(jù)。（2）探明了渠道表面太陽輻射的時空分布規(guī)律，陰、陽坡太陽輻射時長、輻射量存在顯著差異，其表面溫度由外界氣溫和太陽輻射量共同決定，計算的陰、陽坡日均溫差最大為3.5 ℃，導(dǎo)致了陰坡早于陽坡15 d 凍結(jié)，晚于陽坡1 d 融化，凍深發(fā)展差異明顯，最大相差31 cm，渠道陰陽坡溫度場分布不對稱且不同步。（3）凍結(jié)期渠道陰坡結(jié)冰區(qū)域及量值明顯大于陽坡，“夜凍晝消”的凍融循環(huán)作用使基土淺層含水量明顯增多，而融化期陽坡冰融化速度大于陰坡，水分場呈不對稱、不同步、層狀分布等特點。（4）渠坡板和底板的最大凍脹變形分別發(fā)生在1/4～1/3 坡長位置和偏陽坡位置，陰、陽坡變形隨時間的發(fā)展過程呈不均勻、不同步、不對稱性，計算的最大變形分別為12.48 和8.47 cm。（5）渠道陰陽坡效應(yīng)產(chǎn)生的本質(zhì)在于太陽輻射和晝夜溫差引起的渠道能量收支差異、溫度梯度和凍結(jié)速率變化所致，即熱邊界引起的渠道基土的水-熱-力動態(tài)耦合作用。（6）前期設(shè)計時，宜采用多年氣溫擬合值和太陽輻射參數(shù)確定熱邊界；已建工程凍脹評價時，宜采用日氣溫值和太陽輻射參數(shù)確定熱邊界，并輔以監(jiān)測值對其修正。