燃气涡轮叶片传热特性分析Ⅰ：燃气涡轮叶片内冷通道管网计算 刘明林.doc

中国工程热物理学会 热机气动热力学学术会议论文 编号：112007燃气涡轮叶片传热特性分析：燃气涡轮叶片内冷通道管网计算 刘明林 李涛 薛兴旭 颜培刚 黄洪雁 韩万金 （哈尔滨工业大学458信箱，推进理论与技术研究所，哈尔滨，150001）（电话 E-mail：）摘要：本文通过FORTRAN编程建立了考虑旋转、变截面、摩擦以及加热的管网计算模型。该计算模型可以用来预测涡轮叶片冷却通道中压力、温度以及流量的分配情况，从而迅速的发现局部高温区以及流量分配是否合理，然后改进冷却结构，提出优化方案。文中详细给出了管网计算的发展现状、计算原理、求解方法和计算步骤，并建立了一套适用于光滑冷却通道的阻力系数和换热系数的经验关联式，将带有实验数据和三维计算的MarkII导向叶片进行管网程序的验证，结果吻合较好，从而可以为叶片的热分析、强度以及寿命计算提供原始数据。关键词：燃气涡轮叶片；一元流动；管网计算0 引言 管网计算是根据管内一元流动理论发展起来的一门科学，流体在管内流动时会发生传输与瞬变现象，这种现象存在于很多实际的生产领域，如流控技术和气动液压技术中的动态过程分析；水电站设计和运行中的水击现象研究；天然气、石油的输运过程；航空发动机的通气系统与飞机燃油系统的设计等。因而发展管网计算的数值算法意义是十分重大的。将管网计算应用到航空发动机的通气系统中，国内学者起步较晚，现在许多高校和科研院所都投入了大量的精力对管网计算进行研究。1996年西北工业大学刘松龄等提出的管网计算是采用流量残差修正压力的方法，在收敛性和计算精度上都有较大的改进。2005年南航郭文应用压力修正方法求解动量方程，并与能量方程耦合计算，使方程组线性化，节省了计算时间，提高了计算稳定性。2009年北航陶智等应用压力修正方法、流热耦合计算方法，并考虑了阻力系数、换热系数、物性参数等的修正，在此基础上开发了程序，并应用到有经验关联式的模型中进行验证，最后将经过验证的程序应用到某发动机空气冷却系统初步设计方案的校核和改进中，在此基础之上又发展了流固耦合一体化和非稳态管网计算模型。2010年哈尔滨汽轮机厂冯永志等采用简化的动量方程，对F级燃气轮机转子空气系统的压气机某级抽气管路进行了计算，结果模拟较好。1 管网建模思想一个完整的管路系统包括各种尺寸的管道、各种形式的管道附件，若分别对他们进行编程计算，计算量太大，迭代起来比较慢。根据质量连续方程，串连管路内质量流量相同，因此，可将管网划分成各段独立支路（或称为节流单元）进行处理，从而加快了收敛速度。管网进出口称为边界节点，这些节点处的压强或流量为已知的；除边界节点以外，其它节点均称为内节点，在计算管网内流量分布以及压力分布时，内节点处的压强和流量未知，支路中流体从高压力节点流向低压力节点。如下图一、图二分别为简化的冷却结构和对应的管网模型。图一 某叶片局部冷却结构图二 叶片局部冷却结构管网模型管网分析操作步骤如下1、针对具体冷却结构划分节点，形成管网计算模型。2、对具体问题进行一次三维计算，提取内冷通道壁面温度。3、求出各节点元件的阻力系数、换热系数，计算流通的横截面积、侧面积、距离转子中心的距离等各种参数。4、给定边界条件：入口静压，静温，出口静压。5、为避免每次迭代的误差过大，引入欠松弛因子来保证稳定性和收敛性。6、管网计算完毕之后提取相关数据与三维计算进行对比。2 管网计算原理利用一元流动的连续性方程、动量方程和能量方程来耦合求解压力和温度。2.1 控制方程2.11 动量方程考虑旋转、变截面、摩擦以及温度的动量方程1（可以认为是各种因素线性叠加而得）(1)对上式经过一阶差分得：(2)其中：当时，时。式中、为未知量，摩擦阻力系数根据通道几何参数与Re数由经验公式求出，节点处温度作为已知参数给出。2.12 连续性方程连续性方程组在内部节点上建立，即对于i节点，有质量流量平衡方程：(3)式中，表示节点i到节点j的节流单元的质量流量，两节点之间没有节流单元时。对于给定流量的边界节点，方程变为，表示外部流入i节点的流量。2.13 能量方程能量守恒方程说明两节点之间的总焓变化等于控制体本身动能的变化与从外界吸收的热量之和；通过节流单元形式的能量方程和节点的能量方程就可以把各节点的温度求解出来。公式（4）3和（5）分别是节流单元的两种不同形式的能量方程。(4)其中：其中为定压比热容，分别为各节点所对应的总焓，为总温，为静温，为冷气与内冷通道壁面的换热系数，为内冷通道壁面的温度，为冷气与壁面的换热面积。节流单元流体的平均温度取进出口节点温度的算术平均值。冷气侧换热系数根据经验公式求出，换热系数的经验公式一般为或。对上式进行化简可以得到节流单元能量方程的另外一种形式：(5)节点的能量方程形式为：(6)2.2 管网计算的求解方法在管网计算中如果直接通过解控制方程组来求解压力和温度，会导致方程组非线性化，计算不容易收敛。为了简化求解过程，可以把求解分为压力平衡计算和温度平衡计算两部分考虑，压力与温度计算交替进行，直至计算收敛。本算例利用连续迭代的方法求解代数方程组，因而在计算之前需要对被求解的温度、压力和流量提供一个合理的初场，经过验证发现本管网程序对初场的依赖性很大，为解决这一问题可以将入口和出口边界的物理量均匀的分配到各个节点上，这样就容易使计算进行下去。其中使用连续迭代法的求解步骤如下： 1、 建立控制方程：.如果有n个未知变量，则需要n个方程。 2、猜测方程的解：预估控制方程组的初始解，即初场的给定。3、取代方程右边的，计算新变量, 。4、利用松弛因子Under-relax来求解新变量 ，作为下一次求解的初始值。5、用步骤4中的代替。6、重复步骤35，直到变量变化值达到预定精度为止。3 摩擦阻力系数和换热系数的计算管网分析的主要目的是用来对比不同的冷却结构方案，迅速的舍弃一些不合理的冷却方案，所以在建立管网计算模型时可以首先忽略叶身冷气通道中扰流结构的具体形式，可以依据冷却通道全部为光滑管的情况进行计算，再根据计算结果，在适当位置添加强化换热的扰流结构，以控制流量与改善叶片温度场。在管网计算中发现冷却通道中的阻力系数与换热系数的数值对流量和温度具有很大的影响，这两种系数的计算基本上都是依据经验公式，而经验公式的选取与叶片的冷却结构密切相关，因而合理的针对具体冷却结构来选取经验公式直接影响计算的准确性，下面具体给出了针对圆形断面管道摩擦阻力系数和换热系数的计算。3.1 冷气侧摩擦阻力系数的计算阻力系数的定义式为，经验公式的一般形式为，或给出摩擦阻力系数，两者的换算关系为。不同冷却结构的阻力系数不同，流阻系数的计算有两种方法：一是根据经验公式来计算2；二是对于无经验公式可依的冷却结构，可以根据三维计算来获得。对于圆形管道中的层流流动，摩擦阻力系数可以用数学的推演方法计算出准确表达式，并且已为实验所证实。而对于湍流，摩阻系数的计算只能借助于实验以求得经验或半经验公式。下面给出了几种不同流动状态下摩阻系数的经验公式4。1、层流状态（）(7)2、转捩状态（）该状态下的摩擦阻力系数通常按湍流状态下的布拉休斯公式计算。3、湍流状态（）布拉休斯公式（8）4、湍流状态（任意的）(9)3.2 冷气侧换热系数计算一、迪图斯-贝尔特（Dittus-boelter）换热系数准则式5： (10)其中Nu为努塞尔数，Cr为修正系数，Pr为普朗特数，Re为雷诺数，冷却通道流体雷诺数Re由冷气流量(q)、通道当量直径(Dh)和冷气动力粘性系数(u)确定。式（11）给出了流体雷诺数计算的表达式。(11)普朗特数由冷气定压比热()、粘性系数()和导热系数()确定。如式（12）所示。(12)冷却通道壁面对流换热系数与努塞尔数的关系式为(13)适用条件及注意事项：1、实验验证范围：雷诺数，长径比 l/d60，只适用于水力光滑管。2、此公式是流体在管道中湍流充分发展的关联式，如果用于计算Re=230010000的过渡区中的对流换热，则会得出偏高的结果。3、流体与壁面具有中等温差，即对气体不超过50K，对于水不超过2030K。定性温度为流体平均温度（可以用管道进出口两个截面温度的算术平均值来代替），取管内径为特征长度。二、Gnielinski公式5(14)其中符号说明：为管长；为管道直径；为管内流动的Darcy阻力系数，为温差修正系数。适用条件及注意事项：1、实验验证范围： 说明Gnielinski公式可以应用于过渡区。2、Gnielinski公式是迄今计算准确度最高的一个关联式，综合考虑了温差以及长径比的限制，适用于水力光滑管，对于粗糙管可以做初步的估算，阻力系数按粗糙管数值计算。3.3 物性参数的处理在涡轮叶片冷却通道中，由于流动工质工作在压力和温度不断变化的条件下，而工质的物性参数与温度和压力是密切相关的，为了保证管网计算能够和CFX三维计算相对应，对程序中涉及到工质的一些物性参数做了以下处理。1、计算换热系数和滞止参数时需要用到定压比热容，对于空气当t600K时,可以取比热比，但是当温度超出这一界限时可以根据空气在理想状态下的比定压热容与温度的关系式来求。(15)其中，适用温度范围为2731800K，最大误差为0.72%。2、因为冷气处于高压的作用下，因而密度变化是很明显的，但是为了方便处理数据，当马赫数小于1时在一个节流单元可以近似的认为是常数，根据其与压力和温度的关系式进行计算。3、为了保证程序中燃气和冷气的导热系数、动力粘度与CFX三维计算中的计算方法相对应，可以引用CFX中自带的萨德兰公式(16)(17)其中：，4、叶片的导热系数与材料有关，一般为1774 ，在本算例中叶片材料选用ASTM标准的310不锈钢（0Cr25Ni20），导热系数比较低（约为不锈钢的1/4），其物理属性为，密度=8030 kg/m3，比热Cp=502J/kgk其热传导系数K为温度的线性函数：K=0.0115T+9.91054 算例验证MarkII导向叶片是国外文献中公开的带有实验结果与三维仿真对比的一种典型的冷却结构6，因而以该算例进行验证更具有代表性，只要管网计算结果能够与实验结果很好的对应，那么就可以认为管网程序的可用性，在第二部分中会给出MarkII导向叶片的三维计算过程，这里只进行管网的计算和结果对比分析。4.1 管网计算MarkII导向叶片该冷却结构共有十根冷却管道，为了详细分析每根冷却管道中温度、压力和流量沿径向的分布情况，对每根管道人为分成六个节点（如图五）、五个节流单元，所以该管网共有六十个节点，五十个节流单元。程序中需要给定的边界条件是入口的静压、静温，出口的静压；另外还要给出节点处得当量直径以及节流单元之间的长度和侧面积等几何信息，至于内冷通道壁面的温度可以从三维计算结果中提取。 图五 单个管道管网计算的模型4.2 结果对比与分析符号说明：onedime表示管网计算结果，threedime表示CFX三维计算结果，test 5411 表示MarkII导向叶片的5411实验工况所获得的实验数据，element表示叶片中的管道数目1、压力在各管道节流单元上的变化图六 静压在各管道中的变化图七 总压在各管道的变化从上图可以看出两种形式的压力在各自管道中都成下降趋势，管网计算和三维仿真计算的静压和总压几乎重合在一起，这说明了管网程序中压力的求解是非常准确的。2、温度在各管道节流单元上的变化图八 静温在各管道的变化图九 总温在各管道中的变化管网计算的静温和总温都比三维仿真的结果高，这与管网计算所引用的换热系数关联式相关，这可以尝试增加温差效应对换热准则式进行修正；另一方面在提取内部冷却通道壁面温度时，由于给出的是壁面的平均温度，而不是每个节流单元壁面平均温度，这也会造成较大的误差。3、流量的对比图 图十 流量在各管道中的分配情况图十一 流量在各管道的误差图十给出了三维计算、管网计算和实验各冷却管道中的流量情况，图十一是管网计算对比于实验的流量误差，由图十一可以看出流量在各管道中的误差在10以内，因而能够满足实际工程的需要。因为阻力系数的大小直接影响了冷却管道中流量的大小，管网计算得到的流量与实验和三维计算的吻合说明了阻力系数的经验关联式可以用到光滑管道中。4、各冷却通道换热系数对比图十二 换热系数在各管道的变化图十三 不同准则式所对应的换热系数对比图十二对比了管网计算、三维计算和实验各冷却通道中换热系数的变化，可以发现三维计算得到的换热系数偏高，这是因为在三维计算中换热系数的数值大小与网格的质量和湍流模型的选取有密切的关系，因而可以尝试从改进网格质量和湍流模型来获得比较准确的换热系数；通过管网程序计算的换热系数与实验所得到的结果比较接近，这也说明了所选择的换热系数准则式比较适合光滑直管道换热的计算。图十三对比了管网程序中用两种换热系数准则式计算的换热系数，可以看出两者吻合较好，用Gnielinski公式（hc2）计算的换热系数比Dittus-boelter（hc2）计算的换热系数稍低，说明这两种换热系数准则式可以用在光滑管道中。5 小结1、 通过对比管网程序和三维计算结果，静压和总压吻合情况较好，误差都在5%以内。2、 管网程序中计算得到的温度相对三维计算较高，造成这种现象的原因是冷却通道中的壁面温度是按照三维计算管道中的平均温度来处理的，并没有严格的按照节流单元所对应壁面的平均温度来输入，这可以尝试与其它程序或者接口来自动划分管网模型并且提取所需要的数据，这样就既能保证管网计算的精确性，又提高了计算的效率。3、 相对于实验结果，管网程序中计算的流量误差在10%以内，而且能够明显的看出不同管道流量的变化情况，这种变化的趋势可以给冷却结构设计者作为参考