立管渦激振動流固耦合數值模擬

陳天賜，王瑞金，朱澤飛

(杭州電子科技大學機械工程學院，浙江 杭州 310018)

0 引 言

渦激振動(Vortex Induced Vibration，VIV)是海洋立管結構損傷的主要原因[1]。為了提高海洋資源利用率，衍生出一系列研究VIV的成果[2-8]。研究圓柱結構VIV的主要方法分為實驗方法和數值模擬方法，其中數值模擬方法又分為經驗模型和計算流體力學(Computational Fluid Dynamics，CFD)。大長徑比柔性圓柱結構的VIV響應不同于剛性圓柱體。為了解決CFD全尺度計算耗時及收斂較為困難的問題，部分學者使用二維切片法研究立管的振動響應。王哲等[9]采用OPENFOAM軟件把柔性圓柱簡化成若干的二維切片代替三維流場，計算時變頂端張力的渦激振動響應，數值模擬結果與實驗結果基本吻合。Deng等[10]使用厚切片代替二維薄切片，改進了薄二維切片的局限性，使用RANS湍流模型計算了均勻流和階梯流下圓柱的振動響應，發(fā)現切片厚度為1/50時能滿足三維計算的精度要求。大長徑比圓柱的渦激振動響應的三維效應明顯，切片法無法滿足其計算精度的要求。及春寧等[11]運用浸入邊界法對處于剪切來流條件下的柔性圓柱進行數值模擬，結果表明剪切率的增大使得流動方向呈現多模態(tài)共存現象，其中橫向振動模態(tài)向一階模態(tài)過渡。Li等[12]使用矢量形式內在有限元法(Vector Form Intrinsic Finite Element，VFIFE)耦合了流向與橫向的運動，運用經典范德波方程模擬結構上的脈動流體力，數值模擬結果與實驗結果吻合較好，證明了VFIFE法的有效性。Wang等[13]采用任意拉格朗日-歐拉坐標體系處理耦合面動網格，計算均勻來流和剪切來流下垂直立管的VIV響應。Han等[14]使用雙向流固耦合分析了立管的振動響應，并成功捕捉到相鄰階振型的轉換現象。但是，目前針對脈動流下的渦激振動響應研究甚少，本文采用ANSYS Workbench軟件，使用雙向流固耦合方法分析柔性立管的渦激振動響應，研究脈動流條件下脈動流頻率對立管渦激振動響應的影響。

1 理論分析

1.1 流場求解

先采用有限體積法離散化計算區(qū)域，再線性化N-S方程，從而求解獲得流場信息。N-S方程是由三大守恒定律推導而來，分別是質量、動量以及能量守恒。因為不涉及能量求解，故本文沒有列出能量守恒方程，其余2項守恒控制方程矢量表達如下：

(1)

(2)

式中，ρ為流體的密度，U為速度矢量，μ為動力黏度，p為壓強，f為單位質量流體體積力的矢量表達形式。

1.2 結構求解

立管單元離散化后，其結構響應控制方程如下：

(3)

1.3 流固耦合求解

根據求解方式，流固耦合求解的計算可以分為直接求解法和迭代求解法。直接求解法將流場和固體結構的數據信息放在同一個矩陣下進行求解，耗時長且占據計算資源。迭代求解法是分別在流場和結構求解器中進行求解，通過數據交換界面進行插值運算。本文使用分析軟件ANSYS Workbench進行計算，其中，Fluent模塊用于流場求解，Transient Structural模塊用于結構求解，并使用System Coupling模塊進行數據交換，最終完成耦合計算，計算原理如圖1所示。

圖1 流固耦合原理圖

2 求解參數設置

2.1 物理模型計算

參考文獻[15]的實驗參數，本文模擬了兩端鉸接并受頂端張力作用的PVC材料立管的VIV響應，其材料的物性參數如表1所示。

表1 實驗材料物性參數[15]

進行有預應力作用下的模態(tài)分析，取前四階模態(tài)分析的頻率，并與文獻[15]中實驗測得的第1階頻率進行對比，結果如表2所示。

表2 模態(tài)分析結果 單位：Hz

從表2可以看出，本文第1階有限元分析結果與文獻[15]的誤差不大，在接受范圍內，說明立管模型的建立較為合理，可用于后續(xù)計算。

2.2 求解設置

圖2 二維網格平面圖

圖3 流場整體網格

圖4 立管結構平面網格

2.3 數值模擬方法驗證

選取實驗中的一組實驗數據進行數值模擬方法驗證。頂端張力為60 N，速度為0.15 m/s，實驗與數值模擬的振幅數據均為立管中點的振幅數據，縱坐標的y為立管的橫向位移，得到的模擬結果與實驗結果如圖5所示。

圖5 數值模擬和實驗結果對比

從圖5中可以看出，數值模擬得到的振幅與實驗結果的整體趨勢十分接近，但最大振幅稍有出入，誤差在10%以下，說明仿真結果與實驗結果吻合良好，數值模擬方法正確可行。

2.4 脈動來流速度入口設置

目前，立管渦激振動響應的相關研究多是在均勻來流下進行的，脈動來流下的渦激振動響應研究比較少。由于海浪的影響，流動具有強烈的不穩(wěn)定性，所以不能只考慮定常流動。為此，本文將這種不穩(wěn)定來流簡化為正弦脈動變化形式，研究脈動頻率對渦激振動的影響。脈動來流的具體數學表達式[16]為：

(4)

式中，umax為流速最大值，Am為脈動振幅，f為脈動頻率，t為時間。

本文的數值模擬實驗主要研究脈動頻率對渦激振動的影響，故不考慮其余參數。參數設置如下：umax=0.25 m/s，Am=0.2，f分別取0.5fn，fn，2fn，3fn，其中fn為一階固有頻率，根據本文有限元分析結果，其值為3.18 Hz，頂端張力為260 N。

3 數值模擬結果與分析

3.1 立管振動響應分析

3.1.1 水動力系數分析

阻力系數CD和升力系數CL是描述結構受周期性力的2個無量綱參數，分別代表物體流向與橫向的受力情況。

(5)

(6)

式中，Fx，Fy分別為流向力和橫向力，A為特征長度，本文是指立管的最大迎風橫截面積。

不同脈動頻率下的阻力系數和升力系數的時間歷程曲線如圖6所示。

從圖6(a)可以看出，脈動頻率為固有頻率的2倍時，平均阻力系數最大且波動最為明顯。當脈動頻率低于此值時，脈動頻率越低，平均阻力系數越低。脈動頻率為固有頻率的3倍時，平均阻力系數開始下降，并出現負值，這是由于過高的脈動頻率導致立管在回流作用下產生明顯的回彈效應。不管脈動頻率怎么變化，其阻力系數的波動頻率與脈動頻率的變化趨勢是對應的。從圖6(b)可以看出，升力系數的幅值隨脈動頻率的增大而增大，當脈動頻率增大到固有頻率的2倍及其以上時，升力系數幅值基本不變，唯一區(qū)別是，不同脈動頻率達到穩(wěn)定所需時間的不同。綜上分析可以推斷，存在一種臨界脈動頻率，在此之前，阻力系數與升力系數不斷提高，一旦過了此臨界脈動頻率，阻力系數開始下降，但升力系數幅值卻開始保持穩(wěn)定。

圖6 不同脈動頻率下的水動力系數

3.1.2 振幅響應分析

振幅均方根(Root Mean Square，RMS)反映立管振動的偏離程度。在立管的軸向Z上，脈動頻率分別為3fn，2fn,fn，0.5fn時，不同位置的橫向振幅RMS值如圖7所示。

圖7 不同脈動頻率下橫向振幅RMS值

從圖7可以看出，在立管的軸向長度Z上，振幅均方根最大值出現在中間位置，越靠近兩端，其結果越小。如果流動不是均勻來流而是正弦脈動波動，則會降低振幅響應幅值。在4種脈動頻率中，脈動頻率f為2fn時的振幅結果最大，其次是3fn，fn，0.5fn。在立管軸向長度0.6 m，1.2 m，1.8 m，2.4 m以及3.0 m處分別設置測點1，2，3，4以及5，監(jiān)測橫向振動，橫向無量綱振幅用y/D表示，其中D為立管直徑。4種脈動頻率下，振幅變化曲線如圖8所示。由于測點2和測點4，測點1和測點5呈現對稱趨勢，圖8只給出立管在測點1、測點2以及測點3橫向振幅的時間歷程曲線。

圖8 不同脈動頻率下立管橫向振幅曲線

從圖8可以看出，圖8(a)和圖8(b)中的振幅曲線極其不規(guī)則，并且橫向振幅也較低，其最大振幅比是脈動頻率為0.5fn時的10倍；脈動頻率為2fn或3fn時，振幅可以達到立管直徑的0.6倍。

對不同的振幅曲線進行快速傅里葉變換(Fast Fourier Transform，FFT)分析振動頻率。不同脈動頻率下，橫向與流向振動頻率功率譜密度(Power Spectral Density，PSD)分析結果如圖9所示，左側是橫向振動頻率分析，右側是流向振動頻率分析。

圖9 不同脈動頻率下振幅FFT變換

從圖9(a)可以看出，橫向功率譜密度圖出現多峰現象，脈動頻率與一階固有頻率共存，但脈動頻率占絕大部分；對比流向頻率可以發(fā)現，流向振動主頻率并不是橫向振動頻率的2倍，這和均勻流的結果有差異，也說明分析脈動流條件的必要性。從圖9(b)可以看出，橫向振動頻率也出現多峰值現象，說明此時振動比較復雜且并未穩(wěn)定。對比圖9(c)和圖9(d)發(fā)現，兩者的橫向振動頻率都接近一階固有頻率，但圖9(d)的流向振動頻率要大于圖9(c)。綜上分析可以得出結論，當脈動頻率低于fn時，脈動頻率對橫向振動頻率的影響更為強烈；當脈動頻率大于等于2fn時，脈動頻率對橫向振動頻率的影響開始減弱，流向振動頻率與橫向振動頻率之比并不完全接近2，脈動頻率的增加導致流向的振動頻率一直在增加。

3.2 流場分析

文獻[5]指出，漩渦脫落在不同的振幅響應分支內具有不同的脫落模式，分別是2S，P+S，2P模式等。2S模式下，每過半個周期，從圓柱的兩側釋放1個漩渦；P+S模式下，每過半個周期，圓柱從一側釋放1個漩渦，另一側則釋放1對漩渦；2P模式下，每過半個周期，從圓柱兩側釋放1對漩渦，并且2個漩渦強度大小不一。圖10給出了不同脈動頻率下中間平面在耦合計算第8 s的渦量云圖。

圖10 不同脈動頻率下的渦量云圖

圖10(a)中，由于脈動頻率較小，兩側各存在1個漩渦，且并沒有交替脫落，為較弱的2S模式。圖10(b)為交替脫落的2S模式。從圖10(c)和圖10(d)可以看出，當脈動頻率大于或等于2fn時，圓柱的偏移距初始位置較大，在兩側分別有脫落的1對漩渦，為2P模式。結合3.1.2節(jié)對振幅的分析可知，脈動來流下尾渦脫落模式受到振幅的影響較大。

4 結束語

本文通過模擬脈動來流下的三維立管渦激振動響應，探究了脈動頻率對振動振幅及水動力系數的影響，為立管設計與安全壽命研究提供一定的參考。研究結果發(fā)現，脈動頻率可以降低渦激振動響應；流向振動頻率隨脈動頻率的增大而增大，說明不能忽略脈動流情況下的流向振動。由于計算資源有限，本文只是將海洋波動近似為正弦脈動，與實際的海洋流動有一定出入，后續(xù)將研究更為復雜來流形式下立管的渦激振動響應。