前言：尋找寫作靈感？中文期刊網用心挑選的波浪柱流致振動高保真數值模擬，希望能為您的閱讀和創作帶來靈感，歡迎大家閱讀并分享。
摘要波浪柱具有良好的減阻減振效果,在海洋工程和橋梁工程中應用廣泛。為了研究低雷諾數均勻流下細長柔性波浪柱的流致振動特性,選取了3種不同波高的截面直徑呈正弦變化的波浪柱,分別給予其駐波和行波的初始激勵,對其動力特性和尾流特征進行分析,并與同等條件下的光滑直圓柱進行了對比。研究結果表明,柔性波浪柱在駐波和行波初始激勵下的流致振動特性有很大差別:在駐波激勵下,當波高≥0.2D時,波浪柱振動得到了很好的控制,這主要歸功于波浪柱的特殊結構防止了剪切層與強度較高的渦結構發生摻混,增大了漩渦生成長度;而在行波初始激勵下,波浪柱的控制效果有限,局部雷諾數主導了結構響應、尾流模式和展向變化特征。
關鍵詞低雷諾數;柔性波浪柱;流致振動;駐波;行波;動力特性;尾流特征
波浪柱由于具有良好的減阻和抑制渦激振動的效果,在實際工程中得到越來越廣泛的應用。近年來,針對波浪柱的流致振動特性和減振機理,學者們展開了一系列的研究。Lam等[1-3]采用大渦模擬,對亞臨界雷諾數下波浪柱波長、波幅、入流角對減阻效果的影響進行了分析,發現在合適的波長和波幅條件下,波浪柱在亞臨界雷諾數下減阻效果依然優秀,同時波浪柱減阻效果與偏轉角成反比;Garba-ruk等[4]采用準三維方法將波浪柱的全三維流動結構建模轉化為一系列代表截面的二維流動特性問題,并對預測模型進行了穩定性分析,結果表明,穩定性分析結果與數值模擬結果吻合較好,為以后對波浪柱流場進行準三維模擬提供了可靠的模型框架;Zhang等[5]數值模擬研究了初始激勵下Re=5000時波浪柱的流致振動特性;給圓柱施加了一個橫向正弦形式的激勵,并且固定了激勵的振幅,結果表明,隨著初始激勵頻率的增加,升力系數的分布發生顯著變化,但初始激勵與波浪柱及其流場之間的相互作用機理還待進一步探明;安苗等[6]利用風洞試驗,對波浪柱在高雷諾數下的氣動性能進行了分析,結果表明在高雷諾數下波浪柱的減阻效果有限;趙桂欣等[7]使用了大渦模擬對不同波長波幅組合的有限長波浪柱進行了研究,發現在亞臨界雷諾數下大部分波長波幅組合的波浪柱具有明顯減阻效果。目前研究大都以固定或者彈性支承的剛性波浪柱為對象,然而在實際工程中,柔性波浪柱具有更廣闊的應用前景。本次研究在Re=100的均勻來流下,以截面直徑呈正弦規律變化的柔性波浪柱為研究對象,對其在行波和駐波初始激勵下的尾流特性、動力響應以及背后的耦合效應機制進行了探究。
1計算模型的建立及網格劃分
1.1數值計算模型
均勻流中的波浪柱模型的示意圖如圖1所示。該圓柱直徑沿展向變化由下式確定。其中:Dz為展向相應位置的圓截面直徑;D代表波浪柱的平均直徑,由于波浪柱直徑在展向呈正弦變化,平均直徑D=(Dmax+Dmin)/2,定義波浪柱直徑最大截面所在的展向位置為幾何節點,波浪柱直徑最小截面所在的展向位置為幾何鞍點;A代表波狀表面的波高,共設置了3種波高,分別為0.1D、0.2D、0.3D,展向波長設置為λ=4πD。為了方便對比,在相同來流條件下對直徑為D的光滑圓柱也進行了數值模擬。此處需要說明的是:以下討論中所說的圓柱,只要沒有具體指出是波浪柱還是光滑圓柱,則兩種類型圓柱都包含在內。
1.2數值計算方法
采用Serson等[8]詳細介紹的坐標映射方法并結合式(1)來構建波浪柱。設定均勻流的來流速度為U∞,與之相應的雷諾數定義為Re=U∞D/ν(ν為流體的運動粘度系數),將雷諾數設定為100,使流動保持層流狀態。展向的計算域設置為4πD,與波浪柱一個波長的大小相同。參考Prasanth等[9]的設置,將圓柱的質量比設置為10,并給出圓柱質量比的計算公式m*=4ρc/ρfπD2(ρc和ρf分別是單位長度的圓柱質量和流體密度)。限制圓柱在順流向的運動,只考慮圓柱在橫流向的運動,因為順流向運動的振幅通常僅為橫流向振幅的10%~15%[10]。采用張力梁模型來表征圓柱的動力學響應,張力梁模型的控制方程為ρc?2ξy?t2+b?ξy?t=T?2ξy?z2+F(z,t),(2)其中:b為圓柱的結構阻尼系數,設置為0,允許產生最大的橫向振幅振蕩;T、ξy分別表示作用于兩端的拉力和橫向振幅。由于我們只對鎖定狀態感興趣,結構固定頻率被設置為2πSt,以保證其與漩渦脫落頻率相近,其中St是在相同雷諾數下均勻流流過光滑圓柱的斯特羅哈爾數(Strouhalnumber)。然后,根據Newman等[11]計算出相應的圓柱張力T=m*L2S2t,F(z,t)是沿展向作用在圓柱上的力,由壓力和粘滯應力項的積分計算。三維流動由開源軟件Nektar++中的不可壓縮Navier-Stokes求解器控制,該軟件基于傅立葉譜元法[12]。其中:u是速度;p是壓強。Nektar++中采用的流固耦合方法是(x,y)平面上的Jacobi-Galerkin公式和z方向上的傅立葉展開的混合格式。將(x,y)平面的二維空間域離散為四邊形有限元單元,在每個四邊形有限元單元內,再采用高斯-洛巴托-勒讓德(GLL)插值點關聯的高階拉格朗日多項式作為形狀和權重函數進行離散。關于空間離散化算法的細節可參考文獻[12-13]。時間離散化基于高階剛性穩定分裂格式,允許原始變量在每一個時間步長內獨立處理。類似應用于流體流動控制方程的傅立葉變換,對結構的運動施加了一個在全長范圍內周期性的假設。然后將運動變量和力表示為復傅立葉級數,將拉伸梁模型的偏微分方程解耦為一組常微分方程,用二階Newmark-β方法求解。
1.3網格設置
網格設置如圖2所示。圖2(a)為二維平面基礎網格的示意圖,圓柱被放置在一個C形計算域中,C形計算域由圓柱上游半徑為20D的半圓形計算域和圓柱下游的矩形計算域組成,矩形計算域在流向長度為30D,在橫流向長度為40D。在(x,y)平面上采用了877個基礎網格,同時使用P=6對網格進行差分,如圖2(b)所示,每個基礎網格被分成了(6-1)×(6-1)個高階網格。此外,根據Newman等[11]和Zhu等[14]的研究,將計算域沿展向劃分了64個傅立葉面。計算域的邊界條件設置如下:入口邊界為均勻來流條件,來流速度為(U∞,0,0);在出口邊界上采用紐曼邊界條件,令出口邊界上以及p=0;在頂部和底部邊界上指定(u,v,w)=(U∞,0,0)的遠場邊界條件;在圓柱固壁采用無滑移和不可滲透邊界條件;通過假設無限長波浪柱的展向周期性,在波浪柱的兩端實施了只允許橫向運動的自由邊界條件[11]。此外,Lam等[1]證明了只有單個展向波長的周期邊界條件對于計算波浪柱繞流足夠精確。
1.4網格驗證
為了確保模擬使用的網格具有足夠的精度,并且保證在現有的網格設置下繼續增大網格密度對模擬結果影響較小,需要對測試算例進行網格獨立性驗證。采用在(x,y)平面改變多項式的階數P來改變平面網格的疏密,在展向上通過改變傅立葉面數目來改變z方向網格的疏密。不同網格數目下的測試算例結果對比見表1。測試算例均是在雷諾數Re=100條件均勻來流下的光滑圓柱,所有測試算例光滑圓柱的質量比m*=10。由表1可得,在不同的網格數目下,模擬結果的差異非常小。為了兼顧模擬的精度與計算速度,在主要模擬中使用的網格精度設置為P6Nz64。
2駐波初始激勵下模擬結果分析
2.1動力響應與振幅
橫向振幅(ξy)、阻力系數(Cd)和升力系數(Cl)沿展向的均方根值分布如圖3所示。阻力系數Cd=2Fd/ρfU2∞D,升力系數Cl=2Fl/ρfU2∞D,其中Fd和Fl分別代表圓柱受到的流向和橫向的截面力的分量。為了能夠更加精確地捕捉到升阻力和振幅隨波浪柱波高變化的規律,在此處添加了A=0.15D和A=0.25D兩種工況。由圖3(a)可知,隨著波浪柱波高A的增加,波浪柱的橫向振幅逐漸減小,尤其是在A≥0.2D的工況下,波浪柱的振動幾乎完全被抑制。圖3(b)中的阻力系數均方根分布也顯示出同樣的規律,并且在A≥0.2D的3種工況下,阻力系數相比于光滑圓柱有了明顯降低。而圖3(c)中升力系數的均方根分布顯示,在A<0.2D的兩種工況下,升力系數相比于光滑圓柱反而有所升高,繼續增加波浪柱波高至A≥0.2D后,升力系數才會顯著下降。為了進一步闡明波浪柱波高與振動抑制之間的關系,對橫向振幅及升阻力系數的時空演變規律進行了分析,此處取光滑圓柱、A=0.1D和A=0.2D3種代表工況進行分析。光滑圓柱和波浪柱橫向振幅(ξy)的時空演變云圖見圖4。由圖4可以看到,A=0.1D工況的節點和反節點位置與光滑圓柱相同,節點位于z/D=0、2π、4π,反節點位于z/D=π、3π。同時在任何無量綱時間(tU∞/D)下,A=0.2D工況下波浪柱的橫向振幅均較光滑圓柱大幅降低,而A=0.1D工況下波浪柱抑制振動的效果較差。阻力系數(Cd)的時空演變云圖見圖5。由圖5可以看到,波高A=0.1D工況下阻力系數的時空分布和大小與光滑圓柱的基本相同;而A=0.2D工況下的阻力系數時空分布則與光滑圓柱完全不同,在阻力系數分布中看不到明顯的駐波效應,在z/D=0、4π時阻力系數并未取得最大值,并且在不同無量綱時間處阻力系數最小的點是沿展向變化的,同時A=0.2D波浪柱的最大阻力系數僅有光滑圓柱的一半,說明A=0.2D波浪柱具有良好的減阻效果。波浪柱升力系數(Cl)的時空演變云圖見圖6。由圖6可以看到,A=0.1D工況下波浪柱在展向π<z/D<3π范圍內的升力系數較光滑圓柱有所下降,而在其他展向位置升力系數較光滑圓柱有所升高,這是由于局部雷諾數的增加引起的。而A=0.2D的波浪柱僅在展向的兩側有較高的升力系數,其他位置升力系數非常小。根據波浪柱在駐波初始激勵下對橫向振幅和升阻力的抑制效果,可以將波浪柱分成控制有效(A=0.2D、A=0.3D)和控制失效(A=0.1D)波浪柱兩類。
2.2尾流特性
為了進一步了解不同波高導致波浪柱不同動態響應的內在機制,需對波浪柱和光滑圓柱近尾跡的三維渦模式展開研究。3種代表工況下展向渦量瞬時等值面的透視圖和橫流向俯視圖見圖7。從圖7可以看出,駐波初始激勵下光滑圓柱的展向渦量中存在明顯的交織結構,這與Newman等[11]所觀察到的相同。在A=0.1D工況下也可以大概看到交織的渦結構,并且在幾何節點(波浪柱截面最大展向位置,處在波浪柱的兩端)的展向渦強度較大,而幾何鞍點(波浪柱截面最小展向位置,處在波浪柱的跨中)處的展向渦強度較小,強度較強的渦會在兩端產生較高的截面升力系數,也解釋了圖6(b)的現象。在A=0.2D工況下交織渦結構完全消失,彎曲渦管以交錯的方式在波浪柱上下表面分別脫離,并且在尾流中迅速消失。更重要的是,與光滑圓柱相比,A=0.2D波浪柱剪切層的穩定性較差,并且向上卷動較弱,形成渦旋。渦旋在下游進一步分離,使得渦形成區長度進一步增大。Lin等[15]用數值方法研究了亞臨界流動區域中具有相對較大展向波長的剛性波浪圓柱周圍的流動,也發現了類似的結果。
3行波初始激勵下模擬結果分析
Newman等[11]的研究表明,對于無量綱長度超過4π的光滑直圓柱而言,結構的振動模式最終發展為行波響應。然而,對于長度為4π的光滑直圓柱,經過行波-駐波響應的競爭,最終結構的振動模式發展為駐波形式,這表明了結構的長度是影響結構動力響應模式的關鍵因素之一。對波浪柱在行波初始激勵下的動力響應特性和尾流特征展開研究,賦予結構的行波初始激勵為其中:a是振幅;ω是振動頻率,ω=2πcLz,c=T/ρ表示相速度。
3.1動力響應與振幅
行波初始擾動下光滑圓柱的橫向振幅和升阻力系數云圖見圖8。由圖8(a)所示的光滑圓柱橫向振幅可以看出:隨著時間的演化,結構響應呈現出明顯的行波響應特征,即響應由一端傳遞至另外一端,任意處的結構響應最大幅值一致且同一時刻結構展向不同位置處存在連續相位差,這與Newman等[11]的模擬結果相同,進一步證明了本模擬的準確性。行波初始擾動下3類典型工況結構的橫向振幅時空演化云圖見圖9?？傮w而言,不同工況下結構的最大振幅響應基本一致,相較于光滑圓柱的橫向振幅均有所增加。當A=0.1D時,同樣出現了明顯的行波響應,但結構兩端(z/D=0,4π)和中間(z/D=2π)出現了相對較大振幅響應。隨著波高A的繼續增加,較大振幅響應位置逐漸集中至z/D=π,3π兩位置處,如圖9(b)和圖9(c)所示。該現象呈現出了駐波-行波兩類響應疊加的特征。行波初始擾動下波浪柱的阻/升力系數時空分布見圖10。由圖10可以看到,A=0.1D工況下作用于結構上的流體動力同樣不再是單一的行波響應,而是呈現了明顯的駐波-行波疊加響應。New-man等[11]發現柔性管道在無量綱時間為500~570時結構響應由駐波演化為行波,Facchinetti等[16]的研究同樣表明,行波響應是該類結構的穩定狀態。對于本次研究而言,初始激勵為行波激勵,但行波響應并不是穩定狀態,這與結構的波浪形有關,即結構為一種固定形式(可當作駐波),與行波疊加,從而呈現這一特殊現象。此外,A=0.1D工況下阻力和升力幅值與光滑圓柱相比均有所提升。A=0.2D工況除了呈現與A=0.1D工況一致的駐波-行波疊加響應之外,結構兩端(z/D=0,4π)阻升力明顯大于展向其他位置。這與結構的形式有關,在該工況下結構兩端的直徑為1.4D,跨中(z/D=2π)為0.6D,根據局部雷諾數(Reloc=ρU∞dloc/μ)計算,升阻力與光滑圓柱相比產生了明顯變化(最大阻升力提升40%)。A=0.3D工況除了呈現駐波-行波疊加響應外,流體動力響應的最大特征是失去了對稱性,而該非對稱特征隨時間演化并未消失,而是成為穩定狀態。該現象的出現是由于研究中結構本身為超柔性體,結構自激動力響應呈現強非線性,與展向變化的結構形式以及流體的耦合相互作用,使得整個結構系統呈現出強烈的非線性特征,從而動力響應呈現了分叉特征,當結構響應偏向一側時,響應會朝著該方向繼續發展。此外,由于波浪柱的各工況均呈現出駐波-行波疊加響應,以工況A=0.2D為例做出解釋。A=0.2D工況無量綱時間tU∞/D=0~300的結構橫向振幅時空演化圖見圖11。由圖11(a)可見,在初始行波激勵下,結構呈現行波特征,在tU∞/D=40時行波特征消失,隨后在tU∞/D=90時再次出現行波響應;由圖11(b)可見,隨著時間的演化振幅不斷增大(超過0.6D),且行波響應特征增強;由圖11(c)可見,至tU∞/D=240~260,變為駐波-行波響應疊加穩定狀態。結合圖9(b),駐波-行波響應競爭結束后,結構橫向振幅穩定至0.54D。因此,認為結構橫流向呈現的駐波-行波響應是兩種振動響應模式相互競爭、相互結合的結果。該現象通過渦結構開展進一步討論。
3.2尾流特性
行波初始擾動下4類典型工況展向渦量等值面圖(ω=±0.5)見圖12,圖中黃、藍色分別表示展向正負值渦量。由圖12(a)可以看出,結構展向呈現了規則的行波結構,這與前述以及Newman等[11]的研究結果一致,漩渦從結構兩側交替脫落,與管道結構呈現固定角度。此外,隨著波浪柱波高A值的增加,駐波響應特征出現。由圖12(b)可見,工況A=0.1D時,漩渦不再呈固定角度從管道結構上脫落,位于管道兩端的渦尺度相對更大,至A=0.2D時(見圖12(c)),兩端渦強度顯著大于中間部位。對于工況A=0.3D而言,行波或駐波特征通過渦量等值圖難以識別,結合動力響應分析結構認為是駐波-行波疊加響應的結果。為更加清楚地表明圓柱不同展向位置后側尾流形式的變化規律,提取了各工況的兩個展向位置(z/D=π,2π)的展向渦量云圖見圖13。對于光滑圓柱而言(見圖13(a)和圖13(b)),兩展向位置僅有相位差,尾流形式一致;當A=0.1D時(見圖13(c)和圖13(d)),相對于z/D=π位置處的尾流,z/D=2π處的尾流具有更寬的橫流向渦核間距,該現象意味著位于z/D=2π的結構振幅更大,這與圖8(b)的結果一致;當A=0.2D時(見圖13(e)和圖13(f)),z/D=π位置處的橫流向渦核間距相對更大,認為這與顯著增強的駐波響應有關;對于工況A=0.3D而言,z/D=π位置處的尾流呈現了轉捩特征,而z/D=2π處的尾流出現了流動不分離的特征,認為這與局部雷諾數有關,在該工況下結構節點處的局部雷諾數高達160(對于圓柱繞流而言,層流轉捩發生在雷諾數180左右),因此出現了早期轉捩現象。而該工況時最小局部雷諾數為60(對于圓柱繞流而言,流動失穩發生在雷諾數50左右),因此出現了位于z/D=2π處下游區的穩定流動特征。
4結論
在駐波初始激勵下合理地改變波高,可有效減小波浪柱的振動幅值和作用在其上的阻力。根據波浪柱在駐波初始激勵下對橫向振幅和升阻力的抑制效果,將波浪柱分成控制有效(A=0.2D、A=0.3D)和控制失效(A=0.1D)波浪柱兩類。對于控制失效波浪柱,其動力特性與光滑圓柱基本一致,而對于控制有效波浪柱,其動力響應和受到的流體動力均顯著降低。對于駐波響應具有良好控制效果的波浪柔性柱并未對行波響應起到良好控制效果,即截面尺寸沿展向正弦變化的波浪形圓柱對行波響應的控制效果有限,在行波響應下局部雷諾數主導了結構響應、尾流模式的展向變化特征;波浪形細長柔性柱在行波初始擾動下的動力響應與波高顯著相關,呈現兩類不同振動形式:當A?0.1D時,為典型的行波振動模式;當A>0.1D時,經歷駐波、行波振動競爭后,呈現出了駐波-行波疊加振動模式。
參考文獻:
[1]LamK,LinYF.LareEdgdySimulationofFowlAroundWavyCylindersataSbcriticalRueynoldsNumber[J].InternationalJur-onalofHatandFuidFowell,2008,29(4):1071-1088.
[2]LamK,LinYF.EffectsofWavelengthandAmplitudeofaWavyCylinderinCoss-flowatLrowReynoldsNumbers[J].JournalofFuidMlechanics,2009,620:195-220.
[3]LamK,LinYF,ZouL,etal.InvestigationofTurbulentFowlPastaYawedWavyCylinder[J].JournalofFuidsandSruc-lttures,2010,26(7-8):1078-1097.
[4]GarbarukA,CrouchJD.Quasi-threeDmensionalAinalysisofGlobalIstabilitiesn:OnsetofVortexSeddinhgBehindaWavyCylinder[J].JournalofFuidMlechanics,2011,677:572-588.
[5]ZhangK,KatsuchiH,ZhouD,etal.NumericalSudtyofFowlPastaTansverselryOscillatingWavyCylinderatRe=5000[J].OceanEngineering,2018,169:539-550.
[6]安苗,劉慶寬,孫一飛,等.波浪形圓柱的氣動力特性試驗研究[J].振動與沖擊,2021,40(8):209-215.
[7]趙桂欣,桂洪斌,王曉聰.有限長波浪形圓柱繞流數值模擬[J].哈爾濱工業大學學報,2021,53(6):163-170.
[8]SersonD,MeneghiniJR,SherwinSJ.Velocity-correctionS-chemesfortheIcomnpressibleNavier-StokesEquationsinGeneralCordinateoSystems[J].JournalofComputationalPhysics,2016,316:243-254.
[9]PrasanthTK,MittalSVortex-inducedVbrationsofaCrcu-.iilarClinderatLyowReynoldsNumbers[J].JournalofFuidlMechanics,2008,594:463-491.
[10]EvangelinosC,KarniadakisGE.DynamicsandFowlStruc-turesintheTurbulentWakeofRigidandFexibleClinderslySubjecttoVortex-InducedVbrationsi[J].JournalofFuidMle-chanics,1999,400:91-124.
[11]NewmanDJ,KarniadakisGE.ADirectNumericalSmulationiStudyofFowlPastaFeelryVibratingCable[J].JournalofFluidMechanics,1997,344:95-136.
[12]KarniadakisG,SherwinSS.pectral/hpElementMethodsforComputationalFuidDlynamics[M].OxfordUiversitnyPress,2013.
[13]KarniadakisGE,IsraeliM,OrszagSA.High-orderSlittinpgMethodsfortheIcomnpressibleNavier-StokesEquations[J].JournalofComputationalPhysics,1991,97(2):414-443.
[14]ZhuHongbo,WangR,BaoY,etal.FlowOveraSymmetri-callyCurvedCrcularClinderwiththeFeeSreamiyrtParalleltoThePaneofCrvatureatLluowReynoldsNumber[J].Jour-nalofFuidsandSructureslt,2019,87:23-38.
[15]LinYF,BaiHL,AlamMM,etal.EffectsofLrageSanwisepWavelengthontheWakeofaSnusoidalWiavyCylinder[J].JournalofFuidsandSructureslt,2016,61:392-409.
[16]FacchinettiML,deLnagreE,BiolleyF.Vortex-inducedTav-rellingWavesAonlgaCblea[J].EuropeanJurnalofMoechan-ics-B/Fluids,2003,23(1):199-208.
作者:龐璽源 朱宏博 張春曉 包艷 單位:上海交通大學船舶海洋與建筑工程學院 內蒙古路橋集團有限責任公司