山区大跨窄悬索桥抖振响应时域有限元分析

周 戈，盛光祖，李浩然

（1.武汉地产开发投资集团有限公司，湖北 武汉 430022；2.防灾减灾湖北省重点实验室（三峡大学），湖北 宜昌 443002）

0 引 言

山区大跨度窄悬索桥是一种柔性的高次超静定结构，具有跨越能力大、内力分布均匀、外形轻盈美观等优点，逐渐受到桥梁界的高度关注。但该类桥梁的加劲梁宽度一般设置为2～5 m，随着跨径的增大，整体结构高跨比和宽跨比急剧降低，使得窄悬索桥竖向与横向刚度削弱，导致桥梁整体结构的动力稳定性变差[1]。在风雨、车辆和环境激励下极容易产生较大振幅，降低了行人车辆的舒适度和结构安全性，因此开展窄悬索桥结构的风振响应研究具有重要意义。

针对山区大跨度窄桥面悬索桥的抗风问题，许多学者开展了广泛的研究。李春光等人采用增量分析方法，研究了某在建窄桥面大跨度悬索桥的静风失稳的全过程[2]。长安大学李加武等人提出了提高窄悬索桥抗风能力的措施[3]。研究表明，自然风分解为平均风速（平均风）和脉动风速ui（脉动风）两部分，相对应的风荷载也分为静风荷载和脉动风荷载。其中，平均风是指在固定的一段时间内，风速大小、方向为固定值，不随时间而改变，平均风的周期较长，其性质相当于静力作用，而脉动风是由大气的不规则运动造成，脉动风速大小、方向均随时间按随机规律变化，周期较短，其性质相当于动力作用。目前，对于脉动风随机场模拟主要有以下四种方法：谐波合成法（WAWS）[4]，本征正交分解法（POD）[5]，线性滤波法[6]（AR或ARMA）以及小波变换法。线性滤波法优点是计算量小、效率高、耗时少，是需要满足特定条件的一种模拟方法，在模拟过程中往往需要根据已有的经验判断，模拟精确度不能保证；小波变换法[7]的优势在于模拟非平稳信号，但是模拟结果的精确度不稳定，容易受所选小波基以及小波系数的影响，该方法不建议在模拟平稳随机过程时使用，容易造成不可忽视的误差;谐波合成法理论简单明了易掌握，将随机脉动风速看成一系列余弦波的叠加，并且可以生成无条件稳定的高精度模拟结果。

因此，本文以某山区加劲梁窄悬索桥为工程研究背景，采用谐波合成与FFT转换技术相结合的方法，模拟某山区加劲梁窄悬索桥的三维脉动风场，基于ANSYS大型有限元分析软件中APDL语言建立窄悬索桥风振响应模型，分析山区大跨度窄悬索桥结构的抖振响应。

1 脉动三维风场理论分析及数值建模

1.1 脉动三维风场理论分析

脉动风能量在频域范围中的分布情况通常采用脉动风功率谱函数来表达，脉动风功率谱函数反映了脉动风中不同频率成分对湍流脉动总动能的影响，根据脉动风风向进行分类，可以划分为顺风向脉动风速谱和竖向脉动风速谱。我国《公路桥梁抗风设计规范》（JTG_TF50—2011）中顺风向脉动风速谱采用随高度变化的Kaimal-Simiu谱，其表达式为：

式中：u2*为气流剪切速度，m/s；f为无纲量莫宁（或相似率）坐标,，（z）=uln（z/z），其中 z*0为高度，z0为地面粗糙长度，m；k为卡曼常数，取0.4。

而竖向脉动风速功率谱采用Panofsky谱，表达式为：

一般在湍流风场中，一定空间内两个不同坐标点的风速大小和方向往往会有所差别，也可能是完全不相关的情况[8]。在窄悬索桥结构上某一个模拟点的脉动风荷载达到了最大值时，与该点距离较远的另一个点的脉动风荷载一般不会在同一时间里也达到上限值，因此在一定空间范围之内，与该模拟点距离越远，其脉动风荷载就越不可能同时达到上限值。随机风场并非完全一样，这是因为两个点之间距离越远，脉动风速的相关性就随之越小，其中相干函数的值域范围在（0，1）之间，该值越大，说明两点之间的相关性越强。本文采用Davenport建议的空间相干函数模型来描述桥面水平向两点i与j之间的相关性，表达式为：

式中：Coh（Δij，ω） 为空间 i和 j两点的相干函数；cx为水平向衰减因子。

1.2 谐波合成-FFT

（1）谐波合成法

脉动风速可以看成一个一维n变量零均值的高斯随机过程 ｛（ft）｝ ，这个随机过程包含 f（1t）、f（2t）、f（3t）、...、fn（t）共n个变量。根据Shinozuka理论，脉动风速过程可以用下式表示：

式中：Vi（t）为脉动风速，m/s；N为频率采样个数，取 210；Δω 为圆频率增量，rad/s；ωjk为双索引频率；为随机相位角，rad。

（2）FFT加速技术

虽然谐波合成法能高精度的模拟出脉动风速过程，但是计算量通常较大，耗时长，用时较长主要和以下两点有关：（1）随机脉动风速的谐波项叠加过于复杂，该计算过程需要占用大量的时间；（2）每一个对应的频率项都需要进行互谱密度矩阵的Cholesky分解，而且随着模拟节点数量n增加，功率谱密度矩阵维数也越来越大，所以计算量会按照n2/2的速度增加，为了提高脉动风速过程的计算效率，本文引入FFT转换技术。首先，将式（4）转换为以下表达形式[9]：

式中：Re表示复数的实部；Δt为时间增量，s;p=1，2，...，M×n-1；M=2N；q 为 p/M 的余数，q=1，2，...，M-1；其余量的含义同式（4）。

其中Bji（mΔω）按照下式计算：

1.3 脉动风速值模拟计算结果

根据上述计算原理，本文采用MATLAB数学编程软件，实现脉动风速值模拟的计算机算法。数值模拟时主要参数取值见表1。图1和图2所示为10#模拟点脉动风速时程，图3和图4为10#模拟点脉动风功率谱与目标谱的对比。

表1 谐波合成法模拟脉动风速主要参数

图1 部分模拟点顺风向脉动风速时程

图2 部分模拟点竖风向脉动风速时程

基于谐波合成法（WAWS）与FFT结合的方法，模拟了窄悬索桥主梁上共19个点的顺风向和竖风向脉动风速时程，并展示了其中最具代表性的几个点。从图中可清晰的看到，模拟出来的空间多点脉动风速时程随机性较好。如图1和图2可知，模拟点的顺风向脉动风速在-15～15 m/s范围之间波动，竖风向脉动风速在-6～6 m/s范围之间波动，反映了脉动风随时间随机规律变化的特性。

图3和图4可知，模拟出的自功率谱函数和互功率谱函数与目标谱都有相同的变化趋势，模拟功率谱曲线在目标谱曲线周围有规律的波动。在低频区域（10-3～10-1Hz）范围内，功率谱曲线拟合较好，0.1 Hz以后，两种曲线差距开始发散，整体拟合度较好，说明模拟出来的功率谱精确度高。

图310#点模拟自功率谱与目标谱

图410#点模拟自功率谱与目标谱

2 工程概况及有限元模型

2.1 工程概况

本文以某山区窄悬索桥有工程研究背景，该桥位于湖北省兴山县古洞口水库。某山区加劲梁窄悬索桥整体布置图、加劲梁横截面如图5所示。

图5某山区加劲梁窄悬索桥全桥布置图（单位：cm）

桥梁总跨径185.4 m，主桥为150 m单跨双铰地锚式悬索桥，主桥主梁采用纵横梁格钢结构体系，桥面全宽4.0 m，人行道净宽3.5 m。横梁型钢使用双拼[32a槽钢，中横梁顺桥向间距为2.0 m。横梁两端开孔，用于锚固抗风拉索的上端锚头。双拼[32a横梁上沿顺桥向布置14根纵梁，其中12根[16槽钢以及位于两侧的[40槽钢，与横梁通过焊接连接。桥面板采用10 mm厚花纹钢板。引桥横断面为混凝土小箱梁。全桥共73根吊杆，吊杆间距2.0 m（两端为加强吊杆，间距3.0 m）。桥塔为H形钢筋混凝土结构，截面形状为矩形。

2.2 有限元建模

基于窄悬索桥结构约束特点与设计图纸说明，对悬索桥桥塔底部、加劲梁端部、主缆锚固端和索鞍处约束进行力学简化，某山区窄悬索桥有限元计算分析模型如图6所示。采用ANSYS大型有限元计算分析程序对某山区加劲梁窄悬索桥进行建模，总体坐标系以桥跨方向为X轴，以桥宽方向为Y轴，以桥塔方向为Z轴。根据某山区桥梁结构实际情况进行有限元模拟，全桥共4083个单元。

有限元模型的加劲梁和索塔结构采用空间梁单元BEAM4模拟，BEAM4是一种可以承受拉、压、弯、扭的受力单元，每个节点各有6个自由度。主缆、吊索、抗风系统采用空间杆单元LINK10模拟，LINK10为三维仅受拉或受压单元，能够有效模拟悬索桥索类结构的工作状态。栏杆、索夹等采用节点质量单元MASS21模拟，桥面板采用壳单元SHELL63模拟。

图6 加劲梁窄悬索桥有限元计算分析模型

该模型考虑了主缆外层防腐材料质量的影响，在主缆参数设置中给主缆密度一个大于1的修正系数，主缆、吊杆、抗风主缆和抗风拉索的质量密度按照规范给定容重，根据其有效截面面积进行计算，各类构件材料的弹性模量参考设计图纸取值。

3 全桥抖振时域化分析

3.1 风荷载计算

风荷载包括三个部分：静风力、抖振力和气动自激力，其中抖振力和气动自激力为动荷载，抖振力由脉动风引起，气动自激力由结构气动耦合产生[10]。

（1）抖振力计算

根据Scanlan教授的准定常气动理论，窄悬索桥加劲梁的抖振力可按下式计算：Lb（t）=

式中：ρ为空气密度，kg/m3；U为来流平均风速，m/s；B 为构件截面沿主流方向尺寸，m；C'L、C'D、C'M分别为升力、阻力和升力矩系数曲线斜率；u（t）、w（t）分别为顺风向，m/s和竖风向脉动风速，m/s。

（2）气动自激力计算

大跨加劲梁窄悬索桥主梁的气动自激力可以参考Scanlan教授提出的自激空气动力计算方法，作用在窄悬索桥主梁结构单位长度上受到的气动升力Lae、气动阻力Dae和气动扭矩Mae，可以分别表示为竖向位移h、水平位移p和扭转位移α的函数，采用无量纲气动导数H*i、P*i、A*i（i=1，2，...，6）来表达，计算公式如下：

式中：U、ρ、B 同（10）～（11）；K 为无量纲频率，K=Bω/U；ω 为振动圆频率，rad/s；气动导数 H*i、P*i、A*i（i=1，2，...，6）是无纲量风速 U*=U/fB或者无量纲频率的函数，他们的取值与窄悬索桥加劲梁截面的尺寸有关。

3.2 抖振响应时域化分析

根据脉动风速模拟结果计算出脉动风产生的抖振力，静风荷载根据设计基准风速Vd=29.03 m/s进行计算，气动自激力由主梁结构气动耦合引起，采用ANSYS中MATRIX27矩阵模拟加劲梁气动刚度与气动阻尼。考虑三种荷载的共同作用，基于ANSYS中APDL编程功能实现窄悬索桥结构抖振时域化分析，跨中主梁节点抖振位移响应时程如图7和图8所示。

图7 某山区加劲梁窄悬索桥跨中横向位移抖振时程响应

从图7和图8可以清晰看出，横向与竖向抖振时程响应均有明显上限限值，横向位移振动上限值为16.4 cm，竖向位移振动上限值为8.8 cm，从横向与竖向振动幅值来分析，可以看出横向振动位移较大，说明该桥横向抗弯刚度较小，自振频率较低，需要采取一定的抗风措施来对横向刚度进行加强。

图8某山区加劲梁窄悬索桥跨中竖向位移抖振时程响应

另一方面也能说明抖振为限幅振动，一般不具有灾害性的破坏，但是考虑到发生的频率较高，故也容易引起局部构件疲劳损伤。

4 结 语

本文采用谐波合成与FFT转换技术相结合的方法，模拟了某山区加劲梁窄悬索桥三维脉动风场，考虑平均风引起的静风荷载，脉动风引起的抖振力和主梁气动耦合引起的气动自激力。运用ANSYS大型有限元计算分析程序中APDL编程功能，实现加劲梁窄悬索桥抖振响应时域化分析。计算结果得出，该桥横向抖振位移上限值为16.4 cm，竖向位移振动上限值为8.8 cm，说明某山区横向抗弯刚度更小，出现横向弯曲振型频率会更低，需要采取一定的抗风措施加强横向刚度。另一方面也能看出抖振位移响应时程表现为明显的限幅振动，一般不具有灾害性的破坏，但是考虑到发生的频率较高，故也容易引起局部构件疲劳损伤。