PDF
Carbon Dioxide Corrosion in Steel Pipes
Introduction
A flow mixture of water and carbon dioxide (CO2) passing through a steel pipe can cause significant steel corrosion. Properties such as partial pressure of CO2 and temperature affect the corrosion rate.
This model simulates the corrosion taking place on the steel surface of a pipe carrying CO2 and water mixture. The model is based on several journal papers (Ref. 1–Ref. 3).
Model Definition
The corrosion is investigated at an arbitrary position within a steel pipe through which a flow of dissolved CO2 in water passes. A 1D model is used. No variations along the length of the pipe are considered and the interaction of the mixture with the steel is confined to the boundary layer near the steel surface. The boundary layer thickness is considered to be 50 μm. The model geometry is shown in Figure 1.
Figure 1: The model geometry comprising of the boundary layer adjacent to the steel surface.
All species are assumed to be diluted in water and the mass transport is modeled by diffusion and migration. The Tertiary Current Distribution, Nernst–Planck interface with water-based charge electroneutrality is used in the model. Carbon dioxide hydration, water dissociation, proton reduction reaction, and iron dissolution reaction are accounted for; resulting in four species (in addition to H+ and OH- ions) in the model. The species and diffusion coefficients are tabulated in Table 1
Table 1: Modeled species with their respective diffusion coefficients.
Species
D (m2/s)·10-9
CO2
1.96
H2CO3
2.00
HCO3-
1.11
H+
9.31
OH-
5.26
Fe2+
0.72
The Electrode Surface boundary feature is used to calculate the corrosion potential at the steel surface. The electric potential is set to 0 at the steel surface which will result in a mixed potential condition since the electrolyte potential is left to float. The initial value of the electrolyte potential is set to −0.5 V around the free corrosion potential (Ref. 2).
Fluxes of species converted in the electrochemical reactions, ij/F (Faraday’s constant =  C/mol), are applied on the steel surface where i (SI unit: A/m2) is the current density of j number of electrochemical reactions.
For all species, uniform concentrations of species in chemical equilibrium are used as initial values for concentration (Ref. 1). Concentration of CO2 is set to its initial value at the outer point of the boundary layer.
Equilibrium Reactions
The following equilibrium reactions are present in the electrolyte:
where K1 through K3 are the equilibrium constants at 293.15 K (Ref. 1).
The water dissociation equilibrium reaction is in-built for the Tertiary Current Distribution, Nernst–Planck interface with water-based charge electroneutrality charge transport. The other two reactions are modeled using the Equilibrium Reaction domain node; one for each reaction. The Equilibrium Reaction nodes solve for one additional degree of freedom each, where the additional degree of freedom represents the local reaction rate required in order to fulfill the equilibrium expression. The equilibrium expressions are based on the reaction stoichiometry and equilibrium constant Kk according to
where ci (SI unit: mol/m3) is the concentration of species i, vik the stoichiometric coefficient of species i in reaction k. The activity of a species, a(ci) is given by dividing the concentration with 1 M.
Electrochemical Reactions
The following electrochemical reactions are present at the steel surface:
•
Iron dissolution
•
Proton reduction
Study Settings
The problem is solved with an auxiliary sweep on a stationary solver in order to investigate the impact of important parameters such as partial pressure of CO2 and temperature on corrosion.
Results and Discussion
Figure 2 displays the concentration deviation from the bulk of the four species along the boundary layer at partial pressure of CO2 of 1 bar and 20ºC. The concentration of iron ions is significantly higher at the steel surface due to the dissolution of iron. The deviation of carbon dioxide and bicarbonate ions show considerable hydration of the carbon dioxide. Carbonic acid shows little variation in concentration compared to the bulk throughout the boundary layer.
Figure 2: Deviation in concentration of the species compared to the bulk along the liquid boundary layer.
Figure 3 shows the corrosion rate of the steel surface at different partial pressures of CO2 for operating temperatures ranging from 20ºC to 80ºC. The corrosion rate is directly proportional to the corrosion current (that is, the iron dissolution current, since no other anodic reaction is considered). Increased partial pressure of CO2 and temperature increase the corrosion rate.
Figure 3: Corrosion rate in mm/year for different partial pressures of CO2 and operating temperatures range of 20ºC–80ºC.
Figure 4 and Figure 5 show the change in pH at the steel surface and bulk electrolyte, respectively, for different partial pressures of CO2 and for operating temperatures ranging from 20ºC to 80ºC. The higher pH is observed for lower partial pressure of CO2 at both the electrode surface as well as the bulk electrolyte for all operating temperatures considered in the model. The increase in pH with an increase in temperature is observed to be more significant at the electrode surface when compared to the bulk electrolyte.
Figure 4: pH at the electrode surface for different partial pressures of CO2 and operating temperatures range of 20ºC–80ºC.
Figure 5: pH at the bulk electrolyte for different partial pressures of CO2 and operating temperatures range of 20ºC–80ºC.
References
1. M. Nordsveen, S. Nesic, R. Nyborg, and A. Stangeland, “A Mechanistic Model for Carbon Dioxide Corrosion of Mild Steel in the Presence of Protective Iron Carbonate Films-Part 1: theory and Verification,” Corrosion, vol. 59, no. 5, pp. 443–455, 2003.
2. S. Nesic, J. Postlethwaite, and S. Olsen, “An Electrochemical Model for Prediction of Corrosion of Mild Steel in Aqueous Carbon Dioxide Solutions,” Corrosion, vol. 52, no. 4, pp. 280–294, 1996.
3. A. Kahyarian and S. Nesic, “On the mechanism of carbon dioxide corrosion of mild steel: Experimental investigation and mathematical modeling at elevated pressures and non-ideal solutions,” Corrosion Science, vol. 173, no. 108719, pp. 1–27, 2020.
Application Library path: Corrosion_Module/General_Corrosion/co2_corrosion
Modeling Instructions
From the File menu, choose New.
New
In the New window, click  Model Wizard.
Model Wizard
1
In the Model Wizard window, click  1D.
2
In the Select Physics tree, select Electrochemistry>Tertiary Current Distribution, Nernst-Planck>Tertiary, Water-Based with Electroneutrality (tcd).
3
Click Add.
4
In the Number of species text field, type 4.
5
In the Concentrations table, enter the following settings:
cCO2
cH2CO3
cHCO3
cFe
6
Click  Study.
7
In the Select Study tree, select General Studies>Stationary.
8
Click  Done.
Global Definitions
Load the model parameters from a text file.
Parameters 1
1
In the Model Builder window, under Global Definitions click Parameters 1.
2
In the Settings window for Parameters, locate the Parameters section.
3
Click  Load from File.
4
Browse to the model’s Application Libraries folder and double-click the file co2_corrosion_parameters.txt.
Geometry 1
The geometry consists of a single interval.
Interval 1 (i1)
1
In the Model Builder window, under Component 1 (comp1) right-click Geometry 1 and choose Interval.
2
In the Settings window for Interval, locate the Interval section.
3
In the table, enter the following settings:
Coordinates (m)
0
delta
4
Click  Build Selected.
Tertiary Current Distribution, Nernst-Planck (tcd)
Start defining the physics.
Species Charges 1
1
In the Model Builder window, under Component 1 (comp1)>Tertiary Current Distribution, Nernst-Planck (tcd) click Species Charges 1.
2
In the Settings window for Species Charges, locate the Charge section.
3
In the zcHCO3 text field, type -1.
4
In the zcFe text field, type 2.
Electrolyte 1
1
In the Model Builder window, click Electrolyte 1.
2
In the Settings window for Electrolyte, locate the Diffusion section.
3
In the DcCO2 text field, type DCO2.
4
In the DcH2CO3 text field, type DH2CO3.
5
In the DcHCO3 text field, type DHCO3.
6
In the DcFe text field, type DFe.
7
In the DcH text field, type DH.
8
In the DcOH text field, type DOH.
Initial Values 1
Set the initial values to the concentration of the species in the bulk.
1
In the Model Builder window, click Initial Values 1.
2
In the Settings window for Initial Values, locate the Initial Values section.
3
In the cCO2 text field, type cCO20.
4
In the cH2CO3 text field, type cH2CO30.
5
In the cHCO3 text field, type cHCO30.
6
In the cFe text field, type cFe0.
7
In the phil text field, type phil0.
Electrode Surface 1
Use the Dissolving-Depositing species formulation to estimate the corrosion rate. Also, set fluxes according to the electrochemical reactions at the steel surface.
1
In the Physics toolbar, click  Boundaries and choose Electrode Surface.
2
Select Boundary 1 only.
3
In the Settings window for Electrode Surface, click to expand the Dissolving-Depositing Species section.
4
Click  Add.
5
In the table, enter the following settings:
Species
Density (kg/m^3)
Molar mass (kg/mol)
s1
rho_steel
Mw_steel
Electrode Reaction 1
1
In the Model Builder window, click Electrode Reaction 1.
2
In the Settings window for Electrode Reaction, locate the Stoichiometric Coefficients section.
3
In the n text field, type 2.
4
In the νcFe text field, type -1.
5
In the Stoichiometric coefficients for dissolving-depositing species: table, enter the following settings:
Species
Stoichiometric coefficient (1)
s1
-1
6
Locate the Equilibrium Potential section. In the Eeq,ref(T) text field, type Eeq_ref_Fe.
7
Locate the Electrode Kinetics section. In the i0,ref(T) text field, type i0_ref_Fe.
8
In the αa text field, type alphaa_Fe.
Electrode Surface 1
In the Model Builder window, click Electrode Surface 1.
Electrode Reaction 2
1
In the Physics toolbar, click  Attributes and choose Electrode Reaction.
2
In the Settings window for Electrode Reaction, locate the Equilibrium Potential section.
3
In the Eeq,ref(T) text field, type Eeq_ref_H2.
4
Locate the Electrode Kinetics section. In the i0,ref(T) text field, type i0_ref_H2.
5
In the αa text field, type alphaa_H2.
Equilibrium Reaction 1
Set two equilibrium reactions.
1
In the Physics toolbar, click  Domains and choose Equilibrium Reaction.
2
In the Settings window for Equilibrium Reaction, locate the Domain Selection section.
3
From the Selection list, choose All domains.
4
Locate the Equilibrium Condition section. In the Keq text field, type KCO2H.
5
Locate the Stoichiometric Coefficients section. In the νcCO2 text field, type -1.
6
In the νcH2CO3 text field, type 1.
Equilibrium Reaction 2
1
In the Physics toolbar, click  Domains and choose Equilibrium Reaction.
2
In the Settings window for Equilibrium Reaction, locate the Domain Selection section.
3
From the Selection list, choose All domains.
4
Locate the Equilibrium Condition section. In the Keq text field, type KH2CO3.
5
Locate the Stoichiometric Coefficients section. In the νcH2CO3 text field, type -1.
6
In the νcHCO3 text field, type 1.
7
In the νcH text field, type 1.
Concentration 1
Set bulk concentrations at the rightmost boundary.
1
In the Physics toolbar, click  Boundaries and choose Concentration.
2
Select Boundary 2 only.
3
In the Settings window for Concentration, locate the Concentration section.
4
Select the Species cCO2 check box.
5
Select the Species cFe check box.
6
In the c0,cCO2 text field, type cCO20.
7
In the c0,cFe text field, type cFe0.
Mesh 1
Build a user-defined mesh with a maximum element size in the domain of 1e-6 and at the leftmost boundary 1e-7.
1
In the Model Builder window, under Component 1 (comp1) click Mesh 1.
2
In the Settings window for Mesh, locate the Sequence Type section.
3
From the list, choose User-controlled mesh.
Size
1
In the Model Builder window, under Component 1 (comp1)>Mesh 1 click Size.
2
In the Settings window for Size, locate the Element Size section.
3
Click the Custom button.
4
Locate the Element Size Parameters section. In the Maximum element size text field, type 1e-6.
Size 1
1
In the Model Builder window, right-click Edge 1 and choose Size.
2
In the Settings window for Size, locate the Geometric Entity Selection section.
3
From the Geometric entity level list, choose Boundary.
4
Select Boundary 1 only.
5
Locate the Element Size section. Click the Custom button.
6
Locate the Element Size Parameters section.
7
Select the Maximum element size check box. In the associated text field, type 1e-7.
8
Click  Build Selected.
Study 1
1
In the Model Builder window, click Study 1.
2
In the Settings window for Study, locate the Study Settings section.
3
Clear the Generate default plots check box.
Parametric Sweep
Use parametric sweep to investigate the impact of CO2 partial pressure and temperature.
1
In the Study toolbar, click  Parametric Sweep.
2
In the Settings window for Parametric Sweep, locate the Study Settings section.
3
Click  Add twice.
4
In the table, enter the following settings:
Parameter name
Parameter value list
Parameter unit
PCO2 (Partial pressure of CO2 in bar)
0.01 0.1 range(0.5,0.5,3)
 
T (Operating temperature)
293.15[K] 323.15[K] 353.15[K]
K
5
From the Sweep type list, choose All combinations.
Solution 1 (sol1)
1
In the Study toolbar, click  Show Default Solver.
2
In the Model Builder window, expand the Solution 1 (sol1) node.
3
In the Model Builder window, under Study 1>Solver Configurations>Solution 1 (sol1) click Stationary Solver 1.
4
In the Settings window for Stationary Solver, locate the General section.
5
In the Relative tolerance text field, type 1e-6.
6
In the Study toolbar, click  Compute.
Results
The following steps reproduces the plots from the Results and Discussion section.
Concentrations
1
In the Home toolbar, click  Add Plot Group and choose 1D Plot Group.
2
In the Settings window for 1D Plot Group, type Concentrations in the Label text field.
3
Locate the Data section. From the Dataset list, choose Study 1/Solution 1 (sol1).
4
From the Parameter selection (PCO2) list, choose Last.
5
From the Parameter selection (T) list, choose First.
6
Click to expand the Title section. From the Title type list, choose Manual.
7
In the Title text area, type PCO2=1 bar, T=20<sup>\circ</sup>C.
8
Locate the Plot Settings section.
9
Select the x-axis label check box. In the associated text field, type Distance from steel surface (m).
10
Select the y-axis label check box. In the associated text field, type Concentration deviation (mol/m<sup>3</sup>).
Line Graph 1
1
Right-click Concentrations and choose Line Graph.
2
In the Settings window for Line Graph, locate the Selection section.
3
From the Selection list, choose All domains.
4
Locate the y-Axis Data section. In the Expression text field, type cCO2-cCO20.
5
Click to expand the Title section. From the Title type list, choose None.
6
Locate the x-Axis Data section. From the Parameter list, choose Expression.
7
In the Expression text field, type x.
8
Click to expand the Legends section. Select the Show legends check box.
9
From the Legends list, choose Manual.
10
In the table, enter the following settings:
Legends
CO<sub>2</sub>
Line Graph 2
1
Right-click Line Graph 1 and choose Duplicate.
2
In the Settings window for Line Graph, locate the y-Axis Data section.
3
In the Expression text field, type cH2CO3-cH2CO30.
4
Locate the Legends section. In the table, enter the following settings:
Legends
H<sub>2</sub>CO<sub>3</sub>
Line Graph 3
1
Right-click Line Graph 2 and choose Duplicate.
2
In the Settings window for Line Graph, locate the y-Axis Data section.
3
In the Expression text field, type cHCO3-cHCO30.
4
Locate the Legends section. In the table, enter the following settings:
Legends
HCO<sub>3</sub><sup>-</sup>
Line Graph 4
1
Right-click Line Graph 3 and choose Duplicate.
2
In the Settings window for Line Graph, locate the y-Axis Data section.
3
In the Expression text field, type cFe-cFe0.
4
Locate the Legends section. In the table, enter the following settings:
Legends
Fe<sup>2+</sup>
5
In the Concentrations toolbar, click  Plot.
Corrosion rate
1
In the Home toolbar, click  Add Plot Group and choose 1D Plot Group.
2
In the Settings window for 1D Plot Group, type Corrosion rate in the Label text field.
3
Locate the Data section. From the Dataset list, choose Study 1/Solution 1 (sol1).
4
Locate the Title section. From the Title type list, choose None.
5
Locate the Plot Settings section.
6
Select the x-axis label check box. In the associated text field, type Partial pressure of CO2 (bar).
7
Select the y-axis label check box. In the associated text field, type Corrosion rate (mm/year).
Point Graph 1
1
Right-click Corrosion rate and choose Point Graph.
2
Select Boundary 1 only.
3
In the Settings window for Point Graph, click Replace Expression in the upper-right corner of the y-Axis Data section. From the menu, choose Component 1 (comp1)>Tertiary Current Distribution, Nernst-Planck>Dissolving-depositing species>tcd.vbtot - Total electrode growth velocity - m/s.
4
Locate the y-Axis Data section. From the Unit list, choose mm/yr.
5
Locate the x-Axis Data section. From the Axis source data list, choose PCO2.
6
From the Parameter list, choose Expression.
7
In the Expression text field, type PCO2.
8
Click to expand the Legends section. Select the Show legends check box.
9
From the Legends list, choose Evaluated.
10
In the Legend text field, type T=eval(T) K.
11
In the Corrosion rate toolbar, click  Plot.
pH at surface
1
In the Model Builder window, right-click Corrosion rate and choose Duplicate.
2
In the Settings window for 1D Plot Group, type pH at surface in the Label text field.
3
Locate the Plot Settings section. In the y-axis label text field, type pH at surface.
Point Graph 1
1
In the Model Builder window, expand the pH at surface node, then click Point Graph 1.
2
In the Settings window for Point Graph, locate the y-Axis Data section.
3
In the Expression text field, type tcd.pH.
4
In the pH at surface toolbar, click  Plot.
pH at bulk
1
In the Model Builder window, right-click pH at surface and choose Duplicate.
2
In the Settings window for 1D Plot Group, type pH at bulk in the Label text field.
3
Locate the Plot Settings section. In the y-axis label text field, type pH at bulk.
Point Graph 1
1
In the Model Builder window, expand the pH at bulk node, then click Point Graph 1.
2
In the Settings window for Point Graph, locate the Selection section.
3
Click to select the  Activate Selection toggle button.
4
In the list, select 1.
5
Select Boundary 2 only.
6
In the pH at bulk toolbar, click  Plot.