Discrete dislocation plasticity HELPs understand hydrogen effects in bcc materials

In an attempt to bridge the gap between atomistic and continuum plasticity simulations of hydrogen in iron, we present three dimensional discrete dislocation plasticity simulations incorporating the hydrogen elastic stress and a hydrogen dependent dislocation mobility law. The hydrogen induced stress is incorporated following the formulation derived by Gu and El-Awady (2018) which here we extend to a finite boundary value problem, a microcantilever beam, via the superposition principle. The hydrogen dependent mobility law is based on first principle calculations by Katzarov et al. (2017) and was found to promote dislocation generation and enhance slip planarity at a bulk hydrogen concentration of 0.1 appm; which is typical for bcc materials. The hydrogen elastic stress produced the same behaviour, but only when the bulk concentration was extremely high. In a microcantilever, hydrogen was found to promote dislocation activity which lowered the flow stress and generated more pronounced slip steps on the free surfaces. These observations are consistent with the hydrogen enhanced localized plasticity (HELP) mechanism, and it is concluded that both the hydrogen elastic stress and hydrogen increased dislocation mobility are viable explanations for HELP. However it is the latter that dominates at the low concentrations typically found in bcc metals.


Introduction
Hydrogen can cause premature failure accompanied by a decrease in ductility in metals, which is often referred to as hydrogen embrittlement (HE) ( Ayas et al., 2014 ). The underlying mechanism of this phenomenon has been under debate for the past two decades and several theories have been proposed ( Robertson et al., 2015 ). For non hydride forming material, the prevailing mechanisms are hydrogen enhanced decohesion (HEDE) and hydrogen enhanced localized plasticity (HELP). The HEDE mechanism assumes that dissolved hydrogen reduces the cohesive strength of the material lattice ( Oriani, 1972 ). It is supported by atomistic simulations where cohesive strength decreases with increasing hydrogen content ( Lynch, 2012 ). The HELP mechanism assumes that hydrogen facilitates dislocation activity thereby enhancing plasticity locally ( Birnbaum and Sofronis, 1994 ), which is backed up by experimental observations and theoretical calculations ( Robertson, 2001 ). In recent years, it has become widely accepted that HE is a complex, material and environmental dependent process ( Robertson et al., 2015 ) so that no mechanism applies exclusively. In fact, it has been indicated that the aforementioned mechanisms may act in concert to promote failure ( Novak et al., 2010 ). Barrera et al. (2014) performed continuum level HE simulation considering both mechanisms in welded structures, where microcracks formed at the matrix/carbide interface due to decohesion process followed by localization of plastic flow responsible for the link-up of the microcracks. Most recently, a novel hydrogenenhanced-plasticity mediated decohesion mechanism was proposed for HE in martensitic steels ( Nagao et al., 2018;2012 ). It claims that hydrogen-enhanced mobility of dislocations leads to local stress build-up and hydrogen concentration elevation close to material interfaces where decohesion is promoted. In this mechanism, hydrogen enhanced plasticity is the driving force for HE, which highlights the importance of understanding how hydrogen affects dislocations.

Experimental evidence
Experimental evidence for the HELP mechanism is found in the literature over different length scales. In situ deformation experiments using an environmental cell inside a transmission electron microscope (TEM) ( Ferreira et al., 1998 ) enables direct observation of hydrogen-dislocation interactions involving several dislocation lines. With this technique, Shih et al. (1988) observed hydrogen enhanced dislocation motion in α-Titanium, which was eliminated and recovered when the hydrogen environment was removed and reintroduced. A semi-quantitative relationship between dislocation velocity and hydrogen gas pressure was also extracted. Similar experiments were performed on a variety of materials, including Fe ( Tabata and Birnbaum, 1983 ), Ni ( Robertson and Birnbaum, 1986 ), Al ( Bond et al., 1988 ) and their alloys ( Bond et al., 1989;Rozenak et al., 1990 ). A surprising feature of these studies is that the enhancement in the mobility of dislocations due to hydrogen is independent of crystal structure and is the same for edge, screw, mixed, and partial dislocations ( Robertson, 2001 ). Using a similar technique, Ferreira et al. (1998) observed hydrogen decreased the equilibrium separation in a dislocation pileup in 310S stainless steel and Al; indicating hydrogen has a shielding effect on dislocation interactions. Ferreira et al. (1999) also observed hydrogen suppressed dislocation cross-slip in Al and enhanced slip planarity.
Plasticity is an emergent property due to the collective behavior of a large number of dislocations. Therefore it cannot be fully understood by only studying a small number. Small scale tests such as nanoindentation ( Barnoush and Vehoff, 2010;Nibur et al., 2006 ) and microcantilever bending ( Deng and Barnoush, 2018;Deng et al., 2017;Hajilou et al., 2017 ) in hydrogen environments provide essential new insight. In these experiments, the global load-displacement curve is recorded, as well as dislocation pileups and plastic localisation by microscopy techniques. Nibur et al. (2006) studied hydrogen effects on dislocation activity in austenitic stainless steel with nanoindentation. Hydrogen was observed to facilitate dislocation nucleation and enhance slip planarity. Hajilou et al. (2017) performed in situ bend tests on notched microcantilevers of single crystal Fe-3 wt% Si alloy and observed that hydrogen reduced flow stress and promoted crack propagation. Deng et al. (2017) performed similar experiment on un-notched microcantilevers of single crystalline FeAl and observed that hydrogen reduced the flow stress and enhanced slip localisation. Recently, Deng and Barnoush (2018) investigated hydrogen assisted fracture in single crystalline FeAl notched microcantilevers, which demonstrated that hydrogen leads to global softening and reduced fracture toughness. It should be noted that the decrease in flow stress in these experiments was explained with hydrogen enhanced dislocation generation and the accelerated crack propagation with hydrogen reduced dislocation mobility. Ex situ methods where specimens are loaded in a hydrogen environment and observed at specific locations are also widely adopted, and allow tests to be performed at a larger scale. Four point bend tests on lath martensitic steel were carried out ( Nagao et al., 2012 ), after which the fracture surface was examined with scanning electron microscopy (SEM) and the microstructure immediately beneath the prominent features on the fracture surface was determined by focused-ion beam (FIB) machining and transmission electron microscopy (TEM). Significant plasticity, consistent with the HELP mechanism, was observed beneath hydrogen induced "quasi-cleavage "fracture surface, which inspired the hydrogen-enhanced-plasticity mediated decohesion mechanism ( Nagao et al., 2018 ). A wealth of experimental evidence on the macroscopic scale is also found in the literature, for instance, intensive localized shear was observed on fracture surfaces in hydrogen charged AISI 310S stainless steel ( Ulmer and Altstetter, 1991 ), and hydrogen was observed to change the microvoid failure mode from inter-ligament necking to shear banding in slow strain rate tensile test ( Matsuo et al., 2014 ).

Previous modeling and simulation
The experimental evidence for the HELP mechanism demonstrates two different hydrogen effects; decreasing the equilibrium dislocation spacing while simultaneously increasing the dislocation mobility. To some extent, this can be accounted for by the elastic shielding theory ( Birnbaum and Sofronis, 1994 ). In this theory, hydrogen forms an atmosphere around the dislocations. This generates a first order elastic interaction arising from introducing hydrogen atoms within the dislocation stress field and a second-order interaction resulting from a local change in the elastic modulus caused by solute hydrogen. It was shown by Birnbaum and Sofronis (1994) that these elastic interactions led to a short-range decrease in the elastic force experienced by a dislocation due to another dislocation or an obstacle. In this way, hydrogen decreases dislocation spacing by shielding dislocation interactions and enhances dislocation mobility by shielding elastic stress fields of other obstacles. At very high background concentrations such as c 0 = 1 × 10 4 appm (equivalent with 1 hydrogen atom per 100 solvent atoms) the hydrogen elastic shielding effect can be considerable. Within the same theoretical framework, the experimentally observed enhanced slip planarity was explained by Ferreira et al. (1999) . In fcc Al, the first order elastic interaction dominates, leading to the formation of hydrogen atmospheres primarily around edge components. In this way, the self-energy of an edge dislocation-hydrogen complex is lowered while that of the screw component is less affected. This makes the edge components more favorable when compared to the hydrogen free case, therefore, the proportion of screw component will be decreased and the tendency to cross-slip suppressed with hydrogen. This mechanism is shown to be noticeable with a hydrogen concentration of c 0 > 1 × 10 4 appm. Strictly speaking, the elastic shielding theory is based on an implicit assumption that hydrogen diffuses fast enough that the hydrogen distribution is always in equilibrium, in terms of chemical potential, with a remote reservoir with a uniform background concentration c 0 . In reality, c 0 should be taken as the bulk concentration ( Jagodzinski et al., 20 0 0;Wen et al., 2003 ). It is established that hydrogen diffusivity in bcc metals can be orders of magnitude higher than in fcc metals at ambient temperatures ( Lynch, 2012 ), therefore, this theory should be more applicable to bcc metals ( Gu and El-Awady, 2018 ). However, the low hydrogen solubility in bcc metals gives a small bulk concentration c 0 < 1.0 appm, at which hydrogen elastic shielding is negligible. Even for fcc metals, the concentration used in Birnbaum and Sofronis (1994) is too high. This implies that the hydrogen elastic shielding mechanism, while being viable, might not be the dominant factor in HELP phenomena, as noted by other researchers ( Gavriljuk et al., 2010;Taketomi et al., 2008;Teus et al., 2007 ) and will be revisited in this paper.
Atomistic and quantum mechanics based simulations are increasingly utilized to simulate hydrogen-dislocation interactions, providing fundamental insights and sometimes controversial evidence for the HELP mechanism. Using molecular dynamics, Curtin (2013, 2011) observed hydrogen suppressed dislocation emission from a crack tip which facilitated brittle-cleavage failure in α-Fe. They also studied the effect of hydrogen on a pure edge dislocation in α-Fe and concluded that hydrogen reduces dislocation mobility over a wide range of concentrations ( Song and Curtin, 2014 ), consistent with the solute drag theory and in contradiction with the HELP mechanism. A similar scenario was observed by Bhatia et al. (2014) where hydrogen impeded the motion of pure edge dislocations in α-Fe. According to Song and Curtin (2014) , these results are also applicable to fcc metals due to certain similarities of edge dislocations between the two systems. However, the story is quite different when it comes to screw dislocations. In bcc metals, the mobility of screw dislocations is much lower than for edge. Plasticity is therefore mainly limited by screw dislocations in bcc metals whose mobility is dominated by kink-pair nucleation and migration ( Anderson et al., 2017 ). Using a line tension model, Itakura et al. (2013) studied the effect of hydrogen on screw dislocation mobility in α-Fe using first-principles calculations and revealed a general trend that hydrogen enhances screw mobility by promoting nucleation of kink pairs at a low concentration and decreases the mobility by impeding migration of kink pairs at a high concentration. Based on this work and by performing kinetic Monte Carlo (kMC) simulations, Katzarov et al. (2017) were able to simulate the motion of individual 1/2[111] screw dislocations over long timescales. They provided the first quantitative description of hydrogen effects on screw mobility as a function of bulk hydrogen concentration, temperature and applied stress in α-Fe. According to this study, hydrogen significantly enhances screw mobility up to a bulk concentration c 0 = 5 appm at room temperature, which indicates the HELP mechanism takes place in α-Fe given the small solubility of hydrogen. In their multiscale quantummechanics/molecular mechanics simulation, Zhao and Lu (2011) observed hydrogen lowered the Peierls energy barrier for screw dislocations significantly, indicating enhanced mobility. Comparing this conclusion to their observation of hydrogen decreased mobility in edge dislocations, Bhatia et al. (2014) proposed that hydrogen has a tendency to make the crystal more isotropic in terms of dislocation motion, i.e. increasing screw mobility while decreasing edge mobility ( Itakura et al., 2013 ), which was also hypothesized by Moriya et al. (1979) in an experimental investigation.
Continuum level simulations to investigate HELP have been performed and provided insight into hydrogen induced failure at the macroscopic scale. Sofronis et al. (2001) and Liang et al. (2003) showed that hydrogen could induce shear instability in a plain strain tensile specimen. Ahn et al. (2007) observed hydrogen enhanced void growth and coalescence using a representative material volume approach. Barrera et al. (2016) studied hydrogen effects on plastic behavior of notched tensile plates using a coupled mechanical-diffusion analysis. A similarity of these studies was that the HELP mechanism was incorporated by describing the matrix yield stress as a linearly decreasing function of local hydrogen concentration. Most recently, Yu et al. (2018) adopted a sigmoidal relationship for the hydrogen softening effect to study hydrogen-microvoid interactions and concluded that the results were more comparable to certain experimental observations. Nevertheless, both the linear relation and the sigmoidal relation are attempts to describe the HELP mechanism phenomenologically, due to the lack of quantification of the hydrogen softening effect applicable to the continuum scale, either experimentally or numerically. Although many clues about the physics of hydrogen-dislocation interactions have been revealed by atomistic simulations, they cannot be directly utilized in continuum level studies due to the large gap in spatial and temporal scales.
To bridge this gap, discrete dislocation dynamics (DDD) simulations have proven to be an effective tool. DDD can simulate the aggregate behavior of large dislocation ensembles by simulating the dynamics and interactions of every discrete dislocation segment. DDD has been used to study a wide range of problems assuming isotropic elasticity such as the formation and strength of dislocation junctions ( Capolungo, 2011;Wu et al., 2013 ) and work hardening ( Arsenlis et al., 2007 ). A small number of studies have also been performed using anisotropic elasticity to simulate a Frank-Read source ( Fitzgerald et al., 2012 ) and dislocation loops ( Aubry et al., 2011 ). DDD can be coupled with FEM to simulate a finite domain with mixed boundary conditions; this is usually referred to as discrete dislocation plasticity or DDP. This approach can help understand plasticity due to the collective behavior of a large number of dislocations in a confined volume under various loading scenarios, such as in a micro-beam under bending ( Motz et al., 2008 ) or a micropillar under tension ( Motz et al., 2009 ) or compression ( Ryu et al., 2013 ). In these simulations, the dynamics of individual dislocation segments and their interactions are described in a manner consistent with the physics at the atomistic scale, and the results reflect crystal level plasticity. Within the DDP framework, environmental effects on plasticity can be studied by properly considering the effect of environment on dislocation dynamics or dislocation-obstacle interactions. By incorporating a temperature dependent dislocation mobility law informed by atomistic simulation , Cui et al. (2016) studied the effect of temperature on plasticity in bcc micropillar crystals. Cui et al. (2017) also studied the effect of irradiation on strain bursts in fcc and bcc single crystals by considering the irradiation effects on initial dislocation structures.
Hydrogen was recently incorporated in three dimensional DDD by Gu and El-Awady (2018) , which is a promising method for bridging the gap in the multi-scale simulation of the HELP mechanism. They simulated the role of hydrogen during the shrinking of a dislocation loop and in the arrangement of parallel dislocations. Hydrogen atoms were regarded as point defects, misfitting inclusions in an elastic volume, subjected to the externally applied stress and the dislocation stress fields. The elastic stress field due to hydrogen was derived in three dimensional space, following the approach proposed by Cai et al. (2014) . In this way, the hydrogen elastic energy contribution to dislocation activity was considered. In the case with a high hydrogen diffusivity, e.g. in bcc materials, it was shown that hydrogen tended to compress the dislocation loop in the directions parallel to the Burgers vector and to decrease the spacing of parallel edge dislocations, consistent with hydrogen elastic shielding theory. The work of Gu and El-Awady (2018) is therefore the first three dimensional DDD illustration of the hydrogen elastic shielding mechanism in an infinite medium. DDP simulations of the effect of hydrogen on dislocation plasticity in a confined volume has not been reported in the literature previously. Further, it is noted that the influence of hydrogen in the previous DDD simulations was noticeable only at very high bulk concentrations c 0 ≈ 5 × 10 4 appm and c 0 ≈ 1 × 10 5 appm. This is a consequence of only considering hydrogen elastic shielding and therefore indicates that other hydrogen effects need to be taken into account in order to reproduce the influence of hydrogen on plasticity at the low concentrations observed in bcc materials.

The present work
In this paper, we present three dimensional discrete dislocation plasticity (DDP) simulations incorporating hydrogen in a finite volume with dislocation mobility parameters typical for bcc metals and a physical bulk hydrogen concentration. In addition to hydrogen elastic shielding, hydrogen increased dislocation mobility is also incorporated following the quantitative description introduced by Katzarov et al. (2017) . The simulations were performed using an in-house Matlab based GPU accelerated DDP code; which is a modified version of DDLab developed at LLNL by Arsenlis et al. (2007) .
The micro-cantilever geometry, popular for experimental studies on hydrogen embrittlement, is simulated using the superposition method ( El-Awady et al., 2008;Giessen and Needleman, 1999 ). In contrast to Gu and El-Awady (2018) where dislocation mobilities for edge and screw components were assumed equal, here we assign a higher mobility to edge components than screw, which is more representative of bcc metals.

Peach-Koehler force
The DDD simulations in this work are nodal based ( Arsenlis et al., 2007;Bulatov and Cai, 2006 ), i.e. dislocations are represented by discrete straight segments as shown in Fig. 1 (a), and the dislocation motion is accounted for by evaluating the velocity V k of node k at position X k through a general linear mobility law which describes the dislocation mobility with various modes such as glide, climb and cross-slip.
In general, the force per unit length on a dislocation line segment bounded by nodes k and l is where σ( x ) is the local stress at point x on the dislocation line, b kl is the Burgers vector, and l kl is the unit tangent vector of the dislocation line segment. In the DDD simulation, the non-singular continuum theory of dislocations  is adopted. Using the convention that nodal values are denoted in uppercase, F k is the total nodal force on node k , which is evaluated in three parts where ˜ f i j kl (X k ) is the elastic force due to segment i → j integrated along segment k → l evaluated at node k . This is summed over all segments i → j inside the domain, including the self force due to segment k → l , this is then summed over all nodes l which are connected to node k . ˆ f kl (X k ) is the corrective elastic force on node k integrated along segment k → l obtained using FEM to account for the applied and image stresses. Similarly f c kl (X k ) is the force arising from the dislocation core energy of segment k → l . The PK force due to dislocations arises from the stress fields of straight dislocation segments. Take the configuration in Fig. 1 (a) for instance, the stress tensor at a field point x induced by the dislocation segment 1 → 2 can be expressed as ( Arsenlis et al., 2007 ) ˜ σ 12 where μ is the shear modulus, ν the Poisson's ratio, a the core width and I the second order identity tensor.
and R a = R · R + a 2 . The nodal force at node 3 due to the interaction between segment 4 → 3 with Burgers vector b 43 , and the segment 1 → 2 with stress field ˜ σ 12 is Analytic expressions for Eq. (13) have been given by Arsenlis et al. (2007) .

Dislocation mobility law
A linear dislocation mobility law is adopted. For each segment kl , a drag tensor B kl is determined according to the segment character, the nodal velocity V k at node k is then obtained using where the sum is again over all nodes l connected to node k , and F k is the nodal force determined by Eq.
(2) . To construct a drag tensor suitable for bcc materials , the cases of pure screw and pure edge components are considered first. A pure screw dislocation segment is assigned with isotropic mobility in all directions perpendicular to its line direction, which mimics the so-called "pencil-glide" where B s is a drag coefficient for the motion of screw dislocations and B l B s to allow a node to move rapidly along the line direction. Meanwhile, the mobility of a pure edge segment is anisotropic with respect to glide and climb, where B eg and B ec are the drag coefficients for glide and climb, respectively. The unit vectors are the plane normal n kl = To account for the mobility of mixed segments, an interpolation function was proposed by Cai and Bulatov (2004) , enabling a smooth transition between the two limits The dislocation mobility is inversely proportional to the drag coefficient. In bcc materials, the mobility of a pure edge segment should be greater than a pure screw segment, which is accounted for by assigning B eg < B s . Climb is neglected throughout this work by using a B ec B s . More details regarding the selection of drag coefficients are presented in Section 3.1 . Assuming "pencil glide"for pure screw segments is only valid in bcc metals at high temperatures. Such treatment might exaggerate the cross slip of screw dislocations and a more sophisticated cross-slip law will be implemented in future work. In this paper, however, it is noted that the influence of hydrogen on slip planarity is determined by comparing the proportion of near-screw segments with and without hydrogen; which should remain qualitatively unchanged with a more sophisticated cross-slip law.

Treatment of finite boundary conditions
Dislocations in a finite medium are subjected to the corrective force term ˆ f which appears on the right hand side of Eq. (2) . To evaluate the corrective fields the superposition principle ( Giessen and Needleman, 1995;Weygand et al., 2002 ) is adopted. Denoting the infinite-body fields as ˜ () and the finite-element correction fields with ˆ () , the total stress, strain and displacement fields are expressed as respectively. The following procedure ( El-Awady et al., 2008 ) is adopted to evaluate the image stress field. First, the elastic stress field due to dislocations in an infinite body, ˜ σ, is obtained using Eq. (3) ; the tractions ˜ T = ˜ σ · n on the traction boundaries due to this stress are then calculated, reversed and added to the prescribed traction boundary conditions T ; these modified boundary values, ˆ T = T −˜ T , in addition to any prescribed displacement conditions U on the displacement boundaries are used in an elastic finite element simulation to determine the corrective fields. Finally, the corrective stress field, ˆ σ, is used to evaluate the corrective nodal force ˆ f .

Hydrogen elastic shielding
The hydrogen stress contribution required for three dimensional DDD was recently obtained by Gu and El-Awady (2018) using isotropic linear elasticity theory, based on the work of Cai et al. (2014) . Hydrogen atoms were viewed as point defects of Eshelby type ( Eshelby, 1955 ), and it was assumed that hydrogen diffusion is sufficiently rapid that hydrogen atoms can move with the dislocations. Therefore equilibrium of chemical potential is guaranteed over the domain. This was referred to as the "fast diffusion" scenario which is suitable for bcc materials. Under this assumption, steady-state hydrogen redistribution is established at all times, and is only a function of the dislocation pressure field ˜ where 0 ≤ χ( x ) ≤ 1 is the hydrogen field distribution expressed as the occupancy. Occupancy is related to hydrogen con- where c max is the maximum concentration or equivalently, the number of atomic sites per unit volume where each site can accomodate one hydrogen atom. The volume expansion of the matrix due to a hydrogen atom is V , and χ 0 is the reference occupancy under a uniform pressure P with μ 0 being the reference chemical potential of hydrogen in a stress-free infinite medium. In reality, μ 0 could represent the chemical potential of hydrogen in a remote reservoir and should be viewed as the initial hydrogen potential in the bulk. The hydrogen population as distributed point defects will generate an elastic stress field which is dependent on the spatial distribution of hydrogen occupancy; which is a function of stress. Therefore the hydrogen induced stress is dependent on the dislocation configuration in the domain, as indicated in Eq. (10) . According to Gu and El-Awady (2018) , the hydrogen stress field in three dimensions can be expressed as This is proportional to the last three terms (the none screw terms) of the dislocation stress field in Eq. (3) , indicating that the hydrogen elastic stress field can be evaluated via a similar approach to evaluating the stress generated by a dislocation segment. Following Gu and El-Awady (2018) , the hydrogen induced PK force f is based on virtual hydrogen-related dislocation segments as illustrated in Fig. 1 (b). Due to the fast diffusion assumption for bcc materials, the hydrogen-related dislocation segments are simply the dislocation segments, and the hydrogen induced force at a node at X 3 on the dislocation segment 4 → 3 arising from the hydrogen stress field σ 12 is calculated using and this term needs to be added to the right hand side of Eq. (2) . In a finite volume, the hydrogen stress field σ generates tractions on the surfaces of the body. This will cause the body to deform slightly inducing an additional image stress which will contribute to the force on the dislocation nodes. This corrective stress field ˆ σ which includes the additional hydrogen image stress is obtained using FEM with the modified traction boundary condition ˆ T = T −˜ T −T . The total elastic fields are then the superposition of the corrective FEM fields, ˆ (. ) , the analytic dislocation fields ˜ (. ) and analytic hydrogen fields ( . ) After obtaining ˆ σ and consequently ˆ f , the corrective force which includes the hydrogen image force, the total nodal force on node k can be evaluated where f i j kl is the elastic force due to the hydrogen segment i → j integrated along segment k → l . Finally, it should be noted that the Eshelby inclusion assumption is a simplified treatment which ignores certain aspects (modulus, electronic and magnetic effects) associated with a real hydrogen atom embedded in a material lattice. Nevertheless, such an assumption is widely adopted due to the advantage that explicit analytic solutions exist for the stress/strain fields and elastic energy, which can be readily implemented in the DDD framework ( Cai et al., 2014 ).

Hydrogen dependent mobility law
The influence of hydrogen on dislocation mobility can be directly captured by atomistic simulations. Katzarov et al. (2017) calibrated the effect of hydrogen on the mobility of pure screw dislocations in α-Fe, using kMC based on first principle calculations ( Itakura et al., 2013 ). A quantitative relationship between dislocation mobility and bulk hydrogen concentration were extracted under different temperatures and applied stresses. This work concerns the ambient temperature case, and the corresponding curves are reproduced with permission from Katzarov et al. (2017) in Fig. 2 (a). The hydrogen effects are evaluated under three applied stress levels, and hydrogen enhanced dislocation mobility was observed up to 10 appm in all the cases.
Hydrogen effects on dislocation mobility are incorporated into the DD simulation based on Fig. 2 (a). In order to maintain a linear mobility law, the drag coefficients are assumed to be dependent on hydrogen but independent of applied stress. The curve corresponding to σ app = 100 MPa is adopted throughout in order to produce the highest level of hydrogen sensitivity.
Further, considering the small bulk hydrogen concentration in bcc materials, only the initial part of the curve 0.0 ≤ c H ≤ 5.0 appm was fitted. Plotting this part in Fig. 2 (b), the following linear fitting was obtained, According to this equation, the increase in dislocation mobility due to hydrogen can be as large as 17 times with c H ≥ 3 appm. In this work, however, the maximum increase achieved is less than 10 due to the small initial bulk concentration used.
In DD simulations, the local hydrogen concentration c i j H at each segment is evaluated during each iteration. The drag coefficient for each segment is scaled using Eq. (16) and the nodal velocity is then found using Eq. (5) . The scalar drag coefficients used to construct the drag tensor Eq. (8) for each segment were reduced using where B is the hydrogen free drag coefficient, either B s or B eg . Note that B ec was unchanged to ensure climb did not occur. The local hydrogen concentration c H ( x ) at a field point x is obtained using Eq. (10) . However, the normal components of the non-singular stress field of a straight dislocation segment are all zero on the dislocation line and so Tr (σ) = 0 when the field point is on a dislocation segment. Therefore, in order to evaluate the concentration on a dislocation segment without neglecting the self stress of the segment, the concentration was found at 2 field points at ± a n , either side of the segment midpoint, and then averaged. Where a is the core width in Eq.
(3) and n the glide plane normal of the segment. The kMC investigation into hydrogen effects on screw mobility in Katzarov et al. (2017) was performed considering hydrogen concentration near the dislocation core which was well captured by first principle calculations ( Itakura et al., 2013 ).
Starting from a sensible bulk hydrogen concentration, e.g. c H = 10 appm a very high hydrogen occupancy up to χ = 0 . 05 (equivalent to c = 5 × 10 3 c H ) could build up in the core region due to its high binding energy. Such a high concentration is impossible to achieve within linear elasticity theory. This implies the hydrogen increased dislocation mobility could be viewed as the total effect including the contribution of hydrogen in the core region, whereas hydrogen elastic shielding considers only the weak hydrogen accumulation generated by the remote dislocation elastic stress fields. The local hydrogen concentration obtained in the current DDD simulation is associated with the linear elastic stress field of dislocations which is comparable to the bulk concentration c H used in atomistic simulations. In Katzarov et al. (2017) , a single dislocation line was simulated and the simulation size was 10 0 0 b = 0 . 245 μm . This is comparable to the finite element size used in this work 0.24 μ m i.e. the resolution of the corrective image stress and consequently the local hydrogen concentration is close to the simulation size in Katzarov et al. (2017) . Therefore the local hydrogen concentration obtained in the current DDD simulation should be considered equivalent to the bulk concentration c H in Katzarov et al. (2017) .
Exactly how hydrogen influences edge mobility is currently unclear, since no quantitative description is available and there is no consensus in the literature. Robertson (2001) concluded that hydrogen increases dislocation mobility equally for edge and screw dislocations, while atomistic simulations predict hydrogen reduces edge mobility ( Bhatia et al., 2014;Song and Curtin, 2014 ) and increases screw mobility ( Itakura et al., 2013;Katzarov et al., 2017 ). Experimentally, Xie et al. (2016) observed hydrogen decreased edge dislocation mobility in Al during in situ micropillar cyclic loading tests, and the same author ( Xie, 2018 ) most recently observed hydrogen enhanced screw mobility in Fe using a similar technique. Therefore we tested two cases. Case I assumes hydrogen only increases screw mobility, and case II assumes an equal influence on screw and edge mobility

Geometry for hydrogen informed DDD simulation
With hydrogen informed DDD simulation, the influence of hydrogen on dislocation generation, pile-up and slip planarity is investigated starting from a FR source in an infinite solid subjected to both zero and non-zero out-of-plane PK forces. The findings here will provide rationalization for the phenomena observed in microcantilever bending simulations in Section 4 . In all the simulations, two different hydrogen effects, hydrogen elastic shielding and hydrogen increased dislocation mobility, Fig. 3. Illustration of a pure edge dislocation belonging to the (101)[11 1 ] slip system in local crystallographic coordinate system. The edge dislocation is pinned at both ends, and its line direction is l l = (−1 , 2 , 1) . The vectors in this local coordinate system are marked with a subscript l . A pure edge dislocation belonging to the (101)[11 1 ] slip system is selected as the initial configuration for DDD simulation in an infinite medium. The dislocation lies on the (101) plane and is pinned at both ends, as shown in Fig. 3 . Its behavior is dependent on the applied stress σ l in the local coordinate system. For easier illustration, the model is rotated into a global system, with axes along b l , l l and n l , as shown in Figs. 4 and 5 . The stress tensors in the global coordinate system are readily obtained where R is the rotation matrix from the local (crystal) to the global (dislocation) frame In Fig. 4 , out-of-plane pure shear is applied The PK forces on all dislocation components are in plane, therefore, the initial edge dislocation behaves as an FR source. In Fig. 5 , uniaxial tension along the crystallographic [100] direction in the local coordinate system is applied, which in the global system is expressed as At the beginning, the pure edge dislocation bows out to the left driven by the PK force with a component in the negative x direction. Although this PK force possesses a non-zero component in the positive z direction, it does not cause any out-ofplane motion of the edge segment, since climb is not permitted. During the bowing out, pure screw segments are generated, which move in the y direction and cross slip. Therefore, when the dislocation bows back, the screw segments are unable to annihilate each other unless they first cross slip back onto the same plane. The simple model illustrated in Fig. 4 is used to investigate hydrogen effects on dislocation generation and pile-up, and that in Fig. 5 is used to study hydrogen effects on slip planarity, with the results elaborated in Section 3.2 and 3.3 , respectively. Adopting the material parameters (relevant to Ni), initial loop size and hydrogen concentration values used by Gu and El-Awady (2018) , the same results in terms of PK force magnitude and loop shape change were obtained, which are not repeated here. With the current implementation validated, material parameters relevant to Fe, were used to simulate an FR source. The Burgers vector magnitude is | b | = √ 3 2 a 0 , where a 0 = 2 . 856 Å is the lattice parameter; isotropic linear elasticity is assumed, with a shear modulus μ = 83 GPa and Poisson's ratio ν = 0 . 29 ( Tang and Marian, 2014 ). In the absence of hydrogen, the edge and screw drag coefficients are B eg = 5 × 10 −4 Pa s and B s = 1 × 10 −2 Pa s as used by Wang and Beyerlein (2011) . Climb is disabled in the current work, therefore, only pure screw segments can possibly move out of the glide plane, and a dislocation segment is considered pure screw if the angle θ between its line direction and Burgers vector is θ < 0.01 °. The dislocation core width is selected as a = 10b , as used by Arsenlis et al. (2007) and Ayas et al. (2014) . For DDD simulation in an infinite medium ( Figs. 4 and 5 ), the initial dislocation length is l 0 = 500b in both cases. The applied shear stress in Eq. (21) is τ = 450 MPa, and the uniaxial tension stress is σ = 800 MPa in Eq. (22) . All simulations are for T = 300 K, and the fast diffusion scenario ( Gu and El-Awady, 2018 ) is assumed, given the large diffusivity of hydrogen in bcc materials at this temperature. As for the initial bulk hydrogen concentration, an extremely large value of c 0 = 5 × 10 5 appm is selected to produce a noticeable effect with only hydrogen elastic shielding taken into account; a realistic value of c 0 = 0 . 1 appm is then adopted in the simulations considering hydrogen increased dislocation mobility.

Hydrogen elastic shielding
For Fe, the maximum hydrogen concentration c max in Eq. (10) is ≈ 30% lower than Ni. The lower mobility of screw components in Fe results in an elongated dislocation loop, as shown in Fig. 6 , instead of a circular one in Ni. As the density of screw segments (which have no elastic interaction with hydrogen) is higher in Fe than in Ni the effect of hydrogen shielding is reduced for the same occupancy. Therefore, a higher hydrogen occupancy was required than by Gu and El-Awady (2018) , in order to produce noticeable hydrogen influence.
Starting from an extremely high bulk concentration of c 0 = 5 × 10 5 appm, hydrogen elastic shielding influences the behavior of a FR source, as shown in Fig. 6 . With hydrogen present the dislocation loops generated by the FR source move further, and one extra dislocation loop is generated, indicating that hydrogen shielding can promote dislocation nucleation. In addition, the average spacing between the dislocation loops is slightly decreased due to hydrogen, as also observed by Gu and El-Awady (2018) .
The nodal force contribution due to the dislocation segments and the uniform applied stress ˜ F k + ˆ F k and due to the hydrogen stress F k (scaled × 10) are plotted shortly after the first bow out of the FR source in Fig. 7 (a). The direction of ˜ F k + ˆ F tends to expand the half loop and cause it to pinch off to generate a full loop. The distribution of F k is highly dependent on the character of the segments connected to the node. In general, the values are always much less than the other elastic forces and are largest on the nodes connected to edge dominated segments. F k on nodes on long screw segments are negligible. This is expected as screw dislocations exert zero hydrostatic stress there is no elastic interaction with the hydrogen, e.g. using Eq. (12) shows σ kl = 0 for a pure screw segment. The hydrogen force, F k , along pure edge segments or mixed segments which are edge dominant opposes the other elastic forces, ˜ F k + ˆ F , whereas the directions are the same for mixed segments which are screw dominant, as shown in Fig. 7 (b). This indicates that hydrogen has a tendency to slow down the loop expansion along b . As the anisotropy in the mobility is slightly reduced, this leads to a slightly more circular shape of the half loop. In addition, Fig. 7 (c) shows how the hydrogen force promotes the attraction of the short screw dominant segments. Making it easier for the segments to bow back and pinch off to nucleate a loop. So far, hydrogen elastic shielding has been demonstrated to have a tendency to enhance dislocation nucleation and to decrease the average spacing. However, it should be noted that to achieve this required an extremely high bulk hydrogen concentration, e.g. c 0 = 5 × 10 5 appm. With realistic bulk concentrations for bcc Fe, the elastic hydrogen effect becomes negligible. This was verified using c 0 = 0 . 1 appm. In this case, the FR source behaviour was identical to the hydrogen-free case in Fig. 6 and so is not repeated.

Hydrogen increased mobility
We now consider the effect of hydrogen increased mobility at a realistic bulk concentration of c 0 = 0 . 1 appm. Both case I and case II , as outlined in Eq. (18) , are investigated. The dislocation structures with and without hydrogen after the same simulated time are superimposed in Fig. 8 . In both cases, more loops are generated, in the presence of hydrogen. The total dislocation line length increases with the hydrogen enhanced mobility as shown in Fig. 9 (a). In case I , where hydrogen only enhances screw mobility, the shape of the loop differs from the hydrogen free situation. Whereas in case II , where hydrogen exerts equal influence on both screw and edge mobilities, the shape of the loops is unchanged. Fig. 9 (b) shows how the proportion of screw dominant segments evolves with simulated time. Here, a segment is considered screw dominant if the angle between l and b is less than 5 °. It is noted this criterion is only used during postprocessing to identify screw like  segments. The proportion of screw like segments is reduced in case I but unchanged in case II . This same trend can also be reproduced by ignoring the spatial distribution of hydrogen completely and uniformly increasing the screw mobility (a simplified case I) or uniformly increasing the mobility of all segments (a simplified case II). The same dislocation structure is predicted in case II as without hydrogen, but at an earlier simulated time.
In case I , hydrogen promotes dislocation generation via increased screw mobility and the shape of loops becomes more circular; in case II , hydrogen accelerates dislocation evolution but has little influence on the final dislocation structure in quasi static loading. Case I reduces the anisotropy in the dislocation motion, this was noted by Bhatia et al. (2014) and observed experimentally in Moriya et al. (1979) . Case II , however, finds no direct evidence in the literature except for a very general statement in Robertson (2001) , besides, case II implies hydrogen does not alter the dislocation structure. This is inconsistent with experimental evidence on hydrogen enhanced plasticity. Therefore only case I , i.e. hydrogen increases screw mobility, is used in the following sections.

Hydrogen effects on slip planarity
3.3.1. Hydrogen elastic shielding Ferreira et al. (1999) attempted to explain hydrogen enhanced slip planarity in terms of hydrogen elastic shielding, as reviewed in Section 1.2 . From the energetic point of view, it was shown with linear elasticity theory that hydrogen elastic forces reduced the energy of edge dislocations relative to screw thereby reducing the frequency of cross-slip. It should be noted that this mechanism only occurred at a high concentration, c 0 > 1 × 10 4 appm. Following this idea, we first investigate the role of hydrogen elastic shielding in enhanced slip planarity.
The applied stress was chosen to generate an out-of-plane PK force so that screw segments cross slip. Similar to the case in Section 3.2.1 , an extremely high bulk concentration of c 0 = 5 × 10 5 appm was used. The dislocation structures with and without hydrogen are superimposed in Fig. 10 , and it is clear that out of plane motion is suppressed in the presence of hydrogen.
During each iteration, all segments are checked, and a segment is considered out-of-plane if both its nodes are, | X i · n | > 0. The total dislocation line length and the proportion which is out of plane are both reduced by the hydrogen elastic force. The onset of cross slip is also delayed as shown in Fig. 11 The nodal force on screw like segments is shown in Fig. 12 . For nodes on the slip plane, F k and ˜ F k + ˆ F k act in the same direction, driving the screw like segments apart but in opposition on the edge segments. | ˜ F k + ˆ F k | | F k | but still sufficient to decrease the anisotropy in the loop expansion, decreasing the screw length and hence the amount of cross-slip. This is consistent with the theory in Ferreira et al. (1999) that hydrogen tends to reduce the pure screw population. F k has a small out-of-plane component but is predominantly in plane. Furthermore F k is entirely in plane prior to the initial out-of-plane motion which delays the first cross slip event. The hydrogen elastic force is only noticeable at extremely high bulk concentrations, and so can be neglected for bcc Fe.

Hydrogen increased mobility
The cross slip simulations were repeated with a hydrogen enhanced screw mobility law (case I) and a bulk concentration of c 0 = 0 . 2 appm. Time histories of the proportion of out-of-plane dislocation length are plotted in Fig. 13 (a).
The out-of-plane dislocation motion is delayed and suppressed in the presence of hydrogen. It is therefore hypothesised that hydrogen enhanced slip planarity is due to the hydrogen induced shape change of dislocation loops. The evolution of the proportion of screw like segments is shown in Fig. 13 (b), which decreases in the presence of hydrogen. Indicating that hydrogen enhances slip planarity by reducing the amount of cross slip as the proportion of pure screw segments is reduced.
To provide further support for this hypothesis, three simplified simulations with, B s ≡ 10 B eg , B s ≡ 5 B eg , and B s ≡ 2.5 B eg and no hydrogen, c 0 = 0 , were performed. These ideal cases, shown in Fig. 14 qualitatively resemble the hydrogen enhanced screw mobility, Fig. 13 , and are free of the complexity associated with the hydrogen distribution.
Hydrogen can therefore enhance slip planarity through hydrogen enhanced screw mobility at a realistic bulk concentration of c 0 = 0 . 2 appm for bcc Fe, in contrast to the extremely high value of c 0 = 5 × 10 5 appm required for the hydrogen elastic force. Finally, it should be reemphasized that out-of-plane dislocation motions in this work are due to cross slip of pure screw segments, and climb is not considered. The influence of hydrogen on climb will add further complexity.

Hydrogen effects in a microcantilever
In situ hydrogen charged microcantilever bend tests on single crystal Fe-3% Si alloy and single crystal FeAl alloy were recently performed by Hajilou et al. (2017) and Deng and Barnoush (2018) , respectively. Here we simulate a microcantilever with dimension L x = 12 μm and L y = L z = 3 μm to investigate hydrogen effects in a finite volume. The crystal orientation is aligned with the 100 axis along the model axis as shown in Fig. 15 ; the same orientation used in the experiments by Deng and Barnoush (2018) . However we do not attempt to reproduce the exact experimental data at the present stage, since the microcantilever geometry simulated here is simplified, with a square cross section and without a base or notch.
For simplicity only the (101)[11 1 ] slip system is included in the simulation. The sources were square prismatic loops with eight nodes, 6 of which were fixed, as shown in Fig. 15 . A total of 10 sources were used. The boundary conditions were: with an applied displacement of U = −0 . 169 μm and the remaining surfaces were traction free. Fully integrated linear brick finite elements, 0.24 μ m in size, were used to solve for the corrective elastic fields at every time increment.
Hydrogen elastic shielding and hydrogen enhanced screw mobility are considered with a realistic initial bulk hydrogen concentration for bcc Fe, c 0 = 0 . 1 appm. The sources are subjected to an approximately uniaxial stress, σ xx , which is tensile above and compressive below the neutral plane located at z = 1 . 5 μm . The initial stress state is therefore similar to the slip planarity investigation in an infinite volume discussed in Section 3.3 . The load displacement curve is shown in Fig. 16 which shows a clear softening effect through the hydrogen enhanced screw mobility law and no noticeable effect through hydrogen elastic shielding; consistent with the infinite volume simulations. The hydrogen softening can be rationalised by analyses of the dislocation structure, which is shown in Fig. 17 .  In both cases, dislocation pile-ups form at the neutral plane. However, the number of dislocations in a pile up is higher when hydrogen is included, and the proportion of screw like segments is reduced, as shown in Fig. 18 . It is therefore concluded that hydrogen enhances dislocation activity in bcc Fe and consequently leads to global softening behavior, which is consistent with the HELP mechanism and agrees with the observation of hydrogen reducing the yield stress in microcantilever experiments ( Deng and Barnoush, 2018;Hajilou et al., 2017 ). Our simulations indicate the effect is due to hydrogen increased screw mobility, and the contribution of hydrogen elastic shielding is negligible. Further, Fig. 18 (b) suggests that hydrogen enhances slip planarity due to a reduction in the proportion of screw segments available to cross slip. The simulated dislocation density increases in the presence of hydrogen, this is in contrast to the expectation in some literature ( Nibur et al., 2006 ) that enhanced planarity would lead to hardening. This could be because only one slip system was simulated and so junction formation was not considered. Or it could be due to the existence of a neutral plane which arrests a  large number of dislocations. Multiple slip systems and an additional geometry, e.g. micropillar, with a uniform stress state will be simulated in future work to provide further insight. The component of the dislocation stress, ˜ σ xx , and the total stress, ˜ σ xx + ˆ σ xx , along the beam axis, with and without hydrogen enhanced screw mobility are shown in Fig. 19 . The stress hot spot is due to the pile up of edge dislocations at the neutral plane, and is more pronounced with hydrogen due to the higher dislocation density. As the dislocation field is compressive above the neutral plane and tensile below, the overall stress and reaction force is reduced more when hydrogen is present. The dislocation displacement field, ˜ u , is shown in figure Fig. 20 . When a dislocation segment leaves the material at a free surface it generates a plastic slip step. Hydrogen promotes plastic deformation, and the number of internal and exited segments increases, which increases the magnitude of the ˜ u field and the slip step height respectively.

Conclusions
In this paper, the formulation proposed by Gu and El-Awady (2018) to incorporate the hydrogen elastic force was implemented in DDD to simulate a Frank-Read source and within DDP to simulate a microcantilever. A hydrogen enhanced mobility law for screw dislocations was also implemented; based on hydrogen-kink pair interactions simulated by Katzarov et al. (2017) using kMC.
Hydrogen elastic shielding was found to have a tendency to promote dislocation generation, decrease dislocation spacing and enhance slip planarity. The elastic hydrogen force accelerates the bowing back and annihilation of screw dipoles in FR sources reduces the proportion of screw segments, which limits cross slip. These results are consistent with the HELP mechanism ( Birnbaum and Sofronis, 1994;Ferreira et al., 1999 ). However our simulations demonstrated that this only occurs with a nonphysically high bulk hydrogen concentration, e.g. c 0 = 1 × 10 5 appm in Gu and El-Awady (2018) and c 0 = 5 × 10 5 appm in this work, which is unrealistic for bcc metals. In contrast, a significant hydrogen influence on dislocation activity was observed by implementing a hydrogen concentration dependent mobility law for the screw component of dislocation segments. At realistic bulk hydrogen concentrations, c 0 = 0 . 1 appm and c 0 = 0 . 2 appm, hydrogen was found to promote dislocation generation and enhance slip planarity. Hydrogen reduces the anisotropy in the mobility resulting in more circular loops ( Bhatia et al., 2014;Moriya et al., 1979 ), and this also tends to reduce the screw population, making the bowing back of FR sources easier and reducing cross slip. Both hydrogen elastic shielding and hydrogen increased screw mobility tend to enhance dislocation activity, and the mechanisms of the enhancement are actually quite similar. Therefore, hydrogen increased mobility should be viewed as a complementary explanation for the HELP mechanism, and it contributes a larger part than hydrogen elastic shielding, since it essentially arises from hydrogen-kink pair interactions close to the dislocation core, where very high stress and hydrogen concentrations exist but cannot be captured by linear elasticity theory.
Microcantilever simulations with a hydrogen dependent mobility law showed more pronounced dislocation pile-ups form at the neutral plane, and a reduction in the flow stress. Finally, it should be noted again that dislocation climb was excluded in all the simulations, therefore, hydrogen enhanced slip planarity is a consequence only of hydrogen affected cross slipping in the present work.