golomb.c 4.08 KB
#include <stdio.h>
#include <string.h>
#include <fcntl.h>
#include <errno.h>

#include "fdc_int.h"

// Golomb ruler (CSPLib 006)

#define MMAX 15				// maximum number of marks + 1
#define DMAX (MMAX * (MMAX - 1) / 2)	// number of differences

// upper bound of the variables' domain (XXX: why M^2?)
#define DOMMAX \
  (fixed ? minimum_length[M] : (M * M > MAX_VALUE ? MAX_VALUE : M * M))

// from http://www.research.ibm.com/people/s/shearer/gropt.html
static int minimum_length[] =
  {0,  0,   1,   3,   6,  11,  17,  25,  34,  44,  55,	// 10
      72,  85, 106, 127, 151, 177, 199, 216, 246, 283,	// 20
     333, 356, 372, 425, 480, 492};


int M = 7, D; // actual number of marks and of differences

int main(int argc, char *argv[])
{
  fd_int ruler[MMAX];
  fd_int differences[DMAX];
  int arg;
  int solutions = 0, one_solution = 1;
  int fixed = 0;	     // bound domains by the known minimal ruler length
  int restrict = 0;	     // restrict ruler length
  int len_lb, len_ub;
  int use_label = 0;
  int i, j, d;
  int seed;
  int dev_random;

  fd_init(&argc, &argv);

  for (arg = 1; arg < argc; ++arg)
    if (!strcmp(argv[arg], "--all"))
      one_solution = 0;
    else if (!strcmp(argv[arg], "--fix"))
      fixed = 1;
    else if (!strcmp(argv[arg], "--label"))
      use_label = 1;
    else if (!strncmp(argv[arg], "--length=", sizeof("--length=") - 1))
      {
	char *p = argv[arg] + strlen("--length=");

	len_lb = len_ub = atoi(p);
	if (p = strchr(p, '-'))
	  len_ub = *(p + 1) ? atoi(p + 1) : -1;

	restrict = 1;
      }
    else
      {
	M = atoi(argv[arg]);

	break;
      }

  if (M + 1 > MMAX)
    {
      fprintf(stderr, "marks must be %d or less\n", MMAX - 1);
      return 1;
    }

#ifdef LOCAL_SEARCH
#if 0
  seed = time(0);
#else
  if ((dev_random = open("/dev/urandom", O_RDONLY)) == -1)
    fd__fatal("/dev/urandom: %s\n", strerror(errno));

  read(dev_random, &seed, sizeof(seed));
#endif

  srandom(seed);
  fd__trace("seed = %u\n", seed);
#endif

  if (!fixed)
    fd__trace("length upper bound is %d\n", DOMMAX);

  D = M * (M - 1) / 2;

  ruler[0] = fd_new(0, 0);

  i = 1;
  while (++arg < argc)
    {
      int lb, ub;
      char *p;

      lb = ub = atoi(argv[arg]);
      if (p = strchr(argv[arg], '-'))
	ub = atoi(p + 1);

      ruler[i] = fd_new(lb, ub);

      i++;
    }

  for (; i < M - 1; ++i)
    ruler[i] = fd_new(1, DOMMAX);

  if (restrict)
    {
      if (len_ub == -1)
	len_ub = DOMMAX;

      ruler[M - 1] = fd_new(len_lb, len_ub);
    }
  else
    ruler[M - 1] = fd_new(1, DOMMAX);

  if (use_label)
    fd_label(ruler, M);

  if (!fixed)
    fd_min(ruler[M - 1]);

  for (i = 1; i < M; ++i)
    fd_lt(ruler[i - 1], ruler[i]);

  j = 0;
  for (d = 1; d < M; ++d)
    {
      differences[j++] = ruler[d];

      for (i = d + 1; i < M; ++i)
	{
	  differences[j] = fd_new(1, DOMMAX);

	  fd_var_eq_minus(differences[j++], ruler[i], ruler[i - d]);
	}
    }

  fd_all_different(differences, D);

#if 01
  // from GECODE: ruler[j] - ruler[i] >= sum(1 .. j-i)
  j = 0;
  for (d = 1; d < M; ++d)
    {
      fd_int s = fd_const(d * (d + 1) / 2);

      for (i = d; i < M; ++i)
	fd_ge(differences[j++], s);
    }
#elif 0
  // from GECODE: ruler[j] - ruler[i] >= minimum_length[j - i + 1]
  j = 0;
  for (d = 1; d < M; ++d)
    {
      fd_int grl = fd_const(minimum_length[d + 1]);

      for (i = d; i < M; ++i)
	fd_ge(differences[j++], grl);
    }
#endif

  // discard symetric (wrt the differences) solutions:
  // m_1 - m_0 < m_{M-1} - m_{M-2}
  if (M > 2)
    fd_lt(differences[0], differences[M - 2]);

  while (fd_solve() == FD_OK)
    {
      printf("solution %d:\n", ++solutions);

      for (i = 0; i < M; ++i)
	{
	  fd_print(ruler[i]);
	  putchar(' ');
	}
      putchar('\n');

      j = 0;
      for (d = 1; d < M; ++d)
	{
	  printf("%*s", 2 * d, "");

	  for (i = d; i < M; ++i)
	    {
	      fd_print(differences[j++]);
	      putchar(' ');
	    }
	  putchar('\n');
	}

#if !(defined(LOCAL_SEARCH) || defined(DISTRIBUTED_SOLVER))
      if (one_solution)
#endif
	break;
    }

  if (solutions)
    printf("%d solutions found\n", solutions);
  else
    printf("inconsistent CSP\n");

  fd_end();

  return !solutions;
}