pip.c

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