Andreas Psaltis

Delphine Sinoquet

Alexandre Pagot
email: alexandre.pagot@ifpen.fr

Systematic Optimization Methodology for Simultaneous Process Design and Heat Exchanger Network Synthesis with Surrogate Models -Application to Distillation Column Complex

Distillation units require huge amounts of energy for the separation of the multicomponent mixtures involved in refineries and petrochemical industries. Keeping the same separation standards, energy savings can be achieved in multiple ways. The overall efficiency of the distillation column system is determined from the trade-offs of the Operating Expenditures (OPEX) and Capital investment cost (CAPEX), as there is a strong interaction between the distillation columns and the Heat Exchanger Network (HEN) of the interconnecting streams. Therefore, the highly complex optimization of this trade-off necessitates moving towards process integration solutions. Individually, both distillation columns and HEN are well-defined rigorous problems and numerous methodologies have been developed for finding optimal solutions. The integrative problem however is larger in scale, computational demanding and less robust. In this paper, a systematic optimization methodology for process integration of a multicomponent distillation column complex is presented. The highly nonlinear and time consuming rigorous models of the distillation column are being substituted with simple surrogate models that generate operating responses with adequate accuracy. Moreover, another surrogate model is included in the Mixed Integer Non-Linear Programming (MINLP) formulation to enhance the accuracy of the results by considering the phase changes of the process streams involved in the HEN. The methodology is applied on two case studies of the aromatics separation PARAMAX complex and the results illustrate significant reductions on the Total Annualized Cost. With a scope limited to the benzene and toluene columns, the gain reaches about 15%.

Introduction

The continuous rising prices in energy and stringent environmental regulations have lead Chemical industries to invest towards energy efficient solutions. Particularly, for Separation industries such as Refineries and Petrochemical plants, which are among the most energy intensive, a major concern is to minimize the huge operational demands of the distillation columns which are by far the most preferred unit for separation. The low thermodynamic efficiency of the columns, stemming from the higher temperature profile of the reboiler compared to the condenser, has lead researches on improving the performance of HEN between the streams involved in the distillation column complex, by targeting the minimization of the external heat sources and the associated HEN equipment costs.

Over the last decades, several methodologies have been emerged for solving the HEN synthesis problem and can be classified as sequential or simultaneous approaches. The first approach decomposes the problem in a sequence of objectives. Initially the target is to identify the HEN with the minimum utilities demand. Once the energy requirements are met, the second objective is the reduction of capital cost by minimizing the number of Heat Exchangers. Finally, detailed calculations are performed in order to investigate further reductions on the Heat Exchanger Area cost. The Pinch analysis method, developed by Linnhoff and co-workers [START_REF] Linnhoff | Thermodynamic analysis in the design of process networks[END_REF][START_REF] Linnhoff | The pinch design method for heat exchanger networks[END_REF], is among the most influential methodologies of this approach. The energy profile of the process streams is graphically represented on composite curves and according to the heuristic-based value of the minimum allowable temperature difference (ΔT min ) the pinch point determines the energy requirements from external sources.

The decomposition of objectives is also being addressed in mathematical programming formulations such as the transshipment methodology of [START_REF] Papoulias | A structural optimization approach in process synthesis: Part II -Heat recovery networks[END_REF] and [START_REF] Floudas | Automatic Synthesis of Optimum Heat Exchanger Network Configurations[END_REF]. On the other hand, the simultaneous approach methodologies are characterized by the complex mathematical formulations employed for the optimization of the combined sub-objectives of CAPEX and OPEX within well-defined process restrictions. Moreover, potential heat exchange between process streams is examined in superstructures (Yee and Grossman, 1990;[START_REF] Ciric | Heat Exchanger Network Synthesis without Decomposition[END_REF] with the introduction of binary variables, which are used for the realization of existence or not of a heat exchanger. A review of the main literature for HEN synthesis is presented in [START_REF] Furman | A critical review and annotated bibliography for heat exchanger network synthesis in the 20th century[END_REF] and more recent approaches can be found in Klemes and Kravanja (2013).

Heat integration methodologies demonstrate remarkable energy savings but are limited to the existing process conditions which are defined during Process design, performed in advance. Bounded by these conditions, the results of the sequential procedure are suboptimal compared to Process integration configuration, namely when process design and Heat integration are optimized simultaneously. Process integration takes into account the strong interaction between the process units and the HEN of their interconnection streams and the best scheme corresponds to the optimum trade-off compromise of the conflicting CAPEX and OPEX criteria. Despite the numerous Heat integration methodologies reported in literature, limited publications addressed the holistic optimization problem, mainly because the increased complexity of the problem provides little margins for heuristics and the mathematical formulations require significant computational efforts due to the very large search space and are harder to implement and solve. [START_REF] Duran | Simultaneous Optimization and Heat Integration of Chemical Processes[END_REF] pointed out the difficulties of implementing rigorous models due to the combinatorial nature of the problem and reduced the degrees of freedom to energy optimization. The Heat integration part in their Non-Linear Programming (NLP) formulation is restricted on finding the minimum utilities usage by imposing a set of constraints derived from pinch point analysis. Their work was extended by [START_REF] Lang | Simultaneous optimization and heat integration with process simulators[END_REF] with sequential modular simulators and by [START_REF] Grossmann | A rigorous disjunctive optimization model for simultaneous flowsheet optimization and heat integration[END_REF] where the developed rigorous disjunctive formulation considers isothermal stream cases as well. In contrast, [START_REF] Papalexandri | A decomposition-based approach for process optimization and simultaneous heat integration: Application to an industrial process[END_REF] developed a methodology that takes into account the Total Annualized Cost (TAC). The overall process is decomposed in an inner loop where the Heat integration part corresponding to a fixed process flowsheet is solved and then the information is passed to an outer loop where the associated process flowsheet is optimized. [START_REF] Zhang | A Simultaneous Optimization Strategy for Overall Integration in Refinery Planning Ind[END_REF], following the common practice of employing linear and linearized models in Refinery optimization, took into account in their optimization strategy the interactions of the hydrogen network, the utility and material processing system for maximizing the margin between the total revenue and the operating cost. The simplification of the original problem into a MILP (Mixed Integer Linear Programming) formulation was necessary, according to the authors, because no robust MINLP solver is commercially available for handling refinery scale applications. [START_REF] Ponce-Ortega | Optimal synthesis of heat exchanger networks involving isothermal process streams[END_REF] extended the superstructure of [START_REF] Yee | Simultaneous Optimization Models for Heat Integration -II Heat Exchanger Network Synthesis[END_REF] for HEN retrofit to include operational and structural process modifications and to treat isothermal streams. The modifications can be either operational (e.g. pressure changes) or structural (e.g. addition or replacement of equipment) and both capital and operating costs associated with the retrofit decisions are considered simultaneously by the MINLP methodology.

In order to overcome the numerical complexity involved in the above mentioned deterministic methodologies, probabilistic approaches have been proposed [START_REF] Lazzaretto | Energy, economy and environment as objectives in multicriterion optimization of thermal systems design[END_REF]Gassner and Marechal 2009;Fazlollani and Maréchal 2013) and are characterized by : o Decomposition of the overall optimization procedure into sub-problems which are solved in different interconnected softwares, o Use of commercial simulators as black box rigorous models o Probabilistic algorithms (e.g. evolutionary algorithms) for multi-objective optimization procedure.

The first sub-problem considers technological alternatives, modeled in commercial simulators as a superstructure, for the calculation of the energy and material flows corresponding to a predetermined set of operating conditions. The process output is then processed in an energy and process integration model in different software. The opportunities of heat integration are assessed and the system interactions are optimized with regards to the minimum operating cost. At the same time, binary variables are introduced as utilization factors of each technology defined in the superstructure and the interactions are being evaluated by means of integrated composite curves. The data returned by the technology selection and energy integration steps are further processed in a post calculation Multi-Objective Optimization procedure (MOO). The MOO is based on an evolutionary genetic algorithm that generates a population of points in each iteration. The population represents potential solutions and each solution represents a different trade-off between the optimization objectives. The best points in this population are selected as they approach an optimal solution and displayed on a Pareto front, on which the optimal values of the competing objectives and resulting trade-offs are identified. [START_REF] Ochoa-Estopier | Optimization of Heat-Integrated Crude Oil Distillation Systems. Part II: Heat Exchanger Network Retrofit Model[END_REF] proposed a simultaneous optimization and heat integration framework where the sub-problems are solved in different softwares. Process design is carried out in commercial simulators which serve as black box rigorous models for thermodynamic calculations (SRK), such as Aspen HYSYS or ProII; Heat integration is modeled as LP problem in GAMS and the cost evaluation is performed in Excel and Python. All these softwares are also linked to a derivative-free optimizer for the determination of the decision variables, the values of which define the optimum economic objective. Ochoa-Estopier and Jobson (2015a, 2015b, 2015c) presented a methodology for operational optimization of crude oil distillation units with consideration of retrofit modifications on the associated HEN. The commercial simulators are also used for sampling data generation in order to regress the parameters of the distillation columns which are represented with artificial neural networks. Simulation annealing is used to optimize the operating conditions of the columns and propose structural HEN modifications. In a second level, a nonlinear least square problem is used to enforce the HEN constraints.

Simultaneous employment of different softwares provides an alternative perspective to the conceptual design but at the same time significantly increases the computational load involved in the multiple tasks as a result of the very large search space of operating conditions, which is further increased by the additional time needed for the interconnecting softwares to exchange information. Moreover, a limited number of process engineers, with adequate optimization skills and programming knowledge for all individual softwares involved, can utilize these methodologies for their applications.

The purpose of this paper is to present a novel methodology for the optimization of the Process integration synthesis problem. In the proposed methodology, process unit design, operating conditions and the HEN synthesis of the streams involved are simultaneously considered in an MINLP optimization formulation and the TAC defined from the trade-off between the CAPEX and OPEX criteria, is used as a quantitative index for determining the optimum configuration. A key element of the solution strategy is the replacement of the complex rigorous models needed for process design with simple surrogate models that serve as computationally low-cost substitutes with adequate accuracy. Moreover, phase changes during the HEN synthesis are considered using a modified stagewise superstructure of [START_REF] Yee | Simultaneous Optimization Models for Heat Integration -II Heat Exchanger Network Synthesis[END_REF].

In the following sections, a detailed description of the methodology is presented for energy intensive distillation column applications, which despite their maturity as separation process still offer significant margins for improvements and thus, making them ideal for illustrating the benefits of the proposed solution strategy without loss of generality. The main idea of the solution strategy is to benefit from the different capabilities of three softwares and to carry out sequentially 5 separate tasks as shown in Figure 1. Rigorous models provide all the information needed for an accurate representation of the process of interest but are difficult to implement and solve, especially for large scale systems where the accuracy is often traded with computational efficiency. For this reason, the first step of process representation is performed in commercial simulators for given technology and feasible set of design and operating variables. It should be noted that at this step the existence or not of a heat exchanger is not known and a stream to be heated or cooled is identified with flash units. The simulators are too time-consuming to be used directly in the optimization procedure because of the very large search space and the plant wide nature of the process. Hence, the overall process is partitioned to the consisting elementary units (columns, Batch calculations for the characterization of each elementary unit ) and for each one a series of batch calculations is performed in this stage for various values of the design variables (degrees of freedom) in order to generate a representing set of sampling data needed for the identification of the most influencing variables used to reproduce the response of the process variables associated with the optimization objectives, through surrogate models. The batch calculations have been generated by regular grid with defined sample steps for each variable.

Outline of the solution strategy

Based on the sampling data, the parameterization of the surrogate models (parametric models either linear or nonlinear depending on the modeled function) is performed in the third step, using the powerful regression routines available in MATLAB. The tuned parameters are passed to the MINLP formulation of the final step, where the simultaneous optimization is carried out. In addition, as the solution of NLP systems is highly dependent on the initial conditions, and global optimality cannot be guaranteed, multistart approach for finding the best local solution is also performed in MATLAB. For each set of predetermined initial values on MATLAB, the optimizer returns a local solution and the best one is selected among the set of alternatives stored in MATLAB. Distillation columns are frequently encountered in many chemical processes and they have been extensively studied in the previous decades in order to increase the understanding of the thermodynamic phenomena taking place and to develop models capable to represent their behavior. These models represent the relationship between the basic variables that correspond to the degrees of freedom and the intermediate variables used to describe the behavior of the process and at the same time define the trade-off between the optimization objectives. In this paper, the degrees of freedom of a distillation column are the number of trays in the column (NT su ), the pressure at the top of the column (P su Up ) and the n F column feed stream enthalpies (H n F ,su Feed ). Either directly or via the surrogate models, the value of these variables determines both operating and equipment design criteria and thus the optimal process integration is inferential to the optimal value of these variables. The optimal value of the design variables corresponds to the minimum value of the composite index of TAC that takes into account the trade-off between the OPEX and the annualized CAPEX with an expected payback time of six years.

Methodology for distillation column applications

Process design

Initially, the process synthesis and design are performed in commercial simulators, such as PRO/II and HYSYS. This step focuses on setting the links between the primary units (e.g. drums, columns, reactors) and the interconnecting process streams with the objective of establishing the material balance in respect to the final product specifications and the individual process unit constraints such as maximum solvent recovery and operational requirements of the catalysts. As the material balance shows the amount of resources used for achieving the desirable product rates and purities, generally in separation processes the objective is to investigate its implementation through the selection of the suitable technology, from a set of limited alternatives. Moreover, setting the material balance early on the solution strategy reduces the degrees of freedom and the search space of the optimization procedure and therefore, simplifies the overall process integration synthesis objective without effecting the final solution. For example, if the degrees of freedom are reduced to the minimum, a distillation column can achieve the desired separation purities for given flowrates with a variety of combinations between the number of trays and reboiler's duty.

Partition to elementary units

Once the material balance is fixed, the process flowsheet is partitioned and a series of batch calculations is performed for the characterization of elementary units. In our case, these elementary units are:

-Distillation columns (as shown in Figure 2) -Flash drums which are used o for the calculation of the temperature or/and enthalpy of the interconnecting streams with phase change considerations o the operational bound values of the exchanger Rigorous stage-by-stage (meaning with process simulation software) calculations are not needed; the energy requirement of the column is determined by the reboiler and condenser duties which in turn has significant impact on the HEN synthesis, as a consequence of heat integration. In Table 1 are listed the operating variables for the simulation of the distillation column. The diameter of the column is not an operating variable (in fact it is a key variable for the shell design) but as illustrated in the next section, it is calculated as a function of the reboiler's duty and thus, by convention is considered as such.

Surrogate models

Depending on the purpose and the specific application of the examined study, simulation models can be developed utilizing three distinct approaches:

o Rigorous models: deterministic models used for extensive tray calculations. They illustrate very good performance in terms of accuracy but the rigorous computations are time consuming and thus rigorous models are not preferred for large scale optimization purposes. o Simplified physical approximation models: short-cut models developed from both physical laws and empirical relations. They are more computationally efficient than the rigorous models but are based on some assumptions that induce loss in accuracy and restrict the application to specific conditions. For example, side draw products cannot be considered. o Surrogate models: mathematical or statistical models. These models define the relationship between the degrees of freedom and the variables of interest and are used in cases were computational time is a major concern.

The three modeling approaches for distillation column applications are thoroughly reviewed in the work of [START_REF] Ochoa-Estopier | The use of reduced models for design and optimisation of heat-integrated crude oil distillation systems[END_REF]. In this paper, the surrogate approach is preferred because in the methodology both solution accuracy and time constraints are considered to be equally important.

Distillation Column surrogate model

Having the material balance fixed, the distillation column has 4 degrees of freedom: number of trays (NT su ), pressure at the top of the column (P su Up ), Feed location and the n F column feed streams enthalpy (H n F ,su Feed ). A bottom-up approach was selected in order to identify sequentially the influence of each variable on the response of the model. Keeping constant values of the remaining three variables, we examined the effect of each one within a certain allowable value range dictated by the process standards. From a preliminary sensitivity analysis we concluded that feed location had the least impact on the output and it has been decided to simplify the model by considering that the feed stream is always located at the middle of the column. For each distillation column, the "numerical experiments" in the commercial simulator were performed with process simulation commercial softwares, with respect to equipment and operational constraints, for the following range of the design variables:

NT su min ≤ NT su ≤ 100 0.1 ≤ P su Up ≤ 12

(1)

H n 𝐹𝐹 ,su Fmin ≤ H n 𝐹𝐹 ,su Feed ≤ H n 𝐹𝐹 ,su

Fmax

where NT su min is the minimum NT su corresponding to the separation unit su, H n 𝐹𝐹 ,su Fmin and H n 𝐹𝐹 ,su Fmax are the minimum and maximum allowable values of the Feed streams. These three bounds are also calculated from surrogate models.

For the analysis and construction of the most complicated surrogate model of the distillation column, corresponding to Q su REB , more than 1000 independently and identically distributed initial input points were considered in order to fully investigate the effect of the independent variables on the final response. The model is split in two parts in order to respect the bounds on NT su as shown: The medamodel's response captures with adequate accuracy the actual data and the maximum approximation error corresponds to the edge of the search space as shown in the auxiliary diagrams (Figure 4), namely at the values of NT su close to NT su min . Even for this extreme case (min CAPEXmax OPEX) the absolute error is ~4 MW, which in relative error terms is less than 6 % and remains within the acceptable error margin.

Q su REB =P 1,su Q +P 2,su Q P su REB +P 3,su Q H 1,su Feed +P 4,su Q H 2,
The other complicated surrogate model calculates the column diameter in regards to two design variables and Q su REB : For larger columns (e.g. Xylene Column presented in a case study) the model is allowed to be extended up to the range of 25 m because it will indicate that splitting of the column is necessary at the design stage. From the auxiliary plots, the maximum error is less than 0.15m, which corresponds to 2-3%.

D su = P 1,su D + P 2,su D P su REB + P 3,su D Q su REB + P 4,su D P REB Q su REB + P 5,su D P su REB 2 + P 6,su D Q su REB 2 , su ∊ COL (
The remaining surrogates are modeled as simple cubic functions as follows:

M var = ∑ Par n M var 4 n=1 Var n-1 (4)
where respectively. A typical response of the cubic surrogates is illustrated on Figure 5. From all diagrams, the results are almost identical with the experimental data.

Phase change surrogate model

The mathematical formulations are characterized by extensive computational load, especially in large scale systems where the combinatory nature of the problem is a deterrent factor for applying rigorous MINLP formulations. Instead linear approximations are preferred, not only for the reduced computational time needed to obtain a solution but also due to the convexity that guaranties an optimal solution. The HEN synthesis problem can be easily modeled as MILP and the driving forces are expressed in terms of temperature difference:

Q = FCp(T in -T out ) (5)
The mass flowrate F can be easily calculated from material balance but the Heat capacity (Cp) is a nonlinear term. The majority of the current methodologies consider steams with constant heat capacity. That is, the latent heat involved in the mixed phase is not considered and as a result the phase diagram is a straight line with liquid and vapor phases only. This assumption simplifies the problem and effectively enables implementation of LP formulations but in a separation complex, the streams involved undergo phase transitions and the final solution can often be misleading, especially for the estimation of the heat exchanger area.

A very limited number of publications are dealing with this problem. The disjunctive methodology of Grossman et al (1998) extended the pinch location strategy of [START_REF] Duran | Simultaneous Optimization and Heat Integration of Chemical Processes[END_REF] to both isothermal and non-isothermal streams utilizing Boolean variables but not for HEN applications. [START_REF] Castier | Energy targeting in heat exchanger network synthesis using rigorous physical property calculations[END_REF] highlighted the fact that pinch points can be located within the intervals of Supply and Target Temperatures when Cp is considered as variable. They determined these points with rigorous thermodynamic calculations. Within the Pinch methodology framework also, [START_REF] Liporace | Heat exchanger network synthesis considering change phase streams[END_REF], proposed the decomposition of the phase diagram into 3 sub-sections; one for each phase based on dew and bubble points calculation. With the phase change considerations included, the pinch methodologies are more accurate and realistic but are still bounded from the sequential approach limitation. Almost 20 years later, Ponce-Ortega et al (2008) proposed another disjunction based formulation that extends the work of Yee and Grossmann to both isothermal and non-isothermal phase change transitions. The streams of the resulted MINLP methodology are separated in subsets according to their heat exchange requirements, namely sensible, latent and sensible-latent heat exchange. Therefore for each subset, the energy balance is calculated with three different models and integer constraint realizations. Despite the optimal network design, the methodology is limited from the fact that prior knowledge on the stream's conditions is necessary for the appropriate separation of the streams into subsets and consequently the methodology is not applicable to holistic process integration design where the stream conditions are examined simultaneously. [START_REF] Hasan | Operational modeling of multistream heat exchangers with phase changes[END_REF] modeled the phase change curve with a simple cubic function of ΔΤ with adequate accuracy but the calculations are restricted to predetermined pressure cases.

To satisfy the objectives of the proposed methodology, we need the accuracy of the three models corresponding to each phase and the numerical benefits of a single model. The sudden change in the slope takes place at the bubble and dew points which can be easily calculated, for both hot and cold streams, from the following set of equations:

H i Bub = P 1,i F + P 2,i F P i + P 3,i F P i 2 + P 4,i F P i 3 + P 5,i F P i 4 , i ∊ HP (6.a)

H i Dew = P 6,i F + P 7,i F P i + P 8,i F P i 2 + P 9,i F P i 3 + P 10,i F P i 4 , i ∊ HP (6.b)

The next step is the identification of a mechanism for the slope change, similarly to the integer variables but on a single model. Therefore, in the phase change surrogate model of Eq. 7, we take into advantage the switching properties of sigmoid functions for the transition from one phase to another. More specifically, both sigmoid functions are equal to zero at liquid phase and the Enthalpy-Temperature relationship is restricted to a quadratic function. At the bubble point, the first sigmoid changes value while the other remains zero and effectively changes the slope by contributing another term in the equation. The second sigmoid activates at dew point and causes the final slope change that indicates vapor phase.

T i S = P 13,i F + P 14,i F H i S + P 15,i F H i S 2 + P 16,i F �H i S -H i Bub � -P 15,i F �H i S 2 -H i Bub 2 � 1 + e �-P 11,i F �H i S -H i Bub �� + P 17,i F �H i S -H i Dew � + P 18,i F �H i S -H i Dew � 2 1+e �-P 12,i F �H i S -H i Dew �� , i ∊ HP (7)
Eq. 7 is a continuous function that guaranties, if well-tuned, positive ΔΤ/ΔΗ ratio for all arguments within the characterized range of 0 -500 °C. More importantly, as it can be seen from Figure 6, the surrogate model can also generate the phase change diagram for various pressures (all streams were characterized for 0.1 -20 barg) and thus, provides additional flexibility to the optimization algorithm to locate the most profitable conditions. The phase diagrams of a pure and multicomponent stream are presented in Figure 6. The surrogate model captures with the same efficiency the sharp slope change occurring on bubble and dew point for the pure component stream (left part of Figure 6) and multicomponent stream (right part of Figure 6). From the Parity and Absolute error diagram of Figure 7 we can observe that for both streams, the best fit corresponds to the liquid-mixed phase part and the worst to the vapor phase. In addition to the expected error at the bounds, significant error is observed at the dew points defined at lower pressures.

The absolute error is ~10 °C and corresponds to less than 7%. This is satisfactory, especially if we take into account the margin of estimation error when phase changes are not considered. Moreover, the high error values are located in a certain area where an additional corrective term might solve the issue but streams involved in distillation column complex are in vapor phase very rarely and hence, an instant measure for improving the efficiency is to neglect the vapor phase by removing the second sigmoid from the model.

Parameter tuning

Having identified the structure of the surrogates, it remains to determine the best tuned parameters, from the data generated from the rigorous computer experiments, using an optimization algorithm. For the elementary configurations (Figure 1), PRO/II 8.3 can generate up to 100 sampling points per 1 hour, each one corresponding to a unique value of the design variable set. This means that in one working day, a surplus of information is generated for all surrogate models. The regression is carried out in MATLAB to take into advantage the powerful computational capabilities it offers on dealing such problems. A classical least-square data fitting approach is used and the derived optimization problem is solved by the solver Lsqnonlin for the selection of the best parameters. Lsqnonlin uses the trust region-reflective algorithm, a Newton-based method that demonstrates convergences in less than 10 s for each tuning. With the parameters tuned, the MATLAB interface function wgdx passes the tuned parameters and initial values of the design variables to GAMS to start the optimization procedure.

MINLP formulation

The final step of the proposed methodology is the MINLP formulation for the optimization of the design variables, of each distillation column considered, in order to determine the most cost-effective solution. Process integration is achieved with the simultaneous solution of the three components of the formulation:

a) The surrogate models of the distillation columns b) HEN of the streams connecting the units between them c) Cost functions Moreover, the convention that stream pressure values are determined by the downstream elementary units is also applied. The Cost functions correlate the variables involved in the process integration with economic criteria and constraints and constitute the elements needed for the calculation of the TAC. The Cost functions were taken from [START_REF] Chauvel | Manual of process economic evaluation[END_REF] and proprietary software of IFP Energies Nouvelles.

HEN formulation

A modified stage-wise superstructure of [START_REF] Yee | Simultaneous Optimization Models for Heat Integration -II Heat Exchanger Network Synthesis[END_REF] is employed for the HEN synthesis. The superstructure is partitioned into stages for the investigation of potential heat exchange between pairs of hot and cold process streams as shown in Figure 8. The number of stages is selected arbitrarily but usually the maximum number is defined from the cardinality number of the hot and cold streams set. In general, the optimal configurations correspond to less heat exchangers and a small number of stages helps the algorithm to identify them faster. By convention, hot streams enter the superstructure at the first stage and the cold streams at the last. If the heat demand is not satisfied by process streams, the utilities located at the end of the superstructure cover the energy gap. It should be noted that in the superstructure, only one heat exchanger is allowed per stage for each stream. Figure 9 illustrates the allowable HEN configuration for a reboiler stream which is applicable to the condenser streams as well. A reboiler/condenser stream can exchange heat with process streams (in series for the process stream) and between them in parallel; a utility can be also used if needed. The restriction is imposed because mixed phase streams have very low difference in their Input -Output temperatures, which is practically zero for pure component streams, and heat exchangers in series are not recommended due to pressure drop effects and control issues. As the reboiler and condenser consist in the most energy demanding parts of the distillation column complex, they are the targets of heat integration and a HEN is often judged efficient or not by its ability of reducing the utilities usage for these streams. In the formulation, employment of different types of utilities is also examined for the investigation of the most cost-effective solution. The heaters dominate the OPEX and proper selection among the various types has a profound impact on economics. In this paper, we restrict the selection to Cylindrical and Cabin fuel heaters and high pressure steam (HPS), but if needed, more sophisticated solutions can be included. The Cylindrical heater is the cheapest but has limited range of duty (up to 35 MW); The Cabin is an intermediate alternative when compared to the expensive HPS. On the other hand, Air Cooler is the economically preferred source of cooling and in practice, a trim cooler is added in series for additional cooling of the streams required to achieve Battery Limit (BL) specifications. Trim cooler is not considered in the formulation because is limited to operate at low temperatures and it will increase the complexity of the formulation without any real benefit. Therefore, the cooling capabilities in the formulation are defined by the Air Cooler limits, namely by setting 55 °C as the Battery Limit Temperature. Moreover, we need to take into advantage streams with high temperatures for more 'valuable' heat exchange (e.g. with reboilers) instead of 'wasting' this potential between streams with large temperature difference. Thus, in the formulation, stream splits and heat exchange in parallel for the process streams are not permitted. This assumption also simplifies the algorithm and reduces the already demanding computational load. Based on the above mentioned decisions, the Enthalpy balance for hot and cold streams are formulated as follows:

H i S -H i,k=1 I = 0, i ∊ HP H j S -H j,k=nk I = 0, j ∊ CP H i,k I -H i,k O -� q i,j,k n j j=1 = 0, i ∊ HP, j ∊ CP ∪ REB, k ∊ ST H i,k O -H i,k+1 I = 0, i ∊ HP, k ∊ ST ` H i,k=n k O -H i AI = 0, i ∊ HP H i AI -H i AO -q i A = 0, i ∊ HP ∪ CON H i AO -H i T = 0, i ∊ HP ∪ CON H i S -� � q i,j,k n k k=1 n j j=1 -q i A -H i T =0, i ∊ HP, j ∊ CP∪REB, k∊ ST H j,k O -H j,k I -� q i,j,k n i i=1 = 0, i ∊ HP ∪ CON, j ∊ CP, k ∊ ST H j,k O -H j,k-1 I = 0, j ∊ CP (8) H j T -H j,k=1 O -q j H = 0, j ∊ CP H j T -� � q i,j,k n k k=1 n i i=1 -q j H -H j S = 0, i ∊ HP ∪ CON, j ∊ CP, k ∊ ST 𝐋𝐋 𝐂𝐂 Q su REB = � � q i,j,k n k k=1 n i i=1 + q j H , i ∊ HP ∪ CON, j ∊ REB, k ∊ ST, su ∊ COL 𝐋𝐋 𝐇𝐇 Q su CON = � � q i,j,k n k k=1 n j j=1 + q i A , i ∊ CON, j ∊ CP ∪ REB, k ∊ ST, su ∊ COL Q su CON = � H F,su Feed n F F=1 + Q su REB -H su TOP -� H sd,su Sdraw n sd sd=1 , su ∊ COL
where 𝐋𝐋 𝐇𝐇 and 𝐋𝐋 𝐂𝐂 are the matrices used to denote the links between unit and HEN variables. Sets of binary variables are introduced in the formulation to denote the decisions associated with Heat exchange between process streams (Z i,j,k ), hot stream -cold utility (Z i A ) and Cold stream -Hot utility (Z j H ) and setting the allowable duty range for the corresponding heat exchanger. Moreover, three additional binary variables corresponding to the different heater types are introduced:

Z i,j,k µ ≤ q i,j,k ≤ Z i,j,k M, i ∊ HP ∪ CON, j ∊ CP ∪ REB, k ∊ ST Z j H µ ≤ q j H ≤ Z j H M, j ∊ CP ∪ REB Z i A µ ≤ q A ≤ Z i A M, i ∊ HP ∪ CON Z j Cyl µ ≤ q j Cyl ≤ Z j Cyl M, j ∊ CP ∪ REB Z j Cab µ ≤ q j Cab ≤ Z j Cab M, j ∊ CP ∪ REB (9) Z j HPS µ ≤ q j HPS ≤ Z j HPS M, j ∊ CP ∪ REB q j H = q j Cyl + q j Cab + q j HPS , j ∊ CP ∪ REB Z j Cyl + Z j Cab + Z j HPS = Z j H , j ∊ CP ∪ REB
where μ corresponds the lower bound on the duty and M is the big-M parameter used for relaxing the constraints.

The enthalpy of each stream is also needed for phase change calculations (Eq.7). The pressure variation is not penalized as the cost involved is considered negligible. Binary variables are used to activate or deactivate the following constraints to ensure feasible driving forces for exchangers when they are selected as part of the network structure:

δT i,j,k HOCI =T i,k O -T j,k I +�1-Z i,j,k �M, i∊ HP, j ∊CP, k∊ ST δT i,j,k HICO = T i,k I -T j,k O + �1 -Z i,j,k �M, i ∊ HP, j ∊ CP, k ∊ ST δT j HICO HPS = T HP I -T j T + �1 -Z j HP �M, j ∊ CP ∪ REB δT j HΟCΙ HPS = T HP Ο -T j,k=1 O + �1 -Z j HP �M, j ∊ CP δT j HΟCΙ HPS = T HP Ο -T j 𝑆𝑆 + �1 -Z j HP �M, j ∊ REB δT i,j,k HOCI = T i T -T j,k I + �1 -Z i,j,k �M, i ∊ CON, j ∊ CP , k ∊ ST δT i,j,k HICO = T i S -T j,k O + �1 -Z i,j,k �M, i ∊ CON, j ∊ CP, k ∊ ST δT i,j,k HOCI = T i,k O -T j S + �1 -Z i,j,k �M, i ∊ HP, j ∊ REB, k ∊ ST (10) δT i,j,k HICO = T i,k I -T j T + �1 -Z i,j,k �M, i ∊ HP, j ∊ REB, k ∊ ST δT i,j,k HOCI = T i T -T j S + �1 -Z i,j,k �M, i ∊ CON, j ∊ REB, k ∊ ST δT i,j,k HICO = T i S -T j T + �1 -Z i,j,k �M, i ∊ CON, j ∊ REB, k ∊ ST δT i,j,k HOCO = T i,k O -T j,k O + �1 -z i,j,k �M, i ∊ HP, j ∊ CP, k ∊ ST δT i,j,k HOCO = T i,k O -T j T + �1 -Z i,j,k �M, i ∊ HP, j ∊ REB, k ∊ ST δT i,j,k HOCO = T i T -T j,k O + �1 -Z i,j,k �M, i ∊ CON, j ∊ CP, k ∊ ST δT i,j,k HOCO = T i T -T j T + �1 -Z i,j,k �M, i ∊ CON, j ∊ REB, k ∊ ST δT i HICI_A = T i AI -T Air , i ∊ HP ∪ CON δT i HICI_A ≥ 20, i ∊ HP ∪ CON
Geometric mean of three Temperature points (Ponce-Ortega, 2008) have been proposed to account accurate approximations for phase change exchangers but the implementation requires knowledge of the streams involved in phase change. Despite being accurate for isothermal phase changes only, Chen approximation (1987) is the preferred equation used for LMTD estimation, mainly because of the numerical stability that offers:

LMTD i,j,k Chen = �(δT i,j,k HOCI δT i,j,k HICO )( δT i,j,k HOCI +δT i,j,k HICO 2 � 1 3 � , i ∊ HP ∪ CON, j ∊ CP ∪ REB, k ∊ ST (11)
The Chen approximation is used for the calculation of the Heat exchanger area for process streams (PS) and the area of cold -HPS heater:

Area i,j,k PS = q i,j,k U i,j,k LMTD i,j,k Chen , i ∊ HP ∪ CON, j ∊ CP ∪ REB, k ∊ ST (12) A j HPS = q i HPS U HP LMTD j HPS , j ∊ CP ∪ REB
Moreover, the following equations and constraints are required to ensure monotonic temperature behavior across the stages.

T i S = T i,k=1 I , i ∊ HP T j S = T j,k=𝑛𝑛 𝑘𝑘 I , j ∊ CP T i,k I ≥ T i,k O , i ∊ HP T j,k I ≤ T j,k O , j ∊ CP T i S ≥ T i,k=1 I , i ∊ HP T j S ≤ T j,k=n k I , j ∊ CP (13) T i,k I ≥ T i,k O , k ∊ ST, i ∊ HP T j,k I ≤ T j,k O , k ∊ ST , j ∊ CP T i,k=n k O = T i AI , i ∊ HP T i AI ≥ T i AO , i ∊ HP ∪ CON T i T = T i AO , i ∊ HP T i AI = T i S , i ∊ CON T i AO = T i T , i ∊ CON
Finally, the supply enthalpy of the interconnecting streams is defined by the distillation column surrogate models. For the products, the temperature target and pressure are set by the battery limits specifications. The supply and target values of the reboiler and condenser streams are defined directly from the metamodel outputs:

H j S = 𝐋𝐋 𝐂𝐂 H su BOT , j ∊ CP, su ∊ COL H i S = 𝐋𝐋 𝐇𝐇 H su BOT , i ∊ HP, su ∊ COL H j S = 𝐋𝐋 𝐂𝐂 H su TOP , j ∊ CP, su ∊ COL H i S = 𝐋𝐋 𝐇𝐇 H su TOP , i ∊ HP, su ∊ COL T j S = 𝐋𝐋 𝐂𝐂 T su RI , j ∊ REB, su ∊ COL (14)
T j T = 𝐋𝐋 𝐂𝐂 T su RO , j ∊ REB, su ∊ COL T i S = 𝐋𝐋 𝐇𝐇 T su CI , i ∊ CON, su ∊ COL T i T = 𝐋𝐋 𝐇𝐇 T su CO , i ∊ CON, su ∊ COL
The solution of the MINLP problem offers the means to evaluate the best configuration through the corresponding economic objective TAC and furthermore the optimal design and operating variables are determined directly. In the proposed formulation, retrofit optimization can be also performed when a subset of the variables involved have fixed values. In the following section, both grassroot (full degrees of freedom) and retrofit optimization are presented.

The proposed methodology is applied on an aromatic complex that produces benzene and paraxylene from reformat (effluent of reforming unit). Both are petrochemical intermediates produced at commercial grade with high purities (higher than 99.8 % weight). The complex feed is consisted mainly by aromatic isomers with single ring and carbon numbers ranging from 6 to 10 and some aliphatic compounds that have to be removed. In general more than 40 isomers are present in the complex and the Reformate splitter initiates the separation process by splitting the light (C7 -cut) and heavy reformates. The aliphatic compounds are removed from the light reformate using an aromatic extraction unit and the remaining compounds are fed into a two-column complex for the fractionation of the aromatics from benzene to the heaviest ones. The heavy reformate split is mixed with the bottom stream of the toluene separation unit and fed to another two-column complex for the extraction of the economically valuable paraxylene. Several recycle streams feed the fractionation as shown in Figure 10. The application of the methodology on the aromatic complex has been carried out in 2 steps: initially a first simple case study, with 2 distillation columns for evaluation of different objectives, is examined and finally a scale-up case with 4 columns. Both cases were solved in GAMS 24.4 using the SBB/CONOPT solvers.

Figure 10: Illustration of the aromatic fractionation in an aromatic complex with the 2 case studies.

Case 1: 2 Columns

The first case study involves 9 streams and 2 columns:

• Benzene column (BZC): this column is fed by mono aromatics (S1) and produces at the top (S3) the benzene product at commercial grade. Heavier aromatics are sent to the bottom product (S5). • Toluene column (TOC): fed by the bottom of the benzene column, it separates toluene at the top (S7) and heavier aromatics. The toluene product is send to a catalytic unit while the bottom product (S9) goes to a clay treater for olefin removal.

In addition to the surrogate parameters for each stream and distillation unit, the methodology requires the input and output specifications (for storage or downstream unit requirements) of the complex, namely the supply and target conditions of the entering and exiting streams respectively as shown in Table 2. For the HEN synthesis, the output streams S3,S7 and S9, along with the two condenser streams (S2 and S6), are considered as heat sources for potential heat exchange between the two Feeds (S1 and S5) and the two reboilers (S4 and S8) in links with Figure 12. As mentioned in previous sections, the MINLP formulation is non-convex and the final solution is highly dependent on the initial values of the design variables NT su , P su Up and H n 𝐹𝐹 ,su Feed . In the derivation of the surrogate models, we assumed that the feed tray is located at the middle of the column and as a consequence, different initial values of H n 𝐹𝐹 ,su

Feed have small impact on the final value for given NT su , P su Up . Therefore, the average allowable value of H n 𝐹𝐹 ,su Feed was used in the multistart approach for the 64 cases corresponding to the combinations of the initial value sets of P su Up = {0.1, 4, 8, 12} and NT su = {35, 80} for the two columns. Moreover, in all case studies presented in this paper, the number of stages examined for the HEN synthesis is equal to the cardinality number of the Cold stream set. The selection was made in order to evaluate the capabilities of the formulation on finding the best solution for the full range of possible stream combinations involved in the HEN synthesis and at the same time to evaluate the performance in terms of time.

A screening methodology for identifying the most promising initial conditions is also considered by imposing time constraints of 30 minutes in the optimization software. From the 64 cases, 54 converged to a solution within the time limit and from them only 10 resulted solutions close to the best case of 11,2 M$, as illustrated in Figure 13. In this way, we are able to limit the optimization to a very small of 'good' initial values number (less than 1/6 th of the initial number of cases) which has a significant effect on the total time required for finding the best solution. The two most promising initial values were used for an optimization without time constraints and the best solution is found to be 10,98 M$ (1,96% difference) for approximately 90 minutes using a 2.67 GHz Intel core processor and 16 GB RAM (Figure 11). The corresponding HEN is illustrated in Figure 12. As expected, the algorithm identified that the reboiler streams are the total cost manipulators. The S8 stream has the higher temperature in the system and the only way to be heated is through expensive external heater. The algorithm increases the NT TOC in order to reduce the duty demand and consequently the amount of external heat needed. The trade-off between CAPEX increase and OPEX decrease favors this selection. On the other hand, for the BZC, the algorithm favors to operate the BZC column at the minimum pressure with relatively low NT in order to benefit from the low cooling demands of the condenser and top product and at the same time to increase the reboiler's duty, which is within a temperature interval that permits heat exchange with the TOC condenser (S6). Moreover, the TOC Feed stream S5 is heated with the bottom product of TOC (S9) without utilizing any external energy source for both streams. For this case study, a retrofit optimization is performed on an industrial representation of the distillation column complex. Likewise the first case study all streams leaving the system of four columns have fixed target conditions while internals stream's conditions depend on operating conditions of the columns. For the 19 streams involved in the complex, the cold stream set is consisted by the column feeds and the reboiler streams. All the specifications are listed in Table 4 and the design variables of the industrial application in Table 5.

The optimization is carried out with the conditions of SC case applied on the previous case study, namely the simultaneous consideration of the CAPEX-OPEX trade-offs with the design variables fixed to the corresponding values of the industrial application. Figure 13 illustrates the HEN configuration corresponding to the results of the optimization, listed in Table 6.

In all heat exchangers present in complex, XYC streams are involved in order to take into advantage the higher temperature profile of the XYC. The temperatures of BZC and TOC reboilers (T BZC RI = 144 °C and T TOL RI =170 °C) are significantly lower than the XYC condenser (T XYC CI = 248 °C) and as a consequence their heat demands are covered entirely from the XYC condenser. Neglecting air cooler costs, which are only a minor fraction of the heater costs, we can say that half of the complex operates without providing external energy and thus reducing the cost to the minimum. However, this is not the optimum configuration because the given design conditions yield a δΤ HICI between S13 and S18 of only 19 °C and thus the constraint of 15 °C of δΤ HOCO imposed for unit streams cannot be respected and instead an air cooler is used for cooling the remaining S13 and an expensive heater for S18. This also illustrates the importance of process integration in order increase the savings by selecting the proper design variables that set the conditions for heat exchange between S13 and S18 with a compromise on the CAPEX.

Conclusions

This paper has presented a general methodology for Process integration optimization that simultaneously takes into account the CAPEX-OPEX trade-offs. The methodology finds the optimum process design and operating conditions in addition to the HEN configuration of the interconnecting streams of the units involved. For the HEN synthesis, a continuous monotonic model is included in the formulation in order to accurately take into account the phase changes that might occur during heat exchange between streams.

The sequential solution strategy of the methodology benefits from the capacities in different fields of three softwares and reduces the overall computational time needed in large scale optimization problems. The surrogate models efficiently replace the computational demanding rigorous models of the commercial simulators and curries out the optimization procedure in one software that targets the solution directly without the computationally costly data interchange procedure that usually takes place when an interface of multiple softwares is considered.

For the optimization procedure, the only information required is the surrogate models parameters and the input-output specifications of the process for a given set of initial values of the optimization variables. Following this, the interaction of the process units, represented by the surrogate models, with the HEN of the interconnecting streams enables the formulation to treat the overall process as a 'black box' where neither prior heuristics nor engineering experience is needed in order to generate the optimum solution. The dependence of the proposed MINLP formulation on initial conditions is addressed with a multistart approach that is used as a screening technique for the identification of the most promising initial values and has a significant impact on reducing the overall computational time needed. The methodology is applicable to both grassroot design, when all degrees of freedom are optimized, and retrofit design in the case where a subset of the design variables is fixed to meet the requirements of the application of interest.

Finally, the methodology was applied on a distillation column complex with two case studies. Three different objectives were examined for the optimization problem on the first case study to illustrate that the simultaneous consideration of the CAPEX-OPEX trade-off is the best strategy for minimizing the overall cost. Both CAPEX and OPEX sub-objectives are lower in the grassroot case when compared to a retrofit case. The gain reaches 15% between the study sequential approach (operating conditions optimization, then column design with engineering design rules, and finally exchanger network optimization) and grassroots design (simultaneous optimization of all variables)

Acknowledgment

The license provided by GAMS Software GmbH for evaluating the available MINLP solvers is gratefully acknowledged.

Notation

Figure 1 :

 1 Figure 1: Solution sequence of the proposed methodology.

(

Figure 2 :

 2 Figure 2: Elementary configuration of a distillation column.

Figure 3

 3 Figure 3: Q su REB Surrogate model response for Benzene Column

  3) where P n D ,su D are the n D parameters of the surrogate model. The decision of including an output variable instead of using directly the basic design variables during the model construction was taken for two reasons: a) the model would have been enormous in size b) to avoid numerical difficulties in the optimization algorithm Similarly with Q su REB , the Diameter model shows very good performance and the error corresponds to the higher values of Q su REB .

Figure 4 :

 4 Figure 4: Response of the Diameter surrogate model with auxiliary diagrams

Figure 5 :

 5 Figure 5: Response of the 𝑇𝑇 𝑠𝑠𝑠𝑠 𝐶𝐶𝐶𝐶 surrogate model with auxiliary diagrams

Figure 6 :

 6 Figure 6: Phase change model response for the various Pressure and Temperature cases. Left: Pure component stream. Right: Multicomponent stream.

Figure 7 :

 7 Figure 7: Auxiliary diagrams for the phase change model. Left: Pure component stream. Right: Multicomponent stream.

Figure 8 :

 8 Figure 8: Representation of the stagewise superstructure for the interconnecting streams.

Figure 9 :

 9 Figure 9: Representation of the stagewise superstructure for unit streams

Figure 11 :

 11 Figure 11: Left: Multistart approach results for the 2 Column case with 30 minutes time constraint.Right: Effect of computational time.

  ). The rest is recovered at the top and is mixed with toluene to feed the transalkylation unit.

Figure 13 :

 13 Figure 13: HEN Configuration for the 4-Column case.

  j,k : indicates if a Heat exchanger exists between Hot stream i and Cold stream j at stage k. Z i A : indicates if an AirCooler exists for Hot stream i.

  a Heater exists for Cold stream j. Z j Cyl : indicates if a Cylindrical Heater exists for Cold stream j. High pressure steam is used for heating Cold stream j. Matrices L H : n i x n su Matrix L C : n j x n su Matrix

Table 1 :

 1 Distillation column parameters and variables

	Design Variables	Description

Table 2 :

 2 Input and Output specifications for the 2 Column case

	Stream	Conditions value
	S1	H S	21,127 MW
	S3	T T	55°C
		P	5 barg
	S7	T T	110°C
		P	20 barg
	S9	T T	175°C
		P	20 barg

Table 4 :

 4 Specifications of exiting and entering streams of the 4-column complex

	Stream Specification value
	S1	H S	21.13
	S3	T T	55
		P	5 barg
	S7	T T	110
		P	20 barg
	S9	T T	175
		P	20 barg
	S11	H S	44.47 MW
		H T	48.30 MW
	S14	T T	158
		P	17 barg
	S19	T T	55
		P	5 barg

Table 5 :

 5 Design variables corresponding to the industrial application of the 4-columns case.

	Design variable	value
	P BZC Up (barg)	0.1
	P TOC Up (barg)	0.1
	P XYC Up (barg)	0.8
	P HAC Up (barg)	0.1
	NT BZC	41
	NT TOC	29
	NT XYC	85
	NT HAC	27

Table 6 :

 6 Optimization results of 4-Column case study

	Variable	Grassroot
	TAC (M$)	37.28
	OPEX (M$)	27.14
	CAPEX (M$)	10.14
	H 1,BZC Feed (MW)	21.13
	H 1,TOC Feed (MW)	22.57
	H 1,XYC Feed (MW)	71.31
	H 2,XYC Feed (MW)	40.00
	H 1,HAC Feed (MW)	19.90
	E1 (MW)	37.40
	E2 (MW)	46.18
	E3 (MW)	11.15
	E4 (MW)	3.83

Figure 12: The HEN configuration of the optimal solution for the 2 column case.

The simple case considered, also enables a performance comparison of the proposed simultaneous process integration methodology with a classic sequential approach (SA). Based on heuristics, we concluded that the desirable heat exchange between S6 and S4 is achieved with an increase of TOC pressure to 4 barg. Then, engineering design rules recommended 41 and 39 theoretical trays for BZC and TOC respectively. Using the proposed formulation for the specified values of P su Up and NT su , the sequential approach is performed in two steps:

• HEN optimization according to the basic energy objective of minimizing the utility demand plus the number of heat exchanger present in the configuration which is added as a counterweight for the exchanger's CAPEX. • Detailed calculations for the determination of TAC, based on the previous step results

The first step is representative of the most common approach in industry due to the simplicity and the fast results, while the second step optimizes cost objectives based on the resulted HEN and heuristics determined operating conditions. It is noteworthy that using the proposed formulation for the HEN optimization, the effect of feed stream enthalpy and phase change is also considered, in contrast to the majority of the already existing methodologies. Having the process design variables determined, the TAC objective function of the proposed formulation corresponds to the simultaneous consideration of the sequential approach objectives and can be used for retrofit optimization purposes. Table 3 lists the optimization results of the best solution of the proposed methodology with full degrees of freedom (grassroot), the sequential approach (SA) case and simultaneous case (SC).

All approaches identified that the optimal HEN is consisted by 5 heat exchangers without any change on the H 1,BZC Feed energy level. When all degrees of freedom are considered in the optimization, the maximum level of interaction is exploited and as a result, both CAPEX and OPEX sup-objectives are more cost-effective. This is also highlighted by the fact that heuristics and engineering design rules are able to identify near-optimal pressure values but fail when it comes to the selection of NT. The 41 and 39 trays are plausible choices for reducing the CAPEX but the 5 and 34 extra trays, along with the other optimization variables selected by the optimal case, compensate the additional cost by setting the proper conditions for reducing the values of other capital cost associated variables such as heat exchanger area and thickness. Moreover, the high NT of the TOC is also used for minimizing the demand of the only hot utility needed in the resulting HEN for the heating S9 and has a profound impact on the final OPEX value.

On the other hand, the sequential approach objective is satisfied as the OPEX value targeted is lower than the respective of the simultaneous case (SC) but the 0.03M$ save is overlaid by the 0.72M$ on the annualized CAPEX that corresponds to a 5.47% of the TAC per year. It should be noted that the CAPEX is annualized with a 6 years payback period.

Case 2: 4 Columns

To enlarge the problem we have integrated the 2 remaining columns and the attached streams:

• Xylene column (XYC): this column is used to recover the A8 cut at the top (S14). This cut is sent to the paraxylene extraction unit that has stringent constraints on heavy aromatics containeds. The column is fed by two streams : o The first feed (S10) is the mixture of the TOLC bottom stream and the heavy split of the reformer. o A8 isomerate from isomerization unit (S12). The supplied enthalpy is also fixed.

• Heavy aromatic column (HAC): depending on the targets specified at the design stage, a given fraction of the heaviest aromatics is separated at the bottom of this column and sent to Battery