PyCSP3

# Incremental Solving

Interestingly, one can really pilot the solving process by iteratively adding and/or removing constraints (and also adding/changing the objective), handling a form of incremental solving. To add constraints, we already know that it suffices to call satisfy(). To remove constraints, it suffices to call the function:

unpost()

When this function is called, the last posting operation is discarded: it corresponds to all constraints that were posted by the last call to satisfy(). It is also possible to give the index of the posting operation, and even a second parameter indicating the index of constraint(s) inside the specified posting operation.

Below, we illustrate incremental solving by showing how to enumerate solutions by means of solution-blocking constraints, how to simulate an optimization procedure and how to compute diversified solutions. But first, we need first to import the library PyCSP$^3$:

from pycsp3 import *


## Enumerating Solutions with Solution-Blocking Constraints

For a given CSP $P$, a solution-blocking constraint of $P$ is a constraint that forbids a solution of $P$ (i.e., forbids a complete instantiation of the variables of $P$ corresponding to a solution). An original (but not necessarily efficient) way of enumerating the solutions of $P$ with a solver $S$ (that can, for example, only output a single solution) is to find solutions in sequence with $S$ while posting a new solution-blocking constraint every time a solution is found.

Let us consider the following toy model:

x = VarArray(size=4, dom=range(7))

satisfy(
AllDifferent(x),
Increasing(x),
Sum(x) == 10
);


On the one hand, we can directly compute the number of solutions (and display them) as follows:

if solve(sols=ALL) is SAT:
print(n_solutions())

4


On the other hand, we can enumerate the solutions of this model by successively posting solution-blocking constraints:

cnt = 0
while solve() is SAT:
cnt += 1
print(f"Solution {cnt}: {values(x)}")
satisfy(x != values(x))

Solution 1: [0, 1, 3, 6]
Solution 2: [0, 1, 4, 5]
Solution 3: [0, 2, 3, 5]
Solution 4: [1, 2, 3, 4]


We can display the internal representation of the posted constraints; this way, although a little bit technical, we can see that the presence of the four solution-blocking constraints. Note that ‘le’ and ‘eq’ stands for ‘less than or equal to’ and ‘equal to’.

print(posted())

allDifferent(list:x[])
ordered(list:x[], operator:le)
sum(list:x[], condition:(eq,10))
extension(list:x[], conflicts:(0,1,3,6))
extension(list:x[], conflicts:(0,1,4,5))
extension(list:x[], conflicts:(0,2,3,5))
extension(list:x[], conflicts:(1,2,3,4))


## Simulating an Optimization Procedure

For a given CSP $P$, an independent integer cost function $f$ to be minimized, defined from (the Cartesian product of the domains of) a subset $X$ of variables of $P$ to $\mathbb{Z}$, and a solution $sol$ of $P$ whose cost computed by $f$ is $B$ , a bound-improving constraint of $P$ wrt $f$ and $sol$ is a constraint that forbids all solutions of $P$ with a bound greater than or equal to $B$: it can be written $f(X) < B$. An original (but not necessarily efficient) way of finding an optimal solution of $P$ wrt $f$ with a CSP solver $S$ is to find solutions in sequence with $S$ while posting a new bound-improving constraint every time a solution is found.

We give an illustration with the Prime-looking problem. But, we need first to clear everythinh has beend declared and posted.

clear()


The prime-looking problem is from Martin Gardner: a number is said to be prime-looking if it is composite but not divisible by 2, 3 or 5. We know that the three smallest prime-looking numbers are 49, 77 and 91. Can you find the prime-looking numbers less than 250?

The model is rather simple:

# the number we look for
x = Var(range(250))

# a first divider
d1 = Var(range(2, 250))

# a second divider
d2 = Var(range(2, 250))

satisfy(
x == d1 * d2,
x % 2 != 0,
x % 3 != 0,
x % 5 != 0,
d1 <= d2
);


Let us consider that we have a cost function for this problem; this is simply the variable $x$ (to be maximized). One way of ensuring that we get a better solution after finding a first solution is given by the following piece of code:

if solve() is SAT:
print("The prime-looking number is: ", value(x))
satisfy(x > x.value)
if solve() is SAT:
print("The prime-looking number is: ", value(x))

The prime-looking number is:  49
The prime-looking number is:  77


If we want to find an optimal solution, we can write instead:

while True:
if solve() is not SAT:
break
print("The prime-looking number is: ", value(x))
satisfy(x > x.value)

The prime-looking number is:  77
The prime-looking number is:  91
The prime-looking number is:  119
The prime-looking number is:  121
The prime-looking number is:  133
The prime-looking number is:  161
The prime-looking number is:  203
The prime-looking number is:  209
The prime-looking number is:  217
The prime-looking number is:  221
The prime-looking number is:  247


In some cases, one may be worried of posting many bound-improving constraints, knowing that only the last one is relevant (since it is stronger than the other ones). This can be observed that when displaying all posted constraints.

print(posted())

intension(function:eq(x,mul(d1,d2)))
intension(function:ne(mod(x,2),0))
intension(function:ne(mod(x,3),0))
intension(function:ne(mod(x,5),0))
intension(function:le(d1,d2))
intension(function:gt(x,49))
intension(function:gt(x,77))
intension(function:gt(x,91))
intension(function:gt(x,119))
intension(function:gt(x,121))
intension(function:gt(x,133))
intension(function:gt(x,161))
intension(function:gt(x,203))
intension(function:gt(x,209))
intension(function:gt(x,217))
intension(function:gt(x,221))
intension(function:gt(x,247))


To test another approach, let us first discard all constraints.

unpost(ALL)


We can check that there are no more constraints.

print(posted())

[]


We post again the basic constraints of the model.

satisfy(
x == d1 * d2,
x % 2 != 0,
x % 3 != 0,
x % 5 != 0,
d1 <= d2
);


And, to avoid keeping useless bound-improving constraints, we can store the object (constraint) that was posted previously so as to be able to delete it afterwards. This gives:

objective = None
while True:
if solve() is not SAT:
break
print("The prime-looking number is: ", value(x))
if objective is not None:
objective.delete()
objective = satisfy(x > x.value)

The prime-looking number is:  49
The prime-looking number is:  77
The prime-looking number is:  91
The prime-looking number is:  119
The prime-looking number is:  121
The prime-looking number is:  133
The prime-looking number is:  161
The prime-looking number is:  203
The prime-looking number is:  209
The prime-looking number is:  217
The prime-looking number is:  221
The prime-looking number is:  247


We can check that only one bound-improving constraint is present.

print(posted())

intension(function:eq(x,mul(d1,d2)))
intension(function:ne(mod(x,2),0))
intension(function:ne(mod(x,3),0))
intension(function:ne(mod(x,5),0))
intension(function:le(d1,d2))
intension(function:gt(x,247))


unpost()


We control that it is no more present.

print(posted())

intension(function:eq(x,mul(d1,d2)))
intension(function:ne(mod(x,2),0))
intension(function:ne(mod(x,3),0))
intension(function:ne(mod(x,5),0))
intension(function:le(d1,d2))


As an alternative, it is possible to call the function unpost() that discards the constraint(s) posted at the last call to satisfy().

objective = False
while True:
if solve() is not SAT:
break
print("The prime-looking number is: ", value(x))
if objective:
unpost()
else:
objective = True
satisfy(x > x.value)

The prime-looking number is:  49
The prime-looking number is:  77
The prime-looking number is:  91
The prime-looking number is:  119
The prime-looking number is:  121
The prime-looking number is:  133
The prime-looking number is:  161
The prime-looking number is:  203
The prime-looking number is:  209
The prime-looking number is:  217
The prime-looking number is:  221
The prime-looking number is:  247


## Computing Diversified Solutions

Instead of enumerating solutions in the order “fixed” by the solver, one may want to diversify computed solutions by exploiting some distances. In other words, we may be interested in diverse solutions.

Below, we give a few illustrations. But, we need first to clear everythinh has beend declared and posted.

clear()


Let us consider the following toy model:

n = 8

x = VarArray(size=n, dom=range(5))

satisfy(
Maximum(x) == 4
);


If we want to enumerate 5 solutions while maximizing the Hamming distance between found solutions, we can write this piece of code:

solutions = []
while len(solutions) < 5 and solve() in {SAT, OPTIMUM}:
print("Solution: ", values(x))
solutions.append(values(x))
maximize(
Sum(x[i] != solution[i] for i in range(n) for solution in solutions)
)

Solution:  [0, 0, 0, 0, 0, 0, 0, 4]
Solution:  [1, 1, 1, 1, 1, 1, 4, 0]
Solution:  [2, 2, 2, 2, 2, 4, 1, 1]
Solution:  [3, 3, 3, 3, 4, 2, 2, 2]
Solution:  [4, 4, 4, 4, 3, 3, 3, 3]


Note that the problem instance is initially a CSP, and then becomes a COP because an objective is posted after the first turn of the loop (note also that any new objective overwrites the previous one, if any is present). This is the reason why we check if the solving status is either SAT or OPTIMUM.

clear()


As a second illustration, let us consider the following model:

n = 8

x = VarArray(size=n, dom=range(7))

satisfy(
Sum(x) == 22
);


If we want to enumerate 5 solutions while maximizing the Euclidean distance between found solutions, we can write:

solutions = []
while len(solutions) < 5 and solve() in {SAT,OPTIMUM}:
print("Solution: ", values(x))
solutions.append(values(x))
maximize(
Sum(abs(x[i] - solution[i]) for i in range(n) for solution in solutions)
)

Solution:  [0, 0, 0, 0, 4, 6, 6, 6]
Solution:  [4, 6, 6, 6, 0, 0, 0, 0]
Solution:  [6, 4, 0, 0, 6, 0, 0, 6]
Solution:  [0, 0, 4, 6, 0, 6, 6, 0]
Solution:  [0, 6, 6, 0, 6, 4, 0, 0]