A remarkable result by Ducastelle and Cyrot — Lackmann10 relates the tight-binding local density of states to the local topology. If we describe the density of states in terms of its moments where the pth moment is defined by

Ep„(E) dE

— 1

and recall that by definition

„(-E)=E d(E — Ei)

i

where i labels the eigenvalues, we get 1

Ep E d(E — E,) dE = E Ep = Tr[Hp] where H is the Hamiltonian matrix written on the basis of the eigenvectors. But, the trace of a matrix is invariant with respect to a unitary transformation, that is, change of basis vectors to atomic orbitals i. Therefore,

mt = Tr[Hf] = $>*]* ^mP0

ii

A sum of local moments of the density of states mp^. These diagonal terms of Hp are given by the sum of all chains of length p of the form HjjHjkHki. .. H„j. These in turn can be calculated from the local topol­ogy: a prerequisite for an empirical potential. They consist of all chains of hops along bonds between atoms which start and finish at і (e. g., see Figure 5). By counting the number of such chains, we can build up the local density of states.

Unfortunately, algorithms for rebuilding DOS and deducing the energy using higher moments tend to converge rather slowly, the best being the recursion method.1

The zeroth moment simply tells us how many states there are.

The first moment tells us where the band center is. Taking the band center as the zero of energy, the second moment is as follows:

m2i) = [H2]„ = 5>h ] = E h(rj )2 [3]

./ ./

where h is a two-center hopping integral, which can therefore be written as a pairwise potential.

This result, that the second moment of the tight- binding density of states can be written as a sum of pair potentials, provides the theoretical underpinning for the Finnis-Sinclair (FS) potentials. Referring back to the rectangular band model, we can take the second moment of the local density of states m2i) as a measure of the bandwidth.

This gives the relationship between cohesive energy, bandwidth, and number of neighbors (z,). In

the simplest form W, / that is, the band energy is proportional to the square root of the number of neighbors.

Note that this is only a part of the total energy due to valence bonding. There is also an electrosta­tic interaction between the ions and an exclusion — principle repulsion due to nonorthogonality of the atomic orbitals — it turns out that both of these can be written as a pairwise potential V(r).

The moments principle was laid out in the late 1960s.12 To make a potential, the squared hopping integral is replaced by an empirical pair potential f(rp), which also accounts for the prefactor in eqn [4] and the exact relation between bandwidth and second moment. Once the pairwise potential V(ry) is added, these potentials have come to be known as FS potentials.13

Ecoh = X V (П, )-X. Xfj) N

ijij where Vand f are fitting functions.

Further work14 showed that the square root law held for bands of any shape provided that there was no charge transfer between local DOS and that the Fermi energy in the system was fixed. For bcc, atoms in the second neighbor shell are fairly close, and are normally assumed to have a nonzero hopping integral.

Notice that the first three moments only contain information about the distances to the shells of atoms within the range of the hopping integral. Therefore, a third-moment model with near neighbor hopping could not differentiate between hcp and fcc (in fact, only the fifth moment differentiates these in a near­neighbor hopping model!). This led Pettifor to con­sider a bond energy rather than a band energy, and relate it to Coulson’s definition of chemical bond orders in molecules.15 Generalizing this concept leads to a systematic way of going beyond second moments and generating bond order potentials.

One can investigate the second-moment hypothesis by looking at the density of states of a typical transition metal, niobium, calculated by ab initio pseudopotential plane wave method, Figure 6, and comparing it with the density of states at extremely high pressure. The similarity is striking: as the material is compressed, the band broadens but the structure with five peaks remains unchanged. The s-band is displaced slightly to higher energies at high pressure, but still provides a

low, flat background, which extends from slightly

below the d-band to several electron volts above.

1.10.6.2 Key Points

• In a second-moment approximation, the cohesive (bond) energy is proportional to the square root of the coordination.

• Other contributions to the energy can be written as pairwise potentials.

The use of MD to simulate primary damage forma­tion has become widespread and relatively mature. In addition to the research involving metals discussed above, the approach has also been applied to common

structural ceramics11-14,142 and ceramics ofinterest to the nuclear fuel cycle.15,143,144 However, there are a

number of areas that require further research. Some of these have to do with the most basic aspect of MD simulations, that is, the interatomic potentials that are used. In addition to the Finnis-Sinclair potential for iron that was used as a reference case in this chapter, results from several other iron potentials were mentioned. The choice of potential is never an obvious one, and there have been few studies to systematically compare them. In one of the studies mentioned in Section 1.11.3, the details of how one joins the equilibrium part of the potential to a screened Coulomb potential to account for short — range interactions were shown to significantly influence cascade evolution and defect formation.41 Although a clear difference has been demonstrated, there is no clear path to determining what constitutes the ‘correct’ way to join these potentials.

In the case of iron and other magnetic elements such as chromium, research to address the issue of how magnetism may influence defect formation and behavior has only recently begun.145-147 The effect may be modest in the ballistic phase of the cascade when energies are high, but magnetism must cer­tainly influence the configuration and properties of stable defects. Magnetic effects may also determine critical properties of interstitial clusters such as their migration energy and primary diffusion mechanism, which will strongly influence the nature of radiation damage accumulation. As the standard density func­tional theory fails to fully account for magnetic effects, further developments in electronic structure theory are required in order to provide data for fitting new and more accurate potentials.

The interaction between the atomic and elec­tronic systems has largely been neglected in most of the work discussed above. This may impact the results of MD simulations in at least two ways. First, energetic atoms lose energy in a continuous slowing down process that involves both the elastic collisions MD currently models and electronic excitation and ionization between these elastic collisions with lattice atoms. Because of the energy dependence of elastic scattering cross-sections, neglecting the energy loss between atomic collisions could lead to more diffuse cascades and higher predicted defect survival. The second effect is related to inaccuracies in temper­ature when energy transfer between the electronic and atomic systems (electron-phonon coupling) is neglected. To first order, the atoms remain hotter when energy loss to the electron system is not accounted for. Given the temperature dependence of defect survival and defect clustering discussed above, this clearly has the potential to be significant in any one material. In addition, as electron-phonon coupling varies from one material to another, its neglect may obscure real differences in defect forma­tion between materials.

Finally, the issue of rare events requires more investigation. The need to carry out sufficient simu­lations at a given condition to obtain an accurate estimate for mean behavior was emphasized in the chapter. However, it may be that rare events are also important for the prediction of radiation damage accumulation at longer times or higher doses. If nucleation of extended defects is difficult, which is typically the case at higher temperatures and lower point defect supersaturations, rare events that seed the microstructure with large clusters may largely control the process. One example of a potentially significant rare event is provided by the work of Soneda and coworkers.148 They carried out one hun­dred 50-keV simulations at 600 K to obtain a good statistical description of defect formation at this con­dition. In one of these simulations, 223 stable point defects were created, which was much greater than the average of 130 defects. In addition, a <100> vacancy loop containing 153 vacancies was created. The diameter of the loop was about 2.9 nm, which is large enough to be visible by TEM. The impact of the one-in-a-hundred type events should not be underestimated without further study.

Acknowledgments

The author would like to acknowledge the fruitful collaboration and discussions on cascade damage for many years with Drs. David Bacon and Andrew Calder (University of Liverpool, UK), Lorenzo Malerba (SCK/CEN, Mol, Belgium), and Yuri Osetskiy (ORNL). He was first introduced to MD cascade simulations by Drs. Alan Foreman (deceased) and William Phythian during a short-term assign­ment at the AEA Technology Harwell Laboratory, UK, in the summer of 1994. Various aspects of the author’s research discussed in this chapter were sup­ported by the Division of Materials Sciences and Engineering, and the Office of Fusion Energy Sciences, US Department of Energy and the Office of Nuclear Regulatory Research of the US Nuclear Regulatory Commission at the Oak Ridge National Laboratory under contract DE-AC05-00OR22725 with UT-Battelle, LLC.

Damage accumulation in pure metals during irradi­ation primarily takes place in the formation and evolution of vacancy and SIA-type defects. At tem­peratures higher than recovery stage III, which is the main interest for practical purposes, vacancy clusters normally take the form of voids that result in the change of a volume, that is, swelling. Owing to limitations of space, in the following section we focus only on a description of void evolution.

1.13.5.2.1 Void swelling

The solution obtained from eqns [44] depends on the irradiation temperature. Temperatures below recovery stage II will not be considered here. At temperatures smaller than that corresponding to the recovery stage III, when vacancies are immobile and the interstitials are mobile, the concentration of vacancies will build up. At some irradiation dose, the vacancy concentration will become high enough that mutual recombination of PDs may become the dominant mechanism of the defect loss, thus controlling defect accumulation. In this case,

 

DvCvk2 = Pc(x)fc(x) x2

 

[83]

the dose dependence of PD concentrations can be calculated analytically104	formed. In the following discussion we concentrate on the irradiation doses beyond the transient period, which are of more practical interest. If only voids and edge dislocations are present in the system, and mutual recombination and thermal emission of vacancies from voids and dislocations are both negligible, the balance equations for the concentrations of vacancies and SIAs, Cv and Ci, are given by G — k^DyCy — ZdypdDvCy = 0 G — k2 Di Ci — zfpdDiCi = 0 [93]
At temperatures higher than that corresponding to recovery stage III, both vacancies and SIAs are mobile. Hence, after a certain time of irradiation, called the ‘transient period’, their concentrations reach a steady state. A comprehensive analysis of the time (irradia­tion dose) dependence of PD concentrations for different sink strength can be found in Sizmann.9 The dose dependence of PD concentrations and void swelling obtained by the numerical integration of ME73 is presented in Figure 4. As can be seen, the vacancy supersaturation, (DvCv — DiCi)/DvCVq, becomes positive when the PD concentrations reach steady state and this gives rise to void growth. Also, note that in the transient regime only divacancies are
[94]
Figure 4 Dose dependences of the concentrations of point defects, void swelling, vacancy supersaturation, and void number density calculated in the framework of FP3DM by numerical integration of the master equation, eqn [18]. From Golubov and Ovcharenko.73
dd Zi__ Zv Zd
[97]
Bd
The maximum value of the ratio in the right-hand side of eqn [96] is 1/4, when the sink strengths of voids and dislocation are equal to each other, k/ = Z^Pd. Thus the maximum swelling rate is dS =Bd df max 4 It is easy to show that the swelling rate described by eqn [96] depends only weakly on the variation of the sink strength of voids and dislocations: a differ­ence of an order of magnitude results in a decrease of the swelling rate by a factor of 3 only. To obtain the steady-state swelling rates of ~1% per NRT dpa, which are observed in high-swelling

fcc materials, one would need the bias factor to be about several percent. Data on swelling in electron — irradiated metals resulted in Bd ~ 2 — 4% for the fcc copper24,105,106 (data reported by Glowinski107 were used in Konobeev and Golubov106), ^2% for pure Fe-Cr-Ni alloys,108 and orders of magnitude lower values for bcc metals (e. g., swelling data for molybdenum1 ). Because the electron irradiation produces FPs, it is reasonable to accept these values as estimates of the dislocation bias.

Note that the first attempt to determine Bd by solving the diffusion equations with a drift term deter­mined by the elasticity theory for PD-dislocation interaction as described in Section 1.13.5 showed that the bias is significantly larger than the empirical estimate above. Several works have been devoted to such calculations,96,110-113 which predicted much higher Bd values, for example, ~15% for the bcc iron and ^30% for the fcc copper. With these bias factors, the maximum swelling rates based on Bd/4 should be equal to about 4% and 8% per dpa but such values have never been observed. An attempt to resolve this discrepancy can be found in a recent publication.114

Surprisingly, the steady-state swelling rate of ~1% per NRT dpa has been found in neutron — (and ion-) irradiated materials, for example, in various stainless steels, even though the primary damage in these cases is known to be very different and the void swelling should be described in the framework of the PBM, which gives a rather different description ofthe process. An explanation of this is proposed in Section 1.13.6.

Thanks to fundamental advances, coupled with the development of efficient algorithms and fast compu­ters, the PF technique has become a very powerful and versatile tool for simulating phase transforma­tions and microstructural evolution in materials, as illustrated in this chapter. This technique provides simulation tools that are complementary to atomistic models, such as molecular dynamics and lattice Monte Carlo simulations, and to larger scale approaches, such as finite element models. With some modifications, it can also be employed for mate­rials subjected to irradiation.

In the case of materials subjected to irradiation, specific issues need to be addressed to fully realize the potential of PF modeling. First, a proper descrip­tion of point defects and atom transport requires mobility matrices (or tensors) that capture the kinetic coupling between these different species. In particu­lar, the models reviewed in this chapter do not account for the correlated motion of point defects and atoms, thus leading to unphysical correlation factors in the mobility coefficients. These correlation effects, however, play an essential role in phenomena

400
300
300
200	200
100
100
200
100
200
100
300	300
400
400 0	0
Figure 15 Defect concentration fields corresponding to Figure 14: (a) free vacancies CV, (b) free interstitials Cfnt A + CP, t B (A and B atoms), (c) clustered vacancies Cc;V, and (d) clustered interstitials Cc;A + Cc;B (A and B atoms). Reproduced from Badillo, A.; Bellon, P.; Averback, R. S. to be submitted.

such as irradiation-induced segregation and precipi­tation and are thus required in PFMs aiming for system-specific predictive power. Second, it remains challenging to include in a PFM all the elements of the microstructure relevant to evolution under irra­diation, namely point-defect clusters, dislocations, grain boundaries, and surfaces, although it has been shown here that models handling adequately a subset of these microstructural elements are now becoming available. Third, the numerical integration of the evolution equations is more challenging than for conventional PFMs in the sense that the continuous defect production, as well as the large difference in vacancy and interstitial mobility, usually prevents the use of long integration time steps, even in coarse microstructures. Finally, materials under irradiation constitute nonequilibrium systems that are quite sensitive to the amplitude and the structure of fluc­tuations, in particular the fluctuations resulting from

point defect and point-defect cluster production, and from ballistic mixing of species. A self-consistent and tractable PFM that would include both ther­mal and irradiation-induced fluctuations is still miss­ing. Such a model would be very beneficial for the study of microstructural evolution under irradiation, especially that involving the nucleation of new phases, defect clusters, or gas bubbles see Chapter 1.13, Radiation Damage Theory; Chapter 1.14, Kinetic Monte Carlo Simulations of Irradiation Effects; and Chapter 1.09, Molecular Dynamics.

1.18.3.2.1 Local chemical potential

The thermodynamic state equation defines a chemi­cal potential of species i as the partial derivative of the Gibbs free energy G of the alloy, with respect to the number of atoms of species i, that is, Ni. The resulting chemical potential is a function of the temperature and molar fractions (also called concen­trations) of the alloy components, Ci = Ni/N, N being the total number of atoms. TIP postulates that local chemical potentials depend on local concentrations via the thermodynamic state equa­tion. A chemical potential gradient V/ of species i is then equal to

V/L

kBT

where Ci is the local concentration of species i. In a binary alloy, concentration gradients of the two com­ponents are exactly opposite. The chemical potential gradient of component i is then proportional to the concentration gradient:

V/i kB T
JV = -£ LVp(Xp — xv) [6] p
[10]

Under irradiation, diffusion is controlled by both vacancies and interstitials. The flux of interstitials is also deduced from the atomic fluxes:

Ji = X Jp

p

where F, is called the thermodynamic factor. Fur­thermore, the Gibbs-Duhem relationship58 leads to interdependent chemical potential gradients:

Ck V/k = 0

k = 1 ,r

where the sum runs over the number of species. There­fore, in a binary alloy there is one thermodynamic factor left:

Vm

kBT

where F = Fa = FB. Note, that an alloy at finite temperature contains point defects. They are currently assumed to be at equilibrium with the local alloy com­position, with the local chemical potential equal to zero. When calculating the thermodynamic factor, point defect concentration gradients are neglected. During irradiation, although point defects are not at equilibrium, one assumes that eqn [12] continues to be valid.

Under irradiation, additional driving forces are involved. They correspond to the gradients of vacancy and interstitial chemical potentials, which are usually written in terms of their equilibrium concentrations CVq and Cjeq respectively:

mV = kBT ln(CV/cVq) and mI = kBT ln(Q/Cjeq) [13] leading to an expression of the associated driving force11: —V CV — XVAV CA with xVA

The interstitial driving force has the same form, except that letter V is replaced by letter I. Note, that the equilibrium point defect concentrations may vary with the local alloy composition and stress. Although the variation of the equilibrium vacancy concentration is expected to be mainly chemical, the change of the elastic forces due to a solute redistribution at sinks should not be ignored for the interstitials.11 Due to the lack of experimental data, Wolfer11 introduced the equilibrium vacancy concen­tration as a contribution to a mean vacancy diffusion coefficient expressed in terms of the chemical tracer diffusion coefficients. Composition-dependent tracer diffusion coefficients could then account for the change of equilibrium vacancy concentration, with respect to the local composition.

Within the framework of the TIP, a thermody­namic factor depends on the local value but not on the spatial derivatives of the concentration field. The use of this formalism for continuous RIS models deserves discussion. Indeed, a typical RIS profile covers a few tens of nanometers so that the cell size used to define the local driving forces does not exceed a few lattice parameters. Such a mesoscale
chemical potential is expected to depend not only on the local value, but also on the spatial derivatives of the concentration field. According to Cahn and Hilliard,60 the free-energy model of a nonuniform system can be written as a volume integral of an energy density made up of a homogeneous term plus interface contributions proportional to the squares of concentration gradients. Thus, all contin­uous RIS models that are derived from TIP retain only the homogeneous contribution to the energy density and cannot reproduce interface effects and diffuse-interface microstructures. In particular, an equilibrium segregation profile near a surface is pre­dicted to be flat.

The creation of a stable FP requires that a lattice atom receives an energy greater than Tm, which is the minimum displacement energy. This value has been determined experimentally in many materials by measuring the change in some physical property, such as electrical resistivity or length change, as a function of maximum recoil energy of a target atom. Such experiments are practical only for electron irradiations for which recoil energies can be kept low, but with the irradiation particles still penetrating deeply into, or through, the specimen. Typical values are shown in Table 2.

As a crystal is not homogeneous, the threshold energy depends on the crystallographic direction in which the knock-on atom recoils. The anisotropy of the threshold energy surface has been mapped out in various crystals by measuring the production rate of defects as a function of both the electron energy, near threshold, and the orientation of single crystalline

specimens with respect to the electron beam direc­tion.15’16 The total cross-section for FP production rate is given by the expression

2p я/2

dff(g2; £1) df 2

d02 2p [19]

0 0

v(02, F2; T)dd2

where the subscripts 1 and 2 refer to incoming elec­tron and recoiling ion, respectively, and 0b F1, 02, F2 are the polar and azimuthal angles of the electron beam relative to the crystal axis; 02, f2, are these same angles relative to the beam direction; v is the anisotropic damage function. Near threshold, v = 1 for T> Tn, and 0 for T< Tm. By measuring the production rate for many sample orientations and energies, the damage function can be obtained using eqn [19], although various approximations are required in the deconvolution. The results are illustrated in Figure 10 for Cu.17 It is noteworthy that the minimum threshold energy is located in the vicinity of close-packed directions. This is also true for bcc metals. The anisotropy reflects the basic mechanism of defect production, viz., replacement collision sequences (RCSs), which had been identified by molecular dynamics simulations as early as 1960.1

(a)

The primary knock-on atom in an RCS recoils in the direction of its nearest neighbor, (110) in fcc crystals, and replaces it, with the neighbor recoiling also in the (110) and replacing its neighbor. A vacancy is left at the primary recoil site, and an interstitial is created at the end of the sequence. Replacement sequences are the most efficient way to separate the interstitial far enough from its vacancy, ^2-3 interatomic spacings, for the FP to be stable. While the lengths of these sequences are still debated, it is clear that the mechanism results in both defect pro­duction and atomic mixing. For neutron irradiations, higher energy recoils are numerous, and the average displacement energy, Ed, becomes more relevant for calculations of defect production (see eqn [1]). This value, which can be obtained by averaging over the threshold displacement energy surface, is usually dif­ficult to determine experimentally. A rough estimate, however, can be obtained from, Td ~ 1.4 Tm in fcc metals and 1.6ЕП in bcc metals.19

1.09.1 Introduction

A concept that is fundamental to the foundations of Comprehensive Nuclear Materials is that of microstruc­tural evolution in extreme environments. Given the current interest in nuclear energy, an emphasis on how defects in materials evolve under conditions of high temperature, stress, chemical reactivity, and radiation field presents tremendous scientific and technological challenges, as well as opportunities, across the many relevant disciplines in this important undertaking of our society. In the emerging field of computational science, which may simply be defined as the use of advanced computational capabilities to solve complex problems, the collective contents of Comprehensive Nuclear Materials constitute a set of compelling and specific materials problems that can benefit from science-based solutions, a situation that is becoming increasingly recognized.1-4 In discus­sions among communities that share fundamental scientific capabilities and bottlenecks, multiscale modeling and simulation is receiving attention for its ability to elucidate the underlying mechanisms governing the materials phenomena that are critical to nuclear fission and fusion applications. As illu­strated in Figure 1, molecular dynamics (MD) is an atomistic simulation method that can provide details of atomistic processes in microstructural evolution.

As the method is applicable to a certain range of length and time scales, it needs to be integrated with other computational methods to span the length and time scales of interest to nuclear materials.9

The aim of this chapter is to discuss in elemen­tary terms the key attributes of MD as a principal method of studying the evolution of an assembly of atoms under well-controlled conditions. The intro­ductory section is intended to be helpful to students and nonspecialists. We begin with a definition of MD, followed by a description of the ingredients that go into the simulation, the properties that one can calculate with this approach, and the reasons why the method is unique in computational materials research. We next examine results of case studies obtained using an open-source code to illustrate how one can study the structure and elastic properties of a perfect crystal in equilibrium and the mobility of an edge dislocation. We then return to Figure 1 to provide a perspective on the potential as well as the limitations of MD in multiscale materials modeling and simulation.

Having deduced the functional form of the potential from first principles, it remains to choose the fitting functions and fit their parameters to empirical data. Most papers simply state that ‘the potential was parameterized by fitting to….’ The reality is different.

Firstly, one must decide what functions to use for the various terms. Here, one may be guided by the physics (atomic charge density tails in EAM, square root embedding in FS, Friedel oscillations), by the anticipated usage (short-range potentials will speed up MD, and discontinuities in derivatives may cause spurious behavior), or simply by practicality (Can the potential energy be differentiated to give forces?).

Secondly, one must decide what empirical data to fit. Cohesive energies, elastic moduli, equilibrium lattice parameters, and defect energies are common choices. Accurate ab initio calculations can provide further ‘empirical’ data, notably about relative struc­tural stability, but now increasingly about point defect properties. (It is, of course, possible to calcu­late all the fitting data from ab initio means. Potentials fitted in this way are sometimes referred to as ab initio. While this is pedantically true, the implication that these potentials are ‘better’ than those fitted to exper­imental data is irritating.)

Ab initio MD can also give energy and forces for many different configurations at high temperature. Force matching44 to ab initio data is one of best ways to produce huge amounts of fitting data. There have been many attempts to fully automate this process, but to date, none have produced reliably good poten­tials. This is in part because of the fact that although MD only uses forces (differential of the potential), many essential physical features (barrier heights, structural stability, etc.) do depend on energy, which in MD ultimately comes from integrating the forces. Small systematic errors in the forces, which lead to larger errors in the energies, can then cause major errors in MD predictions. Furthermore, if the poten­tial is being used for kinetic Monte Carlo, the forces are irrelevant.

By using least squares fitting, all the data may be incorporated in the fit, or some data may be fitted exactly and others approximately. However, since the main aim of a potential is transferability to different cases, the stability of the fitting process should be checked. The best way to proceed is to divide the empirical data arbitrarily into groups for fitting and control, to fit using only a part of the data, and then to check the model against the control data. This process can be done many times with different divisions offitting and control. Any parameter whose value is highly sensi­tive to this division should be treated with suspicion.

Structural stability is the most difficult thing to check, since one simply has to check as many struc­tures as possible. In addition to testing the ‘usual suspects,’ fcc, bcc, hcp, A15, o-Ti, MD, or lattice

dynamics can help to check for mechanical instability of trial structures.

Interaction between a 1/2(110){111} edge disloca­tion and a 1/2(110) SIA loop in Ni was first studied by Rodney and Martin,70 and then followed a series of more detailed investigations involving glissile and sessile SIA loops and screw and edge disloca — tions.16,70-74 Reaction R3 was observed with small loops having bL = 1/2(110) intersecting the disloca­tion slip plane.70,75 Glissile clusters were attracted by the dislocation core and absorbed there ather — mally, creating a pair of superjogs. Superjog seg­ments have different structure, as depicted in Figure 14, and to accommodate this a few vacancies are formed, as in the case of the R3 reaction with an SFT (see Figure 12). Each superjog has different mobility, for example, the Lomer-Cottrell segment on the left of Figure 14(b) has high Peierls stress. Therefore, although the jogged dislocation con­tinues gliding under applied stress, it has a signifi­cantly lower velocity because it now experiences higher effective phonon drag.71

An interesting example of an R2 reaction was observed during interaction between a 1/2(110) screw dislocation and a 1/3(111) Frank loop.16 If the loop is not too large and the dislocation is not too fast, the dislocation can absorb the whole loop into a helical turn (see stages in Figure 15). The turn can expand only along the dislocation line and, therefore, if applied stress is maintained, the helix constricts; finally a perfect DL, with the same b as the dislocation, is released as the screw breaks away.

Note that the same mechanism occurs when a screw dislocation intersects a glissile 1/2(110) loop with Burgers vector different from that of the dislocation, that is, absorption into a helical turn and release of a loop with the same Burgers vector as the disloca — tion.76 In the case when the Frank loop is large or the dislocation (either screw or edge) is fast, the dislocation simply shears the loop and creates a step on its surface (see Figure 8 in Rodney16). The probability of this reaction is higher for an edge dislocation, whereas transformation of the loop into a perfect loop is more probable for the screw dislocation.73

The strengthening effect due to dislocation-SIA loop interactions can be significant, especially when a helical turn is formed on a screw dislocation. The total contribution of dislocation-SIA loop interaction to the flow stress under irradiation has not been considered so far because of the large number of possible reactions and the sensitivity of their out­comes in terms of mechanism and strengthening effect to parameters such as interaction geometry, loop size and Burgers vector, strain rate, and temper­ature. An estimate of the flow stress contribution for the case of reaction R2 in Ni was made in Rodney and Martin.75

1.12.4.2.2.1 Body-centered cubic metals

Many atomic-scale simulations of interaction between a 1/2(111){110} dislocation and an SIA loop in bcc metals have been made, particularly using EAM potentials for Fe. They have shown that when the defects come into contact, a new dislocation segment is formed with one of two possible Burgers vectors 1/2(111) and (100), and further evolution depends on features such as loop size, b, position, dislocation character, temperature, and strain rate. Thus, many different outcomes have been observed by atomic scale modeling, but here we can present only a few common examples.

Reaction R3 is most common for small SIA loops, that is, loop absorption on the dislocation line resulting in a pair of superjogs. It occurs when the loop is initially below the dislocation slip plane (tension region of the dislocation strain field) and bL = 1/2(111) is inclined to the slip plane. An important feature of the interaction is the ability of such loops to glide quickly toward the core of the approaching dislocation. SIA loops up to 5 nm in diameter have been simulated.77-80 When small loops (a few tens of SIAs) reach the core, they are fully absorbed athermally creating a double superjog of the size equivalent to the number of SIAs in the

Figure 15 Unfaulting of a 1/3(111) Frank loop by interaction with a screw dislocation in Ni at 300 K and transformation into a helical turn on the dislocation line: (a) initial configuration, (b) first cross-slip, (c) and (d) successive cross-slip events,

(e) configuration at the end of unfaulting, (f) configuration after relaxation and elongation of the helical turn. From Rodney, D. Nucl. Instrum. Meth. Phys. Res. B 2005, 228, 100. Copyright (2005) with permission from Elsevier.

loop. This process does not pin the dislocation and these loops are very weak obstacles to dislocation glide.77 Large loops (more than ^100 SIAs) also glide to make contact with the dislocation line but are not absorbed athermally. Instead, a new disloca­tion segment is formed due to the following energet­ically favorable Burgers vector reaction between the dislocation and loop:

2[111]-1[1T1] = [010] [2]

Five stages of the interaction are presented in Figure 16 for a 5-nm loop containing 331 SIAs. The Burgers vectors are indicated and Figure 16(b) corresponds to the occurrence of the reaction of eqn [2]. The new segment with b = [010] cannot glide in the dislocation slip plane (110) and therefore acts as a strong obstacle to further glide of the dislo­cation. Under increasing applied strain, the dislocation bows out until two long segments of 1/2[111] screw dislocations are formed as shown in Figure 16(c). Figure 16(d) shows the same configuration in [111]

projection. High stress at junctions connecting the dislocation line and loop, the remainder of the origi­nal loop, and the new segment induces the latter to slip down on the (101) plane, and glide of this [010] segment over the loop surface results in the following reaction with the remaining loop:

2[1І1] + [010] = 2[111] [3]

This concludes with the formation of a pair of super­jogs on the original dislocation and results in com­plete absorption of the 5 nm loop. Large loops are strong obstacles in this reaction, stronger than voids with the same number of vacancies.4,80

It should be noted, however, that the interaction just described depends rather sensitively on temper­ature because of the low mobility of screw segments in the bcc metals. Cross-slip of the screw segments is required to allow the [010] segment to glide down. Simulations at T = 100 K showed that although the stages in Figure 16(a)-16(c) still occurred, the [010] segment did not glide under the strain rate imposed in MD and the screw dipole created by the bowing dislocation was annihilated without the loop being transformed according to eqn [3]. The resulting reac­tion was of type R1, for both dislocation and loop were unchanged after dislocation breakaway.78 More details and examples can be found.4,77-82

Competition between reactions R1, R2, and R3 for a 1/2(111){110) edge dislocation and 1/2(111) and (100) SIA loops has been considered in detail.83,84 We cannot describe all the reactions here but some pertinent features are underlined; note that the favorable Burgers vector reaction between a 1/2(111) dislocation and a (100) loop results in a 1/2(111) segment, for example,

І[111]+[Г00]=І[Г11] [4]

Thus, a perfect loop with bL = [100] can be converted into a sessile complex of 1/2 [111] and 1/2 [111] loop segments joined bya [100] dislocation segment.83 Similar conjoined loop complexes were observed in simulations of interactions between two glissile 1/2(111) loops85,86 in bcc Fe. Competition between reaction R3 on one side and R1 and R2 on the other was discussed in Bacon and Osetsky. 2 The earlier conclusion that small loops (< 1 nm) can be absorbed easily and not present strong obstacles to dislocation glide was confirmed. The strength and reaction mechanism for larger loops de­pends on their size and the loading conditions. At low T< 100K, both 1/2(111) and (100) loops are strong obstacles that are not absorbed by 1/2(111) dislocations. As with obstacles that result in the true Orowan mecha­nism, the dislocation unpins by recombination of the screw dipole, and the critical stress is determined by the loop size, similar to eqn [1].At higher Tand/or low strain rate (<107s_1), the mechanism changes from Orowan — like bypassing to complete loop absorption, irrespective of bL. The absorption mechanism for R3 involves propa­gation of the reaction segment over the loop surface. This requires cross-slip ofthe arms ofthe screw dipole drawn out on the pinned dislocation and involves dislocation reactions. Thermally activated glide and/or decomposi­tion of the pinning segment, in turn, depends on loop size, temperature, and bL. Therefore, the obstacle strength of a 1/2(111) loop is, in general, higher than that of a (100) loop because a (100) dislocation segment associated with the former (see eqn [3], Figure 16 and related text) is much less mobile than a 1/2 (111) segment involved in reactions with the latter (eqn [4]).

Much less is known about interaction mechanisms involving screw dislocations in bcc metals. We are aware of only two studies reported to date that con­sidered 1/2(111) and (100) SIA loops.87,88 The main feature in these cases is the ability of a 1/2(111) screw dislocation to absorb a complete or part loop into a temporary helical turn before closing the turn by bowing forward and breaking away. As in the case of fcc metals described above, this provides a power­ful route to reorienting the Burger vectors of differ­ent loops to that of the dislocation. An example of such a reaction is visualized in . Other

reactions observed in these studies87,88 include a reforming of the original glissile loop into a sessile complex of two segments having different b, as described above (reaction R2), and complete restora­tion of the initial loop (reaction R1).

As is clear from the examples discussed above, the obstacle strength associated with different mechan­isms can depend strongly on parameters such as loop size and Burgers vector, interaction geometry, Є and T Unlike the situation revealed for inclusion­like obstacles discussed in Section 1.12.4.1, there does not appear to be a simple correlation between size and strength for loops. Some data related to par­ticular sets of conditions can be found in Terentyev et a/.,83,84,88 as well as a comparison of the obstacle strength of voids and DLs.84 There are cases when loop strengthening is compatible to or even exceeds that of voids containing the same number of point defects, making DLs an important component of radiation-induced hardening.

Successful applications of the PBM have been limited to low irradiation doses (< 1 dpa) and pure metals (e. g., copper). There are two apparent problems preventing the general application of the model at higher doses.

1.13.6.3.1 Swelling saturation at random void arrangement

The PBM predicts a saturation of void size.30 This originates from the mixture of 1D and 3D diffusion- reaction kinetics under cascade damage conditions, the assumption lying at the heart of the model. More specifically, it stems from the fact that the interaction cross-section with a void is proportional to the void radius, R, for 3D migrating vacancies and to R2 for 1D diffusing SIA clusters. As a result, above some critical radius, the latter becomes higher than the former and the net vacancy flux to such voids is negative. In contrast, experiments demonstrate indefinite void growth in the majority of materials and conditions.31- 34 An attempt to resolve this contradiction was under­taken by including thermally activated rotations of the SIA-cluster Burgers vector23,25,27; but it has been shown by Barashev et a/.25 that this is not a solution. Thus, the PBM fails to account for this important and common observation, that is, the indefinite void growth under cascade irradiation. A way ofresolving this issue is discussed in Section 1.13.7.

1.13.6.3.2 Absence of void growth in void lattice

Another problem of the PBM is that it fails to explain swelling saturation at rather low swelling levels (approximately several percent) observed in void lat­tices. In fact, it even predicts an increase in the swelling rate when a random void arrangement is changed to that of a lattice.25 This is because the free channels between voids along close-packed directions, which are formed during void ordering, provide escape routes for 1D migrating SIA clusters to disloca­tions and GBs, thereby allowing 3D migrating vacan­cies to be stored in voids. A possible explanation of the problem is discussed in a forthcoming paper by Golubov eta/.37