利用延遲反饋控制柔性圓柱渦激振動(dòng)數(shù)值研究

鄒 琳，王 程，徐勁力，陶 凡，左紅成

(武漢理工大學(xué) 機(jī)電工程學(xué)院，武漢 430070)

鈍體周圍非對(duì)稱的漩渦脫落在結(jié)構(gòu)上產(chǎn)生周期性的脈動(dòng)壓力，從而引發(fā)渦激振動(dòng)(vortex induced vibration，VIV)現(xiàn)象[1-3]。當(dāng)漩渦脫落頻率接近結(jié)構(gòu)固有頻率時(shí)會(huì)發(fā)生大幅度的振動(dòng)，導(dǎo)致結(jié)構(gòu)疲勞損壞，也可以利用渦激振動(dòng)產(chǎn)生新能源[4-5]。對(duì)于類似于海洋立管的柔性圓柱而言，由于展向具有柔性張力，渦激振動(dòng)通常表現(xiàn)為更為復(fù)雜的駐波響應(yīng)特性、行波響應(yīng)特性和多模態(tài)振動(dòng)特性等。

國(guó)內(nèi)外學(xué)者對(duì)柔性圓柱渦激振動(dòng)展開(kāi)了大量數(shù)值及試驗(yàn)研究。早期數(shù)值研究主要是利用直接數(shù)值模擬方法求解N-S方程，以獲得柔性圓柱的渦激振動(dòng)響應(yīng)，如Newman等[6]通過(guò)數(shù)值模擬方法研究了Re=100和Re=200下柔性圓柱的橫向流振動(dòng)，發(fā)現(xiàn)振動(dòng)模式為行波與駐波混合響應(yīng)模式；Lucor等[7]通過(guò)直接數(shù)值模擬研究了長(zhǎng)徑比超過(guò)500的細(xì)長(zhǎng)圓柱在線型和指數(shù)剪切流流作用下的振動(dòng)響應(yīng)，研究認(rèn)為指數(shù)剪切流速剖面會(huì)產(chǎn)生更大的振動(dòng)頻率范圍，但和實(shí)際情況相比，圓柱單自由度振動(dòng)具有一定的局限性。

Lie等[8]開(kāi)展了剪切流作用下物理模型試驗(yàn)，研究立管長(zhǎng)細(xì)比高達(dá)3 000，分析發(fā)現(xiàn)立管的振動(dòng)響應(yīng)是不規(guī)則的，具有一定的帶寬，并且強(qiáng)度也隨著帶寬的不同有所變化。國(guó)內(nèi)學(xué)者高云等[9]通過(guò)試驗(yàn)研究，深入研究了柔性圓柱體在剪切來(lái)流下渦激振動(dòng)響應(yīng)，發(fā)現(xiàn)當(dāng)折合速度較小時(shí)，響應(yīng)頻率會(huì)全程參與渦激振動(dòng)，隨著折合流速的上升，響應(yīng)頻率則會(huì)間歇性參與渦激振動(dòng)。宋磊建等[10]采用模型試驗(yàn)的方法研究了均勻流下柔性立管的渦激振動(dòng)響應(yīng)特性及渦激力載荷特性，研究發(fā)現(xiàn)均勻流下柔性立管的渦激振動(dòng)響應(yīng)為位移和主導(dǎo)頻率不隨時(shí)間變化的穩(wěn)態(tài)響應(yīng)。

對(duì)于大長(zhǎng)徑比的柔性圓柱體，數(shù)值模擬全尺度的振動(dòng)響應(yīng)存在一定的局限性[11-12]，而試驗(yàn)方法研究成本相對(duì)較高[13-15]，因此有必要建立一種能夠快速預(yù)報(bào)響應(yīng)細(xì)長(zhǎng)圓柱渦激振動(dòng)的經(jīng)驗(yàn)?zāi)Ｐ头椒M細(xì)長(zhǎng)柔性體振動(dòng)響應(yīng)主要特征。尾流振子模型是應(yīng)用較為廣泛的經(jīng)驗(yàn)?zāi)Ｐ?。最早Noack等[16]利用范德波爾振子方程推導(dǎo)出能夠預(yù)測(cè)細(xì)長(zhǎng)鈍體渦脫頻率和角度的經(jīng)驗(yàn)?zāi)Ｐ停脦в袛U(kuò)散項(xiàng)的范德波爾方程能夠準(zhǔn)確預(yù)測(cè)出尾流脫落的定性特征。之后Balasubramanian等[17]在此基礎(chǔ)上利用帶有擴(kuò)散項(xiàng)的范德波爾振子方程，用來(lái)預(yù)測(cè)均勻流和剪切流下直圓柱和錐形圓柱的漩渦脫落。Facchinetti等[18]利用尾流振子方程和結(jié)構(gòu)柔性振動(dòng)方程，研究了三維柔性圓柱的渦激振動(dòng)響應(yīng)，并與試驗(yàn)結(jié)果對(duì)比指出耦合模型能夠準(zhǔn)確預(yù)測(cè)渦激振動(dòng)響應(yīng)的振動(dòng)形態(tài)。高云等[19-20]通過(guò)尾流振子模型，對(duì)線型剪切來(lái)流下兩端鉸接的細(xì)長(zhǎng)柔性圓柱的渦激振動(dòng)響應(yīng)做了詳細(xì)的探究。

綜上所述，當(dāng)前研究主要是對(duì)細(xì)長(zhǎng)柔性圓柱的渦激振動(dòng)響應(yīng)進(jìn)行研究，如何控制三維柔性圓柱渦激振動(dòng)的研究相對(duì)較少。因此，本文借助尾流振子模型，提出利用延遲反饋來(lái)實(shí)現(xiàn)柔性圓柱渦激振動(dòng)主動(dòng)控制，研究延遲時(shí)間τ和延遲增益kd對(duì)柔性圓柱振動(dòng)幅值和振動(dòng)頻率的影響，尋求控制柔性圓柱渦激振動(dòng)的延遲反饋規(guī)律，為海洋立管等柔性鈍體的振動(dòng)控制提供有益探索。

1 數(shù)學(xué)模型

柔性圓柱流固耦合作用下產(chǎn)生振動(dòng)幅值主要受到柔性體張力的影響，本文考慮無(wú)限長(zhǎng)張緊圓柱在受到均勻來(lái)流作用下的渦激振動(dòng)響應(yīng)，忽略結(jié)構(gòu)的彎曲剛度，模型示意如圖1所示。

圖1 受到均勻來(lái)流的無(wú)限長(zhǎng)張緊圓柱Fig.1 An infinite tension cylinder subjected to a stationary uniform flow

圖中,圓柱直徑為D，長(zhǎng)度為L(zhǎng)，兩端受到的軸向張力為N，在均勻來(lái)流Uf作用下產(chǎn)生振動(dòng)，本文只考慮y方向的振動(dòng)位移。參考文獻(xiàn)[21]，本文選取長(zhǎng)徑比L/D=25進(jìn)行討論。在延遲反饋控制的作用下，所選柔性圓柱部分的橫向位移Y可以表示為

(1)

式(1)的右端項(xiàng)P(Z,T)為作用在結(jié)構(gòu)上的流體力，可以表示為P(Z,T)=ρU2DCL0q(Z,T)/4，其中CL0為固定圓柱的參考升力系數(shù)[22]，q(Z,T)是用來(lái)描述尾流運(yùn)動(dòng)的尾跡變量。φ(Z,T)=10-3sin(50πt)為初始位置函數(shù)，圖2計(jì)算了不同初始位置函數(shù)φ(Z,T)下柔性圓柱縱向振動(dòng)位移的均方值，說(shuō)明數(shù)值計(jì)算結(jié)果與φ(Z,T)的選取無(wú)關(guān)。

圖2 不同初始位置函數(shù)φ(Z,T)下柔性圓柱縱向振動(dòng)位移的均方值Fig.2 RMS values of dimensionless displacements at different initial location function φ(Z,T)

尾流的運(yùn)動(dòng)可以用非線性振子來(lái)表示，常見(jiàn)的表示形式為范德波爾振子，寫作

(2)

式中，A和ε為經(jīng)驗(yàn)參數(shù)。引入無(wú)量綱參數(shù)y=Y/D，z=Z/D，t=TΩf，式(1)和式(2)可以寫作如下無(wú)量綱形式

(3)

(4)

2 數(shù)值方法及模型驗(yàn)證

2.1 數(shù)值方法

(5)

式中，k_d為延遲項(xiàng)所對(duì)應(yīng)的索引，將式(5)代入式(3)、式(4)可以得到如下迭代格式

(6)

(7)

(8)

2.2 模型驗(yàn)證

首先對(duì)數(shù)值離散方法的時(shí)間步長(zhǎng)和空間步長(zhǎng)進(jìn)行無(wú)關(guān)性驗(yàn)證。如表1所示，時(shí)間步長(zhǎng)Δt取較大值0.02和0.03時(shí)，對(duì)于較小的空間步長(zhǎng)Δz，數(shù)值計(jì)算結(jié)果發(fā)散。當(dāng)Δt=0.005，Δz=0.05時(shí)，中間節(jié)點(diǎn)處最大振幅為1.849，而Δt=0.01，Δz=0.1時(shí)，最大振幅為1.847，降低約0.1%；當(dāng)Δt=0.01，Δz=0.4時(shí)，最大振動(dòng)幅值為1.785，降低約3.4%。因此在后續(xù)的數(shù)值計(jì)算中，取Δt=0.01，Δz= 0.1，整個(gè)數(shù)值仿真時(shí)間Tend=1 400。

表1 不同時(shí)間和空間步長(zhǎng)下，中間節(jié)點(diǎn)處結(jié)構(gòu)振幅Tab.1 The vibration amplitude of cylinder at node under different space and time steps

為了驗(yàn)證文中所選模型的正確性，取L/D=25，c=4,μ=1.785，St=0.16，與Newman等所用參數(shù)保持一致。不考慮延遲反饋?zhàn)饔孟拢Y(jié)構(gòu)的振動(dòng)位移云圖如圖3(a)所示，橫軸的無(wú)量綱時(shí)間為t/2πSt?？梢园l(fā)現(xiàn)，在時(shí)間段200～240內(nèi)，三維柔性圓柱的振動(dòng)形式為駐波；當(dāng)無(wú)量綱時(shí)間達(dá)到260時(shí)，振動(dòng)形式由駐波逐漸向行波轉(zhuǎn)變，此后一直保持行波振動(dòng)狀態(tài)，這一過(guò)程與Newman等DNS計(jì)算結(jié)果一致(如圖3(b)所示)，說(shuō)明本文的計(jì)算模型是可行的。

圖3 渦激振動(dòng)位移響應(yīng)對(duì)比(L/D=25，c=4,μ=1.785)Fig.3 Comparison of VIV displacement response with time and space (L/D=25，c=4,μ=1.785)

3 結(jié)果和討論

3.1 延遲時(shí)間τ對(duì)振幅和頻率的影響

從圖4(a)～圖4(d)均可以看出，在特定的延遲增益kd下，振幅隨著延遲時(shí)間τ呈周期性變化，且每個(gè)周期內(nèi)振幅最大值(定義為振幅峰值)保持不變，延遲時(shí)間τ并不能改變振幅峰值，其大小只與延遲增益kd有關(guān)，隨著kd的增加，振幅峰值要更加偏離未施加延遲反饋時(shí)的振動(dòng)幅值0.647，當(dāng)kd為0.05，0.10和0.15時(shí)，振動(dòng)峰值分別為0.85，1.18和1.82。從圖4(a)可以看出，當(dāng)延遲時(shí)間τ=0時(shí)，振動(dòng)幅值分別為0.518，0.429和0.363，振動(dòng)幅值隨延遲增益kd的增大而減小，均小于未受控時(shí)的振動(dòng)幅值0.647，最大能夠降低43.9%。觀察圖4(b)～圖4(d)可以發(fā)現(xiàn)，當(dāng)kd符號(hào)相反時(shí)，振動(dòng)幅值變化曲線之間會(huì)相差半個(gè)周期的相位。隨著τ的增加，振幅最小值都大于τ=0時(shí)的振幅值，并且逐漸趨向于振幅峰值，而且隨著kd的增大，趨向的速度越來(lái)越快，如箭頭所示。綜上發(fā)現(xiàn)，不同的延遲時(shí)間τ耦合不同的延遲增益kd，可以實(shí)現(xiàn)不同柔性圓柱振動(dòng)幅值的控制，這將有益于實(shí)現(xiàn)工程實(shí)踐中柔性鈍體的振動(dòng)控制。

如圖5所示，振動(dòng)頻率隨延遲時(shí)間也是呈周期性的變化，振動(dòng)頻率隨著延遲時(shí)間τ的增加逐漸趨向于未施加延遲控制柔性圓柱振動(dòng)頻率f=0.156 4，且通過(guò)圖5(b)～圖5(d)可以發(fā)現(xiàn)，延遲增益kd越大，趨向的速度逐漸越快，如箭頭所示。

3.2 延遲增益kd對(duì)振幅和頻率的影響

圖6(a)和圖6(b)的延遲時(shí)間τ分別是圖4(a)振幅取最大和最小時(shí)所對(duì)應(yīng)的延遲時(shí)間。從圖6(a)中可以看出，當(dāng)kd由0逐漸增加時(shí)，振幅也逐漸增加，振幅只與kd有關(guān)。當(dāng)kd從0逐漸減小，在τ取較小值(3.2，9.6和15.9)時(shí)，延遲增益對(duì)振幅的影響占主導(dǎo)，振幅隨kd減小而減小；當(dāng)τ取值較大的情況下，延遲時(shí)間對(duì)振幅的影響占主導(dǎo)地位，振幅隨kd的減小呈現(xiàn)先減小后增大的趨勢(shì)。在相同負(fù)值kd下，延遲時(shí)間越大振幅越大，這與圖4所示振幅最小值逐漸趨向于峰值的變化趨勢(shì)是一致的，而且kd越大趨向速度越快。由于延遲增益kd由正變負(fù)時(shí)變化曲線會(huì)相差半個(gè)周期的相位，因此圖6(b)中振幅隨kd的變化規(guī)律與圖6(a)中所描述的規(guī)律剛好相反。

圖4 延遲時(shí)間(τ)對(duì)振動(dòng)幅值的影響Fig.4 The effect of delay time (τ) on the vibration amplitude

圖5 延遲時(shí)間(τ)對(duì)振動(dòng)頻率的影響Fig.5 The effect of delay time (τ) on the vibration frequency

圖6 延遲增益kd對(duì)振動(dòng)幅值的影響Fig.6 The effect of delay gain (kd) on the vibration amplitude

為了探究柔性圓柱在延遲反饋?zhàn)饔孟碌哪芰哭D(zhuǎn)化，根據(jù)高云等的研究引入瞬時(shí)能量轉(zhuǎn)化公式

(9)

式中，v(z,t)，y(z,t)和q(z,t)分別為瞬時(shí)速度、位移和尾跡變量。瞬時(shí)能量W為正值表示流場(chǎng)向結(jié)構(gòu)傳遞能量，從而激發(fā)柔性圓柱振動(dòng)；負(fù)值表示結(jié)構(gòu)向流場(chǎng)傳遞能量。

圖7(a)所示為柔性圓柱中間節(jié)點(diǎn)處系統(tǒng)振動(dòng)變量的值，可以發(fā)現(xiàn)未施加延遲反饋時(shí)變量q的最大值為4.821，當(dāng)尾跡變量q增加時(shí)，振動(dòng)位移y也相應(yīng)增加。對(duì)于振動(dòng)幅值較小的情況(τ=6.5)開(kāi)始振動(dòng)時(shí)間要小于振動(dòng)幅值較大的情況(τ=3.2)，而且都要小于未施加延遲反饋時(shí)的開(kāi)始振動(dòng)時(shí)間，說(shuō)明施加延遲反饋控制能夠降低系統(tǒng)的起振時(shí)間。圖7(a)中瞬時(shí)能量W變化曲線表明，W的最小值越小，對(duì)應(yīng)產(chǎn)生的振動(dòng)幅值越大。例如對(duì)于τ=3.2而言，此時(shí)W最小值為-0.07，要小于τ=6.5時(shí)的-0.002，但前者對(duì)應(yīng)的振動(dòng)幅值y更大。結(jié)合瞬時(shí)能量W的含義，負(fù)值越小表明結(jié)構(gòu)傳遞到流場(chǎng)中的能量越多，流場(chǎng)中能量增加會(huì)導(dǎo)致尾跡變量q值增大，從而又會(huì)反過(guò)來(lái)增加結(jié)構(gòu)的振動(dòng)y，因此通過(guò)能量云圖7(b)可以看出，當(dāng)τ=3.2時(shí)，由于能量在結(jié)構(gòu)和流體中相互傳遞，駐波向行波演變的時(shí)間會(huì)增加(白色虛線框?qū)?yīng)的時(shí)間)，對(duì)于抑制振動(dòng)的情況(τ=6.5)，駐波向行波演變時(shí)間會(huì)減小。

上述可以發(fā)現(xiàn)瞬時(shí)能量W最小值對(duì)系統(tǒng)振動(dòng)幅值增加或者減小有很大的影響。而結(jié)構(gòu)振動(dòng)增強(qiáng)的能量來(lái)源是外部延遲反饋控制系統(tǒng)，控制系統(tǒng)的延遲時(shí)間τ決定系統(tǒng)當(dāng)前振動(dòng)狀態(tài)與τ時(shí)間之前振動(dòng)狀態(tài)之間的耦合關(guān)系，這種耦合關(guān)系會(huì)降低或者增加結(jié)構(gòu)振動(dòng)幅值。因此延遲增益kd決定了系統(tǒng)振動(dòng)幅值的上限(振幅峰值)，延遲時(shí)間τ決定系統(tǒng)當(dāng)前振動(dòng)狀態(tài)與之前振動(dòng)狀態(tài)的耦合關(guān)系，結(jié)合圖4(a)的結(jié)論，可以看出延遲時(shí)間τ并不能改變振幅峰值，峰值大小只與延遲增益kd有關(guān)。

圖7 當(dāng)kd=0.15時(shí)，延遲時(shí)間(τ)對(duì)系統(tǒng)的影響Fig.7 The effect of delay time (τ) on VIV system when kd=0.15

圖8(a)和圖8(b)中的延遲時(shí)間分別是圖5(a)頻率取最大和最小時(shí)所對(duì)應(yīng)的延遲時(shí)間。從圖8(a)中可以看出，當(dāng)kd由-0.15逐漸增加到0.15時(shí)，對(duì)于延遲時(shí)間τ=1.4，τ=7.4和τ=13.4而言，振動(dòng)頻率隨kd的增加而增加；當(dāng)τ=19.4和τ=25.6時(shí)，振動(dòng)頻率先保持不變，kd進(jìn)入正值以后頻率逐漸減小，在kd=0.12時(shí)頻率又突然增大，最終值和其他延遲時(shí)間保持相同。對(duì)于圖8(b)而言，振動(dòng)頻率隨著kd增大都呈現(xiàn)逐漸減小的趨勢(shì)。

圖8 延遲增益kd對(duì)振動(dòng)頻率的影響Fig.8 The effect of delay gain (kd) on vibration frequency

從圖9可以發(fā)現(xiàn)，未施加延遲反饋控制時(shí)，振動(dòng)幅值為0.647，行波傳遞周期為6.41，振動(dòng)頻率為0.156(1/6.41)；當(dāng)τ=4.7時(shí)，振動(dòng)幅值增大為0.899，行波傳遞周期為6.71，振動(dòng)頻率減小為0.149；當(dāng)τ=1.4時(shí)，振動(dòng)幅值減小為0.521，行波傳遞周期為5.95，振動(dòng)頻率增大為0.168。相較于未施加延遲的情況，振動(dòng)幅值和振動(dòng)頻率變化呈相反的趨勢(shì)，振動(dòng)幅值增大，振動(dòng)頻率就會(huì)減小。

圖9 振動(dòng)位移和行波頻率時(shí)間的關(guān)系(kd =0.15)Fig.9 The relationship between vibration displacement and frequency of travelling wave (kd =0.15)

3.3 理論解釋

對(duì)于耦合方程式(3)和式(4)，設(shè)方程組的解為諧波行波的形式，寫作

y(z,t)=y0ei(kz-ωt-φ),q(z,t)=q0ei(kz-ωt)

(10)

式中：結(jié)構(gòu)位移y和流體變量q都具有共同角頻率ω以及波數(shù)k；φ為變量之間的相位角，y0和q0為與時(shí)間無(wú)關(guān)的幅值。將式(10)代入式(3)和式(4)中，得到

(11)

式(11)對(duì)角項(xiàng)分別為結(jié)構(gòu)和流體的色散關(guān)系，寫作

(12)

(13)

為了得到實(shí)數(shù)角頻率ω和波數(shù)k，式(11)的行列式為0，由此得到流固耦合系統(tǒng)的色散關(guān)系為

DFS(ω,k;τ,q0)=

DS(ω,k;τ)DF(ω,k;q0)+AMω2=0

(14)

令式(14)實(shí)部和虛部為0，分別得到

-(ω2+c2k2)(-ω2+1)-

(15)

(16)

分析式(15)和式(16)可以發(fā)現(xiàn)，振動(dòng)角頻率ω、波數(shù)k以及尾跡變量q0與延遲參數(shù)kd和τ有關(guān)。首先不考慮延遲參數(shù)，得到角速度ω和q0隨波數(shù)k的變化，如圖10所示。可以發(fā)現(xiàn)，有無(wú)數(shù)個(gè)(ω，k，q0)組合滿足上述兩個(gè)方程，特定的組合與給定的邊界條件有關(guān)。對(duì)于無(wú)限長(zhǎng)柔性圓柱而言，最大的振動(dòng)幅值y0和尾跡變量q0會(huì)發(fā)生在ω=1處，此時(shí)波數(shù)k=0.25，q0=4.89。

確定的尾跡變量q0之后，根據(jù)式(15)和式(16)，ω和k與延遲參數(shù)相關(guān)。從圖11中可以發(fā)現(xiàn)，隨著τ的增加，頻率f=ω/2π逐漸趨向于0.159，這一變化規(guī)律與圖5中的變化規(guī)律一致，同時(shí)波數(shù)k也逐漸趨向于0.25。通過(guò)圖11(c)可以看出，當(dāng)τ較大時(shí)，ω/k的值基本保持為4不變，這與Facchinetti等研究中的結(jié)論ω/k=c相一致。通過(guò)上述分析發(fā)現(xiàn)，當(dāng)延遲時(shí)間τ不斷增加時(shí)，振動(dòng)頻率f和波數(shù)k都趨向于未受控時(shí)的值，分別為0.159和0.25。

圖10 振動(dòng)頻率ω和尾跡變量q0隨波數(shù)k的變化趨勢(shì)Fig.10 The variation of vibration frequency ω and wake variable q0 with wave number k

圖11 頻率ω、波數(shù)k以及ω/k隨延遲時(shí)間τ變化規(guī)律(kd=0.15)Fig.11 The variation of frequency ω,wave number k,and ω/k with delay time τ (kd=0.15)

4 結(jié) 論

本文以三維細(xì)長(zhǎng)柔性圓柱為研究對(duì)象，基于結(jié)構(gòu)梁振動(dòng)模型和尾流振子模型，利用二階中心有限差分方法求解耦合方程，通過(guò)引入時(shí)滯反饋控制，討論了在不同的反饋參數(shù)下(延遲時(shí)間τ和延遲增益kd)柔性圓柱渦激振動(dòng)的幅值和頻率變化情況。主要得出如下結(jié)論：

(1)探究了延遲時(shí)間τ對(duì)系統(tǒng)振幅和頻率的影響，發(fā)現(xiàn)三維柔性圓柱振動(dòng)幅值和頻率隨τ的增加呈周期變化；相反的延遲增益會(huì)使振幅變化曲線相差半個(gè)周期的相位角，而且系統(tǒng)能夠達(dá)到的振幅峰值隨著kd的增大而增大。

(2)探究了在特定τ下kd對(duì)系統(tǒng)振幅和頻率的影響，發(fā)現(xiàn)延遲增益kd決定了系統(tǒng)振動(dòng)幅值的上限；且振動(dòng)頻率和幅值的變化呈相反的趨勢(shì)，振動(dòng)幅值增大，振動(dòng)頻率越小。

(3)通過(guò)理論分析得出流固耦合系統(tǒng)的色散關(guān)系，發(fā)現(xiàn)振動(dòng)頻率ω和波數(shù)k隨著延遲時(shí)間τ的增加最終會(huì)保持不變；振幅峰值只與延遲增益kd有關(guān)，kd越大振幅峰值越大。