大流速地下水作用下多圈凍結孔優(yōu)化布置方法研究

董艷賓，榮傳新，王 彬，楊 凡

（安徽理工大學 土木建筑學院，安徽 淮南 232001）

人工地層凍結法目前是富水軟土層豎井掘砌施工的主要工法，該工法在擬施工的豎井井筒周圍設置1 圈或多圈凍結管，通過凍結管中的低溫冷媒與被凍土體之間不斷進行熱量交換，形成1 道具有一定強度、且具有良好封水性能的凍結壁，從而為井壁的掘砌施工提供1 個穩(wěn)定的施工環(huán)境[1-2]。但隨著施工穿越巖土層中地下水流速的增大，采用傳統(tǒng)的凍結管布置方案時，出現(xiàn)了凍結壁交圈時間增加甚至無法交圈的問題，這導致后期的土體開挖以及井筒的施工無法如期開展，造成了巨大的經濟損失。因此，對大流速地下水作用下的人工立井凍結溫度場的形成規(guī)律展開研究具有重要的工程意義[3-5]。

在試驗研究方面，周曉敏等[6]通過雙管凍結正交模型試驗，研究了常規(guī)鹽水（溫度為-26.0～-30.6℃）凍結工藝中，地下水滲流、孔間距等對飽和砂凍結交圈時間和上下游溫度場發(fā)展規(guī)律的影響。Huang R C[7]進行了滲流條件下的單一凍結管的模型試驗，發(fā)現(xiàn)地下水流會導致凍結管上、下游的凍結面積與無滲流條件相比出現(xiàn)了不同程度的減小。Li Yan Lao[8]通過模型試驗得出了不同流速的地下水作用下，單排三管凍結溫度場的發(fā)展規(guī)律。王朝暉等[9]對液氮（-80 ℃）凍結條件下地下水流速對凍結溫度場的影響開展了研究，其結果表明10 m/d 以上的動水對液氮凍結效果有顯著的影響。Vitel M[10]等設計了1 種與熱力學原理完全一致的水熱耦合數(shù)學模型，用于模擬滲流條件下飽和不可變形的多孔介質的人工地層凍結，完成了高滲流速度條件下三維凍結試驗，為多孔介質相變相關的水熱耦合研究提供了關鍵參數(shù)依據(jù)。Pimentel E 等[11]總結了之前學者的研究成果，并考慮到凍結冷量散失等問題，重新設計了試驗裝置，分別進行了流速為0、1、1.5、2.0、2.1 m/d 的滲流條件下人工地面凍結的大型模型試驗，根據(jù)試驗結果從凍結柱狀體的交圈時間與所需冷量的角度出發(fā)，對凍結壁的交圈形式解進行了討論，為凍結法的數(shù)值計算以及工程應用提供了重要依據(jù)。

由于大型相似模型試驗的限制因素較多，且只能部分還原凍結溫度場的發(fā)展規(guī)律，而數(shù)值計算則可以解決這一問題，因此數(shù)值計算一直是研究水熱耦合問題的重要技術手段。在數(shù)值計算研究方面，楊平等[12]利用多孔介質熱運移理論及達西定律，建立了考慮地下水流作用的單根凍結管凍結峰面發(fā)展的數(shù)學模型，分析了凍結過程中溫度場及地下水流場的變化規(guī)律。高娟[13]以及劉建剛等[14]分別運用有限元的方法對地下水作用下豎井和水平凍結的溫度場形成規(guī)律進行了研究。Vitel M[10,15-16]為了模擬在滲流條件下飽和不可變形多孔介質的人工地面凍結過程，構建了與熱力學一致的水熱數(shù)值模型，并對通常使用的限制性假設進行了簡化，該模型在高滲流速度條件下的三維地面凍結實驗中得到了很好的驗證。黃詩冰等[17-19]通過考慮水/冰相變，開發(fā)了水熱耦合模型來模擬水流對凍結過程的影響，并將該模型與基于COMSOL 多物理場平臺的Nelder-Mead 單純形法相結合，優(yōu)化了圓形隧道周圍凍結管的位置。Ahmed 等[20]結合“蟻群算法”對小流速地下水作用下凍結管的布置位置進行了優(yōu)化設計，從而加快了凍結壁的交圈并使形成的凍結壁的強度更加均勻，為凍結方案的優(yōu)化設計開拓新的思路。

將在上述研究的基礎上，通過對比單管水熱耦合模型與試驗結果，驗證水熱耦合模型的正確性和相關參數(shù)取值的合理性。根據(jù)驗證的水熱耦合模型，采用對淮南礦區(qū)潘一煤礦副井凍結法施工的上游外圈凍結管進行局部加密的方法，來實現(xiàn)對常規(guī)凍結方案的優(yōu)化設計，并通過COMSOL 有限元軟件對優(yōu)化效果進行研究。

1 滲流場和溫度場的耦合模型的建立

1.1 基本假設

在凍結過程中，土體是由土體骨架、水、空氣以及冰構成的多相體系。假設凍結過程中的土體保持為1 種飽和狀態(tài)。考慮研究的重點及主要影響因素，在計算模型過程中做出如下假定：

1）忽略凍結管周邊的熱量損失，即視凍結管內、外壁的溫度相同，且都等同于低溫冷媒溫度。

2）流體作為連續(xù)介質充滿滲流區(qū)的全部空間，初始狀態(tài)為穩(wěn)態(tài)滲流且符合達西定律，不考慮凍結過程中的質量遷移。

3）凍結土層為連續(xù)的、均勻的、各向同性的多孔介質，滿足混合物理論的基本假定。

4）在凍結過程中僅考慮土體基質、水、冰的熱傳導效應及冰、水在相變溫度區(qū)域范圍內發(fā)生，且未凍水含量僅與溫度有關。

1.2 水熱耦合控制方程

凍結過程的溫度變化屬于瞬態(tài)熱傳導問題，根據(jù)基本假定，則在凍結過程中，土體中同時存在冰、水和土體骨架，且三者存在不同的傳熱特征：除了傳導傳熱以外，冰的熱量傳遞還包括相變傳熱；水的熱量傳遞還包括液態(tài)水流動產生的對流傳熱。根據(jù)傳熱及滲流原理，凍結過程考慮相變的溫度場與滲流場控制方程為[5,15]：

式中：t 為時間；T 為控制體溫度；χ 為空隙中未凍水體積含量；ci、cl、cs為冰、水、土骨架的比熱容；λi、λl、λs分別為冰、水、土骨架的導熱系數(shù)；ρi、ρl、ρs為冰、水、土骨架的密度；Lw為單位質量水的相變潛熱值為水的相對速度矢量；Q 為等效熱源強度；αl為水的壓縮率；αs為控制體等效壓縮率；φ 為孔隙率；p 為滲透壓力；q 為流量；K′為等效滲透系數(shù)；Kf、K0為完全凍結區(qū)和未凍區(qū)的滲透系數(shù)；θi、θl、θs為控制體的含冰率、含水率、土骨架含量；g 為重力加速度。

2 水熱耦合模型的驗證

2.1 物理模型試驗的設計

為驗證水熱耦合數(shù)學模型的合理性，采用相似模型試驗進行驗證。為了減小因模型尺寸縮小而帶來的失真影響，并考慮試驗室條件的基礎上，取幾何相似比為C1=3，根據(jù)相似準則，測點布置如圖1。

圖1 測點布置方案Fig.1 Measuring point arrangement

試驗選用粒徑為1 mm 的圓粒砂來模擬多孔介質，砂的物理參數(shù)和熱參數(shù)見表1 和表2。

表1 砂子物理參數(shù)Table 1 Physical parameters of the sand

表2 砂子熱參數(shù)Table 2 Thermal parameters of sand

為了防止外界的環(huán)境溫度對凍結過程產生干擾，對箱體外表面以及凍結管路進行了保溫處理：緊貼箱體外表面布置1 層30 mm 的橡塑保溫板，隨后采用40 mm 的聚氨酯保溫板將整個箱體包裹在其中，接縫處用氯丁膠緊密粘合，外露的凍結管以及凍結干管表面包裹1 層30 mm 的橡塑保溫層。保溫層的布置情況如圖2。

圖2 箱體表面保溫層設置情況Fig.2 Insulation layer setting of case surface

2.2 試驗方法及測試方案

1）進行基礎試驗，測定選用的砂層的導熱系數(shù)、質量比熱、孔隙率、滲透率。

2）試驗設備安裝與調試，標定變頻泵、壓力表、流量計及測溫系統(tǒng)。

3）設置變頻泵的恒壓數(shù)值改變水箱內水壓，調節(jié)進水流量控制閥改變清水流量，控制通過砂層的清水流量，使砂層中滲流水速度達到設計值，由PLC控制變頻泵，保證試驗過程中上下游壓力差恒定。

4）調節(jié)酒精溫度至-32 ℃，開啟酒精泵，通過調節(jié)酒精流量來使酒精供冷量達到設計值，記錄開泵時間。

5）開啟記錄數(shù)據(jù)，間隔5 min 記錄箱體邊界溫度、砂層溫度、進出清水溫度、去回路酒精溫度和環(huán)境溫度。

2.3 數(shù)值計算結果與模型試驗結果對比

依據(jù)COMSOL Multi-physics 中的場方程，利用Porous Media and Subsurface Flow 模 塊 和 Hear Transfer in Porous Media 模塊的瞬態(tài)分析功能，將各個參變量代入方程，并調用借助Matlab 編制的內插函數(shù)與耦合方程，實現(xiàn)了凍結過程中水熱耦合數(shù)學模型的瞬態(tài)求解。將1#～8#測點在流速分別為0 m/d和2 m/d 時的數(shù)值計算結果與模型試驗結果進行對比，兩者溫度隨時間和位置變化對比結果如圖3，兩者的量化比較見表3。

圖3 溫度變化規(guī)律對比分析Fig.3 Comparative analysis of temperature change law

表3 數(shù)值計算結果與試驗結果比較Table 3 Comparison between numerical calculation results and test results

1）由圖3（a）和圖3（b）可知，數(shù)值計算結果與模型試驗結果基本一致。流速為0 m/d 時，1#～8#測點均經歷了快速降溫階段和穩(wěn)定降溫階段；流速為2 m/d 時，上游距離凍結管較遠的1#、2#測點僅有穩(wěn)定降溫階段，而較近的 3#、4#測點和下游 5#、6#、7#、8#測點則都經歷有快速降溫階段和穩(wěn)定降溫階段。

2）由圖3（c）和圖3（d）可知，在流速為 0 m/d 凍結1 500 min 時，對稱點測點1#與8#的溫差為0.2℃、2#與 7#的溫差為 1.0 ℃、3#與 6#的溫差為 0.1℃、4#與5#的溫差為0.7 ℃，對稱測點的溫差絕對值均不大于1.0 ℃，可知測點的溫度變化具有位置對稱性；流速為2 m/d 凍結1 500 min 時，對稱點測點1#與 8#的溫差為 15.1 ℃、2#與 7#的溫差為 10.1 ℃、3#與 6#的溫差為 5.1 ℃、4#與 5#的溫差為 3.4 ℃。這是因為在地下水流的作用下，凍結管上游的一部分冷量會被水流帶到下游，使得上游測點降溫慢，下游測點降溫快；由于水流會抑制凍結管冷量向上游傳遞、并促進向下游傳遞，隨著距離凍結管距離增加，這種抑制和促進作用越明顯。因此關于凍結管對稱的下游測點溫度低于上游測點，且測點距離凍結管越近，兩點的溫度差值越小。

3）由表3 可知，1#～8#測點在凍結 1 000、1 500 min 時，地下水流速分別為 0、2 m/d 時，0 m/d 的數(shù)值計算的平均誤差分別為 0 ℃和-0.4 ℃、2 m/d 的數(shù)值計算的平均誤差為 0.7 ℃和 0.6 ℃，均小于 1.0℃；且僅有個別測點在某些時刻的誤差絕對值超過1 ℃，因此判定數(shù)值計算結果與試驗結果的整體吻合度較高，而個別測點在某些時刻的誤差絕對值超過1 ℃，可能是由于試驗測量誤差導致的。

3 工程應用

3.1 參數(shù)與優(yōu)化方案

淮南礦區(qū)潘一煤礦副井凍結法施工的設計參數(shù)如下：①井筒凈直徑：7.500 m；②擬開挖直徑：8.000 m；③內圈管凍結管布置圈徑：14.00 m；④內圈管數(shù)量：27；⑤內圈管管間距：1.625 m；⑥內圈管凍結管壁溫度：-32 ℃；⑦內圈管凍結管尺寸：φ0.140 m；⑧外圈管凍結管布置圈徑：18.00 m；⑨外圈管數(shù)量：43；⑩外圈管凍結管間距：1.314 m；1○外圈管凍結管壁溫度：-32 ℃；12○外圈管凍結管尺寸：φ0.140 m；13○排間距：2.000 m。其中 200 m 層位處為粉質砂層，該層位地下水流速較大，在按照原設計方案進行凍結法施工過程中，發(fā)現(xiàn)凍結壁發(fā)展速度較慢。

對原凍結方案進行優(yōu)化，優(yōu)化參數(shù)如下：①優(yōu)化凍結管布置圈徑：18.00 m；②優(yōu)化凍結管加密角度：120°；③優(yōu)化凍結管增加數(shù)量：4 根；④優(yōu)化凍結管間距：1.006 m；⑤優(yōu)化凍結管壁溫度：-32 ℃；⑥優(yōu)化凍結管尺寸：φ0.140 m。優(yōu)化前后凍結管的布置方式如圖4。優(yōu)化前后的數(shù)值計算模型如圖5。

圖4 優(yōu)化前后凍結管的布置方式Fig.4 The arrangement of the freezing pipes before and after optimization

圖5 常規(guī)方案與優(yōu)化方案的數(shù)值計算模型Fig.5 Numerical calculation model of conventional scheme and optimization scheme

3.2 優(yōu)化效果

根據(jù)飽和砂土的冰點試驗結果，認為凍結砂土層達到-1 ℃時即形成有效的凍結壁，因此-1 ℃等溫線即為凍結壁的輪廓線，將優(yōu)化前后的不同流速的地下水作用下凍結溫度場在30、60、90 d 的分布規(guī)律繪制成圖，常規(guī)方案與優(yōu)化方案溫度場分析如圖6。

由圖6 通過分析發(fā)現(xiàn)：在流動的地下水作用下，下游區(qū)域的凍結壁的交圈時間早于上游區(qū)域，當凍結時間達到90 d 時，下游區(qū)域的凍結壁厚度明顯大于上游的凍結壁。其主要原因是：地下水在通過上游凍結管作用區(qū)域流向下游區(qū)域時，通過對流傳熱作用將上游區(qū)域的冷量攜帶至下游區(qū)域，在冷量的疊加作用下，下游區(qū)域的凍結速率增加，因此其交圈時間要早于上游區(qū)域；與此同時，部分水流在繞過凍結區(qū)域流向下游的過程中，也會將上游、以及兩側凍結區(qū)域的冷量攜帶至下游區(qū)域，因此下游區(qū)域的凍結壁厚度要大于其他位置。

圖6 常規(guī)方案與優(yōu)化方案溫度場分析Fig.6 Temperature field analysis of conventional scheme and optimization scheme

對凍結方案進行優(yōu)化，即對上游120°的區(qū)域的凍結孔進行加密后，優(yōu)化區(qū)域凍結壁的發(fā)展速度明顯提高，在相同的凍結時間內，優(yōu)化后的凍結壁上游位置的凍結范圍明顯大于優(yōu)化前的凍結范圍，并且地下水流速越大，凍結范圍的差別越明顯。以內、外圈凍結管均形成完整凍結壁為凍結壁完成交圈的判定標準。常規(guī)方案與優(yōu)化方案的凍結壁交圈時間隨流速變化的對比情況如圖7。

圖7 凍結壁交圈時間隨流速變化的關系Fig.7 Relationship between the closure time of frozen wall and the flow rate

由圖7 通過分析可以發(fā)現(xiàn)：2 種方案的凍結壁交圈時間均隨著地下水流速增加而增加，并且在地下水流速達到15 m/d 后，凍結壁交圈時間隨流速增加明顯延長；地下水流速為10、15、20 m/d 時，常規(guī)方案的交圈時間分別為36、49、100 d，優(yōu)化方案的交圈時間分別為29、33、43 d，優(yōu)化方案相較于常規(guī)方案在3 種流速下的凍結壁交圈時間分別提前了7、13、57 d，因此優(yōu)化方案相較于常規(guī)方案可以有效縮短凍結壁交圈時間，并且流速越大，優(yōu)化方案對凍結壁交圈時間的改善越明顯。

當凍結時間達到90 d 和120 d 時，常規(guī)方案和優(yōu)化方案凍結壁的厚度隨流速的變化規(guī)律如圖8。

由圖8 可以發(fā)現(xiàn)：凍結90 d，地下水流為7、13、19 m/d 時，常規(guī)方案的凍結壁厚度分別為6.74、6.05、5.01 m，優(yōu)化方案的凍結壁厚度分別為6.76、6.15、5.34 m，優(yōu)化方案相較于常規(guī)方案的凍結壁厚度分別提高了 0.02、0.10、0.33 m；凍結 120 d，地下水流為9、15、21 m/d 時，常規(guī)方案的凍結壁厚度分別為7.00、6.68、5.63 m，優(yōu)化方案的凍結壁厚度分別為7.04、6.73、6.20 m，優(yōu)化方案相較于常規(guī)方案的凍結壁厚度分別提高了0.04、0.05、0.57 m。因此可知，隨著地下水流速的增加，凍結壁的厚度整體呈減小的趨勢；對應相同的地下水流速，優(yōu)化方案形成的凍結壁的厚度明顯大于常規(guī)方案的凍結壁，并且隨著地下水流速的增加優(yōu)化方案的效果越明顯。

圖8 凍結壁厚度隨地下水流速的變化情況Fig.8 Variation of the thickness of frozen wall with the flow rate of groundwater

綜上，在大流速地下水作用下，通過在對雙圈凍結管布置方案中的上游位置的凍結管進行加密可以有效縮短凍結壁的交圈時間、增大凍結壁的厚度。產生這種效果的主要原因是，在較大流速的地下水作用下，凍結管布置圈的上游位置的凍結效果是影響整個凍結壁交圈時間以及厚度的決定性因素，對上游位置的凍結管進行加密處理后，增加了上流區(qū)域的冷量供應，縮小了相鄰凍結管之間的間距，從而縮短了上游區(qū)域的凍結壁的交圈時間，進而提高了整個凍結壁的凍結效率。

4 結 論

1）通過對上游120°范圍內凍結管進行加密的方式對常規(guī)凍結方案進行了優(yōu)化設計。當?shù)叵滤魉贋?0、15、20 m/d 時，常規(guī)方案的交圈時間分別為36、49、100 d，優(yōu)化方案的交圈時間分別為 29、33、43 d，優(yōu)化方案相較于常規(guī)方案在3 種流速下的凍結壁交圈時間分別提前了7、13、57 d，因此優(yōu)化方案相較于常規(guī)方案可以有效縮短凍結壁交圈時間，并且流速越大，優(yōu)化方案對凍結壁交圈時間的改善越明顯。

2）在較大流速的地下水作用下，凍結管布置圈的上游位置的凍結效果是影響整個凍結壁交圈時間以及厚度的決定性因素，對上游位置的凍結管進行加密處理后，增加了上流區(qū)域的冷量供應，縮小了相鄰凍結管之間的間距，從而縮短了上游區(qū)域的凍結壁的交圈時間，進而提高了整個凍結壁的凍結效率。

3）在對凍結方案進行優(yōu)化設計的時候，凍結管加密區(qū)域的位置以及范圍都會對最終的凍結效果產生顯著的影響，我們將在后期的研究中對該問題進行更深入的研究與分析。