Modeling and Simulation of a Multibed Industrial Hydrotreater with

Jun 7, 2014 - Analysis of the effect of temperature indicated that reactor performance ... reactors.4 The flash program was coupled to a reactor model...
0 downloads 0 Views 506KB Size
Article pubs.acs.org/IECR

Modeling and Simulation of a Multibed Industrial Hydrotreater with Vapor−Liquid Equilibrium Anton Alvarez-Majmutov and Jinwen Chen* CanmetENERGY, Natural Resources Canada, One Oil Patch Drive, Devon, AB T9G 1A8, Canada ABSTRACT: Dealing with vapor−liquid equilibrium (VLE) in hydrotreating reactors is a major obstacle to realistic process simulation. The present study focused on the modeling and simulation of a commercial light-cycle-oil (LCO) hydrotreater with rigorous vapor−liquid equilibrium (VLE) calculations. The hydrotreater was represented as a one-dimensional plug-flow adiabatic reactor with interbed gas quenching. The model accounted for the hydrodesulfurization (HDS) and hydrodearomatization (HDA) reactions. VLE calculations were performed in situ with a calibrated flash calculation program developed in-house specifically for hydrogen−petroleum systems. Significant differences in reactor temperature profiles, hydrotreating conversions, fluid rates, and phase compositions were observed between the simulations with and without VLE. Analysis of the effect of temperature indicated that reactor performance depends heavily on the axial temperature distribution, which is established based on the catalyst bed layout and the selection of bed inlet temperatures. It was observed that, although conversion levels were improved at higher inlet temperatures, there was an increase in quenching gas requirement, vaporization rate, and the formation of high-temperature zones toward the inlet of the reactor. In a similar manner, pressure was found to increase hydrotreating performance markedly but with an extensive heat release. Verification of plug-flow and full-catalyst-wetting conditions during the simulations confirmed that all criteria were satisfied in this study.

1. INTRODUCTION Catalytic hydrotreating is an essential step in petroleum refining operations for producing quality blending components for transportation fuels. It allows for the removal of hydrocarbon contaminants, such as sulfur, nitrogen, and oxygen, and the saturation of aromatic rings and olefins, without much change in the boiling range. Gas-oil hydrotreating, in particular, has received significant attention because of the stringent environmental legislations for diesel fuels. The process is generally carried out in fixed-bed reactors where both gas and liquid phases flow cocurrently through a catalyst bed in the trickleflow regime.1 Accurate design and optimization of commercial units currently demands robust modeling tools that must rely on a more realistic description of the process. However, this can be a challenging task, as there are complex interactions between hydrotreating chemistry and the physical phenomena that take place in this type of three-phase system (gas−liquid−solid). Ideally, all of these events must be coupled into a reactor performance model to achieve the best predicting capability. One of the difficulties related to the modeling of gas-oil hydrotreating is how to deal properly with vapor−liquid equilibrium (VLE). It is a well-established fact that petroleum feedstocks tend to vaporize into the gas phase under hydroprocessing conditions.2−5 This event creates a complex disturbance in the reactor by changing fluid velocities, phase compositions and properties, and consequently flow dynamics. At the catalyst particle level, vaporization also affects reaction kinetics because some of the most volatile reactants in the vapor phase do not have access to the active sites of the catalyst. Conversely, the heaviest and most refractory compounds (e.g., polyaromatics and dibenzothiophenes in gas-oil feeds) are concentrated in the liquid phase, which extends their contact time with the catalyst.6 Additionally, VLE has an important effect on the heat balance of the reactor, as vaporization of Published XXXX by the American Chemical Society

hydrocarbons consumes part of the heat released by chemical reactions.7 This makes it more difficult to predict the temperature rise along the catalyst beds in industrial reactors. Thus, VLE is a highly influential factor for hydrotreating processes that cannot be ignored during process design, modeling, and simulation. Although noticeable progress has been made in the modeling of hydrotreating reactors,8 feedstock vaporization inside the reactors is overlooked for simplicity. As a result, it is implicitly assumed that the actual gas- and liquid-phase flow rates, bulk properties, and compositions are the same as those of the respective streams before being mixed at the entrance of the reactor, which usually is far from reality. Moreover, the change in liquid flow rate due to VLE can define whether a reactor operates under plug-flow and full-catalyst-wetting conditions,4 which is a fundamental aspect to be considered during the formulation of model equations. Most of the time, all of the effects associated with vaporization are lumped into the fitting parameters of the model,9 making it less predictive. Only a few modeling studies10−15 have included flash calculations in the solution of the reactor mass and heat balance equations, but in all cases, the interaction parameters in the equation of state (EoS) were somehow estimated from correlations rather than being properly evaluated from actual VLE data. Experimental results from our prior studies on the VLE behavior of several gas-oil feeds under typical hydrotreating conditions showed that the extent of oil vaporization increases at higher temperatures, lower pressures, and/or higher H2/oil Received: March 10, 2014 Revised: May 12, 2014 Accepted: June 7, 2014

A

dx.doi.org/10.1021/ie501032j | Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Industrial & Engineering Chemistry Research

Article

ratios.2−5 The quantitative information from these studies was used to develop a robust in-house flash calculation program that accurately predicts VLE in gas-oil hydroprocessing reactors.4 The flash program was coupled to a reactor model to analyze the behavior of a commercial single-bed gas-oil hydrotreater.16 Reactor simulations with and without VLE revealed significant differences in hydrotreating performance and axial temperatures, confirming the importance of accounting for VLE in modeling hydrotreating processes. As a follow-up to our research on VLE in trickle-bed hydrotreating reactors, this work aimed at performing modeling and simulation of a more representative industrial-scale hydrotreater configuration, consisting of multiple catalyst beds and interbed gas quenching. The hydrotreating of light cycle oil (LCO) was selected as a case study to illustrate the VLE behavior in a highly exothermic process.



1 dFiG − ki LaL(Ci* − Ci L) = 0 AR dz

(2)

where rj represents the apparent rate of reaction, with j = HDS, HDAPA, HDADA, and HDAMA. The negative sign on ρbrj is for the reactants, whereas the positive sign is for the products. The equilibrium concentration for component i at the vapor−liquid interface in eqs 1 and 2 is defined as

Ci* =

CiG C L K i* CG

(3)

where Ki* is the VLE partition coefficient calculated from the equation of state. A heat balance is also required to describe the adiabatic mode of operation of industrial-scale reactors. Considering the assumption that there are no interface heat transfer limitations, reactor temperature is modeled with the following pseudohomogeneous energy balance equation

2. MODELING APPROACH The hydrotreater was modeled as a one-dimensional plug-flow adiabatic reactor coupled with rigorous VLE calculations. In the reaction scheme, only the hydrodesulfurization (HDS) of organic sulfur compounds and the hydrogenation of aromatics (HDA) were considered, as these reactions are the ones that contribute the most to the total heat release. The reaction− VLE system was described in terms of two parallel groups of chemical components. The first consisted of the species participating in the hydrotreating reactions: five chemical lumps within the oil feed (organic sulfur compounds, monoaromatics, diaromatics, polyaromatics, and naphthenes) and two pure gas compounds (H2 and H2S generated by HDS). The second group was composed of 29 hydrocarbon pseudocomponents and represented the basis for the flash calculations. The model equations were formulated under the following assumptions: • The reactor is operated in plug-flow mode with complete catalyst wetting (to be verified during the simulations). • The system is in the steady-state regime. • The pressure is constant. • VLE is achieved at any position of the catalyst bed. • There is no heat transfer resistance among the gas, liquid, and solid phases. • Liquid−solid external mass transfer can be neglected. • The effect of intraparticle diffusion is implicit in the apparent rate constants. 2.1. Reactor Model Equations. The formulation of the model equations was based on the transport of reactants and products between the vapor and liquid phases in accordance with the two-film theory.17 The continuity equation for component i in the vapor phase can be expressed as −

1 dFi L + ki LaL(Ci* − Ci L) ± ρb rj = 0 AR dz

−(ρG CpGuG + ρL CpLuL) +

dφ dT dT − (HG − HL)G0L G dz dT dz

∑ (−ΔHrj)ρb νjrj = 0 (4)

j

The second term in eq 4 corresponds to the heat absorbed during hydrocarbon vaporization. The third term is the heat generated by the HDS and HDA reactions. Interbed quenching is represented as the mixing of the quenching gas with the hot process fluids coming from a catalyst bed. The quench gas rate ( f Q) required to cool the reacting mixture to a certain bed-inlet temperature (T0) or vice versa can be estimated from the heat balance in the quench zone18

∫T

T0

foutG CpG dT +

out

∫T

T0

foutL CpL dT +

out

∫T

T0

fQ CpQ dT

Q

(5)

=0

The gas flow rate entering the next catalyst bed is adjusted by adding the quenching gas to the catalyst bed gas effluent. Gas composition is also updated by solving the individual mass balances of the vapor-phase components in the quench zone. The boundary conditions for the above ordinary differential equations (ODEs) are z = 0:

FiG = Fi0G , Fi L = Fi0L , T = T0

(6)

This set of coupled ODEs is integrated simultaneously in the axial reactor coordinate using the Runge−Kutta method. The initial values for the molar rates and concentrations of individual components are obtained through a flash calculation at the reactor inlet. The number of catalyst beds and their respective dimensions are specified as discussed below. The integration is performed until reaching the outlet of the catalyst bed, after which quenching gas is injected. The initial conditions for the next catalyst bed are reset after solving the quenching-zone equations. The flash program is used to perform the VLE calculations in situ at every integration step and to estimate the thermophysical properties of both phases based on various correlations.4 Fluid flow rates and phase compositions are thus corrected as functions of reactor temperature, pressure, and H2 /oil ratio. Other model parameters, such as gas−liquid mass transfer coefficients and

(1)

where i represents H2, H2S, sulfur compounds (S), monoaromatics (MA), diaromatics (DA), polyaromatics (PA), and naphthenes (NA). Equation 1 essentially states that the change in molar flow rate of any component in the vapor phase (FiG) depends only on the mass transfer to or from the liquid bulk. The change in flow rate of component i in the liquid phase is a function of mass transfer across the vapor−liquid phase boundary and consumption or generation by chemical reactions B

dx.doi.org/10.1021/ie501032j | Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Industrial & Engineering Chemistry Research

Article

2.3. VLE Modeling. The flash program plugged into the reactor model was specifically designed for oil−hydrogen environments. Thermodynamic modeling of phase equilibrium is based on the Peng−Robinson EoS. For the calculations, the LCO feed is partitioned by boiling point into 29 pseudocomponents using SimDis data. Molecular weights and densities of each pseudocomponent are estimated with boiling-point-based correlations, which were derived from a detailed characterization of eight subfractions from several gas-oil feeds. Critical properties and other parameters were obtained from various correlations that use these three properties as bases. In hydrogen−hydrocarbon systems, the key element of VLE modeling is determining the binary interaction coefficients between hydrogen and the hydrocarbon pseudocomponents. To our best knowledge, no correlations are available for estimating such parameters. For this reason, flash experiments were conducted to obtain VLE data with a variety of gas-oil feeds, covering a wide range of temperatures and pressures.2−5 The interaction coefficients between hydrogen and the pseudocomponents were then obtained by fitting the EoS to the experimental data. In addition, these coefficients were further correlated with the boiling points of the pseudocomponents and the aromatic content of the oil feed as follows

liquid diffusivities, are calculated using correlations reported in the literature.17,19 2.2. Reaction Kinetics. The reaction term in the mass balance equations is described by an apparent kinetic model. The reaction scheme involves the following hydrotreating reactions S + 5H 2 → MA + H 2S

(7)

HDAPA :

PA + 2H 2 ↔ DA

(8)

HDADA :

DA + 2H 2 ↔ MA

(9)

HDAMA :

MA + 3H 2 ↔ NA

(10)

HDS:

For the HDS reaction, it is assumed that dibenzothiophenes (DBTs) are the representative sulfur compounds in the LCO feed.14 Sulfur removal proceeds preferentially through the hydrogenation pathway with H2S and a monoaromatic hydrocarbon as the main products.20 This reaction is considered to be irreversible under typical operating conditions.21 Aromatic compounds are grouped as monoaromatics, diaromatics, and polyaromatics, the last of which are assumed to behave as triaromatics. The various saturation reactions (HDAPA, HDADA, and HDAMA) proceed successively ring by ring and are reversible under process conditions. HDS is modeled with the following Langmuir−Hinshelwood-type rate equation17 −rHDS =

kHDSCSL m1C H2L m2 (1 + K H2SC H2SL)2

di H = A + BBPi

(15)

A = 0.3852 − 0.00241Caromatics%

(16)

B = 0.002245 + 1.96 × 10−5Caromatics%

(17)

Equations 15−17 are valid in the temperature range of 250− 400 °C and pressure range of 5−10 MPa. Figure 1 illustrates

(11)

Hydrogenation kinetics are described with the equations proposed by Chowdhury et al.22 ⎛ ⎞ C −rHDAPA = k HDAPAPH2 n1⎜⎜C PAL − DAL ⎟⎟ K HDAPA ⎠ ⎝

(12)

⎛ ⎞ C −rHDADA = k HDADAPH2 n2⎜⎜C DAL − MAL ⎟⎟ K HDADA ⎠ ⎝

(13)

⎛ C NAL ⎞ ⎟ −rHDAMA = k HDAMAPH2 n3⎜⎜CMAL − K HDAMA ⎟⎠ ⎝

(14)

Figure 1. Experimental and predicted LCO vaporization at different temperatures (P = 7 MPa and H2/oil = 5200 scf/bbl).

The rate parameters used in eqs 11−14 were taken from the literature22,23 and are reported in Table 1.

the experimental and predicted LCO vaporization as a function of temperature. It has been shown that this correlation can reasonably predict VLE for any gas-oil hydrotreating system as long as the total aromatics content is known.

Table 1. Reaction Kinetic Parameters of LCO Hydrotreating Reactions reaction HDS

HDA

parameters

3. DESCRIPTION OF THE PROCESS Figure 2 presents a flow diagram of the hydrotreating unit. The oil feed is first combined with a hydrogen-rich gas stream, and the resulting mixture is heated to the required reaction temperature. The oil-gas stream then passes through a fixedbed reactor packed with hydrotreating catalyst. Depending on the amount of heat release, the reactor might have multiple catalyst beds to inject quenching gas between pairs of beds. The reactor effluent is sent to a high-pressure separator (HPS) where the total liquid product (TLP) is separated from the gases. The gas stream is scrubbed with diethylamine (DEA) to

m1 = 1.6, m2 = 0.56 kHDS = 2.5 × 109 exp(−19384/T), KH2S = 5.0 m3/kmol −ΔHrHDS = −60.3 kJ/(mol of H2) n1 = 0.5, n2 = 0.5, n3 = 1.0 kHDAPA = 2.66 × 105 exp(−15170/T), KHDAPA = 1.0573 × 10−5 exp(8308/T) kHDADA = 8.50 × 102 exp(−12140/T), KHDADA = 5.5660 × 10−11 exp(15741/T) kHDAMA = 6.04 × 102 exp(−12414/T), KHDAMA = 7.4928 × 10−17 exp(24070/T) −ΔHrHDA = −67.0 kJ/(mol of H2) C

dx.doi.org/10.1021/ie501032j | Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Industrial & Engineering Chemistry Research

Article

distribution, first bed = 2.5 m, second bed = 3.5 m, third bed = 6 m; liquid hourly space velocity (LHSV) = 1.5 h−1; feed temperature = 330−340 °C; quench gas temperature = 80 °C; pressure = 5−9 MPa; H2/oil ratio = 5000 scf/bbl. The catalyst bed distribution was pre-established by trial-anderror simulations. The final configuration was selected based on the uniformity of the resulting axial temperature profiles.24 In this case, the reactor was determined to have three catalyst beds of increasing length. The specified H2/oil ratio refers to the total amount of hydrogen fed to the reactor, including quenching gas. Because the recycle gas composition and the ratio between the gas fed to the top of the reactor and the quenching gas were initially unknown, guessed values were carefully selected to start the simulations. Similarly, the treat gas composition was initially assumed to be 100% H2. These variables were recalculated by successive simulations until the hydrogen balance in the gas recycle loop was closed. In most cases, the converged treat gas composition was 96−98 mol % H2 and 4−2 mol % H2S.

4. RESULTS AND DISCUSSION 4.1. Simulations with and without VLE. To demonstrate the impact of VLE on the modeling and simulation of hydroprocessing reactors, simulations were conducted with and without accounting for VLE. For the simulations without VLE, the flash program was simply unplugged from the reactor model. The comparative simulation results are presented in Figures 3−9. Figure 3 presents the axial reactor temperature

Figure 2. Flow diagram of the hydrotreating process.

remove most of the hydrogen sulfide and ammonia generated in the reactor. The resulting hydrogen stream is recompressed and recycled to the reaction system. The feedstock was assumed to be pure LCO even though commercial gas-oil hydrotreater feeds are typically blends of various refinery steams, such as straight-run light gas-oil (LGO), cocker gas-oil, and LCO from the catalytic cracking process. The properties of the LCO feed are listed in Table 2. Table 2. Properties of the Feedstock property density at 15.6 °C, g/cm3 C, wt % H, wt % S, wt % saturates, wt % paraffins naphthenes aromatics, wt % monoaromatics diaromatics polyaromatics SimDis, °C IBP 10 wt % 30 wt % 50 wt % 70 wt % 90 wt % FBP

value 0.9360 88.25 10.54 1.17 15.7 6.8 8.9 83.8 17.7 52.3 13.8

Figure 3. Temperature profiles with and without VLE at T0 = 330 °C, LHSV = 1.5 h−1, P = 5 MPa, and H2/oil = 5000 scf/bbl.

profiles. It can be observed that temperature increases from inlet to outlet of the reactor due to the exothermic hydrogenation reactions. The two sudden temperature drops at 2.5 and 6 m correspond to the interbed quenching zones where the hot bed effluents are mixed with quenching gas. The amount of quenching gas for this set of simulations was adjusted to achieve the same bed inlet temperatures (T0 = 330 °C). When VLE is ignored, the model predicts a significantly higher total temperature rise (ΔT) than that when VLE is considered (83 vs 70 °C). This is because the heat absorbed by vaporization of oil is not considered in the heat balance equation when VLE is ignored. This observation is consistent with that reported by Murali et al.7 The remarkable ΔT values in both cases are due to the extensive hydrogenation of aromatics in the LCO feed (83.8 wt %). Figure 4 shows a comparison of the superficial liquid-phase velocity (uL) profiles. Marked differences can be observed between the two cases. In traditional non-VLE reactor

128.4 225.0 255.3 288.5 326.4 373.9 430.8

LCOs are highly unsaturated, making them highly reactive under hydrogenation conditions. Dilution with less aromatic feeds is therefore required in industrial practice to avoid reactor overheating. However, for the purposes of this study, LCO serves as an ideal example to illustrate the VLE behavior of extremely exothermic reaction systems. The reactor geometry and operating conditions used for the simulations were as follows: total catalyst bed volume = 115.5 m3; catalyst bed diameter = 3.5 m; catalyst bed length D

dx.doi.org/10.1021/ie501032j | Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Industrial & Engineering Chemistry Research

Article

sumption is faster, primarily due to the higher temperature, as discussed in Figure 3. Interbed quenching acts as a make-up source of hydrogen, resulting in a jump in H2/oil ratio at the entrance of the catalyst bed. Higher bed ΔT values, in the case of the simulation without VLE, increase the size of this jump because more quench gas is required to reach the specified inlet temperature. This also reduces the H2/oil ratio at the reactor inlet (∼2560 vs ∼3100 scf/bbl), as quenching gas is withdrawn from the gas stream that is eventually fed through the top of the reactor. Figure 7 shows the HDS and HDA conversion profiles. The simulation without VLE produces a slightly higher HDS

Figure 4. Liquid superficial velocity profiles with and without VLE at T0 = 330 °C, LHSV = 1.5 h−1, P = 5 MPa, and H2/oil = 5000 scf/bbl.

modeling, uL remains practically constant throughout the reactor because all of the hydrocarbons are assumed to be in the liquid phase. With VLE, on the other hand, uL drops at the inlet of the reactor and keeps decreasing along the catalyst bed. Such a pattern is explained by the initial flashing of the feed at the inlet of the reactor followed by the additional vaporization caused by the progressive reactor temperature increase, as illustrated in Figure 5. Also indicated in this figure is that

Figure 7. HDS and HDA conversion profiles with and without VLE at T0 = 330 °C, LHSV = 1.5 h−1, P = 5 MPa, and H2/oil = 5000 scf/bbl.

conversion (94% vs 91%) as a result of higher reactor temperatures. However, the extent of HDA without VLE is lower than that with VLE (46% vs 51%). Such behavior can be explained by analyzing the liquid-phase aromatics concentration profiles, as shown in Figure 8. Whereas PA and DA disappear as

Figure 5. Vaporization profiles T0 = 330 °C, LHSV = 1.5 h−1, P = 5 MPa, and H2/oil = 5000 scf/bbl.

quenching produces a step increase in uL and a drop in oil vaporization at axial positions of 2.5 and 6 m due to the partial condensation of the hydrocarbons in the vapor phase. The behavior of the H2/oil ratio along the reactor is presented in Figure 6. Axial H2/oil ratio profiles usually have inverse shapes with respect to the reactor temperature profiles. The H2/oil ratio tends to decrease along the catalyst bed due to hydrogen consumption. Without VLE, the hydrogen con-

Figure 8. Liquid-phase aromatics concentration profiles with and without VLE at T0 = 330 °C, LHSV = 1.5 h−1, P = 5 MPa, and H2/oil = 5000 scf/bbl.

the hydrogenation reactions progress, MA formation by hydrogenation of DA exceeds their conversion into NA along the first half of the reactor. However, the balance eventually shifts toward NA formation. In the case with VLE, large amounts of aromatics are vaporized into the gas phase, particularly volatile species such as NA and MA. Without the vaporization effect, all hydrocarbons stay in the liquid phase, which is reflected in the observed higher concentration values. This difference, in fact, affects the hydrogenation rates, as higher product concentrations favor the reverse reactions (dehydrogenation), resulting in lower overall HDA conversions.

Figure 6. H2/oil profiles with and without VLE at T0 = 330 °C, LHSV = 1.5 h−1, P = 5 MPa, and H2/oil = 5000 scf/bbl. E

dx.doi.org/10.1021/ie501032j | Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Industrial & Engineering Chemistry Research

Article

temperature profiles. Quenching gas rates were adjusted to have either equal bed inlet temperatures, as in the simulations with T0 set at 330 and 340 °C (Figure 10a), or increasing bed inlet temperatures, as in the case where T0 for the first catalyst bed is 335 °C (Figure 10b). For T0 = 330 °C, it can be observed that reactor temperature is properly distributed over the three catalyst beds, with practically the same ΔT value for each bed (∼23 °C) and, therefore, equal average temperatures. Ideally, a hydrotreating reactor should have equal ΔT values per bed, somewhere between 15 and 30 °C, to make more efficient use of the total catalyst inventory.25 When T0 is increased to 340 °C, the heat release becomes stronger toward the top of the reactor, resulting in a less favorable temperature distribution: The first catalyst bed exhibits a sharp ΔT of ∼27 °C, whereas in the last bed, ΔT decreases to about 20 °C. This type of temperature profile increases the risk of local overheating in the front end of the reactor26 and also affects catalyst performance in the beds below because of the reduced average temperatures. As for the simulation with increasing T0 (Figure 10b), the ΔT values in the first two beds are fairly similar (∼25 °C), but ΔT then decreases to about 20 °C in the last bed. However, in contrast to the case with T0 = 340 °C (Figure 10a), the average temperature actually increases along the reactor, and the exotherm in the first bed is less pronounced. Figure 11 presents the corresponding profiles of the H2/oil ratio. Significant differences can be observed between the

The evolution of the gas-phase composition along the catalyst beds with VLE is illustrated in Figure 9. The simulation

Figure 9. Gas-phase concentration profiles with VLE at T0 = 330 °C, LHSV = 1.5 h−1, P = 5 MPa, and H2/oil = 5000 scf/bbl.

confirms that, under typical hydroprocessing conditions, a variety of hydrocarbon components are present in the gas phase, as opposed to a traditional simulation without VLE. It is observed that the reactor inlet already exhibits noticeable amounts of hydrocarbons, mostly DA and MA, due to the initial flashing of the feed. Products, such as H2S and NA, are progressively released into the gas phase as the reactions move forward. The H2 concentration decreases rapidly in each catalyst bed as a result of hydrogenation reactions. Gas quenching produces a step increase in H2 concentration, which is highly beneficial for maintaining a H2-rich atmosphere throughout the reactor. 4.2. Effect of Temperature. Simulations with VLE were performed at different inlet temperatures to evaluate the effect of this variable on reactor performance. The results are reported in Figures 10−13. Figure 10 presents the resulting

Figure 11. H2/oil profiles at various inlet temperatures with LHSV = 1.5 h−1, P = 5 MPa, and H2/oil = 5000 scf/bbl.

simulations, especially in the first and second catalyst beds. This is essentially related to the amount of quenching gas used in each case, as discussed above. For T0 = 330 and 340 °C, the volume of quenching gas required for maintaining equal bed inlet temperatures reduces the H2/oil ratio at the top of the reactor to about 3100 and 3000 scf/bbl, respectively. For T0 = 335 °C, on the other hand, the H2/oil ratio at the top of reactor is significantly higher (3430 scf/bbl) because less quenching gas is withdrawn from the recycle gas stream. The variation of the step increase in H2/oil ratio between catalyst beds is also determined by the rate of quenching gas. Figure 12 shows the corresponding HDS and HDA conversion profiles. In general, total HDS conversion exceeds 90%, whereas HDA conversion is above 50%. There is a small gain in overall HDS and HDA conversions (4−5% and 6−5%, respectively) when T0 is increased from 330 to either 340 or 335 °C. It is noted that the conversion profiles with T0 = 340 and 335 °C practically overlap despite the differences in temperature distribution. Whereas the average temperature in the first catalyst bed is higher at T0 = 340 °C than at T0 = 335 °C, the overall H2/oil ratio is significantly lower, particularly at

Figure 10. Temperature profiles at various inlet temperatures with LHSV = 1.5 h−1, P = 5 MPa, and H2/oil = 5000 scf/bbl: (a) equal bed inlet temperatures and (b) increasing bed inlet temperatures. F

dx.doi.org/10.1021/ie501032j | Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Industrial & Engineering Chemistry Research

Article

Figure 14. Temperature profiles at various pressures with T0 = 330 °C, LHSV = 1.5 h−1, and H2/oil = 5000 scf/bbl.

Figure 12. HDS and HDA conversion profiles at various inlet temperatures with LHSV = 1.5 h−1, P = 5 MPa, and H2/oil = 5000 scf/ bbl.

generation. High pressures enhance hydrogen concentration in the liquid phase, thus increasing reaction rates. Such behavior is clearly observed in the temperature and conversion profiles presented in Figures 14 and 15, respectively. When pressure is

the reactor inlet (3000 vs 3430 scf/bbl), as seen in Figures 10 and 11. This indicates that the effect of temperature on conversion is counterbalanced by the effect on H2/oil ratio. A decreased H2/oil ratio directly affects catalytic performance because there is less hydrogen available for promoting hydrotreating reactions. In this study, this dependence is properly taken into account in the rate equations by considering the effect of hydrogen concentration in explicit form. Reactor temperature also influences hydrocarbon vaporization profiles, as illustrated in Figure 13. It can be observed

Figure 15. HDS and HDA conversion profiles at various pressures with T0 = 330 °C, LHSV = 1.5 h−1, and H2/oil = 5000 scf/bbl.

increased from 5 to 9 MPa, the total temperature rise increases dramatically from 70 to 100 °C, whereas HDS conversion increases from 91% to 96% and HDA conversion improves markedly from 51% to 66%. In this sense, HDA is particularly favored by higher pressures because increased pressure shifts the hydrogenation/dehydrogenation chemical equilibrium toward hydrogenation reactions. Increasing reactor pressure also reduces oil vaporization, as can be seen in Figure 16. The oil vaporization in the initial flash of the feed drops from 30 wt % at 5 MPa to about 15 wt % at 9 MPa. In the latter case, the

Figure 13. Vaporization profiles at various inlet temperatures with LHSV = 1.5 h−1, P = 5 MPa, and H2/oil = 5000 scf/bbl.

that feed vaporization is a strong function of process temperature. At the top of the reactor, about 30−36 wt % of the LCO is vaporized under the conditions specified, depending on the inlet temperature. Further vaporization takes place as heat is released from the hydrogenation reactions, reaching over 50 wt % of the LCO toward the bottom of the reactor. This is because the hydrocarbon stream becomes more volatile as the aromatics content decreases along the reactor. This is achieved by updating the binary interaction coefficients based on the actual aromatics content, as described in section 2.3. Analogously to the conversion profiles (Figure 12), the vaporization profiles with T0 = 340 and 335 °C also overlap, but only in the first catalyst bed. In this case, the H2/oil ratio also counteracts the effect of temperature, as oil vaporization decreases at lower H2/oil ratios. 4.3. Effect of Pressure. Figures 14−16 show the reactor simulation results at different pressures with VLE taken into account. In general, increasing reactor pressure improves hydrotreating performance, which comes with more heat

Figure 16. Vaporization profiles at various pressures with T0 = 330 °C, LHSV = 1.5 h−1, and H2/oil = 5000 scf/bbl. G

dx.doi.org/10.1021/ie501032j | Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Industrial & Engineering Chemistry Research

Article

robust modeling tools. This work presents the detailed modeling and simulation of an industrial LCO hydrotreater with rigorous VLE calculations. The model consisted of a onedimensional plug-flow adiabatic reactor model coupled with a calibrated flash program. The results led to the following conclusions: • Simulations with and without VLE revealed significant differences in reactor temperature profiles, hydrotreating conversions, fluid flow rates, and phase compositions. Such differences clearly demonstrate the importance of VLE in the formulation of the reactor model equations and estimation of rate parameters. Ignoring VLE consequently leads to an unrealistic representation of hydrotreating processes and reduces model predictive capability as a result of the implicit dependence of rate parameters on VLE. • Whereas increasing the inlet temperature from 330 °C to either 335 or 345 °C improved HDS and HDA to a certain extent (4−5% and 6−5%, respectively), it also increased the quenching gas requirement and initial vaporization rate (from 30 to 36 wt %) and promoted the formation of undesired high-temperature zones at the bottom of the first and second catalyst beds (>360 °C). • Increasing reactor pressure from 5 to 9 MPa improved hydrotreating performance (HDS increased from 91% to 96%, whereas HDA increased from 51% to 66%) and significantly reduced feed vaporization at the reactor inlet from 30 to 15 wt %. Nevertheless, this came with excess heat generation, which was reflected in a dramatic total ΔT increase from 70 to 100 °C. • Analysis of plug-flow and full-catalyst-wetting conditions confirmed that all criteria were satisfied in the simulations, even when vaporization exceeded 50 wt %.

amount of vaporized oil stays below 35 wt % in the entire reactor (three catalysts beds). 4.4. Verification of Plug-Flow and Full-CatalystWetting Conditions. Plug-flow deviations and incomplete catalyst wetting are phenomena that typically occur in tricklebed reactors, especially when feed vaporization is high. Therefore, the assumptions of plug-flow operation and fully wetted catalyst beds were verified using several theoretical criteria reported in the literature.27 The analysis was based on the simulation that had the highest vaporization level, as presented in Figure 13, for which the following conditions were used: T0 = 335 °C, LHSV = 1.5 h−1, P = 5 MPa, and H2/oil = 5000 scf/bbl. Under these conditions, the liquid flow rate decreases the most, and consequently, the plug-flow and fullcatalyst-wetting criteria are more difficult to meet. Figure 17 presents the actual Peclet number (Pe) of the liquid phase and the minimum Peclet number required to meet

Figure 17. Plug-flow verification at T0 = 335 °C, LHSV = 1.5 h−1, P = 5 MPa, and H2/oil = 5000 scf/bbl.



the plug-flow criterion at any position of the reactor. It is observed that the actual Pe value far exceeds the minimum required to meet the plug-flow criterion, confirming that the plug-flow assumption is indeed valid for the conditions used in this study. Figure 18 shows the actual wetting number (W) and

AUTHOR INFORMATION

Corresponding Author

*Tel.: 780-987-8763. Fax: 780-987-5349. E-mail: jichen@nrcan. gc.ca. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS Partial funding for this study was provided by the Canadian Interdepartmental Program of Energy Research and Development (PERD 1.1.3). Comments and suggestions from Dr. Edward Little and Dr. Huali Wang on revising the manuscript are greatly appreciated.



Figure 18. Full-catalyst-wetting verification at T0 = 335 °C, LHSV = 1.5 h−1, P = 5 MPa, and H2/oil = 5000 scf/bbl.

the minimum wetting number required to maintain fullcatalyst-wetting conditions at any position of the reactor. The results clearly indicate that the full-catalyst-wetting criterion is also satisfied, even in the last catalyst bed where vaporization exceeds 50 wt %.

5. CONCLUSIONS Fundamental understanding and accurate prediction of VLE in hydrotreating reactors are important aspects for developing H

NOMENCLATURE A = correlation parameter defined in eq 15 aL = gas−liquid interfacial area, m2/m3 AR = cross-sectional area of the reactor, m2 B = correlation parameter defined in eq 15, 1/°C BP = boiling point, °C C = molar concentration, kmol/m3 Caromatics% = aromatics content in the oil feed, wt % Cp = heat capacity, kJ/kg·K diH = interaction coefficient between hydrogen and hydrocarbon pseudocomponents F = molar flow rate, kmol/s f = mass flow rate, kg/s G0L = superficial liquid mass velocity of the feed, kg/m2·s dx.doi.org/10.1021/ie501032j | Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Industrial & Engineering Chemistry Research



H = phase enthalpy, kJ/kmol KH2S = H2S adsorption equilibrium constant, m3/kmol KHDA = hydrogenation equilibrium constant kHDA = forward hydrogenation reaction rate constant, m3/kg· s kHDS = HDS reaction rate constant, (m3)2.16/kg·kmol1.16·s K* = VLE constant kL = gas−liquid mass transfer coefficient, m/s m = HDS reaction order n = HDA reaction order P = pressure, MPa Pe = Peclet number r = reaction rate, kmol/kg·s T = temperature, °C or K u = superficial velocity, m/s W = wetting number z = reactor axial position, m

Article

REFERENCES

(1) Al-Dahhan, M. H.; Larachi, F.; Dudukovic, M. P.; Laurent, A. High-Pressure Trickle-Bed Reactors. Ind. Eng. Chem. Res. 1997, 36, 3292−3314. (2) Chen, J.; Ring, Z. Vapor−Liquid Phase Equilibrium during Hydrotreating of Light Cycle Oil. Presented at the 2004 AIChE National Spring Meeting, New Orleans, LA, Apr 25−29, 2004. (3) Chen, J.; Jiang, W.; Yang, H.; Hawkins, R.; Ring, Z. Comparative Study of Vapor−Liquid Equilibrium during Hydroprocessing of Different Petroleum Feedstocks. Presented at the 2006 AIChE National Spring Meeting, Orlando, FL, Apr 23−27, 2006. (4) Chen, J.; Wang, N.; Mederos, F.; Ancheyta, J. Vapor−Liquid Equilibrium Study in Trickle-Bed Reactors. Ind. Eng. Chem. Res. 2009, 48, 1096−1106. (5) Munteanu, M. C.; Chen, J. Vapor−Liquid Equilibrium (VLE)Based Modeling for the Prediction of Operating Regimes in a Heavy Gas Oil Hydrotreater. Energy Fuels 2012, 26, 1230−1236. (6) Hoekstra, G. The Effect of Gas-to-Oil Rate in Ultra Low Sulphur Diesel Hydrotreating. Catal. Today 2007, 127, 99−102. (7) Murali, C.; Voolapalli, R. K.; Ravichander, N.; Gokak, D. T.; Choudary, N. V. Trickle Bed Reactor Model to Simulate the Performance of Commercial Diesel Hydrotreating Unit. Fuel 2007, 86, 1176−1184. (8) Mederos, F. S.; Elizalde, I.; Ancheyta, J. Steady-State and Dynamic Reactor Models for Hydrotreatment of Oil Fractions: A Review. Catal. Rev.: Sci. Eng. 2009, 51, 485−607. (9) Kocis, G. R.; Ho, T. C. Effects of Liquid Evaporation on the Performance of Trickle-Bed Reactors. Chem. Eng. Res. Des. 1986, 64, 288−291. (10) Avraam, D. G.; Vasalos, I. A. HdPro: A Mathematical Model of Trickle-bed Reactors for the Catalytic Hydroprocessing of Oil Feedstocks. Catal. Today 2003, 79−80, 275−283. (11) Kumar, H.; Froment, G. F. Mechanistic Kinetic Modeling of the Hydrocracking of Complex Feedstocks, such as Vacuum Gas Oils. Ind. Eng. Chem. Res. 2007, 46, 5881−5897. (12) Pellegrini, L. A.; Gamba, S.; Calemma, V.; Bonomi, S. Modelling of Hydrocracking with Vapor-Liquid Equilibrium. Chem. Eng. Sci. 2008, 63, 4285−4291. (13) Möller, K.; le Grange, P.; Acolla, C. A Two-Phase Model for the Hydrocracking of Fischer−Tropsch-Derived Wax. Ind. Eng. Chem. Res. 2009, 48, 3791−3801. (14) Schweitzer, J.-M.; Lopez-Garcia, C.; Ferre, D. Thermal Runaway Analysis of a Three-Phase Reactor for LCO Hydrotreatment. Chem. Eng. Sci. 2010, 65, 313−321. (15) Lopez Garcia, C.; Hudebine, D.; Schweitzer, J.-M.; Verstraete, J. J.; Ferre, D. In-Depth Modeling of Gas Oil Hydrotreating: From Feedstock Reconstruction to Reactor Stability Analysis. Catal. Today 2010, 150, 279−299. (16) Chen, J.; Mulgundmath, V.; Wang, N. Accounting for Vapor− Liquid Equilibrium in the Modeling and Simulation of a Commercial Hydrotreating Reactor. Ind. Eng. Chem. Res. 2011, 50, 1571−1579. (17) Korsten, H.; Hoffmann, U. Three-Phase Reactor Model for Hydrotreating in Pilot Trickle-Bed Reactors. AIChE J. 1996, 42, 1350−1360. (18) Alvarez, A.; Ancheyta, J. Modeling Residue Hydroprocessing in Multi-Fixed-Bed Reactor System. Appl. Catal. A 2008, 351, 148−158. (19) Goto, S.; Smith, J. M. Trickle-Bed Reactor Performance. Part I. Holdup and Mass Transfer Effects. AIChE J. 1975, 21, 706−713. (20) Girgis, M. J.; Gates, B. C. Reactivities, Reaction Networks, and Kinetics in High-Pressure Catalytic Hydroprocessing. Ind. Eng. Chem. Sci. 1991, 30, 2021−2058. (21) Ali, S. A. Thermodynamics of Hydroprocessing Reactions. In Hydroprocessing of Heavy Oils and Residua; Ancheyta, J., Speight, J. G., Eds.; CRC Press: Boca Raton, FL, 2007; Chapter 4. (22) Chowdhury, R.; Pedernera, E.; Reimert, R. Trickle-Bed Reactor Model for Desulfurization and Dearomatization Diesel. AIChE J. 2002, 48, 126−135. (23) Jaffe, S. B. Kinetics of Heat Release in Petroleum Hydrogenation. Ind. Eng. Chem. Process Des. Dev. 1974, 13, 34−39.

Subscripts

0 = inlet DA = diaromatics G = gas phase H2 = hydrogen H2S = hydrogen sulfide i = component i j = reaction j L = liquid phase MA = monoaromatics NA = naphthenes out = outlet PA = polyaromatics Q = quenching gas S = sulfur Superscripts

* = VLE Greek Symbols

ΔHr = heat of reaction, kJ/(mol of H2) ΔT = temperature rise, °C ν = molar stoichiometric coefficient ρ = phase density, kg/m3 ρb = catalyst bulk density, kg/m3 ϕ = vaporized fraction Abbreviations

DA = diaromatics DEA = dethylamine EoS = equation of state H2/oil = H2-to-oil ratio HDA = hydrodearomatization HDS = hydrodesulfurization HPS = high-pressure separator LCO = light cycle oil LGO = light gas oil LHSV = liquid hourly space velocity MA = monoaromatics NA = naphthenes ODE = ordinary differential equation PA = polyaromatics PDE = partial differential equation S = organic sulfur compound TLP = total liquid product VLE = vapor−liquid equilibrium I

dx.doi.org/10.1021/ie501032j | Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Industrial & Engineering Chemistry Research

Article

(24) Alvarez, A.; Ancheyta, J.; Centeno, G.; Marroquín, G. A Modeling Study of the Effect of Reactor Configuration on the Cycle Length of Heavy Oil Fixed-Bed Hydroprocessing. Fuel 2011, 90, 3551−3560. (25) Gruia, A. Recent Advances in Hydrocracking. In Practical Advances in Petroleum Processing; Hsu, C. S., Robinson, P. R., Eds.; Springer: New York, 2006; Vol. 1, Chapter 8. (26) Alvarez, A.; Ancheyta, J. Transient Behavior of Residual Oil Front-End Hydrodemetallization in a Trickle-Bed Reactor. Chem. Eng. J. 2012, 198, 204−214. (27) Mederos, F. S.; Ancheyta, J.; Chen, J. Review on Criteria to Ensure Ideal Behaviors in Trickle-Bed Reactors. Appl. Catal. A 2009, 355, 1−19.

J

dx.doi.org/10.1021/ie501032j | Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX