平動與轉動受迫諧振圓柱的水動力特性分析

劉名名， 唐國強， 呂 林， 滕 斌

(大連理工大學 海岸和近海工程國家重點實驗室，遼寧 大連 116024)

平動與轉動受迫諧振圓柱的水動力特性分析

劉名名， 唐國強， 呂 林， 滕 斌

(大連理工大學 海岸和近海工程國家重點實驗室，遼寧 大連 116024)

通過求解二維不可壓縮黏性流體的Navier-Stokes方程，對低雷諾數下圓柱同時具有周期性旋轉和橫流向振動的流固耦合問題開展了數值分析研究。重點考察了不同強迫振動頻率(f*([0.5, 2.0])以及轉動與平動之間的相位差(φ=±180°，±120°，±60°及0°)下，圓柱結構的受力特性、鎖定模式及尾跡渦流場演化特征。數值分析結果表明，在周期性轉動和平動的聯合激勵下，流動過程除存在常規的基本鎖定(Primary Lock-in)現象外，還存在非線性準周期鎖定(Quasi-periodic Lock-in)現象；旋轉運動與橫流振動之間的相位差φ對鎖定形式及鎖定的頻率區間有重要影響，對稱的正負相位差可導致高度一致或截然不同的鎖定形式；兩自由度諧振圓柱的水動力系數在鎖定條件下會發生明顯跳躍，總體上隨強迫振動頻率的增大而增大；對于特殊的準周期鎖定情況，尾跡區內存在渦旋結構的多級次非線性分叉行為，尾渦演化過程遵循準周期特性。

圓柱；旋轉振動；橫流向振動；準周期鎖定；尾渦模式

當流體經過非流線型結構時，通常會在其尾流區內產生交替的渦旋脫落，導致作用在結構上的流體作用力表現出周期性特征，這是誘發結構發生渦激振動以及馳振等非線性動力響應(甚至是疲勞破壞)的主要原因之一。圓柱結構在水利工程、海洋工程、化學工業以及航空航天等領域都有十分廣泛的應用。理解和認識圓柱結構與流體的相互作用，特別是多自由度振動圓柱結構的受力以及相應的流動特性，無論是對解決工程實際問題(包括流動控制)，還是認識自然規律，都具有重要的科學研究價值。

對結構施加一定運動規律的強迫振動是研究流固耦合問題的重要方法，這樣允許針對特定的運動形式，揭示結構的受力、流動特性以及流固耦合作用的機理。在以往的研究中，根據圓柱結構受迫振動方式的不同，可劃分為平動(包括橫流向及順流向振動)和轉動(包括單向和雙向)兩種基本形式。對于只具有橫流向單自由度受迫簡諧振動的問題，Williamson等[1]通過物理模型試驗，總結了不同的尾渦脫落模式(如2S、2P及P+S等)及其與圓柱振動幅值和振動頻率之間的依賴關系。在后來的研究工作中，很多學者關注了鎖定現象的發生，以及鎖定區間與圓柱受力及尾渦脫落模式之間的關系[2-6]。橫流向受迫振動的研究成果目前已經成為細長柔性結構渦激振動經驗預報模型的基礎[7-8]。對于圓柱發生順流向強迫振動的問題，人們同樣也針對不同振幅和頻率下的結構受力、鎖定特性以及尾跡模式等開展了大量研究，并取得了豐碩的研究成果[9-12]，這對于認識水流和波浪聯合作用下的結構受力和穩定性等都具有直接的借鑒意義。相比較而言，針對圓柱的受迫旋轉運動問題，相關的研究工作要少一些。通過強迫圓柱的旋轉運動，可改變邊界層內的流動分離特性，從而達到流動控制的目的。Du等[13]發現通過旋轉的方法可以抑制圓柱渦激振動。Bai等[14-15]對旋轉圓柱的流固耦合問題開展了數值模擬研究，揭示了圓柱旋轉速率對圓柱受力及流場特性的影響作用。Lu等[16]則提出了通過圓柱的主動性雙向旋轉來減小圓柱受力的流動控制方法。此外，學者們針對均勻流中圓柱發生受迫雙向旋轉的問題，也給出了圓柱受力和鎖定區間等對旋轉幅值和受迫振動頻率的依賴關系[17-18]。

上述的這些研究工作主要是針對具有單自由度的圓柱受迫振動問題開展的，而對于圓柱在水流中同時具有多個自由度強迫運動的研究工作則相對較少。Blackburn等[19]在橫流向受迫振動圓柱的基礎上引入圓柱旋轉運動，開展了流動控制研究。Nazarinia等通過模型試驗，研究了單向均勻流中，圓柱同時具有橫流向振動與轉動時的受力特性，并觀察到了非線性準周期鎖定的現象。但是，由于受到試驗條件和觀測方法的限制，Nazarinia等[20]沒有詳細給出發生準周期鎖定時的流場演化特征，同時也沒有進一步研究圓柱的振動頻率對其受力和鎖定特性的影響作用。因此，針對同時具有平動和轉動的多自由度圓柱受迫振動問題的研究還有待進一步加強，這一物理過程的背后具有豐富的非線性流固耦合現象，理解和認識其中的物理機理，也是今后開展更為深入的流動控制研究的基礎。本文將通過數值模擬，對同時具有橫流向振動和旋轉振動圓柱的水動力特性開展分析，重點考察圓柱受迫振動頻率以及橫流振動與旋轉振動之間的相位差對圓柱結構受力、鎖定特性以及渦流場演化特征的影響作用。

1 數值模型

1.1 流體運動控制方程

描述均勻不可壓縮黏性牛頓流體運動規律的基本控制方程為連續性方程和Navier-Stokes方程。在任意拉格朗日-歐拉(ALE)觀點下，可寫成如下的無量綱形式

(1)

(2)

式中：xi表示笛卡爾坐標(對于二維問題，規定x1=x，x2=y)；ui表示流體速度在xi方向的分量；p為壓力；t表示時間；Re=UD/υ為雷諾數(υ為流體運動學黏性系數，U為均勻來流速度，D為圓柱直徑)；cj表示ALE觀點下的j方向的網格運動速度分量。

1.2 圓柱強迫運動方程及網格更新方法

設圓柱同時發生周期性的橫流向振動及繞軸心的旋轉運動，其無量綱的平動位移y(t)及角位移θ(t)分別為

(3)

(4)

式中：φ2和φ1分別為旋轉運動與橫流向振動的初始相位，φ=φ2-φ1表示旋轉運動與橫流向振動之間的相位差；f*=fe/fs為無量綱的強迫運動頻率，其中fe表示圓柱的振動頻率，fs為相同雷諾數下固定圓柱的渦脫落頻率；A*=A0/D(A0為圓柱平動位移幅值)和θ*=θ0U/2πfe(θ0為圓柱角位移幅值)分別表示無量綱的平動幅值與轉動幅值。

流體與結構之間的相互作用通過動網格方法進行模擬，相應的網格坐標更新，通過求解如下的Laplace方程來獲得

(5)

式中：Δ=?2/?x2+?2/?y2為Laplace算子；q1和q2分別表示網格發生運動后對應的x和y坐標。邊界條件為已知的圓柱表面計算節點坐標。

1.3 數值求解方法

采用三步有限元方法對流動動量方程進行數值離散，具體形式如下[21]

(6)

(7)

(8)

由于網格更新是在一個完整時步完成后進行的，因此，式(6)～(8)中的網格運動速度cj均采用了n時刻的運動速度。數值計算結果表明，上述的線性化處理方法能夠在保證數值精度的前提下，獲得較高的計算效率。此外，由式(8)可知，要獲得n+1時刻的速度場uin+1，首先需要求得n+1時刻的壓力pn+1。通過對式(8)兩端同時取散度，并應用n+1時刻的連續性方程?uin+1/?xi= 0，可得到如下的壓力泊松方程：

(9)

本文所采用的三步有限元方法是一種高階迎風Taylor-Galerkin格式，并且流速和壓力變量通過投影方法進行解耦，因此在利用標準有限元方法進行空間離散時，速度和壓力可采用同階插值函數。同時，三步有限元方法滿足統一的CFL(Courant-Friedrichs-Lewy)穩定性條件，本文采用如下的動態時間步長：

(10)

式中：Se為網格面積；ue為網格中心點速度；Min(*)表示對所有計算網格取最小值；β=0.2為安全系數。對圓柱受到的壓力和黏性剪切應力沿著圓柱表面進行積分，將合力分別投影到x,y方向上并通過除以0.5ρDU2對其進行無因次化處理，可以得到單位長度的圓柱受到的升力系數和拖曳力系數分別為[22-23]

(11)

(12)

式中：p∞為無窮遠處的參考壓力。

2 計算模型及驗證

2.1 計算域及定解條件

圖1給出了圓柱發生轉動和橫流向兩自由度受迫振動時的示意圖及邊界條件。坐標原點與圓柱圓心重合，整個流體計算域的范圍為90D×145D，圓柱距離來流端60D，距離出口端85D，距離上下邊界分別為45D。

數值模型的邊界條件如下：在入口處指定無量綱速度u=1，垂向速度v=0；側壁采用對稱邊界條件?u/?y=0，v=0，?p/?y=0；出口處指定自由出流邊界條件?u/?x=0，?v/?x=0及相對壓力p=0。圓柱表面施加不可滑移邊界條件：u=-rωsinθ，v=dy/dt+rωcosθ及?p/?n=0，其中，ω=dθ/dt為旋轉角速度。在初始時刻，流場中的速度及壓力分布均設為零。

圖1 計算域及邊界條件

2.2 網格收斂性驗證

表1 網格收斂性驗證結果

2.3 數值模型驗證

本文首先針對圓柱僅發生橫流向受迫振動的問題開展了數值模擬，并與文獻[6]的數值結果進行對比驗證。計算中，Re=185，A*=0.2，f*=0.8、0.9、1.0、1.1、1.12和1.2。圖2給出了本文計算所得流體力系數與文獻[6]結果的對比情況。從圖2可以看出，本文的數值結果與已發表論文結果吻合良好，說明本文所建立的數值分析模型能夠正確模擬圓柱橫向受迫振動流固耦合問題。

為進一步說明本文數值模型的可靠性，對Re=100，f*=2.423，θ*=1.114條件下，受迫旋轉圓柱在均勻流場中的繞流問題開展了數值模擬。圖3給出了本文計算得到的升力系數時變曲線與Choi等[18]數值結果的對比。從圖3可以看出，二者具有良好的一致性，進一步說明了本文模型的正確性。

圖2 流體力系數隨圓柱橫流向受迫振動頻率的變化

圖3 受迫旋轉圓柱的升力系數時間歷程對比

3 數值計算結果與分析

在以下的計算中，將圓柱的旋轉振動幅值和橫流向平移振動幅值分別取為θ*=1.0和A*=0.5(此時圓柱表面的最大切向速度與圓柱的最大平動速度具有相同的數值)，重點考察強迫運動頻率f*(fe/fs)以及旋轉振動與橫流向振動之間的相位差φ的影響作用。計算中，雷諾數Re=100保持不變，相應雷諾數下固定圓柱渦脫落頻率為fs=0.165 1。

3.1 圓柱受力特性及鎖定分析

圖4以相位差φ=-60°，強迫振動頻率f*=0.55、0.80、1.20及2.00為例，給出了升力系數和圓柱橫流向振動位移的時間歷程線，對應的升力系數頻譜分析結果，如圖5所示。對比圖4結果可以看出，隨著圓柱振動頻率f*的增大，升力系數的幅值顯著增加，說明圓柱振動頻率是影響圓柱受力的一個重要因素。從圖4(a)中可以看出，當圓柱的振動頻率f*=0.55時，升力系數的時間歷程線不規則。相應的頻譜分析結果圖5(a)表明，此時的升力頻率成份同時包含了圓柱的振動頻率fe及固定圓柱的渦脫落頻率fs以及由二者之間相互作用所產生的新頻率成份(fs-fe)/2及(3.0fe+fs)/2等。當f*=0.80時，從圖4(b)中可見，此時圓柱所受升力變化基本滿足簡諧振動的形式。圖5(b)的升力系數頻譜分析表明，相應的渦脫落頻率(也即升力主控頻率)與圓柱的振動頻率相等，這意味著圓柱后方的渦旋脫落主要由圓柱的振動所控制，此時尾流場發生鎖定現象。此外，雖然圖5(b)的頻譜中包含了3.0fe頻率成份，但其對升力貢獻很小，因此，此時的升力振蕩主要表現為周期性的簡諧振動的形式。對于f*=1.20的情況，圖4(c)的結果表明，此時的升力系數具有明顯的多頻振動特性，并且每5個圓柱橫流向振動周期對應1個完整的升力系數振蕩周期，即圓柱后方發生一次完整的尾渦脫落周期為圓柱振動周期的5倍。在本文中，將尾渦脫落周期與圓柱振動周期相等的工況定義為基本鎖定(Primary Lock-in)，這與傳統的圓柱渦激振動鎖定的定義是一致的[24]；而將完整的尾渦脫落周期為N(N>1)倍圓柱振動周期的情況定義為準周期鎖定(Quasi-periodic Lock-in)。對應于前述的圖4(c)的情況，準周期鎖定發生在5T的周期上。與圖4(c)的情況相類似，當f*=2.0時，尾流場也存在準周期鎖定現象。如圖4(d)所示，此時每3個橫向(旋轉)振動周期，升力系數出現周期性重復特征，對應的準周期鎖定周期為3T。圖5(c)和(d)分別給出了與圖4(c)和(d)相對應的升力系數的頻譜分析結果。從圖4中可以看出，對于這兩個發生準周期鎖定的情況，相應的渦脫落頻率分別為fe/5和fe/3，即fe/N(N=3, 5，…)，而不是圓柱的振動頻率fe，但升力振蕩的主控頻率均對應圓柱的振蕩頻率。同時，在升力變化以及渦旋脫落的過程中，對應于N倍準周期鎖定，還存在(2M-1)fe/N的次諧波頻率成份(M=1, 2, …)。

根據非線性振動力學的基本理論，僅依靠頻率關系來判斷是否發生鎖定是不全面的。為進一步證實上述的基本鎖定和準周期鎖定現象，圖6給出了與圖4相對應的升力系數與圓柱橫流向振動位移的相圖。從圖6(a)可以看出，由于圓柱升力系數振動的不規則性(未鎖定，參見圖4(a))，升力系數和圓柱振動位移的相圖呈現出雜亂無章的特征，構成典型的非周期性混沌運動狀態。圖6(b)給出了發生基本鎖定時的相圖，從圖中可以看出，由于升力系數的振動周期與圓柱振動周期一致，相圖為嚴格封閉的單環結構，這充分證實了基本鎖定現象的發生。圖5(c)和(d)分別給出了在fe/5和fe/3頻率下，發生準周期鎖定時的相圖。從圖6中可以看出，此時的相圖也表現為封閉的圓環。但是，由于升力振蕩周期分別為圓柱振動周期的5倍和3倍，因此對應的相圖中也同時包含了5個和3個封閉環。

(a) φ=-60°, f*=0.55

(b) φ=-60°, f*=0.80

(c) φ=-60°, f*=1.20

(d) φ=-60°, f*=2.00

(a)

(b)

(c)

(d)

(a) φ=-60°, f*=0.55

(b) φ=-60°, f*=0.80

(c) φ=-60°, f*=1.20

(d) φ=-60°, f*=2.00

Fig.6 Phase diagram of lift coefficient and transverse displacement

另外，對于圖6(b)所示的在f*=0.80情況下的基本鎖定情況，由于相圖的主軸左向傾斜，圓柱橫流向振動位移和升力系數之間的相位差Ψ介于π/2和π之間，這表明流體力起到阻尼的作用。對于一般的渦激振動問題而言，當Ψ∈(π/2, π)時，流體作用力將起到抑制圓柱渦激振動響應的作用。但是，隨著圓柱振動頻率的進一步增大，如圖6(c)所對應的f*=1.20及圖6(d)所對應的f*=2.00的情況，圓柱橫流向振動與升力系數之間的相位差Ψ介于0到π/2之間(相圖主軸右向傾斜)，此時流體向圓柱輸入能量，流體力起到激振力的作用[3]。應該注意的是，對于上述的f*=1.20和2.00的情況，二者均對應準周期鎖定。

本文通過綜合考慮升力系數的頻譜特性以及升力振蕩與圓柱振動之間的相圖來判斷基本鎖定和準周期鎖定的發生與否。表2全面給出了在本文所研究的31個受迫振動頻率(f*=0.5～2.0，頻率間隔0.05)以及7個圓柱轉動與橫向振動之間相位差(φ=±180°, ±120°, ±60°和0°)下，共計217組數值計算所得到的鎖定判別結果。此外，為了便于對比分析，表中也一并列出了圓柱在相同雷諾數(Re=100)、相同振動幅值(A*=0.5)以及相同振動頻率下，僅發生橫流向受迫振動的鎖定區間判別結果。應該說明的是，對于圓柱僅發生橫流向振動的情況下，只存在基本鎖定的現象，相應的鎖定頻率區間范圍為[0.7, 1.2]，這與以往的研究成果是一致的[8]。

表2 不同相位差φ及振動頻率f*下的鎖定區間

Tab.2 Lock-in regime at various phase differencesφand forced oscillatory frequencyf*

f?θ?=0φ(θ?=1)180°120°60°0°-60°-120°-180°0．50×××××13T7T×0．55××××4T×××0．60×××■××■×0．65×■×■■■■■0．70■■■■■■■■0．75■■■■■■■■0．80■■■■■■■■0．85■■■■■×■■0．90■■■■■×■■0．95■■■■■×■■1．00■■■■■××■1．05■■■■■×11T■1．10■×■■■×××1．15■6T■■■××6T1．20■×■■■5T11T×1．25××■■■5T2T×1．30××■■■×2T×1．35××■■■×2T×1．40×5T■■■2T×5T1．45××■■■2T11T×1．50××■■■19T××1．55××■■■×7T×1．60××■■■×××1．65××■■■×××1．70×2T■■■7T3T2T1．75××■■■×3T×1．80×17T■■■×3T17T1．85××■■■×3T×1．90××■■■×3T×1．95××■■■3T3T×2．00××■■■3T××注：×，非鎖定；■，基本鎖定(尾渦脫落周期與圓柱的振動周期T一致)；NT，準周期鎖定(完整的尾渦脫落周期為圓柱振動周期T的N倍)

從表2的結果中可以看出，當引入圓柱的旋轉運動后，不但會誘發準周期鎖定現象，同時，圓柱轉動與橫流向振動之間的相位差φ對鎖定形式以及鎖定的頻率區間都有十分重要的影響作用。表2的分析結果表明：① 當0°≤φ<180°時，鎖定現象主要以基本鎖定形式為主；當-180°<φ<0°時，發生準周期鎖定的概率顯著增加；② 當φ=±180°時，無論是基本鎖定，還是準周期鎖定，正負兩個相位差下所得到的數值計算結果保持嚴格一致；③ 當φ=±60°和±120°時，正負兩種相位差下的鎖定形式和鎖定區間存在較大的差別。

對于φ=60°和120°的情況(正的相位差)，基本鎖定區間從低頻到高頻連續分布，未出現準周期鎖定現象，其中φ=60°的情況，取得最大的基本鎖定區間f*∈[0.6, 2.0]；當φ=-60°和-120°時(負的相位差)，基本鎖定區間主要集中在低頻區，而準周期鎖定主要發生在高頻區域。特別的，對于φ=-60°和-120°這兩種情況，即便圓柱的振動頻率與相同雷諾數下固定圓柱的渦脫落頻率相等，即在f*=fe/fs=1.0的條件下，也未發生鎖定。當φ=-60°時，基本鎖定區間具有最小的頻率范圍，f*∈[0.65, 0.8]。對于同時具有旋轉與平移耦合振動的圓柱結構而言，其鎖定的發生條件以及鎖定形式都與這兩種運動之間的相位差密切相關，特別是在對應的正負相位差下所表現出的顯著性差異，可能是非線性遲滯(hysteresis)效應的體現，其中的物理機理還有待于進一步的深入研究。

(a) 升力系數均方根

(b) 平均拖曳力系數

對于相位差為φ=-120°和-60°的特殊情況，由于它們在低頻區內都存在鎖定區間的明顯變化，因此在相應的頻率范圍內，也都出現了升力幅值的跳躍現象。可見，鎖定形式的轉換是導致升力發生變化的原因之一。在各相位差下，平均拖曳力系數隨振動頻率的變化情況相對比較復雜。從圖7(b)中可以看出，當φ=120°時，平均拖曳力系數在各頻率下幾乎都是最小的，且在整個頻率范圍內變化不大。而與之對應的φ=-120°，其平均拖曳力則一直維持在較高的水平。從圖7(b)中可見，相對于僅有橫流向振動的情形，當φ=±180°，-120°，-60°及0°時，無論是鎖定區間還是非鎖定區間，引入圓柱旋轉運動后，平均拖曳力系數均會有所增大。參照表2的結果可以發現，圖7(b)中有關平均拖曳力系數發生跳躍變化的頻率幾乎都對應著鎖定模式的改變。

3.2 基本鎖定與準周期鎖定下的渦脫落模式

圖8以φ=-60°，f*=0.8的工況為例(對應圖4(b))，給出了在基本鎖定情況下(即渦脫落頻率與圓柱的振動頻率相等)，一個完整渦脫落周期內渦量場的變化情況。圖中，白色表示正渦量，黑色表示負渦量；N0時刻對應圓柱的橫流向振動位移正處于由負轉正的平衡位置，各圖之間的時間間隔為1/4橫向振動周期。

從圖8中可見，在N0時刻，圓柱背流側的下方和上方分別形成正渦A和負渦B，其中正渦A的尺度和強度明顯大于負渦B。由于圓柱旋轉自由度的存在，在圓柱的右下表面可見明顯的剪切層。隨著時間的推移，正渦A在N1時刻發生脫落，并在N2時刻形成新的貼體正渦C。而負渦B在圓柱運動與來流的共同作用下，發生顯著扭曲，同時在圓柱右上表面形成明顯的負渦量剪切層。N2時刻的渦量場與N0時刻構成反對稱結構。在N3時刻，負渦B發生脫落，并形成了新的貼體負渦D。同樣，N3時刻的渦量場與N1時刻也呈反對稱結構。之后，經過T/4時間后，在N4時刻，渦流場又重新恢復到初始時刻N0的狀態，完成一個完整的尾渦脫落過程。從以上的分析中可以看出，對于本例的基本鎖定工況，在一個圓柱諧振周期內，會釋放方向相反的兩個渦(即正渦A和負渦B)，這是典型的2S尾渦脫落模式。值得注意的是，由于本文的圓柱運動還同時涉及到繞軸心的轉動，轉動運動分量的出現會使得圓柱表面出現額外的剪切作用，從而形成附加的剪切層，并影響分離點的位置。在圖8中，由旋轉自由度引起的附加剪切層可明顯見于N2和N4時刻，其剪切作用的強弱與相位差有直接關系，并可影響到圓柱的受力以及鎖定等宏觀特性，這與簡單的單自由度橫流向受迫振動情形是不同的[24]。

圖8 φ=-60°，f*=0.8時一個渦脫落周期的渦量圖

圖9進一步以φ=-60°和f*=1.2條件下的準周期鎖定情況為例(渦脫落周期為圓柱振動周期的5倍，對應圖4(c))，給出了渦流場的演化特性，其中各圖之間的時間間隔為1/2圓柱的振動周期。從圖9中可以看出，在初始時刻N0，圓柱的背流側存在一個由負渦A和正渦B構成的渦對。經過半個圓柱振蕩周期后，在N1時刻，圓柱右上方的負渦A分裂為一個貼體負渦A1和一個獨立負渦A2。由于A2是單獨脫落的負渦，在此將其定義為S-形式。此后的N2和N3階段，圓柱尾流區內未見明顯的渦旋脫落，伴隨著A2不斷向下游輸運，A1和B渦都得到了進一步的發展。在N4時刻，即圓柱完成2個完整的振蕩周期時，上述的貼體負渦A1分裂成一個新的貼體渦A11和一個獨立的分離渦A12，其中后者仍以先前的S-單渦形式發生脫落。與此同時，正渦B也一次性分裂為B1、B2和B3，其中B1以貼體渦的形式存在，而B2和B3以成對的正渦形式同時脫落，本文將其定義為P+形式。應該指出的是，在N4時刻同步出現的S和P渦脫落模式是以往在僅有橫向振動的情況下所未見的，這體現了旋轉自由度的重要作用。在N5時刻，貼體渦B1和A11都得到了進一步的發展，同時B2與B3在對流作用下，他們之間的分離也越發明顯。至N6時刻，貼體渦B1也發生分裂，形成B11和B12，其中B12是以S+的形式脫落，而B11仍附著在圓柱表面。N7和N8是尾流場的進一步發展階段，貼體渦B11和A11得到了充足的拉伸，并具備了進一步分裂的條件。在N9時刻，再次出現了S和P的同步分裂模式。B11分裂為貼體渦B111和單體脫落渦B112，后者構成S+的脫落形式。而A11則一次性分裂為貼體渦A111和一個負的脫落渦對A112與A113，即P+的形式。此后，在經過5個完整圓柱振蕩周期后，流場于N10時刻重新恢復到初始的N0時刻狀態，此時的貼體渦旋B111和A111分別對應N0時刻的B渦和A渦。

圖9 φ=-60°，f*=1.2時完整渦脫落周期的渦量圖

縱觀上述整個渦旋脫落過程，可以發現，N0與N5、N1與N6、N2與N7、N3與N8以及N4與N9，互為反對稱渦量場(渦旋結構)。同時，整個尾跡區的演化過程包括了渦旋的多次二級分裂和三級分裂，這些分裂在本質上是對應著非線性問題中的分叉行為。而分叉的發生正是流動不穩定性的體現，也是層流向湍流轉捩的重要形式之一[25]。

針對上述的5倍振動周期的準周期鎖定情況，相應的尾跡區渦旋分裂過程總結于圖10中。對于初始時刻的正渦A，先后在N1和N4時刻發生兩次二級分裂，之后于N9時刻進一步發生三級分裂。對于正渦B，則是首先在N4時刻發生一次三級分裂，之后于N6和N9時刻分別發生兩次二級分裂。這樣，在一個尾跡演化周期內，由初始的一對渦旋，共分裂生成14個次級渦旋。而就渦旋的脫落過程而言，歷經了{S-}，{S-，P+}，{S+}，{S+，P-}，共形成2正2負4個單渦以及1正1負2個渦對的脫落。本文在此僅給出了以5T為代表的準周期鎖定情況下的尾流場結構演化特征，對于其它具有更低頻率的準周期鎖定工況，其尾渦模式也更為復雜。但就總體特征而言，都具有復雜的多級次渦旋分叉過程。

4 結 論

本文通過高階Taylor-Galerkin迎風有限元方法求解二維不可壓縮Navier-Stokes方程，并結合任意拉格朗日-歐拉動網格方法，針對圓柱在均勻流中，同時發生旋轉運動及橫流平動的兩自由度諧振問題開展了數值模擬研究，重點討論了圓柱諧振頻率f*以及轉動與平動之間的相位差φ對圓柱受力、鎖定特性以及尾渦脫落模式的影響作用。本文研究工作的主要結論如下：

(1) 對于基本鎖定工況，尾渦脫落頻率與圓柱的振動頻率一致；而對于準周期鎖定的情況，完整的渦脫落過程所對應的頻率為fe/N(N=3, 5，…)。此外，隨著流場的演化，升力系數中還存在頻率為(2M-1)fe/N的次諧波成份(M=1, 2, …)。

(2) 圓柱旋轉運動與橫流向振動之間的相位差φ對鎖定區間有重要影響作用。當0°≤φ<180°時，鎖定主要以基本鎖定為主，而-180°<φ<0°時，基本鎖定所對應的頻率范圍比較狹窄，在高頻區廣泛分布準周期鎖定現象；當φ=±180°時，無論是基本鎖定還是準周期鎖定，兩個相位差下的數值計算結果完全一致；而當φ=±120°及±60°時，正負相位差條件下的基本鎖定區間及準周期鎖定的分布均存在較大的差別，體現出流體與結構相互作用的非線性特性。

(3) 圓柱旋轉運動與橫流向振動之間的相位差φ對圓柱的受力也存在較大的影響作用。總體上，在各相位差條件下，升力系數的均方根均隨著圓柱振動頻率的增大而增大；相比于圓柱只發生橫向強迫振動的情況，當引入圓柱旋轉振動后，鎖定區間內的升力系數均方根有不同程度的增大；平均拖曳力系數的特性隨圓柱振動頻率的變化則相對復雜。當φ=120°時，平均拖曳力系數最小，并且隨圓柱振動頻率的增大幾乎不發生變化；相比于圓柱只發生橫向強迫振動的情況，除φ=120°及60°，其他相位差條件下所對應的平均拖曳力系數均增大。

(4) 對于基本鎖定情況，其尾跡區內為典型的2S渦脫落模式；本文以φ=-60°，f*=1.2為例，給出了準周期鎖定下的渦脫落模式。數值計算結果表明，在一個完整的尾渦脫落周期內，會存在渦旋結構的多次二級分裂及三級分裂，這些分裂的過程在本質上對應著非線性問題中的Hopf分叉行為，同時這也是流動不穩定性的具體表現。

[1] WILLIAMSON C H K, ROSHKO A. Vortex formation in the wake of an oscillating cylinder[J]. Journal of Fluid and Structures, 1988, 2(4): 355-381.

[2] LU X Y, DALTON C. Calculation of the timing of vortex formation from an oscillating cylinder[J]. Journal of Fluid and Structures, 1996, 10(5): 527-541.

[3] BLACKBURN H M, HENDERSON R D. A study of two-dimensional flow past an oscillating cylinder[J]. Journal of Fluid Mechanics, 1999, 385: 255-286.

[4] ANAGNOSTOPOULOS P. Numerical study of the flow past a cylinder excited transversely to the incident stream. Part 1: Lock-in zone, hydrodynamic forces and wake geometry[J]. Journal of Fluid and Structures, 2000, 14(6): 819-851.

[5] CARBERRY J, SHERIDAN J, ROCKWELL D. Forces and wake modes of an oscillating cylinder[J]. Journal of Fluid and Structures, 2001, 15: 523-532.

[6] GUILMINEAU E, QUEUTEY P. A numerical simulation of vortex shedding from an oscillating circular cylinder[J]. Journal of Fluid and Structures, 2002, 16(6): 773-794.

[7] SARPKAYA T. A critical review of the intrinsic nature of vortex-induced vibration[J]. Journal of Fluid and Structures, 2004, 19(4): 389-447.

[8] LEONTINI J S, STEWART B E, THOMPSON M C, et al. Wake state and energy transitions of an oscillating cylinder at low Reynolds number[J]. Physics of Fluids, 2006, 18: 067101.

[9] ONGOREN A, ROCKWELL D. Flow structure from an oscillating cylinder. Part 2. Mode competition in the near wake[J]. Journal of Fluid Mechanics, 1988, 191: 225-245.

[10] CETINER O, ROCKWELL D. Streamwise oscillations of a cylinder in a steady current. Part 1. Locked-on states of vortex formation and loading[J]. Journal of Fluid Mechanics, 2001, 427: 1-28.

[11] XU S J, ZHOU Y, WANG M H. A symmetric binary-vortex street behind a longitudinally oscillating cylinder[J]. Journal of Fluid Mechanics, 2006, 556: 27-43.

[12] LEONTINI J S, JACONO D L, THOMPSON M C. Wake state and frequency selection of a streamwise oscillating cylinder[J]. Journal of Fluid Mechanics, 2013, 730: 162-192.

[13] DU L., SUN X. Suppression of vortex-induced vibration using the rotary oscillation of a cylinder[J]. Physics of Fluids, 2015, 27: 023603.

[14] BAI W. Numerical simulation of flow past a rotating and rotary oscillating circular cylinder on unstructured meshes[J]. Coupled System Mechanics, 2013(2): 191-241.

[15] MITTAL S, KUMAR B. Flow past a rotating cylinder[J]. Journal of Fluid Mechanics, 2003, 476: 303-334.

[16] LU L, QIN J M, TENG B, et al. Numerical investigations of life suppression by feedback rotary oscillation of circular cylinder at low Reynolds number[J]. Physics of Fluids, 2001, 23: 033601.

[17] BAEK S J, SUNG H J. Numerical simulation of the flow behind a rotary oscillating circular cylinder[J]. Physics of Fluids, 1998, 10(4): 869-876.

[18] CHOI S, CHOI H, KANG S. Characteristics of flow over a rotationally oscillating cylinder at low Reynolds number[J]. Physics of Fluids, 2002, 14(8): 2767-2777.

[19] BLACKBRUN H M, ELSTON J R, SHERIDAN J. Bluff-body propulsion produced by combined rotary and translation oscillation[J]. Physics of Fluids, 1999, 11(1): 4-6.

[20] NAZARINIA M, JACONO D L, THOMPSON M C, et al. Flow behind a cylinder forced by a combination of oscillatory translational and rotational motions[J]. Physics of Fluids, 2009, 21: 051701.

[21] JIANG C B, KAWAHARA M. A three-step finite element method for unsteady incompressible flows[J]. Computational Mechanics, 1993, 11: 355-370.

[22] 郭曉玲, 唐國強, 劉名名, 等. 低雷諾數下串聯雙圓柱渦激振動機理的數值研究[J]. 振動與沖擊, 2014, 33(4): 60-69.

GUO Xiaoling, TANG Guoqiang, LIU Mingming, et al. Numerical investigation on vortex induced vibration of twin tandem circular cylinders under low Reynolds number[J]. Journal of Vibration and Shock, 2014, 33(4): 60-69.

[23] 劉為民, 谷家揚, 陶延武, 等. 低質量比圓形四立柱渦激運動特性研究[J]. 振動與沖擊, 2015, 34(19): 175-180.

LIU Weimin, GU Jiayang, TAO Yanwu, et al. Vortex induced motion characteristics of four circular columns with a low mass ratio[J]. Journal of Vibration and Shock, 2015, 34 (19): 175-180.

[24] WILLIAMSON C H K, GOVARDHAN R. Vortex-induced vibration[J]. Annual Review of Fluid Mechanics, 2004, 36: 413-455.

[25] 朗道 Л Д, 栗弗席茲 E M. 流體動力學[M]. 5版. 李植 譯, 北京: 高等教育出版社, 2013.

Hydrodynamic characteristics of laminar flow over a circular cylinder having forced rotational and transverse harmonic oscillations

LIU Mingming, TANG Guoqiang, Lü Lin, TENG Bin

(State Key Laboratory of Coastal and Offshore Engineering, Dalian University of Technology, Dalian 1160241, China)

Based on the two-dimensional finite element solution to incompressible viscous Navier-Stokes equations using the Arbitrary Lagrangian-Eulerian (ALE) method, the hydrodynamic characteristics of laminar flow over a circular cylinder having forced rotational and transverse harmonic oscillations were examined. The effects of oscillatory frequency (f*[0.5, 2.0] with an interval of 0.05) and phase difference between the rotation and translation motions (φ=±180°, ±120°, ±60° and 0°) on fluid force, lock-in mode and wake mode were analyzed. The numerical results showed that two distinct lock-in modes are identified under the present two forced oscillation excitations, they are the ordinary primary lock-in and the nonlinear quasi-periodic lock-in; the phase difference between the rotation and translation motions has an important influence on the types of lock-in and the corresponding lock-in frequency intervals, symmetric positive and negative phase differences can cause highly consistent or completely different lock-ins; fluid forces exerted on the circular cylinder generally increase with increase in oscillating frequency, the obvious hydrodynamic coefficient’s jumps are observed under the condition of lock-in; for the specific quasi-periodic lock-in cases, there are multi-stage nonlinear bifurcations of vortex structure in the field of wake flow, the wake evolution process follows quasi-periodic features.

circular cylinder; rotational oscillation; transverse oscillation; quasi-period lock-in; wake mode

國家重點基礎研究發展計劃(973計劃)(2014CB046803);國家自然科學基金(51409035；51279029)

2016-01-12 修改稿收到日期：2016-04-12

劉名名 男，博士生，1986年9月生

呂林 男，副研究員，1976年生

P751

A

10.13465/j.cnki.jvs.2017.11.006