- Email: [email protected]
Melting and shock wave creation in uranium oxide due to Coulomb explosion after a pulsed ionization
Nuclear Instruments and Methods in Physics Research B 358 (2015) 65–71
Contents lists available at ScienceDirect
Nuclear Instruments and Methods in Physics Research B journal homepage: www.elsevier.com/locate/nimb
Melting and shock wave creation in uranium oxide due to Coulomb explosion after a pulsed ionization Zhongyu Li a, Di Chen b, Lin Shao b,⇑ a b
Fundamental Science on Nuclear Safety and Simulation Technology Laboratory, Harbin Engineering University, Harbin 150001, China Department of Nuclear Engineering, Texas A&M University, College Station, TX 77843, USA
a r t i c l e
i n f o
Article history: Received 23 February 2015 Received in revised form 26 April 2015 Accepted 28 April 2015
Keywords: Molecular dynamics simulations Shock wave Ion irradiation Uranium oxide
a b s t r a c t By means of molecular dynamics simulations, we study the effects of pulsed ionization in uranium oxide (UO2), which occurs when UO2 is bombarded with swift ions or fission fragments. A general formula is developed to predict melting radius under various conditions due to electron stripping and Coulomb explosion (CE). A critical density model is suggested in which the melting volume is proportional to ionization period, if the period is above a critical value. The maximum melting radius depends on the time period of structural relaxation above the melting temperature, which increases with increasing initial substrate temperatures due to a lower heat dissipation rate. Furthermore, shock waves are observed to emit from CE core but the kinetic energy wave peak exists only in U sublattices. The absence of kinetic energy waves in O sublattices is explained by their relatively higher thermal vibration which cancels the work done from the compression waves. Ó 2015 Elsevier B.V. All rights reserved.
1. Introduction In numerous previous studies, it has been observed that swift ions can create ion tracks in insulators, semiconductors, and polymers. The fundamental mechanisms of track formation, however, are still subject to some debate [1–5]. Coulomb explosion (CE) and thermal spike model were often compared as competing mechanisms. Bringa and Johnson supported CE mechanism [1], while many other studies favored thermal spike model [2–5]. In the thermal spike model, energetic electrons are created under a high electron stopping power and subsequent electron-lattice coupling leads to melting and lattice disorders [6–8]. In the CE model, swift ions cause strong electron ionization of target atoms, leading to partial electron depletion along an ion track. In semiconducting or insulating substrates, the lack of free electrons makes it difficult to quickly recover the charge neutrality through electron back diffusion. With an ionization period sufficiently long (>a few fs), Coulomb force among the positively charged target atoms lead to CE [9], through the following stages: first, Coulomb potentials quickly convert to kinetic energies of atoms during the ionization pulse. Second, within a factor of 1 ps, kinetic energies are shared among atoms and local thermal equilibrium is reached. This process also leads to local melting. Third, the melting zones are cooled ⇑ Corresponding author. E-mail address: [email protected] (L. Shao). http://dx.doi.org/10.1016/j.nimb.2015.04.077 0168-583X/Ó 2015 Elsevier B.V. All rights reserved.
down and form either amorphous or polycrystalline cores, depending on rates of heat dissipation to the surrounding medium. Swift ion induced ion track formations can be used in microtechnology and nanotechnology for a wide range of applications [10–12]. Ion track formation in uranium oxide (UO2) is particular important for light water reactors, due to the facts that fission fragments from neutron reactions, in a typical kinetic energy of 100 MeV, may lead to structural changes which cause thermal and mechanical property degradation, and influence reactor performance. Experiments on both accelerator based swift ion irradiations and reactor irradiation have observed ion track formation in UO2 [13,14]. Computer simulations have been used to understand CE effects [15–18]. In the present study, we preformed systematic molecular dynamics (MD) simulations, in order to understand CE mechanism. Many previous studies have used MD simulations to understand thermal spike model. Modeling studies on CE phenomenon, however, are very limited. The present study has no intention to evaluate the competing effects between Coulomb explosion model and thermal spike model, but expect to benefit the development of a comprehensive model to possibly combine both effects. 2. Modeling procedure The present study is aimed to understand atomic scale details of CE in UO2, and to identity key parameters which influence the
66
Z. Li et al. / Nuclear Instruments and Methods in Physics Research B 358 (2015) 65–71
microstructure changes after CE. Instead of directly linking the modeling to specific fission production irradiation, CE is introduced by a pulsed ionization. The ionization time period and the substrate temperatures are used as variables for systematic investigation. Another reason to introduce a pulsed ionization is to alleviate issues from molecular dynamics (MD) simulations. MD simulations cannot handle electron subsystems. MD simulations are performed by using the LAMMPS (Large-scale Atomic Molecular Massively Parallel Simulator) code [19]. The interatomic potentials were described by a partly ionic Busing-Ida-type using potential parameters proposed by Yakub et al. [20]. Many interatomic potentials have been developed to describe UO2, as reviewed by Govers et al. [21]. These potentials include the two major types [20]: one has charges fixed by the chemical compositions and the other, the Busing-Ida type [22,23], has fractional effective charges. Introducing effective charges provides an important fitting parameter to match with experimental observation. Most importantly, the potential used in the present study is able to accurately predict melting and existence of a pre-melting structural transition [20]. As to be discussed, the transition from solid to liquid is one important phenomenon introduced by CE. The used potentials include Coulomb interactions and short-range interactions. The short-range interactions consist of O–O van der Waals attractions, and U–O covalent bonding represented by the Morse function and overlap (exchange) repulsion. The unit cell has a dimension of about 22 nm 22 nm 5.5 nm and contains 192,000 atoms. The cell is [1 0 0] oriented. All sides take periodic boundary conditions, including the top and the bottom of the cell. The geometry of the cell is schematically shown in Fig. 1. After thermal relaxation at 300 K for 50 ps, an ionization pulse was introduced in a cylindrical region of radius 20 Å. From 0 to 10 Å, half of the U and O atoms within the ionized region are
stripped by one electron per atom, f.g. U4+ ? U5+, to simulate CE. From 10 to 20 Å, the stripped electrons in the CE region are re-distributed, forming an electron rich cylindrical shell surrounding the electron-stripped core. These electrons are added to local U/O atoms so their Columbic force contributions are included in MD simulations. Without such treatment, charge neutrality cannot be preserved, and the accumulated electric fields from positively charged core are high enough to introduce significant Coulomb force and create abnormal cell vibrations. The approach of redistribution of stripped electrons near the CE core effectively creates charge screening effect to eliminate the long range Coulomb force from the core. At the end of the ionization, these redistributed negative charges, together with positive charges introduced in the CE core, are removed simultaneously under assumption of quick charge recovery. For the present study, ionization period range from 5 fs to 10 fs. The time selections consider the facts that electron ionization, and electron–electron thermalization occur in a typical time scale of 1 fs [2]. Once the ionization pulse ends, the pair Coulomb potential was removed. For the rest of the discussion, the time t = 0 refers to the moment when the ionization starts. The cell is divided into cylindrical slabs and closest atomic separated distance (CSD) of each slab, based on the average value of the first neighbor distance within the slab, is extracted to plot CSD radial distributions. The melting radius can be quantitatively determined since CSD values of melting UO2 are larger than that of amorphous UO2. And both of them are larger than that of crystalline UO2. It needs to point out that a defective crystalline UO2 layer may create high CSD values but they are typically less than that of amorphous UO2 and do not show saturation. For a local molten core surrounded by defective zones, radial CSD values are featured with a step height peak corresponding to the molten core and a decreasing tail corresponding to surrounding defective zone. Therefore, we use half height of
Fig. 1. Kinetic energy distributions in UO2 for an ionization period of (a) 5 fs and (b) 10 fs. Ionization is introduced in the cylindrical region of 10 Å radius. The initial temperature is 300 K. The bottoms are enlarged images showing atom disorder levels. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
67
Z. Li et al. / Nuclear Instruments and Methods in Physics Research B 358 (2015) 65–71
CSD maximum to define the melting radius. This makes it possible to track liquid-to-solid changes as a function of time. The time step for structural relaxation prior to the ionization is 1 fs. The time step during the ionization is 0.1 fs, and that for structural evolution after the ionization is 1 fs. 3. Results and discussion Fig. 1a and b compare snapshots of MD simulations as a function of time, for ionization time period 5 fs and 10 fs, respectively. Color bar represents the kinetic energy of atoms. At time t = 0 s, half of randomly selected atoms within a cylinder of radius r = 10 Å are charged by stripping one electron per atom. Their charge statuses remain until the end of the ionization periods. As shown by the kinetic energy distributions, at the end of the ionization, the CE process does not introduce significant kinetic energies within the core region, due to Coulomb force cancelation from all directions. Instead, the atoms at the track edge gain the highest kinetic energy due to biased force along the radial distributions. At 1 ps, localized energy equilibration is reached and the whole core regions have relatively uniform kinetic energy distributions. This stage is also featured with local melting, with atoms in the core losing short range orders. At longer time t = 20 ps, quenching of the molten core is finished. Kinetic energies are distributed over much wider region. As shown by the enlarged image in the box marked region, there is significant difference in crystalline orders between the short ionization (5 fs) and long ionization (10 fs). For the short ionization, thermal quenching ends with resolidification of the molten core into the defective but crystalline structure. In comparison, the long ionization ends with fully amorphous core. In the enlarged images, orange color refers to U atoms, white color refers to charged U with one electron stripped during the ionization pulse. Green color refers to O atoms and black color refers to electron-stripped O atom during the ionization pulse. Atomic distributions show that both electron-stripped U and O atoms have
3.0
4.0
(a) O-O, τ=5 fs
(d) U-U, τ=5 fs 0 fs 5 fs 1 ps 3 ps 10 ps 20 ps
0 fs 5 fs 1 ps 3.9 3 ps 10 ps 20 ps
2.9
Closest seperation distance (Angstrom)
diffusion from the molten core but O atoms diffuse more due to their relatively small mass. These observations suggest that the ionization period plays an important role to determine whether CE core region will lose or restore crystalline orders after quenching. Although MD simulations are limited by time scales, the simulation time 20 ps is long enough to predict the stable structure formed after prolonged cooling down to room temperature over a much longer time. Fig. 2 compares the U–U interatomic and O–O interatomic CSD values as a function time and under different ionization periods (5 fs, 8 fs, and 10 fs). First of all, CSD values expect to increase upon melting due to density changes. This explains why both U–U and O–O CSD values increase in the CE core region and reaches maximum at time 1.0 ps (marked as blue color) for all ionization conditions. Upon solidification to crystalline structures, CSD values change back to original values. This explains the disappearance of CSD step height like peaks at longer time t P 10 ps for both ionization period 5 fs and 8 fs (Fig. 2a and b). If the molten region is solidified into an amorphous structure, CSD values expect to drop, but to a value which is higher than original crystalline CSD. This is what observed for ionization period 10 fs (Fig. 2c), in which O–O CSD values at longest time t = 20 ps are still higher than original values. However, complexity rises if elemental segregation occurs in re-solidified amorphous region, which leads to reduced U–U CSD values at the core center at t = 20 ps (Fig. 2c). The low U–U CSD values at the core center are not caused by densification. Instead, for all three ionization cases, both U and O atomic densities are lower than the bulk values after CE, caused by significant atom diffusion from the core. Fig. 2b shows that at time 1 ps, the molten core reaches its maximum radius. The O–O CSD curve has a step height like feature with its value in the molten region being roughly a constant, which is characteristic of melting. Another evidence of melting comes from kinetic and potential energy change. Immediately after ionization pulse, exchange between kinetic energies and potential
2.8 2.7 0
20 40 (b) O-O, τ=8 fs
60
0 4.2
20 40 (e) U-U, τ=8 fs
60
0 fs 8 fs 4.0 1 ps 3 ps 10 ps 3.8 20 ps
3.0
2.8
0
20
40
60
3.2
0 fs 8 fs 1 ps 3 ps 10 ps 20 ps 0
4.2
(c) O-O, τ=10 fs
20 40 (f) U-U, τ=10 fs
60
2.8
0
20
40
60
80
0 fs 10 fs 1 ps 3 ps 10 ps 20 ps
0 fs 10 fs 4.0 1 ps 3 ps 10 ps 3.8 20 ps
3.0
80
0
20
40
60
80
Radial position (Angstrom) Fig. 2. Time evolution of radial CSD values of different interatomic systems and ionization periods: (a) 5 fs for O–O, (b) 8 fs for O–O, (c) 10 fs for O–O, (d) 5 fs for U–U, (e) 8 fs for U–U, and (f) 10 fs for U–U. The initial substrate temperature is 300 K. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
68
Z. Li et al. / Nuclear Instruments and Methods in Physics Research B 358 (2015) 65–71
Fig. 3. Kinetic energy distributions in UO2 originally at (a) 823 K and (b) 1223 K. The pulsed ionizations are introduced in the core of 10 Å radius and for a time period of 10 fs. The bottoms are enlarged images showing atom disorder levels.
energies occur. With increasing time, kinetic energies of O atoms decrease and their potential energies increase, and both reach saturations. For a small volume cut from the CE core, we found that the critical time to reach such saturation is about a fraction of 1 ps. Therefore, the time scale 1 ps is long enough to form molten zone. Using O–O CSD curves at 1 ps, melting radii are quantitatively determined by the radial position reaching half height of CSD peak. The radii are 16, 27, and 32 Å for ionization periods 5, 8 and 10 fs, respectively. Obviously, the longer the ionization periods, the larger the melting core radii. A critical energy density model can be used to explain this trend. Assuming that the amount of energy gain from CE is linearly proportion to the ionization period and further assuming that the energy density to cause melting is constant, pffiffiffiffiffiffiffiffiffiffiffiffi then the melting radius can be expressed as r ¼ a t t c , where t is the ionization period and tc is the minimum ionization period required to reach the threshold for melting. The first assumption is supported by comparing the kinetic energy differences as a function of ionization periods and the second assumption comes from the fact that melting is a phase change requiring a critical heat
Closest Separation Distance (Angstrom)
3.2
(a) T=300 K
energy. The best fitting is obtained with tc = 3.2 fs. Any ionization shorter than this will not create melting. As to be discussed, this simple expression does not consider dependence on substrate temperature. It is possible to develop a mode with reduced parameters to include both ionization period and substrate temperature effects. We further investigate the effect of substrate temperature on maximum melting radius. Fig. 3 compares MD simulations for substrate temperatures initially kept at 823 K and 1223 K. In both cases, the CE ionization periods are fixed to 10 fs. Atomic kinetic energy distributions show that dissipation of heat/energy into surrounding medium takes longer time at higher substrate temperature. As shown by the enlarged image from the box marked region, CE molten core ends with amorphous structures after quenching and solidification. Two dash lines in Fig. 3a mark the crystalline-to-amorphous boundaries, as judged from the O–O CSD peak positions. In a comparison, a higher temperature at 1223 K, as shown in Fig. 3b, creates a larger amorphous zone. Fig. 4 compares the time evolution of O–O CSD curves at different initial substrate temperatures. At the end of the ionization
(b) T=823 K
3.1 3.0 2.9 O-O 0.0 fs 10.0 fs 1.0 ps 3.0 ps 10.0 ps 20.0 ps
2.8 2.7 2.6 3.2
(c) T=1023 K
(d) T=1223 K
3.1 3.0 2.9 2.8 2.7 2.6 0
10
20
30
40
50 60 0 10 20 Position (Angstrom)
30
40
50
60
Fig. 4. Time evolution of radial O–O CSD values for ionization period of 10 fs, as a function of time in UO2 initially at temperatures of (a) 300 K, (b) 823 K, (c) 1023 K, and (d) 1223 K.
69
Z. Li et al. / Nuclear Instruments and Methods in Physics Research B 358 (2015) 65–71
function of temperatures is noticeable. Combining the data extracted from Figs. 2 and 4, a general formula can be obtained through best fitting, for prediction of melting radius as a function of substrate temperatures and ionization periods. Fig. 5 compares the modeling obtained data with a smooth surface plot from the equation, given by
pffiffiffiffiffiffiffiffiffiffiffiffi r ¼ ð11:27 þ 0:0031TÞ t tc
Fig. 5. Melting radii as a function of ionization periods and initial substrate temperatures. The shadowed plane is the fitting prediction and the dark cycles are data extracted from CSD curves at fixed time 1 ps.
(time = 10 fs), all CSD curves show a dip-like feature in the core center region. This is not caused by densification. Due to kinetic energy gained from the CE, energetic O atoms are able to move closer and have ultra-low energy collisions, which reduce CSD values. Maximum melting/amorphous radius is reached at 1–3 ps. At longest computation time 20 ps, radius is reduced, and its final value is higher for higher substrate temperatures. Using curves at 1 ps, the maximum melting radii are measured to be 32, 36, 38, and 39 Å for substrate temperatures of 300, 823, 1023, and 1223 K, respectively. Although the radius changes are small, their trend of changes as a
ð1Þ
where r is the melting radius (Å), T is temperature (K), t is ionization period (fs), and tc is minimum ionization period required to induce melting. The above fitting is for initial ionization radius of 10 Å only. The temperature dependence can be understood by the difference in heat dissipation rates. Melting radius is limited by two factors: the amount of energy deposition and the structural relaxation time. Both must be higher than certain critical values. First, the ionization period must be long enough (t > 3.2 fs) to deposit sufficient energy density to cause local melting. Second, even with sufficient energy deposition, melting requires a critical structural relaxation time so exchange of kinetic energy to potential energy and local energy equilibrium can be reached. If substrate temperature is low, the large temperature gradient between CE core and substrate will reduce the structural relaxation time, which most sensitively influences the melting boundary. At an elevated temperature, the heat dissipation is slow and the boundary region will stay longer at high temperature, which promotes melting and increases the melting radius. MD simulations further reveal the creation of shock waves. Fig. 6 plots the kinetic energy of U and O atoms as a function of time. The time t = 10 fs corresponds to the ending of the ionization pulse. At time t = 0.5 ps, the edge of the ionized cylindrical core begin to emit a kinetic energy pulse. This single pulsed pressure wave propagates with increasing time. One very interesting finding is that the shock wave forms only in U sublattice. As shown in Fig. 6, O sublattice does not create such a propagating wave.
2 eV
0 eV
Fig. 6. Cross sectional view of kinetic energies of (a) U atoms and (b) O atoms for ionization period of 10 fs and ionization temperature of 300 K.
70
Z. Li et al. / Nuclear Instruments and Methods in Physics Research B 358 (2015) 65–71
0 10 fs 0.5 ps 1.0 ps
(a) U Wave front
Kinetic energy (eV/atom)
1
0.1
0
20
40
60
(b) O
80 0 10 fs 0.5 ps 1.0 ps
1
0.1
0
20
40
60
80
Radial distance (Angstrom) Fig. 7. Time evolution of radial kinetic energy distributions of (a) U and (b) O atoms, for an ionization period of 10 fs and initial temperature of 300 K.
Fig. 7a and b plot the kinetic energy of U and O atoms, respectively, for the ionization period 10 fs. Immediately at the end of the pulse ionization, kinetic energy of the CE core rises up and peaks at the CE boundary. At time 0.5 ps, a pressure wave propagates to the distance 35 Å, with a peak height of 0.14 eV per atom. At time 1.0 ps, the shock wave peak moves to the distance 61 Å. The calculated wave speed is 5.2 103 m/s. This velocity is close to experimentally measured average longitudinal sound velocity of 5.4 103 m/s at ambient condition [24]. Estimation of pressure associated with local melting is given by DV/V = 3p/4l, where DV is the expansion of the initial volume V due to melting, and l is the shear modulus. Volume expansion upon transition from solid phase to liquid phase was experimentally measured to be 9.6% [25]. The l values range from 70.6 GPa to 93.3 GPa [26–28]. Taking an averaged l value, the estimated pressure is 10.5 GPa. One large difference between Figs. 7a and b is that the pulsed kinetic wave peaks seem missing in O sublattice. The O sublattice system exhibits very small kinetic energy peaks at the positions corresponding to U waves, but their peak heights are almost ignorable. We believe the difference is caused by mass difference. Since U atoms are much heavier than O atoms, the averaged vibration frequency of O atoms is much higher than that of U atoms under the same temperature. Upon compression, if the wave propagation time from one atomic row to its immediate neighbor is longer than the atoms’ vibration period, atoms will have oscillated movements back and forth many times during the time period when the wave front propagates to the next atom row. One consequence of this is that work done by the external compression force is canceled, i.e. atom movement is slower when the atom vibrates in the direction opposite to the wave propagation direction. Using a simple harmonic oscillator to represent an atom of mass m and vibration angular frequency x, its displacement y as a function of t, is described by 2
m
d y dt
2
¼ mx2 y þ F;
06t6s
ð2Þ
where F is the external force which is applied to the system for a time period of s. Solving the equation for y, and obtaining velocity and potential energy by using dy/dt and y, respectively, we obtain the following expression for U total energies
8 < U ¼ ðF sÞ2 when x 1 2m : U 6 2F 2 ðF sÞ2 when x 1 2m 2mx2
ð3Þ
Therefore, for light atoms such as O, their relatively higher vibration frequency x 1 leads to ignorable energy gain when compared with U atoms. The possibility of creating shock waves from a hot core has been discussed in early studies by Ronchi [13]. It was pointed out that surface tracks in neutron irradiated UO2 are caused by the reflection of shock waves from free surface. Under reflection, the atoms at the surface are accelerated to a velocity twice of the initial one [13], which lead to fracture or split off of the free surface, if tracks are shallowly located beneath a free surface. The MD simulations cannot handle electrons. Therefore, the effects of electron–phonon interactions are not included in the present study. Based on the estimation from considering electron free paths, Ronchi suggested that relaxation time for the heat transfer from electron subsystems to lattice atoms is about 10 ps for U atoms [13]. This occurs after maximum melting is reached already, according to Fig. 4. Therefore, the effect of electron–phonon interaction is to increase the temperatures and slower down the heat dissipation of the molten region. According to Eq. (1), this is equivalent to adding surrounding temperature at the liquid–solid interface and consequently increasing the melting radius. Analytical analysis suggested that electron–phonon interactions, purely from electron stopping power of fission fragment of typical energy 100 MeV, contributes to a substrate temperature which can rise up about 2000 K at a time scale about 10 ps. A rough estimation based on Fig. 5 suggests melting radius will increase by about 10 Å if this effect is considered. Although the present study simulates CE process with selected parameters and excludes the heating from electron subsystems, the obtained general formula can be easily modified to consider the realistic conditions through converted ionization levels, time periods and effective substrate temperatures. 4. Summary Pulsed ionization is introduced in UO2 to simulate Coulomb explosion caused by swift ions, in which local and temporal electron depletion under extreme high electronic stopping leads to strong Coulomb repulsive forces among lattice atoms and subsequent high energy deposition. We have found that both the time period of ionization pulses and initial substrate temperatures influence the melting radius of ion tracks. The melting creates amorphous cores after cooling. The radius dependence on ionization period can be explained by a critical energy density model in which the melting volume is proportional to deposited energy if the energy is above a critical value. With increasing substrate temperatures, the melting radius slightly increases. The effect is explained by a longer structural relaxation time at higher temperature, due to a reduced heat dissipation rate when temperature gradient between the melting zone and the surrounding medium is lower. Furthermore, shock waves are created upon melting-induced compression pulses. However, kinetic energy wave fronts are observed only in uranium sublattices. The absence or ignorable magnitude of waves in oxygen sublattices is explained by relatively higher vibration frequencies of oxygen atoms due to their lighter mass. When oxygen atoms’ average vibration speed is higher than the wave speed, the work done by the compression wave is canceled by their back and forth movements. References [1] E.M. Bringa, R.E. Johnson, Phys. Rev. Lett. 88 (2002) 165501. [2] I.M. Lifshitz, I.M. Kaganov, L.V. Tanatarov, J. Nucl. Energy A 12 (1960) 69. [3] M. Toulemonde, W. Assmann, C. Dufour, A. Meftah, F. Studer, C. Trautmann, Mat. Fys. Medd. Kong. Dan. Vid. Selsk. 52 (2006) 263. [4] A. Meftah, F. Brisard, J.M. Costantini, E. Dooryhee, M. Hage-Ali, M. Hervieu, J.P. Stoquert, F. Studer, M. Toulemonde, Phys. Rev. B 49 (1994) 12457.
Z. Li et al. / Nuclear Instruments and Methods in Physics Research B 358 (2015) 65–71 [5] Maik Lang, Jie Lian, Jiaming Zhang, Fuxiang. Zhang, William J. Weber, Christina. Trautmann, Rodney C. Ewing, Phys. Rev. B 79 (2009) 224105. [6] G. Bonfiglioli, A. Ferro, A. Mojoni, J. Appl. Phys. 32 (1961) 2499. [7] G. Szenes, Phys. Rev. B. 51 (1995) 8026. [8] M. Toulemonde, Ch. Dufour, A. Meftah, E. Paumier, Nucl. Instr. Meth. Phys. Res. B 166–167 (2000) 903. [9] R.L. Fleischer, P.B. Price, R.M. Walker, Nuclear Tracks in Solids: Principles and Applications, University of California, Berkeley, 1975. [10] R.L. Fleischer, P.B. Price, R.M. Walker, Rev. Sci. Instrum. 34 (1963) 510. [11] P. Apel, Nucl. Instr. Meth. Phys. Res. B 208 (2003) 11. [12] R.W. DeBlois, C.P. Bean, Rev. Sci. Instrum. 41 (1970) 909. [13] C. Ronchi, J. Appl. Phys. 44 (1973) 3575. [14] N. Ishikawa, T. Sonoda, T. Sawabe, H. Sugai, M. Sataka, Nucl. Instr. Meth. Phys. Res. B 314 (2013) 180. [15] Y. Sasajima, T. Osada, N. Ishikawa, A. Iwase, Nucl. Instr. Meth. Phys. Res. B 314 (2013) 195.
[16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26]
71
V.V. Pisarev, S.V. Starikov, J. Phys. Condens. Matter 26 (2014) 475401. Z. Li, D. Chen, J. Wang, L. Shao, Nucl. Instr. Meth. Phys. Rev. B 307 (2013) 86. Z. Li, D. Chen, L. Shao, J. Appl. Phys. 115 (2014) 143507. S.J. Plimpton, J. Comput. Phys. 117 (1995) 1. E. Yakub, C. Ronchi, D. Staicu, J. Chem. Phys. 127 (2007) 094508. K. Govers, S. Lemehov, M. Hou, M. Verwerft, J. Nucl. Mater. 366 (2007) 161. Y. Ida, Phys. Earth Planet. Inter. 13 (1976) 97. K. Kawamura, Springer Ser. Solid-State Sci. 103 (1992) 88. I.J. Fritz, J. Appl. Phys. 47 (1976) 4353. J.A. Christensen, J. Am. Ceram. Soc. 46 (1963) 607. S.M. Lang, Properties of High-Temperature Ceramics and Cermets, Elasticity and Density at Room Temperature, Monograph 6, National Bureau of Standards, Washington, D.C., 1960. [27] J.B. Wachtman Jr., M.L. Wheat, J. Nucl. Mater. 16 (1965) 39. [28] A. Padel, Ch. de Novion, J. Nucl. Mater. 33 (1969) 40.
Contents lists available at ScienceDirect
Nuclear Instruments and Methods in Physics Research B journal homepage: www.elsevier.com/locate/nimb
Melting and shock wave creation in uranium oxide due to Coulomb explosion after a pulsed ionization Zhongyu Li a, Di Chen b, Lin Shao b,⇑ a b
Fundamental Science on Nuclear Safety and Simulation Technology Laboratory, Harbin Engineering University, Harbin 150001, China Department of Nuclear Engineering, Texas A&M University, College Station, TX 77843, USA
a r t i c l e
i n f o
Article history: Received 23 February 2015 Received in revised form 26 April 2015 Accepted 28 April 2015
Keywords: Molecular dynamics simulations Shock wave Ion irradiation Uranium oxide
a b s t r a c t By means of molecular dynamics simulations, we study the effects of pulsed ionization in uranium oxide (UO2), which occurs when UO2 is bombarded with swift ions or fission fragments. A general formula is developed to predict melting radius under various conditions due to electron stripping and Coulomb explosion (CE). A critical density model is suggested in which the melting volume is proportional to ionization period, if the period is above a critical value. The maximum melting radius depends on the time period of structural relaxation above the melting temperature, which increases with increasing initial substrate temperatures due to a lower heat dissipation rate. Furthermore, shock waves are observed to emit from CE core but the kinetic energy wave peak exists only in U sublattices. The absence of kinetic energy waves in O sublattices is explained by their relatively higher thermal vibration which cancels the work done from the compression waves. Ó 2015 Elsevier B.V. All rights reserved.
1. Introduction In numerous previous studies, it has been observed that swift ions can create ion tracks in insulators, semiconductors, and polymers. The fundamental mechanisms of track formation, however, are still subject to some debate [1–5]. Coulomb explosion (CE) and thermal spike model were often compared as competing mechanisms. Bringa and Johnson supported CE mechanism [1], while many other studies favored thermal spike model [2–5]. In the thermal spike model, energetic electrons are created under a high electron stopping power and subsequent electron-lattice coupling leads to melting and lattice disorders [6–8]. In the CE model, swift ions cause strong electron ionization of target atoms, leading to partial electron depletion along an ion track. In semiconducting or insulating substrates, the lack of free electrons makes it difficult to quickly recover the charge neutrality through electron back diffusion. With an ionization period sufficiently long (>a few fs), Coulomb force among the positively charged target atoms lead to CE [9], through the following stages: first, Coulomb potentials quickly convert to kinetic energies of atoms during the ionization pulse. Second, within a factor of 1 ps, kinetic energies are shared among atoms and local thermal equilibrium is reached. This process also leads to local melting. Third, the melting zones are cooled ⇑ Corresponding author. E-mail address: [email protected] (L. Shao). http://dx.doi.org/10.1016/j.nimb.2015.04.077 0168-583X/Ó 2015 Elsevier B.V. All rights reserved.
down and form either amorphous or polycrystalline cores, depending on rates of heat dissipation to the surrounding medium. Swift ion induced ion track formations can be used in microtechnology and nanotechnology for a wide range of applications [10–12]. Ion track formation in uranium oxide (UO2) is particular important for light water reactors, due to the facts that fission fragments from neutron reactions, in a typical kinetic energy of 100 MeV, may lead to structural changes which cause thermal and mechanical property degradation, and influence reactor performance. Experiments on both accelerator based swift ion irradiations and reactor irradiation have observed ion track formation in UO2 [13,14]. Computer simulations have been used to understand CE effects [15–18]. In the present study, we preformed systematic molecular dynamics (MD) simulations, in order to understand CE mechanism. Many previous studies have used MD simulations to understand thermal spike model. Modeling studies on CE phenomenon, however, are very limited. The present study has no intention to evaluate the competing effects between Coulomb explosion model and thermal spike model, but expect to benefit the development of a comprehensive model to possibly combine both effects. 2. Modeling procedure The present study is aimed to understand atomic scale details of CE in UO2, and to identity key parameters which influence the
66
Z. Li et al. / Nuclear Instruments and Methods in Physics Research B 358 (2015) 65–71
microstructure changes after CE. Instead of directly linking the modeling to specific fission production irradiation, CE is introduced by a pulsed ionization. The ionization time period and the substrate temperatures are used as variables for systematic investigation. Another reason to introduce a pulsed ionization is to alleviate issues from molecular dynamics (MD) simulations. MD simulations cannot handle electron subsystems. MD simulations are performed by using the LAMMPS (Large-scale Atomic Molecular Massively Parallel Simulator) code [19]. The interatomic potentials were described by a partly ionic Busing-Ida-type using potential parameters proposed by Yakub et al. [20]. Many interatomic potentials have been developed to describe UO2, as reviewed by Govers et al. [21]. These potentials include the two major types [20]: one has charges fixed by the chemical compositions and the other, the Busing-Ida type [22,23], has fractional effective charges. Introducing effective charges provides an important fitting parameter to match with experimental observation. Most importantly, the potential used in the present study is able to accurately predict melting and existence of a pre-melting structural transition [20]. As to be discussed, the transition from solid to liquid is one important phenomenon introduced by CE. The used potentials include Coulomb interactions and short-range interactions. The short-range interactions consist of O–O van der Waals attractions, and U–O covalent bonding represented by the Morse function and overlap (exchange) repulsion. The unit cell has a dimension of about 22 nm 22 nm 5.5 nm and contains 192,000 atoms. The cell is [1 0 0] oriented. All sides take periodic boundary conditions, including the top and the bottom of the cell. The geometry of the cell is schematically shown in Fig. 1. After thermal relaxation at 300 K for 50 ps, an ionization pulse was introduced in a cylindrical region of radius 20 Å. From 0 to 10 Å, half of the U and O atoms within the ionized region are
stripped by one electron per atom, f.g. U4+ ? U5+, to simulate CE. From 10 to 20 Å, the stripped electrons in the CE region are re-distributed, forming an electron rich cylindrical shell surrounding the electron-stripped core. These electrons are added to local U/O atoms so their Columbic force contributions are included in MD simulations. Without such treatment, charge neutrality cannot be preserved, and the accumulated electric fields from positively charged core are high enough to introduce significant Coulomb force and create abnormal cell vibrations. The approach of redistribution of stripped electrons near the CE core effectively creates charge screening effect to eliminate the long range Coulomb force from the core. At the end of the ionization, these redistributed negative charges, together with positive charges introduced in the CE core, are removed simultaneously under assumption of quick charge recovery. For the present study, ionization period range from 5 fs to 10 fs. The time selections consider the facts that electron ionization, and electron–electron thermalization occur in a typical time scale of 1 fs [2]. Once the ionization pulse ends, the pair Coulomb potential was removed. For the rest of the discussion, the time t = 0 refers to the moment when the ionization starts. The cell is divided into cylindrical slabs and closest atomic separated distance (CSD) of each slab, based on the average value of the first neighbor distance within the slab, is extracted to plot CSD radial distributions. The melting radius can be quantitatively determined since CSD values of melting UO2 are larger than that of amorphous UO2. And both of them are larger than that of crystalline UO2. It needs to point out that a defective crystalline UO2 layer may create high CSD values but they are typically less than that of amorphous UO2 and do not show saturation. For a local molten core surrounded by defective zones, radial CSD values are featured with a step height peak corresponding to the molten core and a decreasing tail corresponding to surrounding defective zone. Therefore, we use half height of
Fig. 1. Kinetic energy distributions in UO2 for an ionization period of (a) 5 fs and (b) 10 fs. Ionization is introduced in the cylindrical region of 10 Å radius. The initial temperature is 300 K. The bottoms are enlarged images showing atom disorder levels. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
67
Z. Li et al. / Nuclear Instruments and Methods in Physics Research B 358 (2015) 65–71
CSD maximum to define the melting radius. This makes it possible to track liquid-to-solid changes as a function of time. The time step for structural relaxation prior to the ionization is 1 fs. The time step during the ionization is 0.1 fs, and that for structural evolution after the ionization is 1 fs. 3. Results and discussion Fig. 1a and b compare snapshots of MD simulations as a function of time, for ionization time period 5 fs and 10 fs, respectively. Color bar represents the kinetic energy of atoms. At time t = 0 s, half of randomly selected atoms within a cylinder of radius r = 10 Å are charged by stripping one electron per atom. Their charge statuses remain until the end of the ionization periods. As shown by the kinetic energy distributions, at the end of the ionization, the CE process does not introduce significant kinetic energies within the core region, due to Coulomb force cancelation from all directions. Instead, the atoms at the track edge gain the highest kinetic energy due to biased force along the radial distributions. At 1 ps, localized energy equilibration is reached and the whole core regions have relatively uniform kinetic energy distributions. This stage is also featured with local melting, with atoms in the core losing short range orders. At longer time t = 20 ps, quenching of the molten core is finished. Kinetic energies are distributed over much wider region. As shown by the enlarged image in the box marked region, there is significant difference in crystalline orders between the short ionization (5 fs) and long ionization (10 fs). For the short ionization, thermal quenching ends with resolidification of the molten core into the defective but crystalline structure. In comparison, the long ionization ends with fully amorphous core. In the enlarged images, orange color refers to U atoms, white color refers to charged U with one electron stripped during the ionization pulse. Green color refers to O atoms and black color refers to electron-stripped O atom during the ionization pulse. Atomic distributions show that both electron-stripped U and O atoms have
3.0
4.0
(a) O-O, τ=5 fs
(d) U-U, τ=5 fs 0 fs 5 fs 1 ps 3 ps 10 ps 20 ps
0 fs 5 fs 1 ps 3.9 3 ps 10 ps 20 ps
2.9
Closest seperation distance (Angstrom)
diffusion from the molten core but O atoms diffuse more due to their relatively small mass. These observations suggest that the ionization period plays an important role to determine whether CE core region will lose or restore crystalline orders after quenching. Although MD simulations are limited by time scales, the simulation time 20 ps is long enough to predict the stable structure formed after prolonged cooling down to room temperature over a much longer time. Fig. 2 compares the U–U interatomic and O–O interatomic CSD values as a function time and under different ionization periods (5 fs, 8 fs, and 10 fs). First of all, CSD values expect to increase upon melting due to density changes. This explains why both U–U and O–O CSD values increase in the CE core region and reaches maximum at time 1.0 ps (marked as blue color) for all ionization conditions. Upon solidification to crystalline structures, CSD values change back to original values. This explains the disappearance of CSD step height like peaks at longer time t P 10 ps for both ionization period 5 fs and 8 fs (Fig. 2a and b). If the molten region is solidified into an amorphous structure, CSD values expect to drop, but to a value which is higher than original crystalline CSD. This is what observed for ionization period 10 fs (Fig. 2c), in which O–O CSD values at longest time t = 20 ps are still higher than original values. However, complexity rises if elemental segregation occurs in re-solidified amorphous region, which leads to reduced U–U CSD values at the core center at t = 20 ps (Fig. 2c). The low U–U CSD values at the core center are not caused by densification. Instead, for all three ionization cases, both U and O atomic densities are lower than the bulk values after CE, caused by significant atom diffusion from the core. Fig. 2b shows that at time 1 ps, the molten core reaches its maximum radius. The O–O CSD curve has a step height like feature with its value in the molten region being roughly a constant, which is characteristic of melting. Another evidence of melting comes from kinetic and potential energy change. Immediately after ionization pulse, exchange between kinetic energies and potential
2.8 2.7 0
20 40 (b) O-O, τ=8 fs
60
0 4.2
20 40 (e) U-U, τ=8 fs
60
0 fs 8 fs 4.0 1 ps 3 ps 10 ps 3.8 20 ps
3.0
2.8
0
20
40
60
3.2
0 fs 8 fs 1 ps 3 ps 10 ps 20 ps 0
4.2
(c) O-O, τ=10 fs
20 40 (f) U-U, τ=10 fs
60
2.8
0
20
40
60
80
0 fs 10 fs 1 ps 3 ps 10 ps 20 ps
0 fs 10 fs 4.0 1 ps 3 ps 10 ps 3.8 20 ps
3.0
80
0
20
40
60
80
Radial position (Angstrom) Fig. 2. Time evolution of radial CSD values of different interatomic systems and ionization periods: (a) 5 fs for O–O, (b) 8 fs for O–O, (c) 10 fs for O–O, (d) 5 fs for U–U, (e) 8 fs for U–U, and (f) 10 fs for U–U. The initial substrate temperature is 300 K. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
68
Z. Li et al. / Nuclear Instruments and Methods in Physics Research B 358 (2015) 65–71
Fig. 3. Kinetic energy distributions in UO2 originally at (a) 823 K and (b) 1223 K. The pulsed ionizations are introduced in the core of 10 Å radius and for a time period of 10 fs. The bottoms are enlarged images showing atom disorder levels.
energies occur. With increasing time, kinetic energies of O atoms decrease and their potential energies increase, and both reach saturations. For a small volume cut from the CE core, we found that the critical time to reach such saturation is about a fraction of 1 ps. Therefore, the time scale 1 ps is long enough to form molten zone. Using O–O CSD curves at 1 ps, melting radii are quantitatively determined by the radial position reaching half height of CSD peak. The radii are 16, 27, and 32 Å for ionization periods 5, 8 and 10 fs, respectively. Obviously, the longer the ionization periods, the larger the melting core radii. A critical energy density model can be used to explain this trend. Assuming that the amount of energy gain from CE is linearly proportion to the ionization period and further assuming that the energy density to cause melting is constant, pffiffiffiffiffiffiffiffiffiffiffiffi then the melting radius can be expressed as r ¼ a t t c , where t is the ionization period and tc is the minimum ionization period required to reach the threshold for melting. The first assumption is supported by comparing the kinetic energy differences as a function of ionization periods and the second assumption comes from the fact that melting is a phase change requiring a critical heat
Closest Separation Distance (Angstrom)
3.2
(a) T=300 K
energy. The best fitting is obtained with tc = 3.2 fs. Any ionization shorter than this will not create melting. As to be discussed, this simple expression does not consider dependence on substrate temperature. It is possible to develop a mode with reduced parameters to include both ionization period and substrate temperature effects. We further investigate the effect of substrate temperature on maximum melting radius. Fig. 3 compares MD simulations for substrate temperatures initially kept at 823 K and 1223 K. In both cases, the CE ionization periods are fixed to 10 fs. Atomic kinetic energy distributions show that dissipation of heat/energy into surrounding medium takes longer time at higher substrate temperature. As shown by the enlarged image from the box marked region, CE molten core ends with amorphous structures after quenching and solidification. Two dash lines in Fig. 3a mark the crystalline-to-amorphous boundaries, as judged from the O–O CSD peak positions. In a comparison, a higher temperature at 1223 K, as shown in Fig. 3b, creates a larger amorphous zone. Fig. 4 compares the time evolution of O–O CSD curves at different initial substrate temperatures. At the end of the ionization
(b) T=823 K
3.1 3.0 2.9 O-O 0.0 fs 10.0 fs 1.0 ps 3.0 ps 10.0 ps 20.0 ps
2.8 2.7 2.6 3.2
(c) T=1023 K
(d) T=1223 K
3.1 3.0 2.9 2.8 2.7 2.6 0
10
20
30
40
50 60 0 10 20 Position (Angstrom)
30
40
50
60
Fig. 4. Time evolution of radial O–O CSD values for ionization period of 10 fs, as a function of time in UO2 initially at temperatures of (a) 300 K, (b) 823 K, (c) 1023 K, and (d) 1223 K.
69
Z. Li et al. / Nuclear Instruments and Methods in Physics Research B 358 (2015) 65–71
function of temperatures is noticeable. Combining the data extracted from Figs. 2 and 4, a general formula can be obtained through best fitting, for prediction of melting radius as a function of substrate temperatures and ionization periods. Fig. 5 compares the modeling obtained data with a smooth surface plot from the equation, given by
pffiffiffiffiffiffiffiffiffiffiffiffi r ¼ ð11:27 þ 0:0031TÞ t tc
Fig. 5. Melting radii as a function of ionization periods and initial substrate temperatures. The shadowed plane is the fitting prediction and the dark cycles are data extracted from CSD curves at fixed time 1 ps.
(time = 10 fs), all CSD curves show a dip-like feature in the core center region. This is not caused by densification. Due to kinetic energy gained from the CE, energetic O atoms are able to move closer and have ultra-low energy collisions, which reduce CSD values. Maximum melting/amorphous radius is reached at 1–3 ps. At longest computation time 20 ps, radius is reduced, and its final value is higher for higher substrate temperatures. Using curves at 1 ps, the maximum melting radii are measured to be 32, 36, 38, and 39 Å for substrate temperatures of 300, 823, 1023, and 1223 K, respectively. Although the radius changes are small, their trend of changes as a
ð1Þ
where r is the melting radius (Å), T is temperature (K), t is ionization period (fs), and tc is minimum ionization period required to induce melting. The above fitting is for initial ionization radius of 10 Å only. The temperature dependence can be understood by the difference in heat dissipation rates. Melting radius is limited by two factors: the amount of energy deposition and the structural relaxation time. Both must be higher than certain critical values. First, the ionization period must be long enough (t > 3.2 fs) to deposit sufficient energy density to cause local melting. Second, even with sufficient energy deposition, melting requires a critical structural relaxation time so exchange of kinetic energy to potential energy and local energy equilibrium can be reached. If substrate temperature is low, the large temperature gradient between CE core and substrate will reduce the structural relaxation time, which most sensitively influences the melting boundary. At an elevated temperature, the heat dissipation is slow and the boundary region will stay longer at high temperature, which promotes melting and increases the melting radius. MD simulations further reveal the creation of shock waves. Fig. 6 plots the kinetic energy of U and O atoms as a function of time. The time t = 10 fs corresponds to the ending of the ionization pulse. At time t = 0.5 ps, the edge of the ionized cylindrical core begin to emit a kinetic energy pulse. This single pulsed pressure wave propagates with increasing time. One very interesting finding is that the shock wave forms only in U sublattice. As shown in Fig. 6, O sublattice does not create such a propagating wave.
2 eV
0 eV
Fig. 6. Cross sectional view of kinetic energies of (a) U atoms and (b) O atoms for ionization period of 10 fs and ionization temperature of 300 K.
70
Z. Li et al. / Nuclear Instruments and Methods in Physics Research B 358 (2015) 65–71
0 10 fs 0.5 ps 1.0 ps
(a) U Wave front
Kinetic energy (eV/atom)
1
0.1
0
20
40
60
(b) O
80 0 10 fs 0.5 ps 1.0 ps
1
0.1
0
20
40
60
80
Radial distance (Angstrom) Fig. 7. Time evolution of radial kinetic energy distributions of (a) U and (b) O atoms, for an ionization period of 10 fs and initial temperature of 300 K.
Fig. 7a and b plot the kinetic energy of U and O atoms, respectively, for the ionization period 10 fs. Immediately at the end of the pulse ionization, kinetic energy of the CE core rises up and peaks at the CE boundary. At time 0.5 ps, a pressure wave propagates to the distance 35 Å, with a peak height of 0.14 eV per atom. At time 1.0 ps, the shock wave peak moves to the distance 61 Å. The calculated wave speed is 5.2 103 m/s. This velocity is close to experimentally measured average longitudinal sound velocity of 5.4 103 m/s at ambient condition [24]. Estimation of pressure associated with local melting is given by DV/V = 3p/4l, where DV is the expansion of the initial volume V due to melting, and l is the shear modulus. Volume expansion upon transition from solid phase to liquid phase was experimentally measured to be 9.6% [25]. The l values range from 70.6 GPa to 93.3 GPa [26–28]. Taking an averaged l value, the estimated pressure is 10.5 GPa. One large difference between Figs. 7a and b is that the pulsed kinetic wave peaks seem missing in O sublattice. The O sublattice system exhibits very small kinetic energy peaks at the positions corresponding to U waves, but their peak heights are almost ignorable. We believe the difference is caused by mass difference. Since U atoms are much heavier than O atoms, the averaged vibration frequency of O atoms is much higher than that of U atoms under the same temperature. Upon compression, if the wave propagation time from one atomic row to its immediate neighbor is longer than the atoms’ vibration period, atoms will have oscillated movements back and forth many times during the time period when the wave front propagates to the next atom row. One consequence of this is that work done by the external compression force is canceled, i.e. atom movement is slower when the atom vibrates in the direction opposite to the wave propagation direction. Using a simple harmonic oscillator to represent an atom of mass m and vibration angular frequency x, its displacement y as a function of t, is described by 2
m
d y dt
2
¼ mx2 y þ F;
06t6s
ð2Þ
where F is the external force which is applied to the system for a time period of s. Solving the equation for y, and obtaining velocity and potential energy by using dy/dt and y, respectively, we obtain the following expression for U total energies
8 < U ¼ ðF sÞ2 when x 1 2m : U 6 2F 2 ðF sÞ2 when x 1 2m 2mx2
ð3Þ
Therefore, for light atoms such as O, their relatively higher vibration frequency x 1 leads to ignorable energy gain when compared with U atoms. The possibility of creating shock waves from a hot core has been discussed in early studies by Ronchi [13]. It was pointed out that surface tracks in neutron irradiated UO2 are caused by the reflection of shock waves from free surface. Under reflection, the atoms at the surface are accelerated to a velocity twice of the initial one [13], which lead to fracture or split off of the free surface, if tracks are shallowly located beneath a free surface. The MD simulations cannot handle electrons. Therefore, the effects of electron–phonon interactions are not included in the present study. Based on the estimation from considering electron free paths, Ronchi suggested that relaxation time for the heat transfer from electron subsystems to lattice atoms is about 10 ps for U atoms [13]. This occurs after maximum melting is reached already, according to Fig. 4. Therefore, the effect of electron–phonon interaction is to increase the temperatures and slower down the heat dissipation of the molten region. According to Eq. (1), this is equivalent to adding surrounding temperature at the liquid–solid interface and consequently increasing the melting radius. Analytical analysis suggested that electron–phonon interactions, purely from electron stopping power of fission fragment of typical energy 100 MeV, contributes to a substrate temperature which can rise up about 2000 K at a time scale about 10 ps. A rough estimation based on Fig. 5 suggests melting radius will increase by about 10 Å if this effect is considered. Although the present study simulates CE process with selected parameters and excludes the heating from electron subsystems, the obtained general formula can be easily modified to consider the realistic conditions through converted ionization levels, time periods and effective substrate temperatures. 4. Summary Pulsed ionization is introduced in UO2 to simulate Coulomb explosion caused by swift ions, in which local and temporal electron depletion under extreme high electronic stopping leads to strong Coulomb repulsive forces among lattice atoms and subsequent high energy deposition. We have found that both the time period of ionization pulses and initial substrate temperatures influence the melting radius of ion tracks. The melting creates amorphous cores after cooling. The radius dependence on ionization period can be explained by a critical energy density model in which the melting volume is proportional to deposited energy if the energy is above a critical value. With increasing substrate temperatures, the melting radius slightly increases. The effect is explained by a longer structural relaxation time at higher temperature, due to a reduced heat dissipation rate when temperature gradient between the melting zone and the surrounding medium is lower. Furthermore, shock waves are created upon melting-induced compression pulses. However, kinetic energy wave fronts are observed only in uranium sublattices. The absence or ignorable magnitude of waves in oxygen sublattices is explained by relatively higher vibration frequencies of oxygen atoms due to their lighter mass. When oxygen atoms’ average vibration speed is higher than the wave speed, the work done by the compression wave is canceled by their back and forth movements. References [1] E.M. Bringa, R.E. Johnson, Phys. Rev. Lett. 88 (2002) 165501. [2] I.M. Lifshitz, I.M. Kaganov, L.V. Tanatarov, J. Nucl. Energy A 12 (1960) 69. [3] M. Toulemonde, W. Assmann, C. Dufour, A. Meftah, F. Studer, C. Trautmann, Mat. Fys. Medd. Kong. Dan. Vid. Selsk. 52 (2006) 263. [4] A. Meftah, F. Brisard, J.M. Costantini, E. Dooryhee, M. Hage-Ali, M. Hervieu, J.P. Stoquert, F. Studer, M. Toulemonde, Phys. Rev. B 49 (1994) 12457.
Z. Li et al. / Nuclear Instruments and Methods in Physics Research B 358 (2015) 65–71 [5] Maik Lang, Jie Lian, Jiaming Zhang, Fuxiang. Zhang, William J. Weber, Christina. Trautmann, Rodney C. Ewing, Phys. Rev. B 79 (2009) 224105. [6] G. Bonfiglioli, A. Ferro, A. Mojoni, J. Appl. Phys. 32 (1961) 2499. [7] G. Szenes, Phys. Rev. B. 51 (1995) 8026. [8] M. Toulemonde, Ch. Dufour, A. Meftah, E. Paumier, Nucl. Instr. Meth. Phys. Res. B 166–167 (2000) 903. [9] R.L. Fleischer, P.B. Price, R.M. Walker, Nuclear Tracks in Solids: Principles and Applications, University of California, Berkeley, 1975. [10] R.L. Fleischer, P.B. Price, R.M. Walker, Rev. Sci. Instrum. 34 (1963) 510. [11] P. Apel, Nucl. Instr. Meth. Phys. Res. B 208 (2003) 11. [12] R.W. DeBlois, C.P. Bean, Rev. Sci. Instrum. 41 (1970) 909. [13] C. Ronchi, J. Appl. Phys. 44 (1973) 3575. [14] N. Ishikawa, T. Sonoda, T. Sawabe, H. Sugai, M. Sataka, Nucl. Instr. Meth. Phys. Res. B 314 (2013) 180. [15] Y. Sasajima, T. Osada, N. Ishikawa, A. Iwase, Nucl. Instr. Meth. Phys. Res. B 314 (2013) 195.
[16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26]
71
V.V. Pisarev, S.V. Starikov, J. Phys. Condens. Matter 26 (2014) 475401. Z. Li, D. Chen, J. Wang, L. Shao, Nucl. Instr. Meth. Phys. Rev. B 307 (2013) 86. Z. Li, D. Chen, L. Shao, J. Appl. Phys. 115 (2014) 143507. S.J. Plimpton, J. Comput. Phys. 117 (1995) 1. E. Yakub, C. Ronchi, D. Staicu, J. Chem. Phys. 127 (2007) 094508. K. Govers, S. Lemehov, M. Hou, M. Verwerft, J. Nucl. Mater. 366 (2007) 161. Y. Ida, Phys. Earth Planet. Inter. 13 (1976) 97. K. Kawamura, Springer Ser. Solid-State Sci. 103 (1992) 88. I.J. Fritz, J. Appl. Phys. 47 (1976) 4353. J.A. Christensen, J. Am. Ceram. Soc. 46 (1963) 607. S.M. Lang, Properties of High-Temperature Ceramics and Cermets, Elasticity and Density at Room Temperature, Monograph 6, National Bureau of Standards, Washington, D.C., 1960. [27] J.B. Wachtman Jr., M.L. Wheat, J. Nucl. Mater. 16 (1965) 39. [28] A. Padel, Ch. de Novion, J. Nucl. Mater. 33 (1969) 40.