CN114678089A - Method for determining morphology of irradiation bubbles in nuclear material and influence of irradiation bubbles on mechanical and thermal properties - Google Patents

Abstract

The invention discloses a method for determining the appearance of irradiation bubbles in a nuclear material and the influence of the irradiation bubbles on the mechanical and thermal properties, which comprises the steps of collecting physical property parameters of the nuclear material required in the irradiation process, and obtaining elastic energy density by solving an elastic equilibrium equation; constructing a free energy equation by combining the elastic energy density; establishing and solving a rate theoretical model of irradiation bubble nucleation to obtain initial bubble distribution, and using the initial bubble distribution as an initial condition of phase field simulation; constructing and solving a phase field evolution equation on the basis of introducing a polycrystalline structure and considering temperature gradient, interface energy anisotropy and diffusion coefficient anisotropy; carrying out visualization processing on the numerical solution of the phase field evolution equation to obtain the appearance of bubble nucleation growth and grain growth in the nuclear material irradiation process; counting all bubbles to obtain the size distribution and porosity parameters of the irradiated bubbles in the nuclear material; by counting the size distribution and porosity of the irradiation bubbles, the irradiation swelling and heat conduction rule of the nuclear material is obtained, so that the service performance of the nuclear material is predicted.

Description

Method for determining morphology of irradiation bubbles in nuclear material and influence of irradiation bubbles on mechanical and thermal properties

Technical Field

The invention relates to the field of nuclear materials, in particular to a method for determining the appearance of irradiation bubbles in a nuclear material and the influence of the irradiation bubbles on the mechanical and thermal properties.

Background

Nuclear materials are closely related to nuclear reactor safety. During operation of nuclear reactors, nuclear materials are subjected to harsh operating environments, and fission reactions, accompanied by the development of fuel consumption, release a large amount of fission products, including solids and fission gases, wherein the fission gases account for about 25% of the total atomic number of all fission products. Because of the very large migration energy and very low solubility of gas atoms, fission gases generated during nuclear reactions are often trapped in voids or grain boundaries, forming intra-or inter-crystalline irradiation bubbles. A large number of bubbles are continuously gathered on the grain boundary, finally, a special through hole structure is formed on the grain boundary, and meanwhile, the volume of a fission product is larger than that of a substance before fission, so that the volume of a material element can be increased along with the development of burnup, and the change is called radiation swelling. Radiation swelling of the nuclear fuel can cause interaction of the fuel and cladding and stress on the cladding, causing radial deformation and lateral stretching of the cladding tube, ultimately resulting in cladding failure and serious threat to the operational safety of the nuclear reactor. The influence of the solid fission product on the irradiation swelling is well researched; however, the behavior of fission gases is complex, and further research into this area is needed. Therefore, it is very important to study the growth of the irradiation bubble nucleation of the nuclear material and the influence of the irradiation bubble nucleation on the irradiation swelling. In addition, because nuclear materials are often in service conditions of high temperature and extreme irradiation, the generated thermal strain and creep strain can greatly increase the stress in the materials, and therefore, the research on the evolution of irradiation bubbles of the materials under the action of the stress is particularly important.

Disclosure of Invention

The invention aims to provide a phase field simulation method for researching the formation and evolution of intragranular bubbles of a nuclear material after the nuclear material is irradiated and the swelling of the nuclear material, namely a method for determining the appearance of the irradiated bubbles in the nuclear material and the influence of the irradiated bubbles on the mechanical and thermal properties.

In order to achieve the purpose, the invention adopts the following technical scheme:

a method for determining the morphology of irradiation bubbles in a nuclear material and the influence thereof on the mechanical and thermal properties comprises the following steps:

s1, selecting phase field variables describing the evolution process of the nuclear material irradiation bubble morphology; collecting physical parameters of the nuclear material, including the formation energy, diffusion coefficient and gradient coefficient of point defects, and the mobility of bubbles and grain boundaries; and calculating the mobility of the point defect considering the diffusion anisotropy;

s2, calculating the elastic modulus distribution related to the position; combining Van der Waals equation to obtain initial intrinsic strain distribution; establishing and solving an elastic balance equation to obtain the stress-strain distribution and the elastic energy density of the irradiation bubbles;

s3, obtaining a free energy density function of a matrix phase of the nuclear material based on a material thermodynamic theory according to the physical parameters of the nuclear material collected in the step S1, obtaining a free energy density function of a bubble phase in the nuclear material by deducing a Van der Waals equation under the conditions of high temperature and high pressure, and obtaining a total free energy equation by combining interface gradient energy, bubble-crystal boundary interaction energy, polycrystal interaction energy and elastic energy density obtained in the step S2;

s4, establishing a speed theoretical model of irradiation bubble nucleation; obtaining bubble size distribution at the initial stage of irradiation bubble nucleation by solving a rate theoretical model, and taking the size distribution as an initial condition of phase field simulation;

s5, constructing a phase field evolution equation according to the total free energy equation obtained in the step S3 and comprehensively considering factors of temperature gradient, generation of point defects, compounding of the point defects and disappearance of the point defects at the defect trap under the irradiation condition; carrying out numerical solution on a phase field evolution equation by adopting a semi-implicit Fourier spectrum method;

s6, carrying out visualization processing on the numerical solution obtained in the step S5 to obtain the appearance of bubble nucleation growth and grain growth in the nuclear material irradiation process; counting all bubbles by adopting a communicated region marking method to obtain the size distribution and porosity parameters of the irradiated bubbles in the nuclear material; carrying out classified statistics on the intercrystalline bubbles and the intragranular bubbles;

s7, calculating and obtaining the effective thermal conductivity and the radiation swelling rule of the nuclear material based on the bubble distribution and the porosity obtained in the step S6.

In one implementation method, S1 specifically includes:

s1.1, selecting phase field variables for describing the evolution process of the nuclear material irradiation bubble morphology, wherein the phase field variables for expressing the point defect concentration need to be selected: concentration of gas atoms c_gVacancy concentration c_vAnd interstitial atom concentration c_i(ii) a And simultaneously selecting a bubble phase sequence parameter eta for distinguishing a bubble phase from a matrix phase: η is 1.0 when in the bubble phase and 0.0 when in the matrix phase; the parameter phi of multiple crystal sequences needs to be selected_iDenotes the grain structure, wherein the index i denotes the ith grain and within the ith grain is_i＝1.0；

S1.2, collecting physical parameters of the nuclear material, including formation energy, diffusion coefficient, mobility and gradient coefficient of gas atoms, vacancies and interstitial atoms, and mobility of bubbles and grain interfaces;

s1.3 calculating the mobility of the point defect taking into account the diffusion anisotropy, the vacancy mobility M, from the diffusion coefficients of the gas atoms, vacancies and interstitial atoms of S1.2_vThe following expression is given:

wherein D_vVacancy diffusion coefficient in tensor form, V_mIs unit volume of crystal lattice, R is ideal gas constant, T is absolute temperature;

wherein D_v ⁰In order to have an effective diffusion coefficient for the vacancies,

and

the components of vacancy diffusion, surface diffusion and grain boundary diffusion are respectively as follows:

wherein w_b、w_sAnd w_GBWeight coefficients of bulk diffusion, surface diffusion and grain boundary diffusion, I is unit matrix, T_s ^vAnd T_GB ^vThe tensor coefficients of the surface and the grain boundary are respectively expressed as the following expressions:

in one implementation method, S2 specifically includes:

s2.1 calculate the position-dependent elastic modulus distribution: in polycrystalline materials, the modulus of elasticity is related to the orientation of the grains, the modulus of elasticity of each oriented grain being determined by coordinate transformation, the tensor C of the modulus of elasticity_ijkl(r) the expression is as follows:

wherein the content of the first and second substances,

is the modulus of elasticity in the matrix and,

and

respectively, the rotation matrix a of the p-th crystal grain to the reference coordinate_ijIa, jb, kc, ld components in (1), rotation matrix a_ijThe expression is determined by the Euler angle of the corresponding crystal grain, and the rotation torque matrix under the three-dimensional coordinate is as follows:

in the formula, theta, xi and psi are three Euler angles for determining the orientation of crystal grains;

s2.2 initial intrinsic Strain

Consisting of two parts, the first part being the intrinsic strain of the matrix phase

The second part being the intrinsic strain of the bubble phase

r is a position vector;

in which the intrinsic strain of the matrix phase is caused by defect inhomogeneities

Comprises the following steps:

wherein epsilon^g0Is the expansion coefficient of lattice constant, delta, due to the introduction of gas atoms_ijIs a function of the Kronecker function,

is the equilibrium concentration of gas atoms, h (η) ═ η³(6η²-15 η +10) is an interpolation function constructed such that when η is 1.0, h (η) is 1.0 in the bubble phase; when η is 0.0, h (η) is 0.0 in the matrix phase; the intrinsic strain of the bubble is very complex because it depends on the internal pressure of the bubble and the size of the bubble, when the internal pressure of the bubble is greater than 2 gamma_sAt/r, the outward pressure acts on the lattice around the bubble, so the characteristic strain of the bubble is positive, which means that the bubble is in an overpressure state; when the bubble pressure is equal to 2 gamma_sAt/r, the intrinsic strain is zero; a negative intrinsic strain value indicates that the bubble is in an under-pressure state, which causes elastic relaxation of the crystal lattice around the bubble towards the center of the bubble; gamma ray_sIs the surface energy of the bubble, and r is the radius of the bubble;

intrinsic strain of bubble phase

The description is as follows:

wherein, P_gIs the internal pressure of the bubbles, 2 gamma_sPer considered as external hydrostatic pressure, C₁₁The gas pressure inside the bubble is always positive, being the modulus of elasticity in the main direction, while the intrinsic strain inside the bubble may be negative; this definition indicates that intrinsic strain is a function of gas pressure, bubble radius, surface tension and material elastic properties; the gas pressure in the bubbles can be obtained using the van der waals equation:

where Ω is the molecular volume of the corresponding core material, b is the Van der Waals constant, k_BBoltzmann constant;

s2.3, stress-strain distribution is obtained by solving an elastic balance equation which is as follows:

in the formula sigma_ijIs stress, r_jIs the unit displacement in the j direction; epsilon_kl(r) is the total strain equivalent to ε_ij(r); total strain epsilon_ij(r) is the sum of the uniform and non-uniform strains:

in the formula

And δ ε_ij(r) uniform strain and non-uniform strain, respectively; uniform strain

The expression is as follows:

in the formula (I), the compound is shown in the specification,

for additional stress, V is the simulation zone volume; non-uniform strain delta epsilon_ij(r) the expression is as follows:

in the formula u_i(r) and u_j(r)Substituting equation (17) for the displacement components in the i and j directions into elastic equilibrium equation (14) yields:

obtaining a displacement field equation after simplifying and sorting:

in the formula

Is a position-independent elastic modulus, C'_ijkl(r) is the position-dependent modulus of elasticity; the following results are obtained after fourier transform:

in the formula

Is the displacement field of Fourier space, k is the wave vector in Fourier space, k_jIs a component of a wave vector;

is the Green tensor; the subscript k indicates that this partial computation is performed in fourier space; equation (20) is solved using perturbation iteration, and the displacement field in the nth iteration is as follows:

the initial value of the iteration is taken as

Displacement field u in real space_k(r) is prepared from

Performing inverse Fourier transform, obtaining a displacement field by iteratively solving an elastic equilibrium equation, and further calculating to obtain stress-strain distribution;

calculating an elastic energy density from the obtained stress strain, the elastic energy density having the following expression:

in one implementation method, S3 specifically includes:

s3.1, based on the theoretical derivation of thermodynamics, calculating a free energy density function of a matrix phase of the nuclear material according to the following formula:

wherein the condition c is always satisfied when deriving the free energy of the matrix, since it is assumed that the defects occupy only the perfect lattice_g+c_v+c_i+c_m＝1.0；c_mIs the concentration of a perfect lattice, in at.%; e_g ^f、E_v ^f、E_i ^fRespectively the formation energy of gas atoms, vacancies and interstitial atoms in the matrix phase;

s3.2 the free energy density function of the bubble phase of the nuclear material is as follows:

wherein A and B are constants,

is a function of the density of free energy within the bubble in relation to the concentration of gas atoms; the Gibbs free energy expression G (p) of the bubbles is:

wherein p is pressure; v_bIs the volume of the bubble; g (p)₀) Is at the reference state pressure p₀Lower gibbs free energy; the van der waals equation at high pressure is:

p(V_b-nb)＝nk_BT (26)

wherein n is the number of gas atoms; the gibbs free energy from which the non-ideal gas is obtained is:

the gas atoms may occupy lattice or non-lattice sites in the bubble phase, so the concentration of gas atoms in the bubble is defined as:

wherein V_siteIs the volume of a single crystal lattice, from which the pressure:

the free energy density function of the bubble that brings the above results into account is:

wherein, mu₀Is a reference state p₀Chemical potential of (a);

s3.3, introducing an interpolation function and a double-trap potential function to obtain a bulk free energy density function according to the obtained free energy density functions of the matrix phase and the bubble phase:

f_bulk(c_g,c_v,c_i,η,T)＝h(η)f_b(c_g,c_v,c_i,T)+[1-h(η)]f_m(c_g,c_v,c_i,T)+wg(η) (31)

wherein h (η) ═ η³(6η²-15 η +10) is an interpolation function constructed such that when η is 1.0, h (η) is 1.0 in the bubble phase; when η is 0.0, h (η) is 0.0 in the matrix phase; g (η) ═ η²(1-η)²Is a double potential well function, w is the potential well height;

the polycrystalline interaction free energy density function is expressed as follows:

where α, β, γ and δ are image-only parameters.

The total free energy equation is obtained by combining a free energy density function, an elastic energy density, an interface gradient energy, a bubble and grain boundary interaction energy and a polycrystal interaction energy:

in the formula, κ_g、κ_v、κ_i、κ_ηThe gradient term coefficients are respectively gas atom concentration, vacancy concentration, interstitial atom concentration and bubble phase sequence parameters, and the last four terms of the formula represent gradient energy.

In one implementation method, S4 specifically includes:

s4.1, establishing a speed theoretical model of irradiation bubble nucleation growth behavior by combining the physical parameters of the nuclear material obtained in S1: the irradiation of nuclear material generates a large number of gas atoms which tend to accumulate in low density areas in order to reduce the energy of the nuclear material and minimize strain, a process known as nucleation; the change in free energy during nucleation is caused by a spherical-radius core, and this free energy Δ F is represented by the following equation:

wherein Δ f is the difference between the free energies of the bubble and the substrate, r is the bubble radius, and σ is the surface energy; while a critical kernel corresponds to the maximum of this free energy, since any larger kernel will have less energy and will grow spontaneously; thereby obtaining the critical radius R of the nucleus^CAnd activation energy Δ F^C：

Randomly introducing bubble cores in a nucleation stage, and ensuring that the average nucleation rate is matched with the expected nucleation rate; assuming that all points in the matrix can become nucleation sites, the probability of forming a core at one atomic position in one nucleation time interval can be calculated by considering the nucleation rate:

wherein Z is a Zeldovich correction factor, N is the atomic number, beta^*Is the probability of the critical core changing into the supercritical core, and tau is the inoculation time of the core; t is t₀Time intervals for calculating nucleation; the total nucleation probability P is thus obtained:

P＝1-exp(-J^*t₀) (38)

s4.2, solving a formula (34) and a formula (38) to obtain nucleation probabilities corresponding to bubbles with different sizes, namely the size distribution of the bubbles at the initial stage of irradiation bubble nucleation; in the rate theory nucleation algorithm, the bubble radius and the nucleation probability of each nucleation point are calculated, and the bubble size distribution at the initial stage of irradiation bubble nucleation is used as an initial condition to be input into the phase field simulation.

In one implementation method, S5 specifically includes:

s5.1 considering temperature gradient under irradiation condition, generation of point defect, compounding of point defect and disappearance process evolution equation of point defect at defect trap are as follows:

evolution equation of gas atomic concentration:

evolution equation of vacancy atom concentration:

evolution equation of interstitial concentration:

evolution equation of the bubble phase field variable:

evolution equation of the polycrystal phase field variable:

in the formula, t is simulation time; f is total free energy; l is the free interface mobility; m_g、M_v、M_iAtomic mobilities of gas atoms, vacancy atoms, and interstitial atoms, respectively;

thermal fluctuation terms of gas atoms, vacancy atoms, interstitial atoms and bubble phases respectively;

respectively the generation rates of gas atoms and point defects under irradiation conditions;

is the recombination rate of point defect vacancies and interstitial atoms,

is the absorption term of point defects by grain boundaries. The influence of the temperature gradient on bubble migration and grain growth is considered in the evolution equation, wherein Q is heat transfer efficiency, and C is a constant;

and

the production rates of gas atoms, vacancies and interstitial gas atoms under irradiation conditions respectively; the generation rate of gas atoms is:

in the formula (f)_rIs the rate of fission; ran is a random number between 0 and 1; Λ is a constant;

the generation of interstitial atoms and vacancies uses the following expression:

wherein R is₁And R₂Two random numbers uniformly distributed between 0 and 1 randomly generated at each time step and each space lattice point; e is a biased constant that can be varied to represent unequal numbers of vacancies and interstitials created by an ex-situ damage; parameter P_cascRepresenting the probability of a cascade collision per unit time and unit volume, V_GRepresenting empty spaces due to cascading collision eventsA maximum increment; in the formula eta<The condition of 0.8 ensures that the cascade collisions occur only in the matrix phase, and not in the bubble phase;

when the point defect vacancies and interstitial atoms meet, they recombine to form a perfect lattice, expressed as the recombination rate:

in the formula, v_rThe recombination rate is: v is_r＝v_b+η²v_s，v_bAnd v_sThe recombination rates of point defects at bulk and interface, respectively; v. of_b＝4πr_iv(D_i+D_v)/Ω，r_ivIs the composite volume radius; d_iAnd D_vThe diffusion coefficients of interstitial atoms and vacancies, respectively; Ω is the lattice volume;

for the absorption term of the point defects by the grain boundary, the expression is as follows:

wherein

The grain boundary absorption factor represents the strength of absorption of point defects by grain boundaries; phi ═ Σ phi_i ²As a function of the position of the crystal grains, phi is 1.0 in the crystal grains, and phi is less than 1.0 at the crystal boundary,

is the equilibrium concentration of gas atoms, vacancies, interstitial atoms;

s5.2, solving a phase field evolution equation by using a semi-implicit Fourier spectrum method, wherein the solution of the evolution equation of the gas atom concentration, the vacancy concentration and the interstitial atom concentration is shown as the following formula:

the evolution equation for the bubble phase field variable and the polycrystal phase field variable is solved as shown in the following formula:

wherein the superscript n represents the value of the portion in the nth time step, the free energy f is the bulk free energy density, the sum of the polycrystalline free energy density and the elastic free energy density, and k is (k ═ k)₁,k₂) Is the vector coordinate of the Fourier space, and delta t is the simulated time step;

in one implementation method, S6 specifically includes:

s6.1, carrying out visualization processing on the numerical solution result in the step S5 by using visualization software to obtain the appearance of bubble nucleation growth and grain growth in the nuclear material irradiation process:

representing visual variables of an irradiation bubble evolution simulation process:

in the formula (I), the compound is shown in the specification,

is a defined visual variable; the characteristics of the bonding diffusion interface are known: inside the crystal grain

At grain boundaries

In the air bubble

Therefore, the crystal grains, the crystal boundaries and the bubbles are distinguished by distinguishing the visual variable values of each point in the space;

s6.2, counting the size distribution of the bubbles in the output result by using an algorithm of the link region mark, wherein the method comprises the following steps:

in the connected domain marking, scanning from left to right and from top to bottom during the first marking, each effective pixel is set to be a label value, and the judgment rule is as follows:

(1) when the left adjacent pixel and the upper adjacent pixel of the pixel are invalid values, setting a new label value, label +1, for the pixel;

(2) when the left-adjacent pixel or the upper-adjacent pixel of the pixel has a valid value, assigning label of the valid value pixel to the label value of the pixel;

(3) when the left adjacent pixel and the upper adjacent pixel of the pixel are both effective values, selecting the smaller label value to be assigned to the label value of the pixel;

finally, counting the area of the region contained in each label value to obtain the size distribution and porosity parameters of the bubbles;

s6.3, counting the intercrystalline bubbles, and realizing by using the following algorithm:

(1) numbering each bubble on the basis of S6.2 statistics, and recording the position of each bubble;

(2) judging whether a certain bubble is an intergranular bubble, wherein the method comprises the following steps: taking the range larger than the bubble as a calculation region, and calculating a multi-order parameter phi_iThe sum Γ within this range_iFinally, multiply it to obtain

If lambda is 0.0, the bubble is in the grain; if λ ≠ 0.0 then the bubble is at the grain boundary;

(3) sequentially judging whether all bubbles are intercrystalline bubbles, respectively obtaining intercrystalline bubble distribution and intragranular bubble distribution, and obtaining intercrystalline bubble coverage rate by using the area ratio of the intercrystalline bubbles to the total area of a crystal boundary;

in one implementation method, S7 specifically includes:

s7.1 calculating the effective thermal conductivity of the nuclear material, the heat transfer equation needs to be solved:

in the formula k_localFor local thermal conductivity, k is a function of the heat transfer properties in spatial position, when located inside a grain_local＝k_bulk，k_bulkIs the thermal conductivity within the grains; when located on grain boundaries k_local＝k_GB，k_GBIs the thermal conductivity at the grain boundaries; when located in the bubble k_local＝k_bubble，k_bubbleIs the thermal conductivity within the bubble; solving a heat transfer equation (56) by using a finite difference method to obtain an equilibrium temperature field and effective thermal conductivity;

s7.2, combining the size distribution and the porosity of the bubbles calculated in the step S6 to obtain an irradiation swelling rule of the nuclear material irradiated bubbles, wherein the fuel swelling caused by the irradiated bubbles in the nuclear material is estimated by the following formula:

where Δ V is the change in fuel volume, V₀Is the initial fuel volume, V_fThe fuel swelling was calculated in a two-dimensional simulation as the ratio of the bubble area to the total area for the final fuel volume including the bubble-induced swelling.

Compared with the prior art, the invention has the following advantages:

by the method, the simulation precision of the evolution process of the irradiation bubbles under the irradiation condition is greatly improved; the influence of various objective conditions on the evolution of the irradiation bubbles is fully considered, and the application range of the phase field model is expanded;

furthermore, by introducing and solving an elastic equilibrium equation, the internal pressure and the elastic energy of the bubbles are fully considered, so that the simulation process is closer to the actual physical process;

furthermore, the bubble phase free energy density is deduced from the Van der Waals equation, so that the defect of using an empirical function in the traditional phase field simulation is avoided, and the simulation precision of the bubble evolution process is improved;

furthermore, a connected region marking method is adopted to accurately count and classify the irradiation bubbles in evolution, so that the accuracy of a simulation result is improved;

furthermore, the influence rule of factors such as grain size, external stress, anisotropic diffusion coefficient and the like on the nucleation and growth process of the irradiation bubble and the thermal performance of the irradiation bubble is fully considered in the simulation, and the application range of the phase field model is expanded.

Furthermore, the invention is implemented by self programming of Fortran language, and also can be implemented by software such as Matlab and the like, and has good flexibility and expansibility.

Drawings

FIG. 1 is a flow chart of the method of the present invention.

FIG. 2a, FIG. 2b, FIG. 2c are the stress components σ of a single irradiation bubble, respectively_xxDistribution, stress component σ_yyDistribution, stress component σ_xyAnd (4) distribution.

FIG. 3 is a graph of mean bubble diameter as a function of irradiation intensity and its comparison with experimental results.

FIG. 4a, FIG. 4b, and FIG. 4c show the distribution of irradiation bubbles generated by polycrystalline tungsten with a grain size of 181nm, a grain size of 45nm, and a grain size of 32nm, respectively.

FIG. 5 is a plot of porosity of polycrystalline tungsten as a function of grain size.

FIG. 6 is a graph of irradiated bubble density of polycrystalline tungsten as a function of grain size.

Detailed Description

The invention will be described in further detail below with reference to the accompanying drawings and specific embodiments.

The main concept of the invention is as follows:

as shown in fig. 1, obtaining physical parameters of a nuclear material required in an evolution process of irradiation bubbles, solving elastic energy in a coupling manner and constructing a free energy equation; establishing a phase field model of irradiation bubble morphology evolution on the basis of introducing a polycrystalline structure and considering temperature gradient, interface energy anisotropy and diffusion coefficient anisotropy, and establishing a phase field evolution equation; establishing and solving a rate theoretical model of irradiation bubble nucleation to obtain initial bubble distribution, and using the initial bubble distribution as an initial condition of phase field simulation; solving a phase field evolution equation by using a semi-implicit Fourier spectrum method to obtain the influence rule of the crystal grain size, the external stress and the anisotropic diffusion coefficient factors on the irradiation bubble nucleation and growth process and the thermal performance. The irradiation swelling and heat conduction rule of the nuclear material is obtained by counting the parameters such as the size, the density and the porosity of the irradiation bubbles and combining the visualization result, so that the service performance of the nuclear material is predicted.

In order to make the technical solution and advantages of the present invention more clear, the present invention is further described in detail below with reference to the accompanying drawings and examples thereof. It should be understood that the specific embodiments described herein are merely illustrative of the invention and are not limiting of the invention.

(1) Selecting phase field variables and collecting physical parameters of the nuclear material

Selecting a phase field variable for describing the evolution process of the nuclear material irradiation bubble morphology, wherein the phase field variable for expressing the point defect concentration is required to be selected: concentration of gas atoms c_gVacancy concentration c_vAnd interstitial atom concentration c_i(ii) a Simultaneous selection of bubble phase sequence parametersη, for distinguishing between the bubble phase and the matrix phase: η is 1.0 when in the bubble phase and 0.0 when in the matrix phase; the parameter phi of multiple crystal sequences needs to be selected_iDenotes a grain structure in which the index i denotes the ith grain and within the ith grain is phi_i＝1.0；

Collecting physical parameters of the nuclear material, including formation energy, diffusion coefficient, mobility and gradient coefficient of gas atoms, vacancies and interstitial atoms, and mobility of bubbles and grain interfaces;

from the diffusion coefficients of the gas atoms, vacancies and interstitial atoms of S1.2, the mobility of the point defects taking into account the diffusion anisotropy, the vacancy mobility M, is calculated_vThe following expression is given:

wherein D_vVacancy diffusion coefficient in tensor form, V_mIs the unit volume of the crystal lattice, R is the ideal gas constant, and T is the absolute temperature;

wherein D_v ⁰In order to have an effective diffusion coefficient for the vacancies,

and

the components of vacancy diffusion, surface diffusion and grain boundary diffusion are respectively as follows:

wherein w_b、w_sAnd w_GBWeight coefficients of bulk diffusion, surface diffusion and grain boundary diffusion, I is unit matrix, and T is_s ^vAnd T_GB ^vThe tensor coefficients of the surface and the grain boundary are respectively expressed as the following expressions:

(2) solving elastic balance equation to obtain elastic energy density

Calculating a position-dependent elastic modulus distribution: in polycrystalline materials, the modulus of elasticity is related to the orientation of the grains, the modulus of elasticity of each oriented grain being determined by a coordinate transformation, the modulus of elasticity tensor C_ijkl(r) the expression is as follows:

wherein the content of the first and second substances,

is the modulus of elasticity in the matrix and,

and

respectively, the rotation matrix a of the p-th crystal grain to the reference coordinate_ijIa, jb, kc, ld components in (1), rotation matrix a_ijThe expression is determined by the euler angle of the corresponding die,the three-dimensional coordinate lower spin torque matrix is as follows:

in the formula, theta, xi and psi are three Euler angles for determining the orientation of crystal grains;

initial intrinsic strain

Consisting of two parts, the first part being the intrinsic strain of the matrix phase

The second part being the intrinsic strain of the bubble phase

r is a position vector;

in which the intrinsic strain of the matrix phase is caused by defect inhomogeneities

Comprises the following steps:

wherein epsilon^g0Is the expansion coefficient of lattice constant, delta, due to the introduction of gas atoms_ijIs a function of the Kronecker function,

is the equilibrium concentration of gas atoms, h (η) ═ η³(6η²-15 η +10) is an interpolation function constructed such that when η is 1.0, h (η) is 1.0 in the bubble phase; when η is 0.0, h (η) is 0.0 in the matrix phase; the intrinsic strain of the bubble is very complex because of itDepending on the internal pressure of the bubble and the size of the bubble, when the internal pressure of the bubble is greater than 2 gamma_sAt/r, the outward pressure acts on the lattice around the bubble, so the characteristic strain of the bubble is positive, which means that the bubble is in an overpressure condition; when the bubble pressure is equal to 2 gamma_sAt/r, the intrinsic strain is zero; a negative intrinsic strain value indicates that the bubble is in an under-pressure state, which causes elastic relaxation of the crystal lattice around the bubble towards the center of the bubble; gamma ray_sIs the surface energy of the bubble, and r is the radius of the bubble;

intrinsic strain of bubble phase

The description is as follows:

wherein, P_gIs the internal pressure of the bubbles, 2 gamma_sConsidering/r as external hydrostatic pressure, C₁₁The gas pressure inside the bubble is always positive, being the modulus of elasticity in the main direction, while the intrinsic strain inside the bubble may be negative; this definition indicates that intrinsic strain is a function of gas pressure, bubble radius, surface tension and material elastic properties; the gas pressure in the bubbles can be obtained using the van der waals equation:

where Ω is the molecular volume of the corresponding core material, b is the Van der Waals constant, k_BBoltzmann constant;

stress-strain distribution is obtained by solving an elastic equilibrium equation which is as follows:

in the formula sigma_ijIs stress, r_jIs the unit displacement in the j direction; epsilon_kl(r) is the total strain equivalent to ε_ij(r)；ε_ij(r) is the sum of the uniform and non-uniform strains:

in the formula

And δ ε_ij(r) uniform strain and non-uniform strain, respectively; uniform strain

The expression is as follows:

in the formula (I), the compound is shown in the specification,

for additional stress, V is the simulation zone volume; non-uniform strain delta epsilon_ij(r) the expression is as follows:

in the formula u_i(r) and u_j(r) is the displacement component in the i and j directions, and the elastic balance equation (14) is obtained by substituting equation (17):

obtaining a displacement field equation after simplifying and sorting:

in the formula

Is a position-independent elastic modulus, C'_ijkl(r) is the position-dependent modulus of elasticity; the following is obtained after the Fourier transform of the above formula:

in the formula

Is the displacement field of Fourier space, k is the wave vector in Fourier space, k_jIs a component of a wave vector;

is the Green tensor; the subscript k indicates that this partial computation is performed in fourier space; equation (20) is solved using perturbation iteration, and the displacement field in the nth iteration is as follows:

the initial value of the iteration is taken as

Displacement field u in real space_k(r) is prepared from

Performing inverse Fourier transform, obtaining a displacement field by iteratively solving an elastic equilibrium equation, and further calculating to obtain stress-strain distribution;

calculating an elastic energy density from the obtained stress strain, the elastic energy density having the following expression:

(3) constructing the free energy equation

Based on thermodynamic theoretical derivation, calculating a free energy density function of a matrix phase of the nuclear material according to the following formula:

wherein the condition c is always satisfied when deriving the free energy of the matrix, since it is assumed that the defects occupy only the perfect lattice_g+c_v+c_i+c_m＝1.0；c_mIs the concentration of a perfect lattice, in at.%; e_g ^f、E_v ^f、E_i ^fRespectively the formation energy of gas atoms, vacancies and interstitial atoms in the matrix phase;

the free energy density function of the core material bubble phase is as follows:

wherein A and B are constants,

is a function of the density of free energy within the bubble in relation to the concentration of gas atoms; the Gibbs free energy expression G (p) of the bubbles is:

wherein p is pressure; v_bIs the volume of the bubble; g (p)₀) Is at a reference state pressure p₀Lower gibbs free energy; the van der waals equation at high pressure is:

p(V_b-nb)＝nk_BT (26)

wherein n is the number of gas atoms; the gibbs free energy from which the non-ideal gas is obtained is:

the gas atoms in the bubble phase may occupy lattice or non-lattice sites, so the concentration of gas atoms in the bubble is defined as:

wherein V_siteThe volume of the individual crystal lattice, from which the pressure:

the free energy density function of the bubble that brings the above results into account is:

wherein, mu₀Is a reference state p₀Chemical potential of (a);

introducing an interpolation function and a double-trap potential function to obtain a bulk free energy density function according to the obtained free energy density functions of the matrix phase and the bubble phase:

f_bulk(c_g,c_v,c_i,η,T)＝h(η)f_b(c_g,c_v,c_i,T)+[1-h(η)]f_m(c_g,c_v,c_i,T)+wg(η) (31)

wherein h (η) ═ η³(6η²-15 η +10) is an interpolation function constructed such that when η is 1.0, h (η) is 1.0 in the bubble phase; when η is 0.0, h (η) is 0.0 in the matrix phase; g (η) ═ η²(1-η)²Is a double potential well function, w is the potential well height;

the polycrystalline interaction free energy density function is expressed as follows:

where α, β, γ and δ are image-only parameters.

The total free energy equation is obtained by combining a free energy density function, an elastic energy density, an interface gradient energy, a bubble and grain boundary interaction energy and a polycrystal interaction energy:

in the formula, κ_g、κ_v、κ_i、κ_ηThe gradient term coefficients of gas atom concentration, vacancy concentration, interstitial atom concentration and bubble phase sequence parameter are respectively, and the last four terms of the formula represent gradient energy.

(4) Obtaining initial bubble size distribution by combining rate theory

Establishing a speed theoretical model of irradiation bubble nucleation growth behavior: the irradiation of nuclear material generates a large number of gas atoms which tend to accumulate in low density areas in order to reduce the energy of the nuclear material and minimize strain, a process known as nucleation; the change in free energy during nucleation is caused by a spherical-radius core, and this free energy Δ F is represented by the following equation:

wherein Δ f is the difference between the free energy of the bubble and the free energy of the substrate, r is the radius of the bubble, and σ is the surface energy; while a critical kernel corresponds to the maximum of this free energy, since any larger kernel will have less energy and will grow spontaneously; thereby obtaining the critical radius R of the nucleus^CAnd activation energy Δ F^C：

The bubble core is randomly introduced in the nucleation stage, and the average nucleation rate is required to be ensured to be matched with the expected nucleation rate; assuming that all points in the matrix can become nucleation sites, the probability of forming a core at one atomic position in one nucleation time interval can be calculated by considering the nucleation rate:

wherein Z is a Zeldovich correction factor, N is the atomic number, beta^*Is the probability of the critical core changing into the supercritical core, and tau is the inoculation time of the core; t is t₀Time intervals for calculating nucleation; the total nucleation probability P is thus obtained:

P＝1-exp(-J^*t₀) (38)

solving a formula (34) and a formula (38) to obtain nucleation probabilities corresponding to bubbles with different sizes, namely the size distribution of the bubbles at the initial stage of irradiation bubble nucleation; in the rate theory nucleation algorithm, the bubble radius and the nucleation probability of each nucleation point are calculated, and the bubble size distribution at the initial stage of irradiation bubble nucleation is used as an initial condition to be input into the phase field simulation.

(5) Establishing and solving phase field evolution equation

Considering the evolution equation of the temperature gradient, the generation of the point defect, the compounding of the point defect and the disappearance process of the point defect at the defect trap under the irradiation condition as follows:

evolution equation of gas atomic concentration:

evolution equation of vacancy atom concentration:

evolution equation of interstitial concentration:

evolution equation of bubble phase field variable:

evolution equation of the polycrystal phase field variable:

in the formula, t is simulation time; f is total free energy; l is the free interface mobility; m_g、M_v、M_iAtomic mobilities of gas atoms, vacancy atoms, and interstitial atoms, respectively;

thermal fluctuation terms of gas atoms, vacancy atoms, interstitial atoms and bubble phases respectively;

respectively the generation rates of gas atoms and point defects under irradiation conditions;

is the recombination rate of point defect vacancies and interstitial atoms,

is the absorption term of point defects by grain boundaries. The influence of the temperature gradient on bubble migration and grain growth is considered in the evolution equation, wherein Q is heat transfer efficiency, and C is a constant;

and

the production rates of gas atoms, vacancies and interstitial gas atoms under irradiation conditions respectively; the generation rate of gas atoms is:

in the formula (f)_rIs the rate of fission; ran is a random number between 0 and 1; Λ is a constant;

the generation of interstitial atoms and vacancies uses the following expression:

wherein R is₁And R₂Two random numbers uniformly distributed between 0 and 1 randomly generated at each time step and each space lattice point; e is a biased constant that can be varied to represent unequal numbers of vacancies and interstitials created by an off-site damage; parameter P_cascRepresenting the probability of a cascade collision per unit time and unit volume, V_GRepresents the maximum amount of increase in vacancies due to a cascading collision event; in the formula eta<The condition of 0.8 ensures that the cascade collisions occur only in the matrix phase and not in the bubble phase;

when the point defect vacancies and interstitial atoms meet, they recombine to form a perfect lattice, expressed as the recombination rate:

in the formula, v_rThe recombination rate is: v is_r＝v_b+η²v_s，v_bAnd v_sThe rate of recombination of point defects at the bulk and interface, respectively；v_b＝4πr_iv(D_i+D_v)/Ω，r_ivIs the composite volume radius; d_iAnd D_vThe diffusion coefficients of interstitial atoms and vacancies, respectively; Ω is the lattice volume;

for the absorption term of the point defects by the grain boundary, the expression is as follows:

wherein

The grain boundary absorption factor represents the strength of absorption of point defects by grain boundaries; phi ═ Σ phi_i ²As a function of the position of the crystal grains, phi is 1.0 in the crystal grains, and phi is less than 1.0 at the crystal boundary,

is the equilibrium concentration of gas atoms, vacancies, interstitial atoms;

solving the evolution equation of the phase field by using a semi-implicit Fourier spectrum method, and solving the evolution equation of the gas atom concentration, the vacancy concentration and the interstitial atom concentration as shown in the following formula:

solving the evolution equation of the bubble phase field variable and the polycrystal phase field variable as shown in the following formula:

wherein the superscript n represents the value of the part in the nth time step, the free energy f is the sum of the bulk free energy density, the polycrystal free energy density and the elastic free energy density, and k is (k)₁,k₂) Is the vector coordinate of the Fourier space, and delta t is the simulated time step;

(6) visualization process and simulation result statistics

And (3) carrying out visualization processing on the numerical solution result of the phase field evolution equation by using visualization software to obtain the appearance of bubble nucleation growth and grain growth in the nuclear material irradiation process:

representing visual variables of an irradiation bubble evolution simulation process:

in the formula (I), the compound is shown in the specification,

is a defined visual variable; the characteristics of the bonding diffusion interface are known: inside the crystal grain

At grain boundaries

In the air bubble

Therefore, the crystal grains, the crystal boundaries and the bubbles are distinguished by distinguishing the visual variable values of each point in the space;

and (3) counting the bubble size distribution in the output result by using an algorithm of the link region mark, wherein the method comprises the following steps:

in the connected domain marking, scanning from left to right and from top to bottom during the first marking, each effective pixel is set to be a label value, and the judgment rule is as follows:

(1) when the left adjacent pixel and the upper adjacent pixel of the pixel are invalid values, setting a new label value, label +1, for the pixel;

(2) when the left-adjacent pixel or the upper-adjacent pixel of the pixel has a valid value, assigning label of the valid value pixel to the label value of the pixel;

(3) when the left adjacent pixel and the upper adjacent pixel of the pixel are both effective values, selecting the smaller label value to be assigned to the label value of the pixel;

finally, counting the area of the region contained in each label value to obtain the size distribution and porosity parameters of the bubbles;

counting intercrystalline bubbles, and realizing the counting by using the following algorithm:

(1) numbering each bubble on the basis of S6.2 statistics, and recording the position of each bubble;

(2) judging whether a certain bubble is an intergranular bubble, wherein the method comprises the following steps: taking the range larger than the bubble as a calculation region, and calculating a multi-order parameter phi_iThe sum Γ within this range_iFinally, multiply it to obtain

If lambda is 0.0, the bubble is in the grain; if λ ≠ 0.0 then the bubble is at the grain boundary;

(3) sequentially judging whether all bubbles are intercrystalline bubbles, respectively obtaining intercrystalline bubble distribution and intragranular bubble distribution, and obtaining intercrystalline bubble coverage rate by using the area ratio of the intercrystalline bubbles to the total area of a crystal boundary;

(7) calculating effective thermal conductivity and radiation swelling

Calculating the effective thermal conductivity of the nuclear material requires solving the heat transfer equation:

in the formula k_localFor local thermal conductivity, k is a function of the heat transfer properties in spatial position, when located inside a grain_local＝k_bulk，k_bulkIs the thermal conductivity within the grains; when located on grain boundaries k_local＝k_GB，k_GBIs the thermal conductivity at grain boundaries; when located in the bubble k_local＝k_bubble，k_bubbleIs the thermal conductivity within the bubble; solving a heat transfer equation (56) by using a finite difference method to obtain an equilibrium temperature field and effective thermal conductivity;

obtaining the irradiation swelling rule of the nuclear material irradiation bubbles according to the size distribution and the porosity of the bubbles, wherein the fuel swelling caused by the irradiation bubbles in the nuclear material is estimated by the following formula:

where Δ V is the change in fuel volume, V₀Is the initial fuel volume, V_fThe fuel swelling was calculated as the ratio of the bubble area to the total area in the two-dimensional simulation for the final fuel volume including the bubble-induced swelling. The irradiation swelling and heat conduction rule of the nuclear material are obtained through calculation, so that the service performance of the nuclear material is predicted.

Some specific examples are provided below.

Example one

In order to make the objects, technical solutions and advantages of the present invention more apparent, the present invention is described in further detail below with reference to the accompanying drawings and embodiments.

The invention provides a simulation method for an evolution process of nuclear material irradiation bubbles, according to a flow in a specific embodiment, a phase field model is utilized to carry out numerical simulation on the evolution process of the irradiation bubbles of pure tungsten at a temperature of 1273K, and simulation use parameters are shown in a table 1:

TABLE 1 pure tungsten simulation parameters of use

Firstly, elastic energy is solved in a coupling mode, the elastic energy density is obtained by solving a mechanical equilibrium equation, and the expression is as follows:

wherein C is_ijkl(r) is the elastic constant tensor,. epsilon_ij ⁰Is the initial intrinsic strain, which is used to denote the inelastic strained portion,

to uniform strain, δ ε_ij(r) is non-uniform strain.

And then, constructing a free energy equation, deducing based on a thermodynamic theory, and calculating a free energy density function of the matrix phase of the nuclear material according to the following formula:

wherein the condition c is always satisfied when deriving the free energy of the matrix, since it is assumed that the defects occupy only the perfect lattice_g+c_v+c_i+c_m＝1.0；c_mIs the concentration of a perfect lattice, in at.%; k is a radical of_BBoltzmann constant; t is the absolute temperature; e_g ^f、E_v ^f、E_i ^fRespectively the formation energy of gas atoms, vacancies and interstitial atoms in the matrix phase;

the free energy density function of the bubble phase of the core material is as follows:

wherein A and B are constants,

the function of the density of free energy in the gas bubble in relation to the concentration of gas atoms is calculated by the van der Waals equation under high pressure

Wherein V_siteVolume of a single lattice, μ₀And p₀B is the van der waals constant for chemical potential and pressure in the reference state;

introducing an interpolation function and a double-trap potential function to obtain a bulk free energy density function according to the obtained free energy density functions of the matrix phase and the bubble phase:

f_bulk(c_g,c_v,c_i,η,T)＝h(η)f_b(c_g,c_v,c_i,T)+[1-h(η)]f_m(c_g,c_v,c_i,T)+wg(η) (31)

wherein h (η) ═ η³(6η²-15 η +10) is an interpolation function constructed such that when η is 1.0, h (η) is 1.0 in the bubble phase; when η is 0.0, h (η) is 0.0 in the matrix phase; g (η) ═ η²(1-η)²Is a double potential well function, w is the potential well height;

the polycrystalline interaction free energy density function is expressed as follows:

where α, β, γ and δ are image-only parameters.

The total free energy equation is obtained by combining a free energy density function, an elastic energy density, an interface gradient energy, a bubble and grain boundary interaction energy and a polycrystal interaction energy:

in the formula, κ_g、κ_v、κ_i、κ_ηGradient term coefficients of gas atom concentration, vacancy concentration, interstitial atom concentration and bubble phase sequence parameters are respectively, and the last four terms of the formula represent gradient energy;

constructing a phase field evolution equation, wherein the evolution equation of the temperature gradient, the generation of point defects, the compounding of the point defects and the disappearance process of the point defects at the defect trap under the irradiation condition is considered as follows:

evolution equation of gas atomic concentration:

evolution equation of vacancy atom concentration:

evolution equation of interstitial concentration:

evolution equation of the bubble phase field variable:

evolution equation of the polycrystal phase field variable:

in the formula, t is simulation time; f is total free energy; l is the free interface mobility; m_g、M_v、M_iAtomic mobilities of gas atoms, vacancy atoms, and interstitial atoms, respectively;

thermal fluctuation terms of gas atoms, vacancy atoms, interstitial atoms and bubble phases respectively;

respectively the generation rates of gas atoms and point defects under irradiation conditions;

is the recombination rate of point defect vacancies and interstitial atoms,

is the absorption term of point defects by grain boundaries. The influence of the temperature gradient on bubble migration and grain growth is considered in the evolution equation, wherein Q is heat transfer efficiency, and C is a constant;

solving a phase field evolution equation by a semi-implicit Fourier spectrum method, and carrying out visualization processing on numerical results by using Paraview, wherein the stress distribution of a single bubble in a simulation result is shown in fig. 2a, 2b and 2 c: FIG. 2a shows the stress component σ_xxDistribution, FIG. 2b is the stress component σ_yyDistribution, FIG. 2c is the stress component σ_xyAnd (4) distribution. The elastic field distribution generated by a single bubble is obtained through the simulation calculation of the first embodiment, and the correctness of the solution of the coupling elastic energy is ensured.

Example two

In the example, phase field simulation research is performed on the evolution of bubbles under the irradiation condition of the single crystal tungsten on the basis of the first embodiment. Initial bubble distribution was obtained with rate theory: establishing a speed theoretical model of irradiation bubble nucleation growth behavior: the nuclear material generates a large amount of gas atoms during irradiation, and in order to reduce the energy of the system and minimize the strain, the gas atoms tend to gather in a low-density area, and the process of forming bubbles is called a nucleation process; the change in free energy during nucleation is caused by a spherical-radius core, and this free energy Δ F is represented by the following equation:

wherein Δ f is the difference between the free energies of the bubble and the substrate, r is the bubble radius, and σ is the surface energy; while a critical kernel corresponds to the maximum of this free energy, since any larger kernel will have less energy and will grow spontaneously; thereby obtaining the critical radius R of the nucleus^CAnd activation energy Δ F^C：

The bubble core is randomly introduced in the nucleation stage, and the average nucleation rate is required to be ensured to be matched with the expected nucleation rate; assuming that all points in the matrix can become nucleation sites, the probability of forming a core at one atomic position in one nucleation time interval can be calculated by considering the nucleation rate:

wherein Z is a Zeldovich correction factor, N is the atomic number, beta^*Is the probability of the critical core changing into the supercritical core, and tau is the inoculation time of the core; t is t₀Time intervals for calculating nucleation; the total nucleation probability P is thus obtained:

P＝1-exp(-J^*t₀) (38)

solving a formula (34) and a formula (38) to obtain nucleation probabilities corresponding to bubbles with different sizes, namely the size distribution of the bubbles at the initial stage of irradiation bubble nucleation; in a rate theory nucleation algorithm, the bubble radius and the nucleation probability of each nucleation point are required to be calculated, and the bubble size distribution at the initial stage of irradiation bubble nucleation is used as an initial condition to be input into phase field simulation;

constructing a phase field evolution equation, wherein the evolution equation of the temperature gradient, the generation of point defects, the compounding of the point defects and the disappearance process of the point defects at the defect trap under the irradiation condition is considered as follows:

evolution equation of gas atomic concentration:

evolution equation of vacancy atom concentration:

evolution equation of interstitial concentration:

evolution equation of the bubble phase field variable:

evolution equation of the polycrystal phase field variable:

wherein L is the free interface mobility; m_g、M_v、M_iRespectively being gas atoms, vacancy atoms and interstitial atomsAtomic mobility of (a);

thermal fluctuation terms of gas atoms, vacancy atoms, interstitial atoms and bubble phases respectively;

respectively the generation rates of gas atoms and point defects under irradiation conditions;

is the recombination rate of point defect vacancies and interstitial atoms,

is the absorption term of point defects by grain boundaries. The influence of the temperature gradient on bubble migration and grain growth is considered in the evolution equation, wherein Q is heat transfer efficiency, and C is a constant;

and

the production rates of gas atoms, vacancies and interstitial gas atoms under irradiation conditions respectively; the generation rate of gas atoms is:

in the formula (f)_rIs the rate of fission; ran is a random number between 0 and 1; Λ is a constant;

the generation of interstitial atoms and vacancies uses the following expression:

wherein R is₁And R₂Is random at each time step and each space lattice pointTwo random numbers generated that are evenly distributed between 0 and 1; e is a biased constant that can be varied to represent unequal numbers of vacancies and interstitials created by an ex-situ damage; parameter P_cascRepresenting the probability of a cascade collision per unit time and unit volume, V_GRepresents the maximum amount of increase in vacancies due to a cascading collision event; in the formula eta<The condition of 0.8 ensures that the cascade collisions occur only in the matrix phase and not in the bubble phase;

solving a phase field evolution equation by a semi-implicit Fourier spectrum method, performing visualization processing on numerical results by using Paraview, and comparing the average bubble diameter of different irradiation intensities with experimental results as shown in FIG. 3, wherein the average bubble diameter is increased along with the increase of the irradiation intensity. The calculation results of the phase field simulation are basically consistent with the experimental results.

EXAMPLE III

In this embodiment, phase field simulation research is performed on the evolution of bubbles under the irradiation condition of polycrystalline tungsten on the basis of the second embodiment, polycrystalline structures with different sizes are introduced in the second embodiment, and the process evolution equations considering the generation of point defects, the recombination of point defects, the disappearance of point defects at a defect trap, and the like, are obtained as follows:

gas atomic concentration field evolution equation:

evolution equation of vacancy atom concentration field:

evolution equation of interstitial concentration field:

evolution equation of bubble phase sequence parameters:

evolution equation of the polycrystalline order parameter:

wherein L is the free interface mobility; m_g、M_v、M_iAtomic mobilities of gas atoms, vacancy atoms, and interstitial atoms, respectively;

thermal fluctuation terms of gas atoms, vacancy atoms, interstitial atoms and bubble phases respectively;

respectively the generation rates of gas atoms and point defects under irradiation conditions;

is the recombination rate of point defect vacancies and interstitial atoms.

For the absorption term of the point defects by the grain boundary, the expression is as follows:

wherein

The grain boundary absorption factor indicates the strength of absorption of point defects by grain boundaries. Phi ═ Σ phi_i ²As a function of the position of the crystal grains, phi is 1.0 in the crystal grains, and phi is less than 1.0 at the crystal grain boundary.

Solving an evolution equation by a Fourier spectrum method, and performing visualization processing on numerical results by using Paraview, wherein the distribution of irradiation bubbles generated by polycrystalline tungsten with different grain sizes is shown in FIG. 4a, FIG. 4b and FIG. 4 c: both bubble size and density decrease as the grain size becomes smaller. This is because the smaller size grains have larger grain boundary area, and the grain boundary can absorb a large amount of surrounding point defects, so that the small size grains have smaller radiation swelling, which also provides a possible idea for improving the radiation resistance of pure tungsten. The average bubble size for different grain sizes is shown in fig. 5, the bubble density for different grain sizes is shown in fig. 6, and the average bubble size and bubble density variations for different grain sizes are substantially consistent with experimental observations. Through the simulation research of the embodiment, the influence of the grain size on the evolution of the irradiation bubbles is quantitatively analyzed, and a theoretical explanation is provided for experimental phenomena.

Claims (8)

1. A method for determining the appearance of irradiation bubbles in a nuclear material and the influence of the irradiation bubbles on the mechanical and thermal properties is characterized in that: the method comprises the following steps:

s1, selecting a phase field variable describing the nuclear material irradiation bubble shape evolution process; collecting physical parameters of the nuclear material, including the formation energy, diffusion coefficient and gradient coefficient of point defects, and the mobility of bubbles and grain boundaries; and calculating mobility of the point defect considering the diffusion anisotropy;

s2, calculating the elastic modulus distribution related to the position; combining Van der Waals equation to obtain initial intrinsic strain distribution; establishing and solving an elastic balance equation to obtain irradiation bubble stress-strain distribution and elastic energy density;

s3, obtaining a free energy density function of a matrix phase of the nuclear material based on a material thermodynamic theory according to the physical parameters of the nuclear material collected in the step S1, obtaining a free energy density function of a bubble phase in the nuclear material by deducing a Van der Waals equation under the conditions of high temperature and high pressure, and obtaining a total free energy equation by combining interface gradient energy, bubble-crystal boundary interaction energy, polycrystal interaction energy and elastic energy density obtained in the step S2;

s4, establishing a speed theoretical model of irradiation bubble nucleation; obtaining bubble size distribution at the initial stage of irradiation bubble nucleation by solving a rate theoretical model, and taking the size distribution as an initial condition of phase field simulation;

s5, constructing a phase field evolution equation according to the total free energy equation obtained in the step S3 and comprehensively considering factors of temperature gradient, generation of point defects, compounding of the point defects and disappearance of the point defects at a defect trap under an irradiation condition; carrying out numerical solution on a phase field evolution equation by adopting a semi-implicit Fourier spectrum method;

s6, carrying out visualization processing on the numerical solution obtained in the step S5 to obtain the appearance of bubble nucleation growth and grain growth in the nuclear material irradiation process; counting all bubbles by adopting a communicated region marking method to obtain the size distribution and porosity parameters of the irradiated bubbles in the nuclear material; carrying out classified statistics on the intercrystalline bubbles and the intragranular bubbles;

s7, calculating and obtaining the effective thermal conductivity and the radiation swelling rule of the nuclear material based on the bubble distribution and the porosity obtained in the step S6.

2. The method of determining the morphology of irradiated bubbles in nuclear material and their effect on the thermodynamic and thermal properties of nuclear material as claimed in claim 1 wherein: the step S1 includes the following sub-steps:

s1.1, selecting phase field variables for describing the evolution process of the nuclear material irradiation bubble morphology, wherein the phase field variables for expressing the point defect concentration need to be selected: concentration of gas atoms c_gVacancy concentration c_vAnd interstitial atom concentration c_i(ii) a And simultaneously selecting a bubble phase sequence parameter eta for distinguishing a bubble phase from a matrix phase: η is 1.0 when in the bubble phase and 0.0 when in the matrix phase; the parameter phi of the multiple crystal order needs to be selected_iDenotes a grain structure in which the index i denotes the ith grain and within the ith grain is phi_i＝1.0；

S1.2, collecting physical parameters of the nuclear material, including formation energy, diffusion coefficient, mobility and gradient coefficient of gas atoms, vacancies and interstitial atoms, and mobility of bubbles and grain interfaces;

s1.3 calculating the mobility of the point defect taking into account the diffusion anisotropy, the vacancy mobility M, from the diffusion coefficients of the gas atoms, vacancies and interstitial atoms of S1.2_vThe following expression is given:

wherein D_vVacancy diffusion coefficient in tensor form, V_mIs the unit volume of the crystal lattice, R is the ideal gas constant, and T is the absolute temperature;

wherein D_v ⁰In order to have an effective diffusion coefficient for the vacancies,

and

the components of vacancy diffusion, surface diffusion and grain boundary diffusion are respectively as follows:

wherein w_b、w_sAnd w_GBWeight coefficients of bulk diffusion, surface diffusion and grain boundary diffusion, I is unit matrix, and T is_s ^vAnd T_GB ^vThe tensor coefficients of the surface and the grain boundary are respectively expressed as the following expressions:

3. the method of determining the morphology of irradiated bubbles in nuclear material and their effect on the thermodynamic and thermal properties of nuclear material as claimed in claim 1 wherein: the step S2 includes the following sub-steps:

s2.1 calculate the position-dependent elastic modulus distribution: in polycrystalline materials, the modulus of elasticity is related to the orientation of the grains, the modulus of elasticity of each oriented grain being determined by coordinate transformation, the tensor C of the modulus of elasticity_ijkl(r) the expression is as follows:

wherein, the first and the second end of the pipe are connected with each other,

is the modulus of elasticity in the matrix and,

and

respectively, the rotation matrix a of the p-th crystal grain to the reference coordinate_ijIa, jb, kc,ld components, rotation matrix a_ijThe expression is determined by the Euler angle of the corresponding crystal grain, and the rotation torque matrix under the three-dimensional coordinate is as follows:

in the formula, theta, xi and psi are three Euler angles for determining the orientation of crystal grains;

s2.2 initial intrinsic Strain

Consisting of two parts, the first part being the intrinsic strain of the matrix phase

The second part being the intrinsic strain of the bubble phase

r is a position vector;

in which the intrinsic strain of the matrix phase is caused by defect inhomogeneities

Comprises the following steps:

wherein epsilon^g0Is the expansion coefficient of lattice constant, delta, due to the introduction of gas atoms_ijIs a function of the Kronecker function,

is the equilibrium concentration of gas atoms, h (η) ═ η³(6η²-15 η +10) is an interpolation function constructed such that when η is 1.0, h (η) is 1.0 in the bubble phase; when η is 0.0, h (η) is 0.0 in the matrix phase; the intrinsic strain of the bubble is very complex because it depends on the internal pressure of the bubble and the size of the bubble, when the internal pressure of the bubble is greater than 2 gamma_sAt/r, the outward pressure acts on the lattice around the bubble, so the characteristic strain of the bubble is positive, which means that the bubble is in an overpressure condition; when the bubble pressure is equal to 2 gamma_sAt/r, the intrinsic strain is zero; a negative intrinsic strain value indicates that the bubble is in an under-pressure state, which causes elastic relaxation of the crystal lattice around the bubble towards the center of the bubble; gamma ray_sIs the surface energy of the bubble, and r is the radius of the bubble;

intrinsic strain of bubble phase

The description is as follows:

wherein, P_gIs the internal pressure of the bubbles, 2 gamma_sConsidering/r as external hydrostatic pressure, C₁₁The gas pressure inside the bubble is always positive, being the modulus of elasticity in the main direction, while the intrinsic strain inside the bubble may be negative; this definition indicates that intrinsic strain is a function of gas pressure, bubble radius, surface tension and material elastic properties; the gas pressure in the bubbles can be obtained using the van der waals equation:

where Ω is the molecular volume of the corresponding core material, b is the Van der Waals constant, k_BBoltzmann constant;

s2.3, stress-strain distribution is obtained by solving an elastic equilibrium equation which is as follows:

in the formula sigma_ijIs stress, r_jIs the unit displacement in the j direction; epsilon_kl(r) is the total strain equivalent to ε_ij(r); total strain epsilon_ij(r) is the sum of the uniform and non-uniform strains:

in the formula

And δ ε_ij(r) uniform strain and non-uniform strain, respectively; uniform strain

The expression is as follows:

in the formula (I), the compound is shown in the specification,

for additional stress, V is the simulation zone volume; non-uniform strain delta epsilon_ij(r) the expression is as follows:

in the formula u_i(r) and u_j(r) is the displacement component in the i and j directions, and the elastic balance equation (14) is obtained by substituting equation (17):

obtaining a displacement field equation after simplifying and sorting:

in the formula

Is a position-independent elastic modulus, C'_ijkl(r) is a position-dependent modulus of elasticity; the following is obtained after the Fourier transform of the above formula:

in the formula

Is the displacement field of Fourier space, k is the wave vector in Fourier space, k_jIs a component of a wave vector;

is the green tensor; the subscript k indicates that this partial computation is performed in fourier space; equation (20) is solved using perturbation iteration, and the displacement field in the nth iteration is as follows:

the initial value of the iteration is taken as

Displacement field u in real space_k(r) is prepared from

Is obtained by inverse Fourier transformIteratively solving an elastic balance equation to obtain a displacement field, and further calculating to obtain stress-strain distribution;

calculating an elastic energy density from the obtained stress strain, the elastic energy density having the following expression:

4. the method of determining the morphology of irradiated bubbles in nuclear material and their effect on the thermodynamic and thermal properties of nuclear material as claimed in claim 1 wherein: the step S3 includes the following sub-steps:

s3.1, based on the theoretical derivation of thermodynamics, calculating a free energy density function of a matrix phase of the nuclear material according to the following formula:

wherein the condition c is always satisfied when deriving the free energy of the matrix, since it is assumed that the defects occupy only the perfect lattice_g+c_v+c_i+c_m＝1.0；c_mIs the concentration of a perfect lattice, in at.%;

respectively the formation energy of gas atoms, vacancies and interstitial atoms in the matrix phase;

s3.2 the free energy density function of the bubble phase of the nuclear material is as follows:

wherein A and B are constants,

is the free energy density in the bubble related to the gas atom concentrationA degree function; the Gibbs free energy expression G (p) of the bubbles is:

wherein p is pressure; v_bIs the volume of the bubble; g (p)₀) Is at a reference state pressure p₀Lower gibbs free energy; the van der waals equation at high pressure is:

p(V_b-nb)＝nk_BT (26)

wherein n is the number of gas atoms; the gibbs free energy from which the non-ideal gas is obtained is:

the gas atoms may occupy lattice or non-lattice sites in the bubble phase, so the concentration of gas atoms in the bubble is defined as:

wherein V_siteThe volume of the individual crystal lattice, from which the pressure:

the free energy density function of the bubble that brings the above results into account is:

wherein, mu₀Is a reference state p₀Chemical potential of (a);

s3.3, introducing an interpolation function and a double-trap potential function to obtain a bulk free energy density function according to the obtained free energy density functions of the matrix phase and the bubble phase:

f_bulk(c_g,c_v,c_i,η,T)＝h(η)f_b(c_g,c_v,c_i,T)+[1-h(η)]f_m(c_g,c_v,c_i,T)+wg(η)(31)

wherein h (η) ═ η³(6η²-15 η +10) is an interpolation function constructed such that when η is 1.0, h (η) is 1.0 in the bubble phase; when η is 0.0, h (η) is 0.0 in the matrix phase; g (η) ═ η²(1-η)²Is a double potential well function, w is the potential well height;

the polycrystalline interaction free energy density function is expressed as follows:

wherein α, β, γ and δ are image-only parameters;

the total free energy equation is obtained by combining a free energy density function, an elastic energy density, an interface gradient energy, a bubble and grain boundary interaction energy and a polycrystal interaction energy:

in the formula, κ_g、κ_v、κ_i、κ_ηThe gradient term coefficients of gas atom concentration, vacancy concentration, interstitial atom concentration and bubble phase sequence parameter are respectively, and the last four terms of the formula represent gradient energy.

5. The method of determining the morphology of irradiated bubbles in nuclear material and their effect on the thermodynamic and thermal properties of nuclear material as claimed in claim 1 wherein: the step S4 includes the following sub-steps:

s4.1, establishing a speed theoretical model of irradiation bubble nucleation growth behavior by combining the physical parameters of the nuclear material obtained in S1: the irradiation of nuclear material generates a large number of gas atoms which tend to accumulate in low density areas in order to reduce the energy of the nuclear material and minimize strain, a process known as nucleation; the change in free energy during nucleation is caused by a spherical-radius core, and this free energy Δ F is represented by the following equation:

wherein Δ f is the difference between the free energies of the bubble and the substrate, r is the bubble radius, and σ is the surface energy; while a critical kernel corresponds to the maximum of this free energy, since any larger kernel will have less energy and will grow spontaneously; thereby obtaining the critical radius R of the nucleus^CAnd activation energy Δ F^C：

Randomly introducing bubble cores in a nucleation stage, and ensuring that the average nucleation rate is matched with the expected nucleation rate; assuming that all points in the matrix can become nucleation sites, the probability of forming a core at one atomic position in one nucleation time interval can be calculated by considering the nucleation rate:

wherein Z is a Zeldovich correction factor, N is the atomic number, beta^*Is the probability of the critical core changing into the supercritical core, and tau is the inoculation time of the core; t is t₀Time intervals for calculating nucleation; the total nucleation probability P is thus obtained:

P＝1-exp(-J^*t₀) (38)

s4.2, solving a formula (34) -a formula (38) to obtain nucleation probabilities corresponding to bubbles with different sizes, namely, the size distribution of the bubbles at the initial stage of irradiation bubble nucleation; in the rate theory nucleation algorithm, the bubble radius and the nucleation probability of each nucleation point are calculated, and the bubble size distribution at the initial stage of irradiation bubble nucleation is used as an initial condition to be input into the phase field simulation.

6. The method of determining the morphology of irradiated bubbles in nuclear material and their effect on the thermodynamic and thermal properties of nuclear material as claimed in claim 1 wherein: the step S5 includes the following sub-steps:

s5.1 considering the temperature gradient under the irradiation condition, the generation of point defects, the compounding of the point defects and the disappearance process evolution equation of the point defects at the defect trap are as follows:

evolution equation of gas atomic concentration:

evolution equation of vacancy atom concentration:

evolution equation of interstitial concentration:

evolution equation of the bubble phase field variable:

evolution equation of the polycrystal phase field variable:

in the formula, t is simulation time; f is total free energy; l is the free interface mobility; m_g、M_v、M_iAtomic mobilities of gas atoms, vacancy atoms, and interstitial atoms, respectively;

thermal fluctuation terms of gas atoms, vacancy atoms, interstitial atoms and bubble phases respectively;

respectively the generation rates of gas atoms and point defects under irradiation conditions;

is the recombination rate of point defect vacancies and interstitial atoms,

is the absorption term of point defects by grain boundaries. The influence of the temperature gradient on bubble migration and grain growth is considered in the evolution equation, wherein Q is heat transfer efficiency, and C is a constant;

and

the production rates of gas atoms, vacancies and interstitial gas atoms under irradiation conditions respectively; the generation rate of gas atoms is:

in the formula (f)_rIs the rate of fission; ran is a random number between 0 and 1; Λ is a constant;

the generation of interstitial atoms and vacancies uses the following expression:

wherein R is₁And R₂Two random numbers uniformly distributed between 0 and 1 randomly generated at each time step and each space lattice point; e is a biased constant that can be varied to represent unequal numbers of vacancies and interstitials created by an ex-situ damage; parameter P_cascRepresenting the probability of a cascade collision per unit time and unit volume, V_GRepresents the maximum amount of increase in vacancies due to a cascading collision event; in the formula eta<The condition of 0.8 ensures that the cascade collisions occur only in the matrix phase and not in the bubble phase;

when the point defect vacancies and interstitial atoms meet, they recombine to form a perfect lattice, expressed as the recombination rate:

in the formula, v_rThe recombination rate is: v is_r＝v_b+η²v_s，v_bAnd v_sThe recombination rates of point defects at bulk and interface, respectively; v. of_b＝4πr_iv(D_i+D_v)/Ω，r_ivIs the composite volume radius; d_iAnd D_vThe diffusion coefficients of interstitial atoms and vacancies, respectively; Ω is the lattice volume;

for the absorption term of the point defects by the grain boundary, the expression is as follows:

wherein

The grain boundary absorption factor represents the strength of absorption of point defects by grain boundaries; phi ═ Σ phi_i ²As a function of the position of the crystal grains, phi is 1.0 in the crystal grains, and phi is less than 1.0 at the crystal boundary,

is the equilibrium concentration of gas atoms, vacancies, interstitial atoms;

s5.2, solving a phase field evolution equation by using a semi-implicit Fourier spectrum method, wherein the solution of the evolution equation of the gas atom concentration, the vacancy concentration and the interstitial atom concentration is shown as the following formula:

solving the evolution equation of the bubble phase field variable and the polycrystal phase field variable as shown in the following formula:

wherein the superscript n represents the value of the portion in the nth time step, the free energy f is the bulk free energy density, the sum of the polycrystalline free energy density and the elastic free energy density, and k is (k ═ k)₁,k₂) Is the vector coordinate of the fourier space, Δ t is the time step of the simulation.

7. The method of determining the morphology of irradiated bubbles in nuclear material and their effect on the thermodynamic and thermal properties of nuclear material as claimed in claim 1 wherein: the step S6 includes the following sub-steps:

s6.1, carrying out visualization processing on the numerical solution result in the S5 by using visualization software to obtain the appearance of bubble nucleation growth and grain growth in the nuclear material irradiation process:

representing visual variables of an irradiation bubble evolution simulation process:

in the formula (I), the compound is shown in the specification,

is a defined visual variable; the characteristics of the bonding diffusion interface are known: inside the crystal grain

At grain boundaries

In the air bubble

Therefore, the crystal grains, the crystal boundaries and the bubbles are distinguished by distinguishing the visual variable values of each point in the space;

s6.2, counting the size distribution of the bubbles in the output result by using an algorithm of the link region mark, wherein the method comprises the following steps:

in the connected domain marking, scanning from left to right and from top to bottom during the first marking, each effective pixel is set to be a label value, and the judgment rule is as follows:

(1) when the left adjacent pixel and the upper adjacent pixel of the pixel are invalid values, setting a new label value, label +1, for the pixel;

(2) when the left-adjacent pixel or the upper-adjacent pixel of the pixel has a valid value, assigning label of the valid value pixel to the label value of the pixel;

(3) when the left adjacent pixel and the upper adjacent pixel of the pixel are both effective values, selecting the smaller label value to be assigned to the label value of the pixel;

finally, counting the area of the region contained in each label value to obtain the size distribution and porosity parameters of the bubbles;

s6.3, counting the intercrystalline bubbles, and realizing by using the following algorithm:

(1) numbering each bubble on the basis of S6.2 statistics, and recording the position of each bubble;

(2) judging whether a certain bubble is an intergranular bubble, wherein the method comprises the following steps: taking the range larger than the bubble as a calculation region, and calculating a multi-order parameter phi_iThe sum Γ within this range_iFinally, multiply it to obtain

If lambda is 0.0, the bubble is in the grain; if λ ≠ 0.0 then the bubble is at the grain boundary;

(3) and sequentially judging whether all the bubbles are intercrystalline bubbles, respectively obtaining intercrystalline bubble distribution and intragranular bubble distribution, and obtaining the intercrystalline bubble coverage rate by using the area ratio of the intercrystalline bubbles to the total area of the crystal boundary.

8. The method of determining the morphology of irradiated bubbles in nuclear material and their effect on the thermodynamic and thermal properties of nuclear material as claimed in claim 1 wherein: the step S7 includes the following sub-steps:

s7.1 calculating the effective thermal conductivity of the nuclear material, the heat transfer equation needs to be solved:

▽·(k_local▽T)＝0 (56)

in the formula k_localFor local thermal conductivity, k is a function of the heat transfer properties in spatial position, when located inside a grain_local＝k_bulk，k_bulkIs the thermal conductivity within the grains; when located on grain boundaries k_local＝k_GB，k_GBIs the thermal conductivity at the grain boundaries; when located in the bubble k_local＝k_bubble，k_bubbleIs the thermal conductivity within the bubble; solving a heat transfer equation (56) by using a finite difference method to obtain an equilibrium temperature field and effective thermal conductivity;

s7.2, combining the size distribution and the porosity of the bubbles calculated in the step S6 to obtain an irradiation swelling rule of the nuclear material irradiated bubbles, wherein the fuel swelling caused by the irradiated bubbles in the nuclear material is estimated by the following formula:

where Δ V is the change in fuel volume, V₀Is the initial fuel volume, V_fThe fuel swelling was calculated in a two-dimensional simulation as the ratio of the bubble area to the total area for the final fuel volume including the bubble-induced swelling.

Images