Solution methods, systems, storage media and applications of core subchannel coupling model

CN122242175BActive Publication Date: 2026-09-01HUAZHONG UNIV OF SCI & TECH
View PDF 0 Cites 0 Cited by

Patent Information

Application Number
CN202610695087.1
Authority / Receiving Office
CN · China
Patent Type
Patents(China)
Current Assignee / Owner
Filing Date
2026-05-20
Publication Date
2026-09-01
Estimated Expiration
2046-05-20

AI Technical Summary

Technical Problem

[0010]针对现有技术的以上缺陷或改进需求,本发明提供了一种堆芯子通道耦合模型求解方法、系统、存储介质及应用,其目的在于提高堆芯子通道耦合模型全隐式求解的收敛稳定性和计算效率,由此解决现有技术中收敛稳定性不足、固定点预迭代次数适应性有限、牛顿类方法对初值较为敏感以及雅可比矩阵显式构造代价较高的问题

Benefits of technology

1. 本发明融合了固定点迭代与无雅可比牛顿-克雷洛夫方法的各自优势,提出了一种用于堆芯子通道耦合模型的全隐式全局统一求解方法及系统,充分利用固定点迭代早期收敛的快速性和无雅可比牛顿-克雷洛夫方法的强稳定性,从而确保堆芯子通道耦合模型的高效性和强鲁棒性。

✦ Generated by Eureka AI based on patent content.

Smart Images

  • Figure CN122242175B_ABST
    Figure CN122242175B_ABST
Patent Text Reader

Abstract

This invention relates to the field of advanced nuclear reactor thermal-hydraulic analysis and multiphysics coupling calculation technology, specifically to a solution method, system, storage medium, and application for core subchannel coupling models. The method establishes a fully implicit discrete core subchannel coupling model and corresponding global nonlinear residual equations based on a drift flow model. It proposes a globally unified solution method combining fixed-point iteration with a Jacobi-Newton-Krylov solution. During the solution process, the global solution vector is first updated using fixed-point iteration to obtain the predicted solution, and a switching criterion is constructed based on the iteration increments of adjacent fixed points. When the switching condition is met, the fixed-point iteration operator is introduced as a physical preprocessing step into the Jacobi-Newton-Krylov solution framework. The solution is then completed by combining the Jacobi-vector product approximation constructed using finite-difference perturbation, the Krylov subspace iteration method, and a line search strategy, yielding the core thermal-hydraulic parameter distribution and effectively improving the stability and efficiency of the fully implicit solution of the core subchannel coupling model.
Need to check novelty before this filing date? Find Prior Art

Description

Technical Field

[0001] This invention belongs to the field of advanced nuclear reactor thermal-hydraulic analysis and multiphysics coupling calculation, and relates to the solution method, system, storage medium and application of core subchannel coupling model. More specifically, it relates to a fully implicit global unified solution method and system for core subchannel coupling model. Background Technology

[0002] Thermal-hydraulic analysis of nuclear reactor cores is an important part of nuclear power plant safety design, operation analysis, and accident assessment. To obtain the distribution patterns of thermal-hydraulic parameters such as flow rate, pressure, temperature, and void fraction within the core, subchannel analysis methods have become one of the commonly used methods in core thermal-hydraulic analysis due to their ability to balance computational accuracy and efficiency.

[0003] In the analysis of subchannels within the reactor core, it is typically necessary to perform coupled calculations of the mass, momentum, and energy transfer processes between the subchannels based on the two-phase flow governing equations and corresponding constitutive relations. As the complexity of the analysis objects and computational conditions continues to increase, the variable coupling relationships in the reactor core's thermal-hydraulic processes become more complex, and the nonlinear characteristics become more pronounced, thus placing higher demands on the stability, convergence, and computational efficiency of numerical solution methods.

[0004] Existing subchannel solutions often employ explicit or semi-implicit methods. While these methods are relatively mature, they often require increasing the number of computational steps to ensure stability when dealing with highly coupled core thermal-hydraulic problems, resulting in high overall computational costs. In contrast, fully implicit unified solution methods can discretize and solve multiple governing equations and related constitutive relations as a whole, thus showing great promise for application in core subchannel analysis.

[0005] However, fully implicit methods require the global coupling solution of the mass equation, momentum equation, energy equation, and related constitutive relations, resulting in a large and highly nonlinear system of algebraic equations, which significantly increases the difficulty of the solution. How to efficiently and stably solve this type of globally coupled nonlinear equation system has become a key problem in fully implicit subchannel analysis.

[0006] For the aforementioned nonlinear equation system, the fixed-point iterative method is often used as a fundamental solution approach due to its simplicity and ease of construction. This method typically achieves iterative solutions by updating the coupled equation system step by step. However, when the coupling is strong, it often suffers from slow convergence speed and sensitivity to parameter settings. To improve convergence performance, it is usually necessary to set relaxation factors or adopt empirical adjustment strategies, which to some extent increases the parameter dependence during the method's use.

[0007] Newton-like methods are another commonly used approach for solving nonlinear equation systems, exhibiting fast local convergence and theoretically high solution efficiency. However, for large-scale nonlinear equation systems formed by core subchannel models, traditional Newton methods typically require explicit construction, storage, and updating of the Jacobian matrix, leading to significant implementation complexity and storage overhead. To address this, the Jacobian-less Newton-Krylov method approximates the Jacobian-vector product using finite difference methods, avoiding explicit construction of the Jacobian matrix and offering advantages in solving large-scale nonlinear problems.

[0008] Nevertheless, the Jacobi-less Newton-Krylov method still belongs to the Newton class of solution methods, which places high demands on initial values ​​and preprocessing. When the initial value differs significantly from the objective solution or the preprocessing is inadequate, convergence and computational efficiency are greatly affected. To improve this, existing techniques often employ a joint solution strategy of performing several fixed-point pre-iterations before switching to the Jacobi-less Newton-Krylov method, utilizing the fixed-point iterations to obtain better predicted initial values ​​or preprocessing information.

[0009] However, in existing methods, the number of pre-iterations used for predicting the initial solution or as a preprocessing operator in fixed-point iteration is usually a pre-set fixed value, lacking the ability to adaptively adjust for specific computational processes. When the number of pre-iterations is set inappropriately, the fixed-point iteration stage may fail to play its full role, or unnecessary iterative calculations may occur, thus affecting the overall solution efficiency and convergence performance. Therefore, how to establish a reasonable switching mechanism between the fixed-point iteration stage and the non-Jacobi Newton-Krylov solution stage to balance solution stability and computational efficiency remains a pressing technical problem to be solved in this field. Summary of the Invention

[0010] To address the aforementioned deficiencies or improvement needs of existing technologies, this invention provides a solution method, system, storage medium, and application for a core subchannel coupling model. The aim is to improve the convergence stability and computational efficiency of the fully implicit solution of the core subchannel coupling model, thereby solving the problems of insufficient convergence stability, limited adaptability of fixed-point pre-iteration times, sensitivity of Newton-like methods to initial values, and high cost of explicit construction of the Jacobian matrix in existing technologies.

[0011] To achieve the above objectives, according to one aspect of the present invention, a method for solving a core subchannel coupling model is provided, comprising the following steps: S1: Divide the reactor core into sub-channel control volume grids, construct a fully implicit discrete core sub-channel coupling model based on the control equations and constitutive equations of the drift flow model; construct a global solution vector according to the coupling relationship between the variables in the core sub-channel coupling model, and establish the corresponding global nonlinear residual equation; S2: Construct a fixed-point iterative solution process, define a fixed-point iterative update of the core sub-channel coupling model as a fixed-point iterative operator, and use the fixed-point iterative operator to iteratively update the global solution vector to obtain the predicted solution; in the fixed-point iterative solution process, monitor the increment of two adjacent fixed-point iterations, and construct a switching criterion for the fixed-point iterative stage; S3: Determine whether the fixed-point iteration result meets the preset convergence condition; if it does, obtain the global solution of the core sub-channel coupling model; otherwise, perform a switching criterion judgment: If the switching criterion does not meet the preset switching conditions, then continue to execute step S2; If the switching criterion satisfies the preset switching condition, the fixed-point iterative solution is terminated, and the current fixed-point iterative solution is used as the initial prediction value for the Jacobi-less Newton-Krylov solution process. S4: Construct a Jacobi-free Newton-Krylov solution process, introduce the fixed-point iterative operator as a physical preprocessing into the Jacobi-free Newton-Krylov solution framework, and construct the Jacobi-vector product approximation using the finite difference perturbation method. S5: In the process of solving the problem without Jacobi Newton-Krylov, the global nonlinear residual equation is linearized to obtain the corresponding local linear equation. The Krylov subspace iteration method is used, combined with the physical preprocessing and the Jacobi-vector product approximation, to solve the local linear equation and obtain the increment of the global solution vector. S6: Update the global solution vector according to the increment of the global solution vector, and determine whether the updated global solution vector satisfies the global convergence condition; if it satisfies the condition, obtain the global solution of the core subchannel coupling model and output the corresponding core thermal-hydraulic parameter distribution; if it does not satisfy the condition, return to step S5.

[0012] Further, in step S1, the governing equations and constitutive equations of the drift flow model include the mixed energy equation, liquid phase energy equation, mixed continuity equation, axial momentum equation, and transverse momentum equation. These equations are discretized on the sub-channel control volume grid, and each term in the discrete equations is implicitly represented by the unknown quantity to be solved in the current iteration step to establish the fully implicit discrete equation set of the core sub-channel coupling model.

[0013] Further, in step S1, based on the coupling relationship between the variables in the core subchannel coupling model, the mixing enthalpy, axial mass flow rate, lateral mass flow rate, axial pressure difference and lateral pressure corresponding to each subchannel control body are assembled into a global solution vector.

[0014] Further, in step S1, the global solution vector is subjected to dimensionless scaling. The scaled dimensionless global solution vector... satisfy: , in, The original global solution vector, matrix It is a diagonal scaling matrix, whose diagonal elements are derived from the components corresponding to the original global solution vector. The norm is determined.

[0015] Further, in step S2, the fixed-point iteration operator takes the current global solution vector as input, and scans from the core inlet to the core outlet in a progressive manner along the axial direction. In each axial layer, it updates the mixing enthalpy, axial mass flow rate, lateral mass flow rate, axial pressure difference and lateral pressure corresponding to each sub-channel control body in sequence, and outputs the updated global solution vector.

[0016] Furthermore, in step S2, the switching criterion is based on the norm ratio of the iteration increments of two adjacent fixed points. Perform construction and judgment: , in, This represents the current iteration step number of the fixed-point iteration operator. , , The first , , The dimensionless global solution vector output by the next iteration. Represent the L2 norm of a vector; set the switching threshold. and the maximum number of fixed-point iterations ,in ;when and When, the fixed-point iteration operator is used to continue iterative calculation; when or When necessary, switch to the Jacobi-less Newton-Krylov method for subsequent global solution.

[0017] Furthermore, in step S6, let the first... The dimensionless global solution vector of the second-most Jacobian Newton-Krylov iteration is The corresponding increment is A line search strategy is used to determine the step size factor. The dimensionless global solution vector is then updated according to the following formula: , in, That is, the first The dimensionless global solution vector of the second-order Jacobi Newton-Krylov iteration.

[0018] According to another aspect of the present invention, a system for solving a core subchannel coupling model is provided, comprising a processor and a memory, wherein the memory stores a computer program, and the processor executes the computer program to implement the core subchannel coupling model solving method as described in any of the preceding claims.

[0019] According to another aspect of the invention, a computer-readable storage medium is provided having a computer program stored thereon, which, when executed by a processor, implements the core subchannel coupling model solution method as described in any of the preceding claims.

[0020] According to another aspect of the present invention, an application of the core subchannel coupling model solution method as described in any of the preceding claims is provided for solving the core subchannel coupling model under single-phase flow conditions and vapor-liquid two-phase flow conditions in a reactor core.

[0021] In summary, the technical solutions conceived in this invention, compared with the prior art, can achieve the following beneficial effects: 1. This invention integrates the advantages of fixed-point iteration and the Newton-Krylov method without Jacobi, and proposes a fully implicit global unified solution method and system for the core subchannel coupling model. It makes full use of the fast early convergence of fixed-point iteration and the strong stability of the Newton-Krylov method without Jacobi, thereby ensuring the efficiency and strong robustness of the core subchannel coupling model.

[0022] 2. This invention introduces an adaptive switching mechanism based on the relative change of adjacent iteration increments in the fixed-point iteration stage. This mechanism can determine the switching timing according to the convergence characteristics of the fixed-point iteration process, making the connection between the fixed-point iteration and the Jacobi-less Newton-Krylov solution process more reasonable. It avoids the additional computational overhead or decrease in solution effect caused by improper setting of the pre-iteration number, thereby improving the convergence stability of the fully implicit global solution and thus improving the algorithm's adaptability to different computational conditions and overall solution efficiency.

[0023] 3. This invention introduces a fixed-point iterative operator into the Jacobi-free Newton-Krylov solution framework as a physical preprocessor, and uses a finite difference perturbation method to construct the Jacobi-vector product approximation. This eliminates the need for explicit derivation, construction, and storage of the Jacobi matrix, thereby reducing the difficulty of program implementation and memory overhead. This is beneficial for the fully implicit global solution of large-scale heap core sub-channel coupling models. Attached Figure Description

[0024] Figure 1 This is a flowchart of the fully implicit global unified solution method for the core subchannel coupling model according to a preferred embodiment of the present invention.

[0025] Figure 2This is a schematic diagram of the component division in Embodiment 1 of the present invention.

[0026] Figure 3 This is a comparison of the residual convergence curves of the method described in Embodiment 1 of the present invention and the traditional fixed-point iteration method.

[0027] Figure 4 This is a comparison of the number of residual function calls and computation time between the present invention and the traditional fixed-point iteration method in Embodiment 1 of the present invention. Detailed Implementation

[0028] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative and not intended to limit the invention. Furthermore, the technical features involved in the various embodiments of this invention described below can be combined with each other as long as they do not conflict with each other.

[0029] Please see Figure 1 The core subchannel coupling model solution method provided in this embodiment is a fully implicit global unified solution method for core subchannel coupling models. This method first performs subchannel control volume mesh generation on the reactor core to establish the corresponding core subchannel coupling model. Then, based on a solution framework combining fixed-point iteration and Newton-Krylov without Jacobi, the core subchannel coupling model is globally solved, specifically including: The reactor core is meshed to determine the computational domain; a core sub-channel coupling model is established within the computational domain, and the core sub-channel coupling model is solved based on the fully implicit global unified solution method disclosed in this invention; more specifically, this invention includes: S1: Divide the reactor core into sub-channel control volume grids, construct a fully implicit discrete core sub-channel coupling model based on the control equations and constitutive equations of the drift flow model; construct a global solution vector according to the coupling relationship between the variables in the core sub-channel coupling model, and establish the corresponding global nonlinear residual equation; S2: Construct a fixed-point iterative solution process, define one fixed-point iterative update of the core sub-channel coupling model as a fixed-point iterative operator, and use the fixed-point iterative operator to iteratively update the global solution vector to obtain the predicted solution; in the fixed-point iterative solution process, monitor the increment of two adjacent fixed-point iterations, and construct a switching criterion for the fixed-point iterative stage; S3: Determine whether the fixed-point iteration result satisfies the preset convergence condition; if it does, obtain the global solution of the core sub-channel coupling model; otherwise, perform a switching criterion judgment: If the switching criterion does not meet the preset switching conditions, then continue to execute step S2; If the switching criterion satisfies the preset switching condition, the fixed-point iterative solution is terminated, and the current fixed-point iterative solution is used as the initial prediction value for the Jacobi-less Newton-Krylov solution process. S4: Construct a Jacobi-free Newton-Krylov solution process, introduce the fixed-point iterative operator as a physical preprocessing into the Jacobi-free Newton-Krylov solution framework, and construct the Jacobi-vector product approximation using the finite difference perturbation method. S5: In the process of solving the Newton-Krylov equation without Jacobi, the global nonlinear residual equation is linearized to obtain the corresponding Newton linear equation. The Krylov subspace iteration method is used, combined with the physical preprocessing and the Jacobi-vector product approximation, to solve the Newton linear equation and obtain the increment of the global solution vector. S6: Update the global solution vector according to the increment of the global solution vector, and determine whether the updated global solution vector satisfies the global convergence condition; if it satisfies the condition, obtain the global solution of the core subchannel coupling model and output the corresponding core thermal-hydraulic parameter distribution; if it does not satisfy the condition, return to step S5.

[0030] Preferably, the purpose of this embodiment is to organically combine fixed-point iteration and the Jacobi-free Newton-Krylov solution method. During computation, it leverages the advantage of quickly obtaining physically consistent global solutions in the early stages of fixed-point iteration, and switches to the Jacobi-free Newton-Krylov solution when the fixed-point iteration yields poor results, thus satisfying computational robustness and efficiency. The flowchart of this embodiment is shown below. Figure 1 As shown, the specific steps and implementation process are as follows: Step S1: Establish a fully implicit discrete core sub-channel coupling model In this embodiment, the core region is first divided into subchannel control volume grids based on the arrangement of fuel assemblies and the geometric topological relationships between fuel rods. Specifically, a subchannel division method centered on the coolant is adopted, discretizing the core into several axial layers along the axial direction. Adjacent axial layers are connected by axial flow, and adjacent subchannels are coupled by lateral flow and lateral pressure difference, thereby forming the core subchannel control volume grid.

[0031] On the sub-channel control volume grid, a set of core thermal-hydraulic control equations is established using a drift flow model. In this embodiment, the control equations include mixing energy equations, liquid phase energy equations, mixing continuity equations, axial momentum equations, and transverse momentum equations. For single-phase flow conditions, the vapor-liquid two-phase drift flow model can be degenerated into a single-phase form; for two-phase flow conditions, relevant constitutive relations such as vapor content, interphase slip, and drift velocity are retained to describe the vapor-liquid two-phase flow and heat transfer processes.

[0032] Furthermore, the aforementioned governing equations are discretized on the sub-channel control volume grid. During discretization, the unknowns to be solved in the current iteration step are implicitly represented, thereby establishing a fully implicit discrete equation system. The constitutive relations may include, but are not limited to, frictional resistance correlations, interphase slip relations, drift velocity relations, wall heat transfer correlations, phase transition relations, and mixing coefficient correlations. The boundary conditions may include core inlet flow boundary and outlet pressure boundary.

[0033] In this embodiment, the following main unknowns are selected for each sub-channel control body: enthalpy of mixing; axial mass flow rate; lateral mass flow rate; axial pressure difference; and lateral pressure.

[0034] Based on the order of the unknowns in all sub-channel control volumes, all unknowns are assembled into a global solution vector. Accordingly, all discrete control equations are reorganized into a unified nonlinear residual form, and a global nonlinear residual equation is established: , in, This represents the global nonlinear residual vector composed of all sub-channels, all axial layers, and the corresponding conservation equations.

[0035] To mitigate the adverse effects of magnitude differences between different physical quantities on the numerical solution process, this embodiment also performs dimensionless scaling on the global solution vector. Specifically, a diagonal scaling matrix is ​​used. For the original global solution vector Scaling is performed to obtain a dimensionless global solution vector. : , Among them, matrix It is a diagonal scaling matrix, whose diagonal elements can be derived from the components of the original global solution vector. The norm is determined. After scaling, the equivalent global nonlinear residual equation for the dimensionless variables can be obtained.

[0036] Step S2: Construct the fixed-point iterative solution process and fixed-point iterative operator To improve the quality of initial values ​​in the subsequent global unified solution process, a fixed-point iterative solution process is first constructed in this embodiment. A complete fixed-point iterative update of the core sub-channel coupling model is defined as the fixed-point iterative operator. The fixed-point iterative operator uses the current dimensionless global solution vector. As input, output the updated dimensionless global solution vector. ,Right now: , In this embodiment, the fixed-point iterative operator adopts a progressive scanning method along the axial direction, processing each axial layer sequentially from the core inlet to the core outlet. In each axial layer, the mixing enthalpy, axial mass flow rate, lateral mass flow rate, axial pressure difference, and lateral pressure corresponding to each sub-channel control body in that layer are updated based on the discrete control equations, boundary conditions, and lateral coupling relationships between adjacent sub-channels, to obtain the update result of the current axial layer. After completing the update of the current axial layer, the update continues for the next axial layer until the scanning of all axial layers of the entire core is completed, and the dimensionless global solution vector after one fixed-point iteration is output.

[0037] In some implementations, the fixed-point iterative operator updates the unknowns of each axial layer while simultaneously updating the relevant constitutive parameters based on the updated thermal state, so as to reflect the influence of changes in frictional resistance, interphase exchange, drift velocity, two-phase density, vapor content, and mixing characteristics on the coefficients of the coupling equation.

[0038] During the fixed-point iteration process, the increment between two adjacent fixed-point iterations is monitored, and a switching criterion for the fixed-point iteration stage is constructed accordingly. This allows for determining whether to switch to the Jacobi-less Newton-Krylov solution process based on the convergence characteristics of the fixed-point iteration process.

[0039] Step S3: Fixed-point iterative convergence judgment and switching criterion judgment In this embodiment, after each fixed-point iteration, it is first determined whether the current fixed-point iteration result meets the preset convergence condition. If the preset convergence condition is met, the current fixed-point iteration result is output as the global solution of the core sub-channel coupling model; if the preset convergence condition is not met, it is further determined whether to continue the fixed-point iteration or switch to the Jacobi-Newton-Krylov solution process based on the switching criterion.

[0040] Specifically, let the first The dimensionless global solution vector output by the next fixed-point iteration is Then the norm ratio of the increments of two adjacent fixed-point iterations will be... Defined as: , in, This represents the current iteration step number of the fixed-point iteration operator. , , The first , , The dimensionless global solution vector output by the next iteration. Representing vectors Norm.

[0041] In this embodiment, a switching threshold is set. and the maximum number of fixed-point iterations ,in In some implementations, A value of 0.7 to 0.9 is acceptable, with a preferred value of around 0.8. The appropriate value can be selected from 5 to 50 times, depending on the core size, the complexity of the operating conditions, and the calculation requirements.

[0042] The fixed-point iteration operator is used to continue iterative calculation when the following conditions are met: , The fixed-point iteration is terminated when any of the following conditions are met, and the current fixed-point iteration solution is used as the initial value for the Jacobi-less Newton-Krylov solution process: , In some implementations, the preset convergence condition for the fixed-point iteration can be set based on the residual norm, the norm of the difference between two adjacent iteration solution vectors, or a combination thereof.

[0043] Step S4: Construct a Jacobi-free Newton-Krylov solution process If the fixed-point iteration does not reach the preset convergence requirement and the switching condition is met, the current fixed-point iteration result is used as the initial value of the Jacobi-less Newton-Krylov solution process to construct a Jacobi-less Newton-Krylov solution framework for solving the global nonlinear residual equation.

[0044] In this embodiment, the Jacobi-Newton-Krylov method is used to solve the global nonlinear residual equation shown in the following equation: , To improve the efficiency of solving the subsequent Newton linear equations, in this embodiment, the fixed-point iterative operator constructed in step S2 is used. As a physical preprocessing step, a Jacobi-free Newton-Krylov solution framework is introduced. The fixed-point iterative operator reflects the main physical processes in the core subchannel coupling model, such as axial transport, lateral flow, and pressure coupling. Therefore, it can approximate the main coupling relationships in the original Newton linear equations, thereby improving the solution conditions of the linear system and increasing the solution efficiency of the Krylov subspace iterative method.

[0045] In some implementations, to facilitate the use of fixed-point iteration operators in constructing the physical preprocessed residual form, an equivalent preprocessing residual function can be defined: , in, To introduce the equivalent preprocessed residual function after the fixed-point iteration operator, This refers to the fixed-point iteration operator defined in step S2.

[0046] Furthermore, in this embodiment, the Jacobian matrix is ​​not explicitly constructed. Instead, a finite-difference perturbation method is used to construct the Jacobian-vector product approximation during the Krylov subspace iteration process. For any trial vector... The Jacobian-vector product can be approximated as follows: , in, This represents the Jacobian matrix of the residual function at the current iteration point. The perturbation parameters are selected adaptively based on machine accuracy, vector norm, and variable scale to balance approximation accuracy and numerical stability.

[0047] Using the equivalent preprocessing residual function In the implementation method, the preprocessed Jacobian-vector product is approximated using the following formula: , in, This represents the Jacobian matrix corresponding to the equivalent preprocessing residual function at the current iteration point.

[0048] Step S5: Linearization of residual equations and solution of local linear equations In this embodiment, for the first The second, uncompareable Newton-Krylov iteration, at the current iteration point Linearizing the global nonlinear residual equation yields the corresponding Newtonian linear equation: , in, This represents the Jacobian matrix of the residual function at the current iteration point. This represents the increment of the global solution vector to be solved in the current Newton iteration step.

[0049] Because core channel coupling models typically involve large scale, numerous unknowns, strong coupling relationships, and asymmetric Jacobian matrices, explicitly forming and storing the Jacobian matrix incurs high computational and storage costs. Therefore, in this embodiment, the Krylov subspace iteration method, combined with the aforementioned physical preprocessing and the Jacobian-vector product approximation, is used to solve the Newton linear equations, thereby obtaining the increment of the global solution vector. .

[0050] In this embodiment, the Krylov subspace iteration method can be the GMRES method; in other embodiments, iterative methods suitable for asymmetric linear systems, such as BiCGSTAB, can also be used. By introducing the physical preprocessing, the number of Krylov iteration steps can be reduced, the efficiency of solving linear equations can be improved, and the robustness of the overall nonlinear solution process can be enhanced.

[0051] Step S6: Update the global solution vector and determine global convergence. In obtaining the first The increment corresponding to the next Jacobian-Krylov iteration is incomparable. Then, the dimensionless global solution vector is updated. In this embodiment, to improve the stability and robustness of the nonlinear solution process, a line search strategy can be used to determine the step size factor. The dimensionless global solution vector is updated as follows: , Among them, step size factor It can be determined based on the residual descent criterion, Armijo condition, or other commonly used line search criteria.

[0052] After each update, the global nonlinear residual is recalculated. The algorithm then determines whether the global convergence condition is met. This global convergence condition can be set based on the residual norm, the relative change in the solution vector, or a combination of both. For example, the global unified solution can be considered complete when the following condition is met: , in For the first The residual of the second Newton iteration. For the initial residual, To achieve convergence tolerance. Desirable to The specific value can be set according to the core size, the complexity of the operating conditions, and the target calculation accuracy.

[0053] Once the global convergence condition is met, the dimensionless global solution vector is restored to its corresponding physical quantity form, yielding the final global solution of the core sub-channel coupling model, and the corresponding core thermal-hydraulic parameter distribution is output. The thermal-hydraulic parameter distribution may include flow distribution, pressure distribution, mixing enthalpy distribution, transverse flow distribution, and parameter distributions related to two-phase flow in each sub-channel.

[0054] To verify the advantages of this invention in terms of solution efficiency and computational overhead, the following specific embodiments are provided. Embodiment 1 is a problem of symmetrically partitioning 1 / 8 of the entire core into thermal subchannels, with the component partitioning as follows: Figure 2As shown, the axial mesh is divided into 44 layers. During the calculation, the relaxation factor in the fixed-point iteration stage is set to 0.6, and the convergence parameter in the Newton-Krylov solution stage is set to the corresponding preset value.

[0055] Figure 3 The residual convergence curves of the method described in this invention and the traditional fixed-point iterative method are shown in comparison. Figure 3 It can be seen that the residual decreases slowly throughout the solution process using the traditional fixed-point iterative method, requiring approximately several hundred nonlinear iterations to achieve convergence accuracy. The overall convergence curve exhibits a clear slow decay characteristic, remaining at a relatively high residual level even after more than 100 nonlinear iterations, and finally converging around the 700th iteration.

[0056] In contrast, the method described in this invention first pre-adjusts the global solution vector using fixed-point iteration in the initial stage. When the norm ratio of the fixed-point iteration increments meets a preset switching condition, it automatically switches to the Jacobi-free Newton-Krylov solution stage. Figure 3 As can be seen, this invention switches to the Newton-Krylov solution stage after completing a small number of fixed-point pre-iterations. After the switch, the residuals show a significant and rapid decrease, dropping rapidly from approximately 0.05 to the order of 10^(-7) around the 20th nonlinear iteration, demonstrating the quadratic convergence characteristic of the Newton method. The overall number of nonlinear iterations required is significantly less than that of the traditional fixed-point iteration method. This shows that using a fixed number of pre-iterations is insufficient to meet the solution requirements of different examples and computational conditions. This invention introduces an adaptive switching mechanism based on the incremental change relationship between two adjacent fixed-point iterations. This mechanism can determine the switching timing according to the actual convergence characteristics of the fixed-point iteration process, reducing unnecessary iterative calculations while ensuring the solution effect, thereby improving the overall solution efficiency and stability.

[0057] Figure 4 This paper presents a quantitative comparison between the present invention and the traditional fixed-point iterative method in terms of the number of residual function calls and computation time. Figure 4 As shown in the bar chart of (a), the traditional fixed-point iteration method calls the residual function 710 times in this example, while the method described in this invention only requires 233 calls, a reduction of approximately 67% compared to the traditional method. Figure 4 As can be seen from the bar chart in (b), the computation time of the traditional fixed-point iteration method is 840s, while the computation time of the method described in this invention is 395s, which is about 53% shorter than the traditional method.

[0058] The above comparison results show that by combining fixed-point iteration with the Jacobi-less Newton-Krylov method and introducing an adaptive switching mechanism based on the incremental change relationship between two adjacent fixed-point iterations, the present invention can switch to the Newton-Krylov solution stage with faster convergence speed in a timely manner after the initial adjustment of the fixed-point iteration. This significantly reduces the number of residual function calls while greatly reducing the computation time, and achieves efficient global unified solution of the core subchannel coupling model.

[0059] Those skilled in the art will readily understand that the above description is merely a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, and improvements made within the spirit and principles of the present invention should be included within the scope of protection of the present invention.

Claims

1. A method for solving a core subchannel coupling model, characterized in that, Includes the following steps: S1: Divide the reactor core into sub-channel control volume grids, construct a fully implicit discrete core sub-channel coupling model based on the control equations and constitutive equations of the drift flow model; construct a global solution vector according to the coupling relationship between the variables in the core sub-channel coupling model, and establish the corresponding global nonlinear residual equation; S2: Construct a fixed-point iterative solution process, define a fixed-point iterative update of the core sub-channel coupling model as a fixed-point iterative operator, and use the fixed-point iterative operator to iteratively update the global solution vector to obtain the predicted solution; During the fixed-point iterative solution process, the increment of two adjacent fixed-point iterations is monitored, and a switching criterion for the fixed-point iteration stage is constructed. The switching criterion is based on the norm ratio of the iteration increments of two adjacent fixed points. Perform construction and judgment: in, This represents the current iteration step number of the fixed-point iteration operator. , , The first , , The dimensionless global solution vector output by the next iteration. Represent the L2 norm of a vector; set the switching threshold. and the maximum number of fixed-point iterations ,in ;when and When, the fixed-point iteration operator is used to continue iterative calculation; when or When necessary, switch to the Jacobi-less Newton-Krylov method for subsequent global solution; S3: Determine whether the fixed-point iteration result meets the preset convergence condition; if it does, obtain the global solution of the core sub-channel coupling model; otherwise, perform a switching criterion judgment: If the switching criterion does not meet the preset switching conditions, then continue to execute step S2; If the switching criterion satisfies the preset switching condition, the fixed-point iterative solution is terminated, and the current fixed-point iterative solution is used as the initial prediction value for the Jacobi-less Newton-Krylov solution process. S4: Construct a Jacobi-free Newton-Krylov solution process, introduce the fixed-point iterative operator as a physical preprocessing into the Jacobi-free Newton-Krylov solution framework, and construct the Jacobi-vector product approximation using the finite difference perturbation method. S5: In the process of solving the problem without Jacobi Newton-Krylov, the global nonlinear residual equation is linearized to obtain the corresponding local linear equation. The Krylov subspace iteration method is used, combined with the physical preprocessing and the Jacobi-vector product approximation, to solve the local linear equation and obtain the increment of the global solution vector. S6: Update the global solution vector according to the increment of the global solution vector, and determine whether the updated global solution vector satisfies the global convergence condition; if it satisfies the condition, obtain the global solution of the core subchannel coupling model and output the corresponding core thermal-hydraulic parameter distribution; if it does not satisfy the condition, return to step S5.

2. The solution method for the core sub-channel coupling model according to claim 1, characterized in that, In step S1, the governing equations and constitutive equations of the drift flow model include the mixed energy equation, liquid phase energy equation, mixed continuity equation, axial momentum equation, and transverse momentum equation. These equations are discretized on the sub-channel control volume grid, and each term in the discrete equations is implicitly represented by the unknowns to be solved in the current iteration step to establish the fully implicit discrete equation set of the core sub-channel coupling model.

3. The solution method for the core sub-channel coupling model according to claim 2, characterized in that, In step S1, based on the coupling relationship between the variables in the core subchannel coupling model, the mixing enthalpy, axial mass flow rate, lateral mass flow rate, axial pressure difference and lateral pressure corresponding to each subchannel control body are assembled into a global solution vector.

4. The solution method for the core sub-channel coupling model according to claim 3, characterized in that, In step S1, the global solution vector is scaled to be dimensionless. The scaled dimensionless global solution vector is... satisfy: in, The original global solution vector, matrix It is a diagonal scaling matrix, whose diagonal elements are derived from the components corresponding to the original global solution vector. The norm is determined.

5. The method for solving the core sub-channel coupling model according to claim 2, characterized in that, In step S2, the fixed-point iterative operator takes the current global solution vector as input, and scans from the core inlet to the core outlet in a progressive manner along the axial direction. In each axial layer, it updates the mixing enthalpy, axial mass flow rate, lateral mass flow rate, axial pressure difference and lateral pressure corresponding to each sub-channel control body in sequence, and outputs the updated global solution vector.

6. The solution method for the core sub-channel coupling model according to claim 4, characterized in that, In step S6, let the first... The dimensionless global solution vector of the second-incomparable Newton-Krylov iteration is The corresponding increment is A line search strategy is used to determine the step size factor. The dimensionless global solution vector is then updated according to the following formula: in, That is, the first The dimensionless global solution vector of the second-order Jacobi Newton-Krylov iteration.

7. A solution system for a coupled model of a reactor core subchannel, characterized in that, It includes a processor and a memory, wherein the memory stores a computer program, and the processor executes the computer program to implement the solution method for the core subchannel coupling model as described in any one of claims 1 to 6.

8. A computer-readable storage medium having a computer program stored thereon, characterized in that, When the computer program is executed by the processor, it implements the solution method for the core subchannel coupling model as described in any one of claims 1 to 6.

9. The application of the solution method for the core sub-channel coupling model according to any one of claims 1 to 6, characterized in that, This is used for solving the core subchannel coupling model under single-phase flow and vapor-liquid two-phase flow conditions in the reactor core.