CN115828666A - High-precision numerical simulation method applied to lead-based stack SGTR accident analysis - Google Patents

High-precision numerical simulation method applied to lead-based stack SGTR accident analysis Download PDF

Info

Publication number
CN115828666A
CN115828666A CN202211400808.XA CN202211400808A CN115828666A CN 115828666 A CN115828666 A CN 115828666A CN 202211400808 A CN202211400808 A CN 202211400808A CN 115828666 A CN115828666 A CN 115828666A
Authority
CN
China
Prior art keywords
analysis program
equation
energy
modeling
lead
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
CN202211400808.XA
Other languages
Chinese (zh)
Inventor
张斌
杨文强
王仁杰
缪凡
单建强
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Xian Jiaotong University
Original Assignee
Xian Jiaotong University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Xian Jiaotong University filed Critical Xian Jiaotong University
Priority to CN202211400808.XA priority Critical patent/CN115828666A/en
Publication of CN115828666A publication Critical patent/CN115828666A/en
Pending legal-status Critical Current

Links

Images

Abstract

The invention discloses a high-precision numerical simulation method applied to lead-based stack SGTR accident analysis, which utilizes the advantage of accurate simulation of an integrated serious accident analysis program on related problems of a pipeline to model a pipeline part of a system, and adopts a reactor core disintegration analysis program which can accurately analyze the interaction problems of phase change and multi-fluid to model in a region with complex interaction during the accident, and the two programs are coupled through a coupling component to realize the accurate analysis capability on the lead-based stack SGTR accident. The invention can fully utilize respective advantages of the integrated serious accident analysis program and the reactor core disintegration analysis program, fully reflect the geometric characteristics of a modeling object, and refine a modeling area, thereby improving the precision of SGTR accident simulation of the lead-based reactor.

Description

High-precision numerical simulation method applied to lead-based stack SGTR accident analysis
Technical Field
The invention belongs to the technical field of nuclear reactor thermal hydraulic power, and particularly relates to a high-precision numerical simulation method applied to lead-based reactor SGTR accident analysis.
Background
The lead-based stack steam generator heat transfer tube breakage (SGTR) accident is a typical lead-based stack accident, and in a steam generator, because there is only one wall separation between liquid lead and water, and there are large pressure difference and temperature difference between the liquid lead and the water, the large mechanical stress and thermal stress easily cause the reactor steam generator to generate vibration and corrosion effect in the operation process, so the heat transfer tubes of a primary loop of the steam generator have the risk of accident occurrence. When the lead-based reactor has an SGTR accident, high-pressure water on the secondary side can be injected into the primary side, so that the high-pressure water is in direct contact with a low-pressure high-temperature heavy metal coolant, and the water can be quickly evaporated at the moment and is accompanied with the generation of pressure waves, which possibly threatens the safety of the reactor. The following phenomena are mainly involved in the SGTR accident and can adversely affect the reactor: the first is the interaction between the coolants, which causes steam to migrate to the core, causing a positive reaction. The second phenomenon is the generation of pressure waves, which can have an effect on the structural integrity of the pressure vessel. The third phenomenon is the mechanical impact generated by the heavy metal liquid, which can cause the shaking of the liquid lead/lead bismuth pool. Therefore, in the design and safety analysis of lead-based stacks, the analysis of SGTR accidents should be emphasized. Accurate simulation is carried out on the complicated accident process and phenomenon of the lead-based stack SGTR, and then correct overall safety assessment is always difficult to study.
Disclosure of Invention
The invention aims to solve the technical problem of providing a numerical simulation method for comprehensively utilizing the advantages of an integrated serious accident analysis program and a reactor core disintegration analysis program and coupling the advantages of the integrated serious accident analysis program and the reactor core disintegration analysis program for the SGTR accident of a lead-based reactor, fully utilizing the advantages of the integrated serious accident analysis program and the reactor core disintegration analysis program, simultaneously considering the analysis of multiphase and multicomponent complex phenomena and pipeline flowing behaviors, carrying out more accurate modeling calculation for the SGTR accident of the lead-based reactor, and improving the accuracy of the accident analysis.
The invention adopts the following technical scheme:
a high-precision numerical simulation method applied to lead-based reactor SGTR accident analysis is characterized in that an accident area is divided into a pipeline area and an area (generally near a steam generator break or even a reactor core part) with complex interaction, an integrated serious accident analysis program is used for modeling the pipeline area, the integrated serious accident analysis program can establish a control equation for a modeled controller, a semi-implicit method is adopted for solving, meanwhile, the reactor core disintegration analysis program is used for modeling the area with the complex interaction, the reactor core disintegration analysis program can divide the components of a material, a conservation equation is established for the divided components, a four-step method is used for solving the conservation equation, a coupling assembly is used for realizing data transmission between the two, and modeling and analysis of the SGTR accident are completed.
The method comprises the steps of firstly, modeling a pipeline part of a system by using an integrated serious accident analysis program, wherein the integrated serious accident analysis program adopts a mode of dividing control bodies and connecting flow channels to complete modeling of the pipeline, the two control bodies are connected by a flow channel, the retention time of fluid is not considered in the flow channel, and the mass or energy in all the flow channels is not considered, namely the mass and the energy of the fluid are stored in the control bodies. The controller distinguishes 'liquid' and 'gas' parts in the body, representing liquid phase and gas phase, and adopts two fluid models to respectively use mass, momentum and energy conservation laws for the gas phase and the liquid phase, wherein a mass conservation constant differential equation of a component m in a controller i can be expressed by the following formula:
Figure BDA0003934870120000021
in the formula: m i,m -controlling the total mass/kg of component m in body i, the angle m representing the water, the liquid droplet or various gas item components, j-deFlow control body, σ ij -the direction of the flow of the gas,
Figure BDA0003934870120000022
to control body j
Figure BDA00039348701200000210
The fraction of the phase(s) is,
Figure BDA00039348701200000211
is a gas phase or a liquid phase in the flow channel, t-time/s,
Figure BDA0003934870120000023
the density of component m in the defluidizing control volume j and the inflowing control volume d/kg m -3 And d represents the incoming flow control volume,
Figure BDA0003934870120000024
-in the defluidizing control body j
Figure BDA00039348701200000212
Flow velocity of phase/m.s, A j Flow area/m of the defluidizing control body j 2 ,F j The opening of the flow area of the defluidizing control body j,
Figure BDA0003934870120000025
controlling the non-flowing mass source term/kg-s of component m in body i -1
The energy conservation ordinary differential equation for each phase can be expressed by the following equation:
Figure BDA0003934870120000026
in the formula:
Figure BDA0003934870120000027
in the control body i
Figure BDA00039348701200000213
Total internal energy of phase/J;
Figure BDA0003934870120000028
-specific enthalpy/J.kg -1
Figure BDA0003934870120000029
In the control body i
Figure BDA00039348701200000214
Non-flowing energy source term of phase/J.s -1 (ii) a The above equation ignores the gravitational potential energy and volume average kinetic energy terms of the fluid;
neglecting the derivative terms of the gas and liquid phase contributions over time, the conservation of momentum ordinary differential equation for each phase is represented by:
Figure BDA0003934870120000031
in the formula: l is j Flow channel inertial length/m, g gravity acceleration/m · s -2 ,ΔP j The head of the pump/Pa, K * -coefficient of friction of profile and wall, f 2,j -coefficient of momentum exchange between phases, L 2,j -the effective length of the inter-phase force action/m,
Figure BDA0003934870120000032
-the amount of change in velocity after passing through the flow channel,
Figure BDA0003934870120000033
-to
Figure BDA00039348701200000310
Another phase of (2), density
Figure BDA0003934870120000034
The density of the inflow control body is taken, and each item on the right side of the upper formula and the middle formula is respectively the influence of pressure difference, gravity, pump, friction and interphase force and convection.
Furthermore, the integrated serious accident analysis program adopts a semi-implicit method to disperse and solve the control equation. In each time step, the mass of the component m in the control body i at the end of the time step is equal to the sum of the mass of the control body at the last moment, the mass of the inflow and outflow control bodies in the time delta t and the mass of the source item. After adopting the linear prediction assumption, separating the unknown terms of the discrete momentum conservation equation to obtain:
Figure BDA0003934870120000035
k in the above formula 1 ,k 2 ,k 3 All of which are known at the present time,
Figure BDA0003934870120000036
indicating the pressure in the control volume i at the previous time step,
Figure BDA0003934870120000037
indicating the pressure in the control volume k at the previous time step,
Figure BDA0003934870120000038
the velocity of the gas is shown for a new time step,
Figure BDA0003934870120000039
the speed of the new time-step liquid is expressed, the equation is used as an iterative equation for solving the speed, the speed of the liquid and the gas can be solved by simultaneously connecting all control bodies under the condition of the assumed pressure and vacuole share state parameters, the mass and energy increment of the gas and the liquid in each control body can be calculated according to the new time speed, then the pressure, the temperature and the vacuole share state parameters in each control body are calculated according to the state equation, and finally the difference between the predicted value and the value after iteration is judged; and repeating the processes until the predicted value is consistent with the calculation state parameter after iteration, ending the whole iteration, and finishing the solution.
Further, the modeling of the region where the complex interaction occurs is accomplished using a core disassembly analysis program, and the method for modeling the region where the complex interaction occurs by the core disassembly analysis program is as follows: the modeling area comprises the vicinity of a steam generator break and a reactor core part, the modeling is carried out on the area by a grid division mode, because the reactor core disintegration analysis program adopts an R-Z cylindrical coordinate system for modeling, a heat transfer pipe which is mainly considered and has rupture is used as a cylindrical coordinate rotational symmetry center for modeling each part outwards, namely the heat transfer pipe part is used for cylindrical modeling, and other parts are equivalently used for annular modeling according to the principle that the radial effective area is equal.
The method for dividing the material components by the reactor core disintegration analysis program comprises the following steps: in the core disintegration analysis program, different materials are divided into density components and energy components according to the difference of density and temperature, fluid is divided into a gas velocity field and two liquid velocity fields according to the difference of flow behaviors, and each material component is required to be distributed at a corresponding position of a modeling area according to initial conditions during modeling.
The reactor core disintegration analysis program respectively establishes mass, energy and momentum conservation equations according to the divided density components, energy components and velocity fields, and the forms are as follows:
Figure BDA0003934870120000041
Figure BDA0003934870120000042
Figure BDA0003934870120000043
each term in the conservation equation is a local mean within the desired computational grid, subscripts M, q, M represent density, velocity, and energy components, respectively, t is time,
Figure BDA0003934870120000044
and Γ m 、Γ M Respectively the macroscopic density and the mass transfer rate per unit volume resulting from the density component M and the energy component M,
Figure BDA0003934870120000045
is the velocity field vector of the velocity field q, p represents the pressure,
Figure BDA0003934870120000046
is the acceleration of gravity, K qs Is a momentum exchange function between the velocity field q and the structural component, K q′q Is a momentum transfer function between the velocity fields q and q', VM q Is a virtual mass term of the gas phase, H (x) is a step function, Γ qq′ Is the quality transmission rate between q and q', e M And alpha M Is the specific energy and volume fraction of the M component, Q N 、Q M 、Q H Respectively the nuclear heat release rate, the mass and the energy exchange rate caused by the heat transmission,
Figure BDA0003934870120000047
representing the term that exists between the gas and average liquid velocity interface.
Further, the reactor core disintegration analysis program is solved based on a time decomposition method, a four-step solving method is used, and numerical diffusion is reduced by adopting high-order spatial difference. The four-step method comprises the following solving process:
solving conservation equations of mass, momentum and energy in each cell divided during modeling, obtaining transmission terms in grids on the premise of neglecting convection terms, separately processing the transmission in grids and the convection among grids, and updating nuclear heating terms of all fluid and fuel components;
the second step is used for solving the fluid convection term, a source term in a conservation equation is ignored during solving, then the conservation equation is integrated, firstly, a high-order space difference method or a first-order donor element grid difference method is used for estimating the value of each variable at the end of a time step, then the updating of quality and energy is carried out, the derivative of a state equation and the derivative of speed are calculated, and the convection term of the momentum conservation equation is processed by adopting a finite difference method;
the third step is a pressure iteration process, a multivariable Newton Simpson method is adopted to keep synchronization of the speed and the pressure of the time step end, and a matrix solving mode adopted in the pressure iteration is a direct inversion method and an incomplete LU decomposition double conjugate gradient method; the incomplete LU decomposition double-conjugate gradient method is more suitable for processing complex problems, particularly under the condition that the grid number is more than 1000, a Gaussian elimination method is adopted in one-dimensional calculation; the residual in the pressure iteration is derived from the second step estimate, and includes: the total density of the liquid, the density of the steel, the density of the lead bismuth and control particles, the total density of the mixed steam, the steam temperature and the difference value of the grid pressure calculated by using the state equation; the internal state equation iteration is saved by defining the state equation as a function of the grid pressure, the six residual errors are expanded to a first-order Taylor series, and the relation between the six residual errors and the pressure is obtained by using the state equation derivative and the speed derivative obtained in the second step;
and the fourth step is to keep the synchronization of the quality, momentum and energy of the convection item based on a semi-implicit algorithm, and simultaneously update each variable of the interface area convection item by using an explicit donor element method and a first-order donor element method.
Further, the method for realizing the data transmission between the integrated severe accident analysis program and the core disassembly analysis program by using the coupling assembly comprises the following steps: three synchronization events are used for realizing the exchange of time steps, the rollback and the data exchange; the three synchronous events are realized in a way that the first synchronization is to read the initial time step length of the integrated serious accident analysis program and the reactor core disintegration analysis program, exchange the information of the two programs, select the smaller value of the time step length of the two programs and calculate the smaller value as the common time step long time step length of the two programs; the second synchronization is that when the time steps calculated by the two programs are not consistent, the time step values determined by the two programs for the first time are halved, the number is used as a new time step of the two programs, the first step is returned to carry out the calculation again, the step is repeated until the time steps calculated by the two programs are consistent, the third synchronization is mutual data exchange, in the data exchange process, the two programs write data needing to be transmitted in a shared memory firstly, then a signal for completing the writing is sent out, when the two programs receive the signal for completing the writing sent by the other program, the required data are read from the shared memory of the other program, a signal for completing the data reading is sent out, and when the two programs receive the signal for completing the data reading of the other program, the calculation process is returned to the beginning.
Compared with the prior art, the invention has at least the following beneficial effects:
the invention provides a high-precision numerical simulation method applied to lead-based reactor SGTR accident analysis, which integrates respective advantages of an integrated serious accident analysis program and a reactor core disassembly analysis program, divides an accident occurring region into a pipeline region and a region where complex interaction occurs aiming at the lead-based reactor SGTR accident, respectively uses the integrated serious accident analysis program and the reactor core disassembly analysis program for modeling, and uses a coupling component to realize data interaction among the programs. The integrity of a modeling area is ensured to a great extent, the modeling refinement is improved to a great extent, and the accuracy of a calculation result can be ensured.
Furthermore, an integrated serious accident analysis program is used, modeling of the pipeline is realized in a mode that the control body is connected with the flow channel, and mass, energy and momentum conservation equations are respectively established for gas phase and liquid phase in the control body, so that accurate modeling calculation of the pipeline part in SGTR accident analysis is ensured.
Further, modeling the region with complex interaction by using a reactor core disintegration analysis program, establishing a conservation equation of mass, energy and momentum according to density components, energy components and a velocity field, and solving by using a four-step method, thereby ensuring the fine analysis of the accident region.
Furthermore, a coupling assembly is used, data interaction between an integrated serious accident analysis program and a reactor core disintegration analysis program is realized through three synchronous events, and finally modeling analysis on the lead-based reactor SGTR accident is completed, so that the correctness of data exchange is ensured, and a more accurate calculation result is obtained.
In conclusion, the method can give full play to the advantages of an integrated serious accident analysis program and a reactor core disintegration analysis program, ensure the integrity of the SGTR accident modeling area of the lead-based reactor, realize accurate data exchange among the programs by utilizing the coupling assembly, and obtain accurate calculation results for the problems of pipeline flow and the complex interaction of multi-phase and multi-fluid related to accidents. The accuracy of the calculation result can be ensured.
The technical solution of the present invention is further described in detail by the accompanying drawings and embodiments.
Drawings
FIGS. 1a and 1b are longitudinal and transverse cross-sectional views, respectively, of an ALFERD reactor;
FIG. 2 is a schematic diagram of an ALFERD reactor SGTR fault modeling;
FIG. 3 is a schematic view of the inside of the control body in an integrated severe accident analysis program;
FIG. 4 is a schematic diagram of the flow paths in an integrated severe accident analysis program;
FIG. 5 is a flow of solving the integrated severe accident analysis program;
FIG. 6 is a four-step calculation process;
FIG. 7 is a diagram of MpCCI parallel coupling scheme;
FIG. 8 is a schematic diagram of coupled signal synchronization;
FIG. 9 is a block diagram of a module coupling calculation flow;
FIG. 10 is a graph illustrating steam generator pressure changes during a SGTR event in an ALFERD reactor with seven tubes burst.
Detailed Description
The invention relates to a high-precision numerical simulation method applying lead-based reactor SGTR accident analysis, which takes modeling calculation analysis of an ALFERD reactor SGTR accident as an example and explains the analysis steps of the lead-based reactor SGTR accident in the invention:
ALFRED is a 300MWth (120 MWe) integrated pool type lead-cooled fast reactor. All major components are contained within the reactor vessel and the flow of liquid metallic lead within the reactor is shown in figures 1a and 1b, with the flow path of the lead being indicated by the arrows. When an SGTR accident occurs, water with the temperature lower than the saturation of 7K interacts with low-pressure high-temperature molten lead at about 180bar, the calculation simulates the extreme working condition that 7 heat transfer pipes are cracked at the same time, and the change conditions of parameters such as pressure in the steam generator after the SGTR accident occurs are researched.
First, modeling the pipeline part by using an integrated serious accident analysis program, as shown in the left half part of fig. 2, modeling a single heat transfer PIPE to control the bodies PIPE121 to PIPE128 as a flow pipeline of liquid lead, setting PIPE129 and PIPE120 as time-dependent control bodies as boundary conditions to control the flow of the liquid lead, and setting the control bodies PIPE101 to PIPE110 as a pipeline part of water, where the liquid water absorbs heat and gradually vaporizes, adding a thermal structure between the two pipelines for simulating heat exchange between the liquid lead and the water.
The schematic diagram of the control body is shown in fig. 3, the schematic diagram of the flow channel is shown in fig. 4, two control bodies are connected by a flow channel, the residence time of the fluid in the flow channel is not considered, and the mass or the energy in all the flow channels is not considered, namely the mass and the energy of the fluid are stored in the control bodies. The controller distinguishes 'liquid' and 'gas' parts in the body, representing liquid phase and gas phase, and adopts two fluid models to respectively use mass, momentum and energy conservation laws for the gas phase and the liquid phase, wherein a mass conservation constant differential equation of a component m in a controller i can be expressed by the following formula:
Figure BDA0003934870120000081
in the formula: m i,m The total mass/kg of component m in control body i, the angle m representing the composition of the water, the droplets or the various gas items, j-the defluidizing control body, σ ij -the direction of the flow of the gas,
Figure BDA0003934870120000082
to control body j
Figure BDA00039348701200000814
The fraction of the phase(s) is,
Figure BDA00039348701200000815
is a gas phase or a liquid phase in the flow channel, t-time/s,
Figure BDA0003934870120000083
the density of component m in the defluidizing control volume j and the inflowing control volume d/kg m -3 And d represents the incoming flow control volume,
Figure BDA0003934870120000084
-in the defluidizing control body j
Figure BDA00039348701200000816
Flow velocity of phase/m.s, A j Flow area/m of the defluidizing control body j 2 ,F j The opening of the flow area of the defluidizing control body j,
Figure BDA0003934870120000085
controlling the non-flowing mass source term/kg-s of component m in body i -1 . The energy conservation ordinary differential equation for each phase can be expressed by the following equation:
Figure BDA0003934870120000086
in the formula:
Figure BDA0003934870120000087
in the control body i
Figure BDA00039348701200000817
Total internal energy of phase/J;
Figure BDA0003934870120000088
specific enthalpy/J.kg -1
Figure BDA0003934870120000089
In the control body i
Figure BDA00039348701200000818
Non-flowing energy source term/J.s of phase -1 (ii) a The above equation ignores the gravitational potential energy and volume average kinetic energy terms of the fluid;
neglecting the derivative terms of the gas and liquid phase contributions over time, the conservation of momentum ordinary differential equation for each phase is represented by:
Figure BDA00039348701200000810
in the formula: l is j Flow channel inertial length/m, g gravitational acceleration/m.s -2 ,ΔP j The head of the pump/Pa, K * -coefficient of friction of profile and wall, f 2,j -coefficient of momentum exchange between phases, L 2,j -the effective length of the inter-phase force action/m,
Figure BDA00039348701200000811
-the amount of change in velocity after passing through the flow channel,
Figure BDA00039348701200000812
-to
Figure BDA00039348701200000819
Another phase of (2), density
Figure BDA00039348701200000813
The density of the inflow control body is taken, and each item on the right side of the upper formula and the middle formula is respectively the influence of pressure difference, gravity, pump, friction and interphase force and convection.
And the integrated serious accident analysis program adopts a semi-implicit method to disperse and solve the control equation. In each time step, the mass of the component m in the control body i at the end of the time step is equal to the sum of the mass of the control body at the last moment, the mass of the inflow and outflow control body in the time delta t and the mass of the source item. After adopting the linear prediction assumption, separating the unknown terms of the discrete momentum conservation equation to obtain:
Figure BDA0003934870120000091
k in the above formula 1 ,k 2 ,k 3 All of which are known at the present time,
Figure BDA0003934870120000092
indicating the pressure in the control volume i at the previous time step,
Figure BDA0003934870120000093
indicating the pressure in the control volume k at the previous time step,
Figure BDA0003934870120000094
the velocity of the gas is shown for a new time step,
Figure BDA0003934870120000095
indicating the velocity of the liquid at the new time step. The equation is used as an iterative equation for solving the speed, the speed can be solved by simultaneously combining all control bodies under the assumed state parameters of pressure, vacuole share and the like, the mass and energy increment can be calculated according to the new moment speed, then the state parameters of pressure, temperature, vacuole share and the like are calculated according to the state equation, and finally the difference between the predicted values of pressure and the like and the values after iteration is judged. The above processes are repeated until the predicted value is consistent with the post-iteration calculation state parameter, the whole iteration is finished, the solution is completed, and the solution flow of the integrated serious accident analysis program is shown in fig. 5.
Then, a core disassembly analysis program is used to model the peripheral region and the core part where the breach has occurred in the steam generator, as shown in the right half part of fig. 2, because the multiphase multicomponent modeling is based on an R-Z cylindrical coordinate system, the principle of the equivalent radial effective area needs to be followed during modeling, but in this problem, the broken heat transfer tube is an object to be considered in an important way, and is placed at the center to restore the cylindrical geometry, and other components are modeled as a ring shape according to the principle of the equivalent area and are distributed around the broken heat transfer tube. The reactor was divided into 28 radial grids and 50 axial grids, and in fig. 2, the area 1 represents lead, the area 2 water, the heat transfer tube and grid structure is represented by the area 3, the area 4 argon, the area 5 non-calculation area, and the position where the rupture occurred is located at the 11 axial grids. Finally, in order to enable the two programs to smoothly realize data interaction, a control body PIPE100 is added on the side of the integrated serious accident analysis program and is used for receiving boundary information of the reactor core disintegration analysis program and transmitting the boundary information to the integrated serious accident analysis program.
The core disintegration analysis program divides each component into a density component and an energy component according to the difference of density and temperature, and divides the fluid into a gas velocity field and two liquid velocity fields according to the difference of flow behaviors. Modeling requires distribution of material components in appropriate locations based on initial conditions. The program establishes conservation equations of mass, energy and momentum according to the divided density components, energy components and velocity fields respectively, and the conservation equations are in the following forms:
Figure BDA0003934870120000101
Figure BDA0003934870120000102
Figure BDA0003934870120000103
each term in the conservation equation is a local mean value in the required computational grid, subscripts M, q, M represent density, velocity and energy components, respectively, t is time,
Figure BDA0003934870120000104
and Γ m 、Γ M The macroscopic density and the mass transfer rate per unit volume resulting from the density component M and the energy component M respectively,
Figure BDA0003934870120000105
is the velocity field vector of the velocity field q, p represents the pressure,
Figure BDA0003934870120000106
is the acceleration of gravity, K qs Is a momentum exchange function between the velocity field q and the structural component, K q′q Is a momentum transfer function between the velocity fields q and q', VM q Is a virtual mass term of the gas phase, H (x) is a step function, Γ qq′ Is the quality transmission rate between q and q', e M And alpha M Is the specific internal energy and volume fraction of the M component, Q N 、Q M 、Q H Respectively the nuclear heat release rate, the energy exchange rate caused by mass and heat transmission,
Figure BDA0003934870120000107
representing the term that exists between the gas and average liquid velocity interfaces. The specific composition tables are shown in tables 1 to 3.
Table 1 structural Components
Figure BDA0003934870120000108
Figure BDA0003934870120000111
TABLE 2 liquid Components
Figure BDA0003934870120000112
TABLE 2 liquid Components
Figure BDA0003934870120000113
TABLE 3 gas composition
Figure BDA0003934870120000114
The reactor core disintegration analysis program is solved based on a time decomposition method, and numerical diffusion is reduced by using a four-step solving method and adopting high-order spatial difference. The four-step solution flow is shown in fig. 6. Solving a conservation equation of mass, momentum and energy in the cells, obtaining a transmission item in the grids on the premise of neglecting a convection item, separately processing the transmission in the grids and the convection between the grids, and updating nuclear heating items of all fluid and fuel components.
The second step is used for solving the fluid convection term, the source term in the conservation equation is ignored during solving, and then the conservation equation is integrated. Firstly, estimating the value of each variable at the end of a time step by using a high-order space difference method or a first-order donor element grid difference method, then updating the quality and the energy, calculating the derivative of a state equation and the derivative of the velocity, and processing the convection term of a momentum conservation equation by using a finite difference method.
The third step is a pressure iteration process, a multivariable Newton Simpson method is adopted to keep synchronization of the speed and the pressure of the time step end, and a matrix solving mode adopted in the pressure iteration is a direct inversion method and an incomplete LU decomposition double conjugate gradient method; the incomplete LU decomposition double conjugate gradient method is more suitable for processing complex problems, particularly under the condition that the grid number is more than 1000, a Gaussian elimination method is adopted in one-dimensional calculation; the residual in the pressure iteration is derived from the estimate of the second step, and includes: the total density of the liquid, the density of the steel, the density of the lead bismuth and control particles, the total density of the mixed steam, the steam temperature and the difference value of the grid pressure calculated by using the state equation; the internal state equation iteration is saved by defining the state equation as a function of the grid pressure, the six residual errors are expanded to a first-order Taylor series, and the relation between the six residual errors and the pressure is obtained by using the state equation derivative and the speed derivative obtained in the second step;
and the fourth step is to keep the synchronization of the quality, momentum and energy of the convection item based on a semi-implicit algorithm, and simultaneously update each variable of the interface area convection item by using an explicit donor element method and a first-order donor element method. Table 4 gives a summary of the variable updates for the 4-step method, where ρ is (1) 、ρ (2) 、ρ (3) And ρ (4) The updated density rho of the 1 st, 2 nd, 3 rd and 4 th steps (n) Density at the beginning of the time step, p (n+1) Density at the end of the time step, v (1) 、v (2) 、v (3) And v (4) Updated speed, v, of steps 1, 2, 3 and 4, respectively (n) Speed at the beginning of the time step, v (n+1) The velocity at the end of the time step, t is the time, e is the specific internal energy,
TABLE 4 variable update procedure of four-step method
Figure BDA0003934870120000121
Figure BDA0003934870120000131
TABLE 4 continuation four-step method variable updating process
Figure BDA0003934870120000132
In order to ensure that the calculation of the integrated serious accident analysis program and the calculation of the reactor core disintegration analysis program can realize synchronous parallelism, the time step propulsion and the data exchange of the integrated serious accident analysis program and the reactor core disintegration analysis program must be ensured to be strictly unified. In the coupling calculation, the two programs need to make the respective time step and the cycle step number the same, and the time point for data exchange is reasonable, and in order to realize the points, the invention refers to the program scheduling sequence in the mcpci, and as shown in fig. 7, three synchronous events are used to realize the time step exchange, the backoff time step exchange and the data exchange.
The integrated severe accident analysis program and the core disassembly analysis program are substantially the same in calculation flow, and include 4 steps of determining an initial value of a system parameter, predicting or determining a time step size, calculating a thermodynamic power, checking whether convergence and advancing the time step, however, the integrated severe accident analysis program and the core disassembly analysis program have certain differences in the sequence of implementing the steps, and therefore, it is necessary to determine whether the time steps are unified in the coupling program. The specific method comprises the following steps: first, the initial time step of the integrated severe accident analysis program and the core disassembly analysis program is read and the information of the two is exchanged, so that the first synchronization between the two is completed. Then, selecting the smaller value of the time step of the two programs, using the value as the time step of the two programs, and calculating the time step as a unified time step, and then, because the solving methods in the two programs are different, the time steps caused by convergence judgment after the calculation of the two programs are ended are probably not consistent, and then, information needs to be returned, namely, the second synchronization, and the specific operation is as follows: the time step values determined by the two programs for the first time are halved, the values are used as new time step lengths of the two programs, the calculation is returned to the first step, when the two programs are calculated and the converged time step lengths are judged to be consistent, the second step can be considered to be finished, then the third synchronization is carried out, namely, the calculated state parameters of the two programs are updated, and then the mutual data exchange is carried out. A schematic diagram of signal synchronization between programs is shown in fig. 8. This completes three synchronization events in a time-stepped computation, the flow chart of which is shown in fig. 9.
It has been mentioned in the foregoing that, the technical solution of the present invention is applied to model and calculate the SGTR accident of the alfard-reflective stack, and a pressure change curve of the breach position is given, as shown in fig. 10, it can be seen from the graph that the calculated result is consistent with the pressure change trend of the calculated result of the reference program, but since the whole steam generator can be modeled when the technical solution of the present invention is used for modeling, and the reference program only partially models the steam generator, the pressure calculated result obtained by applying the present invention has a slower decrease trend. The technical scheme of the invention can be used for completely modeling the whole problem area, can fully exert the advantages of two programs and avoids the defects of the two programs.
In conclusion, the method established by the invention integrates the respective advantages of the integrated serious accident analysis program and the reactor core disintegration analysis program, aiming at the lead-based reactor SGTR accident, the area where the accident occurs is divided into a pipeline area and an area where complex interaction occurs, meanwhile, the integrated serious accident analysis program which is good at calculating the flow problem in the pipeline is used for modeling calculation of the pipeline area, the reactor core disintegration analysis program which is good at calculating the multi-phase and multi-fluid complex interaction problem is used for modeling calculation of the area where complex interaction occurs, the modeling integrity of the area where the accident occurs is ensured to a great extent, the modeling is improved to a great extent, meanwhile, the coupling assembly is used for realizing accurate data transmission between the two programs, and the accuracy of the calculation result can be ensured.

Claims (9)

1. A high-precision numerical simulation method applied to SGTR accident analysis of a lead-based reactor is characterized in that a region of the lead-based reactor where the SGTR accident occurs is divided into a pipeline region and a region where complex interaction occurs, namely, the vicinity of a steam generator break and a reactor core part, an integrated serious accident analysis program is used for modeling the pipeline region, the integrated serious accident analysis program can establish a control equation for a modeled controller and solve the control equation discretely, the reactor core disassembly analysis program is used for modeling the region where the complex interaction occurs, the reactor core disassembly analysis program divides components of materials, meanwhile, a conservation equation is established and solved, a coupling assembly is used for realizing data transmission between the integrated serious accident analysis program and the reactor core disassembly analysis program, and modeling and analysis of the SGTR accident are completed.
2. The high-precision numerical simulation method applied to SGTR accident analysis of the lead-based stacks as claimed in claim 1, wherein the method for modeling the pipeline area by using an integrated serious accident analysis program is as follows: when modeling calculation is carried out on the pipeline through an integrated serious accident analysis program, the pipeline is simulated through a mode that a control body and a flow channel are connected mutually.
3. The high-precision numerical simulation method applied to SGTR accident analysis of the lead-based stacks as claimed in claim 1, wherein the method for the integrated serious accident analysis program to establish the control equation for the modeled control body is as follows: the control body distinguishes liquid and gas parts, respectively represents a liquid phase and a gas phase, two fluid models are adopted, and a control equation is established for the gas phase and the liquid phase by using mass, momentum and energy conservation laws, wherein a mass conservation constant differential equation of a component m in a control body i is expressed by the following formula:
Figure FDA0003934870110000011
in the formula: m i,m The total mass/kg of component m in control body i, the angle m representing the composition of the water, the droplets or the various gas items, j-the defluidizing control body, σ ij -the direction of the flow of the gas,
Figure FDA0003934870110000012
to control body j
Figure FDA0003934870110000013
The fraction of the phase(s) is,
Figure FDA0003934870110000014
is a gas phase or a liquid phase in the flow channel, t-time/s,
Figure FDA0003934870110000015
the density of component m in the defluidizing control volume j and the inflowing control volume d/kg m -3 And d represents the incoming flow control volume,
Figure FDA0003934870110000016
-in the defluidizing control body j
Figure FDA0003934870110000017
The flow velocity of the phase/m.s,A j flow area/m of the defluidizing control body j 2 ,F j -opening of the flow area of the defluidizing control body j,
Figure FDA0003934870110000018
controlling the non-flowing mass source term/kg-s of component m in body i -1
The energy conservation ordinary differential equation for each phase is expressed by the following equation:
Figure FDA0003934870110000021
in the formula:
Figure FDA0003934870110000022
in the control body i
Figure FDA0003934870110000023
Total internal energy of phase/J;
Figure FDA0003934870110000024
specific enthalpy/J.kg -1
Figure FDA0003934870110000025
In the control body i
Figure FDA0003934870110000026
Non-flowing energy source term/J.s of phase -1 (ii) a The above equation ignores the gravitational potential energy and volume average kinetic energy terms of the fluid;
neglecting the derivative terms of the gas and liquid phase contributions over time, the conservation of momentum ordinary differential equation for each phase is represented by:
Figure FDA0003934870110000027
in the formula: l is j -the length of inertia of the flow channel/m,g-gravitational acceleration/m.s -2 ,ΔP j head/Pa, K of the pump * -coefficient of friction of form and wall, f 2,j -coefficient of momentum exchange between phases, L 2,j -the effective length of the inter-phase force action/m,
Figure FDA0003934870110000028
-the amount of change in velocity after passing through the flow channel,
Figure FDA0003934870110000029
-to
Figure FDA00039348701100000210
Another phase of (2), density
Figure FDA00039348701100000211
The density of the inflow control body is taken, and each item on the right side of the upper formula and the middle formula is respectively the influence of pressure difference, gravity, pump, friction and interphase force and convection.
4. The high-precision numerical simulation method applied to SGTR accident analysis of the lead-based stacks as claimed in claim 1, wherein the method for discretely solving the control equation by an integrated serious accident analysis program is as follows: the integrated serious accident analysis program adopts a semi-implicit method to disperse and solve a control equation, in each time step, the mass of a component m in a control body i at the end of the time step is equal to the mass of the control body at the last moment, and the mass of an inflow and outflow control body and the mass of a source item in delta t time, and after a linear prediction hypothesis is adopted, an unknown item of a discrete momentum conservation equation is separated to obtain:
Figure FDA00039348701100000212
k in the above formula 1 ,k 2 ,k 3 All of which are known at the present time,
Figure FDA00039348701100000213
indicating the pressure in the control volume i at the previous time step,
Figure FDA00039348701100000214
indicating the pressure in the control volume k at the previous time step,
Figure FDA00039348701100000215
the velocity of the gas is shown for a new time step,
Figure FDA00039348701100000216
the speed of the new time-step liquid is expressed, the equation is used as an iterative equation for solving the speed, the speed of the liquid and the gas can be solved by simultaneously connecting all control bodies under the condition of the assumed state parameters of pressure and vacuole share, the mass and energy increment of the gas and the liquid in each control body can be calculated according to the new time speed, then the state parameters of the pressure, the temperature and the vacuole share in each control body are calculated according to the state equation, and finally the difference between the predicted value and the value after iteration is judged; and repeating the processes until the predicted value is consistent with the post-iteration calculation state parameter, finishing the whole iteration, and finishing the solution.
5. The high-precision numerical simulation method applied to SGTR accident analysis of the lead-based reactor as set forth in claim 1, wherein the core disassembly analysis program models the regions where complex interactions occur by: the modeling region comprises the vicinity of a steam generator break and a reactor core part, the region is modeled by a grid division mode, because a reactor core disassembly analysis program adopts an R-Z cylindrical coordinate system for modeling, a heat transfer pipe which is seriously considered and is broken is used as a cylindrical coordinate rotational symmetry center for outward modeling of each part, namely the heat transfer pipe part is used for cylindrical modeling, and other parts are equivalently used for annular modeling according to the principle that the radial effective areas are equal.
6. The high-precision numerical simulation method applied to SGTR accident analysis of the lead-based reactor as claimed in claim 1, wherein the method for dividing the material components by the core disintegration analysis program is as follows: in the reactor core disintegration analysis program, different materials are divided into density components and energy components according to different densities and temperatures, fluid is divided into a gas velocity field and two liquid velocity fields according to different flow behaviors, and each material component needs to be distributed at a corresponding position of a modeling area according to initial conditions during modeling.
7. The high-precision numerical simulation method applied to SGTR accident analysis of the lead-based reactor as claimed in claim 1, wherein the method for establishing the conservation equation by the reactor core disassembly analysis program is as follows: the reactor core disintegration analysis program respectively establishes mass, energy and momentum conservation equations according to the divided density components, energy components and velocity fields, and the forms are as follows:
Figure FDA0003934870110000031
Figure FDA0003934870110000032
Figure FDA0003934870110000033
each term in the conservation equation is a local mean within the desired computational grid, subscripts M, q, M represent density, velocity, and energy components, respectively, t is time,
Figure FDA0003934870110000041
and Γ m 、Γ M Respectively the macroscopic density and the mass transfer rate per unit volume resulting from the density component M and the energy component M,
Figure FDA0003934870110000042
is speedThe velocity field vector of the degree field q, p stands for pressure,
Figure FDA0003934870110000043
is the acceleration of gravity, K qs Is a momentum exchange function between the velocity field q and the structural component, K q′q Is a momentum transfer function between the velocity fields q and q', VM q Is a virtual mass term of the gas phase, H (x) is a step function, Γ qq′ Is the quality transmission rate between q and q', e M And alpha M Is the specific internal energy and volume fraction of the M component, Q N 、Q M 、Q H Respectively the nuclear heat release rate, the mass and the energy exchange rate caused by the heat transmission,
Figure FDA0003934870110000044
representing the term that exists between the gas and average liquid velocity interfaces.
8. The high-precision numerical simulation method applied to SGTR accident analysis of the lead-based reactor as claimed in claim 1, wherein the method for solving the conservation equation by the core disassembly analysis program is as follows: the solution of the reactor core disintegration analysis program to the conservation equation is based on a time decomposition method, and a four-step method is used, wherein the four-step method comprises the following solving flow:
solving conservation equations of mass, momentum and energy in each cell divided during modeling, obtaining a transmission item in a grid on the premise of ignoring a convection item, separately processing the transmission in the grid and the convection among the grids, and updating nuclear heating items of all fluid and fuel components;
the second step is used for solving the fluid convection term, a source term in a conservation equation is ignored during solving, then the conservation equation is integrated, firstly, a high-order space difference method or a first-order donor element grid difference method is used for estimating the value of each variable at the end of a time step, then the updating of quality and energy is carried out, the derivative of a state equation and the derivative of speed are calculated, and the convection term of the momentum conservation equation is processed by adopting a finite difference method;
the third step is a pressure iteration process, a multivariable Newton Simpson method is adopted to keep synchronization of the speed and the pressure of the time step end, and a matrix solving mode adopted in the pressure iteration is a direct inversion method and an incomplete LU decomposition double conjugate gradient method; the incomplete LU decomposition double-conjugate gradient method is more suitable for processing complex problems, particularly under the condition that the grid number is more than 1000, a Gaussian elimination method is adopted in one-dimensional calculation; the residual in the pressure iteration is derived from the second step estimate, and includes: the total density of the liquid, the density of the steel, the density of the lead bismuth and control particles, the total density of the mixed steam, the steam temperature and the difference value of the grid pressure calculated by using the state equation; the internal state equation iteration is saved by defining the state equation as a function of the grid pressure, the six residual errors are expanded to a first-order Taylor series, and the relation between the six residual errors and the pressure is obtained by using the state equation derivative and the speed derivative obtained in the second step;
and the fourth step is to keep the synchronization of the quality, momentum and energy of the convection item based on a semi-implicit algorithm, and simultaneously update each variable of the interface area convection item by using an explicit donor element method and a first-order donor element method.
9. The high-precision numerical simulation method applied to the SGTR fault analysis of the lead-based reactor as claimed in claim 1, wherein the method for realizing the data transmission between the integrated severe fault analysis program and the core disassembly analysis program by using the coupling assembly comprises the following steps: three synchronization events are used for realizing the exchange of time steps, the rollback and the data exchange; the three synchronization events are realized in such a way that the first synchronization is to read the initial time step length of the integrated severe accident analysis program and the reactor core disassembly analysis program, exchange the information of the two programs, select the smaller value of the time step lengths of the two programs and use the value as the common time step long time step length of the two programs for calculation; the second synchronization is that when the time steps calculated by the two programs are not consistent, the time step values determined by the two programs for the first time are halved, the number is used as a new time step of the two programs, the first step is returned to carry out the calculation again, the step is repeated until the time steps calculated by the two programs are consistent, the third synchronization is mutual data exchange, in the data exchange process, the two programs write data needing to be transmitted in a shared memory firstly, then a signal for completing the writing is sent out, when the two programs receive the signal for completing the writing sent by the other program, the required data are read from the shared memory of the other program, a signal for completing the data reading is sent out, and when the two programs receive the signal for completing the data reading of the other program, the calculation process is returned to the beginning.
CN202211400808.XA 2022-11-09 2022-11-09 High-precision numerical simulation method applied to lead-based stack SGTR accident analysis Pending CN115828666A (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202211400808.XA CN115828666A (en) 2022-11-09 2022-11-09 High-precision numerical simulation method applied to lead-based stack SGTR accident analysis

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202211400808.XA CN115828666A (en) 2022-11-09 2022-11-09 High-precision numerical simulation method applied to lead-based stack SGTR accident analysis

Publications (1)

Publication Number Publication Date
CN115828666A true CN115828666A (en) 2023-03-21

Family

ID=85527456

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202211400808.XA Pending CN115828666A (en) 2022-11-09 2022-11-09 High-precision numerical simulation method applied to lead-based stack SGTR accident analysis

Country Status (1)

Country Link
CN (1) CN115828666A (en)

Similar Documents

Publication Publication Date Title
CN105247623B (en) Methods for simulating the flow of fluid in nuclear reactor and for calculating the mechanical deformation of fuel assemblies
Poitevin et al. Tritium breeder blankets design and technologies in Europe: Development status of ITER Test Blanket Modules, test & qualification strategy and roadmap towards DEMO
Angelucci et al. STH-CFD codes coupled calculations applied to HLM loop and pool systems
Toti et al. Coupled system thermal-hydraulic/CFD analysis of a protected loss of flow transient in the MYRRHA reactor
Smith Assessment of CFD codes used in nuclear reactor safety simulations
Toti et al. Extension and application on a pool-type test facility of a system thermal-hydraulic/CFD coupling method for transient flow analyses
CN111125972B (en) Hydraulic load analysis method for water loss accident of break of nuclear power plant
CN115828666A (en) High-precision numerical simulation method applied to lead-based stack SGTR accident analysis
CN105247622B (en) Methods for simulating the flow of a fluid in nuclear reactor and for calculating the mechanical deformation of fuel assemblies
Tentner et al. Modeling of two-phase flow in a BWR fuel assembly with a highly-scalable CFD code
Toti et al. Development and preliminary validation of a STH-CFD coupling method for multiscale thermal-hydraulic simulations of the MYRRHA reactor
Vivaldi A 3D model to solve U-tube steam generator secondary side thermal-hydraulics with coupled primary-to-secondary side heat transfer
Xiaoying et al. Fluid-structure interaction analysis using finite element methods for IRWST of reactor building
Kamali Development of a High-Fidelity Aero-Thermo-Elastic Analysis and Design Capability
Zwijsen et al. Fast reactor multi-scale and multi-physics modelling at NRG
Shin et al. Comparative study of constitutive relations implemented in RELAP5 and TRACE–Part I: Methodology & wall friction
Beccantini et al. Simulation of hydrogen release and combustion in large scale geometries: models and methods
Rockholm et al. Compact Steam Generator Numerical Analysis and Mechanics Design
Bai et al. Long-Term Simulation of Sodium Dynamics During a Large Leakage Sodium-Water Reaction
Soler-Martinez Semi-implicit thermal-hydraulic coupling of advanced subchannel and system codes for pressurized water reactor transient applications
Akin Numerical Simulation of Corium Jet Fragmentation
Corzo et al. Thermal-hydraulic 0D/3D coupling in OpenFOAM: Validation and application in nuclear installations
Archambeau et al. Qualification of code-Saturne for thermal-hydraulics single phase nuclear applications
Gu et al. Code development of Single-phase 2-D pressure wave propagation
Van Antwerpen Modelling a pebble bed high temperature gas-cooled reactor using a System-CFD approach

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination