高速飞机的全机气动力数值分析

　　A 辑第 19 卷增刊 　　　　　　水 动 力 学 研 究 与 进 展 　　　　　　Ser. A , Vol. 19 ,Supplement 　　　2004 年 12 月 　　　　　　　　JOURNAL OF HYDRODYNAMICS　　　　　　　　　　Dec. , 2004 文章编号 :100024874 (2004)增刊20838211 高速飞机的全机气动力数值 分析 定性数据统计分析pdf销售业绩分析模板建筑结构震害分析销售进度分析表京东商城竞争战略分析 Ξ 李向群 , 安亦然 , 王世安 程暮林 , 张锡金 , 孙惠中 (北京大学力学与工程科学系 , 北京 100871) 摘 　要 : 　利用 StarCD 和 FLUENT对一真实飞机模型进行了大量工况下的数值模拟 ,并将计算结果与试验结果进行了对 比。根据计算的结果 ,论证了试验中存在的一些问题。最后将两种成品软件对同型号的飞机模拟结果对比 ,说明成品软件计 算结果的准确性。 关 　键 　词 : 　NS方程 ;全机数值模拟 ;动力影响 ;并行计算 中图分类号 : 　V212. 11 　　　文献标识码 :A Numerical aerodynamic analysis of an high speed airplane LI Xiang2qun , AN Yi2ran , WANG Shi2an CHENG Mu2lin , ZHANG Xi2jin , SUN Hui2zhong (Dept . of Mechanics , Peking University , Beijing 100871 , China) Abstract :　In this paper the real airplane computation was taken for massive flight conditions by using commercial code StarCD and Flu2 ent , and the results were compared with relevant wind tunnel testing. According to the computation and analyzing , some wind tunnel testing defects were found. In order to make sure of the reliability of the code we used the part of cases was redone by another CFD code StarCD , and the consistency of the results was approved. Key words :　NS2equation ; numerical simulation ; parallel computation 1 　引言 在现代飞行器的设计[4 ]中 ,模拟计算和风洞试 验、空间试飞一起构成了获得飞行器气动力数据的 三种手段。采用计算流体力学[3 ]方法可缩短周期、 降低费用 ,特别在初选阶段 ,优化选型需要不断改变 参数、重复计算。对那些目前不能在特定的飞行状 态下进行试验的未来飞行器来说 ,数值模拟方法可 以减少其设计风险 ,并在风洞试验前预先筛选设计 方案 气瓶 现场处置方案 .pdf气瓶 现场处置方案 .doc见习基地管理方案.doc关于群访事件的化解方案建筑工地扬尘治理专项方案下载 。然而 ,在流体力学[7 ]中 ,不论是用欧拉方程或 是 NS方程对飞机进行全机数值模拟计算都是一个Ξ 收稿日期 : 　2004202215 作者简介 : 　李向群 (1978～) ,男 ,湖南衡阳人 ,硕士。 © 1995-2005 Tsinghua Tongfang Optical Disc Co., Ltd. All rights reserved. 比较复杂的问题。在计算机技术充分发展以前 ,全 机模拟很难进行。现在有了大规模并行计算技术和 大型计算机以及针对高速流动研究的计算流体力学 方法 ,对全机数值模拟不但可行 ,而且在技术先进的 国家中占整个气动设计和分析工作的百分之八十到 九十。[1 ,2 ] 本文利用成品软件 ,对我国某一型号战斗机进 行了各种情况下的全机数值模拟。与已有的试验数 据相比表明计算的精度足以实际应用 ,并能根据计 算数据的演变规律发现试验异常。在证实数据可信 的基础上 ,计算了大量目前无法进行风洞试验的工 况 ,为新型飞机性能设计提供了依据。本次第一轮 计算工况超过 300 个 ,是目前国内最大规模的飞机 气动力计算。 2 　飞机几何外形 这次数值模拟的飞机正是自我国自行设计的某 型号战斗机。在进行数值模拟的过程中 ,飞机几何 外形完全保留 ,不对原几何数据作任何简化 ,保留了 机身上所有的配件。 3 　主控方程 对全机进行数值模拟的主控方程可以采用欧拉 方程 ,也可以采用 NS 方程[6 ,7 ] 。对升力和力矩的计 算采用欧拉方程即可提供足够的精度 ,欧拉方程计 算收敛速度相对较快 ,结果与 NS 方程差别不大。 但欧拉方程无法准确计算阻力 ,而 NS 方程对升力 阻力以及力矩均能计算得很好 ,但收敛速度相对较 慢。 对全机进行数值模拟所采用的欧拉方程为 : 连续性方程 :5ρ5 t + ý ·(ρv→) = 0 (1) 动量守恒方程 :55 t (ρv→) + ý ·(ρv→ v→) = - ý p (2) 能量方程 : 55 t (ρE) + ý ·( v→ (ρE + p) ) = - ý ·( h 5 T5 xi ) (3) E = CpT (4) 其中 p 为静压 ,ρE 为内能 , E为比焓 , h 为热传 导系数。 对全机进行数值模拟所采用的雷诺应力平均 NS方程为 : 5ρ5 t + 55 xi (ρui) = 0 (5)55 t (ρui) + 55 xj (ρuiuj) = - 5 p5 xi + 55 xj [μ(5 ui5 xj + 5 uj5 xi - 23δij 5 ul5 x l ) ] +55 xj ( - ρu′i u′j ) - ρu′i u′j = μt ( 5 ui5 xj + 5 uj5 xi ) - 23 (ρκ+μt 5 ui5 x i )δij (6) 相应湍流模式采用 SA 模式。 4 　网格结构 网格由 ICEM ANSYS[8 ]生成。ICEM与多种 CAD 软件具有相容接口 ,可以直接导入 CAD 软件中的几 何外形数据。本飞机几何外形由 CAD 软件 CATIA 生成。ICEM与 CATIA 有专门接口 ,可以保证导入 的过程中对原几何形体不会有任何改动。ICEM 自 动生成非结构化网格的功能十分强大。对飞机这种 十分复杂的几何形体 , ICEM可以生成高质量的非结 构化网格。 计算区域采用半径为 10 倍机身长度的圆球体 , 整个流场采用非结构化网格。 飞机机身网格示意图见图 1 (图 1 与真实几何 差异较大) 。 938李向群等 :高速飞机的全机气动力数值分析 © 1995-2005 Tsinghua Tongfang Optical Disc Co., Ltd. All rights reserved. 图 1 　飞机表面网格示意图 5 　计算过程 采用耦合求解法 ,利用有限体积法离散控制方 程 ,采用二阶迎风差分 格式 pdf格式笔记格式下载页码格式下载公文格式下载简报格式下载 [3 ,6 ] 。计算过程中监控全 机的升阻力以及力矩的变化趋势 ,升阻力以及力矩基 本不再变化后结束计算。具体判据为 :升阻力以及力 矩系数有效数字前两位保持不变。 共计算了马赫数分别为 0. 6、0. 8、0. 9、1. 05、 1. 2、1. 5 ,攻角分别为 - 4 度、0 度、4 度、8 度、12 度 ,侧 滑角为 0 度、3 度、6 度、8 度不带动力影响的情况 ,以 及马赫数 0. 2、0. 6、0. 85、1. 5 ,攻角分别为 0 度、4 度、8 度、10 度、12 度 ,侧滑角为 0 度、3 度、5 度和马赫数 0. 95、1. 05、1. 2、1. 6 ,攻角分别为 - 2 度、0 度、4 度状态 下带动力影响的情况。 来流均为可压缩理想气体。 外边界条件按“远场压力”[9 ] ,[10 ]设置。 升阻力以气流坐标为准[5 ] 。 6 　并行计算 由于需要计算的情况相当多 ,而且对于这种网格 数目巨大的计算 ,采用普通的计算机计算显然不能完 成。因此需要多台计算机并行计算或利用超级计算 机/ 并行机群等进行计算。本计算主要使用北京大学 “科学与工程计算研究中心”的大型 IBM 计算机和力 学系微机机群进行。使用机群并行计算时考虑所使 用的软件存在并行计算时结点越多利用率可能越低 的问题 ,对本项具体计算来说并不是结点数目越多就 越有利。在经过 8 结点、16 结点、32 结点等一系列测 试后 ,我们认为对本课题采用 16 结点计算时综合效 率最高。结点数目再增加 ,计算速度没有明显变化。 采用 16 结点 ,对一种飞行状况下的计算大约需要 24 小时 ,而采用 32 结点对同样的状况进行计算 ,机时仅 降到 20 小时。并行计算过程中需要将整体计算网格 按照一定的方式分割成 16 个区域 ,然后指定参与计 算的结点计算机和并行方式。并行方式对并行计算 的效率影响很大 ,我们采用的是 NetMPI 并行方式。 这种方式高效而且十分稳定 ,比较适合多台计算机之 间的并行计算。网格分区对计算的影响也很大。平 行于某一方向的分区会导致网格子区域疏密程度差 别太大 ,计算过程中则表现为某些结点计算机负载过 大而其他结点计算机处于空闲等待数据状态。我们 分区原则是尽量以飞机为中心 ,向四面散开 ,这样可 以提高并行计算的效率。 7 　无动力状态下计算结果分析 7. 1 　与试验结果的对比 计算结果与风洞试验结果基本一致 ,升力系数曲 线、极曲线、力矩曲线与试验[5 ]所得的曲线趋势完全 一致。 以上曲线均利用最小二乘法拟合而成。升力系 数曲线采用线性拟合 ,力矩曲线和极曲线采用二次抛 物线拟合。从曲线对比看 ,除 Ma = 1. 5 的工况以外 计算与试验符合得相当好 ,关于 Ma = 1. 5 时差异的 原因将在后面详细讨论。 计算与试验之间的误差参见表 1。 048 水 　动 　力 　学 　研 　究 　与 　进 　展 　　　　　　　　　2004 年增刊 © 1995-2005 Tsinghua Tongfang Optical Disc Co., Ltd. All rights reserved. 表 1 　误差表 马赫数 攻角 升力误差 阻力误差 力矩误差 0. 6 0 9. 15 % 6. 25 % 9. 27 % 0. 6 4 0. 29 % - 3. 72 % - 8. 35 % 0. 6 8 0. 19 % 1. 55 % - 0. 30 % 0. 6 12 0. 87 % - 2. 10 % - 2. 13 % 0. 8 0 - 237. 67 % 7. 49 % 25. 00 % 0. 8 4 - 0. 30 % 14. 71 % - 8. 58 % 0. 8 8 1. 14 % 8. 65 % - 0. 02 % 0. 8 12 2. 85 % 6. 84 % - 4. 60 % 0. 95 0 141. 76 % 4. 79 % 21. 31 % 0. 95 4 7. 01 % 14. 34 % 24. 90 % 0. 95 8 4. 27 % 8. 96 % 27. 07 % 0. 95 12 3. 13 % 13. 76 % 22. 43 % 1. 05 0 - 6. 62 % 1. 38 % - 67. 41 % 1. 05 4 - 4. 10 % 14. 24 % 9. 38 % 1. 05 8 - 2. 86 % 5. 79 % 11. 86 % 1. 05 12 - 0. 27 % 13. 45 % 19. 14 % 1. 2 0 60. 81 % 8. 96 % 133. 94 % 1. 2 4 - 0. 29 % 8. 59 % 16. 06 % 1. 2 8 - 2. 07 % 6. 12 % 11. 71 % 1. 2 12 - 3. 41 % 10. 57 % 5. 55 % 1. 5 0 1. 35 % 4. 46 % - 431. 79 % 1. 5 4 6. 80 % 6. 91 % 36. 18 % 1. 5 8 4. 77 % 6. 32 % 25. 12 % 1. 5 12 - 21. 28 % - 20. 88 % - 8. 59 % 图 2 　Ma0. 6 升力曲线图 (左)及 Ma0. 6 极曲线图 (右) 148李向群等 :高速飞机的全机气动力数值分析 © 1995-2005 Tsinghua Tongfang Optical Disc Co., Ltd. All rights reserved. 图 3 　Ma0. 6 力矩曲线图 (左)及 Ma0. 95 升力曲线图 (右) 图 4 　Ma0. 95 极曲线图 (左)及 Ma0. 95 力矩曲线图 (右) 图 5 　Ma1. 5 升力曲线图 (左)及 Ma1. 5 极曲线图 (右) 表 中 误 差 计 算 采 用 公 式 : 误 差 = 计算数值 - 试验数值 试验数值 ×100 % ,负值表示计算数值大 于试验数值 ,正值表示计算数值小于试验数值。 各项误差中 ,0 度攻角误差最大 ,尤其是升力系 数误差。上述公式误差实为相对误差 ,在试验数值很 小时不能反映真实精度。在 0 度攻角下 ,无论是计算 结果还是试验结果 ,升力系数都很小。从这个角度来 看 ,0 度攻角下升力系数的试验值和计算值相差并不 大 ,计算结果也是合理的。为此 ,苏联中央空气动力 研究院有专门 规定 关于下班后关闭电源的规定党章中关于入党时间的规定公务员考核规定下载规定办法文件下载宁波关于闷顶的规定 :苏联中央空气动力研究院关于气 动阻力容许误差的规定说明在小攻角时可以绝对误 差为准。例如该院的亚跨音速大风洞 T - 128 是一座 兼容性很广的有五个专用试验段的大型风洞 ,用以研 248 水 　动 　力 　学 　研 　究 　与 　进 　展 　　　　　　　　　2004 年增刊 © 1995-2005 Tsinghua Tongfang Optical Disc Co., Ltd. All rights reserved. 究飞行器以及其部件在尽可能接近飞行条件下的气 动特性和气弹特性 ,可以进行低速、高速、跨音速测力 及测压的全机和部件试验等等。该风洞是目前世界 上精度较高的风洞。该风洞的测量精度为 : 阻力系数 :σcx = ±0. 0002 升力系数 :σcy = ±0. 003 俯仰力矩 :σmz = ±0. 0015 马赫数 0. 8、攻角 0 度时的升力系数试验值约为 - 0. 005 ,计算值约为 0. 007 ,与风洞测量精度 ±0. 003 为同一量级 ,这时误差完全可以以绝对误差 为准。 图 6 　Ma1. 5 力矩曲线图 7. 2 　效率因子和零升阻力分析 基于位势流动的小扰动理论[7 ] ,在小迎角范围内 (攻角小于 4 度) ,升阻力系数之间的关系可以用效率 因子 e[4 ] 和零升阻力系数 Cd0 来表示 ,即 : Cd = Cd0 + Cl2 πAe 式中 A 为飞机的展弦比 , Cd 为阻力系数 , Cl 为 升力系数。 利用0度和4度攻角下的升阻力系数可以求出效 率因子 e 及零升阻力系数 Cd0。 在正常情况下 , 效率因子 e 应小于 1 ,尤其是对 于展弦比较小的战斗机。 当马赫数为 0. 6 时 ,试验的结果分别为 :0 度攻 角下升力系数约为 - 0. 008 ,阻力系数约为 0. 018 ;4 度攻角下升力系数约为 0. 22 ,阻力系数约为 0. 025 , 代入效率因子求解公式可得 : Cd0试验 = 0. 0201 e试验 = 0. 784 显然 , e < 1 ,试验结果是合理的。 而将对应的计算结果 :0 度攻角下升力系数约为 - 0. 009 ,阻力系数约为 0. 019 ;4 度攻角下升力系数 约为 0. 22 ,阻力系数约为 0. 024 代入可得 : Cd0计算 = 0. 0212 e计算 = 1. 108 e > 1 ,计算结果从这个角度来看似乎存在不合 理性。但是考虑到在 0 到 4 度之间都存在这种关系 , 仅用两个点数据在求解效率因子存在不准确性 ,因此 要得出准确的效率因子 ,应该用多个点进行最小二乘 法拟合求解。 通过计算 1 度、2 度的升阻力系数 ,再利用最小 二乘法和不同角度组合求解效率因子以及零升阻力 系数的结果参见表 2 (马赫数 0. 6) 。 由表 2 可见 ,除去利用 0、4 两个极端值以外 ,在 马赫数为 0. 6 情况下 ,效率因子 e 都是小于 1 的。这 说明 ,计算结果是合理的。 其他马赫数下计算结果所得的效率因子及零升 阻力系数与试验结果所得的效率因子及零升阻力系 数对比参见表 3。 由表 3 可得 ,试验与计算所得的效率因子均小于 1 ,结果是合理的。 7. 3 　有关阻力计算的分析 在常规的气动计算中 ,飞机阻力的计算向来是一 个公认的难题 ,我们不但计算出飞机的阻力 ,还可以 使用商品软件对阻力进行一些定性的分析 : 一般把物体在流体中运动时所受的阻力分为摩 擦阻力和压差阻力 ,前者是指表面切应力的贡献 ,它 与流体的粘性密不可分 ;后者则为表面压强的贡献 , 也与流体的粘性有着一定的关系。在我们的计算中 , 这两部分阻力是可以分别计算出的 : 例如 :当来流为 0. 8 马赫 ,攻角为 0 度时 ,压差阻 力系数约为 0. 01 ,摩擦阻力系数约为 0. 008。因此此 时的阻力系数为压差阻力系数与摩擦阻力系数之和 , 数值大小约为 0. 019。攻角为 4 度时 ,压差阻力系数 约为 0. 022 ,摩擦阻力系数约为 0. 008 ,对应总的阻力 348李向群等 :高速飞机的全机气动力数值分析 © 1995-2005 Tsinghua Tongfang Optical Disc Co., Ltd. All rights reserved. 表 2 　效率因子及零升阻力系数表 (Ma0. 6) (1) 试验/ 计算 角度选择 零升阻力系数 效率因子 计算 0、1、2、4 0. 01346 0. 9676 试验 0、1、2、4 0. 01177 0. 778 计算 1、2、4 0. 0184 0. 8578 试验 1、2、4 0. 0166 0. 773 计算 0、4 0. 0212 1. 108 试验 0、4 0. 0201 0. 784 计算 1、4 0. 0223 0. 816 试验 1、4 0. 021 0. 774 表 3 　效率因子及零升阻力系数表 (2) 马赫数 零升阻力系数 效率因子 试验 0. 6 0. 02012508 0. 784039 计算 0. 6 0. 02014191 0. 853262 试验 0. 8 0. 02543960 0. 817945 计算 0. 8 0. 023813197 0. 626448 试验 0. 95 0. 040103913 0. 642312 计算 0. 95 0. 045373566 0. 575673 试验 1. 05 0. 049332324 0. 673262 计算 1. 05 0. 039621716 0. 453083 试验 1. 2 0. 031340248 0. 599861 计算 1. 2 0. 029488035 0. 551091 试验 1. 5 0. 007074993 0. 41798 计算 1. 5 0. 009026342 0. 412675 0. 03 ;攻角为 8 度时 ,压差阻力系数约为 0. 061 ,摩擦 阻力系数约为 0. 0074 ,对应总的阻力系数约为 0. 068 ;攻角为 12 度时 ,压差阻力系数约为 0. 14 ,摩擦 阻力系数约为 0. 007 ,对应总的阻力系数约为 0. 147。 从中可以看出摩擦阻力与压差阻力的差别。零 攻角时压差阻力与摩擦阻力之间的差别并不是很大 , 量级基本相当。当攻角增大时压差阻力迅速增大而 摩擦阻力相对变化较小 ,此时压差阻力和摩擦阻力相 差可达一个量级。 对于大攻角 ,以及跨音速和超音速的情况 ,我们 所说的压差阻力还应包括尾涡阻力和激波阻力 ,这两 部分阻力与粘性无关 ,即使使用无粘流体进行计算也 应得到这些阻力结果。这里的尾涡阻力也就是指通 常的诱导阻力 ,主要是指低速无粘流有限翼展机翼与 无限翼展机翼 (即翼展阻力为零)相比 ,多了因为三维 效应而引起的这部分阻力。其流动机理是在正升力 情况下 ,下翼面的高压区气流绕翼梢流向上翼面低压 区 ,这就产生了翼尖涡 ;另外 ,上翼面存在自翼梢向翼 根的展向流动 ,下翼面存在翼根向翼梢的展向流动 , 这样上下翼面气流在机翼后缘流出而汇合时 ,这一上 一下的相反的展向流动形成漩涡 ,自机翼后缘向后拖 出 ,此后缘拖出的涡与翼梢涡组成了机翼后缘之后的 尾涡面。这一对机翼翼梢尾涡在机翼附近有一向下 的诱导速度 ,该诱导速度与机翼来流速度构成合速 448 水 　动 　力 　学 　研 　究 　与 　进 　展 　　　　　　　　　2004 年增刊 © 1995-2005 Tsinghua Tongfang Optical Disc Co., Ltd. All rights reserved. 度 ,使得流向机翼的速度改变了方向 ,有效迎角减少。 机翼上的实际升力是垂直于局部的来流速度的 ,此升 力向对于远前方来流就一个向后的分力 ,这就是诱导 阻力。 自各种阻力的大小和贡献在不同的工况下是不 同的 ,必须对具体问题做具体分析。通过我们的计 算 ,对于不同性质的流动 (亚音速 ,跨/ 超音速)我们可 以估计出它们的量级。 当马赫数为 0. 8 时 ,压差阻力中的激波阻力很 小 ,可以近似认为为 0。这样 ,我们可以以 0. 8 马赫 时的压差阻力为 标准 excel标准偏差excel标准偏差函数exl标准差函数国标检验抽样标准表免费下载红头文件格式标准下载 ,在跨音速段和超音速段中 ,认 为压差阻力与 0. 8 马赫时的压差阻力之间的差别由 激波阻力产生 ,以此来估计不同马赫数下激波阻力的 量级。 对跨音速段 ,4 度攻角下 ,马赫数为 0. 8 时 ,压差 阻力系数约为 0. 022。马赫数为 0. 95 时 ,压差阻力系 数约为 0. 039 ,两者之间的差值约为 0. 017。同理 ,当 马赫数为 1. 05 时 ,压差阻力系数约为 0. 056 ,与马赫 数为 0. 8 时压差阻力系数的差值约为 0. 034 ,因此可 以估计出 ,在跨音速段 ,激波阻力系数的量级大约为 10 - 2。 对超音速段 ,4 度攻角下 ,马赫数为 1. 2 时 ,压差 阻力系数约为 0. 056 ,与马赫数为 0. 8 时压差阻力系 数的差值约为 0. 034。当马赫数为 1. 5 时 ,压差阻力 系数约为 0. 05 ,与马赫数为 0. 8 时压差阻力系数的差 值约为 0. 028。因此可以估计出 ,在超音速段 ,激波 阻力系数的量级也大约为 10 - 2。 激波阻力引起的差别可以从跨音速的极曲线看 出 ,马赫数越大激波阻力越大。 7. 4 　对马赫数 1. 5 攻角 12 度的试验与计算的差别 分析 根据试验结果 ,当马赫数为 1. 5 ,攻角为 12 度时 , 升力系数约为 0. 78 ;而计算的结果与试验结果有较 大差异。根据计算结果 ,当马赫数为 1. 5 ,攻角为 12 度时 ,升力系数约为 0. 61。升力系数差别相当大 ,误 差为 - 21. 28 % ,是所有升力系数误差中最大的。 从马赫数 1. 5 的升力曲线来看 ,我们认为 ,计算 的结果是合理的 ,相反 ,试验的结果可能存在问题。 马赫数 1. 5 时升力曲线如图 5 所示。 在失速前 ,升力曲线是线性的。这一结论从其他 马赫数时试验和计算的升力曲线的趋势也可以得出。 从马赫数 1. 5 时的升力曲线趋势来看 ,在 0 度、4 度、8 度的时候 ,计算与试验结果基本是吻合的。并且基本 上位于同一直线上。但是当攻角为 12 时 ,试验值偏 离了由 0 度、4 度、8 度所拟合的直线 ,数值显然严重 偏大 ,并导致了整体拟合的直线斜率增大 ,截距减小。 此时计算的数值则与前三个角度符合很好 ,依然与前 三个角度的升力系数呈明显的线性关系。 因此 ,我们认为 ,马赫数 1. 5 攻角 12 度时的升力 系数计算结果是合理的 ,而试验的结果存在问题。 从马赫数 - 升力系数图中 ,我们也能得出同样的 结论。 当攻角分别为 4、8 度时 ,马赫数 - 升力系数如图 7 所示。 图 7 　a = 4 马赫数 - 升力系数图 (上) 及 a = 8 马赫数 - 升力系数图 (下) 从这两幅图中可以看出升力系数随马赫数的变 化而变化的趋势。在马赫数为 1 时升力系数达到最 大值 ,当马赫数大于 1 以后 ,升力系数随马赫数增大 迅速变小。总体来说 ,计算和试验的结果是符合的。 但是当攻角为 12 度时 ,马赫数 - 升力系数图则 出现了异常 ,异常表现在试验值上。攻角为 12 度马 赫数 - 升力系数如图 8 所示。 很显然 ,当马赫数大于 1 以后 ,计算值升力系数 的变化是符合规律的 ,而试验值升力系数的变化缓 慢 ,与攻角为 0、4 度时的马赫数 - 升力系数曲线不符 合。当马赫数增大到 1. 5 时试验值和计算值出现了 548李向群等 :高速飞机的全机气动力数值分析 © 1995-2005 Tsinghua Tongfang Optical Disc Co., Ltd. All rights reserved. 图 8 　a = 12 马赫数 - 升力系数图 图 9 　动力影响 Ma0. 6 升力曲线图 (上) 及 Ma0. 6 极曲线图 (下) 很大的差异。 经过反复对比分析我们认定问题在于试验值本 身 ,经过与试验单位的接触后答复是 :“对于马赫数 1. 5 , a = 12 的试验数据 ,我们查询了一下吹风情况。 当时吹风模型抖动厉害 ,测得的数据不准。同意你们 的分析。” 至此 ,马赫数 1. 5 ,a = 12 的试验数据与计算之间 存在的问题得到圆满解决。事实说明了计算结果的 准确性 ,并指出了试验中所出现的问题。 图 10 　动力影响 Ma0. 6 力矩曲线图 (上) 及 Ma0. 85 升力曲 线图 (下) 图 11 　动力影响Ma0. 85 极曲线图 (上)及Ma0. 85 力矩曲线 图 (下) 648 水 　动 　力 　学 　研 　究 　与 　进 　展 　　　　　　　　　2004 年增刊 © 1995-2005 Tsinghua Tongfang Optical Disc Co., Ltd. All rights reserved. 图 12 　动力影响Ma0. 85 极曲线图 (上)及Ma0. 85 力矩曲线 图 (下) 8 　带动力状态下计算结果分析 喷气式发动机对飞机气动特性的影响包括进气 和喷流两部份 ,对于进气的影响 ,由于其量级比较小 , 又缺乏相应的实验数据 ,很难给予准确的计算。相对 于进气影响 ,喷流对飞机气动特性的影响是显著的 , 喷流对气动的影响可以分为直接影响和间接影响 ,直 接影响包括发动机推力和法向力对升阻力和力矩的 贡献 ;间接影响主要是喷流对周围流场的诱导作用 , 使飞机的绕流发生一定的改变 ,从而影响飞机的气动 特性。 基于对动力影响的这些基本认识[4 ] ,我们首先采 用的方法是在进气道和尾喷管处 (在无动力影响的计 算中这些地方的边界条件是固壁条件)对方程提相应 的边界条件。 按照这种想法 ,进气道的边界条件 (如给定速度 或压力条件)很难给出 ,但它的影响不大 ,我们的做法 是根据向北航发动机专家的建议 ,考虑在飞行过程中 发动机进气有两个机制 :被动进气和主动吸气 ,因此 我们对不同马赫数假定一个发动机的压力损失因子 , 也就是将发动机的进气道入口向内延伸一段 ,在此给 定压力条件 ,但此压力较环境压力有一个降低 ,以此 作为边界条件 ,发动机进气入口的真实压力、速度和 质量流量由计算获得。 尾喷管喷流对气动特性的影响是比较显著的 ,计 算时需要尾喷管的每个马赫数的质量流量和相应的 总温、总压。计算过程中在尾喷管处采用喷射气流的 质量流量作为边界条件。在无动力影响的基础上增 加上述两个边界条件后采用同样的数值方法进行计 算 ,但是这样计算产生了一些问题 ,例如 :气道入口的 质量流量与尾喷管喷流质量流量的守恒问题。由于 边界条件的给法 ,内外流并没有耦合在一起 ,这样导 致进气道的质量流量与尾喷管质量流量可能不匹配。 另外 ,这样复杂的计算耗时是非常巨大的 ,计算一个 状态需要 7 天左右。对于如何更好的模拟发动机的 作用我们决定留作下一步研究 ,而先按实际的吹风模 型 ,在发动机部位设置为一直通管道 ,气流可以在管 道内自由流动 ,以此来模拟动力影响状态。这样改动 的主要依据是实际风洞吹风试验采用的模型大致也 是采用这种办法 ,在飞机内部添加管道 ,利用管道内 部的控制环控制流量的大小 ,以达到模拟飞机飞行时 的真实状态。 从计算的角度来看 ,这样做了更改后 ,原来复杂 的进气道入口边条件和尾喷质量流量边条件变为相 应的内点 ,而多出来的机身内几何形状虽然复杂 ,但 是物理条件为相对简单的固壁。这样求解方程时相 对说来更加准确 ,耗时相对减少。 部分计算结果曲线带动力与无动力对比如下 : 图 13 　动力影响 Ma1. 5 力矩曲线图 从计算最后的结果来看 ,“带动力影响 (直管)”的 升力系数相对于无动力状态下的升力系数有明显的 增大 ;在小攻角以及超音速段 ,带动力影响的阻力系 数相对于无动力状态下的阻力系数明显减小。在跨 音速的大攻角段 ,阻力系数并没有十分明显的增大或 减小的趋势。 748李向群等 :高速飞机的全机气动力数值分析 © 1995-2005 Tsinghua Tongfang Optical Disc Co., Ltd. All rights reserved. 9 　结束语 本文的中心内容是针对飞行器的空气动力学性 能进行数值分析。随着计算机硬件和数值模拟软件 的不断发展 ,进行这种大规模的数值模拟已经是完全 可行的了。经过算法验证以及与风洞试验对比表明 , 计算结果是可信的 ,可以实际应用的。它不但可以代 替大部风洞试验直接应用 ,有时还可以发现风洞试验 可能存在的问题[1 ,2 ] 。 参 　考 　文 　献 : [1 ] 　陈耀松. 创新与构思 - 力学小议之二 [J ] . 力学与实践 , 2002 ,24(3) :72. [2 ] 　陈耀松 ,陈沪东 . 力学小议[J ] . 力学与实践 ,2001 ,23 (4) : 74275. [3 ] 　JOHN D. ,ANDERSON JR. Computational Fluid Dynamics The Basics with Applications[M]. McGraw2Hill. 19951 [4 ] 凌茂芙 ,他山之石. 中俄、中乌民机气动技术研究 [M] . 北京 :航空工业出版社 ,1996. [5 ] 范洁川 . 国家标准《空气动力学 概念、量和符号》背景材 料和使用说明[M] .北京 :中国标准出版社 ,1999. [6 ] 吴子牛 . 计算流体力学基本原理[M] .北京 :科学出版社. 20011 [7 ] 吴望一 . 流体力学[M] .北京 :北京大学出版社 . 19821 [8 ] 　ICEM CFD Engineering , ICEM CFD Version4. 2. 2 User Guide [M]. U1S1A1ANSYS1INC1 , 2002. [9 ] 　FLUENT GAMBIT2. 0 User Guide[M]. U1S1A FUENT INC1 , 2001. [10 ] CD2ADAPCO , Star-CD3. 0 User Guide [M]. U 1S1A FLUENT INC1 ,19991848 水 　动 　力 　学 　研 　究 　与 　进 　展 　　　　　　　　　2004 年增刊 © 1995-2005 Tsinghua Tongfang Optical Disc Co., Ltd. All rights reserved.