- Email: [email protected]
Solution of linear differential equations in chemical kinetics through flow graph theory approach
ARTICLE IN PRESS
JID: JTICE
[m5G;May 19, 2015;15:55]
Journal of the Taiwan Institute of Chemical Engineers 000 (2015) 1–9
Contents lists available at ScienceDirect
Journal of the Taiwan Institute of Chemical Engineers journal homepage: www.elsevier.com/locate/jtice
Solution of linear differential equations in chemical kinetics through flow graph theory approach Balasubramanian Periyasamy∗ Chemical Engineering, Faculty of Engineering, Computing and Science, Swinburne University of Technology, Sarawak Campus, 93350 Kuching, Sarawak, Malaysia.
a r t i c l e
i n f o
Article history: Received 16 November 2014 Revised 18 April 2015 Accepted 3 May 2015 Available online xxx Keywords: Cracking kinetics Flow graph theory Determinants
a b s t r a c t A flow graph theory is a method for finding the analytical solution of linear differential equations which arise in chemical kinetics through Cramer’s method of determinants. This article presents the applicability of flow graph theory for deriving the analytical solution of kinetic equations which arise in modeling of complex reaction system such as hydrocracking of heavy oils. A discrete lumped model for hydrocracking of heavy oils was developed and analytical solution for the governing model equations was derived using Laplace transforms earlier. In this work, a new method involving flow graph theory was used instead of Laplace transforms. The kinetic equations which describe the performance of a hydrocracker are governed by linear differential equations and a general analytical solution was successfully derived using flow graph theory. The analytical solution obtained through flow graph theory is similar with the reported solution using Laplace transforms for the kinetic equations of hydrocracking of heavy oils. Furthermore, the relative errors between the experimental data and model calculations using analytical solution of the three lump hydrocracker model are reasonable except for few data points. © 2015 Taiwan Institute of Chemical Engineers. Published by Elsevier B.V. All rights reserved.
1. Introduction A flow graph is a pictorial representation of reaction mechanisms involved in a chemically reacting system [1,2]. This approach determines the analytical solution of linear differential equations using Cramer’s method of determinants [1–4]. In chemical kinetics, the linear differential equations arise as a result of mass balances of first order reactions which occur in either a batch reactor or an ideal plug flow reactor [3,4]. Analytical solution for the linear differential equations can be conveniently represented as a constant multiplied with a time dependent exponential function [1–4]. In flow graph theory, the constant is determined as ratio of determinants of the formation and consumption flow graphs. The consumption flow graph can be formulated on the basis of reaction stoichiometry. The formation flow graph is deduced from consumption flow graph by including the feed source terms and neglecting the final product with zero kinetic constant of the target product. This approach eliminates the usage of classical integration, Laplace transforms, and eigenvalue method for finding analytical solution of linear differential equations which arise in chemical kinetics and engineering [2,4]. A general analytical solution for the first order monomolecular irreversible reactions which occur in a batch reactor was derived using flow graph theory by Bhusare and Balasubramanian [3]. Recently,
∗
Tel.: +6082260842; fax: +6082260813. E-mail address: [email protected]
Nurul Amira and Balasubramanian [4] applied the flow graph theory for finding analytical solution of kinetic equations for the two and three species reversible reactions. In the present work, it was decided to demonstrate the applicability of the flow graph theory approach for finding analytical solution of kinetic equations for the cracking reactions. The illustrative example considered here is a discrete lumped model for hydrocracking of heavy petroleum fractions [5,6]. In hydrocracking, the high boiling point petroleum fractions such as vacuum gas oil or vacuum residue undergo carbon–carbon bond cleavage reactions and produce low boiling point petroleum fractions such as liquefied petroleum gas, naphtha and diesel in the presence of hydrogen gas over a bifunctional solid acid catalyst [7]. The kinetic equations for hydrocracking of heavy petroleum fractions were developed by assuming binary cracking kinetics. Earlier, Krishna and Balasubramanian [6] developed a general discrete lumped model for hydrocracking of heavy oils using true boiling point of the hydrocarbons as the basis. In this model, the cracking reactions which occur within the lumps are also considered along with the cracking reactions between the lumps. A general analytical solution of kinetic equations for hydrocracking of heavy oils was derived using Laplace transforms [6]. A flow graph theory is a recently emerging method for finding analytical solution of linear differential equations which arise in chemical kinetics. This approach is not applied for finding solution of kinetic equations of cracking kinetics which arises in hydrocracking of heavy oils. Therefore, this article presents the applicability of flow
http://dx.doi.org/10.1016/j.jtice.2015.05.002 1876-1070/© 2015 Taiwan Institute of Chemical Engineers. Published by Elsevier B.V. All rights reserved.
Please cite this article as: B. Periyasamy, Solution of linear differential equations in chemical kinetics through flow graph theory approach, Journal of the Taiwan Institute of Chemical Engineers (2015), http://dx.doi.org/10.1016/j.jtice.2015.05.002
ARTICLE IN PRESS
JID: JTICE 2
[m5G;May 19, 2015;15:55]
B. Periyasamy / Journal of the Taiwan Institute of Chemical Engineers 000 (2015) 1–9
A
residue (>620 K). In kinetic modeling, it was assumed that a molecule in the heavy lump undergoes a binary cracking reaction and produces two products which may lie either in the light lump or within a reactant lump [6]. Therefore, a general stoichiometry of the cracking reactions can be represented as
C
B
ki, j,r
Lr −−→ Li + L j kBA
kCB
Fig. 1. A flow graph for the reaction A −→ B −→ C.
graph theory for finding solution of kinetic equations for hydrocracking of heavy oils.
L dwLr =2 δr,i, j r,i, j kr,i, j wL j − i, j,r ki, j,r wLr dt
j
N
2. Basis of flow graph theory A flow graph in a signal system represents a network in which nodes are connected with directed edges. Each node in a signal flow graph indicates a system variable and each edge connecting two nodes acts as a signal multiplier. The direction of signal flow is represented by placing an arrow on the edge and transmittance of signal flow is indicated along the edge. The signal flow graph depicts how signals are transmitted from one system to another and provides the relationship between the systems. In kinetics, a node represents a chemical component undergoing a structural transformation. An edge is a directed line segment joining two nodes. Weighting of an edge is a real gain between nodes and indicates the kinetic constant of a reaction. The reactant and product are represented as an input and output nodes, respectively. The input and output nodes have only outgoing and incoming edges, respectively. Furthermore, the mixed node has both outgoing and incoming edges. The concept of flow graph in chemical kinetics is illustrated with the following example. The schematic representation of flow graph for the series reaction is depicted in Fig. 1. The label A in the flow graph represents a reactant and is an input node as a result of having only outgoing edge. The labels pB and pC denote the output nodes because of having only the incoming edges. Furthermore, the labels B and C are the mixed nodes, and have both incoming and outgoing edges. The transmittance shown in Fig. 1 is a kinetic constant, and the symbols kBA and kCB denote the kinetic constants for the nodes A and B, respectively. The detailed description for the flow graph theory in kinetics can be found elsewhere [1–4].
(1)
where, r varies from NL to 1, j vary from 1 to r and i vary from 1 to j, NL is number of lumps considered, Lr is label of the lump r, and ki,j,r is kinetic constant for binary cracking of hydrocarbons in the reactant lump r into two products which lie in the product lumps i and j. The kinetic equations for hydrocracking of heavy oils [6] were developed by assuming first order irreversible cracking kinetics and are given by r
r
j=r i=1
(2)
i=1 j=1
Eq. (2) can also be applied for determining the product distribution in thermal cracking of heavy petroleum fractions [8]. In Eq. (2), r,i, j = 4i j/r2 (r + 1 )2 is exponential form of stoichiometric kernel for the distribution of products in the lumps r and i from the reactant lump j by virtue of a cracking reaction. The compensation factorδr,i, j = r/(i + r ) is included for making the sum of weight fractions added to the two lumps is equal to unity at all instances of time. Furthermore, it was assumed that kr,i,j = ki,r,j implying the symmetry of the kinetic constants included in Eq. (1). The detailed description for the derivation of Eq. (2) can be found elsewhere [6]. For deriving analytical solution, the coefficients in the first term on right-hand side for the formation of products within the lumps are conveniently grouped using the factor α r . Similarly, the coefficients in the second term on right-hand side for the disappearance of hydrocarbons in a reactant lump are grouped using the factor β r . Thus, these two factors are
αr = 2
r
δr, j,r r, j,r kr, j,r , and
(3a)
j=1
βr =
r r
i, j,r ki, j,r
(3b)
i=1 j=1
Therefore, Eq. (2) can be conveniently written as j NL dwLr =2 δr,i, j r,i, j kr,i, j wL j + (αr − βr )wLr dt
(4)
j=r+1 i=1
3. Illustrative example
The matrix-vector form of Eq. (4) can be represented as
In the following, the derivation of analytical solution for the kinetic equations of hydrocracking of heavy oils through flow graph theory approach is presented. In hydrocracking, the discrete lumps are classified on the basis of true boiling point of the hydrocarbons. The typical lumps classification for hydrocracking of heavy oils is: (i) liquefied petroleum gas (< 315 K), naphtha (315–425 K), middle distillates (425–620 K) and ⎡ α1 − β1 2 2j=1 δ1, j,2 1, j,2 k1, j,2 · · · 2 rj=1 δ1, j,r 1, j,r k1, j,r ··· ⎢ r ⎢ 0 α2 − β2 ··· 2 j=1 δ2, j,r 2, j,r k2, j,r ··· ⎢ ⎢ ⎢ . . . . . . . . . . ⎢ . . . . . ⎢ ⎢ R=⎢ 0 0 ··· αr − βr ··· ⎢ ⎢ ⎢ . . . . . ⎢ . . . . . . . . . . ⎢ ⎢ ⎢ 0 0 ··· 0 ··· ⎣ 0
0
···
0
···
w˙ L = RwL
(5)
where,
w˙ L =
dwLNL −1 dwLNL dwL1 dwL2 dwLr , ,··· ,··· , , dt dt dt dt dt
T
,
wL = (wL1 , wL2 , · · · , wLr , · · · , wLNL −1 , wLNL )T and 2
NL −1
δ1, j,NL −1 1, j,NL −1 k1, j,NL −1
2
NL
j=1
2
NL −1
δ2, j,NL −1 2, j,NL −1 k2, j,NL −1
2
NL
j=1
j=1
δ1, j,NL 1, j,NL k1, j,NL
j=1
δ2, j,NL 2, j,NL k2, j,NL
. . . 2
NL −1 j=1
. . .
δr, j,NL −1 r, j,NL −1 kr, j,NL −1
2
. . .
αNL −1 − βNL −1
NL j=1
δr, j,NL r, j,NL kr, j,NL . . .
2
0
NL j=1
δNL −1, j,NL NL −1, j,NL kNL −1, j,NL
⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦
αNL − βNL
Please cite this article as: B. Periyasamy, Solution of linear differential equations in chemical kinetics through flow graph theory approach, Journal of the Taiwan Institute of Chemical Engineers (2015), http://dx.doi.org/10.1016/j.jtice.2015.05.002
ARTICLE IN PRESS
JID: JTICE
[m5G;May 19, 2015;15:55]
B. Periyasamy / Journal of the Taiwan Institute of Chemical Engineers 000 (2015) 1–9
The general analytical solution of Eq. (5) can be written as
k1,1,2
wL = Dε
(6)
where,
1,2,2
⎡D ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ D=⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣
1.1
D1,2
···
D1,r
···
D1,NL −2
D1,NL −1
D1,NL ⎤
0
D2.2
···
D2,r
···
D2,NL −2
D2,NL −1
D2,NL ⎥
. . .
. . .
···
. . .
···
. . .
. . .
0
0
···
Dr,r
···
Dr,NL −2
Dr,NL −1
. . .
. . .
···
. . .
···
. . .
. . .
0
0
···
0
···
DNL −2.NL −2
DNL −2.NL −1
0
0
···
0
···
0
DNL −1.NL −1
0
0
···
0
···
0
0
NL
L1
2,2,2
2L2
k2,2,2
The following formula can be used for finding the superdiagonal elements of D coefficient matrix.
m
i
j=1 δr, j,i r, j,i kr, j,i Di,m βr − αr + αm − βm
i=r+1
for r < m
(8)
The elements of D coefficient matrix can be filled starting from DNL ,NL and ending withD1,1 . The factor (αr − βr ) in the exponential function of Eq. (6) can be determined using the following formula.
αr − βr = 2
r
δr, j,r r, j,r kr, j,r −
j=1
c =
L1 L2
r r
i, j,r ki, j,r
(9)
i=1 j=1
L1 1,1,1 k1,1,1 − 2δ1,1,1 1,1,1 k1,1,1 − β + α 0
fLr (αm , βm ) c(αm , βm )
k1,2,2
Fig. 2. The consumption flow graph of the two lump cracking reactions.
Fig. 2. Here, the labels p1 and p2 represent the formation of final products from L1 and L2 , respectively, with zero kinetic constants. The kinetic constants k1,1,2 and k1,2,2 represent the formation of products in the lump L1 from the reactant lump L2 . The cracking reactions which occur within the lumps L2 and L1 are denoted by the kinetic constants k2,1,2 , k2,2,2 and k1,1,1 , respectively. The stoichiometric kernels (1,1,2 , 1,2,2 , 2,1,2 and 2,2,2 ) and 1,1,1 represent the distribution of products as a result of binary cracking reactions from the reactant lumps L2 and L1 , respectively. Furthermore, the compensation factors δ i,j,2 , and δ 1,1,1 are included for making the sum of weight fractions added to the two products is equal to unity. The consumption determinant is formulated on the basis of consumption flow graph shown in Fig. 2 and is given by
−2 δ1, j,2 1, j,2 k1, j,2 j=1 2 2 2 i, j,2 ki, j,2 − 2 δ2, j,2 2, j,2 k2, j,2 − β + α i=1 j=1
(10)
k1,1,1
L1 −−→ 2L1
In Eq. 12, the first and second terms of diagonal elements represent the disappearance of hydrocarbons in the reactant lump Lr to the products in all possible ways with the positive sign in front of it and the formation of products within the reactant lump Lr with the negative sign in front of it, respectively. The first element in the second column represents the formation of products in the lump L1 from reactant lump L2 by virtue of cracking reactions with the kinetic constant k1,j,2 (j = 1, and 2). This element has a negative sign in front of it as a result of transmittance of the edge outgoing from the lump L2 and incoming to the lump L1 . Furthermore, the second element of the first column is zero as a result of irreversible cracking reaction considered in Eq. (11). The values of (α − β ) included in Eq. (12) are determined by making c = 0. That is,
c =
2 2 i=1 j=1
1, j,2 ki, j,2 − 2
2
δ2, j,2 2, j,2 k2, j,2 − β + α
j=1
×(1,1,1 k1,1,1 − 2δ1,1,1 1,1,1 k1,1,1 − β + α ) = 0 (11)
In Eq. (11), the labels L2 and L1 represent the reactant and product, respectively. The consumption flow graph for Eq. (11) is depicted in
(12)
j=1
k2,2,2 k1,1,2
p1
p2
L2 −−→ 2L2 L2 −−→ 2L1
- +
L2
where, c is determinant calculated from the consumption flow graph and is formulated on the basis of reaction stoichiometry, and fLr is determinant calculated from the formation flow graph of lump Lr and is deduced from the consumption flow graph. The concept of consumption flow graph for the cracking kinetics is illustrated by considering the simplest two lump cracking model. The stoichiometry of the cracking reactions can be represented as
L2 −−→ L1 + L2
Lump L1 - +
2
where, r varies from NL to 1. Proof In the following, the derivation of formulas for determining the elements of D coefficient matrix and the factor included in the exponential function is presented. In flow graph theory, the elements of D coefficient matrix can be determined using the following formula.
Dr,m =
1,1,1
k1,1,1
Lump L2
(7)
2L1
1,1,1
2,2,2
DNL .NL
Dr,m for r = m, and
1,1,2
1,2,2
L2
L2
m=r+1
2
k2,1,2 2,1,2
In Eq. (6), the diagonal elements of D coefficient matrix are determined by using the following formula.
Dr,r = wLr,0 −
k1,2,2
2,1,2 1,1,2
⎥ ⎥ ⎥ . ⎥ . ⎥ . ⎥ Dr,NL ⎥ ⎥ ⎥ and ⎥ . ⎥ . ⎥ . ⎥ ⎥ DNL −2.NL ⎥ ⎥ DNL −1.NL ⎦
ε = [exp(α1 − β1 )t, exp(α2 − β2 )t, · · · , exp(αr − βr )t, · · · , exp(αNL −2 − βNL −2 )t, exp(αNL −1 − βNL −1 )t, exp(αNL − βNL )t]T
Dr,m =
3
(13)
Eq. (13) is a second order polynomial in (α − β ) and the two roots are
α1 − β1 = 2δ1,1,1 1,1,1 k1,1,1 − 1,1,1 k1,1,1 , and
(14a)
Please cite this article as: B. Periyasamy, Solution of linear differential equations in chemical kinetics through flow graph theory approach, Journal of the Taiwan Institute of Chemical Engineers (2015), http://dx.doi.org/10.1016/j.jtice.2015.05.002
ARTICLE IN PRESS
JID: JTICE 4
[m5G;May 19, 2015;15:55]
B. Periyasamy / Journal of the Taiwan Institute of Chemical Engineers 000 (2015) 1–9
α2 − β2 = 2
2
δ2, j,2 2, j,2 k2, j,2 −
j=1
2 2
i, j,2 ki, j,2
(14b)
f2
Similarly, the consumption determinant for NL lump cracking kinetics is presented in Eq. (15). The first and second terms of diagonal elements represent the disappearance of reactant into products in all possible ways with positive sign in front of it and the formation of products within the reactant lump with negative sign in front of it, respectively. The superdiagonal elements denote the formation of cracked products in the light lump Lr from the heavy lump Lj (r < j) with negative sign in front of it as a result of transmittance of the edge outgoing from lump Lj and incoming to lump Lr , The subdiagonal elements are zero as a result of irreversible cracking reactions according to Eq. (1). The consumption determinant expression for NL lump cracking kinetics can be written as
k1,1,2 k1,2,2
1,2,2 2,1,2 1,1,2
k2,1,2 L2
L2
L1
2,2,2
2L2
k2,2,2
2L1
1,1,1
k1,1,1
2,2,2
Lump L2
Lump L1
p1
L1 L2
Lr
c = LNL −2
LNL −1
⎢ β1 − α1 − β + α ⎢ ⎢ ⎢ ⎢ 0 ⎢ ⎢ ⎢ . ⎢ . ⎢ . ⎢ ⎢ ⎢ 0 ⎢ ⎢ ⎢ ⎢ . . ⎢ . ⎢ ⎢ ⎢ 0 ⎢ ⎣ 0
LNL
Fig. 3. The formation flow graph of lump L2 in the two lump cracking reactions.
L2 −2
0
2
Lr
δ1, j,2 1, j,2 k1, j,2
···
−2
j=1
β2 − α2 − β + α
···
−2
i j=1 r
δ1, j,r 1, j,r k1, j,r δ2 j,r 2, j,r k2, j,r
LNL −2 ···
−2
···
−2
N L −2 j=1 N L −2
j=1
δ1, j,NL −2 1, j,NL −2 k1, j,NL −2 δ2, j,NL −2 2, j,NL −2 k2, j,NL −2
j=1
. . .
. . .
. . .
. . .
. . .
0
···
βr − αr − β + α
···
−2
N L −2
δr, j,NL −2 r, j,NL −2 kr, j,NL −2
j=1
. . .
. . .
. . .
. . .
. . .
0
···
0
···
βNL −2 − αNL −2 − β + α
0
···
0
···
0
0
···
0
···
0
LNL −1 −2 −2
N L −1 j=1 N L −1
LNL
δ1, j,NL −1 1, j,NL −1 k1, j,NL −1
−2
δ2, j,NL −1 2, j,NL −1 k2, j,NL −1
−2
NL j=1 NL
j=1
δ1, j,NL 1, j,NL k1, j,NL δ2, j,NL 2, j,NL k2, j,NL
j=1
. . . −2
N L −1
. . .
δr, j,NL −1 r, j,NL −1 kr, j,NL −1
−2
j=1
−2
N L −1
δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 βNL −1 − αNL −1 − β + α
−2 −2
NL j=1 NL j=1
0
j=1
δr, j,r r, j,r kr, j,r −
r r
i, j,r ki, j,r
δr, j,NL r, j,NL kr, j,NL . . .
j=1
r
NL j=1
. . .
αr − βr = 2
1,1,1
- +
A general formula for determining the factors (α r – β r ) can be written on the basis of Eqs. (13 and 16) and is given by L1
1,1,2
1,2,2
2,1,2
c = (β1 − α1 − β + α )(β2 − α2 − β + α ) · · · (βr − αr − β + α ) · · · (βNL −1 − αNL −1 − β + α )(βNL − αNL − β + α ) (16)
⎡
f1
wL1, 0
wL2 , 0
i=1 j=1
δNL −2, j,NL NL −2, j,NL kNL −2, j,NL δNL −1, j,NL NL −1, j,NL kNL −1, j,NL
⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦
(15)
βNL − αNL − β + α
(17)
i=1 j=1
where, r varies from NL to 1. For two lump cracking kinetics, Eq. (13) is a second order polynomial in (α − β ), and for NL lump cracking kinetics, the consumption determinant expression can be represented as nth order polynomial in (α − β ).
The expression for the consumption determinant can be conveniently represented as the product of differences of the factors (αr − βr ). Therefore, a general formula for finding the consumption determinant on the basis of consumption flow graph for the cracking stoichiometry presented in Eq. (1) is
Please cite this article as: B. Periyasamy, Solution of linear differential equations in chemical kinetics through flow graph theory approach, Journal of the Taiwan Institute of Chemical Engineers (2015), http://dx.doi.org/10.1016/j.jtice.2015.05.002
ARTICLE IN PRESS
JID: JTICE
[m5G;May 19, 2015;15:55]
B. Periyasamy / Journal of the Taiwan Institute of Chemical Engineers 000 (2015) 1–9
c(αr , βr ) =
NL
0 β j − α j + αr − βr =
Consequently, (18)
where, r = NL , NL –1, …, 1. The concept of formation flow graph is illustrated by considering again the two lump model. The formation flow graph of lump L2 is deduced from consumption flow graph shown in Fig. 2 with the consideration that the reacting lump of interest is a target one and a new feed source is added. The formation flow graph of lump L2 is depicted in Fig. 3. The labels f1 and f2 denote the feed source terms of lumps L1 and L2 with the feed weight fractions wL1,0 and wL2,0 , respectively. Here, the label p1 is final product from lump L1 with zero kinetic constant and the label p2 from lump L2 is neglected as a result of target species. The formation determinant of lump L2 is formulated on the basis of formation flow graph shown in Fig. 3. This can be represented as
fL2
NL
fNL (αr , βr ) = wLNL ,0
j=1 j=r
L1 1,1,1 k1,1,1 − 2δ1,1,1 1,1,1 k1,1,1 − β + α 0
L = 1 L2
f wL1,0 wL2,0
(βi − αi + αr − βr )
The D coefficient element DNL ,NL can be determined by substituting Eqs. (21) and (18) into Eq. (10). It results,
DNL ,NL =
fNL (αNL , βNL ) c(αNL , βNL ) wLNL ,0
=
NL
NL
(βi − αi + αNL − βNL )
i=1 i=NL
i=1 i=NL
NL
fNL −1 (αr , βr ) = wLNL−1 ,0
L1 L2
fNL =
Lr
LNL −2 LNL −1 LNL
L2
⎢ β1 − α1 − β + α ⎢ ⎢ ⎢ ⎢ 0 ⎢ ⎢ ⎢ ⎢ . . ⎢ . ⎢ ⎢ ⎢ ⎢ 0 ⎢ ⎢ ⎢ ⎢ . . ⎢ . ⎢ ⎢ ⎢ 0 ⎢ ⎣ 0 0
−2
2
δ1, j,2 1, j,2 k1, j,2
···
−2
j=1
NL
+
···
−2
r
(βi − αi + αr − βr )
(βi − αi + αr − βr )
i=1 i=NL −1,NL
×2
NL
δNL −1, j,NL NL −1, j,NL kNL −1, j,NL wLNL ,0
(23)
j=1
δ1, j,r 1, j,r k1, j,r
···
−2
j=1
β2 − α2 − β + α
(22)
i=1 i=NL −1
Lr r
= wLNL ,0
(βi − αi + αNL − βNL )
Similarly, the formation determinant of lump NL –1 is deduced from Eq. (15) by replacing the column NL –1 by using the feed weight fractions. Consequently,
The elements in the second column of Eq. (19) represent the feed weight fractions of the lumps considered. Alternatively, Eq. (19) can be deduced from Eq. (12) by replacing the second column with feed weight fractions according to Cramer’s method of determinants. Similarly, the formation determinant of lump NL is deduced from Eq. (15) by replacing column NL with feed weight fractions according to Cramer’s method of determinants. Thus, the formation determinant of lump NL is given by L1
(21)
i=1 i=r
(19)
⎡
5
δ2, j,r 2, j,r k2, j,r
NL −2 j=1 NL −2
···
−2
j=1
LNL −2
δ1, j,NL −2 1, j,NL −2 k1, j,NL −2 δ2, j,NL −2 2, j,NL −2 k2, j,NL −2
j=1
. . .
. . .
. . .
. . .
0
···
βr − αr − β + α
···
. . . −2
NL −2
δr, j,NL −2 r, j,NL −2 kr, j,NL −2
j=1
. . .
. . .
. . .
. . .
. . .
βNL −2 − αNL −2 − β + α
0
···
0
···
0
···
0
···
0
0
···
0
···
0 LNL −1
−2 −2
NL −1 j=1 NL −1
f
⎤
δ1, j,NL −1 1, j,NL −1 k1, j,NL −1
wL1,0
δ2, j,NL −1 2, j,NL −1 k2, j,NL −1
wL2,0
j=1
. . . −2
NL −1
δr, j,NL −1 r, j,NL −1 kr, j,NL −1
. . . wLr,0
j=1
. . . −2
NL −1 j=1
δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 βNL −1 − αNL −1 − β + α 0
. . . wLN
L −2,0
wLN
L −1,0
⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦
(20)
wLN
L ,0
Please cite this article as: B. Periyasamy, Solution of linear differential equations in chemical kinetics through flow graph theory approach, Journal of the Taiwan Institute of Chemical Engineers (2015), http://dx.doi.org/10.1016/j.jtice.2015.05.002
ARTICLE IN PRESS
JID: JTICE 6
[m5G;May 19, 2015;15:55]
B. Periyasamy / Journal of the Taiwan Institute of Chemical Engineers 000 (2015) 1–9
The elements of D coefficient (DNL −1,NL and DNL −1,NL −1 ) matrix of lump NL –1 are determined by substituting Eqs. (23) and (18) into Eq. (10). They are
DNL −1,NL =
fNL −1 (αNL , βNL ) c(αNL , βNL ) NL
= 0+
=
DNL −1,NL −1
2
(βi − αi + αNL − βNL )2
i=1 i=NL −1,NL
NL
NL i=1 i=NL
NL j=1
δNL −1, j,NL NL −1, j,NL kNL −1, j,NL wLNL ,0
(βi − αi + αNL − βNL )
δNL −1, j,NL NL −1, j,NL kNL −1, j,NL DNL ,NL βNL −1 − αNL −1 + αNL − βNL
j=1
fNL −1 (αNL −1 , βNL −1 ) = = c(αNL −1 , βNL −1 ) NL
= wLNL−1 ,0 −
2
NL
NL
i=1 i=NL −1
i=1 i=NL −1
NL i=1 i=NL −1
NL
(βi − αi + αNL −1 − βNL −1 )
(βi − αi + αNL −1 − βNL −1 )
(βi − αi + αNL −1 − βNL −1 )2
i=1 i=NL −1,NL
+
wLNL −1,0
(24)
NL j=1
δNL −1, j,NL NL −1, j,NL kNL −1, j,NL wLNL ,0
(βi − αi + αNL −1 − βNL −1 )
δNL −1, j,NL NL −1, j,NL kNL −1, j,NL DNL ,NL = wLNL −1,0 − DNL −1,NL βNL −1 − αNL −1 + αNL − βNL
j=1
(25)
The formation determinant of lump NL –2 is deduced from Eq. (15) by replacing the column NL –2 by using the feed weight fractions. Consequently, NL
fNL −2 (αr , βr ) = wLNL −2,0
i=1 i=NL −2
+
NL
NL
(βi − αi + αr − βr ) +
N L −1
NL
δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 wLNL−1 ,0
δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 2
j=1
i=1 i=NL −2,NL −1,NL
N L −1 j=1
i=1 i=NL −2,NL −1
(βi − αi + αr − βr ) 2
+ (βNL −1 − αNL −1 + αr − βr )2
(βi − αi + αr − βr )2
NL
δNL −1, j,NL NL −1, j,NL kNL −1, j,NL wLNL ,0
j=1
δNL −2, j,NL NL −2, j,NL kNL −2, j,NL wLNL ,0
(26)
j=1
The elements of D coefficient (DNL −2,NL ,DNL −2,NL −1 and DNL −2,NL −2 ) matrix of lump NL–2 are determined by substituting Eqs. (26) and (18) into Eq. (10). The resulting expression is given by
DNL −2,NL =
fNL −2 (αNL , βNL ) c(αNL , βNL ) NL
i=1 i=NL −2,NL −1,NL
= NL −1
L −1 L δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 2 Nj=1 δNL −1, j,NL NL −1, j,NL kNL −1, j,NL wLNL ,0 (βi − αi + αNL − βNL ) 2 Nj=1 L + (βNL −1 − αNL −1 + αNL − βNL )2 Nj=1 δNL −2, j,NL NL −2, j,NL kNL −2, j,NL wLNL ,0 NL (βi − αi + αNL − βNL ) i=1 i=NL
L δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 DNL −1,NL + 2 Nj=1 δNL −2, j,NL NL −2, j,NL kNL −2, j,NL DNL ,NL = βNL −2 − αNL −2 + αNL − βNL i NL j=1 δNL −2, j,i NL −2, j,i kNL −2, j,i Di,NL i=NL −2+1 = βNL −2 − αNL −2 + αNL − βNL
2
j=1
(27)
Please cite this article as: B. Periyasamy, Solution of linear differential equations in chemical kinetics through flow graph theory approach, Journal of the Taiwan Institute of Chemical Engineers (2015), http://dx.doi.org/10.1016/j.jtice.2015.05.002
ARTICLE IN PRESS
JID: JTICE
[m5G;May 19, 2015;15:55]
B. Periyasamy / Journal of the Taiwan Institute of Chemical Engineers 000 (2015) 1–9
fNL −2 (αNL −1 , βNL −1 ) c(αNL −1 , βNL −1 )
DNL −2,NL −1 =
NL
i=1 i=NL −2,NL −1 NL
(βi − αi + αNL −1 − βNL −1 )2
δ
=
NL −1, j,NL kNL −1, j,NL wLNL ,0 NL
i=1 i=NL −1
=
=
j=1
δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 wLNL−1 ,0
L −1 δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 (βi − αi + αNL −1 − βNL −1 ) 2 Nj=1
+
i=1 i=N −2,NL −1,NL NL L 2 j=1 NL −1, j,NL
NL −1
7
(βi − αi + αNL −1 − βNL −1 )
L δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 δNL −1, j,NL NL −1, j,NL kNL −1, j,NL DNL ,NL 2 Nj=1 wLNL−1 ,0 − βNL −2 − αNL −2 + αNL −1 − βNL −1 βNL −1 − αNL −1 + αNL − βNL
2
NL −1
2
NL −1
j=1
δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 DNL −1,NL −1 βNL −2 − αNL −2 + αNL −1 − βNL −1
j=1
(28)
fNL −2 (αNL −2 , βNL −2 ) c(αNL −2 , βNL −2 )
DNL −2,NL −2 =
NL
(βi − αi + αNL −2 − βNL −2 )wLNL−2 ,0
i=1 i=NL −2 NL
+
(βi − αi + αNL −2 − βNL −2 )wLNL−2 ,0 2
i=1 i=NL −2,NL −1,NL
NL
2
=
i=1 i=NL −2,NL −1 NL
j=1
= wLNL −2,0 − −
2
NL −1 j=1
j=1
δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 wLNL −1,0
L −1 δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 (βi − αi + αNL −2 − βNL −2 ) 2 Nj=1
δNL −1, j,NL NL −1, j,NL kNL −1, j,NL wLNL ,0 + (βNL −1 − αNL −1 + αNL −2 − βNL −2 )2 NL
i=1 i=NL −2
2
NL −1
NL j=1
δNL −2, j,NL NL −2, j,NL kNL −2, j,NL wLNL ,0
(βi − αi + αNL −2 − βNL −2 )
L 2 Nj=1 δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 δNL −1, j,NL NL −1, j,NL kNL −1, j,NL DNL ,NL wLNL−1 ,0 − βNL −2 − αNL −2 + αNL −1 − βNL −1 βNL −1 − αNL −1 + αNL − βNL
NL −1 j=1
NL L δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 DNL −1,NL + 2 Nj=1 δNL −2, j,NL NL −2, j,NL kNL −2, j,NL DNL ,NL = wLNL −2,0 − DNL −2, j βNL −2 − αNL −2 + αNL − βNL j=N −2+1 L
(29) A general formula for determining the diagonal elements of D coefficient matrix is deduced from Eq. (29) and is given by
Dr,r = wLr ,0 −
NL
Dr,m
(30)
m=r+1
Similarly, a general formula for finding the superdiagonal elements of D coefficient matrix is deduced from Eq. (27). Thus, the formula can be represented as
2
Dr,m =
m
i
δr, j,i r, j,i kr, j,i Di,m βr − αr − βm + αm
i=r+1
j=1
(31)
Eqs. (30 and 31) for finding elements of D coefficient matrix are consistent with Eqs. (7 and 8). In the following, the analytical solution of Eq. (5) through Laplace transforms [6] is briefly described. The characterization equation [6] used for deriving the general analytical solution of Eq. (5) using Laplace transforms is given by
wL j =
NL
D j,m exp[(αm − βm )t]
(32)
m= j
Eq. (32) can be obtained through classical integration of Eq. (5) for the lumps NL and NL-1 . Applying Laplace transforms to Eq. (5) after substitution of Eq. (32) results,
wLr (s ) =
NL NL wLr,0 Dr,m Dr,m − + s + (βr − αr ) s + (βr − αr ) s + (βm − αm ) m=r+1
where,
Dr,m =
2
m
(33)
m=r+1
i
j=1 δr, j,i r, j,i kr, j,i Di,m βr − αr + αm − βm
i=r+1
(34)
Please cite this article as: B. Periyasamy, Solution of linear differential equations in chemical kinetics through flow graph theory approach, Journal of the Taiwan Institute of Chemical Engineers (2015), http://dx.doi.org/10.1016/j.jtice.2015.05.002
ARTICLE IN PRESS
JID: JTICE 8
[m5G;May 19, 2015;15:55]
B. Periyasamy / Journal of the Taiwan Institute of Chemical Engineers 000 (2015) 1–9
1
Lump L1
Lump L2
Lump L3
20
Lump L2
Lump L3
15 % Relative error
Weight fraction
0.8
Lump L1
0.6
0.4 0.2
10 5 0 0
0.2
0.4
0.6
0.8
1
-5 0
0
0.25
0.5
0.75 1 1.25 Space time (h)
1.5
1.75
2
-10
Experimental weight fractions
(a)
(a) 1
Lump L1
Lump L2
Lump L3
30
0.8
Lump L2
Lump L3
20 0.6
% Relative error
Weight fraction
Lump L1
25
0.4 0.2 0
15 10 5 0 -5
0
0.25 0.5 0.75 1 1.25 1.5 1.75 Space time (h)
2
(b) Fig. 4. Performance curves of the three lump model for hydrocracking of heavy oils. (a) 663 K and (b) 723 K. Symbol: experimental data [9], and solid line: model calculations using Eqs. (6–9). Table 1 Kinetic constants of the three lump model for hydrocracking of heavy oils [6]. Kinetic constant
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
Value (h−1 ) T = 663 K
T = 723 K
0.001 0.008 0.045 0.045 0.090 0.459 0.459 0.354 0.354 0.354
0.001 0.147 0.230 0.230 0.408 2.859 2.859 2.252 2.252 2.252
The Laplace domain Eq. (33) was obtained by applying partialfractions expansion [6] and further rearrangement. The time domain Eq. (6) can be obtained by applying inverse Laplace transform to Eq. (33). The detailed description on the derivation of analytical solution for kinetic equations of heavy oils through Laplace transforms can be found elsewhere [6]. The derived general analytical solution for the kinetic equations of hydrocracking of heavy oils through flow graph theory is similar with the solution obtained using Laplace transforms. The analytical solution of the three lump kinetic model was considered in order to illustrate the virtue of this method for kinetic modeling of hydrocracking process. The time dependent behavior of the cracked products in the three lump model was calculated using Eqs. (6–9) with the literature kinetic constants [6] presented in Table 1 at two
0
0.2
0.4
0.6
0.8
1
-10 -15
Experimental weight fractions
(b) Fig. 5. Relative errors plot for the three lump model. (a) 663 K, and (b) 723 K.
reaction temperatures 663 and 723 K. The number of kinetic parameters included in the three lump model is 10. For parameter reduction purpose, the number of kinetic parameters is reduced to 6 from 10 by assuming a constant value for each group of reactions. For example, a constant value is assumed for the cracking reactions with the kinetic constants k3,1,3 , k3,2,3 and k3,3,3 . The detailed information on model parameter reduction can be found elsewhere [6]. The classification of lumps presented here was based on true boiling point of the hydrocarbons and are: 0–315, 315–523 and 523–823 K, respectively, for the lumps L1 , L2 and L3 . The space time of the six observations [9] used is 0.3333, 0.4, 0.5, 0.6667, 1, and 2. The initial weight fraction considered for the three lump model is 0, 0, and 1 for L1 , L2 and L3 , respectively. The calculated weight fractions of the three lumps at 663 and 723 K are depicted in Fig. 4a and b, respectively using Eqs. (6)–(9). Fig. 5 depicts the relative errors between the experimental [9] and model calculated weight fractions of the three lumps considered. The calculated relative errors are reasonable except for few data points of the three lump hydrocracker model presented here. 4. Concluding remarks Analytical solution of linear differential equations can be conveniently written as a sum of constants multiplied with the appropriate time-dependent exponential functions. In flow graph theory, the constants are determined as ratio of determinants of the formation and consumption flow graphs. The constants are conveniently represented in the form of D coefficient matrix. The consumption flow graph can be formulated easily on the basis of stoichiometry of the cracking kinetics. Here, the subdiagonal elements of consumption determinant are zero as a result of irreversible cracking reactions
Please cite this article as: B. Periyasamy, Solution of linear differential equations in chemical kinetics through flow graph theory approach, Journal of the Taiwan Institute of Chemical Engineers (2015), http://dx.doi.org/10.1016/j.jtice.2015.05.002
JID: JTICE
ARTICLE IN PRESS
[m5G;May 19, 2015;15:55]
B. Periyasamy / Journal of the Taiwan Institute of Chemical Engineers 000 (2015) 1–9
considered. Thus, a general expression is defined as the product of all main diagonal elements included in the consumption determinant. The formation determinant for each lump starting from NL to 1 are defined from consumption determinant by replacing the column of the target lump with feed mass fractions. Moreover, the construction of formation flow graph from consumption flow graph is not essential for defining the formation determinant of each target lump. The elements of D coefficient matrix must be filled starting from DNL ,NL and ending withD1,1 . Furthermore, it is sufficient to determine the elements starting from DNL ,NL and ending with DNL −2,NL −2 for deriving a general formula of the elements of D coefficient matrix. Earlier, a general analytical solution of kinetic equations for hydrocracking of heavy oils was derived using Laplace transforms. This involves the transformation of time domain differential equations to Laplace domain equations by using Laplace transforms of the first order derivatives. Furthermore, a characteristic equation is essential in order to carry out this transformation. Afterward, the analytical solution of kinetic equations can be obtained by using inverse Laplace transforms. This involves the alteration of Laplace domain equations to time domain expressions. The challenges involved for the application of Laplace transforms for finding the solution of simultaneous linear differential equations are: (i) partial-fractions expansion of the complex Laplace domain equations, and (ii) mathematical manipulations required for applying inverse Laplace transforms. The flow graph theory method avoids above-mentioned complications for finding the solution of linear differential equations. This approach uses a Cramer’s method of determinants for finding the solution of linear differential equations. The applicability of flow graph theory is successfully demonstrated for determining the solution of kinetic equations for hydrocracking of heavy oils in this work. The analytical solution obtained for the hydrocracker model
9
through flow graph theory is similar with the reported solution using Laplace transforms. Furthermore, the relative errors between the experimental and model calculated data using the analytical solution of the three lump model are reasonable except for few data points. In addition, the flow graph theory is not applicable for determining the solution of model equations which arise in second order and nonlinear kinetics. However, this method may be applied for deriving the kinetic solution of complex reaction systems such as thermal cracking of ethane and propane. In future, the applicability of flow graph theory will be explored for finding the solution of free radical reaction chemistry based models. References [1] Socol M, Baldea I. A new approach of flow graph theory applied in physical chemistry. J Chin Chem Soc 2006;53:773–81. [2] Socol M, Bâldea I. New method of finding the analytical solutions directly on the base on the reaction mechanism. J Math Chem 2009;45:478–87. [3] Bhusare VH, Balasubramanian P. A new paradigm in kinetic modeling of complex reaction systems. Int Rev Model Simul 2010;3(5):1145–52. [4] Nurul Amira SBH, Balasubramanian P. Exact solution for the kinetic equations of first order reversible reaction systems through flow graph theory approach. Ind Eng Chem Res 2013;52(31):10594–600. [5] Balasubramanian P, Pushpavanam S. Model discrimination in hydrocracking of vacuum gas oil using discrete lumped kinetics. Fuel 2008;87(8-9):1660–72. [6] Krishna PC, Balasubramanian P. Analytical solution for discrete lumped kinetic equations in hydrocracking of heavier petroleum fractions. Ind Eng Chem Res 2009;48(14):6608–17. [7] Ali MA, Tasumi T, Masuda T. Development of heavy oil hydrocracking catalyst using amorphous silica-alumina and zeolite as catalyst supports. Appl Catal A 2002;233:77–90. [8] Mohd Shah MSB, Periyasamy B. Binary cracking kinetics for thermal cracking of vacuum residues in Eureka process. Fuel 2014;134:618–27. [9] El-Kady FY. Hydrocracking of vacuum distillate fraction over bifunctional molybdenum-nickel/silica-alumina catalyst. Indian J Technol 1979;17:176–83.
Please cite this article as: B. Periyasamy, Solution of linear differential equations in chemical kinetics through flow graph theory approach, Journal of the Taiwan Institute of Chemical Engineers (2015), http://dx.doi.org/10.1016/j.jtice.2015.05.002
JID: JTICE
[m5G;May 19, 2015;15:55]
Journal of the Taiwan Institute of Chemical Engineers 000 (2015) 1–9
Contents lists available at ScienceDirect
Journal of the Taiwan Institute of Chemical Engineers journal homepage: www.elsevier.com/locate/jtice
Solution of linear differential equations in chemical kinetics through flow graph theory approach Balasubramanian Periyasamy∗ Chemical Engineering, Faculty of Engineering, Computing and Science, Swinburne University of Technology, Sarawak Campus, 93350 Kuching, Sarawak, Malaysia.
a r t i c l e
i n f o
Article history: Received 16 November 2014 Revised 18 April 2015 Accepted 3 May 2015 Available online xxx Keywords: Cracking kinetics Flow graph theory Determinants
a b s t r a c t A flow graph theory is a method for finding the analytical solution of linear differential equations which arise in chemical kinetics through Cramer’s method of determinants. This article presents the applicability of flow graph theory for deriving the analytical solution of kinetic equations which arise in modeling of complex reaction system such as hydrocracking of heavy oils. A discrete lumped model for hydrocracking of heavy oils was developed and analytical solution for the governing model equations was derived using Laplace transforms earlier. In this work, a new method involving flow graph theory was used instead of Laplace transforms. The kinetic equations which describe the performance of a hydrocracker are governed by linear differential equations and a general analytical solution was successfully derived using flow graph theory. The analytical solution obtained through flow graph theory is similar with the reported solution using Laplace transforms for the kinetic equations of hydrocracking of heavy oils. Furthermore, the relative errors between the experimental data and model calculations using analytical solution of the three lump hydrocracker model are reasonable except for few data points. © 2015 Taiwan Institute of Chemical Engineers. Published by Elsevier B.V. All rights reserved.
1. Introduction A flow graph is a pictorial representation of reaction mechanisms involved in a chemically reacting system [1,2]. This approach determines the analytical solution of linear differential equations using Cramer’s method of determinants [1–4]. In chemical kinetics, the linear differential equations arise as a result of mass balances of first order reactions which occur in either a batch reactor or an ideal plug flow reactor [3,4]. Analytical solution for the linear differential equations can be conveniently represented as a constant multiplied with a time dependent exponential function [1–4]. In flow graph theory, the constant is determined as ratio of determinants of the formation and consumption flow graphs. The consumption flow graph can be formulated on the basis of reaction stoichiometry. The formation flow graph is deduced from consumption flow graph by including the feed source terms and neglecting the final product with zero kinetic constant of the target product. This approach eliminates the usage of classical integration, Laplace transforms, and eigenvalue method for finding analytical solution of linear differential equations which arise in chemical kinetics and engineering [2,4]. A general analytical solution for the first order monomolecular irreversible reactions which occur in a batch reactor was derived using flow graph theory by Bhusare and Balasubramanian [3]. Recently,
∗
Tel.: +6082260842; fax: +6082260813. E-mail address: [email protected]
Nurul Amira and Balasubramanian [4] applied the flow graph theory for finding analytical solution of kinetic equations for the two and three species reversible reactions. In the present work, it was decided to demonstrate the applicability of the flow graph theory approach for finding analytical solution of kinetic equations for the cracking reactions. The illustrative example considered here is a discrete lumped model for hydrocracking of heavy petroleum fractions [5,6]. In hydrocracking, the high boiling point petroleum fractions such as vacuum gas oil or vacuum residue undergo carbon–carbon bond cleavage reactions and produce low boiling point petroleum fractions such as liquefied petroleum gas, naphtha and diesel in the presence of hydrogen gas over a bifunctional solid acid catalyst [7]. The kinetic equations for hydrocracking of heavy petroleum fractions were developed by assuming binary cracking kinetics. Earlier, Krishna and Balasubramanian [6] developed a general discrete lumped model for hydrocracking of heavy oils using true boiling point of the hydrocarbons as the basis. In this model, the cracking reactions which occur within the lumps are also considered along with the cracking reactions between the lumps. A general analytical solution of kinetic equations for hydrocracking of heavy oils was derived using Laplace transforms [6]. A flow graph theory is a recently emerging method for finding analytical solution of linear differential equations which arise in chemical kinetics. This approach is not applied for finding solution of kinetic equations of cracking kinetics which arises in hydrocracking of heavy oils. Therefore, this article presents the applicability of flow
http://dx.doi.org/10.1016/j.jtice.2015.05.002 1876-1070/© 2015 Taiwan Institute of Chemical Engineers. Published by Elsevier B.V. All rights reserved.
Please cite this article as: B. Periyasamy, Solution of linear differential equations in chemical kinetics through flow graph theory approach, Journal of the Taiwan Institute of Chemical Engineers (2015), http://dx.doi.org/10.1016/j.jtice.2015.05.002
ARTICLE IN PRESS
JID: JTICE 2
[m5G;May 19, 2015;15:55]
B. Periyasamy / Journal of the Taiwan Institute of Chemical Engineers 000 (2015) 1–9
A
residue (>620 K). In kinetic modeling, it was assumed that a molecule in the heavy lump undergoes a binary cracking reaction and produces two products which may lie either in the light lump or within a reactant lump [6]. Therefore, a general stoichiometry of the cracking reactions can be represented as
C
B
ki, j,r
Lr −−→ Li + L j kBA
kCB
Fig. 1. A flow graph for the reaction A −→ B −→ C.
graph theory for finding solution of kinetic equations for hydrocracking of heavy oils.
L dwLr =2 δr,i, j r,i, j kr,i, j wL j − i, j,r ki, j,r wLr dt
j
N
2. Basis of flow graph theory A flow graph in a signal system represents a network in which nodes are connected with directed edges. Each node in a signal flow graph indicates a system variable and each edge connecting two nodes acts as a signal multiplier. The direction of signal flow is represented by placing an arrow on the edge and transmittance of signal flow is indicated along the edge. The signal flow graph depicts how signals are transmitted from one system to another and provides the relationship between the systems. In kinetics, a node represents a chemical component undergoing a structural transformation. An edge is a directed line segment joining two nodes. Weighting of an edge is a real gain between nodes and indicates the kinetic constant of a reaction. The reactant and product are represented as an input and output nodes, respectively. The input and output nodes have only outgoing and incoming edges, respectively. Furthermore, the mixed node has both outgoing and incoming edges. The concept of flow graph in chemical kinetics is illustrated with the following example. The schematic representation of flow graph for the series reaction is depicted in Fig. 1. The label A in the flow graph represents a reactant and is an input node as a result of having only outgoing edge. The labels pB and pC denote the output nodes because of having only the incoming edges. Furthermore, the labels B and C are the mixed nodes, and have both incoming and outgoing edges. The transmittance shown in Fig. 1 is a kinetic constant, and the symbols kBA and kCB denote the kinetic constants for the nodes A and B, respectively. The detailed description for the flow graph theory in kinetics can be found elsewhere [1–4].
(1)
where, r varies from NL to 1, j vary from 1 to r and i vary from 1 to j, NL is number of lumps considered, Lr is label of the lump r, and ki,j,r is kinetic constant for binary cracking of hydrocarbons in the reactant lump r into two products which lie in the product lumps i and j. The kinetic equations for hydrocracking of heavy oils [6] were developed by assuming first order irreversible cracking kinetics and are given by r
r
j=r i=1
(2)
i=1 j=1
Eq. (2) can also be applied for determining the product distribution in thermal cracking of heavy petroleum fractions [8]. In Eq. (2), r,i, j = 4i j/r2 (r + 1 )2 is exponential form of stoichiometric kernel for the distribution of products in the lumps r and i from the reactant lump j by virtue of a cracking reaction. The compensation factorδr,i, j = r/(i + r ) is included for making the sum of weight fractions added to the two lumps is equal to unity at all instances of time. Furthermore, it was assumed that kr,i,j = ki,r,j implying the symmetry of the kinetic constants included in Eq. (1). The detailed description for the derivation of Eq. (2) can be found elsewhere [6]. For deriving analytical solution, the coefficients in the first term on right-hand side for the formation of products within the lumps are conveniently grouped using the factor α r . Similarly, the coefficients in the second term on right-hand side for the disappearance of hydrocarbons in a reactant lump are grouped using the factor β r . Thus, these two factors are
αr = 2
r
δr, j,r r, j,r kr, j,r , and
(3a)
j=1
βr =
r r
i, j,r ki, j,r
(3b)
i=1 j=1
Therefore, Eq. (2) can be conveniently written as j NL dwLr =2 δr,i, j r,i, j kr,i, j wL j + (αr − βr )wLr dt
(4)
j=r+1 i=1
3. Illustrative example
The matrix-vector form of Eq. (4) can be represented as
In the following, the derivation of analytical solution for the kinetic equations of hydrocracking of heavy oils through flow graph theory approach is presented. In hydrocracking, the discrete lumps are classified on the basis of true boiling point of the hydrocarbons. The typical lumps classification for hydrocracking of heavy oils is: (i) liquefied petroleum gas (< 315 K), naphtha (315–425 K), middle distillates (425–620 K) and ⎡ α1 − β1 2 2j=1 δ1, j,2 1, j,2 k1, j,2 · · · 2 rj=1 δ1, j,r 1, j,r k1, j,r ··· ⎢ r ⎢ 0 α2 − β2 ··· 2 j=1 δ2, j,r 2, j,r k2, j,r ··· ⎢ ⎢ ⎢ . . . . . . . . . . ⎢ . . . . . ⎢ ⎢ R=⎢ 0 0 ··· αr − βr ··· ⎢ ⎢ ⎢ . . . . . ⎢ . . . . . . . . . . ⎢ ⎢ ⎢ 0 0 ··· 0 ··· ⎣ 0
0
···
0
···
w˙ L = RwL
(5)
where,
w˙ L =
dwLNL −1 dwLNL dwL1 dwL2 dwLr , ,··· ,··· , , dt dt dt dt dt
T
,
wL = (wL1 , wL2 , · · · , wLr , · · · , wLNL −1 , wLNL )T and 2
NL −1
δ1, j,NL −1 1, j,NL −1 k1, j,NL −1
2
NL
j=1
2
NL −1
δ2, j,NL −1 2, j,NL −1 k2, j,NL −1
2
NL
j=1
j=1
δ1, j,NL 1, j,NL k1, j,NL
j=1
δ2, j,NL 2, j,NL k2, j,NL
. . . 2
NL −1 j=1
. . .
δr, j,NL −1 r, j,NL −1 kr, j,NL −1
2
. . .
αNL −1 − βNL −1
NL j=1
δr, j,NL r, j,NL kr, j,NL . . .
2
0
NL j=1
δNL −1, j,NL NL −1, j,NL kNL −1, j,NL
⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦
αNL − βNL
Please cite this article as: B. Periyasamy, Solution of linear differential equations in chemical kinetics through flow graph theory approach, Journal of the Taiwan Institute of Chemical Engineers (2015), http://dx.doi.org/10.1016/j.jtice.2015.05.002
ARTICLE IN PRESS
JID: JTICE
[m5G;May 19, 2015;15:55]
B. Periyasamy / Journal of the Taiwan Institute of Chemical Engineers 000 (2015) 1–9
The general analytical solution of Eq. (5) can be written as
k1,1,2
wL = Dε
(6)
where,
1,2,2
⎡D ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ D=⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣
1.1
D1,2
···
D1,r
···
D1,NL −2
D1,NL −1
D1,NL ⎤
0
D2.2
···
D2,r
···
D2,NL −2
D2,NL −1
D2,NL ⎥
. . .
. . .
···
. . .
···
. . .
. . .
0
0
···
Dr,r
···
Dr,NL −2
Dr,NL −1
. . .
. . .
···
. . .
···
. . .
. . .
0
0
···
0
···
DNL −2.NL −2
DNL −2.NL −1
0
0
···
0
···
0
DNL −1.NL −1
0
0
···
0
···
0
0
NL
L1
2,2,2
2L2
k2,2,2
The following formula can be used for finding the superdiagonal elements of D coefficient matrix.
m
i
j=1 δr, j,i r, j,i kr, j,i Di,m βr − αr + αm − βm
i=r+1
for r < m
(8)
The elements of D coefficient matrix can be filled starting from DNL ,NL and ending withD1,1 . The factor (αr − βr ) in the exponential function of Eq. (6) can be determined using the following formula.
αr − βr = 2
r
δr, j,r r, j,r kr, j,r −
j=1
c =
L1 L2
r r
i, j,r ki, j,r
(9)
i=1 j=1
L1 1,1,1 k1,1,1 − 2δ1,1,1 1,1,1 k1,1,1 − β + α 0
fLr (αm , βm ) c(αm , βm )
k1,2,2
Fig. 2. The consumption flow graph of the two lump cracking reactions.
Fig. 2. Here, the labels p1 and p2 represent the formation of final products from L1 and L2 , respectively, with zero kinetic constants. The kinetic constants k1,1,2 and k1,2,2 represent the formation of products in the lump L1 from the reactant lump L2 . The cracking reactions which occur within the lumps L2 and L1 are denoted by the kinetic constants k2,1,2 , k2,2,2 and k1,1,1 , respectively. The stoichiometric kernels (1,1,2 , 1,2,2 , 2,1,2 and 2,2,2 ) and 1,1,1 represent the distribution of products as a result of binary cracking reactions from the reactant lumps L2 and L1 , respectively. Furthermore, the compensation factors δ i,j,2 , and δ 1,1,1 are included for making the sum of weight fractions added to the two products is equal to unity. The consumption determinant is formulated on the basis of consumption flow graph shown in Fig. 2 and is given by
−2 δ1, j,2 1, j,2 k1, j,2 j=1 2 2 2 i, j,2 ki, j,2 − 2 δ2, j,2 2, j,2 k2, j,2 − β + α i=1 j=1
(10)
k1,1,1
L1 −−→ 2L1
In Eq. 12, the first and second terms of diagonal elements represent the disappearance of hydrocarbons in the reactant lump Lr to the products in all possible ways with the positive sign in front of it and the formation of products within the reactant lump Lr with the negative sign in front of it, respectively. The first element in the second column represents the formation of products in the lump L1 from reactant lump L2 by virtue of cracking reactions with the kinetic constant k1,j,2 (j = 1, and 2). This element has a negative sign in front of it as a result of transmittance of the edge outgoing from the lump L2 and incoming to the lump L1 . Furthermore, the second element of the first column is zero as a result of irreversible cracking reaction considered in Eq. (11). The values of (α − β ) included in Eq. (12) are determined by making c = 0. That is,
c =
2 2 i=1 j=1
1, j,2 ki, j,2 − 2
2
δ2, j,2 2, j,2 k2, j,2 − β + α
j=1
×(1,1,1 k1,1,1 − 2δ1,1,1 1,1,1 k1,1,1 − β + α ) = 0 (11)
In Eq. (11), the labels L2 and L1 represent the reactant and product, respectively. The consumption flow graph for Eq. (11) is depicted in
(12)
j=1
k2,2,2 k1,1,2
p1
p2
L2 −−→ 2L2 L2 −−→ 2L1
- +
L2
where, c is determinant calculated from the consumption flow graph and is formulated on the basis of reaction stoichiometry, and fLr is determinant calculated from the formation flow graph of lump Lr and is deduced from the consumption flow graph. The concept of consumption flow graph for the cracking kinetics is illustrated by considering the simplest two lump cracking model. The stoichiometry of the cracking reactions can be represented as
L2 −−→ L1 + L2
Lump L1 - +
2
where, r varies from NL to 1. Proof In the following, the derivation of formulas for determining the elements of D coefficient matrix and the factor included in the exponential function is presented. In flow graph theory, the elements of D coefficient matrix can be determined using the following formula.
Dr,m =
1,1,1
k1,1,1
Lump L2
(7)
2L1
1,1,1
2,2,2
DNL .NL
Dr,m for r = m, and
1,1,2
1,2,2
L2
L2
m=r+1
2
k2,1,2 2,1,2
In Eq. (6), the diagonal elements of D coefficient matrix are determined by using the following formula.
Dr,r = wLr,0 −
k1,2,2
2,1,2 1,1,2
⎥ ⎥ ⎥ . ⎥ . ⎥ . ⎥ Dr,NL ⎥ ⎥ ⎥ and ⎥ . ⎥ . ⎥ . ⎥ ⎥ DNL −2.NL ⎥ ⎥ DNL −1.NL ⎦
ε = [exp(α1 − β1 )t, exp(α2 − β2 )t, · · · , exp(αr − βr )t, · · · , exp(αNL −2 − βNL −2 )t, exp(αNL −1 − βNL −1 )t, exp(αNL − βNL )t]T
Dr,m =
3
(13)
Eq. (13) is a second order polynomial in (α − β ) and the two roots are
α1 − β1 = 2δ1,1,1 1,1,1 k1,1,1 − 1,1,1 k1,1,1 , and
(14a)
Please cite this article as: B. Periyasamy, Solution of linear differential equations in chemical kinetics through flow graph theory approach, Journal of the Taiwan Institute of Chemical Engineers (2015), http://dx.doi.org/10.1016/j.jtice.2015.05.002
ARTICLE IN PRESS
JID: JTICE 4
[m5G;May 19, 2015;15:55]
B. Periyasamy / Journal of the Taiwan Institute of Chemical Engineers 000 (2015) 1–9
α2 − β2 = 2
2
δ2, j,2 2, j,2 k2, j,2 −
j=1
2 2
i, j,2 ki, j,2
(14b)
f2
Similarly, the consumption determinant for NL lump cracking kinetics is presented in Eq. (15). The first and second terms of diagonal elements represent the disappearance of reactant into products in all possible ways with positive sign in front of it and the formation of products within the reactant lump with negative sign in front of it, respectively. The superdiagonal elements denote the formation of cracked products in the light lump Lr from the heavy lump Lj (r < j) with negative sign in front of it as a result of transmittance of the edge outgoing from lump Lj and incoming to lump Lr , The subdiagonal elements are zero as a result of irreversible cracking reactions according to Eq. (1). The consumption determinant expression for NL lump cracking kinetics can be written as
k1,1,2 k1,2,2
1,2,2 2,1,2 1,1,2
k2,1,2 L2
L2
L1
2,2,2
2L2
k2,2,2
2L1
1,1,1
k1,1,1
2,2,2
Lump L2
Lump L1
p1
L1 L2
Lr
c = LNL −2
LNL −1
⎢ β1 − α1 − β + α ⎢ ⎢ ⎢ ⎢ 0 ⎢ ⎢ ⎢ . ⎢ . ⎢ . ⎢ ⎢ ⎢ 0 ⎢ ⎢ ⎢ ⎢ . . ⎢ . ⎢ ⎢ ⎢ 0 ⎢ ⎣ 0
LNL
Fig. 3. The formation flow graph of lump L2 in the two lump cracking reactions.
L2 −2
0
2
Lr
δ1, j,2 1, j,2 k1, j,2
···
−2
j=1
β2 − α2 − β + α
···
−2
i j=1 r
δ1, j,r 1, j,r k1, j,r δ2 j,r 2, j,r k2, j,r
LNL −2 ···
−2
···
−2
N L −2 j=1 N L −2
j=1
δ1, j,NL −2 1, j,NL −2 k1, j,NL −2 δ2, j,NL −2 2, j,NL −2 k2, j,NL −2
j=1
. . .
. . .
. . .
. . .
. . .
0
···
βr − αr − β + α
···
−2
N L −2
δr, j,NL −2 r, j,NL −2 kr, j,NL −2
j=1
. . .
. . .
. . .
. . .
. . .
0
···
0
···
βNL −2 − αNL −2 − β + α
0
···
0
···
0
0
···
0
···
0
LNL −1 −2 −2
N L −1 j=1 N L −1
LNL
δ1, j,NL −1 1, j,NL −1 k1, j,NL −1
−2
δ2, j,NL −1 2, j,NL −1 k2, j,NL −1
−2
NL j=1 NL
j=1
δ1, j,NL 1, j,NL k1, j,NL δ2, j,NL 2, j,NL k2, j,NL
j=1
. . . −2
N L −1
. . .
δr, j,NL −1 r, j,NL −1 kr, j,NL −1
−2
j=1
−2
N L −1
δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 βNL −1 − αNL −1 − β + α
−2 −2
NL j=1 NL j=1
0
j=1
δr, j,r r, j,r kr, j,r −
r r
i, j,r ki, j,r
δr, j,NL r, j,NL kr, j,NL . . .
j=1
r
NL j=1
. . .
αr − βr = 2
1,1,1
- +
A general formula for determining the factors (α r – β r ) can be written on the basis of Eqs. (13 and 16) and is given by L1
1,1,2
1,2,2
2,1,2
c = (β1 − α1 − β + α )(β2 − α2 − β + α ) · · · (βr − αr − β + α ) · · · (βNL −1 − αNL −1 − β + α )(βNL − αNL − β + α ) (16)
⎡
f1
wL1, 0
wL2 , 0
i=1 j=1
δNL −2, j,NL NL −2, j,NL kNL −2, j,NL δNL −1, j,NL NL −1, j,NL kNL −1, j,NL
⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦
(15)
βNL − αNL − β + α
(17)
i=1 j=1
where, r varies from NL to 1. For two lump cracking kinetics, Eq. (13) is a second order polynomial in (α − β ), and for NL lump cracking kinetics, the consumption determinant expression can be represented as nth order polynomial in (α − β ).
The expression for the consumption determinant can be conveniently represented as the product of differences of the factors (αr − βr ). Therefore, a general formula for finding the consumption determinant on the basis of consumption flow graph for the cracking stoichiometry presented in Eq. (1) is
Please cite this article as: B. Periyasamy, Solution of linear differential equations in chemical kinetics through flow graph theory approach, Journal of the Taiwan Institute of Chemical Engineers (2015), http://dx.doi.org/10.1016/j.jtice.2015.05.002
ARTICLE IN PRESS
JID: JTICE
[m5G;May 19, 2015;15:55]
B. Periyasamy / Journal of the Taiwan Institute of Chemical Engineers 000 (2015) 1–9
c(αr , βr ) =
NL
0 β j − α j + αr − βr =
Consequently, (18)
where, r = NL , NL –1, …, 1. The concept of formation flow graph is illustrated by considering again the two lump model. The formation flow graph of lump L2 is deduced from consumption flow graph shown in Fig. 2 with the consideration that the reacting lump of interest is a target one and a new feed source is added. The formation flow graph of lump L2 is depicted in Fig. 3. The labels f1 and f2 denote the feed source terms of lumps L1 and L2 with the feed weight fractions wL1,0 and wL2,0 , respectively. Here, the label p1 is final product from lump L1 with zero kinetic constant and the label p2 from lump L2 is neglected as a result of target species. The formation determinant of lump L2 is formulated on the basis of formation flow graph shown in Fig. 3. This can be represented as
fL2
NL
fNL (αr , βr ) = wLNL ,0
j=1 j=r
L1 1,1,1 k1,1,1 − 2δ1,1,1 1,1,1 k1,1,1 − β + α 0
L = 1 L2
f wL1,0 wL2,0
(βi − αi + αr − βr )
The D coefficient element DNL ,NL can be determined by substituting Eqs. (21) and (18) into Eq. (10). It results,
DNL ,NL =
fNL (αNL , βNL ) c(αNL , βNL ) wLNL ,0
=
NL
NL
(βi − αi + αNL − βNL )
i=1 i=NL
i=1 i=NL
NL
fNL −1 (αr , βr ) = wLNL−1 ,0
L1 L2
fNL =
Lr
LNL −2 LNL −1 LNL
L2
⎢ β1 − α1 − β + α ⎢ ⎢ ⎢ ⎢ 0 ⎢ ⎢ ⎢ ⎢ . . ⎢ . ⎢ ⎢ ⎢ ⎢ 0 ⎢ ⎢ ⎢ ⎢ . . ⎢ . ⎢ ⎢ ⎢ 0 ⎢ ⎣ 0 0
−2
2
δ1, j,2 1, j,2 k1, j,2
···
−2
j=1
NL
+
···
−2
r
(βi − αi + αr − βr )
(βi − αi + αr − βr )
i=1 i=NL −1,NL
×2
NL
δNL −1, j,NL NL −1, j,NL kNL −1, j,NL wLNL ,0
(23)
j=1
δ1, j,r 1, j,r k1, j,r
···
−2
j=1
β2 − α2 − β + α
(22)
i=1 i=NL −1
Lr r
= wLNL ,0
(βi − αi + αNL − βNL )
Similarly, the formation determinant of lump NL –1 is deduced from Eq. (15) by replacing the column NL –1 by using the feed weight fractions. Consequently,
The elements in the second column of Eq. (19) represent the feed weight fractions of the lumps considered. Alternatively, Eq. (19) can be deduced from Eq. (12) by replacing the second column with feed weight fractions according to Cramer’s method of determinants. Similarly, the formation determinant of lump NL is deduced from Eq. (15) by replacing column NL with feed weight fractions according to Cramer’s method of determinants. Thus, the formation determinant of lump NL is given by L1
(21)
i=1 i=r
(19)
⎡
5
δ2, j,r 2, j,r k2, j,r
NL −2 j=1 NL −2
···
−2
j=1
LNL −2
δ1, j,NL −2 1, j,NL −2 k1, j,NL −2 δ2, j,NL −2 2, j,NL −2 k2, j,NL −2
j=1
. . .
. . .
. . .
. . .
0
···
βr − αr − β + α
···
. . . −2
NL −2
δr, j,NL −2 r, j,NL −2 kr, j,NL −2
j=1
. . .
. . .
. . .
. . .
. . .
βNL −2 − αNL −2 − β + α
0
···
0
···
0
···
0
···
0
0
···
0
···
0 LNL −1
−2 −2
NL −1 j=1 NL −1
f
⎤
δ1, j,NL −1 1, j,NL −1 k1, j,NL −1
wL1,0
δ2, j,NL −1 2, j,NL −1 k2, j,NL −1
wL2,0
j=1
. . . −2
NL −1
δr, j,NL −1 r, j,NL −1 kr, j,NL −1
. . . wLr,0
j=1
. . . −2
NL −1 j=1
δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 βNL −1 − αNL −1 − β + α 0
. . . wLN
L −2,0
wLN
L −1,0
⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦
(20)
wLN
L ,0
Please cite this article as: B. Periyasamy, Solution of linear differential equations in chemical kinetics through flow graph theory approach, Journal of the Taiwan Institute of Chemical Engineers (2015), http://dx.doi.org/10.1016/j.jtice.2015.05.002
ARTICLE IN PRESS
JID: JTICE 6
[m5G;May 19, 2015;15:55]
B. Periyasamy / Journal of the Taiwan Institute of Chemical Engineers 000 (2015) 1–9
The elements of D coefficient (DNL −1,NL and DNL −1,NL −1 ) matrix of lump NL –1 are determined by substituting Eqs. (23) and (18) into Eq. (10). They are
DNL −1,NL =
fNL −1 (αNL , βNL ) c(αNL , βNL ) NL
= 0+
=
DNL −1,NL −1
2
(βi − αi + αNL − βNL )2
i=1 i=NL −1,NL
NL
NL i=1 i=NL
NL j=1
δNL −1, j,NL NL −1, j,NL kNL −1, j,NL wLNL ,0
(βi − αi + αNL − βNL )
δNL −1, j,NL NL −1, j,NL kNL −1, j,NL DNL ,NL βNL −1 − αNL −1 + αNL − βNL
j=1
fNL −1 (αNL −1 , βNL −1 ) = = c(αNL −1 , βNL −1 ) NL
= wLNL−1 ,0 −
2
NL
NL
i=1 i=NL −1
i=1 i=NL −1
NL i=1 i=NL −1
NL
(βi − αi + αNL −1 − βNL −1 )
(βi − αi + αNL −1 − βNL −1 )
(βi − αi + αNL −1 − βNL −1 )2
i=1 i=NL −1,NL
+
wLNL −1,0
(24)
NL j=1
δNL −1, j,NL NL −1, j,NL kNL −1, j,NL wLNL ,0
(βi − αi + αNL −1 − βNL −1 )
δNL −1, j,NL NL −1, j,NL kNL −1, j,NL DNL ,NL = wLNL −1,0 − DNL −1,NL βNL −1 − αNL −1 + αNL − βNL
j=1
(25)
The formation determinant of lump NL –2 is deduced from Eq. (15) by replacing the column NL –2 by using the feed weight fractions. Consequently, NL
fNL −2 (αr , βr ) = wLNL −2,0
i=1 i=NL −2
+
NL
NL
(βi − αi + αr − βr ) +
N L −1
NL
δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 wLNL−1 ,0
δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 2
j=1
i=1 i=NL −2,NL −1,NL
N L −1 j=1
i=1 i=NL −2,NL −1
(βi − αi + αr − βr ) 2
+ (βNL −1 − αNL −1 + αr − βr )2
(βi − αi + αr − βr )2
NL
δNL −1, j,NL NL −1, j,NL kNL −1, j,NL wLNL ,0
j=1
δNL −2, j,NL NL −2, j,NL kNL −2, j,NL wLNL ,0
(26)
j=1
The elements of D coefficient (DNL −2,NL ,DNL −2,NL −1 and DNL −2,NL −2 ) matrix of lump NL–2 are determined by substituting Eqs. (26) and (18) into Eq. (10). The resulting expression is given by
DNL −2,NL =
fNL −2 (αNL , βNL ) c(αNL , βNL ) NL
i=1 i=NL −2,NL −1,NL
= NL −1
L −1 L δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 2 Nj=1 δNL −1, j,NL NL −1, j,NL kNL −1, j,NL wLNL ,0 (βi − αi + αNL − βNL ) 2 Nj=1 L + (βNL −1 − αNL −1 + αNL − βNL )2 Nj=1 δNL −2, j,NL NL −2, j,NL kNL −2, j,NL wLNL ,0 NL (βi − αi + αNL − βNL ) i=1 i=NL
L δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 DNL −1,NL + 2 Nj=1 δNL −2, j,NL NL −2, j,NL kNL −2, j,NL DNL ,NL = βNL −2 − αNL −2 + αNL − βNL i NL j=1 δNL −2, j,i NL −2, j,i kNL −2, j,i Di,NL i=NL −2+1 = βNL −2 − αNL −2 + αNL − βNL
2
j=1
(27)
Please cite this article as: B. Periyasamy, Solution of linear differential equations in chemical kinetics through flow graph theory approach, Journal of the Taiwan Institute of Chemical Engineers (2015), http://dx.doi.org/10.1016/j.jtice.2015.05.002
ARTICLE IN PRESS
JID: JTICE
[m5G;May 19, 2015;15:55]
B. Periyasamy / Journal of the Taiwan Institute of Chemical Engineers 000 (2015) 1–9
fNL −2 (αNL −1 , βNL −1 ) c(αNL −1 , βNL −1 )
DNL −2,NL −1 =
NL
i=1 i=NL −2,NL −1 NL
(βi − αi + αNL −1 − βNL −1 )2
δ
=
NL −1, j,NL kNL −1, j,NL wLNL ,0 NL
i=1 i=NL −1
=
=
j=1
δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 wLNL−1 ,0
L −1 δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 (βi − αi + αNL −1 − βNL −1 ) 2 Nj=1
+
i=1 i=N −2,NL −1,NL NL L 2 j=1 NL −1, j,NL
NL −1
7
(βi − αi + αNL −1 − βNL −1 )
L δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 δNL −1, j,NL NL −1, j,NL kNL −1, j,NL DNL ,NL 2 Nj=1 wLNL−1 ,0 − βNL −2 − αNL −2 + αNL −1 − βNL −1 βNL −1 − αNL −1 + αNL − βNL
2
NL −1
2
NL −1
j=1
δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 DNL −1,NL −1 βNL −2 − αNL −2 + αNL −1 − βNL −1
j=1
(28)
fNL −2 (αNL −2 , βNL −2 ) c(αNL −2 , βNL −2 )
DNL −2,NL −2 =
NL
(βi − αi + αNL −2 − βNL −2 )wLNL−2 ,0
i=1 i=NL −2 NL
+
(βi − αi + αNL −2 − βNL −2 )wLNL−2 ,0 2
i=1 i=NL −2,NL −1,NL
NL
2
=
i=1 i=NL −2,NL −1 NL
j=1
= wLNL −2,0 − −
2
NL −1 j=1
j=1
δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 wLNL −1,0
L −1 δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 (βi − αi + αNL −2 − βNL −2 ) 2 Nj=1
δNL −1, j,NL NL −1, j,NL kNL −1, j,NL wLNL ,0 + (βNL −1 − αNL −1 + αNL −2 − βNL −2 )2 NL
i=1 i=NL −2
2
NL −1
NL j=1
δNL −2, j,NL NL −2, j,NL kNL −2, j,NL wLNL ,0
(βi − αi + αNL −2 − βNL −2 )
L 2 Nj=1 δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 δNL −1, j,NL NL −1, j,NL kNL −1, j,NL DNL ,NL wLNL−1 ,0 − βNL −2 − αNL −2 + αNL −1 − βNL −1 βNL −1 − αNL −1 + αNL − βNL
NL −1 j=1
NL L δNL −2, j,NL −1 NL −2, j,NL −1 kNL −2, j,NL −1 DNL −1,NL + 2 Nj=1 δNL −2, j,NL NL −2, j,NL kNL −2, j,NL DNL ,NL = wLNL −2,0 − DNL −2, j βNL −2 − αNL −2 + αNL − βNL j=N −2+1 L
(29) A general formula for determining the diagonal elements of D coefficient matrix is deduced from Eq. (29) and is given by
Dr,r = wLr ,0 −
NL
Dr,m
(30)
m=r+1
Similarly, a general formula for finding the superdiagonal elements of D coefficient matrix is deduced from Eq. (27). Thus, the formula can be represented as
2
Dr,m =
m
i
δr, j,i r, j,i kr, j,i Di,m βr − αr − βm + αm
i=r+1
j=1
(31)
Eqs. (30 and 31) for finding elements of D coefficient matrix are consistent with Eqs. (7 and 8). In the following, the analytical solution of Eq. (5) through Laplace transforms [6] is briefly described. The characterization equation [6] used for deriving the general analytical solution of Eq. (5) using Laplace transforms is given by
wL j =
NL
D j,m exp[(αm − βm )t]
(32)
m= j
Eq. (32) can be obtained through classical integration of Eq. (5) for the lumps NL and NL-1 . Applying Laplace transforms to Eq. (5) after substitution of Eq. (32) results,
wLr (s ) =
NL NL wLr,0 Dr,m Dr,m − + s + (βr − αr ) s + (βr − αr ) s + (βm − αm ) m=r+1
where,
Dr,m =
2
m
(33)
m=r+1
i
j=1 δr, j,i r, j,i kr, j,i Di,m βr − αr + αm − βm
i=r+1
(34)
Please cite this article as: B. Periyasamy, Solution of linear differential equations in chemical kinetics through flow graph theory approach, Journal of the Taiwan Institute of Chemical Engineers (2015), http://dx.doi.org/10.1016/j.jtice.2015.05.002
ARTICLE IN PRESS
JID: JTICE 8
[m5G;May 19, 2015;15:55]
B. Periyasamy / Journal of the Taiwan Institute of Chemical Engineers 000 (2015) 1–9
1
Lump L1
Lump L2
Lump L3
20
Lump L2
Lump L3
15 % Relative error
Weight fraction
0.8
Lump L1
0.6
0.4 0.2
10 5 0 0
0.2
0.4
0.6
0.8
1
-5 0
0
0.25
0.5
0.75 1 1.25 Space time (h)
1.5
1.75
2
-10
Experimental weight fractions
(a)
(a) 1
Lump L1
Lump L2
Lump L3
30
0.8
Lump L2
Lump L3
20 0.6
% Relative error
Weight fraction
Lump L1
25
0.4 0.2 0
15 10 5 0 -5
0
0.25 0.5 0.75 1 1.25 1.5 1.75 Space time (h)
2
(b) Fig. 4. Performance curves of the three lump model for hydrocracking of heavy oils. (a) 663 K and (b) 723 K. Symbol: experimental data [9], and solid line: model calculations using Eqs. (6–9). Table 1 Kinetic constants of the three lump model for hydrocracking of heavy oils [6]. Kinetic constant
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
Value (h−1 ) T = 663 K
T = 723 K
0.001 0.008 0.045 0.045 0.090 0.459 0.459 0.354 0.354 0.354
0.001 0.147 0.230 0.230 0.408 2.859 2.859 2.252 2.252 2.252
The Laplace domain Eq. (33) was obtained by applying partialfractions expansion [6] and further rearrangement. The time domain Eq. (6) can be obtained by applying inverse Laplace transform to Eq. (33). The detailed description on the derivation of analytical solution for kinetic equations of heavy oils through Laplace transforms can be found elsewhere [6]. The derived general analytical solution for the kinetic equations of hydrocracking of heavy oils through flow graph theory is similar with the solution obtained using Laplace transforms. The analytical solution of the three lump kinetic model was considered in order to illustrate the virtue of this method for kinetic modeling of hydrocracking process. The time dependent behavior of the cracked products in the three lump model was calculated using Eqs. (6–9) with the literature kinetic constants [6] presented in Table 1 at two
0
0.2
0.4
0.6
0.8
1
-10 -15
Experimental weight fractions
(b) Fig. 5. Relative errors plot for the three lump model. (a) 663 K, and (b) 723 K.
reaction temperatures 663 and 723 K. The number of kinetic parameters included in the three lump model is 10. For parameter reduction purpose, the number of kinetic parameters is reduced to 6 from 10 by assuming a constant value for each group of reactions. For example, a constant value is assumed for the cracking reactions with the kinetic constants k3,1,3 , k3,2,3 and k3,3,3 . The detailed information on model parameter reduction can be found elsewhere [6]. The classification of lumps presented here was based on true boiling point of the hydrocarbons and are: 0–315, 315–523 and 523–823 K, respectively, for the lumps L1 , L2 and L3 . The space time of the six observations [9] used is 0.3333, 0.4, 0.5, 0.6667, 1, and 2. The initial weight fraction considered for the three lump model is 0, 0, and 1 for L1 , L2 and L3 , respectively. The calculated weight fractions of the three lumps at 663 and 723 K are depicted in Fig. 4a and b, respectively using Eqs. (6)–(9). Fig. 5 depicts the relative errors between the experimental [9] and model calculated weight fractions of the three lumps considered. The calculated relative errors are reasonable except for few data points of the three lump hydrocracker model presented here. 4. Concluding remarks Analytical solution of linear differential equations can be conveniently written as a sum of constants multiplied with the appropriate time-dependent exponential functions. In flow graph theory, the constants are determined as ratio of determinants of the formation and consumption flow graphs. The constants are conveniently represented in the form of D coefficient matrix. The consumption flow graph can be formulated easily on the basis of stoichiometry of the cracking kinetics. Here, the subdiagonal elements of consumption determinant are zero as a result of irreversible cracking reactions
Please cite this article as: B. Periyasamy, Solution of linear differential equations in chemical kinetics through flow graph theory approach, Journal of the Taiwan Institute of Chemical Engineers (2015), http://dx.doi.org/10.1016/j.jtice.2015.05.002
JID: JTICE
ARTICLE IN PRESS
[m5G;May 19, 2015;15:55]
B. Periyasamy / Journal of the Taiwan Institute of Chemical Engineers 000 (2015) 1–9
considered. Thus, a general expression is defined as the product of all main diagonal elements included in the consumption determinant. The formation determinant for each lump starting from NL to 1 are defined from consumption determinant by replacing the column of the target lump with feed mass fractions. Moreover, the construction of formation flow graph from consumption flow graph is not essential for defining the formation determinant of each target lump. The elements of D coefficient matrix must be filled starting from DNL ,NL and ending withD1,1 . Furthermore, it is sufficient to determine the elements starting from DNL ,NL and ending with DNL −2,NL −2 for deriving a general formula of the elements of D coefficient matrix. Earlier, a general analytical solution of kinetic equations for hydrocracking of heavy oils was derived using Laplace transforms. This involves the transformation of time domain differential equations to Laplace domain equations by using Laplace transforms of the first order derivatives. Furthermore, a characteristic equation is essential in order to carry out this transformation. Afterward, the analytical solution of kinetic equations can be obtained by using inverse Laplace transforms. This involves the alteration of Laplace domain equations to time domain expressions. The challenges involved for the application of Laplace transforms for finding the solution of simultaneous linear differential equations are: (i) partial-fractions expansion of the complex Laplace domain equations, and (ii) mathematical manipulations required for applying inverse Laplace transforms. The flow graph theory method avoids above-mentioned complications for finding the solution of linear differential equations. This approach uses a Cramer’s method of determinants for finding the solution of linear differential equations. The applicability of flow graph theory is successfully demonstrated for determining the solution of kinetic equations for hydrocracking of heavy oils in this work. The analytical solution obtained for the hydrocracker model
9
through flow graph theory is similar with the reported solution using Laplace transforms. Furthermore, the relative errors between the experimental and model calculated data using the analytical solution of the three lump model are reasonable except for few data points. In addition, the flow graph theory is not applicable for determining the solution of model equations which arise in second order and nonlinear kinetics. However, this method may be applied for deriving the kinetic solution of complex reaction systems such as thermal cracking of ethane and propane. In future, the applicability of flow graph theory will be explored for finding the solution of free radical reaction chemistry based models. References [1] Socol M, Baldea I. A new approach of flow graph theory applied in physical chemistry. J Chin Chem Soc 2006;53:773–81. [2] Socol M, Bâldea I. New method of finding the analytical solutions directly on the base on the reaction mechanism. J Math Chem 2009;45:478–87. [3] Bhusare VH, Balasubramanian P. A new paradigm in kinetic modeling of complex reaction systems. Int Rev Model Simul 2010;3(5):1145–52. [4] Nurul Amira SBH, Balasubramanian P. Exact solution for the kinetic equations of first order reversible reaction systems through flow graph theory approach. Ind Eng Chem Res 2013;52(31):10594–600. [5] Balasubramanian P, Pushpavanam S. Model discrimination in hydrocracking of vacuum gas oil using discrete lumped kinetics. Fuel 2008;87(8-9):1660–72. [6] Krishna PC, Balasubramanian P. Analytical solution for discrete lumped kinetic equations in hydrocracking of heavier petroleum fractions. Ind Eng Chem Res 2009;48(14):6608–17. [7] Ali MA, Tasumi T, Masuda T. Development of heavy oil hydrocracking catalyst using amorphous silica-alumina and zeolite as catalyst supports. Appl Catal A 2002;233:77–90. [8] Mohd Shah MSB, Periyasamy B. Binary cracking kinetics for thermal cracking of vacuum residues in Eureka process. Fuel 2014;134:618–27. [9] El-Kady FY. Hydrocracking of vacuum distillate fraction over bifunctional molybdenum-nickel/silica-alumina catalyst. Indian J Technol 1979;17:176–83.
Please cite this article as: B. Periyasamy, Solution of linear differential equations in chemical kinetics through flow graph theory approach, Journal of the Taiwan Institute of Chemical Engineers (2015), http://dx.doi.org/10.1016/j.jtice.2015.05.002