- Email: [email protected]
A calculation methodology proposed for liquid droplet impingement erosion
Nuclear Engineering and Design 242 (2012) 157–163
Contents lists available at SciVerse ScienceDirect
Nuclear Engineering and Design journal homepage: www.elsevier.com/locate/nucengdes
A calculation methodology proposed for liquid droplet impingement erosion Rui Li a,∗ , Michitsugu Mori b,c , Hisashi Ninokata a a
Research Laboratory for Nuclear Reactors, Tokyo Institute of Technology, 2-12-1-N1-5, Ookayama, Meguro-ku, Tokyo 152-8550, Japan Research and Development Centre, Tokyo Electric Power Company, 4-1, Egasaki-cho, Tsurumi-ku, Kanagawa 230-8510, Japan c School of Science and Technology, Meiji University, 1-1-1, Higashi-Mita, Tama-ku, Kawasaki-shi, Kanagawa 214-8571, Japan b
a r t i c l e
i n f o
Article history: Received 6 May 2011 Received in revised form 3 October 2011 Accepted 5 October 2011
a b s t r a c t Bent pipe wall thinning has been often found at the elbow of the drain line and the high-pressure secondary feed-water bent pipe in nuclear reactors. Liquid droplet impingement (LDI) erosion could be regarded as one of the major causes and is a significant issue of the thermal hydraulics and structural integrity in aging and life extension for nuclear power plant safety. In this paper a computational methodology is established for simulation of LDI erosion using computational fluid dynamics (CFD) simulation and theoretical calculation. Two-phase flow numerical simulations are conducted for standard elbow geometry, typically with the pipe diameter of 170 mm. This computational fluid model is built up by incompressible Reynolds Averaged Navier–Stoke equations using standard k–ε turbulence model and the SIMPLE algorithm, and the numerical droplet model adopts the Lagrangian approach. The turbulence damping in vapor–droplets flow is theoretically analyzed by a damping function on the energy spectrum basis of single phase flow. Locally, a droplet impact angle function is employed to determine the overall erosion rate. Finally, the overall and local investigations are combined to purpose a general methodology of LDI erosion prediction procedure, which has been complemented into CFD code. Based on our more physical computational results, comparison with an available accident data was made to prove that our methodology could be an appropriate way to simulate and predict the bent pipe wall thinning phenomena. © 2011 Elsevier B.V. All rights reserved.
1. Introduction In the area of safety analysis for power engineering industry, the bent pipe wall thinning has been often found in nuclear power plants as well as the fossil power plants. For example, JSME annual report records one accident that the rupture pipe is located as Fig. 1 shows, the small hole is also depicted (JSME, 2007). Liquid droplet impingement (LDI) erosion is one of the major causes thus evaluating LDI erosion is an important issue of the thermal hydraulics and structural integrity in aging and life extension for nuclear power plants safety. In general, major parameters that influence LDI erosion are categorized into four groups: flow pattern, chemical conditions, droplet impingement properties and material properties (Okada and Uchida, 2009). The droplet erosion mechanism is not similar with particulate solid sand erosion (Springer, 1976). From the point of flow pattern view, erosion rate is dependent on a number of factors including the droplet size, impact velocity,
∗ Corresponding author. Tel.: +81 3 5734 3062; fax: +81 3 5734 3056. E-mail address: [email protected] (R. Li). 0029-5493/$ – see front matter © 2011 Elsevier B.V. All rights reserved. doi:10.1016/j.nucengdes.2011.10.004
impact frequency, and liquid and gas density and viscosity. As many of these values are unknown for field situations, it is not easy to predict the rate of droplet erosion. It should also be noted that control of many of these factors in laboratory-based tests is problematical, for instance, the droplet velocity, is difficult to be controlled and measured. Therefore, a great deal of care should be required when extrapolating (from low velocity to high velocity) laboratory test results to practical conditions. Hence computational fluid dynamics (CFD) can be used as a valid tool to investigate the LDI phenomena. There are few equations used to predict liquid droplet erosion, although there are many equations developed to predict solid particle erosion in literatures (Tabakoff and Kotwal, 1979; Oka et al., 2005; Oka and Yoshida, 2005). The erosion equations are generally based on particle impact speed and angle, particle and material properties. However, in case of flow in a bent pipe with 90 degree, various flow situations and fluid properties influence particle impact speed and angle. Therefore CFD can be used to predict particle impact speed and angle which can be employed in erosion models to predict erosion. Our model to predict erosion is based on solid particle Edwards (2000) erosion model, and then developed with our considerations of liquid droplets compressibility and
158
R. Li et al. / Nuclear Engineering and Design 242 (2012) 157–163
Nomenclature a b Ad An Ai A an c C Cd dd d (v) E () F f (˛) g H() k l m M p t Re R R r U ud
v Vd W x y Greek ˛ ε
constant in turbulence theory constant in turbulence theory droplet cross-sectional area (m2 ) amplitude in velocity spectrum of fluid (m) constant wall face area (m2 ) amplitude in velocity spectrum of droplets (m) a constant in turbulence theory erosion formulation constant drag force coefficient droplet diameter (m) droplet impact velocity function spectrum of turbulent kinetic energy (m3 s−2 ) force (Pa) impingement angle function gravity acceleration (m s−2 ) damping function of carrier turbulence due to droplets turbulent kinetic energy (m2 s−2 ) macroscopic length scale of turbulence (m) droplets mass (kg) momentum loss or gain (kg m s−1 ) pressure (Pa) time (s) Reynolds number normalized relative velocity between two phases (m s−1 ) erosion rate (kg m−2 s−1 ) radius (m) mean velocity (m s−1 ) droplet velocity (m s−1 ) velocity (m s−1 ) droplet volume (m3 ) mass concentration of droplets coordinate direction (m) coordinate direction (m)
ıtd
volume fraction, droplet impingement angle dissipation rate of turbulent kinetic energy (m3 s−2 ) Kolmogorov microscale of turbulence (m) dynamic viscosity (kg m−1 s−1 ) wave number (m−1 ) Kolmogorov constant density (kg m−3 ) time step (s)
Indices am c d dr n o tp
virtual mass force carrier fluid droplets drag force new value old value two phase
deformation. In this paper, we establish an innovative methodology of LDI erosion rate for evaluating wall thinning phenomena. We proposed a calculation methodology to provide useful information for industry applications. The current paper originalities are (1) applying a droplet–vapor two-way coupling calculation system and a vapor turbulence damping analysis to LDI phenomena, (2)
implementing an impingement angle function into codes for LDI erosion rate calculation in a full scale bent pipe, and (3) making a comparison with an actual accident erosion data to validate our proposed methodology. 2. Computational models and methods 2.1. Outline In this section we introduce our calculation methodology including two-way vapor–droplet coupling computational method, calculation with carrier turbulence damping theory calculation, local impact force function, erosion semi-empirical formulation and finally the procedure to implement codes using a CFD tool and a set of subroutines. The standard elbow bent pipe computational domain and mesh are also given. 2.2. Two way coupling vapor–droplet method Two-phase flow characteristics could be predicted by CFD, a momentum transfer consideration called two-way coupling between carrier fluid and transported droplets is taken into account. The droplet governing equation is expressed in Lagrangian description: (und − uod )
= Fndr + Fnp + Fnam + Fng , ıtd 1 Fdr = − Cd c Ad |ud − u|(ud − u), 2 Fp = −Vd ∇ p, d(ud − u) 1 Fam = − c Vd , 2 dt Fg = md g, md
(1)
where md is the droplet mass; u and ud are the continuous phase velocity and droplet velocity, respectively; n and o denote “new” and “old” values spanning the integration time step ıtd , respectively. In the droplet Lagrangian description, Fdr , Fp , Fam and Fg are the drag force, pressure force, virtual mass force and the gravity force, respectively. ∇ p is the pressure gradient in the carrier fluid, Vd is the droplet volume, Ad is the droplet cross-sectional area, Cd is the drag coefficient, and g is the gravitational acceleration. The lift force is not considered due to the fact that the typical particle diameter in LDI phenomenon is around 10 m rather than submicro meter (Saffman, 1965). As for the carrier gas, its momentum equation is: ∂ui ∂ ∂p + (uj ui − ij ) = − , ∂t ∂xj ∂xi
(2)
where ui is the carrier fluid velocity component in direction xi and the term ij represents the stress tensor components and need closure. In this current study we use the standard k–ε model which has become the workhorse of practical engineering flow calculations for pipe flows (Launder and Spalding, 1972). Here this calculation system is called “one way coupling” system. In two way coupling system, the momentum exchanges between carrier fluid and droplet are obtained from the integration of the Lagrangian equations. The vapor–droplet coupled model accounts for the fact that each droplet may exchange momentum with the carrier fluid, and that a loss from a droplet is a gain for the fluid and vice versa. From the average rate of momentum loss or gain of a single droplet due to drag and other effects over a single time increment ıtd , the total momentum loss or gain of a group of droplets traversing a given cell is expressed as follows: Md =
all droplets
md
(und − uod ) ıtd
− Fnp
− Fng
,
(3)
R. Li et al. / Nuclear Engineering and Design 242 (2012) 157–163
159
Fig. 1. The small hole as a result of pipe thinning and the schematic drawing (JSME, 2007).
which is inserted as a “source” in the momentum equation of the carrier phase. So the coupling between the fluid and the droplet is achieved via the sink term Md in the momentum equation (2) of the carrier fluid: ∂ ∂p ∂ui + (uj ui − ij ) = − + Md , ∂t ∂xj ∂xi
(4)
2.3. Damps of turbulent kinetic energy due to small droplets Since the turbulence of carrier fluid might be attenuated due to the involved droplets, the impinged droplets erosion rate might be affected by the turbulence of carrier fluid. It is applicable to add the turbulence kinetic energy damping of carrier fluid consideration for the LDI rate evaluation (Li and Ninokata, 2011). The turbulence of carrier fluid contains a wide range of eddy sizes. The strain rate of large or “energy-containing” eddies is proportional to the average strain rate of the mean flow. The strain rate of small eddies is larger compared to that of the mean flow and of the large eddies, so that no permanent anisotropy can be induced at small scales. The spectrum of turbulent kinetic energy, expressed as a function of the wave number , could consist of four sub ranges which are defined as follows (Tennekes, 1972): (1) Larger scale subrange :
(2) Transition subrange :
−3/4
E0 () = cε
c c
E1 () = c(a2 + b),
17/4 4
,
0 ≤ ≤ 0 = 10,
0 ≤ ≤ l =
(5)
to compute the damping of turbulent kinetic energy of the carrier due to the droplets:
H() = exp −
(6)
R2 −5/3 d ,
(9)
0
where dd is the droplet diameter, d is droplet density, W is the mass concentration of the droplets, and R is the normalized relative velocity between carrier and droplets. They are given as: ˛d d , ˛d d + ˛c c
W=
R =
1+
a 2 n
An
an −2 · An
(10)
1/2
1
1 + Uc2 2 (dd d /18c )
,
(11)
where An is the amplitude in velocity spectrum of fluid, an is the amplitude in velocity spectrum of droplets, and Uc is the carrier mean velocity. Hence when the vapor flow with a dilute droplet phase is simulated, the damping of turbulent kinetic energy of the carrier fluid due to the droplets should be computed. Then the revised kinetic energy values are used for erosion evaluation. The formalism is applied through the following procedure: Step 1 (:). The standard k–ε model is utilized to the carrier to compute the integral value of the single pure phase turbulent kinetic energy neglecting the effect of droplets:
1 , l
36Wc ε−1/3 dd d
∞
k=
E(v)dv,
(12)
0
(3) Inertial subrange :
E2 () = cε2/3 −5/3 ,
(4) Equilibrium subrange :
=
1 ≤ ,
l =
E3 () = cε2/3 −5/3
1 1 ≤ ≤ = , l
(7)
3 4/3 exp − () , 2
Step 2 (:). The theoretical distributions of the single pure phase turbulent spectra are computed by self-developed original program, as Section 2.1 discussed in four intervals. E(v) (0 ≤ v ≤ ∞),
(8)
The above four sub ranges are for the single pure phase turbulent energy spectrum. As it was pointed out that the dispersed phase may have an influence to the continuous phase, the influence is considered: the dispersion of droplets due to the continuous fluid phase and, vice versa, the impact of the droplets upon the turbulence intensity of the fluid. According the theoretical analysis, the suppression of turbulent kinetic energy occurs in case of small particles. Bottoni and Ninokata (1998) proposed a damping function
(13)
Step 3 (:). The damping function of carrier turbulence H() is calculated by Eqs. (9)–(11), and then the two-phase flow spectra of turbulent kinetic energy are computed by: E(v)tp = E(v)H(v),
(14)
Step 4 (:). The two-phase flow spectra of turbulent kinetic energy are integrated numerically over the full range of wave numbers to obtain the integral value:
ktp =
∞
E(v)tp dv. 0
(15)
160
R. Li et al. / Nuclear Engineering and Design 242 (2012) 157–163
Fig. 2. Schematic diagrams and bent section mesh of a standard elbow.
The value ktp should replace the single phase value k cell by cell in the computational domain of the carrier fluid.
be done. The solid particle (sand) impingement angle function is (Okita and Shirazi, 2010):
Step 5 (:). The standard k–ε model was applied again with the revised kinetic energy ktp using two-way coupling carrier-droplets computational system to evaluate the LDI erosion rate which is our interest and final target.
f (˛) = 5.4˛ − 10.1˛2 + 10.9˛3 − 6.3˛4 + 1.42˛5 ,
2.4. Overall erosion model In the erosion prediction model, the wall mass loss and penetration rates are calculated. When a droplet impinges the wall, the mass loss distributes uniformly over the computational cell of the wall. With both sufficiently small grid spacing and a large number of droplets, approximation errors induced by these assumptions still cannot be eliminated. Once all droplet trajectories have been computed and all wall impingement data is gathered, the total mass loss for all impingements can be compiled to generate a local penetration rate for each cell which lies on the surface of the geometry. Through the subroutines, the droplet impingement information is gathered. The computed erosion rate is the removed mass per area and time in the unit kg m−2 s−1 , which is defined by (Edwards, 2000)
m C(d )f (˛)vd(v) d d
0◦ < ˛ ≤ 90◦ , (17)
Unlike solid particle erosion the scratch stress caused by tangential momentum might be negligibly small, in case of LDI erosion. Therefore, the scratch stress is much less important due to the liquid mobility and deformation. The normal component of impact velocity is major factor for the LDI erosion damage whilst the tangential component has little effect. The effect of droplet impingement angle on LDI erosion was investigated by the authors (Li and Ninokata, 2011) using the volume of fluid (VOF) computational system, which is a two-phase Eulerian–Eulerian approach. Based on our previous conclusion the impingement angle function f(˛) can be fitted as an approximation shown in Eq. (18), which is also implemented into codes to support the calculation of the overall erosion rate. f (˛) = 1.3 × 10−1 + 1.5 × 10−2 ˛ − 6.4 × 10−5 ˛2 ,
0◦ < ˛ ≤ 90◦ . (18)
Nimpinged
Rerosion =
n=1
Aface
,
(16)
where N is the number of droplets which impinge onto the wall, Aface is the area of the cell face at the wall. C (dd ) is a function of the droplet diameter function and impinged wall density (usually stainless steel), 1.18 × 10−11 without unit, v is the droplet velocity before impinging the wall. d(v) is the droplet impact velocity function, 1.73 without unit. C(dd ) and d(v) are determined from the target wall surface materials. f(˛) is the impingement angle function, ˛ is droplet impingement angle with the wall face. 2.5. Impingement angle function LDI erosion usually occurs in the inner surface of a bent pipe, with the droplets carried by the carrier fluid impinging on the surface with different impingement angles. In the Edward’s model the sand particles are assumed as impinged materials. However, in case of liquid droplets due to their compressibility and deformation after impingements, some reasonable modification should
2.6. Computational domain of a standard bent pipe The photograph of a full scale bent pipe with the flow direction is shown in Fig. 1 where the length unit is mm. It can be seen in Fig. 2, in order to use the uniform inlet condition and fully developed outflow condition, the computational domain is extended in both the upstream and downstream parts, and the two extended parts are called the pre-extended and after-extended region, respectively. In our numerical simulations, the pre-extended region is in length of 20 pipe diameters and the after-extended region is up in length of 10 pipe diameters. The mesh system is also shown in Fig. 2, where the first mesh adjacent to inner wall surface is also shown. In the boundary layer region we employ a much denser grid, whilst in the central region a coarse grid is adopted to save the computer resource. Through a network of subroutine to implement the two-way coupling method and erosion formulation including the impingement angle function, the governing equations along with the
R. Li et al. / Nuclear Engineering and Design 242 (2012) 157–163
161
Table 1 Main assumptions for the discrete phase model. Parameter (unit)
Value
Droplet density (kg m−3 ) Uniform diameter (m) Droplet number density (m−3 ) Total mass flow rate (kg s−1 ) Droplet volume fraction Droplet velocity (m s−1 )
998 10 8.78 × 1012 4.59 0.0001 400
3. Results 3.1. Discussion
Fig. 3. The vapor inlet velocity power-law distribution.
boundary conditions are solved using the commercial finitevolume based solver FLUENT 6.3. The SIMPLE algorithm is adopted to deal with the linkage between velocity and pressure. Following assumptions are made in the numerical simulation: (1) the fluid flow is in steady-state; (2) the effect of the gravitational force is negligible; (3) the fluid thermo-physical properties are constant; and (4) the Mach number is less than 0.3, that is the compressibility effect can be the carrier fluid is incompressible.
2.7. Boundary conditions and initial conditions The carrier fluid is vapor. The required conditions are described for the three different parts, which are inlet, wall and outlet, respectively. Following descriptions are the boundary conditions. (1) The implementation of pipe wall boundary is treated as the normal no-slip wall and the wall surface is supposed to be smooth. (2) The vapor density is 0.554 kg m−3 at room pressure (1.01 × 105 Pa) and room temperature (298 K). (3) The outlet is prescribed as outflow. The liquid water density is 998.0 kg m−3 , which is around 1800 times larger than vapor, which means that this large density ratio will bring strong effects on the droplets trajectories towards to impinge the wall. The vapor is the fully developed turbulent flow in the bent pipe, where the vapor velocity inlet profile usually presents a power-law distribution (Sugiyama, 2003). For instance if the inlet vapor mean velocity is 280 m s−1 , the corresponding distribution is shown in Fig. 3, where r is the radius from wall to the center of the pipe. In our calculations we assume that the initial droplet velocity has the same value as the initial mean vapor velocity, and in this paper we briefly refer them as droplet velocity. The impact of droplets on rigid surfaces may produce a wide variety of consequences, according to the size, velocity and material of the impacting elements and the nature of the surface. For example, droplets may adhere, bounce or shatter, and the liquid deposited from the wall surface may retain its droplet form or merge into a liquid film. In our numerical simulations, the treatment of wall and droplet interaction is simplified by the perfect rebound function. Due to the high Reynolds number, we use the standard wall functions. In the practical operating conditions, the droplet parcel should consist of various droplets with different diameters; in this research we define uniform diameter distribution. The droplet number density is assumed to be uniformly at the inlet.
From the above sections three major considerations are taken into account for LDI erosion prediction: (1) two-way coupling, (2) carrier turbulence damping, and (3) droplet impingement angle function. The main target of this current work is the droplet erosion, which is highly dependent on droplet velocity profile as shown in Eq. (16). Since this study is from the point view of fluid dynamics, it is important to discuss the details of each consideration and compare with those results without the above considerations. Therefore we can know that three considerations are absolutely necessary, thus it is approved that our methodology might be correct to predict the LDI erosion. It has to be confirmed that the calculated droplet velocity profiles are more physical using the models of the two-way coupling modeling and turbulence damping function. The erosion prediction with droplet impact angle function also should be made a comparison with the erosion prediction with solid particle impact angle function. Some discussions are addressed in this section. The main assumptions for the discrete droplets are listed in Table 1, which are typical values in LDI phenomena. In our calculations we assume that the initial droplet velocity has the same value as the initial mean vapor velocity, we briefly refer them as droplet velocity. 3.2. Two way coupling In order not to be influenced by the elbow geometry, the carrier and droplet normalized velocity profile are post-processed in a section ahead located in the length of 5 pipe diameters before the elbow. Figs. 4 and 5 show that the vapor and droplet velocity profiles calculated from one way coupling system and two way coupling system, respectively. Please note that unlike the vapor
Fig. 4. Carrier velocity profile using one way coupling and two way coupling system.
162
R. Li et al. / Nuclear Engineering and Design 242 (2012) 157–163
Fig. 5. Droplet velocity profile using one way coupling and two way coupling system.
Fig. 7. Carrier turbulence kinetic energy profiles between the results with damping function and without damping function.
velocity profile can be extracted directly from the calculation mesh value, the droplet velocity profile was a statistic value based on the sampling of Lagrangian tracking, so there might be some discontinues in the distributions (ANSYS, 2002). Due to the huge amount droplets in the carrier vapor, two way coupling calculation carrier velocity profile becomes ‘more laminar’ than single phase flow. As a result the droplet velocity profile also became “more sharp” with higher velocity in the core. Based on the erosion calculation equation (Eq. (16)), we can make a prediction that the one way coupling erosion result should be lower than two way coupling results. Obviously the two way coupling calculations are more physical and reasonable. 3.3. Carrier turbulence damping For the carrier turbulence damping, it is clear that carrier velocity profiles are the same as Fig. 6 shows regardless of using carrier turbulence damping function or not. Using the damping function, the turbulence kinetic energy is damped as Fig. 7 shows. The droplet velocity profile is shown in Fig. 8 and the enlarged droplet velocity profile is shown in Fig. 9. The turbulence kinetic energy is damped due to the huge amount of droplets using the
Fig. 6. Carrier velocity profile between the results with damping function and without damping function.
Fig. 8. Droplet velocity profile between the results with damping function and without damping function.
Fig. 9. Droplet velocity profile between the results with damping function and without damping function (enlarged view).
R. Li et al. / Nuclear Engineering and Design 242 (2012) 157–163 Table 2 Calculation results comparison.
163
4. Conclusions
Case
Impingement angle function
Erosion rate (kg m−2 s−1 )
Erosion rate deviation (%)
1 2
For solid particles (Eq. (17)) For liquid droplets (Eq. (18))
1.56 × 10−6 1.39 × 10−6
100 −10.9
damping function, the droplets velocity are influenced which is also damped (please see the enlarged droplet velocity profile in Fig. 9). Due to the erosion rate highly dependent on droplet velocity, the droplet velocity damp is not a neglected factor to the whole erosion. We could easily predict that the erosion rate will be damped. 3.4. Droplet impingement angle function After the analysis of droplets velocity profile comparison among different considerations, the next step is to investigate the droplet impingement angle function effect on erosion rate. From Table 2 we can see that the erosion rate will be less than 11% using the impingement impact angle function. We can see that it is applicable for the droplet impingement angle function consideration for the LDI rate evaluation. 3.5. Parametric investigations Since the droplet velocity is the key issue for LDI erosion prediction, the velocity parameter survey calculations were carried out, then the comparison with a real plant accident data was performed to analyze our current erosion model. As shown in Fig. 1 considering the rupture hole area and the operation times, the actual erosion rate can be easily reckoned as 4.58 × 10−7 kg m−2 s−1 (JSME, 2007). The trend of overall erosion rate increases as the velocity increases, which means that the logarithmic erosion rate is highly proportional to the velocity. Ejecting with the same droplets number (N in Eq. (16)), as the droplet velocity increase, the droplet will have more momentum thus the droplets number which impinges onto the wall surface (Nimpinged ) also increase, so the velocity exponent is over 2 larger than the d(v) = 1.73 in Eq. (16). This result is a combination of droplet velocity effect and impinged droplet numbers onto the wall surface. From velocity parameter survey shown in Fig. 10, we can make a prediction that droplet velocity is around 455 m s−1 , which might be reasonable for BWR piping system.
Fig. 10. Overall erosion rate vs droplet velocity.
Computational fluid dynamics calculations are conducted to investigate the damping characteristics of carrier fluid turbulence kinetic energy and erosion rate in a full-scale bent pipe. The carrier turbulent kinetic energy in vapor–droplet flow is calculated by a damping function on the energy spectrum basis of single phase flow. The inter-phase effects between the dispersed phase and the continuous phase are accounted for by an Eulerian–Lagrangian approach using two-way coupling calculation system. The major conclusions are drawn as follows: 1. Considering the carrier turbulence damping due to the suspended liquid droplets, a general LDI erosion prediction procedure for bent pipe geometry has been performed to supplement a commercial available CFD code with a new impact angle function. 2. Three computational results using major models (two way coupling, turbulence damping and droplet impingement angle function) are compared by droplet velocity profile and erosion rate to improve that the models make the calculations more physical. 3. The droplet velocity parameter calculations effect is evaluated based on our current methodology and calculation model. Compared to a real nuclear power plant accident data, we made the velocity prediction which may be useful for the industry engineers and annual maintenances. Acknowledgement The authors are grateful to the anonymous reviewers’ comments of his/her knowledge of algorithms and CFD codes, and valuable work to contribute to this paper. References Bottoni, M., Ninokata, H., 1998. An algorithm for damping of turbulence in particulate flow linked to the fluid-dynamic code COMM IX-M. Journal of Nuclear Science and Technology 35, 101–112. Edwards, J.K., McLaury, B.S., Shirazi, S.A., 2000. Evaluation of alternative pipe bend fittings in erosive service. In: Proceedings of ASME FEDSM‘00: ASME 2000 Fluids Engineering Division Summer Meeting, Boston, June 2000. FLUENT Version 6.3 Help Documentation Methodology. ANSYS Group, 2002. Japan Society of Mechanical Engineers Report, 2007. Understanding of wall thinning phenomena due to the liquid droplet impingent in a bent pipe (in Japanese). Launder, B.E., Spalding, D.B., 1972. Lectures in Mathematical Models of Turbulence. Academic Press, London, England. Li, R., Ninokata, H., 2011. A numerical study of impact force caused by liquid droplet impingement onto a rigid wall. Progress in Nuclear Energy 53, 881– 885. Oka, Y.I., Okamura, K., Yoshida, T., 2005. Practical estimation of erosion damage caused by solid particle impact, Part 1: Effects of impact parameters on a predictive equation. Wear 259, 95–101. Oka, Y.I., Yoshida, T., 2005. Practical estimation of erosion damage caused by solid particle impact, Part 2: Mechanical properties of materials directly associated with erosion damage. Wear 259, 102–109. Okada, H., Uchida, S., 2009. Evaluation on local wall thinning of piping due to liquid droplet impingement by coupled analysis of corrosion and flow dynamics. In: Proceedings of The 13th International Topical Meeting on Nuclear Reactor Thermal Hydraulics (NURETH-13), Kanazawa City, Ishikawa Prefecture, Japan, September 27–October 2, 2009. Okita, R., Shirazi, S.A., 2010. Experimental and CFD investigations to evaluate the effects of fluid viscosity and particle size on erosion damage in oil and gas production equipment. In: Proceedings of the ASME 2010 3rd Joint USEuropean Fluids Engineering Summer Meeting and 8th International Conference on Nanochannels, Microchannels, and Minichannels, Montreal, Canada, August 1–5, 2010. Springer, G.S., 1976. Erosion by Liquid Impact. John Wiley and Sons, New York, NY. Saffman, P.G., 1965. The lift on a small sphere in a slow shear flow. The Journal of Fluid Mechanics 22, 385–400. Sugiyama, H., 2003. Introductions to Mechanical Engineering: Fluid Mechanics. Morikita Publication Co., Ltd, pp. 73–74. Tabakoff, W., Kotwal, R., 1979. Erosion study of different material by coal ash particles. Wear 52, 161–170. Tennekes, H., 1972. A First Course in Turbulence. MIT Press, Cambridge, MA.
Contents lists available at SciVerse ScienceDirect
Nuclear Engineering and Design journal homepage: www.elsevier.com/locate/nucengdes
A calculation methodology proposed for liquid droplet impingement erosion Rui Li a,∗ , Michitsugu Mori b,c , Hisashi Ninokata a a
Research Laboratory for Nuclear Reactors, Tokyo Institute of Technology, 2-12-1-N1-5, Ookayama, Meguro-ku, Tokyo 152-8550, Japan Research and Development Centre, Tokyo Electric Power Company, 4-1, Egasaki-cho, Tsurumi-ku, Kanagawa 230-8510, Japan c School of Science and Technology, Meiji University, 1-1-1, Higashi-Mita, Tama-ku, Kawasaki-shi, Kanagawa 214-8571, Japan b
a r t i c l e
i n f o
Article history: Received 6 May 2011 Received in revised form 3 October 2011 Accepted 5 October 2011
a b s t r a c t Bent pipe wall thinning has been often found at the elbow of the drain line and the high-pressure secondary feed-water bent pipe in nuclear reactors. Liquid droplet impingement (LDI) erosion could be regarded as one of the major causes and is a significant issue of the thermal hydraulics and structural integrity in aging and life extension for nuclear power plant safety. In this paper a computational methodology is established for simulation of LDI erosion using computational fluid dynamics (CFD) simulation and theoretical calculation. Two-phase flow numerical simulations are conducted for standard elbow geometry, typically with the pipe diameter of 170 mm. This computational fluid model is built up by incompressible Reynolds Averaged Navier–Stoke equations using standard k–ε turbulence model and the SIMPLE algorithm, and the numerical droplet model adopts the Lagrangian approach. The turbulence damping in vapor–droplets flow is theoretically analyzed by a damping function on the energy spectrum basis of single phase flow. Locally, a droplet impact angle function is employed to determine the overall erosion rate. Finally, the overall and local investigations are combined to purpose a general methodology of LDI erosion prediction procedure, which has been complemented into CFD code. Based on our more physical computational results, comparison with an available accident data was made to prove that our methodology could be an appropriate way to simulate and predict the bent pipe wall thinning phenomena. © 2011 Elsevier B.V. All rights reserved.
1. Introduction In the area of safety analysis for power engineering industry, the bent pipe wall thinning has been often found in nuclear power plants as well as the fossil power plants. For example, JSME annual report records one accident that the rupture pipe is located as Fig. 1 shows, the small hole is also depicted (JSME, 2007). Liquid droplet impingement (LDI) erosion is one of the major causes thus evaluating LDI erosion is an important issue of the thermal hydraulics and structural integrity in aging and life extension for nuclear power plants safety. In general, major parameters that influence LDI erosion are categorized into four groups: flow pattern, chemical conditions, droplet impingement properties and material properties (Okada and Uchida, 2009). The droplet erosion mechanism is not similar with particulate solid sand erosion (Springer, 1976). From the point of flow pattern view, erosion rate is dependent on a number of factors including the droplet size, impact velocity,
∗ Corresponding author. Tel.: +81 3 5734 3062; fax: +81 3 5734 3056. E-mail address: [email protected] (R. Li). 0029-5493/$ – see front matter © 2011 Elsevier B.V. All rights reserved. doi:10.1016/j.nucengdes.2011.10.004
impact frequency, and liquid and gas density and viscosity. As many of these values are unknown for field situations, it is not easy to predict the rate of droplet erosion. It should also be noted that control of many of these factors in laboratory-based tests is problematical, for instance, the droplet velocity, is difficult to be controlled and measured. Therefore, a great deal of care should be required when extrapolating (from low velocity to high velocity) laboratory test results to practical conditions. Hence computational fluid dynamics (CFD) can be used as a valid tool to investigate the LDI phenomena. There are few equations used to predict liquid droplet erosion, although there are many equations developed to predict solid particle erosion in literatures (Tabakoff and Kotwal, 1979; Oka et al., 2005; Oka and Yoshida, 2005). The erosion equations are generally based on particle impact speed and angle, particle and material properties. However, in case of flow in a bent pipe with 90 degree, various flow situations and fluid properties influence particle impact speed and angle. Therefore CFD can be used to predict particle impact speed and angle which can be employed in erosion models to predict erosion. Our model to predict erosion is based on solid particle Edwards (2000) erosion model, and then developed with our considerations of liquid droplets compressibility and
158
R. Li et al. / Nuclear Engineering and Design 242 (2012) 157–163
Nomenclature a b Ad An Ai A an c C Cd dd d (v) E () F f (˛) g H() k l m M p t Re R R r U ud
v Vd W x y Greek ˛ ε
constant in turbulence theory constant in turbulence theory droplet cross-sectional area (m2 ) amplitude in velocity spectrum of fluid (m) constant wall face area (m2 ) amplitude in velocity spectrum of droplets (m) a constant in turbulence theory erosion formulation constant drag force coefficient droplet diameter (m) droplet impact velocity function spectrum of turbulent kinetic energy (m3 s−2 ) force (Pa) impingement angle function gravity acceleration (m s−2 ) damping function of carrier turbulence due to droplets turbulent kinetic energy (m2 s−2 ) macroscopic length scale of turbulence (m) droplets mass (kg) momentum loss or gain (kg m s−1 ) pressure (Pa) time (s) Reynolds number normalized relative velocity between two phases (m s−1 ) erosion rate (kg m−2 s−1 ) radius (m) mean velocity (m s−1 ) droplet velocity (m s−1 ) velocity (m s−1 ) droplet volume (m3 ) mass concentration of droplets coordinate direction (m) coordinate direction (m)
ıtd
volume fraction, droplet impingement angle dissipation rate of turbulent kinetic energy (m3 s−2 ) Kolmogorov microscale of turbulence (m) dynamic viscosity (kg m−1 s−1 ) wave number (m−1 ) Kolmogorov constant density (kg m−3 ) time step (s)
Indices am c d dr n o tp
virtual mass force carrier fluid droplets drag force new value old value two phase
deformation. In this paper, we establish an innovative methodology of LDI erosion rate for evaluating wall thinning phenomena. We proposed a calculation methodology to provide useful information for industry applications. The current paper originalities are (1) applying a droplet–vapor two-way coupling calculation system and a vapor turbulence damping analysis to LDI phenomena, (2)
implementing an impingement angle function into codes for LDI erosion rate calculation in a full scale bent pipe, and (3) making a comparison with an actual accident erosion data to validate our proposed methodology. 2. Computational models and methods 2.1. Outline In this section we introduce our calculation methodology including two-way vapor–droplet coupling computational method, calculation with carrier turbulence damping theory calculation, local impact force function, erosion semi-empirical formulation and finally the procedure to implement codes using a CFD tool and a set of subroutines. The standard elbow bent pipe computational domain and mesh are also given. 2.2. Two way coupling vapor–droplet method Two-phase flow characteristics could be predicted by CFD, a momentum transfer consideration called two-way coupling between carrier fluid and transported droplets is taken into account. The droplet governing equation is expressed in Lagrangian description: (und − uod )
= Fndr + Fnp + Fnam + Fng , ıtd 1 Fdr = − Cd c Ad |ud − u|(ud − u), 2 Fp = −Vd ∇ p, d(ud − u) 1 Fam = − c Vd , 2 dt Fg = md g, md
(1)
where md is the droplet mass; u and ud are the continuous phase velocity and droplet velocity, respectively; n and o denote “new” and “old” values spanning the integration time step ıtd , respectively. In the droplet Lagrangian description, Fdr , Fp , Fam and Fg are the drag force, pressure force, virtual mass force and the gravity force, respectively. ∇ p is the pressure gradient in the carrier fluid, Vd is the droplet volume, Ad is the droplet cross-sectional area, Cd is the drag coefficient, and g is the gravitational acceleration. The lift force is not considered due to the fact that the typical particle diameter in LDI phenomenon is around 10 m rather than submicro meter (Saffman, 1965). As for the carrier gas, its momentum equation is: ∂ui ∂ ∂p + (uj ui − ij ) = − , ∂t ∂xj ∂xi
(2)
where ui is the carrier fluid velocity component in direction xi and the term ij represents the stress tensor components and need closure. In this current study we use the standard k–ε model which has become the workhorse of practical engineering flow calculations for pipe flows (Launder and Spalding, 1972). Here this calculation system is called “one way coupling” system. In two way coupling system, the momentum exchanges between carrier fluid and droplet are obtained from the integration of the Lagrangian equations. The vapor–droplet coupled model accounts for the fact that each droplet may exchange momentum with the carrier fluid, and that a loss from a droplet is a gain for the fluid and vice versa. From the average rate of momentum loss or gain of a single droplet due to drag and other effects over a single time increment ıtd , the total momentum loss or gain of a group of droplets traversing a given cell is expressed as follows: Md =
all droplets
md
(und − uod ) ıtd
− Fnp
− Fng
,
(3)
R. Li et al. / Nuclear Engineering and Design 242 (2012) 157–163
159
Fig. 1. The small hole as a result of pipe thinning and the schematic drawing (JSME, 2007).
which is inserted as a “source” in the momentum equation of the carrier phase. So the coupling between the fluid and the droplet is achieved via the sink term Md in the momentum equation (2) of the carrier fluid: ∂ ∂p ∂ui + (uj ui − ij ) = − + Md , ∂t ∂xj ∂xi
(4)
2.3. Damps of turbulent kinetic energy due to small droplets Since the turbulence of carrier fluid might be attenuated due to the involved droplets, the impinged droplets erosion rate might be affected by the turbulence of carrier fluid. It is applicable to add the turbulence kinetic energy damping of carrier fluid consideration for the LDI rate evaluation (Li and Ninokata, 2011). The turbulence of carrier fluid contains a wide range of eddy sizes. The strain rate of large or “energy-containing” eddies is proportional to the average strain rate of the mean flow. The strain rate of small eddies is larger compared to that of the mean flow and of the large eddies, so that no permanent anisotropy can be induced at small scales. The spectrum of turbulent kinetic energy, expressed as a function of the wave number , could consist of four sub ranges which are defined as follows (Tennekes, 1972): (1) Larger scale subrange :
(2) Transition subrange :
−3/4
E0 () = cε
c c
E1 () = c(a2 + b),
17/4 4
,
0 ≤ ≤ 0 = 10,
0 ≤ ≤ l =
(5)
to compute the damping of turbulent kinetic energy of the carrier due to the droplets:
H() = exp −
(6)
R2 −5/3 d ,
(9)
0
where dd is the droplet diameter, d is droplet density, W is the mass concentration of the droplets, and R is the normalized relative velocity between carrier and droplets. They are given as: ˛d d , ˛d d + ˛c c
W=
R =
1+
a 2 n
An
an −2 · An
(10)
1/2
1
1 + Uc2 2 (dd d /18c )
,
(11)
where An is the amplitude in velocity spectrum of fluid, an is the amplitude in velocity spectrum of droplets, and Uc is the carrier mean velocity. Hence when the vapor flow with a dilute droplet phase is simulated, the damping of turbulent kinetic energy of the carrier fluid due to the droplets should be computed. Then the revised kinetic energy values are used for erosion evaluation. The formalism is applied through the following procedure: Step 1 (:). The standard k–ε model is utilized to the carrier to compute the integral value of the single pure phase turbulent kinetic energy neglecting the effect of droplets:
1 , l
36Wc ε−1/3 dd d
∞
k=
E(v)dv,
(12)
0
(3) Inertial subrange :
E2 () = cε2/3 −5/3 ,
(4) Equilibrium subrange :
=
1 ≤ ,
l =
E3 () = cε2/3 −5/3
1 1 ≤ ≤ = , l
(7)
3 4/3 exp − () , 2
Step 2 (:). The theoretical distributions of the single pure phase turbulent spectra are computed by self-developed original program, as Section 2.1 discussed in four intervals. E(v) (0 ≤ v ≤ ∞),
(8)
The above four sub ranges are for the single pure phase turbulent energy spectrum. As it was pointed out that the dispersed phase may have an influence to the continuous phase, the influence is considered: the dispersion of droplets due to the continuous fluid phase and, vice versa, the impact of the droplets upon the turbulence intensity of the fluid. According the theoretical analysis, the suppression of turbulent kinetic energy occurs in case of small particles. Bottoni and Ninokata (1998) proposed a damping function
(13)
Step 3 (:). The damping function of carrier turbulence H() is calculated by Eqs. (9)–(11), and then the two-phase flow spectra of turbulent kinetic energy are computed by: E(v)tp = E(v)H(v),
(14)
Step 4 (:). The two-phase flow spectra of turbulent kinetic energy are integrated numerically over the full range of wave numbers to obtain the integral value:
ktp =
∞
E(v)tp dv. 0
(15)
160
R. Li et al. / Nuclear Engineering and Design 242 (2012) 157–163
Fig. 2. Schematic diagrams and bent section mesh of a standard elbow.
The value ktp should replace the single phase value k cell by cell in the computational domain of the carrier fluid.
be done. The solid particle (sand) impingement angle function is (Okita and Shirazi, 2010):
Step 5 (:). The standard k–ε model was applied again with the revised kinetic energy ktp using two-way coupling carrier-droplets computational system to evaluate the LDI erosion rate which is our interest and final target.
f (˛) = 5.4˛ − 10.1˛2 + 10.9˛3 − 6.3˛4 + 1.42˛5 ,
2.4. Overall erosion model In the erosion prediction model, the wall mass loss and penetration rates are calculated. When a droplet impinges the wall, the mass loss distributes uniformly over the computational cell of the wall. With both sufficiently small grid spacing and a large number of droplets, approximation errors induced by these assumptions still cannot be eliminated. Once all droplet trajectories have been computed and all wall impingement data is gathered, the total mass loss for all impingements can be compiled to generate a local penetration rate for each cell which lies on the surface of the geometry. Through the subroutines, the droplet impingement information is gathered. The computed erosion rate is the removed mass per area and time in the unit kg m−2 s−1 , which is defined by (Edwards, 2000)
m C(d )f (˛)vd(v) d d
0◦ < ˛ ≤ 90◦ , (17)
Unlike solid particle erosion the scratch stress caused by tangential momentum might be negligibly small, in case of LDI erosion. Therefore, the scratch stress is much less important due to the liquid mobility and deformation. The normal component of impact velocity is major factor for the LDI erosion damage whilst the tangential component has little effect. The effect of droplet impingement angle on LDI erosion was investigated by the authors (Li and Ninokata, 2011) using the volume of fluid (VOF) computational system, which is a two-phase Eulerian–Eulerian approach. Based on our previous conclusion the impingement angle function f(˛) can be fitted as an approximation shown in Eq. (18), which is also implemented into codes to support the calculation of the overall erosion rate. f (˛) = 1.3 × 10−1 + 1.5 × 10−2 ˛ − 6.4 × 10−5 ˛2 ,
0◦ < ˛ ≤ 90◦ . (18)
Nimpinged
Rerosion =
n=1
Aface
,
(16)
where N is the number of droplets which impinge onto the wall, Aface is the area of the cell face at the wall. C (dd ) is a function of the droplet diameter function and impinged wall density (usually stainless steel), 1.18 × 10−11 without unit, v is the droplet velocity before impinging the wall. d(v) is the droplet impact velocity function, 1.73 without unit. C(dd ) and d(v) are determined from the target wall surface materials. f(˛) is the impingement angle function, ˛ is droplet impingement angle with the wall face. 2.5. Impingement angle function LDI erosion usually occurs in the inner surface of a bent pipe, with the droplets carried by the carrier fluid impinging on the surface with different impingement angles. In the Edward’s model the sand particles are assumed as impinged materials. However, in case of liquid droplets due to their compressibility and deformation after impingements, some reasonable modification should
2.6. Computational domain of a standard bent pipe The photograph of a full scale bent pipe with the flow direction is shown in Fig. 1 where the length unit is mm. It can be seen in Fig. 2, in order to use the uniform inlet condition and fully developed outflow condition, the computational domain is extended in both the upstream and downstream parts, and the two extended parts are called the pre-extended and after-extended region, respectively. In our numerical simulations, the pre-extended region is in length of 20 pipe diameters and the after-extended region is up in length of 10 pipe diameters. The mesh system is also shown in Fig. 2, where the first mesh adjacent to inner wall surface is also shown. In the boundary layer region we employ a much denser grid, whilst in the central region a coarse grid is adopted to save the computer resource. Through a network of subroutine to implement the two-way coupling method and erosion formulation including the impingement angle function, the governing equations along with the
R. Li et al. / Nuclear Engineering and Design 242 (2012) 157–163
161
Table 1 Main assumptions for the discrete phase model. Parameter (unit)
Value
Droplet density (kg m−3 ) Uniform diameter (m) Droplet number density (m−3 ) Total mass flow rate (kg s−1 ) Droplet volume fraction Droplet velocity (m s−1 )
998 10 8.78 × 1012 4.59 0.0001 400
3. Results 3.1. Discussion
Fig. 3. The vapor inlet velocity power-law distribution.
boundary conditions are solved using the commercial finitevolume based solver FLUENT 6.3. The SIMPLE algorithm is adopted to deal with the linkage between velocity and pressure. Following assumptions are made in the numerical simulation: (1) the fluid flow is in steady-state; (2) the effect of the gravitational force is negligible; (3) the fluid thermo-physical properties are constant; and (4) the Mach number is less than 0.3, that is the compressibility effect can be the carrier fluid is incompressible.
2.7. Boundary conditions and initial conditions The carrier fluid is vapor. The required conditions are described for the three different parts, which are inlet, wall and outlet, respectively. Following descriptions are the boundary conditions. (1) The implementation of pipe wall boundary is treated as the normal no-slip wall and the wall surface is supposed to be smooth. (2) The vapor density is 0.554 kg m−3 at room pressure (1.01 × 105 Pa) and room temperature (298 K). (3) The outlet is prescribed as outflow. The liquid water density is 998.0 kg m−3 , which is around 1800 times larger than vapor, which means that this large density ratio will bring strong effects on the droplets trajectories towards to impinge the wall. The vapor is the fully developed turbulent flow in the bent pipe, where the vapor velocity inlet profile usually presents a power-law distribution (Sugiyama, 2003). For instance if the inlet vapor mean velocity is 280 m s−1 , the corresponding distribution is shown in Fig. 3, where r is the radius from wall to the center of the pipe. In our calculations we assume that the initial droplet velocity has the same value as the initial mean vapor velocity, and in this paper we briefly refer them as droplet velocity. The impact of droplets on rigid surfaces may produce a wide variety of consequences, according to the size, velocity and material of the impacting elements and the nature of the surface. For example, droplets may adhere, bounce or shatter, and the liquid deposited from the wall surface may retain its droplet form or merge into a liquid film. In our numerical simulations, the treatment of wall and droplet interaction is simplified by the perfect rebound function. Due to the high Reynolds number, we use the standard wall functions. In the practical operating conditions, the droplet parcel should consist of various droplets with different diameters; in this research we define uniform diameter distribution. The droplet number density is assumed to be uniformly at the inlet.
From the above sections three major considerations are taken into account for LDI erosion prediction: (1) two-way coupling, (2) carrier turbulence damping, and (3) droplet impingement angle function. The main target of this current work is the droplet erosion, which is highly dependent on droplet velocity profile as shown in Eq. (16). Since this study is from the point view of fluid dynamics, it is important to discuss the details of each consideration and compare with those results without the above considerations. Therefore we can know that three considerations are absolutely necessary, thus it is approved that our methodology might be correct to predict the LDI erosion. It has to be confirmed that the calculated droplet velocity profiles are more physical using the models of the two-way coupling modeling and turbulence damping function. The erosion prediction with droplet impact angle function also should be made a comparison with the erosion prediction with solid particle impact angle function. Some discussions are addressed in this section. The main assumptions for the discrete droplets are listed in Table 1, which are typical values in LDI phenomena. In our calculations we assume that the initial droplet velocity has the same value as the initial mean vapor velocity, we briefly refer them as droplet velocity. 3.2. Two way coupling In order not to be influenced by the elbow geometry, the carrier and droplet normalized velocity profile are post-processed in a section ahead located in the length of 5 pipe diameters before the elbow. Figs. 4 and 5 show that the vapor and droplet velocity profiles calculated from one way coupling system and two way coupling system, respectively. Please note that unlike the vapor
Fig. 4. Carrier velocity profile using one way coupling and two way coupling system.
162
R. Li et al. / Nuclear Engineering and Design 242 (2012) 157–163
Fig. 5. Droplet velocity profile using one way coupling and two way coupling system.
Fig. 7. Carrier turbulence kinetic energy profiles between the results with damping function and without damping function.
velocity profile can be extracted directly from the calculation mesh value, the droplet velocity profile was a statistic value based on the sampling of Lagrangian tracking, so there might be some discontinues in the distributions (ANSYS, 2002). Due to the huge amount droplets in the carrier vapor, two way coupling calculation carrier velocity profile becomes ‘more laminar’ than single phase flow. As a result the droplet velocity profile also became “more sharp” with higher velocity in the core. Based on the erosion calculation equation (Eq. (16)), we can make a prediction that the one way coupling erosion result should be lower than two way coupling results. Obviously the two way coupling calculations are more physical and reasonable. 3.3. Carrier turbulence damping For the carrier turbulence damping, it is clear that carrier velocity profiles are the same as Fig. 6 shows regardless of using carrier turbulence damping function or not. Using the damping function, the turbulence kinetic energy is damped as Fig. 7 shows. The droplet velocity profile is shown in Fig. 8 and the enlarged droplet velocity profile is shown in Fig. 9. The turbulence kinetic energy is damped due to the huge amount of droplets using the
Fig. 6. Carrier velocity profile between the results with damping function and without damping function.
Fig. 8. Droplet velocity profile between the results with damping function and without damping function.
Fig. 9. Droplet velocity profile between the results with damping function and without damping function (enlarged view).
R. Li et al. / Nuclear Engineering and Design 242 (2012) 157–163 Table 2 Calculation results comparison.
163
4. Conclusions
Case
Impingement angle function
Erosion rate (kg m−2 s−1 )
Erosion rate deviation (%)
1 2
For solid particles (Eq. (17)) For liquid droplets (Eq. (18))
1.56 × 10−6 1.39 × 10−6
100 −10.9
damping function, the droplets velocity are influenced which is also damped (please see the enlarged droplet velocity profile in Fig. 9). Due to the erosion rate highly dependent on droplet velocity, the droplet velocity damp is not a neglected factor to the whole erosion. We could easily predict that the erosion rate will be damped. 3.4. Droplet impingement angle function After the analysis of droplets velocity profile comparison among different considerations, the next step is to investigate the droplet impingement angle function effect on erosion rate. From Table 2 we can see that the erosion rate will be less than 11% using the impingement impact angle function. We can see that it is applicable for the droplet impingement angle function consideration for the LDI rate evaluation. 3.5. Parametric investigations Since the droplet velocity is the key issue for LDI erosion prediction, the velocity parameter survey calculations were carried out, then the comparison with a real plant accident data was performed to analyze our current erosion model. As shown in Fig. 1 considering the rupture hole area and the operation times, the actual erosion rate can be easily reckoned as 4.58 × 10−7 kg m−2 s−1 (JSME, 2007). The trend of overall erosion rate increases as the velocity increases, which means that the logarithmic erosion rate is highly proportional to the velocity. Ejecting with the same droplets number (N in Eq. (16)), as the droplet velocity increase, the droplet will have more momentum thus the droplets number which impinges onto the wall surface (Nimpinged ) also increase, so the velocity exponent is over 2 larger than the d(v) = 1.73 in Eq. (16). This result is a combination of droplet velocity effect and impinged droplet numbers onto the wall surface. From velocity parameter survey shown in Fig. 10, we can make a prediction that droplet velocity is around 455 m s−1 , which might be reasonable for BWR piping system.
Fig. 10. Overall erosion rate vs droplet velocity.
Computational fluid dynamics calculations are conducted to investigate the damping characteristics of carrier fluid turbulence kinetic energy and erosion rate in a full-scale bent pipe. The carrier turbulent kinetic energy in vapor–droplet flow is calculated by a damping function on the energy spectrum basis of single phase flow. The inter-phase effects between the dispersed phase and the continuous phase are accounted for by an Eulerian–Lagrangian approach using two-way coupling calculation system. The major conclusions are drawn as follows: 1. Considering the carrier turbulence damping due to the suspended liquid droplets, a general LDI erosion prediction procedure for bent pipe geometry has been performed to supplement a commercial available CFD code with a new impact angle function. 2. Three computational results using major models (two way coupling, turbulence damping and droplet impingement angle function) are compared by droplet velocity profile and erosion rate to improve that the models make the calculations more physical. 3. The droplet velocity parameter calculations effect is evaluated based on our current methodology and calculation model. Compared to a real nuclear power plant accident data, we made the velocity prediction which may be useful for the industry engineers and annual maintenances. Acknowledgement The authors are grateful to the anonymous reviewers’ comments of his/her knowledge of algorithms and CFD codes, and valuable work to contribute to this paper. References Bottoni, M., Ninokata, H., 1998. An algorithm for damping of turbulence in particulate flow linked to the fluid-dynamic code COMM IX-M. Journal of Nuclear Science and Technology 35, 101–112. Edwards, J.K., McLaury, B.S., Shirazi, S.A., 2000. Evaluation of alternative pipe bend fittings in erosive service. In: Proceedings of ASME FEDSM‘00: ASME 2000 Fluids Engineering Division Summer Meeting, Boston, June 2000. FLUENT Version 6.3 Help Documentation Methodology. ANSYS Group, 2002. Japan Society of Mechanical Engineers Report, 2007. Understanding of wall thinning phenomena due to the liquid droplet impingent in a bent pipe (in Japanese). Launder, B.E., Spalding, D.B., 1972. Lectures in Mathematical Models of Turbulence. Academic Press, London, England. Li, R., Ninokata, H., 2011. A numerical study of impact force caused by liquid droplet impingement onto a rigid wall. Progress in Nuclear Energy 53, 881– 885. Oka, Y.I., Okamura, K., Yoshida, T., 2005. Practical estimation of erosion damage caused by solid particle impact, Part 1: Effects of impact parameters on a predictive equation. Wear 259, 95–101. Oka, Y.I., Yoshida, T., 2005. Practical estimation of erosion damage caused by solid particle impact, Part 2: Mechanical properties of materials directly associated with erosion damage. Wear 259, 102–109. Okada, H., Uchida, S., 2009. Evaluation on local wall thinning of piping due to liquid droplet impingement by coupled analysis of corrosion and flow dynamics. In: Proceedings of The 13th International Topical Meeting on Nuclear Reactor Thermal Hydraulics (NURETH-13), Kanazawa City, Ishikawa Prefecture, Japan, September 27–October 2, 2009. Okita, R., Shirazi, S.A., 2010. Experimental and CFD investigations to evaluate the effects of fluid viscosity and particle size on erosion damage in oil and gas production equipment. In: Proceedings of the ASME 2010 3rd Joint USEuropean Fluids Engineering Summer Meeting and 8th International Conference on Nanochannels, Microchannels, and Minichannels, Montreal, Canada, August 1–5, 2010. Springer, G.S., 1976. Erosion by Liquid Impact. John Wiley and Sons, New York, NY. Saffman, P.G., 1965. The lift on a small sphere in a slow shear flow. The Journal of Fluid Mechanics 22, 385–400. Sugiyama, H., 2003. Introductions to Mechanical Engineering: Fluid Mechanics. Morikita Publication Co., Ltd, pp. 73–74. Tabakoff, W., Kotwal, R., 1979. Erosion study of different material by coal ash particles. Wear 52, 161–170. Tennekes, H., 1972. A First Course in Turbulence. MIT Press, Cambridge, MA.