Some aspects of Monte Carlo method application to nuclear reactors analysis

Some aspects of Monte Carlo method application to nuclear reactors analysis

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,...

635KB Sizes 0 Downloads 19 Views

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