SOME ASPECTS OF MONTE CARLO M E T H O D APPLICATION TO N U C L E A R REACTORS ANALYSIS E. A. GoMiN, L. V. M~xioRov and M. S. YC[)KEVlCH Kurchatov lAg, Kurchatov Square. Moscow. 123182, U.S.S.R.
ABSTRACT - The MCU project is d e s c r i b e d . Its a i m is to guarantee the h i g h p r e c i s i o n of c a l c u l a t i o n of reactor parameters. The project includes creation of MonteCarlo codes, their verification and validation, development of c o n s t a n t b a s e a n d m e t h o d s of a p p l i c a t i o n Monte-Carlo calculations for i n c r e a s i n g the precision of e n g i n e e r m e t h o d s for p r e d i c t i n g of the m a i n reactor parameters.
i.
INTRODUCTION
For a long t i m e the MCU project is d e v e l o p e d in K u r c h a t o v IAE { M a i o r o v a n d Yudkevich, 1988). T h e a i m of the pro.ject is to p r o v i d e the high precision of c a l c u l a t i o n of r e a c t o r p a r a m e t e r s . The project includes creation of Monte-Carlo codes, their verification and validation, development of c o n s t a n t b a s e a n d the m e t h o d s of a p p l i c a t i o n of Monte-Carlo calculations for i n c r e a s i n g the precision of e n g i n e e r m e t h o d s for p r e d i c t i n g the main reactor parameters. A s h o r t r e v i e w of s o m e of these works is presented here. T h e d e t a i l s a r e g i v e n in r e f e r e n c e s .
2.
2.1
General
MCU
CODE.
description.
M a i n f e a t u r e s of M C U c o d e a r e as follows. Code architecture, based on module principle, a l l o w s to compute functionals for transport equation solution with the help of a wide set of different algorithms, mutually-changeable p r o g r a m m o d u l e s for i d e n t i c a l p u r p o s e s with different constants libraries. The choice of d i f f e r e n t o p t i o n s by user depends on his purposes, precision requirements and planned job solution expenses. Constant supplement is b a s e d on r e c o m m e n d e d neutron data. Code permits to use d e t a i l e d d e s c r i p t i o n of c r o s s s e c t i o n s e n e r g e t i c dependence in fully resolved resonances region. Collision modeling in t h e r m a l region is made t a k i n g into a c c o u n t c h e m i c a l b o u n d s , n u c l e i thermal moving and coherent effects.
211
212
E . A . Go'~ux el al.
It is p o s s i b l e to m o d e l a w i d e c l a s s of t h r e e d i m e n s i o n a l systems on the b a s e of a l ~ o r i t h m s either specialized or ~eneral. Calculation of some s t a n d a r d set of the p r a c t i c a l l y important functionals and these for check of i n g e n i o u s e n g i n e e r i n ~ m e t h o d s m a y be d o n e . T h e last p o r t a b l e v e r s i o n of c o d e - M C U - 2 . 0 is written in FORTRAN-77 by M C U - t e a m w h o s e m e m b e r s a r e the authors of this report and M.Gurevich, S.Marin, A.Romanov S o m e m o d u l e s of o t h e r authors are included in MCU e i t h e r [see r e f e r e n c e s ) . Now MCU is used in computers with operating s y s t e m s SVM, OS (IBM), U N I X , VMS and on different personal computers. Previous versions MCU-I.I and MCU-I.2 are already used for a long time, v e r i f i e d on t h o u s a n d s of practical problems and they have shown good fidelity. There exists some publications on MCU's comparison with experiments data and other codes calculations (main references a r e g i v e n in Maiorov and Yudkevich, 1988). 2.2
Code
architecture.
T h e u s a g e of M o n t e - C a r l o m e t h o d s for i n t e g r a l e q u a t i o n s o l u t i o n p e r m i t s to apply wide spectrum of algorithms based on different probabilities densities f u n c t i o n s for p a r t i c l e s trajectories and different tallies of v a l u e s of i n t e r e s t on t h e trajectories. The code structure permits to s i m u l a t e w a l k i n g of several various types of particles. The particle walking f r o m the c o d e p o i n t of view is a random sequence of abstract events: particle birth from external source, collision, flight between two p o i n t s of c o l l i s i o n , death, birth of a different type of secondary particles, w e i g h t c h a n g i n g etc. Concretization of p h y s i c a l m e a n i n g a n d e v e n t i n t e r p r e t a t i o n algorithms are r e a l i z e d by m o d u l e s f r o m c o d e l i b r a r y . M o d u l e of M C U is a s e t of s e g m e n t s , where each has its functional destination and interfaces strictly determined by the c o d e a r c h i t e c t u r e . W i t h r e g a r d to t h e i r destination MCU module types are c l a s s i f i e d so: C - c o n t r o l (governing), P - physical, Ggeometric, R - registration, S - s o u r c e s . A s a r u l e m o d u l e s of every type have alternative versions according to t h e i r d i f f e r e n t possible algorithms of m o d e l i n g . The code has the preprocessor governing the process of compilation of the e x e c u t a b l e p r o g r a m (EP) o p t i o n f r o m the m o d u l e s , chosen by the u s e r . T h e c o m p i l a t i o n means the realization of t h e f o r m u l a : EP=C+P+G+R+S where C,P,G,R,S modules. 2.3
Constants
are
the
alternative
(i)
names
of
corresponding
types
of
libraries.
All o v e r the w o r l d t h e r e e x i s t s a t e n d e n c y to u s e e v a l u a t e d d a t a files in Monte-Carlo calculations. But u n t i l n o w t h e s e f i l e s a r e not the u n i q u e and most preferable (for e x a m p l e , in u n r e s o l v e d resonances region} source of microscopic c r o s s s e c t i o n s for s u c h c a l c u l a t i o n s . W h e n we s t a r t e d o u r work on M C U pro.ject there existed no Soviet verified neutron data files, although well verified ABBN many group library ( A b a g y a n et al, 1981) for fast r e a c t o r s c a l c u l a t i o n s a l r e a d y e x i s t e d . So we c o n c e n t r a t e d our efforts on c r e a t i o n of c o n s t a n t s u p p o r t for l o w e n e r g y r e g i o n . At the same time A B B N l i b r a r y has b e e n a d a p t e d a n d u p d a t e d . N o w M C U c o n s t a n t s base consists of s e v e r a l l i b r a r i e s . E a c h l i b r a r y is o r i e n t a t e d towards a concrete class of p h y s i c a l p r o b l e m s or e n e r g e t i c region. ABBN/TR
library
is b a s e d
on ABBN
library.
t h a n its o r i g i n a l v e r s i o n . A B B N / T R h a s nuclei information had been updated.
ABBN/TR
information
is
more
about
89
comprehensive, nuclei.
Some
A p p t i c a t i o n to n u c t c a r r e a c t o r anat,,s,s
LIPAR
library
and
code
CROSS
connected
with
it
"_1_
(AbaKyan
et
al,
1984)
has
been developed for c r o s s s e c t i o n d e s c r i p t i o n and calculation in resolved resonances region. LIPAR contains Soviet evaluations based on information from different publications and neutron data files. Now there is information a b o u t 94 i s o t o p e s in LIPAR library. CROSS code allows to calculate c r o s s s e c t i o n s in e v e r y g i v e n e n e r g e t i c point. Oneand multilevel Breit-Wigner and Adler-Adler algorithms a l l o w to take into account resonances interference and cross section temperature dependence. CROSS time consumption p e r m i t s to u s e it in M C U e v e n w h e n m o d e l i n g c o l l i s i o n . It e r a s e s the q u e s t i o n a b o u t c r o s s s e c t i o n i n t e r p o l a t i o n precision. Thermal
constants
library
KORT
(Abagyan
and
Yudkevich,
1981,
1989)
describes n e u t r o n i n t e r a c t i o n w i t h n u c l e i for e n e r g i e s l e s s t h a n 5eV. Now K O R T c o n t a i n s d a t a a b o u t more, t h a n 90 m a t e r i a l s , including their phonon spectra, if t h e y are known. N E W R A S ( L i m a n and Maiorov, 1981) and STEN c o d e s h a v e b e e n w r i t t e n for scattering law and thermal cross sections calculations in incoherent Gauss approximation. Advanced STEN code a l g o r i t h m of s c a t t e r i n g law S(q,6 ) c a l c u l a t i o n is based on recurrent relationship S ( q + ~ q , 6 )=~S(q,6-6'
)S(~q,E'
)dE'
{2)
E L C O H R c o d e (Liman, 1982) is used for elastic coherent cross section calculation. TERMAC code (Gomin and Maiorov,1982) generates multi group work library TEPCON using KORT and LIPAR libraries. CROSS, NEWRAS and E L C O H R c o d e s a r e u s e d in T E R M A C for p o i n t c r o s s section calculation. The n u m b e r of groups, their boundaries and inner group spectrum may be arbitrary and p r o v i d e d by the user. NEDAM
library
describes
neutron-nuclei
interaction
in
energy
region
0 . 1 - 1 5 . 2 MeV. It is f o r m e d by C M O K c o d e ( M a r k o v s k y a n d Borisov, 1989) on the b a s e of e v a l u a t e d n e u t r o n d a t a files. N o w we u s e library formed from ENDF/B-4 and ENDL-83 files. NEDAM describes cross sections and other characteristics in f o r m s i m i l a r to i n i t i a l files. But NEDAM's format is o r i e n t e d on its u s a g e in M o n t e - C a r l o collision modules. 2.4
Physical
modules.
M C U i n c l u d e s a set of physical modules which use different constants libraries and different neutron-nuclei interaction models for different energy regions. Combined physical module (SOFISM) is u s e d in M C U for this purpose. Every SOFISM option consists of four alternative submodules t h r e e for n e u t r o n s a n d o n e for p h o t o n s c o l l i s i o n m o d e l i n g . F o r n e u t r o n s the whole energetic region is divided into three unintersected subregions (conditionally fast, resonance, thermal) with boundaries (Et,Ei_l). Modelin~
in e v e r y
subregion
is m a d e
by
special
submodules
of
Pi
(i=i,2,3)
type, w h i c h h a v e a l t e r n a t i v e versions. Submodules of different types may w o r k in the s a m e e n e r g e t i c s u b r e g i o n s . User himself defines a concrete submodule for his .iob (E i and Pi ) in i n p u t d a t a . P4 s u b m o d u l e is used for modeling is
the
photons'
collisions.
Operating
variant
of
physical
(OP)
module
OP=PI+Pz+Pg+P4 For most submodules continuous energy modeling is used approximation is p e r m i t t e d too. F o r r e s o n a n c e cross sections the f o l l o w i n ~ c o d e s a n d c o n s t a n t s libraries are applicable:
but group description CROSS with
214
E . A . Go~,.~ et aL
LIPAR l i b r a r y for c a l c u l a t i o n in e v e r y e n e r g y point w h e n it is necessary, BLOCK5 {author - Rogov A.D.) w i t h library, w h i c h c o n t a i n s resonance cross s e c t i o n s t a b u l a t e d in a set of e n e r g e t i c points (~6000 point for every r e a c t i o n ) for v a r i o u s t e m p e r a t u r e s . The most i m p o r t a n t s u b m o d u l e s of S O F I S M are d e s c r i b e d t h e m (MOFITTG2) uses g r o u p a p p r o x i m a t i o n . In b r a c k e t s given. FIMTOEN
(P3).
Energy
region
10-s-100 eV.
Main
below. Only one of submodule types are
features
are as follows.
The
K O R T d a t a l i b r a r y is used. Nuclei concentrations and temperatures may d e p e n d c o n t i n u o u s l y on s p a t i a l c o o r d i n a t e s , not I/V cross s e c t i o n s may be c a l c u l a t e d by e x t e r n a l codes (CROSS, for example) for every energy point during continous energy change modeling. The inelastic scattering is m o d e l i n g in ideal gas a p p r o x i m a t i o n or in i n c o h e r e n t Gauss approximation. The c o r r e l a t i o n b e t w e e n e n e r g y c h a n g e and s c a t t e r i n g angle c o s i n e is taken into account. Coherent elastic cross sections are used when it is necessary. MOFITTG2 TEPCON
(P9).
library
FIMBROEN
(P2).
Permits
to c a l c u l a t e
in e n e r g y
in multi
group
transport
Energy
region
0-10.5
r e g i o n 0-i eV on the base
of
approximation. MeV.
It
uses
ABBN/TR
and
LIPAR
libraries. All cross sections in every energy region, not including r e s o l v e d r e s o n a n c e s region, are taken from A B B N / T R library. Main features are d e s c r i b e d below. N u c l e i concentrations and temperatures may depend c o n t i n u o u s l y on s p a t i a l c o o r d i n a t e s . C r o s s s e c t i o n r e s o n a n c e s t r u c t u r e may be taken into account by any of three formalisms: Bondarenko's s e l f - s h i e l d i n g f - f a c t o r s , s u b g r o u p s , r e s o n a n c e cross s e c t i o n c a l c u l a t i o n in every energy point in resolved resonances region. Elastic scattering a n i s o t h r o p y is taken into a c c o u n t by the set of equiprobable bins of s c a t t e r i n g angle cosines in cms, scattering is sampled by kinematic formulas. FORTUN
(Pl),
( M a r k o v s k y and Borisov,
1989).
Energy
region
0.1-15.2
MeV.
Sub m o d u l e uses N E D A M library, M a i n features: scattering is sampled by kinematic formulas with continuous change of energy; energetic t r a n s m i s s i o n s for n o n e l a s t i c r e a c t i o n s are described by scattering laws g i v e n in initial n e u t r o n d a t a files. 2.5 G e o m e t r i c
modules.
MCU code p e r m i t s to m o d e l a r b i t r a r y c o m p l e x three d i m e n s i o n a l systems. This p o s s i b i l i t y is p r o v i d e d e d by a library of changeable (see architecture d e s c r i p t i o n ) g e o m e t r i c modules, w h i c h are d e s i g n e d for g e o m e t r y d e s c r i p t i o n and m o d e l i n g of trajectories between particle interactions. Geometric m o d u l e s l i b r a r y has 11 s p e c i a l i z e d g e o m e t r i c m o d u l e s and 2 u n i v e r s a l onesSCG ( a u t h o r - G u r e v i c h M.I.) and GDL (Ostashenko, 1989). Specialized m o d u l e s are d e s i g n e d to model s y s t e m s of l i m i t e d classes (reactors WER, R B M K etc) and more e f f e c t i v e for them, than universal ones. Widening of s p e c i a l i z e d m o d u l e s set is too c o n s u m p t i v e . So now we aim at universal m o d u l e s d e v e l o p m e n t . B o t h SCG and GDL are based on combinatorial method, w h i c h p e r m i t s to c o n s t r u c t d i f f e r e n t g e o m e t r i c s y s t e m s with the help of Boolean operations with primitives. They constitute a special library of p r i m i t i v e s . GDL m o d u l e uses h i e r a r c h i c a l p r i n c i p l e for c o m p l e x s y s t e m s w i t h one time d e s c r i p t i o n of r e p e a t i n g e l e m e n t s . This allows to s i m p l i f y system d e s c r i p t i o n and d i m i n i s h the required memory (with increasing of time c o n s u m p t i o n ) in c o m p a r i s o n w i t h SCG.
Application to nuclear reactor analysis
2.6
Control
T h i s type of by the m e a n s
21
modules. m o d u l e s p r o v i d e s the p o s s i b i l i t y to s o l v e of d i f f e r e n t alternative al~orithms.
transport
equation
M a i n c o m m o n s e g m e n t of d i f f e r e n t c o n t r o l m o d u l e s is a banker. The banker p e r m i t s n o t o n l y to place into memory and sample different types of particles, but to d i f f e r t h e m in p r i o r i t y in q u e u e s . It is convenient for nonanalo~ue m e t h o d s e.g. for R u s s i a n r o u l e t t e a n d for s p l i t t i n g . Critical problems for t h r e e d i m e n s i o n a l systems are solved by using the m e t h o d of the f o r c e d r e g u l a t i o n of the n u m b e r of p a r t i c l e s in generation. The efficiency of these algorithms are analyzed in (Zolotukhin and Maiorov, 1984). In particular there were derived formulae for the systematic bias of evaluation of functionals in such algorithms. Theoretical a n a l y s i s of b i a s m e a n i n g l e a d s to recommendation (Zolotukhin and Maiorov, 1984. p.31 ) that it is desirable to renormalize particle n u m b e r not in e v e r y g e n e r a t i o n , but a f t e r few o n e s m o d e l i n g . T h i s idea has b e e n indeperudently p r o p o s e d a n d realized in MONK code (Brissenden and Garlick, 1986) too. In a p p e n d i x we consider the connection between the ~'esults of the b o t h p u b l i c a t i o n s for the b i a s e s in the e s t i m a t i o n of Kef f, Besides traditional problems, MCU allows to solve critical cases for asymptotic s y s t e m s w h e n l e a k a g e in s o m e d i r e c t i o n s is ~iven by buckling v e c t o r B. T h e m e t h o d is b a s e d on b o u n d a r y c o n d i t i o n s for neutron flux on the f a c e s of c e l l s , w r i t t e n in the form: ¢(r+l where
i
k
are
k
,V)=¢(r,V)exp(
translation
vectors
of
!BI k ) the
(3)
lattice.
It is r e a l i z e d by the u s a g e of c o m p l e x w e i g h t s (Gelbard and Pego, 1977, Maiorov and Zolotukhin, 1984). F o r a s y m p t o t i c lattices with elementary cell in p a r a l l e l e p i p e d f o r m t h e s e w e i g h t s are: W:exp(~(B1n1+Bznz+Bgng) where number
n
is the k k and
number
of
intersections
(4)
) by
neutron
of
cell
Such approach has been used for criticality calculations material parameter. Exponential and critical experiments were s y s t e m s , w h i c h in XY p l a n e a r e d e s c r i b e d precisely and along t a k e n as h o m o g e n e o u s w i t h g i v e n m e a n i n g of l e a k a g e B z. 2.7
faces
with
Bu=(BIk).
Registration
with given a n a l y z e d for Z axis are
modules.
HCU calculates as a standard r e a c t i o n s r a t e s in s p a t i a l a n d fraction - linear functionals, few g r o u p s d i f f u s i v e constants
set the following parameters: flux and energetic regions defined by the user, migration lengths, neutron lifetime, s e t s of for r e a c t o r c e l l s , ~er etc. Some of this
values are calculated with weight of importance function. The last is calculated as sample meaning of number of descendants after some generations for n e u t r o n s w h i c h contribute to c a l c u l a t e d functional value. MCU registration modules architecture p e r m i t s to the u s e r h i m s e l f to widen calculated v a l u e s set. In t h i s c a s e the u s e r is a u t o m a t i c a l l y provided by m e m o r y for r e g i s t r a t i o n , automatical calculation of statistical errors, conservation of i n f o r m a t i o n on r u n s i n t e r r u p t i o n s , possibility to realize one's own estimators etc.
216
£. A. Go',,>~ et u/.
3.
MCU
VERIFICATION.
MCU has been verified with the help of special tests by comparing calculation results with those obtained by other codes and with e x p e r i m e n t a l d a t a . A s s e m b l i e s w i t h d i f f e r e n t fuel a n d d i f f e r e n t moderators in p a r t i c u l a r (ENDF-202), (BAW-1810,1984), (TIC, 1985) and many other a s s e m b l i e s h a v e b e e n e x a m i n e d . It has been compiled CLAD library. CLAD c o n t a i n s a s s e m b l i e s d e s c r i p t i o n s on M C U input language. There are ~150 assicnments of different experimental assemblies. They include data about reactors of practically all types which differ in fissile nuclei, m o d e r a t o r s a n d n e u t r o n s p e c t r u m . C r i t i c a l a s s e m b l i e s s u f f i c i e n t l y s i m p l e in ~ e o m e t r y a r e u s e d for v e r i f i c a t i o n of c o n s t a n t s . G e o m e t r y a n d m a t e r i a l s of a s s e m b l i e s a r e d e s c r i b e d as p r e c i s e l y as author's documents permit. The r e v i e w of c a l c u l a t i o n r e s u l t s c o m p a r i s o n w i t h e x p e r i m e n t a l d a t a is g i v e n in (Maiorov and Y u d k e v i c h , 1988). Briefly of
MCU
thermal
predicted
verification reactors
with
with
precision
of a s s e m b l i e s w i t h of fast reactors unsystematic.
results
may
be s u m m e d
U z35-
fuel
with
"0.5%
( experimental
up as
errors
Pu-fuel exceeds experimental is predicted better than
M C U is w i d e l y u s e d n o w in the U S S R S o v i e t r e a c t o r s of p r a c t i c a l l y all c a s k s etc.
follows:
moderators
one i%
H20, ~0.5%);
had
or
of d o u b l e
contained algorithm
heterogeneity
was
to 104 m i c r o e l e m e n t s of was developed (Bryzgalov
4.
MCU
analyzed fuel. et al,
for
been
done
sphere
For this 1989).
APPLICATION FOR DEVELOPMENT ENGINEERING METHODS.
is Keff
on i - 1 . 5 % , criticality and discrepancies are
parameters of and shipping
when
the
R e c e n t l y m a n y c a l c u l a t i o n s w e r e m a d e to a n a l y z e P W R / V V E R l a t t i c e s ( B e k k e r et al, 1990) and for t i g h t l a t t i c e s ( P s h e n i n et al, 1989). Effect
C
computed
for the calculation of types, in-pool storage
L a r g e s e r i e s of c a l c u l a t i o n s of v a p o u r e f f e c t of the C h e r n o b y l a c c i d e n t w e r e a n a l y z e d .
criticality D20
OF
elements
purpose
the
causes
with
Gd
which special
INGENIOUS
T h e i n g e n i o u s m e t h o d s of n e u t r o n fields computations used in the USSR p e r m i t to t a k e into account neutron spectrum change in heterogeneous reactors cells without complication of algorithms compared to homogenization traditional methods or Galanin-Feinberg method. A short r e v i e w of t h e s e m e t h o d s see in ( M a i o r o v a n d Y u d k e v i c h , 1988). M C U is used to a n a l y z e precision of engineering reactor codes, based on similar m e t h o d s , u s i n g s o m e k i n d of b a l a n c e e q u a t i o n s . One of the aims of MCU a p p l i c a t i o n is to v e r i f y a c c u r a c y of n e w e n g i n e e r i n g c o d e s c o m p a r i n g their r e s u l t s w i t h M o n t e - C a r l o c a l c u l a t i o n s of c o m p l e x a s s e m b l i e s . At the same time the c o e f f i c i e n t s of b a l a n c e e q u a t i o n s m a y be calculated during the M o n t e - C a r l o run w i t h o u t c a r r y i n g out of a n y a d d i t i o n a l e r r o r s in f i n d i n g of c e l l s p a r a m e t e r s , w h i c h are u s e d for e n g i n e e r i n g c o d e s . In t h i s r e p o r t our t a s k is c o n f i n e d to c o n s i d e r i n g s o m e i d e a s of t h i s a p p r o a c h only. The i d e a spectrum
of n e w m e t h o d s is b a s e d on the main hypothesis, that neutron in i d e n t i c a l c e l l s p l a c e d into d i f f e r e n t environment inside the
Application to nuclear reactor anatw~is
217
r e a c t o r m a y be w e l l a p p r o x i m a t e d by a s u p e r p o s i t i o n of p r e s c r i b e d n u m b e r of identical unknown arbitrary functions. Let us propose, that this number e q u a l s N * G + K . H e r e N is the n u m b e r of f a c e s of a c e l l (or their arbitrary elements in g e n e r a l c a s e ) , G - the n u m b e r of e n e r g y Rroups and K is an additional parameter. All the n u m b e r s G, K (and N in g e n e r a l c a s e a r e free parameters of a p p r o x i m a t i o n , w h i c h m a y be c h o s e n by user. From main hypothesis m a y be e a s i l y d e r i v e d the f o l l o w i n g relations: N
In----~ L
nn , ' , 4 - n
e'~nQ
n' =1
(5)
where
I n ,, n - are
G-dimensional
vectors
of
group
currents
and
fluxes
a v e r a g e d o v e r the c e l l f a c e s w i t h n u m b e r s n. Q and R are E-dimensional vectors. Their components are some arbitrary m o m e n t s of f i s s i o n density a n d a n y g i v e n k i n d of r e a c t o r r a t e R (for example, number of secondary fission neutrons) averaged over the volume of the cell. The weight f u n c t i o n s of t h e s e m o m e n t s m a y be c h o s e n by the user. Matrices
~
nn''
U
n'
V p W '
'
depend
on
the
chosen
set
of
functions
to
approximate
the flux i n s i d e the c e l l s . F r o m the m a i n h y p o t h e s i s f o l l o w s , t h a t all t h e s e matrices for identical ceils in different environment are invariable parameters of c e l l s a n d m a y be c a l c u l a t e d in principle before the full reactor calculation from relations (5). For c a l c u l a t i o n of m a t r i c e s it is sufficient to solve MaG*N+K problems every one on conditions when main hypothesis is valid. matrices m a y be f o u n d by m i n i m i z i n g functionals, constructed on the of r e l a t i o n s {5) connecting currents and fluxes calculated in problems.
model Then base model
All the m a t r i c e s for g r o u p s of identical cells may be obtained as a b y p r o d u c t of full s c a l e M o n t e - C a r l o r e a c t o r or m a c r o c e l l calculations. In r e a c t o r or m a c r o c e l l calculation different model problems are calculation of f l u x e s a n d c u r r e n t s on s u r f a c e s of i d e n t i c a l c e l l s p l a c e d in different environment when neutrons fly in cell with different energies (within different energy groups). The most attractive feature of Monte-Carlo m e t h o d s is the p o s s i b i l i t y to s o l v e a set of m o d e l p r o b l e m s d u r i n g o n e c o d e run ( M a i o r o v et el, 1969). T h e full set of r e a c t o r e q u a t i o n s m a y be w r i t t e n by a d d i t i o n to relations (5) ( w i t h n o w k n o w n c o e f f i c i e n t s ) of b o u n d a r y conditions. These are the r e l a t i o n s of c o n t i n u i t y of c u r r e n t s a n d f l u x e s on the i n n e r faces of the adjacent cells and the boundary conditions on e x t e r n a l s u r f a c e of reactor or m a c r o c e l l . The last necessary equation is c r i t i c a l i t y condition:
Q=
This
set
of
equations, general
eff
equations because
case.
It
is
are
of
(6)
+
only
their
formally
coefficients
important
to
note,
similar i
nn'
that
to u s u a l
etc
have
a crude
no
response physical
interpretation
matrices sense of
as r e s p o n s e m a t r i c e s m a y g i v e the w r o n g r e s u l t s ( M a i o r o v et el, 1969). r e a s o n is p r a c t i c a l impossibility to g u e s s t h e c o r r e c t n e u t r o n s p e c t r u m f a c e s of cells. It is possible to achieve much better accuracy
in
Lnn ' The on of
218
E.A.
calculation
of
L
from
GOMIN et al.
solutions
of
representative
set
of
model
Similar approach was described in ( M a i o r o v a n d Y u d k e v i c h , 1988) K=I a n d d e v e l o p e d in ( S h e v e l e v , 1989) for s p e c i a l c a s e L = N * G a n d c h o i c e of a p p r o x i m a t i o n for f i s s i o n d e n s i t y d i s t r i b u t i o n . Now many modifications w h i c h m a y be a n a l y z e d above.
of b a l a n c e e q u a t i o n s are using Monte-Carlo method
described in from viewpoint
problems.
a
for case special
literature, considered
The MCU code was used (Gomin and Maiorov, 1988, G o v o d k o v , 19891 to analyze the a c c u r a c y of the i n g e n i o u s reactor equations, based on the described a b o v e or s i m i l a r ideas. It w a s c o n s i d e r e d a p r o b l e m of c a l c u l a t i o n of flux inside macrocells of V V E R a n d R B M K r e a c t o r s , where large perturbation of s p e c t r u m t a k e s p l a c e n e a r c o n t r o l rods, w a t e r h o l e s a n d w a t e r g a p s outside the m a c r o c e l l s . F o r two d i m e n s i o n a l problems the balance equations were investigated , w h i c h a r e s i m i l a r to t r a d i t i o n a l finite difference diffusive equations, though having higher precision. For the
two d i m e n s i o n a l case when only
lattices two f i r s t
with symmetry cells azimuthal harmonics
it of
^
and
matrices
(for
K=I),
L
is p o s s i b l e to c o n s i d e r I and 0 are important n
have
°.
(Maiorov
and
the
i
property:
Yudkevich,
~
o° : Io-o
I"
Then
n
from
(5)
it
follows
1988):
J=~(O - Qf) (7)
I -J=~(0.-0)
where
J,O
- current
and
flux
averaged
along
all
the
face
of
cell:
J = N -I ~ In O matrices ~ and ~ and vectors a c e l l of e v e r y kind.
=
f are
N -I ! again
0n independent
from
environment
for
^
^
F o r the c a l c u l a t i o n of m a t r i c e s ; a n d 6 it is s u f f i c i e n t to problems (G - g r o u p s n u m b e r ) w h e n m a i n h y p o t h e s i s is v a l i d m a t r i c e s ; a n d 6 m a y be f o u n d by m i n i m i z i n g and fluxes calculated in m o d e l p r o b l e m s are number):
solve MaG and Q=0.
functionals in which connected (m model
model Then
currents problem
m
(8) ~=l Here
6 - are
the
=I
n
statistical
n
errors
of
m
HC
calculations.
m
T h e v e c t o r f is c a l c u l a t e d usually p r o b l e m (e.g. for a c e l l i n s i d e an by u s i n ~ the f i r s t r e l a t i o n (7).
after solving of an additional model i n f i n i t e l a t t i c e or i n s i d e a macrocell)
Application to nuclear reactor anab.~is
21 ~-;
The simplest balance equations for two d i m e n s i o n a l r e a c t o r w i t h cell h a v i n g N s u r f a c e s m a y be w r i t t e n on b a s e of (5) a f t e r e x c l u d i n ~ the equations of continuity of c u r r e n t s a n d f l u x e s on n e i g h b o r s i d e s of c e l l s as M
Jo =N-1 "~. (8-~+60')-'(P~lJn+Qnf.-O-'Jno o
(9)
Q0fo )
n=1
where number
^-
On I =
~
-I~-1 n - n
, and
n - relative
number
of
the
cell
around
cell
with
0.
The similar equations f o l l o w of the t h e o r y of t h e s u r f a c e h a r m o n i c s either (Laletin and Eishin, 1980), w h e r e it is r e c o m m e n d e d to u s e a n o t h e r p a i r of the a n g u l a r m o m e n t s to c a l c u l a t e c o e f f i c i e n t s of (9) (in particular the s e c o n d m o m e n t , so c a l l e d 'level' i n s t e a d of f l u x ) . A n o t h e r d i f f e r e n c e is a s p e c i a l set of m o d e l p r o b l e m s u s e d in ( L a l e t i n a n d E l s h i n , 1980). [n ( G o m i n and M a i o r o v , 1988, O o r o d k o v , 1989) the c o e f f i c i e n t s of equations (9) w e r e c a l c u l a t e d in p a r t i c u l a r by M C U code. The comparison was made b e t w e e n r e s u l t s of M C U and calculations based on equation (9). The a c c u r a c y w a s a c h i e v e d m u c h b e t t e r t h a n by t r a d i t i o n a l methods. The effects of a c c u r a c y of c a l c u l a t i o n of c o e f f i c i e n t s of e q u a t i o n (9) by simplified engineer methods and other approximations were investigated either. In particular the p o s s i b i l i t y to use for r e l a t i o n s (?) a n o t h e r p a i r of a n g u l a r m o m e n t s of flux, t h a n c u r r e n t a n d s c a l a r flux w a s c o n s i d e r e d ('levels'). T h e r e is a n o t h e r i n t e r e s t i n g application of i d e a s c o n s i d e r e d above. The a i m is to r e d u c e v a r i a n c e of c a l c u l a t i o n of flux distribution for large assemblies. In t h i s case the response matrices may be calculated as byproduct of full s c a l e M o n t e - C a r l o reactor calculation (see (8)). After t h a t flux d i s t r i b u t i o n m a y be f o u n d f r o m e q u a t i o n s s i m i l a r to (5,6) or m o r e s i m p l e o n e s , for e x a m p l e (9). S o l u t i o n of m o d e l p r o b l e m s a n d p r e p a r i n ~ of engineering codes constants are a c c e l e r a t e d when combined statistical and deterministic methods are used. In particular VEPS code (Gomin and Maiorov, 1984) for first collisions probabilities (FCP) calculation in complex two and three dimensional s y s t e m s a n d P E R S T (Gomin, 1984) code for integral transport equation solution in t h e r m a l r e g i o n by F C P m e t h o d h a v e b e e n w r i t t e n . F C P in V E P S c o d e are c a l c u l a t e d by n u m e r i c a l integration in unit cube of many dimensional integrals. T h e m e t h o d is a n a l o g o u s to MC, but it uses special nonrandom uniformly distributed meshes, in p a r t i c u l a r [Sobol, 1969) which provides I/N ~-6 law c o n v e r g e n c e of r e s u l t s ( h e r e 5 is a n y s m a l l n u m b e r ) . On the o t h e r h a n d s i m i l a r i t y of the a l ~ o r i t h m to MC p e r m i t s to use for FCP calculations all M C U ~ e o m e t r i c m o d u l e s , achieving enormous gain in time consumption w i t h o u t loss in ~ e n e r a l i t y . P E R S T c o d e m a y be u s e d in p r i n c i p l e as part of M C U in t h e r m a l r e ~ i o n for d e c r e a s e of t i m e e x p e n d i t u r e .
APPENDIX
In [ Z o l o t u k h i n formulae for
and the
Maiorov, 1984) and ( B r i s s e n d e n and Garlick, 1986) the K biases are worked out for continuous and eff finite-discrete models accordingly. It is i n t e r e s t i n ~ to note, t h a t in b o t h publications authors proposed similar initial relations as a base for estimation of bias. It is f o r m u l a (171 in notations of (Brissenden and Garlick, 1986)
Dn + ]
:
FD n /(t,D n )-e n /M
(A.ll
~.-()
E. A . G o s , u x et al.
Here
D n is a r o w v e c t o r h a v i n ~ as its c o m p o n e n t s the expected number of neutrons absorbed in the ~th b o x of p h a s e s p a c e in g e n e r a t i o n n. M - n u m b e r of p a r t i c l e s in o n e ~ e n e r a t l o n , F - the k e r n e l of t r a n s p o r t e q u a t i o n , u - a v e c t o r ot n u m b e r ot f i s s i o n n e u t r o n s per one absorption, e the vector, n which depends on some variance - covarlance matrix, defined in a r t i c l e . The equivalent relation is the f o u n d a t i o n for conclusions a b o u t the b i a s of K in ( Z o l o t u k h i n and Maiorov, 1984 p.51 (13.3)I It was confirmed in elf ' ' ° (Brissenden and Garlick, 1986) by the series solution for Keff-b~as of {A.I ). B u t t h e n f o r b i a s b a s e d on
followin~ estimators
the authors (A.I).
In ( Z o l o t u k h i n and Maiorov. 1 9 8 4 ) it w a s perturbation theory. Comparing for l a r g e
proposed
the
different
used the formalism of n (A.]) with unperturbed
final
the usual equation
D= FD/I,,D~ it
is
easy
to
derive
6Kef r/Ke where
F+
is
the
rf
=
a vector
¢a.~)
simple
-
expression
for
bias
of
Kerr:
(A.3)
M1(r+e)/(r+D)
of
importance
function
and
e
- a value
of e
n
for
larRe
n.
(F+e)
The expected value m a y be c o m p u t e d by M o n t e - C a r l o c o d e if it permits to e s t i m a t e t h e n u m b e r of d e s c e n d a n t s for neutrons in future generations. This possibility Js i n c o r p o r a t e d in MCU to calculate functionals with weight of an i m p o r t a n c e function. Brissenden v a l u e of estimators usage,
and
Garlick
(r+e). of
however
showed,
that
value
of
bias
of
Kef f,
which
may be calculated as the difference of Keff. Their formula is more convenient its
equivalence
with
(A.3)
was
not
proved
by
depends
two to
on
variance practical
ourself.
REFERENCES
Aba~yan
L.P. , B a z a z y a n t s
Constants
for
Reactors
N.O. , N i k o l a e v and
Shieldin~
M.N.,
Tsybulya
Calculation.
A.M.
(1981)
Energoizdat.
Group
Moscow.
In
Russian. Aba~yan
L.P. , T e b i n
V.V..
library
for
sections
Kurchatov Aba~yan
cross
Institute L.P.,
Tekhniki
Series:
M.S.
M.S.
{1984)
calculation
IAE-4009/5.
Yudkevich
{VANT)
¥udkevich
In
(1981)
Yadernye
in
CROSS
- code
resolved
Kort
library.
Constanty.
V
Voprosz
l(40),p.
39.
Atomnoy In
Pshenin V. , Lazarenko A., Spectral Codes Verification Kernenergie, Voi.33, No 8,
Alekseev N., Vaneev using Numerical and pp 335-341.
Brissenden
R.J.,
Biases
A.R.
constants re~ion.
Russian.
B e k k e r R., Results of Benchmarks.
Garlick
and
resonances
{1986).
in
the
estimation
Nauki
i
Russian. Yu. (1990). Experimental of
K
elf
and
Applicution to nuclear r e a c t o r a n a h ~ i s
its
error
b~
Monte-Car[o
methods.
Ann.
Nuc[.
221
Energy,
Vol .
IJ,
No.Z,
Brvzgalov V. [. , G o m t n E . A . , G u r e v t c h M.t. , Kaminsky A.5. , 5ubbotln 5.5, Tebln V.V. (1989) C a l c u [ a t Lon and Experimenta[ Research of Double Heterogeneity Fuel Elements of H['GR. ~ANI'. 5erles Eiztka i lekhnlka Yadernvkh
Reaktorov
(FiTYR),
E N D F - ' ~ O Z . !.'HERMAI. REA('TUR
V
Z,pp.
Computational
Methods
Nuclear
in
Gomin
E.A.
Calcu[at
~ons
( 1985 )
In
June
1974
Monte-Carlo Proceedln~s
I'omputatiorl of of the Topical
Engineering-.
%.2,
p.7-9.
Directional Meetings on
USA,
Virginia.
E.A.,
P['obabii i t ies
L.V.
(1985)
Ca[culation~
[AE-I~O7/5.
FiTYR,
in
VEPS
Code
Three
V 2,np
( 1988 ) Monte-Carlo Reactors Calcu|ation
24-33.
In
G,~revich
M.I.
11990)
Systems.
M.[.
Po[ycells.
First
First
La[etin
Experiment
N.I. , for
institute
Program Complex for 'Fhermalization Re~ion.
for
First
Collisions
Systems.
kurchatov
Elshin
Preclsion VANT
UC-78
Urania
A.V.
Algorithms for 155. in Russian.
Probabilities in
[n
Calculations
Russian.
Probabillties
[AE-5173/5.
Benchmark.
Heterogeneous
Method Usage for Modern Methods
IAE-SIZZ/5.
Collisions
institute
BAW-1810,DOE/ET/34212-41.
Equations
Institute
Modern Mesh Kner~. 68, p.
Collisions
Kurchatov
(19901
kurchatov
Crltlcal
Funct ionals
Russian.
S.S.(1989) Tile Limit Error of of Neutuon Flux i n RBMK. A t o m n .
Gurevich
Neutrons kurchatov
Dimensiona[
Gorodko~ CalcuLation
Compact
Thermal
Systems.
In R u s s i a n .
Gom~n E.A. 0 Maiorov L.V. Analysis of Heterogeneous
Series
['or
L.V. ( I!)82 ) T E R M A C the GrOHp Cross Sections in V 5(27),pp 70-73. In Russian.
Maiorov
Institute
Code
Dimensional
Russian.
Gomin E.A. , Maiorov Calculation of Neutron \'ANT Series FiTYR,
Gomln
PERST
Three
in
[AE-42US/5.
and
Russlan.
23-25.
A~)F [l
for
[n
CuMP[LAII()N.
BENCHMARK
PeSo R. , (1979) Coet ficients.
Ge[bard E.M. , D~t'L'usion
44-50.
Calculations
for
Russian.
Gadoiinia:
Nuclear
Model
Development
1984.
(1980) Reactors.
The
Derivation
Kurchatov
of
Institute
Finite-difference [AE-3280/5.
[n
R,~ssian. Liman G.F. { 1982) Codes CALCST and Elastic Coherent Scatterin~ Cross
ELCOHR for Calculation Section for Thermal
Series:
In
Liman
FiTYR,
V
5(27),
O.F. , Maiorov
L.V.
DD.
74-76.
(1981)
The
of Total and Neutrons. VANT.
Russian. NEWRAS
Pro,~ram
for
Computation
of
the
222
E.A.
Slow
Re,it i-on
32-~{,d.
[n
Scatt~rin~
o["
Kurchatov
Institute
Maiorov
L.V. ,
[hermai
Evaluation
by
Markovskv
,
Borisov
H( 2L > ,
pp.
A.A.
Files.
V.D..
S.V.
(1989)
CaLculation.
Ya.V.(1989). [.M.
(t969) Nauka.
CoLlected
gbsorbers
Private Many Moscow.
Works.
in
Moscow.
[n
Module Micro
,
In
Moscow.
[n
Kuznetsov
Vienna,
Parameters Russian.
FORTUN-85
for Based
Russian.
E.A.
A.S., Lazarenko Rod Lattices. CRP
V~.ERs.
Russian.
Constants
[AE-4774/5.
H.S.
Romanov in F u e l
Physical
Crittcailty
with
[nstitute
,
kapitonova
A.P. (t989). on Safe Core
IAEA.
Graphics to R e p r e s e n t Reactor" Series: FiTY[4, V 2 , T>P. Z 6 - 3 0
-\ssembiies in
Russian.
communi,~atton. Dimensional
[n
(1984)
Constants
Physical
Method
Yudkevich
Machine VANT,
Neutron
Reactors
( 1989 )
V.E., \lekseev N.I., Results for Gadolinium
Ost~ishenko Monte-Carlo
~ [~To; Clusters.
Ener~oatomizdat.
kurchatov
Sidorenko
Iruknano~ u.}a. Ca/cutatton of
Ener~oatomizdat.
Monte-Ca~'Io
Biirn.'-tbte
TIC
([.9£S)
Method.
wi th
Function~.
M.S.
Calculations.
Mana.~ement
Sobo[
Ftl'~R
Russian.
V.G. ( 1984 )
Data
V.V.
In
Monte-CaL'lo
by
A.V., Osipov CaLcuLation
Sheve[ev
kANT.Series:
A.D. , Method to
Zolotukhin
D.i'. ,
EvaLuated
Pshenin
1AE-420~/5. Yudkevich
Calculations on
Frank-kamenetsky Monte--Carlo
Reactors
L.V. ,
Maiorov
Sections.
Russian.
Maiorov L.V.. The A:>piication
in
(_'ross
O,,~.x ,'t ~d.
~uadrature
Formulas
Russian. Akademiai
Ktado,
%o[.
I,~udapest.
and
Haar's
in