Kinetics of CO₂ hydrogenation to methanol on Cu/ZnO/ZrO₂ based on an extensive dataset

Abstract

The kinetics of CO₂ hydrogenation to methanol over a self-developed Cu/ZnO/ZrO₂ (CZZ) catalyst was studied in a wide range of process conditions. Experiments were performed at industrially relevant pressures (30–60 bar) and temperatures (190–250 °C), as well as H₂ to CO₂ ratios between 1 and 6, addressing the use of hydrogen from renewable energy sources and the use of CO₂ as a C1 raw material in Power-to-X technologies. The CZZ catalyst has shown improved performance and higher stability in CO₂ hydrogenation to methanol in comparison to other Cu/ZnO-based catalysts. A mathematical description of the kinetics is crucial to enable model-based design for the industrial implementation of this catalyst. Therefore, a lumped 6-parameter kinetic model was developed to fit the experimental data, resulting in one of the predictive models with the broadest validity range (experimental database of 500 points) for the CZZ system. This new kinetic model is compared to state-of-the-art literature models with more parameters, and our model performs equally well or even better in terms of sensitivity to process parameters and extrapolability.

---

1 Introduction

Methanol is a key intermediate for the chemical industry, being present in the production of various value-added chemicals such as formaldehyde, olefins and dimethyl ether. It is also gaining relevance in the energy transition towards renewable sources and regarding net-zero emission goals, due to the combination of renewable H₂ and CO₂ obtained via carbon capture.^1–3 Methanol is traditionally produced from fossil-based syngas, composed of CO, H₂ and small amounts of CO₂. However, it can also be produced from CO₂/H₂ feeds. The reactions for methanol formation are given in eqn (1) and (2), with the reverse water–gas shift (rWGS) as a parallel reaction (eqn (3)).

CO + 2 H₂ ⇌ CH₃OH; ΔH⁰_rxn,298K = −90.5 kJ mol⁻¹ (1)

CO₂ + 3 H₂ ⇌ CH₃OH + H₂O; ΔH⁰_rxn,298K = −49.3 kJ mol⁻¹ (2)

CO₂ + H₂ ⇌ CO + H₂O; ΔH⁰_rxn,298K = 41.2 kJ mol⁻¹ (3)

The state-of-the-art catalyst Cu/ZnO/Al₂O₃ (CZA) has been extensively studied and optimized for CO hydrogenation to methanol.⁴ Nevertheless, for CO₂ hydrogenation to methanol, this catalyst exhibits limited performance due to the higher occurrence of rWGS and enhanced formation of water, reducing the number of active sites over time.^5,6 On the other hand, Cu/ZnO/ZrO₂ was shown to be more adequate for CO₂ hydrogenation due to the lower affinity of ZrO₂ to water and enhanced CO₂ adsorption.^7–9

Besides an efficient catalyst, kinetic models are essential to optimize reactor and process design. Especially in the case of CO₂ hydrogenation to methanol, innovative designs and process concepts are desired, to overcome the disadvantages of the lower thermodynamic equilibrium conversions achieved in comparison to CO hydrogenation. In a previous publication from our group,¹⁰ a kinetic model proposed previously by the authors¹¹ was used to carry out a techno-economic analysis of a process concept with three reactors in series with inter-stage condensation. In Bagwan et al.,¹² the model proposed by Portha et al.¹³ for CZA was used to simulated isothermal and adiabatic methanol reactors, analyzing the effect of H₂:CO₂ ratio, initial temperature and pressure. For the results to be significant, the models need to be robust and developed in conditions close to the intended operation. Nyári et al.¹⁴ showed that the choice of kinetic model for a techno-economic analysis can lead to a 10% difference in the calculated levelized costs of CO₂-based methanol.

Most reported kinetic models for methanol synthesis were validated for the commercial catalyst CZA; some examples are discussed in more detail in Section 2. Regarding models developed for CZZ, Portha et al.¹³ carried out 10 experiments and refitted the kinetic parameters (pre-exponential factor and activation energy) for CO₂ hydrogenation and reverse water–gas shift from Graaf et al.¹⁵ The adsorption parameters and the kinetic parameters for CO hydrogenation were retrieved from the original model. A similar modeling approach was followed by Marcos et al.,¹⁶ using 45 experimental points for each of the two investigated ZrO₂ polymorphs. Khawaja and Usman¹⁷ used 31 experimental points reported in Arena et al.¹⁸ to test the integral and the differential methods for a power law model and the kinetic models from Graaf et al.¹⁵ and Park et al.¹⁹ The authors concluded that the integral method, for all models, showed smaller deviations. Beyond that, the model from Graaf et al.¹⁵ with refitted parameters provided more accurate results. Dong et al.²⁰ reported 20 experiments measured at atmospheric pressure and refitted the models from Kubota et al.²¹ and from Vanden Bussche and Froment²² with some assumptions reported by the former. It is important to highlight that, in the work from Dong et al.,²⁰ no direct measurements of methanol are reported; instead, the methanol outlet was calculated via carbon balance. Furthermore, the reported methanol outlet concentrations are implausible regarding thermodynamic equilibrium, which predicts only traces of methanol for atmospheric pressure.

In the present work, we significantly expand the validation ranges of these models, providing 500 experimental points. The dataset provided by this publication includes experiments up to 60 bar and several H₂/CO₂ ratios. In Fig. 1, we summarize the validation ranges of each model. The H₂/CO₂ ratio is given in terms of stoichiometric number (SN), as defined in eqn (4), to enable further comparison with experimental data including CO in the feed.

Fig. 1 Training ranges of published models for CZZ. *45 points for each catalyst.

The objective of this work is to propose a kinetic model for methanol synthesis on Cu/ZnO/ZrO₂ and compare it to state-of-the-art models from literature in terms of sensitivity to operating variables and prediction of external datasets.

2 Theory

In this section, an overview of kinetic models for methanol synthesis over the commercial catalyst CZA is presented.

2.1 Graaf et al.

The model proposed by Graaf et al.¹⁵ includes reaction rate equations for CO hydrogenation, CO₂ hydrogenation and reverse water–gas shift reaction. The adsorption model considers two sites: one for CO and CO₂ and one for H₂ and H₂O. H₂ is assumed to adsorb dissociatively. The authors performed a screening, testing each elementary reaction as a potential rate-determining step. These are presented in the SI. The chosen model is the one that provided the lowest deviation for the formation rates of methanol and water. The derived reaction rate equations are given in eqn (5)–(7):

The model contains 12 parameters. The kinetic constants k_j and the adsorption constants K_j are given respectively as Arrhenius and Van't Hoff equations:

The equilibrium terms η are given by eqn (10)–(12):

The K⁰_p(T) parameters correspond to the equilibrium constants. In this work, they are calculated using the NASA 7-coefficient polynomials.^23,24

2.2 Nestler et al.

The kinetic model proposed by Nestler et al.²⁵ is based on the work of Henkel.²⁶ In this work, Graaf's model is simplified by eliminating CO direct hydrogenation, since, according to the calculations presented in the thesis, more than 99% of the produced methanol comes from CO₂. This conclusion is in agreement with theoretical calculations and experiments for the Cu/ZnO system,^27,28 which state that methanol is mostly produced from CO₂. This simplification reduces the model from 12 to 10 parameters, and one more parameter is removed by eliminating the temperature dependence in the CO₂ adsorption. In the work of Nestler et al.,²⁵ the parameters from eqn (6) and (7) are refitted with literature data from Park et al.¹⁹

2.3 Slotboom et al.

The mechanism presented in Slotboom et al.²⁹ is based on the work of Vanden Bussche and Froment²² with mechanistic updates from Grabow and Mavrikakis.³⁰ Instead of having one active site as in Vanden Bussche and Froment,²² the authors propose a mechanism with three active sites, based on Seidel et al.:³¹ one site for oxidized surface centers, one for reduced surface centers and one for hydrogen. It is assumed that some oxygenated intermediates are adsorbed in two sites of the same nature.

The authors propose two models: in a more extended model, with 10 parameters, the adsorption terms are all taken into account, and their temperature dependence is considered. This model is further simplified to 6 parameters by assuming that the hydrogen sites are always occupied and that the adsorption isotherms are linearly dependent. The 6-parameter model is shown in eqn (13) and (14). Similarly to Graaf, the authors decide on the rate-determining steps by testing all elementary steps. These are given in the SI.

2.4 Lacerda de Oliveira Campos et al.

Lacerda de Oliveira Campos et al.¹¹ proposed three lumped kinetic models, based on a microkinetic model developed by the author³² using theoretical calculations reported by Studt et al.^27,33 The models assume three sites – pure Cu, Cu/Zn and one site for hydrogen and water –, however, in the simplified model for CO₂-containing feeds, the Cu site is assumed to be always occupied, reducing the number of parameters in the model to 6. The term ϕ_Zn accounts for the surface changes in the catalyst, depending on the inlet gas composition. The model for CO₂-containing feeds is given in eqn (15) and (16).

3 Experiments

The kinetic studies were performed in an experimental setup developed for conversion of CO_X and H₂ at pressures up to 150 bar. The dimensions of the stainless steel fixed-bed reactor are 440 mm in length and 17.4 mm in inner diameter. The reactor also contains an inner concentric tube (3.175 mm external diameter) for the axial temperature measurements and is heated by four independent circuits, to reduce temperature differences along the length. Hydrogen (H₂, 99.999%), nitrogen (N₂, 99.999%), carbon monoxide (CO, 99.97%) and carbon dioxide (CO₂, 99.7%) are supplied by Air Liquide Germany GmbH. H₂, N₂ and CO supplies are regulated via mass flow controllers (MFCs, Bronkhorst High Tech). CO₂ is cooled down to 0 °C by a cryostat (Huber Minichiller 280 Olé) and pressurized in a high-pressure liquid chromatography (HPLC) pump (WADose-Lite-HP, type 10-SS-ITU-C, Flusys GmbH). The pressure of CO₂ is regulated through a back-pressure regulator (Equilibar), and the flow is dosed via a CoriFlow MFC (Bronkhorst High Tech). The pressure in the system is in addition controlled through a back-pressure regulator (Equilibar). Reactants (via bypass) and products are analyzed with a Fourier transform infrared spectrometer (FTIR, Gasmet CX4000) coupled to a hydrogen analyzer (H₂-/TCD-Process Gas Analyzer CONTHOS). A schematic flowchart of the setup is shown in Fig. 2.

Fig. 2 Schematic flowchart of the reactor setup.

The Cu/ZnO/ZrO₂ (CZZ) catalyst used is prepared by continuous co-precipitation followed by aging, filtering and calcination. Details about the synthesis process and the setup are available in previous works.^34,35 The obtained material has a molar ratio of 62.7% Cu, 29.5% Zn and 7.8% Zr after calcination. Manufacturing details and characterization of the catalyst are available elsewhere.³⁶ This catalyst serves as a reference catalyst at the Institute of Catalysis Research and Technology, due to the reproducibility of the synthesis method and scalability of its production.^37,38

For the experimental run reported here, 1.0 g of CZZ (250–500 μm) was physically mixed with 10.0 g of silicon carbide (SiC, Mineraliengrosshandel Hausen GmbH), to avoid hot spots. This mixture was filled into the reactor, forming a catalytic bed of 2.6 cm length. Pure SiC was placed at the top and bottom of the catalytic bed. The catalyst was reduced at atmospheric pressure using the following procedure: a volume flow of 600 mL_N min⁻¹ (5% v/v of H₂ in N₂) was fed to the reactor, while heating from 90 °C to 120 °C at a rate of 10 °C h⁻¹ took place. With the same flow, heating proceeded until the temperature of 200 °C at a rate of 7.5 °C h⁻¹, and the final temperature was held for 1 hour. Afterwards, the H₂ content was changed to 50% v/v, and the temperature was increased to 220 °C with a ramp of 10 °C h⁻¹. These conditions were maintained for another hour, after which the reactor was purged with N₂ and cooled to 180 °C.

Kinetic experiments were carried out at temperatures between 190 and 250 °C, pressures of 30, 50 and 60 bar and space velocities (SV) of 48 and 72 Nm³ kg_cat⁻¹ h⁻¹. To evaluate the sensitivity of methanol formation and reverse water–gas shift with respect to H₂ and CO₂, the nominal gas composition was varied in a manner that each component (CO₂ or H₂) was held constant at a time, and the other was varied. The N₂ fraction was adjusted accordingly to maintain the total volume flow. Additionally to this variation, experiments with stoichiometric H₂ : CO₂ ratio for the methanol synthesis (3 : 1) without dilution in N₂ were also performed. The used gas compositions are given in Table 1. Detailed experimental data are provided in the SI.

Table 1 Nominal gas compositions used in the experiments (% v/v)

CO2	H2	N2
20	20	60
20	30	50
20	40	40
20	50	30
20	70	10
20	80	0
20	60	20
30	60	10
25	60	15
15	60	25
10	60	30
25	75	0

---

Prior to the kinetic measurements reported in this paper, the reactor was operated for 50 hours at 30 bar, temperatures between 210 and 250 °C and all the gas compositions given in Table 1. The duration of this period was determined based on reference measurements reported in the SI (Fig. S1). After this, each condition was kept for at least 40 minutes, to ensure that the changes in the operating parameters were completed.

4 Modeling

First, the carbon balance for all the data points was calculated as given in eqn (17):

In which the percentages correspond to the volume fractions obtained in the FTIR. The division by 1 + 0.02·%CH₃OH_out accounts for the stoichiometry of the methanol synthesis. All points used for the modeling are within a ±4% deviation. The 500 experimental points were randomized and divided into 5 groups – or folds – containing 100 points each. In the so-called cross-validation method, each group is removed from the complete set at a time, with the remaining points being used for parameter estimation. The removed points are used for validation.³⁹

The reactor was modeled as a plug-flow reactor (PFR) operating under integral conditions. To justify the assumption of plug-flow, we have followed the criteria reported in Raja et al.⁴⁰ to exclude radial convection and axial diffusion in comparison to axial convection. Internal mass transfer limitations were neglected based on the Weisz–Prater criterion.^41,42 According to Levenspiel,⁴³ film diffusion resistance is unlikely to affect the reaction rate in a gas/porous catalyst system, hence it was neglected in our study. Pressure drop along the bed length was calculated based on Ergun's equation⁴⁴ for the worst-case scenario – highest flow velocity –, and the obtained value was lower than 1 mbar, therefore it was neglected.

The mass balance for each component i is given by eqn (18):

Here, ν_ij corresponds to the stoichiometric coefficient of component i in reaction j, and Ṅ corresponds to the total molar flow.

The variation in the number of moles due to the stoichiometry of CO₂ hydrogenation is given in eqn (19):

Fugacities were calculating using the Peng–Robinson equation of state,⁴⁵ following the methodology described in a previous publication from our research group.⁴⁶ The binary interaction parameters k_ij were retrieved from the works of Meng et al.,^47,48 and an effective hydrogen acentric factor ω = −0.05 was used.⁴⁹

The objective function for parameter estimation corresponds to minimizing the normalized squared errors of the carbon-containing compounds, shown in eqn (20):

Three recent models from literature were re-fitted using this objective function: Nestler et al.²⁵ (eqn (6) and (7)), the simplified model from Slotboom et al.²⁹ (eqn (13) and (14)) and the 6-parameter model from Lacerda de Oliveira Campos et al.¹¹ (eqn (15) and (16)). In the latter, the zinc coverage term ϕ_Zn was lumped into the other parameters, since all experimental points in this work would have the same Zn coverage (0.1). Experimental and theoretical studies about the dynamics of Cu/Zn-based catalysts can be found in the literature.^9,36,50,51 Nevertheless, in lumped kinetic models, which are intended for process optimization and design, the focus is on understanding the overall reaction kinetics. In these situations, the common approach is to incorporate these structural-activity dynamics into the existing parameters.

To account for the activity loss of the catalyst with respect to time, an activity term (eqn (21)) was included in all re-fittings, as done in Rodrigues Niquini et al.⁵² As already discussed in our previous publication,⁵² this term should not be extrapolated to industrial catalyst lifetimes and is intended solely to improve the quality of the kinetic parameters.

In this equation, k_d corresponds to the deactivation constant, ToS corresponds to the time on stream, and t₀ is a reference time, in which the activity a_CZZ is equal to 1. In this work, we chose t₀ = 50 h, corresponding to the conditioning time. To monitor the catalytic activity over time, reference measurements were carried out along the experimental campaign. These reference measurements – not used for the model construction – are shown in the SI, along with the model predictions with respect to time on stream. In these plots, it is possible to observe that the measurements for the initial hours follow a different trend, reason why they were neglected in the kinetic model. Furthermore, the decay in methanol and CO formation over time follows the trend predicted by the model, indicating that the inclusion of an activity term is appropriate.

Parameter estimation was carried out in Matlab R2021b, using the built-in fminsearch function, with a tolerance of X equal to 10⁻⁴ and function tolerance equal to 10⁻⁶. Molar fractions and molar flow were integrated using the ode15s function, with absolute and relative tolerances set to 10⁻⁸. Pre-exponential factors, depending on their order of magnitude, were estimated in logarithmic form, to improve numerical sensitivity. To avoid divisions by zero in the reactor inlet, initial molar fractions for all components were set as at least 1 × 10⁻⁷ (0.1 ppm). Different sets of parameters were given as initial guesses, as an attempt to find the global minimum.

5 Results and discussion

5.1 Parameter estimation and new proposed kinetic model

The parameters for the three models, including the deactivation term, were successfully estimated. In the case of Lacerda de Oliveira Campos et al.,¹¹ it was observed that the term K₃ was orders of magnitude greater than 1, meaning, in practice, that the sites for H₂ and H₂O are always occupied. For the model, it translates to K₃ being lumped into the other parameters. Hence, we propose a new kinetic model, containing 6 parameters, including one to describe the catalyst deactivation over time. The reaction rates are given in eqn (22) and (23).

This modified model is marginally more accurate than the original formulation, with χ² values of, respectively, 4.99 and 5.05 for the best sets. A more detailed comparison between these two formulations, as well as the estimated parameters and the parity plots for the original one, are given in the SI.

The parameters estimated for the new proposed model are given in Table 2. As shown by the confidence intervals, all parameters are statistically significant. In Fig. 3, the parity plots for CH₃OH and CO are shown. A slight underestimation is observed for points with higher methanol formation; still, all 500 points are within ±20%. For CO, 455 points are within ±10% deviation, and 492 are within ±20%. The predictions for CO₂ are all within the ±10% range for this work and the other tested models.

Table 2 Estimated parameters and confidence intervals for the kinetic model proposed in this work (eqn (22) and (23)) *R in J mol⁻¹ K⁻¹

Parameter		Value ± 95% CI		Unit
k CO2hyd	ln A	12.760 ± 0.047		mol s^−1 kgcat^−1 bar^−2
	B	−13 095 ± 27
k rWGS	ln A	24.574 ± 0.053		mol s^−1 kgcat^−1 bar^−1.5
	B	−16 278 ± 31
K 2	A	0.1056 ± 0.0035	K 2 = 0.1056	bar^−1.5
k d		(5.673 ± 0.078) × 10^−4	k d = 5.673 × 10^−4	h^−1

---

Fig. 3 Parity plots for (a) CH₃OH and (b) CO for the model proposed in this work. The dashed lines represent a ±20% deviation.

In Tables 3 and 4, the refitted parameters for Slotboom and Nestler models are presented, respectively. For both models, all parameters are statistically different from zero. The parity plots are available in the SI. For the refitted Slotboom model, 409 out of the 500 points are within the ±10% range for methanol and 497 are within the ±20% range, performing similarly to this work. For CO, only 284 points are within the ±10% range. For Nestler, 492 out of the 500 experimental points are within the ±10% range for methanol, being the most accurate of the investigated models for this substance. On the other hand, only 301 points are within the ±10% range for CO.

Table 3 Estimated parameters and confidence intervals for the kinetic model proposed in Slotboom et al.²⁹ (eqn (13) and (14)) *R in J mol⁻¹ K⁻¹

Parameter		Value ± 95% CI		Unit
k CO2hyd	ln A	33.83 ± 0.15		mol s^−1 kgcat^−1 bar^−1
	B	−17 986 ± 48
k rWGS	ln A	41.932 ± 0.045		mol s^−1 kgcat^−1 bar^−0.5
	B	−22 872 ± 14
k H2	A	8.80 ± 0.31	k H2 = 8.80	bar^0.5
k H2O/9	A	950 ± 47	k H2O/9 = 950	—
k d		(8.24 ± 0.95) × 10^−4	k d = 8.24 × 10^−4	h^−1

---

Table 4 Estimated parameters and confidence intervals for the kinetic model proposed in Nestler et al.²⁵ (eqn (6) and (7)) *R in J mol⁻¹ K⁻¹

Parameter		Value ± 95% CI		Unit
k CO2hyd	ln A	3.640 ± 0.055		mol s^−1 kgcat^−1 bar^−1
	B	−5596 ± 22
k rWGS	ln A	25.493 ± 0.019		mol s^−1 kgcat^−1 bar^−0.5
	B	−16 171 ± 13
K CO	ln A	−8.09 ± 0.13		bar^−1
	B	15 180 ± 100
K CO2	ln A	18.173 ± 0.091	K CO2 = 7.8066 × 10^7	bar^−1
K H2O	ln A	−17.845 ± 0.092		bar^−0.5
	B	9795 ± 52
k d		(6.592 ± 0.049) × 10^−4	k d = 6.592 × 10^−4	h^−1

---

The total χ², as well as the errors MSE and MRE for each model, are given in Table 5. The present model has the lowest χ² (4.99) among the tested models. Regarding the deviations for the carbon-containing components, this work clearly outperforms the literature models for CO prediction. For CO₂, the performance of the models is similar. For methanol, the deviations have the same order of magnitude for all models, with Nestler showing the lowest errors.

Table 5 Calculated errors for the present work and the two refitted literature models

			MSE 10^3			MRE 10^2
	Parameters	χ ^2	CO	CO2	CH3OH	CO	CO2	CH3OH
This work	6	4.99	4.12	0.29	5.58	4.92	1.44	6.23
Nestler	10	8.31	14.79	0.24	1.60	9.93	1.36	3.05
Slotboom	7	12.30	18.66	0.21	5.73	10.75	1.28	6.13

---

In Fig. 4, the effect of temperature on CH₃OH and CO formation is analyzed at 30 and 60 bar, as well as the sensitivity of the models. At 30 bar, methanol formation increases with temperature until it approaches thermodynamic equilibrium. At 60 bar, the thermodynamic limitation is less visible. On the other hand, CO formation rises exponentially at the investigated range, following the expected Arrhenius behavior. All models capture the trends well, with slight differences at higher temperatures.

Fig. 4 Variations in temperature for the following conditions: y_CO2,in = 0.18, y_H2,in = 0.59, SV = 48 Nm³ kg_cat⁻¹ h⁻¹, (a and b) p = 31 bar; (c and d) p = 50 bar. Dashed line represents thermodynamic equilibrium.

In Fig. 5, the experiments' and models' behavior for different pressures is shown. Methanol formation increases slightly for higher pressures, and CO formation stays approximately constant. These results show that, from a kinetic point of view, it is beneficial to carry out methanol synthesis at higher pressures, since selectivity to methanol is enhanced. Nevertheless, this increase in methanol formation is in a smaller order of magnitude than thermodynamic equilibrium predictions (dashed line in Fig. 5c). The trends are well captured by the models, with Nestler slightly overestimating the effect of high pressures on methanol formation.

Fig. 5 Variations in pressure for the following conditions: y_CO2,in = 0.18, y_H2,in = 0.59, SV = 48 Nm³ kg_cat⁻¹ h⁻¹, (a and b) T = 200 °C; (c and d) T = 240 °C. Dashed line represents thermodynamic equilibrium.

Fig. 6 shows the effect of the partial pressure of H₂ on the experiments and the models. For methanol formation, an increase in the partial pressure of H₂ leads to an almost linear increase in methanol formation. This result matches recent literature for mechanistic insights, which states that the apparent reaction order of methanol formation concerning H₂ is close to 1 on Cu-based catalysts.^53,54 All models capture this trend well. Regarding CO formation, the partial pressure of H₂ has low sensitivity. The model proposed in this work is the one that better captures this trend, while Nestler shows a slight increase and Slotboom, a slight decrease. This experimental result also matches mechanistic findings.

Fig. 6 Variations in H₂ partial pressure for the following conditions: y_CO2,in = 0.18, T = 210 °C, SV = 72 Nm³ kg_cat⁻¹ h⁻¹, (a and b) p = 31 bar, (c and d) p = 62 bar. Dashed line represents thermodynamic equilibrium.

In Fig. 7, we show the behavior for varying CO₂ partial pressures. CH₃OH and CO formation both present a slight increase for higher CO₂ partial pressures. This observation matches recent mechanistic findings,⁵³ which predict CO₂ as an adsorption species, hindering the reaction rate. This work and the model from Nestler predict this behavior well, as they both have CO₂ in the denominators. Slotboom, nevertheless, neglects CO₂ adsorption in the model by omitting formate (HCOO) from the site coverage.

Fig. 7 Variations in CO₂ partial pressure for the following conditions: y_H2,in = 0.59, T = 210 °C, SV = 72 Nm³ kg_cat⁻¹ h⁻¹, (a and b) p = 31 bar, (c and d) p = 62 bar.

In summary, all three models presented in this work show good results for methanol and CO predictions on CO₂-rich feeds, especially for the sensitivity analyses. The model from Slotboom slightly overestimates the sensitivity with respect to CO₂ (see Fig. 7), which is a hint that including CO₂ as an adsorption species could improve the model, corroborating experimental evidence. Comparing Nestler and this work, the former presents slightly better prediction results for CH₃OH and the latter, for CO. Nevertheless, the present model contains 4 parameters fewer than Nestler, facilitating its implementation in commercial software for plant simulations, such as Aspen Plus and gPROMS.

5.2 Validation with external datasets

The extrapolability of the model was verified by carrying out simulations for external datasets. Even for CO₂-rich applications, it is advantageous for the model to predict experiments with CO in feed, since it is present in setups with recycle or multi-stage reactors. For this, we simulated previous experiments from our group, most of which are reported in Rodrigues Niquini et al.⁵² The experimental conditions are summarized in Table 6, and a complete list is available in the SI. These experiments were carried out on a self-developed Cu/ZnO/ZrO₂ catalyst with a similar composition to the one used in this work.

Table 6 Summary of experimental conditions from the external datasets simulated

Parameter	Rodrigues Niquini et al.^52	Arena et al.^18
SV/Nm^3 kgcat^−1 h^−1	8.9–17.7	9.6–60.0
SN	0.76–3.04	2
CO2 : COX	0.12–1	1
Temperature/°C	190–250	180–240
Pressure/bar	31–61	30–50

---

Parity plots for each model are shown in Fig. 8. The model proposed in this work presents a good performance for lower methanol formations, and slightly underestimates higher methanol yields. Slotboom presents a similar behavior, with a more severe underestimation. However, the refitted model from Nestler predicts almost no reaction for CO-containing feeds. This occurs probably due to mechanistic assumptions for the adsorption term. In Nestler's model, CO* and are the dominant carbon-containing surface species, resulting in a high sensitivity of the model to the fugacity of CO. The model from Slotboom and the one proposed in this work propose more complex adsorption mechanisms, resulting in a lower sensitivity to CO partial pressure.

Fig. 8 Parity plots for (a) CH₃OH and (b) CO for the dataset published in Rodrigues Niquini et al.⁵² The dashed lines represent a ±20% deviation.

We also simulated the experiments from Arena et al.¹⁸ used in the study developed by Khawaja and Usman.¹⁷ These experiments were carried out on a Cu/ZnO/ZrO₂ prepared via reverse co-precipitation, and a summary of the experimental conditions is also given in Table 6. In this case, the activity of the catalyst was assumed to be always equal to 1, since no data for time on stream was given. The parity plots are shown in Fig. 9. For methanol, all models exhibit a good performance, with most of the points within or close to the ±20% range. For CO, the models from Nestler and Slotboom underestimate some of the values, while all points estimated by our model are within or close to the ±20% range.

Fig. 9 Parity plots for (a) CH₃OH and (b) CO for the dataset published in Arena et al.¹⁸ The dashed lines represent a ±20% deviation.

Numerical results for the simulation of these data sets are given in Table 7. For the experiments from Arena et al.,¹⁸ carried out with CO₂ as only carbon source, our model presents the lowest χ² (2.94), followed by Nestler (5.58) and Slotboom (5.93). Similarly to our experimental points, the model from Nestler has a lower error for methanol, and our model predicts CO better. For our previous experiments,⁵² from which most contain CO in the feed, the present model has again the lower χ² (28.8). As already discussed in the parity plots, the refitted model from Nestler presents a poor performance for CO-containing feed, probably because K_CO was overestimated.

Table 7 Calculated errors for the simulation of external datasets

			MSE 10^3			MRE 10^2
Experiments	Model	χ ^2	CO	CO2	CH3OH	CO	CO2	CH3OH
Arena et al.^18	This work	2.94	36.1	0.11	61.8	14.7	0.82	20.7
	Nestler	5.58	138.3	0.18	47.4	30.2	0.88	18.0
	Slotboom	5.93	110.8	0.35	86.6	26.6	1.25	24.1
Rodrigues Niquini et al.^52	This work	28.8	16.3	0.93	84.4	6.83	2.28	22.7
	Nestler	233.8	63.8	1.04	761.3	14.5	2.60	82.2
	Slotboom	36.0	53.0	1.19	72.9	12.8	2.66	21.1

---

Conclusions

The reaction kinetics of CO₂ hydrogenation to methanol was investigated on a Cu/ZnO/ZrO₂ catalyst under process conditions relevant for scale-up applications in the context of Power-to-X technologies. A new 6-parameter kinetic model was developed based on 500 steady-state experiments, which, to the best of our knowledge, is the broadest validity range for this catalyst composition. As a contribution to other scientists in this field, the experimental data are available in the SI. The presented model was compared with state-of-the-art models from the literature. All models provide a reasonable description of the reaction kinetics, whereat our model has the lowest weighted squared error for the carbon-containing compounds, even with fewer parameters. Our model is able to simulate appropriate sensitivity to variations in pressure, temperature and gas composition. The analyzed models were also used to simulated two external datasets, one of them containing also experiments with CO in the feed. Particularly under these conditions, the accuracy of the model presented here is significantly better, indicating its potential for plant simulation involving recycle. The developed model contributes to the understanding of the reaction kinetics for CO₂ hydrogenation to methanol on CZZ catalyst. The equations can easily be implemented into commercial software for model-based optimization and process scale-up. In this work, our measurements were performed in the kinetic regime, since our goal is the modeling of intrinsic kinetics. For this reason, heat and mass transfer limitations are absent. In real industrial applications, such effects will be present and have to be modeled accordingly. For this purpose, the model presented can be coupled with heat and mass transfer models, quantitatively predicting the effect of inlet compositions, operating conditions and changes in catalytic activity over time, leading to improved product quality and higher yields in larger scale production.

Author contributions

Gabriela Rodrigues Niquini: conceptualization, methodology, software, formal analysis, writing – original draft, visualization; Karla Herrera Delgado: conceptualization, resources, writing – review & editing, supervision, project administration, funding acquisition; Stephan Pitter: resources, writing – review & editing, supervision, funding acquisition; Jörg Sauer: conceptualization, writing – review & editing, supervision, project administration, funding acquisition.

Conflicts of interest

There are no conflicts to declare.

Data availability

The experimental data supporting this article have been included as part of the supplementary information (SI). Supplementary information is available. See DOI: https://doi.org/10.1039/d5re00330j.

Acknowledgements

We thank Dr. Lucas Warmuth for the catalyst synthesis, Siegbert Johnsen for the construction of the reactor setup and Tomas Vergara for the support during the experimental work and the fruitful discussions. The authors acknowledge the German Federal Ministry for Economic Affairs and Energy (BMWE) for the funding of the 3D-PROCESS project (reference number 03EN2065E). Support from the Helmholtz Association is also acknowledged.

Notes and references

1. J. Artz, T. E. Müller, K. Thenert, J. Kleinekorte, R. Meys, A. Sternberg, A. Bardow and W. Leitner, Chem. Rev., 2018, 118, 434–504.
2. S. G. Nnabuife, A. K. Hamzat, J. Whidborne, B. Kuang and K. W. Jenkins, Int. J. Hydrogen Energy, 2025, 107, 218–240.
3. H. Chu, Z. Huang, Z. Zhang, X. Yan, B. Qiu and N. Xu, Sep. Purif. Technol., 2024, 343, 127153.
4. A. Beck, M. A. Newton, L. G. A. van de Water and J. A. van Bokhoven, Chem. Rev., 2024, 124, 4543–4678.
5. S. A. Vali, J. Moral-Vico, X. Font and A. Sánchez, Catal. Lett., 2024, 154, 3157–3173.
6. B. Liang, J. Ma, X. Su, C. Yang, H. Duan, H. Zhou, S. Deng, L. Li and Y. Huang, Ind. Eng. Chem. Res., 2019, 58, 9030–9037.
7. F. Arena, K. Barbera, G. Italiano, G. Bonura, L. Spadaro and F. Frusteri, J. Catal., 2007, 249, 185–194.
8. S. Wild, S. Polierer, T. A. Zevaco, D. Guse, M. Kind, S. Pitter, K. Herrera Delgado and J. Sauer, RSC Adv., 2021, 11, 2556–2564.
9. L. Warmuth, M. Steurer, D. Schild, A. Zimina, J.-D. Grunwaldt and S. Pitter, ACS Appl. Mater. Interfaces, 2024, 16, 8813–8821.
10. B. Lacerda de Oliveira Campos, K. John, P. Beeskow, K. Herrera Delgado, S. Pitter, N. Dahmen and J. Sauer, Processes, 2022, 10, 1535.
11. B. Lacerda de Oliveira Campos, K. Herrera Delgado, S. Pitter and J. Sauer, Ind. Eng. Chem. Res., 2021, 60, 15074–15086.
12. F. M. Bagwan, P. Dongapure, A. A. Kulkarni and S. N. Vasireddy, Chem. Eng. Res. Des., 2024, 205, 79–90.
13. J.-F. Portha, K. Parkhomenko, K. Kobl, A.-C. Roger, S. Arab, J.-M. Commenge and L. Falk, Ind. Eng. Chem. Res., 2017, 56, 13133–13145.
14. J. Nyári, D. Izbassarov, Á. I. Toldy, V. Vuorinen and A. Santasalo-Aarnio, Energy Convers. Manage., 2022, 271, 116200.
15. G. Graaf, E. Stamhuis and A. Beenackers, Chem. Eng. Sci., 1988, 43, 3185–3195.
16. F. C. Marcos, F. M. Cavalcanti, D. D. Petrolini, L. Lin, L. E. Betancourt, S. D. Senanayake, J. A. Rodriguez, J. M. Assaf, R. Giudici and E. M. Assaf, Chem. Eng. J., 2022, 427, 130947.
17. S. Khawaja and M. R. Usman, Int. J. Chem. Kinet., 2024, 56, 469–481.
18. F. Arena, G. Mezzatesta, G. Zafarana, G. Trunfio, F. Frusteri and L. Spadaro, J. Catal., 2013, 300, 141–151.
19. N. Park, M.-J. Park, Y.-J. Lee, K.-S. Ha and K.-W. Jun, Fuel Process. Technol., 2014, 125, 139–147.
20. M. Dong, J. Ning, H. Liu, J. Xiong, J. Yang, Z. Huang, Y. Liang and J. Lu, Int. J. Green Energy, 2024, 21, 3573–3587.
21. T. Kubota, I. Hayakawa, H. Mabuse, K. Mori, K. Ushikoshi, T. Watanabe and M. Saito, Appl. Organomet. Chem., 2001, 15, 121–126.
22. K. Bussche and G. Froment, J. Catal., 1996, 161, 1–10.
23. A. Burcat and B. Ruscic, Third Millennium Ideal Gas and Condensed Phase Thermochemical Database for Combustion with Updates from Active Thermochemical Tables, Argonne national laboratory technical report, 2005.
24. E. Goos, A. Burcat and B. Ruscic, Extended Thermochemical Database for Combustion in the Third Millennium, Including Updates from Active Thermochemical Tables, 2021, https://respecth.elte.hu/burcat/THERM.DAT, Accessed on 06-Feb-2025.
25. F. Nestler, A. Schütze, M. Ouda, M. Hadrich, A. Schaadt, S. Bajohr and T. Kolb, Chem. Eng. J., 2020, 394, 124881.
26. T. A. Henkel, PhD thesis, Technische Universität München, 2011.
27. F. Studt, M. Behrens, E. L. Kunkes, N. Thomas, S. Zander, A. Tarasov, J. Schumann, E. Frei, J. B. Varley, F. Abild-Pedersen, J. K. Nørskov and R. Schlögl, ChemCatChem, 2015, 7, 1105–1111.
28. N. D. Nielsen, A. D. Jensen and J. M. Christensen, J. Catal., 2021, 393, 324–334.
29. Y. Slotboom, M. Bos, J. Pieper, V. Vrieswijk, B. Likozar, S. Kersten and D. Brilman, Chem. Eng. J., 2020, 389, 124181.
30. L. C. Grabow and M. Mavrikakis, ACS Catal., 2011, 1, 365–384.
31. C. Seidel, A. Jörke, B. Vollbrecht, A. Seidel-Morgenstern and A. Kienle, Chem. Eng. Sci., 2018, 175, 130–138.
32. B. Lacerda de Oliveira Campos, K. Herrera Delgado, S. Wild, F. Studt, S. Pitter and J. Sauer, React. Chem. Eng., 2021, 6, 868–887.
33. F. Studt, M. Behrens and F. Abild-Pedersen, Catal. Lett., 2014, 144, 1973–1977.
34. S. Polierer, D. Guse, S. Wild, K. Herrera Delgado, T. N. Otto, T. A. Zevaco, M. Kind, J. Sauer, F. Studt and S. Pitter, Catalysts, 2020, 10, 816.
35. D. Guse, S. Polierer, S. Wild, S. Pitter and M. Kind, Chem. Ing. Tech., 2022, 94, 314–327.
36. M. L. Schulte, V. Truttmann, D. E. Doronkin, L. Baumgarten, A. Nicolai, D. A. M. Beltran, F. J. Summ, C. Kiener, L. Warmuth, S. Pitter, E. Saraçi and J.-D. Grunwaldt, Angew. Chem., Int. Ed., 2025, 64, e202423281.
37. D. Guse, L. Warmuth, M. Herfet, K. Adolf, T. A. Zevaco, S. Pitter and M. Kind, Catalysts, 2024, 14, 517.
38. L. Warmuth, T. A. Zevaco and S. Pitter, Inorg. Chem. Commun., 2025, 172, 113753.
39. G. James, D. Witten, T. Hastie and R. Tibshirani, An Introduction to Statistical Learning with Applications in R, Springer, 2021.
40. L. L. Raja, R. J. Kee, O. Deutschmann, J. Warnatz and L. D. Schmidt, Catal. Today, 2000, 59, 47–60.
41. P. Weisz and C. Prater, Adv. Catal., 1954, 6, 143–196.
42. M. A. Vannice, Kinetics of Catalytic Reactions, Springer, 2005.
43. O. Levenspiel, Chemical Reaction Engineering, John Wiley, Sons, 1999.
44. H. S. Fogler, Elements of chemical reaction engineering, Prentice Hall, 4th edn, 2006.
45. D.-Y. Peng and D. B. Robinson, Ind. Eng. Chem. Fundam., 1976, 15, 59–64.
46. N. Delgado Otalvaro, M. Kaiser, K. Herrera Delgado, S. Wild, J. Sauer and H. Freund, React. Chem. Eng., 2020, 5, 949–960.
47. L. Meng and Y.-Y. Duan, Fluid Phase Equilib., 2005, 238, 229–238.
48. L. Meng, Y.-Y. Duan and X.-D. Wang, Fluid Phase Equilib., 2007, 260, 354–358.
49. U. K. Deiters, Fluid Phase Equilib., 2013, 352, 93–96.
50. A. Prašnikar, D. L. Jurković and B. Likozar, Appl. Catal., B, 2021, 292, 120190.
51. S. Kuld, M. Thorhauge, H. Falsig, C. F. Elkjær, S. Helveg, I. Chorkendorff and J. Sehested, Science, 2016, 352, 969–974.
52. G. Rodrigues Niquini, B. Lacerda de Oliveira Campos, K. Herrera Delgado, S. Pitter and J. Sauer, Chem. Eng. J., 2024, 480, 147968.
53. T. C. Lin and A. Bhan, J. Catal., 2024, 429, 115214.
54. A. Karelovic, G. Galdames, J. C. Medina, C. Yévenes, Y. Barra and R. Jiménez, J. Catal., 2019, 369, 415–426.

---