add isl_{map,set}_add
[isl.git] / isl_coalesce.c
blobad602c2120e65b443bf6d10a118de709942ba45c
1 /*
2 * Copyright 2008-2009 Katholieke Universiteit Leuven
4 * Use of this software is governed by the GNU LGPLv2.1 license
6 * Written by Sven Verdoolaege, K.U.Leuven, Departement
7 * Computerwetenschappen, Celestijnenlaan 200A, B-3001 Leuven, Belgium
8 */
10 #include "isl_map_private.h"
11 #include "isl_seq.h"
12 #include "isl_tab.h"
14 #define STATUS_ERROR -1
15 #define STATUS_REDUNDANT 1
16 #define STATUS_VALID 2
17 #define STATUS_SEPARATE 3
18 #define STATUS_CUT 4
19 #define STATUS_ADJ_EQ 5
20 #define STATUS_ADJ_INEQ 6
22 static int status_in(isl_int *ineq, struct isl_tab *tab)
24 enum isl_ineq_type type = isl_tab_ineq_type(tab, ineq);
25 switch (type) {
26 case isl_ineq_error: return STATUS_ERROR;
27 case isl_ineq_redundant: return STATUS_VALID;
28 case isl_ineq_separate: return STATUS_SEPARATE;
29 case isl_ineq_cut: return STATUS_CUT;
30 case isl_ineq_adj_eq: return STATUS_ADJ_EQ;
31 case isl_ineq_adj_ineq: return STATUS_ADJ_INEQ;
35 /* Compute the position of the equalities of basic map "i"
36 * with respect to basic map "j".
37 * The resulting array has twice as many entries as the number
38 * of equalities corresponding to the two inequalties to which
39 * each equality corresponds.
41 static int *eq_status_in(struct isl_map *map, int i, int j,
42 struct isl_tab **tabs)
44 int k, l;
45 int *eq = isl_calloc_array(map->ctx, int, 2 * map->p[i]->n_eq);
46 unsigned dim;
48 dim = isl_basic_map_total_dim(map->p[i]);
49 for (k = 0; k < map->p[i]->n_eq; ++k) {
50 for (l = 0; l < 2; ++l) {
51 isl_seq_neg(map->p[i]->eq[k], map->p[i]->eq[k], 1+dim);
52 eq[2 * k + l] = status_in(map->p[i]->eq[k], tabs[j]);
53 if (eq[2 * k + l] == STATUS_ERROR)
54 goto error;
56 if (eq[2 * k] == STATUS_SEPARATE ||
57 eq[2 * k + 1] == STATUS_SEPARATE)
58 break;
61 return eq;
62 error:
63 free(eq);
64 return NULL;
67 /* Compute the position of the inequalities of basic map "i"
68 * with respect to basic map "j".
70 static int *ineq_status_in(struct isl_map *map, int i, int j,
71 struct isl_tab **tabs)
73 int k;
74 unsigned n_eq = map->p[i]->n_eq;
75 int *ineq = isl_calloc_array(map->ctx, int, map->p[i]->n_ineq);
77 for (k = 0; k < map->p[i]->n_ineq; ++k) {
78 if (isl_tab_is_redundant(tabs[i], n_eq + k)) {
79 ineq[k] = STATUS_REDUNDANT;
80 continue;
82 ineq[k] = status_in(map->p[i]->ineq[k], tabs[j]);
83 if (ineq[k] == STATUS_ERROR)
84 goto error;
85 if (ineq[k] == STATUS_SEPARATE)
86 break;
89 return ineq;
90 error:
91 free(ineq);
92 return NULL;
95 static int any(int *con, unsigned len, int status)
97 int i;
99 for (i = 0; i < len ; ++i)
100 if (con[i] == status)
101 return 1;
102 return 0;
105 static int count(int *con, unsigned len, int status)
107 int i;
108 int c = 0;
110 for (i = 0; i < len ; ++i)
111 if (con[i] == status)
112 c++;
113 return c;
116 static int all(int *con, unsigned len, int status)
118 int i;
120 for (i = 0; i < len ; ++i) {
121 if (con[i] == STATUS_REDUNDANT)
122 continue;
123 if (con[i] != status)
124 return 0;
126 return 1;
129 static void drop(struct isl_map *map, int i, struct isl_tab **tabs)
131 isl_basic_map_free(map->p[i]);
132 isl_tab_free(tabs[i]);
134 if (i != map->n - 1) {
135 map->p[i] = map->p[map->n - 1];
136 tabs[i] = tabs[map->n - 1];
138 tabs[map->n - 1] = NULL;
139 map->n--;
142 /* Replace the pair of basic maps i and j but the basic map bounded
143 * by the valid constraints in both basic maps.
145 static int fuse(struct isl_map *map, int i, int j, struct isl_tab **tabs,
146 int *ineq_i, int *ineq_j)
148 int k, l;
149 struct isl_basic_map *fused = NULL;
150 struct isl_tab *fused_tab = NULL;
151 unsigned total = isl_basic_map_total_dim(map->p[i]);
153 fused = isl_basic_map_alloc_dim(isl_dim_copy(map->p[i]->dim),
154 map->p[i]->n_div,
155 map->p[i]->n_eq + map->p[j]->n_eq,
156 map->p[i]->n_ineq + map->p[j]->n_ineq);
157 if (!fused)
158 goto error;
160 for (k = 0; k < map->p[i]->n_eq; ++k) {
161 int l = isl_basic_map_alloc_equality(fused);
162 isl_seq_cpy(fused->eq[l], map->p[i]->eq[k], 1 + total);
165 for (k = 0; k < map->p[j]->n_eq; ++k) {
166 int l = isl_basic_map_alloc_equality(fused);
167 isl_seq_cpy(fused->eq[l], map->p[j]->eq[k], 1 + total);
170 for (k = 0; k < map->p[i]->n_ineq; ++k) {
171 if (ineq_i[k] != STATUS_VALID)
172 continue;
173 l = isl_basic_map_alloc_inequality(fused);
174 isl_seq_cpy(fused->ineq[l], map->p[i]->ineq[k], 1 + total);
177 for (k = 0; k < map->p[j]->n_ineq; ++k) {
178 if (ineq_j[k] != STATUS_VALID)
179 continue;
180 l = isl_basic_map_alloc_inequality(fused);
181 isl_seq_cpy(fused->ineq[l], map->p[j]->ineq[k], 1 + total);
184 for (k = 0; k < map->p[i]->n_div; ++k) {
185 int l = isl_basic_map_alloc_div(fused);
186 isl_seq_cpy(fused->div[l], map->p[i]->div[k], 1 + 1 + total);
189 fused = isl_basic_map_gauss(fused, NULL);
190 ISL_F_SET(fused, ISL_BASIC_MAP_FINAL);
191 if (ISL_F_ISSET(map->p[i], ISL_BASIC_MAP_RATIONAL) &&
192 ISL_F_ISSET(map->p[j], ISL_BASIC_MAP_RATIONAL))
193 ISL_F_SET(fused, ISL_BASIC_MAP_RATIONAL);
195 fused_tab = isl_tab_from_basic_map(fused);
196 if (isl_tab_detect_redundant(fused_tab) < 0)
197 goto error;
199 isl_basic_map_free(map->p[i]);
200 map->p[i] = fused;
201 isl_tab_free(tabs[i]);
202 tabs[i] = fused_tab;
203 drop(map, j, tabs);
205 return 1;
206 error:
207 isl_tab_free(fused_tab);
208 isl_basic_map_free(fused);
209 return -1;
212 /* Given a pair of basic maps i and j such that all constraints are either
213 * "valid" or "cut", check if the facets corresponding to the "cut"
214 * constraints of i lie entirely within basic map j.
215 * If so, replace the pair by the basic map consisting of the valid
216 * constraints in both basic maps.
218 * To see that we are not introducing any extra points, call the
219 * two basic maps A and B and the resulting map U and let x
220 * be an element of U \setminus ( A \cup B ).
221 * Then there is a pair of cut constraints c_1 and c_2 in A and B such that x
222 * violates them. Let X be the intersection of U with the opposites
223 * of these constraints. Then x \in X.
224 * The facet corresponding to c_1 contains the corresponding facet of A.
225 * This facet is entirely contained in B, so c_2 is valid on the facet.
226 * However, since it is also (part of) a facet of X, -c_2 is also valid
227 * on the facet. This means c_2 is saturated on the facet, so c_1 and
228 * c_2 must be opposites of each other, but then x could not violate
229 * both of them.
231 static int check_facets(struct isl_map *map, int i, int j,
232 struct isl_tab **tabs, int *ineq_i, int *ineq_j)
234 int k, l;
235 struct isl_tab_undo *snap;
236 unsigned n_eq = map->p[i]->n_eq;
238 snap = isl_tab_snap(tabs[i]);
240 for (k = 0; k < map->p[i]->n_ineq; ++k) {
241 if (ineq_i[k] != STATUS_CUT)
242 continue;
243 tabs[i] = isl_tab_select_facet(tabs[i], n_eq + k);
244 for (l = 0; l < map->p[j]->n_ineq; ++l) {
245 int stat;
246 if (ineq_j[l] != STATUS_CUT)
247 continue;
248 stat = status_in(map->p[j]->ineq[l], tabs[i]);
249 if (stat != STATUS_VALID)
250 break;
252 if (isl_tab_rollback(tabs[i], snap) < 0)
253 return -1;
254 if (l < map->p[j]->n_ineq)
255 break;
258 if (k < map->p[i]->n_ineq)
259 /* BAD CUT PAIR */
260 return 0;
261 return fuse(map, i, j, tabs, ineq_i, ineq_j);
264 /* Both basic maps have at least one inequality with and adjacent
265 * (but opposite) inequality in the other basic map.
266 * Check that there are no cut constraints and that there is only
267 * a single pair of adjacent inequalities.
268 * If so, we can replace the pair by a single basic map described
269 * by all but the pair of adjacent inequalities.
270 * Any additional points introduced lie strictly between the two
271 * adjacent hyperplanes and can therefore be integral.
273 * ____ _____
274 * / ||\ / \
275 * / || \ / \
276 * \ || \ => \ \
277 * \ || / \ /
278 * \___||_/ \_____/
280 * The test for a single pair of adjancent inequalities is important
281 * for avoiding the combination of two basic maps like the following
283 * /|
284 * / |
285 * /__|
286 * _____
287 * | |
288 * | |
289 * |___|
291 static int check_adj_ineq(struct isl_map *map, int i, int j,
292 struct isl_tab **tabs, int *ineq_i, int *ineq_j)
294 int changed = 0;
296 if (any(ineq_i, map->p[i]->n_ineq, STATUS_CUT) ||
297 any(ineq_j, map->p[j]->n_ineq, STATUS_CUT))
298 /* ADJ INEQ CUT */
300 else if (count(ineq_i, map->p[i]->n_ineq, STATUS_ADJ_INEQ) == 1 &&
301 count(ineq_j, map->p[j]->n_ineq, STATUS_ADJ_INEQ) == 1)
302 changed = fuse(map, i, j, tabs, ineq_i, ineq_j);
303 /* else ADJ INEQ TOO MANY */
305 return changed;
308 /* Check if basic map "i" contains the basic map represented
309 * by the tableau "tab".
311 static int contains(struct isl_map *map, int i, int *ineq_i,
312 struct isl_tab *tab)
314 int k, l;
315 unsigned dim;
317 dim = isl_basic_map_total_dim(map->p[i]);
318 for (k = 0; k < map->p[i]->n_eq; ++k) {
319 for (l = 0; l < 2; ++l) {
320 int stat;
321 isl_seq_neg(map->p[i]->eq[k], map->p[i]->eq[k], 1+dim);
322 stat = status_in(map->p[i]->eq[k], tab);
323 if (stat != STATUS_VALID)
324 return 0;
328 for (k = 0; k < map->p[i]->n_ineq; ++k) {
329 int stat;
330 if (ineq_i[k] == STATUS_REDUNDANT)
331 continue;
332 stat = status_in(map->p[i]->ineq[k], tab);
333 if (stat != STATUS_VALID)
334 return 0;
336 return 1;
339 /* At least one of the basic maps has an equality that is adjacent
340 * to inequality. Make sure that only one of the basic maps has
341 * such an equality and that the other basic map has exactly one
342 * inequality adjacent to an equality.
343 * We call the basic map that has the inequality "i" and the basic
344 * map that has the equality "j".
345 * If "i" has any "cut" inequality, then relaxing the inequality
346 * by one would not result in a basic map that contains the other
347 * basic map.
348 * Otherwise, we relax the constraint, compute the corresponding
349 * facet and check whether it is included in the other basic map.
350 * If so, we know that relaxing the constraint extend the basic
351 * map with exactly the other basic map (we already know that this
352 * other basic map is included in the extension, because there
353 * were no "cut" inequalities in "i") and we can replace the
354 * two basic maps by thie extension.
355 * ____ _____
356 * / || / |
357 * / || / |
358 * \ || => \ |
359 * \ || \ |
360 * \___|| \____|
362 static int check_adj_eq(struct isl_map *map, int i, int j,
363 struct isl_tab **tabs, int *eq_i, int *ineq_i, int *eq_j, int *ineq_j)
365 int changed = 0;
366 int super;
367 int k;
368 struct isl_tab_undo *snap, *snap2;
369 unsigned n_eq = map->p[i]->n_eq;
371 if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_ADJ_INEQ) &&
372 any(eq_j, 2 * map->p[j]->n_eq, STATUS_ADJ_INEQ))
373 /* ADJ EQ TOO MANY */
374 return 0;
376 if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_ADJ_INEQ))
377 return check_adj_eq(map, j, i, tabs,
378 eq_j, ineq_j, eq_i, ineq_i);
380 /* j has an equality adjacent to an inequality in i */
382 if (any(ineq_i, map->p[i]->n_ineq, STATUS_CUT))
383 /* ADJ EQ CUT */
384 return 0;
385 if (count(eq_j, 2 * map->p[j]->n_eq, STATUS_ADJ_INEQ) != 1 ||
386 count(ineq_i, map->p[i]->n_ineq, STATUS_ADJ_EQ) != 1 ||
387 any(ineq_j, map->p[j]->n_ineq, STATUS_ADJ_EQ) ||
388 any(ineq_i, map->p[i]->n_ineq, STATUS_ADJ_INEQ) ||
389 any(ineq_j, map->p[j]->n_ineq, STATUS_ADJ_INEQ))
390 /* ADJ EQ TOO MANY */
391 return 0;
393 for (k = 0; k < map->p[i]->n_ineq ; ++k)
394 if (ineq_i[k] == STATUS_ADJ_EQ)
395 break;
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 tabs[i] = isl_tab_select_facet(tabs[i], n_eq + k);
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 /* Check if the union of the given pair of basic maps
420 * can be represented by a single basic map.
421 * If so, replace the pair by the single basic map and return 1.
422 * Otherwise, return 0;
424 * We first check the effect of each constraint of one basic map
425 * on the other basic map.
426 * The constraint may be
427 * redundant the constraint is redundant in its own
428 * basic map and should be ignore and removed
429 * in the end
430 * valid all (integer) points of the other basic map
431 * satisfy the constraint
432 * separate no (integer) point of the other basic map
433 * satisfies the constraint
434 * cut some but not all points of the other basic map
435 * satisfy the constraint
436 * adj_eq the given constraint is adjacent (on the outside)
437 * to an equality of the other basic map
438 * adj_ineq the given constraint is adjacent (on the outside)
439 * to an inequality of the other basic map
441 * We consider four cases in which we can replace the pair by a single
442 * basic map. We ignore all "redundant" constraints.
444 * 1. all constraints of one basic map are valid
445 * => the other basic map is a subset and can be removed
447 * 2. all constraints of both basic maps are either "valid" or "cut"
448 * and the facets corresponding to the "cut" constraints
449 * of one of the basic maps lies entirely inside the other basic map
450 * => the pair can be replaced by a basic map consisting
451 * of the valid constraints in both basic maps
453 * 3. there is a single pair of adjacent inequalities
454 * (all other constraints are "valid")
455 * => the pair can be replaced by a basic map consisting
456 * of the valid constraints in both basic maps
458 * 4. there is a single adjacent pair of an inequality and an equality,
459 * the other constraints of the basic map containing the inequality are
460 * "valid". Moreover, if the inequality the basic map is relaxed
461 * and then turned into an equality, then resulting facet lies
462 * entirely inside the other basic map
463 * => the pair can be replaced by the basic map containing
464 * the inequality, with the inequality relaxed.
466 * Throughout the computation, we maintain a collection of tableaus
467 * corresponding to the basic maps. When the basic maps are dropped
468 * or combined, the tableaus are modified accordingly.
470 static int coalesce_pair(struct isl_map *map, int i, int j,
471 struct isl_tab **tabs)
473 int changed = 0;
474 int *eq_i = NULL;
475 int *eq_j = NULL;
476 int *ineq_i = NULL;
477 int *ineq_j = NULL;
479 eq_i = eq_status_in(map, i, j, tabs);
480 if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_ERROR))
481 goto error;
482 if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_SEPARATE))
483 goto done;
485 eq_j = eq_status_in(map, j, i, tabs);
486 if (any(eq_j, 2 * map->p[j]->n_eq, STATUS_ERROR))
487 goto error;
488 if (any(eq_j, 2 * map->p[j]->n_eq, STATUS_SEPARATE))
489 goto done;
491 ineq_i = ineq_status_in(map, i, j, tabs);
492 if (any(ineq_i, map->p[i]->n_ineq, STATUS_ERROR))
493 goto error;
494 if (any(ineq_i, map->p[i]->n_ineq, STATUS_SEPARATE))
495 goto done;
497 ineq_j = ineq_status_in(map, j, i, tabs);
498 if (any(ineq_j, map->p[j]->n_ineq, STATUS_ERROR))
499 goto error;
500 if (any(ineq_j, map->p[j]->n_ineq, STATUS_SEPARATE))
501 goto done;
503 if (all(eq_i, 2 * map->p[i]->n_eq, STATUS_VALID) &&
504 all(ineq_i, map->p[i]->n_ineq, STATUS_VALID)) {
505 drop(map, j, tabs);
506 changed = 1;
507 } else if (all(eq_j, 2 * map->p[j]->n_eq, STATUS_VALID) &&
508 all(ineq_j, map->p[j]->n_ineq, STATUS_VALID)) {
509 drop(map, i, tabs);
510 changed = 1;
511 } else if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_CUT) ||
512 any(eq_j, 2 * map->p[j]->n_eq, STATUS_CUT)) {
513 /* BAD CUT */
514 } else if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_ADJ_EQ) ||
515 any(eq_j, 2 * map->p[j]->n_eq, STATUS_ADJ_EQ)) {
516 /* ADJ EQ PAIR */
517 } else if (any(eq_i, 2 * map->p[i]->n_eq, STATUS_ADJ_INEQ) ||
518 any(eq_j, 2 * map->p[j]->n_eq, STATUS_ADJ_INEQ)) {
519 changed = check_adj_eq(map, i, j, tabs,
520 eq_i, ineq_i, eq_j, ineq_j);
521 } else if (any(ineq_i, map->p[i]->n_ineq, STATUS_ADJ_EQ) ||
522 any(ineq_j, map->p[j]->n_ineq, STATUS_ADJ_EQ)) {
523 /* Can't happen */
524 /* BAD ADJ INEQ */
525 } else if (any(ineq_i, map->p[i]->n_ineq, STATUS_ADJ_INEQ) ||
526 any(ineq_j, map->p[j]->n_ineq, STATUS_ADJ_INEQ)) {
527 changed = check_adj_ineq(map, i, j, tabs, ineq_i, ineq_j);
528 } else
529 changed = check_facets(map, i, j, tabs, ineq_i, ineq_j);
531 done:
532 free(eq_i);
533 free(eq_j);
534 free(ineq_i);
535 free(ineq_j);
536 return changed;
537 error:
538 free(eq_i);
539 free(eq_j);
540 free(ineq_i);
541 free(ineq_j);
542 return -1;
545 static struct isl_map *coalesce(struct isl_map *map, struct isl_tab **tabs)
547 int i, j;
549 for (i = 0; i < map->n - 1; ++i)
550 for (j = i + 1; j < map->n; ++j) {
551 int changed;
552 changed = coalesce_pair(map, i, j, tabs);
553 if (changed < 0)
554 goto error;
555 if (changed)
556 return coalesce(map, tabs);
558 return map;
559 error:
560 isl_map_free(map);
561 return NULL;
564 /* For each pair of basic maps in the map, check if the union of the two
565 * can be represented by a single basic map.
566 * If so, replace the pair by the single basic map and start over.
568 struct isl_map *isl_map_coalesce(struct isl_map *map)
570 int i;
571 unsigned n;
572 struct isl_tab **tabs = NULL;
574 if (!map)
575 return NULL;
577 if (map->n <= 1)
578 return map;
580 map = isl_map_align_divs(map);
582 tabs = isl_calloc_array(map->ctx, struct isl_tab *, map->n);
583 if (!tabs)
584 goto error;
586 n = map->n;
587 for (i = 0; i < map->n; ++i) {
588 tabs[i] = isl_tab_from_basic_map(map->p[i]);
589 if (!tabs[i])
590 goto error;
591 if (!ISL_F_ISSET(map->p[i], ISL_BASIC_MAP_NO_IMPLICIT))
592 tabs[i] = isl_tab_detect_implicit_equalities(tabs[i]);
593 if (!ISL_F_ISSET(map->p[i], ISL_BASIC_MAP_NO_REDUNDANT))
594 if (isl_tab_detect_redundant(tabs[i]) < 0)
595 goto error;
597 for (i = map->n - 1; i >= 0; --i)
598 if (tabs[i]->empty)
599 drop(map, i, tabs);
601 map = coalesce(map, tabs);
603 if (map)
604 for (i = 0; i < map->n; ++i) {
605 map->p[i] = isl_basic_map_update_from_tab(map->p[i],
606 tabs[i]);
607 map->p[i] = isl_basic_map_finalize(map->p[i]);
608 if (!map->p[i])
609 goto error;
610 ISL_F_SET(map->p[i], ISL_BASIC_MAP_NO_IMPLICIT);
611 ISL_F_SET(map->p[i], ISL_BASIC_MAP_NO_REDUNDANT);
614 for (i = 0; i < n; ++i)
615 isl_tab_free(tabs[i]);
617 free(tabs);
619 return map;
620 error:
621 if (tabs)
622 for (i = 0; i < n; ++i)
623 isl_tab_free(tabs[i]);
624 free(tabs);
625 return NULL;
628 /* For each pair of basic sets in the set, check if the union of the two
629 * can be represented by a single basic set.
630 * If so, replace the pair by the single basic set and start over.
632 struct isl_set *isl_set_coalesce(struct isl_set *set)
634 return (struct isl_set *)isl_map_coalesce((struct isl_map *)set);