沈 杰,白 旭
(江蘇科技大學 船舶與海洋工程學院,江蘇 鎮(zhèn)江 212003)
寒區(qū)海域冬季氣候寒冷,海氣交換強烈,濕度很大,大部分時間相對濕度都在95%以上,表現為多霧、濃霧等極端氣候,與其他海域相比,空氣中的水汽凝結以及海浪拍擊等極易在船舶結構表面上發(fā)生凝霜、冰釉等結冰現象[1],這些新增重量會降低船舶干舷高度,甚至引起船身傾斜,造成船舶傾覆[2]。如1965年1月在白令海上10艘蘇聯船舶由于船體嚴重結冰造成結構失穩(wěn)進而傾覆[3],圖1為寒區(qū)船舶上的結冰現象。
圖 1 寒區(qū)船舶結冰現象[4]Fig. 1 Ship icing in cold regions
寒區(qū)船舶上層建筑結冰多為大氣結冰,目前關于大氣結冰的研究多集中風力機結冰和飛機結冰。如李巖等[5]在自行設計的冰風洞中以繞軸旋轉圓柱為對象研究旋轉模型的結冰問題,在該實驗中分析了轉速、結冰時間和圓柱直徑對圓柱結冰形狀的影響規(guī)律。陳維建[6]提出一種適用于各種結冰氣象條件下飛機機翼結冰過程的數值模擬方法,并以NACA0012翼型為模型進行霜冰、明冰的數值計算,將計算結果與試驗數據進行對比,驗證了該方法的有效性。顧聲龍[7]將NURLS-S809翼型作為研究對象,對在不同環(huán)境參數下,風力機葉片表面所結霜冰的形狀和結冰量進行研究,結果顯示結冰厚度和結冰量與來流風速的關系較大,隨著風速的增大,結冰量也增加,兩者之間基本呈線性關系。吳俊杰等[8]為了確定機翼結冰過程中過冷水滴的運動軌跡,以NACA0012翼型為模型,分析在溫度、時間步長和風速等因素不變的條件下,水滴直徑和來流攻角對水滴軌跡的影響。結果表明,隨著過冷水滴直徑的增大,水滴在翼面的分布范圍越廣;來流攻角不同時,過冷水滴在翼面的撞擊區(qū)域有很大不同。杜雁霞等[9]基于液/固相變的基本理論,對飛機結冰過程中的液/固相變傳熱特性進行了研究,建立了水膜和冰層生長模型,將該模型的計算結果與冰風洞的實驗結果進行對比,驗證了模型的可靠性。之后采用該模型對風速、溫度等各參數對冰層生長及速率的影響進行分析,結果顯示冰層生長及速率不但與風速、溫度等參數有關還與液、固相區(qū)的內部傳熱特性有關,但這些研究都是關于較高風速下結構表面的結冰。
由于寒區(qū)船舶上層建筑存在大量的桿件結構,本文采用Fluent和FENSAP-ICE相結合的方法對桿件結構結冰開展數值模擬,選取風速作為敏感參數,分析風速變化對結冰厚度和結冰量的影響。
寒區(qū)結構結冰按冰的類型可分為霜冰、明冰、濕雪和混合冰[10],圖2為霜冰、明冰與風速和溫度的關系。Magne[12]為了確定挪威西海岸Brosviks?ta山的大氣結冰條件,在山上4個位置處各放置一個探測當地風速、溫度、液態(tài)水濃度等影響結冰參數的裝置,通過該實驗確定了該地區(qū)大氣結冰的氣象條件范圍:1)氣溫在-15 ℃~0 ℃范圍內,當氣溫高于0 ℃時不結冰,當氣溫低于-15 ℃時,空氣中水滴會直接結晶從空氣中掉落而不與結構接觸;2)風速不超過4級(7.9 m/s);3)空氣中液態(tài)水濃度為0.05~0.25 g/m3。
圖 2 明冰、霜冰與風速和溫度的關系[11]Fig. 2 Relationship between glaze and rime and wind speed and temperature[11]
圖3 為霜冰的形成過程簡圖。寒區(qū)空氣中的冷水滴與結構表面接觸后還未擴散就完全凝結形成的冰為霜冰[13]。由于水滴幾乎立刻凝結,霜冰的每個冰晶顆粒之間存在著一定量的空氣并維持半球狀,使其具有粗糙的外形和不透明的性質,而且密度和強度都比明冰較低[14]。
圖 3 霜冰的形成過程簡圖Fig. 3 The formation of rime
由于寒區(qū)船舶結構結冰的特殊性,通常對其的研究都是通過現場試驗進行的,而受到我國地理位置的影響,現場試驗是非常耗時且成本極高的[15]。因此,本文結合Fluent和Fensap-ice軟件進行桿件結構的結冰數值模擬研究,計算流程如圖4所示。分為如下步驟:1)使用Fluent計算結構周圍的流場分布,以得到結構周圍空氣流速ua等參數;2)使用式(1)和式(2)求解水滴的運動軌跡和流場中水滴容積系數的分布;3)使用公式計算得到結構表面控制體內收集到的水滴質量;4)使用式(3)和式(4)計算各控制體內凍結冰的質量;5)使用式(5)得到各控制體內冰的厚度;6)將網格向前推進到霜冰的高度處,并重復上述過程直到給定的截止時間。
式中,cf為空氣比熱,為空氣/冰交界面溫度,Levap為蒸發(fā)潛熱,Tice,rec為 冰面溫度,T為遠場空氣溫度。
孟繁鑫等[16]通過自建的引射式結冰風速研究了不同條件下靜止圓柱表面的霜冰和明冰分布,本文選取該試驗的一組霜冰結冰條件(見表1)進行數值模擬計算,并與試驗結果對比分析。
圖 4 結冰計算流程Fig. 4 Icing calculation process
圖5 中的x軸坐標表示圓柱表面一點與X軸負半軸之間的夾角(如圖5中φ),沿順時針方向為正。圖6中本文的冰形計算結果和文獻中的試驗結果變化趨勢相同,數值也在同一個數量級,兩者之間的偏差最大為3.5%,可證明本文采用的結冰數值模擬方法是可靠的。
表 1 計算條件Tab. 1 Calculation conditions
圖 5 x軸坐標定義Fig. 5 X-axis coordinate definition
圖 6 計算結果與文獻結果對比Fig. 6 Comparison of calculation results with literature results
由于本文選取的風向與桿件結構軸向方向垂直,可簡化為二維模型進行數值計算,如圖7所示。圖8為本文計算域模型,其中inlet為速度入口邊界條件,outlet為壓力出口邊界條件,上下邊界為壁面邊界條件。結冰數值模擬條件由1.1節(jié)中的結冰氣象條件確定(詳細計算參數見表2),由圖2可知本文的結冰類型應為霜冰,此外由于不同風速下所結霜冰的密度不同,本文參考Macklin[17]根據冰風洞試驗得到的霜冰密度經驗公式(見公式9)確定霜冰密度。
圖 7 流場示意圖Fig. 7 Flow field diagram
圖 8 計算域模型Fig. 8 Computational domain model
表 2 計算工況Tab. 2 Calculation conditions
上式的應用條件為:
圖 9 不同風速下圓柱表面的局部收集系數分布Fig. 9 Distribution of icing thickness on a cylindrical surface at different wind speeds
圖 10 不同風速下圓柱表面的結冰厚度分布Fig. 10 Distribution of icing thickness on a cylindrical surface at different wind speeds
圖 11 不同風速下圓柱表面的結冰量Fig. 11 Amount of icing on a cylindrical surface at different wind speeds
圖9 ~圖11分別為風速1.0 m/s,2.0 m/s,3.0 m/s,4.0 m/s,5.0 m/s和7.0 m/s時,桿件表面局部水收集系數、冰形和結冰量分布圖。撞擊到桿件表面的水滴,沿桿件表面的分布是不均勻的,使用式(8)可得到結構表面控制體內的局部水收集系數,進而可確定結構表面的水滴撞擊范圍,最終可確定結冰范圍,如圖9中1~7 m/s風速下桿件表面局部水收集系數的分布情況與圖10中結冰的分布情況是一致的。此外,從圖10可看出:1)隨著風速從1 m/s增加到7 m/s,結構表面的冰厚逐漸增加,因為風速增大使單位時間內撞擊到結構上的過冷水滴數量增多,同時風速增大也會帶走更多結冰釋放的潛熱,使結冰更為迅速;2)隨著風速從1 m/s增加到7 m/s,冰厚最大值的位置逐漸向左移動,因為水滴速度隨來流速度的增加而增加,在相同時間內水滴更快的到達結構表面,水滴運動軌跡的偏轉更小,更多的水滴撞擊到0°位置處,冰厚最大值也隨之向0°方向移動;3)風速為7 m/s時,桿件表面的冰厚最大值為16.59 mm僅比風速為1 m/s時增加了31.3%,因為1 m/s時霜冰的密度僅為7 m/s時的22.8%,當質量一定時,冰的密度越小體積就會越大,因此雖然1 m/s時結冰量僅為7 m/s的18.1%,但冰厚并未遠小于7 m/s時的冰厚。從圖11可看出,隨著風速從1 m/s增加7 m/s,桿件表面的結冰量近似線性增長。
本文以寒區(qū)船舶桿件結構結冰為背景,采用Fluent和Fensap-ice相結合的方法對桿件結構結冰開展數值模擬,選取風速作為敏感參數,分析了風速變化對結冰厚度和結冰量的影響,得出如下結論:
1)隨著風速從1 m/s增加到7 m/s,桿件表面冰厚最大值的位置逐漸向左移動,因為水滴速度隨來流速度的增加而增加,在位移相同時水滴更快的到達結構表面,水滴運動軌跡的偏轉更小,更多的水滴撞擊到0°位置處,冰厚最大值也隨之向0°方向移動。
2)風速為7 m/s時桿件表面的冰厚并未遠大于1 m/s時的冰厚。因為風速為1 m/s時桿件表面所結霜冰的密度僅為7 m/s時的22.8%,當質量一定時,冰的密度越小體積就會越大,因此雖然1 m/s時結冰量僅為7 m/s的18.1%,但冰厚并未遠小于7 m/s時的冰厚。
3)在風速為1~7 m/s區(qū)間內,隨著風速的增加,桿件表面的結冰量近似的呈線性方式增加。