趙立財, 殷珂, 王瑞*, 鄭蕓昕, 任翔
(1.中鐵十九局集團第三工程有限公司, 沈陽 110136; 2.長安大學建筑工程學院, 西安 710061)
季節性凍土約占中國國土總面積的53.5%,隨著中國基礎交通建設的高速發展,季凍區隧道數量與日俱增。隧道凍害頻發造成了大量經濟損失與資源浪費,成為困擾寒區隧道工程的一大難題。隧道凍害與溫度場密切相關,探討氣象條件對溫度場的影響是揭示隧道凍害機理的有效途徑,也是優化凍害防控措施的基礎[1]。
Jun等[2]和Li等[3]實測發現隧道不同進深溫度隨時間變化規律類似,不同時刻隧道出入口附近溫差較大,隨進深增加溫差漸趨減小。冬季寒冷氣象條件下,不僅隧道兩端圍巖會發生負溫凍結,洞內也可能出現較廣泛的負溫分布[4-5]。惡劣氣象條件下鋪設保溫層后隧道圍巖仍有較廣泛的負溫區,傳統被動防寒措施無法完全消除隧道凍害[6-7]。
數值計算方法兼具高效率、高精度的優點,目前已廣泛應用于寒區隧道溫度場的研究中。王仁遠等[8]、Tao等[9]和Tan等[10]通過數值方法分析了自然狀態下的隧道溫度場,發現風溫和風速對隧道入口附近溫度分布有較大影響。Yu等[11]建立三維瞬態傳熱模型分析了冬、夏典型氣象條件下溫度場的分布及演變規律。Wang等[12]建立隨機有限元模型,分析了隨機環境及巖性條件下的隧道溫度場。楊波等[13]和高焱等[14]對比分析了自然狀態和列車開行狀態下的溫度場,認為自然風是隧道溫度場分布狀態的主要影響因素。王開運等[15]分析了圍巖及氣象條件對隧道溫度場的影響強度,發現氣象條件的影響更為顯著。鄭余朝等[16]正交分析了隧道溫度對氣象因素的敏感性,并討論了風向夾角交替變化情況下的隧道溫度場。
氣象條件會顯著影響隧道溫度場的分布及演化規律。采用數值方法分析氣象因素對隧道溫度場的影響時,過度抽象的計算參數無疑會導致模擬結果與實際情況產生較大偏差。此外,持時及風向作為重要的氣象參數,其對隧道溫度場的影響規律研究尚不多見。有鑒于此,現依托季凍區某隧道工程,在統計分析工程區近41年氣象數據的基礎上凝練氣象特征參數,采用三維有限元手段,通過正交試驗分析風溫、風速、持時及風向4個氣象因素對季凍區超長隧道溫度場的影響規律。
新建沈陽-佳木斯鐵路線的新賓隧道工程位于遼寧省新賓縣。工程區屬溫帶大陸性季風氣候,平均海拔492 m。新賓隧道為季凍區隧道,起訖里程DK123+425—DK133+600,全長約10 175 m,走向約為北偏東70°。
根據新賓氣象觀測站(54353)近50年的觀測數據顯示,新賓地區全年最冷月份為1月,最冷月平均溫-15 ℃,工程區環境較惡劣,隧道易在冬季發生凍害。圖1為新賓地區1980—2020年日平均溫變化。

圖1 1980—2020年日平均溫變化規律Fig.1 Variation law of daily average temperature from 1980 to 2020
新賓縣1980—2020年最冷月份風向玫瑰圖如圖2所示。為便于討論不同氣象條件下的隧道溫度場,將16方位風向簡化為圖2中的4個風向,可以看出風向1為新賓地區最冷月份的主導風向。

圖2 1980—2020年最冷月風向玫瑰圖Fig.2 Wind roses in the coldest months from 1980 to 2020
當前以數值方法分析隧道溫度場時,氣象參數的選取常忽略其隨機性的特點,導致數值計算與工程實際結合不夠充分。本文研究基于統計思想提取新賓地區1980—2020年最冷月平均溫度及主導風向風速、持續時間的頻率分布如圖3所示。

圖3 1980—2020年最冷月氣象參數頻率分布Fig.3 Frequency distribution of meteorological parameters in the coldest month from 1980 to 2020
從圖3可以看出工程區氣象參數頻率分布有明顯隨機性,經K-S檢驗發現最冷月平均溫度服從正態分布,即x~N(-15.06,2.332);主導風向風速及持時均基本服從對數正態分布,即風速lnx~N(0.109,0.599 2),持時lnx~N(0.828,0.741 2),檢驗結果如表1所示。

表1 檢驗結果Table 1 Test results
基于新賓地區氣象數據的分布特征,設置隧址區最冷月基本氣象參數。風溫服從正態分布,取其期望值-15 ℃為基本氣象參數;主導風向的風速及持時服從對數正態分布,取右側保證率為50%的截斷值為基本氣象參數(圖3),同理確定風向2、3、4的風速及持時參數。新賓地區最冷月基本氣象參數如表2所示。

表2 基本氣象參數Table 2 Basic meteorological parameters
隧道空氣流動遵循質量、動量及能量守恒定律,可通過控制方程對其進行數學描述,分別稱為連續性方程、動量方程及能量方程,其通用表達形式[17]為

(1)
式(1)中:φ為u、v、w、T等通用求解變量;ρ為空氣密度;V為速度矢量;t為時間變量;Г為廣義擴散系數;S為廣義源項;div為向量場的散度;grad為向量梯度。
本文流場選用RANS方法的SST(shear stress transport)k-ω兩方程湍流模型描述,其在近壁面處保留原始k-ε模型,主流區則應用k-ω模型,在保證近壁面區域計算精度的前提下縮減了計算成本,k方程和ω方程為
(2)
(3)
式中:Gk為平均速度梯度湍流動能;Gω為ω方程的值;Гk、Гω為有效擴散項;Yk、Yω為發散項;Sk和Sω為源項;uj為各方向瞬時速度分量;xj為各方向空間坐標。各參數解釋及取值方法詳見參考文獻[18]。
圍巖固體傳熱控制方程為
(4)
式(4)中:Ceq為等效體積熱容;Ts為溫度;λeq為等效導熱系數;Qe為熱源產生的熱量。
分析隧道溫度場時通過構建比熱容C和導熱系數λ隨溫度變化函數考慮含水圍巖水冰相變,假設水冰相變發生在溫度Tm附近一定范圍內(Tm±ΔT),其相關表達式為
(5)
(6)
式中:Cf為凍結圍巖比熱容;Cu為未凍圍巖比熱容;λf為凍結圍巖導熱系數;λu為未凍圍巖導熱系數;L為相變潛熱;T為水(冰)溫度;ΔT為水冰相變溫度區間半徑;Tm為水冰相變溫度區間的中溫度,取為0 ℃。
FLUENT數值計算軟件在流-固耦合傳熱及非線性溫度場分析方面功能強大,采用該軟件建立“空氣-襯砌-圍巖”三維耦合傳熱模型,基于非穩態計算結果討論隧道溫度場分布特征及其對不同氣象因素的敏感度。
模型邊界尺寸取隧道等效直徑3~5倍可保證溫度場計算精度[19],考慮新賓隧道實際情況取模型尺寸為70 m×70 m×10 000 m,模型由二階實體單元建立并在流體區劃分邊界層。模型下邊界設為熱流邊界,熱流密度為0.06 W/m2[20];暴露于隧道外的襯砌面設為熱對流邊界,對流換熱系數取15 W/(m2·℃)[21];流體入口為速度入口,選用Components指令直接指定各方向速度分量;流體出口為壓力出口,工作條件為標準大氣壓;其余邊界為絕熱邊界。
參考文獻[22]以地溫梯度法確定圍巖初始溫度,以當地年平均溫度5.3 ℃為恒溫層溫度[23],取恒溫層深度60 m;地溫梯度為2 ℃/100 m[20],隧道平均埋深約150 m,可設置圍巖初始溫度為7 ℃。模型網格劃分及邊界設置如圖4所示,材料相關參數如表3所示。

圖4 有限元模型Fig.4 Finite element model
通風條件對隧道溫度場分布有顯著影響,且工程區風向會經常性改變。為分析不同通風狀態下的隧道溫度場,以主導風向為基礎結合后續風向設置3種風向組合:風向1+風向2(組合1);風向1+風向3(組合2);風向1+風向4(組合3)。
參考新賓地區最冷月基本氣象參數設置典型工況氣象條件,取主導風向持時為60 h,后續風向持時為30 h,通風90 h為一個周期;各風向風速按表2取值,假設各風向中心線為該風向實際來流方向,風溫為-15 ℃。
典型工況風向組合1通風60 d后入口進深300 m范圍內隧道溫度云圖如圖5所示??梢钥闯?較小進深范圍內隧道南北兩側溫度分布有明顯差異,南側溫度低于北側且低溫區域分布更廣,隨進深增加差異漸趨減小,溫度有整體升高趨勢。
工程中常以隧道二襯底面是否出現負溫作為抗凍設防的判據[14],同時考慮到隧道南北兩側溫度分布存在差異,分別在隧道南北兩側拱墻2.2 m高處初支與二襯之間設置如圖5所示的溫度測線。典型工況風向組合1通風60 d后沿測線溫度及南北溫差分布如圖6所示。

圖6 典型工況風向組合1通風60 d后溫度分布Fig.6 Temperature distribution after 60 days of ventilation in wind direction combination 1 under typical working conditions
從圖6可以看出,隧道南側溫度以隧道兩端為起點隨進深逐漸升高;北側溫度隨進深增加呈升高-降低-升高趨勢,兩端進深10~40 m范圍內出現了溫度下降段。原因是風向組合1條件下隧道南側為主要受風側,出入口端部溫度更低;隨進深增加(10~40 m)冷空氣逐漸擴散,北側溫度隨之降低;隨進深進一步增加冷空氣逐漸擴散均勻,南北側溫度均持續升高。隧道南北溫度場存在差異可能會導致凍脹力分布不均勻,對隧道的抗凍設防造成影響。隧道南北溫差以隧道兩端為起點隨進深先增大后減小,在進深8 m附近溫差達到峰值。
工程區自然環境狀態一定程度上決定了隧道是否發生凍害及凍害的強度,開展隧道溫度場相關討論時宜考慮不同氣象條件。典型工況3種風向組合條件下通風60 d后隧道出入口負溫長度及南北溫差最大處的溫度如圖7所示。

圖7 典型工況隧道入口端溫度場特征值Fig.7 Characteristic values of temperature field at tunnel inlet under typical working conditions
由圖7(a)可知,不同風向組合條件下隧道入口南北兩側負溫長度基本相等,可將其均值定義為入口負溫長度,風向組合3隧道入口負溫長度最大為
2 680 m。各風向組合隧道出口負溫長度均較小,風向組合1南側負溫長度為240 m,其余風向組合南北兩側負溫長度均不超過25 m。隧道入口為抗凍設防重點,出口凍害相對較弱。
由圖7(b)可知,不同風向組合條件下隧道入口南側溫度均低于北側,風向組合1條件下南北兩側溫度最高且溫差最大為7.2 ℃,風向組合3條件下溫度最低且溫差最小為2.5 ℃。隧道出口在風向組合1條件下溫差最大為5.7 ℃,且各風向組合條件下隧道出口溫差均小于隧道入口。
正交試驗設計是多因素多水平試驗中廣泛應用的一種設計方法,依據數理統計學與正交性原理,從全因素水平中選取具有代表性的點開展局部試驗,試驗具備均勻分散、整齊可比的特點。通過正交試驗設計,可在實現全因素分析的同時大量減少試驗次數,降低試驗成本。
采用基于正交試驗的方法分析不同氣象因素對隧道溫度場影響的敏感性,將正交試驗設計與數值模擬相結合,每模擬一個工況的溫度場相當于一次試驗。溫度場模擬仍采用主導風向結合后續風向的組合通風法,每個工況設置3種風向組合。影響隧道溫度場的氣象因素主要有風溫、主導風向風速及持時、后續風向風速及持時,試驗共設置此5個因素,每個因素在典型工況基礎上波動-25%、0、25%和50%的4個水平。正交試驗選用L16(45)正交表僅需開展16次試驗,而全面分析法則需45即1 024次試驗,采用正交試驗法可極大提高試驗經濟性。
由于典型工況條件下不同風向組合隧道出口北側負溫長度發展均不明顯,以隧道入口負溫長度及南北最大溫差、隧道出口南側負溫長度及南北最大溫差為溫度場評價指標。取不同風向相組合條件下的最大值為試驗結果,分別討論隧道出入口附近溫度場對氣象因素的敏感性。
選用L16(45)正交表開展試驗,將影響因素各水平試驗結果的最大差值作為極差,利用極差分析法分析各因素對隧道溫度場的影響強度,隧道通風60 d后各影響因素極差如圖8所示。
從圖8可以看出,風溫對入口負溫長度影響最大,主導風向風速的影響稍小于風溫,后續風向風速影響相對較小,其余因素影響可忽略。對于出口南側負溫長度,風溫同樣是影響最大的因素,主導風向持時、后續風向風速及持時的影響程度較接近,主導風向風速影響最小。
風溫對入口溫差影響最大,主導風向風速的影響弱于風溫,其余因素影響可忽略。對于出口溫差,風溫的影響最大,后續風向持時的影響次之,后續風向風速影響最小。
分析發現風溫對隧道溫度場影響最大,其余因素影響存在差異,討論季凍區隧道溫度場時宜差異化考慮不同氣象因素。溫度場各評價指標隨因素水平變化如圖9所示。

圖9 評價指標隨因素水平變化規律Fig.9 Variation law of evaluation index with factor level
從圖9(a)可以看出,隧道入口負溫長度隨風溫降低、主導及后續風向風速的增加逐漸增大,變化速率從大到小依次為風溫、主導風向風速、后續風向風速。出口南側負溫長度隨風溫降低、后續風向風速及持時的增加逐漸增大,隨主導風向風速的增加逐漸減小。
從圖9(b)可以看出,隧道入口溫差隨風溫降低逐漸增大,隨主導風向風速增加逐漸減小。隧道出口溫差隨風溫降低、后續風向持時的增加逐漸增大,隨主導風向風速及持時的增加逐漸減小;后續風向風速增加出口溫差變化無明顯規律,原因是其他因素干擾導致溫差波動,同時說明該因素對出口溫差影響較小。
依托新賓隧道工程,根據工程區近41年氣象數據總結典型氣象條件參數,通過數值方法分析了不同通風狀態下隧道溫度場分布的基本特征。在此基礎通過正交試驗分析了氣象因素對隧道溫度場影響的敏感性,得出以下結論。
(1)工程區風溫服從正態分布,風速及持續時間服從對數正態分布。由于風向與隧道縱向的夾角使隧道南北兩側存在較大溫差,溫差隨進深先增大后減小;典型氣象條件下,隧道入口為抗凍設防重點。
(2)隧道入口負溫長度及溫差主要受風溫和主導風向風速影響;出口南側負溫長度及溫差受風溫、主導風向及后續風向持時的影響較大。
(3)隧道入口負溫長度隨風溫降低、主導及后續風向風速的增加逐漸增大,溫差隨風溫降低逐漸增大;出口南側負溫長度隨風溫降低、后續風向及風速的增加逐漸增大,溫差隨風溫降低、后續風向持時的增加逐漸增大。