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(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 (k != index[h]-1 &&
1041                     isl_seq_eq(bmap->ineq[k]+1,
1042                                bmap->ineq[index[h]-1]+1, total))
1043                         break;
1044         return h;
1045 }
1046
1047 static struct isl_basic_map *remove_duplicate_constraints(
1048         struct isl_basic_map *bmap, int *progress)
1049 {
1050         unsigned int size;
1051         int *index;
1052         int k, l, h;
1053         int bits;
1054         unsigned total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
1055         isl_int sum;
1056
1057         if (bmap->n_ineq <= 1)
1058                 return bmap;
1059
1060         size = round_up(4 * (bmap->n_ineq+1) / 3 - 1);
1061         bits = ffs(size) - 1;
1062         index = isl_alloc_array(ctx, int, size);
1063         memset(index, 0, size * sizeof(int));
1064         if (!index)
1065                 return bmap;
1066
1067         index[isl_seq_hash(bmap->ineq[0]+1, total, bits)] = 1;
1068         for (k = 1; k < bmap->n_ineq; ++k) {
1069                 h = hash_index(index, size, bits, bmap, k);
1070                 if (!index[h]) {
1071                         index[h] = k+1;
1072                         continue;
1073                 }
1074                 if (progress)
1075                         *progress = 1;
1076                 l = index[h] - 1;
1077                 if (isl_int_lt(bmap->ineq[k][0], bmap->ineq[l][0]))
1078                         swap_inequality(bmap, k, l);
1079                 isl_basic_map_drop_inequality(bmap, k);
1080                 --k;
1081         }
1082         isl_int_init(sum);
1083         for (k = 0; k < bmap->n_ineq-1; ++k) {
1084                 isl_seq_neg(bmap->ineq[k]+1, bmap->ineq[k]+1, total);
1085                 h = hash_index(index, size, bits, bmap, k);
1086                 isl_seq_neg(bmap->ineq[k]+1, bmap->ineq[k]+1, total);
1087                 if (!index[h])
1088                         continue;
1089                 l = index[h] - 1;
1090                 isl_int_add(sum, bmap->ineq[k][0], bmap->ineq[l][0]);
1091                 if (isl_int_is_pos(sum))
1092                         continue;
1093                 if (isl_int_is_zero(sum)) {
1094                         /* We need to break out of the loop after these
1095                          * changes since the contents of the hash
1096                          * will no longer be valid.
1097                          * Plus, we probably we want to regauss first.
1098                          */
1099                         isl_basic_map_drop_inequality(bmap, l);
1100                         isl_basic_map_inequality_to_equality(bmap, k);
1101                 } else
1102                         bmap = isl_basic_map_set_to_empty(bmap);
1103                 break;
1104         }
1105         isl_int_clear(sum);
1106
1107         free(index);
1108         return bmap;
1109 }
1110
1111 static void compute_elimination_index(struct isl_basic_map *bmap, int *elim)
1112 {
1113         int d, i;
1114         unsigned total;
1115
1116         total = bmap->nparam + bmap->n_in + bmap->n_out;
1117         for (d = 0; d < total; ++d)
1118                 elim[d] = -1;
1119         for (d = total - 1, i = 0; d >= 0 && i < bmap->n_eq; ++i) {
1120                 for (; d >= 0; --d) {
1121                         if (isl_int_is_zero(bmap->eq[i][1+d]))
1122                                 continue;
1123                         elim[d] = i;
1124                         break;
1125                 }
1126         }
1127 }
1128
1129 static int reduced_using_equalities(isl_int *c, struct isl_vec *v,
1130         struct isl_basic_map *bmap, int *elim)
1131 {
1132         int d, i;
1133         int copied = 0;
1134         unsigned total;
1135
1136         total = bmap->nparam + bmap->n_in + bmap->n_out;
1137         for (d = total - 1; d >= 0; --d) {
1138                 if (isl_int_is_zero(c[1+d]))
1139                         continue;
1140                 if (elim[d] == -1)
1141                         continue;
1142                 if (!copied) {
1143                         isl_seq_cpy(v->block.data, c, 1 + total);
1144                         copied = 1;
1145                 }
1146                 isl_seq_elim(v->block.data, bmap->eq[elim[d]],
1147                                 1 + d, 1 + total, NULL);
1148         }
1149         return copied;
1150 }
1151
1152 /* Quick check to see if two basic maps are disjoint.
1153  * In particular, we reduce the equalities and inequalities of
1154  * one basic map in the context of the equalities of the other
1155  * basic map and check if we get a contradiction.
1156  */
1157 int isl_basic_map_fast_is_disjoint(struct isl_basic_map *bmap1,
1158         struct isl_basic_map *bmap2)
1159 {
1160         struct isl_vec *v = NULL;
1161         int *elim = NULL;
1162         unsigned total;
1163         int d, i;
1164
1165         if (!bmap1 || !bmap2)
1166                 return -1;
1167         isl_assert(bmap1->ctx, bmap1->nparam == bmap2->nparam, return -1);
1168         isl_assert(bmap1->ctx, bmap1->n_in == bmap2->n_in, return -1);
1169         isl_assert(bmap1->ctx, bmap1->n_out == bmap2->n_out, return -1);
1170         if (bmap1->n_div || bmap2->n_div)
1171                 return 0;
1172         if (!bmap1->n_eq && !bmap2->n_eq)
1173                 return 0;
1174
1175         total = bmap1->nparam + bmap1->n_in + bmap1->n_out;
1176         if (total == 0)
1177                 return 0;
1178         v = isl_vec_alloc(bmap1->ctx, 1 + total);
1179         if (!v)
1180                 goto error;
1181         elim = isl_alloc_array(bmap1->ctx, int, total);
1182         if (!elim)
1183                 goto error;
1184         compute_elimination_index(bmap1, elim);
1185         for (i = 0; i < bmap2->n_eq; ++i) {
1186                 int reduced;
1187                 reduced = reduced_using_equalities(bmap2->eq[i], v, bmap1, elim);
1188                 if (reduced && !isl_int_is_zero(v->block.data[0]) &&
1189                     isl_seq_first_non_zero(v->block.data + 1, total) == -1)
1190                         goto disjoint;
1191         }
1192         for (i = 0; i < bmap2->n_ineq; ++i) {
1193                 int reduced;
1194                 reduced = reduced_using_equalities(bmap2->ineq[i], v, bmap1, elim);
1195                 if (reduced && isl_int_is_neg(v->block.data[0]) &&
1196                     isl_seq_first_non_zero(v->block.data + 1, total) == -1)
1197                         goto disjoint;
1198         }
1199         compute_elimination_index(bmap2, elim);
1200         for (i = 0; i < bmap1->n_ineq; ++i) {
1201                 int reduced;
1202                 reduced = reduced_using_equalities(bmap1->ineq[i], v, bmap2, elim);
1203                 if (reduced && isl_int_is_neg(v->block.data[0]) &&
1204                     isl_seq_first_non_zero(v->block.data + 1, total) == -1)
1205                         goto disjoint;
1206         }
1207         isl_vec_free(bmap1->ctx, v);
1208         free(elim);
1209         return 0;
1210 disjoint:
1211         isl_vec_free(bmap1->ctx, v);
1212         free(elim);
1213         return 1;
1214 error:
1215         isl_vec_free(bmap1->ctx, v);
1216         free(elim);
1217         return -1;
1218 }
1219
1220 int isl_basic_set_fast_is_disjoint(struct isl_basic_set *bset1,
1221         struct isl_basic_set *bset2)
1222 {
1223         return isl_basic_map_fast_is_disjoint((struct isl_basic_map *)bset1,
1224                                               (struct isl_basic_map *)bset2);
1225 }
1226
1227 int isl_map_fast_is_disjoint(struct isl_map *map1, struct isl_map *map2)
1228 {
1229         int i, j;
1230
1231         if (!map1 || !map2)
1232                 return -1;
1233
1234         if (isl_map_fast_is_equal(map1, map2))
1235                 return 0;
1236
1237         for (i = 0; i < map1->n; ++i) {
1238                 for (j = 0; j < map2->n; ++j) {
1239                         int d = isl_basic_map_fast_is_disjoint(map1->p[i],
1240                                                                map2->p[j]);
1241                         if (d != 1)
1242                                 return d;
1243                 }
1244         }
1245         return 1;
1246 }
1247
1248 int isl_set_fast_is_disjoint(struct isl_set *set1, struct isl_set *set2)
1249 {
1250         return isl_map_fast_is_disjoint((struct isl_map *)set1,
1251                                         (struct isl_map *)set2);
1252 }
1253
1254 static struct isl_basic_map *remove_duplicate_divs(
1255         struct isl_basic_map *bmap, int *progress)
1256 {
1257         unsigned int size;
1258         int *index;
1259         int k, l, h;
1260         int bits;
1261         struct isl_blk eq;
1262         unsigned total_var = bmap->nparam + bmap->n_in + bmap->n_out;
1263         unsigned total = total_var + bmap->n_div;
1264         struct isl_ctx *ctx;
1265
1266         if (bmap->n_div <= 1)
1267                 return bmap;
1268
1269         ctx = bmap->ctx;
1270         for (k = bmap->n_div - 1; k >= 0; --k)
1271                 if (!isl_int_is_zero(bmap->div[k][0]))
1272                         break;
1273         if (k <= 0)
1274                 return bmap;
1275
1276         size = round_up(4 * bmap->n_div / 3 - 1);
1277         bits = ffs(size) - 1;
1278         index = isl_alloc_array(ctx, int, size);
1279         memset(index, 0, size * sizeof(int));
1280         if (!index)
1281                 return bmap;
1282         eq = isl_blk_alloc(ctx, 1+total);
1283         if (isl_blk_is_error(eq))
1284                 goto out;
1285
1286         isl_seq_clr(eq.data, 1+total);
1287         index[isl_seq_hash(bmap->div[k], 2+total, bits)] = k + 1;
1288         for (--k; k >= 0; --k) {
1289                 u_int32_t hash;
1290
1291                 if (isl_int_is_zero(bmap->div[k][0]))
1292                         continue;
1293
1294                 hash = isl_seq_hash(bmap->div[k], 2+total, bits);
1295                 for (h = hash; index[h]; h = (h+1) % size)
1296                         if (isl_seq_eq(bmap->div[k],
1297                                        bmap->div[index[h]-1], 2+total))
1298                                 break;
1299                 if (index[h]) {
1300                         *progress = 1;
1301                         l = index[h] - 1;
1302                         isl_int_set_si(eq.data[1+total_var+k], -1);
1303                         isl_int_set_si(eq.data[1+total_var+l], 1);
1304                         eliminate_div(bmap, eq.data, l);
1305                         isl_int_set_si(eq.data[1+total_var+k], 0);
1306                         isl_int_set_si(eq.data[1+total_var+l], 0);
1307                 }
1308                 index[h] = k+1;
1309         }
1310
1311         isl_blk_free(ctx, eq);
1312 out:
1313         free(index);
1314         return bmap;
1315 }
1316
1317 /* Elimininate divs based on equalities
1318  */
1319 static struct isl_basic_map *eliminate_divs_eq(
1320                 struct isl_basic_map *bmap, int *progress)
1321 {
1322         int d;
1323         int i;
1324         int modified = 0;
1325         unsigned off;
1326
1327         if (!bmap)
1328                 return NULL;
1329
1330         off = 1 + bmap->nparam + bmap->n_in + bmap->n_out;
1331
1332         for (d = bmap->n_div - 1; d >= 0 ; --d) {
1333                 for (i = 0; i < bmap->n_eq; ++i) {
1334                         if (!isl_int_is_one(bmap->eq[i][off + d]) &&
1335                             !isl_int_is_negone(bmap->eq[i][off + d]))
1336                                 continue;
1337                         modified = 1;
1338                         *progress = 1;
1339                         eliminate_div(bmap, bmap->eq[i], d);
1340                         isl_basic_map_drop_equality(bmap, i);
1341                         break;
1342                 }
1343         }
1344         if (modified)
1345                 return eliminate_divs_eq(bmap, progress);
1346         return bmap;
1347 }
1348
1349 static struct isl_basic_map *normalize_constraints(struct isl_basic_map *bmap)
1350 {
1351         int i;
1352         isl_int gcd;
1353         unsigned total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
1354
1355         isl_int_init(gcd);
1356         for (i = bmap->n_eq - 1; i >= 0; --i) {
1357                 isl_seq_gcd(bmap->eq[i]+1, total, &gcd);
1358                 if (isl_int_is_zero(gcd)) {
1359                         if (!isl_int_is_zero(bmap->eq[i][0])) {
1360                                 bmap = isl_basic_map_set_to_empty(bmap);
1361                                 break;
1362                         }
1363                         isl_basic_map_drop_equality(bmap, i);
1364                         continue;
1365                 }
1366                 if (F_ISSET(bmap, ISL_BASIC_MAP_RATIONAL))
1367                         isl_int_gcd(gcd, gcd, bmap->eq[i][0]);
1368                 if (isl_int_is_one(gcd))
1369                         continue;
1370                 if (!isl_int_is_divisible_by(bmap->eq[i][0], gcd)) {
1371                         bmap = isl_basic_map_set_to_empty(bmap);
1372                         break;
1373                 }
1374                 isl_seq_scale_down(bmap->eq[i], bmap->eq[i], gcd, 1+total);
1375         }
1376
1377         for (i = bmap->n_ineq - 1; i >= 0; --i) {
1378                 isl_seq_gcd(bmap->ineq[i]+1, total, &gcd);
1379                 if (isl_int_is_zero(gcd)) {
1380                         if (isl_int_is_neg(bmap->ineq[i][0])) {
1381                                 bmap = isl_basic_map_set_to_empty(bmap);
1382                                 break;
1383                         }
1384                         isl_basic_map_drop_inequality(bmap, i);
1385                         continue;
1386                 }
1387                 if (F_ISSET(bmap, ISL_BASIC_MAP_RATIONAL))
1388                         isl_int_gcd(gcd, gcd, bmap->ineq[i][0]);
1389                 if (isl_int_is_one(gcd))
1390                         continue;
1391                 isl_int_fdiv_q(bmap->ineq[i][0], bmap->ineq[i][0], gcd);
1392                 isl_seq_scale_down(bmap->ineq[i]+1, bmap->ineq[i]+1, gcd, total);
1393         }
1394         isl_int_clear(gcd);
1395
1396         return bmap;
1397 }
1398
1399
1400 /* Eliminate the specified variables from the constraints using
1401  * Fourier-Motzkin.  The variables themselves are not removed.
1402  */
1403 struct isl_basic_map *isl_basic_map_eliminate_vars(
1404         struct isl_basic_map *bmap, int pos, unsigned n)
1405 {
1406         int d;
1407         int i, j, k;
1408         unsigned total;
1409         unsigned extra;
1410
1411         if (!bmap)
1412                 return NULL;
1413         total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
1414         extra = bmap->extra - bmap->n_div;
1415
1416         bmap = isl_basic_map_cow(bmap);
1417         for (d = pos + n - 1; d >= pos; --d) {
1418                 int n_lower, n_upper;
1419                 if (!bmap)
1420                         return NULL;
1421                 for (i = 0; i < bmap->n_eq; ++i) {
1422                         if (isl_int_is_zero(bmap->eq[i][1+d]))
1423                                 continue;
1424                         eliminate_var_using_equality(bmap, d, bmap->eq[i], NULL);
1425                         isl_basic_map_drop_equality(bmap, i);
1426                         break;
1427                 }
1428                 if (i < bmap->n_eq)
1429                         continue;
1430                 n_lower = 0;
1431                 n_upper = 0;
1432                 for (i = 0; i < bmap->n_ineq; ++i) {
1433                         if (isl_int_is_pos(bmap->ineq[i][1+d]))
1434                                 n_lower++;
1435                         else if (isl_int_is_neg(bmap->ineq[i][1+d]))
1436                                 n_upper++;
1437                 }
1438                 bmap = isl_basic_map_extend(bmap,
1439                                 bmap->nparam, bmap->n_in, bmap->n_out, 0,
1440                                 0, n_lower * n_upper);
1441                 for (i = bmap->n_ineq - 1; i >= 0; --i) {
1442                         int last;
1443                         if (isl_int_is_zero(bmap->ineq[i][1+d]))
1444                                 continue;
1445                         last = -1;
1446                         for (j = 0; j < i; ++j) {
1447                                 if (isl_int_is_zero(bmap->ineq[j][1+d]))
1448                                         continue;
1449                                 last = j;
1450                                 if (isl_int_sgn(bmap->ineq[i][1+d]) ==
1451                                     isl_int_sgn(bmap->ineq[j][1+d]))
1452                                         continue;
1453                                 k = isl_basic_map_alloc_inequality(bmap);
1454                                 if (k < 0)
1455                                         goto error;
1456                                 isl_seq_cpy(bmap->ineq[k], bmap->ineq[i],
1457                                                 1+total);
1458                                 isl_seq_clr(bmap->ineq[k]+1+total, extra);
1459                                 isl_seq_elim(bmap->ineq[k], bmap->ineq[j],
1460                                                 1+d, 1+total, NULL);
1461                         }
1462                         isl_basic_map_drop_inequality(bmap, i);
1463                         i = last + 1;
1464                 }
1465                 if (n_lower > 0 && n_upper > 0) {
1466                         bmap = normalize_constraints(bmap);
1467                         bmap = remove_duplicate_constraints(bmap, NULL);
1468                         bmap = isl_basic_map_gauss(bmap, NULL);
1469                         bmap = isl_basic_map_convex_hull(bmap);
1470                         if (!bmap)
1471                                 goto error;
1472                         if (F_ISSET(bmap, ISL_BASIC_MAP_EMPTY))
1473                                 break;
1474                 }
1475         }
1476         F_CLR(bmap, ISL_BASIC_MAP_NORMALIZED);
1477         return bmap;
1478 error:
1479         isl_basic_map_free(bmap);
1480         return NULL;
1481 }
1482
1483 struct isl_basic_set *isl_basic_set_eliminate_vars(
1484         struct isl_basic_set *bset, unsigned pos, unsigned n)
1485 {
1486         return (struct isl_basic_set *)isl_basic_map_eliminate_vars(
1487                         (struct isl_basic_map *)bset, pos, n);
1488 }
1489
1490 /* Eliminate the specified n dimensions starting at first from the
1491  * constraints using Fourier-Motzkin, The dimensions themselves
1492  * are not removed.
1493  */
1494 struct isl_set *isl_set_eliminate_dims(struct isl_set *set,
1495         unsigned first, unsigned n)
1496 {
1497         int i;
1498
1499         if (n == 0)
1500                 return set;
1501
1502         set = isl_set_cow(set);
1503         if (!set)
1504                 return NULL;
1505         isl_assert(set->ctx, first+n <= set->dim, goto error);
1506
1507         for (i = 0; i < set->n; ++i) {
1508                 set->p[i] = isl_basic_set_eliminate_vars(set->p[i],
1509                                                         set->nparam + first, n);
1510                 if (!set->p[i])
1511                         goto error;
1512         }
1513         return set;
1514 error:
1515         isl_set_free(set);
1516         return NULL;
1517 }
1518
1519 /* Project out n dimensions starting at first using Fourier-Motzkin */
1520 struct isl_set *isl_set_remove_dims(struct isl_set *set,
1521         unsigned first, unsigned n)
1522 {
1523         set = isl_set_eliminate_dims(set, first, n);
1524         set = isl_set_drop_dims(set, first, n);
1525         return set;
1526 }
1527
1528 /* Project out n inputs starting at first using Fourier-Motzkin */
1529 struct isl_map *isl_map_remove_inputs(struct isl_map *map,
1530         unsigned first, unsigned n)
1531 {
1532         int i;
1533
1534         if (n == 0)
1535                 return map;
1536
1537         map = isl_map_cow(map);
1538         if (!map)
1539                 return NULL;
1540         isl_assert(map->ctx, first+n <= map->n_in, goto error);
1541
1542         for (i = 0; i < map->n; ++i) {
1543                 map->p[i] = isl_basic_map_eliminate_vars(map->p[i],
1544                                                         map->nparam + first, n);
1545                 if (!map->p[i])
1546                         goto error;
1547         }
1548         map = isl_map_drop_inputs(map, first, n);
1549         return map;
1550 error:
1551         isl_map_free(map);
1552         return NULL;
1553 }
1554
1555 /* Project out n dimensions starting at first using Fourier-Motzkin */
1556 struct isl_basic_set *isl_basic_set_remove_dims(struct isl_basic_set *bset,
1557         unsigned first, unsigned n)
1558 {
1559         bset = isl_basic_set_eliminate_vars(bset, bset->nparam + first, n);
1560         bset = isl_basic_set_drop_dims(bset, first, n);
1561         return bset;
1562 }
1563
1564 /* Elimininate divs based on inequalities
1565  */
1566 static struct isl_basic_map *eliminate_divs_ineq(
1567                 struct isl_basic_map *bmap, int *progress)
1568 {
1569         int d;
1570         int i;
1571         unsigned off;
1572         struct isl_ctx *ctx;
1573
1574         if (!bmap)
1575                 return NULL;
1576
1577         ctx = bmap->ctx;
1578         off = 1 + bmap->nparam + bmap->n_in + bmap->n_out;
1579
1580         for (d = bmap->n_div - 1; d >= 0 ; --d) {
1581                 for (i = 0; i < bmap->n_eq; ++i)
1582                         if (!isl_int_is_zero(bmap->eq[i][off + d]))
1583                                 break;
1584                 if (i < bmap->n_eq)
1585                         continue;
1586                 for (i = 0; i < bmap->n_ineq; ++i)
1587                         if (isl_int_abs_gt(bmap->ineq[i][off + d], ctx->one))
1588                                 break;
1589                 if (i < bmap->n_ineq)
1590                         continue;
1591                 *progress = 1;
1592                 bmap = isl_basic_map_eliminate_vars(bmap, (off-1)+d, 1);
1593                 if (F_ISSET(bmap, ISL_BASIC_MAP_EMPTY))
1594                         break;
1595                 bmap = isl_basic_map_drop_div(bmap, d);
1596                 if (!bmap)
1597                         break;
1598         }
1599         return bmap;
1600 }
1601
1602 struct isl_basic_map *isl_basic_map_simplify(struct isl_basic_map *bmap)
1603 {
1604         int progress = 1;
1605         if (!bmap)
1606                 return NULL;
1607         while (progress) {
1608                 progress = 0;
1609                 bmap = normalize_constraints(bmap);
1610                 bmap = eliminate_divs_eq(bmap, &progress);
1611                 bmap = eliminate_divs_ineq(bmap, &progress);
1612                 bmap = isl_basic_map_gauss(bmap, &progress);
1613                 bmap = remove_duplicate_divs(bmap, &progress);
1614                 bmap = remove_duplicate_constraints(bmap, &progress);
1615         }
1616         return bmap;
1617 }
1618
1619 struct isl_basic_set *isl_basic_set_simplify(struct isl_basic_set *bset)
1620 {
1621         return (struct isl_basic_set *)
1622                 isl_basic_map_simplify((struct isl_basic_map *)bset);
1623 }
1624
1625 static void dump_term(struct isl_basic_map *bmap,
1626                         isl_int c, int pos, FILE *out)
1627 {
1628         unsigned in = bmap->n_in;
1629         unsigned dim = bmap->n_in + bmap->n_out;
1630         if (!pos)
1631                 isl_int_print(out, c, 0);
1632         else {
1633                 if (!isl_int_is_one(c))
1634                         isl_int_print(out, c, 0);
1635                 if (pos < 1 + bmap->nparam)
1636                         fprintf(out, "p%d", pos - 1);
1637                 else if (pos < 1 + bmap->nparam + in)
1638                         fprintf(out, "i%d", pos - 1 - bmap->nparam);
1639                 else if (pos < 1 + bmap->nparam + dim)
1640                         fprintf(out, "o%d", pos - 1 - bmap->nparam - in);
1641                 else
1642                         fprintf(out, "e%d", pos - 1 - bmap->nparam - dim);
1643         }
1644 }
1645
1646 static void dump_constraint_sign(struct isl_basic_map *bmap, isl_int *c,
1647                                 int sign, FILE *out)
1648 {
1649         int i;
1650         int first;
1651         unsigned len = 1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
1652         isl_int v;
1653
1654         isl_int_init(v);
1655         for (i = 0, first = 1; i < len; ++i) {
1656                 if (isl_int_sgn(c[i]) * sign <= 0)
1657                         continue;
1658                 if (!first)
1659                         fprintf(out, " + ");
1660                 first = 0;
1661                 isl_int_abs(v, c[i]);
1662                 dump_term(bmap, v, i, out);
1663         }
1664         isl_int_clear(v);
1665         if (first)
1666                 fprintf(out, "0");
1667 }
1668
1669 static void dump_constraint(struct isl_basic_map *bmap, isl_int *c,
1670                                 const char *op, FILE *out, int indent)
1671 {
1672         int i;
1673
1674         fprintf(out, "%*s", indent, "");
1675
1676         dump_constraint_sign(bmap, c, 1, out);
1677         fprintf(out, " %s ", op);
1678         dump_constraint_sign(bmap, c, -1, out);
1679
1680         fprintf(out, "\n");
1681
1682         for (i = bmap->n_div; i < bmap->extra; ++i) {
1683                 if (isl_int_is_zero(c[1+bmap->nparam+bmap->n_in+bmap->n_out+i]))
1684                         continue;
1685                 fprintf(out, "%*s", indent, "");
1686                 fprintf(out, "ERROR: unused div coefficient not zero\n");
1687         }
1688 }
1689
1690 static void dump_constraints(struct isl_basic_map *bmap,
1691                                 isl_int **c, unsigned n,
1692                                 const char *op, FILE *out, int indent)
1693 {
1694         int i;
1695
1696         for (i = 0; i < n; ++i)
1697                 dump_constraint(bmap, c[i], op, out, indent);
1698 }
1699
1700 static void dump_affine(struct isl_basic_map *bmap, isl_int *exp, FILE *out)
1701 {
1702         int j;
1703         int first = 1;
1704         unsigned total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
1705
1706         for (j = 0; j < 1 + total; ++j) {
1707                 if (isl_int_is_zero(exp[j]))
1708                         continue;
1709                 if (!first && isl_int_is_pos(exp[j]))
1710                         fprintf(out, "+");
1711                 dump_term(bmap, exp[j], j, out);
1712                 first = 0;
1713         }
1714 }
1715
1716 static void dump(struct isl_basic_map *bmap, FILE *out, int indent)
1717 {
1718         int i;
1719
1720         dump_constraints(bmap, bmap->eq, bmap->n_eq, "=", out, indent);
1721         dump_constraints(bmap, bmap->ineq, bmap->n_ineq, ">=", out, indent);
1722
1723         for (i = 0; i < bmap->n_div; ++i) {
1724                 fprintf(out, "%*s", indent, "");
1725                 fprintf(out, "e%d = [(", i);
1726                 dump_affine(bmap, bmap->div[i]+1, out);
1727                 fprintf(out, ")/");
1728                 isl_int_print(out, bmap->div[i][0], 0);
1729                 fprintf(out, "]\n");
1730         }
1731 }
1732
1733 void isl_basic_set_dump(struct isl_basic_set *bset, FILE *out, int indent)
1734 {
1735         if (!bset) {
1736                 fprintf(out, "null basic set\n");
1737                 return;
1738         }
1739
1740         fprintf(out, "%*s", indent, "");
1741         fprintf(out, "ref: %d, nparam: %d, dim: %d, extra: %d, flags: %x\n",
1742                         bset->ref, bset->nparam, bset->dim, bset->extra,
1743                         bset->flags);
1744         dump((struct isl_basic_map *)bset, out, indent);
1745 }
1746
1747 void isl_basic_map_dump(struct isl_basic_map *bmap, FILE *out, int indent)
1748 {
1749         if (!bmap) {
1750                 fprintf(out, "null basic map\n");
1751                 return;
1752         }
1753
1754         fprintf(out, "%*s", indent, "");
1755         fprintf(out, "ref: %d, nparam: %d, in: %d, out: %d, extra: %d, flags: %x\n",
1756                 bmap->ref,
1757                 bmap->nparam, bmap->n_in, bmap->n_out, bmap->extra, bmap->flags);
1758         dump(bmap, out, indent);
1759 }
1760
1761 int isl_inequality_negate(struct isl_basic_map *bmap, unsigned pos)
1762 {
1763         unsigned total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
1764         if (!bmap)
1765                 return -1;
1766         isl_assert(bmap->ctx, pos < bmap->n_ineq, return -1);
1767         isl_seq_neg(bmap->ineq[pos], bmap->ineq[pos], 1 + total);
1768         isl_int_sub_ui(bmap->ineq[pos][0], bmap->ineq[pos][0], 1);
1769         F_CLR(bmap, ISL_BASIC_MAP_NORMALIZED);
1770         return 0;
1771 }
1772
1773 struct isl_set *isl_set_alloc(struct isl_ctx *ctx,
1774                 unsigned nparam, unsigned dim, int n, unsigned flags)
1775 {
1776         struct isl_set *set;
1777
1778         isl_assert(ctx, n >= 0, return NULL);
1779         set = isl_alloc(ctx, struct isl_set,
1780                         sizeof(struct isl_set) +
1781                         n * sizeof(struct isl_basic_set *));
1782         if (!set)
1783                 return NULL;
1784
1785         set->ctx = ctx;
1786         isl_ctx_ref(ctx);
1787         set->ref = 1;
1788         set->size = n;
1789         set->n = 0;
1790         set->nparam = nparam;
1791         set->zero = 0;
1792         set->dim = dim;
1793         set->flags = flags;
1794         return set;
1795 }
1796
1797 struct isl_set *isl_set_dup(struct isl_set *set)
1798 {
1799         int i;
1800         struct isl_set *dup;
1801
1802         if (!set)
1803                 return NULL;
1804
1805         dup = isl_set_alloc(set->ctx, set->nparam, set->dim, set->n, set->flags);
1806         if (!dup)
1807                 return NULL;
1808         for (i = 0; i < set->n; ++i)
1809                 dup = isl_set_add(dup,
1810                                   isl_basic_set_dup(set->p[i]));
1811         return dup;
1812 }
1813
1814 struct isl_set *isl_set_from_basic_set(struct isl_basic_set *bset)
1815 {
1816         struct isl_set *set;
1817
1818         if (!bset)
1819                 return NULL;
1820
1821         set = isl_set_alloc(bset->ctx, bset->nparam, bset->dim, 1, ISL_MAP_DISJOINT);
1822         if (!set) {
1823                 isl_basic_set_free(bset);
1824                 return NULL;
1825         }
1826         return isl_set_add(set, bset);
1827 }
1828
1829 struct isl_map *isl_map_from_basic_map(struct isl_basic_map *bmap)
1830 {
1831         struct isl_map *map;
1832
1833         if (!bmap)
1834                 return NULL;
1835
1836         map = isl_map_alloc(bmap->ctx, bmap->nparam, bmap->n_in, bmap->n_out, 1,
1837                                 ISL_MAP_DISJOINT);
1838         if (!map) {
1839                 isl_basic_map_free(bmap);
1840                 return NULL;
1841         }
1842         return isl_map_add(map, bmap);
1843 }
1844
1845 struct isl_set *isl_set_add(struct isl_set *set, struct isl_basic_set *bset)
1846 {
1847         if (!bset || !set)
1848                 goto error;
1849         isl_assert(set->ctx, set->nparam == bset->nparam, goto error);
1850         isl_assert(set->ctx, set->dim == bset->dim, goto error);
1851         isl_assert(set->ctx, set->n < set->size, goto error);
1852         set->p[set->n] = bset;
1853         set->n++;
1854         return set;
1855 error:
1856         if (set)
1857                 isl_set_free(set);
1858         if (bset)
1859                 isl_basic_set_free(bset);
1860         return NULL;
1861 }
1862
1863 void isl_set_free(struct isl_set *set)
1864 {
1865         int i;
1866
1867         if (!set)
1868                 return;
1869
1870         if (--set->ref > 0)
1871                 return;
1872
1873         isl_ctx_deref(set->ctx);
1874         for (i = 0; i < set->n; ++i)
1875                 isl_basic_set_free(set->p[i]);
1876         free(set);
1877 }
1878
1879 void isl_set_dump(struct isl_set *set, FILE *out, int indent)
1880 {
1881         int i;
1882
1883         if (!set) {
1884                 fprintf(out, "null set\n");
1885                 return;
1886         }
1887
1888         fprintf(out, "%*s", indent, "");
1889         fprintf(out, "ref: %d, n: %d, nparam: %d, dim: %d, flags: %x\n",
1890                         set->ref, set->n, set->nparam, set->dim, set->flags);
1891         for (i = 0; i < set->n; ++i) {
1892                 fprintf(out, "%*s", indent, "");
1893                 fprintf(out, "basic set %d:\n", i);
1894                 isl_basic_set_dump(set->p[i], out, indent+4);
1895         }
1896 }
1897
1898 void isl_map_dump(struct isl_map *map, FILE *out, int indent)
1899 {
1900         int i;
1901
1902         if (!map) {
1903                 fprintf(out, "null map\n");
1904                 return;
1905         }
1906
1907         fprintf(out, "%*s", indent, "");
1908         fprintf(out, "ref: %d, n: %d, nparam: %d, in: %d, out: %d, flags: %x\n",
1909                         map->ref, map->n, map->nparam, map->n_in, map->n_out,
1910                         map->flags);
1911         for (i = 0; i < map->n; ++i) {
1912                 fprintf(out, "%*s", indent, "");
1913                 fprintf(out, "basic map %d:\n", i);
1914                 isl_basic_map_dump(map->p[i], out, indent+4);
1915         }
1916 }
1917
1918 struct isl_basic_map *isl_basic_map_intersect_domain(
1919                 struct isl_basic_map *bmap, struct isl_basic_set *bset)
1920 {
1921         struct isl_basic_map *bmap_domain;
1922
1923         if (!bmap || !bset)
1924                 goto error;
1925
1926         isl_assert(bset->ctx, bset->dim == bmap->n_in, goto error);
1927         isl_assert(bset->ctx, bset->nparam == bmap->nparam, goto error);
1928
1929         bmap = isl_basic_map_extend(bmap,
1930                         bset->nparam, bmap->n_in, bmap->n_out,
1931                         bset->n_div, bset->n_eq, bset->n_ineq);
1932         if (!bmap)
1933                 goto error;
1934         bmap_domain = isl_basic_map_from_basic_set(bset,
1935                                                 bset->dim, 0);
1936         bmap = add_constraints(bmap, bmap_domain, 0, 0);
1937
1938         bmap = isl_basic_map_simplify(bmap);
1939         return isl_basic_map_finalize(bmap);
1940 error:
1941         isl_basic_map_free(bmap);
1942         isl_basic_set_free(bset);
1943         return NULL;
1944 }
1945
1946 struct isl_basic_map *isl_basic_map_intersect_range(
1947                 struct isl_basic_map *bmap, struct isl_basic_set *bset)
1948 {
1949         struct isl_basic_map *bmap_range;
1950
1951         if (!bmap || !bset)
1952                 goto error;
1953
1954         isl_assert(bset->ctx, bset->dim == bmap->n_out, goto error);
1955         isl_assert(bset->ctx, bset->nparam == bmap->nparam, goto error);
1956
1957         bmap = isl_basic_map_extend(bmap,
1958                         bset->nparam, bmap->n_in, bmap->n_out,
1959                         bset->n_div, bset->n_eq, bset->n_ineq);
1960         if (!bmap)
1961                 goto error;
1962         bmap_range = isl_basic_map_from_basic_set(bset,
1963                                                 0, bset->dim);
1964         bmap = add_constraints(bmap, bmap_range, 0, 0);
1965
1966         bmap = isl_basic_map_simplify(bmap);
1967         return isl_basic_map_finalize(bmap);
1968 error:
1969         isl_basic_map_free(bmap);
1970         isl_basic_set_free(bset);
1971         return NULL;
1972 }
1973
1974 struct isl_basic_map *isl_basic_map_intersect(
1975                 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
1976 {
1977         if (!bmap1 || !bmap2)
1978                 goto error;
1979
1980         isl_assert(bmap1->ctx, bmap1->nparam == bmap2->nparam, goto error);
1981         isl_assert(bmap1->ctx, bmap1->n_in == bmap2->n_in, goto error);
1982         isl_assert(bmap1->ctx, bmap1->n_out == bmap2->n_out, goto error);
1983
1984         bmap1 = isl_basic_map_extend(bmap1,
1985                         bmap1->nparam, bmap1->n_in, bmap1->n_out,
1986                         bmap2->n_div, bmap2->n_eq, bmap2->n_ineq);
1987         if (!bmap1)
1988                 goto error;
1989         bmap1 = add_constraints(bmap1, bmap2, 0, 0);
1990
1991         bmap1 = isl_basic_map_simplify(bmap1);
1992         return isl_basic_map_finalize(bmap1);
1993 error:
1994         isl_basic_map_free(bmap1);
1995         isl_basic_map_free(bmap2);
1996         return NULL;
1997 }
1998
1999 struct isl_basic_set *isl_basic_set_intersect(
2000                 struct isl_basic_set *bset1, struct isl_basic_set *bset2)
2001 {
2002         return (struct isl_basic_set *)
2003                 isl_basic_map_intersect(
2004                         (struct isl_basic_map *)bset1,
2005                         (struct isl_basic_map *)bset2);
2006 }
2007
2008 struct isl_map *isl_map_intersect(struct isl_map *map1, struct isl_map *map2)
2009 {
2010         unsigned flags = 0;
2011         struct isl_map *result;
2012         int i, j;
2013
2014         if (!map1 || !map2)
2015                 goto error;
2016
2017         if (F_ISSET(map1, ISL_MAP_DISJOINT) &&
2018             F_ISSET(map2, ISL_MAP_DISJOINT))
2019                 FL_SET(flags, ISL_MAP_DISJOINT);
2020
2021         result = isl_map_alloc(map1->ctx, map1->nparam, map1->n_in, map1->n_out,
2022                                 map1->n * map2->n, flags);
2023         if (!result)
2024                 goto error;
2025         for (i = 0; i < map1->n; ++i)
2026                 for (j = 0; j < map2->n; ++j) {
2027                         struct isl_basic_map *part;
2028                         part = isl_basic_map_intersect(
2029                                     isl_basic_map_copy(map1->p[i]),
2030                                     isl_basic_map_copy(map2->p[j]));
2031                         if (isl_basic_map_is_empty(part))
2032                                 isl_basic_map_free(part);
2033                         else
2034                                 result = isl_map_add(result, part);
2035                         if (!result)
2036                                 goto error;
2037                 }
2038         isl_map_free(map1);
2039         isl_map_free(map2);
2040         return result;
2041 error:
2042         isl_map_free(map1);
2043         isl_map_free(map2);
2044         return NULL;
2045 }
2046
2047 struct isl_set *isl_set_intersect(struct isl_set *set1, struct isl_set *set2)
2048 {
2049         return (struct isl_set *)
2050                 isl_map_intersect((struct isl_map *)set1,
2051                                   (struct isl_map *)set2);
2052 }
2053
2054 struct isl_basic_map *isl_basic_map_reverse(struct isl_basic_map *bmap)
2055 {
2056         struct isl_basic_set *bset;
2057         unsigned in;
2058
2059         if (!bmap)
2060                 return NULL;
2061         bmap = isl_basic_map_cow(bmap);
2062         if (!bmap)
2063                 return NULL;
2064         in = bmap->n_in;
2065         bset = isl_basic_set_from_basic_map(bmap);
2066         bset = isl_basic_set_swap_vars(bset, in);
2067         return isl_basic_map_from_basic_set(bset, bset->dim-in, in);
2068 }
2069
2070 /* Turn final n dimensions into existentially quantified variables.
2071  */
2072 struct isl_basic_set *isl_basic_set_project_out(
2073                 struct isl_basic_set *bset, unsigned n, unsigned flags)
2074 {
2075         int i;
2076         size_t row_size;
2077         isl_int **new_div;
2078         isl_int *old;
2079
2080         if (!bset)
2081                 return NULL;
2082
2083         isl_assert(bset->ctx, n <= bset->dim, goto error);
2084
2085         if (n == 0)
2086                 return bset;
2087
2088         bset = isl_basic_set_cow(bset);
2089
2090         row_size = 1 + bset->nparam + bset->dim + bset->extra;
2091         old = bset->block2.data;
2092         bset->block2 = isl_blk_extend(bset->ctx, bset->block2,
2093                                         (bset->extra + n) * (1 + row_size));
2094         if (!bset->block2.data)
2095                 goto error;
2096         new_div = isl_alloc_array(ctx, isl_int *, bset->extra + n);
2097         if (!new_div)
2098                 goto error;
2099         for (i = 0; i < n; ++i) {
2100                 new_div[i] = bset->block2.data +
2101                                 (bset->extra + i) * (1 + row_size);
2102                 isl_seq_clr(new_div[i], 1 + row_size);
2103         }
2104         for (i = 0; i < bset->extra; ++i)
2105                 new_div[n + i] = bset->block2.data + (bset->div[i] - old);
2106         free(bset->div);
2107         bset->div = new_div;
2108         bset->n_div += n;
2109         bset->extra += n;
2110         bset->dim -= n;
2111         bset = isl_basic_set_simplify(bset);
2112         return isl_basic_set_finalize(bset);
2113 error:
2114         isl_basic_set_free(bset);
2115         return NULL;
2116 }
2117
2118 struct isl_basic_map *isl_basic_map_apply_range(
2119                 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
2120 {
2121         struct isl_basic_set *bset;
2122         unsigned n_in, n_out;
2123
2124         if (!bmap1 || !bmap2)
2125                 goto error;
2126
2127         isl_assert(bmap->ctx, bmap1->n_out == bmap2->n_in, goto error);
2128         isl_assert(bmap->ctx, bmap1->nparam == bmap2->nparam, goto error);
2129
2130         n_in = bmap1->n_in;
2131         n_out = bmap2->n_out;
2132
2133         bmap2 = isl_basic_map_reverse(bmap2);
2134         if (!bmap2)
2135                 goto error;
2136         bmap1 = isl_basic_map_extend(bmap1, bmap1->nparam,
2137                         bmap1->n_in + bmap2->n_in, bmap1->n_out,
2138                         bmap2->extra, bmap2->n_eq, bmap2->n_ineq);
2139         if (!bmap1)
2140                 goto error;
2141         bmap1 = add_constraints(bmap1, bmap2, bmap1->n_in - bmap2->n_in, 0);
2142         bmap1 = isl_basic_map_simplify(bmap1);
2143         bset = isl_basic_set_from_basic_map(bmap1);
2144         bset = isl_basic_set_project_out(bset,
2145                                                 bset->dim - (n_in + n_out), 0);
2146         return isl_basic_map_from_basic_set(bset, n_in, n_out);
2147 error:
2148         isl_basic_map_free(bmap1);
2149         isl_basic_map_free(bmap2);
2150         return NULL;
2151 }
2152
2153 struct isl_basic_set *isl_basic_set_apply(
2154                 struct isl_basic_set *bset, struct isl_basic_map *bmap)
2155 {
2156         if (!bset || !bmap)
2157                 goto error;
2158
2159         isl_assert(bset->ctx, bset->dim == bmap->n_in, goto error);
2160
2161         return (struct isl_basic_set *)
2162                 isl_basic_map_apply_range((struct isl_basic_map *)bset, bmap);
2163 error:
2164         isl_basic_set_free(bset);
2165         isl_basic_map_free(bmap);
2166         return NULL;
2167 }
2168
2169 struct isl_basic_map *isl_basic_map_apply_domain(
2170                 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
2171 {
2172         if (!bmap1 || !bmap2)
2173                 goto error;
2174
2175         isl_assert(ctx, bmap1->n_in == bmap2->n_in, goto error);
2176         isl_assert(ctx, bmap1->nparam == bmap2->nparam, goto error);
2177
2178         bmap1 = isl_basic_map_reverse(bmap1);
2179         bmap1 = isl_basic_map_apply_range(bmap1, bmap2);
2180         return isl_basic_map_reverse(bmap1);
2181 error:
2182         isl_basic_map_free(bmap1);
2183         isl_basic_map_free(bmap2);
2184         return NULL;
2185 }
2186
2187 static struct isl_basic_map *var_equal(struct isl_ctx *ctx,
2188                 struct isl_basic_map *bmap, unsigned pos)
2189 {
2190         int i;
2191         i = isl_basic_map_alloc_equality(bmap);
2192         if (i < 0)
2193                 goto error;
2194         isl_seq_clr(bmap->eq[i],
2195                     1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->extra);
2196         isl_int_set_si(bmap->eq[i][1+bmap->nparam+pos], -1);
2197         isl_int_set_si(bmap->eq[i][1+bmap->nparam+bmap->n_in+pos], 1);
2198         return isl_basic_map_finalize(bmap);
2199 error:
2200         isl_basic_map_free(bmap);
2201         return NULL;
2202 }
2203
2204 static struct isl_basic_map *var_more(struct isl_ctx *ctx,
2205                 struct isl_basic_map *bmap, unsigned pos)
2206 {
2207         int i;
2208         i = isl_basic_map_alloc_inequality(bmap);
2209         if (i < 0)
2210                 goto error;
2211         isl_seq_clr(bmap->ineq[i],
2212                     1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->extra);
2213         isl_int_set_si(bmap->ineq[i][0], -1);
2214         isl_int_set_si(bmap->ineq[i][1+bmap->nparam+pos], -1);
2215         isl_int_set_si(bmap->ineq[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_less(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 struct isl_basic_map *isl_basic_map_equal(struct isl_ctx *ctx,
2241                 unsigned nparam, unsigned in, unsigned out, unsigned n_equal)
2242 {
2243         int i;
2244         struct isl_basic_map *bmap;
2245         bmap = isl_basic_map_alloc(ctx, nparam, in, out, 0, n_equal, 0);
2246         if (!bmap)
2247                 return NULL;
2248         for (i = 0; i < n_equal && bmap; ++i)
2249                 bmap = var_equal(ctx, bmap, i);
2250         return isl_basic_map_finalize(bmap);
2251 }
2252
2253 struct isl_basic_map *isl_basic_map_less_at(struct isl_ctx *ctx,
2254                 unsigned nparam, unsigned in, unsigned out, unsigned pos)
2255 {
2256         int i;
2257         struct isl_basic_map *bmap;
2258         bmap = isl_basic_map_alloc(ctx, nparam, in, out, 0, pos, 1);
2259         if (!bmap)
2260                 return NULL;
2261         for (i = 0; i < pos && bmap; ++i)
2262                 bmap = var_equal(ctx, bmap, i);
2263         if (bmap)
2264                 bmap = var_less(ctx, bmap, pos);
2265         return isl_basic_map_finalize(bmap);
2266 }
2267
2268 struct isl_basic_map *isl_basic_map_more_at(struct isl_ctx *ctx,
2269                 unsigned nparam, unsigned in, unsigned out, unsigned pos)
2270 {
2271         int i;
2272         struct isl_basic_map *bmap;
2273         bmap = isl_basic_map_alloc(ctx, nparam, in, out, 0, pos, 1);
2274         if (!bmap)
2275                 return NULL;
2276         for (i = 0; i < pos && bmap; ++i)
2277                 bmap = var_equal(ctx, bmap, i);
2278         if (bmap)
2279                 bmap = var_more(ctx, bmap, pos);
2280         return isl_basic_map_finalize(bmap);
2281 }
2282
2283 struct isl_basic_map *isl_basic_map_from_basic_set(
2284                 struct isl_basic_set *bset, unsigned n_in, unsigned n_out)
2285 {
2286         struct isl_basic_map *bmap;
2287
2288         bset = isl_basic_set_cow(bset);
2289         if (!bset)
2290                 return NULL;
2291
2292         isl_assert(ctx, bset->dim == n_in + n_out, goto error);
2293         bmap = (struct isl_basic_map *) bset;
2294         bmap->n_in = n_in;
2295         bmap->n_out = n_out;
2296         return isl_basic_map_finalize(bmap);
2297 error:
2298         isl_basic_set_free(bset);
2299         return NULL;
2300 }
2301
2302 struct isl_basic_set *isl_basic_set_from_basic_map(struct isl_basic_map *bmap)
2303 {
2304         if (!bmap)
2305                 goto error;
2306         if (bmap->n_in == 0)
2307                 return (struct isl_basic_set *)bmap;
2308         bmap = isl_basic_map_cow(bmap);
2309         if (!bmap)
2310                 goto error;
2311         bmap->n_out += bmap->n_in;
2312         bmap->n_in = 0;
2313         bmap = isl_basic_map_finalize(bmap);
2314         return (struct isl_basic_set *)bmap;
2315 error:
2316         return NULL;
2317 }
2318
2319 /* For a div d = floor(f/m), add the constraints
2320  *
2321  *              f - m d >= 0
2322  *              -(f-(n-1)) + m d >= 0
2323  *
2324  * Note that the second constraint is the negation of
2325  *
2326  *              f - m d >= n
2327  */
2328 static int add_div_constraints(struct isl_basic_map *bmap, unsigned div)
2329 {
2330         int i, j;
2331         unsigned div_pos = 1 + bmap->nparam + bmap->n_in + bmap->n_out + div;
2332         unsigned total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->extra;
2333
2334         i = isl_basic_map_alloc_inequality(bmap);
2335         if (i < 0)
2336                 return -1;
2337         isl_seq_cpy(bmap->ineq[i], bmap->div[div]+1, 1+total);
2338         isl_int_neg(bmap->ineq[i][div_pos], bmap->div[div][0]);
2339
2340         j = isl_basic_map_alloc_inequality(bmap);
2341         if (j < 0)
2342                 return -1;
2343         isl_seq_neg(bmap->ineq[j], bmap->ineq[i], 1 + total);
2344         isl_int_add(bmap->ineq[j][0], bmap->ineq[j][0], bmap->ineq[j][div_pos]);
2345         isl_int_sub_ui(bmap->ineq[j][0], bmap->ineq[j][0], 1);
2346         return j;
2347 }
2348
2349 struct isl_basic_set *isl_basic_map_underlying_set(
2350                 struct isl_basic_map *bmap)
2351 {
2352         if (!bmap)
2353                 goto error;
2354         if (bmap->nparam == 0 && bmap->n_in == 0 && bmap->n_div == 0)
2355                 return (struct isl_basic_set *)bmap;
2356         bmap = isl_basic_map_cow(bmap);
2357         if (!bmap)
2358                 goto error;
2359         bmap->n_out += bmap->nparam + bmap->n_in + bmap->n_div;
2360         bmap->nparam = 0;
2361         bmap->n_in = 0;
2362         bmap->extra -= bmap->n_div;
2363         bmap->n_div = 0;
2364         bmap = isl_basic_map_finalize(bmap);
2365         return (struct isl_basic_set *)bmap;
2366 error:
2367         return NULL;
2368 }
2369
2370 struct isl_basic_map *isl_basic_map_overlying_set(
2371         struct isl_basic_set *bset, struct isl_basic_map *like)
2372 {
2373         struct isl_basic_map *bmap;
2374         struct isl_ctx *ctx;
2375         unsigned total;
2376         int i, k;
2377
2378         if (!bset || !like)
2379                 goto error;
2380         ctx = bset->ctx;
2381         isl_assert(ctx, bset->dim ==
2382                         like->nparam + like->n_in + like->n_out + like->n_div,
2383                         goto error);
2384         if (like->nparam == 0 && like->n_in == 0 && like->n_div == 0) {
2385                 isl_basic_map_free(like);
2386                 return (struct isl_basic_map *)bset;
2387         }
2388         bset = isl_basic_set_cow(bset);
2389         if (!bset)
2390                 goto error;
2391         total = bset->dim + bset->extra;
2392         bmap = (struct isl_basic_map *)bset;
2393         bmap->nparam = like->nparam;
2394         bmap->n_in = like->n_in;
2395         bmap->n_out = like->n_out;
2396         bmap->extra += like->n_div;
2397         if (bmap->extra) {
2398                 unsigned ltotal;
2399                 ltotal = total - bmap->extra + like->extra;
2400                 if (ltotal > total)
2401                         ltotal = total;
2402                 bmap->block2 = isl_blk_extend(ctx, bmap->block2,
2403                                         bmap->extra * (1 + 1 + total));
2404                 if (isl_blk_is_error(bmap->block2))
2405                         goto error;
2406                 bmap->div = isl_realloc_array(ctx, bmap->div, isl_int *,
2407                                                 bmap->extra);
2408                 if (!bmap->div)
2409                         goto error;
2410                 bmap = isl_basic_map_extend(bmap, bmap->nparam,
2411                         bmap->n_in, bmap->n_out, 0, 0, 2 * like->n_div);
2412                 for (i = 0; i < like->n_div; ++i) {
2413                         k = isl_basic_map_alloc_div(bmap);
2414                         if (k < 0)
2415                                 goto error;
2416                         isl_seq_cpy(bmap->div[k], like->div[i], 1 + 1 + ltotal);
2417                         isl_seq_clr(bmap->div[k]+1+1+ltotal, total - ltotal);
2418                         if (add_div_constraints(bmap, k) < 0)
2419                                 goto error;
2420                 }
2421         }
2422         isl_basic_map_free(like);
2423         bmap = isl_basic_map_finalize(bmap);
2424         return bmap;
2425 error:
2426         isl_basic_map_free(like);
2427         isl_basic_set_free(bset);
2428         return NULL;
2429 }
2430
2431 struct isl_basic_set *isl_basic_set_from_underlying_set(
2432         struct isl_basic_set *bset, struct isl_basic_set *like)
2433 {
2434         return (struct isl_basic_set *)
2435                 isl_basic_map_overlying_set(bset, (struct isl_basic_map *)like);
2436 }
2437
2438 struct isl_set *isl_set_from_underlying_set(
2439         struct isl_set *set, struct isl_basic_set *like)
2440 {
2441         int i;
2442
2443         if (!set || !like)
2444                 goto error;
2445         isl_assert(set->ctx, set->dim == like->nparam + like->dim + like->n_div,
2446                     goto error);
2447         if (like->nparam == 0 && like->n_div == 0) {
2448                 isl_basic_set_free(like);
2449                 return set;
2450         }
2451         set = isl_set_cow(set);
2452         if (!set)
2453                 goto error;
2454         for (i = 0; i < set->n; ++i) {
2455                 set->p[i] = isl_basic_set_from_underlying_set(set->p[i],
2456                                                       isl_basic_set_copy(like));
2457                 if (!set->p[i])
2458                         goto error;
2459         }
2460         set->nparam = like->nparam;
2461         set->dim = like->dim;
2462         isl_basic_set_free(like);
2463         return set;
2464 error:
2465         isl_basic_set_free(like);
2466         isl_set_free(set);
2467         return NULL;
2468 }
2469
2470 struct isl_set *isl_map_underlying_set(struct isl_map *map)
2471 {
2472         int i;
2473
2474         map = isl_map_align_divs(map);
2475         map = isl_map_cow(map);
2476         if (!map)
2477                 return NULL;
2478
2479         for (i = 0; i < map->n; ++i) {
2480                 map->p[i] = (struct isl_basic_map *)
2481                                 isl_basic_map_underlying_set(map->p[i]);
2482                 if (!map->p[i])
2483                         goto error;
2484         }
2485         if (map->n == 0)
2486                 map->n_out += map->nparam + map->n_in;
2487         else
2488                 map->n_out = map->p[0]->n_out;
2489         map->nparam = 0;
2490         map->n_in = 0;
2491         return (struct isl_set *)map;
2492 error:
2493         isl_map_free(map);
2494         return NULL;
2495 }
2496
2497 struct isl_set *isl_set_to_underlying_set(struct isl_set *set)
2498 {
2499         return (struct isl_set *)isl_map_underlying_set((struct isl_map *)set);
2500 }
2501
2502 struct isl_basic_set *isl_basic_map_domain(struct isl_basic_map *bmap)
2503 {
2504         struct isl_basic_set *domain;
2505         unsigned n_out;
2506         if (!bmap)
2507                 return NULL;
2508         n_out = bmap->n_out;
2509         domain = isl_basic_set_from_basic_map(bmap);
2510         return isl_basic_set_project_out(domain, n_out, 0);
2511 }
2512
2513 struct isl_basic_set *isl_basic_map_range(struct isl_basic_map *bmap)
2514 {
2515         return isl_basic_map_domain(isl_basic_map_reverse(bmap));
2516 }
2517
2518 struct isl_set *isl_map_range(struct isl_map *map)
2519 {
2520         int i;
2521         struct isl_set *set;
2522
2523         if (!map)
2524                 goto error;
2525         map = isl_map_cow(map);
2526         if (!map)
2527                 goto error;
2528
2529         set = (struct isl_set *) map;
2530         set->zero = 0;
2531         for (i = 0; i < map->n; ++i) {
2532                 set->p[i] = isl_basic_map_range(map->p[i]);
2533                 if (!set->p[i])
2534                         goto error;
2535         }
2536         F_CLR(set, ISL_MAP_DISJOINT);
2537         F_CLR(set, ISL_SET_NORMALIZED);
2538         return set;
2539 error:
2540         isl_map_free(map);
2541         return NULL;
2542 }
2543
2544 struct isl_map *isl_map_from_set(struct isl_set *set,
2545                 unsigned n_in, unsigned n_out)
2546 {
2547         int i;
2548         struct isl_map *map = NULL;
2549
2550         if (!set)
2551                 return NULL;
2552         isl_assert(set->ctx, set->dim == n_in + n_out, goto error);
2553         set = isl_set_cow(set);
2554         if (!set)
2555                 return NULL;
2556         map = (struct isl_map *)set;
2557         for (i = 0; i < set->n; ++i) {
2558                 map->p[i] = isl_basic_map_from_basic_set(
2559                                 set->p[i], n_in, n_out);
2560                 if (!map->p[i])
2561                         goto error;
2562         }
2563         map->n_in = n_in;
2564         map->n_out = n_out;
2565         return map;
2566 error:
2567         isl_set_free(set);
2568         return NULL;
2569 }
2570
2571 struct isl_set *isl_set_from_map(struct isl_map *map)
2572 {
2573         int i;
2574         struct isl_set *set = NULL;
2575
2576         if (!map)
2577                 return NULL;
2578         map = isl_map_cow(map);
2579         if (!map)
2580                 return NULL;
2581         map->n_out += map->n_in;
2582         map->n_in = 0;
2583         set = (struct isl_set *)map;
2584         for (i = 0; i < map->n; ++i) {
2585                 set->p[i] = isl_basic_set_from_basic_map(map->p[i]);
2586                 if (!set->p[i])
2587                         goto error;
2588         }
2589         return set;
2590 error:
2591         isl_map_free(map);
2592         return NULL;
2593 }
2594
2595 struct isl_map *isl_map_alloc(struct isl_ctx *ctx,
2596                 unsigned nparam, unsigned in, unsigned out, int n,
2597                 unsigned flags)
2598 {
2599         struct isl_map *map;
2600
2601         map = isl_alloc(ctx, struct isl_map,
2602                         sizeof(struct isl_map) +
2603                         n * sizeof(struct isl_basic_map *));
2604         if (!map)
2605                 return NULL;
2606
2607         map->ctx = ctx;
2608         isl_ctx_ref(ctx);
2609         map->ref = 1;
2610         map->size = n;
2611         map->n = 0;
2612         map->nparam = nparam;
2613         map->n_in = in;
2614         map->n_out = out;
2615         map->flags = flags;
2616         return map;
2617 }
2618
2619 struct isl_basic_map *isl_basic_map_empty(struct isl_ctx *ctx,
2620                 unsigned nparam, unsigned in, unsigned out)
2621 {
2622         struct isl_basic_map *bmap;
2623         bmap = isl_basic_map_alloc(ctx, nparam, in, out, 0, 1, 0);
2624         bmap = isl_basic_map_set_to_empty(bmap);
2625         return bmap;
2626 }
2627
2628 struct isl_basic_set *isl_basic_set_empty(struct isl_ctx *ctx,
2629                 unsigned nparam, unsigned dim)
2630 {
2631         struct isl_basic_set *bset;
2632         bset = isl_basic_set_alloc(ctx, nparam, dim, 0, 1, 0);
2633         bset = isl_basic_set_set_to_empty(bset);
2634         return bset;
2635 }
2636
2637 struct isl_basic_map *isl_basic_map_universe(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, 0, 0);
2642         return bmap;
2643 }
2644
2645 struct isl_basic_set *isl_basic_set_universe(struct isl_ctx *ctx,
2646                 unsigned nparam, unsigned dim)
2647 {
2648         struct isl_basic_set *bset;
2649         bset = isl_basic_set_alloc(ctx, nparam, dim, 0, 0, 0);
2650         return bset;
2651 }
2652
2653 struct isl_map *isl_map_empty(struct isl_ctx *ctx,
2654                 unsigned nparam, unsigned in, unsigned out)
2655 {
2656         return isl_map_alloc(ctx, nparam, in, out, 0, ISL_MAP_DISJOINT);
2657 }
2658
2659 struct isl_set *isl_set_empty(struct isl_ctx *ctx,
2660                 unsigned nparam, unsigned dim)
2661 {
2662         return isl_set_alloc(ctx, nparam, dim, 0, ISL_MAP_DISJOINT);
2663 }
2664
2665 struct isl_map *isl_map_dup(struct isl_map *map)
2666 {
2667         int i;
2668         struct isl_map *dup;
2669
2670         if (!map)
2671                 return NULL;
2672         dup = isl_map_alloc(map->ctx, map->nparam, map->n_in, map->n_out, map->n,
2673                                 map->flags);
2674         for (i = 0; i < map->n; ++i)
2675                 dup = isl_map_add(dup,
2676                                   isl_basic_map_dup(map->p[i]));
2677         return dup;
2678 }
2679
2680 struct isl_map *isl_map_add(struct isl_map *map, struct isl_basic_map *bmap)
2681 {
2682         if (!bmap || !map)
2683                 goto error;
2684         isl_assert(map->ctx, map->nparam == bmap->nparam, goto error);
2685         isl_assert(map->ctx, map->n_in == bmap->n_in, goto error);
2686         isl_assert(map->ctx, map->n_out == bmap->n_out, goto error);
2687         isl_assert(map->ctx, map->n < map->size, goto error);
2688         map->p[map->n] = bmap;
2689         map->n++;
2690         F_CLR(map, ISL_MAP_NORMALIZED);
2691         return map;
2692 error:
2693         if (map)
2694                 isl_map_free(map);
2695         if (bmap)
2696                 isl_basic_map_free(bmap);
2697         return NULL;
2698 }
2699
2700 void isl_map_free(struct isl_map *map)
2701 {
2702         int i;
2703
2704         if (!map)
2705                 return;
2706
2707         if (--map->ref > 0)
2708                 return;
2709
2710         isl_ctx_deref(map->ctx);
2711         for (i = 0; i < map->n; ++i)
2712                 isl_basic_map_free(map->p[i]);
2713         free(map);
2714 }
2715
2716 struct isl_map *isl_map_extend(struct isl_map *base,
2717                 unsigned nparam, unsigned n_in, unsigned n_out)
2718 {
2719         int i;
2720
2721         base = isl_map_cow(base);
2722         if (!base)
2723                 return NULL;
2724
2725         isl_assert(base->ctx, base->nparam <= nparam, goto error);
2726         isl_assert(base->ctx, base->n_in <= n_in, goto error);
2727         isl_assert(base->ctx, base->n_out <= n_out, goto error);
2728         base->nparam = nparam;
2729         base->n_in = n_in;
2730         base->n_out = n_out;
2731         for (i = 0; i < base->n; ++i) {
2732                 base->p[i] = isl_basic_map_extend(base->p[i],
2733                                 nparam, n_in, n_out, 0, 0, 0);
2734                 if (!base->p[i])
2735                         goto error;
2736         }
2737         return base;
2738 error:
2739         isl_map_free(base);
2740         return NULL;
2741 }
2742
2743 struct isl_set *isl_set_extend(struct isl_set *base,
2744                 unsigned nparam, unsigned dim)
2745 {
2746         return (struct isl_set *)isl_map_extend((struct isl_map *)base,
2747                                                         nparam, 0, dim);
2748 }
2749
2750 static struct isl_basic_map *isl_basic_map_fix_var(struct isl_basic_map *bmap,
2751                 unsigned var, int value)
2752 {
2753         int j;
2754
2755         bmap = isl_basic_map_cow(bmap);
2756         bmap = isl_basic_map_extend(bmap,
2757                         bmap->nparam, bmap->n_in, bmap->n_out, 0, 1, 0);
2758         j = isl_basic_map_alloc_equality(bmap);
2759         if (j < 0)
2760                 goto error;
2761         isl_seq_clr(bmap->eq[j],
2762                     1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->extra);
2763         isl_int_set_si(bmap->eq[j][1 + var], -1);
2764         isl_int_set_si(bmap->eq[j][0], value);
2765         bmap = isl_basic_map_simplify(bmap);
2766         return isl_basic_map_finalize(bmap);
2767 error:
2768         isl_basic_map_free(bmap);
2769         return NULL;
2770 }
2771
2772 struct isl_basic_map *isl_basic_map_fix_input_si(struct isl_basic_map *bmap,
2773                 unsigned input, int value)
2774 {
2775         if (!bmap)
2776                 return NULL;
2777         isl_assert(bmap->ctx, input < bmap->n_in, goto error);
2778         return isl_basic_map_fix_var(bmap, bmap->nparam + input, value);
2779 error:
2780         isl_basic_map_free(bmap);
2781         return NULL;
2782 }
2783
2784 struct isl_basic_set *isl_basic_set_fix_dim_si(struct isl_basic_set *bset,
2785                 unsigned dim, int value)
2786 {
2787         if (!bset)
2788                 return NULL;
2789         isl_assert(bset->ctx, dim < bset->dim, goto error);
2790         return (struct isl_basic_set *)
2791                 isl_basic_map_fix_var((struct isl_basic_map *)bset,
2792                                                 bset->nparam + dim, value);
2793 error:
2794         isl_basic_set_free(bset);
2795         return NULL;
2796 }
2797
2798 struct isl_map *isl_map_fix_input_si(struct isl_map *map,
2799                 unsigned input, int value)
2800 {
2801         int i;
2802
2803         map = isl_map_cow(map);
2804         if (!map)
2805                 return NULL;
2806
2807         isl_assert(ctx, input < map->n_in, goto error);
2808         for (i = 0; i < map->n; ++i) {
2809                 map->p[i] = isl_basic_map_fix_input_si(map->p[i], input, value);
2810                 if (!map->p[i])
2811                         goto error;
2812         }
2813         F_CLR(map, ISL_MAP_NORMALIZED);
2814         return map;
2815 error:
2816         isl_map_free(map);
2817         return NULL;
2818 }
2819
2820 struct isl_set *isl_set_fix_dim_si(struct isl_set *set, unsigned dim, int value)
2821 {
2822         int i;
2823
2824         set = isl_set_cow(set);
2825         if (!set)
2826                 return NULL;
2827
2828         isl_assert(set->ctx, dim < set->dim, goto error);
2829         for (i = 0; i < set->n; ++i) {
2830                 set->p[i] = isl_basic_set_fix_dim_si(set->p[i], dim, value);
2831                 if (!set->p[i])
2832                         goto error;
2833         }
2834         return set;
2835 error:
2836         isl_set_free(set);
2837         return NULL;
2838 }
2839
2840 struct isl_basic_set *isl_basic_set_lower_bound_dim(struct isl_basic_set *bset,
2841         unsigned dim, isl_int value)
2842 {
2843         int j;
2844
2845         bset = isl_basic_set_cow(bset);
2846         bset = isl_basic_set_extend(bset, bset->nparam, bset->dim, 0, 0, 1);
2847         j = isl_basic_set_alloc_inequality(bset);
2848         if (j < 0)
2849                 goto error;
2850         isl_seq_clr(bset->ineq[j], 1 + bset->nparam + bset->dim + bset->extra);
2851         isl_int_set_si(bset->ineq[j][1 + bset->nparam + dim], 1);
2852         isl_int_neg(bset->ineq[j][0], value);
2853         bset = isl_basic_set_simplify(bset);
2854         return isl_basic_set_finalize(bset);
2855 error:
2856         isl_basic_set_free(bset);
2857         return NULL;
2858 }
2859
2860 struct isl_set *isl_set_lower_bound_dim(struct isl_set *set, unsigned dim,
2861                                         isl_int value)
2862 {
2863         int i;
2864
2865         set = isl_set_cow(set);
2866         if (!set)
2867                 return NULL;
2868
2869         isl_assert(set->ctx, dim < set->dim, goto error);
2870         for (i = 0; i < set->n; ++i) {
2871                 set->p[i] = isl_basic_set_lower_bound_dim(set->p[i], dim, value);
2872                 if (!set->p[i])
2873                         goto error;
2874         }
2875         return set;
2876 error:
2877         isl_set_free(set);
2878         return NULL;
2879 }
2880
2881 struct isl_map *isl_map_reverse(struct isl_map *map)
2882 {
2883         int i;
2884         unsigned t;
2885
2886         map = isl_map_cow(map);
2887         if (!map)
2888                 return NULL;
2889
2890         t = map->n_in;
2891         map->n_in = map->n_out;
2892         map->n_out = t;
2893         for (i = 0; i < map->n; ++i) {
2894                 map->p[i] = isl_basic_map_reverse(map->p[i]);
2895                 if (!map->p[i])
2896                         goto error;
2897         }
2898         F_CLR(map, ISL_MAP_NORMALIZED);
2899         return map;
2900 error:
2901         isl_map_free(map);
2902         return NULL;
2903 }
2904
2905 struct isl_map *isl_basic_map_lexmax(
2906                 struct isl_basic_map *bmap, struct isl_basic_set *dom,
2907                 struct isl_set **empty)
2908 {
2909         return isl_pip_basic_map_lexmax(bmap, dom, empty);
2910 }
2911
2912 struct isl_map *isl_basic_map_lexmin(
2913                 struct isl_basic_map *bmap, struct isl_basic_set *dom,
2914                 struct isl_set **empty)
2915 {
2916         return isl_pip_basic_map_lexmin(bmap, dom, empty);
2917 }
2918
2919 struct isl_set *isl_basic_set_lexmin(struct isl_basic_set *bset)
2920 {
2921         struct isl_basic_map *bmap = NULL;
2922         struct isl_basic_set *dom = NULL;
2923         struct isl_map *min;
2924
2925         if (!bset)
2926                 goto error;
2927         bmap = isl_basic_map_from_basic_set(bset, 0, bset->dim);
2928         if (!bmap)
2929                 goto error;
2930         dom = isl_basic_set_alloc(bmap->ctx, bmap->nparam, 0, 0, 0, 0);
2931         if (!dom)
2932                 goto error;
2933         min = isl_basic_map_lexmin(bmap, dom, NULL);
2934         return isl_map_range(min);
2935 error:
2936         isl_basic_map_free(bmap);
2937         return NULL;
2938 }
2939
2940 struct isl_map *isl_basic_map_compute_divs(struct isl_basic_map *bmap)
2941 {
2942         if (!bmap)
2943                 return NULL;
2944         if (bmap->n_div == 0)
2945                 return isl_map_from_basic_map(bmap);
2946         return isl_pip_basic_map_compute_divs(bmap);
2947 }
2948
2949 struct isl_map *isl_map_compute_divs(struct isl_map *map)
2950 {
2951         int i;
2952         struct isl_map *res;
2953
2954         if (!map)
2955                 return NULL;
2956         if (map->n == 0)
2957                 return map;
2958         res = isl_basic_map_compute_divs(isl_basic_map_copy(map->p[0]));
2959         for (i = 1 ; i < map->n; ++i) {
2960                 struct isl_map *r2;
2961                 r2 = isl_basic_map_compute_divs(isl_basic_map_copy(map->p[i]));
2962                 if (F_ISSET(map, ISL_MAP_DISJOINT))
2963                         res = isl_map_union_disjoint(res, r2);
2964                 else
2965                         res = isl_map_union(res, r2);
2966         }
2967         isl_map_free(map);
2968
2969         return res;
2970 }
2971
2972 struct isl_set *isl_basic_set_compute_divs(struct isl_basic_set *bset)
2973 {
2974         return (struct isl_set *)
2975                 isl_basic_map_compute_divs((struct isl_basic_map *)bset);
2976 }
2977
2978 struct isl_set *isl_set_compute_divs(struct isl_set *set)
2979 {
2980         return (struct isl_set *)
2981                 isl_map_compute_divs((struct isl_map *)set);
2982 }
2983
2984 struct isl_set *isl_map_domain(struct isl_map *map)
2985 {
2986         int i;
2987         struct isl_set *set;
2988
2989         if (!map)
2990                 goto error;
2991
2992         map = isl_map_cow(map);
2993         if (!map)
2994                 return NULL;
2995
2996         set = (struct isl_set *)map;
2997         set->dim = map->n_in;
2998         set->zero = 0;
2999         for (i = 0; i < map->n; ++i) {
3000                 set->p[i] = isl_basic_map_domain(map->p[i]);
3001                 if (!set->p[i])
3002                         goto error;
3003         }
3004         F_CLR(set, ISL_MAP_DISJOINT);
3005         F_CLR(set, ISL_SET_NORMALIZED);
3006         return set;
3007 error:
3008         isl_map_free(map);
3009         return NULL;
3010 }
3011
3012 struct isl_map *isl_map_union_disjoint(
3013                         struct isl_map *map1, struct isl_map *map2)
3014 {
3015         int i;
3016         unsigned flags = 0;
3017         struct isl_map *map = NULL;
3018
3019         if (!map1 || !map2)
3020                 goto error;
3021
3022         if (map1->n == 0) {
3023                 isl_map_free(map1);
3024                 return map2;
3025         }
3026         if (map2->n == 0) {
3027                 isl_map_free(map2);
3028                 return map1;
3029         }
3030
3031         isl_assert(map1->ctx, map1->nparam == map2->nparam, goto error);
3032         isl_assert(map1->ctx, map1->n_in == map2->n_in, goto error);
3033         isl_assert(map1->ctx, map1->n_out == map2->n_out, goto error);
3034
3035         if (F_ISSET(map1, ISL_MAP_DISJOINT) &&
3036             F_ISSET(map2, ISL_MAP_DISJOINT))
3037                 FL_SET(flags, ISL_MAP_DISJOINT);
3038
3039         map = isl_map_alloc(map1->ctx, map1->nparam, map1->n_in, map1->n_out,
3040                                 map1->n + map2->n, flags);
3041         if (!map)
3042                 goto error;
3043         for (i = 0; i < map1->n; ++i) {
3044                 map = isl_map_add(map,
3045                                   isl_basic_map_copy(map1->p[i]));
3046                 if (!map)
3047                         goto error;
3048         }
3049         for (i = 0; i < map2->n; ++i) {
3050                 map = isl_map_add(map,
3051                                   isl_basic_map_copy(map2->p[i]));
3052                 if (!map)
3053                         goto error;
3054         }
3055         isl_map_free(map1);
3056         isl_map_free(map2);
3057         return map;
3058 error:
3059         isl_map_free(map);
3060         isl_map_free(map1);
3061         isl_map_free(map2);
3062         return NULL;
3063 }
3064
3065 struct isl_map *isl_map_union(struct isl_map *map1, struct isl_map *map2)
3066 {
3067         map1 = isl_map_union_disjoint(map1, map2);
3068         if (!map1)
3069                 return NULL;
3070         if (map1->n > 1)
3071                 F_CLR(map1, ISL_MAP_DISJOINT);
3072         return map1;
3073 }
3074
3075 struct isl_set *isl_set_union_disjoint(
3076                         struct isl_set *set1, struct isl_set *set2)
3077 {
3078         return (struct isl_set *)
3079                 isl_map_union_disjoint(
3080                         (struct isl_map *)set1, (struct isl_map *)set2);
3081 }
3082
3083 struct isl_set *isl_set_union(struct isl_set *set1, struct isl_set *set2)
3084 {
3085         return (struct isl_set *)
3086                 isl_map_union((struct isl_map *)set1, (struct isl_map *)set2);
3087 }
3088
3089 struct isl_map *isl_map_intersect_range(
3090                 struct isl_map *map, struct isl_set *set)
3091 {
3092         unsigned flags = 0;
3093         struct isl_map *result;
3094         int i, j;
3095
3096         if (!map || !set)
3097                 goto error;
3098
3099         if (F_ISSET(map, ISL_MAP_DISJOINT) &&
3100             F_ISSET(set, ISL_MAP_DISJOINT))
3101                 FL_SET(flags, ISL_MAP_DISJOINT);
3102
3103         result = isl_map_alloc(map->ctx, map->nparam, map->n_in, map->n_out,
3104                                 map->n * set->n, flags);
3105         if (!result)
3106                 goto error;
3107         for (i = 0; i < map->n; ++i)
3108                 for (j = 0; j < set->n; ++j) {
3109                         result = isl_map_add(result,
3110                             isl_basic_map_intersect_range(
3111                                 isl_basic_map_copy(map->p[i]),
3112                                 isl_basic_set_copy(set->p[j])));
3113                         if (!result)
3114                                 goto error;
3115                 }
3116         isl_map_free(map);
3117         isl_set_free(set);
3118         return result;
3119 error:
3120         isl_map_free(map);
3121         isl_set_free(set);
3122         return NULL;
3123 }
3124
3125 struct isl_map *isl_map_intersect_domain(
3126                 struct isl_map *map, struct isl_set *set)
3127 {
3128         return isl_map_reverse(
3129                 isl_map_intersect_range(isl_map_reverse(map), set));
3130 }
3131
3132 struct isl_map *isl_map_apply_domain(
3133                 struct isl_map *map1, struct isl_map *map2)
3134 {
3135         if (!map1 || !map2)
3136                 goto error;
3137         map1 = isl_map_reverse(map1);
3138         map1 = isl_map_apply_range(map1, map2);
3139         return isl_map_reverse(map1);
3140 error:
3141         isl_map_free(map1);
3142         isl_map_free(map2);
3143         return NULL;
3144 }
3145
3146 struct isl_map *isl_map_apply_range(
3147                 struct isl_map *map1, struct isl_map *map2)
3148 {
3149         struct isl_map *result;
3150         int i, j;
3151
3152         if (!map1 || !map2)
3153                 goto error;
3154
3155         isl_assert(map1->ctx, map1->nparam == map2->nparam, goto error);
3156         isl_assert(map1->ctx, map1->n_out == map2->n_in, goto error);
3157
3158         result = isl_map_alloc(map1->ctx, map1->nparam, map1->n_in, map2->n_out,
3159                                 map1->n * map2->n, 0);
3160         if (!result)
3161                 goto error;
3162         for (i = 0; i < map1->n; ++i)
3163                 for (j = 0; j < map2->n; ++j) {
3164                         result = isl_map_add(result,
3165                             isl_basic_map_apply_range(
3166                                 isl_basic_map_copy(map1->p[i]),
3167                                 isl_basic_map_copy(map2->p[j])));
3168                         if (!result)
3169                                 goto error;
3170                 }
3171         isl_map_free(map1);
3172         isl_map_free(map2);
3173         if (result && result->n <= 1)
3174                 F_SET(result, ISL_MAP_DISJOINT);
3175         return result;
3176 error:
3177         isl_map_free(map1);
3178         isl_map_free(map2);
3179         return NULL;
3180 }
3181
3182 /*
3183  * returns range - domain
3184  */
3185 struct isl_basic_set *isl_basic_map_deltas(struct isl_basic_map *bmap)
3186 {
3187         struct isl_basic_set *bset;
3188         unsigned dim;
3189         int i;
3190
3191         if (!bmap)
3192                 goto error;
3193         isl_assert(bmap->ctx, bmap->n_in == bmap->n_out, goto error);
3194         dim = bmap->n_in;
3195         bset = isl_basic_set_from_basic_map(bmap);
3196         bset = isl_basic_set_extend(bset, bset->nparam, 3*dim, 0,
3197                                         dim, 0);
3198         bset = isl_basic_set_swap_vars(bset, 2*dim);
3199         for (i = 0; i < dim; ++i) {
3200                 int j = isl_basic_map_alloc_equality(
3201                                             (struct isl_basic_map *)bset);
3202                 if (j < 0)
3203                         goto error;
3204                 isl_seq_clr(bset->eq[j],
3205                             1 + bset->nparam + bset->dim + bset->extra);
3206                 isl_int_set_si(bset->eq[j][1+bset->nparam+i], 1);
3207                 isl_int_set_si(bset->eq[j][1+bset->nparam+dim+i], 1);
3208                 isl_int_set_si(bset->eq[j][1+bset->nparam+2*dim+i], -1);
3209         }
3210         return isl_basic_set_project_out(bset, 2*dim, 0);
3211 error:
3212         isl_basic_map_free(bmap);
3213         return NULL;
3214 }
3215
3216 /*
3217  * returns range - domain
3218  */
3219 struct isl_set *isl_map_deltas(struct isl_map *map)
3220 {
3221         int i;
3222         struct isl_set *result;
3223
3224         if (!map)
3225                 return NULL;
3226
3227         isl_assert(map->ctx, map->n_in == map->n_out, goto error);
3228         result = isl_set_alloc(map->ctx, map->nparam, map->n_in, map->n, map->flags);
3229         if (!result)
3230                 goto error;
3231         for (i = 0; i < map->n; ++i)
3232                 result = isl_set_add(result,
3233                           isl_basic_map_deltas(isl_basic_map_copy(map->p[i])));
3234         isl_map_free(map);
3235         return result;
3236 error:
3237         isl_map_free(map);
3238         return NULL;
3239 }
3240
3241 /* If the only constraints a div d=floor(f/m)
3242  * appears in are its two defining constraints
3243  *
3244  *      f - m d >=0
3245  *      -(f - (m - 1)) + m d >= 0
3246  *
3247  * then it can safely be removed.
3248  */
3249 static int div_is_redundant(struct isl_basic_map *bmap, int div)
3250 {
3251         int i;
3252         unsigned pos = 1 + bmap->nparam + bmap->n_in + bmap->n_out + div;
3253
3254         for (i = 0; i < bmap->n_eq; ++i)
3255                 if (!isl_int_is_zero(bmap->eq[i][pos]))
3256                         return 0;
3257
3258         for (i = 0; i < bmap->n_ineq; ++i) {
3259                 if (isl_int_is_zero(bmap->ineq[i][pos]))
3260                         continue;
3261                 if (isl_int_eq(bmap->ineq[i][pos], bmap->div[div][0])) {
3262                         int neg;
3263                         isl_int_sub(bmap->div[div][1],
3264                                         bmap->div[div][1], bmap->div[div][0]);
3265                         isl_int_add_ui(bmap->div[div][1], bmap->div[div][1], 1);
3266                         neg = isl_seq_is_neg(bmap->ineq[i], bmap->div[div]+1, pos);
3267                         isl_int_sub_ui(bmap->div[div][1], bmap->div[div][1], 1);
3268                         isl_int_add(bmap->div[div][1],
3269                                         bmap->div[div][1], bmap->div[div][0]);
3270                         if (!neg)
3271                                 return 0;
3272                         if (isl_seq_first_non_zero(bmap->ineq[i]+pos+1,
3273                                                     bmap->n_div-div-1) != -1)
3274                                 return 0;
3275                 } else if (isl_int_abs_eq(bmap->ineq[i][pos], bmap->div[div][0])) {
3276                         if (!isl_seq_eq(bmap->ineq[i], bmap->div[div]+1, pos))
3277                                 return 0;
3278                         if (isl_seq_first_non_zero(bmap->ineq[i]+pos+1,
3279                                                     bmap->n_div-div-1) != -1)
3280                                 return 0;
3281                 } else
3282                         return 0;
3283         }
3284
3285         for (i = 0; i < bmap->n_div; ++i)
3286                 if (!isl_int_is_zero(bmap->div[i][1+pos]))
3287                         return 0;
3288
3289         return 1;
3290 }
3291
3292 /*
3293  * Remove divs that don't occur in any of the constraints or other divs.
3294  * These can arise when dropping some of the variables in a quast
3295  * returned by piplib.
3296  */
3297 static struct isl_basic_map *remove_redundant_divs(struct isl_basic_map *bmap)
3298 {
3299         int i;
3300
3301         if (!bmap)
3302                 return NULL;
3303
3304         for (i = bmap->n_div-1; i >= 0; --i) {
3305                 if (!div_is_redundant(bmap, i))
3306                         continue;
3307                 bmap = isl_basic_map_drop_div(bmap, i);
3308         }
3309         return bmap;
3310 }
3311
3312 struct isl_basic_map *isl_basic_map_finalize(struct isl_basic_map *bmap)
3313 {
3314         bmap = remove_redundant_divs(bmap);
3315         if (!bmap)
3316                 return NULL;
3317         F_SET(bmap, ISL_BASIC_SET_FINAL);
3318         return bmap;
3319 }
3320
3321 struct isl_basic_set *isl_basic_set_finalize(struct isl_basic_set *bset)
3322 {
3323         return (struct isl_basic_set *)
3324                 isl_basic_map_finalize((struct isl_basic_map *)bset);
3325 }
3326
3327 struct isl_set *isl_set_finalize(struct isl_set *set)
3328 {
3329         int i;
3330
3331         if (!set)
3332                 return NULL;
3333         for (i = 0; i < set->n; ++i) {
3334                 set->p[i] = isl_basic_set_finalize(set->p[i]);
3335                 if (!set->p[i])
3336                         goto error;
3337         }
3338         return set;
3339 error:
3340         isl_set_free(set);
3341         return NULL;
3342 }
3343
3344 struct isl_map *isl_map_finalize(struct isl_map *map)
3345 {
3346         int i;
3347
3348         if (!map)
3349                 return NULL;
3350         for (i = 0; i < map->n; ++i) {
3351                 map->p[i] = isl_basic_map_finalize(map->p[i]);
3352                 if (!map->p[i])
3353                         goto error;
3354         }
3355         F_CLR(map, ISL_MAP_NORMALIZED);
3356         return map;
3357 error:
3358         isl_map_free(map);
3359         return NULL;
3360 }
3361
3362 struct isl_basic_map *isl_basic_map_identity(struct isl_ctx *ctx,
3363                 unsigned nparam, unsigned dim)
3364 {
3365         struct isl_basic_map *bmap;
3366         int i;
3367
3368         bmap = isl_basic_map_alloc(ctx, nparam, dim, dim, 0, dim, 0);
3369         if (!bmap)
3370                 goto error;
3371
3372         for (i = 0; i < dim; ++i) {
3373                 int j = isl_basic_map_alloc_equality(bmap);
3374                 if (j < 0)
3375                         goto error;
3376                 isl_seq_clr(bmap->eq[j],
3377                     1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->extra);
3378                 isl_int_set_si(bmap->eq[j][1+nparam+i], 1);
3379                 isl_int_set_si(bmap->eq[j][1+nparam+dim+i], -1);
3380         }
3381         return isl_basic_map_finalize(bmap);
3382 error:
3383         isl_basic_map_free(bmap);
3384         return NULL;
3385 }
3386
3387 struct isl_map *isl_map_identity(struct isl_ctx *ctx,
3388                 unsigned nparam, unsigned dim)
3389 {
3390         struct isl_map *map = isl_map_alloc(ctx, nparam, dim, dim, 1,
3391                                                 ISL_MAP_DISJOINT);
3392         if (!map)
3393                 goto error;
3394         map = isl_map_add(map,
3395                 isl_basic_map_identity(ctx, nparam, dim));
3396         return map;
3397 error:
3398         isl_map_free(map);
3399         return NULL;
3400 }
3401
3402 int isl_set_is_equal(struct isl_set *set1, struct isl_set *set2)
3403 {
3404         return isl_map_is_equal((struct isl_map *)set1, (struct isl_map *)set2);
3405 }
3406
3407 int isl_set_is_subset(struct isl_set *set1, struct isl_set *set2)
3408 {
3409         return isl_map_is_subset(
3410                         (struct isl_map *)set1, (struct isl_map *)set2);
3411 }
3412
3413 int isl_basic_map_is_subset(
3414                 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
3415 {
3416         int is_subset;
3417         struct isl_map *map1;
3418         struct isl_map *map2;
3419
3420         if (!bmap1 || !bmap2)
3421                 return -1;
3422
3423         map1 = isl_map_from_basic_map(isl_basic_map_copy(bmap1));
3424         map2 = isl_map_from_basic_map(isl_basic_map_copy(bmap2));
3425
3426         is_subset = isl_map_is_subset(map1, map2);
3427
3428         isl_map_free(map1);
3429         isl_map_free(map2);
3430
3431         return is_subset;
3432 }
3433
3434 int isl_basic_map_is_equal(
3435                 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
3436 {
3437         int is_subset;
3438
3439         if (!bmap1 || !bmap2)
3440                 return -1;
3441         is_subset = isl_basic_map_is_subset(bmap1, bmap2);
3442         if (is_subset != 1)
3443                 return is_subset;
3444         is_subset = isl_basic_map_is_subset(bmap2, bmap1);
3445         return is_subset;
3446 }
3447
3448 int isl_basic_set_is_equal(
3449                 struct isl_basic_set *bset1, struct isl_basic_set *bset2)
3450 {
3451         return isl_basic_map_is_equal(
3452                 (struct isl_basic_map *)bset1, (struct isl_basic_map *)bset2);
3453 }
3454
3455 int isl_map_is_empty(struct isl_map *map)
3456 {
3457         int i;
3458         int is_empty;
3459
3460         if (!map)
3461                 return -1;
3462         for (i = 0; i < map->n; ++i) {
3463                 is_empty = isl_basic_map_is_empty(map->p[i]);
3464                 if (is_empty < 0)
3465                         return -1;
3466                 if (!is_empty)
3467                         return 0;
3468         }
3469         return 1;
3470 }
3471
3472 int isl_set_is_empty(struct isl_set *set)
3473 {
3474         return isl_map_is_empty((struct isl_map *)set);
3475 }
3476
3477 int isl_map_is_subset(struct isl_map *map1, struct isl_map *map2)
3478 {
3479         int i;
3480         int is_subset = 0;
3481         struct isl_map *diff;
3482
3483         if (!map1 || !map2)
3484                 return -1;
3485
3486         if (isl_map_is_empty(map1))
3487                 return 1;
3488
3489         if (isl_map_is_empty(map2))
3490                 return 0;
3491
3492         diff = isl_map_subtract(isl_map_copy(map1), isl_map_copy(map2));
3493         if (!diff)
3494                 return -1;
3495
3496         is_subset = isl_map_is_empty(diff);
3497         isl_map_free(diff);
3498
3499         return is_subset;
3500 }
3501
3502 int isl_map_is_equal(struct isl_map *map1, struct isl_map *map2)
3503 {
3504         int is_subset;
3505
3506         if (!map1 || !map2)
3507                 return -1;
3508         is_subset = isl_map_is_subset(map1, map2);
3509         if (is_subset != 1)
3510                 return is_subset;
3511         is_subset = isl_map_is_subset(map2, map1);
3512         return is_subset;
3513 }
3514
3515 int isl_basic_map_is_strict_subset(
3516                 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
3517 {
3518         int is_subset;
3519
3520         if (!bmap1 || !bmap2)
3521                 return -1;
3522         is_subset = isl_basic_map_is_subset(bmap1, bmap2);
3523         if (is_subset != 1)
3524                 return is_subset;
3525         is_subset = isl_basic_map_is_subset(bmap2, bmap1);
3526         if (is_subset == -1)
3527                 return is_subset;
3528         return !is_subset;
3529 }
3530
3531 static int basic_map_contains(struct isl_basic_map *bmap, struct isl_vec *vec)
3532 {
3533         int i;
3534         unsigned total;
3535         isl_int s;
3536
3537         total = 1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
3538         if (total != vec->size)
3539                 return -1;
3540
3541         isl_int_init(s);
3542
3543         for (i = 0; i < bmap->n_eq; ++i) {
3544                 isl_seq_inner_product(vec->block.data, bmap->eq[i], total, &s);
3545                 if (!isl_int_is_zero(s)) {
3546                         isl_int_clear(s);
3547                         return 0;
3548                 }
3549         }
3550
3551         for (i = 0; i < bmap->n_ineq; ++i) {
3552                 isl_seq_inner_product(vec->block.data, bmap->ineq[i], total, &s);
3553                 if (isl_int_is_neg(s)) {
3554                         isl_int_clear(s);
3555                         return 0;
3556                 }
3557         }
3558
3559         isl_int_clear(s);
3560
3561         return 1;
3562 }
3563
3564 int isl_basic_map_is_universe(struct isl_basic_map *bmap)
3565 {
3566         if (!bmap)
3567                 return -1;
3568         return bmap->n_eq == 0 && bmap->n_ineq == 0;
3569 }
3570
3571 int isl_basic_map_is_empty(struct isl_basic_map *bmap)
3572 {
3573         struct isl_basic_set *bset = NULL;
3574         struct isl_vec *sample = NULL;
3575         int empty;
3576         unsigned total;
3577
3578         if (!bmap)
3579                 return -1;
3580
3581         if (F_ISSET(bmap, ISL_BASIC_MAP_EMPTY))
3582                 return 1;
3583
3584         if (F_ISSET(bmap, ISL_BASIC_MAP_RATIONAL)) {
3585                 struct isl_basic_map *copy = isl_basic_map_copy(bmap);
3586                 copy = isl_basic_map_convex_hull(copy);
3587                 empty = F_ISSET(copy, ISL_BASIC_MAP_EMPTY);
3588                 isl_basic_map_free(copy);
3589                 return empty;
3590         }
3591
3592         total = 1 + bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
3593         if (bmap->sample && bmap->sample->size == total) {
3594                 int contains = basic_map_contains(bmap, bmap->sample);
3595                 if (contains < 0)
3596                         return -1;
3597                 if (contains)
3598                         return 0;
3599         }
3600         bset = isl_basic_map_underlying_set(isl_basic_map_copy(bmap));
3601         if (!bset)
3602                 return -1;
3603         sample = isl_basic_set_sample(bset);
3604         if (!sample)
3605                 return -1;
3606         empty = sample->size == 0;
3607         if (bmap->sample)
3608                 isl_vec_free(bmap->ctx, bmap->sample);
3609         bmap->sample = sample;
3610
3611         return empty;
3612 }
3613
3614 struct isl_map *isl_basic_map_union(
3615         struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
3616 {
3617         struct isl_map *map;
3618         if (!bmap1 || !bmap2)
3619                 return NULL;
3620
3621         isl_assert(bmap1->ctx, bmap1->nparam == bmap2->nparam, goto error);
3622         isl_assert(bmap1->ctx, bmap1->n_in == bmap2->n_in, goto error);
3623         isl_assert(bmap1->ctx, bmap1->n_out == bmap2->n_out, goto error);
3624
3625         map = isl_map_alloc(bmap1->ctx, bmap1->nparam,
3626                                 bmap1->n_in, bmap1->n_out, 2, 0);
3627         if (!map)
3628                 goto error;
3629         map = isl_map_add(map, bmap1);
3630         map = isl_map_add(map, bmap2);
3631         return map;
3632 error:
3633         isl_basic_map_free(bmap1);
3634         isl_basic_map_free(bmap2);
3635         return NULL;
3636 }
3637
3638 struct isl_set *isl_basic_set_union(
3639                 struct isl_basic_set *bset1, struct isl_basic_set *bset2)
3640 {
3641         return (struct isl_set *)isl_basic_map_union(
3642                                             (struct isl_basic_map *)bset1,
3643                                             (struct isl_basic_map *)bset2);
3644 }
3645
3646 /* Order divs such that any div only depends on previous divs */
3647 static struct isl_basic_map *order_divs(struct isl_basic_map *bmap)
3648 {
3649         int i;
3650         unsigned off = bmap->nparam + bmap->n_in + bmap->n_out;
3651
3652         for (i = 0; i < bmap->n_div; ++i) {
3653                 int pos;
3654                 pos = isl_seq_first_non_zero(bmap->div[i]+1+1+off+i,
3655                                                             bmap->n_div-i);
3656                 if (pos == -1)
3657                         continue;
3658                 swap_div(bmap, i, pos);
3659                 --i;
3660         }
3661         return bmap;
3662 }
3663
3664 static int find_div(struct isl_basic_map *dst,
3665                         struct isl_basic_map *src, unsigned div)
3666 {
3667         int i;
3668
3669         unsigned total = src->nparam + src->n_in + src->n_out + src->n_div;
3670
3671         for (i = 0; i < dst->n_div; ++i)
3672                 if (isl_seq_eq(dst->div[i], src->div[div], 1+1+total) &&
3673                     isl_seq_first_non_zero(dst->div[i]+1+1+total,
3674                                                 dst->n_div - src->n_div) == -1)
3675                         return i;
3676         return -1;
3677 }
3678
3679 struct isl_basic_map *isl_basic_map_align_divs(
3680                 struct isl_basic_map *dst, struct isl_basic_map *src)
3681 {
3682         int i;
3683         unsigned total = src->nparam + src->n_in + src->n_out + src->n_div;
3684
3685         if (src->n_div == 0)
3686                 return dst;
3687
3688         src = order_divs(src);
3689         dst = isl_basic_map_extend(dst, dst->nparam, dst->n_in,
3690                         dst->n_out, src->n_div, 0, 2 * src->n_div);
3691         if (!dst)
3692                 return NULL;
3693         for (i = 0; i < src->n_div; ++i) {
3694                 int j = find_div(dst, src, i);
3695                 if (j < 0) {
3696                         j = isl_basic_map_alloc_div(dst);
3697                         if (j < 0)
3698                                 goto error;
3699                         isl_seq_cpy(dst->div[j], src->div[i], 1+1+total);
3700                         isl_seq_clr(dst->div[j]+1+1+total,
3701                                                     dst->extra - src->n_div);
3702                         if (add_div_constraints(dst, j) < 0)
3703                                 goto error;
3704                 }
3705                 if (j != i)
3706                         swap_div(dst, i, j);
3707         }
3708         return dst;
3709 error:
3710         isl_basic_map_free(dst);
3711         return NULL;
3712 }
3713
3714 struct isl_map *isl_map_align_divs(struct isl_map *map)
3715 {
3716         int i;
3717
3718         map = isl_map_compute_divs(map);
3719         map = isl_map_cow(map);
3720         if (!map)
3721                 return NULL;
3722
3723         for (i = 1; i < map->n; ++i)
3724                 map->p[0] = isl_basic_map_align_divs(map->p[0], map->p[i]);
3725         for (i = 1; i < map->n; ++i)
3726                 map->p[i] = isl_basic_map_align_divs(map->p[i], map->p[0]);
3727
3728         F_CLR(map, ISL_MAP_NORMALIZED);
3729         return map;
3730 }
3731
3732 static struct isl_map *add_cut_constraint(struct isl_map *dst,
3733                 struct isl_basic_map *src, isl_int *c,
3734                 unsigned len, int oppose)
3735 {
3736         struct isl_basic_map *copy = NULL;
3737         int is_empty;
3738         int k;
3739         unsigned total;
3740
3741         copy = isl_basic_map_copy(src);
3742         copy = isl_basic_map_cow(copy);
3743         if (!copy)
3744                 goto error;
3745         copy = isl_basic_map_extend(copy,
3746                 copy->nparam, copy->n_in, copy->n_out, 0, 0, 1);
3747         k = isl_basic_map_alloc_inequality(copy);
3748         if (k < 0)
3749                 goto error;
3750         if (oppose)
3751                 isl_seq_neg(copy->ineq[k], c, len);
3752         else
3753                 isl_seq_cpy(copy->ineq[k], c, len);
3754         total = 1 + copy->nparam + copy->n_in + copy->n_out + copy->extra;
3755         isl_seq_clr(copy->ineq[k]+len, total - len);
3756         isl_inequality_negate(copy, k);
3757         copy = isl_basic_map_simplify(copy);
3758         copy = isl_basic_map_finalize(copy);
3759         is_empty = isl_basic_map_is_empty(copy);
3760         if (is_empty < 0)
3761                 goto error;
3762         if (!is_empty)
3763                 dst = isl_map_add(dst, copy);
3764         else
3765                 isl_basic_map_free(copy);
3766         return dst;
3767 error:
3768         isl_basic_map_free(copy);
3769         isl_map_free(dst);
3770         return NULL;
3771 }
3772
3773 static struct isl_map *subtract(struct isl_map *map, struct isl_basic_map *bmap)
3774 {
3775         int i, j, k;
3776         unsigned flags = 0;
3777         struct isl_map *rest = NULL;
3778         unsigned max;
3779         unsigned total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
3780
3781         assert(bmap);
3782
3783         if (!map)
3784                 goto error;
3785
3786         if (F_ISSET(map, ISL_MAP_DISJOINT))
3787                 FL_SET(flags, ISL_MAP_DISJOINT);
3788
3789         max = map->n * (2 * bmap->n_eq + bmap->n_ineq);
3790         rest = isl_map_alloc(map->ctx, map->nparam, map->n_in, map->n_out,
3791                                 max, flags);
3792         if (!rest)
3793                 goto error;
3794
3795         for (i = 0; i < map->n; ++i) {
3796                 map->p[i] = isl_basic_map_align_divs(map->p[i], bmap);
3797                 if (!map->p[i])
3798                         goto error;
3799         }
3800
3801         for (j = 0; j < map->n; ++j)
3802                 map->p[j] = isl_basic_map_cow(map->p[j]);
3803
3804         for (i = 0; i < bmap->n_eq; ++i) {
3805                 for (j = 0; j < map->n; ++j) {
3806                         rest = add_cut_constraint(rest,
3807                                 map->p[j], bmap->eq[i], 1+total, 0);
3808                         if (!rest)
3809                                 goto error;
3810
3811                         rest = add_cut_constraint(rest,
3812                                 map->p[j], bmap->eq[i], 1+total, 1);
3813                         if (!rest)
3814                                 goto error;
3815
3816                         map->p[j] = isl_basic_map_extend(map->p[j],
3817                                 map->p[j]->nparam, map->p[j]->n_in,
3818                                 map->p[j]->n_out, 0, 1, 0);
3819                         if (!map->p[j])
3820                                 goto error;
3821                         k = isl_basic_map_alloc_equality(map->p[j]);
3822                         if (k < 0)
3823                                 goto error;
3824                         isl_seq_cpy(map->p[j]->eq[k], bmap->eq[i], 1+total);
3825                         isl_seq_clr(map->p[j]->eq[k]+1+total,
3826                                         map->p[j]->extra - bmap->n_div);
3827                 }
3828         }
3829
3830         for (i = 0; i < bmap->n_ineq; ++i) {
3831                 for (j = 0; j < map->n; ++j) {
3832                         rest = add_cut_constraint(rest,
3833                                 map->p[j], bmap->ineq[i], 1+total, 0);
3834                         if (!rest)
3835                                 goto error;
3836
3837                         map->p[j] = isl_basic_map_extend(map->p[j],
3838                                 map->p[j]->nparam, map->p[j]->n_in,
3839                                 map->p[j]->n_out, 0, 0, 1);
3840                         if (!map->p[j])
3841                                 goto error;
3842                         k = isl_basic_map_alloc_inequality(map->p[j]);
3843                         if (k < 0)
3844                                 goto error;
3845                         isl_seq_cpy(map->p[j]->ineq[k], bmap->ineq[i], 1+total);
3846                         isl_seq_clr(map->p[j]->ineq[k]+1+total,
3847                                         map->p[j]->extra - bmap->n_div);
3848                 }
3849         }
3850
3851         isl_map_free(map);
3852         return rest;
3853 error:
3854         isl_map_free(map);
3855         isl_map_free(rest);
3856         return NULL;
3857 }
3858
3859 struct isl_map *isl_map_subtract(struct isl_map *map1, struct isl_map *map2)
3860 {
3861         int i;
3862         if (!map1 || !map2)
3863                 goto error;
3864
3865         isl_assert(map1->ctx, map1->nparam == map2->nparam, goto error);
3866         isl_assert(map1->ctx, map1->n_in == map2->n_in, goto error);
3867         isl_assert(map1->ctx, map1->n_out == map2->n_out, goto error);
3868
3869         if (isl_map_is_empty(map2)) {
3870                 isl_map_free(map2);
3871                 return map1;
3872         }
3873
3874         map1 = isl_map_compute_divs(map1);
3875         map2 = isl_map_compute_divs(map2);
3876         if (!map1 || !map2)
3877                 goto error;
3878
3879         for (i = 0; map1 && i < map2->n; ++i)
3880                 map1 = subtract(map1, map2->p[i]);
3881
3882         isl_map_free(map2);
3883         return map1;
3884 error:
3885         isl_map_free(map1);
3886         isl_map_free(map2);
3887         return NULL;
3888 }
3889
3890 struct isl_set *isl_set_subtract(struct isl_set *set1, struct isl_set *set2)
3891 {
3892         return (struct isl_set *)
3893                 isl_map_subtract(
3894                         (struct isl_map *)set1, (struct isl_map *)set2);
3895 }
3896
3897 struct isl_set *isl_set_apply(struct isl_set *set, struct isl_map *map)
3898 {
3899         if (!set || !map)
3900                 goto error;
3901         isl_assert(set->ctx, set->dim == map->n_in, goto error);
3902         map = isl_map_intersect_domain(map, set);
3903         set = isl_map_range(map);
3904         return set;
3905 error:
3906         isl_set_free(set);
3907         isl_map_free(map);
3908         return NULL;
3909 }
3910
3911 /* There is no need to cow as removing empty parts doesn't change
3912  * the meaning of the set.
3913  */
3914 struct isl_map *isl_map_remove_empty_parts(struct isl_map *map)
3915 {
3916         int i;
3917
3918         if (!map)
3919                 return NULL;
3920
3921         for (i = map->n-1; i >= 0; --i) {
3922                 if (!F_ISSET(map->p[i], ISL_BASIC_MAP_EMPTY))
3923                         continue;
3924                 isl_basic_map_free(map->p[i]);
3925                 if (i != map->n-1) {
3926                         F_CLR(map, ISL_MAP_NORMALIZED);
3927                         map->p[i] = map->p[map->n-1];
3928                 }
3929                 map->n--;
3930         }
3931
3932         return map;
3933 }
3934
3935 struct isl_set *isl_set_remove_empty_parts(struct isl_set *set)
3936 {
3937         return (struct isl_set *)
3938                 isl_map_remove_empty_parts((struct isl_map *)set);
3939 }
3940
3941 struct isl_basic_set *isl_set_copy_basic_set(struct isl_set *set)
3942 {
3943         struct isl_basic_set *bset;
3944         if (!set || set->n == 0)
3945                 return NULL;
3946         bset = set->p[set->n-1];
3947         isl_assert(set->ctx, F_ISSET(bset, ISL_BASIC_SET_FINAL), return NULL);
3948         return isl_basic_set_copy(bset);
3949 }
3950
3951 struct isl_set *isl_set_drop_basic_set(struct isl_set *set,
3952                                                 struct isl_basic_set *bset)
3953 {
3954         int i;
3955
3956         if (!set || !bset)
3957                 goto error;
3958         for (i = set->n-1; i >= 0; --i) {
3959                 if (set->p[i] != bset)
3960                         continue;
3961                 set = isl_set_cow(set);
3962                 if (!set)
3963                         goto error;
3964                 isl_basic_set_free(set->p[i]);
3965                 if (i != set->n-1) {
3966                         F_CLR(set, ISL_SET_NORMALIZED);
3967                         set->p[i] = set->p[set->n-1];
3968                 }
3969                 set->n--;
3970                 return set;
3971         }
3972         isl_basic_set_free(bset);
3973         return set;
3974 error:
3975         isl_set_free(set);
3976         isl_basic_set_free(bset);
3977         return NULL;
3978 }
3979
3980 /* Given two _disjoint_ basic sets bset1 and bset2, check whether
3981  * for any common value of the parameters and dimensions preceding dim
3982  * in both basic sets, the values of dimension pos in bset1 are
3983  * smaller or larger then those in bset2.
3984  *
3985  * Returns
3986  *       1 if bset1 follows bset2
3987  *      -1 if bset1 precedes bset2
3988  *       0 if bset1 and bset2 are incomparable
3989  *      -2 if some error occurred.
3990  */
3991 int isl_basic_set_compare_at(struct isl_basic_set *bset1,
3992         struct isl_basic_set *bset2, int pos)
3993 {
3994         struct isl_basic_map *bmap1 = NULL;
3995         struct isl_basic_map *bmap2 = NULL;
3996         struct isl_ctx *ctx;
3997         struct isl_vec *obj;
3998         unsigned total;
3999         isl_int num, den;
4000         enum isl_lp_result res;
4001         int cmp;
4002
4003         if (!bset1 || !bset2)
4004                 return -2;
4005
4006         bmap1 = isl_basic_map_from_basic_set(isl_basic_set_copy(bset1),
4007                                                 pos, bset1->dim - pos);
4008         bmap2 = isl_basic_map_from_basic_set(isl_basic_set_copy(bset2),
4009                                                 pos, bset2->dim - pos);
4010         if (!bmap1 || !bmap2)
4011                 goto error;
4012         bmap1 = isl_basic_map_extend(bmap1, bmap1->nparam,
4013                         bmap1->n_in, bmap1->n_out + bmap2->n_out,
4014                         bmap2->n_div, bmap2->n_eq, bmap2->n_ineq);
4015         if (!bmap1)
4016                 goto error;
4017         total = bmap1->nparam + bmap1->n_in + bmap1->n_out + bmap1->n_div;
4018         ctx = bmap1->ctx;
4019         obj = isl_vec_alloc(ctx, total);
4020         isl_seq_clr(obj->block.data, total);
4021         isl_int_set_si(obj->block.data[bmap1->nparam+bmap1->n_in], 1);
4022         isl_int_set_si(obj->block.data[bmap1->nparam+bmap1->n_in+
4023                                         bmap1->n_out-bmap2->n_out], -1);
4024         if (!obj)
4025                 goto error;
4026         bmap1 = add_constraints(bmap1, bmap2, 0, bmap1->n_out - bmap2->n_out);
4027         isl_int_init(num);
4028         isl_int_init(den);
4029         res = isl_solve_lp(bmap1, 0, obj->block.data, ctx->one, &num, &den);
4030         if (res == isl_lp_empty)
4031                 cmp = 0;
4032         else if (res == isl_lp_ok && isl_int_is_pos(num))
4033                 cmp = 1;
4034         else if ((res == isl_lp_ok && isl_int_is_neg(num)) ||
4035                   res == isl_lp_unbounded)
4036                 cmp = -1;
4037         else
4038                 cmp = -2;
4039         isl_int_clear(num);
4040         isl_int_clear(den);
4041         isl_basic_map_free(bmap1);
4042         isl_vec_free(ctx, obj);
4043         return cmp;
4044 error:
4045         isl_basic_map_free(bmap1);
4046         isl_basic_map_free(bmap2);
4047         return -2;
4048 }
4049
4050 static int isl_basic_map_fast_has_fixed_var(struct isl_basic_map *bmap,
4051         unsigned pos, isl_int *val)
4052 {
4053         int i;
4054         int d;
4055         unsigned total;
4056
4057         if (!bmap)
4058                 return -1;
4059         total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
4060         for (i = 0, d = total-1; i < bmap->n_eq && d+1 > pos; ++i) {
4061                 for (; d+1 > pos; --d)
4062                         if (!isl_int_is_zero(bmap->eq[i][1+d]))
4063                                 break;
4064                 if (d != pos)
4065                         continue;
4066                 if (isl_seq_first_non_zero(bmap->eq[i]+1, d) != -1)
4067                         return 0;
4068                 if (isl_seq_first_non_zero(bmap->eq[i]+1+d+1, total-d-1) != -1)
4069                         return 0;
4070                 if (!isl_int_is_one(bmap->eq[i][1+d]))
4071                         return 0;
4072                 if (val)
4073                         isl_int_neg(*val, bmap->eq[i][0]);
4074                 return 1;
4075         }
4076         return 0;
4077 }
4078
4079 static int isl_map_fast_has_fixed_var(struct isl_map *map,
4080         unsigned pos, isl_int *val)
4081 {
4082         int i;
4083         isl_int v;
4084         isl_int tmp;
4085         int fixed;
4086
4087         if (!map)
4088                 return -1;
4089         if (map->n == 0)
4090                 return 0;
4091         if (map->n == 1)
4092                 return isl_basic_map_fast_has_fixed_var(map->p[0], pos, val);
4093         isl_int_init(v);
4094         isl_int_init(tmp);
4095         fixed = isl_basic_map_fast_has_fixed_var(map->p[0], pos, &v);
4096         for (i = 1; fixed == 1 && i < map->n; ++i) {
4097                 fixed = isl_basic_map_fast_has_fixed_var(map->p[i], pos, &tmp);
4098                 if (fixed == 1 && isl_int_ne(tmp, v))
4099                         fixed = 0;
4100         }
4101         if (val)
4102                 isl_int_set(*val, v);
4103         isl_int_clear(tmp);
4104         isl_int_clear(v);
4105         return fixed;
4106 }
4107
4108 static int isl_set_fast_has_fixed_var(struct isl_set *set, unsigned pos,
4109         isl_int *val)
4110 {
4111         return isl_map_fast_has_fixed_var((struct isl_map *)set, pos, val);
4112 }
4113
4114 /* Check if dimension dim has fixed value and if so and if val is not NULL,
4115  * then return this fixed value in *val.
4116  */
4117 int isl_set_fast_dim_is_fixed(struct isl_set *set, unsigned dim, isl_int *val)
4118 {
4119         return isl_set_fast_has_fixed_var(set, set->nparam + dim, val);
4120 }
4121
4122 /* Check if input variable in has fixed value and if so and if val is not NULL,
4123  * then return this fixed value in *val.
4124  */
4125 int isl_map_fast_input_is_fixed(struct isl_map *map, unsigned in, isl_int *val)
4126 {
4127         return isl_map_fast_has_fixed_var(map, map->nparam + in, val);
4128 }
4129
4130 /* Check if dimension dim has an (obvious) fixed lower bound and if so
4131  * and if val is not NULL, then return this lower bound in *val.
4132  */
4133 int isl_basic_set_fast_dim_has_fixed_lower_bound(struct isl_basic_set *bset,
4134         unsigned dim, isl_int *val)
4135 {
4136         int i, i_eq = -1, i_ineq = -1;
4137         isl_int *c;
4138         unsigned total;
4139
4140         if (!bset)
4141                 return -1;
4142         total = bset->nparam + bset->dim + bset->n_div;
4143         for (i = 0; i < bset->n_eq; ++i) {
4144                 if (isl_int_is_zero(bset->eq[i][1+bset->nparam+dim]))
4145                         continue;
4146                 if (i_eq != -1)
4147                         return 0;
4148                 i_eq = i;
4149         }
4150         for (i = 0; i < bset->n_ineq; ++i) {
4151                 if (!isl_int_is_pos(bset->ineq[i][1+bset->nparam+dim]))
4152                         continue;
4153                 if (i_eq != -1 || i_ineq != -1)
4154                         return 0;
4155                 i_ineq = i;
4156         }
4157         if (i_eq == -1 && i_ineq == -1)
4158                 return 0;
4159         c = i_eq != -1 ? bset->eq[i_eq] : bset->ineq[i_ineq];
4160         /* The coefficient should always be one due to normalization. */
4161         if (!isl_int_is_one(c[1+bset->nparam+dim]))
4162                 return 0;
4163         if (isl_seq_first_non_zero(c+1, bset->nparam+dim) != -1)
4164                 return 0;
4165         if (isl_seq_first_non_zero(c+1+bset->nparam+dim+1,
4166                                         total - bset->nparam - dim - 1) != -1)
4167                 return 0;
4168         if (val)
4169                 isl_int_neg(*val, c[0]);
4170         return 1;
4171 }
4172
4173 int isl_set_fast_dim_has_fixed_lower_bound(struct isl_set *set,
4174         unsigned dim, isl_int *val)
4175 {
4176         int i;
4177         isl_int v;
4178         isl_int tmp;
4179         int fixed;
4180
4181         if (!set)
4182                 return -1;
4183         if (set->n == 0)
4184                 return 0;
4185         if (set->n == 1)
4186                 return isl_basic_set_fast_dim_has_fixed_lower_bound(set->p[0],
4187                                                                 dim, val);
4188         isl_int_init(v);
4189         isl_int_init(tmp);
4190         fixed = isl_basic_set_fast_dim_has_fixed_lower_bound(set->p[0],
4191                                                                 dim, &v);
4192         for (i = 1; fixed == 1 && i < set->n; ++i) {
4193                 fixed = isl_basic_set_fast_dim_has_fixed_lower_bound(set->p[i],
4194                                                                 dim, &tmp);
4195                 if (fixed == 1 && isl_int_ne(tmp, v))
4196                         fixed = 0;
4197         }
4198         if (val)
4199                 isl_int_set(*val, v);
4200         isl_int_clear(tmp);
4201         isl_int_clear(v);
4202         return fixed;
4203 }
4204
4205 /* Remove all information from bset that is redundant in the context
4206  * of context.  In particular, equalities that are linear combinations
4207  * of those in context are remobed.  Then the inequalities that are
4208  * redundant in the context of the equalities and inequalities of
4209  * context are removed.
4210  */
4211 static struct isl_basic_set *uset_gist(struct isl_basic_set *bset,
4212         struct isl_basic_set *context)
4213 {
4214         int i;
4215         isl_int opt;
4216         struct isl_basic_set *combined;
4217         struct isl_ctx *ctx;
4218         enum isl_lp_result res = isl_lp_ok;
4219
4220         if (!bset || !context)
4221                 goto error;
4222
4223         if (context->n_eq > 0) {
4224                 struct isl_mat *T;
4225                 struct isl_mat *T2;
4226                 struct isl_ctx *ctx = context->ctx;
4227                 context = isl_basic_set_remove_equalities(context, &T, &T2);
4228                 if (!context)
4229                         goto error;
4230                 bset = isl_basic_set_preimage(ctx, bset, T);
4231                 bset = uset_gist(bset, context);
4232                 bset = isl_basic_set_preimage(ctx, bset, T2);
4233                 return bset;
4234         }
4235         combined = isl_basic_set_extend(isl_basic_set_copy(bset),
4236                         0, bset->dim, 0, context->n_eq, context->n_ineq);
4237         context = set_add_constraints(combined, context, 0);
4238         if (!context)
4239                 goto error;
4240         ctx = context->ctx;
4241         isl_int_init(opt);
4242         for (i = bset->n_ineq-1; i >= 0; --i) {
4243                 set_swap_inequality(context, i, context->n_ineq-1);
4244                 context->n_ineq--;
4245                 res = isl_solve_lp((struct isl_basic_map *)context, 0,
4246                         context->ineq[context->n_ineq]+1, ctx->one, &opt, NULL);
4247                 context->n_ineq++;
4248                 set_swap_inequality(context, i, context->n_ineq-1);
4249                 if (res == isl_lp_unbounded)
4250                         continue;
4251                 if (res == isl_lp_error)
4252                         break;
4253                 if (res == isl_lp_empty) {
4254                         bset = isl_basic_set_set_to_empty(bset);
4255                         break;
4256                 }
4257                 isl_int_add(opt, opt, context->ineq[i][0]);
4258                 if (!isl_int_is_neg(opt)) {
4259                         isl_basic_set_drop_inequality(context, i);
4260                         isl_basic_set_drop_inequality(bset, i);
4261                 }
4262         }
4263         isl_int_clear(opt);
4264         if (res == isl_lp_error)
4265                 goto error;
4266         isl_basic_set_free(context);
4267         return bset;
4268 error:
4269         isl_basic_set_free(context);
4270         isl_basic_set_free(bset);
4271         return NULL;
4272 }
4273
4274 struct isl_basic_map *isl_basic_map_gist(struct isl_basic_map *bmap,
4275         struct isl_basic_map *context)
4276 {
4277         struct isl_basic_set *bset;
4278
4279         if (!bmap || !context)
4280                 goto error;
4281
4282         context = isl_basic_map_align_divs(context, bmap);
4283         bmap = isl_basic_map_align_divs(bmap, context);
4284
4285         bset = uset_gist(isl_basic_map_underlying_set(isl_basic_map_copy(bmap)),
4286                          isl_basic_map_underlying_set(context));
4287
4288         return isl_basic_map_overlying_set(bset, bmap);
4289 error:
4290         isl_basic_map_free(bmap);
4291         isl_basic_map_free(context);
4292         return NULL;
4293 }
4294
4295 struct isl_map *isl_map_gist(struct isl_map *map, struct isl_basic_map *context)
4296 {
4297         int i;
4298
4299         map = isl_map_cow(map);
4300         if (!map || !context)
4301                 return NULL;
4302         isl_assert(map->ctx, map->nparam == context->nparam, goto error);
4303         isl_assert(map->ctx, map->n_in == context->n_in, goto error);
4304         isl_assert(map->ctx, map->n_out == context->n_out, goto error);
4305         for (i = 0; i < map->n; ++i)
4306                 context = isl_basic_map_align_divs(context, map->p[i]);
4307         for (i = 0; i < map->n; ++i) {
4308                 map->p[i] = isl_basic_map_gist(map->p[i],
4309                                                 isl_basic_map_copy(context));
4310                 if (!map->p[i])
4311                         goto error;
4312         }
4313         isl_basic_map_free(context);
4314         F_CLR(map, ISL_MAP_NORMALIZED);
4315         return map;
4316 error:
4317         isl_map_free(map);
4318         isl_basic_map_free(context);
4319         return NULL;
4320 }
4321
4322 struct isl_set *isl_set_gist(struct isl_set *set, struct isl_basic_set *context)
4323 {
4324         return (struct isl_set *)isl_map_gist((struct isl_map *)set,
4325                                         (struct isl_basic_map *)context);
4326 }
4327
4328 struct constraint {
4329         unsigned        size;
4330         isl_int         *c;
4331 };
4332
4333 static int qsort_constraint_cmp(const void *p1, const void *p2)
4334 {
4335         const struct constraint *c1 = (const struct constraint *)p1;
4336         const struct constraint *c2 = (const struct constraint *)p2;
4337         unsigned size = isl_min(c1->size, c2->size);
4338         return isl_seq_cmp(c1->c, c2->c, size);
4339 }
4340
4341 static struct isl_basic_map *isl_basic_map_sort_constraints(
4342         struct isl_basic_map *bmap)
4343 {
4344         int i;
4345         struct constraint *c;
4346         unsigned total;
4347
4348         if (!bmap)
4349                 return NULL;
4350         total = bmap->nparam + bmap->n_in + bmap->n_out + bmap->n_div;
4351         c = isl_alloc_array(bmap->ctx, struct constraint, bmap->n_ineq);
4352         if (!c)
4353                 goto error;
4354         for (i = 0; i < bmap->n_ineq; ++i) {
4355                 c[i].size = total;
4356                 c[i].c = bmap->ineq[i];
4357         }
4358         qsort(c, bmap->n_ineq, sizeof(struct constraint), qsort_constraint_cmp);
4359         for (i = 0; i < bmap->n_ineq; ++i)
4360                 bmap->ineq[i] = c[i].c;
4361         free(c);
4362         return bmap;
4363 error:
4364         isl_basic_map_free(bmap);
4365         return NULL;
4366 }
4367
4368 struct isl_basic_map *isl_basic_map_normalize(struct isl_basic_map *bmap)
4369 {
4370         if (!bmap)
4371                 return NULL;
4372         if (F_ISSET(bmap, ISL_BASIC_MAP_NORMALIZED))
4373                 return bmap;
4374         bmap = isl_basic_map_convex_hull(bmap);
4375         bmap = isl_basic_map_sort_constraints(bmap);
4376         F_SET(bmap, ISL_BASIC_MAP_NORMALIZED);
4377         return bmap;
4378 }
4379
4380 static int isl_basic_map_fast_cmp(const struct isl_basic_map *bmap1,
4381         const struct isl_basic_map *bmap2)
4382 {
4383         int i, cmp;
4384         unsigned total;
4385
4386         if (bmap1 == bmap2)
4387                 return 0;
4388         if (bmap1->nparam != bmap2->nparam)
4389                 return bmap1->nparam - bmap2->nparam;
4390         if (bmap1->n_in != bmap2->n_in)
4391                 return bmap1->n_in - bmap2->n_in;
4392         if (bmap1->n_out != bmap2->n_out)
4393                 return bmap1->n_out - bmap2->n_out;
4394         if (F_ISSET(bmap1, ISL_BASIC_MAP_EMPTY) &&
4395             F_ISSET(bmap2, ISL_BASIC_MAP_EMPTY))
4396                 return 0;
4397         if (F_ISSET(bmap1, ISL_BASIC_MAP_EMPTY))
4398                 return 1;
4399         if (F_ISSET(bmap2, ISL_BASIC_MAP_EMPTY))
4400                 return -1;
4401         if (bmap1->n_eq != bmap2->n_eq)
4402                 return bmap1->n_eq - bmap2->n_eq;
4403         if (bmap1->n_ineq != bmap2->n_ineq)
4404                 return bmap1->n_ineq - bmap2->n_ineq;
4405         if (bmap1->n_div != bmap2->n_div)
4406                 return bmap1->n_div - bmap2->n_div;
4407         total = bmap1->nparam + bmap1->n_in + bmap1->n_out + bmap1->n_div;
4408         for (i = 0; i < bmap1->n_eq; ++i) {
4409                 cmp = isl_seq_cmp(bmap1->eq[i], bmap2->eq[i], 1+total);
4410                 if (cmp)
4411                         return cmp;
4412         }
4413         for (i = 0; i < bmap1->n_ineq; ++i) {
4414                 cmp = isl_seq_cmp(bmap1->ineq[i], bmap2->ineq[i], 1+total);
4415                 if (cmp)
4416                         return cmp;
4417         }
4418         for (i = 0; i < bmap1->n_div; ++i) {
4419                 cmp = isl_seq_cmp(bmap1->div[i], bmap2->div[i], 1+1+total);
4420                 if (cmp)
4421                         return cmp;
4422         }
4423         return 0;
4424 }
4425
4426 static int isl_basic_map_fast_is_equal(struct isl_basic_map *bmap1,
4427         struct isl_basic_map *bmap2)
4428 {
4429         return isl_basic_map_fast_cmp(bmap1, bmap2) == 0;
4430 }
4431
4432 static int qsort_bmap_cmp(const void *p1, const void *p2)
4433 {
4434         const struct isl_basic_map *bmap1 = *(const struct isl_basic_map **)p1;
4435         const struct isl_basic_map *bmap2 = *(const struct isl_basic_map **)p2;
4436
4437         return isl_basic_map_fast_cmp(bmap1, bmap2);
4438 }
4439
4440 /* We normalize in place, but if anything goes wrong we need
4441  * to return NULL, so we need to make sure we don't change the
4442  * meaning of any possible other copies of map.
4443  */
4444 struct isl_map *isl_map_normalize(struct isl_map *map)
4445 {
4446         int i, j;
4447         struct isl_basic_map *bmap;
4448
4449         if (!map)
4450                 return NULL;
4451         if (F_ISSET(map, ISL_MAP_NORMALIZED))
4452                 return map;
4453         for (i = 0; i < map->n; ++i) {
4454                 bmap = isl_basic_map_normalize(isl_basic_map_copy(map->p[i]));
4455                 if (!bmap)
4456                         goto error;
4457                 isl_basic_map_free(map->p[i]);
4458                 map->p[i] = bmap;
4459         }
4460         qsort(map->p, map->n, sizeof(struct isl_basic_map *), qsort_bmap_cmp);
4461         F_SET(map, ISL_MAP_NORMALIZED);
4462         map = isl_map_remove_empty_parts(map);
4463         if (!map)
4464                 return NULL;
4465         for (i = map->n - 1; i >= 1; --i) {
4466                 if (!isl_basic_map_fast_is_equal(map->p[i-1], map->p[i]))
4467                         continue;
4468                 isl_basic_map_free(map->p[i-1]);
4469                 for (j = i; j < map->n; ++j)
4470                         map->p[j-1] = map->p[j];
4471                 map->n--;
4472         }
4473         return map;
4474 error:
4475         isl_map_free(map);
4476         return NULL;
4477
4478 }
4479
4480 struct isl_set *isl_set_normalize(struct isl_set *set)
4481 {
4482         return (struct isl_set *)isl_map_normalize((struct isl_map *)set);
4483 }
4484
4485 int isl_map_fast_is_equal(struct isl_map *map1, struct isl_map *map2)
4486 {
4487         int i;
4488         int equal;
4489
4490         if (!map1 || !map2)
4491                 return -1;
4492
4493         if (map1 == map2)
4494                 return 1;
4495         if (map1->nparam != map2->nparam)
4496                 return 0;
4497         if (map1->n_in != map2->n_in)
4498                 return 0;
4499         if (map1->n_out != map2->n_out)
4500                 return 0;
4501
4502         map1 = isl_map_copy(map1);
4503         map2 = isl_map_copy(map2);
4504         map1 = isl_map_normalize(map1);
4505         map2 = isl_map_normalize(map2);
4506         if (!map1 || !map2)
4507                 goto error;
4508         equal = map1->n == map2->n;
4509         for (i = 0; equal && i < map1->n; ++i) {
4510                 equal = isl_basic_map_fast_is_equal(map1->p[i], map2->p[i]);
4511                 if (equal < 0)
4512                         goto error;
4513         }
4514         isl_map_free(map1);
4515         isl_map_free(map2);
4516         return equal;
4517 error:
4518         isl_map_free(map1);
4519         isl_map_free(map2);
4520         return -1;
4521 }
4522
4523 int isl_set_fast_is_equal(struct isl_set *set1, struct isl_set *set2)
4524 {
4525         return isl_map_fast_is_equal((struct isl_map *)set1,
4526                                                 (struct isl_map *)set2);
4527 }
4528
4529 /* Return an interval that ranges from min to max (inclusive)
4530  */
4531 struct isl_basic_set *isl_basic_set_interval(struct isl_ctx *ctx,
4532         isl_int min, isl_int max)
4533 {
4534         int k;
4535         struct isl_basic_set *bset = NULL;
4536
4537         bset = isl_basic_set_alloc(ctx, 0, 1, 0, 0, 2);
4538         if (!bset)
4539                 goto error;
4540
4541         k = isl_basic_set_alloc_inequality(bset);
4542         if (k < 0)
4543                 goto error;
4544         isl_int_set_si(bset->ineq[k][1], 1);
4545         isl_int_neg(bset->ineq[k][0], min);
4546
4547         k = isl_basic_set_alloc_inequality(bset);
4548         if (k < 0)
4549                 goto error;
4550         isl_int_set_si(bset->ineq[k][1], -1);
4551         isl_int_set(bset->ineq[k][0], max);
4552
4553         return bset;
4554 error:
4555         isl_basic_set_free(bset);
4556         return NULL;
4557 }
4558
4559 /* Return the Cartesian product of the basic sets in list (in the given order).
4560  */
4561 struct isl_basic_set *isl_basic_set_product(struct isl_basic_set_list *list)
4562 {
4563         int i;
4564         unsigned dim;
4565         unsigned nparam;
4566         unsigned extra;
4567         unsigned n_eq;
4568         unsigned n_ineq;
4569         struct isl_basic_set *product = NULL;
4570
4571         if (!list)
4572                 goto error;
4573         isl_assert(list->ctx, list->n > 0, goto error);
4574         isl_assert(list->ctx, list->p[0], goto error);
4575         nparam = list->p[0]->nparam;
4576         dim = list->p[0]->dim;
4577         extra = list->p[0]->n_div;
4578         n_eq = list->p[0]->n_eq;
4579         n_ineq = list->p[0]->n_ineq;
4580         for (i = 1; i < list->n; ++i) {
4581                 isl_assert(list->ctx, list->p[i], goto error);
4582                 isl_assert(list->ctx, nparam == list->p[i]->nparam, goto error);
4583                 dim += list->p[i]->dim;
4584                 extra += list->p[i]->n_div;
4585                 n_eq += list->p[i]->n_eq;
4586                 n_ineq += list->p[i]->n_ineq;
4587         }
4588         product = isl_basic_set_alloc(list->ctx, nparam, dim, extra,
4589                                         n_eq, n_ineq);
4590         if (!product)
4591                 goto error;
4592         dim = 0;
4593         for (i = 0; i < list->n; ++i) {
4594                 set_add_constraints(product,
4595                                         isl_basic_set_copy(list->p[i]), dim);
4596                 dim += list->p[i]->dim;
4597         }
4598         isl_basic_set_list_free(list);
4599         return product;
4600 error:
4601         isl_basic_set_free(product);
4602         isl_basic_set_list_free(list);
4603         return NULL;
4604 }