Probabilistic model for pitting of copper canisters

Copper canisters for the final disposal of spent nuclear fuel are expected to be subject to general (uniform) corrosion and a limited degree of localised damage. Although the canister is expected to be exposed to predominantly active conditions, it is nevertheless prudent to assess the probability of pitting corrosion in passive conditions. A probabilistic pitting model is developed for aerobic conditions in water‐saturated bentonite. While aerobic, saturated conditions may be unlikely to occur if saturation is slow or O2 consumption is fast, these conditions have been assumed as a starting point because of the availability of data. The model samples from environmental conditions use a Monte Carlo approach to predict active or passive conditions. Under passive conditions, the model compares corrosion potentials to breakdown and repassivation potentials to determine pit initiations and deaths, respectively. Pit depth distributions can then be calculated. For the assumed conditions, the model predicts that active conditions dominate 90% of the time while the remaining passive conditions can lead to pitting during the first few years primarily controlled by the availability of oxygen.


| INTRODUCTION
Copper canisters in a planned KBS-3 repository at the Forsmark site in Sweden are expected to be subject to general (uniform) corrosion and a limited degree of localised damage. [1][2][3] Based on the results of both short-term laboratory experiments and longer term full-scale in situ tests at the Äspö Hard Rock Laboratory, this is expected to take the form of surface roughening rather than discrete pitting corrosion. [2] The nature of the near-surface environment is expected to result in an active copper surface rather than the passive state required for pit initiation and growth.
Although the canister is expected to be exposed to predominantly active conditions, it is nevertheless prudent to assess the probability of pitting corrosion. Pitting is an inherently stochastic process, [4] which makes a probabilistic approach appropriate. [5][6][7] In addition, there are other sources of uncertainty and variability to take into account, such as the variability in the canister surface environment from one deposition hole to another and uncertainty in the conditions that produce a passive surface.
It is well-understood that the repository environment evolves with time and that this evolution impacts the nature and time dependence of different corrosion processes. [8] Pitting of copper is usually associated with aerobic conditions. [2] Nonuniform wetting of the canister surface during the evolution of the repository saturation may also facilitate the spatial separation of anodic and cathodic processes, leading to localised damage. The repository environment will start out as aerobic and unsaturated and will eventually become saturated and anaerobic. The existence, or not, of either aerobic, saturated or of anaerobic, unsaturated conditions will depend on the relative rates of saturation and evolution of the redox conditions. This paper describes the development of a probabilistic pitting model for aerobic, saturated conditions, and the application to a particular set of environmental conditions. It is acknowledged that aerobic, saturated conditions may be unlikely to occur, however, these conditions were chosen because of the availability of input data for the model. Having developed the methodology, the intention is to assess the probability of pitting under the other sets of environmental conditions in future studies. Extended development and analysis is available in Briggs et al. [9] 2 | BACKGROUND Pitting requires that the surface be in the passive condition. [1,10] There is a general consensus that the passivation of copper is promoted by increasing pH and/or [ ] support active dissolution. [11][12][13] Increasing temperature promotes active dissolution. [13] While there is no formal definition of what constitutes a passive surface, the distinction between active (A) and passive (P) behaviour can be usefully distinguished using cyclic voltammetry. [11][12][13] There are, however, intermediate cases with characteristics of both A and P behaviour, which represent the border between the two conditions. Thus, the boundary between A and P behaviour is not necessarily well-defined and represents another source of uncertainty in the analysis.
Passive surfaces can be subject to localised film breakdown leading to pit initiation and growth. Film breakdown is characterized by the rapid increase in current during the forward (anodic-going) potential scan at a potential referred to as the breakdown potential (E B ). With increasing polarisation, the current increases as the newly initiated pit propagates. If the potential scan is reversed, the current measured on the reverse scan exceeds that on the forward scan (referred to as positive current hysteresis) and is a characteristic of pit growth. At some point, however, the current decreases to that measured during the forward scan indicating that the propagating pit has ceased to grow and that the surface has repassivated at a potential referred to as the repassivation potential (E RP ).
The following criteria for the initiation and growth of pits are used for this study: 1. The surface must be passive for pit initiation to be possible.
2. Pit initiation occurs if the corrosion potential E CORR exceeds the value of the film breakdown potential E B . 3. Once a pit is initiated, propagation occurs as long as E CORR exceeds the value of the repassivation potential E RP .
The use of the criterion E CORR > E B for pit initiation implies that the model predicts the initiation (and subsequent growth) of stable, as opposed to metastable pits. Furthermore, for both pit initiation and repassivation, it is implicitly assumed that the E B and E RP data used in the analysis account for different sources of variability associated with the stochastic nature of pitting and the experimental methods used to measure them. It should also be noted that the approach described above in Steps 1-3 is inherently conservative, as it presumes that pit repassivation can only be achieved by an appropriate shift in potential, which is only triggered by a depletion of oxygen in the system or a change in temperature. In reality, there are many physical conditions that may lead to repassivation: local depletions of oxygen in the vicinity of the pit; the production of corrosion products that may block the pit and associated movement of species; the effect of increasing spatial separation of the anode and cathode regions of the pit as it deepens, to name a few. The net effect of these processes would indicate that pit depths will be overestimated within the current approach, a point further discussed below.

| INPUT DATA
Conditions for the passivation of copper have been examined in previous studies. [2,11,12] More recently, an extensive study relevant to a repository was performed by the University of Western Ontario (UWO), with results summarised by Qin et al. [13] In contrast to the wide range of environmental conditions in the UWO study, the passivation studies of Cong et al. [11] and Imai et al. [12] were limited to a single temperature. In addition, the pH was neither systematically varied nor indicated in the study of Imai et al. [12] Qin et al. [13] considered the effects of temperature, pH, and of the concentrations of Cl − , SO 4 2− , and HCO 3 − on the passivation of copper. Cyclic voltammetry was used to distinguish active and passive behaviour resulting in A/P maps. A total of 722 voltammograms were recorded, of which 458 were categorised as displaying active behaviour and 264 passive behaviour. Single salt as well as binary and ternary salt mixtures were considered, with the compositions and concentrations selected to represent those expected in the bentonite pore water during the early saturated phase in the repository.
The temperature range was selected to cover the range expected during the thermal transient, and the pH was varied from near-neutral to a value sufficiently high to induce passivation in a given salt solution.
There are a number of sources of critical potentials for copper in relevant environments. [1,2] However, many of these studies were of limited scope with a relatively small number of E B and E RP values having been reported. The two most extensive studies are those of Cong et al. [11] and Qin et al. [13] Qin et al. [13] determined E B and E RP values from the subset of cyclic voltammograms that exhibited passive behaviour. Of the 264 voltammograms that displayed passivity, a total of 254 E B values and 203 E RP values were derived, with not all experiments exhibiting clear breakdown or repassivation phenomena.
If the initiation of a stable pit is predicted, pit growth is assumed to continue until such time that the pit repassivates, as indicated by a decrease in E CORR to the value of E RP . The model, therefore, will predict a pit growth duration from which an estimate of the depth of penetration can be determined using an empirical pit growth expression. Typically, the rate of pit (depth wise) growth is time-dependent, with the pit depth D at a time t given by D = At n where the value of the time exponent n is <1 and a value of approximately 0.5 is commonly reported. There is some mechanistic basis for a timedependent pit growth rate, for example, diffusion control or limitation by an iR drop down an ever-deepening pit would both be consistent with a value of n < 1.
For the current model, values for the constant A and the time exponent n are taken from the underground corrosion study of Denison and Romanoff. [14] This study involved the burial of, among others, copper samples in a range of soils for periods up to 14 years. [15] Pit growth expressions were selected from the soil exhibiting the deepest pits; soil #65 (Chino silt loam, pH 8, good aeration).
The value of the E CORR is used to determine (a) if pit initiation is possible and (b) whether repassivation of an actively growing pit occurs. In general, the value of E CORR is expected to be a function of a number of environmental variables and to vary with time as the repository environment evolves. Detailed reactive transport models have been developed to predict the evolution of the corrosion behaviour and of E CORR of copper canisters under aerobic conditions. [16][17][18] However, using such computational complex calculations within a probabilistic model that will be run 10s of thousands of times is not feasible and a simpler expression for E CORR is required.
King et al. [19] described a mixed-potential model for copper in O 2 -containing Cl − solutions that was applicable to a wide range of mass transport conditions. Copper anodically dissolves as Cu(I) in the form of copper(I) chloro complex ions, supported by the cathodic reduction of O 2 . For the low rates of mass transport characteristic of saturated buffer material, E CORR is a function of temperature and of the concentrations of dissolved O 2 and Cl − and is given by the following equation [19] : where E a 0 is the standard potential of the anodic reaction, k a and k −2 are rate constants for different steps in the anodic reaction, n a and n c are the number of electrons for the anodic and cathodic reactions, respectively, D is the diffusivity of the indicated species, and F, R, and T are the Faraday and gas constants and the absolution temperature, respectively. Parameter values are given in References [9,19]. Because of the dependence of E CORR on [O 2 ], we also need to define the evolution of the dissolved oxygen concentration in the repository. Canister surface temperatures were calculated for each of the 6,916 deposition holes at the Forsmark site using the thermal code described by Hedin. [20] Calculations were performed using thermal properties appropriate for saturated buffer and were carried out to a maximum time of 10,000 years. The use of thermal properties for saturated buffer results in lower temperatures than would be expected for unsaturated buffer, which has a lower thermal conductivity. As decreasing temperature favours passivity rather than active behaviour, the current analysis will tend to be conservative in the sense that there is a higher probability of passivation and localised film breakdown. For each Monte Carlo realization, the canister surface temperature-time profiles were sampled using a discrete uniform distribution. A representative set of 100 temperature profiles sampled by the probabilistic model are shown in Figure 1.
The pH of the bentonite pore water will be determined by calcite equilibrium, ion-exchange, and acid-base reactions with the montmorillonite surface, and the composition of the groundwater. A triangular distribution of pore water pH was used, with minimum, peak, and maximum values of 7, 7.4, and 10, respectively. The lower value (pH 7) reflects the lowest measured groundwater pH value at the Forsmark site, while the maximum value (pH 10) is from a measurement from a bentonite slurry activated using sodium bicarbonate to convert clay to sodium form. [21] Although such a high pH would not be expected for the pore water of compacted bentonite, this upper value was selected for the distribution to be conservative as alkaline conditions promote the passivation of copper. The peak value is defined as pH 7.4 to (a) reflect a typical value for compacted bentonite pore water and (b) not give the higher pH values too much weight. The triangular distribution was sampled at random for each Monte Carlo realization.
The distribution of chloride concentrations is selected based on the experimental data for pore water [Cl − ] from the Prototype Repository. [22] This distribution should not be considered as representative for all deposition holes at Forsmark, but more as an illustration of a single possible distribution.
The pore water concentrations of dissolved sulphate and bicarbonate are assumed to be controlled by the solubility of the respective solid phases. [8,23] A number of solubility-controlling solid phases exhibit retrograde solubility, with decreasing solubility with increasing temperature. The geochemical software PHREEQC version 3.4 [24] was used with the prhreeqc.dat thermodynamic database to estimate the equilibrium SO 4  . The pH was not defined but instead was calculated within the PHREEQC simulation. This may result in a discrepancy between the pH calculated in the solubility calculations and the sampled pH for each realization. However, any such discrepancy is expected to have a relatively minor effect given the overall simplification of the pore water chemistry in the model. The concentrations of sulphate and bicarbonate as a function of temperature from the PHREEQC calculations are shown in Figure 2.
The time-dependent dissolved O 2 concentration is required to predict the value of E CORR (Equation 1), where the [O 2 ] refers to the concentration in the buffer pore water, not that at the canister surface. [19] As such, it is taken to be represented by the overall evolution of the redox conditions within the repository. An exponential decrease in oxygen was assumed.
where [O 2 ] init is the initial dissolved oxygen concentration assumed to be 2.4 × 10 −4 mol/L, equivalent to that in an aerated saline pore water solution at 25°C. [25] Complete consumption (defined here as 0.1% remaining) was assumed to occur in between 5 weeks and 5 years. | 311 This range of times is representative of the estimated times to establish anoxic conditions in different experimental studies in granitic rock [26] and in full-scale tests. [27] Assuming a uniform distribution of α results in an oxygen concentration with the larger portion of the distribution tending toward the lower bound of 5 weeks (see Figure 3), considered to be more representative of the majority of experimental observations than the upper bound of 5 years. [26,27]

| MODEL IMPLEMENTATION
The following terminology is used to describe the execution of the probabilistic model. A realization is the term used to represent a single cycle of the overall model. A simulation represents a complete probabilistic assessment and comprises 1,000,000 individual realizations. A simulation of 1,000,000 realizations represents an average of 145 realizations for each of the 6,916 deposition holes. Multiple realizations for a given deposition hole (characterized by a specific temperature profile) account for the variability and uncertainty in the pore water [Cl − ] and pH and of the rate of O 2 consumption. A single realization represents a series of time increments from emplacement through to 100 years (or some other time, defined by the user). At the beginning of the realization, a single temperature profile is selected at random from the collection of 6,916 profiles, one for each deposition hole. In addition, values of the pore water pH and [Cl − ] are selected from the corresponding distributions and remain constant for the entire realization. The time dependence of the [O 2 ] evolution is also defined by selecting a value for the associated time constant α.
For each time step, the following processes are performed: • Determine if the surface is active or passive.
• If the surface is in the active condition, this observation, and the corresponding time, are documented, the time is incremented, and the above steps are repeated. • If the surface is passive but the pit initiation criterion (E B ) is not met, then this observation (and the associated time) is documented, the time is incremented, and the above steps are repeated. • If a pit is predicted to initiate, the extent of pit propagation is estimated based on the period for which E CORR exceeds E RP .
-Pit propagation continues until the surface repassivates. -When incrementing the time following a period of pit propagation, the time increment corresponds to the predicted period of pit propagation. -By repeating these steps, we allow for the possibility of another period of pit propagation should the various criteria be met. • Execution stops once either the corrosion allowance is exceeded (and the corresponding "failure time" is recorded) or the time reaches the maximum value specified for the realization (typically 100 years). Statistical analysis and model implementation was conducted using the base software package MathWorks MATLAB [28] with the addition of the Statistics and Machine Learning Toolbox. The model underwent several stages of development using different classification and regression methods, pit growth expressions, and time discretizations and full details are available in Briggs et al. [9] This paper focuses on a single method of analysis to determine active/passive behaviour as well as E b and E rp predictions.
The input data for the determination of active or passive conditions comprises a set of 722 experimental measurements as a function of different values of T, pH, ], and [HCO 3 − ]. Since the desired outcome is binary (active or passive), the use of classification methods to determine membership within a given group is appropriate. To condition the input data and minimize error due to large deviations in the relative magnitude between data (e.g., variations on the order of 10 −4 for the HCO 3 − concentrations while the variation in F I G U R E 3 Oxygen consumption times for a 1,000 realization Monte Carlo simulation temperature data are on the order of 10 1 -10 2 ), the data were normalized by dividing each value by the standard deviation of the respective data set. The use of a "bagged tree" machine learning algorithm was investigated. Bootstrap aggregating, or bagging, is an algorithm used in statistical classification and regression and is designed to reduce variance and minimize overfitting. The bagging method generates many training sets from the original data (with replacement) and is then used in a decision tree method to apply a classification to each new data set and finally averages the result. [29] The classification tree, or decision tree, uses properties of an observation to move toward decisions about the value of an observation. The building of the classification tree is done many times, or trained through machine learning, to create a tree with the highest accuracy for a given data set. In this way, bagging is an ensemble of methods to obtain higher predictive performance.
In general, when training a machine learning algorithm to a defined data set, overfitting can be a concern. This is analogous to the use of a high-order polynomial to fit a data set comprising one dependent variable. The fit can be so "good" that the generated polynomial goes through each data point, but this does not necessarily represent a satisfactory fit to the data. As a result, the model would accurately predict the existing data, but due to overfitting may not best represent the predictions of new data sets. Therefore, a method of cross-validation was applied during the training of the bagged tree algorithms. Specifically, nfold cross-validation was used where n = 5. In n-fold cross-validation, the data set is randomly sampled into n groups and the selected algorithm is trained on n − 1 groups with the final group being used as validation for testing the algorithm. This process is then repeated n times with each group being used as the validation data set in turn. As a result, the algorithms are resistant to overfitting as its accuracy is based on testing against new data. The final trained algorithm then uses all data groups to train before any predictions are made by the model. Using fivefold cross-validation, the bagged tree algorithm resulted in a prediction accuracy of 89.2% of correctly assigning A or P behaviour to a given set of environmental data.
Similarly, the input data for the source of critical potentials is the study of Qin et al. [13]  ]). A regression model was fitted separately to the E B and E RP data. As a starting point, the MathWorks MATLAB Regression Learner was used to test many possible regression methods. The method with the lowest root mean square error was selected for use in the full probabilistic model. Specifically, the regression tree algorithm with least squares gradient boosting [30] was selected. As with the A/P classification methods, an n-fold crossvalidation (with n = 5) was used for the implementation of regression methods.

| RESULTS
The full probabilistic analysis comprised a simulation with 1,000,000 Monte Carlo realizations over a 100-year period. Figure 4 shows the time dependence of the number of realizations exhibiting active behaviour (in blue), passive behaviour without film breakdown (in black), and passive behaviour with pit growth (in red). At any given time, approximately 90% of the realizations predict that the surface of the canister will be in the active condition based on the sampled pH and ]. This prediction is consistent with our general understanding that copper will corrode actively in the presence of compacted bentonite. [2] For the remaining 10% of the realizations, the canister surface is predicted to be passive. Film breakdown leading to pitting occurs early during the evolution of the repository environment, at a time when there is sufficient O 2 present to shift E CORR to a value more positive than E B . Pitting becomes infrequent at times beyond 2-3 years as the O 2 is consumed and E CORR shifts to more negative values. Then passive conditions become marginally more favoured, possibly because of the decrease in [SO 4 2− ] with the increasing canister temperature during this period. Sulphate exhibits a F I G U R E 4 Time dependence of the number of active, passive, and pitting realizations for the full probabilistic analysis complex effect on the active-passive behaviour of copper. [13] Sulphate will be present both in the groundwater and as gypsum, which is a common accessory mineral in commercial bentonites. Electrochemically, sulphate has been reported to promote both active behaviour, [11,12] as well as passivation. [13] This tendency toward more passive conditions is not a consequence of an increase in [HCO 3 − ] or decrease in temperature, both of which promote passivation, [13] as both parameters are evolving in the opposite direction during this prethermal peak period.
The consequences for pitting of the canister are shown in Figure 5 in the form of a histogram of predicted pit depths. The maximum pit depth for any of the 1,000,000 realizations is of the order of 1 mm, with the majority of pits being <100 μm deep.
The duration of the individual pitting events also varied significantly. The shortest pitting period was 0.0661 years (24 days) and the longest period was 6.92 years. However, there were relatively few long-lived pits, with the median and mean pit lifetimes being 0.20 and 0.41 years, respectively.
In terms of pit growth, the vast majority of pits (87%) are predicted to be <100 μm in depth. Only two of the total of over 124,000 predicted pitting events resulted in pit depths >1 mm. The deepest pit was predicted to be 1.22 mm deep, with a median and mean depth of only 39 and 53 μm, respectively.
In using an electrochemical approach to predicting the pitting behaviour of copper canisters, there are at least two questions related to the difference in the physical size of the electrodes used for the measurements and the canisters themselves.
The first refers to how many pits will there be on a given canister and what will be the distribution of pit depths. The main aim of the current model was to determine the effects of the variability and uncertainty in the environmental conditions in the repository on the pitting behaviour of the canisters. Accordingly, each realization represents the passivation and localised corrosion behaviour in a single deposition hole, in which the environmental conditions are implicitly assumed to be spatially uniform. Therefore, the current model does not provide any information regarding the number of pits that might occur on a given canister or the distribution of pit depths.
However, a full probabilistic analysis comprising 1,000,000 realizations does provide an indication of the maximum pit depth on any canister in the repository. By repeating the analysis many times, those combinations of conditions that lead to passivation, pit initiation, and an extended period of pit growth are inevitably captured.
The second scaling question is whether the use of measurements made on electrodes of the order of a few square centimeters in size can be used to predict the pitting behaviour of a canister with a surface area of approximately 17 m 2 . For passive alloys, such as stainless steel, it has been shown that the value of E B shifts to more negative values with increasing electrode size. [31] This trend is associated with the increasing number of suitable initiation sites (e.g., MnS inclusions) with increasing surface area. However, in the case of copper, the first requirement for pitting is that the surface be in the passive state. The results in Figure 4 indicate that pit initiation is predicted to occur in the majority of cases in which the surface is predicted to be passive. In other words, the probability of pitting is determined primarily by whether the surface is passive rather than by film breakdown. Thus, more negative values of E B due to the use of larger electrodes, or because of any other change in the experimental procedure such as the use of different potential scan rates or of different methods for pretreating the electrode, would not have significantly increased the probability of pit initiation. As passivation, unlike pit initiation, is not a stochastic process, the size of the electrode should have only limited influence on the experimental conditions found to induce passivation. It is concluded, therefore, that the difference in size of the electrodes and canister has no significant impact on the predicted pitting behaviour.

| CONCLUSIONS
A methodology has been developed for the probabilistic assessment of the pitting of copper canisters. The model has been applied to the assessment of pitting under F I G U R E 5 Histogram of pit depths for the full probabilistic analysis. The majority of pits are predicted to be <100 μm (10 −4 m) in depth aerobic, saturated environmental conditions. Such conditions may or may not exist during the evolution of the repository environment, depending upon the relative rates of saturation and redox evolution. Nevertheless, the techniques developed will be useful for the future assessment of localised corrosion under aerobic, unsaturated conditions, which will certainly exist in the repository.
The model is based on well-established electrochemical criteria for pitting corrosion and accounts for variability and uncertainty in both the pitting process and also the relevant repository environmental conditions. The MATLAB software has proven to be suitable for model development and implementation, where also the built-in library of statistical methods was useful for selecting and comparing different classification and regression techniques for input data. The results of the simulations suggest that: • The canister surface is more likely to be in the active, rather than the passive, condition. • The primary factors determining active versus passive behaviour are the Cl − concentration and pH of the pore water. • The canister temperature and the pore water concentrations of sulphate and bicarbonate ions seem to have a minor effect on the active-passive behaviour. • For passive conditions, pit initiation is possible but only at relatively early times when O 2 is present to sufficiently ennoble E CORR . • The duration of pitting varies from a few 10s of days to a few years.
Based on the pit growth expressions used in the analysis, the maximum pit depth is of the order of 1 mm. However, the vast majority (87%) of pits are <100 μm deep.