pip.c

   1 #include <assert.h>
   2 #include <string.h>
   3 #include "isl_set.h"
   4 #include "isl_tab.h"
   5 #include "isl_map_private.h"
   6 #include "isl_sample.h"
   7 #include "isl_scan.h"
   8 #include "isl_seq.h"
   9 #include "isl_ilp.h"
  10
  11 /* The input of this program is the same as that of the "example" program
  12  * from the PipLib distribution, except that the "big parameter column"
  13  * should always be -1.
  14  *
  15  * Context constraints in PolyLib format
  16  * -1
  17  * Problem constraints in PolyLib format
  18  * Optional list of options
  19  *
  20  * The options are
  21  *      Maximize        compute maximum instead of minimum
  22  *      Rational        compute rational optimum instead of integer optimum
  23  *      Urs_parms       don't assume parameters are non-negative
  24  *      Urs_unknowns    don't assume unknowns are non-negative
  25  */
  26
  27 struct pip_options {
  28         struct isl_options      *isl;
  29         unsigned                 verify;
  30 };
  31
  32 struct isl_arg pip_options_arg[] = {
  33 ISL_ARG_CHILD(struct pip_options, isl, "isl", isl_options_arg)
  34 ISL_ARG_BOOL(struct pip_options, verify, 'T', "verify", 0)
  35 };
  36
  37 ISL_ARG_DEF(pip_options, struct pip_options, pip_options_arg)
  38
  39 static struct isl_basic_set *to_parameter_domain(struct isl_basic_set *context)
  40 {
  41         struct isl_dim *param_dim;
  42         struct isl_basic_set *model;
  43
  44         param_dim = isl_dim_set_alloc(context->ctx,
  45                                         isl_basic_set_n_dim(context), 0);
  46         model = isl_basic_set_empty(param_dim);
  47         context = isl_basic_set_from_underlying_set(context, model);
  48
  49         return context;
  50 }
  51
  52 static isl_basic_set *construct_cube(isl_basic_set *context)
  53 {
  54         int i;
  55         unsigned dim;
  56         int range;
  57         isl_int m, M;
  58         struct isl_basic_set *cube;
  59         struct isl_basic_set *interval;
  60         struct isl_basic_set_list *list;
  61
  62         dim = isl_basic_set_total_dim(context);
  63         if (dim == 0) {
  64                 struct isl_dim *dims;
  65                 dims = isl_dim_set_alloc(context->ctx, 0, 0);
  66                 return isl_basic_set_universe(dims);
  67         }
  68
  69         if (dim >= 8)
  70                 range = 4;
  71         else if (dim >= 5)
  72                 range = 6;
  73         else
  74                 range = 30;
  75
  76         isl_int_init(m);
  77         isl_int_init(M);
  78
  79         isl_int_set_si(m, -range);
  80         isl_int_set_si(M, range);
  81
  82         interval = isl_basic_set_interval(context->ctx, m, M);
  83         list = isl_basic_set_list_alloc(context->ctx, dim);
  84         for (i = 0; i < dim; ++i)
  85                 list = isl_basic_set_list_add(list, isl_basic_set_copy(interval));
  86         isl_basic_set_free(interval);
  87         cube = isl_basic_set_product(list);
  88
  89         isl_int_clear(m);
  90         isl_int_clear(M);
  91
  92         return cube;
  93 }
  94
  95 isl_basic_set *plug_in_parameters(isl_basic_set *bset, struct isl_vec *params)
  96 {
  97         int i;
  98
  99         for (i = 0; i < params->size - 1; ++i)
 100                 bset = isl_basic_set_fix(bset,
 101                                          isl_dim_param, i, params->el[1 + i]);
 102
 103         bset = isl_basic_set_remove(bset, isl_dim_param, 0, params->size - 1);
 104
 105         isl_vec_free(params);
 106
 107         return bset;
 108 }
 109
 110 isl_set *set_plug_in_parameters(isl_set *set, struct isl_vec *params)
 111 {
 112         int i;
 113
 114         for (i = 0; i < params->size - 1; ++i)
 115                 set = isl_set_fix(set, isl_dim_param, i, params->el[1 + i]);
 116
 117         set = isl_set_remove(set, isl_dim_param, 0, params->size - 1);
 118
 119         isl_vec_free(params);
 120
 121         return set;
 122 }
 123
 124 /* Compute the lexicographically minimal (or maximal if max is set)
 125  * element of bset for the given values of the parameters, by
 126  * successively solving an ilp problem in each direction.
 127  */
 128 struct isl_vec *opt_at(struct isl_basic_set *bset,
 129         struct isl_vec *params, int max)
 130 {
 131         unsigned dim;
 132         struct isl_vec *opt;
 133         struct isl_vec *obj;
 134         int i;
 135
 136         dim = isl_basic_set_dim(bset, isl_dim_set);
 137
 138         bset = plug_in_parameters(bset, params);
 139
 140         if (isl_basic_set_fast_is_empty(bset)) {
 141                 opt = isl_vec_alloc(bset->ctx, 0);
 142                 isl_basic_set_free(bset);
 143                 return opt;
 144         }
 145
 146         opt = isl_vec_alloc(bset->ctx, 1 + dim);
 147         assert(opt);
 148
 149         obj = isl_vec_alloc(bset->ctx, 1 + dim);
 150         assert(obj);
 151
 152         isl_int_set_si(opt->el[0], 1);
 153         isl_int_set_si(obj->el[0], 0);
 154
 155         for (i = 0; i < dim; ++i) {
 156                 enum isl_lp_result res;
 157
 158                 isl_seq_clr(obj->el + 1, dim);
 159                 isl_int_set_si(obj->el[1 + i], 1);
 160                 res = isl_basic_set_solve_ilp(bset, max, obj->el,
 161                                                 &opt->el[1 + i], NULL);
 162                 if (res == isl_lp_empty)
 163                         goto empty;
 164                 assert(res == isl_lp_ok);
 165                 bset = isl_basic_set_fix(bset, isl_dim_set, i, opt->el[1 + i]);
 166         }
 167
 168         isl_basic_set_free(bset);
 169         isl_vec_free(obj);
 170
 171         return opt;
 172 empty:
 173         isl_vec_free(opt);
 174         opt = isl_vec_alloc(bset->ctx, 0);
 175         isl_basic_set_free(bset);
 176         isl_vec_free(obj);
 177
 178         return opt;
 179 }
 180
 181 struct isl_scan_pip {
 182         struct isl_scan_callback callback;
 183         isl_basic_set *bset;
 184         isl_set *sol;
 185         isl_set *empty;
 186         int stride;
 187         int n;
 188         int max;
 189 };
 190
 191 /* Check if the "manually" computed optimum of bset at the "sample"
 192  * values of the parameters agrees with the solution of pilp problem
 193  * represented by the pair (sol, empty).
 194  * In particular, if there is no solution for this value of the parameters,
 195  * then it should be an element of the parameter domain "empty".
 196  * Otherwise, the optimal solution, should be equal to the result of
 197  * plugging in the value of the parameters in "sol".
 198  */
 199 static int scan_one(struct isl_scan_callback *callback,
 200         __isl_take isl_vec *sample)
 201 {
 202         struct isl_scan_pip *sp = (struct isl_scan_pip *)callback;
 203         struct isl_vec *opt;
 204
 205         opt = opt_at(isl_basic_set_copy(sp->bset), isl_vec_copy(sample), sp->max);
 206         assert(opt);
 207
 208         if (opt->size == 0) {
 209                 isl_set *sample_set;
 210                 sample_set = isl_set_from_basic_set(
 211                                 isl_basic_set_from_vec(sample));
 212                 assert(isl_set_is_subset(sample_set, sp->empty));
 213                 isl_set_free(sample_set);
 214                 isl_vec_free(opt);
 215         } else {
 216                 isl_set *sol;
 217                 isl_set *opt_set;
 218                 opt_set = isl_set_from_basic_set(isl_basic_set_from_vec(opt));
 219                 sol = set_plug_in_parameters(isl_set_copy(sp->sol), sample);
 220                 assert(isl_set_is_equal(opt_set, sol));
 221                 isl_set_free(sol);
 222                 isl_set_free(opt_set);
 223         }
 224
 225         if (!(sp->n++ % sp->stride)) {
 226                 printf("o");
 227                 fflush(stdout);
 228         }
 229
 230         return 0;
 231 }
 232
 233 struct isl_scan_count {
 234         struct isl_scan_callback callback;
 235         int n;
 236 };
 237
 238 static int count_one(struct isl_scan_callback *callback,
 239         __isl_take isl_vec *sample)
 240 {
 241         struct isl_scan_count *sc = (struct isl_scan_count *)callback;
 242         isl_vec_free(sample);
 243
 244         sc->n++;
 245
 246         return 0;
 247 }
 248
 249 static void check_solution(isl_basic_set *bset, isl_basic_set *context,
 250         isl_set *sol, isl_set *empty, int max)
 251 {
 252         isl_basic_set *cube;
 253         struct isl_scan_count sc;
 254         struct isl_scan_pip sp;
 255         int i;
 256         int r;
 257
 258         context = isl_basic_set_underlying_set(context);
 259
 260         cube = construct_cube(context);
 261         context = isl_basic_set_intersect(context, cube);
 262
 263         sc.callback.add = count_one;
 264         sc.n = 0;
 265
 266         r = isl_basic_set_scan(isl_basic_set_copy(context), &sc.callback);
 267         assert (r == 0);
 268
 269         sp.callback.add = scan_one;
 270         sp.bset = bset;
 271         sp.sol = sol;
 272         sp.empty = empty;
 273         sp.n = 0;
 274         sp.stride = sc.n > 70 ? 1 + (sc.n + 1)/70 : 1;
 275         sp.max = max;
 276
 277         for (i = 0; i < sc.n; i += sp.stride)
 278                 printf(".");
 279         printf("\r");
 280         fflush(stdout);
 281
 282         r = isl_basic_set_scan(context, &sp.callback);
 283         assert(r == 0);
 284
 285         printf("\n");
 286
 287         isl_basic_set_free(bset);
 288 }
 289
 290 int main(int argc, char **argv)
 291 {
 292         struct isl_ctx *ctx;
 293         struct isl_basic_set *context, *bset, *copy, *context_copy;
 294         struct isl_set *set;
 295         struct isl_set *empty;
 296         int neg_one;
 297         char s[1024];
 298         int urs_parms = 0;
 299         int urs_unknowns = 0;
 300         int max = 0;
 301         int rational = 0;
 302         int n;
 303         struct pip_options *options;
 304
 305         options = pip_options_new_with_defaults();
 306         assert(options);
 307         argc = pip_options_parse(options, argc, argv);
 308
 309         ctx = isl_ctx_alloc_with_options(options->isl);
 310         options->isl = NULL;
 311
 312         context = isl_basic_set_read_from_file(ctx, stdin, 0, ISL_FORMAT_POLYLIB);
 313         assert(context);
 314         n = fscanf(stdin, "%d", &neg_one);
 315         assert(n == 1);
 316         assert(neg_one == -1);
 317         bset = isl_basic_set_read_from_file(ctx, stdin,
 318                 isl_basic_set_dim(context, isl_dim_set), ISL_FORMAT_POLYLIB);
 319
 320         while (fgets(s, sizeof(s), stdin)) {
 321                 if (strncasecmp(s, "Maximize", 8) == 0)
 322                         max = 1;
 323                 if (strncasecmp(s, "Rational", 8) == 0) {
 324                         rational = 1;
 325                         bset = isl_basic_set_set_rational(bset);
 326                 }
 327                 if (strncasecmp(s, "Urs_parms", 9) == 0)
 328                         urs_parms = 1;
 329                 if (strncasecmp(s, "Urs_unknowns", 12) == 0)
 330                         urs_unknowns = 1;
 331         }
 332         if (!urs_parms)
 333                 context = isl_basic_set_intersect(context,
 334                 isl_basic_set_positive_orthant(isl_basic_set_get_dim(context)));
 335         context = to_parameter_domain(context);
 336         if (!urs_unknowns)
 337                 bset = isl_basic_set_intersect(bset,
 338                 isl_basic_set_positive_orthant(isl_basic_set_get_dim(bset)));
 339
 340         if (options->verify) {
 341                 copy = isl_basic_set_copy(bset);
 342                 context_copy = isl_basic_set_copy(context);
 343         }
 344
 345         if (max)
 346                 set = isl_basic_set_partial_lexmax(bset, context, &empty);
 347         else
 348                 set = isl_basic_set_partial_lexmin(bset, context, &empty);
 349
 350         if (options->verify) {
 351                 assert(!rational);
 352                 check_solution(copy, context_copy, set, empty, max);
 353         } else {
 354                 isl_set_dump(set, stdout, 0);
 355                 fprintf(stdout, "no solution:\n");
 356                 isl_set_dump(empty, stdout, 4);
 357         }
 358
 359         isl_set_free(set);
 360         isl_set_free(empty);
 361         isl_ctx_free(ctx);
 362         free(options);
 363
 364         return 0;
 365 }