本发明涉及核研究，尤其涉及一种稳定高效的反应堆物理热工耦合模拟的全局耦合迭代方法。
背景技术：
1、目前，针对反应堆物理热工耦合模拟的需求，市面上可选取的耦合求解方法主要有picard迭代方法和jacobian-free matrix krylov(jfnk)方法。
2、picard迭代方法基于解耦的理念，将多物理场耦合系统按照其物理意义分解为若干子物理场。在此方法中，各子物理场通过相互独立的求解器交替求解，并通过求解器间的参数传递来实现迭代，直至所有物理场在一个时间步内全部收敛。因此，picard迭代方法在不同子物理场之间添加了一个外循环和收敛判据。
3、jfnk方法是一种全局耦合迭代方法，是newton方法、krylov子空间方法以及jacobian-free方法的有机结合。该方法首先将耦合系统整合为一个非线性方程组，并采用newton方法进行局部线性化，然后利用krylov子空间方法求解线性化方程组。同时，结合jacobian-free方法，避免了jacobian矩阵的显式构建和存储。因此，jfnk方法同步联立求解所有物理量，在稳定性和收敛特性上更有优势。
4、picard迭代方法面临稳定性差、迭代收敛速度缓慢，甚至迭代发散的缺陷。特别是对于耦合系统的耦合非常强烈或者部分物理量对其他物理量的变化十分敏感的情况，picard迭代可能会出现发散现象。此外，该方法仅具备一阶的收敛速度，且“过解”问题会使得计算效率低。最后，picard方法采用解耦的思想，这样处理弱化了各物理场间的耦合作用，仅能用于少数几个物理场的耦合系统。
5、jfnk方法在理论分析和实际应用中都极具挑战：①残差函数的构建具有一定的挑战性，特别是对于高保真物理-热工耦合模拟，其残差函数的构建甚至可能无法实现；②不同物理场的物理量在数量级上存在极大地差异；③jfnk方法的收敛性对newton步初值的选取敏感，若初值选取不合理可能导致求解发散；④jfnk方法中差分步长的选取往往依赖于经验公式，这对计算精度和效率影响较大。
6、picard迭代方法的不稳定性是由多普勒反馈影响导致的，这使得该方法不能保证无条件收敛。其次，在初始迭代步中，对单个物理场进行过多迭代求解收敛，即对一个不精确问题过度精确求解，从而引起“过解”问题。
7、jfnk方法需要在全局耦合框架下以残差函数的形式描述各物理场，而成熟的单物理场程序并未直接提供残差函数，因此，需要根据各子物理场的离散关系式重新构建。其次，jacobian矩阵的行表示不同的物理方程，列表示不同的物理量，这导致jacobian矩阵的元素在数值上存在显著差异。随后，由于jfnk方法具有局部二阶收敛特性，这使得newton步初值的选取需要在其局部收敛区间之内。最后，jfnk方法采用了jacobian-free技术，使用差分步长近似jacobian矩阵和基向量的乘积，这使得差分步长的选取会影响计算精度和收敛性能。此外，jfnk方法还要求所有物理场能且只能采用同一个差分步长。
8、对于picard迭代方法，虽然采用松弛因子可以提高算法的稳定性，但松弛因子的选取具有强烈的经验性质，且伴随耦合系统复杂性的增加，确定适当的松弛因子更加困难。其次，通过非线性交替技术(anm)和局部收敛技术(pct)等可有效缓解“过解”问题，然而这些方法过于死板，不一定适用于各种多物理场问题。
9、对于jfnk方法，非线性预处理基于不动点迭代可利用现有成熟的单物理场程序构建等价的残差函数，但每次krylov子空间迭代都需要调用各物理场程序，这大大降低了计算效率。其次，尽管使用scaling技术可以降低jacobian矩阵条件数，但scaling因子的选取通常依赖于经验和物理背景，且条件数的改善效果有限。
10、总体来说，尽管市面上常见的picard迭代方法和jfnk方法都经过了一定程度的探索和优化，但它们仍存在一定的缺点和局限性。
技术实现思路
1、针对反应堆物理-热工耦合模拟，本发明提供了一种反应堆物理热工耦合模拟的全局耦合迭代方法，即mfnk方法。该方法将耦合系统整合为一个非线性残差方程，并同步联立求解所有耦合物理量，从而确保了较高的稳定性和收敛性能。
2、本发明提供一种反应堆物理热工耦合模拟的全局耦合迭代方法，其包括如下步骤：
3、步骤1：根据用户提供的输入卡，读取几何和材料信息，并使用构造实体几何方法构建中子物理场模型和热工水力模型。其中，中子物理场模型为精细化中子输运模型，热工水力模型为单通道模型和圆柱导热模型；
4、步骤2：中子物理模型和热工水力模型采用了不同的网格划分策略，建立空间网格映射关系，实现不同物理场间耦合参数的传递；
5、步骤3：利用mfnk方法计算反应堆物理-热工耦合；
6、步骤4：通过均匀化处理，建立完全等价的低阶ml-gcmfd模型、单通道模型和圆柱导热模型，从而将高阶非线性残差系统转换为等价的低阶非线性残差系统；
7、步骤5：在求解低阶非线性残差方程前，使用eisenstat-walker策略计算强制因子ηk，并使用混合差分步长公式计算ε；
8、步骤6：采用非精确newton方法求解低阶非线性残差方程，从而得到newton局部线性化方程；
9、步骤7：使用krylov子空间方法求解局部线性化方程，采用mgs-gmres方法执行内迭代求解，直至满足非精确newton方程迭代求解的收敛条件，获得功率密度的变化量
10、步骤8：更新功率密度随后只顺序求解一次各物理场，获得新的功率密度
11、步骤9：使用功率密度计算均方根rms_pd，并与设定的耦合收敛条件进行对比；
12、若满足收敛条件，则耦合系统求解完成；否则，返回步骤4开启下一个newton步的计算；
13、步骤10：基于收敛的功率密度，重新执行一次热工计算，获得收敛的热工参数信息。
14、在具体实施方式中，步骤3：利用mfnk方法计算反应堆物理-热工耦合；
15、具体地，中子物理场和热工水力场的计算过程表示为：
16、
17、式中，为功率密度，可由中子通量计算得到，上标k为newton迭代步；fth和fn_transport为概念性的求解器，包含实际求解器和对计算结果进行处理的所有附加操作；
18、当式(1-1)和(1-2)收敛时，联立简化为一个整体系统：
19、
20、进而，构建的物理-热工耦合系统的高阶非线性残差方程：
21、
22、其中，mfnk方法选取功率密度为耦合物理量，迭代初始值是通过一次picard迭代计算得到；
23、在具体实施方式中，步骤4：采用ml-gcmfd耦合加速方法，在等价的低阶系统中捕捉热工反馈的影响；
24、具体地，通过均匀化处理，建立完全等价的低阶ml-gcmfd模型、单通道模型和圆柱导热模型，并将与中子输运完全等价的ml-gcmfd计算作为中子物理场的求解器；由此建立等价的低阶非线性残差方程，形式如下：
25、
26、低阶非线性残差方程采用了非线性预处理技术，对于线性系统其相应的非线性系统为：
27、
28、将系数矩阵a分解为a＝m-n，进而线性系统可转变为：
29、
30、式中，代表一次picard计算，可表示为：
31、
32、随后，使用m-1作为左非线性预处理子，对式(5)进行同解变换，可得：
33、
34、结合式(7)和(8)，并用代表使用非线性预处理后新的残差方程：
35、
36、在具体实施方式中，步骤5：在求解低阶非线性残差方程前，完成下述参数的准备工作：使用eisenstat-walker策略计算强制因子ηk和使用混合差分步长计算公式获取ε；
37、具体地，所述混合差分步长计算公式如下：
38、
39、其中，为基向量，εmach为计算机的计算精度，sk为第k次newton步的增量，η为“黑箱”计算时的相对偏差；
40、在具体实施方式中，步骤6中，对低阶非线性残差方程进行泰勒展开并忽略高阶项，如下：
41、
42、进一步简化处理，获得newton局部线性化方程，如下所示：
43、
44、在具体实施方式中，步骤7：采用mgs-gmres方法执行内迭代求解，直至满足非精确newton方程迭代求解的收敛条件：
45、
46、计算得到功率密度的变化量此过程中，一方面，jacobian矩阵的构建和存储困难，另一方面，mgs-gmres求解时仅需计算jacobian矩阵和基向量的乘积。为此，mfnk方法中使用了jacobian-free方法，即通过有限差分来近似矩阵向量积，如下所示：
47、
48、其中，为krylov基向量，ε为混合差分步长；
49、在具体实施方式中，步骤8中：基于gmres计算得到的功率密度的变化量通过来更新功率密度，随后将更新后的功率密度作为输入，随后只顺序求解一次各物理场，获得新的功率密度
50、在具体实施方式中，步骤9中，耦合系统设定的收敛标准的计算公式如下：
51、
52、其中，n为残差方程的维度，atol为用户设定的收敛标准。
53、本发明由此提供基于所述方法的反应堆物理热工耦合模拟的全局耦合迭代系统，其包括相应于步骤1至步骤10的模块。
54、本发明也提供一种电子设备，其包括：至少一个处理器；以及，与所述至少一个处理器通信连接的存储器；
55、其中，所述存储器存储有可被所述至少一个处理器执行的指令，所述指令被设置为用于执行上述的方法或执行所述的系统。
56、本发明还提供一种计算机可读存储介质，其存储有计算机可执行指令，所述计算机可执行指令用于执行所述的方法流程或执行所述的反应堆物理热工耦合模拟的全局耦合迭代系统。
57、本发明的mfnk方法通过以下四个策略来消除jfnk方法的缺陷：
58、①不同于jfnk方法采用中子通量和温度作为耦合物理量，mfnk方法选择功率密度为中间求解变量，这样处理不仅降低了非线性方程的计算规模，还由于功率密度作为单一耦合物理量，使得jacobian矩阵的条件数较小，收敛更快，并且无需依赖scaling技术；
59、②采用非线性预处理技术，实施简单，能够实现“黑箱”耦合，即直接利用已有成熟的单物理场求解器来构建残差方程。同时，结合ml-gcmfd耦合加速方法，通过完全等价的低阶系统来捕捉热工反馈的影响，从而巧妙地避免了krylov子空间迭代时频繁调用中子输运计算带来的计算负担，显著提高了计算效率。
60、③当使用ml-gcmfd耦合加速方法时，需对每个newton步需执行一次非线性交替技术的picard计算，这样处理显著降低了mfnk方法对newton步初值的敏感性，增强了稳定性；
61、④mfnk方法采用了一种混合差分步长计算公式，结合了经验差分步长和最佳差分步长各自的优势，进一步提高了计算效率和收敛速度。
62、本发明可以达到的效果在于：采用mfnk方法能够实现对全堆芯高保真物理-热工耦合的稳定、快速、精确和高分辨率模拟。进一步地，本发明还得到以下更深层次的效果：基于本发明提出的mfnk方法，可以进一步扩展应用到反应堆多物理场耦合模拟领域，例如增加燃料性能的耦合等。
1.一种反应堆物理热工耦合模拟的全局耦合迭代方法，其特征在于，包括如下步骤：
2.如权利要求1所述的方法，其特征在于，步骤3：利用mfnk方法计算反应堆物理-热工耦合；
3.如权利要求1所述的方法，其特征在于，步骤4：采用ml-gcmfd耦合加速方法，在等价的低阶系统中捕捉热工反馈的影响；
4.如权利要求1所述的方法，其特征在于，步骤5：在求解低阶非线性残差方程前，完成下述参数的准备工作：使用eisenstat-walker策略计算强制因子ηk和使用混合差分步长计算公式获取ε；
5.如权利要求1所述的方法，其特征在于，步骤6中，对低阶非线性残差方程进行泰勒展开并忽略高阶项，如下：
6.如权利要求1所述的方法，其特征在于，步骤7：采用mgs-gmres方法执行内迭代求解，直至满足非精确newton方程迭代求解的收敛条件：
7.如权利要求1所述的方法，其特征在于，步骤8中：基于gmres计算得到的功率密度的变化量通过来更新功率密度，随后将更新后的功率密度作为输入，随后只顺序求解一次各物理场，获得新的功率密度
8.如权利要求1所述的方法，其特征在于，步骤9中，耦合系统设定的收敛标准的计算公式如下：
9.基于如权利要求1至8任一项所述方法的反应堆物理热工耦合模拟的全局耦合迭代系统，其包括相应于步骤1至步骤10的模块。
10.一种电子设备，其特征在于，包括：至少一个处理器；以及，与所述至少一个处理器通信连接的存储器；
