#### Abstract

The computational fluid dynamics modeling of solid particles hydrodynamic based on the Lagrangian framework for diluted solid-gas flow through 90° gas pipeline bend is carried out to discover the effect of particles size distribution on particles flow pattern and their erosive effect on the bend. Particles size distribution has been obtained experimentally by measuring the sizes of solid particles that are flowing through the gas pipelines of Aghajari gas booster station. Also the erosion rate at the outer wall of the bend is predicted. The pipeline bend under study has a pipe diameter of 56 inches and ratios of the bend radius of the curvature to the pipeline diameter of 1.5. For the validation of computational model, firstly, the computational modeling is performed for a published experimental solid-gas flow data. The computational results include radial gas velocity and radial particle velocity profiles on planes which are at different angles through the bend. The comparison between the predicted numerical results and similar experimental data proves that the predictions of the computational model are acceptable. Finally, the particles' size distributions on each plane through the bend and the erosion rate on the outer wall of the bend have been obtained. The maximum rate of erosion is found to be 3.2 nm/s, occurring between 40 and 65° of the bend.

#### 1. Introduction

In the oil and gas industry, black powder (BP) is the brief name that is used to describe the black materials found inside most of the gas pipelines worldwide. Black powder can be found in several forms, such as wet with a tar-like appearance or dry in the form of a very fine powder [1–5]. It is composed of different forms of iron sulfide (FeS), iron oxides (Fe_{3}O_{4}, FeOOH), and iron carbonate (FeCO_{3}), mechanically mixed or chemically combined with any number of contaminants, such as salts, sand, liquid hydrocarbons, and metal debris [2]. Once BP exists and is moving with the flow, it can represent a serious threat to the integrity of the gas pipelines by eroding compressor components and pipeline control valves, plugging metering instrumentation and filters and reducing the accuracy of the in-line inspection. Also, BP could have major adverse effects on customers by contaminating the customers’ sales gas supply leading to interruptions of the customers’ operations and/or poor quality of products in which the sales gas is used as feedstock [3].

The required fluid velocity has been determined [6, 7] to entrain and carry away BP in liquid and gas pipelines, respectively. These two studies concluded that the velocity required to move BP particles in gas pipelines is independent of particle size and ranges from 10.4 ft per second (fps) to 13.6 fps for 8′′ and 30′′ pipelines, respectively. In liquid pipelines, the water velocity required depends on the equivalent particle size, up to a size of about 5.0 millimeters, after which it depends only on the pipe diameter.

The effect of the drag coefficient and inlet conditions (inlet velocity profile) of solid particles on the particle tracks calculations in vertical and horizontal ducts are studied [8] using the commercial computational fluid dynamics (CFDs) package, CFX 4.4. They found that the drag coefficient needs to be reduced by as much as 35% of the standard value to achieve good agreement with the corresponding experimental data in case of a vertical channel flow. On the other hand, for a horizontal channel flow, it needs to be reduced only by 20% to achieve similar agreement. Regarding the velocity inlet conditions, it was reported [8] that the vertical turbulent flow seems to be insensitive to the inlet conditions, while for a horizontal flow, it is found to be strongly dependent on inlet conditions.

CFD simulations have been performed [9] on a diluted particulate turbulent flow in a 90° duct bend with a radius of curvature equal to a 1.5 duct (225 mm) hydraulic diameter. As in previous works [8], simulations were performed using CFX 4.4, using the differential Reynolds stress model (DRSM) with fully developed inlet conditions to solve the turbulent flow in the bend and also used the same test facility to produce the experimental data used in validating the simulations. In another work [10], the author used different solid size distributions rather than a single uniform particle size and also made use of a modified shear-slip lift force formula, which is consistent with experimental data. From these studies [8, 9], it was concluded that the DRSM did not capture the correct pressure gradient effects within the bend. Also, it was found that even the finer particles (66 micron) experienced a gas-solid segregation due to the centrifugal effect. This segregation was characterized by a local drop in particle concentration near the inner wall and was well reflected in predictions where the averaged velocity profiles discontinued in the locality. The experimental part of the study [10] is reported in more details [11].

CFD-based erosion modeling can be applied to predict erosion in many complex geometries. To assess the viability and accuracy, a comparison between computed and measured particle velocities and erosion, in both water and air flows, in a direct impact test section, was performed [12]. It was found that for a sand/water flow in a direct impact test, not like sand/air, the particle impact velocity is much lower than the velocity of the slurry jet and varies by a wide range. Also, they found that, among the erosion models tested, the Erosion-Corrosion Research Center erosion model and Oka’s et al. erosion model [13, 14] are more accurate within the scope of their work. To evaluate the performance of elbows and plugged tees geometries under erosive service conditions and using the experimental data [15] to validate the simulation results, a procedure was developed [16] to predict erosion in standard elbows, long-radius elbows, and plugged tees. This procedure is implemented into the CFD code CFX 4.2. The relative erosion severity between plugged tees and elbows under diluted gas/liquid-solid flow conditions has been studied, computationally and experimentally. In this study [17], it was shown that the relative erosion severity is greatly affected by the type of carrier fluid (liquid or gas) properties, as simulations showed that water-sand flows in plugged tees cause more erosion than in elbows, while air-sand flows erosion in plugged tees is found to be two orders of magnitude less than the erosion in standard elbows. In another study [18], it was shown that the longer the radius of curvature of the elbow, the less the erosion it experiences due to solid particle impact. Not only the radius of curvature that affects the erosion but also the elbow (bend) orientation was found to have a large effect on particle motion and, therefore, on erosion rates as shown in another study [19].

In the current study, the reported problem in Aghajari gas booster station, consisting of plugging of the compressor’s filters and the erosion of facilities in gas pipelines due to particles existence, will be investigated from particles hydrodynamic point of view.

Solid particulates that are flowing inside gas pipelines (BP) of Aghajari gas booster station have been analyzed. For evaluating the effect of particles size on particles motion and their erosive influence, CFD modeling based on the Lagrangian framework is performed for a 90° gas pipeline bend. Particles size distribution is considered in the modeling by Rosin-Rammler distribution function.

#### 2. Geometry and Flow Conditions

The considered geometry of gas pipeline is a 90° angled bend with ratios of the bend radius of the curvature to the pipeline diameter of 1.5, which is used for detailed modeling of particles motion which is associated with particles size distribution. This CFD modeling is performed based on the Lagrangian framework.

Particles size distribution has been obtained experimentally by measuring the size and the relevant mass of solid particles that are flowing through the gas pipelines of Aghajari gas booster station by the use of woven wire test sieve (WWTS). Particles size and mass distribution are given in Table 1.

Particles are collected at the sampling point of 56-inche-diameter pipe, which is the primary inlet pipeline to Aghajari gas station. After 500 hours, 300 kg of particles was obtained indicating that the mass flow rate of particles is 0.6 kg/hr. The measured density of particles is 2303 kg/m^{3}. The stream of main inlet pipe is distributed between seven compressors from which one of them works at the normal condition. Therefore, the gas flow rate at normal condition which is used for the modeling purpose is 600 SMMCF/H. The pressure and temperature of supplied gas in main inlet pipe are 80 Barg and 40°C, respectively. The physical conditions of flowing gas in the 56-inche-diameter pipe is given in Table 2.

#### 3. Mathematical Model

The commercial CFD software FLUENT 6.3 is used to solve the Reynolds-averaged Navier-Stokes (RANS) equations for continuous gas phase. For considering the effect of particles size distribution on particles motion and particles trajectory, the Lagrangian framework for modeling diluted solid-gas flow is used. Flowing particles in the main gas pipeline of Aghajari station have been gathered and analyzed by woven wire test sieve to determine the size and mass distribution of particles (Table 1). The Rosin-Rammler distribution function is used to specify the fraction of particles with specific sizes. The mass fraction of particles of diameter greater than is given by where is the size constant and is the size distribution parameter.

Particles are considered as solid spheres, which are injected in to the computational domain via surface injection model. In this model, there are eight cells at the inlet boundary so that ten particles are injected from each one at every injection time, at velocity of 20 (m·s^{−1}) and injection angle of zero relative to the surface normal. According to the accomplished measurement, the total mass flow rate of particles is 0.00017 (kg·s^{−1}). It is considered that when a particle collides a wall surface, it retains all of its normal or tangential momentum after the rebound (an elastic collision).

The trajectory of a discrete phase particle is predicted by integrating the force balance on the particle. This force balance equates the particle inertia with the forces acting on the particle and can be written (for the direction in Cartesian coordinates) as follows: The drag force imposed on the particles is given by (15): where The particles path is computed by integrating (5) as follows: The dispersion of particles due to gas phase turbulence is calculated by the Discrete Random Walk Model.

Turbulent dispersion of particles is predicted by integrating the trajectory equations for individual particles, using the instantaneous fluid velocity along the particle path during the integration as follows: Prediction of particle dispersion makes use of the concept of the integral time scale, , which describes the time spent in turbulent motion along the particle path. For small particles that move with the fluid, the integral time becomes the fluid Lagrangian integral time, . This time scale can be approximated as Each eddy is characterized by a Gaussian distributed random velocity fluctuation, , , and a time scale, where is a normally distributed random number, and the remainder of the right-hand side is the local RMS value of the velocity fluctuations.

Since the kinetic energy of turbulence is known at each point in the flow, these values of the RMS fluctuating components can be defined (assuming isotropic flow) as The characteristic lifetime of the eddy is defined as where is a uniform random number between 0 and 1, and is given by (8).

The particle eddy crossing time is defined as where , , and are particle relaxation times, eddy length scale, and the magnitude of relative velocity, respectively. The particle is assumed to interact with the fluid phase eddy over the smaller of the eddy lifetime and the eddy crossing time. When this time is reached, a new value of the instantaneous velocity is obtained by applying a new value of in (9).

Due to the high gas velocity (33 m/s) and high strain rate of fluid near the pipe wall, the realizable model is used for modeling gas phase turbulence. Consider where , are turbulent Prandtl numbers and is turbulent viscosity.

and are turbulent kinetic energy and dissipation rate, respectively. Consider Momentum balance equation for gas phase is The erosion rate at wall boundaries can be evaluated by a new combination of the Tulsa angle dependent model with Huser and Kvernvold model [20, 21]: where , , , , and are the erosion rate, Brinell hardness, particle shape coefficient, particle relative velocity, and a 5th-order polynomial function of impact angle, respectively. The impact angle function is obtained by a fitting operation on the experimental data, as one can see in Figure 1.

**(a)**

**(b)**

The erosion rate is computed assuming reflecting wall boundary condition. For steel material, and (for fully rounded solid particles).

#### 4. Boundary Conditions

According to Table 2, the inlet mass flow rate of gas is set to 3780 (kg·s^{−1}). Solid particles are injected into the domain through eight cells at the inlet boundary at total mass flow rate of 0.00017 (kg·s^{−1}) with particle size distribution that is obtained from (1). Particles are injected normally, to the inlet boundary at velocity of 20 (m·s^{−1}), and their collision with wall boundary is assumed to be elastic. At the outlet boundary, gas static pressure is set to 80 (barg), and particles escape from outlet boundary. particle tracking is accomplished in unsteady mode.

The boundary layer mesh generation is considered at the wall to limit value.

#### 5. Validation of the Mathematical Model

Due to impossibility of measurements on working pipe, the published experimental data for a diluted gas-solid flow through a curved 90° duct bend [11] was used to validate the mathematical model based on the Lagrangian framework. The curved bend is a squared section (15 cm × 15 cm) and has a radius of curvature, , of 1.5 times the duct hydraulic diameter, , (22.5 cm). Gas phase measurements were obtained using a Laser Doppler Anemometer (LDA) at a bulk gas velocity, , of 10 m/s in the absence of solid phase. The solid phase, which is glass spheres with an average diameter of 66 *μ*m, was released into the flow from a fluidized bed. The solids/gas mass loading ratio reached is well below 1%, so as to setup a diluted gas-solid flow regime. The radial velocity profiles of gas and particles are compared with similar measurement data that is obtained from different cross sectional planes through the squared bend (Figure 2).

In Figure 3, the predicted radial distribution of gas velocity is compared with experimental data.

Radial distance, , is computed by (19) where is the curve radius of duct, is the hydraulic diameter of duct, and is the distance of any point on a special cross sectional plane from the origin.

The total number of cells that is generated for the computational domain is 565000. The results of mesh independency check are depicted in Figure 3. In this figure, the effect of cell numbers on gas velocity distribution in radial direction over the 15° cross sectional plane is considered.

As it is shown in Figure 4, by increasing the cross sectional planes angle, the more conformity between predicted profiles and measured profiles is achieved. This can be due to decreasing the radial component of gas velocity.

**(a)**

**(b)**

**(c)**

In Figure 5, the predicted radial distribution of particles velocity is compared with experimental measurements. As can be seen, predicted results show that particles velocity profile does not continue to inner wall. This is due to the radial component of gas velocity that leads to moving the particles toward the outer wall of the bend. At each angle in the pipe, we use a plane that has been divided radially to 20 sections. At each section, the velocity of particles is gathered over time, and the average of values is calculated for the particle velocity at the specified radial position. This procedure of particle velocity calculation may lead to discrepancy between the measured and calculated values. It is obvious that by increasing the number of sections on each plane, the discrepancy between the measured and calculated particle velocities can be decreased.

**(a)**

**(b)**

**(c)**

Figure 6 shows the radial velocity vectors of gas inside the bend.

An erosion model according to (18) was inserted into CFD solver by a user-defined-function. The model was capable of predicting the erosion rate so that the maximum erosion rate was obtained as m/s, which is in good agreement with that of m/s [22]. The erosion rate at the middle of bend, along the symmetry plane, where the maximum erosion takes place, was previously reported [22].

#### 6. Results and Discussion

In this study, the particulates flow with particles size distribution inside a 90° angled bend with ratio of curve radius to pipe diameter of 1.5 and 56 inches of pipe diameter, is considered based on the Lagrangian framework. The particles size distribution is obtained from experimental data (Table 1) and taken into account by Rosin-Rammler distribution function. The mass fraction of particles of diameter greater than is given by . Table 3 explains the relationship between and according to Table 1*.*

The mass flow rate of particles is 0.6 kg/hr which is shared by particles with different diameters according to their mass fractions. In Figure 7, the contours of gas velocity are depicted. It shows that near the inner wall of the bend, maximum gas velocity occurs so that the radial component of gas velocity leads to dropping particles (especially large particles) toward the outer wall. The trajectory of particles in terms of their diameters are shown in Figure 8. The dispersion pattern of solid particles depends on their size and is shown in Figure 9. As can be seen, the larger particles move toward the outer wall of the bend due to radial component of gas velocity. In Figure 10, the size distribution of particles on different cross sectional planes is drawn versus the relative radial distance according to (19). This figure shows that at the outer wall (), the mean diameter of particles is larger than its value at the inner wall () which is in consistence with that mentioned previously about Figure 8. We can see from Figure 10 that at cross sectional plane of 15°, there are some small particles near the outer wall () of the bend. This plane is located at the region, where the radial component of gas velocity starts to increase but still does not reach its final growth. In Figure 11, the mean particle velocity distribution on each plane is depicted. The variation of particles velocity on each plane is close to a straight line. By increasing the angle of the plane, the slope of velocity variation increases. This is because of increasing the gas velocity at the vertical section of the bend.

Contours of erosion rate and numerical drawing of erosion rate in terms of impact angel are, respectively, shown in Figures 12 and 13, describing that the abnormal rate of erosion is related to the angles between 50 and 60°. The maximum erosion rate is obtained as 3.26 nm/s or 0.1 m/year and occurs at the angle of 52°, which is because of impacting solid particles with size of up to 2 micron. This huge rate of erosion would be dangerous for pipeline fittings and should be controlled by, for instance, the removal of particles with size of up to 1 mm.

#### 7. Conclusion

In this study, we have developed a two-phase Lagrangian CFD model to simulate three-dimensional particulates motion in gas pipeline. The effect of particles diameter on its fluidization pattern was considered by Rosin-Rammler distribution function. Analysis of particulates motion in the bend indicates that due to the increasing trend of radial component of gas velocity through the bend, the larger particles are moved toward the outer wall of the bend and they increase the erosion rate at this region. It is found that the erosion rate in this case is very high due to high particle velocity and high particle diameter. The maximum rate of erosion is discovered around angles between 40 and 65°. This study proves that we can use CFD modeling as a powerful tool for assessing particulates motion and their erosion effects inside different industrial instruments.

#### Nomenclature

: | Phase density (kg·m^{−3}) |

: | Phase volume fraction |

: | Velocity vector for each phase (m·s^{−1}) |

: | Stress tensor for each phase (N·m^{−2}) |

: | Inter phase drag coefficient (Kg·m^{−3}·s^{−1}) |

: | Particle diameter (m) |

: | Turbulent kinetic energy of gas phase (m^{2}·s^{−2}) |

: | Turbulent dissipation rate of gas phase (m^{2}·s^{−2}) |

, : | Turbulent Prandtl number |

: | Generation of turbulence kinetic energy due to the mean velocity gradients (Kg.m^{−1}·s^{−3}) |

: | Kinematic viscosity (m^{2}·s^{−1}) |

: | Single particles drag coefficient |

: | Erosion rate (m·s^{−1}) |

: | Brinell hardness |

: | Particle shape coefficient |

: | Impact angle (Deg) |

: | Particle relative velocity (m·s^{−1}). |

#### Conflict of Interests

The authors of the paper do not have a direct or indirect financial relation with the commercial identity mentioned in their paper that might lead to a conflict of interests for any of the authors.