基于弱可压缩模型的直流蒸汽发生器一维两相流计算方法

1.本发明涉及核反应堆热工水力技术领域，具体涉及一种基于弱可压缩模型的直流蒸汽发生器一维两相流计算方法。


背景技术：

2.钠冷快堆(sfr)作为六种第四代核能系统之一，是目前世界研究范围最广泛、技术最成熟、运行经验最丰富的快堆堆型。钠冷快堆一般采用直流式蒸汽发生器(otsg)，蒸汽发生器主要分为蒸发器和过热器两个模块，是钠冷快堆系统的核心设备，因此，对钠-水蒸汽发生器的热工水力理论研究和数值分析具有重大意义。
3.在给水泵断电导致的失给水事故中，给水流量、压力迅速降低，蒸发器内压力随之降低，蒸发器内出现了剧烈的失压闪蒸、流量骤降等强瞬态现象，导致水侧流量沿轴向的变化巨大。水侧发生了复杂了两相流动过程，且可压缩效应明显，基于不可压缩的均相流模型和两流体模型的系统级安全分析程序并不能很好地捕捉到该现象，且计算值与实验值有较大偏差。
4.系统级安全分析程序需要满足计算结果准确、求解过程稳定且快速的特点。在计算失给水事故过程中，基于不可压缩模型的程序由于模型过度简化难以满足计算的准确性，基于完全可压缩模型由于压力导数误差大数值稳定性差，且求解收敛速度远低于不可压缩。因此，不可压缩模型和完全可压缩模型都不适用于失给水工况的求解，亟需一种既能满足计算的准确性，又可以保证求解稳定性和快速性的可压缩计算方法。


技术实现要素：

5.为了克服上述现有技术存在的问题，本发明的目的在于提供一种基于弱可压缩模型的直流蒸汽发生器一维两相流计算方法，该方法是一种能够计算直流蒸汽发生器中可压缩效应明显的瞬态过程的方法，能够实现一维蒸汽发生器水侧两相热工水力问题的求解，同时保证程序计算的准确性、求解的稳定性和快速性。
6.为了达到上述目的，本发明采用如下技术方案：
7.一种基于弱可压缩模型的直流蒸汽发生器一维两相流计算方法，包括如下步骤：
8.步骤一：采用交错网格的方法，对待求解直流蒸汽发生器的蒸发器和过热器进行一维网格节点划分，同时根据蒸汽发生器几何模型确定不同区域网格连接方式；
9.步骤二：获取当前时间步直流蒸汽发生器的入口流量w
in
、入口比焓h
in
和出口压力p
out
作为边界条件，输入当前时间步的焓场h
t
和上一时间步的密度场ρ
t-δt
、压力场p
t-δt
和流量场w
t-δt
，根据h
t
和p
t-δt
计算当前时间步的物性参数；
10.步骤三：计算蒸发器和过热器的流量分布：直流蒸汽发生器内水侧流动属于低马赫数流动，基于此提出了针对直流蒸汽发生器的弱可压缩模型，即假设压力波传递远快于热响应，对流量导数、压力导数显式处理，根据质量守恒方程，得到流量空间分布如式(1)所示，对流量分布wi空间离散，对密度导数显式计算，得到相邻节点流量关系如式(2)所示，从
而根据边界条件获取直流蒸汽发生器轴向流量分布；
[0011][0012][0013]
式中：ρ——流体密度/kg
·
m-3
；w——流体质量流量/kg
·
s-1
；a——流通面积/m2；z——坐标/m；t——时间/s；δt——时间步长/s；i——控制体节点标号；w
it
——当前时间步第i个控制体的质量流量/kg
·
s-1
，——当前时间步第i-1控制体的质量流量/kg
·
s-1
；
[0014]
步骤四：显式计算流量时间导数：根据上一时间步的流量和步骤三得到的流量计算流量时间导数，如式(3)所示；
[0015][0016]
式中w
t
——当前时间步质量流量/kg
·
s-1
；w
t-δt
——上一时间步质量流量/kg
·
s-1
；
[0017]
步骤五：计算蒸发器和过热器的压力分布：根据直流蒸汽发生器弱可压缩模型假设，对压力导数显式处理，根据动量守恒方程，得到压力空间导数如式(4)所示，对压力空间导数空间离散得到相邻节点压力关系如式(5)所示，从而根据压力边界条件计算得到直流蒸汽发生器当前时间步压力p
t
；
[0018][0019][0020]
式中：p——流体压力/pa；g——重力加速度/m
·
s-2
；de——等效直径/m；f——摩擦阻力系数；p
i-1
——第i-1个控制体流体压力/pa；pi——第i个控制体流体压力/pa；δz——控制体轴向高度/m；
[0021]
步骤六：显式计算压力时间导数：根据上一时间步的压力和步骤五得到的压力计算压力时间导数，如式(6)所示；
[0022][0023]
式中p
t
——当前时间步压力/pa；p
t-δt
——上一时间步压力/pa；
[0024]
步骤七：更新物性：利用当前时间步的焓场h
t
和步骤五计算的当前时间步压力p
t
更新物性，用于比焓时间导数计算；
[0025]
步骤八：隐式计算比焓时间导数：根据直流蒸汽发生器弱可压缩模型假设，热响应慢于压力波，需要对比焓方程隐式求解，由能量守恒方程得到比焓时间导数方程，如式(7)所示：
[0026][0027]
式中：h——流体比焓/j
·
kg-1
；hf——水侧对流换热系数/w
·
m-2
·
k-1
；tw——壁面温度/k；t——流体温度/k；
[0028]
步骤九：采用函数迭代法隐式求解步骤八构建的比焓时间导数方程，如满足收敛性要求，则时间步推进，否则更新比焓，返回步骤一在当前时间步继续迭代计算，直至推进到结束时间。
[0029]
和现有技术相比较，本发明具备如下优点：
[0030]
1、本发明提出了一种针对直流蒸汽发生器的弱可压缩模型，可以实现直流蒸汽发生器水侧的可压缩计算。目前国内钠冷快堆蒸汽发生器热工水力计算程序大多数采用了不可压缩模型，在计算类似失给水这种快速降压、降流量的过程中，由于其不可压缩假设难以捕捉到失压闪蒸现象。本发明提出了一种基于弱可压缩模型的直流蒸汽发生器的计算方法，假设热响应慢于压力波，对流量导数和压力导数进行了显式处理，实现了直流蒸汽发生器变化剧烈的瞬态工况的可压缩效应计算；
[0031]
2、本发明提出的弱可压缩模型数值稳定性高。在两相流程序中需要采用经验关系式计算压降，此类经验关系式由稳态数据拟合、修正而来，适用于稳态或变化不太剧烈的瞬态工况计算，对于失给水此类强瞬态工况会有较大偏差，直接考虑流量导数和压力导导致误差较大，造成数值不稳定。本发明的计算方法对流量和压力导数显式处理，通过边界条件和显式导数对流量和压力分布进行了重新计算，有效降低了计算误差，提高了数值稳定性；
[0032]
3、本发明提出的弱可压缩模型求解速度快。对于蒸汽发生器热工水力程序，若采用完全可压缩计算计算需要对三个物理量(压力、流量、比焓)进行隐式求解，三个物理量收敛速度不一致造成了方程组刚性较大，本发明采用的弱可压缩模型只对比焓进行隐式求解，一方面降低了方程组的刚性，另一方面减少了求解方程组的数量，可以将计算速度提升三倍以上；
附图说明
[0033]
图1为本发明基于弱可压缩模型的直流蒸汽发生器的两相流算法框图；
[0034]
图2为直流蒸汽发生器建模与节点划分图；
[0035]
其中
[0036]
1、蒸发器网格，2、过热器网格，3、蒸发器入口，4、蒸发器出口，5、过热器入口，6、过热器出口。
具体实施方式
[0037]
下面结合图1和图2，对本发明提出的基于弱可压缩模型的直流蒸汽发生器两相流计算方法进行详细说明。
[0038]
本发明一种基于弱可压缩模型的直流蒸汽发生器一维两相流计算方法，包括如下步骤：
[0039]
步骤一：采用交错网格的方法，如附图2所示，对待求解直流蒸汽发生器的蒸发器1和过热器2进行一维网格划分，控制体中心(实心点)储存压力和比焓，控制体边界(空心点)
存储流量；同时根据蒸汽发生器几何模型和流动方向确定不同区域网格连接方式，蒸发器二次侧水从蒸发器入口3流入，蒸发器出口4流出，通过接管与过热器入口5连接，最终蒸汽从过热器出口流出6；
[0040]
步骤二：获取当前时间步直流蒸汽发生器的入口流量w
in
、入口比焓h
in
和出口压力p
out
作为边界条件，同时将其传递给对应节点，输入当前时间步的焓场h
t
和上一时间步的密度场ρ
t-δt
、压力场p
t-δt
和流量场w
t-δt
，根据h
t
和p
t-δt
计算当前时间步的动力粘度、比热、热导率等物性参数；
[0041]
步骤三：计算蒸发器和过热器的流量分布：直流蒸汽发生器内水侧流动属于低马赫数流动，基于此提出了针对直流蒸汽发生器的弱可压缩模型，即假设压力波传递远快于热响应，压力波传递迅速，在每个时刻都能建立准稳态压力泊松方程，流量时间导数不再采用完全可压缩下的隐式计算，而对流量导数显式处理，得到流量分布，根据质量守恒方程，得到流量空间分布如式(1)所示，对流量空间离散，对密度导数显式计算，得到相邻节点流量关系如式(2)所示，从而根据入口流量边界条件依次计算从第2到第n+1节点的流量(式(3))，最终获取蒸汽发生器轴向流量分布；
[0042][0043][0044][0045]
式中：ρ——流体密度/kg
·
m-3
；w——流体质量流量/kg
·
s-1
；a——流通面积/m2；z——坐标/m；δz——控制体轴向高度/m；t——时间/s；δt——时间步长/s；i——控制体节点标号；w
it
——当前时间步第i个控制体的质量流量/kg
·
s-1
，——当前时间步第i-1控制体的质量流量/kg
·
s-1
；
[0046]
步骤四：显式计算流量时间导数，用于压力分布的计算：根据上一时间步的压力和步骤三得到的当前时间步的流量计算流量时间导数，如式(4)所示；
[0047][0048]
式中w
t
——当前时间步质量流量/kg
·
s-1
；w
t-δt
——上一时间步质量流量/kg
·
s-1
；
[0049]
步骤五：计算蒸发器和过热器的压力分布：根据提出的直流蒸汽发生器弱可压缩模型，假设压力波传递迅速，在每个时刻都能建立准稳态压力泊松方程，压力时间导数不再采用完全可压缩下的隐式计算，而对流量导数显式处理，得到压力分布，根据动量守恒方程，得到压力空间导数如式(5)式所示，对其空间离散得到相邻节点压力关系如式(6)所示，从而根据出口压力边界依次计算n+1个节点到第2个节点的压力(式(7))，最终得到直流蒸汽发生器当前时间步压力p
t
；
[0050]
[0051][0052][0053]
式中：p——流体压力/pa；g——重力加速度/m
·
s-2
；de——等效直径/m；f——摩擦阻力系数；p
i-1——第i-1个控制体流体压力/pa；pi——第i个控制体流体压力/pa；
[0054]
步骤六：显式计算压力时间导数，用于比焓导数的计算：根据上一时间步长的压力和步骤五得到的压力计算压力时间导数，如式(8)所示；
[0055][0056]
式中p
t
——当前时间步压力/pa；p
t-δt
——上一时间步压力/pa；
[0057]
步骤七：更新物性：压力分布和流量分布计算中使用了上时间步的压力，由于热响应慢于压力波传递，所以利用当前时间步的焓场h
t
和步骤五计算的当前时间步压力p
t
更新物性，用于比焓时间导数计算；
[0058]
步骤八：隐式计算比焓时间导数：根据能量守恒方程得到比焓时间导数方程如式(9)所示，其中密度时间导数采用步骤三中的计算值，压力时间导数采用步骤六中的计算值；
[0059][0060]
式中：h——流体比焓/j
·
kg-1
；hf——水侧对流换热系数/w
·
m-2
·
k-1
；tw——壁面温度/k；t——流体温度/k；
[0061]
步骤九：采用函数迭代方法隐式求解步骤八构建的比焓时间导数方程，如满足收敛性要求，即比焓时间导数方程求解残差小于设定值，则时间步推进，否则更新比焓，返回步骤一在当前时间步继续迭代计算，直至推进到结束时间。

技术特征：
1.一种基于弱可压缩模型的直流蒸汽发生器一维两相流计算方法，其特征在于，包括如下步骤：步骤一：采用交错网格的方法，对待求解直流蒸汽发生器的蒸发器和过热器进行一维网格节点划分，同时根据蒸汽发生器几何模型确定不同区域网格连接方式；步骤二：获取当前时间步直流蒸汽发生器的入口流量w
in
、入口比焓h
in
和出口压力p
out
作为边界条件，输入当前时间步的焓场h
t
和上一时间步的密度场ρ
t-δt
、压力场p
t-δt
和流量场w
t-δt
，根据h
t
和p
t-δt
计算当前时间步的物性参数；步骤三：计算蒸发器和过热器的流量分布：直流蒸汽发生器内水侧流动属于低马赫数流动，基于此提出了针对直流蒸汽发生器的弱可压缩模型，即假设压力波传递远快于热响应，对流量导数、压力导数显式处理，根据质量守恒方程，得到流量空间分布如式(1)所示，对流量分布w
i
空间离散，对密度导数显式计算，得到相邻节点流量关系如式(2)所示，从而根据边界条件获取直流蒸汽发生器轴向流量分布；根据边界条件获取直流蒸汽发生器轴向流量分布；式中：ρ——流体密度/kg
·
m-3
；w——流体质量流量/kg
·
s-1
；a——流通面积/m2；z——坐标/m；t——时间/s；δt——时间步长/s；i——控制体节点标号；——当前时间步第i个控制体的质量流量/kg
·
s-1
，——当前时间步第i-1控制体的质量流量/kg
·
s-1
；步骤四：显式计算流量时间导数：根据上一时间步的流量和步骤三得到的流量计算流量时间导数，如式(3)所示；式中w
t
——当前时间步质量流量/kg
·
s-1
；w
t-δt
——上一时间步质量流量/kg
·
s-1
；步骤五：计算蒸发器和过热器的压力分布：根据直流蒸汽发生器弱可压缩模型假设，对压力导数显式处理，根据动量守恒方程，得到压力空间导数如式(4)所示，对压力空间导数空间离散得到相邻节点压力关系如式(5)所示，从而根据压力边界条件计算得到直流蒸汽发生器当前时间步压力p
t
；；式中：p——流体压力/pa；g——重力加速度/m
·
s-2
；d
e
——等效直径/m；f——摩擦阻力系数；p
i-1
——第i-1个控制体流体压力/pa；p
i
——第i个控制体流体压力/pa；δz——控制体轴向高度/m；步骤六：显式计算压力时间导数：根据上一时间步的压力和步骤五得到的压力计算压力时间导数，如式(6)所示；
式中p
t
——当前时间步压力/pa；p
t-δt
——上一时间步压力/pa；步骤七：更新物性：利用当前时间步的焓场h
t
和步骤五计算的当前时间步压力p
t
更新物性，用于比焓时间导数计算；步骤八：隐式计算比焓时间导数：根据直流蒸汽发生器弱可压缩模型假设，热响应慢于压力波，需要对比焓方程隐式求解，由能量守恒方程得到比焓时间导数方程，如式(7)所示：式中：h——流体比焓/j
·
kg-1
；h
f
——水侧对流换热系数/w
·
m-2
·
k-1
；t
w
——壁面温度/k；t——流体温度/k；步骤九：采用函数迭代法隐式求解步骤八构建的比焓时间导数方程，如满足收敛性要求，则时间步推进，否则更新比焓，返回步骤一在当前时间步继续迭代计算，直至推进到结束时间。2.根据权利要求1所述的一种基于弱可压缩模型的直流蒸汽发生器一维两相流计算方法，其特征在于，步骤九所述的收敛性要求为比焓时间导数方程求解残差小于设定值。

技术总结
本发明公开了一种基于弱可压缩模型的直流蒸汽发生器一维两相流计算方法，旨在通过对物理过程合理假设和求解过程适当简化，实现直流蒸汽发生器可压缩效应明显的强瞬态工况的准确、稳定、快速计算求解。具体步骤如下：1、采用交错网格的方法对蒸汽发生器进行网格节点划分和连接；2、获取边界条件和输入旧时间步物理场，计算物性参数；3、计算蒸汽发生器流量分布；4、显式计算流量时间导数；5、计算蒸汽发生器压力分布；6、显式计算压力时间导数；7、用当前时间步压力更新物性；8、隐式计算比焓时间导数；9、隐式求解焓方程，收敛后时间推进，直至结束时间。束时间。束时间。


技术研发人员：张大林 冯振宇 李新宇 王博 刘雅鹏 田文喜 秋穗正 苏光辉
受保护的技术使用者：西安交通大学
技术研发日：2023.06.29
技术公布日：2023/9/25