Ann. nucl. Energy, Vol. 17, No. 7, pp. 389-392, 1990
0306-4549/90 $3.00+0.00 Copyright © 1990 Pergamon Press plc
Printed in Great Britain. All rights reserved
TECHNICAL NOTE SIMPLE A N A L Y T I C A L FITS TO F I S S I O N - P R O D U C T D E C A Y HEAT M . S. SRIDHARAN Reactor Physics Division, Indira Gandhi Centre for Atomic Research, Kalpakkam 603 102, Tamil Nadu, India
(Received 29 January 1990) Abstract--Several computer codes exist for detailed summation calculation of fission product decay heat. Results o f this calculation are usually fitted to an analytical function. This fit becomes handy when extensive routine estimates are required to be made. Some analytical fits, involving many exponential terms, have been reported in literature. In this paper, it is shown that another type of fit, simple in nature, can be obtained and employed for the decay heat prediction. This is accomplished by dividing the time range into a number of segments. Three case studies involving 47/14/1 time segments are undertaken and analytical fits corresponding to these segments are obtained. The fitting accuracy relative to the summation calculation is examined. Overall deviation (for irradiation and cooling time ranges specified in this paper) is found to be ___1% for the 47 time segments case, + 6 % for the 14 time segments case and + 2 0 % for a single time segment case. These results appear to be satisfactory. The work carried out here, is confined to fast fissions in 239pu.
l. INTRODUCTION In a nuclear reactor, the heat generated by the decay of fission products (FP) forms a major source of power after reactor shutdown. To calculate this heat, several summation codes are used. Usually, for the purpose of routine calculations, the summation codes are used to the extent o f obtaining a fission pulse functionf(t), which is defined as the mean decay energy release rate per unit fission at a time t after an instantaneous burst of fissions. This function, different for different fissionable species, is fitted with a set of exponential functions for each of these species (Tasaka et al., 1983; Goerge et al., 1982). To compute FP decay heat following any irradiation history, f(t) is then integrated over the irradiation history (T, z) and summed for constituent fissionable species of the fuel. Here, T denotes irradiation time and z refers to cooling time. The integrated function nm(T, z) for a fissile nuclide type m, will be a sum of certain numbers of terms involving exponential functions. Thus the fission pulse fro(t) is given by :
or :
H(T,
T) = ~F~H~,(T, r) m
(m over all fissionable nuclides), where Fm is the fission rate in the ruth species. Since the present work is restricted to 239pu, we shall drop the subscript m hereafter. This paper aims at an alternative approach for obtaining H(T,z). Towards this goal, a compact analytical expression is deduced by properly studying the decay characteristics. But then, this expression will have variable coefficients, i.e. the coefficients will be different for different time ranges of T. In the following, we will examine the decay characteristics, deduce a compact expression on their basis and then apply this expression to the results of detailed summation calculation to get the fitting parameters. A comparison of resulting fits vis-d-vis summation calculation will then be made to determine the accuracy of the fits.
fro(t) = ~ ctm,e -am~, I - l,n
where O~mland 2m~are constant parameters, and the integrated decay heat following an irradiation of Tseconds at a constant fission rate of one fission per second and cooling time z is given by :
2. DECAY CHARACTERISTICS AND GOVERNING EQUATIONS
Let us associate an energy factor e~with each fission product nuclide i. ei, then, represents the average decay-energy content of the nuclide such that, on average, this much energy would have been released when the nuclide is ultimately H~,(T, z) = ~ (ct,,t/2,~t) e-~,~'(l --e-A-'r). transformed to a stable one. ei may be called decay heat /=1 capacity or decay heat potential of i. To illustrate this, conThe value of n is equal to 31 in the work of Tasaka et al. sider a simple decay scheme (Fig. la) where a nucleus A ~, (1983). For a mixture of fissionable nuclides, the integral decays to A 4 in steps. E~,E2 and ~3 are the mean decay energies decay heat becomes and e ~, e2 and e3 are the decay heat capacities of the nuclides A l, A 2 and A 3, respectively. H(T, z) = ~Hm(T, z) We then have the following relations : m 389
390
Technical Note
(a)
3. TO DEDUCE A COMPACT ANALYTICAL EXPRESSION Due to energy conservation we can state, in the case of an instantaneous fission event, that : (2b)
dE/dt = -p,
~2' ~ 2
El ' E1
~3' ~3
StoNe
where : E = ~ 8~Ni is the total decay heat capacity and
(2c)
i
p = ~ 2iN~E~is the decay heat release rate.
(2d)
i
(b)
5~
~ , ' J ~ 3
E1
Though equation (2b) is obviously a manifestation of conservation of decay energy, it can be proved by a detailed treatment (see Appendix). It is considered worthwhile to obtain a plot of p t as a function o f E. For this purpose both p and E are obtained as functions of time. To evaluate E ( t ) as defined by equations (1) and (2c) is a very tedious task. Hence a different route is adopted, viz. to find E(t), equation (2b) can be integrated as follows :
g12 + 013 --I
E3 E(t) =
Fig. 1. (a) A simple fission-product decay chain (schematic). (b) A part of fission-product decay chain involving branching (schematic).
83 = E3,
t 2 - 8 3 = E2 (or)82 = E2-F83
=
E2"FE3,
dE = -
p dt.
This can be done by appropriate modification in the summation calculation. For the limit o% one may choose a higher limit, as for example, 1013 s. A plot of p t vs E corresponding to a single fast fission event (in 239pu, for example) can now be obtained as in Fig. 2. In this figure, the curve can be divided properly into different segments so that each segment can be considered to be linear in E, i.e. (3)
81 --82 = E I ( o r ) el = E1@-82 = E1-I-£2-F£3.
tp(t) = # , E + k r
Let us now consider a part of decay chain (Fig. lb), involving branching o f a nuclide A 1 (with mean decay energy ~1 and mean decay heat capacity 80 to A 2 with a transition probability YI2 and to A 3 with a probability g13. We have then components of Et for the two transitions, (81--82)gl2 and (81--83)gl3 :
for an rth segment, where #, and k, are constants for this segment which lies between the time points, say t, and tr+ 1' NOW differentiating equation (3) with respect to t and
1.20
E~ = g~28~--g]282-I-g]381 --gl383 = (g12 + g 1 3 ) 8 1 =
81 - -
~
(g1282 + g 1 3 8 3 )
glk Sk •
k=2,3
0.90
In general, for the ith nuclide species : Ei = 8 i - ~giksk.
c .o
(1)
k
In a nuclear reactor, the changes in fission product inventory arising from the buildup and decay of these nuclides is governed by the following set of differential equations (excluding neutron capture effect) : (dN,] dt) = a,F + ~. 7j,Nj -- 2iN,,
(2a)
>
0.60
~E o. 0.30
J
where : At,. = number of atoms of ith nuclide species, 2~ = decay constant of ith nuclide species, ~,j~= partial decay constant o f the j t h nuclide to the ith nuclide (~j~= gji2j where g/j is the branching fraction for the transmutation process), a~ = fission yield for nuclide i, F = fission rate.
0
I
I
I
I
3
6
9
1;=
E (MeV)
Fig. 2. ~ission-product decay heat as a function of decay heat capacity following an instantaneous pulse of fast neutron fissions of 239pu.
Technical Note using equation (2b), we get : t(dp/dt) = - # , p ,
(4)
where fir = ( l "+'#,).
391
equation (5) since B, = - / / , + 1 [equation (6)]. Once B, is known, other coefficients A, and 6", can be easily obtained using equation (7a). To find H(T, z), the decay heat after a given irradiation time of T seconds (at constant fission rate) and cooling time z seconds, we write :
Equation (4) has the solution :
H ( T, z) = P( T + z)-- P(z).
p(t) = ~(,t-#,,
(5)
where
H(T, T) = H(oo, z ) - - H ( ~ , T+z).
~r = p(tr)t~r r"
Equation (5) holds good for an rth segment bounded by a pair of time points t, and t,+ ~. In what follows, capital T will be used whenever we want to refer to irradiation time. In the case of irradiation over T seconds, at fission rate F(t), the decay heat immediately after T seconds, P(T), is obtained as : P( T) = =
f(t)p(t) dt
?
F(t)p(t) d t +
o
= P(E)+
:/
fi,
F(t)p(t) dt
{m,(t- Z) +Fr,}p(t) dt
[since, for a linearly varying fission rate : F(t) = m , ( t - T,) +F(T,) when 7', ~< t ~< T,+ i, i.e. t lying in rth time segment. We have used F(T,) = FrJ
.'.P(T) = P ( r 3 + m , ;r, tp(t) d t + ( F r - - m , T ~ ) I, r p(t) dt •
= P(T~) + m " 1L -"ArB' --r - o . T
JT~
T. + I] r + (F~ - - m , T J [ A , T ~ q ~ , r,
,
,
Also, H(T, z) = H(oo, z)--H(oo, 7") when T >> T]. 4. RESULTS The present work is confined to fast fission in 239pu only. This work can, of course, be extended to other fissionable materials and a given fuel composition. We have made 3 different case studies involving 47, 14 and 1 time segments and obtained corresponding fitting coefficients A,, B, and (7, of equation (7a). The time points were selected judiciously and the double precision mode was adopted in the computation. In the case study pertaining to 47 time segments, we find that it is enough, if the mantissa of the coefficients A, and B, are up to 8 decimal places (i.e. singie-precision) but it is necessary that the mantissa of C, be in double-precision mode. This is true for all irradiations between 1 and 10 ~° s, the latter being the upper limit in the present study. The resulting fits are able to reproduce the results of detailed calculation (Murthy and Singh, 1976) to which we have fitted within + 1%. Fitting accuracy relative to calculated results is plotted in Fig. 3 for a typical case of 1 yr irradiation. Relative discrepancies of exponential fitting of Tasaka et al. (1983) is also plotted in this figure for comparison. We find that our fits compare better than this. With an increasing irradiation period and decreasing cooling time, the general trend observed is that the minimum precision needed for
where A, and B, are constants for an rth time segment given by: B,=-~,+I A, = ct,lB,
(8)
[An alternative to equation (8) is, however, possible using equation (8) itself, i.e.
0.20
(= -~,)}
D
0.15
.
(6)
X
0.10
For a linearly varying irradiation step with slope m,, knowledge of the parameters P(T,), A, and B, will facilitate calculation of the decay heat at the end of the irradiation step. If more such steps are involved, this equation can be used as a recursive relation to obtain the required decay heat. When m, = 0, i.e. if the rth irradiation step is constant at fission rate/7,, : P( T) = P( T,) + Fr,[A, Tn,]r,
-
0.058
D
0.00 ~L~ [3
o
~ -0.05 ~
-0.10
× Our
(7a)
-0.15 f
B,.
(7b)
If the fission rate is constant at F for all segments included in T, then Fr becomes F in equations (7a) and (7b). We thus (~btain a set of constant parameters for different time segments. The cofficient B, can be found using
work
o JNDC
work
-0.20
where C r = (e(Tr)/FT)--ArTr
[3
h
or
P(T) = Frr[A, T B"+ C,],
0
--
-0.25
D
I t00
I
J
I
I
101 102 103 104 105
I
I
I
106
10 7
10 8
CooLing time (s I Fig. 3. Relative discrepancies (%) between fits and calculated results of the total sensible decay power after a year's irradiation of 239 Pu by fast neutrons.
392
Technical Note
Table 1. Fitting coefficientsa for the case study involving 14 time segments Time range in seconds (lower bound) O.l+lb 0.6+ 1 0.6+2 0.4+3 0.18+4 0.30+4 0.20+5 0.20+6 0.20+7 0.50+7 0.20+8 0.10+9 0.20+9 0.50+9 to 0.10+ 10
A
B
C
0.14627+01 0.39514+02 -0.13811+02 0.12844+02 -0.20204+02 -0.55124+02 -0.20967+02 -0.86445+01 -0.54091+02 -0.23165+02 -0.35612+03 -0.10732+08 0.56448-03 0.56462-05
0.383lOfOO 0.25493 -01 -0.14648+00 0.49025-01 -0.19841+00 -0.41736+00 -0.27584+00 -0.16155+00 -0.33736+00 -0.26865+00 - 0.46063 + 00 -0.10693+01 0.24231+00 0.4472 I + 00
-0.80295+00 -0.39257+02 0.12185+02 -0.10786+02 0. I2328 + 02 0.10152+02 0.10633+02 0.11113+02 0.10689+02 0.10759+02 0.10660+02 0.10616+02 0.10544+02 0.10573+02
5. CONCLUSIONS
In this paper, we have been able to show that simpler analytical fits could be obtained in the place of a single lengthy fit involving many terms. But then, this simplicity is not without penalty. It calls for many such fits for different time segments into which the entire time range is divided. The curve shown in Fig. 2 is quite interesting and it might be worthwhile to probe into it to understand the physics of it from a microscopic viewpoint. Regarding the analytical fits obtained above, they appear to be adequate for routine calculations. Depending upon the accuracy needed, one can choose appropriate fits viz. that correspond to 47/14/l time segments. Since a computer is normally depended upon for decay heat studies pertaining to the design of decay heat removal system, particularly when parametric studies are made, the case study with a large number of time segments (47-segment case study, in our example) can be used with the aid of a computer. author thanks Dr K. P. N. Murthy for his useful comments/criticisms and Shri R. Shankar Singh for his encouragement.
AcknowledgementsThe
a For large irradiation. bReadasO.lE+IorO.lxlO+l.
mantissa of A,, B, and C, decreases. In the case of large irradiation of the order of hundreds of days, it is enough, the parameters A,, B, and C, have mantissa up to five decimal places and the accuracy of the fits is k 0.1% . A similar study has been carried out for the case involving 14 time segments and its characteristics are similar to the above case study except that the fitting discrepancies become larger, viz. & 6%. In the case of large irradiation of the order of hundreds ofdays, the fitting accuracy is f 1% for a cooling time up to lo7 s and +6% for higher cooling times. As a typical example, for large irradiation, we have given the fitting coefficients in Table 1. The decay heat obtained using these fits are in MeV s- ‘. Yet another fit has been obtained by considering the time points (I, and t2) of 1 day and 1000 days (in seconds), respectively. The resulting fit is :
H(T,z) = P(T’ = T+r)-P(T” = F.A[(T+z)B-tS]
= T)
REFERENCES Goerge D. C. et al. (1982) Report LA-9362-M& LANL. Murthy K. P. N. and Shankar Singh R. (1976) Internal Report RRC-15, IGC. Tasaka K. et al. (1983) JNDC Nuclear Data Library of Fission Products, Report JAERI-1287.
APPENDIX In the case of a single fission burst, equation (2a) of the main text reduces to : dNildt = ~Y,,N,-~~N,. Multiplying both sides by &iand taking sum over all i, we get :
MeVs-‘,
c,dt%d
subject to and t, < T+z Q t2. f,
=1ii ~Y,,N,~
ij
= 1 ~y,JQ,-
~kN,E,>
MU++;m
[using equation (l)]. or
since, .&N, = E,
with P(Z) = FA,?] +C, and P(T+r)
= FA,(T+++&,
where A, = -71.25, B, = -0.3589, C, = 10.73 and C2 = 10.69. This is subject to the condition tl < r < t2. This fit is found to have an uncertainty margin of f 10% relative to the summation calculation.
total decay heat capacity and liSi!i = Yik> partial decay constant of the ith nuclide to the kth one. Therefore, (dE/dr) = --p where p = cLiNici is the decay heat release rate.