CN106055850B - 一种获得偏离泡核沸腾型临界热流密度的方法 - Google Patents

Abstract

一种获得偏离泡核沸腾型临界热流密度的方法，步骤如下：1.设置一个初始热流密度q_m，为了获得用于临界热流密度机理模型所需的局部热工水力参数，基于均匀流模型建立管内两相热工水力瞬态求解模型，基于有限差分法来求解两相均匀流偏导，以及采用半隐式差分格式求解离散后的控制方程；2.当步骤1计算所得局部热工水力参数达到稳态时，将步骤1获得的局部热工水力参数作为已知参数赋给偏离泡核沸腾型临界热流密度机理模型，求解临界热流密度q_CHF；3.判断q_m与q_CHF是否一致，若没有收敛，则更新q_m的值，重复步骤1，直至二者相同；本发明方法具有较好的拓展性，对数据依赖性低，能够有效地预测临界热流密度，对反应堆安全分析具有重要意义。

Description

一种获得偏离泡核沸腾型临界热流密度的方法

技术领域

本发明属于反应堆安全分析技术领域，具体涉及一种获得偏离泡核沸腾型临界热流密度的方法。

技术背景

为保证堆芯流通通道的正常运行和完整性，需要对管内流动换热特性进行研究，临界热流密度是反应堆设计的主要热工安全准则之一。临界热流密度是核燃料元件表面发生传热恶化时的热通量，是冷却剂流动沸腾机理发生转变的结果。如果反应堆燃料元件表面发生沸腾临界，将会导致燃料元件表面温度过高从而造成加热壁面烧毁，放射性物质泄露，造成严重的运行事故，直接影响到反应堆的安全运行。因此准确预测临界热流密度对反应堆的设计和安全分析、对反应堆系统综合性能的提高具有重要意义。

目前，临界热流密度的预测方法较多采用适用于不同流体和不同结构流道的各种经验关系式以及临界热流密度查询表。经验关系式是最早使用的预测和研究临界热流密度的一种方法。至今，针对不同的热工水力条件、流型、管道结构尺寸等参数，国内外研究者们已经提出数百种临界热流密度预测经验关系式。但大多数经验关系式是基于一个比较小的数据库开发得到，因此大多数经验关系式适用范围有限，不能扩展至适用范围以外的工况。而临界热流密度查询表可在很宽的参数范围直接查表并插值的方法得到相应参数下的临界热流密度，但其只适用于单一流体，且外推到查询表提供的参数范围之外(例如高流速和高过冷度)时预测精度下降，较经验关系式需要花费更多的计算时间。此外，数据表格的形式使其不便于更新。

基于此，有必要提供一种可避免使用大量的实验数据拟合公式，降低对数据的依赖性，并能够适用于不同的流型和流体工质的更有效率的预测临界热流密度的方法。

发明内容

为了解决上述问题，本发明提供一种获得偏离泡核沸腾型临界热流密度的方法，可避免使用大量的实验数据拟合公式，降低对数据的依赖性，并能够适用于不同的流型和流体工质。

为了达到上述目的，本发明的采用如下技术方案：

一种获得偏离泡核沸腾型临界热流密度的方法，包括如下步骤：

步骤一：给定流道结构尺寸以及边界条件(进口温度、进口流量和出口压力)，设置一个初始热流密度q_m，采用均匀流模型模拟流道内两相瞬态流场，以此获得流通通道中的局部热工水力参数，即压力、流量及焓值，从而为偏离泡核沸腾型临界热流密度计算提供参数；

步骤二：当步骤一中计算所得局部热工水力参数达到稳态时，将步骤一得到的出口处局部热工水力参数作为已知参数赋给偏离泡核沸腾型临界热流密度机理模型，采用该机理模型计算得到临界热流密度q_CHF；；

步骤三：比较步骤一中的初始热流密度q_m与步骤二计算所得临界热流密度q_CHF，，若q_m大于q_CHF，则减小q_m，若q_m小于q_CHF，则增加q_m，采用变化后的热流密度重复步骤一的过程；若二者相同，则得到所求工况下的临界热流密度q_CHF；

步骤一所述的采用均匀流模型模拟流道内两相瞬态流场，采用有限差分法求解两相均匀流偏导，采用半隐式差分格式求解离散方程。

步骤二所述的偏离泡核沸腾型临界热流密度机理模型为微液层蒸干模型，其原理为：在壁面热流密度逐渐缓慢增加过程中，流体达到两相，壁面附近会出现一个拉长的汽块，它由多个小汽泡结合而成；长汽块下裹覆着一层很薄的液体层(微液层)，该长汽块会沿加热壁面移动，此时长汽块会阻碍主流液体对壁面的有效冷却；当长汽块下的微液层蒸干时，该点处的传热会剧烈恶化，从而导致壁温急剧升高，进而发生沸腾临界；求解机理模型的关键公式如下：

q_CHF＝δρ_fh_fgU_b/L_b 公式(1) 式中：

ρ_f——饱和液相密度/kg·m-³；

h_fg——汽化潜热/J·kg-¹；

δ——长汽块下微液层厚度/m；

U_b——长汽块移动速度/m·s-¹；

L_b——长汽块长度/m。

从式(1)中可知，微液层厚度(δ)、汽块移动速度(U_b)和汽块长度(L_b)是求解的关键参数；

其中汽块移动速度U_b通过求解汽块的轴向受力平衡来确定，轴向方向上汽块的受力包括浮力F_b和拖曳力F_d，公式如下：

F_b-F_d＝0 公式(2)

微液层厚度δ通过求解汽块的径向受力平衡来确定，径向方向上的汽块受力平衡包括侧面提升力F_r、蒸发力F_e以及壁面润滑力F_w：

-F_r+F_e+F_w＝0 公式(3)

汽块长度L_b采用Helmholtz临界波长，其公式如下所示：

式中：

ρ_l——液相密度/kg·m^-3；

ρ_g——汽相密度/kg·m^-3；

σ——表面张力/N·m^-1。

本发明和现有技术相比，具有如下优点：

1.采用均匀流模型模拟了两相瞬态流场，可较精确地计算模拟出临界沸腾发生点处的热工水力参数。

2.从微观上对汽泡受力平衡和运动进行研究，得到了较为精细的机理模型，同时也便于后续分析运动条件对模型的影响因素。

附图说明

图1为本发明方法的流程图。

图2为一维空间网格。

图3为汽块受力的示意图。

图4为机理模型计算流程图。

具体实施方式

下面结合附图和具体实施方式对本发明结构进行详细说明。

如图1所示，本发明一种获得偏离泡核沸腾型临界热流密度的方法，具体方法如下：

步骤一：给定流道结构尺寸以及边界条件(进口温度、进口流量和出口压力)，设置一个初始热流密度q_m，采用均匀流模型模拟流道内两相瞬态流场，以此获得流通通道中的局部热工水力参数，即压力、流量及焓值，从而为偏离泡核沸腾型临界热流密度计算提供参数；

步骤二：当步骤一中计算所得局部热工水力参数达到稳态时，将步骤一得到的出口处局部热工水力参数作为已知参数赋给偏离泡核沸腾型临界热流密度机理模型模块，计算得到临界热流密度q_CHF；

步骤三：比较步骤一中的初始热流密度q_m与步骤二计算所得临界热流密度q_CHF，若q_m大于q_CHF，则减小q_m，若q_m小于q_CHF，则增加q_m，采用变化后的热流密度重复步骤一的过程；若二者相同，则得到所求工况下的临界热流密度q_CHF。

下面进行详细说明：

步骤一中：两相瞬态流场采用均匀流模型模拟求解，下面给出具体的模型：

由于所研究的临界热流密度特性只需要提供外部的热工水力参数，对流场模拟的要求不高，因此采用改进的均匀流模型模拟流场。均匀流模型的基本假设为：

①两相具有相等的线速度，即u_g＝u_l＝u_H；

②两相之间处于热力平衡；

③使用合理确定的单相摩擦因数。

运用上述假定于两相流动基本方程，并与单相流体守恒方程类比，定义合适的混合物物性和摩擦因数计算式，可以得到以下守恒方程。

1)质量守恒方程：

式中：

ρ——两相混合物的密度/kg·m^-3；

u——两相混合物的流动速度/m·s^-1；

t——时间/s；

z——长度/m。

考虑两相流体的可压缩性，为了求解方程组，增加如下物性方程：

ρ＝ρ(p,h) 公式(6)

由式(6)求导可得：

式中：

h——两相混合物的焓值/J·kg^-1；

p——压力/Pa。

将式(7)代入质量守恒方程(5)，可化为：

2)动量守恒方程：

上式右边第二项表示摩擦压降，第三项表示重位压降。

3)能量守恒方程：

式中：

P_rh——加热周长/m；

A——流通面积/m²；

q_V——体积释热率/W·m^-3；

①若无内热源，即流体内部没有释热。上式右侧的q_V为零。

②忽略动能变化引起的能量变化，则上式左侧项中的为零。

因此，能量方程可简化为

将式(7)代入上式，可化为

上式中，

步骤一中：上述守恒方程，采用有限差分法求解偏导，采用半隐式差分格式求解离散方程。

在空间上，控制体划分采用交错网格。交错网格中，热力学变量 (压力p，焓值h，密度ρ，空泡份额α等)存储在控制体中心(图2 中虚线位置)，而流体的速度存放在控制体边界处(图2中实线位置)。

在时间上采用半隐式差分格式进行离散。为了加快计算速度，除了连续方程和能量方程中对流项的速度，以及动量方程中的压力梯度项采用隐式格式以外，其他项均采用显示格式。

各守恒方程的离散形式如下：

质量守恒方程：

式中上标n表示旧时层，n+1表示新时层；下标i表示控制体编号。

动量守恒方程：

其中：

时间步长的选取可以违反(u+a)Δt/Δz<1的波速限制，但是要满足 uΔt/Δz<1的流速限制。

能量守恒方程：

其中：

其中P_rh——加热周长/m。

交错网格以及有限差分法的优点是压力梯度项与速度项自然的耦合在一起，并且如果满足Courant法则uΔt/Δz<1则形成一个对角优势矩阵，可保证求解算法的稳定性。所有网格的方程最后可以组成一个压力的三对角矩阵，它可以利用直接消去法求解，也可以利用 Gauss-Seidel迭代法求解。总体计算思路是：①用动量方程求解出速度；②利用上一时层求解压力方程；③用解出的压力回代求解速度和焓值。

步骤二中：微液层蒸干机理的具体模型如下：

微液层蒸干模型假设加热壁面附近产生的小汽泡结合形成大汽块，在汽块下存在非常薄的液相层，称为微液层。汽块移动过程中，当汽块下的液体全部蒸发烧干时，该点处的加热壁面被单相蒸汽覆盖从而导致传热恶化，进而导致沸腾临界发生。求解机理模型的关键公式如下：

q_CHF＝δρ_fh_fgU_b/L_b 公式(19)

式中：

ρ_f——饱和液相密度/kg·m^-3；

h_fg——汽化潜热/J·kg^-1；

δ——汽块下微液层厚度/m；

U_b——汽块移动速度/m·s^-1；

L_b——汽块长度/m。

从(14)式中可知，微液层厚度(δ)、汽块移动速度(U_b)和汽块长度(L_b)是求解的关键参数。

如图3所示，上述三个关键参数可以通过分析汽块受力获得，其中汽块移动速度U_b通过求解汽块的轴向受力平衡来确定。轴向方向上汽块的受力包括浮力F_b和拖曳力F_d，公式如下：

F_b-F_d＝0 公式(20)

其中

式中：

D_b——汽块当量直径/m；

ρ_l——液相的密度/kg·m^-3；

ρ_g——汽相的密度/kg· m^-3；

g——重力加速度/m s^-2；

t——时间/s；

C_d——拖拽系数；

U_bl——汽块中心线所处的径向位置处的主流速度/m s^-1。

联立式(20)-(22)，汽块速度U_b可表达如下：

上式中，拖曳力系数采用如下公式：

汽块当量直径D_b可采用Levy模型：

式中：

σ——表面张力/N·m^-1；

D——管道直径/m；

τ_w——壁面剪切力/kg·m^-1·s^-2。

采用Karman速度分布方程确定当地流速：

式中，

近壁面区的汽块一般处在缓冲区的范围，即对应于Karman速度分布方程的第二个关系式。所以为简化计算，可选用第二个关系式计算当地流速分布，该关系式可以改写为：

因此距离壁面y＝δ+D_B/2处的当地流速可表达为：

汽块长度L_b采用Helmholtz临界波长，其公式如下所示：

微液层厚度δ通过求解汽块的径向受力平衡来确定。径向方向上的汽块受力包括平衡侧面提升力F_r、蒸发力F_e以及壁面润滑力F_w：

-F_r+F_e+F_w＝0 公式(30)

蒸发力是由微液层中液相蒸发产生的小汽泡进入汽块，对汽块冲击而产生的，这个力将推动汽块向远离加热壁面。侧面提升力是由汽块和主流区存在速度差异，使汽块在运动过程中发生旋转而产生，提升力的大小还与汽块在径向方向上的速度梯度有关，它会推动汽块靠近加热壁面。

沸腾临界将要发生时，由微液层蒸发而产生的蒸发力为：

垂直管道内两相湍流流动中的侧面提升力的表达式：

上式中C_l为提升力系数：

其中式(33)中，α为临界沸腾发生点处的空泡份额，式(29)中的速度由式(22)计算获得。

因此，侧面提升力可表达如下：

此外还需要考虑壁面润滑力F_w，壁面润滑力将非常靠近壁面处的汽泡推向中心区，而对离壁面距离较远的汽泡作用力则非常小，其表达式如下：

式中，C_w是壁面润滑系数，如下所示：

上式中，C₁和C₂为经验常数，y_w表示离加热壁面的垂直距离，式 (36)在y_w<5D_b范围内有效。

联立式(30)，(31)，(32)和(35)，微液层厚度δ可表达如下：

求出δ以后，可以通过式(28)计算得到U_bl，然后分别用式(23) 和(29)计算得到U_b和L_b，再将新的U_b和L_b代入式(38)计算得到新的δ，通过迭代计算直到以上几个参数值达到收敛。最后，将迭代计算得到的δ，U_b和L_b代入式(19)来计算临界热流密度值。机理模型的流程图如图4所示。

Claims (3)

1.一种获得偏离泡核沸腾型临界热流密度的方法，其特征在于：包括如下步骤：

步骤一：给定流道结构尺寸以及边界条件即进口温度、进口流量和出口压力，设置一个初始热流密度q_m，采用均匀流模型模拟流道内两相瞬态流场，以此获得流通通道中的局部热工水力参数，即压力、流量及焓值，从而为偏离泡核沸腾型临界热流密度计算提供参数；

步骤二：当步骤一中计算所得局部热工水力参数达到稳态时，将步骤一得到的出口处局部热工水力参数作为已知参数赋给偏离泡核沸腾型临界热流密度机理模型，采用该机理模型计算得到临界热流密度q_CHF；

步骤三：比较步骤一中的初始热流密度q_m与步骤二计算所得临界热流密度q_CHF，若q_m大于q_CHF，则减小q_m，若q_m小于q_CHF，则增加q_m，采用变化后的热流密度重复步骤一的过程；若二者相同，则得到所求工况下的临界热流密度q_CHF。

2.根据权利要求1所述的一种获得偏离泡核沸腾型临界热流密度的方法，其特征在于：步骤一所述的采用均匀流模型模拟流道内两相瞬态流场，采用有限差分法求解两相均匀流偏导，采用半隐式差分格式求解离散方程。

3.根据权利要求1所述的一种获得偏离泡核沸腾型临界热流密度的方法，其特征在于：步骤二所述的偏离泡核沸腾型临界热流密度机理模型为微液层蒸干模型，其原理为：在壁面热流密度逐渐缓慢增加过程中，流体达到两相，壁面附近会出现一个拉长的汽块，它由多个小汽泡结合而成；长汽块下裹覆着一层很薄的液体层即微液层，该长汽块会沿加热壁面移动，此时长汽块会阻碍主流液体对壁面的有效冷却；当长汽块下的微液层蒸干时，该点处的传热会剧烈恶化，从而导致壁温急剧升高，进而发生沸腾临界；求解机理模型的关键公式如下：

q_CHF＝δρ_fh_fgU_b/L_b 公式(1)

式中：

ρ_f——饱和液相密度/kg·m^-3；

h_fg——汽化潜热/J·kg^-1；

δ——长汽块下微液层厚度/m；

U_b——长汽块移动速度/m·s^-1；

L_b——长汽块长度/m；

从式(1)中可知，微液层厚度δ、汽块移动速度U_b和汽块长度L_b是求解的关键参数；

其中汽块移动速度U_b通过求解汽块的轴向受力平衡来确定，轴向方向上汽块的受力包括浮力F_b和拖曳力F_d，公式如下：

F_b-F_d＝0 公式(2)

微液层厚度δ通过求解汽块的径向受力平衡来确定，径向方向上的汽块受力平衡包括侧面提升力F_r、蒸发力F_e以及壁面润滑力F_w：

-F_r+F_e+F_w＝0 公式(3)

汽块长度L_b采用Helmholtz临界波长，其公式如下所示：

式中：

ρ_l——液相密度/kg·m^-3；

ρ_g——汽相密度/kg·m^-3；

σ——表面张力/N·m^-1。