复振型叠加法截断误差及改进

陈华霆 谭平 彭凌云 李志山 周福霖

摘要： 根据复振型叠加理论，对具有非比例阻尼特性的线性体系，分析了复振型截断误差的来源，给出了误差公式。与比例阻尼体系的高阶实振型响应类似，通过反应比表达式可以看出随着输入频率与自振频率比的减小，高阶复振型的响应也趋于荷载的静态响应。基于此推导出了与传统修正方法形式一致的复振型静力修正方法和振型加速度方法，同时根据复振型退化为实振型的条件（实部与虚部成比例关系）建立了实、复振型修正方法的统一。用数值算例对修正方法进行了验证，结果表明：复振型修正方法能够很好地提高振型叠加的计算精度，特别是当结构主要固有频率大于输入卓越频率时修正效果更为显著；还发现附加阻尼较大时，导致部分高阶振型自振频率减小，在一定程度上会降低修正效果。关键词： 线性振动； 动力分析； 复振型叠加； 静力修正； 振型加速度方法

中图分类号： 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〗