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