Comparison between Artiﬁcial Neural Network and Rigorous Mathematical Model in Simulation of Industrial Heavy Naphtha Reforming Process

: In this study, an artiﬁcial neural network (ANN) model was developed and compared with a rigorous mathematical model (RMM) to estimate the performance of an industrial heavy naphtha reforming process. The ANN model, represented by a multilayer feed forward neural network (MFFNN), had (36-10-10-10-34) topology, while the RMM involved solving 34 ordinary differential equations (ODEs) (32 mass balance, 1 heat balance and 1 momentum balance) to predict compo-sitions, temperature, and pressure distributions within the reforming process. All computations and predictions were performed using MATLAB ® software version 2015a. The ANN topology had minimum MSE when the number of hidden layers, number of neurons in the hidden layer, and the number of training epochs were 3, 10, and 100,000, respectively. Extensive error analysis between the experimental data and the predicted values were conducted using the following error functions: coefﬁcient of determination (R 2 ), mean absolute error (MAE), mean relative error (MRE), and mean square error (MSE). The results revealed that the ANN (R 2 = 0.9403, MAE = 0.0062) simulated the industrial heavy naphtha reforming process slightly better than the rigorous mathematical model (R 2 = 0.9318, MAE = 0.007). Moreover, the computational time was obviously reduced from 120 s for the RMM to 18.3 s for the ANN. However, one disadvantage of the ANN model is that it cannot be used to predict the process performance in the internal points of reactors, while the RMM predicted the internal temperatures, pressures and weight fractions very well. Z.M.S.; validation, A.A.A., H.S.M. and formal analysis, Z.M.S.; investigation, data


Introduction
In general, mathematical models can be categorized as deterministic or empirical; the deterministic models are constructed from first-principles equations, whereas empirical models are mathematical functions generalized to fit the data of one or more variables. Modeling catalytic refinery units is critical for designing, optimizing, and controlling tasks, but the mathematical modeling of these units poses several challenges, including the (1) assumptions used to simplify the models, (2) number of lumped components associated with the kinetic modeling, (3) complexity of the mathematical modeling, (4) catalyst activity decay with time, and (5) evaluation of physical properties under severe conditions [1,2].
Petroleum refineries consist of several thermal and catalytic units used to convert and separate petroleum fractions into useful products. Naphtha reforming units convert low-octane number heavy naphtha into a higher-octane number reformate that is the main feedstock for the blending unit to produce gasoline. Industrial naphtha reforming units include three or four reactors adiabatically operated at temperatures ranging from 450 to 520 • C and pressures ranging from 1 to 3.5 MPa [3,4]. Mathematical representations of these units can be done by integration of rigorous models that combine mass, heat, and momentum equations as well as kinetic models that govern the reactions of the hundreds of components within naphtha. In the past few years, several mathematical models have been developed to describe the behavior of naphtha reforming units. Tailleur (2012) used data obtained from commercial and micropilot plants to predict the reaction kinetics and catalyst deactivation parameters. The ordinary differential equations of mass and energy balance were solved using the Runge-Kutta-Fehlberg method in the axial and radial directions. [5] Iranshahi et al. (2014) developed mathematical and kinetic models to represent the continuous catalytic regeneration naphtha reforming process; the kinetic model consists of 32 pseudocomponents and 84 reactions. They obtained acceptable agreement between observed data and simulation results [6]. Elizalde and Ancheyta (2015) proposed the dynamic modeling of the catalytic naphtha reforming reactor using material and heat balances; the reaction network consisted of 20 components plus hydrogen and 53 chemical reactions [7]. Babaqi et al. (2018) simulated the continuous regeneration naphtha reforming process using a reaction network of 36 lumps and 55 reactions. The model has been validated by comparing its results with plant data, in which the mean relative error for the octane number, reformate yield, light gases, hydrogen yield, reactor temperature, and pressure was 1.3%, 2.5%, 0.93%, 0.43%, 1.03%, 2, and 0.6%, respectively [8]. Dong et al. (2018) described a continuous catalytic reforming process using a kinetic model of 27-lump, plug flow reactor model of a 4-zone parallel-series and an empirical catalyst deactivation model. The mean absolute prediction was found to be 0.76%, 0.42%, 0.90%, and 0.50% for paraffins, naphthenes, aromatics, and hydrogen, respectively [9]. Yusuf et al. (2019) used gPROMs ® software for the steady-state and dynamic modeling of industrial catalytic reforming. The 3D representation of 25 profiles (i.e., concentration, temperature, research octane number [RON], and hydrogen yield) with respect to time and reactor height was estimated by solving the partial differential equations governing mass and heat transfer in the process [10]. Shakor et al. (2020) estimated the kinetic parameters for the kinetic model of 32 lumps and 132 reactions by fitting the model predictions with data obtained from the industrial heavy naphtha reforming process. They observed that after 1225 days the catalyst activity decayed to 58.8% of its original activity [11]. Studies concerning the modeling of catalytic heavy naphtha reforming vary greatly in the number of pseudocomponents in the reaction mixture and the number of reactions in the kinetic networks, and hence, in the model predictions, these models are also graded from moderate complexity to very complicated. Pishnamazi et al. (2020) developed a CFD-based simulation model to predict the amount of aromatic capacity, process efficiency, transfer rate, bed temperature, and pressure in case of changes in operating conditions for naphtha reforming units [12].
Ebrahimian and Iranshahi (2020) simulated the thermal coupling of naphtha reforming with propane ammoxidation using a one-dimensional homogenous model for two processes. A genetic algorithm (GA) was applied to optimize the operating conditions of the selected configuration (naphtha-series ammoxidation). The optimum temperature and feed flow rate and the number of tubes in three reactors were selected to be 776.94 K, 2086.2 kmol h −1 , and 395, respectively. [13] Yusuf et al. (2020) estimated the plant performance, temperature, and concentration profiles of the paraffins, naphthenes, and aromatics of semicatalyst regenerative commercial naphtha catalytic reforming using gPROMS software [14].
In a petroleum refinery, the operators need a fast, simple, and accurate methodology to estimate the process predictions. An artificial neural network (ANN) is one of the promising prediction methods that can be used to predict the performance of highly nonlinear operations. ANNs are widely used for the modeling and controlling of complex chemical processes, in which an ANN has been accurately used to simulate distillation columns, [15][16][17] heat exchangers, [18,19] and catalytic reactors [20,21]. A survey of the literature found that very few studies have investigated the application of ANNs in the modeling of heavy naphtha reforming. Sadighi and Mohaddecy (2013) developed a layeredrecurrent artificial neural network to simulate an industrial fixed-bed catalytic reforming unit. They observed that the ANN could simulate the research octane number (RON), flow rate of produced gasoline, and octane barrel level with a mean absolute deviation of 0.238%, 0.813%, and 0.853%, respectively [22]. Elfghi (2016) compared the response surface methodology (RSM) and artificial neural network (ANN) for modeling catalytic naphtha reforming units and optimized the RON of a produced gasoline. The ANN methodology showed a very obvious advantage over RSM. A maximum RON of 98 was obtained at the optimal conditions (T = 521 • C, P = 3.76 MPa, LHSV = 2.02 h −1 ) [23].
The objective of this work was to explore the use of two different models to predict the performance of the heavy naphtha reforming process. The rigorous mathematical model (RMM) was selected as the deterministic model, while the artificial neural network (ANN) was selected as the empirical model.

Process Description and Data Collection
The flow sheet of the semiregenerative heavy naphtha reforming unit which is located in the Al-Doura Refinery (Baghdad, Iraq) is shown in Figure 1. This unit consists of four reactors with four interstage heaters in a series and containing a Pt/Al 2 O 3 catalyst. The reactors' inlet temperature was 470 • C, feed pressure was 2.75 MPa, and hydrogento-hydrocarbons feed ratio was 7 mol/mol. Gas chromatography/mass spectrometry (GC/MS) was used to analyze the samples collected from the feed and reactors effluent.
Catalysts 2021, 11, x FOR PEER REVIEW 3 of 17 [15][16][17] heat exchangers, [18,19] and catalytic reactors [20,21]. A survey of the literature found that very few studies have investigated the application of ANNs in the modeling of heavy naphtha reforming. Sadighi and Mohaddecy (2013) developed a layered-recurrent artificial neural network to simulate an industrial fixed-bed catalytic reforming unit. They observed that the ANN could simulate the research octane number (RON), flow rate of produced gasoline, and octane barrel level with a mean absolute deviation of 0.238%, 0.813%, and 0.853%, respectively [22]. Elfghi (2016) compared the response surface methodology (RSM) and artificial neural network (ANN) for modeling catalytic naphtha reforming units and optimized the RON of a produced gasoline. The ANN methodology showed a very obvious advantage over RSM. A maximum RON of 98 was obtained at the optimal conditions (T = 521 °C, P = 3.76 MPa, LHSV = 2.02 h −1 ) [23]. The objective of this work was to explore the use of two different models to predict the performance of the heavy naphtha reforming process. The rigorous mathematical model (RMM) was selected as the deterministic model, while the artificial neural network (ANN) was selected as the empirical model.

Process Description and Data Collection
The flow sheet of the semiregenerative heavy naphtha reforming unit which is located in the Al-Doura Refinery (Baghdad, Iraq) is shown in Figure 1. This unit consists of four reactors with four interstage heaters in a series and containing a Pt/Al2O3 catalyst. The reactors' inlet temperature was 470 °C, feed pressure was 2.75 MPa, and hydrogento-hydrocarbons feed ratio was 7 mol/mol. Gas chromatography/mass spectrometry (GC/MS) was used to analyze the samples collected from the feed and reactors effluent. Ten datasets were collected from the commercial heavy naphtha reforming process in a long time period (1225 days), starting after using a fresh catalyst in the reactors. Each dataset contains information about the weight fractions of the streams, catalyst loading, and operating conditions for the commercial heavy naphtha reforming unit (see Table 1). Ten datasets were collected from the commercial heavy naphtha reforming process in a long time period (1225 days), starting after using a fresh catalyst in the reactors. Each dataset contains information about the weight fractions of the streams, catalyst loading, and operating conditions for the commercial heavy naphtha reforming unit (see Table 1).

Rigorous Mathematical Model
To simplify the heavy naphtha reforming process, the following simplification assumptions were considered: (1) Plug flow in reactors. (2) Negligible mass and energy transfer in radial direction [24]. (3) First-order homogeneous gas phase reactions [1,2,25]. (4) Steady state adiabatic reactors. (5) Negligible heat loss to atmosphere. The physical and thermodynamic properties of the pseudocomponents was assumed to be equal to the properties of the main component for these pseudocomponents. The physical properties were obtained from reference [26].
Mass and energy changes with respect to catalyst weight are represented in ordinary differential Equations (1) and (2), respectively [27,28]: Catalysts 2021, 11, 1034 5 of 16 The reaction rate was calculated by Equation (3), and a modified Arrhenius Equation (4) was used to calculate reaction rate constant.
where: n represents the number of components, m is the reactions number. The specific heat of components was estimated according to third order polynomial [28]: Heat of reaction is calculated by: The Ergun equation was used to predict the pressure drop [29]: Equation (9) was used to convert the difference in the reactor length to the difference in the catalyst weight.
Time dependent catalyst deactivation was used to represent the rate of catalyst decay with time as shown in Equation (10) [30]: The proposed kinetic model involves 32 pseudo-components and 132 reactions; the pseudo-components are 1 to 11 carbon atoms of normal paraffins, 4 to 11 carbon atoms of isoparaffins, 6 to 11 carbon atoms of naphthenes, and 6 to 11 carbon atoms of aromatics. The kinetic model consists of isomerization, hydrocracking, dehydrogenation, dehydrocyclization, and hydrodealkylation reactions. The activation energies were grouped into nine activation energies according to reaction type. The estimated pre-exponential factors and activation energies are presented in Tables 2 and 3 respectively.

Artificial Neural Network Model
Artificial Neural networks are a series of mathematical algorithms that mimic the processes of the human brain to estimate relationships between massive amounts of data. Based on a set of experimental data, ANN models create nonlinear relationships that relate the independent and dependent variables. A multilayer feed forward neural network is the most commonly used topology in neural networks [31]. The connections between the input and hidden layer are calculated according to values called weights, which represent the strength of the connection between neurons. The output of the jth and kth nodes of the hidden layers are given by the following equations [29]: 0.0111 nP 11 + H 2 → nP 9 + P 2 0.0007 iP 10 → nP 10 0.0155 N 10 + H 2 → nP 10 0.0317 nP 5 → iP 5 0.0073 0.0008 nP 10 + H 2 → nP 7 + nP 3 0.0035 iP 9 → nP 9 0.0149 N 8 → A 8 + 3H 2 0.5233 A 11 + 4H 2 → nP 11 0.0179 nP 10 + H 2 → nP 6 + nP 4 0.0004 iP 9 + H 2 → iP 8 + P 1 0.0047 N 9 + H 2 → nP 9 0.0020 A 11 + 4H 2 → iP 11 0.0501 nP 10 + H 2 → 2nP 5 0.0052 iP 9 + H 2 → iP 7 + P 2 0.0004 N 9 + H 2 → iP 9 0.0182 A 11 + H 2 → A 10 + P 1 0.0038 nP 10 → iP 10 0.0592 iP 9 + H 2 → iP 6 + P 3 0.0008 N 9 + H 2 → N 8 + P 1 0.0026 A 11 + H 2 → A 9 + P 2 0.0435 Table 3. Activation energies of predicted kinetic model. The transfer functions (sigmoid) occur between the layers of the ANN and relate the input and output of the transfer function [32]. Three transfer functions (i.e., tangent, logarithmic, and linear) are most widely employed in neural network models [33]. As with the formation of the output of the hidden nodes, the ANN outputs are produced based on the summation of incoming weighted signals of hidden nodes passing through a specific transfer function f o .
Designing an Artificial Neural Network The ANN model was designed to predict the performance of catalytic heavy naphtha using a set of input parameters. A multilayered feed forward neural network (MFFNN) Catalysts 2021, 11, 1034 7 of 16 was used in this study. The input to the ANN model contained 36 nodes, including feed weight fractions, inlet temperature, inlet pressure, catalyst weight, and dataset time, while the output of the model contained 34 nodes, including product weight fractions, outlet temperature, and outlet pressure. There is a huge difference in the ranges of the ANN input-output data (shown in Table 4), where the weight fractions are within the range of 0-1, while the other variables are much greater. The temperatures, pressures, catalyst weights, and times were normalized within the range from 0-1 according to the Equation (14) [34]: The weight fractions of some components cannot be normalized using Equation (14) because its values may be equal to zero, making the value of the normalized variable equal to infinity as a result of division by zero.
The predicted data of the temperature and pressure were denormalized to their regional values using Equation (15): The tangent sigmoid transfer function for the hidden layer f(X) is given by Equation (16): Different topologies with one, two, and three hidden layer(s) and varying numbers of neurons (1 to 34) in hidden layer(s) were iteratively tested to achieve the optimum topology. The mean square error (MSE) was used to quantify the performance of the neural network topology according to Equation (17): The averages of the MSE values for all neural network topologies are presented in Figure 2. When using only one hidden layer, the levels of the MSE were higher than when using two or three hidden layers (Figure 2a). In addition, the unacceptable values of the MSE were estimated in different values of neurons. Using two hidden layers of an equal number of neurons produced lower values of the MSE than when using only one hidden layer, especially for neurons higher than 9 ( Figure 2b). As shown in Figure 2c, acceptable values of MSE were estimated using three hidden layers of equal numbers of neurons, especially when the number of neurons exceeded 7, in which case the values of the MSE were stable for the lowest level. Therefore, three hidden layers, each comprised of 10 neurons, were selected as the optimum ANN topology to represent the results of the heavy naphtha reforming process. number of neurons produced lower values of the MSE than when using only one hidden layer, especially for neurons higher than 9 ( Figure 2b). As shown in Figure 2c, acceptable values of MSE were estimated using three hidden layers of equal numbers of neurons, especially when the number of neurons exceeded 7, in which case the values of the MSE were stable for the lowest level. Therefore, three hidden layers, each comprised of 10 neurons, were selected as the optimum ANN topology to represent the results of the heavy naphtha reforming process.        Table 5.  Table 5.

Statistical Analyses
Four statistical criteria were estimated to evaluate the performance of the two models developed in this study: correlation coefficient (R 2 ), mean absolute error (MAE), mean relative error (MRE), and mean squared error (MSE). The coefficient of determination (R 2 ) was calculated using Equation (18) [35]: The mean absolute error (MAE) and mean relative error (MRE) were calculated using Equations (19) and (20), respectively.

ANN Model Training
The datasets were divided into three parts: 70% for training, 20% for validation, and 10% for testing, which contained 28, 8, and 4 sub-datasets, respectively. Figure 5 displays a comparison between the real data and the predicted result from the ANN model during the training mode. The MSE was 4.0519 × 10 −5 , 5.0850 × 10 −5 , and 8.9327 × 10 −6 for the weight fractions, temperatures, and pressures, respectively, while the average MSE for all data was 4.1085 × 10 −5 . The values of the MSE for temperature and pressure were lower than for the weight fractions. This was a result of the weight of the objective functions for temperature and pressure being higher than the weight of the weight fractions because the normalized variables for temperature and pressure fell within the range 0-1 while the weight fractions fell within the range from 0-0.26. The very low level of the predicted errors confirms the reliability of the proposed ANN model. fractions, temperatures, and pressures, respectively, while the average MSE for all data was 4.1085 × 10 −5 . The values of the MSE for temperature and pressure were lower than for the weight fractions. This was a result of the weight of the objective functions for temperature and pressure being higher than the weight of the weight fractions because the normalized variables for temperature and pressure fell within the range 0-1 while the weight fractions fell within the range from 0-0.26. The very low level of the predicted errors confirms the reliability of the proposed ANN model.

Comparison between the ANN and RMM Predictions
The ANN model that was validated in the learning mode was applied to estimate the model outputs in the testing mode; four sub-datasets were used in compassion between the actual data and predicted results. The catalyst deactivation term was already present in the two models, which included the ANN model as time normalized from 0-1225 to 0-1 and also added the kinetic model to achieve the predictions of the RMM model. Figures 6-8 show a comparison between the actual and predicted results using the ANN and RMM, while the summary of the predicted errors of two models is illustrated in Table 6. Figure 6a,b represents a comparison between the actual and predicted weight fractions obtained by the RMM and ANN models, respectively. The corresponding values of the R 2 , MAE, and MSE error functions achieved by the ANN model were 0.9403, 0.0062, and 0.002044, respectively, while for the RMM model, they were 0.9318, 0.0070, and 0.0002284, respectively. The MRE does not calculate for the weight fractions because some real values (those that were equal to zero) displayed an infinite relative error. Comparisons between the real and predicted outlet reactor temperatures obtained by the RMM and ANN models are presented in Figure 7a,b, respectively. The corresponding values of the R 2 , MAE, MRE, and MSE error functions achieved by the ANN model were 0.9736, 2.2246, 0.4774 and 21.6587 respectively, while for the RMM model, they were 0.8951, 4.0939, 0.5679 and 37.329, respectively. Figure 8a,b shows a comparison between the actual and predicted outlet reactor pressures obtained by the RMM and ANN models, respectively. The corresponding values of the R 2 , MAE, MRE, and MSE error functions achieved by the ANN model were 0.9467, 0.4010, 1.4129, and 1.1988, respectively, while for the RMM model, they were 0.4859, 1.6248, 6.4959, and 5.303, respectively. For the RMM model, the lower value of R 2 in the case of pressure can be attributed to the fact that the error in pressure was accumulative from the first reactor to last one, in contrast to the temperature, whereas the heat exchangers were used to heat the feed to the required temperature before entering each one of these four reactors, so there would be no accumulative error. From these figures, close mappings between the measured and simulated weight fractions, temperatures, and pressures can be observed for the two models.
respectively. The MRE does not calculate for the weight fractions because some real values (those that were equal to zero) displayed an infinite relative error. Comparisons between the real and predicted outlet reactor temperatures obtained by the RMM and ANN models are presented in Figure 7a,b, respectively. The corresponding values of the R 2 , MAE, MRE, and MSE error functions achieved by the ANN model were 0.9736, 2.2246, 0.4774 and 21.6587 respectively, while for the RMM model, they were 0.8951, 4.0939, 0.5679 and 37.329, respectively. Figure 8a,b shows a comparison between the actual and predicted outlet reactor pressures obtained by the RMM and ANN models, respectively. The corresponding values of the R 2 , MAE, MRE, and MSE error functions achieved by the ANN model were 0.9467, 0.4010, 1.4129, and 1.1988, respectively, while for the RMM model, they were 0.4859, 1.6248, 6.4959, and 5.303, respectively. For the RMM model, the lower value of R 2 in the case of pressure can be attributed to the fact that the error in pressure was accumulative from the first reactor to last one, in contrast to the temperature, whereas the heat exchangers were used to heat the feed to the required temperature before entering each one of these four reactors, so there would be no accumulative error. From these figures, close mappings between the measured and simulated weight fractions, temperatures, and pressures can be observed for the two models.    Both the ANN and RMM models displayed notable errors because predictions made using both models depended on four consecutive calculations to predict the outputs of the four reactors, which led to the accumulation error increasing from the first reactor effluent to that of the last reactor. The training data for the ANN model represented only data for input into the reactor and output from the reactor. Therefore, the ANN model is valid for the data within this range. Using catalyst weights smaller than those used in training the ANN model makes the model display unusual values of the weight fractions.  Both the ANN and RMM models displayed notable errors because predictions made using both models depended on four consecutive calculations to predict the outputs of the four reactors, which led to the accumulation error increasing from the first reactor effluent to that of the last reactor. The training data for the ANN model represented only data for input into the reactor and output from the reactor. Therefore, the ANN model is valid for the data within this range. Using catalyst weights smaller than those used in training the ANN model makes the model display unusual values of the weight fractions.
The heavy naphtha reforming process is nonlinear, and there is a strong interaction between the products' compositions, temperatures, and pressures. The dehydrogenation of paraffins and naphthenes to aromatics entails endothermic reactions that occur quickly; they take place in the first and second reactors, causing a rapid temperature drop in these reactors. The temperature drop through the first two reactors was higher than that of the last two reactors due to the exothermic nature of hydrocracking and the dehydrocyclization reaction involving the paraffins that occurred in the first two reactors. Despite this process being nonlinear, these two models predicted the four reactors' effluent weight fractions, temperatures, and pressures very well, but the ANN model's predictions were slightly better than those of the RMM model. The ANN and RMM models can effectively simulate the complicated chemical reactors, but the ANN model is faster and more accurate than the RMM model, which is in agreement with the conclusions obtained by Elçiçek et al. (2014) [32]. The computational time was 18.3 s for the ANN model and 120 s for the RMM. This great difference in the computation time is due to the fact that solving the ANN model involves substituting input variables in simple algebraic equations to estimate the outputs of the process, while solving the RMM model involves solving 34 ordinary simultaneous differential equations using the fourth-order Runge-Kutta integration method sequentially for each one of the four reactors. The validated ANN and RMM models can be used in the future for accurate simulation of industrial heavy naphtha reforming processes.

Conclusions
In the present study, real data obtained from the heavy naphtha reforming process in a long time period (1225 days) were modeled by both the rigorous mathematical model (RMM) and artificial neural network (ANN) model. The parameters of the ANN model, including the number of hidden layers, number of neurons in the hidden layers, and the number of iterations, were optimized to construct an optimal ANN model having 36-10-10-10-34 topology. To evaluate the goodness of fit, an error analysis was performed using the mean square error (MSE), coefficient of determination (R 2 ), mean relative error (MRE), and mean absolute error (MAE). The ANN model provided a precise and effective prediction of the experimental data with a coefficient of determination (R 2 ) of 0.9403, 0.9736, and 0.9467 for the weight fractions, temperatures, and pressures, respectively. In comparison, for the rigorous mathematical model, the coefficient of determination (R 2 ) was found to be 0.9318, 0.8951, and 0.4859 for the weight fractions, temperatures, and pressures, respectively. The R 2 of the ANN was higher than 0.94 for the weight fractions, temperatures, and pressures, indicating a good fit by the ANN for the testing dataset. All predictions of the error functions yielded lower values for the ANN model than the RMM model, suggesting that the ANN model is the most suitable model to describe the heavy naphtha reforming process. One disadvantage of the ANN model is that it cannot predict the process performance at the intermediate points inside reactors. We conclude that the ANN may be preferable as an alternative approach instead of the RMM to predict the performance of the heavy naphtha reforming process.