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,
14 #include <isl_map_private.h>
15 #include <isl_dim_private.h>
16 #include <isl_morph.h>
17 #include <isl_vertices_private.h>
18 #include <isl_mat_private.h>
24 static __isl_give isl_vertices
*compute_chambers(__isl_take isl_basic_set
*bset
,
25 __isl_take isl_vertices
*vertices
);
27 __isl_give isl_vertices
*isl_vertices_copy(__isl_keep isl_vertices
*vertices
)
36 void isl_vertices_free(__isl_take isl_vertices
*vertices
)
43 if (--vertices
->ref
> 0)
46 for (i
= 0; i
< vertices
->n_vertices
; ++i
) {
47 isl_basic_set_free(vertices
->v
[i
].vertex
);
48 isl_basic_set_free(vertices
->v
[i
].dom
);
52 for (i
= 0; i
< vertices
->n_chambers
; ++i
) {
53 free(vertices
->c
[i
].vertices
);
54 isl_basic_set_free(vertices
->c
[i
].dom
);
58 isl_basic_set_free(vertices
->bset
);
62 struct isl_vertex_list
{
64 struct isl_vertex_list
*next
;
67 static void free_vertex_list(struct isl_vertex_list
*list
)
69 struct isl_vertex_list
*next
;
71 for (; list
; list
= next
) {
73 isl_basic_set_free(list
->v
.vertex
);
74 isl_basic_set_free(list
->v
.dom
);
79 static __isl_give isl_vertices
*vertices_from_list(__isl_keep isl_basic_set
*bset
,
80 int n_vertices
, struct isl_vertex_list
*list
)
83 struct isl_vertex_list
*next
;
84 isl_vertices
*vertices
;
86 vertices
= isl_calloc_type(bset
->ctx
, isl_vertices
);
90 vertices
->bset
= isl_basic_set_copy(bset
);
91 vertices
->v
= isl_alloc_array(bset
->ctx
, struct isl_vertex
, n_vertices
);
94 vertices
->n_vertices
= n_vertices
;
96 for (i
= 0; list
; list
= next
, i
++) {
98 vertices
->v
[i
] = list
->v
;
105 free_vertex_list(list
);
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
)
116 struct isl_vertex_list
*v
= NULL
;
118 if (isl_tab_detect_implicit_equalities(tab
) < 0)
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
);
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
);
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
);
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
;
161 nparam
= isl_basic_set_dim(bset
, isl_dim_param
);
163 vertices
= isl_calloc_type(bset
->ctx
, isl_vertices
);
166 vertices
->bset
= isl_basic_set_copy(bset
);
169 vertices
->n_vertices
= 0;
170 vertices
->n_chambers
= 0;
175 /* Compute the parametric vertices and the chamber decomposition
176 * of the parametric polytope defined using the same constraints
177 * as "bset" in the 0D case.
178 * There is exactly one 0D vertex and a single chamber containing
181 static __isl_give isl_vertices
*vertices_0D(__isl_keep isl_basic_set
*bset
)
183 isl_vertices
*vertices
;
189 nparam
= isl_basic_set_dim(bset
, isl_dim_param
);
191 vertices
= isl_calloc_type(bset
->ctx
, isl_vertices
);
195 vertices
->bset
= isl_basic_set_copy(bset
);
197 vertices
->v
= isl_calloc_array(bset
->ctx
, struct isl_vertex
, 1);
200 vertices
->n_vertices
= 1;
201 vertices
->v
[0].vertex
= isl_basic_set_copy(bset
);
202 if (!vertices
->v
[0].vertex
)
205 vertices
->c
= isl_calloc_array(bset
->ctx
, struct isl_chamber
, 1);
208 vertices
->n_chambers
= 1;
209 vertices
->c
[0].n_vertices
= 1;
210 vertices
->c
[0].vertices
= isl_calloc_array(bset
->ctx
, int, 1);
211 if (!vertices
->c
[0].vertices
)
213 vertices
->c
[0].dom
= isl_basic_set_copy(bset
);
214 if (!vertices
->c
[0].dom
)
219 isl_vertices_free(vertices
);
223 static int isl_mat_rank(__isl_keep isl_mat
*mat
)
228 H
= isl_mat_left_hermite(isl_mat_copy(mat
), 0, NULL
, NULL
);
232 for (col
= 0; col
< H
->n_col
; ++col
) {
233 for (row
= 0; row
< H
->n_row
; ++row
)
234 if (!isl_int_is_zero(H
->row
[row
][col
]))
245 /* Is the row pointed to by "f" linearly independent of the "n" first
248 static int is_independent(__isl_keep isl_mat
*facets
, int n
, isl_int
*f
)
252 if (isl_seq_first_non_zero(f
, facets
->n_col
) < 0)
255 isl_seq_cpy(facets
->row
[n
], f
, facets
->n_col
);
256 facets
->n_row
= n
+ 1;
257 rank
= isl_mat_rank(facets
);
261 return rank
== n
+ 1;
264 /* Check whether we can select constraint "level", given the current selection
265 * reflected by facets in "tab", the rows of "facets" and the earlier
266 * "selected" elements of "selection".
268 * If the constraint is (strictly) redundant in the tableau, selecting it would
269 * result in an empty tableau, so it can't be selected.
270 * If the set variable part of the constraint is not linearly indepedent
271 * of the set variable parts of the already selected constraints,
272 * the constraint cannot be selected.
273 * If selecting the constraint results in an empty tableau, the constraint
274 * cannot be selected.
275 * Finally, if selecting the constraint results in some explicitly
276 * deselected constraints turning into equalities, then the corresponding
277 * vertices have already been generated, so the constraint cannot be selected.
279 static int can_select(__isl_keep isl_basic_set
*bset
, int level
,
280 struct isl_tab
*tab
, __isl_keep isl_mat
*facets
, int selected
,
286 struct isl_tab_undo
*snap
;
288 if (isl_tab_is_redundant(tab
, level
))
291 ovar
= isl_dim_offset(bset
->dim
, isl_dim_set
);
293 indep
= is_independent(facets
, selected
, bset
->ineq
[level
] + 1 + ovar
);
299 snap
= isl_tab_snap(tab
);
300 if (isl_tab_select_facet(tab
, level
) < 0)
304 if (isl_tab_rollback(tab
, snap
) < 0)
309 for (i
= 0; i
< level
; ++i
) {
312 if (selection
[i
] != DESELECTED
)
315 if (isl_tab_is_equality(tab
, i
))
317 else if (isl_tab_is_redundant(tab
, i
))
320 sgn
= isl_tab_sign_of_max(tab
, i
);
324 if (isl_tab_rollback(tab
, snap
) < 0)
333 /* Compute the parametric vertices and the chamber decomposition
334 * of a parametric polytope that is not full-dimensional.
336 * Simply map the parametric polytope to a lower dimensional space
337 * and map the resulting vertices back.
339 static __isl_give isl_vertices
*lower_dim_vertices(
340 __isl_keep isl_basic_set
*bset
)
343 isl_vertices
*vertices
;
345 bset
= isl_basic_set_copy(bset
);
346 morph
= isl_basic_set_full_compression(bset
);
347 bset
= isl_morph_basic_set(isl_morph_copy(morph
), bset
);
349 vertices
= isl_basic_set_compute_vertices(bset
);
350 isl_basic_set_free(bset
);
352 morph
= isl_morph_inverse(morph
);
354 vertices
= isl_morph_vertices(morph
, vertices
);
359 /* Compute the parametric vertices and the chamber decomposition
360 * of the parametric polytope defined using the same constraints
361 * as "bset". "bset" is assumed to have no existentially quantified
364 * The vertices themselves are computed in a fairly simplistic way.
365 * We simply run through all combinations of d constraints,
366 * with d the number of set variables, and check if those d constraints
367 * define a vertex. To avoid the generation of duplicate vertices,
368 * which we may happen if a vertex is defined by more that d constraints,
369 * we make sure we only generate the vertex for the d constraints with
372 * We set up a tableau and keep track of which facets have been
373 * selected. The tableau is marked strict_redundant so that we can be
374 * sure that any constraint that is marked redundant (and that is not
375 * also marked zero) is not an equality.
376 * If a constraint is marked DESELECTED, it means the constraint was
377 * SELECTED before (in combination with the same selection of earlier
378 * constraints). If such a deselected constraint turns out to be an
379 * equality, then any vertex that may still be found with the current
380 * selection has already been generated when the constraint was selected.
381 * A constraint is marked UNSELECTED when there is no way selecting
382 * the constraint could lead to a vertex (in combination with the current
383 * selection of earlier constraints).
385 * The set variable coefficients of the selected constraints are stored
386 * in the facets matrix.
388 __isl_give isl_vertices
*isl_basic_set_compute_vertices(
389 __isl_keep isl_basic_set
*bset
)
397 struct isl_tab_undo
**snap
;
399 struct isl_vertex_list
*list
= NULL
;
401 isl_vertices
*vertices
;
406 if (isl_basic_set_fast_is_empty(bset
))
407 return vertices_empty(bset
);
410 return lower_dim_vertices(bset
);
412 isl_assert(bset
->ctx
, isl_basic_set_dim(bset
, isl_dim_div
) == 0,
415 if (isl_basic_set_dim(bset
, isl_dim_set
) == 0)
416 return vertices_0D(bset
);
418 nvar
= isl_basic_set_dim(bset
, isl_dim_set
);
420 bset
= isl_basic_set_copy(bset
);
421 bset
= isl_basic_set_set_rational(bset
);
425 tab
= isl_tab_from_basic_set(bset
);
428 tab
->strict_redundant
= 1;
431 vertices
= vertices_empty(bset
);
432 isl_basic_set_free(bset
);
437 selection
= isl_alloc_array(bset
->ctx
, int, bset
->n_ineq
);
438 snap
= isl_alloc_array(bset
->ctx
, struct isl_tab_undo
*, bset
->n_ineq
);
439 facets
= isl_mat_alloc(bset
->ctx
, nvar
, nvar
);
440 if (!selection
|| !snap
|| !facets
)
448 if (level
>= bset
->n_ineq
||
449 (!init
&& selection
[level
] != SELECTED
)) {
456 snap
[level
] = isl_tab_snap(tab
);
457 ok
= can_select(bset
, level
, tab
, facets
, selected
,
462 selection
[level
] = SELECTED
;
465 selection
[level
] = UNSELECTED
;
467 selection
[level
] = DESELECTED
;
469 if (isl_tab_rollback(tab
, snap
[level
]) < 0)
472 if (selected
== nvar
) {
473 if (tab
->n_dead
== nvar
) {
474 if (add_vertex(&list
, bset
, tab
) < 0)
485 isl_mat_free(facets
);
491 vertices
= vertices_from_list(bset
, n_vertices
, list
);
493 vertices
= compute_chambers(bset
, vertices
);
497 isl_mat_free(facets
);
501 isl_basic_set_free(bset
);
505 struct isl_chamber_list
{
506 struct isl_chamber c
;
507 struct isl_chamber_list
*next
;
510 static void free_chamber_list(struct isl_chamber_list
*list
)
512 struct isl_chamber_list
*next
;
514 for (; list
; list
= next
) {
516 isl_basic_set_free(list
->c
.dom
);
517 free(list
->c
.vertices
);
522 /* Check whether the basic set "bset" is a superset of the basic set described
523 * by "tab", i.e., check whether all constraints of "bset" are redundant.
525 static int bset_covers_tab(__isl_keep isl_basic_set
*bset
, struct isl_tab
*tab
)
532 for (i
= 0; i
< bset
->n_ineq
; ++i
) {
533 enum isl_ineq_type type
= isl_tab_ineq_type(tab
, bset
->ineq
[i
]);
535 case isl_ineq_error
: return -1;
536 case isl_ineq_redundant
: continue;
544 static __isl_give isl_vertices
*vertices_add_chambers(
545 __isl_take isl_vertices
*vertices
, int n_chambers
,
546 struct isl_chamber_list
*list
)
549 struct isl_chamber_list
*next
;
551 vertices
->c
= isl_alloc_array(vertices
->ctx
, struct isl_chamber
, n_chambers
);
554 vertices
->n_chambers
= n_chambers
;
556 for (i
= 0; list
; list
= next
, i
++) {
558 vertices
->c
[i
] = list
->c
;
564 isl_vertices_free(vertices
);
565 free_chamber_list(list
);
569 /* Can "tab" be intersected with "bset" without resulting in
570 * a lower-dimensional set.
572 static int can_intersect(struct isl_tab
*tab
, __isl_keep isl_basic_set
*bset
)
575 struct isl_tab_undo
*snap
;
577 if (isl_tab_extend_cons(tab
, bset
->n_ineq
) < 0)
580 snap
= isl_tab_snap(tab
);
582 for (i
= 0; i
< bset
->n_ineq
; ++i
) {
583 if (isl_tab_ineq_type(tab
, bset
->ineq
[i
]) == isl_ineq_redundant
)
585 if (isl_tab_add_ineq(tab
, bset
->ineq
[i
]) < 0)
589 if (isl_tab_detect_implicit_equalities(tab
) < 0)
592 if (isl_tab_rollback(tab
, snap
) < 0)
600 static int add_chamber(struct isl_chamber_list
**list
,
601 __isl_keep isl_vertices
*vertices
, struct isl_tab
*tab
, int *selection
)
606 struct isl_tab_undo
*snap
;
607 struct isl_chamber_list
*c
= NULL
;
609 for (i
= 0; i
< vertices
->n_vertices
; ++i
)
613 snap
= isl_tab_snap(tab
);
615 for (i
= 0; i
< tab
->n_con
&& tab
->con
[i
].frozen
; ++i
)
616 tab
->con
[i
].frozen
= 0;
619 if (isl_tab_detect_redundant(tab
) < 0)
622 c
= isl_calloc_type(tab
->mat
->ctx
, struct isl_chamber_list
);
625 c
->c
.vertices
= isl_alloc_array(tab
->mat
->ctx
, int, n_vertices
);
628 c
->c
.dom
= isl_basic_set_from_basic_map(isl_basic_map_copy(tab
->bmap
));
629 c
->c
.dom
= isl_basic_set_set_rational(c
->c
.dom
);
630 c
->c
.dom
= isl_basic_set_cow(c
->c
.dom
);
631 c
->c
.dom
= isl_basic_set_update_from_tab(c
->c
.dom
, tab
);
632 c
->c
.dom
= isl_basic_set_simplify(c
->c
.dom
);
633 c
->c
.dom
= isl_basic_set_finalize(c
->c
.dom
);
637 c
->c
.n_vertices
= n_vertices
;
639 for (i
= 0, j
= 0; i
< vertices
->n_vertices
; ++i
)
641 c
->c
.vertices
[j
] = i
;
648 for (i
= 0; i
< n_frozen
; ++i
)
649 tab
->con
[i
].frozen
= 1;
651 if (isl_tab_rollback(tab
, snap
) < 0)
656 free_chamber_list(c
);
660 struct isl_facet_todo
{
661 struct isl_tab
*tab
; /* A tableau representation of the facet */
662 isl_basic_set
*bset
; /* A normalized basic set representation */
663 isl_vec
*constraint
; /* Constraint pointing to the other side */
664 struct isl_facet_todo
*next
;
667 static void free_todo(struct isl_facet_todo
*todo
)
670 struct isl_facet_todo
*next
= todo
->next
;
672 isl_tab_free(todo
->tab
);
673 isl_basic_set_free(todo
->bset
);
674 isl_vec_free(todo
->constraint
);
681 static struct isl_facet_todo
*create_todo(struct isl_tab
*tab
, int con
)
685 struct isl_tab_undo
*snap
;
686 struct isl_facet_todo
*todo
;
688 snap
= isl_tab_snap(tab
);
690 for (i
= 0; i
< tab
->n_con
&& tab
->con
[i
].frozen
; ++i
)
691 tab
->con
[i
].frozen
= 0;
694 if (isl_tab_detect_redundant(tab
) < 0)
697 todo
= isl_calloc_type(tab
->mat
->ctx
, struct isl_facet_todo
);
701 todo
->constraint
= isl_vec_alloc(tab
->mat
->ctx
, 1 + tab
->n_var
);
702 if (!todo
->constraint
)
704 isl_seq_neg(todo
->constraint
->el
, tab
->bmap
->ineq
[con
], 1 + tab
->n_var
);
705 todo
->bset
= isl_basic_set_from_basic_map(isl_basic_map_copy(tab
->bmap
));
706 todo
->bset
= isl_basic_set_set_rational(todo
->bset
);
707 todo
->bset
= isl_basic_set_cow(todo
->bset
);
708 todo
->bset
= isl_basic_set_update_from_tab(todo
->bset
, tab
);
709 todo
->bset
= isl_basic_set_simplify(todo
->bset
);
710 todo
->bset
= isl_basic_set_sort_constraints(todo
->bset
);
713 ISL_F_SET(todo
->bset
, ISL_BASIC_SET_NORMALIZED
);
714 todo
->tab
= isl_tab_dup(tab
);
718 for (i
= 0; i
< n_frozen
; ++i
)
719 tab
->con
[i
].frozen
= 1;
721 if (isl_tab_rollback(tab
, snap
) < 0)
730 /* Create todo items for all interior facets of the chamber represented
731 * by "tab" and collect them in "next".
733 static int init_todo(struct isl_facet_todo
**next
, struct isl_tab
*tab
)
736 struct isl_tab_undo
*snap
;
737 struct isl_facet_todo
*todo
;
739 snap
= isl_tab_snap(tab
);
741 for (i
= 0; i
< tab
->n_con
; ++i
) {
742 if (tab
->con
[i
].frozen
)
744 if (tab
->con
[i
].is_redundant
)
747 if (isl_tab_select_facet(tab
, i
) < 0)
750 todo
= create_todo(tab
, i
);
757 if (isl_tab_rollback(tab
, snap
) < 0)
764 /* Does the linked list contain a todo item that is the opposite of "todo".
765 * If so, return 1 and remove the opposite todo item.
767 static int has_opposite(struct isl_facet_todo
*todo
,
768 struct isl_facet_todo
**list
)
770 for (; *list
; list
= &(*list
)->next
) {
772 eq
= isl_basic_set_fast_is_equal(todo
->bset
, (*list
)->bset
);
787 /* Create todo items for all interior facets of the chamber represented
788 * by "tab" and collect them in first->next, taking care to cancel
789 * opposite todo items.
791 static int update_todo(struct isl_facet_todo
*first
, struct isl_tab
*tab
)
794 struct isl_tab_undo
*snap
;
795 struct isl_facet_todo
*todo
;
797 snap
= isl_tab_snap(tab
);
799 for (i
= 0; i
< tab
->n_con
; ++i
) {
802 if (tab
->con
[i
].frozen
)
804 if (tab
->con
[i
].is_redundant
)
807 if (isl_tab_select_facet(tab
, i
) < 0)
810 todo
= create_todo(tab
, i
);
814 drop
= has_opposite(todo
, &first
->next
);
821 todo
->next
= first
->next
;
825 if (isl_tab_rollback(tab
, snap
) < 0)
832 /* Compute the chamber decomposition of the parametric polytope respresented
833 * by "bset" given the parametric vertices and their activity domains.
835 * We are only interested in full-dimensional chambers.
836 * Each of these chambers is the intersection of the activity domains of
837 * one or more vertices and the union of all chambers is equal to the
838 * projection of the entire parametric polytope onto the parameter space.
840 * We first create an initial chamber by intersecting as many activity
841 * domains as possible without ending up with an empty or lower-dimensional
842 * set. As a minor optimization, we only consider those activity domains
843 * that contain some arbitrary point.
845 * For each of interior facets of the chamber, we construct a todo item,
846 * containing the facet and a constraint containing the other side of the facet,
847 * for constructing the chamber on the other side.
848 * While their are any todo items left, we pick a todo item and
849 * create the required chamber by intersecting all activity domains
850 * that contain the facet and have a full-dimensional intersection with
851 * the other side of the facet. For each of the interior facets, we
852 * again create todo items, taking care to cancel opposite todo items.
854 static __isl_give isl_vertices
*compute_chambers(__isl_take isl_basic_set
*bset
,
855 __isl_take isl_vertices
*vertices
)
858 isl_vec
*sample
= NULL
;
859 struct isl_tab
*tab
= NULL
;
860 struct isl_tab_undo
*snap
;
862 int *selection
= NULL
;
864 struct isl_chamber_list
*list
= NULL
;
865 struct isl_facet_todo
*todo
= NULL
;
867 if (!bset
|| !vertices
)
870 selection
= isl_alloc_array(vertices
->ctx
, int, vertices
->n_vertices
);
874 nvar
= isl_basic_set_dim(bset
, isl_dim_set
);
875 bset
= isl_basic_set_project_out(bset
, isl_dim_set
, 0, nvar
);
877 tab
= isl_tab_from_basic_set(bset
);
878 for (i
= 0; i
< bset
->n_ineq
; ++i
)
879 if (isl_tab_freeze_constraint(tab
, i
) < 0)
881 if (isl_tab_track_bset(tab
, bset
) < 0)
884 snap
= isl_tab_snap(tab
);
886 sample
= isl_tab_get_sample_value(tab
);
888 for (i
= 0; i
< vertices
->n_vertices
; ++i
) {
889 selection
[i
] = isl_basic_set_contains(vertices
->v
[i
].dom
, sample
);
890 if (selection
[i
] < 0)
894 selection
[i
] = can_intersect(tab
, vertices
->v
[i
].dom
);
895 if (selection
[i
] < 0)
899 if (isl_tab_detect_redundant(tab
) < 0)
902 if (add_chamber(&list
, vertices
, tab
, selection
) < 0)
906 if (init_todo(&todo
, tab
) < 0)
910 struct isl_facet_todo
*next
;
912 if (isl_tab_rollback(tab
, snap
) < 0)
915 if (isl_tab_add_ineq(tab
, todo
->constraint
->el
) < 0)
917 if (isl_tab_freeze_constraint(tab
, tab
->n_con
- 1) < 0)
920 for (i
= 0; i
< vertices
->n_vertices
; ++i
) {
921 selection
[i
] = bset_covers_tab(vertices
->v
[i
].dom
,
923 if (selection
[i
] < 0)
927 selection
[i
] = can_intersect(tab
, vertices
->v
[i
].dom
);
928 if (selection
[i
] < 0)
932 if (isl_tab_detect_redundant(tab
) < 0)
935 if (add_chamber(&list
, vertices
, tab
, selection
) < 0)
939 if (update_todo(todo
, tab
) < 0)
948 isl_vec_free(sample
);
953 vertices
= vertices_add_chambers(vertices
, n_chambers
, list
);
955 for (i
= 0; vertices
&& i
< vertices
->n_vertices
; ++i
) {
956 isl_basic_set_free(vertices
->v
[i
].dom
);
957 vertices
->v
[i
].dom
= NULL
;
962 free_chamber_list(list
);
964 isl_vec_free(sample
);
968 isl_basic_set_free(bset
);
969 isl_vertices_free(vertices
);
973 isl_ctx
*isl_vertex_get_ctx(__isl_keep isl_vertex
*vertex
)
975 return vertex
? isl_vertices_get_ctx(vertex
->vertices
) : NULL
;
978 int isl_vertex_get_id(__isl_keep isl_vertex
*vertex
)
980 return vertex
? vertex
->id
: -1;
983 __isl_give isl_basic_set
*isl_vertex_get_domain(__isl_keep isl_vertex
*vertex
)
985 struct isl_vertex
*v
;
990 v
= &vertex
->vertices
->v
[vertex
->id
];
993 nvar
= isl_basic_set_dim(v
->vertex
, isl_dim_set
);
994 v
->dom
= isl_basic_set_copy(v
->vertex
);
995 v
->dom
= isl_basic_set_project_out(v
->dom
, isl_dim_set
, 0, nvar
);
998 return isl_basic_set_copy(v
->dom
);
1001 __isl_give isl_basic_set
*isl_vertex_get_expr(__isl_keep isl_vertex
*vertex
)
1003 struct isl_vertex
*v
;
1008 v
= &vertex
->vertices
->v
[vertex
->id
];
1010 return isl_basic_set_copy(v
->vertex
);
1013 static __isl_give isl_vertex
*isl_vertex_alloc(__isl_take isl_vertices
*vertices
,
1021 vertex
= isl_alloc_type(vertices
->ctx
, isl_vertex
);
1025 vertex
->vertices
= vertices
;
1030 isl_vertices_free(vertices
);
1034 void isl_vertex_free(__isl_take isl_vertex
*vertex
)
1038 isl_vertices_free(vertex
->vertices
);
1042 __isl_give isl_basic_set
*isl_basic_set_set_integral(__isl_take isl_basic_set
*bset
)
1047 if (!ISL_F_ISSET(bset
, ISL_BASIC_MAP_RATIONAL
))
1050 bset
= isl_basic_set_cow(bset
);
1054 ISL_F_CLR(bset
, ISL_BASIC_MAP_RATIONAL
);
1056 return isl_basic_set_finalize(bset
);
1059 isl_ctx
*isl_cell_get_ctx(__isl_keep isl_cell
*cell
)
1061 return cell
? cell
->dom
->ctx
: NULL
;
1064 __isl_give isl_basic_set
*isl_cell_get_domain(__isl_keep isl_cell
*cell
)
1066 return cell
? isl_basic_set_copy(cell
->dom
) : NULL
;
1069 static __isl_give isl_cell
*isl_cell_alloc(__isl_take isl_vertices
*vertices
,
1070 __isl_take isl_basic_set
*dom
, int id
)
1073 isl_cell
*cell
= NULL
;
1075 if (!vertices
|| !dom
)
1078 cell
= isl_calloc_type(dom
->ctx
, isl_cell
);
1082 cell
->n_vertices
= vertices
->c
[id
].n_vertices
;
1083 cell
->ids
= isl_alloc_array(dom
->ctx
, int, cell
->n_vertices
);
1086 for (i
= 0; i
< cell
->n_vertices
; ++i
)
1087 cell
->ids
[i
] = vertices
->c
[id
].vertices
[i
];
1088 cell
->vertices
= vertices
;
1093 isl_cell_free(cell
);
1094 isl_vertices_free(vertices
);
1095 isl_basic_set_free(dom
);
1099 void isl_cell_free(__isl_take isl_cell
*cell
)
1104 isl_vertices_free(cell
->vertices
);
1106 isl_basic_set_free(cell
->dom
);
1110 /* Create a tableau of the cone obtained by first homogenizing the given
1111 * polytope and then making all inequalities strict by setting the
1112 * constant term to -1.
1114 static struct isl_tab
*tab_for_shifted_cone(__isl_keep isl_basic_set
*bset
)
1118 struct isl_tab
*tab
;
1122 tab
= isl_tab_alloc(bset
->ctx
, bset
->n_ineq
+ 1,
1123 1 + isl_basic_set_total_dim(bset
), 0);
1126 tab
->rational
= ISL_F_ISSET(bset
, ISL_BASIC_SET_RATIONAL
);
1127 if (ISL_F_ISSET(bset
, ISL_BASIC_MAP_EMPTY
)) {
1128 if (isl_tab_mark_empty(tab
) < 0)
1133 c
= isl_vec_alloc(bset
->ctx
, 1 + 1 + isl_basic_set_total_dim(bset
));
1137 isl_int_set_si(c
->el
[0], 0);
1138 for (i
= 0; i
< bset
->n_eq
; ++i
) {
1139 isl_seq_cpy(c
->el
+ 1, bset
->eq
[i
], c
->size
- 1);
1140 if (isl_tab_add_eq(tab
, c
->el
) < 0)
1144 isl_int_set_si(c
->el
[0], -1);
1145 for (i
= 0; i
< bset
->n_ineq
; ++i
) {
1146 isl_seq_cpy(c
->el
+ 1, bset
->ineq
[i
], c
->size
- 1);
1147 if (isl_tab_add_ineq(tab
, c
->el
) < 0)
1155 isl_seq_clr(c
->el
+ 1, c
->size
- 1);
1156 isl_int_set_si(c
->el
[1], 1);
1157 if (isl_tab_add_ineq(tab
, c
->el
) < 0)
1168 /* Compute an interior point of "bset" by selecting an interior
1169 * point in homogeneous space and projecting the point back down.
1171 static __isl_give isl_vec
*isl_basic_set_interior_point(
1172 __isl_keep isl_basic_set
*bset
)
1175 struct isl_tab
*tab
;
1177 tab
= tab_for_shifted_cone(bset
);
1178 vec
= isl_tab_get_sample_value(tab
);
1183 isl_seq_cpy(vec
->el
, vec
->el
+ 1, vec
->size
- 1);
1189 /* Call "fn" on all chambers of the parametric polytope with the shared
1190 * facets of neighboring chambers only appearing in one of the chambers.
1192 * We pick an interior point from one of the chambers and then make
1193 * all constraints that do not satisfy this point strict.
1195 int isl_vertices_foreach_disjoint_cell(__isl_keep isl_vertices
*vertices
,
1196 int (*fn
)(__isl_take isl_cell
*cell
, void *user
), void *user
)
1206 if (vertices
->n_chambers
== 0)
1209 if (vertices
->n_chambers
== 1) {
1210 isl_basic_set
*dom
= isl_basic_set_copy(vertices
->c
[0].dom
);
1211 dom
= isl_basic_set_set_integral(dom
);
1212 cell
= isl_cell_alloc(isl_vertices_copy(vertices
), dom
, 0);
1215 return fn(cell
, user
);
1218 vec
= isl_basic_set_interior_point(vertices
->c
[0].dom
);
1224 for (i
= 0; i
< vertices
->n_chambers
; ++i
) {
1226 isl_basic_set
*dom
= isl_basic_set_copy(vertices
->c
[i
].dom
);
1227 dom
= isl_basic_set_cow(dom
);
1230 for (j
= 0; i
&& j
< dom
->n_ineq
; ++j
) {
1231 isl_seq_inner_product(vec
->el
, dom
->ineq
[j
], vec
->size
,
1233 if (!isl_int_is_neg(v
))
1235 isl_int_sub_ui(dom
->ineq
[j
][0], dom
->ineq
[j
][0], 1);
1237 dom
= isl_basic_set_set_integral(dom
);
1238 cell
= isl_cell_alloc(isl_vertices_copy(vertices
), dom
, i
);
1256 int isl_vertices_foreach_cell(__isl_keep isl_vertices
*vertices
,
1257 int (*fn
)(__isl_take isl_cell
*cell
, void *user
), void *user
)
1265 if (vertices
->n_chambers
== 0)
1268 for (i
= 0; i
< vertices
->n_chambers
; ++i
) {
1270 isl_basic_set
*dom
= isl_basic_set_copy(vertices
->c
[i
].dom
);
1272 cell
= isl_cell_alloc(isl_vertices_copy(vertices
), dom
, i
);
1284 int isl_vertices_foreach_vertex(__isl_keep isl_vertices
*vertices
,
1285 int (*fn
)(__isl_take isl_vertex
*vertex
, void *user
), void *user
)
1293 if (vertices
->n_vertices
== 0)
1296 for (i
= 0; i
< vertices
->n_vertices
; ++i
) {
1299 vertex
= isl_vertex_alloc(isl_vertices_copy(vertices
), i
);
1303 r
= fn(vertex
, user
);
1311 int isl_cell_foreach_vertex(__isl_keep isl_cell
*cell
,
1312 int (*fn
)(__isl_take isl_vertex
*vertex
, void *user
), void *user
)
1320 if (cell
->n_vertices
== 0)
1323 for (i
= 0; i
< cell
->n_vertices
; ++i
) {
1326 vertex
= isl_vertex_alloc(isl_vertices_copy(cell
->vertices
),
1331 r
= fn(vertex
, user
);
1339 isl_ctx
*isl_vertices_get_ctx(__isl_keep isl_vertices
*vertices
)
1341 return vertices
? vertices
->bset
->ctx
: NULL
;
1344 int isl_vertices_get_n_vertices(__isl_keep isl_vertices
*vertices
)
1346 return vertices
? vertices
->n_vertices
: -1;
1349 __isl_give isl_vertices
*isl_morph_vertices(__isl_take isl_morph
*morph
,
1350 __isl_take isl_vertices
*vertices
)
1353 isl_morph
*param_morph
= NULL
;
1355 if (!morph
|| !vertices
)
1358 isl_assert(vertices
->bset
->ctx
, vertices
->ref
== 1, goto error
);
1360 param_morph
= isl_morph_copy(morph
);
1361 param_morph
= isl_morph_remove_dom_dims(param_morph
, isl_dim_set
,
1362 0, isl_morph_dom_dim(morph
, isl_dim_set
));
1363 param_morph
= isl_morph_remove_ran_dims(param_morph
, isl_dim_set
,
1364 0, isl_morph_ran_dim(morph
, isl_dim_set
));
1366 for (i
= 0; i
< vertices
->n_vertices
; ++i
) {
1367 vertices
->v
[i
].dom
= isl_morph_basic_set(
1368 isl_morph_copy(param_morph
), vertices
->v
[i
].dom
);
1369 vertices
->v
[i
].vertex
= isl_morph_basic_set(
1370 isl_morph_copy(morph
), vertices
->v
[i
].vertex
);
1371 if (!vertices
->v
[i
].vertex
)
1375 for (i
= 0; i
< vertices
->n_chambers
; ++i
) {
1376 vertices
->c
[i
].dom
= isl_morph_basic_set(
1377 isl_morph_copy(param_morph
), vertices
->c
[i
].dom
);
1378 if (!vertices
->c
[i
].dom
)
1382 isl_morph_free(param_morph
);
1383 isl_morph_free(morph
);
1386 isl_morph_free(param_morph
);
1387 isl_morph_free(morph
);
1388 isl_vertices_free(vertices
);
1392 /* Construct a simplex isl_cell spanned by the vertices with indices in
1393 * "simplex_ids" and "other_ids" and call "fn" on this isl_cell.
1395 static int call_on_simplex(__isl_keep isl_cell
*cell
,
1396 int *simplex_ids
, int n_simplex
, int *other_ids
, int n_other
,
1397 int (*fn
)(__isl_take isl_cell
*simplex
, void *user
), void *user
)
1401 struct isl_cell
*simplex
;
1403 ctx
= isl_cell_get_ctx(cell
);
1405 simplex
= isl_calloc_type(ctx
, struct isl_cell
);
1408 simplex
->vertices
= isl_vertices_copy(cell
->vertices
);
1409 if (!simplex
->vertices
)
1411 simplex
->dom
= isl_basic_set_copy(cell
->dom
);
1414 simplex
->n_vertices
= n_simplex
+ n_other
;
1415 simplex
->ids
= isl_alloc_array(ctx
, int, simplex
->n_vertices
);
1419 for (i
= 0; i
< n_simplex
; ++i
)
1420 simplex
->ids
[i
] = simplex_ids
[i
];
1421 for (i
= 0; i
< n_other
; ++i
)
1422 simplex
->ids
[n_simplex
+ i
] = other_ids
[i
];
1424 return fn(simplex
, user
);
1426 isl_cell_free(simplex
);
1430 /* Check whether the parametric vertex described by "vertex"
1431 * lies on the facet corresponding to constraint "facet" of "bset".
1432 * The isl_vec "v" is a temporary vector than can be used by this function.
1434 * We eliminate the variables from the facet constraint using the
1435 * equalities defining the vertex and check if the result is identical
1438 * It would probably be better to keep track of the constraints defining
1439 * a vertex during the vertex construction so that we could simply look
1442 static int vertex_on_facet(__isl_keep isl_basic_set
*vertex
,
1443 __isl_keep isl_basic_set
*bset
, int facet
, __isl_keep isl_vec
*v
)
1448 isl_seq_cpy(v
->el
, bset
->ineq
[facet
], v
->size
);
1451 for (i
= 0; i
< vertex
->n_eq
; ++i
) {
1452 int k
= isl_seq_last_non_zero(vertex
->eq
[i
], v
->size
);
1453 isl_seq_elim(v
->el
, vertex
->eq
[i
], k
, v
->size
, &m
);
1457 return isl_seq_first_non_zero(v
->el
, v
->size
) == -1;
1460 /* Triangulate the polytope spanned by the vertices with ids
1461 * in "simplex_ids" and "other_ids" and call "fn" on each of
1462 * the resulting simplices.
1463 * If the input polytope is already a simplex, we simply call "fn".
1464 * Otherwise, we pick a point from "other_ids" and add it to "simplex_ids".
1465 * Then we consider each facet of "bset" that does not contain the point
1466 * we just picked, but does contain some of the other points in "other_ids"
1467 * and call ourselves recursively on the polytope spanned by the new
1468 * "simplex_ids" and those points in "other_ids" that lie on the facet.
1470 static int triangulate(__isl_keep isl_cell
*cell
, __isl_keep isl_vec
*v
,
1471 int *simplex_ids
, int n_simplex
, int *other_ids
, int n_other
,
1472 int (*fn
)(__isl_take isl_cell
*simplex
, void *user
), void *user
)
1477 isl_basic_set
*vertex
;
1478 isl_basic_set
*bset
;
1480 d
= isl_basic_set_dim(cell
->vertices
->bset
, isl_dim_set
);
1481 nparam
= isl_basic_set_dim(cell
->vertices
->bset
, isl_dim_param
);
1483 if (n_simplex
+ n_other
== d
+ 1)
1484 return call_on_simplex(cell
, simplex_ids
, n_simplex
,
1485 other_ids
, n_other
, fn
, user
);
1487 simplex_ids
[n_simplex
] = other_ids
[0];
1488 vertex
= cell
->vertices
->v
[other_ids
[0]].vertex
;
1489 bset
= cell
->vertices
->bset
;
1491 ids
= isl_alloc_array(ctx
, int, n_other
- 1);
1492 for (i
= 0; i
< bset
->n_ineq
; ++i
) {
1493 if (isl_seq_first_non_zero(bset
->ineq
[i
] + 1 + nparam
, d
) == -1)
1495 if (vertex_on_facet(vertex
, bset
, i
, v
))
1498 for (j
= 1, k
= 0; j
< n_other
; ++j
) {
1500 ov
= cell
->vertices
->v
[other_ids
[j
]].vertex
;
1501 if (vertex_on_facet(ov
, bset
, i
, v
))
1502 ids
[k
++] = other_ids
[j
];
1507 if (triangulate(cell
, v
, simplex_ids
, n_simplex
+ 1,
1508 ids
, k
, fn
, user
) < 0)
1519 /* Triangulate the given cell and call "fn" on each of the resulting
1522 int isl_cell_foreach_simplex(__isl_take isl_cell
*cell
,
1523 int (*fn
)(__isl_take isl_cell
*simplex
, void *user
), void *user
)
1529 int *simplex_ids
= NULL
;
1534 d
= isl_basic_set_dim(cell
->vertices
->bset
, isl_dim_set
);
1535 total
= isl_basic_set_total_dim(cell
->vertices
->bset
);
1537 if (cell
->n_vertices
== d
+ 1)
1538 return fn(cell
, user
);
1540 ctx
= isl_cell_get_ctx(cell
);
1541 simplex_ids
= isl_alloc_array(ctx
, int, d
+ 1);
1545 v
= isl_vec_alloc(ctx
, 1 + total
);
1549 r
= triangulate(cell
, v
, simplex_ids
, 0,
1550 cell
->ids
, cell
->n_vertices
, fn
, user
);
1555 isl_cell_free(cell
);
1561 isl_cell_free(cell
);