Nuclear Inst. and Methods in Physics Research B 455 (2019) 1–6
Contents lists available at ScienceDirect
Nuclear Inst. and Methods in Physics Research B journal homepage: www.elsevier.com/locate/nimb
Migration behaviors of helium atoms near tungsten surfaces: A molecular dynamics study Jiechao Cui, Min Li, Baoqin Fu, Qing Hou
T
⁎
Key Lab for Radiation Physics and Technology, Ministry of Education, Institute of Nuclear Science and Technology, Sichuan University, Chengdu 610064, China
ARTICLE INFO
ABSTRACT
Keywords: Molecular dynamics simulations Single helium atom Migration behaviors Orientation surface
Under the fusion environment, the behavior of helium atoms near tungsten (W) surfaces plays a crucial role in surface morphological evolution and near-surface structural evolution of tungsten. In this paper, the migration behaviors of single helium atoms near W (2 1 1) and W (3 1 1) surfaces were investigated using molecular dynamics simulations. The results show that these single atoms can be well-described by the theory of continuous diffusion of particles in a semi-infinite medium, and different types of trap mutations occur for both surfaces. Although the temperature has an impotent impact on the occurrence probabilities of different types of trap mutations and the probabilities of helium atoms escaping from trap sites, there is no a strict relationship between the depth and temperature for trap mutations occurring. For the W (2 1 1) surface, an approximately one-dimensional diffusion may take place along the 〈1 1 1〉 directions when the helium atoms are trapped near 2 layers below the W surface. For the W (3 1 1) surface, the behaviors of trap mutations seem more significant than (2 1 1) and (1 1 1) surfaces at T = 1000 K. To investigate the diffusion and trap mutation processes, the nudged elastic band (NEB) method was also applied. These results can be helpful for understanding of helium release, helium retention, subsurface helium clustering, as well as surface morphological evolution.
1. Introduction The behavior of helium (He) near tungsten (W) surfaces plays a crucial role in surface morphological evolution and near-surface structural evolution of W. For example in ITER, W is chosen as a plasma facing material (PFM) for the divertor and it will suffer from low energy (< 200 eV) and high flux (> 1022–1024 m−2 s−1) irradiation of He plasma. Although the incident energies of He are lower than the displacement energy of W atoms, high-density defects, such as blistering, swelling, embrittlement as well as complex nanostructure [1–7], can still form. This is because that the strong attraction between He–He can cause He accumulation and the formation of small He clusters even in perfect W crystals. He cluster growth will create interstitial W atoms and He bubbles (He-vacancy complexes). Then, the following macroscopic defects mentioned above could be generated. In addition, these defects behaviors depend on He exposure conditions, and can lead to a drastic decline of the physical and mechanical properties which could affect the safe operation of the reactor. Therefore, it has become an exigent task to investigate the subsequent dynamic processes of He near W surfaces to comprehensively understand the complex multiscale behavior and to establish the corresponding models under the fusion reaction condition. As such, the above processes and various other aspects ⁎
of properties, including implantation of He atoms, diffusion and accumulation of He, bubble formation and coalescence, as well as He release, have been the subjects of a number of computational studies. In this article, we aim to clarify migration behaviors of single He atoms near W surfaces using molecular dynamics (MD) method. In fact, many researchers through computational investigations [8–10] have indicated differential He behaviors near W surfaces, including of course migration of He which is affected by surface orientations. By using the MD method, Li [11] et al. reported the retention and distribution of He atoms when they investigated cumulative He bombardments on W surfaces. They further found that the formation and growth of He clusters and the surface evolution of tungsten substrates are influenced by surface orientations and temperatures. Hammond and Wirth [12] also found a pronounced effect of surface orientation on the initial depth of implanted He ions, as well as a difference in reflection and He retention across different surface orientations. In particular, near (1 1 1) and (2 1 1) surfaces, even single He atoms are sufficient to induce the formation of tungsten adatom and substitutional He pairs, while (0 0 1) and (0 1 1) surfaces need two or three He atoms. Even though only one temperature (933 K) was considered in their work, this is a very appealing result, as such He clusters only containing more than 5 atoms could cause the formation of a vacancy and interstitial W atom in the
Corresponding author. E-mail address:
[email protected] (Q. Hou).
https://doi.org/10.1016/j.nimb.2019.06.015 Received 7 August 2018; Received in revised form 1 May 2019; Accepted 12 June 2019 Available online 17 June 2019 0168-583X/ © 2019 Elsevier B.V. All rights reserved.
Nuclear Inst. and Methods in Physics Research B 455 (2019) 1–6
J. Cui, et al.
(311) surface
(211) surface
400K
1000K
2000K
[111]
[233] [211]
[311]
Fig. 1. Merged snapshots of 1000 He atoms near W (2 1 1) and (3 1 1) surfaces at the end of the simulation time for the temperature of T = 400 K, 1000 K and 2000 K, respectively. Colored spheres indicate He atoms with corresponding CSP values and dark blue dots indicate the initial W lattice sites.
bulk [13]. A similar phenomenon, usually called trap mutation, was reported by Perez et al. [14] and Hammond et al. [15] using MD method, Pan et al. [16] based on density functional theory (DFT) calculations, and Barashev et al. [17] using self-evolving atomistic kinetic Monte Carlo (SEAKMC) method. To bridge the MD simulations with other simulation methods, such as kinetic Monte Carlo (KMC) and rate theory (RT), our group also carried out relevant investigations of single He atoms using MD method and found different types of trap mutations occur near (1 1 1) surfaces, while for (0 0 1) and (0 1 1) surfaces, He atoms just quickly escape out. But for the all three surfaces, the behaviors of single He atoms near surfaces can be well-described by the theory of continuous diffusion of particles in a semi-infinite medium [18]. All results above [12,14–18] show a single He may even induce vacancy formation near some W surfaces. In the real environment, these
formed complexes will trap more interstitial He atoms and cause He accumulation, and must affect He release, He retention, subsurface He clustering, as well as surface morphological evolution. Thus, it is necessary to thoroughly research the dependence of the relevant mechanism on the surface orientations. Here, based on MD method, we continued the previous study [18] and focused on the migration behaviors of single He atoms near W (2 1 1) surfaces, where trap mutations can take place according to many reports, and W (3 1 1) surfaces, which were rarely investigated. At different temperatures, the states of the He atoms in the substrates were identified by their centrosymmetry parameters. The rate dependence for the He atoms to escape out of the substrates on the surface orientation and temperature was analyzed and quantified based on the theory for continuous diffusion and thermal desorption. Activation energies of the various He migration pathways were determined by a climbing-image nudged elastic band (NEB) 2
Nuclear Inst. and Methods in Physics Research B 455 (2019) 1–6
J. Cui, et al.
calculation [19,20].
(0 1 1) and (1 1 1) surfaces, the theory of continuous diffusion of particles in a semi-infinite medium is applied here in the same way. As the previous paper has discussed this theory in detail [18], we just explain it briefly here. The He atoms can be considered to be migrating with a diffusion coefficient, D, in the bulk tungsten until arriving at the absorbing layer (tungsten surface, for example) and are absorbed instantly with the probability denoted by R. Then, the number of absorbed He atoms, Na, at time t is
2. Computational methods Using multiple graphics processing units (GPUs) for parallel computing, a MD package MDPSCU (Molecular Dynamics Package of Sichuan University) was applied in all simulations [21]. In this paper, as only a single He atom in one simulation box, two kinds of potentials were applied: a Finnis–Sinclair-type potential obtained by Ackland and Thetford for describing the W-W interaction [22] and a pair-wise potential proposed by Wang et al. for the W-He interaction [23]. Before the MD simulation run, BCC tungsten boxes of size 8.49 α0 × 8.66 α0 × 29.39 α0 for (2 1 1) surface, and size 8.49 α0 × 9.38 α0 × 29.85 α0 for (3 1 1) surface were prepared separately, where α0 = 3.1652 Å is the lattice constant of tungsten. For each box size, 1000 independent replicas were prepared for the purpose of statistical analysis. The MD simulation run comprised two stages. In the first stage, periodic boundary conditions were employed in all three directions and one He atom was inserted at random position in each box. At first, 1000 simulation boxes were quenched to zero temperature. Then, the simulations were performed in the NVT scheme at a given temperature ranging from 400 K to 2000 K, and relaxed for enough time to bring the boxes to thermal equilibrium. The second stage started with periodic boundary conditions in the z-direction canceled. Thus, every simulation box contained two free surfaces in the opposite z-direction, and the He atoms could be treated as a uniform depth distribution initially. Using a time step of 1 fs, this stage lasted 500 ps with system evolving freely. To diminish the uncertainties caused by thermal fluctuation and for the sake of analyzing the evolution information, the simulation boxes were recorded and quenched to zero temperature every 5 ps. Considering the influence of surfaces in z-direction abruptly freed in this stage, simulation boxes in which He atoms have already escaped out of surfaces or been trapped at t = 5 ps were discarded, and only rest simulation boxes, namely surviving simulation boxes, were analyzed in the following section. Accordingly, t = 5 ps was taken as the starting time.
Na (t ) = NB 1
h h
dz
n 0 (z ) (E1 + E2 + E3 + E4 ) 2
(1)
where
E1 =
(1
2R) kerf
h + ( 1) k (2kh 4Dt
(1
2R) k erf
h + ( 1) k + 1 (2kh 4Dt
(1
2R) kerf
h + ( 1) k + 1 (2kh + z ) 4Dt
(1
2R) kerf
h + ( 1) k (2kh + z ) 4Dt
k=0
E2 = k=0
E3 = k=0
E4 = k=0
z)
z)
The Na (t ) is the number of He atoms out of surfaces or trapped near surfaces, namely absorbed He atoms. The NB is the number of surviving simulation boxes. Thus, the Na (t )/ NB is the ratio of absorbed He atoms. The n 0 (z ) = 1/ LZ is the initial distribution of the He atoms, where LZ is the box size in the z-direction. The h = LZ /2 is the distance from the absorbing layer to the center of the box. The time dependence of the ratio of absorbed He atoms for T = 400–2000 K is shown in Fig. 3. The results indicate that while increasing significantly with temperature from 400 K to 1000 K, the ratio of absorbed He atoms changes little with higher temperature T = 1200–2000 K for both surfaces. Especially, the absorbed ratio at T = 1800 K is larger than that at T = 2000 K. Using the value of the diffusion coefficient D obtained by Wang et al. [23], the absorbing probability, R, for W (2 1 1) and (3 1 1) surfaces was extracted by fitting Eq. (1) to the MD simulation data in Fig. 3. The obtained R listed in Table 1 varies around 1 with temperature, and has a maximal value at T = 400 K for both surfaces. By using the NEB method, the potential barrier near the (2 1 1) and (3 1 1) surface, where He atoms escape into the absorbed layer, is calculated. Generally, the potential barrier near the surface is slightly larger than in the bulk. However, there is an exception. When a He atom rounds the trap mutation sites and escapes into the vacuum directly, the potential barrier is smaller. We also investigated the fraction of quick-escaping He atoms in all atoms absorbed by the absorbing layer. The result shown in Fig. 4 indicates that the fraction increases with temperature generally but not monotonously. Similarly for (1 1 1) surfaces at T = 1000–1300 K, the fraction also shows a decreasing trend for (2 1 1) surfaces at T = 1600–1800 K and (3 1 1) surfaces at T = 1400–1600 K. It is very interesting that the fraction for (3 1 1) surfaces is obviously lower than that for (1 1 1) and (2 1 1) surfaces. This phenomenon agrees with the significant aggregation behavior of He atoms near (3 1 1) surfaces shown in Figs. 1 and 2, even at T = 1000 K. The behaviors of the trapped He below about 2 layers of the W (2 1 1) surface are very special, as they are not immobile. This unexpected scene is presented in Fig. 5. In this figure, we analyze the migration process of a trapped He atom thoroughly at T = 400 K. At the beginning, one He atom is trapped below about 2 layers of the W (2 1 1) surface, and pushes a W atom in the second layer under the surface out of its lattice site. The W atom then goes out of the W surface as an adatom. Thus, a pair of one He atom and one W adatom forms. It should be noted that although in a deeper layer, a He atom can occupy the
3. Simulation results According to the previous work on W (0 0 1), (0 1 1) and (1 1 1) surfaces, temperature and surface orientation were indicated to have an important effect on the migration of He atoms [18]. Following the method described above, we performed the MD simulation to explore the roles of temperature and surface orientation on the He behaviors near W (2 1 1) and (3 1 1) surfaces. Fig. 1 shows local merged 1000 He sites near W (2 1 1) and W (3 1 1) surfaces at the end of the simulations (500 ps) for temperatures of T = 400 K, 1000 K, 2000 K, respectively. It is worth noting that only one of the two free surfaces was displayed in Fig. 1. Color relates to the value of centrosymmetry parameter (CSP) [24,25]. He sites and corresponding CSP values indicate He atoms may not escape out of substrate instantly but accumulate near both surfaces. This is due to the formation of substitutional (interstitial) He and stacking (interstitial) W atom, which is also called the trap mutation behavior. Typically at T = 400 K, trap mutation behaviors occur at two sites separately 0.51 a0 and 0.66 a0 below the (2 1 1) surface. It should be noted that He atom does not completely occupy the initial site of the stacked W atom but is just a little below the site when the trap mutation occurs close to the surface. Fig. 2 further validates trap mutation behaviors, which provides snapshots of the depth distributions of He atoms corresponding to Fig. 1, but presents both two free surfaces. The coordinate origin is the center of the boxes. Although the trap mutation occurs near both surfaces, there are some differences. Especially at T = 1000 K, there are few He atoms trapped near (2 1 1) surfaces, similar to (1 1 1) surfaces, while still a huge number near (3 1 1) surfaces. Considering the success of describing He release from W (0 0 1), 3
Nuclear Inst. and Methods in Physics Research B 455 (2019) 1–6
J. Cui, et al.
(a)
(b)
Fig. 2. Histogram of the depth distribution of He atoms at the end of the simulation time at T = 400 K, 1000 K and 2000 K, respectively. (a) (2 1 1) surfaces; (b) (3 1 1) surfaces.
(a)
(b)
Fig. 3. The accumulated number of He atoms absorbed by the absorbing layer as a function of time at T = 400–2000 K. (a) (2 1 1) surfaces; (b) (3 1 1) surfaces. Table 1 The R values obtained by fitting Eq. (1). The diffusion coefficient used here is from Ref. [23]. Temperature (K)
Diffusion coefficient (cm2/ s)
R for W (2 1 1)
R for W (3 1 1)
400 600 800 1000 1200 1400 1600 1800 2000
1.17E−05 4.33E−05 7.16E−05 9.34E−05 1.07E−04 1.19E−04 1.27E−04 1.45E−04 1.59E−04
1.171 0.974 0.924 0.967 0.988 0.962 0.973 0.938 0.888
1.218 1.067 0.963 0.935 1.033 0.977 0.922 0.952 0.858
lattice site of the W atom pushed out, both the He atom and the W atom are not in the lattice site in the present situation. Here the second layer is about 0.39 a0 below the surface, and the most stable site S1 of the trapped He below 2 layers of the W (2 1 1) surface is about 0.66 a0. Besides, there is another main stable site S1' just around the S1 site, about 0.51 a0 below the surface. The potential energy in the S1' site is about 0.1 eV larger than that in the SI site. Even though with the least potential energy in the S1 site, the He atom could not be trapped tightly, but migrates frequently between SI and S1', accompanied with the W adatom moving. The barrier from S1' to SI is about 0.03 eV. Furthermore,
Fig. 4. The percentage of He atoms escaping into vacuum of all the absorbed He atoms. The data of (1 1 1) surfaces are obtained from Ref. [18].
an approximately one-dimensional migration of a He atom may take place along the 〈1 1 1〉 directions, as shown in Fig. 5. It has been demonstrated that one W adatom can migrate along the 〈1 1 1〉 directions on the W (2 1 1) surface when the substrate is prefect. But here when the He atom in the stable site (normally S1') migrates about 0.87 a0. to 4
Nuclear Inst. and Methods in Physics Research B 455 (2019) 1–6
J. Cui, et al.
Fig. 5. At T = 400 K, merged snapshots of simulation boxes during a trapped He atom diffuses along the 〈1 1 1〉 direction near the W (2 1 1) surface. Dark blue spheres: W atoms, green sphere: the initial site of the He atom, red sphere: the last site, gray sphere: stable sites and light blue dot: other He sites. S1' and S1 are two main stable sites.
(1 1 1) surfaces near which only a small fraction of He atoms could quickly escape into vacuum, a kind of trap mutations also occur for the (2 1 1) and (3 1 1) surfaces. In addition, the latter has the most obvious trap mutations phenomenon. These appearances depend on temperature, and the total absorbed probability can be well described by the continuous diffusion theory. For the (2 1 1) surfaces, when a He atom is trapped about 2 layers below the surface, an approximately one-dimensional migration of a He atom may occur along the 〈1 1 1〉 directions subsequently. These results are helpful for understanding the He release of tungsten surfaces and the near surface structural evolution, and also can be useful for guiding models in the methods of larger timespace scales.
Table 2 Ea and A obtained from Eq. (2). Surface orientation
Trapped site
Depth (a0)
Ea (eV)
A (ps−1)
(2 1 1)
S1 S2 S3 S4
0.66 1.21 0.67 1.49
0.15 1.85 0.95 1.28
0.11 24.1 22.1 1.15
(3 1 1)
another stable site, the previous W adatom just returns to its lattice site and a new W atom is pushed out. The migration barrier is about 0.9 eV. Strangely, this value is obviously larger than the barrier when He escape from the surface (0.22 eV obtained from the DFT calculation [16] and 0.15 eV from our MD results). However, when a He atom is trapped in the deeper place, or under other W surfaces ((1 1 1) and (3 1 1), for example), the similar migration behavior does not occur. Thus, to get a more in-depth understanding of this phenomenon, further work of relevant physical processes with DFT seems to be necessary. Finally, we obtain the escaping energy of several typical trapped sites, by proposing an MD simulation method that simulates the thermal desorption experiments. These sites are S1 and S2 (about 1.21 a0 below the surface) for the (2 1 1) surface, as well as S3 and S4 (about 0.67 a0 and 1.49 a0 below the surface separately) for the (3 1 1) surface. For each trap type, 500 independent boxes including a He atom on the trapped site were prepared. In this method, a linearly increasing-temperature (LIT) scheme was applied to heat the boxes with different value β [26]. Thus, the temperature of the simulation boxes at time t is T(t) = T0 + βt, where T0 is the initial temperature. We assume that the follows the Arrhenius equation: escaping frequency Rs Rs (T ) = Aexp ( Esa/ kB T ) , where Esa is the activation energy, and A is the pre-factor. Via gaining the temperature Tm at which the thermal desorption spectrum has a maximum value, the following function can be fitted:
2lnTm
ln =
Esa kB Tm
ln(Esa/ kB A)
Acknowledgements This work was partially supported by the National Natural Science Foundation of China (Grant No. 51501119 and 11775151). References [1] S. Takamura, N. Ohno, D. Nishijima, S. Kajita, Formation of nanostructured tungsten with arborescent shape due to helium plasma irradiation, Plasma Fusion Res. 1 (2006) 051, https://doi.org/10.1585/pfr.1.051. [2] İ. Tanyeli, L. Marot, D. Mathys, M.C.M. van de Sanden, G. De Temmerman, Surface modifications induced by high fluxes of low energy helium ions, Sci. Rep. 5 (2015) 9779, https://doi.org/10.1038/srep09779. [3] A. Al-Ajlony, J.K. Tripathi, A. Hassanein, Low energy helium ion irradiation induced nanostructure formation on tungsten surface, J. Nucl. Mater. 488 (2017) 1, https://doi.org/10.1016/j.jnucmat.2017.02.029. [4] G. De Temmerman, K. Bystrov, R.P. Doerner, L. Marot, G.M. Wright, K.B. Woller, D.G. Whyte, J.J. Zielinski, Helium effects on tungsten under fusion-relevant plasma loading conditions, J. Nucl. Mater. 438 (2013) 578, https://doi.org/10.1016/j. jnucmat.2013.01.012. [5] M.J. Baldwin, R.P. Doerner, Helium induced nanoscopic morphology on tungsten under fusion relevant plasma conditions, Nucl. Fusion 48 (2008) 035001, , https:// doi.org/10.1088/0029-5515/48/3/035001. [6] S. Kajita, W. Sakaguchi, N. Ohno, N. Yoshida, T. Saeki, Formation process of tungsten nanostructure by the exposure to helium plasma under fusion relevant plasma conditions, Nucl. Fusion 49 (2009) 095005, , https://doi.org/10.1088/ 0029-5515/49/9/095005. [7] S.B. Gilliam, S.M. Gidcumb, N.R. Parikh, D.G. Forsythe, B.K. Patnaik, J.D. Hunn, L.L. Snead, G.P. Lamaze, Retention and surface blistering of helium irradiated tungsten as a first wall material, J. Nucl. Mater. 347 (2005) 289, https://doi.org/ 10.1016/j.jnucmat.2005.08.017. [8] F. Ferroni, K.D. Hammond, B.D. Wirth, Sputtering yields of pure and helium-implanted tungsten under fusion-relevant conditions calculated using molecular dynamics, J. Nucl. Mater. 458 (2015) 419, https://doi.org/10.1016/j.jnucmat.2014. 12.090. [9] F. Sefta, N. Juslin, B.D. Wirth, Helium bubble bursting in tungsten, J. Appl. Phys. 114 (2013) 243518, , https://doi.org/10.1063/1.4860315. [10] A.M. Ito, Y. Yoshimoto, S. Saito, A. Takayama, H. Nakamura, Molecular dynamics simulation of a helium bubble bursting on tungsten surfaces, Phys. Scr. T159 (2014)
(2)
Then the values of Esa and corresponding A could then be obtained, and the results are listed in Table 2. 4. Conclusion Using MD simulations, we presented our further efforts to explore the migration behavior of single He atoms near W(2 1 1) and (3 1 1) surfaces at temperatures ranging from 400 K to 2000 K. Similar to W 5
Nuclear Inst. and Methods in Physics Research B 455 (2019) 1–6
J. Cui, et al. 014062, , https://doi.org/10.1088/0031-8949/2014/T159/014062. [11] M. Li, J. Cui, J. Wang, Q. Hou, Effect of cumulative helium bombardments on the surface and structural properties of tungsten, Nucl. Instrum. Methods Phys. Res. Sect. B Beam Interact. Mater. Atoms. 352 (2015) 92, https://doi.org/10.1016/j. nimb.2014.12.023. [12] K.D. Hammond, B.D. Wirth, Crystal orientation effects on helium ion depth distributions and adatom formation processes in plasma-facing tungsten, J. Appl. Phys. 116 (2014) 143301, , https://doi.org/10.1063/1.4897419. [13] Y. You, Clustering of H and He, and their effects on vacancy evolution in tungsten in a fusion environment, Nucl. Fusion 54 (2014) 103007, , https://doi.org/10.1088/ 0029-5515/54/10/103007. [14] D. Perez, L. Sandoval, S. Blondel, B.D. Wirth, B.P. Uberuaga, A.F. Voter, The mobility of small vacancy/helium complexes in tungsten and its impact on retention in fusion-relevant conditions, Sci. Rep. 7 (2017) 2522, https://doi.org/10.1038/ s41598-017-02428-2. [15] K.D. Hammond, S. Blondel, L. Hu, D. Maroudas, B.D. Wirth, Large-scale atomistic simulations of low-energy helium implantation into tungsten single crystals, Acta Mater. 144 (2018) 561, https://doi.org/10.1016/j.actamat.2017.09.061. [16] G.Y. Pan, Y.G. Li, Y.S. Zhang, C.G. Zhang, Z. Zhao, Z. Zeng, First-principles investigation of the orientation influenced He dissolution and diffusion behaviors on W surfaces, RSC Adv. 7 (2017) 25789, https://doi.org/10.1039/C7RA03281A. [17] A.V. Barashev, H. Xu, R.E. Stoller, The behavior of small helium clusters near free surfaces in tungsten, J. Nucl. Mater. 454 (2014) 421, https://doi.org/10.1016/j. jnucmat.2014.08.033. [18] X. Wang, Z. Wu, Q. Hou, Effects of temperature and surface orientation on migration behaviours of helium atoms near tungsten surfaces, J. Nucl. Mater. 465 (2015) 455, https://doi.org/10.1016/j.jnucmat.2015.06.011.
[19] S.A. Trygubenko, D.J. Wales, A doubly nudged elastic band method for finding transition states, J. Chem. Phys. 120 (2004) 2082, https://doi.org/10.1063/1. 1636455. [20] D. Sheppard, P. Xiao, W. Chemelewski, D.D. Johnson, G. Henkelman, A generalized solid-state nudged elastic band method, J. Chem. Phys. 136 (2012) 074103, , https://doi.org/10.1063/1.3684549. [21] Q. Hou, M. Li, Y. Zhou, J. Cui, Z. Cui, J. Wang, Molecular dynamics simulations with many-body potentials on multiple GPUs – The implementation, package and performance, Comput. Phys. Commun. 184 (2013) 2091, https://doi.org/10.1016/ j.cpc.2013.03.026. [22] G.J. Ackland, R. Thetford, An improved N-body semi-empirical model for bodycentred cubic transition metals, Philos. Mag. A. 56 (1987) 15, https://doi.org/10. 1080/01418618708204464. [23] J. Wang, Y. Zhou, M. Li, Q. Hou, Atomistic simulations of helium behavior in tungsten crystals, J. Nucl. Mater. 427 (2012) 290, https://doi.org/10.1016/j. jnucmat.2012.05.020. [24] C.L. Kelchner, S.J. Plimpton, J.C. Hamilton, Dislocation nucleation and defect structure during surface indentation, Phys. Rev. B. 58 (1998) 11085, https://doi. org/10.1103/PhysRevB.58.11085. [25] R.E. Rudd, Void growth in bcc metals simulated with molecular dynamics using the Finnis-Sinclair potential, Philos. Mag. 89 (2009) 3133, https://doi.org/10.1080/ 14786430903222529. [26] J. Cui, Z. Wu, Q. Hou, Estimation of the lifetime of small helium bubbles near tungsten surfaces – A methodological study, Nucl. Instrum. Methods Phys. Res. Sect. B Beam Interact. Mater. Atoms. 383 (2016) 136, https://doi.org/10.1016/j.nimb. 2016.07.001.
6