Advanced Control for Hydrogen Pyrolysis Installations

: Today, hydrogen production plays an important part in the industry due to the increasing use of hydrogen in signiﬁcant domains, such as chemistry, transportation, or energy. In this paper, we aim to design a numerical control solution based on the thermodynamic analysis of the pyrolysis reactions for hydrogen production and to present novel research developments that highlight industrial applications. Beginning with the evaluation of the technological aspects for the pyrolysis chemical process, the paper studies the thermodynamic evaluation of the system equilibrium for the pyrolysis reactions set, to recommend an appropriate automatic control solution for hydrogen pyrolysis installations. The numerical control architecture is organized on two levels, a control level dedicated to key technological parameters, and a supervisory decision level for optimizing the conversion performances of the pyrolysis process. The data employed for modelling, identiﬁcation, control, and optimization tasks, were obtained from an experimental platform. The scientiﬁc results can be implemented on dedicated equipment, to achieve an optimal exploitation of the industrial pyrolysis process.


Introduction
The chemical and petrochemical industry is an area of interest for automatic control, and some of the most representative applications find their place in this domain. The chemical and petrochemical industry has experienced an unparalleled development through the evolution of technological equipment and production lines, as well as through the expansion of production capacity and use of automatic control tools. The pyrolysis installations are plants of high importance in the petrochemical industry. A multi-component set of the pyrolysis reactions is produced mainly by the cracking of hydrocarbons in the presence of steam [1][2][3]. The cracking reactor is the heart of this chemical complex process. For our research interest, hydrogen is the main component in the pyrolysis reactions, being used in significant industries, like chemical, transportation, or energy. Hydrogen production becomes mainstream in the future and the hydrogen energy technology is beginning to penetrate society. Projects of fuel cell vehicles have been carried out many times in Europe, America, and Asia. Fuel cell vehicles were sold, and fuel cell buses also began to run on ordinary roads. The production volume of these fuel cell vehicles will surely grow, and it seems that it will be important to mainstream automobiles in the future [4,5]. For example, more than 100 hydrogen stations have been built in Japan, and it is planned that 160 more will be built in 2020 and about 900 stations will be constructed by 2030 [6][7][8]. Furthermore, the use of hydrogen storage alloy is promoted as an energy stockpile. Hydrogen energy plants using storage alloys have been commercialized in recent years. Hydrogen storage is a technology that enables clean energy on a yearly basis, and it can be regarded as a distinguishing technology. In this way, the development of hydrogen energy technology has made remarkable progress [9][10][11].
The operating and security conditions required by pyrolysis installations for the hydrogen production involves the most recent solutions offered by automatic control tools [12,13]. For this reason, a study of the thermodynamic reactions is needed before proposing a control strategy.
This paper represents the extended version of the authors' work presented at the 2019 IFAC ACD conference in Bologna, Italy [14]. The newly added contributions are related to the extended evaluation of the state of the art in the domain, which justifies the related novelty and innovation. The interest in hydrogen production in the transportation industry has been justified, and the computing methods for process modelling, identification, control, and supervision were detailed. New simulations were added.
After a short introduction, in Section 2, the paper presents aspects related to the thermodynamic equilibrium analysis of the hydrogen pyrolysis reactions set and to the reproduction possibility of the results obtained on industrial scale. Section 3 proposes the appropriate control systems, designed for key parameters of the pyrolysis reactor. In Section 4, the supervisory system is designed, with the aim to optimize the pyrolysis reactor operating point. Finally, conclusions are given in Section 5.

Thermodynamics of the Pyrolysis Chemical Process
The system of pyrolysis reactions is oriented towards obtaining hydrogen at the thermodynamic equilibrium, when the maximum quantity of hydrogen is obtained. The composition of the reaction products is determined at the thermodynamic equilibrium using a numerical procedure. Ethane and water (vapors) are used as raw materials (reagents). Under constant temperature T and pressure P conditions, the system with simultaneous reactions involving a number N of components is characterized by the thermodynamic function-the enthalpy, from Relation (1) [1-3]: where n is a vector with n i the number of moles corresponding to each component; i and µ are the vector of the chemical potentials of each component involved in the reaction. At thermodynamic equilibrium, the enthalpy of the reaction has the minimum value and respects the conservation of the atomic mass constraints (number of atoms).
where a ji are the elements of the atomic matrix representing the number of atoms of the element j in component i of the system; b j is the atomic abundance of the chemical element j in the system (number of atoms of element j originally introduced into the system); and m is the number of chemical elements in the reaction system. For the equilibrium condition, the problem of nonlinear optimization with equality restrictions is formulated as: The system of pyrolysis reactions is described by enthalpy G, given by the following representation: where µ 0 i (T) is the standard chemical potential of component i, and R is the ideal gas constant. We consider the reaction system to obtain hydrogen as main product, in the following conditions: -reaction temperature: T = 1100 • K; -pression: P = 3.8 atm; ratio no.of moles H 2 O no.of moles C 2 H 6 = 4/1 = r; -reaction products: H 2 , CO, CO 2 , CH 4 , C 2 H 2 , O 2 , H 2 O, C 2 H 4 , C 2 H 6 (N = 9); -the atomic matrix in Table 1; the values of the standard chemical potentials in Table 2; The optimization problem for thermodynamic analysis in (3) is reformulated as follows: where I = G RT is the criterion function, x is the vector of component concentrations (x = n), h and g express the constraint functions from (3).
In order to solve the problem (5), a nonlinear programming method is proposed through penalty techniques [15,16], which will compute the optimal solution x*. The penalty technique turns (5) into the equivalent problem: min x∈C adm where F (k) is the flexible tolerance function calculated at step k; D(x (k) ) is the restriction violation (penalty) function.
In the minimization algorithm, F (k) is decreasing positively: The penalty criterion C(x) is defined as follows, using its quadratic form: with δ i the Heaviside operator: and the admissible domain C adm given by: This penalty technique works based on the concept of quasi-admissibility, defined as: For the case where the intermediate search point x(k) is not admissible (the inequalities from Relation (5) are not respected), the algorithm proceeds to minimize C(x) until (12) is met, and then the criterion function I(x) will be minimized respecting C adm .
The algorithm is working to minimize the penalty function C(x) and to reduce the tolerance function F (k) in the search process for the solution x*, ensuring corresponding values, such that C(x) respects (12).
The minimization for criterion functions I(x) and C(x) is performed with the evolutionary direct search algorithm BOX, recommended for solving the problems formulated in (6). The BOX method described in [15,16], constructs a COMPLEX geometric figure in the n-dimensional space, which moves in a controlled sequence towards the solution x*, by reflection, expansion, or contraction operations.
The flexible tolerance function F (k) is calculated at each search level as follows: where m is the number of equality restrictions, x where a is the side of the initial geometrical search figure.
Considering the pyrolysis reaction simulation at the parameters: T = 1100 K, P = 3.8 atm, the results from Table 3, with the hydrogen concentration at thermodynamic equilibrium, were obtained. Based on the results of this thermodynamic analysis for the pyrolysis process, a control solution for the cracking reactor to the industrial plant is proposed.

Nominal Control Level Design
By taking the results of the thermodynamic reaction analysis from Section 2 into account, a control solution for the pyrolysis reactor was proposed, in order to handle two major aspects concerning the operation of the chemical process: to maintain a good ratio between the quantities of reactants, steam and ethane, that are given to the installation, in a proportion 4:1, and to ensure the reaction conditions regarding the temperature and the pressure inside the reactor. Furthermore, we have considered a control structure which involves three closed loops, presented in Figure 1: FRC-1, control loop for steam/ethane flow ratio; TRC-2, control loop for temperature; PRC-3 control loop for pressure. For the control solution, RST polynomial control was recommended to assure the construction of the two-level proposed architecture, offering independent performances in tracking and control [12,[17][18][19].
The design of the polynomials R(z −1 ), S(z −1 ) and T(z −1 ) is performed in two steps. The closed-loop poles are defined by the desired characteristic polynomial, P(z −1 ), as regulation objective, and R(z −1 ) and S(z −1 ) polynomials are computed by the Diophantine equation: where the polynomialsÂ z −1 andB z −1 represent the discrete model of the process. The Equation (15) has a single unique solution for: The polynomial Equation (15) can be solved by means of the linear algebraic system: where M is the Sylvester matrix of dimensions (na + ns) × (na + ns) and where p i are the coefficients of the characteristic polynomial P(z −1 ), s i and r i are the coefficients of the S(z −1 ) and R(z −1 ) unknown polynomials, respectively. The solution of the algebraic Equation (17) is obtained by the inversion of the matrix M: The polynomial T(z −1 ) is determined in order to find the reference trajectory y * z −1 , the output of the trajectory generator, given by the pair (Â n z −1 ,B n z −1 ), and to assure the output of the global control system to reference trajectory, y z −1 = y * z −1 (see Figure 2).

Control of Stream/Ethane Flow Ratio
The first loop ensures the control of the steam flow to the hydrogen reactor, respecting the imposed ratio. Since steam flow respects the ratio imposed with ethane flow, it is thus achieved the control of both parameters. In order to determine the control algorithm for the steam flow which is pumped into the reactor, we have analytically determined the mathematical model based on the flow equations provided by the parameters of the technological pipe, from the experimental platform [18,19].
Let us consider the pipe with Q-flow rate, ∆P-differential pressure across the pipe, α-flow coefficient, S-pipe section, and ρ-fluid density. Considering the short pipe is equivalent to a hydraulic resistance for which the following relation is valid: Using the technological values for the steam flow process, after some calculus operations (linearization and normalization), we will arrive to the non-dimensional linearized model.
For the fluid flow in stationary state, the forces acting on the system are in equilibrium: where ∆P 0 S is the active force of pressure on the liquid and Q 2 0 ρ 2α 2 S 2 S is the reactive force due to the restriction.
In the dynamic state, the difference between both forces is compensated by the inertial force (time variation of quantity of movement): where M represents the mass of the fluid in the pipe and ν is its velocity. Therefore, we have: In general, we work in stationary state and we consider small variations around an operating point, which translates into: From Equation (23) and Relations (24) and (25), we get: By extracting the stationary regime expressed by (26), we get: In this equation, which corresponds to the non-linear knowledge model, ∆p(t) represents the variation ∆(∆P(t)) of the differential pressure.
By ignoring the quadratic term ∆Q 2 (t) in Equation (27), and using the normalized values for Q(t) and p(t) with respect to stationary state, with the notations: the final model becomes: By including the actuator and the flow sensor [12,20] for a sampling period of 0.5 s, the discrete dynamic model of the steam flow process has been obtained: Taking into consideration this model, the controller was designed to ensure the desired performance using the RST polynomial control algorithm, represented in Figure 2.
The control polynomials R, S and T were computed using the poles allocation method [20]. For the ratio control system, we have obtained the control polynomials in (31): To validate the performances in closed loop, the Matlab/Simulink simulation environment was used. The results of the simulation are shown in Figure 3. It can be noticed that the RST controller ensures in closed loop both tracking and regulation performances for the steam/ethane ratio system [12,21].

Control of Reaction Parameters
It is of utmost importance that the values of the temperature and the pressure are maintained within the boundaries of the admissible operating domain so that the process functions properly and the reactor walls are not submitted to any risk of degradation. In this case, the temperature within the median section of the reactor is controlled by the CH 4 combustion flow that heats the reactor [22]. The automatic control solution should be able to provide the possibility of maintaining the temperature value in the (1100-1150 K) interval and the pressure value, in the interval (3.7-4.4 atm). By controlling the temperature and the pressure inside the reactor, the reaction equilibrium state is ensured.
In this case, the proposed solution must attain two tasks: the data acquisition and identification of the dynamic models from the experimental pyrolysis platform and, respectively, the control system design of the major parameters of the process. The physical data acquired from the plant are used for the identification of the mathematical models. After the identification process was completed, by using the Least Squares Recursive Experimental Method (LSREM), the models were validated using the whiteness test [12]. A dedicated software for experimental identification, WINPIM [20], was used.
The model identified and validated for the reactor temperature is: Moreover, the model identified and validated for the reactor pressure is: The controllers were designed to ensure the desired tracking and regulation performances, using RST polynomial control algorithms. The discrete control polynomials R, S, T were computed by poles allocation methods using a model-based control design software WINREG [20].
For the temperature control system, we have obtained the control polynomials in (34): (34) Figure 4 shows the dynamic response of the closed loop temperature control system when a 20% step disturbance is added. For the pressure control system, we have obtained the control polynomials in (35): Similarly, Figure 5 shows the response of the closed loop pressure control system when a 20% step disturbance is added. For all three systems (FRC-1, TRC-2, PRC-3), the R, S polynomials were computed solving the polynomial Diophantine equation and the T polynomial was computed by assuring the tracking performances for the closed loop systems [18,19].
The RST controllers ensure that the desired tracking and regulation performances for the closed loop systems are met.

Robust Control Solution
The main goal is to design the robust control systems for temperature and pressure that can provide a way of preserving the performances obtained in simulation on the physical pyrolysis reactor in the presence of disturbances (ethane quality) or modeling uncertainties [18][19][20].
The robustness design is based on the disturbance-output sensitivity function of closed loop system. The sensitivity function for the polynomial RST systems (computed) defined as: whereÂ z −1 andB z −1 are the two polynomials of the identified process model, and R z −1 and S z −1 are the control polynomials computed earlier. It is desirable for the maximum value of the output sensitivity function to be as low as possible.
For each robust control system, we imposed a usually maximum of 6 dB for the output sensitivity function.
The sensitivity function of the previously designed nominal temperature control system does not meet this requirement, as shown in Figure 6. To improve the controller, a pair of complex auxiliary poles was added to the characteristic polynomial, given by the denominator of the sensitivity function. The new robust control polynomials are: The sensitivity function of the new closed loop system meets the imposed 6 dB limit (1.7 dB), which makes the new controller a robust one, see Figure 7. For the pressure control system, we have used a similar approach, since we encountered a similar situation, see Figure 8. We added the robust correction, by adding a pair of real stable poles to the characteristic polynomial, and obtained the following control polynomials: The new controller is a robust one, as illustrated in Figure 9. To summarize, the two control systems achieve the desired robustness margin while preserving the nominal performances.

Supervisory Level Design
The optimal operating point of the process can be determined by solving an optimization problem [15,16]. The solution of this problem represents the optimal decision for the hydrogen pyrolysis installation. In this case, the objective is to maximize the hydrogen concentration from the reactor output products. In order to solve this problem, we used a number of experimental data collected directly from the experimental platform: steam flow (y 1 ) m 3 /h , ethane flow (y 2 ) m 3 /h , reactor pressure (y 3 ) [atm], reactor temperature (y 4 ) [K] and hydrogen concentration (z(y)) [%]. The data sets are represented in Table A1.
The multivariable nonlinear structure of the model, was evaluated by: represents the parameters vector.
The structure of the model was determined after a sensitivity test of the output model z in relation to the inputs y. The parameter vector was estimated using the Least Squares Method, resulting the following values:â 0 = −7.228894,â 1 = 0.003917,â 2 = 2192.393311,â 3 = 0.082378,â 4 = 0.00003.
Considering the criterium function J(y) = z(y), the optimization problem is: To solve the problem in (34), the BOX algorithm was used, as it is recommended for nonlinear criterion functions with explicit constraints [15,16,20]. The BOX method generates, in the search space, the COMPLEX geometric figure that initially respects the restrictions in (42). The algorithm calculates the values of the criterion function J(y) and advances towards the maximum solution through controlled operations of reflection, expansion and contraction in the admissible domain, without resorting to the calculation of the gradient of the criterion function.
The optimal solution, obtained from a dedicated software (SISCON) [20] is within the imposed limitations: The resulting optimal steam/ethane flow ratio is r * = 1954.595/443.219 = 4.41.
The quantity of hydrogen produced in this caseĴ(y * ) = 43.6575% represents the maximum value from the total of the pyrolysis reaction products. By implementing this procedure, the hydrogen concentration has been increased by 2.14% in relation to the mean value of 41.522% of the data sets, which corresponds to an approximate of 10.14 [m 3 /h] of hydrogen.
The optimal decision (y * ) from this supervisory level is automatically transferred as setpoints to the control systems [23][24][25]. The research area can be extended to the estimation of the optimal combustion regime for the heating reactor, aiming to decrease industrial pollution and to introduce the stochastic optimization approach for an efficient exploitation of the pyrolysis installation [26].

Conclusions
A modern control solution for hydrogen pyrolysis installation, an important part of the petrochemical platform, was designed, based on robust and supervisory control. The RST control structure was recommended to assure the demands imposed by the two-level architecture, assuring the independent performances in tracking and control. A set of dynamic models for process control systems and the robust control polynomial RST algorithms were computed to support the implementation.
An optimal operating point was computed to maximize the hydrogen production, using mathematical nonlinear programming techniques.
The simulation work was supported by dedicated software: WINPIM for identification, WINREG for control design and SISCON for optimization. The experimental data for the computational support were obtained from a petrochemical experimental platform.
The results obtained in this study can be implemented in technological and numerical architecture configuration, using dedicated equipment for process control to obtain the optimal exploitation on an industrial pyrolysis platform.

Conflicts of Interest:
The authors declare no conflict of interest.