雙索非線性振動瞬時相位差及其來源

孫測世，林俊強，鄧正科

(1.重慶交通大學 土木工程學院，重慶 400074；2.湖南城市學院設計研究院有限公司，長沙 410005)

工程中索結構的振動危害極大，一直是人們關注的重點問題[1-3]。從相位來看，各索若處于同步振動狀態(同相)，即拉索的慣性力方向一致時，拉索對錨固結構的作用力瞬間增大數倍，可造成錨固結構嚴重受損；若各索處在異步振動狀態(異相)，則極易引起鄰索間的碰撞[4-6]，導致索結構的使用壽命大大縮減。

目前，對多索結構的研究已表明，參數差異對相鄰拉索的相對振動存在較大影響。孫測世[7]進行了斜拉橋全橋模型試驗，發現當激勵頻率達到一定值時多根索會發生大幅振動，且相鄰索的振動間均存在一定相位差，導致產生碰撞。孟新田[8]建立了斜拉橋的多索-梁動力學簡化模型，分析了雙索間在倍頻關系下拉索的大幅振動問題。趙躍宇等[9]建立了一種雙索質量塊模型，分析了拉索與拉索之間在振動過程中的相互影響，研究發現單根拉索參數的變動對相鄰拉索振動的影響較大。國內外學者對于拉索的相對振動的研究在模態分析中往往以“同相”與“反相”模態體現，Lepidi等[10]采用一個具有雙索的簡化模型研究了索承重橋梁的模態耦合，研究表明拉索存在同相和反相的局部模態。Ahmad等[11]建立了雙索網絡結構，研究了輔助索剛度對雙索網絡平面內自由振動特性的影響，發現雙索網結構在平面內自由振動中存在同相和反相模態。Abdel-Ghaffar等[12]首次研究了拉索振動過程對斜拉橋振動產生的影響，研究表明MECS模態分析中可獲得大量拉索局部模態，其中包括相鄰拉索的同相和反相模態。吳慶雄等[13]對多索-梁結構的研究中發現二索-梁結中也存在同相和反相模態，并且在此基礎上建立了多索梁模型。Zhou等[14]提出并分析了在錨索錨固附近布置有阻尼器的兩根平行拉索系統，研究發現第三模態下兩平行索之間振幅同相，第四模態下兩平行拉索振幅反相。

對其他類似結構的研究中也表明相頻特性隨結構參數變化。最為典型的結構為擺，根據參數的差異，其相位角、相位差會有不同的影響和變化，會呈現同步或異步運動的現象[15-17]；Kapitaniak等[18]通過在平臺上兩個平面彈性擺加水平激勵，發現了穩定的同異步運動狀態，并分析了同異步振動下的分岔情形。

綜上，索的相頻特性與物理參數和激勵密切相關，相鄰拉索會因參數關系改變而出現同步或異步振動。本文對雙水平索質量塊模型展開對索的瞬時相位差研究，揭示不同參數下高階近似項中引起雙水平索非線性振動瞬時相位差的因素及其影響規律。

1 力學模型與多尺度分析

1.1 基本假定

如圖1所示雙水平索與質量塊耦合振動力學模型，其中質量塊M僅考慮沿索軸向的運動X(t)，不考慮質量塊阻尼。其中，跨度為l，跨中垂度為d。索的單位長度質量m，以索的軸向方向為x方向，重力加速度g的方向為y方向。x、y方向對應位移分別用u、v表示。

圖1 雙水平索與橋面耦合力學模型Fig.1 The coupling mechanical model of double horizontal cable and bridge deck

其它基本假設包括：①不考慮索和梁的材料非線性；②不計拉索的抗彎剛度、抗扭剛度以及抗剪剛度；③認為兩根拉索的重力垂度均為拋物線；④拉索的本構關系服從胡克定律并且各點受力均勻；⑤不計橋塔振動對索的影響；⑥只考慮橋面的豎向運動。

1.2 振動方程的建立

根據Hamilton原理，拉索的動力學方程為

(1)

其中

式中：m為單位長度質量;E為拉索彈性模量;A為拉索截面面積;εd為拉索振動中的動應變。y(x)水平索的靜態構型，其表達式可為：

y(x)=4f(1-x)x

(2)

式中，f為懸索垂跨比。

當引入質量塊位移函數X(t)和阻尼項后可得到

(3)

引入如式(4)中的無量綱變換

(4)

(5)

類似地，將質量塊的運動方程進行無量綱化后可得

(6)

1.3 離散化模型

懸索的振動可以表示為各階模態振動的組合，因此可利用分離變量法令

(7)

(8)

式中:v1(x,t)、v2(x,t)分別表示索1和索2的動位移；qn(t)、Bn(t)為索1和索2振動的廣義時間坐標；Φ(x)為索1和索2振型函數，其正對稱模態函數Φ1(x)和反對稱模態函數Φ2(x)如下[19]

(9a)

(9b)

利用Galerkin方法進行離散，可將包含時間二階項和空間二階項無量綱偏微分方程簡化為只含有時二階項的無量綱離散控制方程，將相應的量代入，聯合各式并簡化系數，則可得相應雙水平索質量塊耦合振動方程組如下

(10a)

(10b)

(10c)

式中，各系數如下：

1.4 攝動分析

基于多尺度法對式(10a)～(10c)求近似解，令：

q(t;ε)=εq0(T0,T1,T2)+ε2q1(T0,T1,T2)+

ε3q2(T0,T1,T2)+…,

(11a)

B(t;ε)=εB0(T0,T1,T2)+ε2B1(T0,T1,T2)+

ε3B2(T0,T1,T2)+…,

(11b)

Z(t;ε)=εZ0(T0,T1,T2)+ε2Z1(T0,T1,T2)+

ε3Z2(T0,T1,T2)+…,

(11c)

其中，Tn=εnt(n=0,1,2,3,…)，代入式(10a)～(10c)并按ε的冪次整理得到：

ε：

(12a)

(12b)

(12c)

ε2：

(13a)

(13b)

(13c)

ε3：

(14c)

設式(12a)～(12c)的復數解為

qk1=A1eiω1T0+cc

(15a)

Bk1=A2eiω2T0+cc

(15b)

Z1=A3eiω3T0+cc

(15c)

代入式(13a)～(13c)，令長期項為0，可得特解：

(16a)

(16b)

(16c)

式中，各系數如下：

再將式(15a)～(16c)代入式(14a)～(14c)取雙索和質量塊三階長期項為0可得：

其中：

Tem=-am1F4k-am1F3k-3am4+

為便于求解式(17a)～(17c)，假定質量塊僅提供拉索軸向激勵X(t)=QcosΩt，其中Q為質量塊的無量綱振動幅值，為質量塊的實際激勵幅值Qb除以索長l，即：Q=Qb/l。令Amn=1/2aneiβn(n=1,2)，代入式(17)，并設穩態振動時

經分離實、虛部可得平均方程

(18a)

(18c)

(18d)

由式(18)可知，至此問題轉化為四個方程求解a1、a2、γ1和γ2四個未知數的非線性代數方程組，雖然一般可以通過直接簡化方程，研究非平凡穩態解的穩定性，但是由于簡化過程存在耦合項，直接求解比較困難，因此引入直角坐標的方法，做以下轉換

sinγ1=p1m/a1，cosγ1=q1m/a1

sinγ2=p1n/a2，cosγ2=q1n/a2

再利用數值計算方法求解方程組，最后得雙索響應近似解為

2 不同計算方法的對比驗證

為驗證基于多尺度法的計算結果的正確性，本節利用算例將數值結果與Runge-Kutta法直接積分、有限元法計算結果進行對比。在給定參數m1=m2=136 kg/m、l=440 m、EA=2.98×109N、H1=8 000 kN、Q=0.000 1等基礎上，改變另一根拉索的索力參數。為結果具有一般性，選用垂跨比這一在工程中較常用和直觀的無量綱參數進行分析。按f1=0.009，f2=0.011或0.01，考慮兩個不同垂跨比之比的工況，即f2/f1=0.01/0.009=1.11和f2/f1=0.011/0.009=1.22。

利用有限元軟件ANSYS建立雙索質量塊模型，采用LINK10模擬拉索單元。質量塊采用MPC184單元模擬，不考慮阻尼力與彈性恢復力，僅作為一個連接兩根拉索的物體，具有無窮大的剛度和一定的質量，通過給質量塊一個x方向的動位移，使雙索產生參數振動。同時，采用四階Runge-Kutta編制MATLAB程序直接對常微分方程式(10)進行直接數值積分。

下文將以f2/f1=1.11和1.22，Ω=1.05作為算例進行分析，模擬兩組參數下的相對振動，繪制出相應的時程曲線圖與Runge-Kutta法求得的響應以及非線性有限元分析結果進行對比。為方便比較，將有限元的有量綱的數據轉化為無量綱，對比結果如圖2和圖3所示。可見，在參數條件為f2/f1=1.11，Ω=1.05時，索1與索2處于同步振動；在參數條件為f2/f1=1.22，Ω=1.05時，索1和索2處于異步振動。

比較兩組不同拉索參數下的Runge-Kutta法與多尺度法所得的時程曲線，可見兩組不同參數下兩種求解方法的結果吻合良好。有限元分析得到的時程曲線的上波峰到下波峰的距離均有微小的減小，但二倍頻項與漂移項的影響確有微小增大。這是由于進行理論分析與Runge-Kutta法計算時，均忽略了面內外耦合。在進行非線性有限元分析時，給與了面外微小擾動，振動過程面內外有一定的耦合，分散了系統的振動能量，因此振動幅值會有所減小。綜上可見，無論是Runge-Kutta法還是非線性有限元分析，結果都與多尺度法結果有較好的吻合，驗證了多尺度方法求解結果的可靠性。

(a) 多尺度法與R-K法結果對比

3 高階解對相頻特性的影響

為研究高階解對相頻特性的影響，繼續通過第2章算例分別考察高階近似項中的相移值γ1和γ2、漂移項等因素對瞬時相頻特性的影響，進一步闡述不同參數下雙水平索的相對振動變化規律。

3.1 高階近似項對瞬時相位的影響

圖4為f2/f1=1.11，Ω=1.15時的時程曲線。其中，圖4(a)為僅考慮線性近似解雙索響應時程曲線，圖4(b)為考慮ε的前兩階時雙索的時程曲線圖。橫坐標為無量綱時間變量，縱坐標為無量綱響應幅值。

(a) 未考慮高階近似項

對比圖4(a)、圖4(b)可知，未考慮高階近似項時，兩根索在振動過程中瞬時相位差值變化波動不大，且索1受高階項的影響相對較小；在考慮高階近似項時，受到高階近似項中二倍頻項及漂移項的影響，響應曲線線型會產生較大的變化。故可推斷出兩根拉索在相移值為零(γ=0)時，兩根拉索會產生一定的瞬時相位差異，因此，繪制了雙索在f2/f1=1.11時的復平面圖。

如圖5所示，在未考慮高階近似項時，兩根拉索的復平面圖都為圓形，相同時刻對應點均在外圓的半徑上，相角差值為零。當考慮高階近似項時，索1與索2的復平面圖中心點均產生一定的偏移，索2尤為明顯，從而使復平面圖中對應曲線形狀發生變化。兩根拉索的瞬時相位差也會因為高階項的存在出現差異，如圖6所示。為了便于比較，定義無量綱參數為

(20)

式中，Δp為兩根拉索間響應的瞬時相位差。

(a) 未考慮高階近似項

在未考慮高階近似項時，即在相移值γ1和γ2為零時，其瞬時相位差幅值也為零，如圖6中虛線所示。當考慮高階近似項時，受到高階近似項的影響，瞬時相位差會呈現很大的變化幅度，相位差時程曲線成周期變化，在該參數條件下，最大相位差可達到0.2π，如圖6中實線所示。

為更加直觀的分析在不同參數下的瞬時相位差變化情況，定義新的無量綱參數

(21)

式中，Δpmax為Δp的幅值。

取各頻率下兩根拉索間響應的瞬時相位差幅值pmax，繪制pmax-Ω曲線圖，將考慮高階近似解與僅考慮線性近似解時瞬時相位差隨頻率變化曲線圖進行對比，如圖7所示。

由圖7可知，當Ω<1.05時，雙索處于同步或異步振動(圖7(a)中為異步振動，圖7(b)中為同步振動)，兩根曲線幾乎重合，說明高階近似項對雙索間的瞬時相位差無明顯影響；當Ω>1.05時，相位差逐漸變大；當Ω>1.10后，差異明顯增大，此時高階近似項對相位差的影響不可忽略。圖7(b)中，在特定頻率下，高階近似項甚至能使瞬時相位增大0.8π。

(a) f2/f1=1.22

綜上所述，考慮高階近似項后，拉索間的瞬時相位差會有很大的變化，因此，在對索的相頻特性分析時，不能忽略高階近似項的影響。

3.2 高階近似項中各因素對瞬時相位的影響

在3.1節基礎上，為進一步追溯高階近似項中導致其“不能忽略”的因素。將第2章已經驗證的MATLAB程序進行適當修改，使高階項中的相移值γ以及漂移項的影響剝離出來進行分析。

3.2.1 相移值γ對瞬時相位的影響

由式(19a)、(19b)可知，在雙索的近似解中，無論是漂移項還是二倍頻率項均含有相移值γ，可見相移值γ的存在可能對響應會有一定的影響，因此繪制出f2/f1=1.11，Ω=1.15時的瞬時相位時程曲線圖與瞬時相位差曲線圖，如圖8所示，兩者進行比較，研究相移值γ對相頻特性的影響。

由圖8(a)可知，瞬時相位差時程曲線不僅向右平移而且向下平移。向右平移是因為兩根拉索的相移值γ不為零，響應時程與瞬時相位時程均產生一定的右移。而因為兩根索相移值存在一定的差異，從而導致瞬時相位差值整體加上或減去一個定值，所以瞬時相位差時程曲線會產生一定的下移。故圖8(a)中Δt和Δ(Δp)的數值大小取決于兩根拉索在同一頻率下的相移值大小(Δt為瞬時相位差時程曲線幅值往右偏移的時間差，Δ(Δp)為瞬時相位差時程曲線幅值往下偏移的瞬時相位差的差值)。

(a) γ=0和γ≠0的瞬時相位差

考慮了相移值后，如圖8(b)所示，瞬時相位差時程曲線會向右平移，兩曲線在同一時刻上對應兩點的相位差數值大小就是相移值γ。雙索結構由于兩根拉索的參數差異，相移值也可能會具有較大差異，而且瞬時相位時程曲線線型也具有一定的差異，故瞬時相位差也會受到較大的影響。

對于雙索近似解，如式(19a)、(19b)中含有二倍頻項與漂移項，分別繪制出γ=0和γ≠0時各項系數的變化曲線，如圖9所示。

(a) f2=0.01

由圖9可知，高階近似項各項系數中均含有響應幅值a1、a2。在向上掃頻時，由于響應幅值a1、a2增大，各項系數絕對值也相應的增大。將各項參數進行對比，可見僅第一項與第四項在γ=0和γ≠ 0時較大的激勵頻率Ω區域時出現明顯差異。但第一項與第四項的數值大小相對于第二項與第三項為小量，數值上可以忽略不計。因此可知，γ對高階項引起的相位差幅值的影響可以忽略。

3.2.2 漂移項對瞬時相位的影響

從式(19a)、(19b)中可見，響應幅值是漂移項中的參數，因此響應幅值對瞬時相位可能也有一定的影響。根據文獻[20]表明，瞬時相位的變化與漂移項在響應振幅中所占比重有關。所以漂移項會影響瞬時相位數值變化，可通過漂移項在響應振幅中所占比重對瞬時相位差幅值進行定性研究。由圖9可知，可將第一項和第四項的幅值忽略，令式(19a)、(19b)中二倍頻率項的余弦值為1，振動達到振幅，則可定義新的無量綱參數

(22a)

(22b)

式中:D1、D2為索1、索2響應近似解中的漂移項(此處僅考慮絕對漂移量);Av1、Av2為索1、索2響應近似解中不考慮漂移項后的響應幅值，稱為無漂移幅值。β1與β2分別反映索1、索2的漂移項在響應幅值中所占的比重。

取f2/f1=1.11時，β1與β2隨Ω的變化曲線如圖10所示。

圖10 β-Ω曲線Fig.10 The curves of β-Ω

由圖10可知，β1與β2曲線在Ω<1.05時基本重合，且均約等于0，此時高階解中幾乎無漂移項，因此，圖7中對應頻率范圍內兩根拉索相位差很小；在Ω>1.10時β2明顯大于β1，這是因為索力較小的索，漂移項的占比明顯增大，進而導致圖7中該頻率范圍內的瞬時相位顯著增大。

綜上所述，高階近似項引起的相位差幅值大小受γ的影響不大，而高階項中的漂移項的變化對瞬時相位的影響卻有決定性作用，因此可利用β值的變化情況對瞬時相位差幅值作定性判斷。

4 結 論

(1) 對比分析表明，Runge-Kutta法、非線性有限元法和多尺度法的結果具有較好的一致性，不同參數下，雙索會呈現同步或異步振動的現象。

(2) 高階近似項對拉索間的瞬時相位差有影響，相頻分析時不能忽略，在一定頻率下，甚至能使瞬時相位增大0.8π。進一步的分析表明：高階近似項中相移值γ對瞬時相位差的貢獻可以忽略，而漂移項對瞬時相位差卻有決定性作用。

(3) 雙索系統的瞬時相位差來源于線性解中相移值γ的差異以及高階項中漂移項占比參數β的差異兩個方面。