isl_coalesce.c

   1 /*
   2  * Copyright 2008-2009 Katholieke Universiteit Leuven
   3  * Copyright 2010      INRIA Saclay
   4  *
   5  * Use of this software is governed by the GNU LGPLv2.1 license
   6  *
   7  * Written by Sven Verdoolaege, K.U.Leuven, Departement
   8  * Computerwetenschappen, Celestijnenlaan 200A, B-3001 Leuven, Belgium
   9  * and INRIA Saclay - Ile-de-France, Parc Club Orsay Universite,
  10  * ZAC des vignes, 4 rue Jacques Monod, 91893 Orsay, France
  11  */
  12
  13 #include "isl_map_private.h"
  14 #include "isl_seq.h"
  15 #include "isl_tab.h"
  16
  17 #define STATUS_ERROR            -1
  18 #define STATUS_REDUNDANT         1
  19 #define STATUS_VALID             2
  20 #define STATUS_SEPARATE          3
  21 #define STATUS_CUT               4
  22 #define STATUS_ADJ_EQ            5
  23 #define STATUS_ADJ_INEQ          6
  24
  25 static int status_in(isl_int *ineq, struct isl_tab *tab)
  26 {
  27         enum isl_ineq_type type = isl_tab_ineq_type(tab, ineq);
  28         switch (type) {
  29         case isl_ineq_error:            return STATUS_ERROR;
  30         case isl_ineq_redundant:        return STATUS_VALID;
  31         case isl_ineq_separate:         return STATUS_SEPARATE;
  32         case isl_ineq_cut:              return STATUS_CUT;
  33         case isl_ineq_adj_eq:           return STATUS_ADJ_EQ;
  34         case isl_ineq_adj_ineq:         return STATUS_ADJ_INEQ;
  35         }
  36 }
  37
  38 /* Compute the position of the equalities of basic map "i"
  39  * with respect to basic map "j".
  40  * The resulting array has twice as many entries as the number
  41  * of equalities corresponding to the two inequalties to which
  42  * each equality corresponds.
  43  */
  44 static int *eq_status_in(struct isl_map *map, int i, int j,
  45         struct isl_tab **tabs)
  46 {
  47         int k, l;
  48         int *eq = isl_calloc_array(map->ctx, int, 2 * map->p[i]->n_eq);
  49         unsigned dim;
  50
  51         dim = isl_basic_map_total_dim(map->p[i]);
  52         for (k = 0; k < map->p[i]->n_eq; ++k) {
  53                 for (l = 0; l < 2; ++l) {
  54                         isl_seq_neg(map->p[i]->eq[k], map->p[i]->eq[k], 1+dim);
  55                         eq[2 * k + l] = status_in(map->p[i]->eq[k], tabs[j]);
  56                         if (eq[2 * k + l] == STATUS_ERROR)
  57                                 goto error;
  58                 }
  59                 if (eq[2 * k] == STATUS_SEPARATE ||
  60                     eq[2 * k + 1] == STATUS_SEPARATE)
  61                         break;
  62         }
  63
  64         return eq;
  65 error:
  66         free(eq);
  67         return NULL;
  68 }
  69
  70 /* Compute the position of the inequalities of basic map "i"
  71  * with respect to basic map "j".
  72  */
  73 static int *ineq_status_in(struct isl_map *map, int i, int j,
  74         struct isl_tab **tabs)
  75 {
  76         int k;
  77         unsigned n_eq = map->p[i]->n_eq;
  78         int *ineq = isl_calloc_array(map->ctx, int, map->p[i]->n_ineq);
  79
  80         for (k = 0; k < map->p[i]->n_ineq; ++k) {
  81                 if (isl_tab_is_redundant(tabs[i], n_eq + k)) {
  82                         ineq[k] = STATUS_REDUNDANT;
  83                         continue;
  84                 }
  85                 ineq[k] = status_in(map->p[i]->ineq[k], tabs[j]);
  86                 if (ineq[k] == STATUS_ERROR)
  87                         goto error;
  88                 if (ineq[k] == STATUS_SEPARATE)
  89                         break;
  90         }
  91
  92         return ineq;
  93 error:
  94         free(ineq);
  95         return NULL;
  96 }
  97
  98 static int any(int *con, unsigned len, int status)
  99 {
 100         int i;
 101
 102         for (i = 0; i < len ; ++i)
 103                 if (con[i] == status)
 104                         return 1;
 105         return 0;
 106 }
 107
 108 static int count(int *con, unsigned len, int status)
 109 {
 110         int i;
 111         int c = 0;
 112
 113         for (i = 0; i < len ; ++i)
 114                 if (con[i] == status)
 115                         c++;
 116         return c;
 117 }
 118
 119 static int all(int *con, unsigned len, int status)
 120 {
 121         int i;
 122
 123         for (i = 0; i < len ; ++i) {
 124                 if (con[i] == STATUS_REDUNDANT)
 125                         continue;
 126                 if (con[i] != status)
 127                         return 0;
 128         }
 129         return 1;
 130 }
 131
 132 static void drop(struct isl_map *map, int i, struct isl_tab **tabs)
 133 {
 134         isl_basic_map_free(map->p[i]);
 135         isl_tab_free(tabs[i]);
 136
 137         if (i != map->n - 1) {
 138                 map->p[i] = map->p[map->n - 1];
 139                 tabs[i] = tabs[map->n - 1];
 140         }
 141         tabs[map->n - 1] = NULL;
 142         map->n--;
 143 }
 144
 145 /* Replace the pair of basic maps i and j by the basic map bounded
 146  * by the valid constraints in both basic maps and the constraint
 147  * in extra (if not NULL).
 148  */
 149 static int fuse(struct isl_map *map, int i, int j,
 150         struct isl_tab **tabs, int *eq_i, int *ineq_i, int *eq_j, int *ineq_j,
 151         __isl_keep isl_mat *extra)
 152 {
 153         int k, l;
 154         struct isl_basic_map *fused = NULL;
 155         struct isl_tab *fused_tab = NULL;
 156         unsigned total = isl_basic_map_total_dim(map->p[i]);
 157         unsigned extra_rows = extra ? extra->n_row : 0;
 158
 159         fused = isl_basic_map_alloc_dim(isl_dim_copy(map->p[i]->dim),
 160                         map->p[i]->n_div,
 161                         map->p[i]->n_eq + map->p[j]->n_eq,
 162                         map->p[i]->n_ineq + map->p[j]->n_ineq + extra_rows);
 163         if (!fused)
 164                 goto error;
 165
 166         for (k = 0; k < map->p[i]->n_eq; ++k) {
 167                 if (eq_i && (eq_i[2 * k] != STATUS_VALID ||
 168                              eq_i[2 * k + 1] != STATUS_VALID))
 169                         continue;
 170                 l = isl_basic_map_alloc_equality(fused);
 171                 if (l < 0)
 172                         goto error;
 173                 isl_seq_cpy(fused->eq[l], map->p[i]->eq[k], 1 + total);
 174         }
 175
 176         for (k = 0; k < map->p[j]->n_eq; ++k) {
 177                 if (eq_j && (eq_j[2 * k] != STATUS_VALID ||
 178                              eq_j[2 * k + 1] != STATUS_VALID))
 179                         continue;
 180                 l = isl_basic_map_alloc_equality(fused);
 181                 if (l < 0)
 182                         goto error;
 183                 isl_seq_cpy(fused->eq[l], map->p[j]->eq[k], 1 + total);
 184         }
 185
 186         for (k = 0; k < map->p[i]->n_ineq; ++k) {
 187                 if (ineq_i[k] != STATUS_VALID)
 188                         continue;
 189                 l = isl_basic_map_alloc_inequality(fused);
 190                 if (l < 0)
 191                         goto error;
 192                 isl_seq_cpy(fused->ineq[l], map->p[i]->ineq[k], 1 + total);
 193         }
 194
 195         for (k = 0; k < map->p[j]->n_ineq; ++k) {
 196                 if (ineq_j[k] != STATUS_VALID)
 197                         continue;
 198                 l = isl_basic_map_alloc_inequality(fused);
 199                 if (l < 0)
 200                         goto error;
 201                 isl_seq_cpy(fused->ineq[l], map->p[j]->ineq[k], 1 + total);
 202         }
 203
 204         for (k = 0; k < map->p[i]->n_div; ++k) {
 205                 int l = isl_basic_map_alloc_div(fused);
 206                 if (l < 0)
 207                         goto error;
 208                 isl_seq_cpy(fused->div[l], map->p[i]->div[k], 1 + 1 + total);
 209         }
 210
 211         for (k = 0; k < extra_rows; ++k) {
 212                 l = isl_basic_map_alloc_inequality(fused);
 213                 if (l < 0)
 214                         goto error;
 215                 isl_seq_cpy(fused->ineq[l], extra->row[k], 1 + total);
 216         }
 217
 218         fused = isl_basic_map_gauss(fused, NULL);
 219         ISL_F_SET(fused, ISL_BASIC_MAP_FINAL);
 220         if (ISL_F_ISSET(map->p[i], ISL_BASIC_MAP_RATIONAL) &&
 221             ISL_F_ISSET(map->p[j], ISL_BASIC_MAP_RATIONAL))
 222                 ISL_F_SET(fused, ISL_BASIC_MAP_RATIONAL);
 223
 224         fused_tab = isl_tab_from_basic_map(fused);
 225         if (isl_tab_detect_redundant(fused_tab) < 0)
 226                 goto error;
 227
 228         isl_basic_map_free(map->p[i]);
 229         map->p[i] = fused;
 230         isl_tab_free(tabs[i]);
 231         tabs[i] = fused_tab;
 232         drop(map, j, tabs);
 233
 234         return 1;
 235 error:
 236         isl_tab_free(fused_tab);
 237         isl_basic_map_free(fused);
 238         return -1;
 239 }
 240
 241 /* Given a pair of basic maps i and j such that all constraints are either
 242  * "valid" or "cut", check if the facets corresponding to the "cut"
 243  * constraints of i lie entirely within basic map j.
 244  * If so, replace the pair by the basic map consisting of the valid
 245  * constraints in both basic maps.
 246  *
 247  * To see that we are not introducing any extra points, call the
 248  * two basic maps A and B and the resulting map U and let x
 249  * be an element of U \setminus ( A \cup B ).
 250  * Then there is a pair of cut constraints c_1 and c_2 in A and B such that x
 251  * violates them.  Let X be the intersection of U with the opposites
 252  * of these constraints.  Then x \in X.
 253  * The facet corresponding to c_1 contains the corresponding facet of A.
 254  * This facet is entirely contained in B, so c_2 is valid on the facet.
 255  * However, since it is also (part of) a facet of X, -c_2 is also valid
 256  * on the facet.  This means c_2 is saturated on the facet, so c_1 and
 257  * c_2 must be opposites of each other, but then x could not violate
 258  * both of them.
 259  */
 260 static int check_facets(struct isl_map *map, int i, int j,
 261         struct isl_tab **tabs, int *ineq_i, int *ineq_j)
 262 {
 263         int k, l;
 264         struct isl_tab_undo *snap;
 265         unsigned n_eq = map->p[i]->n_eq;
 266
 267         snap = isl_tab_snap(tabs[i]);
 268
 269         for (k = 0; k < map->p[i]->n_ineq; ++k) {
 270                 if (ineq_i[k] != STATUS_CUT)
 271                         continue;
 272                 tabs[i] = isl_tab_select_facet(tabs[i], n_eq + k);
 273                 for (l = 0; l < map->p[j]->n_ineq; ++l) {
 274                         int stat;
 275                         if (ineq_j[l] != STATUS_CUT)
 276                                 continue;
 277                         stat = status_in(map->p[j]->ineq[l], tabs[i]);
 278                         if (stat != STATUS_VALID)
 279                                 break;
 280                 }
 281                 if (isl_tab_rollback(tabs[i], snap) < 0)
 282                         return -1;
 283                 if (l < map->p[j]->n_ineq)
 284                         break;
 285         }
 286
 287         if (k < map->p[i]->n_ineq)
 288                 /* BAD CUT PAIR */
 289                 return 0;
 290         return fuse(map, i, j, tabs, NULL, ineq_i, NULL, ineq_j, NULL);
 291 }
 292
 293 /* Both basic maps have at least one inequality with and adjacent
 294  * (but opposite) inequality in the other basic map.
 295  * Check that there are no cut constraints and that there is only
 296  * a single pair of adjacent inequalities.
 297  * If so, we can replace the pair by a single basic map described
 298  * by all but the pair of adjacent inequalities.
 299  * Any additional points introduced lie strictly between the two
 300  * adjacent hyperplanes and can therefore be integral.
 301  *
 302  *        ____                    _____
 303  *       /    ||\                /     \
 304  *      /     || \              /       \
 305  *      \     ||  \     =>      \        \
 306  *       \    ||  /              \       /
 307  *        \___||_/                \_____/
 308  *
 309  * The test for a single pair of adjancent inequalities is important
 310  * for avoiding the combination of two basic maps like the following
 311  *
 312  *       /|
 313  *      / |
 314  *     /__|
 315  *         _____
 316  *         |   |
 317  *         |   |
 318  *         |___|
 319  */
 320 static int check_adj_ineq(struct isl_map *map, int i, int j,
 321         struct isl_tab **tabs, int *ineq_i, int *ineq_j)
 322 {
 323         int changed = 0;
 324
 325         if (any(ineq_i, map->p[i]->n_ineq, STATUS_CUT) ||
 326             any(ineq_j, map->p[j]->n_ineq, STATUS_CUT))
 327                 /* ADJ INEQ CUT */
 328                 ;
 329         else if (count(ineq_i, map->p[i]->n_ineq, STATUS_ADJ_INEQ) == 1 &&
 330                  count(ineq_j, map->p[j]->n_ineq, STATUS_ADJ_INEQ) == 1)
 331                 changed = fuse(map, i, j, tabs, NULL, ineq_i, NULL, ineq_j, NULL);
 332         /* else ADJ INEQ TOO MANY */
 333
 334         return changed;
 335 }
 336
 337 /* Check if basic map "i" contains the basic map represented
 338  * by the tableau "tab".
 339  */
 340 static int contains(struct isl_map *map, int i, int *ineq_i,
 341         struct isl_tab *tab)
 342 {
 343         int k, l;
 344         unsigned dim;
 345
 346         dim = isl_basic_map_total_dim(map->p[i]);
 347         for (k = 0; k < map->p[i]->n_eq; ++k) {
 348                 for (l = 0; l < 2; ++l) {
 349                         int stat;
 350                         isl_seq_neg(map->p[i]->eq[k], map->p[i]->eq[k], 1+dim);
 351                         stat = status_in(map->p[i]->eq[k], tab);
 352                         if (stat != STATUS_VALID)
 353                                 return 0;
 354                 }
 355         }
 356
 357         for (k = 0; k < map->p[i]->n_ineq; ++k) {
 358                 int stat;
 359                 if (ineq_i[k] == STATUS_REDUNDANT)
 360                         continue;
 361                 stat = status_in(map->p[i]->ineq[k], tab);
 362                 if (stat != STATUS_VALID)
 363                         return 0;
 364         }
 365         return 1;
 366 }
 367
 368 /* Basic map "i" has an inequality "k" that is adjacent to some equality
 369  * of basic map "j".  All the other inequalities are valid for "j".
 370  * Check if basic map "j" forms an extension of basic map "i".
 371  *
 372  * In particular, we relax constraint "k", compute the corresponding
 373  * facet and check whether it is included in the other basic map.
 374  * If so, we know that relaxing the constraint extends the basic
 375  * map with exactly the other basic map (we already know that this
 376  * other basic map is included in the extension, because there
 377  * were no "cut" inequalities in "i") and we can replace the
 378  * two basic maps by thie extension.
 379  *        ____                    _____
 380  *       /    ||                 /     |
 381  *      /     ||                /      |
 382  *      \     ||        =>      \      |
 383  *       \    ||                 \     |
 384  *        \___||                  \____|
 385  */
 386 static int is_extension(struct isl_map *map, int i, int j, int k,
 387         struct isl_tab **tabs, int *eq_i, int *ineq_i, int *eq_j, int *ineq_j)
 388 {
 389         int changed = 0;
 390         int super;
 391         struct isl_tab_undo *snap, *snap2;
 392         unsigned n_eq = map->p[i]->n_eq;
 393
 394         snap = isl_tab_snap(tabs[i]);
 395         tabs[i] = isl_tab_relax(tabs[i], n_eq + k);
 396         snap2 = isl_tab_snap(tabs[i]);
 397         tabs[i] = isl_tab_select_facet(tabs[i], n_eq + k);
 398         super = contains(map, j, ineq_j, tabs[i]);
 399         if (super) {
 400                 if (isl_tab_rollback(tabs[i], snap2) < 0)
 401                         return -1;
 402                 map->p[i] = isl_basic_map_cow(map->p[i]);
 403                 if (!map->p[i])
 404                         return -1;
 405                 isl_int_add_ui(map->p[i]->ineq[k][0], map->p[i]->ineq[k][0], 1);
 406                 ISL_F_SET(map->p[i], ISL_BASIC_MAP_FINAL);
 407                 drop(map, j, tabs);
 408                 changed = 1;
 409         } else
 410                 if (isl_tab_rollback(tabs[i], snap) < 0)
 411                         return -1;
 412
 413         return changed;
 414 }
 415
 416 /* For each non-redundant constraint in "bmap" (as determined by "tab"),
 417  * wrap the constraint around "bound" such that it includes the whole
 418  * set "set" and append the resulting constraint to "wraps".
 419  * "wraps" is assumed to have been pre-allocated to the appropriate size.
 420  * wraps->n_row is the number of actual wrapped constraints that have
 421  * been added.
 422  * If any of the wrapping problems results in a constraint that is
 423  * identical to "bound", then this means that "set" is unbounded in such
 424  * way that no wrapping is possible.  If this happens then wraps->n_row
 425  * is reset to zero.
 426  */
 427 static int add_wraps(__isl_keep isl_mat *wraps, __isl_keep isl_basic_map *bmap,
 428         struct isl_tab *tab, isl_int *bound, __isl_keep isl_set *set)
 429 {
 430         int l;
 431         int w;
 432         unsigned total = isl_basic_map_total_dim(bmap);
 433
 434         w = wraps->n_row;
 435
 436         for (l = 0; l < bmap->n_ineq; ++l) {
 437                 if (isl_seq_is_neg(bound, bmap->ineq[l], 1 + total))
 438                         continue;
 439                 if (isl_seq_eq(bound, bmap->ineq[l], 1 + total))
 440                         continue;
 441                 if (isl_tab_is_redundant(tab, bmap->n_eq + l))
 442                         continue;
 443
 444                 isl_seq_cpy(wraps->row[w], bound, 1 + total);
 445                 if (!isl_set_wrap_facet(set, wraps->row[w], bmap->ineq[l]))
 446                         return -1;
 447                 if (isl_seq_eq(wraps->row[w], bound, 1 + total))
 448                         goto unbounded;
 449                 ++w;
 450         }
 451         for (l = 0; l < bmap->n_eq; ++l) {
 452                 if (isl_seq_is_neg(bound, bmap->eq[l], 1 + total))
 453                         continue;
 454                 if (isl_seq_eq(bound, bmap->eq[l], 1 + total))
 455                         continue;
 456
 457                 isl_seq_cpy(wraps->row[w], bound, 1 + total);
 458                 isl_seq_neg(wraps->row[w + 1], bmap->eq[l], 1 + total);
 459                 if (!isl_set_wrap_facet(set, wraps->row[w], wraps->row[w + 1]))
 460                         return -1;
 461                 if (isl_seq_eq(wraps->row[w], bound, 1 + total))
 462                         goto unbounded;
 463                 ++w;
 464
 465                 isl_seq_cpy(wraps->row[w], bound, 1 + total);
 466                 if (!isl_set_wrap_facet(set, wraps->row[w], bmap->eq[l]))
 467                         return -1;
 468                 if (isl_seq_eq(wraps->row[w], bound, 1 + total))
 469                         goto unbounded;
 470                 ++w;
 471         }
 472
 473         wraps->n_row = w;
 474         return 0;
 475 unbounded:
 476         wraps->n_row = 0;
 477         return 0;
 478 }
 479
 480 /* Check if the constraints in "wraps" from "first" until the last
 481  * are all valid for the basic set represented by "tab".
 482  * If not, wraps->n_row is set to zero.
 483  */
 484 static int check_wraps(__isl_keep isl_mat *wraps, int first,
 485         struct isl_tab *tab)
 486 {
 487         int i;
 488
 489         for (i = first; i < wraps->n_row; ++i) {
 490                 enum isl_ineq_type type;
 491                 type = isl_tab_ineq_type(tab, wraps->row[i]);
 492                 if (type == isl_ineq_error)
 493                         return -1;
 494                 if (type == isl_ineq_redundant)
 495                         continue;
 496                 wraps->n_row = 0;
 497                 return 0;
 498         }
 499
 500         return 0;
 501 }
 502
 503 /* Return a set that corresponds to the non-redudant constraints
 504  * (as recorded in tab) of bmap.
 505  *
 506  * It's important to remove the redundant constraints as some
 507  * of the other constraints may have been modified after the
 508  * constraints were marked redundant.
 509  * In particular, a constraint may have been relaxed.
 510  * Redundant constraints are ignored when a constraint is relaxed
 511  * and should therefore continue to be ignored ever after.
 512  * Otherwise, the relaxation might be thwarted by some of
 513  * these constraints.
 514  */
 515 static __isl_give isl_set *set_from_updated_bmap(__isl_keep isl_basic_map *bmap,
 516         struct isl_tab *tab)
 517 {
 518         bmap = isl_basic_map_copy(bmap);
 519         bmap = isl_basic_map_cow(bmap);
 520         bmap = isl_basic_map_update_from_tab(bmap, tab);
 521         return isl_set_from_basic_set(isl_basic_map_underlying_set(bmap));
 522 }
 523
 524 /* Given a basic set i with a constraint k that is adjacent to either the
 525  * whole of basic set j or a facet of basic set j, check if we can wrap
 526  * both the facet corresponding to k and the facet of j (or the whole of j)
 527  * around their ridges to include the other set.
 528  * If so, replace the pair of basic sets by their union.
 529  *
 530  * All constraints of i (except k) are assumed to be valid for j.
 531  *
 532  * However, the constraints of j may not be valid for i and so
 533  * we have to check that the wrapping constraints for j are valid for i.
 534  *
 535  * In the case where j has a facet adjacent to i, tab[j] is assumed
 536  * to have been restricted to this facet, so that the non-redundant
 537  * constraints in tab[j] are the ridges of the facet.
 538  * Note that for the purpose of wrapping, it does not matter whether
 539  * we wrap the ridges of i around the whole of j or just around
 540  * the facet since all the other constraints are assumed to be valid for j.
 541  * In practice, we wrap to include the whole of j.
 542  *        ____                    _____
 543  *       /    |                  /     \
 544  *      /     ||                /      |
 545  *      \     ||        =>      \      |
 546  *       \    ||                 \     |
 547  *        \___||                  \____|
 548  *
 549  */
 550 static int can_wrap_in_facet(struct isl_map *map, int i, int j, int k,
 551         struct isl_tab **tabs, int *eq_i, int *ineq_i, int *eq_j, int *ineq_j)
 552 {
 553         int changed = 0;
 554         struct isl_mat *wraps = NULL;
 555         struct isl_set *set_i = NULL;
 556         struct isl_set *set_j = NULL;
 557         struct isl_vec *bound = NULL;
 558         unsigned total = isl_basic_map_total_dim(map->p[i]);
 559         struct isl_tab_undo *snap;
 560         int n;
 561
 562         set_i = set_from_updated_bmap(map->p[i], tabs[i]);
 563         set_j = set_from_updated_bmap(map->p[j], tabs[j]);
 564         wraps = isl_mat_alloc(map->ctx, 2 * (map->p[i]->n_eq + map->p[j]->n_eq) +
 565                                         map->p[i]->n_ineq + map->p[j]->n_ineq,
 566                                         1 + total);
 567         bound = isl_vec_alloc(map->ctx, 1 + total);
 568         if (!set_i || !set_j || !wraps || !bound)
 569                 goto error;
 570
 571         isl_seq_cpy(bound->el, map->p[i]->ineq[k], 1 + total);
 572         isl_int_add_ui(bound->el[0], bound->el[0], 1);
 573
 574         isl_seq_cpy(wraps->row[0], bound->el, 1 + total);
 575         wraps->n_row = 1;
 576
 577         if (add_wraps(wraps, map->p[j], tabs[j], bound->el, set_i) < 0)
 578                 goto error;
 579         if (!wraps->n_row)
 580                 goto unbounded;
 581
 582         snap = isl_tab_snap(tabs[i]);
 583
 584         tabs[i] = isl_tab_select_facet(tabs[i], map->p[i]->n_eq + k);
 585         if (isl_tab_detect_redundant(tabs[i]) < 0)
 586                 goto error;
 587
 588         isl_seq_neg(bound->el, map->p[i]->ineq[k], 1 + total);
 589
 590         n = wraps->n_row;
 591         if (add_wraps(wraps, map->p[i], tabs[i], bound->el, set_j) < 0)
 592                 goto error;
 593
 594         if (isl_tab_rollback(tabs[i], snap) < 0)
 595                 goto error;
 596         if (check_wraps(wraps, n, tabs[i]) < 0)
 597                 goto error;
 598         if (!wraps->n_row)
 599                 goto unbounded;
 600
 601         changed = fuse(map, i, j, tabs, eq_i, ineq_i, eq_j, ineq_j, wraps);
 602
 603 unbounded:
 604         isl_mat_free(wraps);
 605
 606         isl_set_free(set_i);
 607         isl_set_free(set_j);
 608
 609         isl_vec_free(bound);
 610
 611         return changed;
 612 error:
 613         isl_vec_free(bound);
 614         isl_mat_free(wraps);
 615         isl_set_free(set_i);
 616         isl_set_free(set_j);
 617         return -1;
 618 }
 619
 620 /* Given two basic sets i and j such that i has exactly one cut constraint,
 621  * check if we can wrap the corresponding facet around its ridges to include
 622  * the other basic set (and nothing else).
 623  * If so, replace the pair by their union.
 624  *
 625  * We first check if j has a facet adjacent to the cut constraint of i.
 626  * If so, we try to wrap in the facet.
 627  *        ____                    _____
 628  *       / ___|_                 /     \
 629  *      / |    |                /      |
 630  *      \ |    |        =>      \      |
 631  *       \|____|                 \     |
 632  *        \___|                   \____/
 633  */
 634 static int can_wrap_in_set(struct isl_map *map, int i, int j,
 635         struct isl_tab **tabs, int *ineq_i, int *ineq_j)
 636 {
 637         int changed = 0;
 638         int k, l;
 639         unsigned total = isl_basic_map_total_dim(map->p[i]);
 640         struct isl_tab_undo *snap;
 641
 642         for (k = 0; k < map->p[i]->n_ineq; ++k)
 643                 if (ineq_i[k] == STATUS_CUT)
 644                         break;
 645
 646         isl_assert(map->ctx, k < map->p[i]->n_ineq, return -1);
 647
 648         isl_int_add_ui(map->p[i]->ineq[k][0], map->p[i]->ineq[k][0], 1);
 649         for (l = 0; l < map->p[j]->n_ineq; ++l) {
 650                 if (isl_tab_is_redundant(tabs[j], map->p[j]->n_eq + l))
 651                         continue;
 652                 if (isl_seq_eq(map->p[i]->ineq[k],
 653                                 map->p[j]->ineq[l], 1 + total))
 654                         break;
 655         }
 656         isl_int_sub_ui(map->p[i]->ineq[k][0], map->p[i]->ineq[k][0], 1);
 657
 658         if (l >= map->p[j]->n_ineq)
 659                 return 0;
 660
 661         snap = isl_tab_snap(tabs[j]);
 662         tabs[j] = isl_tab_select_facet(tabs[j], map->p[j]->n_eq + l);
 663         if (isl_tab_detect_redundant(tabs[j]) < 0)
 664                 return -1;
 665
 666         changed = can_wrap_in_facet(map, i, j, k, tabs, NULL, ineq_i, NULL, ineq_j);
 667
 668         if (!changed && isl_tab_rollback(tabs[j], snap) < 0)
 669                 return -1;
 670
 671         return changed;
 672 }
 673
 674 /* Check if either i or j has a single cut constraint that can
 675  * be used to wrap in (a facet of) the other basic set.
 676  * if so, replace the pair by their union.
 677  */
 678 static int check_wrap(struct isl_map *map, int i, int j,
 679         struct isl_tab **tabs, int *ineq_i, int *ineq_j)
 680 {
 681         int changed = 0;
 682
 683         if (count(ineq_i, map->p[i]->n_ineq, STATUS_CUT) == 1)
 684                 changed = can_wrap_in_set(map, i, j, tabs, ineq_i, ineq_j);
 685         if (changed)
 686                 return changed;
 687
 688         if (count(ineq_j, map->p[j]->n_ineq, STATUS_CUT) == 1)
 689                 changed = can_wrap_in_set(map, j, i, tabs, ineq_j, ineq_i);
 690         return changed;
 691 }
 692
 693 /* At least one of the basic maps has an equality that is adjacent
 694  * to inequality.  Make sure that only one of the basic maps has
 695  * such an equality and that the other basic map has exactly one
 696  * inequality adjacent to an equality.
 697  * We call the basic map that has the inequality "i" and the basic
 698  * map that has the equality "j".
 699  * If "i" has any "cut" inequality, then relaxing the inequality
 700  * by one would not result in a basic map that contains the other
 701  * basic map.
 702  */
 703 static int check_adj_eq(struct isl_map *map, int i, int j,
 704         struct isl_tab **tabs, int *eq_i, int *ineq_i, int *eq_j, int *ineq_j)
 705 {
 706         int changed = 0;
 707         int k;
 708
 709         if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_ADJ_INEQ) &&
 710             any(eq_j, 2 * map->p[j]->n_eq, STATUS_ADJ_INEQ))
 711                 /* ADJ EQ TOO MANY */
 712                 return 0;
 713
 714         if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_ADJ_INEQ))
 715                 return check_adj_eq(map, j, i, tabs,
 716                                         eq_j, ineq_j, eq_i, ineq_i);
 717
 718         /* j has an equality adjacent to an inequality in i */
 719
 720         if (any(ineq_i, map->p[i]->n_ineq, STATUS_CUT))
 721                 /* ADJ EQ CUT */
 722                 return 0;
 723         if (count(eq_j, 2 * map->p[j]->n_eq, STATUS_ADJ_INEQ) != 1 ||
 724             count(ineq_i, map->p[i]->n_ineq, STATUS_ADJ_EQ) != 1 ||
 725             any(ineq_j, map->p[j]->n_ineq, STATUS_ADJ_EQ) ||
 726             any(ineq_i, map->p[i]->n_ineq, STATUS_ADJ_INEQ) ||
 727             any(ineq_j, map->p[j]->n_ineq, STATUS_ADJ_INEQ))
 728                 /* ADJ EQ TOO MANY */
 729                 return 0;
 730
 731         for (k = 0; k < map->p[i]->n_ineq ; ++k)
 732                 if (ineq_i[k] == STATUS_ADJ_EQ)
 733                         break;
 734
 735         changed = is_extension(map, i, j, k, tabs, eq_i, ineq_i, eq_j, ineq_j);
 736         if (changed)
 737                 return changed;
 738
 739         changed = can_wrap_in_facet(map, i, j, k, tabs, eq_i, ineq_i, eq_j, ineq_j);
 740
 741         return changed;
 742 }
 743
 744 /* Check if the union of the given pair of basic maps
 745  * can be represented by a single basic map.
 746  * If so, replace the pair by the single basic map and return 1.
 747  * Otherwise, return 0;
 748  *
 749  * We first check the effect of each constraint of one basic map
 750  * on the other basic map.
 751  * The constraint may be
 752  *      redundant       the constraint is redundant in its own
 753  *                      basic map and should be ignore and removed
 754  *                      in the end
 755  *      valid           all (integer) points of the other basic map
 756  *                      satisfy the constraint
 757  *      separate        no (integer) point of the other basic map
 758  *                      satisfies the constraint
 759  *      cut             some but not all points of the other basic map
 760  *                      satisfy the constraint
 761  *      adj_eq          the given constraint is adjacent (on the outside)
 762  *                      to an equality of the other basic map
 763  *      adj_ineq        the given constraint is adjacent (on the outside)
 764  *                      to an inequality of the other basic map
 765  *
 766  * We consider six cases in which we can replace the pair by a single
 767  * basic map.  We ignore all "redundant" constraints.
 768  *
 769  *      1. all constraints of one basic map are valid
 770  *              => the other basic map is a subset and can be removed
 771  *
 772  *      2. all constraints of both basic maps are either "valid" or "cut"
 773  *         and the facets corresponding to the "cut" constraints
 774  *         of one of the basic maps lies entirely inside the other basic map
 775  *              => the pair can be replaced by a basic map consisting
 776  *                 of the valid constraints in both basic maps
 777  *
 778  *      3. there is a single pair of adjacent inequalities
 779  *         (all other constraints are "valid")
 780  *              => the pair can be replaced by a basic map consisting
 781  *                 of the valid constraints in both basic maps
 782  *
 783  *      4. there is a single adjacent pair of an inequality and an equality,
 784  *         the other constraints of the basic map containing the inequality are
 785  *         "valid".  Moreover, if the inequality the basic map is relaxed
 786  *         and then turned into an equality, then resulting facet lies
 787  *         entirely inside the other basic map
 788  *              => the pair can be replaced by the basic map containing
 789  *                 the inequality, with the inequality relaxed.
 790  *
 791  *      5. there is a single adjacent pair of an inequality and an equality,
 792  *         the other constraints of the basic map containing the inequality are
 793  *         "valid".  Moreover, the facets corresponding to both
 794  *         the inequality and the equality can be wrapped around their
 795  *         ridges to include the other basic map
 796  *              => the pair can be replaced by a basic map consisting
 797  *                 of the valid constraints in both basic maps together
 798  *                 with all wrapping constraints
 799  *
 800  *      6. one of the basic maps has a single cut constraint and
 801  *         the other basic map has a constraint adjacent to this constraint.
 802  *         Moreover, the facets corresponding to both constraints
 803  *         can be wrapped around their ridges to include the other basic map
 804  *              => the pair can be replaced by a basic map consisting
 805  *                 of the valid constraints in both basic maps together
 806  *                 with all wrapping constraints
 807  *
 808  * Throughout the computation, we maintain a collection of tableaus
 809  * corresponding to the basic maps.  When the basic maps are dropped
 810  * or combined, the tableaus are modified accordingly.
 811  */
 812 static int coalesce_pair(struct isl_map *map, int i, int j,
 813         struct isl_tab **tabs)
 814 {
 815         int changed = 0;
 816         int *eq_i = NULL;
 817         int *eq_j = NULL;
 818         int *ineq_i = NULL;
 819         int *ineq_j = NULL;
 820
 821         eq_i = eq_status_in(map, i, j, tabs);
 822         if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_ERROR))
 823                 goto error;
 824         if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_SEPARATE))
 825                 goto done;
 826
 827         eq_j = eq_status_in(map, j, i, tabs);
 828         if (any(eq_j, 2 * map->p[j]->n_eq, STATUS_ERROR))
 829                 goto error;
 830         if (any(eq_j, 2 * map->p[j]->n_eq, STATUS_SEPARATE))
 831                 goto done;
 832
 833         ineq_i = ineq_status_in(map, i, j, tabs);
 834         if (any(ineq_i, map->p[i]->n_ineq, STATUS_ERROR))
 835                 goto error;
 836         if (any(ineq_i, map->p[i]->n_ineq, STATUS_SEPARATE))
 837                 goto done;
 838
 839         ineq_j = ineq_status_in(map, j, i, tabs);
 840         if (any(ineq_j, map->p[j]->n_ineq, STATUS_ERROR))
 841                 goto error;
 842         if (any(ineq_j, map->p[j]->n_ineq, STATUS_SEPARATE))
 843                 goto done;
 844
 845         if (all(eq_i, 2 * map->p[i]->n_eq, STATUS_VALID) &&
 846             all(ineq_i, map->p[i]->n_ineq, STATUS_VALID)) {
 847                 drop(map, j, tabs);
 848                 changed = 1;
 849         } else if (all(eq_j, 2 * map->p[j]->n_eq, STATUS_VALID) &&
 850                    all(ineq_j, map->p[j]->n_ineq, STATUS_VALID)) {
 851                 drop(map, i, tabs);
 852                 changed = 1;
 853         } else if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_CUT) ||
 854                    any(eq_j, 2 * map->p[j]->n_eq, STATUS_CUT)) {
 855                 /* BAD CUT */
 856         } else if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_ADJ_EQ) ||
 857                    any(eq_j, 2 * map->p[j]->n_eq, STATUS_ADJ_EQ)) {
 858                 /* ADJ EQ PAIR */
 859         } else if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_ADJ_INEQ) ||
 860                    any(eq_j, 2 * map->p[j]->n_eq, STATUS_ADJ_INEQ)) {
 861                 changed = check_adj_eq(map, i, j, tabs,
 862                                         eq_i, ineq_i, eq_j, ineq_j);
 863         } else if (any(ineq_i, map->p[i]->n_ineq, STATUS_ADJ_EQ) ||
 864                    any(ineq_j, map->p[j]->n_ineq, STATUS_ADJ_EQ)) {
 865                 /* Can't happen */
 866                 /* BAD ADJ INEQ */
 867         } else if (any(ineq_i, map->p[i]->n_ineq, STATUS_ADJ_INEQ) ||
 868                    any(ineq_j, map->p[j]->n_ineq, STATUS_ADJ_INEQ)) {
 869                 changed = check_adj_ineq(map, i, j, tabs, ineq_i, ineq_j);
 870         } else {
 871                 changed = check_facets(map, i, j, tabs, ineq_i, ineq_j);
 872                 if (!changed)
 873                         changed = check_wrap(map, i, j, tabs, ineq_i, ineq_j);
 874         }
 875
 876 done:
 877         free(eq_i);
 878         free(eq_j);
 879         free(ineq_i);
 880         free(ineq_j);
 881         return changed;
 882 error:
 883         free(eq_i);
 884         free(eq_j);
 885         free(ineq_i);
 886         free(ineq_j);
 887         return -1;
 888 }
 889
 890 static struct isl_map *coalesce(struct isl_map *map, struct isl_tab **tabs)
 891 {
 892         int i, j;
 893
 894         for (i = map->n - 2; i >= 0; --i)
 895 restart:
 896                 for (j = i + 1; j < map->n; ++j) {
 897                         int changed;
 898                         changed = coalesce_pair(map, i, j, tabs);
 899                         if (changed < 0)
 900                                 goto error;
 901                         if (changed)
 902                                 goto restart;
 903                 }
 904         return map;
 905 error:
 906         isl_map_free(map);
 907         return NULL;
 908 }
 909
 910 /* For each pair of basic maps in the map, check if the union of the two
 911  * can be represented by a single basic map.
 912  * If so, replace the pair by the single basic map and start over.
 913  */
 914 struct isl_map *isl_map_coalesce(struct isl_map *map)
 915 {
 916         int i;
 917         unsigned n;
 918         struct isl_tab **tabs = NULL;
 919
 920         if (!map)
 921                 return NULL;
 922
 923         if (map->n <= 1)
 924                 return map;
 925
 926         map = isl_map_align_divs(map);
 927
 928         tabs = isl_calloc_array(map->ctx, struct isl_tab *, map->n);
 929         if (!tabs)
 930                 goto error;
 931
 932         n = map->n;
 933         for (i = 0; i < map->n; ++i) {
 934                 tabs[i] = isl_tab_from_basic_map(map->p[i]);
 935                 if (!tabs[i])
 936                         goto error;
 937                 if (!ISL_F_ISSET(map->p[i], ISL_BASIC_MAP_NO_IMPLICIT))
 938                         tabs[i] = isl_tab_detect_implicit_equalities(tabs[i]);
 939                 if (!ISL_F_ISSET(map->p[i], ISL_BASIC_MAP_NO_REDUNDANT))
 940                         if (isl_tab_detect_redundant(tabs[i]) < 0)
 941                                 goto error;
 942         }
 943         for (i = map->n - 1; i >= 0; --i)
 944                 if (tabs[i]->empty)
 945                         drop(map, i, tabs);
 946
 947         map = coalesce(map, tabs);
 948
 949         if (map)
 950                 for (i = 0; i < map->n; ++i) {
 951                         map->p[i] = isl_basic_map_update_from_tab(map->p[i],
 952                                                                     tabs[i]);
 953                         map->p[i] = isl_basic_map_finalize(map->p[i]);
 954                         if (!map->p[i])
 955                                 goto error;
 956                         ISL_F_SET(map->p[i], ISL_BASIC_MAP_NO_IMPLICIT);
 957                         ISL_F_SET(map->p[i], ISL_BASIC_MAP_NO_REDUNDANT);
 958                 }
 959
 960         for (i = 0; i < n; ++i)
 961                 isl_tab_free(tabs[i]);
 962
 963         free(tabs);
 964
 965         return map;
 966 error:
 967         if (tabs)
 968                 for (i = 0; i < n; ++i)
 969                         isl_tab_free(tabs[i]);
 970         free(tabs);
 971         return NULL;
 972 }
 973
 974 /* For each pair of basic sets in the set, check if the union of the two
 975  * can be represented by a single basic set.
 976  * If so, replace the pair by the single basic set and start over.
 977  */
 978 struct isl_set *isl_set_coalesce(struct isl_set *set)
 979 {
 980         return (struct isl_set *)isl_map_coalesce((struct isl_map *)set);
 981 }