Groovy Constraint Programming

4578d99b560b2a470e05288ef6766ac2?s=47 paulking
November 21, 2020

Groovy Constraint Programming

High-level languages can be characterised by the programming paradigm(s) they support. Broadly speaking, the three most common paradigms are the imperative (OO & procedural), functional and logic (or constraint) styles. If you are watching the top-ten list of trendy new languages, you probably aren't seeing any that claim support for logic (or constraint) programming as their primary paradigm - yet for certain kinds of problems, it is arguably the preferred style to use. This talk looks at this less-used style. What is it, when might you use it and how.

The examples use the Apache Groovy programming language and other (mostly JVM) languages which support this style of programming using libraries. Nearly all examples should be easily ported to your favourite programming language.

4578d99b560b2a470e05288ef6766ac2?s=128

paulking

November 21, 2020
Tweet

Transcript

  1. 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
  2. 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
  3. 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
  4. Friends of Apache Groovy Open Collective

  5. 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
  6. 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
  7. Constraint/Logic Programming Description • Declarative style of expressing problem •

    Model: variables + domains + constraints • Solver: does all the work
  8. 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
  9. Constraint Programming Picture: http://osullivan.ucc.ie/AAAI_OSullivan.pdf

  10. 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
  11. 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 …
  12. 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]
  13. 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)"
  14. 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
  15. 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
  16. 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)
  17. 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]" } }
  18. 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
  19. 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
  20. 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]
  21. 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
  22. 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
  23. 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.
  24. 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.
  25. 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.
  26. 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
  27. 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
  28. 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
  29. 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
  30. 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
  31. 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
  32. 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
  33. 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!
  34. 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)]
  35. 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]]
  36. 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)
  37. 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]
  38. 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/
  39. @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,
  40. Global Constraint Catalog https://sofdem.github.io/gccat/gccat/sec5.html

  41. 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
  42. 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
  43. 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)
  44. 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
  45. 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
  46. 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
  47. 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)
  48. 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
  49. 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
  50. 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
  51. 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
  52. 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
  53. 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
  54. 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
  55. 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
  56. 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
  57. 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
  58. 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
  59. 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
  60. 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
  61. 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
  62. 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
  63. 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
  64. 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 …
  65. SEND + MORE = MONEY Graphs generated using cpviz and

    choco-cpviz
  66. SEND + MORE = MONEY (Other JVM language examples) Prolog

  67. 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.
  68. 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
  69. 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)
  70. 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
  71. 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
  72. 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
  73. 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.
  74. 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
  75. Dietary restrictions

  76. Dietary restrictions (Excel)

  77. Dietary restrictions (Excel)

  78. Dietary restrictions (Google sheets)

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

  80. 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 ############################################
  81. 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.
  82. 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" }
  83. 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
  84. 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<IntVar>() def select = new SimpleSelect<IntVar>(all, null, new IndomainMin<IntVar>()) 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)
  85. // 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<IntVar>() def select = new SimpleSelect<IntVar>(all, null, new IndomainMin<IntVar>()) 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)
  86. @PlanningSolution @ToString(includeNames = true, includePackage = false) class DietSolution {

    @PlanningEntityCollectionProperty List<Food> foods @ValueRangeProvider(id = "amount") CountableValueRange<Integer> 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()
  87. @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), } } …
  88. 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)
  89. … private static ToIntFunction<Food> 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
  90. 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
  91. 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