Evaluation of the fraction of delayed photoneutrons for TMSR-SF1

NUCLEAR ENERGY SCIENCE AND ENGINEERING

Evaluation of the fraction of delayed photoneutrons for TMSR-SF1

Rui-Min Ji，

Ye Dai，

Gui-Feng Zhu，

Shi-He Yu，

Yang Zou，

Gui-Min Liu

Nuclear Science and Techniques

Vol.28, No.9

Article number 135

Published in print 01 Sep 2017

Available online 21 Aug 2017

DOI：10.1007/s41365-017-0285-9

36501

The 10 MW_th solid fueled thorium molten salt reactor (TMSR-SF1) is a FLiBe-salt-cooled pebble bed reactor to be deployed in 5–10 years, designed by the TMSR group. Due to a large amount of beryllium in the core, the photoneutrons are produced via (γ, n) reactions. Some of them are generated a long time after the fission event, and therefore are considered as delayed neutrons. In this paper, we redefine the effective delayed neutrons into two fractions: the delayed fission neutron fraction and the delayed photoneutron fraction. With some reasonable assumptions, the inner product method and the k-ratio method are adopted for studying the effective delayed photoneutron fraction. In the k-ratio method, the Monte Carlo code MCNP6 is used to evaluate the effective photoneutron fraction as the ratio between the multiplication factors with and without contribution of the delayed neutrons and photoneutrons. In the inner product method, with the Monte Carlo and deterministic codes together, we use the adjoint neutron flux as a weighting function for the neutrons and photoneutrons generated in the core. Results of the two methods agree well with each other, but the k-ratio method requires much more computing time for the same precision.

PhotoneutronEffective delayed neutron fractionTMSR-SF1.

1 Introduction

Fluoride-salt-cooled High-temperature Reactor (FHR) is a new reactor concept that has gained great attention worldwide [1, 2]. Different from the gas cooled reactors, the liquid salt enables operation with higher power density and provides higher heat capacity. With unique features, an FHR contains technologies derived from earlier reactor designs [3]. In recent years, the Thorium Molten Salt Reactor (TMSR) group of Chinese Academy of Sciences (CAS) proposed the 10 MW_th Solid Fueled TMSR (TMSR-SF1) as China’s first FHR to be deployed in 5–10 years [4-6]. In the design, TMSR-SF1 uses coated particle fuel and the beryllium-based salt FLiBe (⁷Li₂BeF₄) as its primary coolant.

Since there is a large amount of beryllium in the core, photoneutrons are generated via (γ, n) reaction. This reaction requests the threshold energy of 1.67 MeV for the photons. In general, the photoneutrons have little contribution to the overall neutron balance. However, like the delayed fission neutrons which are generated from a few milliseconds to a few minutes after the fission event, some photoneutrons are generated from a few minutes to a few days after the fission event. Thus the photoneutrons are considered as delayed neutrons. Since they have an insignificant contribution and the equilibrium state is difficult to reach, some researchers omitted the photoneutrons in experimental analysis. Therefore, there are measurement errors caused by neglecting the photoneutrons [7-12]. In the LITR reactor, which uses beryllium reflectors, it was found that the control rod worth would be underestimated for about 6% if neglecting the photoneutrons in the rod-drop measurement [10]. In VVRSZM type research reactor of KFKI-AEKI [8], which uses beryllium reflectors, the reactivity meter gave false values without taking into account the photoneutrons. In the miniature neutron source reactor using beryllium reflectors, the reactivity meter gave greater values by considering the photoneutrons. [11]

Considering that the photoneutrons may impact certain dynamic characteristics, in this paper we focus on evaluating the effect of photoneutrons on β_eff for the TMSR-SF1. This paper is organized as follows. Section 2 briefly introduces TMSR-SF1 reactor. Section 3 presents the description of the methodology and the codes used. Section 4 presents results and analyses. Finally, the conclusions are presented in Section 5.

2 TMSR-SF1 core configurations

The pebbles in the TMSR-SF1 are 6 cm in diameter. TRI-structural ISOtropic (TRISO) particles are made up of 17.0% enriched UO₂ kernels coated with a low density buffer layer, an inner pyrolytic carbon (iPyC) layer, an intermediate silicon carbide (SiC) layer and an outer pyrolytic carbon (oPyC) layer. Parameters of the fuel element are listed in Table 1. A schematic representation of the TMSR-SF1 reactor is shown in Fig. 1. The TMSR-SF1 configuration consists of an active zone, a reflector zone and the vessel. The active zone includes a cylinder and two circular truncated cones. The radius of the cylinder is sized at Φ 135 cm×180 cm. The circular truncated cones are sized at Φ 135 cm×30 cm. In the active zone, the fuel pebbles are randomly distributed and the coolant is circulating through the pebbles. Thickness of the upper and lower reflectors is 30 cm in minimum. The side reflector is 75 cm thick. In the reflector zone, the graphite blocks host the control rods, neutron detectors, temperature detectors, etc. To protect the control rods and detectors from FLiBe, Hastelloy-N alloy tubes are used in the channels. The nominal temperatures of pebble kernel, layers, matrix, shell, FLiBe coolant and reflector graphite are 1050, 1000, 980, 940, 900 and 880 K, respectively. The FLiBe density is 1.97 g/cm³ and Li-7 enrichment is 99.99 %. The reflector density is 1.825 g/cm³.

TMSR-SF1 fuel element characteristic parameters

Parameters	Values
Pebble diameter (cm)	6.0
Fuel zone diameter (cm)	5.0
U loading per pebble (g)	7.0
²³⁵U enrichment (%)	17
Equivalent boron content in fuel (ppm)	4
Graphite density (g/cm³)	1.73
Equivalent boron content in graphite (ppm)	3
Fuel kernel radius (mm)	0.25
UO₂ density (g/cm³)	10.4
Coating materials	Buffer/iPyC/SiC/oPyC
Layer thickness (mm)	0.095/0.040/0.035/0.040
Layer densities (g/cm³)	1.10/1.90/3.18/1.90

Fig. 1

(Color online) The TMSR-SF1 model

3 Methods

3.1. Theoretical analysis

The time-dependent transport equation for a system without external source can be written as [13]

L Φ = ω G Φ,

(1)

where, Φ is a vector including the neutron flux (φ) and photon flux (ψ), L is the loss operator, G is the generation operator, and ω=1/k is the inverse of the classical multiplication eigenvalue. It can be further written as follows

Φ = (\begin{matrix} ϕ \\ ψ \end{matrix}),

(2)

L = (\begin{matrix} M^{n} \\ 0 \end{matrix} \begin{matrix} 0 \\ M^{γ} + H^{ph} \end{matrix}),

(3)

G = (\begin{matrix} F_{^{p}}^{n} + D_{^{d}}^{n} & H^{ph} \\ F_{^{p}}^{γ} + D_{^{d}}^{γ} & 0 \end{matrix}),

(4)

where, M is the removal operator; and the superscripts n and γ denote neutronic and photonic operators or functions, respectively. So, Mⁿ is the neutron operator and M^γ is the photon removal operator. For a system containing Be, neutrons are generated from fissions (Fⁿ) and (γ, n) reactions (H^ph). The Fⁿ includes a prompt term (F_pⁿ) and a delayed term (D_dⁿ). H^ph, which describes the neutrons produced by Be(γ,n) reaction, includes a prompt term (H_p^ph) and a delayed term (H_d^ph) according to the photon source. The photons are generated directly after the fission events (F_p^γ) or from the decay of the precursors (D_d^γ). Definitions of the operators are listed in Table 2, which are fairly conventional but some terms deserve further clarifications. χⁿ_p,m and χⁿ_d,m are prompt and delayed neutron spectra for fissionable nuclide m, respectively. vⁿ_p,m and vⁿ_d,m are the prompt and delayed neutron yields per fission for nuclide m, respectively. Σ_f,m and σ_f,m are the macroscopic and microscopic fission cross sections for nuclide m, respectively. N_m(r) is the atomic density for nuclide m at position r. χ^γ_p,m is prompt photon spectra for fissionable nuclide m. v^γ_p,m is the prompt photon yield per fission for fissionable nuclide m. Σ_γ,m' and σ_γ,m are the macroscopic and microscopic capture cross sections for nuclide m', respectively. N_m'(r) is the atomic density for nuclide m' at position r. χⁿ_d,m'' is the delayed photon spectra for nuclide m''. λ^γ_d,m'' is the delay constant for nuclide m''. N_m''(r) is the atomic density for nuclide m'' at position r. Σ_γn and σ_γn are the macroscopic and microscopic Be(γ, n) cross section. N_Be(r) is the atomic density for nuclide Be at position r. μ_c is the Compton effect which alters the direction and energy of the incident photons. μ represents the removal of photons due to the photoelectric reactions and pair productions.

Notations of the operators in Eqs. (2 - 6)

Operator	Neutron	Photon
Prompt source	$F_{p}^{n} ϕ \equiv \sum_{m} χ_{p, m}^{n} \iint ν_{p, m}^{n} Σ_{f, m} (r, E^{'}) ϕ (r, t, Ω^{'}, E^{'}) d E^{'} d Ω^{'}$ $Σ_{f, m} (r, E^{'}) = σ_{f, m} (E^{'}) N_{m} (r)$	$\begin{array}{l} F_{p}^{γ} ϕ \equiv \sum_{m} χ_{p, m}^{γ} \iint ν_{p, m}^{γ} Σ_{f, m} (r, E^{'}) ϕ (r, t, Ω^{'}, E^{'}) d E^{'} d Ω^{'} \\ \begin{matrix} + \end{matrix} \sum_{m'} \iint Σ_{γ, m "} (r, E^{'}) ϕ (r, t, Ω^{'}, E^{'}) d E^{'} d Ω^{'} \end{array}$ $Σ_{γ, m'} (r, E^{'}) = σ_{γ, m'} (E^{'}) N_{m'} (r)$
Delayed source	$D_{d}^{n} ϕ \equiv \sum_{m} χ_{d, m}^{n} \iint ν_{d, m}^{n} Σ_{f, m} (r, E^{'}) ϕ (r, t, Ω^{'}, E^{'}) d E^{'} d Ω^{'}$	$D_{d}^{γ} ϕ \equiv \sum_{m^{″}} χ_{d, m^{″}}^{γ} λ_{d, m^{″}}^{γ} N_{m^{″}} (r)$
Removal	$\begin{array}{l} M^{n} ϕ \equiv Ω \cdot \nabla ϕ + Σ_{t} ϕ \\ - \iint Σ_{s} (r, E^{'} \to E, Ω^{'} \to Ω) ϕ (r, t, Ω^{'}, E^{'}) d E^{'} d Ω^{'} \end{array}$	$\begin{array}{l} M^{γ} ψ \equiv Ω \cdot \nabla ψ + μ ψ \\ - \iint μ_{c} (r, E^{'} \to E, Ω^{'} \to Ω) ψ (r, t, Ω^{'}, E^{'}) d E^{'} d Ω^{'} \end{array}$
Photoneutron source	$H^{ph} ψ \equiv χ^{γn} (E) \iint Σ_{γn} (r, E^{'}, Ω^{'}) ψ (r, t, Ω^{'}, E^{'}) d E^{'} d Ω^{'}$ ; $Σ_{γn} (r, E^{'}) = σ_{γn} (E^{'}) N_{Be} (r)$

The neutron sources (Sⁿ) can be written as

S^{n} = （ F_{^{p}}^{n} + D_{^{d}}^{n} ） ϕ + H^{ph} ψ .

(5)

From the Perturbation Theory, the effective delayed neutron fraction is defined as

β_{eff} \equiv \frac{〈 ϕ^{+} S_{d}^{n} 〉}{〈 ϕ^{+} S_{}^{n} 〉},

(6)

where, the symbol <···> denotes integration over the full phase space, φ⁺ is the adjoint neutron flux, S_dⁿ is the delayed term in the neutron source (Sⁿ) which includes the conventional delayed fission term (D_dⁿφ) and the delayed photoneutron term (H_d^phψ) in this case. The total effective delayed neutron fraction (β_eff^total) can be expressed as

β_{eff}^{total} = β_{eff}^{n} + β_{eff}^{ph} = \frac{〈 ϕ^{+} D_{^{d}}^{n} ϕ 〉}{〈 ϕ^{+} S_{}^{n} 〉} + \frac{〈 ϕ^{+} H_{d}^{ph} ψ 〉}{〈 ϕ^{+} S_{}^{n} 〉} .

(7)

The first term, β_effⁿ, is the fraction of effective delayed neutrons generated by the fission products. The second term, β_eff^ph, is the fraction of effective neutrons created by the delayed photons. In Ref. [16], β_eff^ph is represented as the sum of the fraction of effective neutrons created by the prompt photons (β_eff^ph,p) and the fraction of effective neutrons created by the delayed photons (β_eff^ph,d), i.e., β_eff^ph = β_eff^ph,d + β_eff^ph,p. In this paper, β_eff^ph,p is removed from β_eff^ph. The prompt photons are generated within 10⁻¹¹ seconds after the fission event and the photons are transported “instantaneously”, as the photon speed is much greater than the neutron speed.

Since the delayed fission neutrons are generated from a few milliseconds to a few minutes after the fission event, the photoneutrons from the prompt photons are excluded from the delayed neutron source. B. Dionne and N. Hanan evaluated β_eff^ph,d without considering the delayed photons in the calculation code[16]. The upper limit of β_eff^ph,d was represented by twice as β_eff^ph,p and the lower limit was zero, as the total energy produced per fission for the delayed photons was slightly smaller than the prompt ones. In this paper, β_eff^ph is still evaluated by β_eff^ph,p, i.e., β_eff^ph = β_eff^ph,d ≈ β_eff^ph,p based on the approximation that the prompt photons provide about the same photoneutron production rate as the delayed photons in reactor with the precursors in equilibrium. For ²³⁵U, the average yields of the prompt and delayed photons per fission are 7.2727 and 6.7273, respectively [14]. The average energies of the prompt and delayed photons are 6.61 MeV and 6.11 MeV, respectively. The approximation has been verified in BCNT, which is a reactor using beryllium reflector [25].

3.2. The k-ratio method

In 1997, M.M. Bretscher introduced the k-ratio method for a system without photoneutrons [17]. With some approximations, the effective delayed fission neutron fraction can be evaluated by the ratio between the total and the prompt multiplication factors. This is the reason why this method is called the k-ratio method. It has been widely used worldwide [18-21].

For the system with photoneutrons, one can still obtain the exact definition of the delayed fission neutron fraction based on the fact that the photoneutron term is several orders of magnitude smaller than the fission neutron term. Hence, the effective delayed fission neutron fraction (β_effⁿ) can be obtained by MCNP6 [22] calculation with and without delayed neutrons (TOTNU card) while omitting the photoneutrons, which leads to Eq. (9), where the superscript N denotes “MODE N” calculation in MCNP6 and the subscript p represents that generation of the delayed neutron is suppressed.

β_{eff}^{n} \approx 1 - \frac{k_{p}^{N}}{k^{N}} .

(9)

The k-ratio method has been extended to calculate the effective photoneutron fraction (β^ph) [16], which can be performed by coupled neutron-photon transport (MODE N P). (γ, n) reactions are modeled with MPN card. It is worth mentioning that β_eff^ph is evaluated by β_eff^ph,p in this paper, which can be estimated by Eq. (10), where the superscript “N+P” denotes “MODE N P” calculation in MCNP6.

β_{eff}^{ph} = β_{eff}^{total} - β_{eff}^{n} = \frac{k^{N + P} - k_{p}^{N}}{k^{N + P}} - \frac{k^{N} - k_{p}^{N}}{k^{N}_{}} .

(10)

3.3. The inner product method

With the same assumptions, the definition of β_eff^ph can be rewritten as Eq. (11) after spatial and energy discretization,

β_{eff}^{ph} = \frac{〈 ϕ^{+} H_{d}^{ph} 〉}{〈 ϕ^{+} S_{}^{n} 〉} = \frac{\sum_{I} \sum_{n} φ_{_{n, I}}^{+} \sum_{g'} V_{I} Σ_{γ n, g', I} ψ_{_{g', I}}^{}}{\sum_{I} \sum_{n} φ_{_{n, I}}^{+} χ_{n, I} \sum_{n'} V_{I} υ Σ_{f, n', I} φ_{_{n', I}}^{}},

(11)

where, the generation process of the photons is omitted and the photon flux Ψ is used; Σ_γ,n represents the (γ,n) reaction; V_I is the volume of cell; the index I refers to the volume cell of the spatial mesh; the indexes n and g refer to the neutron and photon energy group, respectively.

The deterministic and Monte Carlo codes were used together for the effective delayed photoneutron fraction calculations [25], the flow chart of which is shown in Fig. 2. The solution of the adjoint neutron transport equation was unavailable because of the continuous-energy treatment of nuclear data in Monte Carlo Code [18]. The adjoint neutron flux could be obtained through the adjoint diffusion calculation of CITATION [23]. The multi-group constants were calculated by the lattice code PIJ [24]. With the cut-off energy at 1.855 eV, four group (3 fast and 1 thermal groups, Table 3) constants were generated for the CITATION diffusion calculation based on R-θ-Z geometry.

Group structure in CITATION calculations

Groups	Energy boundaries	Comment
1	10 MeV–749 keV	Fast, the prompt fission neutron
2	749 keV–18 keV	Fast, the delayed fission neutron
3	18 keV–1.855 eV	Fast
4	1.855 eV–0 eV	Thermal

Fig. 2

(Color online) Calculation flow chart

3.4. Codes and nuclear data libraries

TMSR-SF1 uses coated particle fuels, so the double heterogeneity, caused by the distribution of TRISO fuel particles in the pebbles and the pebble distribution in the core, must be modeled correctly. The general-purpose continuous-energy Monte Caro code MCNP6 ver1.0 was used for explicitly representing the geometry of both fuel pebbles and the individual particle kernels to account for the double heterogeneity effect. It is worth mentioning that all MCNP6 calculations presented in this work were performed with ENDF/B-VII.0 cross sections. The deterministic codes including PIJ and CITATION were used for the forward and adjoint neutron flux calculations. The collision probability method code (PIJ) was used to solve the energy and space dependent Boltzmann neutron transport equation numerically by the collision probability method. This method has been widely used for treating the double heterogeneity energy self-shielding effect of the pebble bed reactor [1]. At every level of heterogeneity, an artificial material is generated to give transmission probability of the aggregate of each level’s constituents. The homogenized material replaces the heterogeneous geometry on the next level of heterogeneity from fine to coarse geometric detail. CITATION, the multi-group diffusion code, was used for three-dimensional core calculation using explicit description of space. Table 4 is a summary of the codes used in this work.

Summary of codes used in this work

Items	Deterministic codes		Monte Carlo
	Fuel cell code	Core calculation code
Code	PIJ	CITATION	MCNP6 1.0
Theory	Collision probability	Multi-dimensional diffusion	Monte Carlo
Model	Ball	3-D (R-θ-Z)	3-D
Nuclear data	JENDL-3.3	-	ENDF/B-VII.0
Number of groups	107(Cut-off energy 1.86 eV)	3+1(Fast + Thermal)	Continuous

4 Results

MCNP6 kept high geometric fidelity on the assumption that the pebbles were regularly distributed in the core and the coated particles were regularly distributed in the fuel zone. Meanwhile, a four energy group and 64(θ)×34(R)×85(Z) spatial meshes were modeled by CITATION for TMSR-SF1.

4.1. Neutron flux and adjoint neutron flux

To validate results calculated by CITATION, the neutron flux integrated over energy space calculated by MCNP6 and CITATION were compared. The axial neutron flux in given Z meshes was averaged over θ and R directions and the radial neutron flux in given R meshes was averaged over θ and Z directions. The axial and radial neutron flux distributions are shown in Fig. 3. The neutron flux distributions calculated by the two codes agree well, except the flux in the reflector zone.

Fig. 3

Axial (a) and radial (b) neutron flux distribution calculated by MCNP6 and CITATION.

The normalized axial and radial adjoint neutron fluxes calculated by CITATION are shown in Fig. 4. The adjoint neutron fluxes were integrated over energy space at the cut-off energy of 1.86 eV. The axial adjoint neutron flux in given Z meshes was averaged over θ and R directions and the radial adjoint neutron flux in given R meshes was averaged over θ and Z directions. It can be seen that the thermal neutron had higher weight than the fast neutron in the active zone. The adjoint neutron flux decreased rapidly in the reflector zone as expected, indicating that neutrons in the reflector zone were of much lower importance than the neutrons in the active zone.

Fig. 4

(Color online) Axial (a) and radial (b) adjoint neutron flux distribution calculated by CITATION.

4.2. The delayed photoneutron fraction

4.2.1. The k-ratio method

Using the k-ratio method, the delayed fission neutron and photoneutron fraction of TMSR-SF1 were calculated with MCNP6. As Be(γ, n) reaction rate is small, 10¹⁰ neutron histories were performed. The results are given Table 5.The critical eigenvalues were obtained with a standard deviation (1σ) of about 0.00001. The k-ratio method, which is based on calculating the multiplication factor with and without the contribution of the photoneutrons and delayed fission neutrons, is time-consuming. Also, it introduces relatively large statistical errors for the effective delayed photoneutrons fraction.

Effective delayed neutron fraction (in pcm) calculated by the k-ratio method for TMSR-SF1.

Mode N		Mode N, P
k^N	k_P^N	k_P^N+P
0.98708 ± 0.00001	0.98040 ± 0.00001	0.98712 ± 0.00001
βⁿ_eff = 677 ± 1.4		β_eff^ph =4.0 ± 1.4

4.2.2. The inner product method

With the same meshes, the axial and radial distributions of photon flux integrated over the energy space were calculated by MCNP6, As shown in Fig. 5(a), the photon flux decreased quickly from the core center to the edges. It can be deduced that photons are mostly generated in the process of fission. In Fig. 5(b), there is a small peak at about 140 cm. The main reason is neutron absorption in Hastelloy-N alloy vessel.

Fig. 5

(Color online) Axial (a) and radial (b) photon flux distribution calculated by MCNP6

The calculated photoneutron fraction combining MCNP6 and CITATION is listed in Table 6, based on the methodology explained in Section 3. It is worth mentioning that the photon fluxes are calculated by MCNP6, so there are statistical uncertainties. It is found that the weighted photoneutrons generated in the active zone take majority share, about 96.65%, while the ones outside the active region account for only 3.35%. This coincides with the photon flux distribution.

Effective photoneutron fraction calculated by the inner product method

Weighted fission neutrons in the active zone		4.44E-15 (±2.70E-20)
Weighted photoneutrons in the core		1.43E-19 (±3.62E-22)
β_eff^ph (pcm)		3.22 (±0.01)
Regions with coolant	Coolant volume (m³)	Contribution for β_eff^ph
Active zone with fuel pebbles	0.99	96.65%
Top and bottom reflector zone	0.50	1.26%
Side reflector zone	1.20	2.09%

5 Conclusion

The k-ratio and the inner product methods were used to estimate the effective delayed photoneutron fraction in TMSR-SF1. The Monte Carlo code MCNP6 was used to calculate the photoneutron fraction using the k-ratio method. The inner product methods were performed by Monte Carlo code and deterministic codes together. The neutron flux calculation of CITATION was validated by Monte Carlo code. It is shown that the neutron flux distributions performed by the Monte Carlo code and deterministic codes agree with each other in general. The effective delayed photoneutron fraction calculated by the k-ratio method is about 4.0 pcm (± 1.4 pcm). The result calculated by the inner product method is 3.22 pcm (± 0.01 pcm). The results show that there are good agreements between two methods.

Influence of the effective delayed photoneutron fraction on the dynamic characteristics will be analyzed for TMSR-SF1.

References

T. Allen, S. Ball, E Blandford et al.,

Preliminary fluoride salt-cooled high temperature reactor (FHR) subsystems definition, functional requirement definition and licensing basis event (LBE) identification white paper

. University of California, Berkeley, USA, UCBTH-12-001(2012)