Multi-objective simulation–optimization via kriging surrogate models applied to natural gas liquefaction process design

.


Introduction
Multi-objective optimization is a key to success in engineering design when the decision-making process involves multiple, competing objectives. However, to use optimization techniques in those tasks, models that describe the considered system are required. For complex chemical processes, modular-sequential process simulators (flowsheeting programs) are often resorted as simulated models. These programs provide accurate predictions for a variety of processes by using rigorous thermodynamic and unit operation models and numerical methods. The motivating problem for this work is the design of natural gas liquefaction processes considering the conflicting goals of minimizing both power consumption and heat transfer area utilization. The main challenge that makes the use of modular-sequential process simulators attractive in this application is modeling the multi-stream heat exchangers (MSHEs). That is especially the case because the streams in these heat exchangers undergo phase transitions, and their compositions and flow rates vary as they are design variables. There are attempts in the literature to model the MSHEs with algebraic formulation using Generalized Disjunctive Programming for dealing with phase change [1], superstructure-based approach [2], and nondifferentiable L.F. Santos et al. solutions of the following: The objectives ∶ R ↦ R and constraints ∶ R ↦ R are computationally-intensive, moderately noisy, black-box functions. and are the number of objectives and constraints, respectively, and is the dimension of . Since the noise is neglected for function evaluation, the expected value of these functions are approximated to the simulated value.
A suitable application for the Problem (MOSO) is the optimal design of LNG processes considering conflicting goals. The purpose of such processes is to cool the natural gas stream to −160 • C using refrigeration cycle(s) so it is liquefied at atmospheric pressure. These processes are mandatory for the natural gas safe and efficient commercialization and transportation through long distances, i.e. when the natural gas resources are far apart from the end consumer [6]. This energyintensive part of the LNG supply chain is responsible for 40%-60% of the cost of the final project [7]. For that, energy optimization of LNG processes have been extensively studied in the past [8]. However, the energy-optimum design solution without heat exchanger area limit might not be appropriate [9]. That is because the optimum design of these processes presents a clear trade-off between energy consumption and equipment size [10]. The concern of achieving high energy efficiency and diminished equipment size is further emphasized in offshore processes, where the plant site size is significantly restricted.
Although there is a vast literature on single-objective optimization of natural gas liquefaction processes, the evaluation of optimal trade-off in these processes' design has been timidly addressed. Khan et al. [10] applied multi-objective optimization to the dual-mixed refrigerant LNG process. The BOX optimization method was used to minimize separately the specific compression power and the heat exchanger area. The authors identified that these objectives are competing ones, and they used the non-dominated sorting genetic algorithm version two (NSGA-II [11]) implemented in the ''gamultiobj'' function from MAT-LAB to deal with balancing the objectives. Ghorbani et al. [12] used Genetic Algorithms (GAs) for single and multiple objectives to optimize the propane-precooled mixed-refrigerant (C3MR) LNG process. The competing goals of maximizing exergy efficiency and minimizing total product cost were considered. A trade-off non-dominated solution was determined with 7.25 $/h of product cost and 50.71% exergy efficiency. Song et al. [13] optimized the nitrogen expansion LNG process with carbon dioxide expansion precooling. The authors used a simple GA to individually minimize the specific energy consumption and maximize the liquefaction rate. The multi-objective scenario was considered using the NSGA-II and increased the liquefaction rate from 0.77 to 0.81, while diminishing the energy consumption by 10.1% compared to a base case.
Nguyen et al. [14] provided an interesting comparison of the singlemixed refrigerant (SMR), single-expander, and dual-expander LNG processes considering the trade-off between power consumption and heat exchanger area. A multi-objective GA was employed for the optimization, and different natural gas compositions were considered. Their findings indicate that the SMR is more energetically efficient at the price of having increased heat exchanger area utilization and process complexity. Primabudi et al. [15] optimized the C3MR process considering exergy efficiency and total product cost objectives. The NSGA-II solver was applied and the Pareto feasible solutions are between 0.557 and 0.613 for exergetic efficiency and between 45600 and 52776 $/h for the total cost of the product. Mofid et al. [16] optimized the parallel nitrogen expansion LNG process using a multi-objective particle swarm optimization (PSO) algorithm in two scenarios. At design, the objectives of minimizing energy consumption and heat exchangers size were considered. At operation, the minimization of energy consumption and maximization of liquefaction rate, with fixed equipment size are performed. The non-dominated solutions of the multi-objective PSO approach dominated the base case and GA results.
The literature on multi-objective optimization of natural gas liquefaction processes is entirely based on meta-heuristics, namely NSGA-II and multi-objective PSO. Recent findings have shown that, to this class of simulation-optimization problems, surrogate-based approaches coupled with gradient-based optimization can perform better than stateof-the-art meta-heuristics, such as GA and PSO [17]. Consequently, frameworks based on surrogate models applied to the multi-objective problem in the design of LNG processes can improve the assessment of their competing objectives. This work's objective is to apply a multiobjective optimization framework based on surrogates to the optimal design LNG processes using reliable process simulators and considering competing objectives.
The proposed methodology is an adaptation of our previous work [17] to the multi-objective problem. This framework includes using kriging models to replace the black-box objectives functions and constraints and -constraint method to handle competing objectives. The surrogate optimization subproblems are embedded into a nonlinear programming problem and solved with CONOPT local solver in the General Algebraic Modeling System (GAMS). This paper embodies as case study the optimal design of the SMR and C3MR processes for natural gas liquefaction considering the minimization of both power consumption and heat transfer area utilization. The Pareto Front determined by the present approach is compared with results from NSGA-II, a well-established and highly employed meta-heuristic method in multiobjective optimization of LNG processes. Suitable trade-off solutions of the SMR and C3MR processes are analyzed. These processes are compared regarding dominance in the objective functions space, and their suitability is discussed. Therefore, the novelty of the present work is two-fold. First, to the best of the knowledge of the authors, this is the first application of multi-objective surrogate-based optimization to design natural gas liquefaction processes. Second, an in-depth comparison between the optima results of the SMR and C3MR processes is performed. The comparison assessment does not consider only the single-objective energy optimum of each process as it is usually performed in the literature. Instead, the Pareto solutions of these processes with respect to the competing objectives of power consumption and heat exchanger area utilization are compared and a dominance inflection point is derived.

Kriging-based multi-objective optimization
An extension to the single-objective kriging-based optimization approach in Santos et al. [17] is proposed to solve the multi-objective simulation-optimization problems as in Problem (MOSO). The present framework relies on kriging surrogate modeling, -constraint method, and gradient-based solver.

Kriging model
The data-driven, interpolating model called kriging is defined as the sum of a regression model and an estimated error. The latter is considered to be the realization of a Gaussian process that predicts the deviation between the regression model and the underlying function. It is assumed to have zero mean and covariance given by some spatial correlation function. Kriging models are derived to be the best linear unbiased predictors at unsampled points [18]. In other words, the kriging model interpolates the data, and the departure between regression and actual function is more correlated for points that are closer together in space than those that are far apart.
(1) are available. For smoothness of the surrogate models, the Gaussian correlation function is used as Given the data and correlation model, it is possible to define the matrix ( ) ∈ R × ×( + ) of the correlation function evaluated at the sampled points for each objective and constraint, and the matrix ( , ) ∈ R ×( + ) of the correlation function evaluated at an unsampled point and the data for each objective and constraint. Both and are functions of the correlation parameters ∈ R ×( + ) , which are defined by minimizing the log-likelihood of the model given the data of each function [19].
In the present work, the kriging models assume a constant regression term (ordinary kriging), such that the error model captures the behavior of the underlying function that generates the data. Considering the above-mentioned assumptions and following the derivation of kriging as the best linear unbiased predictor as in Santos et al. [19] for and , these surrogates becomê in whicĥ= −1 ∕ −1 , = 1, … , + are the generalized least square solution of the ordinary kriging for each objective and constraint.

-constraint method
The -constraint method consists of selecting one objective to be optimized ( = ), while the others ( ≠ ) are converted into constraints such that [20] The multi-objective problem in Problem (MOSO) can be reformulated into an -single-objective simulation-optimization ( -SO) one as follows:

( -SO)
Varying the values of yields non-dominated solutions of Problem (MOSO), forming the Pareto Front [21].

Multi-objective surrogate optimization problem
The kriging models defined in Eq. (3) are readily available in mathematical notation and can be implemented in specialized software, algebraic modeling language systems, to be solved with efficient gradient-based optimization tools [17]. Algebraically, these kriging models in Eq. (3) arê in which ∶, = −1 ( −̂). Notice that, the memory-intensive matrix does not need to be stored. Instead, only , , and do. Given the algebraic formulation of kriging model in Eq. (5), it is possible to rewrite the problem in Problem ( -SO) into the -constrained surrogate optimization subproblem ( -SBO) explicitly such as min ∈,̃. .

( -SBO)
The multi-objective surrogate optimization in Problem ( -SBO) is an algebraic approximated substitute of the original black-box optimization in Problem ( -SO). Its solution is a possibly feasible, Pareto optimal candidate of the original problem considering the value. It is worth mentioning that in the nonlinear programming model in Problem ( -SBO), ,̂, and are constant. The resulting surrogate problem has + 1 variables, which arẽand , and + nonlinear constraints. Notice that the variablẽreceives the residual of̂( ) and is only necessary because some algebraic modeling systems, such as GAMS, solve for the minimization or maximization of a variable and not for the residual of an equation directly. Otherwise,̂( ) would be consider as objective function and Problem ( -SBO) resembling the definition of the Problem ( -SO). Therefore, the present krigingbased framework transforms the multi-objective black-box optimization in Problem (MOSO) into simpler -constrained nonlinear programming problems solved sequentially to provide samples that are possible non-dominated solutions to the original problem.

LNG processes
The optimal design of the SMR and C3MR processes, considering the trade-off between energy consumption and equipment size are assessed as case study in the present work. These refrigeration processes use a mixture undergoing a thermodynamic cycle to produce a heat sink in the evaporation stage to liquefy the natural gas stream. The C3MR process adds a refrigeration cycle operating with pure propane that evaporates at different pressure levels precooling the natural gas and mixed refrigerant streams. Fig. 1 illustrates the process flow diagrams of both the SMR and the C3MR.
The mixed-refrigerant cycle employs a compression system with four stages, considering coolers in between compressors and flash separators for recovering the liquid phase. The refrigerant condensate is sub-cooled in the MSHE-1 separately from the vapor phase, which is further cooled, liquefied, and possibly sub-cooled in the MSHE-2. These multi-component streams are mixed back inside the MSHE-1 in the cold pass, as represented in Fig. 1. A four-stage compression system is also considered in the propane cycle. The propane stream is expanded in valves at four pressure levels, and the resulting condensate is separated in the flash drums. The condensate streams are totally evaporated in different passes of the MSHE-100 to pre-cool the natural gas and mixed-refrigerant streams.

Process simulation
The rigorous process simulations are performed using Aspen HYSYS ® V9, employing as thermodynamic package Peng-Robinson that is suitable for nonpolar mixtures like hydrocarbons and nitrogen. It is supplied to the liquefaction process a natural gas stream (NG) at 5000 kPa and 32 • C. In Table 1, the NG composition jointly with some process considerations and other simulation parameters following the literature [22] are presented. The NG mass flow rate is 1 kg/h as a basis of calculation. However, it does not loose of generality, once the results can be reproduced for different natural gas flow rates. That can be achieved by multiplying the extensive variables of refrigerant component mass flow rate, specific power consumption, and specific heat conductance (UA), to be presented later, by the factor of increase/decrease. The mixed-refrigerant refrigerant is composed of nitrogen, methane, ethane, propane, and i-pentane. Each component mass flow rate is a decision variable of the optimization problem. The heavy component i-pentane is considered only for the SMR to account for the cooling task at high temperatures as it increases the bubble point temperature of the mixture [23]. This is not required in the C3MR because that is achieved by the propane precooling cycle. The pressure of discharge and suction of the mixed-refrigerant cycle, and , are also decision variables, while the intermediate pressures come out of the constant pressure ratio assumption given by ( ∕ ) 1∕4 . The temperatures of hot streams leaving an MSHE are equal and set to the expansion temperature (stream 11l in SMR and 9l in C3MR) for MSHE-L.F. Santos et al. 1, which is a design variable. To ensure the LNG requirements at the final pressure as exposed in Table 1, the temperature of −149.7 • C is fixed in the final expansion stage.
The propane cycle in the C3MR has some degrees of freedom that can be satisfied with the following assumptions. Setting discharge pressure (stream 20') to the minimum, while above saturation pressure of propane at 40 • C at stream 1'. This pressure is set to 1398 kPa. The propane streams leave MSHE-100 (3l', 6l', 9l', and 12l') as a saturated vapor at the same temperature of the inlet streams, such that only the heat of vaporization is used as heat sink. The propane streams enter MSHE-100 (2l', 5l', 8l', and 11l') 3 • C below the temperature of the cold-end hot streams (14,15,16,17) to guarantee the minimum temperature driving force. The latter is calculated to achieve a constant temperature step at each pass of MSHE-100 that is 14 − , ∕4. The final temperature of precooling cycle (stream 1 and NG5 in Fig. 1 (ii)) is a decision variable of the process optimization. Given these assumptions and decision variables, the precooling cycle is fully determined and the mass flow rate of propane, as well as the intermediate pressures of the cycle, come from the energy balances performed by the process simulator.

Process optimization
The goal of design in natural gas liquefaction processes is to reduce energy consumption while maintaining small pieces of equipment. Also, a minimum value of temperature difference between hot and cold streams all through the MSHE 1 and 2 has to be assured, once in MSHE-100 it is already specified as elaborated in Section 3.1. The kriging model, which is to be used to replace objectives and constraints, requires some degree of smoothness of the black-box functions, and that may not be the case for this constraint at every configuration of decision variables. In fact, the nonlinear behavior of enthalpy with temperature in mixtures, which is worsened by phase change, can cause a change in the pinch point location of the multi-stream heat exchanger. The change in the pinch location from a stream regime to another, e.g. from vapor to vapor and liquid or liquid phase, can cause an abrupt Table 1 Simulation parameters and process considerations for the LNG processes, based on the available literature [22].  Fig. 2 with a sensitivity analysis around a carefully selected configuration of decision variables, where the non-smooth behavior of this constraint happens.
To overcome this issue, Santos et al. [17] propose to discretize the MSHEs 1 and 2 in = 1, … , × sections to make it more likely L.F. Santos et al. that the temperature driving force constraints are smooth functions. is the number of sections in each of the multi-stream heat exchangers. It is also illustrated in Fig. 2 with = 10, for which it is possible to see that the × functions (solid lines) are less stiff than the one without discretization (dashed line). The optimal process design can be formulated as the following multi-objective optimization problem . .
In Problem (LNG), ( ) is the power consumption of unit that is contained in compressors and pumps set, .
( ) is given by the enthalpy change caused by the compression as if it was isentropic divided by the isentropic efficiency, , where is the temperature that provides equal entropy in the inlet and outlet streams (  For the SMR process, ={MSHE-1, MSHE-2, C-1, C-2, C-3, C-4}, and for the C3MR, ={MSHE-1, MSHE-2, MSHE-100, C-1, C-2, C-3, C-4}. The cold utility used in the coolers {C-1, . . . , C-4} is considered to be at 30 • C. Notice that these calculations, jointly with the equation of state and flash calculations, are performed in the process simulator and these results are readily available to the user. The work consumption and UA values are divided by the mass flow rate of LNG,̇, to generalize the process energy consumption and heat conductance in terms of LNG production. Therefore, the total power consumption and UA are readily available for a different natural gas flow rate by multiplying the resulting LNG production by the respective specific values. ( ) and ( ) are the temperatures of hot and cold composite curves in the th segment of the MSHE-1 and MSHE-2, and is the set of the points from composite curves calculations that belong to section [17]. These values are also readily available to the user in the ''LNG Exchanger'' unit operation of the Aspen HYSYS process simulator, which, to account for the nonlinear behavior of enthalpy with temperature, performs a temperature discretization in the hot and cold streams and assures energy balance in each temperature stage. Notice that the temperature driving force constraints are not required because it is already assured by setting the temperature of the cold streams 3 • C below the cold-end of the hot streams, as explained in Section 3.1.
As discussed in Section 3.1, there are eight decision variables in for each process so = 8.  Table 2 for SMR and C3MR. These values are set to be between 50% below and above a heuristically determined base case, i.e. [ . These bounds are intended to include promising regions while avoiding unstable simulations convergence because of flowsheet configuration that is far beyond the feasible space. The lower bound of is set to 0 to account for refrigerant mixtures without nitrogen. The final precooling and expansion temperatures are bounded heuristically.

Optimization framework
The present multi-objective simulation-optimization framework is presented in Algorithm 1. Also, a schematic illustration is introduced in Fig. 3 to improve the clarity of the presented algorithm as well as the connections between the pieces of software. This framework substitutes the black-box objectives and constraints functions with kriging models and employs the -constraint method to handle the competing objectives. The algebraic approximated optimization subproblems are embedded into a nonlinear programming problem and solved with CONOPT local solver in GAMS. The calculation of and at given is achieved by a function written in MATLAB that connects to the process simulator Aspen HYSYS. This connection is via the MATLAB built-in function ''actxserver'' that exposes the Aspen HYSYS objects and functions to the programming environment. The values of are set to the respective object property (stream or unit operation), and the and values are acquired after converging the simulation.
Step 1 of Algorithm 1 is to generate 0 samples of . For that, the maximization of the sample points' minimum distance in the search space  is performed, employing a Latin Hypercube algorithm. With programming and simulation environments connected, the value of and are calculated in the simulation for each ∈ . The initial data is defined as 0 = [ ].
L.F. Santos et al. Given 0 , , and , step 2 of Algorithm 1 is the single-objective optimization for each objective, solved using the framework for constrained, black-box optimization proposed in Santos et al. [17]. It consists of using the data in 0 to fit kriging models for the objective and constraints kriging-based adaptive sampling procedure ( or regular ones). These models are implemented in a nonlinear programming problem in GAMS and solved with CONOPT gradient-based solver. The solution, which is a promising feasible and decreasing candidate, is evaluated in the simulation and the , , and values are attached to the dataset. The surrogate models are upgraded with the new data and this process is repeated until the surrogate problem optimization fails to provide an improving candidates five times or simulation evaluation budget is reached. Talk about stopping criteria. The reader is invited to read Santos et al. [17] for more information on this algorithm and implementation.
The third step is to use the value of 1 at the solutions of each single-objective problem, 1 ( * 1 ) and 1 ( * 2 ), to bound the vector. The -constraint vector is defined as the equally distributed -vector from 1 ( * 1 ) to 1 ( * 2 ). Fourth step consists of appending the single-objective solutions to the data 2 . The subscript of stands for the number of non-dominated solutions it already stores. In this case, it is two that account for the single-objective solutions. For the other entries of , the -constrained single-objective optimization in Problem ( -SO) is solved using a slight modification of the aforementioned framework for constrained, black-box optimization proposed in Santos et al. [17]. The difference now is that the surrogate optimization problem includes the kriging-based -constrain for the 1 . Each optimization solution is appended to −1 , so that = −1 ∪ [ * , ( * ), ( * )]. The last step of the proposed framework is to eliminate dominated and infeasible solutions from , which contains the non-optimal initial data ( 0 ) and possibly non-improving iterations from the algorithm. The new set, , is the Pareto Front as it only contains non-dominated solutions. In other words, the solutions that are infeasible ( ) > 0 are eliminated. Also, every solution that is dominated by another is eliminated, i.e. if there exists at least one ∈ such that ( ) < ( ) for all , then is dominated and eliminated. Notice that the data in and, consequently, in the non-dominated solutions are given by the rigorous simulation and not by the surrogates models. Therefore, the optimal natural gas liquefaction processes designed with Algorithm 1 would conform with rigorous thermodynamic and unit operation models, instead of surrogate approximations.

Results
The kriging-based multi-objective optimization approach is applied to the optimum design of SMR and C3MR processes for LNG production, considering = 2 objectives of minimizing power consumption and heat exchanger area utilization. The optimization parameters are initial sample size 0 = 10 , function evaluation budget = 20 , number of decision variables = 8, and number of sections into which each MSHE is divided = 10. The number of non-dominated solutions is heuristically determined as the number of objectives times the size of the problem plus the number of objectives for the single optima: = × + = 18. The NSGA-II parameters are set to the default of the algorithm from the MATLAB optimization toolbox (''gamultobj''), except for individuals (10 ) and generations (20) that are selected heuristically and to match the function evaluation budget of 20 for a fair comparison between algorithms. To account for the randomness of the initial sample in both approaches, each algorithm was run three times. Fig. 4 illustrates the Pareto Front delimited by the non-dominated solutions achieved with the present framework as well as the NSGA-II algorithm. For both SMR and C3MR, the present Pareto Fronts dominate the meta-heuristics results. The optimal SMR processes have non-dominated solutions that present energy consumption as low as 1045 kJ/kg-LNG and specific UA value of 212.2 kJ/( • C kg-LNG). Meanwhile, the extreme solutions for NSGA-II are 1222 kJ/kg-LNG of energy consumption and specific UA value of 216.7 kJ/( • C kg-LNG). The nondominated solutions of the C3MR process have energy consumption and specific UA as low as 980.3 kJ/kg-LNG and 266.9 kJ/( • C kg-LNG), while the NSGA-II have 1054 kJ/kg-LNG and 284.7 kJ/( • C kg-LNG). These results show that the present framework was also able to provide broader Pareto Fronts in the single-objectives directions. It can be inferred that the ranking method of dealing with competing objectives in population-based approaches was unable to spread the L.F. Santos et al. non-dominated solutions enough to reach close to the single-optimum solutions. Therefore, bounded objective function methods such as theconstraint with an appropriate single-objective optimization framework of the transformed multi-objective problem have shown to be more efficient. The limitations of the present approach are for many objectives, as the size of the vector grows quickly with the number of objectives.
The shape of these Pareto curves shows that the rate of change in the heat exchanger area with specific power consumption is more pronounced at low values of the latter. That is expected from the physical point-of-view because the heat exchanger area is roughly directly proportional to the inverse o the temperature differences between hot and cold composite curves. And, to achieve low power consumption, these temperature differences have to be very small, which makes its inverse (roughly area) explodes to high values. Based on that, an interesting trade-off solution would be not so small in power consumption because the UA values increase too quickly (high rate of change), but not so small on UA as it requires to increase power consumption too much (low rate of change). Trade-off solutions can be derived numerically from the Pareto curve in various ways [24]. In the case of heat exchanger area and power consumption objectives, a decisionmaker could choose from the non-dominated solutions considering weight factors that take into account their costs and environmental impacts. However, these values are volatile and uncertain, requiring a more in-depth study for a single application. On the other hand, a simple, more robust, and yet good choice of non-dominated solution that balances high energy efficiency and heat exchanger area utilization is the one given by the reference point criterion [25]. This consists of choosing the non-dominated solution with smaller euclidean norm to a reference point in the objective function space. The reference point is here considered to be the minimum value of and power consumption ([min ∈ ( ), min ∈ ( )]), such that where is the trade-off solution (see Eq. (6)), and is the set of non-dominated solutions. For both processes, these solutions are the fifth non-dominated one, which are highlighted in Fig. 4 with the magenta star. For the SMR process, the trade-off solution presents a 24.40% increase in work consumption and a 43.83% increase in UA compared to the energy and UA optima. For the C3MR process, it presents a 16.70% increase in work consumption and a 14.54% increase in UA compared to the energy and UA optima. For a further investigation of the process design point of view of these optimization results, Tables 3 and 4 present the decision variables, objectives, and constraints of the energyoptimum, UA-optimum, and the trade-off solution of SMR and C3MR, respectively.
From the multi-objective optimization results of the SMR process in Table 3, it is evident that the energy-optimum solution presents a greater mass flow rate of the heaviest components to achieve phase L.F. Santos et al. separation so that the liquid phase is the most responsible for the cooling task down to of −27.43 • C in MSHE-1. This helps save power consumption in the compressors. Also, the overall pressure ratio of 10.53 is the smallest among the three non-dominated solutions analyzed in this table. The UA-optimum solution presents an LNG process configuration that exploits high pressure ratio (60) of the refrigeration cycle to increase the temperature difference between hot and cold composite curves in the MSHEs to diminish their heat transfer area. The composition and flow rate of the multi-component refrigerant are also determined by the optimization approach to maintaining the high-temperature driving force throughout the MSHEs, as the power consumption is not been taken care of. The trade-off solution presents an intermediate cycle pressure ratio (42.22) compared to the other two solutions. The amount of heavy component (iso-pentane) in the mixed refrigerant is also in-between values of single-objective optima. The same behavior can be observed in the multi-objective optimization results of the C3MR process in Table 4. The overall cycle pressure ratio goes from 12.63 to 32.13, and finally to 59.57 for energy-optimum, trade-off solution, and UA-optimum, respectively. The detailed results of specific ∈ in Tables 3 and 4 show that the power consumption in the mixed-refrigerant cycle follows the net power consumption trend. On the other hand, the compressors of the propane cycle consume the most energy for the energy-optimum solution, which is counterintuitive. The mixed-refrigerant for this solution is richer in heavier components (ethane and propane) so it requires more propane heat sink and, consequently, more compression power. That pays off with smaller overall pressure ratio as explained above. The detailed results of specific UA ∈ in Tables 3 and 4 present clearly that the main contributors for diminishing or augmenting the overall heat transfer area utilization are the MSHE-1 and MSHE-2. Therefore, the overall UA is mostly dependent on the mixed-refrigerant composition and flow rate, even for the C3MR process.  energy and UA optima. In other words, the non-dominated solutions are not fine-tuned interpolations of single-objective ones. This interesting finding shows the capability of a multi-objective optimization framework to determine complex and unique non-dominated solutions forming a Pareto curve for a decision-maker to choose from. Therefore, the importance of multi-objective optimization approaches to design problems with competing objectives is determining the rate of change of one objective with the other and deriving trade-off solutions that are not mere adjusted interpolations of single-objective ones. The same quality of trade-off analysis would not be possible by simply using single-objective optimization of a function that weights the competing objectives.
Another important use of the multi-objective optimization approach is to serve as a tool for helping to decide among technologies. Fig. 5 illustrates a comparison between the SMR and C3MR liquefaction processes based on specific power consumption and heat exchange area utilization. Given economic and industrial site scenarios, it would be readily possible to determine which is the natural gas liquefaction technology most suited considering the price of electricity, installed power, heat exchanger costs, and available plant site area. For the given natural gas conditions and process considerations, the C3MR solutions dominate the SMR ones up to around 1355 kJ/kg-LNG of specific power consumption and 274.5 kJ/( • C kg-LNG) of specific UA. It means that, if energy consumption is required to be below this value, the C3MR natural gas liquefaction process should be employed. However, if the overall UA value is expected to be smaller than 274.5 kJ/( • C kg-LNG), then the SMR process should be selected. This inflection point that separates the dominance between C3MR and SMR in objectives functions space is drawn in Fig. 5 with the blue point. All-in-all, C3MR processes are more energetically efficient for flexible heat exchanger area utilization, whereas SMR processes are more energetically efficient for strict heat exchanger area utilization. SMR processes can achieve small values of UA that C3MR cannot, while the former can achieve small values of power consumption that the latter cannot.
Comparing the optimization results with the literature is challenging because the different author uses different process considerations such as natural gas pressure, temperature, and composition, compressors and pumps efficiency, pressure drop in heat exchangers, LNG specifications, cooling temperature, and process constraints. Table 5 presents the literature results of natural gas liquefaction processes optimized to minimize the specific power consumption that considered the same parameters and considerations presented in Table 1. These results show the competitiveness of the kriging-based that surpassed the recent literature, except for the SMR proposed by Qyyum et al. [26]. They included hydraulic turbines in the refrigeration cycle instead of Joule-Thomson valves, which increases the energy efficiency but is not conclusive L.F. Santos et al.

Table 5
Comparison with literature results on energy-optimum natural gas liquefaction processes.
Khan et al. [27] Khan et al. [22] Qyyum et al. [26] Ali et al. [28] Ali et al. [29] Majeed et al. [  without an economic assessment. Therefore, these results provide a flavor of the quality of the processes designed in this paper using the kriging-based multi-objective approach.

Conclusions
This paper presented an application of a kriging-based optimization framework to solve constrained, black-box, multi-objective, simulationoptimization problem of the optimal design of SMR and C3MR natural gas liquefaction processes, considering the minimization of both power consumption and heat exchanger area utilization and using a reliable process simulator. The framework is based on fitting kriging models to simulation data to substitute the black-box objectives and constraints with a simple algebraic formulation. The -constraint methodology is employed to reformulate the multi-objective problem into several single-objective ones. The surrogate, single-objective optimization subproblems are solved in GAMS with CONOPT solver to provide a sampling tool of promising feasible and non-dominated solutions to the original problem.
The Pareto Fronts determined by the present approach to the SMR and C3MR processes dominate the ones achieved with well-established multi-objective meta-heuristics of NSGA-II and are also broader in the single-objective directions. The objective functions of non-dominated solutions go as low as 1045 kJ/kg-LNG and specific UA value of 212.2 kJ/( • C kg-LNG) for SMR and 980.3 kJ/kg-LNG and 266.9 kJ/( • C kg-LNG) for C3MR. The trade-off solutions that present the minimum sum of relative objectives were analyzed. These trade-off results balance the objectives by a 24.40% increase in work consumption and a 43.83% increase in UA for the SMR process and a 16.70% increase in work consumption and a 14.54% increase in UA for the C3MR process, compared to the energy and UA optima. The comparison of Pareto solutions of the SMR and C3MR processes showed that the C3MR solutions dominate the SMR ones up to around 1355 kJ/kg-LNG of specific power consumption and 274.5 kJ/( • C kg-LNG) of specific UA. This finding provides a numerical metric for choosing between the C3MR, which is suitable for energetically efficient applications, and the SMR, which is more appropriate for restricted heat exchanger area utilization. A comparison of optimal LNG processes with previous works showed the competitiveness of the kriging-based that surpassed the recent literature.
The importance of multi-objective optimization approaches to design problems with competing objectives is determining the rate of change of one objective with another, deriving complex trade-off solutions that are not mere adjusted interpolations of single-objective ones, and providing a more thorough and robust platform to compare among choices with competing objectives. For the case of natural gas liquefaction process design, it allowed deriving promising SMR and C3MR solutions that balance power consumption and heat transfer area utilization. Also, an inflection point in the objective space is observed. This point separates dominance between C3MR and SMR processes, making it easier for a decision-maker to choose the appropriate technology for specific scenarios.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Data availability
Data will be made available on request.