復振型疊加法截斷誤差及改進

陳華霆 譚平 彭凌云 李志山 周福霖

摘要： 根據復振型疊加理論，對具有非比例阻尼特性的線性體系，分析了復振型截斷誤差的來源，給出了誤差公式。與比例阻尼體系的高階實振型響應類似，通過反應比表達式可以看出隨著輸入頻率與自振頻率比的減小，高階復振型的響應也趨于荷載的靜態響應。基于此推導出了與傳統修正方法形式一致的復振型靜力修正方法和振型加速度方法，同時根據復振型退化為實振型的條件（實部與虛部成比例關系）建立了實、復振型修正方法的統一。用數值算例對修正方法進行了驗證，結果表明：復振型修正方法能夠很好地提高振型疊加的計算精度，特別是當結構主要固有頻率大于輸入卓越頻率時修正效果更為顯著；還發現附加阻尼較大時，導致部分高階振型自振頻率減小，在一定程度上會降低修正效果。關鍵詞： 線性振動； 動力分析； 復振型疊加； 靜力修正； 振型加速度方法

中圖分類號： O321； TU311.3文獻標志碼： A文章編號： 10044523（2017）04055608

DOI：10.16385/j.cnki.issn.10044523.2017.04.005

引言

對于線性振動體系，一般采用振型疊加方法進行動力分析。通常振型疊加法[1]利用較少的無阻尼（實）振型表達結構的位移模式，通過坐標轉換得到規模較小的用振型廣義坐標表示的運動方程。這種方法對于位移響應有較好的近似，但對于內力和應力誤差較大。為得到較精確的結果往往需要計算大量的振型，以確保充分反映了荷載的空間與頻率分布。為了考慮被截斷的高階振型的影響，Williams[2]提出了振型加速度方法。由于公式中含有振型坐標的加速度項，故由此得名。振型加速度方法還有一些等價的形式[35]，這些方法的共同特點是高階振型的響應被等效成荷載的靜態響應，但前提是荷載的最高頻率與振型的自振頻率比至少在0.5以內[6]。當然，為了振型加速度方法的有效性還需要荷載的空間分布向量中含有高階振型的成分。

以上是針對比例阻尼或經典體系提出的振型疊加修正方法。與傳統的振型疊加方法一樣，非比例或非經典阻尼體系的復振型疊加法也存在振型截斷的修正問題。Kulkarni and Ng[7]，Borino and Muscolino[8]，Aveni and Muscolino[9]均討論過非比例阻尼體系振型疊加的修正問題，他們是在無阻尼振型空間中進行的修正，其修正效果受阻尼矩陣非比例特性大小的影響。Traill Nash[10]在復振型空間中推導了復振型加速度方法的表達式，但式中含有荷載關于時間的一階導數項，在實際應用中不方便。在國內，主要有周錫元[1112]研究了非比例阻尼體系的動力計算問題，但對于復振型疊加截斷的修正問題尚未涉及。

本文簡要敘述了復振型位移方法，即截斷的復振型疊加方法，對復振型截斷誤差的來源進行了分析，并嚴格證明了高階復振型的響應可由荷載的靜態響應代替，基于此推導出了復振型靜力修正方法和振型加速度方法。當結構為比例阻尼或無阻尼體系時，復振型修正方法可退化為實振型修正方法，兩者可以統一起來。文中對上述方法進行了算例驗證，結果表明復振型修正方法能夠提高振型疊加的計算精度，適合在實際工程中應用。

1復振型疊加法

離散為n階自由度的線性結構體系在荷載f（t）作用下的運動方程可表達為M+C+Ku=f（t）（1）式中M，C和K分別為n×n階的質量矩陣、阻尼矩陣和剛度矩陣，其中通常M，K為正定對稱陣；u為n×1階的節點位移向量。

若阻尼矩陣不滿足Caughey解耦條件[13]，則傳統的實振型分解方法對于方程（1）不再適用。然而，這可以通過Foss變換[14]在狀態空間下利用復振型對運動方程進行解耦。

利用恒等式M-M=0，將方程（1）寫為狀態方程的形式A+Bv=p（t）（2）式中A=0M

MC， B=-M0

0K，p（t）=0

f（t），v=

u。

第4期陳華霆，等： 復振型疊加法截斷誤差及改進振 動 工 程 學 報第30卷與方程（2）對應的特征值問題為（λA+B）ψ=0，其中λ，ψ是特征值和相應的特征向量。由于系數矩陣A，B為2n×2n階的對稱非正定矩陣，通常特征值與特征向量為復數，對于低阻尼體系成共軛對出現，即

λi，i=-ξiωi±iωi1-ξ2i，（i=1，2，3，…，n）；ψi=λii

i；i=ii

i

式中ωi，ξi分別為第i階固有圓頻率和振型阻尼比，均為實數；i=-1為純虛數單位；復向量i，i共同稱為第i階復振型，上標“—”表示共軛。本文亦主要對低阻尼體系進行探討。

假設體系的特征向量滿足完備性，則狀態向量可表達為特征向量的線性組合，即v（t）=（t）

u（t）=Ψq（t）（3）式中Ψ=[ψ1…ψn1…n]為特征向量矩陣，q（t）=[q1（t），…，qn（t），1（t），…，n（t）]T為廣義坐標向量。利用特征向量關于矩陣A，B的正交性，可得n組解耦的復振型運動方程i（t）-λiqi（t）=Tif（t）ai

q·-i（t）-ii（t）=Tif（t）i （4）式中ai=ψTiAψi，i=TiAi。上述方程為一階常微分方程，可通過經典的數學方法求解。

利用式（3）可求得全部振型參與計算的結構位移與速度向量u（t）=∑ni=1[iqi（t）+ii（t）]=Φq（t）（5a）

（t）=∑ni=1[λiiqi（t）+iii（t）]=ΦΛq（t）（5b）式中Φ=[1…n1…n]為復振型矩陣，Λ=diag（λ1，…，λn，1，…，n）為譜矩陣。

對于自由度數目較大的大型結構體系，通常取前d階振型參與計算，這可以很大程度地降低計算量。因此，結構的位移和速度近似為

ud（t）≈∑dn=1[iqi（t）+ii（t）]=Φdqd（t）（6a）

d（t）≈∑di=1[λiiqi（t）+iii（t）]=ΦdΛdqd（t） （6b）

式中Φd，Λd與qd（t）分別為前d階振型組成的振型矩陣、譜矩陣和廣義坐標向量。上式即為截斷的復振型疊加方法，也可以稱為復振型位移方法，但振型截斷勢必降低計算的精度，故需要對式（6）的截斷誤差進行估計并進行修正。

2振型截斷誤差

振型截斷對應的位移、速度響應由式（6）確定，它只是方程（1）的近似解。存在一個二階微分運動方程與式（6）對應，將其表達為Md+Cd+Kud=fd（t）（7）式中fd（t）不同于f（t），是對應于位移ud（t）的外荷載。

求fd（t）需要利用振型方程（4），為此需要將式（7）表達成狀態方程的形式。對式（6a）求導，并考慮式（6b），得等式Md（t）-Md（t）=MΦda-1dΦTdf（t）≡（t）式中ad=diag（a1，…，ad，1，…，d）。利用上式，可將式（7）寫成狀態方程的形式Ad+Bvd=pd（t）（8）式中pd（t）=（t）

fd（t）， vd=d

ud。

將vd（t）=Ψdqd（t）代入式（8），同時考慮特征方程，可得AΨdd（t）-AΨdΛdqd（t）=pd（t）（9）將式（4）對應的前d階復振型運動方程表達為矩陣形式d（t）-Λdqd（t）=a-1dΨTdp（t）上式兩邊左乘AΨd，與式（9）對比可得pd（t）=AΨda-1dΨTdp（t）（10）上式即為狀態空間下前d階特征向量對應的外荷載，而且它與后n-d階特征向量是正交的，從而不反映高階振型的貢獻。可將其進一步展開為

fd（t）= [M（ΦdΛda-1dΦTd）+

C（Φda-1dΦTd）f（t）]（11）

可見，位移ud（t）僅反映了荷載fd（t）的影響，fd（t）與f（t）的差別是式（6）截斷誤差的根本來源。荷載的截斷誤差可表達為

ed=f（t）-fd（t）=[I- M（ΦdΛda-1dΦTd）-

C（Φda-1dΦTd）]f（t）（12）

式中I為n階單位陣。若取全部振型參與計算，則荷載截斷誤差ed=0，證明如下：

對式（5a）求導，并與式（5b）相減，同時考慮復振型運動方程（4）可得Φa-1ΦT=0（13）同時，利用關系ΨBΨT=-ΨAΨTΛ=-aΛ及式（13）可得M-1=ΦΛa-1ΦT（14）從而，當d=n時，把式（13），（14）代入式（12）可得ed=0。

3高階振型的準靜態響應

通常，荷載f（t）的幅值和空間分布是隨時間變化的。為了討論方便，假設荷載的空間分布不隨時間變化，僅幅值隨時間變化。這樣，荷載可表示為其分布向量r0和幅值函數f（t）乘積的形式，即f（t）=r0f（t）通過傅里葉變換可將荷載f（t）分解為一系列簡諧荷載分量。因此，研究結構體系在單個諧振荷載作用下的響應，并不失一般性。不仿假設時間函數f（t）=sin（ωet），ωe為輸入頻率，則體系的n組振型運動方程為i（t）-λiqi（t）=pisin（ωet）

q·-i（t）-ii（t）=isin（ωet） 式中pi=Tir0ai，i=Tir0i。

上式中兩方程共軛，可取其一進行分析。假設結構在初始時刻處于靜止狀態，則方程的解qi（t）可表達為

qi（t）=-piω2e+λ2i[λisin（ωet）+ωecos（ωet）]+

piωeω2e+λ2ieλit（15）

等式第2項是瞬態解，隨著時間的延續會逐漸衰減掉，故最終只剩下穩態解

qi（t）=-piω2e+λ2i[λisin（ωet）+ωecos（ωet）]（16）

記pi=αi+iβi，αi和βi分別為pi的實部和虛部。在pi作用下的靜態響應為

q0i=-piλi=αiξi-βi1-ξ2iωi+iαi1-ξ2i+βiξiωi

從而，反應比可表達為

Di（t）=qi（t）q0i=DRisin（ωet-θRi）+

iDIisin（ωet-θIi）（17）

式中DRi=[1-（1-2ξ2iγ2i）]2+ξ2iγ2i（1+γ2i）2〖〗（1-γ2i）2+（2ξiγi）2；θRi=tan-1ξiγi（1+γ2i）1-（1-2ξ2i）γ2i，θRi∈（0，π）；DIi=γi1-ξ2i（1-γ2i）2+（2ξiγi）2；θIi=tan-1-1-γ2i2ξiγi，θIi∈π2，3π2，γi=ωeωi為輸入頻率與固有頻率比。

動力放大系數DRi，DIi與相位角θRi，θIi是關于阻尼比ξi、 頻率比γi的函數， 其函數曲線如圖1， 2所

圖1動力放大系數

Fig.1Dynamic magnification factors圖2相位角

Fig.2Phase angles

示。可見，動力放大系數實部、虛部與單自由度體系的位移、速度放大系數非常相近，當頻率比遠小于1時，放大系數實部接近于1，而其虛部接近于0，基本與阻尼無關；當頻率比遠大于1時，兩者都趨于0，也基本不受阻尼影響；當頻率比在1附近時，兩者都達到最大值，并且對阻尼非常敏感。相位角表示響應滯后于荷載的時間，從圖中可以看出，當頻率比接近0時，相位角實部接近于0，與荷載同相位；隨著頻率比的增大，實部相位大小先是隨著阻尼增大而增大，然后是減小；另外，需要指出的是，對于阻尼比為0的情況，實部相位始終為0。除了90°的相位差之外，相位角虛部曲線形狀與單自由度位移相位角一致，在頻率比為1時，所有阻尼比均通過180°，與荷載異相位。

由于隨著頻率比的減小，DRi→1，θRi→0，DIi→0及θIi→90°，因此，對于高階振型的動態反應可用其準靜態響應來代替，即qi（t）=q0isin（ωet）=-piλisin（ωet）（18）4截斷誤差修正〖*2〗4.1靜力修正不妨將結構的位移響應分為兩部分，第1部分為低階振型的貢獻，第2部分為高階振型的貢獻

標準化。可以看出，對于柔性結構，靜力修正效果不明顯，其原因是結構響應由處在地震卓越頻率范圍之內的低階振型控制，以動力響應為主；另外，還可以看出，隨著阻尼的增大，高階振型的貢獻有所增加（β=5時），靜力修正效果尤其是對于第一階振型變差，這在剛性結構中更為顯著，其原因是阻尼增大會導致一部分高階振型的頻率降低，如圖3（a）所示，從而使得利用靜力變形表示高階振型的貢獻變得不合理。相比于柔性結構，剛性結構的高階振型貢獻較大，尤其是對于基底剪力，剛性結構的靜力修正效果相當明顯，但當阻尼較大時，修正效果有所下降，如β=5時的基底剪力的修正效果。因此，振型加速度方法更適用于剛性結構，但當阻尼較大時需要考慮稍多的振型參與分析。〖BT1〗〖STHZ〗

7 結 論本文經過理論推導和實例分析，得出以下結論：（1）復振型疊加分析中振型截斷的誤差來源于動力荷載的截斷，可以通過荷載截斷誤差估計衡量結構響應的計算精度。（2）通過理論推導可得出與實振型表達形式一致的復振型靜力修正方法和振型加速度方法，并論證了兩種方法是相通的。（3）〖JP3〗復振型截斷的修正方法能夠將實振型的修正方法統一起來，可將其作為復振型修正方法的特例。〖JP〗（4）〖JP2〗算例分析表明，當結構主要的固有頻率大于輸入卓越頻率時，復振型修正方法更為有效，如算例中的剛性結構；另外，增加附加阻尼會減小部分高階振型對應的自振頻率，一定程度上降低修正效果。〖JP〗〖HJ*4/9〗〖BT3〗參考文獻：〖WT5”BZ〗〖HT5”SS〗[1] 〖ZK（#〗Clough R W， Joseph Penzien. Dynamic of Structures[M].3rd ed. Berkeley California 94704， USA： Computers and Structures， Inc.， 1995.[2] Williams D. Dynamic loads in aeroplanes under given impulsive loads with particular reference to landing and gust loads on a large flying boat[R]. Great Britain RAE Report. SME3309 and 3316， 1945.[3] Maddox N R. On the number of modes necessary for accurate response and resulting forces in dynamic analyses[J]. Journal of Applied Mechanics， ASME， 1975，42：516—517.[4] Hansteen O E， Bell K. On the accuracy of mode superposition analysis in structural dynamics[J]. Earthquake Engineering and Structural Dynamics， 1979， 7（5）：405—411.[5] 〖JP2〗Cornwell R E， Craig R R， Johnston C P. On the application of the mode acceleration method to structural dynamics problems[J]. Earthquake Engineering and Structural Dynamics， 1983， 11（6）：679—688.〖JP〗[6] Soriano H L， Filho F V. On the modal acceleration method in structural dynamics mode truncation and static correction[J]. Computer and Structure， 1988， 29（5）： 777—782.[7] Kulkarni S H， Ng S F. Inclusion of higher modes in the analysis of nonclassically damped systems[J]. Earthquake Engineering and Structural Dynamics， 1992， 21：543—549.[8] Borino G， Muscolino G. Modesuperposition methods in dynamic analysis of classically and nonclassically damped linear systems[J]. Earthquake Engineering and Structural Dynamics， 1986，14：705—717.[9] D′Aveni A， Muscolino G. Improved dynamic correction method in seismic analysis of both classically and nonclassically damped structures[J]. Earthquake Engineering and Structural Dynamics， 2001， 30：501—517.[10]〖KG*2〗TraillNash R W. Modal methods in the dynamics of systems with nonclassical damping[J]. Earthquake Engineering and Structural Dynamics， 1981， 9：153—169.[11]〖KG*2〗周錫元.一般有阻尼線性體系地震反應的振型分解方法[C].中國地震工程研究進展，北京：地震出版社，1992.〖KG2*2〗Zhou Xiyuan. The Modal Superposition Method of Earthquake Response of General Damped Linear System[C]. Published in Research Progress of Earthquake Engineering in China， Beijing： Seismological Press， 1992.[12]〖KG*2〗周錫元，俞瑞芳.非比例阻尼線性體系基于規范反應譜的CCQC法[J].工程力學，2006，23 （2）：10—17.〖KG2*2〗Zhou Xiyuan， Yu Ruifang. CCQC method for seismic response of nonclassically damped linear system based on code response spectra[J]. Engineering Mechanics， 2006， 23 （2）：10—17.[13]〖KG*2〗Caughey T K， O′Kelley M E J. Classical normal modes in damped linear dynamic systems[J]. Journal of Applied Mechanics， ASME， 1965，32：583—588.[14]〖KG*2〗Foss F K. Coordinates which uncouple the linear dynamic systems[J]. Journal of Applied Mechanics， ASME， 1958， 24：361—364.[15]〖KG*2〗Chopra A K. Dynamics of Structures： Theory and Applications to Earthquake Engineering[M]. 4th ed. USA： Prentice Hall， 2012.〖ZK）〗〖FL）〗〖HT〗〖HJ*3〗〖KH*8/9D〗〖JZ（〗〖WT4HZ〗Truncation error and improvement of complex modal superposition approach〖WT〗CHEN Huating1 ， TANG Ping1 ， PENG Lingyun2 ， LI Zhishan1 ， ZHOU Fulin1〖WT5”BZ〗（1.State Key Laboratory for Seismic Reduction Control & Structural Safety （Cultivation）， Guangzhou University， Guangzhou 510405， China； 2.College of Architecture and Civil Engineering， Beijing University of Technology， Beijing 100124， China）〖WT〗〖JZ）〗〖KH*8/9D〗〖WT5”HZ〗Abstract： Based on the complex modal superposition theory， the original error of complex mode truncation， which is used to analyze the nonclassically damped linear system， is studied and the error formula is obtained. Similar to higher modal responses of classical damping systems， the form of response ratio shows that with the decrease of the ratio of input frequency and natural frequency the displacements of higher complex modal equations can also represented by static displacement of dynamic load. Therefore， because of this fact， static correction method and modal acceleration approach about nonclassical damping systems are derived， whose forms are the same as that of the conventional methods. In addition， the correction methods of both real and complex modal superposition approaches are unified according to the condition of complex modes degenerating to real modes （real part is proportional to imaginary part）. At last， the effectiveness of correction methods is verified through a numerical example， particularly for the case that the major natural frequencies are beyond the input significant frequency range； in addition， it is observed that rising the added damping will reduce， to some degree， the efficiency of the modified method since the added damping decreases some highorder natural frequencies.〖KH2D〗Key words： 〖ZK（〗linear vibration； dynamic analysis； complex modal superposition approach； static correction； modal acceleration method〖ZK）〗〖KH2D〗〖HTH〗作者簡介： 〖ZK（〗〖HTSS〗陳華霆（1988—），男，博士研究生。電話： （020）86395053； Email： sdcht2008@yeah.net〖ZK）〗〖HTH〗通訊作者： 〖ZK（〗〖HTSS〗譚〖KG1〗平（1973—），男，教授，博士生導師。電話： （020）86395007； Email： ptan@foxmail.com〖ZK）〗〖HT〗〖HJ〗〖WT〗〖ST〗〖LM〗