波浪場(chǎng)中浮冰搖蕩運(yùn)動(dòng)特性

朱仁慶， 張曦， 李志富 (江蘇科技大學(xué) 船舶與海洋工程學(xué)院，江蘇 鎮(zhèn)江 212003)

隨著全球氣候變暖，海冰冰厚變薄使得北極資源的開發(fā)以及航線開通成為可能。許多國(guó)家競(jìng)相儲(chǔ)備開發(fā)北極資源和北極航道，準(zhǔn)備隨時(shí)進(jìn)駐北極。因此，合理認(rèn)識(shí)和利用冰緣區(qū)內(nèi)波-冰作用規(guī)律，有利于加強(qiáng)對(duì)冰區(qū)內(nèi)波浪傳播特性的了解和提高極地資源利用率具有重要意義。對(duì)于浮冰與波浪的相互作用，Harms[1-2]通過(guò)試驗(yàn)回歸出恒密度厚度的二維浮冰漂移運(yùn)動(dòng)估算公式。當(dāng)浮冰的橫搖周期略大于波浪周期時(shí)，浮冰具有最大漂移速度。Lever等[3]在水池試驗(yàn)中研究了冰山自身形狀對(duì)其在波浪中運(yùn)動(dòng)的影響。Huang等[4]研究發(fā)現(xiàn)浮冰漂移速度隨著橫搖周期與波浪周期的比值增大而增大。Bennetts等[5]在試驗(yàn)中發(fā)現(xiàn)浮冰的存在阻礙了波浪傳播，且隨著波陡的增大，浮冰對(duì)波浪的阻礙作用增強(qiáng)。Yiew等[6]發(fā)現(xiàn)當(dāng)波陡較大時(shí)，波浪沖洗現(xiàn)象對(duì)浮冰的運(yùn)動(dòng)響應(yīng)有明顯抑制作用。郭春雨等[7]試驗(yàn)研究了浮冰在波浪中的縱向運(yùn)動(dòng)，同時(shí)發(fā)現(xiàn)了波浪沖洗現(xiàn)象。波浪在冰區(qū)的傳播比例隨著波長(zhǎng)和波頻率的減小而增大，且波浪傳播比例增大也增強(qiáng)了其對(duì)浮冰的破壞能力[8]。

由于模型試驗(yàn)具有成本過(guò)高、存在設(shè)備干擾的缺點(diǎn)，數(shù)值模擬方法被越來(lái)越多的用于研究波-冰作用。其中，勢(shì)流分析方法因其成本低、計(jì)算速度快等優(yōu)點(diǎn)受到了學(xué)者的青睞。Shen等[9]通過(guò)理論研究了反射波、阻力系數(shù)和附加質(zhì)量影響浮體漂移運(yùn)動(dòng)的的因素，如波的反射、附加質(zhì)量等。Meylan[10]研究發(fā)現(xiàn)浮冰運(yùn)動(dòng)很大程度上依賴于浮冰形狀，波浪散射主要取決于入射波方向。Bennetts等[11]研究了非均質(zhì)浮冰對(duì)波浪的散射。Smith等[12]對(duì)均質(zhì)浮冰對(duì)波浪散射現(xiàn)象求解進(jìn)行了研究，并在此方法的基礎(chǔ)上發(fā)展了非均質(zhì)浮冰對(duì)波浪散射的求解。Meylan等[13]通過(guò)改進(jìn)莫里森方程研究大波長(zhǎng)規(guī)則波中浮冰的縱蕩運(yùn)動(dòng)，與試驗(yàn)結(jié)果具有較好的一致性。Toffoli等[14]提出了波浪海冰相互作用的試驗(yàn)?zāi)Ｐ蛠?lái)驗(yàn)證冰緣區(qū)波浪衰減的理論模型。由于波-冰作用中固有的強(qiáng)非線性現(xiàn)象無(wú)法通過(guò)勢(shì)流理論精確計(jì)算出，因此計(jì)算流體力學(xué)動(dòng)態(tài)模擬技術(shù)的迅速發(fā)展研究波-冰作用起到了重要的作用。Bai等[15]研究表明粘流軟件比勢(shì)流軟件更適合于模擬浮冰運(yùn)動(dòng)。Huang等[16]通過(guò)粘流軟件模擬了波浪沖洗冰體表面和波浪散射現(xiàn)象。

在浮冰與波浪的相互作用中，浮冰在波浪作用下會(huì)發(fā)生非形變運(yùn)動(dòng)。與此同時(shí)，波浪場(chǎng)也會(huì)由于浮冰的運(yùn)動(dòng)而發(fā)生改變。目前關(guān)于波-冰作用的研究對(duì)于相互影響已經(jīng)有了充分了解，但對(duì)波浪沖洗現(xiàn)象對(duì)浮冰運(yùn)動(dòng)的影響研究較少。波浪沖洗冰體表面現(xiàn)象雖有分析，但是很少對(duì)波浪沖洗現(xiàn)象、浮冰運(yùn)動(dòng)以及波浪參數(shù)之間定性定量關(guān)系的研究。本文通過(guò)粘流模型數(shù)值模擬了波浪中浮冰的運(yùn)動(dòng)，利用頻譜分析研究了浮冰運(yùn)動(dòng)響應(yīng)。分析了波浪沖洗現(xiàn)象對(duì)浮冰運(yùn)動(dòng)的影響，對(duì)比在有/無(wú)波浪沖洗現(xiàn)象工況下的浮冰運(yùn)動(dòng)。并研究了不同波浪參數(shù)下波浪沖洗現(xiàn)象對(duì)浮冰垂蕩、縱搖運(yùn)動(dòng)的影響。

1 波-冰作用數(shù)值計(jì)算方法

1.1 流體運(yùn)動(dòng)控制方程

浮冰區(qū)域流體運(yùn)動(dòng)滿足連續(xù)方程和動(dòng)量方程?？紤]到湍流流動(dòng)，通常采用時(shí)均法將流體運(yùn)動(dòng)控制方程中的各物理量分解為2部分，分別為時(shí)間平均值和相對(duì)于平均值的脈動(dòng)值，得到雷諾時(shí)均方程組，其表達(dá)式為：

(1)

(2)

Pk-ρ(ε-ε0)+Sk

(3)

(4)

1.2 數(shù)值造波理論

本文通過(guò)速度入口造波方式生成波浪，通過(guò)在入口處和出口處分別通過(guò)松弛區(qū)域消波技術(shù)和阻尼消波方法實(shí)現(xiàn)消波從而減少壁面反射。本文選取斯托克斯五階波浪[18-19]，使模擬的波浪更符合實(shí)際的規(guī)則波。在造波邊界處的速度和波面瞬時(shí)升高滿足以下條件：

x方向速度：

(5)

z方向速度：

(6)

波面瞬時(shí)升高：

(7)

其中各系數(shù)分別為：

式中：ω、d、k分別是圓頻率、水深和波數(shù)；c=cosh(kd)；λ為常系數(shù)，其余各項(xiàng)系數(shù)參考文獻(xiàn)[20]。

入口消波區(qū)的消波原理是通過(guò)在動(dòng)量方程中添加源項(xiàng)qφ強(qiáng)迫波浪在速度入口附近轉(zhuǎn)換成入射波：

qφ=-γρ(φ-φ*)

(8)

式中：γ為強(qiáng)迫系數(shù)；φ為數(shù)值模擬中的動(dòng)量輸運(yùn)方程的結(jié)果；φ*為原動(dòng)量方程的理論解。冰水相互作用的數(shù)值水池模型如圖1所示。

圖1 冰水相互作用數(shù)值水池Fig.1 Diagram of the numerical wave tank

出口消波區(qū)的原理是通過(guò)在Z軸方向上添加阻尼完成消波作用，其中阻尼項(xiàng)為：

(9)

式中：xsd為起始消波位置；xed為結(jié)束消波位置；f1、f2、nd為阻尼消波的參數(shù)；w為Z方向的速度。

1.3 浮冰運(yùn)動(dòng)方程

在隨浮冰平移和旋轉(zhuǎn)的運(yùn)動(dòng)坐標(biāo)系中，浮冰的運(yùn)動(dòng)方程為：

(10)

(11)

式中：m是浮冰質(zhì)量；V是浮冰運(yùn)動(dòng)速度；Ω是浮冰旋轉(zhuǎn)角速度；I是浮冰轉(zhuǎn)動(dòng)慣量；M是浮冰受到的合外力矩；f是浮冰受到的的合外力。

2 波-冰作用計(jì)算模型設(shè)置

2.1 波浪與計(jì)算模型驗(yàn)證

2.1.1 波浪驗(yàn)證

根據(jù)極地冰緣區(qū)出現(xiàn)的波浪高度和頻率分布規(guī)律和波浪參數(shù)與浮冰特征尺寸的比例關(guān)系，在相同的波浪高度下，對(duì)浮冰在一系列波長(zhǎng)(波陡偏大)波浪中的運(yùn)動(dòng)進(jìn)行了數(shù)字模擬和分析研究，工況1～5波高恒定為0.04 m，并且入射波長(zhǎng)為0.6、0.7、0.8、0.9、1.0 m。

液面波高探測(cè)點(diǎn)位置選取為圓柱所在位置，其中選取工況5的波浪來(lái)記錄探測(cè)點(diǎn)處液面升高過(guò)程，將斯托克斯五階波理論解與其監(jiān)測(cè)到的波高時(shí)歷曲線數(shù)據(jù)進(jìn)行比較。對(duì)于波陡小的波浪，造波效果與其設(shè)置的波高H方向的網(wǎng)格尺寸Δz有重要的關(guān)聯(lián)，而能否生成目標(biāo)波浪決定了計(jì)算結(jié)果的準(zhǔn)確與否。因此需要再進(jìn)行收斂性分析，這里分別選取3組不同大小的自由液面處波高方向網(wǎng)格和3組不同時(shí)間步進(jìn)行生成波浪進(jìn)行分析，Δz分別為H/5、H/10及H/15。波長(zhǎng)方向的網(wǎng)格尺寸Δx是波高方向網(wǎng)格尺寸Δz的4倍。網(wǎng)格與時(shí)間步設(shè)置如表1所示。

表1 不同網(wǎng)格尺寸與時(shí)間步設(shè)置Table 1 Different mesh size and time step settings

如圖2所示，數(shù)值水池生成的波浪在頻率、波高等波浪參數(shù)在網(wǎng)格B上的誤差最小，與斯托克斯五階波理論解吻合較好。因此，波高方向上網(wǎng)格尺寸Δz為H/10被選取進(jìn)行數(shù)值模擬。將該套網(wǎng)格在3種不同時(shí)間步下進(jìn)行了計(jì)算，如圖3所示，數(shù)值模擬時(shí)間步Δt為0.002 s的工況的結(jié)果收斂。因此，出于效率和資源上的綜合考慮，波高方向網(wǎng)格尺寸Δz為H/10，時(shí)間步長(zhǎng)Δt為0.002 s被選入進(jìn)行后面的數(shù)值模擬。

圖2 波高方向不同網(wǎng)格尺寸Fig.2 Dimensions of different mesh numbers in the wave height direction

圖3 不同時(shí)間步工況與理論值對(duì)比Fig.3 Comparison of wave conditions and theoretical solution values at different time steps

2.1.2 模型驗(yàn)證

本文選取文獻(xiàn)[6]中單個(gè)浮冰與規(guī)則波相互作用試驗(yàn)研究，提取試驗(yàn)中浮冰的垂蕩幅值時(shí)歷曲線來(lái)驗(yàn)證數(shù)值模擬計(jì)算結(jié)果。浮冰模型為圓柱體，半徑為0.2 m，厚度為0.015 m，密度為636 kg/m3，浮冰邊緣附有薄壁結(jié)構(gòu)，薄壁結(jié)構(gòu)高為0.05 m，厚度為0.025 m，在浮冰周圍施加彈簧力以防止浮冰產(chǎn)生大幅度漂移。在保證無(wú)池壁反射的前提下，適當(dāng)縮減了水池的長(zhǎng)度和寬度，數(shù)值模型中水池寬2 m，長(zhǎng)7 m，水深0.83 m。入射波頻率為1.25 Hz，波幅為0.085 m。浮冰中心距離造波入口4 m處。計(jì)算模型中網(wǎng)格劃分遵循上述網(wǎng)格參數(shù)設(shè)置，在浮冰周圍加密以保證計(jì)算精度，網(wǎng)格劃分圖如圖4所示，其中對(duì)于k-ε模型需要滿足30≤y+≤300，一般以接近30為佳，因此y+取值為30，第1層網(wǎng)格厚度為0.000 47 m。

圖4 重疊網(wǎng)格劃分Fig.4 The arrangement of the mesh

由圖5可知，試驗(yàn)值與數(shù)值模擬中的浮冰垂蕩運(yùn)動(dòng)基本一致，浮冰垂蕩運(yùn)動(dòng)峰值以及到達(dá)相鄰峰值所需的周期大致吻合，由于偏差較小，數(shù)值模擬結(jié)果相較于實(shí)驗(yàn)結(jié)果沒有明顯偏高或偏低趨勢(shì)，因此可以將偏差視為數(shù)值計(jì)算或?qū)嶒?yàn)測(cè)量精度的判斷依據(jù)，說(shuō)明本文采用的數(shù)值模擬方法是合理和準(zhǔn)確的。

圖5 浮冰垂蕩運(yùn)動(dòng)對(duì)比Fig.5 The comparison of ice floe heave motion

2.2 三維模型和網(wǎng)格劃分

如圖6所示，數(shù)值水池尺寸為4.5 m×2 m×0.8 m，水深0.5 m。數(shù)值水池的前后邊界設(shè)置為對(duì)稱邊界條件；左側(cè)邊界、上側(cè)邊界和下側(cè)邊界設(shè)置為速度入口邊界條件；右側(cè)邊界設(shè)置為壓力出口邊界條件。入口以及出口消波區(qū)區(qū)域大小均取為一個(gè)波長(zhǎng)。

圖6 波-冰作用的三維數(shù)值水池Fig.6 Three-dimensional numerical tank of wave-ice interaction

將浮冰視為剛體，由于實(shí)際海況上浮冰互相之間經(jīng)常會(huì)受到碰撞和波浪的沖洗，浮冰邊角的磨損使其形狀近似于圓柱，如圖7所示。將浮冰模型進(jìn)行處理便于計(jì)算，故簡(jiǎn)化為圓柱體和帶一定高度薄壁的圓柱體2種形狀，浮冰中心與入口處的距離是1.5 m。其中圓柱參數(shù)：直徑0.5 m；厚度0.05 m；吃水0.03 m；密度600 kg/m3。為對(duì)無(wú)波浪沖洗現(xiàn)象進(jìn)行研究，將圓柱浮冰周圍增設(shè)高為0.05 m薄壁，研究表明薄壁對(duì)浮冰運(yùn)動(dòng)的影響可忽略不計(jì)[21]。

圖7 2種浮冰模型Fig.7 The geometry of two kind of ice floes

本文通過(guò)重疊網(wǎng)格技術(shù)模擬浮冰運(yùn)動(dòng)，在浮體運(yùn)動(dòng)區(qū)域以及自由液面區(qū)域進(jìn)行網(wǎng)格加密，如圖8為重疊網(wǎng)格。

圖8 重疊網(wǎng)格Fig.8 Overset mesh

3 波-冰作用計(jì)算結(jié)果及分析

波浪的參數(shù)對(duì)浮冰運(yùn)動(dòng)響應(yīng)有著關(guān)鍵的影響。為了進(jìn)一步研究浮冰與入射波的相互作用，通過(guò)頻譜分析研究了不同波長(zhǎng)下浮冰垂蕩η3、縱搖η5運(yùn)動(dòng)響應(yīng)。

3.1 波浪沖洗現(xiàn)象

浮冰在波浪作用下會(huì)發(fā)生大幅的搖蕩運(yùn)動(dòng)，海水沖洗浮冰表面的現(xiàn)象頻頻出現(xiàn)。由數(shù)值模擬可以看出，波長(zhǎng)增加，波浪沖洗現(xiàn)象會(huì)隨之減弱。圖9為在不同波長(zhǎng)條件下一個(gè)周期的浮冰表面的上浪現(xiàn)象。圖9(a)是波長(zhǎng)為0.6 m時(shí)波浪沖洗狀態(tài)，浮冰表面上浪嚴(yán)重；圖9(b)是波長(zhǎng)為1.0 m時(shí)的波浪沖洗狀態(tài)，浮冰表面未完全被水體覆蓋，與圖9(a)相比，其所受的波浪沖洗現(xiàn)象較弱。

圖9 不同波長(zhǎng)的波浪沖洗冰體表面過(guò)程Fig.9 Overwash process of waves with different wavelength

圖10為λ=0.6 m時(shí)的浮冰搖蕩運(yùn)動(dòng)頻譜分析圖。在相同的波浪激勵(lì)頻率作用下，沒有受到波浪沖洗的浮冰在運(yùn)動(dòng)幅值上要大于受到波浪沖洗的浮冰。產(chǎn)生上述差異的主要原因是在波浪與浮冰相互作用中存在波浪沖洗冰體現(xiàn)象，一方面堆積在冰體表面的水阻礙浮冰的運(yùn)動(dòng)；另一方面冰體表面的水體會(huì)產(chǎn)生碰撞繼而引起水體破碎，使波浪產(chǎn)生一定程度的衰減，使得浮冰接受的波浪能量較少。

圖10 波長(zhǎng)λ=0.6 m的浮冰運(yùn)動(dòng)幅值對(duì)比Fig.10 Comparison of the motion amplitudes of the two ic floes, the wave condition was λ=0.6 m

圖11為不同波長(zhǎng)下的波浪沖洗現(xiàn)象對(duì)浮冰的垂蕩和縱搖的影響。浮冰的垂蕩與縱搖隨著波長(zhǎng)減小而減小。當(dāng)波長(zhǎng)較短時(shí)，波陡較大，波浪沖洗冰體表面現(xiàn)象較嚴(yán)重。由于堆積在冰體表面的水體較多以及波浪能量衰減較大，一定程度上阻礙了浮冰運(yùn)動(dòng)，使其運(yùn)動(dòng)幅度減弱。

圖11 無(wú)護(hù)欄浮冰在不同波長(zhǎng)下的運(yùn)動(dòng)響應(yīng)對(duì)比Fig.11 Comparison of the motion response of floating ice without guardrail at different wavelengths

由此可見，在給定波高情況下，波長(zhǎng)越小，即波陡越大，波浪沖洗現(xiàn)象越嚴(yán)重。浮冰運(yùn)動(dòng)在波浪沖洗現(xiàn)象作用下受到了一定的抑制，波浪沖洗現(xiàn)象越強(qiáng)，其對(duì)浮冰的抑制作用越強(qiáng)。

3.2 波浪散射現(xiàn)象

通過(guò)模擬加薄壁的浮冰在波浪中的運(yùn)動(dòng)來(lái)研究波-冰作用。通過(guò)數(shù)值模擬可以看出，短波長(zhǎng)下波-冰作用更加劇烈。浮冰的存在一定程度上阻礙了波浪的運(yùn)動(dòng)，波浪場(chǎng)存在較明顯的散射現(xiàn)象。

圖12(a)對(duì)應(yīng)波長(zhǎng)較短的工況，波浪場(chǎng)產(chǎn)生了明顯的波浪散射現(xiàn)象，浮冰前側(cè)和后側(cè)的波浪液面形狀由于散射的作用發(fā)生了明顯的擾動(dòng)現(xiàn)象。圖12(b)對(duì)應(yīng)波長(zhǎng)較長(zhǎng)的工況，波浪散射現(xiàn)象較弱，浮冰周圍波形沒有明顯變化。

圖12 有護(hù)欄浮冰在不同波長(zhǎng)下與波浪相互作用Fig.12 Floes with guardrails interact with waves at different wavelengths

如圖13所示，在波浪激勵(lì)頻率下，圖13(a)中波長(zhǎng)越小，浮冰的垂蕩運(yùn)動(dòng)幅值也越小。當(dāng)波長(zhǎng)較短時(shí)，浮冰大大阻礙了波浪的傳遞，同時(shí)也對(duì)波浪場(chǎng)造成了一定散射。由于波浪的散射消耗了一定的波能，使得浮冰的運(yùn)動(dòng)響應(yīng)和受到的波浪力也隨之減弱。隨著波長(zhǎng)增加，浮冰相對(duì)波長(zhǎng)的尺度減小，進(jìn)而浮冰對(duì)波浪的作用也減小，浮冰接收的波浪能量變大，因而波浪激勵(lì)頻率下的浮冰運(yùn)動(dòng)也隨之劇烈。圖13(b)中，浮冰的縱搖運(yùn)動(dòng)隨著波長(zhǎng)的減小愈發(fā)劇烈。然而當(dāng)波長(zhǎng)變化至0.8 m時(shí)，波-冰特征長(zhǎng)度比值為1.6，浮冰的縱搖運(yùn)動(dòng)開始減弱。這是因?yàn)椴ɡ松⑸涞某潭入S著波陡的增加愈發(fā)嚴(yán)重，波能耗散也愈發(fā)嚴(yán)重，從而使浮冰受到的作用力和運(yùn)動(dòng)響應(yīng)減弱。

圖13 有護(hù)欄浮冰在不同波長(zhǎng)下的運(yùn)動(dòng)對(duì)比Fig.13 Comparison of the movement amplitude of ice floes with guardrails at different wavelengths

綜上所述，隨著波陡增大，波浪的散射現(xiàn)象越來(lái)越明顯。波浪散射現(xiàn)象對(duì)浮冰垂蕩具有一定程度的抑制作用。研究發(fā)現(xiàn)波浪散射現(xiàn)象越強(qiáng)，浮冰垂蕩受到的抑制作用也越強(qiáng)。其中，對(duì)于浮冰的縱搖運(yùn)動(dòng)，隨著波陡的增大，縱搖運(yùn)動(dòng)幅值增大，當(dāng)波陡增大到一定值時(shí)，浮冰縱搖運(yùn)動(dòng)隨著波陡的減小而減小。

4 結(jié)論

1)通過(guò)對(duì)比有沖洗/無(wú)沖洗浮冰的運(yùn)動(dòng)發(fā)現(xiàn)，波浪沖洗現(xiàn)象對(duì)浮冰的運(yùn)動(dòng)具有明顯的抑制作用。

2)隨著波陡變大，波浪沖洗現(xiàn)象和波浪散射現(xiàn)象增強(qiáng)，波浪沖洗現(xiàn)象對(duì)浮冰運(yùn)動(dòng)的抑制作用增強(qiáng)。同樣地，波浪散射現(xiàn)象對(duì)浮冰垂蕩運(yùn)動(dòng)具有抑制作用，波長(zhǎng)越小，浮冰垂蕩運(yùn)動(dòng)受到的抑制作用越強(qiáng)。

3)波陡越大，浮冰的縱搖運(yùn)動(dòng)越大，當(dāng)波陡增大到一定值時(shí)，由于波浪能量耗散嚴(yán)重，浮冰縱搖運(yùn)動(dòng)幅值隨波陡增大呈現(xiàn)出下降趨勢(shì)。

本文計(jì)算方法和分析結(jié)論對(duì)于進(jìn)一步研究浮冰在波浪下的運(yùn)動(dòng)以及波浪中浮冰與結(jié)構(gòu)物的相互作用具有參考價(jià)值。