© T. Firmansyah et al., published by IFP Energies nouvelles, 2018

Nomenclature

E : Modulus of elasticity

I _y : Area moment of inertia

k : Constant

L : Beam length

ṁ: Mass transfer (kg/s)

S : Source term

α : Volume fraction

T : Temperature (K)

t : Physical time (s)

t _r : Fluid residence time inside radiant coil (s)

W : Energy (J)

w : Uniform load per unit length

ρ : Density (kg/m³)

τ : Normalized time

μ : Viscosity (kg/m.s)

μ ₀ : Reference viscosity (kg/m.s)

k _eff : Effective thermal conductivity (W/m.K)

Subscripts

n : Mode of vibration

l: Liquid phase

v: Vapor phase

1 Introduction

A process heater is a direct-fired heat exchanger that utilizes the hot flue gases of combustion to increase the temperature of a feed flowing through coils of tubes aligned throughout the heater [1–3].

Process heaters are widely used throughout the chemical processing and hydrocarbon industries for heating crude oil in the petroleum refining and petrochemical sectors. Most of the unit operations in these plants require fired heaters and furnaces [4].

The roles served by process heaters in chemical and refinery plants vary from simple heating or providing sensible heat and raising the temperature of the feed charge to heating and partial evaporation of the feed charge, where equilibrium is established between the unvaporised liquid and the vapor. The feed charge leaves the furnace in the form of a partially evaporated liquid in equilibrium.

Direct fired heaters consist mainly of two distinct heating sections: a radiant section, also called as a combustion chamber or firebox, and a convection section followed by the stack [4]. The radiant tubes, either horizontal or vertical, are usually located along the walls in the radiant section of the heater and receive radiant heat direct from the burners.

Rather than entering the radiant section directly, the feed charge enters the coil inlet in the convection section where it is preheated before entering to the radiant section. The purpose of the convection section is to gain as much heat as possible from the combustion gases leaving the radiant section. This is done through removal of heat from the flue gas to preheat the contents of the tubes and thus the temperature of the flue gas exiting the stack is significantly reduced. A third section, known as a shield or shock section, separates the two major heating sections. It contains those tubes close to the radiation section that shield the remaining convection section tubes from direct radiation. Fired heaters operated in crude oil distillation plants are large and complex items. There are two main types of fired heaters; Cabin or Box heater and Vertical Cylindrical heater [5–9]. The Cabin heater is preferred for large heat duty applications (approximately 20 MW and more) [7] and contains horizontally oriented tubes (or tube banks) in both radiant and convection parts of heater.

The Vertical Cylindrical heater, on the other hand, is preferred for small and/or medium heat duty applications (below 20 MW). In the cylindrical-type furnace, the radiation section is in the shape of a cylinder with a vertical axis, and the burners are located on the floor at the base of the cylinder. The heat exchange area covers the vertical walls and therefore exhibits circular symmetry with respect to the heating assembly. In the radiant section, the tubes may be in a circular pattern around the walls of the fire box or they may be in a cross or octagonal design which will expose them to firing from both sides. The shield and convection tubes are normally horizontal. The present study concerns with this type of heater, namely splitter reboiler.

Crude oil heated in a fired heater flows inside tube coil in both the convective and the radiant sections. While the flow of fluid in coil of convective section is one-phase liquid flow, the fluid flow inside radiant tube coil is accompanied by boiling and evaporation process, i.e. occurrence of two-phase flow [10].

Two-phase flow of heated oil inside heater radiant coil is accompanied by boiling and increasing of oil evaporation. Thus, the simultaneous and continuous changes of two-phase flow regime are a result of this heating process. The orientation of tube coil significantly influences the nature and characteristics of changes of these flow regimes.

2 Description of problem

The presence of dual phase flow in any piping systems potentially causes a Flow Induced Vibration (FIV) problem. This has been a source of maintenance requirement for heaters in many engineering fields including petroleum pipelines, power and processing plants, and heat exchangers because of the high costs involved in structural damage caused by the fluid-structure interaction. Detail review of research works in two-phase FIV in piping systems can be found in reference [11] and more recently in reference [12].

In general, there are three different types of FIV [12]:

1. “Momentum fluctuation” is caused by the density difference between two phases and large FIV is generated due to the change in flow direction and the impact force of two-phase mixture incident on piping component structures such as elbows and T-junctions.

2. “Thermal-hydraulic vibration is associated with phase change”, induced from the nature of two-phase flow which involves phase change due to the energy transfer through interfacial boundary and/or energy generation within the phase, e.g. boiling and condensation.

3. “Bubble-induced vibration” is due to the dynamics of various shapes and sizes of bubbles that induce sloshing, fluctuation and disturbances within the flow fields.

The above mechanisms occur in FIV phenomena along with the fluctuating characteristics of two-phase flow, namely momentum, pressure and void fraction fluctuations.

More pronounced structural vibration tends to occur in slug flow regime compared to other flow regimes. When a slug unit consisting of slug bubble and liquid passes through the pipe turning segments such as elbow and bend, force acting on the structure boundary significantly fluctuates due to major difference in momentum/mixture flux density, and thus induces structural vibration. This can be exacerbated further when the frequency of momentum/pressure fluctuations matches with the piping natural frequency leading to resonance mode.

Past studies in internal two-phase FIV revealed that the contribution of momentum flux fluctuation is significant only in low frequency range, which is between 0–50 Hz [12]. So, it is important that natural frequency of the piping structure should be set much higher than the two-phase characteristics frequency range (0–50 Hz) to avoid resonance.

The present study discusses frequent vibration problems encountered in a splitter reboiler in a naphtha hydrotreating unit, which caused fretting-wear damage at some of the top return bends of the radiant tube coil. Close inspection of the heater tubes confirmed lateral/sideway vibration visible by naked eyes at these locations. The fretting-wear damage was found exactly at contact areas between the top return bends and the hanger supports, which likely indicated constant rubbing between them during vibration.

A splitter reboiler is a vertical cylindrical process heater, which has common sections as shown in Figure 1.

3 Model description

A root-cause analysis was carried out through a combination of process simulation and Computational Fluid Dynamics (CFD). The CFD analysis emphasized on internal flow inside the radiant tube coil being a possible source of vibration. It is based on a hypothesis that the presence of vapor pockets following boiling and evaporation process within the primary liquid phase leads to a slug flow regime and thus causes FIV problem in the radiant tube coil.

Figure 2 illustrates a schematic of stages conducted in the root-cause analysis. Prior to developing a CFD model a process simulation was first carried out to calculate a vapor to liquid ratio inside the radiant tube. This was then used in CFD as boundary condition in the multiphase simulations involving both primary liquid and secondary vapor phases. All simulations in CFD were performed as transient calculations, i.e. time-dependent simulations. Results showed flow pattern of liquid and vapor phases inside the radiant tube and consequently the amount of fluctuation of pressure generated in time. Several transient simulations were carried out to see the impact of different vapor to liquid ratio on pressure oscillations inside the radiant tube. This includes simulation scenarios with varying operating flow rates.

Any simulated pressure/velocity fluctuation would lead to an inherent vibration problem if its frequency is close to the natural frequency of the radiant tube coil. An empirical formula for calculating natural frequency of tube vibration was used to determine whether the resulting frequency of pressure fluctuation inside the radiant tube is within the acceptable limit to avoid physical vibration problem.

3.1 Model assumptions and limitations

The current model has the following assumptions and limitations:

• Other sources of vibration such as combustion instability could also play a part in causing radiant tube oscillation [13, 14]. However, detailed firebox combustion modeling was outside the scope of the current work and thus, combustion instability was not considered.

• Heat flux is assumed constant across circumference of the radiant tubes, i.e. heat flux of the tube side facing the heater wall is considered the same as that facing the burner. No fouling nor coke deposition is assumed anywhere inside the radiant tubes.

• Since the radiant tube coils are oriented symmetrically around the heater vertical axis, as illustrated in a top view in Figure 3, the heat flux is thus assumed to be radially symmetrical. With this assumption, only one pass of the two radiant tube coils was modeled and for further simplification, the radiant tube coil was stretched as if the heater wall was in plane, as shown graphically in Figure 4.

  Fig. 3. Sketch of top view of radiant tube coils.

  Fig. 4. Model geometry of radiant tube coil.

• In the absence of detailed profile of local heat flux distribution along the radiant section, an approximate heat flux profile was adopted using a typical profile of vertical tube heater as reported in [15].

• Since transient two-phase flow simulations require high computational resources and thus long computational time, the radiant tube coil was represented in CFD as a 2D model by further assuming that the flow field is symmetrical about the mid-plane cross-section of the tube. This reduced the number of required computational cells significantly and hence hugely reduced computational resources and time.

3.2 Model governing equations

Considering that the flow consists of a mixture of liquid and vapor it must be represented by the governing equations for a two-phase flow system. Volume-Of-Fluid (VOF) multiphase model can suitably describe slug flow regime [16].

VOF formulation is based on condition that the two phases are not interpenetrating. In each control volume, i.e. computational cell, the volume fractions of the two phases sum to unity. The variables and properties in any given cell are either purely representative of one of the phases, or representative of a mixture of the phases, depending upon the volume fraction values. So, for the liquid volume fraction in the cell, α _l, the following three conditions are possible:

• α _l = 0: the cell is empty of the liquid;

• α _l = 1: the cell is full of the liquid;

• 0 < α _l < 1: the cell contains interface between the liquid and the vapor.

The governing continuity equation below solves the volume fraction of the liquid phase.(1)where ṁ_lv is the mass transfer from liquid phase, l, to vapor phase, v, and ṁ_lv vice versa. The source term on the right hand side of equation (1), S _αq , by default is zero or can be specified for setting up mass source for liquid phase.

The volume fraction of vapor phase is computed based on the following constraint:(2)

The density in each cell is calculated based on the presence of the component phases in each control volume. For liquid and vapor system, the density is given by(3)

A single momentum equation is solved throughout the domain, and the resulting velocity field is shared among the phases. The momentum equation, described below, is dependent on the volume fractions of both liquid and vapor phases through the properties density, ρ, and viscosity, μ.(4)

The energy equation is also shared among the phases, as described below:(5)

The energy, W, and temperature, T, are treated as mass-averaged variables in the VOF model. Equation (6) below represents the energy:(6)where W _l and W _v are based on the specific heat of liquid and vapor phases respectively and the shared temperature. The properties density, ρ, and effective thermal conductivity, k _eff, are shared by the phases. The source term, S _h, contains contributions from radiation, as well as any other volumetric heat sources.

3.3 CFD model

3.3.1 Geometry and computational mesh

A two-dimensional CFD geometry of the radiant tube coil covering the entire coil pass from inlet to outlet of the radiant section is shown in Figure 5.

Total number of computational cells were about 150 K, all in the form of structured rectangular mesh. Sufficient near-wall boundary layers were implemented to capture any large wall-bounded variable gradients. Figure 6 shows a zoom-in of structured mesh with boundary layers.

3.3.2 Model boundary conditions and material properties

From process simulation calculations, the vapor fraction was set up accordingly as an inlet boundary condition in the CFD model. At the inlet of CFD model, the vapor fraction was specified as zero, i.e. all incoming flow was in the form of liquid. The inlet velocity of liquid was specified using the volumetric flow rate of liquid and the inlet area. All inlet parameters used in the model are listed in Table 1.

Fluid and thermal properties of both the liquid and vapor phases are listed in Table 2.

The radiant tube coil is made of steel with the following material properties, as tabulated in Table 3.

The heat source from burners was modeled as a heat flux distribution along the longitudinal span of the radiant tube coil. This adopted a typical heat flux profile of vertical tube heater as reported in [15] taking into account the design maximum heat flux given by the manufacturer. Figure 7 shows the heat flux variation along the height of the radiant tube.

The maximum heat flux density is around 51.8 K kcal/m³ h at about 4.5 m height.

The mass transfer effects during VOF simulations were modeled with the evaporation-condensation mechanism with a defined saturation temperature taken from process simulation results.

All CFD simulations in this study were performed using a Finite Volume Method (FVM) based CFD software, ANSYS Fluent 16.2. An implicit formulation of volume fraction was implemented for faster and stable transient calculations. High order spatial and temporal discretization schemes were used to achieve good solution accuracy.

Five sets of CFD transient simulations were executed covering different flow rates starting from normal operating flow rate, as listed in Table 4.

4 Vibration analysis

The resulting pressure fluctuations from CFD simulations in a time wave form were converted into frequency form using Fast Fourier Transform (FFT) method. The dominant frequency was then compared with the natural frequency of the radiant tube coil.

The natural frequency of the radiant tube coil was estimated using simply supported beam theory [17]. The formula of natural frequency for simply supported beam is given by equation (7):(7)where k _n is a constant, n mode of vibration, E modulus of elasticity, I _y area moment of inertia, w uniform load per unit length including the tube weight and L beam length. The constant k _n varies with type of supports and mode of vibration.

Using steel as material of beam and taking into account the tube weight, Table 5 lists all the required variables for estimating equation (7).

For hinge supports at both ends of the beam, the values of k for different mode of vibration [17] and the estimated natural frequency modes are shown in Table 6.

5 Simulation results

5.1 Process simulations

Calculations of vapor fractions at three locations of the radiant coil were done based on thermocouple measurements and estimated pressures. A combination of process simulation software and spreadsheet calculations was employed to estimate the vapor fractions. The charge was characterized using distillation and density data. Then, using process data for the same period of time, based on API 530 code [18], the fluid temperatures inside the coils have been estimated using the tube metal temperatures given by the thermocouples TI₁, TI₂ and TI₃, as shown in Figure 8. Pressure profile has been assumed to be linearly decreasing from the pump charge to the heater outlet. The vapor fractions at discrete locations can then be calculated using the characterized feed, estimated temperatures and pressures. Table 7 summarizes the estimated vapor fractions and pressures at the three thermocouple locations.

The inlet of radiant coil is at the same side as the thermocouple TI₃.

5.2 CFD simulations

5.2.1 Normal operating flow rate

Because of heating and boiling process, a two-phase flow is generated downstream of the radiant tube inlet. Figure 9 shows steady state vapor phase volume fraction after 120 s, which clearly indicates the liquid undergoing evaporation process towards the end of the first radiant tube downstream of the inlet (far right). The blue color represents full liquid phase while red color represents full vapor phase. Following this evaporation process, slug liquid type of flow also appeared in the second radiant tube. As a result of continuous vaporization, the rest of radiant tube coil flow is dominated by a dual phase characteristics.

All calculated pressures across all of the top return bends were recorded at each time step during the course of transient CFD simulation. However, for clarity purpose, Figure 10 displays the computed pressures vs. normalized time, as defined in equation (8), at three representative locations of the top return bends only, i.e. top bends #1, #4 and #8. The numbering of these top return bend locations is shown in Figure 8. Clearly the dual phase flow inside the radiant coil has generated significant pressure fluctuations in all locations of the top return bends. This behavior of flow instability potentially leads to flow induced vibration issues if the resulting dominant frequency of this pressure fluctuation is close to the natural frequency of the radiant tube system.(8)

where τ denotes the normalized time, t is physical time (seconds) and t _r residence time (seconds) of fluid inside the radiant tube coil.

5.2.2 Fast Fourier Transform analysis

Assessment of possible flow induced vibration was conducted by converting the pressure fluctuation in the time wave form to dominant frequency analysis through FFT approach. The Power Spectral Density vs. Frequency plot was then generated as shown in Figure 11 from the resulting pressure fluctuation recorded at the top return bend #4. The dominant frequency at this location was found at 13.2 Hz. Table 8 lists all the dominant frequencies from FFT analysis of dual phase flow pulses at all the top return bends.

As evident, the dominant frequencies mostly stay in the range of 12.5–14.3 Hz.

5.2.3 Resonant vibration prediction

Resonance mode would occur if frequency of the flow induced vibration is within 20% of the calculated natural frequency of the tube coil [17, 19]. Table 9 lists percentage differences between the resulting dominant frequencies and the natural/fundamental frequency of the tube, i.e. f ₁ = 11.3 Hz.

Most of dominant frequencies were on the higher side of the natural frequency, but still within the band of 20% deviation. Therefore, resonance mode most certainly occurred at the top return bend numbers 1–4, 6, and 7. This is in line with at least two locations of severe grooving found during heater tube inspection, i.e. top return bend numbers 4 and 6. So CFD results confirm that resonant vibrations were caused by dual phase flow pattern inducing pressure fluctuation at majority of the top return bends of the radiant tube coil.

5.2.4 Possible routes to restrain resonant vibrations

There are two main possible routes to reduce or remove resonant vibrations. Table 10 lists these two alternate routes with possible modifications of related parameters respectively.

Ideally the radiant tube coil should have been designed to have the natural frequency far from the frequency of flow induced vibration. A detailed review of flow induced vibration study [12] mentions that frequencies of flow induced vibrations usually fall within a range of 0–50 Hz and therefore recommends that the natural frequency of the structure must be designed above 50 Hz.

Increasing natural frequency of the radiant tube coil can be achieved by modifying the tube dimensions in terms of the length and/or the diameter. This option, however, requires major overhaul of the radiant section of the heater, which is not always feasible. A more straightforward option is by installing an additional support, which will significantly increase the natural frequency of the radiant tube coil.

Dual phase flow inside the radiant tube coil has been clearly demonstrated by the current simulations as the source of the flow instability and thus induced resonant vibration. Therefore, modifying this dual phase flow regime by means of reducing vapor fraction can surely reduce flow instability and thus avoid resonant vibration. Impacts of modifying flow rates and installing additional support are further discussed below.

Increasing flow rates certainly increased the overall velocity inside the tube and vice-versa for decreasing flow rates. Figure 12 clearly demonstrates the impacts of varying flow rates on the velocity pattern. Red color represents the upper scale of velocity (10 m/s) whereas the dark blue is zero velocity.

More importantly, increasing flow rate has also led to an increase in slug formation. Figure 13 shows contour of vapor fraction for the three different flow rates. Clearly by decreasing the flow rate the slug formation also decreases and this has direct consequence to the level of pressure fluctuation, i.e. frequency of the induced vibration.

Table 11 shows the change in dominant frequencies of the pressure fluctuations following change of flow rates.

Clearly the increasing flow rates tend to increase the frequency of the pressure fluctuations as a direct consequence to the increasing appearance of slugs. Table 11 also showed that lowering flow rates have more significant impacts on changing the flow induced frequencies away from the natural frequency of the radiant coil. While increasing flow rate by 20% above the normal operating flow rate still produced resonant frequencies within 20% band of the natural frequency at locations 1, 7 and 8, lowering flow rate even by 10% below the normal operating flow rate has removed resonant vibration completely at all locations. This is clearly demonstrated in Figure 14 where the lowering flow rates have moved the flow induced frequencies away from the 20% band of natural frequency deviation represented by the light red shaded band.

So, while it is recommended to check the current capacity utilization compared to design, the CFD simulation predicts that a flow decrease by about 10% below the normal operating flow would keep the tubes away from vibration phenomenon.

Installing additional support on the radiant tube coil shall definitely increase its natural frequency. Courant and Hilbert [20] showed that the optimum locations of the supports should be at the nodal points of higher vibration mode without the supports. For simply supported beam, the optimum location of additional support is always halfway of the beam span [21]. Placing an intermediate support with sufficient stiffness at the optimum location can double the natural frequency of simply supported beam [21].

Figure 15 depicts the possible maximum increase in natural frequency of the radiant tube coil by installing additional support at the mid-height of the radiant tube and thus avoiding the resonance mode. In that scenario, all frequencies of the flow induced fluctuations will be significantly below the 20% band of the new natural frequency deviation represented by the yellow shaded band.

6 Conclusion

A combined study of process, CFD and vibration analysis was conducted to find the root-cause of vibration problem in a splitter reboiler of a naphtha hydrotreating unit. The numerical model was developed with VOF multiphase method in ANSYS Fluent software. CFD simulations have shown formation of slug liquid flow along vapor pockets and thus inducing pressure fluctuations inside the tube coils in the radiant section. FFT was employed to determine the dominant frequencies of the predicted pressure fluctuations. Some of the resulting dominant frequencies were within 20% band of the estimated natural frequency of the tube, which could lead to resonance mode. Two locations from CFD prediction using the operating flow rate matched with the locations of severe grooving as observed during a site inspection.

Decreasing the flow rate has been demonstrated to be more practically implementable scenario rather than increasing the flow rate, in order to steer the operating conditions out of the vibrational mode. Vibration of the tubes can also be avoided by adding support rings at intermediate locations (e.g. midway) thereby modifying the natural frequency of the tube.

Acknowledgments

The authors gratefully acknowledge the support of ADNOC Refining Research Centre Division (RCD) and ADNOC Refining management. This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sector.

References

All Tables

All Figures