Analytical Solution for Discrete Lumped Kinetic Equations in

6608

Ind. Eng. Chem. Res. 2009, 48, 6608–6617

Analytical Solution for Discrete Lumped Kinetic Equations in Hydrocracking of Heavier Petroleum Fractions P. C. Krishna and P. Balasubramanian* Department of Chemical Engineering, Indian Institute of Technology Guwahati, Guwahati 781 039, India

This paper provides a general analytical solution for the full stoichiometry based discrete lumped kinetic model. The reduced stoichiometry and carbon number basis are the subset of the full stoichiometry model. The exponential form of stoichiometric kernels is considered to study the product distribution for the full and reduced stoichiometry models. Here, full stoichiometry favors for the formation of heavier lumps, while the reduced stoichiometry favors for the formation of lighter lumps. The exact solution for the discrete lumped model is validated with the experimental data cited in the literature. The model validation reveals that the full stoichiometry is the best performing model for hydrocracking of heavier petroleum fractions. 1. Introduction Hydrocracking is aimed at the conversion of high boiling point petroleum fractions such as vacuum gas oil into the low boiling point fractions such as gases, naphtha, kerosene, diesel, and gas oil in the presence of hydrogen gas. Normally, the hydrocracking process is carried out in a high pressure trickle bed reactor. Kinetic modeling plays vital role in hydrocracking to predict the performance of the process in terms of conversion of the feedstock and selectivity for the products and, also, to optimize the process variables to meet the market demand. Several researchers have proposed different kinetic models for hydrocracking,1-35 and the models were broadly classified into two categories namely, lumped kinetic modeling and fundamental kinetic modeling methods. The lumped kinetic modeling approach is further classified into continuous lumping,1-4 and discrete lumping1,4-24 methods. In the former, the product distributioninahydrocrackerisgovernedbyanintegral-differential equation, and the latter is governed by ordinary differential equations. Structure oriented lumping modeling25-27 and single event kinetic modeling28-35 methods belong to the fundamental kinetic modeling of the hydrocracking process. A discrete lumped kinetic model is applied to predict the performance of the hydrocracker unit in a refinery. In discrete lumping, the feedstock and products are classified into lumps based on either the true boiling point or molecular weight range of the hydrocarbons. The discrete lumped models published in the literature were based on the component identity.5-17,19,20 That is, the stoichiometry of the macroscopic cracking reactions was proposed based on the molecules present in the heavier lumps which undergo cracking reactions, and it yields the product in the lighter lumps. Recently, discrete lumped kinetic models for hydrocracking of vacuum gas oil were developed from continuous kinetics using carbon number and true boiling point (TBP) as the basis.4 Binary scission was considered to model the hydrocracking process. In the former, fixing the first product uniquely determines the second one and the two products are formed independently in the latter basis. Stoichiometric kernels (random scission, symmetric, and exponential form) were proposed to identify the nature of the cracking reactions and to study the product distribution in a hydrocracker. Here, the TBP basis discrete lumped model with exponential form of stoichiometric kernel exhibits the best performing model * To whom correspondence should be addressed. Tel.: +91-3612582263. Fax: +91-361-2690762. E-mail: [email protected].

for hydrocracking of vacuum gas oil as compared with the carbon number basis model. In the earlier studies, the so-called full stoichiometry (cracking reactions within the lumps that are considered) is not included. Also, the analytical solution to determine the weight fraction of the cracked products is not derived. Therefore, in this paper, a general TBP basis discrete lumped model for hydrocracking is presented to account for the cracking reactions occurring within the lumps. The general kinetic model is also applicable for the reduced stoichiometry (cracking reactions within the lumps that are not considered) and the carbon number basis. The exponential form of stoichiometric kernels is considered to study the product distribution for the full and reduced stoichiometry models. Thus, the objectives are (i) to derive a general exact solution for the governing mass balance equations, (ii) to validate the performance of the full stoichiometry model with the exponential form of stoichiometric kernel using the data cited in the literature, and (iii) to compare the reliability of the full stoichiometry model with the reduced stoichiometry model. 2. Discrete Lumping Approach In the following, the derivation for the TBP basis discrete lumped model and the analytical solution for the governing mass balance equations are presented. 2.1. Stoichiometry. A general full stoichiometry of the hydrocracking reactions is represented as ki,j,r

cr f ci + cj

(1)

where r varies from 1 to NL, i and j vary from 1 to r, NL is the number of lumps considered, cr is the molar concentration of the hydrocarbons present in the lump r, and ki,j,r represents the kinetic constant for cracking of reactant lump r into two products, lumps i and j. In hydrocracking, the discrete lumps are defined based upon wide range of true boiling point of the hydrocarbons. For example, the typical hydrocracker product fractions defined in terms of true boiling point range33 are liquefied petroleum gas (LPG) (620 K). Here, in general, the high boiling point hydrocarbons present in the residue lump react to yield the two products in the low boiling point lumps such as LPG, naphtha, and middle distillates. However, in reality, there

10.1021/ie900178m CCC: $40.75  2009 American Chemical Society Published on Web 06/09/2009

Ind. Eng. Chem. Res., Vol. 48, No. 14, 2009

is a possibility of either the formation of two products in the residue lump itself or one product in the residue lump and the other product in the low boiling range lumps as a result of the structure of hydrocarbons present in the residue lump, and also, a high reaction temperature is required to crack the heavier hydrocarbons into the low boiling point fractions. Similarly, the cracking reactions can also occur within the middle distillates and naphtha boiling point range lumps. Therefore, the indices i and j vary from 1 to r in eq 1, which accounts for the cracking reactions that occur within the lumps. The aforementioned stoichiometry can also be applied for reduced stoichiometry4,18 and carbon number basis4,18 models. In the former, r varies from 2 to NL and i and j vary from 1 to r - 1. For the latter case, r varies from 2 to NL, i varies from 1 to r - 1, and j ) r - i. In the present work, it is assumed that the cracking is a binary process in which only two products are formed in each cracking event. Also, the assumption include the cracking reactions to be of first order and irreversible. Furthermore, coke formation kinetics, isomerization reactions on acid sites, and the reactions which occur on metal sites such as (de)hydrogenation are neglected. 2.2. Number of Kinetic Constants. The number of kinetic constants Nk involved in the three models is calculated using the following formulas: Full stoichiometry model18 Nk )

2M(2M + 1)(2M + 2) 6

when NLis even and NL ) 2M (2)

Nk ) (2M + 1)(2M + 2)(2M + 3) 6

2M(2M + 1)(2M - 1) 6 2M(2M + 1)(2M + 2) 6

kNL-1,j,NL

cNL f cNL-1 + cj where j varies from 1 to NL - 1 kNL-2,j,NL

cNL f cNL-2 + cj where j varies from 1 to NL - 2 and so on k1,1,NL

cNL f 2c1

Thus, the number of kinetic constants for all possible reactions from the reactant lump NL is reduced to the value of NL from (NL + 1)NL/2. Furthermore, the aforementioned parameter reduction method is also applicable for the reactant lumps from NL - 1 to 1. The formulas to determine the number of kinetic constants for the full stoichiometry model are given by 2M(2M + 1) 2

when NL is even and NL ) 2M (9)

Nk )

when NL is even and NL ) 2M

when NL is odd and NL ) 2M + 1

(8b)

(3)

(4) Nk )

(8a)

In eq 8a, the kinetic constants kNL,j,NL are assumed to be constant values and, thus, the number of kinetic constants Nk is reduced to 1 from NL. Similarly, the kinetic constants for the reactant lump NL, which yield the products in the lumps NL - 1 to 1 for the following stoichiometry, are assumed as constant values.

Nk )

Reduced stoichiometry model Nk )

cNL f cNL + cj

whenNL is odd and NL ) 2M + 1

6609

kNL,j,NL

(2M + 1)(2M + 2) when NL is odd and NL ) 2M + 1 2 (10)

In the same analogy, the formulas to determine the number of kinetic constants for the reduced stoichiometry model are given by

(5) Nk )

Carbon number basis model18

2M(2M - 1) when NL is even and NL ) 2M 2 (11)

Nk ) M + M2 when NL is even and NL ) 2M + 1 (6) Nk ) M + M2 when NL is odd andNL ) 2M + 1

(7)

The number of kinetic constants involved in the abovementioned three stoichiometry models is detailed in Table 1 for NL varying from 3 to 10. Obviously, the number of kinetic constants Nk increases with an increase in number of lumps NL considered. 2.3. Kinetic Parameter Reduction. The following assumptions4 are made to reduce the number of kinetic constants Nk for the full and reduced stoichoimetry models. In the full stoichiometry model, the number of kinetic constants is NL when the reactant lump NL yielding the first product in the reactant lump NL and the second product in the lump j (here, j varies from 1 to NL). That is, the stoichiometry of the cracking reactions is represented as

Table 1. Total Number of Kinetic Constants Involved in the Discrete Lumped Model for the Three Stoichiometry Cases Considered Nk reduced stoichiometry

full stoichiometry

NL

carbon number basis

case 1a

case 2b

case 1a

case 2b

3 4 5 6 7 8 9 10

2 4 6 9 12 16 20 25

4 10 20 35 56 84 120 165

3 6 10 15 21 28 36 45

10 20 35 56 84 120 165 220

6 10 15 21 28 36 45 55

a Nk value is obtained for the full and reduced stoichiometry by using eqs 2 and 3 and 4 and 5, respectively. b Nk value is obtained for the full and reduced stoichiometry by using eqs 9 and 10 and 11 and 12, respectively.

6610


Nk )

2M(2M + 1) 2

when NL is odd and NL ) 2M + 1 (12)

Table 1 presents the number of kinetic constants for the full and reduced stoichiometry models using eqs 9 and 10 and 11 and 12, respectively. 2.4. Mass Balance Equations. The molar concentration of the hydrocarbons cr in the lump r is governed by the following ordinary differential equations when first order irreversible reaction kinetics is considered. N

j

L dcr )2 dt j)r

∑ ∑ Ω(r, i, j)k

r

r,i,jcj

-

i)1

1 to NL. The exponential form of kernels for the full and reduced stoichiometry models and the symmetric kernel for the carbon number basis model are given in Appendix A. 2.5. Analytical Solution. The analytical expressions are proposed to determine the weight fraction of lumps by using Laplace transforms. The proof of analytical solution using Laplace transforms36 is detailed in Appendix B. The analytical solution for the weight fraction of lumps can also be obtained by recursively integrating1 eq 15 starting r from NL to 1. Therefore, a general exact solution for the weight fraction of lump r is represented as

r

∑ ∑ Ω(i, j, r)k

NL

i,j,rcr

wr )

i)1 j)1

(13) In eq 13, Ω represents the exponential form of stoichiometric kernel for the full and reduced stoichiometry models. Here, full stoichiometry favors the formation of heavier lumps while the reduced stoichiometry favors the formation of lighter lumps. The stoichiometric kernel exhibits symmetric for the carbon number basis model. Also, Ω(r,i,j) represents the distribution of the products in the lumps r and i, from the reactant lump j by virtue of cracking reaction. Furthermore, the sum of molar concentration of the hydrocarbons present in all the lumps starting from 1 to NL is not unity and increases with time. That NL NL ci * constant for all time instances. However, ∑i)1 ici is, ∑i)1 ) constant for the carbon number basis model and it is not a constant for the full and reduced stoichiometry models. In hydrocracking, the experimental data is available in the form of weight fraction of the discrete lumps. The weight fraction of the hydrocarbons wr in the lump r from molar concentration of the hydrocarbons cr is calculated using the relationship

( )

j

∑ ∑δ i)1

-

(16)

m

Rm ) 2

∑δ

m,j,mΩ(m, j, m)km,j,m

(17)

j)1

m

∑ ∑ Ω(i, j, m)k

βm )

(18)

i,j,m

i)1 j)1

r

r,i,jΩ(r, i, j)kr,i,jwj

exp[(Rm - βm)t]

where r varies from NL to 1. In eq 16, D is the time independent NL × NL coefficient matrix and is a function of initial weight fraction of the lumps and the kinetic constant values. The D coefficients can vary in the following two ways, namely, (i) the choice of initial weight fraction of the lumps with the fixed kinetic constants and (ii) the choice of kinetic constants with the fixed initial weight fraction of the lumps. Here, the term Rm describes the sum of kinetic constants for cracking of lump m which yield the products in the lumps m and j (here, j varies from 1 to m). Similarly, βm represents the sum of kinetic constants for the cracking of lump m which give products in the lumps ranging from 1 to m. The expressions for the terms Rm and βm are represented as

(14)

In eq 14, F0 is the density of the mixture and Mr is the molecular weight of the hydrocarbons in the lump r. Here, the density of the mixture and the molecular weight of the lumps are considered constant. The weight fraction of the hydrocarbons in the lump r is obtained by substituting eq 14 in eq 13, and the resulting expression is given by N

r,m

m

Mr c wr ) F0 r

L dwr )2 dt j)r

∑D

m)r

r

∑ ∑ Ω(i, j, r)k

i,j,rwr

i)1 j)1

(15) In eq 15, δr,i,j ) Mr/Mj and it is reasonably approximated as r/(i + r) to normalize the weight fraction of the products added to the lumps r and i from the reactant lump j by virtue of cracking reaction. The coefficient δr,i,j are conveniently represented using three indices because the value of coefficients (δr,i,j and δi,r,j) for the weight fraction of the products added to the lumps r and i from the reactant lump j is r/(i + r) and i/(i + r), respectively. Here, the sum of the two indices (δr,i,j and δi,r,j) is equal to unity. Also, the value of the coefficient δi,i,j is 0.5 for the weight fraction of two products added to the lump i from the reactant lump j. Thus, the sum of weight fraction of all the lumps obtained from the mass balance eq 15 is equal to unity at all time instances. In eqs 13 and 15, it is assumed that kr,i,j ) ki,r,j implying the symmetry of the kinetic constants. Also, the index r varies from

[

]

The matrix form of the D coefficients is given by · · · D1,NL · · · D2,NL · · · D3,NL l · · · DNL,NL

D1,1 D1,2 D1,3 D2,2 D2,3 0 D3,3 0 0 l l l 0 0 0

(19)

The elements of D coefficient matrix are determined by using the following formulas. The expression to determine the main diagonal elements in the D coefficient matrix is represented as NL

Dr,r ) wr,0 -

∑D

r,j

(20)

j)r+1

The formula to determine the superdiagonal elements in the D coefficient matrix is given by 2 Dr,m )

m

j

j)r+1

i)1

∑ ∑δ

r,i,jΩ(r, i, j)kr,i,jDj,m

(βr - βm) + (Rm - Rr)

The D coefficient matrix must be filled starting from ending with D1,1.

(21) NL,NLa

nd


Finally, the subdiagonal elements in the D coefficient matrix are expressed as Dr,m ) 0

(22)

2 w1 ) w1,0 +

where r varies from m + 1 to NL and m varies from 1 to NL 1. 3. Illustrative Example The applicability of the proposed analytical solution is demonstrated with NL ) 3 by considering the full and reduced stoichiometry models with the exponential form of the stoichiometric kernels. The discrete lumped models are validated with two sets of experimental data available in the literature.5,37 The full stoichiometry of the three lump kinetic model is given by k3,1,3

k2,1,2

k3,2,3

k2,2,2

k3,3,3

k1,1,2

∑δ

(

(β1 - β2) + (R2 - R1) 3

2

c2 f 2c1

(23)

c3 f c1 + c2 k2,2,3

c3 f 2c2 k1,1,3

c3 f 2c1 In eq 23, it is assumed that ki,j,r ) kj,i,r implying the symmetry of the kinetic constants. The definition of lumps considered in the model validation using El-Kady37 experimental data for lump 1, lump 2, and lump 3 in terms of true boiling point range are 0-315, 315-523, and 523-823 K, respectively. Here, in general, the hydrocarbons present in the heavier lump 3 react to yield the products in lighter lumps 2 and 1. However, in reality, there is a possibility of either the formation of two products in the heavier lump 3 or one product in the heavier lump 3 and the other product in the lighter lumps as a result of the structure of hydrocarbons present in the heavier lump, and also, a high reaction temperature is required to crack the heavier hydrocarbons into the low boiling point fractions. Similarly, the cracking reactions can also occur within the lighter lumps. Thus, in this work, the so-called full stoichiometry model is considered to study the product distribution in a hydrocracker with the minimum number of kinetic constants. The analytical expression for the weight fraction of the three lumps is derived using eqs 16-18, 20, and 21 and the details of the derivations are given in Appendix C. The final analytical expressions are shown here. The weight fraction of lump 3 is given by w3 ) w3,0 exp[(R3 - β3)t] 3 ) 1δ3,j,3Ω(3,j,3)k3,j,3

and β3 ) ∑i

(24) 3 ) 1 ∑j )

The weight fraction of lump 2 is represented as 3

2 w2 )

2,j,3Ω(2, j, 3)k2,j,3w3,0

(β2 - β3) + (R3 - R2) (exp[(R2 - β2)t] - 1) +

(

2

2

∑δ

3

3

2

∑δ

∑δ

2

1,j,2Ω(1, j, 2)k1,j,2

j)1

)

2,j,3Ω(2, j, 3)k2,j,3w3,0

j)1

(β2 - β3) + (R3 - R2)

)

+

1,j,3Ω(1, j, 3)k1,j,3w3,0

j)1

(β1 - β3) + (R3 - R1) (exp[(R3 - β3)t] - 1)

k2,1,3

2∑j

∑δ j)1

w2,0 -

k1,1,1

c3 f c2 + c3 c2 f 2c2

where R3 ) 3 1Ω(i,j,3)ki,j,3

1,j,2Ω(1, j, 2)k1,j,2

j)1

c3 f c1 + c3 c2 f c1 + c2 c1 f 2c1

c3 f 2c3

6611

2

∑δ

2,j,3Ω(2, j, 3)k2,j,3w3,0

j)1

(exp[(R3 - β3)t] (β2 - β3) + (R3 - R2) exp[(R2 - β2)t]) + w2,0exp[(R2 - β2)t]

(25)

2 2 2 where R2 ) 2∑j)1 δ2,j,2Ω(2,j,2)k2,j,2 and β2 ) ∑i)1 ∑j)1 Ω(i,j,2)ki,j,2. Furthermore, the weight fraction of lump 1 is given by

(26)

where R1 ) β1 ) k1,1,1. Here, the sum of weight fraction of the three lumps from eqs 24-26 is equal to the sum of initial weight fraction of the three lumps as a result of the values of the coefficient δr,i,j. That is, the sum of weight fraction of the three lumps is always equal to unity. The analytical weight fraction expressions for the reduced stoichiometry model is obtained by substituting Ri ) 0 (i ) 1, 2, and 3) in eqs 24-26, and for the carbon number basis model, the weight fraction expressions are deduced by considering the value of kinetic constants k2,1,3 and k1,1,2. 3.1. Parameter Estimation The proposed analytical solution for the full and reduced stoichiometry models is validated with the experimental data of El-Kady37 and Zhorov et al.5 available in the literature. The objective function considered in the parameter optimization is to minimize the residual sum of squares error (RSSE) between the experimental data and the model calculated values. It is represented as N

f(k) )

NL

∑ ∑ (w

exp(i, j)

- wmodel(i, j))2

(27)

i)1 j)1

In eq 27, N is the number of observation, wexp is the experimental weight fraction data, and wmodel is the model calculated weight fraction values. Here, unweighed regression is used for parameter optimization. That is, equal weighing is applied to all responses. The parameter optimization method is implemented using fmincon routine available in MATLAB software package. The number of observations N and number of lumps NL considered from El-Kady37 experimental data are 6 and 3, respectively. The number of kinetic parameters for the full stoichiometry model is 6, and the initial weight fraction is 0, 0, and 1 for lump 1, lump 2, and lump 3, respectively. The space time (in hours) for the six observation is 0.3333, 0.4, 0.5, 0.6667, 1, and 2. The kinetic constants involved in the full stoichiometry model are estimated through unweighed regression at four reaction temperatures (in kelvin) 663, 683, 703, and 723. The estimated kinetic parameters are shown in Table 2. It is observed that the kinetic constants for the cracking reactions from lump 3 increase with an increase in temperature from 663 to 723 K. The value of kinetic constants is doubled at 723 K. The estimated kinetic constants k3,1,3, k3,2,3, and k3,3,3 are significant at four reaction temperatures, and this confirms that reactions have occurred in lump 3. Here, the exponential form of the

6612


Table 2. Estimated Kinetic Constants for the Four Reaction Temperatures using El-Kady37 Experimental Data temperature (K) kinetic constant (h-1)

663

683

703

723

663

full stoichiometry k1,1,1 k1,1,2 k2,1,2 k2,2,2 k1,1,3 k2,1,3 k2,2,3 k3,1,3 k3,2,3 k3,3,3 RSSE

0.001 0.008 0.045 0.045 0.090 0.459 0.459 0.354 0.354 0.354 0.001

0.001 0.112 0.183 0.183 0.201 0.772 0.772 0.639 0.639 0.639 0.005

0.001 0.191 0.289 0.289 0.296 1.341 1.341 1.100 1.100 1.100 0.015

stoichiometric kernel favors the formation of heavier products in lump 3. The estimated kinetic constants (k2,1,3 and k2,2,3) for the formation of the products in the low boiling point range lumps are high as compared with the kinetic constants (k3,1,3, k3,2,3, and k3,3,3) for the formation of heavier products in lump 3. At low reaction temperature, the high boiling point hydro-

683

703

723

reduced stoichiometry 0.001 0.147 0.230 0.230 0.408 2.859 2.859 2.252 2.252 2.252 0.004

0.001

0.001

0.001

0.001

0.001 0.305 0.305

0.001 0.526 0.526

0.001 0.906 0.906

0.001 1.910 1.910

0.002

0.007

0.016

0.011

carbons are formed in lump 3, and subsequently, it undergo cracking reactions to yield products in the low boiling point range lumps. The possibility of formation of products in lump 1 directly from reactant lump 3 is low. High reaction temperature (723 K) favors the formation of lighter products in lumps 1 and 2 from lump 3.

Figure 1. Parity diagrams of the three lumps obtained for the full and reduced stoichiometry models with the exponential form of stoichiometric kernels using the experimental data of El-Kady.37 (a and d) Lump 1, (b and e) lump 2, and (c and f) lump 3. (a-c) Full stoichiometry model. (d-f) Reduced stoichiometry model. Thin solid line represents the trend line.


6613

Figure 2. Weight fraction versus space time plots of the three lumps obtained for the full and reduced stoichiometry models in hydrocracking using the experimental data of El-Kady.37 (a and d) Lump 1, (b and e) lump 2, and (c and f) lump 3. (a-c) Full stoichiometry model. (d-f) Reduced stoichiometry model. (symbol) El-Kady experimental data.37 (solid line) Model calculated values obtained by using the exact solution eqs 16-18, 20, and 21.

The values of kinetic constants (k2,1,2 and k2,2,2) for the reactions from lump 2 increase with an increase in temperature from 663 to 703 K, and theses kinetic constants suddenly decrease with an increase in temperature from 703 to 723 K. This is due to the limitations of the discrete lumped model. The estimated kinetic constant (k1,1,1) for the reaction in the lump 1 is 0.001 at all the four temperatures. It signifies that there is no further reaction occurring in lump 1 because the boiling point range of the lump 1 is less than 315 K, and the lighter products are formed from the heavier hydrocarbons present in lumps 2 and 3. The kinetic parameters are also estimated for the reduced stoichiometry model and they are shown in Table 2. It is observed that the values of kinetic constants (k2,1,3 and k2,2,3) increases with an increase in temperature from 663 to 723 K. The estimated kinetic constants for k1,1,2 and k1,1,3 are 0.001 at all the four temperatures. Here, the lighter products are formed in lumps 1 and 2 from the heavier lump 3 and the reduced stoichiometry with exponential form of stoichiometric kernel favors the formation of products in the lighter lumps. The residual sum of squares errors for the reduced stoichiometry model is high as compared with the full stoichiometry model. The parity diagrams of the weight fraction of lumps for the full and reduced stoichiometry models are shown in Figure 1. The slope and intercept values of the parity diagrams for the three lumps are also given in Figure 1. For the ideal case, the slope and intercept values should be equal to 1 and 0, respectively.

For the full stoichiometry model, the slope and intercept values of the parity diagrams are close to the ideal case. The weight fraction of the lumps versus space time plots for the three lumps at four reaction temperatures are depicted in Figure 2 for the full and reduced stoichiometry models. In this model validation, the full stoichiometry with the exponential form of stoichiometric kernel shows best fit with the data as compared with the reduced stoichiometry model as a result of low RSSE. For better parameter estimates, the number of observation should be greater than or equal to the number of parameters to be estimated through regression. In the following, the performance of the discrete lumped model is described with the limited experimental data. That is, the number of observation is less than the number parameters to be estimated. The full stoichiometry model is validated with the experimental data of Zhorov et al.5 Here, the definition of lumps considered in terms of true boiling point range for lump 1, lump 2, and lump 3 are 0-315, 315-623, and 623-873 K, respectively. The initial weight fractions of the lumps considered are 0, 0.3, and 0.7. The number of observations N considered here is 3, and the space time (in hours) for the three observation is 0.25, 0.5, and 1.0. The estimated kinetic parameters for the three reaction temperatures (in kelvin) at 673, 698, and 723 are presented in Table 3. The value of kinetic constants for the cracking reactions from the lump 3 increases with an increase in temperature from 673 to 723 K. Meanwhile, there is a sudden increase in the value of kinetic constants at 723 K, and this

6614


Table 3. Estimated Kinetic Constants for the Three Reaction Temperatures using Zhorov et al.5 Experimental Data temperature (K) kinetic constant (h-1)

673

698

723

673

full stoichiometry k1,1,1 k1,1,2 k2,1,2 k2,2,2 k1,1,3 k2,1,3 k2,2,3 k3,1,3 k3,2,3 k3,3,3 RSSE

0.001 0.001 0.001 0.001 0.126 0.577 0.577 0.586 0.586 0.586 2.5 × 10-4

0.001 0.001 0.001 0.001 0.458 1.519 1.519 1.444 1.444 1.444 5.1 × 10-3

high temperature enhances the formation of the products in the lighter lumps from the heavier lump 3. Moreover, the estimated kinetic constants for the cracking reactions from lump 2 and lump 1 at three reaction temperatures are constant as a result of availability of limited experimental data. Furthermore, the estimated kinetic constants for the reduced stoichiometry model are given in Table 3. The parity diagrams of the three lumps obtained for the full and reduced stoichiometry models are depicted in Figure 3. Here also, the residual sum of squares for

698

723

reduced stoichiometry 0.001 0.001 0.001 0.001 0.710 7.766 7.766 6.076 6.076 6.076 5 × 10-3

0.001

0.001

0.001

0.001 0.427 0.427

0.001 1.100 1.100

0.001 5.288 5.288

7.3 × 10-4

8.1 × 10-3

1.1 × 10-2

the full stoiochiometry model is low as compared with the reduced stoichiometry model. Again, the full stoichiometry model shows best fit with the experimental data of Zhorov et al.5 4. Conclusions A general full stoichiometry based discrete lumped model is presented for hydrocracking of heavier petroleum fractions using

Figure 3. Parity diagrams of the three lumps obtained for the full and reduced stoichiometry models with the exponential form of stoichiometric kernels using the experimental data of Zhorov et al.5 (a and d) Lump 1, (b and e) lump 2, and (c and f) lump 3. (a-c) Full stoichiometry model. (d-f) Reduced stoichiometry model. Thin solid line represents the trend line.


true boiling point as the basis. The reduced stoichiometry and carbon number basis are the subset of full stoichiometry model. The exponential form of stoichiometric kernels is considered for the full and reduced stoichiometry models. The kernel is symmetric for the carbon number basis model. The product distribution in a hydrocracker is governed by the ordinary differential equations. In order to reduce the computational time for parameter optimization and for better parameter estimates, the analytical solution is proposed for the system of mass balance equations. The analytical solution provides an explicit description of the model performance without the usage of a numerical method to integrate the initial value ordinary differential equations thereby reducing the computational time. Nowadays, with the advent of use of computer facilities, it is easier to solve the system of differential equations using the mathematical software packages available without any constraints. This is the main advantage of the numerical methods. However, the numerical solutions are always an approximation and this may cause errors during the estimation of derivatives in parameter optimization. Also, the convergence of the numerical method is slow as compared with the analytical one. This is the limitation of the numerical solution. The analytical solution provides better parameter estimates through regression when dealing with large number of kinetic parameters that describe the model. Therefore, a general analytical solution is proposed to determine the weight fraction of the hydrocracker products. The regression is performed using the analytical solution for full and reduced stoichiometry models with the two sets of experimental data5,37 available in the literature. The parameter estimation results reveal that the full stoichiometry model shows best fit with the experimental data as compared with the reduced stoichiometry model.

Ω(i, r - i, r) )

6i(r - i) r(r2 - 1)

6615

(A.5)

where r varies from 2 to NL and i varies from 1 to r - 1. Appendix B: Derivation of Analytical Solution In this appendix, the proof of analytical solution for the weight fraction of lump r is derived using Laplace transforms. The rearranged form of weight fraction distribution eq 15 to derive the analytical solution is given by dwr )2 dt

NL

j

j)r+1

i)1

∑ ∑δ

r,i,jΩ(r, i, j)kr,i,jwj

- βrwr + Rrwr (B.1)

where r varies from 1 to NL r

βr )

r

∑ ∑ Ω(i, j, r)k

(B.2)

i,j,r

i)1 j)1

and r

Rr ) 2

∑δ

r,j,rΩ(r, j, r)kr,j,r

(B.3)

j)1

In eq B.1, a general weight fraction expression wj for the formation of hydrocarbons in the lump r from the heavier lump j is obtained by classical integration of eq 15 for the lumps NL and NL - 1. Thus, the weight fraction wj of the lump j follow the characteristic equation NL

Acknowledgment

∑D

wj )

The fruitful comments from the anonymous reviewers are gratefully acknowledged. Appendix A: Stoichiometric Kernels The normalization and symmetric conditions used to obtain the stoichiometric kernels4 are as follows: r

dwr )2 dt

NL

r,i,jΩ(r, i, j)kr,i,j

j)r+1

i)1

(A.2)

The exponential form of the stoichiometric kernel for the full stoichiometry model is given by Ω(i, j, r) )

4ij r (r + 1)2 2

4(r - i)(r - j) r2(r - 1)2

wr(s) )

j,m

m)j

exp[(Rm - βm)t] - βrwr + Rrwr (B.5)

NL

j

∑ ∑δ

r,i,jΩ(r, i, j)kr,i,j

j)r+1

i)1

Dj,m m - Rm)

∑ s + (β m)j

s + (βr - Rr)

(B.6)

The rearranged form of Laplace domain eq B.6 is given by wr(s) )

(A.4)

where r varies from 2 to NL and i and j vary from 1 to r - 1. The symmetric form of the stoichiometric kernel for the carbon number basis model4 is given by

wr,0 + s + (βr - Rr) NL

2

(A.3)

where r varies from 1 to NL and i and j vary from 1 to r. The exponential form of the stoichiometric kernel for the reduced stoichiometry model4 is represented as Ω(i, j, r) )

∑D

Applying Laplace transform36 to eq B.5 results in the following expression

where r varies from 1 to NL. Ω(i, j, r) ) Ω(j, i, r)

(B.4)

NL

j

∑ ∑δ

(A.1)

i)1 j)1

exp[(Rm - βm)t]

Substituting eq B.4 in eq B.1, the resulting weight fraction distribution equation is given by

r

∑ ∑ Ω(i, j, r) ) 1

j,m

m)j

wr,0 s + (βr - Rr)

NL

∑

m)r+1

Dr,m + s + (βr - Rr) NL

∑

m)r+1

where

Dr,m s + (βm - Rm)

(B.7)

6616

Ind. Eng. Chem. Res., Vol. 48, No. 14, 2009 m

j

∑ ∑δ

2

r,i,jΩ(r, i, j)kr,i,jDj,m

j)r+1

Dr,m )

i)1

(B.8)

[(βr - βm) + (Rm - Rr)]

The analytical solution for the weight fraction of lump r can then be obtained by applying inverse Laplace transform in eq B.7. The resulting analytical solution is given by

coefficients D2,2 and D2,3 are substituted in eq C.3 to get the weight fraction expression for the lump 2. Furthermore, the expression for the weight fraction of lump 1 is given by w1 ) D1,1 + D1,2 exp[(R2 - β2)t] + D1,3 exp[(R3 - β3)t] (C.6) where

NL

wr(t) ) (wr,0 -

∑

D1,1 ) w1,0 - D1,2 - D1,3

Dr,m) exp[(Rr - βr)t] +

(C.7)

m)r+1 NL

∑

2

Dr,m exp[(Rm - βm)t]

(B.9)

2

m)r+1

D1,2 )

This reduces to

∑D

r,m

1,i,2Ω(1, i, 2)k1,i,2D2,2

i)1

(β1 - β2) + (R2 - R1)

(C.8)

and

NL

wr(t) )

∑δ

exp[(Rm - βm)t]

(B.10)

3

m)r

2 D1,3 )

where

j

∑ ∑δ

1,i,jΩ(1, i, j)k1,i,jDj,3

j)2 i)1

(β1 - β3) + (R3 - R1)

(C.9)

NL

Dr,r ) wr,0 -

∑

Dr,m

(B.11)

m)r+1

In this paper, it is true that a general analytical solution for the weight fraction of lumps shown in eq B.10 is identical to that of the characteristic eq B.4. In order to derive the analytical solution for the mass balance eq 15 using Laplace transforms, eq B.4 is considered in the above-mentioned derivation. The analytical solution for the system of mass balance equations for the first-order irreversible kinetics can also be derived without using Laplace transforms. However, the derivation of an exact solution is illustrated here using Laplace transforms. Appendix C: Analytical Expressions for the Three Lump Kinetic Model The weight fraction expression for lump 3 determined using eqs 16-18 is represented as w3 ) D3,3 exp[(R3 - β3)t]

(C.1)

The coefficient D3,3 is deduced from eq 20, and it becomes D3,3 ) w3,0

(C.2)

Hence, eq C.1 can be rewritten using eq C.2 to get the weight fraction expression for lump 3. Similarly, the weight fraction of lump 2 is expressed as w2 ) D2,2 exp[(R2 - β2)t] + D2,3 exp[(R3 - β3)t] (C.3) The D coefficients are calculated from eqs 20 and 21, and the corresponding expressions are as follows: D2,2 ) w2,0 - D2,3

(C.4)

3

2 D2,3 )

∑δ

2,i,3Ω(2, i, 3)k2,i,3D3,3

i)1

(β2 - β3) + (R3 - R2)

(C.5)

Substituting eq C.2 in eq C.5 results in the expression for D2,3, and it is used in C.4 to get the coefficient D2,2. Afterward, the

Literature Cited (1) Kehlen, H.; Ra¨tzsch, M. T.; Bergmann, J. Continuous kinetics for first order degradation reactions in polydisperse mixtures. Chem. Eng. Sci. 1988, 43 (3), 609–616. (2) Browarzik, D.; Kehlen, H. Hydrocracking process of n-alkanes by continuous kinetics. Chem. Eng. Sci. 1994, 49 (6), 923–926. (3) Laxminarasimhan, C. S.; Verma, R. P.; Ramachandran, P. A. Continuous lumping model for simulation of hydrocracking. AIChE J. 1996, 42 (9), 2645–2653. (4) Balasubramanian, P.; Pushpavanam, S. Model discrimination in hydrocracking of vacuum gas oil using discrete lumped kinetics. Fuel 2008, 87 (8-9), 1660–1672. (5) Zhorov, Y. M.; Panchenkov, G. M.; Tatarintseva, G. M.; Kumin, S. T.; Zen kovskii, S. M. Chemical scheme and structure of mathematical description of hydrocracking. Int. Chem. Eng. 1971, 11 (2), 256–258. (6) Stangeland, B. E.; Kittrell, J. R. Jet fuel selectivity in hydrocracking. Ind. Eng. Chem. Process Des. DeV. 1972, 11 (1), 15–20. (7) Stangeland, B. E. A kinetic model for the prediction of hydrocracker yields. Ind. Eng. Chem. Process Des. DeV. 1974, 13 (1), 71–76. (8) Ko¨seoglu, O. R.; Phillips, C. R. Kinetics and product yield distribution in the CoO-MoO3/Al2O3 catalysed hydrocracking of athabasca bitumen. Fuel 1988, 67 (10), 1411–1416. (9) Krishna, R.; Saxena, A. K. Use of an axial-dispersion model for kinetic description of hydrocracking. Chem. Eng. Sci. 1989, 44 (3), 703– 712. (10) Mohanty, S.; Saraf, D. N.; Kunzru, D. Modeling of a hydrocracking reactor. Fuel Proc. Tech. 1991, 29 (1-2), 1–17. (11) Longstaff, D. C.; Deo, M. D.; Hanson, F. V. Hydrotreatment of bitumen from the whiterocks oil sand deposit. Fuel 1994, 73 (9), 1523– 1530. (12) Ayasse, A. R.; Nagaishi, H.; Chan, E. W.; Gray, M. R. Lumped kinetics of hydrocracking of bitumen. Fuel 1997, 76 (11), 1025–1033. (13) Korre, S. C.; Klein, M. T.; Quann, R. J. Hydrocracking of polynuclear aromatic hydrocarbons. Development of rate law through inhibition studies. Ind. Eng. Chem. Res. 1997, 36 (6), 2041–2050. (14) Callejas, M. A.; Martıńez, M. T. Hydrocracking of maya residue. Kinetics and product yield distributions. Ind. Eng. Chem. Res. 1999, 38 (9), 3285–3289. (15) Pacheco, M. A.; Dassori, C. G. Hydrocracking: An improved kinetic model and reactor modeling. Chem. Eng. Commun. 2002, 189 (12), 1684– 1704. (16) Reinhardt, J.; Balfanz, U.; Dimmig, T. New kinetic method models VGO hydrocracking. Oil Gas J. 2002, 100 (18), 88. (17) Sańchez, S.; Rodrı´guez, M. A.; Ancheyta, J. Kinetic model for moderate hydrocracking of heavy oils. Ind. Eng. Chem. Res. 2005, 44 (25), 9409–9413. (18) Balasubramanian, P.; Bettina, S. J.; Pushpavanam, S.; Balaraman, K. S. Kinetic parameter estimation in hydrocracking using a combination of genetic algorithm and sequential quardratic programming. Ind. Eng. Chem. Res. 2003, 42 (20), 4723–4731.

Ind. Eng. Chem. Res., Vol. 48, No. 14, 2009 (19) Bhutani, N.; Rangaiah, G. P.; Ray, A. K. First principle data based and hybrid modeling and optimization of an industrial hydrocracking unit. Ind. Eng. Chem. Res. 2006, 45 (23), 7807–7816. (20) Sańchez, S.; Ancheyta, J. Effect of pressure on the kinetics of moderate hydrocracking of maya crude oil. Energy Fuels 2007, 21 (2), 653–661. (21) Fernandes, F. A. N.; Teles, U. M. Modeling and optimization of fischertropsch products in hydrocracking. Fuel Proc. Tech. 2007, 88, 207–214. (22) Pellegrini, L.; Bonomi, S.; Gamba, S.; Calemma, V.; Molinari, D. The all components hydrocracker model. Chem. Eng. Sci. 2007, 62 (1820), 5013–5020. (23) Pellegrini, L.; Gamba, S.; Calemma, V.; Bonomi, S. Modeling of hydrocracking with vapor-liquid equilibrium. Chem. Eng. Sci. 2008, 63 (17), 4285–4291. (24) Castano, P.; Arandes, J. M.; Pawelec, B.; Olazar, M.; Bilbao, J. Kinetic modeling for assessing the product distribution in toluene hydrocracking on a Pt/ HZSM-5 catalyst. Ind. Eng. Chem. Res. 2008, 47 (4), 1043–1050. (25) Quann, R. J.; Jaffe, S. B. Structure-oriented lumping: describing the chemistry of complex hydrocarbon mixtures. Ind. Eng. Chem. Res. 1992, 31 (11), 2483–2497. (26) Quann, R. J.; Jaffe, S. B. Building useful models of complex reaction systems in petroleum refining. Chem. Eng. Sci. 1995, 51 (10), 1615–1631. (27) Jaffe, S. B.; Freund, H.; Olmstead, W. N. Extension of structureoriented lumping to vacuum residue. Ind. Eng. Chem. Res. 2005, 44 (26), 9840–9852. (28) Baltanas, M. A.; Van Raemdonck, K. K.; Froment, G. F.; Mohedas, S. R. Fundamental kinetic modeling of hydroisomerization and hydrocracking on noble metal-loaded faujasite. 1. rate parameters for hydroisomerization. Ind. Eng. Chem. Res. 1989, 28 (7), 899–910. (29) Svoboda, G. D.; Vynckier, E.; Debrabandere, B.; Froment, G. F. Single-event rate parameters for paraffin hydrocracking on a Pt/US-Y zeolite. Ind. Eng. Chem. Res. 1995, 34 (11), 3793–3800.

6617

(30) Martens, G. G.; Marin, G. B.; Martens, J. A.; Jacobs, P. A.; Baron, G. V. A fundamental kinetic model for hydrocracking of C8 to C12 alkanes on Pt/US-Y zeolites. J. Catal. 2000, 195 (2), 253–267. (31) Martens, G. G.; Thybaut, J. W.; Marin, G. B. Single-event rate parameters for the hydrocracking of cycloalkanes on Pt/US-Y zeolites. Ind. Eng. Chem. Res. 2001, 40 (8), 1832–1844. (32) Thybaut, J. W.; Marin, G. B.; Baron, G. V.; Jacobs, P. A.; Martens, J. A. Alkene protonation enthalpy determination from fundamental kinetic modeling of alkane hydroconversion on Pt/H-(US)-Y zeolite. J. Catal. 2001, 202 (2), 324–339. (33) Martens, G. G.; Marin, G. B. Kinetics for hydrocracking based on structural classes: model development and application. AIChE J. 2001, 47 (7), 1607–1622. (34) Kumar, H.; Froment, G. F. A generalized mechanistic kinetic model for the hydroisomerization and hydrocracking of long chain paraffins. Ind. Eng. Chem. Res. 2007, 46 (12), 4075–4090. (35) Kumar, H.; Froment, G. F. Mechanistic kinetic modeling of the hydrocracking of complex feedstocks, such as vacuum gas oils. Ind. Eng. Chem. Res. 2007, 46 (18), 5881–5897. (36) Eykholt, G. R. Analytical solution for networks of irreversible firstorder reactions. Water Res. 1999, 33 (3), 814–826. (37) El-Kady, F. Y. Hydrocracking of Vacuum Distillate Fraction over Bifunctional Molybdenum-Nickel/Silica-Alumina Catalyst. Indian J. Technol. 1979, 17, 176–183.

ReceiVed for reView February 2, 2009 ReVised manuscript receiVed April 13, 2009 Accepted May 13, 2009 IE900178M

Analytical Solution for Discrete Lumped Kinetic Equations in

Recommend Documents