Algorithm: Multi-cut L-shaped Method Problem Sizes Two-stage Stochastic Programming Model Small Deterministic Problem Model 10 scenarios 100 scenarios 1,000 scenarios # of Constraints 1,369 13,080 130,170 1,301,070 # of Variables 3,937 37,248 370,338 3,701,240 # of Non-zeros 8,910 85,451 850,271 8,498,429 Two-stage Stochastic Programming Model Deterministic Full Problem Model 10 scenarios 100 scenarios 1,000 scenarios # of Constraints 6,373 61,284 610,374 6,101,280 # of Variables 19,225 182,496 1,815,816 18,149,077 # of Non-zeros 41,899 402,267 4,004,697 40,028,872 Note: Problems with red statistical data are not able to be solved by DWS Page 25
Algorithm: Multi-cut L-shaped Method Two-stage SP Model Master problem Scenario sub- problems y 1 y Scenario 2 sub-problems Master y problem y x S Page 26
Algorithm: Multi-cut L-shaped Method Standard L-shaped Method Solve master problem to get a lower bound (LB) cuts Add cut Solve the subproblem to get an upper bound (UB) Yes No STOP UB – LB < Tol ? Page 27
Algorithm: Multi-cut L-shaped Method Expected Recourse Function The expected recourse function Q(x) is convex and piecewise linear Each optimality cut supports Q(x) from below Page 28
Algorithm: Multi-cut L-shaped Method Multi-cut L-shaped Method Solve master problem to get a lower bound (LB) cuts Add cut Solve the subproblem to get an upper bound (UB) Yes No UB – LB < Tol ? STOP Page 29
Algorithm: Multi-cut L-shaped Method Example 220 Standard L-Shaped Upper_bound Standard L-Shaped Lower_bound 210 Multi-cut L-Shaped Upper_bound Multi-cut L-Shaped Lower_bound 200 Cost ($MM) 190 180 170 160 150 140 1 21 41 61 81 101 121 141 161 181 Iterations Page 30
Background Optimal Design of Responsive Process Supply Chains You, Grossmann (2008) Objective : design supply chains under responsive and economic criteria with consideration of inventory management and demand uncertainty
Background Problem Statement W here? W hat? W hen? DCs Customers Suppliers Plants Production Network Network Structure Costs and prices Operational Plan Production and transportation time Max: Net present value Demand information Max: Responsiveness Production Schedule
Example Production Network of Polystyrene Resins Three types of plants: Plant I : Ethylene + Benzene Styrene (1 products) Plant II : Styrene Solid Polystyrene (SPS) (3 products) Plant III : Styrene Expandable Polystyrene (EPS) (2 products) Basic Production Network Multi Product Single Product Multi Product Source: Data Courtesy Nova Chemical Inc. http://www.novachem.com/
Distribution Center Customer Location Possible Plant Site Supplier Location Location Map Example
Example Potential Network Superstructure WA IL NV Plant Site MI SPS CA Ethylene II Styrene I III TX Benzene AZ TX EPS SPS OK II TX Ethylene GA I III Styrene Styrene NC EPS Benzene Plant Site TX Plant Site CA MS FL Ethylene PA Styrene I LA OH Benzene III EPS Plant Site LA IA MA AL MN Suppliers Plant Sites Distribution Centers Customers
Model & Algorithm Responsiveness - Lead Time � Lead Time: The time of a supply chain network to respond to customer demands and preferences in the worst case Responsiveness Lead Time Lead Time is a measure of responsiveness in SCs
Model & Algorithm Lead Time for A Linear Supply Chain • A supply chain network = ∑ Linear supply chains � Assume information transfer instantaneously Suppliers Plants Distribution Centers Customers Information
Model & Algorithm Lead Time for Deterministic Demand Transportation Delay = Transportation Time Production Delay = Residence Time (single product plants) Lead Time
Example Lead Time of SCN • Lead time of a supply chain network (deterministic demand) � The longest lead time for all the paths in the network ( worst case ) � Example: A simple SC with all process are dedicated Path 1 7.7 days Path 2 8.0 days For Path 1: 2 + 1.5 + 0.5 + 1.2 + 1.8+ 0.7 = 7.7 days For Path 2: 2 + 1.5 + 0.2 + 2.6 + 1.2 + 0.5 = 8.0 days Lead Time = max {7.7, 8.0} = 8.0 days
Inventory (Safety Stock) Lead Time under Demand Uncertainty Model & Algorithm
Example Expected Lead Time of SCN • Expected Lead time of a supply chain network (uncertain demand) � The longest expected lead time for all the paths in the network ( worst case ) � Example: A simple SC with all process are dedicated P 1 =20% Path 1 2.1 days P 2 =20% Path 2 2.0 days P For Path 1: (2 + 1.5 + 0.5 + 1.2 + 1.8) × 20% + 0.7 = 2.1 days For Path 2: (2 + 1.5 + 0.2 + 2.6 + 1.2) × 20% + 0.5 = 2.0 days Expected Lead Time = max {2.1, 2.0} = 2.1 days Safety Stock
Model & Algorithm Stock-out Probability (P) • Chance constraint for stockout probability � Integrate lead time, inventory management, demand uncertainty Chance constraint Generalized Disjunctive Programming Target Demand Target Demand MINLP d L d M d U Safety Stock Safety Stock
Model & Algorithm Objective Functions • Responsiveness � Measured by expected lead time • Economics � Measured by net present value (NPV) Sales income Purchase cost Operating cost Transport cost Investment cost Inventory cost
Model & Algorithm Bi-criterion Multiperiod MINLP Formulation NPV Choose Discrete (0-1), continuous variables Pareto Curve • Objective Function: � Max : Net Present Value Bi-criterion � Min : Expected Lead time • Constraints: � Expected Lead Time Network structure constraints Suppliers – plant sites Relationship � Cyclic scheduling constraints Plant sites – Distribution Center Assignment constraint Input and output relationship of a plant Sequence constraint Distribution Center – Customers Demand constraint Cost constraint � Operation planning constraints Production constraint Production constraint Cost constraint � Capacity constraint Probabilistic constraints Mass balance constraint Chance constraint for stock out Demand constraint (reformulations) Upper bound constraint
Model & Algorithm Procedure for Pareto Optimal Curve NPV Maximize: NPV – ε · Lead Time ( ε = 0.001 ) Pareto Curve Highest NPV Impossible! Lowest NPV Minimize: Lead Time Lead Time Longest Lead Time Shortest Lead Time
Example Case Study WA IL NV Plant Site MI SPS CA Ethylene II Styrene I III TX Benzene AZ TX EPS SPS OK II TX Ethylene GA I III Styrene Styrene NC EPS Benzene Plant Site TX Plant Site CA MS FL Ethylene PA Styrene I LA OH Benzene III EPS Plant Site LA IA MA AL Possible Plant Site Possible Plant Site Possible Plant Site Distribution Center Distribution Center Distribution Center MN Supplier Location Supplier Location Supplier Location Customer Location Customer Location Customer Location Suppliers Plant Sites Distribution Centers Customers • Problem Size: • Solution Time: � Solver: GAMS/BARON � # of Discrete Variables: 215 � Direct Solution: > 2 weeks � # of Continuous Variables: 8126 � Proposed Algorithm: ~ 4 hours � # of Constraints: 14617
Example Pareto Curves – with and without safety stock 750 Best Choice 700 650 600 NPV (M$) 550 More Responsive 500 450 400 with safety stock 350 without safety stock 300 1.5 2 2.5 3 3.5 4 4.5 5 5.5 Expected Lead Time (day)
Example Safety Stock Levels - Expected Lead Time 200 EPS in DC2 SPS in DC2 EPS in DC1 SPS in DC1 150 Safety Stock (10^4 T) More inventory, 100 more responsive 50 Responsiveness 0 1.51 2.17 2.83 3.48 4.14 4.8 Expceted Lead Time (day)
Optimal Network Structure 750 700 Pareto Curve 650 (B) (C) NPV (M$) 600 550 500 (A) 450 400 1.5 2 2.5 3 3.5 4 4.5 5 5.5 Expected Lead Time (day)
Example Optimal Network Structure – (A) Shortest Expected Lead Time = 1.5 day NPV = $489.39 MM WA IL NV Plant Site MI SPS CA II Styrene III TX AZ TX EPS SPS OK II TX Ethylene GA I III Styrene Styrene OH EPS Benzene Plant Site TX Plant Site CA MS FL Ethylene PA Styrene I LA Benzene III NC EPS Plant Site LA IA MA AL MN Suppliers Plant Sites Distribution Centers Customers
Example Optimal Network Structure – (B) Expected Lead Time = 2.96 days NPV = $644.46 MM WA IL NV Plant Site MI SPS CA II Styrene III TX AZ TX EPS SPS OK II TX Ethylene GA I III Styrene Styrene OH EPS Benzene Plant Site TX Plant Site CA MS FL Ethylene PA Styrene I LA Benzene III NC EPS Plant Site LA IA MA AL MN Suppliers Plant Sites Distribution Centers Customers
Example Optimal Network Structure – (C) Longest Expected Lead Time = 5.0 day NPV = $690 MM WA IL NV Plant Site MI SPS CA II Styrene III TX AZ TX EPS SPS OK II TX Ethylene GA I III Styrene Styrene OH EPS Benzene Plant Site TX Plant Site CA MS FL Ethylene PA Styrene I LA Benzene III NC EPS Plant Site LA IA MA AL MN Suppliers Plant Sites Distribution Centers Customers
Simultaneous Tactical Planning and Enterprise Optimization Production Scheduling Goal: Improve the asset utilization of geographically distributed assets and reduce cost to serve by improving enterprise wide tactical production planning. Production Plant Customer Multi-scale optimization: temporal and spatial integration 53
Production Planning for Parallel Batch Reactors Erdirik, Grossmann (2006) A Reaction 1 F1 STORAGE Materials: INTERMEDIATE STORAGE � Raw materials, Intermediates, Finished products F2 � Unit ratios (lbs of needed material per lb of material B STORAGE Reaction 2 produced) F3 C Reaction 3 Production Site: STORAGE F4 � Reactors: � Products it can produce � Batch sizes for each product � Batch process time for each product (hr) � Operating costs ($/hr) for each material � Sequence dependent change-over times /costs � => Lost capacity (hrs per transition for each material pair) � Time the reactor is available during a given month (hrs) due date due date due date due date due date due date Customers: � Monthly forecasted demands for desired products � Price paid for each product week 1 week 1 week 2 week 2 week t week t 54
Problem Statement Problem Statement DETERMINE THE PRODUCTION PLAN: DETERMINE THE PRODUCTION PLAN: � Production quantities � Inventory levels � Number of batches of each product � Assignments of products to available processing equipment � Sequence of production in each processing equipment OBJECTIVE: OBJECTIVE: To Maximize Profit . Profit = Sales – Costs Costs =Operating Costs + Inventory Costs +Transition Costs 55
Approaches to Planning and Scheduling Decomposition Simultaneous Planning and Scheduling Decomposition Simultaneous Planning and Scheduling Sequential Hierarchical Approach Detailed scheduling over the entire horizon Planning Planning months, years days, weeks Scheduling Scheduling Challenges: Challenges: � Different models / different time � Very Large Scale Problem scales � Solution times quickly intractable � Mismatches between the levels Goal: Planning model that integrates major aspects of scheduling 56
57 Results for Detailed MILP Scheduling Model: 4 reactors,6 products (1 week)
MILP Detailed Scheduling Model MILP Detailed Scheduling Model Objective Function: ( ) ∑∑ ∑∑∑∑ ∑∑ ∑∑∑∑∑ � = ⋅ − ⋅ − ⋅ + + − ⋅ + + ˆ ofit CP S COP XB CINV INV INVFIN INVINT CTRA Z Z Z Pr ( ) i t i t i t i m l t i t i t i t i t i k i k m l t i k m l t i k m l t , , , , , , , , , , , , , , , , , , , , , , , i t i m l t i t i k m l t Assignment constraints and Processing times: ∑ ≤ ∈ ∀ ∈ ∩ ∀ ∀ W i IM m l L m L t m t 1 ( ), ( ( ) ( )), , i m l t , , , i ∑ ∑ ≥ ∀ ∈ ∩ ≠ ∀ ∀ W W l L m L t l N t m t ( ( ) ( )), ( ), , i m l t i m l t , , , , , , i IM m ∈ i IM m ∈ ( ) ( ) = ⋅ ∀ ∈ ∀ ∈ ∩ ∀ ∀ PT BT W i IM m l L m L t m t ( ), ( ( ) ( )), , i m l t i m i m l t , , , , , , , = ⋅ ⋅ ∀ ∈ ∀ ∈ ∩ ∀ ∀ X R BT W i IM m l L m L t m t ( ), ( ( ) ( )), , i m l t i m i m i m l t , , , , , , , , Detailed timing constraints and sequence dependent change : ≥ + − Z W W 1 ∀ ∈ ∀ ∈ ≠ ∀ ∈ ∩ ≠ ∀ ∀ i IM m i IM m i i l L m L t l Nt m t ( ), ' ( ), ' , ( ( ) ( )), , , + i i m l t i m l t i m l t , ', , , , , , ', , 1, ∑∑ = τ ⋅ TR Z ∀ ∈ ∩ ≠ ∀ ∀ l L m L t l Nt m t ( ( ) ( )), , , m l t i i i i m l t , , , ' , ', , , i i ' � ≥ + − Z W W 1 ∀ ∈ ≠ ∀ ∈ ∩ ∀ ∀ ≠ i i IM m i i l L m L t m t t Ht , ' ( ), ' , ( ( ) ( )), , , + i i m l t i m l t i m t , ', , , , , , ', ,'1', 1 ∑ ∑∑ ∑∑ = + + τ ⋅ + − + τ ⋅ ˆ Te Ts PT Z TRT TX Z ( ) ∀ m l t , , m l t m l t i m l t i k i k m l t m l t m l t i i i i m l t , , , , , , , , , , , , , , , , , ' , ', , , i i k i i ' = + TRT TRT TRT 1 2 ∀ m l t , , m l t m l t m l t , , , , , , = TX TRT 1 ∀ m l t , , m l t m l t , , , , ≤ ⋅ TRT UPPER Y 1 ∀ m l t , , + m l t m l t , , , 1, ≤ ⋅ − TRT UPPER Y 2 (1 ) ∀ m l t , , + m l t m l t , , , 1,
MILP Detailed Scheduling Model MILP Detailed Scheduling Model Mass and Inventory Balances: = + FIN TRA X INVP INVINT ∈ ∀ ∈ ∩ ∀ ∀ i IFINT l L m L t m t , ( ( ) ( )), , i i m l t i m l t i m l t , , , , , , , , , ∑ ∑ ∑ = + + TRA INT INVINT INVP AA AA ∈ ∀ ∈ ∩ ∀ ∈ ∀ ∀ ≠ ∀ i IFINT l L m L t l L t m m m t , ( ( ) ( )), ' ( ), , ' , i i m l t i m l t i m l m l t i m l m l t , , , , , , , , , , ', , , , ', ', > ∈ ≠ ∈ l l l L m m m l L m ' , ' ( ) ' ' ( ') ∑ ∑ ∑ = + + + FIN INT X INVP INVP AA AA ∈ ∀ ∈ ∩ ∀ ∈ ∀ ∀ ≠ ∀ i IFINT l L m L t l L t m m m t , ( ( ) ( )), ' ( ), , ' , i i m l t i m l t i m l t i m l m l t i m l m l t , , , , , , , , , , , , , ', , , , ', ', > ∈ ≠ ∈ l l l L m m m l L m ' , ' ( ) ' ' ( ') ∑ + ≤ + ⋅ − Ts PT Ts BigW YY (1 ) ∀ ∈ ∩ ∀ ∈ ∩ ∀ ∀ ≠ ∀ l L m L t l L m L t m m m t ( ( ) ( )), ' ( ( ') ( )), , ' , m l t i m l t m l t t l l m m t , , , , , ', ', , ', , ', i ∑ ≤ + + ⋅ Ts Ts PT BigW YY ( ) ∀ ∈ ∩ ∀ ∈ ∩ ∀ ∀ ≠ ∀ l L m L t l L m L t m m m t ( ( ) ( )), ' ( ( ') ( )), , ' , m l t m l t i m l t t l l m m t ', ', , , , , , , ', , ', i ≤ ⋅ AA UBOUND YY ( ) ∈ ∀ ∈ ∩ ∀ ∈ ∩ ∀ ∀ ≠ ∀ i IFINT l L m L t l L m L t m m m t , ( ( ) ( )), ' ( ( ') ( )), , ' , i i m l m l t i l l m m t , , , ', ', , ', , ', ∑ ∑ ≤ ⋅ AA UBOUND W ( ) ∈ ∀ ∈ ∩ ∀ ∀ i IFINT l L m L t m t , ( ( ) ( )), , i i m l m l t i i m l t , , , ', ', , , , ≠ ∈ ∩ m m l L m L t ' ' ( ( ') ( )) ∑ ≤ ⋅ AA UBOUND W ( ) ∈ ∀ ∈ ∩ ∀ ∀ i IFINT l L m L t m t , ( ( ) ( )), , i i m l m l t i i m l t , , , ', ', , , , ∈ ∩ l L m L t ' ( ( ) ( )) ∑ ≤ ⋅ AA UPBOUND W ( ) i ∈ IFINT ∀ ∈ l L m ∩ L t ∀ l ∈ L m ∩ L t ∀ m ∀ m ≠ m ∀ t , ( ( ) ( )), ' ( ( ') ( )), , ' , i m l m l t i m i m l m l t i , , , ', ', , ', , , ', ', ∈ i ENDINT i i ' ( ', ) ∑ ∑ ∑ ∑ INT + INT = + INT INV INVP INVC INV ∈ ∀ i IFINT t , i i t − i m l t i m l t i t , 1 , , , , , , , ∈ ∩ ∈ ∩ m l L m L t m l L m L t ( ( ) ( )) ( ( ) ( )) ∑ ∑ + = + FIN FIN FIN INV INVP S INV ∈ ∀ i IFINT t , − i i t i m l t i t i t , 1 , , , , , ∈ ∩ m l L m L t ( ( ) ( )) ∑ ∑ ∑ ∑ = α ⋅ + + X INVC AA AA ( ) ∈ ∀ ∈ ∩ ∀ ∀ i IE l L m L t m t , ( ( ) ( )), , i m l t i i i m l t i m l m l t i m l m l t i , , , , ' ', , , , , '', , , , ', ', , , < ≠ ∈ ∩ i l l m m l L m L t ' '' ' ' ( ( ') ( )) l ∈ L m ∩ L t ' ( ( ) ( )) ∑ ∑ FIN + = + FIN INV X S INV ∈ ∪ ∀ i IE IF t ( ), − i i i t i m l t i t i t , 1 , , , , , ∈ ∩ m l L m L t ( ( ) ( ))
Proposed MILP Planning Models Replace the detailed timing constraints by: Model A. (Relaxed Planning Model) � Constraints that underestimate the sequence dependent changeover times � Weak upper bounds (Optimistic Profit) Model B. (Detailed Planning Model) � Sequencing constraints for accounting for transitions rigorously (Traveling salesman constraints) � Tight upper bounds (Realistic estimate Profit) 60
Generic Form of Proposed MILP Planning Models � Mass Balances on State Nodes P j t , SET ∈ ∈ I I P F P F S P S j t , FP FP i m t i m t , , , , I ∈ ∈ SETO I C S j � Time Balance Constraints on Equipment Reactor R Reactor R Changeover Changeover Product 1 time Product 2 time …………………….. Available time for R � Objective Function 61
Key Variables for Model YP :the assignment of products to units at each time period i m t , , NB :number of each batches of each product on each unit at each period imt FP :amount of material processed by each task imt Products: A, B, C, D, E, F Reactor 1 or Reactor 2 or Reactor 3 A A A B D B E B R1 F B F B B C D R2 F B F E C D R3 T3 T1 T2 = = YP YP 1 1 A reactor time , 1, 1 B reactor time , 2, 2 = = NB NB 3 2 A reactor time B reactor time , 1, 1 , 2, 2 62
Proposed Model B (Detailed Planning) Sequence dependent changeovers (traveling salesman constraints): � Changeovers within each time period: 1. Generate a cyclic schedule where total transition time is minimized. KEY VARIABLE: ZP ' :becomes 1 if product i is after product i’ on unit m at time period t, zero otherwise ii mt P1 P1, P2, P3, P4, P5 ZP P1, P2, M, T = 1 P5 P4 P2 ZP P2, P3, M, T = 1 P4 P3 2. Break the cycle at the pair with the maximum transition time to obtain the sequence. KEY VARIABLE: ZZP :becomes 1 if the link between products i and i’ is to be broken, zero otherwise ii mt ' P1 => P4 → P5 → P1 → P2 → P3 P5 P2 P4 P3 ? ZZP P4, P3, M, T 63
MILP Model According to the location of the link to be broken: P2, P3, P4, P5, P1 P2, P3, P4, P5, P1 ZZP P1, P2, M, T = 1 ZZP P1, P2, M, T = 1 P1 P1 P3, P4, P5, P1, P2 P3, P4, P5, P1, P2 ZZP P2, P3, M, T = 1 ZZP P2, P3, M, T = 1 The sequence with the minimum total transition time is the P4 P4 P2 P2 P4, P5, P1, P2, P3 P4, P5, P1, P2, P3 ZZP P3, P4, M, T = 1 ZZP P3, P4, M, T = 1 optimal sequence within time period t. P5, P1, P2, P3, P4 P5, P1, P2, P3, P4 ZZP P4, P5, M, T = 1 ZZP P4, P5, M, T = 1 P4 P4 P3 P3 P1, P2, P3, P4, P5 P1, P2, P3, P4, P5 ZZP P5, P1, M, T = 1 ZZP P5, P1, M, T = 1 ∑ = ∀ YP ZP i m t , , imt ii mt ' i ' ∑ = ∀ YP ZP i m t ', , i mt ii mt ' ' i ∧ ∧ ∧ ∧ ⇔ ⇔ Generate the cycle and break the cycle to find the [ [ ] ] ≠ ¬ ≠ ¬ ∀ ∀ YP YP YP YP ZP ZP i m t i m t , , , , optimum sequence where transition times are minimized. imt imt i mt i mt iimt iimt ' ' i i ' ' i i ≥ ≥ ∀ ∀ YP YP ZP ZP i m t i m t , , , , imt imt i i m t i i m t , , , , , , + + ≤ ≤ ∀ ∀ ≠ ≠ ZP ZP YP YP i i i i i m t i m t 1 1 , ' , ' , , , , i i m t i i m t i m t i m t , , , , , , ', ', , , ∑ ∑ ≥ ≥ − − ∀ ∀ ZP ZP YP YP YP YP i m t i m t , , , , i i m t i i m t i m t i m t i m t i m t , , , , , , , , , , ', ', , , ≠ ≠ i i i i ' ' ∑∑ = ∀ ZZP m t 1 , ii mt ' i i ' ≤ ∀ ZZP ZP i i m t , ', , ii mt ii mt ' ' Having determining the sequence, we can determine the total transition time within each week. 64
Limitation: Large Problems The proposed planning models may be expensive to solve for long term horizons. ROLLING HORIZON APPROACH : ROLLING HORIZON APPROACH : Problem 1 Problem 1 Model B Model A Model A Model A Model A Week 1 Problem 2 Problem 2 Fixed Model B Model A Model A Model A Week 2 Problem 3 Problem 3 Fixed Fixed Model B Model A Model A Week 3 � The detailed planning period (Model B) moves as the model is solved in time. � Future planning periods include only underestimations for transition times. 65 *Ref. Dimitriadis et al, 1997
EXAMPLE: 5 Products, 2 Reactors, 1 Week Method Number of Number of Number of Time Solution binary continuous Equations (CPUs) ($) variables variables Relaxed Planning 20 49 67 0.046 1,680,960.0 A % 6.484 Difference Detailed Planning 140 207 335 0.296 1,571,960.0 B Scheduling 594 2961 2537 150 1,571,960.0 Obj Function Items ($) Relaxed Planning Detailed Planning Scheduling Detailed Planning and Sales 2,652,800 2,440,000 2,440,000 Operating Costs 971,840 868,000 868,000 Scheduling are Identical! Transition Costs 0 40 40 Inventory Costs 0 0 0 Gantt Chart: 66
EXAMPLE 2 - 15 Products, 6 Reactors, 48 Weeks � Determine the plan for 15 products, 6 reactors plant so as to maximize profit. A Reaction 1 B Reaction 2 • 15 Products, A,B,C,D,E,F,G,H,J,K,L,M,N,O,P R1 C Reaction 3 • B, G and N are produced in 2 stages. D Reaction 4 • 6 Reactors, R1,R2,R3,R4,R5,R6 E Reaction 5 R2 • End time of the week is defined as due dates Reaction 6 F • Demands are lower bounds Reaction 7 G R3 H Reaction 8 J Reaction 9 K Reaction 10 R4 L Reaction 11 M Reaction 12 R5 N Reaction 13 Relaxed planning yields 21% overestimation of profit O Reaction 14 P Reaction 15 R6 number of number of binary continuous number of time solution method variables variables equations (CPU s) ($) relaxed planning (A) 2,592 5,905 9,361 362 224,731,683 rolling horizon (RH) 10,092 25,798 28,171 11,656 184,765,965 rolling horizon (RH**) 1,950 25,798 28,171 4,554 182,169,267 67
Decomposable MILP Problems Complicating Constraints Complicating Variables x 1 x 2 x 3 y x 1 x 2 x 3 A D 1 D 1 A D 2 D 2 D 3 D 3 ∑ T + T T c x a y c x max max i i complicating complicating = = i n 1,.. st Ax b constraints variables + = = st Ay D x d i n 1,.. = = D x d i n i i i 1,.. i i i ≥ ≥ = y x i n 0, 0, 1,.. ∈ = = ≥ i x X x x i n x { , 1,.. , 0} i i Benders decomposition Note: can reformulate by defining Lagrangean decomposition = y y + Complicating constraints i i 1 and apply Lagrangean decomposition 68
Lagrangean Relaxation (Fisher, 1985) � MILP optimization problems can often be modeled as problems with complicating constraints. � The complicating constraints are added to the objective function (i.e. dualized) with a penalty term (Lagrangean multiplier) proportional to the amount of violation of the dualized constraints. � The Lagrangean problem is easier to solve (eg. can be decomposed) than the original problem and provides an upper bound to a maximization problem. 69
Lagrangean Relaxation = Z cx max ≤ s t Ax b . . (IP) ≤ Dx e ∈ n Assume integers only x Z + Easily extended cont. vars. Ax ≤ b Assume that is complicating constraint = + − Z u cx u b Ax ( ) max ( ) LR ≤ Dx e ∈ n x Z + ≥ where u Lagrange multipliers 0 70
Lagrangean Relaxation = max = + − Z cx Z u cx u b Ax ( ) max ( ) LR ≤ Ax b ≤ Dx e Complicating Constraint ≤ Dx e ∈ n x Z + ∈ n x Z ≥ where u 0 + This is a relaxation of original problem because: Ax ≤ b i) removing the constraint relaxes the original feasible space, − Ax ≥ ≥ b Z LR u Z ( ) 0 ( ) ii) always holds as in the original space since ≥ u 0 . and Lagrange multiplier is always ⇒ ≥ Z u Z ( ) Lagrangean Relaxation Yields Upper Bound LR 71
Lagrangean Relaxation = Z cx max Original problem: ≤ s t Ax b . . ≤ Dx e ∈ n x Z + Z LR u ( 1 ) = + − Z u cx u b Ax ( ) max ( ) Relaxed problem: LR Z LR u ≤ ( 2 ) Dx e ∈ n x Z Z LR u ( 3 ) + Z dual D = Z Z u Lagrangean dual: min ( ) gap Z D LR ≥ u 0 72
Graphical Interpretation = + − Z u cx u b Ax ( ) max ( ) Relaxed problem: LR ≤ Dx e ∈ n x Z + = Z Z u Lagrangean dual: min ( ) D LR ≥ u 0 { } = + − Z cx u b Ax min max ( ) Combine Relaxed and D ≥ ≥ u x 0 0 Lagrangean Dual Problems: ≤ Dx e ∈ n x Z + 73
Graphical Interpretation ′ { } = Z cx max = + − D Z cx u b Ax min max ( ) D ≤ ≥ Ax b ≥ u x 0 0 ≤ Dx e ∈ ∈ n ≤ x Conv x Z Dx e ( , ) Nice Proof + ∈ n x Z Frangioni (2005) ≥ + x 0 Optimization of Lagrange multipliers (dual) can be interpreted as optimizing the primal objective function on the intersection of the convex hull of non- { } ∈ n ≤ , x x Z Dx e complicating constraints set and the LP relaxation of the + { } ∈ ≤ , n x A b x Z x relaxed constraints set . + 74
Graphical Interpretation ′ = Z cx max D ≤ Ax b cx ∈ ∈ ≤ n x Conv x Z Dx e ( , ) + ≥ x 0 Z LP Z D Z dual gap Conv { } ∈ n ≤ , x x Z Dx e + { } ≤ x Dx e { } ≤ x Ax b 75
Theorem Lagrangean relaxation yields a bound at least as tight as LP relaxation ≤ ≤ ≤ Z P Z Z u Z ( ) ( ) D LR LP cx Z LP Z D Z Conv { } ∈ ≤ , n x x Z Dx e + { } ≤ x Dx e { } ≤ x Ax b 76
Lagrangean Decomposition (Guignard & Kim, 1987) � Lagrangean Decomposition is a special case of Lagrangean Relaxation. � Define variables for each set of constrain, add constraints equating different variables (new complicating constraints) to the objective function with some penalty terms. ′ = = + − Z cx Z v cx v y x = max Z cx ( ) max ( ) max LD ≤ Ax b ≤ ≤ Ax b s t Ax b . . ≤ Dy e ≤ ≤ Dx e Dy e = x y ∈ n ∈ x Z + n x Z + New complicating ∈ n x Z ∈ constraints n y Z + + ∈ n y Z + Dualize x = y 77
Lagrangean Decomposition = + − Z v cx v y x ( ) max ( ) LD ≤ Ax b ≤ Dy e ∈ n x Z + ∈ n y Z + Subproblem 2 Subproblem 1 = − = max Z v c v x Z v vy ( ) max ( ) ( ) LD LD 1 2 ≤ ≤ Ax b Dy e ∈ ∈ n n x Z y Z + + ( ) = + Lagrangean Z Z v Z v min ( ) ( ) LD LD LD dual 1 2 ≥ v 0 78
Notes � Lagrangean decomposition is different from other possible relaxations because every constraint in the original problem appears in one of the subproblems. Graphically: The optimization of Lagrangean multipliers can be interpreted as optimizing the primal objective function on the intersection of the convex hulls of constraint sets. Subproblem 2 Subproblem 1 79
Graphical Interpretation? Subproblem 2 cx Subproblem 1 Z LP Z LR Z LD Z Conv { } ∈ n ≤ , x x Z Dx e + { } Conv { } ≤ x Dx e { } ∈ ≤ , n x A b x Z x ≤ x Ax b + Note: Z LR , Z LD refer to dual solutions 80
Theorem � The bound predicted by “Lagrangean decomposition” is at least as tight as the one provided by “Lagrangean relaxation” (Guignard and Kim, 1987) � For a maximization problem ≤ ≤ ≤ Z P Z Z Z ( ) LD LR LP Solution of Dual Problem Z LR or Z LD Piecewise linear => Non-differentiable minimum u or ν 81
How to iterate on multipliers u? Assuming Dx ≤ d is a bounded polyhedron (polytope) with extreme points = k x k K 1,2... , then + − ≤ ∈ = + − k k cx u b Ax Dx d x X cx u b Ax max{ ( ) , } max { ( )} = x k K 1,... => Dual problem + − = η η ≥ + − = k k k k K cx u b Ax cx u b Ax k K min max{ ( )} min{ ( ), 1,.. } ≥ ≥ = u k u 0 0 1,.. Cutting plane approach η min subgradient η ≥ + − = k k K n = no. extreme points s t cx u b Ax k K . . ( ), 1,.. n iteration n ≥ η ∈ u R 1 0, Note: x k generated from max{cx + u k (b-Ax) subproblems 82
Subgradient Optimization Approach = − k k s b Ax ( ) Subgradient + = + μ k k k u u s Steepest descent search 1 Update formula for multipliers (Fisher, 1985) + = 2 + α − − − k k LB k k k u 1 u Z Z b Ax b Ax ( )( ) / k LD α ∈ where [0,2] k Note: Can also use bundle methods for nondifferentiable optimization Lemarechal, Nemirovski, Nesterov (1995) 83
Solution of Langrangean Decomposition 1. Iterative search in multilpliers of dual 2. Perform branch and bound search Select MaxI, ε , a k where LP relaxation is replaced by Set UB = + ∞ , LB= - ∞ Lagrangean relaxation/decomposition to Solve (RP’) to find v 0 a) Obtain tighter bound b) Decompose MILP For k = 1..K Typically in Stochastic Programming Caroe and Schultz (1999) Solve (P1) and (P2): Goel and Grossmann (2006) Upper Bound Obtain Z LD Tarhan and Grossmann (2008) Remarks Solve (P) with fixed binaries Lower Bound 1. Methods can be extended to NLP, MINLP or use heuristics: Obtain Z LB 2. Size of dual gap depends greatly on how problems are decomposed | Z LD - Z LB |< ε ? YES Return Z LB & 3. From experience gap often decreases with k = k+1 or k=MaxI? Current Solution problem size. NO Update u k Notes: Heuristic due to dual gap 84 Obtaining Lower Bound might be tricky
Multisite Distribution Network Jackson, Grossmann,Wassick, Hoffman (2002) North America SI TE D SI TE A Latin America SI TE E Europe SI TE F SI TE B Africa/MidEast SI TE G SI TE C � Objective: Develop model and effective solution strategy for large-scale multiperiod planning with Nonlinear Process Models
Multisite Distribution Model •Develop Multisite Model to determine: 1)What products to manufacture in each site 2)What sites will supply the products for each market 3)Production and inventory plan for each site � Objective: Maximize Net Present Value • Challenges/Optimization Bottlenecks: Large-Scale NLP –Interconnections between time periods & sites/markets � Apply Lagrangean Decompostion Method
Spatial Decomposition = PR M PR M − PR M PR M PROFIT SCost SALES PCost PROD , , , , max * * S S S S ( ) + λ − PR M PR M PR M PROD SALES , , , S S S Market M SI TE S S M SITE CONSTRAINT S : Market CONSTRAINT S : ( ) ( ) PR M ≤ PR M ≤ f PROD , f SALES , 0 0 S S ( ) ( ) − + λ − λ PR M PR M PR M PR M PR M PR M PR M PR M PCost PROD PROD SCost SALES SALES , , , , , , , , max max S S S S S S S S Site SUBPROBLEM for all S (NLP) Market SUBPROBLEM for all M (LP)
Temporal Decomposition PR INV , − S t 1 •Decompose at each time period Market M SITE S •Duplicate variables for PR INV , Inventories for each time S t period •Apply Langrangean Market M SITE S Decomposition Algorithm PR INV , + S t 1
Multisite Distribution Model - Spatial •3 Multi-Plant Sites, 3 Geographic Markets •Solved with GAMS/Conopt2 # Time Variables/ Optimal Solution Full Space Lagrangean % Within Full Periods Solution Time Solution Time Optimal Constraints Profit (million-$) Solution (months) (CPU sec) (CPU sec) 2 3345 / 2848 164 52 10 10% 4 6689 / 5698 326 478 127 11% 6 10033 /8548 497 1605 279 9% 8 13377/11398 666 2350 550 9%
Multisite Distribution Model - Temporal • 3 Multi-Plant Sites, 3 Geographic Markets • Solved with GAMS/Conopt2 # Time Variables/ Optimal Solution Full Space Lagrangean % Within Full Periods Solution Time Solution Time Optimal Solution Constraints Profit (million-$) (months) (CPU sec) (CPU sec) 3 5230 / 5005 116.05 395 97 2.2 6 9973 / 8551 236.53 2013 138 2.3 12 19945 /17101 474.18 10254 278 2.2 Temporal much smaller gap! Reason: material balances not violated at each time period
Stochastic Inventory System Replenishment Inventory Level Lead Time Reorder Point Order placed Safety Stock Time • Inventory System under Demand Uncertainty � Total Inventory = Working Inventory (WI) + Safety Stock (SS) � Estimate WI with Economic Order Quantity (EOQ) model
Economic Order Quantity Model Replenishment Average Order quantity Inventory ( Q ) Inventory Level (Q/2) Constant Demand Rate = D Time � F = Fixed ordering cost for each replenishment � h = Unit inventory holding cost
Economic Order Quantity Model D Q FD 2 = ⋅ + ⋅ ⇒ = F Q h Q Total Cost * h 2 Order cost Holding cost Economic Order Quantity (EOQ) Annual Cost Total Cost Curve Holding Cost Curve Minimum Total Cost Order Cost Curve Order quantity Q Optimal Order Quantity (Q*)
( Q,r ) Inventory Policy Order Inventory Level placed Replenishment Order Quantity ( Q ) Reorder Point ( r ) Time Lead Time � When inventory level falls to r , order a quantity of Q � Reorder Point ( r ) = Demand over Lead Time
Stochastic Inventory Model Stochastic Inventory = Working Inventory (EOQ) + Safety Stock Inventory Level Reorder Point ( r ) Safety Stock Time Place Receive order order Lead Time Reorder Point = Expected Demand over Lead Time + Safety Stock
Safety Stock Level (Service Level) Lead time = L Safety Stock
Risk-Pooling Effect* • Single retailer: • Decentralized system: � Each retailer maintains its own inventory � Demand at each retailer is • Centralized system: � All retailers share common inventory � Integrated demand * GD Eppen, “ Effect of centralization in a multi-location newsboy problem ” , Management Science, 1979, 25(5), 498
Supply Chain Design with Stochastic Inventory Management • Given: A potential supply chain You, Grossmann (2008) � Including fixed suppliers, retailers and potential DC locations � Each retailer has uncertain demand, using ( Q , r ) policy � Assume all DCs have identical lead time L (lumped to one supplier) Suppliers Distribution Retailers Centers
Problem Statement • Objective: (Minimize Cost) � Total cost = DC installation cost + transportation cost + fixed order cost + working inventory cost + safety stock cost • Major Decisions (Network + Inventory) � Network: number of DCs and their locations, assignments between retailers and DCs (single sourcing), shipping amounts � Inventory: number of replenishment, reorder point, order quantity, neglect inventories in retailers supplier retailer DC Supplier Distribution Centers Retailers
v(x)= g + ax Working inventory cost transportation cost EOQ cost Annual EOQ cost at a DC: ordering cost
Recommend
More recommend
