Kinetic Modeling of Coke Formation and Deactivation in the Catalytic

Tapan Bera , Joris W. Thybaut , and Guy B. Marin. Industrial ... Gladys Jiménez-García , Roberto Quintana-Solórzano , and Rafael Maya-Yescas. Industri...
2 downloads 0 Views 236KB Size
14

Ind. Eng. Chem. Res. 2003, 42, 14-25

Kinetic Modeling of Coke Formation and Deactivation in the Catalytic Cracking of Vacuum Gas Oil Tarek M. Moustafa† and Gilbert F. Froment* Department of Chemical Engineering, Texas A&M University, College Station, Texas 77843-3122

The kinetics of coke formation in the catalytic cracking of vacuum gas oil are expressed in terms of elementary steps of carbocation chemistry, along the lines used for the main reactions of this process. The reaction network is generated by a computer, without any lumping of components and elementary steps. The rate coefficients depend on the type of elementary step and the nature of the reacting and produced carbocations or molecules, but the single-event concept accounts for the other structural aspects. The application of this concept, a number of plausible assumptions, and thermodynamic constraints tremendously reduce the number of rate parameters to be determined from experimental data. The effect of coke deposition on the activity of the catalyst is accounted for by means of deactivation functions written in terms of the coke content of the catalyst, not contact time. The kinetic model is applied in the simulation of a commercial riser reactor for the catalytic cracking of three types of vacuum gas oil. The results are presented in terms of the usual commercial fractions but also in terms of more detailed fractions, which provide further insight in the process and the catalyst deactivation. 1. Introduction Modern catalytic cracking of vacuum gas oil (VGO), a major oil-refining process, is carried out in riser reactors, a technology that permits continuous transfer of the catalyst to the regenerator, where the coke deposited by the cracking is removed by controlled combustion. Feedstocks such as VGO comprise a large number of components, each of which generates a complex reaction network. By way of example, the catalytic cracking network of n-dodecane contains 201 n- and i-paraffins and 1123 olefins and that of n-pentylcyclohexane, a C11 mononaphthene, 2697 naphthenic and acyclic components. The number of species in such networks significantly increases with the carbon number. Because of this complexity (and also because of very incomplete feedstock analysis), the modeling of processes dealing with heavy oil fractions has been traditionally based upon the characterization of feedstock and products in terms of “lumps”, fractions defined partly by their overall chemical nature and partly by their boiling range.1 From the late 1980s onward, schemes reflecting more the chemistry of the process and based upon individual molecules were developed. Liguras and Allen2 described VGO conversion in terms of a relatively large number of pseudocomponents, most of which are lumps in their own way, however. Klein et al.3,4 generated these pseudocomponents from analytical characteristics using Monte Carlo simulation. Instead, Quann and Jaffe5,6 expressed the chemical transformations in terms of typical structures of the molecules, without completely eliminating lumps. Retaining lumps in the reaction scheme and associated kinetic model inevitably leads to rate parameters that depend on the feedstock composition, so that additional experimental work is needed for each new feedstock. * Corresponding author. E-mail: [email protected]. Tel: (979) 845-0406. Fax: (979) 845-6446. † On leave of absence from Chemical Engineering Department, Cairo University, Cairo, Egypt.

The present investigation follows the lines developed by Froment and co-workers.7-16 The reaction scheme is developed in terms of elementary steps of carbenium and carbonium ion chemistry, like (de)protonation, hydride and methyl shift, PCP-branching isomerization, hydride transfer, β-scission, protolytic scission, (de)alkylation, and cyclization shown in Figure 1. The network becomes gigantic, but it is completely generated by a computer.11 The molecules and ions are characterized by vectors, while the elementary steps of their transformation are performed by means of Boolean relation matrices. When the chemical analysis is incomplete and comprises lumps, all possible components belonging to the definition of these lumps are considered in the network generation. Because the number of types of elementary steps is much smaller than the number of molecules and the number of reactions between these, the approach limits the number of rate parameters. It is clear, however, that the rate coefficient of a given type of elementary step, e.g., cracking by β-scission, also depends not only on the structural aspects of the involved species, like the nature and number of C atoms on the carbon atom carrying the positive charge (primary, secondary, and tertiary ions), but also on its vicinity, while the chain length may also play a role. Clearly, considering all of these aspects still leads to an unrealistic number of rate parameters, unless they are accounted for through a fundamental modeling of the parameters. The effects of the structural aspects on the rate coefficient of an elementary step were expressed in terms of transition-state chemistry, statistical thermodynamics, and introduction of the single-event concept. The modeling of the rate coefficient and the introduction of thermodynamic constraints led to a tractable number of parameters. Until now the coke formation in catalytic cracking was modeled in an empirical way by relating it to the contact time between the reacting gases and the catalyst.1 Instead, in the present paper the coke formation is linked to the steps of the “main ” process, thus becoming an extension of the approach outlined above. After

10.1021/ie0204538 CCC: $25.00 © 2003 American Chemical Society Published on Web 11/28/2002

Ind. Eng. Chem. Res., Vol. 42, No. 1, 2003 15

Figure 1. Elementary steps of carbocation chemistry.

Figure 2. Elementary steps of coke formation in catalytic cracking.

“coke” is defined, the reaction network that generates it is developed in terms of the elementary steps of carbocation chemistry in the same way as it was done for the “main” process. This is the only way to reduce the number of rate parameters appearing in the modeling of coke formation, the more so that many steps are of the same type as those encountered in the “main” process. The effect of the coke content of the catalyst on its deactivation is accounted for through a number of deactivation functions related to the nature of the elementary steps. The parameters of the model are determined from published experimental data, and the kinetic equations are applied in the simulation of a commercial riser reactor for the catalytic cracking of VGOs. 2. Coke Formation in Terms of Elementary Steps Coke is defined as a carbenium ion with a size and structure preventing its desorption, thus permanently covering the active site(s) it was formed upon. The literature agrees that coke mainly consists of polyaromatic moieties.17 In the present work on the catalytic cracking of VGO, coke was considered to consist of five or more aromatic rings and of a total number of C atoms exceeding 40. The types of elementary steps of carbenium/carbonium ion chemistry shown in Figure 1 that play an important role in coke formation are mainly alkylation, cyclization, (de)protonation, and hydride transfer, more specific examples of which are given in Figure 2. Cyclization does not increase the size of the carbenium ion but leads to ring formation, a prerequisite for

Figure 3. Reaction scheme for coke formation in the catalytic cracking of VGO.

coke formation. Transformation into aromatic rings is possible by a succession of deprotonations and hydridetransfer steps. A simplified representation of the coking scheme is given in Figure 3. The molecules considered in the model developed for the catalytic cracking of VGO are nparaffins (individually up to C40) and further per C number i-paraffins; n-olefins; i-olefins; mono-, di-, tri-, and tetra-ring naphthenes, aromatics, and aromatic olefins; and mono-, di-, and tri-ring naphthenoaromatics. It is clear from the preceding that coke, or the intermediate given that name, is not a dead species and may interact with other species of the reaction mixture. This is in agreement with more global experimental observations. 3. Modeling the Rate Coefficient of Elementary Steps The rate coefficient of an elementary step of carbenium ion chemistry largely depends on the nature of the involved carbenium ions. When methylcarbenium and primary carbenium ions, which are much less stable than secondary and tertiary ions, are disregarded, this would reduce the number of rate coefficients of an elementary step like PCP-branching isomerization, e.g., to only four [kPCP(s;s), kPCP(s;t), kPCP(t;s), and kPCP(t;t)], were it not that other effects of structure also had to be accounted for in the modeling of the rate coefficient, thus greatly enhancing their number.

16

Ind. Eng. Chem. Res., Vol. 42, No. 1, 2003

To reduce the number of independent rate parameters, Froment and co-workers7-16 factored out the structure effect from the change of standard entropy associated with the transformation of a reactant into a product through an activated complex. From transitionstate theory, the rate coefficient can be written as

k′ )

( ) (

kBT ∆S° exp h R

q

exp -

∆H° RT

)

q

(1)

According to statistical thermodynamics, the standard entropy of a component consists of translational, vibrational, and rotational contributions. The latter, considered as representative for the structure, each consist of an intrinsic term S ˆ ° and a term related to the symmetry number of the species, σ:

ˆ rot° - R ln σ Srot° ) S

(2)

Accounting, in addition, for chirality leads to a global symmetry number, σgl, defined by

σgl ) σ/2n

(3)

in which n is the number of chiral centers in the species. The change in the standard entropy due to symmetry changes associated with the formation of the activated complex is then given by

(4)

Expliciting ∆S°q in eq 1 by means of eq 4 leads to

( ) ( ) (

)

kBT σgl,r ∆S ˆ° ∆H°q exp exp h σgl,* R RT

(5)

or

k′ ) nek

(6)

The rate coefficient of an elementary step, k′, is now written as a multiple of that of a “single event”, k. The number of single events, ne, is the ratio of the global symmetry numbers of the reactant and the activated complex:

ne ) σgl,r/σgl,*

(7)

The single-event concept reduces the number of frequency factors, and by extension the number of rate coefficients, of the above elementary steps to four, whatever the structure of the ions involved in this step. The effect of the latter is accounted for by ne. The calculation of the symmetry numbers requires the knowledge of the structure of the species. Modern quantum chemical packages are of great help in this. The temperature dependence of the single-event rate coefficient is written in terms of the Arrhenius equation:

k)A ˜ exp(-Ea/RT)

4. Thermodynamic Constraints on the Parameters A further reduction of the number of rate parameters is necessary in the steps involving olefins because the number of olefins and cyclic olefins rapidly increases with carbon number. Therefore, it was assumed that the rate coefficient for the protonation of olefins/cyclic olefins does not depend on the olefin. This is justified by the results of quantum chemical calculations showing that the structure of the activated complex is close to that of the reacting olefin/cyclic olefin. Further reduction of the number of rate coefficients is based on thermodynamic constraints. By way of example, two such constraints and their consequences are discussed here.7 4.1. Constraints on the Rate Coefficients for Deprotonation. The equilibrium constant for the isomerization of olefin O1 into olefin O2 proceeding over a secondary carbenium ion can be written as

Kisom(O1SO2) )

kPr(O1;s) kDep(s;O2) kDep(s;O1) kPr(O2;s)

(8)

A ˜ is a modified frequency factor because the structure

(9)

Introducing the above-mentioned assumption concerning the protonation rate coefficients leads to

kDep(m;Oj) ) kDep(m;Or) Kisom(OrSOj)

σgl,r ∆Ssymm° ) R ln σgl,*

k′ )

effect on ∆S° has been accounted for by factoring out the number of single events, ne.

(10)

where m represents a secondary or tertiary carbenium ion and Or a reference olefin with the double bond in such a position that secondary as well as tertiary ions can be formed by protonation. From the relationship (10), it follows that within a given carbon number only two rate coefficients are independent: kDep(s;Or) and kDep(t;Or). 4.2. Constraints on the Rate Coefficients for Isomerization. From the analysis of equilibria in a reaction sequence consisting of olefin protonation, carbenium ion isomerization, and deprotonation, Baltanas et al.7 derived the relationship

kisom(t;s) kisom(s:t)

)

kPr(s) kDep(t;Or) kPr(t) kDep(s;Or)

(11)

which reduces the number of independent rate coefficients for any type of isomerization from four to three. In addition, it follows from eq 11 that the ratio kDep(t;Or)/kDep(s;Or) has to be identical for all reference olefins, regardless of their carbon number because the other two ratios in eq 11 are independent of the olefin. Consequently, for a given carbon number, there are only two independent single-event rate coefficients for the deprotonation and methyl shift elementary steps. 5. Rate Equations for Cyclization and Alkylation at the Single Event Level The rate of coke formation is dealt with in the same way as that of the main reactions, i.e., in terms of single events. Consider cyclization of the paraffinic side chain of a specific di-ring aromatic molecule, diAro, through an elementary step of a secondary carbenium ion, R+(s), formed, e.g., by hydride transfer between diAro and paraffinic carbenium ions and disappearing by depro-

Ind. Eng. Chem. Res., Vol. 42, No. 1, 2003 17

tonation and hydride transfer with paraffinic ions. The concentration of R+(s) follows from a pseudo-steadystate balance between the rates of formation and disappearance

CR+(s) )

(ne)htkht(s) CB+PdiAroi

(12)

npar

(ne)avg depkDep(s) + kht

(ne)htPP ∑ j)1

j

∑CR j

(13)

+

rcycl ) (ne)cyclkcycl(s) CR+(s)

(14)

where kcycl(s) is the single-event rate coefficient for the cyclization of R+(s) and (ne)cycl the number of single events. Eliminating the inaccessible carbenium ion concentration in eq 14 in favor of gas-phase partial pressures by means of eq 12 leads to

(ne)cyclkcycl(s) (ne)htkht(s) CB+PdiAroi npar

(ne)avg depkDep(s) + kht

(ne)htPP ∑ j)1

(15)

j

If the source of the carbenium ion were the protonation of a diaromatic molecule having an olefinic side chain, the rate of cyclization would be written as

[

rcycl ) (ne)cyclkcycl(s) kPr(s) (1 - CB+)

(ne)PrK(OjSOr) PAr ∑ j)1 npar

(ne)avg depkDep(s) + kht

(ne)htPP ∑ j)1

j

]

i

(16)

where K(OjSOr) is the equilibrium constant for the isomerization between an olefin Oj and a reference Or having the same carbon number. An alkylation elementary step leading to coke requires a large poly-ring aromatic molecule and a secondary or tertiary carbenium ion, R+(m). The latter can be generated, e.g., by hydride transfer, and disappears by deprotonation and hydride transfer. A pseudo-steadystate balance on the corresponding rates yields the concentration CR+(m) that enters into the rate equation of the alkylation step

ralk ) (ne)alkkalk(m) PArCR+(m)

(ne)avg depkDep + kht

(ne)htPPj ∑ j)1

When the source of the carbenium ion is the protonation of a polyaromatic molecule with an olefinic side chain, the rate of alkylation can be written as

[

∑(n )

(ne)alkkalk(m) PArkPr(m) (1 - CB+)

e PrK(OjSOr)

j)1

]

PAri

npar

∑(n )

(ne)avg depkDep + kht

j)1

e htPPj

(19)

6. Rate Equations for Coke Formation in the Catalytic Cracking of VGO

The cyclization rate equation is written as

rcycl )

(18)

npar

ralk )

where (ne)avg dep is the average of ne for the deprotonation steps of R+(s), PPj the partial pressure of a paraffin or a saturated side chain of a naphthene or an aromatic component, Pj, and CB+ the concentration of Brønsted acid sites covered with carbenium ions, relative to the total concentration of Brønsted sites:

C B+ )

ralk )

(ne)alkkalk(m) PAr[(ne)htkht(m) CB+PkArj]

In the present approach, the generation of the VGO catalytic cracking network does not consider any lump. Even with incomplete feedstock characterization, it starts from all possible individual components, present as such in the feed or possibly contained in the analytical lumps, and considers all of their possible elementary steps. Consequently, the number of intermediates and products in the catalytic cracking of a VGO is extremely large (5 160 748, up to C40). In the simulation stage, however, lumping has to be introduced, to limit the number of continuity equations to be integrated and to come to a tractable output. Choosing molecular lumps on the basis of the analysis permits their assignment of initial partial pressures. Using GC-MS analysis, e.g., HVGO was decomposed into some 130 components and lumps.12 The complete network is thus reduced to a network containing the individual n-paraffins, the iparaffins, n-olefins, and i-olefins per C number, mono-, di-, tri-, and tetranaphthenes, mono-, di-, tri-, and tetraaromatics, and mono-, di-, and trinaphthenoaromatics, all per C number up to C40. The total number of these components and lumps amounts to 670. They are linked by a network of 45 000 “global” reactions, the rates of which are constructed from those of the involved elementary steps. By way of example, the rate equation for the cyclization of side chains of a lump consisting of tri-ring aromatic molecules, triAroL, leading to a naphtheno triring aromatic lump can be written as

rcycl )

∑s (ne)cyclkcycl(s) CR

+(s)

1

+

∑t (ne)cyclkcycl(t) CR

+(t)

1

(20) or after elimination of CR+(m) by means of a pseudosteady-state balance on rates of formation and disappearance:

rcycl ) (LC)cycl(s) kcycl(s) Fcycl(s) PtriAroL + (LC)cycl(t) kcycl(t) Fcycl(t) PtriAroL (21)

(17)

LC(m), the so-called lumping coefficient, is given by where PAr represents the partial pressure of the polyaromatic molecule. Expliciting CR+(m) leads to

LC(m) )

(ne)avg ht(ne)cyclyL(triAro)L ∑ m

(22)

18

Ind. Eng. Chem. Res., Vol. 42, No. 1, 2003

The summation includes all of the cyclization elementary steps converting tri-ring aromatic carbenium ions of type m (either secondary or tertiary), generated out of the components of the (triAro)L lump. (ne)avg ht is the average number of single events of hydride transfer transforming the tri-ring aromatic component into the reactant carbenium ion, (ne)cycl is the number of single events of the cyclization elementary step, and yL(triAro)L is the mole fraction of a component of the lump (triAro)L. If the composition of this lump does not correspond to equilibrium, the mole fraction can be approximated by the reciprocal of the number of (triAro)L components. The lumping coefficient depends on the topology of the network of elementary steps and the choice of the lumps but does not contain the single-event rate coefficients. Therefore, the lumping coefficients for the various global reactions between lumps do not have to be calculated online, only once and for all for the generated reaction network. When the carbenium ion involved in the cyclization is formed by hydride transfer from a tri-ring aromatic molecule and disappears by deprotonation and hydride transfer with paraffins, Fcycl, is given by

Fcycl )

[

khtCB+ npar

(ne)avg depkDep + kht

(ne)htPP ∑ j)1

]

(23)

j

CB+ is the sum of all of the carbenium ion concentrations and PPj is the partial pressure of the jth molecular components of the tri-ring aromatic lump with which the tri-ring aromatic carbenium ions react by hydride transfer. An example of coke formation by alkylation in a VGO is the transformation of a lump consisting of tri-ring aromatic molecules with various side chains, triAroL, by means of secondary and tertiary di-ring aromatic carbenium ions generated out of a di-ring aromatic lump, diAroL. This leads to a penta-ring aromatic carbenium ion lump (viz. Figure 2) whose components will not desorb and that will be categorized as coke. The rate of alkylation can be written as

ralk )

(ne)alkkalk(s:s)CR (s) PtriAro + ∑ s;s (ne)alkkalk(t;t)CR (t) PtriAro ∑ t;t

(24)

ralk ) LCalk(s;s)Falk(s;s) kalk(s;s) PdiAroPtriAro + LCalk(t;t) Falk(t;t) kalk(t;t) PdiAroPtriAro (25) in which

∑ (ne)ht(ne)alkyL(diAroL)

molecule, by means of hydride transfer, and disappears by deprotonation and hydride transfer with paraffins, Falk equals

Falk )

[

khtCB+ npar

(ne)avg depkDep + kht

(26)

m;m

with yL(diAroL) the mole fraction of a specific component within the diaromatic lump. When the carbenium ion of the alkylation step is formed out of a di-ring aromatic

(ne)htPP ∑ j)1

]

(27)

j

and when the carbenium ion is generated from protonation in the olefinic side chain of a di-ring aromatic molecule

[

K(OjSOr) ∑ j)1

kPr(1 - CB+)

(ne)avg depkDep + kht

or, after elimination of the inaccessible carbenium ion concentration, CR+, by means of a pseudo-steady-state balance on its rate of formation and disappearance:

LCalk(m;m) )

Figure 5. Lumping coefficient vs carbon number for the (s;s) alkylation of a tetra-ring aromatic lump with a C20 tetra-ring aromatic.

Falk )

+

+

Figure 4. Lumping coefficient vs carbon number for the (s;s) alkylation of a tri-ring aromatic lump with a C20 tetra-ring aromatic carbenium ion.

npar

(ne)htPP ∑ j)1

]

(28)

j

Figure 4 shows the evolution of the lumping coefficients with the carbon number for the alkylation of a tri-ring aromatic species with a C20 tetra-ring aromatic carbenium ion. The lumping coefficient decreases with the carbon number because the fraction of the molecules that can undergo this specific alkylation rapidly decreases with the carbon number. The same behavior is observed for tetra-ring aromatic lumps, as shown in Figure 5. Alkylation leads to a number of lumps with a different carbon number as well as a different structure, depending upon the reacting species, and thus substantially increases the number of components of the reaction network. Because those corresponding to “coke” do not desorb and do not participate in further elementary steps, their identity does not matter in the set of continuity equations used in the reactor simulation, so that one lump suffices to represent coke.

Ind. Eng. Chem. Res., Vol. 42, No. 1, 2003 19

Figure 6. Distribution of the mono-ring aromatics in VGO-A. Table 1. Deactivation Functions for the Various Types of Elementary Steps exp(-R1Cc)

exp(-R2Cc)

exp(-R3Cc)

protonation protolytic scission hydride transfer

alkylation

deprotonation methyl shift PCP branching cyclization β-scission

Because the various elementary steps are deactivated to a different extent as a result of coke formation, it is not sufficient to have one deactivation function in the kinetic model. Instead, three deactivation functions were introduced, corresponding each to a number of elementary steps, and these are given in Table 3. The first deactivation function is assigned to elementary steps which can be considered as bimolecular and involving a molecule and a carbenium ion, as in hydride transfer, or a molecule and a proton of the surface, as in hydride abstraction, followed by protolytic scission of the carbonium ion containing a pentavalent carbon atom. In this sequence the hydride abstraction is the rate-determining step. The second deactivation function relates to the alkylation steps, which are also bimolecular and involve an aromatic molecule and an aromatic carbenium ion. This is the ultimate step yielding coke. The third deactivation function relates to a number of monomolecular elementary steps of carbenium ions such as methyl shift, PCP branching, deprotonation, β-scission, and cyclization. 8. Feedstock Composition

7. Effect of Coke on the Catalyst Activity The effect of coke formation is accounted for by multiplying the rates of the elementary steps by a deactivation function.18,19 It is expressed in terms of the coke content of the catalyst rather than the contact time, as is frequently done. The approach requires a kinetic equation for the coke formation, as developed in the present work. The exponential deactivation function exp(-RCC) is an approximation of the combined effect of site coverage and pore blockage. This functional form is generally observed in deactivation by coke formation but, in particular, in catalytic cracking.19-21

The kinetic model used here for the simulation of the catalytic cracking reactor requires a more detailed characterization of the feedstock than what is generally available and which is generally limited to geographical origin, ASTM boiling range, API gravity, Conradson carbon, UOP K factor, aniline point, and sulfur content. What is preferably required for an optimal application of the present approach is an analysis per C number of the hydrocarbon families of the VGO. There are very few examples of this in the literature. Table 2 shows such an analysis for a partially hydrogenated VGO, HVGO, obtained from the application of the GC-MS

Table 2. Composition of HVGO (wt %) no. of C

nPAR

14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 sum

0.002 0.016 0.077 0.167 0.365 0.729 1.01 1.43 1.56 1.86 2.07 1.65 1.27 0.9 0.538 0.343 0.183 0.115 0.04 14.3

iPAR

MNA

0.002 0.002 0.039 0.152 0.333 0.57 0.91 1.27 1.6 1.72 1.45 1.36 1.11 0.723 0.426 0.268 0.163 0.07 12.2

0.004 0.032 0.128 0.339 0.614 1.03 1.47 1.92 2.38 2.79 2.59 2.38 2.07 1.55 1.15 0.712 0.405 0.236 21.8

DNA

0.021 0.131 0.312 0.613 1.05 1.38 1.77 2.22 2.46 2.27 2.11 1.93 1.52 1.14 0.608 1.248 0.015 19.8

TNA

0.022 0.109 0.247 0.541 0.65 0.873 1.24 1.37 1.09 0.897 0.615 0.403 0.163 0.09 0.092 0.125 8.52

QNA

0.018 0.196 0.447 0.723 0.793 0.873 0.693 0.533 0.439 0.395 0.235 0.118 0.082 0.06 0.047 5.65

MAR

DAR

0.827 2.33 2.71 2.32 1.61 0.718 0.246 0.086 0.012

0.047 0.359 1.01 0.95 0.915 0.335 0.168 0.113 0.02

0.323 0.459 0.126 0.078 0.047 0.01

3.917

1.045

10.89

TAR

NMA

NDA

NTA

0.413 0.18 0.5 0.47 0.23 0.063 0.039 0.023

0.413

1.15

0.355

sum 1.287 2.871 4.242 4.195 3.703 2.537 3.263 4.917 6.605 8.186 10.17 11.1 9.583 8.456 7.02 4.969 3.34 1.943 1.083 0.533 100

Table 3. Analysis of Lagomedio VGO hydrocarbon

formula

content (wt %)

hydrocarbon

formula

content (wt %)

paraffins 1-ring naphthenes 2-ring naphthenes 3-ring naphthenes 4-ring naphthenes 5-ring naphthenes alkylbenzenes benzocycloparaffins

CnH2n+2 CnH2n CnH2n-2 CnH2n-4 CnH2n-6 CnH2n-8 CnH2n-6 CnH2n-8

11.0 15.2 11.0 6.5 2.8 0.2 4.0 5.5

naphthalenes naphthocycloparaffins/biphenyls triaromatics tetraaromatics thiophenes benzothiophenes dibenzothiophenes

CnH2n-12 CnH2n-14 CnH2n-18 CnH2n-24 CnH2n-4S CnH2n-10S CnH2n-16S

2.8 10.1 9.4 3.2 0.6 8.3 9.3

20

Ind. Eng. Chem. Res., Vol. 42, No. 1, 2003

Figure 7. Experimental and predicted distributions in HVGO fractions: (a) n- and i-paraffins; (b) naphthenes; (c) aromatics.

Ind. Eng. Chem. Res., Vol. 42, No. 1, 2003 21 Table 4. Values for the Single-Event Rate Coefficients at 923 K and for the Deactivation Constantsa no. 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 a

k (823 K)

component/lump

10-7

2.8887 × 0.0021 0.0251 20986 0.0267 0.0140 0.0004 0.0009 0.0103 0.0003 0.0301 7.3582 × 10-6 1.4716 × 10-6 4.1206 × 10-5 1.8353 × 10-6 0.0182 0.0560 0.2180 0.5390 0.0035 0.0913 0.0028 0.2656 2.2875 × 10-9 7.1056 × 10-8 1.0678 × 10-6 3.3825 × 10-5 0.0002 1.1183 × 10-5 2.3095 × 10-5 0.0007 4.0409 × 10-7 1.1440 2.3455 × 10-7 1.0288 × 10-7 0.0003 0.0017 0.0108 2.0163 × 10-5

acyclic

naphthenic

aromatic

elementary step

units

hydride transfer protonation (s) protonation (t) deprotonation (p) deprotonation (s) deprotonation (t) PCP branching (s,s) β-scission (s,s) β-scission (s,t) β-scission (t,s) β-scission (t,t) protolytic scission (p) protolytic scission (s) protolytic scission (t) hydride transfer (ring) protonation(r) (s) protonation(r) (t) deprotonation(r) (s) deprotonation(r) (t) β-scission(r) (s,s) β-scission(r) (s,t) β-scission(r) (t,s) β-scission(r) (t,t) protolytic scission(r) (p) protolytic scission(r) (s) protolytic scission(r) (t) aromatization (s) aromatization (t) dealkylation (s) dealkylation (t) cyclization disproportionation (r) protonation/deprotonation alkylation (s,s) alkylation (t,t) protolytic scission (p) protolytic scission (s) protolytic scission (t) hydride transfer

kmol/(kg of catalyst)‚h‚kPa kmol/(kg of catalyst)‚h‚kPa kmol/(kg of catalyst)‚h‚kPa kmol/(kg of catalyst)‚h kmol/(kg of catalyst)‚h kmol/(kg of catalyst)‚h kmol/(kg of catalyst)‚h kmol/(kg of catalyst)‚h kmol/(kg of catalyst)‚h kmol/(kg of catalyst)‚h kmol/(kg of catalyst)‚h kmol/(kg of catalyst)‚h‚kPa kmol/(kg of catalyst)‚h‚kPa kmol/(kg of catalyst)‚h‚kPa kmol/(kg of catalyst)‚h‚kPa kmol/(kg of catalyst)‚h‚kPa kmol/(kg of catalyst)‚h‚kPa kmol/(kg of catalyst)‚h kmol/(kg of catalyst)‚h kmol/(kg of catalyst)‚h kmol/(kg of catalyst)‚h kmol/(kg of catalyst)‚h kmol/(kg of catalyst)‚h kmol/(kg of catalyst)‚h‚kPa kmol/(kg of catalyst)‚h‚kPa kmol/(kg of catalyst)‚h‚kPa kmol/(kg of catalyst)‚h kmol/(kg of catalyst)‚h kmol/(kg of catalyst)‚h kmol/(kg of catalyst)‚h kmol/(kg of catalyst)‚h kmol/(kg of catalyst)‚h kmol/(kg of catalyst)‚h‚kPa kmol/(kg of catalyst)‚h‚kPa kmol/(kg of catalyst)‚h‚kPa kmol/(kg of catalyst)‚h‚kPa kmol/(kg of catalyst)‚h‚kPa kmol/(kg of catalyst)‚h‚kPa kmol/(kg of catalyst)‚h‚kPa

Deactivation constants: R1 ) 5.0; R2 ) 30; R3 ) 50 (kg of catalyst/kg of coke).

technique to samples of saturated and unsaturated fractions obtained by separation of the HVGO on a silica gel column with n-hexane as a mobile phase, a rather straightforward technique.12 Struppe et al.,25 Ramaswamy et al.,26 and Altgelt and Boduszynski27 also reported detailed information and analyses per carbon number for VGOs of various origins. Nevertheless, when such a detailed composition is not available, it has to be constructed from the available incomplete information by means of a model. The twoparameter Gaussian distribution was found to provide a good fit for the distribution of various hydrocarbon families (n-paraffins, mono-ring naphthenes, di-ring aromatics, etc.) in VGO as a function of the carbon number, nC. The mass fraction of a component with carbon number nC(i) belonging to a hydrocarbon family j is given by

xij )

Rj

x2πσj

e-0.5[nC(i)-µj/σj]

2

(29)

where µj and σj are parameters representing the mean and the variance of the distribution and Rj is the total content of the hydrocarbon family in the feedstock. The above equation was tested against detailed product distributions of various hydrocarbon families in VGOs. Both the mean, µj, and the variance, σj, were obtained by fitting the experimental compositions by means of an optimization routine, e.g., Marquardt.29

An example is given in Figure 6 which compares the actual and predicted distributions per C number of monoaromatics in a VGO.27 Figure 7 illustrates that even for a HVGO, the detailed composition of which is given in Table 1, the fit of the distributions is quite satisfactory, as confirmed by parity plots and Kolmogorov-Smirnov and Shapiro-Wilk statistical tests at the 95% confidence level. With less complete analyses, the task goes beyond the fitting of existing curves. An example is given in Table 3, which is an analysis by Sheppard et al.23 of Lagomedio VGO with a boiling range between 310 and 510 °C. This analysis does not distinguish between n- and i-paraffins, and only their global content is given. Naphthenes and aromatics are split per number of rings by means of high-resolution mass spectroscopy, but not per C number, so that the latter distributions, in particular the location of their maximum, have to be estimated while the global contents of Table 3 have to be satisfied. Additional insight, like boiling points and ranges, is helpful in this case. Figure 7 shows the distribution curves developed for the various hydrocarbon families of this Lagomedio VGO using eq 29. 9. Values for the Rate Parameters The single-event approach, assumptions, and thermodynamic constraints reduce the number of rate coefficients for the complete cracking network, including coke formation, to 39. These were mainly estimated on

22

Ind. Eng. Chem. Res., Vol. 42, No. 1, 2003

10. Simulation of an Industrial Riser Reactor The riser reactor was modeled using a one-dimensional heterogeneous model with plug flow and slip between the gas and solid phases.30 The steady-state continuity and energy equations for the components and lumps can be written as

usg

dCi

)

dz

(mGOcpGO + mscps)

dT

∑j rjFs(1 - (z))φC

(30)

)

dz

∑i ∑j rij(-∆H)ijFs(1 - (z))ΩφC

(31)

and for the coke content of the catalyst

ms

Figure 8. Predicted distribution for the various hydrocarbon families in Lagomedio VGO: (a) n- and i-paraffins; (b) naphthenes; (c) aromatics.

the basis of cracking data on the above-mentioned Lagomedio VGO23 and Maya VGO,24 for which a detailed cracking product analysis, required for a serious test of the cracking model, was available. Initial values for the rate parameters were taken from Dewachtere et al.12 These were obtained by applying the single-event approach to the cracking of selected key components but on a first-generation catalyst. The parameter estimation procedure involved both the Rosenbrock-Storey28 and Marquardt29 routines. The objective function consisted of a sum of squares of residuals of liquified petroleum gas (LPG), gasoline, light cycle oil (LCO), HFO, and coke yields but also the yields of the PIONA components of gasoline, per C number. The deactivation parameters were estimated simultaneously. Table 4 shows the values of the single-event rate coefficients and of the deactivation constants at 823 K for a Davison FX equilibrium catalyst. The single-event coefficients exhibit trends that are well established in carbenium ion chemistry from work on single components.22 The set leads to a good fit of the experimental data on the abovementioned VGO, thus confirming the validity of the single-event kinetic model for the prediction of the conversion and the detailed product distribution.

dCC dz

)

∑l rClFs(1 - (z))ΩφC

(32)

with the following boundary conditions: at z ) 0, Ci ) (Ci)0; Cc ) (Cc)0; T ) T0. The rj are the rate equations for the “global” reactions between the 670 components and lumps, constructed from the rates of the elementary steps in which these components and lumps are involved, including coke formation. φC is the appropriate deactivation function, depending upon the type of elementary step. The feedstock composition chosen for the major part of the simulations is the Lagomedio VGO of Table 2, with the distribution per C number shown in Figure 8.23 The liquid feed rate was 218 900 kg/h and the catalystto-oil ratio 8.5. The evolution of the gasoline, LPG, LCO, HFO, and coke yields along the riser is shown in Figure 9 for a bottom temperature of the riser of 823 K. The adiabatic temperature drop was on the order of 25 Κ. The profiles and exit values of the yields of commercial fractions shown in Figure 9a are in good agreement with the industrial trends and values.31,32 Figure 9b also shows the evolution of a somewhat detailed composition of the gasoline fraction. The normal paraffin content is very low and the olefins content high, in agreement with results of commercial catalytic cracking. Even this information can be further detailed if desired. Parts c and d of Figure 9 show the evolution of the composition of the LCO and HFO fraction, which could also be further detailed. The exit conversion [100 - (LCO wt % yield) - (HFO wt % yield)] amounts to 67.04. Increasing the bottom temperature of the riser to 843 K and to 876 K increases the conversion to 70.62 and to 76.10, while the corresponding gasoline yields are 45.05 and 45.39 wt %, respectively. The gasoline yield clearly levels off, while the LPG yield rises to 17.58 and 22.02, indicating what is called overcracking. Figure 10 shows the evolution of the deactivation functions. These level off from a certain conversion onward, like the coking yield of Figure 9a. This is also a well-known feature. It is obvious from the simulations that, for the operating conditions selected here, the major evolution occurs in the bottom half of the riser. Case 2 of Table 5 reports the exit yields and the gasoline composition from a simulation in which φ1 ) φ2 ) φ3 ) 1, meaning no deactivating effect of coke. The exit conversion obviously increases. A comparison of the various yields at the intermediate conversion of 67.04%,

Ind. Eng. Chem. Res., Vol. 42, No. 1, 2003 23

Figure 9. Catalytic cracking of Lagomedio VGO. Riser simulation. Evolution of the wt % yields of families and lumps along the adiabatic riser. Catalyst-to-oil weightt ratio ) 8.5. Bottom temperature of the riser: 823 K. Exit conversion: 67.04%. (a) Evolution of the main commercial fractions along the riser. (b) Composition of the gasoline fraction. (c) Yields of LCO naphthenes. (d) Yields of LCO aromatics: 1, mono-ring aromatics; 2, di-ring aromatics; 3, tri-ring aromatics. (e) Yields of HFO naphthenes: 1, mono-ring naphthenes; 2, di-ring naphthenes; 3, tri-ring naphthenes; 4, tetra-ring naphthenes. (f) Yields of HFO aromatics: 1, mono-ring aromatics; 2, di-ring aromatics; 3, tri-ring aromatics; 4, tetra-ring aromatics.

Figure 10. Catalytic cracking of Lagomedio VGO. Evolution of the deactivation functions along the adiabatic riser. Conditions of Figure 9.

corresponding to the exit conversion of the Lagomedio VGO cracking of Figure 9a-f, does not reveal significant differences, despite selective deactivation of the elementary steps, evidenced by Figure 10, but that is because the fractions are complex mixtures whose components all undergo the various types of elementary steps. Simulations were also carried out for an inlet pressure of 3 bar. It is seen from case 3 of Table 5 how the conversion increases by some 3% absolute, while the gasoline yield drops by 1% and the LPG and coke yields increase by some 2% absolute. Case 4 of Table 5 shows exit data for the cracking of HVGO of Table 2 and Figure 7. The coke yield is much lower than that simulated for Lagomedio, which is not surprising given the low aromatics content of the feedstock (17.8 wt %, including the naphthenoaromatics

24

Ind. Eng. Chem. Res., Vol. 42, No. 1, 2003

Table 5. Exit Conversion, wt % Yield of Commercial Fractions, and Gasoline Composition for Catalytic Cracking of VGOa conversion, % yield, wt % fuel gas LPG gasoline LCO HFO coke gasoline composition, wt % n-paraffins i-paraffins aromatics naphthenes olefins

case 1

case 2

case 3

case 4

case 5

67.04

71.49

70.24

63.54

68.94

2.79 15.11 44.34 15.51 17.45 4.81

3.27 17.35 45.65 14.59 13.92 5.23

1.88 18.01 44.13 13.38 16.37 6.24

3.66 14.23 45.43 22.31 14.15 0.22

2.52 14.72 44.97 12.06 19.01 6.74

3.84 29.57 35.45 15.59 16.55

4.19 31.44 33.96 16.06 14.34

3.03 30.33 36.51 11.36 18.78

7.41 47.78 12.35 21.54 10.92

2.92 22.78 44.51 13.46 16.33

a Case 1: base case, Lagomedio VGO; T bottom of riser, 823 K; total pressure, 1 bar; deactivation by coke (viz. Figure 9). Case 2: Lagomedio VGO; no deactivating effect of coke. Case 3: Lagomedio VGO; total pressure, 3 bar. Case 4: HVGO; operating conditions as in case 1. Case 5: Maya VGO; operating conditions as in case 1.

vs 53.2%). For equal conversion, the gasoline yield from a Maya VGO with boiling range of 343-499 °C, comprising 10.1 wt % paraffins, 28.2 wt % naphthenes, and 61.7 wt % aromatics and naphthenoaromatics, was 2% lower than that of Lagomedio, but the coke was multiplied by a factor of 1.4, as can be seen from case 5 in Table 5. 11. Conclusion It is clear that the very fundamental and detailed single-event model presented here can be applied to the simulation of industrial catalytic cracking of complex feedstocks such as VGO. Through its fundamental nature, the basic rate coefficients of the model are invariant with respect to the feedstock composition, a tremendous benefit for refinery operation where the feedstock is frequently modified. Evidently, the model provides a wealth of detailed information on the effect of the operating conditions and feedstock composition on the exit product distribution and should thus be a powerful tool for guiding the refinery catalytic cracking. Evidently, more work remains to be done. The model does not consider the transformation of the sulfur- and nitrogen-containing components and the effect of their decomposition products on the catalyst activity. Also, today’s catalytic cracking incorporates substantial fractions of atmospheric residue to the cracker feed. To account for this, the behavior of very different molecules has to be accounted for. Finally, heavy feedstocks also contain vanadium and nickel complexes, which undergo transformations before depositing the metal on the catalyst and thus deactivating it and/or affecting the selectivities. It is believed that enough information is available in the literature, be it fragmented, to also enable the incorporation of these aspects and thus lead to a comprehensive model that would be sufficiently realistic to guide the operation of tomorrow’s cracking, even online. Nomenclature CC ) coke content of catalyst, kg of coke/kg of catalyst cpg, cps ) specific heat of gas or solid, kJ/kg‚K CR+ ) concentration of carbenium ion

Ea ) activation energy for chemical transformation, kJ/ kmol h ) Planck’s constant, kmol‚m2/s ∆Hr ) reaction enthalpy, kJ/kmol ∆H° ) difference of standard enthalpy between the reactant and activated complex, kJ/kmol k, k′ ) rate coefficient of respectively a single event or elementary step, kmol/(kg of catalyst)‚h‚kPa; h-1; kPa-1 h-1 kB ) Boltzmann constant, kmol‚m2/s2‚K LC ) lumping coefficient ms ) total mass flow rate of solid, kg/s ne ) number of single events in an elementary step p ) partial pressure, kPa rj ) rate of reaction j, kmol/(kg of catalyst)‚s rC ) rate of coke formation, kg of coke/(kg of catalyst)‚s ∆S° ) difference of standard entropy between the reactant and activated complex, kJ/kmol‚K T ) absolute temperature, K usg ) superficial gas flow velocity, mf3 /mr2‚s z ) axial distance in the riser, mr Greek Letters (z) ) void fraction in the riser, mf3/mr3 Fs ) specific mass of the solid, kg of catalyst/m3 of catalyst σgl,r, σgl ) global symmetry number for respectively the reactant and activated complex Ω ) cross-sectional area of the riser, mr2 Subscripts ht ) hydride transfer Pr ) protonation Proto ) protolytic scission s ) secondary carbenium ion t ) tertiary carbenium ion β ) β-scission

Literature Cited (1) Jacob, S. M.; Gross, B.; Voltz, S. E.; Weeckman, V. W. A Lumping and Reaction Scheme for Catalytic Cracking. AIChE J. 1976, 22, 701. (2) Liguras, D. K.; Allen, D. T. Structural Models for Catalytic Cracking: Reactions of Simulated Oil Mixtures. Ind Eng. Chem Res. 1989, 28, 674. (3) Klein, M. T.; Neurock, M.; Nigam, A.; Libanati, C. Monte Carlo Modeling of Complex Reaction Systems: An Asphaltene Example. In Chemical Reactions in Complex MixturessThe Mobil Workshop; Sapre, A. V., Krambeck, F. J., Eds.; Van Nostrand Reinhold: New York, 1991; pp 126-142. (4) Watson, B. A.; Klein, M. T.; Harding, R. H. Mechanistic Modeling of n-heptane Cracking on ZSM-5. Ind Eng. Chem. Res. 1996, 35, 1506. (5) Quann, R. J.; Jaffe, S. B. Structure Oriented Lumping Describing the Chemistry of Complex Hydrocarbon Mixtures. Ind Eng. Chem. Res. 1992, 31, 2483. (6) Quann, R. J.; Jaffe, S. B. Building Useful Models of Complex Reaction Systems in Petroleum Refining. Chem. Eng. Sci. 1996, 31, 1615. (7) Baltanas, M.; Van Raemdonck, K.; Froment, G. F.; Mohedas, R. S. Fundamental Kinetic Modeling of Hydroisomerization and Hydrocracking on Noble Metal Loaded Faujasites. Ind Eng. Chem. Res. 1989, 28, 899. (8) Vynckier, E.; Froment, G. F. Modeling of the Kinetics of Complex Processes based upon Elementary Steps. In Kinetic and Thermodynamic Lumping of Multicomponent Mixtures; Astarita, G., Sandler, S. I., Eds.; Elsevier: Amsterdam, The Netherlands, 1991. (9) Feng, W.; Vynckier, E.; Froment, G. F. Single Event Kinetics of Catalytic Cracking. Ind. Eng. Chem. Res. 1993, 32, 2997. (10) Svoboda, G. D.; Vynckier, E.; Debrabandere, B.; Froment, G. F. Single Event Rate Parameters for Paraffins Hydrocracking on Pt/US-Y Zeolite. Ind. Eng. Chem. Res. 1995, 34, 3793. (11) Baltanas, M.; Froment, G. F. Computer Generation of Reaction Networks and Calculation of Product Distributions in

Ind. Eng. Chem. Res., Vol. 42, No. 1, 2003 25 the Hydroisomerization and Hydrocracking of Paraffins on Ptcontaining Bifunctional Catalysts. Comput. Chem. Eng. 1985, 9, 71. (12) Dewachtere, N.; Santaella, F.; Froment, G. F. Application of a Single Event Kinetic Model in the Simulation of an Industrial Riser Reactor for the Catalytic Cracking of Vacuum Gas Oil. Chem. Eng. Sci. 1999, 54, 3653. (13) Froment, G F. Kinetic Modeling of Acid-Catalyzed Oil Refining Processes. Catal. Today 1999, 52, 153. (14) Park, T.-Y.; Froment, G. F. Kinetic Modeling of the MTO ProcesssI. Model Formulation. Ind. Eng. Chem. Res. 2001, 40, 4172. (15) Park, T.-Y.; Froment, G. F. Kinetic Modeling of the MTO ProcesssII. Experimental Results Model Discrimination and Parameter Estimation. Ind Eng. Chem. Res. 2001, 40, 4187. (16) Froment, G. F. Kinetic Modeling of Complex Catalytic Reactions. Rev. Inst. Fr. Pe´ t. 1991, 40, 491. (17) Guisnet, M.; Magnoux, P.; Martin, D. Roles of Acidity and Pore Structure in the Deactivation of Zeolites by Carbonaceous Deposits. In Catalyst Deactivation 1997; Bartholomew, C. H., Fuentes, G. A., Eds.; Studies in Surface Science and Catalysis 111; Elsevier: New York, 1997; p 1. (18) Froment, G. F. The Modeling of Catalyst Deactivation by Coke Formation. In Catalyst Deactivation 1991; Bartholomew, C., Butt, J. B., Eds.; Elsevier: New York, 1991; pp 53-83. (19) Froment, G. F. Coking and Deactivation in Complex Catalytic Processes. Proceedings of the XIth International Conference on Chemical Reactors, Institute of Catalysis, Novosibirsk, Jun 18-21, 1996; Abstracts, Part 1, pp 24-45. (20) Froment, G. F.; De Meyer, J.; Derouane, E. G. Deactivation of Zeolite Catalysts by Coke Formation. J. Catal. 1990, 124, 391. (21) Beirnaert, H. C.; Vermeulen, R.; Froment, G. F. A Recycle Electrobalance Reactor for the Study of Catalyst Deactivation by Coke Formation. Stud. Surf. Sci. Catal. 1994, 88, 97-112. (22) Martens, J. A.; Jacobs, P. A. Conceptual Background for the Conversion of Hydrocarbons on Heterogeneous Acid Catalysts. In Theoretical Aspects of Heterogeneous Catalysis; Moffat, J. B., Ed.; Van Nostrand Reinhold: New York, 1990; p 52.

(23) Sheppard, C. M.; Green, J. B.; Vanderveen, J. W. Relating Feedstock Composition to Product Slate and Composition in Catalytic Cracking. 4. Feedstocks derived from Lagomedio Crude. Energy Fuels 1998, 12, 320-328. (24) Green, J. B.; Zagula, E. J.; Reynolds, J. W.; Young, L. L.; McWilliams, T. B.; Green, J. A. Relating Feedstock Composition to Product Slate and Composition in Catalytic Cracking. 3. Feedstocks derived from Maya, a Mexican Crude. Energy Fuels 1997, 11, 46-60. (25) Struppe, H. G.; Janke, H.; Deutsch, K.; Grunov, S. The general physicochemical characteristics and group-structure hydrocarbon composition of vacuum gas oil (350-540 °C) of industrial West Siberia crude oil. Pet. Chem. 1988, 1, 1-11. (26) Ramaswamy, V.; Singh, D. S.; Krishna, R. Characterization of Vacuum Gas Oil from North Gujarat crude mix. Indian J. Technol. 1989, 27, 85-88. (27) Algelt, K. H.; Boduszynski, M. M. Composition and Analysis of Heavy Petroleum Fractions, 1st ed.; Marcel Dekker: New York, 1994. (28) Rosenbrock, H. H.; Storey, C. Computational Techniques for Chemical Engineers; Pergamon Press: New York, 1996. (29) Marquardt, D. An Algorithm for Least-Squares Estimation of Nonlinear Parameters. J. Soc. Ind. Appl. Math. 1963, 11, 431. (30) Froment, G. F.; Bischoff, K. B. Chemical Reactor Analysis and Design, 2nd ed.; John Wiley: New York, 1990. (31) Young, G. W. Realistic Assessment of FCC catalyst performance in the Laboratory. In Fluid Catalytic Cracking; Science and Methodology; Magee, J. S., Mitchell, M. M., Jr., Eds.; Studies in Surface Science and Catalysis 76; Elsevier: New York, 1993. (32) Decroocq, D. Catalytic Cracking of Heavy Petroleum Fractions; Technip: Paris, 1984.

Received for review June 17, 2002 Revised manuscript received October 7, 2002 Accepted October 10, 2002 IE0204538