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
1434         bmap = isl_basic_map_cow(bmap);
1435         for (d = pos + n - 1; d >= pos; --d) {
1436                 int n_lower, n_upper;
1437                 if (!bmap)
1438                         return NULL;
1439                 for (i = 0; i < bmap->n_eq; ++i) {
1440                         if (isl_int_is_zero(bmap->eq[i][1+d]))
1441                                 continue;
1442                         eliminate_var_using_equality(bmap, d, bmap->eq[i], NULL);
1443                         isl_basic_map_drop_equality(bmap, i);
1444                         break;
1445                 }
1446                 if (i < bmap->n_eq)
1447                         continue;
1448                 n_lower = 0;
1449                 n_upper = 0;
1450                 for (i = 0; i < bmap->n_ineq; ++i) {
1451                         if (isl_int_is_pos(bmap->ineq[i][1+d]))
1452                                 n_lower++;
1453                         else if (isl_int_is_neg(bmap->ineq[i][1+d]))
1454                                 n_upper++;
1455                 }
1456                 bmap = isl_basic_map_extend(bmap,
1457                                 bmap->nparam, bmap->n_in, bmap->n_out, 0,
1458                                 0, n_lower * n_upper);
1459                 extra = bmap->extra - bmap->n_div;
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, isl_basic_set_copy(set->p[i]));
1829         return dup;
1830 }
1831
1832 struct isl_set *isl_set_from_basic_set(struct isl_basic_set *bset)
1833 {
1834         struct isl_set *set;
1835
1836         if (!bset)
1837                 return NULL;
1838
1839         set = isl_set_alloc(bset->ctx, bset->nparam, bset->dim, 1, ISL_MAP_DISJOINT);
1840         if (!set) {
1841                 isl_basic_set_free(bset);
1842                 return NULL;
1843         }
1844         return isl_set_add(set, bset);
1845 }
1846
1847 struct isl_map *isl_map_from_basic_map(struct isl_basic_map *bmap)
1848 {
1849         struct isl_map *map;
1850
1851         if (!bmap)
1852                 return NULL;
1853
1854         map = isl_map_alloc(bmap->ctx, bmap->nparam, bmap->n_in, bmap->n_out, 1,
1855                                 ISL_MAP_DISJOINT);
1856         if (!map) {
1857                 isl_basic_map_free(bmap);
1858                 return NULL;
1859         }
1860         return isl_map_add(map, bmap);
1861 }
1862
1863 struct isl_set *isl_set_add(struct isl_set *set, struct isl_basic_set *bset)
1864 {
1865         if (!bset || !set)
1866                 goto error;
1867         isl_assert(set->ctx, set->nparam == bset->nparam, goto error);
1868         isl_assert(set->ctx, set->dim == bset->dim, goto error);
1869         isl_assert(set->ctx, set->n < set->size, goto error);
1870         set->p[set->n] = bset;
1871         set->n++;
1872         return set;
1873 error:
1874         if (set)
1875                 isl_set_free(set);
1876         if (bset)
1877                 isl_basic_set_free(bset);
1878         return NULL;
1879 }
1880
1881 void isl_set_free(struct isl_set *set)
1882 {
1883         int i;
1884
1885         if (!set)
1886                 return;
1887
1888         if (--set->ref > 0)
1889                 return;
1890
1891         isl_ctx_deref(set->ctx);
1892         for (i = 0; i < set->n; ++i)
1893                 isl_basic_set_free(set->p[i]);
1894         free(set);
1895 }
1896
1897 void isl_set_dump(struct isl_set *set, FILE *out, int indent)
1898 {
1899         int i;
1900
1901         if (!set) {
1902                 fprintf(out, "null set\n");
1903                 return;
1904         }
1905
1906         fprintf(out, "%*s", indent, "");
1907         fprintf(out, "ref: %d, n: %d, nparam: %d, dim: %d, flags: %x\n",
1908                         set->ref, set->n, set->nparam, set->dim, set->flags);
1909         for (i = 0; i < set->n; ++i) {
1910                 fprintf(out, "%*s", indent, "");
1911                 fprintf(out, "basic set %d:\n", i);
1912                 isl_basic_set_dump(set->p[i], out, indent+4);
1913         }
1914 }
1915
1916 void isl_map_dump(struct isl_map *map, FILE *out, int indent)
1917 {
1918         int i;
1919
1920         if (!map) {
1921                 fprintf(out, "null map\n");
1922                 return;
1923         }
1924
1925         fprintf(out, "%*s", indent, "");
1926         fprintf(out, "ref: %d, n: %d, nparam: %d, in: %d, out: %d, flags: %x\n",
1927                         map->ref, map->n, map->nparam, map->n_in, map->n_out,
1928                         map->flags);
1929         for (i = 0; i < map->n; ++i) {
1930                 fprintf(out, "%*s", indent, "");
1931                 fprintf(out, "basic map %d:\n", i);
1932                 isl_basic_map_dump(map->p[i], out, indent+4);
1933         }
1934 }
1935
1936 struct isl_basic_map *isl_basic_map_intersect_domain(
1937                 struct isl_basic_map *bmap, struct isl_basic_set *bset)
1938 {
1939         struct isl_basic_map *bmap_domain;
1940
1941         if (!bmap || !bset)
1942                 goto error;
1943
1944         isl_assert(bset->ctx, bset->dim == bmap->n_in, goto error);
1945         isl_assert(bset->ctx, bset->nparam == bmap->nparam, goto error);
1946
1947         bmap = isl_basic_map_extend(bmap,
1948                         bset->nparam, bmap->n_in, bmap->n_out,
1949                         bset->n_div, bset->n_eq, bset->n_ineq);
1950         if (!bmap)
1951                 goto error;
1952         bmap_domain = isl_basic_map_from_basic_set(bset,
1953                                                 bset->dim, 0);
1954         bmap = add_constraints(bmap, bmap_domain, 0, 0);
1955
1956         bmap = isl_basic_map_simplify(bmap);
1957         return isl_basic_map_finalize(bmap);
1958 error:
1959         isl_basic_map_free(bmap);
1960         isl_basic_set_free(bset);
1961         return NULL;
1962 }
1963
1964 struct isl_basic_map *isl_basic_map_intersect_range(
1965                 struct isl_basic_map *bmap, struct isl_basic_set *bset)
1966 {
1967         struct isl_basic_map *bmap_range;
1968
1969         if (!bmap || !bset)
1970                 goto error;
1971
1972         isl_assert(bset->ctx, bset->dim == bmap->n_out, goto error);
1973         isl_assert(bset->ctx, bset->nparam == bmap->nparam, goto error);
1974
1975         bmap = isl_basic_map_extend(bmap,
1976                         bset->nparam, bmap->n_in, bmap->n_out,
1977                         bset->n_div, bset->n_eq, bset->n_ineq);
1978         if (!bmap)
1979                 goto error;
1980         bmap_range = isl_basic_map_from_basic_set(bset,
1981                                                 0, bset->dim);
1982         bmap = add_constraints(bmap, bmap_range, 0, 0);
1983
1984         bmap = isl_basic_map_simplify(bmap);
1985         return isl_basic_map_finalize(bmap);
1986 error:
1987         isl_basic_map_free(bmap);
1988         isl_basic_set_free(bset);
1989         return NULL;
1990 }
1991
1992 struct isl_basic_map *isl_basic_map_intersect(
1993                 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
1994 {
1995         if (!bmap1 || !bmap2)
1996                 goto error;
1997
1998         isl_assert(bmap1->ctx, bmap1->nparam == bmap2->nparam, goto error);
1999         isl_assert(bmap1->ctx, bmap1->n_in == bmap2->n_in, goto error);
2000         isl_assert(bmap1->ctx, bmap1->n_out == bmap2->n_out, goto error);
2001
2002         bmap1 = isl_basic_map_extend(bmap1,
2003                         bmap1->nparam, bmap1->n_in, bmap1->n_out,
2004                         bmap2->n_div, bmap2->n_eq, bmap2->n_ineq);
2005         if (!bmap1)
2006                 goto error;
2007         bmap1 = add_constraints(bmap1, bmap2, 0, 0);
2008
2009         bmap1 = isl_basic_map_simplify(bmap1);
2010         return isl_basic_map_finalize(bmap1);
2011 error:
2012         isl_basic_map_free(bmap1);
2013         isl_basic_map_free(bmap2);
2014         return NULL;
2015 }
2016
2017 struct isl_basic_set *isl_basic_set_intersect(
2018                 struct isl_basic_set *bset1, struct isl_basic_set *bset2)
2019 {
2020         return (struct isl_basic_set *)
2021                 isl_basic_map_intersect(
2022                         (struct isl_basic_map *)bset1,
2023                         (struct isl_basic_map *)bset2);
2024 }
2025
2026 struct isl_map *isl_map_intersect(struct isl_map *map1, struct isl_map *map2)
2027 {
2028         unsigned flags = 0;
2029         struct isl_map *result;
2030         int i, j;
2031
2032         if (!map1 || !map2)
2033                 goto error;
2034
2035         if (F_ISSET(map1, ISL_MAP_DISJOINT) &&
2036             F_ISSET(map2, ISL_MAP_DISJOINT))
2037                 FL_SET(flags, ISL_MAP_DISJOINT);
2038
2039         result = isl_map_alloc(map1->ctx, map1->nparam, map1->n_in, map1->n_out,
2040                                 map1->n * map2->n, flags);
2041         if (!result)
2042                 goto error;
2043         for (i = 0; i < map1->n; ++i)
2044                 for (j = 0; j < map2->n; ++j) {
2045                         struct isl_basic_map *part;
2046                         part = isl_basic_map_intersect(
2047                                     isl_basic_map_copy(map1->p[i]),
2048                                     isl_basic_map_copy(map2->p[j]));
2049                         if (isl_basic_map_is_empty(part))
2050                                 isl_basic_map_free(part);
2051                         else
2052                                 result = isl_map_add(result, part);
2053                         if (!result)
2054                                 goto error;
2055                 }
2056         isl_map_free(map1);
2057         isl_map_free(map2);
2058         return result;
2059 error:
2060         isl_map_free(map1);
2061         isl_map_free(map2);
2062         return NULL;
2063 }
2064
2065 struct isl_set *isl_set_intersect(struct isl_set *set1, struct isl_set *set2)
2066 {
2067         return (struct isl_set *)
2068                 isl_map_intersect((struct isl_map *)set1,
2069                                   (struct isl_map *)set2);
2070 }
2071
2072 struct isl_basic_map *isl_basic_map_reverse(struct isl_basic_map *bmap)
2073 {
2074         struct isl_basic_set *bset;
2075         unsigned in;
2076
2077         if (!bmap)
2078                 return NULL;
2079         bmap = isl_basic_map_cow(bmap);
2080         if (!bmap)
2081                 return NULL;
2082         in = bmap->n_in;
2083         bset = isl_basic_set_from_basic_map(bmap);
2084         bset = isl_basic_set_swap_vars(bset, in);
2085         return isl_basic_map_from_basic_set(bset, bset->dim-in, in);
2086 }
2087
2088 /* Turn final n dimensions into existentially quantified variables.
2089  */
2090 struct isl_basic_set *isl_basic_set_project_out(
2091                 struct isl_basic_set *bset, unsigned n, unsigned flags)
2092 {
2093         int i;
2094         size_t row_size;
2095         isl_int **new_div;
2096         isl_int *old;
2097
2098         if (!bset)
2099                 return NULL;
2100
2101         isl_assert(bset->ctx, n <= bset->dim, goto error);
2102
2103         if (n == 0)
2104                 return bset;
2105
2106         bset = isl_basic_set_cow(bset);
2107
2108         row_size = 1 + bset->nparam + bset->dim + bset->extra;
2109         old = bset->block2.data;
2110         bset->block2 = isl_blk_extend(bset->ctx, bset->block2,
2111                                         (bset->extra + n) * (1 + row_size));
2112         if (!bset->block2.data)
2113                 goto error;
2114         new_div = isl_alloc_array(ctx, isl_int *, bset->extra + n);
2115         if (!new_div)
2116                 goto error;
2117         for (i = 0; i < n; ++i) {
2118                 new_div[i] = bset->block2.data +
2119                                 (bset->extra + i) * (1 + row_size);
2120                 isl_seq_clr(new_div[i], 1 + row_size);
2121         }
2122         for (i = 0; i < bset->extra; ++i)
2123                 new_div[n + i] = bset->block2.data + (bset->div[i] - old);
2124         free(bset->div);
2125         bset->div = new_div;
2126         bset->n_div += n;
2127         bset->extra += n;
2128         bset->dim -= n;
2129         bset = isl_basic_set_simplify(bset);
2130         return isl_basic_set_finalize(bset);
2131 error:
2132         isl_basic_set_free(bset);
2133         return NULL;
2134 }
2135
2136 struct isl_basic_map *isl_basic_map_apply_range(
2137                 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
2138 {
2139         struct isl_basic_set *bset;
2140         unsigned n_in, n_out;
2141
2142         if (!bmap1 || !bmap2)
2143                 goto error;
2144
2145         isl_assert(bmap->ctx, bmap1->n_out == bmap2->n_in, goto error);
2146         isl_assert(bmap->ctx, bmap1->nparam == bmap2->nparam, goto error);
2147
2148         n_in = bmap1->n_in;
2149         n_out = bmap2->n_out;
2150
2151         bmap2 = isl_basic_map_reverse(bmap2);
2152         if (!bmap2)
2153                 goto error;
2154         bmap1 = isl_basic_map_extend(bmap1, bmap1->nparam,
2155                         bmap1->n_in + bmap2->n_in, bmap1->n_out,
2156                         bmap2->extra, bmap2->n_eq, bmap2->n_ineq);
2157         if (!bmap1)
2158                 goto error;
2159         bmap1 = add_constraints(bmap1, bmap2, bmap1->n_in - bmap2->n_in, 0);
2160         bmap1 = isl_basic_map_simplify(bmap1);
2161         bset = isl_basic_set_from_basic_map(bmap1);
2162         bset = isl_basic_set_project_out(bset,
2163                                                 bset->dim - (n_in + n_out), 0);
2164         return isl_basic_map_from_basic_set(bset, n_in, n_out);
2165 error:
2166         isl_basic_map_free(bmap1);
2167         isl_basic_map_free(bmap2);
2168         return NULL;
2169 }
2170
2171 struct isl_basic_set *isl_basic_set_apply(
2172                 struct isl_basic_set *bset, struct isl_basic_map *bmap)
2173 {
2174         if (!bset || !bmap)
2175                 goto error;
2176
2177         isl_assert(bset->ctx, bset->dim == bmap->n_in, goto error);
2178
2179         return (struct isl_basic_set *)
2180                 isl_basic_map_apply_range((struct isl_basic_map *)bset, bmap);
2181 error:
2182         isl_basic_set_free(bset);
2183         isl_basic_map_free(bmap);
2184         return NULL;
2185 }
2186
2187 struct isl_basic_map *isl_basic_map_apply_domain(
2188                 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
2189 {
2190         if (!bmap1 || !bmap2)
2191                 goto error;
2192
2193         isl_assert(ctx, bmap1->n_in == bmap2->n_in, goto error);
2194         isl_assert(ctx, bmap1->nparam == bmap2->nparam, goto error);
2195
2196         bmap1 = isl_basic_map_reverse(bmap1);
2197         bmap1 = isl_basic_map_apply_range(bmap1, bmap2);
2198         return isl_basic_map_reverse(bmap1);
2199 error:
2200         isl_basic_map_free(bmap1);
2201         isl_basic_map_free(bmap2);
2202         return NULL;
2203 }
2204
2205 static struct isl_basic_map *var_equal(struct isl_ctx *ctx,
2206                 struct isl_basic_map *bmap, unsigned pos)
2207 {
2208         int i;
2209         i = isl_basic_map_alloc_equality(bmap);
2210         if (i < 0)
2211                 goto error;
2212         isl_seq_clr(bmap->eq[i],
2213                     1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->extra);
2214         isl_int_set_si(bmap->eq[i][1+bmap->nparam+pos], -1);
2215         isl_int_set_si(bmap->eq[i][1+bmap->nparam+bmap->n_in+pos], 1);
2216         return isl_basic_map_finalize(bmap);
2217 error:
2218         isl_basic_map_free(bmap);
2219         return NULL;
2220 }
2221
2222 static struct isl_basic_map *var_more(struct isl_ctx *ctx,
2223                 struct isl_basic_map *bmap, unsigned pos)
2224 {
2225         int i;
2226         i = isl_basic_map_alloc_inequality(bmap);
2227         if (i < 0)
2228                 goto error;
2229         isl_seq_clr(bmap->ineq[i],
2230                     1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->extra);
2231         isl_int_set_si(bmap->ineq[i][0], -1);
2232         isl_int_set_si(bmap->ineq[i][1+bmap->nparam+pos], -1);
2233         isl_int_set_si(bmap->ineq[i][1+bmap->nparam+bmap->n_in+pos], 1);
2234         return isl_basic_map_finalize(bmap);
2235 error:
2236         isl_basic_map_free(bmap);
2237         return NULL;
2238 }
2239
2240 static struct isl_basic_map *var_less(struct isl_ctx *ctx,
2241                 struct isl_basic_map *bmap, unsigned pos)
2242 {
2243         int i;
2244         i = isl_basic_map_alloc_inequality(bmap);
2245         if (i < 0)
2246                 goto error;
2247         isl_seq_clr(bmap->ineq[i],
2248                     1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->extra);
2249         isl_int_set_si(bmap->ineq[i][0], -1);
2250         isl_int_set_si(bmap->ineq[i][1+bmap->nparam+pos], 1);
2251         isl_int_set_si(bmap->ineq[i][1+bmap->nparam+bmap->n_in+pos], -1);
2252         return isl_basic_map_finalize(bmap);
2253 error:
2254         isl_basic_map_free(bmap);
2255         return NULL;
2256 }
2257
2258 struct isl_basic_map *isl_basic_map_equal(struct isl_ctx *ctx,
2259                 unsigned nparam, unsigned in, unsigned out, unsigned n_equal)
2260 {
2261         int i;
2262         struct isl_basic_map *bmap;
2263         bmap = isl_basic_map_alloc(ctx, nparam, in, out, 0, n_equal, 0);
2264         if (!bmap)
2265                 return NULL;
2266         for (i = 0; i < n_equal && bmap; ++i)
2267                 bmap = var_equal(ctx, bmap, i);
2268         return isl_basic_map_finalize(bmap);
2269 }
2270
2271 struct isl_basic_map *isl_basic_map_less_at(struct isl_ctx *ctx,
2272                 unsigned nparam, unsigned in, unsigned out, unsigned pos)
2273 {
2274         int i;
2275         struct isl_basic_map *bmap;
2276         bmap = isl_basic_map_alloc(ctx, nparam, in, out, 0, pos, 1);
2277         if (!bmap)
2278                 return NULL;
2279         for (i = 0; i < pos && bmap; ++i)
2280                 bmap = var_equal(ctx, bmap, i);
2281         if (bmap)
2282                 bmap = var_less(ctx, bmap, pos);
2283         return isl_basic_map_finalize(bmap);
2284 }
2285
2286 struct isl_basic_map *isl_basic_map_more_at(struct isl_ctx *ctx,
2287                 unsigned nparam, unsigned in, unsigned out, unsigned pos)
2288 {
2289         int i;
2290         struct isl_basic_map *bmap;
2291         bmap = isl_basic_map_alloc(ctx, nparam, in, out, 0, pos, 1);
2292         if (!bmap)
2293                 return NULL;
2294         for (i = 0; i < pos && bmap; ++i)
2295                 bmap = var_equal(ctx, bmap, i);
2296         if (bmap)
2297                 bmap = var_more(ctx, bmap, pos);
2298         return isl_basic_map_finalize(bmap);
2299 }
2300
2301 struct isl_basic_map *isl_basic_map_from_basic_set(
2302                 struct isl_basic_set *bset, unsigned n_in, unsigned n_out)
2303 {
2304         struct isl_basic_map *bmap;
2305
2306         bset = isl_basic_set_cow(bset);
2307         if (!bset)
2308                 return NULL;
2309
2310         isl_assert(ctx, bset->dim == n_in + n_out, goto error);
2311         bmap = (struct isl_basic_map *) bset;
2312         bmap->n_in = n_in;
2313         bmap->n_out = n_out;
2314         return isl_basic_map_finalize(bmap);
2315 error:
2316         isl_basic_set_free(bset);
2317         return NULL;
2318 }
2319
2320 struct isl_basic_set *isl_basic_set_from_basic_map(struct isl_basic_map *bmap)
2321 {
2322         if (!bmap)
2323                 goto error;
2324         if (bmap->n_in == 0)
2325                 return (struct isl_basic_set *)bmap;
2326         bmap = isl_basic_map_cow(bmap);
2327         if (!bmap)
2328                 goto error;
2329         bmap->n_out += bmap->n_in;
2330         bmap->n_in = 0;
2331         bmap = isl_basic_map_finalize(bmap);
2332         return (struct isl_basic_set *)bmap;
2333 error:
2334         return NULL;
2335 }
2336
2337 /* For a div d = floor(f/m), add the constraints
2338  *
2339  *              f - m d >= 0
2340  *              -(f-(n-1)) + m d >= 0
2341  *
2342  * Note that the second constraint is the negation of
2343  *
2344  *              f - m d >= n
2345  */
2346 static int add_div_constraints(struct isl_basic_map *bmap, unsigned div)
2347 {
2348         int i, j;
2349         unsigned div_pos = 1 + bmap->nparam + bmap->n_in + bmap->n_out + div;
2350         unsigned total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->extra;
2351
2352         i = isl_basic_map_alloc_inequality(bmap);
2353         if (i < 0)
2354                 return -1;
2355         isl_seq_cpy(bmap->ineq[i], bmap->div[div]+1, 1+total);
2356         isl_int_neg(bmap->ineq[i][div_pos], bmap->div[div][0]);
2357
2358         j = isl_basic_map_alloc_inequality(bmap);
2359         if (j < 0)
2360                 return -1;
2361         isl_seq_neg(bmap->ineq[j], bmap->ineq[i], 1 + total);
2362         isl_int_add(bmap->ineq[j][0], bmap->ineq[j][0], bmap->ineq[j][div_pos]);
2363         isl_int_sub_ui(bmap->ineq[j][0], bmap->ineq[j][0], 1);
2364         return j;
2365 }
2366
2367 struct isl_basic_set *isl_basic_map_underlying_set(
2368                 struct isl_basic_map *bmap)
2369 {
2370         if (!bmap)
2371                 goto error;
2372         if (bmap->nparam == 0 && bmap->n_in == 0 && bmap->n_div == 0)
2373                 return (struct isl_basic_set *)bmap;
2374         bmap = isl_basic_map_cow(bmap);
2375         if (!bmap)
2376                 goto error;
2377         bmap->n_out += bmap->nparam + bmap->n_in + bmap->n_div;
2378         bmap->nparam = 0;
2379         bmap->n_in = 0;
2380         bmap->extra -= bmap->n_div;
2381         bmap->n_div = 0;
2382         bmap = isl_basic_map_finalize(bmap);
2383         return (struct isl_basic_set *)bmap;
2384 error:
2385         return NULL;
2386 }
2387
2388 struct isl_basic_map *isl_basic_map_overlying_set(
2389         struct isl_basic_set *bset, struct isl_basic_map *like)
2390 {
2391         struct isl_basic_map *bmap;
2392         struct isl_ctx *ctx;
2393         unsigned total;
2394         int i, k;
2395
2396         if (!bset || !like)
2397                 goto error;
2398         ctx = bset->ctx;
2399         isl_assert(ctx, bset->dim ==
2400                         like->nparam + like->n_in + like->n_out + like->n_div,
2401                         goto error);
2402         if (like->nparam == 0 && like->n_in == 0 && like->n_div == 0) {
2403                 isl_basic_map_free(like);
2404                 return (struct isl_basic_map *)bset;
2405         }
2406         bset = isl_basic_set_cow(bset);
2407         if (!bset)
2408                 goto error;
2409         total = bset->dim + bset->extra;
2410         bmap = (struct isl_basic_map *)bset;
2411         bmap->nparam = like->nparam;
2412         bmap->n_in = like->n_in;
2413         bmap->n_out = like->n_out;
2414         bmap->extra += like->n_div;
2415         if (bmap->extra) {
2416                 unsigned ltotal;
2417                 ltotal = total - bmap->extra + like->extra;
2418                 if (ltotal > total)
2419                         ltotal = total;
2420                 bmap->block2 = isl_blk_extend(ctx, bmap->block2,
2421                                         bmap->extra * (1 + 1 + total));
2422                 if (isl_blk_is_error(bmap->block2))
2423                         goto error;
2424                 bmap->div = isl_realloc_array(ctx, bmap->div, isl_int *,
2425                                                 bmap->extra);
2426                 if (!bmap->div)
2427                         goto error;
2428                 bmap = isl_basic_map_extend(bmap, bmap->nparam,
2429                         bmap->n_in, bmap->n_out, 0, 0, 2 * like->n_div);
2430                 for (i = 0; i < like->n_div; ++i) {
2431                         k = isl_basic_map_alloc_div(bmap);
2432                         if (k < 0)
2433                                 goto error;
2434                         isl_seq_cpy(bmap->div[k], like->div[i], 1 + 1 + ltotal);
2435                         isl_seq_clr(bmap->div[k]+1+1+ltotal, total - ltotal);
2436                         if (add_div_constraints(bmap, k) < 0)
2437                                 goto error;
2438                 }
2439         }
2440         isl_basic_map_free(like);
2441         bmap = isl_basic_map_finalize(bmap);
2442         return bmap;
2443 error:
2444         isl_basic_map_free(like);
2445         isl_basic_set_free(bset);
2446         return NULL;
2447 }
2448
2449 struct isl_basic_set *isl_basic_set_from_underlying_set(
2450         struct isl_basic_set *bset, struct isl_basic_set *like)
2451 {
2452         return (struct isl_basic_set *)
2453                 isl_basic_map_overlying_set(bset, (struct isl_basic_map *)like);
2454 }
2455
2456 struct isl_set *isl_set_from_underlying_set(
2457         struct isl_set *set, struct isl_basic_set *like)
2458 {
2459         int i;
2460
2461         if (!set || !like)
2462                 goto error;
2463         isl_assert(set->ctx, set->dim == like->nparam + like->dim + like->n_div,
2464                     goto error);
2465         if (like->nparam == 0 && like->n_div == 0) {
2466                 isl_basic_set_free(like);
2467                 return set;
2468         }
2469         set = isl_set_cow(set);
2470         if (!set)
2471                 goto error;
2472         for (i = 0; i < set->n; ++i) {
2473                 set->p[i] = isl_basic_set_from_underlying_set(set->p[i],
2474                                                       isl_basic_set_copy(like));
2475                 if (!set->p[i])
2476                         goto error;
2477         }
2478         set->nparam = like->nparam;
2479         set->dim = like->dim;
2480         isl_basic_set_free(like);
2481         return set;
2482 error:
2483         isl_basic_set_free(like);
2484         isl_set_free(set);
2485         return NULL;
2486 }
2487
2488 struct isl_set *isl_map_underlying_set(struct isl_map *map)
2489 {
2490         int i;
2491
2492         map = isl_map_align_divs(map);
2493         map = isl_map_cow(map);
2494         if (!map)
2495                 return NULL;
2496
2497         for (i = 0; i < map->n; ++i) {
2498                 map->p[i] = (struct isl_basic_map *)
2499                                 isl_basic_map_underlying_set(map->p[i]);
2500                 if (!map->p[i])
2501                         goto error;
2502         }
2503         if (map->n == 0)
2504                 map->n_out += map->nparam + map->n_in;
2505         else
2506                 map->n_out = map->p[0]->n_out;
2507         map->nparam = 0;
2508         map->n_in = 0;
2509         return (struct isl_set *)map;
2510 error:
2511         isl_map_free(map);
2512         return NULL;
2513 }
2514
2515 struct isl_set *isl_set_to_underlying_set(struct isl_set *set)
2516 {
2517         return (struct isl_set *)isl_map_underlying_set((struct isl_map *)set);
2518 }
2519
2520 struct isl_basic_set *isl_basic_map_domain(struct isl_basic_map *bmap)
2521 {
2522         struct isl_basic_set *domain;
2523         unsigned n_out;
2524         if (!bmap)
2525                 return NULL;
2526         n_out = bmap->n_out;
2527         domain = isl_basic_set_from_basic_map(bmap);
2528         return isl_basic_set_project_out(domain, n_out, 0);
2529 }
2530
2531 struct isl_basic_set *isl_basic_map_range(struct isl_basic_map *bmap)
2532 {
2533         return isl_basic_map_domain(isl_basic_map_reverse(bmap));
2534 }
2535
2536 struct isl_set *isl_map_range(struct isl_map *map)
2537 {
2538         int i;
2539         struct isl_set *set;
2540
2541         if (!map)
2542                 goto error;
2543         map = isl_map_cow(map);
2544         if (!map)
2545                 goto error;
2546
2547         set = (struct isl_set *) map;
2548         set->zero = 0;
2549         for (i = 0; i < map->n; ++i) {
2550                 set->p[i] = isl_basic_map_range(map->p[i]);
2551                 if (!set->p[i])
2552                         goto error;
2553         }
2554         F_CLR(set, ISL_MAP_DISJOINT);
2555         F_CLR(set, ISL_SET_NORMALIZED);
2556         return set;
2557 error:
2558         isl_map_free(map);
2559         return NULL;
2560 }
2561
2562 struct isl_map *isl_map_from_set(struct isl_set *set,
2563                 unsigned n_in, unsigned n_out)
2564 {
2565         int i;
2566         struct isl_map *map = NULL;
2567
2568         if (!set)
2569                 return NULL;
2570         isl_assert(set->ctx, set->dim == n_in + n_out, goto error);
2571         set = isl_set_cow(set);
2572         if (!set)
2573                 return NULL;
2574         map = (struct isl_map *)set;
2575         for (i = 0; i < set->n; ++i) {
2576                 map->p[i] = isl_basic_map_from_basic_set(
2577                                 set->p[i], n_in, n_out);
2578                 if (!map->p[i])
2579                         goto error;
2580         }
2581         map->n_in = n_in;
2582         map->n_out = n_out;
2583         return map;
2584 error:
2585         isl_set_free(set);
2586         return NULL;
2587 }
2588
2589 struct isl_set *isl_set_from_map(struct isl_map *map)
2590 {
2591         int i;
2592         struct isl_set *set = NULL;
2593
2594         if (!map)
2595                 return NULL;
2596         map = isl_map_cow(map);
2597         if (!map)
2598                 return NULL;
2599         map->n_out += map->n_in;
2600         map->n_in = 0;
2601         set = (struct isl_set *)map;
2602         for (i = 0; i < map->n; ++i) {
2603                 set->p[i] = isl_basic_set_from_basic_map(map->p[i]);
2604                 if (!set->p[i])
2605                         goto error;
2606         }
2607         return set;
2608 error:
2609         isl_map_free(map);
2610         return NULL;
2611 }
2612
2613 struct isl_map *isl_map_alloc(struct isl_ctx *ctx,
2614                 unsigned nparam, unsigned in, unsigned out, int n,
2615                 unsigned flags)
2616 {
2617         struct isl_map *map;
2618
2619         map = isl_alloc(ctx, struct isl_map,
2620                         sizeof(struct isl_map) +
2621                         n * sizeof(struct isl_basic_map *));
2622         if (!map)
2623                 return NULL;
2624
2625         map->ctx = ctx;
2626         isl_ctx_ref(ctx);
2627         map->ref = 1;
2628         map->size = n;
2629         map->n = 0;
2630         map->nparam = nparam;
2631         map->n_in = in;
2632         map->n_out = out;
2633         map->flags = flags;
2634         return map;
2635 }
2636
2637 struct isl_basic_map *isl_basic_map_empty(struct isl_ctx *ctx,
2638                 unsigned nparam, unsigned in, unsigned out)
2639 {
2640         struct isl_basic_map *bmap;
2641         bmap = isl_basic_map_alloc(ctx, nparam, in, out, 0, 1, 0);
2642         bmap = isl_basic_map_set_to_empty(bmap);
2643         return bmap;
2644 }
2645
2646 struct isl_basic_set *isl_basic_set_empty(struct isl_ctx *ctx,
2647                 unsigned nparam, unsigned dim)
2648 {
2649         struct isl_basic_set *bset;
2650         bset = isl_basic_set_alloc(ctx, nparam, dim, 0, 1, 0);
2651         bset = isl_basic_set_set_to_empty(bset);
2652         return bset;
2653 }
2654
2655 struct isl_basic_map *isl_basic_map_universe(struct isl_ctx *ctx,
2656                 unsigned nparam, unsigned in, unsigned out)
2657 {
2658         struct isl_basic_map *bmap;
2659         bmap = isl_basic_map_alloc(ctx, nparam, in, out, 0, 0, 0);
2660         return bmap;
2661 }
2662
2663 struct isl_basic_set *isl_basic_set_universe(struct isl_ctx *ctx,
2664                 unsigned nparam, unsigned dim)
2665 {
2666         struct isl_basic_set *bset;
2667         bset = isl_basic_set_alloc(ctx, nparam, dim, 0, 0, 0);
2668         return bset;
2669 }
2670
2671 struct isl_map *isl_map_empty(struct isl_ctx *ctx,
2672                 unsigned nparam, unsigned in, unsigned out)
2673 {
2674         return isl_map_alloc(ctx, nparam, in, out, 0, ISL_MAP_DISJOINT);
2675 }
2676
2677 struct isl_set *isl_set_empty(struct isl_ctx *ctx,
2678                 unsigned nparam, unsigned dim)
2679 {
2680         return isl_set_alloc(ctx, nparam, dim, 0, ISL_MAP_DISJOINT);
2681 }
2682
2683 struct isl_map *isl_map_dup(struct isl_map *map)
2684 {
2685         int i;
2686         struct isl_map *dup;
2687
2688         if (!map)
2689                 return NULL;
2690         dup = isl_map_alloc(map->ctx, map->nparam, map->n_in, map->n_out, map->n,
2691                                 map->flags);
2692         for (i = 0; i < map->n; ++i)
2693                 dup = isl_map_add(dup, isl_basic_map_copy(map->p[i]));
2694         return dup;
2695 }
2696
2697 struct isl_map *isl_map_add(struct isl_map *map, struct isl_basic_map *bmap)
2698 {
2699         if (!bmap || !map)
2700                 goto error;
2701         isl_assert(map->ctx, map->nparam == bmap->nparam, goto error);
2702         isl_assert(map->ctx, map->n_in == bmap->n_in, goto error);
2703         isl_assert(map->ctx, map->n_out == bmap->n_out, goto error);
2704         isl_assert(map->ctx, map->n < map->size, goto error);
2705         map->p[map->n] = bmap;
2706         map->n++;
2707         F_CLR(map, ISL_MAP_NORMALIZED);
2708         return map;
2709 error:
2710         if (map)
2711                 isl_map_free(map);
2712         if (bmap)
2713                 isl_basic_map_free(bmap);
2714         return NULL;
2715 }
2716
2717 void isl_map_free(struct isl_map *map)
2718 {
2719         int i;
2720
2721         if (!map)
2722                 return;
2723
2724         if (--map->ref > 0)
2725                 return;
2726
2727         isl_ctx_deref(map->ctx);
2728         for (i = 0; i < map->n; ++i)
2729                 isl_basic_map_free(map->p[i]);
2730         free(map);
2731 }
2732
2733 struct isl_map *isl_map_extend(struct isl_map *base,
2734                 unsigned nparam, unsigned n_in, unsigned n_out)
2735 {
2736         int i;
2737
2738         base = isl_map_cow(base);
2739         if (!base)
2740                 return NULL;
2741
2742         isl_assert(base->ctx, base->nparam <= nparam, goto error);
2743         isl_assert(base->ctx, base->n_in <= n_in, goto error);
2744         isl_assert(base->ctx, base->n_out <= n_out, goto error);
2745         base->nparam = nparam;
2746         base->n_in = n_in;
2747         base->n_out = n_out;
2748         for (i = 0; i < base->n; ++i) {
2749                 base->p[i] = isl_basic_map_extend(base->p[i],
2750                                 nparam, n_in, n_out, 0, 0, 0);
2751                 if (!base->p[i])
2752                         goto error;
2753         }
2754         return base;
2755 error:
2756         isl_map_free(base);
2757         return NULL;
2758 }
2759
2760 struct isl_set *isl_set_extend(struct isl_set *base,
2761                 unsigned nparam, unsigned dim)
2762 {
2763         return (struct isl_set *)isl_map_extend((struct isl_map *)base,
2764                                                         nparam, 0, dim);
2765 }
2766
2767 static struct isl_basic_map *isl_basic_map_fix_var(struct isl_basic_map *bmap,
2768                 unsigned var, int value)
2769 {
2770         int j;
2771
2772         bmap = isl_basic_map_cow(bmap);
2773         bmap = isl_basic_map_extend(bmap,
2774                         bmap->nparam, bmap->n_in, bmap->n_out, 0, 1, 0);
2775         j = isl_basic_map_alloc_equality(bmap);
2776         if (j < 0)
2777                 goto error;
2778         isl_seq_clr(bmap->eq[j],
2779                     1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->extra);
2780         isl_int_set_si(bmap->eq[j][1 + var], -1);
2781         isl_int_set_si(bmap->eq[j][0], value);
2782         bmap = isl_basic_map_simplify(bmap);
2783         return isl_basic_map_finalize(bmap);
2784 error:
2785         isl_basic_map_free(bmap);
2786         return NULL;
2787 }
2788
2789 struct isl_basic_map *isl_basic_map_fix_input_si(struct isl_basic_map *bmap,
2790                 unsigned input, int value)
2791 {
2792         if (!bmap)
2793                 return NULL;
2794         isl_assert(bmap->ctx, input < bmap->n_in, goto error);
2795         return isl_basic_map_fix_var(bmap, bmap->nparam + input, value);
2796 error:
2797         isl_basic_map_free(bmap);
2798         return NULL;
2799 }
2800
2801 struct isl_basic_set *isl_basic_set_fix_dim_si(struct isl_basic_set *bset,
2802                 unsigned dim, int value)
2803 {
2804         if (!bset)
2805                 return NULL;
2806         isl_assert(bset->ctx, dim < bset->dim, goto error);
2807         return (struct isl_basic_set *)
2808                 isl_basic_map_fix_var((struct isl_basic_map *)bset,
2809                                                 bset->nparam + dim, value);
2810 error:
2811         isl_basic_set_free(bset);
2812         return NULL;
2813 }
2814
2815 struct isl_map *isl_map_fix_input_si(struct isl_map *map,
2816                 unsigned input, int value)
2817 {
2818         int i;
2819
2820         map = isl_map_cow(map);
2821         if (!map)
2822                 return NULL;
2823
2824         isl_assert(ctx, input < map->n_in, goto error);
2825         for (i = 0; i < map->n; ++i) {
2826                 map->p[i] = isl_basic_map_fix_input_si(map->p[i], input, value);
2827                 if (!map->p[i])
2828                         goto error;
2829         }
2830         F_CLR(map, ISL_MAP_NORMALIZED);
2831         return map;
2832 error:
2833         isl_map_free(map);
2834         return NULL;
2835 }
2836
2837 struct isl_set *isl_set_fix_dim_si(struct isl_set *set, unsigned dim, int value)
2838 {
2839         int i;
2840
2841         set = isl_set_cow(set);
2842         if (!set)
2843                 return NULL;
2844
2845         isl_assert(set->ctx, dim < set->dim, goto error);
2846         for (i = 0; i < set->n; ++i) {
2847                 set->p[i] = isl_basic_set_fix_dim_si(set->p[i], dim, value);
2848                 if (!set->p[i])
2849                         goto error;
2850         }
2851         return set;
2852 error:
2853         isl_set_free(set);
2854         return NULL;
2855 }
2856
2857 struct isl_basic_set *isl_basic_set_lower_bound_dim(struct isl_basic_set *bset,
2858         unsigned dim, isl_int value)
2859 {
2860         int j;
2861
2862         bset = isl_basic_set_cow(bset);
2863         bset = isl_basic_set_extend(bset, bset->nparam, bset->dim, 0, 0, 1);
2864         j = isl_basic_set_alloc_inequality(bset);
2865         if (j < 0)
2866                 goto error;
2867         isl_seq_clr(bset->ineq[j], 1 + bset->nparam + bset->dim + bset->extra);
2868         isl_int_set_si(bset->ineq[j][1 + bset->nparam + dim], 1);
2869         isl_int_neg(bset->ineq[j][0], value);
2870         bset = isl_basic_set_simplify(bset);
2871         return isl_basic_set_finalize(bset);
2872 error:
2873         isl_basic_set_free(bset);
2874         return NULL;
2875 }
2876
2877 struct isl_set *isl_set_lower_bound_dim(struct isl_set *set, unsigned dim,
2878                                         isl_int value)
2879 {
2880         int i;
2881
2882         set = isl_set_cow(set);
2883         if (!set)
2884                 return NULL;
2885
2886         isl_assert(set->ctx, dim < set->dim, goto error);
2887         for (i = 0; i < set->n; ++i) {
2888                 set->p[i] = isl_basic_set_lower_bound_dim(set->p[i], dim, value);
2889                 if (!set->p[i])
2890                         goto error;
2891         }
2892         return set;
2893 error:
2894         isl_set_free(set);
2895         return NULL;
2896 }
2897
2898 struct isl_map *isl_map_reverse(struct isl_map *map)
2899 {
2900         int i;
2901         unsigned t;
2902
2903         map = isl_map_cow(map);
2904         if (!map)
2905                 return NULL;
2906
2907         t = map->n_in;
2908         map->n_in = map->n_out;
2909         map->n_out = t;
2910         for (i = 0; i < map->n; ++i) {
2911                 map->p[i] = isl_basic_map_reverse(map->p[i]);
2912                 if (!map->p[i])
2913                         goto error;
2914         }
2915         F_CLR(map, ISL_MAP_NORMALIZED);
2916         return map;
2917 error:
2918         isl_map_free(map);
2919         return NULL;
2920 }
2921
2922 struct isl_map *isl_basic_map_lexmax(
2923                 struct isl_basic_map *bmap, struct isl_basic_set *dom,
2924                 struct isl_set **empty)
2925 {
2926         return isl_pip_basic_map_lexmax(bmap, dom, empty);
2927 }
2928
2929 struct isl_map *isl_basic_map_lexmin(
2930                 struct isl_basic_map *bmap, struct isl_basic_set *dom,
2931                 struct isl_set **empty)
2932 {
2933         return isl_pip_basic_map_lexmin(bmap, dom, empty);
2934 }
2935
2936 struct isl_set *isl_basic_set_lexmin(struct isl_basic_set *bset)
2937 {
2938         struct isl_basic_map *bmap = NULL;
2939         struct isl_basic_set *dom = NULL;
2940         struct isl_map *min;
2941
2942         if (!bset)
2943                 goto error;
2944         bmap = isl_basic_map_from_basic_set(bset, 0, bset->dim);
2945         if (!bmap)
2946                 goto error;
2947         dom = isl_basic_set_alloc(bmap->ctx, bmap->nparam, 0, 0, 0, 0);
2948         if (!dom)
2949                 goto error;
2950         min = isl_basic_map_lexmin(bmap, dom, NULL);
2951         return isl_map_range(min);
2952 error:
2953         isl_basic_map_free(bmap);
2954         return NULL;
2955 }
2956
2957 struct isl_map *isl_basic_map_compute_divs(struct isl_basic_map *bmap)
2958 {
2959         if (!bmap)
2960                 return NULL;
2961         if (bmap->n_div == 0)
2962                 return isl_map_from_basic_map(bmap);
2963         return isl_pip_basic_map_compute_divs(bmap);
2964 }
2965
2966 struct isl_map *isl_map_compute_divs(struct isl_map *map)
2967 {
2968         int i;
2969         struct isl_map *res;
2970
2971         if (!map)
2972                 return NULL;
2973         if (map->n == 0)
2974                 return map;
2975         res = isl_basic_map_compute_divs(isl_basic_map_copy(map->p[0]));
2976         for (i = 1 ; i < map->n; ++i) {
2977                 struct isl_map *r2;
2978                 r2 = isl_basic_map_compute_divs(isl_basic_map_copy(map->p[i]));
2979                 if (F_ISSET(map, ISL_MAP_DISJOINT))
2980                         res = isl_map_union_disjoint(res, r2);
2981                 else
2982                         res = isl_map_union(res, r2);
2983         }
2984         isl_map_free(map);
2985
2986         return res;
2987 }
2988
2989 struct isl_set *isl_basic_set_compute_divs(struct isl_basic_set *bset)
2990 {
2991         return (struct isl_set *)
2992                 isl_basic_map_compute_divs((struct isl_basic_map *)bset);
2993 }
2994
2995 struct isl_set *isl_set_compute_divs(struct isl_set *set)
2996 {
2997         return (struct isl_set *)
2998                 isl_map_compute_divs((struct isl_map *)set);
2999 }
3000
3001 struct isl_set *isl_map_domain(struct isl_map *map)
3002 {
3003         int i;
3004         struct isl_set *set;
3005
3006         if (!map)
3007                 goto error;
3008
3009         map = isl_map_cow(map);
3010         if (!map)
3011                 return NULL;
3012
3013         set = (struct isl_set *)map;
3014         set->dim = map->n_in;
3015         set->zero = 0;
3016         for (i = 0; i < map->n; ++i) {
3017                 set->p[i] = isl_basic_map_domain(map->p[i]);
3018                 if (!set->p[i])
3019                         goto error;
3020         }
3021         F_CLR(set, ISL_MAP_DISJOINT);
3022         F_CLR(set, ISL_SET_NORMALIZED);
3023         return set;
3024 error:
3025         isl_map_free(map);
3026         return NULL;
3027 }
3028
3029 struct isl_map *isl_map_union_disjoint(
3030                         struct isl_map *map1, struct isl_map *map2)
3031 {
3032         int i;
3033         unsigned flags = 0;
3034         struct isl_map *map = NULL;
3035
3036         if (!map1 || !map2)
3037                 goto error;
3038
3039         if (map1->n == 0) {
3040                 isl_map_free(map1);
3041                 return map2;
3042         }
3043         if (map2->n == 0) {
3044                 isl_map_free(map2);
3045                 return map1;
3046         }
3047
3048         isl_assert(map1->ctx, map1->nparam == map2->nparam, goto error);
3049         isl_assert(map1->ctx, map1->n_in == map2->n_in, goto error);
3050         isl_assert(map1->ctx, map1->n_out == map2->n_out, goto error);
3051
3052         if (F_ISSET(map1, ISL_MAP_DISJOINT) &&
3053             F_ISSET(map2, ISL_MAP_DISJOINT))
3054                 FL_SET(flags, ISL_MAP_DISJOINT);
3055
3056         map = isl_map_alloc(map1->ctx, map1->nparam, map1->n_in, map1->n_out,
3057                                 map1->n + map2->n, flags);
3058         if (!map)
3059                 goto error;
3060         for (i = 0; i < map1->n; ++i) {
3061                 map = isl_map_add(map,
3062                                   isl_basic_map_copy(map1->p[i]));
3063                 if (!map)
3064                         goto error;
3065         }
3066         for (i = 0; i < map2->n; ++i) {
3067                 map = isl_map_add(map,
3068                                   isl_basic_map_copy(map2->p[i]));
3069                 if (!map)
3070                         goto error;
3071         }
3072         isl_map_free(map1);
3073         isl_map_free(map2);
3074         return map;
3075 error:
3076         isl_map_free(map);
3077         isl_map_free(map1);
3078         isl_map_free(map2);
3079         return NULL;
3080 }
3081
3082 struct isl_map *isl_map_union(struct isl_map *map1, struct isl_map *map2)
3083 {
3084         map1 = isl_map_union_disjoint(map1, map2);
3085         if (!map1)
3086                 return NULL;
3087         if (map1->n > 1)
3088                 F_CLR(map1, ISL_MAP_DISJOINT);
3089         return map1;
3090 }
3091
3092 struct isl_set *isl_set_union_disjoint(
3093                         struct isl_set *set1, struct isl_set *set2)
3094 {
3095         return (struct isl_set *)
3096                 isl_map_union_disjoint(
3097                         (struct isl_map *)set1, (struct isl_map *)set2);
3098 }
3099
3100 struct isl_set *isl_set_union(struct isl_set *set1, struct isl_set *set2)
3101 {
3102         return (struct isl_set *)
3103                 isl_map_union((struct isl_map *)set1, (struct isl_map *)set2);
3104 }
3105
3106 struct isl_map *isl_map_intersect_range(
3107                 struct isl_map *map, struct isl_set *set)
3108 {
3109         unsigned flags = 0;
3110         struct isl_map *result;
3111         int i, j;
3112
3113         if (!map || !set)
3114                 goto error;
3115
3116         if (F_ISSET(map, ISL_MAP_DISJOINT) &&
3117             F_ISSET(set, ISL_MAP_DISJOINT))
3118                 FL_SET(flags, ISL_MAP_DISJOINT);
3119
3120         result = isl_map_alloc(map->ctx, map->nparam, map->n_in, map->n_out,
3121                                 map->n * set->n, flags);
3122         if (!result)
3123                 goto error;
3124         for (i = 0; i < map->n; ++i)
3125                 for (j = 0; j < set->n; ++j) {
3126                         result = isl_map_add(result,
3127                             isl_basic_map_intersect_range(
3128                                 isl_basic_map_copy(map->p[i]),
3129                                 isl_basic_set_copy(set->p[j])));
3130                         if (!result)
3131                                 goto error;
3132                 }
3133         isl_map_free(map);
3134         isl_set_free(set);
3135         return result;
3136 error:
3137         isl_map_free(map);
3138         isl_set_free(set);
3139         return NULL;
3140 }
3141
3142 struct isl_map *isl_map_intersect_domain(
3143                 struct isl_map *map, struct isl_set *set)
3144 {
3145         return isl_map_reverse(
3146                 isl_map_intersect_range(isl_map_reverse(map), set));
3147 }
3148
3149 struct isl_map *isl_map_apply_domain(
3150                 struct isl_map *map1, struct isl_map *map2)
3151 {
3152         if (!map1 || !map2)
3153                 goto error;
3154         map1 = isl_map_reverse(map1);
3155         map1 = isl_map_apply_range(map1, map2);
3156         return isl_map_reverse(map1);
3157 error:
3158         isl_map_free(map1);
3159         isl_map_free(map2);
3160         return NULL;
3161 }
3162
3163 struct isl_map *isl_map_apply_range(
3164                 struct isl_map *map1, struct isl_map *map2)
3165 {
3166         struct isl_map *result;
3167         int i, j;
3168
3169         if (!map1 || !map2)
3170                 goto error;
3171
3172         isl_assert(map1->ctx, map1->nparam == map2->nparam, goto error);
3173         isl_assert(map1->ctx, map1->n_out == map2->n_in, goto error);
3174
3175         result = isl_map_alloc(map1->ctx, map1->nparam, map1->n_in, map2->n_out,
3176                                 map1->n * map2->n, 0);
3177         if (!result)
3178                 goto error;
3179         for (i = 0; i < map1->n; ++i)
3180                 for (j = 0; j < map2->n; ++j) {
3181                         result = isl_map_add(result,
3182                             isl_basic_map_apply_range(
3183                                 isl_basic_map_copy(map1->p[i]),
3184                                 isl_basic_map_copy(map2->p[j])));
3185                         if (!result)
3186                                 goto error;
3187                 }
3188         isl_map_free(map1);
3189         isl_map_free(map2);
3190         if (result && result->n <= 1)
3191                 F_SET(result, ISL_MAP_DISJOINT);
3192         return result;
3193 error:
3194         isl_map_free(map1);
3195         isl_map_free(map2);
3196         return NULL;
3197 }
3198
3199 /*
3200  * returns range - domain
3201  */
3202 struct isl_basic_set *isl_basic_map_deltas(struct isl_basic_map *bmap)
3203 {
3204         struct isl_basic_set *bset;
3205         unsigned dim;
3206         int i;
3207
3208         if (!bmap)
3209                 goto error;
3210         isl_assert(bmap->ctx, bmap->n_in == bmap->n_out, goto error);
3211         dim = bmap->n_in;
3212         bset = isl_basic_set_from_basic_map(bmap);
3213         bset = isl_basic_set_extend(bset, bset->nparam, 3*dim, 0,
3214                                         dim, 0);
3215         bset = isl_basic_set_swap_vars(bset, 2*dim);
3216         for (i = 0; i < dim; ++i) {
3217                 int j = isl_basic_map_alloc_equality(
3218                                             (struct isl_basic_map *)bset);
3219                 if (j < 0)
3220                         goto error;
3221                 isl_seq_clr(bset->eq[j],
3222                             1 + bset->nparam + bset->dim + bset->extra);
3223                 isl_int_set_si(bset->eq[j][1+bset->nparam+i], 1);
3224                 isl_int_set_si(bset->eq[j][1+bset->nparam+dim+i], 1);
3225                 isl_int_set_si(bset->eq[j][1+bset->nparam+2*dim+i], -1);
3226         }
3227         return isl_basic_set_project_out(bset, 2*dim, 0);
3228 error:
3229         isl_basic_map_free(bmap);
3230         return NULL;
3231 }
3232
3233 /*
3234  * returns range - domain
3235  */
3236 struct isl_set *isl_map_deltas(struct isl_map *map)
3237 {
3238         int i;
3239         struct isl_set *result;
3240
3241         if (!map)
3242                 return NULL;
3243
3244         isl_assert(map->ctx, map->n_in == map->n_out, goto error);
3245         result = isl_set_alloc(map->ctx, map->nparam, map->n_in, map->n, map->flags);
3246         if (!result)
3247                 goto error;
3248         for (i = 0; i < map->n; ++i)
3249                 result = isl_set_add(result,
3250                           isl_basic_map_deltas(isl_basic_map_copy(map->p[i])));
3251         isl_map_free(map);
3252         return result;
3253 error:
3254         isl_map_free(map);
3255         return NULL;
3256 }
3257
3258 /* If the only constraints a div d=floor(f/m)
3259  * appears in are its two defining constraints
3260  *
3261  *      f - m d >=0
3262  *      -(f - (m - 1)) + m d >= 0
3263  *
3264  * then it can safely be removed.
3265  */
3266 static int div_is_redundant(struct isl_basic_map *bmap, int div)
3267 {
3268         int i;
3269         unsigned pos = 1 + bmap->nparam + bmap->n_in + bmap->n_out + div;
3270
3271         for (i = 0; i < bmap->n_eq; ++i)
3272                 if (!isl_int_is_zero(bmap->eq[i][pos]))
3273                         return 0;
3274
3275         for (i = 0; i < bmap->n_ineq; ++i) {
3276                 if (isl_int_is_zero(bmap->ineq[i][pos]))
3277                         continue;
3278                 if (isl_int_eq(bmap->ineq[i][pos], bmap->div[div][0])) {
3279                         int neg;
3280                         isl_int_sub(bmap->div[div][1],
3281                                         bmap->div[div][1], bmap->div[div][0]);
3282                         isl_int_add_ui(bmap->div[div][1], bmap->div[div][1], 1);
3283                         neg = isl_seq_is_neg(bmap->ineq[i], bmap->div[div]+1, pos);
3284                         isl_int_sub_ui(bmap->div[div][1], bmap->div[div][1], 1);
3285                         isl_int_add(bmap->div[div][1],
3286                                         bmap->div[div][1], bmap->div[div][0]);
3287                         if (!neg)
3288                                 return 0;
3289                         if (isl_seq_first_non_zero(bmap->ineq[i]+pos+1,
3290                                                     bmap->n_div-div-1) != -1)
3291                                 return 0;
3292                 } else if (isl_int_abs_eq(bmap->ineq[i][pos], bmap->div[div][0])) {
3293                         if (!isl_seq_eq(bmap->ineq[i], bmap->div[div]+1, pos))
3294                                 return 0;
3295                         if (isl_seq_first_non_zero(bmap->ineq[i]+pos+1,
3296                                                     bmap->n_div-div-1) != -1)
3297                                 return 0;
3298                 } else
3299                         return 0;
3300         }
3301
3302         for (i = 0; i < bmap->n_div; ++i)
3303                 if (!isl_int_is_zero(bmap->div[i][1+pos]))
3304                         return 0;
3305
3306         return 1;
3307 }
3308
3309 /*
3310  * Remove divs that don't occur in any of the constraints or other divs.
3311  * These can arise when dropping some of the variables in a quast
3312  * returned by piplib.
3313  */
3314 static struct isl_basic_map *remove_redundant_divs(struct isl_basic_map *bmap)
3315 {
3316         int i;
3317
3318         if (!bmap)
3319                 return NULL;
3320
3321         for (i = bmap->n_div-1; i >= 0; --i) {
3322                 if (!div_is_redundant(bmap, i))
3323                         continue;
3324                 bmap = isl_basic_map_drop_div(bmap, i);
3325         }
3326         return bmap;
3327 }
3328
3329 struct isl_basic_map *isl_basic_map_finalize(struct isl_basic_map *bmap)
3330 {
3331         bmap = remove_redundant_divs(bmap);
3332         if (!bmap)
3333                 return NULL;
3334         F_SET(bmap, ISL_BASIC_SET_FINAL);
3335         return bmap;
3336 }
3337
3338 struct isl_basic_set *isl_basic_set_finalize(struct isl_basic_set *bset)
3339 {
3340         return (struct isl_basic_set *)
3341                 isl_basic_map_finalize((struct isl_basic_map *)bset);
3342 }
3343
3344 struct isl_set *isl_set_finalize(struct isl_set *set)
3345 {
3346         int i;
3347
3348         if (!set)
3349                 return NULL;
3350         for (i = 0; i < set->n; ++i) {
3351                 set->p[i] = isl_basic_set_finalize(set->p[i]);
3352                 if (!set->p[i])
3353                         goto error;
3354         }
3355         return set;
3356 error:
3357         isl_set_free(set);
3358         return NULL;
3359 }
3360
3361 struct isl_map *isl_map_finalize(struct isl_map *map)
3362 {
3363         int i;
3364
3365         if (!map)
3366                 return NULL;
3367         for (i = 0; i < map->n; ++i) {
3368                 map->p[i] = isl_basic_map_finalize(map->p[i]);
3369                 if (!map->p[i])
3370                         goto error;
3371         }
3372         F_CLR(map, ISL_MAP_NORMALIZED);
3373         return map;
3374 error:
3375         isl_map_free(map);
3376         return NULL;
3377 }
3378
3379 struct isl_basic_map *isl_basic_map_identity(struct isl_ctx *ctx,
3380                 unsigned nparam, unsigned dim)
3381 {
3382         struct isl_basic_map *bmap;
3383         int i;
3384
3385         bmap = isl_basic_map_alloc(ctx, nparam, dim, dim, 0, dim, 0);
3386         if (!bmap)
3387                 goto error;
3388
3389         for (i = 0; i < dim; ++i) {
3390                 int j = isl_basic_map_alloc_equality(bmap);
3391                 if (j < 0)
3392                         goto error;
3393                 isl_seq_clr(bmap->eq[j],
3394                     1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->extra);
3395                 isl_int_set_si(bmap->eq[j][1+nparam+i], 1);
3396                 isl_int_set_si(bmap->eq[j][1+nparam+dim+i], -1);
3397         }
3398         return isl_basic_map_finalize(bmap);
3399 error:
3400         isl_basic_map_free(bmap);
3401         return NULL;
3402 }
3403
3404 struct isl_map *isl_map_identity(struct isl_ctx *ctx,
3405                 unsigned nparam, unsigned dim)
3406 {
3407         struct isl_map *map = isl_map_alloc(ctx, nparam, dim, dim, 1,
3408                                                 ISL_MAP_DISJOINT);
3409         if (!map)
3410                 goto error;
3411         map = isl_map_add(map,
3412                 isl_basic_map_identity(ctx, nparam, dim));
3413         return map;
3414 error:
3415         isl_map_free(map);
3416         return NULL;
3417 }
3418
3419 int isl_set_is_equal(struct isl_set *set1, struct isl_set *set2)
3420 {
3421         return isl_map_is_equal((struct isl_map *)set1, (struct isl_map *)set2);
3422 }
3423
3424 int isl_set_is_subset(struct isl_set *set1, struct isl_set *set2)
3425 {
3426         return isl_map_is_subset(
3427                         (struct isl_map *)set1, (struct isl_map *)set2);
3428 }
3429
3430 int isl_basic_map_is_subset(
3431                 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
3432 {
3433         int is_subset;
3434         struct isl_map *map1;
3435         struct isl_map *map2;
3436
3437         if (!bmap1 || !bmap2)
3438                 return -1;
3439
3440         map1 = isl_map_from_basic_map(isl_basic_map_copy(bmap1));
3441         map2 = isl_map_from_basic_map(isl_basic_map_copy(bmap2));
3442
3443         is_subset = isl_map_is_subset(map1, map2);
3444
3445         isl_map_free(map1);
3446         isl_map_free(map2);
3447
3448         return is_subset;
3449 }
3450
3451 int isl_basic_map_is_equal(
3452                 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
3453 {
3454         int is_subset;
3455
3456         if (!bmap1 || !bmap2)
3457                 return -1;
3458         is_subset = isl_basic_map_is_subset(bmap1, bmap2);
3459         if (is_subset != 1)
3460                 return is_subset;
3461         is_subset = isl_basic_map_is_subset(bmap2, bmap1);
3462         return is_subset;
3463 }
3464
3465 int isl_basic_set_is_equal(
3466                 struct isl_basic_set *bset1, struct isl_basic_set *bset2)
3467 {
3468         return isl_basic_map_is_equal(
3469                 (struct isl_basic_map *)bset1, (struct isl_basic_map *)bset2);
3470 }
3471
3472 int isl_map_is_empty(struct isl_map *map)
3473 {
3474         int i;
3475         int is_empty;
3476
3477         if (!map)
3478                 return -1;
3479         for (i = 0; i < map->n; ++i) {
3480                 is_empty = isl_basic_map_is_empty(map->p[i]);
3481                 if (is_empty < 0)
3482                         return -1;
3483                 if (!is_empty)
3484                         return 0;
3485         }
3486         return 1;
3487 }
3488
3489 int isl_set_is_empty(struct isl_set *set)
3490 {
3491         return isl_map_is_empty((struct isl_map *)set);
3492 }
3493
3494 int isl_map_is_subset(struct isl_map *map1, struct isl_map *map2)
3495 {
3496         int i;
3497         int is_subset = 0;
3498         struct isl_map *diff;
3499
3500         if (!map1 || !map2)
3501                 return -1;
3502
3503         if (isl_map_is_empty(map1))
3504                 return 1;
3505
3506         if (isl_map_is_empty(map2))
3507                 return 0;
3508
3509         diff = isl_map_subtract(isl_map_copy(map1), isl_map_copy(map2));
3510         if (!diff)
3511                 return -1;
3512
3513         is_subset = isl_map_is_empty(diff);
3514         isl_map_free(diff);
3515
3516         return is_subset;
3517 }
3518
3519 int isl_map_is_equal(struct isl_map *map1, struct isl_map *map2)
3520 {
3521         int is_subset;
3522
3523         if (!map1 || !map2)
3524                 return -1;
3525         is_subset = isl_map_is_subset(map1, map2);
3526         if (is_subset != 1)
3527                 return is_subset;
3528         is_subset = isl_map_is_subset(map2, map1);
3529         return is_subset;
3530 }
3531
3532 int isl_basic_map_is_strict_subset(
3533                 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
3534 {
3535         int is_subset;
3536
3537         if (!bmap1 || !bmap2)
3538                 return -1;
3539         is_subset = isl_basic_map_is_subset(bmap1, bmap2);
3540         if (is_subset != 1)
3541                 return is_subset;
3542         is_subset = isl_basic_map_is_subset(bmap2, bmap1);
3543         if (is_subset == -1)
3544                 return is_subset;
3545         return !is_subset;
3546 }
3547
3548 static int basic_map_contains(struct isl_basic_map *bmap, struct isl_vec *vec)
3549 {
3550         int i;
3551         unsigned total;
3552         isl_int s;
3553
3554         total = 1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
3555         if (total != vec->size)
3556                 return -1;
3557
3558         isl_int_init(s);
3559
3560         for (i = 0; i < bmap->n_eq; ++i) {
3561                 isl_seq_inner_product(vec->block.data, bmap->eq[i], total, &s);
3562                 if (!isl_int_is_zero(s)) {
3563                         isl_int_clear(s);
3564                         return 0;
3565                 }
3566         }
3567
3568         for (i = 0; i < bmap->n_ineq; ++i) {
3569                 isl_seq_inner_product(vec->block.data, bmap->ineq[i], total, &s);
3570                 if (isl_int_is_neg(s)) {
3571                         isl_int_clear(s);
3572                         return 0;
3573                 }
3574         }
3575
3576         isl_int_clear(s);
3577
3578         return 1;
3579 }
3580
3581 int isl_basic_map_is_universe(struct isl_basic_map *bmap)
3582 {
3583         if (!bmap)
3584                 return -1;
3585         return bmap->n_eq == 0 && bmap->n_ineq == 0;
3586 }
3587
3588 int isl_basic_map_is_empty(struct isl_basic_map *bmap)
3589 {
3590         struct isl_basic_set *bset = NULL;
3591         struct isl_vec *sample = NULL;
3592         int empty;
3593         unsigned total;
3594
3595         if (!bmap)
3596                 return -1;
3597
3598         if (F_ISSET(bmap, ISL_BASIC_MAP_EMPTY))
3599                 return 1;
3600
3601         if (F_ISSET(bmap, ISL_BASIC_MAP_RATIONAL)) {
3602                 struct isl_basic_map *copy = isl_basic_map_copy(bmap);
3603                 copy = isl_basic_map_convex_hull(copy);
3604                 empty = F_ISSET(copy, ISL_BASIC_MAP_EMPTY);
3605                 isl_basic_map_free(copy);
3606                 return empty;
3607         }
3608
3609         total = 1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
3610         if (bmap->sample && bmap->sample->size == total) {
3611                 int contains = basic_map_contains(bmap, bmap->sample);
3612                 if (contains < 0)
3613                         return -1;
3614                 if (contains)
3615                         return 0;
3616         }
3617         bset = isl_basic_map_underlying_set(isl_basic_map_copy(bmap));
3618         if (!bset)
3619                 return -1;
3620         sample = isl_basic_set_sample(bset);
3621         if (!sample)
3622                 return -1;
3623         empty = sample->size == 0;
3624         if (bmap->sample)
3625                 isl_vec_free(bmap->ctx, bmap->sample);
3626         bmap->sample = sample;
3627
3628         return empty;
3629 }
3630
3631 struct isl_map *isl_basic_map_union(
3632         struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
3633 {
3634         struct isl_map *map;
3635         if (!bmap1 || !bmap2)
3636                 return NULL;
3637
3638         isl_assert(bmap1->ctx, bmap1->nparam == bmap2->nparam, goto error);
3639         isl_assert(bmap1->ctx, bmap1->n_in == bmap2->n_in, goto error);
3640         isl_assert(bmap1->ctx, bmap1->n_out == bmap2->n_out, goto error);
3641
3642         map = isl_map_alloc(bmap1->ctx, bmap1->nparam,
3643                                 bmap1->n_in, bmap1->n_out, 2, 0);
3644         if (!map)
3645                 goto error;
3646         map = isl_map_add(map, bmap1);
3647         map = isl_map_add(map, bmap2);
3648         return map;
3649 error:
3650         isl_basic_map_free(bmap1);
3651         isl_basic_map_free(bmap2);
3652         return NULL;
3653 }
3654
3655 struct isl_set *isl_basic_set_union(
3656                 struct isl_basic_set *bset1, struct isl_basic_set *bset2)
3657 {
3658         return (struct isl_set *)isl_basic_map_union(
3659                                             (struct isl_basic_map *)bset1,
3660                                             (struct isl_basic_map *)bset2);
3661 }
3662
3663 /* Order divs such that any div only depends on previous divs */
3664 static struct isl_basic_map *order_divs(struct isl_basic_map *bmap)
3665 {
3666         int i;
3667         unsigned off = bmap->nparam + bmap->n_in + bmap->n_out;
3668
3669         for (i = 0; i < bmap->n_div; ++i) {
3670                 int pos;
3671                 pos = isl_seq_first_non_zero(bmap->div[i]+1+1+off+i,
3672                                                             bmap->n_div-i);
3673                 if (pos == -1)
3674                         continue;
3675                 swap_div(bmap, i, pos);
3676                 --i;
3677         }
3678         return bmap;
3679 }
3680
3681 static int find_div(struct isl_basic_map *dst,
3682                         struct isl_basic_map *src, unsigned div)
3683 {
3684         int i;
3685
3686         unsigned total = src->nparam + src->n_in + src->n_out + src->n_div;
3687
3688         for (i = 0; i < dst->n_div; ++i)
3689                 if (isl_seq_eq(dst->div[i], src->div[div], 1+1+total) &&
3690                     isl_seq_first_non_zero(dst->div[i]+1+1+total,
3691                                                 dst->n_div - src->n_div) == -1)
3692                         return i;
3693         return -1;
3694 }
3695
3696 struct isl_basic_map *isl_basic_map_align_divs(
3697                 struct isl_basic_map *dst, struct isl_basic_map *src)
3698 {
3699         int i;
3700         unsigned total = src->nparam + src->n_in + src->n_out + src->n_div;
3701
3702         if (!dst || !src)
3703                 goto error;
3704
3705         if (src->n_div == 0)
3706                 return dst;
3707
3708         src = order_divs(src);
3709         dst = isl_basic_map_extend(dst, dst->nparam, dst->n_in,
3710                         dst->n_out, src->n_div, 0, 2 * src->n_div);
3711         if (!dst)
3712                 return NULL;
3713         for (i = 0; i < src->n_div; ++i) {
3714                 int j = find_div(dst, src, i);
3715                 if (j < 0) {
3716                         j = isl_basic_map_alloc_div(dst);
3717                         if (j < 0)
3718                                 goto error;
3719                         isl_seq_cpy(dst->div[j], src->div[i], 1+1+total);
3720                         isl_seq_clr(dst->div[j]+1+1+total,
3721                                                     dst->extra - src->n_div);
3722                         if (add_div_constraints(dst, j) < 0)
3723                                 goto error;
3724                 }
3725                 if (j != i)
3726                         swap_div(dst, i, j);
3727         }
3728         return dst;
3729 error:
3730         isl_basic_map_free(dst);
3731         return NULL;
3732 }
3733
3734 struct isl_basic_set *isl_basic_set_align_divs(
3735                 struct isl_basic_set *dst, struct isl_basic_set *src)
3736 {
3737         return (struct isl_basic_set *)isl_basic_map_align_divs(
3738                 (struct isl_basic_map *)dst, (struct isl_basic_map *)src);
3739 }
3740
3741 struct isl_map *isl_map_align_divs(struct isl_map *map)
3742 {
3743         int i;
3744
3745         map = isl_map_compute_divs(map);
3746         map = isl_map_cow(map);
3747         if (!map)
3748                 return NULL;
3749
3750         for (i = 1; i < map->n; ++i)
3751                 map->p[0] = isl_basic_map_align_divs(map->p[0], map->p[i]);
3752         for (i = 1; i < map->n; ++i)
3753                 map->p[i] = isl_basic_map_align_divs(map->p[i], map->p[0]);
3754
3755         F_CLR(map, ISL_MAP_NORMALIZED);
3756         return map;
3757 }
3758
3759 static struct isl_map *add_cut_constraint(struct isl_map *dst,
3760                 struct isl_basic_map *src, isl_int *c,
3761                 unsigned len, int oppose)
3762 {
3763         struct isl_basic_map *copy = NULL;
3764         int is_empty;
3765         int k;
3766         unsigned total;
3767
3768         copy = isl_basic_map_copy(src);
3769         copy = isl_basic_map_cow(copy);
3770         if (!copy)
3771                 goto error;
3772         copy = isl_basic_map_extend(copy,
3773                 copy->nparam, copy->n_in, copy->n_out, 0, 0, 1);
3774         k = isl_basic_map_alloc_inequality(copy);
3775         if (k < 0)
3776                 goto error;
3777         if (oppose)
3778                 isl_seq_neg(copy->ineq[k], c, len);
3779         else
3780                 isl_seq_cpy(copy->ineq[k], c, len);
3781         total = 1 + copy->nparam + copy->n_in + copy->n_out + copy->extra;
3782         isl_seq_clr(copy->ineq[k]+len, total - len);
3783         isl_inequality_negate(copy, k);
3784         copy = isl_basic_map_simplify(copy);
3785         copy = isl_basic_map_finalize(copy);
3786         is_empty = isl_basic_map_is_empty(copy);
3787         if (is_empty < 0)
3788                 goto error;
3789         if (!is_empty)
3790                 dst = isl_map_add(dst, copy);
3791         else
3792                 isl_basic_map_free(copy);
3793         return dst;
3794 error:
3795         isl_basic_map_free(copy);
3796         isl_map_free(dst);
3797         return NULL;
3798 }
3799
3800 static struct isl_map *subtract(struct isl_map *map, struct isl_basic_map *bmap)
3801 {
3802         int i, j, k;
3803         unsigned flags = 0;
3804         struct isl_map *rest = NULL;
3805         unsigned max;
3806         unsigned total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
3807
3808         assert(bmap);
3809
3810         if (!map)
3811                 goto error;
3812
3813         if (F_ISSET(map, ISL_MAP_DISJOINT))
3814                 FL_SET(flags, ISL_MAP_DISJOINT);
3815
3816         max = map->n * (2 * bmap->n_eq + bmap->n_ineq);
3817         rest = isl_map_alloc(map->ctx, map->nparam, map->n_in, map->n_out,
3818                                 max, flags);
3819         if (!rest)
3820                 goto error;
3821
3822         for (i = 0; i < map->n; ++i) {
3823                 map->p[i] = isl_basic_map_align_divs(map->p[i], bmap);
3824                 if (!map->p[i])
3825                         goto error;
3826         }
3827
3828         for (j = 0; j < map->n; ++j)
3829                 map->p[j] = isl_basic_map_cow(map->p[j]);
3830
3831         for (i = 0; i < bmap->n_eq; ++i) {
3832                 for (j = 0; j < map->n; ++j) {
3833                         rest = add_cut_constraint(rest,
3834                                 map->p[j], bmap->eq[i], 1+total, 0);
3835                         if (!rest)
3836                                 goto error;
3837
3838                         rest = add_cut_constraint(rest,
3839                                 map->p[j], bmap->eq[i], 1+total, 1);
3840                         if (!rest)
3841                                 goto error;
3842
3843                         map->p[j] = isl_basic_map_extend(map->p[j],
3844                                 map->p[j]->nparam, map->p[j]->n_in,
3845                                 map->p[j]->n_out, 0, 1, 0);
3846                         if (!map->p[j])
3847                                 goto error;
3848                         k = isl_basic_map_alloc_equality(map->p[j]);
3849                         if (k < 0)
3850                                 goto error;
3851                         isl_seq_cpy(map->p[j]->eq[k], bmap->eq[i], 1+total);
3852                         isl_seq_clr(map->p[j]->eq[k]+1+total,
3853                                         map->p[j]->extra - bmap->n_div);
3854                 }
3855         }
3856
3857         for (i = 0; i < bmap->n_ineq; ++i) {
3858                 for (j = 0; j < map->n; ++j) {
3859                         rest = add_cut_constraint(rest,
3860                                 map->p[j], bmap->ineq[i], 1+total, 0);
3861                         if (!rest)
3862                                 goto error;
3863
3864                         map->p[j] = isl_basic_map_extend(map->p[j],
3865                                 map->p[j]->nparam, map->p[j]->n_in,
3866                                 map->p[j]->n_out, 0, 0, 1);
3867                         if (!map->p[j])
3868                                 goto error;
3869                         k = isl_basic_map_alloc_inequality(map->p[j]);
3870                         if (k < 0)
3871                                 goto error;
3872                         isl_seq_cpy(map->p[j]->ineq[k], bmap->ineq[i], 1+total);
3873                         isl_seq_clr(map->p[j]->ineq[k]+1+total,
3874                                         map->p[j]->extra - bmap->n_div);
3875                 }
3876         }
3877
3878         isl_map_free(map);
3879         return rest;
3880 error:
3881         isl_map_free(map);
3882         isl_map_free(rest);
3883         return NULL;
3884 }
3885
3886 struct isl_map *isl_map_subtract(struct isl_map *map1, struct isl_map *map2)
3887 {
3888         int i;
3889         if (!map1 || !map2)
3890                 goto error;
3891
3892         isl_assert(map1->ctx, map1->nparam == map2->nparam, goto error);
3893         isl_assert(map1->ctx, map1->n_in == map2->n_in, goto error);
3894         isl_assert(map1->ctx, map1->n_out == map2->n_out, goto error);
3895
3896         if (isl_map_is_empty(map2)) {
3897                 isl_map_free(map2);
3898                 return map1;
3899         }
3900
3901         map1 = isl_map_compute_divs(map1);
3902         map2 = isl_map_compute_divs(map2);
3903         if (!map1 || !map2)
3904                 goto error;
3905
3906         for (i = 0; map1 && i < map2->n; ++i)
3907                 map1 = subtract(map1, map2->p[i]);
3908
3909         isl_map_free(map2);
3910         return map1;
3911 error:
3912         isl_map_free(map1);
3913         isl_map_free(map2);
3914         return NULL;
3915 }
3916
3917 struct isl_set *isl_set_subtract(struct isl_set *set1, struct isl_set *set2)
3918 {
3919         return (struct isl_set *)
3920                 isl_map_subtract(
3921                         (struct isl_map *)set1, (struct isl_map *)set2);
3922 }
3923
3924 struct isl_set *isl_set_apply(struct isl_set *set, struct isl_map *map)
3925 {
3926         if (!set || !map)
3927                 goto error;
3928         isl_assert(set->ctx, set->dim == map->n_in, goto error);
3929         map = isl_map_intersect_domain(map, set);
3930         set = isl_map_range(map);
3931         return set;
3932 error:
3933         isl_set_free(set);
3934         isl_map_free(map);
3935         return NULL;
3936 }
3937
3938 /* There is no need to cow as removing empty parts doesn't change
3939  * the meaning of the set.
3940  */
3941 struct isl_map *isl_map_remove_empty_parts(struct isl_map *map)
3942 {
3943         int i;
3944
3945         if (!map)
3946                 return NULL;
3947
3948         for (i = map->n-1; i >= 0; --i) {
3949                 if (!F_ISSET(map->p[i], ISL_BASIC_MAP_EMPTY))
3950                         continue;
3951                 isl_basic_map_free(map->p[i]);
3952                 if (i != map->n-1) {
3953                         F_CLR(map, ISL_MAP_NORMALIZED);
3954                         map->p[i] = map->p[map->n-1];
3955                 }
3956                 map->n--;
3957         }
3958
3959         return map;
3960 }
3961
3962 struct isl_set *isl_set_remove_empty_parts(struct isl_set *set)
3963 {
3964         return (struct isl_set *)
3965                 isl_map_remove_empty_parts((struct isl_map *)set);
3966 }
3967
3968 struct isl_basic_set *isl_set_copy_basic_set(struct isl_set *set)
3969 {
3970         struct isl_basic_set *bset;
3971         if (!set || set->n == 0)
3972                 return NULL;
3973         bset = set->p[set->n-1];
3974         isl_assert(set->ctx, F_ISSET(bset, ISL_BASIC_SET_FINAL), return NULL);
3975         return isl_basic_set_copy(bset);
3976 }
3977
3978 struct isl_set *isl_set_drop_basic_set(struct isl_set *set,
3979                                                 struct isl_basic_set *bset)
3980 {
3981         int i;
3982
3983         if (!set || !bset)
3984                 goto error;
3985         for (i = set->n-1; i >= 0; --i) {
3986                 if (set->p[i] != bset)
3987                         continue;
3988                 set = isl_set_cow(set);
3989                 if (!set)
3990                         goto error;
3991                 isl_basic_set_free(set->p[i]);
3992                 if (i != set->n-1) {
3993                         F_CLR(set, ISL_SET_NORMALIZED);
3994                         set->p[i] = set->p[set->n-1];
3995                 }
3996                 set->n--;
3997                 return set;
3998         }
3999         isl_basic_set_free(bset);
4000         return set;
4001 error:
4002         isl_set_free(set);
4003         isl_basic_set_free(bset);
4004         return NULL;
4005 }
4006
4007 /* Given two _disjoint_ basic sets bset1 and bset2, check whether
4008  * for any common value of the parameters and dimensions preceding dim
4009  * in both basic sets, the values of dimension pos in bset1 are
4010  * smaller or larger then those in bset2.
4011  *
4012  * Returns
4013  *       1 if bset1 follows bset2
4014  *      -1 if bset1 precedes bset2
4015  *       0 if bset1 and bset2 are incomparable
4016  *      -2 if some error occurred.
4017  */
4018 int isl_basic_set_compare_at(struct isl_basic_set *bset1,
4019         struct isl_basic_set *bset2, int pos)
4020 {
4021         struct isl_basic_map *bmap1 = NULL;
4022         struct isl_basic_map *bmap2 = NULL;
4023         struct isl_ctx *ctx;
4024         struct isl_vec *obj;
4025         unsigned total;
4026         isl_int num, den;
4027         enum isl_lp_result res;
4028         int cmp;
4029
4030         if (!bset1 || !bset2)
4031                 return -2;
4032
4033         bmap1 = isl_basic_map_from_basic_set(isl_basic_set_copy(bset1),
4034                                                 pos, bset1->dim - pos);
4035         bmap2 = isl_basic_map_from_basic_set(isl_basic_set_copy(bset2),
4036                                                 pos, bset2->dim - pos);
4037         if (!bmap1 || !bmap2)
4038                 goto error;
4039         bmap1 = isl_basic_map_extend(bmap1, bmap1->nparam,
4040                         bmap1->n_in, bmap1->n_out + bmap2->n_out,
4041                         bmap2->n_div, bmap2->n_eq, bmap2->n_ineq);
4042         if (!bmap1)
4043                 goto error;
4044         total = bmap1->nparam + bmap1->n_in + bmap1->n_out + bmap1->n_div;
4045         ctx = bmap1->ctx;
4046         obj = isl_vec_alloc(ctx, total);
4047         isl_seq_clr(obj->block.data, total);
4048         isl_int_set_si(obj->block.data[bmap1->nparam+bmap1->n_in], 1);
4049         isl_int_set_si(obj->block.data[bmap1->nparam+bmap1->n_in+
4050                                         bmap1->n_out-bmap2->n_out], -1);
4051         if (!obj)
4052                 goto error;
4053         bmap1 = add_constraints(bmap1, bmap2, 0, bmap1->n_out - bmap2->n_out);
4054         isl_int_init(num);
4055         isl_int_init(den);
4056         res = isl_solve_lp(bmap1, 0, obj->block.data, ctx->one, &num, &den);
4057         if (res == isl_lp_empty)
4058                 cmp = 0;
4059         else if (res == isl_lp_ok && isl_int_is_pos(num))
4060                 cmp = 1;
4061         else if ((res == isl_lp_ok && isl_int_is_neg(num)) ||
4062                   res == isl_lp_unbounded)
4063                 cmp = -1;
4064         else
4065                 cmp = -2;
4066         isl_int_clear(num);
4067         isl_int_clear(den);
4068         isl_basic_map_free(bmap1);
4069         isl_vec_free(ctx, obj);
4070         return cmp;
4071 error:
4072         isl_basic_map_free(bmap1);
4073         isl_basic_map_free(bmap2);
4074         return -2;
4075 }
4076
4077 static int isl_basic_map_fast_has_fixed_var(struct isl_basic_map *bmap,
4078         unsigned pos, isl_int *val)
4079 {
4080         int i;
4081         int d;
4082         unsigned total;
4083
4084         if (!bmap)
4085                 return -1;
4086         total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
4087         for (i = 0, d = total-1; i < bmap->n_eq && d+1 > pos; ++i) {
4088                 for (; d+1 > pos; --d)
4089                         if (!isl_int_is_zero(bmap->eq[i][1+d]))
4090                                 break;
4091                 if (d != pos)
4092                         continue;
4093                 if (isl_seq_first_non_zero(bmap->eq[i]+1, d) != -1)
4094                         return 0;
4095                 if (isl_seq_first_non_zero(bmap->eq[i]+1+d+1, total-d-1) != -1)
4096                         return 0;
4097                 if (!isl_int_is_one(bmap->eq[i][1+d]))
4098                         return 0;
4099                 if (val)
4100                         isl_int_neg(*val, bmap->eq[i][0]);
4101                 return 1;
4102         }
4103         return 0;
4104 }
4105
4106 static int isl_map_fast_has_fixed_var(struct isl_map *map,
4107         unsigned pos, isl_int *val)
4108 {
4109         int i;
4110         isl_int v;
4111         isl_int tmp;
4112         int fixed;
4113
4114         if (!map)
4115                 return -1;
4116         if (map->n == 0)
4117                 return 0;
4118         if (map->n == 1)
4119                 return isl_basic_map_fast_has_fixed_var(map->p[0], pos, val);
4120         isl_int_init(v);
4121         isl_int_init(tmp);
4122         fixed = isl_basic_map_fast_has_fixed_var(map->p[0], pos, &v);
4123         for (i = 1; fixed == 1 && i < map->n; ++i) {
4124                 fixed = isl_basic_map_fast_has_fixed_var(map->p[i], pos, &tmp);
4125                 if (fixed == 1 && isl_int_ne(tmp, v))
4126                         fixed = 0;
4127         }
4128         if (val)
4129                 isl_int_set(*val, v);
4130         isl_int_clear(tmp);
4131         isl_int_clear(v);
4132         return fixed;
4133 }
4134
4135 static int isl_set_fast_has_fixed_var(struct isl_set *set, unsigned pos,
4136         isl_int *val)
4137 {
4138         return isl_map_fast_has_fixed_var((struct isl_map *)set, pos, val);
4139 }
4140
4141 /* Check if dimension dim has fixed value and if so and if val is not NULL,
4142  * then return this fixed value in *val.
4143  */
4144 int isl_set_fast_dim_is_fixed(struct isl_set *set, unsigned dim, isl_int *val)
4145 {
4146         return isl_set_fast_has_fixed_var(set, set->nparam + dim, val);
4147 }
4148
4149 /* Check if input variable in has fixed value and if so and if val is not NULL,
4150  * then return this fixed value in *val.
4151  */
4152 int isl_map_fast_input_is_fixed(struct isl_map *map, unsigned in, isl_int *val)
4153 {
4154         return isl_map_fast_has_fixed_var(map, map->nparam + in, val);
4155 }
4156
4157 /* Check if dimension dim has an (obvious) fixed lower bound and if so
4158  * and if val is not NULL, then return this lower bound in *val.
4159  */
4160 int isl_basic_set_fast_dim_has_fixed_lower_bound(struct isl_basic_set *bset,
4161         unsigned dim, isl_int *val)
4162 {
4163         int i, i_eq = -1, i_ineq = -1;
4164         isl_int *c;
4165         unsigned total;
4166
4167         if (!bset)
4168                 return -1;
4169         total = bset->nparam + bset->dim + bset->n_div;
4170         for (i = 0; i < bset->n_eq; ++i) {
4171                 if (isl_int_is_zero(bset->eq[i][1+bset->nparam+dim]))
4172                         continue;
4173                 if (i_eq != -1)
4174                         return 0;
4175                 i_eq = i;
4176         }
4177         for (i = 0; i < bset->n_ineq; ++i) {
4178                 if (!isl_int_is_pos(bset->ineq[i][1+bset->nparam+dim]))
4179                         continue;
4180                 if (i_eq != -1 || i_ineq != -1)
4181                         return 0;
4182                 i_ineq = i;
4183         }
4184         if (i_eq == -1 && i_ineq == -1)
4185                 return 0;
4186         c = i_eq != -1 ? bset->eq[i_eq] : bset->ineq[i_ineq];
4187         /* The coefficient should always be one due to normalization. */
4188         if (!isl_int_is_one(c[1+bset->nparam+dim]))
4189                 return 0;
4190         if (isl_seq_first_non_zero(c+1, bset->nparam+dim) != -1)
4191                 return 0;
4192         if (isl_seq_first_non_zero(c+1+bset->nparam+dim+1,
4193                                         total - bset->nparam - dim - 1) != -1)
4194                 return 0;
4195         if (val)
4196                 isl_int_neg(*val, c[0]);
4197         return 1;
4198 }
4199
4200 int isl_set_fast_dim_has_fixed_lower_bound(struct isl_set *set,
4201         unsigned dim, isl_int *val)
4202 {
4203         int i;
4204         isl_int v;
4205         isl_int tmp;
4206         int fixed;
4207
4208         if (!set)
4209                 return -1;
4210         if (set->n == 0)
4211                 return 0;
4212         if (set->n == 1)
4213                 return isl_basic_set_fast_dim_has_fixed_lower_bound(set->p[0],
4214                                                                 dim, val);
4215         isl_int_init(v);
4216         isl_int_init(tmp);
4217         fixed = isl_basic_set_fast_dim_has_fixed_lower_bound(set->p[0],
4218                                                                 dim, &v);
4219         for (i = 1; fixed == 1 && i < set->n; ++i) {
4220                 fixed = isl_basic_set_fast_dim_has_fixed_lower_bound(set->p[i],
4221                                                                 dim, &tmp);
4222                 if (fixed == 1 && isl_int_ne(tmp, v))
4223                         fixed = 0;
4224         }
4225         if (val)
4226                 isl_int_set(*val, v);
4227         isl_int_clear(tmp);
4228         isl_int_clear(v);
4229         return fixed;
4230 }
4231
4232 static struct isl_basic_set *isl_basic_set_reduce_using_equalities(
4233         struct isl_basic_set *bset, struct isl_basic_set *context)
4234 {
4235         int i;
4236         int *elim;
4237
4238         if (!bset || !context)
4239                 goto error;
4240
4241         bset = isl_basic_set_cow(bset);
4242         if (!bset)
4243                 goto error;
4244
4245         elim = isl_alloc_array(ctx, int, bset->dim);
4246         if (!elim)
4247                 goto error;
4248         set_compute_elimination_index(context, elim);
4249         for (i = 0; i < bset->n_eq; ++i)
4250                 set_reduced_using_equalities(bset->eq[i], bset->eq[i],
4251                                                         context, elim);
4252         for (i = 0; i < bset->n_ineq; ++i)
4253                 set_reduced_using_equalities(bset->ineq[i], bset->ineq[i],
4254                                                         context, elim);
4255         isl_basic_set_free(context);
4256         free(elim);
4257         bset = isl_basic_set_simplify(bset);
4258         bset = isl_basic_set_finalize(bset);
4259         return bset;
4260 error:
4261         isl_basic_set_free(bset);
4262         isl_basic_set_free(context);
4263         return NULL;
4264 }
4265
4266 static struct isl_basic_set *remove_shifted_constraints(
4267         struct isl_basic_set *bset, struct isl_basic_set *context)
4268 {
4269         unsigned int size;
4270         isl_int ***index;
4271         int bits;
4272         int k, h, l;
4273
4274         if (!bset)
4275                 return NULL;
4276
4277         size = round_up(4 * (context->n_ineq+1) / 3 - 1);
4278         bits = ffs(size) - 1;
4279         index = isl_alloc_array(ctx, isl_int **, size);
4280         memset(index, 0, size * sizeof(isl_int **));
4281         if (!index)
4282                 return bset;
4283
4284         for (k = 0; k < context->n_ineq; ++k) {
4285                 h = set_hash_index(index, size, bits, context, k);
4286                 index[h] = &context->ineq[k];
4287         }
4288         for (k = 0; k < bset->n_ineq; ++k) {
4289                 h = set_hash_index(index, size, bits, bset, k);
4290                 if (!index[h])
4291                         continue;
4292                 l = index[h] - &context->ineq[0];
4293                 if (isl_int_lt(bset->ineq[k][0], context->ineq[l][0]))
4294                         continue;
4295                 bset = isl_basic_set_cow(bset);
4296                 if (!bset)
4297                         goto error;
4298                 isl_basic_set_drop_inequality(bset, k);
4299                 --k;
4300         }
4301         free(index);
4302         return bset;
4303 error:
4304         free(index);
4305         return bset;
4306 }
4307
4308 /* Remove all information from bset that is redundant in the context
4309  * of context.  In particular, equalities that are linear combinations
4310  * of those in context are remobed.  Then the inequalities that are
4311  * redundant in the context of the equalities and inequalities of
4312  * context are removed.
4313  */
4314 static struct isl_basic_set *uset_gist(struct isl_basic_set *bset,
4315         struct isl_basic_set *context)
4316 {
4317         int i;
4318         isl_int opt;
4319         struct isl_basic_set *combined;
4320         struct isl_ctx *ctx;
4321
4322         if (!bset || !context)
4323                 goto error;
4324
4325         if (context->n_eq > 0) {
4326                 struct isl_mat *T;
4327                 struct isl_mat *T2;
4328                 struct isl_ctx *ctx = context->ctx;
4329                 struct isl_basic_set *reduced_context;
4330                 reduced_context = isl_basic_set_remove_equalities(
4331                                         isl_basic_set_copy(context), &T, &T2);
4332                 if (!reduced_context)
4333                         goto error;
4334                 bset = isl_basic_set_preimage(ctx, bset, T);
4335                 bset = uset_gist(bset, reduced_context);
4336                 bset = isl_basic_set_preimage(ctx, bset, T2);
4337                 bset = isl_basic_set_reduce_using_equalities(bset, context);
4338                 return bset;
4339         }
4340         if (!context->n_ineq)
4341                 goto done;
4342         bset = remove_shifted_constraints(bset, context);
4343         combined = isl_basic_set_extend(isl_basic_set_copy(bset),
4344                         0, bset->dim, 0, context->n_eq, context->n_ineq);
4345         context = set_add_constraints(combined, context, 0);
4346         if (!context)
4347                 goto error;
4348         ctx = context->ctx;
4349         isl_int_init(opt);
4350         for (i = bset->n_ineq-1; i >= 0; --i) {
4351                 int redundant;
4352                 set_swap_inequality(context, i, context->n_ineq-1);
4353                 context->n_ineq--;
4354                 redundant = isl_basic_set_constraint_is_redundant(&context,
4355                                 context->ineq[context->n_ineq], &opt, NULL);
4356                 if (redundant == -1) {
4357                         isl_int_clear(opt);
4358                         goto error;
4359                 }
4360                 if (F_ISSET(context, ISL_BASIC_MAP_EMPTY)) {
4361                         bset = isl_basic_set_set_to_empty(bset);
4362                         break;
4363                 }
4364                 context->n_ineq++;
4365                 set_swap_inequality(context, i, context->n_ineq-1);
4366                 if (redundant) {
4367                         isl_basic_set_drop_inequality(context, i);
4368                         isl_basic_set_drop_inequality(bset, i);
4369                 }
4370         }
4371         isl_int_clear(opt);
4372 done:
4373         isl_basic_set_free(context);
4374         return bset;
4375 error:
4376         isl_basic_set_free(context);
4377         isl_basic_set_free(bset);
4378         return NULL;
4379 }
4380
4381 struct isl_basic_map *isl_basic_map_gist(struct isl_basic_map *bmap,
4382         struct isl_basic_map *context)
4383 {
4384         struct isl_basic_set *bset;
4385
4386         if (!bmap || !context)
4387                 goto error;
4388
4389         context = isl_basic_map_align_divs(context, bmap);
4390         bmap = isl_basic_map_align_divs(bmap, context);
4391
4392         bset = uset_gist(isl_basic_map_underlying_set(isl_basic_map_copy(bmap)),
4393                          isl_basic_map_underlying_set(context));
4394
4395         return isl_basic_map_overlying_set(bset, bmap);
4396 error:
4397         isl_basic_map_free(bmap);
4398         isl_basic_map_free(context);
4399         return NULL;
4400 }
4401
4402 struct isl_map *isl_map_gist(struct isl_map *map, struct isl_basic_map *context)
4403 {
4404         int i;
4405
4406         map = isl_map_cow(map);
4407         if (!map || !context)
4408                 return NULL;
4409         isl_assert(map->ctx, map->nparam == context->nparam, goto error);
4410         isl_assert(map->ctx, map->n_in == context->n_in, goto error);
4411         isl_assert(map->ctx, map->n_out == context->n_out, goto error);
4412         for (i = 0; i < map->n; ++i)
4413                 context = isl_basic_map_align_divs(context, map->p[i]);
4414         for (i = 0; i < map->n; ++i) {
4415                 map->p[i] = isl_basic_map_gist(map->p[i],
4416                                                 isl_basic_map_copy(context));
4417                 if (!map->p[i])
4418                         goto error;
4419         }
4420         isl_basic_map_free(context);
4421         F_CLR(map, ISL_MAP_NORMALIZED);
4422         return map;
4423 error:
4424         isl_map_free(map);
4425         isl_basic_map_free(context);
4426         return NULL;
4427 }
4428
4429 struct isl_set *isl_set_gist(struct isl_set *set, struct isl_basic_set *context)
4430 {
4431         return (struct isl_set *)isl_map_gist((struct isl_map *)set,
4432                                         (struct isl_basic_map *)context);
4433 }
4434
4435 struct constraint {
4436         unsigned        size;
4437         isl_int         *c;
4438 };
4439
4440 static int qsort_constraint_cmp(const void *p1, const void *p2)
4441 {
4442         const struct constraint *c1 = (const struct constraint *)p1;
4443         const struct constraint *c2 = (const struct constraint *)p2;
4444         unsigned size = isl_min(c1->size, c2->size);
4445         return isl_seq_cmp(c1->c, c2->c, size);
4446 }
4447
4448 static struct isl_basic_map *isl_basic_map_sort_constraints(
4449         struct isl_basic_map *bmap)
4450 {
4451         int i;
4452         struct constraint *c;
4453         unsigned total;
4454
4455         if (!bmap)
4456                 return NULL;
4457         total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
4458         c = isl_alloc_array(bmap->ctx, struct constraint, bmap->n_ineq);
4459         if (!c)
4460                 goto error;
4461         for (i = 0; i < bmap->n_ineq; ++i) {
4462                 c[i].size = total;
4463                 c[i].c = bmap->ineq[i];
4464         }
4465         qsort(c, bmap->n_ineq, sizeof(struct constraint), qsort_constraint_cmp);
4466         for (i = 0; i < bmap->n_ineq; ++i)
4467                 bmap->ineq[i] = c[i].c;
4468         free(c);
4469         return bmap;
4470 error:
4471         isl_basic_map_free(bmap);
4472         return NULL;
4473 }
4474
4475 struct isl_basic_map *isl_basic_map_normalize(struct isl_basic_map *bmap)
4476 {
4477         if (!bmap)
4478                 return NULL;
4479         if (F_ISSET(bmap, ISL_BASIC_MAP_NORMALIZED))
4480                 return bmap;
4481         bmap = isl_basic_map_convex_hull(bmap);
4482         bmap = isl_basic_map_sort_constraints(bmap);
4483         F_SET(bmap, ISL_BASIC_MAP_NORMALIZED);
4484         return bmap;
4485 }
4486
4487 static int isl_basic_map_fast_cmp(const struct isl_basic_map *bmap1,
4488         const struct isl_basic_map *bmap2)
4489 {
4490         int i, cmp;
4491         unsigned total;
4492
4493         if (bmap1 == bmap2)
4494                 return 0;
4495         if (bmap1->nparam != bmap2->nparam)
4496                 return bmap1->nparam - bmap2->nparam;
4497         if (bmap1->n_in != bmap2->n_in)
4498                 return bmap1->n_in - bmap2->n_in;
4499         if (bmap1->n_out != bmap2->n_out)
4500                 return bmap1->n_out - bmap2->n_out;
4501         if (F_ISSET(bmap1, ISL_BASIC_MAP_EMPTY) &&
4502             F_ISSET(bmap2, ISL_BASIC_MAP_EMPTY))
4503                 return 0;
4504         if (F_ISSET(bmap1, ISL_BASIC_MAP_EMPTY))
4505                 return 1;
4506         if (F_ISSET(bmap2, ISL_BASIC_MAP_EMPTY))
4507                 return -1;
4508         if (bmap1->n_eq != bmap2->n_eq)
4509                 return bmap1->n_eq - bmap2->n_eq;
4510         if (bmap1->n_ineq != bmap2->n_ineq)
4511                 return bmap1->n_ineq - bmap2->n_ineq;
4512         if (bmap1->n_div != bmap2->n_div)
4513                 return bmap1->n_div - bmap2->n_div;
4514         total = bmap1->nparam + bmap1->n_in + bmap1->n_out + bmap1->n_div;
4515         for (i = 0; i < bmap1->n_eq; ++i) {
4516                 cmp = isl_seq_cmp(bmap1->eq[i], bmap2->eq[i], 1+total);
4517                 if (cmp)
4518                         return cmp;
4519         }
4520         for (i = 0; i < bmap1->n_ineq; ++i) {
4521                 cmp = isl_seq_cmp(bmap1->ineq[i], bmap2->ineq[i], 1+total);
4522                 if (cmp)
4523                         return cmp;
4524         }
4525         for (i = 0; i < bmap1->n_div; ++i) {
4526                 cmp = isl_seq_cmp(bmap1->div[i], bmap2->div[i], 1+1+total);
4527                 if (cmp)
4528                         return cmp;
4529         }
4530         return 0;
4531 }
4532
4533 static int isl_basic_map_fast_is_equal(struct isl_basic_map *bmap1,
4534         struct isl_basic_map *bmap2)
4535 {
4536         return isl_basic_map_fast_cmp(bmap1, bmap2) == 0;
4537 }
4538
4539 static int qsort_bmap_cmp(const void *p1, const void *p2)
4540 {
4541         const struct isl_basic_map *bmap1 = *(const struct isl_basic_map **)p1;
4542         const struct isl_basic_map *bmap2 = *(const struct isl_basic_map **)p2;
4543
4544         return isl_basic_map_fast_cmp(bmap1, bmap2);
4545 }
4546
4547 /* We normalize in place, but if anything goes wrong we need
4548  * to return NULL, so we need to make sure we don't change the
4549  * meaning of any possible other copies of map.
4550  */
4551 struct isl_map *isl_map_normalize(struct isl_map *map)
4552 {
4553         int i, j;
4554         struct isl_basic_map *bmap;
4555
4556         if (!map)
4557                 return NULL;
4558         if (F_ISSET(map, ISL_MAP_NORMALIZED))
4559                 return map;
4560         for (i = 0; i < map->n; ++i) {
4561                 bmap = isl_basic_map_normalize(isl_basic_map_copy(map->p[i]));
4562                 if (!bmap)
4563                         goto error;
4564                 isl_basic_map_free(map->p[i]);
4565                 map->p[i] = bmap;
4566         }
4567         qsort(map->p, map->n, sizeof(struct isl_basic_map *), qsort_bmap_cmp);
4568         F_SET(map, ISL_MAP_NORMALIZED);
4569         map = isl_map_remove_empty_parts(map);
4570         if (!map)
4571                 return NULL;
4572         for (i = map->n - 1; i >= 1; --i) {
4573                 if (!isl_basic_map_fast_is_equal(map->p[i-1], map->p[i]))
4574                         continue;
4575                 isl_basic_map_free(map->p[i-1]);
4576                 for (j = i; j < map->n; ++j)
4577                         map->p[j-1] = map->p[j];
4578                 map->n--;
4579         }
4580         return map;
4581 error:
4582         isl_map_free(map);
4583         return NULL;
4584
4585 }
4586
4587 struct isl_set *isl_set_normalize(struct isl_set *set)
4588 {
4589         return (struct isl_set *)isl_map_normalize((struct isl_map *)set);
4590 }
4591
4592 int isl_map_fast_is_equal(struct isl_map *map1, struct isl_map *map2)
4593 {
4594         int i;
4595         int equal;
4596
4597         if (!map1 || !map2)
4598                 return -1;
4599
4600         if (map1 == map2)
4601                 return 1;
4602         if (map1->nparam != map2->nparam)
4603                 return 0;
4604         if (map1->n_in != map2->n_in)
4605                 return 0;
4606         if (map1->n_out != map2->n_out)
4607                 return 0;
4608
4609         map1 = isl_map_copy(map1);
4610         map2 = isl_map_copy(map2);
4611         map1 = isl_map_normalize(map1);
4612         map2 = isl_map_normalize(map2);
4613         if (!map1 || !map2)
4614                 goto error;
4615         equal = map1->n == map2->n;
4616         for (i = 0; equal && i < map1->n; ++i) {
4617                 equal = isl_basic_map_fast_is_equal(map1->p[i], map2->p[i]);
4618                 if (equal < 0)
4619                         goto error;
4620         }
4621         isl_map_free(map1);
4622         isl_map_free(map2);
4623         return equal;
4624 error:
4625         isl_map_free(map1);
4626         isl_map_free(map2);
4627         return -1;
4628 }
4629
4630 int isl_set_fast_is_equal(struct isl_set *set1, struct isl_set *set2)
4631 {
4632         return isl_map_fast_is_equal((struct isl_map *)set1,
4633                                                 (struct isl_map *)set2);
4634 }
4635
4636 /* Return an interval that ranges from min to max (inclusive)
4637  */
4638 struct isl_basic_set *isl_basic_set_interval(struct isl_ctx *ctx,
4639         isl_int min, isl_int max)
4640 {
4641         int k;
4642         struct isl_basic_set *bset = NULL;
4643
4644         bset = isl_basic_set_alloc(ctx, 0, 1, 0, 0, 2);
4645         if (!bset)
4646                 goto error;
4647
4648         k = isl_basic_set_alloc_inequality(bset);
4649         if (k < 0)
4650                 goto error;
4651         isl_int_set_si(bset->ineq[k][1], 1);
4652         isl_int_neg(bset->ineq[k][0], min);
4653
4654         k = isl_basic_set_alloc_inequality(bset);
4655         if (k < 0)
4656                 goto error;
4657         isl_int_set_si(bset->ineq[k][1], -1);
4658         isl_int_set(bset->ineq[k][0], max);
4659
4660         return bset;
4661 error:
4662         isl_basic_set_free(bset);
4663         return NULL;
4664 }
4665
4666 /* Return the Cartesian product of the basic sets in list (in the given order).
4667  */
4668 struct isl_basic_set *isl_basic_set_product(struct isl_basic_set_list *list)
4669 {
4670         int i;
4671         unsigned dim;
4672         unsigned nparam;
4673         unsigned extra;
4674         unsigned n_eq;
4675         unsigned n_ineq;
4676         struct isl_basic_set *product = NULL;
4677
4678         if (!list)
4679                 goto error;
4680         isl_assert(list->ctx, list->n > 0, goto error);
4681         isl_assert(list->ctx, list->p[0], goto error);
4682         nparam = list->p[0]->nparam;
4683         dim = list->p[0]->dim;
4684         extra = list->p[0]->n_div;
4685         n_eq = list->p[0]->n_eq;
4686         n_ineq = list->p[0]->n_ineq;
4687         for (i = 1; i < list->n; ++i) {
4688                 isl_assert(list->ctx, list->p[i], goto error);
4689                 isl_assert(list->ctx, nparam == list->p[i]->nparam, goto error);
4690                 dim += list->p[i]->dim;
4691                 extra += list->p[i]->n_div;
4692                 n_eq += list->p[i]->n_eq;
4693                 n_ineq += list->p[i]->n_ineq;
4694         }
4695         product = isl_basic_set_alloc(list->ctx, nparam, dim, extra,
4696                                         n_eq, n_ineq);
4697         if (!product)
4698                 goto error;
4699         dim = 0;
4700         for (i = 0; i < list->n; ++i) {
4701                 set_add_constraints(product,
4702                                         isl_basic_set_copy(list->p[i]), dim);
4703                 dim += list->p[i]->dim;
4704         }
4705         isl_basic_set_list_free(list);
4706         return product;
4707 error:
4708         isl_basic_set_free(product);
4709         isl_basic_set_list_free(list);
4710         return NULL;
4711 }