大型轴流风机叶片的气动弹性数值分析研究.pdf

机械工程学报 JOURNAL OF MECHANICAL ENGINEERING 第 45卷第 11 期 2009 年 11 月 Vol.45No.11 Nov.20 09 DOI10.3901/JME.2009.11.133 大 型 轴 流 风 机 叶 片 的 气 动 弹 性 数 值 分 析 研 究 * 毛军 1 杨立国 2 郗艳红 1 1. 北京交通大学土木建筑工程学院北京100044； 2. 中国建筑科学研究院建筑结构研究所北京100013 摘要在叶轮机械的流固耦合中以弹性体流固耦合问题最具难度。在这类问题中，需要考虑叶片在运转过程中的弹性变形与 振动对流场的影响以及振荡流场的反作用。以北京地铁 18号专 用轴流风机为例，在考虑叶片变形和流场之间相互影响的耦 合状态下，探讨机叶片气固耦合问题的计算方法。采用 CFX 软件进行流场计算、ANSYS 软件进行结构计算，以 MFX- ANSYS/CFX 为数据耦合平台，应用弱耦合分区法对风机叶片的气动弹性进行数值模拟分析。并与不考虑叶片弹性变形的情 况进行对比，发现考虑气动弹性的最大应力几乎是不考虑气动弹性的最大应力的两倍，相应的安全系数相差较大，不考虑叶 片的气动弹性容易高估叶片结构的安全性。 分析考虑弹性变形对叶片结构安全性的影响，说明叶片的气动弹性分析的必要性 和重要性。 关键词流固耦合弱耦合分区法数值模拟气动弹性风机叶片 中图分类号O35 NumericAnalysis of on the Pneumatic Elasticity of LargeAxial-flow Fan Blade MAOJun 1 YANGLiguo 2 XIYanhong 1 1. School of Civil Engineering, Beijing Jiaotong University, Beijing 100044; 2. Institute of Building Structures, ChinaAcademy of Building Research, Beijing 100013 AbstractIn the fluid-solid coupling of turbomachine,the elastor fluid-solid coupling isthe most difficulty problem.In such kind of problem, it is needed to consider the elastic deation of blade in the operating process,the influence of vibration on the flow field, and the reaction of oscillating flow field. The No.18 axial flow fan in Beijing subway is taken as an example, in the coupling state with the consideration of the interaction between blade deation and flow field, the calculation of gas-solid coupling of fan blade is discussed. REPEATABILITY software is adopted to carry out flow field calculation and ANSYS software for structure calculation, MFX-ANSYS/CFS is taken as data coupling plat, and weak coupling zoning is applied for numerical simulation analysis of the pneumatic elasticity of fan blade. And it is compared to the situation that is without elastic deation considered. It is found that the maximum stress with pneumatic elasticity considered is almost two times of the maximum stress without pneumatic elasticity considered, and the corresponding safety coefficients are in great difference. Without the pneumatic elasticity considered, it is easy to overestimate the safety of blade structure. Through analysis of the influence of elastic deation on the safety of blade structure, the necessityand importance of analyzing the pneumatic elasticity of blade is indicated. Key wordsFluid-solid couplingWeakcoupling zoning Numeric simulationPneumatic elasticityFan blade 0前言* 在大型旋转流体机械中，叶片的稳定性问题由 * 国家科技部重大基础研究前期专项研究5和北京交通大 学科技基金6XM6资助项目。5 收到初稿，6 6 收到修改稿 于直接关系到流体机械的安全运行，成为十分重要 的研究课题。叶片的稳定性取决于作用于叶片的气 动荷载以及叶片的振动变形特性。以大型风机、压 气机等为例，工作气体的非定常流动对叶片将施加 一个非定常作用力，使叶片产生振动，而振动着的 叶片反过来又作用于流场，使流场的非定常程度发 生变化。因此，气固耦合振动是这类问题的特色。 200 CA04700 20002200811120090 1 学兔兔 w w w .x u e t u t u .c o m 机械工程学报第 45卷第 11期期 134 从问题的实质看，叶片不应被视为刚体，而是弹性 体，需要考虑弹性变形与振动对流场的影响以及振 荡流场的反作用，因而对这类问题的研究具有相当 的难度。吴文权[1]采用基于离散涡的方法，把叶片 归结为一个单质点的具有阻尼的弹性系统，来研究 流场与叶片振荡的耦合气动弹性问题。而事实上, 对于三维叶片的弹性结构， 当受到非定常力作用时， 叶片上各点的动态响应非常复杂。黄典贵[2]将叶片 与流场作为一个系统，对理想流场中流固耦合作用 下的叶片动态特性进行了计算，但由于忽略粘性， 并假定空气静止， 因此， 与实际情况有一定的差距。 目前，叶轮机械的气固耦合即气动弹性问题成为流 固耦合的一个重要分支[3]，由于在试验上存在较大 困难，主要采用计算流体力学的方法进行分析[4]， 试验间接验证。 本文以北京地铁 18号轴流通风机为 例，对风机叶片的气固耦合计算进行探讨，对其在 正常运转条件下的气动弹性进行数值分析，初步发 现了气动弹性对风机叶片结构强度的影响。旨在为 大型风机的故障诊断和安全性分析提供参考，对风 机叶片的动力学特性设计也具有重要意义。 1耦合计算方法 求解流固耦合问题一般采用两种耦合求解法 强耦合法和分区强耦合或弱耦合求解法，强耦合法 的关键在于构造出控制方程，通过改写流体、结构 控制方程使其成为同一种形式，然后对控制方程直 接求解。直接法概念清楚，它主要用于流体结构 耦合场的理论分析。利用这种方法求解问题时必须 进行某些适当的简化，因而难以保证计算的准确程 度。此外，要求流场、结构网格在网格上一致，而 在实际问题中，流场分析和结构分析数值计算所要 求的网格往往相差甚远。因此该方法并不适用于实 际工程应用中。 分区求解法分为强耦合分区法和弱耦合分区 法。强耦合分区求解法是利用导出的耦合项完全隐 含的隐式公式， 在每一个时间步求解非线性方程组， 计算出全场变量的值。弱耦合分区求解法是在每一 时间步内分别对计算流体动力学Computational fluid dynamics,CFD 方 程 和 计 算 结 构 动 力 学 Computational structural dynamics, CSD方程依次 求解。即通过 CFD 求解 N-S 方程来模拟气动力， 通过耦合高精度的有限元结构分析即 CSD 来模拟 结构动力响应进行联合求解。结构场和流场的计算 结果通过搭建中间数据交换平台彼此交换信息，从 而实现两个场的耦合求解。此方法优点在于在求解 过程中FD 计算和SD 计算模块相互独立，保持 了 CFD 和 CSD 计算模块的完整性，充分发挥在各 自领域的长处，实现数值模拟的模块化。这种计算 方法相对成熟，结果也有较高的精度，但同时计算 过程相当复杂。这种方法也是目前用得最多的流固 耦合气动弹性分析方法[4]。 这里采用弱耦合分区法， 在考虑叶片变形和流场之间相互影响的耦合状态 下，对风机叶片的气动弹性进行数值模拟分析。 2计算模型的建立 2.1流固耦合平台 CFX 流场计算和 ANSYS 结构计算之间的耦合 平台是 MFX-ANSYS/CFX。 通过该平台把流场中的 压力计算结果和结构场中的位移计算结果进行交 换，进行耦合计算。即在流体区域内求解非定常流 体方程，在固体的区域内求解固体方程，在每一个 时间步长后，将求解流体方程得到的气动力传递给 固体，将求解固体方程得到的叶片位移作为边界条 件传递给流体，然后进行下一步迭代。 MFX-ANSYS/CFX 耦合平台是 ANSYS 公司 2005 年针对CFD 软件CFX 和 CSD软件 ANSYS 新 开发出来的耦合平台，嵌在 ANSYS 软件中，真正 实现了 CFD/CSD 之间的耦合计算。 2.2流场模型 本文采用流固耦合的数值方法研究地铁风机 叶片的变形和流场之间的相互影响。因此准确地模 拟移动边界和移动网格下的非定常流场对于流固耦 合数值研究非常关键。图 1 给出了风机的主要结构 尺寸，叶片的旋转直径为 1 800 mm，叶片偏长，超 过 500 mm一般叶片长 400 mm；叶片偏薄，由叶 片根部到顶部逐渐变薄， 中线最小厚度不超过 6mm。 图 1风机的几何结构和尺寸 2.2.1几何模型及网格划分 风机的动叶片为三维叶片，从根部到叶顶的安 装角和工作角都是变化的， 共 14 片； 静叶片为二维 叶片，共片，如图所示。由于结构的几何对称 CC 92 学兔兔 w w w .x u e t u t u .c o m 月 2009 年 11 月毛军等大型轴流风机叶片的气动弹性数值分析研究 135 性，并考虑到耦合计算的计算量大，这里对几何计 算模型进行了简化， 计算中转子采用 1/14 的叶道几 何，即选取一个流道进行计算。为了减小计算时进 口及出口位置对叶轮内部流态的影响，适当将计算 域的进口边界及出口边界向进口及出口延伸。 图 2风机的三维动、静叶片和轮毂及风筒 采用专业三维建模软件建模，三维非结构化网 格。相对于结构化网格，非结构化网格计算过程比 较复杂，但局部加密比较容易，对不规则空间适应 能力强，易于显示流场的细微结构。高质量的网格 非常重要，应该避免过密或过疏的网格[5]。 由于叶片非常薄，叶片边缘的尺寸相对非常 小，在划分网格过程中，网格尺寸大小的均匀过渡 很困难。应用专业网格工具 ICEM 进行网格划分， 采用面和体结合生成的方法生成网格。先对叶片边 缘小尺寸面进行面网格划分，然后生成体网格。该 方法在保证叶片边缘小尺寸面的网格质量的同时， 也很好地控制体了网格的数量和质量，叶片表面的 面网格数为 22 188，计算域的体网格数为 311 387， 如图 3、4 所示。 图 3叶片几何模型和表面网格 图计算域体网格划分 2.2.2数学模型 对于湍流的模拟，目前在工程中普遍采用雷诺 平 均 Reynolds-averagedNavier-Stokes equations, RANS方法[6]，其基本思想是建立经验公式或方程 使雷诺方程封闭，从而可对其求解。本文采用在工 程计算中应用较为广泛的基于 RANS 模拟的标准 k-ε模型， 标准 k-ε模型是基于湍动动能k和耗散率 ε两输运方程的半经验模型，控制方程式为 i i kku tx ρρ t k jkj k G xx ρ ε σ 1 i i u tx ρ ερ ε 2 12 t k jj CGC xxkk εε ε εεε ρ σ 2 式中ρ为气体的密度；为气体的动力粘度系数； t为紊流粘性系数；Gk为由于平均速度梯度引起的 湍动动能生成项， 常数C1ε1.44,C2ε 1.92，C0.09， σ k1.0，σε1.3，其取值一般通过试验结果确定。 边界条件由该风机在额定转速980 r/min运行 的测试数据给出。气流进口一般要求给出叶轮进 口的速度、压力、密度或相应的相容条件。本文采 用压力入流条件，定义进口静压力相对大气压力为 零，即没有附加的压力作用。气流出口取质量出 口边界条件，整个风道的流量为72 kg/s，取一个流 道计算即为5.14kg/s。 根据实际工况建立的计算模型。为了考虑旋转 的影响，将转动参考系设在于叶轮一起旋转的流体 区域上，计算域绕x轴旋转，转速– 102.625 rad/s。 风道内壁、轮毂、叶片均取无滑移壁面边界条件。 轮毂、叶片均设置移动壁面，旋转速度为0rad/s相 对于坐标系；风道内壁设置移动壁面，设置旋转速 度为102.65 rad/s。定义叶片表面为流固耦合交界 面，所传输的变量为压强。计算域两侧取轴对称周 期性边界条件。 计算条件是初始流场为定常流场；在后面的 耦合计算中，流场为非定常流场。 2.3叶片结构的力学模型 2.3.1几何模型 用ANSYS中提供的实体单元对风机叶片进行 模拟，如图5所示。叶片与轮毂的连接处为位移约 束，其他的面定义为流固耦合交界面。转动中心在 风机的轴线处。叶片的材料为铸铝合金。通过试验 测定的有关力学性能参数为弹性模量63 5 4 E. 学兔兔 w w w .x u e t u t u .c o m 机械工程学报第 45卷第 11期期 136 GPa；泊松比0.313ν；质量密度ρ2 627 kg/m3。 图 5有限元计算的叶片网格 2.3.2数学模型 结构的动力响应由以下运动方程式描述[7] t mu ku cuF3 式中 m 为质量矩阵，k 为刚度矩阵，c 为阻尼矩阵。 Ft表示作用在结构上由流体施加的荷载，表示表 面压力和切应力。对于膜结构，切应力通常不可以 忽略；然而对于其他结构，由于表面所受流体施加 的切应力量值非常小，且结构本身对切应力反映不 敏感，因此可以忽略切应力影响，而只考虑垂直于 表面的流体压力作用。结构的动力响应由结构的位 移 u，速度 u和加速度 u描述。值得指出的是，流 体阻尼包含在式3的右边项中，不包括在阻尼矩 阵c中。 2.4耦合数据的传递模型 基于弱耦合分区方法进行气动弹性数值模拟 的时间步迭代推进计算流程如图6所示，图6中竖 向表示时间步推进过程，水平向代表当前状态，箭 头指向表示计算流程。 图 6时间步迭代推进算法示意图 图6中列出了两种时间步迭代推进算法[7]，它们 的区别在于迭代计算的开始端， 按照图6a的时间步 迭代推进步骤，计算从CFD模块开始，经过一个时 间步长t， 将计算得到的结构风荷载p返回给CSD 模块做瞬态过程分析；再经过一个时间步长t， SD模块计算得到的结构运动速度返回给FD模 块，FD网格发生变形，在新的状态求解流场，如 此反复迭代直至达到设定的收敛标准。 图6b的迭代 计算从CSD模块开始，经过一个时间步长t，将 计算得到的结构运动速度 u 返回给CSD计算模块。 作为数值迭代计算方法， 当时间步长t取得足够小 时， 这两种时间步迭代推进方法所取得的结果一致。 3耦合计算关键问题的处理 3.1CFD 模型和 CSD 模型变量单位的统一 在CFX模型中，变量单位可以指定，但在 ANSYS中，不能指定变量的单位，只要求变量单 位统一。 所以在CFX和ANSYS中单位一定要统一。 3.2耦合面的定义 CFD/CSD耦合计算非常关键的一步就是耦合 面的定义，要分别在CFD软件和CSD软件中定义 相对应的耦合面。 这里CFX中用边界名称定义交接 面，ANSYS中用数值符号定义交接面，在ANSYS 中定义交界面的对应情况。 3.3耦合时间步长的确定 耦合时间步长的给定非常重要。时间步长给的 太小，将极大地增加计算时间。而耦合时间过大， 不能正确反映气动弹性分析的真实情况。对于耦合 时间步长，本文采取等于瞬态动力计算综合时间步 长Integrated time step，ITS办法。瞬态动力计算综 合时间步长通过结构的自振频率来给定，一般为所 关心模态最高频率20倍的倒数，即1 20Sf。其 中S是瞬态动力计算综合时间步长，f是自振频率。 4计算结果及分析 4.1叶片的自振模态计算 通常，在进行结构瞬态动力分析之前要进行结 构自振模态分析。通过模态分析来定义荷载的瞬态 动力计算综合时间步长。下表是叶片的自振模态频 率的模拟计算结果。 表模型前 6 阶自振频率的计算结果 模态阶数 n 123456 频率 f/Hz 30.25067.087182.230278.440299.550475.990 通过模态计算，选取一阶模态为参考，耦合计 算时间步长取为0.001 5 s在本文中耦合时间不长， 迭代100步，共计0.15 s的时间。 4.2叶片的气动力分析 图7和图8给出了定常流不考虑耦合叶片的 叶盆和叶背压力分布图。 从图7、8中可以看出， 叶 盆以正压力为主，压力分布很不均匀，变化范围很 大，最大压力主要出现在叶片顶部前端；叶背主要 CC C 学兔兔 w w w .x u e t u t u .c o m 月 2009 年 11 月毛军等大型轴流风机叶片的气动弹性数值分析研究 137 以负压为主，同时也出现明显的正压区，压力分布 也很不均匀，变化范围很大，最大负压主要出现在 叶片顶部前端。因此，由于流动的复杂性和气固 耦合的作用，叶片表面的压力分布非常复杂，极不 均匀。 图 7定常流下叶盆气动压力分布图 图 8定常流下叶背气动压力分布图 图9、10给出了耦合第100步时叶片叶盆和叶 背压力分布图。从图9、10中可以看出，叶片叶背 表面的气动压力发生了明显的变化，压力的最大值 明显地减小，最小值明显增大。同时，叶盆压力的 分布规律也发生了明显的变化。 图耦合第步叶盆气动压力分布图 图 10耦合第 100 步叶背气动压力分布图 4.3叶片结构的响应分析 从初始条件开始计算， 得到0～0.15 s耦合时间 段内的位移响应曲线，如图11～14所示。无论 顺风向x方向或横风向y方向，叶片结构位移响 应曲线频率都比较单一。随着时间的推进，由于阻 尼的作用，位移幅值逐渐减小。结构在顺风向和横 风向的位移的响应频率相同，均为39.23 Hz，接近 叶片的一阶固有频率，但有一定差距，并不产生 共振。 图 11叶片顶部一点 x方向位移时程响应 图 12叶片顶部一点 y方向位移时程响应 图 3叶片中部一点 x方向位移时程响应 9100 1 学兔兔 w w w .x u e t u t u .c o m 机械工程学报第 45卷第 11期期 138 图 14叶片中部一点 y方向位移时程响应 4.4叶片的应力和位移分析 图11～14的气动弹性计算结果表明， 在不出现 共振的情况下，叶片顶部的最大位移绝对值出现 在耦合的第6步，即第0.009 s6 0.001 5 s。此时， 叶片结构的最大应力为25.8 MPa，出现在叶片根 部，是最不利的情况，如图15、16所示。故对此进 行进一步分析。 图 15耦合第 6 步叶背应力图 图 16耦合第 6 步叶盆应力图 这里， 叶片材料的许用应力150～240 MPa，许 用安全系数[ ]5n [8]。 在考虑气动弹性影响的情况下，最大应力为 25.8 MPa，计算得到的最小安全系数为 min 150 5.81[ ] 25.8 nn≥ 在不考虑气动弹性影响的情况下，最大应力为 13.4 MPa[9]，计算得到的最小安全系数为 min 150 11.2[ ] 13.4 nn≥ 说明该型地铁专用轴流风机在额定转速下运 转时，不论是否考虑气动弹性，叶片在运转过程中 都能满足强度要求。但考虑气动弹性的最大应力几 乎是不考虑气动弹性的最大应力的两倍。在考虑气 动弹性影响的情况下，安全系数已经接近许用安全 系数，一旦由于偶发因素，安全系数可能超过许用 安全系数而导致叶片断裂。虽然本文的数值模拟没 有经过试验的严格验证， 但计算分析与图17所示的 叶片从根部断裂的情况比较吻合。因此，从某种意 义上说，得到了间接的证明，结果是有一定可信 度的。 图 17某地铁车站该型风机的叶片断裂情况 从叶片根部断裂 可见，在不考虑气动弹性影响的情况下，安全 系数远大于许用安全系数， 这与实际情况相差较大， 易夸大叶片结构的安全性。因此，在风机的结构和 运行工况设计中，有必要考虑气动弹性的影响，对 风机叶片特别是长而薄的叶片进行气动弹性分析。 5结论 本文阐述了基于MFX-ANSYS/CFX耦合平台， 应用弱耦合分区法进行风机叶片气动弹性数值模拟 分析的过程和方法。计算表明以下几点。 1应用基于MFX-ANSYS/CFX耦合平台的弱 耦合分区法可以有效地计算风机叶片的气固耦合特 性。在计算过程中，耦合面的定义和耦合时间步长 的确定非常关键。 2对于长而薄的风机叶片，气动弹性对叶片 表面的压力分布规律和压力大小以及结构的位移响 应都有明显的影响，结构在顺风向和横风向的位移 的响应频率相同。 3对于长而薄的风机叶片，考虑气动弹性的 最大应力几乎是不考虑气动弹性的最大应力的两 倍，相应的安全系数相差较大。不考虑叶片的气动 学兔兔 w w w .x u e t u t u .c o m 月 2009 年 11 月毛军等大型轴流风机叶片的气动弹性数值分析研究 139 弹性易高估叶片结构的安全性，因此，在校核轴流 风机叶片结构的强度时应考虑气动弹性的作用。 在这里，仅对风机叶片的一个流道的流场进行 了计算，未对全通道流场进行计算，后者对计算机 硬件和软件有更高的要求。 参考文献 [1] 吴文权. 叶栅气动弹性离散涡数值仿真 II 数值试验[J]. 工程热物理学报，1994，154372-377. WU Wenquan. Numerical simulation for aero-elasticity in turbo-machines with vortexIn II, numerical tests[J]. Journal of Engineering Thermophysics, 1994, 154372-377. [2] 黄典贵. 理想流场中流固耦合作用下叶片的动态特性 研究[J]. 气轮机技术，1998，404235-238. HUANG Diangui. Study on the dynamic characteristics for fluid-solid coupling problems in theoretical flow fields[J]. Turbine Technology,1998, 404235-238. [3] 邢景棠，周盛，崔尔杰，等. 流固耦合力学概述[J]. 力 学进展，1997，27119-38. XING Jingtang, ZHOU Sheng, CUI Erjie, et al. A survey on the fluid-solid interaction mechanics[J]. Advances in Mechanics,1997,27119-38. [4] 杨伟. 基于 RANS 的结构风荷载和响应的数值模拟研 究[D]. 上海同济大学，2004. YANG Wei. Numerical simulation of wind load on atructures based on RANS[D]. ShanghaiTongji Univer- sity，2004. [5] 田铖，张欢，由世俊，等. 利用 FLUENT 软件模拟地 铁专用轴流风机的内部流场一对称翼叶片轴流风 机[J]. 流体机械，2003，311113-15. TIAN Cheng, ZHANG Huan, YOU Shijun, et al. CFD simulation of axial-flow fan in subway by FLUNT software1axial-flow fan with symmetrical blade[J]. Fluid Machinery, 2003,311113-15. [6] 王福军. 计算流体动力学分析CFD 软件原理与应用 [M]. 北京清华大学出版社，2004. WANG Fujun. Computational fluid dynamic analysis- principle and application of CFD software[M]. Beijing Tsinghua University Press, 2004. [7] LI Zhenyin. Parallel computations of 3D unsteady com- pressible euler equationswith structural couping[D]. AmericaPurdue University, 2002. [8] 商景泰. 通风机手册[M]. 北京 机械工业出版社， 1996. SHANG Jingtai. Manual for ventilator[M]. Beijing China Machine Press,1996. [9] 毛军，薛琳，王丹，等. 关于地铁轴流风机叶片发生断 裂和出现裂纹的分析[J]. 中国安全科学学报，20035 75-80. MAO Jun, XUE Lin, WANG Dan, et al. Analysis of the rupture and flaw happened to vanes of axial fans in subway[J].China Safety Science Journal,2003575-80. 作者简介毛军，男，1966 年出生，博士，副教授。主要研究方向为流 体力学与流体机械以及安全技术与工程。主要科研成果包括地铁 18 号轴流风机叶片失效分析及对策、 地铁通风系统的优化设计、地铁列车 新型离心式通风风机、地铁主通风系统节能运行与变风量控制方案； 地 铁火灾过程的非线性规律研究、北京地铁火灾应急预案、强侧风对列车 运行安全性影响等。曾主持和参加国家、铁道部、北京市科委、北京交 通大学等资助的科技基金项目和企业委托的科研项目。 E-mailjunmao2003 学兔兔 w w w .x u e t u t u .c o m