幹貨 | OpenFOAM渦激振動仿真詳解

渦激振動的定義

當雷諾數超過40時,流體在圓柱體上流動會發生渦脫落現象,從而在圓柱體上施加周期性的作用力。如果圓柱體受到彈性支撐,則該力的橫向分量就會導致圓柱體的橫向振動,這一現象稱為渦激振動現象。

彈性安裝圓柱形結構的渦激振動現象在許多工程應用中經常碰到,如換熱器管、海洋立管、海上平臺結構、海底管道、海上風力渦輪機柱基以及輸電線路等等。當渦脫落的頻率接近系統的固有頻率時,系統會發生大幅度橫向振動,這就是鎖定或同步現象。鎖定現象發生時,流體在圓柱結構上施加大幅值周期性流體力,會導致這些柔性結構發生疲勞失效,造成一定的風險隱患。

近幾十年,科學傢們對渦激振動現象展開瞭大量的研究。在進行設計時,難免要與已有的文獻研究結果進行比較。本文將從筆者自身經驗出發,分享如何使用OpenFOAM模擬渦激振動,以及如何取得與文獻相符合的結果的一些個人經驗,供初學者參考。

OpenFOAM具體設置

OpenFOAM版本為esi的v2006版本,使用其動網格功能來實現渦激振動計算。本文所用到的算例在gitee上可供下載:

https://gitee.com/chen-jingle219/vortex_induced_vibration/blob/master/README.md

  • 流場網格

使用OpenFOAM自帶的動網格技術(dynamicMesh),用到的網格是在ICEM畫的結構網格。如下圖所示,在圓柱附近設置瞭一個很大范圍的圈,包裹著圓柱。正常的圓柱繞流算例不會那麼復雜,但為瞭應對渦激振動在特定約化速度Ur下幅值較大的情況,特地留瞭充足的網格變形范圍。如果處理不當,很容易造成局部網格變得極端小,從而計算發散。這裡,圓柱直徑D=0.0016m(沿用一些文獻的尺寸),網格范圍為30D X 50D。

  • 0:邊界條件以及初始條件

OpenFOAM邊界條件以及初始條件在case/0文件夾中面統一設置,可以在我的gitee中看到。值得一提的是,要使用動網格,需要加pointDisplacement 文件。然後在U文件中,圓柱表面的邊界類型由 fixedValue 改為 movingWallVelocity:

其它邊界條件以及初始條件與圓柱繞流算例基本相同。

  • constant: 動網格設置

constant裡面除瞭接下來著重講的dynamicMeshDict以外,還有常見的transportPropertiesDict和turbulencePropertiesDict兩個Dict文件。將運動粘度nu設置為水的粘度1e-6,湍流模型設置為層流。對於dynamicMeshDict,主要講解需要改動的部分:

(1)運動邊界以及網格變動范圍

如上圖,WALL表示圓柱面對應的邊界,如果自己的算例是別的名字,那麼這裡就得改。innerDistance代表以運動邊界中心為圓心,網格變形的最小距離;outerDistance代表網格變形最大的距離。這兩者之間就是一個網格變形環。這兩個距離需要根據自己網格的實際尺寸來設置。如果想要全局變形,那麼可以將innerDistance設置為0,outerDistance設置成超過網格最大尺寸的一個數。

(2)質量、慣量以及密度

質量由質量比、圓柱長度、圓柱直徑和流體密度計算得到,公式在圖裡。質量比可以理解為固體的質量和其排開水的質量的比值,代表其在水裡有多“重”。momentOfInertia為固體在這個慣性系裡面三個方向的慣量,因為本算例不涉及轉動,所以這個無關緊要,但不能設為0,否則會報錯。密度設置為水的密度1000kg/m3。

(3)約束constraints

Constraints可以有約束的意思,且約束的程度更強,有限定的意思。在這裡,constraints就是給圓柱的運動方向進行設定。如圖所示,yline表示其隻能在y方向上運動,zAxis則是不可繞z軸轉動。這樣,圓柱隻能沿y軸進行平移。其它的constraints類型也在圖中給出:fixedPoint可以把圓柱固定在某一點處;zPlane可以讓圓柱在z平面移動,也就是x、y兩個方向都可以運動。

(4)約束restraints

Restraints約束的程度較弱,在這裡,是給圓柱運動施加一定的阻力。Restraints有很多種,適合渦激振動的就是圖中所示的linearSpring,也就是一個彈簧-阻尼系統。這個彈簧有定點(anchor)和動點(refAttachmentPt)兩個end points。計算時,anchor始終固定,而refAttachmentPt則是動點的一個參考點,用來計算彈簧長度的變化量的。Stiffness、damping和restLength顧名思義分別是剛度、阻尼和彈簧原長,根據需要進行改動。至於應該怎麼給?這些都是根據約化速度Ur來定的,在任何一篇渦激振動的文章裡頭都有對應的公式。

(5) system/fvSolution

fvScheme基本與圓柱繞流的一致。fvSolution則需要加上pcorr的矩陣求解器:

其它基本與pimpleFoam算例裡的一致。

運行-後處理-使用腳本提取圓柱位移信息

  • 運行

按自己機器的核數並行計算即可。需要註意,提取圓柱運動的位移,需要在運行時加上> log.pimpleFoam,這樣,輸出的信息就能存在log文件裡面。

  • 後處理

使用cat命令來提取圓柱的位移信息:

最後,圓柱的位移曲線就會存在plotfile.dat裡面,就可以用matlab處理或者用tecplot進行查看。

計算結果

Ur=3~10共8組數據,與文獻結果符合較好,說明之前的設置是比較好的。每個人算的結果都可能會有一點不一樣,這與網格質量、數值格式、矩陣求解器設置等等相關。

分享一些經驗

  • 盡量用現實的變量進行計算而非無量綱數。寫求解器計算時都會根據無量綱數設置一些參數,如將U、D設為1,那麼nu就為Re的倒數;但在OpenFOAM裡不推薦這樣做。一方面,不可壓求解器pimpleFoam等壓力方程的p實際上是p/rho,在量綱上就會比較亂;另一方面,設計不同算例時,換算也容易出錯。建議物理常數如nu、m、rho等都按照實際對應的物理量來設置,再根據公式來換算其他物理量。最後再對結果進行無量綱化,減小出錯率。
  • 直接改剛度,一個算例算N遍。在計算渦激振動時,一般都是固定瞭Re,改變圓柱剛度來改變Ur的。這意味著,流場的情況基本不變,隻需要改變圓柱的剛度即可。在計算不同Ur值時,筆者是用瞭1個算例來算N次,期間不斷修改dynamicMeshDict裡面的stiffness。本次的算例繼承瞭上一個算例的流場,節省瞭大量拉渦-形成卡門渦街-相對穩定的運算時間,提高實驗效率。

赞(0)