CFD方法的大跨度桥梁抖振荷载

刘婷婷，张文首

（大连理工大学 工业装备结构分析国家重点实验室，辽宁 大连116024）

风荷载是大跨度桥梁的主要环境荷载，其中，抖振荷载作为一种长期作用的随机荷载，会使桥梁结构局部构件和连接部位处于交变应力作用，影响桥梁疲劳寿命，是大跨度桥梁抗风设计的重点问题。

目前的抖振荷载计算主要依靠风洞试验方法。在一定缩比的桥梁节段模型上进行定常流风洞试验，测定抖振力三分力气动系数，然后带入有限元模型（通常为鱼骨头模型）进行分析［1－5］。但是，风洞试验方法的周期长、费用高、结果的可视性较差，并且，无法针对桥梁的精细有限元模型进行直接分析，制约了大跨度桥梁抖振分析方法的发展。计算流体力学（Computational Fluid Dynamic，CFD）的发展使数值模拟方法可以替代物理的风洞试验，能够节约试验成本、进行重复试验，弥补风洞试验所测数据不足［6－8］。

采用CFD数值模拟替代传统的风洞试验，改进了大跨度桥梁抖振荷载计算。在香港青马大桥精细有限元模型上，模拟桥板表面压力分布，根据CFD模拟结果确定节点抖振力，由此求解桥板局部应力响应时程。CFD模拟能够获得桥板表面任何位置压力，弥补了风洞试验安装测压传感器有限的缺点。而且，可以针对桥梁精细有限元模型，直接确定节点抖振力，将风致抖振分析范围推进到局部构件，为桥梁风振疲劳研究提供了有力的工具，促进了桥梁健康监测及局部抗风设计的发展。

1 抖振力

1．1 传统抖振力

根据拟定常假定和片条原理，长度为L桥板节段弹性中心的阻力Dbfe，升力Lbfe及弯矩Mbfe可以表达为（如图1）［9］：

式中，ρ是空气密度，U是平均风速，B为桥板宽度，u（t）和w（t）为水平及竖直方向脉动风速时程。CD、CL和C M是与攻角α相关的无量纲三分力系数，主要通过风洞试验在桥板弹性中心测得，C′D＝dCD／dα、C′L＝ dCL／dα、C′M＝ dCM／dα。χDbu，χDbw，χMbu，χMbw，χLbu，χLbw是气动导纳函数，当长度L较小时可看作常数1。

图1 桥板弹性中心及节点抖振力

1．2 桥板压力分布的CFD模拟

实质上，采用以上传统方法计算抖振力时，就是将作用在桥板上的压力合成为弹性中心的合力，因此，在风洞试验中只关注弹性中心的气动参数，而忽略了桥板的压力分布及其产生的局部作用。另外，随着有限元建模技术的发展，已经可以建立更精细的桥梁有限元模型，以适应复杂分析的需要。例如，为了研究交通荷载引起的大跨度悬索桥疲劳问题，文献10建立了青马大桥的大型有限元模型，采用梁单元模拟桥板的桁架构型，并建立桥板节点的局部模型；针对青马大桥健康监测系统——WASHMS（Wind and Structure Health Monitoring System），更准确的青马大桥有限元模型也被实现，并应用到风振引起的桥梁疲劳问题的预测［11－12］。与鱼骨头模型不同，采用精细有限元模型进行抖振计算时，需要的是单元节点的集中抖振力，这是单一采用风洞试验方法无法获得的。

桥板风压分布是计算单元节点抖振力关键，常用方法是制作布置了压力传感器的桥板节段模型，通过风洞试验获得。CFD技术的发展使数值模拟方法可以替代物理的风洞试验，并能够节约试验成本、弥补风洞试验测量数据的不足。计算流体力学CFD的思想是：把时间及空间上连续的物理量（速度或压力），用有限个离散点上的变量集合代替，通过流动基本方程建立关于这些离散点上变量关系的代数方程组，然后进行求解获得近似值。通过这种对流体流动的数值模拟，可以得到复杂问题的流场内各个位置上的基本物理量（速度、压力、温度等）的分布，甚至这些物理量随时间变化的情况。

对于桥板压力分布，CFD数值模拟的对象是桥板周围流场，可以采用二维定常不可压缩控制方程［13－14］，考虑湍流效应，引入标准k－ε双方程湍流模型。

其中，u和v分别是流场内沿水平轴和竖直轴的速度，ρ是空气密度，p为压力，k和ε分别是是湍流动能和湍流动能耗散率；有效粘性系数υeff＝υ＋υ1，其中，υ为运动粘性系数，是湍流粘性系数；系数为：Cμ＝0．09、C1＝1．44、C2＝1．92、σk＝1．0以及σε＝1．3。

1．3 节点抖振力

若通过CFD数值模拟已知桥板压力分布，可将压力转换为作用在桥板表面梁单元的节点力：

其中，˜Fl1和˜Fl2是第l个单元两端节点力，nl是CFD

式中，p k1和p k2为计算网格节点处压力，l k是网格长度。

将式3计算的单元节点力由单元坐标转换为p－h－α坐标后，叠加各相邻单元节点力，就可将桥板压力等效为各单元节点力，其中第j个节点力为：是作用的节点力总数。则节点气动参数可以表示为：

式中，是CFD模拟中的平均风速。

最后，式1可以改写为节点抖振力形式：

2 精细有限元模型的结构运动方程

在全局坐标系x－y－z中，在抖振力作用下，桥梁运动方程为：

其中，a（t）＝ ｛u1，w1，…，um，wm｝T是2m维脉动风速向量，m为桥板节段总数；Pbf为6N×2m气动参数矩阵，由每个桥板节段的气动参数矩阵Q ibf组成。

大跨度桥梁的有限元模型的自由度较大，计算量十分庞大，通常采用模态叠加法进行分析。节点位移X（t）可表示为：模拟计算中将此单元划分的网格数量；和是在第k个网格末端的压力合力（如图1），其表达为：

式中，Φ＝ ［Φ1，Φ2，…，ΦNΦ］是6N×NΦ模态矩阵，q（t）＝ ｛q1（t），q2（t），…，q NΦ（t）｝T是广义位移向量，NΦ是计算中采用的模态总数。式7的模态运动方程为：

模态运动方程（10）可采用直接积分法求解广义位移时程q（t），然后由式9确定节点位移。如果已知某单元模态应力Γr，此单元应力为：

3 算例研究

3．1 青马大桥精细有限元模型

香港青马大桥是连接市区和大屿山机场的主要通路，是世界目前最长的公铁两用悬索桥之一，主跨达1 377 m。青马大桥三维精细有限元模型采用商业软件MSC NASTRAN建立（图2）。有限元模型的几何构型与真实桥梁一致，采用三维2节点梁单元模拟加劲桁架，主索及吊杆用圆截面梁单元模拟，桥面采用正交各向异性板单元，用等效刚度方法考虑加劲肋对桥面刚度影响。2个桥塔是多室空间结构，采用空间梁单元进行划分，并对其刚度进行换算。模型边界条件完全符合实际情况。对该模型进行了动力特性分析，与实测结果符合很好［6］。

图2 青马大桥三维精细有限元模型和主梁3种典型截面

3．2 CFD数值模拟

青马桥的桥板压力分布计算采用商用软件ANSYS FLOTRAN进行计算。为了模拟定常流风洞试验，分析域和边界条件设定如图3所示。其中，流体分析域采用二维流体单元FLUID141。在入口边界上，来流均匀分布，湍流强度小于1%。青马桥主梁主要有3种截面形状，对应有限元模型分别是主跨12节点截面、青衣塔10节点截面和8节点截面（图2）。其中，12节点桥板截面为主梁的主要构型，计算网格划分如图5所示，另外2个截面网格划分与其相似。

图3 计算域与边界条件

图4 12节点桥板截面CFD计算模型

3．3 计算参数

将主梁桥板沿长度方向分为120个18 m长的节段，根据CFD数值模拟获得的桥板压力分布，由式3－6计算桥梁精细有限元模型节点抖振力。其中，桥板水平和竖直方向脉动风速时程采用谐波合成法进行模拟［15］。风速谱分别采用Simiu顺向风谱和Lumley－Pnofsky竖向风谱，即：

图5 青马大桥主跨中点脉动风速时程（平均风速U＝18 m／s）

式中，u＊为摩擦速度，该文为1．15 m／s；n为频率，单位是Hz，f（z）＝nz／U为折减频率。主梁离海平面平均高度z＝60 m，模拟的时程采样频率为S50 Hz，频率上限为25πHz，频率段数为214，时长为10 min；Davenport形式的相关函数中的衰减因子取为16。在本算例中，分别对强风和一般风速下的桥梁应力响应进行分析，取主梁的平均风速分别为40 m／s、18 m／s和8 m／s，图5是平均风速为18 m／s时，主跨中点顺风方向u（t）和竖直方向w（t）的风速时程。

对青马桥局部应力进行数值分析时，采用模态叠加法，用MSC．Nastran对桥梁精细有限元模型进行动力分析，主要考虑前80阶位移模态和应力模态，模态阻尼比为1%。采用Newmark－β法求解广义位移q（t），其中β＝0．25。

3．4 应力响应结果及分析

在青马大桥精细有限元模型中，共有15 904个梁单元模拟主梁桁架。通过对所有梁单元应力计算结果分析比较后，发现青马大桥局部应力较大的位置主要集中在：桥墩 M2（139．5 m）、马湾侧桥塔（495 m）、以及青衣侧桥塔（1 872 m）（如图6所示）。图7所示为马湾侧桥塔的应力较大的杆件单元，主要是主梁桁架的下部顺桥向、主跨一侧的水平杆件，以及下部横向杆件。其他两个位置的重点单元部位与马湾侧桥塔相似，但杆件应力值与重点杆件的数量均不及马湾侧桥塔位置。

图6 青马大桥危险部位

图7 马湾侧桥塔重点单元位置

表1是马湾侧桥塔的主梁桁架重点杆件，单元中部的应力标准差计算值，包括了梁单元截面上部和下部。从应力标准差结果可以看出，不同的平均风速下，单元应力随平均风速增加而增大，并呈现平方增长趋势，因此，在强风或台风时，更加需要关注局部杆件的应力变化。

外侧桁架下部的顺桥向水平杆件单元（单元34111）的单元应力最大。在平均风速40 m／s，即通过青马桥监测数据推算的120 a一遇最大风速，可达6．53 MPa，因此，在抖振力的长期作用下，更容易产生疲劳损伤，在监测维护时，需要对这类杆件进行重点监护，图8为单元34111在平均风速18 m／s时，截面下部单元应力时程。

该文计算的是梁单元中部应力值，在进行疲劳分析时，还需要考虑杆件连接部位应力集中问题。通过杆件单元应力计算结果，可确定危险杆件位置，然后参考钢结构疲劳规范获得应力集中系数，或者对连接部位建立更具体的局部有限元模型，进行热点应力分析。

图8 单元34111截面下部单元应力时程

4 结论

采用CFD数值模拟方法替代传统风洞试验，获得大跨度桥梁主梁表面压力分布，确定节点抖振力，在香港青马大桥精细有限元模型上，求解了桥板局部应力响应。在不同风速下，对计算整个主梁桁架单元应力，确定了青马桥应力较大、需要关注的重点位置和重点杆件单元，这些单元均处于桥塔与主梁连接部位。这些关键杆件的具体位置及单元应力响应的确定，为桥梁抖振引起的疲劳分析以及结构健康监测设计奠定了基础。

［1］项海帆．现代桥梁抗风理论与实践［M］．北京：人民交通出版社，2005．

［2］SUN D K，XU Y L，KO J M．Fully coupled buffeting analysis of long－span cable－supported bridges：formulation［J］．Journal of Sound and Vibration，1999，228：569－588．

［3］XU Y L，KO J M，ZHANG W S．Vibration studies of Tsing Ma suspension bridge［J］．Journal of Bridge Engineering，1997，2（4）：149－156．

［4］ZHU L D，XU Y L．Buffeting response of long－span cable－supported bridges under skew winds．Part 1：theory［J］．Journal of Sound and Vibrations，2005，281：647－673．

［5］ZHU L D，XU Y L．Buffeting response of long－span cable－supported bridges under skew winds．Part 2：case study［J］．Journal of Sound and Vibrations，2005，281：675－697．

［6］SHINICHI KURODA．Numerical simulation of flow around a box girder of a long span suspension bridge［J］．Journal of Wind Engineering and Industrial Aerodynamics，1997（67／68）：239－252．

［7］SHUJI SHIRAI， TOSHIO UEDA． Aerodynamic simulation by CFD on flat box girder of super－long－span suspension bridge［J］．Journal of Wind Engineering and Industrial Aerodynamics，2003，91：279－290．

［8］SHIGERU WATANABEA，KOICHIRO FUMOTO．Aerodynamic study of slotted box girder using computational fluid dynamics［J］．Journal of Wind Engineering and Industrial Aerodynamics，2008，96：1885－1894．

［9］XU Y L，GUO W W，CHEN J，et al．Shum and H．Xia，Dynamic response of suspension bridge to typhoon and trains．II：numerical results［J］．Journal of Structural Engineering，2007，133（1）：12－21．

［10］LI Z X，ZHOU T Q，CHAN T H T，et al，Multi－scale numerical analysis on dynamic response and local damage in long－span bridges ［J］． Engineering Structures，2007（29）：1507－1524．

［11］LIU T T，XU Y L，ZHANG W S，et al．Buffetinginduced stresses in a long suspension bridge：structural health monitoring oriented stress analysis［J］．Wind and Structure，2009，12（6）：479－504．

［12］XU Y L，LIU T T，ZHANG W S．Buffeting－induced fatigue damage assessment of a long suspension bridge［J］．International Journal of Fatigue，2009，31：575－586．

［13］约翰D．安德森．计算流体力学基础及其应用［M］．吴颂平，刘赵淼．译．北京：机械工业出版社，2007．

［14］BLAZEK J．Computational Fluid Dynamics Principles and Applications［M］．Oxford：Elsevier Science Ltd．

［15］CAO Y H，XIANG H F，ZHOU Y．Simulation of stochastic wind velocity field on long－span bridges［J］．Journal of Engineering Mechanics，2000，126（1）：1－6．