基于CFD仿真的平板非线性气动力系统特征研究

祝志文 袁涛+陈政清+邓燕华��

摘 要：为研究平板气动力系统的非线性特征，基于非定常雷诺时均NavierStokes （RANS）方程和SST kω湍流模型，数值模拟了在单位位移激励下平板非定常运动的绕流场，获得了作用在平板上的气动力时程，并基于Volterra理论开展了平板非线性气动力系统识别.研究表明，本文建立的平板非线性气动力模型能对一定頻率带宽和一定幅值范围的激励产生合理的响应；在本文研究的强迫运动位移幅值和频率范围内，平板非线性气动力模型响应没有表现出对振动幅值和频率的明显相关性，且其气动力的非线性效应并不明显，因而可以认为小攻角下的平板绕流属于气动力弱非线性系统.本文研究证明了CFD模拟在桥梁主梁气动力系统识别上的明显优势.

关键词：气动力非线性；平板；Volterra理论；CFD

中图分类号：U448.213 文献标志码：A

不断突破桥梁跨径记录的跨江跨海工程需求和新材料、新结构形式的应用，可能会使得大跨度桥梁的频率和阻尼不断降低，因而导致结构或其构件对风的敏感程度大大增加.合理评价大跨度桥梁气动性能，准确预测桥梁风致振动稳定性，对保障大跨度桥梁设计合理和运营安全具有深远意义.然而，现有桥梁抗风研究涉及的风荷载，均是建立在气动力的线性模型之上.以气弹分析广泛采用Scanlan自激力模型为例[1]，对二自由度桥梁断面而言，其气动自激力表达为位移、速度与颤振导数的线形组合形式，也即气弹自激力是通过线性叠加多项而得到，其中每一项与位移或速度有关，实际反映的是不同广义位移对应的相对攻角所产生的效应，这是一种典型的线性颤振理论.因此，Scanlan这一桥梁颤振理论也是基于线性和小扰动理论，没有考虑桥梁气动外形对气流场可能产生的大扰动，以及大跨度桥梁大幅运动对气流场的强烈干扰，如流动分离、再附和涡系的相互作用.这些流动现象是桥梁气流场的非线性行为，将导致气动力响应呈现非线性特性，这在Scanlan的早期研究中就有阐述，后续研究者也有类似报导[2].流动非线性对大跨度桥梁气弹响应和颤振稳定性的影响目前并不清楚.因此，对气动外形明显钝化的桥梁主梁，或大幅度振动的桥梁结构，研究气动力的非线性，有助于给出大跨度桥梁风致振动稳定性的合理评价.

3.2 计算域网格及流动条件

包含平板主梁模型的计算域及分区如图4所示.采用圆形计算域，计算域外圆到平面中心的距离均为16B，对应的模型堵塞度不到0.2%，因此无需考虑模型堵塞度带来的影响.采用计算域分区划分以控制网格的正交性和网格缩放比.在平板模型外第一个分区为椭圆形“刚性网格区”， CFD模拟模型运动时，该区域网格与模型刚性固定，在每一时间步上与模型同步运动.该刚性网格区外边界为椭圆，对该区域作四边形结构网格划分.计算域的绝大部分区域为“静止网格区”，该区域外边界是计算域外边界，内边界为离开刚性网格区外椭圆一定距离，且包围刚性网格区的圆形.静止网格区采用四边形单元剖分，从内到外采用合适的网格放大比例.静止网格区和刚性网格区在CFD模拟过程中一直使用计算开始时的网格系统.在“静止网格区”和“刚性网格区”之间为“动网格区”，“动网格区”采用三角形单元剖分.在每一时间步上，“动网格区”根据桥梁断面的运动位置并由设定的网格系统质量要求可能需重新进行网格划分.紧靠桥梁断面的区域流场变量变化剧烈，特别是断面迎风侧和断面法向，因此网格划分适应流场变量的变化程度，并沿各个方向采用适度的网格放大率，实现与动网格区域网格的平顺过渡，如图5所示.计算域的网格总数为96 226.

边界条件设置为：计算域入口水平均匀速度边界，湍流度为零；下游出口施加流动出口边界条件；模型表面使用无滑移壁面条件；初始场采用入口速度初始化.数值计算采用非定常二阶隐式格式，采用速度压力解耦的SIMPLE算法，二阶格式离散压力方程，动量、湍动能和湍流耗散率方程均采用二阶迎风格式，数值模拟的Re=3.08×105，本文所有数值模拟均基于CFD专用程序Fluent 6.3.26开展.

3.3 单位采样激励的实现

在CFD模拟中，连续型Volterra级数无法使用，需采用离散型的Volterra级数形式.对小攻角下平板模型，因其有较好的气动外形，气动力中的非线性可能并不显著.对此类弱非线性气动力系统，本文尝试识别其一阶核.

本文考虑桥梁竖向运动，计算由此产生的非线性气动力.根据本文公式（4）～（5）定义的离散时间域内的单位脉冲和二维单位脉冲信号，编制UDF程序实现单位位移激励.图6为单位采样函数，按照定义，作为位移量的u0的值应该取1，也即为很短时间内强迫桥梁断面的竖向位移量.但当赋予单位采样函数实际物理意义时，本文将u0取值为0.002B.为了适应数值计算，定义在τ时刻给予计算模型单位采样激励，而在τ+Δt时刻回到初始位置，其中Δt为非定常CFD计算时间步长，本文为0.005 s.由于单位位移激励作用的时间极短，这种极短的位移激励在现有的风洞试验条件下几乎无法实现，但这在CFD模拟中可以方便地实现，这表明了CFD模拟的巨大优势.

系统的一阶核是时间的一维函数，一阶核表示非线性系统的线性子系统响应，它可以反映一定的非线性行为，与纯粹的线性系统响应是不一样的.由公式（10）可知，在识别系统的一阶核时，并不需要延时的激励作用.CFD模拟时，强迫平板模型在很短的时间快速向上运动（为正），此时，模型将受到为负的气动升力作用.

图7为平板模型非线性气动升力系数的一阶核，前3个时间步出现波动，这是突加模型运动导致的气动力振荡.图8为平板模型非线性气动扭矩系数的一阶核，也在前两个时间步出现波动.这种剧烈的波动与CFD模拟的时间步长和运动幅值有较大关系，过小的时间步长和过大的运动幅值会加剧这种波动，可能引起CFD模拟的数值误差.因此，识别核时需要尽可能减少这种剧烈波动，本文选取了最佳的时间步长和运动幅值.从图7和图8时程的局部放大可知，一阶核很快就衰减，说明某一时间步的模型运动，仅对后续系统的气动力产生有限时间步的影响.

3.5 气动力模型仿真

设平板非线性气动力系统属于时不变系统.对弱非线性气动力系统，如不考虑二阶及以上核，在获得了非线性系统的一阶核后，就可建立对应的非线性气动力系统模型.该系统对于任何输入u[t]的响应w[t]，都可通过卷积公式（2）求得.为验证本文非线性气动力模型的合理性，可对气动力模型进入输入激励，观察模型的气动力响应.本文竖向运动激励采用正弦位移输入[12]：

u（t）=h0sin （2πfht） （18）

式中h0和fh分别为正弦位移输入的幅值和频率.

为考虑平板模型运动幅值和频率的非线性特征，本文位移输入的幅值考虑大幅和小幅情况，对应的幅值分别为1×10-2B和5×10-2B，本文分别称之为小幅运动和大幅运动状态；模型运动的频率分别取1 Hz和4 Hz，也即对应低频运动和高频运动状态，运动幅值和频率的组合分别称之为小幅低频运动和大幅高频运动.显然，后者可能比前者更容易激发出气动力的非线性特征.

从小幅低频运动和大幅高频运动下的气动力输出对比可见，基于一阶核建立的平板非线性气动力模型，对竖弯单频激励和一定幅值范围的激励均具有较好的适用性.另外，在本文研究的强迫运动位移幅值和频率范围内，没有表现出对振动幅值和频率的相关性，且其气动力的非线性效应并不明显，这也证明，0°攻角下的平板绕流，属于气动力弱非线性系统.

4 结 论

本文尝试借助CFD模拟技术通过Volterra理论开展了平板非线性气动力系统识别，可得下述结论：

1）模拟时间极短的单位位移激励作用下的平板模型運动，这种复杂位移激励在现有的风洞试验条件下几乎无法实现，但在CFD模拟中可方便地实现，这表明了CFD模拟在桥梁气动力系统识别上的巨大优势.

2）平板气动力系统的非线性模型，可基于Volterra理论并借助CFD模拟手段得到识别，本文识别了平板气动力系统的非线性一阶核，由此建立的非线性模型能对竖弯单频激励和一定幅值范围的激励均具有较好的适用性；

3）在本文研究的强迫运动位移幅值和频率范围内，平板非线性气动力模型没有表现出对振动幅值和频率的相关性，且其气动力的非线性效应并不明显，也即0o攻角下的平板绕流，属于气动力弱非线性系统.

对钝化的桥梁断面和主梁大幅度运动情况，有待开展进一步的研究，以揭示实际桥梁主梁气动力系统的非线性特征.

参考文献

[1] SCANLAN R H， TOMKO J J. Airfoil and bridge deck flutter derivatives [J]. Journal of Engineering Mechanics， ASCE， 1971， 97（6）：1171-1737.

[2] CHEN Z Q， YU X D，YANG G ， et al. Windinduced selfexcited loads on bridges [J]. Journal of Structural Engineering， 2005， 131（12）：1783-1793.

[3] ZHU Z W. LES prediction of aerodynamics and coherence analysis of fluctuating pressure on box girders of longspan bridges[J]. Computers and Fluids， 2015， 110： 169-180.

[4] ZHU Z W， GU M， CHEN Z Q. Wind tunnel and CFD study on identification of flutter derivatives of a longspan selfanchored suspension bridge[J]. ComputerAided Civil and Infrastructure Engineering， 2007， 22：541-554.

[5] ZHU Z W， GU M. Identification of flutter derivatives of bridge decks using CFDbased discretetime aerodynamic models [J].Wind and Structures， An International Journal， 2014，18（3）： 215-233.

[6] VOLTERRA V. Theory of functionals and of integral and integrodifferential equations [M]. New York： Dover Pubns， 2005：63-105.

[7] SILVA W A. Discretetime linear and nonlinear aerodynamic impulse responses for efficient CFD analyses [D]. Williamsburg， Virginia：College of William & Mary， 1997： 25-68.

[8] MENTER F R. Twoequation eddyviscosity turbulence models for engineering applications [J].AIAA Journal， 1994， 32（8）：269-289.

[9] 祝志文，夏昌，邓燕华.钝体矩形断面绕流场机理与主分量分析[J].湖南大学学报：自然科学版，2012，39（7）：7-13.

ZHU Zhiwen，XIA Chang， DENG Yanhua. Flow mechanisms around bluff rectangular cylinders and its principle component analysis[J]. Journal of Hunan University： Natural Science， 2012，39（7）：7-13.（In Chinese）

[10]ZHU Z W， CHEN Z Q， GU M. CFD based simulations of flutter characteristics of ideal thin plates with and without central slot [J]. Wind and Structures， 2009， 12（1）：1-19.

[11]GU M， ZHANG R X， XIANG H F. Identification of flutter derivatives of bridge decks [J]. Journal of Wind Engineering and Industrial Aerodynamics， 2000， 84（2）：151-162.

[12]祝志文，陈政清，顾明.薄平板气动导数的数值仿真技术[J].湖南大学学报：自然科学版，2005，32（5）：11-15.