These effects have been thoroughly analysed by Parks et al. o m calculating and measuring spectra, average fission and absorption cross-sections, 239Pu to 235U fission ratio, effective multiplication constants, and temperature coefficient of reactivity for bare, homogeneous, critical graphite-235U cubes.

These calculations have been performed at different temperatures for various C/U ratios, using both crystal and free gas models.

The difference in neutron spectra between the two models is shown in Figs. 6.1 and 6.2 for 300 К and 1200 K.

Fig. 6.2. Thermal neutron spectra at 1200 К for six different C/235U atom ratios. (From
Parks etal, Nucl. Sci. Fngng»1)

The spectra are normalized to give the same number of total absorptions below

1.0 eV. At room temperature (300 K) the effect of chemical binding is important for all loadings for which there is a significant number of thermal neutrons. At very high 235U loadings spectral effects of binding are less pronounced because of the lack of thermal neutrons.

At very low 235U loading again the effect of binding must disappear because the spectrum will show a Maxwellian distribution, independently from the scattering model.

At 1200К the effects of chemical binding are negligible.

The differences in spectrum averaged cross-sections is shown in Table 6.2. One can see the effect of the hardening of the spectrum due to chemical binding effects (the Boron cross-section averaged over a Maxwellian distribution at 300 К is of 669 barns.)

Table 6.3 shows the percentage difference in effective multiplication constants.

100[(/Ccrystai kgas )/к gas] *

For heavy 235U loadings the thermal spectrum is unimportant and the difference is small. At the limit of zero absorption this difference must also vanish. A peak in the percentage difference appears around the atom ratio C/235U = 10,000.

Table 6.4 gives the effect of crystal binding on the temperature coefficient

1 dk к дТ’

In this case of bare cubes the change in leakage accounts for most of the negative temperature coefficient. The calculated temperature coefficients are considerably smaller for the crystal model than for the free gas model.

A very important problem in core design consists in obtaining the desired power distribution, and in keeping it sufficiently constant during burn-up. Usually the limits imposed on reactor materials refer to temperature and neutron dose on the different fuel and moderator components. The reactor designer must try to operate every part of the core as near as possible to these limits in order to obtain the maximum performance. As it has been discussed in the first chapter of this book, in HTRs the problem of thermal limitations is not so stringent as in other reactor types.

If we assume for simplicity a well-defined temperature limitation, the best perfor­mance will be obtained with an isothermalization of the core. This would require a radially flat power distribution. In the axial direction a constant temperature of the
coated particles at the centre of the fuel elements would require an exponential variation of the power distribution. This point requires a deeper investigation.

The gas temperature at the height 2 of a coolant channel is given by

Tg(z)=Tg0 + -^-f Q(z’)dz’ (10.4)

cpm Jo

where Tg0 = gas inlet temperature,

cp = specific heat of coolant gas, m = coolant mass flow per coolant channel,

Q(z) = power per unit length in the direction corresponding to the coolant channel under consideration.

The fuel temperature is given by

Tf(z)=Tg(z)+Q(z)R (10.5)

where R is the heat resistance between the coolant and the fuel.

Equation (10.5) is valid either for the centre of the fuel or for its surface, according to the value chosen for R.

From (10.4) and (10.5) we have

If we require an axially constant fuel temperature we impose

(10.7)

from eqns. (10.6) and (10.7) we have

fT + lT^Q =° dz Rcpm

whose solution is

Q = Qo e-,,/Rvi,)z; (10.8)

an axially constant fuel temperature requires an exponential axial power distribution given by (10.8). Changing the value R one can obtain either a constant surface temperature or a constant fuel central temperature.

Whatever the chosen power and temperature distribution, the reactor designer is confronted with the problem of finding a fissile and fertile material distribution which satisfies the requirement at the beginning of the core life (or at the beginning of a reloading period) and as a function of burn-up. In the case of fresh fuel it is rather simple to fulfil such a requirement iterating on the fissile concentration by means of a diffusion code. Usually a continuous variation of the fuel composition is not possible and a certain number of axial and radial zones must be chosen. This implies step variations in concentration and hence in temperature. In the case of radial power flattening for example, the natural decrease of flux and power towards the core edge must be compensated by step variations in the fissile concentration.

A procedure to determine the fissile loading distribution iterating with the help of a diffusion code is described in ref. 1.

It is more difficult to obtain a power distribution which remains invariant with burn-up. In this case besides varying the fissile one has also to vary fertile material concentration. This can be done on a trial-and-error basis, but a mathematical treatment is possible and has been used in the General Atomic GASP code.<2) The power density as function of time and space is given by

P(r, t) = 2/(r, t)<f>{r, t)

and if the power density is constant with burn-up

where superscript Fe denotes fertile material, Ft fissile material, В bred fuel. The first term on the right-hand side represents the increase of 2/ due to breeding. (This equation neglects time delay due to decay of 233Pa or 239Np and losses due to neutron absorption in those isotopes.) The second term represents the decrease of 2/ due to neutron capture and fission in the fissile material.

Equation (10.11) is rearranged as

From (10.10) and (10.12) we obtain

ФІГ, 0)

■ФІГ, 0)J A(r, t’) dt’

2/F’ (r, t) = 2,F’ (r, 0) [ 1 — ф (r, 0) Jo A (r, t’) dt ’ J

It is further assumed that A{r, t) is separable in space and time

A{r, t) = A(r,0)A(t),

,(r s фіг, 0)

I ) r-t J

l-S(r)J A(t)dt

2/F‘ (r, t) = 2/’ (r, 0) [l — Sir) Jo‘ A it) dt]

where

S(r) = ф(г,0)А(г,0).

If the quantity S(r) in eqns. (10.16) and (10.17) is made independent of r the flux and power distributions will be space-invariant with burn-up. Therefore the GASP code makes the quantity Sir) in every region equal to the core average value S. This criterion
must be satisfied together with the condition that the initial power distribution corresponds to the specifications. The fulfilment of this second condition is program­med in GASP using the method of ref. 1.

These considerations did not take into account recharge operations, so that this method can only be applied to the interval between two subsequent reloadings.

The isothermalization of a reactor is a costly operation. Radial flattening requires a higher fissile concentration near the radial core boundary, and an axially exponential power distribution requires a shift of the fissile material toward the gas inlet side. Both operations increase the neutron leakage. Besides they require the manufacturing of many types of fuel elements and a complicated fuel-management scheme.

This means that the reactor designer should weigh the advantages of isothermalizing the core against its disadvantages. Besides one must consider that it is impossible to isothermalize both the central fuel temperature and the surface termperature, and that limits are imposed for both temperatures. An axially exponential power distribution will also increase the maximum temperature gradient across the fuel which would be lower in the case of axially constant power.

In the case of high-temperature reactors although an isothermalization is not strictly necessary, some form of power flattening is required, especially in cases where a very high gas outlet temperature is required. It must be noted that apart from power flattening, a temperature flattening can be sometimes obtained by varying the coolant flow through the fuel elements (gagging). Particular care should be taken to analyse the influence of control rods and burnable poisons on the power distribution.(3)

Another way of eliminating the angular dependence of the Boltzmann equation consists in expanding the angular dependent variables in terms of spherical harmonics:

ф(г, Е,П) = І Z, Фіт(г< E)Pim(Cl) (4.9)

|=0 m = —I

and in an isotropic medium

ї,(Е’-»Е;П’-»П) = 2 s,(E’-»E)P,(cos0o) (4.10)

I

considering that in an isotropic medium %s is only dependent on the angle во between П and O’.

Substituting (4.9) and (4.10) in (4.4), multiplying each term of the resulting equation by Pim (fl) , conjugate complex of Pim (fl), and integrating over В one gets an equation for each component фіт from which the angular dependence has disappeared (see ref. 1, p. 226). Here again the value chosen for l determines the order of the approximation, which is usually indicated by the symbol P,. Usually an odd value of l is chosen because these approximations give normally better results.

The space dependence in the resulting set of differential equations can be treated either with numerical or with analytical methods. Analytical solutions are possible for various simple geometries (e. g. slab, sphere, cylinder).

A case frequent in reactor applications is that of cylindrical geometry with many concentric regions of different properties. In this case the components фіт satisfying the differential equations are a combination of Bessel functions whose coefficients are determined by the boundary conditions: continuity at each inner boundary and reflective condition or vacuum condition at the outer surface.

These types of analytical solutions can lead to very fast computer calculations, but especially in the higher-order approximations numerical difficulties are easily encoun­tered. Furthermore, only a few simple geometries can be treated with these methods (see ref. 8). For this reason analytical solutions of the Pi equations are nowadays

seldom used for reactor calculations, while on the contrary spherical harmonic expansions are often used in connection with numerical methods (e. g. expansion of anisotropic scattering in S„ methods). The first terms of the spherical harmonics expansion of the neutron flux have a well-known physical meaning.

Defining the total neutron flux

it can be easily verified (see ref. 1) that the components of the vector J are:

Mr, E) = y

jx(r, E) = ^у2~112(фи- — фи), (4.14)

/Дг, Е) = —у і2-,,2(<К-і+</’п).

Also the first terms of the expansion (4.10) have a physical significance.

Integrating (4.10) over all directions ft we have

Xs(E’^E) = 4tts0(E’^E), (4.15)

in general defining

XAE’^E)=^Ts,(E’^E)

eqn. (4.10) becomes

2s(E’-> E, fio) = І 2.1 (£’ -» E)P, (цо)

i-o 47Г

where цо = cos во — If we define the average cosine of the scattering angle

Usually the source is isotropic so that

Some isotopes have important resonances at low energy. Typical examples are 239Pu, 240Pu, MIPu and the fission products 135Xe, 149Sm, 103Rh. These resonances occur in an energy range where up-scattering starts to play an important role, so that many of the approximate methods we have seen in this chapter cannot be used.

Doppler broadening of such resonances has usually a very small effect. In most cases these resonances are treated with multi-group methods using a fine group structure in the resonance region. Geometrical effects are treated either with self-shielding factors in zero-dimensional codes or with space dependent multi-group transport calculations (see Chapter 8). A proper treatment of these low-lying resonances is important especially in the calculation of temperature coefficients.

The effect of moderator temperature on the thermal neutron spectrum has been discussed in Chapter 7.

The thermal neutron spectrum in HTRs is usually so different from a Maxwellian distribution that it is no longer possible to speak of a “neutron temperature”. The model developed by Parks to treat the crystal binding in graphite is normally used at low temperature, while the free gas model is adequate above 1000-1200K.

The use of the free gas model at low temperature can lead to errors in the temperature coefficient of the order of 20-50%. The thermal neutron spectrum acts on the factors tj, / and P, h of eqn. (11.2).

The term (Мт))(дт)ІдТ) depends upon the nuclear properties of the fuel and its sign is determined by the percentage of the various fissile isotopes present in the reactor.

The term (1 lf)(df/dT) represents the competition between thermal absorption in fuel and in other materials. Its contribution is negative if the absorption in the fissile isotopes drops with energy more quickly than the absorption in the remaining isotopes. It is then clear that fissile isotopes having pronounced thermal resonances like 239Pu (0.296 eV) are very important in determining this contribution. The increase in temperature shifts the spectrum towards these resonances, giving a negative contribution to the temperature coefficient. This trend is then reversed at very high temperatures. Doppler broadening of such resonances is generally a second-order effect. Uranium-235 has always a negative contribution. Uranium-233 has a positive contribution at very high temperatures, while the sign of its contribution at operating temperature depends on the core design and composition.

Among the fission products a considerable negative contribution at very high temperature is given by l03Rh because of its resonance at 1.26 eV.

A positive contribution is usually given by the isotopes: U9Sm, l5lSm, 155Gd, 157Gd, ll3Cd, 153Eu and especially 135Xe whose absorption cross-section drops strongly with increasing energy.

The effect of 135Xe on temperature coefficient is particularly important. As a consequence of load following the l35Xe concentration changes strongly with time and the moderator temperature coefficient tends to become positive, or less negative, when Xe reaches its maximum concentration.

Let us consider the diffusion equation in its multi-group form (4.40) without external sources,

where N is the total number of groups.

Finite difference methods can be used for the solution of this equation in one-, two- or three-dimensional geometries.

It is not the intention of this book to give a full treatment of all possible numerical methods available to solve the multi-group diffusion equation. We will limit ourselves to the description of a few typical cases (for further details see refs. 4 and 5). It is convenient to split eqn. (5.2) into the two following equations:

DiV^i — ї’іф, +Q, =0, (5.3)

Q, = £ (2^+77-AW. £/")</>-. (5-4)

n = 1 К eff

where (5.3) has taken the form of the one-group equation and Q is the source which in the multi-group cases relates group і to all other groups n.

In the finite difference method the reactor has to be subdivided in a high number of regions whose boundaries are defined by constant values of one of the coordinates (e. g. in three-dimensional XYZ geometry those boundaries are planes of constant X, constant У or constant Z, while in two-dimensional RZ geometry the boundaries are lines of constant R and constant Z). The parameters of the diffusion equation have to be constant within each region, but can change from one region to another. The neutron flux is calculated at mesh points which are normally put at the edge of each region. It is also possible to place the mesh point at the centre of each region (see ref. 7). We will here consider the first possibility only. Let us first consider a one-dimensional slab geometry, with a total of M mesh points.

Equation (5.3) takes the form

D, ^+Q, = 0.

ax

Л„

Let us integrate this equation over the interval corresponding to mesh point xm:

‘(Q, — Х,,ф,) dx = 0;

we obtain for each mesh point x„

The expression will be formally more complicated if the constants D and X, vary from one mesh to the next. These M equations can be written in matrix notation. For the ith energy group we have:

Ai Фі = qi (i group index)

where A is a matrix and Ф and q two vectors. In the monoenergetical case the problem is limited to the solution of this system of equations.

In a similar way it is possible to treat a two-dimensional case. Let us consider a Cartesian XY geometry (Fig. 5.1).

Avt

Equation (5.3) takes the form

again we integrate over the area corresponding to the mesh (s, t) +1 xs, Уі +1 у»

—7— r*+————

2 2 obtaining the five-point difference equation:

Also this equation can be written in matrix notation

А, Фі = q„ (5.7)

It is this case necessary to order the two-dimensional array of the фі(х!, у,) into a one­dimensional vector Ф, sweeping one row after the other.

To complete the system (5.7) one has to add to eqns. (5.5) or (5.6) the boundary conditions for the points lying at the boundaries. These boundaries can be either external or internal (the so-called “non-diffusion regions” representing, for example, control rods).

This is one method of obtaining the finite difference equations, but it is also possible to use Taylor series and variational methods.

It is possible to use a nine-point difference equation in order to express V2ф but in this case difficulties arise when the coefficient D changes from mesh to mesh.

What has been written for XY geometry can be easily extended to RZ or R6 geometries.

The Rd geometry presents some particular problems. The centre point has usually to be omitted because it represents all meshes Є for R = 0. Besides a cyclic boundary condition may appear

ф (r, в = Ф(г, в)

where n is a positive integer.

A geometry peculiar to HTRs is the 60-degree uniform triangular mesh which is very convenient for reactors with hexagonal fuel elements. In this geometry it is difficult to represent a circular outer reflector boundary which has to be approximated, with the consequence of azimuthal variation of the reflector thickness. This can be compensated to some degree by artificially changing the density in some of the reflector blocks. The following seven-point difference equation can be obtained (see ref. 13).

(6D +1h 2lt)ф0 — 2 <t»D = -2h 2Qo■

i = ]

As in the previous cases this formula takes a formally more complicated form if the point considered lies at the corner between three blocks with different composition.

Up to now we have only dealt with eqn. (5.3) which treats only one energy group. Usually the so-called source iteration method is used to solve multi-group problems.

The source Qi consists of two parts, a fission source Qfl and a scattering source Qsi

Qi = Qfi + Qsi,

Qfi = г — 2 Xif’n IfnQn,

*v n = l

Qsi = X

П 7й /

The matrix eqn. (5.7) can then also be written in the form:

АФ = S<h +7- F<£, к

here A is a supermatrix including the At of eqn. (5.7) for all group i, and the same applies to the matrices S and F. SO represents the scattering source (slowing-down and up-scattering) and (l/k)FO represents the fission source. The source iteration tech­nique consists of an iteration on the fission source (l/k)FO.

Defining the matrix

M = A — S

we have

МФ = r FS. к

The iterative procedure used to solve this equation consists of choosing a first guess for the fluxes Ф at all groups and mesh points and for the eigenvalue k. Let us denote these first guesses by Ф0) and kw. Having these values we can calculate a source (l/k<1))FO<I) and iterate according to:

МФ<3) = pL, ^Ф15’0 (5.8)

where s is the iteration index. Those iterations are called outer iterations.

If up-scattering is present the inversion of matrix M implicit in this process may present some difficulties and the process can be improved subdividing S into upper and lower triangular matrices U and L.(,6)The matrix U contains the up-scatter components only and L contains the down-scatter components only. We have then the following iteration procedure:

АФ",= L<J><5)+ 17Ф<5_|> + ТТГГТ) РФ’**11

к

M = A-L

мФи,= Lfo<s“,)+p^Fo(s-,).

The lower triangular matrix L contains only terms related to energy groups whose flux has already been calculated, i. e. they are known terms, so that the inversion of matrix M can be performed with well-established procedures (see ref. 17). This matrix has dimensions of the order of the number of groups times the number of mesh points. In the case of one-dimensional problems the inversion of matrix M is performed by direct methods (Gauss elimination), while for two — or three-dimensional problems iterative methods are necessary (e. g. over-relaxation method). These iterations are usually called inner iterations. The outer iterations are performed according to eqns. (5.8) or (5.8a) and the eigenvalue kis) is given by the Rayleigh quotient (see ref. 17)

с, [РФ^РФ^1

[РФ<!)]ТРФ<5_1)

where [РФ<5)]Т is the row vector transpose of FФ<5). This expression is physically justified by the fact that к is the ratio between the fission sources of two subsequent generations of neutrons. In a non-discretized form it would indeed take the form:

The iteration process is continued until convergence is achieved. The value of к given by the above expression is the overall ксЯ of the reactor. In analogy to this overall value it is also possible to obtain a local value k„ for each mesh point p,

{РФ(% ■„

p {F ф1*-”}

where (РФ), is the component of the vector РФ corresponding to mesh point p.

Usually the overall kctt converges more rapidly than the flux shape which is related to the local value kp.

Because of this reason two convergence criteria are usually used. The calculation is stopped when

1,(0 _ r.(«-l)

•V eff »v eff

kcS(s)

and

I,(S) _ K(S)

Jv max л mm

where fcmax and kmh are respectively the maximum and minimum values of the set kpu).

Of course, the number of iterations required depends on the guesses which have been used. For this reason a considerable amount of computer time can be saved if good guesses can be found, for example performing in sequence series of calculations where only minor changes occur in the reactor (e. g. successive time steps of a burn-up calculation, or power distributions at different control-rod positions).

Sometimes, instead of iterating on kc«, criticality searches are performed in which other parameters like fuel or poison concentration, or reactor dimensions are varied until a specified kc« is obtained.

An important example is the search on control-rod insertion, as this is closely analogous to the way a reactor is operated/6* These options are available in many diffusion or transport theory codes. When the flux distribution in an undercritical system with a source present is required (e. g. reactor start-up) a so-called source calculation is performed in which no eigenvalue is computed, but the method of solution is rather similar to the one described above.

From any set of reactor core calculations (e. g. spectrum calculations followed by space-dependent diffusion or transport calculations) one gets the energy and space distribution of the neutron flux together with the eigenvalue fccff. This is all what a reactor designer needs. On the other hand, a better understanding of what is happening in the reactor can be obtained if more detailed information is given about where, at what energy and in which isotopes neutrons have been absorbed in the reactor.

Once the energy and space distribution of the neutron flux is known a neutron balance can be easily obtained. These balances are usually calculated for any given core region or for the whole core (e. g. in the case of a zero-dimensional calculation). Starting from a fission neutron it is possible to calculate the probability this neutron has of leaking out of the system (or of the considered region) or of being absorbed in the various isotopes. A typical example is given in Table 8.4 for a thorium fuel cycle.

From a balance of this type it is possible to evaluate the relative importance of the isotopes present in the system. It is also possible to calculate the conversion factor as a ratio between the production rate and loss rate of fissile isotopes (in the production rate one must consider the decay of isotopes like 233Pa and 239Np).

In the case of few group diffusion or transport calculations where only macroscopic cross-sections are present in the computer code, a distinction between the different isotopes is not possible. The balance is then limited to the calculation of the total neutron production, absorption and leakage over the whole energy range and for each individual energy group. This regional and group-wise balance calculations are also the best way to obtain the region and energy-dependent bucklings from a space-dependent diffusion calculation. For each group and region fission sources, slowing-down sources and total absorption are known. The difference gives ИВ2ф and knowing the average

flux it is possible to obtain B2 which can be used for buckling iterations in spectrum calculation codes.

Another type of balance consists of resolving the eigenvalue ксЯ in its various components. This break down of ксЯ was in the past used for calculational purposes. The modern development of computer methods has made this unnecessary, but still desirable for a better understanding of core physics. A first break down of ксЯ separates the multiplication factor for an infinite medium from the fast and thermal non­leakage probabilities Pf and Pth,

fee* = k«P/P, h. (8.24)

A further breakdown of k* is given by the classical four-factor formula

к» = 67]fp. (8.25)

This is historically based on a two-group formalism for natural or low enriched uranium reactors with a well thermalized neutron spectrum. The following assumptions are more or less explicitly included in this break down.

(a) No up-scattering is present between thermal and fast group.

(b) All fissions in the fissile elements occur in the thermal group, the only fast fissions taking place in 238U.

(c) All fast absorptions are resonance absorptions in the fertile material.

In this way one can define:

e fast-fission factor (238U fissions),

7] number of fission neutrons per absorption in fissile materials,

/ ratio of thermal absorptions in the fissile material over total thermal absorptions (thermal utilization), p resonance escape probability.

As none of the above assumptions is strictly valid for HTRs the four-factor formula can only be used with some slight change in the definitions. Up-scattering from thermal to fast group must be neglected (some up-scattering is present in HTRs even with a group boundary around 2 eV). This implies also a change in the definition of the slowing-down
cross-section 21^2 from fast to thermal group, in order not to modify the neutron balance. The exact two group equation is (index 1 means fast, 2 means thermal)

DtB,2(f>i + 2ai</>i + 2 і-*гф і = — (v2/i</>i + Г’Х/г^г) + 2г^і</>2,

•Ceff

D2B2 Ф2 + 2a2</>2 + 2 г—1</>2 = 2l-.2</> 1-

We impose now two fictitious transfer cross-sections 2*^2 and 2Li

2г— 1 = 0,

2i-.2 ф 1 = 2і-.г ф 1 — 2г-.і</> 1.

The first of these equations represents the condition of no up-scattering, while the second defines 2І-г in such a way as not to alter the number of neutrons exchanged between the two groups.

As ф 1 and ф2 are known it is possible to obtain 2L2>

21,2 = 2 „2-22-.^

Ф1

and the two group equations take the following form:

(DiBi2 + 2o1 + 2i—г)</>і ~ T” (V2/1 ф 1 + v2/2</>2),

•Ceff

(D2B2 + 2аг)</>2 = 2l-.2</>l-

The above defined parameters can now be calculated according to the following expressions:

this p could be called “fast absorption escape probability” because it does not only consider resonance absorptions in fertile materials, but all epithermal absorptions.

/ can only be calculated if microscopic cross-sections are present in the computer code

where Nu’233 means atomic concentrations of 233U, etc.

Other parameters which are sometimes calculated to give a better appreciation of the situation to the designer are

D2

la2

P.

T 2., + SU

M2 = L2 + t

L2 is related to the total length travelled by the thermal neutrons and т is related to the total length travelled by fast neutrons. This can be seen from eqn. (4.37) considering that S„i + 21^2 is the removal cross-section from group 1.

The quantity M2 is therefore related to the total length travelled by a neutron in the course of its life.

M is called migration length and gives a feeling of the distance at which a perturbation introduced in a certain point of the core can still be felt. The quantity Vr is called slowing-down length.

If one performs a condensation to one energy group only, the diffusion equation takes the form

(8.39)

In this case

so that the diffusion equation becomes

V2<f> +

In the case of a critical uniform reactor with energy-independent boundary conditions the flux satisfies the wave equation

V2<f> +В2ф =0

so that the eigenvalue В2 (geometrical buckling) must be

is called material buckling. In the case of a just critical reactor BM2 = B2.

Some computer codes for spectrum calculation and cross-section averaging, after having performed the cross-section calculation in the desired number of groups, repeat the condensation to two groups in order to produce the following parameters: ксП, Pf,

Pth, є, tj, /, p, /с», L t, M2, Bm. This is especially useful in the calculation of temperature coefficients, where the break down of kcS in its components is very often used.

As explained in § 12.1 135Xe and 149Sm are the only fission products who need to be treated explicitly in core dynamics.

The fission product chain involving l35Xe is:

Fission

/3 6.7 h

Xe •

9.2h

Cs

For 135I and,35Xe eqn. (9.1) takes the following form (the,35I absorption cross-section is negligibly small):

where Xf = macroscopic fission cross-section,

o-xe = l35Xe microscopic capture cross-section,

I(r, t) = l35I concentration,

Xe(r, t) = l35Xe concentration,

>1, yXe = fission yields of l35I and l35Xe,

Ai, Axe = decay constants of l35l and 135Xe.

Samarium-149 is a stable isotope part of a rather complicated chain (see §9.1).

Most of these isotopes have a direct fission yield which for simplicity is not shown here. This chain is very important for burn-up calculations but does not require full treatment in the dynamic calculations. The concentrations of l48MPm, ‘“’Pm and l48Sm can be considered constant during the transient, so that the chain can be broken at this point. The half-life of l49Nd is short (~ 1.7 h) in comparison with that of l49Pm (—47 h) so that also l49Nd does not need to be explicitly treated: its yield will be simply added to the one of l49Pm. The,49Sm equations take then the form

-—РдГ~Г-’ f~ = yu, Pm2/<Mr, t) — A’"pml49Pm(r, () + (Я«“рП11®мРт(г, t)4>(r, t)

Ol

+ o-‘«Pml48Pm<f>(r, f), (12.25)

-—^ = A’"pml49Pm(r, f) — i;“’Sm’49Sm(r, t)4>{r, t) + (r«Sm148Sm(r. t)<f>(r, t).

01

The initial condition will have to be calculated taking the complete fission product chain into account. For a fresh core, however, the concentrations of l48MPm, l48Pm and l48Sm are zero and this simplifies the equation to the form which is most commonly encountered in textbooks. This simplified form can also be used after burn-up, but it can easily lead to errors of 10-20% in the l49Sm concentration. There is, in any case, no difficulty in treating the complete system (12.25) on a digital computer.

For the treatment of Xe and Sm transients (long time dynamics), beside eqns. (12.24) and (12.25) a further equation is needed, relating the poison concentration to the neutron flux. This could, for example, be eqn. (12.1) but, as already mentioned, it is not necessary to treat in this case the delayed neutrons, and the neutron flux can be

considered as quasi-stationary. The multi-group diffusion equation [see eqn. (4.40)]

N N

— D, V2$,(r, f) + 2,чфі(г, t) = 2 2»o, k~i<Mr, f) + Л’і 2 Скї/кфкіг, t) (12.26)

к = 1 к = 1

is normally used. (In many computer codes the treatment is limited to one energy group.) The poison concentration appears in the coefficients 2,,- of eqn. (12.26).

In general all the coefficients of (12.26) are temperature-dependent and the heat- transfer equations will have to be solved together with eqns. (12.24), (12.25) and (12.26). This system of equations is usually programmed in the computer codes used to study the reactor operation (long time dynamics) and the spatial stability. The temperatures can in most of these cases be treated as quasi-stationary.

Although these equations are space-dependent, preliminary assessment calculations are often performed in the zero-dimensional approximation.

In this case all variables are supposed to be separable in space and time. Taking as an example eqn. (12.24) we have

ф(г, t) = <MO/o(r), Xe(r, t) = x(t)/„(r), I(r, t) = i(t)f0(r).

Substituting in (12.24)

A complication is given by the last term of the Xe equation where, because of the product of flux and Xe concentration, the term /o(r) cannot be eliminated.

As an approximation /o(r) is substituted by its mean value

The definition (12.27) can be such that /0 = 1.0, in which case the term fair) will completely disappear from eqn. (12.28). (This way of simplifying eqn. (12.28) is mathematically incorrect, but it can be proved that the error made here is small.)

In this space-independent case eqn. (12.26) is usually limited to one energy group and the leakage term DV2<f> is substituted by a buckling term — В2ф. The Xe absorption appears in the term 2,<f>. This term after substitution with (12.27) would include/02(r) so that the space-dependence would not disappear from this equation. This difficulty is avoided using in this equation an average Xe concentration weighted over the product of the flux ф and the importance ф‘ (see §4.13).

Assuming ф’ = ф we have

From eqn. (12.28) with /0(г) = /0 = 1 one can see that the equilibrium 135Xe concentra­tion is given by

л — = + (12.29)

This shows the dependence of the 135Xe concentration on the flux level, and its saturation for high flux values [this saturation value is given by the limit for ф -><*> of (12.29)].

Analytical solutions of eqn. (12.28) can be found for the cases in which ф remains constant with time after an initial step variation.

It is possible to see that if the power level, and hence the flux of the reactor, is suddenly reduced the loss of l35Xe due to absorption decreases. Consequently the 135Xe concentration builds up to a maximum from the decaying l35I which has been previously formed. Ultimately the radioactive decay of l35Xe takes over and the total l35Xe concentration drops off, after having gone through a maximum. The reactor must have enough reactivity invested in control rods to be able to “override” this Xe peak, otherwise the reactor will be shut-down by the poisoning effect of Xe. For analogous reasons a minimum in the Xe concentration occurs after a power increase. (See also § 12.16 and Fig. 12.8.)

Similar phenomena appear for 149Sm but they are not so important for the reactor operation because of the smaller absorption cross-section of 149Sm and the longer decay constant of l49Pm. Since 149Sm is a stable nuclide it does not go through a maximum after shut-down but tends toward an equilibrium concentration which does not depend on the flux level.

Typical equilibrium Xe concentrations range for high temperature reactors between 1 x Ю14 and 4x 1013 atoms/cm2 and its reactivity worth is of the order of 2-3.5%.

For 149Sm at average burn-up the equilibrium concentration ranges from 2 x 10’5 to Зх 1014 atoms/cm2 and its reactivity worth is of the order of 0.5-0.7%.

The zero-dimensional methods we have now analysed are very often used for preliminary calculations, but an accurate value of the excess reactivity needed for Xe override can in most cases only be obtained with a space-dependent calculation.

In the design of nuclear reactors neutron physics calculations are needed in order to determine and optimize a certain number of quantities like:

Reactivity.

Spatial distribution of power and temperature.

Maximum local values of power density, temperature, fast neutron dose, burn-up.

Quantities related to safety and control (temperature coefficient, control-rod require­ment and effectiveness, load-following capability).

Costs (capital investment, fuel-cycle costs, etc.).

There are various ways of obtaining these data. The calculation methods evolved with time because of the increase in the amount and detail of the available experimental data and of the computational effort which it was possible to invest in the problem. This last point has been very strongly influenced by the introduction of increasingly faster computers. The computer cost for a given calculation has decreased by almost a factor of 10 in the last decade. The modern tendency, now almost universally followed, is to base the calculations on very general basic cross-section data for all nuclides involved. The most recent libraries can describe with rather good accuracy all interesting nuclear reactions in the energy range significant in reactor physics. Starting from these data it is possible to calculate any reactor system with an accuracy which is in great part dependent on the computational effort involved and the complication of the core and fuel geometry.

In the past the tendency was more oriented towards direct measurements of the quantities which were of interest. This involved rather expensive measurements centred on particular core arrangements. Simplified calculation methods were used and the results were adapted to the measured values through correlation of some rather arbitrarily chosen parameters. In this way it was possible to interpolate between the experimental data. This procedure was expensive and rather unreliable for extrapola­tion to situations after a high burn-up where it is difficult to obtain good experimental data. Furthermore, a consistent optimization of the various parameters which can be varied by the reactor designer is only possible with the more accurate modern approach. The older correlation methods (which have also been called the macroscopic approach) evolved slowly towards the present ones (microscopic approach) with the improvement of the computer methods and of the data. This was done going to a more detailed description of the energy dependence of flux and cross-sections (increase in the number of groups) and increasing the detail of the space-dependent calculations (use of transport theory or Monte Carlo instead of diffusion, use of two — or three­dimensional calculations). This evolution towards the use of more basic data and more

detailed calculations is not yet completed and we still observe a tendency towards methods which are numerically more complicated, but physically simpler because they more exactly simulate the phenomena occurring in the reactor.

The work in the physics of high-temperature reactors has been started rather late compared to the work for the other thermal and fast reactors. Besides, as we will analyse in detail further in this book, there are many fundamental differences between the physics of HTRs and of the conventional thermal reactors. This meant that a straightforward application of the “recipes” used at that time for thermal reactor calculations was not possible. For this reason the HTR physicists were among the first to rely on detailed computer calculations, and the Dragon and AVR reactors have been among the first reactors to be built without a previous series of reactor physics experiments on similar lattices. This was, of course, facilitated by the use of highly enriched uranium, so that the reactivity calculations did not need the accuracy necessary to reach criticality with natural uranium.

The existence of good theoretical methods does not eliminate the need of checking their accuracy and reliability. Because of this reason series of measurements have been performed also for high-temperature reactors. This point will be discussed later; we need only note that the measurements have been used to check the accuracy of the calculation methods and of the data and to point out the areas where improvement was needed, but not for correlation purposes.

The treatment of energy and space-dependent neutron behaviour in the complicated geometry of the reactor core with space dependent burn-up, composition and tempera­ture and with variably inserted control rods is a tremendous task even for the best computers now available. This means that it is necessary to introduce simplifications, e. g. restrictions in the treatment of energy dependence (low number of groups) and of space dependence (low number of regions and mesh points), but this is not possible beyond a certain limit without an unacceptable loss in accuracy.

A separation of the space and energy dependence is not possible on the whole system, but it may be still sufficiently accurate in smaller domains. This procedure results in a series of zero dimensional spectrum calculations performed for different parts of the reactor, then combined with a space dependent calculation in which the treatment of energy dependence is greatly simplified (few group calculations). It is difficult to justify this method from a strictly mathematical point of view, but normally it works well and it is widely accepted as a standard method for reactor calculations. According to the complexity of the problem the designer has to choose the treatment of energy dependence (number of energy groups) for the space dependent calculation. This method implies the presence of large regions with rather small spectral changes, but this is usually the case for the normal size of power reactors.

The problem is complicated by the fact that beside a macroscopic neutron spectrum variation between the various reactor regions, we have a microscopic variation of the neutron spectrum across the fuel elements and the moderator in any reactor which is not completely homogeneous. Although high-temperature reactors are much more homogeneous than most reactor types, this problem cannot be neglected. Here again the mathematical treatment can range from energy-dependent multi-dimensional cell calculations to zero-dimensional treatments where the fine structure of the neutron flux is taken into account multiplying all cross-sections by energy-dependent self-shielding factors defined for each material in such a way as to give, in the homogenized calculations, the same reaction rates as in the real cell. The accuracy and level of sophistication required by any reactor calculation depends on the design stage considered.

As with any engineering problem, the process of designing a nuclear reactor goes normally through various stages of increasing sophistication. Because of the very high number of parameters which have to be defined, one must first perform a set of survey calculations covering a rather broad range of possible cases. Economic evaluations performed in conjunction with these surveys will then define the optimum value of these parameters. This is usually an iterative procedure because some data will have to be given a trial value which can only be checked and corrected after having performed the calculation.

This problem will be dealt with in detail in the chapter about optimization. We shall start here with a brief analysis of the methods and data.

7.1. General considerations

As we have seen in Chapter 3, the neutron reactions in the energy range between a few eV up to about 100 keV take place in rather sharp resonances. Although resonances occur for all heavy nuclei and involve scattering, fission and absorption cross-sections, only the resonance absorptions of 238U and 232Th are very important in HTRs, as in most thermal reactors.

Resonance absorption in some fission product or structural material can be treated in many computer codes, but it has no great influence on the reactor behaviour. Fission resonances are of no practical importance for thermal reactors (with possibly one exception: the low-lying Pu resonances).

Resonance absorption in 232Th and 238U is on the contrary one of the most important phenomena in HTR reactor physics. The amount of absorption in these fertile materials determines the core reactivity and the conversion factor, hence influencing the whole burn-up behaviour of the reactor. Besides the Doppler effect in these resonances is of primary importance in determining the temperature coefficient of the core, and consequently the dynamic behaviour of the reactor and its safety characteristics. The transport theory methods we have seen up to now are theoretically perfectly suited to treat the resonance absorption in any reactor.

In practice the very strong energy dependence of cross-sections in this energy range implies the treatment of a very high number of energy points (up to 14,700 discrete neutron velocities are considered in the General Atomic GAR codeU)). Furthermore, the fuel geometry is rather complicated in most reactors, and in HTRs the use of coated particle fuel creates a double heterogeneity both on a microscopic (coated particle) and on a macroscopic (fuel compact) level. Because of the high cross-sections at the resonance peaks diffusion treatment of the space dependence is not usually possible and transport or Monte Carlo methods have to be adopted.

The integral form of the Boltzmann equation has proved to be the best transport tool for these problems, so that collision probability methods are normally used for refined calculations.

Following these calculation schemes it is possible to obtain accurate results taking into account overlapping and mutual shadowing of different resonances and complex geometries. Unfortunately this is a very lengthy procedure necessitating computation costs which are too high for routine or survey calculations. The number of assumptions, simplifications and cook-book type recipes which have been developed in order to perform quick resonance absorption calculations has made of this field one of the most inextricable domains of reactor theory.

Here is given a short survey of the approximations which are most commonly made.

Of course not all of these approximations need or can be made at the same time.

1. The resonance energy lies below the fission spectrum so that as neutron source only slowing down needs to be considered.

2. With the exception of the treatment of the lowest Pu resonances, thermal motion of the moderator (and hence up-scattering) can be neglected.

3. Overlapping and mutual shadowing of resonances is usually not very important in HTRs. It is then possible to assume that between resonances the flux has its asymptotic 1 IE behaviour, and each resonance can be treated separately.

4. The asymptotic ЦЕ behaviour can be assumed also at resonance energy for the outside moderator.

5. The resonance can be considered as being so narrow that any collision, even with heavy absorber atoms, slows down the neutron outside the energy range of the resonance. (Narrow resonance approximation.)

6. In broad resonances the slowing down due to scattering with heavy absorber nuclei can be neglected. (Infinite mass approximation.)

7. In many calculational procedures the geometry is simplified to the treatment of only two regions (fuel and moderator).

8. Collision probabilities can be calculated assuming flat flux in those two regions.

9. A simple rational formula can be used to calculate collision probabilities (Wigner rational approximation).

10. Considering that a neutron can avoid being absorbed in resonances either escaping from the fuel region or being slowed down by scattering in the fuel region, one can simulate a heterogeneous fuel-moderator geometry with a homogeneous system. This equivalent system has the properties of the fuel region to which a fictitious scattering cross-section has been added in order to substitute geometrical escape with slowing down outside the resonance energy. (Equivalence relation.)

11. For any moderator an equivalent amount of hydrogen can be found so that any mixture absorber-moderator can be substituted with an equivalent hydrogen-absorber mixture. In this way only the data for this equivalent mixture need being stored on a computer.

Let us consider in more detail various cases and approximations.