pip.c

   1 /*
   2  * Copyright 2008-2009 Katholieke Universiteit Leuven
   3  *
   4  * Use of this software is governed by the MIT license
   5  *
   6  * Written by Sven Verdoolaege, K.U.Leuven, Departement
   7  * Computerwetenschappen, Celestijnenlaan 200A, B-3001 Leuven, Belgium
   8  */
   9
  10 #include <assert.h>
  11 #include <string.h>
  12 #include <isl_map_private.h>
  13 #include <isl/aff.h>
  14 #include <isl/set.h>
  15 #include "isl_sample.h"
  16 #include "isl_scan.h"
  17 #include <isl_seq.h>
  18 #include <isl_ilp_private.h>
  19 #include <isl/printer.h>
  20 #include <isl_point_private.h>
  21 #include <isl_vec_private.h>
  22 #include <isl/options.h>
  23 #include <isl_config.h>
  24
  25 /* The input of this program is the same as that of the "example" program
  26  * from the PipLib distribution, except that the "big parameter column"
  27  * should always be -1.
  28  *
  29  * Context constraints in PolyLib format
  30  * -1
  31  * Problem constraints in PolyLib format
  32  * Optional list of options
  33  *
  34  * The options are
  35  *      Maximize        compute maximum instead of minimum
  36  *      Rational        compute rational optimum instead of integer optimum
  37  *      Urs_parms       don't assume parameters are non-negative
  38  *      Urs_unknowns    don't assume unknowns are non-negative
  39  */
  40
  41 struct options {
  42         struct isl_options      *isl;
  43         unsigned                 verify;
  44         unsigned                 format;
  45 };
  46
  47 #define FORMAT_SET      0
  48 #define FORMAT_AFF      1
  49
  50 struct isl_arg_choice pip_format[] = {
  51         {"set",         FORMAT_SET},
  52         {"affine",      FORMAT_AFF},
  53         {0}
  54 };
  55
  56 ISL_ARGS_START(struct options, options_args)
  57 ISL_ARG_CHILD(struct options, isl, "isl", &isl_options_args, "isl options")
  58 ISL_ARG_BOOL(struct options, verify, 'T', "verify", 0, NULL)
  59 ISL_ARG_CHOICE(struct options, format, 0, "format",
  60         pip_format, FORMAT_SET, "output format")
  61 ISL_ARGS_END
  62
  63 ISL_ARG_DEF(options, struct options, options_args)
  64
  65 static __isl_give isl_basic_set *set_bounds(__isl_take isl_basic_set *bset)
  66 {
  67         isl_size nparam;
  68         int i, r;
  69         isl_point *pt, *pt2;
  70         isl_basic_set *box;
  71
  72         nparam = isl_basic_set_dim(bset, isl_dim_param);
  73         if (nparam < 0)
  74                 return isl_basic_set_free(bset);
  75         r = nparam >= 8 ? 4 : nparam >= 5 ? 6 : 30;
  76
  77         pt = isl_basic_set_sample_point(isl_basic_set_copy(bset));
  78         pt2 = isl_point_copy(pt);
  79
  80         for (i = 0; i < nparam; ++i) {
  81                 pt = isl_point_add_ui(pt, isl_dim_param, i, r);
  82                 pt2 = isl_point_sub_ui(pt2, isl_dim_param, i, r);
  83         }
  84
  85         box = isl_basic_set_box_from_points(pt, pt2);
  86
  87         return isl_basic_set_intersect(bset, box);
  88 }
  89
  90 static __isl_give isl_basic_set *to_parameter_domain(
  91         __isl_take isl_basic_set *context)
  92 {
  93         isl_size dim;
  94
  95         dim = isl_basic_set_dim(context, isl_dim_set);
  96         if (dim < 0)
  97                 return isl_basic_set_free(context);
  98         context = isl_basic_set_move_dims(context, isl_dim_param, 0,
  99                     isl_dim_set, 0, dim);
 100         context = isl_basic_set_params(context);
 101         return context;
 102 }
 103
 104 /* If "context" has more parameters than "bset", then reinterpret
 105  * the last dimensions of "bset" as parameters.
 106  */
 107 static __isl_give isl_basic_set *move_parameters(__isl_take isl_basic_set *bset,
 108         __isl_keep isl_basic_set *context)
 109 {
 110         isl_size nparam, nparam_bset, dim;
 111
 112         nparam = isl_basic_set_dim(context, isl_dim_param);
 113         nparam_bset = isl_basic_set_dim(bset, isl_dim_param);
 114         if (nparam < 0 || nparam_bset < 0)
 115                 return isl_basic_set_free(bset);
 116         if (nparam == nparam_bset)
 117                 return bset;
 118         dim = isl_basic_set_dim(bset, isl_dim_set);
 119         if (dim < 0)
 120                 return isl_basic_set_free(bset);
 121         bset = isl_basic_set_move_dims(bset, isl_dim_param, 0,
 122                                             isl_dim_set, dim - nparam, nparam);
 123         return bset;
 124 }
 125
 126 /* Plug in the initial values of "params" for the parameters in "bset" and
 127  * return the result.  The remaining entries in "params", if any,
 128  * correspond to the existentially quantified variables in the description
 129  * of the original context and can be ignored.
 130  */
 131 static __isl_give isl_basic_set *plug_in_parameters(
 132         __isl_take isl_basic_set *bset, __isl_take isl_vec *params)
 133 {
 134         int i;
 135         isl_size n;
 136
 137         n = isl_basic_set_dim(bset, isl_dim_param);
 138         if (n < 0)
 139                 bset = isl_basic_set_free(bset);
 140         for (i = 0; i < n; ++i)
 141                 bset = isl_basic_set_fix(bset,
 142                                          isl_dim_param, i, params->el[1 + i]);
 143
 144         bset = isl_basic_set_remove_dims(bset, isl_dim_param, 0, n);
 145
 146         isl_vec_free(params);
 147
 148         return bset;
 149 }
 150
 151 /* Plug in the initial values of "params" for the parameters in "set" and
 152  * return the result.  The remaining entries in "params", if any,
 153  * correspond to the existentially quantified variables in the description
 154  * of the original context and can be ignored.
 155  */
 156 static __isl_give isl_set *set_plug_in_parameters(__isl_take isl_set *set,
 157         __isl_take isl_vec *params)
 158 {
 159         int i;
 160         isl_size n;
 161
 162         n = isl_set_dim(set, isl_dim_param);
 163         if (n < 0)
 164                 set = isl_set_free(set);
 165         for (i = 0; i < n; ++i)
 166                 set = isl_set_fix(set, isl_dim_param, i, params->el[1 + i]);
 167
 168         set = isl_set_remove_dims(set, isl_dim_param, 0, n);
 169
 170         isl_vec_free(params);
 171
 172         return set;
 173 }
 174
 175 /* Compute the lexicographically minimal (or maximal if max is set)
 176  * element of bset for the given values of the parameters, by
 177  * successively solving an ilp problem in each direction.
 178  */
 179 static __isl_give isl_vec *opt_at(__isl_take isl_basic_set *bset,
 180         __isl_take isl_vec *params, int max)
 181 {
 182         isl_size dim;
 183         isl_ctx *ctx;
 184         struct isl_vec *opt;
 185         struct isl_vec *obj;
 186         int i;
 187
 188         dim = isl_basic_set_dim(bset, isl_dim_set);
 189         if (dim < 0)
 190                 goto error;
 191
 192         bset = plug_in_parameters(bset, params);
 193
 194         ctx = isl_basic_set_get_ctx(bset);
 195         if (isl_basic_set_plain_is_empty(bset)) {
 196                 opt = isl_vec_alloc(ctx, 0);
 197                 isl_basic_set_free(bset);
 198                 return opt;
 199         }
 200
 201         opt = isl_vec_alloc(ctx, 1 + dim);
 202         assert(opt);
 203
 204         obj = isl_vec_alloc(ctx, 1 + dim);
 205         assert(obj);
 206
 207         isl_int_set_si(opt->el[0], 1);
 208         isl_int_set_si(obj->el[0], 0);
 209
 210         for (i = 0; i < dim; ++i) {
 211                 enum isl_lp_result res;
 212
 213                 isl_seq_clr(obj->el + 1, dim);
 214                 isl_int_set_si(obj->el[1 + i], 1);
 215                 res = isl_basic_set_solve_ilp(bset, max, obj->el,
 216                                                 &opt->el[1 + i], NULL);
 217                 if (res == isl_lp_empty)
 218                         goto empty;
 219                 assert(res == isl_lp_ok);
 220                 bset = isl_basic_set_fix(bset, isl_dim_set, i, opt->el[1 + i]);
 221         }
 222
 223         isl_basic_set_free(bset);
 224         isl_vec_free(obj);
 225
 226         return opt;
 227 error:
 228         isl_basic_set_free(bset);
 229         isl_vec_free(params);
 230         return NULL;
 231 empty:
 232         isl_vec_free(opt);
 233         opt = isl_vec_alloc(ctx, 0);
 234         isl_basic_set_free(bset);
 235         isl_vec_free(obj);
 236
 237         return opt;
 238 }
 239
 240 struct isl_scan_pip {
 241         struct isl_scan_callback callback;
 242         isl_basic_set *bset;
 243         isl_set *sol;
 244         isl_set *empty;
 245         int stride;
 246         int n;
 247         int max;
 248 };
 249
 250 /* Check if the "manually" computed optimum of bset at the "sample"
 251  * values of the parameters agrees with the solution of pilp problem
 252  * represented by the pair (sol, empty).
 253  * In particular, if there is no solution for this value of the parameters,
 254  * then it should be an element of the parameter domain "empty".
 255  * Otherwise, the optimal solution, should be equal to the result of
 256  * plugging in the value of the parameters in "sol".
 257  */
 258 static isl_stat scan_one(struct isl_scan_callback *callback,
 259         __isl_take isl_vec *sample)
 260 {
 261         struct isl_scan_pip *sp = (struct isl_scan_pip *)callback;
 262         struct isl_vec *opt;
 263
 264         sp->n--;
 265
 266         opt = opt_at(isl_basic_set_copy(sp->bset), isl_vec_copy(sample), sp->max);
 267         assert(opt);
 268
 269         if (opt->size == 0) {
 270                 isl_point *sample_pnt;
 271                 sample_pnt = isl_point_alloc(isl_set_get_space(sp->empty), sample);
 272                 assert(isl_set_contains_point(sp->empty, sample_pnt));
 273                 isl_point_free(sample_pnt);
 274                 isl_vec_free(opt);
 275         } else {
 276                 isl_set *sol;
 277                 isl_set *opt_set;
 278                 opt_set = isl_set_from_basic_set(isl_basic_set_from_vec(opt));
 279                 sol = set_plug_in_parameters(isl_set_copy(sp->sol), sample);
 280                 assert(isl_set_is_equal(opt_set, sol));
 281                 isl_set_free(sol);
 282                 isl_set_free(opt_set);
 283         }
 284
 285         if (!(sp->n % sp->stride)) {
 286                 printf("o");
 287                 fflush(stdout);
 288         }
 289
 290         return sp->n >= 1 ? isl_stat_ok : isl_stat_error;
 291 }
 292
 293 static void check_solution(isl_basic_set *bset, isl_basic_set *context,
 294         isl_set *sol, isl_set *empty, int max)
 295 {
 296         struct isl_scan_pip sp;
 297         isl_int count, count_max;
 298         int i, n;
 299         int r;
 300
 301         context = set_bounds(context);
 302         context = isl_basic_set_underlying_set(context);
 303
 304         isl_int_init(count);
 305         isl_int_init(count_max);
 306
 307         isl_int_set_si(count_max, 2000);
 308         r = isl_basic_set_count_upto(context, count_max, &count);
 309         assert(r >= 0);
 310         n = isl_int_get_si(count);
 311
 312         isl_int_clear(count_max);
 313         isl_int_clear(count);
 314
 315         sp.callback.add = scan_one;
 316         sp.bset = bset;
 317         sp.sol = sol;
 318         sp.empty = empty;
 319         sp.n = n;
 320         sp.stride = n > 70 ? 1 + (n + 1)/70 : 1;
 321         sp.max = max;
 322
 323         for (i = 0; i < n; i += sp.stride)
 324                 printf(".");
 325         printf("\r");
 326         fflush(stdout);
 327
 328         isl_basic_set_scan(context, &sp.callback);
 329
 330         printf("\n");
 331
 332         isl_basic_set_free(bset);
 333 }
 334
 335 int main(int argc, char **argv)
 336 {
 337         struct isl_ctx *ctx;
 338         struct isl_basic_set *context, *bset, *copy, *context_copy;
 339         struct isl_set *set = NULL;
 340         struct isl_set *empty;
 341         isl_pw_multi_aff *pma = NULL;
 342         int neg_one;
 343         char s[1024];
 344         int urs_parms = 0;
 345         int urs_unknowns = 0;
 346         int max = 0;
 347         int rational = 0;
 348         int n;
 349         struct options *options;
 350
 351         options = options_new_with_defaults();
 352         assert(options);
 353         argc = options_parse(options, argc, argv, ISL_ARG_ALL);
 354
 355         ctx = isl_ctx_alloc_with_options(&options_args, options);
 356
 357         context = isl_basic_set_read_from_file(ctx, stdin);
 358         assert(context);
 359         n = fscanf(stdin, "%d", &neg_one);
 360         assert(n == 1);
 361         assert(neg_one == -1);
 362         bset = isl_basic_set_read_from_file(ctx, stdin);
 363
 364         while (fgets(s, sizeof(s), stdin)) {
 365                 if (strncasecmp(s, "Maximize", 8) == 0)
 366                         max = 1;
 367                 if (strncasecmp(s, "Rational", 8) == 0) {
 368                         rational = 1;
 369                         bset = isl_basic_set_set_rational(bset);
 370                 }
 371                 if (strncasecmp(s, "Urs_parms", 9) == 0)
 372                         urs_parms = 1;
 373                 if (strncasecmp(s, "Urs_unknowns", 12) == 0)
 374                         urs_unknowns = 1;
 375         }
 376         if (!urs_parms)
 377                 context = isl_basic_set_intersect(context,
 378                 isl_basic_set_positive_orthant(isl_basic_set_get_space(context)));
 379         context = to_parameter_domain(context);
 380         bset = move_parameters(bset, context);
 381         if (!urs_unknowns)
 382                 bset = isl_basic_set_intersect(bset,
 383                 isl_basic_set_positive_orthant(isl_basic_set_get_space(bset)));
 384
 385         if (options->verify) {
 386                 copy = isl_basic_set_copy(bset);
 387                 context_copy = isl_basic_set_copy(context);
 388         }
 389
 390         if (options->format == FORMAT_AFF) {
 391                 if (max)
 392                         pma = isl_basic_set_partial_lexmax_pw_multi_aff(bset,
 393                                                                 context, &empty);
 394                 else
 395                         pma = isl_basic_set_partial_lexmin_pw_multi_aff(bset,
 396                                                                 context, &empty);
 397         } else {
 398                 if (max)
 399                         set = isl_basic_set_partial_lexmax(bset,
 400                                                                 context, &empty);
 401                 else
 402                         set = isl_basic_set_partial_lexmin(bset,
 403                                                                 context, &empty);
 404         }
 405
 406         if (options->verify) {
 407                 assert(!rational);
 408                 if (options->format == FORMAT_AFF)
 409                         set = isl_set_from_pw_multi_aff(pma);
 410                 check_solution(copy, context_copy, set, empty, max);
 411                 isl_set_free(set);
 412         } else {
 413                 isl_printer *p;
 414                 p = isl_printer_to_file(ctx, stdout);
 415                 if (options->format == FORMAT_AFF)
 416                         p = isl_printer_print_pw_multi_aff(p, pma);
 417                 else
 418                         p = isl_printer_print_set(p, set);
 419                 p = isl_printer_end_line(p);
 420                 p = isl_printer_print_str(p, "no solution: ");
 421                 p = isl_printer_print_set(p, empty);
 422                 p = isl_printer_end_line(p);
 423                 isl_printer_free(p);
 424                 isl_set_free(set);
 425                 isl_pw_multi_aff_free(pma);
 426         }
 427
 428         isl_set_free(empty);
 429         isl_ctx_free(ctx);
 430
 431         return 0;
 432 }