Article Cite This: Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX
pubs.acs.org/IECR
Use of Distributed Threshold Fouling Model at Local Scale in Numerical Simulations Anthony Chambon,† Zoe ́ Anxionnaz-Minvielle,*,‡ Jean-François Fourmigue,́ ‡ Nathalie Guintrand,† Aureĺ ien Davailles̀ ,§ and Fred́ eŕ ic Ducros‡ †
TOTAL Research & Technology Gonfreville, BP 27, F-76700 Harfleur, France CEA, LITEN, DTBH, Laboratory of Heat Exchangers and Reactors, Université Grenoble-Alpes, F-38000 Grenoble, France § TechnipFMC, 6-8 Allée de l’Arche, F-92973 Paris La Défense, France Downloaded via BUFFALO STATE on August 2, 2019 at 15:03:20 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.
‡
ABSTRACT: Numerical simulation is studied for shell and tube heat exchangers incorporating the usage of the local Ebert and Panchal threshold fouling model on the overall fouling prediction. The overall fouling rate obtained by averaging local fouling rates along the whole heat exchanger was found to be consistent with the fouling rate calculated with overall quantities, at least for the shell-and-tube heat exchanger tested, under the operating conditions studied (299.5 °C film temperature, 4.4 Pa wall shear stress) and for the period of linear fouling. The deposition rate, which is minimally affected by the spatial variations of the thermohydraulic properties, controls the fouling rate. Heterogeneities are caused by variations of the local wall shear stress which varies in a large proportion (2.5 to 11.7 Pa) but whose influence on the fouling rate is low. The wide range of local shear stress recorded is the consequence of the flow establishment in the tubes.
1. INTRODUCTION Fouling of heat exchangers is the buildup of deposit layers on the heat transfer surfaces. Crude oil fouling in refineries causes several operating, financial, and safety issues. Before entering the atmospheric column, crude oil is heated in the preheat train by recovering heat from product streams. However, the thermal performance is reduced by fouling due to an increase of the thermal resistance leading to additional fuel consumption and CO2 emissions. Moreover, the gradual decrease in the tube cross-sectional area due to the fouling layer growth requires more pumping power to maintain the throughput. All these consequences are the cause of numerous additional expenditures. According to Coletti and Hewitt,1 the cost of the fouling in the U.S. refineries is estimated at USD 3.6 billion of which USD 2.26 billion is only linked to preheat train inefficiency. It also represents a release of 88 million tons of carbon dioxide (2.5% of the worldwide emissions in 2009). The U.S. Department of Energy2 considers that the energy required for this process could be reduced by 15% with a better management of fouling. Energy spent for crude oil fractionation (the most energy-consuming step) is estimated at 192 TWh a year for all the American refineries. A 15% savings thus represents 6,340 bbl a day or 1% of throughput of the biggest U.S. refinery. Although most of the phenomena responsible for crude oil fouling remain unknown, several efforts were accomplished to © XXXX American Chemical Society
understand it. Several pieces of literature suggest links between fouling of crude oil at the hot end of the preheat train to the precipitation of asphaltenes,3,4 a high molecular weight and sulfur-rich component of crude oil. These molecules will then diffuse toward the heat transfer surfaces. Then, fouling precursors undergo chemical transformations and stick to the wall.1 At the same time, the crude oil flow acts on the fouling layer by removing it or by inhibiting deposition phenomena.5 Thus, the net fouling rate is defined as an equilibrium between a deposition and the elimination rate of the deposit.6 The crude type, the high temperatures levels, and the flow velocity seem to be the most influential parameters on the fouling rate in the preheat train.7 Phenomenological fouling models for crude oil are hindered by the complex composition and the great variety of the crudes combined with a relative lack of knowledge of the exact fouling mechanisms involved. To overcome these difficulties, threshold fouling models were revealed as an easy-to-use tool to quantify fouling from readily measurable operating variables (flow velocity and temperature). However, adjustment of the model constants is very dependent on the nature of the crude and the Received: Revised: Accepted: Published: A
April 29, 2019 June 18, 2019 June 27, 2019 June 27, 2019 DOI: 10.1021/acs.iecr.9b02328 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX
Article
Industrial & Engineering Chemistry Research experimental device8 and makes sense only in a global way (say at the component scale, i.e., single heat exchanger scale), since it relies on global performance measures; but combined with a thermohydraulic numerical simulation, threshold fouling models provide efficient predicting tools.9−12 Called dynamic threshold fouling models, they allow for simulatation of the overall behavior of a shell-and-tube refinery heat exchanger from a distributed threshold fouling model. Temperature and pressure profiles, thickness of the deposits, and fouling resistances can be monitored all along the heat exchanger and during a fouling-cleaning cycle. The threshold fouling models allow the prediction of fouling rates from measurable operating parameters such as flow velocity and temperature. When implemented in a thermohydraulic numerical simulation, local fouling rates are assessed from discretized operating parameters. Distributed models are a promising way for the HEX design and maintenance improvement.9−12 However, the use of such modelings in a predictive manner encounters some specific theoretical difficulties, coming along with the needs to define local physical models, i.e., models able to describe the physics at the cell scale in a CFD framework. To that respect, a local version of semiempirical fouling models should be assessed, which is one of the objectives of this work through the study of the effect of local use of the Panchal et al.13 threshold fouling model on the fouling prediction when used in a shell and tube heat exchanger numerical simulation. For that purpose, a single CFD numerical simulation of the heat exchanger in clean conditions is performed and compared to a set of experimental data. Results in terms of a posteriori predicted fouling rates are compared according to the averaging method (space discretization or numerical and experimental averages). The discretization also allows for quantification of the relative importance of the fouling mechanisms (deposition and removal terms) and operating variables. To validate the CFD thermohydraulic simulation and compare the predicted fouling rates, experimental data in clean conditions are used. As a consequence, the experimental setup is briefly described in the first section. Then, the CFD simulation of the experimental heat exchanger and its validation are detailed. Finally, the prediction of the averaged local fouling rates obtained with the CFD simulation is compared with the prediction of the overall fouling rate obtained from averaged input parameters (both numerical and experimental operating parameters).
Figure 1. Picture of the test section, a small-scale shell-and-tube heat exchanger.
to-fluid temperature gradients are mimicked. Crude oil used for fouling runs comes from the Black Sea area. Adjustable parameters of the rig are inlet residue temperature and crude oil velocity. Crude oil and residue circuits are pressurized with nitrogen at the same pressure as the industrial process. Chemical analyses (density, asphaltene content, Conradson residue, sulfur content, and S-value) are regularly performed to monitor any change in fluid composition and avoid fouling materials depletion. 2.2. Data Reduction. 2.2.1. Fouling Rate Correlation. Predicted fouling rate is assessed with the Panchal et al.13 correlation in eq 1: dR f dt
exp
ji E zyz = αRe βPr −0.33 expjjj− − γτ j RT zzz f{ k
(1) 6
This equation is based on the Kern and Seaton formulation, i.e., it assumes the fouling rate as a difference between a deposition rate (first term on the RHS) and a removal rate (second term). Re, Pr, and Tf [K] refer respectively to the Reynolds number, the Prandtl number, and the film temperature of crude oil. The activation energy of fouling reactions is denoted E [J/mol], R [J/mol/K] stands for the gas constant, and τ [Pa] is the wall shear stress. To compare the fouling rates predicted from experimental and numerical operating data, the model parameters used for this work are α = 23.4 m2 K/J β = −0.8 γ = 6.88 × 10−12 m2 K/J/Pa Ea = 75.1 kJ/mol They have been optimized from several fouling tests performed with a specific crude oil type and a heat exchanger configuration where fouling occurs on both the tube-side and the shell-side.15 2.2.2. Operating Conditions and Fluid Properties. In the following sections, only a single set of clean operating conditions is considered and summarized in Table 1.
2. EXPERIMENTAL DATA 2.1. Experimental Setup. The setup used to get data has been described in previous works.1,14,15 The reader may refer to these references for additional information. The rig aims to reproduce the operating conditions (temperature and pressure) encountered at the end of refineries preheat train where fouling is promoted by high wall temperatures. As an industrial process, crude oil is heated by the residue (atmospheric tower bottom) through the test section illustrated in Figure 1. The test section is a shell-and-tube heat exchanger which reproduces, at a small scale, those operated in refineries. It has been designed to be a scaled-down reproduction (Ø0.45 × 2.70 m) of the warmest heat exchanger of the refinery preheat train (same tube size, material, heat flux, and heat transfer coefficients). Thus, the fluids circulating in the test section are identical to refineries and flow with the same scheme (i.e., residue in the shell side and crude oil in the tubes). Hence, wall-
Table 1. Operating Conditions of a Fouling Run inlet temperature [°C] outlet temperature [°C] mass flow rate [kg/s] heat duty [kW] pressure drop [mbar]
crude
residue
293.9 ± 0.2 300.2 ± 0.2 2.343 ± 0.012 43.6 ± 3.4 565 ± 4.5
314.3 ± 0.2 309.7 ± 0.2 3.200 ± 0.012 43.6 ± 4.5 928 ± 4.5
Reynolds and Prandtl numbers are estimated at bulk temperature Tb [K], which is the mean temperature between inlet Ti [K] and outlet To [K] temperature, eq 2: Tb = 0.5(Ti + To) B
(2) DOI: 10.1021/acs.iecr.9b02328 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX
Article
Industrial & Engineering Chemistry Research
For the studied fouling run, experimental tube-side data reduction leads to Re = 57,239 Pr = 12.2 Tf = 299.5 °C 2.2.2.2. Shear Stress. Wall shear stress τ [Pa] in a circular pipe is calculated with eq 6
Physical properties of the crude oil and the residue are assumed to be constant along the heat exchanger (Table 2). They are extrapolated with PRO-II software from properties assessed at ambient temperature. Table 2. Physical Properties of Fluids and Carbon Steel Used in the CFD Study Tb [°C] density [kg/m3] specific heat [J/kg/K] conductivity [W/m/K] viscosity [Pa s]
crude
residue
carbon steel
297.1 675 2980 9.18 × 10−2 3.75 × 10−4
312.0 730 2925 8.94 × 10−2 7.23 × 10−4
7800 470 45
τ = 0.5fρu 2
where f is the Fanning friction factor, ρ [kg/m ] is the density of oil, and u [m/s] is the crude velocity inside the tubes. Tube roughness is assumed negligible, so the Blasius correlation (2 × 103 < Re < 105) is used to find the Fanning friction factor according to eq 7: f = 0.0791Re−0.25
2.2.2.1. Film Temperature. Film temperature Tf [K] is the average of the wall temperature Tw [K] and the bulk temperature Tb [K] according to the definition proposed by Ebert and Panchal16 in eq 3. This definition of film temperature has been chosen to be consistent with Ebert and Panchal16 and Panchal et al.13 works. Tf = 0.45Tb + 0.55Tw
Q hi A i
(3)
3. DETAILS OF THE CFD MODEL 3.1. Solving Method. The Ebert and Panchal correlation (see eq 1) is strongly nonlinear with respect to its main arguments (Re, Pr, and Tf) that are evaluated through their average values within the whole heat exchanger for a given operating point. This property makes the use of this model using local values (such as local Reynolds, temperatures, ...) questionable. Indeed, the multiplication of their average values is different from the average of their multiplication and may affect the final result of the fouling rate: our goal is to evaluate the impact of what can be considered as a theoretical shortcoming. So, we will perform a thermohydraulic simulation to study the influence of operating data discretization and averaging method on predicted fouling rates. Flow in the test section is turbulent on both sides of the heat exchanger: residue Reynolds number (calculated in the crossflow section of the shell) is 9,800 and
(4)
2
where Ai [m ] is the internal heat transfer area of the tube bundle, and hi [W/m2/K] is the local convective heat transfer coefficient inside the tubes derived from the Dittus-Boelter correlation, eq 5: hi d i = 0.0243Re 0.8Pr 0.4 λc
(7)
Shear stress for the fouling run is τ = 4.4 Pa. These operating data, Reynolds number, Prandtl number, film temperature, and shear stress, are used to predict the experimental fouling rate, denoted (dRf/dt)exp below, from clean conditions.
Wall temperature is calculated from a heat balance carried out on the tube bundle of the test section, eq 4 Tw = Tb +
(6) 3
(5)
Here, λc [W/m/K] refers to the crude oil thermal conductivity. Tubes of the tube bundle have an internal diameter of di = 19.86 mm.
Figure 2. Overall heat exchanger geometry and close-up view showing the sealing strips. C
DOI: 10.1021/acs.iecr.9b02328 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX
Article
Industrial & Engineering Chemistry Research crude oil Reynolds number is 57,200. A turbulence model is required. The k-ε RNG (ReNormalization Group) turbulence model is chosen to simulate crude oil and residue flow since Yang et al.17 and Pal et al.18 demonstrated its suitability for hydrocarbons flowing in a shell-and-tube heat exchanger. Since a wall function is used to simulate the boundary layer, the mesh refinement close to the walls is not necessary. It allows the reduction of the total cell number of the mesh. In addition to continuity and momentum conservation equations, the equation of energy conservation is also solved (conduction and forced convection). Simulation is carried out at steady state with FLUENT R17.2. 3.2. Heat Exchanger Geometry. According to the TEMA (Tubular Exchanger Manufacturers Association, Inc.) heat exchanger designation, the test section is a BES-type shell-andtube heat exchanger. Some simplifications were made to build the FLUENT geometry (Figure 2): • Dead zones at the end of the heat exchanger are removed (shell head over the bundle floating head); • 8 mm-diameter tie rods in the bundle are also removed; • Baffles and tube-sheets are assumed to have no thickness; • Tube-to-baffle and baffle-to-shell clearances are neglected. Nevertheless, some details which could affect the fluid flow are integrated in the 3D-geometry, such as • The impingement plate in front of the shell inlet; • The four sealing strips around the tube bundle. Inlet and outlet piping are long enough to prevent numerical instabilities due to entrance and exit effects. Shell and head volumes are separated by a thin empty volume. This geometric trick allows for thermally insulating the shell side from the tube side while simplifying the task of the mesher since it is no longer necessary to connect meshes at head and shell interface. 3.3. Meshing. The mesh has been sized and optimized according to velocity field plot obtained with several grids tested to reach the best compromise between accuracy and computational effort required. It is a hybrid mesh containing both hexahedral and tetrahedral cells. The tube bundle and most of the shell volume are meshed with hexahedral cells, whereas front and rear head, impingement plate, and areas close to nozzle connections are meshed with tetrahedral cells. This choice was imposed by the complex geometry for which the sweeping method cannot be applied. 3.3.1. Cell Size. To capture the fluid motion between the tubes of the tube bundle, a minimum of three or four cells in the gap is required. Since the gap is 6.35 mm wide, the typical dimension of a cell must be 1.6 mm. All the other cell sizes of the mesh come from this typical dimension to avoid jumps in cell size and get a consistent mesh. 3.3.2. Axial Length of Mesh Cells. To connect hexahedral and tetrahedral mesh and avoid deformed cells, tetrahedron dimensions must be similar in each space direction. The typical cell size (1.6 mm) imposes the axial cell size at the tips of the tubes. After several tests of meshes, a cell length of 5 mm in crossflow volumes has been chosen. To ensure a continuous junction of 2 mm and 5 mm mesh areas, buffer zones are implemented where the cell size increases gradually (Figure 3). 3.3.3. Tetrahedral Mesh. Impingement plate and residue outlet nozzle areas are meshed with tetrahedral cells having a 1.6 mm typical size to ensure the junction with tetrahedral adjacent mesh. Front and rear head are meshed with 15 mm tetrahedrons (Figure 4).
Figure 3. Close-up of the grid in a tube.
Figure 4. Close-up of the header mesh with tetrahedrons.
3.3.4. Mesh Size and Wall Function. Operating conditions of Table 1 and dimensions of the first cell adjacent to the tube wall lead to dimensionless wall distance y+ = 44 and y+ = 131, respectively, for the inner and outer side of tubes. These values match with FLUENT User’s Guide recommendation (30 < y+ < 300) for the ‘enhanced wall f unction’. 3.3.5. Grid Independence Considerations. The current mesh reaches a total cell number of 14 million. Within the tube: as previously mentioned, we make use of a standard mesh sizing in the normal direction to the wall, as well as in the flow direction since flow properties vary slowly (the flow to be described is attached boundary layer over a flat tube type). The situation in the intertube space is slightly different since the mesh chosen is cruder: this is of second order effect in the present study since we only focus on fouling occurring inside the tube, the outer part of the tube being only required to keep the thermal balance correct. These statements are also supported by the comparison between experiments and numerical predictions, as described in section 4.1. 3.4. Case Setup. 3.4.1. Boundary Conditions. Test section inlet nozzles (crude oil and residue) are set to velocity inlet, and inlet temperature is constant (see values in Table 1). Outputs nozzles are set to pressure outlet with zero pressure as reference. 3.4.2. Initial Conditions. Calculation is initialized with the “standard” option. The level of turbulence is initialized (at inlets and in fluid volumes) so that the inequality 8 is verified. D
DOI: 10.1021/acs.iecr.9b02328 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX
Article
Industrial & Engineering Chemistry Research Cμ
k2
