Line data Source code
1 : /*
2 : * Copyright 2010 INRIA Saclay
3 : *
4 : * Use of this software is governed by the MIT license
5 : *
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 : */
10 :
11 : #include <isl_map_private.h>
12 : #include <isl_aff_private.h>
13 : #include <isl/set.h>
14 : #include <isl_seq.h>
15 : #include <isl_tab.h>
16 : #include <isl_space_private.h>
17 : #include <isl_morph.h>
18 : #include <isl_vertices_private.h>
19 : #include <isl_mat_private.h>
20 : #include <isl_vec_private.h>
21 :
22 : #define SELECTED 1
23 : #define DESELECTED -1
24 : #define UNSELECTED 0
25 :
26 : static __isl_give isl_vertices *compute_chambers(__isl_take isl_basic_set *bset,
27 : __isl_take isl_vertices *vertices);
28 :
29 0 : __isl_give isl_vertices *isl_vertices_copy(__isl_keep isl_vertices *vertices)
30 : {
31 0 : if (!vertices)
32 0 : return NULL;
33 :
34 0 : vertices->ref++;
35 0 : return vertices;
36 : }
37 :
38 0 : __isl_null isl_vertices *isl_vertices_free(__isl_take isl_vertices *vertices)
39 : {
40 : int i;
41 :
42 0 : if (!vertices)
43 0 : return NULL;
44 :
45 0 : if (--vertices->ref > 0)
46 0 : return NULL;
47 :
48 0 : for (i = 0; i < vertices->n_vertices; ++i) {
49 0 : isl_basic_set_free(vertices->v[i].vertex);
50 0 : isl_basic_set_free(vertices->v[i].dom);
51 : }
52 0 : free(vertices->v);
53 :
54 0 : for (i = 0; i < vertices->n_chambers; ++i) {
55 0 : free(vertices->c[i].vertices);
56 0 : isl_basic_set_free(vertices->c[i].dom);
57 : }
58 0 : free(vertices->c);
59 :
60 0 : isl_basic_set_free(vertices->bset);
61 0 : free(vertices);
62 :
63 0 : return NULL;
64 : }
65 :
66 : struct isl_vertex_list {
67 : struct isl_vertex v;
68 : struct isl_vertex_list *next;
69 : };
70 :
71 0 : static struct isl_vertex_list *free_vertex_list(struct isl_vertex_list *list)
72 : {
73 : struct isl_vertex_list *next;
74 :
75 0 : for (; list; list = next) {
76 0 : next = list->next;
77 0 : isl_basic_set_free(list->v.vertex);
78 0 : isl_basic_set_free(list->v.dom);
79 0 : free(list);
80 : }
81 :
82 0 : return NULL;
83 : }
84 :
85 0 : static __isl_give isl_vertices *vertices_from_list(__isl_keep isl_basic_set *bset,
86 : int n_vertices, struct isl_vertex_list *list)
87 : {
88 : int i;
89 : struct isl_vertex_list *next;
90 : isl_vertices *vertices;
91 :
92 0 : vertices = isl_calloc_type(bset->ctx, isl_vertices);
93 0 : if (!vertices)
94 0 : goto error;
95 0 : vertices->ref = 1;
96 0 : vertices->bset = isl_basic_set_copy(bset);
97 0 : vertices->v = isl_alloc_array(bset->ctx, struct isl_vertex, n_vertices);
98 0 : if (n_vertices && !vertices->v)
99 0 : goto error;
100 0 : vertices->n_vertices = n_vertices;
101 :
102 0 : for (i = 0; list; list = next, i++) {
103 0 : next = list->next;
104 0 : vertices->v[i] = list->v;
105 0 : free(list);
106 : }
107 :
108 0 : return vertices;
109 : error:
110 0 : isl_vertices_free(vertices);
111 0 : free_vertex_list(list);
112 0 : return NULL;
113 : }
114 :
115 : /* Prepend a vertex to the linked list "list" based on the equalities in "tab".
116 : * Return isl_bool_true if the vertex was actually added and
117 : * isl_bool_false otherwise.
118 : * In particular, vertices with a lower-dimensional activity domain are
119 : * not added to the list because they would not be included in any chamber.
120 : * Return isl_bool_error on error.
121 : */
122 0 : static isl_bool add_vertex(struct isl_vertex_list **list,
123 : __isl_keep isl_basic_set *bset, struct isl_tab *tab)
124 : {
125 : unsigned nvar;
126 0 : struct isl_vertex_list *v = NULL;
127 :
128 0 : if (isl_tab_detect_implicit_equalities(tab) < 0)
129 0 : return isl_bool_error;
130 :
131 0 : nvar = isl_basic_set_dim(bset, isl_dim_set);
132 :
133 0 : v = isl_calloc_type(tab->mat->ctx, struct isl_vertex_list);
134 0 : if (!v)
135 0 : goto error;
136 :
137 0 : v->v.vertex = isl_basic_set_copy(bset);
138 0 : v->v.vertex = isl_basic_set_cow(v->v.vertex);
139 0 : v->v.vertex = isl_basic_set_update_from_tab(v->v.vertex, tab);
140 0 : v->v.vertex = isl_basic_set_simplify(v->v.vertex);
141 0 : v->v.vertex = isl_basic_set_finalize(v->v.vertex);
142 0 : if (!v->v.vertex)
143 0 : goto error;
144 0 : isl_assert(bset->ctx, v->v.vertex->n_eq >= nvar, goto error);
145 0 : v->v.dom = isl_basic_set_copy(v->v.vertex);
146 0 : v->v.dom = isl_basic_set_params(v->v.dom);
147 0 : if (!v->v.dom)
148 0 : goto error;
149 :
150 0 : if (v->v.dom->n_eq > 0) {
151 0 : free_vertex_list(v);
152 0 : return isl_bool_false;
153 : }
154 :
155 0 : v->next = *list;
156 0 : *list = v;
157 :
158 0 : return isl_bool_true;
159 : error:
160 0 : free_vertex_list(v);
161 0 : return isl_bool_error;
162 : }
163 :
164 : /* Compute the parametric vertices and the chamber decomposition
165 : * of an empty parametric polytope.
166 : */
167 0 : static __isl_give isl_vertices *vertices_empty(__isl_keep isl_basic_set *bset)
168 : {
169 : isl_vertices *vertices;
170 :
171 0 : if (!bset)
172 0 : return NULL;
173 :
174 0 : vertices = isl_calloc_type(bset->ctx, isl_vertices);
175 0 : if (!vertices)
176 0 : return NULL;
177 0 : vertices->bset = isl_basic_set_copy(bset);
178 0 : vertices->ref = 1;
179 :
180 0 : vertices->n_vertices = 0;
181 0 : vertices->n_chambers = 0;
182 :
183 0 : return vertices;
184 : }
185 :
186 : /* Compute the parametric vertices and the chamber decomposition
187 : * of the parametric polytope defined using the same constraints
188 : * as "bset" in the 0D case.
189 : * There is exactly one 0D vertex and a single chamber containing
190 : * the vertex.
191 : */
192 0 : static __isl_give isl_vertices *vertices_0D(__isl_keep isl_basic_set *bset)
193 : {
194 : isl_vertices *vertices;
195 :
196 0 : if (!bset)
197 0 : return NULL;
198 :
199 0 : vertices = isl_calloc_type(bset->ctx, isl_vertices);
200 0 : if (!vertices)
201 0 : return NULL;
202 0 : vertices->ref = 1;
203 0 : vertices->bset = isl_basic_set_copy(bset);
204 :
205 0 : vertices->v = isl_calloc_array(bset->ctx, struct isl_vertex, 1);
206 0 : if (!vertices->v)
207 0 : goto error;
208 0 : vertices->n_vertices = 1;
209 0 : vertices->v[0].vertex = isl_basic_set_copy(bset);
210 0 : vertices->v[0].dom = isl_basic_set_params(isl_basic_set_copy(bset));
211 0 : if (!vertices->v[0].vertex || !vertices->v[0].dom)
212 : goto error;
213 :
214 0 : vertices->c = isl_calloc_array(bset->ctx, struct isl_chamber, 1);
215 0 : if (!vertices->c)
216 0 : goto error;
217 0 : vertices->n_chambers = 1;
218 0 : vertices->c[0].n_vertices = 1;
219 0 : vertices->c[0].vertices = isl_calloc_array(bset->ctx, int, 1);
220 0 : if (!vertices->c[0].vertices)
221 0 : goto error;
222 0 : vertices->c[0].dom = isl_basic_set_copy(vertices->v[0].dom);
223 0 : if (!vertices->c[0].dom)
224 0 : goto error;
225 :
226 0 : return vertices;
227 : error:
228 0 : isl_vertices_free(vertices);
229 0 : return NULL;
230 : }
231 :
232 : /* Is the row pointed to by "f" linearly independent of the "n" first
233 : * rows in "facets"?
234 : */
235 0 : static int is_independent(__isl_keep isl_mat *facets, int n, isl_int *f)
236 : {
237 : int rank;
238 :
239 0 : if (isl_seq_first_non_zero(f, facets->n_col) < 0)
240 0 : return 0;
241 :
242 0 : isl_seq_cpy(facets->row[n], f, facets->n_col);
243 0 : facets->n_row = n + 1;
244 0 : rank = isl_mat_rank(facets);
245 0 : if (rank < 0)
246 0 : return -1;
247 :
248 0 : return rank == n + 1;
249 : }
250 :
251 : /* Check whether we can select constraint "level", given the current selection
252 : * reflected by facets in "tab", the rows of "facets" and the earlier
253 : * "selected" elements of "selection".
254 : *
255 : * If the constraint is (strictly) redundant in the tableau, selecting it would
256 : * result in an empty tableau, so it can't be selected.
257 : * If the set variable part of the constraint is not linearly independent
258 : * of the set variable parts of the already selected constraints,
259 : * the constraint cannot be selected.
260 : * If selecting the constraint results in an empty tableau, the constraint
261 : * cannot be selected.
262 : * Finally, if selecting the constraint results in some explicitly
263 : * deselected constraints turning into equalities, then the corresponding
264 : * vertices have already been generated, so the constraint cannot be selected.
265 : */
266 0 : static int can_select(__isl_keep isl_basic_set *bset, int level,
267 : struct isl_tab *tab, __isl_keep isl_mat *facets, int selected,
268 : int *selection)
269 : {
270 : int i;
271 : int indep;
272 : unsigned ovar;
273 : struct isl_tab_undo *snap;
274 :
275 0 : if (isl_tab_is_redundant(tab, level))
276 0 : return 0;
277 :
278 0 : ovar = isl_space_offset(bset->dim, isl_dim_set);
279 :
280 0 : indep = is_independent(facets, selected, bset->ineq[level] + 1 + ovar);
281 0 : if (indep < 0)
282 0 : return -1;
283 0 : if (!indep)
284 0 : return 0;
285 :
286 0 : snap = isl_tab_snap(tab);
287 0 : if (isl_tab_select_facet(tab, level) < 0)
288 0 : return -1;
289 :
290 0 : if (tab->empty) {
291 0 : if (isl_tab_rollback(tab, snap) < 0)
292 0 : return -1;
293 0 : return 0;
294 : }
295 :
296 0 : for (i = 0; i < level; ++i) {
297 : int sgn;
298 :
299 0 : if (selection[i] != DESELECTED)
300 0 : continue;
301 :
302 0 : if (isl_tab_is_equality(tab, i))
303 0 : sgn = 0;
304 0 : else if (isl_tab_is_redundant(tab, i))
305 0 : sgn = 1;
306 : else
307 0 : sgn = isl_tab_sign_of_max(tab, i);
308 0 : if (sgn < -1)
309 0 : return -1;
310 0 : if (sgn <= 0) {
311 0 : if (isl_tab_rollback(tab, snap) < 0)
312 0 : return -1;
313 0 : return 0;
314 : }
315 : }
316 :
317 0 : return 1;
318 : }
319 :
320 : /* Compute the parametric vertices and the chamber decomposition
321 : * of a parametric polytope that is not full-dimensional.
322 : *
323 : * Simply map the parametric polytope to a lower dimensional space
324 : * and map the resulting vertices back.
325 : */
326 0 : static __isl_give isl_vertices *lower_dim_vertices(
327 : __isl_keep isl_basic_set *bset)
328 : {
329 : isl_morph *morph;
330 : isl_vertices *vertices;
331 :
332 0 : bset = isl_basic_set_copy(bset);
333 0 : morph = isl_basic_set_full_compression(bset);
334 0 : bset = isl_morph_basic_set(isl_morph_copy(morph), bset);
335 :
336 0 : vertices = isl_basic_set_compute_vertices(bset);
337 0 : isl_basic_set_free(bset);
338 :
339 0 : morph = isl_morph_inverse(morph);
340 :
341 0 : vertices = isl_morph_vertices(morph, vertices);
342 :
343 0 : return vertices;
344 : }
345 :
346 : /* Compute the parametric vertices and the chamber decomposition
347 : * of the parametric polytope defined using the same constraints
348 : * as "bset". "bset" is assumed to have no existentially quantified
349 : * variables.
350 : *
351 : * The vertices themselves are computed in a fairly simplistic way.
352 : * We simply run through all combinations of d constraints,
353 : * with d the number of set variables, and check if those d constraints
354 : * define a vertex. To avoid the generation of duplicate vertices,
355 : * which we may happen if a vertex is defined by more that d constraints,
356 : * we make sure we only generate the vertex for the d constraints with
357 : * smallest index.
358 : *
359 : * We set up a tableau and keep track of which facets have been
360 : * selected. The tableau is marked strict_redundant so that we can be
361 : * sure that any constraint that is marked redundant (and that is not
362 : * also marked zero) is not an equality.
363 : * If a constraint is marked DESELECTED, it means the constraint was
364 : * SELECTED before (in combination with the same selection of earlier
365 : * constraints). If such a deselected constraint turns out to be an
366 : * equality, then any vertex that may still be found with the current
367 : * selection has already been generated when the constraint was selected.
368 : * A constraint is marked UNSELECTED when there is no way selecting
369 : * the constraint could lead to a vertex (in combination with the current
370 : * selection of earlier constraints).
371 : *
372 : * The set variable coefficients of the selected constraints are stored
373 : * in the facets matrix.
374 : */
375 0 : __isl_give isl_vertices *isl_basic_set_compute_vertices(
376 : __isl_keep isl_basic_set *bset)
377 : {
378 : struct isl_tab *tab;
379 : int level;
380 : int init;
381 : unsigned nvar;
382 0 : int *selection = NULL;
383 : int selected;
384 0 : struct isl_tab_undo **snap = NULL;
385 0 : isl_mat *facets = NULL;
386 0 : struct isl_vertex_list *list = NULL;
387 0 : int n_vertices = 0;
388 : isl_vertices *vertices;
389 :
390 0 : if (!bset)
391 0 : return NULL;
392 :
393 0 : if (isl_basic_set_plain_is_empty(bset))
394 0 : return vertices_empty(bset);
395 :
396 0 : if (bset->n_eq != 0)
397 0 : return lower_dim_vertices(bset);
398 :
399 0 : isl_assert(bset->ctx, isl_basic_set_dim(bset, isl_dim_div) == 0,
400 : return NULL);
401 :
402 0 : if (isl_basic_set_dim(bset, isl_dim_set) == 0)
403 0 : return vertices_0D(bset);
404 :
405 0 : nvar = isl_basic_set_dim(bset, isl_dim_set);
406 :
407 0 : bset = isl_basic_set_copy(bset);
408 0 : bset = isl_basic_set_set_rational(bset);
409 0 : if (!bset)
410 0 : return NULL;
411 :
412 0 : tab = isl_tab_from_basic_set(bset, 0);
413 0 : if (!tab)
414 0 : goto error;
415 0 : tab->strict_redundant = 1;
416 :
417 0 : if (tab->empty) {
418 0 : vertices = vertices_empty(bset);
419 0 : isl_basic_set_free(bset);
420 0 : isl_tab_free(tab);
421 0 : return vertices;
422 : }
423 :
424 0 : selection = isl_alloc_array(bset->ctx, int, bset->n_ineq);
425 0 : snap = isl_alloc_array(bset->ctx, struct isl_tab_undo *, bset->n_ineq);
426 0 : facets = isl_mat_alloc(bset->ctx, nvar, nvar);
427 0 : if ((bset->n_ineq && (!selection || !snap)) || !facets)
428 : goto error;
429 :
430 0 : level = 0;
431 0 : init = 1;
432 0 : selected = 0;
433 :
434 0 : while (level >= 0) {
435 0 : if (level >= bset->n_ineq ||
436 0 : (!init && selection[level] != SELECTED)) {
437 0 : --level;
438 0 : init = 0;
439 0 : continue;
440 : }
441 0 : if (init) {
442 : int ok;
443 0 : snap[level] = isl_tab_snap(tab);
444 0 : ok = can_select(bset, level, tab, facets, selected,
445 : selection);
446 0 : if (ok < 0)
447 0 : goto error;
448 0 : if (ok) {
449 0 : selection[level] = SELECTED;
450 0 : selected++;
451 : } else
452 0 : selection[level] = UNSELECTED;
453 : } else {
454 0 : selection[level] = DESELECTED;
455 0 : selected--;
456 0 : if (isl_tab_rollback(tab, snap[level]) < 0)
457 0 : goto error;
458 : }
459 0 : if (selected == nvar) {
460 0 : if (tab->n_dead == nvar) {
461 0 : isl_bool added = add_vertex(&list, bset, tab);
462 0 : if (added < 0)
463 0 : goto error;
464 0 : if (added)
465 0 : n_vertices++;
466 : }
467 0 : init = 0;
468 0 : continue;
469 : }
470 0 : ++level;
471 0 : init = 1;
472 : }
473 :
474 0 : isl_mat_free(facets);
475 0 : free(selection);
476 0 : free(snap);
477 :
478 0 : isl_tab_free(tab);
479 :
480 0 : vertices = vertices_from_list(bset, n_vertices, list);
481 :
482 0 : vertices = compute_chambers(bset, vertices);
483 :
484 0 : return vertices;
485 : error:
486 0 : free_vertex_list(list);
487 0 : isl_mat_free(facets);
488 0 : free(selection);
489 0 : free(snap);
490 0 : isl_tab_free(tab);
491 0 : isl_basic_set_free(bset);
492 0 : return NULL;
493 : }
494 :
495 : struct isl_chamber_list {
496 : struct isl_chamber c;
497 : struct isl_chamber_list *next;
498 : };
499 :
500 0 : static void free_chamber_list(struct isl_chamber_list *list)
501 : {
502 : struct isl_chamber_list *next;
503 :
504 0 : for (; list; list = next) {
505 0 : next = list->next;
506 0 : isl_basic_set_free(list->c.dom);
507 0 : free(list->c.vertices);
508 0 : free(list);
509 : }
510 0 : }
511 :
512 : /* Check whether the basic set "bset" is a superset of the basic set described
513 : * by "tab", i.e., check whether all constraints of "bset" are redundant.
514 : */
515 0 : static isl_bool bset_covers_tab(__isl_keep isl_basic_set *bset,
516 : struct isl_tab *tab)
517 : {
518 : int i;
519 :
520 0 : if (!bset || !tab)
521 0 : return isl_bool_error;
522 :
523 0 : for (i = 0; i < bset->n_ineq; ++i) {
524 0 : enum isl_ineq_type type = isl_tab_ineq_type(tab, bset->ineq[i]);
525 0 : switch (type) {
526 0 : case isl_ineq_error: return isl_bool_error;
527 0 : case isl_ineq_redundant: continue;
528 0 : default: return isl_bool_false;
529 : }
530 : }
531 :
532 0 : return isl_bool_true;
533 : }
534 :
535 0 : static __isl_give isl_vertices *vertices_add_chambers(
536 : __isl_take isl_vertices *vertices, int n_chambers,
537 : struct isl_chamber_list *list)
538 : {
539 : int i;
540 : isl_ctx *ctx;
541 : struct isl_chamber_list *next;
542 :
543 0 : ctx = isl_vertices_get_ctx(vertices);
544 0 : vertices->c = isl_alloc_array(ctx, struct isl_chamber, n_chambers);
545 0 : if (!vertices->c)
546 0 : goto error;
547 0 : vertices->n_chambers = n_chambers;
548 :
549 0 : for (i = 0; list; list = next, i++) {
550 0 : next = list->next;
551 0 : vertices->c[i] = list->c;
552 0 : free(list);
553 : }
554 :
555 0 : return vertices;
556 : error:
557 0 : isl_vertices_free(vertices);
558 0 : free_chamber_list(list);
559 0 : return NULL;
560 : }
561 :
562 : /* Can "tab" be intersected with "bset" without resulting in
563 : * a lower-dimensional set.
564 : * "bset" itself is assumed to be full-dimensional.
565 : */
566 0 : static isl_bool can_intersect(struct isl_tab *tab,
567 : __isl_keep isl_basic_set *bset)
568 : {
569 : int i;
570 : struct isl_tab_undo *snap;
571 :
572 0 : if (bset->n_eq > 0)
573 0 : isl_die(isl_basic_set_get_ctx(bset), isl_error_internal,
574 : "expecting full-dimensional input",
575 : return isl_bool_error);
576 :
577 0 : if (isl_tab_extend_cons(tab, bset->n_ineq) < 0)
578 0 : return isl_bool_error;
579 :
580 0 : snap = isl_tab_snap(tab);
581 :
582 0 : for (i = 0; i < bset->n_ineq; ++i) {
583 0 : if (isl_tab_ineq_type(tab, bset->ineq[i]) == isl_ineq_redundant)
584 0 : continue;
585 0 : if (isl_tab_add_ineq(tab, bset->ineq[i]) < 0)
586 0 : return isl_bool_error;
587 : }
588 :
589 0 : if (isl_tab_detect_implicit_equalities(tab) < 0)
590 0 : return isl_bool_error;
591 0 : if (tab->n_dead) {
592 0 : if (isl_tab_rollback(tab, snap) < 0)
593 0 : return isl_bool_error;
594 0 : return isl_bool_false;
595 : }
596 :
597 0 : return isl_bool_true;
598 : }
599 :
600 0 : static int add_chamber(struct isl_chamber_list **list,
601 : __isl_keep isl_vertices *vertices, struct isl_tab *tab, int *selection)
602 : {
603 : int n_frozen;
604 : int i, j;
605 0 : int n_vertices = 0;
606 : struct isl_tab_undo *snap;
607 0 : struct isl_chamber_list *c = NULL;
608 :
609 0 : for (i = 0; i < vertices->n_vertices; ++i)
610 0 : if (selection[i])
611 0 : n_vertices++;
612 :
613 0 : snap = isl_tab_snap(tab);
614 :
615 0 : for (i = 0; i < tab->n_con && tab->con[i].frozen; ++i)
616 0 : tab->con[i].frozen = 0;
617 0 : n_frozen = i;
618 :
619 0 : if (isl_tab_detect_redundant(tab) < 0)
620 0 : return -1;
621 :
622 0 : c = isl_calloc_type(tab->mat->ctx, struct isl_chamber_list);
623 0 : if (!c)
624 0 : goto error;
625 0 : c->c.vertices = isl_alloc_array(tab->mat->ctx, int, n_vertices);
626 0 : if (n_vertices && !c->c.vertices)
627 0 : goto error;
628 0 : c->c.dom = isl_basic_set_copy(isl_tab_peek_bset(tab));
629 0 : c->c.dom = isl_basic_set_set_rational(c->c.dom);
630 0 : c->c.dom = isl_basic_set_cow(c->c.dom);
631 0 : c->c.dom = isl_basic_set_update_from_tab(c->c.dom, tab);
632 0 : c->c.dom = isl_basic_set_simplify(c->c.dom);
633 0 : c->c.dom = isl_basic_set_finalize(c->c.dom);
634 0 : if (!c->c.dom)
635 0 : goto error;
636 :
637 0 : c->c.n_vertices = n_vertices;
638 :
639 0 : for (i = 0, j = 0; i < vertices->n_vertices; ++i)
640 0 : if (selection[i]) {
641 0 : c->c.vertices[j] = i;
642 0 : j++;
643 : }
644 :
645 0 : c->next = *list;
646 0 : *list = c;
647 :
648 0 : for (i = 0; i < n_frozen; ++i)
649 0 : tab->con[i].frozen = 1;
650 :
651 0 : if (isl_tab_rollback(tab, snap) < 0)
652 0 : return -1;
653 :
654 0 : return 0;
655 : error:
656 0 : free_chamber_list(c);
657 0 : return -1;
658 : }
659 :
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;
665 : };
666 :
667 0 : static void free_todo(struct isl_facet_todo *todo)
668 : {
669 0 : while (todo) {
670 0 : struct isl_facet_todo *next = todo->next;
671 :
672 0 : isl_tab_free(todo->tab);
673 0 : isl_basic_set_free(todo->bset);
674 0 : isl_vec_free(todo->constraint);
675 0 : free(todo);
676 :
677 0 : todo = next;
678 : }
679 0 : }
680 :
681 0 : static struct isl_facet_todo *create_todo(struct isl_tab *tab, int con)
682 : {
683 : int i;
684 : int n_frozen;
685 : struct isl_tab_undo *snap;
686 : struct isl_facet_todo *todo;
687 :
688 0 : snap = isl_tab_snap(tab);
689 :
690 0 : for (i = 0; i < tab->n_con && tab->con[i].frozen; ++i)
691 0 : tab->con[i].frozen = 0;
692 0 : n_frozen = i;
693 :
694 0 : if (isl_tab_detect_redundant(tab) < 0)
695 0 : return NULL;
696 :
697 0 : todo = isl_calloc_type(tab->mat->ctx, struct isl_facet_todo);
698 0 : if (!todo)
699 0 : return NULL;
700 :
701 0 : todo->constraint = isl_vec_alloc(tab->mat->ctx, 1 + tab->n_var);
702 0 : if (!todo->constraint)
703 0 : goto error;
704 0 : isl_seq_neg(todo->constraint->el, tab->bmap->ineq[con], 1 + tab->n_var);
705 0 : todo->bset = isl_basic_set_copy(isl_tab_peek_bset(tab));
706 0 : todo->bset = isl_basic_set_set_rational(todo->bset);
707 0 : todo->bset = isl_basic_set_cow(todo->bset);
708 0 : todo->bset = isl_basic_set_update_from_tab(todo->bset, tab);
709 0 : todo->bset = isl_basic_set_simplify(todo->bset);
710 0 : todo->bset = isl_basic_set_sort_constraints(todo->bset);
711 0 : if (!todo->bset)
712 0 : goto error;
713 0 : ISL_F_SET(todo->bset, ISL_BASIC_SET_NORMALIZED);
714 0 : todo->tab = isl_tab_dup(tab);
715 0 : if (!todo->tab)
716 0 : goto error;
717 :
718 0 : for (i = 0; i < n_frozen; ++i)
719 0 : tab->con[i].frozen = 1;
720 :
721 0 : if (isl_tab_rollback(tab, snap) < 0)
722 0 : goto error;
723 :
724 0 : return todo;
725 : error:
726 0 : free_todo(todo);
727 0 : return NULL;
728 : }
729 :
730 : /* Create todo items for all interior facets of the chamber represented
731 : * by "tab" and collect them in "next".
732 : */
733 0 : static int init_todo(struct isl_facet_todo **next, struct isl_tab *tab)
734 : {
735 : int i;
736 : struct isl_tab_undo *snap;
737 : struct isl_facet_todo *todo;
738 :
739 0 : snap = isl_tab_snap(tab);
740 :
741 0 : for (i = 0; i < tab->n_con; ++i) {
742 0 : if (tab->con[i].frozen)
743 0 : continue;
744 0 : if (tab->con[i].is_redundant)
745 0 : continue;
746 :
747 0 : if (isl_tab_select_facet(tab, i) < 0)
748 0 : return -1;
749 :
750 0 : todo = create_todo(tab, i);
751 0 : if (!todo)
752 0 : return -1;
753 :
754 0 : todo->next = *next;
755 0 : *next = todo;
756 :
757 0 : if (isl_tab_rollback(tab, snap) < 0)
758 0 : return -1;
759 : }
760 :
761 0 : return 0;
762 : }
763 :
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.
766 : */
767 0 : static int has_opposite(struct isl_facet_todo *todo,
768 : struct isl_facet_todo **list)
769 : {
770 0 : for (; *list; list = &(*list)->next) {
771 : int eq;
772 0 : eq = isl_basic_set_plain_is_equal(todo->bset, (*list)->bset);
773 0 : if (eq < 0)
774 0 : return -1;
775 0 : if (!eq)
776 0 : continue;
777 0 : todo = *list;
778 0 : *list = todo->next;
779 0 : todo->next = NULL;
780 0 : free_todo(todo);
781 0 : return 1;
782 : }
783 :
784 0 : return 0;
785 : }
786 :
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.
790 : */
791 0 : static int update_todo(struct isl_facet_todo *first, struct isl_tab *tab)
792 : {
793 : int i;
794 : struct isl_tab_undo *snap;
795 : struct isl_facet_todo *todo;
796 :
797 0 : snap = isl_tab_snap(tab);
798 :
799 0 : for (i = 0; i < tab->n_con; ++i) {
800 : int drop;
801 :
802 0 : if (tab->con[i].frozen)
803 0 : continue;
804 0 : if (tab->con[i].is_redundant)
805 0 : continue;
806 :
807 0 : if (isl_tab_select_facet(tab, i) < 0)
808 0 : return -1;
809 :
810 0 : todo = create_todo(tab, i);
811 0 : if (!todo)
812 0 : return -1;
813 :
814 0 : drop = has_opposite(todo, &first->next);
815 0 : if (drop < 0)
816 0 : return -1;
817 :
818 0 : if (drop)
819 0 : free_todo(todo);
820 : else {
821 0 : todo->next = first->next;
822 0 : first->next = todo;
823 : }
824 :
825 0 : if (isl_tab_rollback(tab, snap) < 0)
826 0 : return -1;
827 : }
828 :
829 0 : return 0;
830 : }
831 :
832 : /* Compute the chamber decomposition of the parametric polytope respresented
833 : * by "bset" given the parametric vertices and their activity domains.
834 : *
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.
839 : *
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.
844 : *
845 : * For each of the 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.
853 : */
854 0 : static __isl_give isl_vertices *compute_chambers(__isl_take isl_basic_set *bset,
855 : __isl_take isl_vertices *vertices)
856 : {
857 : int i;
858 : isl_ctx *ctx;
859 0 : isl_vec *sample = NULL;
860 0 : struct isl_tab *tab = NULL;
861 : struct isl_tab_undo *snap;
862 0 : int *selection = NULL;
863 0 : int n_chambers = 0;
864 0 : struct isl_chamber_list *list = NULL;
865 0 : struct isl_facet_todo *todo = NULL;
866 :
867 0 : if (!bset || !vertices)
868 : goto error;
869 :
870 0 : ctx = isl_vertices_get_ctx(vertices);
871 0 : selection = isl_alloc_array(ctx, int, vertices->n_vertices);
872 0 : if (vertices->n_vertices && !selection)
873 0 : goto error;
874 :
875 0 : bset = isl_basic_set_params(bset);
876 :
877 0 : tab = isl_tab_from_basic_set(bset, 1);
878 0 : if (!tab)
879 0 : goto error;
880 0 : for (i = 0; i < bset->n_ineq; ++i)
881 0 : if (isl_tab_freeze_constraint(tab, i) < 0)
882 0 : goto error;
883 0 : isl_basic_set_free(bset);
884 :
885 0 : snap = isl_tab_snap(tab);
886 :
887 0 : sample = isl_tab_get_sample_value(tab);
888 :
889 0 : for (i = 0; i < vertices->n_vertices; ++i) {
890 0 : selection[i] = isl_basic_set_contains(vertices->v[i].dom, sample);
891 0 : if (selection[i] < 0)
892 0 : goto error;
893 0 : if (!selection[i])
894 0 : continue;
895 0 : selection[i] = can_intersect(tab, vertices->v[i].dom);
896 0 : if (selection[i] < 0)
897 0 : goto error;
898 : }
899 :
900 0 : if (isl_tab_detect_redundant(tab) < 0)
901 0 : goto error;
902 :
903 0 : if (add_chamber(&list, vertices, tab, selection) < 0)
904 0 : goto error;
905 0 : n_chambers++;
906 :
907 0 : if (init_todo(&todo, tab) < 0)
908 0 : goto error;
909 :
910 0 : while (todo) {
911 : struct isl_facet_todo *next;
912 :
913 0 : if (isl_tab_rollback(tab, snap) < 0)
914 0 : goto error;
915 :
916 0 : if (isl_tab_add_ineq(tab, todo->constraint->el) < 0)
917 0 : goto error;
918 0 : if (isl_tab_freeze_constraint(tab, tab->n_con - 1) < 0)
919 0 : goto error;
920 :
921 0 : for (i = 0; i < vertices->n_vertices; ++i) {
922 0 : selection[i] = bset_covers_tab(vertices->v[i].dom,
923 0 : todo->tab);
924 0 : if (selection[i] < 0)
925 0 : goto error;
926 0 : if (!selection[i])
927 0 : continue;
928 0 : selection[i] = can_intersect(tab, vertices->v[i].dom);
929 0 : if (selection[i] < 0)
930 0 : goto error;
931 : }
932 :
933 0 : if (isl_tab_detect_redundant(tab) < 0)
934 0 : goto error;
935 :
936 0 : if (add_chamber(&list, vertices, tab, selection) < 0)
937 0 : goto error;
938 0 : n_chambers++;
939 :
940 0 : if (update_todo(todo, tab) < 0)
941 0 : goto error;
942 :
943 0 : next = todo->next;
944 0 : todo->next = NULL;
945 0 : free_todo(todo);
946 0 : todo = next;
947 : }
948 :
949 0 : isl_vec_free(sample);
950 :
951 0 : isl_tab_free(tab);
952 0 : free(selection);
953 :
954 0 : vertices = vertices_add_chambers(vertices, n_chambers, list);
955 :
956 0 : for (i = 0; vertices && i < vertices->n_vertices; ++i) {
957 0 : isl_basic_set_free(vertices->v[i].dom);
958 0 : vertices->v[i].dom = NULL;
959 : }
960 :
961 0 : return vertices;
962 : error:
963 0 : free_chamber_list(list);
964 0 : free_todo(todo);
965 0 : isl_vec_free(sample);
966 0 : isl_tab_free(tab);
967 0 : free(selection);
968 0 : if (!tab)
969 0 : isl_basic_set_free(bset);
970 0 : isl_vertices_free(vertices);
971 0 : return NULL;
972 : }
973 :
974 0 : isl_ctx *isl_vertex_get_ctx(__isl_keep isl_vertex *vertex)
975 : {
976 0 : return vertex ? isl_vertices_get_ctx(vertex->vertices) : NULL;
977 : }
978 :
979 0 : int isl_vertex_get_id(__isl_keep isl_vertex *vertex)
980 : {
981 0 : return vertex ? vertex->id : -1;
982 : }
983 :
984 : /* Return the activity domain of the vertex "vertex".
985 : */
986 0 : __isl_give isl_basic_set *isl_vertex_get_domain(__isl_keep isl_vertex *vertex)
987 : {
988 : struct isl_vertex *v;
989 :
990 0 : if (!vertex)
991 0 : return NULL;
992 :
993 0 : v = &vertex->vertices->v[vertex->id];
994 0 : if (!v->dom) {
995 0 : v->dom = isl_basic_set_copy(v->vertex);
996 0 : v->dom = isl_basic_set_params(v->dom);
997 0 : v->dom = isl_basic_set_set_integral(v->dom);
998 : }
999 :
1000 0 : return isl_basic_set_copy(v->dom);
1001 : }
1002 :
1003 : /* Return a multiple quasi-affine expression describing the vertex "vertex"
1004 : * in terms of the parameters,
1005 : */
1006 0 : __isl_give isl_multi_aff *isl_vertex_get_expr(__isl_keep isl_vertex *vertex)
1007 : {
1008 : struct isl_vertex *v;
1009 : isl_basic_set *bset;
1010 :
1011 0 : if (!vertex)
1012 0 : return NULL;
1013 :
1014 0 : v = &vertex->vertices->v[vertex->id];
1015 :
1016 0 : bset = isl_basic_set_copy(v->vertex);
1017 0 : return isl_multi_aff_from_basic_set_equalities(bset);
1018 : }
1019 :
1020 0 : static __isl_give isl_vertex *isl_vertex_alloc(__isl_take isl_vertices *vertices,
1021 : int id)
1022 : {
1023 : isl_ctx *ctx;
1024 : isl_vertex *vertex;
1025 :
1026 0 : if (!vertices)
1027 0 : return NULL;
1028 :
1029 0 : ctx = isl_vertices_get_ctx(vertices);
1030 0 : vertex = isl_alloc_type(ctx, isl_vertex);
1031 0 : if (!vertex)
1032 0 : goto error;
1033 :
1034 0 : vertex->vertices = vertices;
1035 0 : vertex->id = id;
1036 :
1037 0 : return vertex;
1038 : error:
1039 0 : isl_vertices_free(vertices);
1040 0 : return NULL;
1041 : }
1042 :
1043 0 : void isl_vertex_free(__isl_take isl_vertex *vertex)
1044 : {
1045 0 : if (!vertex)
1046 0 : return;
1047 0 : isl_vertices_free(vertex->vertices);
1048 0 : free(vertex);
1049 : }
1050 :
1051 0 : isl_ctx *isl_cell_get_ctx(__isl_keep isl_cell *cell)
1052 : {
1053 0 : return cell ? cell->dom->ctx : NULL;
1054 : }
1055 :
1056 0 : __isl_give isl_basic_set *isl_cell_get_domain(__isl_keep isl_cell *cell)
1057 : {
1058 0 : return cell ? isl_basic_set_copy(cell->dom) : NULL;
1059 : }
1060 :
1061 0 : static __isl_give isl_cell *isl_cell_alloc(__isl_take isl_vertices *vertices,
1062 : __isl_take isl_basic_set *dom, int id)
1063 : {
1064 : int i;
1065 0 : isl_cell *cell = NULL;
1066 :
1067 0 : if (!vertices || !dom)
1068 : goto error;
1069 :
1070 0 : cell = isl_calloc_type(dom->ctx, isl_cell);
1071 0 : if (!cell)
1072 0 : goto error;
1073 :
1074 0 : cell->n_vertices = vertices->c[id].n_vertices;
1075 0 : cell->ids = isl_alloc_array(dom->ctx, int, cell->n_vertices);
1076 0 : if (cell->n_vertices && !cell->ids)
1077 0 : goto error;
1078 0 : for (i = 0; i < cell->n_vertices; ++i)
1079 0 : cell->ids[i] = vertices->c[id].vertices[i];
1080 0 : cell->vertices = vertices;
1081 0 : cell->dom = dom;
1082 :
1083 0 : return cell;
1084 : error:
1085 0 : isl_cell_free(cell);
1086 0 : isl_vertices_free(vertices);
1087 0 : isl_basic_set_free(dom);
1088 0 : return NULL;
1089 : }
1090 :
1091 0 : void isl_cell_free(__isl_take isl_cell *cell)
1092 : {
1093 0 : if (!cell)
1094 0 : return;
1095 :
1096 0 : isl_vertices_free(cell->vertices);
1097 0 : free(cell->ids);
1098 0 : isl_basic_set_free(cell->dom);
1099 0 : free(cell);
1100 : }
1101 :
1102 : /* Create a tableau of the cone obtained by first homogenizing the given
1103 : * polytope and then making all inequalities strict by setting the
1104 : * constant term to -1.
1105 : */
1106 0 : static struct isl_tab *tab_for_shifted_cone(__isl_keep isl_basic_set *bset)
1107 : {
1108 : int i;
1109 0 : isl_vec *c = NULL;
1110 : struct isl_tab *tab;
1111 :
1112 0 : if (!bset)
1113 0 : return NULL;
1114 0 : tab = isl_tab_alloc(bset->ctx, bset->n_eq + bset->n_ineq + 1,
1115 0 : 1 + isl_basic_set_total_dim(bset), 0);
1116 0 : if (!tab)
1117 0 : return NULL;
1118 0 : tab->rational = ISL_F_ISSET(bset, ISL_BASIC_SET_RATIONAL);
1119 0 : if (ISL_F_ISSET(bset, ISL_BASIC_MAP_EMPTY)) {
1120 0 : if (isl_tab_mark_empty(tab) < 0)
1121 0 : goto error;
1122 0 : return tab;
1123 : }
1124 :
1125 0 : c = isl_vec_alloc(bset->ctx, 1 + 1 + isl_basic_set_total_dim(bset));
1126 0 : if (!c)
1127 0 : goto error;
1128 :
1129 0 : isl_int_set_si(c->el[0], 0);
1130 0 : for (i = 0; i < bset->n_eq; ++i) {
1131 0 : isl_seq_cpy(c->el + 1, bset->eq[i], c->size - 1);
1132 0 : if (isl_tab_add_eq(tab, c->el) < 0)
1133 0 : goto error;
1134 : }
1135 :
1136 0 : isl_int_set_si(c->el[0], -1);
1137 0 : for (i = 0; i < bset->n_ineq; ++i) {
1138 0 : isl_seq_cpy(c->el + 1, bset->ineq[i], c->size - 1);
1139 0 : if (isl_tab_add_ineq(tab, c->el) < 0)
1140 0 : goto error;
1141 0 : if (tab->empty) {
1142 0 : isl_vec_free(c);
1143 0 : return tab;
1144 : }
1145 : }
1146 :
1147 0 : isl_seq_clr(c->el + 1, c->size - 1);
1148 0 : isl_int_set_si(c->el[1], 1);
1149 0 : if (isl_tab_add_ineq(tab, c->el) < 0)
1150 0 : goto error;
1151 :
1152 0 : isl_vec_free(c);
1153 0 : return tab;
1154 : error:
1155 0 : isl_vec_free(c);
1156 0 : isl_tab_free(tab);
1157 0 : return NULL;
1158 : }
1159 :
1160 : /* Compute an interior point of "bset" by selecting an interior
1161 : * point in homogeneous space and projecting the point back down.
1162 : */
1163 0 : static __isl_give isl_vec *isl_basic_set_interior_point(
1164 : __isl_keep isl_basic_set *bset)
1165 : {
1166 : isl_vec *vec;
1167 : struct isl_tab *tab;
1168 :
1169 0 : tab = tab_for_shifted_cone(bset);
1170 0 : vec = isl_tab_get_sample_value(tab);
1171 0 : isl_tab_free(tab);
1172 0 : if (!vec)
1173 0 : return NULL;
1174 :
1175 0 : isl_seq_cpy(vec->el, vec->el + 1, vec->size - 1);
1176 0 : vec->size--;
1177 :
1178 0 : return vec;
1179 : }
1180 :
1181 : /* Call "fn" on all chambers of the parametric polytope with the shared
1182 : * facets of neighboring chambers only appearing in one of the chambers.
1183 : *
1184 : * We pick an interior point from one of the chambers and then make
1185 : * all constraints that do not satisfy this point strict.
1186 : * For constraints that saturate the interior point, the sign
1187 : * of the first non-zero coefficient is used to determine which
1188 : * of the two (internal) constraints should be tightened.
1189 : */
1190 0 : isl_stat isl_vertices_foreach_disjoint_cell(__isl_keep isl_vertices *vertices,
1191 : isl_stat (*fn)(__isl_take isl_cell *cell, void *user), void *user)
1192 : {
1193 : int i;
1194 : isl_vec *vec;
1195 : isl_cell *cell;
1196 :
1197 0 : if (!vertices)
1198 0 : return isl_stat_error;
1199 :
1200 0 : if (vertices->n_chambers == 0)
1201 0 : return isl_stat_ok;
1202 :
1203 0 : if (vertices->n_chambers == 1) {
1204 0 : isl_basic_set *dom = isl_basic_set_copy(vertices->c[0].dom);
1205 0 : dom = isl_basic_set_set_integral(dom);
1206 0 : cell = isl_cell_alloc(isl_vertices_copy(vertices), dom, 0);
1207 0 : if (!cell)
1208 0 : return isl_stat_error;
1209 0 : return fn(cell, user);
1210 : }
1211 :
1212 0 : vec = isl_basic_set_interior_point(vertices->c[0].dom);
1213 0 : if (!vec)
1214 0 : return isl_stat_error;
1215 :
1216 0 : for (i = 0; i < vertices->n_chambers; ++i) {
1217 : int r;
1218 0 : isl_basic_set *dom = isl_basic_set_copy(vertices->c[i].dom);
1219 0 : if (i)
1220 0 : dom = isl_basic_set_tighten_outward(dom, vec);
1221 0 : dom = isl_basic_set_set_integral(dom);
1222 0 : cell = isl_cell_alloc(isl_vertices_copy(vertices), dom, i);
1223 0 : if (!cell)
1224 0 : goto error;
1225 0 : r = fn(cell, user);
1226 0 : if (r < 0)
1227 0 : goto error;
1228 : }
1229 :
1230 0 : isl_vec_free(vec);
1231 :
1232 0 : return isl_stat_ok;
1233 : error:
1234 0 : isl_vec_free(vec);
1235 0 : return isl_stat_error;
1236 : }
1237 :
1238 0 : isl_stat isl_vertices_foreach_cell(__isl_keep isl_vertices *vertices,
1239 : isl_stat (*fn)(__isl_take isl_cell *cell, void *user), void *user)
1240 : {
1241 : int i;
1242 : isl_cell *cell;
1243 :
1244 0 : if (!vertices)
1245 0 : return isl_stat_error;
1246 :
1247 0 : if (vertices->n_chambers == 0)
1248 0 : return isl_stat_ok;
1249 :
1250 0 : for (i = 0; i < vertices->n_chambers; ++i) {
1251 : isl_stat r;
1252 0 : isl_basic_set *dom = isl_basic_set_copy(vertices->c[i].dom);
1253 :
1254 0 : cell = isl_cell_alloc(isl_vertices_copy(vertices), dom, i);
1255 0 : if (!cell)
1256 0 : return isl_stat_error;
1257 :
1258 0 : r = fn(cell, user);
1259 0 : if (r < 0)
1260 0 : return isl_stat_error;
1261 : }
1262 :
1263 0 : return isl_stat_ok;
1264 : }
1265 :
1266 0 : isl_stat isl_vertices_foreach_vertex(__isl_keep isl_vertices *vertices,
1267 : isl_stat (*fn)(__isl_take isl_vertex *vertex, void *user), void *user)
1268 : {
1269 : int i;
1270 : isl_vertex *vertex;
1271 :
1272 0 : if (!vertices)
1273 0 : return isl_stat_error;
1274 :
1275 0 : if (vertices->n_vertices == 0)
1276 0 : return isl_stat_ok;
1277 :
1278 0 : for (i = 0; i < vertices->n_vertices; ++i) {
1279 : isl_stat r;
1280 :
1281 0 : vertex = isl_vertex_alloc(isl_vertices_copy(vertices), i);
1282 0 : if (!vertex)
1283 0 : return isl_stat_error;
1284 :
1285 0 : r = fn(vertex, user);
1286 0 : if (r < 0)
1287 0 : return isl_stat_error;
1288 : }
1289 :
1290 0 : return isl_stat_ok;
1291 : }
1292 :
1293 0 : isl_stat isl_cell_foreach_vertex(__isl_keep isl_cell *cell,
1294 : isl_stat (*fn)(__isl_take isl_vertex *vertex, void *user), void *user)
1295 : {
1296 : int i;
1297 : isl_vertex *vertex;
1298 :
1299 0 : if (!cell)
1300 0 : return isl_stat_error;
1301 :
1302 0 : if (cell->n_vertices == 0)
1303 0 : return isl_stat_ok;
1304 :
1305 0 : for (i = 0; i < cell->n_vertices; ++i) {
1306 : isl_stat r;
1307 :
1308 0 : vertex = isl_vertex_alloc(isl_vertices_copy(cell->vertices),
1309 0 : cell->ids[i]);
1310 0 : if (!vertex)
1311 0 : return isl_stat_error;
1312 :
1313 0 : r = fn(vertex, user);
1314 0 : if (r < 0)
1315 0 : return isl_stat_error;
1316 : }
1317 :
1318 0 : return isl_stat_ok;
1319 : }
1320 :
1321 0 : isl_ctx *isl_vertices_get_ctx(__isl_keep isl_vertices *vertices)
1322 : {
1323 0 : return vertices ? vertices->bset->ctx : NULL;
1324 : }
1325 :
1326 0 : int isl_vertices_get_n_vertices(__isl_keep isl_vertices *vertices)
1327 : {
1328 0 : return vertices ? vertices->n_vertices : -1;
1329 : }
1330 :
1331 0 : __isl_give isl_vertices *isl_morph_vertices(__isl_take isl_morph *morph,
1332 : __isl_take isl_vertices *vertices)
1333 : {
1334 : int i;
1335 0 : isl_morph *param_morph = NULL;
1336 :
1337 0 : if (!morph || !vertices)
1338 : goto error;
1339 :
1340 0 : isl_assert(vertices->bset->ctx, vertices->ref == 1, goto error);
1341 :
1342 0 : param_morph = isl_morph_copy(morph);
1343 0 : param_morph = isl_morph_dom_params(param_morph);
1344 0 : param_morph = isl_morph_ran_params(param_morph);
1345 :
1346 0 : for (i = 0; i < vertices->n_vertices; ++i) {
1347 0 : vertices->v[i].dom = isl_morph_basic_set(
1348 0 : isl_morph_copy(param_morph), vertices->v[i].dom);
1349 0 : vertices->v[i].vertex = isl_morph_basic_set(
1350 0 : isl_morph_copy(morph), vertices->v[i].vertex);
1351 0 : if (!vertices->v[i].vertex)
1352 0 : goto error;
1353 : }
1354 :
1355 0 : for (i = 0; i < vertices->n_chambers; ++i) {
1356 0 : vertices->c[i].dom = isl_morph_basic_set(
1357 0 : isl_morph_copy(param_morph), vertices->c[i].dom);
1358 0 : if (!vertices->c[i].dom)
1359 0 : goto error;
1360 : }
1361 :
1362 0 : isl_morph_free(param_morph);
1363 0 : isl_morph_free(morph);
1364 0 : return vertices;
1365 : error:
1366 0 : isl_morph_free(param_morph);
1367 0 : isl_morph_free(morph);
1368 0 : isl_vertices_free(vertices);
1369 0 : return NULL;
1370 : }
1371 :
1372 : /* Construct a simplex isl_cell spanned by the vertices with indices in
1373 : * "simplex_ids" and "other_ids" and call "fn" on this isl_cell.
1374 : */
1375 0 : static isl_stat call_on_simplex(__isl_keep isl_cell *cell,
1376 : int *simplex_ids, int n_simplex, int *other_ids, int n_other,
1377 : isl_stat (*fn)(__isl_take isl_cell *simplex, void *user), void *user)
1378 : {
1379 : int i;
1380 : isl_ctx *ctx;
1381 : struct isl_cell *simplex;
1382 :
1383 0 : ctx = isl_cell_get_ctx(cell);
1384 :
1385 0 : simplex = isl_calloc_type(ctx, struct isl_cell);
1386 0 : if (!simplex)
1387 0 : return isl_stat_error;
1388 0 : simplex->vertices = isl_vertices_copy(cell->vertices);
1389 0 : if (!simplex->vertices)
1390 0 : goto error;
1391 0 : simplex->dom = isl_basic_set_copy(cell->dom);
1392 0 : if (!simplex->dom)
1393 0 : goto error;
1394 0 : simplex->n_vertices = n_simplex + n_other;
1395 0 : simplex->ids = isl_alloc_array(ctx, int, simplex->n_vertices);
1396 0 : if (!simplex->ids)
1397 0 : goto error;
1398 :
1399 0 : for (i = 0; i < n_simplex; ++i)
1400 0 : simplex->ids[i] = simplex_ids[i];
1401 0 : for (i = 0; i < n_other; ++i)
1402 0 : simplex->ids[n_simplex + i] = other_ids[i];
1403 :
1404 0 : return fn(simplex, user);
1405 : error:
1406 0 : isl_cell_free(simplex);
1407 0 : return isl_stat_error;
1408 : }
1409 :
1410 : /* Check whether the parametric vertex described by "vertex"
1411 : * lies on the facet corresponding to constraint "facet" of "bset".
1412 : * The isl_vec "v" is a temporary vector than can be used by this function.
1413 : *
1414 : * We eliminate the variables from the facet constraint using the
1415 : * equalities defining the vertex and check if the result is identical
1416 : * to zero.
1417 : *
1418 : * It would probably be better to keep track of the constraints defining
1419 : * a vertex during the vertex construction so that we could simply look
1420 : * it up here.
1421 : */
1422 0 : static int vertex_on_facet(__isl_keep isl_basic_set *vertex,
1423 : __isl_keep isl_basic_set *bset, int facet, __isl_keep isl_vec *v)
1424 : {
1425 : int i;
1426 : isl_int m;
1427 :
1428 0 : isl_seq_cpy(v->el, bset->ineq[facet], v->size);
1429 :
1430 0 : isl_int_init(m);
1431 0 : for (i = 0; i < vertex->n_eq; ++i) {
1432 0 : int k = isl_seq_last_non_zero(vertex->eq[i], v->size);
1433 0 : isl_seq_elim(v->el, vertex->eq[i], k, v->size, &m);
1434 : }
1435 0 : isl_int_clear(m);
1436 :
1437 0 : return isl_seq_first_non_zero(v->el, v->size) == -1;
1438 : }
1439 :
1440 : /* Triangulate the polytope spanned by the vertices with ids
1441 : * in "simplex_ids" and "other_ids" and call "fn" on each of
1442 : * the resulting simplices.
1443 : * If the input polytope is already a simplex, we simply call "fn".
1444 : * Otherwise, we pick a point from "other_ids" and add it to "simplex_ids".
1445 : * Then we consider each facet of "bset" that does not contain the point
1446 : * we just picked, but does contain some of the other points in "other_ids"
1447 : * and call ourselves recursively on the polytope spanned by the new
1448 : * "simplex_ids" and those points in "other_ids" that lie on the facet.
1449 : */
1450 0 : static isl_stat triangulate(__isl_keep isl_cell *cell, __isl_keep isl_vec *v,
1451 : int *simplex_ids, int n_simplex, int *other_ids, int n_other,
1452 : isl_stat (*fn)(__isl_take isl_cell *simplex, void *user), void *user)
1453 : {
1454 : int i, j, k;
1455 : int d, nparam;
1456 : int *ids;
1457 : isl_ctx *ctx;
1458 : isl_basic_set *vertex;
1459 : isl_basic_set *bset;
1460 :
1461 0 : ctx = isl_cell_get_ctx(cell);
1462 0 : d = isl_basic_set_dim(cell->vertices->bset, isl_dim_set);
1463 0 : nparam = isl_basic_set_dim(cell->vertices->bset, isl_dim_param);
1464 :
1465 0 : if (n_simplex + n_other == d + 1)
1466 0 : return call_on_simplex(cell, simplex_ids, n_simplex,
1467 : other_ids, n_other, fn, user);
1468 :
1469 0 : simplex_ids[n_simplex] = other_ids[0];
1470 0 : vertex = cell->vertices->v[other_ids[0]].vertex;
1471 0 : bset = cell->vertices->bset;
1472 :
1473 0 : ids = isl_alloc_array(ctx, int, n_other - 1);
1474 0 : if (!ids)
1475 0 : goto error;
1476 0 : for (i = 0; i < bset->n_ineq; ++i) {
1477 0 : if (isl_seq_first_non_zero(bset->ineq[i] + 1 + nparam, d) == -1)
1478 0 : continue;
1479 0 : if (vertex_on_facet(vertex, bset, i, v))
1480 0 : continue;
1481 :
1482 0 : for (j = 1, k = 0; j < n_other; ++j) {
1483 : isl_basic_set *ov;
1484 0 : ov = cell->vertices->v[other_ids[j]].vertex;
1485 0 : if (vertex_on_facet(ov, bset, i, v))
1486 0 : ids[k++] = other_ids[j];
1487 : }
1488 0 : if (k == 0)
1489 0 : continue;
1490 :
1491 0 : if (triangulate(cell, v, simplex_ids, n_simplex + 1,
1492 : ids, k, fn, user) < 0)
1493 0 : goto error;
1494 : }
1495 0 : free(ids);
1496 :
1497 0 : return isl_stat_ok;
1498 : error:
1499 0 : free(ids);
1500 0 : return isl_stat_error;
1501 : }
1502 :
1503 : /* Triangulate the given cell and call "fn" on each of the resulting
1504 : * simplices.
1505 : */
1506 0 : isl_stat isl_cell_foreach_simplex(__isl_take isl_cell *cell,
1507 : isl_stat (*fn)(__isl_take isl_cell *simplex, void *user), void *user)
1508 : {
1509 : int d, total;
1510 : isl_stat r;
1511 : isl_ctx *ctx;
1512 0 : isl_vec *v = NULL;
1513 0 : int *simplex_ids = NULL;
1514 :
1515 0 : if (!cell)
1516 0 : return isl_stat_error;
1517 :
1518 0 : d = isl_basic_set_dim(cell->vertices->bset, isl_dim_set);
1519 0 : total = isl_basic_set_total_dim(cell->vertices->bset);
1520 :
1521 0 : if (cell->n_vertices == d + 1)
1522 0 : return fn(cell, user);
1523 :
1524 0 : ctx = isl_cell_get_ctx(cell);
1525 0 : simplex_ids = isl_alloc_array(ctx, int, d + 1);
1526 0 : if (!simplex_ids)
1527 0 : goto error;
1528 :
1529 0 : v = isl_vec_alloc(ctx, 1 + total);
1530 0 : if (!v)
1531 0 : goto error;
1532 :
1533 0 : r = triangulate(cell, v, simplex_ids, 0,
1534 : cell->ids, cell->n_vertices, fn, user);
1535 :
1536 0 : isl_vec_free(v);
1537 0 : free(simplex_ids);
1538 :
1539 0 : isl_cell_free(cell);
1540 :
1541 0 : return r;
1542 : error:
1543 0 : free(simplex_ids);
1544 0 : isl_vec_free(v);
1545 0 : isl_cell_free(cell);
1546 0 : return isl_stat_error;
1547 : }
|
