export isl_qpolynomial_add_dims
[isl.git] / isl_coalesce.c
blob8f953df1eb96b0fca1a706562b67cbbeb12a2ca7
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"
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
25 static int status_in(isl_int *ineq, struct isl_tab *tab)
27 enum isl_ineq_type type = isl_tab_ineq_type(tab, ineq);
28 switch (type) {
29 default:
30 case isl_ineq_error: return STATUS_ERROR;
31 case isl_ineq_redundant: return STATUS_VALID;
32 case isl_ineq_separate: return STATUS_SEPARATE;
33 case isl_ineq_cut: return STATUS_CUT;
34 case isl_ineq_adj_eq: return STATUS_ADJ_EQ;
35 case isl_ineq_adj_ineq: return STATUS_ADJ_INEQ;
39 /* Compute the position of the equalities of basic map "i"
40 * with respect to basic map "j".
41 * The resulting array has twice as many entries as the number
42 * of equalities corresponding to the two inequalties to which
43 * each equality corresponds.
45 static int *eq_status_in(struct isl_map *map, int i, int j,
46 struct isl_tab **tabs)
48 int k, l;
49 int *eq = isl_calloc_array(map->ctx, int, 2 * map->p[i]->n_eq);
50 unsigned dim;
52 dim = isl_basic_map_total_dim(map->p[i]);
53 for (k = 0; k < map->p[i]->n_eq; ++k) {
54 for (l = 0; l < 2; ++l) {
55 isl_seq_neg(map->p[i]->eq[k], map->p[i]->eq[k], 1+dim);
56 eq[2 * k + l] = status_in(map->p[i]->eq[k], tabs[j]);
57 if (eq[2 * k + l] == STATUS_ERROR)
58 goto error;
60 if (eq[2 * k] == STATUS_SEPARATE ||
61 eq[2 * k + 1] == STATUS_SEPARATE)
62 break;
65 return eq;
66 error:
67 free(eq);
68 return NULL;
71 /* Compute the position of the inequalities of basic map "i"
72 * with respect to basic map "j".
74 static int *ineq_status_in(struct isl_map *map, int i, int j,
75 struct isl_tab **tabs)
77 int k;
78 unsigned n_eq = map->p[i]->n_eq;
79 int *ineq = isl_calloc_array(map->ctx, int, map->p[i]->n_ineq);
81 for (k = 0; k < map->p[i]->n_ineq; ++k) {
82 if (isl_tab_is_redundant(tabs[i], n_eq + k)) {
83 ineq[k] = STATUS_REDUNDANT;
84 continue;
86 ineq[k] = status_in(map->p[i]->ineq[k], tabs[j]);
87 if (ineq[k] == STATUS_ERROR)
88 goto error;
89 if (ineq[k] == STATUS_SEPARATE)
90 break;
93 return ineq;
94 error:
95 free(ineq);
96 return NULL;
99 static int any(int *con, unsigned len, int status)
101 int i;
103 for (i = 0; i < len ; ++i)
104 if (con[i] == status)
105 return 1;
106 return 0;
109 static int count(int *con, unsigned len, int status)
111 int i;
112 int c = 0;
114 for (i = 0; i < len ; ++i)
115 if (con[i] == status)
116 c++;
117 return c;
120 static int all(int *con, unsigned len, int status)
122 int i;
124 for (i = 0; i < len ; ++i) {
125 if (con[i] == STATUS_REDUNDANT)
126 continue;
127 if (con[i] != status)
128 return 0;
130 return 1;
133 static void drop(struct isl_map *map, int i, struct isl_tab **tabs)
135 isl_basic_map_free(map->p[i]);
136 isl_tab_free(tabs[i]);
138 if (i != map->n - 1) {
139 map->p[i] = map->p[map->n - 1];
140 tabs[i] = tabs[map->n - 1];
142 tabs[map->n - 1] = NULL;
143 map->n--;
146 /* Replace the pair of basic maps i and j by the basic map bounded
147 * by the valid constraints in both basic maps and the constraint
148 * in extra (if not NULL).
150 static int fuse(struct isl_map *map, int i, int j,
151 struct isl_tab **tabs, int *eq_i, int *ineq_i, int *eq_j, int *ineq_j,
152 __isl_keep isl_mat *extra)
154 int k, l;
155 struct isl_basic_map *fused = NULL;
156 struct isl_tab *fused_tab = NULL;
157 unsigned total = isl_basic_map_total_dim(map->p[i]);
158 unsigned extra_rows = extra ? extra->n_row : 0;
160 fused = isl_basic_map_alloc_dim(isl_dim_copy(map->p[i]->dim),
161 map->p[i]->n_div,
162 map->p[i]->n_eq + map->p[j]->n_eq,
163 map->p[i]->n_ineq + map->p[j]->n_ineq + extra_rows);
164 if (!fused)
165 goto error;
167 for (k = 0; k < map->p[i]->n_eq; ++k) {
168 if (eq_i && (eq_i[2 * k] != STATUS_VALID ||
169 eq_i[2 * k + 1] != STATUS_VALID))
170 continue;
171 l = isl_basic_map_alloc_equality(fused);
172 if (l < 0)
173 goto error;
174 isl_seq_cpy(fused->eq[l], map->p[i]->eq[k], 1 + total);
177 for (k = 0; k < map->p[j]->n_eq; ++k) {
178 if (eq_j && (eq_j[2 * k] != STATUS_VALID ||
179 eq_j[2 * k + 1] != STATUS_VALID))
180 continue;
181 l = isl_basic_map_alloc_equality(fused);
182 if (l < 0)
183 goto error;
184 isl_seq_cpy(fused->eq[l], map->p[j]->eq[k], 1 + total);
187 for (k = 0; k < map->p[i]->n_ineq; ++k) {
188 if (ineq_i[k] != STATUS_VALID)
189 continue;
190 l = isl_basic_map_alloc_inequality(fused);
191 if (l < 0)
192 goto error;
193 isl_seq_cpy(fused->ineq[l], map->p[i]->ineq[k], 1 + total);
196 for (k = 0; k < map->p[j]->n_ineq; ++k) {
197 if (ineq_j[k] != STATUS_VALID)
198 continue;
199 l = isl_basic_map_alloc_inequality(fused);
200 if (l < 0)
201 goto error;
202 isl_seq_cpy(fused->ineq[l], map->p[j]->ineq[k], 1 + total);
205 for (k = 0; k < map->p[i]->n_div; ++k) {
206 int l = isl_basic_map_alloc_div(fused);
207 if (l < 0)
208 goto error;
209 isl_seq_cpy(fused->div[l], map->p[i]->div[k], 1 + 1 + total);
212 for (k = 0; k < extra_rows; ++k) {
213 l = isl_basic_map_alloc_inequality(fused);
214 if (l < 0)
215 goto error;
216 isl_seq_cpy(fused->ineq[l], extra->row[k], 1 + total);
219 fused = isl_basic_map_gauss(fused, NULL);
220 ISL_F_SET(fused, ISL_BASIC_MAP_FINAL);
221 if (ISL_F_ISSET(map->p[i], ISL_BASIC_MAP_RATIONAL) &&
222 ISL_F_ISSET(map->p[j], ISL_BASIC_MAP_RATIONAL))
223 ISL_F_SET(fused, ISL_BASIC_MAP_RATIONAL);
225 fused_tab = isl_tab_from_basic_map(fused);
226 if (isl_tab_detect_redundant(fused_tab) < 0)
227 goto error;
229 isl_basic_map_free(map->p[i]);
230 map->p[i] = fused;
231 isl_tab_free(tabs[i]);
232 tabs[i] = fused_tab;
233 drop(map, j, tabs);
235 return 1;
236 error:
237 isl_tab_free(fused_tab);
238 isl_basic_map_free(fused);
239 return -1;
242 /* Given a pair of basic maps i and j such that all constraints are either
243 * "valid" or "cut", check if the facets corresponding to the "cut"
244 * constraints of i lie entirely within basic map j.
245 * If so, replace the pair by the basic map consisting of the valid
246 * constraints in both basic maps.
248 * To see that we are not introducing any extra points, call the
249 * two basic maps A and B and the resulting map U and let x
250 * be an element of U \setminus ( A \cup B ).
251 * Then there is a pair of cut constraints c_1 and c_2 in A and B such that x
252 * violates them. Let X be the intersection of U with the opposites
253 * of these constraints. Then x \in X.
254 * The facet corresponding to c_1 contains the corresponding facet of A.
255 * This facet is entirely contained in B, so c_2 is valid on the facet.
256 * However, since it is also (part of) a facet of X, -c_2 is also valid
257 * on the facet. This means c_2 is saturated on the facet, so c_1 and
258 * c_2 must be opposites of each other, but then x could not violate
259 * both of them.
261 static int check_facets(struct isl_map *map, int i, int j,
262 struct isl_tab **tabs, int *ineq_i, int *ineq_j)
264 int k, l;
265 struct isl_tab_undo *snap;
266 unsigned n_eq = map->p[i]->n_eq;
268 snap = isl_tab_snap(tabs[i]);
270 for (k = 0; k < map->p[i]->n_ineq; ++k) {
271 if (ineq_i[k] != STATUS_CUT)
272 continue;
273 if (isl_tab_select_facet(tabs[i], n_eq + k) < 0)
274 return -1;
275 for (l = 0; l < map->p[j]->n_ineq; ++l) {
276 int stat;
277 if (ineq_j[l] != STATUS_CUT)
278 continue;
279 stat = status_in(map->p[j]->ineq[l], tabs[i]);
280 if (stat != STATUS_VALID)
281 break;
283 if (isl_tab_rollback(tabs[i], snap) < 0)
284 return -1;
285 if (l < map->p[j]->n_ineq)
286 break;
289 if (k < map->p[i]->n_ineq)
290 /* BAD CUT PAIR */
291 return 0;
292 return fuse(map, i, j, tabs, NULL, ineq_i, NULL, ineq_j, NULL);
295 /* Both basic maps have at least one inequality with and adjacent
296 * (but opposite) inequality in the other basic map.
297 * Check that there are no cut constraints and that there is only
298 * a single pair of adjacent inequalities.
299 * If so, we can replace the pair by a single basic map described
300 * by all but the pair of adjacent inequalities.
301 * Any additional points introduced lie strictly between the two
302 * adjacent hyperplanes and can therefore be integral.
304 * ____ _____
305 * / ||\ / \
306 * / || \ / \
307 * \ || \ => \ \
308 * \ || / \ /
309 * \___||_/ \_____/
311 * The test for a single pair of adjancent inequalities is important
312 * for avoiding the combination of two basic maps like the following
314 * /|
315 * / |
316 * /__|
317 * _____
318 * | |
319 * | |
320 * |___|
322 static int check_adj_ineq(struct isl_map *map, int i, int j,
323 struct isl_tab **tabs, int *ineq_i, int *ineq_j)
325 int changed = 0;
327 if (any(ineq_i, map->p[i]->n_ineq, STATUS_CUT) ||
328 any(ineq_j, map->p[j]->n_ineq, STATUS_CUT))
329 /* ADJ INEQ CUT */
331 else if (count(ineq_i, map->p[i]->n_ineq, STATUS_ADJ_INEQ) == 1 &&
332 count(ineq_j, map->p[j]->n_ineq, STATUS_ADJ_INEQ) == 1)
333 changed = fuse(map, i, j, tabs, NULL, ineq_i, NULL, ineq_j, NULL);
334 /* else ADJ INEQ TOO MANY */
336 return changed;
339 /* Check if basic map "i" contains the basic map represented
340 * by the tableau "tab".
342 static int contains(struct isl_map *map, int i, int *ineq_i,
343 struct isl_tab *tab)
345 int k, l;
346 unsigned dim;
348 dim = isl_basic_map_total_dim(map->p[i]);
349 for (k = 0; k < map->p[i]->n_eq; ++k) {
350 for (l = 0; l < 2; ++l) {
351 int stat;
352 isl_seq_neg(map->p[i]->eq[k], map->p[i]->eq[k], 1+dim);
353 stat = status_in(map->p[i]->eq[k], tab);
354 if (stat != STATUS_VALID)
355 return 0;
359 for (k = 0; k < map->p[i]->n_ineq; ++k) {
360 int stat;
361 if (ineq_i[k] == STATUS_REDUNDANT)
362 continue;
363 stat = status_in(map->p[i]->ineq[k], tab);
364 if (stat != STATUS_VALID)
365 return 0;
367 return 1;
370 /* Basic map "i" has an inequality "k" that is adjacent to some equality
371 * of basic map "j". All the other inequalities are valid for "j".
372 * Check if basic map "j" forms an extension of basic map "i".
374 * In particular, we relax constraint "k", compute the corresponding
375 * facet and check whether it is included in the other basic map.
376 * If so, we know that relaxing the constraint extends the basic
377 * map with exactly the other basic map (we already know that this
378 * other basic map is included in the extension, because there
379 * were no "cut" inequalities in "i") and we can replace the
380 * two basic maps by thie extension.
381 * ____ _____
382 * / || / |
383 * / || / |
384 * \ || => \ |
385 * \ || \ |
386 * \___|| \____|
388 static int is_extension(struct isl_map *map, int i, int j, int k,
389 struct isl_tab **tabs, int *eq_i, int *ineq_i, int *eq_j, int *ineq_j)
391 int changed = 0;
392 int super;
393 struct isl_tab_undo *snap, *snap2;
394 unsigned n_eq = map->p[i]->n_eq;
396 snap = isl_tab_snap(tabs[i]);
397 tabs[i] = isl_tab_relax(tabs[i], n_eq + k);
398 snap2 = isl_tab_snap(tabs[i]);
399 if (isl_tab_select_facet(tabs[i], n_eq + k) < 0)
400 return -1;
401 super = contains(map, j, ineq_j, tabs[i]);
402 if (super) {
403 if (isl_tab_rollback(tabs[i], snap2) < 0)
404 return -1;
405 map->p[i] = isl_basic_map_cow(map->p[i]);
406 if (!map->p[i])
407 return -1;
408 isl_int_add_ui(map->p[i]->ineq[k][0], map->p[i]->ineq[k][0], 1);
409 ISL_F_SET(map->p[i], ISL_BASIC_MAP_FINAL);
410 drop(map, j, tabs);
411 changed = 1;
412 } else
413 if (isl_tab_rollback(tabs[i], snap) < 0)
414 return -1;
416 return changed;
419 /* For each non-redundant constraint in "bmap" (as determined by "tab"),
420 * wrap the constraint around "bound" such that it includes the whole
421 * set "set" and append the resulting constraint to "wraps".
422 * "wraps" is assumed to have been pre-allocated to the appropriate size.
423 * wraps->n_row is the number of actual wrapped constraints that have
424 * been added.
425 * If any of the wrapping problems results in a constraint that is
426 * identical to "bound", then this means that "set" is unbounded in such
427 * way that no wrapping is possible. If this happens then wraps->n_row
428 * is reset to zero.
430 static int add_wraps(__isl_keep isl_mat *wraps, __isl_keep isl_basic_map *bmap,
431 struct isl_tab *tab, isl_int *bound, __isl_keep isl_set *set)
433 int l;
434 int w;
435 unsigned total = isl_basic_map_total_dim(bmap);
437 w = wraps->n_row;
439 for (l = 0; l < bmap->n_ineq; ++l) {
440 if (isl_seq_is_neg(bound, bmap->ineq[l], 1 + total))
441 continue;
442 if (isl_seq_eq(bound, bmap->ineq[l], 1 + total))
443 continue;
444 if (isl_tab_is_redundant(tab, bmap->n_eq + l))
445 continue;
447 isl_seq_cpy(wraps->row[w], bound, 1 + total);
448 if (!isl_set_wrap_facet(set, wraps->row[w], bmap->ineq[l]))
449 return -1;
450 if (isl_seq_eq(wraps->row[w], bound, 1 + total))
451 goto unbounded;
452 ++w;
454 for (l = 0; l < bmap->n_eq; ++l) {
455 if (isl_seq_is_neg(bound, bmap->eq[l], 1 + total))
456 continue;
457 if (isl_seq_eq(bound, bmap->eq[l], 1 + total))
458 continue;
460 isl_seq_cpy(wraps->row[w], bound, 1 + total);
461 isl_seq_neg(wraps->row[w + 1], bmap->eq[l], 1 + total);
462 if (!isl_set_wrap_facet(set, wraps->row[w], wraps->row[w + 1]))
463 return -1;
464 if (isl_seq_eq(wraps->row[w], bound, 1 + total))
465 goto unbounded;
466 ++w;
468 isl_seq_cpy(wraps->row[w], bound, 1 + total);
469 if (!isl_set_wrap_facet(set, wraps->row[w], bmap->eq[l]))
470 return -1;
471 if (isl_seq_eq(wraps->row[w], bound, 1 + total))
472 goto unbounded;
473 ++w;
476 wraps->n_row = w;
477 return 0;
478 unbounded:
479 wraps->n_row = 0;
480 return 0;
483 /* Check if the constraints in "wraps" from "first" until the last
484 * are all valid for the basic set represented by "tab".
485 * If not, wraps->n_row is set to zero.
487 static int check_wraps(__isl_keep isl_mat *wraps, int first,
488 struct isl_tab *tab)
490 int i;
492 for (i = first; i < wraps->n_row; ++i) {
493 enum isl_ineq_type type;
494 type = isl_tab_ineq_type(tab, wraps->row[i]);
495 if (type == isl_ineq_error)
496 return -1;
497 if (type == isl_ineq_redundant)
498 continue;
499 wraps->n_row = 0;
500 return 0;
503 return 0;
506 /* Return a set that corresponds to the non-redudant constraints
507 * (as recorded in tab) of bmap.
509 * It's important to remove the redundant constraints as some
510 * of the other constraints may have been modified after the
511 * constraints were marked redundant.
512 * In particular, a constraint may have been relaxed.
513 * Redundant constraints are ignored when a constraint is relaxed
514 * and should therefore continue to be ignored ever after.
515 * Otherwise, the relaxation might be thwarted by some of
516 * these constraints.
518 static __isl_give isl_set *set_from_updated_bmap(__isl_keep isl_basic_map *bmap,
519 struct isl_tab *tab)
521 bmap = isl_basic_map_copy(bmap);
522 bmap = isl_basic_map_cow(bmap);
523 bmap = isl_basic_map_update_from_tab(bmap, tab);
524 return isl_set_from_basic_set(isl_basic_map_underlying_set(bmap));
527 /* Given a basic set i with a constraint k that is adjacent to either the
528 * whole of basic set j or a facet of basic set j, check if we can wrap
529 * both the facet corresponding to k and the facet of j (or the whole of j)
530 * around their ridges to include the other set.
531 * If so, replace the pair of basic sets by their union.
533 * All constraints of i (except k) are assumed to be valid for j.
535 * However, the constraints of j may not be valid for i and so
536 * we have to check that the wrapping constraints for j are valid for i.
538 * In the case where j has a facet adjacent to i, tab[j] is assumed
539 * to have been restricted to this facet, so that the non-redundant
540 * constraints in tab[j] are the ridges of the facet.
541 * Note that for the purpose of wrapping, it does not matter whether
542 * we wrap the ridges of i around the whole of j or just around
543 * the facet since all the other constraints are assumed to be valid for j.
544 * In practice, we wrap to include the whole of j.
545 * ____ _____
546 * / | / \
547 * / || / |
548 * \ || => \ |
549 * \ || \ |
550 * \___|| \____|
553 static int can_wrap_in_facet(struct isl_map *map, int i, int j, int k,
554 struct isl_tab **tabs, int *eq_i, int *ineq_i, int *eq_j, int *ineq_j)
556 int changed = 0;
557 struct isl_mat *wraps = NULL;
558 struct isl_set *set_i = NULL;
559 struct isl_set *set_j = NULL;
560 struct isl_vec *bound = NULL;
561 unsigned total = isl_basic_map_total_dim(map->p[i]);
562 struct isl_tab_undo *snap;
563 int n;
565 set_i = set_from_updated_bmap(map->p[i], tabs[i]);
566 set_j = set_from_updated_bmap(map->p[j], tabs[j]);
567 wraps = isl_mat_alloc(map->ctx, 2 * (map->p[i]->n_eq + map->p[j]->n_eq) +
568 map->p[i]->n_ineq + map->p[j]->n_ineq,
569 1 + total);
570 bound = isl_vec_alloc(map->ctx, 1 + total);
571 if (!set_i || !set_j || !wraps || !bound)
572 goto error;
574 isl_seq_cpy(bound->el, map->p[i]->ineq[k], 1 + total);
575 isl_int_add_ui(bound->el[0], bound->el[0], 1);
577 isl_seq_cpy(wraps->row[0], bound->el, 1 + total);
578 wraps->n_row = 1;
580 if (add_wraps(wraps, map->p[j], tabs[j], bound->el, set_i) < 0)
581 goto error;
582 if (!wraps->n_row)
583 goto unbounded;
585 snap = isl_tab_snap(tabs[i]);
587 if (isl_tab_select_facet(tabs[i], map->p[i]->n_eq + k) < 0)
588 goto error;
589 if (isl_tab_detect_redundant(tabs[i]) < 0)
590 goto error;
592 isl_seq_neg(bound->el, map->p[i]->ineq[k], 1 + total);
594 n = wraps->n_row;
595 if (add_wraps(wraps, map->p[i], tabs[i], bound->el, set_j) < 0)
596 goto error;
598 if (isl_tab_rollback(tabs[i], snap) < 0)
599 goto error;
600 if (check_wraps(wraps, n, tabs[i]) < 0)
601 goto error;
602 if (!wraps->n_row)
603 goto unbounded;
605 changed = fuse(map, i, j, tabs, eq_i, ineq_i, eq_j, ineq_j, wraps);
607 unbounded:
608 isl_mat_free(wraps);
610 isl_set_free(set_i);
611 isl_set_free(set_j);
613 isl_vec_free(bound);
615 return changed;
616 error:
617 isl_vec_free(bound);
618 isl_mat_free(wraps);
619 isl_set_free(set_i);
620 isl_set_free(set_j);
621 return -1;
624 /* Set the is_redundant property of the "n" constraints in "cuts",
625 * except "k" to "v".
626 * This is a fairly tricky operation as it bypasses isl_tab.c.
627 * The reason we want to temporarily mark some constraints redundant
628 * is that we want to ignore them in add_wraps.
630 * Initially all cut constraints are non-redundant, but the
631 * selection of a facet right before the call to this function
632 * may have made some of them redundant.
633 * Likewise, the same constraints are marked non-redundant
634 * in the second call to this function, before they are officially
635 * made non-redundant again in the subsequent rollback.
637 static void set_is_redundant(struct isl_tab *tab, unsigned n_eq,
638 int *cuts, int n, int k, int v)
640 int l;
642 for (l = 0; l < n; ++l) {
643 if (l == k)
644 continue;
645 tab->con[n_eq + cuts[l]].is_redundant = v;
649 /* Given a pair of basic maps i and j such that j stick out
650 * of i at n cut constraints, each time by at most one,
651 * try to compute wrapping constraints and replace the two
652 * basic maps by a single basic map.
653 * The other constraints of i are assumed to be valid for j.
655 * The facets of i corresponding to the cut constraints are
656 * wrapped around their ridges, except those ridges determined
657 * by any of the other cut constraints.
658 * The intersections of cut constraints need to be ignored
659 * as the result of wrapping on cur constraint around another
660 * would result in a constraint cutting the union.
661 * In each case, the facets are wrapped to include the union
662 * of the two basic maps.
664 * The pieces of j that lie at an offset of exactly one from
665 * one of the cut constraints of i are wrapped around their edges.
666 * Here, there is no need to ignore intersections because we
667 * are wrapping around the union of the two basic maps.
669 * If any wrapping fails, i.e., if we cannot wrap to touch
670 * the union, then we give up.
671 * Otherwise, the pair of basic maps is replaced by their union.
673 static int wrap_in_facets(struct isl_map *map, int i, int j,
674 int *cuts, int n, struct isl_tab **tabs,
675 int *eq_i, int *ineq_i, int *eq_j, int *ineq_j)
677 int changed = 0;
678 isl_mat *wraps = NULL;
679 isl_set *set = NULL;
680 isl_vec *bound = NULL;
681 unsigned total = isl_basic_map_total_dim(map->p[i]);
682 int max_wrap;
683 int k;
684 struct isl_tab_undo *snap_i, *snap_j;
686 if (isl_tab_extend_cons(tabs[j], 1) < 0)
687 goto error;
689 max_wrap = 2 * (map->p[i]->n_eq + map->p[j]->n_eq) +
690 map->p[i]->n_ineq + map->p[j]->n_ineq;
691 max_wrap *= n;
693 set = isl_set_union(set_from_updated_bmap(map->p[i], tabs[i]),
694 set_from_updated_bmap(map->p[j], tabs[j]));
695 wraps = isl_mat_alloc(map->ctx, max_wrap, 1 + total);
696 bound = isl_vec_alloc(map->ctx, 1 + total);
697 if (!set || !wraps || !bound)
698 goto error;
700 snap_i = isl_tab_snap(tabs[i]);
701 snap_j = isl_tab_snap(tabs[j]);
703 wraps->n_row = 0;
705 for (k = 0; k < n; ++k) {
706 if (isl_tab_select_facet(tabs[i], map->p[i]->n_eq + cuts[k]) < 0)
707 goto error;
708 if (isl_tab_detect_redundant(tabs[i]) < 0)
709 goto error;
710 set_is_redundant(tabs[i], map->p[i]->n_eq, cuts, n, k, 1);
712 isl_seq_neg(bound->el, map->p[i]->ineq[cuts[k]], 1 + total);
713 if (add_wraps(wraps, map->p[i], tabs[i], bound->el, set) < 0)
714 goto error;
716 set_is_redundant(tabs[i], map->p[i]->n_eq, cuts, n, k, 0);
717 if (isl_tab_rollback(tabs[i], snap_i) < 0)
718 goto error;
720 if (!wraps->n_row)
721 break;
723 isl_seq_cpy(bound->el, map->p[i]->ineq[cuts[k]], 1 + total);
724 isl_int_add_ui(bound->el[0], bound->el[0], 1);
725 tabs[j] = isl_tab_add_eq(tabs[j], bound->el);
726 if (isl_tab_detect_redundant(tabs[j]) < 0)
727 goto error;
729 if (!tabs[j]->empty &&
730 add_wraps(wraps, map->p[j], tabs[j], bound->el, set) < 0)
731 goto error;
733 if (isl_tab_rollback(tabs[j], snap_j) < 0)
734 goto error;
736 if (!wraps->n_row)
737 break;
740 if (k == n)
741 changed = fuse(map, i, j, tabs,
742 eq_i, ineq_i, eq_j, ineq_j, wraps);
744 isl_vec_free(bound);
745 isl_mat_free(wraps);
746 isl_set_free(set);
748 return changed;
749 error:
750 isl_vec_free(bound);
751 isl_mat_free(wraps);
752 isl_set_free(set);
753 return -1;
756 /* Given two basic sets i and j such that i has not cut equalities,
757 * check if relaxing all the cut inequalities of i by one turns
758 * them into valid constraint for j and check if we can wrap in
759 * the bits that are sticking out.
760 * If so, replace the pair by their union.
762 * We first check if all relaxed cut inequalities of i are valid for j
763 * and then try to wrap in the intersections of the relaxed cut inequalities
764 * with j.
766 * During this wrapping, we consider the points of j that lie at a distance
767 * of exactly 1 from i. In particular, we ignore the points that lie in
768 * between this lower-dimensional space and the basic map i.
769 * We can therefore only apply this to integer maps.
770 * ____ _____
771 * / ___|_ / \
772 * / | | / |
773 * \ | | => \ |
774 * \|____| \ |
775 * \___| \____/
777 * _____ ______
778 * | ____|_ | \
779 * | | | | |
780 * | | | => | |
781 * |_| | | |
782 * |_____| \______|
784 * _______
785 * | |
786 * | |\ |
787 * | | \ |
788 * | | \ |
789 * | | \|
790 * | | \
791 * | |_____\
792 * | |
793 * |_______|
795 * Wrapping can fail if the result of wrapping one of the facets
796 * around its edges does not produce any new facet constraint.
797 * In particular, this happens when we try to wrap in unbounded sets.
799 * _______________________________________________________________________
801 * | ___
802 * | | |
803 * |_| |_________________________________________________________________
804 * |___|
806 * The following is not an acceptable result of coalescing the above two
807 * sets as it includes extra integer points.
808 * _______________________________________________________________________
810 * |
811 * |
813 * \______________________________________________________________________
815 static int can_wrap_in_set(struct isl_map *map, int i, int j,
816 struct isl_tab **tabs, int *eq_i, int *ineq_i, int *eq_j, int *ineq_j)
818 int changed = 0;
819 int k, m;
820 int n;
821 int *cuts = NULL;
823 if (ISL_F_ISSET(map->p[i], ISL_BASIC_MAP_RATIONAL) ||
824 ISL_F_ISSET(map->p[j], ISL_BASIC_MAP_RATIONAL))
825 return 0;
827 n = count(ineq_i, map->p[i]->n_ineq, STATUS_CUT);
828 if (n == 0)
829 return 0;
831 cuts = isl_alloc_array(map->ctx, int, n);
832 if (!cuts)
833 return -1;
835 for (k = 0, m = 0; m < n; ++k) {
836 enum isl_ineq_type type;
838 if (ineq_i[k] != STATUS_CUT)
839 continue;
841 isl_int_add_ui(map->p[i]->ineq[k][0], map->p[i]->ineq[k][0], 1);
842 type = isl_tab_ineq_type(tabs[j], map->p[i]->ineq[k]);
843 isl_int_sub_ui(map->p[i]->ineq[k][0], map->p[i]->ineq[k][0], 1);
844 if (type == isl_ineq_error)
845 goto error;
846 if (type != isl_ineq_redundant)
847 break;
848 cuts[m] = k;
849 ++m;
852 if (m == n)
853 changed = wrap_in_facets(map, i, j, cuts, n, tabs,
854 eq_i, ineq_i, eq_j, ineq_j);
856 free(cuts);
858 return changed;
859 error:
860 free(cuts);
861 return -1;
864 /* Check if either i or j has a single cut constraint that can
865 * be used to wrap in (a facet of) the other basic set.
866 * if so, replace the pair by their union.
868 static int check_wrap(struct isl_map *map, int i, int j,
869 struct isl_tab **tabs, int *eq_i, int *ineq_i, int *eq_j, int *ineq_j)
871 int changed = 0;
873 if (!any(eq_i, 2 * map->p[i]->n_eq, STATUS_CUT))
874 changed = can_wrap_in_set(map, i, j, tabs,
875 eq_i, ineq_i, eq_j, ineq_j);
876 if (changed)
877 return changed;
879 if (!any(eq_j, 2 * map->p[j]->n_eq, STATUS_CUT))
880 changed = can_wrap_in_set(map, j, i, tabs,
881 eq_j, ineq_j, eq_i, ineq_i);
882 return changed;
885 /* At least one of the basic maps has an equality that is adjacent
886 * to inequality. Make sure that only one of the basic maps has
887 * such an equality and that the other basic map has exactly one
888 * inequality adjacent to an equality.
889 * We call the basic map that has the inequality "i" and the basic
890 * map that has the equality "j".
891 * If "i" has any "cut" (in)equality, then relaxing the inequality
892 * by one would not result in a basic map that contains the other
893 * basic map.
895 static int check_adj_eq(struct isl_map *map, int i, int j,
896 struct isl_tab **tabs, int *eq_i, int *ineq_i, int *eq_j, int *ineq_j)
898 int changed = 0;
899 int k;
901 if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_ADJ_INEQ) &&
902 any(eq_j, 2 * map->p[j]->n_eq, STATUS_ADJ_INEQ))
903 /* ADJ EQ TOO MANY */
904 return 0;
906 if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_ADJ_INEQ))
907 return check_adj_eq(map, j, i, tabs,
908 eq_j, ineq_j, eq_i, ineq_i);
910 /* j has an equality adjacent to an inequality in i */
912 if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_CUT))
913 return 0;
914 if (any(ineq_i, map->p[i]->n_ineq, STATUS_CUT))
915 /* ADJ EQ CUT */
916 return 0;
917 if (count(eq_j, 2 * map->p[j]->n_eq, STATUS_ADJ_INEQ) != 1 ||
918 count(ineq_i, map->p[i]->n_ineq, STATUS_ADJ_EQ) != 1 ||
919 any(ineq_j, map->p[j]->n_ineq, STATUS_ADJ_EQ) ||
920 any(ineq_i, map->p[i]->n_ineq, STATUS_ADJ_INEQ) ||
921 any(ineq_j, map->p[j]->n_ineq, STATUS_ADJ_INEQ))
922 /* ADJ EQ TOO MANY */
923 return 0;
925 for (k = 0; k < map->p[i]->n_ineq ; ++k)
926 if (ineq_i[k] == STATUS_ADJ_EQ)
927 break;
929 changed = is_extension(map, i, j, k, tabs, eq_i, ineq_i, eq_j, ineq_j);
930 if (changed)
931 return changed;
933 changed = can_wrap_in_facet(map, i, j, k, tabs, eq_i, ineq_i, eq_j, ineq_j);
935 return changed;
938 /* Check if the union of the given pair of basic maps
939 * can be represented by a single basic map.
940 * If so, replace the pair by the single basic map and return 1.
941 * Otherwise, return 0;
943 * We first check the effect of each constraint of one basic map
944 * on the other basic map.
945 * The constraint may be
946 * redundant the constraint is redundant in its own
947 * basic map and should be ignore and removed
948 * in the end
949 * valid all (integer) points of the other basic map
950 * satisfy the constraint
951 * separate no (integer) point of the other basic map
952 * satisfies the constraint
953 * cut some but not all points of the other basic map
954 * satisfy the constraint
955 * adj_eq the given constraint is adjacent (on the outside)
956 * to an equality of the other basic map
957 * adj_ineq the given constraint is adjacent (on the outside)
958 * to an inequality of the other basic map
960 * We consider six cases in which we can replace the pair by a single
961 * basic map. We ignore all "redundant" constraints.
963 * 1. all constraints of one basic map are valid
964 * => the other basic map is a subset and can be removed
966 * 2. all constraints of both basic maps are either "valid" or "cut"
967 * and the facets corresponding to the "cut" constraints
968 * of one of the basic maps lies entirely inside the other basic map
969 * => the pair can be replaced by a basic map consisting
970 * of the valid constraints in both basic maps
972 * 3. there is a single pair of adjacent inequalities
973 * (all other constraints are "valid")
974 * => the pair can be replaced by a basic map consisting
975 * of the valid constraints in both basic maps
977 * 4. there is a single adjacent pair of an inequality and an equality,
978 * the other constraints of the basic map containing the inequality are
979 * "valid". Moreover, if the inequality the basic map is relaxed
980 * and then turned into an equality, then resulting facet lies
981 * entirely inside the other basic map
982 * => the pair can be replaced by the basic map containing
983 * the inequality, with the inequality relaxed.
985 * 5. there is a single adjacent pair of an inequality and an equality,
986 * the other constraints of the basic map containing the inequality are
987 * "valid". Moreover, the facets corresponding to both
988 * the inequality and the equality can be wrapped around their
989 * ridges to include the other basic map
990 * => the pair can be replaced by a basic map consisting
991 * of the valid constraints in both basic maps together
992 * with all wrapping constraints
994 * 6. one of the basic maps extends beyond the other by at most one.
995 * Moreover, the facets corresponding to the cut constraints and
996 * the pieces of the other basic map at offset one from these cut
997 * constraints can be wrapped around their ridges to include
998 * the unione of the two basic maps
999 * => the pair can be replaced by a basic map consisting
1000 * of the valid constraints in both basic maps together
1001 * with all wrapping constraints
1003 * Throughout the computation, we maintain a collection of tableaus
1004 * corresponding to the basic maps. When the basic maps are dropped
1005 * or combined, the tableaus are modified accordingly.
1007 static int coalesce_pair(struct isl_map *map, int i, int j,
1008 struct isl_tab **tabs)
1010 int changed = 0;
1011 int *eq_i = NULL;
1012 int *eq_j = NULL;
1013 int *ineq_i = NULL;
1014 int *ineq_j = NULL;
1016 eq_i = eq_status_in(map, i, j, tabs);
1017 if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_ERROR))
1018 goto error;
1019 if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_SEPARATE))
1020 goto done;
1022 eq_j = eq_status_in(map, j, i, tabs);
1023 if (any(eq_j, 2 * map->p[j]->n_eq, STATUS_ERROR))
1024 goto error;
1025 if (any(eq_j, 2 * map->p[j]->n_eq, STATUS_SEPARATE))
1026 goto done;
1028 ineq_i = ineq_status_in(map, i, j, tabs);
1029 if (any(ineq_i, map->p[i]->n_ineq, STATUS_ERROR))
1030 goto error;
1031 if (any(ineq_i, map->p[i]->n_ineq, STATUS_SEPARATE))
1032 goto done;
1034 ineq_j = ineq_status_in(map, j, i, tabs);
1035 if (any(ineq_j, map->p[j]->n_ineq, STATUS_ERROR))
1036 goto error;
1037 if (any(ineq_j, map->p[j]->n_ineq, STATUS_SEPARATE))
1038 goto done;
1040 if (all(eq_i, 2 * map->p[i]->n_eq, STATUS_VALID) &&
1041 all(ineq_i, map->p[i]->n_ineq, STATUS_VALID)) {
1042 drop(map, j, tabs);
1043 changed = 1;
1044 } else if (all(eq_j, 2 * map->p[j]->n_eq, STATUS_VALID) &&
1045 all(ineq_j, map->p[j]->n_ineq, STATUS_VALID)) {
1046 drop(map, i, tabs);
1047 changed = 1;
1048 } else if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_ADJ_EQ) ||
1049 any(eq_j, 2 * map->p[j]->n_eq, STATUS_ADJ_EQ)) {
1050 /* ADJ EQ PAIR */
1051 } else if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_ADJ_INEQ) ||
1052 any(eq_j, 2 * map->p[j]->n_eq, STATUS_ADJ_INEQ)) {
1053 changed = check_adj_eq(map, i, j, tabs,
1054 eq_i, ineq_i, eq_j, ineq_j);
1055 } else if (any(ineq_i, map->p[i]->n_ineq, STATUS_ADJ_EQ) ||
1056 any(ineq_j, map->p[j]->n_ineq, STATUS_ADJ_EQ)) {
1057 /* Can't happen */
1058 /* BAD ADJ INEQ */
1059 } else if (any(ineq_i, map->p[i]->n_ineq, STATUS_ADJ_INEQ) ||
1060 any(ineq_j, map->p[j]->n_ineq, STATUS_ADJ_INEQ)) {
1061 if (!any(eq_i, 2 * map->p[i]->n_eq, STATUS_CUT) &&
1062 !any(eq_j, 2 * map->p[j]->n_eq, STATUS_CUT))
1063 changed = check_adj_ineq(map, i, j, tabs,
1064 ineq_i, ineq_j);
1065 } else {
1066 if (!any(eq_i, 2 * map->p[i]->n_eq, STATUS_CUT) &&
1067 !any(eq_j, 2 * map->p[j]->n_eq, STATUS_CUT))
1068 changed = check_facets(map, i, j, tabs, ineq_i, ineq_j);
1069 if (!changed)
1070 changed = check_wrap(map, i, j, tabs,
1071 eq_i, ineq_i, eq_j, ineq_j);
1074 done:
1075 free(eq_i);
1076 free(eq_j);
1077 free(ineq_i);
1078 free(ineq_j);
1079 return changed;
1080 error:
1081 free(eq_i);
1082 free(eq_j);
1083 free(ineq_i);
1084 free(ineq_j);
1085 return -1;
1088 static struct isl_map *coalesce(struct isl_map *map, struct isl_tab **tabs)
1090 int i, j;
1092 for (i = map->n - 2; i >= 0; --i)
1093 restart:
1094 for (j = i + 1; j < map->n; ++j) {
1095 int changed;
1096 changed = coalesce_pair(map, i, j, tabs);
1097 if (changed < 0)
1098 goto error;
1099 if (changed)
1100 goto restart;
1102 return map;
1103 error:
1104 isl_map_free(map);
1105 return NULL;
1108 /* For each pair of basic maps in the map, check if the union of the two
1109 * can be represented by a single basic map.
1110 * If so, replace the pair by the single basic map and start over.
1112 struct isl_map *isl_map_coalesce(struct isl_map *map)
1114 int i;
1115 unsigned n;
1116 struct isl_tab **tabs = NULL;
1118 if (!map)
1119 return NULL;
1121 if (map->n <= 1)
1122 return map;
1124 map = isl_map_align_divs(map);
1126 tabs = isl_calloc_array(map->ctx, struct isl_tab *, map->n);
1127 if (!tabs)
1128 goto error;
1130 n = map->n;
1131 for (i = 0; i < map->n; ++i) {
1132 tabs[i] = isl_tab_from_basic_map(map->p[i]);
1133 if (!tabs[i])
1134 goto error;
1135 if (!ISL_F_ISSET(map->p[i], ISL_BASIC_MAP_NO_IMPLICIT))
1136 if (isl_tab_detect_implicit_equalities(tabs[i]) < 0)
1137 goto error;
1138 if (!ISL_F_ISSET(map->p[i], ISL_BASIC_MAP_NO_REDUNDANT))
1139 if (isl_tab_detect_redundant(tabs[i]) < 0)
1140 goto error;
1142 for (i = map->n - 1; i >= 0; --i)
1143 if (tabs[i]->empty)
1144 drop(map, i, tabs);
1146 map = coalesce(map, tabs);
1148 if (map)
1149 for (i = 0; i < map->n; ++i) {
1150 map->p[i] = isl_basic_map_update_from_tab(map->p[i],
1151 tabs[i]);
1152 map->p[i] = isl_basic_map_finalize(map->p[i]);
1153 if (!map->p[i])
1154 goto error;
1155 ISL_F_SET(map->p[i], ISL_BASIC_MAP_NO_IMPLICIT);
1156 ISL_F_SET(map->p[i], ISL_BASIC_MAP_NO_REDUNDANT);
1159 for (i = 0; i < n; ++i)
1160 isl_tab_free(tabs[i]);
1162 free(tabs);
1164 return map;
1165 error:
1166 if (tabs)
1167 for (i = 0; i < n; ++i)
1168 isl_tab_free(tabs[i]);
1169 free(tabs);
1170 return NULL;
1173 /* For each pair of basic sets in the set, check if the union of the two
1174 * can be represented by a single basic set.
1175 * If so, replace the pair by the single basic set and start over.
1177 struct isl_set *isl_set_coalesce(struct isl_set *set)
1179 return (struct isl_set *)isl_map_coalesce((struct isl_map *)set);