Main

May 19, 2014

JaCoP 4.1 released (with support for float variables)

JaCoP version 4.1 has been released:
Dear users,

We have just released JaCoP 4.1. Feel free to clone JaCoP repository - https://github.com/radsz/jacop and contribute. In this release we introduce floating-point variables (FloatVar) and constraints on these variables. The methods used in floating point domain of JaCoP are based on floating point intervals and consistency methods build around them. We provide a set of basic arithmetic constraints as well as square root, absolute value, trigonometric constraints (sin, cos,tan, asin, acos, atan), exponential constraint and natural logarithm. Moreover element constraint and equation between integer and float variables makes it possible to build models that mix different JaCoP variables. We have also included experimental implementation of several features that is not consider as final and can change in future.

The floating point domain is fully integrated with JaCoP solver. The same search methods can be used as well as floating point variables can be used as cost functions for minimization. We provide also specialized optimization methods that are specially developed for floating point variables.

The released version provides also flatzinc interpreter for floating point variables. Minizinc models containing float variables can be compiled using standard mzn2fzn compiler and used by JaCoP solver.

This release fixes also few bugs.

Best regards, JaCoP Core Developer Team
The GitHub repo is: github.com/radsz/jacop.

I tested with some MiniZinc models with float variables, and it worked quite nice.

April 16, 2012

JaCoP version 3.2 released (bug fixes and Scala interface)

From JaCoP version 3.2 released :
Dear users, We have just released JaCoP 3.2. This is a release that fixes few bugs as well as provides an interface from Scala to JaCoP. Examples of using Scala are provided in ExamplesScala package. Additional discussion and examples are also available at

JaCoP/Scala and at Hakan Kjellerstand blog at A first look at scalaJaCoP (Scala wrapper for JaCoP)

best regards,

Core Developer Team
The latest version of JaCoP can be downloaded from Sourceforge jacop-solver. Also, see my JaCoP page.

August 09, 2011

A first look at scalaJaCoP (Scala wrapper for JaCoP)

In late May, the JaCoP team announced Scala and JaCoP, a Scala wrapper for JaCoP:
Scala programming language gets more and more acceptance in the community. It compiles to Java byte code and is executed using JVM. This makes it possible to use JaCoP directly in Scala by importing its packages. However, since Scala offers a nice way to define domain specific languages (DSL), we have developed an experimental implementation of the DSL for JaCoP solver.

Our implementation uses several Scala features to define a clear way of using JaCoP solver. We overloaded operators to specify easily arithmetical, logical and set constraints. We also use the package object to gather global constraints and search methods. This together with Scala features, such as implicit conversions, simplifies use of JaCoP and makes it easier to understand.
Some day after I started to test the first (beta) version and have had some interesting discussion with the developers, about bugs and also about syntax and functionality. Some weeks ago, they released a new version, which I personally think is better than the first version.

Programming in a high level programming language like Scala makes constraint programming much easier than using plain Java. With scalaJaCoP the JaCoP team has done a great job extending the world of constraint programming to other realms.

Below I explain some of the concepts I've met during my tests, as well as some gotchas. My public models has been collected in my JaCoP/Scala page.

Update 2011-08-11: Some of the gotchas mentioned below has been fixed.

Documentation

The current release of scalaJaCoP don't include any documentation, just the wrapper code and about 20 examples. Here is a little guide how to learn more about it:
  • JaCoP and the supported constraints are presented in JaCoP Library User’s Guide and JaCoP API
  • The Scala file jacop.Scala (in scalaJaCoP directory) is where the core JaCoP is wrapped.
  • The Scala file globals.scala (in scalaJaCoP directory) contains the global constraints supported by scalaJaCoP which is many (most) of JaCoP's global constraints. If you miss any JaCoP constraint, I think that the JaCoP team will be grateful for feedback.
Also, my JaCoP/Scala models may hopefully be of some use.

Running scalaJaCoP

Running scalaJaCoP is quite easy after installing JaCoP and Scala:
  • Download scalaJaCoP: Go to Sourceforge, select the latest JaCoP version and then download scalaJaCoP.zip. It contains two directories: scalaJaCoP (core files) and ExamplesScala (examples)
  • Compile the two core files: jacop.scala and global.scala:
    scalac scalaJaCoP/*.scala
  • Compile the examples:
    scalac ExamplesScala/*.scala
  • To compile an individual model:
    scalac ExamplesScala/Adder.scala
    scala ExamplesScala.Adder

    The program fsc (fast background compiler) is faster when compiling:
    fsc ExamplesScala/Adder.scala
    scala ExamplesScala.Adder


    Or as I tend to do:
    fsc SendMoreMoney.scala && scala SendMoreMoney
Note: If running into trouble after updating scalaJaCoP, it might be a good idea to remove the class files.

SendMoreMoney.scala

Here is a complete scalaJaCoP model, the standard SEND+MORE=MONEY problem: SendMoreMoney.scala. This examplifies much of what scalaJaCoP has to offer, in terms of high level constraints etc.
import scalaJaCoP._
object SendMoreMoney extends App with jacop {
  // define variables
  val s = new IntVar("s", 0, 9)
  val e = new IntVar("e", 0, 9)
  val n = new IntVar("n", 0, 9)
  val d = new IntVar("d", 0, 9)
  val m = new IntVar("m", 0, 9)
  val o = new IntVar("o", 0, 9)
  val r = new IntVar("r", 0, 9)
  val y = new IntVar("y", 0, 9)
  val all = Array(s,e,n,d,m,o,r,y)

  // constraints
  alldifferent(all)

  // note the # operator

            1000*s + 100*e + 10*n + d +
            1000*m + 100*o + 10*r + e #=
  10000*m + 1000*o + 100*n + 10*e + y

  s #> 0
  m #> 0

  // labeling
  val result =
     satisfyAll(search(all, first_fail, indomain_min), printIt) 

  statistics

  // output callback
  def printIt() {
    println("" +      s.value + e.value + n.value + d.value + " + " + 
                      m.value + o.value + r.value + e.value + " = " +
            m.value + o.value + n.value + e.value + y.value)

  }
}
Here are some things to note:

# for JaCoP's operators

The operators with "#" are operators that will be converted to JaCoP constraints. They are specially "tagged" in this way to distiguish them from Scala's operators. (Using operators without "#" are allowed in the current scalaJaCoP version but are deprecated.)

Output

Output can be done in two ways. Either as in this model by using a callback (which I happen to call printIt) that is called for all solutions. The function is placed last in satisfyAll. This is my preferred way.

The other way is to print according to JaCoP's default way of printing, i.e. just the last solution. Then satisfyAll should not have any callback function, of course. Many of the examples in the scalaJaCoP distribution use this principle.

Statistics

For getting statistics of the solution, just place statistics somewhere in the code. Then scalaJaCoP will print the (predefined) statistics information. Here is the statistics for the SendMoreMoney model:
Search statistics:
==================
Search nodes : 10
Search decisions : 5
Wrong search decisions : 5
Search backtracks : 5
Max search depth : 5
Number solutions : 1

Alternative version

Here is an alternative approach of the same problem (from SendMoreMoney2.scala) using an array for the variables and sum (as scalar product).
import scalaJaCoP._
object SendMoreMoney2 extends App with jacop {
  val all = Array.tabulate(8)(i => new IntVar("x["+i+"]", 0, 9))
  val Array(s,e,n,d,m,o,r,y) = all

  alldifferent(all)

  val values1 = Array(1000,100,10,1)
  val values2 = Array(10000,1000,100,10,1)

    ( sum(List(s,e,n,d), values1) 
     +  
      sum(List(m,o,r,e), values1) 
    )  #=  
  sum(List(m,o,n,e,y), values2)

  s #> 0
  m #> 0
  // ... labeling and output as above
}

Other things in scalaJaCoP

Below are some other things I've found during my tests that may be of some help when programming in scalaJaCoP.

Defining constraint functions

Defining constraints (decompositions) is quite straightforward. Here is an example from SequenceSum.scala to define sequence_sum, to ensure that the sum of each subsequence of length s sums to m:
def sequence_sum(y: Array[IntVar], m: IntVar, s: Int) = {
  val n = y.length
  for(i <- 0 until n - s + 1) {
    sum( Range(i,i+s).map(j => x(j) ).toList) #= m
  }
}

Summing an array

In Diet2 there is an example how to sum an array. Here we see both how to define and arrays, and also using sum constraint over a collection of calculated values.
// define variables
val x = Array.tabulate(n)(i => new IntVar("x["+i+"]", 0, 10))
val cost = sum(Array.tabulate(n)(i => x(i)*price(i)).toList)

// constraint
Array.tabulate(n)(
  nutr => sum(Array.tabulate(n)(i => x(i) * all(nutr)(i)))  #>= limits(nutr)
)

Matrices

In scalaJaCoP, matrices can be defined as List of Lists, e.g. the connection matrix in MapColoring.scala:
val connections =
   List(List(0, 0, 1, 1, 1, 1),
        List(0, 0, 0, 1, 0, 0),
        List(1, 0, 0, 1, 1, 0),
        List(1, 1, 1, 0, 1, 1),
        List(1, 0, 1, 1, 0, 0),
        List(1, 0, 0, 1, 0, 0))
The constraint to make the color different for the neighbours (the object of the model) is then defined as:
  for (c1 <- 0 until num_countries;
       c2 <- 0 until c1) {
    if (connections(c1)(c2)==1) {
      color(c1) #\= color(c2)
    }
  }
Note: the element access in the array (and matrix) are with (), i.e. not brackets ([]) as in Java.

Another note: there is no support for using an IntVar as an index in an IntVar array. For this the element constraint must be used. See below for more comments about element.

In WhoKilledAgatha.scala there is an example of an element constraint on a matrix.

yield

yield is not a JaCoP or Constraint Programming thing but a nice way in Scala to collect values given some condition or loop construct. The possibility to do these kind of constructs is why I like scalaJaCoP over the more complex way in Java.

Here's an example from Sudoku to get the boxes:
for(i <- 0 until reg; j <- 0 until reg) {
  alldistinct(  (for{ r <- i*reg until i*reg+reg;
                      c <- j*reg until j*reg+reg
                   } yield x(r)(c)).toArray
               )
Note that we have to convert the List made by the yield expression to an array (with toArray) since this is what alldictict expect.

Another example of yield is in Minesweeper.

IntVar intervals

This is a thing in JaCoP that has bitten me a few times before. JaCoP has defined the limits of IntVar to the be -10000000..10000000 and if a variable overflows it is considered to fail. This is, for example, what happens in DivisibleBy9Through1.scala for n larger than 8.

element

Element in JaCop is another thing that I have had some problems with before, namely that element is per default 1-base, and for use with 0-based one have to use an extra argument (offset) to the constraint to force it to be 0-based: e.g. element(p(i-1), x, p(i), -1).

Here is a 0-based definition of circuit (from CircuitTest0Based.scala), as well as a "extractor" of the path from the circuit (circuit_path), nice for representation of the solution.
def circuit_me(x: Array[IntVar]) = {

  val len = x.length
  val z = Array.tabulate(len)(i=>
           new IntVar("z("+i+")", 0, len-1))
  alldifferent(x)
  alldifferent(z)
  z(0) #= x(0)
  // then put the orbit of x(0) in z(1..n-1)
  for(i <- 1 until len) {
    // we must use element/4, with an offset
    element(z(i-1), x, z(i), -1)
  }
  z(len-1) #= 0
}

def circuit_path(x: Array[IntVar], p:Array[IntVar]) = {
  alldifferent(p)
  p(0) #= 0 // path always starts with 0
  for(i <- 1 until n) {
    element(p(i-1), x, p(i), -1)
  }
} 
(A 1-based version of the above is in CircuitTest.scala.)

Some other models that use element: Langford.scala (Langford's number problem), and WhoKilledAgatha.scala.

Reification

Unfortunately scalaJaCoP don't have a high level support (i.e. syntactic sugar) for reifications. In AllDifferentExcept0.scala there is an example of how to use reification for decomposing the global constraint alldifferent_except_0. The high level version could maybe be something like:
 // This is NOT valid scalaJaCoP code:
for(i <- 0 until y.length; j <- 0 until i) {
  (y(i) #\= 0 AND y(j) #\= 0) -> (y(i) #\= y(j))
}
As of now, one have to write it something like this, i.e. define a BoolVar with the constraint, and then make the reification.
// alldifferent except 0
def alldifferent_except_0(y: Array[IntVar]) = {
  for(i <- 0 until y.length; j <- 0 until i) {
    val b = new BoolVar("b")
    b <=> AND((y(i) #\= 0), (y(j) #\= 0)); 
    b -> (y(i) #\= y(j))
  }
}
A model that use many refications is WhoKilledAgatha.scala. This was one of the models I struggled with quite a long time before I realized about the special treatment for 0-based element.

More about the #operator

Update 2011-08-11: Please note that these problems with #-operator has been fixed in the latest version of scalaJaCoP.

Another thing to be aware of in scalaJaCoP is when using both a (Scala) integer variable and a IntVar: it may matter if the (Scala) integer variable is to the left or right of the #operator. For example, in Xkcd.scala (Xkcd knapsack problem) there are some examples of this, shown here.

Here x is an IntVar array, defined as:

val x = Array.tabulate(num_prices)(i => new IntVar("x["+i+"]", 0, 10))

The object of the model is to ensure that the number of items in x sums to a certain total, namely 1505, which is defined as a Scala val variable.

val total = 1505

Now, if total is on the left side of the #= operator, Scala cannot figure out how to convert to an IntVar:

// don't work
total #= sum(Array.tabulate(num_prices)(i=> x(i)*price(i)))

This results in errors like error: reassignment to val.

However, there is a couple ways that do work, for example moving the integer variable to the right side, which I tend to do:

sum(Array.tabulate(num_prices)(i=> x(i)*price(i))) #= total

Or explicitly convert total to an IntVar with intToIntVar or (:IntVar)

intToIntVar(total) #= sum(Array.tabulate(num_prices)(i=> x(i)*price(i)))
(total:IntVar) #= sum(Array.tabulate(num_prices)(i=> x(i)*price(i)))


Often, this is not much a problem, but it can be

abs (decomposition with reification)

Update 2011-08-11: The latest version of scalaJaCoP includes a definition of abs, so this comment (and the defined function) is not needed anymore.

In AllIntervals.scala, I needed the abs constraint but I couldn't find any in scalaJaCoP. Here is a brutal version using reifications, hardwired to handle the case of c=abs(a-b):
def my_abs(a:IntVar, b:IntVar, c:IntVar) {
  val b1 = new BoolVar("b1")
  b1 <=> (a #> b)
  b1 -> (c #= a - b)
  ~b1 -> (c #= b - a)
}
(This can most probably be written both more elegant and more general.)

Number of solutions/Timeout

Among the things that where added to the new version of scalaJaCoP was support for timeout and for getting a specific number of solutions, both very useful, especifally for debugging.
  • numberSolutions(n) where n is the number of solutions to print
  • timeOut(s) where s is in seconds

My JaCoP/Scala models

All my public scalaJaCoP models are collected in my JaCoP/Scala page but listed here for easy access. Most of them have also been implemented in another CP systems (see Common constraint programming problems).

April 26, 2011

JaCoP version 3.1 is released

JaCoP version 3.1 has been released. Download information is here. (Sourceforge)

From the announcement:
Dear users,

We have just released JaCoP 3.1. This is a release that fixes few bugs as well as provides a new Binpacking constraint. The Binpacking constraint allows modeling of problems exhibiting bin packing problem structure. The addition of this constraint introduced changes to minizinc/fz library in JaCoP since the definitions of bin_packing.mzn, bin_packin_capa.mzn and bin_packing_load.mzn have been added. The JaCoP guide has been updated to explain how to use Binpacking constraint.

best regards, Core Developer Team
Don't forget that there is a SVN version as well.

February 07, 2011

JaCoP version 3.1 is released

JaCoP version 3.1. has been released. Download from SourceForge.

From the announcement:
Dear users, We have just released JaCoP 3.1. This is a release that fixes few bugs as well as provides a new Binpacking constraint. The Binpacking constraint allows modeling of problems exhibiting bin packing problem structure. The addition of this constraint introduced changes to minizinc/fz library in JaCoP since the definitions of bin_packing.mzn, bin_packin_capa.mzn and bin_packing_load.mzn have been added. The JaCoP guide has been updated to explain how to use Binpacking constraint. best regards, Core Developer Team

November 30, 2010

JaCoP version 3.0 (final) is released

JaCoP version 3.0 (final) is released. Download via Sourceforge. From the release message:

Dear users,

We have just released JaCoP 3.0 final. Since the previous release (RC2) we have fixed one rarely triggering bug in ElementInteger constraint, as well as added code supporting generation of data for CP-viz framework. Here is the list of most important changes since version 2.4.

1) The introduction of Network Flow constraint, which allows efficient modeling of problems with network flow structure. The constraint allows to specify minimum cost flow problem in terms of finite domain variables. It uses network simplex algorithm as an propagator for network-like structure. It is a very expressive constraint where arc flow and arc unit costs can be specified by variables.

2) The introduction of set package forced changes in the design of core interfaces. There are IntVar and SetVar classes which inherit from Var class. This change allowed to refactor and improve set package so it is designed in cleaner and more efficient manner.

3) The introduction of special domain representation SmallDenseDomain, which uses bits within a long primitive type to encode small domains of range not larger than 64.

4) The introduction of Soft-alldifferent as well as Soft-GCC global constraints. The soft-alldifferent constraint implements Variable and Decomposition based violation metrics. The soft-gcc constraint implements value based violation metric. Both are decomposed constraints based on network flow constraint.

5) Examples have been updated by moving multiple solving functionality from main function to test function, so user can easily see what is the best model just by looking at main function. BIBD example has been added. Examples with Set variables have been updated to reflect the changes.

6) A number of bug fixes and changes in flatzinc interface to better comply with minizinc 1.1. We have also added into minizinc predicates networkflow that uses newly introduced JaCoP Network Flow constraint.

7) A number of minor API changes to improve intuitiveness of use (e.g. order of arguments in constructors).

8) The JaCoP guide has been updated to reflect the changes and additions to the newest version.

best regards,

Core Developer Team

Also see:

September 07, 2010

My updated models for JaCoP Version 3.0 (RC1)

As mentioned in JaCoP version 3.0 (RC1) is released , JaCoP has a new version 3.0 (as of writing - Release Candidate 1). It can be downloaded via SourceForge.

I have now changed my JaCoP models so they compile and run with version 3.0.

Changes of the models

In all models some of these changes was done:
  • FDStore is now Store
  • Variable is now IntVar (or Var)
  • FDV is now IntVar (or Var)
  • int [] res = label.getSolutionListener().getSolution(s); is now written as
    Domain [] res = label.getSolutionListener().getSolution(s)
    Printing the result then works as before.
  • Min and Max has changed the order of the parameters: store.impose(new Min(x[0], x)) is now written store.impose(new Min(x, x[0])) (and the same for Max)

JaCoP Models

Here are all my JaCoP models, and their data files. Some of these use the helper class HakankUtil.java.

Note: Most of these models where written some years ago. In the JaCoP distribution, directory ExamplesJaCoP, there many more examples that shows the new features for version 3.0. (Some of these are my models with an updated touch.)

September 01, 2010

JaCoP version 3.0 (RC1) is released

JaCoP version 3.0 (RC1) is released. Download this version from JaCoP's SourceForge page.

From the announcement JaCoP version 3.0 (RC1) is released :

Dear users,

We make available JaCoP 3.0 RC1, which has been used in Minizinc Challenge 2010. We are still working on the newest functionality (1) as well as cleaning API to reduce the need of further changes in the future. We look for your feedback so we can make the final release 3.0 as soon as possible. Here is the list of most important changes since version 2.4.

1) The introduction of Network Flow constraint, which allows efficient modeling of problems with network flow structure. The constraint allows to specify minimum cost flow problem in terms of finite domain variables. It uses network simplex algorithm as an propagator for network-like structure. It is a very expressive constraint where arc flow and arc unit costs can be specified by variables.

2) The introduction of set package forced changes in the design of core interfaces. There are IntVar and SetVar classes which inherit from Var class. This change allowed to refactor and improve set package so it is designed in cleaner and more efficient manner.

3) The introduction of special domain representation SmallDenseDomain, which uses bits within a long primitive type to encode small domains of range not larger than 64.

4) The introduction of Soft-alldifferent as well as Soft-GCC global constraints. The soft-alldifferent constraint implements Variable and Decomposition based violation metrics. The soft-gcc constraint implements value based violation metric. Both are decomposed constraints based on network flow constraint.

5) Examples have been updated by moving multiple solving functionality from main function to test function, so user can easily see what is the best model just by looking at main function. BIBD example has been added. Examples with Set variables have been updated to reflect the changes.

6) A number of bug fixes and changes in flatzinc interface to better comply with minizinc 1.1. We have also added into minizinc predicates networkflow that uses newly introduced JaCoP Network Flow constraint.

7) A number of minor API changes to improve intuitiveness of use (e.g. order of arguments in constructors).

8) The JaCoP guide has been updated to reflect the changes and additions to the newest version.

best regards,

Core Developer Team

JaCoP is now available from SourceForge

Earlier today JaCoP became available from SourceForge: jacop-solver. The summary is


JaCoP is a Java Constraint Programming solver. It provides a significant number of (global) constraints to facilitate efficient modeling of combinatorial problems, as well as modular design of search.

It's great that JaCoP now is available via SVN and to be able to preview features in upcoming versions. There are some interesting things coming up, for example the network flow constraint that was mentioned some months ago at the JaCoP site: Network Flow Constraint .

Here are some more links:

Related: My JaCoP page.

March 03, 2010

Survey of Nonogram Solvers updated: one addition is a JaCoP solver

Survey of Paint-by-Number Puzzle Solvers by Jan Wolter has been updated. One interesting addition to the solver is a JaCoP solver: Ben Raziel's JaCoP Solver:


[from Ben Raziel:]
"""
We use JaCoP (a free java CSP solving library) to do the line solving for us. For searching, we use a modified version of JaCoP's DFS search. I employed your probing approach and tweaked it a little in order to achieve better performance.

The basic idea is to probe more cells, in order to reach a contradiction - which means we don't have to guess our next move (since a contradiction tells us what color the current cell is). The probes are ordered into "levels": each level contains cells where a contradiction is more likely to be found - which minimizes the number of cells that we have to probe until a contradiction is reached.
"""
The solver is not currently publically available, but the author plans to make it available, once this is cleared with the University.

From the summary of this solver

Assessment: Very good at hard puzzles, a bit slow at easy ones.

...

The run times on simple problems are not spectacular. Java is an interpreted language, and thus, naturally, rather slow. It is furthermore also built on top of a general solving library, JaCoP, and those typically add some overhead. But my suspicion is that part of this is also a basic design trade off - the authors were focusing more on solving hard puzzles fast than on solving easy puzzles fast.

But this really works very well on the harder puzzles. It solves every puzzle in the full 2,491 webpbn puzzle data set in under 15 minutes and only the "Lion" puzzle, which most solvers can't solve at all, takes more than 5 minutes. No other solver tested has accomplished this.

Of course, this solver, like pbnsolve, was trained on the webpbn data set, which certainly gives it an advantage when being tested against that data set. In fact the authors had the full dataset available from early in their development cycle, which is more than can be said for webpbn.

There are also other updates, such as general reflections about the state of Nonogram solvers.

November 08, 2009

Update on Nonogram: Jan Wolter's Survey and my own new benchmark

Survey of Paint-by-Number Puzzle Solvers

In Some new models, etc I mentioned the great Survey of Paint-by-Number Puzzle Solvers, created by Jan Wolter (also author of the Nonogram solver pbnsolve).

In this survey he included both Gecode's Nonogram solver written by Mikael Lagerkvist as well as my own Nonogram model (with Gecode/FlatZinc).

Since the last update on the former blog post the following has happeded:
  • both our solver has now the assessment: "An amazingly good solver, especially for a simple demo program", and are placed 4, 5, and 6 of 10 tested systems
  • my Gecode/FlatZinc model has been tested for "Full results"; it came 4 out of 5.
  • my Nonogram model with the Lazy FD solver is now included in the "Sample result", at 6'th place
It seems than Wolter has appreciated constraint programming as a general tool for solving these kind of combinatorial problems, much for its ease of experimentation, e.g. with labeling strategies and (for the MiniZinc models) changing solvers:

From the analysis of Lagerkvist's Gecode model:
This is especially impressive because the faster solvers are large, complex programs custom designed to solve paint-by-number puzzles. This one is a general purpose solver with just a couple hundred custom lines of code to generate the constraints, run the solver, and print the result. Considering that this is a simple application of a general purpose solving tool rather than hugely complex and specialized special purpose solving tool, this is an astonishingly good result.

Getting really first class search performance usually requires a lot of experimentation with different search strategies. This is awkward and slow to do if you have to implement each new strategies from scratch. I suspect that a tool like Gecode lets you try out lots of different strategies with relatively little coding to implement each one. That probably contributes a lot to getting to better solvers faster.
From the analysis of my MiniZinc model:
If you tried turning this into a fully useful tool rather than a technology demo, with input file parsers and such, it would get a lot bigger, but clearly the constraint programming approach has big advantages, achieving good search results with low development cost.

...

These two results [Gecode/FlatZinc and LazyFD] highlight the advantage of a solver-independent modeling language like MiniZinc. You can describe your problem once, and then try out a wide variety of different solvers and heuristics without having to code every single one from scratch. You can benefit from the work of the best and the brightest solver designers. It's hard to imagine that this isn't where the future lies for the development of solvers for applications like this.
And later in the Conclusions:
The two constraint-based systems by Lagerkvist and Kjellerstrand come quite close in performance to the dedicated solvers, although both are more in the category of demonstrations of constraint programming than fully developed solving applications. The power of the underlaying search libraries and the ease of experimentation with alternative search heuristics obviously serves them well. I think it very likely that approaches based on these kinds of methods will ultimately prove the most effective.
I think this is an important lesson: Before starting to write very special tools, first try a general tool like a constraint programming system and see how well it perform.

The Lazy FD solver and the Lion problem

Most of the problems in the Sample Results where solved by some solver within the time limit of 30 minutes. However, one problem stand out as extra hard: The Lion problem. When I tested the MiniZinc's Lazy FD solver on my machine I was very excited that it took just about 2 minutes, and mentioned this to Wolter. He also tested this but on his 64-bit machine it took 80 minutes to solve (and since it was above the time limit this is not in the result table). This is how he describes the Lazy FD solver:
But the remarkable thing is that [the Lazy FD solver] solves almost everything. Actually, it even solved the Lion puzzle that no other solver was able to solve, though, on my computer, it took 80 minutes to do it. Now, I didn't allow a lot of other solvers 80 minutes to run, but that's still pretty impressive. (Note that Kjellerstrand got much faster solving times for the Lion than I did. Lagerkvist also reported that his solver could solve it, but I wasn't able to reproduce that result even after 30 CPU hours. I don't know why.)
After some discussion, we come to the conclusion that the differences was probably due to the fact that I use a 32-bit machine (and the 32-bit version of MiniZinc) with 2 Gb memory, and Wolter use a 64-bit machine with 1 Gb memory.

One should also note that the all other solvers was compiled without optimization, including Gecode/FlatZinc; however the LazyFD does not come with source so it is running optimized. This may be an unfair advantage to the LazyFD solver.

My own benchmark of the Sample Results

The times in the Sample Results is, as mentioned above, for solvers compiled with no optimization. I have now run the same problems on my machine (Linux Mandriva, Intel Dual 3.40GHz, with 2Gb memory), but the solvers uses the standard optimization. All problems was run with a time limit of 10 minutes (compared to Wolters 30 minutes) and searched for 2 solutions, which checks for unique solutions. The last three problems (Karate, Flag, Lion) has multiple solutions, and it is considered a failure if not two where found in the time limit. I should also note that during the benchmark I am using the machine for other things, such as surfing etc.

The problems
I downloaded the problems from Wolter's Webbpbn: Puzzle Export. For copyright reasons I cannot republish these models, but it is easy to download each problem. Select ".DZN" for the MiniZinc files, and "A compiled in C++ format" for Gecode. There is no support for Comet's format, but it's quite easy to convert a .dzn file to Comet's.

The solvers + labeling strategies
Here is a description of each solver and its labeling strategy:
  • fz, "normal" (column_row)
    MiniZinc model with Gecode/FlatZinc. The usual labeling in nonogram_create_automaton2.mzn, i.e. where the columns are labeled before rows:
    solve :: int_search(
          [x[i,j] | j in 1..cols, i in 1..rows], 
          first_fail, 
          indomain_min, 
          complete) 
    satisfy;
    
  • fz, "row_column"
    MiniZinc model with Gecode/FlatZinc. Here the order of labeling is reversed, rows are labeled before columns. Model is nonogram_create_automaton2_row_column.mzn
    solve :: int_search(
          [x[i,j] | i in 1..rows, j in 1..cols], 
          first_fail, 
          indomain_min, 
          complete) 
    satisfy;
    
  • fz, "mixed"
    MiniZinc model with Gecode/FlatZinc: nonogram_create_automaton2_mixed.mzn.
    I have long been satisfied with the "normal" labeling in the MiniZinc model because P200 (the hardest problem I until now has tested) was solved so fast. However, the labeling used in the Comet Nonogram model described in Comet: Nonogram improved: solving problem P200 from 1:30 minutes to about 1 second, and which is also used in the Gecode model, is somewhat more complicated since it base the exacl labeling by comparing the hints for the rows and the column.

    I decided to try this labeling in MiniZinc as well. However, labeling in MiniZinc is not so flexible as in Comet and Gecode. Instead we have to add a dedicated array for the labeling (called labeling):
    array[1..rows*cols] of var 1..2: labeling;
    
    and then copy the element in the grid to that array based on the relation between rows and column hints:
    constraint
          % prepare for the labeling
          if rows*row_rule_len < cols*col_rule_len then
               % label on rows first
               labeling = [x[i,j] | i in 1..rows, j in 1..cols]
          else 
               % label on columns first
               labeling = [x[i,j] | j in 1..cols, i in 1..rows]
          endif
          /\
          % .... 
    
    and last, the search is now based just on this labeling array:
    solve :: int_search(
            labeling,
            first_fail, 
            indomain_min, 
            complete)
    satisfy;
    
  • jacop, "normal"
    MiniZinc normal model with JaCoP/FlatZinc using the same model as for fz "normal".
  • lazy, satisfy
    Model: nonogram_create_automaton2_mixed.mzn. This use the MiniZinc LazyFD solver with the search strategy:
    solve satisfy;
    
    This labeling is recommended by the authors of LazyFD. See MiniZinc: the lazy clause generation solver for more information about this solver.

    Note: The solver in MiniZinc latest official version (1.0.3) don't support set vars. Instead I (and also Jan Wolter) used the latest "Release Of The Day" version (as per 2009-11-02).
  • Comet, normal
    Model: nonogram_regular.co. This is the Comet model I described in Comet: Nonogram improved: solving problem P200 from 1:30 minutes to about 1 second. No changes has been done.
  • Gecode, normal
    This is the Nonogram model distributed with Gecode version 3.2.1. The labeling is much like the one used in the Comet model, as well as fz, "mixed". (In fact the labeling in the Gecode model was inspired by the labeling in the Comet model).


Here is the results. For each model (+ labeling strategy) two values are presented:
  • time (in seconds)
  • number of failures if applicable (the LazyFD solver always returns 0 here).
The result
Model fz
normal
fz
row_column
fz
mixed
jacop
normal
lazy
satisfy
Comet
normal
Gecode
normal
Dancer
(#1)
0.48s
(0)
0.31s
(0)
1.00s
(0)
3.64s
(0)
0.91s
(0)
0.691s
(0)
0.199s
(0)
Cat
(#6)
0.24s
(0)
0.24s
(0)
0.25s
(0)
1.20s
(0)
1.13s
(0)
0.6s
(0)
0.025s
(0)
Skid
(#21)
0.24s
(13)
0.23s
(3)
0.28s
(13)
0.78s
(13)
1.37s
(0)
0.586s
(0)
0.217s
(0)
Bucks
(#27)
0.32s
(3)
0.32s
(9)
0.37s
(3)
0.96s
(3)
2.37s
(0)
1.366s
(5)
0.026s
(2)
Edge
(#23)
0.16s
(25)
0.17s
(25)
0.18s
(25)
0.59s
(25)
0.31s
(0)
0.521s
(43)
0.175s
(25)
Smoke
(#2413)
0.27s
(5)
0.27s
(8)
0.28s
(8)
0.83s
(5)
1.44s
(0)
0.616s
(14)
0.275s
(5)
Knot
(#16)
0.42s
(0)
0.43s
(0)
0.48s
(0)
1.19s
(0)
8.15s
(0)
1.307s
(0)
0.329s
(0)
Swing
(#529)
0.95s
(0)
0.94s
(0)
0.96s
(0)
2.19s
(0)
21.94s
(0)
1.782s
(0)
0.431s
(0)
Mum
(#65)
0.64s
(20)
0.64s
(22)
0.66s
(22)
1.68s
(20)
16.34s
(0)
1.268s
(39)
0.491s
(22)
Tragic
(#1694)
340.32s
(394841)
1.02s
(255)
436.92s
(394841)
--
(198329)
45.97s
(0)
477.39s
(702525)
1.139s
(255)
Merka
(#1611)
--
(361587)
1.44s
(79)
--
(294260)
--
(136351)
80.92s
(0)
1.654s
(46)
0.645s
(13)
Petro
(#436)
2.97s
(1738)
143.09s
(106919)
3.42s
(1738)
7.27s
(1738)
9.86s
(0)
3.103s
(3183)
151.09s
(106919)
M_and_m
(#4645)
1.41s
(89)
601.27s
(122090)
1.59s
(89)
3.43s
(89)
66.98s
(0)
2.215s
(162)
2.797s
(428)
Signed
(#3541)
1.87s
(929)
23.12s
(6484)
28.23s
(6484)
5.75s
(929)
73.02s
(0)
20.369s
(12231)
1.648s
(929)
Light
(#803)
600.47s--
(400660)
--
(621547)
--
(485056)
601.53s--
(171305)
8.64s
(0)
--
(0)
--
(538711)
Forever
(#6574)
4.14s
(17143)
7.86s
(30900)
6.22s
(17143)
12.21s
(17143)
3.27s
(0)
7.5s
(27199)
8.077s
(30900)
Hot
(#2040)
--
(303306)
--
(330461)
--
(248307)
--
(119817)
165.72s
(0)
--
(0)
--
(312532)
Karate
(#6739)
95.78s
(215541)
67.27s
(130934)
133.43s
(215541)
373.02s
(215541)
19.32s
(0)
120.02s
(272706)
80.56s
(170355)
Flag
(#2556)
--
(1686545)
5.69s
(14915)
7.93s
(14915)
--
(243222)
9.29s
(0)
7.28s
(24678)
3.998s
(16531)
Lion
(#2712)
--
(542373)
--
(1124697)
--
(420216)
--
(187215)
115.56s
(0)
--
(0)
--
(869513)

Some conclusions, or rather notes

Here are some conclusions (or notes) about the benchmark.
  • The same solver Gecode/FlatZinc is here compared with three different labelings. There is no single labeling that is better than the other. I initially has some hopes that the "mixed" labeling should take the best labeling from the two simpler row/columns labelings, but this is not really the case. For example for Tragic the row_column strategy is better than "normal" and "mixed". I am, however, somewhat tempted, to use the "row_column" labeling, but the drawback is that "P200" problem (not included in Wolfter's sample problems) takes much longer with this labeling.
  • The same model and labeling but with different solvers is compared: Gecode/FlatZinc is faster than JaCoP/FlatZinc on all the problems. For the easier problems this could be explained by the extra startup time of Java for JaCoP, but that is not the complete explanation for the harder problems. Note: Both Gecode/FlatZinc and JaCoP/FlatZinc has dedicated and fast regular constraints (whereas the LazyFD, and the Comet solvers use a decomposition).
  • The LazyFD solver is the only one that solves all problems (including Lion), but is somewhat slower on the middle problems than most of the others. It emphasizes that this is a very interesting solver.
  • It is also interesting to compare the result of the Comet model and Gecode/FlatZinc "mixed", since they use the same principle of labeling. However there are some differences. First, the MiniZinc model with Gecode/FlatZinc use a dedicated regular constraint, and Comet use my own decomposition of the constraint. For the Merka problem the Comet version outperforms the Gecode/FlatZinc version, otherwise it's about the same time (and number of failures).
  • The Light problem: It is weird that this problem was solved in almost exactly 10 minutes (the timeout is 10 minutes) for Gecode/FlatZinc and JaCoP/FlatZinc. The solutions seems correct but I'm suspicious of this. Update: Christian Schulte got me on the right track. Here is was happened: The first (unique) solution was found pretty quick and was printed, but the solvers could not prove a unique solution so it timed out. JaCoP/FlatZinc actually printed "TIME-OUT" but I didn't observe that. Case closed: They both FAILED on this test. Thanks, Christian.End update
As said above, I can only agree with Jan Wolter in his comment that the ease of experimenting, for example changing labeling, and solver for the FlatZinc solvers, is a very nice feature.

Last word

No benchmark or comparison of (constraint programming) models is really complete without the reference to the article On Benchmarking Constraint Logic Programming Platforms. Response to Fernandez and Hill's "A Comparative Study of Eight Constraint Programming Languages over the Boolean and Finite Domains" by Mark Wallace, Joachim Schimpf, Kish Shen, Warwick Harvey. (Link is to ACM.)

From the Abstract:
... The article analyses some pitfalls in benchmarking, recalling previous published results from benchmarking different kinds of software, and explores some issues in comparative benchmarking of CLP systems.

October 01, 2009

MiniZinc: 151 new Nonogram problem instances (from JaCoP)

The latest versions of JaCoP (download here) includes 151 Nonogram instances (ExamplesJaCoP/nonogramRepository/). I wrote about these and the included Nonogram solver in JaCoP: a request from the developers (Knapsack and Geost) and Nonogram labeling.

When I (beta) tested the new FlatZinc support for this version, I converted these instances to MiniZinc's data file (.dzn) for some tests. It is a fun way of testing a system.

Now all these problem instances have been published. They are to be used with the MiniZinc Nonogram model nonogram_create_automaton2.mzn. See At last 2: A Nonogram solver using regular written in "all MiniZinc" for more about this model.

All problem instances - in MiniZinc's .dzn format - are available in www.hakank.org/minizinc/nonogram_examples, and also packed in the Zip file nonogram_examples.zip. The name of each file corresponds to the files in JaCoP's ExamplesJaCoP/nonogramRepository/; the file data000.nin is here called nonogram_jacop_data000.dzn, etc. (A larger file, nonogram_examples_with_fzn.zip also includes the generated .fzn files. Note: I used mzn2fzn for MiniZinc version ROTD/2009-09-13 for this. See below how to generate these files.)

Batch running with JaCoP/Fz2jacop

Running many FlatZinc models with JaCoP's Fz2jacop is not optimal because of the overhead of the Java startup. Instead I used a Java program, JaCoPFznTest.java (based on a program from Krzysztof Kuchcinski, one of JaCoP's main developers. Thanks, Kris.).

The main call is

fz.main(new String[] {"-s", "-a", m});

which means that statistics should be shown and also require to return all solutions of the problem.

The FlatZinc (.fzn) files was generated from the .dzn files with the following command:

foreach_file 'mzn2fzn -G jacop --data $f -o $f.fzn nonogram_create_automaton2.mzn' '*.dzn'

ExamplesJaCoP.Nonogram vs. JaCoP/Fz2jacop

This time I just compared the run time of the two different approaches for solving Nonograms with JaCoP:
* The Java program ExamplesJaCoP.Nonogram
* Running JaCoP's MiniZinc/FlatZinc solver Fz2jacop

The program instances is, with one exception the same as in the distributed ExamplesJaCoP.Nonogram. Since I wanted both program to solve for all solutions, I excluded instance #83 (see below) since it has too many solutions. However, the P200 problem in included since it is hardwired in ExamplesJaCoP.Nonogram. This means that it it still 151 problems running.

The result:
* ExampleJaCoP.Nonogram took 14.8 seconds
* The Fz2jacop version took 17.8 seconds

I'm quite impressed with both of these results, especially Fz2jacop's As shown in JaCoP: a request from the developers (Knapsack and Geost) and Nonogram labeling, many problems are solved in "0" milliseconds, but still.

Instance #83

As mentioned above, problem instance #83 was excluded since it has too many solutions. Here are two of them. It looks like a person (alien?) standing on a spaceship waving.
                              #           #         
                            #       #               
                      # # # #       # # #           
                          #   # # #   # # #         
                  #                   #     #       
                #         #     # # # # #           
              #     #   #           # # #     #     
              #   #   #             # # #   #       
              # #   #               # # #           
              # #                   #   #           
                #                   #   #           
          # # # # # # # # # # # # # # # # #         
      # # # # # # # # # # # # # # # # # # # # #     
    # # # #     # #     # # #     # #     # # # #   
  # # # # #     # #     # # #     # #     # # # # # 
    # # # #     # #     # # #     # #     # # # #   
      # # # # # # # # # # # # # # # # # # # # #     
          # # # # # # # # # # # # # # # # #         
              # # # # # # # # # # # # #             
                  #       #       #                 
                # #       #       # #               
                #         #         #               
              # #         #         # #             
              #           #           #             
            # # #       # # #       # # #           
----------
                              #           #         
                            #       #               
                      # # # #       # # #           
                          #   # # #   # # #         
                  #                   #     #       
                #         #     # # # # #           
              #     #   #           # # #   #       
              #   #   #             # # #     #     
              # #   #               # # #           
              # #                   #   #           
                #                   #   #           
          # # # # # # # # # # # # # # # # #         
      # # # # # # # # # # # # # # # # # # # # #     
    # # # #     # #     # # #     # #     # # # #   
  # # # # #     # #     # # #     # #     # # # # # 
    # # # #     # #     # # #     # #     # # # #   
      # # # # # # # # # # # # # # # # # # # # #     
          # # # # # # # # # # # # # # # # #         
              # # # # # # # # # # # # #             
                  #       #       #                 
                # #       #       # #               
                #         #         #               
              # #         #         # #             
              #           #           #             
            # # #       # # #       # # #           
----------
This output was generated with the (very new) option {noresult} in mzn_show.pl. (This option don't display the real numeric results.)

September 16, 2009

JaCoP version 2.4 released

From JaCoP's news JaCoP version 2.4 is released:
Dear all,

We are happy to announce the release of a new version of our
Java-based solver JaCoP. This new version 2.4 has a number of new
features in addition to some bug fixes. The most important additions
in this version are:

1. The flatzinc interface that makes it possible to run minizinc
programs using JaCoP. The distribution contains number of different
minizinc examples.

2. Geometrical constraint, geost, based on pruning algorithms
originally proposed by Nicolas Beldiceanu et al. This constraint makes
it possible to define placement problems of non-convex objects in
k-dimensional space.

3. Knapsack constraint, which is based on the work published by Irit
Katriel et al. We extend the original work in number of ways, for example by
making it possible to use non-boolean quantity variables.

4. Set constraints defining typical operation on sets using set
interval variables.

This work would not be possible without help of several people. We
would like to thank Hakan Kjellerstrand for his help in testing
flatzinc parser as well as providing a number of examples in minizinc
format. We would also like to thank Meinolf Sellmann for his comments
on the initial implementation of knapsack constraint which have helped to
improve it further. Marc-Olivier Fleury has implemented most of the
functionality of the geost constraint. Robert Åkemalm has implemented the
first version of set constraint package. Wadeck Follonier has implemented
the first version of the Knapsack constraint. We would like to thank them
for their contributions.

As always feel free to send us ( radoslaw [dot] szymanek [at] gmail [dot] com )
feedback. We are always looking for cooperation to improve JaCoP. If you miss
some functionality of JaCoP we can help you to develop it so it can be done
efficiently and fast.


best regards,
Radoslaw Szymanek and Kris Kuchcinski
The latest version of JaCoP can be downloaded here.

For more information, see:
Also see My JaCoP page.

MiniZinc/FlatZinc support

I have especially tested the FlatZinc solver (Fz2jacop) and it is fast. For example, here is the statistics for nonogram_create_automaton2.mzn with the P200 problem (I wrote about this some days ago here):
Model variables : 629
Model constraints : 50

Search CPU time : 640ms
Search nodes : 1040
Search decisions : 520
Wrong search decisions : 520
Search backtracks : 520
Max search depth : 22
Number solutions : 1

Total CPU time : 1010ms
Note: JaCoP uses a special optimized regular constraint.

The FlatZinc solver has the following options:
$ java JaCoP.fz.Fz2jacop --help
Usage: java Fz2jacop [] .fzn
Options:
    -h, --help
        Print this message.
    -a, --all-solutions
    -t , --time-out 
         - time in second.
    -s, --statistics
    -n , --num-solutions 
         - limit on solution number.
Great work!

August 23, 2009

Java Specification Requests JSR 331: Constraint Programming API

The Java Specification Request JSR 331: Constraint Programming API is a request for consolidating the Java constraint solvers.
This specification defines a Java runtime API for constraint programming. The CP API prescribes a set of fundamental operations used to define and solve constraint satisfaction and optimization problems.

The request

2.1 Please describe the proposed Specification:

This specification defines a Java runtime API for constraint programming. The standardized CP API aims to make CP technology more accessible for business software developers. It intends to create CP Standards that will define unified business interface(s) for modeling and solving real-world decision support problems as Constraint Satisfaction and Optimization Problems. Having unified interfaces will allow commercial application developers to model their problems in such a way that the same model can be tried with different CP Solvers. This will minimize vendor dependence without limiting vendor's innovation. At the same time the CP API will help to bring the latest research results to real-world business applications. The CP API prescribes a set of fundamental operations used to define and solve constraint satisfaction and optimization problems. The sets of operations will allow a user to express its decision support problem in terms of:
- Decision Variables defined as constrained integer, real, and set variables
- Constraints defined on these variables using the standardized unary, binary and global constraints.
Clearly separating the problem definition from problem resolution, the sets of operations will allow a user to apply standardized search strategies to find feasible and optimal solutions of the problem. The CP API will also allow a user to create custom constraints and search strategies in a unified way making them independent from a particular CP solver.

...

2.6 Why isn't this need met by existing specifications?

Today CP is a proven optimization technique and many CP solvers empower real-world business applications in such areas as scheduling, planning, configuration, resource allocation, and real-time decision support. However, the absence of standards still limits the acceptance of CP by the business world.

While dissimilar vendor-specific CP API specifications exist, their syntactic and semantic differences are significant enough to cause costly difficulties for application builders and platform vendors.

The development schedule

2.13 Please describe the anticipated schedule for the development of this specification.

1. Early Draft Review: November 20, 2009

2. Complete the Specification, Public Draft/Final Release: April 15, 2010

3. Maintenance Review: June 30, 2010

2.14 Please describe the anticipated working model for the Expert Group working on developing this specification.

The JSR Specification Lead has already contacted major CP vendors and well-known CP experts and involved them in a preliminary discussion at www.cpstandards.org. We will invite them to join the Expert Group right after the initial JSR Review and acceptance. We plan to have an kick-off meeting during the oncoming major CP conference CP-2009 that will take place at Portugal on 20-24 September 2009. A preliminary Java CP API that will become a basis for the JSR is being discussed by CP experts at the Discussion Forum at www.cpstandards.org. During the first 3 months we plan to finalize the Expert Group that will include the most active contributors to the Early Draft.

Related

Other pages of interest:
* JSR #331 Constraint Programming API - JSR Review Ballot * Public forum for JSR 331.
* Standardization of Constraint Programming (www.cpstandards.org), forum with some discussion about the JSR

Also, see:
* My Choco page
* My JaCoP page
* My constraint programming page

August 16, 2009

JaCoP: a request from the developers (Knapsack and Geost) and Nonogram labeling

Here are some JaCoP related news.

A request from the developers: Knapsack and Geost

This is a request from the JaCoP developers, from the article Geost and Knapsack constraint:
Dear all, We have tested intiial but already optimized implementations of Geost constraint as well as Knapsack constraint. We need testers who will use and test the constraint in their problems. Feel free to contact me by email ( radoslaw [dot] szymanek {at} gmail (dot) com) to get source code for any of those two constraints. Geost constraint has been shown to be even 1000 times faster for problems of size 1000+ objects than other open-source implementation available. Knapsack constraint is also very efficient as it uses the sublinear algorithm presented in the literature. We have extended Knapsack constraint to be able to handle 0..n quantity variables as well knapsack capacity and profit specified as finite domain variables and not just integers. Geost constraint is a result of Master thesis done by Marc-Olivier Fleury. Knapsack constraint is based on the student project completed by Wadeck Follonier. best regards, Radoslaw Szymanek

Nonogram

In the latest version of JaCoP (2.3.2, from May 2009) there are 151 test instances for the Nonogram model. One of these is the well known (infamous) problem P200, which I have written about before, for example in Comet: Nonogram improved: solving problem P200 from 1:30 minutes to about 1 second, where I - aided by Pascal Van Hentenryck - happened to find a labeling that was quite good for the P200 problem (and not very bad for the other instances).

The labeling used in the JaCoP distribution of Nonogram.java is (according to the source) the "Zigzag based variable ordering". I thought it would be fun to compare it with the labeling used in the Comet model (also later used in the Gecode model).

The "new" labeling

The "new" labeling described in Comet: Nonogram improved: solving problem P200 from 1:30 minutes to about 1 second is simply to order the variables according to the length of the row rules and column rules:
if (board.length*row_rules.length <  board[0].length*col_rules.length) {
    // i first
    for (int i = 0; i < row_rules.length; i++) {
	for (int j = 0; j < col_rules.length; j++) {
            vars.add(board[i][j]);
	}
    }
} else {
    // j first
    for (int j = 0; j < col_rules.length; j++) {
        for (int i = 0; i < row_rules.length; i++) {
            vars.add(board[i][j]);
	}
    }
}
Here are the statistics of running the P200 problem with the two different labelings. P200 is the first (and hard coded) problem presented in Nonogram.java.

Original labeling:
Nodes : 6472
Decisions : 3236
Wrong Decisions : 3236
Backtracks : 3236
Max Depth : 38

...

Execution time = 3135 ms
The "New" labeling:
Nodes : 1040
Decisions : 520
Wrong Decisions : 520
Backtracks : 520
Max Depth : 22

...

Execution time = 638 ms
So, for P200 this new labeling was an improvement: nearly 5 times better time and about 6 times fewer backtracks.

Running all instances

But how good is the new labeling overall? As noted above, there are 150 test instances (plus P200). I ran all of them and compared the individual results with the two labelings. Note that the instance 83 is the "alien waving problem" instance which has a lot of solutions; I changed the code so it just finds one solution. For the rest of the instances the solving time is for finding all solutions (which happens to be unique).

All the times below are the reported time from the Nonogram model, excluding time for printing the model and the solution. It should also be noted that the time is from one single run.

Total solving time
Running all the 150 + 1 problem instances, the total solving time was:
Total original: 6102
Total new     : 5462
Total diff    : 640
The total of differences means that the newer labeling has better (total) solving time, about 0.5 second.

Individual instances
For a complete listing of the test instances see nonogram_all.txt. Here we will concentrate on a few cases.

Below is the times of the instances where the differences between the two labelings where larger than 10ms, what I call "dominating" instances:
38: orig:34 new:21 diff:13
51: orig:11 new:35 diff:-24
67: orig:35 new:10 diff:25
84: orig:1222 new:963 diff:259
91: orig:872 new:2800 diff:-1928
p200: orig:3006 new:700 diff:2306
As we can see, there are three problem that dominates the total times: 84, 91, and P200. Note that for one of these (instance #91) the new labeling is actually worse.

Here are all differences grouped depending on which labeling was the best. Last (in parenthesis) is the number of instances in each group:
original_best diffs: 
3 1 1 1 10 2 3 2 2 1 1 3 3 2 2 2 1 1 24 1 3 4 2 2 2 1 1 9 2 1 1928 1 1 1 2 1 1 1 2 3 1 2 1 2 1 1 2 1 3 (49)
new_best diffs     : 
1 3 1 2 1 3 4 3 3 4 1 3 1 1 2 1 13 1 1 2 2 1 3 1 1 25 2 1 3 259 1 3 3 1 2 2 3 3 3 1 1 3 1 1 1 3 1 2306  (48)
same               : 
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 (55)

Summary

To summarise: the "new" labeling is overall quite good, but not better than the original labeling in all cases. It was better for two of the three dominating problem instances, and worse on one. Well, as always, there is no silver bullet.

The dominating instances

I don't know the name of the two dominating instances #84 and #91 so here they are, together with the infamous P200, using the output from JaCoP's Nonogram model.
Problem P200
 00   00     000
0000  0 0    0 0000
 0000 0 00   0     0
0000  0 0 0  0  0  0
 00   0 0 00  000  00000
      0   0 0   0  00  0
     000  0 00000  0  00
   000 00  00    0 00  00
  00 0  0000   0 0 0    0
 00      00 0 00 0     00
 0  0     0 000 00   000
 0  0 00  0000000  000
 0 00     00 0 00000
 000  00 00  0  00
  000   00   0   00
    00000    0    00
  00   00    0     00
 0000   00   0    00
000000   00  000 00
0000000   0000 000  00
 0000000    0000   0000
0000000      0      0000
000000       0     0000
 0000       00      00
  00         0
Problem 84
                 0
     0           00
     00         000
     000       00 0
     0 00     000 0
     00 00   000 00
      00 00 0000 00  00
      000 00000  00   00
       00 00000 000 0 0
       000 0000 0000 00
 0000000000 00 0000 00
  0  0000000 0 000 000
   00 0000000 000 000
    00 000000 00 0000
     00 0000 00 000000000
00000000 000 0 00000  0
  0000000 00  0000  00
   0000000   00   000
    0000000 00 00000
     000000   0000
       000 00000000
         000000000000
      0      000000
       0
        0
Problem 91
      0      0           
      00    00           
      00000000           
     00  00  00          
    00        00         
    00  0  0  00         
    00        00         
   00000    00000        
   00 00000000 00        
   0000 0000 0000        
   000000 0000000        
    0 0000 000 0         
    0000 00 0000         
     0 000000 0          
     00 00 0 00          
      00000000           
00    00000000           
 000    0000             
   000000000000          
     000000000000        
    000       00000000000
   000           00000000
  000                0000
000                    00
00                       

See also

See also My JaCoP page for some JaCoP models and links.

May 09, 2009

Learning Constraint Programming IV: Logical constraints: Who killed Agatha? revisited

Here is a comparison of how different constraint programming systems implements the logical constraints in the problem Who Killed Agatha?, also known as The Dreadsbury Mansion Murder Mystery.. In Learning constraint programming - part II: Modeling with the Element constraint the problem was presented and showed how the systems implements the Element constraints.

Problem formulation from The h1 Tool Suite
Someone in Dreadsbury Mansion killed Aunt Agatha. Agatha, the butler, and Charles live in Dreadsbury Mansion, and are the only ones to live there. A killer always hates, and is no richer than his victim. Charles hates noone that Agatha hates. Agatha hates everybody except the butler. The butler hates everyone not richer than Aunt Agatha. The butler hates everyone whom Agatha hates. Noone hates everyone. Who killed Agatha?
Originally from F. J. Pelletier: Seventy-five problems for testing automatic theorem provers., Journal of Automated Reasoning, 2: 191-216, 1986.

Here we see compare how different constraint programming systems implement the three emphasized conditions in the problem formulation above:
  • the concept of richer
  • Charles hates noone that Agatha hates
  • No one hates everyone
All models defines the concepts of hates and richer as two matrices. The matrix declarations are omitted in the code snippets below.

Models

Here are the different models for the Who killed Agatha? problem. JaCoP and Choco has two version for how to implement the Element constraint, see the link above. Also, there is no Gecode/R version of this problem.

Defining the concept of richer

First we define the concept of richer, which consists of two parts:
  • No one is richer than him-/herself
  • if i is richer than j then j is not richer than i
This is an antisymmetric relation which is explored somewhat more (with an alternative predicate of the concept) in the MiniZinc model antisymmetric.mzn.

The logical concept used here is equivalence (if and only of).

MiniZinc

% define the concept of richer: no one is richer than him-/herself
forall(i in r) (
   richer[i,i] = 0
)

/\ % if i is richer than j then j is not richer than 
forall(i, j in r where i != j) (
   richer[i,j] = 1 <-> richer[j,i] = 0
)

Comet

Note that Comet don't have support for equivalence directly, instead we have to use two implications. Update: equivalence in Comet is written as == (I tested <=> which didn't work). Thanks to Pascal Van Hentenryck for pointing this out.
// no one is richer than him-/herself
forall(i in r)
  m.post(richer[i,i] == 0);

// if i is richer than j then j is not richer than i
forall(i in r, j in r : i != j) {
  /* earlier version: two implications
  m.post(richer[i,j] == 1 => richer[j,i] == 0);
  m.post(richer[j,i] == 0 => richer[i,j] == 1);
  */
  // equivalence
  m.post((richer[i,j] == 1) == (richer[j,i] == 0));
}

JaCoP

//  no one is richer than him-/herself
for(int i = 0; i < n; i++) {
    store.impose(new XeqC(richer[i][i], 0));
}

//  if i is richer than j then j is not richer than i
for(int i = 0; i < n; i++) {
    for(int j = 0; j < n; j++) {
        if (i != j) {
             store.impose(
                          new Eq(
                                 new XeqC(richer[i][j], 1),
                                 new XeqC(richer[j][i], 0)
                                 )
                          );

        }
    }
}

Choco

//   a) no one is richer than him-/herself
for(int i = 0; i < n; i++) {
    m.addConstraint(eq(richer[i][i], 0));
}

//   b) if i is richer than j then j is not richer than i
for(int i = 0; i < n; i++) {
    for(int j = 0; j < n; j++) {
        if (i != j) {
             m.addConstraint(
                             ifOnlyIf(
                                    eq(richer[i][j], 1),
                                    eq(richer[j][i], 0)
                                      )
                             );
        }
    }
}

Gecode

// no one is richer than him-/herself
for(int i = 0; i < n; i++) {
  rel(*this, richer_m(i,i), IRT_EQ, 0, opt.icl());
}

// if i is richer than j then j is not richer than i
for(int i = 0; i < n; i++) {
  for(int j = 0; j < n; j++) {
    if (i != j) {
      post(*this, tt(
              eqv(
                 richer_m(j,i) == 1, // <=>
                 richer_m(i,j) == 0)), 
                 opt.icl());
    }
  }
}

Charles hates noone that Agatha hates

Here is the definitions of the condition Charles hates noone that Agatha hates, which simply mean the implication:
  if Agatha hates X then Charles don't hate X

MiniZinc

When starting to modeling these kind of problems, I tend to follow the order of the conditions, which here means that the Charles part is before the Agatha part. When remodeling in another system the order tends to be fixed (cf the Comet version).
forall(i in r) (
   hates[charles, i] = 0 <- hates[agatha, i] = 1
)

Comet

forall(i in r)
  m.post(hates[agatha, i] == 1 => hates[charles, i] == 0);

JaCoP

for(int i = 0; i < n; i++) {
    store.impose(
                 new IfThen(
                            new XeqC(hates[agatha][i], 1),
                            new XeqC(hates[charles][i], 0)
                            )
                 );
}

Choco

I tend to copy/paste the models for Choco and JaCoP and just change the functions that are different. A consequence of this is that some special feature in one of these two systems is not used.
for(int i = 0; i < n; i++) {
    m.addConstraint(
                    implies(
                            eq(hates[agatha][i], 1),
                            eq(hates[charles][i], 0)
                            )
                    );
}

Gecode

for(int i = 0; i < n; i++) {
   post(*this, tt(
                  imp(
                        hates_m(i,agatha) == 1, 
                        // =>
                        hates_m(i,charles) == 0)), 
                        opt.icl());
}

No one hates everyone

This is the last condition to compare: No one hates everyone. It is implemented by a sum of the number of persons that each person hates, and this sum must be 2 or less. Note that it is possible to hate oneself.

MiniZinc

forall(i in r) (
  sum(j in r) (hates[i,j]) <= 2  
)

Comet

  forall(i in r) 
    m.post(sum(j in r) (hates[i,j]) <= 2);

JaCoP

Note: We could save the XlteqC constraint by restrict the domain of a_sum to 0..2 instead of 0..n (=3) but this explicit use of XlteqC is more clear.
for(int i = 0; i < n; i++) {
    FDV a[] = new FDV[n];
    for (int j = 0; j < n; j++) {
        a[j] = new FDV(store, "a"+j, 0, 1);
        a[j] = hates[i][j];
    }
    FDV a_sum = new FDV(store, "a_sum"+i, 0, n);
    store.impose(new Sum(a, a_sum));
    store.impose(new XlteqC(a_sum, 2));
}

Choco

Note: sum is an operator, which makes the condition somewhat easier to state than in JaCoP.
for(int i = 0; i < n; i++) {
    IntegerVariable a[] = makeIntVarArray("a", n, 0, 1);
    for (int j = 0; j < n; j++) {
        a[j] = hates[i][j];
    }
    m.addConstraint(leq(sum(a), 2));
}

Gecode

In Gecode this condition is quite easy to state by using linear. In order to use this there is a Matrix "view" of the hates matrix hates_m.
Matrix hates_m(hates, n, n);
// ...
for(int i = 0; i < n; i++) {
  linear(*this, hates_m.row(i), IRT_LQ, 2, opt.icl());
}

End comment

The mandatory end comment: There are probably better ways of modeling the problem than shown above, either by changing some details or by model the problem completely different. Maybe this will be done sometime...

May 08, 2009

Learning Constraint Programming III: decomposition of a global constraint: alldifferent_except_0

The series Learning Constraint Programming is a preparation for My talk at SweConsNet Workshop 2009: "Learning Constraint Programming (MiniZinc, JaCoP, Choco, Gecode/R, Comet, Gecode): Some Lessons Learned". Confusingly, the entries is not numbered in any logical order. Sorry about that.

Here are the previous entries:

Introduction

The global constraint alldifferent_except_0 (or the more general variant alldifferent_except_c) is one of my favorite global constraints. It is very handy to use when 0 (or any constant c) is coded as an unknown or irrelevant value. Then we can constraint the rest to be all distinct.

The great Global Constraint Catalog entry alldifferent_except_0) explains this constraint as:
Purpose
Enforce all variables of the collection VARIABLES to take distinct values, except those variables that are assigned to 0.

Example
​(<5,​0,​1,​9,​0,​3>)​
The alldifferent_except_0 constraint holds since all the values (that are different from 0) 5, 1, 9 and 3 are distinct.

Models

I have modeled a decomposition of alldifferent_except_0 in the following models, where the constraint is just tested perhaps combined with some other constraint, e.g. sorted or that there must be at least some zeros:

- MiniZinc alldifferent_except_0.mzn
- Comet: alldifferent_except_0.co
- Gecode/R: all_different_except_0.rb
- Choco: AllDifferentExcept0_test.java
- JaCoP: AllDifferentExcept0_test.java
- Gecode: alldifferent_except_0.cpp

Some models using alldifferent_except_0

And here is some real use of the constraint:

- Nonogram (Comet): nonogram.co. (A faster model using the regular constraint, nonogram_regular.co, is described here and here).
- I wrote about alldifferent_except_0 in Pi Day Sudoku 2009. However, as faster way of solving the problem was found and is described in Solving Pi Day Sudoku 2009 with the global cardinality constraint). Note: the competition is still on, so there is no link to any model.
- Sudoku generate (Comet): sudoku_generate.co
- all paths graph (MiniZinc): all_paths_graph.mzn
- Cube sum (MiniZinc): cube_sum.mzn
- Message sending (MiniZinc): message_sending.mzn

As the first two entries indicates, there may be faster solutions than using (a decomposition) of alldifferent_except_0, but even as a decomposition is a great conceptual tool when modeling a problem.

Implementations

In the implementations below we also see how to define a function (predicate) in the constraint programming systems.

For the Gecode/R model there are different approaches:
- "standard" ("direct") approach where we loop over all different pairs of elements and ensures that if both values are different from 0 then they should be different
- using count
- using global cardinality ("simulated" in Gecode/R, see below)

Also, in some models we use the slighly more general version alldifferent_except_c where c is any constant (e.g. "Pi" in the Pi Day Sudoku puzzle mentioned above.

MiniZinc:

Model: alldifferent_except_0.mzn.
predicate all_different_except_0(array[int] of var int: x) =
   let {
      int: n = length(x)
   }
   in
   forall(i,j in 1..n where i != j) (
        (x[i] > 0 /\ x[j] > 0) 
        -> 
        x[i] != x[j]
   )
;

// usage:
constraint all_different_except_0(x);

Comet:

Model: alldifferent_except_0.co.
function void alldifferent_except_0(Solver m, var{int}[] x) {
  int n = x.getSize();
  forall(i in 1..n, j in i+1..n) {
    m.post(
           x[i] > 0 && x[j] > 0 
           => 
           x[i] != x[j]
           );
  }
}

// usage
exploreall {
  // ...
  alldifferent_except_0(m, x);
}

Gecode/R

Model:all_different_except_0.rb.

When modeling the constraint in Gecode/R, I experimented with different approaches. The reification variant all_different_except_0_reif is actually quite fast.
# The simplest and the fastest implementation 
# using count for 1..max (poor man's global cardinality)
def all_different_except_0
  (1..self[0].domain.max).each{|i|
    self.count(i).must <= 1
  }
end

# global cardinality version using an extra array with the counts
def global_cardinality(xgcc)
  (self[0].domain.max+1).times{|i|
    xgcc[i].must == self.count(i)
  }
end

#
# The standard approach using reification.
#
def all_different_except_0_reif(x)
  n = x.length
  b1_is_an bool_var_matrix(n,n)
  b2_is_an bool_var_matrix(n,n)
  b3_is_an bool_var_matrix(n,n)
  n.times{|i|
    n.times{|j|
      if i != j then 
        x[i].must_not.equal(0, :reify => b1[i,j]) 
        x[i].must_not.equal(0, :reify => b2[i,j]) 
        x[i].must_not.equal(x[j], :reify => b3[i,j])
        (b1[i,j] & b2[i,j]).must.imply(b3[i,j])
      else 
        b1[i,j].must.true
        b2[i,j].must.true
        b3[i,j].must.true
      end
    }
  }
end

    # ...
    # usage:
    x.all_different_except_0
    # all_different_except_0_gcc(x)
    # all_different_except_0_reif(x)

Choco

Model: AllDifferentExcept0_test.java.

Note that here alldifferent_except_0 is derived from the more general version alldifferent_except_c.
//
// decomposition of alldifferent except 0
//
public void allDifferentExcept0(CPModel m, IntegerVariable[] v) {
    allDifferentExceptC(m, v, 0);
}

//
// slightly more general: alldifferent except c
//
public void allDifferentExceptC(CPModel m, IntegerVariable[] v, int c) {
    int len = v.length;
    for(int i = 0; i < v.length; i++) {
        for(int j = i+1; j < v.length; j++) {
            m.addConstraint(ifThenElse(
                                       and(
                                           gt(v[i], c), 
                                           gt(v[j], c)
                                           ), 
                                       neq(v[i], v[j]),
                                       TRUE
                                       )
                            );
        }
    }    
}

// ...
// usage: 

    m.addConstraint(allDifferent(queens));

JaCoP

Model: AllDifferentExcept0_test.java

This is exactly the same approach as the Choco version.
//
// decomposition of alldifferent except 0
//
public void allDifferentExcept0(FDstore m, FDV[] v) {
    allDifferentExceptC(m, v, 0);
} // end allDifferentExcept0


//
// slightly more general: alldifferent except c
//
public void allDifferentExceptC(FDstore m, FDV[] v, int c) {
    int len = v.length;
    for(int i = 0; i < v.length; i++) {
        for(int j = i+1; j < v.length; j++) {
            m.impose(new IfThen(
                                       new And(
                                           new XneqC(v[i], c), 
                                           new XneqC(v[j], c)
                                           ), 
                                       new XneqY(v[i], v[j])
                                       )
                            );
        }
    }        
} // end allDifferentExceptC


        // ...
        // usage:
        allDifferentExcept0(store, x);

Gecode

alldifferent_except_0.cpp

The Gecode version is very succint since it use overloaded boolean operators. Very nice.
void alldifferent_except_0(Space& space, IntVarArray x, IntConLevel icl = ICL_BND) {
  for(int i = 0; i < x.size(); i++) {
    for(int j = i+1; j < x.size(); j++) {
      post(space,
           tt(
           imp(x[i] != 0 && x[j] != 0, 
           // =>
           x[i] != x[j])),
           icl
           );
    }
  }
} // alldifferent_except_0


// ...
// usage:
    alldifferent_except_0(*this, x, opt.icl());

May 05, 2009

Learning constraint programming - part II: Modeling with the Element constraint

As indicated last in Learning constraint programming (languages) - part I here are some findings when implementing Crossword, Word square, and Who killed Agatha?. See links below for the implementations.

Spoiled
The first constraint programming system i learned after constraint logic programming in Prolog was MiniZinc. When implemented the problems below I realized that I have been quite spoiled by using MiniZinc. The way MiniZinc (and also Comet) supports the Element constraint, i.e. access of variable arrays/matrices, is very straightforward in these systems and it doesn't matter whether the array to access is an array of integers or of non-variable variable integers. In the Java (Choco, JaCoP) and C++ systems (Gecode), however, this is another matter. Due to different circumstances I have not implemented these models in Gecode/R.

Element in MiniZinc and Comet
Accessing arrays and matrices in MiniZinc and Comet is simply done by using the [] construct, no matter what the type of the array or the index are (I assume integers and variable integers here). For the other systems we must explicitly use the Element constraint (called nth in Choco).

Crossword

This is a standard constraint programming problem, used as a running example in Apt's great book Principles of Constraint Programming. Here is a formulation of the problem (stated differently than in the book):
Place the words listed below in the following crossword. The '#' means a blocked cell, and the numbers indicate the overlappings of the words.

      1   2   3   4   5
    +---+---+---+---+---+
  1 | 1 |   | 2 |   | 3 |
    +---+---+---+---+---+
  2 | # | # |   | # |   |
    +---+---+---+---+---+
  3 | # | 4 |   | 5 |   |
    +---+---+---+---+---+
  4 | 6 | # | 7 |   |   |
    +---+---+---+---+---+
  5 | 8 |   |   |   |   |
    +---+---+---+---+---+
  6 |   | # | # |   | # |
    +---+---+---+---+---+
                         
We can use the following words
* AFT * ALE * EEL * HEEL * HIKE * HOSES * KEEL * KNOT * LASER * LEE * LINE * SAILS * SHEET * STEER * TIE

Models


MiniZinc: crossword.mzn
Choco: CrossWord.java
JaCoP: CrossWord.java
Gecode/R: (Not implemented)
Comet: crossword.co
Gecode: crossword.cpp

Note: I have seen more general models for solving crossword problems in Choco, JaCoP, Gecode/R, and Gecode with constructions other that the simple Elements used here. Since I wanted to compare the same way of solving the problem using the same Element-construct this may be an unfair comparison between the systems. Well, this is at least a finding how to implement this problem by Element...

Explanation of variables
The matrix A is the individual chars of the words (Comet variant):
int A[1..num_words,1..word_len] = 
  [
   [h, o, s, e, s], //  HOSES
   [l, a, s, e, r], //  LASER
   [s, a, i, l, s], //  SAILS
   [s, h, e, e, t], //  SHEET
   [s, t, e, e, r], //  STEER
   [h, e, e, l, 0], //  HEEL
   [h, i, k, e, 0], //  HIKE
   [k, e, e, l, 0], //  KEEL
   [k, n, o, t, 0], //  KNOT
   [l, i, n, e, 0], //  LINE
   [a, f, t, 0, 0], //  AFT
   [a, l, e, 0, 0], //  ALE
   [e, e, l, 0, 0], //  EEL
   [l, e, e, 0, 0], //  LEE
   [t, i, e, 0, 0]  //  TIE
];
overlapping is the matrix of the overlapping cells)
This is the Comet version:
[
 [1, 3, 2, 1],   //  s
 [1, 5, 3, 1],   //  s 
  
 [4, 2, 2, 3],   //  i
 [4, 3, 5, 1],   //  k
 [4, 4, 3, 3],   //  e
  
 [7, 1, 2, 4],   //  l
 [7, 2, 5, 2],   //  e
 [7, 3, 3, 4],   //  e
  
 [8, 1, 6, 2],   //  l
 [8, 3, 2, 5],   //  s
 [8, 4, 5, 3],   //  e
 [8, 5, 3, 5]    //  r
 ];
E is the variable array of which the word to use for the different overlappings. This is in fact the only variable (array) that is needed in the problem, apart from the utility/convenience variables.

The main constraint for the crossword example in each system is stated thus:

MiniZinc:
forall(i in 1..num_overlapping) (
   A[E[overlapping[i,1]], overlapping[i,2]] =  A[E[overlapping[i,3]], overlapping[i,4]]
)
Comet:
  forall(i in 1..num_overlapping) {
    cp.post(A[E[overlapping[i,1]], overlapping[i,2]] ==  A[E[overlapping[i,3]], overlapping[i,4]], onDomains);
  }
Choco:
Note that Choco has a special Element which support two dimensional arrays (matrix), which we use.
for(int I = 0; I < num_overlapping; I++) {
  IntegerVariable tmp = makeIntVar("tmp" + I, 1, 26);
  M.addConstraint(nth(E[overlapping[I][0]], W[overlapping[I][1]], AX, tmp));
  M.addConstraint(nth(E[overlapping[I][2]], W[overlapping[I][3]], AX, tmp));
}
JaCoP:
Here we had used some trickery by using a transposed version of the word matrix since JaCoP has no special Element constraint for two dimensional arrays.
for (int I = 0; I < num_overlapping; I++) {
   FDV tmp = new FDV(store, "TMP" + I, 0, num_words*word_len);
   store.impose(new Element(E[overlapping[I][0]], words_t[overlapping[I][1]], tmp));
   store.impose(new Element(E[overlapping[I][2]], words_t[overlapping[I][3]], tmp));
}
Gecode:
This is more complicated compared to the two Java systems since in Gecode we use an array (of length rows*cols) to simulate the matrix. (There is a Matrix "view" in Gecode but the indices must be of type integer, not IntVar, so it can not be used.) Also, the constraints plus and mult takes IntVar as argument.

The first overlapped crossing is "expanded" like this (Gecode is 0-based):
   A[E[overlapping[i,0]], overlapping[i,1]] // MiniZinc/Comet
   =>
   a1 = A[ E[I*4+0] * word_len + overlapping[I*4+1]] // Gecode
Here is the complete code. The comments hopefully explains what is going on.

First we define an utility function for accessing the element according to the above principle.
/**
 *
 * Special version of element for an array version of a "matrix" words,
 * E is an integer variable array, C is an array of IntVars for
 * the offset j in the words "matrix".
 *
 * The call 
 *    element_offset(*this, words, E[i], word_len_v, C[j], res, opt.icl());
 *
 * corresponds to:
 *    res = words[E[i], j] --> words[E[i]*word_len+J]
 *
 */
void element_offset(Space& space,
                   IntArgs words,
                   IntVar e,
                   IntVar word_len,
                   IntVar c,
                   IntVar res,
                   IntConLevel icl = ICL_DOM) {

      element(space, words, 
              plus(space, 
                   mult(space, 
                        e, 
                        word_len, icl), 
                   c, icl), 
              res, icl);
}

The function is then used as follows:
for(int I = 0; I < num_overlapping; I++) {
   IntVar e1(*this, 0, num_overlapping*4);
   IntVar e2(*this, 0, num_overlapping*4);

   IntVarArray o(*this, 4, 0, num_overlapping*4);
   for(int J = 0; J < 4; J++) {
     post(*this, o[J] == overlapping[I*4+J], opt.icl());
   }

   element(*this, E, o[0], e1, opt.icl());      // e1 = E[I*4+0]
   element(*this, E, o[2], e2, opt.icl());      // e2 = E[I*4+2]

   IntVar a1(*this, 0, num_words*word_len);
      
   element_offset(*this, A, e1, word_len_v, o[1], a1, opt.icl());
   element_offset(*this, A, e2, word_len_v, o[3], a1, opt.icl());
}
(The same element_offset function is also used in the word_square problem below.) It took quite a time to get the function and temporary variables (and their domains) right. With training (and the element_offset as a skeleton) similiar problems should be easier to implement.

Note: this is not a bashing of Gecode. Gecode is a great system and it happens that for this specific problem, Gecode does not has the appropriate support. I should also mention that it was a long time since I programmed in C++ and am little rusty.

As mentioned earlier, I have been very spoiled by the MiniZinc (and Comet) constructs. Also: I'm a very 'lazy' person (in the Perl sense of the word lazy) and likes the agile programming languages - Perl, Ruby, Python etc - much for their high level constructs.

Word square problem

The word problem is a cousin to the crossword problem, and is described in Wikipedia's Word_square:
A word square is a special case of acrostic. It consists of a set of words, all having the same number of letters as the total number of words (the "order" of the square); when the words are written out in a square grid horizontally, the same set of words can be read vertically.
An example of order 7 found by the Comet model where we see that the first row word (aalborg) is also the first column word.

aalborg
abaiser
lantaca
bittman
osamine
recanes
granese

Models

Here are the models for solving the Word square problem:
MiniZinc: word_square.mzn
Choco: WordSquare.java
JaCoP: WordSquare.java
Gecode/R: (Not implemented it Gecode/R)
Comet: word_square.co
Gecode: word_square.cpp

It is somewhat easier than the crossword problem. As before, E is the array of the index of the words to use, and words is an matrix of the words. Also, these models is an experiment of how to read a file, the word list /usr/dict/words (standard on Unix/Linux systems).

MiniZinc
forall(I, J in 1..word_len) (
  words[E[I], J] = words[E[J],I]
)
Comet
  forall(i in 1..word_len) {
    forall(j in 1..word_len) {    
      cp.post(words[E[i], j] == words[E[j],i], onDomains);
    }
  }
JaCoP
// 
// The overlappings (crossings).
// Note that we use a transposed word matrix for the Element.
//
for(int i = 0; i < word_length ; i++) {
    for(int j = 0; j < word_length ; j++) {
        // Comet: words[E[i], j] ==  words[E[j],i]
        FDV tmp = new FDV(store, "tmp" + i + " " + j, 0, dict_size);
        store.impose(new Element(E[i], words[j], tmp));
        store.impose(new Element(E[j], words[i], tmp));
    }
}
Choco
// Constants for the nth constraint below
IntegerVariable [] C = new IntegerVariable[dict_size];
for (int I = 0; I < word_length; I++) {
    C[I] = makeIntVar("C"+I, I,I);
}

// The overlappings (crossings)
for(int I = 0; I < word_length ; I++) {
    for(int J = 0; J < word_length ; J++) {
        // Comet: words[E[i], j] ==  words[E[j],i]
        IntegerVariable tmp = makeIntVar("tmp" + I + " " + J, 0, dict_size);
        M.addConstraint(nth(E[I], C[J], words, tmp));
        M.addConstraint(nth(E[J], C[I], words, tmp));
    }
}
Gecode
Note that this model use the same function element_offset that was used in the Crossword problem. It took some time to realize that it also could be used here.
// convenience variables for the element constraints below
// since element, plus, and mult wants IntVars.
IntVar word_len_v(*this, word_len, word_len);
IntVarArray C(*this, word_len, 0, word_len-1);
for(int i = 0; i < word_len; i++) {
  rel(*this, C[i], IRT_EQ, i, opt.icl());
}

for(int i = 0; i < word_len; i++) {
  for(int j = 0; j < word_len; j++) {
    // words[E[i], j] ==  words[E[j],i]

    IntVar tmp(*this, 0, num_words);

    // tmp == words[E[i], j] --> words[E[i]*word_len+j]
    element_offset(*this, words, E[i], word_len_v, C[j], tmp, opt.icl());

    // tmp == words[E[j], i]  --> words[E[j]*word_len+i]
    element_offset(*this, words, E[j], word_len_v, C[i], tmp, opt.icl());
  }
}

Who killed Agatha?

This is a standard benchmark for theorem proving, also known as The Dreadsbury Mansion Murder Mystery.

Problem formulation from The h1 Tool Suite
Someone in Dreadsbury Mansion killed Aunt Agatha. Agatha, the butler, and Charles live in Dreadsbury Mansion, and are the only ones to live there. A killer always hates, and is no richer than his victim. Charles hates noone that Agatha hates. Agatha hates everybody except the butler. The butler hates everyone not richer than Aunt Agatha. The butler hates everyone whom Agatha hates. Noone hates everyone. Who killed Agatha?
Originally from F. J. Pelletier: Seventy-five problems for testing automatic theorem provers., Journal of Automated Reasoning, 2: 191-216, 1986.

Models


MiniZinc: who_killed_agatha.mzn
JaCoP : WhoKilledAgatha.java
JaCoP : WhoKilledAgatha_element.java
Choco: WhoKilledAgatha.java
Choco: WhoKilledAgatha_element.java
Comet: who_killed_agatha.co
Gecode: who_killed_agatha.cpp

In Some new Gecode models I wrote about the findings in implemented this problem in Gecode and compared to Comet/MiniZinc.

The models use two 3x3 matrices for representing the two relations hates and richer with 0..1 as domain (i.e. boolean). The Element constraints is used for implementing the condition A killer always hates, and is no richer than his victim. where the_killer is an integer variable; the_victim is in some models replaced by agatha (the integer 0). The interesting thing here is that at least one of the indices are integer variables, which caused the difficulties in the two problems above.

These models also use a lot of boolean constructs. A comparison of how these are implemented in the different CP systems may be described in a future blog post.

MiniZinc
hates[the_killer, the_victim] = 1 /\
richer[the_killer, the_victim] = 0
Comet:
m.post(hates[the_killer, the_victim] == 1);
m.post(richer[the_killer, the_victim] == 0);
Note: In the models below I have simplified the problem by using agatha (defined as the integer 0) instead of the integer variable the_victim. This is not a problem since we know that Agatha is the victim, and is the reason why the Element is easier to use than for Crossword and Word square.

JaCoP variant 1 (no Element):
JaCoP don't have a direct support for the case when the index i (in matrix[i][j]) is an integer variable so the first variant of modeling the condition A killer always hates, and is no richer than his victim. does not use Element at all. In WhoKilledAgatha.java we simply loop over all integers (0..2), check if "this" i equals the_killer and then we can use two integers for accessing the matrices. Also, note the IfThen construct.
for(int i = 0; i < n; i++) {
    store.impose(
                 new IfThen(
                            new XeqC(the_killer, i),
                            new XeqC(hates[i][agatha], 1)
                            )
                 );
    store.impose(
                 new IfThen(
                            new XeqC(the_killer, i),
                            new XeqC(richer[i][agatha], 0)
                            )
                 );
}
This was the first variant I implemented, but then I recall the "trickery" used in Crossword and Word square where the matrices where transposed and Element could be used. The problem with this approach is that all constraints must be rewritten in a way that may be confusing. Come to think of it, maybe the names of the matrices should have been changed to is_hated_by and poorer.

JaCoP variant 2 (transposed matrices, Element)
This method of transposition and using Element is implemented in WhoKilledAgatha_element.java. The constraint is now much simpler:
int shift = -1;
for(int i = 0; i < n; i++) {
    store.impose(new Element(the_killer, hates[agatha], one, shift));
    store.impose(new Element(the_killer, richer[agatha], zero, shift));
}
Note: the Element in JaCoP defaults to start index as 1, but has support for shifting it to 0, by using -1 as the shift parameter. Choco variant 1 (no Element)
I implemented exact the same principle that was used in the two JaCoP model in the two Choco models. The first - no Element - is WhoKilledAgatha.java:

for(int i = 0; i < n; i++) {
    m.addConstraint(implies(
                            eq(the_killer,i),
                            eq(hates[i][agatha], 1)
                            )
                    );
    m.addConstraint(implies(
                            eq(the_killer, i),
                            eq(richer[i][agatha], 0)
                            )
                    );
}
Choco variant 2 (transposed matrices, nth)
Note: one and zero are integer variables since nth cannot handle plain integers as the last argument.
for(int i = 0; i < n; i++) {
   m.addConstraint(nth(the_killer, hates[agatha], one));
   m.addConstraint(nth(the_killer, richer[agatha], zero)
}

Comments

Here we have seen - not surprisingly - that using the Element constraint is quite different depending of which CP system we use and it can be easy or not so easy. It was my explicit intension to see how to solve the same problem as similiar as possible. We should also note that most (if not all) problems can be modeled in many ways, some not using Element at all.

One last comment: The two Java models of Who killed Agatha? took quite a long time to implement. The main reason for that was not the the handling of Element but was due to a bug of confusing the two matrices in one of the conditions. Sigh.

February 21, 2009

JaCoP: Features request for version 2.4

In Submit features requests for JaCoP 2.4 Radoslaw Szymanek asks for feature requests for JaCoP version 2.4:


Dear all,

We are already planning release 2.4. We will include global constraint - bounded Knapsack constraint. Given time we will also look into creating a parser from Flatzinc to jaCoP.

Please let us know of other features you would like to see in future releases of JaCoP. Requests for additional global constraints and other requests which extend modeling power of JaCoP will have a preference. Please contact Radoslw Szymanek to submit your request ( radoslaw { . } szymanek [@] gmail ( dot ) com ).

best regards,

Radoslaw Szymanek

I especially look forward to the FlatZinc parser.


By the way, the News section has a feed (which I have missed): RSS 2.0, Atom.

January 15, 2009

Some models in Gecode/R (Ruby interface to Gecode)

Gecode/R is a great Ruby interface to Gecode (implemented in C++). At last I now have done some modeling in Gecode/R and it's nice.

The models and some other information about Gecode/R can be seen at My Gecode/R page.

Since Ruby is a very high level language, the modelling can be done "high levelish", more so than for example with the Java solvers Choco and JaCoP. And that is a feature I really like.

I'm still learning Gecode/R and have surely missed some stuff. There are not many examples in the package (these are also commented at Examples). Some things have been of great help
* The Sitemap
* The Documentation, especially the RDocs
* The test files in the specs directory.


An example: Survo Puzzle
From Wikipedia's Survo Puzzle


In a Survo puzzle the task is to fill an m * n table by integers 1,2,...,m*n so
that each of these numbers appears only once and their row and column sums are
equal to integers given on the bottom and the right side of the table.
Often some of the integers are given readily in the table in order to
guarantee uniqueness of the solution and/or for making the task easier.

E.g. the puzzle 128/2008 is presented with the following clues:


* * * * * * 30
* * 18 * * * 86
* * * * * * 55
22 11 42 32 27 37

where * marks a digit to find. The number to the right is the row sums which the row must sum to, and the last row is the column sums.

The unique solution of the problem is (by running the program ruby survo_puzzle.rb survo_puzzle_128_2008.txt)


Solution #1
4 1 10 5 3 7 = 30
12 8 18 16 15 17 = 86
6 2 14 11 9 13 = 55
= = = = = =
22 11 42 32 27 37

Statistics:
propagations: 3784
failures: 174
clones: 179
memory: 25740
commits: 461
Number of solutions: 1

The relevant constraint programming code is below (slightly edited). I think it's quite nice.


....
def initialize(clues, rowsums, colsums)
r = rowsums.length # number of rows
c = colsums.length # number of columns
x_is_an int_var_matrix(r, c, 1..r*c) # the solution matrix
x.must_be.distinct # all values in x must be distinct
# values in clues with values > 0 are copied straight off to x
r.times{|i|
c.times{|j|
x[i,j].must == clues[i][j] if clues[i][j] > 0
}
}
r.times{|i| x[i,0..c-1].sum.must == rowsums[i] } # check row sums
c.times{|j| x.transpose[j,0..r-1].sum.must == colsums[j] } # check column sums
branch_on x, :variable => :smallest_size, :value => :min
end

The full program is here: survo_puzzle.rb, and three data files:
survo_puzzle_126_2008.txt
survo_puzzle_127_2008.txt
survo_puzzle_128_2008.txt


The Gecode/R models
Below are the models shown at My Gecode/R page. The selection is mostly for comparison with models implemented in other constraint programming languages, see respectively:
My MiniZinc page (which has a lot more models)
My JaCoP page
My Choco page

The models first written (e.g. diet, least_diff) are very "un-Rubyish" since I wanted to model straight after the MiniZinc models. Maybe I Rubyfi these later.

After a while the modeling went quite easy, and both de Bruijn and Minesweeper was done surprisingly fast. I do, however, miss MiniZinc's sum construct, since it which make some things easier (e.g. the neighbour summing in minesweeper.rb).

The execution time of the models os approximately the same time as the corresponding MiniZinc model with the Gecode/flatzinc solver which is normally quite fast. The big exception of these examples is coins_grid which seems to be slow for the constraint programming systems, but fast with linear programming systems, e.g. with the MiniZinc solvers ECLiPSe/ic and MiniZinc/mip.

References to the problems etc are in the header of the model.

January 05, 2009

Tom Schrijvers: "Monadic Constraint Programming" (in Haskell)

Tom Schrijvers presents in the blog post Monadic Constraint Programming a draft version of the paper Monadic Constraint Programming written by him, Peter Stuckey, and Philip Wadler:


A constraint programming system combines two essential components: a constraint solver and a search engine. The constraint solver reasons about satisfiability of conjunctions of constraints, and the search engine controls the search for solutions by iteratively exploring a disjunctive search tree defined by the constraint program. In this paper we give a monadic definition of constraint programming where the solver is defined as a monad threaded through the monadic search tree. We are then able to define search and search strategies as first class objects that can themselves be built or extended by composable search transformers. Search transformers give a powerful and unifying approach to viewing search in constraint programming, and the resulting constraint programming system is first class and extremely flexible.


Prototype code in Haskell can be downloaded here.