If in a system the sources are isotropic and only the isotropic component of the flux is required, it is possible to write the transport equation in a simpler form. Let us treat the monoenergetic time-independent case.

Let us consider a neutron starting from r’ in the direction ft. The probability that it reaches the point r is given by

= exp -[t(x)]

where

t(x)= X,(r — x’fl) dx’

JO

is the optical distance (probability of having a reaction along the path x), and X, is the total cross-section. The neutron flux in r is given by

ф(г, П) = exp [-r(x)]Q(r-xft, ft) dx

where Q is the source representing fission or scattering. If the sources can be assumed to be independent of ft the equation can be integrated over the angle. Analytical expressions are only possible in particular cases, but in general one obtains an expression of the form

(4.53)

In terms of collision density фХ, we have:

(4.54)

The kernel K(r — r’|) is then the probability of neutron penetration from point r’ to point r. K(]r — r’|)2,(r) is the probability of a neutron starting from point r’ to suffer its first collision in r.

In the case of a multi-group calculation eqn. (4.54) holds for every group, and the source Q,(r) is given for each group / by

Q,(r) = Sj(r) + d>,(r)2s(r) = Si(r)+ ф,(г)Х,(г)с

where Xs IX, = c is the number of secondary neutrons per collision and S,(r) represents the fission and slowing-down source and is therefore dependent on the flux of all groups. It is always possible to calculate numerically the kernels K(r — r’|).

These methods are usually called “collision probability” methods because the number of neutrons colliding in one region is calculated from the number of neutrons produced or scattered in other regions multiplied by the probability that they will make their next collision in the region considered. This implies a discretization of the reactor space in a certain number of regions.

In general collision probabilities depend not only on geometry and composition of the regions considered but also on the spatial flux distribution within each region.

What makes these methods attractive is the fact that, in practice, if the regions are small compared with the mean free path, the collision probabilities are rather independent from the flux distribution so that they can be calculated assuming flat flux in each region. The above-mentioned assumption of isotropic sources refers both to fission and scattered neutrons.

This last condition is not always satisfied, but is sufficiently accurate for graphite because of its relatively high atomic weight. Besides in the thermal energy range crystal bindings and thermal motion of the graphite atoms tend to give a more isotropic scattering. In a more general energy-dependent case we have

This expression is sometimes known as Peierl’s equation. Here we can distinguish a transport kernel

K(r-r’,E)

and a scattering kernel

2,(r£’->£)•

Energy is then discretized in G groups and space in N regions obtaining the following expression:

= 2 PS f Qi*v, + v, 2 slrvl

i=l L к=1 J

where 2,8 = total cross-section region j group g,

Vj = volume region j,

Фf = flux in group g region j

2ЇГ8 = scattering kernel discretized for region i, group к to g,

Q8 = source region і group g,

Pfi = collision probability: probability that a neutron of group g originating in region і makes its next collision in region j.

We use here a superscript to indicate the energy group in order not to generate confusion with the various other subscripts.

In this equation the transport kernel K(r-r’,E) has been discretized in the collision probabilities P%. This implies the assumption that the collision probabilities are not dependent on the flux distribution within each region and can be calculated with the flat flux approximation.

One can then obtain for P% the following expression:

the exponential factor 5.‘(s)ds

0

gives the probability that the neutron does not suffer collisions along the path s. From (4.57) follows the reciprocity relation

SiaV, PS = SSVjPS. (4.58)

In an infinite system

N

2 Pn = і

j-i

because a neutron born in region і will certainly make its next collision system regions, if leakage is excluded. The probabilities P8 are normally calculated by numerical integration of eqn. (4.57). Monte Carlo methods can be used for special cases. Once these collision probabilities are known it is possible to calculate the fluxes ФіК by means of eqn. (4.56). Usually iterative methods are used, e. g. inserting a trial
solution for the fluxes on the right-hand side of eqn. (4.56) and obtaining an improved solution on the left-hand side (power iteration method).

A typical example of this method is given by the THERMOS"21 code where only the thermal energy range is considered in the geometry of a reactor cell. A unit slowing — down source is assumed over the cell.

Another example is given by the PIJ code"31 which allows an r — в geometry upon which can be superimposed rather general distributions of cylindrical “rods” themselves subdivided in r and в. A typical HTR block representation for the PIJ code is given in Fig.

4.2. "4)

Collision probability methods are also used to calculate control rods (e. g. MINOTAUR code.’15’

Some extension of the collision probability method has been developed in order to treat anisotropic sources and scattering,"4’ but this is not very important in graphite­moderated systems. Approximate methods are sometimes used to calculate collision probabilities. For cylindrical geometry the widely used Bonalumi method approximates the fate of neutrons beyond the second boundary crossing."6’

Sometimes for spectrum calculations the reactor volume is discretized in very big regions (MULTICELL method,"7’ and the flat flux approximation is no longer valid: corrections based upon diffusion theory have been developed for this purpose."8’

Collision probability methods are best suited to cases where the required accuracy is not very high and approximate methods can be used to calculate collision probabilities. If a high accuracy is demanded an S„ representation may be more advantageous."4’19’

The averaging of cross-section for the reflector regions can be done with any of the codes for homogeneous spectrum calculations which have been analysed.

A fission source is usually specified in order to provide neutrons to start the calculation. The error due to the introduction of this source can usually be neglected. The effect of the neutrons leaking from the core in the reflector can be taken into account using negative bucklings. Unfortunately, these bucklings are very important in the neutron balance of the reflector and usually they are not known to a very high accuracy. This can strongly influence the neutron spectrum, and negative fluxes can be easily obtained, so that it is usually better to calculate with zero bucklings. These calculations have to be repeated for the various parts of the reflectors, having different composition, density and temperature.

This rather crude method is often sufficiently accurate, but the problems posed by these calculations show that in this case the assumption made in performing space — independent spectrum calculation is no longer valid, since the spectrum, in this case, is not only determined by the characteristics of the region under consideration, but also greatly influenced by the neighbouring regions.

The thermal diffusion length varies from 10 to 30 cm depending on the core loading, and in this range very important spectral changes occur. An example is shown in Fig. 8.2<20) where the thermal spectrum is plotted at different positions. This is a rather extreme example because its fuel loading corresponds to the ratio C/U = 2500 which is very high for an HTR power station, but may occur in experimental reactors. If a high accuracy is needed a one-dimensional diffusion calculation with a high number of groups (e. g. 26 in the Fort St. Vrain case) can be performed including core and reflector. The space-dependent spectrum obtained in this way can be used for the condensation of cross-sections for the whole reflector or for some part of it. Similar methods can be used at other interfaces between regions with strong spectral differences.

In highly simplified calculations it is possible to avoid the explicit treatment of the reflector, increasing the size of the core by an amount (reflector saving) defined in such a way as to give to the bare reactor the same kcя of the actual reflected reactor. This

treatment neglects all spectral changes between core and reflector and therefore gives usually a very inaccurate power distribution.

In the general case Д/с is function of temperature, control-rod position, and any variable describing the external disturbance which can be cause of an accident (with the exception of the loss of coolant, all accidents appear in Д/с).

The temperature feedback is usually separated in fuel and moderator effect, so that

Д/с = І " Q/(T’) dT’ + " ©„(Т’ЫТ’ + Дкс+Дк. (12.16)

JT/o Jt-c

where 0/(T) = fuel temperature coefficient,

©m(T) = moderator temperature coefficient,

Д/Сс = reactivity given by the control rods,

Д/са = reactivity disturbancy causing the accident.

The integrals appear because the temperature coefficients are usually temperature- dependent.

Equation (12.16) is generally used for the zero-dimensional calculations. If space — dependence is taken into account one can still use this equation for each region in which the reactor is subdivided, or, in a more general way, it is possible to consider the temperature-dependence of the parameters of eqn. (12.1). In this case Д/с will not explicitly appear in the equations and the poisons introduced by the control system appear in 2,. In any case it is necessary to calculate the reactor temperature in zero or more space dimensions.

The explicit solution of eqns. (12.1) [or its simplified versions (12.2) or (12.7)] with
(12.16) and the heat-transfer equations is a formidable, if not impossible task. On the other hand, a numerical solution can be easily programmed on a digital computer.

The concept of the high-temperature gas-cooled reactor with a prismatic core construction was first seriously investigated at the Atomic Energy Research Establish­ment, Harwell, early in 1956. At about the same time work started in Western Germany on the so-called pebble-bed version of the system in which the core consisted of a randomly packed bed of spherical fuel elements. The Harwell work resulted eventually in the construction of the 20-megawatt reactor experiment of the OECD High Temperature Reactor Project (DRAGON) at Winfrith, England; the reactor reaching its nominal power in April 1966. Also resulting indirectly from the Harwell concept, the Peach Bottom Reactor of the Philadelphia Electric Company, developed and designed by General Atomic of San Diego, achieved full power operation in May 1967. These two experimental high-temperature reactors were quickly followed by a third when in February 1968, the AVR pebble-bed reactor, built by the Brown Boveri/Krupp Reaktorbau GmbH, for Arbeitsgemeinschaft Versuchsreaktor utility group, came into operation at the Kernforschungsanlage, Jiilich.

Over the subsequent years, these three reactors have been operated with consider­able success, demonstrating in a convincing manner the feasibility of the high- temperature reactor and confirming the favourable characteristics claimed for it. Now a much larger version of the prismatic core HTR has been brought into operation at Fort Saint Vrain near Denver, Colorado. This plant, like the Peach Bottom Reactor, is a product of the General Atomic Company and it is intended to generate 300 megawatts of electric power on the network of the Public Service Company of Colorado. A large pebble-bed reactor of similar power output, the Thorium Hoch-Temperatur Reaktor (THTR), is under construction at Schmehausen in West Germany.

The special features of the high-temperature reactor are its use of helium as a coolant and the manner in which graphite is employed not only as a moderator for the neutrons, but also as the structural material of the reactor assembly, troublesome metal cladding and tubing being entirely eliminated from the core. Fissile and fertile materials are present as oxide or carbide kernels of the coated particles. These particles, up to about 1 mm in diameter with pyrocarbon and silicon carbide coating layers to contain and prevent the escape of fission products, are incorporated into the graphite structure of the fuel elements by consolidating them in a matrix of graphite powder and resin moulded under pressure and then carburized at high temperature to form robust fuel bodies.

The excellent characteristics of the HTR results from this combination of helium, graphite and coated particle fuel. Very high temperatures in excess of 800°C for the mixed core outlet gas and even above 1000°C for hot channels have been sustained over

extended periods of operation in the experimental reactors. This makes it possible not only to run the reactor in conjunction with the most advanced steam power plant, but also to contemplate the direct use of the helium in gas turbines and to consider the application of the very hot helium to carry out important industrial processes, giving the HTR great development potential.

A further merit of the HTR stems from the refractory nature of the graphite and its enormous thermal capacity coupled with negative reactivity temperature coefficients which are readily achievable. This makes the reactor particularly insensitive to power excursion and loss of coolant accidents. The HTR is consequently a very safe type of reactor. It is also clean and accessible, a consequence of the remarkable degree of fission product retention of the consolidated coated particle fuel and the absence of activated material arising from the corrosion and erosion of metal cladding and other primary circuit materials that can be a source of trouble in reactors having liquid coolants.

The physics of the high-temperature reactor is of special interest and significance because of its versatility with regard to fuel cycles. Coated particle fuels can use any of the fissile materials and either uranium-238 or thorium as fertile material. According to the level of fissile enrichment of the fuel kernels, the size and volume loading of the coated particles, the ratios of carbon to fissile and fertile atoms and the degree of heterogeneity of the core are variable to an extent that is impossible in other reactors. Absence of metal cladding greatly improves the neutron economy in comparison with other types of reactor, a factor which helps in the achievement of high conversion ratios. Furthermore, it has been demonstrated beyond all doubt that coated particle fuel can sustain a degree of burn-up which is quite unmatched in any other reactor system. The result of all these factors is a reactor which can be readily adapted to the changing economics of the uranium market. In the same reactor, for example, one can change from a low enriched uranium cycle, with a very high burn-up, to a thorium cycle, possibly sacrificing burn-up for a high conversion ratio to match changing cir­cumstances with respect to the cost and availability of uranium ore, separative work and reprocessing capacity.

This book is concerned with the physics aspects of the high-temperature reactor, a subject on which there has so far been the lack of a systematic and comprehensive text despite the very considerable volume of reports that have been written on the various aspects of the subject. Several years ago the Dragon Project decided to sponsor such a work with the co-operation of the Commission of the European Communities (Directorate-General for industrial and Technological Affairs) and invited Dr. Massimo to undertake the authorship.

Dr. Massimo is well qualified to write on the subject having been a member of the Dragon Project Staff from 1959 to 1963 during which time he was involved in the early physics assessments of the 20-MW Reactor Experiment. He was also attached by the Dragon Project for a short period to the General Atomic Company where he participated in the work on the Peach Bottom Reactor. Subsequently, as a member of the staff of Brown Boveri/Krupp, he was engaged in the design study of the THTR reactor.

In his book, the author describes the state of the art, at the time of writing, in the field of HTR reactor physics. Following a survey of the basic theory of modern reactor physics, the book covers the methods and computer codes employed in high —

temperature reactor calculations by the Dragon Project, and by the various groups working in the subject in European and American national research centres and in industry.

High-temperature reactor physics, having reached a state of maturity, is nevertheless still evolving. Any book, such as this, therefore, cannot maintain indefinitely an up-to — the-minute presentation of the subject. However, the book is aimed primarily at the students of HTR physics who are preparing to enter the field as well as technologists of other disciplines who are working on the system and in this role it should remain a useful and relevant text book for many years.

Winfrith, Dorset

July 1975

L. R. Shepherd

Chief Executive, OECD High Temperature Reactor Project (DRAGON)

In a system without leakage or absorption, fission neutrons are slowed down by successive collisions with moderator atoms until their energy reaches an equilibrium with the thermal motion of the moderator. This means that in such an ideal system the low-energy part of the neutron spectrum is a Maxwellian distribution corresponding to the moderator temperature. In real systems a certain number of neutrons is absorbed or leaks out of the reactor before reaching this equilibrium.

While the slowing-down of fast neutrons is not influenced by the thermal motion of the moderator, below a certain neutron energy this thermal motion is no longer negligible in comparison with the neutron energy. In this lower energy range (thermal energy range) neutrons can gain energy in a collision (up-scattering).

As a Maxwellian distribution goes only asymptotically to zero toward higher energies, the upper limit of the thermal range cannot be exactly defined but it must be set at an energy above which up-scattering becomes negligible.

Usually energies of 1 eV or lower are used, but in HTRs, because of the high temperature, limits between 2 and 4 eV are more appropriate. In case of negligible neutron losses the scattering properties of the moderator are unimportant because the neutron spectrum is always a Maxwellian distribution. The higher the losses, the more the neutron spectrum depends on the competition between loss and moderation, and hence on the scattering properties of the moderator. At high neutron energies all moderators can be treated as monoatomic gases, but at low neutron energy the molecular and crystal forces binding the moderator atoms start playing a role in the moderation process.

Burn-up is usually expressed in megawatt-days per tonne (MWd/t) of fuel, referred to the thermal power the fuel has produced in the reactor.

In high-temperature reactors where the fuel is intimately mixed with a part of the moderator this unit can give rise to some confusion. It is then necessary to specify it as megawatt-days per tonne of heavy metal, excluding in this way the weight of carbon and other light elements present in the fuel.

Because of the high burn-up attained by HTR fuel, sometimes the unit megawatt — days per kilogramme (or gigawatt-days per tonne) is used. Other units often used in HTR calculations are FIFA (Fissions per Initial Fissile Atom) and FIMA (Fissions per Initial Metal Atom), sometime expressed in percent. The burn-up measured in MWd/t or in fima differs only by a constant factor. Using the rough approximation that the fission of 1 g of fuel per day corresponds to 1 MW, fima 1 = 105 MWd/t (or more exactly fima 1 = 0.93 x 105 MWd/t).

The relation between fima and fifa depends upon the ratio of fissile to fertile concentration in the fresh fuel. In order to calculate the relation between the energy produced and the number of fissions one has to consider the following distribution of
the fission energy:

Kinetic energy of fission fragments Kinetic energy of fission neutrons

Prompt у Delayed у

P

Neutrinos

to these values one should add the y-rays produced by n, у reactions in various reactor materials.

Except for the neutrinos one can consider that all this energy is converted into thermal energy within reactor core, reflector and shielding. While the kinetic energy of the fission fragments is completely converted into thermal energy within the coated particles, other items, like the kinetic energy of the fission neutrons, are more uniformly spread over the core. One can in general consider that 92% of the power is produced in the fuel and 8% in the moderator. From the above data follows that 1 joule = 3.2 x 1010 fissions (this quantity is weakly dependent on core size).

4.1. The neutron transport equation

The behaviour of neutrons in any medium can be described by the transport
equation, first used by Boltzmann for the description of the behaviour of gas molecules.

Neutrons can indeed be thought as forming a kind of a gas whose treatment is simplified by the fact that neutrons do not interact with each other. The general time, position and angle-dependent form of the transport equation is

mean:

neutron velocity corresponding to energy E, neutron angular density,

total neutron cross-sections (generally function of r and E), neutron source, space coordinate,

unit vector in the direction of the neutron motion,

energy,

time,

scattering cross-section from E’, ft’ into E, ft.

This equation represents simply a neutron balance in a volume element dV for the neutrons having energy between E and E + dE and flight direction in the solid angle dfl around ft.

The first term on the right side is the leakage out of dV, the second term represents the loss due to absorption and scattering, the third gives the sources due to scattering from other directions and energies and the fourth is the source term (including fission and external sources). We do not repeat here the well-known difinitions of neutron angular density, cross-sections, etc. (see refs. 1 and 2).

Defining the neutron angular flux

Ф(г, E, ft, t) = vN(r, E, ft, t),

eqn. (4.1) becomes

— дф(г’ f; a—~ = — ftVt/r(r, £, ft, t) — Х, ф(г, E, ft, О v at

+ J Xs (E’ -* E, ft’ -* Г1)ф(г, E’, ft’, 1) dCT + S(r, E, ft, t).

(4.2)

The boundary conditions for this equation are continuity at the interface between two media;

at the outer free surface ф{г, Е,£1, t) = 0 for directions entering the system; ф(г, £, ft, t) = 0 for energies greater than the maximum source or fission neutron energy. For reactors this limit corresponds to 10-15 MeV.

The transport equation (4.2) without external sources takes the form:
j дф(г, Е,(1, t) = _nv^( E () — х, ф(г, E, ft, t)

V ot

+ JXs (E’ -*E, ft’-* П)ф(г, E’, ft’, t) dE’ d(l’

+ J Х/(Е’)ф(г, E’, ft’, t)v(E’) dE’ dft’, (4.3)

here the last term represents the fission source which is supposed to be isotropic;

x(E) is the fission spectrum, probability that a fission neutron is generated with energy E

£ x(E)dE=l;

v(E) is the average number of neutrons generated by a fission induced by a neutron having energy E.

If the leakage and absorption terms are equal to the production terms (fission and slowing down) the time derivative of the left-hand side of eqn. (4.3) vanishes.

In this case the reactor is said to be critical and the neutron population does not change with time.

For a given material composition all terms on the right-hand side of eqn. (4.3) are fixed except the first one. This means that criticality can only be reached for a particular value of the term

ftVt//(r, E, ft)

which represents the neutron leakage out of the systems and depends on the reactor geometry. In this way the critical geometry is defined. For a given reactor criticality can be obtained changing X, by means of insertion or extraction of neutron absorbers. Static calculations of non-critical systems are usually performed artificially dividing the term v(E) of eqn. (4.3) by a factor keff so that the equation is always satisfied without considering the time dependence of the neutron population.

Again we take the limit a0-> 1 of the first integral of eqn. (7.13) while the second and third integrals are treated with the NR approximation obtaining

7.3. Computer methods for resonance calculations in heterogeneous geometry

Various computer codes have been developed for detailed resonance-absorption calculations. Equation (7.14) is numerically solved in the General Atomic ZUT code<5) which is incorporated in the GGC-IV and GGC-V’6,7’ codes for neutron spectra calculation. The ZUT code deals only with resolved resonances while the unresolved energy range is covered by the TUZ code’51 (see § 7.12). The more sophisticated General Atomic GAROL code’8’ considers also two space regions but takes into account overlapping and mutual shadowing of resonances solving the following equations in any specified number of energy points:

VoX, o<f>o = (l-Po)Vo2j^-f ‘фо(Е’) dE’

і = 1 1 QJ і J E С

M ж t <*E/a( / T? t

+ P, V, 2 t—^ Ф.(Я’) dE’ + V0(l — Po)Qo(E) + ViP, Q,(E)

(7.17)

and another just like it but with subscripts 1 and 0 reversed. Nki is the atomic density of isotope і in region k, asi(E’) is the microscopic scattering cross-section of isotope i, Qk(E) is the source in region к, M is the total number of isotopes considered.

The collision probability method is also used in the General Atomic MICROX code<9) which solves the neutron slowing-down equations on a detailed energy grid for a two-region lattice cell.

In both GAROL and MICROX the detailed shape of the Doppler broadened cross-sections is stored on data tapes in an ultrafine grid of discrete energy mesh points covering the resolved resonance energy range in 9000-15,000 points. Because of storage limitation this treatment has to be limited to the most important nuclides.

Other codes are not limited to two regions. Even if physically only two regions are usually present in HTRs the use of more regions allows a better approximation because collision probabilities are usually calculated assuming flat sources in each region.

For example, the Ispra collision probability code PETARD’10’ can be used in typical cases with 2200 groups, 15 regions, 2 compositions and 50 nuclides. The French RACOLE code allows up to 50 concentrical regions.’10 Four regions are used in the UKAEA SDR code with a total of 120,000 energy points.“2)

As we have seen, the aim of all optimizations is to minimize the cost of the energy produced by the plant. This cost includes investment cost, fuel costs, operating expenses, taxes, etc.

The reactor physicist is usually only concerned with the fuel-cycle-cost calculation, although the influence of core parameters on other investment costs must be taken into account in the optimization (e. g. the influence of the core power density on the pressure vessel size).

Expenditures and returns are strongly time dependent. For each segment of fuel expenses start long before this fuel produces energy, and continue after discharge from the reactor. Revenues are obtained from the sale of energy and in some cases from the sale of discharged fuel. Besides these costs usually change during the reactor lifetime.

The fuel cycle cost is the cost at which energy must be sold in order to recover all the expenses related to the fuel consumption and ownership. The correct value of this cost makes the total present worth of the incomes equal to the total present worth of the expenditures.

The present worth P of a payment R is defined as

P = R(l + iTy (10.18)

where і is the interest rate and у is the number of years after which the payment is made. This formula is a recognition of the fact that an interest has to be paid on money.

The factor

F = (l + iTv (10.19)

is known as “present worth factor” or “discount factor”. A discussion of the present worth method can be found in ref. 15. If interests are compounded more frequently than annually, expressions (10.18) and (10.19) are still valid if і is the interest per time period and у is the number of periods.

In order to compare fuel cycles involving slight time displacements (e. g. the difference between a cooling time of 100 or 120 days for spent fuel) the method of continuous discounting can be used."6’ The time period is divided in n intervals and the

Fig. 10.2. Time sequence of typical fuel segment (from GA-A10593).

fictitious annual interest g is introduced so that

•■—КГ

The new interest rate per period is not і In because continuous compounding would then make the annual interest greater than i. The expression for continuous compound­ing is then obtained if n -»30:

or

(1 + /) [ — e 8t (10.20)

so that

g = In (1 + і) (10.21)

and expression (10.19) is substituted by

F= e

With the present worth method all expenditures and revenues are referred to an arbitrary reference time, which is often taken to be the reactor start-up.

A difficulty may be posed by the definition of the interest rate to be used in this method. This interest must take into account taxes and the fraction of the investment represented by stocks and by bonds, each with their own required rates of return.

Another powerful tool for solving particular reactor problems is given by the Monte Carlo method, which replaces a deterministic problem by an analogous one consisting of a game to be played many times. These numerical calculations which make use of statistical variables can be used to simulate the physical phenomena actually happening in a reactor where the neutrons have a stochastic behaviour. A high number of neutron

histories are simulated on the computer, every interesting event is recorded and statistical averages are calculated. Random numbers are used to decide at which place collisions occur and of which type they are (absorption, scattering, etc.). Although only a series of numbers can be defined as random, one speaks usually of random numbers. They are usually generated by arithmetic procedures in the range between 0 and 1.

Apart from the actual simulation of the physical phenomena the Monte Carlo method can be used to obtain numerical solutions of integral or differential equations in many space dimensions. Seen from this second point of view the Monte Carlo methods can be thought of as a numerical tool for solving the Boltzmann transport equation.

An example of the stochastic methods for a non-stochastic problem can be given by the evaluation of the area under the curve of Fig. 4.3. N random points can be chosen within the dotted rectangle and if M is the number of points falling below the curve the area can be obtained after having calculated the ratio M/N. The accuracy of the calculation will of course depend on number of points chosen. Equally well known is the example of the dice. The problem consists of determining whether a dice presents any irregularity in its geometry or mass distribution. The normal solution consists in performing a number of physical measurements on the dice, the Monte Carlo solution consist of playing many times with it in order to determine whether there is any irregularity in the distribution of the results of the game.

The result of a Monte Carlo calculation is, of course, a statistical variable and its accuracy is roughly inversely proportional to the square root of the number of events contributing to it (e. g. number of neutrons undergoing a collision in a certain region). This means that the number of variables which can be calculated by this method should be relatively small in order not to increase enormously the computational effort required (e. g. it is relatively easy to calculate the total neutron flux in a region, but it becomes very difficult to have its accurate space and energy distribution within this region). The square root law means also that in order to increase the accuracy of a calculation by a factor of 10 one must increase the number of samples by a factor of 100 so that it can be very expansive to increase the accuracy beyond a certain limit.

On the other hand, with Monte Carlo methods the transition to a higher number of dimensions implies only a linear increase in computing effort. In the above-mentioned example of the calculation of an area under a curve a numerical integration would have been much more efficient, but for a similar problem in a five — or six-dimensional system the Monte Carlo method would have been more advantageous.

This means that it is relatively easy to calculate with Monte Carlo methods complicated three-dimensional geometries for which S„ codes could prove impossible for most computers now available. Usually the particles which are followed by these computer calculations represent more than one neutron. The number of neutrons represented is called weight of the particle. In this way if, for example, a collision has a 20% absorption and 80% scattering probability the weight after collision is reduced by the factor 0.8, while the rest is recorded as absorbed.

Monte Carlo methods can be very inefficient in some cases and it is then necessary to modify the game in order to reduce the variance or error. For example, in many cases there are neutrons that because of their position, energy and direction can contribute more than others to the final answers. A typical example is given by the study of neutron attenuation through a shield in which, without special variance reduction techniques, an enormous number of neutrons should be followed across the shield in order to have a few neutrons at the end of it. One of the methods of variance reduction consists of changing the weight of the particles, so that regions of high interest are studied by means of many particles having little weight while regions of lower interest are studied by few particles having high weight. Weights can be changed as particles cross given space or energy boundaries. An increase of weight is performed by a game called “Russian Roulette” in which the particle has a probability p of surviving and 1 — p of being killed. The surviving particles are followed with a weight increased by the factor 1 Ip as imposed by the conservation of the number of neutrons. The opposite game, consisting in a decrease in weight and subsequent increase of the number of particles, is called splitting. The Russian Roulette can also be used in cases in which it is not of interest to follow particles with very low weight. Monte Carlo methods are also used successfully to calculate the effect of localized perturbations in a reactor. In that case the unperturbed case is first calculated and then weight standards are chosen so as to enhance greatly the diffusion of particles into the perturbing region. In the field of high-temperature-reactor design Monte Carlo methods, while unsuitable for complete reactor calculations, can be used to treat a number of detailed problems. Among these we can mention resonance absorption calculations, neutron streaming through big holes, shielding. In this last field the Monte Carlo method is very often the only possibility of obtaining sufficiently accurate results. The criticism which is often made to Monte Carlo calculations is that they provide results, but give very little information on the physical phenomena determining them. (For details see ref. 22.)