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