isl_map.c

   1 #include <string.h>
   2 #include <strings.h>
   3 #include "isl_ctx.h"
   4 #include "isl_blk.h"
   5 #include "isl_equalities.h"
   6 #include "isl_list.h"
   7 #include "isl_lp.h"
   8 #include "isl_seq.h"
   9 #include "isl_set.h"
  10 #include "isl_map.h"
  11 #include "isl_map_private.h"
  12 #include "isl_map_piplib.h"
  13 #include "isl_sample.h"
  14 #include "isl_vec.h"
  15
  16 static struct isl_basic_map *basic_map_init(struct isl_ctx *ctx,
  17                 struct isl_basic_map *bmap,
  18                 unsigned nparam, unsigned n_in, unsigned n_out, unsigned extra,
  19                 unsigned n_eq, unsigned n_ineq)
  20 {
  21         int i;
  22         size_t row_size = 1 + nparam + n_in + n_out + extra;
  23
  24         bmap->block = isl_blk_alloc(ctx, (n_eq + n_ineq) * row_size);
  25         if (isl_blk_is_error(bmap->block)) {
  26                 free(bmap);
  27                 return NULL;
  28         }
  29
  30         bmap->eq = isl_alloc_array(ctx, isl_int *, n_eq + n_ineq);
  31         if (!bmap->eq) {
  32                 isl_blk_free(ctx, bmap->block);
  33                 free(bmap);
  34                 return NULL;
  35         }
  36
  37         if (extra == 0) {
  38                 bmap->block2 = isl_blk_empty();
  39                 bmap->div = NULL;
  40         } else {
  41                 bmap->block2 = isl_blk_alloc(ctx, extra * (1 + row_size));
  42                 if (isl_blk_is_error(bmap->block2)) {
  43                         free(bmap->eq);
  44                         isl_blk_free(ctx, bmap->block);
  45                         free(bmap);
  46                         return NULL;
  47                 }
  48
  49                 bmap->div = isl_alloc_array(ctx, isl_int *, extra);
  50                 if (!bmap->div) {
  51                         isl_blk_free(ctx, bmap->block2);
  52                         free(bmap->eq);
  53                         isl_blk_free(ctx, bmap->block);
  54                         free(bmap);
  55                         return NULL;
  56                 }
  57         }
  58
  59         for (i = 0; i < n_eq + n_ineq; ++i)
  60                 bmap->eq[i] = bmap->block.data + i * row_size;
  61
  62         for (i = 0; i < extra; ++i)
  63                 bmap->div[i] = bmap->block2.data + i * (1 + row_size);
  64
  65         bmap->ctx = ctx;
  66         isl_ctx_ref(ctx);
  67         bmap->ref = 1;
  68         bmap->flags = 0;
  69         bmap->c_size = n_eq + n_ineq;
  70         bmap->ineq = bmap->eq + n_eq;
  71         bmap->nparam = nparam;
  72         bmap->n_in = n_in;
  73         bmap->n_out = n_out;
  74         bmap->extra = extra;
  75         bmap->n_eq = 0;
  76         bmap->n_ineq = 0;
  77         bmap->n_div = 0;
  78         bmap->sample = NULL;
  79
  80         return bmap;
  81 }
  82
  83 struct isl_basic_set *isl_basic_set_alloc(struct isl_ctx *ctx,
  84                 unsigned nparam, unsigned dim, unsigned extra,
  85                 unsigned n_eq, unsigned n_ineq)
  86 {
  87         struct isl_basic_map *bmap;
  88         bmap = isl_basic_map_alloc(ctx, nparam, 0, dim, extra, n_eq, n_ineq);
  89         return (struct isl_basic_set *)bmap;
  90 }
  91
  92 struct isl_basic_map *isl_basic_map_alloc(struct isl_ctx *ctx,
  93                 unsigned nparam, unsigned in, unsigned out, unsigned extra,
  94                 unsigned n_eq, unsigned n_ineq)
  95 {
  96         struct isl_basic_map *bmap;
  97
  98         bmap = isl_alloc_type(ctx, struct isl_basic_map);
  99         if (!bmap)
 100                 return NULL;
 101
 102         return basic_map_init(ctx, bmap,
 103                                nparam, in, out, extra, n_eq, n_ineq);
 104 }
 105
 106 static void dup_constraints(
 107                 struct isl_basic_map *dst, struct isl_basic_map *src)
 108 {
 109         int i;
 110         unsigned total = src->nparam + src->n_in + src->n_out + src->n_div;
 111
 112         for (i = 0; i < src->n_eq; ++i) {
 113                 int j = isl_basic_map_alloc_equality(dst);
 114                 isl_seq_cpy(dst->eq[j], src->eq[i], 1+total);
 115         }
 116
 117         for (i = 0; i < src->n_ineq; ++i) {
 118                 int j = isl_basic_map_alloc_inequality(dst);
 119                 isl_seq_cpy(dst->ineq[j], src->ineq[i], 1+total);
 120         }
 121
 122         for (i = 0; i < src->n_div; ++i) {
 123                 int j = isl_basic_map_alloc_div(dst);
 124                 isl_seq_cpy(dst->div[j], src->div[i], 1+1+total);
 125         }
 126         F_SET(dst, ISL_BASIC_SET_FINAL);
 127 }
 128
 129 struct isl_basic_map *isl_basic_map_dup(struct isl_basic_map *bmap)
 130 {
 131         struct isl_basic_map *dup;
 132
 133         if (!bmap)
 134                 return NULL;
 135         dup = isl_basic_map_alloc(bmap->ctx, bmap->nparam,
 136                         bmap->n_in, bmap->n_out,
 137                         bmap->n_div, bmap->n_eq, bmap->n_ineq);
 138         if (!dup)
 139                 return NULL;
 140         dup->flags = bmap->flags;
 141         dup_constraints(dup, bmap);
 142         dup->sample = isl_vec_copy(bmap->ctx, bmap->sample);
 143         return dup;
 144 }
 145
 146 struct isl_basic_set *isl_basic_set_dup(struct isl_basic_set *bset)
 147 {
 148         struct isl_basic_map *dup;
 149
 150         dup = isl_basic_map_dup((struct isl_basic_map *)bset);
 151         return (struct isl_basic_set *)dup;
 152 }
 153
 154 struct isl_basic_set *isl_basic_set_copy(struct isl_basic_set *bset)
 155 {
 156         if (!bset)
 157                 return NULL;
 158
 159         if (F_ISSET(bset, ISL_BASIC_SET_FINAL)) {
 160                 bset->ref++;
 161                 return bset;
 162         }
 163         return isl_basic_set_dup(bset);
 164 }
 165
 166 struct isl_set *isl_set_copy(struct isl_set *set)
 167 {
 168         if (!set)
 169                 return NULL;
 170
 171         set->ref++;
 172         return set;
 173 }
 174
 175 struct isl_basic_map *isl_basic_map_copy(struct isl_basic_map *bmap)
 176 {
 177         if (!bmap)
 178                 return NULL;
 179
 180         if (F_ISSET(bmap, ISL_BASIC_SET_FINAL)) {
 181                 bmap->ref++;
 182                 return bmap;
 183         }
 184         return isl_basic_map_dup(bmap);
 185 }
 186
 187 struct isl_map *isl_map_copy(struct isl_map *map)
 188 {
 189         if (!map)
 190                 return NULL;
 191
 192         map->ref++;
 193         return map;
 194 }
 195
 196 void isl_basic_map_free(struct isl_basic_map *bmap)
 197 {
 198         if (!bmap)
 199                 return;
 200
 201         if (--bmap->ref > 0)
 202                 return;
 203
 204         isl_ctx_deref(bmap->ctx);
 205         free(bmap->div);
 206         isl_blk_free(bmap->ctx, bmap->block2);
 207         free(bmap->eq);
 208         isl_blk_free(bmap->ctx, bmap->block);
 209         isl_vec_free(bmap->ctx, bmap->sample);
 210         free(bmap);
 211 }
 212
 213 void isl_basic_set_free(struct isl_basic_set *bset)
 214 {
 215         isl_basic_map_free((struct isl_basic_map *)bset);
 216 }
 217
 218 int isl_basic_map_alloc_equality(struct isl_basic_map *bmap)
 219 {
 220         struct isl_ctx *ctx;
 221         if (!bmap)
 222                 return -1;
 223         ctx = bmap->ctx;
 224         isl_assert(ctx, bmap->n_eq + bmap->n_ineq < bmap->c_size, return -1);
 225         isl_assert(ctx, bmap->eq + bmap->n_eq <= bmap->ineq, return -1);
 226         if (bmap->eq + bmap->n_eq == bmap->ineq) {
 227                 isl_int *t;
 228                 int j = isl_basic_map_alloc_inequality(bmap);
 229                 if (j < 0)
 230                         return -1;
 231                 t = bmap->ineq[0];
 232                 bmap->ineq[0] = bmap->ineq[j];
 233                 bmap->ineq[j] = t;
 234                 bmap->n_ineq--;
 235                 bmap->ineq++;
 236         } else
 237                 isl_seq_clr(bmap->eq[bmap->n_eq] +
 238                       1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div,
 239                       bmap->extra - bmap->n_div);
 240         return bmap->n_eq++;
 241 }
 242
 243 int isl_basic_set_alloc_equality(struct isl_basic_set *bset)
 244 {
 245         return isl_basic_map_alloc_equality((struct isl_basic_map *)bset);
 246 }
 247
 248 int isl_basic_map_free_equality(struct isl_basic_map *bmap, unsigned n)
 249 {
 250         if (!bmap)
 251                 return -1;
 252         isl_assert(bmap->ctx, n <= bmap->n_eq, return -1);
 253         bmap->n_eq -= n;
 254         return 0;
 255 }
 256
 257 int isl_basic_map_drop_equality(struct isl_basic_map *bmap, unsigned pos)
 258 {
 259         isl_int *t;
 260         if (!bmap)
 261                 return -1;
 262         isl_assert(bmap->ctx, pos < bmap->n_eq, return -1);
 263
 264         if (pos != bmap->n_eq - 1) {
 265                 t = bmap->eq[pos];
 266                 bmap->eq[pos] = bmap->eq[bmap->n_eq - 1];
 267                 bmap->eq[bmap->n_eq - 1] = t;
 268         }
 269         bmap->n_eq--;
 270         return 0;
 271 }
 272
 273 int isl_basic_set_drop_equality(struct isl_basic_set *bset, unsigned pos)
 274 {
 275         return isl_basic_map_drop_equality((struct isl_basic_map *)bset, pos);
 276 }
 277
 278 void isl_basic_map_inequality_to_equality(
 279                 struct isl_basic_map *bmap, unsigned pos)
 280 {
 281         isl_int *t;
 282
 283         t = bmap->ineq[pos];
 284         bmap->ineq[pos] = bmap->ineq[0];
 285         bmap->ineq[0] = bmap->eq[bmap->n_eq];
 286         bmap->eq[bmap->n_eq] = t;
 287         bmap->n_eq++;
 288         bmap->n_ineq--;
 289         bmap->ineq++;
 290         F_CLR(bmap, ISL_BASIC_MAP_NORMALIZED);
 291 }
 292
 293 int isl_basic_map_alloc_inequality(struct isl_basic_map *bmap)
 294 {
 295         struct isl_ctx *ctx;
 296         if (!bmap)
 297                 return -1;
 298         ctx = bmap->ctx;
 299         isl_assert(ctx, (bmap->ineq - bmap->eq) + bmap->n_ineq < bmap->c_size,
 300                         return -1);
 301         F_CLR(bmap, ISL_BASIC_MAP_NO_IMPLICIT);
 302         F_CLR(bmap, ISL_BASIC_MAP_NO_REDUNDANT);
 303         F_CLR(bmap, ISL_BASIC_MAP_NORMALIZED);
 304         isl_seq_clr(bmap->ineq[bmap->n_ineq] +
 305                       1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div,
 306                       bmap->extra - bmap->n_div);
 307         return bmap->n_ineq++;
 308 }
 309
 310 int isl_basic_set_alloc_inequality(struct isl_basic_set *bset)
 311 {
 312         return isl_basic_map_alloc_inequality((struct isl_basic_map *)bset);
 313 }
 314
 315 int isl_basic_map_free_inequality(struct isl_basic_map *bmap, unsigned n)
 316 {
 317         if (!bmap)
 318                 return -1;
 319         isl_assert(bmap->ctx, n <= bmap->n_ineq, return -1);
 320         bmap->n_ineq -= n;
 321         return 0;
 322 }
 323
 324 int isl_basic_set_free_inequality(struct isl_basic_set *bset, unsigned n)
 325 {
 326         return isl_basic_map_free_inequality((struct isl_basic_map *)bset, n);
 327 }
 328
 329 int isl_basic_map_drop_inequality(struct isl_basic_map *bmap, unsigned pos)
 330 {
 331         isl_int *t;
 332         if (!bmap)
 333                 return -1;
 334         isl_assert(bmap->ctx, pos < bmap->n_ineq, return -1);
 335
 336         if (pos != bmap->n_ineq - 1) {
 337                 t = bmap->ineq[pos];
 338                 bmap->ineq[pos] = bmap->ineq[bmap->n_ineq - 1];
 339                 bmap->ineq[bmap->n_ineq - 1] = t;
 340                 F_CLR(bmap, ISL_BASIC_MAP_NORMALIZED);
 341         }
 342         bmap->n_ineq--;
 343         return 0;
 344 }
 345
 346 int isl_basic_set_drop_inequality(struct isl_basic_set *bset, unsigned pos)
 347 {
 348         return isl_basic_map_drop_inequality((struct isl_basic_map *)bset, pos);
 349 }
 350
 351 int isl_basic_map_alloc_div(struct isl_basic_map *bmap)
 352 {
 353         if (!bmap)
 354                 return -1;
 355         isl_assert(bmap->ctx, bmap->n_div < bmap->extra, return -1);
 356         return bmap->n_div++;
 357 }
 358
 359 int isl_basic_map_free_div(struct isl_basic_map *bmap, unsigned n)
 360 {
 361         if (!bmap)
 362                 return -1;
 363         isl_assert(bmap->ctx, n <= bmap->n_div, return -1);
 364         bmap->n_div -= n;
 365         return 0;
 366 }
 367
 368 /* Copy constraint from src to dst, putting the vars of src at offset
 369  * dim_off in dst and the divs of src at offset div_off in dst.
 370  * If both sets are actually map, then dim_off applies to the input
 371  * variables.
 372  */
 373 static void copy_constraint(struct isl_basic_map *dst_map, isl_int *dst,
 374                             struct isl_basic_map *src_map, isl_int *src,
 375                             unsigned in_off, unsigned out_off, unsigned div_off)
 376 {
 377         isl_int_set(dst[0], src[0]);
 378         isl_seq_cpy(dst+1, src+1, isl_min(dst_map->nparam, src_map->nparam));
 379         if (dst_map->nparam > src_map->nparam)
 380                 isl_seq_clr(dst+1+src_map->nparam,
 381                                 dst_map->nparam - src_map->nparam);
 382         isl_seq_clr(dst+1+dst_map->nparam, in_off);
 383         isl_seq_cpy(dst+1+dst_map->nparam+in_off,
 384                     src+1+src_map->nparam,
 385                     isl_min(dst_map->n_in-in_off, src_map->n_in));
 386         if (dst_map->n_in-in_off > src_map->n_in)
 387                 isl_seq_clr(dst+1+dst_map->nparam+in_off+src_map->n_in,
 388                                 dst_map->n_in - in_off - src_map->n_in);
 389         isl_seq_clr(dst+1+dst_map->nparam+dst_map->n_in, out_off);
 390         isl_seq_cpy(dst+1+dst_map->nparam+dst_map->n_in+out_off,
 391                     src+1+src_map->nparam+src_map->n_in,
 392                     isl_min(dst_map->n_out-out_off, src_map->n_out));
 393         if (dst_map->n_out-out_off > src_map->n_out)
 394                 isl_seq_clr(dst+1+dst_map->nparam+dst_map->n_in+out_off+src_map->n_out,
 395                                 dst_map->n_out - out_off - src_map->n_out);
 396         isl_seq_clr(dst+1+dst_map->nparam+dst_map->n_in+dst_map->n_out, div_off);
 397         isl_seq_cpy(dst+1+dst_map->nparam+dst_map->n_in+dst_map->n_out+div_off,
 398                     src+1+src_map->nparam+src_map->n_in+src_map->n_out,
 399                     isl_min(dst_map->extra-div_off, src_map->n_div));
 400         if (dst_map->n_div-div_off > src_map->n_div)
 401                 isl_seq_clr(dst+1+dst_map->nparam+dst_map->n_in+dst_map->n_out+
 402                                 div_off+src_map->n_div,
 403                                 dst_map->n_div - div_off - src_map->n_div);
 404 }
 405
 406 static void copy_div(struct isl_basic_map *dst_map, isl_int *dst,
 407                      struct isl_basic_map *src_map, isl_int *src,
 408                      unsigned in_off, unsigned out_off, unsigned div_off)
 409 {
 410         isl_int_set(dst[0], src[0]);
 411         copy_constraint(dst_map, dst+1, src_map, src+1, in_off, out_off, div_off);
 412 }
 413
 414 static struct isl_basic_map *add_constraints(struct isl_basic_map *bmap1,
 415                 struct isl_basic_map *bmap2, unsigned i_pos, unsigned o_pos)
 416 {
 417         int i;
 418         unsigned div_off;
 419
 420         if (!bmap1 || !bmap2)
 421                 goto error;
 422
 423         div_off = bmap1->n_div;
 424
 425         for (i = 0; i < bmap2->n_eq; ++i) {
 426                 int i1 = isl_basic_map_alloc_equality(bmap1);
 427                 if (i1 < 0)
 428                         goto error;
 429                 copy_constraint(bmap1, bmap1->eq[i1], bmap2, bmap2->eq[i],
 430                                 i_pos, o_pos, div_off);
 431         }
 432
 433         for (i = 0; i < bmap2->n_ineq; ++i) {
 434                 int i1 = isl_basic_map_alloc_inequality(bmap1);
 435                 if (i1 < 0)
 436                         goto error;
 437                 copy_constraint(bmap1, bmap1->ineq[i1], bmap2, bmap2->ineq[i],
 438                                 i_pos, o_pos, div_off);
 439         }
 440
 441         for (i = 0; i < bmap2->n_div; ++i) {
 442                 int i1 = isl_basic_map_alloc_div(bmap1);
 443                 if (i1 < 0)
 444                         goto error;
 445                 copy_div(bmap1, bmap1->div[i1], bmap2, bmap2->div[i],
 446                          i_pos, o_pos, div_off);
 447         }
 448
 449         isl_basic_map_free(bmap2);
 450
 451         return bmap1;
 452
 453 error:
 454         isl_basic_map_free(bmap1);
 455         isl_basic_map_free(bmap2);
 456         return NULL;
 457 }
 458
 459 static struct isl_basic_set *set_add_constraints(struct isl_basic_set *bset1,
 460                 struct isl_basic_set *bset2, unsigned pos)
 461 {
 462         return (struct isl_basic_set *)
 463                 add_constraints((struct isl_basic_map *)bset1,
 464                                 (struct isl_basic_map *)bset2, 0, pos);
 465 }
 466
 467 struct isl_basic_map *isl_basic_map_extend(struct isl_basic_map *base,
 468                 unsigned nparam, unsigned n_in, unsigned n_out, unsigned extra,
 469                 unsigned n_eq, unsigned n_ineq)
 470 {
 471         struct isl_basic_map *ext;
 472         unsigned flags;
 473         int dims_ok;
 474
 475         base = isl_basic_map_cow(base);
 476         if (!base)
 477                 goto error;
 478
 479         dims_ok = base && base->nparam == nparam &&
 480                   base->n_in == n_in && base->n_out == n_out &&
 481                   base->extra >= base->n_div + extra;
 482
 483         if (dims_ok && n_eq == 0 && n_ineq == 0)
 484                 return base;
 485
 486         isl_assert(base->ctx, base->nparam <= nparam, goto error);
 487         isl_assert(base->ctx, base->n_in <= n_in, goto error);
 488         isl_assert(base->ctx, base->n_out <= n_out, goto error);
 489         extra += base->extra;
 490         n_eq += base->n_eq;
 491         n_ineq += base->n_ineq;
 492
 493         ext = isl_basic_map_alloc(base->ctx, nparam, n_in, n_out,
 494                                         extra, n_eq, n_ineq);
 495         if (!base)
 496                 return ext;
 497         if (!ext)
 498                 goto error;
 499
 500         flags = base->flags;
 501         ext = add_constraints(ext, base, 0, 0);
 502         if (ext) {
 503                 ext->flags = flags;
 504                 F_CLR(ext, ISL_BASIC_SET_FINAL);
 505         }
 506
 507         return ext;
 508
 509 error:
 510         isl_basic_map_free(base);
 511         return NULL;
 512 }
 513
 514 struct isl_basic_set *isl_basic_set_extend(struct isl_basic_set *base,
 515                 unsigned nparam, unsigned dim, unsigned extra,
 516                 unsigned n_eq, unsigned n_ineq)
 517 {
 518         return (struct isl_basic_set *)
 519                 isl_basic_map_extend((struct isl_basic_map *)base,
 520                                         nparam, 0, dim, extra, n_eq, n_ineq);
 521 }
 522
 523 struct isl_basic_set *isl_basic_set_cow(struct isl_basic_set *bset)
 524 {
 525         return (struct isl_basic_set *)
 526                 isl_basic_map_cow((struct isl_basic_map *)bset);
 527 }
 528
 529 struct isl_basic_map *isl_basic_map_cow(struct isl_basic_map *bmap)
 530 {
 531         if (!bmap)
 532                 return NULL;
 533
 534         if (bmap->ref > 1) {
 535                 bmap->ref--;
 536                 bmap = isl_basic_map_dup(bmap);
 537         }
 538         F_CLR(bmap, ISL_BASIC_SET_FINAL);
 539         return bmap;
 540 }
 541
 542 struct isl_set *isl_set_cow(struct isl_set *set)
 543 {
 544         if (!set)
 545                 return NULL;
 546
 547         if (set->ref == 1)
 548                 return set;
 549         set->ref--;
 550         return isl_set_dup(set);
 551 }
 552
 553 struct isl_map *isl_map_cow(struct isl_map *map)
 554 {
 555         if (!map)
 556                 return NULL;
 557
 558         if (map->ref == 1)
 559                 return map;
 560         map->ref--;
 561         return isl_map_dup(map);
 562 }
 563
 564 static void swap_vars(struct isl_blk blk, isl_int *a,
 565                         unsigned a_len, unsigned b_len)
 566 {
 567         isl_seq_cpy(blk.data, a+a_len, b_len);
 568         isl_seq_cpy(blk.data+b_len, a, a_len);
 569         isl_seq_cpy(a, blk.data, b_len+a_len);
 570 }
 571
 572 struct isl_basic_set *isl_basic_set_swap_vars(
 573                 struct isl_basic_set *bset, unsigned n)
 574 {
 575         int i;
 576         struct isl_blk blk;
 577
 578         if (!bset)
 579                 goto error;
 580
 581         isl_assert(bset->ctx, n <= bset->dim, goto error);
 582
 583         if (n == bset->dim)
 584                 return bset;
 585
 586         bset = isl_basic_set_cow(bset);
 587         if (!bset)
 588                 return NULL;
 589
 590         blk = isl_blk_alloc(bset->ctx, bset->dim);
 591         if (isl_blk_is_error(blk))
 592                 goto error;
 593
 594         for (i = 0; i < bset->n_eq; ++i)
 595                 swap_vars(blk,
 596                           bset->eq[i]+1+bset->nparam, n, bset->dim - n);
 597
 598         for (i = 0; i < bset->n_ineq; ++i)
 599                 swap_vars(blk,
 600                           bset->ineq[i]+1+bset->nparam, n, bset->dim - n);
 601
 602         for (i = 0; i < bset->n_div; ++i)
 603                 swap_vars(blk,
 604                           bset->div[i]+1+1+bset->nparam, n, bset->dim - n);
 605
 606         isl_blk_free(bset->ctx, blk);
 607
 608         F_CLR(bset, ISL_BASIC_SET_NORMALIZED);
 609         return bset;
 610
 611 error:
 612         isl_basic_set_free(bset);
 613         return NULL;
 614 }
 615
 616 struct isl_set *isl_set_swap_vars(struct isl_set *set, unsigned n)
 617 {
 618         int i;
 619         set = isl_set_cow(set);
 620         if (!set)
 621                 return NULL;
 622
 623         for (i = 0; i < set->n; ++i) {
 624                 set->p[i] = isl_basic_set_swap_vars(set->p[i], n);
 625                 if (!set->p[i]) {
 626                         isl_set_free(set);
 627                         return NULL;
 628                 }
 629         }
 630         F_CLR(set, ISL_SET_NORMALIZED);
 631         return set;
 632 }
 633
 634 static void constraint_drop_vars(isl_int *c, unsigned n, unsigned rem)
 635 {
 636         isl_seq_cpy(c, c + n, rem);
 637         isl_seq_clr(c + rem, n);
 638 }
 639
 640 struct isl_basic_set *isl_basic_set_drop_dims(
 641                 struct isl_basic_set *bset, unsigned first, unsigned n)
 642 {
 643         int i;
 644
 645         if (!bset)
 646                 goto error;
 647
 648         isl_assert(bset->ctx, first + n <= bset->dim, goto error);
 649
 650         if (n == 0)
 651                 return bset;
 652
 653         bset = isl_basic_set_cow(bset);
 654         if (!bset)
 655                 return NULL;
 656
 657         for (i = 0; i < bset->n_eq; ++i)
 658                 constraint_drop_vars(bset->eq[i]+1+bset->nparam+first, n,
 659                                      (bset->dim-first-n)+bset->extra);
 660
 661         for (i = 0; i < bset->n_ineq; ++i)
 662                 constraint_drop_vars(bset->ineq[i]+1+bset->nparam+first, n,
 663                                      (bset->dim-first-n)+bset->extra);
 664
 665         for (i = 0; i < bset->n_div; ++i)
 666                 constraint_drop_vars(bset->div[i]+1+1+bset->nparam+first, n,
 667                                      (bset->dim-first-n)+bset->extra);
 668
 669         bset->dim -= n;
 670         bset->extra += n;
 671
 672         F_CLR(bset, ISL_BASIC_SET_NORMALIZED);
 673         bset = isl_basic_set_simplify(bset);
 674         return isl_basic_set_finalize(bset);
 675 error:
 676         isl_basic_set_free(bset);
 677         return NULL;
 678 }
 679
 680 static struct isl_set *isl_set_drop_dims(
 681                 struct isl_set *set, unsigned first, unsigned n)
 682 {
 683         int i;
 684
 685         if (!set)
 686                 goto error;
 687
 688         isl_assert(set->ctx, first + n <= set->dim, goto error);
 689
 690         if (n == 0)
 691                 return set;
 692         set = isl_set_cow(set);
 693         if (!set)
 694                 goto error;
 695
 696         for (i = 0; i < set->n; ++i) {
 697                 set->p[i] = isl_basic_set_drop_dims(set->p[i], first, n);
 698                 if (!set->p[i])
 699                         goto error;
 700         }
 701         set->dim -= n;
 702
 703         F_CLR(set, ISL_SET_NORMALIZED);
 704         return set;
 705 error:
 706         isl_set_free(set);
 707         return NULL;
 708 }
 709
 710 struct isl_basic_map *isl_basic_map_drop_inputs(
 711                 struct isl_basic_map *bmap, unsigned first, unsigned n)
 712 {
 713         int i;
 714
 715         if (!bmap)
 716                 goto error;
 717
 718         isl_assert(bmap->ctx, first + n <= bmap->n_in, goto error);
 719
 720         if (n == 0)
 721                 return bmap;
 722
 723         bmap = isl_basic_map_cow(bmap);
 724         if (!bmap)
 725                 return NULL;
 726
 727         for (i = 0; i < bmap->n_eq; ++i)
 728                 constraint_drop_vars(bmap->eq[i]+1+bmap->nparam+first, n,
 729                                  (bmap->n_in-first-n)+bmap->n_out+bmap->extra);
 730
 731         for (i = 0; i < bmap->n_ineq; ++i)
 732                 constraint_drop_vars(bmap->ineq[i]+1+bmap->nparam+first, n,
 733                                  (bmap->n_in-first-n)+bmap->n_out+bmap->extra);
 734
 735         for (i = 0; i < bmap->n_div; ++i)
 736                 constraint_drop_vars(bmap->div[i]+1+1+bmap->nparam+first, n,
 737                                  (bmap->n_in-first-n)+bmap->n_out+bmap->extra);
 738
 739         bmap->n_in -= n;
 740         bmap->extra += n;
 741
 742         F_CLR(bmap, ISL_BASIC_MAP_NORMALIZED);
 743         bmap = isl_basic_map_simplify(bmap);
 744         return isl_basic_map_finalize(bmap);
 745 error:
 746         isl_basic_map_free(bmap);
 747         return NULL;
 748 }
 749
 750 static struct isl_map *isl_map_drop_inputs(
 751                 struct isl_map *map, unsigned first, unsigned n)
 752 {
 753         int i;
 754
 755         if (!map)
 756                 goto error;
 757
 758         isl_assert(map->ctx, first + n <= map->n_in, goto error);
 759
 760         if (n == 0)
 761                 return map;
 762         map = isl_map_cow(map);
 763         if (!map)
 764                 goto error;
 765
 766         for (i = 0; i < map->n; ++i) {
 767                 map->p[i] = isl_basic_map_drop_inputs(map->p[i], first, n);
 768                 if (!map->p[i])
 769                         goto error;
 770         }
 771         map->n_in -= n;
 772         F_CLR(map, ISL_MAP_NORMALIZED);
 773
 774         return map;
 775 error:
 776         isl_map_free(map);
 777         return NULL;
 778 }
 779
 780 /*
 781  * We don't cow, as the div is assumed to be redundant.
 782  */
 783 static struct isl_basic_map *isl_basic_map_drop_div(
 784                 struct isl_basic_map *bmap, unsigned div)
 785 {
 786         int i;
 787         unsigned pos;
 788
 789         if (!bmap)
 790                 goto error;
 791
 792         pos = 1 + bmap->nparam + bmap->n_in + bmap->n_out + div;
 793
 794         isl_assert(bmap->ctx, div < bmap->n_div, goto error);
 795
 796         for (i = 0; i < bmap->n_eq; ++i)
 797                 constraint_drop_vars(bmap->eq[i]+pos, 1, bmap->extra-div-1);
 798
 799         for (i = 0; i < bmap->n_ineq; ++i) {
 800                 if (!isl_int_is_zero(bmap->ineq[i][pos])) {
 801                         isl_basic_map_drop_inequality(bmap, i);
 802                         --i;
 803                         continue;
 804                 }
 805                 constraint_drop_vars(bmap->ineq[i]+pos, 1, bmap->extra-div-1);
 806         }
 807
 808         for (i = 0; i < bmap->n_div; ++i)
 809                 constraint_drop_vars(bmap->div[i]+1+pos, 1, bmap->extra-div-1);
 810
 811         if (div != bmap->n_div - 1) {
 812                 int j;
 813                 isl_int *t = bmap->div[div];
 814
 815                 for (j = div; j < bmap->n_div - 1; ++j)
 816                         bmap->div[j] = bmap->div[j+1];
 817
 818                 bmap->div[bmap->n_div - 1] = t;
 819         }
 820         F_CLR(bmap, ISL_BASIC_MAP_NORMALIZED);
 821         isl_basic_map_free_div(bmap, 1);
 822
 823         return bmap;
 824 error:
 825         isl_basic_map_free(bmap);
 826         return NULL;
 827 }
 828
 829 static void eliminate_div(struct isl_basic_map *bmap, isl_int *eq, unsigned div)
 830 {
 831         int i;
 832         unsigned pos = 1 + bmap->nparam + bmap->n_in + bmap->n_out + div;
 833         unsigned len;
 834         len = 1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
 835
 836         for (i = 0; i < bmap->n_eq; ++i)
 837                 if (bmap->eq[i] != eq)
 838                         isl_seq_elim(bmap->eq[i], eq, pos, len, NULL);
 839
 840         for (i = 0; i < bmap->n_ineq; ++i)
 841                 isl_seq_elim(bmap->ineq[i], eq, pos, len, NULL);
 842
 843         /* We need to be careful about circular definitions,
 844          * so for now we just remove the definitions of other divs that
 845          * depend on this div and (possibly) recompute them later.
 846          */
 847         for (i = 0; i < bmap->n_div; ++i)
 848                 if (!isl_int_is_zero(bmap->div[i][0]) &&
 849                     !isl_int_is_zero(bmap->div[i][1 + pos]))
 850                         isl_seq_clr(bmap->div[i], 1 + len);
 851
 852         isl_basic_map_drop_div(bmap, div);
 853 }
 854
 855 struct isl_basic_map *isl_basic_map_set_to_empty(struct isl_basic_map *bmap)
 856 {
 857         int i = 0;
 858         unsigned total;
 859         if (!bmap)
 860                 goto error;
 861         total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
 862         isl_basic_map_free_div(bmap, bmap->n_div);
 863         isl_basic_map_free_inequality(bmap, bmap->n_ineq);
 864         if (bmap->n_eq > 0)
 865                 isl_basic_map_free_equality(bmap, bmap->n_eq-1);
 866         else {
 867                 isl_basic_map_alloc_equality(bmap);
 868                 if (i < 0)
 869                         goto error;
 870         }
 871         isl_int_set_si(bmap->eq[i][0], 1);
 872         isl_seq_clr(bmap->eq[i]+1, total);
 873         F_SET(bmap, ISL_BASIC_MAP_EMPTY);
 874         return isl_basic_map_finalize(bmap);
 875 error:
 876         isl_basic_map_free(bmap);
 877         return NULL;
 878 }
 879
 880 struct isl_basic_set *isl_basic_set_set_to_empty(struct isl_basic_set *bset)
 881 {
 882         return (struct isl_basic_set *)
 883                 isl_basic_map_set_to_empty((struct isl_basic_map *)bset);
 884 }
 885
 886 static void swap_equality(struct isl_basic_map *bmap, int a, int b)
 887 {
 888         isl_int *t = bmap->eq[a];
 889         bmap->eq[a] = bmap->eq[b];
 890         bmap->eq[b] = t;
 891 }
 892
 893 static void swap_inequality(struct isl_basic_map *bmap, int a, int b)
 894 {
 895         if (a != b) {
 896                 isl_int *t = bmap->ineq[a];
 897                 bmap->ineq[a] = bmap->ineq[b];
 898                 bmap->ineq[b] = t;
 899         }
 900 }
 901
 902 static void set_swap_inequality(struct isl_basic_set *bset, int a, int b)
 903 {
 904         swap_inequality((struct isl_basic_map *)bset, a, b);
 905 }
 906
 907 static void swap_div(struct isl_basic_map *bmap, int a, int b)
 908 {
 909         int i;
 910         unsigned off = bmap->nparam + bmap->n_in + bmap->n_out;
 911         isl_int *t = bmap->div[a];
 912         bmap->div[a] = bmap->div[b];
 913         bmap->div[b] = t;
 914
 915         for (i = 0; i < bmap->n_eq; ++i)
 916                 isl_int_swap(bmap->eq[i][1+off+a], bmap->eq[i][1+off+b]);
 917
 918         for (i = 0; i < bmap->n_ineq; ++i)
 919                 isl_int_swap(bmap->ineq[i][1+off+a], bmap->ineq[i][1+off+b]);
 920
 921         for (i = 0; i < bmap->n_div; ++i)
 922                 isl_int_swap(bmap->div[i][1+1+off+a], bmap->div[i][1+1+off+b]);
 923         F_CLR(bmap, ISL_BASIC_MAP_NORMALIZED);
 924 }
 925
 926 static void eliminate_var_using_equality(struct isl_basic_map *bmap,
 927         unsigned pos, isl_int *eq, int *progress)
 928 {
 929         unsigned total;
 930         int k;
 931
 932         total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
 933         for (k = 0; k < bmap->n_eq; ++k) {
 934                 if (bmap->eq[k] == eq)
 935                         continue;
 936                 if (isl_int_is_zero(bmap->eq[k][1+pos]))
 937                         continue;
 938                 if (progress)
 939                         *progress = 1;
 940                 isl_seq_elim(bmap->eq[k], eq, 1+pos, 1+total, NULL);
 941         }
 942
 943         for (k = 0; k < bmap->n_ineq; ++k) {
 944                 if (isl_int_is_zero(bmap->ineq[k][1+pos]))
 945                         continue;
 946                 if (progress)
 947                         *progress = 1;
 948                 isl_seq_elim(bmap->ineq[k], eq, 1+pos, 1+total, NULL);
 949                 F_CLR(bmap, ISL_BASIC_MAP_NORMALIZED);
 950         }
 951
 952         for (k = 0; k < bmap->n_div; ++k) {
 953                 if (isl_int_is_zero(bmap->div[k][0]))
 954                         continue;
 955                 if (isl_int_is_zero(bmap->div[k][1+1+pos]))
 956                         continue;
 957                 if (progress)
 958                         *progress = 1;
 959                 isl_seq_elim(bmap->div[k]+1, eq,
 960                                 1+pos, 1+total, &bmap->div[k][0]);
 961                 F_CLR(bmap, ISL_BASIC_MAP_NORMALIZED);
 962         }
 963 }
 964
 965
 966 struct isl_basic_map *isl_basic_map_gauss(
 967         struct isl_basic_map *bmap, int *progress)
 968 {
 969         int k;
 970         int done;
 971         int last_var;
 972         unsigned total_var;
 973         unsigned total;
 974
 975         if (!bmap)
 976                 return NULL;
 977
 978         total_var = bmap->nparam + bmap->n_in + bmap->n_out;
 979         total = total_var + bmap->n_div;
 980
 981         last_var = total - 1;
 982         for (done = 0; done < bmap->n_eq; ++done) {
 983                 for (; last_var >= 0; --last_var) {
 984                         for (k = done; k < bmap->n_eq; ++k)
 985                                 if (!isl_int_is_zero(bmap->eq[k][1+last_var]))
 986                                         break;
 987                         if (k < bmap->n_eq)
 988                                 break;
 989                 }
 990                 if (last_var < 0)
 991                         break;
 992                 if (k != done)
 993                         swap_equality(bmap, k, done);
 994                 if (isl_int_is_neg(bmap->eq[done][1+last_var]))
 995                         isl_seq_neg(bmap->eq[done], bmap->eq[done], 1+total);
 996
 997                 eliminate_var_using_equality(bmap, last_var, bmap->eq[done],
 998                                                 progress);
 999
1000                 if (last_var >= total_var &&
1001                     isl_int_is_zero(bmap->div[last_var - total_var][0])) {
1002                         unsigned div = last_var - total_var;
1003                         isl_seq_neg(bmap->div[div]+1, bmap->eq[done], 1+total);
1004                         isl_int_set_si(bmap->div[div][1+1+last_var], 0);
1005                         isl_int_set(bmap->div[div][0],
1006                                     bmap->eq[done][1+last_var]);
1007                         F_CLR(bmap, ISL_BASIC_MAP_NORMALIZED);
1008                 }
1009         }
1010         if (done == bmap->n_eq)
1011                 return bmap;
1012         for (k = done; k < bmap->n_eq; ++k) {
1013                 if (isl_int_is_zero(bmap->eq[k][0]))
1014                         continue;
1015                 return isl_basic_map_set_to_empty(bmap);
1016         }
1017         isl_basic_map_free_equality(bmap, bmap->n_eq-done);
1018         return bmap;
1019 }
1020
1021 struct isl_basic_set *isl_basic_set_gauss(
1022         struct isl_basic_set *bset, int *progress)
1023 {
1024         return (struct isl_basic_set*)isl_basic_map_gauss(
1025                         (struct isl_basic_map *)bset, progress);
1026 }
1027
1028 static unsigned int round_up(unsigned int v)
1029 {
1030         int old_v = v;
1031
1032         while (v) {
1033                 old_v = v;
1034                 v ^= v & -v;
1035         }
1036         return old_v << 1;
1037 }
1038
1039 static int hash_index(isl_int ***index, unsigned int size, int bits,
1040                         struct isl_basic_map *bmap, int k)
1041 {
1042         int h;
1043         unsigned total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
1044         u_int32_t hash = isl_seq_hash(bmap->ineq[k]+1, total, bits);
1045         for (h = hash; index[h]; h = (h+1) % size)
1046                 if (&bmap->ineq[k] != index[h] &&
1047                     isl_seq_eq(bmap->ineq[k]+1, index[h][0]+1, total))
1048                         break;
1049         return h;
1050 }
1051
1052 static int set_hash_index(isl_int ***index, unsigned int size, int bits,
1053                           struct isl_basic_set *bset, int k)
1054 {
1055         return hash_index(index, size, bits, (struct isl_basic_map *)bset, k);
1056 }
1057
1058 static struct isl_basic_map *remove_duplicate_constraints(
1059         struct isl_basic_map *bmap, int *progress)
1060 {
1061         unsigned int size;
1062         isl_int ***index;
1063         int k, l, h;
1064         int bits;
1065         unsigned total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
1066         isl_int sum;
1067
1068         if (bmap->n_ineq <= 1)
1069                 return bmap;
1070
1071         size = round_up(4 * (bmap->n_ineq+1) / 3 - 1);
1072         bits = ffs(size) - 1;
1073         index = isl_alloc_array(ctx, isl_int **, size);
1074         memset(index, 0, size * sizeof(isl_int **));
1075         if (!index)
1076                 return bmap;
1077
1078         index[isl_seq_hash(bmap->ineq[0]+1, total, bits)] = &bmap->ineq[0];
1079         for (k = 1; k < bmap->n_ineq; ++k) {
1080                 h = hash_index(index, size, bits, bmap, k);
1081                 if (!index[h]) {
1082                         index[h] = &bmap->ineq[k];
1083                         continue;
1084                 }
1085                 if (progress)
1086                         *progress = 1;
1087                 l = index[h] - &bmap->ineq[0];
1088                 if (isl_int_lt(bmap->ineq[k][0], bmap->ineq[l][0]))
1089                         swap_inequality(bmap, k, l);
1090                 isl_basic_map_drop_inequality(bmap, k);
1091                 --k;
1092         }
1093         isl_int_init(sum);
1094         for (k = 0; k < bmap->n_ineq-1; ++k) {
1095                 isl_seq_neg(bmap->ineq[k]+1, bmap->ineq[k]+1, total);
1096                 h = hash_index(index, size, bits, bmap, k);
1097                 isl_seq_neg(bmap->ineq[k]+1, bmap->ineq[k]+1, total);
1098                 if (!index[h])
1099                         continue;
1100                 l = index[h] - &bmap->ineq[0];
1101                 isl_int_add(sum, bmap->ineq[k][0], bmap->ineq[l][0]);
1102                 if (isl_int_is_pos(sum))
1103                         continue;
1104                 if (isl_int_is_zero(sum)) {
1105                         /* We need to break out of the loop after these
1106                          * changes since the contents of the hash
1107                          * will no longer be valid.
1108                          * Plus, we probably we want to regauss first.
1109                          */
1110                         isl_basic_map_drop_inequality(bmap, l);
1111                         isl_basic_map_inequality_to_equality(bmap, k);
1112                 } else
1113                         bmap = isl_basic_map_set_to_empty(bmap);
1114                 break;
1115         }
1116         isl_int_clear(sum);
1117
1118         free(index);
1119         return bmap;
1120 }
1121
1122 static void compute_elimination_index(struct isl_basic_map *bmap, int *elim)
1123 {
1124         int d, i;
1125         unsigned total;
1126
1127         total = bmap->nparam + bmap->n_in + bmap->n_out;
1128         for (d = 0; d < total; ++d)
1129                 elim[d] = -1;
1130         for (d = total - 1, i = 0; d >= 0 && i < bmap->n_eq; ++i) {
1131                 for (; d >= 0; --d) {
1132                         if (isl_int_is_zero(bmap->eq[i][1+d]))
1133                                 continue;
1134                         elim[d] = i;
1135                         break;
1136                 }
1137         }
1138 }
1139
1140 static void set_compute_elimination_index(struct isl_basic_set *bset, int *elim)
1141 {
1142         return compute_elimination_index((struct isl_basic_map *)bset, elim);
1143 }
1144
1145 static int reduced_using_equalities(isl_int *dst, isl_int *src,
1146         struct isl_basic_map *bmap, int *elim)
1147 {
1148         int d, i;
1149         int copied = 0;
1150         unsigned total;
1151
1152         total = bmap->nparam + bmap->n_in + bmap->n_out;
1153         for (d = total - 1; d >= 0; --d) {
1154                 if (isl_int_is_zero(src[1+d]))
1155                         continue;
1156                 if (elim[d] == -1)
1157                         continue;
1158                 if (!copied) {
1159                         isl_seq_cpy(dst, src, 1 + total);
1160                         copied = 1;
1161                 }
1162                 isl_seq_elim(dst, bmap->eq[elim[d]], 1 + d, 1 + total, NULL);
1163         }
1164         return copied;
1165 }
1166
1167 static int set_reduced_using_equalities(isl_int *dst, isl_int *src,
1168         struct isl_basic_set *bset, int *elim)
1169 {
1170         return reduced_using_equalities(dst, src,
1171                                         (struct isl_basic_map *)bset, elim);
1172 }
1173
1174 /* Quick check to see if two basic maps are disjoint.
1175  * In particular, we reduce the equalities and inequalities of
1176  * one basic map in the context of the equalities of the other
1177  * basic map and check if we get a contradiction.
1178  */
1179 int isl_basic_map_fast_is_disjoint(struct isl_basic_map *bmap1,
1180         struct isl_basic_map *bmap2)
1181 {
1182         struct isl_vec *v = NULL;
1183         int *elim = NULL;
1184         unsigned total;
1185         int d, i;
1186
1187         if (!bmap1 || !bmap2)
1188                 return -1;
1189         isl_assert(bmap1->ctx, bmap1->nparam == bmap2->nparam, return -1);
1190         isl_assert(bmap1->ctx, bmap1->n_in == bmap2->n_in, return -1);
1191         isl_assert(bmap1->ctx, bmap1->n_out == bmap2->n_out, return -1);
1192         if (bmap1->n_div || bmap2->n_div)
1193                 return 0;
1194         if (!bmap1->n_eq && !bmap2->n_eq)
1195                 return 0;
1196
1197         total = bmap1->nparam + bmap1->n_in + bmap1->n_out;
1198         if (total == 0)
1199                 return 0;
1200         v = isl_vec_alloc(bmap1->ctx, 1 + total);
1201         if (!v)
1202                 goto error;
1203         elim = isl_alloc_array(bmap1->ctx, int, total);
1204         if (!elim)
1205                 goto error;
1206         compute_elimination_index(bmap1, elim);
1207         for (i = 0; i < bmap2->n_eq; ++i) {
1208                 int reduced;
1209                 reduced = reduced_using_equalities(v->block.data, bmap2->eq[i],
1210                                                         bmap1, elim);
1211                 if (reduced && !isl_int_is_zero(v->block.data[0]) &&
1212                     isl_seq_first_non_zero(v->block.data + 1, total) == -1)
1213                         goto disjoint;
1214         }
1215         for (i = 0; i < bmap2->n_ineq; ++i) {
1216                 int reduced;
1217                 reduced = reduced_using_equalities(v->block.data,
1218                                                 bmap2->ineq[i], bmap1, elim);
1219                 if (reduced && isl_int_is_neg(v->block.data[0]) &&
1220                     isl_seq_first_non_zero(v->block.data + 1, total) == -1)
1221                         goto disjoint;
1222         }
1223         compute_elimination_index(bmap2, elim);
1224         for (i = 0; i < bmap1->n_ineq; ++i) {
1225                 int reduced;
1226                 reduced = reduced_using_equalities(v->block.data,
1227                                                 bmap1->ineq[i], bmap2, elim);
1228                 if (reduced && isl_int_is_neg(v->block.data[0]) &&
1229                     isl_seq_first_non_zero(v->block.data + 1, total) == -1)
1230                         goto disjoint;
1231         }
1232         isl_vec_free(bmap1->ctx, v);
1233         free(elim);
1234         return 0;
1235 disjoint:
1236         isl_vec_free(bmap1->ctx, v);
1237         free(elim);
1238         return 1;
1239 error:
1240         isl_vec_free(bmap1->ctx, v);
1241         free(elim);
1242         return -1;
1243 }
1244
1245 int isl_basic_set_fast_is_disjoint(struct isl_basic_set *bset1,
1246         struct isl_basic_set *bset2)
1247 {
1248         return isl_basic_map_fast_is_disjoint((struct isl_basic_map *)bset1,
1249                                               (struct isl_basic_map *)bset2);
1250 }
1251
1252 int isl_map_fast_is_disjoint(struct isl_map *map1, struct isl_map *map2)
1253 {
1254         int i, j;
1255
1256         if (!map1 || !map2)
1257                 return -1;
1258
1259         if (isl_map_fast_is_equal(map1, map2))
1260                 return 0;
1261
1262         for (i = 0; i < map1->n; ++i) {
1263                 for (j = 0; j < map2->n; ++j) {
1264                         int d = isl_basic_map_fast_is_disjoint(map1->p[i],
1265                                                                map2->p[j]);
1266                         if (d != 1)
1267                                 return d;
1268                 }
1269         }
1270         return 1;
1271 }
1272
1273 int isl_set_fast_is_disjoint(struct isl_set *set1, struct isl_set *set2)
1274 {
1275         return isl_map_fast_is_disjoint((struct isl_map *)set1,
1276                                         (struct isl_map *)set2);
1277 }
1278
1279 static struct isl_basic_map *remove_duplicate_divs(
1280         struct isl_basic_map *bmap, int *progress)
1281 {
1282         unsigned int size;
1283         int *index;
1284         int k, l, h;
1285         int bits;
1286         struct isl_blk eq;
1287         unsigned total_var = bmap->nparam + bmap->n_in + bmap->n_out;
1288         unsigned total = total_var + bmap->n_div;
1289         struct isl_ctx *ctx;
1290
1291         if (bmap->n_div <= 1)
1292                 return bmap;
1293
1294         ctx = bmap->ctx;
1295         for (k = bmap->n_div - 1; k >= 0; --k)
1296                 if (!isl_int_is_zero(bmap->div[k][0]))
1297                         break;
1298         if (k <= 0)
1299                 return bmap;
1300
1301         size = round_up(4 * bmap->n_div / 3 - 1);
1302         bits = ffs(size) - 1;
1303         index = isl_alloc_array(ctx, int, size);
1304         memset(index, 0, size * sizeof(int));
1305         if (!index)
1306                 return bmap;
1307         eq = isl_blk_alloc(ctx, 1+total);
1308         if (isl_blk_is_error(eq))
1309                 goto out;
1310
1311         isl_seq_clr(eq.data, 1+total);
1312         index[isl_seq_hash(bmap->div[k], 2+total, bits)] = k + 1;
1313         for (--k; k >= 0; --k) {
1314                 u_int32_t hash;
1315
1316                 if (isl_int_is_zero(bmap->div[k][0]))
1317                         continue;
1318
1319                 hash = isl_seq_hash(bmap->div[k], 2+total, bits);
1320                 for (h = hash; index[h]; h = (h+1) % size)
1321                         if (isl_seq_eq(bmap->div[k],
1322                                        bmap->div[index[h]-1], 2+total))
1323                                 break;
1324                 if (index[h]) {
1325                         *progress = 1;
1326                         l = index[h] - 1;
1327                         isl_int_set_si(eq.data[1+total_var+k], -1);
1328                         isl_int_set_si(eq.data[1+total_var+l], 1);
1329                         eliminate_div(bmap, eq.data, l);
1330                         isl_int_set_si(eq.data[1+total_var+k], 0);
1331                         isl_int_set_si(eq.data[1+total_var+l], 0);
1332                 }
1333                 index[h] = k+1;
1334         }
1335
1336         isl_blk_free(ctx, eq);
1337 out:
1338         free(index);
1339         return bmap;
1340 }
1341
1342 /* Elimininate divs based on equalities
1343  */
1344 static struct isl_basic_map *eliminate_divs_eq(
1345                 struct isl_basic_map *bmap, int *progress)
1346 {
1347         int d;
1348         int i;
1349         int modified = 0;
1350         unsigned off;
1351
1352         if (!bmap)
1353                 return NULL;
1354
1355         off = 1 + bmap->nparam + bmap->n_in + bmap->n_out;
1356
1357         for (d = bmap->n_div - 1; d >= 0 ; --d) {
1358                 for (i = 0; i < bmap->n_eq; ++i) {
1359                         if (!isl_int_is_one(bmap->eq[i][off + d]) &&
1360                             !isl_int_is_negone(bmap->eq[i][off + d]))
1361                                 continue;
1362                         modified = 1;
1363                         *progress = 1;
1364                         eliminate_div(bmap, bmap->eq[i], d);
1365                         isl_basic_map_drop_equality(bmap, i);
1366                         break;
1367                 }
1368         }
1369         if (modified)
1370                 return eliminate_divs_eq(bmap, progress);
1371         return bmap;
1372 }
1373
1374 static struct isl_basic_map *normalize_constraints(struct isl_basic_map *bmap)
1375 {
1376         int i;
1377         isl_int gcd;
1378         unsigned total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
1379
1380         isl_int_init(gcd);
1381         for (i = bmap->n_eq - 1; i >= 0; --i) {
1382                 isl_seq_gcd(bmap->eq[i]+1, total, &gcd);
1383                 if (isl_int_is_zero(gcd)) {
1384                         if (!isl_int_is_zero(bmap->eq[i][0])) {
1385                                 bmap = isl_basic_map_set_to_empty(bmap);
1386                                 break;
1387                         }
1388                         isl_basic_map_drop_equality(bmap, i);
1389                         continue;
1390                 }
1391                 if (F_ISSET(bmap, ISL_BASIC_MAP_RATIONAL))
1392                         isl_int_gcd(gcd, gcd, bmap->eq[i][0]);
1393                 if (isl_int_is_one(gcd))
1394                         continue;
1395                 if (!isl_int_is_divisible_by(bmap->eq[i][0], gcd)) {
1396                         bmap = isl_basic_map_set_to_empty(bmap);
1397                         break;
1398                 }
1399                 isl_seq_scale_down(bmap->eq[i], bmap->eq[i], gcd, 1+total);
1400         }
1401
1402         for (i = bmap->n_ineq - 1; i >= 0; --i) {
1403                 isl_seq_gcd(bmap->ineq[i]+1, total, &gcd);
1404                 if (isl_int_is_zero(gcd)) {
1405                         if (isl_int_is_neg(bmap->ineq[i][0])) {
1406                                 bmap = isl_basic_map_set_to_empty(bmap);
1407                                 break;
1408                         }
1409                         isl_basic_map_drop_inequality(bmap, i);
1410                         continue;
1411                 }
1412                 if (F_ISSET(bmap, ISL_BASIC_MAP_RATIONAL))
1413                         isl_int_gcd(gcd, gcd, bmap->ineq[i][0]);
1414                 if (isl_int_is_one(gcd))
1415                         continue;
1416                 isl_int_fdiv_q(bmap->ineq[i][0], bmap->ineq[i][0], gcd);
1417                 isl_seq_scale_down(bmap->ineq[i]+1, bmap->ineq[i]+1, gcd, total);
1418         }
1419         isl_int_clear(gcd);
1420
1421         return bmap;
1422 }
1423
1424
1425 /* Eliminate the specified variables from the constraints using
1426  * Fourier-Motzkin.  The variables themselves are not removed.
1427  */
1428 struct isl_basic_map *isl_basic_map_eliminate_vars(
1429         struct isl_basic_map *bmap, int pos, unsigned n)
1430 {
1431         int d;
1432         int i, j, k;
1433         unsigned total;
1434
1435         if (!bmap)
1436                 return NULL;
1437         total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
1438
1439         bmap = isl_basic_map_cow(bmap);
1440         for (d = pos + n - 1; d >= pos; --d) {
1441                 int n_lower, n_upper;
1442                 if (!bmap)
1443                         return NULL;
1444                 for (i = 0; i < bmap->n_eq; ++i) {
1445                         if (isl_int_is_zero(bmap->eq[i][1+d]))
1446                                 continue;
1447                         eliminate_var_using_equality(bmap, d, bmap->eq[i], NULL);
1448                         isl_basic_map_drop_equality(bmap, i);
1449                         break;
1450                 }
1451                 if (i < bmap->n_eq)
1452                         continue;
1453                 n_lower = 0;
1454                 n_upper = 0;
1455                 for (i = 0; i < bmap->n_ineq; ++i) {
1456                         if (isl_int_is_pos(bmap->ineq[i][1+d]))
1457                                 n_lower++;
1458                         else if (isl_int_is_neg(bmap->ineq[i][1+d]))
1459                                 n_upper++;
1460                 }
1461                 bmap = isl_basic_map_extend(bmap,
1462                                 bmap->nparam, bmap->n_in, bmap->n_out, 0,
1463                                 0, n_lower * n_upper);
1464                 for (i = bmap->n_ineq - 1; i >= 0; --i) {
1465                         int last;
1466                         if (isl_int_is_zero(bmap->ineq[i][1+d]))
1467                                 continue;
1468                         last = -1;
1469                         for (j = 0; j < i; ++j) {
1470                                 if (isl_int_is_zero(bmap->ineq[j][1+d]))
1471                                         continue;
1472                                 last = j;
1473                                 if (isl_int_sgn(bmap->ineq[i][1+d]) ==
1474                                     isl_int_sgn(bmap->ineq[j][1+d]))
1475                                         continue;
1476                                 k = isl_basic_map_alloc_inequality(bmap);
1477                                 if (k < 0)
1478                                         goto error;
1479                                 isl_seq_cpy(bmap->ineq[k], bmap->ineq[i],
1480                                                 1+total);
1481                                 isl_seq_elim(bmap->ineq[k], bmap->ineq[j],
1482                                                 1+d, 1+total, NULL);
1483                         }
1484                         isl_basic_map_drop_inequality(bmap, i);
1485                         i = last + 1;
1486                 }
1487                 if (n_lower > 0 && n_upper > 0) {
1488                         bmap = normalize_constraints(bmap);
1489                         bmap = remove_duplicate_constraints(bmap, NULL);
1490                         bmap = isl_basic_map_gauss(bmap, NULL);
1491                         bmap = isl_basic_map_convex_hull(bmap);
1492                         if (!bmap)
1493                                 goto error;
1494                         if (F_ISSET(bmap, ISL_BASIC_MAP_EMPTY))
1495                                 break;
1496                 }
1497         }
1498         F_CLR(bmap, ISL_BASIC_MAP_NORMALIZED);
1499         return bmap;
1500 error:
1501         isl_basic_map_free(bmap);
1502         return NULL;
1503 }
1504
1505 struct isl_basic_set *isl_basic_set_eliminate_vars(
1506         struct isl_basic_set *bset, unsigned pos, unsigned n)
1507 {
1508         return (struct isl_basic_set *)isl_basic_map_eliminate_vars(
1509                         (struct isl_basic_map *)bset, pos, n);
1510 }
1511
1512 /* Eliminate the specified n dimensions starting at first from the
1513  * constraints using Fourier-Motzkin, The dimensions themselves
1514  * are not removed.
1515  */
1516 struct isl_set *isl_set_eliminate_dims(struct isl_set *set,
1517         unsigned first, unsigned n)
1518 {
1519         int i;
1520
1521         if (n == 0)
1522                 return set;
1523
1524         set = isl_set_cow(set);
1525         if (!set)
1526                 return NULL;
1527         isl_assert(set->ctx, first+n <= set->dim, goto error);
1528
1529         for (i = 0; i < set->n; ++i) {
1530                 set->p[i] = isl_basic_set_eliminate_vars(set->p[i],
1531                                                         set->nparam + first, n);
1532                 if (!set->p[i])
1533                         goto error;
1534         }
1535         return set;
1536 error:
1537         isl_set_free(set);
1538         return NULL;
1539 }
1540
1541 /* Project out n dimensions starting at first using Fourier-Motzkin */
1542 struct isl_set *isl_set_remove_dims(struct isl_set *set,
1543         unsigned first, unsigned n)
1544 {
1545         set = isl_set_eliminate_dims(set, first, n);
1546         set = isl_set_drop_dims(set, first, n);
1547         return set;
1548 }
1549
1550 /* Project out n inputs starting at first using Fourier-Motzkin */
1551 struct isl_map *isl_map_remove_inputs(struct isl_map *map,
1552         unsigned first, unsigned n)
1553 {
1554         int i;
1555
1556         if (n == 0)
1557                 return map;
1558
1559         map = isl_map_cow(map);
1560         if (!map)
1561                 return NULL;
1562         isl_assert(map->ctx, first+n <= map->n_in, goto error);
1563
1564         for (i = 0; i < map->n; ++i) {
1565                 map->p[i] = isl_basic_map_eliminate_vars(map->p[i],
1566                                                         map->nparam + first, n);
1567                 if (!map->p[i])
1568                         goto error;
1569         }
1570         map = isl_map_drop_inputs(map, first, n);
1571         return map;
1572 error:
1573         isl_map_free(map);
1574         return NULL;
1575 }
1576
1577 /* Project out n dimensions starting at first using Fourier-Motzkin */
1578 struct isl_basic_set *isl_basic_set_remove_dims(struct isl_basic_set *bset,
1579         unsigned first, unsigned n)
1580 {
1581         bset = isl_basic_set_eliminate_vars(bset, bset->nparam + first, n);
1582         bset = isl_basic_set_drop_dims(bset, first, n);
1583         return bset;
1584 }
1585
1586 /* Elimininate divs based on inequalities
1587  */
1588 static struct isl_basic_map *eliminate_divs_ineq(
1589                 struct isl_basic_map *bmap, int *progress)
1590 {
1591         int d;
1592         int i;
1593         unsigned off;
1594         struct isl_ctx *ctx;
1595
1596         if (!bmap)
1597                 return NULL;
1598
1599         ctx = bmap->ctx;
1600         off = 1 + bmap->nparam + bmap->n_in + bmap->n_out;
1601
1602         for (d = bmap->n_div - 1; d >= 0 ; --d) {
1603                 for (i = 0; i < bmap->n_eq; ++i)
1604                         if (!isl_int_is_zero(bmap->eq[i][off + d]))
1605                                 break;
1606                 if (i < bmap->n_eq)
1607                         continue;
1608                 for (i = 0; i < bmap->n_ineq; ++i)
1609                         if (isl_int_abs_gt(bmap->ineq[i][off + d], ctx->one))
1610                                 break;
1611                 if (i < bmap->n_ineq)
1612                         continue;
1613                 *progress = 1;
1614                 bmap = isl_basic_map_eliminate_vars(bmap, (off-1)+d, 1);
1615                 if (F_ISSET(bmap, ISL_BASIC_MAP_EMPTY))
1616                         break;
1617                 bmap = isl_basic_map_drop_div(bmap, d);
1618                 if (!bmap)
1619                         break;
1620         }
1621         return bmap;
1622 }
1623
1624 struct isl_basic_map *isl_basic_map_simplify(struct isl_basic_map *bmap)
1625 {
1626         int progress = 1;
1627         if (!bmap)
1628                 return NULL;
1629         while (progress) {
1630                 progress = 0;
1631                 bmap = normalize_constraints(bmap);
1632                 bmap = eliminate_divs_eq(bmap, &progress);
1633                 bmap = eliminate_divs_ineq(bmap, &progress);
1634                 bmap = isl_basic_map_gauss(bmap, &progress);
1635                 bmap = remove_duplicate_divs(bmap, &progress);
1636                 bmap = remove_duplicate_constraints(bmap, &progress);
1637         }
1638         return bmap;
1639 }
1640
1641 struct isl_basic_set *isl_basic_set_simplify(struct isl_basic_set *bset)
1642 {
1643         return (struct isl_basic_set *)
1644                 isl_basic_map_simplify((struct isl_basic_map *)bset);
1645 }
1646
1647 static void dump_term(struct isl_basic_map *bmap,
1648                         isl_int c, int pos, FILE *out)
1649 {
1650         unsigned in = bmap->n_in;
1651         unsigned dim = bmap->n_in + bmap->n_out;
1652         if (!pos)
1653                 isl_int_print(out, c, 0);
1654         else {
1655                 if (!isl_int_is_one(c))
1656                         isl_int_print(out, c, 0);
1657                 if (pos < 1 + bmap->nparam)
1658                         fprintf(out, "p%d", pos - 1);
1659                 else if (pos < 1 + bmap->nparam + in)
1660                         fprintf(out, "i%d", pos - 1 - bmap->nparam);
1661                 else if (pos < 1 + bmap->nparam + dim)
1662                         fprintf(out, "o%d", pos - 1 - bmap->nparam - in);
1663                 else
1664                         fprintf(out, "e%d", pos - 1 - bmap->nparam - dim);
1665         }
1666 }
1667
1668 static void dump_constraint_sign(struct isl_basic_map *bmap, isl_int *c,
1669                                 int sign, FILE *out)
1670 {
1671         int i;
1672         int first;
1673         unsigned len = 1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
1674         isl_int v;
1675
1676         isl_int_init(v);
1677         for (i = 0, first = 1; i < len; ++i) {
1678                 if (isl_int_sgn(c[i]) * sign <= 0)
1679                         continue;
1680                 if (!first)
1681                         fprintf(out, " + ");
1682                 first = 0;
1683                 isl_int_abs(v, c[i]);
1684                 dump_term(bmap, v, i, out);
1685         }
1686         isl_int_clear(v);
1687         if (first)
1688                 fprintf(out, "0");
1689 }
1690
1691 static void dump_constraint(struct isl_basic_map *bmap, isl_int *c,
1692                                 const char *op, FILE *out, int indent)
1693 {
1694         int i;
1695
1696         fprintf(out, "%*s", indent, "");
1697
1698         dump_constraint_sign(bmap, c, 1, out);
1699         fprintf(out, " %s ", op);
1700         dump_constraint_sign(bmap, c, -1, out);
1701
1702         fprintf(out, "\n");
1703
1704         for (i = bmap->n_div; i < bmap->extra; ++i) {
1705                 if (isl_int_is_zero(c[1+bmap->nparam+bmap->n_in+bmap->n_out+i]))
1706                         continue;
1707                 fprintf(out, "%*s", indent, "");
1708                 fprintf(out, "ERROR: unused div coefficient not zero\n");
1709         }
1710 }
1711
1712 static void dump_constraints(struct isl_basic_map *bmap,
1713                                 isl_int **c, unsigned n,
1714                                 const char *op, FILE *out, int indent)
1715 {
1716         int i;
1717
1718         for (i = 0; i < n; ++i)
1719                 dump_constraint(bmap, c[i], op, out, indent);
1720 }
1721
1722 static void dump_affine(struct isl_basic_map *bmap, isl_int *exp, FILE *out)
1723 {
1724         int j;
1725         int first = 1;
1726         unsigned total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
1727
1728         for (j = 0; j < 1 + total; ++j) {
1729                 if (isl_int_is_zero(exp[j]))
1730                         continue;
1731                 if (!first && isl_int_is_pos(exp[j]))
1732                         fprintf(out, "+");
1733                 dump_term(bmap, exp[j], j, out);
1734                 first = 0;
1735         }
1736 }
1737
1738 static void dump(struct isl_basic_map *bmap, FILE *out, int indent)
1739 {
1740         int i;
1741
1742         dump_constraints(bmap, bmap->eq, bmap->n_eq, "=", out, indent);
1743         dump_constraints(bmap, bmap->ineq, bmap->n_ineq, ">=", out, indent);
1744
1745         for (i = 0; i < bmap->n_div; ++i) {
1746                 fprintf(out, "%*s", indent, "");
1747                 fprintf(out, "e%d = [(", i);
1748                 dump_affine(bmap, bmap->div[i]+1, out);
1749                 fprintf(out, ")/");
1750                 isl_int_print(out, bmap->div[i][0], 0);
1751                 fprintf(out, "]\n");
1752         }
1753 }
1754
1755 void isl_basic_set_dump(struct isl_basic_set *bset, FILE *out, int indent)
1756 {
1757         if (!bset) {
1758                 fprintf(out, "null basic set\n");
1759                 return;
1760         }
1761
1762         fprintf(out, "%*s", indent, "");
1763         fprintf(out, "ref: %d, nparam: %d, dim: %d, extra: %d, flags: %x\n",
1764                         bset->ref, bset->nparam, bset->dim, bset->extra,
1765                         bset->flags);
1766         dump((struct isl_basic_map *)bset, out, indent);
1767 }
1768
1769 void isl_basic_map_dump(struct isl_basic_map *bmap, FILE *out, int indent)
1770 {
1771         if (!bmap) {
1772                 fprintf(out, "null basic map\n");
1773                 return;
1774         }
1775
1776         fprintf(out, "%*s", indent, "");
1777         fprintf(out, "ref: %d, nparam: %d, in: %d, out: %d, extra: %d, flags: %x\n",
1778                 bmap->ref,
1779                 bmap->nparam, bmap->n_in, bmap->n_out, bmap->extra, bmap->flags);
1780         dump(bmap, out, indent);
1781 }
1782
1783 int isl_inequality_negate(struct isl_basic_map *bmap, unsigned pos)
1784 {
1785         unsigned total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
1786         if (!bmap)
1787                 return -1;
1788         isl_assert(bmap->ctx, pos < bmap->n_ineq, return -1);
1789         isl_seq_neg(bmap->ineq[pos], bmap->ineq[pos], 1 + total);
1790         isl_int_sub_ui(bmap->ineq[pos][0], bmap->ineq[pos][0], 1);
1791         F_CLR(bmap, ISL_BASIC_MAP_NORMALIZED);
1792         return 0;
1793 }
1794
1795 struct isl_set *isl_set_alloc(struct isl_ctx *ctx,
1796                 unsigned nparam, unsigned dim, int n, unsigned flags)
1797 {
1798         struct isl_set *set;
1799
1800         isl_assert(ctx, n >= 0, return NULL);
1801         set = isl_alloc(ctx, struct isl_set,
1802                         sizeof(struct isl_set) +
1803                         n * sizeof(struct isl_basic_set *));
1804         if (!set)
1805                 return NULL;
1806
1807         set->ctx = ctx;
1808         isl_ctx_ref(ctx);
1809         set->ref = 1;
1810         set->size = n;
1811         set->n = 0;
1812         set->nparam = nparam;
1813         set->zero = 0;
1814         set->dim = dim;
1815         set->flags = flags;
1816         return set;
1817 }
1818
1819 struct isl_set *isl_set_dup(struct isl_set *set)
1820 {
1821         int i;
1822         struct isl_set *dup;
1823
1824         if (!set)
1825                 return NULL;
1826
1827         dup = isl_set_alloc(set->ctx, set->nparam, set->dim, set->n, set->flags);
1828         if (!dup)
1829                 return NULL;
1830         for (i = 0; i < set->n; ++i)
1831                 dup = isl_set_add(dup, isl_basic_set_copy(set->p[i]));
1832         return dup;
1833 }
1834
1835 struct isl_set *isl_set_from_basic_set(struct isl_basic_set *bset)
1836 {
1837         struct isl_set *set;
1838
1839         if (!bset)
1840                 return NULL;
1841
1842         set = isl_set_alloc(bset->ctx, bset->nparam, bset->dim, 1, ISL_MAP_DISJOINT);
1843         if (!set) {
1844                 isl_basic_set_free(bset);
1845                 return NULL;
1846         }
1847         return isl_set_add(set, bset);
1848 }
1849
1850 struct isl_map *isl_map_from_basic_map(struct isl_basic_map *bmap)
1851 {
1852         struct isl_map *map;
1853
1854         if (!bmap)
1855                 return NULL;
1856
1857         map = isl_map_alloc(bmap->ctx, bmap->nparam, bmap->n_in, bmap->n_out, 1,
1858                                 ISL_MAP_DISJOINT);
1859         if (!map) {
1860                 isl_basic_map_free(bmap);
1861                 return NULL;
1862         }
1863         return isl_map_add(map, bmap);
1864 }
1865
1866 struct isl_set *isl_set_add(struct isl_set *set, struct isl_basic_set *bset)
1867 {
1868         if (!bset || !set)
1869                 goto error;
1870         isl_assert(set->ctx, set->nparam == bset->nparam, goto error);
1871         isl_assert(set->ctx, set->dim == bset->dim, goto error);
1872         isl_assert(set->ctx, set->n < set->size, goto error);
1873         set->p[set->n] = bset;
1874         set->n++;
1875         return set;
1876 error:
1877         if (set)
1878                 isl_set_free(set);
1879         if (bset)
1880                 isl_basic_set_free(bset);
1881         return NULL;
1882 }
1883
1884 void isl_set_free(struct isl_set *set)
1885 {
1886         int i;
1887
1888         if (!set)
1889                 return;
1890
1891         if (--set->ref > 0)
1892                 return;
1893
1894         isl_ctx_deref(set->ctx);
1895         for (i = 0; i < set->n; ++i)
1896                 isl_basic_set_free(set->p[i]);
1897         free(set);
1898 }
1899
1900 void isl_set_dump(struct isl_set *set, FILE *out, int indent)
1901 {
1902         int i;
1903
1904         if (!set) {
1905                 fprintf(out, "null set\n");
1906                 return;
1907         }
1908
1909         fprintf(out, "%*s", indent, "");
1910         fprintf(out, "ref: %d, n: %d, nparam: %d, dim: %d, flags: %x\n",
1911                         set->ref, set->n, set->nparam, set->dim, set->flags);
1912         for (i = 0; i < set->n; ++i) {
1913                 fprintf(out, "%*s", indent, "");
1914                 fprintf(out, "basic set %d:\n", i);
1915                 isl_basic_set_dump(set->p[i], out, indent+4);
1916         }
1917 }
1918
1919 void isl_map_dump(struct isl_map *map, FILE *out, int indent)
1920 {
1921         int i;
1922
1923         if (!map) {
1924                 fprintf(out, "null map\n");
1925                 return;
1926         }
1927
1928         fprintf(out, "%*s", indent, "");
1929         fprintf(out, "ref: %d, n: %d, nparam: %d, in: %d, out: %d, flags: %x\n",
1930                         map->ref, map->n, map->nparam, map->n_in, map->n_out,
1931                         map->flags);
1932         for (i = 0; i < map->n; ++i) {
1933                 fprintf(out, "%*s", indent, "");
1934                 fprintf(out, "basic map %d:\n", i);
1935                 isl_basic_map_dump(map->p[i], out, indent+4);
1936         }
1937 }
1938
1939 struct isl_basic_map *isl_basic_map_intersect_domain(
1940                 struct isl_basic_map *bmap, struct isl_basic_set *bset)
1941 {
1942         struct isl_basic_map *bmap_domain;
1943
1944         if (!bmap || !bset)
1945                 goto error;
1946
1947         isl_assert(bset->ctx, bset->dim == bmap->n_in, goto error);
1948         isl_assert(bset->ctx, bset->nparam == bmap->nparam, goto error);
1949
1950         bmap = isl_basic_map_extend(bmap,
1951                         bset->nparam, bmap->n_in, bmap->n_out,
1952                         bset->n_div, bset->n_eq, bset->n_ineq);
1953         if (!bmap)
1954                 goto error;
1955         bmap_domain = isl_basic_map_from_basic_set(bset,
1956                                                 bset->dim, 0);
1957         bmap = add_constraints(bmap, bmap_domain, 0, 0);
1958
1959         bmap = isl_basic_map_simplify(bmap);
1960         return isl_basic_map_finalize(bmap);
1961 error:
1962         isl_basic_map_free(bmap);
1963         isl_basic_set_free(bset);
1964         return NULL;
1965 }
1966
1967 struct isl_basic_map *isl_basic_map_intersect_range(
1968                 struct isl_basic_map *bmap, struct isl_basic_set *bset)
1969 {
1970         struct isl_basic_map *bmap_range;
1971
1972         if (!bmap || !bset)
1973                 goto error;
1974
1975         isl_assert(bset->ctx, bset->dim == bmap->n_out, goto error);
1976         isl_assert(bset->ctx, bset->nparam == bmap->nparam, goto error);
1977
1978         bmap = isl_basic_map_extend(bmap,
1979                         bset->nparam, bmap->n_in, bmap->n_out,
1980                         bset->n_div, bset->n_eq, bset->n_ineq);
1981         if (!bmap)
1982                 goto error;
1983         bmap_range = isl_basic_map_from_basic_set(bset,
1984                                                 0, bset->dim);
1985         bmap = add_constraints(bmap, bmap_range, 0, 0);
1986
1987         bmap = isl_basic_map_simplify(bmap);
1988         return isl_basic_map_finalize(bmap);
1989 error:
1990         isl_basic_map_free(bmap);
1991         isl_basic_set_free(bset);
1992         return NULL;
1993 }
1994
1995 struct isl_basic_map *isl_basic_map_intersect(
1996                 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
1997 {
1998         if (!bmap1 || !bmap2)
1999                 goto error;
2000
2001         isl_assert(bmap1->ctx, bmap1->nparam == bmap2->nparam, goto error);
2002         isl_assert(bmap1->ctx, bmap1->n_in == bmap2->n_in, goto error);
2003         isl_assert(bmap1->ctx, bmap1->n_out == bmap2->n_out, goto error);
2004
2005         bmap1 = isl_basic_map_extend(bmap1,
2006                         bmap1->nparam, bmap1->n_in, bmap1->n_out,
2007                         bmap2->n_div, bmap2->n_eq, bmap2->n_ineq);
2008         if (!bmap1)
2009                 goto error;
2010         bmap1 = add_constraints(bmap1, bmap2, 0, 0);
2011
2012         bmap1 = isl_basic_map_simplify(bmap1);
2013         return isl_basic_map_finalize(bmap1);
2014 error:
2015         isl_basic_map_free(bmap1);
2016         isl_basic_map_free(bmap2);
2017         return NULL;
2018 }
2019
2020 struct isl_basic_set *isl_basic_set_intersect(
2021                 struct isl_basic_set *bset1, struct isl_basic_set *bset2)
2022 {
2023         return (struct isl_basic_set *)
2024                 isl_basic_map_intersect(
2025                         (struct isl_basic_map *)bset1,
2026                         (struct isl_basic_map *)bset2);
2027 }
2028
2029 struct isl_map *isl_map_intersect(struct isl_map *map1, struct isl_map *map2)
2030 {
2031         unsigned flags = 0;
2032         struct isl_map *result;
2033         int i, j;
2034
2035         if (!map1 || !map2)
2036                 goto error;
2037
2038         if (F_ISSET(map1, ISL_MAP_DISJOINT) &&
2039             F_ISSET(map2, ISL_MAP_DISJOINT))
2040                 FL_SET(flags, ISL_MAP_DISJOINT);
2041
2042         result = isl_map_alloc(map1->ctx, map1->nparam, map1->n_in, map1->n_out,
2043                                 map1->n * map2->n, flags);
2044         if (!result)
2045                 goto error;
2046         for (i = 0; i < map1->n; ++i)
2047                 for (j = 0; j < map2->n; ++j) {
2048                         struct isl_basic_map *part;
2049                         part = isl_basic_map_intersect(
2050                                     isl_basic_map_copy(map1->p[i]),
2051                                     isl_basic_map_copy(map2->p[j]));
2052                         if (isl_basic_map_is_empty(part))
2053                                 isl_basic_map_free(part);
2054                         else
2055                                 result = isl_map_add(result, part);
2056                         if (!result)
2057                                 goto error;
2058                 }
2059         isl_map_free(map1);
2060         isl_map_free(map2);
2061         return result;
2062 error:
2063         isl_map_free(map1);
2064         isl_map_free(map2);
2065         return NULL;
2066 }
2067
2068 struct isl_set *isl_set_intersect(struct isl_set *set1, struct isl_set *set2)
2069 {
2070         return (struct isl_set *)
2071                 isl_map_intersect((struct isl_map *)set1,
2072                                   (struct isl_map *)set2);
2073 }
2074
2075 struct isl_basic_map *isl_basic_map_reverse(struct isl_basic_map *bmap)
2076 {
2077         struct isl_basic_set *bset;
2078         unsigned in;
2079
2080         if (!bmap)
2081                 return NULL;
2082         bmap = isl_basic_map_cow(bmap);
2083         if (!bmap)
2084                 return NULL;
2085         in = bmap->n_in;
2086         bset = isl_basic_set_from_basic_map(bmap);
2087         bset = isl_basic_set_swap_vars(bset, in);
2088         return isl_basic_map_from_basic_set(bset, bset->dim-in, in);
2089 }
2090
2091 /* Turn final n dimensions into existentially quantified variables.
2092  */
2093 struct isl_basic_set *isl_basic_set_project_out(
2094                 struct isl_basic_set *bset, unsigned n, unsigned flags)
2095 {
2096         int i;
2097         size_t row_size;
2098         isl_int **new_div;
2099         isl_int *old;
2100
2101         if (!bset)
2102                 return NULL;
2103
2104         isl_assert(bset->ctx, n <= bset->dim, goto error);
2105
2106         if (n == 0)
2107                 return bset;
2108
2109         bset = isl_basic_set_cow(bset);
2110
2111         row_size = 1 + bset->nparam + bset->dim + bset->extra;
2112         old = bset->block2.data;
2113         bset->block2 = isl_blk_extend(bset->ctx, bset->block2,
2114                                         (bset->extra + n) * (1 + row_size));
2115         if (!bset->block2.data)
2116                 goto error;
2117         new_div = isl_alloc_array(ctx, isl_int *, bset->extra + n);
2118         if (!new_div)
2119                 goto error;
2120         for (i = 0; i < n; ++i) {
2121                 new_div[i] = bset->block2.data +
2122                                 (bset->extra + i) * (1 + row_size);
2123                 isl_seq_clr(new_div[i], 1 + row_size);
2124         }
2125         for (i = 0; i < bset->extra; ++i)
2126                 new_div[n + i] = bset->block2.data + (bset->div[i] - old);
2127         free(bset->div);
2128         bset->div = new_div;
2129         bset->n_div += n;
2130         bset->extra += n;
2131         bset->dim -= n;
2132         bset = isl_basic_set_simplify(bset);
2133         return isl_basic_set_finalize(bset);
2134 error:
2135         isl_basic_set_free(bset);
2136         return NULL;
2137 }
2138
2139 struct isl_basic_map *isl_basic_map_apply_range(
2140                 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
2141 {
2142         struct isl_basic_set *bset;
2143         unsigned n_in, n_out;
2144
2145         if (!bmap1 || !bmap2)
2146                 goto error;
2147
2148         isl_assert(bmap->ctx, bmap1->n_out == bmap2->n_in, goto error);
2149         isl_assert(bmap->ctx, bmap1->nparam == bmap2->nparam, goto error);
2150
2151         n_in = bmap1->n_in;
2152         n_out = bmap2->n_out;
2153
2154         bmap2 = isl_basic_map_reverse(bmap2);
2155         if (!bmap2)
2156                 goto error;
2157         bmap1 = isl_basic_map_extend(bmap1, bmap1->nparam,
2158                         bmap1->n_in + bmap2->n_in, bmap1->n_out,
2159                         bmap2->extra, bmap2->n_eq, bmap2->n_ineq);
2160         if (!bmap1)
2161                 goto error;
2162         bmap1 = add_constraints(bmap1, bmap2, bmap1->n_in - bmap2->n_in, 0);
2163         bmap1 = isl_basic_map_simplify(bmap1);
2164         bset = isl_basic_set_from_basic_map(bmap1);
2165         bset = isl_basic_set_project_out(bset,
2166                                                 bset->dim - (n_in + n_out), 0);
2167         return isl_basic_map_from_basic_set(bset, n_in, n_out);
2168 error:
2169         isl_basic_map_free(bmap1);
2170         isl_basic_map_free(bmap2);
2171         return NULL;
2172 }
2173
2174 struct isl_basic_set *isl_basic_set_apply(
2175                 struct isl_basic_set *bset, struct isl_basic_map *bmap)
2176 {
2177         if (!bset || !bmap)
2178                 goto error;
2179
2180         isl_assert(bset->ctx, bset->dim == bmap->n_in, goto error);
2181
2182         return (struct isl_basic_set *)
2183                 isl_basic_map_apply_range((struct isl_basic_map *)bset, bmap);
2184 error:
2185         isl_basic_set_free(bset);
2186         isl_basic_map_free(bmap);
2187         return NULL;
2188 }
2189
2190 struct isl_basic_map *isl_basic_map_apply_domain(
2191                 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
2192 {
2193         if (!bmap1 || !bmap2)
2194                 goto error;
2195
2196         isl_assert(ctx, bmap1->n_in == bmap2->n_in, goto error);
2197         isl_assert(ctx, bmap1->nparam == bmap2->nparam, goto error);
2198
2199         bmap1 = isl_basic_map_reverse(bmap1);
2200         bmap1 = isl_basic_map_apply_range(bmap1, bmap2);
2201         return isl_basic_map_reverse(bmap1);
2202 error:
2203         isl_basic_map_free(bmap1);
2204         isl_basic_map_free(bmap2);
2205         return NULL;
2206 }
2207
2208 static struct isl_basic_map *var_equal(struct isl_ctx *ctx,
2209                 struct isl_basic_map *bmap, unsigned pos)
2210 {
2211         int i;
2212         i = isl_basic_map_alloc_equality(bmap);
2213         if (i < 0)
2214                 goto error;
2215         isl_seq_clr(bmap->eq[i],
2216                     1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div);
2217         isl_int_set_si(bmap->eq[i][1+bmap->nparam+pos], -1);
2218         isl_int_set_si(bmap->eq[i][1+bmap->nparam+bmap->n_in+pos], 1);
2219         return isl_basic_map_finalize(bmap);
2220 error:
2221         isl_basic_map_free(bmap);
2222         return NULL;
2223 }
2224
2225 static struct isl_basic_map *var_more(struct isl_ctx *ctx,
2226                 struct isl_basic_map *bmap, unsigned pos)
2227 {
2228         int i;
2229         i = isl_basic_map_alloc_inequality(bmap);
2230         if (i < 0)
2231                 goto error;
2232         isl_seq_clr(bmap->ineq[i],
2233                     1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div);
2234         isl_int_set_si(bmap->ineq[i][0], -1);
2235         isl_int_set_si(bmap->ineq[i][1+bmap->nparam+pos], -1);
2236         isl_int_set_si(bmap->ineq[i][1+bmap->nparam+bmap->n_in+pos], 1);
2237         return isl_basic_map_finalize(bmap);
2238 error:
2239         isl_basic_map_free(bmap);
2240         return NULL;
2241 }
2242
2243 static struct isl_basic_map *var_less(struct isl_ctx *ctx,
2244                 struct isl_basic_map *bmap, unsigned pos)
2245 {
2246         int i;
2247         i = isl_basic_map_alloc_inequality(bmap);
2248         if (i < 0)
2249                 goto error;
2250         isl_seq_clr(bmap->ineq[i],
2251                     1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div);
2252         isl_int_set_si(bmap->ineq[i][0], -1);
2253         isl_int_set_si(bmap->ineq[i][1+bmap->nparam+pos], 1);
2254         isl_int_set_si(bmap->ineq[i][1+bmap->nparam+bmap->n_in+pos], -1);
2255         return isl_basic_map_finalize(bmap);
2256 error:
2257         isl_basic_map_free(bmap);
2258         return NULL;
2259 }
2260
2261 struct isl_basic_map *isl_basic_map_equal(struct isl_ctx *ctx,
2262                 unsigned nparam, unsigned in, unsigned out, unsigned n_equal)
2263 {
2264         int i;
2265         struct isl_basic_map *bmap;
2266         bmap = isl_basic_map_alloc(ctx, nparam, in, out, 0, n_equal, 0);
2267         if (!bmap)
2268                 return NULL;
2269         for (i = 0; i < n_equal && bmap; ++i)
2270                 bmap = var_equal(ctx, bmap, i);
2271         return isl_basic_map_finalize(bmap);
2272 }
2273
2274 struct isl_basic_map *isl_basic_map_less_at(struct isl_ctx *ctx,
2275                 unsigned nparam, unsigned in, unsigned out, unsigned pos)
2276 {
2277         int i;
2278         struct isl_basic_map *bmap;
2279         bmap = isl_basic_map_alloc(ctx, nparam, in, out, 0, pos, 1);
2280         if (!bmap)
2281                 return NULL;
2282         for (i = 0; i < pos && bmap; ++i)
2283                 bmap = var_equal(ctx, bmap, i);
2284         if (bmap)
2285                 bmap = var_less(ctx, bmap, pos);
2286         return isl_basic_map_finalize(bmap);
2287 }
2288
2289 struct isl_basic_map *isl_basic_map_more_at(struct isl_ctx *ctx,
2290                 unsigned nparam, unsigned in, unsigned out, unsigned pos)
2291 {
2292         int i;
2293         struct isl_basic_map *bmap;
2294         bmap = isl_basic_map_alloc(ctx, nparam, in, out, 0, pos, 1);
2295         if (!bmap)
2296                 return NULL;
2297         for (i = 0; i < pos && bmap; ++i)
2298                 bmap = var_equal(ctx, bmap, i);
2299         if (bmap)
2300                 bmap = var_more(ctx, bmap, pos);
2301         return isl_basic_map_finalize(bmap);
2302 }
2303
2304 struct isl_basic_map *isl_basic_map_from_basic_set(
2305                 struct isl_basic_set *bset, unsigned n_in, unsigned n_out)
2306 {
2307         struct isl_basic_map *bmap;
2308
2309         bset = isl_basic_set_cow(bset);
2310         if (!bset)
2311                 return NULL;
2312
2313         isl_assert(ctx, bset->dim == n_in + n_out, goto error);
2314         bmap = (struct isl_basic_map *) bset;
2315         bmap->n_in = n_in;
2316         bmap->n_out = n_out;
2317         return isl_basic_map_finalize(bmap);
2318 error:
2319         isl_basic_set_free(bset);
2320         return NULL;
2321 }
2322
2323 struct isl_basic_set *isl_basic_set_from_basic_map(struct isl_basic_map *bmap)
2324 {
2325         if (!bmap)
2326                 goto error;
2327         if (bmap->n_in == 0)
2328                 return (struct isl_basic_set *)bmap;
2329         bmap = isl_basic_map_cow(bmap);
2330         if (!bmap)
2331                 goto error;
2332         bmap->n_out += bmap->n_in;
2333         bmap->n_in = 0;
2334         bmap = isl_basic_map_finalize(bmap);
2335         return (struct isl_basic_set *)bmap;
2336 error:
2337         return NULL;
2338 }
2339
2340 /* For a div d = floor(f/m), add the constraints
2341  *
2342  *              f - m d >= 0
2343  *              -(f-(n-1)) + m d >= 0
2344  *
2345  * Note that the second constraint is the negation of
2346  *
2347  *              f - m d >= n
2348  */
2349 static int add_div_constraints(struct isl_basic_map *bmap, unsigned div)
2350 {
2351         int i, j;
2352         unsigned div_pos = 1 + bmap->nparam + bmap->n_in + bmap->n_out + div;
2353         unsigned total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
2354
2355         i = isl_basic_map_alloc_inequality(bmap);
2356         if (i < 0)
2357                 return -1;
2358         isl_seq_cpy(bmap->ineq[i], bmap->div[div]+1, 1+total);
2359         isl_int_neg(bmap->ineq[i][div_pos], bmap->div[div][0]);
2360
2361         j = isl_basic_map_alloc_inequality(bmap);
2362         if (j < 0)
2363                 return -1;
2364         isl_seq_neg(bmap->ineq[j], bmap->ineq[i], 1 + total);
2365         isl_int_add(bmap->ineq[j][0], bmap->ineq[j][0], bmap->ineq[j][div_pos]);
2366         isl_int_sub_ui(bmap->ineq[j][0], bmap->ineq[j][0], 1);
2367         return j;
2368 }
2369
2370 struct isl_basic_set *isl_basic_map_underlying_set(
2371                 struct isl_basic_map *bmap)
2372 {
2373         if (!bmap)
2374                 goto error;
2375         if (bmap->nparam == 0 && bmap->n_in == 0 && bmap->n_div == 0)
2376                 return (struct isl_basic_set *)bmap;
2377         bmap = isl_basic_map_cow(bmap);
2378         if (!bmap)
2379                 goto error;
2380         bmap->n_out += bmap->nparam + bmap->n_in + bmap->n_div;
2381         bmap->nparam = 0;
2382         bmap->n_in = 0;
2383         bmap->extra -= bmap->n_div;
2384         bmap->n_div = 0;
2385         bmap = isl_basic_map_finalize(bmap);
2386         return (struct isl_basic_set *)bmap;
2387 error:
2388         return NULL;
2389 }
2390
2391 struct isl_basic_map *isl_basic_map_overlying_set(
2392         struct isl_basic_set *bset, struct isl_basic_map *like)
2393 {
2394         struct isl_basic_map *bmap;
2395         struct isl_ctx *ctx;
2396         unsigned total;
2397         int i, k;
2398
2399         if (!bset || !like)
2400                 goto error;
2401         ctx = bset->ctx;
2402         isl_assert(ctx, bset->dim ==
2403                         like->nparam + like->n_in + like->n_out + like->n_div,
2404                         goto error);
2405         if (like->nparam == 0 && like->n_in == 0 && like->n_div == 0) {
2406                 isl_basic_map_free(like);
2407                 return (struct isl_basic_map *)bset;
2408         }
2409         bset = isl_basic_set_cow(bset);
2410         if (!bset)
2411                 goto error;
2412         total = bset->dim + bset->extra;
2413         bmap = (struct isl_basic_map *)bset;
2414         bmap->nparam = like->nparam;
2415         bmap->n_in = like->n_in;
2416         bmap->n_out = like->n_out;
2417         bmap->extra += like->n_div;
2418         if (bmap->extra) {
2419                 unsigned ltotal;
2420                 ltotal = total - bmap->extra + like->extra;
2421                 if (ltotal > total)
2422                         ltotal = total;
2423                 bmap->block2 = isl_blk_extend(ctx, bmap->block2,
2424                                         bmap->extra * (1 + 1 + total));
2425                 if (isl_blk_is_error(bmap->block2))
2426                         goto error;
2427                 bmap->div = isl_realloc_array(ctx, bmap->div, isl_int *,
2428                                                 bmap->extra);
2429                 if (!bmap->div)
2430                         goto error;
2431                 bmap = isl_basic_map_extend(bmap, bmap->nparam,
2432                         bmap->n_in, bmap->n_out, 0, 0, 2 * like->n_div);
2433                 for (i = 0; i < like->n_div; ++i) {
2434                         k = isl_basic_map_alloc_div(bmap);
2435                         if (k < 0)
2436                                 goto error;
2437                         isl_seq_cpy(bmap->div[k], like->div[i], 1 + 1 + ltotal);
2438                         isl_seq_clr(bmap->div[k]+1+1+ltotal, total - ltotal);
2439                         if (add_div_constraints(bmap, k) < 0)
2440                                 goto error;
2441                 }
2442         }
2443         isl_basic_map_free(like);
2444         bmap = isl_basic_map_finalize(bmap);
2445         return bmap;
2446 error:
2447         isl_basic_map_free(like);
2448         isl_basic_set_free(bset);
2449         return NULL;
2450 }
2451
2452 struct isl_basic_set *isl_basic_set_from_underlying_set(
2453         struct isl_basic_set *bset, struct isl_basic_set *like)
2454 {
2455         return (struct isl_basic_set *)
2456                 isl_basic_map_overlying_set(bset, (struct isl_basic_map *)like);
2457 }
2458
2459 struct isl_set *isl_set_from_underlying_set(
2460         struct isl_set *set, struct isl_basic_set *like)
2461 {
2462         int i;
2463
2464         if (!set || !like)
2465                 goto error;
2466         isl_assert(set->ctx, set->dim == like->nparam + like->dim + like->n_div,
2467                     goto error);
2468         if (like->nparam == 0 && like->n_div == 0) {
2469                 isl_basic_set_free(like);
2470                 return set;
2471         }
2472         set = isl_set_cow(set);
2473         if (!set)
2474                 goto error;
2475         for (i = 0; i < set->n; ++i) {
2476                 set->p[i] = isl_basic_set_from_underlying_set(set->p[i],
2477                                                       isl_basic_set_copy(like));
2478                 if (!set->p[i])
2479                         goto error;
2480         }
2481         set->nparam = like->nparam;
2482         set->dim = like->dim;
2483         isl_basic_set_free(like);
2484         return set;
2485 error:
2486         isl_basic_set_free(like);
2487         isl_set_free(set);
2488         return NULL;
2489 }
2490
2491 struct isl_set *isl_map_underlying_set(struct isl_map *map)
2492 {
2493         int i;
2494
2495         map = isl_map_align_divs(map);
2496         map = isl_map_cow(map);
2497         if (!map)
2498                 return NULL;
2499
2500         for (i = 0; i < map->n; ++i) {
2501                 map->p[i] = (struct isl_basic_map *)
2502                                 isl_basic_map_underlying_set(map->p[i]);
2503                 if (!map->p[i])
2504                         goto error;
2505         }
2506         if (map->n == 0)
2507                 map->n_out += map->nparam + map->n_in;
2508         else
2509                 map->n_out = map->p[0]->n_out;
2510         map->nparam = 0;
2511         map->n_in = 0;
2512         return (struct isl_set *)map;
2513 error:
2514         isl_map_free(map);
2515         return NULL;
2516 }
2517
2518 struct isl_set *isl_set_to_underlying_set(struct isl_set *set)
2519 {
2520         return (struct isl_set *)isl_map_underlying_set((struct isl_map *)set);
2521 }
2522
2523 struct isl_basic_set *isl_basic_map_domain(struct isl_basic_map *bmap)
2524 {
2525         struct isl_basic_set *domain;
2526         unsigned n_out;
2527         if (!bmap)
2528                 return NULL;
2529         n_out = bmap->n_out;
2530         domain = isl_basic_set_from_basic_map(bmap);
2531         return isl_basic_set_project_out(domain, n_out, 0);
2532 }
2533
2534 struct isl_basic_set *isl_basic_map_range(struct isl_basic_map *bmap)
2535 {
2536         return isl_basic_map_domain(isl_basic_map_reverse(bmap));
2537 }
2538
2539 struct isl_set *isl_map_range(struct isl_map *map)
2540 {
2541         int i;
2542         struct isl_set *set;
2543
2544         if (!map)
2545                 goto error;
2546         map = isl_map_cow(map);
2547         if (!map)
2548                 goto error;
2549
2550         set = (struct isl_set *) map;
2551         set->zero = 0;
2552         for (i = 0; i < map->n; ++i) {
2553                 set->p[i] = isl_basic_map_range(map->p[i]);
2554                 if (!set->p[i])
2555                         goto error;
2556         }
2557         F_CLR(set, ISL_MAP_DISJOINT);
2558         F_CLR(set, ISL_SET_NORMALIZED);
2559         return set;
2560 error:
2561         isl_map_free(map);
2562         return NULL;
2563 }
2564
2565 struct isl_map *isl_map_from_set(struct isl_set *set,
2566                 unsigned n_in, unsigned n_out)
2567 {
2568         int i;
2569         struct isl_map *map = NULL;
2570
2571         if (!set)
2572                 return NULL;
2573         isl_assert(set->ctx, set->dim == n_in + n_out, goto error);
2574         set = isl_set_cow(set);
2575         if (!set)
2576                 return NULL;
2577         map = (struct isl_map *)set;
2578         for (i = 0; i < set->n; ++i) {
2579                 map->p[i] = isl_basic_map_from_basic_set(
2580                                 set->p[i], n_in, n_out);
2581                 if (!map->p[i])
2582                         goto error;
2583         }
2584         map->n_in = n_in;
2585         map->n_out = n_out;
2586         return map;
2587 error:
2588         isl_set_free(set);
2589         return NULL;
2590 }
2591
2592 struct isl_set *isl_set_from_map(struct isl_map *map)
2593 {
2594         int i;
2595         struct isl_set *set = NULL;
2596
2597         if (!map)
2598                 return NULL;
2599         map = isl_map_cow(map);
2600         if (!map)
2601                 return NULL;
2602         map->n_out += map->n_in;
2603         map->n_in = 0;
2604         set = (struct isl_set *)map;
2605         for (i = 0; i < map->n; ++i) {
2606                 set->p[i] = isl_basic_set_from_basic_map(map->p[i]);
2607                 if (!set->p[i])
2608                         goto error;
2609         }
2610         return set;
2611 error:
2612         isl_map_free(map);
2613         return NULL;
2614 }
2615
2616 struct isl_map *isl_map_alloc(struct isl_ctx *ctx,
2617                 unsigned nparam, unsigned in, unsigned out, int n,
2618                 unsigned flags)
2619 {
2620         struct isl_map *map;
2621
2622         map = isl_alloc(ctx, struct isl_map,
2623                         sizeof(struct isl_map) +
2624                         n * sizeof(struct isl_basic_map *));
2625         if (!map)
2626                 return NULL;
2627
2628         map->ctx = ctx;
2629         isl_ctx_ref(ctx);
2630         map->ref = 1;
2631         map->size = n;
2632         map->n = 0;
2633         map->nparam = nparam;
2634         map->n_in = in;
2635         map->n_out = out;
2636         map->flags = flags;
2637         return map;
2638 }
2639
2640 struct isl_basic_map *isl_basic_map_empty(struct isl_ctx *ctx,
2641                 unsigned nparam, unsigned in, unsigned out)
2642 {
2643         struct isl_basic_map *bmap;
2644         bmap = isl_basic_map_alloc(ctx, nparam, in, out, 0, 1, 0);
2645         bmap = isl_basic_map_set_to_empty(bmap);
2646         return bmap;
2647 }
2648
2649 struct isl_basic_set *isl_basic_set_empty(struct isl_ctx *ctx,
2650                 unsigned nparam, unsigned dim)
2651 {
2652         struct isl_basic_set *bset;
2653         bset = isl_basic_set_alloc(ctx, nparam, dim, 0, 1, 0);
2654         bset = isl_basic_set_set_to_empty(bset);
2655         return bset;
2656 }
2657
2658 struct isl_basic_map *isl_basic_map_universe(struct isl_ctx *ctx,
2659                 unsigned nparam, unsigned in, unsigned out)
2660 {
2661         struct isl_basic_map *bmap;
2662         bmap = isl_basic_map_alloc(ctx, nparam, in, out, 0, 0, 0);
2663         return bmap;
2664 }
2665
2666 struct isl_basic_set *isl_basic_set_universe(struct isl_ctx *ctx,
2667                 unsigned nparam, unsigned dim)
2668 {
2669         struct isl_basic_set *bset;
2670         bset = isl_basic_set_alloc(ctx, nparam, dim, 0, 0, 0);
2671         return bset;
2672 }
2673
2674 struct isl_map *isl_map_empty(struct isl_ctx *ctx,
2675                 unsigned nparam, unsigned in, unsigned out)
2676 {
2677         return isl_map_alloc(ctx, nparam, in, out, 0, ISL_MAP_DISJOINT);
2678 }
2679
2680 struct isl_set *isl_set_empty(struct isl_ctx *ctx,
2681                 unsigned nparam, unsigned dim)
2682 {
2683         return isl_set_alloc(ctx, nparam, dim, 0, ISL_MAP_DISJOINT);
2684 }
2685
2686 struct isl_map *isl_map_dup(struct isl_map *map)
2687 {
2688         int i;
2689         struct isl_map *dup;
2690
2691         if (!map)
2692                 return NULL;
2693         dup = isl_map_alloc(map->ctx, map->nparam, map->n_in, map->n_out, map->n,
2694                                 map->flags);
2695         for (i = 0; i < map->n; ++i)
2696                 dup = isl_map_add(dup, isl_basic_map_copy(map->p[i]));
2697         return dup;
2698 }
2699
2700 struct isl_map *isl_map_add(struct isl_map *map, struct isl_basic_map *bmap)
2701 {
2702         if (!bmap || !map)
2703                 goto error;
2704         isl_assert(map->ctx, map->nparam == bmap->nparam, goto error);
2705         isl_assert(map->ctx, map->n_in == bmap->n_in, goto error);
2706         isl_assert(map->ctx, map->n_out == bmap->n_out, goto error);
2707         isl_assert(map->ctx, map->n < map->size, goto error);
2708         map->p[map->n] = bmap;
2709         map->n++;
2710         F_CLR(map, ISL_MAP_NORMALIZED);
2711         return map;
2712 error:
2713         if (map)
2714                 isl_map_free(map);
2715         if (bmap)
2716                 isl_basic_map_free(bmap);
2717         return NULL;
2718 }
2719
2720 void isl_map_free(struct isl_map *map)
2721 {
2722         int i;
2723
2724         if (!map)
2725                 return;
2726
2727         if (--map->ref > 0)
2728                 return;
2729
2730         isl_ctx_deref(map->ctx);
2731         for (i = 0; i < map->n; ++i)
2732                 isl_basic_map_free(map->p[i]);
2733         free(map);
2734 }
2735
2736 struct isl_map *isl_map_extend(struct isl_map *base,
2737                 unsigned nparam, unsigned n_in, unsigned n_out)
2738 {
2739         int i;
2740
2741         base = isl_map_cow(base);
2742         if (!base)
2743                 return NULL;
2744
2745         isl_assert(base->ctx, base->nparam <= nparam, goto error);
2746         isl_assert(base->ctx, base->n_in <= n_in, goto error);
2747         isl_assert(base->ctx, base->n_out <= n_out, goto error);
2748         base->nparam = nparam;
2749         base->n_in = n_in;
2750         base->n_out = n_out;
2751         for (i = 0; i < base->n; ++i) {
2752                 base->p[i] = isl_basic_map_extend(base->p[i],
2753                                 nparam, n_in, n_out, 0, 0, 0);
2754                 if (!base->p[i])
2755                         goto error;
2756         }
2757         return base;
2758 error:
2759         isl_map_free(base);
2760         return NULL;
2761 }
2762
2763 struct isl_set *isl_set_extend(struct isl_set *base,
2764                 unsigned nparam, unsigned dim)
2765 {
2766         return (struct isl_set *)isl_map_extend((struct isl_map *)base,
2767                                                         nparam, 0, dim);
2768 }
2769
2770 static struct isl_basic_map *isl_basic_map_fix_var(struct isl_basic_map *bmap,
2771                 unsigned var, int value)
2772 {
2773         int j;
2774
2775         bmap = isl_basic_map_cow(bmap);
2776         bmap = isl_basic_map_extend(bmap,
2777                         bmap->nparam, bmap->n_in, bmap->n_out, 0, 1, 0);
2778         j = isl_basic_map_alloc_equality(bmap);
2779         if (j < 0)
2780                 goto error;
2781         isl_seq_clr(bmap->eq[j],
2782                     1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div);
2783         isl_int_set_si(bmap->eq[j][1 + var], -1);
2784         isl_int_set_si(bmap->eq[j][0], value);
2785         bmap = isl_basic_map_simplify(bmap);
2786         return isl_basic_map_finalize(bmap);
2787 error:
2788         isl_basic_map_free(bmap);
2789         return NULL;
2790 }
2791
2792 struct isl_basic_map *isl_basic_map_fix_input_si(struct isl_basic_map *bmap,
2793                 unsigned input, int value)
2794 {
2795         if (!bmap)
2796                 return NULL;
2797         isl_assert(bmap->ctx, input < bmap->n_in, goto error);
2798         return isl_basic_map_fix_var(bmap, bmap->nparam + input, value);
2799 error:
2800         isl_basic_map_free(bmap);
2801         return NULL;
2802 }
2803
2804 struct isl_basic_set *isl_basic_set_fix_dim_si(struct isl_basic_set *bset,
2805                 unsigned dim, int value)
2806 {
2807         if (!bset)
2808                 return NULL;
2809         isl_assert(bset->ctx, dim < bset->dim, goto error);
2810         return (struct isl_basic_set *)
2811                 isl_basic_map_fix_var((struct isl_basic_map *)bset,
2812                                                 bset->nparam + dim, value);
2813 error:
2814         isl_basic_set_free(bset);
2815         return NULL;
2816 }
2817
2818 struct isl_map *isl_map_fix_input_si(struct isl_map *map,
2819                 unsigned input, int value)
2820 {
2821         int i;
2822
2823         map = isl_map_cow(map);
2824         if (!map)
2825                 return NULL;
2826
2827         isl_assert(ctx, input < map->n_in, goto error);
2828         for (i = 0; i < map->n; ++i) {
2829                 map->p[i] = isl_basic_map_fix_input_si(map->p[i], input, value);
2830                 if (!map->p[i])
2831                         goto error;
2832         }
2833         F_CLR(map, ISL_MAP_NORMALIZED);
2834         return map;
2835 error:
2836         isl_map_free(map);
2837         return NULL;
2838 }
2839
2840 struct isl_set *isl_set_fix_dim_si(struct isl_set *set, unsigned dim, int value)
2841 {
2842         int i;
2843
2844         set = isl_set_cow(set);
2845         if (!set)
2846                 return NULL;
2847
2848         isl_assert(set->ctx, dim < set->dim, goto error);
2849         for (i = 0; i < set->n; ++i) {
2850                 set->p[i] = isl_basic_set_fix_dim_si(set->p[i], dim, value);
2851                 if (!set->p[i])
2852                         goto error;
2853         }
2854         return set;
2855 error:
2856         isl_set_free(set);
2857         return NULL;
2858 }
2859
2860 struct isl_basic_set *isl_basic_set_lower_bound_dim(struct isl_basic_set *bset,
2861         unsigned dim, isl_int value)
2862 {
2863         int j;
2864
2865         bset = isl_basic_set_cow(bset);
2866         bset = isl_basic_set_extend(bset, bset->nparam, bset->dim, 0, 0, 1);
2867         j = isl_basic_set_alloc_inequality(bset);
2868         if (j < 0)
2869                 goto error;
2870         isl_seq_clr(bset->ineq[j], 1 + bset->nparam + bset->dim + bset->n_div);
2871         isl_int_set_si(bset->ineq[j][1 + bset->nparam + dim], 1);
2872         isl_int_neg(bset->ineq[j][0], value);
2873         bset = isl_basic_set_simplify(bset);
2874         return isl_basic_set_finalize(bset);
2875 error:
2876         isl_basic_set_free(bset);
2877         return NULL;
2878 }
2879
2880 struct isl_set *isl_set_lower_bound_dim(struct isl_set *set, unsigned dim,
2881                                         isl_int value)
2882 {
2883         int i;
2884
2885         set = isl_set_cow(set);
2886         if (!set)
2887                 return NULL;
2888
2889         isl_assert(set->ctx, dim < set->dim, goto error);
2890         for (i = 0; i < set->n; ++i) {
2891                 set->p[i] = isl_basic_set_lower_bound_dim(set->p[i], dim, value);
2892                 if (!set->p[i])
2893                         goto error;
2894         }
2895         return set;
2896 error:
2897         isl_set_free(set);
2898         return NULL;
2899 }
2900
2901 struct isl_map *isl_map_reverse(struct isl_map *map)
2902 {
2903         int i;
2904         unsigned t;
2905
2906         map = isl_map_cow(map);
2907         if (!map)
2908                 return NULL;
2909
2910         t = map->n_in;
2911         map->n_in = map->n_out;
2912         map->n_out = t;
2913         for (i = 0; i < map->n; ++i) {
2914                 map->p[i] = isl_basic_map_reverse(map->p[i]);
2915                 if (!map->p[i])
2916                         goto error;
2917         }
2918         F_CLR(map, ISL_MAP_NORMALIZED);
2919         return map;
2920 error:
2921         isl_map_free(map);
2922         return NULL;
2923 }
2924
2925 struct isl_map *isl_basic_map_lexmax(
2926                 struct isl_basic_map *bmap, struct isl_basic_set *dom,
2927                 struct isl_set **empty)
2928 {
2929         return isl_pip_basic_map_lexmax(bmap, dom, empty);
2930 }
2931
2932 struct isl_map *isl_basic_map_lexmin(
2933                 struct isl_basic_map *bmap, struct isl_basic_set *dom,
2934                 struct isl_set **empty)
2935 {
2936         return isl_pip_basic_map_lexmin(bmap, dom, empty);
2937 }
2938
2939 struct isl_set *isl_basic_set_lexmin(struct isl_basic_set *bset)
2940 {
2941         struct isl_basic_map *bmap = NULL;
2942         struct isl_basic_set *dom = NULL;
2943         struct isl_map *min;
2944
2945         if (!bset)
2946                 goto error;
2947         bmap = isl_basic_map_from_basic_set(bset, 0, bset->dim);
2948         if (!bmap)
2949                 goto error;
2950         dom = isl_basic_set_alloc(bmap->ctx, bmap->nparam, 0, 0, 0, 0);
2951         if (!dom)
2952                 goto error;
2953         min = isl_basic_map_lexmin(bmap, dom, NULL);
2954         return isl_map_range(min);
2955 error:
2956         isl_basic_map_free(bmap);
2957         return NULL;
2958 }
2959
2960 struct isl_map *isl_basic_map_compute_divs(struct isl_basic_map *bmap)
2961 {
2962         if (!bmap)
2963                 return NULL;
2964         if (bmap->n_div == 0)
2965                 return isl_map_from_basic_map(bmap);
2966         return isl_pip_basic_map_compute_divs(bmap);
2967 }
2968
2969 struct isl_map *isl_map_compute_divs(struct isl_map *map)
2970 {
2971         int i;
2972         struct isl_map *res;
2973
2974         if (!map)
2975                 return NULL;
2976         if (map->n == 0)
2977                 return map;
2978         res = isl_basic_map_compute_divs(isl_basic_map_copy(map->p[0]));
2979         for (i = 1 ; i < map->n; ++i) {
2980                 struct isl_map *r2;
2981                 r2 = isl_basic_map_compute_divs(isl_basic_map_copy(map->p[i]));
2982                 if (F_ISSET(map, ISL_MAP_DISJOINT))
2983                         res = isl_map_union_disjoint(res, r2);
2984                 else
2985                         res = isl_map_union(res, r2);
2986         }
2987         isl_map_free(map);
2988
2989         return res;
2990 }
2991
2992 struct isl_set *isl_basic_set_compute_divs(struct isl_basic_set *bset)
2993 {
2994         return (struct isl_set *)
2995                 isl_basic_map_compute_divs((struct isl_basic_map *)bset);
2996 }
2997
2998 struct isl_set *isl_set_compute_divs(struct isl_set *set)
2999 {
3000         return (struct isl_set *)
3001                 isl_map_compute_divs((struct isl_map *)set);
3002 }
3003
3004 struct isl_set *isl_map_domain(struct isl_map *map)
3005 {
3006         int i;
3007         struct isl_set *set;
3008
3009         if (!map)
3010                 goto error;
3011
3012         map = isl_map_cow(map);
3013         if (!map)
3014                 return NULL;
3015
3016         set = (struct isl_set *)map;
3017         set->dim = map->n_in;
3018         set->zero = 0;
3019         for (i = 0; i < map->n; ++i) {
3020                 set->p[i] = isl_basic_map_domain(map->p[i]);
3021                 if (!set->p[i])
3022                         goto error;
3023         }
3024         F_CLR(set, ISL_MAP_DISJOINT);
3025         F_CLR(set, ISL_SET_NORMALIZED);
3026         return set;
3027 error:
3028         isl_map_free(map);
3029         return NULL;
3030 }
3031
3032 struct isl_map *isl_map_union_disjoint(
3033                         struct isl_map *map1, struct isl_map *map2)
3034 {
3035         int i;
3036         unsigned flags = 0;
3037         struct isl_map *map = NULL;
3038
3039         if (!map1 || !map2)
3040                 goto error;
3041
3042         if (map1->n == 0) {
3043                 isl_map_free(map1);
3044                 return map2;
3045         }
3046         if (map2->n == 0) {
3047                 isl_map_free(map2);
3048                 return map1;
3049         }
3050
3051         isl_assert(map1->ctx, map1->nparam == map2->nparam, goto error);
3052         isl_assert(map1->ctx, map1->n_in == map2->n_in, goto error);
3053         isl_assert(map1->ctx, map1->n_out == map2->n_out, goto error);
3054
3055         if (F_ISSET(map1, ISL_MAP_DISJOINT) &&
3056             F_ISSET(map2, ISL_MAP_DISJOINT))
3057                 FL_SET(flags, ISL_MAP_DISJOINT);
3058
3059         map = isl_map_alloc(map1->ctx, map1->nparam, map1->n_in, map1->n_out,
3060                                 map1->n + map2->n, flags);
3061         if (!map)
3062                 goto error;
3063         for (i = 0; i < map1->n; ++i) {
3064                 map = isl_map_add(map,
3065                                   isl_basic_map_copy(map1->p[i]));
3066                 if (!map)
3067                         goto error;
3068         }
3069         for (i = 0; i < map2->n; ++i) {
3070                 map = isl_map_add(map,
3071                                   isl_basic_map_copy(map2->p[i]));
3072                 if (!map)
3073                         goto error;
3074         }
3075         isl_map_free(map1);
3076         isl_map_free(map2);
3077         return map;
3078 error:
3079         isl_map_free(map);
3080         isl_map_free(map1);
3081         isl_map_free(map2);
3082         return NULL;
3083 }
3084
3085 struct isl_map *isl_map_union(struct isl_map *map1, struct isl_map *map2)
3086 {
3087         map1 = isl_map_union_disjoint(map1, map2);
3088         if (!map1)
3089                 return NULL;
3090         if (map1->n > 1)
3091                 F_CLR(map1, ISL_MAP_DISJOINT);
3092         return map1;
3093 }
3094
3095 struct isl_set *isl_set_union_disjoint(
3096                         struct isl_set *set1, struct isl_set *set2)
3097 {
3098         return (struct isl_set *)
3099                 isl_map_union_disjoint(
3100                         (struct isl_map *)set1, (struct isl_map *)set2);
3101 }
3102
3103 struct isl_set *isl_set_union(struct isl_set *set1, struct isl_set *set2)
3104 {
3105         return (struct isl_set *)
3106                 isl_map_union((struct isl_map *)set1, (struct isl_map *)set2);
3107 }
3108
3109 struct isl_map *isl_map_intersect_range(
3110                 struct isl_map *map, struct isl_set *set)
3111 {
3112         unsigned flags = 0;
3113         struct isl_map *result;
3114         int i, j;
3115
3116         if (!map || !set)
3117                 goto error;
3118
3119         if (F_ISSET(map, ISL_MAP_DISJOINT) &&
3120             F_ISSET(set, ISL_MAP_DISJOINT))
3121                 FL_SET(flags, ISL_MAP_DISJOINT);
3122
3123         result = isl_map_alloc(map->ctx, map->nparam, map->n_in, map->n_out,
3124                                 map->n * set->n, flags);
3125         if (!result)
3126                 goto error;
3127         for (i = 0; i < map->n; ++i)
3128                 for (j = 0; j < set->n; ++j) {
3129                         result = isl_map_add(result,
3130                             isl_basic_map_intersect_range(
3131                                 isl_basic_map_copy(map->p[i]),
3132                                 isl_basic_set_copy(set->p[j])));
3133                         if (!result)
3134                                 goto error;
3135                 }
3136         isl_map_free(map);
3137         isl_set_free(set);
3138         return result;
3139 error:
3140         isl_map_free(map);
3141         isl_set_free(set);
3142         return NULL;
3143 }
3144
3145 struct isl_map *isl_map_intersect_domain(
3146                 struct isl_map *map, struct isl_set *set)
3147 {
3148         return isl_map_reverse(
3149                 isl_map_intersect_range(isl_map_reverse(map), set));
3150 }
3151
3152 struct isl_map *isl_map_apply_domain(
3153                 struct isl_map *map1, struct isl_map *map2)
3154 {
3155         if (!map1 || !map2)
3156                 goto error;
3157         map1 = isl_map_reverse(map1);
3158         map1 = isl_map_apply_range(map1, map2);
3159         return isl_map_reverse(map1);
3160 error:
3161         isl_map_free(map1);
3162         isl_map_free(map2);
3163         return NULL;
3164 }
3165
3166 struct isl_map *isl_map_apply_range(
3167                 struct isl_map *map1, struct isl_map *map2)
3168 {
3169         struct isl_map *result;
3170         int i, j;
3171
3172         if (!map1 || !map2)
3173                 goto error;
3174
3175         isl_assert(map1->ctx, map1->nparam == map2->nparam, goto error);
3176         isl_assert(map1->ctx, map1->n_out == map2->n_in, goto error);
3177
3178         result = isl_map_alloc(map1->ctx, map1->nparam, map1->n_in, map2->n_out,
3179                                 map1->n * map2->n, 0);
3180         if (!result)
3181                 goto error;
3182         for (i = 0; i < map1->n; ++i)
3183                 for (j = 0; j < map2->n; ++j) {
3184                         result = isl_map_add(result,
3185                             isl_basic_map_apply_range(
3186                                 isl_basic_map_copy(map1->p[i]),
3187                                 isl_basic_map_copy(map2->p[j])));
3188                         if (!result)
3189                                 goto error;
3190                 }
3191         isl_map_free(map1);
3192         isl_map_free(map2);
3193         if (result && result->n <= 1)
3194                 F_SET(result, ISL_MAP_DISJOINT);
3195         return result;
3196 error:
3197         isl_map_free(map1);
3198         isl_map_free(map2);
3199         return NULL;
3200 }
3201
3202 /*
3203  * returns range - domain
3204  */
3205 struct isl_basic_set *isl_basic_map_deltas(struct isl_basic_map *bmap)
3206 {
3207         struct isl_basic_set *bset;
3208         unsigned dim;
3209         int i;
3210
3211         if (!bmap)
3212                 goto error;
3213         isl_assert(bmap->ctx, bmap->n_in == bmap->n_out, goto error);
3214         dim = bmap->n_in;
3215         bset = isl_basic_set_from_basic_map(bmap);
3216         bset = isl_basic_set_extend(bset, bset->nparam, 3*dim, 0,
3217                                         dim, 0);
3218         bset = isl_basic_set_swap_vars(bset, 2*dim);
3219         for (i = 0; i < dim; ++i) {
3220                 int j = isl_basic_map_alloc_equality(
3221                                             (struct isl_basic_map *)bset);
3222                 if (j < 0)
3223                         goto error;
3224                 isl_seq_clr(bset->eq[j],
3225                             1 + bset->nparam + bset->dim + bset->n_div);
3226                 isl_int_set_si(bset->eq[j][1+bset->nparam+i], 1);
3227                 isl_int_set_si(bset->eq[j][1+bset->nparam+dim+i], 1);
3228                 isl_int_set_si(bset->eq[j][1+bset->nparam+2*dim+i], -1);
3229         }
3230         return isl_basic_set_project_out(bset, 2*dim, 0);
3231 error:
3232         isl_basic_map_free(bmap);
3233         return NULL;
3234 }
3235
3236 /*
3237  * returns range - domain
3238  */
3239 struct isl_set *isl_map_deltas(struct isl_map *map)
3240 {
3241         int i;
3242         struct isl_set *result;
3243
3244         if (!map)
3245                 return NULL;
3246
3247         isl_assert(map->ctx, map->n_in == map->n_out, goto error);
3248         result = isl_set_alloc(map->ctx, map->nparam, map->n_in, map->n, map->flags);
3249         if (!result)
3250                 goto error;
3251         for (i = 0; i < map->n; ++i)
3252                 result = isl_set_add(result,
3253                           isl_basic_map_deltas(isl_basic_map_copy(map->p[i])));
3254         isl_map_free(map);
3255         return result;
3256 error:
3257         isl_map_free(map);
3258         return NULL;
3259 }
3260
3261 /* If the only constraints a div d=floor(f/m)
3262  * appears in are its two defining constraints
3263  *
3264  *      f - m d >=0
3265  *      -(f - (m - 1)) + m d >= 0
3266  *
3267  * then it can safely be removed.
3268  */
3269 static int div_is_redundant(struct isl_basic_map *bmap, int div)
3270 {
3271         int i;
3272         unsigned pos = 1 + bmap->nparam + bmap->n_in + bmap->n_out + div;
3273
3274         for (i = 0; i < bmap->n_eq; ++i)
3275                 if (!isl_int_is_zero(bmap->eq[i][pos]))
3276                         return 0;
3277
3278         for (i = 0; i < bmap->n_ineq; ++i) {
3279                 if (isl_int_is_zero(bmap->ineq[i][pos]))
3280                         continue;
3281                 if (isl_int_eq(bmap->ineq[i][pos], bmap->div[div][0])) {
3282                         int neg;
3283                         isl_int_sub(bmap->div[div][1],
3284                                         bmap->div[div][1], bmap->div[div][0]);
3285                         isl_int_add_ui(bmap->div[div][1], bmap->div[div][1], 1);
3286                         neg = isl_seq_is_neg(bmap->ineq[i], bmap->div[div]+1, pos);
3287                         isl_int_sub_ui(bmap->div[div][1], bmap->div[div][1], 1);
3288                         isl_int_add(bmap->div[div][1],
3289                                         bmap->div[div][1], bmap->div[div][0]);
3290                         if (!neg)
3291                                 return 0;
3292                         if (isl_seq_first_non_zero(bmap->ineq[i]+pos+1,
3293                                                     bmap->n_div-div-1) != -1)
3294                                 return 0;
3295                 } else if (isl_int_abs_eq(bmap->ineq[i][pos], bmap->div[div][0])) {
3296                         if (!isl_seq_eq(bmap->ineq[i], bmap->div[div]+1, pos))
3297                                 return 0;
3298                         if (isl_seq_first_non_zero(bmap->ineq[i]+pos+1,
3299                                                     bmap->n_div-div-1) != -1)
3300                                 return 0;
3301                 } else
3302                         return 0;
3303         }
3304
3305         for (i = 0; i < bmap->n_div; ++i)
3306                 if (!isl_int_is_zero(bmap->div[i][1+pos]))
3307                         return 0;
3308
3309         return 1;
3310 }
3311
3312 /*
3313  * Remove divs that don't occur in any of the constraints or other divs.
3314  * These can arise when dropping some of the variables in a quast
3315  * returned by piplib.
3316  */
3317 static struct isl_basic_map *remove_redundant_divs(struct isl_basic_map *bmap)
3318 {
3319         int i;
3320
3321         if (!bmap)
3322                 return NULL;
3323
3324         for (i = bmap->n_div-1; i >= 0; --i) {
3325                 if (!div_is_redundant(bmap, i))
3326                         continue;
3327                 bmap = isl_basic_map_drop_div(bmap, i);
3328         }
3329         return bmap;
3330 }
3331
3332 struct isl_basic_map *isl_basic_map_finalize(struct isl_basic_map *bmap)
3333 {
3334         bmap = remove_redundant_divs(bmap);
3335         if (!bmap)
3336                 return NULL;
3337         F_SET(bmap, ISL_BASIC_SET_FINAL);
3338         return bmap;
3339 }
3340
3341 struct isl_basic_set *isl_basic_set_finalize(struct isl_basic_set *bset)
3342 {
3343         return (struct isl_basic_set *)
3344                 isl_basic_map_finalize((struct isl_basic_map *)bset);
3345 }
3346
3347 struct isl_set *isl_set_finalize(struct isl_set *set)
3348 {
3349         int i;
3350
3351         if (!set)
3352                 return NULL;
3353         for (i = 0; i < set->n; ++i) {
3354                 set->p[i] = isl_basic_set_finalize(set->p[i]);
3355                 if (!set->p[i])
3356                         goto error;
3357         }
3358         return set;
3359 error:
3360         isl_set_free(set);
3361         return NULL;
3362 }
3363
3364 struct isl_map *isl_map_finalize(struct isl_map *map)
3365 {
3366         int i;
3367
3368         if (!map)
3369                 return NULL;
3370         for (i = 0; i < map->n; ++i) {
3371                 map->p[i] = isl_basic_map_finalize(map->p[i]);
3372                 if (!map->p[i])
3373                         goto error;
3374         }
3375         F_CLR(map, ISL_MAP_NORMALIZED);
3376         return map;
3377 error:
3378         isl_map_free(map);
3379         return NULL;
3380 }
3381
3382 struct isl_basic_map *isl_basic_map_identity(struct isl_ctx *ctx,
3383                 unsigned nparam, unsigned dim)
3384 {
3385         struct isl_basic_map *bmap;
3386         int i;
3387
3388         bmap = isl_basic_map_alloc(ctx, nparam, dim, dim, 0, dim, 0);
3389         if (!bmap)
3390                 goto error;
3391
3392         for (i = 0; i < dim; ++i) {
3393                 int j = isl_basic_map_alloc_equality(bmap);
3394                 if (j < 0)
3395                         goto error;
3396                 isl_seq_clr(bmap->eq[j],
3397                     1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div);
3398                 isl_int_set_si(bmap->eq[j][1+nparam+i], 1);
3399                 isl_int_set_si(bmap->eq[j][1+nparam+dim+i], -1);
3400         }
3401         return isl_basic_map_finalize(bmap);
3402 error:
3403         isl_basic_map_free(bmap);
3404         return NULL;
3405 }
3406
3407 struct isl_map *isl_map_identity(struct isl_ctx *ctx,
3408                 unsigned nparam, unsigned dim)
3409 {
3410         struct isl_map *map = isl_map_alloc(ctx, nparam, dim, dim, 1,
3411                                                 ISL_MAP_DISJOINT);
3412         if (!map)
3413                 goto error;
3414         map = isl_map_add(map,
3415                 isl_basic_map_identity(ctx, nparam, dim));
3416         return map;
3417 error:
3418         isl_map_free(map);
3419         return NULL;
3420 }
3421
3422 int isl_set_is_equal(struct isl_set *set1, struct isl_set *set2)
3423 {
3424         return isl_map_is_equal((struct isl_map *)set1, (struct isl_map *)set2);
3425 }
3426
3427 int isl_set_is_subset(struct isl_set *set1, struct isl_set *set2)
3428 {
3429         return isl_map_is_subset(
3430                         (struct isl_map *)set1, (struct isl_map *)set2);
3431 }
3432
3433 int isl_basic_map_is_subset(
3434                 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
3435 {
3436         int is_subset;
3437         struct isl_map *map1;
3438         struct isl_map *map2;
3439
3440         if (!bmap1 || !bmap2)
3441                 return -1;
3442
3443         map1 = isl_map_from_basic_map(isl_basic_map_copy(bmap1));
3444         map2 = isl_map_from_basic_map(isl_basic_map_copy(bmap2));
3445
3446         is_subset = isl_map_is_subset(map1, map2);
3447
3448         isl_map_free(map1);
3449         isl_map_free(map2);
3450
3451         return is_subset;
3452 }
3453
3454 int isl_basic_map_is_equal(
3455                 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
3456 {
3457         int is_subset;
3458
3459         if (!bmap1 || !bmap2)
3460                 return -1;
3461         is_subset = isl_basic_map_is_subset(bmap1, bmap2);
3462         if (is_subset != 1)
3463                 return is_subset;
3464         is_subset = isl_basic_map_is_subset(bmap2, bmap1);
3465         return is_subset;
3466 }
3467
3468 int isl_basic_set_is_equal(
3469                 struct isl_basic_set *bset1, struct isl_basic_set *bset2)
3470 {
3471         return isl_basic_map_is_equal(
3472                 (struct isl_basic_map *)bset1, (struct isl_basic_map *)bset2);
3473 }
3474
3475 int isl_map_is_empty(struct isl_map *map)
3476 {
3477         int i;
3478         int is_empty;
3479
3480         if (!map)
3481                 return -1;
3482         for (i = 0; i < map->n; ++i) {
3483                 is_empty = isl_basic_map_is_empty(map->p[i]);
3484                 if (is_empty < 0)
3485                         return -1;
3486                 if (!is_empty)
3487                         return 0;
3488         }
3489         return 1;
3490 }
3491
3492 int isl_set_is_empty(struct isl_set *set)
3493 {
3494         return isl_map_is_empty((struct isl_map *)set);
3495 }
3496
3497 int isl_map_is_subset(struct isl_map *map1, struct isl_map *map2)
3498 {
3499         int i;
3500         int is_subset = 0;
3501         struct isl_map *diff;
3502
3503         if (!map1 || !map2)
3504                 return -1;
3505
3506         if (isl_map_is_empty(map1))
3507                 return 1;
3508
3509         if (isl_map_is_empty(map2))
3510                 return 0;
3511
3512         diff = isl_map_subtract(isl_map_copy(map1), isl_map_copy(map2));
3513         if (!diff)
3514                 return -1;
3515
3516         is_subset = isl_map_is_empty(diff);
3517         isl_map_free(diff);
3518
3519         return is_subset;
3520 }
3521
3522 int isl_map_is_equal(struct isl_map *map1, struct isl_map *map2)
3523 {
3524         int is_subset;
3525
3526         if (!map1 || !map2)
3527                 return -1;
3528         is_subset = isl_map_is_subset(map1, map2);
3529         if (is_subset != 1)
3530                 return is_subset;
3531         is_subset = isl_map_is_subset(map2, map1);
3532         return is_subset;
3533 }
3534
3535 int isl_basic_map_is_strict_subset(
3536                 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
3537 {
3538         int is_subset;
3539
3540         if (!bmap1 || !bmap2)
3541                 return -1;
3542         is_subset = isl_basic_map_is_subset(bmap1, bmap2);
3543         if (is_subset != 1)
3544                 return is_subset;
3545         is_subset = isl_basic_map_is_subset(bmap2, bmap1);
3546         if (is_subset == -1)
3547                 return is_subset;
3548         return !is_subset;
3549 }
3550
3551 static int basic_map_contains(struct isl_basic_map *bmap, struct isl_vec *vec)
3552 {
3553         int i;
3554         unsigned total;
3555         isl_int s;
3556
3557         total = 1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
3558         if (total != vec->size)
3559                 return -1;
3560
3561         isl_int_init(s);
3562
3563         for (i = 0; i < bmap->n_eq; ++i) {
3564                 isl_seq_inner_product(vec->block.data, bmap->eq[i], total, &s);
3565                 if (!isl_int_is_zero(s)) {
3566                         isl_int_clear(s);
3567                         return 0;
3568                 }
3569         }
3570
3571         for (i = 0; i < bmap->n_ineq; ++i) {
3572                 isl_seq_inner_product(vec->block.data, bmap->ineq[i], total, &s);
3573                 if (isl_int_is_neg(s)) {
3574                         isl_int_clear(s);
3575                         return 0;
3576                 }
3577         }
3578
3579         isl_int_clear(s);
3580
3581         return 1;
3582 }
3583
3584 int isl_basic_map_is_universe(struct isl_basic_map *bmap)
3585 {
3586         if (!bmap)
3587                 return -1;
3588         return bmap->n_eq == 0 && bmap->n_ineq == 0;
3589 }
3590
3591 int isl_basic_map_is_empty(struct isl_basic_map *bmap)
3592 {
3593         struct isl_basic_set *bset = NULL;
3594         struct isl_vec *sample = NULL;
3595         int empty;
3596         unsigned total;
3597
3598         if (!bmap)
3599                 return -1;
3600
3601         if (F_ISSET(bmap, ISL_BASIC_MAP_EMPTY))
3602                 return 1;
3603
3604         if (F_ISSET(bmap, ISL_BASIC_MAP_RATIONAL)) {
3605                 struct isl_basic_map *copy = isl_basic_map_copy(bmap);
3606                 copy = isl_basic_map_convex_hull(copy);
3607                 empty = F_ISSET(copy, ISL_BASIC_MAP_EMPTY);
3608                 isl_basic_map_free(copy);
3609                 return empty;
3610         }
3611
3612         total = 1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
3613         if (bmap->sample && bmap->sample->size == total) {
3614                 int contains = basic_map_contains(bmap, bmap->sample);
3615                 if (contains < 0)
3616                         return -1;
3617                 if (contains)
3618                         return 0;
3619         }
3620         bset = isl_basic_map_underlying_set(isl_basic_map_copy(bmap));
3621         if (!bset)
3622                 return -1;
3623         sample = isl_basic_set_sample(bset);
3624         if (!sample)
3625                 return -1;
3626         empty = sample->size == 0;
3627         if (bmap->sample)
3628                 isl_vec_free(bmap->ctx, bmap->sample);
3629         bmap->sample = sample;
3630
3631         return empty;
3632 }
3633
3634 struct isl_map *isl_basic_map_union(
3635         struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
3636 {
3637         struct isl_map *map;
3638         if (!bmap1 || !bmap2)
3639                 return NULL;
3640
3641         isl_assert(bmap1->ctx, bmap1->nparam == bmap2->nparam, goto error);
3642         isl_assert(bmap1->ctx, bmap1->n_in == bmap2->n_in, goto error);
3643         isl_assert(bmap1->ctx, bmap1->n_out == bmap2->n_out, goto error);
3644
3645         map = isl_map_alloc(bmap1->ctx, bmap1->nparam,
3646                                 bmap1->n_in, bmap1->n_out, 2, 0);
3647         if (!map)
3648                 goto error;
3649         map = isl_map_add(map, bmap1);
3650         map = isl_map_add(map, bmap2);
3651         return map;
3652 error:
3653         isl_basic_map_free(bmap1);
3654         isl_basic_map_free(bmap2);
3655         return NULL;
3656 }
3657
3658 struct isl_set *isl_basic_set_union(
3659                 struct isl_basic_set *bset1, struct isl_basic_set *bset2)
3660 {
3661         return (struct isl_set *)isl_basic_map_union(
3662                                             (struct isl_basic_map *)bset1,
3663                                             (struct isl_basic_map *)bset2);
3664 }
3665
3666 /* Order divs such that any div only depends on previous divs */
3667 static struct isl_basic_map *order_divs(struct isl_basic_map *bmap)
3668 {
3669         int i;
3670         unsigned off = bmap->nparam + bmap->n_in + bmap->n_out;
3671
3672         for (i = 0; i < bmap->n_div; ++i) {
3673                 int pos;
3674                 pos = isl_seq_first_non_zero(bmap->div[i]+1+1+off+i,
3675                                                             bmap->n_div-i);
3676                 if (pos == -1)
3677                         continue;
3678                 swap_div(bmap, i, pos);
3679                 --i;
3680         }
3681         return bmap;
3682 }
3683
3684 static int find_div(struct isl_basic_map *dst,
3685                         struct isl_basic_map *src, unsigned div)
3686 {
3687         int i;
3688
3689         unsigned total = src->nparam + src->n_in + src->n_out + src->n_div;
3690
3691         for (i = 0; i < dst->n_div; ++i)
3692                 if (isl_seq_eq(dst->div[i], src->div[div], 1+1+total) &&
3693                     isl_seq_first_non_zero(dst->div[i]+1+1+total,
3694                                                 dst->n_div - src->n_div) == -1)
3695                         return i;
3696         return -1;
3697 }
3698
3699 struct isl_basic_map *isl_basic_map_align_divs(
3700                 struct isl_basic_map *dst, struct isl_basic_map *src)
3701 {
3702         int i;
3703         unsigned total = src->nparam + src->n_in + src->n_out + src->n_div;
3704
3705         if (!dst || !src)
3706                 goto error;
3707
3708         if (src->n_div == 0)
3709                 return dst;
3710
3711         src = order_divs(src);
3712         dst = isl_basic_map_extend(dst, dst->nparam, dst->n_in,
3713                         dst->n_out, src->n_div, 0, 2 * src->n_div);
3714         if (!dst)
3715                 return NULL;
3716         for (i = 0; i < src->n_div; ++i) {
3717                 int j = find_div(dst, src, i);
3718                 if (j < 0) {
3719                         j = isl_basic_map_alloc_div(dst);
3720                         if (j < 0)
3721                                 goto error;
3722                         isl_seq_cpy(dst->div[j], src->div[i], 1+1+total);
3723                         isl_seq_clr(dst->div[j]+1+1+total,
3724                                                     dst->n_div - src->n_div);
3725                         if (add_div_constraints(dst, j) < 0)
3726                                 goto error;
3727                 }
3728                 if (j != i)
3729                         swap_div(dst, i, j);
3730         }
3731         return dst;
3732 error:
3733         isl_basic_map_free(dst);
3734         return NULL;
3735 }
3736
3737 struct isl_basic_set *isl_basic_set_align_divs(
3738                 struct isl_basic_set *dst, struct isl_basic_set *src)
3739 {
3740         return (struct isl_basic_set *)isl_basic_map_align_divs(
3741                 (struct isl_basic_map *)dst, (struct isl_basic_map *)src);
3742 }
3743
3744 struct isl_map *isl_map_align_divs(struct isl_map *map)
3745 {
3746         int i;
3747
3748         map = isl_map_compute_divs(map);
3749         map = isl_map_cow(map);
3750         if (!map)
3751                 return NULL;
3752
3753         for (i = 1; i < map->n; ++i)
3754                 map->p[0] = isl_basic_map_align_divs(map->p[0], map->p[i]);
3755         for (i = 1; i < map->n; ++i)
3756                 map->p[i] = isl_basic_map_align_divs(map->p[i], map->p[0]);
3757
3758         F_CLR(map, ISL_MAP_NORMALIZED);
3759         return map;
3760 }
3761
3762 static struct isl_map *add_cut_constraint(struct isl_map *dst,
3763                 struct isl_basic_map *src, isl_int *c,
3764                 unsigned len, int oppose)
3765 {
3766         struct isl_basic_map *copy = NULL;
3767         int is_empty;
3768         int k;
3769         unsigned total;
3770
3771         copy = isl_basic_map_copy(src);
3772         copy = isl_basic_map_cow(copy);
3773         if (!copy)
3774                 goto error;
3775         copy = isl_basic_map_extend(copy,
3776                 copy->nparam, copy->n_in, copy->n_out, 0, 0, 1);
3777         k = isl_basic_map_alloc_inequality(copy);
3778         if (k < 0)
3779                 goto error;
3780         if (oppose)
3781                 isl_seq_neg(copy->ineq[k], c, len);
3782         else
3783                 isl_seq_cpy(copy->ineq[k], c, len);
3784         total = 1 + copy->nparam + copy->n_in + copy->n_out + copy->n_div;
3785         isl_seq_clr(copy->ineq[k]+len, total - len);
3786         isl_inequality_negate(copy, k);
3787         copy = isl_basic_map_simplify(copy);
3788         copy = isl_basic_map_finalize(copy);
3789         is_empty = isl_basic_map_is_empty(copy);
3790         if (is_empty < 0)
3791                 goto error;
3792         if (!is_empty)
3793                 dst = isl_map_add(dst, copy);
3794         else
3795                 isl_basic_map_free(copy);
3796         return dst;
3797 error:
3798         isl_basic_map_free(copy);
3799         isl_map_free(dst);
3800         return NULL;
3801 }
3802
3803 static struct isl_map *subtract(struct isl_map *map, struct isl_basic_map *bmap)
3804 {
3805         int i, j, k;
3806         unsigned flags = 0;
3807         struct isl_map *rest = NULL;
3808         unsigned max;
3809         unsigned total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
3810
3811         assert(bmap);
3812
3813         if (!map)
3814                 goto error;
3815
3816         if (F_ISSET(map, ISL_MAP_DISJOINT))
3817                 FL_SET(flags, ISL_MAP_DISJOINT);
3818
3819         max = map->n * (2 * bmap->n_eq + bmap->n_ineq);
3820         rest = isl_map_alloc(map->ctx, map->nparam, map->n_in, map->n_out,
3821                                 max, flags);
3822         if (!rest)
3823                 goto error;
3824
3825         for (i = 0; i < map->n; ++i) {
3826                 map->p[i] = isl_basic_map_align_divs(map->p[i], bmap);
3827                 if (!map->p[i])
3828                         goto error;
3829         }
3830
3831         for (j = 0; j < map->n; ++j)
3832                 map->p[j] = isl_basic_map_cow(map->p[j]);
3833
3834         for (i = 0; i < bmap->n_eq; ++i) {
3835                 for (j = 0; j < map->n; ++j) {
3836                         rest = add_cut_constraint(rest,
3837                                 map->p[j], bmap->eq[i], 1+total, 0);
3838                         if (!rest)
3839                                 goto error;
3840
3841                         rest = add_cut_constraint(rest,
3842                                 map->p[j], bmap->eq[i], 1+total, 1);
3843                         if (!rest)
3844                                 goto error;
3845
3846                         map->p[j] = isl_basic_map_extend(map->p[j],
3847                                 map->p[j]->nparam, map->p[j]->n_in,
3848                                 map->p[j]->n_out, 0, 1, 0);
3849                         if (!map->p[j])
3850                                 goto error;
3851                         k = isl_basic_map_alloc_equality(map->p[j]);
3852                         if (k < 0)
3853                                 goto error;
3854                         isl_seq_cpy(map->p[j]->eq[k], bmap->eq[i], 1+total);
3855                         isl_seq_clr(map->p[j]->eq[k]+1+total,
3856                                         map->p[j]->n_div - bmap->n_div);
3857                 }
3858         }
3859
3860         for (i = 0; i < bmap->n_ineq; ++i) {
3861                 for (j = 0; j < map->n; ++j) {
3862                         rest = add_cut_constraint(rest,
3863                                 map->p[j], bmap->ineq[i], 1+total, 0);
3864                         if (!rest)
3865                                 goto error;
3866
3867                         map->p[j] = isl_basic_map_extend(map->p[j],
3868                                 map->p[j]->nparam, map->p[j]->n_in,
3869                                 map->p[j]->n_out, 0, 0, 1);
3870                         if (!map->p[j])
3871                                 goto error;
3872                         k = isl_basic_map_alloc_inequality(map->p[j]);
3873                         if (k < 0)
3874                                 goto error;
3875                         isl_seq_cpy(map->p[j]->ineq[k], bmap->ineq[i], 1+total);
3876                         isl_seq_clr(map->p[j]->ineq[k]+1+total,
3877                                         map->p[j]->n_div - bmap->n_div);
3878                 }
3879         }
3880
3881         isl_map_free(map);
3882         return rest;
3883 error:
3884         isl_map_free(map);
3885         isl_map_free(rest);
3886         return NULL;
3887 }
3888
3889 struct isl_map *isl_map_subtract(struct isl_map *map1, struct isl_map *map2)
3890 {
3891         int i;
3892         if (!map1 || !map2)
3893                 goto error;
3894
3895         isl_assert(map1->ctx, map1->nparam == map2->nparam, goto error);
3896         isl_assert(map1->ctx, map1->n_in == map2->n_in, goto error);
3897         isl_assert(map1->ctx, map1->n_out == map2->n_out, goto error);
3898
3899         if (isl_map_is_empty(map2)) {
3900                 isl_map_free(map2);
3901                 return map1;
3902         }
3903
3904         map1 = isl_map_compute_divs(map1);
3905         map2 = isl_map_compute_divs(map2);
3906         if (!map1 || !map2)
3907                 goto error;
3908
3909         for (i = 0; map1 && i < map2->n; ++i)
3910                 map1 = subtract(map1, map2->p[i]);
3911
3912         isl_map_free(map2);
3913         return map1;
3914 error:
3915         isl_map_free(map1);
3916         isl_map_free(map2);
3917         return NULL;
3918 }
3919
3920 struct isl_set *isl_set_subtract(struct isl_set *set1, struct isl_set *set2)
3921 {
3922         return (struct isl_set *)
3923                 isl_map_subtract(
3924                         (struct isl_map *)set1, (struct isl_map *)set2);
3925 }
3926
3927 struct isl_set *isl_set_apply(struct isl_set *set, struct isl_map *map)
3928 {
3929         if (!set || !map)
3930                 goto error;
3931         isl_assert(set->ctx, set->dim == map->n_in, goto error);
3932         map = isl_map_intersect_domain(map, set);
3933         set = isl_map_range(map);
3934         return set;
3935 error:
3936         isl_set_free(set);
3937         isl_map_free(map);
3938         return NULL;
3939 }
3940
3941 /* There is no need to cow as removing empty parts doesn't change
3942  * the meaning of the set.
3943  */
3944 struct isl_map *isl_map_remove_empty_parts(struct isl_map *map)
3945 {
3946         int i;
3947
3948         if (!map)
3949                 return NULL;
3950
3951         for (i = map->n-1; i >= 0; --i) {
3952                 if (!F_ISSET(map->p[i], ISL_BASIC_MAP_EMPTY))
3953                         continue;
3954                 isl_basic_map_free(map->p[i]);
3955                 if (i != map->n-1) {
3956                         F_CLR(map, ISL_MAP_NORMALIZED);
3957                         map->p[i] = map->p[map->n-1];
3958                 }
3959                 map->n--;
3960         }
3961
3962         return map;
3963 }
3964
3965 struct isl_set *isl_set_remove_empty_parts(struct isl_set *set)
3966 {
3967         return (struct isl_set *)
3968                 isl_map_remove_empty_parts((struct isl_map *)set);
3969 }
3970
3971 struct isl_basic_set *isl_set_copy_basic_set(struct isl_set *set)
3972 {
3973         struct isl_basic_set *bset;
3974         if (!set || set->n == 0)
3975                 return NULL;
3976         bset = set->p[set->n-1];
3977         isl_assert(set->ctx, F_ISSET(bset, ISL_BASIC_SET_FINAL), return NULL);
3978         return isl_basic_set_copy(bset);
3979 }
3980
3981 struct isl_set *isl_set_drop_basic_set(struct isl_set *set,
3982                                                 struct isl_basic_set *bset)
3983 {
3984         int i;
3985
3986         if (!set || !bset)
3987                 goto error;
3988         for (i = set->n-1; i >= 0; --i) {
3989                 if (set->p[i] != bset)
3990                         continue;
3991                 set = isl_set_cow(set);
3992                 if (!set)
3993                         goto error;
3994                 isl_basic_set_free(set->p[i]);
3995                 if (i != set->n-1) {
3996                         F_CLR(set, ISL_SET_NORMALIZED);
3997                         set->p[i] = set->p[set->n-1];
3998                 }
3999                 set->n--;
4000                 return set;
4001         }
4002         isl_basic_set_free(bset);
4003         return set;
4004 error:
4005         isl_set_free(set);
4006         isl_basic_set_free(bset);
4007         return NULL;
4008 }
4009
4010 /* Given two _disjoint_ basic sets bset1 and bset2, check whether
4011  * for any common value of the parameters and dimensions preceding dim
4012  * in both basic sets, the values of dimension pos in bset1 are
4013  * smaller or larger then those in bset2.
4014  *
4015  * Returns
4016  *       1 if bset1 follows bset2
4017  *      -1 if bset1 precedes bset2
4018  *       0 if bset1 and bset2 are incomparable
4019  *      -2 if some error occurred.
4020  */
4021 int isl_basic_set_compare_at(struct isl_basic_set *bset1,
4022         struct isl_basic_set *bset2, int pos)
4023 {
4024         struct isl_basic_map *bmap1 = NULL;
4025         struct isl_basic_map *bmap2 = NULL;
4026         struct isl_ctx *ctx;
4027         struct isl_vec *obj;
4028         unsigned total;
4029         isl_int num, den;
4030         enum isl_lp_result res;
4031         int cmp;
4032
4033         if (!bset1 || !bset2)
4034                 return -2;
4035
4036         bmap1 = isl_basic_map_from_basic_set(isl_basic_set_copy(bset1),
4037                                                 pos, bset1->dim - pos);
4038         bmap2 = isl_basic_map_from_basic_set(isl_basic_set_copy(bset2),
4039                                                 pos, bset2->dim - pos);
4040         if (!bmap1 || !bmap2)
4041                 goto error;
4042         bmap1 = isl_basic_map_extend(bmap1, bmap1->nparam,
4043                         bmap1->n_in, bmap1->n_out + bmap2->n_out,
4044                         bmap2->n_div, bmap2->n_eq, bmap2->n_ineq);
4045         if (!bmap1)
4046                 goto error;
4047         total = bmap1->nparam + bmap1->n_in + bmap1->n_out + bmap1->n_div;
4048         ctx = bmap1->ctx;
4049         obj = isl_vec_alloc(ctx, total);
4050         isl_seq_clr(obj->block.data, total);
4051         isl_int_set_si(obj->block.data[bmap1->nparam+bmap1->n_in], 1);
4052         isl_int_set_si(obj->block.data[bmap1->nparam+bmap1->n_in+
4053                                         bmap1->n_out-bmap2->n_out], -1);
4054         if (!obj)
4055                 goto error;
4056         bmap1 = add_constraints(bmap1, bmap2, 0, bmap1->n_out - bmap2->n_out);
4057         isl_int_init(num);
4058         isl_int_init(den);
4059         res = isl_solve_lp(bmap1, 0, obj->block.data, ctx->one, &num, &den);
4060         if (res == isl_lp_empty)
4061                 cmp = 0;
4062         else if (res == isl_lp_ok && isl_int_is_pos(num))
4063                 cmp = 1;
4064         else if ((res == isl_lp_ok && isl_int_is_neg(num)) ||
4065                   res == isl_lp_unbounded)
4066                 cmp = -1;
4067         else
4068                 cmp = -2;
4069         isl_int_clear(num);
4070         isl_int_clear(den);
4071         isl_basic_map_free(bmap1);
4072         isl_vec_free(ctx, obj);
4073         return cmp;
4074 error:
4075         isl_basic_map_free(bmap1);
4076         isl_basic_map_free(bmap2);
4077         return -2;
4078 }
4079
4080 static int isl_basic_map_fast_has_fixed_var(struct isl_basic_map *bmap,
4081         unsigned pos, isl_int *val)
4082 {
4083         int i;
4084         int d;
4085         unsigned total;
4086
4087         if (!bmap)
4088                 return -1;
4089         total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
4090         for (i = 0, d = total-1; i < bmap->n_eq && d+1 > pos; ++i) {
4091                 for (; d+1 > pos; --d)
4092                         if (!isl_int_is_zero(bmap->eq[i][1+d]))
4093                                 break;
4094                 if (d != pos)
4095                         continue;
4096                 if (isl_seq_first_non_zero(bmap->eq[i]+1, d) != -1)
4097                         return 0;
4098                 if (isl_seq_first_non_zero(bmap->eq[i]+1+d+1, total-d-1) != -1)
4099                         return 0;
4100                 if (!isl_int_is_one(bmap->eq[i][1+d]))
4101                         return 0;
4102                 if (val)
4103                         isl_int_neg(*val, bmap->eq[i][0]);
4104                 return 1;
4105         }
4106         return 0;
4107 }
4108
4109 static int isl_map_fast_has_fixed_var(struct isl_map *map,
4110         unsigned pos, isl_int *val)
4111 {
4112         int i;
4113         isl_int v;
4114         isl_int tmp;
4115         int fixed;
4116
4117         if (!map)
4118                 return -1;
4119         if (map->n == 0)
4120                 return 0;
4121         if (map->n == 1)
4122                 return isl_basic_map_fast_has_fixed_var(map->p[0], pos, val);
4123         isl_int_init(v);
4124         isl_int_init(tmp);
4125         fixed = isl_basic_map_fast_has_fixed_var(map->p[0], pos, &v);
4126         for (i = 1; fixed == 1 && i < map->n; ++i) {
4127                 fixed = isl_basic_map_fast_has_fixed_var(map->p[i], pos, &tmp);
4128                 if (fixed == 1 && isl_int_ne(tmp, v))
4129                         fixed = 0;
4130         }
4131         if (val)
4132                 isl_int_set(*val, v);
4133         isl_int_clear(tmp);
4134         isl_int_clear(v);
4135         return fixed;
4136 }
4137
4138 static int isl_set_fast_has_fixed_var(struct isl_set *set, unsigned pos,
4139         isl_int *val)
4140 {
4141         return isl_map_fast_has_fixed_var((struct isl_map *)set, pos, val);
4142 }
4143
4144 /* Check if dimension dim has fixed value and if so and if val is not NULL,
4145  * then return this fixed value in *val.
4146  */
4147 int isl_set_fast_dim_is_fixed(struct isl_set *set, unsigned dim, isl_int *val)
4148 {
4149         return isl_set_fast_has_fixed_var(set, set->nparam + dim, val);
4150 }
4151
4152 /* Check if input variable in has fixed value and if so and if val is not NULL,
4153  * then return this fixed value in *val.
4154  */
4155 int isl_map_fast_input_is_fixed(struct isl_map *map, unsigned in, isl_int *val)
4156 {
4157         return isl_map_fast_has_fixed_var(map, map->nparam + in, val);
4158 }
4159
4160 /* Check if dimension dim has an (obvious) fixed lower bound and if so
4161  * and if val is not NULL, then return this lower bound in *val.
4162  */
4163 int isl_basic_set_fast_dim_has_fixed_lower_bound(struct isl_basic_set *bset,
4164         unsigned dim, isl_int *val)
4165 {
4166         int i, i_eq = -1, i_ineq = -1;
4167         isl_int *c;
4168         unsigned total;
4169
4170         if (!bset)
4171                 return -1;
4172         total = bset->nparam + bset->dim + bset->n_div;
4173         for (i = 0; i < bset->n_eq; ++i) {
4174                 if (isl_int_is_zero(bset->eq[i][1+bset->nparam+dim]))
4175                         continue;
4176                 if (i_eq != -1)
4177                         return 0;
4178                 i_eq = i;
4179         }
4180         for (i = 0; i < bset->n_ineq; ++i) {
4181                 if (!isl_int_is_pos(bset->ineq[i][1+bset->nparam+dim]))
4182                         continue;
4183                 if (i_eq != -1 || i_ineq != -1)
4184                         return 0;
4185                 i_ineq = i;
4186         }
4187         if (i_eq == -1 && i_ineq == -1)
4188                 return 0;
4189         c = i_eq != -1 ? bset->eq[i_eq] : bset->ineq[i_ineq];
4190         /* The coefficient should always be one due to normalization. */
4191         if (!isl_int_is_one(c[1+bset->nparam+dim]))
4192                 return 0;
4193         if (isl_seq_first_non_zero(c+1, bset->nparam+dim) != -1)
4194                 return 0;
4195         if (isl_seq_first_non_zero(c+1+bset->nparam+dim+1,
4196                                         total - bset->nparam - dim - 1) != -1)
4197                 return 0;
4198         if (val)
4199                 isl_int_neg(*val, c[0]);
4200         return 1;
4201 }
4202
4203 int isl_set_fast_dim_has_fixed_lower_bound(struct isl_set *set,
4204         unsigned dim, isl_int *val)
4205 {
4206         int i;
4207         isl_int v;
4208         isl_int tmp;
4209         int fixed;
4210
4211         if (!set)
4212                 return -1;
4213         if (set->n == 0)
4214                 return 0;
4215         if (set->n == 1)
4216                 return isl_basic_set_fast_dim_has_fixed_lower_bound(set->p[0],
4217                                                                 dim, val);
4218         isl_int_init(v);
4219         isl_int_init(tmp);
4220         fixed = isl_basic_set_fast_dim_has_fixed_lower_bound(set->p[0],
4221                                                                 dim, &v);
4222         for (i = 1; fixed == 1 && i < set->n; ++i) {
4223                 fixed = isl_basic_set_fast_dim_has_fixed_lower_bound(set->p[i],
4224                                                                 dim, &tmp);
4225                 if (fixed == 1 && isl_int_ne(tmp, v))
4226                         fixed = 0;
4227         }
4228         if (val)
4229                 isl_int_set(*val, v);
4230         isl_int_clear(tmp);
4231         isl_int_clear(v);
4232         return fixed;
4233 }
4234
4235 static struct isl_basic_set *isl_basic_set_reduce_using_equalities(
4236         struct isl_basic_set *bset, struct isl_basic_set *context)
4237 {
4238         int i;
4239         int *elim;
4240
4241         if (!bset || !context)
4242                 goto error;
4243
4244         bset = isl_basic_set_cow(bset);
4245         if (!bset)
4246                 goto error;
4247
4248         elim = isl_alloc_array(ctx, int, bset->dim);
4249         if (!elim)
4250                 goto error;
4251         set_compute_elimination_index(context, elim);
4252         for (i = 0; i < bset->n_eq; ++i)
4253                 set_reduced_using_equalities(bset->eq[i], bset->eq[i],
4254                                                         context, elim);
4255         for (i = 0; i < bset->n_ineq; ++i)
4256                 set_reduced_using_equalities(bset->ineq[i], bset->ineq[i],
4257                                                         context, elim);
4258         isl_basic_set_free(context);
4259         free(elim);
4260         bset = isl_basic_set_simplify(bset);
4261         bset = isl_basic_set_finalize(bset);
4262         return bset;
4263 error:
4264         isl_basic_set_free(bset);
4265         isl_basic_set_free(context);
4266         return NULL;
4267 }
4268
4269 static struct isl_basic_set *remove_shifted_constraints(
4270         struct isl_basic_set *bset, struct isl_basic_set *context)
4271 {
4272         unsigned int size;
4273         isl_int ***index;
4274         int bits;
4275         int k, h, l;
4276
4277         if (!bset)
4278                 return NULL;
4279
4280         size = round_up(4 * (context->n_ineq+1) / 3 - 1);
4281         bits = ffs(size) - 1;
4282         index = isl_alloc_array(ctx, isl_int **, size);
4283         memset(index, 0, size * sizeof(isl_int **));
4284         if (!index)
4285                 return bset;
4286
4287         for (k = 0; k < context->n_ineq; ++k) {
4288                 h = set_hash_index(index, size, bits, context, k);
4289                 index[h] = &context->ineq[k];
4290         }
4291         for (k = 0; k < bset->n_ineq; ++k) {
4292                 h = set_hash_index(index, size, bits, bset, k);
4293                 if (!index[h])
4294                         continue;
4295                 l = index[h] - &context->ineq[0];
4296                 if (isl_int_lt(bset->ineq[k][0], context->ineq[l][0]))
4297                         continue;
4298                 bset = isl_basic_set_cow(bset);
4299                 if (!bset)
4300                         goto error;
4301                 isl_basic_set_drop_inequality(bset, k);
4302                 --k;
4303         }
4304         free(index);
4305         return bset;
4306 error:
4307         free(index);
4308         return bset;
4309 }
4310
4311 /* Remove all information from bset that is redundant in the context
4312  * of context.  In particular, equalities that are linear combinations
4313  * of those in context are remobed.  Then the inequalities that are
4314  * redundant in the context of the equalities and inequalities of
4315  * context are removed.
4316  */
4317 static struct isl_basic_set *uset_gist(struct isl_basic_set *bset,
4318         struct isl_basic_set *context)
4319 {
4320         int i;
4321         isl_int opt;
4322         struct isl_basic_set *combined;
4323         struct isl_ctx *ctx;
4324
4325         if (!bset || !context)
4326                 goto error;
4327
4328         if (context->n_eq > 0) {
4329                 struct isl_mat *T;
4330                 struct isl_mat *T2;
4331                 struct isl_ctx *ctx = context->ctx;
4332                 struct isl_basic_set *reduced_context;
4333                 reduced_context = isl_basic_set_remove_equalities(
4334                                         isl_basic_set_copy(context), &T, &T2);
4335                 if (!reduced_context)
4336                         goto error;
4337                 bset = isl_basic_set_preimage(ctx, bset, T);
4338                 bset = uset_gist(bset, reduced_context);
4339                 bset = isl_basic_set_preimage(ctx, bset, T2);
4340                 bset = isl_basic_set_reduce_using_equalities(bset, context);
4341                 return bset;
4342         }
4343         if (!context->n_ineq)
4344                 goto done;
4345         bset = remove_shifted_constraints(bset, context);
4346         combined = isl_basic_set_extend(isl_basic_set_copy(bset),
4347                         0, bset->dim, 0, context->n_eq, context->n_ineq);
4348         context = set_add_constraints(combined, context, 0);
4349         if (!context)
4350                 goto error;
4351         ctx = context->ctx;
4352         isl_int_init(opt);
4353         for (i = bset->n_ineq-1; i >= 0; --i) {
4354                 int redundant;
4355                 set_swap_inequality(context, i, context->n_ineq-1);
4356                 context->n_ineq--;
4357                 redundant = isl_basic_set_constraint_is_redundant(&context,
4358                                 context->ineq[context->n_ineq], &opt, NULL);
4359                 if (redundant == -1) {
4360                         isl_int_clear(opt);
4361                         goto error;
4362                 }
4363                 if (F_ISSET(context, ISL_BASIC_MAP_EMPTY)) {
4364                         bset = isl_basic_set_set_to_empty(bset);
4365                         break;
4366                 }
4367                 context->n_ineq++;
4368                 set_swap_inequality(context, i, context->n_ineq-1);
4369                 if (redundant) {
4370                         isl_basic_set_drop_inequality(context, i);
4371                         isl_basic_set_drop_inequality(bset, i);
4372                 }
4373         }
4374         isl_int_clear(opt);
4375 done:
4376         isl_basic_set_free(context);
4377         return bset;
4378 error:
4379         isl_basic_set_free(context);
4380         isl_basic_set_free(bset);
4381         return NULL;
4382 }
4383
4384 struct isl_basic_map *isl_basic_map_gist(struct isl_basic_map *bmap,
4385         struct isl_basic_map *context)
4386 {
4387         struct isl_basic_set *bset;
4388
4389         if (!bmap || !context)
4390                 goto error;
4391
4392         context = isl_basic_map_align_divs(context, bmap);
4393         bmap = isl_basic_map_align_divs(bmap, context);
4394
4395         bset = uset_gist(isl_basic_map_underlying_set(isl_basic_map_copy(bmap)),
4396                          isl_basic_map_underlying_set(context));
4397
4398         return isl_basic_map_overlying_set(bset, bmap);
4399 error:
4400         isl_basic_map_free(bmap);
4401         isl_basic_map_free(context);
4402         return NULL;
4403 }
4404
4405 struct isl_map *isl_map_gist(struct isl_map *map, struct isl_basic_map *context)
4406 {
4407         int i;
4408
4409         map = isl_map_cow(map);
4410         if (!map || !context)
4411                 return NULL;
4412         isl_assert(map->ctx, map->nparam == context->nparam, goto error);
4413         isl_assert(map->ctx, map->n_in == context->n_in, goto error);
4414         isl_assert(map->ctx, map->n_out == context->n_out, goto error);
4415         for (i = 0; i < map->n; ++i)
4416                 context = isl_basic_map_align_divs(context, map->p[i]);
4417         for (i = 0; i < map->n; ++i) {
4418                 map->p[i] = isl_basic_map_gist(map->p[i],
4419                                                 isl_basic_map_copy(context));
4420                 if (!map->p[i])
4421                         goto error;
4422         }
4423         isl_basic_map_free(context);
4424         F_CLR(map, ISL_MAP_NORMALIZED);
4425         return map;
4426 error:
4427         isl_map_free(map);
4428         isl_basic_map_free(context);
4429         return NULL;
4430 }
4431
4432 struct isl_set *isl_set_gist(struct isl_set *set, struct isl_basic_set *context)
4433 {
4434         return (struct isl_set *)isl_map_gist((struct isl_map *)set,
4435                                         (struct isl_basic_map *)context);
4436 }
4437
4438 struct constraint {
4439         unsigned        size;
4440         isl_int         *c;
4441 };
4442
4443 static int qsort_constraint_cmp(const void *p1, const void *p2)
4444 {
4445         const struct constraint *c1 = (const struct constraint *)p1;
4446         const struct constraint *c2 = (const struct constraint *)p2;
4447         unsigned size = isl_min(c1->size, c2->size);
4448         return isl_seq_cmp(c1->c, c2->c, size);
4449 }
4450
4451 static struct isl_basic_map *isl_basic_map_sort_constraints(
4452         struct isl_basic_map *bmap)
4453 {
4454         int i;
4455         struct constraint *c;
4456         unsigned total;
4457
4458         if (!bmap)
4459                 return NULL;
4460         total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
4461         c = isl_alloc_array(bmap->ctx, struct constraint, bmap->n_ineq);
4462         if (!c)
4463                 goto error;
4464         for (i = 0; i < bmap->n_ineq; ++i) {
4465                 c[i].size = total;
4466                 c[i].c = bmap->ineq[i];
4467         }
4468         qsort(c, bmap->n_ineq, sizeof(struct constraint), qsort_constraint_cmp);
4469         for (i = 0; i < bmap->n_ineq; ++i)
4470                 bmap->ineq[i] = c[i].c;
4471         free(c);
4472         return bmap;
4473 error:
4474         isl_basic_map_free(bmap);
4475         return NULL;
4476 }
4477
4478 struct isl_basic_map *isl_basic_map_normalize(struct isl_basic_map *bmap)
4479 {
4480         if (!bmap)
4481                 return NULL;
4482         if (F_ISSET(bmap, ISL_BASIC_MAP_NORMALIZED))
4483                 return bmap;
4484         bmap = isl_basic_map_convex_hull(bmap);
4485         bmap = isl_basic_map_sort_constraints(bmap);
4486         F_SET(bmap, ISL_BASIC_MAP_NORMALIZED);
4487         return bmap;
4488 }
4489
4490 static int isl_basic_map_fast_cmp(const struct isl_basic_map *bmap1,
4491         const struct isl_basic_map *bmap2)
4492 {
4493         int i, cmp;
4494         unsigned total;
4495
4496         if (bmap1 == bmap2)
4497                 return 0;
4498         if (bmap1->nparam != bmap2->nparam)
4499                 return bmap1->nparam - bmap2->nparam;
4500         if (bmap1->n_in != bmap2->n_in)
4501                 return bmap1->n_in - bmap2->n_in;
4502         if (bmap1->n_out != bmap2->n_out)
4503                 return bmap1->n_out - bmap2->n_out;
4504         if (F_ISSET(bmap1, ISL_BASIC_MAP_EMPTY) &&
4505             F_ISSET(bmap2, ISL_BASIC_MAP_EMPTY))
4506                 return 0;
4507         if (F_ISSET(bmap1, ISL_BASIC_MAP_EMPTY))
4508                 return 1;
4509         if (F_ISSET(bmap2, ISL_BASIC_MAP_EMPTY))
4510                 return -1;
4511         if (bmap1->n_eq != bmap2->n_eq)
4512                 return bmap1->n_eq - bmap2->n_eq;
4513         if (bmap1->n_ineq != bmap2->n_ineq)
4514                 return bmap1->n_ineq - bmap2->n_ineq;
4515         if (bmap1->n_div != bmap2->n_div)
4516                 return bmap1->n_div - bmap2->n_div;
4517         total = bmap1->nparam + bmap1->n_in + bmap1->n_out + bmap1->n_div;
4518         for (i = 0; i < bmap1->n_eq; ++i) {
4519                 cmp = isl_seq_cmp(bmap1->eq[i], bmap2->eq[i], 1+total);
4520                 if (cmp)
4521                         return cmp;
4522         }
4523         for (i = 0; i < bmap1->n_ineq; ++i) {
4524                 cmp = isl_seq_cmp(bmap1->ineq[i], bmap2->ineq[i], 1+total);
4525                 if (cmp)
4526                         return cmp;
4527         }
4528         for (i = 0; i < bmap1->n_div; ++i) {
4529                 cmp = isl_seq_cmp(bmap1->div[i], bmap2->div[i], 1+1+total);
4530                 if (cmp)
4531                         return cmp;
4532         }
4533         return 0;
4534 }
4535
4536 static int isl_basic_map_fast_is_equal(struct isl_basic_map *bmap1,
4537         struct isl_basic_map *bmap2)
4538 {
4539         return isl_basic_map_fast_cmp(bmap1, bmap2) == 0;
4540 }
4541
4542 static int qsort_bmap_cmp(const void *p1, const void *p2)
4543 {
4544         const struct isl_basic_map *bmap1 = *(const struct isl_basic_map **)p1;
4545         const struct isl_basic_map *bmap2 = *(const struct isl_basic_map **)p2;
4546
4547         return isl_basic_map_fast_cmp(bmap1, bmap2);
4548 }
4549
4550 /* We normalize in place, but if anything goes wrong we need
4551  * to return NULL, so we need to make sure we don't change the
4552  * meaning of any possible other copies of map.
4553  */
4554 struct isl_map *isl_map_normalize(struct isl_map *map)
4555 {
4556         int i, j;
4557         struct isl_basic_map *bmap;
4558
4559         if (!map)
4560                 return NULL;
4561         if (F_ISSET(map, ISL_MAP_NORMALIZED))
4562                 return map;
4563         for (i = 0; i < map->n; ++i) {
4564                 bmap = isl_basic_map_normalize(isl_basic_map_copy(map->p[i]));
4565                 if (!bmap)
4566                         goto error;
4567                 isl_basic_map_free(map->p[i]);
4568                 map->p[i] = bmap;
4569         }
4570         qsort(map->p, map->n, sizeof(struct isl_basic_map *), qsort_bmap_cmp);
4571         F_SET(map, ISL_MAP_NORMALIZED);
4572         map = isl_map_remove_empty_parts(map);
4573         if (!map)
4574                 return NULL;
4575         for (i = map->n - 1; i >= 1; --i) {
4576                 if (!isl_basic_map_fast_is_equal(map->p[i-1], map->p[i]))
4577                         continue;
4578                 isl_basic_map_free(map->p[i-1]);
4579                 for (j = i; j < map->n; ++j)
4580                         map->p[j-1] = map->p[j];
4581                 map->n--;
4582         }
4583         return map;
4584 error:
4585         isl_map_free(map);
4586         return NULL;
4587
4588 }
4589
4590 struct isl_set *isl_set_normalize(struct isl_set *set)
4591 {
4592         return (struct isl_set *)isl_map_normalize((struct isl_map *)set);
4593 }
4594
4595 int isl_map_fast_is_equal(struct isl_map *map1, struct isl_map *map2)
4596 {
4597         int i;
4598         int equal;
4599
4600         if (!map1 || !map2)
4601                 return -1;
4602
4603         if (map1 == map2)
4604                 return 1;
4605         if (map1->nparam != map2->nparam)
4606                 return 0;
4607         if (map1->n_in != map2->n_in)
4608                 return 0;
4609         if (map1->n_out != map2->n_out)
4610                 return 0;
4611
4612         map1 = isl_map_copy(map1);
4613         map2 = isl_map_copy(map2);
4614         map1 = isl_map_normalize(map1);
4615         map2 = isl_map_normalize(map2);
4616         if (!map1 || !map2)
4617                 goto error;
4618         equal = map1->n == map2->n;
4619         for (i = 0; equal && i < map1->n; ++i) {
4620                 equal = isl_basic_map_fast_is_equal(map1->p[i], map2->p[i]);
4621                 if (equal < 0)
4622                         goto error;
4623         }
4624         isl_map_free(map1);
4625         isl_map_free(map2);
4626         return equal;
4627 error:
4628         isl_map_free(map1);
4629         isl_map_free(map2);
4630         return -1;
4631 }
4632
4633 int isl_set_fast_is_equal(struct isl_set *set1, struct isl_set *set2)
4634 {
4635         return isl_map_fast_is_equal((struct isl_map *)set1,
4636                                                 (struct isl_map *)set2);
4637 }
4638
4639 /* Return an interval that ranges from min to max (inclusive)
4640  */
4641 struct isl_basic_set *isl_basic_set_interval(struct isl_ctx *ctx,
4642         isl_int min, isl_int max)
4643 {
4644         int k;
4645         struct isl_basic_set *bset = NULL;
4646
4647         bset = isl_basic_set_alloc(ctx, 0, 1, 0, 0, 2);
4648         if (!bset)
4649                 goto error;
4650
4651         k = isl_basic_set_alloc_inequality(bset);
4652         if (k < 0)
4653                 goto error;
4654         isl_int_set_si(bset->ineq[k][1], 1);
4655         isl_int_neg(bset->ineq[k][0], min);
4656
4657         k = isl_basic_set_alloc_inequality(bset);
4658         if (k < 0)
4659                 goto error;
4660         isl_int_set_si(bset->ineq[k][1], -1);
4661         isl_int_set(bset->ineq[k][0], max);
4662
4663         return bset;
4664 error:
4665         isl_basic_set_free(bset);
4666         return NULL;
4667 }
4668
4669 /* Return the Cartesian product of the basic sets in list (in the given order).
4670  */
4671 struct isl_basic_set *isl_basic_set_product(struct isl_basic_set_list *list)
4672 {
4673         int i;
4674         unsigned dim;
4675         unsigned nparam;
4676         unsigned extra;
4677         unsigned n_eq;
4678         unsigned n_ineq;
4679         struct isl_basic_set *product = NULL;
4680
4681         if (!list)
4682                 goto error;
4683         isl_assert(list->ctx, list->n > 0, goto error);
4684         isl_assert(list->ctx, list->p[0], goto error);
4685         nparam = list->p[0]->nparam;
4686         dim = list->p[0]->dim;
4687         extra = list->p[0]->n_div;
4688         n_eq = list->p[0]->n_eq;
4689         n_ineq = list->p[0]->n_ineq;
4690         for (i = 1; i < list->n; ++i) {
4691                 isl_assert(list->ctx, list->p[i], goto error);
4692                 isl_assert(list->ctx, nparam == list->p[i]->nparam, goto error);
4693                 dim += list->p[i]->dim;
4694                 extra += list->p[i]->n_div;
4695                 n_eq += list->p[i]->n_eq;
4696                 n_ineq += list->p[i]->n_ineq;
4697         }
4698         product = isl_basic_set_alloc(list->ctx, nparam, dim, extra,
4699                                         n_eq, n_ineq);
4700         if (!product)
4701                 goto error;
4702         dim = 0;
4703         for (i = 0; i < list->n; ++i) {
4704                 set_add_constraints(product,
4705                                         isl_basic_set_copy(list->p[i]), dim);
4706                 dim += list->p[i]->dim;
4707         }
4708         isl_basic_set_list_free(list);
4709         return product;
4710 error:
4711         isl_basic_set_free(product);
4712         isl_basic_set_list_free(list);
4713         return NULL;
4714 }