Given the short time scale and small volume associated with atomic displacement cascades, it is not currently possible to directly observe their behavior by any avail­able experimental method. Some of their characteristics have been inferred by experimental techniques that can examine the fine microstructural features that form after low doses of irradiation. The experimental work that provides the best estimate of stable Frenkel pair produc­tion involves cryogenic irradiation and subsequent annealing while measuring a parameter such as electri­cal resistivity.26,27 Less direct experimental measure­ments include small angle neutron scattering,28 X-ray scattering, positron annihilation spectroscopy, 0 and field ion microscopy. 1 More broadly, transmission elec­tron microscopy (TEM) has been used to characterize the small point defect clusters such as microvoids, dislocation loops, and stacking fault tetrahedra that are formed as the cascade collapses.32-36

The primary tool for investigating radiation dam­age formation in displacement cascades has been computer simulation using MD, which is a computa­tionally intensive method for modeling atomic sys­tems on the time and length scales appropriate to displacement cascades. The method was pioneered by Vineyard and coworkers at Brookhaven National Laboratory,37 and much of the early work on atomistic simulations is collected in a review by Beeler.38 Other methods, such as those based on the BCA,2 , 1 have also been used to study displacement cascades. The binary collision models are well suited for very high — energy events, which require that the interatomic potential accurately simulate only close encounters between pairs of atoms. This method requires sub­stantially less computer time than MD but provides less detailed information about lower energy colli­sions where many-body effects become important. In addition, in-cascade recombination and clustering can only be treated parametrically in the BCA. When the necessary parameters have been calibrated using the results of an appropriate database of MD cascade results, the BCA codes have been shown to reproduce the results of MD simulations reasonably well.39,40

A detailed description of the MD method is given in Chapter 1.09, Molecular Dynamics, and

will not be repeated here. Briefly, the method relies on obtaining a sufficiently accurate analytical inter­atomic potential function that describes the energy of the atomic system and the forces on each atom as a function of its position relative to the other atoms in the system. This function must account for both attractive and repulsive forces to obtain the appropri­ate stable lattice configuration. Specific values for the adjustable coefficients in the function are obtained by ensuring that the interatomic potential leads to reasonable agreement with measured material para­meters such as the lattice parameter, lattice cohesive energy, single crystal elastic constants, melting tem­perature, and point defect formation energies. The process of developing and fitting interatomic poten­tials is the subject of Chapter 1.10, Interatomic Potential Development. One unique aspect arises when using MD and an empirical potential to investi­gate radiation damage, viz. the distance of closest approach for highly energetic atoms is much smaller than that obtained in any equilibrium condition. Most potentials are developed to describe equilibrium con­ditions and must be modified or ‘stiffened’ to account for these short-range interactions. Chapter 1.10, Interatomic Potential Development, discusses a common approach in which a screened Coulomb potential is joined to the equilibrium potential for this purpose. However, as Malerba points out,41 criti­cal aspects of cascade behavior can be sensitive to the details of this joining process.

When this interatomic potential has been derived, the total energy of the system of atoms being simulated can be calculated by summing over all the atoms. The forces on the atoms are obtained from the gradient of the interatomic potential. These forces can be used to calculate the atom’s accelerations according to Newton’s second law, the familiar F = ma (force = mass x acceleration), and the equations of motion for the atoms can be solved by numerical integration using a suitably small time step. At the end of the time step, the forces are recalculated for the new atomic positions and this process is repeated as long as necessary to reach the time or state of interest. For energetic PKA, the initial time step may range from ^1 to 10 x 10~18 s, with the maximum time step limited to ~1—10 x 10~15s to maintain acceptable numerical accuracy in the integration. As a result, MD cascade simulations are typically not run for times longer than 10—100 ps. With periodic boundary conditions, the size of the simulation cell needs to be
large enough to prevent the cascade from interacting with periodic images of itself. Higher energy events therefore require a larger number of atoms in the cell. Typical MD cascade energies and the approxi­mate number of atoms required in the simulation are listed in Table 1. With periodic boundaries, it is important that the cell size be large enough to avoid cascade self-interaction. For a given energy, this size depends on the material and, for a given material, on the interatomic potential used. Different inter­atomic potentials may predict significantly different cascade volumes, even though little variation is even­tually found in the number of stable Frenkel pair.42 Using a modest number of processors on a modern parallel computer, the clock time required to com­plete a high-energy simulation with several million atoms is generally less than 48 h. Longer-term evolu­tion of the cascade-produced defect structure can be carried out using Monte Carlo (MC) methods as discussed in Chapter 1.14, Kinetic Monte Carlo Simulations of Irradiation Effects.

The process of conducting a cascade simulation requires two steps. First, a block of atoms of the desired size is thermally equilibrated. This permits the lattice thermal vibrations (phonon waves) to be established for the simulated temperature and typi­cally requires a simulation time of approximately 10 ps. This equilibrated atom block can be saved and used as the starting point for several subsequent cas­cade simulations. Subsequently, the cascade simula­tions are initiated by giving one ofthe atoms a defined amount of kinetic energy, EMD, in a specified direc­tion. Statistical variability can be introduced by either
further equilibration of the starting block, choosing a different PKA or PKA direction, or some combination of these. The number of simulations required at any one condition to obtain a good statistical description of defect production is not large. Typically, only about 8-10 simulations are required to obtain a small stan­dard error about the mean number of defects pro­duced; the scatter in defect clustering parameters is larger. This topic will be discussed further below when the results are presented. Most of the cascade simulations discussed below were generated using a [135] PKA direction to minimize directional effects such as channeling and directions with particularly low or high displacement thresholds. The objective has been to determine mean behavior, and investigations of the effect of PKA direction generally indicate that mean values obtained from [135] cascades are repre­sentative of the average defect production expected in cascades greater than about 1 keV.43 A stronger influ­ence of PKA direction can be observed at lower ener­gies as discussed in Stoller and coworkers.4 ,

In the course of the simulation, some procedure must be applied to determine which of the atoms should be characterized as being in a defect state for the purpose of visualization and analysis. One approach is to search the volume of a Wigner-Seitz cell, which is centered on one of the original, perfect lattice sites. An empty cell indicates the presence of a vacancy and a cell containing more than one atom indicates an interstitial-type defect. A more simple geometric criterion has been used to identify defects in most of the results presented below. A sphere with a radius equal to 30% of the iron lattice parameter is

centered on the perfect lattice sites, and a search similar to that just described for the Wigner-Seitz cell is carried out. Any atom that is not within such a sphere is identified as part of an interstitial defect and each empty sphere identifies the location of a vacancy. The diameter of the effective sphere is slightly less than the spacing of the two atoms in a dumbbell interstitial (see below). A comparison of the effective sphere and Wigner-Seitz cell approaches found no significant difference in the number of stable point defects identified at the end of cascade simulation, and the effective sphere method is faster computationally. The drawback to this approach is that the number of defects identi­fied by the algorithm must be corrected to account for the nature of the interstitial defect that is formed. In order to minimize the lattice strain energy, most interstitials are found in the dumbbell configuration; the energy is reduced by distributing the distortion over multiple lattice sites. In this case, the single interstitial appears to be composed of two interstitials separated by a vacancy. In other cases, the interstitial configuration is extended further, as in the case of the crowdion in which an interstitial may be visualized as three displaced atoms and two empty lattice sites. These interstitial configurations are illustrated in Figure 4, which uses the convention adopted throughout this chapter, that is, vacancies are displayed as red spheres and interstitials as green spheres. A simple postprocessing code was used to determine the true number of point defects, which are reported below.

Most MD codes describe only the elastic colli­sions between atoms; they do not account for energy

loss mechanisms such as electronic excitation and ionization. Thus, the initial kinetic energy, EMD, given to the simulated PKA in MD simulations is more analogous to Td in eqn [2] than it is to the PKA energy, which is the total kinetic energy of the recoil in an actual collision. Using the values of EMD in Table 1 as a basis, the corresponding EPKA and nNRT for iron, and the ratio of the damage energy to the PKA energy, have been calculated using the procedure described in Norgett and coworkers.19 and the recommended 40 eV displacement thresh­old.16 These values are also listed in Table 1, along with the neutron energy that would yield EPKA as the average recoil energy in iron. This is one-half of the maximum energy given by eqn [1]. As mentioned above, the difference between the MD cascade energy, or damage energy, and the PKA energy increases as the PKA energy increases. Discussions of cascade energy in the literature on MD cascade simulations are not consistent with respect to the use of the term PKA energy. The third curve in Figure 3 shows the calculated number of Frenkel pair predicted by the NRT model if the PKA energy is used in eqn [2] rather than the damage energy. The difference between the two sets of NRT values is substantial and is a measure of the ambiguity associated with being vague in the use of terminology. It is recom­mended that the MD cascade energy should not be referred to as the PKA energy. For the purpose of comparing MD results to the NRT model, the MD cascade energy should be considered as approxi­mately equal to the damage energy (Td in eqn [2]).

In reality, energetic atoms lose energy continu­ously by a combination of electronic and nuclear reactions, and the typical MD simulation effectively deletes the electronic component at time zero. The effects of continuous energy loss on defect produc­tion have been investigated in the past using a damp­ing term to slowly remove kinetic energy.46 The related issues of how this extracted energy heats the electron system and the effects of electron-phonon coupling on local temperature have also been exam — ined.47-50 More recently, computational and algorith­mic advances have enabled these phenomena to be investigated with higher fidelity.51 Some of the work just referenced has shown that accounting for the electronic system has a modest quantitative effect on defect formation in displacement cascades. For exam­ple, Gao and coworkers found a systematic increase in defect formation as they increased the effective electron-phonon coupling in 2, 5, and 10 keV cascade simulations in iron,50 and a similar effect was reported
by Finnis and coworkers.47 However, the primary physical mechanisms of defect formation that are the focus of this chapter can be understood in the absence of these effects.