晏致濤，李正良，楊振華

(1.重慶大學 土木工程學院，重慶 400045;2.山地城鎮(zhèn)建設與新技術教育部重點實驗室(重慶大學)，重慶 400045)

覆冰導線在風荷載作用下可能作低頻大幅自激振動，這種現(xiàn)象稱為舞動。在2009年底至2010年初，河南、河北、湖北、湖南、山東、遼寧等13個省輸電線路發(fā)生大面積導線舞動現(xiàn)象。舞動發(fā)生風速一般較低，幅值會達到垂度的數(shù)倍，所以舞動形成后對輸電線路造成極大的破壞作用，導致電線相互碰撞、絕緣子損壞、金具破壞，嚴重時甚至使鐵塔倒塌，發(fā)生大面積電力供應癱瘓［1］。

舞動發(fā)生時，導線會產生豎直、水平和扭轉等復雜耦合振動。最早的輸電線舞動分析理論是20世紀30年代 DenHartog［2］提出的垂直運動舞動理論，隨后，Nigol提出了扭轉舞動理論［3，4］，Haaker［5］用單(扭轉)自由度模型進行了進一步分析，這幾種舞動理論均為基于單自由度體系的解析分析。之后，一些學者進行了2自由度的導線舞動模型分析，如包括豎向和扭轉運動的導線舞動的2自由度模型［6，7］及考慮了水平運動和仰俯運動的共同作用的2自由度模型［8］。在實際的野外觀測中，導線的舞動是包含水平、豎直以及扭轉等各個方向的振動?；谶@種現(xiàn)象，Yu［9-11］建立了包含2個平動和1個扭轉自由度的3自由度舞動分析模型。Desai［12］等人采用有限元的思想，提出了用3結點等參索單元建立結點3自由度模型研究導線的舞動。

基于3自由度的舞動分析模型得到了極大的發(fā)展，許多研究者結合具體工程實例進行了分析［14-18］。該類分析多數(shù)基于傳統(tǒng)的索單元理論［19］。由于索單元理論一般均忽略抗彎剛度、抗扭剛度和抗剪剛度，多數(shù)學者均添加一個抗扭剛度來表征輸電線的抗扭能力，采用McConnell［20］提出的耦合系數(shù) BT考慮軸向和扭轉的耦合剛度，在舞動分析的時候大都將該耦合系數(shù)取值為零。然而，Luongo［21，22］認為采用上述耦合系數(shù)簡單地考慮軸向與扭轉的耦合剛度是不合理的，并通過解析方法說明抗彎剛度對扭轉剛度的影響不可忽略，全面考慮抗彎、抗拉剛度以及抗扭剛度是非常有必要的。

在Zhu［23］提出的空間曲梁插值函數(shù)的基礎上，筆者考慮了導線運動過程中的扭轉非線性以及覆冰產生的偏心慣性作用，提出了結點6自由度覆冰輸電線舞動分析模型，建立了模擬覆冰導線舞動的基于更新Lagrange格式的有限元方程，并進行了求解。由于考慮了抗彎剛度、抗拉剛度與抗扭剛度的影響，該模型能更加準確地描述輸電線的各自由度耦合舞動行為，有助于闡明導線舞動的發(fā)生機理。

1 覆冰導線有限元模型

輸電線分析模型及坐標系如圖1(a)所示，Kx、Ky、Kz分別表示輸電線兩端絕緣子和輸電塔在x，y，z方向上的邊簡化界條件。圖1(b)為輸電線的斷面形式，由于覆冰后，導線截面的質心移至P點。圖1(c)為空間坐標系下的3結點覆冰輸電線索梁單元。每個單元是位于局部隨轉坐標系，每個結點有6個自由度——x，y，z方向的平動和轉動。

圖1 輸電線及3結點覆冰導線模型Fig.1 Model of transmission lines and three-node iced cable

1.1 位移插值函數(shù)

3結點6自由度覆冰輸電線的有限元模型位移插值函數(shù)可以參考曲梁理論［23］，按多項式構造，可得t時刻單元位移向量:

式中，κ為單元的曲率，下劃線為含曲率的項。若曲率為零，則該插值函數(shù)對應普通直梁單元。根據(jù)圖1(c)中坐標系及符號表示，令結點位移ue={u11u21u31θ11θ21θ31u12u22u32θ12θ22θ32u13u23u33θ13θ23θ33}T，單元內部任一點位移為u={u1u2u3θ1θ2θ3}T，則式(1)可以寫成位移插值的矩陣形式，有:

式中，A為式(1)多項式插值的系數(shù)矩陣，Tr可由位移插值函數(shù)代入各結點位移求得，顯然，3結點6自由度單元模型的型函數(shù)為ATr。

1.2 曲率—位移關系

在圖1(c)中隨轉坐標系下，忽略索的剪切變形，考慮軸向、兩個方向的彎曲以及扭轉變形，其應變-位移關系為:

式中:下標中“，”表示求導;ε為軸向應變;φ1為局部坐標系下u1方向彎曲曲率;φ2為局部坐標系下u2方向彎曲應變;η為局部坐標系下u3方向扭轉應變;下劃線部分表示應變的非線性部分。從應變可見，扭轉應變與軸向變形也有關。

2 覆冰導線舞動的有限元方程

根據(jù)Hamilton原理，并考慮輸電線的阻尼力，可以得到運動方程:

式中tM、tK、tC為t時刻的質量、剛度與阻尼矩陣，t+ΔtF 為t+Δt時刻的荷載向量。

2.1 單元質量矩陣

圖1(b)所示的覆冰輸電線上的一個橫截面在局部坐標系下x，y，z方向的質量密度矩陣為μ，則空間3結點覆冰輸電線單元的一致質量矩陣為:

式中，上標“e”表示單元，Ar為截面面積;Sy為截面對y軸的靜矩;Iy為截面對y軸的慣性矩;Iz為截面對z軸的慣性矩;J為截面對形心的極慣性矩。質量密度矩陣μ中，Sy和Sz以及極慣性矩J可以沿輸電線長度s方向改變。

2.2 單元剛度矩陣

在局部坐標系下建立單元剛度矩陣Ke:

式中，tE、tG分別為t時刻的彈性模量和剪切模量，矩陣A括號里的數(shù)字表示行向量。

式中，ΔT、ΔM1、ΔM2、ΔM3分別為軸力及力矩的線性增量，可以通過線性剛度表達式(8)乘以變形增量求得。從式中可見，軸力、彎矩以及扭矩是相互耦合的，并且扭矩是非線性的。

式中，g為重力加速度，Sz為覆冰輸電線對z軸的靜矩;diag為取對角陣符號。

2.3 單元阻尼矩陣

采用瑞利阻尼，表達式為:

式中，αv、βv為瑞利阻尼系數(shù)。

2.4 單元荷載向量

作用在導線上的引起舞動的空氣動力與作用在其他結構上的常見荷載不同。氣動力依賴于覆冰導線的幾何非線性和風攻角α。通過覆冰導線的風洞試驗，測量得到的對應于攻角α的阻力系數(shù)C1(α)、升力系數(shù)C2(α)和扭矩系數(shù)C3(α)，進而計算得到升力F1、阻力F2和扭矩M3。單元荷載向量Fe可表示為:

式中，F(xiàn)1、F2、F3分別為局部坐標系三個方向上的結點力;M1、M2、M3分別為局部坐標系三個方向上的彎矩。ρair為空氣密度;d為裸導線的直徑;Vrel為相對風速，一般取Vrel≈Vw，Vw為作用在導線上的平均風速。

2.5 坐標轉換矩陣

上述單元質量矩陣、單元剛度矩陣以及荷載列陣等是在局部隨動坐標下推導出來的，在求解整體方程時，需要將它轉換到整體坐標。從整體坐標矩陣至局部坐標矩陣的轉換矩陣可以分兩步完成，第一步將全局坐標系的x坐標軸旋轉至輸電線的弦長方向，y方向旋轉至輸電線的豎向，z方向旋轉至輸電線的側向;第二步將輸電線各個結點x，y，z方向旋轉成圖1(c)中的1、2、3方向，中間結點可以不旋轉。

2.6 邊界條件

在一跨導線檔距兩端邊界處，既有遠跨導線連接，又有絕緣子連接。邊界處的剛度系為相應邊界剛度與絕緣子剛度的疊加。其中x方向邊界剛度系數(shù)為x方向剛度系數(shù)為絕緣子x方向擺動剛度系數(shù)與遠跨導線x方向剛度系數(shù)相疊加。y方向剛度系數(shù)為絕緣子y方向擺動剛度系數(shù)，一般假定為無窮剛性。z方向剛度系數(shù)為絕緣子z方向擺動剛度系數(shù)與遠跨導線z方向剛度系數(shù)相疊加。上述剛度組裝疊加過程一直進行到耐張塔為止，即耐張塔處結點是鉸接的。

3 動力平衡方程的求解

求解覆冰輸電線非線性舞動動態(tài)問題的Newmark-β法的步驟歸納如下:

(1)計算Newmark-β法中a0-a7等8個積分參數(shù)。根據(jù)新的位形計算覆冰輸電線單元的內力tT、tM1、tM2及tM3，根據(jù)式(8)～ 式(11)確定t時刻結構的剛度矩陣tK和新的有效剛度矩陣tK^:

(2)形成有效載荷向量

(3)求解位移增量

(4)在每一時間步內進行動力平衡迭代(Newton-Raphson平衡迭代)，求解第i次迭代位移增量的修正量，并計算修正后的位移增量:

檢查收斂性，若收斂，即得到本時間步最終的位移增量t+Δtu，進而求得 t+Δt時刻的位移、速度以及加速度。

4 算例

算例1 為了說明結點6自由度模型分析輸電線的正確性，編制程序，分析了文獻［19］中的一根輸電線。該輸電線為一小垂跨比輸電線，輸電線質量0.821 7 kg/m，剛度為 EI為 11.07 Nm2，輸電線的軸力為15 860 N，輸電線長為13.385 m。表1表示文獻［19］中所列出 Irvine 方法［5］、Desai方法［12］、試驗值與本文理論分析該輸電線的模態(tài)分析結果比較。結果表明，對于上該輸電線，由于垂跨比較小，上述4種理論均與試驗吻合較好。由于曲梁理論推導沒有小垂跨比假設，因此，結點6自由度模型可以適用于更廣泛的輸電線分析。

表1 不同理論模型分析頻率(Hz)Tab.1 Mode analysis with different compute model(Hz)

表2 單根導線的物理參數(shù)Tab.2 Physical parameters of single conductor

算例2 上述算例表示結點6自由度模型計算線性模態(tài)是正確的，為了考察其非線性舞動特性，這里采用6自由度模型計算文獻［12］的一個經典覆冰導線舞動算例。算例物理參數(shù)如表2，氣動參數(shù)可以參考原文獻。該導線是覆冰形狀為D型截面，有野外舞動觀測的數(shù)據(jù)，現(xiàn)有眾多的數(shù)值分析文獻大都利用該野外觀測數(shù)據(jù)進行驗證，遺憾的是，野外原始數(shù)據(jù)中缺少了兩個關鍵的截面參數(shù)Sy、Sz。文獻［12］假設了這兩個參數(shù)，分析結果與野外數(shù)據(jù)吻合較好，本文取值按表2。利用結點6自由度單元進行該導線的舞動分析，編制計算程序，采用Newmark-β法解動力方程。得到初始攻角為100、風速4.1 m/s時的導線舞動軌跡，如圖2～圖4所示。

圖2 本文計算結果與文獻［12］比較Fig.2 Comparison between reference［12］and this paper

將計算結果與文獻［12］中的計算結果［圖2(a)］進行對比。［圖2(a)］的中心點位置為(0.04，1.38)，圖2b)的中心點位置為(0.003，0.005)，這是由于本文圖示位置是基于平衡位置;計算y方向的最大振幅為0.79 m，在側向的最大振幅為0.02 m，幅值平面的夾角為 0.8°，計算結果有效。

圖3～圖4分別為輸電線舞動時中點的豎向位移及扭轉位移時程曲線，計算初值取豎向振幅為0.001m，計算步長取0.001 s，大約在2 500 s振動趨于穩(wěn)定。由于扭轉剛度較低，該導線產生較大的扭轉角，屬于扭轉大變形。傳統(tǒng)的3自由度舞動分析模型在處理扭轉角時沒有考慮抗彎剛度的影響，并且扭轉變形是基于小變形、線彈性假設推導的。從圖4可以看出，上述假設并不成立，需要考慮扭轉非線性的影響。由于本例在4.1m/s風速下發(fā)生的是頻率比1∶1∶1(側向頻率:豎向頻率:扭轉頻率)的共振，因此，在舞動發(fā)生時，豎向振動、側向振動以及扭轉角振動的鎖定頻率是一致的，振幅是相互耦聯(lián)的。從時程分析計算可知，發(fā)生舞動與初始條件無關，是一種自激振動。

圖5 不同抗彎剛度對舞動豎向位移的影響Fig.5 Vertical galloping displacement of different bend stiffness

圖6 不同抗彎剛度對舞動側向位移的影響Fig.6 Lateral galloping displacement of different bend stiffness

圖7 不同抗彎剛度對舞動豎向位移的影響Fig.7 Torsional galloping angle of different bend stiffness

圖5～圖7表示假設該導線不同的截面抗彎剛度采用結點6自由度模型分析輸電線舞動動態(tài)振幅和靜態(tài)偏移。導線參數(shù)和氣動力參數(shù)仍然如表2，只是原文獻缺乏導線的抗彎剛度數(shù)據(jù)，這里根據(jù)原文獻的截面形式，假使導線截面均勻分布，確定導線橫截面上各個方向的抗彎剛度EI為294 N·m2，圖中橫坐標n表示EI的倍數(shù)。實心點曲線表示舞動的動態(tài)振幅，空心點曲線表示舞動的靜態(tài)偏移。圖5表明，隨抗彎剛度增加，豎向舞動振幅及靜態(tài)偏移均增大;圖6表明，隨抗彎剛度增加，側向舞動振幅及靜態(tài)偏移均減小;圖7表明，隨抗彎剛度增加，豎向舞動振幅減小，但靜態(tài)偏移變化不大。可見，由于覆冰輸電線抗彎剛度與抗扭剛度屬同一個數(shù)量級，其對輸電線舞動的扭轉有較大影響，又因為根據(jù)式(13)所推導的準定常氣動力產生了與平動相關的耦聯(lián)項，因此對輸電線舞動的平動振幅也有較大影響。

5 結論

傳統(tǒng)的輸電線3自由度舞動分析模型無法考慮扭轉與平動的耦合作用及扭轉運動的非線性。利用3結點6自由度曲梁理論求解輸電線的舞動問題可以合理地反映上述關系。基于6自由度曲梁的應變—位移關系，建立了結點包含具有3個平動自由度和3個轉角自由度的輸電線舞動分析有限元模型?？紤]覆冰導線所受空氣動力的非線性和導線大幅運動的幾何非線性，利用虛功原理建立基于更新Lagrange格式的覆冰導線非線性運動方程。采用 Newmark時間積分和Newton-Raphson非線性迭代法求解有限元方程。算例分析表明，輸電線的抗彎截面模量對輸電線的平動和扭轉有較大的影響?，F(xiàn)有的舞動機理如垂直運動激發(fā)理論、扭轉激發(fā)理論以及慣性偏心激發(fā)理論實際上是基于單自由度、2自由度或3自由度的理論，是6自由度耦合運動的特例或近似。

［1］郭應龍，李國興，尤傳永.輸電線路舞動［M］.北京:中國電力出版社，2003.

［2］ Den Hartog J P.Mechanical vibration［M］.New York:Mcgraw-Hill，1956.

［3］Nigol O，Buchan P G.Conductor galloping partⅠ-den hartogmechanism ［J］. IEEE Transactions on Power Apparatus and System，1981，100(2):699-707.

［4］Nigol O，Buchan P G.Conductor galloping partⅡ-torsional mechanism［J］.IEEE Transactions on Power Apparatus and System，1981，100(2):708-720.

［5］Haaker T I，Vanoudheusden B W.One-degree-of-freedom rotational galloping understrong wind conditions［J］.International Journal of Non-linear Mechanics，1997，32(5):803-814.

［6］ Byun G S，Egbert R I.2-degree-of-freedom analysis of power-line galloping by describing function methods［J］.Electric Power Systems Research，1991，21(3):187-193.

［7］Yu P，Shah A H，Popplewell N.Inertially coupled galloping of iced conductors［J］.Journal of Applied Mechanicstransactions of the ASME，1992，59(1):140-145.

［8］ Jones K F.Coupled vertical and horizontal galloping［J］.Journal of Engineering Mechanics，1992，118(1):92-107.

［9］Yu P，Desai Y M，Shah A H.3-degree-of-freedom model for galloping-1. formulation［J］. Journal of Engineering Mechanics-ASCE，1993，119(12):2404-2425.

［10］Yu P，Desai Y M，Popplewell N.3-degree-of-freedom model for galloping-2. solutions［J］. JournalofEngineering Mechanics-ASCE，1993，119(12):2426-2448.

［11］ Yu P，Popplewell N，Shah A H.Instability trends of inertially coupled galloping:periodic vibrations［J］.Journal of Sound and Vibration，1995，183(4):679-691.

［12］ Desai Y M，Yu P，Popplewell N，et al.Finite element modeling of transmission line galloping［J］.Computers and Structures，1995，57(3):407-420.

［13］ Desai Y M，Yu P，Shan A H，et al.Perturbation-based finite element analyses of transmission lines galloping［J］.Journal of Sound and Vibration，1996，191(4):469-489.

［14］樊社新，何國金，廖小平，等.結冰導線舞動機制分析［J］.中國電機工程學報，2006(14):131-133.

［15］范欽珊，官 飛，趙坤民，等.覆冰導線舞動的機理分析及動態(tài)模擬［J］.清華大學學報(自然科學版)，1995(2):34-40.

［16］胡 杰，郭應龍，惲俐麗.三分裂輸電導線舞動的有限元分析與計算［J］.水利電力機械，1996(4):9-12.

［17］劉小會，嚴 波，張宏雁，等.覆冰導線舞動非線性數(shù)值模擬方法［J］.應用數(shù)學和力學，2009，30(4)，:457-468.

［18］Wang J W，Lilien J L.Overhead electrical transmission line galloping-a full multi-span 3-DOF model，some application and design recommendations［J］.IEEE Transaction on Power Delivery，1998，13(3):909-915.

［19］Irvine H M，Cable structures［M］.New York:MIT Press，1981.

［20］ Mc Connell K G，Chang C N.A study of the axial-torsional coupling effect on a sagged transmission line［J］.Experimental Mechanics，1986，26(4):324-329.

［21］Angelo L，Zulli D ，Piccardo G.Analytical and numerical approaches tononlineargallopingofinternally resonant suspended cables［J］.Journal of Sound and Vibration，2008，(315):375-393.

［22］Angelo L，Zulli D ，Piccardo G.On the effect of twist angle on nonlinear galloping of suspended cables［J］.Computers and Structures，2009，87(15-16):1-12.

［23］Zhu Z H，Meguid S A.Vibration analysis of a new curved beam element［J］.Journal of Sound and Vibration，2008，309(1-2):86-95.