Forced periodic operations of a chemical reactor for methanol synthesis – The search for the best scenario based on Nonlinear Frequency Response Method. Part I Single input modulations

Abstract In this two-part paper a comprehensive study of the potential to improve performance criteria of a methanol synthesis reactor through forced periodical operations is presented. The study uses the Nonlinear Frequency Response method, a powerful analytical and approximate tool which gives an answer whether and under which conditions certain periodic operation would lead to improvement of process performance. To demonstrate the method, isothermal and isobaric methanol synthesis in a lab-scale CSTR is considered. In Part I, the analysis is performed for single input modulations. Partial pressures of each reactant in the feed stream and the total inlet volumetric flow-rate are considered as possible modulated inputs. The results show that modulations of single inputs essentially do not provide potential for significant improvements. In Part II, the study will be extended to analysis of periodic operations with simultaneous modulations of two inputs and conditions offering significant performance enhancements will be identified.


Introduction
Operating processes in forced periodic mode is one way of Process Intensification (PI), which represents a set of innovative principles in process or equipment design leading to significant process improvements (Van Gerven and Stankiewicz, 2009). The term forced periodic operation refers to a case when one or more inputs of a system are periodically modulated around their corresponding steady-state value(s) (Petkovska and Seidel-Morgenstern, 2013). In chemical engineering, the standard way to design and operate continuous processes is based on the optimal steady-state design and a control system which keeps all outputs as close as possible to theirs optimal steady-state values. Nevertheless, it is a known fact that perturbing the system periodically can sometimes result in better performance than the optimal steady-state operation (Silveston and Hudgins, 2013). Forced periodic operations of chemical reactors have been of interest for many research groups worldwide, involved in numerous theoretical, numerical or experimental studies (Armstrong and Teixeira, 2020;Bailey 1973;Douglas and Rippin 1966;Douglas 1967;Douglas 1972;Renken 1972;Schadlich et al. 1983;Silveston 1987;1998;Silveston and Hudgins, 2013;Sterman and Ydstie 1990a, 1990b, 1991Chen et al., 1994). These studies showed that the time-average indicators of chemical reactor performance, such as conversion, selectivity, production rates, productivity, could be improved by implementing forced periodic operations.
Although the physical reasons can be different, it could be said that the process improvement owing to forced periodic operations is always a consequence of process nonlinearity (Petkovska and Seidel-Morgenstern, 2013). Also, it is important to know that the resulting performances could be improved, deteriorated or unchanged, in comparison to the steady-state performance (Douglas and Rippin, 1966). Considering that there are many ways to apply forced periodic operations, e.g.: different inputs which can be modulated, different shapes of the modulated input(s), different values of the forcing parameters (amplitude(s), frequency, phase difference, etc.), it is a challenging task to find the mode which would lead to the highest improvement (Parulekar, 2003;Silveston et al., 1995). It is therefore of economic importance to carry out theoretical studies for assessing the effects of forced periodic operations of chemical processes, before any experimental studies (Chen et al., 1994). In our previous work (Marković et al, 2008;Nikolić, 2016;Petkovska and Seidel-Morgenstern, 2013;Petkovska et al., 2018) we introduced the Nonlinear Frequency Response (NFR) method as a reliable analytical tool for evaluating possible improvements and finding the best forcing parameters (Nikolić, 2016;Živković et al., 2020b).
The Nonlinear Frequency Response (NFR) method is based on the analysis of the frequency response of weakly nonlinear systems (Nikolić, 2016;Petkovska and Seidel-Morgenstern, 2013;Petkovska et al., 2018). It is a general, mathematically based theoretical and approximate method which answers the following questions (Nikolić, 2016;Petkovska and Seidel-Morgenstern, 2013;Petkovska et al., 2018): • Can the process performance be improved by periodic input modulations or not?
 Which input(s) should be periodically modulated in order to improve the process performance?
• Which forcing parameters (amplitude(s) and frequency of the input modulations, as well as the phase difference in the case of simultaneous modulation of two inputs) should be used?
• What would be the extent of the possible improvement?
Also recently, the NFR method was further advanced into the so-called computer-aided Nonlinear Frequency Response (cNFR) method, by developing a user friendly software application for implementing the NFR method, making its application much easier (Živković et al., 2020a). Another recent development was establishing and implementing a new methodology for optimizing forced periodic operations, combining the cNFR approach and multi-objective optimization the one step optimization (Živković et al., 2020b). The main advantage of this approach is that using the NFR approach leads to objective functions which are defined as algebraic expressions of all optimization parameters, which drastically shortens the needed computing time. Also, the optimization is performed in a single step, meaning that all optimization parameters: the steady state point around which the forced periodic operation should be performed and the forcing parameters (frequency, amplitudes and phase difference) are determined at the same time.
In this paper, the NFR method is used for analysing the potential of forced periodical operations of a chemical reactor in which the catalytic methanol synthesis from synthesis gas occurs.
Methanol is an important basic chemical which is produced on large scale in chemical industry and used as starting material for production of paraffins, olefins and other organic chemicals, as well as fuel or fuel additives (Fiedler et al., 2000). The important fact is that methanol can be used as an energy carrier (Olah, 2004). Excess electrical energy from renewable resources (wind or solar) can be used to generate hydrogen which, combined with CO and/or CO 2 (from biomass or agriculture waste streams) in the presence of a suitable catalyst, can be converted into methanol as a medium for chemical energy storage (Larsten and Sonderberg, 2013;Martin, 2016;Olah, 2005;Raeuchle et al., 2016;. However, methanol production for energy storage using renewable resources deals with unavoidable fluctuations in the supply of hydrogen, CO and/or CO 2 . In some cases, these fluctuations could lead to improvement of methanol production. Previous experimental investigations showed that significant improvement in methanol production could be achieved by forced periodic operations (Chanchlani et al., 1992(Chanchlani et al., , 1994Silveston, 1987). Experimental results presented in (Chanchlani et al., 1992(Chanchlani et al., , 1994 showed that the improvement of about 35% increase of methanol production relative to steady state is possible (when H 2 and CO 2 in the feed stream are periodically modulated) for an isothermal packed bed reactor when Cu/ZnO and Cu/ZnO/Al 2 O 3 were used as catalyst at 225 o C and 2.86 MPa.
In this two-part manuscript we use the NFR method to perform a systematic search for the best periodic process which would intensify the process of methanol synthesis. In Part I of this work, the NFR analysis is performed for cases of single input modulations. In Part II, the cases of simultaneous modulations of two inputs will be analysed.

Nonlinear frequency response method for single input modulations
By definition, frequency response is the quasi-stationary response of a stable system to a periodic (sinusoidal or co-sinusoidal) input modulation around its steady-state value (Douglas, 1972).
Frequency response is obtained when the transient response becomes negligible (theoretically for infinite time). For linear systems, frequency response is a periodic function of the same shape and frequency as the input function, but with different amplitude, with a phase shift and the mean value which is equal to the steady-state one.
On the other hand, frequency response of a nonlinear system is a complex periodic function. For a weakly nonlinear system (Weiner and Spina, 1980), if the input x is modulated in a cosinewave form, with amplitude and frequency ω, around a steady-state value x s : after long enough (theoretically infinite) time, the output of the system would contain the basic harmonic (y I ) which has the same frequency as the input modulation, a non-periodic (the socalled DC) component (y DC ) and an infinite number of higher harmonics (y II , y III ,…) (Douglas, 1972;Weiner and Spina, 1980): where y(t) represents the output, y s its steady-state value, while B i and φ i are the amplitude and the phase shift of the i-th harmonic of the output, respectively.
In order to evaluate a forced periodic operation around a steady-state point, only the timeaverage value of the periodic steady-state response is of interest. Using equation (2) it is easily concluded that the DC (non-periodic) component of the frequency response equals the difference between the time-average and the steady-state value of the output. Using the concept of higher order FRFs, the DC component can be written as the following infinite series (Weiner and Spina, 1980): For weakly nonlinear systems, the significance of different terms in equation (3) decreases with the increase of the corresponding FRF order. As a consequence, the DC component can be approximated with its dominant term, which is proportional to the asymmetrical second order function and the square of the input amplitude (Marković et al., 2008): Equation (4) is the foundation of the NFR method for evaluating periodic operations with one modulated input. The sign of y,x,x (2) ( , − ) determines whether the periodic operation would be superior to the corresponding steady-state one, while its magnitude determines the possible improvement.
In this paper the focus is on the production of methanol from syngas (a mixture of CO, CO 2 and H 2 ) using a commercial Cu/ZnO/Al 2 O 3 catalyst. The overall reaction mechanism assumes the reactions of CO and CO 2 hydrogenation: СО 2 + 3Н 2 ⇄ СН 3 ОН + Н 2 О (6) and the reverse water-gas shift reaction (RWGS): (see e.g. Graaf et al. (1988)).

Kinetic model
The NFR analysis presented in this work is based on a mathematical model of the reactor incorporating a reaction kinetic model of methanol synthesis presented in Seidel et al. ( , 2020, which showed reasonable agreement with steady state and dynamic experimental data from the Ph.D. Thesis of Vollbrecht (2007). The model is based on a Langmuir-Hinshelwood mechanism which implies three main steps: adsorption of the reactants on the catalyst surface, reaction of the adsorbed species and desorption of the reaction products. Adsorption and desorption are assumed to be in equilibrium. Further, the model assumes three different active centres on the catalytic surface, i.e.
The fraction of the reduced surface centres was denoted by . Following Ovesen et al. (1997), changes in the catalyst morphology due to the oxidizing influence of CO 2 and H 2 O and the reducing influence of CO and H 2 were also taken into account and modelled with the following dynamic equation: In equation (8) it is assumed that the maximal value of the fraction of the reduced centres is limited to . In the current study the maximal value = 0.9 was used ). The equilibrium constants were fitted separately to the steady state data , whereas the dynamic rate constants k 1 + , k 2 + were fitted to dynamic data .
Finally, also an ideal gas phase was assumed. Catalyst deactivation, and further side reactions were neglected. With all of these assumptions the following lumped reaction rate expression were obtained :  For the reaction of CO hydrogenation (Eq. (5))  For the reaction of CO 2 hydrogenation (Eq. (6)) 2 = 2 2 ( 2 2 2 − 1 2 3 2 2 ) * 2 ⊗ 4 (12)  For the reverse water-gas shift reaction (RWGS) (Eq. (7)) The corresponding relative amounts of free active surface centres are given with the following expressions : The reaction rate constants k j were determined based on the modified Arrhenius equation: with the reference temperature T ref =523.15 K Vollbrecht, 2007) and j=1,2,3 corresponding to CO hydrogenation, CO 2 hydrogenation and reversed water-gas shift reaction, respectively.
The equilibrium constants of the chemical reactions defined by equations (5-7), as functions of temperature (Vollbrecht, 2007), are given in Appendix A.
The values of the kinetic parameters used in this paper are given in Table 1. These values are somewhat different than the ones reported in our previous publications , as they have been refitted to the experimental data of Vollbrecht (2007), by using ≤ 0.9 in the constraint set of the nonlinear least squares problem. In Table 1, the specific amount of surface centres q sat is also given.

Application of the NFR method to evaluate the potential of single input forced periodic operation of methanol synthesis reactor
In this Chapter, the NFR method is applied for the analysis and evaluation of possible improvement of methanol production, for forced periodic operations with single input modulations. The analysis is performed for a laboratory-scale uniformly mixed reactor such as the Micro-Berty reactor, which was used for kinetic measurements (Vollbrecht, 2007) on which the kinetic model used in this study is based. Thus, the theoretical results present below will also serve as a basis for a planned later experimental validation using this same reactor type.

Mathematical model
The mathematical model of the catalytic reactor for methanol synthesis is based on the following assumptions (the assumptions listed here are consistent with the assumptions used for the kinetic model presented above):  The reaction occurs in an isothermal and isobaric CSTR,  The gas phase is ideal in the range of operation parameters,  The adsorption equilibrium between the solid and the fluid phase exists,  The adsorption processes follow the Langmuir-Hinshelwood mechanism with the maximal adsorption capacity q sat (Table 1),  The catalyst deactivation can be neglected,  The reaction mechanism is defined with equations (5-7) and all other reactions can be neglected.
The Micro-Berty reactor can be modelled as a continuous stirred tank reactor (CSTR). As stated above, the reactor system in which methanol synthesis occurs is established for the case when total pressure (p tot ) is constant Considering that during methanol synthesis the total number of moles is decreasing and that the total pressure is held constant, the volumetric outlet flow-rate is also changing.
The mathematical model of the analysed system can be described with the following equations:  material balances for each component i  total material balance for the case when total pressure is held constant  the equation describing the catalyst dynamics The outlet volumetric flow-rate (̇) is evaluated based on the total material balance (Eq. (20)) which can be reformulated as follows: The adsorption equilibrium is described with the competitive adsorption Langmuir isotherm. The The dimensionless frequency (ω) (Tab.2) is defined based on the steady-state residence time ( 0, ) calculated using the steady-state inlet volumetric flow-rate (̇0 , ) and volume of the reactor (V G ), as follows: For applying the NFR analysis, all nonlinear terms in the mathematical model need to be given in the polynomial form or expanded in Taylor series around a previously established steady-state point (Petkovska and Seidel-Morgenstern, 2013;Nikolić, 2016;Petkovska et al., 2018).
Therefore, the nonlinear terms (reaction rate expressions ) from the mathematical model ) are replaced by their Taylor series expansions, which are given in Appendix C.
After incorporating the dimensionless variables (Tab. 2) in the mathematical model equations ) + ⋯ The auxiliary coefficients q, Q, , U, s, S as well as the auxiliary parameters E 1 -E 4 used in this dimensionless mathematical model (Eqs. (23-25) are given Appendices C and D.

Inputs, outputs and frequency response functions (FRFs)
The inputs which can be modulated for the analysed system are:  partial pressure of CO 2 in the feed stream,  partial pressure of CO in the feed stream,  partial pressure of H 2 in the feed stream and  total volumetric flow-rate of the feed stream.
The outputs of the analysed system are the following variables:  the partial pressures of all components in the outlet stream,  the fraction of reduced active surface centres of the catalyst in the reactor (which in fact represents the state of the catalyst in the reactor) and  the volumetric flow-rate of the outlet stream.
The vectors of inputs X and outputs Y in the dimensionless form are defined, as follows: The FRFs which correlate an output y (y=1,…,7) with a modulated input x (x=1,…,4) will be denoted as G-functions. For implementation of the NFR method for evaluating the potential forced periodic operations, it was necessary to derive:  the first order frequency response functions marked as , (1) ( ),  the asymmetrical second order frequency response functions marked as , , (2) ( , − ).
The G-FRFs were derived by implementing a standard derivation procedure which was given in our previous publications (Petkovska and Seidel-Morgenstern, 2013;Nikolić, 2016;.

3. Derivation of the FRFs
The periodic modulation of input , defined as a dimensionless inlet partial pressure of CO 2 , CO or H 2 (for x=1, 2 or 3) or dimensionless flow-rate (x=4), with a forcing frequency ω and forcing amplitude , in the shape of a co-sinusoidal function of frequency, is defined as follows: In the cases when the partial pressure of one of the reactants is the modulated input, the partial pressure of the inert (N 2 ) is adjusted in order to assure isobaric conditions (constant total pressure) in the reactor.
For the general case, when input X x is periodically modulated, the output Y y , based on the Volterra series (Volterra, 1959) can be written in the following way: The solution of this matrix equation gives the matrix of all first order FRFs for all combinations of outputs and inputs: By collecting the non-periodic terms with (( 2 ) 2 0 ), a set of linear algebraic equations defining the asymmetrical second order G-FRFs is obtained, which are again written in the matrix form and given with Eq. (31) The definitions of the coefficients used in equations (29) and (30) are given in Appendix E.
The definitions of the coefficients used in equation (31) and (32) are given in Appendix F.

Identification and evaluation of regions of possible improvement
The main goal of implementing forced periodic operations is to improve the reactor performance, e.g. through increase of methanol production, conversion or yield. All these performance criteria can be evaluated based on the time-average outlet molar flow-rate of methanol, which has been chosen as the main indicator of possible improvement that should be maximized.
The methanol molar flow-rate can be evaluated from the methanol partial pressure and volumetric flow-rate of the outlet stream: It is convenient to use the dimensionless molar flow-rate of methanol, which is defined as a relative deviation from its steady-state value, in an analogous way as the dimensionless partial pressures ( The non-periodic (DC) component of the outlet molar flow-rate of methanol, which is the measure of improvement of methanol production, can be evaluated in the following way: The outlet molar flow-rate of methanol is an additional output which is of interest, which can be associated to additional sets of FRFs, which will be denoted as H-functions. If one of the inputs X x is modulated in a co-sinusoidal way, the DC component of outlet molar flow-rate of methanol can be approximately evaluated using the corresponding H ASO FRF: Based on equations (35 and 36) it is relatively easy to derive a relation between the H ASO FRFs and the previously derived G-FRFs, corresponding to the methanol partial pressure and the outlet flow-rate. The asymmetrical second order H-FRF is: ( )), = 1, 2, 3 or 4 Based on the NFR method, the mean (time-average) value of the outlet molar flow rate of methanol for co-sinusoidal modulation of input X x , can be approximately calculated using the following expression: is the outlet molar flow-rate of methanol in steady state, while (̇3 ) is the time average value of the outlet molar flows-rate over an integer number of periods P: Based on the mean value of the methanol outlet molar flow rate, several performance indicators were defined. One of them is the normalized methanol production rate per unit mass of catalyst for the periodic operation (PO): Other two performance indicators analysed are yield of methanol based of total carbon: and yield of methanol based on hydrogen: It should be noticed that for single input modulations, the mean values of the molar flow-rates of the reactants in the feed stream are identical to their steady-state values. Using this fact and equation (38), the yields defined in equations (42) and (43) can be evaluated based on their steady-state values and function 1, , ( , − ): and where and are the yields of methanol based on total carbon and based on hydrogen, respectively, corresponding to the chosen steady-state point.
Based on the sign of ASO H-FRF 1, , (2) ( , − ), it is possible to predict whether the improvement owing to periodic modulation of the input X x is possible at all, or not (Petkovska and Seidel-Morgenstern, 2013;Nikolić, 2016;Petkovska et al., 2018). The improvement can be achieved only if 1, , ( , − ) is positive.

Simulation results and discussion
In this Section, the simulation results based on the NFR analysis, for periodically operated isothermal, isobaric, lab-scale Micro-Berty reactor are given, for cases of single input modulations of the reactant partial pressures in the feed stream, or its volumetric flow-rate. The analysis was performed for a lab-scale reactor of the volume of the reaction mixture (i.e. gas phase) V G =10.3 ml and with a mass of catalyst =0.00395 kg.

Choosing the optimal steady-state for analysis
The first step in the analysis of forced periodic operation is to determine the optimal steady state, around which the system inputs should be modulated.
The optimal steady-state was chosen based on multi-objective optimization with two objective functions: normalized outlet molar flow-rate of methanol (mmol/(min kg cat )) and yield of methanol based on total carbon, which both need to be maximized. The multi-objective optimization problem was solved using ɛ-Constraint method (Haimes et al., 1971). The variables optimized were the mole fractions of all reactants (CO 2 , CO and H 2 ) in feed stream and the reactor temperature. The values of reactor pressure, the flow-rate of the feed stream and the mole fraction of the inert (N 2 ) were fixed. The optimization was performed in the range of validity of the kinetic model (Vollbrecht, 2007 (mole fractions of CO 2 and CO between 0 and 1, mole fraction of H 2 between 0.5 and 1 and temperature between 473 and 533 K).
The resulting Pareto front with the marked selected optimal steady-state point is given in Appendix G (Figure G.1). An overview of the optimization results for that selected steady-state point is given in Table 3.

Results for single input modulations around their optimal steady-state values
The simulation results of NFR analysis for single input modulations around optimal steady-state are presented here. The asymmetrical second order H-FRFs which correlate the outlet molar flow-rate of methanol to the 4 inputs related to the feed reactor steam (partial pressures of CO 2 , CO, H 2 and the volumetric flow-rate), are given in Figure 1, as a function of dimensionless forcing frequency. Based on results presented in Fig.1, it can be concluded that:  For single input modulations of CO 2 , CO and inlet volumetric flow-rate around the optimal stead-state, the H-ASO FRFs which correlate the outlet molar flow-rate of methanol to modulated inputs (H (2) 1,1,1 (ω,-ω), H (2) 1,2,2 (ω,-ω), H (2) 1,4,4 (ω,-ω)), are negative and tend to zero for high forcing frequencies. Consequently, periodic modulations of these inputs cannot improve the process of methanol synthesis.
 For single input modulation of H 2 partial pressure around its optimal steady-state value, the corresponding H-ASO FRF which correlates the outlet molar flow-rate of methanol to the modulated input (H (2) 1,3,3 (ω,-ω)) is positive for forcing frequencies higher than 0.55 and also tends to zero for high forcing frequencies. The maximal value of this ASO FRF of 0.049 is obtained for dimensionless forcing frequency ω≈1. The maximal possible increase of the normalized outlet molar flow rate of methanol corresponding to this is 0.13% (when the highest possible forcing amplitude is used), which is practically insignificant.

Analysis of maximal possible improvement for single input modulations
The results of the previous section show that, under the defined conditions, it is not possible to improve the reactor performance by periodic modulation of partial pressures of CO and CO 2 or the volumetric flow-rate around the chosen optimal steady-state point, as corresponding ASO H-  A short overview of these cases is given in Table 4. In this table, Table 4 show that the possible improvement is very small (the highest improvement is obtained for partial pressure of CO in the feed stream as the modulated input, and it is possible in the whole frequency range).
Nevertheless, despite the fact that some improvement is possible in comparison to the corresponding steady-states around which the four inputs should be modulated, all performance criteria corresponding to the periodic operations with modulation of the inlet partial pressures of CO 2 , CO, H 2 or the volumetric flow-rate of the feed stream, presented in Table 3, are by far worse than the performance criteria corresponding to the optimal steady-state defined in Section 5.1 (Table 3). So, the conclusion of this analysis is that none of these four periodic operations with single input modulations is acceptable.

Conclusions
The goal of this two-part manuscript is to present the results of a comprehensive study of the potential of using forced periodic operations in order to improve the performance of a chemical reactor for isothermal and isobaric methanol synthesis from syngas. Four potential forced periodic inputs are considered: partial pressures of all reactants (CO 2 , CO and H 2 ) in the feed stream and its total volumetric flow-rate. Because of the quantity of the obtained results, the manuscript needed to be split into two parts.
In this first part of the manuscript, only the results for single input modulations are shown. Here are the most important results of this analysis:  Periodic modulations of the partial pressures of CO 2 , CO and volumetric flow-rate of the feed stream, around the optimal steady state, always result with performance deterioration, instead of performance improvement, in the whole frequency range.
Therefore, such periodic operations are unacceptable.
 Periodic modulation of the partial pressure of H 2 around the optimal steady state could lead to improvement of the reactor performances for some forcing frequencies with the maximal possible improvement of 0.13% which is practically insignificant.
 It is possible to find some cases for which some limited improvement can be achieved with periodic modulations of the analysed inputs, if the inputs would be modulated around some other steady-state points. Nevertheless, all these cases correspond to reactor performances that are much worse than for the optimal steady-state process. Accordingly, these periodic operations are also unacceptable.
 Analysis of periodic operations with simultaneous modulation of two inputs is the next logical step, as it is a well-known fact that such operations have high potential for improvement (Petkovska and Seidel-Morgenstern, 2013;Felischak et al., 2021), owing to the cross-effect between the two modulated inputs, which can be easily adjusted by adjusting the phase difference between the two inputs. This analysis will be presented in the second part of our manuscript.
It is important to point out that frequency response functions for the four inputs, derived and presented here, are necessary for the analysis of the periodic operations with simultaneous modulation of two inputs, i.e. the analysis and results presented in Part II would not be possible without the results presented in Part I.
In both parts of our manuscript, the analysis was performed by using the nonlinear frequency response analysis. Even for a complex case, such as the reactor for methanol synthesis, with four potential modulated inputs and a large number of outputs, the NFR method was proven as a very useful and efficient tool for evaluating whether the reactor performance could be improved by using forced periodic operations, or not.

Appendix
Appendix A -Temperature dependence of equilibrium constants of chemical reactions [Vollbrecht, 2007] log (