Report

Combining Linear Programming Based Decomposition Techniques with Constraint Programming Menkes van den Briel Member of Research Staff NICTA and ANU [email protected] CP-based column generation Application Reference Urban transit crew management T.H. Yunes., A.V. Moura, C.C. de Souza. Solving very large crew scheduling problems to optimality. Proceedings of ACM symposium on Applied Computing, pages 446-451, 2000. T.H. Yunes., A.V. Moura, C.C. de Souza. Hybrid column generation approaches for urban transit crew management problems. Transportation Science 39(2):273-288, 2005. Travelling tournament K. Easton, G.L. Nemhauser, and M.A. Trick. Solving the travelling tournament problem: A combined integer programming and constraint programming approach. Proceedings of Practice and Theory of Automated Timetabling, volume 2740 of Lecture Notes in Computer Science, pages 100-112. Springer, 2002. Two-dimensional bin packing D. Pisinger, M. Sigurd. Using decomposition techniques and constraint programming for solving the two-dimensional bin-packing problem. Journal on Computing 19(1):36-51, 2007. Graph coloring S. Gualandi. Enhancing constraint programming-based column generation for integer programs. PhD thesis, Politechnico di Milano, 2008. Constrained cutting stock T. Fahle, M. Sellmann. Cost based filtering for the constrained knapsack problem. Annals of Operations Research 115(1):73-93, 2002. Employee timetabling S. Demassey, G. Pesant, L.M. Rousseau. A cost-regular based hybrid column generation approach. Constraints 11(4):315-333, 2006. Wireless mesh networks A. Capone, G. Carello, I. Filippini, S. Gualandi, F. Malucelli. Solving a resource allocation problem in wireless mess networks: A comparison between a CP-based and a classical column generation. Networks 55(3):221-233, 2010. Multi-machine scheduling R. Sadykov, L.A. Wolsey. Integer programming and constraint programming in solving a multimachine assignment scheduling problem with deadlines and release dates. Journal on Computing 18(2):209-217, 2006. CP-based column generation Application Reference Airline crew assignment U. Junker, S.E. Karisch, N. Kohl, B. Vaaben, T. Fahle, M. Sellmann. A framework for constraint programming based column generation. Proceedings of Principles and Practice of Constraint Programming, volume 1713 of Lecture Notes in Computer Science, pages 261-274, 1999. T. Fahle, U. Junker, S.E. Karisch, N. Kohl, M. Sellmann, B. Vaaben. Constraint programming based column generation for crew assignment. Journal of Hueristics 8(1):59-81, 2002. M. Sellmann, K. Zervoudakis, P. Stamatopoulos, T. Fahle. Crew assignment via constraint programming: integrating column generation and heuristic tree search. Annals of Operations Research 115(1):207-225, 2002. Vehicle routing with time windows L.M. Rousseau. Stabilization issues for constraint programming based column generation. Proceedings of Integration of AI and OR techniques in CP for Combinatorial Optimization, volume 3011 of Lecture notes in Computer Science, pages 402-408. Springer, 2004. L.M. Rousseau, M. Gendreau, G. Pesant, F. Focacci. Solving VRPTWs with constraint programming based column generation. Annals of Operations Research 130(1):199-216, 2004. CP-based Benders decomposition Application Reference Parallel machine scheduling V. Jain, I.E. Grossmann. Algorithms for hybrid MILP/CP models for a class of optimization problems. INFORMS Journal on Computing 13(4):258-276, 2001. Polypropylene batch scheduling C. Timpe. Solving planning and scheduling problems with combined integer and constraint programming. OR Spectrum 24(4):431-448, 2002. Call center scheduling T. Benoist, E. Gaudin, B. Rottembourg. Constraint programming contribution to Benders decomposition: A case study. Principles and Practice of Constraint Programming, volume 2470 of Lecture Notes in Computer Science, pages 603-617. Springer, 2002. Multi-machine scheduling J.N. Hooker. A hybrid method for planning and scheduling. Principles and Practice of Constraint Programming, volume 3258 of Lecture Notes in Computer Science, pages 305-316. Springer, 2004. J.N. Hooker. Planning and scheduling to minimize tardiness. Principles and Practice of Constraint Programming, volume 3709 of Lecture Notes in Computer Science, pages 314-327. Springer, 2005. CP versus IP CP IP Variables Finite domain Continuous, Binary, Integer Constraints Symbolic: alldifferent cumulative Linear, algebraic: (+, –, *, =, ≤, ≥) Inference Constraint propagation LP relaxation Local Feasible Global Optimal CP versus IP • “MILP is very efficient when the relaxation is tight and models have a structure that can be effectively exploited” • “CP works better for highly constrained discrete optimization problems where expressiveness of MILP is a major limitation” • “From the work that has been performed, it is not clear whether a general integration strategy will always perform better than either CP or an MILP approach by itself. This is especially true for the cases where one of these methods is a very good tool to solve the problem at hand. However, it is usually possible to enhance the performance of one approach by borrowing some ideas from the other” – Source: Jain and Grossmann, 2001 Outline • • • • • Background Introduction Dantzig Wolfe decomposition Benders decomposition Conclusions What is your background? • Have implemented Benders and/or Dantzig Wolfe decomposition • Have heard about Benders and/or Dantzig Wolfe decomposition • Have seen Bender and/or Dances with Wolves Things to take away • A better understanding of how to combine linear programming based decomposition techniques with constraint programming • A better understanding of column generation, Dantzig Wolfe decomposition and Benders decomposition • A whole lot of Python code with example implementations Helpful installations • Python 2.6.x or 2.7.x – “Python is a programming language that lets you work more quickly and integrate your systems more effectively” – http://www.python.org/getit/ • Gurobi (Python interface) – “The state-of-the-art solver for linear programming (LP), quadratic and quadratically constrained programming (QP and QCP), and mixed-integer programming (MILP, MIQP, and MIQCP)” – http://www.gurobi.com/products/gurobi-optimizer/try-for-yourself • NetworkX – “NetworkX is a Python language software package for the creation, manipulation, and study of the structure, dynamics, and functions of complex networks” – http://networkx.lanl.gov/download.html Abbreviations • • • • • • • • • Artificial Intelligence (AI) Constraint Programming (CP) Constraint Satisfaction Problem (CSP) Integer Programming (IP) Linear Programming (LP) Mixed Integer Programming (MIP) Mixed Integer Linear Programming (MILP) Mathematical Programming (MP) Operations Research (OR) Outline • • • • • Background Introduction Dantzig Wolfe decomposition Benders decomposition Conclusions What is decomposition? • “Decomposition in computer science, also known as factoring, refers to the process by which a complex problem or system is broken down into parts that are easier to conceive, understand, program, and maintain” – Source: http://en.wikipedia.org/wiki/Decomposition_(computer_science) • Decomposition in linear programming is a technique for solving linear programming problems where the constraints (or variables) of the problem can be divided into two groups, one group of “easy” constraints and another of “hard” constraints “easy” versus “hard” constraints • Referring to the constraints as “easy” and “hard” may be a bit misleading – The “hard” constraints need not be very difficult in themselves, but they can complicate the linear program making the overall problem more difficult to solve – When the “hard” constraints are removed from the problem, then more efficient techniques could be applied to solve the resulting linear program Example G = (N, A), source s, sink t • Shortest path problem (P) Min (i,j)A cijxij s.t. 1 for i = s Source j:(i,j)A xij – j:(j,i)A xji = 0 for iN – {s, t} Flow -1 for i = t Sink xij {0, 1} • Resource constrained shortest path problem (NP-complete) Min (i,j)A cijxij s.t. 1 for i = s Source j:(i,j)A xij – j:(j,i)A xji = 0 for iN – {s, t} Flow -1 for i = t Sink (i,j)A dijxij ≤ C Capacity xij {0, 1} Example • Assignment problem (P) Max i=1,…, m, j=1,…,n cijxij s.t. j=1,…,n xij = 1 for 1 ≤ i ≤ m i=1,…,m xij = 1 for 1 ≤ j ≤ n xij {0, 1} m jobs, n machines Job Machine • Generalized assignment problem (NP-complete) Max i=1,…, m, j=1,…,n cijxij s.t. j=1,…,n xij = 1 for 1 ≤ i ≤ m Job i=1,…,m dijxij ≤ Cj for 1 ≤ j ≤ n Capacity xij {0, 1} Example • Consider developing a strategic corporate plan for several production facilities. Each facility has its own capacity and production constraints, but decisions are linked together at the corporate level by budgetary considerations Common constraints Facility 1 Independent constraints Facility 2 Facility n “easy” versus “hard” variables • Referring to the variables as “easy” and “hard” may be a bit misleading – The “hard” variables need not be very difficult in themselves, but they can complicate the linear program making the overall problem more difficult to solve – When the “hard” variables are removed from the problem, then more efficient techniques could be applied to solve the resulting linear program Example • Capacitated facility location problem (NP-complete) Min i=1,…,n,j=1,…,m cijxij + j=1,…,m fjyj s.t. i=1,…,m xij ≥ 1 for j = 1,…, n Demand j=1,…,n dixij ≤ Ciyi for i = 1,…, m Roll xij ≤ yi for i = 1,…, m j = 1,…, n Flow impl. xij ≥ 0 m facilities, n customers yi {0, 1} Example Common variables • Consider solving a multi period scheduling problem. Each period has its own set of variables but is linked together through resource consumption variables Period 1 Period 2 Period n Independent variables Outline • • • • • Background Introduction Dantzig Wolfe decomposition Benders decomposition Conclusions Background • Primal Min cx s.t. Ax ≥ b x≥0 [y] • Dual Max yTb s.t. yTA ≤ c y≥0 [x] Background • Primal Min cx s.t. Ax ≥ b x≥0 • Dual Max bTy s.t. ATy ≤ cT y≥0 [y] x c cx A Ax [x] y b bT bTy AT ATy cT Travelling salesman • G = (N, A), cost cij x y 0 20 19 1 1 1 2 17 15 3 14 6 4 12 12 5 12 3 6 9 8 7 15 20 8 19 11 9 7 5 7 0 2 4 6 3 9 1 8 5 Travelling salesman • G = (N, A), cost cij x y 0 20 19 1 1 1 2 17 15 3 14 6 4 12 12 5 12 3 6 9 8 7 15 20 8 19 11 9 7 5 7 0 2 4 6 3 9 1 8 5 Cost 60.78 Travelling salesman • Variables xij is 1 if arc (i, j) is on the shortest tour, 0 otherwise • Formulation Min (i,j)A cijxij s.t. i:(i,j)A xij = 1 for j N j:(i,j)A xij = 1 for i N i,jS:(i,j)A xij ≤ |S| – 1 for S N xij {0, 1} Inflow Outflow Subtour Travelling salesman • Variables xij is 1 if arc (i, j) is on the shortest tour, 0 otherwise • Formulation Min (i,j)A cijxij s.t. i:(i,j)A xij = 1 j:(i,j)A xij = 1 xij {0, 1} for j N for i N Inflow Outflow Example code Travelling salesman • G = (N, A), cost cij x y 0 20 19 1 1 1 2 17 15 3 14 6 4 12 12 5 12 3 6 9 8 7 15 20 8 19 11 9 7 5 7 0 2 4 6 3 9 1 8 5 Subtour 0, 2, 7 Travelling salesman • G = (N, A), cost cij x y 0 20 19 1 1 1 2 17 15 3 14 6 4 12 12 5 12 3 6 9 8 7 15 20 8 19 11 9 7 5 7 0 2 4 6 3 9 1 8 5 Subtour 0, 8, 1, 9 Travelling salesman • G = (N, A), cost cij x y 0 20 19 1 1 1 2 17 15 3 14 6 4 12 12 5 12 3 6 9 8 7 15 20 8 19 11 9 7 5 7 0 2 4 6 3 9 1 8 5 Subtour 0, 8, 2, 7 Travelling salesman • G = (N, A), cost cij x y 0 20 19 1 1 1 2 17 15 3 14 6 4 12 12 5 12 3 6 9 8 7 15 20 8 19 11 9 7 5 7 0 2 4 6 3 9 1 8 5 Cost 79.98 Travelling salesman • G = (N, A), cost cij x y 0 20 19 1 1 1 2 17 15 3 14 6 4 12 12 5 12 3 6 9 8 7 15 20 8 19 11 9 7 5 7 0 2 4 6 3 9 1 8 5 Cost 60.78 LPs with many constraints • The number of constraints that are tight (or active) is at most equal to the number of variables, so even with many constraints (possibly exponential many) only a small subset will be tight in the optimal solution A Active Non-active Row generation in the primal… x c cx A Ax b … is column generation in the dual y bT bTy AT ATy cT …and vice versa x y c cx A Ax bT bTy AT ATy cT b Column generation in the primal = Row generation in the dual Resource constrained shortest path • G = (N, A), source s, sink t, for each (i, j) A, cost cij, resource demand dij, and resource capacity C Capacity = 14 2 1,10 1 1,1 4 2,3 1,7 1,2 10,1 10,3 5,7 3 i 12,3 cij, dij 6 2,2 5 j Source: Desrosiers and Lübbecke, 2005 Resource constrained shortest path • G = (N, A), source s, sink t, for each (i, j) A, cost cij, resource demand dij, and resource capacity C Capacity = 14 2 1,10 1 1,1 4 2,3 1,7 1,2 10,1 10,3 5,7 3 i 12,3 cij, dij 6 2,2 5 j Cost 13 Demand 13 Resource constrained shortest path • Variables xij is 1 if arc (i, j) is on the shortest path, 0 otherwise • Formulation Min (i,j)A cijxij s.t. j:(i,j)A xij – j:(j,i)A xji = (i,j)A dijxij ≤ C xij {0, 1} 1 for i = s Source 0 for iN – {s, t} Flow -1 for i = t Sink Capacity Example code Resource constrained shortest path • Variables k is 1 if path k is the shortest path, 0 otherwise • Formulation Min kK ckk s.t. kK k = 1 kK dkk ≤ C k ≥ 0 Convex Capacity Arc versus path • Arc variables 2 • Path variables 4 1 2 6 4 6 1 3 5 3 5 2 4 2 4 1 6 3 5 1 6 3 5 Example code Revised Simplex method • Min cx s.t. Ax ≥ b x≥0 Add slack variables • Min z = cx s.t. Ax = b x≥0 • Let x be a basic feasible solution, such that x = (xB, xN) where xB is the vector of basic variables and xN is the vector of non-basic variables Revised Simplex method • Min z = cx s.t. Ax = b x≥0 • Min z = cBxB + cNxN s.t. BxB + ANxN = b xB, xN ≥ 0 x = (xB, xN), c = (cB, cN), A = (B, AN) Rearrange • Min z = cBxB + cNxN s.t. xB = B-1b – B-1ANxN xB, xN ≥ 0 Revised Simplex method • Min z = cBxB + cNxN s.t. xB = B-1b – B-1ANxN xB, xN ≥ 0 Substitute • Min z = cBB-1b + (cN – cBB-1AN)xN s.t. xB = B-1b – B-1ANxN xB, xN ≥ 0 Revised Simplex method • Min z = cBB-1b + (cN – cBB-1AN)xN s.t. xB = B-1b – B-1ANxN xB, xN ≥ 0 • At the end of each iteration we have – – – – – Current value of non-basic variables xN = 0 Current objective function value z = cBB-1b Current value of basic variables xB = B-1b Objective coefficients of basic variables 0 Objective coefficients of non-basic variables (cN – cBB-1AN) are the so-called reduced costs – With a minimization objective we want non-basic variables with negative reduced costs Revised Simplex method • Simplex algorithm 1. Select new basic variable (xN to enter the basis) 2. Select new non-basic variable (xB to exit the basis) 3. Update data structures Revised Simplex method • Simplex algorithm xS = b (slack variables equal rhs) x\S = 0 (non-slack variables equal 0) while minj{(cj – cBB-1Aj)} < 0 1. Select new basic variable j : (cj – cBB-1Aj) < 0 2. Select new non-basic variable j’ by increasing xj as much as possible 3. Update data structures by swapping columns between matrix B and matrix AN Example • Min z = – x1 – 2x2 s.t. – 2x1 + x2 ≥ 2 – x1 + 2x2 ≥ 7 x1 ≥ 7 x1 , x2 ≥ 0 • Min z = – x1 – 2x2 s.t. – 2x1 + x2 + x3 = 2 – x1 + 2x2 + x4 = 7 x1 + x5 = 7 x1 , x2 , x3 , x4 , x 5 ≥ 0 Add slack variables Example • Simplex method bsc x1 x2 x3 • Revised Simplex method x4 x5 rhs x2 bsc x3 x4 x5 rhs -z -1 -2 0 0 0 0 -2 -z 0 0 0 0 x3 -2 1 1 0 0 2 1 x3 1 0 0 2 x4 -1 2 0 1 0 7 2 x4 0 1 0 7 x5 1 0 0 0 1 3 0 x5 0 0 1 3 bsc x1 x2 x3 x4 x5 rhs x1 bsc x3 x4 x5 rhs -z -5 0 2 0 0 4 -5 -z 2 0 0 4 x2 -2 1 1 0 0 2 -2 x3 1 0 0 2 x4 3 0 -2 1 0 3 3 x4 -2 1 0 3 x5 1 0 0 0 1 3 1 x5 0 0 1 3 Example • Simplex method bsc x1 x2 x3 • Revised Simplex method x4 x5 rhs x3 -z 0 0 -3/4 5/3 0 9 -3/4 x2 0 1 -1/3 2/3 0 4 x1 1 0 -2/3 1/3 0 x5 0 0 2/3 -1/3 1 bsc x1 x2 x3 x4 x5 bsc -z x3 x4 x5 rhs 2 0 0 9 -1/3 x2 -1/3 2/3 0 4 1 -2/3 x1 -2/3 1/3 0 1 2 2/3 x5 2/3 -1/3 1 2 rhs bsc x3 x4 x5 rhs -z 0 0 0 1 2 13 -z 0 0 0 13 x2 0 1 0 1/2 1/2 5 x2 0 1/2 1/2 5 x1 1 0 0 0 1 3 x1 0 0 1 3 x3 0 0 1 -1/2 3/2 3 x3 1 -1/2 3/2 3 Column generation • Simplex algorithm xS = b (slack variables equal rhs) x\S = 0 (non-slack variables equal 0) while minj{(cj – cBB-1Aj)} < 0 1. Select new basic variable j : (cj – cBB-1Aj) < 0 2. Select new non-basic variable j’ by increasing xj as much as possible 3. Update data structures by swapping columns between matrix B and matrix AN In column generation, rather than checking the reduced cost for each variable, a subproblem is solved to find a variable with negative reduced cost LPs with many variables • The number of basic (non-zero) variables is at most equal to the number of constraints, so even with many variables (possibly exponential many) only a small subset will be in the optimal solution A xB xN Column generation • (cN – cBB-1AN) < 0 Substitute • (cN – yTAN) < 0 Column generation • (cN – yTAN) < 0 • Primal Min cx s.t. Ax ≥ b x≥0 • Dual Max yTb s.t. yTA ≤ c y≥0 x Column with negative reduced cost c cx A Ax Row with violated rhs b y bT bTy AT ATy cT Resource constrained shortest path • Variables k is 1 if path k is the shortest path, 0 otherwise • Formulation Min kK ckk s.t. kK k = 1 kK dkk ≤ C k {0, 1} Convex Capacity Resource constrained shortest path • Primal Min kK ckk s.t. kK k = 1 kK dkk ≤ C k ≥ 0 [] [] • Dual Max + C s.t. + dk ≤ ck = free ≤0 Need to find a path for which ck – – dk < 0 Implicitly search all paths by optimizing Min (i,j)A (cij – dij) s.t. Source, Flow, Sink [k] Resource constrained shortest path • G = (N, A), source s, sink t, for each (i, j) A, cost cij, resource demand dij, and resource capacity C Capacity = 14 1 2 1 1 4 2 1 1 6 10 10 5 3 i 2 12 (cij – dij) 5 j Resource constrained shortest path • Master Min kK ckk s.t. kK k = 1 kK dkk ≤ C k ≥ 0 • Subproblem Min (i,j)A (cij – dij)xij s.t. j:(i,j)A xij – j:(j,i)A xji = Convex Capacity 1 for i = s Source 0 for iN – {s, t} Flow -1 for i = t Sink • Add variable to master if (i,j)A (cij – dij)xij – < 0 Example code Cutting stock • Roll width W, m orders of di rolls of length li, i = 1,…, m 100 11 x 4x 4x 2x di 12 31 36 45 li Cutting stock • Roll width W, m orders of di rolls of length li, i = 1,…, m 12 12 31 12 12 12 12 36 12 12 36 36 36 12 31 12 12 31 31 45 45 11 x 4x 4x 2x di 12 31 36 45 li 96 96 98 100 100 Rolls 5 Cutting stock • Variables xik is the number of times order i is cut from roll k yk is 1 if roll k is used, 0 otherwise • Formulation Min k=1,…,K yk s.t. k=1,…,K xik ≥ di i=1,…,n lixik ≤ Wyk xik ≥ 0 and integer yk {0, 1} for i = 1,…, n for k = 1,…, K Demand Roll Example code Cutting stock • Variables k is the number of times cutting pattern k is used • Formulation Min kK k s.t. kK aikk ≥ di for i = 1,…, m k ≥ 0 and integer Demand Cutting stock • Cutting pattern variables k 12 12 36 36 aik [2, 0, 2, 0] k 12 aik 12 [2, 1, 0, 1] 31 11 x 4x 4x 2x 45 12 31 36 45 Cutting stock • Primal Min kK k s.t. kK aikk ≥ di k ≥ 0 [i] • Dual Max i=1,…,n dii s.t. i=1,…,n aiki ≤ 1 i ≥ 0 Need to find a cutting pattern for which 1 – i=1,…,n aiki < 0 Implicitly search all cutting patterns by optimizing Max i=1,…,n aii s.t. i=1,…,n liai ≤ W ai ≥ 0 and integer [k] Cutting stock • m items with value i and weight li, i = 1,…, m, maximum allowed weight W 100lbs 12 12 36 36 $0.125, 12lbs 0.125 12 $0.33, 31lbs 0.33 31 $0.50, 36lbs 0.50 36 $0.50, 45lbs 0.50 45 i li Cutting stock • Master Min kK k s.t. kK aikk ≥ di for i = 1,…, m k ≥ 0 • Subproblem Max i=1,…,m aii s.t. i=1,…,m liai ≤ W ai ≥ 0 and integer • Add variable to master if 1 – aii < 0 Demand Example code Generalized assignment • n jobs, m machines, cost cij, demand dij, capacity Ci Job 1 1 2 1 17, 8 23, 15 2 21, 15 16, 7 3 22, 14 21, 23 4 18, 23 16, 22 5 24, 8 17, 11 2 1 36 2 34 i Cj 3 4 5 j cij, dij Generalized assignment • n jobs, m machines, cost cij, demand dij, capacity Ci Job 1 1 2 1 17, 8 23, 15 2 21, 15 16, 7 3 22, 14 21, 23 4 18, 23 16, 22 5 24, 8 17, 11 2 1 36 30 2 34 29 3 4 Cost 95 5 j cij, dij i Cj Generalized assignment • Variables xij is 1 if job j is assigned to machine i, 0 otherwise • Formulation Max i=1,…,m,j=1,…,n cijxij s.t. i=1,…,m xij = 1 j=1,…,n dijxij ≤ Ci xij {0, 1} for 1 ≤ j ≤ n Job for 1 ≤ i ≤ m Capacity Example code Generalized assignment • Variables ik is 1 if machine i has job assignment k, 0 otherwise • Formulation Max i=1,…,m,k=1,…,Ki cikik s.t. i=1,…,m,k=1,…,Ki aijkik = 1 k=1,…,Ki ik = 1 ik {0, 1} for 1 ≤ j ≤ n Job for 1 ≤ i ≤ m Convexity Generalized assignment • Job assignment variables 1 2 ik 1 1 ik 3 4 5 aijk [1, 0, 1, 0, 1] 2 2 1 3 4 5 aijk [0, 1, 0, 1, 0] 2 Generalized assignment • Formulation Max i=1,…,m,k=1,…,Ki cikik s.t. i=1,…,m,k=1,…,Ki aijkik = 1 k=1,…,Ki ik = 1 ik {0, 1} for 1 ≤ j ≤ n Job for 1 ≤ i ≤ m Convexity Common constraints Machine 1 Independent constraints Machine 2 Machine n Generalized assignment • Primal Max i=1,…,m,k=1,…,Ki cikik s.t. i=1,…,m,k=1,…,Ki aijkik = 1 k=1,…,Ki ik = 1 ik ≥ 0 • Dual Min j=1,…,n j + i=1,…,m i s.t. j=1,…,n aijkj + i ≥ cik j = free i = free Need to find a cutting pattern for which j=1,…,n (cik – aijkj ) – i > 0 for i = 1,…,m Implicitly search all cutting patterns by optimizing Max j=1,…,n (cij – aijj ) s.t. j=1,…,n dijaij ≤ Ci aij ≥ 0 and integer Generalized assignment • n items with value j and weight dij, j = 1,…, n, maximum allowed weight W Job 36lbs $44.00, 8lbs 1 1 1 44, 8 2 55, 15 3 51, 14 4 52, 23 5 55, 8 Job 40, 15 2 37, 7 $51.00, 14lbs 3 43, 23 $52.00, 23lbs 4 34, 22 5 41, 11 $55.00, 8lbs 1 36 3 4 2 5 1 2 1 $55.00, 15lbs 2 2 1 3 4 5 2 34 Generalized assignment • Master Max i=1,…,m,k=1,…,Ki cikik s.t. i=1,…,m,k=1,…,Ki aijkik = 1 k=1,…,Ki ik = 1 ik {0, 1} for 1 ≤ j ≤ n Job for 1 ≤ i ≤ m Convexity • Subproblem (for each machine i) Max j=1,…,n (cij – aijj ) s.t. j=1,…,n dijaij ≤ Ci aij ≥ 0 and integer • Add variable to master if j=1,…,n (cij – aijj ) – i > 0 Example code History of column generation 1958: A suggested computation for maximal multicommodity network flows L.R. Ford and D.R. Fulkerson 1960: Decomposition principle for linear programs G.B. Dantzig and P. Wolfe “Credit is due to Ford and Fulkerson for their proposal for solving multicommodity network problems as it served to inspire the present development.” 1961: A linear programming approach to the cutting-stock problem P.C. Gilmore and R.E. Gomory 1963: A linear programming approach to the cutting-stock problem–Part II P.C. Gilmore and R.E. Gomory 1969: A column generation algorithm for a ship scheduling problem L.E. Appelgren Solving integer programs by column generation 1984: Routing with time windows by column generation Y. Dumas, F. Soumis and M. Desrochers 1998: Branch-and-price: column generation for solving huge integer programs C. Barnhart, E.L. Johnson, G.L. Nemhauser, M.W.P. Savelsbergh and P.H. Vance 2000: On Dantzig-Wolfe decomposition in integer programming and ways to perform branching in a branch-and-price algorithm F. Vanderbeck 2005: A primer in column generation J. Desrosiers and M.E. Lubbecke 2011: Branching in branch-and-price: a generic scheme F. Vanderbeck CP-based column generation 1999: A framework for constraint programming based column generation U. Junker, S.E. Karisch, N. Kohl, B. Vaaben, T. Fahle and M. Sellmann 2000: Solving very large crew scheduling problems to optimality T.H. Yunes, A.V. Moura and C.C. de Souza CP-based column generation Application Reference CP used to solve subproblem CP used within Branch-and-Price Urban transit crew management T.H. Yunes., A.V. Moura, C.C. de Souza. 2000. Y Y T.H. Yunes., A.V. Moura, C.C. de Souza. 2005. Y Y Travelling tournament K. Easton, G.L. Nemhauser, and M.A. Trick. 2002. Y Y Two-dimensional bin packing D. Pisinger, M. Sigurd. 2007. Y Y Graph coloring S. Gualandi. 2008. Y Y Constrained cutting stock T. Fahle, M. Sellmann. 2002. Y N Employee timetabling S. Demassey, G. Pesant, L.M. Rousseau. 2006. Y Y Wireless mesh networks A. Capone, G. Carello, I. Filippini, S. Gualandi, F. Malucelli. 2010. Y N Multi-machine scheduling R. Sadykov, L.A. Wolsey. 2006. Y N Source: Gualandi and Malucelli, 2009 CP-based column generation Application Reference CP used to solve subproblem CP used within Branch-and-Price Airline crew assignment U. Junker, S.E. Karisch, N. Kohl, B. Vaaben, T. Fahle, M. Sellmann. 1999. Y N T. Fahle, U. Junker, S.E. Karisch, N. Kohl, M. Sellmann, B. Vaaben. 2002. Y N M. Sellmann, K. Zervoudakis, P. Stamatopoulos, T. Fahle. 2002. Y N L.M. Rousseau. 2004. Y N L.M. Rousseau, M. Gendreau, G. Pesant, F. Focacci. 2004. Y Y Vehicle routing with time windows Source: Gualandi and Malucelli, 2009 CP-based column generation • Typical implementation Dual information Master Linear programming Subproblem New columns Constraint programming Outline • • • • • Background Introduction Dantzig Wolfe decomposition Benders decomposition Conclusions Two-stage optimization Stage 1 Solution values Stage 2 Benders decomposition Solution values Stage 1 Stage 2 Benders cuts Benders decomposition “Learn from ones mistakes” Distinguish primary variables from secondary variables Search over primary variables (master problem) For each trial value of primary variables, solve problem over secondary variables (subproblem) If solution is suboptimal/infeasible, find out why and design a constraint that rules out not only this solution but a large class of solutions that are suboptimal/infeasible for the same reason (Benders cut) Add Benders cut to the master problem and resolve Solution values Master Subproblem Benders cuts Capacitated facility location • m facilities, n customers, cost cij, demand dj, capacity Ci, fixed cost fi Cust 1 2 3 1 2 4 5 2 3 3 4 3 4 1 2 4 5 2 1 5 7 6 3 10, 3 10, 4 10, 4 Ci, fi 1 6 2 4 3 8 4 7 5 5 j dj 1 2 3 i cij Capacitated facility location • m facilities, n customers, cost cij, demand dj, capacity Ci, fixed cost fi Cust 1 2 3 1 2 4 5 2 3 3 4 3 4 1 2 4 5 2 1 5 7 6 3 10, 3 10, 4 10, 4 6 2 4 3 8 4 7 5 5 j dj 1 2 3 Cost 21.29 Ci, fi 1 i cij Capacitated facility location • Variables xij fraction of demand supplied by facility i to cusomter j yi is 1 if facility i is open, 0 otherwise • Formulation Min i=1,…,n,j=1,…,m cijxij + j=1,…,m fjyj s.t. i=1,…,m xij ≥ 1 for j = 1,…, n j=1,…,n dixij ≤ Ciyi for i = 1,…, m xij ≤ yi for i = 1,…, m j = 1,…, n xij ≥ 0 yi {0, 1} Demand Roll Flow Example code Benders decomposition • Min cx + dy s.t. Ax ≥ b Px + Qy ≥ r x ≥ 0 and integer y≥0 Master • Min cx + s.t. Ax ≥ b x ≥ 0 and integer ≥0 • Solution values Benders cuts Subproblem Min dy s.t. Qy ≥ r – Px y≥0 What if the subproblem is infeasible? Benders decomposition • Primal, dual possibilities Dual Primal Optimal Unbounded Infeasible Optimal Yes No No Unbounded No No Yes Infeasible No Yes Yes Benders decomposition • Min cx + dy s.t. Ax ≥ b Px + Qy ≥ r x ≥ 0 and integer y≥0 Master • Min cx + s.t. Ax ≥ b optimality cuts feasibility cuts x ≥ 0 and integer ≥0 • Solution values Benders cuts Subproblem Min dy s.t. Qy ≥ r – Px y≥0 Benders decomposition • Min dy s.t. Qy ≥ r – Px y≥0 [u] • Max uT(r – Px) s.t. uTQ ≤ d u≥0 • Optimal • Optimality cut ≥ ukT(r – Px) • Infeasible • Infeasibility cut vkT(r – Px) ≤ 0 [y] Benders decomposition • Min cx + dy s.t. Ax ≥ b Px + Qy ≥ r x ≥ 0 and integer y≥0 Master • Min cx + s.t. Ax ≥ b ≥ ukT(r – Px) vkT(r – Px) ≤ 0 x ≥ 0 and integer ≥0 • Solution values Benders cuts Subproblem Max dy s.t. Qy ≤ r – Px y≥0 Benders decomposition START Solve master problem yes Terminate? no Solve sub problem Add optimality cut yes Is optimal? END Add feasibility cut no Capacitated facility location • Variables xij fraction of demand supplied by facility i to cusomter j yi is 1 if facility i is open, 0 otherwise • Formulation Min i=1,…,n,j=1,…,m cijxij + j=1,…,m fjyj s.t. i=1,…,m xij ≥ 1 for j = 1,…, n j=1,…,n dixij ≤ Ciyi for i = 1,…, m xij ≤ yi for i = 1,…, m j = 1,…, n xij ≥ 0 yi {0, 1} Demand Roll Flow Capacitated facility location • Master Min j=1,…,m fjyj + s.t. optimality cuts feasibility cuts yi {0, 1} ≥0 • Subproblem Min i=1,…,n,j=1,…,m cijxij s.t. i=1,…,m xij ≥ 1 for j = 1,…, n j=1,…,n dixij ≤ Ciyi for i = 1,…, m xij ≤ yi for i = 1,…, m j = 1,…, n xij ≥ 0 Demand Roll Flow Capacitated facility location • Subproblem primal Min i=1,…,n,j=1,…,m cijxij s.t. i=1,…,m xij ≥ 1 j=1,…,n dixij ≤ Ciyi xij ≤ yi xij ≥ 0 [j] [i] [ij] • Subproblem dual Max j=1,…,m j + i=1,…,n Ciyii + i=1,…,n,j=1,…,m yiij s.t. j + dii + ij ≥ 1 [xij] j ≥ 0 i ≤ 0 ij ≤ 0 Capacitated facility location • Master Min j=1,…,m fjyj + s.t. ≥ j=1,…,m j + i=1,…,n Ciiyi + i=1,…,n,j=1,…,m ij yi j=1,…,m j + i=1,…,n Ciiyi + i=1,…,n,j=1,…,m ij yi ≤ 0 yi {0, 1} ≥0 Example code Benders decomposition for stochastic prog. Scenario 1 Master Scenario 2 Scenario 3 Capacitated facility location • m facilities, n customers, cost cij, demand dj, capacity Ci, fixed cost fi Cust 1 2 3 1 2 4 5 2 3 3 4 3 4 1 2 4 5 2 1 5 7 6 3 10, 3 10, 4 10, 4 Ci, fi 1 6 5 4 2 4 3 2 3 8 7 6 4 7 6 5 5 5 4 3 j dj 1 2 3 i cij Example code CP-based Benders decomposition • Typical implementation(?) Solution values Master Constraint programming Subproblem Benders cuts Linear programming CP-based Benders decomposition • Recent developments Solution values Master Integer programming Subproblem Benders cuts Constraint programming CP-based Benders decomposition Application Reference Master problem Subproblem Parallel machine scheduling V. Jain, I.E. Grossmann. 2001. MILP CP Polypropylene batch scheduling C. Timpe. 2002. MILP CP Call center scheduling T. Benoist, E. Gaudin, B. Rottembourg. 2002. CP LP Multi-machine scheduling J.N. Hooker. 2004. MILP CP J.N. Hooker. 2005. MILP CP Source: Hooker, 2006 Nested Benders decomposition • Nested Benders decomposition – When the subproblem is decomposed into a master and subproblem Master Master Forward pass Solve master problems Sub Backward pass Solve subproblems and add Benders cuts Sub Master Master Sub Sub Outline • • • • • Introduction Background Dantzig Wolfe decomposition Benders decomposition Conclusions Why use decomposition? • Many real-world systems contain loosely connected components, and as a result, the corresponding mathematical models present a certain structure that can be exploited • It may be your only choice when solving the model without decomposition is impossible, because it is too large (memory error or timeout) When is decomposition likely most effective? • When you have either complicating constraints or complicating variables Dantzig Wolfe decomposition Benders decomposition Further reading • Column Generation – Guy Desaulniers, Jacques Desrosiers, Marius M. Solomon • Decomposition Techniques in Mathematical Programming – Antonio J. Conejo, Enrique Castillo, Roberto Minguez and Raquel Garcia-Bertrand • Linear Programming and Network Flows – Mokhtar S. Bazaraa, John J. Jarvis, Hanif D. Sherali From imagination to impact