FEM formulation for mass diffusion through UMATHT subroutine

A26585 - Unrestricted

Report

FEM formula on for mass diffusion through UMATHT subrou ne

Theory and implementa on

Author(s)

Antonio Alvaro, Philippe Mainçon, Vidar Osen

SINTEF Materials and Chemistry

Document History

VERSION DATE VERSION DESCRIPTION 1.0 16/12/2014 First dra send to QA 1.1 20/12/2014 Final version

Contents

1 Introduc on 4

1.1 Transport of hydrogen in the bulk of materials . . . 4

1.1.1 Diffusible hydrogen . . . 4

1.1.2 Trapping . . . 6

1.1.3 Local Equilibrium . . . 7

2 FEM formula on for mass diffusion 8 2.1 Assump ons . . . 8

2.2 Total form . . . 8

2.3 Weak form . . . 8

2.4 Incremental form . . . 9

2.5 Discrete form . . . 9

3 Implementa on in UMATHT and UMAT 11 3.1 Comparison with Krom et al. . . 12

3.2 Deriva on of stress gradients and equivalent plas c strain . . . 12

4 Units and orders of magnitude 13 4.1 Concentra on conversion factors . . . 13

5 Usage of the user subrou ne UMATHT 15 5.1 Geometry and element connec vity . . . 15

5.2 Input constants . . . 15

1 Introduc on

This report summarizes the work done during the first year (2014) of ROP project toward building a modelling framework for describing mass diffusion in the bulk of materials. This framework has, within this project, the particular scope of simulating the transport stage of atomic hydrogen in iron or, more specifically, toward the site where degradation, and therefore Hydrogen Embrittlement (HE) occurs.

This introductory section is therefore devoted to this topic, with focus on the definitions of the the different material/mechanical contributions to the determination of hydrogen distribution in bulk iron.

Nevertheless, despite the clear application, the model framework is of general validity, i.e. it can be applied and used to describe mass diffusion for any solvent/solute system, given the correct parameter/inputs.

1.1 Transport of hydrogen in the bulk of materials

Once hydrogen is absorbed into the metal, its distribution is traditionally divided in two main contributions.

The first contribution is the hydrogen at the Normal Interstitial Lattice Sites (NILS), and indicated withCL; the second contribution is given by hydrogen atoms trapped at the material’s “imperfection” such as dislocations, grain boundaries, vacancies, inclusions or precipitates. Such contribution is indicated with C_T. These two contributions are locally in equilibrium [12], and the total final hydrogen distribution is obtained by their sum.

1.1.1 Diffusible hydrogen

In the particular case of steel, iron lattice basic structure exists in two forms: body-centered cubic (bcc) ferrite (Fe-α) and face-centered cubic (fcc) austenite (Fe-γ). The property of the steels, in terms of solubility and diffusivity, are therefore inferred from the geometrical disposition of the atoms within the lattice structures and the consequent size of the different interstitial sites. The main “engine” for hydrogen mobility is given by the concentration gradients which is generated, for instance, due to hydrogen adsorption from external sources.

Fick’s first law relates the diffusive flux to the concentration field, by postulating that the flux goes from regions of high concentration to regions of low concentration, with a magnitude that is proportional to the concentration gradient (spatial derivative) [5]. In one (spatial) dimension:

Fick’s first law: J =−D (∂C_L

∂x )

(1) where J is the hydrogen flux, C_Lis the diffusible hydrogen, Dis the diffusion coefficient or diffusivity, which is proportional to the velocity of the diffusing particle, andxis the position. The negative sign arises because diffusion occurs in the direction opposite to the direction of increasing concentration. It has to be noted that Fick’s first law, as written in Eq. 1, is consistent only for isotropic medium whose structure and diffusion properties are the same relatively to all directions. Fick’s first law applies to steady state systems, where concentration keeps constant. In many cases of diffusion, the concentration within the diffusion volume however changes with time. With some calculations based on the first Fick’s law and mass balance, Fick’s second law in the case of 3-dimensional diffusion is obtained:

Fick‘s second law: ∂C_L

∂t =−D∇²C_L (2)

wheretis time and∇is thedeloperator. In a 3D system with perpendicular coordinates (x, y, z), this is a Cartesian coordinate systemR³,delis defined as:

∇=i ∂

∂x+j ∂

∂y+k ∂

∂z (3)

where i, j and k are the unit-vectors in the direction of the respective coordinate (the standard basis inR³).

Eq. 1 and 2 defines just the basic laws to describe hydrogen flux within an unstressed perfect lattice structure.

In reality, the correct way to look at hydrogen diffusion through the steel is not as a straight line trajectory. As showed by Jiang and Carter [7] through first principle Density Functional Theory (DFT) calculations, H-atom diffuses through bcc Fe not via a straight trajectory line, but rather jumps from one t-site (t stays for tetrahedral) to another neighboring t-site by a curve path distorted toward octahedral sites.

Figure 1: Interstitial sites present in a BCC crystal structure [2].

As mentioned before, the small size of hydrogen allows higher solubility and greater mobility than other elements. However, hydrogen atoms are still larger than the interstitial sites: r_H = 0.53 Å, to be compared withri = 0.19Å for o-sites in austenitic steel (the largest among the interstitial sites in steels). This induces a distortion of the host lattice and the resulting stress and displacement fields interact with other defects. The

“strength” of this type of interaction is quantified by the hydrogen’s partial molar volume, i.e. the unconstrained volume dilatation of the metal containing one mole of hydrogen. An isochore introduction of hydrogen into a lattice creates a hydrostatic compression stress [3, 4]. Thereby, the mean hydrostatic stress affects the hydrogen solubility in the host metal, and hydrostatic stress gradients affect hydrogen diffusion [10]. Therefore, a more complete description for hydrogen diffusion should consider both for hydrogen concentrations and hydrostatic stress gradients. Abaqus [1] provides the following mass diffusion equation already embedded with the finite element code:

∂CL

∂t =D∇²CL+D

( VH

R(T −T^Z) )

∇CL∇p+D

( VH

R(T −T^Z) )

CL∇²p (4) In conclusion, relation 4 indicates that the energy necessary to introduce a hydrogen atom in the lattice increases together with the decrease of hydrogen concentration gradients and is decreased by dilatational hydro- static stress: areas in front of cracks are therefore zones of strong accumulation for diffusible hydrogen.

1.1.2 Trapping

Hydrogen distribution is determined not only by the concentration gradients and lattice dilatation due to hydro- static stress as described by 4, but it is recognized that also trapping is a potent mechanism for H segregation [12], [16], [13]. Traps can be seen as local perturbations of the lattice structure: typical traps are dislocations, vacancies, grain boundaries, phase boundaries, inclusions and precipitates. Trapping reduces the amount of mobile hydrogen inferring a decrease of the apparent diffusivity and increasing the local solubility favoring segregation. The ability of traps to “hold” hydrogen atoms is associated with hydrogen binding energy and activation energy for hydrogen release (∆E_T) and, based on the values of these terms, trapping sites are usually categorized as reversible and irreversible traps. Reversible traps concern binding energies usually lower than 60/70kJ/mol[17] and hydrogen can typically be released by simple tempering. On the contrary, for irrevers- ible traps normal tempering result uneffective since the energy barrier to be overcome in order for hydrogen to regain mobility is higher. A typical example of irreversible trap sites is the interface between non-metallic inclusion or precipitates jump.

Of special interest is hydrogen trapping related to dislocation which has been studied extensively in the last decades [16],[6],[3], [4], [14]. In a nutshell, dislocations, which in metals are associated with the plastic de- formation development, can be imagined as moving traps, having therefore a great impact on the total hydrogen distribution in the material. Kumnick and Johnson [9] firstly calculated the variation of the amount and the densities of trapping sites in steel at different cold-working levels. Their results are graphically presented (see Fig. 2) by Sofronis and McMeeking [15], who also proposed the following mathematical fitting:

Log(N_T) = 23.26−2.33exp(−5.5ϵ_p) (5) whereN_T is the trap density andϵ_pis the equivalent plastic strain: trap densities are independent on temperat- ure, increase sharply with deformation at low deformation levels and more gradually with further deformation, reaching saturation at plastic strain greater than 80%.

Figure 2: Relation between number of traps and plastic strain fitted from Kumnick and Johnson experimental results [15].

1.1.3 Local Equilibrium

As summarized in the previous sections, two types of hydrogen concentrations can be distinguished: hydrogen in NILS, indicated byC_Land hydrogen in traps,C_T. Oriani [12] proposed the principle that for rapid trap filling kinetics, these two “populations” are locally in equilibrium. Such principle is described by the relation

C_T = K^αN_βN^T

LCL

1 +_βN^K

LCL

(6) whereN_T andN_Lare the available sites for hydrogen in traps and in lattice sites, respectively. N_Lis calculated through:

NL= NA

VM

(7) whereN_Ais the Avogadro number andV_M is the molar volume of iron. αandβare parameters which indicates the site occupancies for traps and lattice sites, respectively. The parameterαis taken as 1 whileβis assigned to 6 under assumption of tetrahedral site occupancy. KT represents the equilibrium constant between lattice and trap sites:

K_T =exp

(−∆E_T RT

)

(8)

∆E_T is the trap binding energy, always negative, and assigned according to the type of traps under consideration [9]. The inherent meaning of Oriani’s equilibrium principle is that, at local level, the two populations of hydrogen affect one another: variations of lattice hydrogen modifies the reversibly trapped hydrogen and viceversa. It is important also to point out that such theory considers exclusively reversible traps.

2 FEM formula on for mass diffusion

The FEM formulation of hydrogen diffusion here implemented has been derived from by [8] which has been obtained through Eq. 1. to Eq. 8 . The formulation is revisited here.

2.1 Assump ons The system is modeled by

C˙L+ ˙CT +∇ ·J = 0 (9) with

J = −DL∇CL−DpCL∇p (10) C_T

NT −CT

= C_L

NL−CL

KT (11)

K_T = e^−∆RTET (12)

D_p = D_LV_H

RT (13)

2.2 Total form

By multiplying denominator and divisor byCL, (11) can be rewritten as CT = KTNTCL

N_L+C_L(K_T −1) (14)

AssumingNLandKT not to vary with time

∂C_L

∂t +∂C_T

∂t = ∂C_L

∂t +∂C_T

∂C_L

∂C_L

∂t + ∂C_T

∂N_T

∂N_T

∂t (15)

= [

1 + NTNLKT

(N_L+C_L(K_T −1))² ]∂CL

∂t +

[ KTCL

NL+CL(KT −1) ]∂NT

∂t (16) so finaly

0 =f (

CL,C˙L

)

+∇ ·J(CL) (17)

with

f (

CL,C˙L

)

= (

1 + NTNLKT

(N_L+C_L(K_T −1))² )

C˙L+ KTCLN˙T

N_L+C_L(K_T −1) (18) This is a single scalar non-linear differential equation inCL,NT andp.CT has been eliminated using (14).

In the following we will assume NT andp to be known (weak connection to another solution). We further assume all other coefficients to be known and constant in time and space.

2.3 Weak form

The weak form of the differential equation is obtained by requiring that for any “virtual variation”δC_Lof the unknown fieldCL, equation (17) multiplied byδCLand integrated overV, must be verified

0 =∫

VδC_Lf dV + ∫

VδC_L∇ ·J dV (19)

=∫

VδC_Lf dV + ∫

V∇ ·(

δC_LJ)

dV −∫

V∇δC_L·J dV (20)

=∫

VδC_Lf dV + ∫

Sn·(

δC_LJ)

dS −∫

V∇δC_L·J dV (21)

0 =∫

VδCLf dV + ∫

SδCLϕ dS −∫

V∇δCL·J dV (22)

with

ϕ = n·J (23)

wherenis the outward pointing normal toSandϕis the hydrogen flow out throughS.

The integration by parts has allowed to introduce boundary conditions: at any point ofS, eitherϕmust be known (non-essential boundary condition) orδC_Lmust be zero (essential boudary condition).

2.4 Incremental form

The incremental form ofE = 0inC_LandC˙_Lis 0 =E+ ∂E

∂CL

∆C_L+ ∂E

∂C˙L

∆ ˙C_L (24)

and hence the incremental form of (22) is (assumingϕto be a known boundary condition) 0 =

∫

VδC_Lf dV +

∫

SδC_Lϕ dS −

∫

V∇δC_L·J dV (25) +

∫

VδCL

∂f

∂CL

∆CLdV −

∫

V∇δCL· ∂J

∂CL

∆CLdV (26)

+

∫

VδC_L ∂f

∂C˙_L∆ ˙C_LdV (27)

or

0 =

∫

VδC_L [(

1 + NTNLKT

(N_L+C_L(K_T −1))² )

C˙_L+ KTCLN˙T

N_L+C_L(K_T −1) ]

dV +

∫

SδCLϕ dS

−

∫

V∇δCL·[

−DL∇CL−DpCL∇p] dV +

∫

VδC_L

[−2 (K_T −1)N_TN_LK_TC˙_L

(NL+CL(KT −1))³ + K_TN_LN˙_T (NL+CL(KT −1))²

]

∆C_LdV

−

∫

V∇δC_L·[

−D_L∇ −Dp∇p]

∆C_LdV +

∫

VδC_L [

1 + N_TN_LK_T (NL+CL(KT −1))²

]

∆ ˙C_LdV (28)

Since (28) must hold for anyδCL, it must hold for anyδCLand we obtain:

0 =

∫

VSC

[(

1 + NTNLKT

(N_L+C_L(K_T −1))² )

C˙L+ KTCLN˙T

N_L+C_L(K_T −1) ]

dV +

∫

SS_C ϕ dS

−

∫

VG^T_C ·[

−DL∇CL−DpCL∇p] dV

+

∫

VS_C

[−2 (K_T −1)N_TN_LK_TC˙_L

(N_L+C_L(K_T −1))³ + K_TN_LN˙_T (N_L+C_L(K_T −1))²

]

S_C dV ·∆C_L

−

∫

VG^T_C ·[−DL]GC dV ·∆CL

−

∫

VG^T_C ·[

−Dp∇p]

SC dV ·∆CL

+

∫

VS_C [

1 + N_TN_LK_T (NL+CL(KT −1))²

]

S_C dV ·∆C˙_L (33)

This can be rewritten

K·∆C_L+M·∆C˙_L=F (34) with

K =

∫

VS_C

[−2 (K_T −1)N_TN_LK_TC˙_L

(NL+CL(KT −1))³ + K_TN_LN˙_T (NL+CL(KT −1))²

] S_CdV +

∫

VG^T_C·[D_L]G_CdV +

∫

VG^T_C·[ D_p∇p]

S_CdV (35)

M =

∫

VSC

[

1 + NTNLKT

(N_L+C_L(K_T −1))² ]

SCdV (36)

F = −

∫

VSC

[(

1 + NTNLKT

(N_L+C_L(K_T −1))² )

C˙L+ KTCLN˙T

N_L+C_L(K_T −1) ]

dV

−

∫

SSCϕ dS

−

∫

VG^T_C·[

DL∇CL+DpCL∇p]

dV (37)

Note that (37) does not contain second order derivatives if eitherCLorp, thanks to the sequence weak form- partial integration-incremental form. The alternative sequence incremental form-weak form-partial integration yields the same formulation, except for second derivatives inF.

For the special case in which there are no traps (hydrogen is only dissolved in the lattice), it is found by settingN_T = ˙N_T = 0, or equivalently,K_T = 0in the above expressions.

3 Implementa on in UMATHT and UMAT

As mentioned in the previous sections, Abaqus’ option for mass transfer analysis cannot be used because it does not allow to specify the non-linear diffusivity and solubility due to dislocation traps. However, an alternative way can be found by taking advantage of the similarity between Fourier’s and Fick’s law (i.e. solute concen- tration can be treated as temperature): the above formulation is implemented in Abaqus by using UMATHT, a user definable subroutine for heat transfer.

The following table shows the relation between the inputs and outputs to MATHT and the present formulation.

Function Variable Heat Transfer Mass Diffusion

UMATHT U U(θ, t) CL+CT

DUDT ^∂U_∂θ ^∂(C_∂C^L+CT)

L = 1 + ^NTNLKT

(NL+CL(KT−1))²

DUDG ^∂U

∂∇θ

∂(CL+CT)

∂∇CL = 0

TEMP θ CL

DTEMP ∆θ ∆C_L

DTEMDX ∇θ ∇CL

FLUX f(

θ,∇θ)

J(

CL,∇CL

)=−DL∇CL−DpCL∇p

DFDT ^∂f_{∂θ ∂C}^∂J

L =−D_p∇p

DFDG ^∂f

∂∇θ

∂J

∂∇CL =−DL

UMAT RPL ∆U 0

DDSDDT ^∂∆σ_∂θ ^∂∆σ_∂C

L

DRPLDE ^∂∆U_{∂ϵ N} ^KTCL

L+CL(KT−1)

∂NT

∂ϵ

DRPLDT ^∂∆U_∂θ ^{−2(KT−1)NTNLKT}

(NL+CL(KT−1))³

Table 1: Correspondence with UMATHT and UMAT

Thinking of a future extension to coupled analysis (diffusion and continuum nechanics), the absence of a matrix containing the term

K_TC_L N_L+C_L(K_T −1)

implies a weak connection between the analyses, which may be a source of convergence problem.

Following [11], UMATHT allows an exact and (thought not Newton-optimal) solution of the uncoupled hydrogen diffusion problem. Moriconi [11] further uses this UMATHT implementation in a coupled diffusion- mechanical analysis. Apparently, Abaqus limits what mechanical material models can be used in a coupled analysis, so it was necessary to implement UMAT to provide the relevant material plasticicty model. So the coupled analysis was implemented without implementing UEL (User defined ELement subroutine). Note that this does not allow to model a strongly coupled problem (due to the above-mentioned missing term), for that it seems necessary to implement a complete element (UEL).

3.1 Comparison with Krom et al.

For largeK_T, one can introduce the approximationK_T −1≈K_T. Then it is possible to show that

∂C_T

∂C_L = C_T

(

1−^C_N^T_T)

C_L (38)

∂C_T

∂N_T = C_T

N_T (39)

so that the differential equation they solve is the correct one (even though they offer no proof for it in that publication).

3.2 Deriva on of stress gradients and equivalent plas c strain

In our solution, the gradient of the pressure in one element is computed as follow:

The pressure at the Gauss points (center) of all neighboring elements is considered and a linear interpolation of the pressure as a function of position in the deformed mesh is introduced. The slope of this linear relation is use as a pressure gradient.This is just one of many possible approaches to finding the pressure gradient. All of them have in common that they are imperfect (computing slopes on numerical data is only apparently easy).

Nevertheless, it is not known which approach Abaqus is used to compute the gradient.

The pressure is given by Abaqus as the third stress invariant for the integration point. The stress invariants can be read by callingGETVRM()and specify'SINV'as the value for the parameterVAR.

The plastic strain is given by Abaqus as the seventh plastic strain value for the integration point. The strain values can be read by callingGETVRM()and specify'PE'as the value for the parameterVAR.

4 Units and orders of magnitude

In order to build a model of the most general validity, particular attention needs to be put into the use of the unit system as well as in the list of the constants. The modelling framework is therefore built buy using con- stants which will eventually allow the study of mass diffusion of systems which are different than hydrogen/iron.

The unit system used is MPa meaning that the basic units are : mass in [tonne], length in [mm], time in [s], temperature in [K], pressure in [MPa], energy in [mJ] and so forth. Table 2 summarises the units for the quantities described in the previous section along with their expected order of magnitude for hydrogen/iron system.

Symbol Order of magnitude Unit

C_L mol·mm⁻³

CT mol·mm⁻³

NL ≈10⁻⁰⁴ mol·mm⁻³

N_T 0.16603·10^{NT1−NT2·eNT3·ϵp}* mol·mm⁻³

N_T1 −8.74* mol·mm⁻³

N_T2 2.33*

NT3 −5.5*

J mol·m⁻²·s⁻¹

∆E_T ≈10⁶ mol⁻¹·tonne·mm²·s⁻² R 8.314·10³ mol⁻¹·K⁻¹·tonne·

mm²·s⁻²

T 293. K

∇C_L mol·mm⁻⁴

D_L ≈10⁻⁶ mm²·s⁻¹

Dp ≈10⁻⁹ tonne⁻¹·mm³·s

p tonne·mm⁻¹·

s⁻²(M P a)

∇p tonne·mm⁻²·s⁻²

Table 2: Units and orders of magnitude (*: constants obtained form Sofronis and McMeecking [15] interpolation of Kumnick and Johnson[9] experimental data)

4.1 Concentra on conversion factors

C[

at_H·m⁻_{F e}³]

=C[wppm]·n_A[

atoms·mol⁻¹]

·( z_H[

tonne·mol⁻¹]

·ρ_{F e}[

tonne·mm⁻³])₋1

(42) C[

mol·mm⁻³]

=C[

atH ·m⁻_{F e}³]

·10⁻⁹[

m³·mm⁻³]

·( nA

[atoms·mol⁻¹])₋₁

(43) C[appm] =C[

at_H ·m⁻_{F e}³]

·z_{F e}[

tonne·mm⁻³]

·( n_A[

atoms·mol⁻¹]

·ρ_{F e}[

tonne·mm⁻³])₋1

(44)

The constants that are necessary to perform the conversions described above (N.B. the constants defined below and the consequently calculated conversion factors are relative to hydrogen/iron system) are the following:

V_M : M olar V olume: 7.119·10⁻⁶ m³·mol⁻¹ = 7.119·10³ mm³·mol⁻¹ which indicates the volume occupied by one mole of substance at a given temperature and pressure;

VH : P artial M olar V olume : 2·10⁻⁶ m³ ·mol⁻¹ = 2·10³ mm³ ·mol⁻¹ it is a thermodynamic quantity which indicates the unconstrained volume dilatation of a metal (iron in this case) due to the absorption/introdution of one mole of a solute (hydrogen);

z_H : Hydrogen M olar M ass: 1.008g·mol⁻¹ = 1.008·10⁻⁶tonne·mol⁻¹mass of an atom of hydrogen;

zF e: Iron M olar M ass: 55.845g·mol⁻¹ = 55.845·10⁻⁶tonne·mol⁻¹mass of an atom of pure iron;

nA : Avogadro N umber : 6.023·10²³ atoms·mol⁻¹ it defines the number of atoms per mole for any substance;

ρF e: Iron Density : 7.8747tonne·m⁻³= 7.8747·10⁻⁹ tonne·mm⁻³.

Finally, Table 3 reports conversion factors as calculated by formulas in the group of Eq. (44) (unit of the value to convert “enters” from the column and the converted results “exit” from the row):

CONC wppm appm atH ·m⁻_{F e}³ mol·mm⁻³ wppm 1 55.402 4.70·10³⁰ 7.81·10⁻³

appm 0.01805 1 8.493·10²⁸ 1.41·10⁻¹⁰ at_H·m⁻_{F e}³ 2.12·10⁻³¹ 1.17·10⁻²⁹ 1 1.66·10⁻³³ mol·mm⁻³ 128 7.092·10⁹ 6.023·10³² 1

Table 3: Concetration conversion factors

5 Usage of the user subrou ne UMATHT

5.1 Geometry and element connec vity

UMATHT reads the geometry (node positions) and element connectivity from a filetopology.dat. This file needs to be generated before the calculation can start. A separate program,ReadTopologyhas been created to produce thetopology.datfile. This program requires that this information is available in a*.filfile.

Instruct Abaqus to produce a.filfile by inserting adding code like this in the input file:

**

** OUTPUT REQUESTS

**

*Restart, write, frequency=0

**

*Output, History

*Element Output

**

*El File, Position=Integration Point COORD

**

** FIELD OUTPUT: F-Output-1

**

*Output, field, variable=PRESELECT

*Output, history, frequency=0

*End Step

To produce the*.filfile, run Abaqus with the commanddatacheck:

abaqus datacheck interactive user="UMATHT.OBJ" job="myjobname" input="myjobname.inp"

The filetopology.datcan then be produced by running abaqus readTopology myjobname

5.2 Input constants

To specify a material that uses the formulation in UMATHT, specify the material like this in the input file of Abaqus:

*User Material, constants=10, type=THERMAL

It should be followed by two lines containing the a list of 10 parameters; 8 on the first line and 2 on the last line:

1. D : Mass diffusivity

8. EB: Trap binding energy

9. Material temperature in Kelvin (used to calculateK_T)

10. R: Universal gas constant (with units consistent withE_Band material temperature).

Example:

*USER MATERIAL, constants=10, type=THERMAL

** DL DP VM NT_C0 NT_C1 NT_C2 NT_C3 EB

3.99E-6, 3.4629E-9, 7.116E+3, 0.16603, -8.740, 2.330, -5.5, -28.8E+6

** T R

293.15, 8.3144621E+3

** in front of the line means that the line is commented while *indicates an Abaqus command line with appro- priate keywords.

References

[1] ABAQUS. Version 6.11. analysis user’s manual. ABAQUS Inc., 2012.

[2] W.D. Callister and D.G. Rethwisch. Materials science and engineering: an introduction. John Wiley and Sons New York, 1997.

[3] J.P. Chateau, D. Delafosse, and T. Magnin. Numerical simulations of hydrogen-dislocation interac- tions in fcc stainless steels.: part i: hydrogen-dislocation interactions in bulk crystals. Acta materialia, 50(6):1507–1522, 2002.

[4] J.P. Chateau, D. Delafosse, and T. Magnin. Numerical simulations of hydrogen–dislocation interactions in fcc stainless steels.: part ii: hydrogen effects on crack tip plasticity at a stress corrosion crack. Acta materialia, 50(6):1523–1538, 2002.

[5] J. Crank. The mathematics of diffusion. Oxford University Press, USA, 1979.

[6] D. Delafosse and T. Magnin. Hydrogen induced plasticity in stress corrosion cracking of engineering systems. Engineering Fracture Mechanics, 68(6):693–729, 2001.

[7] D.E. Jiang and E.A. Carter. Diffusion of interstitial hydrogen into and through bcc fe from first principles.

Physical Review B, 70(6):64102, 2004.

[8] A.H.M. Krom, R.W.J. Koers, and A. Bakker. Hydrogen transport near a blunting crack tip.Journal of the Mechanics and Physics of Solids, 47(4):971–992, 1999.

[9] A.J. Kumnick and H.H. Johnson. Deep trapping states for hydrogen in deformed iron. Acta Metallurgica, 28(1):33–39, 1980.

[10] J.C.M. Li, R.A. Oriani, and L.S. Darken. The thermodynamics of stressed solids.Zeitschrift fur Physikalis- che Chemie, 49(3-5):271–290, 1966.

[11] C. Moriconi. Modelisation de la propagation de fissure de fatigue assistee par l’hydrogene gazeux dans le materiaux metalliques. These de doctorat, Ecole Nationale Superieure de Mecanique et d’Aerotechnique, 2012.

[12] R.A. Oriani. The diffusion and trapping of hydrogen in steel. Acta Metallurgica, 18(1):147–157, 1970.

[13] G.M. Pressouyre. Trap theory of hydrogen embrittlement. Acta Metallurgica, 28(7):895–911, 1980.

[14] I.M. Robertson. The effect of hydrogen on dislocation dynamics. Engineering Fracture Mechanics, 68(6):671–692, 2001.

[15] P. Sofronis and R.M. McMeeking. Numerical analysis of hydrogen transport near a blunting crack tip.

Journal of the Mechanics and Physics of Solids, 37(3):317–350, 1989.

Technology for a be er society

www.sintef.no