畸形波作用下浮體的載荷與響應(yīng)研究新方法

陳旭東，朱仁慶，楊 帆，紀(jì)仁偉

(江蘇科技大學(xué) 船舶與海洋工程學(xué)院，江蘇 鎮(zhèn)江 212003)

0 引 言

畸形波（freak wave）是一種波高巨大、波峰異常尖瘦的不規(guī)則波[1–2]，具有極強(qiáng)的非線性。Klinting等[3]認(rèn)為，自然海況下，最大波高分別大于有效波高與前后相鄰波高的2倍，而小于其波峰高度65%的波浪可稱為畸形波。畸形波能量巨大，對(duì)海上平臺(tái)及船舶構(gòu)成很大的威脅，因而引起海洋工程界的特別重視。近些年來(lái)國(guó)外學(xué)者通過(guò)模型試驗(yàn)與基于勢(shì)流理論的時(shí)域模擬方法對(duì)海洋結(jié)構(gòu)物如半潛平臺(tái)、FPSO遭遇畸形波過(guò)程進(jìn)行研究，得到了平臺(tái)在畸形波作用下運(yùn)動(dòng)響應(yīng)的基本規(guī)律[4–6]。Rudman等[7]基于SPH方法，對(duì)畸形波砰擊浪向角和預(yù)張力對(duì)張力腿平臺(tái)運(yùn)動(dòng)的影響進(jìn)行深入研究，并且對(duì)每一條張力腿的最大張力給出合理預(yù)測(cè)。Zhao等[8]基于自主研發(fā)的CIP模型開(kāi)展了極端波浪條件下浮體兩自由度響應(yīng)的模擬研究，初步驗(yàn)證了CIP方法建立的波浪水槽對(duì)該類問(wèn)題的適用性。

本文以Ansys Workbench為計(jì)算平臺(tái)，重點(diǎn)研究浮體的中橫剖面切片模型在畸形波作用下的載荷與運(yùn)動(dòng)響應(yīng)問(wèn)題，這樣浮體結(jié)構(gòu)的六自由度運(yùn)動(dòng)就可以簡(jiǎn)化為三自由度運(yùn)動(dòng)，即橫搖、垂蕩和橫蕩。載荷問(wèn)題主要包括上浪水位和砰擊壓強(qiáng)，本文提出“單元貼片法”，便捷有效地解決了追蹤浮式結(jié)構(gòu)物表面“動(dòng)點(diǎn)”的難題。

1 理論基礎(chǔ)與控制方程

1.1 理論基礎(chǔ)

按照Longuet-Higgins[9]提出的經(jīng)典海浪模型，固定點(diǎn)的波面表達(dá)式為：

式中：為組成波（Component waves）個(gè)數(shù)；和分別為第i 個(gè)組成波的波數(shù)和角頻率；為第i個(gè)組成波的相位；且滿足個(gè)余弦波在特定位置， 特定時(shí)刻的峰值疊加的可表示為[10]：

為第i個(gè)組成波波幅，其表達(dá)式為：

式中：為波浪頻譜，是頻率的函數(shù)。本文采用JONSWAP譜描述隨機(jī)海浪特性，其特征參數(shù)有：有義波高，譜峰周期，關(guān)于該譜的理論描述詳見(jiàn)文獻(xiàn)[10]。聯(lián)立式(1)～式(3)可得不規(guī)則波波面表達(dá)式為：

1.2 控制方程

連續(xù)性方程：

N-S方程：

式中：表示，方向；表示，方向的速度；為流體的密度；為動(dòng)力粘性系數(shù)；為重力加速度；為壓力。

1.3 自由液面追蹤原理

自由液面的追蹤采用VOF法，具體來(lái)說(shuō)，采用流體體積函數(shù)法來(lái)標(biāo)記自由液面，其表達(dá)式為：

式中：定義為流體單元內(nèi)第相流體所占體積與該單元的體積比。若，則表示該單元內(nèi)全部為第相流體；若，則表示該單元內(nèi)沒(méi)有第相流體；若，則表示該單元為交界面單元。文中表示空氣相，表示水相。

2 單元貼片法

本文在建立波物相互作用的二維數(shù)值模型時(shí)，將數(shù)值波浪水池與浮體共同定義一個(gè)“薄片”厚度（基于切片理論），使流場(chǎng)與浮體生成真實(shí)的交界面。由于波浪在寬度上與浮體一致，波物運(yùn)動(dòng)在任意縱剖面上的投影情況完全相同，因此可以看作是二維運(yùn)動(dòng)。這樣，在寬度方向上，可以用一個(gè)緊貼運(yùn)動(dòng)浮體的面元（簡(jiǎn)稱動(dòng)面元）來(lái)代替該面元的中心點(diǎn)，如圖1（a）所示，同時(shí)在流場(chǎng)中同一位置處挖出這樣一個(gè)面元，如圖1（b）所示，然后在耦合求解器中將二者以流固交界面的形式連接在一起，則通過(guò)追蹤流場(chǎng)中的該動(dòng)面元上的信息即可知道運(yùn)動(dòng)著浮體上的某個(gè)點(diǎn)上的信息。該動(dòng)面元功能類似于物理試驗(yàn)中測(cè)量結(jié)構(gòu)應(yīng)變的“應(yīng)變片”，故而本文將該面元稱作“單元貼片”。

圖 1 單元貼片法Fig. 1 Patch element method

3 算 例

3.1 數(shù)值模型

數(shù)值水池尺寸如下：長(zhǎng)40 m，高3 m，水深2 m，尾部10 m為消波區(qū)域，如圖2所示。結(jié)構(gòu)物初始時(shí)刻處于正浮狀態(tài)，與造波板的水平距離為，吃水為，干舷為。

圖 2 數(shù)值模型示意圖Fig. 2 Schematic diagram of the numerical model

流場(chǎng)網(wǎng)格劃分采用結(jié)構(gòu)網(wǎng)格和非結(jié)構(gòu)網(wǎng)格相結(jié)合的方法，浮體周圍采用三角形非結(jié)構(gòu)網(wǎng)格，通過(guò)動(dòng)網(wǎng)格技術(shù)中的擴(kuò)散光順?lè)▉?lái)更新網(wǎng)格，以精確模擬浮體的運(yùn)動(dòng)。外圍區(qū)域采用規(guī)則的四邊形結(jié)構(gòu)網(wǎng)格進(jìn)行劃分，如圖2所示，流場(chǎng)網(wǎng)格數(shù)量為118 256。

3.2 邊界條件及求解設(shè)置

水池頂部為pressure-inlet壓力入口邊界條件，底部及尾壁均為wall無(wú)滑移固壁條件，左端為piston運(yùn)動(dòng)邊界，用以模擬平推式造波板的運(yùn)動(dòng)，其運(yùn)動(dòng)規(guī)律為：

式中：為第i個(gè)組成波與造波板運(yùn)動(dòng)速度之間的水力傳遞函數(shù)，其表達(dá)式為

將式（4）代入式（9）得到推板造畸形波的速度公式：

在Fluent中采用非定常分離隱式求解器求解，壓力方程選用加權(quán)體積力格式（Body Force Weighted），壓力速度耦合方式采用PISO算法，體積分?jǐn)?shù)法為幾何重構(gòu)（Geo-Reconstruct），壓力參考值為一個(gè)標(biāo)準(zhǔn)大氣壓。流場(chǎng)采用 RNG k-e模式，動(dòng)量方程中的瞬態(tài)項(xiàng)采用一階迎風(fēng)格式，時(shí)間步長(zhǎng)為0.01 s。

3.3 計(jì)算分析

3.3.1 畸形波的生成

采用JONSWAP譜[10]描述隨機(jī)海浪特性，譜峰周期= 1.8 s，有義波高=0.012 7 m，得到譜形曲線如圖4所示。

圖 3 浮體周圍流場(chǎng)網(wǎng)格分布Fig. 3 Mesh distribution of flow field around floating body

圖 4 目標(biāo)譜Fig. 4 Shape of target wave spectrum

設(shè)在目標(biāo)譜的譜峰頻率兩側(cè)分別略去能量峰值的0.1%和5%，由此可得頻譜范圍=2.318～9.891 rad/s。取組成波數(shù)目，各組成波波幅總和=0.09 m，設(shè)計(jì)聚焦時(shí)間與聚焦位置分別為=15 s，=10 m，計(jì)算結(jié)果如圖5所示。

3.3.2 畸形波作用下不同橫截面浮體的運(yùn)動(dòng)響應(yīng)分析

本節(jié)針對(duì)不同橫截面浮體在畸形波中的運(yùn)動(dòng)響應(yīng)進(jìn)行計(jì)算分析，圖6（a）為無(wú)甲板室浮體，6（b）為有甲板室浮體，浮體初始位置=10 m。

圖 5 畸形波Fig. 5 Freak wave

圖 6 不同橫截面的浮體Fig. 6 Floating bodies with different cross section

圖7給出了2種工況下浮體運(yùn)動(dòng)響應(yīng)時(shí)歷。由圖可知，兩者的垂蕩運(yùn)動(dòng)變化在遭遇畸形波時(shí)刻及附近比較一致。甲板上浪水體沖擊甲板室會(huì)造成浮體的橫蕩和橫搖運(yùn)動(dòng)明顯不同于無(wú)甲板室時(shí)的情況：有甲板室的浮體在遭遇畸形波之前甚至大幅度反向漂移，畸形波襲擊過(guò)后，才逐漸向水池尾端漂移，并于t=25 s時(shí)達(dá)到與無(wú)甲板室浮體相同的橫蕩值；對(duì)于橫搖運(yùn)動(dòng)，在遭遇畸形波之前，2種工況的橫搖角變化幾乎一致，之后約2 s內(nèi)，無(wú)甲板室浮體的正向橫搖角幅值大于有甲板室浮體的正向橫搖角幅值，反向橫搖角幅值則呈現(xiàn)相反的情況，約t=17 s后，有甲板室的橫搖角幅值無(wú)論正反向均大于無(wú)甲板室的橫搖角幅值，并且在受畸形波襲擊后，較長(zhǎng)時(shí)間處于劇烈的橫搖運(yùn)動(dòng)中，不容易恢復(fù)正浮狀態(tài)，既體現(xiàn)了很強(qiáng)的波物非線性作用的特征，又毫無(wú)規(guī)律可言。

圖8給出了不同時(shí)刻上述2種工況下畸形波與浮體相互作用時(shí)的流場(chǎng)形態(tài)、浮體運(yùn)動(dòng)姿態(tài)及浮體結(jié)構(gòu)馮米斯應(yīng)力云圖。無(wú)甲板室浮體在遭遇畸形波襲擊后，上浪水體會(huì)很快經(jīng)甲板從另一側(cè)流出，對(duì)浮體背浪一側(cè)的流場(chǎng)造成擾動(dòng)，從而進(jìn)一步影響了浮體的運(yùn)動(dòng)狀態(tài)；而有甲板室浮體在遭遇畸形波襲擊后，由于甲板室壁的阻礙作用，上浪水體并未沖過(guò)甲板室，而是瞬間幾乎全部被擋回，對(duì)浮體迎浪一側(cè)的流場(chǎng)形成強(qiáng)烈干擾，甚至產(chǎn)生了氣泡和漩渦，這對(duì)浮體結(jié)構(gòu)大為不利。由于浮體結(jié)構(gòu)剛性足夠，因而可以承受畸形波浪的砰擊而幾乎不會(huì)產(chǎn)生變形。對(duì)于有甲板室的浮體而言，受到上浪水體的的猛烈砰擊后，甲板室與甲板連接處會(huì)發(fā)生應(yīng)力集中，極易造成結(jié)構(gòu)破壞，在設(shè)計(jì)初期就應(yīng)注意采取結(jié)構(gòu)加強(qiáng)措施。

圖 7 不同橫截面的浮體的運(yùn)動(dòng)響應(yīng)時(shí)歷Fig. 7 Time history of the motion response of floating bodies with different cross section

3.3.3 畸形波作用下抑制浮體漂移對(duì)甲板上浪的影響

為了探討畸形波作用下抑制浮體漂移對(duì)甲板上浪的影響，對(duì)3.2.2節(jié)中浮體（圖6（b））的橫蕩運(yùn)動(dòng)作了彈簧力抑制處理，使其在方向幾乎無(wú)漂移，即只保留了垂蕩與橫搖2個(gè)自由度。本節(jié)主要研究的是水動(dòng)力和運(yùn)動(dòng)響應(yīng)，因此在浮體迎浪一側(cè)共設(shè)置了4個(gè)監(jiān)測(cè)點(diǎn)，，，，其中，，用于監(jiān)測(cè)甲板上浪水位高度，間距為0.1 m；則監(jiān)測(cè)上浪水體對(duì)甲板室壁的砰擊壓強(qiáng)，與甲板垂直距離為0.01 m。浮體初始位置=10 m。

圖 8 畸形波與不同橫截面的浮體相互作用的云圖Fig. 8 Nephogram of interaction between freak wave and floating bodieswith different cross section

圖9為浮體在自由浮動(dòng)和限制漂移2種情況下的水動(dòng)力監(jiān)測(cè)點(diǎn)數(shù)值對(duì)比。從9（a）～9（c）中可以發(fā)現(xiàn)，在同一時(shí)刻，在遭遇畸形波的前后一段時(shí)間內(nèi)，P1，P2，P3處的上浪水位依次遞增，大約在第21 s后，P1，P2，P3處的上浪水位又依次遞減，由此可以推斷出甲板上浪的整個(gè)過(guò)程：伴隨著畸形波浪到達(dá)，水體涌上甲板，當(dāng)沖到甲板室壁時(shí)，水體順著甲板室壁面迅速向上爬升，產(chǎn)生瞬時(shí)飛濺現(xiàn)象，接著水體翻卷回流。當(dāng)畸形波列繞過(guò)浮體后，浮體搖蕩仍在繼續(xù)，且仍有少量水體涌上甲板，但其動(dòng)能很大程度上減小，因此在受到甲板室阻擋時(shí)幾乎沒(méi)有發(fā)生飛濺現(xiàn)象，因而此時(shí)上浪水位變化也較為緩和。圖9（d）中P0處監(jiān)測(cè)到的上浪水體對(duì)甲板室壁的砰擊壓強(qiáng)時(shí)歷證實(shí)了上述推斷的合理性。

自由浮動(dòng)情況下的甲板上浪水位整體上低于抑制漂移情況下甲板上浪水位，在14～19 s這段時(shí)間內(nèi)顯得尤為明顯，對(duì)于P1和P2，僅在約23 s時(shí)，抑制漂移情況下的上浪水位被自由浮動(dòng)情況反超，對(duì)于P3，這種水位反超現(xiàn)象提前約3 s。這主要是因?yàn)橐种聘◇w漂移后，浮體幾乎無(wú)法順浪前移，快速來(lái)襲的波浪遇到浮體壁面阻擋瞬時(shí)發(fā)生瞬時(shí)飛濺，后續(xù)的波浪在此瞬時(shí)堆積，一定程度上提升了涌上甲板的水體高度。而自由浮動(dòng)情況下，由于浮體具有順波浪方向的動(dòng)能，這在一定程度上緩解了來(lái)襲波浪的沖擊，因而涌上甲板的水體高度也相應(yīng)地降低。2種工況下浮體的運(yùn)動(dòng)響應(yīng)時(shí)歷如圖10所示。

圖 9 自由運(yùn)動(dòng)與抑制漂移兩種情況下的監(jiān)測(cè)點(diǎn)數(shù)值對(duì)比Fig. 9 Comparison between free-floating and x-fixed on monitoring point data

由圖10可知，自由浮動(dòng)和限制漂移2種約束情況下，浮體的垂蕩在畸形波浪來(lái)襲之前幾乎一致，在波浪聚焦開(kāi)始的一段時(shí)間內(nèi)出現(xiàn)偏差，且自由浮動(dòng)情況下的垂蕩值大于抑制漂移情況下的垂蕩值，約第20 s過(guò)后2種情況下的垂蕩值又持平。抑制漂移后，浮體的橫搖運(yùn)動(dòng)比自由浮動(dòng)時(shí)的情況更加劇烈，畸形波來(lái)襲前后，浮體從右傾約6°的極值短時(shí)間內(nèi)搖擺到左傾約12°的極值，擺幅可達(dá)約18°，而對(duì)于自由浮動(dòng)的浮體，相同時(shí)間范圍內(nèi)的擺幅約為13°。

圖 10 自由運(yùn)動(dòng)與抑制漂移兩種約束情況下的浮體運(yùn)動(dòng)響應(yīng)時(shí)歷Fig. 10 Time history of the motion response of floating bodies under different constraint coditions, respectively as free-floating and x-fixed

3.3.4 浮體在水池中的放置地點(diǎn)對(duì)甲板上浪的影響

為了探究浮體在水池中的放置地點(diǎn)對(duì)甲板上浪的影響，本節(jié)將3.2.2節(jié)中浮體（圖6（b））分別放置于距水池左端10 m、9.5 m和10 m處，且均抑制掉其橫蕩運(yùn)動(dòng)，使其在方向幾乎無(wú)漂移，同時(shí)只保留了垂蕩與橫搖2個(gè)自由度。

圖11給出了浮體上的水動(dòng)力監(jiān)測(cè)點(diǎn)數(shù)據(jù)時(shí)歷。比較同一放置地點(diǎn)、不同水位監(jiān)測(cè)點(diǎn)的數(shù)據(jù)可知，在設(shè)計(jì)波能聚焦時(shí)刻附近，離甲板邊緣相對(duì)較遠(yuǎn)的水位監(jiān)測(cè)點(diǎn)的上浪水位也相對(duì)較高，這主要是因?yàn)椴ɡ擞可霞装搴笤诩装迳隙虝r(shí)間堆積，遇到甲板室壁瞬間翻卷回流，一部分水體越過(guò)甲板室壁繼續(xù)向前流動(dòng)，甲板室前的上浪水體則通常呈現(xiàn)前低后高的近似梯形或三角形分布[11]；比較不同放置地點(diǎn)、同一水位監(jiān)測(cè)點(diǎn)的數(shù)據(jù)可知，在波能聚焦時(shí)刻附近，距造波板相對(duì)較近的同一水位監(jiān)測(cè)點(diǎn)的上浪水位相對(duì)較高。同時(shí)觀察P0砰擊壓強(qiáng)時(shí)歷可發(fā)現(xiàn)，浮體位于=9.5 m時(shí)，P0約在12.68 s時(shí)刻瞬間達(dá)到5 093.91 Pa的峰值，波浪動(dòng)能異常巨大，之后又分別在約14.2 s和15.80 s時(shí)刻出現(xiàn)了954.13 Pa和854.75 Pa的大峰值砰擊壓強(qiáng)，體現(xiàn)了波浪砰擊的隨機(jī)性與強(qiáng)非線性特點(diǎn)。相比之下，浮體在=10.5 m時(shí)的P0砰擊壓強(qiáng)相對(duì)最小，僅在約14.66 s時(shí)達(dá)到529.05 Pa的大峰值。3種工況下浮體的運(yùn)動(dòng)響應(yīng)時(shí)歷如圖12所示。

圖 11 浮體不同放置位置情況下的監(jiān)測(cè)點(diǎn)數(shù)值對(duì)比Fig. 11 Comparison among different positions of floating body on monitoring point data

圖 12 不同放置地點(diǎn)情況下的浮體運(yùn)動(dòng)響應(yīng)時(shí)歷Fig. 12 Time history of the motion response of floating bodies located at different positions

由圖12可知，浮體在=9.5 m處時(shí)正向橫搖幅值最大，隨著放置位置的后移，正向橫搖幅值先減小后增大，而浮體在=10.5 m處時(shí)的反向橫搖幅值可達(dá)到18°，但是都顯示出隨機(jī)性特點(diǎn)。這表明，從浮體穩(wěn)性考慮，放置于波能聚焦處未必是最危險(xiǎn)的狀態(tài)，偏離波能聚焦處的橫搖力矩大一些，因而橫搖運(yùn)動(dòng)也相對(duì)劇烈。

4 結(jié) 語(yǔ)

1）基于線性波疊加原理，采用能量聚焦方法生成了畸形波。

2）基于流固耦合理論與切片原理，對(duì)不同橫截面浮體在畸形波中的運(yùn)動(dòng)響應(yīng)進(jìn)行計(jì)算分析，結(jié)果表明，相較于無(wú)甲板室浮體而言，有甲板室浮體在提高浮體重心高度的同時(shí)，會(huì)顯著改變上浪水體的運(yùn)動(dòng)形態(tài)，對(duì)浮體所處的流場(chǎng)，特別是迎浪一側(cè)造成強(qiáng)烈干擾，進(jìn)而使得浮體運(yùn)動(dòng)更加非線性化。

3）提出“單元貼片”研究方法，分別對(duì)自由浮動(dòng)及抑制漂移2種情況下浮體與畸形波相互作用的載荷與響應(yīng)問(wèn)題進(jìn)行了對(duì)比研究，發(fā)現(xiàn)甲板上浪水體對(duì)迎浪側(cè)的流場(chǎng)形成強(qiáng)烈的擾動(dòng)，并不同程度地伴隨著氣泡和漩渦的產(chǎn)生，而背浪一側(cè)則顯得相對(duì)平靜一些，而上浪水體的劇烈砰擊使得上層建筑與甲板連接處的結(jié)構(gòu)馮米斯應(yīng)力值瞬間增加，容易使此處材料發(fā)生屈服失效甚至破壞，故應(yīng)注意采取適當(dāng)?shù)慕Y(jié)構(gòu)加強(qiáng)措施。同時(shí)，抑制浮體漂移后，相應(yīng)的上浪水位高度、上浪水體砰擊壓強(qiáng)及浮體的橫搖運(yùn)動(dòng)也增加，但垂蕩運(yùn)動(dòng)在波浪聚焦處開(kāi)始減小?；谏鲜龇椒?，對(duì)浮體在水池中不同位置時(shí)的甲板上浪情況進(jìn)行數(shù)值模擬，結(jié)果表明，浮體處于同一位置時(shí)，在設(shè)計(jì)波能聚焦時(shí)刻附近，離甲板邊緣相對(duì)較遠(yuǎn)的水位監(jiān)測(cè)點(diǎn)的上浪水位相對(duì)較高；浮體處于不同位置時(shí)，在波能聚焦時(shí)刻附近，距造波板相對(duì)較近的同一水位監(jiān)測(cè)點(diǎn)的上浪水位相對(duì)較高。

[1] 楊冠聲, 董艷秋, 陳學(xué)闖. 畸形波(freak wave)[J]. 海洋工程,2002, (4): 105–108.YANG Guan-sheng, DONG Yan-qiu, CHEN Xue-chuang.Freak wave[J]. The Ocean Engineering, 2002, (4): 105–108.

[2] KIMURA A, OHTA T. Probability of the freak wave appearance in a 3-dimensional sea condition[C]// Proceedings of the Coastal Engineering Conference, 1994, 23–28.

[3] KLINTING P, SAND S. Analysis of prototype freak waves.Conf. on nearshore hydrodynamics[C]// ASCE. 1987:618–632.

[4] GUNTHER F C, CHRISTIAN E S, KATJA S. Freak wave impact on semisubmersible time-domain analysisof motions and forces [C]// Proceeding of The Thirteenth (2003)International Offshore and Polar Engineerin Conference. 2003,449–456.

[5] GUNTHER F C, MARCO K. Influence of the bow shape on loads in high and steep waves [C]// Proceedings of the ASME 2011 30th International Conference on Ocean, Offshore and Arctic Engineering. 2010: OMAE 2010–2014.

[6] GUNTHER F C, SURESH R, MARCO K. Time domain comparison with experiments for ship motions and structure loads on a container ship in abnormal waves [C]// Proceedings of the ASME 2011 30th International Conference on Ocean,Offshore and Arctic Engineering. 2011: OMAE 2011–50316.

[7] RUDMAN M, CLCARY P W. Rogue wave impact on a tension leg platform: the effect of wave incidence angle and mooring line tension[J]. Ocean Engineering, 2013(61):123–138.

[8] ZHAO Xi-zeng, HU Chang-hong. Numerical and experimental study on a 2-D floating body under extreme wavc conditions[J]. Applied Ocean Research, 2012, 35: 1–13.

[9] LONGUET-HIGGINS M S, COKELET E D. The deformation of steep surface waves on water. I. A numerical method of computation[J]. Proceedings of the Royal Society of London.A. Mathematical and Physical Sciences, 1976, 350(1660):1–26.

[10] 趙艷. 強(qiáng)非線性波與海洋浮式結(jié)構(gòu)物的相互作用[D]. 鎮(zhèn)江：江蘇科技大學(xué), 2014.

[11] 劉利琴, 王賓, 沈文君. 甲板上浪船舶的橫搖運(yùn)動(dòng)[J]. 振動(dòng).測(cè)試與診斷, 2012(S1): 14–19, 146.LIU Li-qin, WANG Bin, SHEN Wen-jun. The rolling motion of the ship of green water[J]. Journal of Vibration,Measurement & Diagnosis, 2012(S1): 14–19, 146.