A crack tip in α-Fe presents a preferential trap site for hydrogen, and sufficient concentration of hydrogen can change the incipient crack tip deformation response, causing a transition from a ductile to a brittle failure mechanism for inherently ductile alloys. In this work, the effect of hydrogen segregation around the crack tip on deformation in α-Fe was examined using atomistic simulations and the continuum based Rice-Thompson criterion for various modes of fracture (I, II, and III). The presence of a hydrogen rich region ahead of the crack tip was found to cause a decrease in the critical stress intensity factor required for incipient deformation for various crack orientations and modes of fracture examined here. Furthermore, the triaxial stress state ahead of the crack tip was found to play a crucial role in determining the effect of hydrogen on the deformation behavior. Overall, the segregation of hydrogen atoms around the crack tip enhanced both dislocation emission and cleavage behavior suggesting that hydrogen has a dual role during the deformation in α-Fe.

Hydrogen embrittlement (HE) is a phenomenon that affects both the physical and chemical properties of several intrinsically ductile metals.1–4 In the past, various mechanisms in the literature including hydrogen-enhanced decohesion,5–9 hydrogen-enhanced local plasticity,10–12 adsorption-induced hydrogen-enhanced dislocation emission,13,14 and cleavage15–17 have been developed to explain the HE phenomenon in iron (α-Fe) as well as other crystalline materials. However, these theories are often divergent, resulting in a lack of agreement on the underlying mechanism. For instance, first principles18 and atomistic19,20 calculations of hydrogen (H) segregation along grain boundaries have reported a significant increase in the susceptibility for decohesion of grain boundaries. On the other hand, in situ TEM findings10,21,22 suggest that sufficient concentration of H around dislocations can result in enhanced dislocation mobility. Furthermore, the theoretical12,23 and experimental10,21,22 evidence shows that the presence of H atoms can alter the elastic interactions between the dislocations and microstructural obstacles (shielding effect), thereby decreasing the resistance to dislocation glide. Moreover, based on electrochemical nano-indentation experiments,24,25 it was proposed that the presence of H reduces the defect formation energy, thereby resulting in a decrease in the stress required for dislocation nucleation.26 Overall, the enhanced dislocation activity leads to an increased tendency to form dislocation pileups that could result in void nucleation.27 

In general, H permeates the microstructure, either through rapid diffusion along grain boundary networks or by ingress through cracks, and migrates along with the dislocations emitted from the crack tip. Thus, it becomes difficult to experimentally examine the effect of H on microstructural features individually, resulting in a lack of comprehensive understanding of the HE mechanism. This has motivated a great deal of modelling efforts, particularly molecular dynamics19,27–30 and quantum mechanics-based31–33 simulations that can accurately describe the interatomic interaction and provide atomic insights into the HE mechanism. These studies include hydrogen-dislocation core interactions,30,32 grain boundary (GB) segregation,19,20,34–36 and the effect of H on crack tip events.28,37 The crack tip presents a preferential trap site for H due to the high hydrostatic stresses ahead of the crack tip. Therefore, the accumulation of sufficient H concentration ahead of the crack tip has been reported to alter the crack tip deformation mechanism.38,39 Furthermore, during cyclic loading, the presence of a H rich environment drastically reduces the threshold stress intensity factor,40–42 resulting in a decrease of the functional life of structural alloys. Hence, a comprehensive understanding of the crack tip deformation during HE phenomena is critical to the development of an accurate life prediction tool.43,44

In this work, we examine the influence of H on the incipient crack tip events in α-Fe using atomistic simulations coupled with the displacement field corresponding to a semi-infinite crack subjected to various levels of applied stress intensity factor. A range of crack orientations in α-Fe was investigated to understand the effect of cleavage/slip plane orientation on the deformation behavior while subjecting to different modes of fracture. This research article is organized in the following manner. Section II describes the computational methods for studying the incipient crack tip response and the introduction of equilibrium H concentrations ahead of the crack tip. Then, we discuss the effect of H on the crack tip deformation by quantifying the critical stress intensity factor and analyzing the deformation behavior (cleavage or dislocation emission) in light of small scale yielding. For instance, we observed that the presence of H atoms ahead of the crack tip caused a decrease in the critical stress intensity factor across the various crack orientations and modes of fracture. However, the cracks subjected to Mode I deformation observed the largest decrease in the critical stress intensity factor required for the incipient crack tip events. This clearly highlights the role of stress triaxiality (in general, the stress-state) in the HE phenomena. Finally, the dual role of H in the crack tip deformation was observed by examining the incipient crack tip deformation across various crack orientations and modes of fracture. Overall, the present study provides novel insights into the hydrogen embrittlement mechanism in α-Fe and its effect on crack tip deformation.

A parallel molecular dynamics code (large-scale atomic/molecular massively parallel simulator, LAMMPS45) with a semi-empirical embedded atom method (EAM)46 potential was used to study the role of H in the incipient crack tip deformation for various orientations in α-Fe (refer to Table I for details). However, it is essential to understand the limitations of atomistic simulations for studying various material phenomena/properties.47–51 For instance, as discussed in the study by Adlakha et al.,52 for a small-scale yielding, the domain of K-dominance (K is the stress intensity factor) around a crack tip is a fraction of the crack length, i.e., in the atomistic simulations, if the crack is on the order of nano-meters then the K-dominant region will be on the order of a few angstroms. Thus, if a dislocation is emitted from such a small crack, it will be pushed into a non-K region of the specimen and in that region its final position is controlled by the entire specimen geometry and loading (non K-dominated region). Furthermore, the new emitted dislocation will exert image stresses on the crack tip which will influence subsequent growth behavior.53,54 Therefore, no events beyond the first emission can possibly be properly/accurately captured through such simulations with K based displacement boundary conditions. Hence, in the present work, only the incipient crack tip deformation event is discussed as opposed to the subsequent events leading up to the crack growth.

TABLE I.

The different crystallographic orientations employed in the Modes I, II, and III specimen.

XYZ
Orientation-A 100 010 001 
Orientation-B 11¯0 110 001 
Orientation-C 100 011 01¯1 
Orientation-D 111 1¯1¯2 11¯0 
XYZ
Orientation-A 100 010 001 
Orientation-B 11¯0 110 001 
Orientation-C 100 011 01¯1 
Orientation-D 111 1¯1¯2 11¯0 

In this work, we employed the EAM potential of Ramasubramaniam et al.55 to describe the Fe-H system, which is based on the Fe EAM potential of Hepburn and Ackland.56 The Fe-H potential was parameterized using an extensive database of energies and configurations from density functional theory (DFT) calculations of the dissolution and diffusion of H in bulk α-Fe, the binding of H to free surfaces, vacancies, and dislocations as well as other cross interactions between H and Fe. Moreover, the binding and the formation energies corresponding to multiple H-segregations to bulk α-Fe are consistent with the values predicted using ab initio calculations and experimentation.57 This potential has been used to study H-dislocation core interactions,30,32,37 GB segregation,19,20 and subsequent decohesion19 as well as the effect of H on generalized stacking fault energies58 and surface energies.55,59 In particular, Solanki et al.19 have used the Fe-H potential and showed that the local atomic arrangements within the grain boundary region and the resulting structural units have a significant influence on the magnitude of the binding energies of H atoms, which was further verified by Wang et al. in their recent effort.35 Similarly, Bhatia et al.30 have shown that the dislocation interaction with a H vacancy complex has a significant influence on obstacle strength, which was further validated by Tehranchi et al.60 in their recent work. Overall, this potential has been shown to be in good agreement with DFT results. Therefore, in this work the aforementioned potential was deemed appropriate for studying the role of H in the crack tip deformation.

Molecular statics (0 K) fracture simulations were performed using the displacement field corresponding to a semi-infinite crack subjected to the appropriate stress intensity factor, i.e., KI, KII, and KIII at infinity (Fig. 1). These simulations were used to quantify the critical stress intensity factor for the incipient crack tip deformation and to corroborate with the Rice-Thompson (RT) predictions61 (see  Appendix A). The minimum dimension for the entire specimen at equilibrium was approximately 200 Å in radius and was 40 Å thick along the Z direction (approximately 300 000 atoms). An atomistically sharp through crack of length 40 Å was introduced by deleting atoms within a small region (Fig. 1). A free boundary condition was prescribed along the X and Y direction (crack plane normal), while the Z direction (crack front direction) was modelled using periodic boundary conditions. The atoms within 10 Å of the cylindrical surface were defined as boundary atoms (grey region in Fig. 1). The cylindrical specimens were loaded by prescribing displacements along the X(u), Y(v), and Z(w) directions using the following boundary conditions to model the various modes of fracture.

  • Mode I
    uIr,θ=KI1+ϑEr2πcosθ224ϑ+2sin2θ2,vIr,θ=KI1+ϑEr2πsinθ244ϑ+2cos2θ2.
    (1)
  • Mode II
    uIIr,θ=KII1+ϑEr2πsinθ244ϑ+2cos2θ2,vIIr,θ=KII1+ϑEr2πcosθ224ϑ+2sin2θ2,
    (2)
  • Mode III
    wIIIr,θ=4KIII1+ϑEr2πsinθ2,
    (3)

where KI,KII, and KIII are the stress intensity factors for Mode I, II, and III loading, respectively, ϑ is the Poisson's ratio (ϑ= 0.3), E is the Young's modulus, and the equations are in polar coordinates with r and θ corresponding to the position of atoms within the boundary region with origin defined at the crack tip as shown in Fig. 1. The displacement fields were varied by incrementally changing the stress intensity factor (ΔK = 1.5 × 10−3 MPa m);62,63 each incremental displacement was followed by a nonlinear conjugate gradient energy minimization process.

FIG. 1.

A schematic representation of the fracture specimen that was subjected to various modes of fracture by displacing boundary atoms (grey) within the atomistic framework to study the effect of H on the incipient crack tip deformation.

FIG. 1.

A schematic representation of the fracture specimen that was subjected to various modes of fracture by displacing boundary atoms (grey) within the atomistic framework to study the effect of H on the incipient crack tip deformation.

Close modal

The influence of H on the incipient crack tip deformation behavior was examined by introducing H atoms around the crack tip. The H occupancy at a tetrahedral site (θi) is dependent on the H binding energy (Ebi) and the temperature (T), which can be expressed in the following manner:64 

θi1θi=θbulk1θbulkexpEbikBT,
(4)

where θbulk is the atomic fraction of H for the whole system, and kB is the Boltzmann constant. Note that the H atomic fraction of approximately 4 × 10−4 was chosen based on previous experimental works (see Refs. 35, 64, and 65) The binding energy indicates the energetic preference of a H atom to occupy an interstitial site in the vicinity of the crack tip as compared to a bulk lattice site. The binding energy for H to occupy a tetrahedral site near the crack tip can be approximated by evaluating the elastic energy contribution σmiΔV. Here, σmi is the atomic hydrostatic stress around site i and ΔV is the atomic misfit volume due to the presence of a H atom. Hydrogen atoms were initially placed around the crack tip based on the occupation probability [Eq. (4)]. Furthermore, the Monte-Carlo (MC) method was employed to obtain a realistic H distribution around the GB at a finite temperature (300 K) for 2 ns. In order to study incipient crack tip events, the atomic structure obtained at the end of the MC method was deformed using the methodology described earlier (refer Sec. III A).

In this section, we examine the incipient crack tip event along four different crack orientations as tabulated in Table I, under the Mode I, Mode II, and Mode III loading conditions. Note that since the K-dominant region will be on the order of a few angstroms, no crack tip events beyond the first event (emission versus cleave) were quantified; see more details in Refs. 51 and 52. Furthermore, the observed atomistic behaviors are then compared with the theoretical prediction based on the work of Rice-Thompson (RT)61 (see  Appendix A). In order to understand the effect of H on the atomic scale crack tip deformation, a nominal concentration of H (4 × 10−4) was introduced ahead of the crack tip. The incipient crack tip deformation mechanism (cleavage/dislocation nucleation) and the corresponding stress intensity factor KIc were compiled and compared for the two different environments (H free and H rich) across a wide range of crack orientations and loading conditions. The findings presented in this work aid in developing a better understanding of the hydrogen embrittlement phenomena in α-Fe. Finally, we highlight the role of local stress-state along with the H segregation in the incipient crack tip behavior under various conditions.

1. Orientation-A—(010)[001]

The influence of H on the incipient events for the 010001 crack orientation was examined here. Based on the surface energies of various closed packed planes (Table IV), the stress intensity factor for cleavage along the {011} planes, KI{011}=0.83MPam was found to be higher than the KIc predicted along the (010) plane [refer to Fig. 7(a) and Table V]. Based on atomistic simulations, we further observed that the crack advances through cleavage along the {011} planes aligned at an angle of 45° to the crack plane for both the H free and charged scenarios. However, previous experimental studies on Fe-3%Si for this crack orientation have observed crack propagation along the (010) plane.66,67 Such variation in comparisons could be attributed to the empirical nature of interatomic potential and the defect free conditions of atomistic modeling (pure Fe vs. Fe-3%Si).62,68 Nonetheless, both the modeling and experimental findings66,67 observe a cleavage/brittle response for the 010001 crack orientation under Mode I loading conditions. Interestingly, the presence of H (c0=4×104) ahead of the crack tip was found to decrease the stress intensity factor from 0.59MPam (no H) [Fig. 2(a)] to 0.48MPam (with H) allowing the crack to extend in a cleavage manner along the {011} plane [Fig. 2(b)]. The decrease in the stress intensity factor for cleavage advance KIc due to the presence of H atoms can be attributed to the decrease in surface energies for the closed packed planes.69 In the case of Mode II loading, the RT criteria predict the nucleation of a mixed dislocation along the {112}⟨111⟩ slip system [refer to Fig. 7(b) and Table V for details]. However, using atomistic simulations we observed the nucleation of an edge dislocation along the {010}⟨100⟩ slip system KII=1.04MPam [Fig. 2(c)]. Additionally, a BCC to FCC phase transformation70 was also observed ahead of the crack tip that is a known artifact of the high stress state achieved ahead of the crack tip in atomistic simulations.69 In the presence of H ahead of the crack tip, the stress intensity factor required for the emission of the edge dislocation decreases to KIIH=0.94MPam [Fig. 2(d)]. In the case of Mode III loading, the nucleation of a screw dislocation along the {110}[001] slip system KIII=0.44MPam was observed [Fig. 2(d)] for the H free scenario ahead of the crack tip. The dislocation nucleation along the {110}[001] slip system is not in agreement with the RT criteria predictions [refer to Table V and Fig. 7(c)]. On the other hand, in the presence of H (c0=4×104) ahead of the crack tip, the nucleation of a mixed character dislocation on the {110}⟨111⟩ slip system was observed at a stress intensity factor of 0.25MPam. In summary, the stress intensity factor required for both cleavage and dislocation emission was found to decrease in the presence of H for the different loading conditions examined for the 010001 crack orientation. This clearly highlights the dual role of the effect of H on the crack tip deformation, i.e., dislocation emission or crack extension through cleavage.

FIG. 2.

(a) The incipient cleavage crack advance along the {110} plane for the 010001 orientation subjected to Mode I loading conditions. (b) The presence of a H rich environment ahead of the crack tip reduces the critical stress intensity factor required for cleavage. The atoms were colored according to atomic volume estimated by the Voronoi tessellation. The crack tip response under Mode II (c) and (d) and Mode III (e) and (f) was examined for both with and without H environments. (c)–(f) The atoms were colored according to the local atomic shear strain invariant (γs).

FIG. 2.

(a) The incipient cleavage crack advance along the {110} plane for the 010001 orientation subjected to Mode I loading conditions. (b) The presence of a H rich environment ahead of the crack tip reduces the critical stress intensity factor required for cleavage. The atoms were colored according to atomic volume estimated by the Voronoi tessellation. The crack tip response under Mode II (c) and (d) and Mode III (e) and (f) was examined for both with and without H environments. (c)–(f) The atoms were colored according to the local atomic shear strain invariant (γs).

Close modal

2. Orientation-B—(110)[001]

Here, we investigate the effect of H on the incipient crack tip deformation for the 110001 crack orientation (refer to Table I for details) subjected to different modes of fracture. Under Mode I loading conditions, the RT criteria predict a cleavage advance of the crack along the (110) plane at a stress intensity factor of 0.87MPam [refer to Fig. 7(d) and Table V]. The onset of cleavage advance along the (110) plane was observed at a stress intensity factor of 1.01MPam and 0.83MPam for the H free and H charged scenarios, respectively. In the case of Mode II loading, the nucleation of a mixed dislocation along the {110}⟨111⟩ slip system was observed for both the H free KII=0.58MPam and H charged KIIH=0.42MPam environments. This incipient crack tip deformation behavior agreed with the RT criteria predictions [refer to Fig. 7(e) and Table V] and with the introduction of hydrogen we observed a decrease in the critical stress intensity factor.13 Lastly, for the Mode III loading conditions, we observed nucleation of a mixed character dislocation along the 11¯0[111] slip system for both scenarios KIII=0.4MPamandKIIIH=0.37MPam in agreement with theoretical predictions [refer to Fig. 7(f) and Table V].

3. Orientation-C—(011)[01¯1]

In this section, we will discuss the effect of H on the 01101¯1 crack orientation (refer to Table I). Based on the RT criteria and linear elastic fracture mechanics (LEFM) theory, the crack advance (cleavage) along the (011) plane was predicted to occur at a critical stress intensity factor, KIc=0.87MPam for Mode I loading conditions (refer  Appendix A). However, using atomistic simulations, an edge dislocation KI=1.00MPam along the 21111¯1¯ slip system was observed that was in agreement with previous findings71 [Fig. 3(a)]. Interestingly, the introduction of H atoms ahead of the crack tip alters the incipient crack tip deformation mechanism (edge dislocation  cleavage advance at KIH=0.89MPam) [Fig. 3(b)]. The transition in the incipient deformation mechanism ahead of the crack tip can be attributed to a reduction in the surface energy of the (011) plane due to hydrogen, which is in agreement with both the theoretical and experirmental observations.38,69 Next, in the case of Mode II loading, an edge dislocation along the 2111¯11 slip system for both the H free KII=0.37MPam and H charged KIIH=0.38MPam specimens was emitted from the crack tip [Figs. 3(c) and 3(d)], which is in disagreement with the theoretical predictions [refer to Fig. 7(h) and Table V]. In this case, there was no noticeable effect of hydrogen on the incipient crack tip deformation. For the Mode III loading, the emission of a mixed dislocation on the {110}⟨111⟩ slip system for both cases KIII=0.56MPamandKIIIH=0.72MPam was observed, which was found to be in agreement with RT criteria predicitons [refer Fig. 7(i) and Table V]. Furthermore, the location of the emitted dislocation at a stress intensity factor of 0.72MPam was compared in different environments (with and without H) in Figs. 3(e) and 3(f), to highlight the large distance traversed by the dislocation in a H free enviroment as compared to the presence of H. A plausible explanation for the observed behavior could be that the presence of H for this crack orientation alters the elastic strain fields, resulting in a significant increase in the resistance for dislocation nucleation and propogation. However, this is contrary to the hypothesis proposed by Lynch;13 therefore, further investigations are needed to verify this hypothesis with more focus on the kinetic aspect of hydrogen.

FIG. 3.

The effect of H on the crack tip deformation for the 01101¯1 orientation was examined under the influence of different modes of fracture. (a) The emission of an edge dislocation along the 21111¯1¯ slip system was observed for Mode I loading conditions in the absence of H (atoms were colored according to the local shear strain invariant). (b) The presence of a H rich environment ahead of the crack was found to induce a cleavage response (atoms were colored based on the atomically derived volumetric strain). (c) In the presence of H atoms, an edge dislocation along the 21111¯1¯ slip system was emitted ahead of the crack tip under Modes II and (d) shows the dislocation line and the Burger's vector direction using the DXA algorithm.72 (e) and (f) In the case of Mode III loading, the nucleation of a mixed dislocation was observed for both the environments, i.e., (e) no H and (f) with H. (c)–(e) The atoms were colored according to the local atomic shear strain invariant (γs).

FIG. 3.

The effect of H on the crack tip deformation for the 01101¯1 orientation was examined under the influence of different modes of fracture. (a) The emission of an edge dislocation along the 21111¯1¯ slip system was observed for Mode I loading conditions in the absence of H (atoms were colored according to the local shear strain invariant). (b) The presence of a H rich environment ahead of the crack was found to induce a cleavage response (atoms were colored based on the atomically derived volumetric strain). (c) In the presence of H atoms, an edge dislocation along the 21111¯1¯ slip system was emitted ahead of the crack tip under Modes II and (d) shows the dislocation line and the Burger's vector direction using the DXA algorithm.72 (e) and (f) In the case of Mode III loading, the nucleation of a mixed dislocation was observed for both the environments, i.e., (e) no H and (f) with H. (c)–(e) The atoms were colored according to the local atomic shear strain invariant (γs).

Close modal

4. Orientation-D—(1¯1¯2)[11¯0]

In this section, we will discuss our findings on the effect of H on the crack tip deformation mechanism for the 1¯1¯211¯0 orientation (refer to Table I). Under the Mode I loading conditions, based on LEFM considerations, crack growth was predicted to occur along the 1¯1¯2 plane (KI(1¯1¯2)=0.92MPam) (for details refer to  Appendix A). Additionally, the stress intensity factor for cleavage advance along the 1¯1¯1 plane KI1¯1¯1=0.97MPam closely matched the minimum stress intensity factor along the 1¯1¯2 plane [refer to Fig. 7(j)]. Interestingly, the presence of H atoms ahead of the crack tip had no effect on the observed crack advance along the 1¯1¯2 plane at KI=0.94MPam. In the case of Mode II loading, an edge dislocation emitted from the crack tip along the 1¯1¯2 plane at KII=0.37MPam for a H free environment [Fig. 4(a)]. This was consistent with the RT criteria predictions for incipient crack tip deformation [refer to Fig. 7(k) and Table V]. In the presence of H atoms ahead of the crack tip, an edge dislocation was emitted along the 112 plane from the crack tip at a significantly lower critical stress intensity factor KIIH=0.19MPam compared to the H free scenario KII=0.37MPam and was found to be in agreement with previous work58 [Fig. 4(b)]. In the case of Mode III loading, a mixed dislocation along the 01¯11¯11 slip system was emitted at a KIII of 0.26MPam and 0.36MPam for the H free and charged scenarios [refer to Figs. 4(c) and 4(d)]. The nucleation of a mixed dislocation was in qualitative agreement with the RT criteria [refer to Fig. 7(k) and Table V]. However, once again the presence of H ahead of the crack tip results in an increase in the stress intensity factor required for dislocation emission.

FIG. 4.

The effect of H on the crack tip deformation for the 1¯1¯211¯0 orientation subjected to different modes of fracture, i.e., Mode II (a) and (b) and Mode III (c) and (d). (a) The emission of an edge dislocation along the 1¯1¯2111 slip system was observed for Mode II loading conditions in the absence of H. (b) The presence of a H rich environment ahead of the crack caused dislocation emission along the 1121¯11 slip system for Mode II. (c) and (d) In the case of Mode III, dislocation emission along the 01¯11¯11 slip system was observed both in the (c) absence and (d) presence of H. The atoms were colored according to the local atomic shear strain invariant (γs).

FIG. 4.

The effect of H on the crack tip deformation for the 1¯1¯211¯0 orientation subjected to different modes of fracture, i.e., Mode II (a) and (b) and Mode III (c) and (d). (a) The emission of an edge dislocation along the 1¯1¯2111 slip system was observed for Mode II loading conditions in the absence of H. (b) The presence of a H rich environment ahead of the crack caused dislocation emission along the 1121¯11 slip system for Mode II. (c) and (d) In the case of Mode III, dislocation emission along the 01¯11¯11 slip system was observed both in the (c) absence and (d) presence of H. The atoms were colored according to the local atomic shear strain invariant (γs).

Close modal

The RT criteria61 are based on the Peierls-Nabarro concept73,74 that the shear stress along the slip plane is a periodic function of the shear displacement. This identifies the unstable stacking fault energy as the relevant material property necessary to estimate the stress intensity factor required for dislocation emission ahead of the crack tip. However, this criterion does not take into account the tensile softening effect along the slip plane that results in the decrease of the energy required for dislocation nucleation.75 Furthermore, the RT criteria fail to account for the lattice trapping,76 surface correction to the bulk unstable stacking fault energy77,78 and anisotropic elastic effects.79 Therefore, in this work we have focused on comparing the predicted crack tip deformation mechanism with the atomistic observations instead of the magnitude of the stress intensity factor required for the incipient deformation ahead of the crack. Interestingly, we observed dislocations with a ⟨100⟩ Burger's vector that require overcoming a large energy barrier for the dislocation nucleation (refer to Table II and  Appendix A). However, there have been similar observations using the quasi-continuum framework during the study of crack tip deformation in α-Fe.62 In the case of crack orientation A, crack advance occurred along the {110} plane (Table II) instead of the {100} plane.66,67 The lower surface energy along the {110} surface for the EAM potential instead of the most active {100} plane experimentally (refer to Table IV) is a plausible explanation. In summary, the RT criteria were found to have mixed success in predicting atomistic observations presented here (refer to Table II and  Appendix A). Similar discrepancies have been noted in previous studies.62,68,80,81

TABLE II.

The effect of H on the incipient crack tip deformation mechanism for the various crack orientations under different modes of fracture.

OrientationFracture modeWithout HWith H
Cleavage along the (110) plane Cleavage along the 11¯0 plane 
II Edge dislocation along the 010001 slip system Edge dislocation along the 010001 slip system 
III Screw dislocation along the 11¯0100 slip system Mixed dislocation along the 11011¯1 slip system 
Cleavage along the 110 plane Cleavage along the 110 plane 
II Mixed dislocation along the 1101¯11 slip system Mixed dislocation along the 1101¯11 slip system 
III Mixed dislocation along the 11¯0111 slip system Mixed dislocation along the 1¯10111 slip system 
Edge dislocation along the 21111¯1¯ slip system Cleavage along the 011 plane 
II Edge dislocation 2111¯11 slip system Edge dislocation 2111¯11 slip system 
III Mixed dislocation along the 01111¯1 slip system Mixed dislocation along the 01111¯1 slip system 
Cleavage along the 1¯1¯2 plane Cleavage along the 1¯1¯2 plane 
II Edge dislocation along the 1¯1¯2111 slip system Edge dislocation along the 1121¯11 slip system 
III Mixed dislocation along the 01¯11¯11 slip system Mixed dislocation along the 01¯11¯11 slip system 
OrientationFracture modeWithout HWith H
Cleavage along the (110) plane Cleavage along the 11¯0 plane 
II Edge dislocation along the 010001 slip system Edge dislocation along the 010001 slip system 
III Screw dislocation along the 11¯0100 slip system Mixed dislocation along the 11011¯1 slip system 
Cleavage along the 110 plane Cleavage along the 110 plane 
II Mixed dislocation along the 1101¯11 slip system Mixed dislocation along the 1101¯11 slip system 
III Mixed dislocation along the 11¯0111 slip system Mixed dislocation along the 1¯10111 slip system 
Edge dislocation along the 21111¯1¯ slip system Cleavage along the 011 plane 
II Edge dislocation 2111¯11 slip system Edge dislocation 2111¯11 slip system 
III Mixed dislocation along the 01111¯1 slip system Mixed dislocation along the 01111¯1 slip system 
Cleavage along the 1¯1¯2 plane Cleavage along the 1¯1¯2 plane 
II Edge dislocation along the 1¯1¯2111 slip system Edge dislocation along the 1121¯11 slip system 
III Mixed dislocation along the 01¯11¯11 slip system Mixed dislocation along the 01¯11¯11 slip system 

The ability of H atoms to segregate around the crack tip and induce a sudden failure has been linked with the tensile hydrostatic stress field ahead of the crack tip.5 However, experimental observations39,82–88 have been unable to reach an agreement on the exact role of the stress state in the H induced cracking phenomena. For instance, the role of the tensile hydrostatic stress field in inducing failure has been questioned due to the observation of H cracking in the absence of a hydrostatic stress field.87–89 Therefore, to provide novel insights we have examined the atomic scale effects of H on the crack tip deformation mechanism under various modes of loading conditions. Due to the length and time scale limitations of atomistic simulations, we have focused on understanding the effect of H on the incipient inelastic crack tip deformation in α-Fe. For the crack orientations examined here, it was found that the presence of H around the crack tip generally decreases the stress intensity factor required for the initation of inelastic deformation (Fig. 5). In the case of Mode I loading, the presence of H atoms ahead of the crack tip was found to effect the largest decrease in the stress intensity factor required for brittle crack advance (cleavage) (refer to Fig. 5 and Table II). Furthermore, the drop observed in the stress intensity factor due to the presence of H reduces from Mode I to Mode III loading for the various crack orientations examined here (Fig. 5). These findings clearly highlight the role of stress state ahead of the crack tip during H embrittlement. Therefore, these atomic-scale insights provide a better understanding of the role of H in crack tip deformation in α-Fe.

FIG. 5.

A summary of the stress intensity factor for the incipient deformation of the various crack orientations subjected to Mode I, II, and III loading conditions. The black and red markers distinguish between the no H and H rich c0=4×104 scenarios ahead of the crack tip. Note: the filled red markers indicate a decrease in stress intensity factor for incipient deformation due to the presence of H. The partially filled red marker indicates a change in the deformation mechanism due to the presence of H ahead of the crack tip. In the case of orientation D, the KI was unchanged under both environments.

FIG. 5.

A summary of the stress intensity factor for the incipient deformation of the various crack orientations subjected to Mode I, II, and III loading conditions. The black and red markers distinguish between the no H and H rich c0=4×104 scenarios ahead of the crack tip. Note: the filled red markers indicate a decrease in stress intensity factor for incipient deformation due to the presence of H. The partially filled red marker indicates a change in the deformation mechanism due to the presence of H ahead of the crack tip. In the case of orientation D, the KI was unchanged under both environments.

Close modal

The following conclusions can be drawn based on studying the influence of H on the incipient crack tip events in α-Fe:

  1. In the case of Mode I loading, the presence of H ahead of the crack tip was found to significantly decrease the stress intensity factor required for crack propagation (brittle/cleavage response). Furthermore, in the case of the 01101¯1 crack orientation the presence of H was found to cause a transition in the incipient deformation mechanism, i.e., an edge dislocation (no H) to cleavage advance (with H).

  2. In the case of Mode II and III loading, H enhanced dislocation emission was observed for the different crack orientations examined here. Thus, providing merit to the notion that the presence of H along the slip planes ahead of the crack can decrease the interatomic interactions resulting in a decrease in the energy required for dislocation nucleation.13 

  3. Despite the debate regarding the role of stress state in H induced failure, this work provides atomistic insights that H embrittlement occurs across the various modes of fracture examined here. However, in the case of Mode III, the presence of H was found to have a reduced impact on the stress intensity factor required for dislocation emission.

  4. Finally, in this study, we have successfully demonstrated the co-existence of both characteristics of H embrittlement: (a) the segregation of sufficient H atoms ahead of the crack tip induces dislocation emission;13 and (b) the presence of H promotes a brittle crack response by causing a decrease in the stress intensity factor under Mode I loading conditions.

The authors gratefully acknowledge support from the Office of Naval Research under Contract N00014-16-1-2174.

The Rice-Thompson (RT) criteria61 and the linear elastic fracture mechanics (LEFM) theory90 were employed to theoretically examine the critical stress intensity factor and the mechanism (cleavage/dislocation character) of the incipient crack tip deformation for the various crack orientations (Table I) subjected to different modes of fracture (I, II and III). The RT criteria are dependent on the relative angle of the slip plane to the crack plane (χ) and the angle between normal to the crack tip and the slip direction (ϕ) (Fig. 6). The critical stress intensity factor for dislocation emission ahead of the crack tip subjected to the various modes of facture can be expressed in the following manner:

KIc=1cos2χ/2sinχ/2cosϕ2μ1ϑcos2ϕ+1ϑsin2ϕγusKIIc=1cosχ/213sin2χ/2cosϕ2μ1ϑcos2ϕ+1ϑsin2ϕγusKIIIc=1cosχ/22μ1ϑcos2ϕ+1ϑsin2ϕγus,
(A1)

here, μ is the shear modulus and γus is the unstable stacking fault energy for the various possible slip systems listed in Table III, calculated using the EAM potential described in Sec. II.

FIG. 6.

A schematic representation of the slip plane (green) inclined at an angle, χ w.r.t the crack plane and ϕ defines the inclination between the crack normal and the Burger's vector.

FIG. 6.

A schematic representation of the slip plane (green) inclined at an angle, χ w.r.t the crack plane and ϕ defines the inclination between the crack normal and the Burger's vector.

Close modal
TABLE III.

The unstable stacking fault energies (γus) along the various slip systems in α-Fe.

{100}⟨010⟩{110}⟨001⟩{110}⟨111⟩{112}⟨111⟩{123}⟨111⟩
γus (J/m21.80 1.36 0.71 0.81 0.79 
{100}⟨010⟩{110}⟨001⟩{110}⟨111⟩{112}⟨111⟩{123}⟨111⟩
γus (J/m21.80 1.36 0.71 0.81 0.79 

Based on LEFM theory,90 the critical stress intensity factor required for cleavage crack growth under Mode I loading can be estimated by the following expression:

KIc=2Eγs1ϑ2C12+C22C1=34cosω/2+14sin3ω/2C2=14sinω/2+sin3ω/2,
(A2)

where ω is the angle between the crack plane and the cleavage plane and γs is the surface energy. The surface energy along the various possible cleavage planes is listed in Table IV.

TABLE IV.

The surface energies along the various cleavage planes in α-Fe, calculated using the discussed EAM potential.55 

{100}{110}{111}{112}
Surface energy, γs (J/m21.75 1.65 1.85 1.96 
{100}{110}{111}{112}
Surface energy, γs (J/m21.75 1.65 1.85 1.96 

Using the theoretical relationships described earlier, the two energetically favorable slip/cleavage plane orientations ahead of the crack tip are schematically illustrated for the various crack orientations and the different modes of fracture in Fig. 7. Finally, the theoretically predicted incipient crack tip deformation for the various crack orientations examined here is presented in Table V.

FIG. 7.

(a) Illustration of the two energetically favorable slip/cleavage plane orientations ahead of the crack tip are schematically presented for the various crack orientations along with the different modes of fracture obtained using theoretical relationships.

FIG. 7.

(a) Illustration of the two energetically favorable slip/cleavage plane orientations ahead of the crack tip are schematically presented for the various crack orientations along with the different modes of fracture obtained using theoretical relationships.

Close modal
TABLE V.

A summary of the theoretical predictions for the incipient crack tip event for the various crack orientations under different modes of fracture.

Crack orientationFracture modeTheoretical predictions
(010)[001] Cleavage along the (010) planes, KIc=0.73MPam 
II Mixed dislocation along the 1¯21111¯ slip system, KIIc=0.49MPam 
III Mixed dislocation along the 1101¯11 slip system, KIIIc=0.29MPam 
(110)[001] Cleavage along the (011) planes, KIc=0.87MPam 
II Mixed dislocation along the (110)[111] slip system, KIIc=0.47MPam 
III Mixed dislocation along the 01111¯1 slip system, KIIIc=0.44MPam 
01101¯1 Cleavage along the (011) plane, KIc=0.87MPam 
II Edge dislocation on the (011)[100] slip system, KIIc=0.56MPam 
III Mixed dislocation along the 01111¯1 slip system, KIIIc=0.41MPam 
1¯1¯211¯0 Cleavage along the 1¯1¯2 plane, KIc=0.92MPam 
II Edge dislocation along the 1¯1¯2111 slip system, KIIc=0.43MPam 
III Mixed dislocation along the 01¯11¯11 slip system, KIIIc=0.43MPam 
Crack orientationFracture modeTheoretical predictions
(010)[001] Cleavage along the (010) planes, KIc=0.73MPam 
II Mixed dislocation along the 1¯21111¯ slip system, KIIc=0.49MPam 
III Mixed dislocation along the 1101¯11 slip system, KIIIc=0.29MPam 
(110)[001] Cleavage along the (011) planes, KIc=0.87MPam 
II Mixed dislocation along the (110)[111] slip system, KIIc=0.47MPam 
III Mixed dislocation along the 01111¯1 slip system, KIIIc=0.44MPam 
01101¯1 Cleavage along the (011) plane, KIc=0.87MPam 
II Edge dislocation on the (011)[100] slip system, KIIc=0.56MPam 
III Mixed dislocation along the 01111¯1 slip system, KIIIc=0.41MPam 
1¯1¯211¯0 Cleavage along the 1¯1¯2 plane, KIc=0.92MPam 
II Edge dislocation along the 1¯1¯2111 slip system, KIIc=0.43MPam 
III Mixed dislocation along the 01¯11¯11 slip system, KIIIc=0.43MPam 
1.
H. C.
Rogers
,
Science
159
,
1057
(
1968
).
2.
R. D.
Merrick
,
Mater. Perform.
28
,
2
(
1989
).
3.
National Research Council and National Academy of Engineering
,
The Hydrogen Economy: Opportunities, Costs, Barriers, and R&D Needs
(
National Academies Press
,
2004
).
4.
P. R.
Rhodes
,
L. A.
Skogsberg
, and
R. N.
Tuttle
,
Corrosion
63
,
63
(
2007
).
5.
A. R.
Troiano
,
Trans. ASM
52
,
54
(
1960
).
6.
R. A.
Oriani
,
Ber. Bunsengesellschaft Phys. Chem.
76
,
848
(
1972
).
7.
R. A.
Oriani
and
P. H.
Josephic
,
Acta Metall.
22
,
1065
(
1974
).
8.
W. W.
Gerberich
,
P.
Marsh
,
J.
Hoehn
,
S.
Venkataraman
,
H.
Huang
,
T.
Magnin
, and
J. M.
Gras
, in
Corrosion-Deformation Interactions (CDI'92)
(
Les Edition de Physique
,
Les Ulis, France
,
1993
), p.
325
.
9.
W. W.
Gerberich
,
P.
Marsh
,
J.
Hoehn
,
S.
Venkataraman
, and
H.
Huang
, in
Corrosion-Deformation Interactions (CDI'92)
, edited by
T.
Magnin
and
J. M.
Gras
(
Les Editions de Physique
,
Les Ulis
,
1993
), p.
325
.
10.
H. K.
Birnbaum
and
P.
Sofronis
,
Mater. Sci. Eng., A
176
,
191
(
1994
).
11.
D. P.
Abraham
and
C. J.
Altstetter
,
Metall. Mater. Trans. A
26
,
2859
(
1995
).
12.
P.
Sofronis
,
Y.
Liang
, and
N.
Aravas
,
Eur. J. Mech.-A/Solids
20
,
857
(
2001
).
13.
S. P.
Lynch
,
Acta Metall.
36
,
2639
(
1988
).
14.
Y.
Liang
and
P.
Sofronis
,
J. Mech. Phys. Solids
51
,
1509
(
2003
).
15.
S.
Gahr
,
M. L.
Grossbeck
, and
H. K.
Birnbaum
,
Acta Metall.
25
,
125
(
1977
).
16.
D. S.
Shih
,
I. M.
Robertson
, and
H. K.
Birnbaum
,
Acta Metall.
36
,
111
(
1988
).
17.
Y.
Oda
and
H.
Noguchi
,
Int. J. Fract.
132
,
99
(
2005
).
18.
M.
Yamaguchi
,
K.-I.
Ebihara
,
M.
Itakura
,
T.
Kadoyoshi
,
T.
Suzudo
, and
H.
Kaburaki
,
Metall. Mater. Trans. A
42
,
330
(
2011
).
19.
K. N.
Solanki
,
M. A.
Tschopp
,
M. A.
Bhatia
, and
N. R.
Rhodes
,
Metall. Mater. Trans. A
44
,
1365
(
2013
).
20.
M.
Rajagopalan
,
M. A.
Tschopp
, and
K. N.
Solanki
,
JOM
66
,
129
(
2014
).
21.
P.
Ferreira
,
I.
Robertson
, and
H.
Birnbaum
,
Acta Mater.
46
,
1749
(
1998
).
22.
S.
Wang
,
M. L.
Martin
,
P.
Sofronis
,
S.
Ohnuki
,
N.
Hashimoto
, and
I. M.
Robertson
,
Acta Mater.
69
,
275
(
2014
).
23.
P.
Sofronis
and
H. K.
Birnbaum
,
J. Mech. Phys. Solids
43
,
49
(
1995
).
24.
A.
Barnoush
and
H.
Vehoff
,
Acta Mater.
58
,
5274
(
2010
).
25.
A.
Barnoush
,
N.
Kheradmand
, and
T.
Hajilou
,
Scr. Mater.
108
,
76
(
2015
).
26.
K. N.
Solanki
,
D. K.
Ward
, and
D. J.
Bammann
,
Metall. Mater. Trans. A
42
,
340
(
2011
).
27.
I.
Adlakha
and
K. N.
Solanki
,
Proc. R. Soc. A
472
,
20150617
(
2016
).
28.
R.
Matsumoto
,
S.
Taketomi
,
S.
Matsumoto
, and
N.
Miyazaki
,
Int. J. Hydrogen Energy
34
,
9576
(
2009
).
29.
M.
Rajagopalan
,
M. A.
Bhatia
,
M. A.
Tschopp
,
D. J.
Srolovitz
, and
K. N.
Solanki
,
Acta Mater.
73
,
312
(
2014
).
30.
M. A.
Bhatia
,
S.
Groh
, and
K. N.
Solanki
,
J. Appl. Phys.
116
,
064302
(
2014
).
31.
M.
Yamaguchi
,
Metall. Mater. Trans. A
42
,
319
(
2011
).
32.
Y.
Zhao
and
G.
Lu
,
Modell. Simul. Mater. Sci. Eng.
19
,
065004
(
2011
).
33.
M.
Itakura
,
H.
Kaburaki
, and
M.
Yamaguchi
,
Acta Mater.
60
,
3698
(
2012
).
34.
I.
Adlakha
and
K. N.
Solanki
,
Acta Mater.
118
,
64
(
2016
).
35.
S.
Wang
,
M. L.
Martin
,
I. M.
Robertson
, and
P.
Sofronis
,
Acta Mater.
107
,
279
(
2016
).
36.
M.
Rajagopalan
,
I.
Adlakha
,
M. A.
Tschopp
, and
K. N.
Solanki
,
JOM
69
,
1398
(
2017
).
37.
S.
Taketomi
,
R.
Matsumoto
, and
N.
Miyazaki
,
Int. J. Mech. Sci.
52
,
334
(
2010
).
38.
H.
Vehoff
and
W.
Rothe
,
Acta Metall.
31
,
1781
(
1983
).
39.
T.
Tabata
and
H. K.
Birnbaum
,
Scr. Metall.
17
,
947
(
1983
).
40.
R. P.
Gangloff
,
Mater. Sci. Eng., A
103
,
157
(
1988
).
41.
Y.
Murakami
,
T.
Kanezaki
,
Y.
Mine
, and
S.
Matsuoka
,
Metall. Mater. Trans. A
39
,
1327
(
2008
).
42.
D. P.
Williams
and
H. G.
Nelson
,
Metall. Trans.
1
,
63
(
1970
).
43.
L.
Hagn
,
Mater. Sci. Eng., A
103
,
193
(
1988
).
44.
R.
Wang
,
Int. J. Fatigue
30
,
1376
(
2008
).
45.
S.
Plimpton
,
J. Comput. Phys.
117
,
1
(
1995
).
46.
M. S.
Daw
and
M. I.
Baskes
,
Phys. Rev. B
29
,
6443
(
1984
).
47.
J. A.
Zimmerman
,
H.
Gao
, and
F. F.
Abraham
,
Modell. Simul. Mater. Sci. Eng.
8
,
103
(
2000
).
48.
J. A.
Zimmerman
,
E. B.
Webb
 III
,
J. J.
Hoyt
,
R. E.
Jones
,
P. A.
Klein
, and
D. J.
Bammann
,
Modell. Simul. Mater. Sci. Eng.
12
,
S319
(
2004
).
49.
E.
Bitzek
,
J. R.
Kermode
, and
P.
Gumbsch
,
Int. J. Fract.
191
,
13
(
2015
).
50.
M. F.
Horstemeyer
,
Practical Aspects of Computational Chemistry
(
Springer
,
Dordrecht
,
2009
), pp.
87
135
.
51.
P.
Chowdhury
and
H.
Sehitoglu
,
Fatigue Fract. Eng. Mater. Struct.
39
,
652
(
2016
).
52.
I.
Adlakha
,
M. A.
Bhatia
,
M. A.
Tschopp
, and
K. N.
Solanki
,
Philos. Mag.
94
,
3445
(
2014
).
53.
S. M.
Ohr
,
S.‐J.
Chang
, and
R.
Thomson
,
J. Appl. Phys.
57
,
1839
(
1985
).
54.
G.
Schoeck
,
Philos. Mag. A
63
,
111
(
1991
).
55.
A.
Ramasubramaniam
,
M.
Itakura
, and
E. A.
Carter
,
Phys. Rev. B
79
,
174101
(
2009
).
56.
D. J.
Hepburn
and
G. J.
Ackland
,
Phys. Rev. B
78
,
165115
(
2008
).
57.
E.
Hayward
and
C.
Deo
,
J. Phys.: Condens. Matter
23
,
425402
(
2011
).
58.
S.
Taketomi
,
R.
Matsumoto
, and
N.
Miyazaki
,
Acta Mater.
56
,
3761
(
2008
).
59.
D.
Jiang
and
E. A.
Carter
,
Surf. Sci.
547
,
85
(
2003
).
60.
A.
Tehranchi
,
X.
Zhang
,
G.
Lu
, and
W.
Curtin
,
Modell. Simul. Mater. Sci. Eng.
25
,
025001
(
2017
).
61.
J. R.
Rice
,
J. Mech. Phys. Solids
40
,
239
(
1992
).
62.
I. R.
Vatne
,
A.
Stukowski
,
C.
Thaulow
,
E.
Østby
, and
J.
Marian
,
Mater. Sci. Eng., A
560
,
306
(
2013
).
63.
W.-S.
Ko
,
J. B.
Jeon
,
C.-H.
Lee
,
J.-K.
Lee
, and
B.-J.
Lee
,
Modell. Simul. Mater. Sci. Eng.
21
,
025012
(
2013
).
64.
R. A.
Oriani
,
Acta Metall.
18
,
147
(
1970
).
65.
H.
Vehoff
and
H.-K.
Klameth
,
Acta Metall.
33
,
955
(
1985
).
66.
E.
Taki
,
Y.
Kawakami
,
M.
Otsu
, and
K.
Takashima
,
J. Solid Mech. Mater. Eng.
1
,
779
(
2007
).
67.
M.
Tanaka
,
E.
Tarleton
, and
S. G.
Roberts
,
Acta Mater.
56
,
5123
(
2008
).
68.
J. J.
Möller
and
E.
Bitzek
,
Modell. Simul. Mater. Sci. Eng.
22
,
045002
(
2014
).
69.
D. E.
Jiang
and
E. A.
Carter
,
Acta Mater.
52
,
4801
(
2004
).
70.
Z.
Nishiyama
,
Sci. Rep. Tohoku Imp. Univ.
23
,
637
(
1934
).
71.
A.
Machov
and
G. E.
Beltz
,
Mater. Sci. Eng., A
387–389
,
414
(
2004
).
72.
A.
Stukowski
and
K.
Albe
,
Modell. Simul. Mater. Sci. Eng.
18
,
085001
(
2010
).
73.
R.
Peierls
,
Proc. Phys. Soc.
52
,
34
(
1940
).
74.
F. R. N.
Nabarro
,
Proc. Phys. Soc.
59
,
256
(
1947
).
75.
Y.
Sun
,
G. E.
Beltz
, and
J. R.
Rice
,
Mater. Sci. Eng., A
170
,
67
(
1993
).
76.
R.
Thomson
,
C.
Hsieh
, and
V.
Rana
,
J. Appl. Phys.
42
,
3154
(
1971
).
77.
J.
Knap
and
K.
Sieradzki
,
Phys. Rev. Lett.
82
,
1700
(
1999
).
78.
K.
Sieradzki
and
R. C.
Cammarata
,
Phys. Rev. Lett.
73
,
1049
(
1994
).
79.
Y.
Sun
and
G. E.
Beltz
,
J. Mech. Phys. Solids
42
,
1905
(
1994
).
80.
J.-S.
Wang
and
P. M.
Anderson
,
Acta Metall. Mater.
39
,
779
(
1991
).
81.
J.-S.
Wang
and
S. D.
Mesarovic
,
Acta Metall. Mater.
43
,
3837
(
1995
).
82.
C. S.
John
and
W. W.
Gerberich
,
Metall. Trans.
4
,
589
(
1973
).
83.
X.
Chen
,
T.
Foecke
,
M.
Lii
,
Y.
Katz
, and
W. W.
Gerberich
,
Eng. Fract. Mech.
35
,
997
(
1990
).
84.
L.-J.
Qiao
,
W.-Y.
Chu
, and
C.-M.
Hsiao
,
Scr. Metall.
21
,
7
(
1987
).
85.
M.
Takemoto
,
Corrosion
42
,
585
(
1986
).
86.
S.-H.
Chen
,
Y.
Katz
, and
W. W.
Gerberich
,
Scr. Metall. Mater.
24
,
1125
(
1990
).
87.
W.-Y.
Chu
,
C.-M.
Hsiao
,
S.-Y.
Ju
, and
C.
Wang
,
Corrosion
38
,
446
(
1982
).
88.
W.-Y.
Chu
,
C.-M.
Hsiao
, and
B.-J.
Xu
,
Metall. Trans. A
17
,
711
(
1986
).
89.
K. F.
McGuinn
and
M.
Aballe
,
Br. Corros. J.
17
,
18
(
1982
).
90.
T. L.
Anderson
and
T. L.
Anderson
,
Fracture Mechanics: Fundamentals and Applications
(
CRC Press
,
2005
).