add isl_set_plain_cmp
[isl.git] / isl_coalesce.c
blobb99c7b53d2b8b1a791f4aeab576e2769720836dc
1 /*
2 * Copyright 2008-2009 Katholieke Universiteit Leuven
3 * Copyright 2010 INRIA Saclay
5 * Use of this software is governed by the GNU LGPLv2.1 license
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
13 #include "isl_map_private.h"
14 #include <isl/seq.h>
15 #include "isl_tab.h"
16 #include <isl_mat_private.h>
18 #define STATUS_ERROR -1
19 #define STATUS_REDUNDANT 1
20 #define STATUS_VALID 2
21 #define STATUS_SEPARATE 3
22 #define STATUS_CUT 4
23 #define STATUS_ADJ_EQ 5
24 #define STATUS_ADJ_INEQ 6
26 static int status_in(isl_int *ineq, struct isl_tab *tab)
28 enum isl_ineq_type type = isl_tab_ineq_type(tab, ineq);
29 switch (type) {
30 default:
31 case isl_ineq_error: return STATUS_ERROR;
32 case isl_ineq_redundant: return STATUS_VALID;
33 case isl_ineq_separate: return STATUS_SEPARATE;
34 case isl_ineq_cut: return STATUS_CUT;
35 case isl_ineq_adj_eq: return STATUS_ADJ_EQ;
36 case isl_ineq_adj_ineq: return STATUS_ADJ_INEQ;
40 /* Compute the position of the equalities of basic map "i"
41 * with respect to basic map "j".
42 * The resulting array has twice as many entries as the number
43 * of equalities corresponding to the two inequalties to which
44 * each equality corresponds.
46 static int *eq_status_in(struct isl_map *map, int i, int j,
47 struct isl_tab **tabs)
49 int k, l;
50 int *eq = isl_calloc_array(map->ctx, int, 2 * map->p[i]->n_eq);
51 unsigned dim;
53 dim = isl_basic_map_total_dim(map->p[i]);
54 for (k = 0; k < map->p[i]->n_eq; ++k) {
55 for (l = 0; l < 2; ++l) {
56 isl_seq_neg(map->p[i]->eq[k], map->p[i]->eq[k], 1+dim);
57 eq[2 * k + l] = status_in(map->p[i]->eq[k], tabs[j]);
58 if (eq[2 * k + l] == STATUS_ERROR)
59 goto error;
61 if (eq[2 * k] == STATUS_SEPARATE ||
62 eq[2 * k + 1] == STATUS_SEPARATE)
63 break;
66 return eq;
67 error:
68 free(eq);
69 return NULL;
72 /* Compute the position of the inequalities of basic map "i"
73 * with respect to basic map "j".
75 static int *ineq_status_in(struct isl_map *map, int i, int j,
76 struct isl_tab **tabs)
78 int k;
79 unsigned n_eq = map->p[i]->n_eq;
80 int *ineq = isl_calloc_array(map->ctx, int, map->p[i]->n_ineq);
82 for (k = 0; k < map->p[i]->n_ineq; ++k) {
83 if (isl_tab_is_redundant(tabs[i], n_eq + k)) {
84 ineq[k] = STATUS_REDUNDANT;
85 continue;
87 ineq[k] = status_in(map->p[i]->ineq[k], tabs[j]);
88 if (ineq[k] == STATUS_ERROR)
89 goto error;
90 if (ineq[k] == STATUS_SEPARATE)
91 break;
94 return ineq;
95 error:
96 free(ineq);
97 return NULL;
100 static int any(int *con, unsigned len, int status)
102 int i;
104 for (i = 0; i < len ; ++i)
105 if (con[i] == status)
106 return 1;
107 return 0;
110 static int count(int *con, unsigned len, int status)
112 int i;
113 int c = 0;
115 for (i = 0; i < len ; ++i)
116 if (con[i] == status)
117 c++;
118 return c;
121 static int all(int *con, unsigned len, int status)
123 int i;
125 for (i = 0; i < len ; ++i) {
126 if (con[i] == STATUS_REDUNDANT)
127 continue;
128 if (con[i] != status)
129 return 0;
131 return 1;
134 static void drop(struct isl_map *map, int i, struct isl_tab **tabs)
136 isl_basic_map_free(map->p[i]);
137 isl_tab_free(tabs[i]);
139 if (i != map->n - 1) {
140 map->p[i] = map->p[map->n - 1];
141 tabs[i] = tabs[map->n - 1];
143 tabs[map->n - 1] = NULL;
144 map->n--;
147 /* Replace the pair of basic maps i and j by the basic map bounded
148 * by the valid constraints in both basic maps and the constraint
149 * in extra (if not NULL).
151 static int fuse(struct isl_map *map, int i, int j,
152 struct isl_tab **tabs, int *eq_i, int *ineq_i, int *eq_j, int *ineq_j,
153 __isl_keep isl_mat *extra)
155 int k, l;
156 struct isl_basic_map *fused = NULL;
157 struct isl_tab *fused_tab = NULL;
158 unsigned total = isl_basic_map_total_dim(map->p[i]);
159 unsigned extra_rows = extra ? extra->n_row : 0;
161 fused = isl_basic_map_alloc_space(isl_space_copy(map->p[i]->dim),
162 map->p[i]->n_div,
163 map->p[i]->n_eq + map->p[j]->n_eq,
164 map->p[i]->n_ineq + map->p[j]->n_ineq + extra_rows);
165 if (!fused)
166 goto error;
168 for (k = 0; k < map->p[i]->n_eq; ++k) {
169 if (eq_i && (eq_i[2 * k] != STATUS_VALID ||
170 eq_i[2 * k + 1] != STATUS_VALID))
171 continue;
172 l = isl_basic_map_alloc_equality(fused);
173 if (l < 0)
174 goto error;
175 isl_seq_cpy(fused->eq[l], map->p[i]->eq[k], 1 + total);
178 for (k = 0; k < map->p[j]->n_eq; ++k) {
179 if (eq_j && (eq_j[2 * k] != STATUS_VALID ||
180 eq_j[2 * k + 1] != STATUS_VALID))
181 continue;
182 l = isl_basic_map_alloc_equality(fused);
183 if (l < 0)
184 goto error;
185 isl_seq_cpy(fused->eq[l], map->p[j]->eq[k], 1 + total);
188 for (k = 0; k < map->p[i]->n_ineq; ++k) {
189 if (ineq_i[k] != STATUS_VALID)
190 continue;
191 l = isl_basic_map_alloc_inequality(fused);
192 if (l < 0)
193 goto error;
194 isl_seq_cpy(fused->ineq[l], map->p[i]->ineq[k], 1 + total);
197 for (k = 0; k < map->p[j]->n_ineq; ++k) {
198 if (ineq_j[k] != STATUS_VALID)
199 continue;
200 l = isl_basic_map_alloc_inequality(fused);
201 if (l < 0)
202 goto error;
203 isl_seq_cpy(fused->ineq[l], map->p[j]->ineq[k], 1 + total);
206 for (k = 0; k < map->p[i]->n_div; ++k) {
207 int l = isl_basic_map_alloc_div(fused);
208 if (l < 0)
209 goto error;
210 isl_seq_cpy(fused->div[l], map->p[i]->div[k], 1 + 1 + total);
213 for (k = 0; k < extra_rows; ++k) {
214 l = isl_basic_map_alloc_inequality(fused);
215 if (l < 0)
216 goto error;
217 isl_seq_cpy(fused->ineq[l], extra->row[k], 1 + total);
220 fused = isl_basic_map_gauss(fused, NULL);
221 ISL_F_SET(fused, ISL_BASIC_MAP_FINAL);
222 if (ISL_F_ISSET(map->p[i], ISL_BASIC_MAP_RATIONAL) &&
223 ISL_F_ISSET(map->p[j], ISL_BASIC_MAP_RATIONAL))
224 ISL_F_SET(fused, ISL_BASIC_MAP_RATIONAL);
226 fused_tab = isl_tab_from_basic_map(fused);
227 if (isl_tab_detect_redundant(fused_tab) < 0)
228 goto error;
230 isl_basic_map_free(map->p[i]);
231 map->p[i] = fused;
232 isl_tab_free(tabs[i]);
233 tabs[i] = fused_tab;
234 drop(map, j, tabs);
236 return 1;
237 error:
238 isl_tab_free(fused_tab);
239 isl_basic_map_free(fused);
240 return -1;
243 /* Given a pair of basic maps i and j such that all constraints are either
244 * "valid" or "cut", check if the facets corresponding to the "cut"
245 * constraints of i lie entirely within basic map j.
246 * If so, replace the pair by the basic map consisting of the valid
247 * constraints in both basic maps.
249 * To see that we are not introducing any extra points, call the
250 * two basic maps A and B and the resulting map U and let x
251 * be an element of U \setminus ( A \cup B ).
252 * Then there is a pair of cut constraints c_1 and c_2 in A and B such that x
253 * violates them. Let X be the intersection of U with the opposites
254 * of these constraints. Then x \in X.
255 * The facet corresponding to c_1 contains the corresponding facet of A.
256 * This facet is entirely contained in B, so c_2 is valid on the facet.
257 * However, since it is also (part of) a facet of X, -c_2 is also valid
258 * on the facet. This means c_2 is saturated on the facet, so c_1 and
259 * c_2 must be opposites of each other, but then x could not violate
260 * both of them.
262 static int check_facets(struct isl_map *map, int i, int j,
263 struct isl_tab **tabs, int *ineq_i, int *ineq_j)
265 int k, l;
266 struct isl_tab_undo *snap;
267 unsigned n_eq = map->p[i]->n_eq;
269 snap = isl_tab_snap(tabs[i]);
271 for (k = 0; k < map->p[i]->n_ineq; ++k) {
272 if (ineq_i[k] != STATUS_CUT)
273 continue;
274 if (isl_tab_select_facet(tabs[i], n_eq + k) < 0)
275 return -1;
276 for (l = 0; l < map->p[j]->n_ineq; ++l) {
277 int stat;
278 if (ineq_j[l] != STATUS_CUT)
279 continue;
280 stat = status_in(map->p[j]->ineq[l], tabs[i]);
281 if (stat != STATUS_VALID)
282 break;
284 if (isl_tab_rollback(tabs[i], snap) < 0)
285 return -1;
286 if (l < map->p[j]->n_ineq)
287 break;
290 if (k < map->p[i]->n_ineq)
291 /* BAD CUT PAIR */
292 return 0;
293 return fuse(map, i, j, tabs, NULL, ineq_i, NULL, ineq_j, NULL);
296 /* Both basic maps have at least one inequality with and adjacent
297 * (but opposite) inequality in the other basic map.
298 * Check that there are no cut constraints and that there is only
299 * a single pair of adjacent inequalities.
300 * If so, we can replace the pair by a single basic map described
301 * by all but the pair of adjacent inequalities.
302 * Any additional points introduced lie strictly between the two
303 * adjacent hyperplanes and can therefore be integral.
305 * ____ _____
306 * / ||\ / \
307 * / || \ / \
308 * \ || \ => \ \
309 * \ || / \ /
310 * \___||_/ \_____/
312 * The test for a single pair of adjancent inequalities is important
313 * for avoiding the combination of two basic maps like the following
315 * /|
316 * / |
317 * /__|
318 * _____
319 * | |
320 * | |
321 * |___|
323 static int check_adj_ineq(struct isl_map *map, int i, int j,
324 struct isl_tab **tabs, int *ineq_i, int *ineq_j)
326 int changed = 0;
328 if (any(ineq_i, map->p[i]->n_ineq, STATUS_CUT) ||
329 any(ineq_j, map->p[j]->n_ineq, STATUS_CUT))
330 /* ADJ INEQ CUT */
332 else if (count(ineq_i, map->p[i]->n_ineq, STATUS_ADJ_INEQ) == 1 &&
333 count(ineq_j, map->p[j]->n_ineq, STATUS_ADJ_INEQ) == 1)
334 changed = fuse(map, i, j, tabs, NULL, ineq_i, NULL, ineq_j, NULL);
335 /* else ADJ INEQ TOO MANY */
337 return changed;
340 /* Check if basic map "i" contains the basic map represented
341 * by the tableau "tab".
343 static int contains(struct isl_map *map, int i, int *ineq_i,
344 struct isl_tab *tab)
346 int k, l;
347 unsigned dim;
349 dim = isl_basic_map_total_dim(map->p[i]);
350 for (k = 0; k < map->p[i]->n_eq; ++k) {
351 for (l = 0; l < 2; ++l) {
352 int stat;
353 isl_seq_neg(map->p[i]->eq[k], map->p[i]->eq[k], 1+dim);
354 stat = status_in(map->p[i]->eq[k], tab);
355 if (stat != STATUS_VALID)
356 return 0;
360 for (k = 0; k < map->p[i]->n_ineq; ++k) {
361 int stat;
362 if (ineq_i[k] == STATUS_REDUNDANT)
363 continue;
364 stat = status_in(map->p[i]->ineq[k], tab);
365 if (stat != STATUS_VALID)
366 return 0;
368 return 1;
371 /* Basic map "i" has an inequality "k" that is adjacent to some equality
372 * of basic map "j". All the other inequalities are valid for "j".
373 * Check if basic map "j" forms an extension of basic map "i".
375 * In particular, we relax constraint "k", compute the corresponding
376 * facet and check whether it is included in the other basic map.
377 * If so, we know that relaxing the constraint extends the basic
378 * map with exactly the other basic map (we already know that this
379 * other basic map is included in the extension, because there
380 * were no "cut" inequalities in "i") and we can replace the
381 * two basic maps by thie extension.
382 * ____ _____
383 * / || / |
384 * / || / |
385 * \ || => \ |
386 * \ || \ |
387 * \___|| \____|
389 static int is_extension(struct isl_map *map, int i, int j, int k,
390 struct isl_tab **tabs, int *eq_i, int *ineq_i, int *eq_j, int *ineq_j)
392 int changed = 0;
393 int super;
394 struct isl_tab_undo *snap, *snap2;
395 unsigned n_eq = map->p[i]->n_eq;
397 snap = isl_tab_snap(tabs[i]);
398 tabs[i] = isl_tab_relax(tabs[i], n_eq + k);
399 snap2 = isl_tab_snap(tabs[i]);
400 if (isl_tab_select_facet(tabs[i], n_eq + k) < 0)
401 return -1;
402 super = contains(map, j, ineq_j, tabs[i]);
403 if (super) {
404 if (isl_tab_rollback(tabs[i], snap2) < 0)
405 return -1;
406 map->p[i] = isl_basic_map_cow(map->p[i]);
407 if (!map->p[i])
408 return -1;
409 isl_int_add_ui(map->p[i]->ineq[k][0], map->p[i]->ineq[k][0], 1);
410 ISL_F_SET(map->p[i], ISL_BASIC_MAP_FINAL);
411 drop(map, j, tabs);
412 changed = 1;
413 } else
414 if (isl_tab_rollback(tabs[i], snap) < 0)
415 return -1;
417 return changed;
420 /* For each non-redundant constraint in "bmap" (as determined by "tab"),
421 * wrap the constraint around "bound" such that it includes the whole
422 * set "set" and append the resulting constraint to "wraps".
423 * "wraps" is assumed to have been pre-allocated to the appropriate size.
424 * wraps->n_row is the number of actual wrapped constraints that have
425 * been added.
426 * If any of the wrapping problems results in a constraint that is
427 * identical to "bound", then this means that "set" is unbounded in such
428 * way that no wrapping is possible. If this happens then wraps->n_row
429 * is reset to zero.
431 static int add_wraps(__isl_keep isl_mat *wraps, __isl_keep isl_basic_map *bmap,
432 struct isl_tab *tab, isl_int *bound, __isl_keep isl_set *set)
434 int l;
435 int w;
436 unsigned total = isl_basic_map_total_dim(bmap);
438 w = wraps->n_row;
440 for (l = 0; l < bmap->n_ineq; ++l) {
441 if (isl_seq_is_neg(bound, bmap->ineq[l], 1 + total))
442 continue;
443 if (isl_seq_eq(bound, bmap->ineq[l], 1 + total))
444 continue;
445 if (isl_tab_is_redundant(tab, bmap->n_eq + l))
446 continue;
448 isl_seq_cpy(wraps->row[w], bound, 1 + total);
449 if (!isl_set_wrap_facet(set, wraps->row[w], bmap->ineq[l]))
450 return -1;
451 if (isl_seq_eq(wraps->row[w], bound, 1 + total))
452 goto unbounded;
453 ++w;
455 for (l = 0; l < bmap->n_eq; ++l) {
456 if (isl_seq_is_neg(bound, bmap->eq[l], 1 + total))
457 continue;
458 if (isl_seq_eq(bound, bmap->eq[l], 1 + total))
459 continue;
461 isl_seq_cpy(wraps->row[w], bound, 1 + total);
462 isl_seq_neg(wraps->row[w + 1], bmap->eq[l], 1 + total);
463 if (!isl_set_wrap_facet(set, wraps->row[w], wraps->row[w + 1]))
464 return -1;
465 if (isl_seq_eq(wraps->row[w], bound, 1 + total))
466 goto unbounded;
467 ++w;
469 isl_seq_cpy(wraps->row[w], bound, 1 + total);
470 if (!isl_set_wrap_facet(set, wraps->row[w], bmap->eq[l]))
471 return -1;
472 if (isl_seq_eq(wraps->row[w], bound, 1 + total))
473 goto unbounded;
474 ++w;
477 wraps->n_row = w;
478 return 0;
479 unbounded:
480 wraps->n_row = 0;
481 return 0;
484 /* Check if the constraints in "wraps" from "first" until the last
485 * are all valid for the basic set represented by "tab".
486 * If not, wraps->n_row is set to zero.
488 static int check_wraps(__isl_keep isl_mat *wraps, int first,
489 struct isl_tab *tab)
491 int i;
493 for (i = first; i < wraps->n_row; ++i) {
494 enum isl_ineq_type type;
495 type = isl_tab_ineq_type(tab, wraps->row[i]);
496 if (type == isl_ineq_error)
497 return -1;
498 if (type == isl_ineq_redundant)
499 continue;
500 wraps->n_row = 0;
501 return 0;
504 return 0;
507 /* Return a set that corresponds to the non-redudant constraints
508 * (as recorded in tab) of bmap.
510 * It's important to remove the redundant constraints as some
511 * of the other constraints may have been modified after the
512 * constraints were marked redundant.
513 * In particular, a constraint may have been relaxed.
514 * Redundant constraints are ignored when a constraint is relaxed
515 * and should therefore continue to be ignored ever after.
516 * Otherwise, the relaxation might be thwarted by some of
517 * these constraints.
519 static __isl_give isl_set *set_from_updated_bmap(__isl_keep isl_basic_map *bmap,
520 struct isl_tab *tab)
522 bmap = isl_basic_map_copy(bmap);
523 bmap = isl_basic_map_cow(bmap);
524 bmap = isl_basic_map_update_from_tab(bmap, tab);
525 return isl_set_from_basic_set(isl_basic_map_underlying_set(bmap));
528 /* Given a basic set i with a constraint k that is adjacent to either the
529 * whole of basic set j or a facet of basic set j, check if we can wrap
530 * both the facet corresponding to k and the facet of j (or the whole of j)
531 * around their ridges to include the other set.
532 * If so, replace the pair of basic sets by their union.
534 * All constraints of i (except k) are assumed to be valid for j.
536 * However, the constraints of j may not be valid for i and so
537 * we have to check that the wrapping constraints for j are valid for i.
539 * In the case where j has a facet adjacent to i, tab[j] is assumed
540 * to have been restricted to this facet, so that the non-redundant
541 * constraints in tab[j] are the ridges of the facet.
542 * Note that for the purpose of wrapping, it does not matter whether
543 * we wrap the ridges of i around the whole of j or just around
544 * the facet since all the other constraints are assumed to be valid for j.
545 * In practice, we wrap to include the whole of j.
546 * ____ _____
547 * / | / \
548 * / || / |
549 * \ || => \ |
550 * \ || \ |
551 * \___|| \____|
554 static int can_wrap_in_facet(struct isl_map *map, int i, int j, int k,
555 struct isl_tab **tabs, int *eq_i, int *ineq_i, int *eq_j, int *ineq_j)
557 int changed = 0;
558 struct isl_mat *wraps = NULL;
559 struct isl_set *set_i = NULL;
560 struct isl_set *set_j = NULL;
561 struct isl_vec *bound = NULL;
562 unsigned total = isl_basic_map_total_dim(map->p[i]);
563 struct isl_tab_undo *snap;
564 int n;
566 set_i = set_from_updated_bmap(map->p[i], tabs[i]);
567 set_j = set_from_updated_bmap(map->p[j], tabs[j]);
568 wraps = isl_mat_alloc(map->ctx, 2 * (map->p[i]->n_eq + map->p[j]->n_eq) +
569 map->p[i]->n_ineq + map->p[j]->n_ineq,
570 1 + total);
571 bound = isl_vec_alloc(map->ctx, 1 + total);
572 if (!set_i || !set_j || !wraps || !bound)
573 goto error;
575 isl_seq_cpy(bound->el, map->p[i]->ineq[k], 1 + total);
576 isl_int_add_ui(bound->el[0], bound->el[0], 1);
578 isl_seq_cpy(wraps->row[0], bound->el, 1 + total);
579 wraps->n_row = 1;
581 if (add_wraps(wraps, map->p[j], tabs[j], bound->el, set_i) < 0)
582 goto error;
583 if (!wraps->n_row)
584 goto unbounded;
586 snap = isl_tab_snap(tabs[i]);
588 if (isl_tab_select_facet(tabs[i], map->p[i]->n_eq + k) < 0)
589 goto error;
590 if (isl_tab_detect_redundant(tabs[i]) < 0)
591 goto error;
593 isl_seq_neg(bound->el, map->p[i]->ineq[k], 1 + total);
595 n = wraps->n_row;
596 if (add_wraps(wraps, map->p[i], tabs[i], bound->el, set_j) < 0)
597 goto error;
599 if (isl_tab_rollback(tabs[i], snap) < 0)
600 goto error;
601 if (check_wraps(wraps, n, tabs[i]) < 0)
602 goto error;
603 if (!wraps->n_row)
604 goto unbounded;
606 changed = fuse(map, i, j, tabs, eq_i, ineq_i, eq_j, ineq_j, wraps);
608 unbounded:
609 isl_mat_free(wraps);
611 isl_set_free(set_i);
612 isl_set_free(set_j);
614 isl_vec_free(bound);
616 return changed;
617 error:
618 isl_vec_free(bound);
619 isl_mat_free(wraps);
620 isl_set_free(set_i);
621 isl_set_free(set_j);
622 return -1;
625 /* Set the is_redundant property of the "n" constraints in "cuts",
626 * except "k" to "v".
627 * This is a fairly tricky operation as it bypasses isl_tab.c.
628 * The reason we want to temporarily mark some constraints redundant
629 * is that we want to ignore them in add_wraps.
631 * Initially all cut constraints are non-redundant, but the
632 * selection of a facet right before the call to this function
633 * may have made some of them redundant.
634 * Likewise, the same constraints are marked non-redundant
635 * in the second call to this function, before they are officially
636 * made non-redundant again in the subsequent rollback.
638 static void set_is_redundant(struct isl_tab *tab, unsigned n_eq,
639 int *cuts, int n, int k, int v)
641 int l;
643 for (l = 0; l < n; ++l) {
644 if (l == k)
645 continue;
646 tab->con[n_eq + cuts[l]].is_redundant = v;
650 /* Given a pair of basic maps i and j such that j sticks out
651 * of i at n cut constraints, each time by at most one,
652 * try to compute wrapping constraints and replace the two
653 * basic maps by a single basic map.
654 * The other constraints of i are assumed to be valid for j.
656 * The facets of i corresponding to the cut constraints are
657 * wrapped around their ridges, except those ridges determined
658 * by any of the other cut constraints.
659 * The intersections of cut constraints need to be ignored
660 * as the result of wrapping one cut constraint around another
661 * would result in a constraint cutting the union.
662 * In each case, the facets are wrapped to include the union
663 * of the two basic maps.
665 * The pieces of j that lie at an offset of exactly one from
666 * one of the cut constraints of i are wrapped around their edges.
667 * Here, there is no need to ignore intersections because we
668 * are wrapping around the union of the two basic maps.
670 * If any wrapping fails, i.e., if we cannot wrap to touch
671 * the union, then we give up.
672 * Otherwise, the pair of basic maps is replaced by their union.
674 static int wrap_in_facets(struct isl_map *map, int i, int j,
675 int *cuts, int n, struct isl_tab **tabs,
676 int *eq_i, int *ineq_i, int *eq_j, int *ineq_j)
678 int changed = 0;
679 isl_mat *wraps = NULL;
680 isl_set *set = NULL;
681 isl_vec *bound = NULL;
682 unsigned total = isl_basic_map_total_dim(map->p[i]);
683 int max_wrap;
684 int k;
685 struct isl_tab_undo *snap_i, *snap_j;
687 if (isl_tab_extend_cons(tabs[j], 1) < 0)
688 goto error;
690 max_wrap = 2 * (map->p[i]->n_eq + map->p[j]->n_eq) +
691 map->p[i]->n_ineq + map->p[j]->n_ineq;
692 max_wrap *= n;
694 set = isl_set_union(set_from_updated_bmap(map->p[i], tabs[i]),
695 set_from_updated_bmap(map->p[j], tabs[j]));
696 wraps = isl_mat_alloc(map->ctx, max_wrap, 1 + total);
697 bound = isl_vec_alloc(map->ctx, 1 + total);
698 if (!set || !wraps || !bound)
699 goto error;
701 snap_i = isl_tab_snap(tabs[i]);
702 snap_j = isl_tab_snap(tabs[j]);
704 wraps->n_row = 0;
706 for (k = 0; k < n; ++k) {
707 if (isl_tab_select_facet(tabs[i], map->p[i]->n_eq + cuts[k]) < 0)
708 goto error;
709 if (isl_tab_detect_redundant(tabs[i]) < 0)
710 goto error;
711 set_is_redundant(tabs[i], map->p[i]->n_eq, cuts, n, k, 1);
713 isl_seq_neg(bound->el, map->p[i]->ineq[cuts[k]], 1 + total);
714 if (add_wraps(wraps, map->p[i], tabs[i], bound->el, set) < 0)
715 goto error;
717 set_is_redundant(tabs[i], map->p[i]->n_eq, cuts, n, k, 0);
718 if (isl_tab_rollback(tabs[i], snap_i) < 0)
719 goto error;
721 if (!wraps->n_row)
722 break;
724 isl_seq_cpy(bound->el, map->p[i]->ineq[cuts[k]], 1 + total);
725 isl_int_add_ui(bound->el[0], bound->el[0], 1);
726 if (isl_tab_add_eq(tabs[j], bound->el) < 0)
727 goto error;
728 if (isl_tab_detect_redundant(tabs[j]) < 0)
729 goto error;
731 if (!tabs[j]->empty &&
732 add_wraps(wraps, map->p[j], tabs[j], bound->el, set) < 0)
733 goto error;
735 if (isl_tab_rollback(tabs[j], snap_j) < 0)
736 goto error;
738 if (!wraps->n_row)
739 break;
742 if (k == n)
743 changed = fuse(map, i, j, tabs,
744 eq_i, ineq_i, eq_j, ineq_j, wraps);
746 isl_vec_free(bound);
747 isl_mat_free(wraps);
748 isl_set_free(set);
750 return changed;
751 error:
752 isl_vec_free(bound);
753 isl_mat_free(wraps);
754 isl_set_free(set);
755 return -1;
758 /* Given two basic sets i and j such that i has no cut equalities,
759 * check if relaxing all the cut inequalities of i by one turns
760 * them into valid constraint for j and check if we can wrap in
761 * the bits that are sticking out.
762 * If so, replace the pair by their union.
764 * We first check if all relaxed cut inequalities of i are valid for j
765 * and then try to wrap in the intersections of the relaxed cut inequalities
766 * with j.
768 * During this wrapping, we consider the points of j that lie at a distance
769 * of exactly 1 from i. In particular, we ignore the points that lie in
770 * between this lower-dimensional space and the basic map i.
771 * We can therefore only apply this to integer maps.
772 * ____ _____
773 * / ___|_ / \
774 * / | | / |
775 * \ | | => \ |
776 * \|____| \ |
777 * \___| \____/
779 * _____ ______
780 * | ____|_ | \
781 * | | | | |
782 * | | | => | |
783 * |_| | | |
784 * |_____| \______|
786 * _______
787 * | |
788 * | |\ |
789 * | | \ |
790 * | | \ |
791 * | | \|
792 * | | \
793 * | |_____\
794 * | |
795 * |_______|
797 * Wrapping can fail if the result of wrapping one of the facets
798 * around its edges does not produce any new facet constraint.
799 * In particular, this happens when we try to wrap in unbounded sets.
801 * _______________________________________________________________________
803 * | ___
804 * | | |
805 * |_| |_________________________________________________________________
806 * |___|
808 * The following is not an acceptable result of coalescing the above two
809 * sets as it includes extra integer points.
810 * _______________________________________________________________________
812 * |
813 * |
815 * \______________________________________________________________________
817 static int can_wrap_in_set(struct isl_map *map, int i, int j,
818 struct isl_tab **tabs, int *eq_i, int *ineq_i, int *eq_j, int *ineq_j)
820 int changed = 0;
821 int k, m;
822 int n;
823 int *cuts = NULL;
825 if (ISL_F_ISSET(map->p[i], ISL_BASIC_MAP_RATIONAL) ||
826 ISL_F_ISSET(map->p[j], ISL_BASIC_MAP_RATIONAL))
827 return 0;
829 n = count(ineq_i, map->p[i]->n_ineq, STATUS_CUT);
830 if (n == 0)
831 return 0;
833 cuts = isl_alloc_array(map->ctx, int, n);
834 if (!cuts)
835 return -1;
837 for (k = 0, m = 0; m < n; ++k) {
838 enum isl_ineq_type type;
840 if (ineq_i[k] != STATUS_CUT)
841 continue;
843 isl_int_add_ui(map->p[i]->ineq[k][0], map->p[i]->ineq[k][0], 1);
844 type = isl_tab_ineq_type(tabs[j], map->p[i]->ineq[k]);
845 isl_int_sub_ui(map->p[i]->ineq[k][0], map->p[i]->ineq[k][0], 1);
846 if (type == isl_ineq_error)
847 goto error;
848 if (type != isl_ineq_redundant)
849 break;
850 cuts[m] = k;
851 ++m;
854 if (m == n)
855 changed = wrap_in_facets(map, i, j, cuts, n, tabs,
856 eq_i, ineq_i, eq_j, ineq_j);
858 free(cuts);
860 return changed;
861 error:
862 free(cuts);
863 return -1;
866 /* Check if either i or j has a single cut constraint that can
867 * be used to wrap in (a facet of) the other basic set.
868 * if so, replace the pair by their union.
870 static int check_wrap(struct isl_map *map, int i, int j,
871 struct isl_tab **tabs, int *eq_i, int *ineq_i, int *eq_j, int *ineq_j)
873 int changed = 0;
875 if (!any(eq_i, 2 * map->p[i]->n_eq, STATUS_CUT))
876 changed = can_wrap_in_set(map, i, j, tabs,
877 eq_i, ineq_i, eq_j, ineq_j);
878 if (changed)
879 return changed;
881 if (!any(eq_j, 2 * map->p[j]->n_eq, STATUS_CUT))
882 changed = can_wrap_in_set(map, j, i, tabs,
883 eq_j, ineq_j, eq_i, ineq_i);
884 return changed;
887 /* At least one of the basic maps has an equality that is adjacent
888 * to inequality. Make sure that only one of the basic maps has
889 * such an equality and that the other basic map has exactly one
890 * inequality adjacent to an equality.
891 * We call the basic map that has the inequality "i" and the basic
892 * map that has the equality "j".
893 * If "i" has any "cut" (in)equality, then relaxing the inequality
894 * by one would not result in a basic map that contains the other
895 * basic map.
897 static int check_adj_eq(struct isl_map *map, int i, int j,
898 struct isl_tab **tabs, int *eq_i, int *ineq_i, int *eq_j, int *ineq_j)
900 int changed = 0;
901 int k;
903 if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_ADJ_INEQ) &&
904 any(eq_j, 2 * map->p[j]->n_eq, STATUS_ADJ_INEQ))
905 /* ADJ EQ TOO MANY */
906 return 0;
908 if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_ADJ_INEQ))
909 return check_adj_eq(map, j, i, tabs,
910 eq_j, ineq_j, eq_i, ineq_i);
912 /* j has an equality adjacent to an inequality in i */
914 if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_CUT))
915 return 0;
916 if (any(ineq_i, map->p[i]->n_ineq, STATUS_CUT))
917 /* ADJ EQ CUT */
918 return 0;
919 if (count(ineq_i, map->p[i]->n_ineq, STATUS_ADJ_EQ) != 1 ||
920 any(ineq_j, map->p[j]->n_ineq, STATUS_ADJ_EQ) ||
921 any(ineq_i, map->p[i]->n_ineq, STATUS_ADJ_INEQ) ||
922 any(ineq_j, map->p[j]->n_ineq, STATUS_ADJ_INEQ))
923 /* ADJ EQ TOO MANY */
924 return 0;
926 for (k = 0; k < map->p[i]->n_ineq ; ++k)
927 if (ineq_i[k] == STATUS_ADJ_EQ)
928 break;
930 changed = is_extension(map, i, j, k, tabs, eq_i, ineq_i, eq_j, ineq_j);
931 if (changed)
932 return changed;
934 if (count(eq_j, 2 * map->p[j]->n_eq, STATUS_ADJ_INEQ) != 1)
935 return 0;
937 changed = can_wrap_in_facet(map, i, j, k, tabs, eq_i, ineq_i, eq_j, ineq_j);
939 return changed;
942 /* The two basic maps lie on adjacent hyperplanes. In particular,
943 * basic map "i" has an equality that lies parallel to basic map "j".
944 * Check if we can wrap the facets around the parallel hyperplanes
945 * to include the other set.
947 * We perform basically the same operations as can_wrap_in_facet,
948 * except that we don't need to select a facet of one of the sets.
950 * \\ \\
951 * \\ => \\
952 * \ \|
954 * We only allow one equality of "i" to be adjacent to an equality of "j"
955 * to avoid coalescing
957 * [m, n] -> { [x, y] -> [x, 1 + y] : x >= 1 and y >= 1 and
958 * x <= 10 and y <= 10;
959 * [x, y] -> [1 + x, y] : x >= 1 and x <= 20 and
960 * y >= 5 and y <= 15 }
962 * to
964 * [m, n] -> { [x, y] -> [x2, y2] : x >= 1 and 10y2 <= 20 - x + 10y and
965 * 4y2 >= 5 + 3y and 5y2 <= 15 + 4y and
966 * y2 <= 1 + x + y - x2 and y2 >= y and
967 * y2 >= 1 + x + y - x2 }
969 static int check_eq_adj_eq(struct isl_map *map, int i, int j,
970 struct isl_tab **tabs, int *eq_i, int *ineq_i, int *eq_j, int *ineq_j)
972 int k;
973 int changed = 0;
974 struct isl_mat *wraps = NULL;
975 struct isl_set *set_i = NULL;
976 struct isl_set *set_j = NULL;
977 struct isl_vec *bound = NULL;
978 unsigned total = isl_basic_map_total_dim(map->p[i]);
980 if (count(eq_i, 2 * map->p[i]->n_eq, STATUS_ADJ_EQ) != 1)
981 return 0;
983 for (k = 0; k < 2 * map->p[i]->n_eq ; ++k)
984 if (eq_i[k] == STATUS_ADJ_EQ)
985 break;
987 set_i = set_from_updated_bmap(map->p[i], tabs[i]);
988 set_j = set_from_updated_bmap(map->p[j], tabs[j]);
989 wraps = isl_mat_alloc(map->ctx, 2 * (map->p[i]->n_eq + map->p[j]->n_eq) +
990 map->p[i]->n_ineq + map->p[j]->n_ineq,
991 1 + total);
992 bound = isl_vec_alloc(map->ctx, 1 + total);
993 if (!set_i || !set_j || !wraps || !bound)
994 goto error;
996 if (k % 2 == 0)
997 isl_seq_neg(bound->el, map->p[i]->eq[k / 2], 1 + total);
998 else
999 isl_seq_cpy(bound->el, map->p[i]->eq[k / 2], 1 + total);
1000 isl_int_add_ui(bound->el[0], bound->el[0], 1);
1002 isl_seq_cpy(wraps->row[0], bound->el, 1 + total);
1003 wraps->n_row = 1;
1005 if (add_wraps(wraps, map->p[j], tabs[j], bound->el, set_i) < 0)
1006 goto error;
1007 if (!wraps->n_row)
1008 goto unbounded;
1010 isl_int_sub_ui(bound->el[0], bound->el[0], 1);
1011 isl_seq_neg(bound->el, bound->el, 1 + total);
1013 isl_seq_cpy(wraps->row[wraps->n_row], bound->el, 1 + total);
1014 wraps->n_row++;
1016 if (add_wraps(wraps, map->p[i], tabs[i], bound->el, set_j) < 0)
1017 goto error;
1018 if (!wraps->n_row)
1019 goto unbounded;
1021 changed = fuse(map, i, j, tabs, eq_i, ineq_i, eq_j, ineq_j, wraps);
1023 if (0) {
1024 error: changed = -1;
1026 unbounded:
1028 isl_mat_free(wraps);
1029 isl_set_free(set_i);
1030 isl_set_free(set_j);
1031 isl_vec_free(bound);
1033 return changed;
1036 /* Check if the union of the given pair of basic maps
1037 * can be represented by a single basic map.
1038 * If so, replace the pair by the single basic map and return 1.
1039 * Otherwise, return 0;
1041 * We first check the effect of each constraint of one basic map
1042 * on the other basic map.
1043 * The constraint may be
1044 * redundant the constraint is redundant in its own
1045 * basic map and should be ignore and removed
1046 * in the end
1047 * valid all (integer) points of the other basic map
1048 * satisfy the constraint
1049 * separate no (integer) point of the other basic map
1050 * satisfies the constraint
1051 * cut some but not all points of the other basic map
1052 * satisfy the constraint
1053 * adj_eq the given constraint is adjacent (on the outside)
1054 * to an equality of the other basic map
1055 * adj_ineq the given constraint is adjacent (on the outside)
1056 * to an inequality of the other basic map
1058 * We consider seven cases in which we can replace the pair by a single
1059 * basic map. We ignore all "redundant" constraints.
1061 * 1. all constraints of one basic map are valid
1062 * => the other basic map is a subset and can be removed
1064 * 2. all constraints of both basic maps are either "valid" or "cut"
1065 * and the facets corresponding to the "cut" constraints
1066 * of one of the basic maps lies entirely inside the other basic map
1067 * => the pair can be replaced by a basic map consisting
1068 * of the valid constraints in both basic maps
1070 * 3. there is a single pair of adjacent inequalities
1071 * (all other constraints are "valid")
1072 * => the pair can be replaced by a basic map consisting
1073 * of the valid constraints in both basic maps
1075 * 4. there is a single adjacent pair of an inequality and an equality,
1076 * the other constraints of the basic map containing the inequality are
1077 * "valid". Moreover, if the inequality the basic map is relaxed
1078 * and then turned into an equality, then resulting facet lies
1079 * entirely inside the other basic map
1080 * => the pair can be replaced by the basic map containing
1081 * the inequality, with the inequality relaxed.
1083 * 5. there is a single adjacent pair of an inequality and an equality,
1084 * the other constraints of the basic map containing the inequality are
1085 * "valid". Moreover, the facets corresponding to both
1086 * the inequality and the equality can be wrapped around their
1087 * ridges to include the other basic map
1088 * => the pair can be replaced by a basic map consisting
1089 * of the valid constraints in both basic maps together
1090 * with all wrapping constraints
1092 * 6. one of the basic maps extends beyond the other by at most one.
1093 * Moreover, the facets corresponding to the cut constraints and
1094 * the pieces of the other basic map at offset one from these cut
1095 * constraints can be wrapped around their ridges to include
1096 * the union of the two basic maps
1097 * => the pair can be replaced by a basic map consisting
1098 * of the valid constraints in both basic maps together
1099 * with all wrapping constraints
1101 * 7. the two basic maps live in adjacent hyperplanes. In principle
1102 * such sets can always be combined through wrapping, but we impose
1103 * that there is only one such pair, to avoid overeager coalescing.
1105 * Throughout the computation, we maintain a collection of tableaus
1106 * corresponding to the basic maps. When the basic maps are dropped
1107 * or combined, the tableaus are modified accordingly.
1109 static int coalesce_pair(struct isl_map *map, int i, int j,
1110 struct isl_tab **tabs)
1112 int changed = 0;
1113 int *eq_i = NULL;
1114 int *eq_j = NULL;
1115 int *ineq_i = NULL;
1116 int *ineq_j = NULL;
1118 eq_i = eq_status_in(map, i, j, tabs);
1119 if (!eq_i)
1120 goto error;
1121 if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_ERROR))
1122 goto error;
1123 if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_SEPARATE))
1124 goto done;
1126 eq_j = eq_status_in(map, j, i, tabs);
1127 if (!eq_j)
1128 goto error;
1129 if (any(eq_j, 2 * map->p[j]->n_eq, STATUS_ERROR))
1130 goto error;
1131 if (any(eq_j, 2 * map->p[j]->n_eq, STATUS_SEPARATE))
1132 goto done;
1134 ineq_i = ineq_status_in(map, i, j, tabs);
1135 if (!ineq_i)
1136 goto error;
1137 if (any(ineq_i, map->p[i]->n_ineq, STATUS_ERROR))
1138 goto error;
1139 if (any(ineq_i, map->p[i]->n_ineq, STATUS_SEPARATE))
1140 goto done;
1142 ineq_j = ineq_status_in(map, j, i, tabs);
1143 if (!ineq_j)
1144 goto error;
1145 if (any(ineq_j, map->p[j]->n_ineq, STATUS_ERROR))
1146 goto error;
1147 if (any(ineq_j, map->p[j]->n_ineq, STATUS_SEPARATE))
1148 goto done;
1150 if (all(eq_i, 2 * map->p[i]->n_eq, STATUS_VALID) &&
1151 all(ineq_i, map->p[i]->n_ineq, STATUS_VALID)) {
1152 drop(map, j, tabs);
1153 changed = 1;
1154 } else if (all(eq_j, 2 * map->p[j]->n_eq, STATUS_VALID) &&
1155 all(ineq_j, map->p[j]->n_ineq, STATUS_VALID)) {
1156 drop(map, i, tabs);
1157 changed = 1;
1158 } else if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_ADJ_EQ)) {
1159 changed = check_eq_adj_eq(map, i, j, tabs,
1160 eq_i, ineq_i, eq_j, ineq_j);
1161 } else if (any(eq_j, 2 * map->p[j]->n_eq, STATUS_ADJ_EQ)) {
1162 changed = check_eq_adj_eq(map, j, i, tabs,
1163 eq_j, ineq_j, eq_i, ineq_i);
1164 } else if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_ADJ_INEQ) ||
1165 any(eq_j, 2 * map->p[j]->n_eq, STATUS_ADJ_INEQ)) {
1166 changed = check_adj_eq(map, i, j, tabs,
1167 eq_i, ineq_i, eq_j, ineq_j);
1168 } else if (any(ineq_i, map->p[i]->n_ineq, STATUS_ADJ_EQ) ||
1169 any(ineq_j, map->p[j]->n_ineq, STATUS_ADJ_EQ)) {
1170 /* Can't happen */
1171 /* BAD ADJ INEQ */
1172 } else if (any(ineq_i, map->p[i]->n_ineq, STATUS_ADJ_INEQ) ||
1173 any(ineq_j, map->p[j]->n_ineq, STATUS_ADJ_INEQ)) {
1174 if (!any(eq_i, 2 * map->p[i]->n_eq, STATUS_CUT) &&
1175 !any(eq_j, 2 * map->p[j]->n_eq, STATUS_CUT))
1176 changed = check_adj_ineq(map, i, j, tabs,
1177 ineq_i, ineq_j);
1178 } else {
1179 if (!any(eq_i, 2 * map->p[i]->n_eq, STATUS_CUT) &&
1180 !any(eq_j, 2 * map->p[j]->n_eq, STATUS_CUT))
1181 changed = check_facets(map, i, j, tabs, ineq_i, ineq_j);
1182 if (!changed)
1183 changed = check_wrap(map, i, j, tabs,
1184 eq_i, ineq_i, eq_j, ineq_j);
1187 done:
1188 free(eq_i);
1189 free(eq_j);
1190 free(ineq_i);
1191 free(ineq_j);
1192 return changed;
1193 error:
1194 free(eq_i);
1195 free(eq_j);
1196 free(ineq_i);
1197 free(ineq_j);
1198 return -1;
1201 static struct isl_map *coalesce(struct isl_map *map, struct isl_tab **tabs)
1203 int i, j;
1205 for (i = map->n - 2; i >= 0; --i)
1206 restart:
1207 for (j = i + 1; j < map->n; ++j) {
1208 int changed;
1209 changed = coalesce_pair(map, i, j, tabs);
1210 if (changed < 0)
1211 goto error;
1212 if (changed)
1213 goto restart;
1215 return map;
1216 error:
1217 isl_map_free(map);
1218 return NULL;
1221 /* For each pair of basic maps in the map, check if the union of the two
1222 * can be represented by a single basic map.
1223 * If so, replace the pair by the single basic map and start over.
1225 struct isl_map *isl_map_coalesce(struct isl_map *map)
1227 int i;
1228 unsigned n;
1229 struct isl_tab **tabs = NULL;
1231 if (!map)
1232 return NULL;
1234 if (map->n <= 1)
1235 return map;
1237 map = isl_map_align_divs(map);
1239 tabs = isl_calloc_array(map->ctx, struct isl_tab *, map->n);
1240 if (!tabs)
1241 goto error;
1243 n = map->n;
1244 for (i = 0; i < map->n; ++i) {
1245 tabs[i] = isl_tab_from_basic_map(map->p[i]);
1246 if (!tabs[i])
1247 goto error;
1248 if (!ISL_F_ISSET(map->p[i], ISL_BASIC_MAP_NO_IMPLICIT))
1249 if (isl_tab_detect_implicit_equalities(tabs[i]) < 0)
1250 goto error;
1251 if (!ISL_F_ISSET(map->p[i], ISL_BASIC_MAP_NO_REDUNDANT))
1252 if (isl_tab_detect_redundant(tabs[i]) < 0)
1253 goto error;
1255 for (i = map->n - 1; i >= 0; --i)
1256 if (tabs[i]->empty)
1257 drop(map, i, tabs);
1259 map = coalesce(map, tabs);
1261 if (map)
1262 for (i = 0; i < map->n; ++i) {
1263 map->p[i] = isl_basic_map_update_from_tab(map->p[i],
1264 tabs[i]);
1265 map->p[i] = isl_basic_map_finalize(map->p[i]);
1266 if (!map->p[i])
1267 goto error;
1268 ISL_F_SET(map->p[i], ISL_BASIC_MAP_NO_IMPLICIT);
1269 ISL_F_SET(map->p[i], ISL_BASIC_MAP_NO_REDUNDANT);
1272 for (i = 0; i < n; ++i)
1273 isl_tab_free(tabs[i]);
1275 free(tabs);
1277 return map;
1278 error:
1279 if (tabs)
1280 for (i = 0; i < n; ++i)
1281 isl_tab_free(tabs[i]);
1282 free(tabs);
1283 isl_map_free(map);
1284 return NULL;
1287 /* For each pair of basic sets in the set, check if the union of the two
1288 * can be represented by a single basic set.
1289 * If so, replace the pair by the single basic set and start over.
1291 struct isl_set *isl_set_coalesce(struct isl_set *set)
1293 return (struct isl_set *)isl_map_coalesce((struct isl_map *)set);