A peridynamic mechano-chemical damage model for stress-assisted corrosion

Electrochimica Acta 323 (2019) 134795

Contents lists available at ScienceDirect

Electrochimica Acta journal homepage: www.elsevier.com/locate/electacta

A peridynamic mechano-chemical damage model for stress-assisted corrosion Siavash Jafarzadeh a, Ziguang Chen b, **, Shumin Li a, Florin Bobaru a, * a b

Department of Mechanical and Materials Engineering, University of Nebraska-Lincoln, Lincoln, NE, 68588-0526, USA Department of Mechanics, Huazhong University of Science and Technology, Wuhan, 430074, China

a r t i c l e i n f o

a b s t r a c t

Article history: Received 7 June 2019 Received in revised form 27 August 2019 Accepted 30 August 2019 Available online 31 August 2019

We introduce a mechano-chemical coupled peridynamic model for stress-assisted corrosion. In the model, the activation-controlled anodic dissolution is a function of nonlocal dilatation of material near the corrosion front. We use a new experimental design to measure anodic dissolution rate for copper under activation-controlled and applied tensile stress. The experimental results are used to calibrate and then validate the new model. A case for stress-assisted corrosion in a complex geometry with nonuniform stress distribution is included to demonstrate the model's capabilities. © 2019 Elsevier Ltd. All rights reserved.

Keywords: Corrosion Damage Stress-dependent corrosion Peridynamics Stress concentration

1. Introduction The influence of stress on corrosion rates and damage evolution in metallic parts of machines and structures that carry mechanical load is being actively investigated [1e11]. The enormous cost of failure of corrosion-affected systems, as well as providing maintenance and implementing corrosion prevention methods, require new computational tools that can predict the time evolution of corrosion-induced damage in the materials. The effect of stress on corrosion rate is not completely understood. It can differ among corrosion systems and different loading conditions. Elastic and plastic strains, compressive or tensile, can affect any of the chemical reactions (anodic, cathodic, passivation, etc.) present during the corrosion process [2,3,5,6,9e13]. Nearly all experimental investigations suggest that tensile elastic stress is likely to increase the anodic dissolution rate [1,2,4,6e9,12]. This behavior is usually explained via a mechano-chemical theory by Gutman [14]. While analytical models based on this theory or similar approaches are useful for overall estimation of stressassisted corrosion in the elastic regime, computational corrosion models that consider stress dependency in their governing

* Corresponding author. ** Corresponding author. E-mail addresses: [email protected] (Z. Chen), [email protected] (F. Bobaru). https://doi.org/10.1016/j.electacta.2019.134795 0013-4686/© 2019 Elsevier Ltd. All rights reserved.

equations have the additional advantage of detailed predictions for corrosion damage evolution in various types of corrosion such as pitting, crevice, galvanic, etc., in arbitrary geometries [15e18]. The focus of this study is on modeling mechano-chemical effect of stress on the anodic dissolution rate. This phenomenon is a part of stress corrosion cracking, which involves crack initiation and growth under mechanical loading in a corrosive environment [19,20]. Only a few computational models have been developed to date to address stress-assisted corrosion with non-uniform stress distribution. Based on mechano-chemistry analysis, Sarkar and Aquino [21] defined the current density being dependent on hydrostatic pressure and surface curvature. In their model, stress distribution and curvature of the front are inputs for computing current density and the motion of the dissolution front [21]. This model is used for a simple geometry with known analytical stress field and curvature. In real cases with complex morphologies, calculating local curvature at dissolution front can be a challenge. Wang and Han [16] developed a coupled model that computes, at each time step, stress distribution via the finite element method (FEM), then calculates current density distribution based on Gutman's theory [14], and simulates the motion of the corrosion front via cellular automata (CA). Wang and Han [16] used their model to simulate growth and coalescence of two pits under uniaxial tension via a 2D simulation. A similar coupled FEM-CA model proposed by Fatoba et al. [17] is

2

S. Jafarzadeh et al. / Electrochimica Acta 323 (2019) 134795

used for 2D simulation of pitting corrosion under tensile stress. In that study, the CA parameters are optimized for a particular experiment. While CA models are shown to be useful for qualitative conclusions and can reproduce some of the realistic features of localized corrosion, the model parameters and even the corrosion time need to be adjusted by optimization for each actual corrosion case [15]. The reason is that the CA parameters are not electrochemo-physical quantities, but they are usually probabilities that correlate with heuristic rules [22]. This limit the prediction capability of CA models [23e25]. Another drawback of the FEM-CA coupled models is the moving boundary process, which, in an FEM analysis requires remeshing, hence, significantly increases the computational cost of simulations. Phase-field theory is also used to model stress-assisted corrosion [18]. In the phase-field model by Lin et al. [18] the corrosion rate depends on the strain energy on the surface. The model is used for 2D and 3D simulations of mechanochemical corrosion in metallic composites. However, the results and the model's theoretical basis are not validated against experimental observations. In this paper, we introduce a computational mechano-chemical coupled model for stress-assisted corrosion using the peridynamic (PD) theory. The peridynamic corrosion model has been recently shown to be a reliable and flexible tool for simulation of corrosion damage [26e30]. To date, it has been used to simulate the damage evolution in pitting and intergranular corrosion in various alloys [27e30]. The model considers corrosion as a type of damage, which facilitates its coupling to mechanical analysis. Another advantage of the PD corrosion model, is that it can capture the experimentally observed subsurface damage at the corrosion front, which has a gradual change in composition and mechanical integrity compared with the intact bulk substrate [31e35]. In addition, the corrosion evolution occurs with an autonomous phase change process that eliminates the computationally expensive process for explicit tracking of the corrosion front and domain re-meshing [26,28,36]. We also introduce an experimental procedure for calibration of the model parameter that specify the dependency of the anodic reaction rate to tensile elastic stress. We calibrate the model for copper and verify it against several stress levels. We use the model for a more complex geometry to demonstrate the capability of our stress-dependent corrosion computational model in arbitrary domains. Advantages of the new model compared with existing ones are also discussed. The novelty of this work consists in coupling the existing peridynamic model for corrosion damage with the mechanical peridynamic model. The new model is capable of simulating degradation caused by the combined action of mechano-chemical phenomena, as the corrosion rate is directly influenced by the stress/strain in the body. Also new is the validation of the mechanochemical coupled model for stress-assisted corrosion against experiments we conducted.

a segment in 1D, a disk in 2D and a sphere in 3D) centered at x, is called the horizon region of x, and is denoted by Hx. Its radius is called the horizon size (or simply, the horizon) and often denoted by d. Fig. 1 shows a schematic PD body and the horizon for a generic point. Objects that transfer the pair-wise information between points are called bonds. Representing physical phenomena as such, leads to governing equations with integrals instead of spatial derivatives, which are common in the classical local theories. As a result, in PD models the domain or the fields are not required to be smooth, not even continuous. Hence, PD can easily model phenomena in bodies with evolving discontinuities, like crack propagation in fracture problems and damage propagation in corrosion problems [37,38]. Dynamic brittle and ductile fracture [39e42], thermally-driven cracks [43], fatigue cracking [44], failure of composites [45,46], damage of porous and granular materials [47e49], pitting corrosion damage [28,30], intergranular corrosion damage [29], and stress-corrosion cracking [50,51] are among the various applications of PD models.

2.2. The peridynamic model for elastic behavior and brittle fracture According to bond-based PD [52,53], equilibrium equation is:

ð

f ðb y  y; b x  x; xÞdVbx þ b ¼ 0

(1)

Hx

where f is the micro-force density that the neighbor point b x exerts on x; y and b y are respectively the position vectors of points x and b x in deformed configuration, and b is the body force on point x. Fig. 2 plots schematics of a peridynamic body at reference and deformed configurations. Eq. (1) states that the summation of forces on a generic point x, is zero in the equilibrium state. For derivation details see Refs. [37,52]. The micro-force density for a microelastic solid is expressed as [53]:

f ðb y  y; b x  x; xÞ ¼ mðx; b x ÞcðE; dÞsð b y  y; b x  xÞnð b y  yÞ

(2)

where s is the bond strain (or relative deformation):

2. Model description We briefly describe the peridynamic (PD) models for elasticity and corrosion damage. We then introduce the new coupled mechano-chemical model for stress-dependent corrosion. 2.1. The peridynamic theory Peridynamic is an extension of continuum mechanics that describes the behavior of continuous bodies, bodies with evolving discontinuities, and discrete bodies, within a unified nonlocal framework [37,38]. In this theory, physical phenomena at a generic material point x are described via interactions of that point with other points in its neighborhood. This neighborhood (which can be

Fig. 1. Schematic of a peridynamic body, with a generic point x, and its horizon (domain of nonlocal interactions).

S. Jafarzadeh et al. / Electrochimica Acta 323 (2019) 134795

3

b are respectively, the displacement vectors for points x and b Fig. 2. Schematic of a peridynamic body in the reference and deformed configurations. u and u x.

s¼

y  yk  k b x  xk kb x  xk kb

(3)

vC ¼ vt

n is the unit vector of the deformed bond:

n¼

b yy b  yk ky

(4)

c is called the micromodulus and is calculated using the Young modulus E:

c¼

9E

(5)

pd3

Eq. (5) gives the relationship for 2D plane stress, under the assumption of a micromodulus that is constant over the horizon region. For other options, these relationships are shown in several papers [38,54,55]. m in Eq. (2) is a history-dependent binary scalar-valued quantity [53]:

8 <

(6)

A bond is “broken” when its strain exceeds a limit which is calculated from the material's critical fracture energy [38]. Accordingly, damage is defined for each point as a scalar field ranging in [0, 1]:

dðxÞ ¼ 1 

Hx

Hx

b C   C k d; b d dVbx ; b k x  xk2

(8)

where C ¼ Cðx; tÞ is the concentration of metal at point x and time t, b ¼ Cð b and C x ; tÞ. Likewise, d ¼ dðx; tÞ is the damage at point x and time t, and b d ¼ dð b x ;tÞ. k is the damage-dependent micro-diffusivity of x b x diffusion bonds. Similar to the PD model for elasticity, in which mechanical bonds transmit pairwise forces between neighboring points and nodes store the displacement information (see Eq. (1)), diffusion bonds in nonlocal mass transfer transmit mass fluxes between points and nodes carry mass concentration information. k is defined for each bond as a function of damage, meaning that its value depends on the phases at the two ends of the bond [30]:

8 > >  b  < kL ðDÞ k d; d ¼ 0 > > : kdiss

; dðx; tÞ ¼ 1 &dð b x ; tÞ ¼ 1 ; dðx; tÞ < 1&dð b x ; tÞ < 1 ; ½dðx; tÞ or dð b x ; tÞ < 1&½dðx; tÞ or dð b x ; tÞ ¼ 1 (9)

if x b x bond is intact if x b x bond is broken

mðx; bx Þ ¼ 1 :0

ð

ð

mðx; bx ÞdVbx ð

(7) dVbx

Hx

In the discretized version, nodal damage is the ratio of the broken bonds to the total number of bonds at a node. A material described by this model behaves as linear elastic up to failure, and the failure is brittle.

2.3. The peridynamic corrosion damage model The PD corrosion damage model [26,30] assumes a damage value of one (d ¼ 1) for the liquid or the corroded phase where electrolyte is present, and damage less than one (0  d < 1) for the solid phase. The corrosion process is formulated based on a damage-dependent nonlocal diffusion equation [30]:

According to Eq. (9), if both points are in the liquid phase, the k ¼ kL is calculated based on the diffusivity of the electrolyte, D [56,57]:

kL ¼

4D

pd2

(10)

Similar to Eq. (5), Eq. (10) gives the relationship for 2D and assuming that the microdiffusivity is constant over the horizon. Relationships for 1D and 3D are shown elsewhere [26]. Note that different kernels can be used in Eqs. (1) and (8), which result in different forms for the micro-modulus and micro-diffusivity values given in Eqs. (5) and (10) [58,59]. Eq. (9) states that no diffusion takes place if both points belong to the solid phase. However, if one point is solid and the other point is liquid, the bond is an interfacial bond which carries the anodic dissolution micro-flux. For such bonds, k ¼ kdiss is called the microdissolvability [30]. Activation controlled anodic current density (denoted by i) scales linearly with kdiss . Any corrosion behavior can then be captured in this model by defining a proper formula for kdiss corresponding to the anodic current density. Fig. 3 shows schematics of the metal-electrolyte domain in PD corrosion models. While Eq. (8) models mass transfer via anodic dissolution at the interface and diffusion in the electrolyte, a concentrationdependent damage model, describes the damage propagation during the corrosion process [26,30]:

4

S. Jafarzadeh et al. / Electrochimica Acta 323 (2019) 134795

Fig. 3. Schematic of the metal-electrolyte domain in peridynamic corrosion models [30].

8 > > > > > > >
1

solid  Cðx; tÞ dðx; tÞ ¼ > C solid  Csat > > > > > 0 > :

  sm Vm is ¼ iu exp RT

; Cðx; tÞ  Csat ; Csat < Cðx; tÞ < Csolid

(11)

; Cðx; tÞ ¼ Csolid

Eq. (11) implies that when concentration in the solid drops below the saturation concentration of metal ions in the electrolyte (Csat ) by dissolution, the solid phase changes to liquid (we set d ¼ 1). The damage for intact solid (full-concentration of Csolid ) is zero. In this model, the solid region within the d distance of the metal/ electrolyte is the dissolving region (see Fig. 3) and possesses a concentration value between Csat and Csolid . Eq. (11) provides a relationship that determines the damage for this region according to the concentration. This partially damaged dissolving solid region with the gradual change in concentration, is in fact, the diffuse corrosion layer [31] observed experimentally near the corrosion front, where, over a thickness of several micrometers, gradual changes in composition and degraded mechanical properties take place [31e35]. The concentration-dependent damage in Eq. (11) translates into mechanical damage by removing mechanical bonds in the corroded region. This is done by calculating mðx; b x Þ for all mechanical bonds such that damage evaluated from Eq. (7) returns the same value as that obtained from Eq. (11). The mðx; b x Þ values are found via a stochastic procedure based on d and b d. Details are given elsewhere [26,30]. Since corrosion is approached as mechanical damage in the model, the PD corrosion model can be easily coupled with PD mechanical and fracture models in order to simulate stress-assisted corrosion (presented here) and stress corrosion cracking [50,51].

where is and iu are respectively the anodic current densities in stressed and unstressed material. sm is the mean stress or the hydrostatic pressure, Vm is the molar volume, R is the gas constant, and T is the absolute temperature. The hydrostatic pressure in PD can be expressed in terms of nonlocal dilatation q and bulk modulus K [60]. Then, Eq. (12) can also be written as:

  Vm K q is ¼ iu exp RT

In this section, we introduce the new mechano-chemical peridynamic model for stress-assisted corrosion. According to Gutman's theory [14], under elastic deformation, the anodic current density is exponentially related to hydrostatic stress:

(13)

While this theoretical relationship gives reasonable approximations for changes in current density in some experiments [4,6], there are cases for which the theory significantly overestimates the stress influence [2]. A possible reason for this could be that the theory, being a simplified continuum model for homogenous material and smooth surfaces, it ignores local factors such as microscale heterogeneities, surface roughness, and other complex smallscale effects, which may play an important role. In order to obtain a model with a wider applicability, one can replace the theoretical parameters with empirical constants to be calibrated to experimental data [8]. mK Accordingly, we replace the scalar term VRT with an unknown constant g, found by calibration to experiments, and determine the stress dependency of anodic dissolution. g can be calibrated for each material from experimental tests. We use, therefore, the following constitutive relationship for kdiss in our stress-assisted corrosion PD model for the elastic deformation regime:

ksdiss ðxL ; xS Þ ¼ kudiss exp½gqðxS Þ 2.4. A coupled mechano-chemical peridynamic damage model

(12)

(14)

where xL and xS are respectively the liquid end and the solid end of the interfacial bond, ksdiss is the micro-dissolvality for the stressed material. kudiss is the micro-dissolvability for unstressed material and is calibrated to iu . This equation expresses that the microdissolvability of an interfacial bond is exponentially related to dilatation of the solid end of that bond. The nonlocal dilatation for a bond-based material in 3D is [38]:

S. Jafarzadeh et al. / Electrochimica Acta 323 (2019) 134795

ð

qðxÞ ¼ 3

Hx

y  yk  k b x  xkdVbx kb ð

(15)

x  xkdVbx kb

5

To validate our model, we first perform a novel experiment to observe and characterize the activation-controlled uniform corrosion of copper in NaCl solution under elastic deformation. We use this experimental results in section 4.1 to obtain stress-dependency coefficient in our model for corrosion of copper under loading.

Hx

Nonlocal dilatation is three times the average extension of mechanical bonds, and is analogous to volumetric strain in the local theory (which is three times the mean strain). In this model, dilatation is calculated for points inside the dissolving solid region. Since points in this region are located near the solid surface, their horizon does not cover only solid points (see Fig. 3). As a result, they suffer from a nonlocal phenomenon called the peridynamic surface effect [61], which means that points near the surface are slightly softer that the material in the bulk because of their incomplete horizon. To reduce this effect one can use the PD surface correction methods [61,62]. Note that Eq. (15) is for the 3D case. Nonlocal dilatation formula in 2D is similar to Eq. (15), and only differs in terms of the scalar coefficient [63]. In this study however, we introduce an alternative formula to approximate dilatation for 2D cases, that reduces the surface effect and is cheaper to compute [63].

ksdiss ðxL ; xS Þ ¼ kudiss expfga½smax ðxS Þ þ smin ðxS Þg

(16)

smax and smin are respectively the maximum and the minimum bond-strains among the mechanical bonds connected to the point xS . The constant a is 1 for plane strain and 12 for plane stress conditions. Derivation of this formula is provided in Appendix A. To simulate stress-assisted corrosion numerically for arbitrary geometries and boundary conditions, we use the algorithm shown in Fig. 4. The domain geometry, metal/electrolyte subdomains, D, d, kudiss and g are given as inputs. The details of the numerical implementation are provided in Appendix B. If the corrosion time steps are too small, and the computation of equilibrium at each time step is time consuming, one can solve the equilibrium and update bond stretches every several steps of corrosion to reduce the computational cost.

3. Stress-assisted corrosion in copper In this part, we introduce a new experimental design for stressassisted corrosion. The novelty in this experiment is to maintain the activation controlled uniform corrosion regime under application of stress. 3.1. Experimental procedure A commercial Cu plate (ThyssenKrupp Materials NA, Inc.) with a purity of 99.94% was used in stress-corrosion (SC) tests. The metallographic image via an optical microscope (Leica DM2700 M) in Fig. 5a indicates an isotropic microstructure with the grain size no larger than 50 mm on the observed direction. Dynamic polarization curve (scanning rate: 1.0 mV s1; from 2.0 to þ1.0 V) of Cu plate in 1M NaCl solution at room temperature was recorded by an electrochemical analyzer (CHI6062E, CH Instruments) using platinum foil and saturated calomel electrode (SCE) as the counter and the reference electrodes, respectively. The curve in Fig. 5b suggests that the Cu sample passivates in 1 M NaCl when the applied potential is between 0.35 and 0.1 V. To induce a uniform corrosion and prevent the pitting corrosion of Cu specimens during the SC tests, a positive potential of þ0.4 V was applied. Therefore, the recorded current divided by the corroding area can represent the true corrosion rate. If the corrosion is not uniform, the area includes passive surfaces, and consequently, the calculated current density will not reveal the actual anodic dissolution rate. Rectangular samples (80 mm  10 mm, Fig. 6a) were cut from the 0.635 mm thick Cu plate for the SC tests. Sample surfaces were finely ground to 800 grit with silicon carbide abrasive papers, followed by rinsing with ethanol to remove the debris and contaminants. To accurately measure the corrosion rate, only a circular area

Fig. 4. Graphical illustration of the coupled stress-assisted corrosion solver.

6

S. Jafarzadeh et al. / Electrochimica Acta 323 (2019) 134795

Fig. 5. a) Microstructure of Cu specimen etched by HCl/FeCl3/H2O/Glycerol solution; b) the dynamic polarization curve of a Cu specimen in 1 M NaCl solution.

(with diameter of 3.2 mm) at the center of the specimen was exposed to corrosive solutions. The rest of the immersed area was covered by Scotch tape to avoid corrosion. Before recording the corrosion current, the exposed surface was electrochemically cleaned by applying a bias of þ0.4 V (vs SCE) for 10 s.

In order to investigate the effect of tensile stress on the uniform corrosion of Cu, a spring-loaded stressing frame (Fig. 6b) was used to apply uniaxial tension [64]. Desired tensile stresses were applied through the compression of the loading spring, and the displacement of the loading spring was controlled by turning of nut at the

Fig. 6. a) Sketch of the Cu specimen. b) Mechanical drawing and c) Photo of the setup with a spring-loaded stressing frame used for stress corrosion experiment.

S. Jafarzadeh et al. / Electrochimica Acta 323 (2019) 134795

right side. The spring rate is 130.4 N/mm, which was experimentally calibrated through a material test system (MTS model 318.10). The applied tensile stresses in SC tests are lower than the yield strength (258.6 MPa) to avoid the plastic deformation in specimens. Meanwhile, 100 mL 1 M NaCl solution was circulated (60 mL/min) via two pumps to promote the activation-controlled corrosion mode by avoiding the local accumulation of Cu ions and domination of transport kinetics [65]. The input tube (inner diameter of 2.4 mm) for the flowing solution was placed over the corrosion spot to maximally flushing away the corrosion-produced species. The corrosion of Cu in NaCl solution was driven by applying a constant bias (þ0.4 V, vs SCE) using an electrochemical analyzer (CHI6062E, CH Instruments), and the corrosion current was simultaneously recorded. Three-electrode configuration was used, including Pt counter electrode, Cu working electrode and SCE reference electrode. 3.2. Experimental results and discussion First, we verify whether the designed set up omits the transport controlled effect and reveals the activation controlled process. This is performed by comparing the current densities of two experiments: one with flowing solution, and one with static solution (both in absence of stress). As shown in Fig. 7, the corrosion current density of Cu specimen in flowing NaCl solution is significantly different with the corrosion in static solution. In the first ~5 s, the recorded current includes both capacitive current and the current for breaking up the initial passivation surface. The corrosion current density for the Cu specimen immersed in the static solution decreases dramatically (up to ~20 s), indicating the diffusioncontrolled corrosion. However, the current density for the corrosion in flowing solution maintain a stable status at much higher level, suggesting that the corrosion is mainly controlled by activation mode. The motion of solution can prevent the accumulation of metal ions and consequently minimize the transport-control effect. This procedure enables us to estimate, quantitatively, the tensile stress effect on the activation-controlled corrosion. The experiment with flowing flow was repeated several times the standard deviation of about 4% was observed. We also needed to verify whether the desired uniform corrosion is achieved with the choice of þ0.4 V potential and flowing 1M NaCl

Fig. 7. Evolution of corrosion current densities in static and flowing solutions, under the same applied potential (þ0.4 V vs SCE).

7

solution. The uniform corrosion of Cu specimen has been experimentally verified by 3D surface profile analysis. A Keyence laser scanning microscope (VK-X200k) was used to measure the 3D surface profile of corrosion spot. Fig. 8 shows the profile of a corrosion spot (corroded in flowing solution for 300 s under the bias of þ0.4 V vs SCE, in absence of stress). The exposed circular region is almost uniformly corroded without noticed localization, with an average depth of 12 mm. 2D cross-sectional profiles (curves in Fig. 8) show that the central area is rougher than the edge of the corrosion spot, which may be induced by the variation in flowing rate from the edge to the center. In order to quantitatively identify the effect of elastic stress on corrosion current, the corrosion potential was set to the fixed value of þ0.4 V (vs SCE). Five different tensile stresses of 0, 87, 130, 174, 200 MPa were applied to the Cu samples for while recording corrosion currents. The experimental results are presented in Fig. 9. Fig. 9 shows that elastic tensile stress increases the anodic dissolution rate. This supports the mechano-chemical theories that suggest elastic stress can produce a shift in the chemical potential of electrodes, and increase the corrosion rates [14,21]. Note that during the first 3 s, the capacitive current, caused by ions accumulating in front of the electrodes, is dominant. Afterwards, the Faraday current, generated by the electrochemical reaction, increases and becomes dominant. Therefore, curves in both Figs. 7 and 9, exhibit “wave hollow” phenomenon over the first couple of seconds, as the transition between the capacitive and the Faraday currents takes place. The measurements presented in Fig. 9 are used in the next

Fig. 8. Uniform corrosion of Cu sample analyzed by the 3D profile measurement. The exposed circle was corroded for 300 s under the potential of þ0.4 V vs SCE, in the circulating NaCl solution. The photo of the corroded sample is shown next to the 3D profile.

8

S. Jafarzadeh et al. / Electrochimica Acta 323 (2019) 134795

Fig. 9. Effect of tensile elastic stress on activation controlled corrosion rate in Copper exposed to 1M NaCl flowing solution and under þ0.4 V(vs SCE) applied potential.

domain, discretized with uniform grid spacing of Dx ¼ 0:5 mm. The horizon size is 2 mm which is close to the typical thickness of the diffuse corrosion layer [31e35]. No-flux boundary conditions on sides and bottom of the domain is considered for corrosion model. On the top surface a d-thick layer of nodes is the initial electrolyte. Since the flowing electrolyte flashes the dissolved metal ions, we consider zero concentration in the whole electrolyte phase, and do not solve the diffusion in it for this activation controlled case. This acts as a boundary condition on top of the dissolving solid. Details of boundary conditions in PD corrosion models are found elsewhere [30]. In this study, for the simulations that involve stress, we apply displacement controlled boundary condition on the sides, such that the desired stress is induced in the solid domain. The amount of applied displacement to produce the desired stress is easily calculated from Hook's law. Fig. 10 shows schematics of the 2D domain and initial and boundary conditions for stress-assisted corrosion simulations. Concentration of intact solid is calculated from density and molar mass. For the tested material Csolid ¼ 141 M. From these simulations current density is calculated according to Ref. [30]:

Pn

i¼1 Cðx; t

iðtÞ ¼ zF section to calibrate and validate the introduced PD mechanochemical model for stress-assisted corrosion discussed in previous section. In order to have one value that represents the experimental current density generated from the electrochemical reaction for each stress level, the measured current densities plotted in Fig. 9, are averaged for values in the 30 s interval between 10 and 40 s of the measurement. The current densities for different tensile stresses are calculated and shown in Table 1. 4. Simulation results and discussion

4.1. Calibration of model parameters Two parameters in our model need to be calibrated: microdissolvability of unstressed material (kudiss ), and the stressdependency coefficient (g). This is done via two trial simulations and two measured current densities: one with stress and one for the unstressed case. For our simulations in this study we use a 2D 300  100 mm2

(17)

where z is the charge number for metal ions, F is the Faraday's constant, n is the total number of nodes, DVi is the volume of grid cell i, Dt is the time step, and A is the area of the corroding surface. To find kudiss we used its linear correlation with the current trial density. We perform a trial simulation with some ku; and diss calculate its current density itrial using Eq. (17). Given the experiu mentally measured current density for the unstressed material iu , we use the following relationship to calculate kudiss :

kudiss

In this part, we first use two of the measured current densities to calibrate our 2D model for stress-assisted corrosion in copper. Then, we use the calibrated model to predict the remaining measured corrosion rates in various stress levels. In the end we use the model for corrosion of a tapered plate with multiple holes under uniaxial tension, to show the applicability of the introduced computational model in general geometries and boundary conditions.

P þ DtÞDVi  ni¼1 Cðx; t  DtÞDVi Að2DtÞ

ku;trial diss

¼

iu itrial u

(18)

With this described computational set up, and knowing iu ¼ 1.45 mA mm2 from experiment, we calibrated kudiss to be 1.05  102 m1 s1 for copper. We also need to calibrate g. Since kdiss linearly correlates with i, given Eq. (16) one can write:

is ¼ iu expfag½smax ðxÞ þ smin ðxÞg;

(19)

For some dissolving solid point x on the surface. For certain applied stress s with a known current density is, we perform a trial simulation with a trial gtrial and a ¼ 0:5 in plane stress, which results in itrial current density. Then one can write: s

is trial is

¼

expfag½smax ðxÞ þ smin ðxÞg expfagtrial ½smax ðxÞ þ smin ðxÞg

¼ expfaðg  gtrial Þ½smax ðxÞ þ smin ðxÞg

(20)

Combining Eqs. (19) and (20) gives the relationship for g: Table 1 Current densities for each tensile stress level (averaged between 10 and 40 s in Fig. 9). 2

Stress/MPa

Current density/mA mm

0 87 130 174 200

1.45 1.53 1.59 1.74 1.82

For more accurate measurements, one can use Scanning electrochemical microscopy (SECM) [66] to locally measure the surface reactivity under various stresses.

. 3  is log itrial s g ¼ gtrial 41   trial . 5 iu log is 2

(21)

Using a trial simulation with the stress value 174 MPa and knowing the corresponding measured current density iu ¼ 1.74 mA mm2, we calibrated g to be 275.28 for copper. If the influence of stress is considered to be independent of the nature of the electrolyte and the potential, the calibrated model can be used to simulate stress assisted corrosion for copper under any elastic stress level.

S. Jafarzadeh et al. / Electrochimica Acta 323 (2019) 134795

9

Fig. 10. The 2D model set up, and the initial and boundary conditions for stress assisted corrosion. Here, u is the axial displacement, L is the length of the domain, and s is the desired tensile stress.

4.2. Model prediction In this part, we use the calibrated model to predict the corrosion rate for the stress levels: 87, 130, and 200 MPa. Using the same discretization and boundary conditions explained in previous section, five simulations were performed using the five tensile stress levels 0, 87, 130, 174, 200 MPa. Time step is Dt ¼ 0:1 s and for simulations with non-zero stress, a static solver updates the bond stretches after each 10 s of corrosion time. Current densities are then obtained using Eq. (17). Fig. 11 show the current densities from simulations in comparison to the measured values. We observe that the model successfully predicts the increase in current density with a reasonable accuracy. Trivially, the simulation matches the experiments for the stresses 0 and 174 MPa, since the model parameters are calibrated to these two tests. Fig. 12 shows evolution of the “uniform” stress-assisted corrosion damage under uniaxial tension of s ¼ 130 MPa at t ¼ 5, 10 and 15 min. As observed, the corroding surface in our simulation shows a micro-level roughness, which makes our simulations closer to real world observations. The under lying mechanism for this feature of the model is the stochastic procedure employed for breaking mechanical bonds to update the damage according to the values obtained from Eq. (11). Note that the main characteristic features in simulations like the current density are not affected by this stochasticity and remain deterministic [30]. The stochastic procedure, introduces a micro-level randomness in the model that accounts for the heterogeneities that are too small for the grid resolution to

Fig. 12. Evolution of uniform corrosion damage in peridynamic simulation of uniform stress-assisted corrosion.

include explicitly. Since after each (or a few) time steps of corrosion, the equilibrium is solved and bond stretches are updated for the new geometry, deformation is affected by the surface roughness, which is the case in actual corrosion. This is one of the unique advantages of our model for stress-assisted corrosion. Explicit representation of microstructural heterogeneities is possible only via micro-mechanical approaches [67] which could be highly expensive to compute, since they require extra fine grids or coupling with large scale models in multiscale frameworks [68,69]. The simulations presented in this part aimed for model verification and used a simple geometry. This model can be used for any desired geometry with complex stress distribution. 4.3. Example for a non-uniform stress distribution case

Fig. 11. Predicted current densities from simulation results and the measured values from experiments.

In this part, we use the calibrated model to simulate stressassisted corrosion for a trapezoidal plate with three inline holes, under uniaxial tension. In this computational test the edges of the holes are exposed to electrolyte. Fig. 13 shows the description of the geometry and boundary conditions in this problem. The left end is fixed in horizontal direction, and at the right end, a tensile loading is applied. Each hole experiences a different

10

S. Jafarzadeh et al. / Electrochimica Acta 323 (2019) 134795

Fig. 13. The 2D domain and the initial and boundary conditions for simulation of stress-assisted corrosion in a plate with three inline holes. Corrosion happens only inside the holes.

magnitude of a non-uniform stress distribution because of the tapered shape of the plate. This type of 2D setup is similar to actual cases where a metal plate with holes for fasteners undergoes crevice and galvanic corrosion inside the holes [70e73]. The simulation is performed using the same discretization and model parameters used in the uniform corrosion of copper in the previous section. Fig. 14 shows the corrosion damage progression in time, and the corresponding distribution of dilatation: As observed, concentration of negative hydrostatic pressure (tensile stress) leads to faster corrosion rate, while the positive pressure (compressive stress) slows down the dissolution. As a result, holes corrode faster in vertical direction. We also note that the severity of corrosion varies among the holes, since the plate is tapered, and consequently, the stress concentration is higher for the hole near the narrower end. Again, we observe the effect of the stochastic process employed in updating mechanical damage reflected in the slightly rough and non-symmetric corrosion front inside the holes. This is similar to the actual corrosion damage observed in experiments. This model can be easily extended to 3D for simulation of more

practical conditions. An advantage of PD corrosion model is that it could easily simulate stress corrosion cracking when a fracture criterion (a PD mechanical damage model) is used after solving the equilibrium equation [33,50]. Note that the stress-assisted corrosion model presented here is limited to elastic deformations (we assumed that only hydrostatic pressure can influence anodic dissolution). In order to extend the model to stress-assisted corrosion in the elasto-plastic regime, one needs to adopt a chemo-mechanical theory that accounts for the influence of plastic deformation on anodic dissolution. Further improvements can be made via incorporating the role played by cathodic reactions as well as the presence of corrosion products and passive films. To simulate circumstances where the mass-transfer effects are not negligible, PD diffusion or advection-diffusion equations could be solved in the liquid domain, to find the concentration of chemical species in the electrolyte and possibly impose their influence on the local corrosion rate (via modifying micro-dissolvability).

Fig. 14. Corrosion damage time evolution (top to bottom) in peridynamic simulation of stress-assisted corrosion in a tapered plate with holes. Corrosion happens only inside the holes: dilatation/hydrostatic pressure distribution (left pictures), and profiles of the corroded sample (right pictures), at different time steps, showing that highly stressed regions corrode faster.

S. Jafarzadeh et al. / Electrochimica Acta 323 (2019) 134795

5. Conclusions In this study, we introduced a mechano-chemical peridynamic corrosion damage model to simulate stress-assisted corrosion. This model considered an exponential dependency of anodic dissolution rate with nonlocal dilatation. A constant in the exponent determined the sensitivity of anodic current density to applied stress needs to be calibrated to the measured current density in an experimental test under application of some elastic stress. As an example, we provided a special stress-corrosion test on copper that produced conditions for activation-controlled uniform corrosion under applied tensile stress. The model parameters were calibrated using experimental data from these tests. We then used the calibrated model to predict the corrosion rate at other elastic stress levels. A simulation of stress-assisted corrosion for a more complex geometry was presented to show the capability and generality of the model. The new coupled peridynamic stress-assisted corrosion model showed the following advantages: propagation of damage was autonomous and explicit tracking of the corrosion front was not necessary; the diffuse corrosion layer with degraded mechanical properties was automatically captured as part of the solution method; micro-level surface roughness was reproduced; dependency of corrosion rate to elastic deformation was validated against experiments; simulation of stress corrosion cracking could be readily achieved with a fracture criteria added to the model.

11

horizon size is very small compared with the problem geometrical dimensions. As a result, nonlocal dilatation can be approximated by the classical volumetric strain:

q y ε1 þ ε2 þ ε3

A.1

where ε1 , ε2 ; and ε3 are the classical principal strains. In 2D:

qyaðε1 þ ε2 Þ

A.2

where a ¼ 1 for plane-strain condition, since ε3 ¼ 0. Under planestress conditions,

ε3 ¼ 

n 1n

ðε1 þ ε2 Þ;

A.3

where n is the Poisson's ratio. In 2D bond-based peridynamics n ¼ 13 , for plane-stress [38]. Consequently, by substituting A.3 in A.1, the constant a in A.2 for plane-stress in obtained as 12. Since d is very small compared with the domain's length scales, and assuming small elastic deformations, the disk-shaped horizon for each point undergoes an approximately homogeneous deformation and changes into an ellipse. Note that, for the points near discontinuities like cracks or the corrosion front, the horizon region changes from a disk-sector to an ellipse-sector. Homogeneous deformation of a generic point's horizon region is illustrated in Fig. A1. The principal strains ε1 and ε2 correspond to the directions of the long and short axes of the ellipse, respectively.

Fig. A.1. Horizon of a genetic point x before and after a homogeneous deformation and its principal axes.

Acknowledgments This work has been supported by the ONR project #N00014-151-2034 “SCC: the Importance of Damage Evolution in the Layer Affected by Corrosion” (program manager William Nickerson), by the AFOSR MURI Center for Materials Failure Prediction through Peridynamics (program managers Jaimie Tiley, David Stargel, Ali Sayir, Fariba Fahroo, and James Fillerup), and by a University of Nebraska-Lincoln System Science Award. The computations have been performed utilizing the Holland Computing Center of the University of Nebraska, which receives support from the Nebraska Research Initiative. Z.C. acknowledges the support from the Natural Science Foundation of China (No. 11802098) and the Natural Science Foundation of Hubei Province (No. 2018CFB111). Appendix A. Nonlocal dilatation approximation in 2D For derivation of Eq. (16) as the alternative formula for nonlocal dilatation in 2D, we refer to the classical local notion of dilatation. Note that the PD theory for elasticity and diffusion problems converges to the local classical theory as the horizon goes to zero [57e59,74]. In this study, similar to most other PD models, the

Then one can write:

ε1 ðxÞ y

y 1  yk  k b x 1  xk kb ¼ maxfsðx; b x Þg ¼ smax ðxÞ x 1  xk kb bx

A.4

ε2 ðxÞ y

x 2  xk y 2  yk  k b kb ¼ minfsðx; b x Þg ¼ smin ðxÞ x 2  xk kb bx

A.5

where b y 1 , and b y 2 are the positions of two points in the family of y along the principal directions 1 and 2, respectively. Position vectors b x 1 , and b x 2 are the points in the reference configuration that are mapped by the deformation into b y 1 and b y 2 , respectively. The scalars smax and smin are, respectively, the maximum and the minimum bond strains of all mechanical bonds at y. Using,Eqs. A.4 and A.5 2D nonlocal dilatation can be approximated by:

qya½smax ðxS Þ þ smin ðxS Þ

A.6

By comparing Eq. (A.6) with Eq. (15), we observe that the integral of strains over the horizon is replaced with finding the

12

S. Jafarzadeh et al. / Electrochimica Acta 323 (2019) 134795

maximum and the minimum of strains, which is less sensitive to surface effect and faster to compute. Note that if strains vary significantly over the horizon region scale (e.g. if horizon size is large compared with the domain dimensions, or in cases with large, localized deformations), and the assumption of homogenous deformation is no longer valid, nonlocal dilatation may be significantly different from the local dilatation. In such cases, the information provided by the maximum and minimum bond-strains is not sufficient to obtain the local or nonlocal dilatation values.

(B.3) for displacement vectors (u). Details are available in Refs. [44,77]. The discretized PD damage-dependent diffusion (Eq. (8)), with the meshfree method with one-point Gaussian quadrature, for a generic point xi and at time tn , is:

C_ni ¼

X n C nj  C ni kij   DVij xj  xi 2 j

The superscript n associate quantities with nth time step. knij for each bond is calculated from damage at its ends:

Appendix B. Discretization and numerical methods For spatial discretization, we use uniform grid spacing in both directions of the Cartesian coordinate system (see Fig. B1).

B.4

knij



dni ; dnj



¼

8 > > > <

kL ðDÞ 0

> u > > : kdiss exp½gqðxS Þ

; dni ¼ 1&dnj ¼ 1 ; dni < 1&dnj < 1 i h i h n ; di or dnj < 1& dni or dnj ¼ 1 B.5

qðxS Þ in Eq. (B.5) is calculated using Eq. (16), for the solid end of xi xj bond. Damage is computed from the discretized concentrationdependent damage equation:

dni ¼

8 > > > > > > > > <

1

Csolid  C ni > Csolid  Csat > > > > > > 0 > :

; if C ni  Csat ; if Csat < C ni < Csolid

B.6

; if C ni ¼ Csolid

We use the forward Euler method to update the concentration values for time step n þ 1: n C nþ1 ¼ C ni þ C_ i Dt i

Fig. B.1. The uniform grid spacing in 2D. In this study d ¼ 4:02  Dx.

To discretize bond-based PD mechanical equilibrium (Eq. (1)) we use a meshfree method with one-point Gaussian quadrature [44,75]:

      yj  yi   xj  xi  yj  yi   mij c DVij þ bi ¼ 0 xj  xi  yj  yi

X j

B.1

The subscripts i and j are used to associate quantities with the current node x, and with the family node b x respectively, in the discretized domain. DVij is the partial volume of the grid cell associated with node xj and covered by the horizon of node xi . For grid cells fully covered inside horizon DVij ¼ Dx2 . For grid cells on the horizon edge that are not fully located inside the circle (see Fig. B1), DVij is approximated via an efficient algorithm explained elsewhere [54,76]. mij in Eq. (B.1) is calculated for xi xj bond from Eq. (6). The discrete version of Eq. (7) for damage is:

P j di ¼ 1  P

mij

j1

B.2

Let uj ¼ yj  xj , and ui ¼ yi  xi be the displacements associated with xi and xj nodes, respectively (See Fig. 2). Eq. (B.1), in terms of displacements, is:

X

mij c

j

    uj  ui þ xj  xi   xj  xi    xj  xi 

u  ui þ xj  xi  j  uj  ui þ xj  xi  DVij þ bi

¼0 B.3 We use the nonlinear conjugate gradient method to solve Eq.

B.7

After updating concentrations, damage values are calculated using Eq. (B.6), and mij are updated using Eq. (B.2) and the stochastic procedure described in Refs. [26,30]. mij values are passed on to the mechanical solver when static solver is called. More details on the numerical solution to the PD corrosion model is available in Ref. [30]. Note that nonlocal computations in 3D tend to be expensive. For such problems, significant efficiencies can be gained via GPU-based or parallel computations (see Ref. [30]), or faster numerical methods such as the boundary-adapted spectral method (see Ref. [78]). References [1] S.J. Kim, H.G. Jung, K.Y. Kim, Effect of tensile stress in elastic and plastic range on hydrogen permeation of high-strength steel in sour environment, Electrochim. Acta 78 (2012) 139e146. [2] L.Y. Xu, Y.F. Cheng, An experimental investigation of corrosion of X100 pipeline steel under uniaxial elastic stress in a near-neutral pH solution, Corros. Sci. 59 (2012) 103e109. [3] L.Y. Xu, Y.F. Cheng, Corrosion of X100 pipeline steel under plastic strain in a neutral pH bicarbonate solution, Corros. Sci. 64 (2012) 145e152. [4] X.H. Wang, X.H. Tang, L.W. Wang, C. Wang, W.Q. Zhou, Synergistic effect of stray current and stress on corrosion of API X65 steel, J. Nat. Gas Sci. Eng. 21 (2014) 474e480. [5] Y. Zheng, Y. Li, J. Chen, Z. Zou, Effects of tensile and compressive deformation on corrosion behaviour of a MgeZn alloy, Corros. Sci. 90 (2015) 445e450. [6] M.Y. Bao, C.Q. Ren, M.Y. Lei, X. Wang, A. Singh, X.Y. Guo, Electrochemical behavior of tensile stressed P110 steel in CO2 environment, Corros. Sci. 112 (2016) 585e595. [7] T. Wu, M. Yan, J. Xu, Y. Liu, C. Sun, W. Ke, Mechano-chemical effect of pipeline steel in microbiological corrosion, Corros. Sci. 108 (2016) 160e168. [8] H.-Q. Yang, Q. Zhang, S.-S. Tu, Y. Wang, Y.-M. Li, Y. Huang, Effects of inhomogeneous elastic stress on corrosion behaviour of Q235 steel in 3.5% NaCl solution using a novel multi-channel electrode technique, Corros. Sci. 110 (2016) 1e14.

S. Jafarzadeh et al. / Electrochimica Acta 323 (2019) 134795 [9] S.J. Kim, Effect of the elastic tensile load on the electrochemical corrosion behavior and diffusible hydrogen content of ferritic steel in acidic environment, Int. J. Hydrogen Energy 42 (2017) 19367e19375. [10] O. Kuwazuru, K. Ode, M. Yamada, A.J. Kassab, E. Divo, Experimental and boundary element method study on the effect of stress on the polarization curve of cast aluminum alloy in sodium chloride solution, Corros. Sci. 132 (2018) 136e145. [11] X.H. Zhao, Y.R. Feng, S.W. Tang, J.X. Zhang, Electrochemical corrosion behavior of 15Cr-6Ni-2Mo stainless steel with/without stress under the coexistence of CO2 and H2S, Int. J. Electrochem. Sci. 13 (2018) 6296e6309. [12] G. Wu, P.M. Singh, Effect of elastic and plastic strain on general corrosion and metastable pitting of steels, in: CORROSION 2018, NACE International, 2018. [13] D.H. Xia, J.Q. Wang, Z.B. Qin, Z.M. Gao, Z. Wu, J.H. Wang, L.X. Yang, W.B. Hu, J.L. Luo, Sulfur induced corrosion (SIC) mechanism of steam generator (SG) tubing at micro scale: a critical review, Mater. Chem. Phys. 233 (2019) 133e140. [14] E.M. Gutman, Mechanochemistry of Solid Surfaces, World Scientific Publishing Company, 1994. €n, A.M. Zimer, E.C. Pereira, L.H. Mascaro, O.M. Bruno, B. De [15] P. Van der Weee Baets, Modeling pitting corrosion by means of a 3D discrete stochastic model, Corros. Sci. 82 (2014) 133e144. [16] H. Wang, E.-H. Han, Computational simulation of corrosion pit interactions under mechanochemical effects using a cellular automaton/finite element model, Corros. Sci. 103 (2016) 305e311. [17] O.O. Fatoba, R. Leiva-Garcia, S.V. Lishchuk, N.O. Larrosa, R. Akid, Simulation of stress-assisted localised corrosion using a cellular automaton finite element approach, Corros. Sci. 137 (2018) 83e97. [18] C. Lin, H. Ruan, S.-Q. Shi, Phase field study of mechanico-electrochemical corrosion, Electrochim. Acta 310 (2019) 240e255. [19] P. Marcus, Corrosion Mechanisms in Theory and Practice, CRC Press, 2011. [20] D.-H. Xia, Y. Behnamian, J.-L. Luo, Factors influencing sulfur induced corrosion on the secondary side in pressurized water reactors (PWRs), J. Electrochem. Soc. 166 (2019) C49eC64. [21] S. Sarkar, W. Aquino, Changes in electrodic reaction rates due to elastic stress and stress-induced surface patterns, Electrochim. Acta 111 (2013) 814e822. [22] S. Jafarzadeh, Z. Chen, F. Bobaru, Computational modeling of pitting corrosion, Corros. Rev. (2019). https://doi.org/10.1515/corrrev-2019-0049. , D. Fe ron, J. Badiali, [23] D. Di Caprio, C. Vautrin-Ul, J. Stafiej, J. Saunier, A. Chausse Morphology of corroded surfaces: contribution of cellular automaton modelling, Corros. Sci. 53 (2011) 418e425. [24] J. Stafiej, D. di Caprio, L. Bartosik, Corrosion-passivation processes in a cellular automata based simulation study, J. Supercomput. 65 (2013) 697e709. [25] B.P. Rusyn, R.V. Tors'ka, A.Y. Pokhmurs'kyi, Modeling of the evolution of corrosion pitting with the use of cellular automata, Mater. Sci. 50 (2015) 706e713. [26] Z.G. Chen, F. Bobaru, Peridynamic modeling of pitting corrosion damage, J. Mech. Phys. Solids 78 (2015) 352e381. [27] Z. Chen, G. Zhang, F. Bobaru, The influence of passive film damage on pitting corrosion, J. Electrochem. Soc. 163 (2016) C19eC24. [28] S. Jafarzadeh, Z. Chen, F. Bobaru, Peridynamic modeling of repassivation in pitting corrosion of stainless steel, Corrosion 74 (2018) 393e414. [29] S. Jafarzadeh, Z. Chen, F. Bobaru, Peridynamic modeling of intergranular corrosion damage, J. Electrochem. Soc. 165 (2018) C362eC374. [30] S. Jafarzadeh, Z. Chen, J. Zhao, F. Bobaru, Pitting, lacy covers, and pit merger in stainless steel: 3D peridynamic models, Corros. Sci. 150 (2019) 17e31. [31] S. Li, Z. Chen, F. Wang, B. Cui, L. Tan, F. Bobaru, Analysis of corrosion-induced diffusion layer in ZK60A magnesium alloy, J. Electrochem. Soc. 163 (2016) C784eC790. [32] N. Badwe, X. Chen, D.K. Schreiber, M.J. Olszta, N.R. Overman, E.K. Karasz, A.Y. Tse, S.M. Bruemmer, K. Sieradzki, Decoupling the role of stress and corrosion in the intergranular cracking of noble-metal alloys, Nat. Mater. 17 (2018) 887e893. [33] S. Li, Z. Chen, L. Tan, F. Bobaru, Corrosion-induced embrittlement in ZK60A Mg alloy, Mater. Sci. Eng. A 713 (2018) 7e17. [34] R. Vallabhaneni, T.J. Stannard, C.S. Kaira, N. Chawla, 3D X-ray microtomography and mechanical characterization of corrosion - induced damage in 7075 aluminium (Al) alloys, Corros. Sci. 139 (2018) 97e113. [35] D. Yavas, P. Mishra, A. Alshehri, P. Shrotriya, K.R. Hebert, A.F. Bastawros, Nanoindentation study of corrosion-induced grain boundary degradation in a pipeline steel, Electrochem. Commun. 88 (2018) 88e92. [36] B. Marzban, A Multiphysics Computational Framework for Understanding Cell and Microtissue Morphogenesis, University of Rhode Island, 2018. [37] S.A. Silling, R. Lehoucq, Peridynamic theory of solid mechanics, Adv. Appl. Mech. 44 (2010) 73e168. [38] F. Bobaru, J.T. Foster, P.H. Geubelle, S.A. Silling, Handbook of Peridynamic Modeling, CRC Press, 2016. [39] Y.D. Ha, F. Bobaru, Characteristics of dynamic brittle fracture captured with peridynamics, Eng. Fract. Mech. 78 (2011) 1156e1168. [40] D. Kamensky, M. Behzadinasab, J.T. Foster, Y. Bazilevs, Peridynamic modeling of frictional contact, J. Peridynamics and Nonlocal Model. (2019) 1e15. [41] J. Mehrmashhadi, Y. Tang, X. Zhao, Z. Xu, J.J. Pan, Q. Van Le, F. Bobaru, The effect of solder joint microstructure on the drop test failureda peridynamic analysis, IEEE Trans. Compon. Packag. Manuf. Technol. 9 (2019) 58e71. [42] M. Behzadinasab, J.T. Foster, The third Sandia Fracture Challenge: peridynamic blind prediction of ductile fracture characterization in additively

13

manufactured metal, Int. J. Fract. 218 (2019) 97e109. [43] Z. Xu, G. Zhang, Z. Chen, F. Bobaru, Elastic vortices and thermally-driven cracks in brittle materials with peridynamics, Int. J. Fract. 209 (2018) 203e222. [44] G. Zhang, Q. Le, A. Loghin, A. Subramaniyan, F. Bobaru, Validation of a peridynamic model for fatigue cracking, Eng. Fract. Mech. 162 (2016) 76e94. [45] F. Bobaru, J. Mehrmashadi, Z. Chen, S. Niazi, Intraply fracture in fiberreinforced composites: a peridynamic analysis, in: ASC 33rd Annual Technical Conference & 18th US-Japan Conference on Composite Materials, Seattle, 2018, p. 9. [46] J. Mehrmashhadi, Z. Chen, J. Zhao, F. Bobaru, A Stochastically Homogenized peridynamic model for Intraply fracture in Fiber-Reinforced composites, Compos. Sci. Technol. 182 (2019) 107770. [47] Z. Chen, S. Niazi, G. Zhang, F. Bobaru, Peridynamic Functionally Graded and Porous Materials: Modeling Fracture and Damage, Handbook of Nonlocal Continuum Mechanics for Materials and Structures, 2017, pp. 1e35. [48] M. Behzadinasab, T.J. Vogler, A.M. Peterson, R. Rahman, J.T. Foster, Peridynamics modeling of a shock wave perturbation decay experiment in granular materials with intra-granular fracture, J. Dyn. Behav. Mater. 4 (2018) 529e542. [49] Z. Chen, S. Niazi, F. Bobaru, A Peridynamic Model For Damage And Fracture In Porous Materials, Int. J. Rock Mech. Min. Sci. 122 (2019) 104059. [50] Z. Chen, S. Jafarzadeh, S. Li, F. Bobaru, Q. Qian, Peridynamic mechano-chemical modeling of stress corrosion cracking, in: IRF2018: 6th International Conference Integrity-Reliability-Failure, INEGI/FEUP, 2018, pp. 22e26. [51] S. Jafarzadeh, Z. Chen, J. Zhao, F. Bobaru, 3D Peridynamic Models for Pitting Corrosion and Stress Corrosion Cracking, Meeting Abstracts, 2019. MA201901 999. [52] S.A. Silling, Reformulation of elasticity theory for discontinuities and longrange forces, J. Mech. Phys. Solids 48 (2000) 175e209. [53] S.A. Silling, Introduction to peridynamics, in: F. Bobaru, J.T. Foster, P.H. Geubelle, S.A. Silling (Eds.), Handbook of Peridynamic Modeling, CRC Press, 2016, pp. 63e98. [54] F. Bobaru, G. Zhang, Why do cracks branch? A peridynamic investigation of dynamic brittle fracture, Int. J. Fract. 196 (2015) 59e98. [55] G. Zhang, G.A. Gazonas, F. Bobaru, Supershear damage propagation and subRayleigh crack growth from edge-on impact: a peridynamic analysis, Int. J. Impact Eng. 113 (2018) 73e87. [56] F. Bobaru, M. Duangpanya, A peridynamic formulation for transient heat conduction in bodies with evolving discontinuities, J. Comput. Phys. 231 (2012) 2764e2785. [57] J. Zhao, Z. Chen, J. Mehrmashhadi, F. Bobaru, Construction of a peridynamic model for transient advection-diffusion problems, Int. J. Heat Mass Transf. 126 (2018) 1253e1266. [58] Z. Chen, F. Bobaru, Selecting the kernel in a peridynamic formulation: a study for transient heat diffusion, Comput. Phys. Commun. 197 (2015) 51e60. [59] Z. Chen, D. Bakenhus, F. Bobaru, A constructive peridynamic kernel for elasticity, Comput. Methods Appl. Mech. Eng. 311 (2016) 356e373. [60] S.A. Silling, M. Epton, O. Weckner, J. Xu, E. Askari, Peridynamic states and constitutive modeling, J. Elast. 88 (2007) 151e184. [61] Q. Le, F. Bobaru, Surface corrections for peridynamic models in elasticity and fracture, Comput. Mech. 61 (2018) 499e518. [62] E. Madenci, M. Dorduncu, A. Barut, N. Phan, Weak form of peridynamics for nonlocal essential and natural boundary conditions, Comput. Methods Appl. Mech. Eng. 337 (2018) 598e631. [63] Q. Le, W. Chan, J. Schwartz, A two-dimensional ordinary, state-based peridynamic model for linearly elastic solids, Int. J. Numer. Methods Eng. 98 (2014) 547e561. [64] ASTM Standard G49-85 (reapproved 2011): Standard Practice for Preparation and Use of Direct Tension Stress-Corrosion Test Specimens. [65] T.P. Hoar, J.M. West, Mechano-chemical anodic dissolution, Nature 181 (1958), 835-835. [66] D.H. Xia, J.Q. Wang, Z. Wu, Z.B. Qin, L.K. Xu, W.B. Hu, Y. Behnamian, J.L. Luo, Sensing corrosion within an artificial defect in organic coating using SECM, Sens. Actuators B Chem. 280 (2019) 235e242. [67] K. Makarian, S. Santhanam, Utility of 2D finite element simulations for predicting effective thermomechanical properties of particle-reinforced composites, in: ASME 2018 International Mechanical Engineering Congress and Exposition, American Society of Mechanical Engineers, 2018. V009T012A047. [68] K. Zare Rami, Y.-R. Kim, M. Khedmati, G. Nsengiyumva, H. Alanazi, Two-way linked multiscale method integrated with nanomechanical tests and cohesive zone fracture to model highly heterogeneous binding materials, J. Eng. Mech. 144 (2018), 04018095. [69] T. You, Y.-R. Kim, K. Zare Rami, D.N. Little, Multiscale modeling of asphaltic pavements: comparison with field performance and parametric analysis of design variables, J. Transport. Eng., Part B: Pavements 144 (2018), 04018012. [70] D.C. Hansen, Metal corrosion in the human body: the ultimate bio-corrosion scenario, Electrochem. Soc. Interface 17 (2008) 31. [71] U.S. Department of Transportation, Federal aviation administeration, 2018. Advisory Circular No.43e4B. [72] W. Faes, S. Lecompte, Z.Y. Ahmed, J. Van Bael, R. Salenbien, K. Verbeken, M. De Paepe, Corrosion and corrosion prevention in heat exchangers, Corros. Rev. 37 (2019) 131e155. [73] G.S. Frankel, Corrosion Degradation of Coated Aluminum Alloy Systems through Galvanic Interactions, The Ohio State University, 2017.

14

S. Jafarzadeh et al. / Electrochimica Acta 323 (2019) 134795

[74] P.K. Jha, R. Lipton, Numerical convergence of nonlinear nonlocal continuum models to local elastodynamics, Int. J. Numer. Methods Eng. 114 (2018) 1389e1410. [75] S.A. Silling, E. Askari, A meshfree method based on the peridynamic model of solid mechanics, Comput. Struct. 83 (2005) 1526e1535. [76] F. Bobaru, Y.D. Ha, Adaptive refinement and multiscale modeling in 2d peridynamics, Int. J. Multiscale Com. 9 (2011) 635e659.

[77] Q. Van Le, Relationship between Microstructure and Mechanical Properties in Bi 2 Sr 2 CaCu 2 Ox Round Wires Using PeridynamicSimulation, North Carolina State University, 2014. [78] S. Jafarzadeh, A. Larios, F. Bobaru, Efficient Solutions for Nonlocal Diffusion Problems via Boundary-Adapted Spectral Methods, 2019 arXiv preprint arXiv: 1905.03875.