1. Introduction

The controllable fusion energy is presently considered as one of the most promising green energy, which will be competent to solve the energy crisis in the future because of the abundant fuels, deuterium (D) and tritium (T) in oceans. On the other hand, uncontrollable fusion energy is applied in the design of fusion weapons. They perform with hydrogen (H) isotopes as the fuels. It is widely known that the T is radioactive and it will decay into He isotope (³He) as the final products with a half-life of about 12.346 years [1]. Once the T atoms penetrate into materials or nuclear fuels, they will raise He, H, and H-He synergic effects in materials[1-5]. The hydrogen and their isotope atom may greatly influence materials [6]. It could threaten the safety and operating life of T relevant facilities, such as first wall and T storage vessels [7-9]. He atoms are hardly soluble in metals, and they can migrate easily therein. Wide investigation has shown that interstitial He atoms can be deeply trapped at pre-existing traps, such as vacancies, voids, dislocations, and grain boundaries (GBs), to nucleate bubbles [10-12]. The accumulation of He atoms in materials may result in severe degradation of macroscopic properties of T storing materials, such as the volume swelling, work hardening, and intra- and inter-granular embrittlement[13,14]. These phenomena will lead to the failure of T storage materials.

T atoms penetrate into T storage materials at extremely high atomic concentrations. It is necessary for H storage materials to release T atoms as completely as possible. At the same time, the T decay produced He atoms should not be released from the materials. As a consequence, the He retaining property is of great importance for materials in the application in T storage. Pd and its alloys are commonly used for H isotope’s processing and storage due to their high ability to retain the ³He generated by T decay[15]. The primary He behaviors have been widely studied with experimental and computational methods in Pd. The segregation of He atoms have been widely observed and they tend to bind to each other to form He clusters. These He clusters could grow up by continually absorbing He atoms and activating defects including self-interstitial atoms, dislocation loops, and stacking fault in materials [16-18]. The presence of H isotopes may greatly influence the primary behaviors of He in Pd crystal, a face-centered crystal (FCC). It is of great meaning to study the effects of H isotopes on the He behaviors in Pd. As the development of high performance computers increases, computational methods have become more and more popular in material science society. Ab initio calculations have been widely performed to calculate primary defect behaviors from the point of electron character [19, 20]. Molecular dynamics method is used to model the clustering behaviors of defects [21]. The accuracy of the interatomic potential plays a key role in atom behaviors investigation.

In this work, a set of ab initio calculation-based interatomic potentials have been fitted to describe the interaction between Pd, H, and He atoms. The molecular dynamics (MD) method is utilized to study the effects of H on He behaviors in order to analogically reveal the T impacts on He behaviors in Pd. We show the relevant fitting techniques, potential testing results, and all parameters of relevant models in the second part of the work. The potential fitting results are given and discussed in Sect. 3 in details. Meanwhile, the computational results on H effects on He behaviors have also been presented and discussed in this part. Finally, we clearly conclude the effects of H on He behaviors in Pd in Sect. 4.

2. Models and methods

In order to investigate the effects of H atoms on the He behaviors, a set of interaction potentials should be available to MD method. The interactions between Pd atoms are described with the modified analytical embedding atomic (MAEAM) potential developed by Hu et al [22]. As is referenced, the L-J potentials are applied to describe the He-He, H-H, and H-He interactions in Pd [23]. Both the Pd-He and Pd-H interactions are described with the Morse potential [24], as follows:

where the parameters ${\varnothing }_0$，$\text{α}$, and ${\text{r}}_0$ are adjustable parameters which need to be accurately fitted. The least square method is a commonly used method for the development of interatomic potentials, and the method is also applied to fit the Pd-H potential. The formula can be expressed as follows:

In Eq. (1), the parameter U is the object function and the adjustable parameters should be the best one if U obtains the smallest value. F_i are reference values obtained from ab initio calculation or experimental results. The ${\text{ω}}_i$ represents the weight value of each equation used in fitting processes and $f_i\left(\lambda \right)$ is the real value of our functions. It is well known that the solution of a group of equations with three variables needs no less than three equations. As a consequence, we have taken six H defects into consideration in Pd, and the configurations of H-participated defects are shown in Fig. 1. Since the position occupation is of great importance for impurities in host lattice, the formation energy of H atoms at tetrahedral and octahedral interstitial sites are consequently taken as important target values. All formation energy is constructed and calculated with ab initio calculation method. We have introduced the generalized gradient approximation (GGA) for the treating of the exchange-correlation energy [25]. The widely used projector augmented wave pseudopotentials (PAW) is utilized in the calculations [26]. Since the size of the supercell may affect the results of the defect formation energy, ab initio calculations are performed with supercells at different sizes and the last size of the supercell is set as 3a₀×3a₀×3a₀, which contains about 108 atoms. The cut-off energy is set as 500 eV, and the k-point mesh is tested and set as 5×5×5 at last. These results are taken as reference values in the fitting of Pd-H interatomic potential. The quenching optimization method is used to search the best parameters.

Fig. 1Full-screen View

The configuration of H clusters in Pd at size of 2 to 4 corresponding to (a), (b), (c) and (d), where the blue and orange balls represents Pd and H atoms.

We have investigated the diffusion and releasing behaviors of He atoms in H contained Pd crystal with MD method implemented in a MOLDY codes [27]. The computational models are briefly shown in Fig. 2 and 3 for He diffusion and releasing in Pd. The size of the box is set as 10a₀×10a₀×10a₀, which contains more than 4000 atoms in total, and a₀ is the lattice constants of Pd at computing temperatures. The diffusion simulations are performed at the temperature of 300, 600, 900, and 1200 K with 3D periodic boundary conditions. The box sizes of He releasing are set as 12a₀×12a₀×25a₀. The periodic boundary conditions are applied in the direction of x and y, while a surface is introduced in z direction. A void is introduced in Pd at the depth of 3a₀ with the radius equaling to one nanometer, and then the void is filled with He atoms at given density to construct the He bubble releasing models in Pd membrane. Randomly distributed H atoms are introduced into Pd at the atomic concentrations of 10%, and the releasing temperature is set to 300 K. The time steps equal to 1 fs and all simulations are performed in volume and temperature constant NVT ensemble.

Fig. 2Full-screen View

The He diffusion models in pure Pd and H contained Pd, where the blue, green and red balls represent Pd, He, and H atoms respectively. The He atom is initially sited at octahedral interstitial site in Pd.

Fig. 3Full-screen View

The He releasing model from the Pd surface. The green and red balls are host Pd and He atoms respectively. The depth of the He bubble from the surface to bulk is indicated with the variable h.

3. Results and discussion

In the first column of Table 1, $E_f^T$, $E_f^o$ represent the formation energy of tetrahedral and octahedral interstitial hydrogen atoms, respectively. $E_f^{{\text{H}}_2},$ $E_f^{{\text{H}}_3}$, $E_f^{{\text{H}}_4}$, $E_f^{{\text{H}}_5}$ mean the formation energy of H₂, H₃, H₄, and H₅ clusters, respectively.

A set of accurate interatomic potential is necessary to investigate H effects on He behaviors in Pd. We have fitted the Pd-H interatomic potential based on ab initio calculations. The formation energy of H defect is calculated using the ab initio and the results have shown that the configurations of H clusters have high geometry symmetries. As shown in Fig. 1(a), H atoms in the H₂ cluster stay at a pair of second nearest neighbor octahedral interstitial sites. H₃ cluster has a configuration of equilateral triangle with the normal vector pointed to <111> direction in Pd. The four H atoms in H₄ cluster construct a regular tetrahedron. The configuration of H₅ clusters is also clearly shown in Fig. 1 (d). Further calculations have shown that a single H atom favorably stay at the octahedral interstitial sites in Pd with the formation energy of about -3.46 eV, which is a little lower than the formation energy of tetrahedral interstitial H (-3.41 eV). All ab initio results are listed in the last column within Table 1 and specific parameters are listed in Table 2. The occupation order of H in Pd is very important, which is taken into consideration in the Pd-H potential fitting processes. The fitting results have also been listed in the second column of Table 1. It can be seen from the Table that present results agree well with the target DFT values.

We have carefully analyzed the mean square displacement (MSD) of the He atom in Pd doped with H at the concentration of 15% (at.) under temperatures lower than 1000 K. It is found that the He atom diffuses little in Pd with 15% H contained at 300, 600, and 900 K. He atom keeps on vibrating nearby its initial sites until the temperature increases to 1200 K. In order to study the effects of H content on the diffusion range of He atom, the predominant attentions have been focused on the diffusing trajectory of He atom in models of pure Pd and Pd with H at the concentrations of 5% and 15%. Comparisons are made and the results are shown in Fig. 4. It can be clearly seen from Fig. 4 that He atoms in pure Pd diffuses quickly and the trajectory pervades the whole box. The He atom in Pd doping with 5% (at.) H atoms diffuse much less than that in pure Pd, which has the smallest diffusing trajectory. With the increase of H concentration, it is difficult for He atom to diffuse in Pd, which indicates that the presence of H in Pd efficiently obstacles the migration of He atoms. It can be seen from the formation energy of H defects that the introduction of H atoms in Pd greatly decreases the total energy of the system, which leads to the stabilization of the system. As shown in Fig. 5, a He atom can be trapped by H atoms on or near its diffusion paths. Furthermore, it is noted from the ab initio calculation results that H atoms are strongly bond to Pd atoms, and they severely constrict the migration passageway of He atoms, which efficiently elevates the migration barriers of He atom. Once the He atoms meet with H atoms, its diffusion will be transiently restrained, which can be evidently seen from the MSD curve.

Fig. 4Full-screen View

The diffusion range of a single interstitial He atom in Pd with H at atomic concentrations of 0.0% (a), 5% (b) and15% (c).

Fig. 5Full-screen View

The diffusing MSD curve of a He atom in Pd with H at the concentration of 9.0% at. as the function of time. The insert figure is the local microstructure of models corresponding to the time segments shown with dash lines.

The He releasing behaviors is of great significance for Pd applied as the He retaining or T storage material. Contrastive investigation on He releasing have been performed in pure and in H doped Pd. The helium releasing processes have been shown in Fig. 6. Subfigures in both rows orderly represent the snapshots of helium releasing levels at the moments of 0, 0.5, 5.0, 10, 200, 250, and 500 ps from left to right. It can be seen that He atoms begin to release at 0.5 ps from the surface of pure Pd membrane. The surface starts to crack with a few He atoms released at 5 ps. As the computational time increases, more and more He atoms are released. He atoms are particularly released at 200 ps, and most of the helium atoms are released at 500 ps. Analysis on the microstructure of the surface have been made in the processes of He release and it is found that the surface quickly deforms at 0.5 ps due to the high press from the near-surface He bubbles. Then a brittle fracture of the surface has been observed accompanying with the explosive release of helium atoms, and the accelerating release of helium bubble occurs at 200 ps. The radius of the chasm remains to be smaller than a few lattice constants. In fact, the rate of helium release depends both on the amount density of helium atoms in the bubble and the depth of the bubble. Furthermore, the temperature greatly influences the releasing rate of helium bubble in pure Pd membrane.

Fig. 6Full-screen View

The releasing of He in H free and H containing Pd with surface, sub figures in each line represents configuration snapshots at time of 0, 0.5, 5.0, 10, 200, 250, and 500 ps in the order from left to right.

Further study is done to investigate the effects of bubble composition on the releasing behavior of He from Pd film surface, and the results indicate that the pure He bubble releases from the surface, but the H-He mixed bubble doesn’t release while the depths remains equal. The mixing of H atoms with He into the bubble efficiently low down the intra-press of gas bubble. It results in the fact that the intra-press of the mixing bubble is not high enough to crack the surface of Pd membrane. H atoms are introduced into the Pd membrane at the atomic concentration of 15%, and the release of He bubbles has also been simulated at 300 K. Compared with the pure Pd membrane, the He releases slightly slower. The surface begins to deform at 5.0 ps and no He atoms are released from the surface. The accelerating release of He atom delays about 50 ps in H-doped Pd membrane. It can also be seen from Fig. 6, that He atoms release from the surface of pure Pd membrane in the form of clusters. On the contrary, the He releases from H-doped Pd film one by one. Microstructure analysis shows that the size of the crack gap of H-doped Pd membrane is smaller than that in pure Pd film. The reason is that the introduction of H in Pd membrane slightly decreases the fragility of the host material which leads to separate release of He atoms from surface. The H atoms decrease the total energy of the system and increase the releasing barrier of near-surface He atoms. The strong bonding of H and Pd atoms increases the surface tension and minimizes the releasing passageway of He atom, which causes the slower release of He atoms from H-doped Pd membrane.

4. Conclusion

The Pd-H interatomic potential function is accurately fitted based on ab initio calculations with the Morse formula. The effects of H on the He behaviors in Pd have been be carefully investigated with the interatomic potential. Both the effects of H on He diffusing and releasing behaviors are simulated with MD method. It is noted that the ab initio based interatomic potential accurately describes the interaction between Pd atoms and H atoms. The introduction of H atoms in Pd greatly slows down the diffusion of He in Pd and the release of He bubble from the H-doped membrane can also be slightly hindered.