A discontinuous galerkin flux reconstruction method for solving thermal-hydraulic equations of nuclear reactors

Abstract

The application discloses a discontinuous Galerkin flux reconstruction method for solving a nuclear reactor thermal hydraulic equation, 1, pre-preparation; 2, spatial discretization; 3, format derivation; 4, define numerical flux; 5, time discretization; 6, posteriori error estimation. When the linear diffusion term in the control equation is processed, compared with the local discontinuous Galerkin (LDG) method, the temporary variable does not need to be introduced, the complexity of the numerical format is significantly reduced, the calculation time is saved, the growth rate of the calculation efficiency is positively correlated with the grid number and the piecewise polynomial degree, compared with the interior penalty discontinuous Galerkin (IPDG) method, the experience penalty parameter does not need to be introduced, and the method is more robust and physically accurate. In the case that the linear diffusion term is discretized by using the central flux and the odd-degree polynomial, the method has higher precision than the LDG method. The compatibility of the Dirichlet boundary condition of the method is stronger than that of the LDG method, the error asymptotic behavior of the solution shows convergence, and the method is suitable for long-time simulation of a well-posed problem which needs to strictly satisfy the physical boundary condition in engineering.

Description

Technical Field

[0001] This invention belongs to the field of nuclear reactor thermal-hydraulic analysis technology, specifically relating to a discontinuous Galerkin flux reconstruction method for solving nuclear reactor thermal-hydraulic equations. Background Technology

[0002] The thermal-hydraulic characteristics of nuclear reactors, including the coolant flow behavior in the primary loop and the temperature distribution of fuel elements, can be predicted using numerical simulation techniques. This numerical simulation relies on the numerical solution of partial differential equations in thermal-hydraulic systems, such as the mass, energy, and momentum conservation equations for fluids and the thermal conductivity equations for solids, all of which can be reduced to convection-diffusion equations. Traditional high-precision numerical discretization methods, such as the Discontinuous Galerkin (DG) method and higher-order finite volume methods, suffer from high computational complexity, high memory consumption, and insufficient boundary condition compatibility when solving convection-diffusion equations. The differential form of flux reconstruction (FR) based on the DG method offers higher computational efficiency but suffers from poor stability. The aforementioned DG and FR methods typically only discretize the first-order spatial derivative. For the discretization of the second spatial derivative, introducing temporary variables (such as Local Discontinuous Galerkin, LDG) will increase memory usage and computation time; introducing penalty terms will cause the stability and accuracy of the scheme to depend on empirical parameters, reducing the robustness and accuracy of the scheme (such as Interior Penalty Discontinuous Galerkin, IPDG). Summary of the Invention

[0003] To overcome the problems of existing numerical schemes, this invention proposes a discontinuous Galerkin flux reconstruction method (DG-FR) for solving the thermal-hydraulic equations of nuclear reactors. This method combines the high-order accuracy advantage of the DG method in handling convection terms with the FR method in simplifying the discretization of linear diffusion terms. It does not require the introduction of temporary variables or penalty terms, significantly reduces the complexity of the numerical scheme, saves computation time, and has good compatibility with Dirichlet boundary conditions.

[0004] To achieve the above objectives, the present invention adopts the following technical solution: A method for reconstructing discontinuous Galerkin fluxes in solving the thermal-hydraulic equations of a nuclear reactor includes the following steps: Step 1, Preliminary Preparation: First, define the boundary value problem to be solved, including the governing equations, the range of values ​​for spatial and temporal independent variables, initial conditions, and boundary conditions; then, perform spatial mesh generation, dividing the core solution domain into several discrete computational cells, and define the function spaces where the numerical solution and experimental function reside; define the numerical solution space as a piecewise polynomial space, which allows for jumps in the numerical solution at cell interfaces, and the experimental function space and the numerical solution space remain consistent; The governing equations solved by the discontinuous Galerkin flux reconstruction method are: (1) This is the mass conservation equation for the reactor core coolant, where... It is the density of the coolant, kg·m -3 ; It is the coolant velocity, m·s -1 ; It is time, s; (2) This is the energy conservation equation for the reactor core coolant, where... It is the total energy per unit mass of coolant, J·kg -1 ; It is the thermal conductivity of the coolant, W·m -1 ·K -1 ; It is the coolant temperature, in K; It is the gravitational force acting on a unit mass of coolant, N·kg -1 ; It is the static pressure of the coolant, in Pa; It is the viscous stress tensor, Pa; (3) This is the momentum conservation equation for the reactor core coolant; (4) This is the thermal conductivity equation for the fuel pellets, where... It is the density of fuel pellets, kg·m -3 ; It is the specific heat capacity at constant volume of fuel pellets, J·kg -1 ·K -1 ; It is the temperature of the fuel pellets, in K; It is the thermal conductivity of the fuel pellet, W·m -1 ·K -1 ; It is the power (W·kg) that generates heat per unit mass of fuel pellets. -1 The equation can still describe the heat conduction problem of the fuel cladding even after removing the source term.

[0005] The spatial independent variable takes values ​​within the specific solution region of the reactor core, while the temporal independent variable takes values ​​within the transient analysis duration.

[0006] Initial conditions include the initial velocity, pressure, spatial distribution of temperature of the coolant, and the initial temperature distribution of the fuel pellets and cladding.

[0007] Boundary conditions include thermal boundary conditions and velocity boundary conditions for the coolant, fuel pellets, and cladding.

[0008] Step 2, Spatial Discretization: Perform Galerkin variation on the defined control equations within the partitioned discrete computational units, that is, multiply both sides of the control equations by the experimental function and integrate within the discrete computational units. Use integration by parts to process the integrals of the convection and diffusion terms in the control equations to obtain the continuous weak form of the control equations; then replace the exact solution space with the numerical solution space and introduce numerical flux at the interface of the discrete computational units to obtain the semi-discrete weak form of the control equations. Step 3, Scheme Derivation: Introduce affine transformation to map discrete computational units to standard units, and then use Legendre polynomials to define basis functions for the numerical solution space; expand the numerical solution into a linear combination of basis functions, take the experimental function as the basis function, and substitute them into the semi-discrete weak form of the control equations; use the orthogonality of Legendre polynomials to simplify the integrals in the control equations. Step 4: Define Numerical Fluxes. Adjacent cells within the computational domain communicate only through numerical fluxes at their interfaces, representing weak coupling. In the governing equations, the numerical flux for the convection term is defined using conventional methods, such as upwind flux. This ensures the mobility of the discrete convection term and prevents non-physical oscillations under strong convection conditions with coolant. The numerical flux for the diffusion term, however, requires the introduction of the Flux Reconstruction (FR) method. First, Radau polynomials are introduced to reconstruct the numerical solutions on the computational cells, ensuring continuity of the reconstructed solutions at the cell interfaces. Then, the gradient of the reconstructed numerical solutions is calculated to define the numerical flux for the diffusion term. This numerical flux employs a central flux or alternating flux, similar to the LDG method, which ensures the stability of the discrete diffusion term. Step 5, Time Discretization: A third-order Strong Stability Preserving Runge-Kutta (SSPRK) explicit scheme is introduced to discretize the semi-discrete weak form of the governing equations, resulting in a fully discrete scheme. Solving this fully discrete scheme yields a high-precision numerical solution to the boundary value problem, namely the discrete velocity, pressure, and temperature fields of the core coolant, as well as the discrete temperature fields of the fuel pellets and cladding. The time step and spatial step in the fully discrete scheme must satisfy the CFL condition to ensure the stability of the explicit time-layer propagation. In actual calculations… The numerical solution of the time layer is obtained by approximating the initial conditions with the best square. The integrals in the numerical scheme are all calculated using the Gauss-Legendre quadrature formula. Step 6, Posterior error estimation: Define the root mean square error of equidistant solution points to verify the accuracy of the method, obtain the absolute error and convergence order of the numerical format, and thus determine whether the calculated numerical solution can accurately reflect the physical information of the reactor core.

[0009] Compared with the prior art, the present invention has the following advantages: When dealing with the linear diffusion term in the thermal-hydraulic equations of nuclear reactors, compared to the LDG method, the DG-FR method proposed in this invention does not require the introduction of temporary variables, significantly reducing the complexity of the numerical scheme and saving computation time. Furthermore, the growth rate of computational efficiency is positively correlated with the number of grid cells and the degree of the piecewise polynomials (see Tables 10 and 11). Compared to the IPDG method, the DG-FR method does not require the introduction of empirical penalty parameters, making it more robust and physically accurate. When the linear diffusion term is discretized using the central flux plus an odd-degree polynomial, the DG-FR method exhibits a higher convergence order than the LDG method (see Table 9). In addition, the DG-FR method shows stronger compatibility with Dirichlet boundary conditions than the LDG method, and its solution's asymptotic error behavior demonstrates convergence (see Table 9). Figures 2-5 This indicates that the DG-FR method is suitable for long-term simulation of well-determined problems in engineering that require strict satisfaction of physical boundary conditions. Attached Figure Description

[0010] Figure 1 This is a flowchart of the discontinuous Galerkin flux reconstruction method.

[0011] Figures 2-5 This is an error evolution graph of the viscous Burgers equation. Detailed Implementation

[0012] The present invention will now be described in detail with reference to the accompanying drawings and examples: like Figure 1 As shown, this invention proposes a method for reconstructing discontinuous Galerkin flux by solving the thermal-hydraulic equations of a nuclear reactor, comprising the following six steps: 1. Preparatory work; 2. Spatial discretization; 3. Scheme derivation; 4. Definition of numerical flux; 5. Temporal discretization; 6. Posterior error estimation. The specific implementation of the DG-FR method of this invention is illustrated using the conservation-type viscous Burgers equation as an example. This equation contains nonlinear convection terms and linear diffusion terms. After introducing pressure gradient and gravity, it can be generalized into a non-conservative one-dimensional fluid momentum conservation equation, namely, the one-dimensional form of equation (3), which can be used to describe the flow conditions of the core coolant. The specific boundary value problem is as follows: (5) In the formula It is the fluid velocity, m·s -1 ; It is a one-dimensional spatial coordinate, m; yes Domain; It is the maximum simulation duration, in seconds; It is the fluid velocity at time 0, in m·s. -1 ; yes Fluid velocity at [location], m·s -1 ; yes Fluid velocity at [location], m·s -1 ; It is space as The time is Fluid velocity at [location], m·s -1 .

[0013] The discontinuous Galerkin flux reconstruction method includes the following steps: Step 1, Preliminary Preparations: The boundary value problem is Equation (5). Mesh generation is first defined... One grid node: These nodes will Dissected Units: , yes The left endpoint; yes The right endpoint. Then the numerical solution space is defined as a piecewise polynomial space: ; Step 2, Spatial Discretization: Let... By performing Galerkin variation on the governing equations, replacing the exact solution space with a numerical solution space, and introducing numerical flux at the element interfaces, a semi-discrete weak form of the viscous Burgers equations is obtained. (6) In the formula It is the numerical flux of the convective term; It is the numerical flux of the diffusion term; It is the numerical solution of the fluid velocity, in m·s -1 ; It is the experimental function; yes exist The value at; yes exist The value at; yes exist The value at; yes exist The value at; yes exist The value at; yes exist The value at; Step 3, Scheme Derivation: Using Affine Transformation Discrete computing units Mapping to standard cells ,in yes midpoint, yes The length of the equation is then defined using Legendre polynomials. The numerical solution is expanded into a linear combination of the basis functions, and the experimental function is taken as the basis function. These are substituted into equation (6), and the integral is simplified using the orthogonality of Legendre polynomials. Further derivation yields: (7) In the formula yes Legendre polynomials of degree 1; yes Legendre polynomials of degree 1; yes exist The value at; yes exist The value at the location; unit Numerical solutions within Expanded to: In the formula It is in the numerical solution expansion Polynomial basis functions; yes exist The coefficient at time; yes exist The coefficient at time; yes right The derivative; yes right The derivative; yes right The derivative; Step 4: Define numerical flux: The upwind flux can be selected, as shown in the following expression: (8) In the formula yes exist The value at; yes exist The value at; This requires introducing the concept of flux reconstruction (FR) for definition. First, we need to define... within Reconstruction is performed to obtain a globally continuous numerical solution for fluid velocity. ,m·s -1 : (9) In the formula yes Right-hand Radau polynomial. yes Second-left Radau polynomials are defined as follows: (10) In the formula yes The second-order Legendre polynomial; then, with the help of This will complete the task. The numerical flux can be defined using a central flux or an alternating flux similar to the LDG method, as defined below: Central flux: (11) Alternating flux: (12) In the formula yes exist The value at; yes exist The value at; yes exist The value at; yes exist The value at; yes exist The value at; yes exist The value at; yes exist The value at; yes exist The value at; Step 5, Time Discretization: A third-order strong stability-preserving Runge-Kutta (SSPRK) explicit scheme is introduced to discretize (7) in time. Specifically, take... Substituting into (7), we get: (13) This formula is the first step of SSPRK. yes Time layer ; yes Time layer ; yes exist The value at; yes exist The value at; yes exist The value at; yes exist The value at the location; let (13) be denoted as Then the remaining two steps of SSPRK are: (14) This yields the fully discrete scheme of the viscous Burgers equation. In the equation... It is the time step; yes Time layer ; yes Time layer ; It is the first step of SPRK. ; It is the second step of SSPRK. ; yes The vector formed; yes The vector formed; yes The vector formed; It is an operator defined by (13); yes Acting on The result; yes Acting on The result; yes Acting on The result; and The CFL condition needs to be satisfied to ensure the stability of explicit time-level propagation. In actual calculations, Time layer By performing optimal square approximation on the initial conditions, the integrals in the numerical scheme can all be calculated using the Gauss-Legendre quadrature formula. Step 6, Posterior Error Estimation: Define the following root mean square error to verify the algorithm's accuracy: (15) In the formula It is the number of discrete solution points selected in each time layer, and these solution points are evenly distributed; It is the first Discrete solution points; It is space as The time is Precise fluid velocity at [location], m·s -1 ; It is space as , Time layer and spatial step size Numerical fluid velocity, m·s -1 ; yes Time layer and spatial step size The root mean square error; this error definition is only applicable to verifying simple boundary value problems with known analytical solutions. If the analytical solution to the problem is unknown, then the following root mean square error needs to be defined: (16) It can also be used to verify the accuracy of algorithms, where... It is space as , Time layer and spatial step size Numerical fluid velocity, m·s -1 .

[0014] After completing the numerical solution of the viscous Burgers equation, we can obtain the appropriate value using equation (15). The following error record table is sometimes available: Table 1. Central Flux Error Record Table (k=1) Table 2. Central Flux Error Record Table (k=2) Table 3. Central Flux Error Record Table (k=3) Table 4. Central Flux Error Record Table (k=4) Table 5 Alternating Flux Error Record Table (k=1) Table 6 Alternating Flux Error Record Table (k=2) Table 7 Alternating Flux Error Record Table (k=3) Table 8 Alternating Flux Error Record Table (k=4) By summarizing the contents of Tables 1-8, we can obtain the following summary table of convergence orders: Table 9 Summary of Convergence Orders As shown in Table 9, when the linear diffusion term is discretized using the central flux plus an odd-degree polynomial, the DG-FR method has a higher convergence order than the LDG method.

[0015] By comparison The runtime is used to compare the computational efficiency of the LDG and DG-FR methods, with the time unit being seconds. Each comparison group maintains the same spatiotemporal step size, number of grids, polynomial degree, and numerical flux type for both numerical formats, ensuring that the runtime depends only on the numerical format type. The computation time record table is as follows: Table 10 Record of Central Flux Calculation Duration Table 11 Record of Alternating Flux Calculation Duration As shown in Tables 10 and 11, compared with the LDG method, the DG-FR method proposed in this invention can significantly reduce the computation time of the numerical format, and the growth rate of computational efficiency is positively correlated with the number of grids and the degree of piecewise polynomials.

[0016] Figures 2-5 This is the error evolution graph for the viscous Burgers equation. It can be seen that when the solution time is long... Figure 3 The error of one tends to a small positive number, while the others tend to 0.

[0017] It is normal for the error to tend to zero because numerical methods themselves have numerical dissipation, which reduces the amplitude of the numerical solution; while the equations contain physical dissipation terms, which reduce the amplitude of the analytical solution. The physical dissipation and numerical dissipation are in the same direction and can work together to gradually reduce the error over time, eventually tending to zero. Figure 3 The error does not tend to 0 because the LDG scheme using alternating fluxes is not completely compatible with the Dirichlet boundary conditions, and the non-homogeneity of the error equation causes the error to tend to a non-zero steady-state value.

[0018] The comparison shows that in the DG-FR scheme, regardless of whether the diffusion term is discretized using the central flux or the alternating flux, the scheme is compatible with the Dirichlet boundary conditions. This indicates that the DG-FR method is superior to the LDG method in terms of boundary condition compatibility.

[0019] In summary, the advantages of this invention are that, when dealing with linear diffusion terms in the governing equations, compared to the Local Discontinuous Galerkin (LDG) method, it eliminates the need for temporary variables, significantly reducing the complexity of the numerical scheme and saving computation time. Furthermore, the growth rate of computational efficiency is positively correlated with the number of grid cells and the degree of the piecewise polynomials. Compared to the Interior Penalty Discontinuous Galerkin (IPDG) method, it eliminates the need for empirical penalty parameters, resulting in greater robustness and physical accuracy. When the linear diffusion terms are discretized using a central flux plus an odd-degree polynomial, this invention achieves higher accuracy than the LDG method. Moreover, this invention exhibits stronger compatibility with Dirichlet boundary conditions than the LDG method, and its solution's asymptotic error behavior demonstrates convergence, making it suitable for long-term simulations of well-posed problems in engineering that require strict adherence to physical boundary conditions. This invention provides theoretical and algorithmic support for the transition from low-order approximations to high-order accurate simulations in nuclear reactor thermal-hydraulic analysis.

Claims

1. A method for reconstructing discontinuous Galerkin fluxes in solving the thermal-hydraulic equations of a nuclear reactor, characterized in that: Includes the following steps: Step 1, Preliminary Preparation: First, define the boundary value problem to be solved, including the governing equations, the range of values ​​for spatial and temporal independent variables, initial conditions, and boundary conditions; then, perform spatial mesh generation, dividing the core solution domain into several discrete computational cells, and define the function spaces where the numerical solution and experimental function reside; define the numerical solution space as a piecewise polynomial space, which allows for jumps in the numerical solution at cell interfaces, and the experimental function space and the numerical solution space remain consistent; Step 2, Spatial Discretization: Perform Galerkin variation on the defined control equations within the partitioned discrete computational units, that is, multiply both sides of the control equations by the experimental function and integrate within the discrete computational units. Use integration by parts to process the integrals of the convection and diffusion terms in the control equations to obtain the continuous weak form of the control equations; then replace the exact solution space with the numerical solution space and introduce numerical flux at the interface of the discrete computational units to obtain the semi-discrete weak form of the control equations. Step 3, Scheme Derivation: Introduce affine transformation to map discrete computational units to standard units, and then use Legendre polynomials to define basis functions for the numerical solution space; expand the numerical solution into a linear combination of basis functions, take the experimental function as the basis function, and substitute them into the semi-discrete weak form of the control equations; use the orthogonality of Legendre polynomials to simplify the integrals in the control equations. Step 4: Define numerical flux: Adjacent cells within the computational domain communicate only through numerical fluxes at their interfaces, representing a weak coupling. In the governing equations, the numerical flux of the convection term ensures the mobility of the discrete convection term, preventing non-physical oscillations under strong convection conditions with coolant. The numerical flux of the diffusion term requires the introduction of the flux reconstruction method (FR) for definition. First, Radau polynomials are introduced to reconstruct the numerical solutions on the computational cells, ensuring the continuity of the reconstructed solutions at the cell interfaces. Then, the gradient of the reconstructed numerical solutions is calculated to define the numerical flux of the diffusion term. This numerical flux employs a central flux or alternating flux similar to the LDG method, which ensures the stability of the discrete diffusion term. Step 5, Time Discretization: A third-order Strong Stability Preserving Runge-Kutta (SSPRK) explicit scheme is introduced to discretize the semi-discrete weak form of the governing equations, resulting in a fully discrete scheme. Solving this fully discrete scheme yields a high-precision numerical solution to the boundary value problem, namely the discrete velocity, pressure, and temperature fields of the core coolant, as well as the discrete temperature fields of the fuel pellets and cladding. The time and space steps in the fully discrete scheme must satisfy the CFL condition to ensure the stability of the explicit time-layer propagation. In actual calculations… The numerical solution of the time layer is obtained by approximating the initial conditions with the best square. The integrals in the numerical scheme are all calculated using the Gauss-Legendre quadrature formula. Step 6, Posterior error estimation: Define the root mean square error of equidistant solution points to verify the accuracy of the method, obtain the absolute error and convergence order of the numerical format, and thus determine whether the calculated numerical solution can accurately reflect the physical information of the reactor core.

2. The method for reconstructing discontinuous Galerkin flux by solving the thermal-hydraulic equations of a nuclear reactor according to claim 1, characterized in that: In step 1, the governing equations solved by the discontinuous Galerkin flux reconstruction method are: (1) This is the mass conservation equation for the reactor core coolant, where... It is the density of the coolant, kg·m -3 ; It is the coolant velocity, m·s -1 ; It is time, s; (2) This is the energy conservation equation for the reactor core coolant, where... It is the total energy per unit mass of coolant, J·kg -1 ; It is the thermal conductivity of the coolant, W·m -1 ·K -1 ; It is the coolant temperature, in K; It is the gravitational force acting on a unit mass of coolant, N·kg -1 ; It is the static pressure of the coolant, in Pa; It is the viscous stress tensor, Pa; (3) This is the momentum conservation equation for the reactor core coolant; (4) This is the thermal conductivity equation for the fuel pellets, where... It is the density of fuel pellets, kg·m -3 ; It is the specific heat capacity at constant volume of fuel pellets, J·kg -1 ·K -1 ; It is the temperature of the fuel pellets, in K; It is the thermal conductivity of the fuel pellet, W·m -1 ·K -1 ; It is the power (W·kg) that generates heat per unit mass of fuel pellets. -1 The equation can still describe the heat conduction problem of the fuel cladding even after removing the source term.

3. The method for reconstructing discontinuous Galerkin flux by solving the thermal-hydraulic equations of a nuclear reactor according to claim 1, characterized in that: In step 1, the range of values ​​for the spatial independent variable is the specific solution region within the reactor core, and the range of values ​​for the temporal independent variable is the duration of the transient analysis.

4. The method for reconstructing discontinuous Galerkin flux by solving the thermal-hydraulic equations of a nuclear reactor according to claim 1, characterized in that: In step 1, the initial conditions include the initial velocity, pressure, spatial distribution of temperature of the coolant, and the initial temperature distribution of the fuel pellets and cladding.

5. The method for reconstructing discontinuous Galerkin flux by solving the thermal-hydraulic equations of a nuclear reactor according to claim 1, characterized in that: In step 1, the boundary conditions include thermal boundary conditions and velocity boundary conditions for the coolant, fuel pellets, and cladding.

Images