際際滷

際際滷Share a Scribd company logo

Simultaneousnonlinear two dimensional modeling of tubular reactor of hydrogen production unit of Shazand Arak refinery

Arash Nasiri
Arash Nasiri

This document presents mathematical models for simulating a packed tubular reactor for methane steam reforming. It develops two-dimensional partial differential equation models to generate radial and axial plots of component concentrations and temperature over time. Both steady-state and transient regimes are considered. The models assume ideal gas behavior and include mass and energy balance equations coupled through a reaction rate constant. Results are presented from solving the PDEs numerically using MATLAB, showing temperature, concentration and production rate profiles along the reactor under different operating conditions. In conclusion, the models capture the endothermic nature of the reforming reactions and how temperature initially decreases but then rises sharply in the reactor due to high heat flux.

1 of 7
Download to read offline
Abstract This paper developed two mathematical models of a packed tubular reactor loaded with reforming catalysts which results in a simultaneous, nonlinear, two dimensional PDE in cylindrical coordinates. It generates 2-D radial and axial plots of components concentration and temperature. Here both regimes of steady state and transient flow are assumed. We also ignored the presence of WGS reaction and analyzed the consequences by comparing the results with experimental data. KeywordsCatalytic reactor, Reformer modeling, methane Steam reforming, Shazand Arak refinery, Hydrogen production unit. I. INTRODUCTION Ethan steam reforming is one of the most common methods to produce hydrogen at industrial scale nowadays. It is a large scale operation carries out in 240 rows of plug reactors pass through a firebox of reformer. In this process the reformer is charged with a mixture of CH4 and H2O at a molar ratio from 1:3 to 1:4 [i ]. Although the lower steam ratio could lead to a reduction in the life of reformer tubes, economical restrictions limit our choices. At unit 17 this ratio is 1:3.9[ii ]. Reforming is extremely endothermic and is limited by thermodynamic equilibrium. Therefore, the development of catalysts significantly increases the conversion of the reactants. The temperature inside the tubes evolves from 545 C to 860 C and the pressure alters from 24.5 barg to 21.30 barg. Indeed precise simulations requires vast knowledge and kinetic and thermodynamic data of the process, Xu & Froment [iii ,iv ] derived twenty-one sets of three rate equations. They also examined the accuracy and credit of Langmuir-Hinshelwood rate equation. Rostrup-Nielsen [v ] investigated the details of catalysts deactivation and principles of coke formation. While taking all of the aspects of the above methods into account leads us to theoretical insolvable PDEs, logical results will appear by making reasonable assumptions to reduce the complexity of intrinsic kinetics. In order to halt cracking reaction and to prevent the formation of hot spot and overheated bands, reformer is charged with two types of catalysts [i ,ii ]. While their shapes are similar and their density values are close, their MS student at Sharif University of technology, International campus (E-Mail: arash.nasiris@gmail.com). MS student at Babol Noshirvani University of Technology (E-Mail: ahosseini.che.eng@gmail.com). composition differs. One is a Nickel oxide dispersed on calcium aluminate ceramic support catalyst, Katalco57-4Q and the other is a Nickel oxide dispersed on calcium aluminate ceramic support promoted by alkali catalyst, Katalco25-4Q. II. A BRIEF DESCRIPTION OF THE REFORMER SECTION The treated feed vapor is mixed with a controlled quantity of superheated process steam and is preheated in the convection section of the reformer furnace (firebox). The hot feed plus steam mixture is distributed over the catalyst tubes of the reformer, where the hydrocarbons in the feed gas are converted to hydrogen, carbon monoxide and carbon dioxide in the presence of steam over the nickel catalyst (Ni/Ca-Al2O3). The produced gas, which leaves the reformer, is essentially a mixture of hydrogen, carbon monoxide, carbon dioxide, methane and steam, which will make the equilibrium to approach a certain degree. The corresponding (dry gas) methane concentration is generally referred to as methane slip at the outlet temperature. (1). CnHm+nH2O nCO + (n+m/2) H2 (2)WGS reaction CO+H2O CO2 +H2 Besides the above-mentioned reactions there are a number of side reactions, which are not desirable: (3)Boudouard reaction: 2CO C+CO2 (4)CO reduction: CO+H2 C+H2O (5)Methane Cracking CH4 C+2H2 These reactions are suppressed by applying an excess of steam, so that eventually formed carbon will be removed by the reversed CO reduction reaction (4). A low steam to carbon ratio can lead to carbon deposition and thus catalyst damage. The overall heat effect of the steam reforming reactions is strongly endothermic i.e. heat has to be supplied externally to achieve the required conversion. This heat is provided by combustion of PSA purge gas as the priority fuel and natural gas as the make-up fuel. III. MATHEMATICAL MODEL AND EQUATIONS DERIVATION I. ASSUMPTIONS AND VARIABLES The following assumptions were considered: 1. Unsteady State operation Simultaneous, nonlinear, 2D modeling of tubular reactor of CH4 reforming in H2 production unit of Shazand Arak refinery Arash Nasiri Savadkouhi, Seyed Adel Hoseini M
2. Non-iso thermal condition. 3. Plug flow on reaction zone. 4. Reaction progress in circumferential direction is ignored. Since the maximum ratio of the lateral perimeter to height is 2() (5.65)/ (1200) =0.0295 5. All the reactants are gaseous and follow the ideal gas equation. Since CH4 is non-polar and the temperature is also high enough. 6. The First law of Fick is utilized with the constant diffusion coefficient. 7. The conduction heat rates evaluated from Fouriers law 8. Heat transfer by radiation is neglected. Since the temperature inside is less than 1000 C [vi ] 9. Pressure drop through the tubes is ignored. While Preal=3.2 bar 10. We considered that the reformer is loaded by a uniform catalyst. 11. WGS reaction is ignored. 12. Although the true shape of catalyst is shown below to simplify the calculations we consider one catalyst as four simple cylindrical catalysts. Fig 1. True shape of catalyst and the assumed shape. II.Governing Equations We start with the derivation of the PDEs. The physical system is a convectiondiffusionreaction system, for example, a tubular reactor, which is modeled in cylindrical coordinates, (r, z) (we assume angular symmetry so that the third cylindrical coordinate, 慮, is not required). The PDE system has two dependent variables, ca and Tk, which physically are the reacting fluid concentration, ca, and the fluid temperature, Tk. Thus we require two PDEs. These two dependent variables are functions of three independent variables: r, the radial coordinate; z, the axial coordinate; and t, time, as illustrated by Figure 2. The solution will be ca(r, z, t) and Tk(r, z, t) in numerical form as a function of r, z, t. Fig 2.Represention of a convection-diffusion-reaction system in a plug reactor Fig 3. Precise view of the reactor charged with catalysts. The mass balance for an incremental element of length z is: (1) Division by 2rrz and minor rearrangement gives: (2) Or in the limit r 0, z 0, We may briefly discuss the physical meaning of each term in this PDE. For the differential volume in r and z, As mentioned before we now assume Ficks first law for the flux with a mass diffusivity, Dc (1) We obtain Or expanding the radial group, (2) Equation (2) is the required material balance for ca(r, z, t) and would be sufficient for the calculation of ca(r, z, t)
except we consider the case of an exothermic reaction that liberates heat so that system is not isothermal. That is, the reacting fluid temperature varies as a function of r, z, and t so that a second PDE is required, an energy balance, for the computation of the temperature Tk (r, z, t). Further, Eq. (2) is coupled to the temperature through the reaction rate coefficient kr as it changes dramatically by changes in temperature. Also note that the volumetric rate of reaction in Eq. (2),kr.ca 2 , is second order; that is, the rate of reaction is proportional to the square of the concentration ca .This kinetic term is therefore a source of non-linearity in the model. The thermal energy balance for the incremental section is Divided by 2r.r.z. gives Or in the limit r 0, z 0, Briefly, the physical meaning of each term in this PDE for the differential volume in r and z is: Also, there is one technical detail about this energy balance which we may consider. The derivative in t is based on the specific heat Cp reflecting the fluid enthalpy. More generally, an energy balance is based on the specific heat Cv reflecting the fluid internal energy. For a liquid, for which pressure effects are negligible, the two are essentially the same, although our system is gaseous, we assumed that P is constant therefore we use Cp instead. If we now assume Fouriers first law for the flux with a thermal conductivity, k, (3) We obtain (after its division by Cp) Or expanding the radial group, (4) Where Dt =k/Cp. Equation (4) is the required energy balance for T(r, z, t). The reaction rate constant, kr, is given by (5) Note from Eq. (5) that kr is also nonlinearly related to the temperature Tk as eE/RTk . This Arrhenius temperature dependency is a strong nonlinearity, and explains why chemical reaction rates are generally strongly dependent on temperature; of course, this depends, to some extent, on the magnitude of the activation energy, E. The variables and parameters of equations (1)(5) and the subsequent auxiliary conditions are summarized in Table 1 (in cgs units). Equation (2) requires one initial condition (IC) in t (since it is the first order in t), two boundary conditions (BCs) in r (since it is the second order in r), and one BC in z (since it is the first order in z). (6) (7) (8) (9) Equation (4) also requires one IC in t, two BCs in r, and one BC in z. (10) (11) (12) (13) The solution to equations. (1) (13) gives ca(r, z, t), Tk(r, z,
t) as a function of the independent variables r, z, t. Now assuming steady state flow, the requisite equations for the model can be gained presuming that the reaction is the first order one with respect to CH4. One of the difficulties of simplifications is the presence of rate coefficient of reaction in equations since it follows LHHW model which is LangmuirHinshelwoodHougen Watson (LHHW) kinetic; it was proposed with nine parameters, where the surface decomposition of methane is assumed as the rate determining step (RDS), and while all other reaction steps are set as reversible [vii ]. Fortunately kinetic scheme of Hougan-Watson can be easily replaced by the first-order kinetics. Besides while k is generally considered as a function of temperature, here it is assumed to be constant. Although Agnelli et.al [viii ] disclosed three different values of k, KA, KB at various temperatures for diverse rate models, they stipulated that regression techniques can be utilized to yield values of k. They also approved the accuracy of the first order reactions with respect to CH4. The overall reactions considering the enthalpies are: CH4+H2O CO + 3 H2 H298 k=206 KJ/mol CO+H2O CO2 + H2 H298 k= -41 KJ/mol CH4+2H2O CO2 + 4 H2 H298 k=165 KJ/mol Rate equations utilized here are gained from the results of Xu & Froment researches. The mechanism of adsorption and desorption contains thirteen steps while only three of them limit the rate of reaction. (14) (15) (16) (17) Table 1. variables and parameters of the 2D PDE model, unsteady state Variable/Parameter Symbol Reactant concentration Ca (r, z, t) Temperature Tk (r, z, t) Reaction rate constant ki Time t Radial position r Axial position z Mass diffusion flux qm Energy conduction flux qh Reactants and Products Entering concentration Cae Entering temperature Tke Initial concentration Cao Initial temperature Tko average wall temperature Tw Reactor radius R0 Reactor length Z1 Linear fluid velocity Mass diffusivity Dc Thermal diffusivity Dt=k/Cp Fluid density Fluid specific heat Cp Heat of reaction (948 K) H Specific rate constant K0 Activation energy(948 K) E Gas constant R Thermal conductivity k Heat transfer coefficient h Table 1. continue Units/Value None gmol/cm3 None K None cm3 /(gmol.s) None s None cm None cm None gmol/(cm2 .s) None cal/(cm2 .s) None References CH4 H2 None Calculated based on [ii ] 0.3 0.3 None [ii ] 818.15K - 0 gmol/cm3 [ii ] 818.15K Calculated based on [ii ] 1123K [ii ] 5.65cm [ii ] 1200cm Calculated based on [ii ] cm2 /s 2.88 0.881 Calculated based on[ix ] cm2 /s 11.8 4379.4 [x ] None [ii ] 0.00676 g/cm3 [ii ] 0.62 cal/g.K [vi ] 53536 cal/(gmol) calculated based on[xi ] 4.798 cm3 /(gmol.s) [vi ] 57383.9 cal/(gmol) - 1.987 cal/(gmol.K) calculated based on[x ] 0.000408(cal.cm)/(s.cm2 .K) calculated based on[x ] 0300..(cal)/(s.cm2 .K)
(18) (19) IV. RESULTS Down below the results of solving the PDEs with using MATLAB for both the steady and unsteady flow is illustrated. Fig 3. T gradient at first second. Qnormal operation=8977 Fig 4. T gradient through the length of reformer Qnormal operation=8977 Table 2. variables and parameters of the steady state model Names Variables Methane A Water B Carbon Monoxide C Hydrogen D Entrance Temperature T Cross sectional area of reactor S Specific heat of component i Cpi Reaction rate expression rA Catalyst bulk density 痢b Heating inside reactor q Partial pressure of component i Pi Total pressure PT Molar flow rate component i ni Total molar flow rate nT Equilibrium constant K(T) Constant [viii ] KA Constant [viii ] KB Reactor length z Table 2. continue Values Units 1 Molar rate of x/molar rate of methane 3.51 0 .00983 818 Kelvin 380 cm2 None cal/(g mole. kelvin) None g mole/ (volume of catalyst.second) 0.860 gram catalyst/ volume of catalyst None cal/(second. cm) None KPa 2449.4 KPa None g mole/ sec None g mole/ sec 0.45 - None - None KPa-1 1200 cm
Fig 5. T gradient through the reformer Qmin=175 J.m-1.min-1 Fig 6. Products production rates Fig 7. Reactants consumption rates Fig 8.surface plot of methane concentration at t= 200s Fig 9. Concentration and temperature gradient at t =150s and t=200s Fig 9. Concentration and temperature changes at t=200s and t=150s Fig 10. Surface plot of methane temperature at t=200s
Fig 11. Concentration and temperature gradient at t =50s and t=100s V. RESULTS AND CONCLUSION We can perceive the following notes from the plots: Since the reaction is extremely endothermic, at the entrance of the reactor the temperature may decrease instantly but it will rise radically since the heat flux is extremely high. This can be understood by looking into figure 3. Figure 4 predicts the outlet temperature of the reactor. The model estimates its order value 917 C while the true outlet temperature is 860 C. This 57 C difference highlights the truth that ignoring WGS reaction, was a reasonable assumption since it only affects the models results from less than ten percent in comparison with real data. Figure 5 illustrates the outlet temperature of the reactor when the furnaces are working by its 2% of normal operational load. As it is illustrated, if the reformer works at this ratio, the outlet temperature will reach to 818 C which means that all the generated heat by furnace is used to proceed the reaction, while the reformer works normally by a ratio fifty one times larger in order to maximize the production. Moreover, the increase of heat flux helps halting WGS reaction. Figure 6 estimates the outlet methane to be approximately zero while based on unit PFD, the ratio of (rate of methane in outlet stream)/ (rate of methane in inlet stream) is 0.19, the difference however can be related to omitting the probable reactions in the first assumptions. Figures 8 and 10 show clearly the radial variation of concentration and temperature along the axial length of the reactor at time t =200, when a steady-state profile is going to be reached. The radial variation of concentration is relatively small however the temperature Tk (figure 10) shows dramatic increase through the length of reactor. There is a significant difference in the results of temperature gradient graphs based on two different steady state and transient flow models. While the steady state model predicts that the components may reach the final temperature of 980 C when it leaves reformer, the transient flow model predicts that the components reach the final temperature less than 50 seconds. The reason for this noticeable difference is related to the ignorance of convection mechanism of heat transfer in steady state flow. The same conclusion can be derived for concentration plots. We can assume that the concentration would be slightly different from those of being plotted. Since the term diffusivity may not always be appropriate (particularly for a high-Reynolds-number situation) in the sense that the radial dispersion in a reactor or similar physical system is probably not due solely to molecular diffusion (which is generally a small effect) as much as stream splitting and flow around the internal obstructions known as eddy diffusion. Thus the term dispersion coefficient might be more appropriate. Also, dispersion coefficients are frequently measured experimentally since they typically reflect the effects of complex internal flow patterns. Eventually it seems that the presented models can generate reasonable results. Although some magnificent reactions ignored on purpose, the results approved that through the reaction conditions, prepared by reforming, the main progressive reaction is methane reforming. REFERENCES i Handbook of Industrial Catalysts, FUNDAMENTAL AND APPLIED CATALYSIS. M. V. Twigg, M. S. Spencer, Johnson Matthey. ii Shazand Arak Refinery Expansion and Upgrading Project, Unit 17 operation Manual and its PFD. iii Xu J, Froment G.F. Methane Steam Reforming, Methanation and Water-Gas Shift: a. Intrinsic Kinetics. AIChE J. 1989;35(88) iv Xu J, Froment G.F. Methane Steam Reforming, M. and WGS: b. Diffusional limitations and Reactor Simulation. AIChE J.1989;35(97) v Rostrup-Nielsen J.R. Catalysis Science and Technoogy (Volume 5). Edited by Anderson J.R. and Boudart M. Springer-Verlag, 1984 vi Denis Tschumperle, Mark Nikolardot. Fiber cooling modelisation during draw using CFD vii Y.-J. Wu, Ping Li, J.-G. Yu, A.F. Cunha, A.E. Rodrigues, Sorption-enhanced steam reforming of ethanol on NiMgAl multifunctional materials: Experimental and numerical investigation, Chemical Engineering Journal, 2013, 231, 36 viii Agnelli, M.E, Ponzi, E,N, and Yeramian,A.A, 1987. Catalytic Deactivation on Methane Steam Reforming Catalysts. 2 Kinetic Study. ix W.J.Massman. A review of Molecular Diffusivities of H2O, CO2, CH4, CO, O3, SO2, NH3, N2O, NO and NO2 in Air, O2 and N2 near STP x Perrys Chemical Engineers Handbook xi D.W.Allen & E.R.Gerhard and M.R.Likins,Jr. Kinetics of the Methane-Steam Reaction.

More Related Content

Simultaneousnonlinear two dimensional modeling of tubular reactor of hydrogen production unit of Shazand Arak refinery

  • 1. Abstract This paper developed two mathematical models of a packed tubular reactor loaded with reforming catalysts which results in a simultaneous, nonlinear, two dimensional PDE in cylindrical coordinates. It generates 2-D radial and axial plots of components concentration and temperature. Here both regimes of steady state and transient flow are assumed. We also ignored the presence of WGS reaction and analyzed the consequences by comparing the results with experimental data. KeywordsCatalytic reactor, Reformer modeling, methane Steam reforming, Shazand Arak refinery, Hydrogen production unit. I. INTRODUCTION Ethan steam reforming is one of the most common methods to produce hydrogen at industrial scale nowadays. It is a large scale operation carries out in 240 rows of plug reactors pass through a firebox of reformer. In this process the reformer is charged with a mixture of CH4 and H2O at a molar ratio from 1:3 to 1:4 [i ]. Although the lower steam ratio could lead to a reduction in the life of reformer tubes, economical restrictions limit our choices. At unit 17 this ratio is 1:3.9[ii ]. Reforming is extremely endothermic and is limited by thermodynamic equilibrium. Therefore, the development of catalysts significantly increases the conversion of the reactants. The temperature inside the tubes evolves from 545 C to 860 C and the pressure alters from 24.5 barg to 21.30 barg. Indeed precise simulations requires vast knowledge and kinetic and thermodynamic data of the process, Xu & Froment [iii ,iv ] derived twenty-one sets of three rate equations. They also examined the accuracy and credit of Langmuir-Hinshelwood rate equation. Rostrup-Nielsen [v ] investigated the details of catalysts deactivation and principles of coke formation. While taking all of the aspects of the above methods into account leads us to theoretical insolvable PDEs, logical results will appear by making reasonable assumptions to reduce the complexity of intrinsic kinetics. In order to halt cracking reaction and to prevent the formation of hot spot and overheated bands, reformer is charged with two types of catalysts [i ,ii ]. While their shapes are similar and their density values are close, their MS student at Sharif University of technology, International campus (E-Mail: arash.nasiris@gmail.com). MS student at Babol Noshirvani University of Technology (E-Mail: ahosseini.che.eng@gmail.com). composition differs. One is a Nickel oxide dispersed on calcium aluminate ceramic support catalyst, Katalco57-4Q and the other is a Nickel oxide dispersed on calcium aluminate ceramic support promoted by alkali catalyst, Katalco25-4Q. II. A BRIEF DESCRIPTION OF THE REFORMER SECTION The treated feed vapor is mixed with a controlled quantity of superheated process steam and is preheated in the convection section of the reformer furnace (firebox). The hot feed plus steam mixture is distributed over the catalyst tubes of the reformer, where the hydrocarbons in the feed gas are converted to hydrogen, carbon monoxide and carbon dioxide in the presence of steam over the nickel catalyst (Ni/Ca-Al2O3). The produced gas, which leaves the reformer, is essentially a mixture of hydrogen, carbon monoxide, carbon dioxide, methane and steam, which will make the equilibrium to approach a certain degree. The corresponding (dry gas) methane concentration is generally referred to as methane slip at the outlet temperature. (1). CnHm+nH2O nCO + (n+m/2) H2 (2)WGS reaction CO+H2O CO2 +H2 Besides the above-mentioned reactions there are a number of side reactions, which are not desirable: (3)Boudouard reaction: 2CO C+CO2 (4)CO reduction: CO+H2 C+H2O (5)Methane Cracking CH4 C+2H2 These reactions are suppressed by applying an excess of steam, so that eventually formed carbon will be removed by the reversed CO reduction reaction (4). A low steam to carbon ratio can lead to carbon deposition and thus catalyst damage. The overall heat effect of the steam reforming reactions is strongly endothermic i.e. heat has to be supplied externally to achieve the required conversion. This heat is provided by combustion of PSA purge gas as the priority fuel and natural gas as the make-up fuel. III. MATHEMATICAL MODEL AND EQUATIONS DERIVATION I. ASSUMPTIONS AND VARIABLES The following assumptions were considered: 1. Unsteady State operation Simultaneous, nonlinear, 2D modeling of tubular reactor of CH4 reforming in H2 production unit of Shazand Arak refinery Arash Nasiri Savadkouhi, Seyed Adel Hoseini M
  • 2. 2. Non-iso thermal condition. 3. Plug flow on reaction zone. 4. Reaction progress in circumferential direction is ignored. Since the maximum ratio of the lateral perimeter to height is 2() (5.65)/ (1200) =0.0295 5. All the reactants are gaseous and follow the ideal gas equation. Since CH4 is non-polar and the temperature is also high enough. 6. The First law of Fick is utilized with the constant diffusion coefficient. 7. The conduction heat rates evaluated from Fouriers law 8. Heat transfer by radiation is neglected. Since the temperature inside is less than 1000 C [vi ] 9. Pressure drop through the tubes is ignored. While Preal=3.2 bar 10. We considered that the reformer is loaded by a uniform catalyst. 11. WGS reaction is ignored. 12. Although the true shape of catalyst is shown below to simplify the calculations we consider one catalyst as four simple cylindrical catalysts. Fig 1. True shape of catalyst and the assumed shape. II.Governing Equations We start with the derivation of the PDEs. The physical system is a convectiondiffusionreaction system, for example, a tubular reactor, which is modeled in cylindrical coordinates, (r, z) (we assume angular symmetry so that the third cylindrical coordinate, 慮, is not required). The PDE system has two dependent variables, ca and Tk, which physically are the reacting fluid concentration, ca, and the fluid temperature, Tk. Thus we require two PDEs. These two dependent variables are functions of three independent variables: r, the radial coordinate; z, the axial coordinate; and t, time, as illustrated by Figure 2. The solution will be ca(r, z, t) and Tk(r, z, t) in numerical form as a function of r, z, t. Fig 2.Represention of a convection-diffusion-reaction system in a plug reactor Fig 3. Precise view of the reactor charged with catalysts. The mass balance for an incremental element of length z is: (1) Division by 2rrz and minor rearrangement gives: (2) Or in the limit r 0, z 0, We may briefly discuss the physical meaning of each term in this PDE. For the differential volume in r and z, As mentioned before we now assume Ficks first law for the flux with a mass diffusivity, Dc (1) We obtain Or expanding the radial group, (2) Equation (2) is the required material balance for ca(r, z, t) and would be sufficient for the calculation of ca(r, z, t)
  • 3. except we consider the case of an exothermic reaction that liberates heat so that system is not isothermal. That is, the reacting fluid temperature varies as a function of r, z, and t so that a second PDE is required, an energy balance, for the computation of the temperature Tk (r, z, t). Further, Eq. (2) is coupled to the temperature through the reaction rate coefficient kr as it changes dramatically by changes in temperature. Also note that the volumetric rate of reaction in Eq. (2),kr.ca 2 , is second order; that is, the rate of reaction is proportional to the square of the concentration ca .This kinetic term is therefore a source of non-linearity in the model. The thermal energy balance for the incremental section is Divided by 2r.r.z. gives Or in the limit r 0, z 0, Briefly, the physical meaning of each term in this PDE for the differential volume in r and z is: Also, there is one technical detail about this energy balance which we may consider. The derivative in t is based on the specific heat Cp reflecting the fluid enthalpy. More generally, an energy balance is based on the specific heat Cv reflecting the fluid internal energy. For a liquid, for which pressure effects are negligible, the two are essentially the same, although our system is gaseous, we assumed that P is constant therefore we use Cp instead. If we now assume Fouriers first law for the flux with a thermal conductivity, k, (3) We obtain (after its division by Cp) Or expanding the radial group, (4) Where Dt =k/Cp. Equation (4) is the required energy balance for T(r, z, t). The reaction rate constant, kr, is given by (5) Note from Eq. (5) that kr is also nonlinearly related to the temperature Tk as eE/RTk . This Arrhenius temperature dependency is a strong nonlinearity, and explains why chemical reaction rates are generally strongly dependent on temperature; of course, this depends, to some extent, on the magnitude of the activation energy, E. The variables and parameters of equations (1)(5) and the subsequent auxiliary conditions are summarized in Table 1 (in cgs units). Equation (2) requires one initial condition (IC) in t (since it is the first order in t), two boundary conditions (BCs) in r (since it is the second order in r), and one BC in z (since it is the first order in z). (6) (7) (8) (9) Equation (4) also requires one IC in t, two BCs in r, and one BC in z. (10) (11) (12) (13) The solution to equations. (1) (13) gives ca(r, z, t), Tk(r, z,
  • 4. t) as a function of the independent variables r, z, t. Now assuming steady state flow, the requisite equations for the model can be gained presuming that the reaction is the first order one with respect to CH4. One of the difficulties of simplifications is the presence of rate coefficient of reaction in equations since it follows LHHW model which is LangmuirHinshelwoodHougen Watson (LHHW) kinetic; it was proposed with nine parameters, where the surface decomposition of methane is assumed as the rate determining step (RDS), and while all other reaction steps are set as reversible [vii ]. Fortunately kinetic scheme of Hougan-Watson can be easily replaced by the first-order kinetics. Besides while k is generally considered as a function of temperature, here it is assumed to be constant. Although Agnelli et.al [viii ] disclosed three different values of k, KA, KB at various temperatures for diverse rate models, they stipulated that regression techniques can be utilized to yield values of k. They also approved the accuracy of the first order reactions with respect to CH4. The overall reactions considering the enthalpies are: CH4+H2O CO + 3 H2 H298 k=206 KJ/mol CO+H2O CO2 + H2 H298 k= -41 KJ/mol CH4+2H2O CO2 + 4 H2 H298 k=165 KJ/mol Rate equations utilized here are gained from the results of Xu & Froment researches. The mechanism of adsorption and desorption contains thirteen steps while only three of them limit the rate of reaction. (14) (15) (16) (17) Table 1. variables and parameters of the 2D PDE model, unsteady state Variable/Parameter Symbol Reactant concentration Ca (r, z, t) Temperature Tk (r, z, t) Reaction rate constant ki Time t Radial position r Axial position z Mass diffusion flux qm Energy conduction flux qh Reactants and Products Entering concentration Cae Entering temperature Tke Initial concentration Cao Initial temperature Tko average wall temperature Tw Reactor radius R0 Reactor length Z1 Linear fluid velocity Mass diffusivity Dc Thermal diffusivity Dt=k/Cp Fluid density Fluid specific heat Cp Heat of reaction (948 K) H Specific rate constant K0 Activation energy(948 K) E Gas constant R Thermal conductivity k Heat transfer coefficient h Table 1. continue Units/Value None gmol/cm3 None K None cm3 /(gmol.s) None s None cm None cm None gmol/(cm2 .s) None cal/(cm2 .s) None References CH4 H2 None Calculated based on [ii ] 0.3 0.3 None [ii ] 818.15K - 0 gmol/cm3 [ii ] 818.15K Calculated based on [ii ] 1123K [ii ] 5.65cm [ii ] 1200cm Calculated based on [ii ] cm2 /s 2.88 0.881 Calculated based on[ix ] cm2 /s 11.8 4379.4 [x ] None [ii ] 0.00676 g/cm3 [ii ] 0.62 cal/g.K [vi ] 53536 cal/(gmol) calculated based on[xi ] 4.798 cm3 /(gmol.s) [vi ] 57383.9 cal/(gmol) - 1.987 cal/(gmol.K) calculated based on[x ] 0.000408(cal.cm)/(s.cm2 .K) calculated based on[x ] 0300..(cal)/(s.cm2 .K)
  • 5. (18) (19) IV. RESULTS Down below the results of solving the PDEs with using MATLAB for both the steady and unsteady flow is illustrated. Fig 3. T gradient at first second. Qnormal operation=8977 Fig 4. T gradient through the length of reformer Qnormal operation=8977 Table 2. variables and parameters of the steady state model Names Variables Methane A Water B Carbon Monoxide C Hydrogen D Entrance Temperature T Cross sectional area of reactor S Specific heat of component i Cpi Reaction rate expression rA Catalyst bulk density 痢b Heating inside reactor q Partial pressure of component i Pi Total pressure PT Molar flow rate component i ni Total molar flow rate nT Equilibrium constant K(T) Constant [viii ] KA Constant [viii ] KB Reactor length z Table 2. continue Values Units 1 Molar rate of x/molar rate of methane 3.51 0 .00983 818 Kelvin 380 cm2 None cal/(g mole. kelvin) None g mole/ (volume of catalyst.second) 0.860 gram catalyst/ volume of catalyst None cal/(second. cm) None KPa 2449.4 KPa None g mole/ sec None g mole/ sec 0.45 - None - None KPa-1 1200 cm
  • 6. Fig 5. T gradient through the reformer Qmin=175 J.m-1.min-1 Fig 6. Products production rates Fig 7. Reactants consumption rates Fig 8.surface plot of methane concentration at t= 200s Fig 9. Concentration and temperature gradient at t =150s and t=200s Fig 9. Concentration and temperature changes at t=200s and t=150s Fig 10. Surface plot of methane temperature at t=200s
  • 7. Fig 11. Concentration and temperature gradient at t =50s and t=100s V. RESULTS AND CONCLUSION We can perceive the following notes from the plots: Since the reaction is extremely endothermic, at the entrance of the reactor the temperature may decrease instantly but it will rise radically since the heat flux is extremely high. This can be understood by looking into figure 3. Figure 4 predicts the outlet temperature of the reactor. The model estimates its order value 917 C while the true outlet temperature is 860 C. This 57 C difference highlights the truth that ignoring WGS reaction, was a reasonable assumption since it only affects the models results from less than ten percent in comparison with real data. Figure 5 illustrates the outlet temperature of the reactor when the furnaces are working by its 2% of normal operational load. As it is illustrated, if the reformer works at this ratio, the outlet temperature will reach to 818 C which means that all the generated heat by furnace is used to proceed the reaction, while the reformer works normally by a ratio fifty one times larger in order to maximize the production. Moreover, the increase of heat flux helps halting WGS reaction. Figure 6 estimates the outlet methane to be approximately zero while based on unit PFD, the ratio of (rate of methane in outlet stream)/ (rate of methane in inlet stream) is 0.19, the difference however can be related to omitting the probable reactions in the first assumptions. Figures 8 and 10 show clearly the radial variation of concentration and temperature along the axial length of the reactor at time t =200, when a steady-state profile is going to be reached. The radial variation of concentration is relatively small however the temperature Tk (figure 10) shows dramatic increase through the length of reactor. There is a significant difference in the results of temperature gradient graphs based on two different steady state and transient flow models. While the steady state model predicts that the components may reach the final temperature of 980 C when it leaves reformer, the transient flow model predicts that the components reach the final temperature less than 50 seconds. The reason for this noticeable difference is related to the ignorance of convection mechanism of heat transfer in steady state flow. The same conclusion can be derived for concentration plots. We can assume that the concentration would be slightly different from those of being plotted. Since the term diffusivity may not always be appropriate (particularly for a high-Reynolds-number situation) in the sense that the radial dispersion in a reactor or similar physical system is probably not due solely to molecular diffusion (which is generally a small effect) as much as stream splitting and flow around the internal obstructions known as eddy diffusion. Thus the term dispersion coefficient might be more appropriate. Also, dispersion coefficients are frequently measured experimentally since they typically reflect the effects of complex internal flow patterns. Eventually it seems that the presented models can generate reasonable results. Although some magnificent reactions ignored on purpose, the results approved that through the reaction conditions, prepared by reforming, the main progressive reaction is methane reforming. REFERENCES i Handbook of Industrial Catalysts, FUNDAMENTAL AND APPLIED CATALYSIS. M. V. Twigg, M. S. Spencer, Johnson Matthey. ii Shazand Arak Refinery Expansion and Upgrading Project, Unit 17 operation Manual and its PFD. iii Xu J, Froment G.F. Methane Steam Reforming, Methanation and Water-Gas Shift: a. Intrinsic Kinetics. AIChE J. 1989;35(88) iv Xu J, Froment G.F. Methane Steam Reforming, M. and WGS: b. Diffusional limitations and Reactor Simulation. AIChE J.1989;35(97) v Rostrup-Nielsen J.R. Catalysis Science and Technoogy (Volume 5). Edited by Anderson J.R. and Boudart M. Springer-Verlag, 1984 vi Denis Tschumperle, Mark Nikolardot. Fiber cooling modelisation during draw using CFD vii Y.-J. Wu, Ping Li, J.-G. Yu, A.F. Cunha, A.E. Rodrigues, Sorption-enhanced steam reforming of ethanol on NiMgAl multifunctional materials: Experimental and numerical investigation, Chemical Engineering Journal, 2013, 231, 36 viii Agnelli, M.E, Ponzi, E,N, and Yeramian,A.A, 1987. Catalytic Deactivation on Methane Steam Reforming Catalysts. 2 Kinetic Study. ix W.J.Massman. A review of Molecular Diffusivities of H2O, CO2, CH4, CO, O3, SO2, NH3, N2O, NO and NO2 in Air, O2 and N2 near STP x Perrys Chemical Engineers Handbook xi D.W.Allen & E.R.Gerhard and M.R.Likins,Jr. Kinetics of the Methane-Steam Reaction.