Progress and Opportunities in Computational Modeling of Localized Corrosion

Most metals exist in natural environments in the form of compounds, such as oxides, carbonates, and sulfides. Through smelting technology, people consume energy to extract metal materials with excellent performance from these compounds, and people utilize metallic materials in all aspects of life. The extracted metal materials are comparatively unstable with higher energy than their natural compound. Corrosion is a process in which a chemical reaction occurs between a metal and its environment (usually the natural environment) and slowly returns to its stable state from an unstable state. In this process, the structural materials may deteriorate, which can cause severe harm to society’s development.

First of all, corrosion will directly lead to a decline in the life of metal materials, resulting in tremendous economic losses. The annual economic loss due to corrosion in China reached 2.1 trillion yuan (2014), accounting for about 3.34% of the gross domestic product (GDP) of that year.¹  The annual economic loss caused by corrosion in the world is about 2.5 trillion US dollars (2013).²  The impact of corrosion is not just measurable as a share of GDP. Corrosion causes enormous waste of resources and economic losses, pollutes the environment seriously, and even leads to fire and explosion, resulting in heavy casualties. A typical instance is the leakage of radioactive water from the Fukushima nuclear power plant in 2013 due to corrosion.³  It will take 40 y to completely repair the environmental contamination caused by the incident. On November 22, 2013, the Donghuang oil pipeline of Sinopec in Qingdao leaked and exploded, resulting in 62 deaths, 136 injuries, and a direct economic loss of nearly 752 million yuan. The subsequent investigation found that the direct cause of the accident was the corrosion at the junction of the oil pipeline and the drainage canal, and the pipeline was thinned and cracked, leading to crude oil leakage.⁴  On August 1, 2014, a gas explosion occurred in Kaohsiung, Taiwan due to the corrosion of old pipelines.⁵  In addition, corrosion will also affect the appearance and use of civil architecture, scenic spots, and daily necessities, and poses a serious threat to human health and travel safety. Corrosion also has a serious impact on cultural relics, historical sites, and national security. With the rapid development of China’s infrastructure, especially the continuous advancement of the national marine strategy, China has paid increasing attention to the problem of material corrosion. In fact, corrosion has become a stumbling block for the implementation of China’s “One Belt, One Road” strategy, and the regions along the “21st Century Maritime Silk Road” such as Southeast Asia and the coast of the Indian Ocean are facing serious corrosion problems. In a recently published report, Wu summarized and analyzed in detail the intense impact of corrosion on China’s economy, environment, and national defense.⁶  In May 2015, the “Made in China 2025” issued by the State Council of China clearly pointed out 10 key development areas, many of which involved corrosion research, such as new materials, marine engineering equipment and high-tech ships, energy-saving and environmental protection equipment and resource recycling, etc.⁷ 

Benefiting from the rapid development of various high-precision material characterization methods, a large number of experimental observations have been made on the evolution of various corrosion damages at different scales. Due to the complexity of the corrosion process, the current cognizance of the corrosion evolution process and the simulation and prediction of corrosion damage are still challenging. The National Research Council of the United States published “Research Opportunities in the Field of Corrosion Science and Engineering” in 2011. The book summarizes four major scientific research opportunities related to corrosion that can also be said to be four major challenges:⁸  (1) development of low-cost and environmentally friendly corrosion-resistant materials and coatings; (2) high-fidelity modeling of corrosion evolution in service environments; (3) quantitatively link accelerated corrosion tests under laboratory conditions with long-term corrosion observations in service environments; and (4) accurately predict the remaining service life. Among the four challenges, items (2) and (4) are in connection with corrosion modeling directly. The establishment of corrosion models will also greatly accelerate the development of corrosion-resistant materials and coatings. Additionally, they will provide insight into the viability and advantages of various corrosion mitigation methods.

In this paper, several modeling methods of localized corrosion in the last 10 y are introduced, including cellular automata (CA), finite element method (FEM), phase field (PF) model, and peridynamics (PD) model. The strengths and weaknesses of each method are analyzed in this present review. Additionally, several difficulties that need to be solved are pointed out, providing guidance for further research on localized corrosion simulation.

The establishment of corrosion models is required to combine the understanding of the corrosion mechanism, the knowledge of the materials in the specific application environment, and the calibration of the model parameters with the corrosion test data. The purpose of corrosion modeling is to predict the corrosion evolution process of materials in a specific application environment and the resulting changes in other properties and to reasonably estimate the expected corrosion behavior of materials. Over the past two decades, researchers’ progress in corrosion models has made the prediction of corrosion damage evolution a realistic possibility in the near future. In recent years, the experimental characterization of corrosion at different scales has deepened the understanding of corrosion mechanisms. At the same time, with the further enhancement of computing power, modeling research, and numerical prediction of corrosion have received more and more attention and development. The following introduces some important localized corrosion models and their corresponding results.

CA, proposed by mathematicians Stanislaw M. Ulam and John von Neumann, is a mathematical model to describe complex phenomena in nature.⁹  CA technique has been used in simulation studies in various fields, including public transport systems,¹⁰  the spread of forest fires,¹¹  and the growth of cells.¹²  Because the CA model is on the basis of the localized reaction rules of the system, it can reflect the influence of complex physical and chemical systems at different scales, including the molecular or atomic scale, and the properties of these systems can be qualitatively described in the macroscopic scale. It is more intuitive and convenient to model the corrosion evolution process using this technology. In recent years, a growing number of scholars have used CA to conduct microscopic modeling of corrosion and realize the numerical simulation of the metal corrosion process. The application of CA in the corrosion field has experienced three stages: simple CA model, complex CA model, and three-dimensional (3D) CA model. We first introduce the mathematical basis of CA, and then summarize the application and research status of CA method in the field of corrosion.

The mathematical definition of CA is:¹³  Σ is the set of cellular states, K = |Σ| is the number of states, N is the scale of CA, i = 0, 1, … , N − 1 is the serial number of each cell, the state of the cell i at the moment t is marked as ⁠, is neighbor of cell i, represents the scale of the neighbor cells, is renewal rule of cellular, is the configuration of CA. Then CA {Σ^N, φ} can be regarded as a mapping on configuration space, that is φ: Σ^N → Σ^N, where and the rules of CA are also often represented by f. CA is mainly composed of elements such as cells, cell space, cell neighbor, and cell rules. In the CA model, each cell has a specific state at a given time, and the evolutions of the cells and the updating of the internal states follow certain rules. A cell is the most basic unit of the CA model. As a state point, a cell is discretely distributed in one-dimensional or multidimensional space. Each cell is in a specific state at every moment, which can be binary {0,1} or discrete set form {S₀, S₁, … , S_n}. The cell space is a collection of distribution points of cells on a spatial grid. In one-dimensional space, there is only one way to divide cells. However, in high-dimensional space, there can be multiple ways to divide cells. In common two-dimensional (2D) space, the main types of cell division are triangular cells, quadrangular cells, and hexagonal cells.

When using a CA model for simulation, the simulated objects are limited, and the selected model space is bounded. This means that the evolution rules at the boundary cannot be defined based on the processing methods used in other regions of the cellular space. Therefore, special processing must be done at the boundary of the cell space. One possible solution is to extend the neighbor of the cell at the boundary. The most common boundary conditions include periodic boundary conditions, fixed boundary conditions, adiabatic boundary conditions, and mapped boundary conditions.¹³  During the evolution of a CA model, certain rules must be followed, and these rules are based on the local range of the cellular space. This means that the state of a cell at the next moment is determined by the state of the cell itself and its neighbors at that moment. Therefore, we must specify which cells belong to the neighbor by defining the neighbor rules of the cell. In general, the neighborhood of a one-dimensional CA model is determined by a radius. For 2D and high-dimensional CA, the neighbors are related to the cellular spatial structure of the model. For example, in the case of a two-dimensional quadrilateral space structure, the most commonly used cellular neighbors include the Von Neumann type, Moore type, and extended Moore type.¹³  Cellular rules are the core aspect of the CA model that must be formulated according to the physical nature of the research object. The evolution rules are applicable to every cell in the cellular space, making them global in nature.

Stafiej, et al.,²²  established a CA model of the metal corrosion-passivation process based on the anode and cathode reactions of space separation and simulated the evolution process of corrosion pit morphology when the protective layer of the metal surface was locally destroyed or the protective layer was not covered on the metal surface. Reis, et al.,²³  defined metal cells in active and passivated states, acidic, neutral, and alkaline solution cells, and then described the formation and dissolution process of the passivated layer on the metal surface. Pidaparti, et al.,^24-26  established two stages of pitting initiation and propagation to simulate multipit corrosion. In the initiation stage, it is assumed that all cells do not corrode, and the chaos degree of cells increases with each time step. When the chaos degree reaches a certain value, the cells become corrosive cells and enter the pitting propagation stage, with corresponding rules for cellular transformation. The changes in the cellular state depend on the surrounding cellular and the pH, dissolution parameters, temperature, and potential difference between the metal and the solution. Taleb, et al.,²⁷  and Vautrin-Ul, et al.,^28-29  established probabilistic CA models based on electrochemical and diffusion mechanisms and simulated the generation of the “island” phenomenon in the corrosion process by replacing the actual electrochemical reaction with probabilistic events in the CA model. This deviation between the electrochemical reaction and Faraday’s law in real solutions is explained. Saunier, et al.,^30-31  simulated the growth process of the corrosion product layer at the solid/solution interface in a metal corrosion system with the diffusion process as the control factor and further explored the diffusion behavior of cations in the corrosion product layer. The thickness of the corrosion product layer, the concentration distribution of the cation, and the intensity and position of the diffusion front were compared under different reaction rates. Fatoba, et al.,³²  used CA to describe the corrosion components and utilized finite element analysis (FEA) with a stress concentration effect to simulate the damage caused by the interaction between corrosion and mechanics, which occurs in the initial environment-induced cracking. This provides a physical framework for modeling damage evolution under the combination of the corrosive environment and mechanical load, which is applicable to stress corrosion cracking (SCC) and corrosion fatigue mechanisms.

To present the corrosion microstructure and initiation process more intuitively, researchers have used the CA method to establish a 3D model to simulate the corrosion process of metal. Van der Weeën, et al.,³³  developed a 3D model that considered the role of corrosive chloride ions in solution. The model focused on mass transport, IR drop, pit initiation, metal dissolution, and cathodic protection, but ignored special corrosion behaviors such as passivation and bimetallic corrosion. Pérez-Brokate, et al.,^34-35  established a 3D CA model based on the electrochemical half reaction of spatial separation in the corrosion process, ion diffusion, acid-base neutralization of the solution, and passivation characteristics of the oxide layer, and studied the influence of pitting morphology and solution acidity on the pitting process. Zenkri, et al.,³⁶  explored the evolution process of metal corrosion damage under three different initial conditions using a 3D random CA model. The model considered the local acidification of the solution, and the research revealed that the autonomous spatial separation of the anode and cathode regions is not only related to the metal matrix and surface passivation film but also to the local acidification of the solution. Compared with 2D simulation, CA simulation in 3D space can describe the morphology of pits and the evolution process of corrosion damage more intuitively, but the amount of calculation considerably increases. Xuefeng, et al.,³⁷  presented a 3D CA model that efficiently simulates multipit corrosion, which considers the breakdown of the passive layer that plays a crucial role in the corrosion process. The results can contribute to simulating the corrosion process and interactions of multipit corrosion covered by a protective coating.

The FEM is a popular numerical technique utilized to solve problems in engineering and mathematical physics that involve behaviors described by differential equations. These differential equations can describe a wide range of physical phenomena, from electrical and mechanical systems to chemical and fluid flow problems. The FEM uses different types of discretization methods to divide the domain of interest into smaller elements and approximate the differential equations with numerical model equations to provide numerical solutions for a given set of boundary and initial conditions. It is important to note that the selection of the appropriate governing equation is a critical step in any FEM-based computational study. In electrochemical and corrosion systems, there are two main modeling approaches using the FEM: those that use the Nernst-Planck equation (N-P) and those that use the Laplace equation.^38-40  In this context, we will first introduce the mathematical basis of these two approaches and then discuss their application to local corrosion simulation.

This assumption simply states that any bulk material (the electrolyte in this case) cannot sustain a net charge. Application of this assumption provides the (n +1)th equation needed and, as a result, the electrolyte potential as well as the concentration of each charged species are obtained.

Using the N-P equation provides complete transient descriptions of the distributions of potential, current density, and concentrations of important species. However, this approach is computationally expensive due to its modeling complexity and execution time. It is important to consider the difference in time scales when modeling both fast steps (such as Faradaic reaction at the interface and ion migration under electrostatic force, etc.) and slow steps (such as species diffusion under a concentration gradient). These varying time scales, along with highly nonlinear electrochemical kinetics as boundary conditions, require highly refined spatial meshing at the relevant boundaries and smaller time steps during the transient study. Furthermore, N-P equation-based modeling generally assumes electroneutrality, and the selection of reference ions may result in discrepancies in the modeling results.

Both approaches use boundary conditions to describe the electrochemical kinetics of reactions, which are generally highly nonlinear. These nonlinear boundary conditions can be described mathematically in some cases, such as via the Butler-Volmer equation for cathodic/anodic reactions under either charge-transfer control, mass-transfer control, or mixed control. However, in many cases, the polarization behavior measured does not follow any prescribed law, such as in a system containing an active-passive transition. In these and other cases, a precise description of the kinetics is required via numerical fits to the polarization behavior for both anodic and cathodic reactions. Regardless of the means by which the electrochemical kinetics are described, they serve as the most critical aspect of the problem statement.

Engelhardt, et al.,⁴⁵  presented a simplified method to estimate pit propagation rates in steel in dilute sodium chloride solutions. They assumed that if the rate of an electrode reaction depends only on the potential, the pit growth rate depends only on the concentration of those species that determine the potential distribution near the metal within the cavity. They compared analytical expressions for calculating the propagation rates of cylindrical and hemispherical pits and explored the influence of aggressive anions on the propagation rate using experimental data. Laycock, et al.,⁴⁶  proposed an accurate model, incorporating a minimum cation concentration for active dissolution, to explore the observation of porous covers made of metal and oxide over the pit mouth and the resultant pit growth morphology for stainless steel. The lace-like porous structure of the pit cover in stainless steel arises from the strong dependence of the dissolution rate on dissolved metal ion concentration. Later, they developed a mathematical model to simulate the propagation of pitting in stainless steel.⁴²  The FEM was used to simulate 2D pit growth, which showed good agreement with experimental results concerning pit morphology and propagation stability. In situ synchrotron x-ray radiography was used to verify the simulation results.⁴⁷ 

Galvele and Gravano^48-49  built a pitting model based on the assumption that metal ions hydrolyze inside the pits and that corrosion products are transported by diffusion. They found that the dominant reason for passivity breakdown at the initial stages of pit growth is localized acidification due to metal ion hydrolysis. Malki, et al.,⁵⁰  developed a corrosion pitting model for the polarization mode and made more detailed calculations to compare the behavior of ferritic and austenitic stainless steel. They obtained the formulation of critical conditions for pit stabilization in terms of pit geometry and applied potential. Amri, et al.,⁵¹  proposed a numerical model of the steady-state behavior of single pitting, which provided an interesting view of the behavior of the chemistry and electrochemistry inside the active pit.

The influence of mass transport can also play an important role in the development of pitting. Alkire, et al.,^52-53  investigated the mass transfer between the interior surface of small rectangular cavities and the solution flowing past the cavity opening using both theoretical and experimental techniques. This approach can be a powerful tool for better understanding phenomena involving cavities, such as polishing, corrosion, rinsing, and roughness. Recently, Srinivasan, et al.,^54-55  examined the flux from a one-dimensional artificial pit electrode corroding under a salt film using experimental and modeling techniques. Their work provided new insight into the origin of the dependence of the measured repassivation potential with pitting depth and contributed toward a quantitative framework relating the various critical influences controlling pitting.

The passive film formed on a metal surface can have a direct impact on the size and life of the pitting. The formation of the passive film depends on the critical potential or current density criterion. Marshall, et al.,⁵⁶  adopted the FEM to simulate and compare the critical pit radius of stainless steel predicted by thermodynamic and kinetic repassivation criteria. They utilized experimental electrochemical boundary conditions to obtain the kinetics of active pits and incorporated geometric and environmental parameters, such as the shape and size of the pit, solution concentration, and thickness of the water layer, to evaluate their impact on the repassivation criteria of pitting. The results showed that parameters, such as lower chloride concentration, smaller cathode diameter, and thinner water layers, lowered the stability of a pit and increased the probability of repassivation, independent of pit shape.

Xu and Pickering⁵⁷  developed a computational model to calculate the distribution of potential and current in the electrolyte phase on the electrode surface for a system. A boundary variation and a trial and error technique were incorporated into the FEM to evaluate the values of the critical distance into a crevice for pure iron in buffered ammoniacal (pH = 9.7) and acetic acid (pH = 4.6) solutions. The results were in good agreement with the available literature values.⁵⁸  The presented model can also be used in the prediction of IR-induced crevice corrosion under open-circuit conditions and in the prediction of precracked samples, critical temperature of SCC, design of anodic protection systems, and design of some corrosion-resistant alloys and coatings.

On the basis of crevice corrosion simulation of single stainless steel and single titanium, Ding, et al.,⁶³  conducted crevice corrosion simulation of 304 stainless steel and titanium overlapping, which is based on a multiphysical model. This work provides a lot of favorable information for the study of crevice corrosion mechanism that is difficult to obtain in experiments, which performs a significant function in the prevention and control of crevice corrosion.

Galvanic corrosion refers to a kind of electrochemical corrosion, sometimes also known as contact corrosion, which occurs when metal components form corrosive batteries in an electrolyte solution and generate galvanic current, which speeds up the dissolution of metals with lower corrosion potential and slows down the dissolution of metals with higher corrosion potential. A new numerical model has been developed to capture the evolution of cross-sectional microstructure by incorporating a moving mesh technique into COMSOL^† MultiPhysics.⁶⁴  The model used an arbitrary Lagrangian Eulerian (ALE) method to determine the moving boundary of the corroding phase (α phase) and velocity of the interface. It is found that a configuration with a continuous network of β phase around the α phase provides the best corrosion resistance. The corrosion rate is found to be higher during the initial stages due to increase in the β phase fraction. However, corrosion is halted after the α phase is preferentially dissolved in the electrolyte solution, exposing the β phase which is electrochemically less active. This trend is obtained from both the numerical model and scanning vibrating electrode technique experiments.

In addition to the study of microgalvanic corrosion, the FEM can also contribute to the overall corrosion distribution of many engineering structures. Deshpande^39,70  presented a Laplace equation-based numerical model to predict the corrosion rate of a galvanic couple, which was capable of tracking the moving boundary of the corroding components. The corrosion rates simulated are in good agreement with those tested from the experimental data. In addition, the FEM model can provide a view on the performance of sacrificial metallic coatings which protect the underlying materials substrate.^71-76 

SCC refers to corrosion damage that can promote crack initiation and propagation through continuous dissolution in the anode area of the crack tip under a corrosive environment and tensile stress. Currently, the simulation of stress corrosion processes using FEM is not widely used due to the environmental impact.^77-79  Some studies instead concentrate on analyzing the pitting or stress distribution of the entire structure,^80-81  while others focus on examining the correlation between various influencing factors such as chloride concentration, temperature, pH, cation species, and SCC.^82-83 

PF is another important method for simulating moveable interface moveable problems. In the PF method, the free energy of the system is assumed to be a function of the field variables and their gradients of each material point. The derivation of the governing equation of the PF variables follows thermodynamics, that is, the dynamic simulated process follows the minimization of the system’s free energy. The position of the phase interface is implicitly solved by the control equation, and it is not necessary to apply boundary conditions on the interface, which also contributes to the application of this method to the morphological evolution problem. The PF model has been used for many years to simulate the evolution of interfaces between phases in various problems such as solidification,⁸⁴  microstructural evolution,⁸⁵  phase transitions in ferroelectric⁸⁶  and ferromagnetic⁸⁷  materials, etc. PF models have been recently used in simulating corrosion through modeling the evolution of the metal/electrolyte interface.⁸⁸  We first introduce the mathematical basis of PF model, and then summarize the application and research status of the PF model in the field of corrosion.

A PF approach for modeling corrosion allows for the implicit approximation of interface evolution by evaluating the distribution of an auxiliary field variable φ over the entire system.⁸⁸  In contrast to sharp interface models (SIMs), which rely on explicit interfaces, a PF model assumes constant values in the bulk of each phase with φ = 0 for the electrolyte and φ = 1 for the metal. The PF variable φ varies continuously across the diffuse interface with a finite thickness l, allowing for a more accurate representation of the pitting boundary. Thus, in a PF model, the auxiliary PF variable φ is utilized to track pitting interface evolution, with a thickness of l representing the diffuse pitting boundary.

where L is the interface kinetics parameter and M is the diffusion mobility for mass transport. The two equations are known as the Allen-Cahn and Cahn-Hilliard equations, respectively. It is worth noting that in practice, only one of the gradient terms (∇c′ or ∇φ) is needed to approximate the energy contribution from the diffuse interface. For simplicity, the concentration gradient energy coefficient α_c is assumed to be 0 in this case. Although Equations (10) and (11) represent the general form of a PF model, the definition of the local free energy density f(c′, φ) is what differentiates various models. By using an appropriate local free energy density f(c′, φ) with the proper assumptions, one can simulate the propagation of pitting corrosion by solving Equations (10) and (11), respectively.

Brewick⁹²  proposed a computational modeling framework to study the effect of crystallographic orientation on corrosion pits by relying on corrosion potentials and elastic constants. However, the model is constrained by assumptions and is currently limited to 2D analysis. In addition, the current model does not contain potential distribution, which also affects the accuracy and fidelity of the model to a certain extent. A 3D PF model was established to study the influence of substrate orientation on the morphology of crystal pits, which shows that the morphology of crystallographic pit changes significantly with the orientation of the matrix, and has the same symmetry as the orientation of the matrix.⁹³  This is the first 3D PF model of crystallographic pits, which is helpful to predict complex shapes of pits, thus promoting the forefront of pitting modeling.

The PF model currently lacks standardized commercial software, leaving researchers to write programs from scratch, requiring a significant amount of computation. As a result, the development of efficient PF calculation methods at both spatial and temporal scales is of utmost importance. Gao, et al.,⁹⁴  proposed a pitting PF model which adopted spatiotemporal adaptive FEM to solve two main challenges. The first challenge is the need for an intensely small time step size. To solve this problem, they combined the Rosenbrock-Euler exponential integrator with the Crank-Nicolson scheme for time discretization. In addition, the adaptive time stepping formula is derived by using the characteristic that the corrosion interface velocity decreases with time. The second one is the requirement for very fine mesh generation. To solve this problem, they propose a simple and effective adaptive mesh generation strategy, which greatly reduces the computational cost. In future work, this method may be extended to other types of corrosion models.

Xiao, et al.,⁹⁵  proposed a quantitative PF to study crevice corrosion. The overpotential through the Butler-Volmer equation was considered as well as the potential distribution and chemical reaction between different kinds of ions in the electrolyte. The system is the same as the system studied by Sharland, et al.^43-44  It should be noted that Sharland’s work was on the basis of steady-state corrosion kinetics, while Xiao’s model is time dependent. They considered six kinds of ions and related chemical reactions. In addition, some physicochemical properties related to corrosion under different metal potentials including overpotential, pH, and corrosion rate, were also examined. The 2D and 3D PF corrosion simulation can be realized by adjusting the modeling area of ion species in the quantitative PF model. This investigation provides a sound framework for further study on 2D and 3D local corrosion.

A PF model was presented to simulate the SCC initiated from surface defects and corrosion pits.⁹⁸  In order to couple the localized corrosion phenomenon with mechanical stress, the film rupture-dissolution-repassivation mechanism was adopted, and the SCC growth rate is considered as a power function of the stress intensity factor by determining the proportionality constant between the velocity of the interface and the interface kinetics parameter. The proposed PF model could validate the availability of the empirical model in reproducing the growth rate of SCC.

Nguyen, et al.,⁹⁹  developed a PF framework for simulating the fracture growth caused by stress corrosion incorporating the impact of material anisotropy. The classical phase transition model for dissolution of material is coupled with the mechanics in an acceptable way, which provides an efficient tool for exploring the competition between electrochemical and mechanical contributing to fracture. The model can also be extended to heterogeneous materials or polycrystalline systems by adopting appropriate conditions.

In order to establish a new PF model for investigating mechanical-electrochemical corrosion, Allen-Cahn type equation for governing the phase transformation from metal to ions in a liquid electrolyte, NP equations for diffusion, Poisson’s equation for electrostatic field distribution, and mechanical equilibrium equation for elastic energy density evaluation were integrated.¹⁰⁰  The numerical results demonstrated that the effect of stress concentration at the tip of a corrosion pit promotes a higher corrosion rate resulting in an accelerated failure of a metallic structure.

Based on the work of Mai and Soghrati,⁹⁸  a new theoretical and numerical framework was presented to simulate mechanical-assisted corrosion in elastoplastic materials, which is capable of capturing the process of pitting and SCC, as well as pit-to-crack transition.¹⁰¹  FEM and implicit integration method were adopted to solve the electrochemical-mechanical coupling equation. The transition from pits to cracks is the natural result of the model and it lays a foundation for mechanical evaluation of engineering materials and structures. Future developments may involve extending the model to other cracking mechanisms such as cathodic reaction-driven cracking. Lin and Ruan¹⁰²  proposed the PF model of mechano-chemical-coupled SCC. To identify quantitatively the critical condition for pit-to-crack transition, the relative-rate parameter difference between SCC and mere corrosion is defined so as to identify quantitatively the critical condition for pit-to-crack transition. When the relative-rate parameter is greater than 1, the pit-crack transition occurs, which characterizes the critical condition where stress-induced degradation occurs faster than electrochemical dissolution. Moreover, the relationship between the stress intensity factor and cracking velocity is exponential, indicating that this is an autocatalytic process caused by stress and corrosion acceleration. The results also show that the SCC was more significant when a metallic component was applied stress. In this case, the effect of the initial pit depth or surface damage was critical in the meanwhile.

The appearance of the PD theory has brought new vigor to the simulation of local corrosion. PD theory was first proposed by Silling¹⁰³  and applied to solid mechanics, which successfully solved the singularity of classical continuum mechanics when solving discontinuity problems. After more than 20 y of development, the PD method has been successfully used to simulate material crack growth^104-107  and heat transfer,^108-111  and to study multiscale¹¹²  and multiphysical field coupling problems.¹¹³  With the indepth study of PD theory by researchers, PD was used to simulate local corrosion for the first time in 2015,¹¹⁴  and relatively satisfactory results were obtained, which laid the foundation for PD’s application in corrosion. In this review, the research progress of PD methods in local corrosion (pitting corrosion, crevices corrosion, IGC, galvanic corrosion, and stress cracking) is reviewed.

The effect of PD simulated repassivation and recrystallization on pitting corrosion of 304 stainless steel in 0.1 mol/L NaCl solution did not fully consider the electrochemical reaction on the electrode surface covered by the salt film, resulting in smaller pitting and lace-cover in 3D than the experimental results.^117-118  Jafarzadeh, et al.,¹¹⁹  improved the model by using the classical pitting theory,¹²⁰  and simulated the pitting growth under coverage of the salt film. Compared with the original model,¹¹⁴  the results of this model are closer to the real condition.

Jafarzadeh, et al.,¹²¹  introduced a new crevice corrosion damage PD model and verified it with experiments in the literature. The electrochemical mechanism is considered as a simple pattern that ion concentration is dependent on corrosion rate. This simplification makes it feasible to define local dissolved microfluxes at the solid/liquid interface. By modifying the PD formula to accommodate the event horizon of any shape, calculations can be made efficiently for areas with extreme aspect ratios, such as those observed in cracks. The model was verified by the experimental data of bolt washers. The simulation predicted the site and depth of a deep corrosion trench formed at certain distances from the crevice opening. The experimental crevice corrosion process is reproduced by far-field boundary conditions. The results also show that the two factors considered (diffusion-driven transport of dissolved metal ions and dissolution dependent on their concentration) are effective in simulating the evolution of crevice corrosion.

IGC is a kind of local corrosion occurring along the grain boundaries of materials. As the anodic reaction rate of grain boundaries is larger than that of grains, grain boundaries are always corroded preferentially than grains. According to the heterogeneous material corrosion model,¹¹⁴  grains can be regarded as the parent phase and grain boundaries as the second phase when studying IGC.¹²²  The model can also simulate the process from single IGC to the whole grain dissolution, that is, the corrosion of polycrystalline materials. The microparticles or precipitated phases of IMC on grain boundaries may be the second phase. Both grains and grain boundaries dissolve, and the crystal orientation can affect the corrosion of polycrystalline materials. It is worth noting that there are differences between IGC and galvanic corrosion here. IGC does not consider galvanic coupling and mainly focuses on the influence of grain size, orientation, and other factors on localized corrosion.

Research has shown that FEM-based or PF models of corrosion in galvanic couples often require an artificial “step-down” in the geometry at the material interface of the couple in order to produce reasonable results.¹²³  This step is necessary for models based on partial differential equations (PDEs) due to difficulties in assigning appropriate boundary conditions at the galvanic couple interface and properly initializing the motion of the corrosion front. However, the introduction of such artificial geometrical modifications reduces the generality of the model, and when mechanical loadings are applied to a galvanic couple system, the stress profile near the interface may differ from the actual one, particularly in the early stages of the corrosion process. Additionally, if there is a strong stress dependency on the corrosion rate, the probability of obtaining incorrect results for coupled problems may be further increased. In contrast, the PD model presented here does not require any artificial changes in geometry at the interface to accurately predict corrosion progression as observed in experiments.

The CA method has been successfully applied in the simulation of many electrochemical problems. It is a powerful tool to study the complex physical and chemical processes such as pitting corrosion, formation and dissolution of passivated film, material diffusion, and formation of corrosion products. It is capable of simulating the mesoscopic corrosion damage evolution process. It should be pointed out that the simulation of the corrosion process by using CA also has its limitations, which are mainly manifested as the limitations of cellular state and cellular transformation rules. Although some basic rules of electrochemical reactions are introduced into the CA model, they are only expressed qualitatively, and the model cannot be calibrated by experimental observation. The results obtained by the model, especially the dynamic evolution process, have no specific time scale. Therefore, the CA model can be called a semideterministic model.

The principle of FEM is to solve the governing equation of a corrosion electric field. Its calculation is relatively high in spatial scale and time scale but accompanied by relatively weak characterization capability of microscopic characteristics. FEM is suitable for solving nonlinear and nonuniform corrosive medium problems, so it has certain advantages for the simulation of corrosion development and the dynamics process of evolution. However, the cost is expensive and prohibitive in 3D, and it is difficult to deal with infinite and semi-infinite domain problems.

PF is a flexible framework for the prediction and simulation of corrosion due to the autonomous movable interface. In the PF model, the phase transition and mass transport process are characterized by two distinguished but coupled PDEs, which significantly increase the calculation burden and algorithm complexity.¹²⁹  Although the PF model introduces a length scale into the model through the thickness of the diffusion layer (the PF transition region between phases is the input of the problem), they are still local models and impose classical boundary conditions on the set of PDEs. The diffusion interface in the PF model leads to the transition area of the corrosion front. Moreover, the damage evolution in the PF model is strongly dependent on the specific choice of the energy functional.^130-132 

As a nonlocal theory, PD has a mature application in the aspects of material fracture and heat transfer and can simulate the corrosion damage evolution of metal materials. The model can successfully simulate the local corrosion process through surface diffusion, ion transport, electrochemical reaction at the metal/solution interface, and free movement of the solid/liquid interface. The PD corrosion damage model can simulate the damage of the secondary surface layer and the formation of lacy cover in pitting corrosion, the growth of pitting corrosion under the salt film, simplifying the autocatalytic effect in crevice corrosion, IGC, potential distribution at the interface of corroded metal/electrolyte solution, and crack growth at pitting under stress. PD corrosion damage model can be calibrated by other numerical simulation results and experimental results, but it is still in the developing stage.

Table 1 compares the models reviewed in this study in terms of governing equation, numerical method, main applications, advantages, and disadvantages.

Table 1.

Comparison of Cmputational Models for Localized Corrosion¹¹⁵ 

View large

View Large

The surface potential and current density distribution of the mobile interface can be influenced by metal corrosion particles, which can change the corrosion rate. The dynamic tracking of the corrosion interface is therefore crucial to build an accurate simulation model and is currently a major research challenge.¹³³  In recent years, a dynamic interface model based on FEM has been proposed and used to simulate the evolution process of the active metal dissolution interface^134-135  but the SIM, which needs to track the metal interface at every step, requires a lot of remapping in the whole simulation process. To ensure that each step of geometric evolution can meet the computational requirements, several methods have been introduced to arbitrary modeling, including the level set method (LSM)^136-137  and the ALE¹³⁸ . However, maintaining a high-quality grid throughout the simulation process remains a challenge, especially for complex corrosion modeling problems

With the continuous development of computer technology, some complex dynamic interface models have also been successfully implemented in the research of tracking corrosion interfaces, such as (i) FEM, some researchers have successfully used the finite volume method to simulate metal pitting.¹³⁹  This method is especially suitable for dealing with mathematical problems with spatial discontinuities. The location of the pitting interface is determined according to the ion concentration in each volume element, thus eliminating the difficulty of creating a new grid.¹⁴⁰  (ii) Meshfree method, which avoids the process of grid generation in numerical calculation and uses some arbitrarily distributed coordinate points to construct interpolation function discrete governing equation to track the expansion of complex interface, is a feasible method to simulate the problem of corrosion moving interface.¹⁴¹  (iii) Extended/generalized FEM (E/GFEM), which is a major improvement of the traditional FEM, uses discontinuous form function to characterize the discontinuity in the domain of domain, and is also a mesh-independent FEM.^142-143  There have been relevant studies on the simulation of pitting corrosion by using it in combination with LSM.¹⁴⁴ 

Although there are more applications of the method for tracking moving interfaces, most of the methods are still in the theoretical stage. In terms of engineering applications and on a larger scale, popular models like PF and PD are not supported by definite commercial software and most researchers still need to write programs line by line. Implementation of these methods requires a certain reserve of professional knowledge, relevant mathematical basis, and strong computer programming ability, making it difficult for quite a few researchers and restricting their application to a large extent. Furthermore, the computational cost of implementing both models is prohibitively expensive, especially in 3D. However, recent advancements have considerably improved this situation.^145-147  A fast convolution-based method has been proposed for efficient discretization of PD/nonlocal models, reducing memory allocation by several orders of magnitude. The result shows that PD problems that would have required years of computations with existing discretization methods can now be solved in a matter of days with the proposed method, making fast computation of fracture and damage with high accuracy possible.

In the first place, computational modeling of localized corrosion has advanced to the point that major uncertainties in predictions are due to the lack of experimental measures of input parameters. Therefore, it is highly desirable to adapt emerging multiphysics programs for modeling corrosion processes, allowing corrosion scientists and engineers to access guidelines for material design more easily. Several mechanisms contribute to structural corrosion, including processes at the megascale, material deformation, electrochemical reactions, mass transport, and fluid dynamics.^148-149 

Due to the interdependency between these mechanisms, it is highly desirable to use a multiphysics approach to model and simulate the corrosion process within a feasible timeframe. The problem is multidisciplinary, with a large number of variables and controlling factors involved such as material composition, strengthening mechanisms, microstructure, loading modes and frequency, machining and residual stresses, environmental composition, cyclic effects of environmental variables, etc. Expertise in many scientific fields is required for the analysis and control of corrosion phenomena in structural alloys.

Moreover, corrosion of materials is often not solely caused by a single physical driver. It is often accompanied by various influences such as stress, flow, thermal, and magnetic fields. The existence of multiple physical influences makes it difficult to model accurately.¹⁵¹  Taking aluminum alloy materials used in seaplanes, ships, and offshore photovoltaic supports as examples, their service environment is mainly the seawater/atmosphere interface, which is a typical environment with multifield and multiphase couplings.^152-154  The corrosion of these structural materials under the influence of temperature, pH, dissolved oxygen, sea/atmosphere alternations, and waves is a key factor affecting the long-term stable use of marine engineering and weapons equipment. However, marine tests are complicated, time-consuming, and costly. Therefore, it is particularly important to use computational models to simulate the local corrosion damage of metal materials under multifield and multiphase couplings.

• In this paper, the characteristics and the application of four numerical simulation methods in localized corrosion, including CA, FEM, PF model, and PD model, are reviewed. The review shows that the CA is a powerful tool for studying pitting corrosion and passivation film formation and complex microprocesses, such as dissolution and substance diffusion, widely used in the microscopic mechanism research of metal corrosion. FEM is appropriate for solving nonlinear and nonuniform corrosive medium problems, so it has certain advantages for the simulation of corrosion development and the dynamics process of evolution. As another important method to simulate the problem of the variable interface, the PF method does not require the application of boundary conditions at the interface, which promotes the application in the problem of corrosion morphology evolution. At present, PD, as a nonlocal theory, can successfully simulate the corrosion damage evolution of metal, which can be calibrated by other numerical simulation results and experimental results. However, the corrosion of metal materials is often affected by multiple influencing factors and multiple physics fields, so the research on computational models on the evolution of local corrosion damage is still in the development stage.

• In summary, numerical simulation has become one of the important means to study metal corrosion and protection, and with the rapid development of software science and computer technology, the operation speed and accuracy of numerical simulation will be greatly improved, numerical simulation in the field of corrosion application will be more prominent advantages. According to the different corrosion systems and the problem to be solved of selecting the appropriate simulation method, one can make full use of its strengths and circumvent its weaknesses.

This research effort was financially supported by Southern Marine Science and Engineering Guangdong Laboratory (Zhuhai) [Fund, D5110200598]. Generous thanks are extended to the Fund Monitors for their encouragement and support of this exhaustive study.