張會霞，李 師，鄒昶方

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

0 引言

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

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

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

1 無關(guān)性分析

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

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

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

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

2 橫搖激勵下晃蕩壓力特點

液體固有頻率計算公式：

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

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

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

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

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

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

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

圖4 fr=0.81 Hz 點P 處壓力時程曲線Fig.4 Pressure time history curve at point P with fr=0.81 Hz

3 復合激勵下液艙晃蕩壓力特點

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

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

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

表2 復合激勵工況Tab.2 Compound excitation working condition

3.1 數(shù)值模型驗證

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

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

圖6 Constantin L 垂蕩激勵實驗自由液面變化Fig.6 Free surface change in Constantin L heaving excitation experiment

3.2 復合激勵下晃蕩壓力時域特征

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

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

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

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

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

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

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

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

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

3.3 復合激勵下晃蕩壓力頻域特征

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

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

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

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

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

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

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

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

4 結(jié)語

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

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

2）橫搖與垂蕩復合激勵中，橫搖激勵頻率起主導作用，當fr=fh時，垂蕩激勵會減小晃蕩壓力，fh>fr時垂蕩激勵放大晃蕩壓力的作用會更加明顯，但效果沒有改變橫搖激勵頻率顯著。

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

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