Surface irregularities as sources for corrosion fatigue

SURFACE IRREGULARITIES AS SOURCES FOR

CORROSION FATIGUE

Andrey P. Jivkov

Solid Mechanics, Malmö University, 205 06 Malmö, Sweden

Abstract

Corrosion fatigue crack nucleation from surface irregularity is modelled as a mov-ing boundary value problem. The model is based on material dissolution proportional to the surface stretch. Dissolution and re-passivation processes are forming the ge-ometry of the crack tip, thus creating con-ditions for strain concentration. No crack growth criterion is used. The interaction between the electrochemical processes and the deformation of the crack tip region is incorporated in continuum mechanical theory. Elastic-perfectly plastic materials under low frequency cyclic load are con-sidered. The model simulates how cracks form and grow in a single continuous process. The resulting natural variation of lengths of the formed cracks makes them grow with different rates. One crack after another falls into a wake behind a larger crack and the crack tip load of the smaller decreases leading to its arrest.

Keywords: C rack evolution, Fatigue cor-rosion, Moving boundary, Finite elements.

1. Introduction

Crack nucleation and short crack growth are important phases in fatigue corrosion (FC) and stress corrosion cracking (SCC). Pitting has been found to be a major mechanism for crack initiation. Experi-ments have shown that as a rule the largest pits are responsible for the emergence of cracks [1,2]. The formation of pits on a smooth surface is less understood, but is usually assigned to the interaction between surface defects and electrochemical reac-tions. While this assertion locates the ori-gin of pitting, it cannot answer the ques-tion about the size of the forming pit. It is

believed that the synergy between the strain energy, introduced into the solid by mechanical loading, and the potential of the environment-solid interface, is size determining for the formed pits. Such a concept has been used in studying mor-phology evolution of solid surfaces via surface diffusion of atoms [3,4]. In this case the interface potential coincides with the chemical potential of the surface. Un-der corrosion conditions, apart from the surface diffusion, there is another atomic-level mechanism, namely direct addition and removal of particles from and to the environment. The tendency of this mecha-nism to lower the system free energy is represented by an electrostatic addition to the chemical potential, which makes the interface potential electrochemical [5]. Both mechanisms have been taken into account in a recent work on chemical etch-ing [6]. The work shows a very good agreement between experimentally meas-ured surface morphology evolution and theoretically obtained wavelength for the surface roughness. The approach of [6], as a bi-product, predicts with reasonable ac-curacy pit sizes, reported in another recent experimental study dealing with corrosion pit formation under fatigue loading [7]. Since the surface diffusion mechanism is relevant for processes involving high stresses and temperatures and small size scales, the addition-removal mechanism is perhaps size determining in pit formation. For a developed corrosion crack, however, the surface diffusion may become increas-ingly important as will be discussed later in this paper.

Although the above energy approach may be successful in determining pit size, its form is not sufficient to explain crack formation from a present pit. The exis-tence of passive surface films, controlling

the access of the solution to a fresh metal surface, and corresponding repassivation processes have to be taken into account [8]. It is the interaction between mechanical loading, deforming or rupturing the pas-sive film, and repassivation, which defines the balance necessary for strain and disso-lution localization in the crack tip region [9]. A strong support for the key role of this interaction is provided by recent ex-perimental studies [10,11], showing that active loading is an essential prerequisite for continuing SCC. While the chemical composition, structure and growth of pas-sive films are now well understood [12], there are very few works that analyse the effect of deformation on films properties. This is probably because of the thinness of the film, typically 1-3 nm. An early pro-posal states that the deformed film rup-tures when a critical tensile strain is reached and uncovers bare metal that sub-sequently corrodes [13]. This view re-quires brittleness of the film, a feature that is not always observed. What could be said with confidence for non-brittle films is that the conductive properties of the film are significantly modified by surface deformation [14,15] in a direction facilitat-ing dissolution for tension and in the op-posite direction for compression. The ki-netics of charge transfer at the environ-ment-solid interface was found unaltered by deformation [15]. This means that a tensile surface deformation increases the dissolved volume but without changing the kinetics of the dissolution reaction. Such an observation allows one to think of the deformation-dissolution interaction in an unified way from mechanical point of view – the surface strain leads to increased material loss irrespective of whether the increase is due to film rupture or higher film conductance.

It must be further recognised that the ob-served corrosion cracks do not have atomically sharp tips and hence possess intrinsic length scale – the size of the dis-solution process region or crack width. A model for strain-driven corrosion growth of cracks with such realistic geometries has been proposed in [16]. The physical reasoning behind the model is shortly de-scribed in the next section. In mechanical terms, the model leads to a moving bound-ary value problem. It has been later sug-gested that strain-driven dissolution is responsible for the pit development and

short crack growth [17]. The moving boundary approach of [16] has been used in [17] to study crack nucleation and com-petition from a single or groups of surface pits in an elastic material. This work ex-tends the study to elastic-perfectly plastic materials and shows results characteristic for them.

2. Strain-driven dissolution

A solid body with a surface in contact with a corrosive environment is considered. A passive film is assumed to cover the sur-face at rest. The rate of dissolution in the presence of intact film and at rest is about three orders of magnitude lower than the rate for bare metal [14,15] and will be ne-glected in the proposed model. Hence, for fully covered surface the dissolution is not allowed. Upon load application the film deforms together with the bulk and as a measure of this deformation is chosen the strain parallel to the surface, defined as ε = (ds - dS) / ds, where dS and ds are the lengths of a differential surface arc ele-ment before and after the deformation, respectively. The deformation causes en-hanced dissolution of metal ions from the bulk to the environment, which is viewed as surface evolution in normal direction. The model proposed in [16,17] states that the rate of the surface evolution, R, is pro-portional to the surface strain, ε, via

(

f

)

d H ε

R

R= ε ε− , (1) Where Rd is the dissolution rate obtained

for an unstrained bare surface, the strain εf

is a characteristic property of the film and

H(x) is the Heaviside step function, i.e. H(x) = 1 if x > 0, and H(x) = 0 otherwise.

In cases of brittle films, the characteristic strain εf, has been identified as the strain at

film rupture, given in the literature as εf =

8÷10×10-4 [18]. The physical reasoning for (1) in such cases lies in smearing dis-crete fracture instances along the surface. If the strain at a number of surface points reaches the rupture strain, εf, the film

breaks at these points and a corresponding number of film fragments separated by gaps of bare metal are created. The film fragments are supposed to relax. The un-protected gaps are exposed to the corro-sive solution and the metal dissolves,

ad-vancing that portion of the surface. As-suming a continuous formulation of this process leads to (1). In cases of non-brittle films, the strain εf, may be attributed a

meaning of threshold in accordance with the experimental results [14], where a stress level was found below which no changes in the film conductance were ob-served and the film was as protective as in the absence of stress. The reported stress level was 20% of the yield strength of the considered material (steel 316L), which gives a rough estimate of the characteristic strain in (1) as εf = 2.5×10-4. Note, that

irrespective of the physical process, the evolution rate at load application is limited from above by the bare surface dissolution rate, because of the chosen strain measure with the property ε < 1. Although the lin-ear relation (1) has no experimental sup-port, it is believed to be a good first ap-proximation representing the general trend in the deformation-dissolution relation and sufficient to show the contribution of the mechanics to crack formation stage.

When the film is deformed to either rup-ture or changed conductivity, a repassiva-tion process starts, which in short time restores the film properties of an intact film at rest [19], i.e. interrupts the active dissolution. This is the main reason for the need of active loading as asserted in [10,11]. The repassivation time may vary from a couple of seconds to a couple of minutes, depending on the solid-environment composition. The repassiva-tion process causes an exponential de-crease of the dissolution rate, Rd, and

hence the surface evolution rate, R, with

time [20]. In order to avoid an explicit time dependence of the results, a simpli-fied integral approach is adopted here, where the load is thought of as applied in cycles of rectangular shape, i.e. quasi-statically for each cycle. The cycle period is supposed sufficiently large to allow full restoration of the film before the next cy-cle is applied. Given peak duration of a load cycle, T0, one can integrate the

par-ticular repassivation law (e.g. the one given in [20]) over T0 to obtain the

aver-age surface advance per cycle as shown in [17]. In this way the rate equation (1) is translated to a finite difference equation, where Rd should be understood as the

ad-vance that would be obtained for time T0 if

the dissolution-repassivation processes start from a bare metal surface, and R as

the apparent surface advance resulting from the initial presence of intact film and subsequently introduced surface strain, ε.

3. Problems formulation and solution Consider a plane body occupying the re-gion 0 ≤X1 ≤B and X2 ≤ B with respect

to a fixed coordinate system (X1, X2). The

body surface at X1 = 0, X2 ≤ B is in

con-tact with a corrosive environment and ini-tially covered by a passive film. The rup-ture strain of the passive film is chosen as εf = 0.001. On the contact surface, a single

pit or several competing pits in form of waves are taken. Fig. 1 illustrates the ge-ometry for the general case. The width of every pit is chosen as W = 10-3B and the

depth D = 10-4B. The material of the body

is supposed elastic-perfectly plastic with elastic modulus E = 206 GPa, Poisson’s

ratio ν = 0.3, and yield strength σy. The

body is considered to be in plane strain. The components of the displacement and the traction vectors are denoted by U1 and

U2, and T1 and T2, respectively. The

boundary conditions of the quasi-static problem of one cycle are given as: T1 = T2

= 0 and (1) along X1 = 0 (including pits);

U1 = T2 = 0 along X1 = B; T1 = 0, U2 = u

along X2 = B; T1 = 0, U2 = -u along X2 = -B,

where u is the magnitude of the applied

displacement during the cycle peak. It is chosen as u = εf B, which ensures a

con-stant strain field in a rectangular body without surface flaws, which is exactly on the threshold to break the protective film of the flat surface.

Figure 1: Wavy surface geometry.

Note, that the problem stated above is a moving boundary value problem. The evo-lution rate (1) constitutes an additional

W W X1 X2 D E, ν, σy εf

boundary condition, which compels the boundary to evolve, such that the redis-tributed stresses are maintained in an in-stantaneous equilibrium with the applied loads. The equilibrium part of the problem at each load cycle is solved using a com-mercial FEM-based software [21]. Con-stant strain triangular finite elements are used in the analysis. An in-house proce-dure for surface evolution tracking and geometry re-meshing is used. The evolv-ing parts of the surface are described usevolv-ing cubic B-splines. These are reshaped after each cycle by displacing the surface nodes according to (1) and re-meshed to accu-rately describe the new shape. The value of the dissolution rate used in the calcula-tions is rd = 1. The re-meshing process is

based on a length criterion – minimal edge size, emin; and on a curvature criterion –

maximal angle between two consecutive edges, γmax. These two parameters

deter-mine entirely the geometry of the forming corrosion crack, an observation that will be an important topic of the discussion section. In the simulations, the values emin

= 5×10-6B and γmax = 0.15 are used. The

interior is re-meshed by a Delauney-type triangulation procedure. The residual stresses and plastic strains, obtained after a load cycle closing, have been neglected when defining the boundary value prob-lem for the next cycle.

4. Results

Simulations have been performed for a series of material yield strengths, such that σ∞_/σ

y = 0, 0.025, 0.05, 0.075, 0.1, 0.15,

0.2, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5, where the remotely applied stress is σ∞ = Eεf. For

edge cracks, the load interval σ∞/σy < 0.32

defines small scale yielding (SSY) condi-tions according to ASTM convention. Cracks evolutions for σ∞/σy = 0.2 are

shown in Fig. 2 for the cases of single (a), double (b), and triple (c) pits. In picture (a), the profiles of the initial pit, the surface at crack incubation (to be shortly defined), and the crack at extension 1.5W are

pre-sented. The crack length is denoted by a,

and the crack width in the tip region by 2ρ. The incubation length is ainc. In picture (b)

the initial and incubation profiles are shown as well as the geometry at arrest of one of the cracks. The arrest length is aarr.

Picture (c) shows the profiles at arrests of the first and the second cracks, together

with the initial pit. All lengths depicted in the pictures are found nearly independent of the material yield strength. The cracks are found to attain width of 2ρ ≈ 15emin,

independently of the initial pit width, W.

This result is clarified in the discussion.

Figure 2: Cracks evolving from single (a), double (b) and triple (c) surface pits.

Fig. 3 shows the dependence of the nor-malised crack growth rate, R / Rd, on the

square root of the normalised crack exten-sion, a / W, for three values of the yield

strength, σ∞/σy = 0, 0.15, 0.3, and for a

crack emerging from a single pit. Note, that the first value represents elastic mate-rial. The plots help to unambiguously de-fine the crack incubation point as the onset of (approximately) linear relation between growth rate and square root of crack length. This is also the point where the evolving pit surface attains a curvature at

2ρ a = 1.5W ainc≈ 10ρ W

a

2W ainc≈ 10ρ aarr≈ 2W a ≈ 4W

b

4W

c

aarr≈ 4W a ≈ 8W aarr < 2W

the tip (approximately) equal to the one of the subsequently propagating crack. The cracks in all studied cases maintained nearly constant width, 2ρ, during the first stage of propagation shown in Fig. 3.

Figure 3: Crack growth rates vs. square root of crack extension.

At the same time, the crack growth rate after incubation shows strong dependence on the material yield strength. This sug-gests that shape changes of the crack tip region occur. To demonstrate the extent of these changes an appropriate scaling factor is used. Consider the elastic solution for a deep slender notch with circular tip region of radius ρ, e.g. [22]. The maximum strain at the tip is given by:

ρ

π

ε

E K_I tip 2 = , (2) where KI is the stress intensity factor for a

crack of same length as the notch. In the present contexts

ε

tip = R / Rd,

σ

∞ = E

ε

f.

Further, edge cracks are of concern. Using these observations, define a dimensionless shape factor as

(

)

a R R f d

ρ

ε

φ

243 . 2 = . (3) The dependence of the shape factor on the yield strength is presented in Fig. 4. Note from (2), that

φ

= 1 if the notch has exactly circular tip shape and the material is elas-tic. Note also that for a given yield strength,

φ

is proportional to the slope of the corresponding line in Fig. 3 represent-ing crack propagation stage.

Figure 4: Shape factor vs. yield strength. The results shown hereto in Figs. 2a, 3, 4 represent the cracks behaviour for short extensions. What is meant by short will be clarified next. If the crack growth is fol-lowed for sufficiently long time, the inter-action between the corrosion and the large plastic deformations cause blunting of the tip region. This blunting continues until the crack tip region attains a shape that is maintained self-similar in the continuation of the propagation. The self-similar growth is characterised by an approxi-mately constant rate. The overall growth rate behaviour is illustrated in Fig. 5 for the material with

σ

∞/

σ

y = 0.2.

Figure 5: Extended growth behaviour. The blunting and the stage of self-similar growth were observed for a non-linear elastic material [17], but were not quanti-tatively characterised there. For the cases of elastic-perfectly plastic materials stud-ied here, the onset of blunting was found to be approximately determined by the relation (a /

ρ

) (

σ

∞2/E

σ

y) = 0.025. If short

crack growth is defined as the KI

-controlled stage, i.e. where growth rate is proportional to square root of crack length,

0,4 0,6 0,8 1,0 0,000 0,005 0,010 0,015 0,020 σ∝ / σy = 0 σ∝ / σ_y = 0.15 σ∝ / σ_y = 0.3 Propagation Incubation R / R d (a / W)1/2 0,0 0,5 1,0 1,5 2,0 2,5 3,0 0,00 0,01 0,02 0,03 0,04 0,05 0,06

Short crack growth

Blu n ting Sel f-s imilar g rowth In cubatio n R / R d (a / W)1/2 0,0 0,1 0,2 0,3 0,4 0,5 1,0 1,2 1,4 1,6 1,8 2,0 2,2 SSY limit

φ

σ∝ / σ_y

then the extend of this growth is depend-ing both on the material yield strength and on the local geometry via the corrosion crack width. In terms of crack opening displacement,

δ

, for a corresponding ideal crack in a respective material, this relation translates to

δ

≈ 0.05

ρ

.

Finally, Fig. 6 shows crack extensions as functions of number of cycles (or time) for the cases of double (a) and triple (b) pits in a material with

σ

∞/

σ

y = 0.2. The crack

length is scaled with the initial pit size, and the number of cycles is scaled with ratio between the accepted in calculations dissolution rate and the electrochemical one. For comparison, the extension of a single crack in the same material is shown with dashed line. As already mentioned in connection with Fig. 2, the lengths scales shown are found independent of the mate-rial yield strength. Time scale is the only one affected by the material properties.

Figure 6: Cracks extensions vs. time for double (a) and triple (b) initial pits.

5. Discussion

The present study aimed at characterising the influence of the material yield strength on the initiation, short propagation and

competition of fatigue corrosion cracks. The moving boundary model of FC is in-trinsically connected to length and time scales, which determine various stages of evolution and competition. The time di-mension is introduced via the electro-chemically determined dissolution rate, Rd,

and is further influenced by the material properties or load level, i.e. by the ratio

σ

∞_/

_σ

y. This influence is outside the scope

of the present work. Therefore, the discus-sion does not address questions about when various events happen in time. The controlling length parameters are of major concern. Different length scales were found to determine crack evolution and crack competition.

Three stages of corrosion crack evolution are observed in Fig. 5. The approximately linear relation between crack growth rate and square root of crack extension shows that the stress and strain fields surrounding the crack tip are KI-controlled. The

inter-val under KI-control is called short crack

growth. The incubation stage is then de-fined as preceding the short crack growth. During incubation, localisation of strains and dissolution at the tip occur leading to the fully formed crack. A stage, character-ised by large deformations in the tip re-gion succeeds the short crack growth. In this non-linear regime, crack tip blunting and self-similar growth are then observed. The major parameter controlling evolution is the crack width, 2

ρ

. It governs the incu-bation length and the short crack growth limit via the relations ainc ≈ 10

ρ

and

δ

≈

0.05

ρ

, respectively, as visible from the results presented in Figs. 3, 4, 5. An addi-tional parameter is the curvature at the tip, which may lead to about double as high crack growth rate for very high loads or very soft plastic materials, as for linear elastic material, Fig. 4. For the latter case, the crack tip was found to attain nearly circular shape at incubation and to main-tain this shape in the propagation stage. Hence, the plasticity has an effect of sharpening the crack tip region, which increases with decreasing material yield strength.

Crack competition, opposite to evolution, is primarily controlled by the distance be-tween the existing pits. As visible in Figs. 2, 6 the competition is characterised by initial repulsion between the growing cracks, which leads to bended crack paths

0,0 0,1 0,2 0,3 0,4 0,5 0 2 4 6 8 10 Single crack a Arrest a_arr / W a / W a / W N (r_d / R_d ) 0,0 0,1 0,2 0,3 0,4 0,5 0 2 4 6 8 10 Single crack b Arrest a_arr / W a / W a_arr / W a / W N ( r_d / R_d )

and lower crack growth rates than the one for a single crack. The inhomogeneity of the mesh used in the FEM analysis gives priority to one of the cracks to initially grow slightly faster. In reality, there may be various microstructural causes that would produce similar priority. The rate difference increases with increased crack extension and leads to eventual arrest of the slower crack. At arrest, the continuing crack returns to its straight propagation as a Mode I crack, i.e. perpendicularly to the remotely applied load. For all the materi-als studied, the arrested crack has an ex-tension approximately equal to the dis-tance to the continually propagating crack. The latter has already grown to a double as big extension. From these observations, one may conclude, that the propagating cracks arrest all shorter cracks that fall into the sectors of right angle centred at the propagating crack tips. This result may be used to statistically simulate distribu-tions of arrested cracks along a solid sur-face [17].

In the study, the crack width and the crack tip radius of curvature were entirely de-termined by the criteria for tracking the evolving surface. Thus the control parame-ters, emin and

γ

max, set the limits for

ρ

and

for accuracy of shape description. This was confirmed by additional simulations with varying control parameters and same initial geometry as in Fig. 1. The results show independence of the crack width as well as the other lengths discussed in the previous paragraph on the initial pit ge-ometry. What the initial geometry does affect is the time scale during the incuba-tion stage [17]. Hence, the forming crack tends to minimise its width as much as the controlling parameters allow. This is rec-ognised as a valuable feature of the mov-ing boundary approach, because it leads to cracks attaining realistic (in the sense of corrosion cracks) geometrical shapes. An example of a transdendritic fatigue corro-sion crack in a pressure vessel component is given in Fig. 7.

In reality the width of the forming corro-sion crack may be set by materials micro-structure or by additional physical proc-esses. In the former case, the grain bound-ary thickness might be the determining length parameter, which is typical for in-tergranular corrosion cracking. When the corrosion cracking is not of this type as in

Fig. 7, the latter case should be considered. As reviewed in the introduction, a candi-date for a physical process, not included in the present model, is the diffusion of at-oms along the solid surface. While the interaction between the surface deforma-tion and the dissoludeforma-tion and passivadeforma-tion processes tends to minimise the forming crack width, the surface diffusion tends to operate in the opposite direction by bring-ing material particles from less to more strained surface regions. One may only speculate at this point, that the balance between these two tendencies determines the size of the nucleated crack, as the role of the surface diffusion increases with the increasing surface strains in the tip region. Further, one may speculate that the same balance is in reality not allowing the cor-rosion cracks to blunt and grow self-similarly, as was the result of the present work after the short crack growth.

Figure 7: Near-tip region of FC crack. (With permission from Vattenfall, Sweden)

6. Conclusions

The proposed moving boundary value ap-proach for modelling fatigue corrosion and its numerical solution allow for studying various evolution and competition parame-ters with sufficient accuracy. Corrosion crack evolution is primarily controlled by the width of the formed crack, 2

ρ

. The crack extension at incubation, given by ainc

≈ 10

ρ

, is independent of the initial geome-try and of the material properties. Growth under KI-control follows the incubation

stage until the crack extension reaches the value a ≈ 0.025

ρ

(E

σ

y/

σ

∞2). For larger

crack extensions the used model predicts crack tip blunting and self-similar growth

under constant rate. The distance between the cracks controls crack competition. Cracks arrest when they fall into a shadow zone of a larger propagating crack. This zone is approximately the right angular sector centred at the growing crack tip. Further work is required to examine the physics determining the size of the formed corrosion cracks.

Acknowledgements

The financial support of the Swedish Cen-tre for Nuclear Technology (SKC) and the Knowledge Foundation (KKS) is highly appreciated.

References

[1] Chen, G.S., Wan, K.-C., Gao, M., Wei, R.P., Flournoy, T.H. (1996), Mat. Sci. Eng., A219, p.126.

[2] Rokhlin, S.I., Kim, J.Y., Nagy, H., Zoofan, B. (1999), Eng. Fract. Mech., 62, No.4/5, p.425.

[3] Asaro, R.J., Tiller, W.A. (1972), Metall. Trans., 3, p.1789.

[4] Freund, L.B. (1995), Int. J. Solids Struct., 32, No.6/7, p.911.

[5] Linford, R.G. (1978), Chem. Rev., 78, No.2, p.81.

[6] Kim, K.-S., Hurtado, J.A., Tan, H. (1999), Phys. Rev. Lett., 83, No.19, p.3872. [7] Nakai, Y., Shimizu, Y., Fujiwara, S., Ogawa, T. (2001), In: Advances in Frac-ture Research, Proceedings of the 10th In-ternational Conference on Fracture (on CD), Ravi-Chandar, K., Karihaloo, B.L., Kishi, T., Ritchie, R.O., Yokobori Jr., A.R. and Yokobori, T (Eds), Pergamon. [8] Scully, J.C. (1979), In: Mechanisms of Environment Sensitive Cracking of Mate-rials, Proceedings of Metals Society Con-ference, pp.1-17, Swann, P.R., Ford, F.P. and Westwood, A.R.C. (Eds), The Metals Society, London.

[9] Parkins, R.N. (1992), JOM-J. Min. Met. Mat. S., 44, No.12, p.12.

[10] Kussmaul, K., Blind, D., Läpple, V. (1997), Nucl. Eng. Des., 168, No.1-3, p.53. [11] Heldt J., Seifert, H.P. (1999), In: Ninth International Symposium on Envi-ronmental Degradation of Materials in Nuclear Power Systems – Water Reactors, pp.901-909, Ford, F.P., Bruemmer, S.M. and Was, G.S. (Eds), The Minerals, Met-als & MateriMet-als Society (TMS).

[12] Olsson, C.-O.A., Landolt, D. (2003), Electochim. Acta, in press.

[13] Ford, F.P. (1996), Corrosion, 52, No.5, p.375.

[14] Yang, Q., Luo, J.L. (2001), Electro-chim. Acta, 46, p.851.

[15] Vignal, V., Oltra, R., Verneau, M., Coudreuse, L. (2001), Mat. Sci. Eng., A303, p.173.

[16] Jivkov, A.P., Ståhle, P. (2002), Eng. Fract. Mech., 69, No.18, p.2095.

[17] Jivkov, A.P. (2002), Research Report MUMAT2002:5, Malmö University, to appear in Theor. Appl. Fract. Mech.

[18] Engelhardt, G.R., Macdonald, D.D., Urquidi-Macdonald, M. (1999), Corros. Sci., 41, No.12, p.2267.

[19] Macdonald, D.D. (1999), Pure Appl. Chem., 71, No.6, p.951.

[20] Engelhardt, G., Urquidi-Macdonald, M., Macdonald, D.D. (1997), Corros. Sci.,

39, No.3, p.419.

[21] ABAQUS User’s Manual, Version 6.3 (2002), Hibbitt, Karlsson & Sorensen Inc.

[22] Tada, H., Paris, P.C ., Irwin, G.R. (2000), The stress analysis of cracks handbook, 3d Ed. ASME Press, New York.