Powder Technology 334 (2018) 9–26
Contents lists available at ScienceDirect
Powder Technology journal homepage: www.elsevier.com/locate/powtec
An investigation of erosion prediction for 15° to 90° elbows by numerical simulation of gas-solid flow M.R. Banakermani, Hamid Naderan ⁎, Majid Saffar-Avval Mechanical Engineering Dept., Amirkabir University of Technology, Tehran, Iran
a r t i c l e
i n f o
Article history: Received 18 August 2017 Received in revised form 26 January 2018 Accepted 14 April 2018 Available online 17 April 2018 Keywords: Erosion Elbow Gas-solid OpenFOAM Eulerian-Lagrangian Particle impact frequency
a b s t r a c t Erosion is a common problem in various parts of the pipeline industry. In this study, computational fluid dynamics is employed for analysis of pipeline erosion due to gas transported particles. Total removed volume and maximum erosion rate in different elbow geometries from 15° to 90°elbows, for a range of mass loadings are studied. Particle tracking is achieved using combination of Eulerian and Lagrangian methods. In the first stage, gas-solid flow is validated for both straight-pipe and elbow geometries. After that, the erosion model is validated with available data for a 90°, 76.2 mm diameter standard elbow. In the final stage, simulations for a range of elbow angles are performed for two different flow orientations: horizontal inlet and outlet flow directions (H-H flow) and the cases in which inlet flows are vertical and outlet flows either lie in the horizontal plane or make some angles with it (V-H flow). In addition to erosion rate, some important particle-wall impact-related variables such as impact speed, impact angle and impact frequency are presented for elbows with different configurations and angles. Results show that, in general, for a fixed inlet condition and bend geometry, maximum erosion rate in the V-H configuration is greater than that for the H-H orientation. However, total annual eroded volume in the H-H configuration is greater than that for the V-H case. In addition, in both V-H and H-H cases, for the range of mass loadings investigated, the maximum erosion rate increases steadily when the elbow angle increases from 15° to 90° but the rate of total eroded volume remains relatively constant for each value of sand rate. Results of this study are helpful in optimal selection of pipeline elbows for erosion protection. © 2018 Published by Elsevier B.V.
1. Introduction Solid particle erosion is a common problem in different industries including the oil and gas industries. Pipes, fittings and rotary machinery in gas transport pipelines are subject to the impact of particles that are carried over with gas. Although particle filters and cyclone separators are widely used in pipeline network, the presence of debris particles in pipes is unavoidable. In addition to the sand coming from wells, particles are constantly produced due to some unwanted and not completely understood reactions that take place inside the pipeline. These particles are called black powder and occur in wide ranges of size and composition as reported in the literature [1,2]. Erosion may cause rupture and failure of pipes, so in order to maintain safe conditions, estimation of the potential damage of particles on wall material should be properly regarded during the design process. One of the most sensitive locations, prone to erosion is where the flow changes direction in pipe bends. Bends with different angles are widely used in pipeline industry. Different forming processes are employed to produce field cold bends, factory made bends or elbows, mitered bends and elasto-plastic bends. Due to the inertia of particles, they do ⁎ Corresponding author. E-mail address:
[email protected]. (H. Naderan).
https://doi.org/10.1016/j.powtec.2018.04.033 0032-5910/© 2018 Published by Elsevier B.V.
not follow the flow streamline in bends and hit the wall, causing damage to the surface. Different investigators have proposed many erosion equations and models. Some models are empirical, some mechanistic and some semi-mechanistic. A variety of parameters is reported to affect the amount of mass being eroded. Some parameters are related to the properties of particles such as impact speed, impact angle, size, density, hardness and shape of particles. Another class of parameters are related to the target wall material properties such as ductility, density and hardness. The final set of parameters are those that are related to the flow characteristics such as fluid density, viscosity, turbulent fluctuations and flow pattern in the considered geometry. Parsi et al. [3] reviewed these parameters and some of the erosion equations which were proposed in the literature. One of the earliest erosion equations is the model proposed by Finnie et al. [4]. This model calculates the volume removed by a single particle's cutting action when impinging a ductile wall material. At small impact angles, the model works fine but for near normal impacts, the model does not predict the removed volume accurately [4,5]. Hutchings [6,7] used a critical plastic strain limit to theoretically develop an equation for erosion of metals during normal impact by spherical particles. He also performed some tests on ductile metals, under impact with spherical and sharp particles at oblique angles, and described the mechanism of material removal. It was found
10
M.R. Banakermani et al. / Powder Technology 334 (2018) 9–26
that material removal from work hardened metals, due to the concentration of deformation in the surface layer, is easier compared to annealed metals for which the impact energy is spread through a large volume during collision [6]. Bitter [8,9] described the erosion mechanism as the effect of simultaneous action of cutting wear and repeated deformation. Other researchers employed equations proposed by Bitter and Hutchings and calibrated their constants. Huang et al. [10] developed a mechanistic erosion model and calibrated the constants used in their model by testing on some pure materials. In the model, it was assumed that the particle's tangential velocity component was responsible for cutting removal, while deformation damage removal was related to normal velocity component. Jafari et al. [11] modified Huang et al.'s [10] model by considering the effect of fluid turbulence on particle impact velocity in the erosion equation and investigated erosion rates in horizontal pipes. In a recent work, Arabnejad et al. [12] developed a semi-mechanistic erosion model by performing several direct impingement tests on different target materials including alloy materials such as stainless steel and carbon steel. Much experimental and numerical research has been performed on erosion in 90° elbows. Chen et al. [13] examined the effect of particlewall rebound models on erosion pattern in a 90°bend and tested their model by performing some experimental erosion tests. Pereira et al. [14], numerically investigated the effect of some parameters such as number of computational particles, wall surface roughness models and wall restitution and friction coefficients on erosion pattern and the amount of eroded material in a 90° elbow. They found that while the erosion pattern predicted by all of the studied empirical erosion equations were similar, the model proposed by Oka et al. [15,16] provided the closest estimation with experimental data for the erosion penetration rate [14]. However, it should be noted that Oka et al. [15,16] model is developed under experiments with high speed flow of gases and doesn't seem to be suitable for gas pipeline industries. Duarte et al. [17], in their numerical research, found that increase in mass loading could cause a reduction of the penetration ratio due to cushioning effect. So they recommended that inter-particle collisions must be always taken into account even at low mass loadings [17]. In the majority of numerical research, using various erosion equations, a combination of semi-elliptic and “V”-shaped erosion scars in 90°elbows was reported. However, recently Solnordal et al. [18], conducted tests with a surface profiler and obtained a detailed map of erosion depth in a standard 90° elbow; The “V”-shaped scars in their metal loss patterns were not present and just the semi-elliptic part was observed [18]. They also found that the accuracy of numerical modeling in predicting mass loss pattern and maximum erosion depth was related to the ability of modeling rough wall collisions. In another work, Vieira et al. [19] used non-invasive ultrasonic technology to measure metal loss at 16 different locations in a standard stainless steel elbow caused by air-sand flow. Their tests were conducted with different values for air velocities and sand mass loadings. They also performed CFD simulations to compare the maximum erosion rate and its location with their experimental measurements for the standard elbow. The erosion model used in CFD simulations was generated from the results of PIV measurements of the slip velocity between sand and air in a direct impact geometry. Elbows with different angles exist in the pipeline industry for connecting two straight pipes together and changing the flow direction. Although many numerical and experimental works were performed in the literature on 90°elbows, to the best of the author's knowledge, there has been no such research on erosion in elbows with angles other than 90°. In this study one of the most recent erosion models, Arabnejad et al. [12] model, is used for investigating erosion patterns in elbows with different angles. Since every part of the numerical simulation may influence the results of surface erosion, validation is done for different parts of the simulation. The validation process consists of three stages; in the first stage, numerical simulations are performed to validate the developed
solver in OpenFOAM 2.3.1 for gas-solid flows. For this purpose, experimental LDV measurements in a straight pipe from Tsuji et al. [20] are used. In the second stage, gas-solid flow pattern is validated with LDV measurements of Kliafas et al. [21] in a 90°curved duct. In the final stage, experimental measurements of Vieira et al. [19] are selected to validate metal loss pattern and maximum erosion rate in a 90°standard elbow [19]. After the validation phase, elbows with different angles, from 15° to 90°, are simulated with OpenFOAM, using URANS method to solve for Navier-Stokes equations, RNG k-epsilon model for turbulence closure and Lagrangian particle tracking. One-way coupling for interaction of fluid phase and dispersed phase is considered. Particle-wall collisions are simulated by Sommerfeld's shadow model for wall roughness [22] and collisions between particles are modeled by Sommerfeld's stochastic approach [23]. Also the MOB dispersion model, introduced by Moissette, Oesterle and Boulet is applied to regenerate fluid field velocity fluctuations [24]. Amani and Nobari tuned the coefficients of the MOB model for simulation of evaporating sprays and the coefficients are set according to their work [25]. Results of the CFD simulation, give a more in depth insight about the underlying influential events that lead to material removal, such as impact angle, impact speed and impact frequency of the particles. Erosion mass loss patterns, maximum erosion rate and total annual volume loss due to erosion are calculated for elbows of different angles and several inlet flow velocities and sand mass loadings. 2. Mathematical model To model two-phase gas-solid flow, Eulerian formulation is employed for the fluid flow and Lagrangian particle tracking is used for analysis of particle motion. Both one-way and two-way coupling between the fluid and solid phases are used in validation cases and erosion calculation simulations. Sections 2.1 and 2.2 describe the governing equations for the fluid and solid phases respectively. 2.1. Fluid phase modeling The governing partial differential equations of fluid phase are the unsteady, incompressible, turbulent and isothermal flow equations, represented using tensor notation as follows: ∂ ðui Þ ¼ 0 ∂xi
ρf
" !# ∂ ui u j ∂ðui Þ ∂p ∂ ∂ui ∂u j þ ρf ¼− þ ðμ þ μ t Þ þ þ Suip þ ρ f gi ∂t ∂xi ∂x j ∂x j ∂xi ∂x j
ð1Þ
ð2Þ
In Eq. (2), Suip represents the momentum source term due to the presence of particles. When simulating using the one-way coupling assumption, this term is neglected. In case of two-way coupling, its magnitude equals the force applied by fluid on particles but in the opposite direction, according to Newton's third law. OpenFOAM is employed for numerical solution of Eqs. (1) and (2), using finite volume method and unstructured three-dimensional grids. The PIMPLE algorithm in OpenFOAM, which is the combination of SIMPLE and PISO algorithms, is used for coupling of pressure and velocity fields. During the simulations, to close the equation set and therefore to model μt, standard k − ε and RNG k − ε turbulence models are selected for straight pipe and elbow geometries, respectively. OpenFOAM has many schemes for numerical discretization of each term of the equations, but the combination used in the solver must be consistent and optimum. The schemes used in the current work are provided in Table 1.
M.R. Banakermani et al. / Powder Technology 334 (2018) 9–26
11
Table 1 Different discretization schemes (used in the fvSchemes file of project). Gradient terms discretization ∇p Gauss linear ! Gauss linear ∇U Time derivative discretization ∂ Euler method
Divergence terms discretization ! Bounded Gauss upwind ∇:ðφ U Þ ∇. (φk) Bounded Gauss upwind ∇. (φε) Bounded Gauss upwind Laplacian terms discretization 2 ∇ (φ) Gauss linear corrected
∂t
2.2. Particle phase modeling The fluid exerts several forces and a torque on the particles. Newton's second law governs the motion of the particles. The equations for linear and angular motion of particles are mp
Ip
dU p ¼ F g þ F D þ F SL þ F RL dt
ð3Þ
dωp ¼ TR dt
ð4Þ Fig. 2. Cross section view of grid used for validation of Tsuji et al.'s experiment.
Fg is the sum of gravity and buoyancy forces exerted on a particle: ρf F g ¼ mp g 1− ρp
! ð5Þ
FD is the drag force which is the most important force acting on particles and is calculated by assuming a spherical shape [26]. FD ¼
3ρ f mp C D U f −U p U f −U p 4ρp Dp
ð6Þ
The drag coefficient is calculated from the relation developed for spherical particles, as follows [26]: ( CD ¼
; 24 Rep −1 1 þ 0:15 Re0:687 p
Rep ≤ 1000
0:44;
Rep N 1000
ð7Þ
Rep is the particle Reynolds number and is computed by Eq. (8) based on the velocity of the particle relative to the fluid: ð8Þ Rep ¼ ρ f U f −U p Dp =μ
FSL is the force exerted on particle due to the presence of shear flow around it and is calculated by the relation developed by Saffman [27]: F SL ¼ 1:615Dp μ f Res0:5 C sl U f −U p ω f
ð9Þ
In Eq. (9), ωf = ∇ × Uf, is the vorticity of fluid, Res ¼ ρ f Dp 2 jω f j=μ f is by definition shear flow Reynolds number and, Csl is the coefficient of slip-shear lift force. Mei [28] and Crowe [29] proposed the following relation for a broad range of particle Reynolds numbers:
C sl ¼
( pffiffiffiffi pffiffiffiffi Re 1−0:331 β e− =10 þ 0:331 β; Rep ≤ 40 p pffiffiffiffiffiffiffiffiffiffi 0:0524 β Rep ; Rep N 40
ð10Þ
where β = 0.5Res/Rep. FRL is the slip-rotational lift force exerted on particle due to the difference between the fluid and particle rotational speed which is
Fig. 1. Particle-particle collision takes place on a collision cylinder locus [23].
12
M.R. Banakermani et al. / Powder Technology 334 (2018) 9–26
Fig. 3. Normalized axial velocity profiles of gas-solid flow on vertical symmetry plane of the pipe Ug = 10 m/s, (a) m = 0.4, (b) m = 1.1.
Fig. 4. Fluid turbulence intensity profiles on vertical symmetry plane of the pipe. Ug = 10 m/s, (a) m = 0.4, (b) m = 1.1.
Fig. 5. Particle velocity fluctuations intensity profile on vertical symmetry plane of the pipe. Ug = 10 m/s, (a) m = 0.4, (b) m = 1.1.
M.R. Banakermani et al. / Powder Technology 334 (2018) 9–26
13
Fig. 6. Computational grid used for validation with Kliafas et al.'s experiment. (a) side view, (b) cross section view.
Ω = 0.5∇ × Uf − ωp. Rubinow and Keller proposed Eq. (11) for this force [30]: F RL ¼
π 3 Rep Dp ρ f C Ω U f −U p 8 ReR rl
ð11Þ
where by definition, ReR ¼ ρ f Dp 2 jΩj=μ f and Clr is the slip rotational lift coefficient, which based on Oesterle and Bui Dinh's work for Rep b 2000 is [31].
C rl ¼ 0:45 þ
ReR 0:3 0:4 −0:45 e−0:05684 ReR Rep Rep
ð12Þ
In order to consider the torque acting on a rotating particle, Eq. (13) from Rubinow and Keller [30] is used: T¼
ρf 2
5 Dp C R jΩjΩ 2
ð13Þ
in which CR is the rotation coefficient and is computed by from Eq. (14): 8 > > <
64π ; ReR CR ¼ 12:9 128:4 > > ; : pffiffiffiffiffiffiffi þ ReR ReR
Rep ≤ 32 32 b Rep b 1000
Fig. 7. Normalized streamwise velocity profiles of gas-solid flow on vertical symmetry plane of the duct. (a) Duct inlet 0°, (b) 30°, (c) 45°, (d) duct outlet 90°.
ð14Þ
14
M.R. Banakermani et al. / Powder Technology 334 (2018) 9–26
Fig. 8. Computational grid used for validation study of erosion model with Vieira et al. experiment. Left: cross section view, Right: side view.
Table 2 Constants of Arabnejad's erosion model. Dp (μm)
C
K
εw(kg/ms2)
Utsh(m/s)
Pw(Pa)
n
fs
150 300
0.0125 0.0125
0.4 0.4
1.4 × 1011 1.4 × 1011
8 2
14.97 14.97
2.49 2.39
0.5 1
2.2.1. Particle-particle interaction Collisions between particles are modeled by the stochastic method developed by Sommerfeld [23]. In this approach, for every real particle, a
fictitious particle is generated at every Lagrangian time step, based on local diameter and velocity distribution of particles. After that, a collision probability is computed using the kinetic theory of gases, which is compared to a Gaussian random number. If the probability number is less than the random number, a collision will occur on a collision cylinder (Fig. 1) between the real and the fictitious particle and the post collision properties of the real particle are computed. It should be noted that in the simulations, restitution and friction coefficients for inter-particle collision calculations are assumed to be constant and equal to 0.9 and 0.4 for sand particles [32], respectively. For plastic particles, by assuming less elasticity, these values are set to 0.85 and 0.5, respectively.
Fig. 9. Calculated metal loss pattern corresponding to test #4 of Table 3.
Table 3 Comparison of Arabnejad et al. erosion model results with Vieira's UT Data and Vieira's CFD prediction in V-H configuration for flow direction in standard elbow. Test number
Gas velocity
Sand size
Sand rate
Max. UT data Vieira's data
CFD Max. Er. Rate pred. by Vieira
Pred. over data
CFD Max. Er. Rate pred. Arabnejad model
Pred. over data
#
(m/s)
(μm)
kg/day
mm/year
mm/year
(−)
mm/year
(−)
1 2 3 4 5 6 7 8 9
15 15 15 15 15 11 11 27 27
300 300 300 300 150 300 150 300 150
154 192 452 103 237 288 254 256 206
22.1 19.3 58.5 14.7 13.2 16.9 6.5 129.6 54
30.3 37.5 87.4 20.1 19.3 19.9 7.5 234.6 69.3
1.37 1.94 1.49 1.36 1.46 1.17 1.15 1.81 1.28
28.01 34.96 79.93 19.13 13.98 23.13 8.08 254.44 55.2
1.26 1.81 1.37 1.30 1.06 1.37 1.24 1.96 1.02
M.R. Banakermani et al. / Powder Technology 334 (2018) 9–26
15
Fig. 10. Side view of computational grid used to study erosion in elbows with different angles.
2.2.2. Particle-wall interaction Particle-wall interactions are simulated by considering wall roughness model proposed by Sommerfeld [22]. In this approach, the impact angle of particle with wall is composed of two parts. The first part is the angle derived from particle tracking, and the second part is calculated from a stochastic distribution. θ01 ¼ θ1 þ Δγξ
ð15Þ
In Eq. (15), ξ is a random number with zero mean and unit variance and Δγ is called the roughness angle parameter, which depends on the wall roughness structure and particle diameter [22].
2.3. Erosion model A recently developed erosion model by Arabnejad et al. [12] is used for computing the removed mass from target walls. Arabnejad et al.
_ sand ¼ 192 kg=day for V-H installation of elbows. (Note: different scale for each elbow angle). Fig. 11. Predicted metal loss pattern obtained at Ugas = 15 m/s, dp = 300 μm, m
16
M.R. Banakermani et al. / Powder Technology 334 (2018) 9–26
[12], proposed a semi-mechanistic erosion model, based on the results of direct impingement tests on several target materials. They assumed that material removal is due to two mechanisms, deformation and cutting. The volume swept by the cutting mechanism is computed from Eq. (16): 8 mp U p n sinθ½2Kcosθ−sinθ > > > ; θ b tan−1 K < 2K 2 P w Volc ¼ > mp U p n ðcosθÞ2 > > : ; θ ≥ tan−1 K 2P w
ð16Þ
In Eq. (16), Pw is the plastic flow stress of target wall material and is assumed to be equal to its Vickers hardness index, Hv (GPa). K is an empirical constant and n is the velocity exponent which is reported to be 2.41 for most of the tested materials [12]. mpand Up are particle's mass and particle's impact speed. The deformation damage removal is calculated from Eq. (17), developed by Bitter [9]:
VolD ¼
2 1 mp U p sinθ−U tsh 2 εw
Finally, total erosion ratio is obtained: kg CVolc þ VolD ER ¼ F s ρw kg mp
ð18Þ
where Fs in Eq. (18), is the sharpness factor, which is 1.0 for particles with sharp corners, and 0.25 for fully rounded particles, ρw is the wall material density and C is the cutting erosion coefficient, introduced to account for the difference between the swept volume predicted by the ideal impact and the actual impact. Therefore, the erosion rate can be found from Eq. (19):
Er ¼
Δh ERG_ sand ¼ t ρw Acell
ð19Þ
In Eq. (19), Δh is the erosion depth, t is the simulation time, ρwall is the density of wall material, Acell is the face area of the computational cell which is hit by particles and G_ is the rate of sand hitting the sand
ð17Þ
where in Eq. (17), Utsh is the threshold velocity, that is the minimum velocity above which deformation damage removal is not negligible and εw(kg/ms2) is the deformation wear factor, depending on the target wall material hardness.
cell. Changing the units to (mm/year), which is more applicable in pipeline industries is done by multiplying by a constant factor. Integrating the erosion rate over the entire surface, gives the total eroded volume over a year: ΔV r mm3 =year ¼
Z Er ðmm=yearÞdA
_ sand ¼ 192 kg=day for V-H installation of elbows. Fig. 12. Predicted metal loss pattern obtained at Ugas = 15 m/s, dp = 300 μm, m
ð20Þ
M.R. Banakermani et al. / Powder Technology 334 (2018) 9–26
3. Validation of gas-solid flow and erosion model Validation of gas-solid flow is performed in both a straight pipe and a curved duct geometry by comparing the simulation results with Tsuji et al.'s [20] and Kliafas et al.'s [21] experimental data. To validate the erosion model, Vieira et al.'s [19] experimental data are used.
3.1. Gas-solid flow validation in a straight pipe geometry For validation of the developed two-phase solver in OpenFOAM 2.3.1, Tsuji et al.'s [20] experimental data is used. In their work, LDV measurements are performed in a 30 mm horizontal pipe for standard air and plastic particles with an average diameter of 210 μm and density of 1000 kg/m3. To simulate the experiment, half of the pipe is modeled with respect to the vertical symmetry plane of the pipe. Cross section view of the computational grid used is shown in Fig. 2. In the first stage, to ensure fully developed flow, single-phase gas flow is solved in a three-meter long pipe with k − ε turbulence model. The obtained velocity, turbulence kinetic energy and turbulence energy dissipation profiles at the outlet and the calculated pressure gradient along the pipe are used as inlet condition for simulation of two-phase flow in a one-meter long pipe. The two-way coupling assumption is used here, and solid particles are injected from random positions at the inlet patch with size distribution given by Tsuji et al. [20]. Interparticle collisions are modeled using Sommerfeld's stochastic collision
17
method and particle-wall impacts are modeled by considering restitution coefficient of 0.9 and friction coefficient of 0.1. Fig. 3 compares non-dimensional mean velocity profiles of gas and solid phases on the vertical symmetry plane of the pipe, for two different particle mass loadings, m (kgsolid/kggas), with the experimental results of Tsuji et al. [20]. Fig. 4 shows the fluid turbulence intensity and Fig. 5 is a plot of the intensity of velocity fluctuations of particles, both compared to Tsuji et al.'s [20] measurements. 3.2. Gas-solid flow validation in a curved duct geometry Validation of gas-solid flow is performed also in a curved duct geometry, using Kliafas and Holt's [21] experimental data. They reported LDV measurements inside a square section 0.1m × 0.1m curved duct with R/D = 1.76 with air as the carrier fluid and sand particles as dispersed phase, where R is the mean radius of curvature and D is the hydraulic diameter of the cross section, which is equal to 100 mm. The computational grid used for the simulation is shown in Fig. 6. Air at standard conditions enters the duct with a uniform velocity of 52.19 m/s. Size and density of particles are 50 μm and 2990 kg/m3, respectively. Since the mass-loading of particles in this experiment was low, gassolid flow simulation is performed using one-way approach with k − ε turbulence model. Sommerfeld's stochastic model is used for particleparticle collisions and Sommerfeld's shadow model is used for simulating wall roughness with 5 ≤ Δγ ≤ 10. In Fig. 7, the normalized
_ sand ¼ 192 kg=day for V-H installation of elbows. Fig. 13. Contours of average particle impact speed at Ugas = 15 m/s, dp = 300 μm, m
18
M.R. Banakermani et al. / Powder Technology 334 (2018) 9–26
streamwise velocity is plotted against the normalized radial coordinate, r∗, which is defined in Eq. (21). ri and ro are radius of curvature of the inner and outer walls: r ¼
r−ri r o −r i
ð21Þ
3.3. Validation of erosion model To evaluate the accuracy of erosion models in estimating the amount and pattern of metal loss, Vieira et al.'s [19] experimental data are used. They measured metal loss at sixteen different locations of a 3 in. diameter (D = 76.2 mm) stainless steel standard 90°elbow with R/D = 1, by ultrasonic transducers mounted on its surface. Sand particles with diameters of 150 μm and 300 μm were used as erodent particles. The range of airflow velocity was from 11 m/s to 27 m/s and several sand mass loadings were tested. The computational grid used for validation study is shown in Fig. 8. RNG k-epsilon model to solve for turbulence, One-way coupling between the phases, stochastic collision between particles and wall roughness angle parameter of 5 ≤ Δγ ≤ 10 are considered in the simulations. Furthermore, restitution coefficient of 0.8 and friction coefficient of 0.2 were found to give the best fit to the experimental results. The constants
of erosion model for the pipe wall material and properties of particles are set according to Table 2. In OpenFOAM every computational parcel is a set of particles which is represented by nParticle. OpenFOAM's patchInjection model is used to inject parcels to the computational domain. To determine the amount of mass (or the number of parcels) being injected in each time step, duration of injection, mass of parcels added to the system during injection (which is set according to sand flow rates) and number of parcels added in every second should be specified. In each second, using OpenFOAM's patchInjection model, nearly 5 × 105 computational parcels were injected from randomly-changing positions on the inlet patch of the pipe with a uniform velocity equal to the bulk flow velocity for each test. This resulted in nParticle of approximately 1.06 in the test cases. Injections were continued until the amount of erosion rate did not change with time. In Fig. 9, the metal loss pattern obtained in one of the test cases is illustrated. Other simulation cases with different inlet conditions also had similar patterns and only different erosion rate magnitudes were observed. By inspecting the erosion profile on the outer wall of the elbow, the maximum metal loss occurs at around 55°curvature angle (V-H configuration). The angle is measured from inlet side to the flow outlet, such that 0° is the inlet of the elbow and 90°corresponds to the outlet. In Table 3, maximum metal losses obtained from the current simulations are compared to the numerical simulations and experimental measurements of Vieira et al. [19].
_ sand ¼ 192 kg=day for V-H installation of elbows. Fig. 14. Contours of average particle impact angle at Ugas = 15 m/s, dp = 300 μm, m
M.R. Banakermani et al. / Powder Technology 334 (2018) 9–26
_ sand ¼ 192 kg=day for V-H installation of elbows. Fig. 15. Contours of particle impact frequency at Ugas = 15 m/s, dp = 300 μm, m
_ sand ¼ 192 kg=day for H-H installation of elbows. Air-sand flow goes from left to right. Fig. 16. Predicted metal loss pattern obtained at Ugas = 15 m/s, dp = 300 μm, m
19
20
M.R. Banakermani et al. / Powder Technology 334 (2018) 9–26
_ sand ¼ 192 kg=day for H-H installation of elbows. Air-sand flow goes from left to right. Fig. 17. Contours of average particle impact speed at Ugas = 15 m/s, dp = 300 μm, m
As can be seen from Table 3, a relatively close agreement is obtained between current numerical predictions and the experimental data. 4. Erosion in elbows with different angles In order to compare erosion results for different elbow angles, several elbows with angles from 15°to 90° are numerically simulated at certain fixed inlet gas-solid flow conditions, equal to Vieira et al.'s [19]
experimental conditions. Both vertical inlet-horizontal outlet (V-H) and completely horizontal (H-H) flow directions are modeled for each elbow angle with inlet conditions 1 to 4 given in Table 3. Throughout this paper, where both inlet and outlet flow directions are horizontal, the configuration is called H-H arrangement. For the cases in which inlet direction is vertical and outlet direction either lies in the horizontal plane or makes different angles with it, the configuration is called V-H arrangement. The cross section of computational grids is the same as
_ sand ¼ 192 kg=day for H-H installation of elbows. Air-sand flow goes from left to right. Fig. 18. Contours of average particle impact angle at Ugas = 15 m/s, dp = 300 μm, m
M.R. Banakermani et al. / Powder Technology 334 (2018) 9–26
21
_ sand ¼ 192 kg=day for H-H installation of elbows. Air-sand flow goes from left to right. Fig. 19. Contours of particle impact frequency at Ugas = 15 m/s, dp = 300 μm, m
in Fig. 8 and Fig. 10 show the side view of computational grids for each elbow angle. The grid cells are hexahedral for better numerical stability and accuracy [14].
All settings used for each case are completely identical to the settings used for validation of 90° elbow in Vieira's experiment (Table 3, Tests #1 to #4), and only the geometry and in H-H cases, the flow direction is
Fig. 20. Erosion rate (a), impact angle (b), impact frequency (c) and impact speed (d), presented at the symmetry plane of the elbow from the bend inlet to the bend outlet for V-H flow installation.
22
M.R. Banakermani et al. / Powder Technology 334 (2018) 9–26
Fig. 21. Erosion rate (a), impact angle (b), impact frequency (c) and impact speed (d), presented at the symmetry plane of the outlet pipe from the bend outlet to 1 m further downstream for V-H flow installation.
Fig. 22. Erosion rate (a), impact angle (b), impact frequency (c) and impact speed (d), presented close to the symmetry plane of the elbow from the bend inlet to the bend outlet for H-H flow installation.
M.R. Banakermani et al. / Powder Technology 334 (2018) 9–26
23
Fig. 23. Erosion rate (a), impact angle (b), impact frequency (c) and impact speed (d), presented close to the symmetry plane of the outlet pipe from the bend outlet to 1 m further downstream for H-H flow installation.
changed. Fig. 11 shows the contours of metal loss in V-H flow direction from the side view, for one of the inlet flow conditions. In Fig. 12, to make the comparison easier, the contours of Fig. 11 are plotted in a different view and one scale is selected for all of the bend angles. Also to illustrate the process of material removal in different cases, in the same manner as Fig. 12, contours of important variables regarding the
interaction of particles and wall, including particle's impact speed, impact angle and impact frequency are represented in Figs. 13–15 respectively. Similarly, in Figs. 16–19, erosion rate and the above mentioned variables are depicted for the H-H flow direction. It should be noted that the impact speed and impact angle contours are obtained by taking the average on all the particles that impact each cell face. Since different
Fig. 24. Maximum erosion rate for different elbow angles and sand mass loadings in V-H and H-H elbow configurations.
24
M.R. Banakermani et al. / Powder Technology 334 (2018) 9–26
Fig. 25. Total eroded volume for different elbow angles and sand mass loadings in V-H and H-H elbow configurations.
geometries with different wall's cell face areas are studied, the reported impact frequency is defined as the number of particles hit the wall, per unit time per unit area. Considering Figs. 12, 16, for both V-H and H-H configurations, the location of maximum erosion is very close to the symmetry plane of the elbows. Therefore, to examine more carefully the variation of erosion rate and the influence of impact-related variables on it, these values are plotted on the outer wall of the elbows, along the symmetry plane for V-H configurations, and very close to it for H-H ones. The data are extracted both from the bend portion, starting from 0° to the bend angle (Figs. 20, 22) and from the outlet pipe's wall, starting from the bend outlet to 1 m further downstream (Figs. 21, 23). For more quantitative details, the maximum erosion rate and total annual eroded volume, which are two important parameters in the analysis of surface erosion, are plotted in Figs. 24, 25 respectively for each elbow angle. 5. Results and discussion As can be seen in Figs. 12, 16, an elliptic pattern for metal loss is formed in all elbow angles, but by increasing elbow angle, an additional “V” shape pattern also appears. The “V” shape pattern becomes more visible by increasing the elbow angle from 60° to 90°. However, Solnordal et al. [18] showed that for V-H arrangement of 90° elbows
the “V” shape scars are not present in reality. Therefore, appearance of the “V” shape scars specially in V-H configurations, might be an artefact of modeling, although because of careful modeling of particle-particle collisions and wall effects, the scars in this work are not as evident as in most other author's simulations. In addition, these scars in H-H arrangements are much more evident than the ones in V-H arrangements. This could be due to the effect of gravity on trajectories of particles that causes a focusing of secondary particle-wall impacts in these regions. It is clear in Fig. 19 that the upper branch of “V” (the secondary impact region) is struck more by the particles than the lower one, which indicates that after the primary impact on the lower wall, particles move upward along the upper branch of the produced “V” scar. By comparing Fig. 12 with Fig. 15 and Fig. 16 with Fig. 19, it is evident that regardless of flow orientation and bend angle, the pattern of erosion rate is similar to the pattern of particle impact frequency. So this variable seems to be the most important one in determining material removal pattern in both H-H and V-H flows. Considering Figs. 20, 22, it is clear that erosion rate and other impactrelated variables follow approximately a similar trend from the start to the end of the bend for all the elbow angles. Additionally, it is found from the results that the point of maximum erosion rate on the elbow's wall is around 55° in the V-H flow (Fig. 20a) and around 50° in the H-H flow (Fig. 22a). Further investigating the erosion rate and other variables along the outlet pipe's wall, in Figs. 21 and 23, shows that the
Fig. 26. Normalized maximum erosion rate in elbows with bend angles from 15° to 90°, for different sand mass flow rates: (1) 154 kg/day, (2) 192 kg/day, (3) 452 kg/day, (4) 103 kg/day.
M.R. Banakermani et al. / Powder Technology 334 (2018) 9–26
values along the pipe's length are spread out. Besides, it is understood from Figs. 21a, 23a that how much the outlet pipe is being influenced by erosion. Both in V-H and H-H flows, for 75° and 90° elbows, the maximum erosion rate on the outlet pipe's wall is negligible in comparison with other bend angles. Considering the outlet pipe's wall, in V-H flow from 60° to 15°, the amount of maximum erosion rate decreases and the location of maximum moves further downstream along the flow direction. In H-H flow, a similar trend is observed but the point of maximum erosion rate for 60° elbow seems to be moved slightly upstream (on the bend portion) and 45° and 30° elbows have higher maximum erosion rates than 60°elbow on the outlet pipe's wall. Examining the area that is struck by the particles, in geometries with lower bend angles, more surface area is affected by cutting action. Due to the ductility of wall material in this study, cutting removal is more important than the deformation mechanism. Therefore, a larger area is affected by erosion for cases with smaller bend angles. However, according to Eq. (16), the amount of metal loss depends on the particle's impact angle and as the elbow angle decreases, the particles hit the wall with smaller angles, creating small erosion rates. It can be said that flow orientation is influential on both the maximum erosion rate and the total removed volume. From Fig. 24, it is evident that generally for the V-H flow orientations the maximum erosion rates in every elbow angle are higher than the ones for the H-H flow orientation. The reason for such behavior can be associated to the effect of gravity on the distribution of particles. From comparing Fig. 13 to Fig. 17 (also comparing Fig. 20d to Fig. 22), it can be seen that compared to the H-H arrangement, in the V-H flow configurations particles hit the wall with higher speeds. This is due to the effect of gravity that in H-H flow configurations causes more particles to shift to near the wall, where the fluid velocity is lower. But in the V-H flow configurations, the particles mainly travel near the center of the pipe where the velocities are higher. Comparing Fig. 15 to Fig. 19, also clarifies that in the V-H flow, the elliptic main eroded area, is being hit with a higher impact frequency than H-H flow. In H-H configurations, gravity moves the particles downward and they hit a larger area that causes lower impact frequency (or particle focusing) at the elliptic main eroded area, which is around the maximum erosion location. On the other hand, the total annual eroded volume is slightly higher in the H-H flow orientations, which means that more surface is influenced by the erosion phenomena and integration over a larger area would give a slightly more overall volume removal rate. In Figs. 26, 27, for every value of sand rate, the maximum erosion rate and the total annual eroded volume for every elbow angle are normalized by the value of the maximum erosion rate and the total annual eroded volume in the 90°elbow, respectively. As it is clear from Fig. 26, all the curves virtually coincide with each other, which means that regardless of the value of sand rate, the maximum erosion rate for every elbow angle is
25
approximately a particular fraction of the maximum erosion rate in 90°elbow. On the other hand, as can be seen in Fig. 27, in comparison to Fig. 26 the variation is very slight. This is because the rate of total eroded volume is proportional to the kinetic energy of the particles entering the pipe. Since the sand flow rate and sand mass loading is the same for all the bend angles in each test condition, the total removed volume over a year remains almost constant for both flow configurations and the curves follow a similar trend. However, in H-H arrangements, generally the total annual eroded volume has a maximum for the 30° elbow. 6. Conclusion In this work, a two-phase particulate flow solver is developed in OpenFOAM 2.3.1 which combines Eulerian fluid flow simulation with Lagrangian particle tracking. In the first stage, the gas-solid solver is validated with available experimental LDV data in both a straight pipe and a curved duct geometry. Then, a recently developed semi-mechanistic erosion model is validated by comparing the simulation results with the Vieira et al.'s [19] experimental UT data in a 90°standard 76.2 mm diameter elbow. After which, the same elbow diameter and flow conditions are considered for analyzing erosion in elbows with angles ranging from 15°to 90° in vertical to horizontal (V-H) and horizontal to horizontal (H-H) flow arrangements. Erosion rate and important impactrelated variables such as particle impact speed, impact frequency and impact angle were investigated for these cases. It was observed that the pattern of metal loss and particle-wall impact frequency are very similar and the location of maximum erosion rate on elbow's wall is about 55° for the V-H configurations and about 50° for the H-H ones. Furthermore, it was found that in elbows with angles from 15° to 45°, only an elliptic metal loss pattern is created, but by increasing the elbow angle, a “V” shape scar appears also, which intensifies as the elbow angle is increased. This scar, which its intensity in the current work is low compared to most of other author's simulations, especially for the V-H arrangements is speculated to be an artefact of modeling, but for H-H arrangements secondary particle-wall impacts also could be effective. Comparing the results, it is evident that the maximum erosion rate for the V-H flow orientation is always larger than that for the H-H arrangement. However, integrating the erosion rates over the total wall's area, the rate of total volume removal is the opposite. By normalizing the obtained results for maximum erosion rate with the values of it in 90° elbow, it was observed that regardless of the value of sand flow rate, the maximum erosion rate in a specific elbow angle could be expressed as a particular fraction of the value corresponding to 90° elbow and this fraction is a function of elbow angle. Doing the same for the rate of total volume removal gives an approximately constant value, which is nearly independent of the elbow angle.
Fig. 27. Normalized total eroded volume in elbows with bend angles from 15°to 90°, for different sand mass flow rates: (1) 154 kg/day, (2) 192 kg/day, (3) 452 kg/day, (4) 103 kg/day.
26
M.R. Banakermani et al. / Powder Technology 334 (2018) 9–26
List of symbols Symbol
Units (SI) 2
References Description
Acell C CD Clr Cls CR Dp Er ER
m − − − − − m m/s
Fg FRL Fs FSL g G_
N N N N m/s2 kg/s
Sum of gravity and buoyancy forces, see Eq. (5) Slip-rotational lift force, see Eq. (11) Particle sharpness (shape) factor, see Eq. (18) Saffman lift force, see Eq. (9) Gravitational acceleration of earth Rate of sand hitting a cell of the wall, see Eq. (19)
Δh Hv Ip k K
m GPa kgm2 m2/s2 −
Erosion depth, see Eq. (19) Vickers hardness Particle moment of inertia Turbulence kinetic energy, see Table 1. Ratio of vertical to horizontal projected contact area, see Eq. (16) Particle mass loading
sand
kg wall kg solid
Face area of a computational cell on the wall, see Eq. (19) Cutting erosion factor, see Eq. (18) Coefficient of drag force, see Eq. (7) Slip-rotational lift coefficient, see Eq. (12) coefficient of slip-shear lift force, see Eq. (10) Rotation coefficient, see Eq. (14) Particle diameter Erosion rate, see Eq. (19) Erosion ratio, see Eq. (18)
m
kg solid kg gas
mp _ sand m n p Pw r r∗ R ri ro Rep ReR
kg kg/s − Pa Pa m − m m m − −
Particle mass Sand mass flow rate through the elbows Ratio of contact area to removed area, see Eq. (16) Fluid pressure, see Eq. (2) Wall material flow pressure, see Eq. (16) Elbow radial coordinate direction Normalized radial coordinate direction, see Eq. (21) Elbow mean radius of curvature Radius of curvature of the inner wall of elbow, see Eq. (16) Radius of curvature of the outer wall of elbow, see Eq. (16) Particle Reynolds number, see Eq. (8)
Res
−
Shear flow Reynolds number, Res ¼ ρ f Dp 2 jω f j=μ f
Suip
kg ms2
t T Uf Up Utsh VolC VolD ∇Vr
s N. m m/s m/s m/s m3 m3 m3/s
Source term of particles in fluid momentum equation, see Eq. (2) Time Torque acting on a rotating particle, see Eq. (13) Carrier fluid velocity Particle velocity Threshold velocity of deformation erosion, see Eq. (17) Removed volume by cutting mechanism, see Eq. (16) Removed volume by deformation mechanism, see Eq. (17) Rate of total volume removal, see Eq. (20)
Greek symbols β − Δγ Deg. ε m2/s3 kg εw 2 θ θ1′
ms
Deg. Deg.
μ μt ξ
Pa. s Pa. s −
ρf ρp ρw φ ωf ωp Ω
kg/m3 kg/m3 kg/m3 − 1/s 1/s 1/s
Rotational Reynolds number, ReR ¼ ρ f Dp 2 jΩj=μ f
Ratio of Reynolds numbers, β = 0.5Res/Rep, see Eq. (10) Wall roughness angle parameter, see Eq. (15) Rate of dissipation of turbulence energy, see Table 1 Deformation wear factor Particle-wall impact angle Particle-wall impact angle obtained by wall roughness modeling, see Eq. (15) Carrier fluid viscosity Fluid turbulence viscosity, see Eq. (2) Random number with zero mean and unit variance, see Eq. (15) Carrier fluid density Particle's density Wall material density The variable being solved, see Table 1 Carrier fluid vorticity, ωf = ∇ × Uf Particle angular velocity Difference between the fluid and particle rotational speed, Ω = 0.5∇ × Uf − ωp.
[1] A. Sherik, S. Zaidi, E. Tuzan, J. Perez, Black Powder in Gas Transmission Systems, CORROSION 2008, NACE International, 2008. [2] M. Saremi, M. Kazemi, The effect of black powder composition on the erosion of compressor's impeller in gas transmission line, Adv. Mater. Res. (2011) 1514–1518. [3] M. Parsi, K. Najmi, F. Najafifard, S. Hassani, B.S. McLaury, S.A. Shirazi, A comprehensive review of solid particle erosion modeling for oil and gas wells and pipelines applications, J. Nat. Gas Sci. Eng. 21 (2014) 850–873. [4] I. Finnie, The mechanism of erosion of ductile metals, 3rd US National Congress of Applied Mechanics, 1958. [5] I. Finnie, D. McFadden, On the velocity dependence of the erosion of ductile metals by solid particles at low angles of incidence, Wear 48 (1978) 181–190. [6] I. Hutchings, R. Winter, J. Field, Solid particle erosion of metals: the removal of surface material by spherical projectiles, Proceedings of the Royal Society of London A: Mathematical, Physical and Engineering Sciences, The Royal Society 1976, pp. 379–392. [7] I. Hutchings, A model for the erosion of metals by spherical particles at normal incidence, Wear 70 (1981) 269–281. [8] J. Bitter, A study of erosion phenomena: part II, Wear 6 (1963) 169–190. [9] J. Bitter, A study of erosion phenomena part I, Wear 6 (1963) 5–21. [10] C. Huang, S. Chiovelli, P. Minev, J. Luo, K. Nandakumar, A comprehensive phenomenological model for erosion of materials in jet flow, Powder Technol. 187 (2008) 273–279. [11] M. Jafari, Z. Mansoori, M.S. Avval, G. Ahmadi, A. Ebadi, Modeling and numerical investigation of erosion rate for turbulent two-phase gas–solid flow in horizontal pipes, Powder Technol. 267 (2014) 362–370. [12] H. Arabnejad, A. Mansouri, S. Shirazi, B. McLaury, Development of mechanistic erosion equation for solid particles, Wear 332 (2015) 1044–1050. [13] X. Chen, B.S. McLaury, S.A. Shirazi, Application and experimental validation of a computational fluid dynamics (CFD)-based erosion prediction model in elbows and plugged tees, Comput. Fluids 33 (2004) 1251–1272. [14] G.C. Pereira, F.J. de Souza, D.A. de Moro Martins, Numerical prediction of the erosion due to particles in elbows, Powder Technol. 261 (2014) 105–117. [15] Y.I. Oka, K. Okamura, T. Yoshida, Practical estimation of erosion damage caused by solid particle impact: part 1: effects of impact parameters on a predictive equation, Wear 259 (2005) 95–101. [16] Y. Oka, T. Yoshida, Practical estimation of erosion damage caused by solid particle impact: part 2: mechanical properties of materials directly associated with erosion damage, Wear 259 (2005) 102–109. [17] C.A.R. Duarte, F.J. de Souza, V.F. dos Santos, Numerical investigation of mass loading effects on elbow erosion, Powder Technol. 283 (2015) 593–606. [18] C.B. Solnordal, C.Y. Wong, J. Boulanger, An experimental and numerical analysis of erosion caused by sand pneumatically conveyed through a standard pipe elbow, Wear 336 (2015) 43–57. [19] R.E. Vieira, A. Mansouri, B.S. McLaury, S.A. Shirazi, Experimental and computational study of erosion in elbows due to sand particles in air flow, Powder Technol. 288 (2016) 339–353. [20] Y. Tsuji, Y. Morikawa, LDV measurements of an air–-solid two-phase flow in a horizontal pipe, J. Fluid Mech. 120 (1982) 385–409. [21] Y. Kliafas, M. Holt, LDV measurements of a turbulent air-solid two-phase flow in a 90 bend, Exp. Fluids 5 (1987) 73–85. [22] M. Sommerfeld, N. Huber, Experimental analysis and modelling of particle-wall collisions, Int. J. Multiphase Flow 25 (1999) 1457–1489. [23] M. Sommerfeld, Validation of a stochastic Lagrangian modelling approach for interparticle collisions in homogeneous isotropic turbulence, Int. J. Multiphase Flow 27 (2001) 1829–1858. [24] S. Moissette, B. Oesterlé, P. Boulet, Temperature fluctuations of discrete particles in a homogeneous turbulent flow: a Lagrangian model, Int. J. Heat Fluid Flow 22 (2001) 220–226. [25] E. Amani, M. Nobari, Systematic tuning of dispersion models for simulation of evaporating sprays, Int. J. Multiphase Flow 48 (2013) 11–31. [26] L. Schiller, Z. Naumann, A drag coefficient correlation, Ver. Dtsch. Ing. (1935) 77–318. [27] P. Saffman, The lift on a small sphere in a slow shear flow, J. Fluid Mech. 22 (1965) 385–400. [28] R. Mei, An approximate expression for the shear lift force on a spherical particle at finite Reynolds number, Int. J. Multiphase Flow 18 (1992) 145–147. [29] C. Crowe, M. Sommerfeld, Y. Tsuji, Fundamentals of Gas-Particle and Gas-Droplet Flows, CRC press Boca Raton, 1998. [30] S. Rubinow, J.B. Keller, The transverse force on a spinning sphere moving in a viscous fluid, J. Fluid Mech. 11 (1961) 447–459. [31] B. Oesterle, T.B. Dinh, Experiments on the lift of a spinning sphere in a range of intermediate Reynolds numbers, Exp. Fluids 25 (1998) 16–22. [32] C.A.R. Duarte, F.J. de Souza, R. de Vasconcelos Salvo, V.F. dos Santos, The role of interparticle collisions on elbow erosion, Int. J. Multiphase Flow 89 (2017) 1–22.