1 /* 2 * Copyright 2005-2007 Universiteit Leiden 3 * Copyright 2008-2009 Katholieke Universiteit Leuven 4 * Copyright 2010 INRIA Saclay 5 * 6 * Use of this software is governed by the MIT license 7 * 8 * Written by Sven Verdoolaege, Leiden Institute of Advanced Computer Science, 9 * Universiteit Leiden, Niels Bohrweg 1, 2333 CA Leiden, The Netherlands 10 * and K.U.Leuven, Departement Computerwetenschappen, Celestijnenlaan 200A, 11 * B-3001 Leuven, Belgium 12 * and INRIA Saclay - Ile-de-France, Parc Club Orsay Universite, 13 * ZAC des vignes, 4 rue Jacques Monod, 91893 Orsay, France 14 */ 15 16 #include <isl_map_private.h> 17 #include <isl_factorization.h> 18 #include <isl_space_private.h> 19 #include <isl_mat_private.h> 20 21 /* Return the isl_ctx to which "f" belongs. 22 */ 23 isl_ctx *isl_factorizer_get_ctx(__isl_keep isl_factorizer *f) 24 { 25 if (!f) 26 return NULL; 27 return isl_basic_set_get_ctx(f->bset); 28 } 29 30 static __isl_give isl_factorizer *isl_factorizer_alloc( 31 __isl_keep isl_basic_set *bset, __isl_take isl_morph *morph, 32 int n_group) 33 { 34 isl_factorizer *f = NULL; 35 int *len = NULL; 36 37 if (!morph) 38 return NULL; 39 40 if (n_group > 0) { 41 len = isl_alloc_array(morph->dom->ctx, int, n_group); 42 if (!len) 43 goto error; 44 } 45 46 f = isl_alloc_type(morph->dom->ctx, struct isl_factorizer); 47 if (!f) 48 goto error; 49 50 f->bset = isl_basic_set_copy(bset); 51 f->morph = morph; 52 f->n_group = n_group; 53 f->len = len; 54 55 return f; 56 error: 57 free(len); 58 isl_morph_free(morph); 59 return NULL; 60 } 61 62 __isl_null isl_factorizer *isl_factorizer_free(__isl_take isl_factorizer *f) 63 { 64 if (!f) 65 return NULL; 66 67 isl_basic_set_free(f->bset); 68 isl_morph_free(f->morph); 69 free(f->len); 70 free(f); 71 return NULL; 72 } 73 74 void isl_factorizer_dump(__isl_take isl_factorizer *f) 75 { 76 int i; 77 78 if (!f) 79 return; 80 81 isl_morph_print_internal(f->morph, stderr); 82 fprintf(stderr, "["); 83 for (i = 0; i < f->n_group; ++i) { 84 if (i) 85 fprintf(stderr, ", "); 86 fprintf(stderr, "%d", f->len[i]); 87 } 88 fprintf(stderr, "]\n"); 89 } 90 91 __isl_give isl_factorizer *isl_factorizer_identity(__isl_keep isl_basic_set *bset) 92 { 93 return isl_factorizer_alloc(bset, isl_morph_identity(bset), 0); 94 } 95 96 __isl_give isl_factorizer *isl_factorizer_groups(__isl_keep isl_basic_set *bset, 97 __isl_take isl_mat *Q, __isl_take isl_mat *U, int n, int *len) 98 { 99 int i; 100 isl_size nvar, off; 101 isl_space *space; 102 isl_basic_set *dom; 103 isl_basic_set *ran; 104 isl_morph *morph; 105 isl_factorizer *f; 106 isl_mat *id; 107 108 nvar = isl_basic_set_dim(bset, isl_dim_set); 109 off = isl_basic_set_var_offset(bset, isl_dim_set); 110 if (nvar < 0 || off < 0 || !Q || !U) 111 goto error; 112 113 id = isl_mat_identity(bset->ctx, 1 + off); 114 Q = isl_mat_diagonal(isl_mat_copy(id), Q); 115 U = isl_mat_diagonal(id, U); 116 117 space = isl_basic_set_get_space(bset); 118 dom = isl_basic_set_universe(isl_space_copy(space)); 119 space = isl_space_drop_dims(space, isl_dim_set, 0, nvar); 120 space = isl_space_add_dims(space, isl_dim_set, nvar); 121 ran = isl_basic_set_universe(space); 122 morph = isl_morph_alloc(dom, ran, Q, U); 123 f = isl_factorizer_alloc(bset, morph, n); 124 if (!f) 125 return NULL; 126 for (i = 0; i < n; ++i) 127 f->len[i] = len[i]; 128 return f; 129 error: 130 isl_mat_free(Q); 131 isl_mat_free(U); 132 return NULL; 133 } 134 135 struct isl_factor_groups { 136 int *pos; /* for each column: row position of pivot */ 137 int *group; /* group to which a column belongs */ 138 int *cnt; /* number of columns in the group */ 139 int *rowgroup; /* group to which a constraint belongs */ 140 }; 141 142 /* Initialize isl_factor_groups structure: find pivot row positions, 143 * each column initially belongs to its own group and the groups 144 * of the constraints are still unknown. 145 */ 146 static int init_groups(struct isl_factor_groups *g, __isl_keep isl_mat *H) 147 { 148 int i, j; 149 150 if (!H) 151 return -1; 152 153 g->pos = isl_alloc_array(H->ctx, int, H->n_col); 154 g->group = isl_alloc_array(H->ctx, int, H->n_col); 155 g->cnt = isl_alloc_array(H->ctx, int, H->n_col); 156 g->rowgroup = isl_alloc_array(H->ctx, int, H->n_row); 157 158 if (!g->pos || !g->group || !g->cnt || !g->rowgroup) 159 return -1; 160 161 for (i = 0; i < H->n_row; ++i) 162 g->rowgroup[i] = -1; 163 for (i = 0, j = 0; i < H->n_col; ++i) { 164 for ( ; j < H->n_row; ++j) 165 if (!isl_int_is_zero(H->row[j][i])) 166 break; 167 g->pos[i] = j; 168 } 169 for (i = 0; i < H->n_col; ++i) { 170 g->group[i] = i; 171 g->cnt[i] = 1; 172 } 173 174 return 0; 175 } 176 177 /* Update group[k] to the group column k belongs to. 178 * When merging two groups, only the group of the current 179 * group leader is changed. Here we change the group of 180 * the other members to also point to the group that the 181 * old group leader now points to. 182 */ 183 static void update_group(struct isl_factor_groups *g, int k) 184 { 185 int p = g->group[k]; 186 while (g->cnt[p] == 0) 187 p = g->group[p]; 188 g->group[k] = p; 189 } 190 191 /* Merge group i with all groups of the subsequent columns 192 * with non-zero coefficients in row j of H. 193 * (The previous columns are all zero; otherwise we would have handled 194 * the row before.) 195 */ 196 static int update_group_i_with_row_j(struct isl_factor_groups *g, int i, int j, 197 __isl_keep isl_mat *H) 198 { 199 int k; 200 201 g->rowgroup[j] = g->group[i]; 202 for (k = i + 1; k < H->n_col && j >= g->pos[k]; ++k) { 203 update_group(g, k); 204 update_group(g, i); 205 if (g->group[k] != g->group[i] && 206 !isl_int_is_zero(H->row[j][k])) { 207 isl_assert(H->ctx, g->cnt[g->group[k]] != 0, return -1); 208 isl_assert(H->ctx, g->cnt[g->group[i]] != 0, return -1); 209 if (g->group[i] < g->group[k]) { 210 g->cnt[g->group[i]] += g->cnt[g->group[k]]; 211 g->cnt[g->group[k]] = 0; 212 g->group[g->group[k]] = g->group[i]; 213 } else { 214 g->cnt[g->group[k]] += g->cnt[g->group[i]]; 215 g->cnt[g->group[i]] = 0; 216 g->group[g->group[i]] = g->group[k]; 217 } 218 } 219 } 220 221 return 0; 222 } 223 224 /* Update the group information based on the constraint matrix. 225 */ 226 static int update_groups(struct isl_factor_groups *g, __isl_keep isl_mat *H) 227 { 228 int i, j; 229 230 for (i = 0; i < H->n_col && g->cnt[0] < H->n_col; ++i) { 231 if (g->pos[i] == H->n_row) 232 continue; /* A line direction */ 233 if (g->rowgroup[g->pos[i]] == -1) 234 g->rowgroup[g->pos[i]] = i; 235 for (j = g->pos[i] + 1; j < H->n_row; ++j) { 236 if (isl_int_is_zero(H->row[j][i])) 237 continue; 238 if (g->rowgroup[j] != -1) 239 continue; 240 if (update_group_i_with_row_j(g, i, j, H) < 0) 241 return -1; 242 } 243 } 244 for (i = 1; i < H->n_col; ++i) 245 update_group(g, i); 246 247 return 0; 248 } 249 250 static void clear_groups(struct isl_factor_groups *g) 251 { 252 if (!g) 253 return; 254 free(g->pos); 255 free(g->group); 256 free(g->cnt); 257 free(g->rowgroup); 258 } 259 260 /* Determine if the set variables of the basic set can be factorized and 261 * return the results in an isl_factorizer. 262 * 263 * The algorithm works by first computing the Hermite normal form 264 * and then grouping columns linked by one or more constraints together, 265 * where a constraints "links" two or more columns if the constraint 266 * has nonzero coefficients in the columns. 267 */ 268 __isl_give isl_factorizer *isl_basic_set_factorizer( 269 __isl_keep isl_basic_set *bset) 270 { 271 int i, j, n, done; 272 isl_mat *H, *U, *Q; 273 isl_size nvar, first; 274 struct isl_factor_groups g = { 0 }; 275 isl_factorizer *f; 276 277 nvar = isl_basic_set_dim(bset, isl_dim_set); 278 first = isl_basic_set_var_offset(bset, isl_dim_set); 279 if (nvar < 0 || first < 0 || isl_basic_set_check_no_locals(bset) < 0) 280 return NULL; 281 282 if (nvar <= 1) 283 return isl_factorizer_identity(bset); 284 285 H = isl_mat_alloc(bset->ctx, bset->n_eq + bset->n_ineq, nvar); 286 if (!H) 287 return NULL; 288 isl_mat_sub_copy(bset->ctx, H->row, bset->eq, bset->n_eq, 289 0, 1 + first, nvar); 290 isl_mat_sub_copy(bset->ctx, H->row + bset->n_eq, bset->ineq, bset->n_ineq, 291 0, 1 + first, nvar); 292 H = isl_mat_left_hermite(H, 0, &U, &Q); 293 294 if (init_groups(&g, H) < 0) 295 goto error; 296 if (update_groups(&g, H) < 0) 297 goto error; 298 299 if (g.cnt[0] == nvar) { 300 isl_mat_free(H); 301 isl_mat_free(U); 302 isl_mat_free(Q); 303 clear_groups(&g); 304 305 return isl_factorizer_identity(bset); 306 } 307 308 done = 0; 309 n = 0; 310 while (done != nvar) { 311 int group = g.group[done]; 312 for (i = 1; i < g.cnt[group]; ++i) { 313 if (g.group[done + i] == group) 314 continue; 315 for (j = done + g.cnt[group]; j < nvar; ++j) 316 if (g.group[j] == group) 317 break; 318 if (j == nvar) 319 isl_die(bset->ctx, isl_error_internal, 320 "internal error", goto error); 321 g.group[j] = g.group[done + i]; 322 Q = isl_mat_swap_rows(Q, done + i, j); 323 U = isl_mat_swap_cols(U, done + i, j); 324 } 325 done += g.cnt[group]; 326 g.pos[n++] = g.cnt[group]; 327 } 328 329 f = isl_factorizer_groups(bset, Q, U, n, g.pos); 330 331 isl_mat_free(H); 332 clear_groups(&g); 333 334 return f; 335 error: 336 isl_mat_free(H); 337 isl_mat_free(U); 338 isl_mat_free(Q); 339 clear_groups(&g); 340 return NULL; 341 } 342 343 /* Given the factorizer "f" of a basic set, 344 * call "test" on each resulting factor as long as each call succeeds. 345 */ 346 __isl_give isl_bool isl_factorizer_every_factor_basic_set( 347 __isl_keep isl_factorizer *f, 348 isl_bool (*test)(__isl_keep isl_basic_set *bset, void *user), 349 void *user) 350 { 351 int i, n; 352 isl_bool every = isl_bool_true; 353 isl_size nparam, nvar; 354 isl_basic_set *bset; 355 356 if (!f) 357 return isl_bool_error; 358 nparam = isl_basic_set_dim(f->bset, isl_dim_param); 359 nvar = isl_basic_set_dim(f->bset, isl_dim_set); 360 if (nparam < 0 || nvar < 0) 361 return isl_bool_error; 362 363 bset = isl_basic_set_copy(f->bset); 364 bset = isl_morph_basic_set(isl_morph_copy(f->morph), bset); 365 366 for (i = 0, n = 0; i < f->n_group; ++i) { 367 isl_basic_set *factor; 368 369 factor = isl_basic_set_copy(bset); 370 factor = isl_basic_set_drop_constraints_involving(factor, 371 nparam + n + f->len[i], nvar - n - f->len[i]); 372 factor = isl_basic_set_drop_constraints_involving(factor, 373 nparam, n); 374 factor = isl_basic_set_drop(factor, isl_dim_set, 375 n + f->len[i], nvar - n - f->len[i]); 376 factor = isl_basic_set_drop(factor, isl_dim_set, 0, n); 377 every = test(factor, user); 378 isl_basic_set_free(factor); 379 380 if (every < 0 || !every) 381 break; 382 383 n += f->len[i]; 384 } 385 386 isl_basic_set_free(bset); 387 388 return every; 389 } 390