CFD仿真分析（共8篇）
CFD仿真分析 篇1
引言
近年来,能源问题的加剧和日益严格的排放法规促使涡轮增压技术在内燃机领域得到了广泛的应用。涡轮增压能使内燃机在动力性、经济性、比质量和排放等性能指标上有很大的提高和改善,它是实现现代内燃机高效率、低污染的重要技术措施之一。目前,涡轮增压器正向高增压、高转速、高效率以及更宽广的运行范围发展,这对涡轮增压器的性能提出了更高的要求,而涡轮增压器性能提高的关键是提高压气机和涡轮的气动效率。压气机和涡轮叶片的内部流动很复杂,理论难以精确分析其流场,CFD仿真为研究其内部流场分布和各种损失的原因提供了重要的手段。涡轮和压气机的仿真计算和内部流场有一定的相似性,本文以某压气机为例,利用FIRE仿真软件计算压气机的效率和流量,并分析其内部流场。这能极大的缩短压气机的设计周期,对于改善压气机性能以及设计新压气机具有重要意义。
1、计算模型及网格划分
1.1 建立几何模型
采用逆向工程,利用三维激光扫描技术测量得到了叶轮表面和压气机内表面的点云数据,并使用CAD软件Solidworks建立了叶轮和蜗壳内部流动区域几何模型(图1、图2)。其中,压气机内部流动区域的入口和出口段进行了加长,这样的边界流动比较均匀,有助于边界条件的设定和计算收敛。
1.2 网格划分
压气机和涡轮均分成三部分建模:入口部分,转子部分和出口部分。转子部分与入口部分和出口部分之间的耦合采用的是FIRE中的多重参考系模型MultiReferenceFrame,MRF)。转子部分的网格在计算时保持静止,在惯性坐标系中以作用的哥氏力和离心力进行定常计算,而入口和出口区域则是在惯性坐标系里进行定常计算。在3个部分的交界面处交换惯性坐标系下的流体参数,保证了交界面的连续性,达到了用定常计算来研究非定常计算的目的。计算区域的网格使用FIRE软件中的Fame Advanced Hybrid模块生成,包括六面体网格和四面体网格。入口部分最大网格尺寸设为2 mm,最小网格尺寸设为1 mm。转子部分最大网格尺寸设为0.8 mm,最小网格尺寸设为0.1 mm,从而使得叶轮的叶尖间隙处至少有4层网格,保证了计算精度。出口部分最大网格尺寸设为2mm,最小网格尺寸设为0.5 mm。最后,入口部分生成55548个网格,转子部分生成594777个网格,出口部分生成92383个网格。计算区域的总网格数为742708个,如图3所示。
1.3 边界条件设置
进口设置为标准大气的条件。取总压为101000 Pa,温度为30℃,该条件跟试验很接近,计算结果可以和试验数据直接对比,比较方便。出口设置静压,计算不同的压比只改变出口的静压。计算时从低压比不断向高压比推进,高压比的计算利用低压比的结果,这样不仅提高了计算的稳定性,还大大加快了计算速率。
2、计算结果与分析
流量和效率是压气机最重要的参数。计算时设定一定的转速,固定进口边界条件,不断改变出口静压,得到不同压比下的压气机的流量、效率。然后改变转速,最终得到压气机的脉谱图。由于计算的点较多,而且不同转速下,流动有一定的相似。在低压比时,压气机易堵塞,效率偏低;在中等压比时效率较高;在高压比时容易喘振。本文选择叶片转速为9万转时的流场进行分析,分别取低压比(压比1.1)、中压比(压比1.24)和高压比(压比1.33)三个点与实验结果对比,如表1所示。该压气机的叶轮半径较小,9万转的转速是压气机工作时较低的转速,所以压比偏低,压比1.33就接近喘振边界了。
由上表1看出,仿真得到的效率偏高,而在中低压比时流量偏小,接近喘振时流量偏大。初步分析有几点原因:
1.网格数不够,受设备和计算时间的限制,网格数不能过多;网格质量缺陷,蜗壳有一部分和压气机叶轮形状都比较复杂,局部网格质量不高,导致计算精度下降。2.在几何模型的建模过程中,由于经过了三维激光扫描和CAD软件的处理,最终得到的几何模型在形状上不可能与真实的计算区域完全一致。3.边界条件的差异,实际中壁面不绝热;更重要的是壁面摩擦的状况跟实际是有差距的。4.实验数据本身也存在误差。
三个状态下的壁面附近的压力场分布如图4所示。(下文为叙述方便,压比为1.11状态称为A,压比为1.25的状态称为B,压比为1.33的状态称为C。)
从三个状态下的静压分布可以看出,压气机进气由于弯折,进入叶轮前压力不均匀。刚要进入叶片区域前有一个压力降低的过程,这是由于该段的流通截面有一定的渐缩,压力能向动能转化造成的。这段加速过程有利于减小进气的不均匀。在叶片区域气流压力增压较多,进入无叶扩压器区域后,压力接近最大值,对于B、C,在蜗壳内也有一定的压力增加。但是对于A,压力场似乎比较反常。最大压力出现在蜗壳的较小的流通截面处,说明在该处气流不通畅,流速低,动能较多的转为压力能。这说明该点过于偏离设计流量,蜗壳的流通效果很不好。
进一步分析A效率不高的原因。整个计算域的马赫数没有超过1,说明没有产生激波而造成压气机堵塞。图5、图6为分别为某长叶片吸力面和压力面的静压分布。吸力面入口叶尖有一高压区,而压力面则相反有一低压区,说明气流在进口的叶尖处有强烈的撞击并造成了较强的气流分离,这里会有较大能量损失。从图7过压气机轴线的切面看出,该处有一块高压区,迫使气流向外壁流动,而外壁急剧转弯,气流产生强烈分离,这样有效的流通截面就减小了,流动阻力大大增加。该扩压器的宽度仅为2.5 mm,不太适合大流量,但是较窄的扩压器适应更高的压比。
对于A和B,它们的流量都要比实验值低,主要原因是叶尖漏气,如图8所示,叶片顶端聚压气机外罩的间隙虽然只有0.4 mm,但是由于本压气机很小,到叶片出口处叶高仅有2.5 mm,会有相当一部分气体从间隙流过。
对于C,比较反常的是它的流量比实验值偏大。从图9看出,在高压比下,流体的回流明显,而且回流首先在从叶轮出的壁面开始向进口扩展,在涡轮叶片进口处卷起大团漩涡。该状态即将达到压气机的喘振边界,气流的分离以及入流与回流的掺混都很厉害。此时,湍流模型的精度会大大降低,计算精度受到挑战。再者,仿真的进气是很稳定的,与实验有些波动的进气有差距,也会导致接近喘振时流量偏高。
仿真得到的效率都偏高,也意味着能量损失小,出口温度偏低。整个流动过程,叶轮中的损失应该占大部分,在叶轮处温升的应该很可观的,但是以状态A为例,从叶片进口到叶片出口,平均温升仅6℃左右,说明在这一区域仿真的精度还不够高。从图10看出,高温区主要在叶尖吸力面一侧,由于该处存在较强的流量分离和能量损失。而沿叶片轴线方向上,温升较小,说明摩擦损失不多。还有选择旋转区域时,没有对静止的外壁作特殊处理。如果外壁跟着叶片一起旋转,虽然外壁对气流扰动不大,但是气流与外壁的相对速度减小了,从而减小了摩擦损失。
3、结论
离心式压气机基于AVL FIRE的三维流动CFD分析可以得到压气机内部的流场信息,并可以从中找出模型内存在较大流动损失区域和造成损失的原因,从而进行有针对性的优化。同时也可以得到压气机的流量特性图,这对压气机的前期评估和与发动机的匹配预测起到非常重要的指导作用。
压气机在低压比时效率较低,但对于本文来说,并不是超音速导致的堵塞,因为计算转速较低,流场还没达到音速。而是在较高流量时压气机的局部损失较大。在高压比时,压气机的倒流很厉害,计算精度难保证。中压比时由于碰撞损失和气流分离减少,效率和计算的精度都较高。
压气机的内部流动比较复杂,网格要求很高,加上逆压梯度,对计算的精度带来巨大的挑战。下一步将继续提高计算精度,将尝试不同的湍流模型,更加细密的网格和更加合理的边界条件。
参考文献
[1]周龙保.内燃机学[M].机械工业出版社,2005.
[2]肖昕,李云清.车用涡轮增压器蜗壳内三维流场模拟分析[J].汽车技术,2011.
[3]车用小型离心式压气机结构参数优化及流量特性的仿真计算[D]西安交通大学硕士学位论文,2010.
[4]最新增压器设计、生产新工艺新技术与故障诊断排除及质量检验标准规范实用手册[M].北方工业出版社,2007.
CFD仿真分析 篇2
(南京航空航天大学机械结构力学与控制国家重点实验室， 江苏 南京 210016)
引 言
直升机乘坐舒适性要求的提高，对减小直升机振动水平提出了更高的要求。直升机的振动主要来自于旋翼系统及其与机身之间的相互耦合。旋翼/机身耦合系统的振动响应分析，是进行减振处理，降低全机振动水平的基础。
已有许多国外学者对直升机振动响应分析进行了研究。早期的工作大多是分别求出旋翼和机身的阻抗，令两者的阻抗在桨毂处匹配，进而求得机身的振动响应[1～3]。Warmbrodt和Friedmann最先建立了旋翼/机身耦合系统的非线性气动弹性响应分析模型，考虑了机身的6个刚体自由度[4]。Fledel通过有限元方法，对旋翼/机身耦合系统的振动响应进行了研究，考虑了旋翼与机体的运动和惯性载荷之间的相互耦合[5]。之后，Friedmann和Chopra等分别建立了旋翼/弹性机身的非线性气动弹性方程，采用梁单元建立了机身的三维有限元模型，对机身典型位置的振动响应进行了研究[6～9]。
准确的机身结构模型和旋翼流场分析，是直升机振动响应分析的关键。上述研究均采用准定常或半经验的二维气动力公式，机身则采用刚体模型或基于梁单元的简化模型，其分析精度还有待进一步提高。
近年来，随着计算机技术和计算流体力学(CFD)的迅速发展，基于非定常Euler/N-S方程的CFD方法，成为目前最有效的旋翼流场分析方法[10～12]。这为旋翼/机身耦合系统分析精度的提高，提供了新的思路。目前CFD方法已经用于求解单独旋翼的气动弹性问题，但对于旋翼/机身耦合系统的响应分析，由于问题难度更大，还未见到相关文献。
本文将CFD方法引入到旋翼/机身耦合系统的振动响应分析中，通过动态嵌套网格与动网格技术来实现桨叶网格的刚体运动和弹性变形，利用CSD软件建立精细的机身三维有限元模型，然后通过带配平的松耦合迭代方法求解系统响应。该方法能准确模拟机身结构和旋翼流场的特性，具有较高的分析精度和工程适用性。以某型两片桨叶跷跷板式直升机为算例，对悬停及定常前飞状态下，机身典型位置的振动响应进行了分析，将计算结果与实验值进行了比较，并考察了机身弹性对系统响应分析的影响。
1 计算方法及理论
1.1 旋翼/机身耦合动力学方程
旋翼为两片桨叶跷跷板式构型，桨叶采用15自由度的中等变形梁单元[13]。理论气动力采用准定常气动力模型和Drees线性入流模型。旋翼运动方程与机身方程通过桨毂中心点耦合在一起，旋翼载荷通过桨毂传到机身上，而机身的运动也会影响到旋翼载荷。
在常规的旋翼气动弹性响应分析中，仅考虑桨毂力常值部分对全机配平方程的影响，在定常前飞时，则假定机身固定不动。但桨毂力的谐波部分会使机身产生刚体运动，进而与旋翼运动发生耦合。因此，在机身振动水平分析中，还必须考虑机身的刚体运动自由度。
系统的动力学方程通过哈密顿变分原理得到
(1)
式中 δU为弹性能变分，δT为动能变分，δW为外力虚功。对这些能量有贡献的是桨叶和机身，各能量的变分可表示为：
(2)
(3)
(4)
下标b，F分别表示桨叶和机身，Nb表示桨叶片数。
此次分析中，通过MSC/Nastran软件建立机身的有限元模型，并根据模态试验结果，对有限元模型进行修正，以保证模型有足够的分析精度。然后将机身方程与旋翼方程组装起来，得到旋翼/机身耦合系统的动力学方程。
机身有限元模型如图1所示，机身主承力结构采用梁单元和壳单元，其他设备质量则采用集中质量单元分布到对应位置处。机身前10阶弹性模态的计算值和实验结果见表1，可以看出，有限元模型的固有模态计算值与实验结果吻合良好。
图1 机身有限元模型
表1 机身弹性模态
为减小计算量，对桨叶运动方程和机身弹性运动方程进行了模态缩聚。为避免出现矩阵奇异，将机身刚体运动自由度分开来，放到方程右端。惯性坐标系下，得到旋翼-刚体/弹性机身耦合系统的非线性动力学方程为：
(5)
式中M，C，K为系统的质量、阻尼和刚度矩阵；FL和FNL表示线性和非线性力向量；p是模态缩聚后的广义坐标。下标b，G，r和e分别表示与桨叶、桨毂万向铰、机身刚体运动和机身弹性模态相对应的自由度。旋翼/机身耦合系统方程是一个非线性，周期系数常微分方程组，本文采用时间有限元方法进行求解。机身刚体运动则由配平方程来求解，并通过与式(5)迭代，得到耦合方程的解。
1.2 CFD分析模型
旋翼流场采用非定常Euler/N-S方程进行求解，并通过动网格与动态嵌套网格方法来实现桨叶的弹性变形和刚体运动。为节省计算时间，本文控制方程采用非定常Euler方程。计算中采用双时间步进行时间推进，每个旋转周期进行360个物理时间步的计算，每个物理时间步内进行20个伪时间步的迭代计算，对流通量采用Roe格式计算。
图2(a)是旋翼系统的嵌套网格装配示意图，其中背景网格由7个网格块组成，总网格单元数为1 749 000。单片桨叶网格如图2(b)所示，桨叶采用C-H型结构网格，由8个子网格块组成，共171 900个网格单元。两组桨叶网格和背景网格组成整个嵌套网格系统，网格总单元数为2 092 800。
图2 旋翼系统网格图
1.3 分析流程
直升机的配平计算是旋翼气动弹性分析中不可缺少地重要环节，本文采用自由飞行条件下的6自由度平衡方程，来计算直升机的操纵输入(包括总距θ75和周期变距θ1c,θ1s)、飞行姿态角(俯仰角αs和滚转角φs)和尾桨总距θT。
旋翼-刚体/弹性机身耦合系统的响应，采用带配平的松耦合迭代方法进行求解，CFD和CSD数据每旋转一周交换一次，具体分析步骤如下：
5．重复步骤3和4直到系统的响应收敛。此时机身的振动响应即为旋翼/机身耦合系统的非线性振动响应。
2 结果与讨论
本文以某两片桨叶跷跷板式直升机为算例，直升机的具体参数如表2所示。每片桨叶划分为6个具有15自由度的中等变形梁单元，各个单元的截面参数见表3。其中，li/R是无因次化的单元长度，mi/m0是无因次化的桨叶线密度，该算例中取参考线密度m0=2.63 kg/m。EIz/m0Ω2R4，EIy/m0Ω2R4和GJ/m0Ω2R4分别为无因次化的摆振、挥舞和扭转刚度。
表2 直升机参数
表3 桨叶单元参数
2.1 悬停状态
在悬停状态下，采用本文所提出的旋翼-刚体/弹性机身耦合系统分析方法，对旋翼轴根部和机身前端点(如图1所示)的振动响应进行分析，并将计算结果与实验值进行比较。
图3是旋翼轴根部振动加速度功率谱密度曲线的计算值和实验值的对比。可以看出，对于2倍通过频率(2/rev)和4倍通过频率(4/rev)的振动响应，计算值与实验值吻合良好。试验中有一定的通过频率响应成分，而计算中没有，这主要是由旋翼各片桨叶的气动和质量不平衡造成的，而在计算中没有考虑这些因素。此外，实验中还有较大的52～54 Hz频率响应，这是由发动机(额定转速3 200 r/min)振动造成的。旋翼轴根部靠近发动机，因此其受发动机振动的影响较大，从图3(c)中可以看出，当主旋翼所引起的响应较小时，发动机振动和噪声的影响将明显增大。
图4是机身前端点处的加速度功率谱密度曲线，由于此点远离发动机，因此其振动受发动机的影响明显减小。
图3 旋翼轴根部振动自功率谱密度曲线
图4 机身前端点振动自功率谱密度曲线
图5 机体纵向振动随前飞速度的变化
2.2 前飞状态
在定常前飞状态下，对机身前端点和旋翼轴根部的振动响应进行计算，研究前飞速度对机身振动水平的影响，并比较考虑和不考虑机身弹性运动对分析结果的影响。
图5～7是机体的纵向、横向和垂向加速度响应幅值随前飞速度的变化曲线。从图中可以看出，当直升机由悬停转入前飞状态时，机体振动水平会随着前飞速度的增加而迅速增大，并在前进比μ=0.05左右时到达一个局部峰值；当前进比μ在0.05～0.2范围内时，机体振动水平则基本保持不变；当μ>0.2后，机体的振动水平又会随着前飞速度的增加而迅速增大。
图6 机体横向振动随前飞速度的变化
图7 机体垂向振动随前飞速度的变化
对比图5～7可以看出，机体的纵横向振动响应明显高于垂向振动响应，这说明旋翼的横向力扰动相对较大。在前进比μ<0.2时，各个自由度均以2/rev的频率响应分量为主，4/rev频率的响应相对较小。而当μ>0.2时，4/rev频率的振动响应所占比重则逐渐增大。
对比仅考虑机身刚体运动和综合考虑刚体/弹性运动两种分析模型，可以看出，对旋翼轴根部位置，两者的分析结果误差较小；而对机身前端点，两者的误差则较大。这是因为旋翼轴弹性相对较小，因此旋翼轴根部振动主要以刚体运动为主；而机身前端点远离旋翼，且结构弹性较大，因而其振动中弹性运动所占的分量比较大。由此可见，在进行直升机振动水平分析时，机身的弹性影响是不能忽略的。
3 结 论
本文提出了一种直升机旋翼-刚体/弹性机身耦合系统振动响应分析方法，综合考虑了桨叶的几何非线性，机身的弹性和刚体运动等。桨叶采用15自由度中等变形梁单元，通过基于动网格和动态嵌套网格技术的非定常Euler/N-S方程，与准定常气动力模型相结合的方法来求解弹性旋翼流场，利用CSD软件建立精细的机身三维有限元模型，最后用松耦合迭代方法对系统的响应进行分析。该方法能准确模拟机身结构动力学和旋翼流场的特性，具有很好的分析精度和工程适用性，可广泛应用于直升机振动水平分析和减振设计。
以某小型直升机为例，通过MSC/Patran软件建立精细的机身三维有限元模型，并根据模态试验对模型进行修正。计算了悬停状态下机身典型位置的振动响应，计算结果与实验值吻合良好。研究了前飞速度对机身振动水平的影响，结果表明，在小速度和大速度前飞时，机身振动水平随前飞速度的增加而增大，在中等速度段，机体振动水平基本保持不变。
分析中采用了仅考虑机身刚体运动和综合考虑刚体/弹性运动两种模型，结果表明，在旋翼轴根部，机身振动以刚体运动为主，而在机身前端处的机体振动响应中，结构弹性的影响非常大。因此，在进行直升机振动水平分析时，必须考虑机身结构弹性的影响。
参考文献：
[1] Gerstenberger W, Wood E R. Analysis of helicopter aeroelastic characteristics in high speed flight[J]. AIAA Journal, 1963, 1(10): 2 366—2 381.
[2] Staley J A, Sciarra J J. Coupled rotor/airframe vibration prediction methods[R]. NASA SP-352, 1974, 81—90.
[3] Hsu T K, Peters D A. Coupled rotor/airframe vibration analysis by a combined harmonic-balance impedance-matching method [J]. Journal of the American Helicopter Society, 1982, 27(1): 25—34.
[4] Warmbrodt Friedmann. Formulation of coupled rotor/fuselage equations of motion [J]. Vertica, 1979, 3: 245—271.
[5] Fledel. Coupled rotor/airframe vibration analysis[D]. Aerospace Engineering, University of Maryland, 1989.
[6] Thiem Chiu, Peretz P Friedmann. A coupled helicopter rotor/fuselage aeroelastic response model for ACSR[A]. AIAA Conference[C]. New Orleans,USA.1995, 574—600.
[7] Richard C Cribbs, Peretz P Friedmann, Thiem Chiu. Coupled helicopter rotor/flexible fuselage aeroelastic model for control of structural response[J]. AIAA Journal, 2000, 38(10):1 777—1 788.
[8] Senthilvel Vellaichamy, Inderjit Chopra. Aeroelastic response of helicopters with flexible fuselage modeling[A]. AIAA Conference[C]. Dallas,USA.1992： 2 015—2 026.
[9] Hyeonsoo Yeo, Inderjit Chopra. Coupled rotor/fuselage vibration analysis for teetering rotor and test data comparison[J]. Journal of Aircraft, 2001, 38(1):111—121.
[10] Klausdieter Pahlke, Berend G， van der Wall. Chimera simulations of multibladed rotors in high-speed forward flight with weak fluid-structure-coupling [J]. Aerospace Science and Technology, 2005, 9: 379—389.
[11] Mark Potsdam, Hyeonsoo Yeo, Wayne Johnson. Rotor airloads prediction using loose aerodynamic/structural coupling [J]. Journal of Aircraft, 2006，43(3)：732—742.
CFD仿真分析 篇3
通过三维CFD软件FIRE计算某型汽油机的铸铁排气歧管内部的气体流动情况。该计算为瞬态计算, 在排气歧管的进出口设置质量流量边界, 其周期性的流量及温度边界数据来源于一维热力学软件BOOST的计算结果。
为了评估排气歧管在各工况下的气流流动情况, 我们选择了三个典型工况点进行评估, 即额定工况 (4 850 r/min) 、倒拖工况 (4 850 r/min) 以及怠速工况 (700 r/min) 。对于每个工况点, 都计算四个发动机循环, 以保证计算的收敛性。
1 CFD分析方法
1.1 计算模型
一般来说, 排气歧管的计算都只是分析歧管本身。但是由于对一个CFD分析来说, 边界条件的准确性直接决定了整个计算的精确性, 而如果只是排气歧管本身作为计算模型的话, 我们很难在歧管入口给出准确的气流分布情况。因此, 为了保证进入排气歧管的气流与实际情况相符, 在排气歧管的计算模型前加上排气道的模型, 则整个计算域的模模型型见见图图11。。
1.2 热力学边界条件
一维热力学计算的主要目的就是给三维计算提供精确的边界条件, 边界条件的准确与否直接决定了三维计算的准确性, 因此一维热力学计算是整个排气歧管CFD分析中十分重要的一步。
图2所示的是在BOOST中搭建的一维热力学计算模型, 其中测点CY1, CY2, CY3, CY4表示四个缸排气道的进口, 测点MP12, MP13, MP14, MP15分别是排气歧管四个缸支管的进口, 而测点MP6则表示排气歧管的出口。
1.3 计算网格的生成
我们采用FIRE自带的半自动网格生成器FAME ADVANCED HYBRID进行CFD网格的划分。该网格划分工具能够生成以六面体为主的网格 (六面体约占85%以上) , 加上少量的金字塔网格、四面体网格。
为了满足计算的精度, 同时又要尽量减少计算的时间, 对于不同的工况划分不同数目的网格模型。其中对于气体流量较高的工况 (额定工况以及倒拖工况) 划分更精细的网格, 网格数约为36万;而对于气体流量较小的工况 (怠速工况) , 则采用更粗的网格模型, 网格数仅为7.3万。
1.4 计算设置描述
由于排气歧管主要用于将缸内燃烧的废气排出气缸, 其受到的热负荷较大, 因此, 排气歧管主要考查其承受高热负荷工作的能力, 所以排气歧管的分析需要考虑整个发动机循环。而整个发动机循环过程中, 排气歧管内部的气流每时每刻都不同, 因此, 排气歧管的计算是瞬态的。
对应于瞬态计算模式, 其余计算设置为:设定管内空气流动为可压缩粘性湍流流动, 不考虑损失项, 将空气视为理想气体。湍流模型选择kzeta-f方程, 因为此方程计算精度更高。离散格式选择迎风离散, 而压力与速度耦合算法选择SILMPLE。
2 计算结果分析
2.1 额定工况
图3为排气歧管换热系数分布云图。从图中可以看到, 局部换热系数主要在613~1 590 W/ (m2·K) 之间。其中换热系数最高的部位分布在排气歧管三、四缸之后与排气总管结合的部位。经过加权平均, 额定工况下排气歧管近壁表面的平均换热系数为1 020 W/ (m2·K) 定工况来说, 较高的换热系数有利于排气歧管的散热, 而此排气歧管的换热系数高达1 020 W/ (m2·K) , 是非常理想的数值。
图4为排气歧管温度分布云图。从图中可以看到, 近壁面气流的温度在757~930℃之间。其中温度最高的部位出现在三四缸支管的结合部分 (图中圆圈所示) 。而排气歧管近壁面的平均温度为875℃, 这是铸铁排气歧管可以承受的范围。
2.2 倒拖工况
图5、图6分别是倒拖工况下排气歧管换热系数分布云图以及温度分布云图。由于倒拖工况即发动机不点火, 通过电机拖动发动机达到额定工况转速 (4 850 r/min) 的一种评价工况。因此, 尽管其转速与额定工况相当, 但是由于发动机没有点火, 所以其热负荷明显比额定工况低得多。
从图中可以看到, 倒拖工况下排气歧管的局部换热系数主要在243~421 W/ (m2·K) 之间, 而其平均换热系数为333 W/ (m2·K) 。近壁面气流的温度在270~335℃之间, 其平均温度为316℃。
2.3 怠速工况
图7、图8分别是怠速工况下排气歧管换热系数分布云图以及温度分布云图。从图中可以看到, 怠速工况下排气歧管的局部换热系数仅在35~73W/ (m2·K) 之间, 而其平均换热系数只有53 W/ (m2·K) , 明显比额定工况和倒拖工况小。但在怠速工况下, 排气歧管近壁面气流的温度在281~376℃之间, 其平均温度为343℃。虽然比起额定工况小很多, 但是却高于倒拖工况。这是因为虽然怠速工况转速比倒拖工况小很多, 但是怠速工况下发动机仍然点火燃烧, 而倒拖工况是由电机拖动的, 发动机并没有点火。尽管其转速比怠速工况大得多, 其热负荷却比怠速工况低。
3 结论
通过对排气歧管的内流场进行数值分析, 分别得到了额定工况、倒拖工况以及怠速工况, 歧管内部的换热系数以及温度分布。其中额定工况下换热系数和温度均最高, 但是数值都在一个合理的范围内。而倒拖工况的换热系数高于怠速工况, 其温度略低于怠速工况。
这种应用CFD技术对排气歧管进行研究的方法, 不仅能够获得歧管内部的宏观流动特性, 而且能够获得内部流场的大量微观信息, 为排气歧管的优化设计和改进提供有利依据。
摘要：首先通过一维热力学软件BOOST计算排气歧管的边界条件, 然后通过三维CFD软件FIRE计算排气歧管内部的气体流动情况。该计算为瞬态计算, 选择了额定工况、倒拖工况以及怠速工况三个典型工况点。通过计算, 得到了排气歧管的换热系数分布以及温度分布。
关键词：CFD,排气歧管,换热系数,温度
参考文献
[1]王福军.计算流体动力学分析[M].北京:清华大学出版社, 2004.
[2]郭立新, 韩颖, 惠涵, 等.CFD-FE耦合计算分析某汽油机排气歧管热负荷[J].现代车用动力, 2009, (2) .
CFD仿真分析 篇4
当前进口清扫车的价格非常高昂, 有些型号不能适合我国的路面情况。所以在清扫车的吸尘原理及结构方面有待深入分析以提高吸尘效果。气力输送系统主要由风机、吸嘴、过渡吸筒、沉降室等部分组成。
目前国内对气力输送整体进行仿真分析的工作还很少, 本研究针对全系统仿真做一些尝试性工作。本文利用CFD软件, 对洗扫车的气力输送部分进行仿真分析, 包括风机、吸嘴、沉降室的建模、网格划分及后处理。分析风机内部流场, 以直观的方式表示出气力输送系统的速度场、压力场。
2 风机的仿真分析
2.1 建模
本文采用Solidworks对风机流域建立实体模型, 风机可以分为风叶和蜗壳两部分。其中为了分析方便, 将风叶末端放大为一个圆面, 便于边界设定, 风叶中心也做类似处理。图1为叶轮造型图, 图2为风机整体造型。
2.2 网格划分
将风机总体分为蜗壳和风叶, 并分别命名为“static”和“moving”;选择风机入口面为“velocity-in”。出口面命名为“pressure out”, 风叶壁面命名为“impeller”, 风叶大端端面命名为“interface2-2”, 风叶小端端面命名为“interface2-1”。与之对应, 蜗壳与风叶大端端面重合的面命名为“interface1-2”, 另一个面命名为“interface1-1”。这样就把模型的重要边界独立出来, 方便后续FLUENT设置边界条件。风机实物图如图3所示, 具体网格划分如图4所示。
2.3 风机内部流场分析
为了便于观察风机内部流场, 创建两个面:沿着风叶中间的面和吸风口中间平面。现逐一对各平面结果进行分析, 如图5、6所示针对离心风机进行整体分析。风机静压从吸风口到叶轮中心逐渐降低, 经叶轮旋转带动逐渐向蜗壳外部增高。随着气体流出蜗壳, 流动损失, 静压减小。并且静压最小中心靠近蜗壳最大处。
动压不能用压力表直接读出。如图7、8, 其值是和流动速度有关, 从风叶中心到风叶边缘速度逐渐增高, 动压增大, 风叶边缘达到最大, 蜗壳内动压由于气体相互作用而减小, 在出风口处逐渐趋于稳定, 在截面变化的地方也会受到影响。
3 吸嘴及沉降室的仿真分析
3.1 建立吸嘴及沉降室几何模型
吸嘴是一个封闭的异型吸气罩, 由于风机高速抽取垃圾箱内的空气, 使之形成较强的负压, 从而在吸嘴处形成高速补充气流, 以气力输送的方式, 把吸嘴四周的尘土吸入垃圾箱, 实现吸扫功能。图9为为吸嘴结构图;图10为吸嘴、沉降室整体模型。
3.2 网格划分
本文应用三种划分网格方法并行进行生成网格, 一种是四面体 (Tetrahedrons) 网格划分方法, 即Patch Conforming Method。一种是对吸嘴面进行细化也就是局部网格控制, 即Refinement命令。最后一种是Body Sizing网格划分方法, 本文综合以上三种方法进行网格划分, 如图10为沉降室实物图, 图11所示为吸嘴到沉降室整体的网格模型, 单元尺寸为30 mm。
3.3 吸嘴及沉降室的流场分析
气体进入吸尘口并由连接管进入沉降室的过程中, 满足动量守恒和质量守恒。由于空气在吸尘口中的流动属于高雷诺数的湍流流动, 而且不属于强旋流, 本文控制方程选择标准κ-ε方程模型。
整体吸嘴料箱模型计算中参数设置如下:
速度入口:34 m/s (根据风机模拟结果)
出口条件:outflow
收敛残差:0.001
湍流模型:Standardκ-ε模型
求解器:Pressure-Based Solver
运行环境:在入口处设置参考压力101 325 Pa (忽略重力)
当风机达到工作转速后, 吸嘴底面由于内部的负压以及和地面距离10mm的封闭空间, 在大气压的作用下, 空气从外部通过缝隙进入到吸嘴内部, 从而带动灰尘颗粒进入到沉降室, 进而达到清扫目的。
图12所示的静压图表示压力从吸嘴外部逐渐减小, 外部最大, 到吸嘴过渡接口处压力减小到一稳定值。图13显示的是速度矢量图, 由图可知吸筒部分气流速度相对较高。在沉降室内部, 气流流经遮罩板时, 由于遮罩板的形状使得气流向沉降室后部流动, 形成涡流, 再被吸入到出气口。图14的速度流线图很好地描述了气体的流动过程。
3.4 新建截面的的计算分析
为了直观显示吸嘴内部压力变化及速度情况, 新建界面进行描述。图15显示的是截面为X=1.45m的速度云图, 速度从吸嘴外部到吸筒变化较大, 吸筒转折处速度达到最大, 经过遮罩板缓冲后速度逐渐减小, 形成涡流后气流到出风口速度逐渐增大, 从而形成循环过程。
3结语
利用CFD对清扫车气力输送系统进行仿真分析得到以下结论:
(1) 通过仿真分析可以发现结构上的缺陷, 提出改进措施, 从而达到节能、提高吸尘效率的效果。
(2) 通过FLUENT对风机流场数值模拟结果的后处理, 可以显示风机内部各种流动参数的云图、矢量图和流体流动状态, 方便进行设计和参数调整。
(3) 通过仿真分析可以直观了解从吸嘴到料箱理论上的压力及速度变化情况, 以便于设计人员进行适当的改进。
参考文献
[1]杨华, 等.基于CFD紊流计算的离心泵叶型优化设计[J].扬州大学学报 (自然科学版) .2007 (03) .
[2]席光, 等.叶轮机械气动优化设计中的近似模型方法及其应用[J].西安交通大学学报.2007 (02) .
[3]杨策, 胡良军, 何永玲.叶片前缘后掠对离心压气机性能的影响[J].北京理工大学学报.2006 (11) .
[4]陈波, 高学林, 袁新.基于NURBS的叶片全三维气动优化设计[J].工程热物理学报.2006 (05) .
[5]杨春朝, 等.基于流场模拟的真空清扫车吸尘口数设计[J].中南大学学报 (自然科学版) .2012 (09) .
[6]姜兆文, 成凯, 耿宇明.吸扫式扫路车吸嘴流场性能分析[J].专用汽车.2012 (06) .
CFD仿真分析 篇5
关键词：火炬塔,CFD,数值仿真,安全评估
0 引言
随着经济的发展,我国开始承办越来越多的国际大型运动会,如2008年北京奥运会、2010年广州亚运会、2011年深圳世界大学生运动会等,即将举办的还有东亚运动会、世界青年运动会等。大型综合运动会属于人员密集活动,安全问题至关重要,需要进行专门评估,火炬塔的安全评估就是其中之一。大型运动会的火炬塔一般体型巨大,火炬燃烧所产生的热量也较高,如2008年北京奥运会,主火炬塔高30m,最大燃气流量约7000m3/h,燃烧功率约70MW。火炬塔需要在运动会期间连续燃烧,一般在火炬塔附近会存在一些设备或结构,在燃烧期间是否会对周围设备或结构造成影响,是火炬塔运行期间必须要考虑的问题。如火炬塔对其周围钢结构的影响、对电子显示屏的影响等等。目前,在国内举行的大型运动会,如奥运会、亚运会等,均采用本文提供的方法进行了安全评估。而在国外,未见有开展该方面研究的报道。本文提出了大型运动会火炬塔安全评估方法,以某大型运动会主火炬塔为例,采用CFD仿真方法,对火炬塔进行了安全分析。
1 分析方法与仿真工具
1.1 分析方法
火炬塔燃烧的热量主要以辐射和对流的形式向周围传递,周围物体吸收热量后,温度会逐步升高,物体的温升应在其安全范围之内。评价火炬塔燃烧对周围物体的影响最直接的方法就是进行试验研究。但由于客观原因,试验研究往往难以实现。另一种方法就是数值分析法。因火炬燃烧及其热量传递可归结为计算流体力学和计算传热学的问题,因此可借助CFD软件进行分析[1,2,3]。根据问题特点,提出如下分析步骤:
(1)确定安全判定标准:根据分析对象的特点,建立安全标准,如温度标准、应力标准等,本分析主要涉及到物体的温升,因此主要建立温度判定标准;
(2)确定计算场景:根据分析内容,确定相应计算场景,一般考虑最不利场景,如火炬塔燃料流量最大的燃烧工况;
(3)CFD仿真计算:利用CFD前处理器仿真模型,设置控制方程及边界条件,划分网格,选择收敛方法和收敛准则,进行模拟计算;
(4)结果分析和安全判定:通过后处理程序对计算结果进行处理,将结果和判定标准进行对比,以确定模拟对象的安全性;
(5)提出防护措施:若存在不安全的情况,则根据仿真结果,提出相应的防护措施,之后再进行仿真计算,直至满足安全要求[4]。
1.2 仿真工具
目前常用的CFD软件有PHOENICS、FLUENT、CFX、STAR-CD等[5,6,7,8]。理论上,各软件均可完成分析,但从易用性和准确性方面考虑,本文采用PHOENICS。
PHOENICS由英国CHAM 公司开发,是世界上第一套计算流体力学与计算传热学商用软件。可用于求解可压缩或不可压缩、单相或多相流体的稳态或非稳态流动,已广泛应用于航空航天、船舶、汽车、安全、暖通空调、环境、能源动力、化工等各个领域。安全分析是PHOENICS软件一个十分重要的应用领域,其在此方面的应用于已有十五年之久的历史。其可用于通风/排烟分析、可燃、毒性气体的泄漏分析、污染物扩散分析、火灾模拟等方面[9]。
1.3 控制方程
PHOENICS所求解的控制方程的通用形式如式(1)所示[10]:
其中:φ为求解变量,可表示温度、速度、压力、焓等,ρ为密度,t为时间,
PHOENICS采用有限体积法对控制方程进行离散处理,基本离散格式如式(2)所示:
其中,a为控制方程离散后的系数,P为控制体节点,角标W,E,N,S,H,L,T分别代表离散节点周围的八个节点及时间。在PHOENICS中,边界条件是按源项来处理的。离散格式可选择一阶迎风、混合格式、QUICK格式等,默认情况下采用QUICK格式。PHOENICS采用交错网格法进行控制方程的离散,采用SIMPLEST算法进行流场计算。
2 分析内容与判定标准
该大型运动会火炬塔周围主要有LED显示屏、威亚及看台。看台顶部设有ETFE(乙烯-四氟乙烯共聚物)膜材料。数值仿真的主要内容就是分析火炬燃烧能否对各物体造成不利影响,具体分析如下。
(1)火炬对LED显示屏的影响:
LED显示屏距离火炬46m,采用铝合金材料制成,显示屏元件在温度超过70℃时会出现故障。因此,需要模拟火炬燃烧对LED显示屏的影响,以判断其是否处于安全的温度范围内。
(2)火炬对威亚的影响:
在火炬上方15m处设有威亚,威亚为钢丝制作,用于人员高空表演。钢材的机械性能随温度的不同而有变化,当温度升高时,钢材的屈服强度、抗拉强度和弹性模量的总趋势是降低的,但在150℃以下时变化不大。因此,保守认为威亚的温度不能超过150℃。
(3)火炬对看台顶棚膜材料的影响:
看台距火炬最近51m,看台屋顶的ETFE膜材料在275℃时会熔化,当超过80℃后,可能产生不可恢复的变形,会影响建筑效果。因此,在火炬的各种燃烧工况下,ETFE膜材料的温度不应超过80℃。
3 计算场景
计算场景设置主要考虑最不利情况,火炬塔采用天然气作为燃料,燃烧功率按其燃料最大流量确定,同时考虑了风速影响,其中2m/s为当地平均风速,其他风速值主要考虑燃烧火焰在风的作用下产生偏离后距分析物体最近的情况。共设置5个计算场景,见表1。
4 分析过程及结果
4.1 模型建立
火炬塔总高26m,火炬盆直径6m,高4.8m,因主要考虑火炬塔燃烧对LED显示屏、威亚及看台膜材料的影响,同时考虑计算成本,因此采用各个物体单独计算的思路,分别用PHOENICS虚拟现实前处理器VR-Editor建立了三个计算模型,见图1所示。
4.2 计算条件
环境温度为32℃,将火灾视为质量源和热源的结合体,天然气流量2500m3/h,天然气热值44MJ/m3,燃烧效率为100%。火炬在运动会期间会持续燃烧,因此,采用稳态计算。辐射模型采用IMMERSOL模型,湍流模型采用k-epsilon模型,物体吸收系数为1。采用非均匀网格,火炬燃烧区域局部加密,最大网格尺寸0.5m,最小网格尺寸0.2m。
4.3 计算结果
各场景计算结果见图2-4所示,经分析可得如下结论:
(1)当风速为平均风速2m/s时,LED显示屏最高温度为53℃;当风速为最大风速12m/s时,LED显示屏最高温度为64℃。因此,LED显示屏所能达到的最高温度为68℃。
(2)当风速为平均风速2m/s时,威亚最高温度为63℃;当风速为4m/s时,威亚最高温度为142℃。因此,威亚所能达到的最高温度为142℃。
(3)当风速为6m/s时,膜结构最高温度为44℃,此温度为膜结构所能达到的最高温度。
因此,根据安全判据可知LED显示屏、威亚和看台膜结构在火炬连续燃烧期间是安全的。
5 结语
(1)大型运动会期间,主火炬塔一直处于燃烧状态,由于其燃烧产生的热量巨大,可能会对周围物体产生不利影响,为保证运动会的顺利进行,应对其进行风险评估,以便采取措施,消除安全隐患。
(2)按照本文所提供的方法,采用CFD数值仿真的方法可以较好地完成对火炬塔的安全评估。
(3)CFD仿真计算应考虑最不利场景,并应重点考虑室外风速的影响,在较大风速作用下,火炬塔火焰会发生较大偏转,更容易对周围物体产生不利影响。
(4)在无法进行试验的情况下,采用CFD仿真分析不失为一种好的选择。但条件允许时,宜优先采用试验和数值仿真相结合的方法,这样,试验可以为数值仿真提供边界条件且可以验证仿真结果,同时仿真分析可以为试验提供指导并且可以对试验所不能涵盖的工况进行模拟分析。
参考文献
[1]李万平.计算流体力学[M].武汉:华中科技大学出版社,2004
[2]吴玉剑,潘旭海.障碍物地形条件下重气泄漏扩散实验的CFD模拟验证[J].中国安全生产科学技术,2010,6(3):13-17WUYu-jian,PAN Xu-hai.Simulation and verification ofCFD on dispersion of heavy gas leakage in obstacle terrain[J].Journal of Safety Science and Technology,2010,6(3):13-17
[3]唐保金.燃气管道泄漏及扩散规律的研究[D].济南:山东建筑大学,2009.
[4]郭铜修.基于风险评价技术的企业安全评估体系应用研究[J].中国安全生产科学技术,2008,4(5):19-23GUO Tong-xiu.Applying study on enterprise safety as-sessment system based on risk evaluation techniques[J].Journal of Safety Science and Technology,2008,4(5):19-23
[5]姚征,陈康民.CFD通用软件综述[J].上海理工大学学报,2002,24(2):139-144YAO Zheng,CHEN Kang-ming.Review on the commer-cial CFD software[J].University of Shanghai for Scienceand Technology,2002,24(2):139-144
[6]唐高胤,张礼敬,付道兴,等.南京长江隧道火灾数值模拟[J].中国安全生产科学技术,2011,7(12):74-79TANG Gao-yin,ZHANG Li-jing,FU Dao-xing,et al.The numerical simulation on the Nanjing YANGTZE rivertunnel fire[J].Journal of Safety Science and Technolo-gy,2011,7(12):74-79
[7]张培红,刘岩,宋波.空气幕的不同送风角度对深埋地铁火灾烟气控制数值分析[J].中国安全生产科学技术,2009,5(1):21-25ZHANG Pei-hong,LIU Yan,SONG Bo.Numerical simula-tion on the smoke control effect by the supplying angle ofair curtain in deep buried metro station fire[J].Journal ofSafety Science and Technology,2009,5(1):21-25
[8]Stephen M O,Douglqs J.C.An Updated InternationalSurvey of Computer Models for Fire and Smoke[J].Jour-nal of Fire Protection Enginerring,2003,13(5):87-109
[9]Ludwing J C,Malin M R,Spalding D B.PHOENICSdocumentation[R].Concentration Heat and MomentumLimited,2005
CFD仿真分析 篇6
1 竖直U型地埋管传热分析
U型地埋管换热器与土壤之间的传热过程由以下几部分组成: (1) U型管内循环液体与U型管内壁间对流换热; (2) U型管内外壁之间的导热; (3) U型管外壁与钻井内回填材料之间的传热; (4) 回填材料与周围的土壤进行热交换; (5) 土壤的传热。在系统运行时, U型管中的循环液是靠地能循环泵强制循环的, 因此循环流体与管壁之间的热交换方式是强制热传导换热, 管内流体主要属于对流换热;地能循环泵停止运行时, U型管中的循环液基本是静止的, 此时流体与管壁之间和流体内部是自然传导换热。土壤和回填土是以固体物质为主, 含有一定水分的多孔介质。因此土壤和回填土中的传热包括热传导与热对流, 在传热的同时发生水分迁移, 是传热和传质相复合的过程。
此处我们分析的U型管内循环液是定常不可压缩的, 其密度为常数。为得到较好的换热能力, 一般要使其处于紊流状态。因此在进行CFD模拟时, 选用了k一ε双方程模型来解得水管中水的流动和传热的控制方程, 并将能量方程与回填土和土壤中的传热过程耦合起来, 另外循环液受到物理守恒定律的支配, 即流体应满足质量守恒方程、动量守恒方程和能量守恒方程。由于被系统采用的循环液黏性较小, 其影响仅限于贴近U型内管壁面的薄层内, 从而形成一个边界层, 因此分析管内流体时还应考虑边界层的影响。
2 等直径竖直U型地埋管模型仿真
考虑到土壤成份的多样性, 岩层结构的复杂性及地埋管换热过程的特殊性, 在保证模型精度符合工程要求的前提下, 作以下假设以简化模型:
(1) 地埋管换热器周围土壤均匀性, 即认为不同深度处土壤成份均匀且温度具有一致; (2) U型地埋管两支管具有对称性, 即假定地埋管没有受到安装和填埋所造成的偏离挤压变形等因素影响, 认为两支管沿钻孔轴线对称分布, 且位于钻井的中心部位; (3) 由于U型埋管长径比很大, 相对管长而言管径可以忽略 (此处管长50m, 管径32mm) , 可认为的循环液的温度在径向上没有变化; (4) 假定岩土及循环液的热物性是定值, 不随运行时间和温度场的变化而变化; (5) 忽略U型管外管壁和回填材料、回填材料和土壤之间的接触热阻, 即认为回填材料与U型管外管壁、土壤之间接触良好; (6) 忽略地下水迁移的影响; (7) 钻孔间距足够大, 忽略孔与孔之间的传热影响
2.1 等直径U型管的模型建立
针对本系统的地下部分, 采用FLUENT软件的前处理模块建立几何三维模型, 模型包括的几何体有U型地埋管内的流体、U型地埋管、回填材料和土壤。考虑到地下部分的对称性、节省计算时间等因素, 因此将模型简化为一半。模型的主要参数为:钻孔深度50m, 受到影响的土壤半径取2m, U型地埋管的外径为32mm, 壁厚为3mm, 管材为高密度PE管, U型地埋管两支管中心距为48mm, 钻孔直径为127mm。
2.2 模型的网格划分
网格划分是一个关键步骤, 网格质量对Fluent计算精度和计算效率有直接的影响。由于模型中的U型管是一个细长结构, 相对于U型管回填土加上土壤的体积所占体积大得多, 为了保证网格数量容易控制和较好的网格质量, 所以选择结构化网格来划分网格。所建立的模型共有4个体模型组成, 需分别对其划分网格。网格划分的原则是:在直管部分的管壁和管内流体, 因其温度和速度变化较缓慢, 所以网格可画得疏一些;在弯头部分, 由于速度变化较剧烈, 所以这部分管壁和管内流体应该密一点, 体各体模型的网格划分如如图1～3所示。
2.3 边界条件的定义
边界条件包括流动变量和热变量在边界处的值。在这一环节中, 将U型管内的循环液定义为“FLUID (流体) ”, 材料为warter (水) , 将U型管、回填材料和土壤区域均定义为“SOLID (固体) ”, 固体区域仅需要输入材料类型, U型管为高密度PE, 回填料为膨润土和细沙混合物。然后进行边界条件类型的定义:U型地埋管进口inlet的类型设为“VELOCITY INLET (速度入口) ”, 之后设定进口流速及水温;出口outlet设为“OUTFLOW (自由出流) ”型;径向方向上距离U型地埋管最远处及模型底部认为是无限远处, 定义为同一类型的“WALL (壁面) ”;模型顶部除进出口外还有U型管、回填材料和土壤的顶面, 将它们定义为另一类型的“WALL”;建立的模型为对称模型的一半, 所以要将对称面sym定义为“SYMMERY (对称) ”型, 表明对称面上各参数的梯度都为零。由于是三维非稳态传热问题, 故计算模型的时间特征选择“Unsteady (非稳态) ”, 选择能量方程“Energy Equation”;考虑到管内流体呈紊流时, 换热效果最好, 系统运行时的流速能保证处于紊流状态, 因此设定紊流模型, 控制方程为“k-epsilon (2 eqn) ”;设置进口inlet的速度和温度;还有其它壁面的温度等相关参数。所定义管内循环液、U型管、回填材料及土壤的物性参数见表1。
2.4 模拟计算
在定义了所有边界条件后, 对系统连续运行48小时进行模拟仿真。系统启动时, 测得管内流体温度为21.6℃, 流量为9.5m3/h, 单口井内管内流体流速折算为0.19m/s, 且认定整个测试过程基本不变。为了能较准确的模拟出水温度, 应该舍去系统运行前期所测得的数据, 因为前期U型管内的水和土壤未完全达到热平衡, 且热泵运行也不稳定, 根据实测数据去掉系统运行前6个小时的数据, 得到实测进水温度的平均值为14.3℃, 以此作为模拟进水温度。
本文对土壤温度作了简化处理, 认为U型管内的液体温度即为土壤的初始温度。因为在系统运行前, 让U型管内的循环液与地埋管周围的回填材料和土壤进行充分换热, 达到热平衡, 所以可以认为两者的温度一致。对整个钻井深度实测的管内液体的平均值约21.8℃, 因此我们将土壤初始温度设为21.8℃。
由图7和8可以发现:出水温度的模拟值与实测值的偏差在前几个小时偏差较大, 但是随着运行时间增加偏差缩小。模拟值与实测值的最大相对误差为13% (偏差为1.9℃) , 最小偏差为1.8%, 平均偏差为4.6%, 这在工程允许范围内的。
3 结语
(1) 国内外对地源热泵的研究在理论方面已经比较成熟, 为了使得模拟结果更贴近实际, 近期在国外又引入了守恒方程、模糊理论等方法, 这对于开展研究有一定的导向作用, 但对于一些企业, 由于科研人员质量的原因, 要想以理论为切入点推广地源热泵技术的话还是有一定难度的。
(2) 打实验井虽然能得到比较合理的换热参数, 但是由于所需的场地比较大, 且耗费了人力和物力, 因此企业从经济利益出发为了节约成本, 也不利于广泛推广地源热泵技术。
(3) CFD软件模拟有其特有的优点, 只要建模、网格划分已经边界条件正确就能得到比较贴近实际的第一手数据, 因此对于企业来说还是有比较大的空间的。
参考文献
[1]余延顺, 马最良, 姚杨.土壤蓄冷与土壤耦合热泵集成系统的模拟研究[J].暖通空调, 2005, 35 (6) :1-5.
[2]杨卫波, 施明恒.地源热泵中U型埋管传热过程的数值模拟[J].东南大学学报, 2007, 37 (1) :78-83.
[3]吴玉庭, 顾中煊, 马重芳, 等.U型管传热量影响因素的数值模拟研究[J].工程热物理学报, 2007, 28 (1) :116-118.
[4]袁艳平, 雷波, 余南阳, 等.地源热泵地埋管换热器传热研究 (1) :综述[J].暖通空调, 2008, 38 (4) :25-32.
[5]Hwang Yujin, Lee JK, Jeong YM, etal.Cooling performanceof a vertical ground-coupled heat pump system installed ina school building[J].Renewable Energy, 2008, 34:578-582.
[6]杨卫波, 施明恒, 陈振乾.土壤源热泵夏季运行特性的实验研究[J].太阳能学报, 2007, 28 (9) :1012-1016.
[7]Kara Y A.Ex periment per formance evaluat ion of aclosed-loop vertical ground source heat pump in theheating mode using energy analysis method[J].InternationalJournal of Energy Research, 2007, 31:1504-1516.
[8]王华军, 赵军, 沈亮.地源热泵系统长期运行特性的实验研究[J].华北电力大学学报, 2007, 34 (2) :52-54.
天然气泄漏扩散规律的CFD仿真 篇7
关键词：天然气,泄漏,扩散规律,CFD,模拟
CFD方法,又称计算流体动力学[1],是自计算机出现之后,迅速发展起来的一种应用广泛的计算方法,其基本思想是把原来在时间域及空间域上连续的物理量的场,如速度场和压力场,用一系列有限个离散点上的变量值的集合来代替,通过一定的原则和方式建立起关于这些离散点上场变量之间关系的代数方程组,然后求解代数方程组获得场变量的近似值[2]。
针对管道泄漏的模拟分析,可以将天然气的泄漏过程分为: 小孔泄漏和管线断裂两种情况,前者泄漏口径相对管线直径较小,其水力直径一般不应超过管线直径的30% ,后者适用于泄漏口径较大的情况[14]。
1泄漏理论
传统的泄漏范围确定一般是利用半理论半经验公式,应用fluent软件,成功的对天然气两种泄漏方式: 小口径泄漏和管线断裂分别选取若干工况进行了模拟[3]。
管道的泄漏理论目前被广泛采用的主要有Levenspie,Crow等提出的孔隙模型和管道模型[15]。空隙模型只适用于相对管道尺寸泄漏口径很小的情况,该模型将管道视为一内部压力恒定的大容器,将气体出流视为一稳定过程。因该模型为考虑因紧急切断装置动作而形成的不稳定状态,故其出流处的流速和状态参数均比实际偏大。管道模型只适用于管道完全破裂的情况。天然气管道的泄漏事故规模相对于管道尺寸来说一般是大规模( 大孔径) 泄漏事故,但相对于整个扩散范围来说,天然气的泄漏属于大气淹没射流,均可看作小孔射流。
为了计算泄漏速度将天然气的泄漏过程进行简化:
( 1 ) 将天然气作为理想气体来考虑;
( 2) 假设其泄漏过程为小孔喷流方式;
( 3) 气体从小孔中喷出的过程看做是绝热过程。
利用伯努利方程和绝热方程,泄漏速度v0的计算公式为:
式中:——— 流速系数,表示实际流速与理论流速之比,一般为0. 97 ~ 0. 98
k———绝热指数,即定压热容与定容热容之比,对于天然气一般取k = 1. 3
R———气体常数,8. 314 J / ( mol·K)
T1——— 气体温度, K
P0——— 环境压力, Pa
P1———气体泄漏前压力,也就是管道内介质的压力,Pa
2泄漏模型建立
由于计算空间比较大,为了保证网格精度以及计算的快速收敛,考虑对模型进行以下简化:
首先,虽然天然气在大气中的扩散是三维的,但由于本文并不考虑地形起伏作用对扩散的影响,可认为扩散在空间中各个二维截面上是相同的,因此考虑建立二维模型[5]。
实际泄漏过程中,泄漏率会随着泄漏的进行而逐步降低, 并最终与大气压达到平衡[6],本文中,考虑到实际中泄漏时间并不唯一,因此采用固定泄漏速度的方法,检测泄漏点周围参考点处的天然气浓度随泄漏时间的变化,可以将模拟结果推广到不同泄漏时间、不同泄漏量的多种情况。
在以上简化基础上,取泄漏周围100 m × 100 m的区域为计算区域。
( 1) 网格划分
由于在泄漏孔出口附近有较高的压力及速度梯度[7],因此划分较密的网格,以保证计算精度及收敛速度。
( 2) 边界条件及初始条件设定
首先启动非定常求解器,选用标准k - ε 模型来考虑湍流的影响,选择组分输运模型,在材料面板中加载天然气,在操作压力选项板中勾选gravity选项,将重力加速度设为y方向 - 9. 8,操作温度为288. 16 K,天然气和空气均为理想气体, 操作密度为空气密度1. 225[8]。设置完成后进入边界条件设置面板:
速度入口边界为已知速度: 104m / s。
湍流边界条件: 湍流边界条件的设定是模型解算的重要部分,在此进行介绍:
( 1) 湍流强度
定义: 速度波动的均方根与平均速度的比值,小于1% 为低湍流强度,高于10% 为高湍流强度[9]。
计算公式:
式中: I———湍流强度
Re———雷诺数
( 2) 湍流尺度及水力直径
湍流尺度 ( turbulence length) 通常计算方式:
式中: L———特征尺度,可认为是水力直径
0. 07———因数,基于充分发展的湍流管流中的混合长度的最大值[10]
在species面板中,设定天然气的入口体积分数为1。
压力出口条件,设定为一个大气压 。
模拟扩散时间为120 s。
3模拟结果分析
针对无风和有风两种工况进行模拟,两种工况的参数见表1。
3.1无风工况模拟结果分析
工况一运行120 s后计算完成,分别给出扩散0. 05 s、60和120 s时天然气浓度云图,结果见图1 ~ 图3 ( 图例中数字表示天然气浓度,单位为1,下同) 。
由天然气扩散过程可知,无风时,天然气扩散时,浓度范围主要分布在泄漏口上方垂直于地面的空间内,且向周围区域扩散并不明显[12]。
3.2有风工况模拟结果分析
工况二的模拟分析将入口尺寸改为0. 2 m,左侧边界由出口边界改为入口边界且速度为112 m/s,其他参数设置同工况一,分别给出扩散0. 05 s、60 s和120 s时天然气浓度云图,结果见图4 ~ 图6。
由天然气扩散过程可知,有风时,天然气扩散时,浓度范围具有明显的定向性,大都分布在风的下风向,且主要分布在泄漏口上方垂直于地面的空间内,且向下风向区域扩散明显, 风速越大,定向性越明显,分布规律越显著。
3.3模拟结果对比分析
( 1) 管道输送压力的影响
相同的管道在不同的管道压力下扩散的范围也不相同,一般来说管道压力上升则扩散范围增加。
( 2) 风向及风速的影响
风向主要影响泄漏气云的扩散方向,大部分泄漏气体总是分布在下风向。
风速影响泄漏气云的扩散速度和被空气稀释的速度,风速越大,则大气的湍流越强,空气的稀释作用就越强。一般情况下当风速为1 ~ 5 m/s时,有利于泄漏气云的扩散,危险区域较大; 若风速再大,则泄漏气体在地面的浓度变稀。若无风天, 则泄漏气体以泄漏源为中心向四周扩散。
( 3) 泄漏孔径的影响
泄漏孔径增大,使得扩散范围增加。
4结论
CFD仿真分析 篇8
关键词：抗性消声器,计算流体动力学,压力损失,数值仿真
1 引言
抗性消声器因其结构简单、制造容易而广泛应用到发动机排气噪声的治理中,但在实现降噪的同时,排气动力学性能也会发生一定的变化,为了保证发动机正常工作,插入消声器的压力损失必须控制在允许范围内。因此有必要研究影响抗性消声器压力损失的因素,并给出这些因素对压力损失的影响规律。目前研究消声器压力损失的传统方法一般采用经验公式法和试验法,显然这两种方法不能满足复杂结构消声器的设计要求。本文借助于计算流体动力学(CFD)技术,仿真了抗性消声器内部气体流场分布,分析了影响内部流场的因素,计算了各因素影响下的压力损失,为复杂消声器的设计提供了一定的理论基础。
2 抗性消声器的仿真分析
以图1所示的抗性消声器[1]为分析对象,对其进行仿真分析。该消声器为中间连通管的双级膨胀腔消声器,左端为入口,右端为出口。
采用标准k-ε双方程模型[1]对该消声器内部流场进行数学建模,并用下式表示:
其中,ρ为空气密度,准为通用变量,u为速度矢量,Γ为准相对应的广义扩散系数,S为广义源项。
模型计算边界条件的确定:消声器进口气体的流速和温度为速度入口边界条件;出口背压为压力出口边界条件;气体在消声器内部流动,故消声器壳体是固体边界,假设固壁绝热,无摩擦,在无滑移条件下,内表面流体的速度等于壁面的速度且为零。对于消声器废气的物性参数通过温度插值多项式来实现。该插值多项式可表示为:
视流过消声器的气体为不可压理想流体,并考虑到气体的实际流动性,暂按表1的参数设定排气废气的物性参数。
按照以上的分析与设置,采用三维数值仿真软件Fluent对该消声器的内部流动进行仿真计算,当进口流速设定为20m/s,温度为38℃时其计算结果如图2至图5所示。
由图可以看出,在消声器中心轴线处气流速度较大,内部速度变化规律为从壁面处开始沿径向呈抛物线形。内部压力变化较大处发生在内插管两端部,该处流速与压力都急剧下降,从而产生较大漩涡,导致压力损失。另外,气体漩涡的存在必然会导致再生噪声,故在消声器的设计中应考虑漩涡气流的影响。
为了更清楚地说明问题,对仿真计算出的压力损失与实测压力损失进行比较,如表2所示。从表中数据差异可以看出,两者差别较小,具有一致性,从而表明了用CFD方法进行消声器压力损失的分析是可行的。
3 改变结构后的消声器压力损失计算及分析
在验证了仿真算法可靠性的基础上,本文对算例消声器的结构做了几种改变,分别研究了内插管、中间挡板位置以及进口空气流速对消声器压力损失的影响。消声器速度入口边界条件均为:0~100m/s(步长为10m/s),11个速度下分别计算,进气温度为38℃。出口背压为压力出口边界条件,且为一个标准大气压;消声器废气的物性参数由式(2)和表1计算得到。
3.1 内插管对消声器压力损失的影响
将算例消声器结构中的内插管去掉,其余结构保持不变,在Fluent软件中计算两种结构消声器在进口空气流速为0~100m/s时的压力损失值。计算结果如图6所示。
图6中曲线1代表算例消声器压力损失随入口流速的变化而变化的趋势,曲线2代表算例消声器去掉内插管后压力损失随入口流速的变化而变化的趋势。因此,两种结构消声器的压力损失规律均与流体流速成正比。由于内插管安装造成内插管出口流速比同位置没有内插管的流速大且压降较小,因此无内插管消声器的压力损失远大于带内插管消声器的压力损失。
3.2 中间挡板位置对消声器压力损失的影响
本文中的算例消声器是一个双膨胀腔消声器(如图1),它的第一膨胀腔长度为600mm,第二膨胀腔长度为400mm。为了研究探讨中间挡板位置对消声器压力损失的影响,在这里对算例消声器的结构做如下两种改变:将算例消声器的第一膨胀腔长度减小为500mm,此时第二膨胀腔长度为500mm,为了方便描述,将此种结构消声器命名为类型1消声器;将算例消声器的第一膨胀腔长度增加为700mm,此时第二膨胀腔长度为300mm。为了方便描述,将此种结构消声器命名为类型2消声器。在Fluent软件中计算类型1消声器和类型2消声器在进口空气流速为0~100m/s时的压力损失值,计算结果如图7所示。
图7中曲线1、2、3分别代表算例消声器、类型1消声器、类型2消声器的压力损失随入口流速的变化规律。该规律依然是压力损失随流速的增加而变大。
如图7,改变结构后的两种类型消声器的压力损失均略大于算例消声器的压力损失,而类型1与类型2消声器的压力损失相差较小,这表明对于双级膨胀腔消声器的压力损失来说,中间挡板位置对其影响较小。类型1消声器的挡板位于膨胀腔的中间位置,类型2消声器的挡板使第一、第二膨胀腔的长度差别较大,而算例消声器挡板两边的膨胀腔长度差别不是很大,算例消声器的压力损失也是最小的。由此可以看出,当双级膨胀腔消声器的挡板位于中间位置或挡板使两边膨胀腔长度差别较大时,消声器压力损失较大。
4 结论
本文采用计算流体力学方法仿真分析了某一算例消声器的内部流场,并计算得到了该算例消声器的压力损失,将数值计算结果与试验结果相比较,可以看出二者吻合良好,这说明采用计算流体动力学方法研究消声器的压力损失是可行的。在此基础上,对算例消声器的结构进行了一些改变,研究探讨消声器基本结构和进口空气流速对消声器压力损失的影响,得到如下结论:
(1)内插管对消声器的压力损失有比较大的影响,无内插管消声器增加内插管可以使消声器的压力损失有较大幅度的降低。
(2)中间挡板的位置对双级膨胀腔消声器的压力损失有一定影响,但影响比较小,当双级膨胀腔消声器的挡板位于中间位置或挡板使两边膨胀腔长度差别较大时,消声器压力损失较大。故工程设计中可通过合理安排双级膨胀腔消声器的挡板位置达到减小消声器压力损失的目的。
(3)一般情况下,随着进口气体流速的增加,消声器的压力损失变大。
参考文献
[1]郭小林.排气冷却消声器性能研究及热应力分析[D].哈尔滨工程大学,2008.
[2]张智辉,陈军,王树宗,等.消声器内部流场的数值模拟[J].振动与冲击,2006,25(6):21-24.
[3]朱红均,林元华,谢龙汉.Fluent12流体分析及工程仿真[M].北京:清华大学出版社,2011.