CN112699620B - 基于计算流体力学的反应堆堆芯热工水力特性分析方法 - Google Patents
基于计算流体力学的反应堆堆芯热工水力特性分析方法 Download PDFInfo
- Publication number
- CN112699620B CN112699620B CN202110026314.9A CN202110026314A CN112699620B CN 112699620 B CN112699620 B CN 112699620B CN 202110026314 A CN202110026314 A CN 202110026314A CN 112699620 B CN112699620 B CN 112699620B
- Authority
- CN
- China
- Prior art keywords
- fuel rod
- node
- temperature
- heat exchange
- reactor core
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
- G06F30/28—Design optimisation, verification or simulation using fluid dynamics, e.g. using Navier-Stokes equations or computational fluid dynamics [CFD]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/10—Geometric CAD
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2113/00—Details relating to the application field
- G06F2113/08—Fluids
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
- G06F2119/08—Thermal analysis or thermal optimisation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
- G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Mathematical Physics (AREA)
- Pure & Applied Mathematics (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Computational Mathematics (AREA)
- Geometry (AREA)
- General Engineering & Computer Science (AREA)
- Computer Hardware Design (AREA)
- Evolutionary Computation (AREA)
- Algebra (AREA)
- Data Mining & Analysis (AREA)
- Fluid Mechanics (AREA)
- Computing Systems (AREA)
- Operations Research (AREA)
- Databases & Information Systems (AREA)
- Software Systems (AREA)
- Monitoring And Testing Of Nuclear Reactors (AREA)
Abstract
本发明公开了一种基于计算流体力学的反应堆堆芯热工水力特性分析方法,该方法主要内容如下:1、建立冷却剂流动换热特性分析模型;2、建立燃料棒热工特性分析模型;3、建立流固耦合换热模型;4、建立流固耦合求解方法。本方法通过对反应堆堆芯内复杂棒束结构实施简化处理,可以精确计算反应堆堆芯燃料棒温度分布和冷却剂流场和温度场,同时显著降低计算资源消耗,具有建模简单、计算速度快的特点,为在工程实际中运用计算流体力学数值模拟方法对反应堆全堆芯进行热工水力特性分析提供了一种新的计算方法。
Description
技术领域
本发明属于反应堆安全分析技术领域,具体涉及到正常运行及事故工况的反应堆功率变化条件下,反应堆堆芯热工水力特性分析方法。
背景技术
反应堆堆芯是核电厂一回路***的核心部件,而堆芯燃料包壳是防止反应堆内的放射性物质向外泄露的第一道屏障,因此反应堆堆芯结构的完整性——尤其是燃料包壳的完整性对反应堆安全来说至关重要,在反应堆堆芯热工水力分析研究中,堆芯燃料棒及冷却剂温度分布一直是研究的重点。而在堆芯中核燃料及控制装置以棒束结构分布,造成堆芯内流动及传热特性的复杂化,对堆芯建模及热工水力数值模拟分析造成了较大的困难。
目前应用于反应堆热工水力分析程序中,对燃料棒温度分布计算的方法一般有三种:1.一维***程序——如RELAP5程序、TRACE程序通过简单的热构件计算反应堆燃料的平均温度分布;2.子通道分析程序——如COBRA程序对燃料棒逐个建模,并通过数值求解得到燃料棒温度分布;3.计算流体动力学(CFD)模拟——如ANSYS FLUENT软件一般采用对燃料棒进行建模和网格划分,再逐个燃料棒求解温度分布。其中,一维***分析程序无法获得堆芯燃料棒准确的温度分布,因此对功率分布不均匀的瞬态计算过程模拟效果很差;子通道程序建模复杂,而且在全堆模拟中,边界条件需要人为指定,对反应堆堆芯入口不规律瞬变过程模拟能力不足;CFD软件对燃料棒精确建模工作量大,网格量大,导致计算速度慢。
发明内容
本发明要解决的技术问题在于,提供一种低计算资源消耗、高计算精度的反应堆堆芯热工水力特性分析方法。
本发明解决该技术问题所采用技术方案为:提出一种基于计算流体力学的反应堆堆芯热工水力特性分析方法。针对反应堆堆芯内燃料组件复杂的棒束结构造成堆芯内部流动及传热特性复杂化的问题,本发明基于计算流体力学,建立冷却剂流动换热特性分析模型、燃料棒热工特性分析模型,建立流固耦合换热模型,通过流固耦合求解方法计算获得堆芯燃料棒温度分布和冷却剂流场和温度场,实现反应堆堆芯热工水力特性分析。
一种基于计算流体力学的反应堆堆芯热工水力特性分析方法,其特征在于:具体包括如下步骤:
步骤1:建立冷却剂流动换热特性分析模型
(1)建立反应堆堆芯的计算模型
建立反应堆堆芯的几何模型,将反应堆堆芯的几何模型导入网格生成软件中划分控制体,形成反应堆堆芯的网格模型,将反应堆堆芯的网格模型导入计算流体力学软件中,形成反应堆堆芯的计算模型;
(2)建立冷却剂流动换热特性分析模型
为了将反应堆堆芯内复杂棒束结构的影响简化为对冷却剂的流动阻力,提出如下假设:①堆芯的流动阻力使用包含粘性损失、惯性损失的压降描述:
式中:是沿m坐标方向的动量方程中的压降,m和n是x,y,z坐标方向,Dmn和Cmn分别是粘性损失项和惯性损失项的系数矩阵,μ是冷却剂的动力粘性系数/kg·m-1·s-1,ρ是冷却剂的密度/kg·m-3,|v|是速度的模/m·s-1,vn是沿n坐标方向的速度/m·s-1;②计算域中燃料棒与定位格架固体部件体积不变,位置无移动;③孔隙率是各向同性的;基于上述假设,对流体控制方程中的动量方程和能量方程中添加源项:在动量方程中添加包含粘性损失、惯性损失的压降作为动量源项,在能量方程中添加燃料棒表面的对流换热量作为能量源项;从而建立了冷却剂流动换热特性分析模型;
步骤2:建立燃料棒热工特性分析模型
燃料棒热工特性分析模型中考虑燃料棒表面对流换热和燃料棒内部导热,包括表面对流换热模型和内部导热模型,分别建立表面对流换热模型和内部导热模型,步骤如下:
(1)建立表面对流换热模型
提出如下假设:①反应堆燃料棒及控制棒结构均匀地分布在堆芯中;②单个控制体内燃料棒几何结构不发生剧烈变化;③单个控制体内所有燃料棒具有相同的热工状态;④单个控制体内所有燃料棒具有相同的几何特性且已知;
基于上述假设,建立反应堆堆芯内的换热面积与反应堆堆芯体积的线性关系式:
式中:AT,cell为控制体内的换热面积/m2,Vcell为控制体体积/m3,AT,core为反应堆堆芯内的换热面积/m2,Vcore为反应堆堆芯体积/m3;
将该线性关系式应用于步骤1中建立的反应堆堆芯的计算模型中,用控制体的体积表示控制体内的换热面积;最终将控制体内的换热面积代入到表面对流换热公式中,计算得到燃料棒表面的对流换热量:
QT=hTAcell(Tr-Tf)
式中:QT为燃料棒表面的对流换热量/W,hT为燃料棒表面的对流换热系数/W·m-2·K-1,Acell为控制体内的换热面积/m2,Tf为流体主流温度/K,Tr为燃料棒表面温度/K;
通过上述步骤2的步骤(1)建立了表面对流换热模型;
(2)内部导热模型建立步骤
根据燃料棒的结构特点将燃料棒沿径向划分N个计算节点,其中外侧两个计算节点分别位于包壳内外两侧边界,其它计算节点沿燃料棒芯块中心向外侧依次布置,且有一个计算节点在芯块边界上,在包壳内外两侧各有一个计算节点;忽略燃料棒的轴向导热,沿燃料棒径向建立温度方程:
式中:ρi为计算节点i处材料密度/kg·m-3,CP,i为计算节点i处材料比热容/kJ·kg-1·K-1,Vi为计算节点i的等效控制体体积/m3,Ti为计算节点i处的温度/K,Qi-1,i为从计算节点i-1传导到计算节点i的热量/kW,Qi+1,i为计算节点i+1传导到计算节点i的热量/kW,Qi为计算节点i处单位体积释热率/kW·m-3;
燃料棒内部计算节点间传导的热量通过计算节点间温差得到:
Qi-1,i=Ki-1,i(Ti-1-Ti)
Qi+1,i=Ki+1,i(Ti+1-Ti)
式中:Ti-1为计算节点i-1处的温度/K,Ti+1为计算节点i+1处的温度/K,Ki-1,i为计算节点i-1和i之间的导热系数/kW·K-1,等同于Ki,i-1,Ki+1,i为计算节点i+1和i之间的导热系数/kW·K-1,等同于Ki,i+1,均为热导率和燃料棒几何参数的函数:
式中:ki是计算节点i处的热导率/kW·m-1·K-1,ki-1是计算节点i-1处的热导率/kW·m-1·K-1,ki+1是计算节点i+1处的热导率/kW·m-1·K-1,ri是计算节点i处的芯块直径/m,ri-1是计算节点i处的芯块直径/m,ΔXj为控制体轴向高度/m,为计算节点i的体积平均半径/m,为计算节点i-1的体积平均半径/m,为计算节点i+1的体积平均半径/m,其表达式为:
根据燃料棒的导热特点,建立如下边界条件:①燃料棒中心为对称边界,因此对于计算节点1的控制方程有K0,1=Ki-1,i=0;②在燃料芯块与包壳之间存在间隙导热,即在i=N-1处Qi,i-1=Hgap(Ti-Ti-1),其中Hgap为考虑换热面积的间隙导热率/kW·K-1,Qi,i-1为从计算节点i传导到计算节点i-1的热量/kW,Ti为计算节点i处的温度/K,Ti-1为计算节点i-1处的温度/K;③在燃料棒外表面,即在i=N处,Ki+1,i=H,Ti+1=Tf,则有Qi+1,i=H(Tf-Ti),其中,Ki+1,i为计算节点i+1和i之间的导热系数/kW·K-1,H为考虑换热面积的燃料棒表面的对流换热系数/kW·K-1,Qi+1,i为计算节点i+1传导到计算节点i的热量/kW,Tf为流体主流温度/K,Ti+1为计算节点i+1处的温度/K,Ti为计算节点i处的温度/K;
将温度方程在各个计算节点上化为差分方程:
式中:Ti n为计算节点i处在第n个时间步下的温度/K,Ti n-1为计算节点i处在第n-1个时间步下的温度/K,Ti n+1为计算节点i处在第n+1个时间步下的温度/K;
将差分方程变形为:
联立各个计算节点的差分方程形成方程组,其中方程组的系数矩阵为三对角矩阵,通过高斯迭代法求解得到燃料棒内各计算节点温度;
通过上述步骤2的步骤(2)建立了内部导热模型;
步骤3:建立流固耦合换热模型
冷却剂与燃料棒之间的换热系数由Dittus-Bolter关系式计算:
Nu=0.023Re0.8Pr0.4
式中:Nu为努塞尔数,Re为雷诺数,Pr为普朗特数;Nu和Re的表达式为:
式中:k为冷却剂热导率/W·m-1·K-1,L为特征长度/m,ρ为冷却剂密度/kg·m-3,μ为冷却剂的动力粘性系数/kg·m-1·s-1,U为冷却剂流速/m·s-1;
由于k、ρ、μ、U均为冷却剂热物性和流动状态的参数,因此冷却剂流场和温度场决定了燃料棒表面的对流换热系数,进而影响燃料棒表面的对流换热量;同时,燃料棒表面的对流换热量又影响冷却剂流场和温度场,二者存在耦合关系;因此在燃料棒表面的对流换热量的求解过程中,需对冷却剂流动换热特性分析模型和燃料棒热工特性分析模型之间进行流固耦合;
具体耦合过程为:通过燃料棒热工特性分析模型获得燃料棒温度分布的计算结果;根据冷却剂流场分布,通过上述Dittus-Bolter公式获得燃料棒表面的对流换热系数;将燃料棒表面的对流换热系数、冷却剂温度、燃料棒温度代入步骤2中的表面对流换热公式中,获得燃料棒表面的对流换热量;
通过上述步骤3建立了流固耦合换热模型;
步骤4:进行流固耦合求解
在计算开始时对步骤1中构建的反应堆堆芯的计算模型进行初始化,假定初始燃料棒温度分布、冷却剂流场和温度场;在计算过程中,燃料棒热工特性分析模型通过TN-1时刻内的堆芯冷却剂流场和温度场获得燃料棒表面的对流换热系数,并以此作为边界条件完成TN时刻的计算,得到TN时刻燃料棒温度分布的计算结果;冷却剂流动换热特性分析模型根据燃料棒温度分布的计算结果,获得燃料棒表面的对流换热量,并将其添加为步骤1中的能量方程的源项进行求解,得到TN时刻的冷却剂流场和温度场,自此完成TN-1时刻至TN时刻的求解;
步骤5:判断冷却剂流速、温度、燃料棒内部温度的残差是否满足收敛判据——即小于预设的残差值,若不满足则重复步骤4,直至收敛;最终获得在反应堆堆芯内燃料棒温度分布、冷却剂流场和温度场的计算结果,实现反应堆堆芯热工水力特性分析。
本发明具有以下优点和效果:
1.将堆芯内燃料棒和控制棒等复杂棒束结构简化,从而避免了逐个对燃料棒建立复杂几何模型;
2.克服了一维***分析程序无法获得堆芯燃料棒准确的温度分布、对功率分布不均匀瞬态计算过程响应模拟较差的问题;
3.对堆芯堆芯内复杂流动状态进行简化,通过冷却剂的流动阻力反映堆芯区域中燃料棒、控制棒等固体结构的影响;
4.通过燃料棒热工特性分析模型与冷却剂流动换热特性分析模型之间的数据交换,形成流固耦合换热模型,能够在分析冷却剂流动换热特性的同时获得堆芯内燃料棒温度分布;
5.计算资源消耗少、计算速度快且计算结果精度高。
本发明已经通过实践证明,该方法能够准确获得堆芯内燃料棒温度分布、冷却剂流场和温度场的计算结果,本发明中提出的基于计算流体力学的反应堆堆芯热工水力特性分析方法可应用于三维计算流体力学数值模拟中,在显著降低计算资源消耗的同时实现反应堆堆芯燃料棒温度分布的精确计算,具有建模简单、计算速度快的特点。
附图说明
图1为反应堆堆芯燃料棒及燃料组件示意图。
图2为燃料棒热工特性分析模型计算节点划分示意图。
图3为流固耦合换热模型框图。
图4为反应堆堆芯热工水力特性分析方法框图。
具体实施方式
为了对本发明的技术特征、目的和效果有更加清楚的理解,现结合附图详细说明本发明的具体实施方式。
本发明提出了一种基于计算流体力学的反应堆堆芯热工水力特性分析方法,框架如图4所示。该方法具体实施方式如下:
步骤1:建立冷却剂流动换热特性分析模型
(1)建立反应堆堆芯计算模型
通过几何建模软件建立反应堆堆芯的几何模型,在建模时对反应堆堆芯整体建立几何模型,而不需要对单个燃料组件或单个燃料棒进行精细建模。将反应堆堆芯的几何模型导入网格生成软件——如ICEM CFD中,对几何模型沿径向和轴向划分控制体,形成反应堆堆芯的网格模型。将反应堆堆芯的网格模型导入计算流体力学软件——如Fluent中,建立反应堆堆芯的计算模型。
(2)建立冷却剂流动换热特性分析模型
为了将如图1所示的反应堆堆芯内复杂棒束结构的影响简化为对冷却剂的流动阻力,提出如下假设:①堆芯的流动阻力用包含粘性损失、惯性损失的压降描述:
式中:是沿m坐标方向的动量方程中的压降,m和n是x,y,z坐标方向,Dmn和Cmn分别是粘性损失项和惯性损失项的系数矩阵,μ是冷却剂的动力粘性系数/kg·m-1·s-1,ρ是冷却剂的密度/kg·m-3,|v|是速度的模/m·s-1,vn是沿n坐标方向的速度/m·s-1;②计算域中燃料棒与定位格架固体部件体积不变,位置无移动;③孔隙率是各向同性的;基于上述假设,对流体控制方程中的动量方程和能量方程中添加源项:在动量方程中添加包含粘性损失、惯性损失的压降作为动量源项,在能量方程中添加燃料棒表面的对流换热量作为能量源项,对流换热量的具体值通过燃料棒热工特性分析模型计算;从而建立了冷却剂流动换热特性分析模型。
步骤2:建立燃料棒热工特性分析模型
燃料棒热工特性分析模型中考虑燃料棒表面对流换热和燃料棒内部导热,包括表面对流换热模型和内部导热模型,分别建立表面对流换热模型和内部导热模型,步骤如下:
(1)建立表面对流换热模型
提出如下假设:①反应堆燃料棒及控制棒结构均匀地分布在堆芯中;②单个控制体内燃料棒几何结构不发生剧烈变化;③单个控制体内所有燃料棒具有相同的热工状态;④单个控制体内所有燃料棒具有相同的几何特性且已知。
基于上述假设,反应堆堆芯内的换热面积与反应堆堆芯体积是线性相关的,建立反应堆堆芯内的换热面积与反应堆堆芯体积的线性关系式,即:
式中,AT,cell为控制体内的换热面积/m2,Vcell为控制体体积/m3,AT,core为反应堆堆芯内的换热面积/m2,Vcore为反应堆堆芯体积/m3。
由于反应堆堆芯内的换热面积和反应堆堆芯体积是已知的,将该关系式应用于步骤1中建立的反应堆堆芯的计算模型中,用控制体的体积表示控制体内的换热面积,即:
最终将控制体内的换热面积代入到表面对流换热公式中,计算得到燃料棒表面的对流换热量:
QT=hTAcell(Tr-Tf)
式中:QT为燃料棒表面的对流换热量/W,hT为燃料棒表面的对流换热系数/W·m-2·K-1,Acell为控制体内的换热面积/m2,Tf为流体主流温度/K,Tr为燃料棒表面温度/K。
通过上述步骤2的步骤(1)建立了表面对流换热模型。
(2)内部导热模型建立步骤
根据燃料棒的结构特点将燃料棒沿径向划分多个计算节点,其中外侧两个计算节点分别位于包壳内外两侧边界,其它计算节点沿燃料棒芯块中心向外侧依次布置,且有一个计算节点在芯块边界上,在包壳内外两侧各有一个计算节点;忽略燃料棒的轴向导热,沿燃料棒径向建立温度方程:
式中:ρi为计算节点i处材料密度/kg·m-3,CP,i为计算节点i处材料比热容/kJ·kg-1·K-1,Vi为计算节点i的等效控制体体积/m3,Ti为计算节点i处的温度/K,Qi-1,i为从计算节点i-1传导到计算节点i的热量/kW,Qi+1,i为计算节点i+1传导到计算节点i的热量/kW,Qi为计算节点i处单位体积释热率/kW·m-3。
燃料棒内部计算节点间传导的热量通过计算节点间温差得到:
Qi-1,i=Ki-1,i(Ti-1-Ti)
Qi+1,i=Ki+1,i(Ti+1-Ti)
式中:Ti-1为计算节点i-1处的温度/K,Ti+1为计算节点i+1处的温度/K,Ki-1,i为计算节点i-1和i之间的导热系数/kW·K-1,等同于Ki,i-1,Ki+1,i为计算节点i+1和i之间的导热系数/kW·K-1,等同于Ki,i+1,均为热导率和燃料棒几何参数的函数:
式中:ki是计算节点i处的热导率/kW·m-1·K-1,ki-1是计算节点i-1处的热导率/kW·m-1·K-1,ki+1是计算节点i+1处的热导率/kW·m-1·K-1,ri是计算节点i处的芯块直径/m,ri-1是计算节点i处的芯块直径/m,ΔXj为控制体轴向高度/m,为计算节点i的体积平均半径/m,为计算节点i-1的体积平均半径/m,为计算节点i+1的体积平均半径/m,其表达式为:
根据燃料棒的导热特点,建立如下边界条件:①燃料棒中心为对称边界,因此对于计算节点1的控制方程有K0,1=Ki-1,i=0;②在燃料芯块与包壳之间存在间隙导热,即在i=N-1处Qi,i-1=Hgap(Ti-Ti-1),其中Hgap为考虑换热面积的间隙导热率/kW·K-1,Qi,i-1为从计算节点i传导到计算节点i-1的热量/kW,Ti为计算节点i处的温度/K,Ti-1为计算节点i-1处的温度/K;③在燃料棒外表面,即在i=N处,
Ki+1,i=H,Ti+1=Tf,则有Qi+1,i=H(Tf-Ti),其中,Ki+1,i为计算节点i+1和i之间的导热系数/kW·K-1,H为考虑换热面积的燃料棒表面的对流换热系数/kW·K-1,Qi+1,i为计算节点i+1传导到计算节点i的热量/kW,Tf为流体主流温度/K,Ti+1为计算节点i+1处的温度/K,Ti为计算节点i处的温度/K。
将温度方程在各个计算节点上化为差分方程:
式中:Ti n为计算节点i处在第n个时间步下的温度/K,Ti n-1为计算节点i处在第n-1个时间步下的温度/K,Ti n+1为计算节点i处在第n+1个时间步下的温度/K。
将差分方程变形为:
根据边界条件可知,i=1时:
A1=0
对每个计算节点建立能量方程后联立求解,系数矩阵为:
该系数矩阵为三对角矩阵,可由高斯迭代法求解得到燃料棒内各计算节点温度。
步骤3:建立流固耦合换热模型
冷却剂与燃料棒之间的换热系数由Dittus-Bolter关系式计算:
Nu=0.023Re0.8Pr0.4
式中:Nu为努塞尔数,Re为雷诺数,Pr为普朗特数;Nu和Re的表达式为:
式中:k为冷却剂热导率/W·m-1·K-1,L为特征长度/m,ρ为冷却剂密度/kg·m-3,μ为冷却剂的动力粘性系数/kg·m-1·s-1,U为冷却剂流速/m·s-1。
由于k、ρ、μ、U均为冷却剂热物性和流动状态的参数,因此冷却剂流场和温度场决定了燃料棒表面的对流换热系数,进而影响燃料棒表面的对流换热量;同时,燃料棒表面的对流换热量又影响冷却剂流场和温度场,二者存在耦合关系;因此在燃料棒表面的对流换热量的求解过程中,需对冷却剂流动换热特性分析模型和燃料棒热工特性分析模型之间进行流固耦合。
具体耦合过程为:通过燃料棒热工特性分析模型获得燃料棒温度分布的计算结果;根据冷却剂流场分布,通过上述Dittus-Bolter公式获得燃料棒表面的对流换热系数;将燃料棒表面的对流换热系数、冷却剂温度、燃料棒温度代入步骤2中的表面对流换热公式中,获得燃料棒表面的对流换热量。
通过上述步骤3建立了流固耦合换热模型。
步骤4:进行流固耦合求解
流固耦合求解方法如图3所示。在计算开始时对步骤1中构建的反应堆堆芯的计算模型进行初始化,假定初始燃料棒温度分布、冷却剂流场和温度场;在计算过程中,燃料棒热工特性分析模型通过TK-1时刻内的堆芯冷却剂流场和温度场获得燃料棒表面的对流换热系数,并以此作为边界条件完成TK时刻的计算,得到TK时刻燃料棒温度分布的计算结果;冷却剂流动换热特性分析模型根据燃料棒温度分布的计算结果,获得燃料棒表面的对流换热量,并将其添加为步骤1中的能量方程的源项进行求解,得到TK时刻的冷却剂流场和温度场,自此完成TK-1时刻至TK时刻的求解。
步骤5:判断冷却剂流速、温度、燃料棒内部温度的残差是否满足收敛判据——即小于预设的残差值,若不满足则重复步骤4,直至收敛;最终获得在反应堆堆芯内燃料棒温度分布、冷却剂流场和温度场的计算结果,实现反应堆堆芯热工水力特性分析。
Claims (1)
1.一种基于计算流体力学的反应堆堆芯热工水力特性分析方法,其特征在于:具体包括如下步骤:
步骤1:建立冷却剂流动换热特性分析模型
(1)建立反应堆堆芯的计算模型
建立反应堆堆芯的几何模型,将反应堆堆芯的几何模型导入网格生成软件中划分控制体,形成反应堆堆芯的网格模型,将反应堆堆芯的网格模型导入计算流体力学软件中,形成反应堆堆芯的计算模型;
(2)建立冷却剂流动换热特性分析模型
为了将反应堆堆芯内复杂棒束结构的影响简化为对冷却剂的流动阻力,提出如下假设:①堆芯的流动阻力使用包含粘性损失、惯性损失的压降描述:
式中:▽pm是沿m坐标方向的动量方程中的压降,m和n是x,y,z坐标方向,Dmn和Cmn分别是粘性损失项和惯性损失项的系数矩阵,μ是冷却剂的动力粘性系数/kg·m-1·s-1,ρ是冷却剂的密度/kg·m-3,|v|是速度的模/m·s-1,vn是沿n坐标方向的速度/m·s-1;②计算域中燃料棒与定位格架固体部件体积不变,位置无移动;③孔隙率是各向同性的;基于上述假设,对流体控制方程中的动量方程和能量方程中添加源项:在动量方程中添加包含粘性损失、惯性损失的压降作为动量源项,在能量方程中添加燃料棒表面的对流换热量作为能量源项;从而建立了冷却剂流动换热特性分析模型;
步骤2:建立燃料棒热工特性分析模型
燃料棒热工特性分析模型中考虑燃料棒表面对流换热和燃料棒内部导热,包括表面对流换热模型和内部导热模型,分别建立表面对流换热模型和内部导热模型,步骤如下:
(1)建立表面对流换热模型
提出如下假设:①反应堆燃料棒及控制棒结构均匀地分布在堆芯中;②单个控制体内燃料棒几何结构不发生剧烈变化;③单个控制体内所有燃料棒具有相同的热工状态;④单个控制体内所有燃料棒具有相同的几何特性且已知;
基于上述假设,建立反应堆堆芯内的换热面积与反应堆堆芯体积的线性关系式:
式中:AT,cell为控制体内的换热面积/m2,Vcell为控制体体积/m3,AT,core为反应堆堆芯内的换热面积/m2,Vcore为反应堆堆芯体积/m3;
将该线性关系式应用于步骤1中建立的反应堆堆芯的计算模型中,用控制体的体积表示控制体内的换热面积;最终将控制体内的换热面积代入到表面对流换热公式中,计算得到燃料棒表面的对流换热量:
QT=hTAcell(Tr-Tf)
式中:QT为燃料棒表面的对流换热量/W,hT为燃料棒表面的对流换热系数/W·m-2·K-1,Acell为控制体内的换热面积/m2,Tf为流体主流温度/K,Tr为燃料棒表面温度/K;
通过上述步骤2的步骤(1)建立了表面对流换热模型;
(2)内部导热模型建立步骤
根据燃料棒的结构特点将燃料棒沿径向划分N个计算节点,其中外侧两个计算节点分别位于包壳内外两侧边界,其它计算节点沿燃料棒芯块中心向外侧依次布置,且有一个计算节点在芯块边界上,在包壳内外两侧各有一个计算节点;忽略燃料棒的轴向导热,沿燃料棒径向建立温度方程:
式中:ρi为计算节点i处材料密度/kg·m-3,CP,i为计算节点i处材料比热容/kJ·kg-1·K-1,Vi为计算节点i的等效控制体体积/m3,Ti为计算节点i处的温度/K,Qi-1,i为从计算节点i-1传导到计算节点i的热量/kW,Qi+1,i为计算节点i+1传导到计算节点i的热量/kW,Qi为计算节点i处单位体积释热率/kW·m-3;
燃料棒内部计算节点间传导的热量通过计算节点间温差得到:
Qi-1,i=Ki-1,i(Ti-1-Ti)
Qi+1,i=Ki+1,i(Ti+1-Ti)
式中:Ti-1为计算节点i-1处的温度/K,Ti+1为计算节点i+1处的温度/K,Ki-1,i为计算节点i-1和i之间的导热系数/kW·K-1,等同于Ki,i-1,Ki+1,i为计算节点i+1和i之间的导热系数/kW·K-1,等同于Ki,i+1,均为热导率和燃料棒几何参数的函数:
式中:ki是计算节点i处的热导率/kW·m-1·K-1,ki-1是计算节点i-1处的热导率/kW·m-1·K-1,ki+1是计算节点i+1处的热导率/kW·m-1·K-1,ri是计算节点i处的芯块直径/m,ri-1是计算节点i-1处的芯块直径/m,ΔXj为控制体轴向高度/m,为计算节点i的体积平均半径/m,为计算节点i-1的体积平均半径/m,为计算节点i+1的体积平均半径/m,其表达式为:
根据燃料棒的导热特点,建立如下边界条件:①燃料棒中心为对称边界,因此对于计算节点1的控制方程有K0,1=Ki-1,i=0;②在燃料芯块与包壳之间存在间隙导热,即在i=N-1处Qi,i-1=Hgap(Ti-Ti-1),其中Hgap为考虑换热面积的间隙导热率/kW·K-1,Qi,i-1为从计算节点i传导到计算节点i-1的热量/kW,Ti为计算节点i处的温度/K,Ti-1为计算节点i-1处的温度/K;③在燃料棒外表面,即在i=N处,Ki+1,i=H,Ti+1=Tf,则有Qi+1,i=H(Tf-Ti),其中,Ki+1,i为计算节点i+1和i之间的导热系数/kW·K-1,H为考虑换热面积的燃料棒表面的对流换热系数/kW·K-1,Qi+1,i为计算节点i+1传导到计算节点i的热量/kW,Tf为流体主流温度/K,Ti+1为计算节点i+1处的温度/K,Ti为计算节点i处的温度/K;
将温度方程在各个计算节点上化为差分方程:
式中:Ti n为计算节点i处在第n个时间步下的温度/K,Ti n-1为计算节点i处在第n-1个时间步下的温度/K,Ti n+1为计算节点i处在第n+1个时间步下的温度/K;
将差分方程变形为:
联立各个计算节点的差分方程形成方程组,其中方程组的系数矩阵为三对角矩阵,通过高斯迭代法求解得到燃料棒内各计算节点温度;
通过上述步骤2的步骤(2)建立了内部导热模型;
步骤3:建立流固耦合换热模型
冷却剂与燃料棒之间的换热系数由Dittus-Bolter关系式计算:
Nu=0.023Re0.8Pr0.4
式中:Nu为努塞尔数,Re为雷诺数,Pr为普朗特数;Nu和Re的表达式为:
式中:k为冷却剂热导率/W·m-1·K-1,L为特征长度/m,ρ为冷却剂密度/kg·m-3,μ为冷却剂的动力粘性系数/kg·m-1·s-1,U为冷却剂流速/m·s-1;
由于k、ρ、μ、U均为冷却剂热物性和流动状态的参数,因此冷却剂流场和温度场决定了燃料棒表面的对流换热系数,进而影响燃料棒表面的对流换热量;同时,燃料棒表面的对流换热量又影响冷却剂流场和温度场,二者存在耦合关系;因此在燃料棒表面的对流换热量的求解过程中,需对冷却剂流动换热特性分析模型和燃料棒热工特性分析模型之间进行流固耦合;
具体耦合过程为:通过燃料棒热工特性分析模型获得燃料棒温度分布的计算结果;根据冷却剂流场分布,通过上述Dittus-Bolter公式获得燃料棒表面的对流换热系数;将燃料棒表面的对流换热系数、冷却剂温度、燃料棒温度代入步骤2中的表面对流换热公式中,获得燃料棒表面的对流换热量;
通过上述步骤3建立了流固耦合换热模型;
步骤4:进行流固耦合求解
在计算开始时对步骤1中构建的反应堆堆芯的计算模型进行初始化,假定初始燃料棒温度分布、冷却剂流场和温度场;在计算过程中,燃料棒热工特性分析模型通过TN-1时刻内的堆芯冷却剂流场和温度场获得燃料棒表面的对流换热系数,并以此作为边界条件完成TN时刻的计算,得到TN时刻燃料棒温度分布的计算结果;冷却剂流动换热特性分析模型根据燃料棒温度分布的计算结果,获得燃料棒表面的对流换热量,并将其添加为步骤1中的能量方程的源项进行求解,得到TN时刻的冷却剂流场和温度场,自此完成TN-1时刻至TN时刻的求解;
步骤5:判断冷却剂流速、温度、燃料棒内部温度的残差是否满足收敛判据——即小于预设的残差值,若不满足则重复步骤4,直至收敛;最终获得在反应堆堆芯内燃料棒温度分布、冷却剂流场和温度场的计算结果,实现反应堆堆芯热工水力特性分析。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202110026314.9A CN112699620B (zh) | 2021-01-08 | 2021-01-08 | 基于计算流体力学的反应堆堆芯热工水力特性分析方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202110026314.9A CN112699620B (zh) | 2021-01-08 | 2021-01-08 | 基于计算流体力学的反应堆堆芯热工水力特性分析方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN112699620A CN112699620A (zh) | 2021-04-23 |
CN112699620B true CN112699620B (zh) | 2022-10-28 |
Family
ID=75513539
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202110026314.9A Active CN112699620B (zh) | 2021-01-08 | 2021-01-08 | 基于计算流体力学的反应堆堆芯热工水力特性分析方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN112699620B (zh) |
Families Citing this family (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN113793711B (zh) * | 2021-09-13 | 2022-12-02 | 西安交通大学 | 一种锂冷核反应堆与斯特林发电机耦合传热特性分析方法 |
CN114036871B (zh) * | 2021-11-26 | 2023-10-20 | 中国核动力研究设计院 | 基于瞬态分析的堆内棒状结构湍流激振分析方法及装置 |
CN114282460B (zh) * | 2021-12-31 | 2023-05-16 | 西安交通大学 | 一种铅铋快堆堆芯热工水力特性分析方法 |
CN114444413B (zh) * | 2022-01-21 | 2023-10-24 | 西安交通大学 | 一种板状燃料堆芯亚通道级三维热工水力分析方法 |
CN114757123B (zh) * | 2022-04-20 | 2023-06-20 | 西安交通大学 | 一种用于板形核燃料堆芯的跨维度流固耦合分析方法 |
CN115186603B (zh) * | 2022-06-22 | 2023-04-18 | 西安交通大学 | 一种计算燃料棒与夹持机构相互作用的方法 |
CN115577583B (zh) * | 2022-09-06 | 2023-07-04 | 西安交通大学 | 一种铅基反应堆绕丝定位燃料棒流致振动的分析方法 |
Citations (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JP2009156745A (ja) * | 2007-12-27 | 2009-07-16 | Hitachi-Ge Nuclear Energy Ltd | 炉心性能計算装置 |
WO2015017599A1 (en) * | 2013-08-01 | 2015-02-05 | Newcomb Eric William | Power demand management using thermal hydraulic generator |
CN108875213A (zh) * | 2018-06-19 | 2018-11-23 | 哈尔滨工程大学 | 一种反应堆堆芯热工水力多尺度分析方法 |
CN109299494A (zh) * | 2018-07-24 | 2019-02-01 | 哈尔滨工程大学 | 一种反应堆堆芯热工水力多尺度耦合计算的数据重构方法 |
CN109903870A (zh) * | 2019-03-15 | 2019-06-18 | 西安交通大学 | 一种核动力***跨维度耦合模拟方法 |
CN110598304A (zh) * | 2019-09-06 | 2019-12-20 | 西安交通大学 | 一种空间核电推进***球床反应堆物理热工耦合分析方法 |
CN111651944A (zh) * | 2020-06-05 | 2020-09-11 | 中国原子能科学研究院 | 核反应堆计算*** |
Family Cites Families (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110543704B (zh) * | 2019-08-19 | 2021-05-28 | 西安交通大学 | 一种在反应堆堆芯流场计算中考虑局部结构影响的修正方法 |
-
2021
- 2021-01-08 CN CN202110026314.9A patent/CN112699620B/zh active Active
Patent Citations (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JP2009156745A (ja) * | 2007-12-27 | 2009-07-16 | Hitachi-Ge Nuclear Energy Ltd | 炉心性能計算装置 |
WO2015017599A1 (en) * | 2013-08-01 | 2015-02-05 | Newcomb Eric William | Power demand management using thermal hydraulic generator |
CN108875213A (zh) * | 2018-06-19 | 2018-11-23 | 哈尔滨工程大学 | 一种反应堆堆芯热工水力多尺度分析方法 |
CN109299494A (zh) * | 2018-07-24 | 2019-02-01 | 哈尔滨工程大学 | 一种反应堆堆芯热工水力多尺度耦合计算的数据重构方法 |
CN109903870A (zh) * | 2019-03-15 | 2019-06-18 | 西安交通大学 | 一种核动力***跨维度耦合模拟方法 |
CN110598304A (zh) * | 2019-09-06 | 2019-12-20 | 西安交通大学 | 一种空间核电推进***球床反应堆物理热工耦合分析方法 |
CN111651944A (zh) * | 2020-06-05 | 2020-09-11 | 中国原子能科学研究院 | 核反应堆计算*** |
Non-Patent Citations (5)
Title |
---|
Thermal hydraulic simulation of reactor of HTR-PM based on thermal-fluid network and SIMPLE algorithm;Yangping Zhou,et al.;《Progress in Nuclear Energy》;20130131;第83-93页 * |
Thermo-hydraulic Analysis of Compact Heat Exchanger for a Simple Recuperated sCO2 Brayton Cycle;Pandey,V.,et al.;《renewable and sustainable energy reviews》;20200912;第1-43页 * |
一体化小型压水堆一回路热工水力特性数值模拟研究;陈焕栋;《中国论文全文数据库 工程科技II辑》;20190515;C034-58 * |
某型空间堆堆芯热工水力特性数值分析;任董国 等;《核科学与工程》;20201031;第742-750页 * |
空间核反应堆电源热工水力特性研究综述;代智文 等;《原子能科学技术》;20190731;第1296-1309页 * |
Also Published As
Publication number | Publication date |
---|---|
CN112699620A (zh) | 2021-04-23 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN112699620B (zh) | 基于计算流体力学的反应堆堆芯热工水力特性分析方法 | |
CN110598324B (zh) | 一种核反应堆弥散型板型燃料元件堆芯流固耦合计算方法 | |
CN113094947B (zh) | 一种核反应堆堆芯核热耦合分析方法 | |
CN108170924B (zh) | 一种用于核电厂蒸汽发生器传热管堵流工况模型建立方法 | |
CN109903870B (zh) | 一种核动力***跨维度耦合模拟方法 | |
CN114282460B (zh) | 一种铅铋快堆堆芯热工水力特性分析方法 | |
CN111291494B (zh) | 用于核反应堆triso燃料颗粒的多尺度多物理场耦合模拟方法 | |
CN114444413B (zh) | 一种板状燃料堆芯亚通道级三维热工水力分析方法 | |
CN115270660B (zh) | 空间热离子反应堆瞬态行为多尺度多物理场耦合分析方法 | |
JP6404904B2 (ja) | 原子炉容器内部での流体の流れのシミュレーション方法、及び、原子炉炉心の部品の機械的変形の計算方法、及び、関連するコンピュータープログラム製品 | |
CN111274748A (zh) | 池式钠冷快堆非能动余热排出***跨维度耦合计算方法 | |
CN104331539A (zh) | 一种核电站管道热分层效应疲劳评价方法及*** | |
Kim et al. | Development of thermal-hydraulic analysis methodology for multiple modules of water-cooled breeder blanket in fusion DEMO reactor | |
CN116070424A (zh) | 一种铅铋快堆铅池热分层降阶分析方法 | |
Gu et al. | Verification of a HC-PK-CFD coupled program based a benchmark on beam trip transients for XADS reactor | |
Hernandez et al. | Development of a CFD-based model to simulate loss of flow transients in a small lead-cooled reactor | |
Zadeh et al. | Effect of 3-D moderator flow configurations on the reactivity of CANDU nuclear reactors | |
CN114757123A (zh) | 一种用于板形核燃料堆芯的跨维度流固耦合分析方法 | |
CN105247622B (zh) | 核反应堆中模拟流体流动及计算燃料组件机械形变的方法 | |
Wardhani et al. | Analysis of coolant flow distribution to the reactor core of modified TRIGA Bandung with plate-type fuel | |
Feng et al. | Development of a multi-physics coupling system based on ICoCo interface and its validation on NEA-OECD core transient benchmark | |
Zhu et al. | A CFD Method for Thermal-Hydraulic Calculation of Pebble Bed HTGRs | |
CN115221813B (zh) | 一种核-热-力跨维度耦合方法 | |
CN116504431A (zh) | 一种钠冷快堆堆芯核热耦合方法 | |
CN116484764A (zh) | 一种钠冷快堆多维度耦合计算方法 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |