何廣華， 荊芃霖， 王威， 張婧文， 金瑞佳， 劉濤

(1.哈爾濱工業(yè)大學 機電工程學院，黑龍江 哈爾濱 150001； 2.哈爾濱工業(yè)大學(威海) 海洋工程學院，山東， 威海 264209； 3.山東船舶技術研究院，山東 威海 264209； 4.交通運輸部天津水運工程科學研究院， 天津 300450)

在陸地資源開發(fā)利用已趨于飽和的今天，海洋資源作為地球自然資源的主要組成部分得到了世界各國的青睞，海洋工程結構物的開發(fā)及應用的競爭日趨激烈。因海洋結構物及其組成部分之間多存在著窄縫，如：靠泊的船體與碼頭、并排的船與船、浮式生產儲油卸油裝置(floating production storage and offloading，F(xiàn)PSO)的月池及箱式超大型海洋結構物各模塊間隙等。在入射波浪的持續(xù)砰擊下，間隙內的波面抬高呈往復運動。在某些特定波頻率下，其振幅要遠大于入射波幅，對結構物的作用力也隨之顯著增大，且力的變化極為劇烈。這種現(xiàn)象被稱為窄縫共振。

窄縫共振關系到海洋工程結構物的正常生產作業(yè)安全，對結構使用壽命也存在影響。研究窄縫共振問題，不僅要深入分析現(xiàn)象，更要明晰其產生的機理，對于多浮體海洋結構物及FPSO等系統(tǒng)的設計具有積極指導意義，同時也對規(guī)避海域作業(yè)時窄縫共振的不利影響具有重要的工程應用價值。

線性理論應用到求解兩箱體間的窄縫共振問題[1-8]，學者們基于線性理論進行了變參數(shù)研究和共振模態(tài)研究?；诎虢馕龇椒ê蛣萘骼碚摰姆椒m可高效地對雙箱共振問題進行求解，但因忽略了流體的粘性作用，共振區(qū)間內預報雙箱間窄縫處的波面抬高常常偏差較大，且多為固定雙箱，未考慮雙箱在入射波作用下的運動響應。

隨著計算機技術的不斷發(fā)展計算流體力學(computational fluid dynamics，CFD)方法更多的應用到窄縫共振問題的研究中[9-13]。CFD方法考慮流體粘性，可更直觀地展示流場細節(jié)，便于對精細流場進行分析。Saitoh等[14]做了不同窄縫寬度與浮體吃水的固定雙箱模型實驗，對不同入射波浪作用下窄縫共振波高進行了研究，發(fā)現(xiàn)窄縫內流體的共振波數(shù)、共振波高受箱體吃水深度及箱體間距影響。Iwata等[15]隨后采用與Saitoh相似的實驗方法，研究了3個毗鄰排放箱體間隙內水體的共振現(xiàn)象。物理模型實驗的方法更加直觀，與數(shù)值模擬相結合的方法能更加深入地對窄縫共振現(xiàn)象進行分析。

由上述文獻研究可見，對于雙箱共振問題的研究大多基于雙箱固定不動的假定[6-15]，與多浮體海洋結構物實際運動響應的情況相差較大，仍需進一步研究。本文對自由垂蕩雙箱在不同頻率入射波激勵下的運動響應及雙箱間窄縫處水體共振問題，采用實驗和數(shù)值相結合的方法研究，深入分析了窄縫水體共振機理，給出了該類物理實驗改進建議。自由垂蕩雙箱的窄縫共振問題新穎，涉及入射波砰擊、箱體自由垂蕩響應、窄縫處水體活塞運動之間的耦合作用，現(xiàn)象復雜且影響因素多，本文對于進一步研究多浮體結構物垂蕩及窄縫水體共振問題具有積極作用。

1 實驗模型及設置

本模型實驗水槽長68 m，寬1.0 m，深1.0 m，水槽前后兩端均設有消波裝置。兩箱體尺寸完全一致，如圖1所示，箱體寬度B=0.5 m, 雙箱間窄縫寬度Bg=0.05 m, 水深h=0.5 m, 箱體吃水d=0.153 m。

實驗所用浪高儀為TS-DWG型電容式數(shù)字信號浪高儀，絕對誤差≤0.1 mm，浪高儀放置于雙箱間窄縫的中心，如圖1所示。浪高儀的采樣頻率設為50 Hz，采用電磁式非接觸六分量傳感器LIBERTY 2408系統(tǒng)來測定箱體位移。實驗方箱由壁厚10 mm的PVC硬板制成，方箱沿水槽方向寬度為967 mm，箱體剛度及水密性滿足實驗要求。水槽內壁粘貼有機玻璃板，兩箱體沿固定在玻璃板上的導軌做垂蕩運動，該導軌由耐腐蝕、低噪音的SGR10內置雙軸心滾輪直線結構組成如圖1所示。

圖1 實驗設置Fig.1 Schematic of experimental setup

實驗中對箱體位移和窄縫處的波面抬高進行實時監(jiān)測。圖2為實驗中的經典工況kh=2.40(共振點)時窄縫處波面抬高和箱體A垂蕩運動的時歷曲線，k為波數(shù)、H0為入射波高、hg為窄縫處波面抬高、yA為箱體A的垂向位移。波面抬高周期性較好，波面穩(wěn)定，波峰處較為平緩。箱體A垂蕩運動呈周期性變化，運動達到相對穩(wěn)定狀態(tài)。圖3為共振頻率kh=2.40條件下，一個波周期內實驗的照片，圖中長實線為靜水面位置，短虛線為窄縫內波面，窄縫內水體的運動幅值大、變化劇烈。

2 數(shù)值模型

2.1 有限體積法

STAR-CCM+使用有限體積法求解雷諾平均納維斯托克斯(Reynolds-averaged Navier Stokes，RANS)方程，將流體域離散化為有限個控制體積，采用壓力耦合方程組的半隱式方法(semi-implicit method for pressure-linked equations，SIMPLE)實現(xiàn)壓力和速度之間的隱式求解。積分形式的控制方程為：

質量守恒方程：

(1)

動量守恒方程：

(2)

空間守恒方程：

(3)

式中：v是流體速度矢量；vb表示控制體積表面的速度；n是控制體積曲面的單位法向量；控制體積的面積為S；體積為V；流體密度為ρ；時間為t；壓力為p；流體應力以張量T的形式表示；I是單位張量；物體力的矢量表示為b。

通過流體體積 (volume of fluid，VOF)多相模型引入體積分數(shù)來計算出多相流中的自由液面位置；使用高分辨率交界面捕捉(high resolution interface capturing，HRIC)離散格式來精確捕捉自由液面的流體變化。

圖2 實驗時歷曲線Fig.2 Temporal data of the experiment

圖3 一個周期內實驗照片F(xiàn)ig.3 Photos of this experiment in one period

2.2 DFBI模型

動態(tài)流體相互作用(dynamic fluid body interaction，DFBI)模型可用來求解剛性體的流固動態(tài)耦合，使用6自由度(6-DOF)求解器計算結構物的6自由度運動，6-DOF求解器采用二階精度的梯形方案，可以實現(xiàn)6個自由度同時求解。在總體坐標系下，物體的平動方程為：

(4)

式中：m是物體質量；f是作用在物體上的合力；vc是質心的速度。

在隨體坐標系下物體的轉動方程：

(5)

式中：M是慣性矩張量；ω是剛體的角速度；Mn是作用在物體上的合力矩。

轉動慣量的張量展開式為：

(6)

2.3 計算域布置和邊界條件

計算域的邊界條件設置如圖4所示，箱體A和B尺寸與實驗箱體完全一致，流域速度入口邊界通過場函數(shù)設定自由液面位置和一階波的速度，壓力出口邊界給出了自由面條件，出口處的壓力被設定為波浪的靜水壓力。在計算域的入口和出口參照實驗條件使用力消波(wave forcing)進行消波處理。

圖4 計算域和邊界條件Fig.4 Computational domain and boundary condition

2.4 網(wǎng)格劃分

網(wǎng)格劃分采用了切割體、拉伸體、網(wǎng)格重疊等技術組合，保證網(wǎng)格質量。具體網(wǎng)格結構如圖5所示，箱體周圍網(wǎng)格逐步加密(圖5左圖)以精確捕捉作用在箱體上的波浪力和箱體運動響應，遠離箱體區(qū)域使用漸疏網(wǎng)格，減弱遠方邊界處的波浪反射并減少網(wǎng)格數(shù)量(圖5右圖)，同時自由液面附近網(wǎng)格適當加密捕捉波面分布。在箱體周圍的網(wǎng)格為隨體運動網(wǎng)格，該網(wǎng)格在背景網(wǎng)格上重疊平動，計算過程中各網(wǎng)格形狀無變化，精確預報箱體的運動響應，計算效率較高。

圖5 箱體周圍網(wǎng)格(左)及漸疏網(wǎng)格(右)Fig.5 Grid around the bodies (left) and gradually thinning grid (right)

3 結果與討論

3.1 單箱垂蕩運動收斂性驗證

本文基于該數(shù)值水池對單個箱體在波浪中的運動響應進行數(shù)值模擬，數(shù)值模型的收斂性驗證見文獻[16]。本節(jié)進一步研究了無因次波數(shù)kb在0.4～1.2的典型工況(b=0.5B)。為對比分析，本文數(shù)值波浪水池和箱體參數(shù)與Rodríguez的實驗[17]保持一致，所得箱體的響應幅值算子(RAO)如圖6所示：kb=0.4,0.5,0.8,1.0,1.2工況與實驗值對比吻合極好，相對誤差僅在5%左右；kb=0.6,0.7這2個工況的相對誤差稍有增大；數(shù)值方法與實驗方法對于共振頻率的預報一致。綜合考慮，可認為本文所用數(shù)值模型在求解箱體垂蕩運動響應問題上有較好的可靠性與準確性。

圖6 單箱RAO對比驗證Fig.6 Comparison verification of RAO of single box

3.2 雙箱自由垂蕩運動數(shù)值結果與實驗結果對比

針對固定雙箱問題，文獻[15]中給出了計算共振頻率的詳細推導過程和計算公式：

(7)

式中：ωg為共振頻率；g為重力加速度。由式(7)結合色散公式可計算得出，本文中固定雙箱發(fā)生共振時的無因次波數(shù)kh為2.27，文獻[14]中通過物理模型實驗得到的共振波數(shù)為kh=2.154，二者之間相對誤差為5.37%。

因2個箱體均考慮了垂向的自由度，雙箱共振問題變得更加復雜；不僅存在窄縫處水體共振，還涉及窄縫處水體和箱體A、B垂蕩響應之間相互耦合問題。本文采用實驗方法進行研究；重復性良好，實驗結果取平均值。實驗中的入射波頻率kh包括 2.00、2.20、2.30、2.35、2.40、2.45、2.50、2.60和2.80，實驗結果見圖7、8。數(shù)值模擬工況與實驗保持一致，入射波高均為H0=0.024 m，在kh在2.0～2.8內均勻取點，并在共振點附近進行加密取點來精確尋找共振頻率。圖7為雙箱間窄縫處的波面抬高，將本文實驗及數(shù)值結果與固定雙箱實驗[14]相對比分析，發(fā)現(xiàn)：箱體可自由垂蕩運動時，窄縫處波面抬高的共振點、共振區(qū)間和共振點處的幅值均較雙箱固定條件下有較大的差異；固定雙箱共振點幅值要遠大于自由垂蕩雙箱，為實驗結果的1.85倍，數(shù)值結果的2.17倍；自由垂蕩雙箱的實驗和數(shù)值結果雖存在一定的差異，但2條曲線的變化趨勢一致，兩者對共振頻率的預報基本一致。雙箱間窄縫處的水體有2種共振形式：活塞式和晃蕩式。雙箱間窄縫寬度較大時共振為晃蕩式，寬度較小時為活塞式[18]；本文雙箱的窄縫寬度較小，水體共振形式表現(xiàn)為活塞式共振。

圖7 窄縫處波面抬高Fig.7 Wave elevation at narrow gap

自由垂蕩雙箱在共振時的窄縫波面抬高均小于固定雙箱情況，因波面抬高是影響箱體所受波浪力的重要因素，因此，若在多模塊海洋結構物的共振問題研究中，若不考慮結構物的垂蕩響應，會在錯誤預報共振區(qū)間及共振頻率的同時，高估雙箱共振現(xiàn)象的危害。

圖8為箱體A和箱體B在波浪砰擊下的垂蕩運動的RAO，數(shù)值結果與實驗結果二者變化規(guī)律相同，隨kh增加，二者在數(shù)值上也趨于一致，尤其kh在2.5～2.8內吻合較好。表1給出了在差異較大的區(qū)間內，不同kh條件下實驗與數(shù)值結果的比值，易知該區(qū)間內二者比值十分接近，故而認為實驗中存在一個固定比例的阻尼力或能量耗散，該現(xiàn)象將在下一小節(jié)中進行詳細分析。

圖8 箱體A和箱體B的RAOFig.8 RAO of BOX-A and BOX-B

表1 實驗值與數(shù)值結果之比Table 1 Ratio of experimental data and numerical data

圖9為箱體A和箱體B分別在水平和垂直方向的無量綱波浪力，波浪力的無量綱化為：

式中：ρ為流體密度；e為計算域厚度；H0為入射波高；FAx、FAy、FBx、FBy分別為箱體A和箱體B所受水平波浪力和垂直波浪力。

由圖9可知：箱體A在垂向和水平向所受波浪力均分別大于箱體B；箱體B所受垂直波浪力隨kh的增大而趨于0，即入射波長越小，箱體B在背浪側所受的波浪砰擊影響越小；當波長增大時，箱體A和箱體B在水平方向所受波浪力趨于一定值；雙箱水平波浪力和垂直波浪力的變化趨勢分別一致。

圖9 雙箱所受波浪力Fig.9 Wave force on twin boxes

為對自由垂蕩運動雙箱的水動力特性進一步分析，圖10給出了雙箱各參數(shù)無因次化后的時歷曲線，其中yA為箱體A垂向位移，yB為為箱體B垂向位移。由圖10可知：箱體B所受水平波浪力與窄縫處波面抬高的時歷曲線相位相同，窄縫處的波面抬高是影響箱體B水平波浪力的主要因素；箱體A所受水平波浪力則由于箱體自由垂蕩運動與窄縫處的波面抬高的耦合影響，窄縫處的波面抬高不再是其主要影響因素，這點與固定雙箱條件下有著明顯的不同[12]；箱體A和箱體B在垂直方向所受波浪力則主要受其本身的垂向位移影響，其時歷曲線相位相反。

圖10 共振點處各參數(shù)無因次量Fig.10 Dimensionless parameters at resonance point

4 誤差分析

由3.2節(jié)可知，本實驗與數(shù)值模擬的結果之間存在一定誤差，而數(shù)值模型已經過系統(tǒng)的有效性與可靠性驗證[16]，因此基于數(shù)值結果對模型實驗進行誤差分析。因實驗中對每個工況都進行了多次實驗，且重復性良好，因此排除偶然結果的影響。由表1可知，在一定波數(shù)下，實驗值與數(shù)值解之比較為穩(wěn)定，而箱體A和箱體B垂蕩RAO要小于數(shù)值結果，窄縫處的波面抬高的實驗值要大于數(shù)值結果。綜上可知：共振點處雙箱運動幅值越大，窄縫處的波面抬高越低，符合能量守恒原理。

數(shù)值模擬與實驗結果存在偏差有諸多因素，如：導軌存在對流場的影響、導軌安裝(摩擦與偏向)的影響、箱體與水池邊壁間縫隙的影響等。下面對幾個關鍵因素進行分析與討論。

4.1 導軌安裝誤差的影響

安裝誤差見示意圖11，是夸張化的導軌安裝誤差示意圖。如圖11所示，安裝誤差可能使得箱體在導軌中運動時存在較小阻尼，而且導軌與滑輪本身又存在些許的摩擦。由于水槽寬度僅為1.0 m，人員的作業(yè)環(huán)境和作業(yè)條件較差，在導軌安裝過程中容易產生偏差；在導軌偏向上的任何一絲偏差，對結果會產生較大的影響。

4.2 箱體-水槽的間隙影響

在安裝箱體時，為保證水槽側壁的安全，在箱體與水槽之間留有33 mm的縫隙，如圖12所示，圖12為實驗裝置的俯視圖。在波浪的砰擊下，會有一部分的波浪從該縫隙中穿過。為探究該縫隙的存在對于實驗的具體影響，再次通過數(shù)值方法1∶1模擬了3個典型工況(kh=2.20、2.40和2.60)下，窄縫存在時的雙箱算例，結果如圖13、14所示。由于實驗中池壁與箱體的間距較小，此時相當于把水槽沿長度方向分成了3個分水槽，3個分水槽在結構物周圍有一定的相互干涉作用，但總體來說：箱體-水槽的間隙對計算結果影響較小。

圖12 箱體與水池側壁相對位置示意圖(俯視圖)Fig.12 Signal of the relative position of box and wall of wave tank (vertical view)

圖13 箱體與水槽壁間隙對垂蕩RAO的影響Fig.13 Effect of the slit between the boxes and the wave tank wall on heave RAO of boxes

4.3 導軌存在對流場的影響

導軌會引起其周邊流場的劇烈變化，從而影響雙箱的運動響應。本節(jié)1∶1還原了有導軌時的三維雙箱共振算例，仍然選取典型工況kh=2.20，2.40，2.60，其結果如圖15所示。由圖15可知，導軌的存在對雙箱共振特性的影響十分明顯，尤其是窄縫處的波面抬高和箱體B垂蕩RAO，有導軌時窄縫處波面抬高的RAO與箱體B垂蕩RAO均介于實驗值與數(shù)值結果之間，可見導軌的存在對流場的干涉程度較大，對垂蕩雙箱的共振響應影響較大。

圖14 箱體與水槽壁間隙對窄縫處波面抬高的影響Fig.14 Effect of the slit between the boxes and the wave tank wall on wave elevation at narrow gap

圖15 導軌對垂蕩和窄縫處波面抬高的影響Fig.15 Effect of the guide on heave RAO of boxes and wave elevation at narrow gap

5 結論

1) 本文建立的二維粘流數(shù)值波浪水池可正確求解波浪中單箱和雙箱的垂蕩運動問題，且具有良好的穩(wěn)定性與準確性。

2) 在研究雙箱共振問題時，如果忽略兩毗鄰結構的垂蕩影響，假設為固定不動條件，會在預報共振區(qū)間與共振頻率時出現(xiàn)較大偏差。

3) 考慮箱體的垂蕩運動后，各箱體在窄縫共振發(fā)生時的水動力會發(fā)生變化。共振點處，箱體A的水平波浪力受到窄縫處波面抬高和箱體垂蕩運動的耦合影響；箱體B的水平波浪力受窄縫處的波面抬高影響；箱體A所受力的幅值都大于箱體B。

4) 本文所用實驗方法能準確預報可垂蕩雙箱的窄縫共振頻率；有利于自由垂蕩雙箱在不同入射波頻率下窄縫處波面抬高、箱體垂蕩位移等參數(shù)的變化趨勢，對于研究箱間水體的窄縫共振問題具有積極意義。