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