The effect of hydrogen on the evolution of intergranular cracking: a cross-scale study using first-principles and cohesive finite element methods

Abstract

A combination of first-principles and cohesive finite element (CFE) cross-scale calculations is performed to examine the effect of hydrogen on the intergranular cracking of aluminum. First-principles calculations based on density functional theory are investigated to assess the adsorption and diffusion of H along the Al (111) surface and the cohesive energies of the grain boundaries (GBs) with different H concentrations. CFE calculations are used to simulate the evolution of intergranular cracking induced by hydrogen segregation in GBs. To combine first-principles calculations and CFE modelling, the GB cohesive energies are used as inputs for the fracture energies in the CFE calculations. The results show that H atoms diffuse into the interstitial sites in the bulk material and segregated to the GBs by overcoming the energy barriers. The cohesive energies of the GBs decrease linearly with increasing of H concentration. The application of these cross-scale approaches is very efficient for investigating the evolution of hydrogen induced intergranular cracking.

---

1. Introduction

Intergranular cracking is the primary failure mode for stress corrosion cracking (SCC) of high-strength Al alloys under atmospheric conditions. It has been widely suggested that hydrogen embrittlement (HE) induced by the diffusion of hydrogen increases a material's susceptibility to SCC.^1,2 Sample preparation, processing, and the service environment can result in residual hydrogen in materials. F. Y. Guo et al.³ studied the dissociation adsorption of a H₂O molecule on an Al (111) surface using first-principle calculations based on density functional theory (DFT). The method of transition states search (TSS) was adopted to investigate the energy barriers of H₂O dissociation on Al (111) surfaces. The results show that the hydrogen atom dissociating from H₂O needs to overcome a 248.32 kJ mol⁻¹ of energy barrier on clean Al (111) surface, while the dissociating energy barrier decreases to 128.53 kJ mol⁻¹ on the O pre-adsorbed Al (111) surface with the aid of O atom. H and OH dissociated from H₂O can be adsorbed on surface by the form of chemical bond separately.

The mechanisms for SCC in metals can be classified as either (i) anodic or (ii) cathodic, e.g. hydrogen evolution, adsorption, diffusion, or embrittlement. Stress-assisted localized dissolution was initially considered the mechanism for intergranular stress corrosion cracking, and this perspective persisted until the effect of hydrogen on the fracture process was identified.⁴ When hydrogen is adsorbed from the external environment onto the material, it preferentially diffuses along grain boundaries (GBs), eventually causing embrittlement. H atoms can diffuse to depths of several millimetres toward regions of high stress, become trapped at gain boundaries, and induce intergranular fractures without any contribution from the corrosion process.⁵

The behaviour of HE leads to a degradation of mechanical properties^6–8 such as fracture resistance, elongation to failure, and transition from ductility to brittleness, all of which are induced by interactions with hydrogen. The three fundamental mechanisms for HE, namely, hydride-induced embrittlement,¹¹ hydrogen-enhanced decohesion,¹² and hydrogen-enhanced localized plasticity^12,13 have been well documented.^9–13 Experimental and theoretical efforts^14–17 using conventional mechanical and in situ tests have been conducted to understand the mechanism for HE. One approach involves measuring the differences in mechanical properties between a hydrogen-charged sample and an uncharged hydrogen sample.^18,19 However, providing a control to guarantee the reproducibility of the results is difficult for such tests. In situ testing techniques such as environmental transmission electron microscopy (TEM) can address the dislocation nucleation, multiplication, and movement, all of which can significantly influence the mechanical properties of the metals. Most metals such as Fe, Al, Ni, and their alloys have been subjected to in situ TEM tests to investigate the effect of hydrogen gas introduction to the TEM chamber on the mobility of dislocations observed in the experiments.¹³ However, in situ tests^20,21 have stringent requirements for experimental environments and conditions that limit their range of applications.

Because of the difficulties of experimental testing for HE, calculation and simulation methods have been proposed. Finite element (FE) method and damage mechanics model have been applied to SCC.^22,23 It is also necessary to investigate the role of GBs damage in intergranular cracking behavior. The cohesive element approach can be used to explicitly account for the damage to GBs and has been reported in a number of cases in the literatures.^24–28 Igor Simonovski and Leon Cizelj²⁸ used a cohesive element approach to model intergranular cracking at the grain level in which cohesive elements with zero physical thickness were directly inserted between the adjacent grains. The fracture energy is one of the most parameters for cohesive finite elements (CFE) calculations, which usually estimated by some approximate algorithms.^27,28 The process for the calculation of fracture energy is complex with an obvious fluctuation, while the GB fracture energies can be defined as the cohesive energies which can be deduced by surface energies and GB energies. Many experimental and theoretical^29–33 efforts have been made to determine the surface and GB energies. Vasily V. Bulatov et al.²⁹ presented and justified a concise hypothesis on the topography of the functional space of interface energies and constructed a closed-form function that quantitatively describes energy variations in the 5-space of macroscopic parameters defining grain boundary geometry. This function is universal for the crystallography class of face-centered cubic metals. L. Vitos et al.³⁰ have used DFT to establish a database of surface energies for low index surfaces of 60 metals in the periodic table. The results show that the surface energy of Al (111), Al (100) and Al (110) is 1.199, 1.347 and 1.271 J m⁻². The work of W. R. Tyson et al.³³ shows that the surface energy of Al (111) is 1.143 J m⁻².

This paper presents a novel method that combines first-principles and cohesive finite element (CFE) calculations to avoid the inaccuracy of conventional mechanical tests and the complexity of in situ tests. This method can be used in conjunction with experimental tests. The primary objective of this work is to understand the effect of hydrogen on the evolution on intergranular crack initiation and propagation in a cross-scale process from an atomic to a macroscopic scale. We initially report a first-principle study on the adsorption and diffusion of H atoms fracture energy is obtained by first-principle calculation through grain boundary energy and surface energy. Finally, a 3D finite element model of a polycrystalline material containing GBs with cohesive elements is then presented to investigate the effect of H on the intergranular cracking of a high-strength aluminium alloy.

2. Computational models

Fig. 1 shows a schematic of the entire cross-scale process. H atoms are adsorbed onto the surface from the outer environment, diffuse into the bulk and segregate to the GBs to decrease the cohesive energy, which can be obtained by first-principles calculations. The weakened GBs easily crack when stressed. The evolution of intergranular cracking is studied using CFE calculations.

Fig. 1 The schematic of the entire cross-scale process for the effect of H on intergranular cracking.

2.1 The first-principles calculation

All calculations of atomic models presented in this section are performed using the MedeA-VASP 5.3 program,^34,35 a fast and highly reliable electronic structure method based on density functional theory (DFT).³⁶ The calculations are conducted with a plane-wave basis via the projector-augmented wave method.³⁷ The GGA-PBE exchange–correlation function³⁸ is used to describe the interaction. The convergence for the electronic iterations is 10⁻⁵ eV and is attained by using a normal (blocked Davidson) algorithm. Periodical boundary conditions are used, creating an infinite, periodically repeating system. A spin-polarized magnetic calculation is performed to account for the scalar spin-up and spin-down magnetic moments for each atom. The post-processing of the results involved structural and charge density constructs using VESTA.³⁹

2.1.1 Adsorption and transition states search calculations. The calculations for the adsorption and transition states search of the H atoms are conducted on 6-layer slabs of an Al (111) surface with a 15 Å vacuum gap. A (3 × 3), (2 × 3) and (2 × 2) supercell is used for the adsorption calculations. Free movement of the adsorbates is permitted on the four upper layers, and the bottom two layers are fixed. The slab models are calculated based on a 4 × 4 × 1, 6 × 4 × 1, and 6 × 6 × 1 Monkhorst–Pack grids⁴⁰ for the (3 × 3), (2 × 3) and (2 × 2) supercell, respectively. For the TSS calculations of the H atoms from surface to bulk, a (3 × 3) supercell is adopted, and a nudged elastic band is used to map the minimum energy path between the initial and final systems with three intermediate images and a spring constant of 5 eV Å⁻². The initial images are created using linear interpolation. Transition states (TS) are acquired solely for the highest saddle point and optimized to the best extent. The image closest to a saddle point is allowed to ascend to the saddle point if the greatest force on an atom was smaller than 0.5 eV Å⁻¹. The “reaction co-ordinate” is defined as the interval coefficient for linear interpolation. The cutoff energy for the adsorption and TSS calculations is 300 eV. The robust and efficient Bader scheme is used to divide the slab model into atomic regions.⁴¹ The MedeA-VASP includes a property module; this approach is based on a classic algorithm proposed by Bader and implemented by Henkelman et al.⁴²

2.1.2 Cohesive energy calculations. The GB strength can be measured in terms of its cohesive energy. The cohesive energy is defined as the reversible free energy change for the formation free surfaces from the GB. Approximating the free energy as the heat of formation from VASP, the cohesive energy can be calculated based on the following equation:

E^GB_coh = γ¹_surf + γ²_surf − γ_GB (1)

where γ¹_surf and γ²_surf are the energies of the two surfaces that form the GB and γ_GB is the GB energy.

The surface energy can be obtained from:

γ_surf = (E_surf − nE^bulk_{Al atom})/2A (2)

where E_surf is the total energy of the surface, E^bulk_{Al atom} is the energy of a single Al atom in the bulk material, n is the number of atoms of the whole slab surface model, and A is the GB area.

The GB energy can be calculated as follows:

γ_GB = (E_GB − n₁E^bulk_{Al atom} − n₂E^bulk_{Al atom})/2A (3)

where E_GB is the total energy of the GB, and n₁ and n₂ are the number of atoms on the two surfaces.

We assume that H atoms bond with each other to form H₂ molecules that escape from GBs through cracks when the intergranular cracks are produced and have no effect on the fresh surfaces. The equation for the hydrogen induced decohesion for GBs is

where m is the number of H atoms, E_H2 is the energy of the H₂ molecule, and γ^H_GB is the energy of the GB with H atoms, which is:

γ^H_GB = γ_GB + (E^H_GB − E_GB)/A (5)

where E^H_GB is the total energy of the GB with H atoms.

Substituting eqn (2)–(5) into eqn (1), the cohesive energy of the GB with H atoms is:

E^H_coh = E¹_surf + E²_surf + mE_H2 + E_GB − 2E^H_GB/2A (6)

Each term can be calculated directly using VASP software.

More than 10 types of atomic models of GBs are simulated to calculate the cohesive energies based on their surface and GBs energies. Different concentrations of H atoms are introduced into GBs with different twist angles and misorientations to quantify the decrease in GB cohesive energy according to eqn (6).

2.2 The CFE calculation

2.2.1 The CFE model. A 3D finite element polycrystalline microstructure model is generated using a 3D Voronoi tessellation based on Neper,⁴³ as shown in Fig. 2. Cohesive elements with zero physical thickness are used for the GBs. The model contained 300 grains and 1387 grain boundaries. Linear tetrahedron elements (ABAQUS type C3D4) are used for meshing the grains and linear triangular elements (ABAQUS type COH3D6) are used for meshing the grain boundaries. The crystallographic orientation of the grains is orthotropic and randomly distributed.

Fig. 2 The 3D polycrystalline model with 300 grains generated by Neper.

2.2.2 Loads and boundary conditions. A displacement of 9 × 10⁻⁴ mm is applied on the right surface along the x direction. The nodes of the left surface are constrained in the x, y, and z directions.

2.2.3 Material parameters. The following anisotropic elastic properties for a 2024 Al alloy are used for the CFE model:⁴⁴

E₁ = 92.1 GPa, E₂ = E₃ = 86.7 GPa, μ₁₂ = 0.34, μ₁₃ = 0.32, μ₂₃ = 0.33,

where E_i (i, j = 1, 2, 3) are the Young's moduli and μ_ij (i, j = 1, 2, 3) are the Poisson ratios.

The remaining Poisson's ratios can be calculated using:

For an anisotropic elasticity material, the relationship between stress and strain is

where σ_i and ε_i are normal stress and strain, τ_ij and γ_ij are the shear stress and shear strain, and [C_ij] is the stiffness matrix calculated based on the principle of mechanics of materials. The results are listed in Table 1.

Table 1 Elastic constants of CFE calculations for grains

Constant	C 11	C 22	C 33	C 12	C 13	C 23	C 44	C 55	C 66
Value (Gpa)	91.9	86.52	86.52	43.89	43.90	41.63	41.63	43.9	43.89

---

A cohesive zone approach^45–47 with damage initiation and evolution as implemented in ABAQUS is used as a constitutive law for the GBs. The 3D cohesive elements are embedded into the grain boundaries to represent all possible fracture paths. The constitutive response of the cohesive elements is based on the traction–separation law that assumes initial linear elastic behaviour followed by the initiation and evolution of damage (see Fig. 3).

Fig. 3 Example of traction–separation response (not to scale).

The bilinear traction–separation laws represent the relationship between traction T and crack opening displacement δ, where the indices n, s, and t denote the normal and two orthogonal shear directions of the cohesive elements, respectively. T⁰_n (T⁰_s, T⁰_t) refers to the maximum traction that the cohesive element can bear before failure, δ₀ is the critical separation at which the effective traction T reaches the maximum traction T₀, δ_f is the separation at failure (assumed to be 20δ₀) and the area under the curve G_c is the fracture energy. The tractions of a cohesive element T_n, T_s, and T_t are given by the following equation:

where K_nn, K_ss, and K_tt are the stiffnesses in the normal and two shear directions, respectively, and K_ss = K_tt = K_nn/2(1 + μ). ε_n, ε_s, and ε_t are defined as δ_n/t₀, δ_s/t₀, and δ_t/t₀, respectively, where t₀ represents the constitutive thickness of a cohesive element. t₀ is mostly different from the geometric thickness, which is typically close or equal to zero. Based on previous calculation,²⁸ a value of t₀ = 0.001 mm is used in this study and we assume that the specific deformation at the damage initialization point of ε⁰_n is assumed to be 0.001. This results in δ⁰_n = ε⁰_nt₀ = 10⁻⁶ mm and δ^f_n = 20δ⁰_n = 2 × 10⁻⁵ mm, which are reasonable values based on estimations of the yield stress and fracture energy of a 2024 Al alloy. The fracture energy G_c is equal to the cohesive energy, E_coh, calculated using first principles. Note that there are only three independent values of G_c, T₀, and δ because G_c = (T₀ × δ_f)/2 and K = T₀/δ₀. Viscous regularization is applied to improve the convergence and a value of μ = 0.001 is used.²⁸

3. Results and discussion

3.1 The adsorption and dissolution of H atoms on Al (111) surface

The surface energy and work function of a clean 6-layer-Al (111) surface are 1.06 J m⁻² and 4.12 eV, respectively, which agree well with the corresponding experimental values of 1.14 J m⁻² (ref. 33) and 4.24 eV.⁴⁸Fig. 4 shows the four highly symmetric adsorption sites (top, bridge, fcc, and hcp) that are considered. The adsorption and solution energies are calculated according to the following expression:

E_ads/sol = E_{Al (111)} + nE_H − E_{Al (111)–nH} (10)

where E_{Al (111)–nH} is the total energy of the adsorbed structure of H and Al (111) surface, E_{Al (111)} is the total energy of the clean Al (111), n denotes the number of the H atoms, and E_H is the energy of a single H atom. Thus, a higher energy indicates that interaction between the H atoms and the surface is more stable.

Fig. 4 The optimized six-layer Al (111) surface: (a) top view, four highly symmetric adsorption sites are considered: the top site, bridge site, fcc site, and hcp site. L and S indicate the long axis and short axis of repeating unit; (b) side view of Al atoms in different layers.

The adsorption structures for the four symmetric adsorption sites are calculated for the adsorption of a single H atom onto the Al (111) surface with a H atom coverage of 1/9 ML. A single H atom can only be adsorbed onto the fcc site and top sites, and if the initial site of the H atom is an hcp or bridge, it will spontaneously migrate to fcc site. The geometric parameters, adsorption energies and charge transfer are shown in Table 2. For the fcc-site adsorption, the average distance between the Al and H atoms is 1.92 Å, which is larger than the distance for top-site adsorption, 1.62 Å. A shorter distance between Al atom and H atom implies a stronger interaction, however, it is not reasonable to determine the adsorption strength only based on the distance because the H atom can bond with one Al atom on the surface and three Al atoms at fcc site. The fcc-site adsorption energy is 2.03 eV is larger than the top-site adsorption energy of 1.85 eV indicating that fcc sites are more stable adsorption site than top site on Al (111) surface at the coverage of 1/9 ML.

Table 2 Geometric parameters, solution energy and charge transfer of single H atom along Al (111) surface at an H coverage of 1/9 ML

Stable structures	Average distance (Å) Al–H	Solution energy (eV)	Charge transfer Al → H
H(fcc)	1.92	2.03	0.91
H(top)	1.62	1.85	0.61
H(TIS1-1)	1.87	1.52	1.20
H(TIS1-2)	1.91	1.62	1.06
H(OIS1)	2.07	1.54	1.29
H(TIS2)	1.87	1.62	1.12
H(OIS2)	2.05	1.54	1.32

---

H atoms can diffuse into interstitial sites in the bulk through the surface. There are five types of subsurface and bulk interstitial sites, including tetrahedral interstice sites (TIS) and octahedral interstice sites (OIS), as shown in Fig. 5. The geometric parameters, solution energies, and charge transfer are listed in Table 2. Compared with the adsorption of a single H atom on the surface, the dissolutions in bulk occur with a lower energy and higher charge transfer, as illustrated in Fig. 6. The adsorption energy increases with the charge transfer when H atom is on the surface, while the solution energy exhibits an inverse relationship with the charge transfer when H atom is in the interstitial sites. It can be inferred that the relationships are primarily due to the change in number of neighbouring Al atoms around the H atom.

Fig. 5 The interstitial sites of subsurface and bulk of Al (111) surface: (a) TIS_1-1, the first type of tetrahedral interstice site of subsurface; (b) TIS_1-2, the second type of tetrahedral interstice site of subsurface; (c) OIS₁, the octahedral interstice site of subsurface; (d) TIS₂, the tetrahedral interstice site of bulk, (e) OIS₂, the octahedral interstice site of bulk.

Fig. 6 The adsorption energy and solution energy of single H atom along Al (111) surface and corresponding charge transfer with different sites at an H coverage of 1/9 ML.

With the coverage of H atoms increases to 2/9 ML, five stable adsorption structures are achieved, as shown in Table 3. Two H atoms are adsorbed on a six-layer (3 × 3) Al (111) surface. H atoms can still only be adsorbed onto fcc site and top sites. Due to the differences in adsorption energy, the site stability ranking is as follow:

H(fcc) − H(fcc) > H(fcc) − H(top) − S > H(fcc) − H(top) − L > H(top) − H(top) − L > H(top) − H(top) − S.

Table 3 Stable structures of H atoms adsorption on Al (111) surface at an H coverage of 2/9 ML

Stable structures	Top view	Side view
H(fcc)–H(fcc)
H(fcc)–H(top)–S
H(fcc)–H(top)–L
H(top)–H(top)–S
H(top)–H(top)–L

---

The ranking further illustrates that the fcc site is a more stable site than the top site. The geometric parameters, adsorption energies and charge transfer are shown in Table 4. The adsorption energies are consistent with the sum of the adsorption of single H atom and the charge transfer remains almost the same as for single fcc site or top site adsorption compared with the previous calculation listed in Table 2. To investigate the effect of H coverage on adsorption behaviour, a smaller surface (2 × 3) supercell is constructed with an H coverage of 1/3 ML. The stable adsorption structures are shown in Table 5. Compared with the geometric parameters, adsorption energies and charge transfer for an H coverage of 2/9 ML, the adsorption structures for an H coverage of 1/3 ML are similar (Table 4), except for the formation of H₂. Comparing the H coverage of 2/9 ML with that of 1/3 ML, the adsorption energies of two H atoms on Al (111) surface basically remain the same for the same adsorption behaviour, which is equal to the sum of single H atom adsorption on the same site. The reason is that the distance between H–H (see Table 4) is much larger than the sum of radius of H atom (∼1.6 Å) so that there is no interaction between H atoms. We infer that the coverage of H has no obvious effect on the adsorption behaviour, when H adsorb on Al (111) surface by the form of atom. H₂ produces when the H–H distance is small enough for H₂ to form at the initial structure H(fcc)–H(bridge) or H(hcp)–H(bridge) at the H coverage of 1/3 ML. The phenomenon is not observed at an H coverage of 2/9 ML indicating that the formation of H₂ is related with to the coverage of H atoms. We infer that the formation of H₂ will become easier as the H coverage increases. The generation of H₂ can cause hydrogen blistering theory,⁴⁹ which is an important factor for HE.

Table 4 Geometric parameters, adsorption energies and charge transfer of H atoms on Al (111) surface at an H coverage of 2/9 and 1/3 ML

Coverage (ML)	Stable structures	Average distance (Å)			Adsorption energy (eV)	Charge transfer
		Al–H(fcc)	Al–H(top)	H–H		H(fcc)	H(top)
2/9	H(fcc)–H(fcc)	1.92	—	2.71	4.07	0.90	—
	H(fcc)–H(top)–S	1.90	1.62	2.38	3.90	0.90	0.62
	H(fcc)–H(top)–L	1.93	1.62	3.45	3.89	0.90	0.60
	H(top)–H(top)–S	—	1.62	2.85	3.70	0.58	0.58
	H(top)–H(top)–L	—	1.62	4.94	3.77	0.60	0.60
1/3	H(fcc)–H(fcc)	1.93	—	2.68	4.20	0.90	—
	H(fcc)–H(top)–S	1.90	1.62	2.34	3.93	0.90	0.61
	H(fcc)–H(top)–L	1.93	1.63	3.35	3.92	0.88	0.62
	H(top)–H(top)–S	—	1.62	2.85	3.70	—	0.57
	H(top)–H(top)–L	—	1.62	4.91	3.77	—	0.59
	H2 formation	—	—	0.75	4.52	—	—

---

Table 5 Stable structures of H atoms adsorption on Al (111) surface at an H coverage of 1/3 ML

Stable structures	Top view	Side view
H(fcc)–H(fcc)
H(fcc)–H(top)–S
H(fcc)–H(top)–L
H(top)–H(top)–S
H(top)–H(top)–L
H2 formation

---

3.2 Diffusion pathways for hydrogen from surface to bulk

The diffusion properties of H atoms provides a good reference for H atoms segregation at GBs. The diffusion pathways of H atom are investigated using a TSS method at an H coverage of 1/9 ML to gain an insight into the transition from surface to bulk diffusion. There are multiple possible pathways for diffusion from fcc, the stable adsorption site, to bulk OIS₂, the stable solution site. As shown in Fig. 7a, the H atom jumps from the top site to the fcc site on the surface with a diffusion energy barrier of 7.02 kJ mol⁻¹, which is higher than that of H diffusion along the top–bridge–fcc pathway (3.66 kJ mol⁻¹) or the top–hcp–bridge pathway (3.71 kJ mol⁻¹). In contrast, the energy barriers of the fcc–top, fcc–bridge and fcc–hcp pathways, which are 17.79, 8.24 and 13.73 kJ mol⁻¹ respectively, are higher. This indicates that the H atom at the top site prefers to jump through the bridge or hcp site to the fcc site. However, diffusion from fcc site to other sites on the surface is difficult.

Fig. 7 The diffusion pathways of H atom from Al (111) surface to bulk: (a) pathways and energy barrier of all calculated pathways; (b) the TSS result of optimal pathway of H atom from surface fcc site to bulk TIS₂ site.

For the diffusion from surface to subsurface, we considered the diffusion path from the bridge, fcc and hcp sites on the surface to the TIS_1-1, TIS_1-2 and OIS₁ sites on the subsurface. The diffusion energy barriers for the H atom jumping bridge to TIS_1-1, bridge to OIS₁, fcc to TIS_1-1, fcc to OIS₁, hcp to TIS_1-1 and hcp to OIS₁ are 47.98, 43.28, 49.88 and 54.15, 40.15 and 50.22 kJ mol⁻¹, respectively. Thus, it is much more difficult for the H atom to jump from surface to subsurface than to different sites on the surface. We infer that this higher energy barrier is due to the change of the number of neighboring Al atoms around H atom. The migrations of the H atom within the subsurface are also investigated. The energy barriers between the TISs, TIS_1-1 and TIS_1-2, are 14.58 kJ mol⁻¹ and 24.28 kJ mol⁻¹, which indicates that the two TISs are not equally site influenced by the surface. The energy barriers of TIS_1-1–OIS₁ and TIS_1-2–OIS₁ are 7.54 kJ mol⁻¹ and 13.06 kJ mol⁻¹, respectively.

For the diffusion of the H atom from subsurface to bulk, four pathways, namely, TIS_1-1–TIS₂, TIS_1-2–OIS₂, OIS₁–TIS₂ and OIS₁–OIS₂, were considered. The energy barriers for the four pathways are 15.73, 24.41, 17.44 and 16.58 kJ mol⁻¹, respectively. For the migration of H the atom within the bulk material, TIS₂–OIS₂, the energy barrier is 16.75 kJ mol⁻¹.

In conclusion, the optimal diffusion pathway of H from the surface to the bulk is fcc–hcp–TIS_1-1–TIS₂, as shown in Fig. 7b. The H atom at the fcc site jumps to the subsurface TIS_1-1 site through the hcp site and finally reaches the TIS₂ bulk site. The transition state (TS) of the fcc–hcp path occurs at the hcp site and the energy barrier is 13.73 kJ mol⁻¹. The TS for the hcp–TIS_1-1 pathway occurs at the position close to the TIS_1-1 site with a much higher energy battier of 40.15 kJ mol⁻¹. The diffusion energy barrier for the TIS_1-1–TIS₂ path decreases to 15.73 kJ mol⁻¹ and the TS occurs at the 3/4 site between the TIS_1-1 and the TIS₂ sites.

3.3 The effect of hydrogen on intergranular cracking of Al alloys

H atoms prefer to diffuse from an outer environment to the bulk along GBs⁵⁰ where the energy barrier for diffusion is much lower than that at a complete surface. The ability to resist intergranular cracking is related to the strength of the GBs, which can be quantified using the cohesive energy. The cohesive energy is determined by the surface energy and GB energy. The surface energies of the (111), (100) and (110) planes are 1.06 J m⁻², 1.25 J m⁻² and 1.28 J m⁻², respectively, and are consistent with previous calculation and experimental results.⁵¹ The energies of the more than 10 types of GBs generated by [110], [100] and [111] twist angles with different misorientations are shown in Table 6 and are consistent with previous work.^52,53 The energies of the GB with a [110] twist angle indicate that there is no obvious fluctuation with the increasing of misorientation but that there is a greater energy gap between different twist angles. The cohesive energy of the GBs have been calculated using eqn (1)–(3) with the first-principle method. The cohesive energy varies little for the GBs with the same twist angle because the same surface energy and similar GBs energies. The cohesive energies of GB [001]/22.6°, [011]/38.9°, and [111]/60° with different H concentrations have been calculated using eqn (4)–(6). The cohesive energy decreases linearly with the increasing H atoms concentrations (see Fig. 8), which indicates a hydrogen induced decohesion of HE mechanism. Comparing the slope of the lines, the effect of H on the decohesion for the three GBs is different. The cohesive energies of GBs with the same H concentrations are determined by their initial values and its descending rate.

Table 6 GBs energies of Al with different misorientations of [110], [111] and [100] twist angle

Twist axis	Misorientation °	γ GB J m^−2	γ surf J m^−2	E coh J m^−2
[110]	9.45	0.69	1.28	1.87
[110]	10.02	0.65		1.91
[110]	13.7	0.66		1.90
[110]	14.42	0.59		1.97
[110]	19.47	0.66		1.90
[110]	22	0.67		1.89
[110]	29.5	0.69		1.87
[110]	31.59	0.60		1.96
[110]	35.26	0.63		1.93
[110]	44.71	0.63		1.93
[111]	21.79	0.25	1.06	1.87
[111]	30	0.18		1.94
[100]	22.6	0.35	1.25	2.15
[100]	28.07	0.41		2.09
[100]	36.87	0.25		2.25

---

Fig. 8 The cohesive energies of GBs [001]/22.6°, [011]/38.9°, and [111]/60° as a function of H concentrations.

3.4 The CFE simulation of intergranular cracking

To further investigate the hydrogen induced decohesion for GBs, the evolution of intergranular cracks is explored using a CFE method. It begins to appear slight stress concentration in GBs instead of cracks, when the cohesive energy decreases to about 1.5 J m⁻². So we set eight linearly decreasing values of the cohesive energies less than 1.5 are adopted to simulate the intergranular cracking behaviour of the 2024 Al alloy described in Section 2.2. The cohesive energies and other input parameters according to the traction–separation law for the CFE calculation in Table 7.

Table 7 The input parameters of GBs for CFE calculations

G c = Ecoh (mJ mm^−2)	T ^0n (MPa)	K nn (MPa)	K ss = Ktt (MPa)
1.406 × 10^−3	140.6	140.66 × 10^6	52.86 × 10^6
1.336 × 10^−3	133.6	133.6 × 10^6	50.23 × 10^6
1.265 × 10^−3	126.5	126.5 × 10^6	47.56 × 10^6
1.195 × 10^−3	119.5	119.5 × 10^6	44.92 × 10^6
1.125 × 10^−3	112.5	112.5 × 10^6	42.29 × 10^6
1.055 × 10^−3	105.5	105.5 × 10^6	39.66 × 10^6
0.984 × 10^−3	98.4	98.4 × 10^6	36.99 × 10^6
0.914 × 10^−3	91.4	91.4 × 10^6	34.36 × 10^6

---

Fig. 9 shows the Mises stress nephogram of 3D polycrystal model with different fracture energies. The stress is inhomogeneously distributed in the model because the grains are randomly orientated. The microcracks usually initiates at a triple line (marked by white arrows) where several grains come together on the surfaces or in the interior of the material. Although there is a stress concentration (marked by red box) at the GBs, no obvious cracks are found on the surface for a fracture energy of 1.406 J mol⁻¹. When the fracture energy decreases to 1.265 J mol⁻¹, the stress concentrations expend from the GBs to the interior of the grains, producing microcracks. The stress concentrations are relaxed because of the crack propagation when the fracture energy is 1.125 J mol⁻¹ and multiple cracks are observed on the surface. The cracks propagate along GBs, which may lead to the fracture of the polycrystal as the fracture energy decreases further.

Fig. 9 The CFE calculation results of Mises stress nephogram of the polycrystalline models under the uniaxial tension with fracture energy of (a) 1.406 J m⁻²; (b) 1.265 J m⁻²; (c) 1.215 J m⁻² and (d) 0.984 J m⁻².

Fig. 10 shows the GB damage evolution defined as scalar stiffness degradation (SDEG) according to different fracture energies. The red cohesive elements (SDEG = 1) have failed, suggesting that cracks initiated. To quantify the degree of the damage, the number of failed cohesive elements is determined.

Fig. 10 The CFE calculation results of SDEG nephogram of the polycrystalline models under the uniaxial tension with fracture energy of (a) 1.406 J m⁻²; (b) 1.265 J m⁻²; (c) 1.215 J m⁻² and (d) 0.984 J m⁻².

Fig. 11 shows that the percentage of failed elements increases with decreasing fracture energy and that the descending rate firstly increases and then decreases. The polycrystal fractures completely when the descending rate of damage approaches zero. Fig. 12 shows the relationship between the percentage of failed elements and the deformation displacement of polycrystal with different fracture energies. The percentage of failed elements increases with the increasing of deformation displacement. For equivalent deformation displacement, the polycrystal with a lower fracture energy fractures more easily.

Fig. 11 The percentage of failed cohesive elements (D_dam) as a function of cohesive energies (E_coh) under uniaxial tension.

Fig. 12 The relationship between the percentage of failed cohesive elements and deformation displacement of the polycrystalline models under the uniaxial tension with different cohesive energies.

The results of the CFE calculations reveal that the initiation and propagation of cracks is closely associated with the fracture energy, which is sensitive to the H concentrations segregated in the GBs. The simulation is a good representation for the hydrogen-induced decohesion of GBs HE mechanism.

4. Conclusions

The approach presented here combines first-principles and CFE calculations to investigate the effect of H atoms on intergranular cracking from an atomic scale to a grain level. The effects of H atoms on intergranular embrittlement from the adsorption and diffusion of H atoms to intergranular cracking induced by H segregation was studied. The simulations and calculations investigated the intergranular embrittlement mechanism from an energy standpoint. The new approach realizes the cross-scale calculation via the calculation of cohesive energy. H atoms can be adsorbed on fcc site and top site. Fcc site is the most stable H-adsorbed site on Al (111) surface at the coverages of 1/9, 2/9 and 1/3 ML, and top site takes second place. H₂ can be generated when the coverage of H atoms increases to 1/3 ML. The optimal diffusion pathway from the most stable surface fcc site to the TIS site which is more stable than OIS in bulk, is fcc–hcp–TIS_1-1–TIS₂. The migration of H atom is easier on surface than that in the interior and H atom must overcome a higher energy barrier to move from surface to subsurface. Subsequently, the energy barrier for bulk migration of the H is lower. The segregation of H atoms at GBs weakens the strength of GBs. The cohesive energies of GBs decrease linearly with increasing of H concentrations. Thus CFE calculations effectively simulate the evolution of intergranular cracking, which is a good representation for HE of hydrogen induced decohesion for GBs.

Acknowledgements

This work is supported by the National Natural Science Foundation of China (No. 51222106), the Fundamental Research Funds for the Central Universities (No. 230201306500002) and National Basic Research Program of China (973 Program) (No. 2014CB643300).

References

1. P. V. Petroyiannis, A. T. Kermanidis, P. Papanikos and Sp. G. Pantelakis, Theor. Appl. Fract. Mech., 2004, 41, 173–183.
2. R. G. Song, W. Dietzel, B. J. Zhang, W. J. Liu, M. K. Tseng and A. Atrens, Acta Mater., 2004, 52, 4727–4743.
3. F. Y. Guo, C. G. Long, J. Zhang, Z. Zhang, C. H. Liu and K. Yu, Appl. Surf. Sci., 2015, 324, 584–589.
4. R. J. Gest and A. R. Troiano, Corrosion, 1974, 30, 274–279.
5. E. Pouillier, A.-F. Gourgues, D. Tanguy and E. P. Busso, Int. J. Plast., 2012, 34, 139–153.
6. H. Matsui, H. Kimura and S. Moriya, J. Mater. Sci. Eng., 1979, 40, 207–216.
7. G. Bertolino, G. Meyer and J. Perez Ipiña, J. Alloys Compd., 2002, 330–332, 408–413.
8. H. Ji, I. I.-J. Park, S.-M. Lee and Y.-K. Lee, J. Alloys Compd., 2014, 598, 205–212.
9. A. Troiano, Trans. ASM, 1960, 52, 54–80.
10. J. Hirth, Metall. Trans. A, 1980, 11, 861–890.
11. D. F. Teter, I. M. Robertson and H. K. Birnbaum, Acta Mater., 2001, 49, 4313–4323.
12. M. Koyama, C. C. Tasan, E. Akiyama, K. Tsuzaki and D. Raabe, Acta Mater., 2014, 70, 174–187.
13. P. J. Ferreira, I. M. Robertson and H. K. Birnbaum, Acta Mater., 1998, 46, 1749–1757.
14. M. Niwa, T. Shikama and A. Yonezu, J. Mater. Sci. Eng. A, 2015, 624, 52–61.
15. J. Song and W. A. Curtin, Acta Mater., 2011, 59, 1557–1569.
16. Z. Guo, M. Zhao, C. Li, S. Chen and L. Rong, J. Mater. Sci. Eng. A, 2012, 555, 77–84.
17. O. M. Alyousif and R. Nishimura, Corros. Sci., 2012, 52, 7–13.
18. Y. Wang, X. Wang, J. Gong, L. Shen and W. Dong, Int. J. Hydrogen Energy, 2014, 39, 13909–13918.
19. L. Zhang, Z. Li, J. Zheng, Y. Zhao, P. Xu, C. Zhou and X. Li, Int. J. Hydrogen Energy, 2013, 38, 8208–8214.
20. A. Barnoush and H. Vehoff, Acta Mater., 2010, 58, 5274–5285.
21. A. Barnoush and H. Vehoff, Corros. Sci., 2008, 50, 259–267.
22. S. G. Xu and Y. L. Zhao, Eng. Failure Anal., 2013, 34, 24–34.
23. M. Kikuchi, Y. Wada, Y. Shimizu and Y. L. Li, Int. J. Pressure Vessels Piping, 2012, 90–91, 2–8.
24. T. Luther and C. Konke, Eng. Fract. Mech., 2009, 76, 2332–2343.
25. P. Zhang, M. Karimpounr, D. Balint, J. Lin and D. Farrugia, Comput. Mater. Sci., 2012, 64, 84–89.
26. I. Simonovski and L. Cizelj, Comput. Mater. Sci., 2011, 50, 1606–1618.
27. J. Zhang, Z. H. Chen and C. F. Dong, J. Mater. Eng. Perform., 2015, 24, 4908–4918.
28. I. Simonovski and L. Cizelj, Eng. Fract. Mech., 2013, 110, 364–377.
29. V. V. Bulatov, B. W. Reed and M. Kumar, Surf. Sci., 2014, 65, 161–175.
30. L. Vitos, A. V. Ruban, H. L. Skriver and J. Kollár, Surf. Sci., 1998, 411, 186–202.
31. P. J. Feibelman, Phys. Rev. B: Condens. Matter Mater. Phys., 1996, 53, 13740–13746.
32. M. Yamamoto, C. T. Chan and K. M. Ho, Phys. Rev. B: Condens. Matter Mater. Phys., 1994, 50, 7932–7939.
33. W. R. Tyson and W. A. Miller, Surf. Sci., 1977, 62, 267–276.
34. G. Kresse and J. Hafner, Phys. Rev. B: Condens. Matter Mater. Phys., 1993, 47, 558–561.
35. G. Kresse and J. Furthmüller, J. Phys. Chem. B, 1996, 54, 11169–11186.
36. W. Kohn and L. J. Sham, Phys. Rev. A, 1965, 140, 1133–1138.
37. G. Kresse and D. Joubert, Phys. Rev. B: Condens. Matter Mater. Phys., 1999, 59, 1758–1775.
38. J. P. Perdew and K. Burke, Phys. Rev. Lett., 1996, 77, 3865–3868.
39. K. Momma and F. Izumi, J. Appl. Crystallogr., 2008, 41, 653–658.
40. H. Monkhorst and J. Pack, Phys. Rev. B: Solid State, 1976, 13, 5188–5192.
41. R. F. W. Bader, Atoms in Molecules-A Quantum Theory, Oxford University Press, New York, 1990.
42. G. Henkelman, A. Arnaldsson and H. Jónsson, Comput. Mater. Sci., 2006, 36, 354–360.
43. R. Quey, P. R. Dawson and F. Barbee, Comput. Meth. Appl. Mech. Eng., 2011, 200, 1729–1745.
44. M. N. Krivosheina, S. V. Kobenko and M. A. Kozlova, J. Phys.: Conf. Ser., 2009, 181, 012083.
45. T. Luther and C. Konle, Eng. Fract. Mech., 2009, 76, 2332–2343.
46. P. Zhang, M. Karimpounr, D. Balint, J. Lin and D. Farrugia, Comput. Mater. Sci., 2012, 64, 84–89.
47. I. Simonovski and L. Cizelj, Comput. Mater. Sci., 2011, 50, 1606–1618.
48. H. B. Michaelson, J. Appl. Phys., 1977, 48, 4729–4733.
49. J. B. Condon and T. Schober, J. Nucl. Mater., 1993, 207, 1–24.
50. A. Oudriss, J. Creus, J. Bouhattate, C. Savall, B. Peraudeau and X. Feaugas, Scr. Mater., 2012, 66, 37–40.
51. L. Vitos, A. V. Ruban, H. L. Skriver and J. Kollár, Surf. Sci., 1998, 411, 186–202.
52. V. V. Bulatov, B. W. Reed and M. Kumar, Acta Mater., 2014, 65, 151–175.
53. D. L. Olmsted, S. M. Foiles and E. A. Holm, Acta Mater., 2009, 57, 3694–3703.

---