Study of Industrial Naphtha Catalytic Reforming Reactions via Modelling and Simulation

Abstract

Steady state and dynamic modelling and simulation of catalytic reforming unit of Kaduna Refining & Petrochemical Company, NNPC (Nigeria) was carried to find out the behaviour of the reactions under both steady and unsteady state conditions. The basic model together with kinetic and thermodynamic parameters and properties were taken from the literature but is developed in gPROMs (an equation oriented modelling software) model building platform for the first time rather than in MATLAB or other modelling platform used by other researchers in the past. The simulation was performed using gPROMs and the predictions were validated against those available in the literature. The validated model was then used to monitor the behaviour of the temperature, concentrations of paraffins, naphthenes and aromatics with respect to both time and height of the reactor of the industrial refinery of Nigeria. Hydrogen yield, Research octane number (RON) and temperature profiles are also reported. The components behave similarly in terms of reactions in the reactors but the time to attain quasi-steady state is different. The results are in good agreement with the industrial plant data.

1. Introduction

Naphtha is a low octane cut from atmospheric distillation of crude oil from the topping plant usually in the range between 30 °C–200 °C and constitutes typically 15–30% by weight of the crude oil [1]. This cut contains some impurities such as sulphur, oxygen and nitrogen compounds. This low run naphtha, sometimes referred to as straight run, requires reforming to a more valuable higher run naphtha called reformate which is blended into a pool of gasoline for use as a motor spirit. SR naphtha can be Light naphtha a fraction with boiling point ranging from 30 °C to 90 °C, containing the C₅ and C₆ hydrocarbons or Heavy naphtha which is the fraction with boiling range from 90 °C to 200 °C. The medium naphtha is a fraction of this heavy cut that boils below 150 °C and includes mostly C₇–C₉ hydrocarbons [2]. Prior to reforming, the SR is treated in a hydrotreating unit. Reforming is a form of reaction where the components paraffins, naphthenes and aromatics are rearranged to produce higher gasolines. These components typically contain 40–70 wt% paraffins, 20–50 wt% naphthenes, 5–20 wt% aromatics [3]. The various reforming reactions are dehydrogenation reactions of the paraffins to naphthenes and dehydrogenation of the naphthenes to aromatics, dehyrocylization of paraffins to naphthenes, isomerization of normal paraffins to iso-paraffins, cracking of paraffins and naphthenes and hydrodealkylation of aromatics. Modelling and simulation of catalytic reforming has been considered by different researchers to determine the plant performance or product distribution behaviour. Hu carried out [4] steady state modelling of semi regenerative catalytic reformer with 17 lumps and monitored the plant performance on stream. Modelling and parameter estimation was carried out by Stijepovic with benchmarking with industrial data and a good prediction of concentration and hydrogen was obtained [5]. Optimization of catalyst distribution was carried by [6] to determine the optimal catalyst distribution in the reactors. The dynamics of a process tells a lot about the behaviour of the system and the effect of perturbation due to changes in operational conditions such as temperature. It gives an indication on the transient behaviour and the time for the system to reach steady state. This information helps whenever there is need for plant scale up or process modification. Few works have been carried out on the dynamics of the reaction [7,8,9,10]. The modelling and simulation works carried out by [11,12] (steady state) and [10] (dynamic) use the same reactor configurations, same feed throughput and feed stock properties. It requires rigorous modelling method to model and simulate dynamically series of reactors concurrently together having three or four heaters and reactors in series. It is a difficult task and in [10], one reactor was simulated as a first contribution In this work, first a modelling and simulation of catalytic reforming unit by [11,12] was performed on gPROMs software and was compared with that performed on MATLAB ODE45 by the earlier authors. The validated model is then used to study the steady state and dynamic behaviour of an industrial commercial plant of Kaduna Refining and Petrochemical Company, using four reactors was performed sequentially using the gPROMs software (Version 5.0.0 Process System Enterprise Limited, London, UK). gPROMs software was found capable and robust in performing complex dynamic modelling and simulations.

2. Process Description and Process Variables

Naphtha catalytic reforming reactions include Dehydrogenation, Dehydrocyclization, Isomerization and Hydrocracking reactions which are performed on both the metallic and acidic parts of the catalyst. A simple network reaction is given in Figure 1.

M represents the metal catalyst site, A represents the Acid catalyst site. The main reactions involved in reforming processes are catalysed essentially either by the acid support or the metal functions.

Dehydrogenation	Metallic function
Dehydrocyclization	Metallic + acidic functions
Isomerization	Metallic + acidic functions
Hydrogenolysis	Metallic function
Hydrocracking	Metallic + acidic function

A Naphtha catalytic reforming unit (CRU) typically is made up of a preheat system, that is, a long vertical heat exchanger with large surface area to allow for feed preheating prior to further heating in the charge heater usually three or four heaters and reactors in series and a low pressure separator. The hydrogen is flashed and separated in the separator where part of it is recycled to the meet the feed again continuously before it enters the first heater, while the heavier liquid (reformate) is sent to a stabilization column. Reformate is obtained as a bottom product from the stabilization column while the liquefied petroleum gases (LPG) are obtained from the top of the stabilization after condensation in a trim cooler and condenser. It is important to mention that the flashed hydrogen from the low pressure separator goes through a drying section, a dryer consisting some alumina balls in order to reduce the amount of water vapour in the gas to a maximum of 25 ppm before recycling back to the feed. Also the bottom of the low pressure separator, liquid reformate prior to stabilization is re-contacted with hydrogen from a dry gas compressor at high pressure (high pressure separator) in order to strip more of hydrogen if still contained in the reformate. This high pressure separator otherwise called re-contactor is the source of hydrogen to other units in refinery where it is distributed at high pressure. These reforming reactions globally are endothermic in nature, multiple charge heaters are used to raise the temperature of the reaction to the desired value (400–500 °C) [13]. During the reaction in the catalytic bed in the reactor, one of the most significant reaction which is the dehydrogenation of naphthenes to aromatics, is very fast and more endothermic. This results in a high decrease of temperature in the reactor. The temperature of the exit stream from the first reactor is raised and sent to the next reactor. The rate of reaction of the feed decreases whenever it moves through the reactors in series. The reaction becomes less endothermic and the differential temperature along the reactor falls and the quantity of heat required to reheat the stream also decrease. The treated naphtha fed to the reactor is combined with a gas stream containing about 60 to 90 mole % hydrogen whose pressure depends on compressor capacity. Thus, total reactor charge is the naphtha feed plus the hydrogen recycle gas stream. The reactor feed must be raised to the proper temperature for the reforming reactions to occur when the charge contacts the catalyst. As shown in Figure 2 total reactor charge is heated, at first by exchange with effluent from the last reactor and is finally brought up to the first reactor inlet temperature in the first charge heater which is also called the first inter heater. The reactor effluent-to-feed exchanger recovers the heat from the reactor effluent and provides it to the reactor feed. Thus, it is one key to energy conservation in a catalytic reformer. The reactor effluent which may be as high as 750 to 790 K, must be cooled to 315 to 320 K for flash separation of hydrogen from reformate. This heat exchange is carried out in several banks of heat exchangers, arranged for parallel flow. Some of the newer units use one large single-pass vertical exchanger for this purpose. Catalytic reforming is a gas phase process [13]. After passing through the reactor effluent-to-feed exchanger and the charge heater, the total reactor charge is 100% vapour at the appropriate reaction temperature and is ready to contact the reforming catalyst. The most commonly used catalyst type is platinum on alumina support. Recent developments in catalysis have produced some bi-metallic and tri-metallic catalysts, such as Pt-Re or Pt-Ir on alumina support [14]. The flow scheme in Figure 2 shows the four reactors in series.

The hydrogen circulation to the naphtha charge is called the hydrogen recycle. The hydrogen recycle is maintained in order of control the rates of deactivation of the catalyst. Excess hydrogen from the separator is sent to feed gas or to hydrogen consuming units in the refinery such as hydro-treaters and petrochemical plants. The separator liquid, comprised mostly of the desired reformate product but also containing hydrogen, methane, ethane, propane and butanes is pumped to the reformate stabilizer. The stabilized reformate is sent to storage for gasoline blending while the C₃ and C₄ are sent to LPG storage tanks. The reformer that is modelled in this work is a semi regenerative type of reformer consisting of four reactors in series. When the catalyst gets deactivated after time and has reached end of run (EOR) where the quality of the product is not obtained, the plant is shut down and the deactivated catalyst is regenerated and this is called a catalyst cycle. A cycle may be terminated for a number of reasons but the most common is poor catalyst performance due to deactivation of the catalyst by coke formation on its surface as a result of severity of operation mainly high temperature. This usually leads to a loss of reformate yield. The regeneration procedure restores the catalyst to near fresh-catalyst condition. The temperature, pressure, space velocity, H₂/Hc ratio and feed quality are the variables affecting rate of the reaction and the sensitivity of the product distribution to these variables are studied in the result section of this work. A summary of the effect of these variables are shown in

Table 1. Effect of variables on product distribution, coke and Research octane number (RON).

Increased		RON	Reformate Yield	H2 Yield	Coke Deposit
Pressure
Temperature
Space velocity
H2/Hc ratio
Naphtha Quality	A + 0.85 N
	End boiling point
	Initial boiling point

[15].

Process Classification

Naphtha catalytic reforming process is typically classified based on the method and procedure of operation and regeneration of catalyst. They are semi catalyst regeneration, cyclic catalyst regeneration and continuous catalyst regeneration. Different refineries have different processes but most widely used is the semi catalyst regeneration, then continuous catalyst regeneration and the less commonly used is cyclic catalyst regeneration. Nowadays, more catalytic reformers are developed and designed for continuous catalyst regeneration and the previous semi catalyst regeneration units are being revamped to operate in continuous regeneration mode. A semi catalyst regeneration process typically consists of three or four reactors in series in a fixed bed catalyst system and operates within a time line where the catalyst gets deactivated, a cycle length and it is regenerated. Usually a period of 6 months to one year. The advantage of this process type is the simplicity of the mode of operation but one major drawback is that the plant has to be shut down for catalyst regeneration after end of run (EOR). Another drawback is the application of high pressure and hydrogen hydrocarbon (HHR) to avoid quick catalyst deactivation. The cyclic catalyst regeneration unit is similar to the semi catalyst regeneration mode but it has an extra swing reactor which is used as a spare reactor whenever the catalyst of one of the reactors is deactivated. It serves as a reactor for the unit while it is regenerated and this prevents the total shut down of the plant. The drawback of this process is the complexity and challenges of reactor swing from one to the other although it has the advantage of use of lower pressure than the semi catalyst regeneration mode. The continuous catalyst regeneration mode is a moving bed mode where the catalyst is continuously regenerated. It does not require shut down for catalyst regeneration and this implies less down time for the unit and more production time. The major challenge of this process is the complexity of mechanism of engaging and disengaging of the product from the last reactor to the regenerator. The hydrogen from the catalyst to be regenerated has to be displaced with nitrogen before it is sent to the regenerator to avoid explosion since the regenerator operates at a high temperature of about 600–700 °C. The catalyst after cooling must be swept of nitrogen before is sent to the reactor. An important aspect in this type of process is disharmony at times between the catalyst suppliers and the technology owners where some catalysts used may not be suitable for the technology of the process. The major operational draw backs of the previous modes have solved in this mode of operation. It uses low pressure and low hydrogen hydrocarbon ratio and gives higher hydrogen purity and higher yield of aromatics. In this work, a semi catalyst regeneration mode of CRU, in Kaduna refining and petrochemical company is modelled.

3. Kinetic Models

In this work, various models have been developed for both steady state and dynamic process of the catalytic reforming reaction with different lumps to investigate the behaviour of the reaction and products distribution. The reaction steps and equations of the model are given in this section.

Table 2. Table of reaction steps and rate constants of the Model (adapted from [10]).

Reaction Steps of the Reactions	Rate Constants (kg/h·kg)−1
Dehydrogenation reactions of paraffins (P) to naphthenes (N)
1. $P_{11}\to N_{11}$ + $H_2$	0.00	k1
2. $P_{10}\to N_{10}$ + $H_2$	2.54	k7
3. $P_9\to N_9$ + $H_2$	1.81	k13
4. $P_8\to N_8$ + $H_2$	1.33	k18
5. $P_7\to N_7$ + $H_2$	0.58	k23
6. $P_6\to N_6$ + $H_2$	0.00	k27
7. $P_6\to$ $MCP$ + $H_2$	0.00	k28
Hydrocracking reactions of paraffins (P)
8. $P_{11}$ + $H_2\to$ $P_{10}$ + $P_1$	0.00	k2
9. $P_{11}$ + $H_2\to$ $P_9$ + $P_2$	0.00	k3
10. $P_{11}$ + $H_2\to$ $P_8$ + $P_3$	0.00	k4
11. $P_{11}$ + $H_2\to$ $P_7$ + $P_4$	0.00	k5
12. $P_{11}$ + $H_2\to$ $P_6$ + $P_5$	0.00	k6
13. $P_{10}$ + $H_2\to$ $P_9$ + $P_1$	0.49	k8
14. $P_{10}$ + $H_2\to$ $P_8$ + $P_2$	0.63	k9
15. $P_{10}$ + $H_2\to$ $P_7$ + $P_3$	1.09	k10
16. $P_{10}$ + $H_2\to$ $P_6$ + $P_4$	0.89	k11
17. $P_{10}$ + $H_2\to$ $2P_5$	1.24	k12
18. $P_9$ + $H_2\to$ $P_8$ + $P_1$	0.30	k14
19. $P_9$ + $H_2\to$ $P_7$ + $P_2$	0.39	k15
20. $P_9$ + $H_2\to$ $P_6$ + $P_3$	0.68	k16
21. $P_9$ + $H_2\to$ $P_5$ + $P_4$	0.55	k17
22. $P_8$ + $H_2\to$ $P_7$ + $P_1$	0.19	k19
23. $P_8$ + $H_2\to$ $P_6$ + $P_2$	0.25	k20
24. $P_8$ + $H_2\to$ $P_5$ + $P_3$	0.43	k21
25. $P_8$ + $H_2\to$ $2P_4$	0.35	k22
26. $P_7$ + $H_2\to$ $P_6$ + $P_1$	0.14	k24
27. $P_7$ + $H_2\to$ $P_5$ + $P_2$	0.18	k25
28. $P_7$ + $H_2\to$ $P_4$ + $P_3$	0.32	k26
29. $P_6$ + $H_2\to$ $P_5$ + $P_1$	0.14	k27
30. $P_6$ + $H_2\to$ $P_4$ + $P_2$	0.18	k29
31. $P_6$ + $H_2\to$ $2P_3$	0.27	k31
32. $P_5$ + $H_2\to$ $P_4$ + $P_1$	0.12	k32
33. $P_5$ + $H_2\to$ $P_3$ + $P_2$	0.15	k33
Dehydrogenation reactions of Naphthenes(N)
34. $N_{11}\to A_{11}$ + $3H_2$	0.00	k35
35. $N_{10}\to A_{10}$ + $3H_2$	24.5	k40
36. $N_9\to A_9$ + $3H_2$	24.5	k44
37. $N_8\to A_8+3H_2$	21.5	k49
38. $N_7\to A_7+3H_2$	9.03	k52
39. $N_6\to A_6+3H_2$	4.02	k54
40. $N_{11}$ + $H_2\to$ $P_{11}$	0.00	k34
41. $N_{10}$ + $H_2\to$ $P_{10}$	0.54	k39
42. $N_9$ + $H_2\to$ $P_9$	0.54	k44
43. $N_8$ + $H_2\to$ $P_8$	0.47	k48
44. $N_7+H_2\to P_7$	0.20	k51
45. $N_6+H_2\to P_6$	1.48	k53
46. $N_6\to MCP$	0.00	k55
47. $MCP+H_2\to P_6$	0.00	k56
48. $MCP\to N_6$	0.00	k57
Hydrocracking reactions of Naphthenes(N)
49. $N_{11}$ + $H_2\to$ $N_{10}$ + $P_1$	0.00	k36
50. $N_{11}$ + $H_2\to$ $N_9$ + $P_2$	0.00	k37
51. $N_{11}$ + $H_2\to$ $N_8$ + $P_3$	0.00	k38
52. $N_{10}$ + $H_2\to$ $N_9$ + $P_1$	1.84	k41
53. $N_{10}$ + $H_2\to$ $N_8$ + $P_2$	1.34	k42
54. $N_{10}$ + $H_2\to$ $N_7$ + $P_3$	0.80	k43
55. $N_9$ + $H_2\to$ $N_8$ + $P_1$	1.27	k45
56. $N_9$ + $H_2\to$ $N_7$ + $P_2$	1.27	k46
57. $N_8+H_2\to N_7+P_1$	0.09	k47
Hydrodealkylation reactions of Aromatics(A)
58. $A_{11}+H_2\to A_{10}+P_1$	0.00	k59
59. $A_{11}+H_2\to A_9+P_2$	0.00	k60
60. $A_{10}+H_2\to A_9+P_1$	0.06	k62
61. $A_{10}+H_2\to A_8+P_2$	0.06	k63
62. $A_{10}+H_2\to A_7+P_3$	0.00	k64
63. $A_9+H_2\to A_8+P_1$	0.05	k66
64. $A_9+H_2\to A_7+P_2$	0.05	k67
65. $A_8+H_2\to A_7+P_1$	0.01	k69
66. $A_{11}+4H_2\to P_{11}$	0.00	k58
67. $A_{10}+4H_2\to P_{10}$	0.00	k61
68. $A_9+4H_2\to P_9$	0.16	k65
69. $A_8+4H_2\to P_8$	0.16	k68
70. $A_7+4H_2\to P_7$	0.16	k70
71. $A_6+3H_2\to N_6$	0.45	k71

shows the reaction steps and rate constants for the various reactions used in the model.

The kinetic rate equations for the components are given in Equations (1)–(25) where $SV$ is space velocity and $\frac{1}{sv}$ is residence time.

$$\begin{array}{cc}\hfill \frac{d{P}_1}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}=& K_2{P}_{11}+K_8{P}_{10}+K_{14}{P}_9+K_{19}{P}_8+K_{24}{P}_7+K_{29}{P}_6+K_{32}{P}_5+K_{36}{N}_{11}\hfill \\ \hfill & +K_{41}{N}_{10}+K_{45}{N}_9+K_{50}{N}_8+K_{50}{A}_{11}+K_{62}{A}_{10}+K_{66}{A}_9+K_{99}{A}_8\hfill \end{array}$$

(1)

$$\begin{array}{cc}\hfill \frac{d{P}_2}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}=K_3{P}_{11}+=& K_9{P}_{10}+K_{15}{P}_9+K_{20}{P}_8+K_{25}{P}_7+K_{30}{P}_6+K_{33}{P}_5+K_{37}{N}_{11}+\hfill \\ \hfill & K_{42}{N}_{30}+K_{47}{N}_9+K_{60}{A}_{11}+K_{62}{A}_{10}+K_{67}{A}_9\hfill \end{array}$$

(2)

$$\begin{array}{cc}\hfill \frac{d{P}_3}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}==& K_4{P}_{11}+K_{10}{P}_{10}+K_{16}{P}_9+K_{21}{P}_8+K_{26}{P}_7+2K_{31}{P}_6+K_{33}{P}_5+K_{38}{N}_{11}\hfill \\ \hfill & +K_{43}{N}_{10}+K_{64}{A}_{10}\hfill \end{array}$$

(3)

$$\frac{d{P}_4}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}=K_5{P}_{11}+K_{11}{P}_{10}+K_{17}{P}_9+2K_{22}{P}_8+K_{25}{P}_7+K_{10}{P}_6+K_{32}{P}_5$$

(4)

$$\frac{d{P}_5}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}=K_6{P}_{11}+2K_{12}{P}_{10}+K_{17}{P}_{19}+K_{21}{P}_8+K_{24}{P}_7+K_{29}{P}_6+(K_{32}+K_{33})P_5$$

(5)

$$\begin{array}{cc}\hfill \frac{d{P}_6}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}==& K_6{P}_{11}+K_{11}{P}_{10}+K_{16}{P}_9+K_{20}{P}_8+K_{25}{P}_7+K_{53}{N}_6+K_{96}MPC\hfill \\ \hfill & -(K_{27}+K_{28}+K_{29}+K_{30}+K_{31})P_6\hfill \end{array}$$

(6)

$$\begin{array}{cc}\hfill \frac{d{P}_7}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}==& K_5{P}_{11}+K_{10}{P}_{10}+K_{15}{P}_9+K_{19}{P}_8+K_{51}{N}_7+K_{70}{A}_7\hfill \\ \hfill & -(K_{23}+K_{24}+K_{25}+K_{26}+K_{27})P_7\hfill \end{array}$$

(7)

$$\begin{array}{cc}\hfill \frac{d{P}_8}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}==& K_4{P}_{11}+K_9{P}_{10}+K_{14}{P}_9+K_{48}{N}_8+K_{68}{A}_8\hfill \\ \hfill & -(K_{18}+K_{19}+K_{20}+K_{21}+K_{22})P_8\hfill \end{array}$$

(8)

$$\frac{d{P}_9}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}=K_3{P}_{11}+K_8{P}_{10}+K_{44}{N}_9+K_{66}{A}_9-(K_{13}+K_{14}+K_{15}+K_{16}+K_{17})P_9$$

(9)

$$\frac{d{P}_{10}}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}=K_2{P}_{11}+K_{39}{N}_{10}+K_{61}{A}_{10}-(K_7+K_8+K_9+K_{10}+K_{11}+K_{12})P_{10}$$

(10)

$$\frac{d{P}_{11}}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}=K_{34}{N}_{11}+K_{58}{A}_{11}-(K_1+K_2+K_3+K_4+K_5+K_6)P_1$$

(11)

$$\frac{dMCP}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}=K_{28}{P}_6+K_{55}{N}_6-(K_{56}+K_{57})$$

(12)

$$\frac{d{N}_6}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}=K_{27}{P}_6+K_{57}MCP+K_{71}{A}_6-(K_{53}+K_{54}+K_{55})N_6$$

(13)

$$\frac{d{N}_7}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}=K_{23}{P}_7+K_{43}{N}_{10}+K_{47}{N}_{10}+K_{50}{N}_9-(K_{51}+K_{52})N_7$$

(14)

$$\frac{d{N}_8}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}=K_{18}{P}_8+K_{38}{N}_{11}+K_{42}{N}_{10}+K_{46}{N}_9-(K_{48}+K_{49}+K_{50})N_8$$

(15)

$$\frac{d{N}_9}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}=K_{13}{P}_9+K_{37}{N}_{11}+K_{41}{N}_{10}-(K_{44}+K_{45}+K_{46}+K_{47})N_9$$

(16)

$$\frac{d{N}_{10}}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}=K_7{P}_{10}+K_{36}{N}_{11}-(K_{39}+K_{40}+K_{41}+K_{42}+K_{43})N_{10}$$

(17)

$$\frac{d{N}_{11}}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}=K_1{P}_{11}-(K_{34}+K_{35}+K_{36}+K_{37}+K_{38})N_{11}$$

(18)

$$\frac{d{A}_6}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}=K_{49}{N}_8+K_{63}{A}_{10}+K_{66}{A}_{10}-(K_{68}+K_{69})A_8$$

(19)

$$\frac{d{A}_7}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}=K_{49}{N}_8+K_{63}{A}_{10}+K_{66}{A}_{10}-(K_{68}+K_{69})A_8$$

(20)

$$\frac{d{A}_8}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}=K_{49}{N}_8+K_{63}{A}_{10}+K_{66}{A}_{10}-(K_{68}+K_{69})A_8$$

(21)

$$\frac{d{A}_9}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}=K_{45}{N}_9+K_{60}{A}_{11}+K_{62}{A}_{10}-(K_{65}+K_{66}+K_{67})A_9$$

(22)

$$\frac{d{A}_{10}}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}=K_{40}{N}_{10}+K_{59}{A}_{11}-(K_{61}+K_{62}+K_{63}+K_{64})A_{10}$$

(23)

$$\frac{d{A}_{11}}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}=K_{35}{N}_{11}-(K_{58}+K_{59}+K_{60})A_{11}$$

(24)

$$\begin{array}{cc}\hfill \frac{d{H}_2}{d(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$SV$}\right.)}==& a_1{P}_{11}+a_2{P}_{10}+a_3{P}_9++a_4{P}_8+a_5{P}_7+a_6{P}_6+a_7{P}_5+b_1{N}_{11}+b_2{N}_{10}\hfill \\ \hfill & +b_3{N}_9+b_4{N}_8+b_5{N}_7+b_6{N}_6+c_1{A}_{11}+c_2{A}_{10}+c_3{A}_9+c_4{A}_8+c_5{A}_7\hfill \\ \hfill & +c_6{A}_6\hfill \end{array}$$

(25)

where a₁ = $k_1-(k_2+k_3+k_4+k_5+k_6)$, a₂ = $k_7-(k_8+k_9+k_{10}+k_{11}+k_{12}$, a₃ = $k_{13}-(k_{14}+k_{15}+k_{16}+k_{17})$, a₄ = $k_{18}-(k_{19}+k_{20}+k_{21}+k_{22})$, a₅ = $k_{23}-(k_{24}+k_{25}+k_{26})$, a₆ = $k_{27}+k_{28}-(k_{29}+k_{30}+k_{31})$, a₇ = $-(k_{32}+k_{33})$, b₁ = $3k_{35}-(k_{34}+k_{36}+k_{37}+k_{38})$, b₂ = $3k_{40}-(k_{39}+k_{41}+k_{42}+k_{43}$), b₃ = $3k_{45}-(k_{44}+k_{46}+k_{47})$, b₄ = $3k_{49}-(k_{48}+k_{50})$, b₅ = $3k_{52}-k_{51}$, b₆ = $3k_{54}-k_{53}$, c₁ = ($4k_{58}+k_{59}+k_{60}$), c₂ = $4k_{61}+k_{62}+k_{63}+k_{64}$, c₃ = $4k_{65}+k_{66}+k_{67}$, c₄ = $4k_{68}+k_{69}$, c₅ = $4k_{70}$, c₆ = $3k_{71}$.

$$K_i=K_i^0\left[\frac{E_{Ai}}{R}\left(\frac{1}{T_0}-\frac{1}{T}\right)\right]{\left(\frac{P}{P_0}\right)}^{a_k}$$

(26)

where K_i is the kinetic constants for the reactions, T₀ and P₀ are the reference temperatures and pressures. E_Ai (activation energy) and w (exponential values of pressure) are taken from [10] and are given in

Table 3. Activation Energies for the reactions and Exponential values of Pressure (adapted from [10]).

Activation Energy of Reactions
Reactions	Activation Energy(kcal/mol)
Dehydrocyclization of paraffins ($P_n\to$ $N_n$)	45.00
Hydrocracking of paraffins ($P_n\to$ $P_{n-i}+P_i$)	55.00
Dehydrogenation of naphthenes ($N_n\to$ $A_n$)	30.00
Hydrodealkylation of naphthenes ($N_n\to$ $N_{n-i}+P_i$)	55.00
Ring opening of naphthenes ($N_n\to$ $P_n$)	45.00
Hydrodealkylation of aromatics ($A_n\to$ $A_{n-i}+P_i$)	40.00
Ring opening of aromatics ($A_n\to$ $P_n$)	45.00
Hydrogenation of aromatics ($A_n\to$ $N_n$ )	30.00
Exponential values of Pressure
Reaction	w
Other reactions	0
Hydrocracking of paraffins	0.433
Dehydrogenation/hydrogenation of aromatics	0.0
Hydrodealkylation of aromatic & naphthenes	0.5
Dehydrocyclization of paraffins	−0.7

. k_i, E_Ai, w, T₀ and P₀ are also obtained from [10].

4. Mathematical Model

The kinetic rate equations from Equations (1)–(25) which determine the rate of the reaction of each component is used in the mathematical modelling Equations of (27)–(32) to determine the concentration, temperature and pressure distributions of the components describing the system. The rates of the equations are defined as $r_j$ in Equations (27) and (28) and $r_j{\gamma }_i$ summations of the rate equations and their stoichiometric coefficients. These model equations, (27)–(32), in cooperated with Equations (1)–(25) are solved simultaneously on gPROMs. Equations (27) and (28) are equations representing the dynamic behaviour of the system and change with time dt. At steady state, the system does not change with time and dt becomes zero. So the equations can be used for both steady state and dynamic behaviour depending on the application.

$$-\frac{dF}{dW}+{{\displaystyle \sum }}_{j=1}^{Nr}\left(r_j{\gamma }_i\right)=\frac{\epsilon d\left(ci\right)}{\rho bdt}.$$

(27)

$$\frac{dT}{dt}=\frac{-({{\displaystyle \sum }}_{i=1}^{NC}{F}_{i}C{P}_i)\frac{dT}{dW}{{\displaystyle \sum }}_{j=1}^{Nr}\left(r_j{\gamma }_i\right)\left(-\Delta H_{Rj}\right)}{Cpcat+\frac{\epsilon }{\rho b}{{\displaystyle \sum }}_{i=1}^{NC}{C}_{i}C{P}_i}$$

(28)

$$\frac{dP}{dz}=\left[1.75\times {10}^{-5}\frac{\left(1-\epsilon \right)}{{\epsilon }^3}\frac{G^2}{\rho d_P{g}_C}+1.5\times {10}^{-5}\frac{\left(1-\epsilon \right)}{{\epsilon }^3}\frac{G\mu }{\rho d_P^2{g}_C}\right]\times \rho A$$

(29)

$$\Delta H_R={\displaystyle \sum }{v}_P{H}_{fP}-{\displaystyle \sum }{v}_r{H}_{fr}$$

(30)

$$H_{ri}=H_{ri}^0+{{\displaystyle \int }}_{298K}^{T}CPdT$$

(31)

$$CP=A1+B1T+C1T^2+D1T^3$$

(32)

5. RON Model

The estimation of the Research octane number, RON of hydrocarbons for the feed and products can be done using different methods of prediction. RON can be calculated for each pure component by using a polynomial equation that correlates to the normal boiling point [16]. Equation (33) represents the model equation for the RON estimation.

RON = Aa + BbT + CcT² + DdT³ + EeT⁴

(33)

where T = TBP × 0.01, TBP represents the normal boiling point (°C) and Aa, Bb, Cc, Dd and Ee, are coefficients. The RON of a hydrocarbon mixture is calculated by assuming that the mixture consists of paraffins, naphthenic hydrocarbons and aromatics. The equation is expressed as the sum of the RON for each pure component multiplied by the volume fractions of the components. Equation (27) is the material balance equation where $\rho b$ is the bulk density and $\epsilon$ the bed voidage. Equation (28) is the heat balance of the reaction to determine the temperature behaviour where $-\Delta H_{Rj}$ is the heat of the reaction of jth component and Cpcat the heat capacity of the catalyst. Equation (29) is the pressure profile of the reaction where Dp is particle diameter, G the mass flux, $\epsilon$ the bed voidage, $\rho$ the gas density and A the reactor cross sectional area. $H_{ri}^0$, the standard heat of formation of the components. The constants A, B, C, D used for the heat capacity, Cp, in Equations (31) and (32) are taken from [16]. These equations are solved simultaneously using gPROMs to determine the behaviour of the system. gPROMs is a robust mathematical software that can solve the equations for all the four reactors in series dynamically given the behaviour of the paraffins, naphthenes and aromatics as well as the dynamics of the temperature in the reactors. All the components behaviours were determined.

6. Model Validation

The model validation was carried out by modelling and simulating the commercial data of [11,12] on gPROMs to ascertain the capability and ruggedness of the mathematical software in complex modelling. The modelling of [11] was performed using MATLAB software and the result is compared with that obtained with gPROMs. The configuration of the commercial reformer and Feed Properties are given in

Table 4. Configuration of Commercial Catalytic Reformer [11].

Reactor	Height (m)	Diameter (m)	Catalyst (Kg)	WHSV (h−1)
1	4.902	2.438	9130	16
2	5.410	2.819	13,820	10.6
3	6.452	2.971	22,820	6.4
4	8.208	3.505	42,580	3.4

;

Table 5. Feed stock Properties of Commercial Catalytic Reformer [11].

Property
MW	104.8
Specific gravity	0.7406
IBP	88
10%	101
90%	155
EBP	180
Total paraffins	59.11
Total naphthenes	20.01
Total aromatics	20.88

respectively.

The commercial reformer of [11] was simulated with the properties of Table 4 and Table 5 respectively with four semi regenerative reactors in series. The reformer’s through put is 30 MBPD at inlet temperature of 495 °C, pressure of 10.5 kg/cm² using MATLAB mathematical tool on ODE45.

The simulation was performed on gPROMs software to validate its capability in modelling this problem. The result obtained from gPROMs is compared with that obtained with MATLAB from [11] and are shown in

Table 6. Comparison between results from MATLAB [11] and gPROMs.

Components	Actual	Simulated [11]	Simulated gPROMs	Absolute Difference [11]	Absolute Difference gPROMs
P5	4.83	5.37	5.94	0.21	1.11
P6	15.08	15.17	13.18	0.25	1.9
P7	11.2	11.32	11.31	0.28	0.11
P8	4.79	4.44	5.52	0.18	0.73
P9	2.3	2.20	2.21	0.10	0.09
P10	0.15	0.16	0.33	0.01	0.18
P11	0.01	0.01	0.00	0.00	0.01
MCP	1.25	1.22	1.23	0.11	0.02
N6	0.19	0.07	0.23	0.04	0.04
N7	0.38	0.35	0.50	0.03	0.12
N8	0.59	0.56	0.36	0.03	0.23
N9	0.14	0.15	0.11	0.01	0.03
N10	0.02	0.02	0.023	0.00	0.003
N11	0.00	0.00	0.00	0.00	0
A6	5.43	5.68	5.36	0.25	0.07
A7	15.03	14.96	13.62	0.07	0.41
A8	18.63	19.22	17.78	0.61	0.85
A9	13.68	13.07	12.15	0.61	1.53
A10	4.72	4.56	4.34	0.16	0.38
A11	1.58	1.51	1.32	0.07	0.26

. The result showed a good and comparable result between the actual and that simulated with gPROMs. This is an indication of the capability of the mathematical tool in solving complex model equations and hence can be used as a validation for modelling other catalytic reformers using same kinetic models.

7. Results and Discussion

The model was simulated based on the kinetic data from [11,12] and validated with Kaduna Refining & Petrochemical Company, KRPC Commercial Catalytic Reformer data from Nigeria.

Table 7. Configuration of KRPC Commercial Catalytic Reformer.

Reactor Number	Height (m)	Diameter (m)	Catalyst (Kg)	WHSV (h−1)
1	5.63	1.9	9572	5.56
2	5.83	2.1	12119	4.39
3	6.51	2.3	16231	3.28
4	7.26	2.7	24938	2.13

and

Table 8. Feed Properties KRPC Commercial Catalytic Reformer.

FEED STOCK PROPERTIES
MW	103.7
SG	0.76
IBP	88
5%	97
10%	102
15%	110
30%	118
50%	131
75%	147
95%	152
EBP	170
RVP	0.97
SULFUR(PPM)	0.32

show the Configuration and Feed Properties of KRPC Reformer. The results show a good agreement as shown in Figure 3 and

Table 9. Relative errors between Industrial and Simulated results of KRPC plant.

COMP	Simulated	Industrial	Error
P5	4.22	4.81	0.595
P6	7.16	7.07	0.089
P7	8.22	9.35	1.135
P8	7.20	7.01	0.194
P9	5.84	3.94	1.902
P10	4.70	1.59	3.111
MCP	0.42	0.42	0
N6	2.76	3.09	0.326
N7	2.34	2.49	0.152
N8	0.31	0.68	0.365
N9	0.20	0.34	0.136
N10	0.19	0.12	0.077
A6	4.29	3.53	0.764
A7	16.02	16.57	0.550
A8	18.71	19.72	1.009
A9	12.27	14.87	2.597
A10	4.86	4.21	0.646

.

7.1. Steady State Simulation Results

The time derivatives of the model equation are set to zero to simulate the process at steady state. The modelling results from the model equations of (27)–(32) gives the concentration and temperature distribution only across the reactor height. In this work, four reactors are simulated with the reactor configurations given in Table 4. The total reactor height is 25.232 m which is total sum of the heights of the four reactors as in Table 4 and is denoted as height of reactor in Section 7.1 therefore the height in Figure 4, Figure 5, Figure 6, Figure 7, Figure 8 and Figure 9 refer to reactors total height of 25.232 m. In Section 7.2, the different reactors heights for the four reactors were used differently to show variation with both the height and time as shown in Figure 10, Figure 11, Figure 12, Figure 13 and Figure 14.

Figure 4 shows the concentration profile of paraffins along the reactor height. P₅, P₆, P₇, P₈, P₉, P₁₀ are paraffins with hydrocarbon number 5, 6, 7, 8, 9, 10 respectively. Normal Paraffin components basically undergo three major reactions during naphtha reforming. They undergo hydrocracking to lighter paraffins that is, methane, ethane, propane and butane as in shown in reactions 8 to 33 in Table 2, isomerization to iso-paraffins which are not modelled in this work, it is a slow reaction with slow reaction rate, dehydrocyclization to naphthenes, a very slow reaction which leads to a decrease in the paraffins. These reactions are reactions 1–7 from Table 2. Dehydrocyclization reaction becomes easier as the molecular weight of the paraffins increases as in P₈, P₉ and P₁₀ while P₇ shows little increase. P₅ and P₆ increase due to hydrocracking. The main effects of hydrocracking are decrease of paraffins (C₅+) in the reformate, decrease in hydrogen production, increase of LPG production and hydrogenolysis. The isomerization reactions are fast, slightly exothermic and do not affect the number of carbon atoms. The thermodynamic equilibrium of iso-paraffins to paraffins depends mainly on the temperature and pressure has no effect. The paraffins isomerization results in a slight increase of the octane number. These reactions are promoted by the acidic function of the catalyst support. The paraffin dehydrocyclization step becomes easier as the molecular weight of the paraffin increases, from Table 2 we can see that the rates increase from 0.0, 0.58, 1.33,1.81, 2.54 from P₆ to P₁₁ however the tendency of paraffins to hydrocrack increases concurrently. Kinetically, the rate of dehydrocyclization increases with low pressure and high temperature. To sum up, the dehydrocyclization of P₆ paraffins to benzene is more difficult than that of C₇ paraffin to toluene, which itself is more difficult than that of C₈ paraffin to xylenes. Accordingly, the most suitable fraction to feed a reforming process is the C₇–C₁₀ fraction.

Figure 5 shows the behaviour of the naphthenes along the reactor height. N₆, N₇, N₈, N₉, N₁₀ are naphthenes with hydrocarbon number 6, 7, 8, 9, 10 respectively. The fastest reaction is aromatization that is, dehydrogenation of naphthenes, resulting in the large temperature drop due to its highly endothermic nature. The drastic reduction is due to the aromatization to aromatics. Thermodynamically the reaction is highly endothermic and is favoured by high temperature and low pressure. In addition, the higher the number of carbon atoms the higher the aromatics production at equilibrium from N₈, N₉, N₁₀. This can be seen from Table 2 in equations (34)–(39) where the reaction rates increase from 4.02, 9.03, 21.5, 24.5 and 24.5 respectively for N₆, N₇, N₈, N₉ and N₁₀. From a kinetic point of view, the rate of reaction increases with temperature. The naphthenes also undergo hydrocracking to lighter hydrocarbons leading to their decrease as shown in reactions 49–57 in Table 2.

Figure 6 shows the behaviour of aromatics along the reactor height. A₆, A₇, A₈, A₉, A₁₀ are aromatics with hydrocarbon number 5, 6, 7, 8, 9, 10 respectively. The dehydrogenation of naphthenes as seen from Figure 5 increases the aromatics shown in Figure 6. The sharp increase is as a result of drastic aromatization of the naphthenes favoured by the metallic sites of the catalyst. The rate of benzene formation is lower due to lower carbon number of N₆ while A₇, A₈, A₉ increase rapidly due to higher carbon number of N₇, N₈ and N₉. This can be seen from reactions 58–71 in Table 2.

Figure 7 shows the temperature profile of the reactions along the reactor height. The sharp drop in first and second reactor is due to the more endothermic dehydrogenation reaction of naphthenes to aromatics. Typically, dehydrogenation and isomerization reactions take place in the first reactor, dehydrogenation, isomerization, dehydrogenation and cracking in the second followed by dehydrogenation and cracking in the third and fourth reactors.

Figure 8 shows the increase in hydrogen yield along the height of the reactor. Figure 9 shows the increase in RON along the reactor height. This is due to the increase in the aromatics along the reactor height. The RON increases from reactor one to four due to increase in aromatization reaction which gives a higher research octane number, a parameter for antiknock in the gasoline engine. Straight run naphtha has a low RON and cannot be used as gasoline and hence the reforming reaction of the components to give a gasoline with higher research octane number.

7.2. Dynamic Simulation of the Four Reactors

The modelling results from the model Equations of (27)–(32) gives the concentration and temperature distributions only across the reactor height and the time simultaneously where their behaviours are predicted with changing both time height. In this work, four reactors are simulated with the reactor configurations given in Table 5.

Figure 10, Figure 11, Figure 12 and Figure 13 shows the 3D profiles of the dynamics of the temperatures in the four reactors respectively. It is a combination of the change along both the reactor height and time. For the variation with time, the interesting aspect could be seeing clearly in Figure 14 that the exit temperatures falls and rise with time before reaching a final temperature value unlike the steady state where it shows only a drop in temperature without seeing such interesting aspect as in Figure 7. From Figure 14, temperature decreases along the height simultaneously with time. The dynamic simulation was performed over 10 h but the time to reach steady state was 0.25 h maximum, hence the representation in Figure 14 stops at 0.4 h for clarity purpose. The drop in temperature in first reactor is mainly due to the high endothermic reaction of naphthenes conversion to aromatics. Dehydrogenation is the main reaction in this reactor, the fastest reaction among all the reactions and reaches a steady state at about 0.05 h. This can also be attributed to the lowest residence time of 0.179 h for the reactions in the first reactor. There are also some endothermic reactions in the second reactor where dehydrogenation and isomerization reactions continue, thereby increasing slightly the temperature with time. The temperature here reaches steady state at 0.12 h which is longer than the first value. In reactor three there is lesser temperature drop due to exothermic reactions of isomerization and hydrocracking, reaching a steady state at 0.20 h. Reactor four has the highest exothermic reactions and one can see a sharp rise in temperature. It reaches a steady state at 0.25 h.

Figure 15, Figure 16, Figure 17, Figure 18, Figure 19 and Figure 20 give the dynamic behaviours of the paraffins. Unlike in Figure 4 where it shows the steady state behaviour of the paraffins changing with the reactor height, in these Figures, the time to reach steady state in all the reactors where studied. The dynamics were run for 10 h but for clarity purpose it stops at 1 h. It takes 0.18 h in reactor 1 to attain steady state, 0.20 h to attain steady state in reactor 2, 0.22 h to attain steady state for reactor 3 while it takes 0.28 h to attain steady state in reactor 4.

Figure 21, Figure 22, Figure 23, Figure 24 and Figure 25 give the dynamic behaviours of the naphthenes. Unlike in Figure 5 where it shows the steady state behaviour of the naphthenes changing with the reactor height, in these Figures, the time to reach steady state in all the reactors where studied. The dynamics were run for 10 h but for clarity purpose it stops at 1 h. It takes 0.18 h in reactor 1 to attain steady state, 0.20 h to attain steady state in reactor 2, 0.22 h to attain steady state for reactor 3 while it takes 0.28 h to attain steady state in reactor 4.

Figure 26, Figure 27, Figure 28, Figure 29 and Figure 30 give the dynamic behaviours of the aromatics. Unlike in Figure 6 where it shows the steady state behaviour of the aromatics changing with the reactor height, in these Figures, the time to reach steady state in all the reactors where studied. The dynamics were run for 10 h but for clarity purpose it stops at 1 h. It takes 0.18 h in reactor 1 to attain steady state, 0.20 h to attain steady state in reactor 2, 0.22 h to attain steady state for reactor 3 while it takes 0.28 h to attain steady state in reactor 4.

8. Conclusions

In this paper, the catalytic naphtha reforming process of Kaduna refining and petrochemical company, Kaduna were studied via steady state and dynamic modelling using gPROMs software. The basic model of the process is taken from the literature and was developed in gPROMs platform and was validated using literature data. Note, while [10] simulated one of four reactors (in series) this work simulated all 4 reactors dynamically. The validated model is then used for Kaduna refining to study the behaviour of the concentrations of the paraffins, naphthenes and aromatics components along the reactor height and, the temperature profile of the reactors. The simulated result from gPROMs was in good agreement with industrial plant data of KRPC. The sensitivity of the reaction to the research octane number (RON) was also presented with increase from the first to the fourth reactor due to increase in aromatization. The dynamic models in form of partial differential equations were also solved numerically on gPROMs software and the concentration and temperature profiles for the 25 components determined. The time to attain quasi steady state in the reactor was determined. It is observed that the first reactor reaches the quasi steady state fastest compared to the second, third and fourth reactors due to the type of the reactions in the four different reactors. The 3D representation of the dynamics has also been presented. gPROMs for the first time is used to model and simulate naphtha catalytic reforming reactions and was found to be suitable and robust.