双向流固耦合的圆柱式桥墩绕流的三维数值模拟

任 辉，周亦唐，李 睿

(昆明理工大学建筑工程学院，云南昆明650000)

桥梁是道路跨越河流的主要形式，它起到保障公路运输畅通和宣泄洪水的作用［1］。桥墩建成后，桥位水域阻力增加，过水断面减小，从而对周围的水动力环境产生一定影响。当流体绕过物体时，在一定的情况下物体后部会发生漩涡脱落现象，特别是由此产生的周期性交替的作用力，会诱发结构体振动，导致疲劳破坏。为了提高桥梁抗水毁的能力，研究了桥墩与桥墩周围三维绕流流场的相互耦合作用，对桥梁的安全具有非常重要的意义。

张兆顺，等［2］指出结构体绕流问题呈现出明显的三维效应。结构体绕流不但涉及流动分离、涡旋的生成和脱落，还包括旋涡的相互干扰等许多基本的问题，掌握流体流过结构体障碍物的特性对于许多工程设计问题非常重要，而且对许多实际工程问题也有非常重要的指导作用。笔者以云南省漾濞(跃进)至云龙(诺邓)二级公路1合同跃进段K0+570大桥为研究的工程背景，通过系统的计算机仿真模拟，分析在一定水深、流速、流向等条件下桥墩周围流场的分布情况，为建立桥梁抗水毁能力的评价体系和方法提供参考依据。

1 漩涡脱落的形成与机理

桥墩周围的流场结构主要包括墩前水面涌波、桥墩迎水面向下水流和尺度很大的漩涡。因圆柱表面存在边界层，在柱体上发生的一个重要现象时边界层分离(图1)，边界层外的流体可近似认为是理想流体，因此没有能量损失。边界层外流体由于没有动能损失，仍在向前流动，这样在此区域形成流体的旋转运动，亦即成漩涡，使流体不再贴着柱体表面流动，而是从柱体表面分离出去，并形成向下游延展的自由剪切层。漩涡是柱体左右两侧交替地、周期性地发生的。当一个漩涡向下游泄放，即自柱体脱落并向下游移动时，它对柱体的影响及相应的升力也随之减小直到消失，而下一个漩涡又从对面一侧发生，并产生同前一个相反方向的升力。因此，每一“对”漩涡具有相反的升力，并共同构成一个垂直于流向的交变力的周期，引起结构振动［3-4］。

图1 圆柱体表面边界层分离Fig.1 Boundary layer separation of cylinder surface

2 数学模型

对不可压缩黏性流体，在直角坐标系下，其运动规律可用连续性方程和动量方程来描述［5］。

1)连续性方程(质量守恒方程)。在流场中，流体不断的流进或流出控制体，控制体中所包含的质量可能随时间变化。但是无论流动如何复杂，就控制体而言，必须满足质量守恒定律，即单位时间内通过控制面流入控制体的质量之和等于单位时间内控制体中质量的增量：

式中：Vx，Vy，Vz分别为流体在3个方向上的速度，m/s。

2)运动方程(动量守恒方程)。运动着的流体在任意时刻和任意空间点处的流体质点在各外力作用下，将会以一定方式改变其运动学参数。将动量守恒定律用于运动着的黏性流体质点上，即可得到诸流动参数间的特定关系，其数学表达形式即为Navier-Stokes方程(简称N-S方程)：

式中：ρ为流体密度，kg/m3;μ为流体动力黏性系数，Pa·s;f为单位质量流体体积力，N。

3 工程算例

3.1 模型建立及边界条件

选用漾濞(跃进)—云龙(诺邓)二级公路1合同跃进段上一座40～20 m预应力T形连续梁桥为例。前10孔采用结构简支桥面连续，后30孔采用结构先简支后结构连续，全桥均为5孔一联，共8联，设9道伸缩缝;下部结构，桥墩采用独柱桥墩，桥台采用重力式桥台(图2)。对其圆形独柱式桥墩在三维的空间里进行流体绕流与双向流固耦合的数值模拟。

图2 全桥立面展示(单位：cm)Fig.2 Elevation of the whole bridge

为更好地求解耦合方程，必须合理地定义好流场的初始条件和边界条件。笔者参考已公开发表试验结果的圆柱模型［6-8］，定义如下。

1)结构初始条件：在刚性柱壁面边界条件为流固耦合面处X，Z方向无滑移，Y方向设为Free，并设置弹性约束条件E=300 N/m3(允许Y方向产生位移)。将圆柱表面与流体接触的面设为流固耦合面，即流体计算结果与结构计算结果相互耦合的联系。

2)流体的初始条件：流场参考压强1 atm，水温25℃。流场均匀来流速度为u=0.26 m/s。

3)远场条件：流场远处(流场出口)，平均静态压强P=0 Pa。流场两壁面采用无滑移壁面。圆柱周围的流场利用ANSYS CFX计算流体力学软件求得，计算采用有限体积法及SIMPLE算法，离散格式采用Second Order Upwind二阶迎风格式。

3.2 模型尺寸及网格划分

上游入口须取得足够远以保证流场充分发展，下游出口离桥墩的距离也应取得足够大以保证流场达到稳定［9］。参考已公开发表的模型试验确定模型尺寸如下所述，其中外部流域为200 cm×100 cm，水深为20 cm，流速为0.26 m/s。桥墩直径 D为10 cm。建模时流体区域用长方体模拟，桥墩用圆柱体模拟。桥墩布置于距上游5D处，出口边界选在距离圆柱中心15D处。流体左右两边界距桥墩几何中心(视桥墩为轴对称结构)相等。物理模型和计算网格的划分如图3、图4，网格划分时对桥墩周围和流体交界面附近网格进行加密(采用ICEM进行网格划分，加密区域采用分块实现)，以保证计算精度。

图3 计算的物理模型(单位：cm)Fig.3 Physical model obtained by the calculation

图4 计算网格划分Fig.4 Meshes division obtained by calculation

4 结果讨论

模型数值计算是在CFX求解器下进行，数值计算中时间步长取为0.05 s。图5为由于圆柱桥墩受到流体作用，使得桥墩质点位于z=H/2处质点坐标为(45，0，10)和 z=4H/5 处质点坐标为(45，0，16)的断面位移时程曲线在22 s左右位移呈现周期性变化，图中Ys表示圆柱在Y方向的位移。如图6，在t=28及30 s时，流场运动是相同的，而t=29 s时的流场与t=28 s时的具有明显的对称性。即物体后部两侧产生了周期性漩涡脱落。

图5 位移时程曲线Fig.5 The displacement time history curve

图6 质点位于z=4H/5断面流场速度云图Fig.6 Velocity nephograph when the particle in the plane of z=4H/5

图7为Re=260时圆柱升力系数CL和阻力系数CD随时间的变化曲线图。文中，Z方向平均的升力和阻力系数分别定义为［10］：

式中：FL，FD分别为作用于圆柱上总的升力和阻力。

Strouhal数是描述圆柱绕流的一个重要参数，其定义为［10］：式中：St为Strouhal数，表征旋涡的脱离情况;fs为圆柱涡脱频率;D为发生涡街柱体的迎风面宽度(柱体迎面宽度即等于圆柱体直径);U为来流速度。

图7 Re=260时圆柱升力系数CL和阻力系数CD曲线Fig.7 Curve of lift coefficient CLand drag coefficient CDof circular cylinder at Re=260

由图7阻力系数与升力系数时程曲线可以看出，升力系数时程曲线在32 s左右后显现较为稳定的周期性。对升力时程曲线进行频谱变换得到涡脱的主频率为0.439 1 Hz，并由此计算得到Strouhal数为 0.168 8。

由图8可见，在流动折转较大的区域，流动速度及压力增加，其最大值出现在桥墩迎水面两侧圆心角约70°～－90°范围内。由图9可见在流体正对方向(即前驻点)Cp为1，随着向圆柱两侧扩展，Cp值减小，在60°左右处达到最小，之后又快速增大到一个较稳定值并在圆柱背后形成平坦的分布，由于流体力学将Cp开始平坦分布式的位置作为漩涡的脱离点，这与之前所述的边界层分离即漩涡脱落现象较为吻合，其中Cp为平均压力系数;θ为圆柱表面测点所对应的角度。

图8 桥墩附近动压力分布Fig.8 Dynamic pressure distribution in the vicinity of pier

图9 平均压力系数周向分布Fig.9 Distribution of average pressure coefficients

5 结论

1)采用ANSYS CFX对桥墩附近流场进行模拟，对Re=260的低雷诺数圆柱绕流，采用层流模型并结合合理的网格划分，模拟圆柱表面周向压力系数分布，预测分离点出现的位置，计算了阻力系数、升力系数及Strouhal数。数值计算的收敛曲线大概在700步(30 s)左右后显示稳定的周期性。

2)Re=260的低雷诺数圆柱绕流，由于圆柱表面涡脱在分离交界区的相互作用圆柱升力系数和阻力系数出现周期性的振动。

［1］孙中堂.山区桥梁浅基础水毁分析与防治［J］.山东交通科技，2002(3)：20-21.Sun Zhongtang.Analysis and prevention of shallow foundation of bridge in mountainous area［J］.Journal of Shandong Transportation Science and Technology，2002(3)：20-21.

［2］张兆顺，崔桂香，许春晓.湍流理论与模型［M］.北京：清华大学出版社，2005.Zhang Zhaoshun，Cui Guixiang，Xu Chunxiao.Theory and Modeling of Turbulence［M］.Beijing：Tsinghua University Press，2005.

［3］张阿漫，戴绍仕.流固耦合动力学［M］.北京：国防工业出版社，2002.Zhang A’man，Dai Shaoshi.Fluid Structure Interaction Dynamics［M］.Beijing：National Defense Industry Press，2010.

［4］高东光.公路与桥梁水毁防治技术［M］.北京：人民交通出版社，2002.Gao Dongguang.Innovative Technology for Highway and Bridge Flood Damage Prevention［M］.Beijing：China Communications Press，2002.

［5］江春波，张永良，丁则平.计算流体力学［M］.北京：中国电力出版社，2007.Jiang Chunbo，Zhang Yongliang，Ding Zeping.Computational Fluid Dynamics［M］.Beijing：China Electric Power Press，2007.

［6］Dong S，Karniadakis G E.DNS of flow past a stationary and oscillating cylinder at Re=10 000［J］.Journal of Fluids and Structures，2005，20(4)：519-531.

［7］Lam K，Wang F K，So R M C.Three-dimensional nature of vortices in the near wake of a wavy cylinder［J］.Journal of Fluids and Structures，2004，19(6)：815-833.

［8］Briaud J L，Chen H C，Wang J.Pier and Contraction Scour in Cohesive Soils［R］//National Cooperative Highway Project.Washington D C：Transportation Research Board，2004.

［9］凌建明，林小平，赵鸿铎.圆柱形桥墩附近三维流场及河床局部冲刷分析［J］同济大学学报：自然科学版，2007，35(5)：582-586.Ling Jianming，Lin Xiaoping，Zhao Hongduo.Analysis of three-dimensional flow field and local scour of riverbed around cylindrical pier［J］.Journal of Tongji University：Natural Science，2007，35(5)：582-586.

［10］白桦，李加武，夏勇.低雷诺数圆柱绕流数值模拟及控制措施［J］.建筑科学与工程学报，2010，27(4)：39-42.Bai Hua，Li Jiawu，Xia Yong.Numerical simulation and control measures of flow around circular cylinders at low Reynolds number［J］.Journal of Architecture and Civil Engineering，2010，27(4)：39-42.