- Email: [email protected]
On the group method for the approximate solution of a hierarchy of rate equations describing nucleation and growth kinetics
Journal of Nuclear Materials 59 (1976) 175-182 @ North-Holl~d Publishing Company
ON THE GROUP METHOD FOR THE APPROXIMATE SOLUTION OF A HIERARCHY RATE EQUATIONS DESCRIBING NUCLEATION AND GROWTH KINETICS
OF
M.R. HAYNS L%eoreticafPhysics Division, AERE Hatwell, Oxfordshire OX1 I ORA, UK Received 6 October 1975
An approximate scheme devised by Kiritani for the solution of a large number of rate equations describing the nucleation and growth kinetics of point defect aggregates in irradiated materials is discussed in the light of serious objections raised by Koiwa as to the validity of such a scheme. By means of calculations on the particular problem of the nucleation and growth of interstitial dislocation loops it is shown that the grouping method is a useful approximation and suggestions are made as to the best choice of parameters for its application. Since it is an approximate scheme and is being applied to very complex systems, which can only yield to numerical methods, it is pointed out that great care is required to ensure that the method is not used inappropriately. Un schdma approximatif propose’par Kiritani pour re’soudre un grand nombre d’e’quations de vitesse decrivant les cine’tiques de germination et de croissance d’agregats de de’fauts ponctuels dans les ma&iaux irradie’s est discutia la Iumikre d’objections s&ieuses soulevees par Koiwa au sujet de la validite’d’un tel schima. A l’aide de calculs sur le probleme particulier de la ge~ination et de la croissance de boucles de dislocations interstitielles, il est montrg que ia m6 thode des groupes est une approximation utile et des suggestions sont faites B propos du meilleur choix des paramktres en vue de l’application de cette mithode. Puisqu’il s’agit d’un sche’ma approximatif et qu’il doit 5tre applique ‘a des systemes t&s complexes, ce qui peut conduire seulementi des me’thodes num&iques, on insiste sur ie fait qu’il faut utiliser cette me’thode avec grande precaution pour que celle-ci soit utilisCe ‘a propos. Ein NLherungsschema, das von Kiritani zur Losung einer grossen Zahl von Geschwindigkeitsgleichungen zur Beschreibung der Keimbildungs- und Wachstumskinetik von Punktdefektaggregaten in bestrahltem Material angegeben wurde, wird im Lichte der von Koiwa erhobencn gravierenden Einwande beziiglich der Giiltigkeit eines solchen Schemas diskutiert. Mit Hilfe von Rechnungen zum besonderen Problem der Keimbildung und des Wachstums von ZwischengitterVersetzungsringen wird gezeigt, dass die Gruppenmethode eine niitzliche Naherung ist. Zu deren Anwendung werden Annahmen hber die beste Auswahl der Parameter gemacht. Da das Schema Gherungsweise gilt und fir sehr komplexe Systeme Verwendung findet, die nur numerisch behandelt werden kBnnen, wird darauf hingeweisen, dass grosse Sorgfab erforderlich ist, urn sicher zu gehen, dass die Methode nur in geeigneter Weise angewandt wird.
it is necessary to consider clusters containing perhaps as many as lo6 defects. Obviously the full set of equations cannot realistically be solved and Kiritani constructed an approximate scheme whereby whole groups of cluster sizes were spanned by a single equation, with the rate constants chosen as being the average for the clusters covered by the group. In this way a scheme could be devised which would follow the growth of quite large clusters, at least in an approximate way. More recently Koiwa f2] has analysed the grouping approximation by means of comparisons between the analytical and numerical results for a
1. Introduction In a recent publication Kiritani [I] has discussed an approximate scheme for the solution of hierarchy of rate equations describing the nucleation and growth of vacancy clusters in irradiated metals. Such a scheme involves equations which describe the rate of accumulation and loss of single point defects by defect aggregates, consequently one rate equation is required to describe the history of a particular cluster size. In order to extend the calculation to clusters large enough to correlate with experimental results, 175
176
M.R. Hayns j A hierarchy
very simple model system which is ameanable to analytic solution. He raises serious doubts as to the validity of the grouping method. If one ever hopes to gain a quantitative understanding of the kinetics of nucleation and growth of defect aggregates in irradiated materials then methods similar to those described by Kiritani would seem to offer the only general course to follow. It is, of course, true that exact analytical solutions can be obtained in a few very special cases [3], but in general the complex nature of the systems involved in irradiation damage in structural materials precludes such treatments. Consequently it is important to establish the applicability of the grouping method in relation to problems of practical interest. In the past we have used the hierarchy method to investigate the nucleation and early growth of rare gas bubbles in operating nuclear fuels [4] and to study the problem of interstitial loops in M316 stainless steel [5], with particular emphasis on the void swelling problem. The above-mentioned work used a hierarchy of single equations, but recently we have included in the scheme a grouping method similar to that described by Kiritani. Since we are unavoidably concerned with the very complex systems found in realistic applications, we must turn to numerical methods and in particular the hierarchy plus grouping has proved to be extremely useful. Therefore it is of prime importance to establish the conditions under which such methods are realistic and where difficulties may be found. It is the purpose of this present work to establish the validity of these methods and to respond to the objections raised by Koiwa.
2. Rate theory model of interstitial and growth
loop nucleation
As Koiwa rightly points out, a detailed, mathematically sound, discussion of the very complicated equations describing even simple nucleation and growth problems is extremely difficult. Consequently we shall confine our comments to the numerical SO]Ution of a particular set of rate equations which are typical of those found in irradiation damage processes. In particular we shall consider the nucleation and growth kinetics of interstitial dislocation loops in M3 16 steel. This system was chosen simply because
of rate equations it forms the basis of a much more detailed study which will be presented elswhere [6]. Previously [5] we considered a similar problem, but used only a single group size hierarchy and we ignored the effects of the thermal emission of vacancies. In the present work both a group method and the thermal emission of vacancies have been included. The system of equations may be written
+ L,CzCv - DiPDOCi
2Ani + 1
Most of the quantities appearing in these equations have been defined previously [5] t. However, the definitions are given again in table 1 for convenience. In order to keep the system of equations relatively simple we have ignored the thermal emission of vacancies from all sinks other than the interstitial dislocation loops. The equations may be briefly interpreted as follows; the first two equations represent the balance of vacancies and interstitials between production, 7 We should note that the sink strength associated with the interstitial loop used in the calculations is a & where n is the number of interstitials in the loop. This sink strength was not found simply by counting the number of edge sites available, thereby giving a purely rate controlled growth, but by solving the diffusion to an effective sphere, the rate controlling process then appears as a geometrical factor which scales the fluxes into the loop according to the ratio of the areas of the effective sphere and the loop edge. The derivation of this sink strength has been given in detail previously [ 71 whilst a more complete analysis will be presented shortly [8].
MR. Hayns /A hierarchy of
rate
equations
Table 1 Definition of the quantities appearing in equations 1 I_____
QUWity
Definition
C”
Vacancy concentration t Interstitial concentration Concentration of interstitial loops containing n atoms Production rate of vacancies and interstitials (assumed equal) Recombination rate Rate constant for the production of loops containing n + 1 atoms by the addition of single atoms to loops of n atoms.
ci cn x K1 RN
Rn = 0.78fiVi Ln
Vi = attempt frequency; Em’ = interstitial migration energy; k = Boltzman constant; 7’ = temperature. Rate constant for the production of loops containing n - 1 atoms by vacancy absorption. L, = 0.78&v,
En
exp(-Emi/k73
exp(-~mv/k~
E ’ = vacancy migration energy. R?te constant for the production of loops containing II + 1 atoms by vacancy emission. E, = 0.78&v,
exp[-EmfjkT]
exp { -478~ + Felt a’/kT}
Emf = vacancy formation energy; 7s~ = stacking fault energy; Fe1 = the elastic energy of the loop [7] and stress effects have been ignored.
Di
D P> $i
Interstitial diffusion coefficient Vacancy diffusion coefficient. Concentration of network dislocations Width of the jth group, given by (rz._ 1 + 1)/M Number of equations defining sing ie loop sizes
.-
t all concen~ations
---
are fractional.
at rate K, and losses by recombination, vacancy absorption by loops, the loss of interstitials to existing loops and by the nucleation of new di-interstiti~ primary loops and the respective losses to an existing network dislocation structure of density pDo Cm/Cm3. Also included in the vacancy equation is the contribution by thermal emission from the existing interstitial loops, with the rate Ei. The remaining equation is representative of the hierarchy descfibing the individual cluster sizes. The terms on the right hand side are to be interpreted as the rate of production of loops of size nj by accumulation of interstitials, at rate Rj_ 1 and thermal emission of vacancies at rate Ei_1 from the next smallest loop size, the losses from the current size by further accumulation of interstitials and thermal emission, the loss by vacancy absorption and finally the backward reaction caused by vacancy absorption in the next highest loop size. The weighting factors in this equation correspond to the grouping method adopted and Ani is .the width of theith group, given by nj - nj_ 1. The coefficients
correspond exactly to these used by both Kiritani [l] and Koiwa [2]. A detailed discussion of the physical basis of these equations is not appropriate here, it has been given previously [.5] and, in relation to the thermal emission of vacancies will be discussed elsewhere [6]. The equations are based upon the assumption of homegeneous-binary nucleation, with such features as thermal dissociation and interstitial trapping not represented in the current equations. Equations (1) are, of course, solved numerically, in our case by means of a sophisticated subroutine based on Gear’s method, which is available as part of the HarweIl subroutine library and has been described in detail 191.
3. The results of the grouping method The grouping method, by collecting together loops of a size range into a single rate equation, the rate
M.R. Hayris /A hierarchy of rate equations
178
constants in which correspond to some average sink strength for the group, can be set up in a variety of ways. Kiritani adopted the scheme whereby the first 100 equations described only single loop sizes, whilst above this number the cluster sizes were obtained from the algorithm, Ani < O.O5(nj_l + 1). That is the size of cluster increases by 5% of the previous cluster size and as more and more equations are considered the width of the cluster can become very large indeed (i.e. Anj = lo4 at nj = 106). Koiwa, on the other hand chose several values of the number of un-grouped equations (30 and 70) whilst his Anj values were fixed at values ranging from 3 to 9. It is one of the purposes of this note to outline a method of specifying the appropriate choice of Anj for a given number of un-grouped equations which represents the optimum value within the general approximate scheme of the grouping method. We have used a’prescription for obtaining the width Anj as follows AHj
=(nj_l
+
1)/M
with M = 20 this expression reduces to the Kiritani equation. We shall show that simply by choosing M = number of un-grouped equations many of the objections raised by Koiwa can be eliminated. In fig. 1 we show the loop size distribution calculated using eqs. (1) but with no grouping. For our current purposes the solution of up to 200 equations is adequate and unless stated otherwise this is true for all of the calculations presented here. All of the results are for similar conditions, that is 723 K, network dislocation density of log cm/cm3 and a damage rate of 10S6 dpa/sec. Therefore these results are typical for solution treated M3 16 steel during fast neutron irradiation; the remaining materials parameters required are given in table 2. Fig. 1 is simply the loop size distribution at the longest time we could obtain using only 200 single equations of 0.54 X 10P4 d.p.a. (54 seconds at a dose rate of lop6 dpa/sec). In fact this is a relatively long time since nucleation * is completed very quickly, in about 10h5 d.p.a. in these calculations. The effects of introducing the group method are summarized in figs. 2 and 3. In fig. 2 we show loop * As previously [5] we are defining the end of nucleation as when the quantity (KzCf - LzC,C~)/C,(Z is less than 10e4.
R,C,
+ KICv)
L
8
12
16 20
2L
28
NO
OF
ATOMS/LOOP
32
36 Lo
LL
LE
52
56
60
Fig. 1. Distribution of loop sizes for parameters defined in tables 1 and 2 with no grouping, after a dose of 0.54 X 10e4 d.p.a. (54 seconds).
size distribution functions at rather early times (1 X lop5 dpa) for values ofM= 1, lo,20 and 50, with the hierarchy constrained to contain 10 single equations. The discontinuities in the distribution function at the cross over oetween single and grouped equations reported by Koiwa are apparent for M = 1, 20, and 50, whilst for N = 10 no discontinuity can be shown on the scale of this graph. We have connected the points with continuous lines, even though the functions are not defined except at the discrete positions, simply to highlight the differences between the curves. In fig. 3 similar results are given, but for a later time (3 X 1O-5 d.p.a.) when the peak of the distribution has passed through the single-group boundary. Here the divergence from the smooth curve obtained from the single equation
179
M.R. Hayns /A hierarchy of rate equations
Table 2 Parameters used for the calculations Quantity
Value
Interstitial diffusion coefficient, Di
10m3 exp(-0.2
Vacancy diffusion coefficient D,
0.6 exp(-1.4
Lattice parameter, a
3.63 X 10e8 cm
Network dislocation density pD”
lo9 cm/cm3
Temperature, T
723 K
Defect production rate, K
lo+ a/a.s/sec
eV/kr) cm’/sec eV/H’) cm2/sec
Vacancy formation energy Emf
1.6 eV
Stacking fault energy
9.4 X 1012 eV/cm’
Effective modulus b1 = p/(1 - u)
4.0 X 10” dynes/cm2
u-
10
9-o-
00-
-
M=l
0
M 10
I
M-20
&
M:50
=
r
0 M:lO
11
lo-
Full line smgle
A M:20 ..M=l
equations only
% x
Full lone single equations only
9-
70-
; ;
z
60-
b
8-
t Y 8
50-
< g F
LO-
7t i Y 8 w
Y E
30-
6-
5-
.?? 0
i J20-
=
L-
d II IO2
I
1
I I I 23L5670
I
II
I
NO
OF
’ 9
1”’
10
11
12
13
IL 15
16
ATOMS/LOOP
Fig. 2. Distributions of loop sizes after lo-’ d.p.a. for various values of M, the group width parameter. Appropriate values of An, the group width are given in table 3. The points for the grouped equations are shown at the mid point of the group. The value for M = 1 at n = 10 is out of scale (see text). All calculations contained hierarchies with 10 single equations.
2 !V l-
I
I
I
I
I
2
L
6
8
10
I
I
I
,
12 1L
16
NO
ATOMS
OF
I,,
1.3 20
I
I
22
2L 26
28
30 32 34
LOOP
Fig. 3. Distributions of loop sizes after 3 X lo-’ d.p.a. for similar conditions to fig. 2.
M. R. Ha.yns /A hierarchy of rate equations
180
solution when M = 20 is clearly shown, whilst the appropriate concentration for M = 1 is well out of the scale of this diagram. However, the curve for M= 10 does not show any discontinuity. In general the pattern noted by Koiwa of a shift in the peak of the distribution function to lower loop sizes whilst also showing a broadening is clearly seen in these figures. At 3 X 10m5 dpa the shift is in fact by one group, that is the peak occurs at 13.6 atoms/loop for M = 10 whilst for the single equations it lies at 16 atoms/loop. The reason for the discontinuities in the distribution function is that by taking too large a value of hj for the first group the number of loops being reduced in size per unit time by backward reactions into the last un-grouped equation is large, therefore artificially increasing the concentration of this loop size. The reverse situation is also true when the value of it is chosen to be very large, causing the width of the first grouped equation to be small. Therefore there is an optimum value of M, for which the transition is as smooth as possible. It is readily shown that for Ani_ 1 - Ani to be a minimum at the single-group boundary then M is simply the number of single equations. Therefore with IO single equations, M = 10 and ~~~=ll,l,~~~=l.l.Wehavecollectedtogether in table 3 values of group widths for IO and 100 single equations for various values of M Since the number of single equations is of necessity integral, we see that for 10 equations the smoothest transition gives the first group width as 10% wider than the last single equation. As the number of singles increases, this transition becomes smoother, at 100, there is only a 1% change. Kiritani’s scheme was an M value of 20, which gives a width of 5.05 to the first group. Our results are all for the rather coarse case of only IO single equations, however, even here the appropriate value of M = IO gives no discontinuities in the
distribution function which can be seen of the curves in figs. 2 and 3, in fact the numerical results indicate that such deviations are of the order of 0.1%. Consequently, it is suggested that the particular scheme adopted by Kiritani would not in fact lead to significant discontinuities in the distribution function. The results given by Koiwa using constant values of An of between 3 and 9 for 30 and 70 single equations would give rise to discontinuities since they are far from the optimum values of I .03 and 1.014 respectively, and these are indeed clearly shown in his results. It is important to emphasize that the accuracy of the grouping method is relative since any grouping will disturb the agreement with the exact result, which, following Koiwa we can assume to be given by the result of calculation using only single equations. However we are discussing an approximate scheme and as such some deviations from the ideal are inevitable, it is simply necessary to understand the inadequacies properly in order to make use of the scheme. In fact the details of the distribution function are not usually of prime importance in analysing irradiation damage experiments. Firstly because the experimental dif~culties ensure that any size distributions given are far coarser than the scale of deviations from the ideal as found by the grouping method (see for example the results of Okamoto and Harkness [ 121). Secondly the parameters of interest are not usually sensitive to the distribution function and thirdly such theoretical models as described here are inevitably based on approximations of the basic physical mechanisms which involves the use of sometimes very badly known quantities such as formation and activation energies, quantities which can introduce into the model a large degree of quantitative uncertainty. Furthermore, the models are based on current thinking as to the most important mechanisms involved,
Table 3 Values of A,, the width of the first group for 10 and 100 single equations for various values of M _. .___~_.. _._-._.- ______ ^ __~_ Number of single equations
1
5
10
20
50
..--.. 100
10
11
2.2
1.1
0.55
0.22
0.11
100
101
20.2
10.1
5.05
2.02
1.01
..___~_.._^.+M
1An
MR. Hayns /A hie~rch~ of rate equations
and others are not even included. In particular in the results presented here we have ignored interstitial trapping as an important mechanism in loop growth. Therefore the results can only be useful for correlations where it is absolutely certain that such processes are not important. When attempting to describe theoretically the complex interactions in a real material under irradiation we must be pragmatic and try to gain an understanding of the relationship between the various physical mechanics which eventually will lead to a fuller picture of the damage process and aid in the interpretation of the large number of experiments currently being performed. Therefore we need a suitable vehicle in which we can represent the physical mechanisms and investigate their ~ter-relation~ip. Such a vehicle is the rate theory coupled with the grouping method. The approximations inherent in it are understood reasonably well and until more sophisticated methods are produced it represents one of the very few ways of gaining understanding of the behaviour of these systems.
181
As an example of the type of data which can usefully be predicted by such a model we consider briefly the prediction of the total interstitial loop concentration. This is of special interest in the void swelling problem, where it has been shown that swelling is very sensitive to the loop concentration [ 1O]. Since experimental measurements of loop density at a sufficient range of temperatures and network densities are not available, the appropriate parameters for use in the void swelling calculations must be obtained by extrapolation f7]_ Obviously a finer de~nition of the appropriate parameters for a given set of conditions is necessary to enable realistic predictions of the void swelling to be made. In fig. 4 we show the total loop concentration as a function of dose for the same conditions as used previously i.e. 723K, pDo = log cm/ cm3,K=10-6d pa/ sec. It is clear that no si~i~c~t differences have been introduced by means of the grouping approximation and that even those values ofM which give rise to large discontinuities in the distribution function predict values for the total loop concentration which are in good agreement with the
r
1.4
x x
Single
O
t-4 = 10
A
M = 20
I
M= 50
I.4
M. 100
I
I
10-6
10-S
DOSE Id Fig. 4.
p.a.
1
equations
only
182
M.li. Haps
f A hierarchy
exact, single equation result. There is, of course a simple explanation and this is that the time scale for nucleation is short and therefore the equilibrium loop concentration is obtained whilst the loops are still relatively small and the eventual concentration is therefore determined by the first few rate equations. In this particular case the number of atoms/loop at the peak of the distribution is 6 at the end of nucleation. Consequently the grouping has little effect on the nucleation stage and only affects the growth kinetics. Therefore, in Kiritani’s case where 100 single equations were used there is good reason to believe that the loop concentrations predicted by his calculations will not be affected by any possible inadequacies in the grouping method. Since he also chose a value of A4 which is close to the optimum we would also expect that the kinetics of growth of the vacancy aggegates are also described reasonably accurately. We can conclude from this discussion that for each problem which is tackled using the group method it is necessary to first obtain prelimina~ information as to the number of single equations to be used to adequately describe the nucleation stage. Once this is established the simple rule given here for the definition of group widths will then ensure that the best level of approximation within the method is obtained. If the particular problem is such that significant nucleation is still continuiung when large complexes also exist then the transition from the grouped to un-grouped equations is more complicated since the continued production of small aggregates continually disturbs the first few grouped equations. This situation requires especial care to ensure that large errors are not introduced. Such a situation can be found during the nucleation of rare gas bubbles in nuclear fuels where nucleation contributes well into the growth stage, and we have considered this problem in detail and have shown by extensive calculations that the grouping method still represents a useful technique [ll].
4. Conclusions By means of a specific example of interstitial dislocation looos we have investieated the erouo method
ofme equations devised by Kiritani [l] for the solution of a large number of rate equations describing the nucleation and growth kinetics of point defect aggregates in irradiated materials. It has been shown that by making a suitable choice of group widths the discontinuities found by Koiwa [2] in the size distribution function can be minimized and that such overall quantities as the total concentration of interstitial loops are very insensitive to the details of the distribution function, as long as sufficient single equations are used to defin adequately the nucleation stage. It is concluded that the validity of the group metlod is established as a useful approximate scheme. As such, care will always be necessary in its application to ensure that the approximations inherent within it are acceptable to the problem in hand. If this is done then the method offers an extremely useful vehicle for the investigation of the complex physical interactions which govern the nucleation and growth kinetics of defect aggregates in irradiated materials.
References
[II M. Kiritani, J. Phys. Sot. Japan 3.5 (1973) 95. N. Yoshida and M. Kiritani, J. Phys. Sot. Japan 35 (1973) 1418. M. Kiritani and N. Yoshida, J. Phys. Sot. Japan 35 (1973) 83. VI M. Koiwa, J. Phys. Sot. Japan 37 (1974) 1532. 131 C. Wagner, Zeit-Fur Elektrochem 65 (1961) 581. [41 M.R. Hayns and Bullough, IAEA conference on ~ermodynamics of Reactor Fuel Elements, (1974) ISM/190. 151 M.R. Hayns, J. Nucl. Mat. 56 (1975) 267. M.R. Hayns and R. Perrin, Proc. Harwell Symposium on Void Swelling, Sept. 1974. [61 M.R. Hayns and E. Savino, to be published. 171 A.D. Brailsford and R. Bullough, J. Nucl. Mat. 44 (1972) 12. [81 R. Bullough and M.R. Hayns, to be published. [91 C.W. Gear, The automatic integration of stiff ordinary differential equations, Proc. of I.F.I.P. Congress 1968. [lOI R. Bullough, B. Eyre and K. Krishan, to appear in Proc. Roy. Sot. Harwell Res. Rep. (1975), AERE R79S2. 1111 M.H. Wood and M.R. Hayns, Harwell Res. Rep. (1975) AERE R8012. f121 P.R. Okamoto and S.D. Harkness, 3. Nucl. Mat. 48 (1973) 204.
ON THE GROUP METHOD FOR THE APPROXIMATE SOLUTION OF A HIERARCHY RATE EQUATIONS DESCRIBING NUCLEATION AND GROWTH KINETICS
OF
M.R. HAYNS L%eoreticafPhysics Division, AERE Hatwell, Oxfordshire OX1 I ORA, UK Received 6 October 1975
An approximate scheme devised by Kiritani for the solution of a large number of rate equations describing the nucleation and growth kinetics of point defect aggregates in irradiated materials is discussed in the light of serious objections raised by Koiwa as to the validity of such a scheme. By means of calculations on the particular problem of the nucleation and growth of interstitial dislocation loops it is shown that the grouping method is a useful approximation and suggestions are made as to the best choice of parameters for its application. Since it is an approximate scheme and is being applied to very complex systems, which can only yield to numerical methods, it is pointed out that great care is required to ensure that the method is not used inappropriately. Un schdma approximatif propose’par Kiritani pour re’soudre un grand nombre d’e’quations de vitesse decrivant les cine’tiques de germination et de croissance d’agregats de de’fauts ponctuels dans les ma&iaux irradie’s est discutia la Iumikre d’objections s&ieuses soulevees par Koiwa au sujet de la validite’d’un tel schima. A l’aide de calculs sur le probleme particulier de la ge~ination et de la croissance de boucles de dislocations interstitielles, il est montrg que ia m6 thode des groupes est une approximation utile et des suggestions sont faites B propos du meilleur choix des paramktres en vue de l’application de cette mithode. Puisqu’il s’agit d’un sche’ma approximatif et qu’il doit 5tre applique ‘a des systemes t&s complexes, ce qui peut conduire seulementi des me’thodes num&iques, on insiste sur ie fait qu’il faut utiliser cette me’thode avec grande precaution pour que celle-ci soit utilisCe ‘a propos. Ein NLherungsschema, das von Kiritani zur Losung einer grossen Zahl von Geschwindigkeitsgleichungen zur Beschreibung der Keimbildungs- und Wachstumskinetik von Punktdefektaggregaten in bestrahltem Material angegeben wurde, wird im Lichte der von Koiwa erhobencn gravierenden Einwande beziiglich der Giiltigkeit eines solchen Schemas diskutiert. Mit Hilfe von Rechnungen zum besonderen Problem der Keimbildung und des Wachstums von ZwischengitterVersetzungsringen wird gezeigt, dass die Gruppenmethode eine niitzliche Naherung ist. Zu deren Anwendung werden Annahmen hber die beste Auswahl der Parameter gemacht. Da das Schema Gherungsweise gilt und fir sehr komplexe Systeme Verwendung findet, die nur numerisch behandelt werden kBnnen, wird darauf hingeweisen, dass grosse Sorgfab erforderlich ist, urn sicher zu gehen, dass die Methode nur in geeigneter Weise angewandt wird.
it is necessary to consider clusters containing perhaps as many as lo6 defects. Obviously the full set of equations cannot realistically be solved and Kiritani constructed an approximate scheme whereby whole groups of cluster sizes were spanned by a single equation, with the rate constants chosen as being the average for the clusters covered by the group. In this way a scheme could be devised which would follow the growth of quite large clusters, at least in an approximate way. More recently Koiwa f2] has analysed the grouping approximation by means of comparisons between the analytical and numerical results for a
1. Introduction In a recent publication Kiritani [I] has discussed an approximate scheme for the solution of hierarchy of rate equations describing the nucleation and growth of vacancy clusters in irradiated metals. Such a scheme involves equations which describe the rate of accumulation and loss of single point defects by defect aggregates, consequently one rate equation is required to describe the history of a particular cluster size. In order to extend the calculation to clusters large enough to correlate with experimental results, 175
176
M.R. Hayns j A hierarchy
very simple model system which is ameanable to analytic solution. He raises serious doubts as to the validity of the grouping method. If one ever hopes to gain a quantitative understanding of the kinetics of nucleation and growth of defect aggregates in irradiated materials then methods similar to those described by Kiritani would seem to offer the only general course to follow. It is, of course, true that exact analytical solutions can be obtained in a few very special cases [3], but in general the complex nature of the systems involved in irradiation damage in structural materials precludes such treatments. Consequently it is important to establish the applicability of the grouping method in relation to problems of practical interest. In the past we have used the hierarchy method to investigate the nucleation and early growth of rare gas bubbles in operating nuclear fuels [4] and to study the problem of interstitial loops in M316 stainless steel [5], with particular emphasis on the void swelling problem. The above-mentioned work used a hierarchy of single equations, but recently we have included in the scheme a grouping method similar to that described by Kiritani. Since we are unavoidably concerned with the very complex systems found in realistic applications, we must turn to numerical methods and in particular the hierarchy plus grouping has proved to be extremely useful. Therefore it is of prime importance to establish the conditions under which such methods are realistic and where difficulties may be found. It is the purpose of this present work to establish the validity of these methods and to respond to the objections raised by Koiwa.
2. Rate theory model of interstitial and growth
loop nucleation
As Koiwa rightly points out, a detailed, mathematically sound, discussion of the very complicated equations describing even simple nucleation and growth problems is extremely difficult. Consequently we shall confine our comments to the numerical SO]Ution of a particular set of rate equations which are typical of those found in irradiation damage processes. In particular we shall consider the nucleation and growth kinetics of interstitial dislocation loops in M3 16 steel. This system was chosen simply because
of rate equations it forms the basis of a much more detailed study which will be presented elswhere [6]. Previously [5] we considered a similar problem, but used only a single group size hierarchy and we ignored the effects of the thermal emission of vacancies. In the present work both a group method and the thermal emission of vacancies have been included. The system of equations may be written
+ L,CzCv - DiPDOCi
2Ani + 1
Most of the quantities appearing in these equations have been defined previously [5] t. However, the definitions are given again in table 1 for convenience. In order to keep the system of equations relatively simple we have ignored the thermal emission of vacancies from all sinks other than the interstitial dislocation loops. The equations may be briefly interpreted as follows; the first two equations represent the balance of vacancies and interstitials between production, 7 We should note that the sink strength associated with the interstitial loop used in the calculations is a & where n is the number of interstitials in the loop. This sink strength was not found simply by counting the number of edge sites available, thereby giving a purely rate controlled growth, but by solving the diffusion to an effective sphere, the rate controlling process then appears as a geometrical factor which scales the fluxes into the loop according to the ratio of the areas of the effective sphere and the loop edge. The derivation of this sink strength has been given in detail previously [ 71 whilst a more complete analysis will be presented shortly [8].
MR. Hayns /A hierarchy of
rate
equations
Table 1 Definition of the quantities appearing in equations 1 I_____
QUWity
Definition
C”
Vacancy concentration t Interstitial concentration Concentration of interstitial loops containing n atoms Production rate of vacancies and interstitials (assumed equal) Recombination rate Rate constant for the production of loops containing n + 1 atoms by the addition of single atoms to loops of n atoms.
ci cn x K1 RN
Rn = 0.78fiVi Ln
Vi = attempt frequency; Em’ = interstitial migration energy; k = Boltzman constant; 7’ = temperature. Rate constant for the production of loops containing n - 1 atoms by vacancy absorption. L, = 0.78&v,
En
exp(-Emi/k73
exp(-~mv/k~
E ’ = vacancy migration energy. R?te constant for the production of loops containing II + 1 atoms by vacancy emission. E, = 0.78&v,
exp[-EmfjkT]
exp { -478~ + Felt a’/kT}
Emf = vacancy formation energy; 7s~ = stacking fault energy; Fe1 = the elastic energy of the loop [7] and stress effects have been ignored.
Di
D P> $i
Interstitial diffusion coefficient Vacancy diffusion coefficient. Concentration of network dislocations Width of the jth group, given by (rz._ 1 + 1)/M Number of equations defining sing ie loop sizes
.-
t all concen~ations
---
are fractional.
at rate K, and losses by recombination, vacancy absorption by loops, the loss of interstitials to existing loops and by the nucleation of new di-interstiti~ primary loops and the respective losses to an existing network dislocation structure of density pDo Cm/Cm3. Also included in the vacancy equation is the contribution by thermal emission from the existing interstitial loops, with the rate Ei. The remaining equation is representative of the hierarchy descfibing the individual cluster sizes. The terms on the right hand side are to be interpreted as the rate of production of loops of size nj by accumulation of interstitials, at rate Rj_ 1 and thermal emission of vacancies at rate Ei_1 from the next smallest loop size, the losses from the current size by further accumulation of interstitials and thermal emission, the loss by vacancy absorption and finally the backward reaction caused by vacancy absorption in the next highest loop size. The weighting factors in this equation correspond to the grouping method adopted and Ani is .the width of theith group, given by nj - nj_ 1. The coefficients
correspond exactly to these used by both Kiritani [l] and Koiwa [2]. A detailed discussion of the physical basis of these equations is not appropriate here, it has been given previously [.5] and, in relation to the thermal emission of vacancies will be discussed elsewhere [6]. The equations are based upon the assumption of homegeneous-binary nucleation, with such features as thermal dissociation and interstitial trapping not represented in the current equations. Equations (1) are, of course, solved numerically, in our case by means of a sophisticated subroutine based on Gear’s method, which is available as part of the HarweIl subroutine library and has been described in detail 191.
3. The results of the grouping method The grouping method, by collecting together loops of a size range into a single rate equation, the rate
M.R. Hayris /A hierarchy of rate equations
178
constants in which correspond to some average sink strength for the group, can be set up in a variety of ways. Kiritani adopted the scheme whereby the first 100 equations described only single loop sizes, whilst above this number the cluster sizes were obtained from the algorithm, Ani < O.O5(nj_l + 1). That is the size of cluster increases by 5% of the previous cluster size and as more and more equations are considered the width of the cluster can become very large indeed (i.e. Anj = lo4 at nj = 106). Koiwa, on the other hand chose several values of the number of un-grouped equations (30 and 70) whilst his Anj values were fixed at values ranging from 3 to 9. It is one of the purposes of this note to outline a method of specifying the appropriate choice of Anj for a given number of un-grouped equations which represents the optimum value within the general approximate scheme of the grouping method. We have used a’prescription for obtaining the width Anj as follows AHj
=(nj_l
+
1)/M
with M = 20 this expression reduces to the Kiritani equation. We shall show that simply by choosing M = number of un-grouped equations many of the objections raised by Koiwa can be eliminated. In fig. 1 we show the loop size distribution calculated using eqs. (1) but with no grouping. For our current purposes the solution of up to 200 equations is adequate and unless stated otherwise this is true for all of the calculations presented here. All of the results are for similar conditions, that is 723 K, network dislocation density of log cm/cm3 and a damage rate of 10S6 dpa/sec. Therefore these results are typical for solution treated M3 16 steel during fast neutron irradiation; the remaining materials parameters required are given in table 2. Fig. 1 is simply the loop size distribution at the longest time we could obtain using only 200 single equations of 0.54 X 10P4 d.p.a. (54 seconds at a dose rate of lop6 dpa/sec). In fact this is a relatively long time since nucleation * is completed very quickly, in about 10h5 d.p.a. in these calculations. The effects of introducing the group method are summarized in figs. 2 and 3. In fig. 2 we show loop * As previously [5] we are defining the end of nucleation as when the quantity (KzCf - LzC,C~)/C,(Z is less than 10e4.
R,C,
+ KICv)
L
8
12
16 20
2L
28
NO
OF
ATOMS/LOOP
32
36 Lo
LL
LE
52
56
60
Fig. 1. Distribution of loop sizes for parameters defined in tables 1 and 2 with no grouping, after a dose of 0.54 X 10e4 d.p.a. (54 seconds).
size distribution functions at rather early times (1 X lop5 dpa) for values ofM= 1, lo,20 and 50, with the hierarchy constrained to contain 10 single equations. The discontinuities in the distribution function at the cross over oetween single and grouped equations reported by Koiwa are apparent for M = 1, 20, and 50, whilst for N = 10 no discontinuity can be shown on the scale of this graph. We have connected the points with continuous lines, even though the functions are not defined except at the discrete positions, simply to highlight the differences between the curves. In fig. 3 similar results are given, but for a later time (3 X 1O-5 d.p.a.) when the peak of the distribution has passed through the single-group boundary. Here the divergence from the smooth curve obtained from the single equation
179
M.R. Hayns /A hierarchy of rate equations
Table 2 Parameters used for the calculations Quantity
Value
Interstitial diffusion coefficient, Di
10m3 exp(-0.2
Vacancy diffusion coefficient D,
0.6 exp(-1.4
Lattice parameter, a
3.63 X 10e8 cm
Network dislocation density pD”
lo9 cm/cm3
Temperature, T
723 K
Defect production rate, K
lo+ a/a.s/sec
eV/kr) cm’/sec eV/H’) cm2/sec
Vacancy formation energy Emf
1.6 eV
Stacking fault energy
9.4 X 1012 eV/cm’
Effective modulus b1 = p/(1 - u)
4.0 X 10” dynes/cm2
u-
10
9-o-
00-
-
M=l
0
M 10
I
M-20
&
M:50
=
r
0 M:lO
11
lo-
Full line smgle
A M:20 ..M=l
equations only
% x
Full lone single equations only
9-
70-
; ;
z
60-
b
8-
t Y 8
50-
< g F
LO-
7t i Y 8 w
Y E
30-
6-
5-
.?? 0
i J20-
=
L-
d II IO2
I
1
I I I 23L5670
I
II
I
NO
OF
’ 9
1”’
10
11
12
13
IL 15
16
ATOMS/LOOP
Fig. 2. Distributions of loop sizes after lo-’ d.p.a. for various values of M, the group width parameter. Appropriate values of An, the group width are given in table 3. The points for the grouped equations are shown at the mid point of the group. The value for M = 1 at n = 10 is out of scale (see text). All calculations contained hierarchies with 10 single equations.
2 !V l-
I
I
I
I
I
2
L
6
8
10
I
I
I
,
12 1L
16
NO
ATOMS
OF
I,,
1.3 20
I
I
22
2L 26
28
30 32 34
LOOP
Fig. 3. Distributions of loop sizes after 3 X lo-’ d.p.a. for similar conditions to fig. 2.
M. R. Ha.yns /A hierarchy of rate equations
180
solution when M = 20 is clearly shown, whilst the appropriate concentration for M = 1 is well out of the scale of this diagram. However, the curve for M= 10 does not show any discontinuity. In general the pattern noted by Koiwa of a shift in the peak of the distribution function to lower loop sizes whilst also showing a broadening is clearly seen in these figures. At 3 X 10m5 dpa the shift is in fact by one group, that is the peak occurs at 13.6 atoms/loop for M = 10 whilst for the single equations it lies at 16 atoms/loop. The reason for the discontinuities in the distribution function is that by taking too large a value of hj for the first group the number of loops being reduced in size per unit time by backward reactions into the last un-grouped equation is large, therefore artificially increasing the concentration of this loop size. The reverse situation is also true when the value of it is chosen to be very large, causing the width of the first grouped equation to be small. Therefore there is an optimum value of M, for which the transition is as smooth as possible. It is readily shown that for Ani_ 1 - Ani to be a minimum at the single-group boundary then M is simply the number of single equations. Therefore with IO single equations, M = 10 and ~~~=ll,l,~~~=l.l.Wehavecollectedtogether in table 3 values of group widths for IO and 100 single equations for various values of M Since the number of single equations is of necessity integral, we see that for 10 equations the smoothest transition gives the first group width as 10% wider than the last single equation. As the number of singles increases, this transition becomes smoother, at 100, there is only a 1% change. Kiritani’s scheme was an M value of 20, which gives a width of 5.05 to the first group. Our results are all for the rather coarse case of only IO single equations, however, even here the appropriate value of M = IO gives no discontinuities in the
distribution function which can be seen of the curves in figs. 2 and 3, in fact the numerical results indicate that such deviations are of the order of 0.1%. Consequently, it is suggested that the particular scheme adopted by Kiritani would not in fact lead to significant discontinuities in the distribution function. The results given by Koiwa using constant values of An of between 3 and 9 for 30 and 70 single equations would give rise to discontinuities since they are far from the optimum values of I .03 and 1.014 respectively, and these are indeed clearly shown in his results. It is important to emphasize that the accuracy of the grouping method is relative since any grouping will disturb the agreement with the exact result, which, following Koiwa we can assume to be given by the result of calculation using only single equations. However we are discussing an approximate scheme and as such some deviations from the ideal are inevitable, it is simply necessary to understand the inadequacies properly in order to make use of the scheme. In fact the details of the distribution function are not usually of prime importance in analysing irradiation damage experiments. Firstly because the experimental dif~culties ensure that any size distributions given are far coarser than the scale of deviations from the ideal as found by the grouping method (see for example the results of Okamoto and Harkness [ 121). Secondly the parameters of interest are not usually sensitive to the distribution function and thirdly such theoretical models as described here are inevitably based on approximations of the basic physical mechanisms which involves the use of sometimes very badly known quantities such as formation and activation energies, quantities which can introduce into the model a large degree of quantitative uncertainty. Furthermore, the models are based on current thinking as to the most important mechanisms involved,
Table 3 Values of A,, the width of the first group for 10 and 100 single equations for various values of M _. .___~_.. _._-._.- ______ ^ __~_ Number of single equations
1
5
10
20
50
..--.. 100
10
11
2.2
1.1
0.55
0.22
0.11
100
101
20.2
10.1
5.05
2.02
1.01
..___~_.._^.+M
1An
MR. Hayns /A hie~rch~ of rate equations
and others are not even included. In particular in the results presented here we have ignored interstitial trapping as an important mechanism in loop growth. Therefore the results can only be useful for correlations where it is absolutely certain that such processes are not important. When attempting to describe theoretically the complex interactions in a real material under irradiation we must be pragmatic and try to gain an understanding of the relationship between the various physical mechanics which eventually will lead to a fuller picture of the damage process and aid in the interpretation of the large number of experiments currently being performed. Therefore we need a suitable vehicle in which we can represent the physical mechanisms and investigate their ~ter-relation~ip. Such a vehicle is the rate theory coupled with the grouping method. The approximations inherent in it are understood reasonably well and until more sophisticated methods are produced it represents one of the very few ways of gaining understanding of the behaviour of these systems.
181
As an example of the type of data which can usefully be predicted by such a model we consider briefly the prediction of the total interstitial loop concentration. This is of special interest in the void swelling problem, where it has been shown that swelling is very sensitive to the loop concentration [ 1O]. Since experimental measurements of loop density at a sufficient range of temperatures and network densities are not available, the appropriate parameters for use in the void swelling calculations must be obtained by extrapolation f7]_ Obviously a finer de~nition of the appropriate parameters for a given set of conditions is necessary to enable realistic predictions of the void swelling to be made. In fig. 4 we show the total loop concentration as a function of dose for the same conditions as used previously i.e. 723K, pDo = log cm/ cm3,K=10-6d pa/ sec. It is clear that no si~i~c~t differences have been introduced by means of the grouping approximation and that even those values ofM which give rise to large discontinuities in the distribution function predict values for the total loop concentration which are in good agreement with the
r
1.4
x x
Single
O
t-4 = 10
A
M = 20
I
M= 50
I.4
M. 100
I
I
10-6
10-S
DOSE Id Fig. 4.
p.a.
1
equations
only
182
M.li. Haps
f A hierarchy
exact, single equation result. There is, of course a simple explanation and this is that the time scale for nucleation is short and therefore the equilibrium loop concentration is obtained whilst the loops are still relatively small and the eventual concentration is therefore determined by the first few rate equations. In this particular case the number of atoms/loop at the peak of the distribution is 6 at the end of nucleation. Consequently the grouping has little effect on the nucleation stage and only affects the growth kinetics. Therefore, in Kiritani’s case where 100 single equations were used there is good reason to believe that the loop concentrations predicted by his calculations will not be affected by any possible inadequacies in the grouping method. Since he also chose a value of A4 which is close to the optimum we would also expect that the kinetics of growth of the vacancy aggegates are also described reasonably accurately. We can conclude from this discussion that for each problem which is tackled using the group method it is necessary to first obtain prelimina~ information as to the number of single equations to be used to adequately describe the nucleation stage. Once this is established the simple rule given here for the definition of group widths will then ensure that the best level of approximation within the method is obtained. If the particular problem is such that significant nucleation is still continuiung when large complexes also exist then the transition from the grouped to un-grouped equations is more complicated since the continued production of small aggregates continually disturbs the first few grouped equations. This situation requires especial care to ensure that large errors are not introduced. Such a situation can be found during the nucleation of rare gas bubbles in nuclear fuels where nucleation contributes well into the growth stage, and we have considered this problem in detail and have shown by extensive calculations that the grouping method still represents a useful technique [ll].
4. Conclusions By means of a specific example of interstitial dislocation looos we have investieated the erouo method
ofme equations devised by Kiritani [l] for the solution of a large number of rate equations describing the nucleation and growth kinetics of point defect aggregates in irradiated materials. It has been shown that by making a suitable choice of group widths the discontinuities found by Koiwa [2] in the size distribution function can be minimized and that such overall quantities as the total concentration of interstitial loops are very insensitive to the details of the distribution function, as long as sufficient single equations are used to defin adequately the nucleation stage. It is concluded that the validity of the group metlod is established as a useful approximate scheme. As such, care will always be necessary in its application to ensure that the approximations inherent within it are acceptable to the problem in hand. If this is done then the method offers an extremely useful vehicle for the investigation of the complex physical interactions which govern the nucleation and growth kinetics of defect aggregates in irradiated materials.
References
[II M. Kiritani, J. Phys. Sot. Japan 3.5 (1973) 95. N. Yoshida and M. Kiritani, J. Phys. Sot. Japan 35 (1973) 1418. M. Kiritani and N. Yoshida, J. Phys. Sot. Japan 35 (1973) 83. VI M. Koiwa, J. Phys. Sot. Japan 37 (1974) 1532. 131 C. Wagner, Zeit-Fur Elektrochem 65 (1961) 581. [41 M.R. Hayns and Bullough, IAEA conference on ~ermodynamics of Reactor Fuel Elements, (1974) ISM/190. 151 M.R. Hayns, J. Nucl. Mat. 56 (1975) 267. M.R. Hayns and R. Perrin, Proc. Harwell Symposium on Void Swelling, Sept. 1974. [61 M.R. Hayns and E. Savino, to be published. 171 A.D. Brailsford and R. Bullough, J. Nucl. Mat. 44 (1972) 12. [81 R. Bullough and M.R. Hayns, to be published. [91 C.W. Gear, The automatic integration of stiff ordinary differential equations, Proc. of I.F.I.P. Congress 1968. [lOI R. Bullough, B. Eyre and K. Krishan, to appear in Proc. Roy. Sot. Harwell Res. Rep. (1975), AERE R79S2. 1111 M.H. Wood and M.R. Hayns, Harwell Res. Rep. (1975) AERE R8012. f121 P.R. Okamoto and S.D. Harkness, 3. Nucl. Mat. 48 (1973) 204.