- Email: [email protected]
Modelling of erosion–corrosion in practical geometries
Corrosion Science 51 (2009) 769–775
Contents lists available at ScienceDirect
Corrosion Science journal homepage: www.elsevier.com/locate/corsci
Modelling of erosion–corrosion in practical geometries Cian Davis *, Patrick Frawley Department of Mechanical and Aeronautical Engineering, University of Limerick, Ireland
a r t i c l e
i n f o
Article history: Received 25 August 2008 Accepted 23 December 2008 Available online 17 January 2009 Keywords: B. CFD B. Wear prediction C. Erosion–corrosion C. Sherwood number
a b s t r a c t The limiting factor in certain instances of erosion–corrosion of steel is the presence of dissolved oxygen in the solution and the transfer of this oxygen to the reacting surface. Computational Fluid Dynamics (CFD) can be used to calculate oxygen diffusion throughout the flow and its transfer to the reacting surface. This was used in a computational model to calculate wear and validated against experimental results, for the first time, of erosion–corrosion wear in a contracting–expanding geometry. It was found that in order to correctly predict erosion–corrosion wear, Sherwood number independent grids were required providing a new metric to evaluate turbulent erosion–corrosion modelling. The predicted wear profile matched very closely with experimental results and overall matching was very good. Downstream of the flow expansion, erosion–corrosion wear was under-predicted. The disparity is due to detached flow for most of this sector where under-predicted radial velocities decrease transfer of oxygen to the reacting surface. This under-prediction is apparent in the downstream section due to the larger relative magnitude of the radial velocity in this sector. Ó 2009 Elsevier Ltd. All rights reserved.
1. Introduction Flowing fluids can have significant effects on rates of corrosion in a number of ways including accelerated mass transfer of reactants and corrosion products [1], fluid shear [2] or impingement of solid particles in the fluid which disrupt protective layers [3]. The ability to model and predict wear processes allows service intervals to be better timed so as to reduce unnecessary checks while not being subject to costly downtime due to equipment failure. Also, modelling ability can be applied at design stage to reduce the susceptibility of parts to wear. Due to the complexity of erosion–corrosion and significant number of inter-related variables, many attempts have been made to identify the controlling factors or relate the process to some other easily measurable quantity. Postlethwaite et al. [4] stated that erosion–corrosion of carbon steel piping in saline-sand slurry was controlled by mass transfer of oxygen to the reacting surface and that the role of the solid particles was to prevent formation of a complete rust film. Lotz and Postlethwaite [5] suggested this increased mass transport of oxygen was due to increased turbulent transfer across the mass-transfer boundary layer. Postlethwaite and Lotz [6] found a significant effect of surface roughness on oxygen mass-transfer in a similar sand-based slurry. Following on from Nesic and Postlethwaite [7,8], Keating and Nesic [9] adapted the work of Lotz [10] to formulate a model for erosion–corrosion. * Corresponding author. Tel.: +353 868045259; fax: +353 61213529. E-mail address: [email protected] (C. Davis). 0010-938X/$ - see front matter Ó 2009 Elsevier Ltd. All rights reserved. doi:10.1016/j.corsci.2008.12.025
This paper will apply the model of Keating and Nesic [9] to an contracting–expanding pipe geometry in order to prove the ability of the model in practical geometries. The contracting–expanding pipe geometry was chosen as detailed experimental data of both wear and flow parameters are available and complex flow characteristics are evident. Also, the optimum grid density will be determined for accurate predictions. 2. Modelling 2.1. Fluid modelling Modelling of the fluid is undertaken using Computational Fluid Dynamics (CFD). All fluid modelling was undertaken with FLUENT v6.3 code, a commercial CFD package. FLUENT uses a finite volume method in its CFD code. This is a well established and well utilised method, found in many CFD codes. The governing fluid equations are applied to a fixed zone in space, known as a control volume. These equations are also known as the conservation equations as they describe the conservation of mass, momentum and energy over the control volume. Detailed explanations of the model are provided by Versteeg and Malalasekera [11], Patankar [12] and Ferziger and Peric [13]. A derivation of the governing equations in an Eulerian formulation is give by Drew [14]. 2.1.1. Turbulence modelling The k e model was first proposed by Launder and Spalding [15] and has been subject of various improvements and modifications since. Yakhot and Orszag [16] developed an improvement to
770
C. Davis, P. Frawley / Corrosion Science 51 (2009) 769–775
the k e model by utilising renormalisation group (RNG) theory. In this RNG k e model, the constants that were derived empirically for the standard k e model are calculated using RNG theory. Also, additional terms and functions are added to the transport equations for k and e. By using this RNG-based approach, the model does not have to be modified to suit the type of flow being modeled and gives better predictions for separated flow or flows in curved geometries while only a small amount of extra computational effort is required – typically 10–15% [17]. 2.1.2. Wall effects There are two methods that can be utilised to resolve the near wall velocity profile numerically. The rigorous approach is to fully resolve the viscosity-affected region. This requires at least ten cells within the Rey < 200 region [18]. At high Reynolds numbers, the viscous sub-layer is so thin, it is difficult to resolve fully [13]. In a large geometry, addition of these extra grid points can impose a significant extra computational effort. This second approach is to utilise a wall function approach. This negates the need for excessive grid points at the boundary and is computationally less expensive. Wall functions are approximations based on simple flows and difficulties may arise in separated flows. Calculation of species transport brings with it the issue of the mass-transfer boundary layer, which is approximately an order of magnitude less than the velocity boundary layer [19]. The FLUENT CFD code includes a hybrid near-wall modelling method that utilises both approaches described above. In the case of a model with a fine mesh at the wall ðyþ 1Þ, the viscous sublayer is fully resolved. However, a model with a fine mesh throughout the model may not be desirable or practicable, for example at areas of negligible velocity or little interest. In these cases, (typically, yþ > 5) a coarse mesh can be applied and a blended wall function is used [20]. Therefore, areas of separation or interest can be easily investigated by increasing resolution of the mesh in that area, without the requirement of a fine mesh throughout the model, or even at all boundaries. This procedure provides a compromise between modelling error and computational effort and is employed to model near wall effects in this study. 2.2. Erosion–corrosion modelling Following on from Nesic and Postlethwaite [7,8], Keating and Nesic [9] adapted the work of Lotz [10] to formulate a model for erosion–corrosion. þþ
ð1Þ
The flux of the corrosive species (oxygen, in this case and denoted by J O2 ) is related to the mass transfer coefficient km , the wall concentration of oxygen, C wO2 and the bulk concentration of oxygen, C bO by 2
J O2 ¼ km ðC bO C wO2 Þ
ð2Þ
2
In the model C wO2 is set to zero as all the oxygen diffusing to the wall is consumed by the corrosion reaction. According to Eq. 1 two moles of Fe react with every mole of O2 . Therefore, the equation for the flux of Fe is given by
ð3Þ
2
2 1
The flux of Fe, given in kmol m s can be converted to a corrosion rate by use of the molar mass of Fe, M Fe (kg/kmol) and its density qFe (kg/m3). The result is an equation giving corrosion rate in mm/year.
CR ¼
2km C bO M Fe 2
qFe
J O2 ¼
D ðC fcO C wO2 Þ 2 Dy
ð5Þ
where Dy is the distance of the first cell centre from the wall, D is the mass diffusion coefficient (m2 s1) and C fcO is the concentration 2 at that point, obtained in CFD by solving the species transport equation. By assuming the concentration at the wall is zero and then combining Eqs. 2 and 5, km is given by as
km ¼
D C fcO2 Dy C b O
ð6Þ
2
Combining Eqs. 4 and 6 gives a erosion–corrosion model similar to that developed by Stack et al. [21].
CR ¼
2DC fcO M Fe 2
qFe Dy
24 60 60 365 103
ð7Þ
By using the assumption of Postlethwaite et al. [4] that the protective layer is removed by the combined action of the fluid and solid phase, the erosion–corrosion process can be assumed to be under the control of oxygen mass-transfer to the corroding surface. Both C fcO and Dy in Eq. 7 can be calculated with the aid of CFD. The 2 rate of erosion–corrosion wear can then be calculated using Eq. 7. The Sherwood number is a ratio of length-scale to diffusive boundary layer thickness and is defined as:
Sh ¼
km d D
ð8Þ
Eq. 4 can be written in terms of Sherwood number in a new correlation for erosion–corrosion wear
CR ¼
2DShC bO M Fe 2
qFe d
24 60 60 365 103
ð9Þ
The only flow-dependent variable in Eq. 9 is the Sherwood number itself. By definition, therefore, accurate erosion–corrosion predictions require a grid which is Sherwood number independent. 3. Model validation 3.1. Fluid model
2Fe ! 2Fe þ 4e O2 þ 2H2 O þ 4e ! 4OH
J Fe ¼ 2km C bO
The mass transfer coefficient km can be calculated in CFD in the following manner [19]. If the first cell from the wall is within the mass transfer boundary layer, the mass transfer between the two is governed by purely diffusive effects and Eq. 2 can be rewritten as
24 60 60 365 103
ð4Þ
Before CFD can be used to predict erosion–corrosion, flow characteristics must be validated for a flow regime similar to intended setup. Experimental investigations of turbulent flow through a sudden expansion have been conducted using laser Doppler anemometry (LDA) by numerous authors [22–24]. Difficulties in accurately measuring close to the wall were often reported. The near-wall layer is of particular interest in this erosion–corrosion model as it relies on calculating mass-transfer across the boundary layer. Founti and Klipfel [25] undertook a study in a similar geometry using diesel oil. Detailed experimental results of axial and radial velocity as well as turbulent measurements were provided at numerous locations downstream of the expansion. Results provided sufficient resolution in the near-wall region and detailed fluid parameters were provided. These results were selected for comparison with the numerical study. 3.1.1. Experimental description Founti and Klipfel conducted their experiments in a closed loop system with a vertically mounted working section. Upstream, the pipe has an internal diameter, d, of 25.5 mm and downstream,
C. Davis, P. Frawley / Corrosion Science 51 (2009) 769–775
771
where C l is given as 0.0845 from RNG theory and us is the friction velocity given by
us ¼
sw qc
0:5 ð13Þ
The relationship between wall shear stress and bulk inlet velocity is given as
sw ¼
q ub d 0:25 f ¼ 0:079 c
the diameter, d, was twice as large at 51 mm. The working section incorporated an expansion with a ratio of 2:1 and fully developed flow in the direction of gravity. A schematic diagram and indication of notation is given in Fig. 1. A mix of diesel oils was used in order that the refractive index of the working fluid matched the walls of the test section. Experimental parameters are given in Table 1. Measurements were taken using LDA and the error in results was given as 1% and 2% for the mean and turbulent measurements respectively. 3.1.2. Computational model In the computational modelling of the problem, steady, incompressible, isothermal, turbulent flow of diesel in a sudden expansion was assumed. Since the flow did not include a swirl component, the flow was modelled as 2D axisymmetric. This allowed significant resolution in areas of interest while maintaining a computationally efficient model. At the inlet to the computational domain, the axial velocity, turbulent kinetic energy and turbulent dissipation rate must be specified. In order to achieve the best approximation of experimental conditions, a fully turbulent flow was specified based on the axial peak velocity, up , at the centreline. The axial velocity at the inlet is related to the distance from the pipe centre line by a power law suggested by Prandtl [26] for fully developed turbulent pipe flow.
ð10Þ
Where d is the radius of the downstream pipe and r is the radial distance from the centreline. By integrating Eq. 10, an expression relating the bulk velocity, ub , to the peak velocity, up , can be obtained
up ¼ 1:224ub
ð11Þ
As an axisymmetric model was used, no tangential velocities are apparent in the model. A zero radial velocity was assumed at the inlet. For fully-developed flow, a linear variation from the wall to a free stream value at the centreline is assumed for the turbulent kinetic energy. The near-wall value of k is given by
u2s ffi knw ¼ pffiffiffiffiffi Cl
ð12Þ
Table 1 Experimental properties and operating conditions [25]. Model properties Density Viscosity Mass flow rate Upstream Reynolds number Downstream Reynolds number Bulk velocity (upstream of expansion) Bulk velocity (downstream of expansion)
ð14Þ
where f is the friction factor and is estimated for smooth pipes using the Blasius equation [27] as
Fig. 1. Diagram of experimental setup [25].
r 1=7 ux ¼ up 1 d
f qc u2b 2
830 kg m3 0.0043201 kg m1 s1 1.7304 kg s1 20,000 10,004 4.082 m s1 1.021 m s1
l
ð15Þ
Eq. 15 is valid for Reynold’s numbers between 4000 and 105 . Turbulent kinetic energy at the centreline us given by
kp ¼ 0:002u2p
ð16Þ
By combining Eqs. 12 and 16, a linear relationship for turbulent kinetic energy is given as
r k ¼ kp þ ðknw kp Þ d
ð17Þ
For fully-developed flow, the turbulent dissipation rate is given as 3=2
e¼
C 3=4 l k l
ð18Þ
where l is the mixing length and is set as jðd rÞ in the near-wall region and 0:09d in the outer layer [27]. By applying the experimental conditions given in Table 1 to Eqs. 10, 17 and 18, inlet profiles are generated for the models, shown in Fig. 2. In order to apply these profiles, a custom subroutine was used. 3.1.3. Fluid simulation results Overall, very good matching of experimental and computational data is obtained. A comparison of experimental and computational results for axial velocity between 4 mm and 50 mm downstream of the expansion are given in Fig. 3 and radial velocity in Fig. 4. Slight over-prediction of axial velocity is noticeable for values of r=0:5d of greater than 0.6 at x ¼ 75mm and x ¼ 100mm and for values of r=0:5d of greater than 0.2 at x ¼ 200mm. Overall, there is poor agreement between experimental and computational radial velocity results. This deviation can be accounted for by two main factors. First, the magnitude of the radial velocity is at least an order of magnitude less than the corresponding magnitude of the axial velocity. Second, the magnitude of the fluctuating component of the radial velocity is approximately double the magnitude of the mean radial velocity component. This results in a 200% turbulence intensity. Reynolds averaged turbulence models are unable to solve such high turbulence intensities [11]. Other studies have also reported difficulties in correctly predicting radial velocities in a similar flow setup [7,19,28]. Due to the small magnitude of the radial component relative to the axial component, the deficiencies in modelling the radial velocity do not have a significant effect on the overall fluid solution. 3.2. Erosion–corrosion model By using CFD to provide inputs for the equations derived above, an erosion–corrosion model was developed using solution variables that are currently solved by the CFD solver, such as oxygen concentration and cell wall distance. The computational model
772
C. Davis, P. Frawley / Corrosion Science 51 (2009) 769–775
Radial position (r/0.5 )
1
0.8
0.6
0.4
0.2
0 0
1
2
3
4
5
0
0.1
0.2
0.3
0
0.4
Turbulent kinetic energy (m 2/s2)
Axial Velocity (m/s)
1000 2000 3000 4000 5000 6000
Turbulent dissipation rate (m2/s3)
Fig. 2. Flow profiles at inlet.
1 x = 4mm
x = 6mm
x = 9mm
x = 50mm
Radial position (r/0.5d)
0.8
0.6
0.4
0.2
0 -1 0
1
2
3
4
5
6 -1 0
1
2
3
4
5
6 -1 0
1
2
3
4
5
6 -1 0
1
2
3
4
5
6
Axial velocity (m/s) Fig. 3. Comparison of experimental ðÞ and computational results (–) of axial velocity.
1 x = 9mm
x = 25mm
0.1 0.2 0.3 -0.3 -0.2 -0.1
0
x = 50mm
x = 75mm
Radial position (r/0.5d)
0.8
0.6
0.4
0.2
0 -0.2 -0.1
0
0.1
-0.2
-0.1
0
0.1
0.2 0 0.05 0.1 0.15 0.2 0.25 0.3
Radial velocity (m/s) Fig. 4. Comparison of experimental ðÞ and computational results (–) of radial velocity.
was validated by comparing the results with the experimental results of Lotz and Postlethwaite [5]. 3.2.1. Experimental description A saline solution with 2% concentration of sand was pumped through a vertically-mounted contracting–expanding nozzle made of carbon-steel experimental specimens (Fig. 5) with flow opposing gravity. The sections were separated with acrylic to provide electri-
cal insulation. This setup allowed investigation of erosion–corrosion in each section. The contracting–expanding geometry allowed investigation of the effects of detaching and reattaching flow. Experimental conditions are given in Table 2. 3.2.2. Model verification Verification of the model is especially important in the current context and is achieved by comparing the results against the
773
C. Davis, P. Frawley / Corrosion Science 51 (2009) 769–775 Table 3 Erosion–corrosion model validation grid sizes. Grid
Number of cells
Solution time (hours)
Full resolution of laminar sublayer Intermediate mesh Wall function mesh
368,900 60,870 8,349
52 15 4
Fig. 5. Example of test section in flow loop test rig [5].
experimental work of Lotz and Postlethwaite [5]. In the opinion of Keating [19], accurate calculation of wall oxygen concentrations in CFD required full resolution of the mass transfer boundary layer. Levich [29] developed an expression to approximate the width of the mass-transfer boundary layer relative to the fluid boundary layer
13 1 dmt 3 df Sc
ð19Þ
where dmt is the width of the mass-transfer boundary layer, df is the width of the fluid boundary layer and Sc is the Schmidt number, which relates momentum diffusivity to mass diffusivity and is given by
Sc ¼
m D
¼
l qD
ð20Þ
For the current system, the Schmidt number is 548, resulting in a mass-transfer boundary layer approximately 0.07 times the width of the fluid boundary layer. On this basis, full resolution of the mass transfer boundary layer would require yþ 0:1. Making a more precise comparison of the fluid and mass-transfer boundary layers is difficult, since the thickness of the masstransfer is arbitrary [30]. It is essentially another way to write the mass-transfer co-efficient. The best way to ensure proper resolution of the mass-transfer boundary layer is to verify the grid is Sherwood number independent. Resolution of the viscous sublayer is accomplished in CFD using so-called low Reynolds number (LRN) models. Applying a mesh with a yþ 0:1 throughout the computational domain, often results in an excessive computational requirement. In most cases, wall meshes are not this fine and the boundary layer is approximated with the use of wall functions. The RNG k e model offers enhanced wall functions, which cover the full range of modelling approaches [18]. For a value of Rey < 200, the laminar sublayer is fully resolved, where Rey is the turbulent Reynolds number given by:
pffiffiffi
qy k Rey ¼ l
ð21Þ
where y is the distance from the cell centre to the wall and k is the turbulent kinetic energy. The model offers accurate matching with small values of yþ and reasonable representation where yþ falls inside the wall buffer region ð3 < yþ < 10Þ.
In order to assess the performance of each boundary layer resolution method, three grids were constructed with varying nearwall mesh densities. The first mesh is designed to accurately resolve the laminar sublayer, while the third will use wall functions. The near-wall density of the second grid falls between these two method and is designed to evaluate model performance where optimum resolution of the laminar sublayer is not achieved. Details of grid size are given in Table 3. The three grids are shown in Fig. 6. 3.2.3. Computational setup Steady, incompressible, isothermal, turbulent flow of saline solution in a sudden contraction–expansion geometry was assumed. Since the flow did not include a swirl component, the flow was modelled as 2D axisymmetric. The same inlet profiles as detailed previously were used for k, e and axial velocity (Fig. 2). The experimental setup given by Lotz and Postlethwaite [5] specify that the water was saturated with oxygen. Lewis [31] provides tables of oxygen saturation in water for varying temperature, pressure and chloride concentration. At a standard pressure of 10 kPa and a temperature of 30 °C, the oxygen saturation is given as 7.4 mg/L. Correcting for a 3% saline concentration gives an oxygen saturation concentration of 6.23 mg/L. This model assumes that any oxygen reaching the wall is consumed by the corrosion reaction. On this basis, the concentration of oxygen at the wall is set to zero. The infinite dilution diffusion coefficient can be estimated using a model for solutes in aqueous solutions proposed by Hayduk and Minhas [32].
D ¼ 1:25 1012 V 0:19 0:292 T 1:52 l A where
is given as:
9:58 ¼ 1:12 VA
Model properties 340,000 170,000 13.2 m s1 3.3 m s1 2% by volume sand 430 lm Carbon steel (AISI MT-1015) 30 °C 3% by weight
ð23Þ
and V A is the solute molar volume at its normal boiling point (given as 25.6 cm3 mol1 for oxygen), T is temperature and l is the viscosity of water in centipoise. This results in a diffusivity coefficient of 1:83 109 m2 s1 for the system outlined. The wear rate in the experiments of Lotz and Postlethwaite [5] is calculated by measuring the change in weight of each experimental section over the experiment. There are twenty-seven experimental data points versus between 251 and 1904 data points in the computational model over the experimental section,
Table 2 Experimental properties and operating conditions given by Lotz and Postlethwaite [5].
Reynolds number (narrow section) Reynolds number (wide section) Bulk velocity (narrow section) Bulk velocity (wide section) Particle concentration Average particle diameter Experimental material Temperature NaCl concentration
ð22Þ
Fig. 6. Erosion–corrosion model validation grids.
774
C. Davis, P. Frawley / Corrosion Science 51 (2009) 769–775
depending on the grid size. In order to compare the data from the two sources, the corrosion rate for each computational data point is calculated using Eq. 7. This wear depth is used to calculate the volume of wear around the section in question. By multiplying this volume by the density, the mass moved can be calculated. As the width and radius of the computational section is known, the area can easily be calculated. By dividing the mass of material removed by the area of the section, the wear rate can be calculated over the same sections and in the same units as used by Lotz and Postlethwaite [5]
Wear rate ¼
pwððr CRÞ2 r2 Þq qðCR2 2rCRÞ ¼ 2r 2prw
ð24Þ
where w is the width of the experimental section, CR is the corrosion wear depth rate defined by Eq. 7, r is the pipe radius and q is the density of steel. In order to confirm that the contribution of the solid phase to wear is only to remove the corrosion product, a separate model was run to calculate solid-particle erosion wear. The discretephase model (DPM) in FLUENT tracks particles and impacts through the computational domain. This was combined with the erosion model developed by Oka et al. [33] and supported by experimental results from O’Mahony et al. [34] for the erosion profile of steel being impacted by sand. This gives an accurate computational prediction of solid-particle erosion. 4. Results and discussion The results of the computational study detailed in the previous section were reviewed in order to evaluate the accuracy of the erosion–corrosion model. Additionally, the results provided by each grid size were evaluated to quantify the effect of boundary layer resolution method on the prediction accuracy. The rate of erosion–corrosion was calculated as according to Eq. 24 for each mesh. The data are plotted against the experimental results of Lotz and Postlethwaite [5] in Fig. 7. The predicted erosion–corrosion rates are compared with experimental results over three sectors. The first three experimental sections make up the first sector up until the flow constriction. The second sector makes up the flow constriction consisting of the next eleven experimental sections. The third sector is the flow expansion after the constriction consisting of the final thirteen
experimental sections. They are numbered as sectors 1, 2 and 3, respectively in Fig. 7 to aid identification. It is immediately apparent from Fig. 7 that the wear predicted by the models with the two coarsest grids (60k cells and 9k cells) is greatly under-predicted. Wear in these two models is approximately an order of magnitude lower than both experimental results and results from the model with the finest grid. It is apparent, therefore, that full resolution of the mass transfer boundary layer is required for prediction of erosion–corrosion wear using the current model. Neither wall functions nor partial resolution of either the mass transfer or fluid boundary layers correctly predict transfer of oxygen to the reacting surface. On this basis, subsequent comparison of predicted erosion–corrosion rates with experimental results will only be conducted using results from the model with the finest grid. In order to verify that the mesh used for verification was grid independent, a fourth mesh with 560k cells was constructed and the governing equations were solved. The maximum Sherwood number differed by only 7%. On this basis, the 368k grid was regarded as grid independent. In the contracted region, wear is under-predicted by the computational model. The point of highest wear is seen on the cusp of contraction. This is to be expected in the area of highest velocity both from the point of view increased oxygen transport across the boundary layer but also, in the case of a second solid phase, increased impacts on a sharp corner resulting in extra solid-particle erosion not considered by the current erosion–corrosion model. The results of the erosion prediction seem to confirm this. The erosion wear in this section is almost two orders of magnitude greater than in other sections. Detailed analysis of the erosion results shows that the wear is confined to the first element just on the cusp with no wear on any other element in the experimental section despite the very fine mesh in the area. A shortcoming of the experimental results is the following. The experimental results provide an average over the experimental section whereas the computational results can provide a much finer resolution, matching reality more closely. Reasonable matching of experimental and computational results are seen in the bulk of the constricted flow area. In future models where sharp corners are involved, the erosion model will be required to improve predictions. Downstream of the flow constriction where the flow once again expands, the wear rate is greatly under-predicted by the computa-
Position in test cell (m) 0
0.1
0.2
0.3
Experimental Data CFD erosion-corrosion prediction (368k cells) CFD erosion-corrosion prediction (60k cells) CFD erosion-corrosion prediction (8k cells) CFD erosion prediction
100
Weight loss rate (g m -2 h-1)
0.4
75 Sector
Sector
Sector
1
2
3
50
25
0
-6
-5
-4
-3
-2
-1
0
1
2
Position in test cell (diameters relative to expansion) Fig. 7. Comparision of CFD model with experimental data of Lotz and Postlethwaite [5].
3
4
C. Davis, P. Frawley / Corrosion Science 51 (2009) 769–775
tional model. The experimental area downstream of the expansion completely encompasses the recirculation region. As noted in Section 3.1.3, the RNG k e model cannot accurately model radial velocities, which are usually greatly under-predicted. In the recirculation region, the axial velocity drops to an almost negligible level, increasing the relative influence of the radial velocity. Since radial velocities are under-predicted by CFD, transfer of oxygen to the boundary layer will be reduced thereby decreasing the predicted wear rate. The greater influence of radial velocity in this region results in a greater under-prediction than in other areas. The erosion rate in this section is negligible. 5. Conclusions Overall, the erosion–corrosion model offers good matching between experimental and computational results as long as the mass and fluid boundary layers are fully resolved. The assumption that the solid-phase serves only to remove the erosion–corrosion product seems valid as the predicted erosion wear is generally negligible and cannot make up the shortfall in wear rate results. However, areas of flow constriction and detachment create problems for the model. In areas of flow constriction, velocities increase substantially and the assumption that the erosion component of the flow only removes the protective layer of corrosion products fails. In these cases, erosion wear becomes a significant factor. In the case of flow detachment, the restrictions of the RNG k e model, which lead to correctly predicted radial velocities and decrease the transport of oxygen to the reacting surface. As a result, erosion–corrosion wear is under-predicted. Alternative turbulence models, such as RSM or LES, may provide better prediction of these radial velocities and result in improved predictions. Due to the significant extra computational effort required for these models, it was not possible to incorporate this investigation into the existing study. A RSM investigation would require approximately three times the effort. A LES investigation would require a 3D grid. Even modelling a quarter of the section would require a minimum of twenty cells in the tangential direction, increasing the size of the model to 7.4 million cells. Consideration should be given to the data provided by Lotz and Postlethwaite [5]. While bulk velocity measurements were provided, neither mass flow nor density measurements were given. While data regarding oxygen saturation was taken from Lewis [31], data provided by Lotz and Postlethwaite [5] would have been preferable. The manner in which wear was calculated is effectively an average over the experimental specimen. A higher resolution of data across the specimen may not only show better agreement, but crucially highlight where the model fails in predicting wear. Finally, a tabulation of experimental data would have provided much improved accuracy over the graphical form of data provided by Lotz and Postlethwaite [5]. References [1] D.C. Silverman, Rotating cylinder electrode – geometry relationships for prediction of velocity-sensitive corrosion, Corrosion 44 (1) (1988) 42–49. [2] D.C. Silverman, Rotating cylinder electrode for velocity sensitivity training, Corrosion 40 (5) (1984) 220–226. [3] J. Postlethwaite, E.B. Tinker, M.W. Hawrylak, Erosion–corrosion in slurry pipelines, Corrosion 30 (8) (1974) 285–290.
775
[4] J. Postlethwaite, M.H. Dobbin, K. Bergevin, The role of oxygen mass transfer in the erosion–corrosion of slurry pipelines, Corrosion 42 (9) (1986) 514– 521. [5] U. Lotz, J. Postlethwaite, Erosion–corrosion in disturbed two phase liquid/ particle flow, Corrosion Science 30 (1990) 95–106. [6] J. Postlethwaite, U. Lotz, Mass transfer at erosion–corrosion roughened surfaces, The Canadian Journal of Chemical Engineering 66 (February) (1988) 75–78. [7] S. Nesic, J. Postlethwaite, Hydrodynamics of disturbed flow and erosion– corrosion. Part i – single-phase flow study, The Canadian Journal of Chemical Engineering 69 (June) (1991) 704–710. [8] S. Nesic, J. Postlethwaite, Hydrodynamics of disturbed flow and erosion– corrosion. Part ii – two-phase flow study, The Canadian Journal of Chemical Engineering 69 (June) (1991) 704–710. [9] A. Keating, S. Nesic, Prediction of two-phase erosion–corrosion in bends, in: Second International Conference on CFD in the Minerals and Process Industries, Melbourne, 1999, pp. 229–236. [10] U. Lotz, Velocity effects in flow induced corrosion, Corrosion (27) (1990) 22. [11] H.K. Versteeg, W. Malalasekera, An introduction to computational fluid dynamics the finite volume method, first ed., Longman Scientific and Technical, 1995. [12] S.V. Patankar, Numerical Heat Transfer and Fluid Flow, first ed., Taylor and Francis, 1980. [13] J.H. Ferziger, M. Peric, Computational Methods for Fluid Dynamics, third ed., Springer, 2002. [14] D.A. Drew, Mathematical modeling of two-phase flows, Annual Review of Fluid Mechanics 15 (1983) 261–291. [15] B.E. Launder, D.B. Spalding, Lectures in Mathematical Models of Turbulence, Academic Press, London, England, 1972. [16] V. Yakhot, S.A. Orszag, Renormalization group analysis of turbulence: I. Basic theory, Journal of Scientific Computing 1 (1986) 1–51. [17] D. Choudhury, Introduction to the renormalization group method and turbulence modeling, Technical report, Fluent Inc. Technical Memorandum TM-107, 1993. [18] Fluent Inc., Fluent User’s Guide Manual. Software Release Version 6.3, 2007. [19] A. Keating, A Model for the Investigation of Two-phase Erosion–Corrosion in Complex Geometries, Mesc, University of Queensland, 1999. [20] B. Kader, Temperature and concentration profiles in fully turbulent boundary layers, International Journal of Heat and Mass Transfer 24 (9) (1981) 1541– 1544. [21] M.M. Stack, N. Corlett, S. Turgoose, Some recent advances in the development of theoretical approaches for the construction of erosion–corrosion maps in aqueous conditions, Wear 233–235 (1999) 535–541. [22] L. Khezzar, J.H. Whitelaw, M. Yiannekis, An experimental study of round sudden expansion flows, in: Proc. 5th Symp. on Turbulent Shear Flows, Cornell University, Ithaca, NY, 1985, pp. 5–25. [23] R.P. Durret, W.H. Stevenson, H.D. Thompson, Radial and axial turbulent flow measurements with an ldv in an axisymmetric sudden expansion air flow, ASME Journal of Fluids and Engineering 110 (1988) 367–372. [24] M. Stieglmeier, C. Tropea, N. Weiser, W. Nitsche, Experimental investigation of the flow through axisymmetric expansions, ASME Journal of Fluids and Engineering 111 (1989) 464–471. [25] M.A. Founti, A. Klipfel, Experimental and computational investigations of nearly dense two-phase sudden expansion flows, Experimental Thermal and Fluid Science 17 (1998) 27–36. [26] F.M. White, Fluid Mechanics, fourth ed., McGraw-Hill, 1999. [27] F.M. White, Viscous Fluid Flow, second ed., McGraw-Hill, 1991. [28] A.P. O’Mahony, Numerical Simulation of Slurry Erosion using Lagrangian and Eulerian Techniques, Ph.d., University of Limerick, 2006. [29] V.G. Levich, Physicochemical Hydrodynamics, translated from the 2nd russian ed. edition., Prentice-Hall, 1962. [30] V. Gold, K.L. Loenino, A.D. McNaught, P. Sehmi, Compendium of Chemical Terminology: IUPAC Recommendations, Blackwell Scientific Publications, 1987. [31] M.E. Lewis, Dissolved Oxygen. U.S. Geological Survey Techniques of WaterResources Investigations, 2006. Book 9, Chapter A6, Section 6.2, Accessed 20th November, 2008, http://pubs.water.usgs.gov/twri9A6/. [32] W. Hayduk, B.S. Minhas, Correlations for prediction of molecular diffusivities in liquids, Canadian Journal of Chemical Engineering 60 (1982) 295–299. [33] Y. Oka, H. Ohnogi, T. Hosokawa, M. Matsumura, The impact angle dependence of erosion damage caused by solid particle impact, Wear 203–204 (1997) 573– 579. [34] A.P. O’Mahony, C.A. Marsh, A.J. Niven, P.J. Frawley, Determination of empirical erosion constants using a sand-blast type erosion rig, in: 7th International Symposium on Fluid Control, Measurement and Visualization, Sorrento, Italy, 2003.
Contents lists available at ScienceDirect
Corrosion Science journal homepage: www.elsevier.com/locate/corsci
Modelling of erosion–corrosion in practical geometries Cian Davis *, Patrick Frawley Department of Mechanical and Aeronautical Engineering, University of Limerick, Ireland
a r t i c l e
i n f o
Article history: Received 25 August 2008 Accepted 23 December 2008 Available online 17 January 2009 Keywords: B. CFD B. Wear prediction C. Erosion–corrosion C. Sherwood number
a b s t r a c t The limiting factor in certain instances of erosion–corrosion of steel is the presence of dissolved oxygen in the solution and the transfer of this oxygen to the reacting surface. Computational Fluid Dynamics (CFD) can be used to calculate oxygen diffusion throughout the flow and its transfer to the reacting surface. This was used in a computational model to calculate wear and validated against experimental results, for the first time, of erosion–corrosion wear in a contracting–expanding geometry. It was found that in order to correctly predict erosion–corrosion wear, Sherwood number independent grids were required providing a new metric to evaluate turbulent erosion–corrosion modelling. The predicted wear profile matched very closely with experimental results and overall matching was very good. Downstream of the flow expansion, erosion–corrosion wear was under-predicted. The disparity is due to detached flow for most of this sector where under-predicted radial velocities decrease transfer of oxygen to the reacting surface. This under-prediction is apparent in the downstream section due to the larger relative magnitude of the radial velocity in this sector. Ó 2009 Elsevier Ltd. All rights reserved.
1. Introduction Flowing fluids can have significant effects on rates of corrosion in a number of ways including accelerated mass transfer of reactants and corrosion products [1], fluid shear [2] or impingement of solid particles in the fluid which disrupt protective layers [3]. The ability to model and predict wear processes allows service intervals to be better timed so as to reduce unnecessary checks while not being subject to costly downtime due to equipment failure. Also, modelling ability can be applied at design stage to reduce the susceptibility of parts to wear. Due to the complexity of erosion–corrosion and significant number of inter-related variables, many attempts have been made to identify the controlling factors or relate the process to some other easily measurable quantity. Postlethwaite et al. [4] stated that erosion–corrosion of carbon steel piping in saline-sand slurry was controlled by mass transfer of oxygen to the reacting surface and that the role of the solid particles was to prevent formation of a complete rust film. Lotz and Postlethwaite [5] suggested this increased mass transport of oxygen was due to increased turbulent transfer across the mass-transfer boundary layer. Postlethwaite and Lotz [6] found a significant effect of surface roughness on oxygen mass-transfer in a similar sand-based slurry. Following on from Nesic and Postlethwaite [7,8], Keating and Nesic [9] adapted the work of Lotz [10] to formulate a model for erosion–corrosion. * Corresponding author. Tel.: +353 868045259; fax: +353 61213529. E-mail address: [email protected] (C. Davis). 0010-938X/$ - see front matter Ó 2009 Elsevier Ltd. All rights reserved. doi:10.1016/j.corsci.2008.12.025
This paper will apply the model of Keating and Nesic [9] to an contracting–expanding pipe geometry in order to prove the ability of the model in practical geometries. The contracting–expanding pipe geometry was chosen as detailed experimental data of both wear and flow parameters are available and complex flow characteristics are evident. Also, the optimum grid density will be determined for accurate predictions. 2. Modelling 2.1. Fluid modelling Modelling of the fluid is undertaken using Computational Fluid Dynamics (CFD). All fluid modelling was undertaken with FLUENT v6.3 code, a commercial CFD package. FLUENT uses a finite volume method in its CFD code. This is a well established and well utilised method, found in many CFD codes. The governing fluid equations are applied to a fixed zone in space, known as a control volume. These equations are also known as the conservation equations as they describe the conservation of mass, momentum and energy over the control volume. Detailed explanations of the model are provided by Versteeg and Malalasekera [11], Patankar [12] and Ferziger and Peric [13]. A derivation of the governing equations in an Eulerian formulation is give by Drew [14]. 2.1.1. Turbulence modelling The k e model was first proposed by Launder and Spalding [15] and has been subject of various improvements and modifications since. Yakhot and Orszag [16] developed an improvement to
770
C. Davis, P. Frawley / Corrosion Science 51 (2009) 769–775
the k e model by utilising renormalisation group (RNG) theory. In this RNG k e model, the constants that were derived empirically for the standard k e model are calculated using RNG theory. Also, additional terms and functions are added to the transport equations for k and e. By using this RNG-based approach, the model does not have to be modified to suit the type of flow being modeled and gives better predictions for separated flow or flows in curved geometries while only a small amount of extra computational effort is required – typically 10–15% [17]. 2.1.2. Wall effects There are two methods that can be utilised to resolve the near wall velocity profile numerically. The rigorous approach is to fully resolve the viscosity-affected region. This requires at least ten cells within the Rey < 200 region [18]. At high Reynolds numbers, the viscous sub-layer is so thin, it is difficult to resolve fully [13]. In a large geometry, addition of these extra grid points can impose a significant extra computational effort. This second approach is to utilise a wall function approach. This negates the need for excessive grid points at the boundary and is computationally less expensive. Wall functions are approximations based on simple flows and difficulties may arise in separated flows. Calculation of species transport brings with it the issue of the mass-transfer boundary layer, which is approximately an order of magnitude less than the velocity boundary layer [19]. The FLUENT CFD code includes a hybrid near-wall modelling method that utilises both approaches described above. In the case of a model with a fine mesh at the wall ðyþ 1Þ, the viscous sublayer is fully resolved. However, a model with a fine mesh throughout the model may not be desirable or practicable, for example at areas of negligible velocity or little interest. In these cases, (typically, yþ > 5) a coarse mesh can be applied and a blended wall function is used [20]. Therefore, areas of separation or interest can be easily investigated by increasing resolution of the mesh in that area, without the requirement of a fine mesh throughout the model, or even at all boundaries. This procedure provides a compromise between modelling error and computational effort and is employed to model near wall effects in this study. 2.2. Erosion–corrosion modelling Following on from Nesic and Postlethwaite [7,8], Keating and Nesic [9] adapted the work of Lotz [10] to formulate a model for erosion–corrosion. þþ
ð1Þ
The flux of the corrosive species (oxygen, in this case and denoted by J O2 ) is related to the mass transfer coefficient km , the wall concentration of oxygen, C wO2 and the bulk concentration of oxygen, C bO by 2
J O2 ¼ km ðC bO C wO2 Þ
ð2Þ
2
In the model C wO2 is set to zero as all the oxygen diffusing to the wall is consumed by the corrosion reaction. According to Eq. 1 two moles of Fe react with every mole of O2 . Therefore, the equation for the flux of Fe is given by
ð3Þ
2
2 1
The flux of Fe, given in kmol m s can be converted to a corrosion rate by use of the molar mass of Fe, M Fe (kg/kmol) and its density qFe (kg/m3). The result is an equation giving corrosion rate in mm/year.
CR ¼
2km C bO M Fe 2
qFe
J O2 ¼
D ðC fcO C wO2 Þ 2 Dy
ð5Þ
where Dy is the distance of the first cell centre from the wall, D is the mass diffusion coefficient (m2 s1) and C fcO is the concentration 2 at that point, obtained in CFD by solving the species transport equation. By assuming the concentration at the wall is zero and then combining Eqs. 2 and 5, km is given by as
km ¼
D C fcO2 Dy C b O
ð6Þ
2
Combining Eqs. 4 and 6 gives a erosion–corrosion model similar to that developed by Stack et al. [21].
CR ¼
2DC fcO M Fe 2
qFe Dy
24 60 60 365 103
ð7Þ
By using the assumption of Postlethwaite et al. [4] that the protective layer is removed by the combined action of the fluid and solid phase, the erosion–corrosion process can be assumed to be under the control of oxygen mass-transfer to the corroding surface. Both C fcO and Dy in Eq. 7 can be calculated with the aid of CFD. The 2 rate of erosion–corrosion wear can then be calculated using Eq. 7. The Sherwood number is a ratio of length-scale to diffusive boundary layer thickness and is defined as:
Sh ¼
km d D
ð8Þ
Eq. 4 can be written in terms of Sherwood number in a new correlation for erosion–corrosion wear
CR ¼
2DShC bO M Fe 2
qFe d
24 60 60 365 103
ð9Þ
The only flow-dependent variable in Eq. 9 is the Sherwood number itself. By definition, therefore, accurate erosion–corrosion predictions require a grid which is Sherwood number independent. 3. Model validation 3.1. Fluid model
2Fe ! 2Fe þ 4e O2 þ 2H2 O þ 4e ! 4OH
J Fe ¼ 2km C bO
The mass transfer coefficient km can be calculated in CFD in the following manner [19]. If the first cell from the wall is within the mass transfer boundary layer, the mass transfer between the two is governed by purely diffusive effects and Eq. 2 can be rewritten as
24 60 60 365 103
ð4Þ
Before CFD can be used to predict erosion–corrosion, flow characteristics must be validated for a flow regime similar to intended setup. Experimental investigations of turbulent flow through a sudden expansion have been conducted using laser Doppler anemometry (LDA) by numerous authors [22–24]. Difficulties in accurately measuring close to the wall were often reported. The near-wall layer is of particular interest in this erosion–corrosion model as it relies on calculating mass-transfer across the boundary layer. Founti and Klipfel [25] undertook a study in a similar geometry using diesel oil. Detailed experimental results of axial and radial velocity as well as turbulent measurements were provided at numerous locations downstream of the expansion. Results provided sufficient resolution in the near-wall region and detailed fluid parameters were provided. These results were selected for comparison with the numerical study. 3.1.1. Experimental description Founti and Klipfel conducted their experiments in a closed loop system with a vertically mounted working section. Upstream, the pipe has an internal diameter, d, of 25.5 mm and downstream,
C. Davis, P. Frawley / Corrosion Science 51 (2009) 769–775
771
where C l is given as 0.0845 from RNG theory and us is the friction velocity given by
us ¼
sw qc
0:5 ð13Þ
The relationship between wall shear stress and bulk inlet velocity is given as
sw ¼
q ub d 0:25 f ¼ 0:079 c
the diameter, d, was twice as large at 51 mm. The working section incorporated an expansion with a ratio of 2:1 and fully developed flow in the direction of gravity. A schematic diagram and indication of notation is given in Fig. 1. A mix of diesel oils was used in order that the refractive index of the working fluid matched the walls of the test section. Experimental parameters are given in Table 1. Measurements were taken using LDA and the error in results was given as 1% and 2% for the mean and turbulent measurements respectively. 3.1.2. Computational model In the computational modelling of the problem, steady, incompressible, isothermal, turbulent flow of diesel in a sudden expansion was assumed. Since the flow did not include a swirl component, the flow was modelled as 2D axisymmetric. This allowed significant resolution in areas of interest while maintaining a computationally efficient model. At the inlet to the computational domain, the axial velocity, turbulent kinetic energy and turbulent dissipation rate must be specified. In order to achieve the best approximation of experimental conditions, a fully turbulent flow was specified based on the axial peak velocity, up , at the centreline. The axial velocity at the inlet is related to the distance from the pipe centre line by a power law suggested by Prandtl [26] for fully developed turbulent pipe flow.
ð10Þ
Where d is the radius of the downstream pipe and r is the radial distance from the centreline. By integrating Eq. 10, an expression relating the bulk velocity, ub , to the peak velocity, up , can be obtained
up ¼ 1:224ub
ð11Þ
As an axisymmetric model was used, no tangential velocities are apparent in the model. A zero radial velocity was assumed at the inlet. For fully-developed flow, a linear variation from the wall to a free stream value at the centreline is assumed for the turbulent kinetic energy. The near-wall value of k is given by
u2s ffi knw ¼ pffiffiffiffiffi Cl
ð12Þ
Table 1 Experimental properties and operating conditions [25]. Model properties Density Viscosity Mass flow rate Upstream Reynolds number Downstream Reynolds number Bulk velocity (upstream of expansion) Bulk velocity (downstream of expansion)
ð14Þ
where f is the friction factor and is estimated for smooth pipes using the Blasius equation [27] as
Fig. 1. Diagram of experimental setup [25].
r 1=7 ux ¼ up 1 d
f qc u2b 2
830 kg m3 0.0043201 kg m1 s1 1.7304 kg s1 20,000 10,004 4.082 m s1 1.021 m s1
l
ð15Þ
Eq. 15 is valid for Reynold’s numbers between 4000 and 105 . Turbulent kinetic energy at the centreline us given by
kp ¼ 0:002u2p
ð16Þ
By combining Eqs. 12 and 16, a linear relationship for turbulent kinetic energy is given as
r k ¼ kp þ ðknw kp Þ d
ð17Þ
For fully-developed flow, the turbulent dissipation rate is given as 3=2
e¼
C 3=4 l k l
ð18Þ
where l is the mixing length and is set as jðd rÞ in the near-wall region and 0:09d in the outer layer [27]. By applying the experimental conditions given in Table 1 to Eqs. 10, 17 and 18, inlet profiles are generated for the models, shown in Fig. 2. In order to apply these profiles, a custom subroutine was used. 3.1.3. Fluid simulation results Overall, very good matching of experimental and computational data is obtained. A comparison of experimental and computational results for axial velocity between 4 mm and 50 mm downstream of the expansion are given in Fig. 3 and radial velocity in Fig. 4. Slight over-prediction of axial velocity is noticeable for values of r=0:5d of greater than 0.6 at x ¼ 75mm and x ¼ 100mm and for values of r=0:5d of greater than 0.2 at x ¼ 200mm. Overall, there is poor agreement between experimental and computational radial velocity results. This deviation can be accounted for by two main factors. First, the magnitude of the radial velocity is at least an order of magnitude less than the corresponding magnitude of the axial velocity. Second, the magnitude of the fluctuating component of the radial velocity is approximately double the magnitude of the mean radial velocity component. This results in a 200% turbulence intensity. Reynolds averaged turbulence models are unable to solve such high turbulence intensities [11]. Other studies have also reported difficulties in correctly predicting radial velocities in a similar flow setup [7,19,28]. Due to the small magnitude of the radial component relative to the axial component, the deficiencies in modelling the radial velocity do not have a significant effect on the overall fluid solution. 3.2. Erosion–corrosion model By using CFD to provide inputs for the equations derived above, an erosion–corrosion model was developed using solution variables that are currently solved by the CFD solver, such as oxygen concentration and cell wall distance. The computational model
772
C. Davis, P. Frawley / Corrosion Science 51 (2009) 769–775
Radial position (r/0.5 )
1
0.8
0.6
0.4
0.2
0 0
1
2
3
4
5
0
0.1
0.2
0.3
0
0.4
Turbulent kinetic energy (m 2/s2)
Axial Velocity (m/s)
1000 2000 3000 4000 5000 6000
Turbulent dissipation rate (m2/s3)
Fig. 2. Flow profiles at inlet.
1 x = 4mm
x = 6mm
x = 9mm
x = 50mm
Radial position (r/0.5d)
0.8
0.6
0.4
0.2
0 -1 0
1
2
3
4
5
6 -1 0
1
2
3
4
5
6 -1 0
1
2
3
4
5
6 -1 0
1
2
3
4
5
6
Axial velocity (m/s) Fig. 3. Comparison of experimental ðÞ and computational results (–) of axial velocity.
1 x = 9mm
x = 25mm
0.1 0.2 0.3 -0.3 -0.2 -0.1
0
x = 50mm
x = 75mm
Radial position (r/0.5d)
0.8
0.6
0.4
0.2
0 -0.2 -0.1
0
0.1
-0.2
-0.1
0
0.1
0.2 0 0.05 0.1 0.15 0.2 0.25 0.3
Radial velocity (m/s) Fig. 4. Comparison of experimental ðÞ and computational results (–) of radial velocity.
was validated by comparing the results with the experimental results of Lotz and Postlethwaite [5]. 3.2.1. Experimental description A saline solution with 2% concentration of sand was pumped through a vertically-mounted contracting–expanding nozzle made of carbon-steel experimental specimens (Fig. 5) with flow opposing gravity. The sections were separated with acrylic to provide electri-
cal insulation. This setup allowed investigation of erosion–corrosion in each section. The contracting–expanding geometry allowed investigation of the effects of detaching and reattaching flow. Experimental conditions are given in Table 2. 3.2.2. Model verification Verification of the model is especially important in the current context and is achieved by comparing the results against the
773
C. Davis, P. Frawley / Corrosion Science 51 (2009) 769–775 Table 3 Erosion–corrosion model validation grid sizes. Grid
Number of cells
Solution time (hours)
Full resolution of laminar sublayer Intermediate mesh Wall function mesh
368,900 60,870 8,349
52 15 4
Fig. 5. Example of test section in flow loop test rig [5].
experimental work of Lotz and Postlethwaite [5]. In the opinion of Keating [19], accurate calculation of wall oxygen concentrations in CFD required full resolution of the mass transfer boundary layer. Levich [29] developed an expression to approximate the width of the mass-transfer boundary layer relative to the fluid boundary layer
13 1 dmt 3 df Sc
ð19Þ
where dmt is the width of the mass-transfer boundary layer, df is the width of the fluid boundary layer and Sc is the Schmidt number, which relates momentum diffusivity to mass diffusivity and is given by
Sc ¼
m D
¼
l qD
ð20Þ
For the current system, the Schmidt number is 548, resulting in a mass-transfer boundary layer approximately 0.07 times the width of the fluid boundary layer. On this basis, full resolution of the mass transfer boundary layer would require yþ 0:1. Making a more precise comparison of the fluid and mass-transfer boundary layers is difficult, since the thickness of the masstransfer is arbitrary [30]. It is essentially another way to write the mass-transfer co-efficient. The best way to ensure proper resolution of the mass-transfer boundary layer is to verify the grid is Sherwood number independent. Resolution of the viscous sublayer is accomplished in CFD using so-called low Reynolds number (LRN) models. Applying a mesh with a yþ 0:1 throughout the computational domain, often results in an excessive computational requirement. In most cases, wall meshes are not this fine and the boundary layer is approximated with the use of wall functions. The RNG k e model offers enhanced wall functions, which cover the full range of modelling approaches [18]. For a value of Rey < 200, the laminar sublayer is fully resolved, where Rey is the turbulent Reynolds number given by:
pffiffiffi
qy k Rey ¼ l
ð21Þ
where y is the distance from the cell centre to the wall and k is the turbulent kinetic energy. The model offers accurate matching with small values of yþ and reasonable representation where yþ falls inside the wall buffer region ð3 < yþ < 10Þ.
In order to assess the performance of each boundary layer resolution method, three grids were constructed with varying nearwall mesh densities. The first mesh is designed to accurately resolve the laminar sublayer, while the third will use wall functions. The near-wall density of the second grid falls between these two method and is designed to evaluate model performance where optimum resolution of the laminar sublayer is not achieved. Details of grid size are given in Table 3. The three grids are shown in Fig. 6. 3.2.3. Computational setup Steady, incompressible, isothermal, turbulent flow of saline solution in a sudden contraction–expansion geometry was assumed. Since the flow did not include a swirl component, the flow was modelled as 2D axisymmetric. The same inlet profiles as detailed previously were used for k, e and axial velocity (Fig. 2). The experimental setup given by Lotz and Postlethwaite [5] specify that the water was saturated with oxygen. Lewis [31] provides tables of oxygen saturation in water for varying temperature, pressure and chloride concentration. At a standard pressure of 10 kPa and a temperature of 30 °C, the oxygen saturation is given as 7.4 mg/L. Correcting for a 3% saline concentration gives an oxygen saturation concentration of 6.23 mg/L. This model assumes that any oxygen reaching the wall is consumed by the corrosion reaction. On this basis, the concentration of oxygen at the wall is set to zero. The infinite dilution diffusion coefficient can be estimated using a model for solutes in aqueous solutions proposed by Hayduk and Minhas [32].
D ¼ 1:25 1012 V 0:19 0:292 T 1:52 l A where
is given as:
9:58 ¼ 1:12 VA
Model properties 340,000 170,000 13.2 m s1 3.3 m s1 2% by volume sand 430 lm Carbon steel (AISI MT-1015) 30 °C 3% by weight
ð23Þ
and V A is the solute molar volume at its normal boiling point (given as 25.6 cm3 mol1 for oxygen), T is temperature and l is the viscosity of water in centipoise. This results in a diffusivity coefficient of 1:83 109 m2 s1 for the system outlined. The wear rate in the experiments of Lotz and Postlethwaite [5] is calculated by measuring the change in weight of each experimental section over the experiment. There are twenty-seven experimental data points versus between 251 and 1904 data points in the computational model over the experimental section,
Table 2 Experimental properties and operating conditions given by Lotz and Postlethwaite [5].
Reynolds number (narrow section) Reynolds number (wide section) Bulk velocity (narrow section) Bulk velocity (wide section) Particle concentration Average particle diameter Experimental material Temperature NaCl concentration
ð22Þ
Fig. 6. Erosion–corrosion model validation grids.
774
C. Davis, P. Frawley / Corrosion Science 51 (2009) 769–775
depending on the grid size. In order to compare the data from the two sources, the corrosion rate for each computational data point is calculated using Eq. 7. This wear depth is used to calculate the volume of wear around the section in question. By multiplying this volume by the density, the mass moved can be calculated. As the width and radius of the computational section is known, the area can easily be calculated. By dividing the mass of material removed by the area of the section, the wear rate can be calculated over the same sections and in the same units as used by Lotz and Postlethwaite [5]
Wear rate ¼
pwððr CRÞ2 r2 Þq qðCR2 2rCRÞ ¼ 2r 2prw
ð24Þ
where w is the width of the experimental section, CR is the corrosion wear depth rate defined by Eq. 7, r is the pipe radius and q is the density of steel. In order to confirm that the contribution of the solid phase to wear is only to remove the corrosion product, a separate model was run to calculate solid-particle erosion wear. The discretephase model (DPM) in FLUENT tracks particles and impacts through the computational domain. This was combined with the erosion model developed by Oka et al. [33] and supported by experimental results from O’Mahony et al. [34] for the erosion profile of steel being impacted by sand. This gives an accurate computational prediction of solid-particle erosion. 4. Results and discussion The results of the computational study detailed in the previous section were reviewed in order to evaluate the accuracy of the erosion–corrosion model. Additionally, the results provided by each grid size were evaluated to quantify the effect of boundary layer resolution method on the prediction accuracy. The rate of erosion–corrosion was calculated as according to Eq. 24 for each mesh. The data are plotted against the experimental results of Lotz and Postlethwaite [5] in Fig. 7. The predicted erosion–corrosion rates are compared with experimental results over three sectors. The first three experimental sections make up the first sector up until the flow constriction. The second sector makes up the flow constriction consisting of the next eleven experimental sections. The third sector is the flow expansion after the constriction consisting of the final thirteen
experimental sections. They are numbered as sectors 1, 2 and 3, respectively in Fig. 7 to aid identification. It is immediately apparent from Fig. 7 that the wear predicted by the models with the two coarsest grids (60k cells and 9k cells) is greatly under-predicted. Wear in these two models is approximately an order of magnitude lower than both experimental results and results from the model with the finest grid. It is apparent, therefore, that full resolution of the mass transfer boundary layer is required for prediction of erosion–corrosion wear using the current model. Neither wall functions nor partial resolution of either the mass transfer or fluid boundary layers correctly predict transfer of oxygen to the reacting surface. On this basis, subsequent comparison of predicted erosion–corrosion rates with experimental results will only be conducted using results from the model with the finest grid. In order to verify that the mesh used for verification was grid independent, a fourth mesh with 560k cells was constructed and the governing equations were solved. The maximum Sherwood number differed by only 7%. On this basis, the 368k grid was regarded as grid independent. In the contracted region, wear is under-predicted by the computational model. The point of highest wear is seen on the cusp of contraction. This is to be expected in the area of highest velocity both from the point of view increased oxygen transport across the boundary layer but also, in the case of a second solid phase, increased impacts on a sharp corner resulting in extra solid-particle erosion not considered by the current erosion–corrosion model. The results of the erosion prediction seem to confirm this. The erosion wear in this section is almost two orders of magnitude greater than in other sections. Detailed analysis of the erosion results shows that the wear is confined to the first element just on the cusp with no wear on any other element in the experimental section despite the very fine mesh in the area. A shortcoming of the experimental results is the following. The experimental results provide an average over the experimental section whereas the computational results can provide a much finer resolution, matching reality more closely. Reasonable matching of experimental and computational results are seen in the bulk of the constricted flow area. In future models where sharp corners are involved, the erosion model will be required to improve predictions. Downstream of the flow constriction where the flow once again expands, the wear rate is greatly under-predicted by the computa-
Position in test cell (m) 0
0.1
0.2
0.3
Experimental Data CFD erosion-corrosion prediction (368k cells) CFD erosion-corrosion prediction (60k cells) CFD erosion-corrosion prediction (8k cells) CFD erosion prediction
100
Weight loss rate (g m -2 h-1)
0.4
75 Sector
Sector
Sector
1
2
3
50
25
0
-6
-5
-4
-3
-2
-1
0
1
2
Position in test cell (diameters relative to expansion) Fig. 7. Comparision of CFD model with experimental data of Lotz and Postlethwaite [5].
3
4
C. Davis, P. Frawley / Corrosion Science 51 (2009) 769–775
tional model. The experimental area downstream of the expansion completely encompasses the recirculation region. As noted in Section 3.1.3, the RNG k e model cannot accurately model radial velocities, which are usually greatly under-predicted. In the recirculation region, the axial velocity drops to an almost negligible level, increasing the relative influence of the radial velocity. Since radial velocities are under-predicted by CFD, transfer of oxygen to the boundary layer will be reduced thereby decreasing the predicted wear rate. The greater influence of radial velocity in this region results in a greater under-prediction than in other areas. The erosion rate in this section is negligible. 5. Conclusions Overall, the erosion–corrosion model offers good matching between experimental and computational results as long as the mass and fluid boundary layers are fully resolved. The assumption that the solid-phase serves only to remove the erosion–corrosion product seems valid as the predicted erosion wear is generally negligible and cannot make up the shortfall in wear rate results. However, areas of flow constriction and detachment create problems for the model. In areas of flow constriction, velocities increase substantially and the assumption that the erosion component of the flow only removes the protective layer of corrosion products fails. In these cases, erosion wear becomes a significant factor. In the case of flow detachment, the restrictions of the RNG k e model, which lead to correctly predicted radial velocities and decrease the transport of oxygen to the reacting surface. As a result, erosion–corrosion wear is under-predicted. Alternative turbulence models, such as RSM or LES, may provide better prediction of these radial velocities and result in improved predictions. Due to the significant extra computational effort required for these models, it was not possible to incorporate this investigation into the existing study. A RSM investigation would require approximately three times the effort. A LES investigation would require a 3D grid. Even modelling a quarter of the section would require a minimum of twenty cells in the tangential direction, increasing the size of the model to 7.4 million cells. Consideration should be given to the data provided by Lotz and Postlethwaite [5]. While bulk velocity measurements were provided, neither mass flow nor density measurements were given. While data regarding oxygen saturation was taken from Lewis [31], data provided by Lotz and Postlethwaite [5] would have been preferable. The manner in which wear was calculated is effectively an average over the experimental specimen. A higher resolution of data across the specimen may not only show better agreement, but crucially highlight where the model fails in predicting wear. Finally, a tabulation of experimental data would have provided much improved accuracy over the graphical form of data provided by Lotz and Postlethwaite [5]. References [1] D.C. Silverman, Rotating cylinder electrode – geometry relationships for prediction of velocity-sensitive corrosion, Corrosion 44 (1) (1988) 42–49. [2] D.C. Silverman, Rotating cylinder electrode for velocity sensitivity training, Corrosion 40 (5) (1984) 220–226. [3] J. Postlethwaite, E.B. Tinker, M.W. Hawrylak, Erosion–corrosion in slurry pipelines, Corrosion 30 (8) (1974) 285–290.
775
[4] J. Postlethwaite, M.H. Dobbin, K. Bergevin, The role of oxygen mass transfer in the erosion–corrosion of slurry pipelines, Corrosion 42 (9) (1986) 514– 521. [5] U. Lotz, J. Postlethwaite, Erosion–corrosion in disturbed two phase liquid/ particle flow, Corrosion Science 30 (1990) 95–106. [6] J. Postlethwaite, U. Lotz, Mass transfer at erosion–corrosion roughened surfaces, The Canadian Journal of Chemical Engineering 66 (February) (1988) 75–78. [7] S. Nesic, J. Postlethwaite, Hydrodynamics of disturbed flow and erosion– corrosion. Part i – single-phase flow study, The Canadian Journal of Chemical Engineering 69 (June) (1991) 704–710. [8] S. Nesic, J. Postlethwaite, Hydrodynamics of disturbed flow and erosion– corrosion. Part ii – two-phase flow study, The Canadian Journal of Chemical Engineering 69 (June) (1991) 704–710. [9] A. Keating, S. Nesic, Prediction of two-phase erosion–corrosion in bends, in: Second International Conference on CFD in the Minerals and Process Industries, Melbourne, 1999, pp. 229–236. [10] U. Lotz, Velocity effects in flow induced corrosion, Corrosion (27) (1990) 22. [11] H.K. Versteeg, W. Malalasekera, An introduction to computational fluid dynamics the finite volume method, first ed., Longman Scientific and Technical, 1995. [12] S.V. Patankar, Numerical Heat Transfer and Fluid Flow, first ed., Taylor and Francis, 1980. [13] J.H. Ferziger, M. Peric, Computational Methods for Fluid Dynamics, third ed., Springer, 2002. [14] D.A. Drew, Mathematical modeling of two-phase flows, Annual Review of Fluid Mechanics 15 (1983) 261–291. [15] B.E. Launder, D.B. Spalding, Lectures in Mathematical Models of Turbulence, Academic Press, London, England, 1972. [16] V. Yakhot, S.A. Orszag, Renormalization group analysis of turbulence: I. Basic theory, Journal of Scientific Computing 1 (1986) 1–51. [17] D. Choudhury, Introduction to the renormalization group method and turbulence modeling, Technical report, Fluent Inc. Technical Memorandum TM-107, 1993. [18] Fluent Inc., Fluent User’s Guide Manual. Software Release Version 6.3, 2007. [19] A. Keating, A Model for the Investigation of Two-phase Erosion–Corrosion in Complex Geometries, Mesc, University of Queensland, 1999. [20] B. Kader, Temperature and concentration profiles in fully turbulent boundary layers, International Journal of Heat and Mass Transfer 24 (9) (1981) 1541– 1544. [21] M.M. Stack, N. Corlett, S. Turgoose, Some recent advances in the development of theoretical approaches for the construction of erosion–corrosion maps in aqueous conditions, Wear 233–235 (1999) 535–541. [22] L. Khezzar, J.H. Whitelaw, M. Yiannekis, An experimental study of round sudden expansion flows, in: Proc. 5th Symp. on Turbulent Shear Flows, Cornell University, Ithaca, NY, 1985, pp. 5–25. [23] R.P. Durret, W.H. Stevenson, H.D. Thompson, Radial and axial turbulent flow measurements with an ldv in an axisymmetric sudden expansion air flow, ASME Journal of Fluids and Engineering 110 (1988) 367–372. [24] M. Stieglmeier, C. Tropea, N. Weiser, W. Nitsche, Experimental investigation of the flow through axisymmetric expansions, ASME Journal of Fluids and Engineering 111 (1989) 464–471. [25] M.A. Founti, A. Klipfel, Experimental and computational investigations of nearly dense two-phase sudden expansion flows, Experimental Thermal and Fluid Science 17 (1998) 27–36. [26] F.M. White, Fluid Mechanics, fourth ed., McGraw-Hill, 1999. [27] F.M. White, Viscous Fluid Flow, second ed., McGraw-Hill, 1991. [28] A.P. O’Mahony, Numerical Simulation of Slurry Erosion using Lagrangian and Eulerian Techniques, Ph.d., University of Limerick, 2006. [29] V.G. Levich, Physicochemical Hydrodynamics, translated from the 2nd russian ed. edition., Prentice-Hall, 1962. [30] V. Gold, K.L. Loenino, A.D. McNaught, P. Sehmi, Compendium of Chemical Terminology: IUPAC Recommendations, Blackwell Scientific Publications, 1987. [31] M.E. Lewis, Dissolved Oxygen. U.S. Geological Survey Techniques of WaterResources Investigations, 2006. Book 9, Chapter A6, Section 6.2, Accessed 20th November, 2008, http://pubs.water.usgs.gov/twri9A6/. [32] W. Hayduk, B.S. Minhas, Correlations for prediction of molecular diffusivities in liquids, Canadian Journal of Chemical Engineering 60 (1982) 295–299. [33] Y. Oka, H. Ohnogi, T. Hosokawa, M. Matsumura, The impact angle dependence of erosion damage caused by solid particle impact, Wear 203–204 (1997) 573– 579. [34] A.P. O’Mahony, C.A. Marsh, A.J. Niven, P.J. Frawley, Determination of empirical erosion constants using a sand-blast type erosion rig, in: 7th International Symposium on Fluid Control, Measurement and Visualization, Sorrento, Italy, 2003.