大跨斜拉桥Rayleigh阻尼系数约束优化解

潘旦光，靳国豪，高莉莉

(北京科技大学 土木系，北京 100083)

我国规范[1]规定，进行斜拉桥、悬索桥、单跨跨径150 m以上等特殊桥梁地震反应分析时可采用时程分析法、多振型反应谱法及功率谱法。用直接积分法进行时程分析时必涉及阻尼矩阵的建立。由于形成阻尼机理复杂，影响因素多[2]，因此无法直接利用构件尺寸、材料性质直接构建阻尼矩阵，而采用构造的方法。由于工程结构阻尼比较易测量，且能反应结构宏观阻尼特性，因此常用阻尼比构建阻尼矩阵。阻尼可分为复阻尼及粘滞阻尼。复阻尼耗能与迫振频率无关，这与结构试验结果相符。阻尼为复数，易于用频域进行分析。对非线性分析只能采用等效线性化方法，无法进行真非线性分析。与复阻尼相比，用时域内粘滞阻尼更方便[3-4]。

目前已有诸多粘滞阻尼矩阵的构造方法[5-7]，其中最具代表性的为Rayleigh阻尼、Caughey阻尼及叠加振型阻尼。董军等[8-9]将Rayleigh阻尼与叠加振型阻尼结合，构造出任意阶模态阻尼比等于精确值的阻尼矩阵。但常用的为质量矩阵与刚度矩阵线性组合的Rayleigh阻尼矩阵,可表达为

C=αM+βK

(1)

式中：α，β分别为质量、刚度比例阻尼系数；M，K分别为质量、刚度矩阵。

计算α，β一般先进行结构模态分析，后指定两阶参考模态(ωi，ωj)的阻尼比等于已知值(ζi*，ζj*)。α，β计算公式为

(2)

由式(1)得各阶模态阻尼比为

ζn=α/2ωn+βωn/2

(3)

进行桥梁设计时,虽事先无法获知模态阻尼比的精确值，但可据已建桥梁阻尼比经验值设各阶模态精确值为常量。如钢筋混凝土桥梁阻尼比为0.05。因此式(3)所得阻尼比仅在两阶参考频率处阻尼比等于精确值，而其它频率阻尼比存在一定误差。由于阻尼比会影响结构动力反应，因此合理指定两阶参考频率成为时程分析关键。结构基频为影响结构动力反应最主要模态。规范[1]中建议第一个参考频率选择基频。而对第二个参考频率选择无明确方法。楼梦麟等[10-11]通过对大跨拱桥计算认为参考频率的选择会显著影响结构反应，应选振型参与系数较大模态作为参考频率；通过对深覆土层模型分析认为计算Rayleigh阻尼系数需考虑输入地震波频谱特性影响。为避免人为选择两阶参考频率所致任意性，刘红石[12]提出的最小二乘法计算α，β可使阻尼比误差在平方意义上最小，但并未考虑各阶模态对动力反应贡献的差异，无法保证对结构反应有显著贡献模态的阻尼比取值合理。杨大彬等[13-14]针对网壳结构地震反应分析，建立基于多参考振型的加权最小二乘法用于计算Rayleigh阻尼系数。潘旦光[15-16]基于振型叠加反应谱理论，详细讨论式(2)中两最优参考频率及荷载空间分布、频谱特性及结构动力特性关系，并建立Rayleigh阻尼系数的无约束优化求解方法。但对桥梁结构而言，与结构所受激励输入方向对应的第1阶整体振动模态为控制结构动力反应的关键模态之一，通常以该阶模态为计算Rayleigh阻尼系数的第1个参考频率。因此，Rayleigh阻尼系数计算实际为合理选取第2个参考频率问题。大跨桥梁参与振动模态数较多，频率密集，不同于普通建筑结构。为此，本文以大跨桥梁为分析对象，以结构峰值位移误差最小为目标，建立求解Rayleigh阻尼系数的优化方法。并将基频阻尼比等于精确值作为约束条件，形成等价于指定最优第2参考频率的约束优化求解方法。并以840 m斜拉桥为例，说明Rayleigh阻尼系数约束优化解的变化规律及计算方法精度。

1 Rayleigh阻尼系数约束优化解

地震输入作用下多自由度体系强迫振动方程[17]为

(4)

通过模态分析获得结构前N阶频率ωn及模态φn(n=1,2,…,N)。由模态叠加反应谱法理论知，基于Rayleigh阻尼与精确阻尼比,分别计算所得第n阶模态对第k自由度最大位移反应ukn为

(5)

Sd(ζn,ωn)=

基于平方和开平方原理，第k自由度位移反应误差为

(7)

为简化计算，将反应谱函数用1阶Taylor级数展开，得

(8)

将式(8)代入式(7)，整理得

(9)

(10)

采用Lagrange乘子法求式(10)约束条件下式(9)最小值，令

(11)

将式(11)分别对X，λ求导，并令相应导数为零得代数方程组为

(12)

2 参数研究

为说明公式的计算精度及影响因素，以图1悬臂梁为例进行讨论。梁长6 m，弹性模量20 GPa，线密度600 kg/m，截面转动惯量1.6×10-3m4。将梁分为10个单元，用集中质量模型，体系共10阶模态，其中前4阶模态分析结果见表1。

图1 悬臂梁计算模型

表1 悬臂梁模态分析结果

设体系基础运动竖向加速度为

(13)

激振频率考虑三种情况：ω=0.7ω1, 6.0ω1，16.0ω1，即分别考虑激振频率小于基频、接近ω2、接近ω3。在任意激振频率下第n阶模态反应谱及导数为

Sd(ζ,ωn)=

(14)

(15)

2.1 Rayleigh阻尼模型计算精度

阻尼为结构固有特性，在线弹性范围内阻尼比与外部荷载无关。而Rayleigh阻尼为近似阻尼，因此为使构建的阻尼矩阵动力反应计算误差小，应使对结构动力反应有显著贡献模态的阻尼比误差越小越好。不同荷载频谱特性不同，Rayleigh阻尼合理参考频率与荷载有关。以含所有模态在精确阻尼比下的模态叠加法计算结果为精确解，记r*。采用Rayleigh阻尼模型计算所得近似解记为r，Rayleigh阻尼模型计算结果相对误差为

(16)

设悬臂梁所有频率精确阻尼比为1%，表2、表3为简谐荷载作用下用Rayleigh阻尼矩阵进行计算所得顶点A的竖向位移及支座B处弯矩相对误差。计算中本文方法的φkn=φn(A)用前4阶模态参与优化计算。为与传统计算方法进行比较，用式(2)计算Rayleigh阻尼系数时指定两阶参考频率考虑三种组合，分别为i=1&j=2,i=1&j=3 ，i=2 &j=3。由表2、表3看出：① 由本文方法所得Rayleigh阻尼矩阵的反应误差均小于或等于传统计算方法，计算精度良好;② 在不同迫振频率下结构动力反应的显著贡献模态不同，当ω/ω1= 0.7，6.0时，结构动力反应由前两阶模态控制。以i=1 &j=2构建阻尼矩阵计算误差最小，而ω/ω1=16.0时，以i=1 &j=3构建阻尼矩阵的计算误差最小。由此可知，计算Rayleigh阻尼系数最优参考频率应据荷载频谱特性作相应调整。

表2 A点位移相对误差(%)

表3 B点弯矩相对误差(%)

2.2 优化参考频率影响因素

为说明本文方法所得Rayleigh阻尼系数特点，不同迫振频率下所得阻尼系数及相应阻尼比见表4。由表4看出，ω/ω1=0.7，6.0时，第2个优化参考频率j=2；而当ω/ω1=16.0时，第2个优化参考频率j=3。结合表2、表3知，本文方法所得Rayleigh阻尼系数为综合考虑结构动力特性及输入荷载频谱的结果，使对结构有显著贡献模态的阻尼取值更合理，可直接获得Rayleigh阻尼系数，无需人为选择参考频率，便于应用。

表4 优化Rayleigh阻尼系数及相应阻尼比

3 工程概况及输入地震波

对简单结构而言，最优两参考频率较易判断；但对大跨桥梁而言，结构显著贡献模态数量较多，人为方法选择合理参考模态较困难。为验证本文方法计算精度，以840 m斜拉桥为例，分析地震作用下优化Rayleigh阻尼系数的合理性及计算精度。大桥为60 m+120 m+480 m+120 m+60 m双塔双索半漂浮体系钢混混合斜拉桥。塔高155.1 m，桥面板以上高103 m。主梁为宽28.5 m、高3 m的单室封闭钢箱梁。塔柱与混凝土主梁用C50混凝土，钢箱梁用Q345qc。塔座、边墩及辅助墩用C40混凝土，承台用C30混凝土。拉索用112根Φj15.24 mm钢绞线。用单脊梁式方法建立有限元模型。塔墩及主梁采用三维梁单元，索用杆单元，建立有限元模型见图2。

图2 斜拉桥有限元模型

有限元模型整体坐标系x方向为顺桥向，y方向为竖向，z方向为横桥向，模型共1 149个节点、1 196个单元。设桥梁各阶模态阻尼比精确值为0.02，以x方向地震输入为例，分析塔顶(A点)、主梁跨中(B点)x方向水平位移及桥墩(C点)处剪力及弯矩在不同阻尼模型下结构动力反应特点与计算精度。其它方向地震反应分析方法类似。

式(9)计算中权重系数wkn涉及第k自由度模态位移φkn，即优化的参考自由度问题。理论上参考自由度可选择任意自由度；但对结构而言，若结构最大位移误差小，则其余自由度动力反应的误差通常也较小。而斜拉桥桥塔最大位移位于桥塔顶端，桥身最大位移位于桥梁跨中，故计算中参考自由度分别取A、B点，即x方向地震输入时分别以φkn=φn(Ax)，φkn=φn(Bx)为权重系数的参考自由度，并比较其对计算结果影响。

基于图2有限元模型，选4条不同类型地震波作为输入，将地震波峰值加速度统一调整为0.1 g，分析不同Rayleigh阻尼系数计算方法引起结构动力反应误差。地震波加速度时程与7个阻尼比(ζ=0.005, 0.01, 0.02, 0.03, 0.05, 0.10, 0.20)的反应谱见图3。实际地震输入的反应谱为不规则曲线，无法建立位移反应谱显式表达式。统计分析结果[18]表明，位移反应谱与阻尼比对数间具有线性关系。即对某确定地震波，自振频率ωn体系的反应谱随阻尼比ζ变化规律可表示为

Sd(ζ,ωn)=an(ωn)+bn(ωn)ln100ζ

(17)

式中：an(ωn)，bn(ωn)为拟合参数，可用最小二乘法计算。

(18)

4 斜拉桥计算结果

4.1 自振特性

对斜拉桥而言，拉索应力对结构动力特性及动力反应有显著影响。为考虑重力对拉索应力影响及拉索几何刚度，斜拉桥模态分析分两步：① 计算重力及初始应力下结构应力；② 进行有预应力的模态分析。斜拉桥主要振动模态见表5。x方向累积振型参与质量见图4。由各模态振型参与质量、模态特征看出，第1阶模态为顺桥向振动第1个显著贡献模态。式(12)计算时顺桥向约束模态取第1阶。

图3 地震波加速度时程及位移反应谱

表5 斜拉桥模态分析结果

图4 斜拉桥累积振型参与质量

4.2 Rayleigh阻尼系数变化规律

由式(9)知，Rayleigh阻尼系数优化计算结果与参与优化计算模态数N有关。顺桥向优化Rayleigh阻尼系数随模态数N的变化见图5。结合图4计算结果知，在累积振型参与质量小于70%时(前128阶模态)，随模态数的增加使累积振型参与质量发生显著变化的模态处Rayleigh阻尼系数亦会发生显著变化。累积振型参与质量超过70%后，部分高阶振型参与质量较大(第538、539阶)，但Rayleigh阻尼系数未在此处发生相应变化。由于振型参与质量系数是用于度量各模态对基底剪力贡献大小的指标，部分高阶局部振型对结构基底剪力有显著影响，但塔顶、主梁位移反应主要由低阶模态控制，而Rayleigh阻尼系数优化计算方程本质为基于位移反应的优化方程，因此当累积振型参与质量超过70%后，Rayleigh阻尼系数优化结果趋于稳定。与高层建筑等结构Rayleigh阻尼系数变化规律有较大不同。

4.3 优化Rayleigh阻尼系数及优化参考频率

不同地震波激励下Rayleigh阻尼系数优化解稳定计算结果见表6。将所得Rayleigh阻尼系数代入式(3)，即得各阶模态阻尼比，除作为约束条件的模态外，第2个等于精确阻尼比对应的固有频率阶数称为优化参考频率j′。由表6看出，① 本文所得优化参考频率不一定位于结构振型参与质量较大固有频率上。因除基频外，桥梁结构无具有显著统治地位模态，最优参考频率为综合考虑多阶模态动力反应结果；② 输入地震波不同最优参考频率也不同，说明本文结构最优参考频率考虑输入地震波频谱特性影响；③ 以主梁跨中位移为优化目标的优化参考频率值小于以塔顶为优化目标的计算结果。此因为主梁位移由少数低阶模态控制(前28阶模态)，而塔顶位移频率成分相对丰富(前52阶模态)。因此，优化方程自动识别出的优化目标位移显著贡献模态有所不同。

图5 顺桥向地震输入Rayleigh阻尼系数优化解

表6 优化Rayleigh阻尼系数及优化参考频率

4.4 Rayleigh阻尼模型反应峰值计算误差

由图4知，当模态个数达539阶时，体系x方向累积振型参与质量超过90%；620阶模态后累积振型参与质量超过99%。为此本文以620阶模态振型分解法计算所得结构动力反应量作为精确解。对传统Rayleigh阻尼计算方法，选i=1 &j=56，i=1 &j=108及i=1 &j=539三种组合方式进行计算。本文方法考虑φkn=φn(Ax)，φkn=φn(Bx)两种情况。塔顶A点、主梁跨中B点在不同地震波输入下水平位移峰值计算误差见表7。El Centro波地震输入作用下塔顶A点、主梁跨中B点水平位移时程及反应谱见图6、图7。

由计算结果看出，① 由误差的平均值看，传统Rayleigh阻尼计算方法与本文方法位移反应误差均小于5%，均满足工程需要。由位移FFT谱知，A、B点水平位移均由第1阶模态反应控制。第1阶模态阻尼比等于精确解时，结构位移反应误差即可控制在合理范围内。② 第2个参考频率越大(如i=1 &j=539)，两参考频率之间模态阻尼比越小于精确解，所得计算结果越大于精确解。本文所用优化参考频率使对结构动力反应有影响部分模态阻尼小于精确解，部分大于精确解，由此达到位移反应误差为最小目的。用本文方法计算所得结果误差均小于传统方法，表明计算精度良好。③ 由精确解的Fourier谱知，A点水平位移主要频率成分小于2 Hz，B点水平位移主要频率成分小于1.2 Hz。约束优化算法中考虑各模态对总位移贡献影响，即表6中φkn=φn(Ax)的第2个参考频率大于φkn=φn(Bx)的原因。

表7 A、B点水平位移反应峰值相对误差(%)

图6 A点水平位移及Fourier谱

不同地震波作用下桥墩C点剪力Fx及弯矩Mz计算误差(Fx，Mz为整体坐标系下计算结果) 见表8。

El Centro波顺桥向地震输入作用下C点剪力Fx与弯矩Mz时程及反应谱见图8、图9。由计算结果看出，① 用传统方法计算的三种组合计算误差平均值均大于5%。且随第2个参考频率增大误差增大。此因对该桥墩剪力、弯矩贡献最大的两阶模态为第28阶(频率1.123 Hz)、第56阶(频率3.058 Hz)。本文所选三种传统Rayleigh阻尼计算组合均使第28阶模态阻尼比小于精确解。第2个参考频率越大该阶模态阻尼比越小，导致时程分析结果较精确解大。② 本文方法虽以位移反应误差为目标函数获得Rayleigh阻尼系数，但仍可兼顾内力反应的频谱特性，使剪力、弯矩计算误差均小于传统方法，计算精度良好。③ 比较φkn=φn(Ax)与φkn=φn(Bx)两组统计结果，以φkn=φn(Ax)为权重函数的计算结果大于精确解，以φkn=φn(Bx)时的结果小于精确解。此因φkn=φn(Bx)时所得第2个参考频率位于第28阶模态附近，而Rayleigh阻尼将使高阶振动由于人为原因造成的阻尼增大而被消除。导致对结构剪力、弯矩有显著影响的第56阶模态反应小于精确解，从而所得内力计算结果小于精确解。而φkn=φn(Ax)优化所得第2阶参考频率位于第28、56阶模态之间，所得Rayleigh阻尼将高估第28阶模态振动同时又低估第56阶模态振动，故误差相对较小且计算结果大于精确解。当结构动力反应小于精确解时，所得计算结果较不安全。因此以塔顶位移为参考自由度进行顺桥向地震激励下Rayleigh阻尼系数的计算更合理。

表8 C点剪力、弯矩峰值相对误差(%)

图8 C点剪力及Fourier谱

5 结 论

本文针对桥梁结构基频模态为结构动力反应显著贡献模态特点，提出Rayleigh阻尼系数的约束优化分析方法。利用悬臂梁简谐振动动力反应，讨论荷载迫振频率对构建合理Rayleigh阻尼影响，结论如下：

(1) Rayleigh阻尼系数约束优化解可综合考虑荷载频谱特性及结构动力特性影响，直接获得Rayleigh阻尼系数，使对结构有显著贡献模态的阻尼比合理。

(2) 大跨斜拉桥塔顶及主梁跨中位移主要由少数低阶模态控制。用传统Rayleigh阻尼计算方法分析时，所得结果计算误差可控制在5%以内；桥墩的剪力、弯矩因其频率成分丰富，部分高阶振型模态对结构的动力反应同样贡献显著；采用人为方法选择两阶参考模态时，选取不当将引起较大计算误差。

(3) 随计算模态个数的增加Rayleigh阻尼系数约束优化解所得α，β趋于稳定。大跨桥梁的计算所用模态个数对应的累积振型参与质量超过70%时，可得α，β的稳定解。

(4) 本文方法所得Rayleigh阻尼系数计算误差在统计意义上小于传统方法。可避免因人为选择参考频率带来的任意性，适合工程结构的计算与分析。考虑顺桥向振动时，以桥身位移为优化目标所得Rayleigh阻尼系数易使结构动力反应小于精确解，安全性较低；以塔顶位移为参考自由度时所得结构动力反应误差较小且大于精确解，因此塔顶位移作为参考自由度更合理。

[1] JTG/T B02-01-2008,公路桥梁抗震设计细则[S].

[2] Jeary A P. Damping in structures [J]. Journal of Wind Engineering and Industrial Aerodynamics,1997,72(1):345-355.

[3] 楼梦麟,潘旦光. 滞后阻尼在土层时域分析中的应用[J]. 同济大学学报, 2004, 32(3): 281-285.

LOU Meng-lin, PAN Dan-guang. Hysteretic damping application in time domain analysis of soil layer[J]. Journal of Tongji University, 2004, 32(3): 281-285.

[4] 何钟怡，廖振鹏，王小华. 关于复阻尼理论的几点注记[J]. 地震工程与工程振动,2002, 22(1):1-6.

HE Zhong-yi, LIAO Zhen-peng, WANG Xiao-hua. Some notes on theory of complex damping[J]. Earthquake Engineering and Engineering Vibration, 2002, 22(1): 1-6.

[5] Luco J E. Anote on classical damping matrices[J]. Earthquake Engineering and Structural Dynamics, 2008, 37(4): 615-626.

[6] Adhikari S. Dampingmodelling using generalized proportional damping [J]. Journal of Sound and Vibration, 2006, 293(1/2): 156-170.

[7] 黄宗明,白绍良,赖明. 结构地震反应时程分析中的阻尼问题评述[J]. 地震工程与工程振动,1996, 16(2): 95-105.

HUANG Zong-ming, BAI Shao-liang, LAI Ming. Review on the damping in earthquake response time-history analysis of structures[J]. Earthquake Engineering and Engineering Vibration, 1996, 16(2): 95-105.

[8] 董军，邓洪洲，王肇民. 结构动力分析阻尼模型研究[J]. 世界地震工程,2000, 16(4): 63-69.

DONG Jun, DENG Hong-zhou, WANG Zhao-min. Studies on the damping models for structural dynamic time history analysis[J]. World Information on Earthquake Engineering, 2000, 16(4): 63-69.

[9] 克拉夫R，彭津J,著.王光远，译.结构动力学(第二版)[M]. 北京：高等教育出版社,2006.

[10] 楼梦麟,张静. 大跨度拱桥地震反应分析中阻尼模型的讨论[J]. 振动与冲击, 2009, 28(5): 22-26.

LOU Meng-lin, ZHANG Jing. Discussion on damping models for seismic response analysis of long-span bridge [J]. Journal of Vibration and Shock, 2009, 28(5): 22-26.

[11] 楼梦麟，邵新刚. 深覆盖土层Rayleigh阻尼矩阵建模问题的讨论[J]. 岩土工程学报, 2013, 35(7):1272-1279.

LOU Meng-lin, SHAO Xin-gang. Discussion on modeling issues of rayleigh damping matrix in soil layer with deep Ddeposit [J]. Chinese Journal of Geotechnical Engineering, 2013, 35(7): 1272-1279.

[12] 刘红石. Rayleigh阻尼比例系数的确定[J]. 噪声与振动控制, 1999(6): 21-22.

LIU Hong-shi. Determination of rayleigh damping scale coefficient[J]. Noise and Vibration Control,1999(6):21-22.

[13] Yang D B, Zhang G Y G, Wu J Z. Computation of rayleigh damping coefficients in seismic time-history analysis of spatial structures[J]. Journal of the International Association for Shell and Spatial Structures, 2010, 51(2): 125-135.

[14] 杨大彬, 张毅刚, 吴金志. 基于多参考振型的Rayleigh阻尼系数计算方法在单层柱面网壳中的应用[J]. 空间结构, 2011,17(3): 8-15.

YANG Da-bin, ZHANG Yi-gang, WU Jin-zhi. Application of the multiple reference modes based computation method of Rayleigh damping coefficients in single-layer cylindrical latticed shell [J]. Spatial Structures, 2011, 17(3): 8-15.

[15] 潘旦光. 直接确定Rayleigh阻尼系数的一种优化方法[J]. 工程力学,2013, 30(9):16-21.

PAN Dan-guang. An optimization method for the direct determination of Rayleigh damping coefficients[J]. Engineering Mechanics, 2013, 30(9): 16-21.

[16] 潘旦光. 地震反应分析中Rayleigh阻尼系数的优化解[J]. 工程力学,2013, 30(11):15-20.

PAN Dan-guang. An optimization solution for Rayleigh damping coefficients in seismic response analysis[J]. Engineering Mechanics, 2013, 30(11):15-20.

[17] Chopra A K. Dynamics ofstructures: theory and applications to earthquake engineering[M]. New Jersey: Englewood Cliffs, Prentice-Hall, 1995.

[18] Newmark N M, Hall W J. Earthquake spectra and design,earthquake engineering research institute[R]. California: Berkeley, 1982:29-37.