Stable and efficient global coupling iteration method for physical and thermal coupling simulation of reactor
The MFNK method addresses instability and inefficiency in reactor physical-thermal coupling by using power density as the coupled quantity and integrating nonlinear preprocessing and ML-gCMFD acceleration, achieving stable and efficient simulations.
Patent Information
- Application Number
- CN202410973151.9
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2024-07-19
- Publication Date
- 2024-12-13
AI Technical Summary
Current methods for reactor physical and thermal coupling simulation, such as the Picard iteration and JFNK method, suffer from instability, slow convergence, and inefficiency, particularly in high-fidelity simulations, due to issues like sensitivity to initial values, empirical parameter selection, and computational burden.
The MFNK method integrates a nonlinear residual equation, uses power density as the coupled physical quantity, employs nonlinear preprocessing, and combines the ML-gCMFD acceleration method to stabilize and accelerate the iteration process, reducing the computational burden and improving convergence.
The MFNK method achieves stable, fast, and accurate high-fidelity simulation of reactor physical-thermal coupling, with improved convergence and efficiency, enabling further extensions to reactor multi-physics field coupling.
Smart Images
Figure SMS_5 
Figure SMS_20 
Figure SMS_27
Abstract
Description
Technical Field
[0001] The invention relates to the field of nuclear research technology, and in particular to a stable and efficient global coupling iteration method for reactor physical and thermal coupling simulation. Background Art
[0002] At present, in response to the needs of reactor physical and thermal coupling simulation, the coupling solution methods available on the market mainly include the Picard iteration method and the Jacobian-Free Matrix Krylov (JFNK) method.
[0003] The Picard iteration method is based on the concept of decoupling, which decomposes the multi-physics coupling system into several sub-physics fields according to their physical meaning. In this method, each sub-physics field is solved alternately by independent solvers, and iteration is achieved by parameter transfer between solvers until all physics fields converge within one time step. Therefore, the Picard iteration method adds an outer loop and convergence criterion between different sub-physics fields.
[0004] The JFNK method is a global coupled iterative method, which is an organic combination of the Newton method, the Krylov subspace method and the Jacobian-Free method. This method first integrates the coupled system into a nonlinear system of equations, uses the Newton method for local linearization, and then uses the Krylov subspace method to solve the linearized system of equations. At the same time, combined with the Jacobian-Free method, the explicit construction and storage of the Jacobian matrix is avoided. Therefore, the JFNK method solves all physical quantities simultaneously, and has more advantages in stability and convergence characteristics.
[0005] The Picard iteration method faces the defects of poor stability, slow iteration convergence speed, and even iteration divergence. Especially for the case where the coupling of the coupled system is very strong or some physical quantities are very sensitive to the changes of other physical quantities, the Picard iteration may diverge. In addition, this method only has a first-order convergence speed, and the "over-solution" problem will make the calculation efficiency low. Finally, the Picard method adopts the idea of decoupling, which weakens the coupling between the physical fields and can only be used for coupled systems of a few physical fields.
[0006] The JFNK method is extremely challenging both in theoretical analysis and in practical application: ① The construction of the residual function is challenging, especially for high-fidelity physical-thermal coupling simulations, where the construction of the residual function may not even be possible; ② The physical quantities of different physical fields differ greatly in magnitude; ③ The convergence of the JFNK method is sensitive to the selection of the initial value of the Newton step. Improper selection of the initial value may lead to divergent solutions; ④ The selection of the differential step size in the JFNK method often relies on empirical formulas, which has a significant impact on the calculation accuracy and efficiency.
[0007] The instability of the Picard iterative method is caused by the influence of Doppler feedback, which makes the method unable to guarantee unconditional convergence. Secondly, in the initial iteration step, too many iterations are performed on a single physical field to converge, that is, an inaccurate problem is solved too accurately, thus causing the "over-solution" problem.
[0008] The JFNK method needs to describe each physical field in the form of a residual function under the global coupling framework, while mature single physical field programs do not directly provide residual functions. Therefore, it is necessary to reconstruct it according to the discrete relations of each sub-physical field. Secondly, the rows of the Jacobian matrix represent different physical equations, and the columns represent different physical quantities, which leads to significant numerical differences in the elements of the Jacobian matrix. Subsequently, since the JFNK method has the characteristics of local second-order convergence, the selection of the initial value of the Newton step needs to be within its local convergence interval. Finally, the JFNK method adopts the Jacobian-Free technology, using the differential step size to approximate the product of the Jacobian matrix and the basis vector, which makes the selection of the differential step size affect the calculation accuracy and convergence performance. In addition, the JFNK method also requires that all physical fields can and can only use the same differential step size.
[0009] For the Picard iterative method, although the use of relaxation factors can improve the stability of the algorithm, the selection of relaxation factors has a strong empirical nature, and with the increase in the complexity of the coupled system, it is more difficult to determine the appropriate relaxation factor. Secondly, the "over-solution" problem can be effectively alleviated by nonlinear alternating techniques (ANM) and local convergence techniques (PCT), but these methods are too rigid and may not be applicable to various multi-physics problems.
[0010] For the JFNK method, nonlinear preprocessing based on fixed point iteration can use the existing mature single physical field program to construct an equivalent residual function, but each Krylov subspace iteration requires calling each physical field program, which greatly reduces the computational efficiency. Secondly, although the use of Scaling technology can reduce the Jacobian matrix condition number, the selection of Scaling factors usually depends on experience and physical background, and the improvement effect of the condition number is limited.
[0011] In general, although the Picard iteration method and JFNK method commonly used in the market have been explored and optimized to a certain extent, they still have certain shortcomings and limitations. Summary of the invention
[0012] For reactor physics-thermal coupling simulation, the present invention provides a global coupling iteration method for reactor physics-thermal coupling simulation, namely, the MFNK method. This method integrates the coupling system into a nonlinear residual equation and simultaneously solves all coupled physical quantities, thereby ensuring high stability and convergence performance.
[0013] The present invention provides a global coupling iteration method for reactor physical thermal coupling simulation, which comprises the following steps:
[0014] Step 1: According to the input card provided by the user, read the geometry and material information, and use the constructive solid geometry method to build the neutron physics model and thermal hydraulic model. Among them, the neutron physics model is a refined neutron transport model, and the thermal hydraulic model is a single channel model and a cylindrical heat conduction model;
[0015] Step 2: The neutron physics model and the thermal hydraulic model use different grid division strategies to establish a spatial grid mapping relationship and realize the transmission of coupling parameters between different physical fields;
[0016] Step 3: Calculate the reactor physics-thermal coupling using the MFNK method;
[0017] Step 4: Through homogenization processing, a completely equivalent low-order ML-gCMFD model, a single-channel model, and a cylindrical thermal conduction model are established, thereby converting the high-order nonlinear residual system into an equivalent low-order nonlinear residual system;
[0018] Step 5: Before solving the low-order nonlinear residual equation, calculate the forcing factor η using the Eisenstat-Walker strategy k , and use the mixed difference step size formula to calculate ε;
[0019] Step 6: Use the inexact Newton method to solve the low-order nonlinear residual equation to obtain the Newton local linearization equation;
[0020] Step 7: Use the Krylov subspace method to solve the local linearized equations, and use the MGS-GMRES method to perform inner iteration until the convergence conditions of the iterative solution of the inexact Newton equation are met to obtain the change in power density.
[0021] Step 8: Update Power Density Then, each physical field is solved only once in sequence to obtain the new power density
[0022] Step 9: Using Power Density Calculate the root mean square RMS_PD and compare it with the set coupling convergence conditions;
[0023] If the convergence condition is met, the coupled system solution is completed; otherwise, return to step 4 to start the calculation of the next Newton step;
[0024] Step 10: Based on the converged power density, re-execute the thermal calculation to obtain converged thermal parameter information.
[0025] In a specific implementation, step 3: using the MFNK method to calculate the reactor physics-thermal coupling;
[0026] Specifically, the calculation process of neutron physics field and thermal hydraulic field is expressed as:
[0027]
[0028] In the formula, is the power density, which can be calculated from the neutron flux, and the superscript k is the Newton iteration step; F TH and F N_transport It is a conceptual solver, which contains the actual solver and all additional operations for processing the calculation results;
[0029] When equations (1-1) and (1-2) converge, the combination is simplified to an overall system:
[0030]
[0031] Furthermore, the high-order nonlinear residual equation of the physical-thermal coupling system is constructed:
[0032]
[0033] Among them, the MFNK method selects power density as the coupled physical quantity and iterates the initial value It is calculated by one Picard iteration;
[0034] In a specific implementation, step 4: using the ML-gCMFD coupled acceleration method to capture the impact of thermal feedback in an equivalent low-order system;
[0035] Specifically, through homogenization processing, completely equivalent low-order ML-gCMFD models, single-channel models and cylindrical thermal conduction models are established, and the ML-gCMFD calculation that is completely equivalent to neutron transport is used as the solver of the neutron physical field; thereby, an equivalent low-order nonlinear residual equation is established, which is as follows:
[0036]
[0037] The low-order nonlinear residual equation uses nonlinear preprocessing technology for linear systems The corresponding nonlinear system is:
[0038]
[0039] Decompose the coefficient matrix A into A = MN, and then the linear system Can be transformed into:
[0040]
[0041] In the formula, represents a Picard calculation and can be expressed as:
[0042]
[0043] Then, use M -1 As a left nonlinear preprocessor, the same solution transformation is performed on equation (5), and we can get:
[0044]
[0045] Combining equations (7) and (8), and using Represents the new residual equation after nonlinear preprocessing:
[0046]
[0047] In a specific implementation, step 5: before solving the low-order nonlinear residual equation, complete the following parameter preparation: use the Eisenstat-Walker strategy to calculate the forcing factor η k and use the mixed difference step calculation formula to obtain ε;
[0048] Specifically, the mixed difference step length calculation formula is as follows:
[0049]
[0050] in, is the basis vector, ε mach is the calculation accuracy of the computer, s k is the increment of the kth Newton step, and η is the relative deviation during the “black box” calculation;
[0051] In a specific implementation, in step 6, the low-order nonlinear residual equation Perform Taylor expansion and ignore higher-order terms as follows:
[0052]
[0053] Further simplifying the process, we obtain the Newton local linearization equation as follows:
[0054]
[0055] In a specific implementation, step 7: using the MGS-GMRES method to perform inner iterative solution until the convergence condition of the iterative solution of the inexact Newton equation is met:
[0056]
[0057] Calculate the change in power density In this process, on the one hand, the construction and storage of the Jacobian matrix are difficult, and on the other hand, MGS-GMRES only needs to calculate the product of the Jacobian matrix and the basis vector when solving. For this reason, the Jacobian-Free method is used in the MFNK method, that is, the matrix-vector product is approximated by finite differences, as shown below:
[0058]
[0059] in, is the Krylov basis vector, ε is the mixed difference step size;
[0060] In a specific implementation, in step 8: based on the change in power density calculated by GMRES pass to update the power density, and then the updated power density As input, each physics field is then solved only once in sequence to obtain the new power density
[0061] In a specific implementation, in step 9, the calculation formula of the convergence standard set by the coupling system is as follows:
[0062]
[0063] Where N is the residual equation dimension, and atol is the convergence criterion set by the user.
[0064] The present invention thus provides a global coupled iterative system for reactor physical and thermal coupled simulation based on the method, which includes modules corresponding to steps 1 to 10.
[0065] The present invention also provides an electronic device, comprising: at least one processor; and a memory communicatively connected to the at least one processor;
[0066] The memory stores instructions executable by the at least one processor, and the instructions are configured to execute the above method or the system.
[0067] The present invention also provides a computer-readable storage medium storing computer-executable instructions, wherein the computer-executable instructions are used to execute the method flow or execute the global coupling iterative system of the reactor physical-thermal coupling simulation.
[0068] The MFNK method of the present invention eliminates the defects of the JFNK method through the following four strategies:
[0069] ① Different from the JFNK method which uses neutron flux and temperature as coupled physical quantities, the MFNK method selects power density as the intermediate solution variable. This not only reduces the computational scale of the nonlinear equation, but also because power density is used as a single coupled physical quantity, the condition number of the Jacobian matrix is smaller, the convergence is faster, and there is no need to rely on Scaling technology;
[0070] ② The nonlinear preprocessing technology is simple to implement and can achieve "black box" coupling, that is, directly use the existing mature single physical field solver to construct the residual equation. At the same time, combined with the ML-gCMFD coupling acceleration method, the influence of thermal feedback is captured through a completely equivalent low-order system, thereby cleverly avoiding the computational burden caused by frequent calls to neutron transport calculations during Krylov subspace iterations, and significantly improving the computational efficiency.
[0071] ③ When using the ML-gCMFD coupled acceleration method, a Picard calculation of the nonlinear alternating technique needs to be performed for each Newton step. This process significantly reduces the sensitivity of the MFNK method to the initial value of the Newton step and enhances stability;
[0072] ④The MFNK method adopts a hybrid difference step calculation formula, which combines the advantages of the empirical difference step and the optimal difference step, further improving the calculation efficiency and convergence speed.
[0073] The effect that can be achieved by the present invention is that the MFNK method can realize stable, fast, accurate and high-resolution simulation of high-fidelity physical-thermal coupling of the entire core. Furthermore, the present invention also obtains the following deeper effect: based on the MFNK method proposed by the present invention, it can be further extended to the field of reactor multi-physics field coupling simulation, such as increasing the coupling of fuel performance, etc. BRIEF DESCRIPTION OF THE DRAWINGS
[0074] Figure 1 Flowchart of the implementation of the MFNK method.
[0075] Figure 2Schematic diagram of the axial and radial structures of the VERA 3D single-element problem.
[0076] Figure 3 Axial relative power distribution under each coupling method.
[0077] Figure 4 VERA 3D single element problem using the convergence history of each coupled method.
[0078] Figure 5 Comparison of computational time for a VERA 3D single element problem using various coupling methods. DETAILED DESCRIPTION
[0079] In order to better understand the technical solution of the present invention, the embodiments of the present invention are described in detail below with reference to the accompanying drawings.
[0080] It should be clear that the described embodiments are only some embodiments of the present invention, not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.
[0081] The terms are explained as follows:
[0082] JFNK: Jacobian-Free Newton Krylov, no Jacobian Newton-Krylov
[0083] ANM: Alternating Nonlinear Method
[0084] PCT: the Partial-Convergence Technique, local convergence technology
[0085] ML-gCMFD: the Multi-Level Generalized equivalence theory based CoarseMesh Finite Difference method, multi-level coarse mesh finite difference method based on generalized equivalence theory
[0086] MGS-GMRES: Generalized Minumum RESidual based on the Modified Gram-Schmit, a generalized minimum residual algorithm based on the modified Gram-Schimit orthogonalization method.
[0087] The present invention proposed an initial MFNK method in the early stage, which combined the inexact Newton method, the Krylov subspace method and the Jacobian-Free method, and adopted the optimal differential step calculation formula proposed by Xu of Purdue University. Although the implementation of the initial MFNK method is relatively simple, and it is superior to the traditional Picard iteration method with relaxation factors in terms of convergence performance and efficiency. However, this method also has certain limitations: ① The neutron transport calculation module needs to be frequently called during the Krylov subspace iteration process, which greatly increases the calculation burden and leads to low calculation efficiency; ② The selection of the initial value of the Newton step still has a large impact on the convergence; ③ Compared with the traditional empirical differential step calculation formula, although the use of the optimal differential step calculation formula can improve the convergence performance, its application has limitations. Especially in the initial Newton iteration step, the optimal differential step has the defect of being too large. Therefore, through further research, an implementation method of the present invention is proposed.
[0088] Example 1: Global coupling iteration method for reactor physical and thermal coupling simulation
[0089] The specific steps are as follows:
[0090] Step 1: According to the input card provided by the user, read the geometry and material information, and use the constructive solid geometry method to build the neutron physics model and thermal hydraulic model. Among them, the neutron physics model is a refined neutron transport model, and the thermal hydraulic model is a single channel model and a cylindrical heat conduction model;
[0091] Step 2: The neutron physics model and the thermal hydraulic model use different grid division strategies to establish a spatial grid mapping relationship and realize the transmission of coupling parameters between different physical fields;
[0092] Step 3: Calculate the reactor physics-thermal coupling using the MFNK method. In this embodiment, for the reactor physics-thermal coupling system, the calculation process of the neutron physical field and the thermal hydraulic field can be expressed as:
[0093]
[0094] In the formula, is the power density, which can be calculated from the neutron flux, and the superscript k is the Newton iteration step. TH and F N_transport It is a conceptual solver that contains the actual solver and all additional operations for processing the calculation results.
[0095] When equations (1-1) and (1-2) converge, they can be simplified into an overall system:
[0096]
[0097] Furthermore, the high-order nonlinear residual equation of the physical-thermal coupling system constructed in this embodiment is:
[0098]
[0099] It can be seen that the MFNK method selects power density as the coupled physical quantity. It should be noted that the initial value of the iteration It is calculated by one Picard iteration;
[0100] Step 4: This embodiment uses the ML-gCMFD coupling acceleration method to capture the impact of thermal feedback in an equivalent low-order system. Specifically, through homogenization, a completely equivalent low-order ML-gCMFD model, a single-channel model, and a cylindrical thermal conduction model are established, and the ML-gCMFD calculation that is completely equivalent to neutron transport is used as the solver of the neutron physics field. At this point, an equivalent low-order nonlinear residual equation is established, which is as follows:
[0101]
[0102] The low-order nonlinear residual equation uses nonlinear preprocessing techniques. Specifically, for the linear system The corresponding nonlinear system is:
[0103]
[0104] Decompose the coefficient matrix A into A = MN, and then the linear system Can be transformed into:
[0105]
[0106] In the formula, represents a Picard calculation and can be expressed as:
[0107]
[0108] Then, use M -1 As a left nonlinear preprocessor, the same solution transformation is performed on equation (5), and we can get:
[0109]
[0110] Combining equations (7) and (8), and using Represents the new residual equation after nonlinear preprocessing:
[0111]
[0112] It can be seen that the format of the new residual equation after nonlinear preprocessing is consistent with the low-order nonlinear residual equation (4). Therefore, the MFNK method in this embodiment applies nonlinear preprocessing technology, which has the following significant improvements and advantages: ① It is simple to implement and can realize "black box" coupling; ② Combined with the ML-gCMFD coupling acceleration method, it cleverly avoids the computational burden caused by a large number of calls to the neutron transport calculation module during the iteration of the Krylov subspace method; ③ The nonlinear preprocessing technology has good accelerated convergence performance;
[0113] Step 5: Before solving the low-order nonlinear residual equation, a series of parameter preparations need to be completed. Specifically, the forced factor η is calculated using the Eisenstat-Walker strategy. k And use the hybrid difference step calculation formula to obtain ε. Considering that the Jacobian-Free method approximates the Jacobian matrix-vector product by using finite differences, the selection of the difference step ε is crucial to the calculation accuracy and convergence performance. The present invention proposes a hybrid difference step calculation formula, which combines the respective advantages of the empirical difference step calculation formula of the NITSOL package and the optimal difference step calculation formula proposed by Xu of Purdue University, ensuring the suitability of the difference step in the entire coupled iterative process. Among them, the hybrid difference step calculation formula is as follows:
[0114]
[0115] in, is the basis vector, ε mach is the calculation accuracy of the computer, s k is the increment of the kth Newton step, and η is the relative deviation during the “black box” calculation;
[0116] Step 6: Use the inexact Newton method to solve the low-order nonlinear residual equation. Perform Taylor expansion and ignore higher-order terms as follows:
[0117]
[0118] Further simplifying the process, we obtain the Newton local linearization equation as follows:
[0119]
[0120] Thus, the complex nonlinear residual equation is converted into a linear system of equations, so that the linear system of equations can be solved by the Krylov subspace method.
[0121] Step 7: Use the Krylov subspace method to solve the Newton local linearization equation. Specifically, the MGS-GMRES method is used to perform inner iteration until the convergence condition of the iterative solution of the inexact Newton equation is met:
[0122]
[0123] Finally, the change in power density is calculated In this process, on the one hand, the construction and storage of the Jacobian matrix are difficult, and on the other hand, MGS-GMRES only needs to calculate the product of the Jacobian matrix and the basis vector when solving. For this reason, the Jacobian-Free method is used in the MFNK method, that is, the matrix-vector product is approximated by finite differences, as shown below:
[0124]
[0125] in, is the Krylov basis vector, ε is the mixed difference step size;
[0126] To better understand the above relationship, refer to Figure 1 First, in order to solve the Newton local linearization equation, that is, formula (12), we use the MGS-GMRES method, which is an iterative method. Therefore, it is necessary to set the convergence condition of the iterative solution, that is, formula (13). When we solve When formula (13) is satisfied, the solution of Newton's local linearization equation is considered complete. According to what I just mentioned, the Krylov subspace method (that is, the MGS-GMRES method) is used to solve it. According to the principle of the MGS-GMRES method, it is necessary to calculate the product of the matrix and the vector in the solution process, that is, This product can be calculated using formula (14).
[0127] Step 8: Change in power density calculated based on GMRES pass to update the power density, and then the updated power density As input, each physics field is then solved only once in sequence to obtain the new power density
[0128] Step 9: Determine the convergence of the MFNK solution. Specifically, use the power density Calculate the root mean square RMS_PD and compare it with the set coupling convergence condition. If the convergence condition is met, the coupled system solution is completed. Otherwise, return to step 4 to start the calculation of the next Newton step. The calculation formula of the convergence standard set for the coupled system is as follows:
[0129]
[0130] Where N is the residual equation Dimension, atol is the convergence standard set by the user;
[0131] Step 10: Based on the converged power density, re-execute the thermal calculation to obtain converged thermal parameter information.
[0132] Example 2: Verification of the method of the present invention
[0133] Based on the VERA#6 benchmark problem provided by the CASL program of the U.S. Department of Energy, a VERA 3D single-element steady-state physics-thermal coupling problem is constructed. Figure 2 The schematic diagram of the structure of the 3D single-cell problem is shown. In the radial direction, the fuel cell includes UO2 fuel pellets, air gaps, cladding and moderators. At the same time, the core active area is divided into 46 layers in the axial direction: 8.2111, 8.2112, 8.2111, 8.21112, 8.2111, 8.2112, 8.2111, 8.065 (30 layers), 7.9212 (5 layers) and 6.6815cm (4 layers). In addition, there is a 30cm reflection layer above and below. Vacuum boundary conditions are used in the axial direction, and total reflection boundary conditions are used in the radial direction. This problem is calculated using one kernel. Table 1 gives the relevant coupling parameters of this problem.
[0134] Table 1. VERA 3D single-element problem coupling parameters
[0135]
[0136] For the VERA 3D single cell problem, the Picard iteration method, the JFNK method and the MFNK method of the present invention are selected for simulation. At the same time, when the root mean square of the power density RMS_PD≤10 -5 The purpose is to test the numerical performance of the MFNK method of the present invention in terms of computational accuracy, stability, convergence performance and computational efficiency.
[0137] Table 2 and Figure 3 The calculation results of the main parameters of each coupling method and the comparison of the axial relative power distribution are given respectively. Figure 4 and Figure 5 The convergence history and computational time of various coupling methods are shown separately.
[0138] Table 2. Calculation results of solving VERA 3D single-element problem using various coupling methods
[0139]
[0140] In terms of calculation accuracy, according to Table 2 and Figure 3 It can be seen that the calculation results of the MFNK method of the present invention are very consistent with those of other coupling methods. At the same time, after taking into account the thermal feedback, the axial relative power distribution of various coupling methods shows good consistency, and the distribution becomes asymmetric, and the relative power peak moves to the lower half of the core, which further confirms the correctness of the MFNK method of the present invention.
[0141] In terms of convergence performance and stability, Figure 4 It can be seen that the MFNK method of the present invention has the fastest convergence speed, followed by the JFNK method, and the Picard method converges the slowest. In addition, the Picard method may have convergence divergence in some cases. Only when the relaxation factor is selected in the range of 0.2 to 0.8 can the convergence of the coupling calculation be guaranteed. In contrast, the MFNK method of the present invention shows higher stability.
[0142] In terms of computational efficiency, Figure 4 As shown in the figure, the calculation time of the MFNK method of the present invention is lower than that of the JFNK method, indicating that the calculation efficiency of the MFNK method is higher. Although in some cases, the calculation time of the MFNK method of the present invention is slightly higher than that of the Picard iterative method, this depends on the selection of the optimal relaxation factor, which is highly empirical, and the relaxation factor has a significant impact on the calculation efficiency of the Picard method.
[0143] In summary, the MFNK method of the present invention has high convergence performance, strong stability and relatively high computational efficiency while ensuring calculation accuracy. It is a stable and efficient global coupled iterative method.
Claims
1. A global coupling iterative method for reactor physical and thermal coupling simulation, characterized in that: The steps include: Step 1: According to the input card provided by the user, read the geometry and material information, and use the constructive solid geometry method to build the neutron physics model and thermal hydraulic model. Among them, the neutron physics model is a refined neutron transport model, and the thermal hydraulic model is a single channel model and a cylindrical heat conduction model; Step 2: The neutron physics model and the thermal hydraulic model use different grid division strategies to establish a spatial grid mapping relationship and realize the transmission of coupling parameters between different physical fields; Step 3: Calculate the reactor physics-thermal coupling using the MFNK method; Step 4: Through homogenization processing, a completely equivalent low-order ML-gCMFD model, a single-channel model, and a cylindrical thermal conduction model are established, thereby converting the high-order nonlinear residual system into an equivalent low-order nonlinear residual system; Step 5: Before solving the low-order nonlinear residual equation, calculate the forcing factor η using the Eisenstat-Walker strategy k , and use the mixed difference step size formula to calculate ε; Step 6: Use the inexact Newton method to solve the low-order nonlinear residual equation to obtain the Newton local linearization equation; Step 7: Use the Krylov subspace method to solve the local linearized equations, and use the MGS-GMRES method to perform inner iteration until the convergence conditions of the iterative solution of the inexact Newton equation are met to obtain the change in power density. Step 8: Update Power Density Then, each physical field is solved only once in sequence to obtain the new power density Step 9: Using Power Density Calculate the root mean square RMS_PD and compare it with the set coupling convergence conditions; If the convergence condition is met, the coupled system solution is completed; otherwise, return to step 4 to start the calculation of the next Newton step; Step 10: Based on the converged power density, re-execute the thermal calculation to obtain converged thermal parameter information.
2. The method according to claim 1, characterized in that Step 3: Calculate the reactor physics-thermal coupling using the MFNK method; Specifically, the calculation process of neutron physics field and thermal hydraulic field is expressed as: In the formula, is the power density, which can be calculated from the neutron flux, and the superscript k is the Newton iteration step; F TH and F N_transport It is a conceptual solver, which contains the actual solver and all additional operations for processing the calculation results; When equations (1-1) and (1-2) converge, the combination is simplified to an overall system: Furthermore, the high-order nonlinear residual equation of the physical-thermal coupling system is constructed: Among them, the MFNK method selects power density as the coupled physical quantity and iterates the initial value It is calculated through a Picard iteration.
3. The method according to claim 1, characterized in that Step 4: Use the ML-gCMFD coupled acceleration method to capture the impact of thermal feedback in an equivalent low-order system; Specifically, through homogenization processing, completely equivalent low-order ML-gCMFD models, single-channel models and cylindrical thermal conduction models are established, and the ML-gCMFD calculation that is completely equivalent to neutron transport is used as the solver of the neutron physical field; thereby, an equivalent low-order nonlinear residual equation is established, which is as follows: The low-order nonlinear residual equation uses nonlinear preprocessing technology for linear systems The corresponding nonlinear system is: Decompose the coefficient matrix A into A = MN, and then the linear system Can be transformed into: In the formula, represents a Picard calculation and can be expressed as: Then, use M -1 As a left nonlinear preprocessor, the same solution transformation is performed on equation (5), and we can get: Combining equations (7) and (8), and using Represents the new residual equation after nonlinear preprocessing:
4. The method according to claim 1, characterized in that Step 5: Before solving the low-order nonlinear residual equation, complete the following parameter preparation: Calculate the forcing factor η using the Eisenstat-Walker strategy k and use the mixed difference step calculation formula to obtain ε; Specifically, the mixed difference step length calculation formula is as follows: in, is the basis vector, ε mach is the calculation accuracy of the computer, s k is the increment of the kth Newton step, and η is the relative deviation during the "black box" calculation.
5. The method according to claim 1, characterized in that In step 6, the low-order nonlinear residual equation Perform Taylor expansion and ignore higher-order terms as follows: Further simplifying the process, we obtain the Newton local linearization equation as follows:
6. The method according to claim 1, characterized in that Step 7: Use the MGS-GMRES method to perform inner iterative solution until the convergence condition of the iterative solution of the inexact Newton equation is met: Calculate the change in power density in, The calculation formula is as follows: in, is the Krylov basis vector, and ε is the mixed difference step size.
7. The method according to claim 1, characterized in that In step 8: The change in power density calculated based on GMRES pass to update the power density, and then the updated power density As input, each physics field is then solved only once in sequence to obtain the new power density 8. The method according to claim 1, characterized in that In step 9, the calculation formula for the convergence criterion set for the coupled system is as follows: Where N is the residual equation dimension, and atol is the convergence criterion set by the user.
9. A global coupled iterative system for reactor physical and thermal coupled simulation based on the method according to any one of claims 1 to 8, comprising modules corresponding to steps 1 to 10.
10. An electronic device, characterized in that: include: at least one processor; and, a memory communicatively coupled to the at least one processor; The memory stores instructions executable by the at least one processor, wherein the instructions are configured to execute the method according to any one of claims 1 to 8 or to execute the system according to claim 9.