A Robust Decomposition Methodology for Synthesis of Flexible Processes with Many Uncertainty Parameters – Application to HEN Synthesis

This contribution presents a new robust decomposition methodology for generating optimal flexible process flow sheets with a large number of uncertain parameters. During the initial steps, first-stage variables are determined by performing mixed-integer nonlinear programming (MINLP) synthesis of a flow sheet at the nominal conditions, and then by exposing the obtained flow sheet sequentially over a set of extreme MINLP scenarios of uncertain parameters. As a result, the sizes of the flow-sheet units gradually increase, and/or new units are added until the required feasibility is achieved. After testing the flexibility of the obtained design, a Monte Carlo stochastic optimization of the second-stage variables is performed using a sampling method in order to obtain an optimum value of the expected objective variable. The advantages of the proposed methodology are the independence of process model sizes from the number of uncertain parameters, the straightforward use of deterministic models for incorporating uncertainty, and relatively simple execution of MINLP synthesis of processes under uncertainty. Thus, it could be used for designing large processes with a large number of uncertain parameters. The methodology is illustrated by synthesis of a flexible Heat Exchanger Network.


Introduction
There are several difficulties connected with the presence of variating or undetermined input parameters, such as identification of uncertainties, availability of data, the types of variables in the system, and especially the availability of appropriate strategies and large computational efforts for solving problems with many uncertain parameters. 1 Synthesis of chemical processes under uncertainty in general involves the following tasks: identifying a flexible process flow-sheet structure, determining the optimal values of design and operating variables, and optimizing the expected value of the objective function so that flexible and lifetime optimal processes will be generated.Such problems are often formulated as two-stage stochastic programs with recourse, in which variables are partitioned into first-and second-stage variables. 2First-stage variables are those related to process topology, i.e., binary variables and sizes of process units that are determined in advance before the values of uncer-tain parameters are realized.Second-stage variables are determined during operation when uncertain parameters have known values, and control variables are adjusted to achieve feasible and optimal solutions.Another approach is by solving a one-stage optimization problem. 3 Both approaches often rely on various sampling-based techniques for evaluating multiple integrals that give the expected values.Such procedures are very expensive to solve because the number of samples greatly increases with the number of uncertain parameters. 4everal methods have been proposed for reducing the number of scenarios; for example, Karuppiah et al. determined a smaller set of scenarios that approximate the optimization problem in a reduced space. 5Novak Pintarič et al. proposed an algorithm for identification of critical scenarios. 6Martin and Martinez applied scenario reduction and sample approximation approaches to problems dealing with formulated products and process design. 7In some studies, the number of trials within Monte Carlo simulations was determined in such a way that sufficiently small errors were obtained. 8 general or standard approach for incorporating uncertainties into design and synthesis of chemical processes has yet to be established; however, many researchers have developed specific algorithms.Bahakim et al. approximated process constraints using Power Series Expansion functions. 9teimel and Engell developed a computer tool that transforms the process superstructure problem into a two-stage stochastic model. 4rocess systems that have drawn a lot of attention regarding uncertainty are Heat Exchanger Networks (HENs).Many authors have used decomposition techniques for solving two-stage stochastic HEN models; for example, a synthesis of flexible HENs was performed by solving multi-period problems by gradually adding extreme periods to the nominal point. 10This procedure was later simplified using a simulation-based method to realize the flexibility tests instead of the model-based method. 11scobar et al. performed a synthesis of flexible and operable HENs by solving a deterministic equivalent through discretization of uncertainty. 12In our previous work, a robust computational methodology for the synthesis and design of flexible Heat Exchanger Networks (HEN) with a large number of uncertain parameters was developed. 13A single period model of HEN was employed to reduce the computational efforts, and promising matches were introduced to enhance solution generation. 14 Gu et  al. proposed an MINLP model for minimizing control action to identify the inactive bypasses and active pairings simultaneously. 15safiade and Fraser presented a new superstructure, named interval based MINLP superstructure, for multiperiod HEN synthesis. 16Isafiade et al. developed a reduced multiperiod HEN model by introducing the stream matches obtained by solving multiperiod HEN models at different minimum approach temperatures and different number of stages. 17ilva et al. applied Particle Swarm Optimization, wherein all variables were optimized simultaneously. 18heng et al. proposed an approach for flexible HEN synthesis under severe operation uncertainty, represented by the probability bounds analysis theory and sampled with a double-loop sampling method. 19Li et al. presented a two-step approach for the synthesis of flexible HENs. 20An example with 11 uncertain parameters was solved by applying a stepwise optimization method also suitable for nonconvex network problems.Li et al. presented a methodology for the synthesis of flexible heat exchanger networks where for large non-convex heat exchanger problems, the degree of flexibility was established as well as the direction of deviations of uncertain parameters using a direction matrix. 21ased on a literature survey, it can be concluded that many approaches have been developed for designing flexible process flow sheets, and flexible HENs in particular; however, most of them deal with problems involving a limited number of uncer-tain parameters and/or a low number of operating periods.The number of uncertain parameters in real industrial problems is large, and operating conditions vary continuously.Therefore, the generation of flexible processes with a large number of uncertainty parameters for the entire lifetime and the entire range of possible uncertain parameters values is necessary.
To accomplish this task, this study aims to develop a robust methodology for considering a large number of uncertain parameters during the synthesis and design of process systems.In order to avoid exponential increase of mathematical models with the number of uncertain parameters, the main purpose of this methodology is to handle a limited number of scenarios through sequential iterations, i.e., within the loop, rather than simultaneously.Only one-scenario, or at most two-scenario models, are solved at each iteration, thus keeping the model size under control.A four-step methodology is described further herein, and illustrated by a case study in two variants.

Problem statement
Synthesis of process flow sheets involves discrete decisions, i.e., optimal selection of process units from the superstructure of alternatives, and continuous decisions, i.e., operating and control variables.The program (P1) mathematically describes a mixed integer nonlinear programming (MINLP) model, which contains, beside discrete and continuous variables, several input parameters, including economic, model, and process parameters.Many of these are subjects of uncertainty and therefore multi-scenario stochastic approaches to process flow-sheet synthesis are required.

{ }
In (P1), the subscript s denotes index of discrete scenarios, S is a set of selected scenarios, E(Z) is an expected value of scalar objective variable Z, p probability of scenarios, c fixed costs of alternatives, f variable cost function, y binary variables representing process topology, d design variables representing capacities and sizes of process units, x s ∈ S operating variables, z control variables, h equality constraints, g inequality constraints, g d design expressions, q uncertain parameters that vary between lower (q LO ) and upper (q UP ) bounds, and A, B, and a are the matrices and vector of constants.
(P1) represents a two-stage stochastic problem with recourse, in which first-stage variables, such as process topology and design variables (y and d), are equal for all scenarios and determined in advance during the process design phase, while second-stage variables, i.e., operating and control variables (x, z), are adjusted recursively during process operation in order to obtain a feasible and optimal solution for each scenario (q s ).The set of selected scenarios S in (P1) is not uniquely defined.
Various existing methods have in common the fact that the number of scenarios increases hugely with the number of uncertain parameters, and thus, the model in (P1) could become unsolvable.The methodology proposed in this paper decomposes the MINLP process synthesis problem under uncertainty (P1) into several steps that are described in the following section.The main advantage is that one-scenario (or at most two-scenario) problems are solved sequentially in the loop, thus, i) avoiding multi-scenario models the sizes of which would increase drastically with the number of uncertain parameters, and ii) providing good initial points for optimizations in subsequent iterations.Another advantage is that deterministic process models can be applied directly for handling uncertainties with no specific preparations, modifications or adjustments, which gives certain robustness to the proposed approach.
The hypothesis is that such an approach would be suitable to apply to large-scale process synthesis problems containing a large number of uncertain parameters.The assumptions are that the ranges of variations of uncertain parameters are known, that uncertain parameters are mutually independent, and that the extreme values of uncertain parameters are critical for feasibility.The latter assumption arose from our experience that, in most chemical process models, certain assumptions about critical extreme points would be appropriate.
A limitation of this approach is that it does not rely on exact theory-based simultaneous procedures, but rather deals with sequentially executed finite numbers of randomly selected scenarios, and therefore, global optimal solutions cannot be guaranteed even for convex problems.However, it can be assumed that even with this limitation, good optimal or near optimal process flow sheets can be obtained by performing a sufficient number of iterations in randomly selected scenarios, the number of which should be kept as low as possible.

Methodology
The main idea of the methodology is to start with an optimal process flow sheet obtained by solving a one-scenario MINLP problem at the nominal values of uncertain parameters.The flow sheet thus obtained is usually inflexible and then sequentially exposed to deviations of uncertain parameters towards their extreme values (vertices) in order to increase its flexibility by: i) enlarging the sizes of the existing process units, and ii) adding those additional process units required for feasibility.In this way, the initial process flow sheet is gradually extended either by increasing the sizes of those units already selected in previous iterations or by adding new process units for assuring feasibility.Through iterations, the initial nominal process flow sheet gradually increases and transforms from inflexible to flexible.The goal is to add only as many additional process units as required for pre-specified deviations of uncertain parameters.The four main steps are as follows (Fig. 1): Step 1 -Generation of initial process flow sheet at nominal conditions.
Step 2 -Determination of first-stage variables for a flexible process flow sheet.Two alternative options were studied: -Option 1 -One-scenario approach, considering only one vertex point at each iteration.-Option 2 -Two-scenario approach, considering simultaneously a nominal point and a vertex point at each iteration.
Step 3 -Determination of flexibility index, and corrections of design variables if needed.
Step 4 -Determination of second-stage variables through stochastic optimization of the flexible process flow sheet derived in Step 2.
The motivation for testing two options in Step 2 originates from the fact that one-scenario problems are easier to solve, but on the other hand, it can be assumed that a two-scenario approach will produce better results due to the presence of both the nominal and the extreme points, which will establish better trade-offs between the first-and second-stage variables during the determination of process topology and process unit sizes.
Step 1 -Generation of initial process flow sheet at nominal conditions The goal of this step is to obtain a process structure at the nominal conditions which, although most probably inflexible, can serve as a good starting point for the generation of a flexible process flow sheet in the second step.Initial process struc-ture is generated by solving the one-scenario MIN-LP problem (P2), in which all uncertain parameters are fixed at their nominal values.
( ) , , 0, 0,1 where superscript N denotes the nominal value.The result of this step is an optimal nominal process flow sheet, which is most likely inflexible for deviations of uncertain parameters from their nominal values.
Step 2 -Determination of first-stage variables The main goal of this step is to determine the first-stage variables, i.e., the topology binary variables and the sizes of process units for feasible operation within specified deviations of uncertain parameters.This step starts with the initial process F i g . 1 -Methodology for synthesis and design of optimal flexible process flow sheets flow sheet generated at the nominal conditions in the previous step.Those process units selected in the previous step are retained in the flow sheet by fixing their binary variables to one, while the sizes of these units are limited downwards by those values obtained at the nominal conditions.
At the beginning of this step, a set of random vertices should be generated (S RSV ).From this point, the process can proceed in two ways: either by a one-scenario approach, where only a vertex point is considered at each iteration, or by a two-scenario approach, where a nominal point and vertex point are considered simultaneously.It is supposed that the approach with two scenarios, from which one is always the nominal point, would produce first-stage variables closer to the optimal stochastic result, because the latter is most probably nearer to the nominal point than to any extreme vertex point.

One-scenario approach
The one-scenario MINLP model is presented as problem (P3).In each subsequent iteration, this problem is solved at a new vertex, with those binary variables that obtained unity values in previous iterations fixed to 1, while the remaining binary variables are optimized either to 1 or 0. Design variables are limited downwards by the optimal values obtained in the previous iteration.In this way, fixed topology options are forced to increase their sizes first, and only then are the new topology alternatives included in the flow sheet when further increases in sizes are insufficient for achieving flexibility at a new vertex.Besides, a solution at each iteration provides a good initialization point for subsequent iterations.
( ) ( , , , ) where s is an index of random vertices within the set S RSV , and y 1 s-1 is a vector of binary variables that obtained unity values at previous iterations and are fixed to 1 at the current iteration.The fixed cost in the objective function accounts only for those binary variables, i.e., process units that are added at ver-tex s, while unity binary variables from previous vertices are not considered.Note that the main idea behind (P3) is to minimize the extension of the current flow sheet necessary to achieve required flexibility at a new vertex.

Two-scenario approach
The two-scenario MINLP model is presented as problem (P4).At each iteration, this problem is solved simultaneously in two scenarios: the first corresponds to the nominal point (q N ), while the second stands for a randomly selected vertex (q s ).The second-stage variables in the objective function are assumed at the nominal conditions, because the stochastic result should be closer to the nominal point than to the extreme vertex point.The size of the problem (P4) would be doubled in comparison to the one of (P3); however, it can be expected that this size would still be manageable in most cases, because exponential growth of the model would be avoided anyway.

Z c y y f d x z h d x z h d x z g d x z By g d x z By d g x z d g x z d d
Ay a d x z x z y y y q q q q q q q − − − The constraints on the left side stand for the nominal point, while those on the right stand for the vertex point.Again, design variables are limited downwards by the values from the previous iteration, and the binary variables of already selected units are fixed to one.
Within the proposed methodology, either problem, (P3) or (P4), is solved successively at generated vertex points until new topology options are added and/or the design variables increase.It is assumed that flexible structure is achieved when the flowsheet structure and the values of design variables stop changing, so usually only a minor part of randomly generated vertices within the set S RSV would be used in this step.The results are the values of the first-stage variables, i.e., process topology and design variables that are most likely flexible for prescribed variations of uncertain parameters.The effectiveness of one-and two-scenario approaches is demonstrated in the Results and discussion section.The process flow sheet obtained in the previous step needs to be tested for flexibility, and flexibility index is determined as that defined by Swaney and Grossmann. 22Flexibility index is evaluated sequentially as a one-scenario NLP problem (P5) over the entire set of randomly selected vertices (S RSV ) for fixed process flow-sheet structure and design variables obtained from the previous step, y opt and d opt .This problem corresponds to a nonlinear problem (NLP), because of fixed binary variables.x z d q q q q q q q d d d In (P5), δ represents a scaled deviation of uncertain parameters from their nominal values, and sgn represents the directions of deviations from the nominal point towards the extreme values, either in negative or positive directions; i.e., sgn values -1 or +1 represent the directions towards lower bound θ LO and upper bound θ UP , respectively.The deviation d is limited by an upper bound in order to prevent unbounded solutions.
The flexibility index (I flex ) is then determined as the minimum of δ s : The value of the flexibility index I flex should be greater than or equal to 1 if the design, y opt and d opt , is to be feasible for the predetermined ranges of uncertain parameters.Based on those deviations, d obtained at randomly selected vertices, the mean value, and the standard deviation can be calculated.If the sample size is large enough, a normal distribution of δ can be assumed, and a confidence level for flexibility index equal to or larger than 1 can be evaluated.
Experience has shown that, in some cases, feasible solutions can be obtained at all testing vertices, but the flexibility index would be slightly lower than 1.A correction problem (P6) can be used to slightly modify the obtained values of design variables in order to achieve a flexibility index of 1.In this model, positive slack variables d sl opt are added to the optimum values of design variables d opt .The slacks are multiplied by a large number M in order to be minimized within the objective function, together with the maximization of d, under the additional constraint that d should be at least 1.The slacks then represent the shortages in the sizes of process units that should be added to optimal values of the previous step in order to reach a flexibility index equal to 1. x z d Step 4 -Determination of second-stage variables Finally, a stochastic Monte Carlo optimization is performed for the process flow sheet and process unit sizes obtained in Step 2 in order to determine the values of the second-stage, i.e., operating and control variables.A single-scenario problem (P7) is solved within a loop for a set of randomly selected points according to distribution functions of uncertain parameters (S MC ) at the fixed values of binary and design variables, y opt and d opt .It should be noted that the set S MC contains inner points from the entire region of uncertain parameters, while the set S RSV used in previous steps contains the extreme values of uncertain parameters, i.e., vertices.As in the previous step, the problem (P7) corresponds to a NLP optimization because binary variables are fixed.rived.If the sample size is large enough, it can be assumed that the objective variable is normally distributed, and its expected value can be expressed as the mean value: MC ( ) where N is the sample size, and E(Z) the expected value.Standard deviation can be evaluated, followed by a determination of the error margin and level of confidence.

Results and discussion
The proposed methodology was tested on a flexible Heat Exchanger Network (HEN) synthesis consisting of 6 hot and 3 cold streams, as defined in Table 1, with their supply temperatures (T s ), target temperatures (T t ), heat capacity flow rates (CF), and heat transfer coefficients (α).
The network is a MINLP model of a multistage superstructure developed by Yee and Grossmann. 23The number of stages within the superstructure was set at six.The objective function was minimum total annual cost (TAC) composed of the utility costs plus annualized investment of heat exchangers.
A HEN system was selected as an example because it has clearly defined first-and second-stage variables, and there are strong interactions and trade-offs between them.First-stage variables were binary variables for selection of heat matches, and areas of heat transfer units.Hot and cold utilities consumption, as well as the intermediate temperatures within the network, were treated as second-stage variables.Altogether, 22 uncertain parameters were defined: -supply temperatures of 9 process streams with ranges of variations ± 10 K from their nominal values, -supply temperatures of 2 utilities with ranges of variations ± 10 K from their nominal values, -heat capacity flow rates of 9 process streams with ranges of variations ± 5 % from their nominal values, and -prices of 2 utilities with ranges of variations ± 5 $ kW -1 y -1 for hot utility and ± 10 $ kW -1 y -1 for cold utility.
To solve the mathematical models defined in the previous section, the General Algebraic Modeling System (GAMS) was used as the solution tool with a DICOPT solver for treating the MINLP problems, and CPLEX for MILP and CONOPT3 for NLP (sub)problems.
Step 1 -Nominal HEN structure In the first step, the optimal HEN topology was derived by solving (P2) at the nominal values of uncertain parameters.It was a threshold problem consisting of 7 process to process matches, 2 heaters, and no cooling (Fig. 2).The total heat transfer area was 1675 m 2 and the total annual cost amounted to 1.207 M$ y -1 .The size of the one-scenario model was 677 constraints and 636 variables, of which 107 were binary.The CPU time per iteration on an average personal computer was around 0.1 s.
The nominal network was incapable of tolerating required deviations of uncertain parameters from their nominal values.Therefore, it was exposed to extreme deviations of uncertain parameters in the second step in order to increase the existing process units and/or add new ones until a flexible structure was obtained.
Step 2 -Generation of flexible HEN structure This step was accomplished in two ways: i) as a one-scenario problem (P3) considering only one vertex point at each iteration, or ii) simultaneously considering the nominal point and vertex point at each iteration, and solving two-scenario MINLP problems (P4).

One-scenario approach
In the one-scenario approach, HEN structure and the areas of heat exchangers stopped changing after approximately 1100 iterations, meaning that Ta b l e 1 -Nominal data for illustrative example of HEN synthesis Cost of heat transfer units ($ y -1 ) = 1800•A 0.65 (A in m 2 ) Price of hot utility = 80 $ kW -1 y -1 Price of cold utility = 20 $ kW -1 y -1 1100 randomly selected vertices were used.The final HEN (Fig. 3) consisted of nine units from the nominal structure (No. 1-9) to which three additional heat exchangers (No. 10-12), five coolers (No. 13-17), and one heater (No. 18) were added.
The area of heat transfer units from step one was increased by 34 %, while the area of newly added units (No. 10-18) was 410 m 2 , giving the total area of the extended network as 2658 m 2 , which corresponds to a 58 % increase.In the two-scenario approach, HEN structure and the areas of heat exchangers stopped changing after 1000 iterations, meaning that 1000 randomly selected vertices were used together with the nominal point.The size of the model was doubled in comparison to that of the one-scenario problem, and consisted of 1300 constraints and 1100 variables, of which 107 were binary.The CPU time for solving the two-scenario problem within one iteration was around 1 s.
Fig. 5 shows how the total areas of both networks in Figs. 3 and 4 increase as the number of iterations increases.It is evident that the two-scenario approach generated HENs with lower total area, and a considerably lower number of iterations, i.e., vertices, were required for stabilizing the HEN area.
Step 3 -Flexibility index of HEN structures The flexibilities of both HEN structures obtained in Step 2 were confirmed by (P6) over 5,000 randomly selected vertices.Fig. 6 presents the gradual increase in the flexibility index towards 1 with the increasing number of vertices used in Step 2.
The flexibility index of the final one-scenario HEN structure was 1.00, with a mean value of scaled deviation d of 1.65 and standard deviation of 0.34.We estimated the probability of the flexibility index being larger than or equal to 1 at 97 %.The flexibility index of the two-scenario structure was 1.00 with the mean value of scaled deviation d of 1.35 and standard deviation of 0.42.It was estimated that the probability of the flexibility index being larger than or equal to 1 was 80 %.The confidence in the results obtained would be sufficient for decision-making in practice.
Step 4 -Stochastic optimization of HEN In this step, a Monte Carlo optimization was performed for a fixed HEN structure and unit sizes as shown in Figs. 3 and 4 over 5,000 randomly selected points of uncertain parameters.Normal distribution was assumed for the supply temperatures and heat capacity flow rates.It was established that all solutions were feasible.The results demonstrate that using a two-scenario approach, a HEN with lower expected TAC was obtained, and fewer iterations were needed to generate a flexible HEN structure than with the one-scenario approach.The hypothesis is thus confirmed, that the presence of nominal and vertex points during determination of first-stage variables will establish better trade-offs and provide better

Conclusions
This study proposed a new robust procedure for generating optimal flexible process flow sheets represented as two-stage optimization problems with recourse, and containing a large number of uncertain parameters.This methodology decomposes the two-stage problem into determination of the first-stage variables followed by flexibility testing and determination of the second-stage variables.First-stage variables are determined by solving one-scenario (or at most two-scenario) MINLP problems starting from the optimal nominal process topology, while increasing the sizes of the existing process units and adding new ones to achieve feasibility over a sufficient set of extreme deviations of uncertain parameters.Second-stage variables are determined with Monte Carlo stochastic optimization by solving one-scenario NLP problems over a set of randomly selected points.The main advantage of this approach is that one-scenario problems are solved sequentially in a loop, thus making the sizes of mathematical models independent of the number of uncertain parameters.
The methodology was demonstrated on the MINLP synthesis of a flexible Heat Exchanger Network containing 22 uncertain parameters.It was confirmed that optimal flexible solutions were obtained with relatively low modeling and computing efforts.The level of confidence achieved would be suitable for quality decision-making in practice.Future efforts should be oriented towards further improvement of this method, and testing on various large-scale process and supply chain systems.

Step 3 -
Determination of flexibility index

≥
Optimal values of operating and control variables are determined for each scenario, and the expected value of the objective function is then de-(P5) ∀s ∈ S RSV (P6) ∀s ∈ S RSV (P7) ∀s ∈ S MC

F i g . 2 -
Nominal HEN design F i g .3 -Optimal HEN design obtained with one-scenario approach F i g .4 -Optimal HEN design obtained with two-scenario approach Two-scenario approach

Fig. 7
Fig. 7 presents how the expected values of TAC increase with the number of vertices used in Step 2 for both HEN structures.The final expected value of TAC for a one-scenario HEN was 1.482 M$ y -1 , and for a two-scenario HEN, 1.364 M$ y -1 .The results demonstrate that using a two-scenario approach, a HEN with lower expected TAC was obtained, and fewer iterations were needed to generate a flexible HEN structure than with the one-scenario approach.The hypothesis is thus confirmed, that the presence of nominal and vertex points during determination of first-stage variables will establish better trade-offs and provide better

F i g . 5 -
Total area of HENs obtained with the one-and two-scenario approaches F i g .6 -Flexibility index of HENs obtained with one-and two-scenario approaches solutions, while keeping the sizes of process models under control and independent of the number of uncertain parameters.
N o m e n c l a t u r e A b b r e v i a t i o n s HEN -Heat Exchanger Network MINLP -Mixed Integer Non-Linear Programming NLP -Non-Linear Programming TAC -Total Annual Cost I n d i c e s s -vertex point S u b -a n d s u p e r s c r i p t s equality constraints I flex -flexibility index M -large penalty scalar N -sample size N RSV -number of randomly selected vertices sgn -sign, direction (+1 or -1)

F i g . 7 -
Expected TAC of HENs obtained by one-and two-scenario approachesT-temperature, K x -vector of operating variables y -vector of binary variables y1 -vector of binary variables with unity value z -vector of control variables Z -scalar objective variableG r e e k l e t t e r s a -individual heat transfer coefficient, kW m -2 K -1 Δ -difference d -scaled deviation of an uncertain parameter q -vector of uncertain parameters