Comparison of statistical evaluation of criticality calculations for reactors VENUS-F and ALFRED

Limitations of correct evaluation of keff in Monte Carlo calculations, claimed in literature, apart from the nuclear data uncertainty, need to be addressed more thoroughly. Respective doubts concern: the proper number of discarded initial cycles, the sufficient number of neutrons in a cycle and the recognition and dealing with the keff bias. Calculations were performed to provide more information on these points with the use of the MCB code, solely for fast cores. We present applied methods and results, such as: calculation results for stability of variance, relation between standard deviation reported by MCNP and this from the dispersion of multiple independent keff values, second order standard deviations obtained from different numbers of grouped results. All obtained results for numbers of discarded initial cycles from 0 to 3000 were analysed leading for interesting conclusions.


Introduction
As stated by F.B. Brown [1] the limitations of correct evaluation of keff to be addressed more thoroughly are:  Discarding a sufficient number of initial cycles (neutron generations) so that the contamination of results is negligible. Performing calculations of keff for sufficient number of neutrons (batch) in each cycle. Recognizing and dealing with the bias in the keff uncertainties.Also Mervin and Maldonado [2] discussed the problem of under-prediction in the uncertainty of keff evaluation by MCNP.To provide some information on these points, when keff and its standard deviation are calculated for a fast subcritical core, calculations were performed with the use of the MCB code [3] for reactors VENUS-F (subcritical) [4,5] and ALFRED [6,7].Methods and results are presented in the following sections.
The scheme of calculations with the use of KCODE was organized in such a way that the results would be convenient for applying the Hartley test of variation [8] and at the same time the calculation will not be excessively long.Respective to the Hartley tables of the test:  the highest number of repetition in a series for which the values in table are evaluated is m = 12, and  the highest number of results (generations) in a series is l = 61.The idea of calculations was to assess the dependence of keff and its standard deviation (SD) on the increasing number of discarded initial cycles (DIC) and constant number of accepted ones (61) and also to compare different ways of the keff (SD) evaluation.
The assumed batch of 10000 neutrons was applied that is regarded as sufficient for producing reliable precision of evaluation.It allowed for performing the calculations for 0, 1, 2, 3, 4, 5, 10, 15, 20, 25, 50, 100, 200, 400, 800, 1500 and 3000 DIC (for 100, 200 and 800 DIC calculations were repeated 2 or 3 times).More discarded cycles resulted in excessive computation time.For reducing the uncertainty of keff and better comparing the standard deviation evaluated by MCB (MCNP), called also the apparent SD [2], with that derived from the dispersion of keff, called the real SD [2], the sets of 12 results were repeated 5 times so the total number of results for keff and its SD for each number of DIC amounted to 60..

Results
For each of the 210 (2 x 5 x 21) rows in Tables 1 and 2, for different number of discarded cycles, the Hartley test was applied answering the question whether all variances from series of m = 12 repetitions are statistically equal.The tested variances taken from MCB code results were based on the least squares solution of a multivariate linear regression of the cycle keff estimates, i.e. collision, absorption and track length, used in MCNP code and reported by T.J. Urbatch et al. [9].Only in one of the 210 cases, both for VENUS and ALFRED, the test answer was negative for the confidence level =0.01 (in the first raw for 100 discarded cycles; keff = 0.96614).This effect, can be regarded as stochastically justified.Thus, the conclusion is that the keff variance produced by MCB code is stable within the 12 repetitions, starting from 0 up to 3000 DIC.Next, for the 60 results (5 sets x 12 repetitions) the two SD were evaluated: one, as the average value of SDs computed for each keff by the MCB code, and the other, from the dispersion of 60 values of keff , in the same way for both reactors.Also the effect of grouping results in subsets was assessed.The subsets of the following multiplicity m can be obtained: 2, 3, 4, 5, 6, 10, 12, 15, 20, 30 and 60.For all of the subsets keff and its both SDs were calculated.The numbers of results n for each of the 3 variables were equal to: 30, 20, 15, 12, 10, 6, 5, 4, 3, 2 and 1, respectively; m x n was always equal to 60.The following definitions of variables are used to present the obtained results, while symbols  are used for SD produced directly by MCNP (apparent) and s for SD resulting from the dispersion of MCNP results (real): ki -keff from MCB for l = 61 generations of 10000 source neutrons kbar = -average keff from one group of results (Eqn. tot -average from all group averages (bar), i.e. the estimated SD of the multiplication factor kbar (Eqn.4) and -the estimated SD of the system multiplication factor ktot, i.e.SD1 (Eqn.5) So, the ktot , tot and for every set of m and n is the same and equal to the one for all 60 results in one group (or 60 "groups" with single result in each).On the way of group calculations we can also evaluate SDs from the dispersion of group averages: of ktot ( ), i.e. the SD from dispersion of kbar (Eqn.6) and of tot ( ) from dispersion of bar (Eqn.7): (6) (7) Another SDs can be obtained from: the dispersion of ki within a group and dispersion of this SD between groups.They are: sj -SD of the single value of ki in the j-th group (Eqn.8) (8) sbar -average of sj from n single group results (Eqn.9) (9) stot -SD of the average of sbar (Eqn.10) (10) One can further cheque the role of such values as (average SD of , from the dispersion of i in groups) (Eqn.11) and further SD of tot -(Eqn.12).
(    Exemplary results of the above described evaluations, for 50 DIC, are presented in Table 1 and on the graphs in Figs 4 and 5. Further, influence of the number of DIC on keff SD ( ) was assessed, both for VENUS and ALFRED.The results are presented in Fig. 6.One can see that the SD is distinctly higher only in case when no initial cycle (generation) is discarded.The difference is more evident for higher number of groups.It is not the case for keff itself (see Fig. 8) when 3 cycles should be discarded as minimum.In Fig. 7 we present values of the same SDs as in Fig. 5, but for all numbers of DIC.This allows for better observation of the effect of grouping on the quality of obtained SD evaluation.It seems obvious from Fig. 7 that for low number of groups (n < 10) the spread of results is very high, yielding both overestimated or underestimated SDs.For more groups (n ≥ 10) with lower number of results (m ≤ 6) the differences are lower and usually the lowest fluctuations between different DICs are observed for 60 "groups" with single result.However, in this case we are not able to evaluate the other SDs shown in the Table 1.In 7 are results of the same SDs but obtained in opposite sequence.In case of SD*, after calculation SDs of in groups, the respective average value was obtained.As the result we can observe more smooth dependence but with the same tendency with regard to the number of DICs.

Conclusions
The following conclusions are only valid for a fast reactor type similar to reactors VENUS and ALFRED used in calculations.
1.With the use of Hartley test of variation it is proven that the keff variance produced by MCB code is stable within 12 independent repetitions, starting from 0 up to 3000 discarded initial cycles.This is thus most probable that the variance will be also stable for all cases of more numerous repetitions and DICs. 2. For 50 DIC, the average SDs, of a single value of keff, obtained from the MCNP evaluation method, -tot (0.000756 for VENUS and 0.000657 for ALFRED) are almost equal to: the same SDs resulting from the dispersion of keff -sbar (0.000748 for VENUS and 0.000623 for ALFRED), as the average from group calculations (without 60 results group), and -((000833 for VENUS and 0.000655 for ALFRED), from all 60 results of keff.3. The respective SD of the average value of keff (from 60 results), called further the first order SD, are: for VENUS (0.000098), stot (0.000097) and (0.000108), and respectively for ALFRED (0.000087, 0.000080 and 0.000085.4. Standard deviations of the above SDs (second order SD) are: (VENUS -0.000009 and ALFRED -0.000007) for the average value of tot and evaluated from the dispersion of 60 i values, as well as: and respectively: VENUS -0.000009 and ALFRED -0,000007 and VENUS -0.000010 and ALFRED -0,000008, both on the basis of group calculations.5. Thus the final result of the VENUS and ALFRED cores evaluation of keff (for 50 DIC) are: 0,966130 ± (0,000098 ± 0,000009) and 1,010099 ± (0,000085 ± 0,000007).This result recalculated for control rods submerging effect, from 40 cm (simulation) to 47 cm (experiment), reads 0.96505, very near to the experimental value: 0.96415±0.00086[5].Values of the first and the second order SD are only slightly dependent on the way of calculation, whether basing on the MCB (MCNP) methodology or on the observation of dispersion of the MCB results for keff and.6.The graphs in Figs. 4 and 5 show certain differences in the spread and distribution of results obtained from group evaluation of SDs: and are more dispersed Energy and Fuels 50 between different number of groups then stot and , in turn distributions of the later are more smooth.This may result from the sequence of operations in calculation.The smooth distributions are obtained when the calculation of SD precedes the operation of averaging, while the dispersed distributions appear in the case of opposite sequence.This effect is also visible for all numbers of DIC in Fig. 7. 7. Results obtained for different number of DIC show that keff SD ( ) is distinctly higher only in case when none of initial cycles was discarded.It is not the case for keff itself when minimum 3 cycles should be discarded.For the dependence of SD on the number of groups it was possible to see the effect of grouping on the quality of evaluated SD.It seems obvious that for low number of groups (n < 10) the spread of results is very high, yielding both overestimated or underestimated SD.For more groups (n ≥ 10) with lower number of results (m ≤ 6) the differences are lower and usually the lowest fluctuations between different DICs are observed for 60 "groups" with single result.However, the grouping is useful for evaluation of the second order SDs.This work was supported by the EC in the frame of FREYA project No 269665 and by Polish Ministry of Higher Education.All calculations were carried in the Academic Computer Centre CYFRONET, AGH.

Fig. 3 .
Fig. 3. Scheme of the calculations of

Table 1 .
Exemplary numerical criticality statistics evaluation of reactors VENUS and ALFRED based on the results from MCB code modelling (for 50 DIC).