isl_qpolynomial_move_dims: plug memory leak and fix reordering of divs
[isl.git] / isl_vertices.c
blob5d2b718edaf45483307108ac0160137b62e6a369
1 /*
2 * Copyright 2010 INRIA Saclay
4 * Use of this software is governed by the GNU LGPLv2.1 license
6 * Written by Sven Verdoolaege, INRIA Saclay - Ile-de-France,
7 * Parc Club Orsay Universite, ZAC des vignes, 4 rue Jacques Monod,
8 * 91893 Orsay, France
9 */
11 #include <isl_set.h>
12 #include <isl_seq.h>
13 #include <isl_tab.h>
14 #include <isl_map_private.h>
15 #include <isl_dim_private.h>
16 #include <isl_morph.h>
17 #include <isl_vertices_private.h>
19 #define SELECTED 1
20 #define DESELECTED -1
21 #define UNSELECTED 0
23 static __isl_give isl_vertices *compute_chambers(__isl_take isl_basic_set *bset,
24 __isl_take isl_vertices *vertices);
26 __isl_give isl_vertices *isl_vertices_copy(__isl_keep isl_vertices *vertices)
28 if (!vertices)
29 return NULL;
31 vertices->ref++;
32 return vertices;
35 void isl_vertices_free(__isl_take isl_vertices *vertices)
37 int i;
39 if (!vertices)
40 return;
42 if (--vertices->ref > 0)
43 return;
45 for (i = 0; i < vertices->n_vertices; ++i) {
46 isl_basic_set_free(vertices->v[i].vertex);
47 isl_basic_set_free(vertices->v[i].dom);
49 free(vertices->v);
51 for (i = 0; i < vertices->n_chambers; ++i) {
52 free(vertices->c[i].vertices);
53 isl_basic_set_free(vertices->c[i].dom);
55 free(vertices->c);
57 isl_ctx_deref(vertices->ctx);
58 free(vertices);
61 struct isl_vertex_list {
62 struct isl_vertex v;
63 struct isl_vertex_list *next;
66 static void free_vertex_list(struct isl_vertex_list *list)
68 struct isl_vertex_list *next;
70 for (; list; list = next) {
71 next = list->next;
72 isl_basic_set_free(list->v.vertex);
73 isl_basic_set_free(list->v.dom);
74 free(list);
78 static __isl_give isl_vertices *vertices_from_list(__isl_keep isl_basic_set *bset,
79 int n_vertices, struct isl_vertex_list *list)
81 int i;
82 struct isl_vertex_list *next;
83 isl_vertices *vertices;
85 vertices = isl_calloc_type(bset->ctx, isl_vertices);
86 if (!vertices)
87 goto error;
88 vertices->ctx = bset->ctx;
89 isl_ctx_ref(bset->ctx);
90 vertices->ref = 1;
91 vertices->v = isl_alloc_array(bset->ctx, struct isl_vertex, n_vertices);
92 if (!vertices->v)
93 goto error;
94 vertices->n_vertices = n_vertices;
96 for (i = 0; list; list = next, i++) {
97 next = list->next;
98 vertices->v[i] = list->v;
99 free(list);
102 return vertices;
103 error:
104 free(vertices);
105 free_vertex_list(list);
106 return NULL;
109 /* Prepend a vertex to the linked list "list" based on the equalities in "tab".
111 static int add_vertex(struct isl_vertex_list **list,
112 __isl_keep isl_basic_set *bset, struct isl_tab *tab)
114 unsigned nvar;
115 unsigned nparam;
116 struct isl_vertex_list *v = NULL;
118 if (isl_tab_detect_implicit_equalities(tab) < 0)
119 return -1;
121 nvar = isl_basic_set_dim(bset, isl_dim_set);
122 nparam = isl_basic_set_dim(bset, isl_dim_param);
124 v = isl_calloc_type(tab->mat->ctx, struct isl_vertex_list);
125 if (!v)
126 goto error;
128 v->v.vertex = isl_basic_set_copy(bset);
129 v->v.vertex = isl_basic_set_cow(v->v.vertex);
130 v->v.vertex = isl_basic_set_update_from_tab(v->v.vertex, tab);
131 v->v.vertex = isl_basic_set_simplify(v->v.vertex);
132 v->v.vertex = isl_basic_set_finalize(v->v.vertex);
133 if (!v->v.vertex)
134 goto error;
135 isl_assert(bset->ctx, v->v.vertex->n_eq >= nvar, goto error);
136 v->v.dom = isl_basic_set_copy(v->v.vertex);
137 v->v.dom = isl_basic_set_project_out(v->v.dom, isl_dim_set, 0, nvar);
138 if (!v->v.dom)
139 goto error;
141 v->next = *list;
142 *list = v;
144 return 0;
145 error:
146 free_vertex_list(v);
147 return -1;
150 /* Compute the parametric vertices and the chamber decomposition
151 * of an empty parametric polytope.
153 static __isl_give isl_vertices *vertices_empty(__isl_keep isl_basic_set *bset)
155 isl_vertices *vertices;
156 unsigned nparam;
158 if (!bset)
159 return NULL;
161 nparam = isl_basic_set_dim(bset, isl_dim_param);
163 vertices = isl_calloc_type(bset->ctx, isl_vertices);
164 if (!vertices)
165 return NULL;
166 vertices->ctx = bset->ctx;
167 isl_ctx_ref(bset->ctx);
168 vertices->ref = 1;
170 vertices->n_vertices = 0;
171 vertices->n_chambers = 0;
173 return vertices;
176 /* Compute the parametric vertices and the chamber decomposition
177 * of the parametric polytope defined using the same constraints
178 * as "bset" in the 0D case.
179 * There is exactly one 0D vertex and a single chamber containing
180 * the vertex.
182 static __isl_give isl_vertices *vertices_0D(__isl_keep isl_basic_set *bset)
184 isl_vertices *vertices;
185 unsigned nparam;
187 if (!bset)
188 return NULL;
190 nparam = isl_basic_set_dim(bset, isl_dim_param);
192 vertices = isl_calloc_type(bset->ctx, isl_vertices);
193 if (!vertices)
194 return NULL;
195 vertices->ctx = bset->ctx;
196 isl_ctx_ref(bset->ctx);
197 vertices->ref = 1;
199 vertices->v = isl_calloc_array(bset->ctx, struct isl_vertex, 1);
200 if (!vertices->v)
201 goto error;
202 vertices->n_vertices = 1;
203 vertices->v[0].vertex = isl_basic_set_copy(bset);
204 if (!vertices->v[0].vertex)
205 goto error;
207 vertices->c = isl_calloc_array(bset->ctx, struct isl_chamber, 1);
208 if (!vertices->c)
209 goto error;
210 vertices->n_chambers = 1;
211 vertices->c[0].n_vertices = 1;
212 vertices->c[0].vertices = isl_calloc_array(bset->ctx, int, 1);
213 if (!vertices->c[0].vertices)
214 goto error;
215 vertices->c[0].dom = isl_basic_set_copy(bset);
216 if (!vertices->c[0].dom)
217 goto error;
219 return vertices;
220 error:
221 isl_vertices_free(vertices);
222 return NULL;
225 static int isl_mat_rank(__isl_keep isl_mat *mat)
227 int row, col;
228 isl_mat *H;
230 H = isl_mat_left_hermite(isl_mat_copy(mat), 0, NULL, NULL);
231 if (!H)
232 return -1;
234 for (col = 0; col < H->n_col; ++col) {
235 for (row = 0; row < H->n_row; ++row)
236 if (!isl_int_is_zero(H->row[row][col]))
237 break;
238 if (row == H->n_row)
239 break;
242 isl_mat_free(H);
244 return col;
247 /* Is the row pointed to by "f" linearly independent of the "n" first
248 * rows in "facets"?
250 static int is_independent(__isl_keep isl_mat *facets, int n, isl_int *f)
252 int rank;
254 if (isl_seq_first_non_zero(f, facets->n_col) < 0)
255 return 0;
257 isl_seq_cpy(facets->row[n], f, facets->n_col);
258 facets->n_row = n + 1;
259 rank = isl_mat_rank(facets);
260 if (rank < 0)
261 return -1;
263 return rank == n + 1;
266 /* Check whether we can select constraint "level", given the current selection
267 * reflected by facets in "tab", the rows of "facets" and the earlier
268 * "selected" elements of "selection".
270 * If the constraint is (strictly) redundant in the tableau, selecting it would
271 * result in an empty tableau, so it can't be selected.
272 * If the set variable part of the constraint is not linearly indepedent
273 * of the set variable parts of the already selected constraints,
274 * the constraint cannot be selected.
275 * If selecting the constraint results in an empty tableau, the constraint
276 * cannot be selected.
277 * Finally, if selecting the constraint results in some explicitly
278 * deselected constraints turning into equalities, then the corresponding
279 * vertices have already been generated, so the constraint cannot be selected.
281 static int can_select(__isl_keep isl_basic_set *bset, int level,
282 struct isl_tab *tab, __isl_keep isl_mat *facets, int selected,
283 int *selection)
285 int i;
286 int indep;
287 unsigned ovar;
288 struct isl_tab_undo *snap;
290 if (isl_tab_is_redundant(tab, level))
291 return 0;
293 ovar = isl_dim_offset(bset->dim, isl_dim_set);
295 indep = is_independent(facets, selected, bset->ineq[level] + 1 + ovar);
296 if (indep < 0)
297 return -1;
298 if (!indep)
299 return 0;
301 snap = isl_tab_snap(tab);
302 if (isl_tab_select_facet(tab, level) < 0)
303 return -1;
305 if (tab->empty) {
306 if (isl_tab_rollback(tab, snap) < 0)
307 return -1;
308 return 0;
311 for (i = 0; i < level; ++i) {
312 int sgn;
314 if (selection[i] != DESELECTED)
315 continue;
317 if (isl_tab_is_equality(tab, i))
318 sgn = 0;
319 else if (isl_tab_is_redundant(tab, i))
320 sgn = 1;
321 else
322 sgn = isl_tab_sign_of_max(tab, i);
323 if (sgn < -1)
324 return -1;
325 if (sgn <= 0) {
326 if (isl_tab_rollback(tab, snap) < 0)
327 return -1;
328 return 0;
332 return 1;
335 /* Compute the parametric vertices and the chamber decomposition
336 * of a parametric polytope that is not full-dimensional.
338 * Simply map the parametric polytope to a lower dimensional space
339 * and map the resulting vertices back.
341 static __isl_give isl_vertices *lower_dim_vertices(
342 __isl_keep isl_basic_set *bset)
344 isl_morph *morph;
345 isl_vertices *vertices;
347 bset = isl_basic_set_copy(bset);
348 morph = isl_basic_set_full_compression(bset);
349 bset = isl_morph_basic_set(isl_morph_copy(morph), bset);
351 vertices = isl_basic_set_compute_vertices(bset);
352 isl_basic_set_free(bset);
354 morph = isl_morph_inverse(morph);
356 vertices = isl_morph_vertices(morph, vertices);
358 return vertices;
361 /* Compute the parametric vertices and the chamber decomposition
362 * of the parametric polytope defined using the same constraints
363 * as "bset". "bset" is assumed to have no existentially quantified
364 * variables.
366 * The vertices themselves are computed in a fairly simplistic way.
367 * We simply run through all combinations of d constraints,
368 * with d the number of set variables, and check if those d constraints
369 * define a vertex. To avoid the generation of duplicate vertices,
370 * which we may happen if a vertex is defined by more that d constraints,
371 * we make sure we only generate the vertex for the d constraints with
372 * smallest index.
374 * We set up a tableau and keep track of which facets have been
375 * selected. The tableau is marked strict_redundant so that we can be
376 * sure that any constraint that is marked redundant (and that is not
377 * also marked zero) is not an equality.
378 * If a constraint is marked DESELECTED, it means the constraint was
379 * SELECTED before (in combination with the same selection of earlier
380 * constraints). If such a deselected constraint turns out to be an
381 * equality, then any vertex that may still be found with the current
382 * selection has already been generated when the constraint was selected.
383 * A constraint is marked UNSELECTED when there is no way selecting
384 * the constraint could lead to a vertex (in combination with the current
385 * selection of earlier constraints).
387 * The set variable coefficients of the selected constraints are stored
388 * in the facets matrix.
390 __isl_give isl_vertices *isl_basic_set_compute_vertices(
391 __isl_keep isl_basic_set *bset)
393 struct isl_tab *tab;
394 int level;
395 int init;
396 unsigned nvar;
397 int *selection;
398 int selected;
399 struct isl_tab_undo **snap;
400 isl_mat *facets;
401 struct isl_vertex_list *list = NULL;
402 int n_vertices = 0;
403 isl_vertices *vertices;
405 if (!bset)
406 return NULL;
408 if (isl_basic_set_fast_is_empty(bset))
409 return vertices_empty(bset);
411 if (bset->n_eq != 0)
412 return lower_dim_vertices(bset);
414 isl_assert(bset->ctx, isl_basic_set_dim(bset, isl_dim_div) == 0,
415 return NULL);
417 if (isl_basic_set_dim(bset, isl_dim_set) == 0)
418 return vertices_0D(bset);
420 nvar = isl_basic_set_dim(bset, isl_dim_set);
422 bset = isl_basic_set_copy(bset);
423 bset = isl_basic_set_set_rational(bset);
424 if (!bset)
425 return NULL;
427 tab = isl_tab_from_basic_set(bset);
428 if (!tab)
429 goto error;
430 tab->strict_redundant = 1;
432 if (tab->empty) {
433 vertices = vertices_empty(bset);
434 isl_basic_set_free(bset);
435 isl_tab_free(tab);
436 return vertices;
439 selection = isl_alloc_array(bset->ctx, int, bset->n_ineq);
440 snap = isl_alloc_array(bset->ctx, struct isl_tab_undo *, bset->n_ineq);
441 facets = isl_mat_alloc(bset->ctx, nvar, nvar);
442 if (!selection || !snap || !facets)
443 goto error;
445 level = 0;
446 init = 1;
447 selected = 0;
449 while (level >= 0) {
450 if (level >= bset->n_ineq ||
451 (!init && selection[level] != SELECTED)) {
452 --level;
453 init = 0;
454 continue;
456 if (init) {
457 int ok;
458 snap[level] = isl_tab_snap(tab);
459 ok = can_select(bset, level, tab, facets, selected,
460 selection);
461 if (ok < 0)
462 goto error;
463 if (ok) {
464 selection[level] = SELECTED;
465 selected++;
466 } else
467 selection[level] = UNSELECTED;
468 } else {
469 selection[level] = DESELECTED;
470 selected--;
471 if (isl_tab_rollback(tab, snap[level]) < 0)
472 goto error;
474 if (selected == nvar) {
475 if (tab->n_dead == nvar) {
476 if (add_vertex(&list, bset, tab) < 0)
477 goto error;
478 n_vertices++;
480 init = 0;
481 continue;
483 ++level;
484 init = 1;
487 isl_mat_free(facets);
488 free(selection);
489 free(snap);
491 isl_tab_free(tab);
493 vertices = vertices_from_list(bset, n_vertices, list);
495 vertices = compute_chambers(bset, vertices);
497 return vertices;
498 error:
499 isl_mat_free(facets);
500 free(selection);
501 free(snap);
502 isl_tab_free(tab);
503 isl_basic_set_free(bset);
504 return NULL;
507 struct isl_chamber_list {
508 struct isl_chamber c;
509 struct isl_chamber_list *next;
512 static void free_chamber_list(struct isl_chamber_list *list)
514 struct isl_chamber_list *next;
516 for (; list; list = next) {
517 next = list->next;
518 isl_basic_set_free(list->c.dom);
519 free(list->c.vertices);
520 free(list);
524 /* Check whether the basic set "bset" is a superset of the basic set described
525 * by "tab", i.e., check whether all constraints of "bset" are redundant.
527 static int bset_covers_tab(__isl_keep isl_basic_set *bset, struct isl_tab *tab)
529 int i;
531 if (!bset || !tab)
532 return -1;
534 for (i = 0; i < bset->n_ineq; ++i) {
535 enum isl_ineq_type type = isl_tab_ineq_type(tab, bset->ineq[i]);
536 switch (type) {
537 case isl_ineq_error: return -1;
538 case isl_ineq_redundant: continue;
539 default: return 0;
543 return 1;
546 static __isl_give isl_vertices *vertices_add_chambers(
547 __isl_take isl_vertices *vertices, int n_chambers,
548 struct isl_chamber_list *list)
550 int i;
551 struct isl_chamber_list *next;
553 vertices->c = isl_alloc_array(vertices->ctx, struct isl_chamber, n_chambers);
554 if (!vertices->c)
555 goto error;
556 vertices->n_chambers = n_chambers;
558 for (i = 0; list; list = next, i++) {
559 next = list->next;
560 vertices->c[i] = list->c;
561 free(list);
564 return vertices;
565 error:
566 isl_vertices_free(vertices);
567 free_chamber_list(list);
568 return NULL;
571 /* Can "tab" be intersected with "bset" without resulting in
572 * a lower-dimensional set.
574 static int can_intersect(struct isl_tab *tab, __isl_keep isl_basic_set *bset)
576 int i;
577 struct isl_tab_undo *snap;
579 if (isl_tab_extend_cons(tab, bset->n_ineq) < 0)
580 return -1;
582 snap = isl_tab_snap(tab);
584 for (i = 0; i < bset->n_ineq; ++i) {
585 if (isl_tab_ineq_type(tab, bset->ineq[i]) == isl_ineq_redundant)
586 continue;
587 if (isl_tab_add_ineq(tab, bset->ineq[i]) < 0)
588 return -1;
591 if (isl_tab_detect_implicit_equalities(tab) < 0)
592 return -1;
593 if (tab->n_dead) {
594 if (isl_tab_rollback(tab, snap) < 0)
595 return -1;
596 return 0;
599 return 1;
602 static int add_chamber(struct isl_chamber_list **list,
603 __isl_keep isl_vertices *vertices, struct isl_tab *tab, int *selection)
605 int n_frozen;
606 int i, j;
607 int n_vertices = 0;
608 struct isl_tab_undo *snap;
609 struct isl_chamber_list *c = NULL;
611 for (i = 0; i < vertices->n_vertices; ++i)
612 if (selection[i])
613 n_vertices++;
615 snap = isl_tab_snap(tab);
617 for (i = 0; i < tab->n_con && tab->con[i].frozen; ++i)
618 tab->con[i].frozen = 0;
619 n_frozen = i;
621 if (isl_tab_detect_redundant(tab) < 0)
622 return -1;
624 c = isl_calloc_type(tab->mat->ctx, struct isl_chamber_list);
625 if (!c)
626 goto error;
627 c->c.vertices = isl_alloc_array(tab->mat->ctx, int, n_vertices);
628 if (!c->c.vertices)
629 goto error;
630 c->c.dom = isl_basic_set_from_basic_map(isl_basic_map_copy(tab->bmap));
631 c->c.dom = isl_basic_set_set_rational(c->c.dom);
632 c->c.dom = isl_basic_set_cow(c->c.dom);
633 c->c.dom = isl_basic_set_update_from_tab(c->c.dom, tab);
634 c->c.dom = isl_basic_set_simplify(c->c.dom);
635 c->c.dom = isl_basic_set_finalize(c->c.dom);
636 if (!c->c.dom)
637 goto error;
639 c->c.n_vertices = n_vertices;
641 for (i = 0, j = 0; i < vertices->n_vertices; ++i)
642 if (selection[i]) {
643 c->c.vertices[j] = i;
644 j++;
647 c->next = *list;
648 *list = c;
650 for (i = 0; i < n_frozen; ++i)
651 tab->con[i].frozen = 1;
653 if (isl_tab_rollback(tab, snap) < 0)
654 return -1;
656 return 0;
657 error:
658 free_chamber_list(c);
659 return -1;
662 struct isl_facet_todo {
663 struct isl_tab *tab; /* A tableau representation of the facet */
664 isl_basic_set *bset; /* A normalized basic set representation */
665 isl_vec *constraint; /* Constraint pointing to the other side */
666 struct isl_facet_todo *next;
669 static void free_todo(struct isl_facet_todo *todo)
671 while (todo) {
672 struct isl_facet_todo *next = todo->next;
674 isl_tab_free(todo->tab);
675 isl_basic_set_free(todo->bset);
676 isl_vec_free(todo->constraint);
677 free(todo);
679 todo = next;
683 static struct isl_facet_todo *create_todo(struct isl_tab *tab, int con)
685 int i;
686 int n_frozen;
687 struct isl_tab_undo *snap;
688 struct isl_facet_todo *todo;
690 snap = isl_tab_snap(tab);
692 for (i = 0; i < tab->n_con && tab->con[i].frozen; ++i)
693 tab->con[i].frozen = 0;
694 n_frozen = i;
696 if (isl_tab_detect_redundant(tab) < 0)
697 return NULL;
699 todo = isl_calloc_type(tab->mat->ctx, struct isl_facet_todo);
700 if (!todo)
701 return NULL;
703 todo->constraint = isl_vec_alloc(tab->mat->ctx, 1 + tab->n_var);
704 if (!todo->constraint)
705 goto error;
706 isl_seq_neg(todo->constraint->el, tab->bmap->ineq[con], 1 + tab->n_var);
707 todo->bset = isl_basic_set_from_basic_map(isl_basic_map_copy(tab->bmap));
708 todo->bset = isl_basic_set_set_rational(todo->bset);
709 todo->bset = isl_basic_set_cow(todo->bset);
710 todo->bset = isl_basic_set_update_from_tab(todo->bset, tab);
711 todo->bset = isl_basic_set_simplify(todo->bset);
712 todo->bset = isl_basic_set_sort_constraints(todo->bset);
713 if (!todo->bset)
714 goto error;
715 ISL_F_SET(todo->bset, ISL_BASIC_SET_NORMALIZED);
716 todo->tab = isl_tab_dup(tab);
717 if (!todo->tab)
718 goto error;
720 for (i = 0; i < n_frozen; ++i)
721 tab->con[i].frozen = 1;
723 if (isl_tab_rollback(tab, snap) < 0)
724 goto error;
726 return todo;
727 error:
728 free_todo(todo);
729 return NULL;
732 /* Create todo items for all interior facets of the chamber represented
733 * by "tab" and collect them in "next".
735 static int init_todo(struct isl_facet_todo **next, struct isl_tab *tab)
737 int i;
738 struct isl_tab_undo *snap;
739 struct isl_facet_todo *todo;
741 snap = isl_tab_snap(tab);
743 for (i = 0; i < tab->n_con; ++i) {
744 if (tab->con[i].frozen)
745 continue;
746 if (tab->con[i].is_redundant)
747 continue;
749 if (isl_tab_select_facet(tab, i) < 0)
750 return -1;
752 todo = create_todo(tab, i);
753 if (!todo)
754 return -1;
756 todo->next = *next;
757 *next = todo;
759 if (isl_tab_rollback(tab, snap) < 0)
760 return -1;
763 return 0;
766 /* Does the linked list contain a todo item that is the opposite of "todo".
767 * If so, return 1 and remove the opposite todo item.
769 static int has_opposite(struct isl_facet_todo *todo,
770 struct isl_facet_todo **list)
772 for (; *list; list = &(*list)->next) {
773 int eq;
774 eq = isl_basic_set_fast_is_equal(todo->bset, (*list)->bset);
775 if (eq < 0)
776 return -1;
777 if (!eq)
778 continue;
779 todo = *list;
780 *list = todo->next;
781 todo->next = NULL;
782 free_todo(todo);
783 return 1;
786 return 0;
789 /* Create todo items for all interior facets of the chamber represented
790 * by "tab" and collect them in first->next, taking care to cancel
791 * opposite todo items.
793 static int update_todo(struct isl_facet_todo *first, struct isl_tab *tab)
795 int i;
796 struct isl_tab_undo *snap;
797 struct isl_facet_todo *todo;
799 snap = isl_tab_snap(tab);
801 for (i = 0; i < tab->n_con; ++i) {
802 int drop;
804 if (tab->con[i].frozen)
805 continue;
806 if (tab->con[i].is_redundant)
807 continue;
809 if (isl_tab_select_facet(tab, i) < 0)
810 return -1;
812 todo = create_todo(tab, i);
813 if (!todo)
814 return -1;
816 drop = has_opposite(todo, &first->next);
817 if (drop < 0)
818 return -1;
820 if (drop)
821 free_todo(todo);
822 else {
823 todo->next = first->next;
824 first->next = todo;
827 if (isl_tab_rollback(tab, snap) < 0)
828 return -1;
831 return 0;
834 /* Compute the chamber decomposition of the parametric polytope respresented
835 * by "bset" given the parametric vertices and their activity domains.
837 * We are only interested in full-dimensional chambers.
838 * Each of these chambers is the intersection of the activity domains of
839 * one or more vertices and the union of all chambers is equal to the
840 * projection of the entire parametric polytope onto the parameter space.
842 * We first create an initial chamber by intersecting as many activity
843 * domains as possible without ending up with an empty or lower-dimensional
844 * set. As a minor optimization, we only consider those activity domains
845 * that contain some arbitrary point.
847 * For each of interior facets of the chamber, we construct a todo item,
848 * containing the facet and a constraint containing the other side of the facet,
849 * for constructing the chamber on the other side.
850 * While their are any todo items left, we pick a todo item and
851 * create the required chamber by intersecting all activity domains
852 * that contain the facet and have a full-dimensional intersection with
853 * the other side of the facet. For each of the interior facets, we
854 * again create todo items, taking care to cancel opposite todo items.
856 static __isl_give isl_vertices *compute_chambers(__isl_take isl_basic_set *bset,
857 __isl_take isl_vertices *vertices)
859 int i;
860 isl_vec *sample = NULL;
861 struct isl_tab *tab = NULL;
862 struct isl_tab_undo *snap;
863 unsigned nvar;
864 int *selection = NULL;
865 int n_chambers = 0;
866 struct isl_chamber_list *list = NULL;
867 struct isl_facet_todo *todo = NULL;
869 if (!bset || !vertices)
870 goto error;
872 selection = isl_alloc_array(vertices->ctx, int, vertices->n_vertices);
873 if (!selection)
874 goto error;
876 nvar = isl_basic_set_dim(bset, isl_dim_set);
877 bset = isl_basic_set_project_out(bset, isl_dim_set, 0, nvar);
879 tab = isl_tab_from_basic_set(bset);
880 for (i = 0; i < bset->n_ineq; ++i)
881 if (isl_tab_freeze_constraint(tab, i) < 0)
882 goto error;
883 if (isl_tab_track_bset(tab, bset) < 0)
884 goto error;
886 snap = isl_tab_snap(tab);
888 sample = isl_tab_get_sample_value(tab);
890 for (i = 0; i < vertices->n_vertices; ++i) {
891 selection[i] = isl_basic_set_contains(vertices->v[i].dom, sample);
892 if (selection[i] < 0)
893 goto error;
894 if (!selection[i])
895 continue;
896 selection[i] = can_intersect(tab, vertices->v[i].dom);
897 if (selection[i] < 0)
898 goto error;
901 if (isl_tab_detect_redundant(tab) < 0)
902 goto error;
904 if (add_chamber(&list, vertices, tab, selection) < 0)
905 goto error;
906 n_chambers++;
908 if (init_todo(&todo, tab) < 0)
909 goto error;
911 while (todo) {
912 struct isl_facet_todo *next;
914 if (isl_tab_rollback(tab, snap) < 0)
915 goto error;
917 if (isl_tab_add_ineq(tab, todo->constraint->el) < 0)
918 goto error;
919 if (isl_tab_freeze_constraint(tab, tab->n_con - 1) < 0)
920 goto error;
922 for (i = 0; i < vertices->n_vertices; ++i) {
923 selection[i] = bset_covers_tab(vertices->v[i].dom,
924 todo->tab);
925 if (selection[i] < 0)
926 goto error;
927 if (!selection[i])
928 continue;
929 selection[i] = can_intersect(tab, vertices->v[i].dom);
930 if (selection[i] < 0)
931 goto error;
934 if (isl_tab_detect_redundant(tab) < 0)
935 goto error;
937 if (add_chamber(&list, vertices, tab, selection) < 0)
938 goto error;
939 n_chambers++;
941 if (update_todo(todo, tab) < 0)
942 goto error;
944 next = todo->next;
945 todo->next = NULL;
946 free_todo(todo);
947 todo = next;
950 isl_vec_free(sample);
952 isl_tab_free(tab);
953 free(selection);
955 vertices = vertices_add_chambers(vertices, n_chambers, list);
957 for (i = 0; vertices && i < vertices->n_vertices; ++i) {
958 isl_basic_set_free(vertices->v[i].dom);
959 vertices->v[i].dom = NULL;
962 return vertices;
963 error:
964 free_chamber_list(list);
965 free_todo(todo);
966 isl_vec_free(sample);
967 isl_tab_free(tab);
968 free(selection);
969 if (!tab)
970 isl_basic_set_free(bset);
971 isl_vertices_free(vertices);
972 return NULL;
975 isl_ctx *isl_vertex_get_ctx(__isl_keep isl_vertex *vertex)
977 return vertex ? vertex->vertices->ctx : NULL;
980 int isl_vertex_get_id(__isl_keep isl_vertex *vertex)
982 return vertex ? vertex->id : -1;
985 __isl_give isl_basic_set *isl_vertex_get_domain(__isl_keep isl_vertex *vertex)
987 struct isl_vertex *v;
989 if (!vertex)
990 return NULL;
992 v = &vertex->vertices->v[vertex->id];
993 if (!v->dom) {
994 unsigned nvar;
995 nvar = isl_basic_set_dim(v->vertex, isl_dim_set);
996 v->dom = isl_basic_set_copy(v->vertex);
997 v->dom = isl_basic_set_project_out(v->dom, isl_dim_set, 0, nvar);
1000 return isl_basic_set_copy(v->dom);
1003 __isl_give isl_basic_set *isl_vertex_get_expr(__isl_keep isl_vertex *vertex)
1005 struct isl_vertex *v;
1007 if (!vertex)
1008 return NULL;
1010 v = &vertex->vertices->v[vertex->id];
1012 return isl_basic_set_copy(v->vertex);
1015 static __isl_give isl_vertex *isl_vertex_alloc(__isl_take isl_vertices *vertices,
1016 int id)
1018 isl_vertex *vertex;
1020 if (!vertices)
1021 return NULL;
1023 vertex = isl_alloc_type(vertices->ctx, isl_vertex);
1024 if (!vertex)
1025 goto error;
1027 vertex->vertices = vertices;
1028 vertex->id = id;
1030 return vertex;
1031 error:
1032 isl_vertices_free(vertices);
1033 return NULL;
1036 void isl_vertex_free(__isl_take isl_vertex *vertex)
1038 if (!vertex)
1039 return;
1040 isl_vertices_free(vertex->vertices);
1041 free(vertex);
1044 __isl_give isl_basic_set *isl_basic_set_set_integral(__isl_take isl_basic_set *bset)
1046 if (!bset)
1047 return NULL;
1049 if (!ISL_F_ISSET(bset, ISL_BASIC_MAP_RATIONAL))
1050 return bset;
1052 bset = isl_basic_set_cow(bset);
1053 if (!bset)
1054 return NULL;
1056 ISL_F_CLR(bset, ISL_BASIC_MAP_RATIONAL);
1058 return isl_basic_set_finalize(bset);
1061 isl_ctx *isl_cell_get_ctx(__isl_keep isl_cell *cell)
1063 return cell ? cell->vertices->ctx : NULL;
1066 __isl_give isl_basic_set *isl_cell_get_domain(__isl_keep isl_cell *cell)
1068 struct isl_chamber *c;
1070 if (!cell)
1071 return NULL;
1073 c = &cell->vertices->c[cell->id];
1075 return isl_basic_set_copy(c->dom);
1078 static __isl_give isl_cell *isl_cell_alloc(__isl_take isl_vertices *vertices,
1079 __isl_take isl_basic_set *dom, int id)
1081 isl_cell *cell;
1083 if (!vertices || !dom)
1084 goto error;
1086 cell = isl_alloc_type(dom->ctx, isl_cell);
1087 if (!cell)
1088 goto error;
1090 cell->vertices = vertices;
1091 cell->dom = dom;
1092 cell->id = id;
1094 return cell;
1095 error:
1096 isl_vertices_free(vertices);
1097 isl_basic_set_free(dom);
1098 return NULL;
1101 void isl_cell_free(__isl_take isl_cell *cell)
1103 if (!cell)
1104 return;
1106 isl_vertices_free(cell->vertices);
1107 isl_basic_set_free(cell->dom);
1108 free(cell);
1111 /* Create a tableau of the cone obtained by first homogenizing the given
1112 * polytope and then making all inequalities strict by setting the
1113 * constant term to -1.
1115 static struct isl_tab *tab_for_shifted_cone(__isl_keep isl_basic_set *bset)
1117 int i;
1118 isl_vec *c = NULL;
1119 struct isl_tab *tab;
1121 if (!bset)
1122 return NULL;
1123 tab = isl_tab_alloc(bset->ctx, bset->n_ineq + 1,
1124 1 + isl_basic_set_total_dim(bset), 0);
1125 if (!tab)
1126 return NULL;
1127 tab->rational = ISL_F_ISSET(bset, ISL_BASIC_SET_RATIONAL);
1128 if (ISL_F_ISSET(bset, ISL_BASIC_MAP_EMPTY)) {
1129 if (isl_tab_mark_empty(tab) < 0)
1130 goto error;
1131 return tab;
1134 c = isl_vec_alloc(bset->ctx, 1 + 1 + isl_basic_set_total_dim(bset));
1135 if (!c)
1136 goto error;
1138 isl_int_set_si(c->el[0], 0);
1139 for (i = 0; i < bset->n_eq; ++i) {
1140 isl_seq_cpy(c->el + 1, bset->eq[i], c->size - 1);
1141 if (isl_tab_add_eq(tab, c->el) < 0)
1142 goto error;
1145 isl_int_set_si(c->el[0], -1);
1146 for (i = 0; i < bset->n_ineq; ++i) {
1147 isl_seq_cpy(c->el + 1, bset->ineq[i], c->size - 1);
1148 if (isl_tab_add_ineq(tab, c->el) < 0)
1149 goto error;
1150 if (tab->empty) {
1151 isl_vec_free(c);
1152 return tab;
1156 isl_seq_clr(c->el + 1, c->size - 1);
1157 isl_int_set_si(c->el[1], 1);
1158 if (isl_tab_add_ineq(tab, c->el) < 0)
1159 goto error;
1161 isl_vec_free(c);
1162 return tab;
1163 error:
1164 isl_vec_free(c);
1165 isl_tab_free(tab);
1166 return NULL;
1169 /* Compute an interior point of "bset" by selecting an interior
1170 * point in homogeneous space and projecting the point back down.
1172 static __isl_give isl_vec *isl_basic_set_interior_point(
1173 __isl_keep isl_basic_set *bset)
1175 isl_vec *vec;
1176 struct isl_tab *tab;
1178 tab = tab_for_shifted_cone(bset);
1179 vec = isl_tab_get_sample_value(tab);
1180 isl_tab_free(tab);
1181 if (!vec)
1182 return NULL;
1184 isl_seq_cpy(vec->el, vec->el + 1, vec->size - 1);
1185 vec->size--;
1187 return vec;
1190 /* Call "fn" on all chambers of the parametric polytope with the shared
1191 * facets of neighboring chambers only appearing in one of the chambers.
1193 * We pick an interior point from one of the chambers and then make
1194 * all constraints that do not satisfy this point strict.
1196 int isl_vertices_foreach_disjoint_cell(__isl_keep isl_vertices *vertices,
1197 int (*fn)(__isl_take isl_cell *cell, void *user), void *user)
1199 int i, j;
1200 isl_vec *vec;
1201 isl_int v;
1202 isl_cell *cell;
1204 if (!vertices)
1205 return -1;
1207 if (vertices->n_chambers == 0)
1208 return 0;
1210 if (vertices->n_chambers == 1) {
1211 isl_basic_set *dom = isl_basic_set_copy(vertices->c[0].dom);
1212 dom = isl_basic_set_set_integral(dom);
1213 cell = isl_cell_alloc(isl_vertices_copy(vertices), dom, 0);
1214 if (!cell)
1215 return -1;
1216 return fn(cell, user);
1219 vec = isl_basic_set_interior_point(vertices->c[0].dom);
1220 if (!vec)
1221 return -1;
1223 isl_int_init(v);
1225 for (i = 0; i < vertices->n_chambers; ++i) {
1226 int r;
1227 isl_basic_set *dom = isl_basic_set_copy(vertices->c[i].dom);
1228 dom = isl_basic_set_cow(dom);
1229 if (!dom)
1230 goto error;
1231 for (j = 0; i && j < dom->n_ineq; ++j) {
1232 isl_seq_inner_product(vec->el, dom->ineq[j], vec->size,
1233 &v);
1234 if (!isl_int_is_neg(v))
1235 continue;
1236 isl_int_sub_ui(dom->ineq[j][0], dom->ineq[j][0], 1);
1238 dom = isl_basic_set_set_integral(dom);
1239 cell = isl_cell_alloc(isl_vertices_copy(vertices), dom, i);
1240 if (!cell)
1241 goto error;
1242 r = fn(cell, user);
1243 if (r < 0)
1244 goto error;
1247 isl_int_clear(v);
1248 isl_vec_free(vec);
1250 return 0;
1251 error:
1252 isl_int_clear(v);
1253 isl_vec_free(vec);
1254 return -1;
1257 int isl_vertices_foreach_cell(__isl_keep isl_vertices *vertices,
1258 int (*fn)(__isl_take isl_cell *cell, void *user), void *user)
1260 int i;
1261 isl_cell *cell;
1263 if (!vertices)
1264 return -1;
1266 if (vertices->n_chambers == 0)
1267 return 0;
1269 for (i = 0; i < vertices->n_chambers; ++i) {
1270 int r;
1271 isl_basic_set *dom = isl_basic_set_copy(vertices->c[i].dom);
1273 cell = isl_cell_alloc(isl_vertices_copy(vertices), dom, i);
1274 if (!cell)
1275 return -1;
1277 r = fn(cell, user);
1278 if (r < 0)
1279 return -1;
1282 return 0;
1285 int isl_vertices_foreach_vertex(__isl_keep isl_vertices *vertices,
1286 int (*fn)(__isl_take isl_vertex *vertex, void *user), void *user)
1288 int i;
1289 isl_vertex *vertex;
1291 if (!vertices)
1292 return -1;
1294 if (vertices->n_vertices == 0)
1295 return 0;
1297 for (i = 0; i < vertices->n_vertices; ++i) {
1298 int r;
1300 vertex = isl_vertex_alloc(isl_vertices_copy(vertices), i);
1301 if (!vertex)
1302 return -1;
1304 r = fn(vertex, user);
1305 if (r < 0)
1306 return -1;
1309 return 0;
1312 int isl_cell_foreach_vertex(__isl_keep isl_cell *cell,
1313 int (*fn)(__isl_take isl_vertex *vertex, void *user), void *user)
1315 int i;
1316 isl_vertex *vertex;
1317 struct isl_chamber *c;
1319 if (!cell)
1320 return -1;
1322 c = &cell->vertices->c[cell->id];
1324 if (c->n_vertices == 0)
1325 return 0;
1327 for (i = 0; i < c->n_vertices; ++i) {
1328 int r;
1330 vertex = isl_vertex_alloc(isl_vertices_copy(cell->vertices),
1331 c->vertices[i]);
1332 if (!vertex)
1333 return -1;
1335 r = fn(vertex, user);
1336 if (r < 0)
1337 return -1;
1340 return 0;
1343 isl_ctx *isl_vertices_get_ctx(__isl_keep isl_vertices *vertices)
1345 return vertices ? vertices->ctx : NULL;
1348 int isl_vertices_get_n_vertices(__isl_keep isl_vertices *vertices)
1350 return vertices ? vertices->n_vertices : -1;
1353 __isl_give isl_vertices *isl_morph_vertices(__isl_take isl_morph *morph,
1354 __isl_take isl_vertices *vertices)
1356 int i;
1357 isl_morph *param_morph = NULL;
1359 if (!morph || !vertices)
1360 goto error;
1362 isl_assert(vertices->ctx, vertices->ref == 1, goto error);
1364 param_morph = isl_morph_copy(morph);
1365 param_morph = isl_morph_remove_dom_dims(param_morph, isl_dim_set,
1366 0, isl_morph_dom_dim(morph, isl_dim_set));
1367 param_morph = isl_morph_remove_ran_dims(param_morph, isl_dim_set,
1368 0, isl_morph_ran_dim(morph, isl_dim_set));
1370 for (i = 0; i < vertices->n_vertices; ++i) {
1371 vertices->v[i].dom = isl_morph_basic_set(
1372 isl_morph_copy(param_morph), vertices->v[i].dom);
1373 vertices->v[i].vertex = isl_morph_basic_set(
1374 isl_morph_copy(morph), vertices->v[i].vertex);
1375 if (!vertices->v[i].vertex)
1376 goto error;
1379 for (i = 0; i < vertices->n_chambers; ++i) {
1380 vertices->c[i].dom = isl_morph_basic_set(
1381 isl_morph_copy(param_morph), vertices->c[i].dom);
1382 if (!vertices->c[i].dom)
1383 goto error;
1386 isl_morph_free(param_morph);
1387 isl_morph_free(morph);
1388 return vertices;
1389 error:
1390 isl_morph_free(param_morph);
1391 isl_morph_free(morph);
1392 isl_vertices_free(vertices);
1393 return NULL;