Comparison of Kinetic-based and Artificial Neural Network Modeling Methods for a Pilot Scale Vacuum Gas Oil Hydrocracking Reactor

An artificial neural network (ANN) and kinetic-based models for a pilot scale vacuum gas oil (VGO) hydrocracking plant are presented in this paper. Reported experimental data in the literature were used to develop, train, and check these models. The proposed models are capable of predicting the yield of all main hydrocracking products including dry gas, light naphtha, heavy naphtha, kerosene, diesel, and unconverted VGO (residue). Results showed that kinetic-based and artificial neural models have specific capabilities to predict yield of hydrocracking products. The former is able to accurately predict the yield of lighter products, i.e. light naphtha, heavy naphtha and kerosene. However, ANN model is capable of predicting yields of diesel and residue with higher precision. The comparison shows that the ANN model is superior to the kinetic-base models. © 2013 BCREC UNDIP. All rights reserved


Introduction
Petroleum refining is in a significant transition period as the industry has moved into the 21st century and the demand for petroleum products has shown a sharp growth in recent years, especially with the entry of China into the automobile market.It means that the demand for transportation fuels will show a steady growth in the next decade, contributing to petroleum product demand patterns that can only be fulfilled by the inclusion of heavier feedstock into refinery operations [1].This heavy feedstock can be converted to lighter ones using thermal and/or catalytic processing in the absence or presence of hydrogen pressure [2].Hydrocracking is one of the most versatile of all petroleumrefining processes which is the most attractive process for production of clean transportation fuels, for example, diesel with high quality [3].
Needed for all industrial processes, the optimal operation is required to guarantee the profitability which can be achieved by process models.These models are used to predict the product yields and qualities, and they are useful for sensitivity analysis, process optimization, and control, design of new plants and selection of suitable hydrocracking catalysts.However, the complexity of hydrocracking feed makes it extremely difficult to characterize and describe its kinetic at a molecular level which can be solved by considering the partition of the species into a few equivalent classes, the socalled lumps or lumping technique, and then assume each class is an independent entity.As it is mentioned in literatures, developing discretelumped kinetic models (e.g., power-law model) for complex catalytic reactions is a common approach as it can give basic information for catalyst screening, reactor design and optimization.In this field, many works were reported in literatures [4][5][6][7][8][9][10][11][12][13][14].Because performing pilot tests for vacuum gas oil (VGO) hydrocracking is expensive, covering all possible operating conditions is not economic.So, it is tried to do experiments for some points, and then find the appropriate region by using the gathered data.It can be concluded that a mathematical model can be really beneficial in this respect, but the accuracy of the developed model should be essentially investigated.This matter is crucial when the aim of performing pilot plant tests is to find the optimum point for a new catalyst which can be possibly loaded into an industrial reactor.
Besides of kinetic-based hydrocracking models which are classified as deterministic or first principal models, the use of an artificial neural network (ANN) which is a 'black box' model can be also beneficial, especially when the former approach cannot describe appropriately a system.In particular, neural networks are nonlinear and they learn (or train) by examples.The user of a neural network gathers representative data, and then invokes training algorithms to learn the structure of data [15].Applying of ANN has been previously applied successfully for modeling of the performance of industrial HDS reactors and delayed coking process which their natures are not far from the hydrocracking process [16][17][18].Furthermore, due to its ability to model the complex and nonlinear problems, the ANN can be a useful approach to model the complex behavior between input and output in the catalytic processes, such as catalyticdielectric barrier discharge plasma reactors [19][20][21].
The present study was aimed at investigating the predictability of kinetic-based and artificial neural network (ANN) models for a pilot scale VGO hydrocracker reactor.This investigation can be significant because it discusses about the usage of mathematical models to find the other operating regions from the existing data.Additionally, the application of black box models such as neural network has been always questionable especially when the number of experimental data is limited or scarce.In this paper without any prejudging about the ability of these kinds of models, it is tried to show the advantages of deterministic (kineticbased) and black box (ANN) models for a momentous refining process.
To do the experiments, the under study reactor bed was loaded with 33 g (50 cm 3 ) of a dual functional amorphous catalyst.The mesh of the catalyst particles and internal diameter of the reactor were 10-20 (0.83-1.65 mm) and about to 1 inch, respectively.

Hydrocracking Reaction Network
This work applied six lumps, i.e. unconverted VGO or residue, diesel, kerosene, light naphtha, heavy naphtha and gas to predict all valuable products of the pilot plant reactor.Figure 1 depicts the reaction pathways associated with this strategy.Note that if all reaction pathways were considered, the model including thirty kinetic parameters (frequency factors and apparent activation energies) which should be estimated using experimental data and it was laborious.Some judgments were normally welcome to reduce the model complexity, without scarifying the accuracy [22].
In Figure 1, k, F, D, K, HN, LN, and G demonstrate rate constant, VGO feed, diesel, kerosene, heavy naphtha, light naphtha and gas respectively.The combination of these nominators, show the path of hydrocracking reaction.For example, kFD represents the rate constant for conversion VGO feed to diesel.

Kinetic Expressions
Mathematical models for the VGO hydrocracking process, in a trickle-bed regime, can be very complex due to the many microscopic and macroscopic effects occurring inside the [24].So, some assumptions were introduced to simplify the model as follows: (1).Hydrocracking was a first order hydrocracking reaction [25].Since hydrogen was present in excess, the rate of hydrocracking can be supposed to be independent of the hydrogen concentration; (2).The pilot reactor operated under isothermal conditions; (3).Hydrogen feed was pure; (4).The feed and products were in the liquid phase; (5).The operation of the pilot unit was steady state; (6).Catalyst activity did not change with time.Hence simulation was only valid for the start of run.
For each reaction, a kinetic expression (R) was formulated as a function of mass concentration (C) and kinetic parameters ( k0, E ).Based on these assumptions, the kinetic constants of proposed model were as the following: Vacuum gas oil or Feed (F): (1) Note that J in Eq. (1) represents diesel (D), kerosene (K), heavy naphtha (HN), light naphtha (LN) and gas (G) lumps.
Now, the reaction rates (R) can be formulated as the following: Vacuum gas oil reaction (RF):

Mass Balance Equations
Plug flow for fixed-bed reactors was assumed in many reported models in literature, consisting of a set of ordinary differential equation (ODEs) with defined boundary conditions.In this paper, to model the hydrocracking reactor, a cell network approach was used which its accuracy was confirmed for trickle bed reactors [26].As shown in Figure 2, the hydrocracking section from the inlet to the outlet was visualized to be divided into a number of well-mixed cells (N=200 for plug flow regime) along the longitude direction.Mixing only occurred within each cell and back mixing was not accounted for the adjacent cells.The main advantage of the aforesaid approach was that the solution of ordinary differential equations was converted into the solution of a set of algebraic equations for a hydrocracking process which needed less time, and it was mathematically stable.

Bulletin of Chemical Reaction Engineering & Catalysis, 8 (2), 2013, 128
Copyright © 2013, BCREC, ISSN 1978-2993 The overall mass balance equations for all lumps are expressed as (12) where is the effectiveness factors for the hydrocracking reactants which is supposed to import the effects of the mass and pore diffusion resistance in the model.The effectiveness factor in the trickle bed regime for a spherical catalyst for an external contacting efficiency of 0.5 and Thiele modulus of 50 was calculated ranging from 0.83 to 0.8 [27].In this research, this factor for the hydrocracking reaction was considered to be equal to 0.8.Also, The " -" is for reactant (feed or VGO), and the "+" sign is for the products.
reactor; Yj is product yield and ρ is the density of stream through reactor.The catalyst volume fraction for this work was about 0.264.
The only remained unknown variable is the density of the stream inside the reactor ρ which can be calculated as follows (17) where, ρj is the density of each lump.
The reaction and mass balance expressions according to Eq. 1 to Eq. 17 were solved simultaneously to evaluate the product yields (Yi) by using Aspen Custom Modeler (ACM) programming environment (AspenTech, 2004).

Parameter Estimation
For the parameter estimation two methods were used as follows:

Un-weighted Method
In this method, the sum of squared error,SQE1, as given below, was minimized whilst all weight function (Wj) were one.(18) where Nt, Ykj meas , and Ykj pred were the number of test runs, measured yield and the predicted one, respectively.

Weighted Method
Before minimizing Eq. 18, the weight functions (W) were determined by minimizing the following expression: (19) where Wj in Eq. 19 was the weight coefficient of lumps, which played a crucial role to have an

Bulletin of Chemical Reaction Engineering & Catalysis, 8 (2), 2013, 129
Copyright © 2013, BCREC, ISSN 1978-2993 evenly distributed error among the predicted yield for the lumps with higher yield like diesel and the lumps with lower yield like light naphtha [23].The subscript ref in Eq. 19 refers to the lump with the lowest yield which was light naphtha in this work.
At first, in order to estimate weight parameters, the objective function presented in Eq. 19 was minimized by solver tool in Excel package by using Newton search method.Then Eq. 18 was minimized by applying these weights and sequencing NL2Sol and Nelder-mead algorithm, which were available in Aspen Custom Modeler software.To promote convergence from a poor initial point, a trust, the approximate region was found with NL2Sol; then to fine tune the parameters; Nelder-Mead simplex method was used.
To compare the simulated and measured product values, the mean square error were expressed in Eq. 20. (20)

Development of Artificial Neural Network Model
The artificial neural network (ANN) was developed by using the Neural Network Toolbox (newlrn function), MATLAB 2010a (The Mathworks).Because there were only 12 sets of data for each lump, a simple architecture was selected for the ANN.Therefore, a layered-recurrent neural network consisting of 2 neurons in the input layer, 4 neurons in hidden layer and five neurons in output layer was built.The inputs of the ANN were temperature and LHSV, and the outputs of that were the yield of hydrocracking products including gas, light naphtha, heavy naphtha, kerosene and diesel.To meet the mass balance criterion, residue was calculated from the mass balance equation.
quired coefficients of the designed network were limited to 53 coefficients, less than the number of training data (60 points for 5 selected lumps).Training of the ANN was carried out with the function 'trainlm' which uses Levenberg-Marquardt optimization method to estimate weights and biases.All parameters of ANN were taken as the defaults of MATLAB.Detailed of ANN algorithm used for the hydrocracking model is presented in Table 3. Training was performed until finding the minimum MSE % between predicted and measured yields of the hydrocracking products.

Results and Discussion
From sixteen sets of test runs in four levels of temperature and LHSV, twelve data were selected for the estimating (or training) of the kinetic parameters (or ANN), and four of them were preserved for evaluating the ability of the prediction.It was tried that the estimating/training and validating data could cover all levels of LHSVs and temperatures.The arrangement of these data is presented in Table 4.
At first, it was assumed that all reaction pathways presented in Figure 1 were activated by the catalyst.Therefore, thirty essential kinetic parameters (activation energy and frequency factors) of the kinetic network (called complete network) were estimated by using twelve estimating data (Table 4) which was equal to 72 observations.At this step all weights in Eq. 18 were assumed to be equal to one.Table 5 presents the estimated values of apparent activation energies and frequency factors.Additionally, rate constants for all reactions and the ratio of them to the highest one (kLNG) were calculated in the average operating temperature (about to 400 0 C).
Data in Table 5 revealed that the rate constant of reactions in the average temperature for kFHN, kFG, kDHN, kDLN, kDG, kLN, kKG and kHNG were significantly lower compared to highest value (kLNG).It is concluded that these paths can be ignored from the kinetic network in Figure 1.Therefore the main reactions involved in the hydrocracking reactor can be shown in Figure 3, called reduced kinetic net-work.The low values of kFHN and kFG and the higher value of kFD and kFK can show that the tendency of the amorphous hydrocracking catalyst to produce middle distillate from the VGO is relatively higher than other products.Moreover, from the negligible rates of kDHN, kDLN and kDG it can be found that hydrocracking of diesel to naphtha and gas was not considerable.
After estimating parameters for the complete kinetic network, to compare the simulated and measured yields, the MSE% was calculated that its average for all lumps was about 10.91%.Also, its value for each lump was calculated and presented in Table 6 under the name of complete network.
Although activation energies for the hydrocracking process are strongly related to the type of feed and catalyst, the estimated values in this work were acceptably comparable with the previous studies.As revealed by Table 5, the apparent activation energy of VGO hydrocracking to diesel and light naphtha were about to 16.6 kcal/mol and 31 kcal/mol (EFD and EFLN), respectively.The reported ones by Aboul-Ghiet [4] for hydrocracking of VGO to middle distillate and naphtha were 13-17.5 kcal/mol and 22-24 kcal/mol, respectively, not far from this research.Also, In this work, the estimated value for the apparent activation energy of VGO to kerosene (EFK) is about to 40.27 kcal/mol which was close to the reported value for hydrocracking of VGO to middle distillate in Sanchez study (39.5 kcal/mol) [28].Furthermore, the activation energy of catalytic cracking of naphtha to gas, reported by Ancheyta et al. [10] was 9-9.92 kcal/mol, close to the reported one for the light naphtha to gas in this work (ELNG).Also, the reported activation energy by the Ancheyta (ENG=9.92kcal/mol) is between the estimated values of light naphtha (ELNG=7.66)and heavy naphtha (EHNG=15.06) in this work.
In order to reduce the number of kinetic parameters involved in the model, negligible constants (Table 5) can be omitted during parameter estimation.After eliminating the least possible reactions, the reduced reaction network is depicted in Figure 3. Now, fourteen remained kinetic parameters should be re-estimated by the measured data.4. The arrangement of data used for developing and testing the models The resulted kinetic parameters for this strategy, called reduced kinetic network, are presented in Table 7.
Comparison between the simulated and measured yields showed that the MSE% of the reduced network decreased to 9.92% which was better than the complete network.Also, the MSE% of all lumps tabulated in Table 6 (under the name of reduced network) confirmed that applying this strategy was efficient to decrease the simulation error which was similar to the results of the previous work [13].
The next try for parameter estimation was done by using the factors presented in Table 8 which were estimated from minimizing of Eq. 19.After substituting these weight factors in Eq. 18, the MSE% of the simulation in comparison to the measured data was 11.06%.But, from Table 6 it can be found that this strategy was really beneficial to decrease the deviations of light lumps i.e., light naphtha, heavy naphtha which can be the strategic products of hydrocracking process [23].The resulted kinetic parameters for this model are presented in Tables 9.
Frequency Factor k0 [m 3 .h - .m  6. MSE% of hydrocracking models after estimating the parameters Finally, the described procedure for developing the artificial neural network was followed to train the model for 12 selected test runs.The fitting results showed that the MSE% of ANN model was 0.275%, considerably lower than the kinetic-based models.In Table 6, the MSE% values for all lumps were presented for training data under the name of ANN network.It can be found that the deviation of simulated values from the measured data has been appreciably decreased for all products.
After training the ANN and estimating the parameters for the kinetic models, these models were used to predict the product yields at the operating conditions selected for the prediction (Table 4).For the prediction data, the calculated deviation (MSE %) between the real and the predicted values for the gas, light naphtha, heavy naphtha, kerosene, diesel and residue in the product stream of the reactor are presented in Table 10.From this table it can be found that ANN model has lower average MSE% than the kinetic models and all kinetic models have high deviation to predict the diesel lump.However it is obvious that the weightedreduced kinetic model is really trustable to predict the yield of light naphtha, heavy naphtha and kerosene as it was mentioned before.Without any exception, all models were unable to predict the reside yield.
Among the kinetic-based models, because all of those have approximately the same deviation for diesel and residue yields, the reduced-weighted approach can be more favorable because of its high accuracy for the mentioned lumps.
The parity plots for the hydrocracking products predicted by the kinetic models and ANN models are presented in Figures 4 to 9. In Figure 4, it can be found that the ANN model appear to be fairly better than the reduced-weighted model to predict the yield of the gas product but the accuracy of its prediction is the same as other kinetic approaches.In contrary, Figures 5 to 7 which are the predicted yields for the light naphtha, heavy naphtha and kerosene confirm the accuracy of the reducedweighted model to predict the related lumps even much better than the ANN model.It is observed form Figures 8 and 9 that although the predicted yields of the diesel and residue are more satisfactory for the ANN model, the results of both ANN and kinetic-based models cannot be trustable for these lumps.
From the above discussion, it is clear that all kinetic models could not predict the yield of diesel with a high accuracy.Moreover, the accuracy of residue prediction was not acceptable for all models including ANN.
In Figure 10 the yield of diesel with temperature and LHSV has been shown, respectively.Figure 10 demonstrates that the yield of diesel decreases with ascending LHSV.This behavior can be observed at 400 o C and 410 o C in the LHSV region from 0.5 to 1.5 h -1 .But yield of diesel is different at the highest temperature (420 o C) and the LHSV upper than 1.5 h -1 , in the operating tem-    In Figure 11 the variations of residue yield versus temperature and LHSV has been shown.It can be found that the trend of residue yield at 420 o C, especially at LHSV of 2 h -1 , was different from the other operating conditions.Additionally at T=410 o C the yield of this product increased sharply form LHSV=0.5 h -1 to LHSV=1 h -1 which was not followed by the other points.It is supposed that due to limitation of test runs to only 12 points, the described uneven trends affected the predictability of models especially for the ANN model which had the lowest MSE% for the training data (Table 6).It can be concluded that the limitation of training points has caused memorizing instead of learning.
Furthermore, during the modeling process the reactor was assumed as an ideal plug one.As a rule for examining the importance of axial mixing, a useful parameter is the ratio L/dp, where L is the length of the bed and dp is the diameter of the catalyst particle.The ratio D/dp, where D is the internal diameter of the reactor, is also frequently used.Thus, a widely accepted empirical criterion is used to design experimental setup or to determine if axial mixing can be neglected and is given as the following [29]

L
According to the presented data for the understudy reactor, the length of the reactor bed was about (10 cm).Therefore the length of the bed to catalyst particle diameter was ranging from 59.8 to 118.5 which implied the possibility of the axialdispersion phenomenon thorough the catalytic bed.So, it can be concluded that the reason for some deviations in yield prediction is possibly because of axial-dispersion phenomena through the catalytic bed which was not included in the model.

Conclusions
Two different approaches i.e., kinetic-based and artificial neural network models were used to predict the product yields of a pilot scale vacuum gas oil hydrocracker.A layered-recurrent neural network consisting of 2 neurons in the input layer, 4 neurons in the hidden layer, and 5 neurons in the output layer were utilized for building the ANN model.The results showed that the MSE% of ANN model for predicting the yield of all products including gas, light naphtha, heavy naphtha, kerosene, diesel and residue was about 12.73% whilst MSE% of the most accurate kinetic-base model, called reduced strategy, was about 17.29%.Furthermore, it was confirmed that the proposed ANN model was more reliable to predict the yield of diesel whose abnormal variations in the high LHSV and temperature regions was observed.It was supposed that the reason for some deviations in yield prediction was possibly due to axial-dispersion   (♦) T=380 0 C, (■)T=400 0 C, (▲)T=410 0 C and (•)T=420 0 C [22]

Figure 1 .
Figure 1.The complete six-lump kinetic model where ε' is catalyst volume fraction; Vb is catalyst bed volume; Fm is feed mass flow rate ; Cj is lumps mass concentration; V is volume flow rate through

Figure 2 .
Figure 2. Schematic representation of series mixed cells neurons LHSV (h -1 ) and Temperature (K) Type of network Layered-recurrent Algorithm used for training Levenberg-Marquardt

Figure 4 .Figure 5 .
Figure 4. Parity plot for the predicted yields of gas product calculated by reduced-weighted and ANN models :

Figure 6 .Figure 7 .Figure 8 .Figure 9 .
Figure 6.Parity plot for the predicted yields of HN product calculated by reduced-weighted and ANN models

Table 1 .
Properties of fresh vacuum gas oil

Table 2 .
Average properties of hydrocracking products

Table 3 .
Detailed of ANN algorithm used for the hydrocracking model

Table 5 .
3cat-1] Estimated kinetic parameters for the complete kinetic network

Table 7 .
3cat-1] Estimated kinetic parameters for the reduced kinetic network

Table 8 .
Estimated factors for weighted estimation

Table 9 .
Estimated kinetic parameters for the reduced-weighted kinetic network

Table 10 .
MSE% of the prediction for the hydrocracking models