Constraint Programming My View Hakan Kjellerstrand http www

Constraint Programming: My View Hakan Kjellerstrand http: //www. hakank. org/index_eng. html My Constraint Programming Blog: http: //www. hakank. org/constraint_programming_blog/

Overview - Presentation of me - What is Constraint Programming - why is it so fascinating? - Key concepts in CP - Examples of modelling with code (G 12 Mini. Zinc and or-tools) - Questions

Hakan Kjellerstrand (hakank@gmail. com) - http: //www. hakank. org/index_eng. html (on the web since 1995) - Software Developer at Bokus (http: //www. bokus. com/), Swedish on-line book store. 1997 -2001, 2005 -today. - Software developer since 1996, mostly web/internet/dot. com - Technical Writer and support, Swedish word processor Cicero 1987– 1994 - Malmoe General Hospital, EDB department etc 1982 -1987 - Education: Political Science, Philosophy and Computer Science Search for “Hakan Kjellerstrand”; that's me.

Constraint Programming - Pure private interest, “private researcher” - Checked out CLP in 1998 and later 2001 - 2007 autumn: Mathematical Programming/Operations Research - 2008 February: Found G 12 Mini. Zinc (blogged at Swedish blog) - Late 2008: Started “My Constraint Programming Blog” - Have now tested about 18 CP systems, 15 of them in public (bug reports, suggestions, examples) - At least one academic paper in the pipeline (as co-author) Other interests: - programming languages - machine learning/data mining, symbolic regression - magic (though not performing any more)

Constraint Programming systems Mini. Zinc : http: //hakank. org/minizinc/zinc. html Choco : http: /hakank. org/choco/ Ja. Co. P : http: //hakank. org/Ja. Co. P/Scala : http: //hakank. org/jacop_scala. html Gecode/R : http: //hakank. org/gecode_r/ Comet : http: //hakank. org/comet/ Gecode : http: //hakank. org/gecode/ ECLi. PSe CLP : http: //hakank. org/eclipse// Tailor/Essence' : http: //hakank. org/tailor/ SICStus Prolog : http: //hakank. org/sicstus/ or-tools/Python : http: //hakank. org/or-tools/#python or-tools/Java : http: //hakank. org/or-tools/#java or-tools/C# : http: //hakank. org/or-tools/#csharp Answer Set Programming : http: //hakank. org/asp/ 861 models 39 models 20 models 18 models 38 models 27 models 166 models 162 models 176 models 29 models 150 models 202 models 36 models 107 models 84 models Common constraint models: http: //hakank. org/common_cp_models/

What is CP used for? In industry CP is used for (among much other things): - scheduling and resource allocation - staff rostering - complex configuration problems - DNA sequencing - vehicle routing - analysis of analog circuits -. . . - in general: combinatorial (optimization)

What do I do with CP? - Understanding CP/OR concepts * bin packing, seating problems, TSP, . . . - Understanding/experimenting with math/CS concepts * clique, hitting set, graph related, . . . - Exhausting/checking combinatorial structures * 17 x 17 problem, Graeco-Latin (magic) squares, . . . - Modelling grid/pencil games * Sudoku, Ken, Killer. Sudoku, Strimko, Rogo, . . . - Recreational math/CS puzzles * 8809 = 6 problem, alphametic problems, . . . - Planning problems (optimal steps for a plan) * 15 -puzzle, M 12, etc.

My main interests and goals with CP Main interest in CP: - learning CP - modelling - testing CP systems - global constraints (decompositions) Method: - Explore (testing a lot of CP systems) - Exploit (modelling) - Explain (blogging) In general: - trying to make CP more known - solve specific problems

What is CP? Some key concepts Techniques: - Constraints (CSP, constraint satisfaction problem) - Constraint store - Search tree - Combinatorial optimization - Decision variables (vs. parameters) - Domains (often finite) - Propagator, special algorithm for a constraint - Fix point Modelling: - Declarative - Global constraints (alldifferent, element, etc) - Decomposition (high level version of a global constraint) - Search strategies, to guide traversing the search tree - Reifications - Bi-direction (multiple flow patterns) - Channelling/dual model - Symmetry breaking

Main principle of CP (simplified) - All variables are declared with a DOMAIN (often finite domain), including arrays, matrices, etc. - The solver PROPAGATES the CONSTRAINTS and variables in the model until: * a FIXPOINT is reached (no changes in the domains) * or assign the variables in the search tree, perhaps via backtracking - If a DOMAIN of a variable is empty then the assignment is illegal and the solver has to backtrack.

Time for some code! G 12 Mini. Zinc for high level concepts: - http: //www. g 12. cs. mu. oz. au/minizinc/ - http: //hakank. org/minizinc/ or-tools (mostly Python and C#): - http: //code. google. com/p/or-tools/ - http: //hakank. org/or-tools/

Sudoku 4 x 4 4 _ 3 1 _ _ _ _ 4 1 _ 2 Constraints: - all rows must have different digits - all columns must have different digits - all boxes (2 x 2) must have different digits Let's see how this can be stated in a CP system.

Sudoku 4 x 4 Constraints: - alldifferent(each ROW) - alldifferent(each COLUMN) - alldifferent(each BLOCK) % each 2 x 2 block Where alldifferent(ARRAY) can be defined as a DECOMPOSITION in Mini. Zinc: predicate alldifferent(array[int] of var int: x) = forall(i, j in 1. . length(x) where i < j) ( x[i] != x[j] ) ;

Sudoku 4 x 4 – Mini. Zinc code % init of hints, search, and output missing int: r = 2; int: n = r*r % declare decision variables array[1. . n, 1. . n] of var 1. . n: x; constraint % rows forall(i in 1. . n) (alldifferent([x[i, j] | j in 1. . n])) / % columns forall(j in 1. . n) (alldifferent([x[i, j] | i in 1. . n])) / % blocks forall(i in 0. . r-1, j in 0. . r-1) ( alldifferent([x[r, c] | r in i*r+1. . i*r+r, c in j*r+1. . j*r+r]) ) ;

Sudoku 4 x 4 – or-tools/Python block_size = 2; line_size = block_size ** 2 line = range(0, line_size); block = range(0, block_size) x = {} for i in line: % declare decision variables for j in line: x[(i, j)] = s. Int. Var(1, line_size) for i in line: % rows s. Add(s. All. Different([x[(i, j)] for j in line])) for j in line: % columns s. Add(s. All. Different([x[(i, j)] for i in line])) for i in block: % blocks for j in block: one_block = [] for di in block: for dj in block: one_block. append(x[(i * block_size + di, j * block_size + dj)]) s. Add(s. All. Different(one_block))

Sudoku 4 x 4 – or-tools/C# int bs = 2; // Block size IEnumerable<int> BLOCK = Enumerable. Range(0, bs); int n = bs * bs; // Grid size IEnumerable<int> RANGE = Enumerable. Range(0, n); Int. Var[, ] x = s. Make. Int. Var. Matrix(n, n, 1, n, "x"); foreach(int i in RANGE) { s. Add((from j in RANGE // Rows select x[i, j]). To. Array(). All. Different()); s. Add((from j in RANGE // Columns select x[j, i]). To. Array(). All. Different()); } foreach(int i in BLOCK) { // Blocks foreach(int j in BLOCK) { s. Add((from di in BLOCK from dj in BLOCK select x[i*bs+di, j*bs+dj] ). To. Array(). All. Different()); } }

Sudoku 4 x 4 – propagation 4 _ _ _ 3 1 _ _ _ 4 1 _ 2 Solution (unique) 4 2 3 1 1 3 2 4 2 3 1 4 4 1 3 2 How does a CP solver reach this solution?

Sudoku 4 x 4 – propagation example (simplified) 4 3 1 1234 1234 1234 4 1 2 1234 Add DOMAINS (1. . 4) to all unknown variables. Hints are FIXED already. Now we will propagate three alldifferent constraints: - alldifferent(ROW) - alldifferent(COLUMN) - alldifferent(BLOCK) This is a very simplified example. Real CP systems use more intelligent propagation.

Sudoku 4 x 4 – simple propagation example 4 3 1 2 1234 1234 4 1 2 1234 Cell (1, 1): Fixed value (4). Cell (1, 2): Reduce: - remove 4 (row, block) - remove 1 (column, block) - remove 3 (block) → Single value: 2

Sudoku 4 x 4 – simple propagation example 4 3 1 2 1234 1 3 1234 4 1 2 1234 Cell (1, 3): Reduce: - remove 4 (row, column) - remove 2 (row) → {1 3}

Sudoku 4 x 4 – simple propagation example 4 3 1 2 1 3 1234 Cell (1, 4): Reduce: - remove 4 (row) - remove 1 (column) - remove 2 (column) → 3 Note: Here we don't go back to fix cell (1, 3). 1234 4 1 2 1234 Cell (2, 1): fixed (3) Cell (2, 2): fixed (1)

Sudoku 4 x 4 – simple propagation example 4 3 1 2 1234 1 3 3 2 1234 4 1 2 1234 Cell (2, 3): Reduce: - remove 3 (row) - remove 1 (row) - remove 4 (column) → 2

Sudoku 4 x 4 – simple propagation example 4 3 1 2 1234 1 3 3 2 4 1 2 1234 4 Cell (2, 4): Reduce: - remove 3 (row) - remove 1 (row, column) - remove 2 (row, column) → 4

Sudoku 4 x 4 – simple propagation example 4 3 1 2 2 1234 1 3 3 2 4 1 2 1234 4 Cell (3, 1): Reduce: - remove 4 (row, column) - remove 1 (row) - remove 3 (column) → 2

Sudoku 4 x 4 – simple propagation example 4 3 1 2 2 1234 1 3 3 2 3 1234 4 1 2 1234 4 Cell (3, 2): Reduce: - remove 1 (row, column, block) - remove 2 (row) - remove 4 (row) → 3

Sudoku 4 x 4 – simple propagation example 4 3 1 2 2 1234 1 3 2 3 1234 4 1 2 1234 Cell (3, 3): Fixed. Cell (3, 4): Fixed. 3 4

Sudoku 4 x 4 – simple propagation example 4 3 1 2 2 1 1 3 3 2 3 1234 4 1 2 1234 4 Cell (4, 1): Reduce - remove 2 (row) - remove 4 (column) - remove 3 (column) → 1

Sudoku 4 x 4 – simple propagation example 4 3 1 2 2 1 1 3 3 2 3 4 4 1 2 1234 4 Cell (4, 2): Reduce - remove 1 (row) - remove 2 (row) - remove 3 (column) → 4

Sudoku 4 x 4 – simple propagation example 4 3 1 2 1 3 3 2 4 Cell (4, 3): Reduce - remove 2 (row, block) - remove 4 (column, block) – - remove 1 (block) → 3 Cell (4, 4): Fixed 2 2 1 3 4 4 1 2 3 Are we finished? No! There is still a variable/cell with no single assignment, i. e. Cell (1, 3).

Sudoku 4 x 4 – simple propagation example 4 3 1 2 1 Cell (1, 3): Reduce - remove 3 (row, block) → 1 3 2 4 And now all variables has been assigned to a single value. 2 1 3 4 4 1 2 3

Sudoku 4 x 4 – simple propagation example 4 2 3 1 1 3 2 4 2 3 1 4 4 1 3 2 … and we got a solution! It is unique – as a Sudoku should be.

Sudoku 4 x 4: Points to take home - Domain reduction is one of the key principles to CP. - Using search heuristics (strategies) influences the order of variable/value selection.

Another model: Minesweeper – in this version – is a simple grid problem: . . 2. 3. 2. . . . 24. 3 1. 34. . . . 3. 3. 3. . - Each number represents how many bombs there are in the nearby cells - The “. ” (dot) represents an unknown cell: either a bomb or empty cell. - Problem: Where are the bombs?

Minesweeper – Mini. Zinc version (the setup) int: X = -1; % the unknowns % >= 0 for number of mines in the neighbourhood array[1. . r, 1. . c] of -1. . 8: game; % decision variables: 0/1 for no bomb/bomb array[1. . r, 1. . c] of var 0. . 1: mines; % problem instance int: r = 6; % rows int: c = 6; % column game = array 2 d(1. . r, 1. . c, [ X, X, 2, X, 3, X, 2, X, X, 2, 4, X, 3, 1, X, 3, 4, X, X, 3, X, X, ]);

Minesweeper – Mini. Zinc version (constraint) % game[1. . n, 1. . n]: the given hints % mines[1. . n, 1. . n]: 0/1 where 1 represent a bomb % X: -1 represents the unknown constraint forall(i in 1. . r, j in 1. . c) ( ( (game[i, j] >= 0 ) → % the hint number must be the number % of all the surrounded bombs game[i, j] = sum(a, b in {-1, 0, 1} where i+a > 0 / j+b > 0 / i+a <= r / j+b <= c ) (mines[i+a, j+b]) ) / % if a hint, then it can't be a bomb (game[i, j] > X -> mines[i, j] = 0) ) ;

Minesweeper – Mini. Zinc version Declarative aspect of Constraint Programming. The ideal: - Just state the requirements (the constraints) - It's now up to the CP solver to find the solution. %. . . % the hint number must be the number % of all the surrounded bombs game[i, j] = sum(a, b in {-1, 0, 1} where i+a > 0 / j+b > 0 / i+a <= r / j+b <= c ) (mines[i+a, j+b]) %. . .

Minesweeper - Solution. . 2. 3. 2. . . . 24. 3 1. 34. . . . 3. 3. 3. . 100001 010110 000010 011100 100011 % 1: Bomb, 0: no bomb

Minesweeper – or-tools/Python mines = {} # define decision variables for i in range(r): for j in range(c): mines[(i, j)] = solver. Int. Var(0, 1) # constraints for i in range(r): for j in range(c): if game[i][j] >= 0: solver. Add( game[i][j] == solver. Sum([mines[i+a, j+b] for a in S for b in S if i + a >=0 and j + b >=0 and i + a < r and j + b < c])) if game[i][j] > X: # This cell cannot be a mine solver. Add(mines[i, j] == 0)

Minesweeper – or-tools/C# // constraints for(int i = 0; i < r; i++) { for(int j = 0; j < c; j++) { if (game[i, j] >= 0) { var tmp = from a in S from b in S where i + a >= 0 && j + b >= 0 && i + a < r && j + b < c select(mines[i+a, j+b]); solver. Add(tmp. To. Array(). Sum() == game[i, j]); } if (game[i, j] > X) { solver. Add(mines[i, j] == 0); } } }

Minesweeper – other implementations Implementations of Minesweeper using the same (as possible) approach: http: //hakank. org/common_cp_models/#minesweeper Answer Set Programming: minesweeper. lp Choco: Mine. Sweeper. java Comet: minesweeper. co ECLi. PSE: minesweeper. ecl Gecode/R: minesweeper. rb Google CP Solver/C#: minesweeper. cs Google CP Solver/Java: Minesweeper. java Google CP Solver/Python: minesweeper. py Ja. Co. P: Mine. Sweeper. java Ja. Co. P/Scala: Minesweeper. scala Mini. Zinc: minesweeper. mzn SICStus: minesweeper. pl Tailor/Essence': minesweeper. eprime Zinc: minesweeper. zinc

Minesweeper – search strategies (Mini. Zinc) To ensure the best/good performance, search strategies (heuristics) are used: - variable selection: which variable to use first_fail, input_order, most_constrained, etc - value selection: which value to try on the selected variable indomain_min, indomain_max, etc - the variables (as array) to work on - first_fail/indomain_min is often good to start with Selecting the best strategy is an art, not an science, and often requires much experimentation. solve : : int_search( [mines[i, j]|i in 1. . r, j in 1. . c], first_fail, indomain_min, complete) satisfy;

Minesweeper – search (or-tools/Python) - variable selection: INT_VAR_DEFAULT, CHOOSE_RANDOM, etc - value selection: INT_VALUE_DEFAULT, ASSIGN_MIN_VALUE - the variables, as array, to work on (here: mines) - INT_VAR_DEFAULT/INT_VALUE_DEFAULT is often good to start with db = solver. Phase([mines[(i, j)] for i in range(r) for j in range(c)], solver. INT_VAR_DEFAULT, solver. INT_VALUE_DEFAULT)

Key concepts in CP modelling Some of the most important - and IMHO the most fascinating concepts in CP regarding the modelling: - Global constraints - Element - Reification (logical relation) - Bi-directedness (multiple flow pattern) - Channelling (dual views) - 1/N/All solutions - Symmetry breaking

Global constraints - “Patterns” in modelling problems - Conceptual tool - Much research in effective propagators - Global Constraint Catalog, 364 global constraints http: //www. emn. fr/z-info/sdemasse/gccat/index. html Some of the most common global constraints: - All different: all elements must be distinct - All different except 0: all elements != 0 must be distinct - Element: indexing with decision variables - Global Cardinality Count: count occurrences of values in array - Decrease/Increase: sortedness - Regular: finite state machine - Cumulative: for scheduling - Circuit: Hamiltonian circuit - Table: Allowed assignments http: //hakank. org/minizinc/#global (~170 decompositions)

Element constraint (z = x[y]) One of the most common – and important - constraint: - z = x[y], where the index (y) is a decision variable Most CP solvers don't have support for this direct syntax. Instead: - z = element(x, y) - z = x. element(y) - element(x, y, z) % Mini. Zinc array[1. . 4] of var 1. . 4: x; var 1. . 4: y; var 1. . 4: z; constraint z = x[y]; % % example: x = [4, 3, 1, 2] y = 2 → z = 3

Element constraint (z = x[y]) // or-tools/C# solver. Add(z == x. Element(y)); % or-tools/Python solver. Add(z == solver. Element(x, y)) // or-tools/Java solver. add. Constraint( solver. make. Equality(z, solver. make. Element(x, y). var()));

Reification: truth values of constraints C => b If constraint C holds then Bool. Var b = 1 else b = 0. Some examples (Mini. Zinc syntax): x and y are var int (domain 1. . 10) b 1 and b 2 are var boolean (x = y) -> (b 1 = 1) (x > 4) -> (y <= 4) b 1 <-> not(b 2) (x = 1 / y = 1) <-> (b 1 = b 2)

Reification – alldifferent_except_0 (Mini. Zinc) The global constraint alldifferent_except_0 is quite useful: All variables that are != 0 must be distinct. Applications: - an Nx. N grid which should contain N distinct values, e. g. one per row and column. Other cells are set to 0. Here is a decomposition using reification in Mini. Zinc: forall(i, j in index_set(x) where i < j) ( (x[i] != 0 / x[j] != 0) -> x[i] != x[j] )

Reification – alldifferent_except_0 (or-tools/C#) In or-tools we have to use Boolean logic. public static void All. Different. Except 0(Sol s, Int. Var[] a) { int n = a. Length; for(int i = 0; i < n; i++) { for(int j = 0; j < i; j++) { // (a[i] != 0 / a[j] != 0) → (a[i] != a[j]) s. Add((a[i] != 0) * (a[j] != 0) <= (a[i] != a[j])); } } }

Reification - differences between or-tools' interfaces or-tools/Python: for i in range(n): for j in range(i): s. Add((a[i] != 0) * (a[j] != 0) <= (a[i] != a[j])) or-tools/Java: for(int i = 0; i < n; i++) { for(int j = 0; j < i; j++) { Int. Var bi = s. make. Is. Different. Cst. Var(a[i], 0); Int. Var bj = s. make. Is. Different. Cst. Var(a[j], 0); Int. Var bij = s. make. Is. Different. Cst. Var(a[i], a[j]); s. add. Constraint( s. make. Less. Or. Equal( s. make. Prod(bi, bj). var(), bij)); } } }

“Bi-direction” of variables Cf Prolog where a variable may be either input, output or both (“multiple flow pattern”). Simple example: - converting a number (variable) TO/FROM its digits (array) given a base [Later note: After the talk Nicolas Beldiceanu mentioned that this concept is better known as “reversibility” in Prolog. ]

“Bi-direction” of variables (or-tools/C#) - converting a number (variable) TO/FROM its digits (array) given a base private static Constraint To. Num(Int. Var[] a, % digits Int. Var num, % number int bbase) { int len = a. Length; Int. Var[] tmp = new Int. Var[len]; for(int i = 0; i < len; i++) { // represent a[i] with exponents in base tmp[i] = (a[i]*(int)Math. Pow(bbase, (len-i-1))). Var(); } % returns “success”/”fail” of the Constraint return tmp. Sum() == num; } // num = 532 a = {5, 3, 2} tmp = [5*100, 3*10, 2*1]

“Bi-direction” of variables (or-tools/C#) - converting a number (variable) TO/FROM its digits (array) given a base int n = 5; Int. Var[] x = solver. Make. Int. Var. Array(n, 0, 9); Int. Var n = solver. Make. Int. Var(0, 99999); Solver. Add(To. Num(x, n)); // FROM number TO array → x = {1, 3, 2, 5, 4} solver. Add(n == 13254); // OR // TO number FROM the digits -> n = 13254 solver. Add(x[0] == 1); solver. Add(x[1] == 3); solver. Add(x[2] == 2); solver. Add(x[3] == 5); solver. Add(x[4] == 4);

“Bi-direction” of variables – “dual view” - converting a number (variable) TO/FROM its digits (array) given a base We can now add constraints on either n or x: a “dual view”. Example: int n = 5; Int. Var[] x = solver. Make. Int. Var. Array(n, 0, 9); Int. Var n = solver. Make. Int. Var(12345, 99999); Solver. Add(To. Num(x, n)); Solver. Add(x. All. Different()); solver. Add(x[3] == 4); solver. Add(n >= 50000); The ability to add constraints quite easily to a model is one of the strengths of CP.

1/N/All solutions Most CP solvers can show any number of solutions. - First solution: e. g. just any valid solution (schedule/seating/etc). - N solutions, e. g. N=2 to check if/ensure that there are > 1 solutions. * Sudoku: must be unique solution (model is wrong if 2 solutions are printed) - All solutions: e. g. as debugging aid * N-queens for N=8 should yield 92 solutions

Symmetry breaking A common technique to reduce the search tree is to use SYMMETRY BREAKING. Some standard constraints for symmetry breaking: - increasing: the array must be sorted - lex_less(x, y): array x is lexicographic less then array y - lex 2(matrix): all rows and columns must be in lexicographic order Special forms: - x[0] < x[k-1]

So, why do I like Constraint Programming? - Excellent tool for certain type of problem solving - Modelling language is often high or very high - Declarative - Reification - Global constraints (efficiency, conceptual tool) - Bi-direction of predicates (“dual view”) - Can be very fast - Can show 1, 2, N, or all solutions But - like any programming paradigm - is no silver bullet. - It takes time to model a problem with much experiments, especially to get search strategies right: art not science - Debugging can be harder than in non-CP programming.

What have I not talked about? Advanced features in or-tools such as - Large Neighbourhood Search - Vehicle Routing - Min/Max-Cost Flow - Graph algorithms - etc. Advanced features in CP in general: - how search heuristics and propagation really works - set variables (can be very handy) - advanced examples - how a CP solver works under the hood

Thank you! - Questions? Thanks to Laurent Perron and the or-tools group for a great system! Hakan Kjellerstrand http: //www. hakank. org/index_eng. html