吳學(xué)敏，黃維平

(中國海洋大學(xué) 山東省海洋工程重點(diǎn)實(shí)驗(yàn)室，山東青島 266100)

渦激振動(dòng)(Vortex Induced Vibration，VIV)可導(dǎo)致深海大柔性立管發(fā)生疲勞破壞。深水立管在深海復(fù)雜環(huán)境中，除與海底、平臺(tái)接觸外無任何支撐。研究發(fā)現(xiàn)，深海立管在漩渦脫落影響下發(fā)生多種振動(dòng)現(xiàn)象:多模態(tài)振動(dòng)、高階模態(tài)振動(dòng)、非對(duì)稱彎曲大變形等［1－6］。正確預(yù)測深水立管渦激振動(dòng)特征對(duì)海洋立管設(shè)計(jì)意義重大。

由于深水立管自身特性及所處復(fù)雜深海環(huán)境等特點(diǎn)，研究深水立管渦激振動(dòng)特征用試驗(yàn)方法不易操作。據(jù)獲得作用于結(jié)構(gòu)的流體力方法不同，預(yù)測深水立管渦激振動(dòng)主要有①半經(jīng)驗(yàn)方法，該法中的尾流振子模型較常見，典型的為基于非線性van der Pol方程的尾流振子模型，該模型因具有自激、自限性特征，故能捕捉到渦激振動(dòng)現(xiàn)象，如:鎖定現(xiàn)象。但在采用該模型時(shí)，需作假設(shè):如當(dāng)鎖定現(xiàn)象發(fā)生時(shí)，假設(shè)漩渦脫落沿立管長度全相關(guān)。②預(yù)測方法(Computational Fluid Dynamics，CFD)，該方法被公認(rèn)為較理想的預(yù)測細(xì)長立管渦激振動(dòng)方法，但由于三維模擬計(jì)算量大等問題，此類模擬較少。而大多數(shù)CFD數(shù)值模擬方法中，均選取立管若干截面，在各個(gè)截面上進(jìn)行二維流體數(shù)值模擬計(jì)算。此二維模擬忽略展向尾流的影響。Yamamoto等［7］用CFD方法模擬細(xì)長圓柱體渦激振動(dòng)，結(jié)果與試驗(yàn)結(jié)果吻合較好;秦延龍等［8］給出較全面的CFD方法研究深水立管渦激振動(dòng)綜述。

目前對(duì)海洋立管的動(dòng)力響應(yīng)分析模型中，大多采用線性化彈性細(xì)長桿模型，忽略了深水立管的幾何非線性效應(yīng)。由于海洋立管渦激振動(dòng)不僅會(huì)出現(xiàn)多模態(tài)參與振動(dòng)現(xiàn)象，且會(huì)出現(xiàn)彎曲大變形現(xiàn)象。因此深水立管的幾何非線性效應(yīng)不可忽略。對(duì)此本文提出考慮深海立管大變形彎曲振動(dòng)模型，并用有限元方法分析均勻流中某大長徑比深海頂張力式立管的渦激振動(dòng)特性。

有研究表明，深海立管的順流向渦激振動(dòng)振幅與橫向振動(dòng)振幅具有相同數(shù)量級(jí)［9］。當(dāng)立管順流向與橫流向固有頻率比與順流向與橫流向激發(fā)頻率比相等時(shí)，立管的渦激振動(dòng)出現(xiàn)雙向共振特征，立管呈“8”字形狀［10－11］。即順流向振動(dòng)對(duì)深水立管疲勞損傷同等重要。本文在所提新模型基礎(chǔ)上對(duì)深海立管進(jìn)行順流向、橫流向自由度渦激振動(dòng)大變形振動(dòng)響應(yīng)特征分析。

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

在笛卡爾坐標(biāo)系下建立運(yùn)動(dòng)方程［12－15］。假設(shè)立管不可伸長，即變形前后弧長s不變。圖1為深水頂張式立管微元段力學(xué)模型。由動(dòng)量、動(dòng)量矩守恒可得長ds微段立管單元平衡方程:

圖1 微元段Fig.1 Microelement

式中:F為截面內(nèi)力;q為立管單位長度分布外力;ρ為立管單位長度質(zhì)量;m為單位長度分布外力矩;ρ=ρtAt+ρiAi為截面力矩，為彎矩與扭矩之和，表示為:

式中:B=EI為截面彎曲剛度;H為截面扭矩。

由式(3)得:

將式(4)代入式(2)并令H=0、m=0，即忽略單元扭矩、外力矩，得:

式中:λ為Lagrange乘子;T(s，t)=r';F'為立管張力;κ2= －r'r?，κ(s，t)為局部曲率。

將式(5)代入式(1)得大撓度細(xì)長梁運(yùn)動(dòng)方程為:

將式(7)寫成頂張力式立管大變形彎曲振動(dòng)微分方程為:

2 數(shù)值計(jì)算

為驗(yàn)證模型式(8)、式(9)，本文開發(fā)出深水立管渦激振動(dòng)分析程序。由于深水立管渦激振動(dòng)的非線性特征，程序采用Newmark－β隱式積分法求解。Newmark－β法為直接積分法，β=0.25時(shí)無條件穩(wěn)定。

設(shè)大變形深海立管頂點(diǎn)為坐標(biāo)原點(diǎn)，u軸為垂直方向坐標(biāo)軸，v軸為水平方向坐標(biāo)軸。將彎曲振動(dòng)方程式(8)、(9)離散后得深海立管大變形彎曲振動(dòng)有限元方程為:

式中:［Mv］為立管橫向質(zhì)量矩陣:

［Mu］為立管豎向慣性系數(shù)矩陣:

［Kv］為立管橫向剛度矩陣:

［Ku］為立管豎向剛度矩陣:

［Cv］為立管橫向阻尼矩陣:

［Cv］= ［Mv］+ ［Kv］

［Cu］為立管豎向阻尼矩陣:

［Cu］= α［Mu］+ β［Ku］

其中:α，β 為瑞雷祖尼系數(shù)，單位 s－1，s，兩者取值為［16］:取體系基頻，取體系第三階自然頻率;{Fu}為立管豎向荷載向量;{Fv}為立管橫向荷載向量;{}為立管橫向加速度向量;{}為立管橫向速度向量;{v}為立管橫向位移向量;{}為立管豎向加速度向量;}為立管豎向速度向量;{u}為立管豎向位移向量;［N］為立管單元插值函數(shù)矩陣;l為單元長度;n為單元數(shù)量。

引入大變形，系統(tǒng)質(zhì)量矩陣、阻尼矩陣、剛度矩陣與節(jié)點(diǎn)曲率相關(guān)。因此在每個(gè)時(shí)間步，有限元方程的質(zhì)量矩陣、剛度矩陣、阻尼矩陣均未變化的。為提高計(jì)算精度，對(duì)每個(gè)時(shí)間步進(jìn)行迭代，時(shí)間步長dt=0.01 s。

3 算例分析

用以上程序?qū)? 500 m水深頂張式立管(TTR)進(jìn)行VIV大變形渦激振動(dòng)響應(yīng)特征分析。該TTR為雙套管(雙屏)生產(chǎn)立管，外套管外徑324 mm，內(nèi)套管外徑222 mm，油管外徑114 mm。基于彎曲剛度等效原則將立管等效為單層管分析模型，等效后內(nèi)徑為292 mm，頂張力系數(shù)取1.4。計(jì)算模型邊界條件一端固定、一端簡支。

3.1 計(jì)算結(jié)果驗(yàn)證

由于Shear7采用基于模態(tài)疊加頻域分析方法，計(jì)算結(jié)果為響應(yīng)的均方根值。因此，本文將計(jì)算結(jié)果統(tǒng)計(jì)分析，給出響應(yīng)均方根值。圖2、圖3為固定－簡支邊界條件的均方根位移與均方根應(yīng)力沿立管軸向分布曲線，虛線為本文分析模型計(jì)算結(jié)果，實(shí)線為商用軟件Shear7分析結(jié)果?？梢钥闯?，本文所提分析模型與Shear7模態(tài)疊加法計(jì)算結(jié)果吻合較好。

圖2 均方根位移曲線(固定－簡支)Fig.2 RMS Displacement along TTR

圖3 均方根應(yīng)力曲線(固定－簡支)Fig.3 RMS Stress along TTR

3.2 動(dòng)力特性分析

圖4為 Re=1.2×105時(shí)水深 Z=300 m、750 m、1350 m處三水平橫截面立管節(jié)點(diǎn)橫向振動(dòng)時(shí)程曲線。并經(jīng)傅里葉變換得圖5三水平截面處橫向振動(dòng)響應(yīng)頻譜。由圖5(a)看出，立管在Z=300 m處有六個(gè)峰值頻率:fv=0.029 9 Hz，0.069 9 Hz，0.119 9 Hz，0.159 9 Hz，0.199 9 Hz，0.259 9 Hz;在 Z=750 m 處有五個(gè)峰值頻率:fv=0.029 9 Hz，0.069 9 Hz，0.119 9 Hz，0.159 9 Hz，0.259 9 Hz;在 Z=1 350 m 處也有五個(gè)峰值頻率:fv=0.029 9 Hz，0.069 9 Hz，0.119 9 Hz，0.159 9 Hz，0.259 9 Hz。可見，立管水平位置不同，其控制模態(tài)亦不同，均發(fā)生多模態(tài)參與現(xiàn)象，參與的模態(tài)數(shù)目也不相同，且參與模態(tài)間頻率間隔較小。研究發(fā)現(xiàn)此時(shí)所對(duì)應(yīng)的頻率峰值均相同。

圖6為Re=1.2×105時(shí)水深 Z=300 m、Z=750 m、1 350 m處三水平橫截面立管表面某點(diǎn)順流向振動(dòng)響應(yīng)頻譜圖。由圖6看出:立管在Z=300 m處有六個(gè)峰值頻率:fv=0.029 9 Hz，0.069 9 Hz，0.119 9 Hz，0.159 9 Hz，0.199 9 Hz，0.259 9 Hz;Z=750 m 處有四個(gè)峰值頻率:fv=0.029 9 Hz，0.119 9 Hz，0.159 9 Hz，0.259 9 Hz;Z=1 350 m處有五個(gè)峰值頻率:fv=0.029 9 Hz，0.069 9 Hz，0.119 9 Hz，0.159 9 Hz，0.259 9 Hz?？梢姡瑯釉诹⒐苋齻€(gè)不同水平位置處，立管順流向振動(dòng)的控制模態(tài)也不同，但對(duì)應(yīng)的頻率峰值相同。均發(fā)生多模態(tài)參與振動(dòng)現(xiàn)象。

圖4 Re=1.2×105時(shí)立管橫向振動(dòng)時(shí)程曲線Fig.4 Time-history curve of cross－ flow vibration when Re=1.2 ×105

圖5 Re=1.2×105時(shí)立管橫向振動(dòng)頻譜圖Fig.5 Spectrogram of cross-flow vibration，when Re=1.2 ×105

圖6 Re=1.2×105時(shí)立管順流向振動(dòng)頻譜圖Fig.6 Spectrogram of in-line vibration，when Re=1.2 ×105

圖7 Re=6 300時(shí)立管橫向振動(dòng)頻譜圖Fig.7 Spectrogram of cross-flow vibration，when Re=6 300

圖8 當(dāng)Re=6 300時(shí)立管順流向振動(dòng)頻譜圖Fig.8 Spectrogram of in-line vibration，when Re=6 300

在Re=1.2×105時(shí)，立管的順流向、橫向振動(dòng)均發(fā)生多階模態(tài)參與振動(dòng)現(xiàn)象，頻率峰值基本相同。

圖7、圖8分別為Re=6 300時(shí)水深Z=300 m、750 m、1350 m處三個(gè)水平橫截面節(jié)點(diǎn)橫向、順流向渦激振動(dòng)響應(yīng)頻譜圖。由二圖可見，低雷諾數(shù)時(shí)立管渦激振動(dòng)以低階單模態(tài)為主，漩渦脫落并未激發(fā)多個(gè)模態(tài)參與振動(dòng)，此時(shí)，橫向頻率峰值為fv=0.024 9 Hz，順流向頻率峰值為fv=0.024 9 Hz，兩向頻率峰值相等。

圖9為兩種雷諾數(shù)時(shí)立管振動(dòng)橫向振幅達(dá)最大的變形圖。Re=6 300時(shí)，立管橫向最大振幅為0.195 5 m，發(fā)生在Z=825 m處;順流向最大振幅0.005 8 m，發(fā)生在Z=795 m處，且立管振動(dòng)以一階單模態(tài)振動(dòng)為主。當(dāng)Re=1.2×105時(shí)，立管橫向最大振幅為0.503 3 m，發(fā)生在Z=1 395 m處;順流向最大振幅發(fā)生在Z=810 m處，且漩渦脫落激起高階模態(tài)參與振動(dòng)。立管渦激振動(dòng)變形并非完全對(duì)稱，且有大變形現(xiàn)象。隨雷諾數(shù)的增大，立管高階模態(tài)被激發(fā)。

圖9 橫向振幅最大時(shí)立管振動(dòng)變形圖Fig.9 Deformation figures of riser，when the max CF amplitude

由數(shù)值模擬結(jié)果知，雷諾數(shù)不同，立管不各位置參與振動(dòng)的模態(tài)也不同。雷諾數(shù)較低時(shí)，立管渦激振動(dòng)以低階單模態(tài)振動(dòng)為主;隨雷諾數(shù)的增大，立管渦激振動(dòng)出現(xiàn)多模態(tài)參與振動(dòng)現(xiàn)象，甚至立管高階模態(tài)也被激發(fā)。渦激振動(dòng)振幅也隨雷諾數(shù)的增大而增大。雷諾數(shù)相同時(shí)，立管發(fā)生渦激振動(dòng)的不同時(shí)刻，其最大振幅發(fā)生的位置亦不同，且立管渦激振動(dòng)均發(fā)生非對(duì)稱大變形現(xiàn)象。

表1 Re=6 300時(shí)立管渦激振動(dòng)響應(yīng)振幅Tab.1 Vortex-induced vibration amplitude of riser when Re=6 300

本文在Re=6 300、Re=1.2×105兩種雷諾數(shù)下，將考慮大變形與不考慮大變形立管渦激振動(dòng)的最大振幅進(jìn)行比較，結(jié)果分別見表1、表2。由表1、表2數(shù)據(jù)知，大變形在低雷諾數(shù)時(shí)對(duì)振幅影響較小，在高模態(tài)數(shù)情況下對(duì)振幅影響增大。

表2 Re=1.2×105時(shí)立管渦激振動(dòng)響應(yīng)振幅Tab.2 Vortex-induced vibration amplitude of riser when Re=1.2 ×105

4 結(jié)論

本文在已有研究基礎(chǔ)上，考慮深水立管因其大長徑比特點(diǎn)，在深海易發(fā)生大變形現(xiàn)象，提出深水頂張力式立管大變形渦激振動(dòng)分析模型，研究立管在不同雷諾數(shù)下的渦激振動(dòng)特性。結(jié)論如下:

(1)深海立管多模態(tài)響應(yīng)與非對(duì)稱彎曲大變形等渦激振動(dòng)現(xiàn)象在低雷諾數(shù)情況下未能體現(xiàn);卻體現(xiàn)在高模態(tài)數(shù)情況下。

(2)考慮大變形的幾何非線性分析模型可用于深水、極深水立管系統(tǒng)渦激振動(dòng)響應(yīng)分析預(yù)測。該模型考慮深水立管在復(fù)雜深海環(huán)境下易發(fā)生大變形現(xiàn)象，將大變形幾何非線性特征引入深水立管渦激振動(dòng)分析中，為尋求更準(zhǔn)確模型分析深水立管方法提供參考。

［1］Modarres-Sadeghi Y，Chasparis F，Triantafyllou M S，et al.Chaotic response is a generic feature of vortex-induced vibrations of flexible risers［J］. Journal of Sound and Vibrations，2011，330(11):2565 －2579.

［2］Wu X D，Ge F，Hong Y S.A view of recent studies on vortex-induced vibrations of long slender cylinders［J］.Journal of Fluid and Structures，2011，28:292 －308.

［3］Williamson C H K，Govardhan R.A brief review of recent results in vortex-induced vibrations［J］.Journal of Wind Engineering and Industrial Aerodynamics，2008，6(6－7):714－735.

［4］Willden R H J，Graham J M R.Multi-modal vortex-induced vibrations if a vertical riser pipe subject to a uniform current profile［J］.European Journal of Mechanics B/Fluids，2004，23(1):209－218.

［5］So R M C，Liu Y，Cui Z X，et a1.Three-dimensional wake effects on flow-induced forces［J］.Journal of Fluids and Structures，2005，20(2):373 －402.

［6］Trim A D，Braaten H， LieH， eta1. Experimental investigation of vortex-induced vibration of long marine risers［J］.Journal of Fluids and Structures，2005，21(3):335 －361.

［7］Yamamoto C T，Meneghin J R.Numerical simulations of vortex-induced vibration on flexible cylinders［J］.Journal of Fluids and Structures，2004，19(4):467－489.

［8］秦延龍，王世澎.海洋立管渦激振動(dòng)計(jì)算方法進(jìn)展［J］.中國海洋平臺(tái)，2008，23(4):14 －17.QIN Yan-long， WANG Shi-peng.Computational methods progress of the vortex-induced vibration of ocean risers［J］.China Offsh Ore Platform，2008，23(4):14－17.

［9］Vandiver J K，Jaiswal V，Jhingran V.Insights on vortexinduced，traveling waves on long risers［J］.Journal of Fluids and Structures，2009，25(4):641 －653.

［10］Dahl J M，Hover F S，Triantafillou M S，et al.Dual resonance in vortex-induced vibrations at subcritical and supercriticalreynolds numbers［J］. JournalofFluid Mechanics，2010，643:395 －424.

［11］Dahl J M，Hover F S，Triantafillou M S.High harmonic forces and predicted vibrations from forced in-line and crossflow cylinder motions［C］.In:International Offshore and Polar Engineering Conference Proceedings，2008:481 －488.

［12］Nordgren R P.On computation of the motion of elastic rods［J］.ASME Journal of Applied Mechanics，1974，41(3):777－780.

［13］Garrett D L.Dynamic analysis of slender rods［J］.Journal of Energy Resources Technology，Transaction of ASME，1982，104(4):302－307.

［14］Pauling J R，Webster W C.A consistent large-amplitude analysis of the coupled response of a TLP and tendon system［J］. Proceedings ofOffshore Mechanics and Arctic Engineering(OMAE)，1986，3:126－133.

［15］Ma W，Webster W C.An analytical approach to cable dynamics:theory and user manual［R］.Sea Grant Project R/OE －26，1994.

［16］克拉夫R，彭津J，著.王光遠(yuǎn)，譯.結(jié)構(gòu)動(dòng)力學(xué)(第二版)［M］.北京:高等教育出版社，2007:183－185.