Abstract
The stochastic nature of pitting corrosion of metallic structures has been widely recognized. It is assumed that this kind of deterioration retains no memory of the past, so only the current state of the damage influences its future development. This characteristic allows pitting corrosion to be categorized as a Markov process. In this paper, two different models of pitting corrosion, developed using Markov chains, are presented. Firstly, a continuous-time, nonhomogeneous linear growth (pure birth) Markov process is used to model external pitting corrosion in underground pipelines. A closed-form solution of the system of Kolmogorov's forward equations is used to describe the transition probability function in a discrete pit depth space. The transition probability function is identified by correlating the stochastic pit depth mean with the empirical deterministic mean. In the second model, the distribution of maximum pit depths in a pitting experiment is successfully modeled after the combination of two stochastic processes: pit initiation and pit growth. Pit generation is modeled as a nonhomogeneous Poisson process, in which induction time is simulated as the realization of a Weibull process. Pit growth is simulated using a nonhomogeneous Markov process. An analytical solution of Kolmogorov's system of equations is also found for the transition probabilities from the first Markov state. Extreme value statistics is employed to find the distribution of maximum pit depths.
1. Introduction
Localized corrosion, specifically pitting corrosion of metals and alloys, constitutes one of the main failure mechanisms of corroding structures such as pressurized containers and pipes. Pits cause failure through perforation of the component wall. In other cases, pits become nucleation centers for cracks [1]. In the oil and gas industry, pitting corrosion is a major problem, especially for transporting pipelines [2].
Pitting corrosion comprises two main processes: pit initiation and stable pit growth (pit repassivation is not considered in this paper). It is accepted that pit initiation can be a consequence of the breakdown of the passive layer caused by random fluctuations in local conditions. This process takes some time, usually called induction (nucleation or initiation) period [3]. Passive layer breakdown, followed by localized metal dissolution, is the most common mechanism of pitting corrosion. However, pitting can also occur as the result of the active dissolution of certain regions of the material at its surface, such as nonmetallic inclusions, which are susceptible and dissolve faster than the rest of the surface [4]; in this case, the pitting induction time is typically shorter.
After a pit has nucleated, it can repassivate (stop growing) immediately or grow and then repassivate. This process is regarded as metastable pitting. If a pit is able to grow indefinitely, it becomes a stable growing pit [5]. Those pits that reach the stable growth regime become part of a pit population that shows a remarkable stochastic behavior [6, 7].
The time evolution of pit depth due to corrosion is often expressed as a power function of time [7–11]: , where is the exposure time, and are parameters related to the corrosion process, and stands for pit initiation time. This same function has been found to govern the pitting corrosion growth in stainless and mild steels and aluminum alloys [11, 12].
Pitting corrosion occurs in a wide range of metals and environments. This fact points out to the universality of this phenomenon and suggests that randomness is an inherent and unavoidable characteristic of this damage over time, so that stochastic models are better suited to describe pitting corrosion than deterministic ones. Localized corrosion cannot be explained without assuming stochastic points of view due to the large scatter in the measurable parameters such as corrosion rate, maximum pit depth, and time to perforation [13]. Many variables of the metal-environment system such as alloy composition and microstructure, and composition of the surrounding media and temperature, are all involved in the pitting process [4]. Such complexity imposes the development of theoretical models and simulation tools for a better understanding of the outcome of the pitting corrosion process. These tools help predict more accurately the time evolution of pit depth in corroding structures as the key factor in structural reliability assessment.
Another important characteristic of the pitting corrosion process that is worth noting is the time and pit-depth dependence of the corrosion rate [6, 7]. It has been established that, for a given pit, the growth rate decreases with time, while for pits with equal lifetimes, the corrosion rate is larger for deeper ones.
Provan and Rodriguez [14] are amongst the first authors to use a nonhomogenous Markov process to model pit depth growth. In their model, the authors divided the space of possible pit depths into discrete, non overlapping states and numerically solved the system of Kolmogorov’s forward equations (1) for the transition probabilities between damage states and . However, Provan and Rodriguez modeled pitting without taking into account the pit generation process and proposed an expression for the intensities of the process that conveyed no physical meaning. They did not discuss the method used to solve the system of equations either. This has made it impossible to reproduce their results (deeper discussion on this topic can be found in [15])
Other authors who made use of Markov chains to model corrosion were Morrison and Worthingham [16]. Making use of a continuous time birth process with linear intensity , they determined the reliability of high-pressure corroding pipelines. For their purpose, these authors divided the space of the load-resistance ratio into discrete states and numerically solved the Kolmogorov’s equations in order to find the intensities of transition between damage states. Afterwards, probability distribution function of the load-resistance ratio was estimated and compared to the distribution obtained from field measurements. Worthingham’s model was further improved by Hong [17], who used an analytical solution to the system of Kolmogorov’s equations for the same homogeneous continuous type of Markov process and derived the process probability transition matrix in order to evaluate the probability of failure. Hong also investigated the influence of pit depth on the load-resistance ratio. The merits and limitations of these two models are discussed in detail in [15, 18].
In recent years, modeling of pitting corrosion with Markov chains has shown new advances. For example, Bolzoni et al. [19] have modeled the first stages of localized corrosion using a continuous-time, three-state Markov process. The Markov states of the metal surface are passivity, meta stability, and stable pit growth. On the other hand, Timashev and coworkers [20] formulated a model based on the use of a continuous-time, discrete-state pure birth homogenous Markov process for stochastically describing the growth of corrosion-caused metal loss. The goal was to assess the conditional probability of pipeline failure and to optimize the maintenance of operating pipelines. In their model, the intensities of the process were calculated by iteratively solving the proposed system of Kolmogorov’s forward equations. The drawbacks of Timashev’s model are discussed in detail in [18].
In the present paper, a review of the Markov models developed by the authors in an attempt to describe pitting corrosion is presented. The first model intends to solve a problem that is crucial in reliability-based pipeline integrity management: the accurate estimation of future pit depth and growth rate distributions from a (single) measured or assumed pit depth distribution. It has been recognized [21] that such estimation can be carried out only if oversimplifications are made, or if additional information, besides the pit depth distribution, is available. In the developed model, it is postulated that in the case of external pitting corrosion in underground pipelines, this additional information can be attained from the available predictive models for pit growth as a function of the soil characteristics [8, 22]. A model for pit growth previously developed by the authors has been used to perform Monte Carlo simulations in order to predict the distribution of maximum pit depths as a function of the pipeline age and the physicochemical characteristics of the soil [9].
A nonhomogenous linear growth pure birth Markov process, with discrete states in continuous time, is used to model external pitting corrosion in underground pipelines. The system of forward Kolmogorov’s equations (1) is analytically solved using the binomial closed-form solution for the transition probabilities between Markov states in a given time interval. This Markov framework is used to predict the time evolution of the pitting depth and rate distributions. The analytical solution becomes available under the assumption that Markov-derived stochastic mean of the pit depth distribution is equal to the deterministic mean of the distribution obtained through Monte Carlo simulations. This assumption is made for different exposure times and different soil classes defined according to soil physicochemical characteristics that are easy to measure in the field. In this way, the transition probabilities are obtained as a function of soil properties for a given time point, and the corrosion rate and future pit depth distributions are predicted. Real-life case studies are presented to illustrate the proposed Markov model. The main advantage of this model resides in its capability of correctly predicting the time evolution of pit depth and rate distributions over time.
One of the main goals of reliability analysis is to estimate the risk of perforation of in-service components produced by the deepest pit. Extreme value statistics is commonly used together with pit depth growth modeling to predict the risk of failure of in-service components and structures [23, 24].
It is recognized [6, 7] that the deepest pits grow at higher rates than the rest of the defects right from the beginning of the corrosion process, so that the maximum pit depths are commonly sampled from an exponential parent distribution that constitutes the right tail of the pit depth distribution of the whole pit population [7]. Besides that, in the corrosion literature [13, 14, 24, 25], it is a well-established fact that the Gumbel extreme value distribution [26, 27] is a good fit to the maximum pit depth distribution obtained by measuring the maximum pit depth on several areas. The cumulative function of the Gumbel distribution, with location parameter and scale parameter , is
It is important to underline that the previously mentioned Markov model has proved successful in modeling the time evolution of the entire pit depth population but failed to correctly describe the evolution of pit-depth extremes. The second model presented in this work focuses on the simulation of the time evolution of extreme pit depths. The model is based on the stochastic description of pitting corrosion, taking into account pit initiation and growth. A nonhomogeneous Poisson process is used to model pit initiation. The distribution of pit nucleation times is simulated using a Weibull process. Pit depth growth is also modeled as a nonhomogeneous Markov process. The system of forward Kolmogorov’s equations (1) is solved analytically for the transition probability from the first state to any th state during a given time interval [28]. From this solution, the cumulative distribution function of pit depth for the one-pit case can be found. This distribution function has an exponential character and corresponds to the parent distribution from which the extremes can be readily drawn. The extreme depths distribution for the -pit case is obtained by the multiplication of the parent cumulative distribution functions of the pits population. This stochastic model is able to predict the time evolution of the probability distribution of maximum pit depths.
2. Stochastic Models of Pitting Corrosion
In this section, two different Markov chain models are presented to describe pitting corrosion. The first one is focused on the description of time evolution of pit depth and rate distributions in operating underground pipelines. The second model deals with the description of maximum pit depths when multiple corrosion pits are taken into account in a laboratory (controlled) pitting corrosion experiment.
In the following, only the definitions relevant to the focus of the presented models are given. The reader is referred to [28, 29] for the formalism theory of Markov processes.
In both models, the possible pit depths constitute the Markov space. The material thickness (along which pits propagate) is divided into discrete Markov states. The corrosion damage (pit) depth, at time , is represented by a discrete random variable such that , with . Furthermore, it is postulated that the probability that the damage that is at the th state at the moment increases by one state in a very small interval of time is expressed as . For a continuous-time, nonhomogenous linear growth Markov process with intensities , the probability that the process in state will be in state at some later time obeys the system of Kolmogorov’s forward equations presented in (1). In this infinitesimal transition scheme, can be interpreted as the jump frequency between the th to the th states during the time interval . This means that the number of states transited by the corrosion pit in a short time interval can be written as
It is not difficult to find out that the functions and have direct physical meaning when Markov processes are used to model pitting corrosion. They are related to the pit growth rate and pit depth, respectively.
2.1. Markov Chain Modeling of Pitting Corrosion Depth and Rate in Underground Pipelines
In [28, page 304], it is shown that for a Markov process defined by the system of (1), the conditional probability of transition from the th state to the th state () in the interval () can be obtained analytically and has the form where is defined by (3).
Equation (4) shows that the increase in pit depth over follows a negative binomial distribution with parameters and . From the transition probability it is possible to estimate the probability distribution function of the pitting corrosion rate over the time interval , when the pit depth is at the th state as
From the distribution function , it is straightforward to derive the pitting rate probability distribution associated with the entire pit population (all possible depths) as
Until this point, we have the transition probabilities from any state to the state in the interval , given by (4). The corrosion rate distribution can also be derived through (6). In principle, if the probability distribution of the corrosion depth at is known, that is, , the pit depth distribution at any future moment in time can be estimated using
For the case of buried pipelines, the initial probability distribution of pit depths can be obtained if the corrosion damage in the pipeline is monitored using in-line inspection (ILI), being the time of the inspection. The values of the probabilities can be estimated from the ratio of the number of corrosion pits with depths in the th state to the total number of observed pits.
It is evident from (4) that the transition probabilities depend on the functions and . At this point, physically sound expression for and should be proposed in order to estimate the evolution of the pit depth and corrosion rate distributions. For that, the knowledge about the pitting corrosion process in soils must be used. It is postulated that the stochastic mean pit depth can be assumed to be equal to the deterministic mean of the pitting process, which can be estimated from the observed evolution of the average pit depth over time as
This equality is true under certain simple assumptions that are explained in [29]. In summary, a sufficient condition for the equality of the stochastic and deterministic means is that for any positive integer , the structure of the process starting from individuals is identical to that of the sum of separated systems each starting from .
Cox and Miller [29] have shown that if the initial damage state is at , so that , then the time-dependent stochastic mean of the linear growth Markov process can be expressed as
If we consider that a power function is an accurate deterministic representation of the pit growth process, one can postulate that the deterministic mean pit depth at time is [8] where and are the pitting proportionality and exponent parameters, respectively, and is the starting time of the pitting corrosion process. In systems where passivity breakdown and/or inclusion dissolution are the prevalent mechanisms for pit initiation, would represent the initiation time of stable pit growth. In the case of underground pipelines, this parameter corresponds to the total elapsed time from pipeline commissioning to coating damage plus the time period when the cathodic protection is still effective preventing or attenuating external pitting corrosion after coating damage.
The increase of the deterministic mean pit depth in the time interval is where can be interpreted as the deterministic intensity (rate) of the process. The pitting rate, obtained by taking the time derivative of (10), obeys the functional form of (11) with and as
If represents the stay time of the corrosion pit in the first state of the chain, then during the time interval between and . If is significantly less than the simulation time span, it is easy to show from (8)–(10) that the value of the function can be approximated as
From this, taking into account (3), it follows that
Note that the intensity of the Markov process is inversely proportional to the exposure time , as it is the case of the deterministic intensity according to (11) and (12).
One can substitute the expression for from (13) into (4) and to show that the probability parameter can be expressed as
Suppose that a pit is in the state at . Let be the average damage rate in the time interval . It can be shown [19] that the instantaneous pitting rate predicted by the stochastic model is
The stochastically predicted instantaneous damage rate agrees with the deterministic rate given by (12). This coincidence is a sign of the adequacy of the proposed Markov pitting-corrosion model.
For the case of underground pipelines, the pitting corrosion damage evolution can be undertaken as follows. The measured or assumed pit depth probability distribution at is used as the initial corrosion damage distribution , where is the time of the pit depth measurement. The value of the probabilities is estimated from the ratio of the number of corrosion pits with depths in the th state to the total number of observed pits. The transition probability function can be identified from (4) and (15) if the parameters and are known.
The estimation of the pitting parameters is possible thanks to the previously developed [8, 9] predictive model for localized corrosion in buried pipelines, which relates the physicochemical conditions of the pipe and soil to parameters and . The model is based on a multivariate nonlinear regression analysis using (10), with the field-measured maximum pit depths as the dependent variable to be predicted and the pipe age and soil-pipe characteristics as the independent ones. The experimental details and corrosion data used to produce the model can be found elsewhere [8, 9]. From this model, the pitting initiation time was estimated to have the values displayed in Table 1 for different soil classes. In Table 1, class “All” corresponds to a general class that includes all the soils found in the field survey carried out in South Mexico [8]. The pitting parameters and in (10) were estimated as linear combinations of soil and pipe characteristics [8] for the different soil classes.
Later on, Caleyo et al. [9] used a Monte Carlo simulation based on (10) and on the obtained linear combinations of environmental variables for the pitting parameters to predict the time evolution of the maximum pit depth probability distribution in underground pipelines. The measured distributions of the considered soil-pipe independent variables were used as inputs in Monte Carlo simulations of the pitting process for each soil class to produce the simulated extreme pit depth distribution for each soil category and exposure time considered. The reader is referred to [9] for details about the used pipe-soil properties distributions and the Monte Carlo simulation framework. The Monte Carlo produced extreme pit depth values for each soil class at different time points were fitted to the generalized extreme value (GEV) distribution [30]: where , , and are the shape, scale, and location parameters of the distribution, respectively.
Figure 1(a) shows the time evolution of the GEV distribution fitted to the Monte Carlo predicted extreme pit depths for a general soil class found in southern Mexico [8]. The inset of Figure 1(a) shows the evolution of the mean and variance of such distributions with time. These average values are used to estimate the pitting proportionality and exponent factors ( and ) associated with the typical (average) values of the predictor variables in this soil category. The time evolution of the mean of the simulated maximum pit depth distributions was fitted to (10) using the corresponding value of given in Table 1. Typical values and , derived following this method, are displayed in Table 2 for the analyzed soil textural classes.
(a)
(b)
From the estimated values of and for the general soil class and the value of shown in the last column of Table 1, it is possible to obtain the function (13) and (15) for typical conditions of the “All” soil category [8, 9].
Figure 1(b) shows the time evolution of and predicted from the parameters associated with the pitting process in soils. Given that is completely determined by the extent of the damage [15], its value is unique for each soil class. Also, the probability is unique in each soil type. However, in Figure 1(b), two different curves of are displayed to illustrate the dependence of on time. They correspond to two distinct values of (see (15)), the time when the initial pit depth distribution is measured. The value of at a given moment in time increases with increasing and decreases when the length of the interval () increases. The increase in with indicates that the mean and variance of the pitting rate decrease as the lifetime of the pitting damage increases, which is in agreement with experimental observations. Also, the form of in Figure 1(b) suggests that, for pits with equal lifetimes, the deeper the pit, the smaller the value of and, therefore, the larger the mean and variance of the pitting rate. This characteristic of the developed model will permit accurately predicting the corrosion rate distribution in a pitting corrosion experiment.
To estimate the evolution of pit depth and pitting rate distributions using the proposed Markov model one can use (7) and (6), respectively. As has been already stated, from the results of an ILI of the pipe and the knowledge of the local soil characteristics, it would be possible to estimate the pitting corrosion damage evolution.
To illustrate the application of the proposed Markov chain model, it is employed in the analysis of an 82 km long operating pipeline, used to transport sweet gas since its commissioning in 1981. This pipeline is coal-tar coated and has a wall thickness of 9.52 mm (0.375′′). It was inspected in 2002 and 2007 using magnetic flux leakage (MFL) ILI. The pit depth distributions present in the pipeline were obtained by calibrating the ILI tools using a methodology described by Caleyo et al. [31]. These distributions are shown in Figure 2 in the form of histograms. The first, grey-shadowed histogram represents the depth distribution of pits located and measured by ILI in 2002, while the hatched histogram represents the depth distribution of pits measured by ILI in mid-2007.
In order to apply the proposed model, the pipe wall thickness was divided in 0.1 mm thick non overlapping Markov states, so that the pitting damage is represented through Markov chains with states ranging from to . Unless otherwise specified, this scheme of discretization of the pipe wall thickness is used hereafter to represent the pitting damage penetration.
The empirical depth distribution observed in 2002 was used as the initial distribution so that years, years, and being the number of pits with measured depth in the th state. The soil characteristics along the pipeline were taken as those of the “All” (generic) soil class and the value of was taken equal to 2.9 years, as indicated in Table 1. From Table 2, a value of 0.780 was assigned to in order to obtain the predicted pit depth distribution ( years) for 2007. The values of can be calculated using the following expression, which is deducted from (4), (7), and (15):
In Figure 2, the Markov chain-predicted pit depth distribution for 2007 is represented with a thick line.
In order to test the degree of coincidence between the empirical distribution (for 2007) and the Markov-model predicted distribution, the two-sample Kolmogorov Smirnov (K-S) and the two-sample Anderson-Darling (A-D) tests were used. This was done by means of a Monte Carlo simulation in which 1000 samples of 50 depth values were generated for each one of the distributions under comparison and then used in pairs for the tests. In the case of the Markov-predicted distribution the data samples were generated as (pseudo) random variates with probability density function (pdf) as displayed in Figure 2 and given by (18). The samples from the empirical data were produced by sampling the experimental 2007 dataset with replacement. The K-S test gave an average value of 0.42, with only 7% of cases where the test failed at 5% significance. By its part, the A-D test produced an average value of 0.39 with an 11% fraction of failed tests. This is an evidence of the validity of the null hypothesis about both samples, modeled and empirical, coming from the same distribution. The good agreement between the empirical pit depth distribution observed in 2007 and the Markov chain-modeled distribution also points to the accuracy of the proposed model.
The corrosion rate distribution associated with the entire population of pits in the 2002–2007 period was also estimated using (19), which was derived from (4), (5), (6), and (15). Consider
The estimated rate distribution is shown in the inset of Figure 2; it can be fitted to a GEV distribution (17) with negative shape parameter , although it is very close to zero. This means that the Weibull and Gumbel subfamilies [27] of the GEV distribution are appropriate to describe the pitting rate in the pipeline over the selected estimation period.
It has to be noticed that the proposed Markov model can be used to predict the progression of other pitting processes. For this to be done, the values of the pitting exponent and starting time must be known for the process. These can be obtained, for example, from the analysis of repeated in-line inspections of the pipeline, from the study of corrosion coupons, or from the analysis of laboratory tests. To show this capability of the model, another example is given, using this time a different experimental data set gathered from two successive MFL-ILIs. The analysis involves a 28 km long pipeline used to transport natural gas since its commissioning in 1985. This pipeline is coal-tar coated, made of API 5L grade X52 steel [32], with a diameter of 457.2 mm (18′′) and a wall thickness of 8.74 mm (0.344′′). The first inspection was conducted in 1996 and the other in 2006. The ILIs were calibrated following procedures described elsewhere [31]. Before using the ILI data, the corrosion defects from both inspections were matched using the following criteria: (i) the matched defects should agree in location, and (ii) the depth of a defect in the second inspection must be equal to or greater than its depth in the first inspection. Only the matched defects were used in the analyses that followed to ensure that only the actual defect depth progression with time was considered (one-by-one if necessary) without including defects that might have nucleated in the interval between the two inspections. The depth distribution of the matched defects of the 1996 inspection was taken as the initial depth distribution fed into the proposed Markov model. The resulting depth distribution was subsequently compared to that of the defects observed in the 2006 inspection.
From the knowledge of the mean depth values from both inspections, an empirical value of the pitting exponent was estimated instead of using the value recommended in Table 2. It was calculated for the specific type of soil in which the pipe under analysis was buried by using the ratio between the empirical pit depth means and , from the ILIs in 1996 and 2006, respectively. Assuming that the pit depth mean complies with (10) and that parameters , and (which depend exclusively on the soil and pipe properties) are the same at the times of both inspections, it is possible to obtain the pitting exponential parameter from where years and years are the times of the inspections, and is the incubation time of the corrosion metal losses in the pipeline. The value of was taken equal to 2.9 years, the average pitting initiation time for a generic soil class (Table 1).
As in the previous example, the wall thickness was divided into (0.1 mm thick) equally spaced Markov states. The observed defect depth was converted to Markov-state units, and the depth distribution was given in terms of the probability for the depth in a state equal or less than at time . The proposed Markov model (18) was applied to estimate the pit depth distribution after a time interval of 10 years, with years, years, and years.
In Figure 3, the final (2006) pit depth distribution of the by-ILI measured and matched 179 defects in this pipeline is shown on a probability density scale in the form of a histogram, together with the pit depth distribution predicted by the Markov model for 2006. Again, the Monte Carlo framework described in the previous example was followed in order to apply the two-sample K-S test. The sample size was 50, and the sampling process was repeated 1000 times. The resulting average value was 0.07, which is enough for not rejecting the null hypothesis that both the empirical and modeled distributions come from the same distribution. From the results of the test and from what it is shown in Figure 3, it can be concluded that the Markov model predicts a pit depth distribution that is close to the experimental one. More important, the model is also capable to correctly describe the form (shape) of the defect depth distribution, as seen also from this figure.
To further explore the capabilities of the model in estimating the time evolution of pit depths, the defects in the initial (1996) inspection were grouped into 6-state (0.6 mm) depth intervals in order to test the model efficacy in predicting the pit depth evolution by depth interval. An example of these grouping and modeling approach is shown in the inset of Figure 3. The chosen initial depth interval is that between 9 and 15 states (0.9 and 1.5 mm). In the inset, a histogram of the matched defect depths in the second ILI together with the distribution (solid line) predicted by the Markov model starting from the chosen initial depth subpopulation is depicted. Here as well, the model is capable of reproducing the experimental results. The same was proven true for the rest of depth intervals in the initial (1996) pit depth population. Therefore, it can be stated that the Markov model is not only accurate in estimating the time evolution of the entire defect depth distribution, but also in estimating the time evolution of defect subpopulations that differ in depths.
Additionally, to check the ability of the proposed Markov model to predict the corrosion rate (CR) distribution, an empirical CR distribution was determined using the data from both inspections based on the observed change in depth of the matched defects over the interval . The empirical rate distribution was then compared with the corrosion rate distribution predicted by the model (19). This comparison can be observed in Figure 4. The two-sample K-S test for 1000 pairs of 50-point CR samples yielded an average value of 0.27 with only 12% of rejections at 5% significance level. Again, the Markov model satisfactorily reproduces the empirical CR distribution.
The reasons behind the foregoing results lie in the ability of the Markov chain model to capture the influence of both the depth and age of the corrosion defects on the deterministic pit depth growth together with its ability to reproduce the stochastic nature of the process, also as a function of pit depth and age. The inset of Figure 4 is a proof of this last statement. It shows the experimental corrosion rate mean and variance dependence on pit depth. To obtain these results, pit depths were grouped into depth intervals of ten states (1 mm) and the corrosion rate values, derived from the comparison of the two ILIs, were averaged within each group. From this inset, it can be concluded that both the corrosion rate mean and variance increase with defect depth, as was previously noted from different experimental evidence [7, 18]. This result is critical to support the criteria that a good corrosion rate model should consider both the age and depth of the to-be-evolved corrosion defects together with the actual shape of the function describing the observed dependence of the corrosion defect depth with time.
2.2. Markov Chain Modeling of Extreme Corrosion Pit Depths
The novelty in this second stochastic model is the consideration of multiple pits in a given corrosion area (called coupon hereafter). Pit initiation is modeled using a nonhomogeneous Poisson process, while the distribution of pit nucleation times is simulated using a Weibull process. The initiation time of each one of the pits produced in a coupon is described as the time to the first failure of a part of the system [33]. So, if pit initiation time is understood as the time to first failure (passive layer breakdown or inclusion dissolution) of an individual part [33], one can assume that it follows a Weibull process with distribution
The pit initiation times are computed as the realizations of a Weibull probability distribution described by (21) with both the scale parameter and the time expressed in days.
After a pit has been generated at time , it is assumed that it starts to grow. The time evolution of pit depth is then modeled using a nonhomogeneous Markov process. The transition probabilities from state to state satisfy the system of forward Kolmogorov’s equations given in (1).
If it is assumed that the pitting damage (depth) is in state 1 at the initial time (after initiation), then, for the Markov process defined by the system of (1), the transition probability from the first state to any th state during the interval , that is, , can be found analytically [28] as where is defined by (3).
From (22) it follows that, for a single pit, the probability that its depth is equal or less than state , after a time increment (), is where is the total number of states in the Markov chain.
It is easy to show [15] that expression (23) can be rewritten as
Equation (24) represents the probability of finding a single pit in a state less than or equal to state after a time interval (). It is worth noting that in (22)–(24) the pit initiation and growth processes are combined by taking into account the pitting initiation time .
As it was stated before, a direct physical meaning can be given to the functions and , which are related to the pit growth rate and pit depth, respectively. Taking into account that the dependence of pit depth on time is a power function (10), the functional dependence of and on time was assumed to be where has dimensions of distance over the th power of time and is less than 1.0.
In order to generalize from the single pit case to the -pit case, the probability that the maximum damage state is less than or equal than a given value for a time interval is estimated under the assumption that all the pits are described by parent distributions of the type given in (24). First, let us consider the simplest situation in which it is assumed that all the pits are generated at the same time . In this case, function (24) represents the cumulative distribution function of a parent population of pits. To find the extreme depth value distribution, the distribution of (24) should be raised to the th power [27]. It can be shown [15, 27] that for a large number of pits () and under a suitable variable transformation, it follows that
Substituting (24) and (25) in expression (27), and after some transformations [15], converges to a Gumbel function given by (2). The involved mathematical formalisms that led to this result can be found in Appendix A of [15].
Summarizing, when the pit initiation times are very short (smaller than the observation time) and/or when it can be considered that they are the same for all the pits, the distribution function for maximum pit depths is obtained by raising (24) to the power of the number of pits in the area of interest. The parameters of the model in this case are the number of pits and the parameters of the function: and (see expression (25)).
Consider now an alternative case, when pits are generated at different times , and expression (24) holds for each one of them. The cumulative distribution functions, , , must be combined in order to estimate the distribution of the deepest pit, under the assumption that pits nucleate and grow independently. In such a case, the probability that the deepest pit is in a state less than or equal to state , at time , can be estimated using the expression
It can be shown that this cumulative function also follows a Gumbel distribution for large [15].
In (28), expression (25) for must be substituted together with the initiation times , which are assumed to follow a Weibull distribution (21). Therefore, function in (28) includes the model parameters , , and to simulate pit initiation, together with parameters and to model pit growth through the function .
At this point, owing to the fact that functions (from (27)) and (from (28)) are extreme value distributions of the Gumbel type, it can be stated that function (24) is the parent distribution that lies in the domain of attraction of the Gumbel distribution for maxima [27]. From (24), it can be observed that is of exponential type; therefore, the parent distribution for pit depth extremes is actually an exponential function. This is in complete agreement with the findings in [7], where it was concluded, after measuring and analyzing all the pits in corroded low carbon steel coupons, that the exponential pit depth distribution fitted to the right tail of the pit depth distribution (of the Normal type, adjusted to the depths of the whole population of pits) is actually the extreme’s parent distribution, which lies in the domain of attraction of the Gumbel extreme value distribution for maxima.
Following this idea, it is soundly to shift the initial Markov state to the depth value that constitutes the starting point of the exponential distribution tail. If we consider that this exponential distribution starts at some depth value , then the Markov state that appears in (24) should be changed to a new variable
Because the pits that contribute to the extreme pit depths values are only those whose depths exceed the threshold , the number of pits that should be taken into account when applying (27) or (28) can be equated to the average number of pits whose depths exceed the threshold value . A more detailed justification of this procedure of variable change for the state in (24) can be found elsewhere [34].
The empirical average number of exceedances over threshold per coupon (to be compared to the model parameter ) can be estimated using (30), which relates with the threshold value and the parameters and of the Gumbel distribution fitted to the observed extreme pit depths. Consider
This expression is readily obtained from the existing relation between the Gumbel and exponential parameters, as has been shown in [7].
To run the proposed Markov model for extremes, (27) (pits initiate all at the same time) or (28) (pits initiate at different times) is computed and fitted to a Gumbel distribution for several times (2). The model parameters (three for the instantaneous pit initiation: , , and , and five for the generation of pits at different times: , , , , and ) are adjusted considering the experimental distributions for the deepest pits. A series of corrosion experiments is considered, each one of them consisting in the exposition of coupons of a corrodible material to a corrosion environment for a given time. After exposure, the depth of the deepest pit in each coupon is measured [6, 7], and the distribution of maximum depths is fitted with a Gumbel distribution (2). If the experiment is carried out for exposure times, and the observed behavior of depth extremes is to be described with the proposed model, the model parameters are adjusted by minimizing a total error function defined as where and are the mean and variance values of the th experimental and estimated extreme value distributions, respectively.
The model parameters are adjusted through minimization of the average value of the error function computed during Monte Carlo simulations, with .
2.2.1. Application of Markov Model for Extremes to Low Carbon Steel Corrosion Experiment
The described Markov model for extreme pit depths was applied to experimental data obtained in pitting corrosion experiments for low-carbon steel reported by Rivas et al. [7]. The details of the experiments can be found elsewhere [7].
In these experiments, groups of 20 coupons of API-5L X52 [32] pipeline steel, with cm2 of exposed area, were immersed in a corroding solution for 1, 3, 7, 15, 21, and 30 days. After the immersion, all the pits with depths greater than 5 μm were measured. The maximum pit depth in each coupon was also recorded. The Gumbel distribution (2) was fitted to the resulting maxima data sets using the Maximum Likelihood Estimator (MLE) [27].
The Gumbel probability density functions fitted to the experimental pit depth maxima for the six exposure times are displayed in Figure 5. They are represented with red solid lines. Table 3 shows the values of the Gumbel location and scale parameters for the fitted distributions.
In their work, Rivas et al. [7] concluded that, as has been previously recognized [4, 35], in low carbon steel pits initiate at the site of MnS inclusions. In the present experiment, the dissolution of these inclusions occurs in a matter of few hours after immersion, and then pits start propagating. The authors also established that the deepest pits of the entire pit population are the oldest ones, meaning that they are the defects that initiate first in time [7]. This conclusion leads to the use of (27) to model the extreme pits using Markov chains. Since the pit initiation time for this experiment has been considered to be very short, the nucleation time in (27) was set to zero.
In order to apply the Markov model, the value of the threshold parameter of (29) should be determined from the empirical corrosion data. Given that represents the depth value from which the exponential tail of the whole pit depth distribution starts, it should coincide with the previously determined [7] threshold pit depth value for the tail of the whole pit depth distribution.
When applying the so-called Peak over Threshold (POT) approach in pitting corrosion experiments, Rivas et al. proposed [7] a simplified method for determining the threshold depth value without the need of measuring the entire pit depth population. The method includes a proposition of the a priori determination of the threshold pit depth value from the fitted Gumbel distribution to the pit depth extremes. It was suggested to take as a depth value for which the Gumbel distribution shows a cumulative probability in the range from 0.00005 to 0.005. This proposition was based on the empirical fact that in the analyzed pitting corrosion experiments [7], the Gumbel distribution for maximum depth starts to rise precisely at the beginning of the exponential tail of the normal distribution of the entire pit depth population. Therefore, as the starting point of the Gumbel distribution, it was proposed to use the 0.5%, 0.05%, or even the 0.005% quantile of the experimental pit depth distributions.
Following this suggestion, and taking into consideration the determined values of the location and scale from the fitted Gumbel distribution, it is possible to determine the quantile, being , 0.0005, or 0.00005. To achieve this, the following equation, obtained through the inversion of expression (1), is used:
Thus, the value of can be used as the threshold () [7]. For the analyzed experiment we are fixing the value of to 0.0005. Substituting this value in (32) and using the estimated location () and scale () parameters for the empirical Gumbel distributions for each one of the exposure times (Table 3), six values of are obtained. These values are plotted in Figure 6 against the exposure times. In this figure, a power function fitted to the empirical threshold values is also displayed. From this curve, the values of to be substituted in (29) are taken. A more detailed proof of the correctness of this choice of can be found in [34].
Once the value of the threshold is determined, the function of (24) is raised to the th power, and the resulting function (27) is fitted to a Gumbel distribution (2). The fitted function is compared with the empirical Gumbel distribution in order to adjust the model parameters. For this particular case, the parameters to take into account are the number of pits , which should approximately match the average experimental number of exceedances over the threshold per coupon, and the parameters and , which define the function . The model parameters are refined by minimizing the error function of (31).
The minimization process gives an adjusted parameter . This value can be approximated to 8 pits. It represents the average number of pits per coupon whose depths exceed the threshold value . This fact can be confirmed by comparing the adjusted value for with the values of the empirical exceedances displayed in the fifth column of Table 3. The exceedances for each exposure period were calculated by substituting the corresponding estimated Gumbel parameters and and the threshold value (all given in Table 3) into (30). The mean value of (from Table 3) equals 8, which is in good agreement with the model prediction.
The probability density functions corresponding to the distributions obtained after raising to the power of the adjusted parameter are shown in Figure 5 for the six exposure times. In this figure, the modeled distributions can be compared with the empirical Gumbel distributions. One can see a good agreement between the modeled and experimental Gumbel distributions. Even for this case, in which the number of pits is small (see Section 2.1 and (27)), both the experimental and modeled distributions agree with the Gumbel model for maxima.
The Gumbel distribution (2) was also fitted to the functions predicted by the model using the Maximum Likelihood Estimator (MLE) [27]. In Table 3, the estimated parameters and for the fitted Gumbel distributions (to the Markov-predicted functions) are displayed together with the corresponding parameters of the empirical Gumbel distributions. The differences between the empirical and modeled parameters do not exceed 8%. The last column of Table 3 shows the average values of the two-sample K-S test performed on 1000 pairs of 50-depth-point samples (one for the empirical and one for the Markov-modeled distributions) for the six exposure periods. The smaller value among them is 0.29, which leads to not reject the null hypothesis that the two samples belong to the same distribution.
The good agreement between the modeled and the empirical Gumbel distributions can also be established from the inset of Figure 5, where the time evolution of the Gumbel mean and variance for the observed data and for the result of the proposed model are compared. These results demonstrate the suitability of the proposed model to describe the initiation and evolution of the maximum pit depths in a pitting corrosion experiment, which is of great importance in many applications such as reliability assessment and risk management.
The advantages of the proposed model compared with previous Markov models (developed by other authors) can also be established. The details regarding this topic can be found in [15].
3. Concluding Remarks
Two Markov chain models have been presented to simulate pitting corrosion. They have been developed and validated using experimental pitting corrosion data. Both models are attractive due to the existence of analytical solutions of the system of Kolmogorov’s forward equations.
The first model describes the time evolution of pit depths that correspond to the general population of defects in underground pipelines. It has been developed using a continuous-time, nonhomogenous linear growth (pure birth) Markov process, under the assumption that the Markov-chain-derived stochastic mean of the pitting damage equals the deterministic, empirical mean of the defect depths. Such an assumption allows the transition probability function to be identified only from the pitting starting time and exponent parameter. Moreover, this assumption requires that the functional form of the stochastic and deterministic instantaneous pitting rates also agree. This supports the idea that the intensities of the transitions in the Markov process are closely related to the pitting damage rate. This model permits predicting the time evolution of pit depth distribution, which is of paramount importance for reliability estimations. It is also able to capture the dependence of the pitting rate on the pit depth and lifetime. This is an advantage of the Markov chain approach over deterministic and other stochastic models for pitting corrosion. It could also be extended to investigate pitting corrosion in environments other than soils, for example, in laboratory experiments.
The second Markov chain pitting corrosion model gives account for the maximum pit depths. In it, pitting corrosion is modeled as the combination of two independent nonhomogeneous in time physical processes, one for pit initiation and one for pit growth. Both processes are combined using well-suited physical and statistical methodologies, such as extreme value statistics, to produce a unified stochastic model of pitting corrosion.
Pit initiation is described by means of a nonhomogeneous Poisson process so that a set of multiple pit nucleation events can be modeled as a Weibull process. Pit growth is modeled as a nonhomogeneous Markov process. Given that the intensity of the process is related to the corrosion rate, its functional form can be proposed from the results of the experimental tests.
Taking into account the experimental evidence that the pit depth parent distribution that leads to the distribution for extremes is the exponential tail of the pit depth population distribution, the solution of the Markov chain is shifted to the threshold pit depth value, where the exponential tail of the pit depth distribution starts. The threshold value is determined from the pit depth value from which the empirical Gumbel distribution for maximum pit depth becomes significant.
Extreme value statistics has been used to show the accuracy of the model describing the experimental observations. The model is capable of predicting not only the correct Gumbel distributions for pitting corrosion maxima in low-carbon steel, but also of estimating the number of extreme pits that has physical sense and that matches the experimental findings.
In order to simulate the whole pitting process, five model parameters are necessary. Two parameters are required to simulate pit initiation as a Weibull process; another two parameters are required to simulate pit growth with the nonhomogeneous Markov process; finally, the number of pits is necessary to combine these two processes. If the assumption that all pits nucleate instantaneously holds, as is the case of the presented experimental example, only three parameters will be necessary to fit the model to the experimental data.
Given the fact that the model parameters and assumptions do not depend on the corroded material, nor on the corrosion environment, the model is suited for different corrosion systems. The model can be easily adapted to describe situations in which the distribution of pit initiation times and the functional form of the time dependence for pit growth differ from those considered in this work.