Tweet
Tweet

Slide 1

Slide 1 text

objectcomputing.com © 2018, Object Computing, Inc. (OCI). All rights reserved. No part of these notes may be reproduced, stored in a retrieval system, or transmitted, in any form or by any means, electronic, mechanical, photocopying, recording, or otherwise, without the prior, written permission of Object Computing, Inc. (OCI) An Introduction to Constraint Programming (with Groovy & other languages) Dr Paul King OCI Groovy Lead V.P. and PMC Chair Apache Groovy @paulk_asert © 2020 Object Computing, Inc. (OCI). All rights reserved. objectcomputing.com

Slide 2

Slide 2 text

Dr Paul King OCI Groovy Lead V.P. and PMC Chair Apache Groovy Slides: https://speakerdeck.com/paulk/groovy-constraint-programming Examples repo: https://github.com/paulk-asert/groovy-constraint-programming ReGinA author: https://www.manning.com/books/groovy-in-action-second-edition Twitter: @paulk_asert

Slide 3

Slide 3 text

Apache Groovy Programming Language • Multi-faceted language • Imperative/OO & functional (+ other styles via libraries) • Dynamic & static • Aligned closely with Java language • Very extensible • 17+ years since inception • ~700M downloads (partial count) • 300+ contributors • 200+ releases • https://www.youtube.com/watch?v=eIGOG- F9ZTw&feature=youtu.be

Slide 4

Slide 4 text

Friends of Apache Groovy Open Collective

Slide 5

Slide 5 text

Computer science algorithms Decision problems Search problems Counting problems Optimization problems • Paradigms: brute force, divide & conquer, search & enumerate, randomized, complexity reduction, recursive, back tracking, graph • Optimization approaches: linear programming, dynamic programming, greedy, heuristics

Slide 6

Slide 6 text

Optimization Overview Optimization: • Finding (optimal) solutions Approaches: • Dynamic programming • Linear programming • Non-linear programming • Constraint programming • Integer programming Aspects: • Combinatorial optimization • Continuous optimization Applications: • Timetabling/scheduling • Vehicle routing • Packing problems • Trading systems • Logistics • Language detection CP Constraint Programming CIP Constraint Integer Programming MINLP Mixed Integer Non-linear Programming MIP Mixed Integer Linear Programming IP Integer Linear Programming SAT Satisfiability LP Linear Programming

Slide 7

Slide 7 text

Constraint/Logic Programming Description • Declarative style of expressing problem • Model: variables + domains + constraints • Solver: does all the work

Slide 8

Slide 8 text

Constraint/Logic Programming Description • Declarative style of expressing problem • Express properties of solution not steps • Model: variables + domains + constraints • Solver: does all the work • Aspects: heuristic search, inference/propagation, symmetry, backtracking • Origins: popularized in logic programming languages such as Prolog but now typically used outside logic programming specific languages • Numerous libraries we’ll look at: Choco, ojalgo and others

Slide 9

Slide 9 text

Constraint Programming Picture: http://osullivan.ucc.ie/AAAI_OSullivan.pdf

Slide 10

Slide 10 text

Constraint Programming Typical domains: • Boolean (SAT problem) • integer • rational • interval (scheduling problems) • linear (approaches to non- linear problems do exist) • finite sets • mixed (2 or more of above) Picture: http://osullivan.ucc.ie/AAAI_OSullivan.pdf

Slide 11

Slide 11 text

Pythagorean triples Three positive integers a, b, c, such that: a2 + b2 = c2 Form right-angled triangles a c b a b c 3 4 5 5 12 13 7 24 25 8 15 17 …

Slide 12

Slide 12 text

Pythagorean triples // imperative for (x in 1..30) for (y in (x+1)..30) // y > x to remove duplicates for (z in (y+1)..30) // z > y for efficiency if (x*x + y*y == z*z) println "($x, $y, $z)" // constraint-based new Model('Pythagorean-triple').with { def x = intVar(1, 30) def y = intVar(1, 30) def z = intVar(1, 30) y.gt(x).post() // y > x to remove duplicates z.gt(y).post() // z > y for efficiency x.mul(x).add(y.mul(y)).eq(z.mul(z)).post() while (solver.solve()) { println "[$x.value, $y.value, $z.value]" } } (3, 4, 5) (5, 12, 13) (6, 8, 10) (7, 24, 25) (8, 15, 17) (9, 12, 15) (10, 24, 26) (12, 16, 20) (15, 20, 25) (18, 24, 30) (20, 21, 29) [3, 4, 5] [6, 8, 10] [5, 12, 13] [9, 12, 15] [8, 15, 17] [7, 24, 25] [12, 16, 20] [10, 24, 26] [15, 20, 25] [20, 21, 29] [18, 24, 30]

Slide 13

Slide 13 text

Pythagorean triples // imperative for (x in 1..30) for (y in (x+1)..30) // y > x to remove duplicates for (z in (y+1)..30) // z > y for efficiency if (x*x + y*y == z*z) println "($x, $y, $z)"

Slide 14

Slide 14 text

Pythagorean triples // imperative for (x in 1..30) for (y in (x+1)..30) // y > x to remove duplicates for (z in (y+1)..30) // z > y for efficiency if (x*x + y*y == z*z) println "($x, $y, $z)" Code explicitly chooses order in which to traverse integer variables

Slide 15

Slide 15 text

Pythagorean triples // imperative for (x in 1..30) for (y in (x+1)..30) // y > x to remove duplicates for (z in (y+1)..30) // z > y for efficiency if (x*x + y*y == z*z) println "($x, $y, $z)" Code has explicitly chosen order for variable values

Slide 16

Slide 16 text

Pythagorean triples // imperative for (x in 1..30) for (y in (x+1)..30) // y > x to remove duplicates for (z in (y+1)..30) // z > y for efficiency if (x*x + y*y == z*z) println "($x, $y, $z)" Manual testing of desired condition (3, 4, 5) (5, 12, 13) (6, 8, 10) (7, 24, 25) (8, 15, 17) (9, 12, 15) (10, 24, 26) (12, 16, 20) (15, 20, 25) (18, 24, 30) (20, 21, 29)

Slide 17

Slide 17 text

Pythagorean triples // constraint-based new Model('Pythagorean-triple').with { def x = intVar(1, 30) def y = intVar(1, 30) def z = intVar(1, 30) y.gt(x).post() // y > x to remove duplicates z.gt(y).post() // z > y for efficiency x.mul(x).add(y.mul(y)).eq(z.mul(z)).post() while (solver.solve()) { println "[$x.value, $y.value, $z.value]" } }

Slide 18

Slide 18 text

Pythagorean triples // constraint-based new Model('Pythagorean-triple').with { def x = intVar(1, 30) def y = intVar(1, 30) def z = intVar(1, 30) y.gt(x).post() // y > x to remove duplicates z.gt(y).post() // z > y for efficiency x.mul(x).add(y.mul(y)).eq(z.mul(z)).post() while (solver.solve()) { println "[$x.value, $y.value, $z.value]" } } Declare variables as unassigned “integer black boxes” with a given range

Slide 19

Slide 19 text

Pythagorean triples // constraint-based new Model('Pythagorean-triple').with { def x = intVar(1, 30) def y = intVar(1, 30) def z = intVar(1, 30) y.gt(x).post() // y > x to remove duplicates z.gt(y).post() // z > y for efficiency x.mul(x).add(y.mul(y)).eq(z.mul(z)).post() while (solver.solve()) { println "[$x.value, $y.value, $z.value]" } } Other relationships between the variables are expressed as constraints

Slide 20

Slide 20 text

Pythagorean triples // constraint-based new Model('Pythagorean-triple').with { def x = intVar(1, 30) def y = intVar(1, 30) def z = intVar(1, 30) y.gt(x).post() // y > x to remove duplicates z.gt(y).post() // z > y for efficiency x.mul(x).add(y.mul(y)).eq(z.mul(z)).post() while (solver.solve()) { println "[$x.value, $y.value, $z.value]" } } Invoke solver. It determines ordering and approach to use. [3, 4, 5] [6, 8, 10] [5, 12, 13] [9, 12, 15] [8, 15, 17] [7, 24, 25] [12, 16, 20] [10, 24, 26] [15, 20, 25] [20, 21, 29] [18, 24, 30]

Slide 21

Slide 21 text

Cryptarithmetic S E N D + M O R E = M O N E Y Replace the letters with decimal digits to make a valid arithmetic sum

Slide 22

Slide 22 text

Cryptarithmetic: Solving by hand S E N D + M O R E = M O N E Y Col 5 4 3 2 1 From column 5, M = 1 since it is the only carry- over possible from the sum of two single digit numbers in column 4 Inspired by: https://en.wikipedia.org/wiki/Verbal_arithmetic

Slide 23

Slide 23 text

Cryptarithmetic: Solving by hand S E N D + 1 O R E = 1 O N E Y Col 5 4 3 2 1 So far: M = 1 Since there is a carry in column 5, O must be less than or equal to M (from column 4). But O cannot be equal to M, so O is less than M. Therefore O = 0.

Slide 24

Slide 24 text

Cryptarithmetic: Solving by hand S E N D + 1 0 R E = 1 0 N E Y Col 5 4 3 2 1 So far: M = 1, O = 0 Since O is 1 less than M, S is either 8 or 9 depending on whether there is a carry in column 4. But if there were a carry in column 4, N would be less than or equal to O (from column 3). This is impossible since O = 0. Therefore there is no carry in column 3 and S = 9.

Slide 25

Slide 25 text

Cryptarithmetic: Solving by hand 9 E N D + 1 0 R E = 1 0 N E Y Col 5 4 3 2 1 So far: M = 1, O = 0, S = 9 If there were no carry in column 3 then E = N, which is impossible. Therefore there is a carry and N = E + 1.

Slide 26

Slide 26 text

Cryptarithmetic: Solving by hand 9 E N D + 1 0 R E = 1 0 N E Y Col 5 4 3 2 1 So far: M = 1, O = 0, S = 9 N = E + 1 If there were no carry in column 2, then ( N + R ) mod 10 = E, and N = E + 1, so ( E + 1 + R ) mod 10 = E which means ( 1 + R ) mod 10 = 0, so R = 9. But S = 9, so there must be a carry in column 2 so R = 8

Slide 27

Slide 27 text

Cryptarithmetic: Solving by hand 9 E N D + 1 0 8 E = 1 0 N E Y Col 5 4 3 2 1 So far: M = 1, O = 0, S = 9, R = 8 N = E + 1 To produce a carry in column 2, we must have D + E = 10 + Y Y is at least 2 so D + E is at least 12

Slide 28

Slide 28 text

Cryptarithmetic: Solving by hand 9 E N D + 1 0 8 E = 1 0 N E Y Col 5 4 3 2 1 So far: M = 1, O = 0, S = 9, R = 8 N = E + 1, D + E >= 12 The only two pairs of available numbers that sum to at least 12 are (5,7) and (6,7) so either E = 7 or D = 7

Slide 29

Slide 29 text

Cryptarithmetic: Solving by hand 9 E N D + 1 0 8 E = 1 0 N E Y Col 5 4 3 2 1 So far: M = 1, O = 0, S = 9, R = 8 N = E + 1, D + E >= 12 Since N = E + 1, E can't be 7 because then N = 8 = R so D = 7

Slide 30

Slide 30 text

Cryptarithmetic: Solving by hand 9 E N 7 + 1 0 8 E = 1 0 N E Y Col 5 4 3 2 1 So far: M = 1, O = 0, S = 9, R = 8, D = 7 N = E + 1, D + E >= 12 E can't be 6 because then N = 7 = D so E = 5 and N = 6

Slide 31

Slide 31 text

Cryptarithmetic: Solving by hand 9 5 6 7 + 1 0 8 5 = 1 0 6 5 Y Col 5 4 3 2 1 So far: M = 1, O = 0, S = 9, R = 8 D = 7, E = 5, N = 6 D + E = 12 so Y = 2

Slide 32

Slide 32 text

Cryptarithmetic: Solving by hand 9 5 6 7 + 1 0 8 5 = 1 0 6 5 2 Solved! M = 1, O = 0, S = 9, R = 8 D = 7, E = 5, N = 6, Y = 2

Slide 33

Slide 33 text

Cryptarithmetic: Solving by hand 9 5 6 7 + 1 0 8 5 = 1 0 6 5 2 Solved! M = 1, O = 0, S = 9, R = 8 D = 7, E = 5, N = 6, Y = 2 S P E N D + L E S S = M O N E Y No solution!

Slide 34

Slide 34 text

Cryptarithmetic: Brute force def solutions(): # letters = ('s', 'e', 'n', 'd', 'm', 'o', 'r', 'y') all_solutions = list() for s in range(1, 10): for e in range(0, 10): for n in range(0, 10): for d in range(0, 10): for m in range(1, 10): for o in range(0, 10): for r in range(0, 10): for y in range(0, 10): if len({s, e, n, d, m, o, r, y}) == 8: send = 1000 * s + 100 * e + 10 * n + d more = 1000 * m + 100 * o + 10 * r + e money = 10000 * m + 1000 * o + 100 * n + 10 * e + y if send + more == money: all_solutions.append((send, more, money)) return all_solutions print(solutions()) S E N D + M O R E = M O N E Y [(9567, 1085, 10652)]

Slide 35

Slide 35 text

Cryptarithmetic: Brute force def solutions() { // letters = ['s', 'e', 'n', 'd', 'm', 'o', 'r', 'y'] def all_solutions = [] for (s in 1..<10) for (e in 0..9) for (n in 0..9) for (d in 0..9) for (m in 1..9) for (o in 0..9) for (r in 0..9) for (y in 0..9) if ([s, e, n, d, m, o, r, y].toSet().size() == 8) { def send = 1000 * s + 100 * e + 10 * n + d def more = 1000 * m + 100 * o + 10 * r + e def money = 10000 * m + 1000 * o + 100 * n + 10 * e + y if (send + more == money) all_solutions.add([send, more, money]) } return all_solutions } print(solutions()) S E N D + M O R E = M O N E Y [[9567, 1085, 10652]]

Slide 36

Slide 36 text

Cryptarithmetic: Brute force from itertools import permutations def solution2(): letters = ('s', 'e', 'n', 'd', 'm', 'o', 'r', 'y') digits = range(10) for perm in permutations(digits, len(letters)): sol = dict(zip(letters, perm)) if sol['s'] == 0 or sol['m'] == 0: continue send = 1000 * sol['s'] + 100 * sol['e'] + 10 * sol['n'] + sol['d'] more = 1000 * sol['m'] + 100 * sol['o'] + 10 * sol['r'] + sol['e'] money = 10000 * sol['m'] + 1000 * sol['o'] + 100 * sol['n'] + 10 * sol['e'] + sol['y'] if send + more == money: return send, more, money print(solution2()) S E N D + M O R E = M O N E Y (9567, 1085, 10652)

Slide 37

Slide 37 text

Cryptarithmetic: Brute force def solution2() { def digits = 0..9 for (p in digits.permutations()) { if (p[-1] < p[-2]) continue def (s, e, n, d, m, o, r, y) = p if (s == 0 || m == 0) continue def send = 1000 * s + 100 * e + 10 * n + d def more = 1000 * m + 100 * o + 10 * r + e def money = 10000 * m + 1000 * o + 100 * n + 10 * e + y if (send + more == money) return [send, more, money] } } print(solution2()) S E N D + M O R E = M O N E Y [9567, 1085, 10652]

Slide 38

Slide 38 text

from csp import Constraint, CSP from typing import Dict, List, Optional class SendMoreMoneyConstraint(Constraint[str, int]): def __init__(self, letters: List[str]) -> None: super().__init__(letters) self.letters: List[str] = letters def satisfied(self, assignment: Dict[str, int]) -> bool: # not a solution if duplicate values if len(set(assignment.values())) < len(assignment): return False # if all vars assigned, check if correct if len(assignment) == len(self.letters): s: int = assignment["S"] e: int = assignment["E"] n: int = assignment["N"] d: int = assignment["D"] m: int = assignment["M"] o: int = assignment["O"] r: int = assignment["R"] y: int = assignment["Y"] send: int = s * 1000 + e * 100 + n * 10 + d more: int = m * 1000 + o * 100 + r * 10 + e money: int = m * 10000 + o * 1000 + n * 100 + e * 10 + y return send + more == money return True … Cryptarithmetic: Constraint programming … letters = ["S", "E", "N", "D", "M", "O", "R", "Y"] possible_digits = {} for letter in letters: possible_digits[letter] = [0, 1, 2, 3, 4, 5, 6, 7, 8, 9] possible_digits["M"] = [1] # can't start with 0 constraint csp: CSP[str, int] = CSP(letters, possible_digits) csp.add_constraint(SendMoreMoneyConstraint(letters)) solution: Optional[Dict[str, int]] = csp.backtracking_search() if solution is None: print("No solution found!") else: print(solution) S E N D + M O R E = M O N E Y {'S': 9, 'E': 5, 'N': 6, 'D': 7, 'M': 1, 'O': 0, 'R': 8, 'Y': 2} Adapted from: https://freecontent.manning.com/constraint-satisfaction-problems-in-python/

Slide 39

Slide 39 text

@Grab('org.choco-solver:choco-solver:4.10.5') import org.chocosolver.solver.Model import org.chocosolver.solver.variables.IntVar def model = new Model("SEND+MORE=MONEY") def S = model.intVar("S", 1, 9) def E = model.intVar("E", 0, 9) def N = model.intVar("N", 0, 9) def D = model.intVar("D", 0, 9) def M = model.intVar("M", 1, 9) def O = model.intVar("O", 0, 9) def R = model.intVar("R", 0, 9) def Y = model.intVar("Y", 0, 9) model.allDifferent(S, E, N, D, M, O, R, Y).post() … Cryptarithmetic: Constraint programming … IntVar[] ALL = [ S, E, N, D, M, O, R, E, M, O, N, E, Y] int[] COEFFS = [ 1000, 100, 10, 1, 1000, 100, 10, 1, -10000, -1000, -100, -10, -1] model.scalar(ALL, COEFFS, "=", 0).post() model.solver.findSolution() S E N D + M O R E = M O N E Y Solution: S=9, E=5, N=6, D=7, M=1, O=0, R=8, Y=2,

Slide 40

Slide 40 text

Global Constraint Catalog https://sofdem.github.io/gccat/gccat/sec5.html

Slide 41

Slide 41 text

Constraint Programming Solver Search Constraint & Variable Store For a particular search strategy For each variable For each value in the domain of the current variable Check if all constraints are satisfied May involve branching, pruning, propagation, arc/path consistency … Constraints have dedicated feasibility/pruning algorithms Propose new constraints Check satisfiability

Slide 42

Slide 42 text

SEND + MORE = MONEY def (S, M) = ['S', 'M'].collect { intVar(it, 1..9) } def (E, N, D, O, R, Y) = ['E', 'N', 'D', 'O', 'R', 'Y'].collect { intVar(it, 0..9) } def C = intVarArray(4, 0..1) // C3 C2 C1 C0 // S E N D // M O R E // ------------- // M O N E Y

Slide 43

Slide 43 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10)

Slide 44

Slide 44 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10) Start by pruning Consider the range of both sides

Slide 45

Slide 45 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10) range(LHS) = 0..1 range(RHS) = 1..9

Slide 46

Slide 46 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10) 1 is the only common value: C[3] = 1 M = 1

Slide 47

Slide 47 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10)

Slide 48

Slide 48 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10) Consider range of both sides to see if a solution is feasible

Slide 49

Slide 49 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10) range(LHS) = range(C[2] + S + M) = range(C[2] + S + 1) = minRange(C[2] + S + 1) .. maxRange(C[2] + S + 1) = 0 + 2 + 1 .. 1 + 9 + 1 = 3 .. 11

Slide 50

Slide 50 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10) range(LHS) = 3..11 range(RHS) = range(O + C[3] * 10) = range(O + 1 * 10) = 0 + 10 .. 9 + 10 = 10 .. 19

Slide 51

Slide 51 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10) range(LHS) = 3..11 range(RHS) = 10..19 overlap = 10..11

Slide 52

Slide 52 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10) Applying new information to LHS 10 <= C[2] + S + 1 <= 11 9 <= C[2] + S <= 10 9 – C[2] <= S <= 10 – C[2] 9 - maxRange(C[2]) <= S <= 10 – minRange(C[2]) 9 – 1 <= S <= 10 – 0 S in 8..10

Slide 53

Slide 53 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10) Applying new information to RHS 10 <= O + 1 * 10 <= 11 0 <= O <= 1 O in 0..1 but 1 is already taken O = 0

Slide 54

Slide 54 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10) Applying new information to RHS 10 <= O + 1 * 10 <= 11 0 <= O <= 1 O in 0..1 but 1 is already taken O = 0

Slide 55

Slide 55 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10) Applying new information to RHS 10 <= O + 1 * 10 <= 11 0 <= O <= 1 O in 0..1 but 1 is already taken O = 0

Slide 56

Slide 56 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10) Reapplying allDifferent

Slide 57

Slide 57 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10) C[1] + E + O = N + C[2] * 10 range(LHS) = 2..10 Applying to RHS 2 <= N + C[2] * 10 <= 10 2 – N <= C[2] * 10 <= 10 - N 2 – 9 <= C[2] * 10 <= 10 – 2 -7 <= C[2] * 10 <= 8 C[2] = 0

Slide 58

Slide 58 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10) Go back to previous constraint now we know C[2] S + 1 = 10 S = 9

Slide 59

Slide 59 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10) And apply allDifferent as previously. Try branching: Guess C[1] = 0

Slide 60

Slide 60 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10) Substitute C[1] = 0 into above constraint: C[1] + E + O = N + C[2] * 10 0 + E + 0 = N + 0 * 10 E = N which can’t be true because of allDifferent, so this isn’t possible and C[1] must be 1

Slide 61

Slide 61 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10) Back to pruning Reconsider the above constraint

Slide 62

Slide 62 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10) Now we have: 1 + E = N, so N > E and we can prune 2 more

Slide 63

Slide 63 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10) C[0] + N + R = E + 10 R = E + 10 – N – C[0] range(RHS) = 2+10-8-1 .. 7+10-3-0 = 3..14

Slide 64

Slide 64 text

SEND + MORE = MONEY allDifferent(S, E, N, D, M, O, R, Y) // C3 C2 C1 C0 C[3] .eq(M) // S E N D (C[2] + S + M).eq(O + C[3] * 10) // M O R E (C[1] + E + O).eq(N + C[2] * 10) // ------------- (C[0] + N + R).eq(E + C[1] * 10) // M O N E Y (D + E).eq(Y + C[0] * 10) Now back to branching …

Slide 65

Slide 65 text

SEND + MORE = MONEY Graphs generated using cpviz and choco-cpviz

Slide 66

Slide 66 text

SEND + MORE = MONEY (Other JVM language examples) Prolog

Slide 67

Slide 67 text

Frobenius numbers Given positive integers a 1 , a 2 , …, a n , find integers that cannot be expressed as a sum: k 1 a 1 + k 2 a 2 + ··· + k n a n where k 1 , k 2 , …, k n are non-negative integers. Largest such number is the Frobenius number for those integers.

Slide 68

Slide 68 text

Frobenius numbers Given positive integers a 1 , a 2 , …, a n , find integers that cannot be expressed as a sum: k 1 a 1 + k 2 a 2 + ··· + k n a n where k 1 , k 2 , …, k n are non-negative integers. Largest such number is the Frobenius number for those integers. Frobenius number for {2, 5} is 3 https://en.wikipedia.org/wiki/Coin_problem

Slide 69

Slide 69 text

Frobenius numbers Given positive integers a 1 , a 2 , …, a n , find integers that cannot be expressed as a sum: k 1 a 1 + k 2 a 2 + ··· + k n a n where k 1 , k 2 , …, k n are non-negative integers. Largest such number is the Frobenius number for those integers. Applications: • Numerical semigroups/monoids in algebraic geometry • Rugby/NFL: impossible scores • Coin/stamp denominations (knapsack-like problem)

Slide 70

Slide 70 text

Frobenius numbers Given positive integers a 1 , a 2 , …, a n , find integers that cannot be expressed as a sum: k 1 a 1 + k 2 a 2 + ··· + k n a n where k 1 , k 2 , …, k n are non-negative integers. Largest such number is the Frobenius number for those integers. Applications: • Numerical semigroups/monoids in algebraic geometry • Rugby/NFL: impossible scores • Coin/stamp denominations (knapsack-like problem) • McNuggets challenge 6 9 20

Slide 71

Slide 71 text

McNuggets challenge A fast-food restaurant sells chicken nuggets in orders of 6, 9, and 20. What is the largest nugget order that cannot be served exactly from this restaurant? 6 9 20

Slide 72

Slide 72 text

McNuggets challenge (Choco) @Grab('org.choco-solver:choco-solver:4.10.5') import org.chocosolver.solver.Model import org.chocosolver.solver.variables.IntVar def limit = 100 int[] coeffs = [6, 9, 20] for (int i : limit..1) { def model = new Model("McNuggets challenge") def packs6 = model.intVar("#6-packs", 0, limit.intdiv(6)) def packs9 = model.intVar("#9-packs", 0, limit.intdiv(9)) def packs20 = model.intVar("#20-packs", 0, limit.intdiv(20)) IntVar[] all = [packs6, packs9, packs20] model.scalar(all, coeffs, "=", i).post() def found = model.solver.findSolution() if (found) { println "$i: $found" } else { println "Not possible to order $i nuggets" break } } … 49: Solution: #6-packs=0, #9-packs=1, #20-packs=2, 48: Solution: #6-packs=8, #9-packs=0, #20-packs=0, 47: Solution: #6-packs=3, #9-packs=1, #20-packs=1, 46: Solution: #6-packs=1, #9-packs=0, #20-packs=2, 45: Solution: #6-packs=6, #9-packs=1, #20-packs=0, 44: Solution: #6-packs=4, #9-packs=0, #20-packs=1, Not possible to order 43 nuggets

Slide 73

Slide 73 text

McNuggets challenge (ojalgo) @Grab('org.ojalgo:ojalgo:48.3.1') import org.ojalgo.optimisation.ExpressionsBasedModel import static org.ojalgo.netio.BasicLogger.debug def model = new ExpressionsBasedModel() def packs6 = model.addVariable('#6-packs').lower(0).integer(true) def packs9 = model.addVariable('#9-packs').lower(0).integer(true) def packs20 = model.addVariable('#20-packs').lower(0).integer(true) def totalNuggets = model.addExpression().weight(1) totalNuggets.set(packs6, 6) totalNuggets.set(packs9, 9) totalNuggets.set(packs20, 20) for (int i : 100..1) { totalNuggets.upper(i) def result = model.maximise() if (Math.round(result.value) < i) { debug("Not possible to order $i nuggets") break } else { debug(result) } } // See: https://www.ojalgo.org/2019/05/the-mcnuggets-challenge/ … OPTIMAL 50.0 @ { 2, 2, 1 } OPTIMAL 49.00000000000001 @ { 0.0, 1.0000000000000004, 2.0 } OPTIMAL 48.0 @ { 8, 0, 0 } OPTIMAL 47.00000000000001 @ { 0.0, 3.0000000000000004, 1.0 } OPTIMAL 45.99999999999999 @ { 0.9999999999999994, 0.0, 2.0 } OPTIMAL 45.0 @ { 0.0, 5.0, 0.0 } OPTIMAL 44.0 @ { 4.0, 0.0, 1.0 } Not possible to order 43 nuggets Alternative approach, rather than trying to solve exactly with Integer variables, find maximal solution with Real variables and check if that was an exact match.

Slide 74

Slide 74 text

Dietary restrictions Minimise cost of diet given: • Must be at least 300 calories • Not more than 10 grams of protein • Not less than 10 grams of carbohydrates • Not less than 8 grams of fat • At least 0.5 units of fish • No more than 1 unit of milk Bread Milk Cheese Potato Fish Yogurt Cost 2.0 3.5 8.0 1.5 11.0 1.0 Protein (g) 4.0 8.0 7.0 1.3 8.0 9.2 Fat (g) 1.0 5.0 9.0 0.1 7.0 1.0 Carbohydrates (g) 15.0 11.7 0.4 22.6 0.0 17.0 Calories 90 120 106 97 130 180 http://support.sas.com/documentation/cdl/en/ormpug/63352/HTML/default/viewer.htm#ormpug_lpsolver_sect018.htm

Slide 75

Slide 75 text

Dietary restrictions

Slide 76

Slide 76 text

Dietary restrictions (Excel)

Slide 77

Slide 77 text

Dietary restrictions (Excel)

Slide 78

Slide 78 text

Dietary restrictions (Google sheets)

Slide 79

Slide 79 text

Dietary restrictions (Integer Linear Programming theory) Image source: https://www.researchgate.net/figure/Example-Integer-Linear-Program-ILP_fig1_215780881

Slide 80

Slide 80 text

Diet problem (ojalgo) def model = new ExpressionsBasedModel() def bread = model.addVariable("Bread").lower(0) def milk = model.addVariable("Milk").lower(0).upper(1) def cheese = model.addVariable("Cheese").lower(0) def potato = model.addVariable("Potato").lower(0) def fish = model.addVariable("Fish").lower(0.5) def yogurt = model.addVariable("Yogurt").lower(0) def cost = model.addExpression("Cost") cost.set(bread, 2.0).set(milk, 3.5).set(cheese, 8.0).set(potato, 1.5).set(fish, 11.0).set(yogurt, 1.0) def protein = model.addExpression("Protein").upper(10) protein.set(bread, 4.0).set(milk, 8.0).set(cheese, 7.0).set(potato, 1.3) .set(fish, 8.0).set(yogurt, 9.2) def fat = model.addExpression("Fat").lower(8) fat.set(bread, 1.0).set(milk, 5.0).set(cheese, 9.0).set(potato, 0.1) .set(fish, 7.0).set(yogurt, 1.0) def carbs = model.addExpression("Carbohydrates").lower(10) carbs.set(bread, 15.0).set(milk, 11.7).set(cheese, 0.4).set(potato, 22.6) .set(fish, 0.0).set(yogurt, 17.0) def calories = model.addExpression("Calories").lower(300) calories.set(bread, 90).set(milk, 120).set(cheese, 106).set(potato, 97) .set(fish, 130).set(yogurt, 180) cost.weight(1.0) def result = model.minimise() // for a variation, see: // https://www.ojalgo.org/2019/05/the-diet-problem/ OPTIMAL 12.08133788110815 @ { 0, 0.05359877488515, 0.44949881665042, 1.86516775720451, 0.5, 0 } ############################################ 0 <= Bread: 0 (2) 0 <= Milk: 0.053599 (3.5) <= 1 0 <= Cheese: 0.449499 (8) 0 <= Potato: 1.865168 (1.5) 0.5 <= Fish: 0.5 (11) 0 <= Yogurt: 0 (1) 10 <= Carbohydrates: 42.959697 8 <= Fat: 8.0 Cost: 12.081338 Protein: 10.0 <= 10 300 <= Calories: 300.0 ############################################

Slide 81

Slide 81 text

Diet problem (Choco) def model = new Model("Diet problem") def unbounded = 100000 // scale quantities by 10, coefficients by 10, products by 100 def bread = model.intVar("Bread", 0, unbounded) def milk = model.intVar("Milk", 0, 10) def cheese = model.intVar("Cheese", 0, unbounded) def potato = model.intVar("Potato", 0, unbounded) def fish = model.intVar("Fish", 5, unbounded) def yogurt = model.intVar("Yogurt", 0, unbounded) IntVar[] all = [bread, milk, cheese, potato, fish, yogurt] def cost = model.intVar("Cost", 0, unbounded) model.scalar(all, [20, 35, 80, 15, 110, 10] as int[], "=", cost).post() def protein = model.intVar("Protein", 0, 1000) model.scalar(all, [40, 80, 70, 13, 80, 92] as int[], "=", protein).post() def fat = model.intVar("Fat", 800, unbounded) model.scalar(all, [10, 50, 90, 1, 70, 10] as int[], "=", fat).post() def carbs = model.intVar("Carbohydrates", 1000, unbounded) model.scalar(all, [150, 117, 4, 226, 0, 170] as int[], "=", carbs).post() def calories = model.intVar("Calories", 30000, unbounded) model.scalar(all, [900, 1200, 1060, 970, 1300, 1800] as int[], "=", calories).post() model.setObjective(Model.MINIMIZE, cost) def found = model.solver.findSolution() if (found) { all.each { println "$it.name: ${it.value / 10}" } [carbs, fat, protein, calories, cost].each { println "$it.name: ${it.value / 100}" } } else { println "No solution" } Bread: 0 Milk: 0 Cheese: 0.5 Potato: 1.9 Fish: 0.5 Yogurt: 0 Carbohydrates: 43.14 Fat: 8.19 Protein: 9.97 Calories: 302.3 Cost: 12.35 Choco supports RealVar but doesn’t have as rich a set of possible constraints for such vars.

Slide 82

Slide 82 text

Diet problem (Choco) def model = new Model("Diet problem") def unbounded = 1000.0d def precision = 0.00001d // scale quantities by 10, coefficients by 10, products by 100 def bread = model.realVar("Bread", 0.0, unbounded, precision) def milk = model.realVar("Milk", 0.0, 1.0, precision) def cheese = model.realVar("Cheese", 0.0, unbounded, precision) def potato = model.realVar("Potato", 0.0, unbounded, precision) def fish = model.realVar("Fish", 0.5, unbounded, precision) def yogurt = model.realVar("Yogurt", 0.0, unbounded, precision) RealVar[] all = [bread, milk, cheese, potato, fish, yogurt] def scalarIbex = { coeffs, var -> def (a, b, c, d, e, f) = coeffs model.realIbexGenericConstraint("$a*{0}+$b*{1}+$c*{2}+$d*{3}+$e*{4}+$f*{5}={6}", [*all, var] as RealVar[]).post(); } def cost = model.realVar("Cost", 0.0, unbounded, precision) scalarIbex([2.0, 3.5, 8.0, 1.5, 11.0, 1.0], cost) def protein = model.realVar("Protein", 0.0, 10.0, precision) scalarIbex([4.0, 8.0, 7.0, 1.3, 8.0, 9.2], protein) def fat = model.realVar("Fat", 8.0, unbounded, precision) scalarIbex([1.0, 5.0, 9.0, 0.1, 7.0, 1.0], fat) def carbs = model.realVar("Carbohydrates", 10.0, unbounded, precision) scalarIbex([15.0, 11.7, 0.4, 22.6, 0.0, 17.0], carbs) def calories = model.realVar("Calories", 300, unbounded, precision) scalarIbex([90, 120, 106, 97, 130, 180], calories) model.setObjective(Model.MINIMIZE, cost) def found = model.solver.findSolution() Bread: 0.025131 .. 0.025137 Milk: 0.000009 .. 0.000010 Cheese: 0.428571 .. 0.428571 Potato: 1.848118 .. 1.848124 Fish: 0.561836 .. 0.561836 Yogurt: 0.000007 .. 0.000010 Carbohydrates: 42.316203 .. 42.316211 Fat: 8.000000 .. 8.000005 Protein: 9.997920 .. 9.997926 Calories: 300.000000 .. 300.000008 Cost: 12.431241 .. 12.431245 Choco does have a plugin (via JNI) for the Ibex C++ constraint processing library which does handle real numbers. def pretty = { var -> def bounds = found.getRealBounds(var) printf "%s: %.6f .. %.6f%n", var.name, *bounds } if (found) { all.each { pretty(it) } [carbs, fat, protein, calories, cost].each { pretty(it) } } else { println "No solution" }

Slide 83

Slide 83 text

Diet problem (Apache Commons Math) import org.apache.commons.math3.optim.linear.* import org.apache.commons.math3.optim.nonlinear.scalar.GoalType import static org.apache.commons.math3.optim.linear.Relationship.* def cost = new LinearObjectiveFunction([2.0, 3.5, 8.0, 1.5, 11.0, 1.0] as double[], 0) static scalar(coeffs, rel, val) { new LinearConstraint(coeffs as double[], rel, val) } def bread_min = scalar([1, 0, 0, 0, 0, 0], GEQ, 0) def milk_min = scalar([0, 1, 0, 0, 0, 0], GEQ, 0) def milk_max = scalar([0, 1, 0, 0, 0, 0], LEQ, 1) def cheese_min = scalar([0, 0, 1, 0, 0, 0], GEQ, 0) def potato_min = scalar([0, 0, 0, 1, 0, 0], GEQ, 0) def fish_min = scalar([0, 0, 0, 0, 1, 0], GEQ, 0.5) def yogurt_min = scalar([0, 0, 0, 0, 0, 1], GEQ, 0) def protein = scalar([4.0, 8.0, 7.0, 1.3, 8.0, 9.2], LEQ, 10) def fat = scalar([1.0, 5.0, 9.0, 0.1, 7.0, 1.0], GEQ, 8) def carbs = scalar([15.0, 11.7, 0.4, 22.6, 0.0, 17.0], GEQ, 10) def calories = scalar([90, 120, 106, 97, 130, 180], GEQ, 300) LinearConstraintSet constraints = [bread_min, milk_min, milk_max, fish_min, cheese_min, potato_min, yogurt_min, protein, fat, carbs, calories] def solution = new SimplexSolver().optimize(cost, constraints, GoalType.MAXIMIZE) if (solution != null) { println "Opt: $solution.value" println solution.point.collect{ sprintf '%.2f', it }.join(', ') } Opt: 12.553783912422674 -0.00, -0.00, 0.41, 1.85, 0.59, 0.00

Slide 84

Slide 84 text

import org.jacop.constraints.LinearInt import org.jacop.core.* import org.jacop.search.* import static org.jacop.core.IntDomain.MaxInt def store = new Store() // scale quantities by 10, coefficients by 10, products by 100 def bread = new IntVar(store, 'bread', 0, MaxInt) def milk = new IntVar(store, 'milk', 0, 10) def cheese = new IntVar(store, 'cheese', 0, MaxInt) def potato = new IntVar(store, 'potato', 0, MaxInt) def fish = new IntVar(store, 'fish', 5, MaxInt) def yogurt = new IntVar(store, 'yogurt', 0, MaxInt) IntVar[] all = [bread, milk, cheese, potato, fish, yogurt] def cost = new IntVar(store, 'Cost', 0, MaxInt) store.impose(new LinearInt(all, [20, 35, 80, 15, 110, 10] as int[], '=', cost)) def protein = new IntVar(store, 'Protein', 0, 1000) store.impose(new LinearInt(all, [40, 80, 70, 13, 80, 92] as int[], '=', protein)) def fat = new IntVar(store, 'Fat', 800, MaxInt) store.impose(new LinearInt(all, [10, 50, 90, 1, 70, 10] as int[], '=', fat)) def carbs = new IntVar(store, 'Carbohydrates', 1000, MaxInt) store.impose(new LinearInt(all, [150, 117, 4, 226, 0, 170] as int[], '=', carbs)) def calories = new IntVar(store, 'Calories', 30000, MaxInt) store.impose(new LinearInt(all, [900, 1200, 1060, 970, 1300, 1800] as int[], '=', calories)) def search = new DepthFirstSearch() def select = new SimpleSelect(all, null, new IndomainMin()) def result = search.labeling(store, select, cost) println result if (result) store.print() Solution cost is 1235 DFS1: DFS([bread = 0, milk = 0, cheese = 5, potato = 19, fish = 5, yogurt = 0], null, null, org.jacop.search.IndomainMin@6 63411de) true *** Store bread = 0 milk = 0 cheese = 5 potato = 19 fish = 5 yogurt = 0 Cost = 1235 Protein = 997 Fat = 819 Carbohydrates = 4314 Calories = 30230 Diet problem (JaCoP)

Slide 85

Slide 85 text

// imports ... def store = new Store() def unbounded = 100000 // scale quantities by 10, coefficients by 10, products by 100 def bread = new IntVar(store, 'bread', 0, unbounded) def milk = new IntVar(store, 'milk', 0, 10) def cheese = new IntVar(store, 'cheese', 0, unbounded) def potato = new IntVar(store, 'potato', 0, unbounded) def fish = new IntVar(store, 'fish', 5, unbounded) def yogurt = new IntVar(store, 'yogurt', 0, unbounded) IntVar[] all = [bread, milk, cheese, potato, fish, yogurt] int[] price = [20, 35, 80, 15, 110, 10] def cost = new IntVar(store, 'Cost', 0, unbounded) def protein = new IntVar(store, 'Protein', 0, 1000) store.impose(new Knapsack([40, 80, 70, 13, 80, 92] as int[], price, all, cost, protein)) def fat = new IntVar(store, 'Fat', 800, unbounded) store.impose(new Knapsack([10, 50, 90, 1, 70, 10] as int[], price, all, cost, fat)) def carbs = new IntVar(store, 'Carbohydrates', 1000, unbounded) // carbs variable can lead to -ve profit which violates knapsack //store.impose(new Knapsack([150, 117, 4, 226, 0, 170] as int[], price, all, cost, carbs)) store.impose(new LinearInt(all, [150, 117, 4, 226, 0, 170] as int[], '=', carbs)) def calories = new IntVar(store, 'Calories', 30000, unbounded) store.impose(new Knapsack([900, 1200, 1060, 970, 1300, 1800] as int[], price, all, cost, calories)) def search = new DepthFirstSearch() def select = new SimpleSelect(all, null, new IndomainMin()) def result = search.labeling(store, select, cost) println result if (result) store.print() Solution cost is 1235 DFS1: DFS([bread = 0, milk = 0, cheese = 5, potato = 19, fish = 5, yogurt = 0], null, null, org.jacop.search.IndomainMin@5 a9d6f02) true *** Store bread = 0 milk = 0 cheese = 5 potato = 19 fish = 5 yogurt = 0 Cost = 1235 Protein = 997 _1 = 0 Fat = 819 _2 = 0 Carbohydrates = 4314 Calories = 30230 _3 = 0 Diet problem (JaCoP with Knapsack)

Slide 86

Slide 86 text

@PlanningSolution @ToString(includeNames = true, includePackage = false) class DietSolution { @PlanningEntityCollectionProperty List foods @ValueRangeProvider(id = "amount") CountableValueRange getAmount() { ValueRangeFactory.createIntValueRange(0, 200) } @PlanningScore HardSoftScore score String pretty() { def sb = new StringBuilder() foods.each { sb << "$it.name: ${it.amount / 100}\n" } for (name in ['fat', 'carbs', 'protein', 'calories', 'cost']) { sb << "Total $name: ${foods.sum{ f -> f."$name" * f.amount / 100 }}\n" } sb << "Score: $score" sb } } Diet problem (OptaPlanner) def unsolved = new DietSolution(foods: [ new Food(name: 'Bread', cost: 2.0, protein: 4.0, fat: 1.0, carbs: 15.0, calories: 90), new Food(name: 'Milk', cost: 3.5, protein: 8.0, fat: 5.0, carbs: 11.7, calories: 120), new Food(name: 'Potato', cost: 1.5, protein: 1.3, fat: 0.1, carbs: 22.6, calories: 97), new Food(name: 'Cheese', cost: 8.0, protein: 7.0, fat: 9.0, carbs: 0.4, calories: 106), new Food(name: 'Fish', cost: 11.0, protein: 8.0, fat: 7.0, carbs: 0.0, calories: 130), new Food(name: 'Yogurt', cost: 1.0, protein: 9.2, fat: 1.0, carbs: 17.0, calories: 180) ]) def factory = SolverFactory.createFromXmlResource("dietSolverConfig.xml") def solver = factory.buildSolver() def solved = solver.solve(unsolved) println solved.pretty()

Slide 87

Slide 87 text

@PlanningEntity @ToString(includePackage = false) class Food { String name @PlanningVariable(valueRangeProviderRefs = "amount") Integer amount // times 100 double cost, protein, fat, carbs, calories } Diet problem (OptaPlanner) class DietConstraintProvider implements ConstraintProvider { @Override Constraint[] defineConstraints(ConstraintFactory factory) { new Constraint[]{ maxField(factory, 'protein', 10), minField(factory, 'fat', 8), minField(factory, 'carbs', 10), minField(factory, 'calories', 300), minFood(factory, 'Fish', 50), maxFood(factory, 'Milk', 100), minCost(factory), } } …

Slide 88

Slide 88 text

private Constraint minField(ConstraintFactory factory, String fieldName, double minAmount) { factory.from(Food).filter(f -> f.amount > 0) .groupBy(sum((ToIntFunction)(f -> (f."$fieldName" * f.amount).toInteger()))) .filter(fs -> fs < minAmount * 100) .penalize("Min $fieldName", ONE_HARD) } private Constraint maxField(ConstraintFactory factory, String fieldName, double maxAmount) { factory.from(Food).filter(f -> f.amount > 0) .groupBy(sum((ToIntFunction)(f -> (f."$fieldName" * f.amount).toInteger()))) .filter(fs -> fs > maxAmount * 100) .penalize("Max $fieldName", ONE_HARD) } private Constraint minFood(ConstraintFactory factory, String foodName, double minAmount) { factory.from(Food) .filter(f -> f.name == foodName && f.amount < minAmount) .penalize("Min $foodName", ONE_HARD) } private Constraint maxFood(ConstraintFactory factory, String foodName, double maxAmount) { factory.from(Food) .filter(f -> f.name == foodName && f.amount > maxAmount) .penalize("Max $foodName", ONE_HARD) } … Diet problem (OptaPlanner)

Slide 89

Slide 89 text

… private static ToIntFunction totalCost = f -> (f.cost * f.amount).toInteger() private static Constraint minCost(ConstraintFactory factory) { factory.from(Food).filter(f -> f.amount > 0).groupBy(sum(totalCost)) .penalize('Min cost', ONE_SOFT, fs -> fs >> 2) } } Diet problem (OptaPlanner) 00:28:33.285 [main] INFO o.o.core.impl.solver.DefaultSolver - Solving started: time spent (103), best score (-6init/-3hard/0soft), environment mode (REPRODUCIBLE), random (JDK with seed 0). 00:28:33.517 [main] INFO o.o.c.i.c.DefaultConstructionHeuristicPhase - Construction Heuristic phase (0) ended: time spent (335), best score (- 1hard/-522soft), score calculation speed (5221/sec), step total (6). 00:28:53.182 [main] INFO o.o.c.i.l.DefaultLocalSearchPhase - Local Search phase (1) ended: time spent (20000), best score (-1hard/-246soft), score calculation speed (87928/sec), step total (2445). 00:28:53.182 [main] INFO o.o.core.impl.solver.DefaultSolver - Solving ended: time spent (20000), best score (-1hard/-246soft), score calculation speed (86511/sec), phase total (2), environment mode (REPRODUCIBLE). Bread: 0.12 Milk: 0.1 Potato: 0.06 Cheese: 0.34 Fish: 0.5 Yogurt: 0.96 Total fat: 8.146 Total carbs: 20.782 Total protein: 16.57 Total calories: 302.46000000000004 Total cost: 9.86 Score: -1hard/-246soft

Slide 90

Slide 90 text

Constraint programming summary • Declarative approach • Involves variables, domains, and constraints • Solvers off sophisticated toolkit for expressing and efficiently solving a wide range of combinatorial & optimization problems • May involve backtracking, constraint propagation and heuristics • Multiple libraries available for most languages

Slide 91

Slide 91 text

References • https://en.wikipedia.org/wiki/Constraint_programming • https://www.cs.upc.edu/~erodri/webpage/cps//theory/cp/intro/slides.pdf • https://www.cse.unsw.edu.au/~tw/brwhkr08.pdf • https://www.analyticsvidhya.com/blog/2017/02/lintroductory-guide-on-linear-programming- explained-in-simple-english/ • Computational Mixed-Integer Programming https://www.birs.ca/cmo-workshops/2018/18w5208/files/GleixnerAmbros.pdf • Integer Linear Programming: Graphical Introduction https://www.youtube.com/watch?v=RhHhy-8sz-4 • Integer Linear Programming: Binary Constraints https://www.youtube.com/watch?v=-3my1TkyFiM https://www.youtube.com/watch?v=B3biWsBLeCw https://www.youtube.com/watch?v=MO8uQnIch6I • Linear Programming: Alternate solutions, Infeasibility, Unboundedness, & Redundancy https://www.youtube.com/watch?v=eMA0LWsRQQ