Advanced Search
RSS LinkedIn Twitter

Journal Archive

Johnson Matthey Technol. Rev., 2022, 66, (2), 137


Discrete Simulation Model of Industrial Natural Gas Primary Reformer in Ammonia Production and Related Evaluation of the Catalyst Performance

Optimising catalyst performance and lifetime

PEER REVIEWED Received 4th March 2021; Revised 20th May 2021; Accepted 26th May 2021; Online 28th May 2021


Article Synopsis

The catalytic steam reforming process of natural gas consumes up to approximately 60% of overall energy used in ammonia production. The optimisation of the reforming catalyst performance can significantly improve the operation of the whole ammonia plant. An online model uses actual process parameters to optimise and reconcile the data of primary reforming products with possibility to predict the catalyst performance. The model uses a combination of commercial simulator and open-source code based on scripts and functions in the form of m-files to calculate various physical properties of reacting gases. The optimisation of steady-state flowsheet, based on real-time plant data from the distributed control system (DCS), is essential for the application of the model at the industrial level. The simplicity of the calculation method used by the model provides the fundamental basis for industrial application in the frame of digitalisation initiative. The principal aim of the optimisation procedure is to change the working curve for methane regarding its equilibrium curve as well as methane outlet molar concentration. This is the critical process parameter in reforming catalyst operation. An industrial top fired primary reformer unit based on Kellogg Inc technology design served for the validation of the model. Calculation procedure is used for continuous online evaluation of the most commercially available primary reformer catalysts. Based on the conducted evaluation, the model can indicate possible recommendations which can mitigate marginal performance and prolong reformer catalyst lifetime.

1. Introduction

In ammonia production, approximately 60% of overall energy consumption relates to the front end of the production process, namely steam methane reforming (SMR) (1). The proper operation of the SMR unit is the primary focus for operators who aim to minimise costs in the whole ammonia plant. Primary reforming is a process in which gases containing hydrocarbon react with steam to generate a product which is a gas with, as much as possible, higher hydrogen content. The reacting conditions are such that the remaining product components also include methane, carbon monoxide, carbon dioxide and some inerts (nitrogen, argon and helium) that may have been present in the feedstock. The basic occurring reactions in the process are highly endothermic (2). As a result, the process operates by passing the mixture of hydrocarbon and steam reactants over a reformer catalyst in a fired furnace.

Regarding such complex equipment, it is important to consider the type of furnace used to transfer heat to the reactants, the catalyst properties such as the activity, life, size and strength together with other operating parameters such as feedstock characteristics, pressure, temperature, volume flows and the desired product composition.

Latham (3) developed a mathematical model of the SMR unit for use in process performance simulations and online monitoring of tube-wall temperatures using a plug flow pattern. According to the model inputs, it is able to calculate temperature profiles for the outer-tube wall, inner-tube wall, furnace gas and process gas. However, it was concluded that to make the current model usable by plant operators, an interface between the plant inputs and the model needs to be created and the model runtime of 4 min may need to be improved. Computational fluid dynamics (CFD) modelling and simulation of SMR reactors and furnaces has been deeply elaborated by Aguirre (4) and the same are applicable for both pilot-scale and full industrial scale furnaces. The author designed the workflow to be executed without the need of an expert user, to be deployed in cloud environment and to be fully or partially used. Model convergence is determined by a difference of standard deviation below 3.0%. In spite of successful trials, a total time of approximately 30 h was required for the simulation and optimisation to achieve convergence. The author suggested further study in the modelling process to speed up the CFD calculation by implementation of the smart-determination of variable numerical computation parameters which could be implemented by the different optimisation schemes. Moreover, Lao et al. (5) demonstrated that CFD software can be employed to create a detailed CFD model of an industrial-scale SMR tube using plug flow. It showed good approximation of the catalyst with the available industrial plant data. Simulation results were very close to the plant data for temperature and species composition. The only drawback of the proposed model was computational time which was at the level of 5 min with the steady state solver. It can be concluded that the CFD modelling technique provides high accuracy, but takes a long time to converge which is impractical for industrial applications. To overcome this drawback Holt et al. (6) created a SMR model in Python based on previous work by Xu and Froment (7, 8). The model was shown to be a reasonable replica of the original work of Xu and Froment (7, 8), however they did not regress the model against a real SMR unit.

The general arrangement of the SMR unit comprises a reformer furnace, reformer tubes and a related catalyst. The general furnace classification according to firing pattern refers to top, side and bottom fired furnace. In all the furnaces, the catalyst is contained in heat resistant alloy reformer tubes, the inside diameter of which ranges typically from 6.0 cm to 20.0 cm and the wall thickness of which is from 0.9 cm to 1.9 cm. Fired lengths in commercial furnaces vary approximately from 2.5 m to 14 m. The most commonly encountered fired lengths, however, are 9.0 m to 12.0 m (1, 2). Firing is usually controlled in such a way that the tube wall temperature is maintained at values that will give reasonable tube lifetime. Different firing control strategies for SMR units have been developed over time using standard and advanced process control approaches based on proportional-integral-derivative controllers or model predictive control techniques (912). Each of them showed different advantages in implementation of control structures to describe the dynamic relationship between the reformer tube wall temperature (process outputs) and manipulated variables and disturbances (process inputs). The main focus of different control strategies is to keep the reformer tube wall temperature at a safe temperature level to protect the reformer tube wall material against mechanical degradation. By design and industry practice, maximum allowable tube wall temperatures will give in-service life of 100,000 h when considering the stress-to-rupture and creep damage properties of the particular alloy used in manufacturing the reformer tube (13).

The discrete model was validated with Kellogg Inc top fired furnace which is characterised by having the burners on the top and firing down. The reformer tubes in such a furnace are installed in parallel rows with the burners between each row.

Catalyst performance has an important effect on hydrocarbon conversion and on the reformer tube wall temperature, which is usually monitored by plant operators, and therefore subject to human error due to the lack of appropriate knowledge about correction actions. Primary reformer catalyst performance is usually estimated in terms of its methane approach to equilibrium (ATE), reformer furnace tube wall temperatures, pressure drop and the presence or absence of hot spots or bands on the reformer tubes. The most important variable for bringing the catalyst performance to the maximum activity is ATE. The ATE represents the difference between the actual value at the exit of the catalyst and the value at which the measured exit gas composition would be at equilibrium (2). The actual methane ATE cannot theoretically be less than zero (a negative number).

Over the years, several applications were developed to make evaluations and recommendations regarding catalyst performance (56, 1416). However, according to literature findings, none of them are designed to continuously (online) evaluate the catalyst performance during SMR operation. Most design calculations involve simply checking designs submitted by the major catalyst suppliers. These submitted requests give specific material balances (both at the inlet and at the exit of the reformer tubes), operating conditions (pressures, temperatures) and furnace configurations (number of tubes and tube dimensions). These data are then entered into the proprietary applications and the results are checked for methane ATE, heat flux versus catalyst size-type and pressure drop. After performing evaluation, the catalyst suppliers submit the report to the operators, which is then used for eventual corrective actions. In most cases the evaluation reports with related recommendations for corrective actions are submitted with a huge time delay and cannot be immediately applied to adjust catalyst performance.

In order to overcome these limitations and bring novelty in this research, the primary aim of this work is the delivery of a sophisticated online model which predicts the performance of reformer catalysts having specific design. The model can simulate the tube side process and provides a detailed profile of the reaction rates, molar gas composition, actual conversion of hydrocarbon feedstock, the pressure and temperature profiles inside of the tubes incrementally. The model can continuously receive real-time plant data from any commercial DCS, which is then reconciled against the model and subsequently used to generate recommendations for the operators. The supplementary novelty in the proposed model is the simple and reliable calculation method with very fast computational routine, which brings to the operator sufficiently understandable recommendations to evaluate the catalyst performance and carry out necessary remediation measures.

Furthermore, additional novelty is the development of appropriate shared memory coupling for communication between steady-state model and MATLAB® via a shared memory area on the host system. The developed communication system with related open structure and user friendly interface enables implementation of the proposed solution to any DCS system which enables time savings in catalyst performance evaluation.

In the frame of a digitalisation initiative and according to the goals of Industry 4.0, the proposed model can bring additional innovative benefits to SMR units to improve productivity and uptime.

A commercial simulator (UniSim® Design R470, Honeywell Inc, USA) was used to build the steady-state flowsheet of the primary reformer furnace, while MATLAB® was used to reconcile the process data from the simulator. MATLAB® uses ‘actxserver’ to create a component object model (COM) automation server that can control the simulator. Through the COM interface, the simulator flowsheet can be opened, data written to and read, saved and closed. The COM interface establishes a two-way communication between the simulator and MATLAB® through shared memory which was built as level-2 S-function. Figure 1 shows the communication structure between steady-state flowsheet in UniSim®, MATLAB® and DCS.

Fig. 1.

The communication structure for online monitoring and optimisation of reformer catalyst

The communication structure for online monitoring and optimisation of reformer catalyst

2. Model Development

2.1. Process Description

The process described herein is based on the Kellogg Inc catalytic high-pressure reforming method for producing ammonia starting with natural gas feed. An ammonia plant steam reforming unit can produce 1360 tonnes per day of liquid ammonia. Figure 2 presents the steady-state flow sheet of the SMR unit build in UniSim® Design R470 with the main process flow designated with the red line.

Fig. 2.

SMR steady-state flowsheet

SMR steady-state flowsheet

Natural gas feed at a pressure of about 32 bar enters the natural gas knock-out drum 120-F for elimination of entrained liquid. The outlet line of 120-F feeds the one-stage centrifugal natural gas feed compressor 102-J driven by back-pressure (40/4 bar) steam turbine 102-JT. Outlet pressure of natural gas is at the level of 42 bar. Hydrogen required for desulfurisation of the natural gas is injected into the paralleled natural gas stream entering the natural gas fired heater 103-B. The outlet temperature of 103-B is 400°C. The heated natural gas stream flows through two reactors in series. The first is the hydrogenator 101-D, which contains a single bed of cobalt-molybdenum catalyst. It converts the organic sulfur compounds to hydrogen sulfide in the presence of the hydrogen injected upstream of 103-B. The natural gas stream next passes into the desulfuriser reactor 102-D, which contains a single bed of zinc-oxide catalyst. In this reactor the hydrogen sulfide is converted to zinc sulfide which remains in the catalyst.

The desulfurised natural gas, plus residual hydrogen, leaves 102-D with a sulfur content of 0.25 ppm and at a temperature of 370°C. The natural gas plus residual hydrogen stream is joined by the process steam in the mixer. The process steam is at a pressure of about 40 bar and a temperature of 392°C. The steam flow is controlled with the steam-to-natural gas (S/NG) molar ratio controller.

The SMR feed gas flows to the mixed feed coil, which is located in the convection section of the SMR furnace. In this coil, the SMR feed is heated to about 510°C. After heating, the SMR feed flows down through ten rows of reformer tubes that are suspended in the radiant box of primary reformer 101-B. Eleven rows of forced draught down fired burners are located in parallel rows to the catalyst tubes, in total 198 burners. They raise the feed temperature to about 790°C at the outlet of the catalyst tubes. In addition, 11 tunnel burners are used to heat the waste gases passing from the radiant to the convection part of the SMR furnace. 520 catalyst tubes with a total length of 10 m and inside diameter of 0.0857 m contain 30 m3 of nickel reformer catalyst. The reformed gas (syngas) then flows to the secondary reformer for further processing.

2.2. Steam Reforming Model

In order to predict the performance of the SMR process, it is necessary to simulate the tube side process and provide a detailed profile of the heat flux, gas composition, carbon forming potential and the pressure inside of the reformer tubes incrementally. The calculations involve solving material and energy balance equations along with reaction kinetic expressions for the nickel catalyst.

The general overall reaction for the steam reforming of any hydrocarbons can be defined as Equation (i) (1, 2):



In this work, steam reforming of the natural gas is described with the following equations, as the methane is the major constituent of the natural gas. Equation (ii) (1, 2):



In parallel with this SMR equilibrium, the water gas shift (WGS) reaction proceeds according to Equation (iii) (1, 2):



Minette et al. (17) in their work stated that the second SMR reaction is often not accounted for assuming it follows directly from combining Equations (i) and (ii). However, the work of Xu and Froment (78) showed that the second SMR reaction expressed by Equation (iv) follows an independent reaction path and must be accounted for in combination with Equations (i) and (ii), as confirmed by the measurements of Minette et al. (17):



As mentioned, the described reactions proceed in indirectly heated reformer tubes filled with nickel-containing reforming catalyst and are controlled to achieve only a partial methane conversion. In a top fired reformer usually up to 65% to 68% conversion based on methane feed can be accomplished, leaving around 10 mol% to 14 mol% methane per dry basis (1, 2).

The overall SMR reaction of methane is endothermic and proceeds with an increase of volume at the elevated pressure of 20 bar to 40 bar and temperatures from 800°C to 1200°C at the exit of the reformer tubes in the presence of metallic nickel catalyst as an active component. Besides pressure and temperature, the S/NG molar ratio has a beneficial effect on the equilibrium methane concentration (18).

Another reason for applying the appropriate (higher) S/NG molar ratio is to prevent carbon deposition on the reforming catalyst. The side effect of carbon deposition is a higher pressure drop and the reduction of catalyst activity. As the rate of endothermic reaction is lowered, this can cause local overheating of the reformer tubes (hot spots and bands) and the premature failure of the tube walls. The carbon formation may occur via Boudouard reaction, methane cracking and carbon monoxide and carbon dioxide reduction. These reactions are reversible with dynamic equilibrium between carbon formation and removal. Under typical steam reforming conditions, Boudouard reaction and carbon monoxide and carbon dioxide reduction cause carbon removal, whilst methane cracking leads to carbon formation in the upper part of the reformer tube (19). Greenfield SMR units based on natural gas regularly use a S/NG molar ratio of around 3.0, while older installations are in the range from 3.5 to 4.0 (1). From the theoretical point of view any S/NG molar ratio which is slightly over 1.0 will prevent cracking, because the rate of carbon removing reactions is faster than the rate of carbon deposition reactions. However, from the practical point of view (catalyst limitations and sufficient quantity of steam for the downstream process step of WGS conversion), the minimum molar ratio which applies at the industrial level is 2.5. To account for all these facts, the model was validated for S/NG molar ratios in the range from 2.0 to 6.0.

The nickel content in relation to the composition and structure of the support differs considerably from one catalyst supplier to another. This is the reason why it is difficult to relate data from industrial plants to laboratory experiments. Reformer simulations frequently use a numerical approach in which the experimental data serves for reaction rate calculations which are described by closed analytical expressions. From the reaction rates perspective, it is possible to calculate the equilibrium gas composition for a given pressure and S/NG molar ratio at different temperatures. On top of this, the equilibrium curve which is defined by the corresponding enthalpy changes versus temperature also presents a useful parameter in the estimation of the catalyst performance. The comparison of the mentioned equilibrium curves with the working curves (working point) and the subsequent operator’s adjustments of the influencing process parameters according to the evaluated recommendations seem a useful tool to improve the catalyst performance.

In order to describe the kinetic conditions which are necessary for the determination of equilibrium methane molar concentration (a measure for the theoretical conversion) and enthalpy change over different nickel catalysts in relation with temperature at different S/NG molar ratios and reforming pressures, the model uses the following reaction rates for the equilibrium Equations (ii) to (iv) (78, 20), Equations (v)(viii):









where r presents the reaction rates for methane, carbon monoxide and carbon dioxide in kmol m–3 s–1; p stands for the species partial pressures (in atm); T is the temperature (in K); while R is the gas constant (in kJ kmol–1 K–1).

Kinetic rate constants ki are given by the general Arrhenius relationship, Equation (ix) (78, 20), where i denotes the number of reactions from Equation (i) to (iii):


The units of k2 and k4 (Equation (ii) and (iv)) are kmol bar0.5 kg–1cat h–1), while the unit of k3 (Equation (iii) is kmol bar–1 kg–1cat h–1).

Table I (20) gives the parameters for the activation energies, Ei, and for the pre-exponential factors, Ai, used in the model, valid for most of the commercial nickel catalysts with either MgAl2O4 or CaAl12O19 support.

Table I

Parameters for the Activation Energies, E i, and for the Pre-Exponential Factors, Ai

Equilibrium reaction Activation energy, Ei
Pre-exponential factor, Ai
Unit Value Unit Value
Reaction no. 2 kJ mol–1 –240.100 kmol bar0.5 kg–1cat h–1 4.22 × 1015
Reaction no. 3 kJ mol–1 –67.130 kmol bar–1 kg–1cat h–1 1.96 × 106
Reaction no. 4 kJ mol–1 –243.900 kmol bar0.5 kg–1cat h–1 1.02 × 1015

Apparent adsorption equilibrium constants Ki in Equation (x) are defined by the general expression given in (78, 20), where i denotes the species in Equations (i), (ii) and (iii) or methane, water, hydrogen and carbon monoxide:


Bi is the pre-exponential factor expressed in bar-1 or unitless, while ΔHi is the absorption enthalpy change expressed in kJ mol–1.

Table II presents the pre-exponential factors and the absorption enthalpy changes for species given in Equation (x), and the same is also valid for most of the commercial nickel catalysts with either MgAl2O4 or CaAl12O19 support.

From Equations (v) to (vii) it can be concluded that the concentration of hydrogen cannot be zero, because dividing with zero would make calculated reaction rates infinite. So, according to this, it is necessary to ensure the minimum content of hydrogen in the natural gas stream to ensure applicability of these equations in the model. From the process side, hydrogen is necessary for two reasons. Firstly, it is important for the removal of organic sulfur compounds present in the natural gas by the cobalt-molybdenum catalyst, as sulfur is a poison for the nickel catalyst (reaction between organic sulfur compounds and hydrogen to give hydrogen sulfide which is subsequently absorbed by zinc oxide bed). Secondly, hydrogen will always keep the nickel catalyst in the reduced state of metallic nickel and hence maintain adequate catalyst activity in the reformer tubes.

Table II

Parameters for the Pre-Exponential Factor, Bi, and for the Absorption Enthalpy Changes ΔHi

Species Pre-exponential factor, Bi
Absorption enthalpy change, ΔHi
Unit Value Unit Value
Methane bar–1 6.65 × 10–4 kJ mol–1 38.280
Water 1.77 × 105 kJ mol–1 –88.680
Hydrogen bar–1 6.12 × 10–9 kJ mol–1 82.900
Carbon monoxide bar–1 8.23 × 10–5 kJ mol–1 70.650

From the general stoichiometry and according to defined reaction rates, the model can calculate the molar flow rates of species i in kmol h–1 in the presence of an adequate quantity of nickel catalyst with the ultimate result of methane and water conversions. The relations used to determine the methane and water conversions are as follows (21, 22), Equations (xi)(xii):


A denotes the catalyst tube cross-sectional area in m2; ρB represents the catalyst bed density in kg m–3; Fi is the molar flow rate of the species methane and water in kmol h–1; while ηi is the effectiveness factor for methane and water.

To account for the variations in reaction rate throughout the catalyst pellet, a parameter called effectiveness factor, η, is defined. This is the ratio of the overall reaction rate in the catalyst pellet and the reaction rate at the external surface of the catalyst pellet. Effectiveness factor is a function of Thiele modulus, Φ, which is related to the catalyst volume and the external surface area of the catalyst pellets. Taking into account reaction rates given by Equations (v)(vii) and following the mechanism given by Xu and Froment (7, 8), the effectiveness factor can be calculated from Equation (xiii):


where p is the partial pressure of the species in bar; r presents the reaction rates for methane, carbon monoxide and carbon dioxide in kmol m–3 s–1; while ξ is the dimensionless intracatalyst coordinate.

Effectiveness factor profiles along the length of the reformer tube are calculated for all key species given in Equations (ii) to (iv) by solving two-point boundary differential equations for the catalyst pellets with the help of scripts and functions in the form of m-files, which was reconciled with the data from the simulator flowsheet.

The algorithm uses the following relationship for calculation of species concentration profiles inside the catalyst layer under reconciled conditions (17), Equations (xiv)(xv):




with the corresponding boundary conditions, Equations (xvi)(xvii):




where ξ is the dimensionless intracatalyst coordinate; De,A is the species effective diffusivity in m3fluid m–1catalyst s–1; p denotes the partial pressure of species in bar; R is the universal gas constant in kJ kmol–1 K–1; T is the bulk fluid temperature in K; h is catalytic layer thickness in m and ρs is the active solid density in kgcatalyst m–3catalyst.

The interfacial (gas-solid) mass and heat transfer limitations are negligible and were not accounted for, because the high volume flow velocity and sufficient turbulence have been assumed which reflects the operation conditions inside of the reformer tubes.

Due to model simplification and minimisation of the computational time the simplest geometry of a slab of catalyst has been assumed, which is a satisfactory assumption for the computational routine required for industrial application. The model has been tested with coating thickness in the range from 10 μm to 50 μm and the best fit with the actual process data was achieved with the catalyst coating of 10 μm.

The species effective diffusivity is determined by Equation (xviii):


where ɛs is the internal void fraction or porosity of the catalyst in m3fluid m–3catalyst; τ denotes the catalyst tortuosity and is the average diffusivity of species A.

The average diffusivity of species is determined by Equation (xix):


where DA is the diffusivity of the reacting species A given by Equation (xx) and S(rp,i) is the void fraction taken by the pores with radii ranging from rp,i to rp,i +1:


where DkA is the Knudsen diffusivity in m3fluid m–1catalyst s–1.

In order to have an appropriate computational speed of effectiveness factor (which is performed by m-file), the actxserver command is used for the interconnection through the COM automation server that controls the simulator. The COM interface establishes a two-way communication between the simulator and MATLAB® through shared memory block, which is built as level-2S-function. The approximation of the catalyst effectiveness factor is determined by correlating the kinetic model results with the plant process data, and the model is validated to get maximum alignment with the actual process data.

Conversions of methane and water are calculated by Equations (xxi)(xxii) (22):





The Ergun equation for the determination of the pressure drop across the plug flow reactor (PFR) is used and solved as an ordinary differential equation (2331), Equation (xxiii):


where ρ denotes the pressure in bar; ρ is the fluid density in kg m–3; v is the fluid velocity in m s–1; dp is the catalyst particle diameter in m; ∈ is the catalyst void fraction and Re is the particle Reynolds number.

The temperature variation of the reacting mixture (natural gas and steam) along the reformer tube is calculated according to the following relationship, Equation (xxiv):


where G is the reacting mixture flow rate in kg h–1; denotes average specific heat of the gas mixture in kJ kg–1 K–1; U is the overall heat transfer coefficient between the reformer tubes and their surrounding in m2 h K kJ–1; Tt,0 is the temperature of the furnace that surrounds the reformer tubes; ΔHi is the enthalpy change in kJ kmol–1; ρB represents the catalyst bed density in kg m–3; ηi is the effectiveness factor for each of the species in reacting mixture and ri is the reaction rates in kmol m–3 s–1.

The reformer catalyst tubes are simulated as PFR in which the flow field is modelled as plug flow, implying that the stream is radially isotropic (without mass or energy gradients). According to this, axial mixing is negligible. As the reactants flow the length of the reformer tube, they are continually consumed, hence, there is an axial variation in the concentration. Since reaction rate is a function of concentration, the reaction rate varies axially. To get the solution for the PFR (axial profiles of compositions, temperature and so forth) the reformer tubes are divided into several sub-volumes. Within each sub-volume, the reaction rate is spatially uniform. A mole balance executes routine calculation procedure in each sub-volume j according to Equation (xxv) (28, 29):



Because the reaction rate is spatially uniform in each sub-volume, the third term reduces to rjdV and at steady state, the above expression reduces to Equation (xxvi):



The firing side (furnace combustion model) was simulated according to the previous work of Zečević and Bolf (32) which is able to calculate adiabatic and real flame temperatures, quality and quantity composition of the waste gases, according to the known composition of the fuel gas and inlet temperatures of fuel and combustion air, with possibility to control all critical process parameters by implementation of proposed gain-scheduled model predictive control.

The basic input requirements for the model are:

  1. Integration information: number of reformer tube segments, minimum step fraction, minimum step length

  2. Tube dimensions: total volume, length and internal diameter of the reformer tube, number of tubes, wall thickness

  3. Tube packing: void fraction

  4. Catalyst data: diameter, sphericity, solid density, solid heat capacity, number of holes, tortuosity, mean pore radius, catalyst characteristic length, catalyst support

  5. Inlet process composition: flow rate, natural gas composition, pressure, temperature

  6. Outside tube wall temperature: measured values

  7. Heat transfer coefficient

  8. Activity coefficient.

2.3. Steady-State Model Assumptions

The aim of the steady-state model is to perform steady-state energy and material balances. In the model development phase, it is important to ensure that the properties of pure species and mixtures are estimated appropriately and this is one of the most important steps that will affect the rest of simulation. The basis for the fluid modelling is the Peng-Robinson equation of state, because the software used (UniSim® Design R470) allows the thermodynamic properties of gaseous hydrocarbons to be reliably predicted. The steady-state model uses modular operations with a non-sequential algorithm. The benefit of this assumption is that the requested information is processed as soon as it is supplied and in parallel the results of any calculation are automatically propagated throughout the flowsheet, both forwards and backwards. This feature enables fast model response and ensures adequate valuable recommendations to the user. Material, energy and composition balances are considered at the same time. Pressure, flow, temperature and composition specifications are equal and uniform at any cross section of the catalyst bed in the reformer tubes. Equality and uniformity assumptions ensure adequate conditions to the model solver to bring simulation results with either specification which simplifies the calculation procedure and speeds up the calculation routine. Due to very large flow rates of process gas the axial diffusion of mass and heat is assumed to be negligible.

The reformer tubes are designed as heterogeneous PFR. The external mass and heat transfer between the catalyst pellet and the reacting gas is negligible, because it is assumed that the internal diffusion is limited (the internal effectiveness factor is small). The interparticle diffusional resistance is considerable because of the heterogenous nature of the catalytic reaction and ensuring appropriate contact between the bulk gas phase and the catalytic active centres. All reformer tubes within the furnace are identical, the overall performance of the reformer tubes is achieved by multiplying with the number of tubes. The reason for this is uniform heat flux pattern inside of the radiation box of SMR furnace.

The model converts all the higher hydrocarbons in natural gas to an equivalent methane before any reformer tube integration takes place to satisfy Xu and Froment’s kinetic mechanism (78). The carbon forming potential applies only to methane. Tube wall temperature readings data must be collected by an adequate measurement technique (for example, optic pyrometry), assumed to be accurate because all other measurement methods are practically unreliable due to high temperatures and fast deterioration of the thermocouples. SMR and WGS reactions are considered as the only reactions in the model to be able to follow the proposed Xu and Froment kinetics (78). Methane reforming reaction with carbon dioxide (dry reforming) is neglected, because of very low intrinsic activity, unfavourable thermodynamics, minimum amounts of ‘pure’ carbon dioxide in the natural gas feed flow and avoidance of carbon deposition which is in the model prevented with higher S/NG molar ratios (33).

2.4. Steady-State Model Limitations

The model is not appropriate or does not give meaningful predictions regarding:

  1. Usage of feedstocks with higher hydrocarbon content such as butane or naphtha. As noted earlier, the model converts these compounds to an equivalent methane before any integrations take place

  2. The model cannot predict tube wall temperatures and these must be assumed or measured in the plant

  3. The model cannot provide any information about hot spot or hot band problems in the furnace.

3. Results and Discussion

The first step of this research is to simulate a chemical reaction system of the SMR operation in an ammonia plant by an appropriate commercial simulator and by application of scripts and functions in the accompanying m-file. The next step is the validation of the simulation results with the literature and process data from a real plant which uses commercially available reforming catalyst. The final and major step is to ensure adequate recommendations for adjustment of the process parameters during SMR operation which will bring the working curve of reforming catalyst performance much closer to the equilibrium curve (working point).

In order to perform model testing, the process parameters from the reference case (ammonia plant in Kutina, Croatia, owned by Petrokemija Plc) were used for the model validation. Table III presents the necessary input variables. This is the basis for the calculation of the equilibrium curves at different S/NG molar ratios, outlet temperatures, pressures and heat load. One part of the process information must be manually entered, while the other part is automatically fed from the DSC system (Table III).

Table III

Input Data for Model Validation

Variable Unit Value

Manual entry
Heated length of reformer tubes m 10.0
Inside diameter of reformer tubes m 0.085
Outside diameter of reformer tubes m 0.110
Number of reformer tubes pcs. 520
Number of rows pcs. 10
Number of arch/tunnel burners pcs. 198/11
Design heat liberation per arch/tunnel burner kW 860/880
Heat duty of reformer MWh 205
Average reformer tube wall temperature °C 930
Catalyst shape Raschig rings with 12 holes
Catalyst dimension 19 × 12 mm/30%; 19 × 16 mm/70%
Catalyst bulk density kg m–3 800
Catalyst porosity 0.51963
Catalyst tortuosity 2.74
Catalyst mean pore radius Å 80
Catalyst characteristic length cm 0.001948
Catalyst active material/support Ni/CaAl12O19
Catalyst quantity tonnes 32.0
Effectiveness factor 0.01
Activity coefficient: start of the run % 130
DSC entry
Natural gas volume flow rate Nm3 h–1 36,500
Hydrogen gas recycle volume flow rate Nm3 h–1 800
MP steam mass flow rate Kg h–1 99,000
Feed gas inlet temperature: to reformer tubes °C 498
Feed gas inlet pressure: to reformer tubes bar 30
Feed gas outlet temperature: exit of reformer tubes °C 790
Feed gas outlet pressure: exit of reformer tubes bar 29.2
Methane outlet molar concentration mol% dry basis 10.30
Natural gas composition: based on in situ gas chromatography
Nitrogen mol% 1.45
Methane mol% 98.36
Carbon dioxide mol% 0.19

As mentioned, the primary aim of this work is to provide process operators with adequate recommendations or advice which will enable them to adjust critical process parameters during SMR operation. With adequate change of process parameters, the working curve (working point) of the reformer catalyst will be brought as much as possible close to the theoretical equilibrium curve with the ultimate goal to improve the reformer catalyst performance and the overall conversion of the SMR endothermic reaction. In order to extract the maximum possible performance of the reformer catalyst, operators will have the possibility to adjust the methane outlet molar concentration and temperature ATE by changing the S/NG molar ratio, heat load or by adjusting the fuel flow distribution to the reformer tubes respectively.

Figures 3 and 4 show the relationship between theoretical equilibrium methane molar concentration (a measure for the theoretical conversion) and temperature, S/NG molar ratios and reforming pressure which is relevant for the endothermic SMR reaction in the industrially interesting range. Simulated results in Figures 3 and 4 show a visible standard pattern of SMR reaction. Namely, at constant reforming pressure, higher S/NG molar ratios and higher temperatures have a beneficial effect on the equilibrium methane molar concentration, but the penalty is a higher energy consumption on account of higher steam mass flow and higher consumption of the fuel gas in the primary reformer furnace. The higher S/NG molar ratio prevents carbon deposition on the catalyst, which may not only increase the pressure drop but also reduce the catalyst activity. On account of this, local reformer tubes overheating (hot spots or bands) will be minimised. In contrast to this, higher temperature will cause higher tube wall temperatures and consequently the shorter lifetime of the reformer tubes due to creep damage and stress as well as rupture.

Fig. 3.

Simulated methane theoretical equilibrium outlet molar concentration versus temperature at various S/NG molar ratios and inlet reforming pressure of 30 bar

Simulated methane theoretical equilibrium outlet molar concentration versus temperature at various S/NG molar ratios and inlet reforming pressure of 30 bar

Fig. 4.

Simulated methane theoretical equilibrium outlet molar concentration versus temperature at various pressures and constant S/NG molar ratio of 3.5

Simulated methane theoretical equilibrium outlet molar concentration versus temperature at various pressures and constant S/NG molar ratio of 3.5

Because the natural gas and steam are at high pressure and the SMR reaction entails an increase in volume, significant savings in compression energy can be achieved if the process operates under elevated pressure. However, thermodynamically, this is unfavourable. Figure 4 shows results of the model test regarding the changes in reforming pressure (on account of the volume increase) versus methane theoretical equilibrium outlet molar concentration. Three different pressures at constant S/NG molar ratio were tested, namely 20 bar, 30 bar and 35 bar, which are common industrial pressure levels.

It is also very important to determine the quantity of heat load to the reformer tubes and to the related catalyst, which will secure sufficient heat for endothermic SMR reaction. Figure 5 shows the relationship between heat load and temperature at different S/NG molar ratios and constant pressure of 30 bar. The major goal is achieving the theoretical methane equilibrium concentrations.

From the simulated results, it can be seen that for higher S/NG molar ratios a higher input of heat load to achieve theoretical methane equilibrium molar concentrations is needed. This effect comes on account of higher feed gas volume flow because of a higher steam mass flow. Although a higher quantity of steam will shift SMR reaction in a favourable direction, consequently it will demand more energy to keep the same methane equilibrium molar concentration, which in the end will have a detrimental effect on the tube wall temperature.

Fig. 5.

Simulated heat loads versus temperature at various S/NG molar ratios and the constant pressure of 30 bar

Simulated heat loads versus temperature at various S/NG molar ratios and the constant pressure of 30 bar

The simulated results are completely in line with the literature findings, which makes the developed model adequately reliable.

According to the simulated results it can be concluded that the best compromise in SMR industrial applications can be achieved when elevated pressure of 30 bar is used (savings in compression energy), medium range of S/NG molar ratio is of 2.8 to 3.5 (energy savings and prevention of carbon deposition) and when the reformer outlet temperature ranges from 780°C to 810°C (which correlates to tube wall temperature up to 930°C: prevention of creep damage and stress to rupture).

With respect to validating the model against real process data, the model was forced to calculate the temperature and pressure drop profile inside of the catalyst tubes with the inlet process gas temperature of 500°C and pressure of 30.5 bar, molar flow of 1530 kmol h–1, average tube wall temperature of 840°C and S/NG molar ratio of 3.5. Figure 6 presents the simulated temperature profile, while Figure 7 shows simulation of the pressure drop along the reformer tubes.

Fig. 6.

Simulated temperature profile along the length of the reformer tube

Simulated temperature profile along the length of the reformer tube

Fig. 7.

Simulated pressure drop profile along the length of the reformer tube

Simulated pressure drop profile along the length of the reformer tube

In the real conditions the outlet reformer tube temperature of the referenced SMR unit at the same inlet process conditions was at the average temperature level of 801°C, while the pressure drop was at the level of 1.58 bar. The difference of 20°C and 0.37 bar give enough validity to be used at the industrial level.

For the same inlet process conditions, the molar flux for all reaction species (methane, carbon dioxide, carbon monoxide, hydrogen and water) as a function of the reformer tube length has been evaluated, and the results are presented in Figure 8. From the simulated results it can be concluded that the production rate of hydrogen increases through the whole reformer tube profile reflecting the methane conversion of 64.64%, methane outlet concentration of 10.30 mol% and ATE of 14°C. From the simulated data it can be seen that the reaction equilibrium is shifted in the direction of reactants approximately in the first 500 mm of the reformer tube length, the same is also visible in Figure 6 through the temperature drop. This can be explained by the endothermic nature of the SMR reaction, because the absorption of heat from the arch burners is mostly implemented in this first part after which the full equilibrium is achieved and the SMR reaction starts to proceed toward reforming products.

Fig. 8.

Simulated molar flux profiles for SMR reaction

Simulated molar flux profiles for SMR reaction

To properly address the evaluation of the catalyst performance in industrial conditions, the model determines the theoretical methane ATE at given process conditions and compares this value with the measured outlet methane molar concentration. Properly designed reformers should, with new catalyst, have methane ATE much lower than 5°C to 10°C. However, plants which have a desulfurisation system often have reformer furnaces operating with methane ATE ranging from 0°C to 3°C. When the evaluation shows this level of methane ATE, this means that the reformer catalyst is of satisfactory performance. Methane ATE levels above 10°C would correspond to marginal performance and would become a factor in discharging the catalyst.

If the process parameters are not optimally adjusted, methane ATE will regularly be in the marginal range. With the appropriate adjustment of process parameters, the reformer catalyst performance can be brought to a satisfactory range below 10°C. As it was mentioned, the temperature at which the exit gas composition would be at the equilibrium is determined in the model by calculating the equilibrium constants from the material balance and determining the corresponding temperature from the correlating equations.

The reference plant in this work is the top fired SMR unit, the primary task of which is the preparation of the synthesis gas for further production of liquid ammonia in the amount of 1360 tonnes per day. According to the data in Table III it can be seen that the SMR unit operates at a pressure of 30 bar, exit temperature is at the level of 790°C and S/NG molar ratio is 3.5. The reformer catalyst was supplied by Clariant, Switzerland, and it has been in operation for one year. All mentioned process parameters result in the outlet methane molar concentration of 10.30 mol% per dry basis. In order to test the catalyst performance against the model, actual process data was fed from the DCS system to the model. The process data was used for the calculation of the theoretical equilibrium curve (the effectiveness factor was additionally reconciled with the developed scripts and functions in appropriate m-file) and the results of the model were compared with the working point. The working point was at the S/NG molar ratio of 3.5 and at the pressure level of 30 bar. Figure 9 shows evaluation results obtained by the model.

Fig. 9.

Evaluation of the reformer catalyst performance during operation of SMR unit

Evaluation of the reformer catalyst performance during operation of SMR unit

According to the theoretical equilibrium curve at given process conditions, methane ATE is at the level of 14°C (790°C minus 776°C), which indicates marginal performance of the reforming catalyst in operation. As the reforming catalyst has only been in operation for one year, excellent activity and no carbon deposition due to higher S/NG molar ratio of 3.5 can be expected from the material balance, the performed evaluation implies that with adequate adjustments of process parameters, catalyst performance can be improved and brought in the satisfactory range below 10°C. Besides that, natural gas used in the current operation has extremely good purity with almost nil content of sulfur compounds. In combination with the sulfur guard beds in the shape of cobalt-molybdenum and zinc oxide catalyst beds, the possibility for sulfur poisoning and eventual sintering is lowered to minimum extent.

The model suggests that the heat load in the reformer furnace is at the adequate value to achieve the theoretical equilibrium methane ATE. According to the performed model evaluation, the main recommendation is to verify the firing conditions inside of the reformer box. After examining flame patterns, tube wall temperatures and the distribution of the fuel volume flow through arch and tunnel burners, it is concluded that there are opportunities for improvement. By adjusting all the mentioned process parameters, the methane ATE was lowered by 6°C, which resulted in the ultimate value of methane ATE of 8°C.

The online model is able to deliver all results and solutions in a time of 30 s after receiving real process data from the DCS system. In comparison with other similar models for example, Latham (3) and Lao et al. (5) achieved computational speed presenting significant advantage for the end users.

This practical example successfully validated the usefulness of the model in evaluating the reformer catalyst during SMR operation. However, the identical-tube model for top-fired SMR furnace developed in this work has potential for further scientific improvements. One of the promising outlooks for further research is modelling of tubes in different radiation environments, expansion to side-fired SMR furnaces and investigation of combustion heat release patterns. A proposal for future work is to develop a fine grid radiative environment model which will be able to group the tubes into several zones taking into consideration outer (refractory wall zone) and inner radiative areas. With this contribution the model will be able to ensure additional recommendations to the users to adjust firing rate with the final goal to reduce energy consumption and subsequently greenhouse gas emission.

4. Conclusion

Continuous evaluation of the primary reformer catalyst characteristics in SMR operation can significantly improve the performance of the overall unit. In order to ensure adequate insight in the performance of the reformer catalyst, a discrete online model for industrial steam-natural gas reformer within ammonia production has been developed. The primary aim of the model is continuous evaluation of the catalyst performance based on the actual process parameters. Commercial simulator and related scripts and functions in the form of m-files were used to develop a steady-state flowsheet. The model has been tested and calculated data reconciled with real process data for the top fired primary reformer designed by Kellogg Inc. A series of differential equations have been used for a very close description of the kinetics, molar flow, temperature and pressure changes along the reformer tubes based on the previous literature models. Reaction rates follow very closely the theoretical model. Methane molar outlet concentration and its ATE were used for the evaluation and optimisation of the reformer catalyst as two main process parameters. The model can calculate theoretical equilibrium curves of methane outlet molar concentration at different temperatures, pressures and steam-to-natural gas molar ratios. According to the calculated theoretical equilibrium curves, the model is able to compare the working point of the reformer catalyst and ensure adequate recommendations to operators to mitigate their potential marginal performance. Use of the model industrial conditions can also prolong the reformer catalyst lifetime. The computational speed is quick enough for application in actual conditions, while the output recommendations are simple enough to be used by working staff. The present paper brings innovation for easier implementation of an online model in an SMR unit, with possibility to ensure adequate recommendations for process parameters adjustment with respect to better catalyst performance. The related calculation method with benefit of a fast communication routine can significantly improve productivity and uptime and subsequently remediate plant performance deviations with possibility to recognise major catalyst efficiencies.



  1. 1.
    M. Appl, “Ammonia: Principles and Industrial Practice”, Wiley-VCH Verlag GmbH, Weinheim, Germany, 1999, 301 pp LINK
  2. 2.
    “Catalyst Handbook”, ed. M. V. Twigg, 2nd Edn., CRC Press, Boca Raton, USA, 1989
  3. 3.
    D. Latham, “Mathematical Modelling of an Industrial Steam Methane Reformer”, Master’s Thesis, Department of Chemical Engineering, Queen’s University, Kingston, Canada, December, 2008, 279 pp LINK
  4. 4.
    A. Aguirre, “Computational Fluid Dynamics Modelling and Simulation of Steam Methane Reforming Reactors and Furnaces”, PhD Thesis, Department of Chemical and Biomolecular Engineering, University of California, Los Angeles, USA, 2017, 223 pp LINK
  5. 5.
    L. Lao, A. Aguirre, A. Tran, Z. Wu, H. Durand and P. D. Christofides, Chem. Eng. Sci., 2016, 148, 78 LINK
  6. 6.
    J. E. Holt, J. K. Kreusser, A. Herritsch and M. Watson, ANZIAM J., 2017, 59, C112 LINK
  7. 7.
    J. Xu and G. F. Froment, AIChE J., 1989, 35, (1), 97 LINK
  8. 8.
    J. Xu and G. F. Froment, AIChE J., 1989, 35, (1), 88 LINK
  9. 9.
    Z. Wu, A. Aguirre, A. Tran, H. Durand, D. Ni and P. D. Christofides, Ind. Eng. Chem. Res., 2017, 56, (20), 6002 LINK
  10. 10.
    L. Sun, “Modelling and MPC for a Primary Gas Reformer”, Masters Thesis, Department of Chemical and Materials Engineering, University of Alberta, Edmonton, Canada, 2013, 88 pp LINK
  11. 11.
    A. Meziou, P. B. Deshpande and I. M. Alatiqi, Int. J. Hydrogen Energy, 1995, 20, (3), 187 LINK
  12. 12.
    I. M. Alatiqi and A. M. Meziou, Comput. Chem. Eng., 1991, 15, (3), 147 LINK
  13. 13.
    C. M. Schillmoller and U. W. van den Bruck, Hydrocarbon Proc., 1984, 63, (12), 55
  14. 14.
    N. A. El Moneim, I. Ismail and M. M. Nasser, Int. J. Novel Res. Dev., 2018, 3, 11
  15. 15.
    D. A. Latham, K. B. McAuley, B. A. Peppley and T. M. Raybold, Fuel Process. Technol., 2011, 92, (8), 1574 LINK
  16. 16.
    J. S. Lee, J. Seo, H. Y. Kim, J. T. Chung and S. S. Yoon, Fuel, 2013, 111, 461 LINK
  17. 17.
    F. Minette, M. Lugo-Pimentel, D. Modroukas, A. W. Davis, R. Gill, M. J. Castaldi and J. De Wilde, Appl. Catal. B: Environ., 2018, 238, 184 LINK
  18. 18.
    J. R. Rostrup-Nielsen, J. Sehested and J. K. Nørskov, Adv. Catal., 2002, 47, 65 LINK
  19. 19.
    A. J. Forester and B. J. Cromarty, ‘Theory and Practice of Steam Reforming’, ICI/Katalco/KTI, UOP, 3rd Annual International Seminar of Hydrogen Plant Operation, 7th–9th June, 1995, Chicago, USA, 1995
  20. 20.
    S. Z. Abbas, V. Dupont and T. Mahmud, Int. J. Hydrogen Energy, 2017, 42, (5), 2889 LINK
  21. 21.
    J. Rostrup-Nielsen, and L. J. Christiansen, “Concepts in Syngas Manufacture”, ed. G. J. Hutchins, Vol. 10, Catalytic Science Series, Imperial College Press, London, UK, 2011
  22. 22.
    S. S. E. H. Elnashaie and F. Uhlig, “Numerical Techniques for Chemical and Biological Engineers Using MATLAB®: A Simple Bifurcation Approach”, Springer Science and Business Media LLC, New York, USA, 2007, 588 pp LINK
  23. 23.
    A. Olivieri and F. Vegliò, Fuel Process. Technol., 2008, 89, (6), 622 LINK
  24. 24.
    F. M. Alhabdan, M. A. Abashar and S. S. E. Elnashaie, Math. Comput. Model., 1992, 16, (11), 77 LINK
  25. 25.
    S. S. E. H. Elnashaie and S. S. Elshishini, “Modelling, Simulation and Optimization of Industrial Fixed Bed Catalytic Reactors”, Topics in Chemical Engineering, Vol. 7, Gordon and Breach Science Publishers SA, Yverdon, Switzerland, 1993
  26. 26.
    E. L. Cussler, “Diffusion: Mass Transfer in Fluid Systems”, 2nd Edn., Cambridge University Press, Cambridge, UK, 1997, 580 pp
  27. 27.
    E. B. Nauman, “Chemical Reactor Design, Optimization, and Scaleup”, 2nd Edn., John Wiley and Sons Inc, Hoboken, USA, 2008, 608 pp LINK
  28. 28.
    G. M. Hampson, Chem. Eng., 1979, 523
  29. 29.
    G. F. Froment, K. B. Bischoff and J. De Wilde, “Chemical Reactor Analysis and Design”, 3rd Edn., John Wiley and Sons, New York, NY, USA, 2010, 606 pp
  30. 30.
    S. S. E. H. Elnashaie and A. Adris, “Fluidized Bed Steam Reformer for Methane”, Proceedings of the IV International Fluidization Conference, Banff, Canada, May, 1989
  31. 31.
    S. S. E. H. Elnashaie, A. M. Adris, M. A. Soliman and A. S. Al-Ubaid, Can. J. Chem. Eng., 1992, 70, (4), 786 LINK
  32. 32.
    N. Zecevic and N. Bolf, Ind. Eng. Chem. Res., 2020, 59, (8), 3458 LINK
  33. 33.
    J.-M. Lavoie, Front. Chem., 2014, 2, 81 LINK



reformer tube cross sectional area; m2


pre-exponential factor; kmol bar5 kg–1cat h–1 or kmol kg–1cat bar–1 h–1


approach to equilibrium; °C or K


pre-exponential factor; bar–1






carbon monoxide


carbon dioxide


computational fluid dynamics


component object model


specific heat capacity of process gas i; kJ kmol–1 K–1


diffusivity of the reacting species A


Knudsen diffusivity; m3fluid m–1catalyst s–1


effective diffusivity; m3fluid m–1catalyst s–1

average diffusivity of species A


distributed control system




activation energy; kJ mol–1


molar flow rate; kmol h–1


mass velocity of the process gas; kg m–2 h–1


catalytic layer thickness; m






adsorption constants for component i; bar–1


equilibrium constants for SR and WGS


rate coefficient for reaction; kmol bar5 kg–1cat h–1 or kmol kg–1cat bar–1


reformer tube heated length; m


mass; kg or tonnes


molecular weight; kg kmol–1


pressure and partial pressure; Pa, bar or atm


partial pressure of component i; bar


plug flow reactor


gas constant; kJ kmol–1 K–1


Reynolds number

rCO, CO2, CH4, H2

rate of formation and disappearance; kmol kg–1cat h–1


rate of chemical reaction; kmol kg–1cat h–1


pore radii


void fraction


steam-to-natural gas molar ratio


steam methane reforming


temperature; °C or K


reformer tube skin temperature; K


time; s, min or h


heat transfer capacity; kJ m–2 h–1 K–1


velocity of the fluid; m s–1


molar rate of conversion; mol%


water gas shift


enthalpy change at 298 K; kJ mol–1


convective heat transfer coefficient in the preheated bed; kJ m–2 h–1 K–1

catalyst bed voidage or porosity; m3fluid m–3cat


process gas thermal conductivity; kJ m–1 h–1 K–1


tube metal thermal conductivity; kJ m–2 h–1 K–1


dynamic viscosity of the fluid; Pa s, N s m–2 or kg m–1 s–1

ηCH4, H2O

constant effectiveness factor


fluid density; kg m–3


catalyst bed density; kg m–3


Thiele modulus


intracatalyst coordinate


catalyst tortuosity

Supplementary Information

MATLAB®, The MathWorks Inc, Natick, USA LINK

UniSim® Design, Version R470, Software for Process Design and Simulation, Honeywell International Inc, Charolotte, USA LINK

The Author

Nenad Zečević joined Petrokemija Plc, Croatia, in 1999 after he obtained his Master’s degree in Chemical Engineering from the University of Zagreb, Faculty of Natural Sciences, Croatia. Having worked as an ammonia plant manager, he is currently a Production & Maintenance Manager with 20 years of operational experience in the fertiliser industry. He is also a PhD candidate from the University of Zagreb, Faculty of Chemical Engineering and Technology, Croatia with the thesis Advancement of Process Control System in Ammonia Production.

Read more from this issue »