Uncertainty Analysis for Process Design and Scheduling Marianthi - PowerPoint PPT Presentation

Surface Reconstruction Ideas Problem definition of surface reconstruction: Given a set of sample points, determine the shape formed by these points Identify points constituting the boundary of the data set � Join boundary points to reconstruct the surface � Analogous to problem of feasibility analysis Determine mathematical representation of the boundary of the feasible region PASI 2005 August 16-25 Iguazu Falls

Improved Feasibility Analysis by α − shapes Given a set of points, determine the shape formed by these points Eliminate maximum possible circles of radius α without eliminating any data point H. Edelsbrunner, 1983 α α → For the α shape degenerates to the 0 α original point set α For α → ∞ the α shape is the convex hull of the original point set (Ken Clarkson http://bell-labs.com/netlib/voronoi/hull.html) PASI 2005 August 16-25 Iguazu Falls

Selection of α value for α − shapes Value of α controls the level of details of the constructed surface. α is a function of sample size ( n ) α is a function of inter-point distance Determination of α value Mandal& Murty, 1997 Construct minimum spanning tree (MST) of sampled data points Evaluate L n = sum of Euclidean distance between points of the MST L n α value = n PASI 2005 August 16-25 Iguazu Falls

Algorithm for Feasibility Analysis Sample the feasible space Construct MST of sampled data to determine α value Points constituting the boundary of the feasible Construct an α shape of region sampled data points Join points by a line in 2-d, by triangle in 3d Obtain polygonal representation of the feasible region PASI 2005 August 16-25 Iguazu Falls

Sampling Technique Implementation of this idea requires sampling of the feasible region Common sampling techniques sample the parameter space based on the distribution of the uncertain parameter Typically, the feasible region constitutes a small fraction of the entire parameter space Uniform sampling of entire parameter space can be expensive Require new technique which samples the feasible region with minimum total function evaluations PASI 2005 August 16-25 Iguazu Falls

Reformulation of Sampling Problem Obtain good sample of feasible region with less function evaluations Sampling problem framed as an optimization problem Objective function V feas (volume of max V feas θ the feasible region) evaluated by ( ) constructing the α -shape using the ≤ subject to: f 0 θ 1 sample points ( ) ≤ f 0 θ 2 � Improved only when the sampled ( ) ≤ f 0 n θ point is feasible θ : sampled parameter value Formulated optimization problem is solved using Genetic Algorithm GA has the inherent property of concentrating around good solutions PASI 2005 August 16-25 Iguazu Falls

Sampling Technique using GA Optimization variables 0 1 1 0 1 0 0 1 1 encoded as a string of bits v 1 v 2 v 3 Strings are appended to chromosome form a chromosome Solution procedure starts with a population of chromosome • Reproduction • Crossover Population of chromosome evolve through : • Mutation Reproduction: identifies good solutions in a population Cannot create makes multiple copies of the good new solution solution eliminates bad solution Crossover creates new solution by swapping portions of chromosomes Number of strings with similarities at certain string position is increased Schema PASI 2005 August 16-25 Iguazu Falls

Performance of the Sampling Technique Simple Random Sampling Sampling using GA 15 10 5 Feasible region bounded by nonconvex constraints θ 2 0 θ 2 -5 -10 -15 -10 -8 -6 -4 -2 0 2 4 θ 1 θ 1 Population size = 20 3000 function evaluations Solution evolved for 2000 ~ 1000 feasible points generations Random Sampling : 950 feasible point generation required 10000 function evaluation PASI 2005 August 16-25 Iguazu Falls

α − shape of the Sampled Data L n ~ 25 α value = n α value of 25 accurately captured the non-convex shape Higher value of α =1000 could not capture non-convexity PASI 2005 August 16-25 Iguazu Falls

Estimation of Feasible Region: α -shape Sample the feasible space by GA formulation 800 feasible points GA : 1800 function evaluations RS : 4000 function evaluations Construct α – shape with the sampled points Determine points forming the boundary of the feasible region Join boundary points with triangle for polygonal representation of feasible region PASI 2005 August 16-25 Iguazu Falls

Determination of Feasible Conditions Determination of conditions inside the feasible surface Point-in-polygon algorithm Ray casting algorithm � Draw semi-infinite ray from point of Jordan Curve Theorem concern � Determine number of times it intersects surface Point inside the surface number of intersections odd Point outside surface number of intersections even/0 Perform point-in-polygon test to check if particular parameter value lies inside feasible region PASI 2005 August 16-25 Iguazu Falls

Performance of α - shape α − shape covers a larger region than convex hull α - shape accurately captures the nonconvex shape The prediction of feasible region by α - shape can be improved by a better sampling technique Point-in-polygon check ~ 0.3 ms α shape could capture ~ 80 % of the feasible region PASI 2005 August 16-25 Iguazu Falls

Capturing Disjoint Feasible Region by α - shape Estimated feasible space Sampled feasible space Temperature Using α – shape it is possible Temperature to capture disjoint feasible regions Methane Methane Oxygen Oxygen PASI 2005 August 16-25 Iguazu Falls

Characterizing the Effects of Uncertainty PASI 2005 August 16-25 Iguazu Falls

Uncertainty Propagation Concentration vs. time plot Specific t=t1 Initial Conditions Not a Point but a Distribution Uncertainty Propagation Dynamic Parameter Variability Model P(A 1 ) t=t2 y =f(A1,A2,…, An) A 1 P(An) Deterministic Parameter Values An PASI 2005 August 16-25 Iguazu Falls

Stochastic Response Surface Method � The outputs are represented as a polynomial chaos expansion (Ghanem and Spanos, 1991) in terms of Hermite polynomials : n + ∑ 1 st order U a a ξ = i 1 0,1 i,1 i 1 = n n n 1 n − 2 nd order ∑ ∑ ∑∑ U a a a ( 2 1) a = + ξ + ξ − + ξξ 2 0,2 i,2 i ii,2 i ij,2 i j i 1 i 1 i 1 j i = = = > � The coefficients of these polynomials are determined through application of an efficient collocation scheme and regression � Direct evaluation of the output pdf’s characteristics (for example, for a single variable second order SRSM approximation) a a 2a Mean = Variance = 2 2 + 0,2 1,2 11,2 PASI 2005 August 16-25 Iguazu Falls

Stochastic Surface Response Method Input Distributions Output Distributions Model Select a set of standard random variables (srvs) and transform inputs Generate a set of in terms of srvs regression points Estimate the coefficients of the output approximation Outputs as a series in srvs with unknown coefficients • Two orders reduction in model runs required compared to Monte Carlo • Output uncertainty expressed as polynomial function of input uncertainty • Direct evaluation of the output pdf’s characteristics PASI 2005 August 16-25 Iguazu Falls

Case Study : Supercritical wet oxidation Reaction A i nom n E a /R UF i � Constant temperature (823K) high OH+H ↔ H2O 1.620E+14 0.00 75 3.16 pressure (246 Bar) oxidation of H 2 H2+OH ↔ H2O+H 1.023E+08 1.60 1660 1.26 and O 2 consisting of 19 reactions H+O2 ↔ HO2 1.481E+12 0.60 0 1.58 and 10 species HO2+HO2 ↔ H2O2+O2 1.866E+12 0.00 775 1.41 H2O2+OH ↔ H2O+HO2 7.826E+12 0.00 670 1.58 � Pre-exponential factors (A i ’s) H2O2+H ↔ HO2+H2 1.686E+12 0.00 1890 2.00 taken to be log-normal random H2O2 ↔ OH+OH 3.000E+14 0.00 24400 3.16 variables. Parameters obtained OH+HO2 ↔ H2O+O2 2.890E+13 0.00 -250 3.16 assuming: H+O2 ↔ OH+O 1.987E+14 0.00 8460 1.16 O+H2 ↔ OH+H 5.117E+04 2.67 3160 1.22 � Computed and literature values 2OH ↔ O+H2O 1.505E+09 1.14 50 1.22 of the multiplicative H2+M ↔ H+H+M 4.575E+19 -1.40 52530 3.0 uncertainty factors (UF) valid H+HO2 ↔ OH+OH 1.686E+14 0.00 440 1.35 for the reaction temperature H+HO2 ↔ H2+O2 4.274E+13 0.00 710 1.35 considered O+HO2 ↔ OH+O2 3.191E+13 0.00 0 1.49 H2O2+H ↔ H2O+OH 1.023E+13 0.00 1800 1.35 � 95% confidence limits provide O+H+M ↔ OH+M 4.711E+18 -1.00 0 10.0 upper and lower bounds. O+O+M ↔ O2+M 1.885E+13 0.00 -900 1.3 H2O2+O ↔ OH+HO2 6.622E+11 0.00 2000 1.35 †Phenix et. al. 1998 PASI 2005 August 16-25 Iguazu Falls

Uncertainty Propagation: Results (Balakrishnan S., P. Georgopoulos,I. Banerjee and M.G. Ierapetritou. AIChE J , 48 2875, 2002) Concentration profiles display time varying distributions � Number of model simulations required by SRSM is orders of � magnitude less than Monte Carlo (723 vs. 15,000) H2 mole fraction vs. time PASI 2005 August 16-25 Iguazu Falls

Design Considering Uncertainty PASI 2005 August 16-25 Iguazu Falls

Process Design Considering Market Demand • Potential Set of Units Flexible Production Plant • Uncertainty in Internal And External Conditions •Increase Plant Flexibility •Minimize Cost Flexibility Trade-off : OPTIMIZATION Cost Objective: Develop a systematic methodology to increase the plant at a minimum cost PASI 2005 August 16-25 Iguazu Falls

Background: Design Under Uncertainty Existing approaches to model uncertainty in design/planning problem depends on nature of model equations * Most models restricted by assumption of convexity * Most models restricted to a rather small number of uncertain parameters Require single model to describe uncertainty propagation irrespective of nature and complexity of the problem *Gal,T. Math.Prog.St. (1984), Jongen,H.T.,Weber,G.W. Ann. Op. Res. (1990), Pistikopoulos and coworkers, PASI 2005 August 16-25 Iguazu Falls

Background : Design under Uncertainty Deterministic Approach : description of uncertainty is provided by specific bounds, or finite number of fixed parameter values Grossmann, Halemane, AIChE (1982); Grossmann, Sargent AIChE (1987) Stochastic Approach : uncertainty described by probability distribution functions Pistikopoulos, Mazzuchi, Comput. Chem. Engg.(1990) Combined multiperiod/ stochastic formulation : combines parametric and stochastic programming approaches to deal with synthesis/planning problems Ierapetritou et al, Comput. Chem. Engg (1996), Hene et al, I&ECR (2002) PASI 2005 August 16-25 Iguazu Falls

Proposed Technique Tabulation technique to map input uncertainty to model output. Query Points ( d,q ) Perform Tabulate Systematic Model Prediction model runs results interpolation F(d,q) High Dimensional Model Representation (HDMR)* technique used to capture the variations of output with changes in the input * Rabitz,H. Alis,O. J. Math. Chem. 25,195(1999) PASI 2005 August 16-25 Iguazu Falls

Design with Parametric Uncertainty:Blackbox Models (Banerjee, I., and M.G. Ierapetritou. Ind. & Eng. Chem. Res , 41 , 6687, 2002) � No Assumptions Regarding System’s Model � Parametric Expression of the Optimal Solution Query Design Feasibility Points Analysis ( d, θ ) Feasible C(d, θ ) Look-up HDMR/ range of SRSM parameters Table Design Model Optimization Prediction Black Box C( θ ) F(d, θ) Model PASI 2005 August 16-25 Iguazu Falls

Feasibility Analysis Feasible region q 2 First order Second order prediction prediction q 1 Sampled points Sampled points for second order for first order approximation approximation PASI 2005 August 16-25 Iguazu Falls

High Dimensional Model Reduction n n = + + + + ∑ ∑ … … … g ( x , x , x ) f f ( x ) f ( x , x ) f ( x , x , , x ) … n n i i ij i j 1 2 0 1 , 2 , n 1 2 = ≤ < ≤ i 1 1 i j n f 0 constant f i (x i ) independent action of variable xi upon the output f ij (x i ,x j ) correlated impact of xi, xj upon the output ... f 1,2,…,n (x 1 ,x 2 ,…,x n ) residual correlated impact � Order of correlation of independent variables diminish rapidly � 2 nd order approximation commonly suffices � Application in complex kinetics modeling (i.e.,atmospheric chemistry, photochemical reaction modeling etc) Evaluation of first order expansion function requires n(s-1) model runs Evaluation of second order expansion requires n(s-1) 2 (n-1)/2 model runs PASI 2005 August 16-25 Iguazu Falls

Nonlinear Multiparametric Problem Optimization problem Feasibility problem Min z Min u Subject to : Subject to : -z- θ 1 + θ 2 2 /2+2 θ 3 3 +d 1 -3d 2 -8 ≤ 0 -z- θ 1 + θ 2 2 /2+2 θ 3 3 +d 1 -3d 2 -8 ≤ u -z- θ 1 /3- θ 2 - θ 3 /3+d 2 +8/3 ≤ 0 -z- θ 1 /3- θ 2 - θ 3 /3+d 2 +8/3 ≤ u z+ θ 1 2 - θ 2 -d 1 + θ 3 - 4 ≤ 0 z+ θ 1 2 - θ 2 -d 1 + θ 3 - 4 ≤ u Where: z is control variable. θ 1 ∈ [0 4]; θ 2 ∈ [0 4]; θ 3 ∈ [0 4] q 1 , q 2 , q 3 are uncertain parameters. d 1 ∈ [1 5]; d 2 ∈ [1 5] d 1 ,d 2 are design variables. PASI 2005 August 16-25 Iguazu Falls

Steps of Proposed algorithm:Feasibility Analysis Step 1: Fix the value of design variable. Determine the feasible region of operation . Constraints bounding Fixed value of θ 1 =2.56 HDMR prediction the feasible region No overprediction; 0.7% underprediction PASI 2005 August 16-25 Iguazu Falls

Steps of Proposed Algorithm: Optimization Problem Actual solution (surface) Step 2 : HDMR prediction Determine the variation (points) of optimal solution with uncertain parameters for the fixed value of design Estimation Error = 3.73 % PASI 2005 August 16-25 Iguazu Falls

Steps of Proposed Algorithm: Design Problem Step 3 : Perform the feasibility analysis and design optimization for different values of design variables. Estimation Error 7% HDMR prediction HDMR prediction Actual solution Uncertainty Uncertainty Error=5.5% Error=0.4% at the mean at the extreme PASI 2005 August 16-25 Iguazu Falls

Process Synthesis Problem Discretize θ in Identify different accordance with design configurations HDMR Solve MINLP at Construct Fix binary variables fixed values of θ look-up table at optimum configuration Update look-up Solve NLP/LP at table fixed values of θ over entire range PASI 2005 August 16-25 Iguazu Falls

Branch and Bound Procedure Binary variables : y1,y2 Relaxed problem Uncertain parameters : θ Y 1 =1 Root Node Y 1 =0 At each node solve NLP/LP Z 1 2 ,Z 2 2 , Z 1 Z 1 2 ,Z 2 2 , Z 1 at fixed values of θ ( θ 1 , θ 2 , θ 3 ) 3 3 θ 1 ,θ 2 ,θ 3 θ 1 ,θ 2 ,θ 3 Node 2 Node 1 Branching criteria: Y 2 =1 Y 2 =0 Choose a node having larger number of better optimal solutions Z 1 2 ,Z 2 2 , Z 1 Z 1 2 ,Z 2 2 , Z 1 3 3 θ 1 ,θ 2 ,θ 3 θ 1 ,θ 2 ,θ 3 Node 4 Node 3 Fathoming criteria: Compare solutions at all θ values. Fathom a node with respect to a particular θ value . PASI 2005 August 16-25 Iguazu Falls

Branch and Bound Procedure (Example) Branching Step : Relaxed problem Compare solutions of Node 1 and Node 2 z 1 2 > z 1 1 ; z 2 2 > z 2 1 ; z 3 2 < z 3 Y 1 =1 Root Node Y 1 =0 1 Selected node for branching: Node 2 Z 1 2 ,Z 2 2 , Z 3 Z 1 1 ,Z 2 1 , Z 3 2 Fathoming Step: 1 θ 1 ,θ 2 ,θ 3 θ 1 ,θ 2 ,θ 3 Compare Node 3 and Node 4 Node 2 Node 1 z 1 4 >z 1 3 ; z 2 4 < z 2 3 ; z 3 4 < z 3 3 Compare Node 3 and Node 1 Y 2 =0 Y 2 =1 z 1 3 > z 1 1 ; z 2 3 > z 2 1 ; z 3 3 < z 3 1 Z 1 4 ,Z 2 4 , Z 2 Z 1 3 , Z 2 3 , Z 3 4 3 Optimal solution θ 1 ,θ 2 ,θ 3 θ 1 ,θ 2 ,θ 3 At θ 1 : z 1 4 [1,1] ; At θ 2 : z 2 3 [1,0] Node 4 Node 3 Fathom Node 1 wrt θ 1 ,θ 2 Branch on Node 1 only for θ 3 PASI 2005 August 16-25 Iguazu Falls

Example Problem with Single Uncertain Parameter B f I 1 Binary variables : 8 O 1 Process Uncertain parameter : 1 1 O 2 I 2 O 4 Process Process I 4 2 4 O 3 I 3 I 7 Process Process ≥ θ P D ( ) 3 7 1 1 I 5 O 5 Process I 8 Process 5 ≥ θ P D ( ) I 6 8 2 2 O 6 Process 6 * Acevedo and Pistikopoulos 1996 PASI 2005 August 16-25 Iguazu Falls

Application of the Proposed Approach Solution of MINLP at discrete θ points * * Step 1: Uncertain range of θ discretized. * * MINLP solved at each discrete * θ value. * * Discretize θ Optimal binary solutions are noted [0,1,0,1,0,1,1,1] [0,1,0,1,1,1,1,1] PASI 2005 August 16-25 Iguazu Falls

Application of the Proposed Approach [0,1,0,1,0,1,1,1] Step 2: Fix binary variable at [0,1,0,1,1,1,1,1] optimal combination. Solve NLP/LP at different θ values over entire range of θ Predict variation of optimal solution for each binary combination over entire range of θ PASI 2005 August 16-25 Iguazu Falls

Application of the Proposed Approach Analyze predicted variation of optimal solution to determine optimal binary configuration and optimal solution Actual solution HDMR prediction Estimation Error = 1.7% [0,1,0,1,0,1,1,1][0,1,0,1,1,1,1] PASI 2005 August 16-25 Iguazu Falls

Error Analysis of Process Synthesis Problem Problem 1(linear): 3 Binary variables ; 1 uncertain parameter Error % 2 2 Problem 2 (nonlinear): 8 Binary variables; 1 Uncertain parameters 1.5 1.5 1 1 Problem 3 (nonlinear): 0.5 0.5 8 Binary variables; 2 Uncertain parameters 0 0 Problem 4(linear): Prob 1 Prob 1 Prob 2 Prob 2 Prob 3 Prob 3 Prob4 Prob4 Prob 5 Prob 5 2 Binary variables; 3 Uncertain parameters Problem 5 (nonlinear): 6 Binary variables; 3 Uncertain parameters PASI 2005 August 16-25 Iguazu Falls

Design Optimization Integrating Market Data Integration of Data Analysis,and Feasibility Quantification at the Process Design Traditionally: Design for the Nominal Point Performance • Limited Flexibility Cluster 3 Cluster 4 • Poor Performance Away Product 2 From the Nominal Cluster 1 Area of Interest Cluster 2 Based on Market Data Product 1 • Increased Plant Flexibility • Better Performance Within the Whole Range of Interest • Larger Profitability Due to the Economy of Scale PASI 2005 August 16-25 Iguazu Falls

Motivation Example Product P2 Product P1 � Produce P 1 and P 2 from A,B,C � Given Demand Data for P 1 and P 2 � MINLP Optimization PASI 2005 August 16-25 Iguazu Falls

Flexibility Plot for Customized Design Development Demand = (9,17) Design configuration = (1,0,0,1,0,0,1,1) (Processes 1, 4, 7, 8) Flexibility Index = 0.21 PASI 2005 August 16-25 Iguazu Falls

Limitations -1 � Underestimation of the Feasible Region � Demand point (6,16) Feasible No New Design Needed PASI 2005 August 16-25 Iguazu Falls

Limitations -2 � Customized Design Development � Demand point (13,20) New Design Required PASI 2005 August 16-25 Iguazu Falls

Moving the System Boundaries Supply chain information Given: - Ranges of product quality Subcooler LPC - Ranges of throughput MAC BAC PPU - Design alternatives MA Reboiler / PL Cond. GOX HPA LOX MPN WN G Given: TA IL - Product specifications Turbine/ of different clients HPC Generator MHx - Client expectations JT Valve RL IC Pump Find: the equipment set that minimizes cost and Find: the set of designs that optimally cover the maximizes flexibility whole space PASI 2005 August 16-25 Iguazu Falls

Required Tools � Accurate description of the feasible space of a process � Data Clustering technique to cluster the demand data into closely packed groups � Development of a unified data analysisprocess optimization framework PASI 2005 August 16-25 Iguazu Falls

Data Analysis (Clustering) � Partition a data set of multi-dimensional vectors into clusters such that patterns within each cluster are more “similar” to each other than to patterns in other clusters. � Quality of Clustering depends on both the similarity measure used by the algorithm and its implementation. PASI 2005 August 16-25 Iguazu Falls

K-Medoid Clustering: PAM � Find representative objects, called medoids, in clusters � PAM (Partitioning Around Medoids) starts from an initial set of medoids and iteratively replaces one of the medoids by one of the non-medoids if it improves the total distance of the resulting clustering. Kaufman and Rousseeuw, Finding Groups In Data (1990) PASI 2005 August 16-25 Iguazu Falls

Design Optimization Integrating Market data Integration with Design Optimization Design/Synthesis with CONVENTIONAL Data Analysis - APPROACH fixed degree of feasibility Clustering (Multi-period Model) Evaluate Cluster For each Centers Simulation for the specific design Cluster center (Black-box models) x 2 Feasibility Check G F Cluster 3 Cluster 1 ΝΟ C DE Check if the B A designs cover Cluster 2 Evaluate Designs’ x the whole 1 Feasible Regions space YES STOP PASI 2005 August 16-25 Iguazu Falls

Air Separation Plant Superstructure Heat Exchange AIR {2 : 0} Distillation Refrigeration {4 : 2} {4 : 4} Compression {4 : 2} PRODUCTS Oxygen/Nitrogen Gaseous/Liquid {main options/suboptions per main option} PASI 2005 August 16-25 Iguazu Falls

Air Separation Case Study:Sample Demand Source Factors for Consumption of Oxygen (tons per ton of product) Product Consumption Ethylene Oxide 1.01 Propylene Oxide 1.26 Sample Plant Capacities (million lb/yr) Ethylene Oxide 150 300 600 Propylene Oxide 200 400 800 PASI 2005 August 16-25 Iguazu Falls

Demand Data Each point represents a different potential customer Pgo (Atm) Fgo (KgMol/Hr) PASI 2005 August 16-25 Iguazu Falls

Results – Iteration 1 (3 Clusters) Simplicial (411.2,50) Approximation (1696.2, 41) Feasible (514,10) Points PASI 2005 August 16-25 Iguazu Falls

Optimality Stage ( 1696.2,35 ) (2056, 41.5) (215.9,10) (1455.4, 37.25) (1865, 40) (343.7,17.4) Lower Cost Lower Cost SAME Similar Flexibility Lower Flexibility Accepted Rejected PASI 2005 August 16-25 Iguazu Falls

Final Design Portfolio Simplicial Approximation Multiple Clients D1 D2 Accurate descriptio A A A of the operability A2 A1 A A3 C1 A1 C3 boundaries Quality D2 y C1 A1 t D1 i D1 l a D1 u Q Trade-offs between Hx cost and flexibility Rf Feasible Points A Dist A A Comp C1 Client Demand B1 Sensitivity to differ A1 B1 D2 increased cost A1 units but higher flexibility Throughput Throughput PASI 2005 August 16-25 Iguazu Falls

Sensitivity to Different Units Change Distillation Option to Increase Throughput Change Refrigeration Option to Increase Quality PASI 2005 August 16-25 Iguazu Falls

Introduction of New Technology Changing Compressor to Larger Capacity DOES NOT Increase Plant Flexibility PASI 2005 August 16-25 Iguazu Falls

Staging Design - Spare Units Initial New Design Design Hx A A Rf A2 A2 Dc C1 C1+A1 Co C1 C1+A1 * 4 1 Quality Inst. Cost 13E6 20.18E6 Oper. Cost 11.E6 12.8E6 Incremental Cost: Distillation: 2.1E6 Compressor: 5.07E6 1233.6 1400 Throughput PASI 2005 August 16-25 Iguazu Falls

Relevance and Importance Manufacturer : Modular-based designs are substantially cheaper than customized designs and can satisfy larger range of demands Customer : Greater flexibility in decision making at design stage as different design alternatives can be considered based on expected demand and economic feasibility PASI 2005 August 16-25 Iguazu Falls

Multi-Period Robust Design Optimization Σ Min Capital and Operating Cost i + Cost { n λ Σ (variance in operating cost) 2 ′ i i i i 2 λ − w (C w C ) ∑ ∑ Robustness i ′ = i 1 i + m n β Σ (Un-met Demands) 2 i i 2 β j 1i 1 w (z ) ∑ ∑ j i = = Subject to: i i i θ = ∀ h (d,x ,u , ) 0 i � Process Constraints i i i θ ≤ ∀ g (d,x ,u , ) 0 i � Un-met Demand Constraints i i i + ≥ ∀ F z F i, j Prod _ j j demand _ j i = periods, λ , β = robustness parameters * Mulvey et al., Robust Optimization of Large Scale Systems, Oper. Res., 1995 PASI 2005 August 16-25 Iguazu Falls

Flexible Module-Based Design Generation Data Analysis - Robust Design Clustering NLP Optimization for (Facility Location) Formulation For each each Cluster Cluster G F Cluster 3 y Cluster 1 D C E Evaluate Designs’ Evaluate Designs’ A B YES Cluster 2 Increase the Feasible Regions Feasible Regions x number of Using Simplicial Using Simplicial Any new clusters by one Approximation Approximation Designs Check if the NO NO designs cover Enforce feasibility Cost the whole For the design Flexibility space YES STOP PASI 2005 August 16-25 Iguazu Falls

Illustrating Case Study F Prod A B k1 F ao C ao x a x b x c k2 x d x e k3 C D k4 k5 α E β � Given rate-constants and demand � Determine CSTR Volume, Input F ao and C ao that minimize overall cost PASI 2005 August 16-25 Iguazu Falls

Demand Data Product E (Mol/Hr) Product B (Mol/Hr) PASI 2005 August 16-25 Iguazu Falls

Results – Iteration 1 (2 Clusters) •V = 49.17 •Capital Cost = $67006.7 •V = 35.83 •Capital Cost = $36924.6 PASI 2005 August 16-25 Iguazu Falls

Iteration 4 (5 Clusters) 11 designs to •V = 44.4 cover demand space Product Ε (Mol/Hr) •Capital Cost = $ 56276.3 •V = 40.1 • V = 26.1 •Capital Cost • Capital Cost = $47821. = $ 24791.9 Product B (Mol/Hr) PASI 2005 August 16-25 Iguazu Falls

Capital Cost vs Feasibility m n j 1 i 1 w(z ) i i i β i i 2 F z F + ≥ ∑∑ j Pr od _ j demand _ j j = = Increasing β drives the objective towards larger Increasing β Infeasibilities designs with higher feasibility Capital Cost PASI 2005 August 16-25 Iguazu Falls

Capital Cost vs Robustness n w (C w C ) λ i i − i i 2 + ′ ∑ ∑ i 1 i ′ = Increasing λ drives the objective towards Increasing λ lower operating cost Operating Cost systems at the expense of the fixed cost Capital Cost PASI 2005 August 16-25 Iguazu Falls

Optimization of Noisy Systems PASI 2005 August 16-25 Iguazu Falls

Optimization of Noisy Functions Inputs x • Programming Model*: Intermediate Outputs 1100 -210 min Z = c T y + F(x, ε ) -220 1000 -230 900 s.t. h(X) = 0 -240 800 -250 g(X) + My ≤ 0 700 -260 600 ∈ ∈ -270 x X, y Y 500 -280 400 -290 300 -300 6 6.5 7 7.5 8 8.5 9 9.5 10 6 6.5 7 7.5 8 8.5 9 9.5 10 •X is a vector of continuous variables, (P, T, Flowrates) Final Outputs •Y is a vector of binary variables, (existence of a 45 particular stream or unit) 40 •The uncertainty ε can propagate or dampen as the 35 process moves forward 30 •Optimality conditions cannot be defined at optima 25 •Conventional algorithms may become trapped in artificial 20 local optima or even fail completely 15 6 6.5 7 7.5 8 8.5 9 9.5 10 *Biegler et. al., Systematic Methods of Chemical Process Design, (1997), 513 PASI 2005 August 16-25 Iguazu Falls

Existing Work DIRECT “DIvided RECTangles in action” (Jones et. al., 1993) x 2 85 x 1 80 •Splits feasible region into hyper-rectangles and samples at center points � global search 75 •Scatter plot created to discover which sample points lie below a prescribed improvement in the objective � local 70 search •If the best point is unsatisfactory, smaller hyper- 65 rectangles are inscribed inside region and sampling at the centers of these new regions continues 60 0.2 0.4 0.6 0.8 1 •Slow to converge, especially if the optimum is along a boundary PASI 2005 August 16-25 Iguazu Falls

Existing Work Multilevel Coordinate Search (Huyer & Neumaier, 1998) •Avoids slow convergence of DIRECT by sampling at boundary points •Newton-based methods/SQP minimize interpolating polynomials to obtain new Irregularly split regions regions for sample points allow larger area to be sampled during local search Implicit Filtering (Choi & Kelley, 1999, Gilmore & Kelley, 1994) σ 2 σ 1 σ 1 •Applies Newton-based methods with step sizes proportional to high-frequency noise, “filtering”, σ 3 or “stepping over” low-frequency noise •Successively decreases the step size as optimum is approached PASI 2005 August 16-25 Iguazu Falls

Existing Work Differential Evolution (Storn & Price, 1995) X i,G =[X i1,G X i2,G …X iD,G ] Gen. G Choose multiplier F ∈ [0,2] Select crossover index CR ∈ [0,1] For i = 1…N: Select Γ (i) = integer ∈ 1…D (ensures at least X 1,G one element from V i,G mixes with X i,G Randomly select integers X 2,G β (j) ∈ U[0,1], j = 1…D r 1 ,r 2 ,r 3 ∈ 1…N, r a ≠ r b ≠ i ( ) ( ) v if (j) CR or j (i) β ≤ = Γ ⎧ ⎫ : ji, G V i,G = X r1,G + F(X r2,G – X r3,G ) X = ⎨ ⎬ ( ) ( ) ⎭ ji, G 1 x if (j) CR and j (i) + β > ≠ Γ : ⎩ ji, G X N,G X 1,G V 1,G X 1,G+1 X 11 V 11 V 11 X 12 V 12 X 12 X 13 V 13 V 13 X 14 V 14 X 14 X 15 V 15 X 15 X 16 V 16 V 16 PASI 2005 August 16-25 Iguazu Falls