jeudi 23 décembre 2010

CP vs MIP

I recently joined n-Side (an optimization company) and yesterday I introduced CP to my colleagues having rather a MIP background. So my yesterday "ennemy" was MIP (OK I agree that this is a bit provocative and the comparisons below are slightly insencere ;-) )







To the question: Why is CP not more famous, I presented this video





and explained that even if every one agrees that MAC is clearly superior over PC, MAC has only 5% market share. CP is the MAC of OR, less used ... but better ;-)

ps: let me now if you find/ use other funy comparisons when you introduce CP.

disclaimer: these comparisons only reflects my strictly personal point of view.

mardi 23 novembre 2010

Go Spurs Go!

I was in Texas in November to attend INFORMS. Precisely, I was in Austin. There were a lot of people: more than 4,000! The hotel was crowded...
Unfortunately, there were only few people attending the 2 sessions about Constraint Programming. I think about 15 people were in the room for each session. However, I found the talks quite interesting and this opinion was share by some other persons. I particularly appreciated the talk given by John Hooker about MDD based solvers. Nevertheless, I think that John should really try to combine some constraints like he started to do with the among constraints instead of claiming that a whole solver can be based on MDD. The use of MDD for combining some constraints sharing only few variables seems to be a very promising idea.
Personnaly, I gave a talk about nurse scheduling, I changed a little bit my CPAIOR talk and it seems that people appreciated it. I also met Ashish Sabharwal who just moved from Cornell to IBM Watson Research Center (near New York city). It seems that some other people are moving there, which has been confirmed by several other persons at IBM. I will investigate more if I have time.

If you want to find more information you can look at Michael Trick's blog where you will find a picture about the Twitters group at Informs (http://mat.tepper.cmu.edu/blog/?p=1258).
On my side, I have also a picture to show. Before, going to Austin I went to San Antonio. I have been really well welcome as you can see on this picture


Everybody can recognize the Spurs because it is written on the Pom-Pom girl dresses...
It was my first NBA game and I really enjoyed it! Tony Parker is a really good player!




jeudi 18 novembre 2010

Vehicle Routing Library

Vincent Furnon, who moved from ILOG to Google, has just released a vehicle routing library.
This means that you can really save money by using some products instead some others. In addition, some nice and clever people will be very happy to help you to develop your application whereas some others will see more and more the life in dark blue.


Here is his presentation about his new development. I really encourage you to test it.


"We just pushed the vehicle routing extension of the CP Solver. It's a 
layer above the CP Solver which lets you model a wide range of vehicle 
routing problems from the Traveling Salesman Problem (and its 
variants, ATSP, TSPTW, ...) to multi-vehicle problems with dimension 
constraints (capacities, time windows) and various "routing" 
constraints (optional nodes, alternate nodes,...). Given it exposes 
the CP variables one can extend the model using the constraints 
available in the CP Solver. 

Concerning the resolution of the problem, most of the "cabling" is
hidden so you just need to call the Solve() method on the RoutingModel
class. However the search can be parametrized using command-line
flags. We are mainly using CP-based local search and large
neighborhood search using routing-specific neighborhoods.
Implementations of Tabu Search and Guided Local Search are available
too and have proven to give good results (especially GLS). 
The basic library is in C++ and a Python-wrapped version is available
too. 
Documentation will be added in the future but you can have a look at
the TSP example (both in C++ and Python) as a starting point. We'll be
adding a CVRPTW example very shortly. "

mercredi 13 octobre 2010

The importance of redundant constraints: Modeling a Bin-Packing

What makes CP really fun is that modeling is an art and you have to understand (or at least have a good feeling) of how the propagation will behave when you write a model.
Detractors of CP will say it's a weakness, Pascal Van Hentenryck will say that "CP is a smart tool for smart people" ;-)
When I wrote 5 years ago my second CP model (first one was n-queens) it was the Balanced Academic Curriculum Problem (BACP).
I started to read the paper "Modelling a balanced academic curriculum problem" (CPAIOR2002) by Brahim Hnich , Zeynep Kiziltan and Toby Walsh. The experiments in this paper show that the problem is really difficult (they cannot prove optimality with an hybridization of CP and MIP).
In reality this problem is not difficult for CP but like many problems, it has bin-packing component and you have to model it correctly. I understood this when I read the very good paper of Paul Shaw: "A Constraint for Bin Packing" (CP2004).

A bin packing constraint implies some sized items that you must place into some capacitated bins. So you have one variable for each item that tells you in which bin you place the item and one variable for each bin telling you how loaded is the bin. A ground instance is bin-packing([1,2,1,3],[3,4,1,5],[4,4,5]) with 5 items placed in bins [1,2,1,3] with sizes [3,4,1,5] so that the load of bins 1,2,3 are [4,4,5].

A model for this constraint in google cp-solver is the following:


def BinPacking(cp, binvars, weights, loadvars):
'''post the load constraint on bins.

constraints forall j: loadvars[j] == sum_i (binvars[i] == j) * weights[i])
'''
for j in range(len(loadvars)):
b = [cp.IsEqualCstVar(binvars[i], j) for i in range(len(binvars))]
cp.Add(cp.ScalProd(b, weights) == loadvars[j])
cp.Add(cp.SumEquality(loadvars, sum(weights)))


For all bins j, you say that the load is equal to the sum of the weights of items placed into that bin j. This is enough to model the semantic of the constraint but not enough for the pruning. It is very important to add the redundant constraints that the sum of the load variables is equal to the sum of the total size of the items. If you don't do that, there is not communication between the bins in the cp store. Consider this example to realize this:

You have 5 items with size (6,5,5,4,4,4,2) to place into three bins with capacity 10. There is clearly a solution: (5,5) in bin one, (4,4,2) in bin two and (6,4) in the third one. Assume now that you have a partial solution where one item of 5 and one item of 4 is placed into bin one. Clearly this partial solution cannot lead to a final solution since you know that the bins must be completely packed at the end (30 units to place and three bins with capa 10). Unfortunately you can see that the first bin cannot be packed completely with this partial solution. Interestingly this inconsistency can only be detected with the redundant constraint "the sum of the load variables = sum of the total size". Indeed, the load of bin 1 is at most 9, and the load of bins 2 and 3 is at most 10 so at most the total load will sum to 29, which is smaller than 30 => fail!

Now if you want to experiment (like the authors of the first BACP model) how difficult BACP is without this redundant constraint, you can comment the redundant constraint in the google-solver cp model available here

If you don't forget the redundant constraint, you'll solve the problem in 68ms and 1000 failures with a default first fail heuristic. In conclusion this problem is trivial for CP and I hope if you read this that you will never forget to post this redundant constraint to model a bin-packing like me 5 years ago ;-)

Note that this problem has been generalized to be more realistic here. I think this new problem is now challenging for CP.

mardi 12 octobre 2010

Improvement or Innovation ?

There is currently a discussion on the blog of Michael Trick about the LDS (Limited Discrepency Search) which has neither any relation with LSD nor with Lucy, Diamond or Sky, but which comes from a request of Matt Ginsberg asking Christian Schulte for changing the licence of Gecode because the use of Limited Discrepency Search for solving some job scheduling problems is patented.

The same story happened when I was at ILOG. So, M. Ginsberg seems to like repeated business :-) I don't know exactly what was the agreement between ILOG and M. Ginsberg. Everybody knows that  ILOG removed the
LDS acronym and name from all the documentation and from the product, but there are some unknown funny stories (because cp is fun) behind that old problem.

- we discovered that we didn't correctly understand the LDS and what we implemented was  not the LDS. Thus, we name it Slice Based Search and it was published under the name Discrepancy-Bounded Depth First Search (J. Christopher Beck and Laurent Perron)
- all that stuff does not really work well. I think that there is no more any CP model of benchmarks using it.
- LDS is clearly outperformed by restarts with randomization.

However, there are some other lessons that are interesting:

First, it is not because it is open source that you have the right to copy everything and put it into your product.

Second, M. Ginsberg is right when he says :"But if a commercial company wants to sell an implementation of LDS for use in job scheduling, it is — in my opinion — not unreasonable for the people who invented LDS to be compensated.". With such an idea I could become rich, especially if I create a company with N. Beldiceanu :-) But how can we implement this idea ? By patents ? I am not sure it is the good solution. We could accept this only when the idea is really an invention, that's mean it is strongly original. In other words, it is an inovation and not only an improvement. Because if it is an improvement we should also compensate the people who introduce the idea before, and so on... But how do we distinguish these two notions? In the comments of Mike's blog we can see that this is the point! T. Walsh proposed the Left First Search : “Thus, in searching the tree we want to minimize right arcs … We call a search algorithm which does this left-first search. At the n+1-th ply of the search we explore all those nodes whose path back to the root includes n right arcs”. Then Toby proposed an exponential algorithms. If you recognize LDS in Toby's words then LDS is an improvement otherwise it could be seen as an inovation...

At the end, I hope that M. Ginsberg made more money by the breakthrough brought by LDS over other search methods than by the patent and the threats associated with it...

lundi 11 octobre 2010

May the force be with you!

I found another vendor of the future headset interface!
It can help you to control the force:
http://company.neurosky.com/products/force-trainer/

Even if they make toys, the company seems to be serious
http://www.neurosky.com/

I should try one of these headsets and keep you inform about my tests :-)

[Edit] If you have an iphone, you can immediatly order one : http://www.plxwave.com/

lundi 4 octobre 2010

CP is fun

CP is fun is the title of this blog.
I needed to find quickly a title when I created the blog and I thought to this sentence. Ten years ago, when I gave talk at Informs, I always finished my talks by this sentence in order to convince people to be more interested in CP. Ten years later I am not sure that I was successfull, except if I consider that Irv Lustig told me once that he really enjoyes when he defines CP models.

When I read Pierre's posts, I really think that CP is fun because we can do a lot of things.
However, one sentence of Pierre drew my attention to optimization languages :
"When you take a close look at the line:
solver.Add(sum([10**(n-i-1)*x[i] for i in range(n)]) == nb)
It is difficult to argue that it is very far from dedicated optimization languages!"

 
I think the association of the concept of objects and libraries kill a lot of languages. Consider C++, do you want to have a kind of Lisp ? Use a library which manipulates lists like lisp! Do you want to use iterations with mapping functions ? Use  the stl or meta-programmation! In addition, you will benefit from all the advantages of the recent compilers and IDE. I think it is really difficult to obtain such a result with a dedicated language.
Thus, I strongly support Pierre.
 
However, this is not the best interface we can imagine, because I think that I found it.
Look at this website http://www.emotiv.com/
Being able to model and to solve some problems with CP while using their products will definitively be a very good way to prove that CP is really fun!

samedi 2 octobre 2010

Open positions

You would like to work on the resolution of problems like green-it, cloud computing or program verification with constraint programming?
You would like to work on industrial problems like some of Google?

We have some open positions at the University of Nice-Sophia Antipolis at the I3S department. We have some grants for master students, phD or postdoctoral stay.
The dates are flexible.

Interested ? Feel free to send me a message : jcregin@gmail.com

vendredi 1 octobre 2010

Multiplication vs Exponent Constraints

In the Dudeney problem we have the constraint:
solver.Add(nb = s*s*s)
and the model that I have shown previously finishes with 92 fails with a labeling on x.
But behind the scene, google solver is smart because it doesn’t post multiplication constraints only. The first multiplication introduces a variable (call it t) with the constraint (t==s^2) (an exponent two constraint). Then the multiplication constraint (nb==t*s) is posted.

Is it important? Yes for the pruning of the domains. And CP is all about pruning!

Let’s consider some domains on our variables to realize this. Assume that s has a domain [1..10] and t has a domain of [25..100]. Then you can prune the interval of s to [5..10] with the power two constraint. With a simple multiplication you can only prune it to [3..10].

If you want to see the effect on the Dudeney problem with only multiplication we can post this such that the solver has no way to know that a power two constraint can be posted:


s1 = solver.IntVar(range(1,9*n+1),'s1')
solver.Add(s1 == s)
solver.Add(nb == s*s1*s)

And now you get 5821 fails!

Can we do better in this case?
Yes if you have a power three constraint ;-) But this is generally not implemented in solvers because it’s not so used.
You can also do better on the Dudeney problem with a labeling on another decision variable but this is another story...

mercredi 29 septembre 2010

Getting all the solutions with google solver

Laurent explained to me how to get all the solutions in a more elegant way than with collectors and assignment. Here is the trick for the Dudeney problem:


solver.NewSearch(db)
while solver.NextSolution():
print nb
solver.EndSearch()


What I still have to discover is:
- How to get the fix point after each post in the python model ?
- How to implement a search (at least variable/value heuristic) in python ?
- How to make LNS and defining a fragment from python ?

mardi 28 septembre 2010

A Martin Gardner problem

Here is a Gardner problem that was posted today on this site
The solution I found manually is really ugly. I would be happy if someone could tell me an elegant way to solve it. In case you want to check your solution here is a google cp-solver model:

  
from constraint_solver import pywrapcp

solver = pywrapcp.Solver('Gardner')
dx = [solver.IntVar(range(10),'x'+str(i)) for i in range(2)]
dy = [dx[1],dx[0]]
x = solver.IntVar(range(1,100),'x'+str(i))
y = solver.IntVar(range(1,100),'x'+str(i))

m = solver.IntVar(range(1,100000),'m')

solver.Add(x == 10*dx[0]+dx[1])
solver.Add(y == 10*dy[0]+dy[1])

solver.Add(x*x - y*y == m*m)

solution = solver.Assignment()
solution.Add(x)
solution.Add(y)
solution.Add(m)
collector = solver.AllSolutionCollector(solution)

solver.Solve(solver.Phase(dx, solver.INT_VAR_DEFAULT,
solver.INT_VALUE_DEFAULT),[collector])

for i in range(collector.solution_count()):
current = collector.solution(i)
x = current.Value(x)
y = current.Value(y)
m = current.Value(m)
print "sol:",x+y+m,"x:",x,"y:",y,"m:",m

print "#fails:",solver.failures()
print "time:",solver.wall_time()

dimanche 26 septembre 2010

Where are CP people located ?

It's really funny that after a few posts on this blog, we can already see a trend.
For those familiar with the CP/CPAIOR conferences, this also confirms where the authors of the publications are mainly coming from. I guess our goal is to contaminate all the map ;-)


US
334
France
73
Canada
44
UK
41
Germany
31
Australia
29
Belgium
25
Sueden
17
Netherland
13
India
12

vendredi 24 septembre 2010

Dudeney number

I discovered yesterday Dudeney Numbers
A Dudeney Numbers is a positive integer that is a perfect cube such that the sum of its decimal digits is equal to the cube root of the number. There are only six Dudeney Numbers and those are very easy to find with CP.
I made my first experience with google cp solver so find these numbers (model below) and must say that I found it very convenient to build CP models in python!
When you take a close look at the line: solver.Add(sum([10**(n-i-1)*x[i] for i in range(n)]) == nb)
It is difficult to argue that it is very far from dedicated optimization languages!

from constraint_solver import pywrapcp



def dudeney(n):

solver = pywrapcp.Solver('Dudeney')

x = [solver.IntVar(range(10),'x'+str(i)) for i in range(n)]

nb = solver.IntVar(range(1,10**n),'nb')

s = solver.IntVar(range(1,9*n+1),'s')



solver.Add(nb == s*s*s)

solver.Add(sum([10**(n-i-1)*x[i] for i in range(n)]) == nb)

solver.Add(sum([x[i] for i in range(n)]) == s)



solution = solver.Assignment()

solution.Add(nb)

collector = solver.AllSolutionCollector(solution)



solver.Solve(solver.Phase(x, solver.INT_VAR_DEFAULT,

solver.INT_VALUE_DEFAULT),[collector])



for i in range(collector.solution_count()):

current = collector.solution(i)

nbsol = current.Value(nb)

print nbsol



print "#fails:",solver.failures()

print "time:",solver.wall_time()



if __name__ == '__main__':

dudeney(6)

lundi 20 septembre 2010

Decomposition of sortedness global constraint

A recent trend of CP that I noticed is the decomposition of global constraints.

See for instance the papers:

· Decomposition of the NValue Constraint. C. Bessiere, G. Katsirelos, N. Narodytska, C.G. Quimper, T. Walsh. Proceedings CP'10, St Andrews, Scotland, pages 114-128.

· Decompositions of All-Different, Global Cardinality and Related Constraints. Christian Bessière, George Katsirelos, Nina Narodytska, Claude-Guy Quimper, and Toby Walsh. In Proceedings of the 21st International Joint Conference on Artificial Intelligence (IJCAI 09), pages 419-424, 2009.

· Decomposing Global Grammar Constraints. Claude-Guy Quimper and Toby Walsh. In Proceedings of the 13th International Conference on Principles and Practice of Constraint Programming (CP 07), pages 590-604, 2007.

You can find these papers on the personal page of Claude Guy

I find these decompositions very useful since it can avoid the error prone task of implementing difficult algorithms for global constraints.

We will make an attempt of decomposing another global constraint: the sortedness constraint (I’m not aware of a decomposition of this one).

The paper describing bound-consistent filtering for this constraint is:

Faster Algorithms for Bound-Consistency of the Sortedness and the Alldifferent Constraint, Kurt Melhorn and Sven Thiel, CP2000.

If you don’t know what the sortedness is about, here is a quick reminder:

sortedness(x[], s[],p[]) takes three f.d . var array as arguments and it enforces that s is a sorted version of x and p represents the permutation of the sort.

For instance sortedness(x[5,2,4,1,3],s[1,2,3,4,5],p[4,2,5,3,1]) is a valid instantiation.

If you don’t have the implementation of sortedness from the paper above at hand, here is what you usually do as a straightforward model:


post(s[i] == x[p[i]])
forall(i in 1..n-1) {
post(x[p[i]] <= x[p[i+1]])
post(s[i] <= s[i+1])
}
post(alldifferent(p))
Unfortunately the above model does not enforce bound consistency.

We will add redundant constraints to strengthen the decomposition and hopefully enforce bound consistency with the decomposition. The redundant constraints are based on the counting sort algorithm. Looking at our example sortedness(x[5,2,4,1,3],s[1,2,3,4,5],p[4,2,5,3,1]) the idea is the following:

How many values in x can be strictly smaller than s[i] (assuming indices starts at 1)? The answer is at most i-1 otherwise value s[i] would be at a greater index. This idea is the basis of the linear time counting sort algorithm. Here is how this idea is translated in CP terms:



minval = min(i in 1..n) x[i].min()
maxval = max(i in 1..n) x[i].max()
#define an array of variable occ with domains {0,...,n} that will represent the number of occurrences of each value
var occ [minval..maxval](0..n)
#gcc to count number of occurrences of each values in the array x
post(gcc(x,occ))
#same for the sorted ones, not sure it is useful
post(gcc(s,occ))
var nbBefore [minval..maxval](0..n) #nbBefore[i] = #{i | x[i]#obviously no values are smaller than minval
post(nbBefore[minval] == 0)
forall(i in minval+1..maxval) {
cp.post(nbBefore[i] == sum(j in minval..i-1) occ[j])
}
forall(i in 1..n) {
#there are strictly less than i values smaller than s[i]
cp.post(nbBefore[s[i]] < i)
}

I still need to prove formally but I think this decomposition achieves bound-consistency for sortedness. Any comments are welcome (if you find it useful, if you find a counter example or a better decomposition or if something is not clear…).

mercredi 15 septembre 2010

Google CP Solver is available!

This is my first scoop!
You can find the source code of the Google CP solver at this address
https://code.google.com/p/or-tools/
you can checkout a read only copy anonymously over HTTP:

svn checkout http://or-tools.googlecode.com/svn/trunk/ or-tools-read-only
If you want to be informed about the future development dont' forget to become a member of the or-tools group:

Now, you can ask for a review :-) or better tell me what you think (not discuss about the fact that you were already aware of this project)
 

Welcome!

Some people asked me to create a blog about constraint programming.
I thought it could be a hard task, because I already need to maintain my personnal webpage, my webpage for my students and the webpages of the 2nd year students at the university. Fortunately, someone (I decided to keep the names secret in this blog) gave me the Google url for managing easily a blog (https://www.blogger.com/start?hl=fr). Thus, I created this blog in 1 minute (plus 1 min due to an internal error) which seems to be better than using wordpress which claims to require only 5 minutes.

So, welcome to this blog where you are going to learn some news about constraint programming and certainly also about some other stuff because I am sure that I will not resist to speak about different thinks.