橫搖-垂蕩復(fù)合激勵(lì)下液艙晃蕩的壓力特性研究

張會(huì)霞，李 師，鄒昶方

(1.江蘇海洋大學(xué) 海洋工程學(xué)院，江蘇 連云港 222005；2.江蘇省海洋資源開發(fā)研究院，江蘇 連云港 222005)

0 引言

隨著能源需求的不斷增加，海上油氣資源開采和運(yùn)輸日益增加，國際市場對(duì)大型液貨船的需求日益增長。液貨船具有寬度較大、裝載率高的特性。在遇到一定范圍的波浪時(shí)，船體運(yùn)動(dòng)會(huì)導(dǎo)致劇烈的液艙晃蕩，并且晃蕩載荷會(huì)對(duì)船體結(jié)構(gòu)和運(yùn)動(dòng)響應(yīng)造成重大安全隱患。

針對(duì)液貨船的液艙晃蕩問題，蔡忠華[1]基于VOF法建立數(shù)值計(jì)算模型，結(jié)果與模型晃蕩試驗(yàn)對(duì)比，驗(yàn)證了數(shù)值方法解決晃蕩問題的可行性，對(duì)含有內(nèi)部構(gòu)件的VLGC 液艙數(shù)值模型研究，發(fā)現(xiàn)內(nèi)部構(gòu)件減小共振頻率下晃蕩沖擊載荷的同時(shí)，改變了艙內(nèi)液體的固有頻率，當(dāng)外激勵(lì)頻率在此頻率附近時(shí)，沖擊載荷反而更大。張秋艷[2]使用ADIDA 軟件模擬縱蕩激勵(lì)下二維矩形液艙晃蕩的數(shù)值模型，其液面波動(dòng)情況與Faltinsen 的解析解對(duì)比，驗(yàn)證了該模型的有效性。對(duì)二維彈性液艙進(jìn)行數(shù)值模擬，分析了不同壁厚、不同頻率對(duì)液體晃蕩的影響。金鑫[3]基于有限差分法對(duì)NS 方程進(jìn)行數(shù)值離散，建立縱蕩與升沉耦合激勵(lì)下流體晃蕩三維數(shù)值模型，基于勢流理論推導(dǎo)出縱蕩與升沉耦合激勵(lì)下理想流體晃蕩解析解，驗(yàn)證了數(shù)值模型。結(jié)果表明當(dāng)激勵(lì)頻率的和值或者差值接近流體系統(tǒng)自然頻率時(shí)，耦合激勵(lì)也可以激發(fā)共振晃蕩波。Zhang[4]采用有限差分法對(duì)垂直激勵(lì)槽內(nèi)的非線性液體晃蕩問題進(jìn)行數(shù)值模擬，通過σ 變換將不規(guī)則液體域映射到一個(gè)矩形區(qū)域，在迭代過程中對(duì)自由液面進(jìn)行預(yù)測，該方法可以精確計(jì)算出垂直晃蕩過程中的自由液面高度與晃蕩壓力。Ida M.Strand[5]對(duì)帶有柔性側(cè)壁的二維矩形液艙進(jìn)行數(shù)值分析，研究晃蕩與柔性艙壁之間的耦合問題，通過分析和數(shù)值計(jì)算，得到了2 種不同膜長度下由搖擺激勵(lì)引起的耦合本征頻率和傳遞函數(shù)，分析結(jié)果與數(shù)值結(jié)果一致。寧德志[6]采用高階邊界元法，模擬了二維剛性矩形液艙在水平與垂直激勵(lì)下的晃蕩運(yùn)動(dòng)，結(jié)果表明在耦合晃蕩問題中，水平激勵(lì)對(duì)波浪的形狀起主導(dǎo)作用。Mashy D.Green[7]基于SPH法，提出改進(jìn)公式，能夠捕獲液艙的晃蕩頻率，并且準(zhǔn)確模擬液艙長時(shí)間晃蕩的液面形狀。Constantin[8–9]對(duì)矩形液艙在垂直簡諧激勵(lì)下的晃蕩運(yùn)動(dòng)進(jìn)行實(shí)驗(yàn)研究，定性和定量解釋了阻力飽和點(diǎn)和耗散效應(yīng)，以及區(qū)分不同晃蕩狀態(tài)的指標(biāo)。建立基于SPH 法的數(shù)值模型和基于彈跳球模型的等效力學(xué)模型，結(jié)果與實(shí)驗(yàn)對(duì)比發(fā)現(xiàn)2 種模型在晃蕩的初始階段均具有較好的預(yù)測能力，但僅有SPH 模型能預(yù)測后續(xù)的流態(tài)。

本文針對(duì)多自由度下液艙晃蕩問題，采用基于VOF 法的多相流模型，建立多組復(fù)合激勵(lì)工況，探究復(fù)合激勵(lì)下液艙晃蕩的壓力變化特點(diǎn)以及頻率對(duì)壓力影響，研究結(jié)果對(duì)艙體結(jié)構(gòu)設(shè)計(jì)具有指導(dǎo)意義。

1 無關(guān)性分析

本文計(jì)算采用二維矩形液艙[10]，液艙的尺寸，橫搖中心與橫搖角大小如圖1 所示，其中P點(diǎn)是25%H時(shí)自由液面高度，也是壓力監(jiān)測點(diǎn)所在位置。

圖1 液艙模型尺寸與橫搖參數(shù)(m)Fig.1 Tank size and rolling parameters (m)

在25%裝載高度下，時(shí)間步大小設(shè)置為0.002 s，橫搖角度為12°，設(shè)置網(wǎng)格數(shù)量分別為1 萬、2.4 萬、3 萬、4 萬、5 萬和6.25 萬，分別計(jì)算在P點(diǎn)的晃蕩載荷變化，見圖2。結(jié)果顯示，隨著網(wǎng)格數(shù)量的增大，計(jì)算精度逐漸提高。當(dāng)網(wǎng)格數(shù)量達(dá)到2.4 萬后，網(wǎng)格數(shù)量的提高對(duì)計(jì)算結(jié)果影響大幅降低，計(jì)算結(jié)果逐漸收斂。

圖2 不同網(wǎng)格數(shù)量下晃蕩壓力計(jì)算結(jié)果Fig.2 Calculation results of sloshing pressure under different grid numbers

2 橫搖激勵(lì)下晃蕩壓力特點(diǎn)

液體固有頻率計(jì)算公式：

式中：fn為n階固有頻率，Hz；ɡ為當(dāng)?shù)刂亓铀俣龋籐表示液艙運(yùn)動(dòng)方向自由液面長度；H表示液艙內(nèi)液體裝載高度。液艙縱搖時(shí)，L等于液艙長度。本文研究包括液艙橫搖狀態(tài)，故L等于液艙寬度。當(dāng)n等于1 時(shí)，即液體的最低固有頻率，此時(shí)的頻率最容易被激發(fā)，且在此頻率附近晃蕩最劇烈，對(duì)艙壁產(chǎn)生的沖擊力也最大。

計(jì)算可得25%H下，矩形液艙的固有頻率為0.81 Hz，理論計(jì)算結(jié)果與實(shí)際頻率有一定的誤差，所以在數(shù)值模擬時(shí)，選擇在固有頻率0.81 Hz 上下增減0.05 Hz，選擇5 個(gè)頻率，工況如表1 所示。為進(jìn)一步確保數(shù)值計(jì)算的精確度，在3 萬網(wǎng)格的基礎(chǔ)上選擇時(shí)間步為0.001 15 s。

表1 25%裝載率橫搖工況Tab.1 Rolling condition of 25% loading rate

使用Fluent 中UDF 驅(qū)動(dòng)液艙模型橫搖轉(zhuǎn)動(dòng)，角速度為：

其中：fr為橫搖激勵(lì)頻率，由于轉(zhuǎn)動(dòng)角度為12°，經(jīng)計(jì)算可得：

由圖3 可知，矩形液艙25%H，橫搖激勵(lì)頻率在0.81 Hz 附近時(shí)，艙壁受到的砰擊載荷最大，所以此時(shí)的共振頻率在0.81 Hz 附近。當(dāng)激勵(lì)頻率小于共振頻率時(shí)，晃蕩載荷隨頻率增加的速率較快，在大于共振頻率時(shí)，晃蕩載荷隨頻率減小的速率較慢。圖4 為單純橫搖激勵(lì)fr=0.81 Hz 時(shí)，晃蕩壓力曲線呈現(xiàn)明顯的雙峰特征，由于液體猛烈砰擊在豎直壁面，所以首峰數(shù)值較大，持續(xù)時(shí)間較短，次峰值較小，持續(xù)時(shí)間較長。

圖3 25%H 橫搖激勵(lì)晃蕩壓力隨頻率分布特點(diǎn)Fig.3 Distribution characteristics of sloshing pressure with rolling excitation under 25%H

圖4 fr=0.81 Hz 點(diǎn)P 處壓力時(shí)程曲線Fig.4 Pressure time history curve at point P with fr=0.81 Hz

3 復(fù)合激勵(lì)下液艙晃蕩壓力特點(diǎn)

基于上述液艙的單純橫搖運(yùn)動(dòng)，對(duì)液艙運(yùn)動(dòng)增加垂蕩激勵(lì)，激勵(lì)幅值為0.1 m。通過UDF 驅(qū)動(dòng)液艙模型的垂蕩運(yùn)動(dòng)，垂向運(yùn)動(dòng)線速度為：

由于垂蕩幅值為0.1 m，經(jīng)過計(jì)算得：

其中，fh為垂蕩激勵(lì)頻率。選擇fh時(shí)，綜合考慮橫搖激勵(lì)頻率與共振頻率，制定復(fù)合激勵(lì)工況如表2 所示。

表2 復(fù)合激勵(lì)工況Tab.2 Compound excitation working condition

3.1 數(shù)值模型驗(yàn)證

針對(duì)這一數(shù)值模型對(duì)解決垂蕩激勵(lì)下液艙晃蕩問題的有效性，選擇Constantin L[9]中的液艙模型，在30%裝載高度，8.3 Hz 垂蕩簡諧激勵(lì)頻率，采用本文數(shù)值模型模擬的自由液面變化如圖5 所示。圖6 為Constantin L實(shí)驗(yàn)中的自由液面變化。可以清晰看出，數(shù)值模擬的自由液面變化與實(shí)驗(yàn)結(jié)果基本吻合，說明此數(shù)值模型可以很好的模擬垂蕩激勵(lì)下的晃蕩問題。

圖5 數(shù)值模擬垂蕩激勵(lì)下自由液面變化Fig.5 Numerical simulation of free surface change under heaving excitation

圖6 Constantin L 垂蕩激勵(lì)實(shí)驗(yàn)自由液面變化Fig.6 Free surface change in Constantin L heaving excitation experiment

3.2 復(fù)合激勵(lì)下晃蕩壓力時(shí)域特征

單純橫搖激勵(lì)下，每個(gè)運(yùn)動(dòng)周期內(nèi)晃蕩壓力變化規(guī)律基本相同。橫搖激勵(lì)在共振頻率附近時(shí)壓力曲線呈現(xiàn)雙峰特征，首峰峰值大且持續(xù)時(shí)間短，次峰相對(duì)較小但持續(xù)時(shí)間較長。復(fù)合激勵(lì)時(shí)，影響晃蕩載荷的頻率因素增多，艙體運(yùn)動(dòng)變得復(fù)雜，此時(shí)壓力曲線會(huì)出現(xiàn)多個(gè)運(yùn)動(dòng)周期組合的大周期。

由圖7(a)可知，當(dāng)fr=0.71 Hz，fh=0.6 Hz，兩者都遠(yuǎn)離一階共振頻率時(shí)，晃蕩壓力數(shù)值較小，結(jié)合氣液分布圖可知此時(shí)晃蕩劇烈程度較低，在單個(gè)運(yùn)動(dòng)周期內(nèi)液面沿著艙壁上下輕微起伏，再由于兩自由度復(fù)合激勵(lì)運(yùn)動(dòng)，其中影響頻率較多，壓力曲線表現(xiàn)出更大的周期性。當(dāng)fr=fh=0.71 Hz 時(shí)，每個(gè)運(yùn)動(dòng)周期內(nèi)艙體運(yùn)動(dòng)變化相同，所以其壓力變化基本一致，由于遠(yuǎn)離一階固有頻率，壓力數(shù)值仍然較小，但由于激勵(lì)頻率相同，壓力曲線的周期性降低。當(dāng)fr=0.71 Hz，fh=0.81 Hz，即fr遠(yuǎn)離共振頻率，但fh在一階共振頻率時(shí)，壓力數(shù)值有較明顯的上升。當(dāng)fr=0.71 Hz，fh=1.13 Hz，即垂蕩激勵(lì)頻率大于一階共振頻率時(shí)，壓力數(shù)值繼續(xù)增大，晃蕩形式變得復(fù)雜，壓力雙峰大小交替周期性變化。與圖3 相比，壓力時(shí)程曲線呈現(xiàn)的周期增多，在fr遠(yuǎn)離共振頻率時(shí)，壓力未出現(xiàn)瞬時(shí)快速增大的現(xiàn)象，fr=0.71 Hz 系列工況的自由液面變化特征如圖7(b)所示，此系列復(fù)合激勵(lì)自由液面變化很小，整個(gè)運(yùn)動(dòng)周期內(nèi)自由液面沿艙壁爬升幅度很低，所以此時(shí)艙壁受到壓力載荷較小。

圖7 工況1～工況9 晃蕩特征Fig.7 Sloshing characteristics under conditions 1～9

圖8(a)表明，fr處于一階共振頻率時(shí)，壓力曲線與單純橫搖激勵(lì)時(shí)類似，壓力載荷曲線呈明顯的雙峰特征，首峰峰值大且持續(xù)時(shí)間短，次峰相對(duì)較小，持續(xù)時(shí)間較長。區(qū)別在于當(dāng)fh=0.68 Hz 遠(yuǎn)離共振頻率時(shí)，壓力曲線呈周期性變化，一開始首峰峰值大且持續(xù)時(shí)間短，次峰峰值較小但時(shí)間相對(duì)較長，兩者差距很大，隨著運(yùn)動(dòng)周期增加，首峰逐漸降低，兩者差距逐漸縮小再增大，差距在0.61～2.19 kPa 之間。當(dāng)fh=fr=0.81 Hz，兩者都處于一階共振頻率時(shí)，艙內(nèi)液體晃蕩較為劇烈，在每個(gè)運(yùn)動(dòng)周期內(nèi)的壓力變化規(guī)律相同，首峰壓力峰值基本穩(wěn)定在3.04 kPa 左右，但相較于fr=0.81 Hz 的其他復(fù)合激勵(lì)工況，此時(shí)的晃蕩壓力峰值略小。fh=0.9 Hz 時(shí)，壓力的變化規(guī)律與fh=0.6 Hz 時(shí)相似，但峰值變化更加明顯，首峰與次峰的差距變化在0.18～3.6 kPa。此系列壓力時(shí)程曲線與圖4fr=0.81 Hz單純橫搖激勵(lì)相比，兩者都有明顯的雙峰特征，但由于復(fù)合激勵(lì)下液艙運(yùn)動(dòng)形式會(huì)周期性變化，所以復(fù)合激勵(lì)下的壓力時(shí)程曲線周期性更多。fr=0.81 Hz 復(fù)合激勵(lì)系列工況的自由液面變化特征如圖8 (b)所示，此系列工況下液體晃蕩劇烈，自由液面沿艙壁爬升較高且出現(xiàn)射流，液體內(nèi)部裹挾了大量氣泡。

圖8 工況19～工況27 晃蕩特征Fig.8 Sloshing characteristics under conditions 19～27

圖9 表明，當(dāng)fr=0.91 Hz，fh=0.6 Hz 時(shí)，激勵(lì)頻率都遠(yuǎn)離一階共振頻率，壓力載荷相對(duì)圖8 有明顯的降低，其壓力峰值基本在1～1.5 kPa 左右，壓力曲線的首峰與次峰呈周期性交替變化，雙峰數(shù)值差距不大，持續(xù)時(shí)間均較短。當(dāng)fh=0.81 Hz 時(shí)，雙峰大小周期性變化特征更加顯著，且此時(shí)轉(zhuǎn)換周期比fh=0.6 Hz 時(shí)更長。當(dāng)fh=fr=0.91 Hz 時(shí)，壓力曲線與單純橫搖激勵(lì)類似，單個(gè)運(yùn)動(dòng)周期內(nèi)壓力載荷較大且持續(xù)時(shí)間較短，此時(shí)由于遠(yuǎn)離一階共振頻率，液體晃蕩相對(duì)較為緩和。當(dāng)fh=1.13 Hz 時(shí)，激勵(lì)頻率均遠(yuǎn)離一階共振頻率，相比于fh=0.6 Hz，壓力載荷首峰與次峰大小依舊呈現(xiàn)周期性變化，此時(shí)壓力大小有一定增大。與圖4 相比，壓力曲線的雙峰特征不明顯且數(shù)值較小。此系列晃蕩的液面變化特征如圖9(b)所示，相比圖8(b)此系列工況液體晃蕩有明顯緩和，但液體沿艙壁爬升仍較高，內(nèi)部裹挾有少量氣泡。

圖9 工況37～工況45 晃蕩特征Fig.9 Sloshing characteristics under conditions 37～45

綜合各工況下壓力時(shí)程曲線與氣液兩相分布圖，可以發(fā)現(xiàn)當(dāng)矩形液艙激勵(lì)頻率fr≠fh時(shí)，其每個(gè)運(yùn)動(dòng)周期內(nèi)的壓力載荷變化會(huì)出現(xiàn)多種類型的壓力載荷變化，多個(gè)運(yùn)動(dòng)周期組合形成的更大周期，這表明復(fù)合激勵(lì)內(nèi)存在更多影響其壓力載荷的頻率。

當(dāng)激勵(lì)頻率fr=fh時(shí)，壓力載荷曲線特征與單純橫搖相同，每個(gè)運(yùn)動(dòng)周期內(nèi)的壓力載荷變化穩(wěn)定，兩者在一階共振頻率0.81 Hz 時(shí)壓力曲線有明顯的雙峰特征，此時(shí)晃蕩最為劇烈，在液體爬升艙壁的極限位置出現(xiàn)射流且液體內(nèi)部裹挾有大量的氣泡。

3.3 復(fù)合激勵(lì)下晃蕩壓力頻域特征

統(tǒng)計(jì)表2 中各工況下的壓力載荷大小，分析其隨頻率變化特點(diǎn)。再以fr=0.81 Hz 為例，對(duì)其各個(gè)工況的壓力時(shí)程曲線進(jìn)行傅里葉分析，研究影響壓力載荷的主要頻率。

圖10 中淺色點(diǎn)代表矩形液艙在25%裝載率時(shí)，單純橫搖激勵(lì)下晃蕩壓力載荷大小，黑色表示表2 各工況下晃蕩壓力載荷大小。可以看出，垂蕩頻率對(duì)晃蕩載荷的影響基本呈V 字變化，拐點(diǎn)出現(xiàn)在fr=fh時(shí)，此時(shí)復(fù)合激勵(lì)下的晃蕩載荷與單純橫搖激勵(lì)相差不大。當(dāng)fr≠fh時(shí)，fh的加入對(duì)晃蕩載荷有增大作用，尤其是當(dāng)fh>fr時(shí)，晃蕩載荷的增加更為顯著，且隨著的增加，晃蕩載荷繼續(xù)呈增加的趨勢。可以看出，圖10(a)和圖10(d)，fr遠(yuǎn)離0.81 Hz 的工況，晃蕩載荷在垂蕩激勵(lì)的放大作用下始終沒有超過圖10(b)和圖10(e)的最小值。同樣的圖10(b)和圖10(e)，在接近一階共振頻率的工況，晃蕩載荷也始終未達(dá)到圖10(c)中的最小值。說明在本文工況中，垂蕩激勵(lì)對(duì)于晃蕩載荷的放大效果有限，起主導(dǎo)作用的始終是橫搖激勵(lì)。

圖10 橫搖與垂蕩復(fù)合激勵(lì)工況下晃蕩壓力變化Fig.10 Variation of sloshing pressure under compound excitation of roll and heave

從圖11 可以看出，各組復(fù)合激勵(lì)下晃蕩載荷隨橫搖激勵(lì)頻率的變化規(guī)律與單純橫搖激勵(lì)時(shí)相同。晃蕩壓力在fr=0.81 Hz 時(shí)最大，隨著橫搖頻率逐漸遠(yuǎn)離共振頻率，復(fù)合激勵(lì)的晃蕩壓力也在逐漸減小。垂蕩激勵(lì)對(duì)晃蕩壓力的放大作用在fr=0.81 Hz，即一階共振頻率時(shí)最為明顯，其效果隨著fr逐漸遠(yuǎn)離一階共振頻率而逐漸降低，放大作用效果也隨著垂蕩頻率的增加而增大。

圖11 fh 不同,晃蕩載荷隨fr 變化情況Fig.11 Variation of sloshing load with fr for different fh

為了進(jìn)一步探究復(fù)合激勵(lì)下頻率對(duì)晃蕩壓力的影響，對(duì)fr=0.81 Hz 系列工況的晃蕩壓力時(shí)程曲線進(jìn)行傅里葉分析，對(duì)其中影響因素較大的頻率極其幅值進(jìn)行統(tǒng)計(jì)，如圖12 所示。

本文研究的復(fù)合激勵(lì)包含橫搖激勵(lì)與垂蕩激勵(lì)，對(duì)壓力時(shí)程曲線進(jìn)行傅里葉變換后可以發(fā)現(xiàn)，對(duì)晃蕩載荷影響最大的為n倍的fr（n為非零正整數(shù)）。如圖12(a)所示，隨著垂蕩激勵(lì)的不斷增大，fr的幅值整體呈下降趨勢，但仍大于圖11(b)中各頻率的幅值，說明此時(shí)橫搖頻率依然是影響晃蕩載荷的主要因素。除了nfr以外，nfh，n(fr+fh)，n(fr?fh)都是影響晃蕩載荷的重要頻率因素。如圖12(b)，3 種頻率的幅值隨著垂蕩激勵(lì)的增大而不斷增大，在欠共振階段幅值增速較慢，在超共振階段，頻率的幅值增速較快，其中增速最快的是fr+fh。

圖12 各主要頻率的幅值變化Fig.12 Amplitude variation of main frequencies

4 結(jié)語

本文基于VOF 法多相流模型模擬橫搖與垂蕩復(fù)合激勵(lì)下二維矩形液艙的晃蕩問題，研究分析了在25%H低裝載下自由液面處壓力載荷的時(shí)域與頻域特征，得到以下結(jié)論：

1）當(dāng)fr≠fh時(shí)，每個(gè)運(yùn)動(dòng)周期內(nèi)的壓力載荷變化會(huì)出現(xiàn)多種類型的壓力載荷變化，多個(gè)運(yùn)動(dòng)周期組合形成的更大周期；當(dāng)fr=fh時(shí)，壓力載荷曲線與單純橫搖相同，每個(gè)運(yùn)動(dòng)周期內(nèi)的壓力載荷變化穩(wěn)定。

2）橫搖與垂蕩復(fù)合激勵(lì)中，橫搖激勵(lì)頻率起主導(dǎo)作用，當(dāng)fr=fh時(shí)，垂蕩激勵(lì)會(huì)減小晃蕩壓力，fh>fr時(shí)垂蕩激勵(lì)放大晃蕩壓力的作用會(huì)更加明顯，但效果沒有改變橫搖激勵(lì)頻率顯著。

3）隨著fh的不斷增大，fr對(duì)晃蕩壓力的影響會(huì)逐漸降低但仍為主要影響因素，fh，fr+fh，fr?fh對(duì)壓力載荷的影響不斷增加，多種頻率的影響使壓力時(shí)程曲線呈現(xiàn)多種周期形式，在這些周期內(nèi)壓力雙峰大小交替性變化。

本文僅研究25%裝載率矩形液艙在橫搖-垂蕩復(fù)合激勵(lì)下晃蕩壓力變化，對(duì)于中高裝載率與其他多種復(fù)合激勵(lì)還需要做進(jìn)一步研究，從而全面探索液艙晃蕩多種運(yùn)動(dòng)下的晃蕩壓力變化規(guī)律，為工程設(shè)計(jì)提供參考與建議。