（计算数学专业论文）地面上两串列立方体的干扰数值模拟.pdf

摘要 摘要 对从高层建筑结构抗风问题中抽象出来的方柱绕流问题进行数值模拟计算 有着非常重要的现实的意义。这一问题抽象出来的经典数学模型n a v i e r s t o k e s 方程的求解非常复杂，伴随着现代计算机技术的蓬勃发展，通过计算机技术来对 复杂数学模型进行数值模拟已成为一种发展趋势。 本文利用c f d 软件f l u e n t 对方柱绕流的流场进行了数值模拟计算，并对结 果进行了详细的分析。本文首先对绕流问题的基本公式、概念、数值模拟模型 和方法进行了了介绍，并着重介绍了k 一占计算模型和s i m p l e 数值方法。本文首 先针对不同网格尺寸采用k 一占模型的三种湍流模式对二维单方柱绕流在雷诺数 2 2 0 0 0 的情况下进行了计算，并将其结果( 包括升阻力系数、s t r o u h a l 数等) 与 实验结果和其他作者的计算结果进行了比较，以证明f l u e n t 软件对不可压缩绕 流流场计算的可行性及不同网格尺寸对模型模拟结果的影响。接着，本文研究 了地面上两串列立方体绕流的干扰问题，分别对间距比为s h = 1 、s h = 2 和s h = 4 三种情况下的流场进行了模拟，并与实验结果进行了比较。 关键词：k 一占模型，d e s 模型，s i m p l e 算法，升阻力系数，斯脱罗哈数 a b s t r a c t n u m e r i c a c o m p u t a t i o no ft h ef l o wa r o u n ds q u a r ec y l i n d e r sd r a w sf r o m i n t e r f e r e n c ee f f e c t so fw i n dl o a d so fac l u s t e ro ft a l lb u i l d i n g s i ti sd i f f i c u l tt os o l v e t h en a v i e r - s t o k e se q u a t i o n sw h i c he x p r e s st h i s p h e n o m e n o n i th a sb e e n a d e v e l o p i n gt e n d e n c yt on u m e r i c a l l ys i m u l a t et h i sc o m p l e xm o d e lw i t ht h eh e l po f c o m p u t e r , b e c a u s eo ft h ec o m p u t e rt e c h n o l o g yi se x p e r i e n c i n gf a s td e v e l o p m e n t i nt h i st h e s i s ，w es t u d yt h ei n c o m p r e s s i b l ef l o wf i e l da r o u n ds q u a r ec y l i n d e r s u s i n gac f d s o f t w a r ef l u e n ta n da n a l y z e dt h er e s u l ti nd e t a i l t h eb a s i cc o n c e p t s ， f o r m u l a s ，m o d e l sa n dn u m e r i c a lm e t h o d sa l ei n t r o d u c e df i r s t 1 1 l ek fm o d e la n dt h e s i m p l em e t h o df o rs o l v i n gi n c o m p r e s s i b l ef l o wf i e l da r ee m p h a s i z e do n w 色f i r s t l y s o l v et h ec l a s s i c a lc a s eo ft h ef l o wp a s tas i n g l es q u a r ec y l i n d e ri nd i f f e r e n t 鲥d sw i t h r e y n o l d sn u m b e ra t2 2 0 0 0t o t e s tt h er e l i a b i l i t yo ft h ef l u e n ts o f t w a r ea n dt o c o m p a r et h ei n f l u e n c eo fn u m e r i c a lr e s u l t sw i t hd i f f e r e n ts i z eo fm e s h a n dt h e nw e c o m p u t ei n t e r f e r e n c ep r o b l e mo ft w os u r f a c e m o u n t e dc u b e si nt a n d e m w ec a r r i e d o nt h es i m u l a t i o nt ot h ef l o wf i e l ds e p a r a t e l yt ot h r e ek i n do fs i t u a t i o nw h i c ha r e s h = 1 、s h = 2a n ds h = 4 ，a n dw ec o m p a r e dt h ec o m p u t e dr e s u l tw i t ht h ee x p e r i m e n t a l r e s t a l t k e yw o r d s ：k fm o d e l ，d e sm e t h o d ，s i m p l em e t h o d ，l i f ta n dd r a gc o e f f i c i e n t s ， s t r o u h a ln u m b e r t t 第一章绪论 1 1 引言 第1 章绪论 1 9 6 5 年英国渡桥电厂两排8 个冷却塔中的后3 个塔在风速并不太大的情况 下发生倒塌，这一事件激起了对建筑群干扰问题的研究。现代城市的发展产生 了密集的建筑群。由于相邻建筑物之间的流场相互干扰，受扰高层建筑的风荷 载大小及性态与其单独存在时相比会有很大的变化，在某些情况下可能会大大 超过其单独存在时的风力，引起抗风设计中的安全问题及建筑物居住者由于建 筑物振动产生的不舒适性问题。工程结构的抗风研究，例如作用于建筑物或结 构等钝体上的气动力预测，主要采用试验的方法在风洞中完成。随着计算机硬 件水平的飞速发展和c f d ( c o m p u t a t i o n a lf l u i dd y n a m i c s ) 技术的不断完善， 出现了与试验相对应的数值模拟方法。数值模拟比风洞实验成本低、所需周期 短、效率高；数值模拟不受模型尺度影响，可以进行全尺度的模拟，克服了风 洞实验中难以满足雷诺数相似的困难，还可以方便地变化各种参数，及早发现 问题。工程结构抗风的数值模拟由于其自身的特点，研究的对象多是钝体绕流 问题。钝体绕流问题广泛存在于桥梁工程、兵器工业工程和现代宇航工程等工 程科学领域中，具有很强的实际意义。例如，海洋平台、斜拉桥、原子能反应 堆中的换热器等都存在类似的流体与结构相互作用的情况。 本章首先阐述数值研究建筑群中建筑物之间相互气动干扰的现实意义，然 后回顾有关此问题的国内外研究现状，最后概要介绍本文的工作。 1 2 课题研究的意义 随着经济的发展和科学技术的进步，近二十年来，国内外建造了大量的重 大工程建筑结构。仅在上海陆家嘴地区，已建和拟建的4 0 0 米以上的结构有5 栋，2 0 0 米以上的超高层建筑有十多栋。计划中的意大利m e s s i n a 大桥的主跨达 3 0 0 0 多米，我国在建的苏通长江大桥是世界第一的斜拉桥，主跨达1 0 8 8 米。我 国奥运会及世博会的申办成功，将建造大量的大跨空间结构。此外，发达国家 第一章绪论 甚至提出了千米高度量级的“空中城市 的概念。强风作用下结构的风荷载和 响应是结构安全性和适用性的控制荷载之一。 在城市规划以及建筑设计中，特别是在城市建设日益向高度化发展的今天， 由于城市建筑物向高层、密集化发展，就会产生“城市街区峡谷”。在“城市 街区峡谷 中产生“城市急流”、“气流死区 ，由此引发一系列的环境问题： 比如在建成的街区内出现风口，或是局部风速过大，严重时将影响到行人以及 附近建筑物：由于风速和风向的改变，在有火灾等紧急情况发生时会出现烟道 效应，加速灾害的传递，增加灾害损失；局部街区可能出现流动迟滞现象，造 成严重的局部空气污染。此种种决定楼宇建设诱发的街区小气候问题将日益受 到人们的关注。有鉴于此，高层建筑群的相互干扰的研究具有很强的实际意义。 高层建筑问气流干扰研究的必要性主要表现在两个方面：一方面是高层建筑本 身的安全问题，另一方面是高层建筑对周围环境风场的影响。高层建筑群对环 境的危害大致是以下两个原因：( 1 ) 建筑物的钝边造成了流动的畸变，流线在拐 角处的密集导致风速的增大，一般认为建筑物迎风面拐角处1 5 米高度的最大 平均风速与上游同高度未受到干扰风速之比大于1 4 的时候，就可能会对行人 造成伤害甚至由于玻璃幕墙等覆盖物的脱落造成行人伤亡；( 2 ) 流动在拐角处从 建筑物表面分离，所形成的旋涡及其脱落，造成拐角风的阵性，形成风速突变， 使行人感觉不舒服。 对于高层建筑物的相互干扰的研究，目前大致有四种方法： ( 1 ) 实测；( 2 ) 风洞实验：( 3 ) 数值模拟；( 4 ) 理论分析 实测试验虽真实可靠，但耗时费钱且使用受限。风洞实验可靠、方便，但 当高层建筑物所处的地形复杂、风的方向速度变化无常时，大型的风洞实验被 需要，这导致耗资巨大，时间周期长。数值模拟比风洞实验成本低、所需周期 短、效率高；数值模拟不受模型尺度影响，可以进行全尺度的模拟，克服了风 洞实验中难以满足雷诺数相似的困难，还可以方便地变化各种参数，及早发现 问题。由于数值方法的居多优点，以及计算机的迅速发展、计算方法的改进， 数值模拟已成为研究高层建筑群的干扰问题的最有潜力的方法。本课题的研究 方法正采用了这种方法。 为了降低强风暴灾害所造成的损失，发达国家进一步加大了研究和开发的 投入。据悉，美国国会最近通过法案，启动了国家减风灾计划( n a t i o n a lw i n d h a z a r d sr e d u c t i o np r o g r a m - - n w h r p ) ，有关抗风研究和开发( t h ew i n dr e l a t e d 2 第一章绪论 r e s e a r c ha n dd e v e l o p m e n t ) 的费用由5 0 0 万美元年增加至2 3 0 0 万美元年。 日本政府也已启动了c o e ( c e n t e ro fe x c e l l e n c e ，2 0 0 3 - 2 0 0 7 ) 计划，增大投 入，开展相关研究。 可靠的湍流模式、先进的计算方法和高性能的计算机，这三者是实现高逼 真度的模拟仿真所需要解决的三个核心问题。美国能源部d o e 与其管辖下的三 大国家实验室在1 9 9 6 年9 月，提出“加速战略计算创新”计划，长远目标在于 通过1 0 0 万亿次级的计算机平台来实现三维、全物理过程、全系统规模的高分 辨率、高逼真度的科学模拟计算。这就意味着将改变过去以实验为主的科学研 究方法，而代之以科学模拟计算为主。创新的科学模拟将成为新世纪促进科技 发展的一项重要工具。本课题的研究是总结和验证部分前人所作的相关研究， 希望能够对我国科学模拟计算发展尽一点微薄之力。 1 3 国内外研究现状 数值研究建筑群中建筑物之间相互气动干扰的模型实质是数值计算钝体绕 流问题，其属于计算流体动力学中建筑钝体空气动力学的研究范畴。下面分别 介绍计算流体动力学、钝体空气动力学以及钝体绕流问题的国内外研究现状。 1 3 1 计算流体动力学的研究现状 c f d ( c 0 m p u t a t i o n a lf l u i dd y n a m i c s ) 即计算流体动力学是2 0 世纪6 0 年代 起伴随计算机技术迅速崛起的学科，经过半个世纪的迅猛发展，这门学科已经 相当成熟。各种c f d 通用性软件包陆续出现，成为商品化软件，为工业界广泛 接受，性能日趋完善，应用范围不断扩大，广泛应用于风工程、热能动力、航 空航天、机械、土木水力、环境化工等诸多领域，风工程领域是c f d 技术应用 的重要领域之一。 各种c f d 通用软件的数学模型的组成都是以n a v i e r - s t o c k e s 方程组与各种 湍流模型为主体，再加上多相流模型、燃烧与化学反应流模型、自由面模型以 及非牛顿流体模型等。大多数附加的模型是在主体方程组上补充一些附加源项、 附加输运方程与关系式。随着应用范围的不断扩大和新方法的出现，新的模型 也在增加。 第一章绪论 c f d 数值求解方法主要有有限差分法、有限单元法、边界元法以及有限分析 法。其中以有限差分法和有限单元法为主。有限差分法从微分方程出发，将计 算区域经过离散处理后，近似的用差分、差商来代替微分、微商，这样微分方 程和边界条件的求解就可以归纳为一个线形代数方程组的数值求解。空调室内 的流场、温度场和浓度场的数值模拟，长期以来几乎全部采用有限差分法。有 限单元法汲取了有限差分法中离散处理的内核，同时又继承了变分计算中选择 试探函数并对求解区域积分的合理方法。在有限单元法中，试探函数的定义和 积分范围不是整个区域，而是从区域中按实际需要划分的单元。经过比较发现 对于边界形状较规则的研究区域如矩形区域，二者模拟效果相同；而对于边界 形状较复杂的区域，有限单元法模拟效果更好。目前大多数的商用c f d 软件都 采用的是有限单元法。 此外，c f d 软件都配有网格生成( 前处理) 与流动显示( 后处理) 模块。网格 生成质量对计算精度与稳定性有很大的影响，因此网格生成能力的强弱是衡量 c f d 软件性能的一个重要因素。网格分为结构性网格和非结构性两大类。目前广 泛采用的是结构性网格。对于较复杂的求解域，构造结构性网格时要根据其拓 扑性质分成若干子域，各子域间采用分区对接或分区重叠技术来实现。非结构 性网格不受求解域的拓扑结构与边界形状限制，构造起来很方便，而且便于生 成自适应网格，能根据流场特征自动调整网格密度，对提高局部区域计算精度 十分有利。但是非结构性网格所需内存量和计算工作量都比结构性网格大很多。 因此，两者结合的复合型网格是网格生成技术的发展方向。目前，f l u e n t 软件 已具有这种功能。 下面具体针对目前市场占有率最高的f l u e n t 软件阐述其应用和主要特点。 c f d 通用软件包的出现与商业化，对c f d 技术在工程中应用的推广起了巨大 的促进作用。1 9 9 8 年，全球市场占有率最高的c f d 软件f l u e n t 正式进入中 国市场，为目前c f d 主流商业软件，其市场占有率达4 0 左右。 f l u e n t 软件设计基于“c f d 计算机软件群的概念 ，针对每一种流动的物 理问题的特点，采用适合于它的数值解法以在计算速度、稳定性和精度等方面 达到最优。 f l u e n t 软件的结构由前处理、求解器及后处理三大模块组成。f l u e n t 软件 中采用g a m b i t 作为专用的前处理软件，使网格可以有多种形状。对二维流动可 以生成三角形和矩形网格；对于三维流动，可以生成四面体、六面体、三角柱 4 第一章绪论 和金字塔等网格：结合具体计算，还可以生成混合网格，其自适应功能，能对 网格进行细分或粗化，或生成不连续网格、可变网格和滑动网格。f l u e n t 软件 采用的二阶上风格式是b a r t h 与j e s p e r s e n 针对非结构网格提出的多维梯度重 构法( m u l t i _ d i m e n s i o n a lg r a d i e n tr e c o n s t r u c t i o n ) ，后来进一步发展，采 用最d , - - 乘法估算梯度，能较好的处理畸变网格的计算。f l u e n t 率先采用非结 构网格使其在技术上处于领先。 f l u e n t 软件的核心部分是n a v i e r s t o c k e s 方程组的求解模块。用压力校正 法作为低速不可压流动的计算方法，包括s i m p l e 、s i m p l e r 、s i m p l e c 、p i s o 等。 采用有限体积法离散方程，其计算精度和稳定性都优于传统编程中使用的有限 差分法。离散格式为对流项二阶迎风插值格式叫u i c k 格式( q u a d r a t i cu p w i n d i n t e r p o l a t i o nf o rc o n v e c t i o nk i n e t i c ss c h e m e ) ，其数值耗散较低，精度 高且构造简单。而对可压缩流动采用耦合法，即连续性方程、动量方程、能量 方程联立求解。 湍流模型是包括f l u e n t 软件在内的c f d 软件的主要组成部分。f l u e n t 软件 配有各种层次的湍流模型，包括有代数模型、一方程模型、二方程模型、湍应 力模型、大涡模拟等。应用最广泛的二方程模型是k 一模型，软件中收录有标 准k 一模型及其几种修正模型。 f l u e n t 软件的后处理模块具有三维显示功能来展现各种流动特性，并能以 动画功能演示非定常过程，从而以直观的形式展示模拟效果，便于进一步的分 析。 1 。3 ，2 钝体空气动力学研究现状 空气绕过钝体时的风场和绕过流线体时有明显的不同，这是因为存在分离 流和剪切层的非定常振动。结构风工程要面对复杂的钝体空气动力学问题，还 要考虑来流紊流所引起的迎风钝体表面的压力脉动以及结构风致振动。 对钝体空气动力学的研究始于1 9 3 3 年，b a i l e y 在英国国家物理实验室( n p l ) 的风洞中率先进行了一座铁路车库的模型试验，并和现场实测的压力系数作了 对比，发现了来流边界层的粗糙高度z 。和房屋高度h 之比的重要影响。1 9 5 8 年m j e n s 据此建议这一无量纲比值作为模型试验的相似常数，称为m j e n s e n 定律。6 0 年代a g d a v e n p o r t 和b j v i c k e r y 发展了边界层风洞技术，进 第一章绪论 行了简单钝体的绕流试验。1 9 7 3 年i csg a r t s h o r e 第一次系统地揭示了不同 比例的矩形钝体在来流紊流风场中的气动现象。他发现随着紊流的加强将减少 边缘分离边界层的曲率半径，使分离流提前与物体再附，从而改变压力分布以 及相应的升力曲线斜率，并且影响振动钝体的驰振现象。此后，bel e e ( 1 9 7 5 ) 和w h m e l b o u r n e ( 1 9 7 5 ) 用流迹显示试验深入地研究了气流在钝体前缘的分离 和再附的规律，进一步明确了紊流强度对流态的影响。 进入8 0 年代以后，钝体空气动力学的研究逐渐和工程应用相结合。主要有 以下几个方面： 1 高层建筑的幕墙动压 由于发生了玻璃幕墙的强风作用下的破坏，世界各国的风载规范试图列入 有关幕墙峰值压力系数的条文。风洞试验显示在建筑物的转角处，由于气流的 分离出现明显的边缘部连续性，最大峰值负压系数可以从常数的3 - 4 增大到8 ，因 而必须引起特别的注意。 2 高层建筑的横风向相应 高层建筑作为钝体绕流将发生边缘处的旋涡脱落并激发横风向的振动。来 流紊流也会对旋涡脱落过程产生重要影响，主要是改变气流的再附，从而影响 响应的大小。与方形断面的建筑相比，圆形断面或削角的方形断面建筑的横风 向响应要小得多，同时驰振的稳定性也能得到改善。 3 桥梁的涡振和抖振响应 形状较钝的桥梁断面会发生涡激振动。在紊流风场中，脉动风压引起的抖 振力将激起桥面的随机振动响应。从流迹显示中可以观察到紊流度对流态( 分 离和再附) 的影响，并改变跨度方向的相关性。 4 大跨悬臂雨篷 大型体育场建筑的兴建提出了大跨悬臂看台雨篷的风载和风振问题。目前 主要通过节段模型和全气弹模型的风洞试验来观察流态，测定压力分布和振动 响应。 由于计算机的迅速发展，计算流体动力学( c f d ) 也从7 0 年代由航空工程逐 渐引入了风工程流域。最初从均匀流场中的圆柱体绕流开始，8 0 年代起，为了 研究紊流场中钝体的绕流，建立了各种紊流模型。利用这些模型进行t - ：维和 三维矩形断面的固定钝体和振动钝体的流场和压力场数值分析，取得了和实验 较一致的结果。 6 第一章绪论 综上所述，由于钝体空气动力学在理论上的困难，目前还处于较初始的阶 段，研究对象还限于简单钝体。对于实际土木工程中复杂外形的振动钝体问题， 目前主要依靠风洞试验来提供设计所需的风载和风振响应值。 1 3 3 钝体绕流问题研究现状 钝体绕流的研究已经有一百多年的历史了。然而钝体绕流中的许多流动现 象和物理机制仍未得到透彻的认识，一直是国内外流体力学学者的研究课题。 钝体绕流中研究最多的是圆柱绕流问题。方柱绕流问题最近亦有很多学者进行 了研究。其它不规则类型的绕流问题由于涉及的参数复杂，研究相对困难。本 文研究了方柱绕流问题。对于方柱绕流问题的国内外研究状况作了如下了解： 国外，k a t s u y ae d a m o t o 和m u t s u t ok a w a h a r a 1 3 用有限元方法对串列放 置的两方柱周围的流场的研究；k i mt 和f l y n n m r 2 用离散涡方法对多物 体周围流场的计算；k a w a i ，h i r o m a s a 和f u j i n a m i ，k i y o s h i 3 对两串列方柱周 围流场的计算；在国内，陈素琴 4 ，5 ，6 对方柱绕流问题作了比较全面的研究， 包括两方柱串列放置、并列放置的流场的计算，以及对低雷诺数和高雷诺数两 种情况下分别交错放置的两方柱干扰的数值计算，并给出了串列、并列以及交 错放置时流场的干扰机理。河海大学的王兴勇 7 用l a t t i c eb o l t z m a n n 方法模 拟了方柱位于流场的边壁和位于流场的中央两种情况下的方柱绕流。y o s h i d a 8 和穆维贤 9 用有限单元法计算了方柱位于流场中央时绕流形成的速度分布、压 力分布以及施加扰动后产生的周期旋涡。哈尔滨工业大学的徐枫 1 0 研究了二 维刚体方柱的涡激振动和非定常绕流，利用计算流体力学软件f l u e n t 6 0 ，引入 雷诺应力模型求解不可压粘性流体n a v i e r s t o k e s 方程，方柱和流场的耦合作用 通过f l u e n t 的动网格技术来实现，进行雷诺数从7 0 0 0 到4 2 0 0 0 范围内的模拟 得到静止方柱的升阻力系数和涡脱频率，其中雷诺数为2 2 0 0 0 结果与试验一致： 同时得到方柱振动的位移和速度时程。哈尔滨工业大学的魏英杰 1 1 应用二阶 全展开e t g 有限元离散格式与大涡模拟相结合的方法，对间距比为1 5 情况下的 并列双方柱绕流进行了数值模拟，给出了两方柱的阻力系数、升力系数以及两 对称点上流向速度随时间的变化历程，并采用谱分析的方法研究了对称边界条 件下并列双方柱绕流的频谱对称性问题结果表明：对称边界条件下，双方柱绕 流运动参量的时域过程虽然是不对称的，但频域过程却是对称的。 7 第一章绪论 1 4 本文的研究内容 有关钝体绕流问题的文章比较多，然而国内外数值研究不同模型钝体绕流 问题的比较的文章相对不多。本人收集了一些相关的论文资料，如刘宇 1 2 等 人采用湍流模型利用g a l e r k i n 有限元方法计算了不同雷诺数下方柱绕流的数值 模拟，得出了不同雷诺数下卡门涡街脱落形态；文章 1 3 、 1 4 、 1 5 3 等采用 r n g 湍流模型对钝体绕流流场的不稳定、非定常进行了数值模拟；m u r a l a m i 1 6 、 1 7 等人使用标准模型对钝体绕流进行了研究，还有 1 8 一 2 6 3 ，但它们都属于 一个模型的研究验证，没有给出不同模型的比较，也没有考虑网格大小对计算 结果的影响。本文主要是对前人提出的相关理论在不同模型及不同网格尺寸下 进行数值模拟验证以及比较各自间的差异性。 同时干扰问题的研究现正成为热点，如文章 2 7 - 3 0 3 都是对干扰问题的研 究。尽管国内外有不少人研究多柱体相互干扰问题，但他们大多是对无限长柱体 进行研究，对地面上柱体的干扰情况研究的较少。r o b e r tj m a r t i n u z z i 和b r i a n h a v e l 3 1 3 在2 0 0 4 年用风洞试验方法对地面上两串列方柱绕流流场进行了实验研 究，但目前用数值模拟方法对此问题研究的还比较少。本文主要是用数值模拟的 方法利用r n g k 一占模型研究地面上两串列方柱体周围的流场特征。 8 第二章方柱绕流问题的计算模型和数值方法 第2 章方柱绕流问题的计算模型和数值方法 2 1 引言 方柱绕流问题的流场特征是方柱周围的边界层和剪切层，这两个特征几乎 完全决定了方柱在流场中的气动力行为。方柱绕流问题中，流场中各处的涡的 尺度并不都是一样的大小，越靠近方柱，涡的尺寸越小。如果在计算中采用能 分辨最小涡尺度的均匀网格，将会带来巨大的计算量。本文利用非均匀网格， 采用网格伸展的网格设计技巧，即当对绕着非流线形体的流场进行分析时，离 非线性体越远，网格分解设计得越粗糙。 本研究采用流体力学软件f l u i n t 对所提出的数学模型使用七一占模型的三种 湍流模式分别进行数值模拟。其中用到的数值方法为s i m p l e 方法。根据文章 3 2 _ 3 6 ，本章对各种模型分别作了介绍。 2 2n a yie r - s t o c k e s 方程 n a v i e r - s t o k e s 方程描述了粘性不可压缩流体的运动。自然界中大量的流体 模型，诸如具有热传导效应的流体动力学模型，磁流体动力学模型，海洋动力学 模型，以及描述像血液流动等管道流的数学模型，其主部均为n a v i e r - s t o k e s 方程。因此，它有很强的物理背景，是当今非线性科学研究中的重点和热点问 题。它的研究现已成为非线性偏微分方程、数值分析和动力系统研究的推动力 量。俄罗斯科学院和著名数学家斯梅尔均把该模型作为本世纪需要解决的重大 问题之一。 在风工程中，流场的描述方程为不可压n a v i e r s t o k e s 方程( 以后称为n s 方程) ，其形式为 d i v = 0( 2 1 ) 旦竺+ 厅露：一 + l 露(229rad g r a d p ) 百+ ” “一 + 一r e 。， 其中历，p ，和r e 分别表示速度向量、压力、时间和雷诺数，上述方程已用特征长 度l 、特征速度u 、流体的密度p 和运动粘性系数v 将各量无量纲化，r e ：一u l 。 9 第二章方柱绕流问题的计算模型和数值方法 2 3k 一占模型 2 3 1 标准k 一占模型的优缺点 最简单的完整湍流模型是两个方程的模型。在f l u e n t 中，标准k s 模型 自从被l a u n d e r a n ds p a l d i n g 提出之后，就变成- r 程流场计算中主要的工具了。适 用范围广、经济、合理的精度，这就是为什么它在工业流场和热交换模拟中有 如此广泛的应用了。它是个半经验的公式，是从实验现象中总结出来的。主要 是基于湍流动能和扩散率。k 方程是个精确方程，e 方程是个由经验公式导出的 方程。标准k 一占模型的方程如下： ( 胱，) = 0 ( 2 3 ) 鲁+ 毒c 叩沪一吉考+ 毒c + 以考+ 等， c 2 4 ， 昙c 肚，+ 毒c 砌沪考l c + 箦，考j + 瓯一胪 c 2 5 ， 昙c ，+ 毒c 胪，2 毒 c + 箦，考 + q 。妻g g 。p i 6 2 c 2 6 ， 式中湍流逐度鸬由f 式确定： 朋= 孵等 ( 2 7 ) c 。是常量。瓯表示由层流速度梯度引起的湍动能生成项。可以定义为 g 。= 一俐- - v - - ，：- , 瓦o u 7 ( 2 8 ) 为了计算方便，b o u s s i n e s q 假设：g i = 从s 2 ，s 为表面张力系数，定义为 s = q 嚣j s o 一 ( 2 9 ) 其中 1 0 第二章方柱绕流问题的计算模型和数值方法 西= 摇+ 等 标准后一s 模型方程中的常数一般取q ：1 4 4 ， 些常量是从试验中得来的。 ( 2 1 0 ) c 2 。= 1 9 2 ，c 。= 0 0 9 ，仃，= 1 3 逐 k s 模型假定流场完全是湍流，分子之间的粘性可以忽略。标准k 一占模型 因而只对完全是湍流的流场有效，并非对任何流动都适用。由于人们已经知道 了k s 模型适用的范围，因此人们对它加以改造，出现了r n g k s 模型和带旋 流修正k 一占模型 2 3 2r n g k 一占模型( r e n o r m a ii z a t i o ng r o u p ) 的优缺点 它是从瞬态n - s 方程中推导出来的，使用了一种叫“r e n o r m a l i z a t i o n g r o u p ( 重整化群) 的数学方法。r n g 的基本思想是通过在任意空间尺度上的一 系列连续的变换，对原本十分复杂的系统或过程实现粗分辨率的描述。r n g k 一占 模型的k 方程与标准k s 模型方程中的k 方程一样，只是e 方程有差别，方 程如下： 考+ 以毒= 毒心毒) + q 妻q g 占一足 ( 2 1 1 ) + 以瓦5 瓦吒面+ 乏q q 占气一足 ( 2 1 1 ) 瓯为湍动能生成项，与标准k - s 模型中的g 定义相同。c l 。= 1 4 2 ，c 2 占= 1 6 8 f 与标准k - 6 方程相比增加了一项r 。， c ，, p r l 3 ( 1 一旦) ， 耻斋等 亿 其中 例孝( s = 厨，s = 争警+ 等】) ( 2 1 3 ) r o 和是模型常数，r o = 4 3 8 = 0 0 1 2 由于r 。项的存在，当刁 玎。时，r 占为负，结果在r a p i d l ys t r a i n e d 流中，r n g 模型产生的湍流粘度要低于标准k 一占模型。因而，r n g 模型相比于标 第二章方柱绕流问题的计算模型和数值方法 准七一s 模型对瞬变流和流线弯曲的影响能作出更好的反应，这也可以解释r n g 模型在某类流动中有很好的表现。 r n g k 一5 模型来源于严格的统计技术。它和标准膏一占模型很相似，但是有以 下改进： ( 1 )r n g 模型在方程中加了一项r 。，有效的改善了精度； ( 2 )标准七一占模型是一种高雷诺数的模型，r n g k 一占理论提供了一个考虑低 雷诺数流动粘性的解析公式。这些公式的效用依靠正确的对待近壁区域。 而且，r n g k s 模型方程中的常数并非用经验方法确定，而是由理论推导出 来的精确值。这些特点使得r n g k s 模型比标准霓一g 模型在更广泛的流动中有更 高的可信度和精度。 2 。3 。3r e ai iz a bie k s 模型的优缺点 它的方程与标准模型相比k 方程是一样的( 常量除外) ，也只有方程有差 别，方程如下； 昙c 脚+ 丢c 盹，2 毒卜+ 箦，考j + 皖一 q m ， 扣，+ 毒c 刚= 毒卜针鹏熹 方程( 2 1 5 ) 中的v 表示运动粘性系数，方程( 2 1 4 ) 中的q 表示由于平 均速度梯度引起的湍动能生成项，它与标准七一s 模型中的q 定义相同。模型中 常量为c 2 = 1 9 ，吒= 1 2 。 像其他的模型一样，漩涡粘度由下式计算： 氓= 口cn 竺 r e a l i z a b l e k f 模型与标准模型r n g k 一占模型的区别在于c 。不再是常 量了，它由下式计算： q 2 不赫 q 1 6 1 2 第二章方柱绕流问题的计算模型和数值方法 u 兰& 岛+ 孬蓼孬移 ( 2 1 7 ) q l ，= q 一2 s u k r a k q 驴= 一啄吨 这里孬盯是在柱坐标下的带有角速度的q 层流旋度，模型常量4 为： 彳。= 4 0 4 ，彳。= 拓c 。s 矽，= l 气c o s - i ( 廊) ，w = s u s j k s k , s ，露= 属瓦， 西= 丢c 善+ 可以看出，是层流应变和旋度的函数。模型常量c 2 、吼和仃。已经为某种 规范流做过优化。模型常量是：c l f = 1 4 4 ，q = 1 9 ，o k = 1 o ，盯。= 1 2 r e a l i z a b l e k s 模型是近期才出现的，i 扫s h i h 提出，比起标准k 一模型来 有两个主要的不同点：一是带旋流修正的k f 模型为湍流粘性增加了一个公式； 二是耗散率e 输运方程不同，为耗散率增加了新的传输方程，这个方程来源于 一个为层流速度波动而作的精确方程，因而被认为更适于表达能谱的传输。 r e a l i z a b l e k s 模型在范围广泛的流动现象中得到验证。 2 4 其他模型 2 4 1s p a i a r t a ii m a r a s 模型 在湍流模型中利用b o u s s i n e s q 逼近，中心问题是怎样计算漩涡粘度。这个 模型被s p a l a r ta n da l l m a r a s 提出，利用解决因湍流动粘滞率而修改的数量方 程。此模型的偏微分方程是： 昙( 聊( p v = g v + 吉秀批讷争( 劬母s 汜1 8 ) 模型的变量中矿是湍流动粘滞率除了近壁区域，这里的g ，是湍流粘度生成 1 3 第二章方柱绕流问题的计算模型和数值方法 的，f v 是被湍流粘度消去，发生在近壁区域。& 是用户定义的，注意到湍流动 能在s p a l a r t a ll m a r a s 没有被计算，但估计雷诺压力时没有被考虑。 2 4 2 标准和s s tk c o 模型 s s t 和标准模型的不同之处是 ( 1 ) 从边界层内部的标准k 一国模型到边界层外部的高雷诺数的k 一占模 型的逐渐转变 ( 2 ) 考虑到湍流剪应力的影响修改了湍流粘性公式 标准k 一彩模型 标准模型是一种经验模型，是基于湍流能量方程和扩散速率方程。 昙( 肚) + 毒( 砌扣毒( r 。考) + 瓯一k + 瓯 ( 2 1 9 ) 瓦0 ( c o 件毒( 删，) _ 考( r 。考) 十瓯一匕+ 瓯 ( 2 2 0 ) 在方程中，g k 是由层流速度梯度而产生的湍流动能。g 。是由缈方程产生 的。k 和f 。表明了后和国的扩散率。k 和y o 是由于扩散产生的湍流。瓯和咒是 用户定义的。其中 l ：+ 丛， o k l ：+ 丝 6 m 这里吼和盯脚是后、缈方程的湍流能量的普朗特数。湍流粘度鸬： ( 2 2 1 ) 鸬爿譬c 口哦c 蒜券跳= 盖) 汜2 2 ， s s t k 一缈模型 s s t 模型更适合对流减压区的计算，另外它还考虑了正交发散项从而使方 程在近壁面和远壁面都适合。它的七方程和标准k 一缈模型是一样的，只有力方 1 4 第二章方柱绕流问题的计算模型和数值方法 程有差别，方程中加了一项正交发散项见： c o ) + 毒( 删。考( r oa o ) + a 口o x j - 匕蝇怛 ( 2 2 3 ) 2 4 3 雷诺应力模型( r s m ) 为了反应r e y n o l d s 应力在各个不同方向之间的运输关系，许多科学工作者 放弃了b o u s s i n e s q 涡粘性假设，设法建立了关于r e y n o l d s 应力的微分输运关 系式，并由此得出r e y n o l d s 应力输运方程模型( r s m ) 。r s m 完全抛弃了湍流粘 性的假设，直接建立以u l u 。为因变量的偏微分方程，因此与两方程模型相比，在 理论上具有更好的通用性和精度。但这套模型求解雷诺应力的所有分量所满足 的微分方程，同时还要求解k 、方程，从而使求解的模型方程数大大增加，对 计算机容量和计算费用的要求也大大增加。目前，r s m 的应用还不如k f 模型 广泛。方程为： 警+ 瓦警= 而妻一弧差+ 伤+ 岛一白 组2 4 ) 尹+ i = 飞吩蔷一吩魄葛+ 伤+ 一白 ( 2 2 4 ) 式中pi i 是雷诺应力再分配项，dl | 是雷诺应力扩散项，l l 是雷诺应力耗 散项，各相关项的封闭模式参阅文献 3 2 。 2 4 4l e s 模型 湍流大涡数值模拟( l e s ) 是有别于直接数值模拟和雷诺平均模式的一种数 值预测湍流的方法。文章 3 7 一 4 2 都是采用l e s 模型进行数值模拟的。l e s 是 把脉动运动在内的湍流瞬时运动量通过某种滤波方法分解成大尺度运动和小尺 度运动两部分。大尺度要通过数值求解运动微分方程直接计算出来，小尺度对 大尺度运动的影响将在运动方程中表现为类似于雷诺应力一样的应力项，该应 力称为亚格子雷诺应力，它们将通过建立模型来模拟。实现大涡数值模拟，首 先要把小尺度脉动过滤掉，然后再导出大尺度运动的控制方程和小尺度运动的 封闭方程。l e s 的控制方程为： 要+ 擎：一三要+ v 黑 ( 2 2 5 ) o t o x i pa x i o x j o x j 1 5 第二章方柱绕流问题的计算模型和数值方法 互：o ( 2 2 6 ) 苏 令“l 甜j = 甜l 甜j + ( 甜l 甜一us “j ) ，则式( 2 2 5 ) 可写作 ( 2 2 7 ) 式中弓= - ( u ，豁- - s t ，甜) 称为亚格子应力，其为不封闭，须构造亚格子应力的封 闭模式。 2 。4 。5 分离涡模型( d e s ) 虽然大涡模拟效果较好，但需要大量的计算资源。这主要是因为它需在时 间与空间上精确求解含能湍流涡，尤其在近墙区域，计算量因所需求解的尺度 加速变小而急剧增加。由此产生了分离涡方法( d e s ) 。d e s 在近墙区域使用雷诺 平均，从而降低对网格的要求。在其他区域则采用l e s ，从而确保模拟精度。d e s 模型中使用r e a l i z a b l e k s 模型作为雷诺平均算法，d y n a m i ck i n e t i ce n e r g y 亚格子模型作为l e s 算法。k 方程的扩散项改为 ，：丝竺 q m i n ( k 3 他s ，c 如) 其中。为0 6 1 ，a 为局部单元的最大边长。 2 5sim pie 算法 s i m p l e 算法：s e m i i m p l i c i tm e t h o df o rp r e s s u r el i n k e de q u a t i o n ，于 1 9 7 2 年，由p a t a n k a r 和s p a l d i n g 最早提出，它现在实际上已经成为许多工程 流动、传热以及反应体系的数值模拟的最重要的方法许多商、 k c f d 软件，女h c f x 与f l u e n t 其核心也都基于s i m p l e 算法。下面以直角坐标下二维层流定常运动 为例介绍该方法。 二维层流定常运动的控制方程为： 1 6 挈器 万一钆p 巫 甄百 第二章方柱绕流问题的计算模型和数值方法 l 眦( p “) + 专d 圳) = 丢( 罢) + 专( 孝) 一卺+ s 。 昙p 圳) + 专( ) = 昙( 罢) + 品( 多) 一号+ s ， 昙d “) 十昙v ) ：o o xd v ( 2 2 8 ) ( 2 2 9 ) ( 2 3 0 ) 这组方程有两个新的特点： ( 1 )动量方程的对流项包含非线性项。 ( 2 )三个方程耦合在一起，没有明显的关于压力的输运方程。 s i m p l e 方法通过迭代技术解决了上面两个问题。 我们使用两套交错网格线，一套为标以大写字母“i 一1 ，i ，i + l ，“j 一1 ，j ， j + l 的实线网格，一套为标以小写字母“i - 1 ，i ，i + 1 ”，“j 一1 ，j ，j + l 的虚线 网格，若标量节点定义在两实线交点，则下标记为( i ，j ) ：x 方向速度u 定义在x 方向虚线，y 方向实线相交点上，标记为( i ，j ) ，y 方向速度v 定义在x 方向实 线，y 方向虚线相交点上，标记为( i ，j ) 通过在标量节点对速度的离散化，我们得到节点( i ，j ) 的离散动量方程： a i , j u = z a 柏材柚+ ( p ，- l ，一p i , j ) 4 ，- ，+ 包， ( 2 3 1 ) a l , j v ，= z a 曲1 ，柏+ (