羅冬冬，朱仁慶

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

0 引 言

近年來，由于地球人口急劇增加和能源的不斷消耗，各國能源危機越發(fā)明顯。同時伴隨陸上石油資源的日益枯竭，使得世界油氣的開發(fā)重點逐步向海洋轉移，尤其是深海區(qū)域。海洋立管在深海海洋平臺中應用最為廣泛，是連接海底的資源與海上作業(yè)平臺的關鍵設備，但也是薄弱易損的構件之一。一旦采油立管在渦激振動作用下發(fā)生斷裂，帶來的不僅僅是經(jīng)濟上的巨大損失，嚴重的是將會造成巨大的環(huán)境污染?？偠灾?，目前對海洋立管進行VIV 問題的研究已成為國內(nèi)外的熱點之一。

經(jīng)過數(shù)十年來對渦激振動的研究，學者們雖然還未完全把握渦激振動的機理，但取得了許多階段性的成果，構成了當今渦激振動研究的基石。無論是實驗研究還是數(shù)值模擬，大多數(shù)的研究主要針對在均勻來流情況下的渦激振動，如張建僑［7］研究質(zhì)量比對柔性細長立管渦激振動的影響。在同一流速下，質(zhì)量比大的立管模型響應位移小，其中空管的渦激振動響應一直處于大振幅的鎖定狀態(tài)下。在相同流速下，質(zhì)量比大的立管模型所激起的模態(tài)更高。王嘉松等［4］對1 根長502.92 m，直徑為0.534 m 的細長立管的進行模擬，指出了當立管發(fā)生大的、非均勻彎曲變形時，多模態(tài)的振動和線性動力響應不應該被忽視。同時表明預張力的增加，能有效的減小響應振幅，從而有效地減小疲勞損傷。但真實海洋來流情況往往比較復雜，不可能出現(xiàn)這種均勻來流的情況，相反更類似于具有一定斜率的剪切流的情況。因此剪切流下的立管渦激振動的模擬更加貼近海洋中的真實情況，符合實際的工程需要。目前，對于剪切流下的立管渦激振動的研究還十分少，尤其在國內(nèi)無論是實驗研究還是數(shù)值模擬都十分稀少。在國外，2003 年ExxonMobil［1］公司在MARINTEK 的深水水池中利用旋轉的架子進行了均勻流和剪切流的一系列的立管渦激振動系列實驗，說明了只有在低階模態(tài)參與的響應中，軸向疲勞損傷十分重要。近年來，Remi Bourguet［3］研究了立管在剪切流下的鎖定情況。結果表明在所有模擬條件中，在渦的發(fā)放與垂向振動局部同步時，在高流速區(qū)域并且長度至少超過圓柱長度的30%，才會發(fā)生鎖定現(xiàn)象。在鎖定情況下，振幅與在相似雷偌數(shù)下受迫振動的剛性圓柱的局部頻率有關。梁勇［12］基于Matlab 編寫了預測渦激振動響應的計算程序，來模擬剪切流下的立管。結果表明約束條件和立管密度對立管的振動響應有著十分重要的影響。

本文利用Workbench 平臺，建立三維流固耦合模型，自編了進口剪切流的程序，采用了k － ω 模型，結合基于薄殼理論的有限元方法，并通過一種新的方法System Coupling 實現(xiàn)流－固耦合交界面的數(shù)據(jù)交換，從而實現(xiàn)剪切流下三維流固耦合模擬。與相關文獻進行對比，研究立管的振動變形響應及尾渦形態(tài)特征等問題。

1 數(shù)值模擬方法概述

1.1 兩場控制方程

1.1.1 流體場控制方程

流場海水可看做粘性不可壓縮的流體，流場滿足以下微分方程:

1)連續(xù)方程

連續(xù)方程是質(zhì)量守恒定律在流體力學中的表現(xiàn)形式，該方程可以表述為:單位時間內(nèi)流入該微元的凈質(zhì)量等于單位時間流出該微元的凈質(zhì)量。對于三維流動，其平面直角坐標形式可表示如下:

式中u，v，w 分別為流體在x，y，z 方向上的速度，m/s。

2)動量方程

動量方程是動量守恒定律在流體力學中的表現(xiàn)形式，其本質(zhì)是滿足牛頓第二定律。該方程可以表述為:對于一給定的流體微元，其具有的動量對時間的導數(shù)等于作用在該微元上的各種力之和。其數(shù)學表達形式為Navier －Stokes 方程(簡稱N －S)方程，對于不可壓縮流體，寫成Cartesian 坐標形式為:

式中:f 為流體在x，y，z 各方向上所受質(zhì)量力，m/s2;p 為流體壓強，Pa;ρ 為密度，kg/m3;ν 為運動粘性系數(shù)，m2/s。

1.1.2 固體場控制方程

對海洋立管的動力響應采用基于三維實體單元的有限元方法進行模擬。對運動方程進行有限元離散后可得到有限元方程:

2 算例分析

2.1 模型的建立

本文模型參數(shù)如表1 所示。

表1 立管模型參數(shù)Tab.1 Parameters of riser model

進口速度為剪切流v=(0.03 × Z +0.2)m/s，Z為立管長度，立管軸線與流場速度入口邊界距離為10 D，與出口邊界距離為30 D，與兩側對稱邊界距離為10 D?？紤]的現(xiàn)有的計算機能力，使用較精密的平面網(wǎng)格體系，如圖1 所示。

圖1 二維網(wǎng)格拓撲圖Fig.1 Two dimensional grid

流場采用k － ω 模型，迭代方法選用具有高穩(wěn)定性的SIMPLEC 算法。選取立管的圓柱面為流固耦合面，受到動網(wǎng)格的限制，時間的離散方式選用了一階隱式(1nd－Order－Implicit)。時間步長為t=0.002 s，計算總時長為6 s。

圖2 結構有限元網(wǎng)格圖Fig.2 Finite element of structure

對于結構模型的網(wǎng)格劃分如圖2 所示，同樣結構中，設定立管的圓柱面為流固耦合面。同時，為了保持迭代耦合的一致性，時間步長與時長同流場相一致。

2.2 海洋立管的運動分析

2.2.1 順流向的運動分析

在順流向方向，立管因為受到拖曳力的影響，發(fā)生彎曲變形，如圖3 所示，在t = 2.57 s 時，彎曲達到最大，立管的最大彎曲變形處發(fā)生在約z/L=0.6 處，最大變形為1.2 D，此后，流場趨于穩(wěn)定，立管變形很小，只是在最大位置處產(chǎn)生小幅振動，振動幅值為0.1 D。如圖4 所示。與文獻［11］對比發(fā)現(xiàn):其運動趨勢與均勻流下立管的運動趨勢相同，并都多出現(xiàn)了順流向的不可忽略的振動，只有最大彎曲變形位置向上移動。該結果說明由于剪切流流速遞增的特性，立管受到水流的拖曳力的中心向上移動，使得立管的最大彎曲變形的位置上移，但順流向上對立管的運動趨勢的影響與均勻流下相同。同時說明了System Coupling 實現(xiàn)流固耦合具有一定的可行性。

圖3 順流向運動效果圖Fig.3 In-line motion of riser

圖4 最大位移處點的位移隨時間的變化圖Fig.4 Scatter diagram of displacement varying with time of point with largest displacement

2.2.2 橫流向的運動分析

從圖5 可看出，在橫流向上，立管的振動在1階、2 階振動模態(tài)之間往復，立管的彎曲變形呈非對稱性，當流體趨于穩(wěn)定時，最終鎖定在2 階振動模態(tài)，同時出現(xiàn)了3 階振動模態(tài)的情況。

而根據(jù)文獻［11］得到立管自振頻率公式為:

式中:md為單位長度隔水管排開海水的質(zhì)量;CA為附加質(zhì)量系數(shù)，一般取1.0;I 為截面慣性矩;T 為預緊力。同時根據(jù)“鎖定”的定義:

圖5 立管不同時刻橫流向的振動效果圖Fig.5 Cross-flow vibration shape of riser at different moment

式中:f 為隔水管實際振動頻率;fn為隔水管理論自振頻率;fst為渦脫落頻率。

由此得到前3 階振動模態(tài)的激發(fā)流速:1 階模態(tài)的激發(fā)流速0.188 m/s;2 階和3 階模態(tài)的激發(fā)流速為0.389 m/s 和0.607 m/s。而本文的最高流速為0.488 9 m/s，最低流速為0.2 m/s，并沒有達到3 階模態(tài)的激發(fā)流速。本文激發(fā)2 階振動模態(tài)的來流約占立管長度的40%，符合上文Remi Bourguet 的指出:在渦的發(fā)放與垂向振動局部同步時，在高流速區(qū)域并且長度至少超過圓柱長度的30%，鎖定才會發(fā)生。所以立管的振動最終鎖定在2 階振動模態(tài)。同時激發(fā)1階振動模態(tài)的來流約占立管長度的60%，對立管的影響也很大，因此可能導致了3 階振動模態(tài)的出現(xiàn)。該現(xiàn)象表明:剪切流增加了立管振動的復雜性，立管出現(xiàn)了比剪切流最大流速下的振動模態(tài)更高1 階的3 階模態(tài)振動，但最終鎖定模態(tài)處于2 階振動模態(tài)。

2.3 海洋立管的流場分析

圖6 為立管在前期和后期不同截面的尾流區(qū)渦的三維特性圖。圖中5 個截面分別對應Z/L=0.1，0.3，0.5，0.7，0.9 的渦量圖。在t = 1.2 時，由于剪切流的關系，在Z/L=0.7 和0.9 的截面上都形成了渦的發(fā)放，而其他3 個界面渦還沒有形成。當最終整個立管的VIV 穩(wěn)定后，渦的發(fā)放形式主要是“2S”形的渦。從圖中可以看出:由于立管的運動和剪切流的作用，渦在立管不同的位置有明顯的區(qū)別，在Z/L=0.1，0.5，0.9 處，立管在橫流向的運動范圍比較小;兩渦在橫流向上的距離明顯小于Z/L=0.3，0.7 處。

在立管渦的發(fā)放穩(wěn)定前，還伴隨其他不同形式的渦，如圖7 所示。圖7 為t = 2.4 s 時，Z/L=0.3、0.5 的截面的渦量圖。圖7(a)為Z/L=0.3 截面渦量圖，此時渦為“P +S”的形式。而到圖7(b)時即Z/L =0.5 截面，渦為“2P”的形式。

圖6 立管不同截面不同時刻的渦量圖Fig.6 Vorticity of riser at different sections for different moment

圖7 t=2.4 s 時渦量圖Fig.7 Vorticity of riser at t=2.4 s

3 結 語

本文基于計算流體力學(CFD)和計算結構力學(CSD)方法，借助商業(yè)軟件Ansys 14.0，參考長徑比為482 的均勻流的立管模擬，實現(xiàn)了VIV 的三維數(shù)值模擬。結果表明，數(shù)值模擬與均勻流有部分相似;大長徑比柔性隔水管發(fā)生VIV 時流場受到明顯擾動，渦的脫落呈現(xiàn)如2P、2S、P+S 等模態(tài);立管最大位移處向上移動，立管在橫流向上出現(xiàn)了大的彎曲變形，出現(xiàn)了1 階、2 階模態(tài)之間的轉變，最終鎖定在2 階模態(tài)。

［1］TOGNARALLI M A．VIV response of a long fexible cylinder in uniform and linearly sheared currents［C］//Offshore Technology Conference，2004．

［2］FARSHIDIANFAR A，ZANGANEH H． A modified wake oscillator model for vortex-induced vibration of circular cylinders for a wide range of mass － damping ratio［J］．Journal of Fluids and Structures，2010．

［3］BAHMANI M H，AKBARI M H． Effects of mass and damping ratios on VIV of a circular cylinder［J］． Ocean Engineering，2010．

［4］REMI B，MICHAEL S T． Lock-in of the vortex-induced vibrations of a long tensioned beam in shear flow［J］．Journal of Fluids and Structures，2011．

［5］WANG Jia-song，ZHAN Lu-lu． Three dimensional numerical simulation of vortex induced vibration for a 500 － m － long marine riser［C］//ISOPE，2012．

［6］LIN Li-ming，ZHONG Xing-fu，WU Ying-xiang． Vortexinduced vibrations of a circular cylinder with different geometric disturbances［C］//ISOPE，2012．

［7］HUANG Zhi-yong，CUI Wei-cheng．Vortex-induced vibration of flexible riser in time domain［D］． Shanghai:Shanghai Jiao Tong University，2008．

［8］ZHANG Jian-qiao，SONG Ji-ning，LV Lin，et al．Analysis on the vortexoinduced vibrations experimental results of flexible risers［J］． China Offshore Platform，2009．

［9］TANG Guo-qiang，TENG Bin，LV Lin，et al．Development and application of an empirical model for vortex-induced vibration of flexible deep water riser［J］．China Offshore Platform，2010．

［10］CAI Jie． Investigation on numerical simulation of the fluidsolidInteraction characteristics for oceanic riser［D］．Shanghai:Shanghai Jiao Tong University，2010．

［11］CAI Jie，YOU Yun-xiang，LI Wei，et al． The VIV characteristics of deep-sea risers with high aspect ratio in a uniform current profile［J］．Chinese Journal of Hydrodynamics，2010．

［12］WANG Cheng-guan，TIAN Zhong-xu，QIAO Xin-qi． Study on three dimensional numerical simulation of VIV characteristics for marine riser［J］． Chinese Journal of Hydrodynamics，2011．

［13］LIANG Yong．Response analysis of flexible riserssubjeted to Vortex-induced vibration in sheared flow［D］． Harbin:Harbin Institute of Technology，2011．