New Phonon Dynamics in a Diamond Crystal during Post-Annealing

A new mesoscopic response in a crystal during post-annealing (PA) after ion implantation has been discovered. PA is an indispensable procedure for activating the implanted dopant center by retrieving the lost crystallinity of the host crystal, which usually elapses over periods longer than a few milliseconds. By using an empirical molecular-dynamics simulation, we examined what actually happened during extended PA for up to one nanosecond in case of sub-keV N2 beam implantation into bulk diamond at room temperature with annealing at 1000 K. While observing the changes in the long-range-order (LRO) parameter, a notable phase transition from the crystalline to the amorphous state (CA transition) occurred synergistically at approximately a few picoseconds (ps). After the CA transition, PA was initiated. A few tens of ps later, a sequence of restoration spikes was observed in LRO parameter with a regular time interval of 70–80 ps. Each sharp peak showed the spontaneous restoration of “amorphous → crystal → amorphous” transitions. This process may trigger a “critical state” that occurs during the slow and steady recovery of the lost crystallinity in the early stages of PA. It was found that annealing not only enhances the vibration amplitude of the individual atoms in the bulk but also causes phonon-assisted restoration in a mesoscopic space.


The Necessity of the Post-Annealing
Following the concept of ion-implantation doping, which was introduced by Shockley, investigations into crystalline defects have continuously improved the functional design of semiconductor devices, as can be seen in the annual reports of the International Technology Roadmap for Semiconductors (ITRS, 2011).The implanted atom (dopant) functions effectively after activation or post-annealing (PA) when lattice vibration mediates the restoration of damaged crystallinity and allocates some of the dopants at the lattice sites.Thus, PA is an indispensable treatment following ion implantation.As a matter of fact, PA does not reset the effect of doping and irradiation completely by an infinitely long period.PA is practically stopped when it reaches an acceptable level.If normal rapid thermal annealing (RTA) is used, PA requires a longer annealing time  a above 1 ms (Kwok et al., 2009;Race, 2010) and a higher annealing temperature T a above approximately 10 3 K.However, studies on methods that could practically replace traditional RTA are still ongoing (Larson, Williams, & Current, 2011).For example, microwave (Thompson et al., 2005) and pulsed-laser beams (Jones, Kuryliw, Murto, Rendon, & Talwar, 2000) have been proposed to shorten  a , and a cluster beam (Tanjo & Naito, 2011b) has been proposed to reduce T a .

The Characterization of the Post-Annealing
Such an acceptable level is characterized in two ways.One criterion for ion implantation is that the density of the activated dopant, measured by sheet resistance (Kwok et al., 2009), should approach a maximum level of 10 17 /cm 3 -10 21 /cm 3 , which is acceptable.For example, if used a B 18 beam with an impact energy of E 0 = 337 eV/ (B atom) into a crystalline Si (c-Si) (Onoda et al., 2010) with a conventional fluence of approximately 10 15 /cm 2 , this density corresponds to a state in which approximately 30% of the dopant is substituted at lattice sites.A simulation of this situation involves counting the number of annihilations between a dopant (at an interstitial site) and a vacancy (or a sink) using kinetic Monte-Carlo (KMC) simulation (Cox & Miller, 1965).Another criterion is to visualize the crystalline phase by capturing TEM snapshots at different times (Tanjo et al., 2011a).In a simulation for this situation, crystalline phases are identified by using snapshots, showing realistic atomic distributions (Park et al., 2005) or radial-distribution function [so-called g(r)] (Marqués, Pelaz, Hernández, Barbolla, & Gilmer, 2001).Although this approach is practical, it is difficult to identify when and how phase transition occurs.

The Numerical Analysis of the Post-Annealing
The interaction between projectiles and materials, such as that in case of ion-implantation doping, was first studied in terms of the primary and elemental processes that transfer energy or momentum mainly to target atoms or electrons (Sigmund, 2006), which occurs in the first stage.In the second stage, transferred energy or momentum is redistributed to the surrounding region in a solid (Nakagawa, 1997).This secondary stage continues for a very long time, in which the many-body interactions of atoms play a role over a mesoscopic region.For a long process that includes a huge span of the realm, a seamless class simulation (ITRS, 2011;Frangi, Cercignani, Mukherjee, & Aluru, 2008) is strongly required.
In the numerical computation of ion-implantation doping, the realm is described by the set of energy E, space X, and time t.Classical molecular dynamics (MD) is an atomistic simulation that solves the equation of motions of every atom by using empirical interatomic potentials (IPs) between atoms involved.Therefore the possible span for an MD simulation is E = 10 3 → 10 −1 eV, X = 10 −11 → 10 −7 m, and t = 10 −17 → 10 −9 s.The computational size of an MD frame defines the higher limits for X max and t max , whereas the actual atomistic view restricts the proper limits for X min and E max , which identify every atom and their noticeable displacements in a condensed material.In addition, the accuracy of an MD simulation depends explicitly on IPs and the time step used to solve differential equations, which define the lower limits for E min and t min , respectively.Herein, we consider the case of ion-implantation doping with a sub keV impact energy (Nakagawa, 2010a), for which the necessary span of the entire PA is E = 10 −2 -10 3 eV, X = 10 −11 -10 −3 m, and t = 10 −17 -10 −3 s for the complete simulation.Therefore, it is evident that an MD cannot deal with the entire process.
In order to examine the huge span of the realm for PA treatment, a typical simulation is composed of MD and KMC simulations (ITRS, 2011; de la Rubia, Soneda, Caturla, & Alonso, 1997), as is used for the first criterion of PA.MD describes the first stage and the earlier part of the secondary stage and focuses on monitoring all of the atoms in the potential field, modified by the chaotic collisions due to dopant implantation.KMC describes the later part of the secondary stage and focuses mainly on tracing the diffusive movement of the dopant on the basis of "stochastic random-walk model for diffusion" (Einstein, 1956).The dopants are assumed to cause stochastic annihilation by chance with vacancies produced by the previous collision stage along the path of a random walk in deformed potential field with many dips (Caturla et al., 2000).

The Information to Be Extracted from the Simulation
The connecting factor between these two methodologies (MD and KMC) is the diffusion constant D, which can be obtained using an MD simulation by averaging the moving distance of the atoms over the entire MD box throughout the long time span (Nakagawa, Iwatani, & Betz, 2006).The factor D is described by the trap-depth (activation energy) of a dip in the potential field and attempt frequency for jumping from a trapped dip to one of the neighboring quasi-stable sites, as expressed by the Arrhenius rule.Typically, attempt frequency is presumed to be constant, i.e., 10 13 /sec (Caturla et al., 2000) or 10 12 /sec (Latz et al., 2012) and corresponds to the frequency of lattice vibration phonons; thus, activation energy is very important for determining the value of D. Many physical studies designed to improve PA have been practically helpful (Larson et al., 2011;Thompson et al., 2005;Jones et al., 2000;Tanjo & Naito, 2011b), and KMC simulations can reproduce some of these successful results using various parameters, including the value of D. However, the KMC method cannot essentially predict how and why a newly proposed heating procedure for PA will be effective in a damaged crystal.The problem is that the KMC method is an optimization method at present (Latz et al., 2012).
In order to examine the mechanism of phase change, a single quantity must be monitored throughout computation for clearly identifying the moment of transition, as is used for the second criterion of PA.There are a few candidates (quantities) for this purpose (Nakagawa, 2009).The range-order parameter of a crystal is one of the best candidates.For this approach, we have developed a crystallographic analysis called the Pixel-Mapping (PM) method (Nakagawa, 2002(Nakagawa, , 2009)), which can be used in conjunction with MD (Nakagawa & Betz, 2001).

Methods
In the present study, we examined the functioning of the trigger prior to the final settlement during PA for the restoration of a damaged crystal by focusing on the long-range-order (LRO) parameter of the crystal.As an example, molecular-nitrogen (N 2 ) implantation into a IIa-type diamond, which can lead to the formation of a nitrogen-vacancy-nitrogen (NV-N) center (Neumann et al., 2009), was investigated using a simulation followed by a subsequent physical experiment (Gaebel et al., 2006).The NV-N center may be a key component of "solid-state quantum computers at room temperature (RT)" (Wootters, 1998) when the intra-pair distance (R N-N ) between the two N atoms is approximately 2.3-2.6 nm (Hanson, Mendoza, Epstein, & Awschalon, 2006).For the physical experiment (Gaebel et al., 2006), a IIa-type diamond was bombarded by an N 2 beam with an impact energy of E 0 = 6-14 keV/(N atom) at RT followed by PA at T a = 1073 K.We have, however, reproduced the NV-N center even with the much lower E 0 of a sub-keV N 2 beam; thus, it can be concluded that annealing is more essential than the impact energy (E 0 ) for producing the NV-N center (Nakagawa, Kanda, & Betz, 2010b).
So far, when an MD simulation has been used to study phonon dynamics, MD extracts the so-called dispersion relation (Henry & Chen, 2008) or examines the thermal transport (Cahill et al., 2003).Both approaches require the presumption of the presence of a global phonon.However, the irradiation effect does not uniformly expand over the entire space of the irradiated crystal; therefore, it is more realistic to assume that a local phonon would be stimulated in a mesoscopic region.When a diamond is excited by a femtosecond laser, a coherent phonon is excited (Ishioka, Hase, Kitajima, & Petek, 2006) and is damped with showing "collapse and revival" within 10 ps (Misochko, Hase, Ishioka, & Kitajima, 2004).While this phonon is an optical coherent phonon, this time-span demonstrates that an atomistic simulation using MD may be ideally suited for the examination of phonon dynamics if searched up to the order of a few nanoseconds (ns), which corresponds to t max of a MD.Therefore, the present study focused on how the heat of PA affects such a mesoscopic phonon effectively, which may mediate the restoration of damaged crystallinity.
PM identifies each NV-N center and the current phase using the value of LRO parameter (Nakagawa, 2002(Nakagawa, , 2009) ) in conjunction with an MD simulation.The MD code traces all the positions of the atoms through high-precision computing.However, in the final comparison of observed results of crystalline defects with experiments, as long as concerning the order or disorder, a coarse-grained view is more significant than a precise one for the identification of the phase (Nakagawa, 2011b).Thus, to investigate the phonon-dynamics observed in an MD supercell, LRO parameter was monitored.In the PM method, LRO parameter is defined by the fraction of appropriate atoms at the correct lattice sites in a ternary crystal with cubic symmetry (Nakagawa, 2007).Namely, CA transition is recognized by the change in this LRO parameter from the crystalline state (LRO = 1) to amorphous state (LRO = 0) (Nakagawa, 2007).A diamond lattice is composed of two FCC lattices.Thus, in the PM method, when LRO = 1, two FCC lattices together form a diamond lattice, whereas when LRO = 0, neither FCC lattice is formed (Nakagawa, 2002(Nakagawa, , 2007)).The amorphous-crystalline (AC) transition is the reverse of this order.

Results
As can be observed in Figure 1, LRO parameter apparently does not change drastically immediately after the formation of the collision cascade (within 10 −15 s), which is produced by the N 2 beam impact with an energy of 200 eV/(N atom) at RT.Because the frequency of lattice vibration is approximately 0.1 ps, the time intervalt for monitoring LRO parameter can be 0.1 ps or above.As is shown in the upper abscissa of Figure 1, different values of t are considered for three different ranges, i.e., t = 0.1 ps, 1 ps, and 10 ps for time t < 3 ps, t < 20 ps, and t ≥ 20 ps, respectively.
In this case, the so-called critical-slowing-down state toward amorphization, which occurred prior to the formation of a new amorphous phase, took approximately 2 ps.After the CA transition the system was thermally stabilized at RT. Then PA at T a = 1000 K was initiated at t = 2.4 ps (indicated by the vertical arrow in Figure 1).The first LRO spike appeared at 60 ps.
The spike represents the instantaneous recovery of crystallinity (jumps in LRO parameter from 0 to 1) and was repeated every 70-80 ps.These restoration spikes also appeared when 300 K < T a <1000 K and with different energies [E 0 = 100, 300, and 500eV/(N atom)].The time interval between restoration spikes was observed to clearly depend more on T a than on E 0 , and it increased as T a decreased.Similar restoration spikes were also obtained when c-Si was bombarded with sub-keV boron clusters (B n ; n ≦ 6) and annealed at T a = 1323 K (Nakagawa, 2011c).
Because of the time limit for MD simulations, LRO parameter was computed up to 1 ns.Following the first restoration spike, the amorphous state remained for a long time, e.g.several tens of ps, and was then followed by the next sharp spike, and subsequently.The period of each spike was less than 10 ps.Such an oscillatory behavior of LRO parameter repeated until 1 ns, as shown in Figure 2. The repetition of the restoration spikes occurred almost regularly, which is a sign of restoration.Many symbols on the ordinate in Figure 2 indicate that CA transition occurred rapidly, as shown clearly in Figure 1, in which time is indicated using a logarithmic scale.The sequence of restoration spikes suggests that a similar change might occur temporally from place-to-place in the period of critical state toward the end of PA; otherwise, only the amorphous state (LRO = 0) would have been continuously observed until the AC transition finally once occurred, as expected from the reverse image of Figure 1.

Conclusions
In conclusion, after the CA transition at a few ps (Nakagawa et al., 2010b;Nakagawa, 2011aNakagawa, , 2012)), the PA treatment caused a sequence of restoration spikes in LRO.These spikes indicate the presence of alternative amorphous and crystalline phases, as if a local phonon carrying the synergistic movement of atoms might have been stimulated.Nonetheless, this phonon mode reacts much more slowly than a coherent optical phonon that exhibits oscillatory damping with a much finer frequency in the order of less than 0.1 ps (Misochko et al., 2004), which is known as the "collapse revival", in the same diamond (Ishioka et al., 2006).In addition, results further suggest the possibility of proposing a coarse viewpoint for the time-scale after the CA transition shown in Figure 1.The restortion spikes indicate the slow alternation between the amorphous and crystalline states prior to complete restoration.
Therefore, we conclude that PA not only enhances vibration amplitude of individual atoms but also prepares the crystal slowly and synergistically for complete restoration.The critical state stage involves the repetition of the spontaneous stimulation of local phonons, leading to temporary restoration.Furthermore, the possibility of a more reasonable class simulation than MD + KMC, with which using a seamless simulation with a coarse viewpoint of both time and space can be developed, can be demonstrated.
Figure 1.CA transition as an effect of ion doping on the LRO parameter up to 300 ps after N 2 implantation [200 eV/(N atom)] into a diamond at the bulk temperature T = RT.The vertical arrow indicates the initiation of heating with T a = 1000 K, whereas the horizontal arrows denote the different monitoring steps (t)