金瑞佳，滕 斌，呂 林

(大連理工大學(xué) 海岸和近海工程國家重點實驗室，遼寧 大連 116024)

應(yīng)用兩次展開法求解系泊柱體慢漂運動方程

金瑞佳，滕 斌，呂 林

(大連理工大學(xué) 海岸和近海工程國家重點實驗室，遼寧 大連 116024)

在深海中系泊的海洋平臺，如Spar平臺，水下部分為帶有系泊的圓柱結(jié)構(gòu)，其水平方向運動響應(yīng)往往具有較低的自振頻率，容易在低頻波浪力(源于非線性的差頻效應(yīng))作用下發(fā)生共振響應(yīng)，使結(jié)構(gòu)發(fā)生大幅水平慢漂。當(dāng)浮體的瞬時位置大幅偏離初始位置時，基于初始平衡位置的攝動展開法會存在較大誤差。針對這一問題，采用兩次展開方法，對大幅慢漂運動開展時域模擬研究。對雙色波作用下自由漂浮圓柱的大幅運動響應(yīng)問題進行數(shù)值分析，并與采用基于初始平衡位置的攝動展開法的計算結(jié)果進行了對比。結(jié)果表明，采用新的兩次展開法可以計算出波浪遭遇頻率的變化和波浪漂移阻尼，而這無法從基于初始平衡位置的常規(guī)攝動展開法中得到，體現(xiàn)了兩次展開法在分析大幅慢漂問題上的優(yōu)勢。

兩次展開；二階時域方法；小波變換；大幅慢漂；波浪漂移阻尼

Abstract：Large drift phenomenon of Spar platform is often observed in deep sea,whose natural frequency in horizontal direction is very low because the difference frequency of the nonlinear wave is close to the natural frequency of the structure.There is a large error if the perturbation expansion based on the initial equilibrium location is adopted when the floater is far away from the initial equilibrium position.For this reason,a twice-expansion time domain method is employed in this paper to study the problem.Using this method,the motion response of a free cylinder is evaluated,and the results are compared with those by perturbation expansion method based on initial equilibrium location.The comparisons show that the twice-expansion method demonstrates the performance on calculating the change of encountering frequency between waves and structures as well as the wave drift damping,while the original method does not.This shows the advantages of the twice-expansion method in solving the large drift problem.

Keywords：twice expansion; second-order time domain method; wavelet transform; large drift motion; wave drift damping

深海中的大型海洋平臺，如Spar平臺，通過錨泊系統(tǒng)定位于海中，水下部分為圓柱結(jié)構(gòu)。Spar平臺系統(tǒng)在水平運動方向通常具有很低的固有頻率，特別容易在非線性波浪差頻的作用下發(fā)生大幅度的水平往復(fù)慢漂運動。當(dāng)平臺發(fā)生大幅慢漂現(xiàn)象時，若將參考系設(shè)置在結(jié)構(gòu)物上，則可以簡化為水流與結(jié)構(gòu)物的相對運動問題。關(guān)于這類問題，Nossen等人[1]在頻域內(nèi)計算了三維固定在波流環(huán)境中受到的一階波浪力，而沒有涉及運動物體的求解計算。Kinoshita和Bao[2,3]通過理論推導(dǎo)獲得了單個圓柱和一列圓柱在波浪和水流中的波浪力和波浪漂移阻尼，但理論推導(dǎo)具有一定的局限性，不能計算任意形狀的物體。Kim[4]和Liu[5]隨后在時域內(nèi)求解了波流與結(jié)構(gòu)物的耦合作用產(chǎn)生的一階運動響應(yīng)，但未計算高階次的問題。上述學(xué)者研究的都是水流流速不變情況下的波浪與結(jié)構(gòu)物耦合問題，而對于浮體結(jié)構(gòu)物的大幅慢漂問題，漂移速度一直在變化，因此相當(dāng)于變速水流與波浪中結(jié)構(gòu)物的相對運動問題，至今相關(guān)研究較少。

同時，針對海洋平臺系泊系統(tǒng)以及波浪的非線性，大量學(xué)者開展了相關(guān)的研究工作。Sarkar和Eatock Taylor[6,7]提出了對波陡和系泊纜索非線性的雙參數(shù)展開方法，并采用半解析的方法計算了有限水深中帶有系泊纜索的截斷圓柱問題。Ormberg和Larsen[8]、Kim[9,10]、Garrett[11]計算了船體與系泊纜索的耦合作用，這些工作主要關(guān)注了系泊纜索的計算方法，而水動力荷載部分則采用比較簡單的線性模型。Low[12-14]同樣認為，深水結(jié)構(gòu)物和浮體的耦合作用中存在大量的非線性影響，并提出一種時域和頻域混合的方法，求解結(jié)構(gòu)物的運動響應(yīng)和系泊系統(tǒng)的張力。Kim等人[15]采用基于頻域攝動展開的Cummins方法計算了波浪荷載下漂浮結(jié)構(gòu)物的水動力荷載。Yang等人[16,17]采用時域的攝動展開法計算水動力荷載。上述的攝動展開均是直接基于初始平衡位置進行，沒有考慮結(jié)構(gòu)物發(fā)生大幅慢漂現(xiàn)象對分析結(jié)果的影響。

深海結(jié)構(gòu)物的大幅度慢漂運動，其幅度甚至?xí)笥诓ㄩL尺度。在此種情況下，若仍采用基于初始平衡位置的攝動展開法，無法準確滿足攝動展開的小振幅假設(shè)，計算會存在較大誤差，而且也無法考慮由于物體運動導(dǎo)致波浪力相位的變化和入射波浪與結(jié)構(gòu)物遭遇頻率的變化。此外，基于初始平衡位置的攝動展開法無法計算結(jié)構(gòu)物發(fā)生大幅慢漂運動時產(chǎn)生的波浪漂移阻尼，這部分阻尼會影響結(jié)構(gòu)物的自振頻率和運動響應(yīng)幅值?；谏鲜鲈?，本文采用了一種在時域內(nèi)求解波浪與結(jié)構(gòu)物相互作用的兩次展開的方法。其中，結(jié)構(gòu)物的運動被分解為兩部分：一是低頻運動；另一個是在低頻運動基礎(chǔ)上的波頻振蕩運動。低頻運動位置稱為瞬時平均位置，而在此平均位置進行波頻振蕩的運動需要通過攝動展開計算得到。此時，這個攝動展開的參數(shù)可以滿足小振幅的基本假設(shè)。由于此瞬時平均位置需要通過對總體運動響應(yīng)進行再次展開獲得，故此方法稱之為兩次展開法。

1 數(shù)學(xué)模型

1.1坐標(biāo)系的定義和運動響應(yīng)的分解

圖 1 坐標(biāo)系及計算域定義Fig.1 Definition sketch of coordinate systems and computation domain

(1)

(2)

H中的元素涉及一階轉(zhuǎn)動分量的乘積，是二階量。

對物體由波頻引起的運動響應(yīng)在平均位置進行攝動展開如下，

ξ(t)=εξ(1)(t)+ε2ξ(2)(t)+…

(3)

α(t)=εα(1)(t)+ε2α(2)(t)+…

(4)

式中：上標(biāo)(1)和(2)分別表示關(guān)于波陡ε=kA的一階和二階分量。

1.2控制方程和邊界條件

由于非線性波浪差頻的影響，浮體會發(fā)生低頻漫漂運動，在大地坐標(biāo)系下的速度勢可以寫為

(5)

式中：φc為物體低頻運動而使周圍流場產(chǎn)生定常擾動速度勢，與入射波浪無關(guān)，為零階勢；φw為大地坐標(biāo)系下的波浪速度勢，可以按照波陡ε展開為一階和二階分量。φc和φw均滿足拉普拉斯方程和海底不可滲透邊界條件。

(6)

式中：φc為單位速度勢。將式(5)代入伯努利方程，

(7)

物體慢漂運動的角頻率相對波頻運動的角頻率是小量，而低頻運動速度對時間變化率和速度平方的計算均會出現(xiàn)角頻率的平方項，可以忽略。因此，伯努利方程可以近似化簡為

(8)

速度勢Φ滿足瞬時物面邊界條件，不可滲透海底條件和動力學(xué)及運動學(xué)瞬時自由水面條件。將物面邊界條件和水面邊界條件分別在低頻平均位置和靜水面做泰勒級數(shù)展開，得到一階和二階的物面和水面邊界條件

(9)

(10)

(11)

(12)

(13)

(14)

為了書寫簡單，式中用Ξ(k)=ξ(k)+α(k)×(x-x0) 表示k階平移和轉(zhuǎn)動的聯(lián)合運動。

上述部分為兩次展開方法中的第一次展開，即基于低頻平均位置做關(guān)于未知量的攝動展開和關(guān)于邊界條件的泰勒級數(shù)展開，不同于以往在初始平衡位置展開。本文推導(dǎo)過程考慮了物體做低頻大幅慢漂運動的影響，因此和基于初始平衡位置的攝動展開方法相比有根本的不同。

1.3積分方程和波浪力

波浪的散射勢和低頻運動產(chǎn)生的定常擾動單位速度勢均需要通過求解積分方程得到，對于波浪散射勢φs，選取Rankine源和它關(guān)于海底的鏡像作為格林函數(shù)：

(15)

式中：R=[(x-x0)2+(y-y0)2+(z-z0)2]1/2,R1=[(x-x0)2+(y-y0)2+(z+z0+2d)2]1/2,x0(x0,y0,z0) 和x(x,y,z)分別表示源點坐標(biāo)和場點坐標(biāo)。

利用格林第二定理，可以得到關(guān)于計算域邊界上速度勢的邊界積分方程：

(16)

邊界S包括淹沒于水中的物體表面SB和從物體到阻尼層外邊界的有限靜水面SF。

而對于定常擾動的單位速度勢φc，同樣可以利用格林第二定理建立物體表面的積分方程來計算其物面上的速度勢。

(17)

自由水面各點的速度勢φc同樣通過式(17)求得，固角系數(shù)選取為1.0。在式(17)中，格林函數(shù)為

(18)

其中:

(19)

通過積分方程求出速度勢函數(shù)以后，作用在瞬時物體濕表面上的壓力可以通過Taylor級數(shù)展開變換到低頻平均物體表面上來表示。按照ε進行分離，力的分量可以按照不同的階數(shù)分解為：

(20)

一階項波浪力的具體表達式為:

(21)

(22)

二階項可以表達為

(23)

(24)

(25)

(26)

(27)

1.4低頻平均位置計算

根據(jù)動量守恒定理，可以分別得到浮體一階和二階運動方程

(28)

(29)

式中：[M]為質(zhì)量矩陣,[B]為粘性阻尼矩陣,[C]為恢復(fù)力矩陣,xc=(xc,yc,zc)是物體的質(zhì)心。E的表達式為:

(30)

為了得到更準確的描述浮體大幅慢漂運動的結(jié)果，本文提出的兩次展開方法中的第一次展開需要在每一時刻的平均位置進行。為此，平均位置的獲得在本方法中顯得至關(guān)重要。根據(jù)模擬得到的運動響應(yīng)的時間序列，使用小波變換方法找到其平均位置[18]。

結(jié)構(gòu)物在大地坐標(biāo)系下運動響應(yīng)的時間序列為χ(nΔt)=(χ1,χ2)，小波函數(shù)為Ψ0(η)。計算在縱蕩和橫蕩兩個方向分別進行，取任意一個方向為例，對于離散序列的連續(xù)小波變換可以定義為具有某一尺度的χ(nΔt)和Ψ0(η)轉(zhuǎn)換形式的卷積

(31)

(32)

為了和連續(xù)小波變換近似，式(31)需要對任意一個尺度做N次卷積，其中N為時間序列點的個數(shù)。通過卷積定理，上式可轉(zhuǎn)換為χ(t)和ψ*傅里葉變換乘積的逆傅里葉變換：

(33)

(34)

再對連續(xù)小波變換進行逆變換，即可得到原始的時間序列

(35)

式中:Cδ為重構(gòu)系數(shù)，為0.776;ψ0(0)為π-1/4。上述部分為兩次展開方法中的第二次展開，通過這次數(shù)值展開求得低頻平均位置，以此作為下一步計算的攝動展開和泰勒級數(shù)展開的平均位置。

由于物體表面網(wǎng)格和自由水面網(wǎng)格跟隨物體一起做水平慢漂運動，當(dāng)物體和水面網(wǎng)格移動到新位置時，當(dāng)前位置的自由水面的速度勢和波面高度需要依據(jù)前一時刻的結(jié)果進行計算。同時，因為網(wǎng)格整體平移，不存在網(wǎng)格變形和拓撲結(jié)構(gòu)的改變，本文采用泰勒展開的方法對當(dāng)前網(wǎng)格上的速度勢和波面高度進行插值計算。

(36)

式中：f表示自由水面各節(jié)點的散射勢φs和波面高度ηs。當(dāng)波浪非線性較強的時候，水線上的點由于波面變化較大而不利于數(shù)值穩(wěn)定，本文利用五點光滑法對節(jié)點的波面和速度勢進行光滑。

2 算例分析與討論

應(yīng)用兩次展開方法，本文計算了截斷圓柱在雙色波作用下的大幅慢漂運動。計算結(jié)果與基于初始平衡位置的攝動展開方法進行比較，后文稱之為直接展開法。文章中選取雙色波的差頻與系統(tǒng)的自振頻率較為接近，以激發(fā)圓柱大幅慢漂運動。

圓柱半徑a=1.0 m，吃水T=2.5 m，計算水深H=10.0 m。結(jié)構(gòu)只在水平方向施加k11=1 300 N/m的線剛度，此時結(jié)構(gòu)在水平方向的自振頻率ω0約為0.3 rad/s。而對于這一系列雙色波，波幅均為0.2 m，角頻率的和頻為6.0 rad/s，沿著x軸正向傳播，不同之處是每組波浪的差頻不同，具體參數(shù)如表 1所示。

表 1 雙色波參數(shù)Tab.1 Parameters of bichromatic waves

對不同Δω下的二階運動響應(yīng)進行傅里葉變換得到不同波浪差頻下的二階運動響應(yīng)幅值譜，如圖 2所示，選取其中幾個工況下的二階運動響應(yīng)對比分別展示結(jié)果，圖 3為Δω=0.20 和0.25 rad/s時二階運動響應(yīng)的對比，圖 4為Δω=0.35 和0.40 rad/s時二階運動響應(yīng)的對比。

圖 2 不同Δω下二階縱蕩響應(yīng)傅里葉變換的幅值Fig.2 Fourier transform amplitude of second order surge motion under different Δω

從圖 2的結(jié)果可以看出，直接展開方法得到二階運動響應(yīng)幅值在結(jié)構(gòu)自振頻率附近達到最大值，而兩次展開方法計算出的二階運動響應(yīng)幅值在小于ω0處達到最大。同時，兩次展開方法計算出的峰值也小于直接展開方法的峰值。發(fā)生這種情況的原因是在兩次展開法計算中考慮了物體大幅慢漂運動，結(jié)果中包含了波浪漂移阻尼這一項，而這一阻尼項十分重要，無法從直接展開中計算得到。因為波浪漂移阻尼的存在不僅影響了結(jié)構(gòu)的自振頻率，同時影響了運動響應(yīng)的幅值，因此在兩次展開法中，結(jié)構(gòu)雙色波差頻在小于ω0的情況下發(fā)生最顯著的共振現(xiàn)象。

相應(yīng)地，系統(tǒng)自振頻率的變化影響了不同雙色波差頻情況下結(jié)構(gòu)物的二階運動響應(yīng)。當(dāng)雙色波差頻較小時，如圖 3所示，即Δω為0.20和0.25 rad/s，比較接近兩次展開法得到的系統(tǒng)自振頻率，故兩次展開法可以計算出更大的二階運動響應(yīng)；同理，當(dāng)雙色波差頻較大時，如圖 4所示，即Δω為0.35和0.40 rad/s，與直接展開法的系統(tǒng)自振頻率接近，故二次展開法計算得到的運動響應(yīng)結(jié)果偏小。因此，當(dāng)結(jié)構(gòu)物發(fā)生大幅慢漂運動時，若仍采用直接展開方法計算結(jié)構(gòu)物的運動響應(yīng)和收到的波浪力，會產(chǎn)生較大的誤差。

圖 3 二階縱蕩運動響應(yīng)對比的時間歷程線Fig.3 Time histories of second-order surge displacement

圖 4 二階運動響應(yīng)對比的時間歷程線Fig.4 Time histories of second-order surge displacement

3 結(jié) 語

本文采用兩次展開的方法來求解非線性波浪作用下結(jié)構(gòu)物慢漂運動。在兩次展開方法中，物體運動由低頻慢漂運動和波頻振蕩運動組成。速度勢分解成波浪引起的一階、二階速度勢和物體低頻運動產(chǎn)生的零階速度勢。波浪引起的各分量根據(jù)波陡進行攝動展開，而物面邊界條件和水面邊界條件在低頻運動的平均物面和靜水面上采用泰勒級數(shù)展開求得，分別建立積分方程來求解波浪的散射勢和物體低頻運動產(chǎn)生的速度勢。初始條件采用二階Stokes波浪場，自由水面條件的時間積分采用四階預(yù)報校正的方法，浮體的運動方程求解采用Newmark-β法。求得運動響應(yīng)結(jié)果后使用小波變換的數(shù)值方法再次展開得到下一時刻攝動展開計算需要的平均位置。

采用本文提出的兩次展開方法計算結(jié)構(gòu)物在波浪中發(fā)生大幅運動時的波浪力和運動響應(yīng)，計算結(jié)果與基于初始平衡平均位置的攝動展開方法相比有較大區(qū)別。文章計算了雙色波與截斷圓柱的耦合作用，比較了兩種方法計算產(chǎn)生的二階運動響應(yīng)，從結(jié)果可以發(fā)現(xiàn)兩次展開方法考慮了攝動展開位置隨著低頻運動的變化，計算結(jié)果中包括了波浪慢漂阻尼的部分，使得結(jié)構(gòu)的自振頻率向低頻偏移，結(jié)構(gòu)的最大響應(yīng)峰值也低于基于初始平衡位置的直接展開方法的結(jié)果。通過波浪與結(jié)構(gòu)物的耦合作用計算，發(fā)現(xiàn)兩次展開方法不僅可以計算出波浪的結(jié)構(gòu)物遭遇頻率的變化，還可以計算出物體大幅慢漂運動時產(chǎn)生的波浪漂移阻尼部分。因此，本位提出的兩次展開方法，可以用于深海中平臺在大幅慢漂問題的時域計算。

[1] NOSSEN J,GRUE J,PALM E.Wave-forces on 3-dimensional floating bodies with small forward speed [J].Journal of Fluid Mechanics，1991,227：135-160.

[2] KINOSHITA T,BAO W.Hydrodynamic forces acting on a circular cylinder oscillating in waves and small current [J].Journal of Marine Science and Technology,1996,1：155-173.

[3] KINOSHITA T,BAO W,SUNAHARA S.The hydrodynamic forces acting on a cylinder array oscillating in waves and current[J].Journal of Marine Science and Technology,1997,2：135-147.

[4] KIM Y,KIM K H,KIM Y.Springing analysis of a seagoing vessel using fully coupled BEM-FEM in the time domain [J].Ocean Engineering,2009,36(11)：785-796.

[5] LIU Z,TENG B,NING D Z,et al.Wave-current interactions with three-dimensional floating bodies [J].Journal of Hydrodynamics,2010,22(2)：229-240.

[6] SARKAR A,EATOCK TAYLOR R.Low-frequency responses of non-linearly moored vessels in random waves：development of a two-scale perturbation method [J].Applied Ocean Research,1998,20(4)：225-236.

[7] SARKER A,EATOCK TAYLOR R.Low-frequency responses of non-linearly moord vessels in random waves：coupled surge,pitch and heave motions [J].Jounal of Fluids and Structures,2001,15：133-150.

[8] ORMBERG H,LARSEN K.Coupled analysis of floater motion and mooring dynamics for a turret-moored ship [J].Applied Ocean Research，1998,20(1-2)：55-67.

[9] KIM M H,KOO B,MERCIER RM,et al.Vessel/mooring/riser coupled dynamic analysis of a turret-moored FPSO compared with OTRC experiment [J].Ocean Engineering,2005,32(14-15)：1 780-1 802.

[10] TAHAR A,KIM M H.Coupled-dynamic analysis of floating structures with polyester mooring lines [J].Ocean Engineering,2008,35(17-18)：1 676-1 685.

[11] GARRETT D L.Coupled analysis of floating production systems [J].Ocean Engineering,2005,32(7)：802-816.

[12] LOW Y M,LANGLEY R S.A hybrid time/frequency domain approach for efficient coupled analysis of vessel/ mooring/ riser dynamics [J].Ocean Engineering,2008,35(5-6)：433-446.

[13] LOW Y M.Prediction of extreme response of floating structures using a hybrid time/frequency domain coupled analysis approach [J].Ocean Engineering,2008,36：1 416-1 428.

[14] LOW Y M.Extending a time/frequency domain hybrid method for riser fatigue anaylsis [J].Applied Ocean Research,2011,33：79-87.

[15] KIM B W,HONG G S,JIN H K,et al.Comparison of linear spring and nonlinear FEM methods in dynamic coupled analysis of floating structure and mooring system [J].Journal of Fluids and Structures,2013,42：205-227.

[16] YANG M D,TENG B,NING D Z,et al.Coupled dynamic analysis for wave interaction with a truss spar and its mooring line/riser system in time domain [J].Ocean Engineering,2012,39：72-87.

[17] YANG M D,TENG B,XIAO L F,et al.Full time-domain nonlinear coupled dynamic analysis of a truss spar and its mooring/riser system in irregular wave [J].Science China-Physics Mechanics & Astronomy,2014,57(1)：152-165.

[18] CHRITOPHER T,GILBERT C P.A practical guide to wavelet analysis [J].Bulletin of American Meteorological Society,1998,79(1)：64-78.

Large drift motion of mooring cylinder by applying twice expansion method

JIN Ruijia,TENG Bin,LYU Lin

(State Key Laboratory of Coastal and Offshore Engineering,Dalian University of Technology,Dalian 116024,China)

P751

A

10.16483/j.issn.1005-9865.2016.05.002

1005-9865(2016)05-0011-09

2015-09-28

國家自然科學(xué)基金資助項目(51379032,51490672)；國家973資助項目(2011CB013703)

金瑞佳(1987-)，男，江蘇武進人，博士研究生，主要從事波浪對海上建筑物作用研究。E-mail ：492472669@qq.com