# Diffusivity in metallic beryllium: The case of the H, C, N and O species

Damien Connétable, Aurélien Riot, François Virot

► To cite this version:

Damien Connétable, Aurélien Riot, François Virot. Diffusivity in metallic beryllium: The case of the H, C, N and O species. Computational Materials Science, 2025, 246, pp.113455. 10.1016/j.commatsci.2024.113455. hal-04753025

HAL Id: hal-04753025

https://hal.science/hal-04753025v1

Submitted on 7 Apr 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

Full length article

# Diffusivity in metallic beryllium: The case of the H, C, N and O species

Damien Connétable<sup>a,*</sup>, Aurélien Riot<sup>a,b</sup>, François Virot<sup>b</sup>

<figure><img src="image_3.png" /></figure>

## ARTICLE INFO

### Keywords:

Beryllium

Diffusion

DFT calculations

## ABSTRACT

The insertion and diffusivity of interstitial species in metallic beryllium are discussed in this work using a multi-scale methodology, coupling first-principles calculations and a multi-site approach. Emphasis is placed on the main interstitial species, i.e., H, C, N and O atoms. The results show that the most stable site is strongly dependent on the nature of the interstitial atom. Indeed, carbon is most stable in octahedral sites and hydrogen tetrahedral sites, while oxygen and nitrogen are most stable in basal tetrahedral sites. From the stable insertion sites and symmetrical saddle points identified, a number of migration pathways were mapped. The diffusion pathways were then completed using nudged elastic band (NEB) calculations. The diffusivity of the atoms shows an isotropic behavior as expected for carbon, which shows a strong anisotropic behavior with faster diffusivity along the basal plane. These theoretical results are in agreement with known experimental data, especially for hydrogen diffusion. The discrepancy between theory and experiment is corrected by the vacancy trapping effect. Finally, this paper theoretically determines the Arrhenius parameters of each diffusing species, in particular carbon and oxygen, for which no data were available.

## 1. Introduction

Beryllium is a lightweight metal with a combination of chemical, physical, mechanical, and nuclear properties that make it suitable for a wide range of high technology applications in the industrial, aerospace, defense, medical, and nuclear fields. For example, due to its stiffness, lightness, and structural stability, beryllium is the structural material of the primary mirror of the James Webb Space Telescope [1]. Beryllium is often alloyed with copper or nickel, resulting in a great combination of properties such as strength, electrical and thermal conductivity [2]. Due to its low atomic number, beryllium can be used, for example, as an X-ray window in the medical field. Beryllium has several applications in the nuclear industry, especially in the field of nuclear fusion. In fact, it could be used as the plasma-facing material in the vacuum vessel of the ITER reactor, mainly due to its compatibility with plasma as well as its good thermal and mechanical properties [3]. In the same field, beryllium pebbles could be part of the preliminary design of the tritium breeding blanket as a neutron multiplier for the nuclear fusion plant [4].

Unfortunately, all these tremendous properties can be compromised by impurities, defects and phase inclusion. The unexpected phase formation in metallic materials, such as oxide, carbide, etc., is controlled by the amount and kinetic rate of the elements that can lead to these compounds. It is therefore important to understand the behavior of impurities in metallic beryllium, in particular their diffusivity, which

plays an important role. Oxygen and carbon are the major impurities in structural and nuclear beryllium, about 1 and 0.1% respectively [5]. These species could then form the well-known beryllium oxide (BeO) and beryllium carbide (Be₂C). In order to purify the plasma atmosphere in the ITER machine, the former is expected to be grown on the first wall of the vacuum vessel. As regards the behavior of hydrogen in metals, it can cause embrittlement of the material. In the case of beryllium, its hydride is stable at low temperatures [6]. In addition, for hydrogen isotopes in beryllium, tritium retention is a major safety issue. Finally, in the ITER facility, nitrogen puffing is being considered to reduce the thermal load on the inner wall of the vacuum vessel.

Thus, for both technological industries and the nuclear fusion field, accurately determining the diffusion of oxygen, carbon, hydrogen and nitrogen in metallic beryllium is key in order to improve the knowledge of the mechanisms governing the formation of phase inclusions as well as their permeation. As far as the diffusivity of these interstitial species is concerned, many works have already covered that of hydrogen, using both experimental [7–9] and theoretical approaches [10–13]. This work is being driven by concerns about tritium retention in plasma-facing components. The energy landscape of interstitial oxygen and carbon (essential in order to determine the diffusion coefficients) has only been studied theoretically [14]. To our knowledge, the behavior of interstitial nitrogen in metallic beryllium has never been investigated

<footnote><sup>a</sup> CIRIMAT, CNRS, INPT, 4 allée Émile Monso, BP 44362, F-31030 Toulouse Cedex 4, France</footnote>

<footnote><sup>b</sup> Institut de Radioprotection et de Sûreté Nucléaire (IRSN), PSN-RES/SAM/LETR, F-13115 Saint-Paul-lez-Durance, France</footnote>

<footnote>* Corresponding author at: CIRIMAT, CNRS, INPT, 4 allée Émile Monso, BP 44362, F-31030 Toulouse Cedex 4, France. E-mail address: damien.connetable@ensiacet.fr (D. Connétable).</footnote>
and a thorough experimental program studying this matter has long been expected [15].

In this paper, the behavior of carbon, hydrogen, oxygen and nitrogen in metallic beryllium was analyzed in terms of insertion and migration energy landscapes. Finally, diffusion coefficients were evaluated. Our results provide valuable insight as to how these impurities behave in beryllium. This paper is organized as follows: after a short presentation of the computational approaches used, the ground state properties of hexagonal closed-packed bare beryllium are briefly discussed. The next Section 4 focuses on characterizing the insertion of the interstitial species. Finally, Section 5 looks at their diffusive behavior and evaluates their diffusion coefficient.

## 2. Methodology

DFT calculations were performed using the Vienna Ab initio Simulation Package, VASP [16]. Simulations were carried out within the Perdew–Burke–Ernzerhof generalized gradient approximation of the exchange and correlation functional (PBE) [17] to the density functional theory [18,19]. Projector augmented wave (PAW) pseudopotentials [20] were used, and a kinetic energy cut-off of 600 eV for the plane-wave basis set was adopted to calculate energies. The pseudo-potential of beryllium includes all of its four electrons, while the pseudo-potentials of H, C, N and O atoms include one, four, five, and six electrons, respectively. Energies calculations were performed on 4 × 4 × 3 supercells, i.e. 96 atoms. The k-mesh grids were then sampled with Γ-centered 8 × 8 × 10 [21]. In order to compare the energies obtained at constant pressure, the computed data at constant volume were corrected using the elasticity theory, these corrections were calculated with the Aneto code [22].

To analyze the energy landscapes and thus identify the different transition states, climbing image NEB calculations [23] were performed using five intermediate positions. The energies of the transition states were further calculated with the same convergence criteria used for the insertion sites.

Additional calculations were performed to compute inter-atomic force constants (IFC). They allowed characterizing the stability of the configurations: the absence of imaginary frequency indicates that the sites are stable, while one (or two) indicates that the configuration is a first- (or second) order transition state. The finite displacement method using the phonopy package [24] was performed on 3 × 3 × 3 supercells, i.e., 54 atoms. The Helmholtz free energies, $F_v$,

$$F_v(T) = k_{BT} \sum_{v=1}^{3N} \int \ln \left[ 2 \sinh \left( \frac{\hbar \omega_{v,q}}{2k_B T} \right) \right] n_v(q) dq. \quad (1)$$

were then calculated. $\omega_{v,q}$ are the frequencies of wave vector $q$ and mode $v$. $n_v(q)$ is the phonon density of the wave vector $q$. Calculations were performed on 20 × 20 × 30 q-mesh grids.

To calculate the elastic tensor ($P$) and relaxation volume ($V$) (in dilute solution), we used the Aneto code developed by Varvenne et al. [25]. It uses the residual stress induced by an impurity when inserted into the perfect crystal, as well as the crystal’s elastic constants.

## 3. Ground-state properties of the Be system

The Be-hcp system belongs to space group 194 (hP2, P63/mmc). Be atoms are located in the 2c sites [26]. The ground state properties calculated with the proposed computational approach are gathered in Table 1, i.e. lattice parameters, $a$ and $c$ (in Å), the cohesive energy, $E_{coh}$ (in eV/atom), and the elastic constants, $C_{ij}$ (in GPa), calculated with a strain-stress method. The bulk modulus, $B_v$, is given by:

$$B_v = \frac{2C_{11} + C_{33} + 2C_{12} + 4C_{13}}{9} \quad (2)$$

<figure><img src="image_4.png" /><figcaption>Fig. 1. Schematic representation of the Be-hcp structure (green balls) with interstitial positions: TE, BT, OC, BO, CN, and BN. Figures were built using VESTA [27].</figcaption></figure>

<caption>Table 1</caption>

Lattice parameters, $a$ and $c$ in Å, and elastic constants, $C_{ij}$ in GPa, of the Be system. The cohesive energy, $E_{coh}$, is in eV per atom. Theoretical values were computed using the PBE functional.

<table><thead><tr><th></th><th>a</th><th>c</th><th>c/a</th><th>B<sub>v</sub></th><th>E<sub>coh</sub></th></tr></thead><tbody><tr><td>Exp. [28]</td><td>2.286</td><td>3.584</td><td>1.568</td><td>121</td><td>3.32</td></tr><tr><td>Theo. [13]</td><td>2.258</td><td>3.549</td><td>1.572</td><td>121</td><td>3.70</td></tr><tr><td>This work</td><td>2.264</td><td>3.566</td><td>1.575</td><td>122</td><td>3.72</td></tr><tr><td></td><td>C<sub>11</sub></td><td>C<sub>33</sub></td><td>C<sub>44</sub></td><td>C<sub>66</sub></td><td>C<sub>12</sub></td><td>C<sub>13</sub></td></tr><tr><td>Exp. [29]</td><td>293.6</td><td>356.7</td><td>166.2</td><td>133.4</td><td>26.8</td><td>14.0</td></tr><tr><td>This work</td><td>319.6</td><td>378.8</td><td>166.1</td><td>145.9</td><td>27.7</td><td>6.5</td></tr></tbody></table>

and $C_{66}$ is equal to $(C_{11}-C_{12})/2$. The cohesive energy $E_{coh}$ is obtained by calculating the difference between the energy of the isolated Be atom ($E_{dft[Be-at]}$) and that of the bulk per unit atom ($E_{dft[Be-hcp]}$):

$$E_{coh} = E_{dft[Be-hcp]} / 2 - E_{dft[Be-at]} \quad (3)$$

The DFT results are in agreement with experimental [28,29] and theoretical values [10,13,30,31]. The theoretical lattice parameters are about 1% lower than those found experimentally, while the cohesive energy is roughly reproduced. Regarding the values of the elastic tensor, it can be noted that there is a 10% difference, at most, between the theoretical values and those obtained experimentally [29]. However, the agreement between theory and experiment is quite good. These theoretical data will be used in the following section for the elastic theory (interaction energy, dipole tensors...).

## 4. Insertion sites

Six possible interstitial positions are reported in the hexagonal close-packed system, hP2. These positions are shown in Fig. 1. The octahedral (labeled OC) and tetrahedral (TE) positions occupy the 2a (0, 0, 0) and the 4f (2/3, 1/3, z) positions, respectively. z is usually about 0.375 in the ideal position. The projection of the OC and TE sites onto the basal plane is the basal octahedral position (denoted BO) in 2b (0, 0, 1/4) and the basal tetrahedral position (BT) in 2d (2/3, 1/3, 1/4). The crowdion configurations BC and CN are located between two neighboring Be atoms in position 6h (x, 2x, 1/4), where $x \approx 0.1666$, and 6g (1/2, 1/2, 1/2) respectively.

To compare the stability of sites, the insertion energy, $E_i$, and zero-point energies of stable positions were calculated. $E_i$ is given by:

$$E_i[Y] = E_{dft[Be+Y]} - E_{dft[Be]} - \mu_{ref}[Y] \quad (4)$$

where $E_{dft[Be+Y]}$ and $E_{dft[Be]}$ are the DFT energies of the system with and without interstitial species, respectively. $\mu_{ref}[Y]$ is the chemical potential of Y (per atom) in its reference state, i.e. the $Y_2$ molecule for
<caption>Table 2</caption>

Summary of computed data to evaluate insertion energies, in eV, of H, C, N, O atoms. The notations *r* and *nr* correspond to the cases of relaxed lattice parameters (constant pressure) or non-relaxed (constant volume) respectively. Energies are in eV. The correction term, *E*<sub>d</sub> (in eV), based on the elasticity theory, is computed using the Aneto code [22]). For the non-relaxed case, the insertion energy is calculated by *E*<sub>i</sub><sup>nr</sup> = *E*<sub>0</sub><sup>nr</sup> + ZPE − *E*<sub>d</sub>. The bold font marks the insertion site identified as the ground state. The labels † and ‡ refer to an unstable configuration with one or two imaginary frequencies, respectively.

<table><thead><tr><th></th><th>OC</th><th>BO</th><th>TE</th><th>BT</th><th>BC</th><th>CN</th></tr></thead><tbody><tr><td rowspan="4">H</td><td>E<sub>0</sub><sup>nr</sup>/E<sub>0</sub><sup>nr</sup></td><td>1.773/1.796</td><td>2.370/-</td><td><strong>1.589/1.623</strong></td><td>1.592/-</td><td>2.319</td></tr><tr><td>ZPE/E<sub>d</sub></td><td>0.005/0.004</td><td>†/-</td><td>0.073/0.007</td><td>†/-</td><td>‡</td></tr><tr><td>E<sub>i</sub><sup>nr</sup>/E<sub>i</sub><sup>nr</sup></td><td>1.778/1.796</td><td>-</td><td>1.663/1.689</td><td>-</td><td>-</td></tr><tr><td>H</td><td>+0.107</td><td>-</td><td>0.0</td><td>-</td><td>-</td></tr><tr><td rowspan="4">C</td><td>E<sub>0</sub><sup>nr</sup>/E<sub>0</sub><sup>nr</sup></td><td><strong>3.536/3.703</strong></td><td>4.756/4.895</td><td>-</td><td>4.293/4.514</td><td>3.894/4.061</td></tr><tr><td>ZPE/E<sub>d</sub></td><td>-0.162/0.040</td><td>†/-</td><td>-</td><td>-0.204/0.124</td><td>†/-</td></tr><tr><td>E<sub>i</sub><sup>nr</sup>/E<sub>i</sub><sup>nr</sup></td><td>3.374/3.501</td><td>-</td><td>-</td><td>4.089/4.184</td><td>-</td></tr><tr><td>H</td><td>0.0</td><td>-</td><td>-</td><td>+0.71</td><td>-</td></tr><tr><td rowspan="4">N</td><td>E<sub>0</sub><sup>nr</sup>/E<sub>0</sub><sup>nr</sup></td><td>1.446/1.600</td><td>3.572/-</td><td>-</td><td><strong>1.117/1.323</strong></td><td>2.155/2.405</td></tr><tr><td>ZPE/E<sub>d</sub></td><td>-0.115/0.031</td><td>‡</td><td>-</td><td>-0.050/0.093</td><td>-0.111/0.048</td></tr><tr><td>E<sub>i</sub><sup>nr</sup>/E<sub>i</sub><sup>nr</sup></td><td>1.331/1.454</td><td>-</td><td>-</td><td>1.067/1.180</td><td>2.044/2.246</td></tr><tr><td>H</td><td>+0.274</td><td>-</td><td>-</td><td>0.0</td><td>1.066</td></tr><tr><td rowspan="4">O</td><td>E<sub>0</sub><sup>nr</sup>/E<sub>0</sub><sup>nr</sup></td><td>-0.266/-0.105</td><td>0.974/-</td><td>-</td><td><strong>-1.666/-1.549</strong></td><td>0.185/-</td></tr><tr><td>ZPE/E<sub>d</sub></td><td>-0.147/0.028</td><td>†/-</td><td>-</td><td>-0.039/0.056</td><td>‡</td></tr><tr><td>E<sub>i</sub><sup>nr</sup>/E<sub>i</sub><sup>nr</sup></td><td>-0.413/-0.280</td><td>-</td><td>-</td><td>-1.705/-1.644</td><td>-</td></tr><tr><td>H</td><td>+1.425</td><td>-</td><td>-</td><td>0.0</td><td>-</td></tr></tbody></table>

H, O and N, and diamond for carbon. The zero-point energies, labeled ZPE, were calculated with the same approach, i.e.:

$$ZPE[Y] = F_v[\text{Be}+\text{Y}] - F_v[\text{Be}] - F_v^{\text{ref}}[Y] \quad (5)$$

where $F_v$ are the Helmholtz free energies computed at 0 K from Eq. (1). For H, C, N and O atoms, $F_v^{\text{ref}}[Y]$ values were extracted from David's paper [32]. Values are listed in Table 2.

As found and described in many metals [32–34], the stable sites vary greatly depending on the chemical nature of the inserted species and its interactions with the atoms of the system. Vibrations were analyzed to identify the configurations that were local minimal. Here, the Be-hcp system shows three different cases.

In the case of carbon, four configurations remain unchanged after atomic relaxation: OC, CN, BT and BO sites, in ascending order of energy. These results are in relatively good agreement with those obtained by Zhang et al. [10], except for the CN sites where our calculations lead to lower energy values. In the study of Zhang et al., the energy difference between BT, BC or BO with respect to the ground state (OC) is 0.72, 0.71 and 1.22 eV respectively. In this work, these values are 0.76 and 1.22 eV respectively (by comparing calculations at constant pressure). The BC sites are found to be unstable in terms of their vibrational properties. By analyzing the frequencies of these configurations, only two types of sites are local minima and can therefore be considered as carbon insertion sites: OC and BT sites. The results also show that the other two configurations, the BO and CN sites, have a single imaginary frequency. They are therefore first-order transition states. This point will be discussed in the next section.

For O and N atoms, four sites are found to be stable: the OC, BO, BT and CN sites. For the O atom, the frequency analysis shows that only the OC and BT sites are stable, as for carbon. For nitrogen, there are three stable sites:OC, BT and CN sites. For both species, the BT site is the most stable configuration. There is good agreement between these results and those obtained by Zhang et al. [10] for oxygen. There is only a small constant “shift” between the values probably due to the DFT value of the energy of the O<sub>2</sub> molecule which is different between the two papers. It is indeed interesting to compare energy differences. Zhang et al. found 1.5, 1.96 and 2.76 eV for OC, CN and BO compared to the BT site (which is also the main stable configuration). In our case, the values are 1.4, 1.85 and 2.64 eV. Similarly to carbon, the BO site should be a first-order transition site for the oxygen atom, whereas CN for oxygen or BO for nitrogen are identified as second-order.

The most debated case remains hydrogen. There are many papers in the literature dealing with hydrogen insertion and there is still some controversy about the most stable site. Previous work suggests that the

<figure><img src="image_5.png" /><figcaption>Fig. 2. Evolution of the energy as a function of the distance to the BT site in the direction of TE. Calculations were performed on 4 × 4 × 3 super-cells.</figcaption></figure>

most stable site for H atoms is either the TE [35] site or the BT site [10–12,30,31,36]. The present results clearly show that five sites are stable after relaxation: the OC, BO, BT, CN and TE sites. A detailed analysis of the positions, regardless of the size of the supercell (3 × 3 × 3, 4 × 4 × 3, 4 × 4 × 4 and 5 × 5 × 5 supercells) shows that the ground state is always the TE site. The final position, z ≈ 0.20, is very close to the BT site. The distance between the two sites is about 0.2 Å. Consequently, the energy difference between the BT and TE sites is very small: it is about 3 meV, in favor of the TE site. A frequency analysis shows that the “ideal” BT site has one imaginary frequency, like the BO site. These results were confirmed using different types of pseudo-potentials for hydrogen and Be with VASP, but also with the Quantum-Espresso code [37]. To complete this discussion, the energy landscape in the vicinity of BT and TE is shown in Fig. 2. The non-parabolic Mexican hat-like energy potential surface with a very small energy barrier is at the root of this debate. In the CN site, hydrogen has two imaginary frequencies, corresponding to a second-order transition state.

In conclusion, H atoms are stable in the pseudo-tetrahedral and octahedral sites. Despite the difficult stability of the TE site, the energy values agree well with those obtained in existing theoretical work: in the range of [1.40 – 1.67] eV for the most stable site, see the work of Ferry et al. [36]. Zero point energies (ZPE) have also been calculated for stable sites. Irrespective of the species introduced, the relative stability of the different sites remains unchanged. Only the data for hydrogen appear in the literature; the calculated ΔZPE is low (less than 100 meV)
<caption>**Table 3**</caption>

Volume of formation, $\Omega_f[X] = V[Be+X] - V[Be]$, and of relaxation, $\Omega_r$ in Å³, calculated using the Aneto code [22], Voronoï[38]/Bader [39] charges (B, in e⁻) of the X interstitial, and elastic dipole tensors, $P_{ij}$ in eV, and relaxation volume tensor, $V_{ij}$ in Å for the stable sites. The numbers in parentheses correspond to the number of electrons of each species used in the pseudo-potential (†). The difference in Voronoïvolume, $\delta V_v = V_v[Be+X] - V_v[Be]$ in $\forall^3$, of each stable site is also given. Values were computed on 4 × 4 × 3 super-cells.

<table><thead><tr><th>atom</th><th>site</th><th>B</th><th>Ω<sub>f</sub> [Å<sup>3</sup>]</th><th>Ω<sub>r</sub> [$\forall^3$]</th><th>δV<sub>v</sub> [$\forall^3$]</th><th>V [Å]</th><th>P [eV]</th></tr></thead><tbody><tr><td rowspan="3">H</td><td rowspan="3">OC</td><td rowspan="3">2.9/1.7 (1†)</td><td rowspan="3">3.4</td><td rowspan="3">3.2</td><td rowspan="3">0.4</td><td>[1.5 0.0 0.0]</td><td>[3.2 0.0 0.0]</td></tr><tr><td>[0.0 1.5 0.0]</td><td>[0.0 3.2 0.0]</td></tr><tr><td>[0.0 0.0 0.4]</td><td>[0.0 0.0 1.0]</td></tr><tr><td rowspan="3">C</td><td rowspan="3">TE</td><td rowspan="3">2.2/1.8 (1†)</td><td rowspan="3">4.3</td><td rowspan="3">4.2</td><td rowspan="3">0.9</td><td>[2.0 0.0 0.0]</td><td>[4.3 0.0 0.0]</td></tr><tr><td>[0.0 2.0 0.0]</td><td>[0.0 4.3 0.0]</td></tr><tr><td>[0.0 0.0 0.4]</td><td>[0.0 0.0 1.0]</td></tr><tr><td rowspan="3">C</td><td rowspan="3">OC</td><td rowspan="3">7.8/4.2 (4†)</td><td rowspan="3">13.5</td><td rowspan="3">13.1</td><td rowspan="3">0.7</td><td>[3.2 0.0 0.0]</td><td>[7.3 0.0 0.0]</td></tr><tr><td>[0.0 3.2 0.0]</td><td>[0.0 7.3 0.0]</td></tr><tr><td>[0.0 0.0 7.2]</td><td>[0.0 0.0 17.3]</td></tr><tr><td rowspan="3">BT</td><td rowspan="3">7.2/4.2 (4†)</td><td rowspan="3">15.6</td><td rowspan="3">15.1</td><td rowspan="3">1.7</td><td rowspan="3">1.7</td><td>[7.7 0.0 0.0]</td><td>[16.8 0.0 0.0]</td></tr><tr><td>[0.0 7.7 0.0]</td><td>[0.0 16.8 0.0]</td></tr><tr><td>[0.0 0.0 0.0]</td><td>[0.0 0.0 0.6]</td></tr><tr><td rowspan="3">N</td><td rowspan="3">OC</td><td rowspan="3">8.0/5.3 (5†)</td><td rowspan="3">12.1</td><td rowspan="3">11.8</td><td rowspan="3">0.5</td><td>[3.8 0.0 0.0]</td><td>[8.5 0.0 0.0]</td></tr><tr><td>[0.0 3.5 0.0]</td><td>[0.0 8.5 0.0]</td></tr><tr><td>[0.0 0.0 5.1]</td><td>[0.0 0.0 12.4]</td></tr><tr><td rowspan="3">BT</td><td rowspan="3">7.8/5.2 (5†)</td><td rowspan="3">14.5</td><td rowspan="3">14.3</td><td rowspan="3">1.2</td><td rowspan="3">1.2</td><td>[7.0 0.0 0.0]</td><td>[15.1 0.0 0.0]</td></tr><tr><td>[0.0 7.0 0.0]</td><td>[0.0 15.1 0.0]</td></tr><tr><td>[0.0 0.0 0.8]</td><td>[0.0 0.0 2.5]</td></tr><tr><td rowspan="3">CN</td><td rowspan="3">7.8/5.4 (5†)</td><td rowspan="3">12.1</td><td rowspan="3">11.9</td><td rowspan="3">1.0</td><td rowspan="3">1.0</td><td>[6.6 0.0 0.0]</td><td>[13.6 0.0 0.0]</td></tr><tr><td>[0.0 1.0 -5.8]</td><td>[0.0 3.3 -12.1]</td></tr><tr><td>[0.0 -5.8 4.6]</td><td>[0.0 -12.1 11.1]</td></tr><tr><td rowspan="3">O</td><td rowspan="3">OC</td><td rowspan="3">8.1/6.4 (6†)</td><td rowspan="3">12.5</td><td rowspan="3">12.3</td><td rowspan="3">0.9</td><td>[3.0 0.0 0.0]</td><td>[6.8 0.0 0.0]</td></tr><tr><td>[0.0 2.6 0.0]</td><td>[0.0 6.8 0.0]</td></tr><tr><td>[0.0 0.0 6.3]</td><td>[0.0 0.0 15.1]</td></tr><tr><td rowspan="3">BT</td><td rowspan="3">8.0/6.3 (6†)</td><td rowspan="3">12.8</td><td rowspan="3">12.7</td><td rowspan="3">1.4</td><td rowspan="3">1.4</td><td>[5.7 0.0 0.0]</td><td>[12.4 0.0 0.0]</td></tr><tr><td>[0.0 5.7 0.0]</td><td>[0.0 12.4 0.0]</td></tr><tr><td>[0.0 0.0 1.7]</td><td>[0.0 0.0 4.4]</td></tr></tbody></table>

compared to the other species considered here (oxygen, nitrogen and carbon).

In order to analyze the interactions between the interstitial species and the Be atoms, additional parameters were calculated: the relaxation volume, $\Omega_f$, calculated directly from the DFT simulations, i.e. the volume difference of the supercell (4 × 4 × 3) with and without interstitial, $\Omega_r$, the relaxation volume calculated using the elasticity theory and the Aneto code [22] (and elastic constants presented above); the difference between the Voronoïvolumes, $\delta V_v = V_v[Be+X] - V_v[Be]$ before and after relaxation. $\Omega_r$ and $\Omega_f$ quantify the global relaxation of the lattice when an atom is inserted into the bulk, and $\delta V_v$ quantifies the local relaxation around the interstitial. We also computed the relaxation volume tensor, $V$, and the elastic dipole, $P$. These data are summarized in Table 3.

First, we can see that the values of $\Omega_f$ and $\Omega_r$ are in excellent agreement. The elasticity theory reproduces the global volume relaxation due to insertion with high accuracy for all species. It can therefore be used here to capture the physics involved in this process. The relaxation volume tensor, $V_{ij}$, quantifies the anisotropy observed during insertion. Therefore, the deformation is mainly felt in the basal plane for almost all configurations, except for OC for O, C, and N atoms. In this case, the insertion induces a strong increase in the c axis. This strong variation of the parameter c can be related to the large difference between the relaxed and non-relaxed energies (the atomic forces were always minimized in the two calculations). It can also be noted that there is no clear correlation between relaxation volumes (global effects), Voronoïvariations (local effects) and site stability. The hypothesis that the smaller the deformation induced by the interstitial, the greater the stability is not confirmed in the case of beryllium. Therefore, “steric” effects cannot be the only explanation for the stability of the sites. The electronic interactions between the interstitials and the lattice are therefore essential to understand the stability of the sites.

Interaction energies with periodic images (calculated with Aneto), $E_d$, were also calculated, see Table 2, or can be deduced from the non-relaxed energies. The correction is generally small, ≈ 120 meV. In the

framework of the elasticity theory, the solid solution must be calculated from the unrelaxed value corrected by $E_d$. In the case of nitrogen, the electronic interactions are large compared to the elastic correction.

To quantify the electronic interactions, the Bader charges were extracted from the DFT calculations, either using the “classic” Bader charge partitioning [39], or based on the Voronoïvolumes proposed by Yu [38]. Depending on the approach, the analysis gives very different results: the charge transfer from the metal to the interstitial is either very strong with the “classic” Bader charge, or small with the Voronoïmethod. A detailed analysis of the results obtained with the classical approach revealed an inconsistency in the calculation. In fact, the Bader charge of Be atoms, with or without interstitials, varies greatly: from 3.5 to 5 electrons, which makes no physical sense. There is a numerical artifact which means that all Be atoms are no longer equivalent. The second Bader charge evaluation, based on Voronoïanalysis, shows little charge transfer between the beryllium lattice and the inserted atom, except for hydrogen (+0.7 electron in OC sites and +0.8 in TE sites). The charge maps are equivalent to those obtained by Zhang et al. [10] and do not add any significant information about the metal and the inserted species.

## 5. Diffusivity

The calculation of the diffusion coefficients required the identification of diffusion pathways and, where necessary, NEB calculations to determine the transition states and energy barriers. The transition state frequencies were systematically calculated.

Regardless of the diffusing species, in the case of the OC-OC straight jumps along the c axis, the atom necessarily passes through the BO saddle point, which has been identified as a first-order transition state, see Table 1. Similarly, in the basal plane, the interstitial species passes through the CN site. However, depending on the nature of the CN site (stable, first or second order transition state), this jump is taken into account or not in the calculation of the diffusion coefficients. All
<figure><img src="image_6.png" /><figcaption>Fig. 3. Diffusion paths (red sticks) of C (top left) and O (top right), H (bottom left) and N (bottom right) atoms in the Be system (green balls).</figcaption></figure>

<caption>Table 4</caption>

Migration energies ($E_{xy}^m$, in eV) of jumps from x to y. Where the transition state symmetry has been identified, the site name is given in brackets; the letters in brackets then correspond to the Wyckoff symmetry notation used (see text). Calculations were performed on 4 × 4 × 3 supercells. The jumps are shown in Fig. 3. To avoid any ambiguity, the jump direction is indicated in brackets if it lies in the basal plane (labeled ab), or along the c axis.

<table><tbody><tr><th colspan="3">Hydrogen</th></tr><tr><th>y \ x</th><th>OC (a)</th><th>TE (f)</th></tr><tr><th>OC (a)</th><td>0.597 (BO, c)</td><td>0.445</td></tr><tr><th>TE (f)</th><td>0.261</td><td>0.003 (BT)</td></tr><tr><th colspan="3">Carbon</th></tr><tr><th>y \ x</th><th>OC (a)</th><th>BT (d)</th></tr><tr><th>OC (a)</th><td>1.195 (BO, c)<br>0.360 (CN, ab)</td><td>0.006</td></tr><tr><th>BT (d)</th><td>0.810</td><td>-</td></tr><tr><th colspan="3">Nitrogen</th></tr><tr><th>y \ x</th><th>OC (a)</th><th>BT (d)</th><th>CN (g)</th></tr><tr><th>OC (a)</th><td>-</td><td>0.924</td><td>0.242</td></tr><tr><th>BT (d)</th><td>0.696</td><td>-</td><td>0.129</td></tr><tr><th>CN (g)</th><td>0.943</td><td>1.163</td><td>-</td></tr><tr><th colspan="3">Oxygen</th></tr><tr><th>y \ x</th><th>OC (a)</th><th>BT (d)</th></tr><tr><th>OC (a)</th><td>1.240 (BO, c)</td><td>1.565</td></tr><tr><th>BT (d)</th><td>0.166</td><td>-</td></tr></tbody></table>

identified diffusion paths for each species are shown in Fig. 3 and the migration energies are summarized in Table 4.

Regarding the migration path of carbon, several paths between OC and BT were considered. Migration along the c axis, corresponding to a jump between two nearest neighbor first OC sites (1NN), is represented by a straight path through the BO site with a barrier energy of 1.20 eV (as previously obtained by Zhang et al. [10]), a second through the CN site (in the basal plane), and a final path between an OC site and a BT site with a lower barrier energy of about 0.8 eV. As mentioned above, CN is a first-order transition state for carbon, so the diffusion process in the xy plane is also considered. The investigated paths from OC to BT result in a very low energy barrier to leave BT, about 10 meV.

This confirms the metastable nature of the BT site. Furthermore, the NEB simulation of the latter path showed that a straight jump between BT sites is not possible, the carbon preferentially going to the OC site first. Thus, from an energetic point of view, the migration landscape of carbon in metallic beryllium amounts to jumps between OC sites, see Fig. 3.

In the case of oxygen, migration between two 1NN BT sites is not considered because of the second-order behavior of the CN site. This statement does not take into account the analysis of Zhang et al. [10] for this migration pathway. To complete the energy landscape of oxygen migration, NEB calculations were performed to investigate the path from BT to OC. The migration energies obtained are close to those previously published [10], in fact the present data are higher by a few tenths of meV. However, contrary to what Zhang et al. found, our calculations identified the BO site as a first order saddle point. Consequently, the present work considers the migration path from the OC site to another along the c-axis.

As seen above, in the case of hydrogen, there are two sites: TE and OC sites. So, the displacement matrix contains three different jumps. There is a jump between two 1NN TE sites, where the transition state is located at the BT site. A NEB calculation was performed to determine the migration path of the TE-OC jump. Finally, for the jumps between two 1NN OC sites, the migration in the basal plane is excluded because the CN site is a second-order transition state. Nevertheless, migration along the c-axis is possible through the BO transition state. The theoretical energy landscape of interstitial hydrogen in metallic beryllium has been studied extensively in the past [10–13]. Unfortunately, all of these works considered the BT site as the ground state. The present work finds that this site is actually a saddle point and identifies the TE site as the most stable insertion site for hydrogen. However, due to the very small gap in terms of distance and energy between BT and TE sites, the interpretation does not change drastically, and our migration energies related to TE→OC→TE agree well with published data without challenging the energy landscape.

The last case not yet studied is nitrogen, for which three stable sites were found. From the frequency analysis presented in the previous section, as well as from the study of possible pathways, three pathways were identified: OC-BT, OC-CN, BT-CN, while the straight CN-CN path was excluded by the NEB simulations (very energetic and saddle point with two imaginary frequencies). The combination of stable sites and the high multiplicity of CN leads to a complex migration path for the interstitial nitrogen atom, see Fig. 3. The calculated migration energies indicate that the energy landscape in the vicinity of CN is flat compared to the other stable sites, meaning that the diffusing nitrogen atom jumps quickly through the OC or BT site.

From calculated frequencies, we can evaluate the jump rates, $\Gamma_{ab}$. $\Gamma_{ab}$ corresponds to the probability of leaving site a to site b, it is expressed by:

$$ \Gamma_{ab}(T) = \frac{k_B T}{h} \frac{Z_{TS}}{Z_a} \exp \left[ - \frac{\Delta E_m^{ab}}{k_B T} \right]. \quad (6) $$

where $\Delta E_m^{ab}$ is the migration energy to go from a to b, $Z_{TS}$ the (vibration) partition function for the transition state and $Z_a$ the partition function for the initial position, i.e. a. Z is related to the Helmholtz free energy by:

$$ F_v(T) = -k_B T \ln Z. \quad (7) $$

For numerical cost reasons, in the case of nitrogen (where transition states exhibit low or no symmetry), the ratio $Z_{TS}/Z_a$ was approximated with the high temperature approximation of the jump frequency, see [40].

Diffusion coefficients were expressed using a multi-state approach [41]. The jump rate notation includes labels of the initial and final sites corresponding to the symmetry notation (see Table 4). In the case
<figure><img src="image_7.png" /><figcaption>Fig. 4. Diffusion coefficients of H and C atoms in the Be system. In the case of hydrogen, experimental data are extracted from Jones [7], Abramov [8], Tazhibaeva [45] and Kizu [9].</figcaption></figure>

where OC and BT sites have to be considered, for oxygen and carbon, according to Wu [42], the diffusion coefficients are given by:

$$D_{x,y}[T] = \frac{a^2 \Gamma_{dd} \Gamma_{ad} + 2\Gamma_{da}\Gamma_{ad} + 3\Gamma_{da}\Gamma_{aa}^x}{2 \frac{\Gamma_{da} + \Gamma_{ad}}{\Gamma_{da} + \Gamma_{ad}}} \quad (8)$$

$$D_z[T] = \frac{c^2 6\Gamma_{dd}\Gamma_{ad} + 3\Gamma_{da}\Gamma_{ad} + 2\Gamma_{da}\Gamma_{aa}^z}{8 \frac{\Gamma_{da} + \Gamma_{ad}}{\Gamma_{da} + \Gamma_{ad}}} \quad (9)$$

In the case of carbon, the probability of diffusing via the site CN, $\Gamma_{dd}^x$ is zero. In the case of oxygen, $\Gamma_{aa}^x$ must be fixed at zero.

In the case where the OC and TE sites are stable, i.e. hydrogen, we can use the formula published by Connétable [43], by eliminating $\Gamma_{ii}$:

$$D_{x,y}[T] = a^2 \frac{\Gamma_{af} \Gamma_{fa}}{2 \frac{\Gamma_{af} + \Gamma_{fa}}{\Gamma_{af} + \Gamma_{fa}}} \quad (10)$$

$$D_z[T] = \frac{c^2 (3 \Gamma_{aa}^z \Gamma_{fa} + 2 \Gamma_{aa}^z \Gamma_{ff}^z + 3 \Gamma_{af} \Gamma_{ff}^z) \Gamma_{fa}}{4 \frac{(2 \Gamma_{af} + \Gamma_{fa})(3 \Gamma_{fa} + 2 \Gamma_{ff}^z)}{(2 \Gamma_{af} + \Gamma_{fa})(3 \Gamma_{fa} + 2 \Gamma_{ff}^z)}} \quad (11)$$

In the case of nitrogen, there are 3 stable sites -OC, BT and CN-, the diffusion coefficient is given by [44]:

$$D_{x,y}[T] = a^2 \frac{A}{B} \quad (12)$$

$$D_z[T] = \frac{c^2 C}{8 B} \quad (13)$$

where

$$A = \Gamma_{ad}(4\Gamma_{da}(\Gamma_{ga} + \Gamma_{gd}) + \Gamma_{dg}(4\Gamma_{ga} + \Gamma_{gd})) + \Gamma_{ag}(\Gamma_{da}(3\Gamma_{ga} + 4\Gamma_{gd}) + \Gamma_{dg}(3\Gamma_{ga} + \Gamma_{gd})) \quad (14)$$

$$B = \Gamma_{ad}(3\Gamma_{gd} + \Gamma_{dg} + \Gamma_{ga}) + \Gamma_{ag}(3\Gamma_{dg} + 3\Gamma_{da} + \Gamma_{gd}) + \Gamma_{da}(\Gamma_{ga} + \Gamma_{gd}) + \Gamma_{dg}\Gamma_{ga} \quad (15)$$

and

$$C = 3\Gamma_{ad}(\Gamma_{da}(\Gamma_{ga} + \Gamma_{gd}) + \Gamma_{dg}(\Gamma_{ga} + \Gamma_{gd})) + 3\Gamma_{ag}\Gamma_{gd}(\Gamma_{da} + \Gamma_{dg}) \quad (16)$$

The diffusion coefficients were therefore calculated from these expressions using the DFT data presented above. The curves are shown in Figs. 4 and 5. The coefficients were then fitted on Arrhenius curves. The parameters are given in Table 5.

It can be seen that the diffusivity is always isotropic, except for carbon which has a strong anisotropic behavior especially at low temperature. At low T (ambient temperature), H diffuses very rapidly in Be, faster than all the other interstitials, while at high temperature, carbon atoms also diffuse rapidly.

There are many experimental papers on hydrogen and its isotope diffusivity [7–9], only one data for N atoms reported by Pemsler [15] (a value of $5 \cdot 10^{-14}$ m²/s at 1025 °C), but none for C and O atoms.

<table><caption>Table 5
Activation energies, $E_a$ in eV, and $D_o$ in m²/s, obtained from an Arrhenius fit.</caption><thead><tr><th rowspan="2">Direction</th><th colspan="2">Hydrogen</th><th colspan="2">Carbon</th></tr><tr><th>E<sub>a</sub></th><th>D<sub>a</sub></th><th>E<sub>a</sub></th><th>D<sub>a</sub></th></tr></thead><tbody><tr><td>x, y</td><td>0.40</td><td>1.3 10<sup>-6</sup></td><td>0.35</td><td>2.2 10<sup>-6</sup></td></tr><tr><td>z</td><td>0.32</td><td>5.3 10<sup>-7</sup></td><td>0.56</td><td>4.1 10<sup>-6</sup></td></tr></tbody><thead><tr><th rowspan="2">Direction</th><th colspan="2">Nitrogen</th><th colspan="2">Oxygen</th></tr><tr><th>E<sub>a</sub></th><th>D<sub>a</sub></th><th>E<sub>a</sub></th><th>D<sub>a</sub></th></tr></thead><tbody><tr><td>x, y</td><td>1.09</td><td>1.2 10<sup>-6</sup></td><td>1.56</td><td>6.2 10<sup>-7</sup></td></tr><tr><td>z</td><td>1.08</td><td>1.2 10<sup>-6</sup></td><td>1.56</td><td>5.8 10<sup>-7</sup></td></tr></tbody></table>

In the case of hydrogen, the theoretical coefficient is found several orders of magnitude faster than the experimental data. This discrepancy can be explained by a very strong interaction of hydrogen and its isotopes with Be vacancies. In fact, Ferry et al. [13] showed a strong hydrogen segregation energy of more than –1.2 eV in a vacancy. Based on this work, the formation of vacancy-hydrogen clusters is also favored which would alter the H diffusion coefficient. In nickel, where the same vacancy-hydrogen interaction is weaker, about –0.3 eV [34,46], a potentially significant effect of the presence of vacancies on hydrogen diffusivity was found as a function of H concentration. The vacancies act as hydrogen traps, slowing down its diffusion. By applying this analysis to the Be system, we can explain the discrepancy found between theoretical and experimental results. To confirm our approach, we used the Oriani model [47], which takes into account the presence of traps in the expression of the diffusion coefficient of H by calculating an effective diffusion coefficient, $D_{eff}$, as [46,48]:

$$D_{\text{eff}}[H, T] = \frac{D_{\text{ideal}}[H, T]}{1 + \sum_{\text{trap}} C_{\text{trap}} \exp(-E_{\text{trap}}/k_{B}T)} \quad (17)$$

where the sum is performed on all the traps present in the system. $E_{trap}$ is then the trapping energy, and $C_{trap}$ the trap concentration in the system. With a value of –1.2 eV as the trapping energy and a trap concentration equal to the vacancy concentration calculated with its formation energy, $E_{vac}^f$, of 0.82 eV [13] ($C_{trap} \approx \exp(-E_{vac}^f/k_B T)$), the effective diffusion coefficient was computed, see Fig. 6. By taking trap effects into account, it is then possible to reconcile theoretical and experimental data.

## 6. Conclusion

In this paper, the insertion and diffusivity of H, C, N and O in the Be system were reexamined. In the case of C and O, results are similar to those found in the literature. These atoms are preferentially located in the OC and BT sites, respectively. For nitrogen, there are no data in the literature, neither theoretical nor experimental. Nitrogen is found to
<figure><img src="image_8.png" /><figcaption>Fig. 5. Diffusion coefficients of N and O atoms in the Be system. In the case of nitrogen, the experimental value is extracted from Pemsler [15].</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>Fig. 6. Diffusion coefficients of hydrogen in the Be system, with and without taking into account the Oriani model, see text. Experimental data are extracted from Jones [7], Abramov [8], Tazhibaeva [45] and Kizu [9].</figcaption></figure>

be preferentially stable in the *BT* sites like oxygen, but an additional stable site was found. In the case of hydrogen, the study carried out using different codes and different system sizes shows that the most stable site is the tetrahedral site, although the difference with the *BT* site is very small at 0 K. The results suggest that this site becomes stable as temperature rises.

The results also show that site stability strongly depends on the interstitial species. Additionally, based on the analysis of elastic effects, which accurately reproduces the effects of relaxation due to the atom insertion, this work highlights that the electronic interactions between the inserted species and beryllium are strong. In particular, a strong charge transfer is observed, which explains why the most stable sites are not necessarily those leading to a deformation of the minimal system.

After having studied the diffusion of the species considered, the values of the diffusion coefficients were calculated. In the case of H, we find that the diffusion coefficient is faster than that which is observed experimentally. These results are consistent with the fact that hydrogen interacts strongly with the intrinsic defects of the Be system, which act as traps for hydrogen. For nitrogen, our results are in relatively good agreement with the only study available in the literature. Finally, for C and O there are no experimental data, so we give Arrhenius parameters and expect an overestimation of the diffusion coefficients. Finally, the diffusion coefficient C is the only one that is strongly anisotropic.

### CRediT authorship contribution statement

**Damien Connétable:** Writing – review & editing, Writing – original draft, Visualization, Validation, Supervision, Methodology, Investigation, Formal analysis, Data curation, Conceptualization. **Aurélien Riot:** Writing – original draft, Validation, Resources, Data curation, Conceptualization. **François Virot:** Writing – review & editing, Writing – original draft, Validation, Methodology, Investigation, Funding acquisition, Data curation, Conceptualization.

### Declaration of competing interest

The authors declare the following financial interests/personal relationships which may be considered as potential competing interests: Virot reports financial support was provided by Euratom Research and Training Programme. If there are other authors, they declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

### Acknowledgments

This work was performed using HPC resources from CALMIP (Grant 2022 and 2023-p0749) and GENCI-TGCC (Grant A0080911483 and A0100911483). This work has been carried out within the framework of the Euro-fusion Consortium, funded by the European Union via the Euratom Research and Training Program (Grant Agreement No 101052200 - EURO-fusion). Views and opinions expressed are however those of the author(s) only and do not necessarily reflect those of the European Union or the European Commission. Neither the European Union nor the European Commission can be held responsible for them.

### Data availability

Data will be made available on request.

### References

[1] P.A. Lightsey, C.B. Atkinson, M.C. Clampin, L.D. Feinberg, James Webb Space Telescope: Large deployable cryogenic telescope in space, Opt. Eng., Bellingham 51 (1) (2012) 011003.

[2] N. Nagel, Beryllium and copper-beryllium alloys, ChemBioEng Reviews 5 (1) (2018) 30–33, http://dx.doi.org/10.1002/cben.201700016, URL https://doi.org/10.1002/cben.201700016.

[3] G. Federici, V. Barabash, R. Doerner, P. Lorenzetto, G. Matthews, A.R. Raffray, Beryllium as a plasma facing material for near-term fusion devices, in: Reference Module in Materials Science and Materials Engineering, 2012.

[4] L. Boccaccine, G. Aiello, J. Aubert, C. Nachmann, T. Barrett, A. Del Nevo, D. Demange, L. Forest, F. Hernandez, P. Norajitra, P. Porempovic, D. Rapisarda, P. Sardain, M. Utili, L. Vala, Objectives and status of eurofusion demo blanket studies, Fusion Eng. Des. 109 (2016) 1199.
[5] A. Goldberg, Beryllium manufacturing processes, 2006.

[6] F. Virot, L. Ferry, Y. Ferro, C. Pardanaud, M. Barrachin, Contribution to a better evaluation of the dust speciation in a case of an accident in iter, Fusion Eng. Des. (2017).

[7] P. Jones, R. Gibson, Hydrogen in beryllium, J. Nucl. Mater. 21 (3) (1967) 353-354, http://dx.doi.org/10.1016/0022-3115(67)90190-0, URL https://www.sciencedirect.com/science/article/pii/0022311567901900.

[8] E. Abramov, M. Riehm, D. Thompson, W. Smeltzer, Deuterium permeation and diffusion in high-purity beryllium, J. Nucl. Mater. 175 (1-2) (1990) 90-95.

[9] K. Kizu, K. Miyazaki, T. Tanabe, Hydrogen permeation and diffusion in beryllium, Fusion Technol. 28 (3P2) (1995) 1205-1210.

[10] P. Zhang, J. Zhao, B. Wen, Retention and diffusion of h, he, o, c impurities in be, J. Nucl. Mater. 423 (1) (2012) 164-169, http://dx.doi.org/10.1016/j.jnucmat.2012.01.027, URL https://www.sciencedirect.com/science/article/pii/S002231152000396.

[11] M.G. Ganchenkova, V.A. Borodin, R.M. Nieminen, Hydrogen in beryllium: Solubility, transport, and trapping, Phys. Rev. B 79 (2009) 134101, http://dx.doi.org/10.1103/PhysRevB.79.134101, URL https://link.aps.org/doi/10.1103/PhysRevB.79.134101.

[12] S. Middleburgh, R. Grimes, Defects and transport processes in beryllium, Acta Mater. 59 (18) (2011) 7095-7103, http://dx.doi.org/10.1016/j.actamat.2011.07.064, URL https://www.sciencedirect.com/science/article/pii/S1359645411005507.

[13] L. Ferry, F. Virot, Y. Ferro, D. Matveev, C. Linsmeier, M. Barrachin, Diffusivity of hydrogen and properties of point defects in beryllium investigated by DFT, J. Nucl. Mater. 524 (2019) 323-329, http://dx.doi.org/10.1016/j.jnucmat.2019.07.016, URL https://www.sciencedirect.com/science/article/pii/S0022311519301278.

[14] P. Zhang, J. Zhao, B. Wen, Vacancy trapping mechanism for multiple hydrogen and helium in be: A first principles study, J. Phys.: Condens. Matter. 24 (2012) 095004.

[15] J. Pemsler, R. Anderson, E. Rapperport, Solubility and Diffusion of Gases in Beryllium, Tech. rep., Textron Inc West Concord MA Nuclear Metals Div. 1963.

[16] G. Kresse, J. Hafner, Ab initio molecular dynamics for liquid metals, Phys. Rev. B 47 (1993) 558-561, http://dx.doi.org/10.1103/PhysRevB.47.558, URL https://link.aps.org/doi/10.1103/PhysRevB.47.558.

[17] J.P. Perdew, K. Burke, M. Ernzerhof, Generalized gradient approximation made simple, Phys. Rev. Lett. 77 (1996) 3865-3868, http://dx.doi.org/10.1103/PhysRevLett.77.3865, URL https://link.aps.org/doi/10.1103/PhysRevLett.77.3865.

[18] P. Hohenberg, W. Kohn, Inhomogeneous electron gas, Phys. Rev. 136 (1964) B864-B871, http://dx.doi.org/10.1103/PhysRev.136.B864.

[19] W. Kohn, L. Sham, Self-consistent equations including exchange and correlation effect, Phys. Rev. 140 (1965) 1133, http://dx.doi.org/10.1103/PhysRev.140.A1133.

[20] G. Kresse, D. Joubert, From ultrasoft pseudopotentials to the projector augmented-wave method, Phys. Rev. B 59 (1999) 1758-1775, http://dx.doi.org/10.1103/PhysRevB.59.1758, URL https://link.aps.org/doi/10.1103/PhysRevB.59.1758.

[21] H.J. Monkhorst, J.D. Pack, Special points for brillouin-zone integrations, Phys. Rev. B 13 (1976) 5188-5192, http://dx.doi.org/10.1103/PhysRevB.13.5188, URL https://link.aps.org/doi/10.1103/PhysRevB.13.5188.

[22] C. Varvenne, E. Clouet, Elastic dipoles of point defects from atomistic simulations, Phys. Rev. B 96 (2017) 224103, http://dx.doi.org/10.1103/PhysRevB.96.224103, URL https://link.aps.org/doi/10.1103/PhysRevB.96.224103.

[23] G. Henkelman, B.P. Uberuaga, H. Jónsson, A climbing image nudged elastic band method for finding saddle points and minimum energy paths, J. Chem. Phys. 113 (22) (2000) 9901-9904, http://dx.doi.org/10.1063/1.1329672, URL arXiv:https://doi.org/10.1063/1.1329672.

[24] A. Togo, F. Oba, I. Tanaka, First-principles calculations of the ferroelastic transition between rutile-type and cacl-type sio₂ at high pressures, Phys. Rev. B 78 (2008) 134106, http://dx.doi.org/10.1103/PhysRevB.78.134106, URL https://link.aps.org/doi/10.1103/PhysRevB.78.134106.

[25] C. Varvenne, F. Bruneval, M.-C. Marinica, E. Clouet, Point defect modeling in materials: Coupling ab initio and elasticity approaches, Phys. Rev. B 88 (2013) 134102, http://dx.doi.org/10.1103/PhysRevB.88.134102, URL https://link.aps.org/doi/10.1103/PhysRevB.88.134102.

[26] M.I. Aroyo, J.M. Perez-Mato, C. Capillas, Bilbao crystallographic server: I, databases and crystallographic computing programs, Z. Kristallogr. - Cryst. Mater. 221 (2009) 15-27, http://dx.doi.org/10.1524/zkri.2006.221.1.15.

[27] K. Momma, F. Izumi, VESTA3 for three-dimensional visualization of crystal, volumetric and morphology data, J. Appl. Crystallogr. 44 (6) (2011) 1272-1276, http://dx.doi.org/10.1107/S0021889811038970.

[28] K.J.H. Mackay, N.A. Hill, Lattice parameter and hardness measurements on high purity beryllium, J. Nucl. Mater. 8 (1963) 263-264, http://dx.doi.org/10.1016/0022-3115(63)90043-6, URL http://www.sciencedirect.com/science/article/pii/0022311563900436.

[29] G.V. Sin'ko, N.A. Smirnov, Relative stability and elastic properties of hcp, bcc, and fcc beryllium under pressure, Phys. Rev. B 71 (2005) 214108, http://dx.doi.org/10.1103/PhysRevB.71.214108, URL https://link.aps.org/doi/10.1103/PhysRevB.71.214108.

[30] A. Allouche, M. Oberkofler, M. Reinelt, C. Linsmeier, Quantum modeling of hydrogen retention in beryllium bulk and vacancies, J. Phys. Chem. C 114 (2010) 3588-3598, http://dx.doi.org/10.1021/jp910806j.

[31] M. Ganchenkova, P. Vladimirov, V. Borodin, Vacancies, interstitials and gas atoms in beryllium: Ab initio study, J. Nucl. Mater. 386-388 (2009) 79.

[32] M. David, A. Prillieux, D. Monceau, D. Connétable, First-principles study of the insertion and diffusion of interstitial atoms (h, c, n and o) in nickel, J. Alloys Compd. 822 (2020) 153555, http://dx.doi.org/10.1016/j.jallcom.2019.153555.

[33] R. Agarwal, D.R. Trinkle, Light-element diffusion in mg using first-principles calculations: Anisotropy and elastodiffusion, Phys. Rev. B 94 (2016) 054106, http://dx.doi.org/10.1103/PhysRevB.94.054106, URL https://link.aps.org/doi/10.1103/PhysRevB.94.054106.

[34] D. Connétable, M. David, Study of vacancy-(h, b, c, n, o) clusters in al using dft and statistical approaches: Consequences on solubility of solutes, J. Alloys Compd. 748 (2018) 12-25, http://dx.doi.org/10.1016/j.jallcom.2018.03.081, URL http://www.sciencedirect.com/science/article/pii/S0925838818309423.

[35] H. Krimmel, M. Fähnle, Properties of hydrogen isotopes in the tokamak plasma-facing material beryllium, J. Nucl. Mater. 231 (1-2) (1996) 159-161.

[36] L. Ferry, Contribution à l'étude du comportement du tritium dans le beryllium (Contexte ITER), (Ph.D. thesis), Aix-Marseille Université, 2017.

[37] P. Giannozzi, S. Baroni, N. Bonini, M. Calandra, R. Car, C. Cavazzoni, D. Ceresoli, G.L. Chiarotti, M. Cococcioni, I. Dabo, A.D. Corso, S. de Gironcoli, S. Fabris, G. Fratesi, R. Gebauer, U. Gerstmann, C. Gougoussis, A. Kokalj, M. Lazzeri, L. Martin-Samos, N. Marzari, F. Mauri, R. Mazzarello, S. Paolini, A. Pasquarello, L. Paulatto, C. Sbraccia, S. Scandolo, G. Sclauzero, A.P. Seitsonen, A. Smogunov, P. Umari, R.M. Wentzcovitch, Quantum espresso: A modular and open-source software project for quantum simulations of materials, J. Phys.: Condens. Matter. 21 (39) (2009) 395502.

[38] M. Yu, D.R. Trinkle, Accurate and efficient algorithm for bader charge integration, J. Chem. Phys. 134 (6) (2011) 064111, http://dx.doi.org/10.1063/1.3553716, URL arXiv:https://pubs.aip.org/aip/jcp/article-pdf/doi/10.1063/1.3553716/14735696/064111_1_online.pdf.

[39] G. Henkelman, A. Arnaldsson, H. Jónsson, A fast and robust algorithm for bader decomposition of charge density, Comput. Mater. Sci. 36 (3) (2006) 354-360, http://dx.doi.org/10.1016/j.commatsci.2005.04.010, URL http://www.sciencedirect.com/science/article/pii/S0927025605001849.

[40] E. Wimmer, W. Wolf, J. Sticht, P. Saxe, C.B. Geller, R. Najafabadi, G.A. Young, Temperature-dependent diffusion coefficients from ab initio computations: Hydrogen, deuterium, and tritium in nickel, Phys. Rev. B 77 (2008) 134305, http://dx.doi.org/10.1103/PhysRevB.77.134305, URL https://link.aps.org/doi/10.1103/PhysRevB.77.134305.

[41] U. Landman, M.F. Shlesinger, Stochastic theory of multistate diffusion in perfect and defective systems. I. Mathematical formalism, Phys. Rev. B 19 (1979) 6207-6219, http://dx.doi.org/10.1103/PhysRevB.19.6207, URL http://link.aps.org/doi/10.1103/PhysRevB.19.6207.

[42] H.H. Wu, P. Wisesa, D.R. Trinkle, Oxygen diffusion in hcp metals from first principles, Phys. Rev. B 94 (2016) 014307, http://dx.doi.org/10.1103/PhysRevB.94.014307, URL https://link.aps.org/doi/10.1103/PhysRevB.94.014307.

[43] D. Connétable, Theoretical study on hydrogen solubility and diffusivity in the γ-tial $L_1$ structure, Int. J. Hydrog. Energy 44 (2019) 12215-12227, http://dx.doi.org/10.1016/j.ijhydene.2019.03.110.

[44] H.H. Wu, D.R. Trinkle, Direct diffusion through interpenetrating networks: Oxygen in titanium, Phys. Rev. Lett. 107 (2011) 045504, http://dx.doi.org/10.1103/PhysRevLett.107.045504, URL http://link.aps.org/doi/10.1103/PhysRevLett.107.045504.

[45] I. Tazhibaeva, V. Shestakov, E. Chikhray, O. Romanenko, A.K. Klepikov, Deuterium permeation through beryllium with surface element composition control, 1994.

[46] Y. Wang, D. Connétable, D. Tanguy, Influence of trap connectivity on h diffusion: Vacancy trapping, Acta Mater. 103 (2016) 334-340, http://dx.doi.org/10.1016/j.actamat.2015.10.018, URL https://www.sciencedirect.com/science/article/pii/S1359645415300185.

[47] R. Oriani, The diffusion and trapping of hydrogen in steel, Acta Metall. 18 (1) (1970) 147-157, http://dx.doi.org/10.1016/0001-6160(70)90078-7, URL https://www.sciencedirect.com/science/article/pii/0001616070900787.

[48] R. Kirchheim, Monte-carlo simulations of interstitial diffusion and trapping—I. One type of traps and dislocations, Acta Metall. 35 (2) (1987) 271-280, http://dx.doi.org/10.1016/0001-6160(87)90235-5, URL https://www.sciencedirect.com/science/article/pii/0001616087902355.