Thermodynamic Optimization of the Ethylene Oligomerization Chemical Process

The use of olefin oligomerization in the synthesis of liquid fuel has broad application prospects in military and civil fields. Here, based on finite time thermodynamics (FTT), an ethylene oligomerization chemical process (EOCP) model with a constant temperature heat source outside the heat exchanger and reactor pipes was established. The process was first optimized with the minimum specific entropy generation rate (SEGR) as the optimization objective, then multi-objective optimization was further performed by utilizing the NSGA-II algorithm with the minimization of the entropy generation rate (EGR) and the maximization of the C10H20 yield as the optimization objectives. The results showed that the point of the minimum EGR was the same as that of SEGR in the Pareto optimal frontier. The solution obtained using the Shannon entropy decision method had the lowest deviation index, the C10H20 yield was reduced by 49.46% compared with the point of reference and the EGR and SEGR were reduced by 59.01% and 18.88%, respectively.


Introduction
Ethylene oligomerization is a reaction that can transform ethylene into nitrogen-free and sulfur-free crude oil, which can be used as clean liquid fuel. It represents a new method for obtaining energy and has good application prospects. However, the popularization and application of this technology are limited by the bottlenecks of high energy consumption and low production rate, which still need further analysis and optimization. The traditional ethylene oligomerization reaction is mainly based on classical thermodynamics and chemical reaction kinetics. Wu et al. [1] calculated the reaction heat, Gibbs free energy and reaction equilibrium constant of an ethylene oligomerization system with a temperature range between 298 K and 700 K, and the influences of temperature and pressure on the equilibrium constant were also studied. Yuan et al. [2] optimized the process conditions of ethylene oligomerization and established the reaction kinetic model. Toch et al. [3][4][5] analyzed the reaction kinetics of ethylene oligomerization and established a one-dimensional heterogeneous industrial reactor model. Chen et al. [6] modeled an olefin oligomerization reactor, which they optimized for the minimum EGR. However, the development of this technology cannot be satisfied by studies only on catalysts and reactors.
In engineering, a chemical process not only includes chemical reactors, but also other components, such as mixers, separators, etc. Thus, the optimization of the whole chemical process is much more reasonable. For a single component, Andresen et al. [7,8] analyzed the minimum EGR of the heat exchange system and studied the mechanical separation process through FTT. Kingston et al. [9] analyzed the EGR of the separation process of air in a low-pressure distillation column. Badescu [10] took ammonia decomposition rate and heat flux as the objective functions and optimized ammonia decomposition reactors by changing the reactor wall temperature, diameter and catalyst particle distribution. Li et al. [11] and heat flux as the objective functions and optimized ammonia decomposition reactors by changing the reactor wall temperature, diameter and catalyst particle distribution. Li et al. [11] and Kong et al. [12] established a steam methane reforming reactor and SO3 decomposition membrane reactor through FTT and optimized them with the minimum EGR. In terms of using the FTT theory [13][14][15][16] for process optimization, Røsjorde et al. [17] studied the process model for the dehydrogenation of propane and optimized the external heat source temperature of the reactor with the minimum EGR as the optimization objective. Kingston and Razzitte [18,19] established a dimethyl ether synthesis process modelwhich included a constant temperature reactor, a compressor with different adiabatic efficiencies and a heat exchanger with a constant temperature difference-solved the minimum EGR and the corresponding optimal operating parameters [18], and further studied the coupling of the two reactors [19]. On the basis of Ref. [6], an EOCP model including a mixer, compressor, heat exchanger and reactor will be established using FTT in this paper; the process will be optimized with the minimum SEGR as the first objective of optimization, and a multi-objective optimization will be further performed utilizing the NSGA-Ⅱ algorithm [20][21][22][23] with minimum EGR and maximum C10H20 yield as the optimization objectives.

Physical Model of EOCP
In industry, raw gas is often stored in a compressed state in gas storage tanks, and the compression of the gas from low pressure to high pressure consumes energy. In order to fully consider the energy consumption of the EOCP, in this study, the feed gas was treated with an ambient temperature (298.15 K) and ambient pressure (0.101 MPa) as the initial state.
The EOCP model is divided into four parts, including the mixer, compressor, heat exchanger and reactor, as shown in Figure 1. According to the reaction equation, the mixed gas used in the process was composed of C2H4, C4H8, C10H20 and N2. According to the Peng-Robinson equation [24], the compression factor of the mixed gas was 0.9948, with a deviation of 0.52% from the ideal gas. Therefore, it could be assumed that the mixed gas was an ideal gas. Firstly, the raw reaction materials C2H4 and C4H8 and the protective gas N2, which does not participate in the reaction, enter the mixer and are combined fully. Secondly, the mixed gas reaches the set temperature and pressure through the action of the heat exchanger and compressor, respectively. Finally, the mixed gas enters the reactor for the chemical reaction. In this process, the working fluid follows the ideal gas equation of state: where R is the molar gas constant in − − ⋅ ⋅ 1 1 J mol K ; p is the pressure in MPa ; n is the amount of substance in mol ; and T is the temperature of the mixed gas in K . The reaction gas flows through each component and complies with the law of mass conservation: According to the reaction equation, the mixed gas used in the process was composed of C 2 H 4 , C 4 H 8 , C 10 H 20 and N 2 . According to the Peng-Robinson equation [24], the compression factor of the mixed gas was 0.9948, with a deviation of 0.52% from the ideal gas. Therefore, it could be assumed that the mixed gas was an ideal gas. Firstly, the raw reaction materials C 2 H 4 and C 4 H 8 and the protective gas N 2 , which does not participate in the reaction, enter the mixer and are combined fully. Secondly, the mixed gas reaches the set temperature and pressure through the action of the heat exchanger and compressor, respectively. Finally, the mixed gas enters the reactor for the chemical reaction. In this process, the working fluid follows the ideal gas equation of state: where R is the molar gas constant in J · mol −1 · K −1 ; p is the pressure in MPa; n is the amount of substance in mol; and T is the temperature of the mixed gas in K. The reaction gas flows through each component and complies with the law of mass conservation: . m out are the inlet and outlet mass flow rates of each component, respectively, in kg/s.

Physical Model of Mixer
It is assumed that the gas mixing process is isothermal and isobaric; thus, the EGR of the mixer can be calculated according to the entropy change rate between the inlet and the outlet, which can be expressed as: where . n k is the molar flow rate of substance k in mol · s −1 and y k is the mole fraction of each component substance.

Physical Model of Compressor
It is assumed that the compression process is irreversible adiabatic compression without a chemical reaction. The irreversibility of the compression process can be reflected by the efficiency of the compressor as follows: where h is the molar enthalpy of the working fluid in J/mol; p in and p out are the inlet and outlet pressures of the compressor in MPa, respectively; T isen out is the outlet temperature of compressor for the ideal reversible process in K; and T out is the outlet temperature of compressor for the real irreversible process in K.
The entropy flow rates at the outlet and inlet of the compressor can be calculated according to the following formula: where C p,k is the constant pressure molar heat capacity of the substance k in J · mol −1 · K −1 . The EGR can be obtained according to the entropy change rate between the inlet and the outlet. The relevant values can be checked according to the chemical manual. The EGR of the compressor is equal to the entropy change rate of the working fluid between the compressor inlet and outlet:

Physical Model of Heat Exchanger
The heat exchange process is irreversible. The heat exchanger model is shown in Figure 2. are the inlet and outlet mass flow rates of each component, respectively, in kg/s.

Physical model of mixer
It is assumed that the gas mixing process is isothermal and isobaric; thus, the EGR of the mixer can be calculated according to the entropy change rate between the inlet and the outlet, which can be expressed as: where  k n is the molar flow rate of substance k in − ⋅ 1 mol s and k y is the mole fraction of each component substance.

Physical Model of Compressor
It is assumed that the compression process is irreversible adiabatic compression without a chemical reaction. The irreversibility of the compression process can be reflected by the efficiency of the compressor as follows: where h is the molar enthalpy of the working fluid in J/mol; in p and out p are the inlet and outlet pressures of the compressor in MPa, respectively; isen out T is the outlet temperature of compressor for the ideal reversible process in K; and out T is the outlet temperature of compressor for the real irreversible process in K.
The entropy flow rates at the outlet and inlet of the compressor can be calculated according to the following formula: 0 , C , i n , i n 5 ( d ) l n 1.01325 10 Pa where ,k p C is the constant pressure molar heat capacity of the substance k in J mol K . The EGR can be obtained according to the entropy change rate between the inlet and the outlet. The relevant values can be checked according to the chemical manual.
The EGR of the compressor is equal to the entropy change rate of the working fluid between the compressor inlet and outlet:

Physical Model of Heat Exchanger
The heat exchange process is irreversible. The heat exchanger model is shown in Fig  (2) There is no chemical reaction in the heat exchange process. The heat exchange process follows the law of energy conservation, which can be expressed as: where q HE = U HE (1/T HE − 1/T HE,a ), i.e., the heat transfer process is assumed to follow the linear phenomenological heat transfer law, and U HE is the heat transfer coefficient in W · K · m −2 ; T HE and T HE,a are the temperatures of the reactant and the heat source, respectively, in K; d HE,i is the inner diameter of the heat exchanger in m; and F k is the molar flow rate of substance k in mol · s −1 .
Then, the local EGR can be expressed by the following formula [25,26]: Finally, the EGR of the heat exchanger can be expressed as: where L HE is the length of the heat exchanger in m.

Physical Model of Reactor
A reactor is a piece of equipment that is used to realize the reaction process. In this study, a multi-tube fixed bed reactor was used. The operating conditions of the reaction tubes in each tube bundle are the same; thus, a single tube is taken as an example in this paper. The model of the reactor is shown in Figure 3. (1) The model of gas in the heat exchange tube is a one-dimensional steady-state model; (2) There is no chemical reaction in the heat exchange process. The heat exchange process follows the law of energy conservation, which can be expressed as: where HE Then, the local EGR can be expressed by the following formula [25,26]: Finally, the EGR of the heat exchanger can be expressed as: where HE L is the length of the heat exchanger in m .

Physical Model of Reactor
A reactor is a piece of equipment that is used to realize the reaction process. In this study, a multi-tube fixed bed reactor was used. The operating conditions of the reaction tubes in each tube bundle are the same; thus, a single tube is taken as an example in this paper. The model of the reactor is shown in Figure 3. Assuming that the reactor model is in a one-dimensional piston flow steady-state mode [6], there is no radial temperature or concentration gradient, and there is also no axial fluid mixing. There are three parts to EGR in the reactor, which are due to heat transfer, viscous flow and the chemical reaction.
The main reactions Ⅰ and Ⅱ that take place in the reactor: Assuming that the reactor model is in a one-dimensional piston flow steady-state mode [6], there is no radial temperature or concentration gradient, and there is also no axial fluid mixing. There are three parts to EGR in the reactor, which are due to heat transfer, viscous flow and the chemical reaction.
The main reactions I and II that take place in the reactor: As the reaction proceeds, the molar flow rates F k for each reaction component satisfy the following mass conservation equations: Entropy 2022, 24, 660 5 of 14 where A c is the cross-sectional area of the reactor in m 2 ; ρ b is the bulk density of the catalytic bed in kg · m 3 ; η j is the effective factor of internal diffusion (j =I, II); and r j is the reaction rate (j =I, II) in mol · kg −1 · s −1 . The reaction rate can be calculated according to the Arrhenius formulas, which can be written as: where p i is the partial pressure; k i is the reaction rate constant; A j refers to the preexponential factor in mol · s −1 · kg −1 · MPa −1 ; K i is the reaction equilibrium constant; and E j is the reaction activation energy in J · mol −1 . The values of A j and E j are given in Table 1, and fit well with the experimental data [6]. Table 1. Values of related parameters in the reaction rate equation [6].

Reaction 1 Reaction 2
The viscous flow process follows the momentum conservation equation [27], which can be expressed as: where ε is the porosity of the catalytic bed; µ m is the viscosity coefficient of the reactant in kg · s −1 · m −1 ; d p is the diameter of the catalyst particles in m; G is the mass flow rate in kg · s −1 · m −2 ; and v m is the average flow rate of the reactants in m/s. According to the theory of nonequilibrium thermodynamics, the local EGR of the reactor is [14,15]: The EGR of the reactor is the integral of the local EGR along the axial direction: where L CR is the length of reactor in m.

Optimization
The EGR of the EOCP is the sum of the EGR of each component: It is desirable for the EGR to be as small as possible while the C 10 H 20 yield is as large as possible; thus, the specific entropy generation rate can be established as follows [16]: S tot ∆F C 10 H 20 (23) In this study, the total molar flow rate at the inlet of the mixer (F T,in ) and the inlet pressure of the reactor (p in ) were chosen as optimization variables. Firstly, the influence of a single variable change on the EGR was analyzed when the other variables and operation parameters were the same as those of the reference process. Secondly, the process performance was optimized with the help of a genetic algorithm, the optimal variable parameter values for the smallest specific entropy generation rate were calculated, and the calculation results of the multi-objective optimization were compared with the reference process. The value ranges of the optimization variables are given by: The mathematical description of multi-objective optimization problem is: The parameters of each component of the reference chemical process are shown in Tables 2-4; the parameters listed in the table are taken from Ref. [6].

Minimization of SEGR
The influence of p in on ∆ . S tot and ∆F C 10 H 20 of EOCP is shown in Figure 4. As can be seen from the figure, with an increase in p in , ∆F C 10 H 20 gradually increased, the rate of increase slowed and the maximum value of ∆F C 10 H 20 was obtained near 4.1 MPa before gradually decreasing; ∆ . S tot increased monotonically. The compressor outlet temperature increased with the increase in p in , and the inlet temperature of the reactor increased under the action of the heat exchanger. The inlet temperature increased with the increase in p in , the reaction rates of exothermic reactions I and II were decreased, and the yield of C 10 H 20 was decreased. When p in was greater than 4.1 MPa, ∆F C 10 H 20 decreased, as the influence of temperature was greater than that of pressure. From 2 MPa to 4.1 MPa, although ∆F C 10 H 20 increased by 12.51%, ∆S tot also increased by 16

Minimization of SEGR
The influence of in p on tot S Δ  and 10 20 C H F Δ of EOCP is shown in Figure 4. As can be seen from the figure, with an increase in in p , 10 20 C H F Δ gradually increased, the rate of increase slowed and the maximum value of 10 20 C H F Δ was obtained near 4.1 MPa before gradually decreasing; tot S Δ  increased monotonically. The compressor outlet temperature increased with the increase in in p , and the inlet temperature of the reactor increased under the action of the heat exchanger. The inlet temperature increased with the increase in in p , the reaction rates of exothermic reactions Ⅰ and Ⅱ were decreased, and the yield of C10H20 was decreased. When in p was greater than 4.1 MPa, 10 20 C H F Δ decreased, as the influence of temperature was greater than that of pressure. From 2 MPa to 4.1 MPa, although 10 20 C H       The optimal variable values for the minimum Spec S Δ  are shown in Table 5. Compared with the reference process, the optimal in p increased by 0.9923 MPa and the optimal T,in   The optimal variable values for the minimum ∆ . S Spec are shown in Table 5. Compared with the reference process, the optimal p in increased by 0.9923 MPa and the optimal F T,in decreased by 0.9 mol/s. It was calculated that ∆

Multi-Objective Optimization
Using the PlatEMO toolbox, we set the population number as 100, the optimization objective number as 2, the optimization variable number as 2 and the genetic algebra as 100, and the Pareto optimal frontier of EOCP based on the objectives of ∆ . S tot minimum and ∆F C 10 H 20 maximum was obtained. LINMAP [28], TOPSIS [29] and Shannon entropy [30] decision-making methods were used, where the weight (w) of ∆ . S tot and ∆F C 10 H 20 was set to 0.5; the optimization results obtained using different decision-making methods were evaluated with the deviation index from Ref. [31].
For the LINMAP decision-making method, the optimal solution i opt is calculated as: where F j is the j-th target value, positive is the positive ideal point and PED i is the Euclidean distance between the i-th feasible solution and the positive ideal point. The LINMAP decision-making method calculates the closest point to the positive ideal point. For the TOPSIS decision-making method, the optimal solution i opt is calculated as: where negative is the negative ideal point and NED i is the Euclidean distance between the i-th feasible solution and the negative ideal point. The TOPSIS decision-making method determines the point farthest from the negative ideal point. For the Shannon entropy decision-making method, F i,j needs to be normalized first: The Shannon entropy of the j-th target can be calculated by the following formula: The optimal solution i opt is calculated as: where w j is the weight of the j-th target. The calculation formula for the deviation index is: S tot is at a minimum and ∆F C 10 H 20 is at a maximum are the least ideal on the Pareto front, as can be seen from Figure 8, which verifies that for the EOCP, the minimum ∆ The distribution of the Pareto front within the variation range of in p and T,in F is shown in Figures 9 and 10. Figure 9 shows that in the Pareto front, the optimization range of pin was 2-4.5 MPa; however, it was mainly distributed in the range of 2-3.5 MPa-that is, the selection of in p in this range would be beneficial to the optimization of EOCP. T,in F was evenly distributed in the whole optimization range, which shows that the selection of T,in F was used to adjust the opposition between decreasing tot S Δ  and increasing 10 20 C H F Δ at the same time. The distribution of the Pareto front within the variation range of p in and F T,in is shown in Figures 9 and 10. Figure 9 shows that in the Pareto front, the optimization range of p in was 2-4.5 MPa; however, it was mainly distributed in the range of 2-3.5 MPa-that is, the selection of p in in this range would be beneficial to the optimization of EOCP. F T,in was evenly distributed in the whole optimization range, which shows that the selection of F T,in  The distribution of the Pareto front within the variation range of in p and T,in F is shown in Figures 9 and 10. Figure 9 shows that in the Pareto front, the optimization range of pin was 2-4.5 MPa; however, it was mainly distributed in the range of 2-3.5 MPa-that is, the selection of in p in this range would be beneficial to the optimization of EOCP. T,in F was evenly distributed in the whole optimization range, which shows that the selection of T,in F was used to adjust the opposition between decreasing tot S Δ  and increasing 10 20 C H F Δ at the same time.   The results of the single-objective optimization, reference point and multi-objective optimization are shown in Table 6. The deviation indexes corresponding to the results of all the optimizations were smaller than the reference point. In the results of the single-objective optimization, the ∆ between the other two decision results. When using the deviation index to evaluate the results, the deviation index of the optimal solution under the Shannon enterprise decision result is the smallest; thus, it is the best optimization scheme. In the results obtained using the Shannon entropy decision method, the value of p in was 2.0315 MPa and the value of F T,in was 0.3909 mol/s. Compared with the reference process, the value of p in was reduced by 32.28% and the value of F T,in was decreased by 60.91%. The final value of ∆ . S tot was reduced by 59.01%, and the value of ∆F C 10 H 20 was decreased by 49.46%. The degree of this reaction was lowered by the Shannon entropy decision method, and the EGR was reduced by the decrease in the C 10 H 20 yield of the process.

Conclusions
An EOCP model with a constant temperature heat source outside the pipe was established in this paper. The process was optimized with the minimum SEGR. In addition, with the aim of making the EGR is as small as possible and the C 10 H 20 yield is as large as possible, the multi-objective optimization of EOCP was carried out using the NSGA-II algorithm, and the Pareto optimal frontier was obtained. The conclusions can be summarized as follows: 1.
The maximum C 10 H 20 yield and the minimum EGR cannot be guaranteed. When only p in was used as the optimization variable, the most economical point was obtained at 4.1 MPa. When only F T,in was used as the optimization variable, the most economical point was taken at 0.119mol/s.

2.
By taking the minimum SEGR as the optimization objective,∆ . S Spec was reduced by 27.7%, ∆S tot was reduced by 89.16% and ∆F C 10 H 20 was reduced by 85.01% compared to the reference process.

3.
A comparison of the decision points indicates that Shannon entropy is the best decision-making method. The corresponding ∆ . S tot of the decision point was reduced by 59.01%, ∆F C 10 H 20 was reduced by 49.46% and ∆S Spec was reduced by 18.88%.  Informed Consent Statement: Not applicable.

Data Availability Statement:
The data that support the findings of this study are available from the corresponding author, Y.Y., upon reasonable request.

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

S
Entropy generation rate, J · K −1 · s −1 Abbreviations isen Ideal reversible process M Mixer C Compressor HE Heat exchanger CR Reactor