核热推进系统分析程序模型与计算方法初步研究

毛晨瑞,吉 宇,孙 俊,郎明刚,石 磊

(清华大学 核能与新能源技术研究院,先进反应堆工程与安全教育部重点实验室,北京 100084)

核热推进(NTP)是利用核能将推进剂加热到高温从喷管膨胀并高速喷出产生推力的空间推进装置,主要由推进剂供应系统、反应堆系统以及喷管等部件组成。推进剂供应系统主要包括推进剂贮箱、涡轮泵、管阀和测控仪表等部件,反应堆系统包括辐射屏蔽、堆芯和控制鼓等部件[1]。与化学火箭发动机相比,NTP具备更高的比冲和与之相当的推力。因此,在时效性要求较高的载人深空探测和轨道机动等任务方面具有巨大的应用潜力。美国和前苏联于20世纪50年代便开始围绕洲际导弹推进、战略导弹助推段拦截和大型有效载荷空间轨道转移等场景开展了NTP系统研究,通过Rover/NERVA、Timberwind/SNTP和SEI等项目,取得了一系列重要成果,具体包括较为完备的NTP系统性能分析方法与程序、NTP反应堆方案、耐高温抗氢蚀核燃料以及关键地面试验设施等,前苏联更是建成了迄今为止成熟度最高的RD-0410核热火箭发动机地面样机[2-4]。但随着化学推进的快速发展,NTP系统丧失了其原有的任务需求,众多项目陆续中止。

近些年来,由于突出的性能优势以及在众多航天任务中的不可替代性,美国重新启动了对NTP系统的研究,美国国家航空航天局(NASA)在2014年开始的核低温推进级(NCPS)项目中以安全、高效、经济作为主要的研发原则[5]。此外,美国国防先进研究计划局(DARPA)在2020年牵头实施的“敏捷的月空间飞行火箭验证项目”(DRACO)目前正快速推进,并计划于2027年实现在轨示范验证[6]。由于NTP工作温度极高,且反应堆系统与推进剂供应系统存在极强的耦合关系,存在实验难度大和研发成本高等挑战,所以针对NTP开发专用的系统设计分析程序是必要的。本文针对NTP系统性能分析需求,提出主要研发内容,并对已有的分析程序与框架进行评估,开展NTP系统瞬态模型以及计算方法研究。

1 NTP系统设计分析概念与内涵

NTP相比传统的反应堆装置存在以下差异:1) 系统为开放的一回路设计,规模小,部件间耦合紧密,动态响应快;2) 系统自身的设计裕量小,尤其是反应堆堆芯;3) 在快速启动、停机和推力调节等典型瞬态过程中,存在系统和堆芯物理热工参数的剧烈变化;4) 受空间应用的限制,无法设计类似于传统反应堆的安全措施和冗余手段。基于以上差异,水堆或高温气冷堆所采用的系统设计与分析方法在NTP中可能会导致较大的偏差,甚至无法获得合理可行的方案。因此有必要对NTP系统的设计与分析方法开展进一步研究。

1.1 NTP系统设计分析思路

根据NTP特点及前期探索,将系统设计和分析活动分为稳态设计点分析与优化、稳态非设计点性能分析与瞬态性能分析3部分,如图1所示[7]。对NTP系统设计而言,首先需要选定装置的一个主要技术状态作为设计和分析的基准,也就是设计点,来进行热力计算。在这个过程中,依据论证的系统构型进行热力循环参数的确定和优化,具体涵盖对要求的性能参数(主要是比冲和推力)、设计极限(如涡轮许用温度限值和各部件的效率等)和设计选择(如泵扬程、涡轮压比和分流比等)之间的关系进行分析,获得NTP系统性能参数的变化规律,作为优化的方向和依据。对空间装置而言,有时还需要根据热力循环参数对系统的质量进行初步估计,实现性能和质量的双目标优化。在这个阶段中,NTP系统可看作是柔性的“橡皮泥”,尺寸和参数都是不定的,有较大的自由度和优化空间。

图1 通用NTP设计程序框架Fig.1 Generic NTP design program framework

经过设计点分析与优化确定初步的热力循环方案后,则可以根据节点参数开展初步的部件和流道设计工作,获得各部件的工作特性。非设计点性能分析是对处于NTP系统运行包线内的所有点进行分析,求解各部件的共同工作状态,获得系统性能参数的变化规律等。此时,可设想完成了一台暂缺详细结构信息的NTP装置的“标准设计”。针对多种“标准设计”开展非设计点性能分析的工作,可进一步优化得到能满足任务需求的NTP系统。非设计点性能分析后,可进一步研究系统和部件参数在各非设计点之间的动态变化过程,即瞬态性能分析。所获得的特性可作为NTP系统启动、停机和推力调节等过程运行方案设计的依据。同时瞬态性能分析工作也可用于研究反应性误引入和叶轮机械断轴等事故场景下的系统参数变化规律。对在动态过程中出现参数超过设计限值的情况,有必要与设计点优化和非设计点性能分析等工作迭代,从而给出更优的系统设计。

1.2 稳态设计分析程序

NTP系统设计点分析与优化以及非设计点性能分析均是基于稳态设计分析程序而开展的。NESS(Nuclear Engine System Simulation)是为了支持美国空间探索计划(SEI)而由美国科学应用国际公司(SAIC)为NASA格伦研究中心(GRC)开发的NTP稳态设计分析程序。该程序结合了西屋公司开发的反应堆模型ENABLER与Aerojet通用公司开发的液体火箭发动机模拟程序ELES[8]。NESS可以设计和确定NTR推进系统组件(包括反应堆)的重量、性能和工作特性,以及发动机子系统参数。NESS程序与NERVA项目成果小型核火箭发动机(SNRE)设计结果进行了对比验证[9]。

除了NTP专用稳态设计程序NESS外,NASA的GRC曾使用通用的NPSS(Numerical Propulsion System Simulation)程序对NTP系统进行了稳态性能分析工作。NPSS被航空航天工业广泛用于叶轮机械、吸气式推进系统、液体火箭发动机、发动机控制系统和系统模型集成的建模[10]。该软件具备初步设计、非设计点性能分析、瞬态性能分析和飞行试验数据关联等能力[11]。

清华大学核能与新能源技术研究院自主开发了稳态性能设计分析软件PANES(Program for Analyzing Nuclear Engine Systems),该软件具备NTP系统循环参数计算、颗粒床反应堆轮廓方案、反射层冷却方案、喷管构型方案设计以及堆芯流动换热特性计算的能力[1,12]。

1.3 瞬态分析程序

表1列出了美国NTP瞬态分析程序情况,针对20世纪开展的NERVA、SNTP以及21世纪开展的NCPS等项目,美国开展了多型NTP系统瞬态分析程序的开发工作。国内的多家单位也开展了NTP系统设计、分析以及应用研究,刘忠恕[13]使用国产Modelica分析平台MWorks建立了NTP系统整机模型,并开展了动态特性研究;钟岩俊[14]使用AMESIM仿真程序开展了小推力核热火箭的循环方案研究;Qi等[15]开发了NTPSS-vPower程序以开展NTP瞬态建模与研究。

表1 美国NTP瞬态分析程序Table 1 NTP transient analysis program of US

随着计算技术发展与仿真技术的发展,瞬态分析软件由专用简化系统模型向通用分析平台及全范围系统模型过渡。而在建模方式上,这些系统分析程序依然采用了“流网-热网”的建模思路。

2 NTP系统瞬态分析程序框架

按照NTP系统设计分析思路,为了形成一体化设计分析能力以提高设计效率,在现有稳态程序PANES的基础上采用“热网-流网”的框架自主开发了NTP系统瞬态分析功能。

2.1 流网模型

以“流阻/感-流容”的形式建立NTP系统分析的流网模型。喷管、阀门、离心泵、涡轮等部件模型可以简化为单一的流阻/感模型。堆芯内流道以及再生冷却通道需要在流阻/感模型的基础上考虑流容模型。在流容模型中考虑流体的可压缩性以及能量输入,将质量连续性方程与能量守恒方程进行简化后如式(1)、(2)所示。在流阻/感模型中将动量守恒方程简化后可获得考虑流体惯性以及阻力系数的质量流量计算方程,如式(3)所示。

(1)

(2)

(3)

其中:p为压力,Pa;t为时间,s;ρ为流体密度,kg/m3;G为质量流量,kg/s;Kβ为体积弹性模量,Pa;a为流体声速,m/s;V为控制体容积,m3;h为流体比焓,J/kg;Qi为流网节点与热网节点间的换热功率,W;Rf为阻力系数,Pa/(kg·s-1)2;A为通流面积,m2;L为管道长度,m;下标1、2分别表示对应控制体进、出口处参数。

图2为流网基础单元示意图,选择流网环节中进口压力、进口比焓以及出口质量流量为流网边界,建立以“流阻/感-流容”(R/I-C)连接顺序的流网基础单元。多个基础单元相连即可描述包含多个节点的任意流网环节。

图2 流网基础单元示意图Fig.2 Schematic diagram of flow network basic unit

2.2 热网模型

与流网模型类似,以“热容-热阻”(R-C)的形式建立燃料组件、喷管壁等热构件的热网模型(图3)。对于热容环节,其连接的热阻环节不局限于图3中给出的形式,可根据建模的精度连接任意热阻环节。式(4)为热阻模型,用于描述固体-流体间对流换热、固体与固体间的一维导热,式(5)为能量守恒下描述热构件温度变化的热容模型。

图3 热网基础单元示意图Fig.3 Schematic diagram of thermal network basic unit

q=A(T2-T1)/Rh=hA(T2-T1)

(4)

(5)

其中:q为热流量,W;T为热构件平均温度,K;Rh为等效热阻,(m2·K)/W;h为对流换热系数,W/(m2·K);A为等效换热面积,m2;m为热构件质量,kg;c为热构件平均比热容,J/(kg·K)。

2.3 瞬态求解算法

通过建立“流网-热网”的物理模型即可获得描述NTP系统的数学模型即常微分方程组(ODE)。求解该方程组,便可得到NTP系统瞬态参数。对该问题,有显式隐式两类求解方式。与显式方法相比,隐式方法通常更加稳定,具有更大的稳定时间步长范围,对初始化要求更低。但是隐式方法中求解非线性方程组需要消耗更多计算资源。同时,由于隐式求解中,各部件关键参数需要耦合求解,提高了程序框架设计以及开发的难度。

基于“流网-热网”框架,可得到描述整个系统的数学模型,并建立系统性能分析的流程,如图4所示。在稳态分析中,可以将ODE中时间导数项置零,此时,稳态求解问题即转化为非线性方程组求根问题,采用牛顿法即可解决。在瞬态分析中,由于隐式方法同样需要求解非线性方程组,因此,在实现隐式方法的同时可以改进稳态分析程序,提高稳态计算速度以实现快速的性能方案设计、分析以及优化。高效的稳态求解也可以为瞬态分析提供可靠的初始值。

图4 NTP系统瞬态分析程序计算流程示意图Fig.4 Schematic diagram of calculation process for transient analysis program in NTP system

3 反应堆动力学模型

反应堆中子动力学模型旨在解决NTP推力调节等瞬态过程中存在一系列反应性引入以及反应性反馈下的堆芯功率计算问题。在动态问题中忽略堆芯内中子密度随空间分布可以引入点堆模型。

3.1 点堆模型动态方程

在时空动力学的基础上忽略空间效应即可获得简化的点堆模型,如式(6)所示,该模型的求解即为ODE的初值问题。由于点堆模型各组缓发中子的时间尺度差别较大,该常微分方程组具备较强的刚性。因此,采用后向欧拉(BE)以及Crank Nicholson(CN)等隐式求解方法可以采用更大的时间步长。

(6)

初值条件:n(t)|t=0=n0,Ci(t)|t=0=Ci0,i=1,2,…,I

其中:n为中子密度,cm-3;ρ为净反应性,$;Ci为先驱核浓度,cm-3;βi为第i组缓发中子份额;β为缓发中子总份额;λi为第i组衰变常数,s-1;Λ为中子每代时间,s;I为缓发中子组数。

3.2 点堆模型求解流程

线性ODE的隐式求解需要在各时间步内递进求解Ax=b,其计算流程如图5所示。其中BE与CN方法的区别体现在如何构造矩阵A和向量b上。在确定各时间步下净反应性ρ(t)后即可采用LU分解对线性方程组进行求解。

图5 点堆动力学方程隐式求解流程图Fig.5 Schematic diagram of implicit solution process for point reactor kinetics model

4 涡轮泵模型

涡轮泵是推进剂供应系统的核心部件,经过喷管再生冷却加热后的推进剂进入涡轮做功驱动离心泵工作以提供驱动压头。涡轮泵瞬态分析模型是系统瞬态分析的基础,涡轮泵模型主要包括离心泵模型、涡轮模型以及转轴模型3个部分。

4.1 离心泵模型

离心泵特性曲线由转矩、压头关于流量系数曲线组成。泵的特性曲线实际上反映的是其稳态性能,在准稳态假设下可认为泵在不同工况下的性能是沿着其特性曲线变化的。式(7)给出了基于流网中“流阻/感”模型而建立的泵流量与增压之间的关系:

(7)

(8)

其中:Δp为扬程,Pa;I为泵内流体的惯性,kg/m4;ρR为流体的参考密度,kg/m3;LR为泵中流体的长度,m;AR为泵的参考流通面积,m2。I与泵的额定参数和泵实际的时间常数相关。

基于上述泵模型可获得瞬态计算流程如图6所示,采用该显式计算方法需要知道当前时刻参数,以及下一时刻的边界参数进口参数、出口压力、泵转速,即可计算所需的泵瞬态参数。

图6 离心泵瞬态计算流程Fig.6 Schematic diagram of transient calculation process for centrifugal pump

4.2 涡轮模型

(9)

(10)

氢涡轮模型具体的计算流程为:根据涡轮进口参数、出口背压、转速以及涡轮工作特性,迭代获得涡轮内气体质量流量和涡轮转矩Tt,同时更新涡轮工作特性,如图7所示。

图7 涡轮模型计算流程图Fig.7 Schematic diagram of turbine model calculation process

4.3 转轴模型

涡轮模型中假设离心泵叶轮、涡轮动叶以及转轴为刚体并绕转轴中心定轴旋转。离心泵和涡轮模型输入参数均包括轴转速,输出参数均包括转矩。通过式(11)刚体定轴转动角动量方程即可将离心泵和涡轮进行联立瞬态求解。

(11)

其中:ω为涡轮泵轴转速,rad/s;Jtp为涡轮泵转动部分转动惯量,kg·m2;Tt为涡轮驱动转矩,N·m;Tp为离心泵阻力力矩,N·m;Tf为摩擦阻力力矩,N·m。

5 结论

本文针对NTP系统研发需求,提出了系统性能设计主要内容应包括稳态设计点性能分析与优化、稳态非设计点性能分析以及瞬态性能分析3个环节。为自主开发具备稳态、瞬态设计分析能力的一体化的NTP系统分析程序,提出了基于“流网-热网”的建模框架,建立了以常微分方程组为基础的系统数学模型。基于上述框架,给出了点堆动力学模型与涡轮泵瞬态模型的瞬态计算分析流程。后续将采用该框架进一步完善系统与各部件模型,并对典型系统方案开展性能分析工作。