Article pubs.acs.org/IECR
Development of Kinetic and Process Models for the Oxidative Desulfurization of Light Fuel, Using Experiments and the Parameter Estimation Technique Amer T. Nawaf, Aysar T. Jarullah,* and Saba A. Gheni Chemical Engineering Department, College of Engineering, Tikrit University, Tikrit, Salahuddin, Iraq
Iqbal M. Mujtaba* Chemical Engineering Division, University of Bradford, Bradford BD7 1DP, United Kingdom ABSTRACT: The oxidative desulfurization (ODS) of light gas oil (LGO) is investigated with an in-house-designed catalyst consisting of cobalt oxide loaded on alumina (γ-Al2O3), in the presence of air as an oxidizing agent, under moderate operating conditions (temperatures of 403−473 K, liquid hourly space velocity (LHSV) = 1−3 h−1, initial concentration = 500−1000 ppm). The Incipient Wetness Impregnation (IWI) method of cobalt oxide over γ-alumina (2% Co3O4/γ-Al2O3) is used for the preparation of the catalyst. The optimal design of experiments is studied to evaluate the effects of several process variables (namely, temperature, LHSV, and concentration of dibenzothiophene) and their optimal values were found to be 473 K, 1 h−1, and 1000 ppm, respectively. For the conversion of dibenzothiophene to sulfone and sulfoxide, the results indicates that the IWI method is suitable to prepare this type of the catalyst. Based on the experiments, mathematical models that represent a threephase reactor for describing the behavior of the ODS process are developed. In order to develop a useful model for simulation, control, design, and scaleup of the oxidation process, accurate evaluation of important process parameters, such as reaction rate parameters, is absolutely essential. For this purpose, the parameter estimation technique available in the general Process Modeling System (gPROMS) software is employed in this work. With the estimated process parameters, further simulations of the process is carried out and the concentration profiles of dibenzothiophene within the reactor are generated.
1. INTRODUCTION Sulfur compounds (mainly benzothiophene (BT), and dibenzothiophene (DBT) and its derivatives) in oil (fuels) are the main source of air pollution, because of the generation of sulfoxides via the combustion process, which leads to acid rain. The traditional mode of sulfur removal in fuels is via catalytic hydrodesulfurization (HDS), commonly known as hydrotreating, and requires modified catalyst and severe operating conditions (temperature, pressure, etc.). This make the hydrotreating process more expensive, compared to other processes.1−3 Thiophene and its aromatic compounds are the main sulfur components found in the oil feedstock. Decreasing the amount of sulfur in fuel has gained significance, because of increasing awareness about the serious consequences of burning sulfur-bearing fuels. The U.S. Environmental Protection Agency (EPA) had instituted new sulfur standards for diesel fuels and gasoline.1 The oxidative desulfurization (ODS) process is regarded as one of the most promising alternative deep desulfurization operations to get ultralow sulfur fuels.1,4,5 The ODS process of sulfur compounds, such as thiophene (Th), BT, DBT, and their compounds, are investigated by employing various solid catalysts, such as Mo-Al2O36,7 and cobalt−aluminum phosphate,8 where the sulfur components are oxidized into their corresponding sulfoxides and sulfones. The greatest advantage of ODS, compared to the HDS process, is that the ODS can be conducted in a liquid phase under moderate conditions. In ODS reactions, the sulfur compound is oxidized by adding oxygen molecules to form the hexavalent sulfur of sulfones.9 © 2015 American Chemical Society
The idea of ODS is actually quite simple. Sulfur compounds are known to be slightly more polar than hydrocarbons of similar structure.10,11 However, oxidized sulfur compounds such as sulfones or sulfoxides are substantially more polar than unoxidized sulfur compounds. This permits the selective removal of sulfur compounds from hydrocarbon via a combination process of selective oxidation and solvent extraction or solid adsorption.12 Before 1980, the most popular oxidants in the study of ODS are nitric acid and nitrogen oxides and are used largely because they have double effects of oxidizing sulfur compounds and nitrating the aromatic compounds to form nitro aromatics with high cetane numbers. However, it has major drawbacks, such as poor selectivity, low yield, and loss of heating value of the treated oil.12,13 Other types of oxidants have also been used, including H2O2/AcOH, H2O2/H2SO4, O3, KMnO4 and BuOOH,14−16 oxygen,17 and O2/aldehyde/cobalt catalysts.18 Three-phase reactors with a fixed bed of catalyst and a cocurrent downflow of gas and liquid are widely utilized in different oil, petrochemical, and chemical industries, besides water industries treating wastewater. Understanding the phenomena that govern the performance of three-phase reactors has played a significant role in the designing of such equipment. The hydrodynamic factors, such as pressure drop, liquid holdup, and Received: Revised: Accepted: Published: 12503
September 4, 2015 November 16, 2015 November 23, 2015 November 23, 2015 DOI: 10.1021/acs.iecr.5b03289 Ind. Eng. Chem. Res. 2015, 54, 12503−12515
Article
Industrial & Engineering Chemistry Research Table 1. List of Chemicals and Materials Used for Catalyst Preparation and Aluminum Oxide Specifications chemicals and materials
purity (%)
cobalt nitrate, Co(NO3)2·6H2O deionized water
99.5
function
manufacturer
active material solvent of active material Aluminum Oxide (γ-Al2O3) Specifications
Alpha Chemical Samarra Company
catalyst
pore volume
bulk density
surface area
particle diameter
particle shape
γ-Al2O3
0.5367 cm3/g
0.671 g/cm3
289 m2/g
1.6 mm
spherical
with continuous stirring, until impregnation of the entire solution is complete. The temperature is kept constant at 373 K, using a bath water. The impregnated γ-alumina is then left to dry overnight in the furnace at 393 K. The objective of this step is to eliminate water. The calcination step is then applied for 5 h in the furnace at 823 K with air. This step converts the metal salts deposited on the γ-Al2O3 into metal oxides, allowing for deposition of active metal oxides on the catalyst support, and, thus, the desired physical and chemical specifications of the catalyst are achieved. The calcination step is conducted at the Fertilizer/Northern Company−Baiji. Figure 1 illustrates the steps and sequence of activities in the catalyst preparation.
catalyst wetting efficiency, together with characterization of reaction kinetics, as well as transport in catalyst particles, are all important and should be considered for developing an accurate model of the process. Plug flow models for the liquid phase with modified external liquid holdup, as well as external contacting catalyst effectiveness parameters, have been suggested by several investigators in the past.19,20 Based on experimental studies with an in-house-designed catalyst, the objective of this study is to develop kinetic models for the ODS process. For this purpose, a full process model available in the public domain is used and the reaction parameters of the model are determined by minimizing (optimization) sum of the squared error between the data obtained experimentally and those predicted by the model. The modeling, simulation, and optimization process of ODS operation are carried out by employing gPROMS software.21
2. EXPERIMENTAL WORK 2.1. Feedstock (Light Gas Oil). Light gas oil (LGO), the feedstock used in this study, has been provided by the North Refineries (Iraq) with the following specifications: specific gravity (Sp gr), 0.851; viscosity at 293 K, 4.9 cSt; flash point, 55 °C; total sulfur, 9.8 ppm; cetane index, 52; and pour point, −39 °C; these specification were determined through testing in the laboratories at North Refineries Company. 2.2. Dibenzothiophene (DBT). The dibenzothiophene (DBT) obtained from Aldrich is chosen to study the reactivity of sulfur in the oxidation reaction. The purity of the sulfur compound is ∼98%. 2.3. Air. Air is used as an oxidizing agent. The oxygen that is contained in air can act as an oxidant to sulfur compounds. 2.4. Catalyst. Chemical compounds are used for catalyst preparation, as follow. 2.4.1. Active Compound Used in the Catalyst Preparation. Specifications of the active compound used in the catalyst preparation are shown in Table 1. 2.4.2. Supported Alumina Oxide (γ-Al2O3). The specification of a commercial spherical-particle alumina oxide-type γ-alumina is also presented in Table 1, which has been used as a carrier in the manufacturing of the catalyst. 2.5. Catalyst Preparation Used in Experimental Work. The cobalt solution (cobalt nitrate supported on alumina (γ-Al2O3)) has been obtained using the Incipient Wetness Impregnation (IWI) method. The preparation procedure is as follows. First, 100 g of the alumina is dried in the furnace at 393 K for 4 h to remove the moisture from alumina before impregnation. Second, 2.1 g of cobalt nitrate is added to 40 cm3 of deionized water (the pore volume of γ-alumina is equal to the volume of deionized water), while the solution is being stirred (using magnetic stirrer) for 1 h at room temperature. The pretreated γ-alumina in the first step is then transferred into a flask under vacuum (using a vacuum pump) for removing gases out of support pores. After that, the solution obtained in the second step is added to γ-alumina at a rate of 15−20 drops/min,
Figure 1. Flowchart of preparation steps.
2.6. Oxidation Operation in a Trickle Bed Reactor. 2.6.1. Apparatus and Procedure. The experiments have been carried out at high temperature and pressure in a trickle bed reactor (TBR) available in Tikrit University (Iraq). TBRs are extensively utilized in oil refineries, and the process flow diagram of this system is presented in Figure 2. Three phases exist in such systems: a solid phase (catalyst bed), a gas phase (air), and a liquid phase (LGO). The continuous oxidation of LGO is carried out in the TBR where the LGO and air are fed in co-current mode. The reactor is made of stainless steel with 12504
DOI: 10.1021/acs.iecr.5b03289 Ind. Eng. Chem. Res. 2015, 54, 12503−12515
Article
Industrial & Engineering Chemistry Research
Figure 2. Process flow diagram of the trickle bed reactor (TBR).
an inside diameter of 1.6 cm and length of 77 cm. Four steel-jacket heaters of equal length are used to control the reactor temperature. The top and bottom parts (30%−35% (by volume) at each end) of the reactor are filled with inert particles to serve as the disengaging part. The inner part (40 vol %) of the reactor contains a packing of cobalt oxide catalyst.22 The LGO is pumped at pressure up to 20 bar with flow rates of 0.0−1.65 L/h. The oxidant gas (air) flows from an air compressor at high pressure and a fixed operation pressure is maintained. The LGO with varying DBT concentration is mixed with air before feeding into the reactor at the desired temperature, allowing DBT to oxidize to sulfones. The outlet from the reactor flows through a heat exchanger to high-pressure gas− liquid separator in order to separate excess air from the treated LGO. The description and specifications of the experimental equipment can be found in the work of Nawaf et al.38 2.6.2. Experimental Runs. The effect of operational parameters on the reactor performance of ODS using cobalt oxide (Co3O4/γ-Al2O3) catalyst is evaluated by varying the temperature between 403 K and 473 K and the liquid hourly space velocity (LHSV) between 1 h−1 and 3 h−1. The concentration of DBT is varied between 500 ppm and 1000 ppm. The oxidation experiments have been conducted in a trickle bed reactor that was packed with 40% catalyst particles under isothermal conditions. The model light gas oil is prepared by adding DBT to hydrotreated light gas oil (containing 2 ppm of DBT) with specified initial concentrations of DBT. The temperature of the LGO feed tank is controlled using a cooling jacket, where the coolant side temperature is maintained below 293 K, to prevent vaporization of the light compounds found in light gas oil. To prevent leaks and to remove any gases and liquid remained from the last run, nitrogen gas is passed through the reactor. LGO mixed with air is then passed through the reactor at a pressure of 2 bar, and the temperature controller is set to the
desired feed temperature. When the air temperature reached the feed injection temperature, the dosing pump was turned on to allow a certain flow of light gas oil, and the temperature is increased at a rate of 293 K per hour until a steady-state temperature is reached. At the end of a run, the LGO dosing pump is turned off, keeping the air gas flow on to backwash any remaining light gas oil. Finally, the air valve is shut off. 2.6.3. Sulfur Measurements (GC-Capillary Chromatography). The DBT concentration in the feed and the product are evaluated according to GC-capillary chromatography. The detailed specifications of the GC-capillary chromatography is given in the work of Nawaf et al.38
3. MATHEMATICAL MODEL OF TRICKLE BED REACTOR (TBR) FOR OXIDATIVE DESULFURIZATION (ODS) REACTION The process model plays a very important role in industry, from operator training, and health and safety, to design, operation, and control.23 Several investigators have suggested that pore diffusion should be considered within the reaction rate constant (which is determined by multiplying the intrinsic rate constant by an effectiveness factor) resulting in a pseudo-homogeneous basic plug flow model, which is adequate for describing the progress of chemical reactions in the liquid phase of a trickle bed reactor (TBR).20,24,25 3.1. Mass Balance Equations. Figure 3 shows a typical TBR with various features (assumption, operation parameters, software used, etc.). The general mass balance for a reactor can be described as ⎡ ⎤ mass disappearance [mass in] = [mass out] + ⎢ ⎥ ⎣ appearance by chemical reaction ⎦ + [accumulation] 12505
(1) DOI: 10.1021/acs.iecr.5b03289 Ind. Eng. Chem. Res. 2015, 54, 12503−12515
Article
Industrial & Engineering Chemistry Research
Figure 3. Required data and available tools for modeling and optimization oxidative desulfurization (ODS) reactions.
experimental observations.26 To appreciate the complexity of the chemical reaction, one should start with an nth-order kinetics:
The input of DBT (FDBT) is defined as the number of moles of DBT per unit time (FDBT = moles/time), and the output of DBT mol/time is expressed as FDBT + dFDBT. The disappearance of DBT is expressed as (−rDBT) dV (in terms of reaction moles per unit time), and the accumulation of dibenzothiophene is zero. FDBT = (FDBT + dFDBT) + ( −rDBT) dV
(2)
dFDBT = d[FDBT0(1 − XDBT)] = −FDBT0 dXDBT
(3)
−rDBT = −
K app = η0ηceK in
∫0
X DBTf
dXDBT −rDBT
(7)
where internal diffusion is described by the catalyst effectiveness factor (η0) and the hydrodynamics are described by the external catalyst wetting efficiency (ηce). The chemical reaction is produced as
(4)
−rDBT = −
The above-mentioned equation accounts for DBT compound in the differential part of catalyst volume (dV). For the catalytic reactor as a whole, the term should be integrated. Now, the feed rate (FDBT0) is fixed, but (−rDBT) certainly depends on the concentration or conversion of materials. τ = C DBT0
(6)
Apparent kinetics are related to the intrinsic kinetics regarding internal diffusion and TBR hydrodynamic influences as follows:27
Since FDBT = CDBTvL, where CDBT is the concentration of DBT (given in terms of moles per unit volume) and vL is the volumetric flow rate (expressed in terms of volume per unit time), we obtain, upon replacement, FDBT0 dXDBT = ( −rDBT) dV
dC DBT n = K appC DBT dt
dc DBT = K inη0ηceC DBT n dt
(8)
Reaction rate constant for ODS reaction (Kin) can be estimated for each reaction utilizing the Arrhenius equation as follows: K in = K 0 e−EA / RT
(9)
K0 is the frequency or pre-exponential factor and (EA) is the activation energy of the reaction. This term fits the experiment well over wide ranges of temperature and is highly proposed from different standpoints as being a very good approximation to the actual temperature dependency. The chemical reaction rate can be expressed as
(5)
3.2. Chemical Reaction Rate. Kinetic models are essential for catalyst testing at a laboratory scale and for comparing various catalysts for a given task, such as ODS. Among the various methods, parameter estimation technique is one of them where a kinetic model is assumed but its parameters are adjusted by comparing the model predictions with the
−rDBT = − 12506
dc DBT = K 0 e−EA / RT η0ηceC DBT n dt
(10)
DOI: 10.1021/acs.iecr.5b03289 Ind. Eng. Chem. Res. 2015, 54, 12503−12515
Article
Industrial & Engineering Chemistry Research If the catalytic reaction of DBT oxidation obeys nth-order kinetics, it can be integrated and one can obtain a final expression: ⎡ ⎤ k 1 ⎢ 1 1 ⎥ = app − n − 1 n − 1 n − 1 ⎢⎣ C DBTf C DBT0 ⎥⎦ LHSV
The critical specific volume of light gas oil (liquid) can be evaluated using the Riazi−Daubert equation:33 VL = 0.285(vc L)1.048
vc L = (7.5214 × 10−3(TmeABP)0.2896 (ρ15.6 )−0.7666 )MwL
(11)
(22)
3.3. Reactor Description. The TBR includes many control variables: mass-transfer coefficients, viscosity and density of the oil, diffusivities, effectiveness factor, and others. These factors are evaluated utilizing the relations presented in the literature as follows. To account for hydrodynamics and other physical impacts, an apparent kinetic constant can be stated: Kapp = Kin f (hydrodynamics) and is rewritten as (note that η0ηceKin is employed instead of Kapp): ⎡ ⎤ 1 ⎢ 1 1 ⎥ = η0ηceK in − n − 1 n − 1 n − 1 ⎢⎣ C DBTf LHSV C DBT0 ⎥⎦
Mean pore radius (rg) is given as
rg =
(12)
⎛ VP ⎜⎛ n + 1 ⎞ K in(C DBT) ⎜ ⎟ SP ⎜⎝⎝ 2 ⎠ Dei
ρp =
ρp ⎞ ⎟ ⎟ ⎠
ReL =
(14)
ReL″ =
GaL =
GaL″ = (16)
ρL uLd p μL (1 − ϵB)
(27)
d p3ρL 2 g μL 2
(28)
d p3ρL 2 g ϵB3 μL 2 (1 − ϵB)3
(29)
⎡ ⎛ dt ⎞2 ⎤ ⎜ ⎟ ⎥ − 2 ⎢ ⎝ d pe ⎠ ⎥ ϵB = 0.38 + 0.073⎢1 + ⎢ ⎛ dt ⎞ 2 ⎥ ⎜ ⎟ ⎢ ⎝ d pe ⎠ ⎥⎦ ⎣
The effective diffusivity is dependent on the Knudsen diffusivity (Dki) and molecular diffusivity (Dmi), which can be evaluated as follows:32,31,30
(30)
where ϵB is the catalyst bed porosity. The equivalent particle diameter (dpe), which is defined as the diameter of a sphere that has the same external surface (or volume) as the actual catalyst particle, is a sufficient particle characteristic, depending on the particle shape and size.
0.5
(18)
d pe = d p = 1.6 mm
(19)
This latter equation is known as the Tyn−Calus correlation. The molar volume of light gas oil (in liters) can be calculated using the following equation: vDBT = 0.285(vcDBT)1.048
(26)
The bed porosity (bed void fraction) for an undiluted catalyst bed may be evaluated utilizing the equations reported by Froment and Bischoff,28 and introduced by Jarullah et al.:35
(17)
⎛ v 0.267T ⎞ Dmi = 8.93 × 10−8⎜⎜ L 0.433 ⎟⎟ ⎝ vDBT μL ⎠
μL
A modified Galileo number is given as
Catalyst porosity (ϵS) can be estimated with the following relationship, based on experiments of total pore volume (Vg):
⎛ T ⎞ Dki = 9700rg ⎜ ⎟ ⎝ Mwi ⎠
ρL uLd p
The Galileo number (GaL) is described as
The effective diffusivity (Dei), where the structure (porosity and tortuosity) of the pore network inside the particle is taken into account through the modeling.29
ϵS = ρp Vg
(25)
A modified Reynolds number is introduced as
(15)
⎡ ⎤ ⎢ ⎥ ϵ 1 ⎥ Dei = S ⎢ ; ⎢ ⎛⎜ 1 ⎞⎟ ⎛⎜ 1 ⎞⎟ ⎥ ⎢⎣ ⎝ Dmi ⎠ + ⎝ Dki ⎠ ⎥⎦
(23)
The Reynolds number (ReL) is given as
ρB 1 − ϵB
Sg
ηce = 1.617ReL 0.146GaL−0.071
Generally, the Thiele modulus (Φ) for nth-order irreversible reaction is given as follows:24 Φ=
2Vg
The external catalyst wetting efficiency of the surface (ηce) is determined at atmosphere pressure utilizing the equation of Al-Dahhan and Dudukovic.34
(13)
n−1
20
The tortuosity factor (; ) of the pore structure, in eq 16, is given by30 ϵS 1 = 1 ; 1 − 2 log(ϵS) (24)
The effectiveness factor (η0) is to be calculated as a function of the Thiele modulus (Φ) with the following correlations employed for spherical particles:28 3(Φ coth Φ − 1) η0 = Φ2
(21)
(31)
The external volume (VP) and the catalyst surface of regular (spherical) shape (SP) can be estimated as
VP =
(20) 12507
4 π (rp)3 3
(32) DOI: 10.1021/acs.iecr.5b03289 Ind. Eng. Chem. Res. 2015, 54, 12503−12515
Article
Industrial & Engineering Chemistry Research SP = 4π (rp)2
In eq 41, Nt is the number of test runs, Cmeas DBT the evaluated product yield, and Cpred the product yield predicted by the DBT model. Using the kinetic parameters reported in the literature (without optimization), the composition of all fractions was estimated via application of model correlations in gPROMS. The comparison between experimental and predicted results using the kinetic parameters published in the literature (without the optimization) is listed in Table 2. As shown in this table, there is a big variation between the estimated and experimental values; therefore, optimization is employed on model variables for minimizing this variation. 4.1. Optimization Problem Formulation for Parameter Evaluation. The parameter evaluation problem formulation is described as follows: Given the reactor configuration, the catalyst, the feedstock, the operation conditions, Optimize the reaction orders of ODS (n1) and reaction constants (k) at various temperatures (403, 443, 473 K, respectively), To minimize the sum of squared errors (OBJ), Subject to constraints on the conversion and linear bounds on all optimization variables. Mathematically, the problem is stated as
(33)
The LGO density (ρL) as a function of process conditions is calculated by the Standing−Katz correlation:20 ρL = ρo + Δρp − ΔρT
(34)
⎛ P ⎞ ⎟ Δρp = (0.167 + 16.181 × 10−0.0425ρ0)⎜ ⎝ 1000 ⎠ ⎛ P ⎞2 ⎟ − 0.01 × (0.299 + 263 × 10−0.0603ρ0)⎜ ⎝ 1000 ⎠ (35)
The temperature employed for equation of the liquid density in this relation is given as ΔρT = [0.0133 + 152.4(ρo + ΔρP )−2.45 ](T − 520) − [8.1 × 10−6 − 0.0622 × 10−0.764(ρ0 +ΔρP )](T − 520)2 (36)
Glaso’s correlation has utilized the following generalized correlation for oil viscosity:36 μL = 3.141 × 1010(T − 460)−3.444 [log10 API]a
(37)
where
min
a = 10.313[log10(T − 460)] − 36.447
(38)
141.5 − 131.5 Sp.gr15.6
f (z , x(z), x̲ (z), u(z), v) = 0 (39)
(C L ≤ C ≤ C U ) (n j L ≤ n j ≤ n j U ) (k i j L ≤ k i j ≤ k i j U )
4. PARAMETER ESTIMATION TECHNIQUES Accurate determination of process parameters is essential to benefit from any model-based activities, such as design, control, scale-up, etc.32,37 Evaluation of reaction rate parameters can be accomplished via parameter estimation technique based on experimental data and model predictions, so that errors between experimental and predicted data are minimized. The experimental data of the ODS process considered in this work were matched against a simple power law kinetic model (eq 8): dc −rDBT = − DBT = K inη0ηceC DBT n dt Calculated yields are estimated via integration of eq 8, where CDBT0 is the feed concentration of DBT:
where f(z, x(z), x(z), u(z), (v) = 0) refers to the mathematical model of the process, z denotes the independent variable, u(z) is the optimization variable, x(z) is the set of all differential and algebraic variables, x(z) represents the derivatives of differential variables, with respect to z, and v refers the fixed parameters. C is the concentration, and CL and CU represent the lower and upper bounds of concentration. The subscripts “L” and “U” respectively denote the lower and upper bounds of the parameters concerned. The optimization solution method utilized by gPROMS is a two-step procedure that is known as the feasible path approach, which has been described in detail in Nawaf et al.38 and Jarullah et al.29,32,35 The values of constant factors employed in the models are listed in Table 3. Note that, to avoid local minima, the solutions are checked by starting with different initial guesses of the parameters. In gPROMS, one can also provide default (initial) values of the parameters with wide bounds (lower and upper bounds).
calc = C DBT f
⎛ ⎞1/(n − 1) C DBT0 n − 1 × LHSV ⎜ ⎟ ⎜C ⎟ n−1 ⎝ DBT0 η0 × ηce × K in × (n − 1) + LHSV ⎠
5. RESULTS AND DISCUSSIONS 5.1. Experimental Results. 5.1.1. Effect of Catalyst Loading on Process Conversion under Operating Conditions. The influence of initial concentration, reaction temperature, LHSV, and catalyst loading on the process conversion is investigated. At temperatures
