白 旭，楊蘇杰，田于逵

（1.江蘇科技大學(xué)船舶與海洋工程學(xué)院，江蘇 鎮(zhèn)江 212003；2.中國(guó)船舶科學(xué)研究中心，江蘇 無錫 214082）

0 引 言

海水飛沫是導(dǎo)致船舶結(jié)冰的主要原因，96%的船舶結(jié)冰與海水飛沫有關(guān)[1]。船舶表面結(jié)冰會(huì)降低船舶的穩(wěn)性和降低結(jié)構(gòu)穩(wěn)定性，影響船上設(shè)備的使用，對(duì)船舶的安全和性能造成了極大影響。研究結(jié)冰特性對(duì)建立有效的船舶結(jié)冰預(yù)報(bào)方法具有重要意義[2]。

現(xiàn)有的結(jié)冰研究多為根據(jù)當(dāng)?shù)貧夂驐l件進(jìn)行的結(jié)冰預(yù)報(bào)，或是利用CFD 軟件對(duì)整體或構(gòu)件進(jìn)行結(jié)冰模擬，使用的公式也是根據(jù)已有數(shù)據(jù)推算的經(jīng)驗(yàn)公式[3]。針對(duì)船舶的結(jié)冰模型有許多缺陷，如忽略了水與冰的換熱、表面粗糙度的變化等[4]。分析及試驗(yàn)表明，冰的生長(zhǎng)規(guī)律和尖端形貌會(huì)影響冰形的預(yù)測(cè)可靠性。由于成核率和穩(wěn)定性的限制，所以冰晶生長(zhǎng)的試驗(yàn)數(shù)據(jù)準(zhǔn)確性差，而采用數(shù)值模擬冰晶生長(zhǎng)的研究較少，且均有局限性。

為了掌握結(jié)冰的演變機(jī)理，開發(fā)船舶結(jié)冰特性的精細(xì)預(yù)測(cè)方法具有重要意義[5]?；诮鹌澅?朗道理論的相場(chǎng)法是模擬凝固組織的一種有效方法，它能夠描述復(fù)雜的晶體形貌。該方法已成為模擬微觀組織演化的一種重要的、通用性極強(qiáng)的方法，成為模擬凝固過程中微觀組織演化的有力工具[6]。

1998 年，Kobayashi[7]首先提出了擴(kuò)散邊界概念，并在此基礎(chǔ)上建立了第一個(gè)相場(chǎng)模型，進(jìn)行了金屬的凝固模擬；Wheeler 等[8]在模擬中將相場(chǎng)參數(shù)和實(shí)際物理量相結(jié)合，建立了WBM 模型；Karma 等[9]建立了KKS模型，并消除了界面厚度的影響，擴(kuò)大了模擬的過冷度范圍，加快了計(jì)算速度；Valeria等[10]提出了一個(gè)模型用來計(jì)算海水結(jié)冰過程中的鹽分流通現(xiàn)象，描述了高鹽度溶液與低鹽度冰相分離的時(shí)間演化過程，但是沒有分析鹽度的影響。

國(guó)內(nèi)學(xué)者對(duì)于相場(chǎng)模型的研究也有一定進(jìn)展。張玉妥等[11]率先使用相場(chǎng)法模擬純物質(zhì)金屬的凝固現(xiàn)象和等軸枝晶的生長(zhǎng)，通過數(shù)值計(jì)算顯示了等軸枝晶的形態(tài)；鄧儒超[6]進(jìn)行了海水冰晶生長(zhǎng)的相場(chǎng)模擬，發(fā)現(xiàn)過冷度和各向異性越大，冰晶的生長(zhǎng)速度越大；陳尚海[5]分析了對(duì)流情況下的海水冰晶生長(zhǎng)凝固過程，枝晶尖端生長(zhǎng)速度隨著流速、擾動(dòng)強(qiáng)度和過冷度的增大而增大；唐楊新[12]給出了適用于海冰演化的多相相場(chǎng)模型的簡(jiǎn)單推導(dǎo)，應(yīng)用Alber-Zhu的相場(chǎng)模型來描述不同種類的海冰之間界面的演化，介紹了海冰增長(zhǎng)的相場(chǎng)模型和不同海冰在構(gòu)型力驅(qū)動(dòng)下相變演化模型的新證明方法。

本文基于相場(chǎng)法，建立船舶結(jié)冰微觀結(jié)構(gòu)特征預(yù)測(cè)的數(shù)值模擬方法，分析不同鹽度條件下結(jié)冰的微觀結(jié)構(gòu)演化特征。海水飛沫引起的結(jié)冰過程是一個(gè)多晶核生長(zhǎng)的復(fù)雜競(jìng)爭(zhēng)過程，它是由振動(dòng)觸發(fā)的非均勻形核開始的，但鹽度對(duì)成核率的影響無從得知，因此在分析時(shí)僅考慮單晶核條件下鹽度對(duì)冰晶生長(zhǎng)的影響。

1 考慮鹽度的相場(chǎng)模型

1.1 純物質(zhì)相場(chǎng)模型

Wheeler 相場(chǎng)模型使用相場(chǎng)和溫度場(chǎng)方程來實(shí)現(xiàn)對(duì)凝固現(xiàn)象的模擬[13]，對(duì)于純物質(zhì)的凝固，相場(chǎng)和溫度場(chǎng)控制方程如下[8]：

相場(chǎng)方程為

式中：ε(θ)為考慮各向異性影響的相場(chǎng)參數(shù)，ε(θ)=εˉ( 1 +γcos(jθ))，γ為各向異性強(qiáng)度，j為各向異性模數(shù)，本文中取j= 6，εˉ為無量綱界面厚度；m為相場(chǎng)遷移率；α為界面自由能；Δ為無量綱過冷度。

溫度場(chǎng)方程為

式中：p'(?)為固相分?jǐn)?shù)p(?)對(duì)?的導(dǎo)數(shù)，p'(?)=30?2(1-?2)；u為無量綱溫度，u=(T-TM)/(TM-T0)。

純物質(zhì)相場(chǎng)模型無法模擬海水的凝固現(xiàn)象，要進(jìn)行海水的凝固模擬以及分析鹽度的影響，需要引入溶質(zhì)的影響。

1.2 考慮鹽度的海水相場(chǎng)模型

研究鹽度對(duì)海水凝固的影響時(shí)，將海水視作二元溶液。在海水凝固過程中，由于溶質(zhì)鹽在水中的溶解度大于在冰中的溶解度，溶質(zhì)會(huì)在固液界面前沿富集，隨著界面推進(jìn)，固液界面將會(huì)向液相排出溶質(zhì)，進(jìn)而影響凝固過程。因此，在進(jìn)行海水凝固模擬時(shí)，需要引入溶質(zhì)再分配的影響。

基于菲克定律，引入溶質(zhì)場(chǎng)方程，以實(shí)現(xiàn)凝固時(shí)計(jì)算域內(nèi)的溶質(zhì)變化：

式中，C取海水中鹽的摩爾分?jǐn)?shù)，D為無量綱擴(kuò)散系數(shù)，Ds為固相擴(kuò)散率，Dl為液相擴(kuò)散率，k0為平衡分配系數(shù)。

計(jì)算域內(nèi)的溶質(zhì)濃度變化會(huì)導(dǎo)致熔點(diǎn)發(fā)生變化，因此引入過冷度的實(shí)時(shí)動(dòng)態(tài)變化來實(shí)現(xiàn)鹽度在凝固過程的影響。

式中，CP為海水比熱容，ΔT為溶液過冷度，mL為液相線斜率，C0為溶液初始濃度，L為單位體積潛熱。

添加溶質(zhì)影響后的相場(chǎng)方程為

溫度場(chǎng)方程為

2 參數(shù)取值和計(jì)算方法

2.1 參數(shù)取值

本文所取物性參數(shù)均按照濃度為3.50%的海水溶液選取，具體取值如表1所示[13]。

表1 海水溶液物性參數(shù)Tab.1 Physical properties of seawater solution

2.2 有限差分法及離散方式

對(duì)于相場(chǎng)方程，空間離散采用中心差分方法，時(shí)間離散采用顯式差分方法。另外，對(duì)于方程中出現(xiàn)的拉普拉斯算子，采用九點(diǎn)差分方法進(jìn)行離散：

計(jì)算時(shí)為了保證計(jì)算結(jié)果收斂，計(jì)算過程的穩(wěn)定性由相場(chǎng)方程決定，時(shí)間步長(zhǎng)需滿足以下條件：

2.3 初始條件和邊界條件

本文計(jì)算域設(shè)置為正方形，網(wǎng)格數(shù)設(shè)置為150×150，網(wǎng)格間距Δx和Δy在無量綱化后均為0.03，時(shí)間步長(zhǎng)Δt= 0.0001 s。

具體晶核位置和條件設(shè)置如公式（11）所示：

相場(chǎng)、溫度場(chǎng)和濃度場(chǎng)均采用Neumann邊界條件，即

2.4 可靠性驗(yàn)證

對(duì)冰晶凝固過程進(jìn)行模擬，并將其和實(shí)際自然條件下生成的冰晶進(jìn)行對(duì)比。使用上文設(shè)置參數(shù)，網(wǎng)格劃分為500×500, 模擬時(shí)間為7500Δt，得到結(jié)果見圖1。通過比較模擬的冰晶形狀和現(xiàn)實(shí)的冰晶形狀[16]，發(fā)現(xiàn)實(shí)驗(yàn)結(jié)果吻合良好，可以進(jìn)行進(jìn)一步模擬計(jì)算。

3 鹽度對(duì)冰晶生長(zhǎng)的影響

世界海水鹽度分布大致為3%～4%[14]。鹽度對(duì)海水凝固的影響，除了會(huì)實(shí)時(shí)影響凝固場(chǎng)內(nèi)的過冷度，還會(huì)直接影響海水的熔點(diǎn)。

取6 個(gè)不同的過冷度作為溫度條件，在每個(gè)溫度下，鹽度從3.00%到3.90%，每隔0.1%取一個(gè)值，計(jì)算不同鹽度的冰晶生長(zhǎng)行為，分析鹽度對(duì)冰晶生長(zhǎng)過程的影響。

3.1 鹽度與溫度關(guān)系

取初始溫度T0為240.85 K、243.85 K、246.85 K、249.85 K、252.85 K、255.85 K，在各溫度條件下引入鹽度引起的熔點(diǎn)變化。

根據(jù)Bodnar[15]的研究，用熔點(diǎn)相比零攝氏度的下降溫度來描述鹽度對(duì)熔點(diǎn)的影響，文中計(jì)算的鹽度條件對(duì)應(yīng)的熔點(diǎn)溫度下降如表2所示。

表2 鹽度對(duì)應(yīng)的熔點(diǎn)下降Tab.2 Salinities corresponding to melting point depressions

三種不同溫度條件下，不同鹽度的過冷度參數(shù)如表3所示。

表3 初始溫度參數(shù)Tab.3 Initial temperature parameters

圖2 顯示了溫度為240.85 K 時(shí)的結(jié)冰圖像，在低過冷度條件下，冰晶以六角冰晶的形式生長(zhǎng)，但是并未生成分枝。圖3 為枝晶尖端速度圖像，尖端速度盡管略有波動(dòng)，但還是逐漸過渡到平穩(wěn)期，這與Ivantsov枝晶生長(zhǎng)理論一致。

3.2 鹽度對(duì)結(jié)冰過程的影響

表4 列舉了三個(gè)溫度條件下枝晶尖端速度的平均值?？傮w來說，枝晶生長(zhǎng)速度隨著溫度降低而增加。鹽度影響下的溶質(zhì)場(chǎng)如圖4 所示。固相中的濃度要低于液相中的濃度，這在固液界面處尤為明顯。這是由于溶質(zhì)在固相中的溶解度低于其在液相中的溶解度，隨著固液界面的推進(jìn)，溶質(zhì)從固相中被排出并富集在固液界面處，與凝固理論相一致。

表4 各溫度尖端速度均值Tab.4 Mean value of tip velocity at different temperatures

圖5 展示了在不同溫度下，鹽度的變化對(duì)尖端速度產(chǎn)生的影響，具體表現(xiàn)為：隨著鹽度的增加，枝晶生長(zhǎng)速度降低。這是因?yàn)辂}度的增長(zhǎng)導(dǎo)致了溶液熔點(diǎn)的降低，在外界溫度不變的情況下，從而使溶液的過冷度減小，減緩了枝晶的生長(zhǎng)。

鹽度變化會(huì)影響枝晶的生長(zhǎng)速度，從圖5 可以看出，溫度為240.85 K 時(shí)，鹽度從3.00%增長(zhǎng)到3.90%，尖端速度從5.280Δx/s變?yōu)?.271Δx/s，降低了19.1%，相比整體速度降低比較明顯。

溫度在240～246 K 之間時(shí)，由于過冷度較大，結(jié)冰速度較快，速度變化量明顯，分別為1.0085Δx/s、1.0985Δx/s、0.9039Δx/s，而當(dāng)溫度在249～255 K 時(shí)，尖端速度的數(shù)值較小，此時(shí)鹽度引起的速度變化量也很小，分別為0.0059Δx/s、0.0041Δx/s、0.0038Δx/s。

溫度從240.85 K 升至246 K 時(shí)，尖端速度隨鹽度變化比例分別為19.1%、33.0%和76.9%，比例逐漸上升，這是因?yàn)殡S著溫度升高，過冷度降低，鹽度引起的過冷度變化占原有過冷度的比例越來越大，導(dǎo)致速度變化越來越明顯。而當(dāng)溫度為249.85K 時(shí)，此時(shí)過冷度進(jìn)一步降低，尖端速度發(fā)生數(shù)量級(jí)的變化，因此鹽度引起的速度改變量與之前相比也顯著減少，減少比例分別為27.3%、29.1%和46.3%。

可以看出，結(jié)冰速度隨鹽度變化而變化，且改變量顯著，如果以海洋結(jié)構(gòu)物作為考慮對(duì)象，那么不管是結(jié)冰尺寸還是結(jié)冰時(shí)間都會(huì)顯著增大，帶來的結(jié)冰量變化會(huì)更加明顯。

4 結(jié) 論

為了提高船舶結(jié)構(gòu)結(jié)冰特性的預(yù)測(cè)精度，本文利用Wheeler 相場(chǎng)模型，并用有限差分法對(duì)不同溫度和鹽度下過冷液體晶體生長(zhǎng)過程進(jìn)行了數(shù)值模擬，提取枝晶尖端速度作為衡量鹽度對(duì)冰晶生長(zhǎng)影響的標(biāo)準(zhǔn)，得到如下結(jié)論：

（1）鹽度對(duì)海水結(jié)冰的影響主要來自于鹽度與熔點(diǎn)之間的聯(lián)系。鹽度越大，海水熔點(diǎn)越低，在相同溫度下導(dǎo)致過冷度越小。因此鹽度越大，海水結(jié)冰越慢；

（2）在凝固過程中，溶質(zhì)從固相析出并被排到液相，被排出的溶質(zhì)在固液界面前富集，進(jìn)一步影響了海水的凝固；

（3）鹽度對(duì)冰晶生長(zhǎng)的影響與過冷度有關(guān)。過冷度在240～246 K 之間時(shí)，結(jié)冰尖端速度較大，鹽度引起的速度變化量分別為1.0085Δx/s、1.0985Δx/s、0.9039Δx/s，而當(dāng)過冷度在249～255 K 時(shí)，此時(shí)尖端速度數(shù)值較小，鹽度引起的尖端速度變化量隨之減小，分別為0.0059Δx/s、0.0041Δx/s、0.0038Δx/s。

本文模擬了鹽度對(duì)微觀尺度海水結(jié)冰過程的影響。傳統(tǒng)的針對(duì)船舶結(jié)構(gòu)結(jié)冰的研究尚有許多不足之處，如《結(jié)構(gòu)大氣結(jié)冰標(biāo)準(zhǔn)》[17]中的結(jié)冰計(jì)算主要依據(jù)經(jīng)驗(yàn)公式，結(jié)冰過程中的許多機(jī)理未得到解釋。通過進(jìn)行微觀尺度的結(jié)冰研究，可以對(duì)結(jié)冰過程中的機(jī)理進(jìn)行分析。下一步將著重于研究液滴撞擊前后的微觀結(jié)冰現(xiàn)象，探究撞擊過程中的能量轉(zhuǎn)化，并將其運(yùn)用于船舶結(jié)構(gòu)宏觀結(jié)冰研究中。