A method for calculating two-dimensional void fraction distribution of fission gas in a fissile aqueous solution system

Abstract

The application discloses a kind of two-dimensional radiolysis gas void fraction distribution calculation methods of fissile aqueous solution system, and relates to nuclear and radiation safety technical field.The method carries out axial and radial partition to the internal cavity region of cylindrical container, and the radiolysis bubble in solution is divided into two groups of bubbles according to its behavior after nucleation, and the nucleation, growth and migration process of radiolysis bubble are considered.Compared with prior art, the method of the application avoids two-phase flow simulation while considering more specific radiolysis gas behavior process, can quickly and accurately analyze the solution system nuclear criticality accident safety, and provides parameter basis for solution system related nuclear criticality safety accident radiation protection design and reactor design.

Description

Technical Field

[0001] This invention belongs to the field of nuclear and radiation safety technology, and specifically relates to a method for calculating the two-dimensional radiolytic gas cavitation distribution in a fissile aqueous solution system. Background Technology

[0002] Homogeneous uranium salt aqueous solutions are commonly found in homogeneous aqueous reactors and spent fuel reprocessing. When a fissile aqueous system reaches criticality, fission products collide with substances in the solution, resulting in radiodecomposition and the formation of radiodecomposition products. Compared to other fission products, radiodecomposition products from fission fragments account for approximately 96%, and the primary products formed by their collisions with water molecules are mainly hydrogen and hydrogen peroxide. Furthermore, the energy rapidly released along the trajectory of fission fragments causes the generation of water vapor; however, if the solution has not reached saturation temperature, the water vapor dissolves rapidly. Therefore, the gas formed along the trajectory of fission fragments is mostly hydrogen.

[0003] Before the hydrogen saturation concentration is reached, the small bubbles formed by the rapid generation and accumulation of hydrogen will dissolve until the saturation concentration is reached. These small bubbles generated along the trajectory of fission fragments will no longer dissolve, and the radiation cracking gas dissolved in the solution will diffuse into the bubbles, causing the bubbles to grow rapidly.

[0004] Increased temperature and rapid bubble growth create inertial pressure, posing a safety threat to the container walls. Simultaneously, these stable bubbles, buoyant and escaping the solution, cause power oscillations. Therefore, modeling and analyzing the behavior of radiolytic bubbles is crucial in the transient analysis of homogeneous uranium salt aqueous solutions. Summary of the Invention

[0005] This invention provides a method for calculating the two-dimensional cavitation distribution of radiolytic gas in a fissile aqueous solution system based on a simplified model of the behavior of radiolytic bubbles in solution.

[0006] This invention is achieved using the following technical solution: A method for calculating the cavitation fraction distribution of a two-dimensional radiolytic gas in a fissile aqueous solution system includes the following steps: (1) Determine the geometric parameters, material properties, mesh generation parameters, and hydrogen production per unit energy of the solution system. The unit is: mol / J turbulence coefficient The unit is m. 2 / s, and bubble-related parameters; divide the system into grids based on the input parameters, and determine the coordinates and dimensions of each grid; determine the system power distribution. The unit is W, and the temperature distribution is... The unit is: ℃ ,in and To control the physical quantities within the body, all of them change with time and space; (2) The radiolysis bubbles in the solution are divided into three types: nucleation bubbles, dissolution bubbles, and growth bubbles. Nucleation bubbles are bubbles generated along the trajectory of fission fragments. Based on whether they grow in the solution, they are divided into two groups of bubbles. Nucleation bubbles and dissolution bubbles are the first group of bubbles, and growth bubbles are the second group of bubbles. The following information about the bubbles needs to be stored using a matrix to reflect the two-dimensional spatial distribution characteristics of the variables: a. The radii of the two sets of bubbles, , The unit is: m ; b. The flow rates of the two sets of bubbles, , The unit is: m / s ; c. The number density of the two groups of bubbles. , The unit is: individual. / m 3 ; d. The mass of gas per unit volume of the two sets of bubbles. , The unit is: mol / m 3 ; e. Control the volume of the two sets of air bubbles inside the body. , The unit is: m 3 ; f. The gas density of the two sets of bubbles, , The unit is: kg / m 3 ; g. Radius of the nucleating bubble The unit is m equilibrium concentration of nucleated bubbles The unit is mol / m 3 First group of bubble dissolution time The unit is seconds (s). h. Equilibrium concentration of growth bubbles The unit is mol / m 3 ; A, the interfacial area density between the growth bubbles and the solution (in m³); mass transfer coefficient between the growth bubbles and the solution (in m³). The unit is m / s; (3) Use the superscript 0 and last to mark the initial value and the value at the previous moment of the relevant parameters; start the calculation from moment 0, the liquid level is The axial position number of the grid where the liquid surface is located Radiolysis gas concentration and , , , All values ​​are set to 0, solution pressure Taken as atmospheric pressure; solution velocity , Set to 0; based on the inner diameter of the cylindrical container. , outer diameter and liquid level and initial density of solution Calculate the total mass of the solution at the initial time. ;Assign the value from the previous time step to the value from the initial time step; (4) According to , , , , Based on the ideal gas law, iterative calculation , ;according to , Calculate the volume of the two sets of radiolysis bubbles , ;according to , and control volume Calculate the void fraction ;make hour, , to be used to calculate the grid position number where the liquid level is located in step (8); (5) Based on solution pressure The solution velocity at the previous moment , Based on the momentum conservation equation and corresponding boundary conditions, the solution flow rate is calculated. , ;according to , and Based on the continuity equation, the density change caused by solution motion is calculated. ;according to and The difference and the solution density at the previous moment Calculate the density change caused by the compressibility of radiolytic bubbles. ;according to and Calculate the change in solution density and according to Calculate the solution density This step only calculates the number. The physical quantities within the following grid region, for The grid area ; (6) According to , , , , Calculation number The isothermal compressibility of the gas-liquid mixture within the following grid region and isobaric expansion coefficient ; (7) According to and The difference , , , Based on the equation of state for gas-liquid mixtures, the number was calculated. Pressure changes within the following grid area Other regions always assume ;based on and the pressure of the previous moment Calculate solution pressure ; (8) Obtained according to step (4) and the results obtained in step (5) Summing the solution masses within each control cell along the axial direction from bottom to top, the total solution mass below the J-th grid layer is calculated. until Obtain the current liquid level position number. ;Calculate the liquid level based on the conservation of the total mass of the solution. ;according to , , , and The void fraction distribution is re-given. Update the cavitation fraction and solution pressure to make , ; (9) Repeat steps (4) to (8) until the maximum relative error of the control body solution pressure before and after the iteration is reached. It meets the convergence requirements and enables coupled calculation of solution flow rate, solution pressure, cavitation fraction, liquid level, and bubble radius; (10) According to and According to the physical quantity at the previous moment , , , and To obtain the redistribution of physical quantities after the liquid level change: , , , and ;according to , and , give and Bubble radius between , , Update the liquid level location number, and make ; (11) Calculate the number Within the following grid area, physical quantities other than c and d in step (2); based on the bubble radius , Solution pressure Herry constant H and surface tension Calculate the equilibrium concentration separately , According to the radius , and pressure Calculate the gas density of the two groups of bubbles respectively. and And based on these parameters and solution density Calculate the flow velocity sequentially and ;according to , , and solution flow rate , Calculate the interphase interface mass transfer coefficient ;according to , Calculate the interfacial area density A; (12) Parameters obtained from the above steps , , , , , A, and input parameters , P, T, construct , , , and The differential equation; according to Determine the liquid level location and provide boundary conditions; calculate according to step (10). , , , and Calculate the current time number The physical quantities within the following grid region, , , , and For numbering In the region and above, all physical quantities are set to 0; (13) By repeating steps (4) to (12), the linear differential equation system in step (12) is repeatedly constructed and solved to achieve coupled computation between the equation systems until And two iterations , , , and Maximum relative error The convergence requirement is met; (14) Use the calculated value at the current moment as the value at the previous moment to update the physical quantity: , , , , , , , , , and Output and Used for reactive feedback calculations; (15) Starting from time zero, advance the time step and repeatedly input the power at different times. , And continue to calculate according to steps (4) to (14) until the time required for transient calculation is reached.

[0007] A further improvement of the present invention is that, in step (1), the geometric parameters of the solution system are a ring-shaped solution region. , outer diameter and liquid level Height of the air cavity region inside the cylindrical container The physical property of the solution material is density. Surface tension Herry constant H, uranium concentration The mesh generation parameters are the axial and radial mesh counts in the solution region and the axial mesh count in the gas cavity region; the bubble-related parameters are the nucleation radius. The dissolution time of the first group of bubbles Radiation gas molecule diffusion rate m 2 / s.

[0008] A further improvement of the present invention is that, in step (3), the total mass of the solution at the initial moment is calculated according to the following relationship:

[0009] In step (4), when the number density of bubbles in each group At that time, the radius of each group of bubbles is equal to the nucleation radius. When the number density of bubbles in each group At that time, based on the ideal gas law, the radii of each group of bubbles are calculated:

[0010] in, The radius of the bubble represents , ; This refers to the mass of a gas per unit volume, with units of: mol / m 3 , Bubble number density, in units of: pcs / m 3 , This is the gas constant, with units of: J / (mol·K) ; Calculate the number according to the following formula. The following grid regions show the bubble volume and cavitation fraction for each group:

[0011]

[0012] .

[0013] A further improvement of the present invention is that, in step (5), the grid number is calculated according to the following equation. Solution flow rate in the following areas:

[0014]

[0015] In this system of equations, the last term on the right-hand side represents the momentum dissipation term resulting from the interaction between the solution and the bubbles, which is related to... Proportional, of which Represents the momentum dissipation coefficient; The boundary conditions for the momentum equation are:

[0016] in, Atmospheric pressure, unit: Pa; Calculate the grid number according to the following equation. The solution density changes in the following regions:

[0017]

[0018]

[0019] Calculate the number based on the following relationship. Solution density within the following grid region: .

[0020] A further improvement of the present invention is that, in step (6), the grid number is calculated according to the following relationship. The isothermal compressibility and isobaric expansion coefficient of gas-liquid mixtures in the following regions:

[0021]

[0022] In the above formula, The isothermal compressibility coefficient of the pure liquid phase is given by: 1 / MPa , The isobaric expansion coefficient of a pure liquid phase is given by: 1 / ℃ , To control the characteristic radius of the air bubble within the body, it is calculated using the following formula: .

[0023] A further improvement of the present invention is that, in step (7), the solution pressure change is calculated based on the gas-liquid mixture state equation according to the following relationship:

[0024] Calculate the number based on the following relationship. Solution pressure within the following grid region: .

[0025] A further improvement of the present invention is that, in step (8), the total mass m of the solution below the J-th grid layer is calculated according to the following relationship:

[0026] Where the subscripts i and j represent the radial and axial grid positions, respectively; M is the total number of axial grids; The liquid level height is determined according to the following formula:

[0027] in, Let be the axial dimension of the j-th layer mesh, in meters (m). The void fraction distribution is re-given based on the following relationship: .

[0028] A further improvement of the present invention is that, in step (10), the distribution of the values ​​of each physical quantity is redetermined according to the following relationship:

[0029] in, Representative except External physical quantities: , , and For the concentration of dissolved radiolytic gases The solution region is a closed system, and it is necessary to ensure the conservation of the amount of dissolved radiolytic gas. The concentration distribution of the dissolved radiolytic gas should be re-given based on the following relationship:

[0030]

[0031] The following relationship is given. and Bubble radius between:

[0032] Among them, if Then j satisfies ;like Then j satisfies .

[0033] A further improvement of this invention is that, in step (11), the internal pressure of the bubble and the concentration of the radiolytic gas satisfy Heryy's law. When the internal pressure of the bubble is in equilibrium with the external pressure, the concentration of the radiolytic gas is called the equilibrium concentration. Based on the bubble radius, the equilibrium concentration is calculated using the following formula:

[0034] In the formula, To balance the concentration, Where is the bubble radius; Based on the ideal gas law, the gas density inside the bubble can be calculated using the following relationship:

[0035] in, The density of the gas inside the bubble is expressed in kg / m³.3 ; The molar mass of the radiolytic gas is given by: kg / mol ; The bubble velocity is described by the limiting upward velocity of the bubble. Based on the balance between drag and buoyancy, the limiting upward velocity of the bubble satisfies the following relationship:

[0036] In the formula, Let g be the bubble velocity, and g be the acceleration due to gravity. The drag coefficient satisfies the following relationship:

[0037]

[0038]

[0039] Where Re is the Reynolds number and Eo is the Etovos number. Dynamic viscosity; The interfacial mass transfer coefficient between the grown bubbles and the solution is calculated using the following formula:

[0040] in, The diffusion rate of gas in solution (m 2 / s), The contact time between the bubble and the surrounding fluid (s); Calculate the interphase interface area density using the following formula:

[0041] Calculate the bubble dissolution time for the second group based on the following formula: .

[0042] A further improvement of the present invention is that, in step (12), the conservation equations for the amount of gaseous substance, the number density of bubbles, and the concentration of radiolytic gas are as follows:

[0043]

[0044]

[0045]

[0046] In the above formula, The mass of gaseous matter within a single nucleating bubble, expressed in units of: mol / indivual; The volume of the control volume is expressed in units of: m 3 ; The percentage of nucleation bubbles that have grown out of the total number of nucleation bubbles, where is The correction factor in this ratio, θ(x) For the Heviside function, when x≤0, θ=0, x>0, θ=1 ; This represents the liquid phase volume fraction, and ; Calculate the mass of gas within a single nucleating bubble using the ideal gas law:

[0047] The boundary conditions for the radiolysis gas concentration at each boundary surface are:

[0048] Since the boundary surface is the liquid surface, radiolytic bubbles escape from the solution only through transport, and the flux of the diffusion term at this boundary surface is 0. Therefore, the liquid surface boundary condition is described by the following relationship: ,

[0049] ,

[0050] The boundary conditions on other boundary surfaces all satisfy the following relationship:

[0051] Where X is the relevant physical quantity, representing , , , , where n is the normal direction of other boundary surfaces.

[0052] Compared with the prior art, the present invention has at least the following beneficial technical effects: 1. The modeling and analysis method applicable to fissile aqueous solution systems described in this invention starts from the actual physical processes of fissile aqueous solution systems, considers the diffusion of dissolved radiolytic gases, bubble nucleation and bubble growth caused by gas-liquid mass transfer, and combines physical reality and reasonable simplification to provide a complete model of the bubble generation, growth and migration process.

[0053] 2. This invention can be used in simulation programs for fissile aqueous solutions, especially for cases where there is a large reactive inertial pressure pulse introduced by an instantaneous criticality. This invention can calculate the magnitude of the inertial pressure and simulate the distribution of radiolytic bubbles.

[0054] 3. Because commonly used radiolytic gas models do not consider the processes of bubble nucleation and growth, and rely on many empirical parameters, simulation programs for fissile aqueous solutions cannot accurately simulate both the first fission peak and broadening simultaneously. Furthermore, CFD two-phase flow simulation is complex, computationally intensive, and inefficient. This invention achieves accurate calculations without requiring a large number of meshes, has a small computational load, and is fast, ensuring both accuracy and computational efficiency. Attached Figure Description

[0055] Figure 1 This is a flowchart illustrating a method for calculating the cavitation fraction distribution of a two-dimensional radiolytic gas in a fissile aqueous solution system, as proposed in this invention.

[0056] Figure 2 This is an attempt to divide the annular cavity region inside the container into grids according to an embodiment of the present invention.

[0057] Figure 3 This is a comparison chart of the power and pressure calculated by the coupling program in the embodiment of the present invention and the experimental measurement results.

[0058] Figure 4 This is a comparison chart of the power and pressure calculated by the coupling program in the embodiment of the present invention and the experimental measurement results.

[0059] Figure 5 This is a cavitation fraction distribution diagram within each grid at 0s, 0.04s, 0.05s, and 0.55s in an embodiment of the present invention. Detailed Implementation

[0060] The present invention will now be described in further detail with reference to the accompanying drawings.

[0061] A type loaded with uranium concentration of The uranyl nitrate solution is labeled with the superscripts 0 and last to indicate the initial and previous values ​​of relevant parameters, respectively, and the initial temperature is... initial solution density In a cylindrical container, the inner diameter of the annular cavity... , outer diameter The method for calculating the cavitation fraction of two-dimensional radiolysis gas in the region includes the following steps: (1) Determine the geometric parameters, material properties, mesh generation parameters, and hydrogen production per unit energy of the solution system. (mol / J), turbidity coefficient (m) 2 / s) and bubble-related parameters include nucleation radius. ( The dissolution time of the first group of bubbles ( and the diffusion rate of radiolytic gas molecules m 2 / s; The grid is divided according to the input parameters, such as... Figure 2 As shown, the total number of grid cells along the axial direction is N, and the total number of grid cells in the radial direction is M. The initial power within the control volume (i,j) is given by neutron transport calculations. The initial temperature of the control body is given through thermal-hydraulic calculations. .

[0062] (2) The radiolysis bubbles are divided into two groups based on whether the nucleation bubbles grow in the solution. The nucleation bubbles and dissolution bubbles are in the first group, and the growing bubbles are in the second group. The following information about the bubbles needs to be stored using a two-dimensional array to reflect the two-dimensional spatial distribution characteristics of the variables: a. The radii of the two sets of bubbles, , The unit is: m ; b. The flow rates of the two sets of bubbles, , The unit is: m / s ; c. The number density of the two groups of bubbles. , The unit is: individual. / m 3 ; d. The mass of gas per unit volume of the two sets of bubbles. , The unit is: mol / m 3 ; e. Control the volume of the two sets of air bubbles inside the body. , The unit is: m 3 ; f. The gas density of the two sets of bubbles, , The unit is: kg / m 3 ; g. Radius of the nucleating bubble The unit is m equilibrium concentration of nucleated bubbles The unit is mol / m 3 First group of bubble dissolution time The unit is seconds (s). h. Equilibrium concentration of growth bubbles The unit is ( mol / m 3 ); the interfacial area density A between the growth bubbles and the solution (in m3); the interfacial mass transfer coefficient between the growth bubbles and the solution. Unit: m / s (3) Starting from time 0, the liquid level is The axial position number of the grid where the liquid surface is located Radiolysis gas concentration and , , , All values ​​are set to 0, solution pressure Taken as atmospheric pressure; solution velocity , Set to 0; based on the inner diameter of the cylindrical container. , outer diameter and liquid level and initial density of solution Calculate the total mass of the solution at the initial time. Use the value from the previous time step as the value from the initial time step. Total mass of solution at initial time The calculation formula is:

[0063] (4) According to , , , , Based on the ideal gas law, iterative calculation , ;according to , Calculate the volume of the two sets of radiolysis bubbles , ;according to , and control volume Calculate the void fraction ;make hour, , to be used to calculate the grid position number where the liquid level is located in step (8); a) When the number density of bubbles in each group At that time, the bubble radius of each group is taken as the nucleation radius. When the number density of bubbles in each group At that time, the radii of each group of bubbles can be calculated based on the ideal gas law:

[0064] in, The bubble radius can represent , ; Mass of gas per unit volume ( mol / m3 ), Bubble number density ( pcs / m 3 ), The gas constant ( J / (mol·K) ).

[0065] b) Calculate the number according to the following formula The following grid regions show the bubble volume and cavitation fraction for each group:

[0066]

[0067]

[0068] (5) Based on solution pressure The solution velocity at the previous moment , Based on the momentum conservation equation and corresponding boundary conditions, the solution flow rate is calculated. , ;according to , and Based on the continuity equation, the density change caused by solution motion is calculated. ;according to and The difference and the solution density at the previous moment Calculate the density change caused by the compressibility of radiolytic bubbles. ;according to and Calculate the change in solution density and according to Calculate the solution density This step only calculates the number. The physical quantities within the following grid region, for The grid area ; a) Ignore pressure changes caused by liquid level changes, grid number The momentum equation for the following region:

[0069]

[0070] In this system of equations, the last term on the right-hand side represents the momentum dissipation term resulting from the interaction between the solution and the bubbles, which is related to... Proportional, of which This represents the momentum dissipation coefficient. The boundary conditions for the momentum equation are:

[0071] in, Atmospheric pressure (Pa).

[0072] b) The grid number can be calculated using the following equation. The solution density changes in the following regions:

[0073]

[0074]

[0075] The number can be calculated based on the following relationship. Solution density within the following grid region:

[0076] (6) According to , , , , Calculation number The isothermal compressibility of the gas-liquid mixture within the following grid region and isobaric expansion coefficient ; Grid number In the following regions, the isothermal compressibility and isobaric expansion coefficient of the gas-liquid mixture are:

[0077]

[0078] In the above formula, The isothermal compressibility coefficient of the pure liquid phase ( 1 / MPa ), The isobaric expansion coefficient of pure liquid phase ( 1 / ℃ ), To control the characteristic radius of the air bubble within the body, it can be calculated using the following formula:

[0079] (7) According to and The difference , , , Based on the equation of state for gas-liquid mixtures, the number was calculated. Pressure changes within the following grid area Other regions always assume ;based on and the pressure of the previous moment Calculate solution pressure ; Equation of state for gas-liquid mixtures:

[0080] The number can be calculated based on the following relationship. Solution pressure within the following grid region:

[0081] (8) Obtained according to step (4) and the results obtained in step (5) Summing the solution masses within each control cell along the axial direction from bottom to top, the total solution mass below the J-th grid layer is calculated. until Obtain the current liquid level position number. ;Calculate the liquid level based on the conservation of the total mass of the solution. ;according to , , , and The void fraction distribution is re-given. Update the cavitation fraction and solution pressure to make , ; a) Total mass m of the solution below the J-th grid layer:

[0082] Where the subscripts i and j represent the radial and axial grid positions, respectively; M is the total number of axial grids.

[0083] b) Liquid level height:

[0084] in, Let be the axial dimension (m) of the j-th layer mesh.

[0085] c) Cavitation fraction:

[0086] (9) Repeat steps (4) to (8) until the maximum relative error of the control body solution pressure before and after the iteration is reached. It meets the convergence requirements and enables coupled calculation of solution flow rate, solution pressure, cavitation fraction, liquid level, and bubble radius; (10) According to and According to the physical quantity at the previous moment , , , and To obtain the redistribution of physical quantities after the liquid level change: , , , and ;according to , and , give and Bubble radius between , , Update the liquid level location number, and make ; a) Redistribute the values ​​of each physical quantity:

[0087] in, Representative except External physical quantities: , , and For the concentration of dissolved radiolytic gases The solution region is a closed system, and it is necessary to ensure the conservation of the amount of dissolved radiolytic gas. The concentration distribution of the dissolved radiolytic gas should be re-given based on the following relationship:

[0088]

[0089] b) Give and Bubble radius between:

[0090] Among them, if Then j satisfies ;like Then j satisfies .

[0091] (11) Calculate the number Within the following grid area, physical quantities other than c and d in step (2); based on the bubble radius , Solution pressure Herry constant H and surface tension Calculate the equilibrium concentration separately , According to the radius , and pressure Calculate the gas density of the two groups of bubbles respectively. and And based on these parameters and solution density Calculate the flow velocity sequentially and ;according to , , , and solution flow rate , Calculate the interphase interface mass transfer coefficient ;according to , Calculate the interfacial area density A; based on , , , , Calculate the dissolution time of the second group of bubbles. ; a) Equilibrium concentration:

[0092] In the formula, To balance the concentration, This represents the bubble radius. , .

[0093] b) Gas density inside the bubble:

[0094] in, The density of the gas inside the bubble (kg / m³) 3 ), The molar mass of the radiolytic gas ( kg / mol ).

[0095] c) The limiting upward velocity of the bubble satisfies the following relationship:

[0096] In the formula, Let g be the bubble velocity, and g be the acceleration due to gravity. The drag coefficient satisfies the following relationship:

[0097]

[0098]

[0099] Where Re is the Reynolds number and Eo is the Etovos number. This refers to dynamic viscosity.

[0100] d) Mass transfer coefficient at the interphase interface between the grown bubbles and the solution:

[0101] in, The diffusion rate of gas in solution (m 2 / s), The bubble time is the contact time (s) between the bubble and the surrounding fluid.

[0102] e) Interfacial area density:

[0103] f) Dissolution time of the second group of bubbles:

[0104] (12) Parameters obtained from the above steps , , , , , A, and input parameters , P, T, construct , , , and The differential equation; according to Determine the liquid level location and provide boundary conditions; calculate according to step (10). , , , and Calculate the current time number The physical quantities within the following grid region, , , , and For numbering In the region and above, all physical quantities are set to 0; a) The conservation equations for the amount of gaseous substance, the number density of bubbles, and the concentration of radiolytic gases are as follows:

[0105]

[0106]

[0107]

[0108]

[0109]

[0110]

[0111]

[0112] In the above formula, The mass of gas within a single nucleating bubble ( mol / indivual); To control the volume of the body ( m 3 ); The percentage of nucleation bubbles that have grown out of the total number of nucleation bubbles, where is The correction factor in this ratio, θ(x) For the Heviside function, when x≤0, θ=0, x>0, θ=1 ; The dissolution time of the first group of bubbles ( ); Liquid volume fraction ( According to the ideal gas law, the mass of gas within a single nucleating bubble satisfies the following relationship:

[0113] b) Boundary conditions: The boundary conditions for the radiolysis gas concentration at each boundary surface are:

[0114] Since the boundary surface is the liquid surface, radiolytic bubbles escape from the solution only through transport, and the flux of the diffusion term at this boundary surface is 0. Therefore, the liquid surface boundary condition can be described by the following relationship: ,

[0115] ,

[0116] The boundary conditions on other boundary surfaces all satisfy the following relationship:

[0117] Where X is the relevant physical quantity, which can represent , , , , where n is the normal direction of other boundary surfaces.

[0118] (13) By repeating steps (4) to (12), the linear differential equation system in step (12) is repeatedly constructed and solved to achieve coupled computation between the equation systems until And two iterations , , , and Maximum relative error The convergence requirement is met; (14) Use the calculated value at the current moment as the value at the previous moment to update the physical quantity: , , , , , , , , , and Output and Used for reactive feedback calculations; (15) Starting from time zero, advance the time step and repeatedly input the power at different times. , And continue to calculate according to steps (4) to (14) until the time required for transient calculation is reached.

[0119] Example An example is provided using the French experimental device SILENE, which has an inner radius of 0.038 m and an outer half-diameter of 0.18 m for its toroidal container. In experiment S2-173 of this series, the loading volume was 41.09 L, the uranium concentration was 71 g / L (92.7 wt%), and the acid concentration was 200 ml / m³. 3 The fuel solution was introduced with a reactivity of 2.98 g / m³ via a control rod. The radial region was divided into N=5 grids, and the axial region into M=26 grids. The number of grids below the initial liquid level was NL-1=10. Due to significant liquid level fluctuations in this experiment, the axial height of the liquid level fluctuation region was set to 0.5 m, with a grid size of 15. The calculation time step was 1Ee⁻⁵ s, the contact time between the bubble and the solution was 1E⁻⁵ s, and the initial solution density was 1163.44 kg / m³. 3 The surface tension is taken as 0.075 N / m. 2The Herry constant is taken as 7.258E-6, the source term coefficient as 5.0E-17, the momentum dissipation coefficient as 100, the hydrogen production rate as 1.35E-7 mol / J, and the mixing coefficient as 0.2m. 2 / s, the diffusion coefficient of radiolytic gas molecules is taken as 7.38E-9m. 2 / s. The isothermal compressibility coefficient of the solution is 4.5E-10. 1 / MPa isobaric expansion coefficient ( 1 / ℃ The following relationship must be satisfied:

[0120]

[0121] The calculation method of this invention is coupled with point pile dynamics and thermal-hydraulic models for calculation. A comparison chart of the power and pressure calculated by the coupled program and the experimental measurements is shown below. Figures 3-4 The cavitation fraction distribution within each grid at 0s, 0.04s, 0.05s, and 0.55s is as follows: Figure 5 As shown.

[0122] Although the present invention has been described in detail above with general descriptions and specific embodiments, modifications or improvements can be made to it, which will be obvious to those skilled in the art. Therefore, all such modifications or improvements made without departing from the spirit of the present invention fall within the scope of protection claimed by the present invention.

Claims

1. A method for calculating the cavitation fraction distribution of a two-dimensional radiolytic gas in a fissile aqueous solution system, characterized in that, Includes the following steps: (1) Determine the geometric parameters, material properties, mesh generation parameters, and hydrogen production per unit energy of the solution system. The unit is: mol / J turbulence coefficient The unit is m. 2 / s, and bubble-related parameters; divide the system into grids based on the input parameters, and determine the coordinates and dimensions of each grid; determine the system power distribution. The unit is W, and the temperature distribution. The unit is: ℃ ,in and To control the physical quantities within the body, all of them change with time and space; (2) The radiolysis bubbles in the solution are divided into three types: nucleation bubbles, dissolution bubbles, and growth bubbles. Nucleation bubbles are bubbles generated along the trajectory of fission fragments. Based on whether they grow in the solution, they are divided into two groups of bubbles. Nucleation bubbles and dissolution bubbles are the first group of bubbles, and growth bubbles are the second group of bubbles. The following information about the bubbles needs to be stored using a matrix to reflect the two-dimensional spatial distribution characteristics of the variables: a. The radii of the two sets of bubbles, , The unit is: m ; b. The flow rates of the two sets of bubbles, , The unit is: m / s ; c. Number density of the two groups of bubbles , The unit is: individual. / m 3 ; d. The mass of gas per unit volume of the two sets of bubbles. , The unit is: mol / m 3 ; e. Controlling the volume of the two sets of air bubbles inside the body. , The unit is: m 3 ; f. The gas density of the two sets of bubbles, , The unit is: kg / m 3 ; g. Radius of the nucleating bubble The unit is m equilibrium concentration of nucleated bubbles The unit is mol / m 3 First group of bubble dissolution time The unit is seconds (s). h. Equilibrium concentration of growth bubbles The unit is mol / m 3 ; A, the interfacial area density between the growth bubbles and the solution (in m³); mass transfer coefficient between the growth bubbles and the solution (in m³). The unit is m / s; (3) Use the superscript 0 and last to mark the initial value and the value at the previous moment of the relevant parameters; start the calculation from moment 0, the liquid level is The axial position number of the grid where the liquid surface is located Radiation gas concentration and , , , All values ​​are set to 0, solution pressure Taken as atmospheric pressure; solution velocity , Set to 0; based on the inner diameter of the cylindrical container. , outer diameter and liquid level and initial density of solution Calculate the total mass of the solution at the initial time. ;Assign the value from the previous time step to the value from the initial time step; (4) According to , , , , Based on the ideal gas law, iterative calculation , ;according to , Calculate the volume of the two sets of radiolysis bubbles , ;according to , and control volume Calculate the void fraction ;make hour, , to be used to calculate the grid position number where the liquid level is located in step (8); (5) Based on solution pressure The solution velocity at the previous moment , Based on the momentum conservation equation and corresponding boundary conditions, the solution flow rate is calculated. , ;according to , and Based on the continuity equation, the density change caused by solution motion is calculated. ;according to and The difference and the solution density at the previous moment Calculate the density change caused by the compressibility of radiolytic bubbles. ;according to and Calculate the change in solution density and according to Calculate the solution density This step only calculates the number. The physical quantities within the following grid region, for The grid area ; (6) According to , , , , Calculation number The isothermal compressibility of the gas-liquid mixture within the following grid region and isobaric expansion coefficient ; (7) According to and The difference , , , Based on the equation of state for gas-liquid mixtures, the number was calculated. Pressure changes within the following grid area Other regions always assume ;based on and the pressure of the previous moment Calculate solution pressure ; (8) Obtained according to step (4) and the results obtained in step (5) Summing the solution masses within each control cell along the axial direction from bottom to top, the total solution mass below the J-th grid layer is calculated. until Obtain the current liquid level position number. ; Calculate the liquid level based on the conservation of the total mass of the solution. ;according to , , , and The void fraction distribution is re-given. ; Update the cavitation fraction and solution pressure to make , ; (9) Repeat steps (4) to (8) until the maximum relative error of the control body solution pressure before and after the iteration is reached. It meets the convergence requirements and enables coupled calculation of solution flow rate, solution pressure, cavitation fraction, liquid level, and bubble radius; (10) According to and According to the physical quantity at the previous moment , , , and To obtain the redistribution of physical quantities after the liquid level change: , , , and ;according to , and , give and Bubble radius between , , Update the liquid level location number, and make ; (11) Calculate the number Within the following grid area, physical quantities other than c and d in step (2); based on the bubble radius , Solution pressure Herry constant H and surface tension Calculate the equilibrium concentration separately , According to the radius , and pressure Calculate the gas density of the two groups of bubbles respectively. and And based on these parameters and solution density Calculate the flow velocity sequentially and ;according to , , and solution flow rate , Calculate the interphase interface mass transfer coefficient ;according to , Calculate the interfacial area density A; (12) Parameters obtained from the above steps , , , , , A, and input parameters , P, T, construct , , , and The differential equation; according to Determine the liquid level location and provide boundary conditions; calculate according to step (10). , , , and Calculate the current time number The physical quantities within the following grid region, , , , and For numbering In the region and above, all physical quantities are set to 0; (13) By repeating steps (4) to (12), the linear differential equation system in step (12) is repeatedly constructed and solved to achieve coupled computation between the equation systems until And two iterations , , , and Maximum relative error The convergence requirement is met; (14) Use the calculated value at the current moment as the value at the previous moment to update the physical quantity: , , , , , , , , , and Output and Used for reactive feedback calculations; (15) Starting from time zero, advance the time step and repeatedly input the power at different times. , And continue to calculate according to steps (4) to (14) until the time required for transient calculation is reached.

2. The method for calculating the cavitation fraction distribution of a two-dimensional radiolytic gas in a fissile aqueous solution system according to claim 1, characterized in that, In step (1), the geometric parameters of the solution system are a ring-shaped solution region. , outer diameter and liquid level Height of the air cavity region inside the cylindrical container The physical property of the solution material is density. Surface tension Herry constant H, uranium concentration The mesh generation parameters are the axial and radial mesh counts in the solution region and the axial mesh count in the gas cavity region; the bubble-related parameters are the nucleation radius. The dissolution time of the first group of bubbles Radiation gas molecule diffusion rate m 2 / s.

3. The method for calculating the cavitation fraction distribution of a two-dimensional radiolytic gas in a fissile aqueous solution system according to claim 1, characterized in that, In step (3), the total mass of the solution at the initial moment is calculated according to the following relationship: In step (4), when the number density of bubbles in each group At that time, the radius of each group of bubbles is equal to the nucleation radius. When the number density of bubbles in each group At that time, based on the ideal gas law, the radii of each group of bubbles are calculated: in, The radius of the bubble represents , ; This refers to the mass of a gas per unit volume, with units of: mol / m 3 , Bubble number density, in units of: pcs / m 3 , This is the gas constant, with units of: J / (mol·K) ; Calculate the number according to the following formula. The following grid regions show the bubble volume and cavitation fraction for each group: 。 4. The method for calculating the cavitation fraction distribution of a two-dimensional radiolytic gas in a fissile aqueous solution system according to claim 1, characterized in that, In step (5), the grid number is calculated according to the following equation. Solution flow rate in the following areas: In this system of equations, the last term on the right-hand side represents the momentum dissipation term resulting from the interaction between the solution and the bubbles, which is related to... Proportional, of which Represents the momentum dissipation coefficient; The boundary conditions for the momentum equation are: in, Atmospheric pressure, unit: Pa; Calculate the grid number using the following equation. The solution density changes in the following regions: Calculate the number based on the following relationship. Solution density within the following grid region: 。 5. The method for calculating the cavitation fraction distribution of a two-dimensional radiolytic gas in a fissile aqueous solution system according to claim 1, characterized in that, In step (6), the grid number is calculated according to the following formula. The isothermal compressibility and isobaric expansion coefficient of gas-liquid mixtures in the following regions: In the above formula, The isothermal compressibility coefficient of the pure liquid phase is given by: 1 / MPa , The isobaric expansion coefficient of a pure liquid phase is given by: 1 / ℃ , To control the characteristic radius of the air bubble within the body, it is calculated using the following formula: 。 6. The method for calculating the cavitation fraction distribution of a two-dimensional radiolytic gas in a fissile aqueous solution system according to claim 1, characterized in that, In step (7), the solution pressure change is calculated based on the equation of state for the gas-liquid mixture according to the following relationship: Calculate the number based on the following relationship. Solution pressure within the following grid region: 。 7. The method for calculating the cavitation fraction distribution of a two-dimensional radiolytic gas in a fissile aqueous solution system according to claim 1, characterized in that, In step (8), the total mass m of the solution below the J-th grid layer is calculated according to the following formula: Where the subscripts i and j represent the radial and axial grid positions, respectively; M is the total number of axial grids; The liquid level height is determined according to the following formula: in, Let be the axial dimension of the j-th layer mesh, in meters (m). The void fraction distribution is re-given based on the following relationship: 。 8. The method for calculating the cavitation fraction distribution of a two-dimensional radiolytic gas in a fissile aqueous solution system according to claim 1, characterized in that, In step (10), the distribution of the values ​​of each physical quantity is given again according to the following relationship: in, Representative except External physical quantities: , , and For the concentration of dissolved radiolytic gases The solution region is a closed system, and it is necessary to ensure the conservation of the amount of dissolved radiolytic gas. The concentration distribution of the dissolved radiolytic gas should be re-given based on the following relationship: The following relationship is given. and Bubble radius between: Among them, if Then j satisfies ;like Then j satisfies .

9. The method for calculating the cavitation fraction distribution of a two-dimensional radiolytic gas in a fissile aqueous solution system according to claim 1, characterized in that, In step (11), the internal pressure of the bubble and the concentration of the radiolytic gas satisfy Heryy's law. When the internal pressure of the bubble is in equilibrium with the external pressure, the concentration of the radiolytic gas is called the equilibrium concentration. The equilibrium concentration is calculated based on the bubble radius using the following formula: In the formula, To balance the concentration, Where is the bubble radius; Based on the ideal gas law, the gas density inside the bubble can be calculated using the following relationship: in, The density of the gas inside the bubble is expressed in kg / m³. 3 ; The molar mass of the radiolytic gas is given by: kg / mol ; The bubble velocity is described by the limiting upward velocity of the bubble. Based on the balance between drag and buoyancy, the limiting upward velocity of the bubble satisfies the following relationship: In the formula, Let g be the bubble velocity, and g be the acceleration due to gravity. The drag coefficient satisfies the following relationship: Where Re is the Reynolds number and Eo is the Etovos number. Dynamic viscosity; The interfacial mass transfer coefficient between the grown bubbles and the solution is calculated using the following formula: in, The diffusion rate of gas in solution (m 2 / s), The contact time between the bubble and the surrounding fluid (s); Calculate the interphase interface area density using the following formula: Calculate the bubble dissolution time for the second group based on the following formula: 。 10. The method for calculating the cavitation fraction distribution of a two-dimensional radiolytic gas in a fissile aqueous solution system according to claim 1, characterized in that, In step (12), the conservation equations for the amount of gaseous substance, the number density of bubbles, and the concentration of radiolytic gas are as follows: In the above formula, The mass of gaseous matter within a single nucleating bubble, expressed in units of: mol / indivual; The volume of the control volume is expressed in units of: m 3 ; The percentage of nucleation bubbles that have grown out of the total number of nucleation bubbles, where is The correction factor in this ratio, θ(x) For the Heviside function, when x≤0, θ=0, x>0, θ=1 ; This represents the liquid phase volume fraction, and ; Calculate the mass of gas within a single nucleating bubble using the ideal gas law: The boundary conditions for the radiolysis gas concentration at each boundary surface are: Since the boundary surface is the liquid surface, radiolytic bubbles escape from the solution only through transport, and the flux of the diffusion term at this boundary surface is 0. Therefore, the liquid surface boundary condition is described by the following relationship: ， ， The boundary conditions on other boundary surfaces all satisfy the following relationship: Where X is the relevant physical quantity, representing , , , , where n is the normal direction of other boundary surfaces.

Images