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