Subscriber access provided by University of South Dakota
Catalysis and Kinetics
Molecular-Level Kinetic Modeling of Lube Base Oil Hydroisomerization Nicholas Evenepoel, Pratyush Agarwal, and Michael T Klein Energy Fuels, Just Accepted Manuscript • DOI: 10.1021/acs.energyfuels.8b02266 • Publication Date (Web): 09 Aug 2018 Downloaded from http://pubs.acs.org on August 15, 2018
Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.
is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.
Page 1 of 15 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Energy & Fuels
Molecular‐Level Kinetic Modeling of Lube Base Oil Hydroisomerization Nicholas Evenepoel 1 , Pratyush Agarwal 1 , Michael T. Klein 1 , 2 1
Department of Chemical and Biomolecular Engineering, University of Delaware, Newark, DE 19716 2 Center for Refining & Petrochemicals, The Research Institute, King Fahd University of Petroleum & Minerals, Dhahran 31261, Saudi Arabia
Abstract A molecular‐level kinetic model was constructed for the hydroisomerization/hydrocracking of a hydrotreated deasphalted oil feed. A kinetic network was developed for the lube base oil production containing 1,105 molecular species and 15,991 reactions grouped into 9 reaction families. The molecular composition of the feedstock was reconstructed by minimizing the difference between experimental data and simulated mixture properties. More specifically, hydrocarbon types, carbon number distribution and sulfur level were matched very accurately. To model the kinetics and reduce the computational load, the Langmuir‐Hinshelwood‐Hougen‐Watson (LHHW) rate law parameters were constrained using the Linear Free‐Energy Relationship (LFER) principles. Using experimental process and product data of the hydroisomerization over a commercial catalyst, the reaction kinetics were optimized based on 141 data points. An excellent agreement was found between experimental and simulated properties of the lube base product of the hydroprocessing at three different temperatures.
1. Introduction Lubricant oils play an important role in several industrial applications. They are used to reduce friction between moving surfaces in contact, transfer heat, and prevent corrosion, among other things.1,2 Additives are often added to the lube base oil to further improve the characteristics or provide additional functionalities. Properties of prime importance for the quality of a lube base oil are the pour point, viscosity, viscosity index, and oxidative and thermal stability.3 The most common lube base oils are mineral oils derived from crude oil fractions. While only specific crude fractions from the high‐ quality Pennsylvania crudes were considered to be a feedstock for lube oil manufacture in the past, nowadays, less desirable deasphalted oil (DAO) and vacuum gas oil (VGO) fractions can also be used for the production of high quality lubricants through upgrading processes. These processes are especially necessary to increase the viscosity index and decrease the pour point of the base oils. For mineral lube base oil production, the DAO and VGO fractions are recovered from a crude oil distillation column and are deasphalted using an extractive precipitation process to remove undesirable asphalts and resins. Often, the product is then hydrotreated to improve its quality and remove impurities naturally present in the feedstock. Chevron developed a processing scheme to produce high quality lubricants from (hydrotreated) DAO and VGO by combining hydrocracking/hydroisomerization operations with catalytic dewaxing. This technology is crucial as it makes mineral lube oils the cheapest, most stable, and most readily available type of lubricants.4‐6 The most important properties of the lubricant, such as the pour point and viscosity index, are highly dependent on the molecular composition and hydrocarbon types in the base oils. Noh et al.3 developed correlations for the viscosity index based on the hydrocarbon type composition. Verdier et al.7 were able to predict the viscosity index of vacuum gas oils even more accurately based on chain length, 1
ACS Paragon Plus Environment
Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 2 of 15
branching position and aromatic content of the individual molecules. The pour point also depends on molecule type and varies even for different isomers.8 Hydrocarbon species present in VGO and DAO include normal paraffins, isoparaffins, olefins, aromatics, and naphthenes, often imbedded with heteroatoms like sulfur and nitrogen. Normal paraffins in the oil improve the viscosity index, but their very high pour point can be problematic. Isoparaffins are the most desired molecules since they produce a low pour point and high viscosity index. Multiring‐naphthenes and aromatics degrade the quality of the lube oil due to their low viscosity index. Severe hydrocracking processes can be used to convert these ringed molecules into more favorable isoparaffins and fewer‐ring naphthenes. Hydroisomerization will convert normal paraffins into isoparaffins and improve the viscosity index. A hydropretreatment will remove most of the olefins and heteroatoms which improves the oxidative and thermal stability of the product.9 Molecular‐level kinetic modeling makes it possible to simulate the molecules and associated reactions in lube oil production to predict the product composition based on (bulk) feed properties, reaction chemistry, and reaction conditions. Tracking the individual species that impact the product lube oil properties can allow for a more accurate prediction of those properties. This makes molecular‐level kinetic modeling an excellent tool to control and optimize hydrotreatment and hydroisomerization processes for lube oil production. In this study, a molecular‐level kinetic model for the hydroisomerization of a hydrotreated vacuum gas oil is developed. The hydroprocessing kinetics are optimized based on experimental data published by Wang et al.10, who report the hydroprocessing over a commercial catalyst with strong acidity at three different reaction temperatures. First, a comprehensive reaction network is generated based on experimental hydroprocessing studies found in literature. Second, the molecular composition of the hydrotreated gas oil feed is calculated by matching experimental data to simulated data. Finally, the kinetics of the different reaction families are extracted and evaluated based on the experimental product composition. An in‐house software, the Kinetic Modeler’s Toolbox (KMT), was used to build and solve the kinetic model.11 While some kinetic models have been developed for the hydrocracking of vacuum gas oil feeds12‐13, this manuscript is unique by modeling the hydroisomerization process of a hydrotreated feed designed for the lube base oil manufacture.
2. Reaction network generation The hydrotreated vacuum gas oil feedstock of the lube oil hydrosiomerization process consists of a complex mixture of organic compounds which undergo several hydrogenation reactions. The overall hydroprocessing chemistry and reaction pathways of this process were developed based on experimental studies in literature.3,4,14 Wang et al.10 analyzed the feed and lube oil products using GC‐ MS. They reported the carbon number distribution of hydrocarbons with carbon numbers up to 46. In the model, organic compounds with a carbon number up to 50 were included for each species type to capture this detail. Table 1 provides an overview of the different organic species in the model. The most prominent hydrocarbon types are the normal paraffins, isoparaffins, naphthenes and aromatics. As it was impossible to include all isomers of each molecule type, a reasonable minimization was made to at least match the level of detail reported by the experimental studies. The isoparaffins with a carbon number smaller than 20 were limited to only one methyl sidechain attached to any carbon atom of the chain. For paraffins with more than 20 carbon atoms in the chain, three monobranched isomers were included where the methyl sidechain is situated at the 2‐position, the Ncarbon/4 position, or the middle of the carbon chain. This resulted in a representative branched isomer distribution which could be further adjusted in future work. Naphthenes and aromatics included hydrocarbons with up to four 5‐ 2
ACS Paragon Plus Environment
Page 3 of 15 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Energy & Fuels
or 6‐membered rings. A vacuum gas oil feed typically contains, besides pure hydrocarbons, sulfur‐ and nitrogen‐containing compounds. Since the feed was hydrotreated prior to hydroisomerization, most of these heteroatom‐containing compounds were presumed to already have been removed from the oil. Wang et al.10 reported only 10.5 and 1 wppm sulfur and nitrogen fractions, respectively, in the feed. The sulfur‐containing compounds were represented by thiophenes, and their derivatives as the sulfur‐carbon bonds of sulfides and mercaptans are more likely to be cracked during the pretreatment.14 Nitrogen‐containing compounds were neglected in this kinetic model due to their very low level in the feed and lack of nitrogen level data of the product. Also, olefins and olefin hydrogenation reactions were excluded in this model as olefins were not present in the feed. Table 1. Species in lube base oil hydroprocessing.
Species type
Sample Species
n‐paraffin
i‐paraffin
Naphthene
Aromatic
Thiophene derivative The included hydrocarbons can undergo various reactions. While the reaction conditions and acidic bifunctional catalyst favored hydroisomerization and hydrocracking reactions, hydrotreatment reactions also happened at lower reaction rates. The reactions can be classified into nine reaction families as shown in Table 2. In the sample reactions in Table 2, R stands for the sidechain which can have various lengths. The reactions take place on the acid and/or metal site of the bifunctional catalyst. Mainly the destructive hydrocracking and hydroisomerization reactions happen on the acid sites. Hydrotreatment reaction types are favored by mild reaction conditions and depend on the metallic active sites. Desulfurization reactions remove the undesired heteroatoms from the lube oil. These reactions cleave either sulfur‐sulfur or carbon‐sulfur bonds. In the case of benzothiophenes, heterocyclic saturation 3
ACS Paragon Plus Environment
Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 4 of 15
occurred before the desulfurization since styrene type molecules are not observed in hydroprocessing.14 Ring saturations were modeled to convert aromatic rings into saturated rings with the same number of carbons on the metal sites. As noted earlier, the number of isomers had to be limited to reduce the computational load of the model without losing accuracy. Ring isomerization transforms cyclohexanes into cyclopentanes reversibly but was limited to mono‐ and di‐ring naphthenes with carbon numbers below 20 in the model. Hydrodealkylation and sidechain cracking reactions were allowed to occur on any naphthenic, aromatic and sulfur‐containing organic compound with a sidechain. Paraffin isomerization reactions were used to convert normal paraffins into isoparaffins which can further undergo cracking to smaller paraffins. The sidechain cracking, paraffin cracking and hydrodealkylation reactions are crucial for hydrocracking processes as they break large molecules into smaller molecules. Table 2. Reaction families for lube base oil hydroprocessing.
Reaction Family
Sample Reactions
Catalyst site
1. Desulfurization
Metal
2. Heterocyclic Saturation
Metal
3. Ring Saturation
Metal
4. Ring Opening
Acid
5. Ring Isomerization
Acid
6. Hydrodealkylation
Acid
7. Sidechain Cracking
Acid
4
ACS Paragon Plus Environment
Page 5 of 15 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Energy & Fuels
8. Paraffin Isomerization
Acid
9. Paraffin Cracking
Acid
With the species and reactions defined, the reaction network was constructed using the INGen (Interactive Network Generator) tool. In INGen, every molecular species has a computational representation in the form of a bond electron matrix and SMILES code. By storing each reaction family as a small reactive matrix subgroup, the library of desired reaction family matrices can easily be applied to all of the molecule bond electron matrices as simple matrix addition operations. The SMILES codes are used to represent unique species in the model and find their properties. Figure 1 shows a simplified computational representation of a paraffin cracking reaction. To minimize the network size while obtaining the desired molecular detail, seeding was used as a model reduction tool as described by Joshi et al.15 Table 3 gives on overview of the obtained network statistics. The final model contained a total of 1,105 species and 15,991 reactions. Table 3. Network statistics generated using the Interactive Network Generator (INGen).
Species Type n‐Paraffin i‐Paraffin Naphthene Monoring Diring Triring Tetraring Aromatic Monoring Diring Triring Tetraring Thiophene derivative H2, H2S, NH3 Total Species
Number 50 205 215 45 93 44 33 507 249 155 70 33 125 3 1 105
Reaction Type Desulfurization Heterocyclic Saturation Ring Saturation Ring Opening Ring Isomerization Hydrodealkylation Sidechain Cracking Paraffin Isomerization Paraffin Cracking
Count 82 43 540 119 31 939 13 750 190 297
Total Reactions
15 991
5
ACS Paragon Plus Environment
Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 6 of 15
Figure 1. Simplified matrix representation of paraffin cracking reaction.
3. Model equations and kinetics The hydroprocessing of hydrocarbon feedstocks is a heterogeneous catalytic process.16 Based on the reaction network, a material balance for every species was generated. Along with the initial conditions of the feed and the reactor, this defined the initial value problem used to solve the kinetic model. The in‐house software KME (Kinetic Modeling Editor) was used to generate this system of material balances.11 The same model equations and kinetic theories were implemented as for an earlier study in the Klein research group on green diesel hydroprocessing.17 In this study, the hydroisomerization process in packed bed reactors with a bifunctional metal/acid catalyst under plug flow conditions was modeled assuming the Langmuir‐Hinshelwood‐Hougen‐Watson (LHHW) kinetic rate law formalism.18,19 For the hydroisomerization on a commercial hydroprocessing catalyst, surface reaction control was assumed to be the limiting reaction step. An explicit hydrogen partial pressure dependence was added to the adsorption denominator. Equation 1 gives the modified LHHW kinetic rate law for reversible reactions used by the KME tool. ∗∏
∏
∏
∗
,
∑
1
1
,
is the surface reaction rate coefficient, is the In this equation, r is the reaction rate, concentration of catalyst sites, is the overall reaction equilibrium constant, is the hydrogen is the adsorption equilibrium pressure dependency factor of the acid or metal catalyst sites and and overall constant, also depending on the catalyst site. The adsorption equilibrium constant equilibrium constant were calculated based on theoretical and empirical correlations. For the overall equilibrium constant, the standard thermodynamic formulation is given by equation 2. Group contribution methods and corresponding state functions were used to calculate the Gibbs free energy and other thermodynamic properties of every molecular species. More specifically, these properties were calculated using Marrero20 and Benson21 group contribution methods. ∆
ln
2
Korre and Klein22 developed a quantitative structure/reactivity relationship (QSRR) to calculate the adsorption equilibrium constant of each reaction. Equation 3 gives a modified version of this relationship used in this study. The adsorption constant of a molecule was defined as a function of its structure and the type of site (acid/metal) the species is adsorbed on. The aromatic ring number, naphthenic ring number, and non‐ring saturated carbon number dependencies were included in the QSRR. The parameters , , , , , , and , had different values for the acid and metal sites. ln
,
,
,
,
,
3
for each reaction i was calculated using the Arrhenius The surface reaction rate coefficient equation as shown by equation 4. The two kinetic parameters in this equation are the Arrhenius constant or pre‐exponential factor A and activation energy E. As in most hydrotreatment processes, tens of thousands of reactions can occur, so the Bell‐Evans‐Polanyi linear free‐energy relationship (LFER) was implemented as an additional method to reduce the number of tunable parameters in the 6
ACS Paragon Plus Environment
Page 7 of 15 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Energy & Fuels
model. 23,24 This theory exploits the experimentally observed linear relationship between reaction rates and certain reaction indices like reaction enthalpy of member reactions in a reaction family. Each reaction family j consists of a homologous set of reactants subject to the same type of reactions. The activation energy for each reaction i was defined by equation 5 with tunable parameters and , . Therefore, the surface reaction rate coefficient calculation for every reaction family now only required three tunable parameters in total, rather than two per individual reaction in a reaction family. As the lube oil hydroisomerization process consisted of 9 reaction families but 15,991 reactions, this was a significant parameter reduction. ln
, ,
,
ln ,
∆
4 5
4. Initial Conditions: Feed composition generation Wang et al.10 experimentally measured feed properties including density, simulated distillation data, hydrocarbon type analysis, carbon number distribution, and sulfur level. To solve the kinetic model, these properties had to be converted to a molecular representation and associated amounts (e.g., wt%). The in‐house software ICG (Initial Condition Generator) was used to find the weight percentage of each individual species in the feed. To minimize the number of tunable parameters, the molecule set was represented as an attribute probability density function (PDF) tree. The PDF tree represents the numerical probability of each structural moiety, like naphthenic rings, aromatic rings, and heteroatoms, present in the molecular set, arranged as a network of parent‐child generations connected by branches. Each member of a generation imposes its probabilities on all molecules that satisfy its conditions. Child generations impose further classifications on members that satisfy their respective conditions. A given molecule can have attributes that satisfy the conditions of multiple members in the same generation of the PDF tree since each independent PDF is normalized to one. This allowed for the merging of child generations that were assumed to behave similarly; most notably, the carbon number distributions of similar molecule functional cores. Figure 2 shows the PDF tree generated for this study. Carbon number distributions of all species satisfying some user‐defined conditions were defined by average and standard deviation parameters assuming histograms with a Gaussian distribution. The overall feed was distributed into n‐paraffin, isoparaffin, aromatic, and naphthenic child members. These members were further classified based on carbon number and ring number. As mentioned in section 2, the only sulfur‐containing compounds considered in the feed were thiophenes and their derivatives. Therefore, aromatics have a sibling classification of ring number and contained sulfur where the underlying carbon number distributions are assumed to be identical. The software juxtaposes the attribute PDF tree and histogram parameters into a simulated molecular composition. A tuning algorithm iteratively optimized the PDF tree values, and therefore molecular composition, by matching the experimental and simulated properties. Bulk properties were calculated from the simulated composition using mixing rules for mixture properties and group contribution methods for individual species as described in section 3.
7
ACS Paragon Plus Environment
Energy & Fuels 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 8 of 15
Figure 2. Probability density function tree for feed reconstruction of lube base oil hydroprocessing.
Figure 3 compares the experimentally obtained carbon number distribution and hydrocarbon type analysis with the simulated values from the reconstructed molecular composition. Table 5 gives the comparison of some other bulk properties. Overall, the composition data, including hydrocarbon type and carbon number distribution, matched the generated molecular set very accurately. The lower density prediction can be attributed to the simple mixing model used for density and the impact of branching of the aromatic and naphthenic species. Although aromatics and naphthenics can have multiple branches rather than the one long chain representation in this model, due to the large number of possible isomers, they were ignored. These branched molecules have higher densities than the single chain but are not included. The simulated true boiling points differed slightly from the experimentally obtained simulated distillation data. A first explanation could be experimental errors and the reliability of the simulated distillation method used in the experimental study. Simulated distillation experiments give only an approximation of the true boiling point ranges for a mixture. This is especially true in this case as a vacuum distillation method, ASTM D 1160, was used to avoid preliminary degrading of long carbon chains before boiling at atmospheric conditions.25 A second reason could be the incapability of the group contribution methods to accurately simulate the boiling points of heavier molecular species. Manually, boiling point data from literature were used to improve the boiling point estimation for chrysene and naphthalene compounds. With these reasons in mind and the good agreement of all other detailed feed properties, it can be concluded that the simulated composition approximates the real molecular composition reasonably well.
8
ACS Paragon Plus Environment
Page 9 of 15
35 Feed composition (wt %)
30
Experimental
25 Simulated
20 15 10 5 0
n‐Paraffins i‐Paraffins Monoring Naph. Diring Naph. Triring Naph. Tetraring Naph. Aromatics C20 C20 C21 C21 C22 C22 C23 C23 C24 C24 C25 C25 C26 C26 C27 C27 C28 C28 C C29 29 C30 C30 C31 C31 C32 C32 C33 C33 C34 C34 C35 C35 C36 C36 C37 C37 C38 C38 C39 C39 C C40 40 C41 C41 C42 C42 C43 C43 C44 C44 C45 C45
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Energy & Fuels
Figure 3. Evaluation of feed reconstruction for lube base oil hydroprocessing: carbon type and number distribution. Table 4. Evaluation of feed reconstruction for lube base oil hydroprocessing: Density, sulfur level and simulated distillation data.
Properties Density at 20 °C ∙ Sulfur level (wppm) ASTM D 1160 distillation (°C) 10 % 30 % 50 % 70 % 90 %
Experimental Simulated 0.8632 0.8405 10.5 10.2 353 493 535 538 547
352 455 492 514 542
5. Kinetic model evaluation Based on the experimental data reported by Wang et al.10, the parameters of the molecular‐level kinetic model were optimized. An objective function of the form given by equation 6 was minimized to reduce the difference in value between each observed and predicted experimental property. is a user defined weight function estimated by the average experimental error in the property method. In cases without sufficient information, an approximation was made based on the order of magnitude of the measurement. 6
Three experiments were done at reaction temperatures of 360, 370 and 380 °C. The space velocity, pressure, hydrogen/oil volumetric ratio and catalyst loading were kept constant at 0.8 h‐1, 9.0 MPa, 800 v/v and 100 mL, respectively. The catalyst was a commercial bifunctional catalyst containing platinum and phosphate on silica‐alumina with a strong acidity.26 For each reaction temperature, the study reported the hydrocarbon type analysis and carbon number distribution of the lube base oil 9
ACS Paragon Plus Environment
Energy & Fuels
product. For the experiment at 380 °C, more detailed product properties were reported including sulfur fractions. The hydroprocessing was followed by a vacuum distillation to refine the product in different distillation cuts based on the boiling point distribution. Only the high overall boiling point cuts (b.p. ≥ 350 °C) are useful as lube base oil and were crucial in evaluating the model optimization. For the low overall boiling point cuts (b.p.
