李可也，袁 理

（中水淮河規(guī)劃設計研究有限公司，合肥 230601）

高海拔高寒地區(qū)冬季氣溫較低， 且白天太陽輻射強烈，近地面的氣溫會迅速上升，日溫差較大。由于高原低溫期持續(xù)時間較長，水庫、河流、渠道內極易結冰。研究表明，冰在氣溫升高時體積膨脹［1］，會對約束其膨脹的水工建筑物產生溫度膨脹力， 又稱靜冰壓力。 較大的靜冰壓力會對水工建筑物產生嚴重危害， 為研究高寒地區(qū)較大日溫差下的靜冰壓力分布規(guī)律， 本文選取冰蓋形狀和約束形式都比較簡單的渡槽進行仿真分析， 通過對其溫度場與應力場的耦合計算，研究渡槽內冰蓋在不同溫升率、不同冰厚、不同約束形式下的靜冰壓力分布。

1 靜冰壓力計算原理

靜冰壓力主要與冰蓋的溫度分布、溫升率、冰蓋厚度等因素有關， 運用溫度加載法對靜冰壓力進行仿真研究時，由于加載速率較慢，因此冰表現(xiàn)出一定的黏彈性［2］，黏彈性材料的溫度變化與材料變形的關系可表示為：

式中ε 為冰的線應變；α 為冰的溫度膨脹系數(shù)，一般認為淡水冰的溫度膨脹系數(shù)為常量；Δθ 是溫度的增量。

對式（1）的時間分量求導，得：

由式（2）可知，冰的應變速率與溫升率呈線性關系，靜冰壓力的大小受溫升率的控制。

靜冰壓力求解過程中首先需要計算出冰蓋的溫度場， 冰蓋內部溫度場的變化主要由冰表面的大氣溫度變化引起的，在有限單元法中，這是一個瞬態(tài)傳熱問題， 根據(jù)能量守恒原理， 瞬態(tài)熱平衡可表達為（以矩陣形式表示）［3］：

式中[K]為熱傳導矩陣，包含導熱系數(shù)、對流系數(shù)及輻射率和形狀系數(shù)；[C]為比熱矩陣，考慮系統(tǒng)內能的增加；{T}是節(jié)點溫度向量；{T·}為溫度對時間的導數(shù)；{Q}為冰點熱流向量，包含熱生率。

計算渡槽內冰蓋溫度場時， 認為淡水冰蓋下表層與水的接觸面的溫度為其融化點0 ℃［4］，冰蓋上表層溫度按第三類邊界條件確定， 即冰層上表面的溫度與大氣溫度相同。第三類邊界條件可表示為：

式中 Tf為熱流介質溫度；Γ 為物體邊界；ζ 為熱換系數(shù)。

求得冰蓋溫度場后， 將有限元模型中各節(jié)點的溫度作為荷載引入應力場計算中， 通過定義材料的熱膨脹系數(shù)α 將溫度場與應力場進行耦合計算，最后得到冰蓋內的應力分布規(guī)律， 作用于邊界上的冰的溫度膨脹力即為冰對水工建筑物的靜冰壓力［4］。

通過原型觀測，計算冰蓋應力場時，做以下兩點假定［5］：

（1）計算初始時，冰蓋內部無溫度膨脹產生的應力，邊界上的靜冰壓力為零。

（2）認為冰蓋邊界與渡槽凍結，仿真計算時，約束類型為四周全約束。

2 直立墻渡槽內的靜冰壓力

為研究高寒地區(qū)靜冰壓力分布規(guī)律， 本文以地處高原地區(qū)平均海拔2300 m 的某縣某直立墻渡槽為例進行了靜冰壓力的仿真研究，據(jù)氣象資料查得，該縣近年來測得的年最低溫度為-20 ℃，最大日溫升達2 h 升高10 ℃。冰蓋各項參數(shù)取值如下：

表1 冰蓋計算參數(shù)

2.1 冰蓋溫度場分布

取直立墻渡槽內2 m×8 m×0.4 m 的冰蓋為研究對象，為取得較精確的結果，將冰蓋沿厚度分為20層，模型共劃分單元32000 個，節(jié)點139041 個，冰蓋有限元模型如圖1。首先計算冰蓋溫度場分布，設冰蓋上表面度為-20 ℃，且在2 h 升高10 ℃，冰蓋下表面溫度為0 ℃，計算初始時，認為冰蓋已形成穩(wěn)定的溫度場。將穩(wěn)態(tài)溫度場作為初始條件進行迭代，最終計算得到溫升2 h 過程的瞬態(tài)溫度場分布， 穩(wěn)態(tài)溫度場云圖和沿程瞬態(tài)溫度場曲線如圖2、圖3。

圖1 冰蓋三維有限元計算模型

圖2 冰蓋穩(wěn)態(tài)溫度場分布云圖

圖3 冰蓋沿程瞬態(tài)溫度場分布曲線

由圖2 和圖3 可看出，渡槽內的穩(wěn)態(tài)溫度場呈線性分布，隨著時間的增加，冰蓋上表面溫度的升高，冰蓋溫度場產生非線性變化，但由于冰的導熱性較差，溫差的升高集中在冰蓋上部1/3以上處，且從表層向下變化程度逐漸減小，瞬態(tài)溫度場在約15 cm 處以下基本與穩(wěn)態(tài)溫度場重合。

2.2 靜冰壓力分布

得到冰蓋的溫度場之后， 將熱單元轉換為結構單元。由于假定冰層內部初始應力為零，這就要求給材料賦予線膨脹系數(shù)，設定材料參考溫度時，將穩(wěn)態(tài)溫度場計算時求出的單元溫度值分布賦予應力場分析時對應的單元。 然后將瞬態(tài)溫度場計算得到的各步溫度作為荷載施加到應力計算的模型中， 經過迭代得出各節(jié)點的應力分量［6］。由于所求結果為壓應力分布，將Y 方向的應力作為輸出結果進行分析，直立墻渡槽內0.4 m 厚冰蓋的靜冰壓力云圖分布如圖4～圖7，冰蓋沿程瞬態(tài)應力場分布曲線如圖8。

圖4 0.5h 時冰蓋的應力場分布云圖

圖5 1h 時冰蓋的應力場分布云圖

圖6 1.5h 時冰蓋的應力場分布云圖

圖7 2h 時冰蓋的應力場分布云圖

圖8 冰蓋沿程瞬態(tài)應力場分布曲線

由圖8 可看出，隨著冰蓋表層溫度的緩慢上升，冰蓋產生的溫度膨脹力受到周圍約束的作用產生靜冰壓力。靜冰壓力沿程分布不均，集中在冰蓋上部1/3以上處，下部趨近于零，表層靜冰壓力最大，靜冰壓力隨著溫度的不斷升高持續(xù)增大，2 h 的最大靜冰壓力值為268.7 kPa，這與溫度場仿真的結果相對應，即溫度變化最大處的靜冰壓力也最大。

2.3 溫升率對靜冰壓力的影響

為了研究溫升率對靜冰壓力的影響，假定大氣溫度在 2 h 內從-20 ℃升高至-15 ℃，其他條件不變， 得到冰蓋的應力場仿真結果，2 h時兩種溫升率下的Y 方向靜冰壓力沿程分布曲線如圖9。

圖9 兩種溫升率下2h 時冰蓋沿程靜冰壓力分布曲線

由圖9 可看出，兩種溫升率下的靜冰壓力均集中于冰蓋1/3 以上處，靜冰壓力最大處同樣位于冰蓋表層。溫升率的大小決定了靜冰壓力的大小，隨著溫升率的增大，靜冰壓力也變大，符合之前的理論分析。

2.4 冰厚對靜冰壓力的影響

為研究不同冰蓋厚度對靜冰壓力的影響， 根據(jù)該縣當?shù)卣鎸嵉谋遒Y料進行假定， 當冰蓋厚度分別為0.3，0.4，0.5，0.6，0.7，0.8 m， 按2 h 溫升10 ℃，其他條件不變時，冰蓋表層最大靜冰壓力值分布如圖10。

圖10 不同厚度冰蓋表層靜冰壓力分布曲線

由圖10 可看出，冰蓋表層靜冰壓力的最大值隨著冰蓋厚度的增大而增大。且由仿真分析可知，不同厚度的冰蓋，其靜冰壓力仍集中作用于冰蓋表層值冰厚15cm處，說明冰蓋表層的溫升產生的靜冰壓力僅對冰蓋上部產生比較大的影響，對冰蓋下部影響很小。

3 結語

（1）采用有限元法對冰蓋的靜冰壓力進行仿真計算具有相當精度， 在進行高寒地區(qū)水工建筑物設計時， 其表層最大靜冰壓力仿真結果可作為估算冰蓋靜冰壓力極值的參考。

（2）冰蓋瞬態(tài)溫度場變化集中于冰蓋上部15cm以上處。

（3）溫升率對冰蓋表面靜冰壓力有較大影響，靜冰壓力隨溫升率增加而增加。

（4）冰蓋厚度對表層最大靜冰壓力值有一定影響，但靜冰壓力仍集中于冰蓋上部15cm 以上處。