The ROCOCO system (Kislitsyna and Nikolayev 2016) is designed to supply constants for Monte-Carlo calculations of both neutron fields and the gamma fields they generate. The initial database of nuclear data used in the system is the Russian national library of evaluated neutron data (ROSFOND) (Nikolayev 2006, Zabrodskaya et al. 2007, RUSFOND 2017). ROCOCO is specific in that it enables the estimator to optimize the level of detail when describing neutron cross-section energy dependences. Cross-sections of key fuel, structural and coolant materials can be described in such detail as the evaluated data permits; cross-sections of secondary nuclides (minor actinides, fission products, etc.) can be described in a 299-group BNAB approximation (Manturov et al. 1996) with regard for the resonance self-shielding by subgroup method or without regard for self-shielding altogether. The energy dependence of gamma rays is described in a 127-group P5 approximation (Koshcheyev et al. 2014). Optimizing the level of detail makes it possible to reduce to a great extent the counting time with no major effect on the result and its error. Where desired, in the process of calculating the energy release in neutron reactions or in gamma-quanta formation matrices, contributions from the decay of radionuclides formed in these reactions (with a half-life of less than three years) can be taken into account. The energy dependence of the elastic scattering anisotropy is described in detail, or in the event of a group or subgroup description of cross-sections, by defining 33 boundaries of 32 equiprobable cosine intervals of the scattering angle. The thermalization effects in calculations of neutron fields are taken into account either in an ideal gas approximation or using 72-group thermalization matrices built based on thermalization files contained (if any) in the ROSFOND library.

It should be noted that the system contains descriptions of detailed dependences of elastic scattering cross-sections and angular distributions on all multi-isotope elements; the relationship between the scattering angle and the energy loss in this case is determined with the use of the energy-dependent effective atomic weight.

The system’s programs are written in the FORTRAN language. The system is easily integrated in Monte-Carlo codes.

Nuclear data, radiation field Monte-Carlo calculation, combination of detailed and group descriptions of cross-section energy dependence

ROCOСO (Russian abbreviation for Calculation and Organization of COmbined COnstants) is a system designed to supply constants for Monte-Carlo calculations of neutron fields and the photon fields they generate in nuclear reactors. It enables detailed description of the neutron cross-section energy dependences (as it is done in MCNP, a commonly known code, (Z-5 Monte-Carlo Team 2003). However, the cross-sections of the nuclides contained in the reactor in low concentrations can be also described with a lower level of detail, that is, in a multi-group approximation, with regard for the resonance self-shielding, by subgroup method or without regard for self-shielding altogether. This naturally increases the rate of computations which is essential in the event of multivariant engineering calculations. Therefore, ROCOCO makes it possible for the user to determine by himself/herself the level of detail with which the energy dependences of cross-sections are described.

The ROCOCO system programs are written in the FORTRAN language. The system operability was tested using various compiler versions (GNU, LAHEY, INTEL, etc.). The system consists of three parts.

The first of these is the COMBINE package of programs designed to form a binary library of initial constants based on evaluated data files in the ENDF6 format (RUSFOND, a Russian library of evaluated neutron data files) and of group constants in the BNAB format (Manturov et al. 1996). The code developers are the users of this system part. In particular, the authors used the COMBINE system based on files of the ROSFOND library and the BNAB-RF group constant library (Koshcheyev et al. 2014) to create a version of the BIBIN-2016 binary library. There is also BIBIN-B7, a version created based on the ENDF/B-VII library. In future, where required, other versions of the binary library can be also formed. The COMBINE system has been described in sufficient detail, so new versions can be created without the authors’ involvement as well.

The second part is used directly in the calculation. The procedures in this block read the computational assignment and an online library is formed on its basis for the required composition of materials and in the format convenient for simulation. In particular, inelastic scattering, thermalization, and fission spectra matrices are presented as accumulated probabilities of processes. It is also here that subgroup parameters of resonance isotopes are generated for the preset temperatures and compositions, while constants with detailed tracking of the neutron energy are recounted for the preset media temperatures. All this data is generated only for the isotopes present in the medium and is recorded into a dedicated exchange file that does not depend on the computation platform (what is meant is a binary or text file of a universal format read similarly by various compilers and in different operational systems, this enabling the user to prepare a library, e.g., on an office computer and to perform calculations on a supercomputer). The procedures in this block are united into an independently executed program module, COLIBRY (Kislitsyna and Nikolayev 2016a).

The third part is intended to be used directly in Monte-Carlo simulation programs and consists of several procedures. The key procedure, COLIBRY, is called only once in the process of a calculation. It refers to the exchange library prepared by the system’s second part, and accommodates them in respective general arrays which are accessible to all operating procedures.

The operating procedures are listed below. The first of them, GEFIS, has the parameters INUC, NGFIS, ig, E, and ves (here and hereafter, the parameters of the procedures are underlined for identification). Evidently, the fission spectrum is defined by which INUC nuclide fission gives birth to the neutron and which NGFIS energy group the neutron that has caused the fission belonged to. The procedure based on these parameters defines the group ig and the energy E of the fission neutron. The parameter ves is equal to unity by default. It is possible to increase artificially the yield of high-energy neutrons, in which case the parameter ves serves to compensate for the fission spectrum distortion. It is assumed in the first generation, when no INUC and NGFIS have yet been defined, that the neutron formed by the fission of the first of the listed neutrons (a priori fissionable) under the action of neutrons in the 50^th group (in which all fissile isotopes have fission cross-sections other than zero). We shall note that there is an extended version of this procedure, called GEFISBETA, which, apart from the functions performed by the GEFIS procedure, forms the vesD array, eight components of which represent probabilities that the emitted fission neutron (in the group ig and with the energy Е) has been formed by the decay of the precursor for the respective delayed neutron group. The GEFISBETA procedure is intended for use in calculating effective fractions of delayed neutrons.

The SIGMA procedure has the parameters SIGMAC, SABS, SFNU, KERM, E, NG, and m. For a neutron with the energy Е, lying in the NG group in the physical zone m, it defines the full macroscopic cross-section SIGMAC, SABS is the macroscopic cross-section of absorption, the “fission neutron formation cross-section” SFNU = n×S_f , and KERM is the kerma value. The group’s energy and number are exchanged among the operating procedures automatically. To determine these values, the SIGMA procedure calculates initially the microscopic constants of all isotopes in the zone m for the energy of the considered neutron. If the neutron has formed in the considered zone, then the cross-sections of all nuclides in this zone are determined for the energy of this neutron. If the neutron enters this area after crossing its boundaries, microconstants are determined only for the nuclides for which they have not been determined in the previous zones on the neutron flight path. Microconstants are entered in shared arrays to be further used in the COLLY procedure and, where required, in recording procedures (other than being part of the ROCOCO system). The SIGMA procedure also processes data for gamma quanta to which groups numbered 301 through 427 are attached (if the calculation assignment requires the gamma quanta paths generated by neutrons to be traced: the number of gamma groups is shown as equal to 127).

The COLLY procedure has the parameters m, NG, E, r, n, etha, NG1, E1, mulab, yi, ngam, and ygam. For the collision of a neutron taking place in the zone m, with the energy Е lying in the group NG, it defines the following type of collision: r = 1 is the elastic scattering; r = 2 is the inelastic scattering; and r = 3 is the absorption. Also n is determined as the number of the nuclide with which the collision took place; the average number of secondary fission neutrons formed in this absorption, E1 and NG1 are the energy and the number of the scattered neutron group; mulab is the scattering angle cosine in the laboratory coordinate system; yi is the set of scattered neutrons (exceeds unity if there is a contribution of the (n, xn) reactions to the inelastic scattering); ngam is the group number of the gamma quanta emitted in the inelastic scattering or capture; and ygam is the number of these quanta (301 ≤ ngam ≤ 427).

Therefore, it is enough for the connection of the ROCOCO system to the reactor Monte-Carlo code that the COLIBRY subprogram (to the inlet of which only the exchange library file is supplied) is called prior to the solution of a particular problem. Then, the GEFIS procedure (defining the fission neutron energy and group) shall be called at the beginning of each selected path. The SIGMA procedure is called prior to tracing the track of the neutron (or the photon) lying within the zone in which the considered neutron or photon was born as the result of a process, or within the zone it entered when crossing the boundary. The COLLY procedure is called for the consideration of each collision.

Therefore, the procedures provide the calling program with all information required for tracing the neutron path and, where required, the path of the gamma quantum born in the neutron reaction. Also values are provided which are required to estimate the breeding ratio and the energy release in terms of range lengths, collisions and absorptions. Since, for each of the neutron track portion lying within the physical zone, the SIGMA procedure determines microconstants for all of the nuclides it contains, the program referring to ROCOCO may include recording procedures that estimate such functionals as the isotope reaction rates, disabled multigroup constants, etc.

The following is shown in the input file for the second block procedures the COLIBRY module includes:

MZONE – the number of physical zones (i.e., of zones differing in composition or temperature);

NIS – the full number of nuclides in at least one of these zones;

the number of the BNAB neutron groups MGRN ≤ 299 – in groups with numbers larger than 227 (i.e., with energies smaller than 4.64 eV), the thermal motion of neutrons is taken into account in calculations; a smaller number of groups (e.g., 298) needs to be specified for taking this into account;

the MGRG gamma group number equal either to zero or to 127 (as stated above, group numbers 301 through 427 are attached to gamma quanta);

NUCLIST – a list of NIS nuclides: four-symbol names of isotopes and two-symbol names of multi-isotope elements;

ID – the array of NIS signs for the level of the cross-section structure description detail for each nuclide: ID = 2 with a detailed description of the cross-section structure; ID = 1 for the use of a 299-group approximation with regard for the intragroup cross-section structure by subgroup method; ID = 0, a 299-group approximation without regard for resonance self-shielding;

THERM – the set of the numbers of thermalization matrices for the BNAB system of constants (if this number is assumed to be equal to zero, then thermalization will be taken into account in a free gas approximation);

TZONE – the array of the physical zone temperatures (MZONE of values) (nuclides with ID = 1 or 2 contained in zones with different temperatures will be considered as different nuclides in the calculation);

RO – a two-dimensional array of nuclear concentrations: concentrations of all nuclides should be shown for each zone.

It should be noted that the above keywords of the assignment parameters can be substituted by analogs adopted in the CONSYST constants supply system (Nikolayev M.N., VANT. Ser.: Yadernye konstanty. 1999, No. 2, p. 148).

In calculations, the cross-sections of “detailed” nuclides outside the resonance region are described with the same level of detail as they are presented in evaluated data files. Detailed energy dependences of the cross-sections, calculated for the temperature of 300 K using the RECONR and BROADER modules of the NJOY 2В code, are used in the region of resolved resonances. For zones with temperatures higher than the indoor temperature, detailed cross-sections are recalculated using the DOPPLER procedure that implements the algorithm similar to that used in the BROADER module (NJOY99.0 2000). In the region of unresolved resonances, the structure of cross-sections is described by subgroup method. Subgroup parameters are computed in the COLIBRY program based on tables of the resonance self-shielding factors for the BNAB-RF constants system (Koshcheyev et al. 2014). Energy-angular distributions of neutrons that have undergone inelastic scattering or have formed as the result of the (n, xn) reactions are described in a 299-group approximation using matrices of probabilities and mean cosines of intergroup transfers from the BNAB-RF constants system. The experience of comparing 299-group calculations with calculations with a detailed description of inelastic interaction processes has shown that errors not only in the breeding ratio but also in group flows turn out to be smaller than the errors caused by the inaccurate knowledge of these processes. As to inelastic scattering, calculations use angular distributions presented as 33 boundaries of 32 equiprobable intervals for the scattering angle cosine in the inertia center system calculated based on initial files with such level of detail for energy with which these distributions are presented in them. These angular distributions for all isotopes were formed using the ACER module (MacFarlane and Kahler AC 2010) of the NJOY code.

The cross-sections of “subgroup” (ID = 1) and group (ID = 0) nuclides are described in a 299-group BNAB approximation. And the resonance structure of cross-sections for “subgroup” nuclides in the region of both unresolved and resolved resonances is described using subgroup parameters obtained based on tables of self-shielding factors interpolated to the temperature of the zone under consideration. In the region of unresolved resonances where two subgroups are normally enough to describe the structure of cross-sections, subgroup parameters are determined using a simplified algorithm. In the region of resolved resonances where two subgroups are often not enough, the universal procedure contained in COLIBRY is used to obtain subgroup parameters based on self-shielding factors (Zhemchugov 2017). As to elastic scattering, a detailed description of the energy dependences of the boundaries of equiprobable intervals are used independent of the ID value.

Energy-angular distributions of neutrons in the thermalization region are described either in a multigroup approximation using the BNAB-RF thermalization matrices, or based on a free gas model using the algorithm described in (Eriksson 1970). As mentioned, the temperature dependence of detailed cross-sections is taken into account using the DOPPLER procedure that reduces the BIBIN-contained detailed cross-sections matching the indoor temperature to the temperatures given in the calculation assignment.

The energy release in conditions of elastic or inelastic scattering or during entrainment or fission is calculated using 299-group BNAB constants. The energy dependence of energy release within each group has the same level of detail as the description of cross-sections. Where the calculation assignment requires the distribution of gamma quanta formed in neutron reactions (MGRG = 127) to be calculated, tables of local energy release are used, or, in the opposite case, tables of full energy release are used which also take into account the energy entrapped by gamma quanta. In both cases, the energy release data takes into account the contribution from the decay of neutron reaction products (with a half-life of up to three years).

Formation of gamma quanta is taken into account for the entrainment, inelastic scattering and fission reactions. 299×127-group matrices of gamma quanta yields from the BNAB-RF are used. These matrices also take into account the yield of gamma quanta in the decay of products from respective reactions. The gamma quanta transport is calculated in a 127-group approximation using macroscopic constants. The anisotropy of both coherent and incoherent scattering is taken into account by specifying probabilities of scattering for each of the three specified angles with respect to each intergroup transfer. The probabilities and the angles have been determined from the condition of five angular momentums being preserved.

We shall note the specific nature of the neutron interaction with nuclei of natural multi-isotope elements. The format of the evaluated data file libraries enables representation of data for natural mixtures of multi-isotope elements, but, in reality, this capability is used to describe the cross-sections of only several materials rarely employed in engineering practice (e.g. such as platinum). ROCOCO enables a detailed description of the cross-section structure for all natural multi-isotope elements. The cross-sections of all isotopes are reduced to a single energy mesh which is thinned by dropping the energy nodes in which the cross-sections of elements can be restored by linear interpolation from neighboring nodes with the preset accuracy (normally 0.1%). Energy dependences of the cross-sections of elements are determined by adding the isotope cross-sections to the weight of their concentrations. Angular distributions of neutrons, elastically scattered on isotopes and presented in the form of expansions in Legendre polynomials, are taken directly from the evaluated data files. The expansion coefficients for all isotopes are reduced to a single energy mesh and are averaged then with the weight of concentrations and the elastic scattering cross-sections. At each point of the energy mesh, detailed dependence of the scattering probability on the angle of scattering in the inertia center system is restored after which 33 boundaries of 32 equiprobable intervals are determined from the scattering angle. In parallel with the cross-sections, the energy dependence of the effective element mass, determined such that the average loss of energy during scattering is preserved, is calculated and stored. As in the event of individual elements, energy-angular distributions of inelastic scattering are described in a multi-group BNAB-RF approximation. The use of the energy dependences of cross-sections for natural elements reduces greatly the number of the nuclides the medium contains, which reduces the cost of the machine time for the calculation. The drift of the results when using this approximation is practically negligible.

To check the serviceability of the ROCOCO system, a program, kinfin, was written to calculate the criticality of infinite media. It enables estimation of breeding ratios for absorptions, k_inf.0, and for collisions, k_inf.1. No path lengths were estimated since, in an infinite homogeneous medium, a flying neutron will definitely experience a collision in this medium, and the estimates for collisions and path lengths will differ only in that the error will be somewhat higher in the latter case (the sample estimation error for the average path length will add to the estimation error for collisions). The multiplication of neutrons thanks to the (n, xn) reactions was taken into account by multiplying the neutron weight by the average multiplicity of the neutrons emitted as the result of inelastic scattering processes. The paths were calculated with the analog consideration of absorptions. The estimates obtained in the several first generations were not taken into account for all generations in the averaging of the results, as a steady-state spectrum of fission neutrons was established in these generations. After the selection among all generations, the results obtained for all generations (except a number of those dropped) were averaged, and the dispersions of each of the estimates (D₀ and D₁) and the coefficient of correlation between them were computed (r). After that, based on the least square method, a combined estimate of the breeding ratio (k_inf) and its dispersion (D) were calculated:

comb = (1 – C)×k_inf.0 + C×k_inf.1,

where

The kinfin program was intended for verification calculations. The verification was done by comparing the calculation results with the data obtained using the MCNP-5 code “tied” to the ROSFOND library. To make it possible to judge on the significance of the discrepancies in the results, media were chosen for the comparison, for which the results of the experimental k_inf determination were known (with some uncertainties). Another measure of the significance of the discrepancies was the difference between the results of the MCNP-5 calculations using the ROSFOND and ENDF/B-7 libraries. The MCNP-5 calculations were performed in the same way as the kinfin calculations, with the analog consideration of absorptions. The results of the k_inf calculation for the four media will be presented below.

The first of these was uranium medium with an enrichment of 5.56%, SCHERZO-5.56 (Chaudat et al. 1974). The fact that the enrichment of 5.56% makes k_inf equal to unity (with an accuracy of 0.25%) was revealed as the result of a computational analysis in a series of experiments performed at the SNEAK, ERMINE, and HARMONIE critical test facilities. The calculation results were obtained on one and the same single-core PC. There were 1020 generations of 2000 paths each selected among in each calculation. The first twenty generations were preliminary, and the further 1000 were active; it was in these that the calculated values were estimated. The results were as follows:

kinfin: k_inf = 0.99750 ± 0.00041 (the counting time is 107 s),

MCNP(RF): k_inf = 0.99778 ± 0.00043 (the counting time is 174 s),

MCNP(B7): k_inf = 0.99288 ± 0.00043 (the counting time is 140 s).

Therefore, the difference between the kinfin and MCNP calculation results, amounting to 0.028%, is much less than the statistical error of this difference equal to ±0.56%. The counting time in the event of kinfin for this problem turned out to be one and a half time less than in the event of MCNP (obviously, thanks largely to the advantage of a multigroup description of the secondary neutron spectra from inelastic scattering and the (n, xn) reactions).

The results of the calculations based on the ROSFOND library agree with the experimentally estimated value of 1.0000 ± 0.0025 within the latter’s error limits.

The result of the calculation based on the ENDF/B-7 library is 0.5% as less as the experimental value, which is twice as high as the experimental error.

A stainless-steel medium with a small addition of high-enrichment uranium was chosen for another test calculation and tested at the KOBRA test facility at IPPE (KBR-9 assembly). The experiment conditions and the procedures for introducing computational corrections for the heterogeneity and the finite dimension of the investigated zone (having undergone an international expert appraisal) are described in detail in (NEA Nucl. Sci. Committee). The content of uranium with an enrichment of 90% was 90.2% by the number of nuclei. The rest of the components were iron, nickel, and chromium with small amounts of carbon, oxygen, silicon, titanium, and manganese impurities (Table 1).

Table 1.Download as CSV 

Composition of the KBR-9 medium (numbers of the nuclei in the volume are presented, barn×cm).

C	3.1453E–4	Mn	8.2518E–4	Si	6.7256E–4	Ti	4.7341E–4	235U	3.7041E–4
O	8.3126E–4	Fe	4.7217E–2	Ni	6.5656E–3	Cr	1.2788E–2	238U	4.0125E–5

The experimental estimate of k_inf for this medium is 1.108 ± 0.009. As in the event of SCHERZO-5.56, 1020 generations of 2000 paths each were selected among; the data for the first twenty paths were dropped. The results are as follows:

kinfin: k_inf =1.10807 ± 0.00051 (the counting time is 83 min),

kinfin (elem.) k_inf =1.10840 ± 0.00041 (the counting time is 64 min),

MCNP(RF): k_inf =1.10841 ± 0.00041 (the counting time is 103 min),

MCNP(B7): k_inf =1.11139 ± 0.00050 (the counting time is 102 min).

The result obtained based on kinfin, with the selection of paths by isotopes, coincided with the MCNP result obtained based also on the ROSFOND library. These results also coincide with the experimental value. Practically the same result was also obtained with the composition of structural materials defined by elements, which provides for a marked reduction in the counting time. The MCNP result, using the ENDF/B-7 library, leads to k_inf being slightly overestimated.

Similar results were obtained also in calculations for the KOBRA-investigated media with k_inf ~ 1 containing nickel (KBR-7), molybdenum (KBR-10), chromium (KBR-15), and zirconium (KBR-16) as the base material. The results of the calculations with the neutron energy traced in detail coincided with those obtained using MCNP-ROSFOND in the limits of statistical errors, but required a 20% smaller time. With the composition of the media having been described by elements, the time gain was 30 to 40% greater.

Two more test problems were considered. The first of them considered an infinite medium with the same composition as the unburnt core of the SVBR fast reactor. It contained 23 nuclides (including five fissile nuclides). The obtained results are as follows:

kinfin: k_inf = 1.2465 ± 0.0005 (the counting time is 20 min),

MCNP: k_inf = 1.2459 ± 0.0005 (the counting time is 21 min).

The counting time decreases to 17 min with the zone composition being described by elements.

The next problem considered the medium with the same composition as the burnt-up (to 11% of heavy atoms) SVBR reactor core containing 169 nuclides. The kinfin calculations described in details 52 nuclides the concentration of which exceeded 0.5%. The rest were described in a group approximation. The results are as follows:

kinfin: k_inf = 1.1436 ± 0.0003 (the counting time is 52 min),

MCNP: k_inf = 1.1437 ± 0.0003 (the counting time is 209 min).

The obtained data shows graphically the advantage of a group description of nuclides contained in miniscule concentrations. This is only natural since one has to search for the cross-sections matching the energy resulting from the selection process for all nuclides independent of their concentration.

Finally, aqueous solutions of high-enrichment uranyl-nitrate and plutonium were calculated. The characteristics of the solutions and the k_eff calculation results are presented in Table 2.

Table 2.Download as CSV 

Results of calculating the criticality of aqueous plutonium and high-enrichment uranium solutions using different codes and libraries of constants

Solution	MMK-RF	MCNP5 (ROSFOND)	MCNP5 (ENDF/B-VII)	Difference from MCNP5 (ROSFOND)
				MMK-RF	ENDF/B-VII
Pu+H2O	0.99388(8)	0.99369(10)	0.99404(10)	0.02%	0.04%
U+H2O	0.90857(8)	0.90880(10)	0.90889(10)	–0.03%	0.01%

No major differences in neutron spectra have been detected (it is only with energies of < 10^–3 eV that the MNCP-5 neutron fluxes turned out to be 2 to 5% smaller than those calculated based on kinfin).

Testing has shown that ROCOCO makes it possible to obtain the same results as the MCNP code with smaller time expenditures thanks to the group description of inelastic scattering indicatrices and with the possibility of a group or subgroup description of the cross-sections of nuclides with a low concentration and a detailed description of the cross-sections of elements (for which there are no evaluated data files).

The connection by V.B. Polevoy of ROCOCO to the MMKFK code (Frank-Kamenetskiy 1978) made it possible to reveal and eliminate a number of drawbacks. There are naturally other drawbacks which will be eliminated as they are revealed. It is to be supposed that connecting ROCOCO to other codes will not lead to major difficulties as well.

It makes sense to use the ROCOCO system in precision codes for engineering calculations of power reactors in which increased speed of response will make it possible to solve problems more complex than criticality problems (solution of space and time problems, estimation of local functionals, etc.), to improve the quality of calculations, and to expand the application scope (Osetskaya et al. 2017, Osetskaya et al. 2017a).

The authors express their profound gratitude to V.N. Koshcheyev, G.B. Lomakin, and V.V. Korobeynikov for their assistance with the study; to V.B. Polevoy who demonstrated the feasibility of connecting the ROCOCO system to engineering codes, and to Ye.V. Zhemchugov who developed and provided the SUBGGROUPS program, one of the code’s key subprograms.