李可也，袁 理

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

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

1 靜冰壓力計算原理

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

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

對式（1）的時間分量求導(dǎo)，得：

由式（2）可知，冰的應(yīng)變速率與溫升率呈線性關(guān)系，靜冰壓力的大小受溫升率的控制。

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

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

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

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

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

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

（1）計算初始時，冰蓋內(nèi)部無溫度膨脹產(chǎn)生的應(yīng)力，邊界上的靜冰壓力為零。

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

2 直立墻渡槽內(nèi)的靜冰壓力

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

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

2.1 冰蓋溫度場分布

取直立墻渡槽內(nèi)2 m×8 m×0.4 m 的冰蓋為研究對象，為取得較精確的結(jié)果，將冰蓋沿厚度分為20層，模型共劃分單元32000 個，節(jié)點139041 個，冰蓋有限元模型如圖1。首先計算冰蓋溫度場分布，設(shè)冰蓋上表面度為-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 可看出，渡槽內(nèi)的穩(wěn)態(tài)溫度場呈線性分布，隨著時間的增加，冰蓋上表面溫度的升高，冰蓋溫度場產(chǎn)生非線性變化，但由于冰的導(dǎo)熱性較差，溫差的升高集中在冰蓋上部1/3以上處，且從表層向下變化程度逐漸減小，瞬態(tài)溫度場在約15 cm 處以下基本與穩(wěn)態(tài)溫度場重合。

2.2 靜冰壓力分布

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

圖4 0.5h 時冰蓋的應(yīng)力場分布云圖

圖5 1h 時冰蓋的應(yīng)力場分布云圖

圖6 1.5h 時冰蓋的應(yīng)力場分布云圖

圖7 2h 時冰蓋的應(yīng)力場分布云圖

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

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

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

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

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

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

2.4 冰厚對靜冰壓力的影響

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

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

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

3 結(jié)語

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

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

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

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