Corrosion Potential Profile Simulation in a Tube under Cathodic Protection

The potential distribution in tubes of a heat exchanger is simulated when applying cathodic polarization to its extremes. The comparison of two methods to achieve this goal is presented: a numeric solution based on boundary elements carried out with the commercial software Beasy-GID and a semianalytical method developed by the authors. The mathematical model, the simplifications considered, and the problem solving are shown. Since both approaches use polarization curves as a boundary condition, experimental polarization curves (voltage versus current density) were determined in the laboratory under flow conditions and cylindrical cell geometry.The results obtained suggest the impossibility of extending the protection along the whole tube length; therefore, other protection methods are considered.


Introduction
The corrosion in heat exchangers in a steam cycle of an electric generating plant represents an important problem for maintenance.The morphology of attack is both generalized and localized; the latter involves drilling shutdowns for the exchange surface and contamination of the circulating fluids.The problem is solved by stopping the operation and sealing the punctured tube, which causes a decrease in heat exchange efficiency.
Condensers are multimetal systems with special conditions fluid dynamics and geometric complexity [1].It is usual to find carbon steel boxes and titanium, stainless steel, or even different copper alloy tubes and plates [2].All these alloys have different electrochemical corrosion characteristics, and due to their electrical contact a galvanic corrosion can be induced [3].In the event of cathodic protection, applied in the box, the plate is cathodically polarized, but its effect may not be seen by the entire tube and a part of it will work as an anode.Therefore, the electrochemical potential across the tube is not constant, and it is important to determine the profile in the longitudinal direction, taking as parameters the potentials applied to the edges.The current and potential distribution in the tubes is theoretically analyzed by authors as Alkire and Mirarefi [4].Astley [5], Verbrugge [6], and Song [7] address the problem considering the one-dimensional Laplace equation as the governing relationship, which is assumed to be linear kinetics; under these conditions, the problem is transformed into a distribution system consisting of an ordinary differential equation (ODE) and various boundary conditions.The above authors validate the models with respect to the results obtained through experimental prototypes.Scully and Hack [3] worked on the same research lines and performed a numerical finite element analysis, as well as a qualitative study from Wagner dimensionless number to determine the influence of plate, tube pair.The data are validated by using pilot scale prototypes.
This issue has been addressed in previous works [8, 9] using a one-dimensional approximation and, as a consequence, the problem is transformed into a system consisting of an (ODE) and several boundary conditions.
In this paper, the resolution of the equations corresponding to current and potential distribution in heat exchanger tubes is studied when a given potential is applied at the edges.We compare the results obtained by using the commercial software Beasy-GID [10] with a method developed by the authors.The method involves an analytical resolution for a simplified model based on the one-dimensional electrochemical approximation proposed by Frumkin [11].
We describe the mathematical model, the simplifications considered, and the resolution of the corresponding mathematical equations.The polarization curve-which serves as boundary condition for the proposed model-is obtained experimentally under flow conditions in a cylindrical cell, to simulate the fluid dynamic conditions within the heat exchanger.

Modeling
In order to analyze the current and potential distribution in a tubular system with electrolyte flow, a mass balance in an infinitesimal control volume  is established [12].For generic species , the accumulation is equal to the input minus the output of the species over a period with a corresponding generation (): In ( 1) that mass balance is stated,  being the flux of the species through the area  and  the molar concentration of the species.
The global balance-for all the species-can be written as follows: Multiplying these terms by the species charges () and the Faraday constant (), we can convert the mass balance into a charge balance: Considering an electrochemical system, with reaction at the interface, it can be assumed that there is no generation in the control volume (within the electrolyte) and then,   is null.Moreover, in a steady state condition the left term of (3) is also zero, resulting in where  is the current density vector.Expanding  in the transport components corresponding to migration, diffusion, and convection (Nernst-Plank development [13]), (4) can be expressed as is the ionic conductivity,   is the diffusivity of species ,  is the electric potential, and V is the advection velocity.If the solution is proposed for an isotropic medium, then ionic conductivity and species diffusivity are homogeneous properties and so their gradient is zero.Furthermore, considering a developed profile in the pipe [14,15] and an incompressible flow, the divergence of the velocity is zero, so (5) can be converted into Due to the electroneutrality in the solution bulk, the contribution corresponding to the third term (i.e., the convective transport) is zero: Considering primary or secondary current distribution [16], the species concentration gradient within the electrolyte can be neglected.Therefore (7) simplifies to The obtained equation is formally the same that corresponds to heat transfer or mass diffusion mathematical models, in steady state.The solutions obtained for these phenomena can be applied to solve the equation of the electric field in an electrochemical system, with appropriate boundary conditions.The procedures used for this solution are analytic, analogue, or digital.
Since the system considered is a tube, it is convenient to express (8) in cylindrical coordinates as follows [17]: Due to symmetry reasons, the derivatives with respect to the polar argument  can be cancelled.Then, the elliptical partial differential equation (PDE) to solve is the following: Frumkin [11] determined experimentally that the potential in an electrolyte varies only in one dimension (along the flow) when the following condition, considered in the cylindrical geometry, verifies In this formula  is the radius of the tube,  is the resistivity of the electrolyte, and  is the proportionality constant derived from the polarization curve (inverse of the linear polarization resistance).Thus if condition (11) is verified, (10) results in Under these conditions, the electrical potential gradient is parallel to the tube walls, and it is related to the current density (also parallel to the pipe walls,   ) flowing in the electrolyte.Applying Ohm's law, Applying Kirchhoff 's law it is possible to establish a relation between the current density normal to the tube wall (  ) and the longitudinal current density, as follows: Combining ( 13) and ( 14) the following equation can be obtained: The electrostatic potential  is linearly related to the electrochemical potential : where  is the chemical potential and  is the number of exchanged electrons.The boundary conditions are expressed usually through linear or semilogarithmic approximations: or The behaviors expressed by the boundary conditions ( 17) and ( 18) can be observed in the experimental curves only in bounded intervals; linear or semilogarithmic conditions cannot be extended to the entire range considered.Finally, for the boundary condition at the pipe edges, both Dirichlet (constant potential at the edges) and Neumann (no potential variation at the edges) conditions [18] can be assumed.In this paper only the first one, that is, the Dirichlet condition, is considered.

Resolution
The ODE for the potential function corresponds to (15) and the boundary condition at the tube walls is obtained from the polarization curve.
Assuming a constant potential at the tube edges, It is important to mention that, under the conditions considered in this paper, the inequality proposed by Frumkin (11) is verified.So as to solve analytically the ODE (15), the boundary condition given by the experimental curve is linearized by the following procedure.For each experimental point, a moving box, containing ten points before and ten points after, is considered.For this set, a simple linear model was fitted using least squares and then, from this model, obtained point by point, the slope of the tangent and its intersection with the vertical axis are obtained.Then, at each point a linear dependency like the following one is obtained.In order to solve the problem, a discretizing approach is performed in the following way: the tube is divided into sections such that the length of each one is decided by a preliminary observation.In fact, shorter sections are taken close to the extremes (where the greatest potential change takes place) and longer ones in the central part of the tube.In each pipe section an edge potential is assigned.The electrochemical properties for each section are considered constant and evaluated as an average of the extreme values.
If the coefficient  of ( 20) is positive, the solution for the potential takes the following form: with Imposing the boundary condition (19), the coefficients  1 and  2 are obtained as follows: Formula ( 21) is used after getting  and  from (20), using the values of  1 and  2 given by ( 22) and (23).The final result is obtained iteratively using for this purpose Microsoft Excel.
For the first iteration the corrosion potential condition is used at the interior points and the potential of ( 19) is imposed at the extreme points.As a closure criterion, the supremum norm of the difference between successive steps must be less than an arbitrary epsilon.

Use of Software Beasy-GID
With the purpose of comparing the results of the previous method, the commercial software Beasy-GID is used.This package is designed to simulate corrosion systems, using the boundary elements method [19].The information required by this software includes the system geometry (which is added as a geometrical model with the GID program), the electrical connections between the different metal components of the system (which are modeled as an electrical circuit at the software interface), and the polarization curve for the involved metals.
The boxes are made of carbon steel, 30 × 30 × 30 cm, sandblasted, and painted using a nonaqueous solvent-based acrylic with a final thickness of 120 microns.
The tube was discontinued and the sections were joined by a nonconductive tube, electrically connecting the sections using an 8 mm 2 cable in order to determine the circulating current.
At certain distances the pipe is performed introducing a zinc wire as the reference electrode.Additionally zinc reference electrodes are incorporated in the boxes (see Figure 1).In order to get the polarization of the system, platinized titanium electrodes 1  diameter and 15 cm long are incorporated in the boxes.The polarization of the system is carried out using a voltage source or controlled current (20 A maximum current and maximum voltages 30 V).The polarization of the tube sheet is chosen close to −1 V, approximately at zinc corrosion potential (the most electronegative brass component).
The electrolyte solution used in the experiment was 0.1 M sodium sulfate, pH 7, with a conductivity of 0.2 S/m.

Determination of the Boundary Conditions
A cell specifically built for this experiment is used (cell shown in Figure 2).Inside it, a half pipe is placed to be tested and it constitutes the working electrode.
In the cell are also placed the counter electrode (platinum platinized) and the reference electrode (saturated calomel electrode, SCE).The electrolyte is recirculated by a centrifugal pump.The flow is measured and controlled through a  rotameter and a ball valve, respectively.The system temperature is controlled by water circulation in a radiator located in the electrolyte reservoir.
The flow used in the experiment is 30 L/min and the electrolyte temperature is maintained at 20 ± 2 ∘ C.
The material tested is admiralty brass, the tube diameter is 0.0254 cm (1  ), and the electrode area is 20 cm 2 .
The working electrode was in contact with the electrolyte for 24 hours in sealed conditions in the cell, before the determination of the curve in flow conditions.The polarization curves are obtained by a linear scanning (5 mV/min, between −600 and +600 mV versus SCE).The experimental curves obtained are used as boundary conditions in the simulation programs.

Results and Discussion
Figure 3 shows the polarization curve obtained for admiralty brass corresponding to a flow condition of 30 L/min in sodium sulfate, pH 7.
The polarization curve in Figure 3 exhibits two plateaus (between −0.04 and −0.01 V versus SCE and above 0.005 V) and an abrupt increase in current in the vicinity of 0 V versus SCE.The response is a consequence of the chemical complexity due to the electrooxidation of the two principal alloying metals (copper and zinc), the influence of microalloying element (tin), and the cathodic reactions that take place: evolution of hydrogen and oxygen reduction [21][22][23].
Figure 4 shows the potential profile for the prototype, the simulation Beasy-GID, and the proposed simulation.
The simulation software Beasy-GID presents potential values which show a more smoothed trend in their evolution than the values corresponding to the experimental prototype.This fact leads to a prediction of the influence of polarization beyond the experimental one, measured with the prototype (approximately 100 diameters ≈ 2.5 m).
The results obtained with our proposed simulation are closer to the values determined with the prototype, but it predicts a somewhat less distance of influence for the protection: 0.15 m (approximately 6 pipe diameters).
It can be observed that the data obtained from the prototype have high variability, possibly due to the conditions of the reference electrodes used, zinc wire, and the turbulence generated at the measurement place.The influence of the imposed potential can be estimated as 20 diameters (0.5 m), which indicates, for a usual 9 m long tube, that the center of it is not reached by the cathodic protection.

Final Considerations
In the first place, the modeling has been presented for the distribution of electrochemical potential in a tube.The resolution of this model was performed using two methods: the first based on commercial software and the second one using a specific program developed by the working group.
Secondly, an experimental cell was built in order to determine current profiles and potentials for cylindrical geometries and electrolyte flow conditions.
The results obtained by the developed program were satisfactory when compared with the experimental results.
Consequently, from the technological viewpoint, both methods can predict the potential profile along the tube, changing boundary conditions, edges potentials and flow conditions, without requiring laborious and costly prototype implementation in a pilot scale.
However, the cathodic protection does not extend beyond 20 pipe diameters and this fact indicates that other methods must be used to generate an effective protection in that area.
Since it is not possible to extend the protection to the entire length of the tube, it is necessary to implement other technological solutions.One industrial possibility is the dosage of corrosion inhibitors.Particularly, in open cooling systems the ferrous sulfate dosage is commonly used [24].
As a closing remark, we recommend combining the effect of cathodic protection, which preserves the vicinity of the tube plate, with the effect of the addition of ferrous sulfate by employing consumable iron anodes [25].

Figure 1 :
Figure 1: Prototype condenser and diagram of the circuit.

Figure 4 :
Figure 4: Plot of potential versus distance to the prototype tube edge: proposed and Beasy-GID simulation.