贮箱内低温推进剂汽化过程的CFD数值仿真

陈 亮 梁国柱 邓新宇 胡 炜

(北京航空航天大学 宇航学院，北京100191)(北京宇航系统工程研究所，北京100076)

液体运载火箭贮箱在加注低温推进剂后至发射前一直处于地面停放状态，贮箱壁面漏热会使推进剂的温度升高并发生汽化，减少推进剂的可用量，并带来安全隐患，因此需要对贮箱内传热传质过程，尤其是汽化过程进行研究.国内外对贮箱内低温推进剂的传热传质过程的仿真采用的方法主要有集总参数法和CFD(Computational Fluid Dynamic)方法.文献［1］对增压贮箱进行了集总参数建模，将贮箱内部分为增压气体、液态推进剂主体部分和气液界面3部分，根据热力学平衡原理建立了相变模型.文献［2］采用集总参数法针对液氧贮箱的增压过程建立了贮箱内气相空间的集总参数模型，认为气枕与液氧表面间的传质过程是由液氧的饱和蒸气压与氧气分压之差决定的.集总参数法在对贮箱的研究中应用比较广泛，但不能反应贮箱中各项参数的空间分布情况.本文研究的贮箱尺寸较大且温度分布对相变有较大影响，因而不适用于集总参数法.

国内外采用CFD方法的研究相对较少.文献

［3］使用有限体积法研究通入热增压气体后液氧贮箱内的物理变化，但忽略了相变过程，文献［4］在此基础上加入了相变模型，认为引起液氧相变的主要原因是增压气体对液氧的加热作用，气体向推进剂的传热是由自然对流引起的.参考文献［5－6］对贮箱内液体主体区域采用二维的CFD方法建模并与气枕的集总参数法模型耦合，计算得到的增压速率与采用集总方法的热力学分析结果相符.文献［7］采用CFD方法建立了液氧贮箱的2维模型，根据热力学平衡原理建立了相变模型.这些研究多采用类似文献［4］的相变模型，气液两相之间换热系数需根据流体的运动状态求解，但在实际情况中，液态推进剂中可能出现大量“气泡”，流动状态复杂，不易求解.本文基于实际问题，根据基本的热力学平衡原理构造相变模型，采用CFD方法对贮箱中的物理过程进行仿真.

1 物理过程分析

本文研究的推进剂贮箱为如图1的高约7m、直径约3 m的液氢贮箱，氢箱与氧箱共底，氢箱后底和氧箱前底间为真空隔热层，贮箱内壁采用铝合金，铝合金壁面外包覆绝热层.贮箱处于地面停放状态，柱段外壁面直接与外界环境接触，两者之间的换热过程为强制对流过程，受外界温度和风速影响.氢箱前底外为火箭舱室，舱内充填氮气，前底与舱室内气体的换热过程可视为自然对流过程.氢箱顶部开口，通过管路与外界大气相通.外界环境以及氧箱(液氧温度明显高于液氢)通过壁面漏热加热贮箱中处于饱和状态的液氢，使液氢发生汽化，这一过程属于非均匀大容器沸腾［8］.在贮箱内部发生的物理过程主要包括推进剂的自然对流过程、推进剂与贮箱内壁面的对流换热过程、贮箱壁的导热过程以及推进剂的相变过程.其中，相变过程包括液体推进剂主体和气枕间的界面上的蒸发和冷凝现象以及液体推进剂主体内部(主要是壁面附近)的汽化现象.

图1 贮箱示意图

2 数学模型

为满足仿真需要，贮箱的数学模型除了包括多相流问题的通用控制方程，还需要建立推进剂相变模型、真实气体模型和固体壁面导热模型.

2.1 流体控制方程

控制方程中考虑流体粘性影响，忽略相变引起的动量变化.

连续方程:

式中，ρ为流体密度;v为速度矢量;Sm为质量源项;表示气液两相间传质过程引起的单位体积的质量交换率.

动量方程:

能量方程:

式中，E为能量;keff为有效热传导系数为有效应力张量;Sh为能量源项，表示两相间传质过程引起的单位体积的能量交换率.

多相流方程

多相流模型选用可以追踪相界面的VOF模型［5］，容积比率的连续方程为

式中，αq，ρq和vq分别为第q相的体积比率、密度和速度矢量;Sαq为第q相的质量源项.

2.2 推进剂相变模型

低温推进剂的相变过程主要出现在加热面(即贮箱内壁面)和气液界面上(包括液体推进剂内部的气泡界面和液体推进剂主体与气枕间的界面).由于贮箱内推进剂近似处于饱和状态，且贮箱压力保持稳定，因此可假设液态推进剂及推进剂蒸汽处于热力学平衡，则气液界面上应满足力平衡条件pg=pf，即气液两相压力相等(不考虑张力影响);热平衡条件Tg=Tf，即气液两相温度相等;相平衡条件pg=pf=psat(T)，即气液界面上压力等于对应温度下的饱和压力［8］.在气液界面上，当推进剂(i表示液氢或者液氧等低温推进剂)的蒸汽分压pig大于气相区域温度对应的饱和蒸汽压psat，i时，发生冷凝，此时气液界面上发生的质量交换和热交换［1］分别为

式中，Δt为时间步长;Vig为气相体积;Mi为推进剂摩尔质量;Til为液态推进剂温度;R为气体常数;Δhi为汽化潜热.

当pig小于psat，i时，发生汽化.此时气液界面上发生的质量交换率和能量交换率分别为

在贮箱壁面附近，壁面漏热会引发推进剂的汽化.在发生汽化前，液态推进剂主体内部在内壁面附近没有气体，因此前述相变模型不适用.假设壁面附近处于热力学平衡状态，根据相平衡原理应满足Tf=Tsat(p)，即液体温度等于对应压力下的饱和温度，根据能量守恒的原理可得壁面附近的汽化过程的质量交换率和能量交换率分别为

式中，ρl为液态推进剂密度;Tsat，i为饱和温度;cpi为液态推进剂的比热容.

连续方程的质量源项Sm由式(5)、式(7)、式(9)决定;能量方程中的能量源项Sh由式(6)、式(8)、式(10)决定.由于相变模型建立在热力学平衡假设的基础之上，因此可用于近似处于饱和状态且压力温度都不发生显著变化的系统中，而对于过热沸腾问题或压力温度发生剧烈变化的系统不适用.

2.3 其他模型

由于氢气在低温时的性质与理想气体有一定差异，故根据试验数据建立真实气体模型如下:

式中，Z为压缩因子，通过对试验数据拟合得到，Z随温度升高而增大，当压力为0.08 MPa，温度为20 K 时，Z 约为0.9.

贮箱固体壁面的传热不仅包括垂直于壁面方向的传热过程，也包括沿壁面方向的传热，能量输运方程如下:

式中，ρw为固体壁面的密度;hw为壁面材料的显焓;kw为热传导系数;Tw为壁面的温度.

3 计算方法

贮箱仿真模型通过Fluent建立.模型网格分为流体区域和固体壁面区域(包括铝合金贮箱壁和绝热层)，靠近固体壁面的流体区域网格进行适当加密.通过试算对比选定节点数为63 818的2维轴对称网格，经验证，不同划分方式和节点数的网格的计算结果相近.

初始条件设置为推进剂处于饱和状态.受重力影响，贮箱内不同深度的压力不同，因此，将液氢的初始温度设置为对应压力下的饱和温度.为简化问题，贮箱固体壁面区域的初始温度设置为与同一深度的液氢的温度相同(实际情况中，在加注液氢前要对贮箱进行预冷，预冷后壁面温度接近液氢温度).

氢箱柱段和前底的外壁面设置为第3类传热边界条件，壁面换热系数根据对流运动状态计算得到.环境压力为0.08MPa，温度为30℃，风速为3 m/s，计算得到柱段外壁面的强制对流换热系数为9.008W/m2/K;氢箱前底外的火箭舱室内充填20℃的氮气，自然对流换热系数为 0.937～1.007 W/m2/K，随外壁面温度变化而变化.氢箱与氧箱的共底设置为第2类传热边界条件，氧箱通过共底向氢箱的平均漏热率根据试验结果设置，为83 W/m2.湍流模型采用适合低雷诺数的RNG k-ε模型，壁面处理方式采用增强壁面函数.求解器选择基于压力的求解器，压力插值方式选用body-force-weighted方案，速度压力耦合方式选用压力和速度耦合的coupled方案.相变模型通过UDF以质量源项和能量源项的形式加入到Fluent模型中.经多次试算比较，时间步长选定为0.04s，实际数值计算的连续性、动量、能量方程的残差分别为 10－4，10－6，10－10量级.

4 计算结果与分析

贮箱的初始压力、温度和两相分布如图2所示.顶部压力为0 Pa，底部为4957 Pa;顶部温度为19.502 K，底部为19.697 K，温度差为 0.195 K，液氢处于饱和状态，呈明显的温度分层现象.随着计算时间的推移，贮箱内推进剂出现对流运动和汽化现象.t=10.0 s和50.0 s时刻贮箱内温度、局部速度、两相和相变分布分别见图3和图4.由于侧壁面和底面漏热，靠近壁面的液氢温度开始升高(如图3a)，贮箱底部开始产生较明显的汽化现象(如图3c、图3d，图3d中黑色表示汽化，白色表示冷凝，灰色表示无相变).汽化产生的气氢带动液氢沿侧壁面和底面向上流动(如图3b).底部温度较高的饱和液氢沿侧壁面向上流动变为“过热”，汽化加剧，上部液氢中出现大量氢气(如图4c);贮箱上部靠中间位置、温度较低的液氢向下流动，变为“过冷”流体，汽化现象逐渐消失，并开始出现冷凝(如图4d).

图2 初始状态的贮箱内各参数分布

图3 t=10 s时刻贮箱内各参数情况分布

t=450 s时刻贮箱内温度、局部速度、两相和相变分布见图5.

图4 t=50 s时刻贮箱内各参数情况分布

图5 t=450 s时刻贮箱内各参数情况分布

随着壁面温度的升高，外壁面与外界环境的换热和内壁面与流体的换热逐渐达到平衡，贮箱内液氢的温度形成相对稳定的中间温度较低，靠近侧壁面和底面的区域温度较高的“V型”温度分布(如图5a).贮箱内液氢的流动状态也随之趋于稳定，形成底部液氢沿侧壁面向上流动，上方液氢从中间向下流动的循环流动(如图5b).随着温度分布和流动状态趋于稳定，贮箱内的汽化现象也趋于稳定，汽化现象主要集中在贮箱上部壁面附近(如图5c，图5d).

仿真中贮箱壁面初始温度设置为与推进剂相同，壁面温度受外界环境传热影响而上升，因此，外壁面换热率减小，内壁面换热率增大，换热率变化曲线如图6a所示，可见内外壁面的换热率趋于相等，300s后基本稳定在约20000W.由于贮箱通大气，箱内压力接近环境大气压，因此贮箱出口流量近似等于液氢单位时间的汽化量(即汽化速率).贮箱出口流量变化曲线如图6b所示，在450 s仿真中，液氢汽化速率在100 s以后开始趋于稳定(100 s以后平均值为23.4 g/s)，结合前文分析，液氢的相变与其温度分布和流动状态密切相关，汽化速率趋于稳定可以间接反映出随着贮箱壁面传热过程逐渐稳定，液氢的温度分布、流动状态和相变情况也逐渐趋于稳定.

图6 换热量和出口流量变化曲线

假设贮箱壁面漏热全部作用于相变过程，推进剂的汽化速率约为44 g/s，与仿真结果相比有出入，这一方面由于壁面漏热中约有15%(与气枕相接触的内壁面对气枕的换热量约为3000 W)用于加热气体并被带出贮箱;另一方面，由于贮箱尺寸较大而仿真时间相对较短，虽然推进剂汽化速率已经趋于稳定，但贮箱内温度分布还继续缓慢变化，一部分壁面漏热加热从上部流入底部变为过冷的液体而尚未直接影响到相变过程.因此，壁面漏热虽是引起汽化的关键因素但并不能直接决定汽化速率，壁面漏热对相变的影响和贮箱内各物理过程间的复杂的相互作用机理有待进一步研究.此外，相变模型建立在热力学平衡假设基础上，根据仿真，贮箱内推进剂的过冷度或过热度在0.1 K以内，接近饱和状态，但局部区域的不平衡状态仍可能对仿真造成一定影响，因此相变模型也需要进一步改进完善，以解决非平衡状态下相变过程的仿真问题.

5 结论

本文根据热力学平衡原理建立了低温推进剂的相变模型，通过对地面停放状态下的液氢贮箱的仿真获得了以下结论:

1)随着贮箱壁面传热过程的稳定，贮箱内部热力学过程趋于稳定.贮箱内液氢出现底部沿侧壁面向上，上部由中间向下的循环流动;液氢温度呈中间温度低，侧壁面和底面附近温度高的“V型”分布.

2)受推进剂对流运动的影响，贮箱内出现汽化的区域主要集中在贮箱上部靠近内壁面处.在450 s的仿真中，液氢的汽化速率趋于稳定，100 s后平均值约为23.4 g/s.

3)影响液氢汽化速率的主要因素的壁面漏热和推进剂对流运动.壁面漏热使靠近壁面的液氢温度升高发生汽化，还影响推进剂的对流过程，是影响汽化速率的关键因素.

References)

［1］Zilliac G，Arif Karabeyoglu M.Modeling of propellant tank pressurization［R］.AIAA-2005-3549，2005

［2］王赞社，顾兆林，冯诗愚，等.低温推进剂贮箱增压过程的传热传质数学模拟［J］.低温工程，2007，160(6):28－37 Wang Zanshe，Gu Zhaolin，Feng Shiyu，et al.Simulation of heat transfer and mass transfer in cryogenic propellant tank pressurization process［J］.Cryogenics，2007，160(6):28－37(in Chinese)

［3］Majumdar A，Steadman T.Numerical modeling of pressurization of a propellant tank［R］.AIAA-99-0879，1999

［4］代予东，赵红轩.运用数学方法模拟推进剂贮箱增压［J］.火箭推进，2003，29(3):34－40 Dai Yudong，Zhao Hongxuan.Simulation of propellant tank pressurization using mathematics method［J］.Journal of Rocket Propulsion，2003，29(3):34－40(in Chinese)

［5］Panzarella C H，Kassemi M.On the validity of purely thermodynamic descriptions of two-phase cryogenic fluid storage［J］.J Fluid Mech，2003，484:41－68

［6］Panzarella C H，Kassemi M.Self-pressurization of large spherical cryogenic tanks in space［J］.Journal of Spacecraft and Rockets，2005，42(2):299－308

［7］Cheng Xianghua，Li Yanzhong，Chen Erfeng，et al.Effect of return inlet on thermal stratication in a rocket tank［J］.Journal of Thermophysics and Heat Transfer，2010，24(1):112－122

［8］徐济鋆.沸腾传热和气液两相流［M］.北京:原子能出版社，2001:215－220 Xu Jijun.Boiling heat transfer and gas-liquid two-phase flow［M］.Beijing:Atomic Energy Press，2001:215－220(in Chinese)