/Users/buildslave/jenkins/workspace/clang-stage2-coverage-R/llvm/tools/polly/lib/External/isl/isl_transitive_closure.c
Line | Count | Source (jump to first uncovered line) |
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_ctx_private.h> |
12 | | #include <isl_map_private.h> |
13 | | #include <isl/map.h> |
14 | | #include <isl_seq.h> |
15 | | #include <isl_space_private.h> |
16 | | #include <isl_lp_private.h> |
17 | | #include <isl/union_map.h> |
18 | | #include <isl_mat_private.h> |
19 | | #include <isl_vec_private.h> |
20 | | #include <isl_options_private.h> |
21 | | #include <isl_tarjan.h> |
22 | | |
23 | | int isl_map_is_transitively_closed(__isl_keep isl_map *map) |
24 | 116 | { |
25 | 116 | isl_map *map2; |
26 | 116 | int closed; |
27 | 116 | |
28 | 116 | map2 = isl_map_apply_range(isl_map_copy(map), isl_map_copy(map)); |
29 | 116 | closed = isl_map_is_subset(map2, map); |
30 | 116 | isl_map_free(map2); |
31 | 116 | |
32 | 116 | return closed; |
33 | 116 | } |
34 | | |
35 | | int isl_union_map_is_transitively_closed(__isl_keep isl_union_map *umap) |
36 | 73 | { |
37 | 73 | isl_union_map *umap2; |
38 | 73 | int closed; |
39 | 73 | |
40 | 73 | umap2 = isl_union_map_apply_range(isl_union_map_copy(umap), |
41 | 73 | isl_union_map_copy(umap)); |
42 | 73 | closed = isl_union_map_is_subset(umap2, umap); |
43 | 73 | isl_union_map_free(umap2); |
44 | 73 | |
45 | 73 | return closed; |
46 | 73 | } |
47 | | |
48 | | /* Given a map that represents a path with the length of the path |
49 | | * encoded as the difference between the last output coordindate |
50 | | * and the last input coordinate, set this length to either |
51 | | * exactly "length" (if "exactly" is set) or at least "length" |
52 | | * (if "exactly" is not set). |
53 | | */ |
54 | | static __isl_give isl_map *set_path_length(__isl_take isl_map *map, |
55 | | int exactly, int length) |
56 | 307 | { |
57 | 307 | isl_space *dim; |
58 | 307 | struct isl_basic_map *bmap; |
59 | 307 | unsigned d; |
60 | 307 | unsigned nparam; |
61 | 307 | int k; |
62 | 307 | isl_int *c; |
63 | 307 | |
64 | 307 | if (!map) |
65 | 0 | return NULL; |
66 | 307 | |
67 | 307 | dim = isl_map_get_space(map); |
68 | 307 | d = isl_space_dim(dim, isl_dim_in); |
69 | 307 | nparam = isl_space_dim(dim, isl_dim_param); |
70 | 307 | bmap = isl_basic_map_alloc_space(dim, 0, 1, 1); |
71 | 307 | if (exactly) { |
72 | 6 | k = isl_basic_map_alloc_equality(bmap); |
73 | 6 | if (k < 0) |
74 | 0 | goto error; |
75 | 6 | c = bmap->eq[k]; |
76 | 301 | } else { |
77 | 301 | k = isl_basic_map_alloc_inequality(bmap); |
78 | 301 | if (k < 0) |
79 | 0 | goto error; |
80 | 301 | c = bmap->ineq[k]; |
81 | 301 | } |
82 | 307 | isl_seq_clr(c, 1 + isl_basic_map_total_dim(bmap)); |
83 | 307 | isl_int_set_si(c[0], -length); |
84 | 307 | isl_int_set_si(c[1 + nparam + d - 1], -1); |
85 | 307 | isl_int_set_si(c[1 + nparam + d + d - 1], 1); |
86 | 307 | |
87 | 307 | bmap = isl_basic_map_finalize(bmap); |
88 | 307 | map = isl_map_intersect(map, isl_map_from_basic_map(bmap)); |
89 | 307 | |
90 | 307 | return map; |
91 | 0 | error: |
92 | 0 | isl_basic_map_free(bmap); |
93 | 0 | isl_map_free(map); |
94 | 0 | return NULL; |
95 | 307 | } |
96 | | |
97 | | /* Check whether the overapproximation of the power of "map" is exactly |
98 | | * the power of "map". Let R be "map" and A_k the overapproximation. |
99 | | * The approximation is exact if |
100 | | * |
101 | | * A_1 = R |
102 | | * A_k = A_{k-1} \circ R k >= 2 |
103 | | * |
104 | | * Since A_k is known to be an overapproximation, we only need to check |
105 | | * |
106 | | * A_1 \subset R |
107 | | * A_k \subset A_{k-1} \circ R k >= 2 |
108 | | * |
109 | | * In practice, "app" has an extra input and output coordinate |
110 | | * to encode the length of the path. So, we first need to add |
111 | | * this coordinate to "map" and set the length of the path to |
112 | | * one. |
113 | | */ |
114 | | static int check_power_exactness(__isl_take isl_map *map, |
115 | | __isl_take isl_map *app) |
116 | 3 | { |
117 | 3 | int exact; |
118 | 3 | isl_map *app_1; |
119 | 3 | isl_map *app_2; |
120 | 3 | |
121 | 3 | map = isl_map_add_dims(map, isl_dim_in, 1); |
122 | 3 | map = isl_map_add_dims(map, isl_dim_out, 1); |
123 | 3 | map = set_path_length(map, 1, 1); |
124 | 3 | |
125 | 3 | app_1 = set_path_length(isl_map_copy(app), 1, 1); |
126 | 3 | |
127 | 3 | exact = isl_map_is_subset(app_1, map); |
128 | 3 | isl_map_free(app_1); |
129 | 3 | |
130 | 3 | if (!exact || exact < 0) { |
131 | 0 | isl_map_free(app); |
132 | 0 | isl_map_free(map); |
133 | 0 | return exact; |
134 | 0 | } |
135 | 3 | |
136 | 3 | app_1 = set_path_length(isl_map_copy(app), 0, 1); |
137 | 3 | app_2 = set_path_length(app, 0, 2); |
138 | 3 | app_1 = isl_map_apply_range(map, app_1); |
139 | 3 | |
140 | 3 | exact = isl_map_is_subset(app_2, app_1); |
141 | 3 | |
142 | 3 | isl_map_free(app_1); |
143 | 3 | isl_map_free(app_2); |
144 | 3 | |
145 | 3 | return exact; |
146 | 3 | } |
147 | | |
148 | | /* Check whether the overapproximation of the power of "map" is exactly |
149 | | * the power of "map", possibly after projecting out the power (if "project" |
150 | | * is set). |
151 | | * |
152 | | * If "project" is set and if "steps" can only result in acyclic paths, |
153 | | * then we check |
154 | | * |
155 | | * A = R \cup (A \circ R) |
156 | | * |
157 | | * where A is the overapproximation with the power projected out, i.e., |
158 | | * an overapproximation of the transitive closure. |
159 | | * More specifically, since A is known to be an overapproximation, we check |
160 | | * |
161 | | * A \subset R \cup (A \circ R) |
162 | | * |
163 | | * Otherwise, we check if the power is exact. |
164 | | * |
165 | | * Note that "app" has an extra input and output coordinate to encode |
166 | | * the length of the part. If we are only interested in the transitive |
167 | | * closure, then we can simply project out these coordinates first. |
168 | | */ |
169 | | static int check_exactness(__isl_take isl_map *map, __isl_take isl_map *app, |
170 | | int project) |
171 | 102 | { |
172 | 102 | isl_map *test; |
173 | 102 | int exact; |
174 | 102 | unsigned d; |
175 | 102 | |
176 | 102 | if (!project) |
177 | 3 | return check_power_exactness(map, app); |
178 | 99 | |
179 | 99 | d = isl_map_dim(map, isl_dim_in); |
180 | 99 | app = set_path_length(app, 0, 1); |
181 | 99 | app = isl_map_project_out(app, isl_dim_in, d, 1); |
182 | 99 | app = isl_map_project_out(app, isl_dim_out, d, 1); |
183 | 99 | |
184 | 99 | app = isl_map_reset_space(app, isl_map_get_space(map)); |
185 | 99 | |
186 | 99 | test = isl_map_apply_range(isl_map_copy(map), isl_map_copy(app)); |
187 | 99 | test = isl_map_union(test, isl_map_copy(map)); |
188 | 99 | |
189 | 99 | exact = isl_map_is_subset(app, test); |
190 | 99 | |
191 | 99 | isl_map_free(app); |
192 | 99 | isl_map_free(test); |
193 | 99 | |
194 | 99 | isl_map_free(map); |
195 | 99 | |
196 | 99 | return exact; |
197 | 99 | } |
198 | | |
199 | | /* |
200 | | * The transitive closure implementation is based on the paper |
201 | | * "Computing the Transitive Closure of a Union of Affine Integer |
202 | | * Tuple Relations" by Anna Beletska, Denis Barthou, Wlodzimierz Bielecki and |
203 | | * Albert Cohen. |
204 | | */ |
205 | | |
206 | | /* Given a set of n offsets v_i (the rows of "steps"), construct a relation |
207 | | * of the given dimension specification (Z^{n+1} -> Z^{n+1}) |
208 | | * that maps an element x to any element that can be reached |
209 | | * by taking a non-negative number of steps along any of |
210 | | * the extended offsets v'_i = [v_i 1]. |
211 | | * That is, construct |
212 | | * |
213 | | * { [x] -> [y] : exists k_i >= 0, y = x + \sum_i k_i v'_i } |
214 | | * |
215 | | * For any element in this relation, the number of steps taken |
216 | | * is equal to the difference in the final coordinates. |
217 | | */ |
218 | | static __isl_give isl_map *path_along_steps(__isl_take isl_space *dim, |
219 | | __isl_keep isl_mat *steps) |
220 | 188 | { |
221 | 188 | int i, j, k; |
222 | 188 | struct isl_basic_map *path = NULL; |
223 | 188 | unsigned d; |
224 | 188 | unsigned n; |
225 | 188 | unsigned nparam; |
226 | 188 | |
227 | 188 | if (!dim || !steps) |
228 | 0 | goto error; |
229 | 188 | |
230 | 188 | d = isl_space_dim(dim, isl_dim_in); |
231 | 188 | n = steps->n_row; |
232 | 188 | nparam = isl_space_dim(dim, isl_dim_param); |
233 | 188 | |
234 | 188 | path = isl_basic_map_alloc_space(isl_space_copy(dim), n, d, n); |
235 | 188 | |
236 | 398 | for (i = 0; i < n; ++i210 ) { |
237 | 210 | k = isl_basic_map_alloc_div(path); |
238 | 210 | if (k < 0) |
239 | 0 | goto error; |
240 | 210 | isl_assert(steps->ctx, i == k, goto error); |
241 | 210 | isl_int_set_si(path->div[k][0], 0); |
242 | 210 | } |
243 | 188 | |
244 | 915 | for (i = 0; 188 i < d; ++i727 ) { |
245 | 727 | k = isl_basic_map_alloc_equality(path); |
246 | 727 | if (k < 0) |
247 | 0 | goto error; |
248 | 727 | isl_seq_clr(path->eq[k], 1 + isl_basic_map_total_dim(path)); |
249 | 727 | isl_int_set_si(path->eq[k][1 + nparam + i], 1); |
250 | 727 | isl_int_set_si(path->eq[k][1 + nparam + d + i], -1); |
251 | 727 | if (i == d - 1) |
252 | 398 | for (j = 0; 188 j < n; ++j210 ) |
253 | 210 | isl_int_set_si(path->eq[k][1 + nparam + 2 * d + j], 1); |
254 | 727 | else |
255 | 1.14k | for (j = 0; 539 j < n; ++j606 ) |
256 | 606 | isl_int_set(path->eq[k][1 + nparam + 2 * d + j], |
257 | 727 | steps->row[j][i]); |
258 | 727 | } |
259 | 188 | |
260 | 398 | for (i = 0; 188 i < n; ++i210 ) { |
261 | 210 | k = isl_basic_map_alloc_inequality(path); |
262 | 210 | if (k < 0) |
263 | 0 | goto error; |
264 | 210 | isl_seq_clr(path->ineq[k], 1 + isl_basic_map_total_dim(path)); |
265 | 210 | isl_int_set_si(path->ineq[k][1 + nparam + 2 * d + i], 1); |
266 | 210 | } |
267 | 188 | |
268 | 188 | isl_space_free(dim); |
269 | 188 | |
270 | 188 | path = isl_basic_map_simplify(path); |
271 | 188 | path = isl_basic_map_finalize(path); |
272 | 188 | return isl_map_from_basic_map(path); |
273 | 0 | error: |
274 | 0 | isl_space_free(dim); |
275 | 0 | isl_basic_map_free(path); |
276 | 0 | return NULL; |
277 | 188 | } |
278 | | |
279 | 74 | #define IMPURE 0 |
280 | 111 | #define PURE_PARAM 1 |
281 | 278 | #define PURE_VAR 2 |
282 | 24 | #define MIXED 3 |
283 | | |
284 | | /* Check whether the parametric constant term of constraint c is never |
285 | | * positive in "bset". |
286 | | */ |
287 | | static isl_bool parametric_constant_never_positive( |
288 | | __isl_keep isl_basic_set *bset, isl_int *c, int *div_purity) |
289 | 10 | { |
290 | 10 | unsigned d; |
291 | 10 | unsigned n_div; |
292 | 10 | unsigned nparam; |
293 | 10 | int i; |
294 | 10 | int k; |
295 | 10 | isl_bool empty; |
296 | 10 | |
297 | 10 | n_div = isl_basic_set_dim(bset, isl_dim_div); |
298 | 10 | d = isl_basic_set_dim(bset, isl_dim_set); |
299 | 10 | nparam = isl_basic_set_dim(bset, isl_dim_param); |
300 | 10 | |
301 | 10 | bset = isl_basic_set_copy(bset); |
302 | 10 | bset = isl_basic_set_cow(bset); |
303 | 10 | bset = isl_basic_set_extend_constraints(bset, 0, 1); |
304 | 10 | k = isl_basic_set_alloc_inequality(bset); |
305 | 10 | if (k < 0) |
306 | 0 | goto error; |
307 | 10 | isl_seq_clr(bset->ineq[k], 1 + isl_basic_set_total_dim(bset)); |
308 | 10 | isl_seq_cpy(bset->ineq[k], c, 1 + nparam); |
309 | 15 | for (i = 0; i < n_div; ++i5 ) { |
310 | 5 | if (div_purity[i] != PURE_PARAM) |
311 | 5 | continue1 ; |
312 | 4 | isl_int_set(bset->ineq[k][1 + nparam + d + i], |
313 | 4 | c[1 + nparam + d + i]); |
314 | 4 | } |
315 | 10 | isl_int_sub_ui(bset->ineq[k][0], bset->ineq[k][0], 1); |
316 | 10 | empty = isl_basic_set_is_empty(bset); |
317 | 10 | isl_basic_set_free(bset); |
318 | 10 | |
319 | 10 | return empty; |
320 | 0 | error: |
321 | 0 | isl_basic_set_free(bset); |
322 | 0 | return isl_bool_error; |
323 | 10 | } |
324 | | |
325 | | /* Return PURE_PARAM if only the coefficients of the parameters are non-zero. |
326 | | * Return PURE_VAR if only the coefficients of the set variables are non-zero. |
327 | | * Return MIXED if only the coefficients of the parameters and the set |
328 | | * variables are non-zero and if moreover the parametric constant |
329 | | * can never attain positive values. |
330 | | * Return IMPURE otherwise. |
331 | | */ |
332 | | static int purity(__isl_keep isl_basic_set *bset, isl_int *c, int *div_purity, |
333 | | int eq) |
334 | 42 | { |
335 | 42 | unsigned d; |
336 | 42 | unsigned n_div; |
337 | 42 | unsigned nparam; |
338 | 42 | isl_bool empty; |
339 | 42 | int i; |
340 | 42 | int p = 0, v = 0; |
341 | 42 | |
342 | 42 | n_div = isl_basic_set_dim(bset, isl_dim_div); |
343 | 42 | d = isl_basic_set_dim(bset, isl_dim_set); |
344 | 42 | nparam = isl_basic_set_dim(bset, isl_dim_param); |
345 | 42 | |
346 | 76 | for (i = 0; i < n_div; ++i34 ) { |
347 | 34 | if (isl_int_is_zero(c[1 + nparam + d + i])) |
348 | 34 | continue29 ; |
349 | 5 | switch (div_purity[i]) { |
350 | 5 | case 4 PURE_PARAM4 : p = 1; break; |
351 | 5 | case 1 PURE_VAR1 : v = 1; break; |
352 | 5 | default: return 0 IMPURE0 ; |
353 | 5 | } |
354 | 5 | } |
355 | 42 | if (!p && isl_seq_first_non_zero(c + 1, nparam) == -138 ) |
356 | 25 | return PURE_VAR; |
357 | 17 | if (!v && isl_seq_first_non_zero(c + 1 + nparam, d) == -1) |
358 | 7 | return PURE_PARAM; |
359 | 10 | |
360 | 10 | empty = parametric_constant_never_positive(bset, c, div_purity); |
361 | 10 | if (eq && empty >= 00 && !empty0 ) { |
362 | 0 | isl_seq_neg(c, c, 1 + nparam + d + n_div); |
363 | 0 | empty = parametric_constant_never_positive(bset, c, div_purity); |
364 | 0 | } |
365 | 10 | |
366 | 10 | return empty < 0 ? -10 : empty ? MIXED1 : IMPURE9 ; |
367 | 10 | } |
368 | | |
369 | | /* Return an array of integers indicating the type of each div in bset. |
370 | | * If the div is (recursively) defined in terms of only the parameters, |
371 | | * then the type is PURE_PARAM. |
372 | | * If the div is (recursively) defined in terms of only the set variables, |
373 | | * then the type is PURE_VAR. |
374 | | * Otherwise, the type is IMPURE. |
375 | | */ |
376 | | static __isl_give int *get_div_purity(__isl_keep isl_basic_set *bset) |
377 | 8 | { |
378 | 8 | int i, j; |
379 | 8 | int *div_purity; |
380 | 8 | unsigned d; |
381 | 8 | unsigned n_div; |
382 | 8 | unsigned nparam; |
383 | 8 | |
384 | 8 | if (!bset) |
385 | 0 | return NULL; |
386 | 8 | |
387 | 8 | n_div = isl_basic_set_dim(bset, isl_dim_div); |
388 | 8 | d = isl_basic_set_dim(bset, isl_dim_set); |
389 | 8 | nparam = isl_basic_set_dim(bset, isl_dim_param); |
390 | 8 | |
391 | 8 | div_purity = isl_alloc_array(bset->ctx, int, n_div); |
392 | 8 | if (n_div && !div_purity2 ) |
393 | 0 | return NULL; |
394 | 8 | |
395 | 11 | for (i = 0; 8 i < bset->n_div; ++i3 ) { |
396 | 3 | int p = 0, v = 0; |
397 | 3 | if (isl_int_is_zero(bset->div[i][0])) { |
398 | 0 | div_purity[i] = IMPURE; |
399 | 0 | continue; |
400 | 0 | } |
401 | 3 | if (isl_seq_first_non_zero(bset->div[i] + 2, nparam) != -1) |
402 | 2 | p = 1; |
403 | 3 | if (isl_seq_first_non_zero(bset->div[i] + 2 + nparam, d) != -1) |
404 | 1 | v = 1; |
405 | 4 | for (j = 0; j < i; ++j1 ) { |
406 | 1 | if (isl_int_is_zero(bset->div[i][2 + nparam + d + j])) |
407 | 1 | continue; |
408 | 0 | switch (div_purity[j]) { |
409 | 0 | case PURE_PARAM: p = 1; break; |
410 | 0 | case PURE_VAR: v = 1; break; |
411 | 0 | default: p = v = 1; break; |
412 | 0 | } |
413 | 0 | } |
414 | 3 | div_purity[i] = v ? p ? 1 IMPURE0 : PURE_VAR1 : PURE_PARAM2 ; |
415 | 3 | } |
416 | 8 | |
417 | 8 | return div_purity; |
418 | 8 | } |
419 | | |
420 | | /* Given a path with the as yet unconstrained length at position "pos", |
421 | | * check if setting the length to zero results in only the identity |
422 | | * mapping. |
423 | | */ |
424 | | static int empty_path_is_identity(__isl_keep isl_basic_map *path, unsigned pos) |
425 | 8 | { |
426 | 8 | isl_basic_map *test = NULL; |
427 | 8 | isl_basic_map *id = NULL; |
428 | 8 | int k; |
429 | 8 | int is_id; |
430 | 8 | |
431 | 8 | test = isl_basic_map_copy(path); |
432 | 8 | test = isl_basic_map_extend_constraints(test, 1, 0); |
433 | 8 | k = isl_basic_map_alloc_equality(test); |
434 | 8 | if (k < 0) |
435 | 0 | goto error; |
436 | 8 | isl_seq_clr(test->eq[k], 1 + isl_basic_map_total_dim(test)); |
437 | 8 | isl_int_set_si(test->eq[k][pos], 1); |
438 | 8 | test = isl_basic_map_gauss(test, NULL); |
439 | 8 | id = isl_basic_map_identity(isl_basic_map_get_space(path)); |
440 | 8 | is_id = isl_basic_map_is_equal(test, id); |
441 | 8 | isl_basic_map_free(test); |
442 | 8 | isl_basic_map_free(id); |
443 | 8 | return is_id; |
444 | 0 | error: |
445 | 0 | isl_basic_map_free(test); |
446 | 0 | return -1; |
447 | 8 | } |
448 | | |
449 | | /* If any of the constraints is found to be impure then this function |
450 | | * sets *impurity to 1. |
451 | | * |
452 | | * If impurity is NULL then we are dealing with a non-parametric set |
453 | | * and so the constraints are obviously PURE_VAR. |
454 | | */ |
455 | | static __isl_give isl_basic_map *add_delta_constraints( |
456 | | __isl_take isl_basic_map *path, |
457 | | __isl_keep isl_basic_set *delta, unsigned off, unsigned nparam, |
458 | | unsigned d, int *div_purity, int eq, int *impurity) |
459 | 30 | { |
460 | 30 | int i, k; |
461 | 30 | int n = eq ? delta->n_eq15 : delta->n_ineq15 ; |
462 | 30 | isl_int **delta_c = eq ? delta->eq15 : delta->ineq15 ; |
463 | 30 | unsigned n_div; |
464 | 30 | |
465 | 30 | n_div = isl_basic_set_dim(delta, isl_dim_div); |
466 | 30 | |
467 | 95 | for (i = 0; i < n; ++i65 ) { |
468 | 65 | isl_int *path_c; |
469 | 65 | int p = PURE_VAR; |
470 | 65 | if (impurity) |
471 | 42 | p = purity(delta, delta_c[i], div_purity, eq); |
472 | 65 | if (p < 0) |
473 | 0 | goto error; |
474 | 65 | if (p != PURE_VAR && p != 17 PURE_PARAM17 && !*impurity10 ) |
475 | 7 | *impurity = 1; |
476 | 65 | if (p == IMPURE) |
477 | 65 | continue9 ; |
478 | 56 | if (eq && p != 23 MIXED23 ) { |
479 | 23 | k = isl_basic_map_alloc_equality(path); |
480 | 23 | if (k < 0) |
481 | 0 | goto error; |
482 | 23 | path_c = path->eq[k]; |
483 | 33 | } else { |
484 | 33 | k = isl_basic_map_alloc_inequality(path); |
485 | 33 | if (k < 0) |
486 | 0 | goto error; |
487 | 33 | path_c = path->ineq[k]; |
488 | 33 | } |
489 | 56 | isl_seq_clr(path_c, 1 + isl_basic_map_total_dim(path)); |
490 | 56 | if (p == PURE_VAR) { |
491 | 48 | isl_seq_cpy(path_c + off, |
492 | 48 | delta_c[i] + 1 + nparam, d); |
493 | 48 | isl_int_set(path_c[off + d], delta_c[i][0]); |
494 | 48 | } else if (8 p == 8 PURE_PARAM8 ) { |
495 | 7 | isl_seq_cpy(path_c, delta_c[i], 1 + nparam); |
496 | 7 | } else { |
497 | 1 | isl_seq_cpy(path_c + off, |
498 | 1 | delta_c[i] + 1 + nparam, d); |
499 | 1 | isl_seq_cpy(path_c, delta_c[i], 1 + nparam); |
500 | 1 | } |
501 | 56 | isl_seq_cpy(path_c + off - n_div, |
502 | 56 | delta_c[i] + 1 + nparam + d, n_div); |
503 | 56 | } |
504 | 30 | |
505 | 30 | return path; |
506 | 0 | error: |
507 | 0 | isl_basic_map_free(path); |
508 | 0 | return NULL; |
509 | 30 | } |
510 | | |
511 | | /* Given a set of offsets "delta", construct a relation of the |
512 | | * given dimension specification (Z^{n+1} -> Z^{n+1}) that |
513 | | * is an overapproximation of the relations that |
514 | | * maps an element x to any element that can be reached |
515 | | * by taking a non-negative number of steps along any of |
516 | | * the elements in "delta". |
517 | | * That is, construct an approximation of |
518 | | * |
519 | | * { [x] -> [y] : exists f \in \delta, k \in Z : |
520 | | * y = x + k [f, 1] and k >= 0 } |
521 | | * |
522 | | * For any element in this relation, the number of steps taken |
523 | | * is equal to the difference in the final coordinates. |
524 | | * |
525 | | * In particular, let delta be defined as |
526 | | * |
527 | | * \delta = [p] -> { [x] : A x + a >= 0 and B p + b >= 0 and |
528 | | * C x + C'p + c >= 0 and |
529 | | * D x + D'p + d >= 0 } |
530 | | * |
531 | | * where the constraints C x + C'p + c >= 0 are such that the parametric |
532 | | * constant term of each constraint j, "C_j x + C'_j p + c_j", |
533 | | * can never attain positive values, then the relation is constructed as |
534 | | * |
535 | | * { [x] -> [y] : exists [f, k] \in Z^{n+1} : y = x + f and |
536 | | * A f + k a >= 0 and B p + b >= 0 and |
537 | | * C f + C'p + c >= 0 and k >= 1 } |
538 | | * union { [x] -> [x] } |
539 | | * |
540 | | * If the zero-length paths happen to correspond exactly to the identity |
541 | | * mapping, then we return |
542 | | * |
543 | | * { [x] -> [y] : exists [f, k] \in Z^{n+1} : y = x + f and |
544 | | * A f + k a >= 0 and B p + b >= 0 and |
545 | | * C f + C'p + c >= 0 and k >= 0 } |
546 | | * |
547 | | * instead. |
548 | | * |
549 | | * Existentially quantified variables in \delta are handled by |
550 | | * classifying them as independent of the parameters, purely |
551 | | * parameter dependent and others. Constraints containing |
552 | | * any of the other existentially quantified variables are removed. |
553 | | * This is safe, but leads to an additional overapproximation. |
554 | | * |
555 | | * If there are any impure constraints, then we also eliminate |
556 | | * the parameters from \delta, resulting in a set |
557 | | * |
558 | | * \delta' = { [x] : E x + e >= 0 } |
559 | | * |
560 | | * and add the constraints |
561 | | * |
562 | | * E f + k e >= 0 |
563 | | * |
564 | | * to the constructed relation. |
565 | | */ |
566 | | static __isl_give isl_map *path_along_delta(__isl_take isl_space *dim, |
567 | | __isl_take isl_basic_set *delta) |
568 | 8 | { |
569 | 8 | isl_basic_map *path = NULL; |
570 | 8 | unsigned d; |
571 | 8 | unsigned n_div; |
572 | 8 | unsigned nparam; |
573 | 8 | unsigned off; |
574 | 8 | int i, k; |
575 | 8 | int is_id; |
576 | 8 | int *div_purity = NULL; |
577 | 8 | int impurity = 0; |
578 | 8 | |
579 | 8 | if (!delta) |
580 | 0 | goto error; |
581 | 8 | n_div = isl_basic_set_dim(delta, isl_dim_div); |
582 | 8 | d = isl_basic_set_dim(delta, isl_dim_set); |
583 | 8 | nparam = isl_basic_set_dim(delta, isl_dim_param); |
584 | 8 | path = isl_basic_map_alloc_space(isl_space_copy(dim), n_div + d + 1, |
585 | 8 | d + 1 + delta->n_eq, delta->n_eq + delta->n_ineq + 1); |
586 | 8 | off = 1 + nparam + 2 * (d + 1) + n_div; |
587 | 8 | |
588 | 39 | for (i = 0; i < n_div + d + 1; ++i31 ) { |
589 | 31 | k = isl_basic_map_alloc_div(path); |
590 | 31 | if (k < 0) |
591 | 0 | goto error; |
592 | 31 | isl_int_set_si(path->div[k][0], 0); |
593 | 31 | } |
594 | 8 | |
595 | 36 | for (i = 0; 8 i < d + 1; ++i28 ) { |
596 | 28 | k = isl_basic_map_alloc_equality(path); |
597 | 28 | if (k < 0) |
598 | 0 | goto error; |
599 | 28 | isl_seq_clr(path->eq[k], 1 + isl_basic_map_total_dim(path)); |
600 | 28 | isl_int_set_si(path->eq[k][1 + nparam + i], 1); |
601 | 28 | isl_int_set_si(path->eq[k][1 + nparam + d + 1 + i], -1); |
602 | 28 | isl_int_set_si(path->eq[k][off + i], 1); |
603 | 28 | } |
604 | 8 | |
605 | 8 | div_purity = get_div_purity(delta); |
606 | 8 | if (n_div && !div_purity2 ) |
607 | 0 | goto error; |
608 | 8 | |
609 | 8 | path = add_delta_constraints(path, delta, off, nparam, d, |
610 | 8 | div_purity, 1, &impurity); |
611 | 8 | path = add_delta_constraints(path, delta, off, nparam, d, |
612 | 8 | div_purity, 0, &impurity); |
613 | 8 | if (impurity) { |
614 | 7 | isl_space *dim = isl_basic_set_get_space(delta); |
615 | 7 | delta = isl_basic_set_project_out(delta, |
616 | 7 | isl_dim_param, 0, nparam); |
617 | 7 | delta = isl_basic_set_add_dims(delta, isl_dim_param, nparam); |
618 | 7 | delta = isl_basic_set_reset_space(delta, dim); |
619 | 7 | if (!delta) |
620 | 0 | goto error; |
621 | 7 | path = isl_basic_map_extend_constraints(path, delta->n_eq, |
622 | 7 | delta->n_ineq + 1); |
623 | 7 | path = add_delta_constraints(path, delta, off, nparam, d, |
624 | 7 | NULL, 1, NULL); |
625 | 7 | path = add_delta_constraints(path, delta, off, nparam, d, |
626 | 7 | NULL, 0, NULL); |
627 | 7 | path = isl_basic_map_gauss(path, NULL); |
628 | 7 | } |
629 | 8 | |
630 | 8 | is_id = empty_path_is_identity(path, off + d); |
631 | 8 | if (is_id < 0) |
632 | 0 | goto error; |
633 | 8 | |
634 | 8 | k = isl_basic_map_alloc_inequality(path); |
635 | 8 | if (k < 0) |
636 | 0 | goto error; |
637 | 8 | isl_seq_clr(path->ineq[k], 1 + isl_basic_map_total_dim(path)); |
638 | 8 | if (!is_id) |
639 | 8 | isl_int_set_si6 (path->ineq[k][0], -1); |
640 | 8 | isl_int_set_si(path->ineq[k][off + d], 1); |
641 | 8 | |
642 | 8 | free(div_purity); |
643 | 8 | isl_basic_set_free(delta); |
644 | 8 | path = isl_basic_map_finalize(path); |
645 | 8 | if (is_id) { |
646 | 2 | isl_space_free(dim); |
647 | 2 | return isl_map_from_basic_map(path); |
648 | 2 | } |
649 | 6 | return isl_basic_map_union(path, isl_basic_map_identity(dim)); |
650 | 0 | error: |
651 | 0 | free(div_purity); |
652 | 0 | isl_space_free(dim); |
653 | 0 | isl_basic_set_free(delta); |
654 | 0 | isl_basic_map_free(path); |
655 | 0 | return NULL; |
656 | 6 | } |
657 | | |
658 | | /* Given a dimension specification Z^{n+1} -> Z^{n+1} and a parameter "param", |
659 | | * construct a map that equates the parameter to the difference |
660 | | * in the final coordinates and imposes that this difference is positive. |
661 | | * That is, construct |
662 | | * |
663 | | * { [x,x_s] -> [y,y_s] : k = y_s - x_s > 0 } |
664 | | */ |
665 | | static __isl_give isl_map *equate_parameter_to_length(__isl_take isl_space *dim, |
666 | | unsigned param) |
667 | 3 | { |
668 | 3 | struct isl_basic_map *bmap; |
669 | 3 | unsigned d; |
670 | 3 | unsigned nparam; |
671 | 3 | int k; |
672 | 3 | |
673 | 3 | d = isl_space_dim(dim, isl_dim_in); |
674 | 3 | nparam = isl_space_dim(dim, isl_dim_param); |
675 | 3 | bmap = isl_basic_map_alloc_space(dim, 0, 1, 1); |
676 | 3 | k = isl_basic_map_alloc_equality(bmap); |
677 | 3 | if (k < 0) |
678 | 0 | goto error; |
679 | 3 | isl_seq_clr(bmap->eq[k], 1 + isl_basic_map_total_dim(bmap)); |
680 | 3 | isl_int_set_si(bmap->eq[k][1 + param], -1); |
681 | 3 | isl_int_set_si(bmap->eq[k][1 + nparam + d - 1], -1); |
682 | 3 | isl_int_set_si(bmap->eq[k][1 + nparam + d + d - 1], 1); |
683 | 3 | |
684 | 3 | k = isl_basic_map_alloc_inequality(bmap); |
685 | 3 | if (k < 0) |
686 | 0 | goto error; |
687 | 3 | isl_seq_clr(bmap->ineq[k], 1 + isl_basic_map_total_dim(bmap)); |
688 | 3 | isl_int_set_si(bmap->ineq[k][1 + param], 1); |
689 | 3 | isl_int_set_si(bmap->ineq[k][0], -1); |
690 | 3 | |
691 | 3 | bmap = isl_basic_map_finalize(bmap); |
692 | 3 | return isl_map_from_basic_map(bmap); |
693 | 0 | error: |
694 | 0 | isl_basic_map_free(bmap); |
695 | 0 | return NULL; |
696 | 3 | } |
697 | | |
698 | | /* Check whether "path" is acyclic, where the last coordinates of domain |
699 | | * and range of path encode the number of steps taken. |
700 | | * That is, check whether |
701 | | * |
702 | | * { d | d = y - x and (x,y) in path } |
703 | | * |
704 | | * does not contain any element with positive last coordinate (positive length) |
705 | | * and zero remaining coordinates (cycle). |
706 | | */ |
707 | | static int is_acyclic(__isl_take isl_map *path) |
708 | 99 | { |
709 | 99 | int i; |
710 | 99 | int acyclic; |
711 | 99 | unsigned dim; |
712 | 99 | struct isl_set *delta; |
713 | 99 | |
714 | 99 | delta = isl_map_deltas(path); |
715 | 99 | dim = isl_set_dim(delta, isl_dim_set); |
716 | 486 | for (i = 0; i < dim; ++i387 ) { |
717 | 387 | if (i == dim -1) |
718 | 99 | delta = isl_set_lower_bound_si(delta, isl_dim_set, i, 1); |
719 | 288 | else |
720 | 288 | delta = isl_set_fix_si(delta, isl_dim_set, i, 0); |
721 | 387 | } |
722 | 99 | |
723 | 99 | acyclic = isl_set_is_empty(delta); |
724 | 99 | isl_set_free(delta); |
725 | 99 | |
726 | 99 | return acyclic; |
727 | 99 | } |
728 | | |
729 | | /* Given a union of basic maps R = \cup_i R_i \subseteq D \times D |
730 | | * and a dimension specification (Z^{n+1} -> Z^{n+1}), |
731 | | * construct a map that is an overapproximation of the map |
732 | | * that takes an element from the space D \times Z to another |
733 | | * element from the same space, such that the first n coordinates of the |
734 | | * difference between them is a sum of differences between images |
735 | | * and pre-images in one of the R_i and such that the last coordinate |
736 | | * is equal to the number of steps taken. |
737 | | * That is, let |
738 | | * |
739 | | * \Delta_i = { y - x | (x, y) in R_i } |
740 | | * |
741 | | * then the constructed map is an overapproximation of |
742 | | * |
743 | | * { (x) -> (x + d) | \exists k_i >= 0, \delta_i \in \Delta_i : |
744 | | * d = (\sum_i k_i \delta_i, \sum_i k_i) } |
745 | | * |
746 | | * The elements of the singleton \Delta_i's are collected as the |
747 | | * rows of the steps matrix. For all these \Delta_i's together, |
748 | | * a single path is constructed. |
749 | | * For each of the other \Delta_i's, we compute an overapproximation |
750 | | * of the paths along elements of \Delta_i. |
751 | | * Since each of these paths performs an addition, composition is |
752 | | * symmetric and we can simply compose all resulting paths in any order. |
753 | | */ |
754 | | static __isl_give isl_map *construct_extended_path(__isl_take isl_space *dim, |
755 | | __isl_keep isl_map *map, int *project) |
756 | 196 | { |
757 | 196 | struct isl_mat *steps = NULL; |
758 | 196 | struct isl_map *path = NULL; |
759 | 196 | unsigned d; |
760 | 196 | int i, j, n; |
761 | 196 | |
762 | 196 | if (!map) |
763 | 0 | goto error; |
764 | 196 | |
765 | 196 | d = isl_map_dim(map, isl_dim_in); |
766 | 196 | |
767 | 196 | path = isl_map_identity(isl_space_copy(dim)); |
768 | 196 | |
769 | 196 | steps = isl_mat_alloc(map->ctx, map->n, d); |
770 | 196 | if (!steps) |
771 | 0 | goto error; |
772 | 196 | |
773 | 196 | n = 0; |
774 | 414 | for (i = 0; i < map->n; ++i218 ) { |
775 | 218 | struct isl_basic_set *delta; |
776 | 218 | |
777 | 218 | delta = isl_basic_map_deltas(isl_basic_map_copy(map->p[i])); |
778 | 218 | |
779 | 827 | for (j = 0; j < d; ++j609 ) { |
780 | 617 | isl_bool fixed; |
781 | 617 | |
782 | 617 | fixed = isl_basic_set_plain_dim_is_fixed(delta, j, |
783 | 617 | &steps->row[n][j]); |
784 | 617 | if (fixed < 0) { |
785 | 0 | isl_basic_set_free(delta); |
786 | 0 | goto error; |
787 | 0 | } |
788 | 617 | if (!fixed) |
789 | 8 | break; |
790 | 617 | } |
791 | 218 | |
792 | 218 | |
793 | 218 | if (j < d) { |
794 | 8 | path = isl_map_apply_range(path, |
795 | 8 | path_along_delta(isl_space_copy(dim), delta)); |
796 | 8 | path = isl_map_coalesce(path); |
797 | 210 | } else { |
798 | 210 | isl_basic_set_free(delta); |
799 | 210 | ++n; |
800 | 210 | } |
801 | 218 | } |
802 | 196 | |
803 | 196 | if (n > 0) { |
804 | 188 | steps->n_row = n; |
805 | 188 | path = isl_map_apply_range(path, |
806 | 188 | path_along_steps(isl_space_copy(dim), steps)); |
807 | 188 | } |
808 | 196 | |
809 | 196 | if (project && *project102 ) { |
810 | 99 | *project = is_acyclic(isl_map_copy(path)); |
811 | 99 | if (*project < 0) |
812 | 0 | goto error; |
813 | 196 | } |
814 | 196 | |
815 | 196 | isl_space_free(dim); |
816 | 196 | isl_mat_free(steps); |
817 | 196 | return path; |
818 | 0 | error: |
819 | 0 | isl_space_free(dim); |
820 | 0 | isl_mat_free(steps); |
821 | 0 | isl_map_free(path); |
822 | 0 | return NULL; |
823 | 196 | } |
824 | | |
825 | | static isl_bool isl_set_overlaps(__isl_keep isl_set *set1, |
826 | | __isl_keep isl_set *set2) |
827 | 410 | { |
828 | 410 | isl_set *i; |
829 | 410 | isl_bool no_overlap; |
830 | 410 | |
831 | 410 | if (!set1 || !set2) |
832 | 0 | return isl_bool_error; |
833 | 410 | |
834 | 410 | if (!isl_space_tuple_is_equal(set1->dim, isl_dim_set, |
835 | 410 | set2->dim, isl_dim_set)) |
836 | 0 | return isl_bool_false; |
837 | 410 | |
838 | 410 | i = isl_set_intersect(isl_set_copy(set1), isl_set_copy(set2)); |
839 | 410 | no_overlap = isl_set_is_empty(i); |
840 | 410 | isl_set_free(i); |
841 | 410 | |
842 | 410 | return isl_bool_not(no_overlap); |
843 | 410 | } |
844 | | |
845 | | /* Given a union of basic maps R = \cup_i R_i \subseteq D \times D |
846 | | * and a dimension specification (Z^{n+1} -> Z^{n+1}), |
847 | | * construct a map that is an overapproximation of the map |
848 | | * that takes an element from the dom R \times Z to an |
849 | | * element from ran R \times Z, such that the first n coordinates of the |
850 | | * difference between them is a sum of differences between images |
851 | | * and pre-images in one of the R_i and such that the last coordinate |
852 | | * is equal to the number of steps taken. |
853 | | * That is, let |
854 | | * |
855 | | * \Delta_i = { y - x | (x, y) in R_i } |
856 | | * |
857 | | * then the constructed map is an overapproximation of |
858 | | * |
859 | | * { (x) -> (x + d) | \exists k_i >= 0, \delta_i \in \Delta_i : |
860 | | * d = (\sum_i k_i \delta_i, \sum_i k_i) and |
861 | | * x in dom R and x + d in ran R and |
862 | | * \sum_i k_i >= 1 } |
863 | | */ |
864 | | static __isl_give isl_map *construct_component(__isl_take isl_space *dim, |
865 | | __isl_keep isl_map *map, int *exact, int project) |
866 | 113 | { |
867 | 113 | struct isl_set *domain = NULL; |
868 | 113 | struct isl_set *range = NULL; |
869 | 113 | struct isl_map *app = NULL; |
870 | 113 | struct isl_map *path = NULL; |
871 | 113 | isl_bool overlaps; |
872 | 113 | |
873 | 113 | domain = isl_map_domain(isl_map_copy(map)); |
874 | 113 | domain = isl_set_coalesce(domain); |
875 | 113 | range = isl_map_range(isl_map_copy(map)); |
876 | 113 | range = isl_set_coalesce(range); |
877 | 113 | overlaps = isl_set_overlaps(domain, range); |
878 | 113 | if (overlaps < 0 || !overlaps) { |
879 | 0 | isl_set_free(domain); |
880 | 0 | isl_set_free(range); |
881 | 0 | isl_space_free(dim); |
882 | 0 |
|
883 | 0 | if (overlaps < 0) |
884 | 0 | map = NULL; |
885 | 0 | map = isl_map_copy(map); |
886 | 0 | map = isl_map_add_dims(map, isl_dim_in, 1); |
887 | 0 | map = isl_map_add_dims(map, isl_dim_out, 1); |
888 | 0 | map = set_path_length(map, 1, 1); |
889 | 0 | return map; |
890 | 0 | } |
891 | 113 | app = isl_map_from_domain_and_range(domain, range); |
892 | 113 | app = isl_map_add_dims(app, isl_dim_in, 1); |
893 | 113 | app = isl_map_add_dims(app, isl_dim_out, 1); |
894 | 113 | |
895 | 113 | path = construct_extended_path(isl_space_copy(dim), map, |
896 | 113 | exact && *exact19 ? &project19 : NULL); |
897 | 113 | app = isl_map_intersect(app, path); |
898 | 113 | |
899 | 113 | if (exact && *exact19 && |
900 | 113 | (*exact = check_exactness(isl_map_copy(map), isl_map_copy(app), |
901 | 19 | project)) < 0) |
902 | 0 | goto error; |
903 | 113 | |
904 | 113 | isl_space_free(dim); |
905 | 113 | app = set_path_length(app, 0, 1); |
906 | 113 | return app; |
907 | 0 | error: |
908 | 0 | isl_space_free(dim); |
909 | 0 | isl_map_free(app); |
910 | 0 | return NULL; |
911 | 113 | } |
912 | | |
913 | | /* Call construct_component and, if "project" is set, project out |
914 | | * the final coordinates. |
915 | | */ |
916 | | static __isl_give isl_map *construct_projected_component( |
917 | | __isl_take isl_space *dim, |
918 | | __isl_keep isl_map *map, int *exact, int project) |
919 | 113 | { |
920 | 113 | isl_map *app; |
921 | 113 | unsigned d; |
922 | 113 | |
923 | 113 | if (!dim) |
924 | 0 | return NULL; |
925 | 113 | d = isl_space_dim(dim, isl_dim_in); |
926 | 113 | |
927 | 113 | app = construct_component(dim, map, exact, project); |
928 | 113 | if (project) { |
929 | 110 | app = isl_map_project_out(app, isl_dim_in, d - 1, 1); |
930 | 110 | app = isl_map_project_out(app, isl_dim_out, d - 1, 1); |
931 | 110 | } |
932 | 113 | return app; |
933 | 113 | } |
934 | | |
935 | | /* Compute an extended version, i.e., with path lengths, of |
936 | | * an overapproximation of the transitive closure of "bmap" |
937 | | * with path lengths greater than or equal to zero and with |
938 | | * domain and range equal to "dom". |
939 | | */ |
940 | | static __isl_give isl_map *q_closure(__isl_take isl_space *dim, |
941 | | __isl_take isl_set *dom, __isl_keep isl_basic_map *bmap, int *exact) |
942 | 83 | { |
943 | 83 | int project = 1; |
944 | 83 | isl_map *path; |
945 | 83 | isl_map *map; |
946 | 83 | isl_map *app; |
947 | 83 | |
948 | 83 | dom = isl_set_add_dims(dom, isl_dim_set, 1); |
949 | 83 | app = isl_map_from_domain_and_range(dom, isl_set_copy(dom)); |
950 | 83 | map = isl_map_from_basic_map(isl_basic_map_copy(bmap)); |
951 | 83 | path = construct_extended_path(dim, map, &project); |
952 | 83 | app = isl_map_intersect(app, path); |
953 | 83 | |
954 | 83 | if ((*exact = check_exactness(map, isl_map_copy(app), project)) < 0) |
955 | 0 | goto error; |
956 | 83 | |
957 | 83 | return app; |
958 | 0 | error: |
959 | 0 | isl_map_free(app); |
960 | 0 | return NULL; |
961 | 83 | } |
962 | | |
963 | | /* Check whether qc has any elements of length at least one |
964 | | * with domain and/or range outside of dom and ran. |
965 | | */ |
966 | | static int has_spurious_elements(__isl_keep isl_map *qc, |
967 | | __isl_keep isl_set *dom, __isl_keep isl_set *ran) |
968 | 83 | { |
969 | 83 | isl_set *s; |
970 | 83 | int subset; |
971 | 83 | unsigned d; |
972 | 83 | |
973 | 83 | if (!qc || !dom || !ran) |
974 | 0 | return -1; |
975 | 83 | |
976 | 83 | d = isl_map_dim(qc, isl_dim_in); |
977 | 83 | |
978 | 83 | qc = isl_map_copy(qc); |
979 | 83 | qc = set_path_length(qc, 0, 1); |
980 | 83 | qc = isl_map_project_out(qc, isl_dim_in, d - 1, 1); |
981 | 83 | qc = isl_map_project_out(qc, isl_dim_out, d - 1, 1); |
982 | 83 | |
983 | 83 | s = isl_map_domain(isl_map_copy(qc)); |
984 | 83 | subset = isl_set_is_subset(s, dom); |
985 | 83 | isl_set_free(s); |
986 | 83 | if (subset < 0) |
987 | 0 | goto error; |
988 | 83 | if (!subset) { |
989 | 76 | isl_map_free(qc); |
990 | 76 | return 1; |
991 | 76 | } |
992 | 7 | |
993 | 7 | s = isl_map_range(qc); |
994 | 7 | subset = isl_set_is_subset(s, ran); |
995 | 7 | isl_set_free(s); |
996 | 7 | |
997 | 7 | return subset < 0 ? -10 : !subset; |
998 | 0 | error: |
999 | 0 | isl_map_free(qc); |
1000 | 0 | return -1; |
1001 | 7 | } |
1002 | | |
1003 | 82 | #define LEFT 2 |
1004 | 82 | #define RIGHT 1 |
1005 | | |
1006 | | /* For each basic map in "map", except i, check whether it combines |
1007 | | * with the transitive closure that is reflexive on C combines |
1008 | | * to the left and to the right. |
1009 | | * |
1010 | | * In particular, if |
1011 | | * |
1012 | | * dom map_j \subseteq C |
1013 | | * |
1014 | | * then right[j] is set to 1. Otherwise, if |
1015 | | * |
1016 | | * ran map_i \cap dom map_j = \emptyset |
1017 | | * |
1018 | | * then right[j] is set to 0. Otherwise, composing to the right |
1019 | | * is impossible. |
1020 | | * |
1021 | | * Similar, for composing to the left, we have if |
1022 | | * |
1023 | | * ran map_j \subseteq C |
1024 | | * |
1025 | | * then left[j] is set to 1. Otherwise, if |
1026 | | * |
1027 | | * dom map_i \cap ran map_j = \emptyset |
1028 | | * |
1029 | | * then left[j] is set to 0. Otherwise, composing to the left |
1030 | | * is impossible. |
1031 | | * |
1032 | | * The return value is or'd with LEFT if composing to the left |
1033 | | * is possible and with RIGHT if composing to the right is possible. |
1034 | | */ |
1035 | | static int composability(__isl_keep isl_set *C, int i, |
1036 | | isl_set **dom, isl_set **ran, int *left, int *right, |
1037 | | __isl_keep isl_map *map) |
1038 | 40 | { |
1039 | 40 | int j; |
1040 | 40 | int ok; |
1041 | 40 | |
1042 | 40 | ok = LEFT | RIGHT; |
1043 | 120 | for (j = 0; j < map->n && ok80 ; ++j80 ) { |
1044 | 80 | isl_bool overlaps, subset; |
1045 | 80 | if (j == i) |
1046 | 40 | continue; |
1047 | 40 | |
1048 | 40 | if (ok & RIGHT) { |
1049 | 40 | if (!dom[j]) |
1050 | 0 | dom[j] = isl_set_from_basic_set( |
1051 | 0 | isl_basic_map_domain( |
1052 | 0 | isl_basic_map_copy(map->p[j]))); |
1053 | 40 | if (!dom[j]) |
1054 | 0 | return -1; |
1055 | 40 | overlaps = isl_set_overlaps(ran[i], dom[j]); |
1056 | 40 | if (overlaps < 0) |
1057 | 0 | return -1; |
1058 | 40 | if (!overlaps) |
1059 | 0 | right[j] = 0; |
1060 | 40 | else { |
1061 | 40 | subset = isl_set_is_subset(dom[j], C); |
1062 | 40 | if (subset < 0) |
1063 | 0 | return -1; |
1064 | 40 | if (subset) |
1065 | 40 | right[j] = 1; |
1066 | 0 | else |
1067 | 0 | ok &= ~RIGHT; |
1068 | 40 | } |
1069 | 40 | } |
1070 | 40 | |
1071 | 40 | if (ok & LEFT) { |
1072 | 40 | if (!ran[j]) |
1073 | 0 | ran[j] = isl_set_from_basic_set( |
1074 | 0 | isl_basic_map_range( |
1075 | 0 | isl_basic_map_copy(map->p[j]))); |
1076 | 40 | if (!ran[j]) |
1077 | 0 | return -1; |
1078 | 40 | overlaps = isl_set_overlaps(dom[i], ran[j]); |
1079 | 40 | if (overlaps < 0) |
1080 | 0 | return -1; |
1081 | 40 | if (!overlaps) |
1082 | 0 | left[j] = 0; |
1083 | 40 | else { |
1084 | 40 | subset = isl_set_is_subset(ran[j], C); |
1085 | 40 | if (subset < 0) |
1086 | 0 | return -1; |
1087 | 40 | if (subset) |
1088 | 40 | left[j] = 1; |
1089 | 0 | else |
1090 | 0 | ok &= ~LEFT; |
1091 | 40 | } |
1092 | 40 | } |
1093 | 40 | } |
1094 | 40 | |
1095 | 40 | return ok; |
1096 | 40 | } |
1097 | | |
1098 | | static __isl_give isl_map *anonymize(__isl_take isl_map *map) |
1099 | 69 | { |
1100 | 69 | map = isl_map_reset(map, isl_dim_in); |
1101 | 69 | map = isl_map_reset(map, isl_dim_out); |
1102 | 69 | return map; |
1103 | 69 | } |
1104 | | |
1105 | | /* Return a map that is a union of the basic maps in "map", except i, |
1106 | | * composed to left and right with qc based on the entries of "left" |
1107 | | * and "right". |
1108 | | */ |
1109 | | static __isl_give isl_map *compose(__isl_keep isl_map *map, int i, |
1110 | | __isl_take isl_map *qc, int *left, int *right) |
1111 | 7 | { |
1112 | 7 | int j; |
1113 | 7 | isl_map *comp; |
1114 | 7 | |
1115 | 7 | comp = isl_map_empty(isl_map_get_space(map)); |
1116 | 21 | for (j = 0; j < map->n; ++j14 ) { |
1117 | 14 | isl_map *map_j; |
1118 | 14 | |
1119 | 14 | if (j == i) |
1120 | 7 | continue; |
1121 | 7 | |
1122 | 7 | map_j = isl_map_from_basic_map(isl_basic_map_copy(map->p[j])); |
1123 | 7 | map_j = anonymize(map_j); |
1124 | 7 | if (left && left[j]) |
1125 | 7 | map_j = isl_map_apply_range(map_j, isl_map_copy(qc)); |
1126 | 7 | if (right && right[j]) |
1127 | 7 | map_j = isl_map_apply_range(isl_map_copy(qc), map_j); |
1128 | 7 | comp = isl_map_union(comp, map_j); |
1129 | 7 | } |
1130 | 7 | |
1131 | 7 | comp = isl_map_compute_divs(comp); |
1132 | 7 | comp = isl_map_coalesce(comp); |
1133 | 7 | |
1134 | 7 | isl_map_free(qc); |
1135 | 7 | |
1136 | 7 | return comp; |
1137 | 7 | } |
1138 | | |
1139 | | /* Compute the transitive closure of "map" incrementally by |
1140 | | * computing |
1141 | | * |
1142 | | * map_i^+ \cup qc^+ |
1143 | | * |
1144 | | * or |
1145 | | * |
1146 | | * map_i^+ \cup ((id \cup map_i^) \circ qc^+) |
1147 | | * |
1148 | | * or |
1149 | | * |
1150 | | * map_i^+ \cup (qc^+ \circ (id \cup map_i^)) |
1151 | | * |
1152 | | * depending on whether left or right are NULL. |
1153 | | */ |
1154 | | static __isl_give isl_map *compute_incremental( |
1155 | | __isl_take isl_space *dim, __isl_keep isl_map *map, |
1156 | | int i, __isl_take isl_map *qc, int *left, int *right, int *exact) |
1157 | 3 | { |
1158 | 3 | isl_map *map_i; |
1159 | 3 | isl_map *tc; |
1160 | 3 | isl_map *rtc = NULL; |
1161 | 3 | |
1162 | 3 | if (!map) |
1163 | 0 | goto error; |
1164 | 3 | isl_assert(map->ctx, left || right, goto error); |
1165 | 3 | |
1166 | 3 | map_i = isl_map_from_basic_map(isl_basic_map_copy(map->p[i])); |
1167 | 3 | tc = construct_projected_component(isl_space_copy(dim), map_i, |
1168 | 3 | exact, 1); |
1169 | 3 | isl_map_free(map_i); |
1170 | 3 | |
1171 | 3 | if (*exact) |
1172 | 3 | qc = isl_map_transitive_closure(qc, exact); |
1173 | 3 | |
1174 | 3 | if (!*exact) { |
1175 | 0 | isl_space_free(dim); |
1176 | 0 | isl_map_free(tc); |
1177 | 0 | isl_map_free(qc); |
1178 | 0 | return isl_map_universe(isl_map_get_space(map)); |
1179 | 0 | } |
1180 | 3 | |
1181 | 3 | if (!left || !right) |
1182 | 0 | rtc = isl_map_union(isl_map_copy(tc), |
1183 | 0 | isl_map_identity(isl_map_get_space(tc))); |
1184 | 3 | if (!right) |
1185 | 0 | qc = isl_map_apply_range(rtc, qc); |
1186 | 3 | if (!left) |
1187 | 0 | qc = isl_map_apply_range(qc, rtc); |
1188 | 3 | qc = isl_map_union(tc, qc); |
1189 | 3 | |
1190 | 3 | isl_space_free(dim); |
1191 | 3 | |
1192 | 3 | return qc; |
1193 | 0 | error: |
1194 | 0 | isl_space_free(dim); |
1195 | 0 | isl_map_free(qc); |
1196 | 0 | return NULL; |
1197 | 3 | } |
1198 | | |
1199 | | /* Given a map "map", try to find a basic map such that |
1200 | | * map^+ can be computed as |
1201 | | * |
1202 | | * map^+ = map_i^+ \cup |
1203 | | * \bigcup_j ((map_i^+ \cup Id_C)^+ \circ map_j \circ (map_i^+ \cup Id_C))^+ |
1204 | | * |
1205 | | * with C the simple hull of the domain and range of the input map. |
1206 | | * map_i^ \cup Id_C is computed by allowing the path lengths to be zero |
1207 | | * and by intersecting domain and range with C. |
1208 | | * Of course, we need to check that this is actually equal to map_i^ \cup Id_C. |
1209 | | * Also, we only use the incremental computation if all the transitive |
1210 | | * closures are exact and if the number of basic maps in the union, |
1211 | | * after computing the integer divisions, is smaller than the number |
1212 | | * of basic maps in the input map. |
1213 | | */ |
1214 | | static int incemental_on_entire_domain(__isl_keep isl_space *dim, |
1215 | | __isl_keep isl_map *map, |
1216 | | isl_set **dom, isl_set **ran, int *left, int *right, |
1217 | | __isl_give isl_map **res) |
1218 | 23 | { |
1219 | 23 | int i; |
1220 | 23 | isl_set *C; |
1221 | 23 | unsigned d; |
1222 | 23 | |
1223 | 23 | *res = NULL; |
1224 | 23 | |
1225 | 23 | C = isl_set_union(isl_map_domain(isl_map_copy(map)), |
1226 | 23 | isl_map_range(isl_map_copy(map))); |
1227 | 23 | C = isl_set_from_basic_set(isl_set_simple_hull(C)); |
1228 | 23 | if (!C) |
1229 | 0 | return -1; |
1230 | 23 | if (C->n != 1) { |
1231 | 0 | isl_set_free(C); |
1232 | 0 | return 0; |
1233 | 0 | } |
1234 | 23 | |
1235 | 23 | d = isl_map_dim(map, isl_dim_in); |
1236 | 23 | |
1237 | 63 | for (i = 0; i < map->n; ++i40 ) { |
1238 | 43 | isl_map *qc; |
1239 | 43 | int exact_i, spurious; |
1240 | 43 | int j; |
1241 | 43 | dom[i] = isl_set_from_basic_set(isl_basic_map_domain( |
1242 | 43 | isl_basic_map_copy(map->p[i]))); |
1243 | 43 | ran[i] = isl_set_from_basic_set(isl_basic_map_range( |
1244 | 43 | isl_basic_map_copy(map->p[i]))); |
1245 | 43 | qc = q_closure(isl_space_copy(dim), isl_set_copy(C), |
1246 | 43 | map->p[i], &exact_i); |
1247 | 43 | if (!qc) |
1248 | 0 | goto error; |
1249 | 43 | if (!exact_i) { |
1250 | 0 | isl_map_free(qc); |
1251 | 0 | continue; |
1252 | 0 | } |
1253 | 43 | spurious = has_spurious_elements(qc, dom[i], ran[i]); |
1254 | 43 | if (spurious) { |
1255 | 38 | isl_map_free(qc); |
1256 | 38 | if (spurious < 0) |
1257 | 0 | goto error; |
1258 | 38 | continue; |
1259 | 38 | } |
1260 | 5 | qc = isl_map_project_out(qc, isl_dim_in, d, 1); |
1261 | 5 | qc = isl_map_project_out(qc, isl_dim_out, d, 1); |
1262 | 5 | qc = isl_map_compute_divs(qc); |
1263 | 15 | for (j = 0; j < map->n; ++j10 ) |
1264 | 10 | left[j] = right[j] = 1; |
1265 | 5 | qc = compose(map, i, qc, left, right); |
1266 | 5 | if (!qc) |
1267 | 0 | goto error; |
1268 | 5 | if (qc->n >= map->n) { |
1269 | 2 | isl_map_free(qc); |
1270 | 2 | continue; |
1271 | 2 | } |
1272 | 3 | *res = compute_incremental(isl_space_copy(dim), map, i, qc, |
1273 | 3 | left, right, &exact_i); |
1274 | 3 | if (!*res) |
1275 | 0 | goto error; |
1276 | 3 | if (exact_i) |
1277 | 3 | break; |
1278 | 0 | isl_map_free(*res); |
1279 | 0 | *res = NULL; |
1280 | 0 | } |
1281 | 23 | |
1282 | 23 | isl_set_free(C); |
1283 | 23 | |
1284 | 23 | return *res != NULL; |
1285 | 23 | error: |
1286 | 0 | isl_set_free(C); |
1287 | 0 | return -1; |
1288 | 23 | } |
1289 | | |
1290 | | /* Try and compute the transitive closure of "map" as |
1291 | | * |
1292 | | * map^+ = map_i^+ \cup |
1293 | | * \bigcup_j ((map_i^+ \cup Id_C)^+ \circ map_j \circ (map_i^+ \cup Id_C))^+ |
1294 | | * |
1295 | | * with C either the simple hull of the domain and range of the entire |
1296 | | * map or the simple hull of domain and range of map_i. |
1297 | | */ |
1298 | | static __isl_give isl_map *incremental_closure(__isl_take isl_space *dim, |
1299 | | __isl_keep isl_map *map, int *exact, int project) |
1300 | 113 | { |
1301 | 113 | int i; |
1302 | 113 | isl_set **dom = NULL; |
1303 | 113 | isl_set **ran = NULL; |
1304 | 113 | int *left = NULL; |
1305 | 113 | int *right = NULL; |
1306 | 113 | isl_set *C; |
1307 | 113 | unsigned d; |
1308 | 113 | isl_map *res = NULL; |
1309 | 113 | |
1310 | 113 | if (!project) |
1311 | 3 | return construct_projected_component(dim, map, exact, project); |
1312 | 110 | |
1313 | 110 | if (!map) |
1314 | 0 | goto error; |
1315 | 110 | if (map->n <= 1) |
1316 | 87 | return construct_projected_component(dim, map, exact, project); |
1317 | 23 | |
1318 | 23 | d = isl_map_dim(map, isl_dim_in); |
1319 | 23 | |
1320 | 23 | dom = isl_calloc_array(map->ctx, isl_set *, map->n); |
1321 | 23 | ran = isl_calloc_array(map->ctx, isl_set *, map->n); |
1322 | 23 | left = isl_calloc_array(map->ctx, int, map->n); |
1323 | 23 | right = isl_calloc_array(map->ctx, int, map->n); |
1324 | 23 | if (!ran || !dom || !left || !right) |
1325 | 0 | goto error; |
1326 | 23 | |
1327 | 23 | if (incemental_on_entire_domain(dim, map, dom, ran, left, right, &res) < 0) |
1328 | 0 | goto error; |
1329 | 23 | |
1330 | 63 | for (i = 0; 23 !res && i < map->n60 ; ++i40 ) { |
1331 | 40 | isl_map *qc; |
1332 | 40 | int exact_i, spurious, comp; |
1333 | 40 | if (!dom[i]) |
1334 | 0 | dom[i] = isl_set_from_basic_set( |
1335 | 0 | isl_basic_map_domain( |
1336 | 0 | isl_basic_map_copy(map->p[i]))); |
1337 | 40 | if (!dom[i]) |
1338 | 0 | goto error; |
1339 | 40 | if (!ran[i]) |
1340 | 0 | ran[i] = isl_set_from_basic_set( |
1341 | 0 | isl_basic_map_range( |
1342 | 0 | isl_basic_map_copy(map->p[i]))); |
1343 | 40 | if (!ran[i]) |
1344 | 0 | goto error; |
1345 | 40 | C = isl_set_union(isl_set_copy(dom[i]), |
1346 | 40 | isl_set_copy(ran[i])); |
1347 | 40 | C = isl_set_from_basic_set(isl_set_simple_hull(C)); |
1348 | 40 | if (!C) |
1349 | 0 | goto error; |
1350 | 40 | if (C->n != 1) { |
1351 | 0 | isl_set_free(C); |
1352 | 0 | continue; |
1353 | 0 | } |
1354 | 40 | comp = composability(C, i, dom, ran, left, right, map); |
1355 | 40 | if (!comp || comp < 0) { |
1356 | 0 | isl_set_free(C); |
1357 | 0 | if (comp < 0) |
1358 | 0 | goto error; |
1359 | 0 | continue; |
1360 | 0 | } |
1361 | 40 | qc = q_closure(isl_space_copy(dim), C, map->p[i], &exact_i); |
1362 | 40 | if (!qc) |
1363 | 0 | goto error; |
1364 | 40 | if (!exact_i) { |
1365 | 0 | isl_map_free(qc); |
1366 | 0 | continue; |
1367 | 0 | } |
1368 | 40 | spurious = has_spurious_elements(qc, dom[i], ran[i]); |
1369 | 40 | if (spurious) { |
1370 | 38 | isl_map_free(qc); |
1371 | 38 | if (spurious < 0) |
1372 | 0 | goto error; |
1373 | 38 | continue; |
1374 | 38 | } |
1375 | 2 | qc = isl_map_project_out(qc, isl_dim_in, d, 1); |
1376 | 2 | qc = isl_map_project_out(qc, isl_dim_out, d, 1); |
1377 | 2 | qc = isl_map_compute_divs(qc); |
1378 | 2 | qc = compose(map, i, qc, (comp & LEFT) ? left : NULL, |
1379 | 2 | (comp & RIGHT) ? right : NULL); |
1380 | 2 | if (!qc) |
1381 | 0 | goto error; |
1382 | 2 | if (qc->n >= map->n) { |
1383 | 2 | isl_map_free(qc); |
1384 | 2 | continue; |
1385 | 2 | } |
1386 | 0 | res = compute_incremental(isl_space_copy(dim), map, i, qc, |
1387 | 0 | (comp & LEFT) ? left : NULL, |
1388 | 0 | (comp & RIGHT) ? right : NULL, &exact_i); |
1389 | 0 | if (!res) |
1390 | 0 | goto error; |
1391 | 0 | if (exact_i) |
1392 | 0 | break; |
1393 | 0 | isl_map_free(res); |
1394 | 0 | res = NULL; |
1395 | 0 | } |
1396 | 23 | |
1397 | 69 | for (i = 0; 23 i < map->n; ++i46 ) { |
1398 | 46 | isl_set_free(dom[i]); |
1399 | 46 | isl_set_free(ran[i]); |
1400 | 46 | } |
1401 | 23 | free(dom); |
1402 | 23 | free(ran); |
1403 | 23 | free(left); |
1404 | 23 | free(right); |
1405 | 23 | |
1406 | 23 | if (res) { |
1407 | 3 | isl_space_free(dim); |
1408 | 3 | return res; |
1409 | 3 | } |
1410 | 20 | |
1411 | 20 | return construct_projected_component(dim, map, exact, project); |
1412 | 0 | error: |
1413 | 0 | if (dom) |
1414 | 0 | for (i = 0; i < map->n; ++i) |
1415 | 0 | isl_set_free(dom[i]); |
1416 | 0 | free(dom); |
1417 | 0 | if (ran) |
1418 | 0 | for (i = 0; i < map->n; ++i) |
1419 | 0 | isl_set_free(ran[i]); |
1420 | 0 | free(ran); |
1421 | 0 | free(left); |
1422 | 0 | free(right); |
1423 | 0 | isl_space_free(dim); |
1424 | 0 | return NULL; |
1425 | 20 | } |
1426 | | |
1427 | | /* Given an array of sets "set", add "dom" at position "pos" |
1428 | | * and search for elements at earlier positions that overlap with "dom". |
1429 | | * If any can be found, then merge all of them, together with "dom", into |
1430 | | * a single set and assign the union to the first in the array, |
1431 | | * which becomes the new group leader for all groups involved in the merge. |
1432 | | * During the search, we only consider group leaders, i.e., those with |
1433 | | * group[i] = i, as the other sets have already been combined |
1434 | | * with one of the group leaders. |
1435 | | */ |
1436 | | static int merge(isl_set **set, int *group, __isl_take isl_set *dom, int pos) |
1437 | 334 | { |
1438 | 334 | int i; |
1439 | 334 | |
1440 | 334 | group[pos] = pos; |
1441 | 334 | set[pos] = isl_set_copy(dom); |
1442 | 334 | |
1443 | 685 | for (i = pos - 1; i >= 0; --i351 ) { |
1444 | 351 | isl_bool o; |
1445 | 351 | |
1446 | 351 | if (group[i] != i) |
1447 | 134 | continue; |
1448 | 217 | |
1449 | 217 | o = isl_set_overlaps(set[i], dom); |
1450 | 217 | if (o < 0) |
1451 | 0 | goto error; |
1452 | 217 | if (!o) |
1453 | 7 | continue; |
1454 | 210 | |
1455 | 210 | set[i] = isl_set_union(set[i], set[group[pos]]); |
1456 | 210 | set[group[pos]] = NULL; |
1457 | 210 | if (!set[i]) |
1458 | 0 | goto error; |
1459 | 210 | group[group[pos]] = i; |
1460 | 210 | group[pos] = i; |
1461 | 210 | } |
1462 | 334 | |
1463 | 334 | isl_set_free(dom); |
1464 | 334 | return 0; |
1465 | 0 | error: |
1466 | 0 | isl_set_free(dom); |
1467 | 0 | return -1; |
1468 | 334 | } |
1469 | | |
1470 | | /* Replace each entry in the n by n grid of maps by the cross product |
1471 | | * with the relation { [i] -> [i + 1] }. |
1472 | | */ |
1473 | | static int add_length(__isl_keep isl_map *map, isl_map ***grid, int n) |
1474 | 1 | { |
1475 | 1 | int i, j, k; |
1476 | 1 | isl_space *dim; |
1477 | 1 | isl_basic_map *bstep; |
1478 | 1 | isl_map *step; |
1479 | 1 | unsigned nparam; |
1480 | 1 | |
1481 | 1 | if (!map) |
1482 | 0 | return -1; |
1483 | 1 | |
1484 | 1 | dim = isl_map_get_space(map); |
1485 | 1 | nparam = isl_space_dim(dim, isl_dim_param); |
1486 | 1 | dim = isl_space_drop_dims(dim, isl_dim_in, 0, isl_space_dim(dim, isl_dim_in)); |
1487 | 1 | dim = isl_space_drop_dims(dim, isl_dim_out, 0, isl_space_dim(dim, isl_dim_out)); |
1488 | 1 | dim = isl_space_add_dims(dim, isl_dim_in, 1); |
1489 | 1 | dim = isl_space_add_dims(dim, isl_dim_out, 1); |
1490 | 1 | bstep = isl_basic_map_alloc_space(dim, 0, 1, 0); |
1491 | 1 | k = isl_basic_map_alloc_equality(bstep); |
1492 | 1 | if (k < 0) { |
1493 | 0 | isl_basic_map_free(bstep); |
1494 | 0 | return -1; |
1495 | 0 | } |
1496 | 1 | isl_seq_clr(bstep->eq[k], 1 + isl_basic_map_total_dim(bstep)); |
1497 | 1 | isl_int_set_si(bstep->eq[k][0], 1); |
1498 | 1 | isl_int_set_si(bstep->eq[k][1 + nparam], 1); |
1499 | 1 | isl_int_set_si(bstep->eq[k][1 + nparam + 1], -1); |
1500 | 1 | bstep = isl_basic_map_finalize(bstep); |
1501 | 1 | step = isl_map_from_basic_map(bstep); |
1502 | 1 | |
1503 | 3 | for (i = 0; i < n; ++i2 ) |
1504 | 6 | for (j = 0; 2 j < n; ++j4 ) |
1505 | 4 | grid[i][j] = isl_map_product(grid[i][j], |
1506 | 4 | isl_map_copy(step)); |
1507 | 1 | |
1508 | 1 | isl_map_free(step); |
1509 | 1 | |
1510 | 1 | return 0; |
1511 | 1 | } |
1512 | | |
1513 | | /* The core of the Floyd-Warshall algorithm. |
1514 | | * Updates the given n x x matrix of relations in place. |
1515 | | * |
1516 | | * The algorithm iterates over all vertices. In each step, the whole |
1517 | | * matrix is updated to include all paths that go to the current vertex, |
1518 | | * possibly stay there a while (including passing through earlier vertices) |
1519 | | * and then come back. At the start of each iteration, the diagonal |
1520 | | * element corresponding to the current vertex is replaced by its |
1521 | | * transitive closure to account for all indirect paths that stay |
1522 | | * in the current vertex. |
1523 | | */ |
1524 | | static void floyd_warshall_iterate(isl_map ***grid, int n, int *exact) |
1525 | 96 | { |
1526 | 96 | int r, p, q; |
1527 | 96 | |
1528 | 195 | for (r = 0; r < n; ++r99 ) { |
1529 | 99 | int r_exact; |
1530 | 99 | grid[r][r] = isl_map_transitive_closure(grid[r][r], |
1531 | 99 | (exact && *exact4 ) ? &r_exact4 : NULL); |
1532 | 99 | if (exact && *exact4 && !r_exact4 ) |
1533 | 0 | *exact = 0; |
1534 | 99 | |
1535 | 204 | for (p = 0; p < n; ++p105 ) |
1536 | 222 | for (q = 0; 105 q < n; ++q117 ) { |
1537 | 117 | isl_map *loop; |
1538 | 117 | if (p == r && q == r105 ) |
1539 | 99 | continue; |
1540 | 18 | loop = isl_map_apply_range( |
1541 | 18 | isl_map_copy(grid[p][r]), |
1542 | 18 | isl_map_copy(grid[r][q])); |
1543 | 18 | grid[p][q] = isl_map_union(grid[p][q], loop); |
1544 | 18 | loop = isl_map_apply_range( |
1545 | 18 | isl_map_copy(grid[p][r]), |
1546 | 18 | isl_map_apply_range( |
1547 | 18 | isl_map_copy(grid[r][r]), |
1548 | 18 | isl_map_copy(grid[r][q]))); |
1549 | 18 | grid[p][q] = isl_map_union(grid[p][q], loop); |
1550 | 18 | grid[p][q] = isl_map_coalesce(grid[p][q]); |
1551 | 18 | } |
1552 | 99 | } |
1553 | 96 | } |
1554 | | |
1555 | | /* Given a partition of the domains and ranges of the basic maps in "map", |
1556 | | * apply the Floyd-Warshall algorithm with the elements in the partition |
1557 | | * as vertices. |
1558 | | * |
1559 | | * In particular, there are "n" elements in the partition and "group" is |
1560 | | * an array of length 2 * map->n with entries in [0,n-1]. |
1561 | | * |
1562 | | * We first construct a matrix of relations based on the partition information, |
1563 | | * apply Floyd-Warshall on this matrix of relations and then take the |
1564 | | * union of all entries in the matrix as the final result. |
1565 | | * |
1566 | | * If we are actually computing the power instead of the transitive closure, |
1567 | | * i.e., when "project" is not set, then the result should have the |
1568 | | * path lengths encoded as the difference between an extra pair of |
1569 | | * coordinates. We therefore apply the nested transitive closures |
1570 | | * to relations that include these lengths. In particular, we replace |
1571 | | * the input relation by the cross product with the unit length relation |
1572 | | * { [i] -> [i + 1] }. |
1573 | | */ |
1574 | | static __isl_give isl_map *floyd_warshall_with_groups(__isl_take isl_space *dim, |
1575 | | __isl_keep isl_map *map, int *exact, int project, int *group, int n) |
1576 | 27 | { |
1577 | 27 | int i, j, k; |
1578 | 27 | isl_map ***grid = NULL; |
1579 | 27 | isl_map *app; |
1580 | 27 | |
1581 | 27 | if (!map) |
1582 | 0 | goto error; |
1583 | 27 | |
1584 | 27 | if (n == 1) { |
1585 | 25 | free(group); |
1586 | 25 | return incremental_closure(dim, map, exact, project); |
1587 | 25 | } |
1588 | 2 | |
1589 | 2 | grid = isl_calloc_array(map->ctx, isl_map **, n); |
1590 | 2 | if (!grid) |
1591 | 0 | goto error; |
1592 | 6 | for (i = 0; 2 i < n; ++i4 ) { |
1593 | 4 | grid[i] = isl_calloc_array(map->ctx, isl_map *, n); |
1594 | 4 | if (!grid[i]) |
1595 | 0 | goto error; |
1596 | 12 | for (j = 0; 4 j < n; ++j8 ) |
1597 | 8 | grid[i][j] = isl_map_empty(isl_map_get_space(map)); |
1598 | 4 | } |
1599 | 2 | |
1600 | 6 | for (k = 0; 2 k < map->n; ++k4 ) { |
1601 | 4 | i = group[2 * k]; |
1602 | 4 | j = group[2 * k + 1]; |
1603 | 4 | grid[i][j] = isl_map_union(grid[i][j], |
1604 | 4 | isl_map_from_basic_map( |
1605 | 4 | isl_basic_map_copy(map->p[k]))); |
1606 | 4 | } |
1607 | 2 | |
1608 | 2 | if (!project && add_length(map, grid, n) < 01 ) |
1609 | 0 | goto error; |
1610 | 2 | |
1611 | 2 | floyd_warshall_iterate(grid, n, exact); |
1612 | 2 | |
1613 | 2 | app = isl_map_empty(isl_map_get_space(grid[0][0])); |
1614 | 2 | |
1615 | 6 | for (i = 0; i < n; ++i4 ) { |
1616 | 12 | for (j = 0; j < n; ++j8 ) |
1617 | 8 | app = isl_map_union(app, grid[i][j]); |
1618 | 4 | free(grid[i]); |
1619 | 4 | } |
1620 | 2 | free(grid); |
1621 | 2 | |
1622 | 2 | free(group); |
1623 | 2 | isl_space_free(dim); |
1624 | 2 | |
1625 | 2 | return app; |
1626 | 0 | error: |
1627 | 0 | if (grid) |
1628 | 0 | for (i = 0; i < n; ++i) { |
1629 | 0 | if (!grid[i]) |
1630 | 0 | continue; |
1631 | 0 | for (j = 0; j < n; ++j) |
1632 | 0 | isl_map_free(grid[i][j]); |
1633 | 0 | free(grid[i]); |
1634 | 0 | } |
1635 | 0 | free(grid); |
1636 | 0 | free(group); |
1637 | 0 | isl_space_free(dim); |
1638 | 0 | return NULL; |
1639 | 2 | } |
1640 | | |
1641 | | /* Partition the domains and ranges of the n basic relations in list |
1642 | | * into disjoint cells. |
1643 | | * |
1644 | | * To find the partition, we simply consider all of the domains |
1645 | | * and ranges in turn and combine those that overlap. |
1646 | | * "set" contains the partition elements and "group" indicates |
1647 | | * to which partition element a given domain or range belongs. |
1648 | | * The domain of basic map i corresponds to element 2 * i in these arrays, |
1649 | | * while the domain corresponds to element 2 * i + 1. |
1650 | | * During the construction group[k] is either equal to k, |
1651 | | * in which case set[k] contains the union of all the domains and |
1652 | | * ranges in the corresponding group, or is equal to some l < k, |
1653 | | * with l another domain or range in the same group. |
1654 | | */ |
1655 | | static int *setup_groups(isl_ctx *ctx, __isl_keep isl_basic_map **list, int n, |
1656 | | isl_set ***set, int *n_group) |
1657 | 121 | { |
1658 | 121 | int i; |
1659 | 121 | int *group = NULL; |
1660 | 121 | int g; |
1661 | 121 | |
1662 | 121 | *set = isl_calloc_array(ctx, isl_set *, 2 * n); |
1663 | 121 | group = isl_alloc_array(ctx, int, 2 * n); |
1664 | 121 | |
1665 | 121 | if (!*set || !group) |
1666 | 0 | goto error; |
1667 | 121 | |
1668 | 288 | for (i = 0; 121 i < n; ++i167 ) { |
1669 | 167 | isl_set *dom; |
1670 | 167 | dom = isl_set_from_basic_set(isl_basic_map_domain( |
1671 | 167 | isl_basic_map_copy(list[i]))); |
1672 | 167 | if (merge(*set, group, dom, 2 * i) < 0) |
1673 | 0 | goto error; |
1674 | 167 | dom = isl_set_from_basic_set(isl_basic_map_range( |
1675 | 167 | isl_basic_map_copy(list[i]))); |
1676 | 167 | if (merge(*set, group, dom, 2 * i + 1) < 0) |
1677 | 0 | goto error; |
1678 | 167 | } |
1679 | 121 | |
1680 | 121 | g = 0; |
1681 | 455 | for (i = 0; i < 2 * n; ++i334 ) |
1682 | 334 | if (group[i] == i) { |
1683 | 124 | if (g != i) { |
1684 | 0 | (*set)[g] = (*set)[i]; |
1685 | 0 | (*set)[i] = NULL; |
1686 | 0 | } |
1687 | 124 | group[i] = g++; |
1688 | 124 | } else |
1689 | 210 | group[i] = group[group[i]]; |
1690 | 121 | |
1691 | 121 | *n_group = g; |
1692 | 121 | |
1693 | 121 | return group; |
1694 | 0 | error: |
1695 | 0 | if (*set) { |
1696 | 0 | for (i = 0; i < 2 * n; ++i) |
1697 | 0 | isl_set_free((*set)[i]); |
1698 | 0 | free(*set); |
1699 | 0 | *set = NULL; |
1700 | 0 | } |
1701 | 0 | free(group); |
1702 | 0 | return NULL; |
1703 | 121 | } |
1704 | | |
1705 | | /* Check if the domains and ranges of the basic maps in "map" can |
1706 | | * be partitioned, and if so, apply Floyd-Warshall on the elements |
1707 | | * of the partition. Note that we also apply this algorithm |
1708 | | * if we want to compute the power, i.e., when "project" is not set. |
1709 | | * However, the results are unlikely to be exact since the recursive |
1710 | | * calls inside the Floyd-Warshall algorithm typically result in |
1711 | | * non-linear path lengths quite quickly. |
1712 | | */ |
1713 | | static __isl_give isl_map *floyd_warshall(__isl_take isl_space *dim, |
1714 | | __isl_keep isl_map *map, int *exact, int project) |
1715 | 115 | { |
1716 | 115 | int i; |
1717 | 115 | isl_set **set = NULL; |
1718 | 115 | int *group = NULL; |
1719 | 115 | int n; |
1720 | 115 | |
1721 | 115 | if (!map) |
1722 | 0 | goto error; |
1723 | 115 | if (map->n <= 1) |
1724 | 88 | return incremental_closure(dim, map, exact, project); |
1725 | 27 | |
1726 | 27 | group = setup_groups(map->ctx, map->p, map->n, &set, &n); |
1727 | 27 | if (!group) |
1728 | 0 | goto error; |
1729 | 27 | |
1730 | 135 | for (i = 0; 27 i < 2 * map->n; ++i108 ) |
1731 | 108 | isl_set_free(set[i]); |
1732 | 27 | |
1733 | 27 | free(set); |
1734 | 27 | |
1735 | 27 | return floyd_warshall_with_groups(dim, map, exact, project, group, n); |
1736 | 0 | error: |
1737 | 0 | isl_space_free(dim); |
1738 | 0 | return NULL; |
1739 | 27 | } |
1740 | | |
1741 | | /* Structure for representing the nodes of the graph of which |
1742 | | * strongly connected components are being computed. |
1743 | | * |
1744 | | * list contains the actual nodes |
1745 | | * check_closed is set if we may have used the fact that |
1746 | | * a pair of basic maps can be interchanged |
1747 | | */ |
1748 | | struct isl_tc_follows_data { |
1749 | | isl_basic_map **list; |
1750 | | int check_closed; |
1751 | | }; |
1752 | | |
1753 | | /* Check whether in the computation of the transitive closure |
1754 | | * "list[i]" (R_1) should follow (or be part of the same component as) |
1755 | | * "list[j]" (R_2). |
1756 | | * |
1757 | | * That is check whether |
1758 | | * |
1759 | | * R_1 \circ R_2 |
1760 | | * |
1761 | | * is a subset of |
1762 | | * |
1763 | | * R_2 \circ R_1 |
1764 | | * |
1765 | | * If so, then there is no reason for R_1 to immediately follow R_2 |
1766 | | * in any path. |
1767 | | * |
1768 | | * *check_closed is set if the subset relation holds while |
1769 | | * R_1 \circ R_2 is not empty. |
1770 | | */ |
1771 | | static isl_bool basic_map_follows(int i, int j, void *user) |
1772 | 479 | { |
1773 | 479 | struct isl_tc_follows_data *data = user; |
1774 | 479 | struct isl_map *map12 = NULL; |
1775 | 479 | struct isl_map *map21 = NULL; |
1776 | 479 | isl_bool subset; |
1777 | 479 | |
1778 | 479 | if (!isl_space_tuple_is_equal(data->list[i]->dim, isl_dim_in, |
1779 | 479 | data->list[j]->dim, isl_dim_out)) |
1780 | 377 | return isl_bool_false; |
1781 | 102 | |
1782 | 102 | map21 = isl_map_from_basic_map( |
1783 | 102 | isl_basic_map_apply_range( |
1784 | 102 | isl_basic_map_copy(data->list[j]), |
1785 | 102 | isl_basic_map_copy(data->list[i]))); |
1786 | 102 | subset = isl_map_is_empty(map21); |
1787 | 102 | if (subset < 0) |
1788 | 0 | goto error; |
1789 | 102 | if (subset) { |
1790 | 3 | isl_map_free(map21); |
1791 | 3 | return isl_bool_false; |
1792 | 3 | } |
1793 | 99 | |
1794 | 99 | if (!isl_space_tuple_is_equal(data->list[i]->dim, isl_dim_in, |
1795 | 99 | data->list[i]->dim, isl_dim_out) || |
1796 | 99 | !isl_space_tuple_is_equal(data->list[j]->dim, isl_dim_in, |
1797 | 99 | data->list[j]->dim, isl_dim_out)) { |
1798 | 0 | isl_map_free(map21); |
1799 | 0 | return isl_bool_true; |
1800 | 0 | } |
1801 | 99 | |
1802 | 99 | map12 = isl_map_from_basic_map( |
1803 | 99 | isl_basic_map_apply_range( |
1804 | 99 | isl_basic_map_copy(data->list[i]), |
1805 | 99 | isl_basic_map_copy(data->list[j]))); |
1806 | 99 | |
1807 | 99 | subset = isl_map_is_subset(map21, map12); |
1808 | 99 | |
1809 | 99 | isl_map_free(map12); |
1810 | 99 | isl_map_free(map21); |
1811 | 99 | |
1812 | 99 | if (subset) |
1813 | 6 | data->check_closed = 1; |
1814 | 99 | |
1815 | 99 | return subset < 0 ? isl_bool_error0 : !subset; |
1816 | 0 | error: |
1817 | 0 | isl_map_free(map21); |
1818 | 0 | return isl_bool_error; |
1819 | 99 | } |
1820 | | |
1821 | | /* Given a union of basic maps R = \cup_i R_i \subseteq D \times D |
1822 | | * and a dimension specification (Z^{n+1} -> Z^{n+1}), |
1823 | | * construct a map that is an overapproximation of the map |
1824 | | * that takes an element from the dom R \times Z to an |
1825 | | * element from ran R \times Z, such that the first n coordinates of the |
1826 | | * difference between them is a sum of differences between images |
1827 | | * and pre-images in one of the R_i and such that the last coordinate |
1828 | | * is equal to the number of steps taken. |
1829 | | * If "project" is set, then these final coordinates are not included, |
1830 | | * i.e., a relation of type Z^n -> Z^n is returned. |
1831 | | * That is, let |
1832 | | * |
1833 | | * \Delta_i = { y - x | (x, y) in R_i } |
1834 | | * |
1835 | | * then the constructed map is an overapproximation of |
1836 | | * |
1837 | | * { (x) -> (x + d) | \exists k_i >= 0, \delta_i \in \Delta_i : |
1838 | | * d = (\sum_i k_i \delta_i, \sum_i k_i) and |
1839 | | * x in dom R and x + d in ran R } |
1840 | | * |
1841 | | * or |
1842 | | * |
1843 | | * { (x) -> (x + d) | \exists k_i >= 0, \delta_i \in \Delta_i : |
1844 | | * d = (\sum_i k_i \delta_i) and |
1845 | | * x in dom R and x + d in ran R } |
1846 | | * |
1847 | | * if "project" is set. |
1848 | | * |
1849 | | * We first split the map into strongly connected components, perform |
1850 | | * the above on each component and then join the results in the correct |
1851 | | * order, at each join also taking in the union of both arguments |
1852 | | * to allow for paths that do not go through one of the two arguments. |
1853 | | */ |
1854 | | static __isl_give isl_map *construct_power_components(__isl_take isl_space *dim, |
1855 | | __isl_keep isl_map *map, int *exact, int project) |
1856 | 111 | { |
1857 | 111 | int i, n, c; |
1858 | 111 | struct isl_map *path = NULL; |
1859 | 111 | struct isl_tc_follows_data data; |
1860 | 111 | struct isl_tarjan_graph *g = NULL; |
1861 | 111 | int *orig_exact; |
1862 | 111 | int local_exact; |
1863 | 111 | |
1864 | 111 | if (!map) |
1865 | 0 | goto error; |
1866 | 111 | if (map->n <= 1) |
1867 | 82 | return floyd_warshall(dim, map, exact, project); |
1868 | 29 | |
1869 | 29 | data.list = map->p; |
1870 | 29 | data.check_closed = 0; |
1871 | 29 | g = isl_tarjan_graph_init(map->ctx, map->n, &basic_map_follows, &data); |
1872 | 29 | if (!g) |
1873 | 0 | goto error; |
1874 | 29 | |
1875 | 29 | orig_exact = exact; |
1876 | 29 | if (data.check_closed && !exact4 ) |
1877 | 0 | exact = &local_exact; |
1878 | 29 | |
1879 | 29 | c = 0; |
1880 | 29 | i = 0; |
1881 | 29 | n = map->n; |
1882 | 29 | if (project) |
1883 | 26 | path = isl_map_empty(isl_map_get_space(map)); |
1884 | 3 | else |
1885 | 3 | path = isl_map_empty(isl_space_copy(dim)); |
1886 | 29 | path = anonymize(path); |
1887 | 62 | while (n) { |
1888 | 33 | struct isl_map *comp; |
1889 | 33 | isl_map *path_comp, *path_comb; |
1890 | 33 | comp = isl_map_alloc_space(isl_map_get_space(map), n, 0); |
1891 | 93 | while (g->order[i] != -1) { |
1892 | 60 | comp = isl_map_add_basic_map(comp, |
1893 | 60 | isl_basic_map_copy(map->p[g->order[i]])); |
1894 | 60 | --n; |
1895 | 60 | ++i; |
1896 | 60 | } |
1897 | 33 | path_comp = floyd_warshall(isl_space_copy(dim), |
1898 | 33 | comp, exact, project); |
1899 | 33 | path_comp = anonymize(path_comp); |
1900 | 33 | path_comb = isl_map_apply_range(isl_map_copy(path), |
1901 | 33 | isl_map_copy(path_comp)); |
1902 | 33 | path = isl_map_union(path, path_comp); |
1903 | 33 | path = isl_map_union(path, path_comb); |
1904 | 33 | isl_map_free(comp); |
1905 | 33 | ++i; |
1906 | 33 | ++c; |
1907 | 33 | } |
1908 | 29 | |
1909 | 29 | if (c > 1 && data.check_closed4 && !*exact4 ) { |
1910 | 0 | int closed; |
1911 | 0 |
|
1912 | 0 | closed = isl_map_is_transitively_closed(path); |
1913 | 0 | if (closed < 0) |
1914 | 0 | goto error; |
1915 | 0 | if (!closed) { |
1916 | 0 | isl_tarjan_graph_free(g); |
1917 | 0 | isl_map_free(path); |
1918 | 0 | return floyd_warshall(dim, map, orig_exact, project); |
1919 | 0 | } |
1920 | 29 | } |
1921 | 29 | |
1922 | 29 | isl_tarjan_graph_free(g); |
1923 | 29 | isl_space_free(dim); |
1924 | 29 | |
1925 | 29 | return path; |
1926 | 0 | error: |
1927 | 0 | isl_tarjan_graph_free(g); |
1928 | 0 | isl_space_free(dim); |
1929 | 0 | isl_map_free(path); |
1930 | 0 | return NULL; |
1931 | 29 | } |
1932 | | |
1933 | | /* Given a union of basic maps R = \cup_i R_i \subseteq D \times D, |
1934 | | * construct a map that is an overapproximation of the map |
1935 | | * that takes an element from the space D to another |
1936 | | * element from the same space, such that the difference between |
1937 | | * them is a strictly positive sum of differences between images |
1938 | | * and pre-images in one of the R_i. |
1939 | | * The number of differences in the sum is equated to parameter "param". |
1940 | | * That is, let |
1941 | | * |
1942 | | * \Delta_i = { y - x | (x, y) in R_i } |
1943 | | * |
1944 | | * then the constructed map is an overapproximation of |
1945 | | * |
1946 | | * { (x) -> (x + d) | \exists k_i >= 0, \delta_i \in \Delta_i : |
1947 | | * d = \sum_i k_i \delta_i and k = \sum_i k_i > 0 } |
1948 | | * or |
1949 | | * |
1950 | | * { (x) -> (x + d) | \exists k_i >= 0, \delta_i \in \Delta_i : |
1951 | | * d = \sum_i k_i \delta_i and \sum_i k_i > 0 } |
1952 | | * |
1953 | | * if "project" is set. |
1954 | | * |
1955 | | * If "project" is not set, then |
1956 | | * we construct an extended mapping with an extra coordinate |
1957 | | * that indicates the number of steps taken. In particular, |
1958 | | * the difference in the last coordinate is equal to the number |
1959 | | * of steps taken to move from a domain element to the corresponding |
1960 | | * image element(s). |
1961 | | */ |
1962 | | static __isl_give isl_map *construct_power(__isl_keep isl_map *map, |
1963 | | int *exact, int project) |
1964 | 111 | { |
1965 | 111 | struct isl_map *app = NULL; |
1966 | 111 | isl_space *dim = NULL; |
1967 | 111 | |
1968 | 111 | if (!map) |
1969 | 0 | return NULL; |
1970 | 111 | |
1971 | 111 | dim = isl_map_get_space(map); |
1972 | 111 | |
1973 | 111 | dim = isl_space_add_dims(dim, isl_dim_in, 1); |
1974 | 111 | dim = isl_space_add_dims(dim, isl_dim_out, 1); |
1975 | 111 | |
1976 | 111 | app = construct_power_components(isl_space_copy(dim), map, |
1977 | 111 | exact, project); |
1978 | 111 | |
1979 | 111 | isl_space_free(dim); |
1980 | 111 | |
1981 | 111 | return app; |
1982 | 111 | } |
1983 | | |
1984 | | /* Compute the positive powers of "map", or an overapproximation. |
1985 | | * If the result is exact, then *exact is set to 1. |
1986 | | * |
1987 | | * If project is set, then we are actually interested in the transitive |
1988 | | * closure, so we can use a more relaxed exactness check. |
1989 | | * The lengths of the paths are also projected out instead of being |
1990 | | * encoded as the difference between an extra pair of final coordinates. |
1991 | | */ |
1992 | | static __isl_give isl_map *map_power(__isl_take isl_map *map, |
1993 | | int *exact, int project) |
1994 | 111 | { |
1995 | 111 | struct isl_map *app = NULL; |
1996 | 111 | |
1997 | 111 | if (exact) |
1998 | 17 | *exact = 1; |
1999 | 111 | |
2000 | 111 | if (!map) |
2001 | 0 | return NULL; |
2002 | 111 | |
2003 | 111 | isl_assert(map->ctx, |
2004 | 111 | isl_map_dim(map, isl_dim_in) == isl_map_dim(map, isl_dim_out), |
2005 | 111 | goto error); |
2006 | 111 | |
2007 | 111 | app = construct_power(map, exact, project); |
2008 | 111 | |
2009 | 111 | isl_map_free(map); |
2010 | 111 | return app; |
2011 | 0 | error: |
2012 | 0 | isl_map_free(map); |
2013 | 0 | isl_map_free(app); |
2014 | 0 | return NULL; |
2015 | 111 | } |
2016 | | |
2017 | | /* Compute the positive powers of "map", or an overapproximation. |
2018 | | * The result maps the exponent to a nested copy of the corresponding power. |
2019 | | * If the result is exact, then *exact is set to 1. |
2020 | | * map_power constructs an extended relation with the path lengths |
2021 | | * encoded as the difference between the final coordinates. |
2022 | | * In the final step, this difference is equated to an extra parameter |
2023 | | * and made positive. The extra coordinates are subsequently projected out |
2024 | | * and the parameter is turned into the domain of the result. |
2025 | | */ |
2026 | | __isl_give isl_map *isl_map_power(__isl_take isl_map *map, int *exact) |
2027 | 3 | { |
2028 | 3 | isl_space *target_dim; |
2029 | 3 | isl_space *dim; |
2030 | 3 | isl_map *diff; |
2031 | 3 | unsigned d; |
2032 | 3 | unsigned param; |
2033 | 3 | |
2034 | 3 | if (!map) |
2035 | 0 | return NULL; |
2036 | 3 | |
2037 | 3 | d = isl_map_dim(map, isl_dim_in); |
2038 | 3 | param = isl_map_dim(map, isl_dim_param); |
2039 | 3 | |
2040 | 3 | map = isl_map_compute_divs(map); |
2041 | 3 | map = isl_map_coalesce(map); |
2042 | 3 | |
2043 | 3 | if (isl_map_plain_is_empty(map)) { |
2044 | 0 | map = isl_map_from_range(isl_map_wrap(map)); |
2045 | 0 | map = isl_map_add_dims(map, isl_dim_in, 1); |
2046 | 0 | map = isl_map_set_dim_name(map, isl_dim_in, 0, "k"); |
2047 | 0 | return map; |
2048 | 0 | } |
2049 | 3 | |
2050 | 3 | target_dim = isl_map_get_space(map); |
2051 | 3 | target_dim = isl_space_from_range(isl_space_wrap(target_dim)); |
2052 | 3 | target_dim = isl_space_add_dims(target_dim, isl_dim_in, 1); |
2053 | 3 | target_dim = isl_space_set_dim_name(target_dim, isl_dim_in, 0, "k"); |
2054 | 3 | |
2055 | 3 | map = map_power(map, exact, 0); |
2056 | 3 | |
2057 | 3 | map = isl_map_add_dims(map, isl_dim_param, 1); |
2058 | 3 | dim = isl_map_get_space(map); |
2059 | 3 | diff = equate_parameter_to_length(dim, param); |
2060 | 3 | map = isl_map_intersect(map, diff); |
2061 | 3 | map = isl_map_project_out(map, isl_dim_in, d, 1); |
2062 | 3 | map = isl_map_project_out(map, isl_dim_out, d, 1); |
2063 | 3 | map = isl_map_from_range(isl_map_wrap(map)); |
2064 | 3 | map = isl_map_move_dims(map, isl_dim_in, 0, isl_dim_param, param, 1); |
2065 | 3 | |
2066 | 3 | map = isl_map_reset_space(map, target_dim); |
2067 | 3 | |
2068 | 3 | return map; |
2069 | 3 | } |
2070 | | |
2071 | | /* Compute a relation that maps each element in the range of the input |
2072 | | * relation to the lengths of all paths composed of edges in the input |
2073 | | * relation that end up in the given range element. |
2074 | | * The result may be an overapproximation, in which case *exact is set to 0. |
2075 | | * The resulting relation is very similar to the power relation. |
2076 | | * The difference are that the domain has been projected out, the |
2077 | | * range has become the domain and the exponent is the range instead |
2078 | | * of a parameter. |
2079 | | */ |
2080 | | __isl_give isl_map *isl_map_reaching_path_lengths(__isl_take isl_map *map, |
2081 | | int *exact) |
2082 | 0 | { |
2083 | 0 | isl_space *dim; |
2084 | 0 | isl_map *diff; |
2085 | 0 | unsigned d; |
2086 | 0 | unsigned param; |
2087 | 0 |
|
2088 | 0 | if (!map) |
2089 | 0 | return NULL; |
2090 | 0 | |
2091 | 0 | d = isl_map_dim(map, isl_dim_in); |
2092 | 0 | param = isl_map_dim(map, isl_dim_param); |
2093 | 0 |
|
2094 | 0 | map = isl_map_compute_divs(map); |
2095 | 0 | map = isl_map_coalesce(map); |
2096 | 0 |
|
2097 | 0 | if (isl_map_plain_is_empty(map)) { |
2098 | 0 | if (exact) |
2099 | 0 | *exact = 1; |
2100 | 0 | map = isl_map_project_out(map, isl_dim_out, 0, d); |
2101 | 0 | map = isl_map_add_dims(map, isl_dim_out, 1); |
2102 | 0 | return map; |
2103 | 0 | } |
2104 | 0 |
|
2105 | 0 | map = map_power(map, exact, 0); |
2106 | 0 |
|
2107 | 0 | map = isl_map_add_dims(map, isl_dim_param, 1); |
2108 | 0 | dim = isl_map_get_space(map); |
2109 | 0 | diff = equate_parameter_to_length(dim, param); |
2110 | 0 | map = isl_map_intersect(map, diff); |
2111 | 0 | map = isl_map_project_out(map, isl_dim_in, 0, d + 1); |
2112 | 0 | map = isl_map_project_out(map, isl_dim_out, d, 1); |
2113 | 0 | map = isl_map_reverse(map); |
2114 | 0 | map = isl_map_move_dims(map, isl_dim_out, 0, isl_dim_param, param, 1); |
2115 | 0 |
|
2116 | 0 | return map; |
2117 | 0 | } |
2118 | | |
2119 | | /* Given a map, compute the smallest superset of this map that is of the form |
2120 | | * |
2121 | | * { i -> j : L <= j - i <= U and exists a_p: j_p - i_p = M_p a_p } |
2122 | | * |
2123 | | * (where p ranges over the (non-parametric) dimensions), |
2124 | | * compute the transitive closure of this map, i.e., |
2125 | | * |
2126 | | * { i -> j : exists k > 0: |
2127 | | * k L <= j - i <= k U and exists a: j_p - i_p = M_p a_p } |
2128 | | * |
2129 | | * and intersect domain and range of this transitive closure with |
2130 | | * the given domain and range. |
2131 | | * |
2132 | | * If with_id is set, then try to include as much of the identity mapping |
2133 | | * as possible, by computing |
2134 | | * |
2135 | | * { i -> j : exists k >= 0: |
2136 | | * k L <= j - i <= k U and exists a: j_p - i_p = M_p a_p } |
2137 | | * |
2138 | | * instead (i.e., allow k = 0). |
2139 | | * |
2140 | | * In practice, we compute the difference set |
2141 | | * |
2142 | | * delta = { j - i | i -> j in map }, |
2143 | | * |
2144 | | * look for stride constraint on the individual dimensions and compute |
2145 | | * (constant) lower and upper bounds for each individual dimension, |
2146 | | * adding a constraint for each bound not equal to infinity. |
2147 | | */ |
2148 | | static __isl_give isl_map *box_closure_on_domain(__isl_take isl_map *map, |
2149 | | __isl_take isl_set *dom, __isl_take isl_set *ran, int with_id) |
2150 | 0 | { |
2151 | 0 | int i; |
2152 | 0 | int k; |
2153 | 0 | unsigned d; |
2154 | 0 | unsigned nparam; |
2155 | 0 | unsigned total; |
2156 | 0 | isl_space *dim; |
2157 | 0 | isl_set *delta; |
2158 | 0 | isl_map *app = NULL; |
2159 | 0 | isl_basic_set *aff = NULL; |
2160 | 0 | isl_basic_map *bmap = NULL; |
2161 | 0 | isl_vec *obj = NULL; |
2162 | 0 | isl_int opt; |
2163 | 0 |
|
2164 | 0 | isl_int_init(opt); |
2165 | 0 |
|
2166 | 0 | delta = isl_map_deltas(isl_map_copy(map)); |
2167 | 0 |
|
2168 | 0 | aff = isl_set_affine_hull(isl_set_copy(delta)); |
2169 | 0 | if (!aff) |
2170 | 0 | goto error; |
2171 | 0 | dim = isl_map_get_space(map); |
2172 | 0 | d = isl_space_dim(dim, isl_dim_in); |
2173 | 0 | nparam = isl_space_dim(dim, isl_dim_param); |
2174 | 0 | total = isl_space_dim(dim, isl_dim_all); |
2175 | 0 | bmap = isl_basic_map_alloc_space(dim, |
2176 | 0 | aff->n_div + 1, aff->n_div, 2 * d + 1); |
2177 | 0 | for (i = 0; i < aff->n_div + 1; ++i) { |
2178 | 0 | k = isl_basic_map_alloc_div(bmap); |
2179 | 0 | if (k < 0) |
2180 | 0 | goto error; |
2181 | 0 | isl_int_set_si(bmap->div[k][0], 0); |
2182 | 0 | } |
2183 | 0 | for (i = 0; i < aff->n_eq; ++i) { |
2184 | 0 | if (!isl_basic_set_eq_is_stride(aff, i)) |
2185 | 0 | continue; |
2186 | 0 | k = isl_basic_map_alloc_equality(bmap); |
2187 | 0 | if (k < 0) |
2188 | 0 | goto error; |
2189 | 0 | isl_seq_clr(bmap->eq[k], 1 + nparam); |
2190 | 0 | isl_seq_cpy(bmap->eq[k] + 1 + nparam + d, |
2191 | 0 | aff->eq[i] + 1 + nparam, d); |
2192 | 0 | isl_seq_neg(bmap->eq[k] + 1 + nparam, |
2193 | 0 | aff->eq[i] + 1 + nparam, d); |
2194 | 0 | isl_seq_cpy(bmap->eq[k] + 1 + nparam + 2 * d, |
2195 | 0 | aff->eq[i] + 1 + nparam + d, aff->n_div); |
2196 | 0 | isl_int_set_si(bmap->eq[k][1 + total + aff->n_div], 0); |
2197 | 0 | } |
2198 | 0 | obj = isl_vec_alloc(map->ctx, 1 + nparam + d); |
2199 | 0 | if (!obj) |
2200 | 0 | goto error; |
2201 | 0 | isl_seq_clr(obj->el, 1 + nparam + d); |
2202 | 0 | for (i = 0; i < d; ++ i) { |
2203 | 0 | enum isl_lp_result res; |
2204 | 0 |
|
2205 | 0 | isl_int_set_si(obj->el[1 + nparam + i], 1); |
2206 | 0 |
|
2207 | 0 | res = isl_set_solve_lp(delta, 0, obj->el, map->ctx->one, &opt, |
2208 | 0 | NULL, NULL); |
2209 | 0 | if (res == isl_lp_error) |
2210 | 0 | goto error; |
2211 | 0 | if (res == isl_lp_ok) { |
2212 | 0 | k = isl_basic_map_alloc_inequality(bmap); |
2213 | 0 | if (k < 0) |
2214 | 0 | goto error; |
2215 | 0 | isl_seq_clr(bmap->ineq[k], |
2216 | 0 | 1 + nparam + 2 * d + bmap->n_div); |
2217 | 0 | isl_int_set_si(bmap->ineq[k][1 + nparam + i], -1); |
2218 | 0 | isl_int_set_si(bmap->ineq[k][1 + nparam + d + i], 1); |
2219 | 0 | isl_int_neg(bmap->ineq[k][1 + nparam + 2 * d + aff->n_div], opt); |
2220 | 0 | } |
2221 | 0 |
|
2222 | 0 | res = isl_set_solve_lp(delta, 1, obj->el, map->ctx->one, &opt, |
2223 | 0 | NULL, NULL); |
2224 | 0 | if (res == isl_lp_error) |
2225 | 0 | goto error; |
2226 | 0 | if (res == isl_lp_ok) { |
2227 | 0 | k = isl_basic_map_alloc_inequality(bmap); |
2228 | 0 | if (k < 0) |
2229 | 0 | goto error; |
2230 | 0 | isl_seq_clr(bmap->ineq[k], |
2231 | 0 | 1 + nparam + 2 * d + bmap->n_div); |
2232 | 0 | isl_int_set_si(bmap->ineq[k][1 + nparam + i], 1); |
2233 | 0 | isl_int_set_si(bmap->ineq[k][1 + nparam + d + i], -1); |
2234 | 0 | isl_int_set(bmap->ineq[k][1 + nparam + 2 * d + aff->n_div], opt); |
2235 | 0 | } |
2236 | 0 |
|
2237 | 0 | isl_int_set_si(obj->el[1 + nparam + i], 0); |
2238 | 0 | } |
2239 | 0 | k = isl_basic_map_alloc_inequality(bmap); |
2240 | 0 | if (k < 0) |
2241 | 0 | goto error; |
2242 | 0 | isl_seq_clr(bmap->ineq[k], |
2243 | 0 | 1 + nparam + 2 * d + bmap->n_div); |
2244 | 0 | if (!with_id) |
2245 | 0 | isl_int_set_si(bmap->ineq[k][0], -1); |
2246 | 0 | isl_int_set_si(bmap->ineq[k][1 + nparam + 2 * d + aff->n_div], 1); |
2247 | 0 |
|
2248 | 0 | app = isl_map_from_domain_and_range(dom, ran); |
2249 | 0 |
|
2250 | 0 | isl_vec_free(obj); |
2251 | 0 | isl_basic_set_free(aff); |
2252 | 0 | isl_map_free(map); |
2253 | 0 | bmap = isl_basic_map_finalize(bmap); |
2254 | 0 | isl_set_free(delta); |
2255 | 0 | isl_int_clear(opt); |
2256 | 0 |
|
2257 | 0 | map = isl_map_from_basic_map(bmap); |
2258 | 0 | map = isl_map_intersect(map, app); |
2259 | 0 |
|
2260 | 0 | return map; |
2261 | 0 | error: |
2262 | 0 | isl_vec_free(obj); |
2263 | 0 | isl_basic_map_free(bmap); |
2264 | 0 | isl_basic_set_free(aff); |
2265 | 0 | isl_set_free(dom); |
2266 | 0 | isl_set_free(ran); |
2267 | 0 | isl_map_free(map); |
2268 | 0 | isl_set_free(delta); |
2269 | 0 | isl_int_clear(opt); |
2270 | 0 | return NULL; |
2271 | 0 | } |
2272 | | |
2273 | | /* Given a map, compute the smallest superset of this map that is of the form |
2274 | | * |
2275 | | * { i -> j : L <= j - i <= U and exists a_p: j_p - i_p = M_p a_p } |
2276 | | * |
2277 | | * (where p ranges over the (non-parametric) dimensions), |
2278 | | * compute the transitive closure of this map, i.e., |
2279 | | * |
2280 | | * { i -> j : exists k > 0: |
2281 | | * k L <= j - i <= k U and exists a: j_p - i_p = M_p a_p } |
2282 | | * |
2283 | | * and intersect domain and range of this transitive closure with |
2284 | | * domain and range of the original map. |
2285 | | */ |
2286 | | static __isl_give isl_map *box_closure(__isl_take isl_map *map) |
2287 | 0 | { |
2288 | 0 | isl_set *domain; |
2289 | 0 | isl_set *range; |
2290 | 0 |
|
2291 | 0 | domain = isl_map_domain(isl_map_copy(map)); |
2292 | 0 | domain = isl_set_coalesce(domain); |
2293 | 0 | range = isl_map_range(isl_map_copy(map)); |
2294 | 0 | range = isl_set_coalesce(range); |
2295 | 0 |
|
2296 | 0 | return box_closure_on_domain(map, domain, range, 0); |
2297 | 0 | } |
2298 | | |
2299 | | /* Given a map, compute the smallest superset of this map that is of the form |
2300 | | * |
2301 | | * { i -> j : L <= j - i <= U and exists a_p: j_p - i_p = M_p a_p } |
2302 | | * |
2303 | | * (where p ranges over the (non-parametric) dimensions), |
2304 | | * compute the transitive and partially reflexive closure of this map, i.e., |
2305 | | * |
2306 | | * { i -> j : exists k >= 0: |
2307 | | * k L <= j - i <= k U and exists a: j_p - i_p = M_p a_p } |
2308 | | * |
2309 | | * and intersect domain and range of this transitive closure with |
2310 | | * the given domain. |
2311 | | */ |
2312 | | static __isl_give isl_map *box_closure_with_identity(__isl_take isl_map *map, |
2313 | | __isl_take isl_set *dom) |
2314 | 0 | { |
2315 | 0 | return box_closure_on_domain(map, dom, isl_set_copy(dom), 1); |
2316 | 0 | } |
2317 | | |
2318 | | /* Check whether app is the transitive closure of map. |
2319 | | * In particular, check that app is acyclic and, if so, |
2320 | | * check that |
2321 | | * |
2322 | | * app \subset (map \cup (map \circ app)) |
2323 | | */ |
2324 | | static int check_exactness_omega(__isl_keep isl_map *map, |
2325 | | __isl_keep isl_map *app) |
2326 | 0 | { |
2327 | 0 | isl_set *delta; |
2328 | 0 | int i; |
2329 | 0 | int is_empty, is_exact; |
2330 | 0 | unsigned d; |
2331 | 0 | isl_map *test; |
2332 | 0 |
|
2333 | 0 | delta = isl_map_deltas(isl_map_copy(app)); |
2334 | 0 | d = isl_set_dim(delta, isl_dim_set); |
2335 | 0 | for (i = 0; i < d; ++i) |
2336 | 0 | delta = isl_set_fix_si(delta, isl_dim_set, i, 0); |
2337 | 0 | is_empty = isl_set_is_empty(delta); |
2338 | 0 | isl_set_free(delta); |
2339 | 0 | if (is_empty < 0) |
2340 | 0 | return -1; |
2341 | 0 | if (!is_empty) |
2342 | 0 | return 0; |
2343 | 0 | |
2344 | 0 | test = isl_map_apply_range(isl_map_copy(app), isl_map_copy(map)); |
2345 | 0 | test = isl_map_union(test, isl_map_copy(map)); |
2346 | 0 | is_exact = isl_map_is_subset(app, test); |
2347 | 0 | isl_map_free(test); |
2348 | 0 |
|
2349 | 0 | return is_exact; |
2350 | 0 | } |
2351 | | |
2352 | | /* Check if basic map M_i can be combined with all the other |
2353 | | * basic maps such that |
2354 | | * |
2355 | | * (\cup_j M_j)^+ |
2356 | | * |
2357 | | * can be computed as |
2358 | | * |
2359 | | * M_i \cup (\cup_{j \ne i} M_i^* \circ M_j \circ M_i^*)^+ |
2360 | | * |
2361 | | * In particular, check if we can compute a compact representation |
2362 | | * of |
2363 | | * |
2364 | | * M_i^* \circ M_j \circ M_i^* |
2365 | | * |
2366 | | * for each j != i. |
2367 | | * Let M_i^? be an extension of M_i^+ that allows paths |
2368 | | * of length zero, i.e., the result of box_closure(., 1). |
2369 | | * The criterion, as proposed by Kelly et al., is that |
2370 | | * id = M_i^? - M_i^+ can be represented as a basic map |
2371 | | * and that |
2372 | | * |
2373 | | * id \circ M_j \circ id = M_j |
2374 | | * |
2375 | | * for each j != i. |
2376 | | * |
2377 | | * If this function returns 1, then tc and qc are set to |
2378 | | * M_i^+ and M_i^?, respectively. |
2379 | | */ |
2380 | | static int can_be_split_off(__isl_keep isl_map *map, int i, |
2381 | | __isl_give isl_map **tc, __isl_give isl_map **qc) |
2382 | 0 | { |
2383 | 0 | isl_map *map_i, *id = NULL; |
2384 | 0 | int j = -1; |
2385 | 0 | isl_set *C; |
2386 | 0 |
|
2387 | 0 | *tc = NULL; |
2388 | 0 | *qc = NULL; |
2389 | 0 |
|
2390 | 0 | C = isl_set_union(isl_map_domain(isl_map_copy(map)), |
2391 | 0 | isl_map_range(isl_map_copy(map))); |
2392 | 0 | C = isl_set_from_basic_set(isl_set_simple_hull(C)); |
2393 | 0 | if (!C) |
2394 | 0 | goto error; |
2395 | 0 | |
2396 | 0 | map_i = isl_map_from_basic_map(isl_basic_map_copy(map->p[i])); |
2397 | 0 | *tc = box_closure(isl_map_copy(map_i)); |
2398 | 0 | *qc = box_closure_with_identity(map_i, C); |
2399 | 0 | id = isl_map_subtract(isl_map_copy(*qc), isl_map_copy(*tc)); |
2400 | 0 |
|
2401 | 0 | if (!id || !*qc) |
2402 | 0 | goto error; |
2403 | 0 | if (id->n != 1 || (*qc)->n != 1) |
2404 | 0 | goto done; |
2405 | 0 | |
2406 | 0 | for (j = 0; j < map->n; ++j) { |
2407 | 0 | isl_map *map_j, *test; |
2408 | 0 | int is_ok; |
2409 | 0 |
|
2410 | 0 | if (i == j) |
2411 | 0 | continue; |
2412 | 0 | map_j = isl_map_from_basic_map( |
2413 | 0 | isl_basic_map_copy(map->p[j])); |
2414 | 0 | test = isl_map_apply_range(isl_map_copy(id), |
2415 | 0 | isl_map_copy(map_j)); |
2416 | 0 | test = isl_map_apply_range(test, isl_map_copy(id)); |
2417 | 0 | is_ok = isl_map_is_equal(test, map_j); |
2418 | 0 | isl_map_free(map_j); |
2419 | 0 | isl_map_free(test); |
2420 | 0 | if (is_ok < 0) |
2421 | 0 | goto error; |
2422 | 0 | if (!is_ok) |
2423 | 0 | break; |
2424 | 0 | } |
2425 | 0 |
|
2426 | 0 | done: |
2427 | 0 | isl_map_free(id); |
2428 | 0 | if (j == map->n) |
2429 | 0 | return 1; |
2430 | 0 | |
2431 | 0 | isl_map_free(*qc); |
2432 | 0 | isl_map_free(*tc); |
2433 | 0 | *qc = NULL; |
2434 | 0 | *tc = NULL; |
2435 | 0 |
|
2436 | 0 | return 0; |
2437 | 0 | error: |
2438 | 0 | isl_map_free(id); |
2439 | 0 | isl_map_free(*qc); |
2440 | 0 | isl_map_free(*tc); |
2441 | 0 | *qc = NULL; |
2442 | 0 | *tc = NULL; |
2443 | 0 | return -1; |
2444 | 0 | } |
2445 | | |
2446 | | static __isl_give isl_map *box_closure_with_check(__isl_take isl_map *map, |
2447 | | int *exact) |
2448 | 0 | { |
2449 | 0 | isl_map *app; |
2450 | 0 |
|
2451 | 0 | app = box_closure(isl_map_copy(map)); |
2452 | 0 | if (exact) |
2453 | 0 | *exact = check_exactness_omega(map, app); |
2454 | 0 |
|
2455 | 0 | isl_map_free(map); |
2456 | 0 | return app; |
2457 | 0 | } |
2458 | | |
2459 | | /* Compute an overapproximation of the transitive closure of "map" |
2460 | | * using a variation of the algorithm from |
2461 | | * "Transitive Closure of Infinite Graphs and its Applications" |
2462 | | * by Kelly et al. |
2463 | | * |
2464 | | * We first check whether we can can split of any basic map M_i and |
2465 | | * compute |
2466 | | * |
2467 | | * (\cup_j M_j)^+ |
2468 | | * |
2469 | | * as |
2470 | | * |
2471 | | * M_i \cup (\cup_{j \ne i} M_i^* \circ M_j \circ M_i^*)^+ |
2472 | | * |
2473 | | * using a recursive call on the remaining map. |
2474 | | * |
2475 | | * If not, we simply call box_closure on the whole map. |
2476 | | */ |
2477 | | static __isl_give isl_map *transitive_closure_omega(__isl_take isl_map *map, |
2478 | | int *exact) |
2479 | 0 | { |
2480 | 0 | int i, j; |
2481 | 0 | int exact_i; |
2482 | 0 | isl_map *app; |
2483 | 0 |
|
2484 | 0 | if (!map) |
2485 | 0 | return NULL; |
2486 | 0 | if (map->n == 1) |
2487 | 0 | return box_closure_with_check(map, exact); |
2488 | 0 | |
2489 | 0 | for (i = 0; i < map->n; ++i) { |
2490 | 0 | int ok; |
2491 | 0 | isl_map *qc, *tc; |
2492 | 0 | ok = can_be_split_off(map, i, &tc, &qc); |
2493 | 0 | if (ok < 0) |
2494 | 0 | goto error; |
2495 | 0 | if (!ok) |
2496 | 0 | continue; |
2497 | 0 | |
2498 | 0 | app = isl_map_alloc_space(isl_map_get_space(map), map->n - 1, 0); |
2499 | 0 |
|
2500 | 0 | for (j = 0; j < map->n; ++j) { |
2501 | 0 | if (j == i) |
2502 | 0 | continue; |
2503 | 0 | app = isl_map_add_basic_map(app, |
2504 | 0 | isl_basic_map_copy(map->p[j])); |
2505 | 0 | } |
2506 | 0 |
|
2507 | 0 | app = isl_map_apply_range(isl_map_copy(qc), app); |
2508 | 0 | app = isl_map_apply_range(app, qc); |
2509 | 0 |
|
2510 | 0 | app = isl_map_union(tc, transitive_closure_omega(app, NULL)); |
2511 | 0 | exact_i = check_exactness_omega(map, app); |
2512 | 0 | if (exact_i == 1) { |
2513 | 0 | if (exact) |
2514 | 0 | *exact = exact_i; |
2515 | 0 | isl_map_free(map); |
2516 | 0 | return app; |
2517 | 0 | } |
2518 | 0 | isl_map_free(app); |
2519 | 0 | if (exact_i < 0) |
2520 | 0 | goto error; |
2521 | 0 | } |
2522 | 0 |
|
2523 | 0 | return box_closure_with_check(map, exact); |
2524 | 0 | error: |
2525 | 0 | isl_map_free(map); |
2526 | 0 | return NULL; |
2527 | 0 | } |
2528 | | |
2529 | | /* Compute the transitive closure of "map", or an overapproximation. |
2530 | | * If the result is exact, then *exact is set to 1. |
2531 | | * Simply use map_power to compute the powers of map, but tell |
2532 | | * it to project out the lengths of the paths instead of equating |
2533 | | * the length to a parameter. |
2534 | | */ |
2535 | | __isl_give isl_map *isl_map_transitive_closure(__isl_take isl_map *map, |
2536 | | int *exact) |
2537 | 116 | { |
2538 | 116 | isl_space *target_dim; |
2539 | 116 | int closed; |
2540 | 116 | |
2541 | 116 | if (!map) |
2542 | 0 | goto error; |
2543 | 116 | |
2544 | 116 | if (map->ctx->opt->closure == ISL_CLOSURE_BOX) |
2545 | 116 | return transitive_closure_omega(map, exact)0 ; |
2546 | 116 | |
2547 | 116 | map = isl_map_compute_divs(map); |
2548 | 116 | map = isl_map_coalesce(map); |
2549 | 116 | closed = isl_map_is_transitively_closed(map); |
2550 | 116 | if (closed < 0) |
2551 | 0 | goto error; |
2552 | 116 | if (closed) { |
2553 | 8 | if (exact) |
2554 | 6 | *exact = 1; |
2555 | 8 | return map; |
2556 | 8 | } |
2557 | 108 | |
2558 | 108 | target_dim = isl_map_get_space(map); |
2559 | 108 | map = map_power(map, exact, 1); |
2560 | 108 | map = isl_map_reset_space(map, target_dim); |
2561 | 108 | |
2562 | 108 | return map; |
2563 | 0 | error: |
2564 | 0 | isl_map_free(map); |
2565 | 0 | return NULL; |
2566 | 108 | } |
2567 | | |
2568 | | static isl_stat inc_count(__isl_take isl_map *map, void *user) |
2569 | 185 | { |
2570 | 185 | int *n = user; |
2571 | 185 | |
2572 | 185 | *n += map->n; |
2573 | 185 | |
2574 | 185 | isl_map_free(map); |
2575 | 185 | |
2576 | 185 | return isl_stat_ok; |
2577 | 185 | } |
2578 | | |
2579 | | static isl_stat collect_basic_map(__isl_take isl_map *map, void *user) |
2580 | 127 | { |
2581 | 127 | int i; |
2582 | 127 | isl_basic_map ***next = user; |
2583 | 127 | |
2584 | 295 | for (i = 0; i < map->n; ++i168 ) { |
2585 | 168 | **next = isl_basic_map_copy(map->p[i]); |
2586 | 168 | if (!**next) |
2587 | 0 | goto error; |
2588 | 168 | (*next)++; |
2589 | 168 | } |
2590 | 127 | |
2591 | 127 | isl_map_free(map); |
2592 | 127 | return isl_stat_ok; |
2593 | 0 | error: |
2594 | 0 | isl_map_free(map); |
2595 | 0 | return isl_stat_error; |
2596 | 127 | } |
2597 | | |
2598 | | /* Perform Floyd-Warshall on the given list of basic relations. |
2599 | | * The basic relations may live in different dimensions, |
2600 | | * but basic relations that get assigned to the diagonal of the |
2601 | | * grid have domains and ranges of the same dimension and so |
2602 | | * the standard algorithm can be used because the nested transitive |
2603 | | * closures are only applied to diagonal elements and because all |
2604 | | * compositions are peformed on relations with compatible domains and ranges. |
2605 | | */ |
2606 | | static __isl_give isl_union_map *union_floyd_warshall_on_list(isl_ctx *ctx, |
2607 | | __isl_keep isl_basic_map **list, int n, int *exact) |
2608 | 94 | { |
2609 | 94 | int i, j, k; |
2610 | 94 | int n_group; |
2611 | 94 | int *group = NULL; |
2612 | 94 | isl_set **set = NULL; |
2613 | 94 | isl_map ***grid = NULL; |
2614 | 94 | isl_union_map *app; |
2615 | 94 | |
2616 | 94 | group = setup_groups(ctx, list, n, &set, &n_group); |
2617 | 94 | if (!group) |
2618 | 0 | goto error; |
2619 | 94 | |
2620 | 94 | grid = isl_calloc_array(ctx, isl_map **, n_group); |
2621 | 94 | if (!grid) |
2622 | 0 | goto error; |
2623 | 189 | for (i = 0; 94 i < n_group; ++i95 ) { |
2624 | 95 | grid[i] = isl_calloc_array(ctx, isl_map *, n_group); |
2625 | 95 | if (!grid[i]) |
2626 | 0 | goto error; |
2627 | 192 | for (j = 0; 95 j < n_group; ++j97 ) { |
2628 | 97 | isl_space *dim1, *dim2, *dim; |
2629 | 97 | dim1 = isl_space_reverse(isl_set_get_space(set[i])); |
2630 | 97 | dim2 = isl_set_get_space(set[j]); |
2631 | 97 | dim = isl_space_join(dim1, dim2); |
2632 | 97 | grid[i][j] = isl_map_empty(dim); |
2633 | 97 | } |
2634 | 95 | } |
2635 | 94 | |
2636 | 207 | for (k = 0; 94 k < n; ++k113 ) { |
2637 | 113 | i = group[2 * k]; |
2638 | 113 | j = group[2 * k + 1]; |
2639 | 113 | grid[i][j] = isl_map_union(grid[i][j], |
2640 | 113 | isl_map_from_basic_map( |
2641 | 113 | isl_basic_map_copy(list[k]))); |
2642 | 113 | } |
2643 | 94 | |
2644 | 94 | floyd_warshall_iterate(grid, n_group, exact); |
2645 | 94 | |
2646 | 94 | app = isl_union_map_empty(isl_map_get_space(grid[0][0])); |
2647 | 94 | |
2648 | 189 | for (i = 0; i < n_group; ++i95 ) { |
2649 | 192 | for (j = 0; j < n_group; ++j97 ) |
2650 | 97 | app = isl_union_map_add_map(app, grid[i][j]); |
2651 | 95 | free(grid[i]); |
2652 | 95 | } |
2653 | 94 | free(grid); |
2654 | 94 | |
2655 | 320 | for (i = 0; i < 2 * n; ++i226 ) |
2656 | 226 | isl_set_free(set[i]); |
2657 | 94 | free(set); |
2658 | 94 | |
2659 | 94 | free(group); |
2660 | 94 | return app; |
2661 | 0 | error: |
2662 | 0 | if (grid) |
2663 | 0 | for (i = 0; i < n_group; ++i) { |
2664 | 0 | if (!grid[i]) |
2665 | 0 | continue; |
2666 | 0 | for (j = 0; j < n_group; ++j) |
2667 | 0 | isl_map_free(grid[i][j]); |
2668 | 0 | free(grid[i]); |
2669 | 0 | } |
2670 | 0 | free(grid); |
2671 | 0 | if (set) { |
2672 | 0 | for (i = 0; i < 2 * n; ++i) |
2673 | 0 | isl_set_free(set[i]); |
2674 | 0 | free(set); |
2675 | 0 | } |
2676 | 0 | free(group); |
2677 | 0 | return NULL; |
2678 | 94 | } |
2679 | | |
2680 | | /* Perform Floyd-Warshall on the given union relation. |
2681 | | * The implementation is very similar to that for non-unions. |
2682 | | * The main difference is that it is applied unconditionally. |
2683 | | * We first extract a list of basic maps from the union map |
2684 | | * and then perform the algorithm on this list. |
2685 | | */ |
2686 | | static __isl_give isl_union_map *union_floyd_warshall( |
2687 | | __isl_take isl_union_map *umap, int *exact) |
2688 | 94 | { |
2689 | 94 | int i, n; |
2690 | 94 | isl_ctx *ctx; |
2691 | 94 | isl_basic_map **list = NULL; |
2692 | 94 | isl_basic_map **next; |
2693 | 94 | isl_union_map *res; |
2694 | 94 | |
2695 | 94 | n = 0; |
2696 | 94 | if (isl_union_map_foreach_map(umap, inc_count, &n) < 0) |
2697 | 0 | goto error; |
2698 | 94 | |
2699 | 94 | ctx = isl_union_map_get_ctx(umap); |
2700 | 94 | list = isl_calloc_array(ctx, isl_basic_map *, n); |
2701 | 94 | if (!list) |
2702 | 0 | goto error; |
2703 | 94 | |
2704 | 94 | next = list; |
2705 | 94 | if (isl_union_map_foreach_map(umap, collect_basic_map, &next) < 0) |
2706 | 0 | goto error; |
2707 | 94 | |
2708 | 94 | res = union_floyd_warshall_on_list(ctx, list, n, exact); |
2709 | 94 | |
2710 | 94 | if (list) { |
2711 | 207 | for (i = 0; i < n; ++i113 ) |
2712 | 113 | isl_basic_map_free(list[i]); |
2713 | 94 | free(list); |
2714 | 94 | } |
2715 | 94 | |
2716 | 94 | isl_union_map_free(umap); |
2717 | 94 | return res; |
2718 | 0 | error: |
2719 | 0 | if (list) { |
2720 | 0 | for (i = 0; i < n; ++i) |
2721 | 0 | isl_basic_map_free(list[i]); |
2722 | 0 | free(list); |
2723 | 0 | } |
2724 | 0 | isl_union_map_free(umap); |
2725 | 0 | return NULL; |
2726 | 94 | } |
2727 | | |
2728 | | /* Decompose the give union relation into strongly connected components. |
2729 | | * The implementation is essentially the same as that of |
2730 | | * construct_power_components with the major difference that all |
2731 | | * operations are performed on union maps. |
2732 | | */ |
2733 | | static __isl_give isl_union_map *union_components( |
2734 | | __isl_take isl_union_map *umap, int *exact) |
2735 | 71 | { |
2736 | 71 | int i; |
2737 | 71 | int n; |
2738 | 71 | isl_ctx *ctx; |
2739 | 71 | isl_basic_map **list = NULL; |
2740 | 71 | isl_basic_map **next; |
2741 | 71 | isl_union_map *path = NULL; |
2742 | 71 | struct isl_tc_follows_data data; |
2743 | 71 | struct isl_tarjan_graph *g = NULL; |
2744 | 71 | int c, l; |
2745 | 71 | int recheck = 0; |
2746 | 71 | |
2747 | 71 | n = 0; |
2748 | 71 | if (isl_union_map_foreach_map(umap, inc_count, &n) < 0) |
2749 | 0 | goto error; |
2750 | 71 | |
2751 | 71 | if (n == 0) |
2752 | 0 | return umap; |
2753 | 71 | if (n <= 1) |
2754 | 58 | return union_floyd_warshall(umap, exact); |
2755 | 13 | |
2756 | 13 | ctx = isl_union_map_get_ctx(umap); |
2757 | 13 | list = isl_calloc_array(ctx, isl_basic_map *, n); |
2758 | 13 | if (!list) |
2759 | 0 | goto error; |
2760 | 13 | |
2761 | 13 | next = list; |
2762 | 13 | if (isl_union_map_foreach_map(umap, collect_basic_map, &next) < 0) |
2763 | 0 | goto error; |
2764 | 13 | |
2765 | 13 | data.list = list; |
2766 | 13 | data.check_closed = 0; |
2767 | 13 | g = isl_tarjan_graph_init(ctx, n, &basic_map_follows, &data); |
2768 | 13 | if (!g) |
2769 | 0 | goto error; |
2770 | 13 | |
2771 | 13 | c = 0; |
2772 | 13 | i = 0; |
2773 | 13 | l = n; |
2774 | 13 | path = isl_union_map_empty(isl_union_map_get_space(umap)); |
2775 | 49 | while (l) { |
2776 | 36 | isl_union_map *comp; |
2777 | 36 | isl_union_map *path_comp, *path_comb; |
2778 | 36 | comp = isl_union_map_empty(isl_union_map_get_space(umap)); |
2779 | 91 | while (g->order[i] != -1) { |
2780 | 55 | comp = isl_union_map_add_map(comp, |
2781 | 55 | isl_map_from_basic_map( |
2782 | 55 | isl_basic_map_copy(list[g->order[i]]))); |
2783 | 55 | --l; |
2784 | 55 | ++i; |
2785 | 55 | } |
2786 | 36 | path_comp = union_floyd_warshall(comp, exact); |
2787 | 36 | path_comb = isl_union_map_apply_range(isl_union_map_copy(path), |
2788 | 36 | isl_union_map_copy(path_comp)); |
2789 | 36 | path = isl_union_map_union(path, path_comp); |
2790 | 36 | path = isl_union_map_union(path, path_comb); |
2791 | 36 | ++i; |
2792 | 36 | ++c; |
2793 | 36 | } |
2794 | 13 | |
2795 | 13 | if (c > 1 && data.check_closed8 && !*exact0 ) { |
2796 | 0 | int closed; |
2797 | 0 |
|
2798 | 0 | closed = isl_union_map_is_transitively_closed(path); |
2799 | 0 | if (closed < 0) |
2800 | 0 | goto error; |
2801 | 0 | recheck = !closed; |
2802 | 0 | } |
2803 | 13 | |
2804 | 13 | isl_tarjan_graph_free(g); |
2805 | 13 | |
2806 | 68 | for (i = 0; i < n; ++i55 ) |
2807 | 55 | isl_basic_map_free(list[i]); |
2808 | 13 | free(list); |
2809 | 13 | |
2810 | 13 | if (recheck) { |
2811 | 0 | isl_union_map_free(path); |
2812 | 0 | return union_floyd_warshall(umap, exact); |
2813 | 0 | } |
2814 | 13 | |
2815 | 13 | isl_union_map_free(umap); |
2816 | 13 | |
2817 | 13 | return path; |
2818 | 0 | error: |
2819 | 0 | isl_tarjan_graph_free(g); |
2820 | 0 | if (list) { |
2821 | 0 | for (i = 0; i < n; ++i) |
2822 | 0 | isl_basic_map_free(list[i]); |
2823 | 0 | free(list); |
2824 | 0 | } |
2825 | 0 | isl_union_map_free(umap); |
2826 | 0 | isl_union_map_free(path); |
2827 | 0 | return NULL; |
2828 | 13 | } |
2829 | | |
2830 | | /* Compute the transitive closure of "umap", or an overapproximation. |
2831 | | * If the result is exact, then *exact is set to 1. |
2832 | | */ |
2833 | | __isl_give isl_union_map *isl_union_map_transitive_closure( |
2834 | | __isl_take isl_union_map *umap, int *exact) |
2835 | 73 | { |
2836 | 73 | int closed; |
2837 | 73 | |
2838 | 73 | if (!umap) |
2839 | 0 | return NULL; |
2840 | 73 | |
2841 | 73 | if (exact) |
2842 | 0 | *exact = 1; |
2843 | 73 | |
2844 | 73 | umap = isl_union_map_compute_divs(umap); |
2845 | 73 | umap = isl_union_map_coalesce(umap); |
2846 | 73 | closed = isl_union_map_is_transitively_closed(umap); |
2847 | 73 | if (closed < 0) |
2848 | 0 | goto error; |
2849 | 73 | if (closed) |
2850 | 2 | return umap; |
2851 | 71 | umap = union_components(umap, exact); |
2852 | 71 | return umap; |
2853 | 0 | error: |
2854 | 0 | isl_union_map_free(umap); |
2855 | 0 | return NULL; |
2856 | 71 | } |
2857 | | |
2858 | | struct isl_union_power { |
2859 | | isl_union_map *pow; |
2860 | | int *exact; |
2861 | | }; |
2862 | | |
2863 | | static isl_stat power(__isl_take isl_map *map, void *user) |
2864 | 0 | { |
2865 | 0 | struct isl_union_power *up = user; |
2866 | 0 |
|
2867 | 0 | map = isl_map_power(map, up->exact); |
2868 | 0 | up->pow = isl_union_map_from_map(map); |
2869 | 0 |
|
2870 | 0 | return isl_stat_error; |
2871 | 0 | } |
2872 | | |
2873 | | /* Construct a map [x] -> [x+1], with parameters prescribed by "dim". |
2874 | | */ |
2875 | | static __isl_give isl_union_map *increment(__isl_take isl_space *dim) |
2876 | 0 | { |
2877 | 0 | int k; |
2878 | 0 | isl_basic_map *bmap; |
2879 | 0 |
|
2880 | 0 | dim = isl_space_add_dims(dim, isl_dim_in, 1); |
2881 | 0 | dim = isl_space_add_dims(dim, isl_dim_out, 1); |
2882 | 0 | bmap = isl_basic_map_alloc_space(dim, 0, 1, 0); |
2883 | 0 | k = isl_basic_map_alloc_equality(bmap); |
2884 | 0 | if (k < 0) |
2885 | 0 | goto error; |
2886 | 0 | isl_seq_clr(bmap->eq[k], isl_basic_map_total_dim(bmap)); |
2887 | 0 | isl_int_set_si(bmap->eq[k][0], 1); |
2888 | 0 | isl_int_set_si(bmap->eq[k][isl_basic_map_offset(bmap, isl_dim_in)], 1); |
2889 | 0 | isl_int_set_si(bmap->eq[k][isl_basic_map_offset(bmap, isl_dim_out)], -1); |
2890 | 0 | return isl_union_map_from_map(isl_map_from_basic_map(bmap)); |
2891 | 0 | error: |
2892 | 0 | isl_basic_map_free(bmap); |
2893 | 0 | return NULL; |
2894 | 0 | } |
2895 | | |
2896 | | /* Construct a map [[x]->[y]] -> [y-x], with parameters prescribed by "dim". |
2897 | | */ |
2898 | | static __isl_give isl_union_map *deltas_map(__isl_take isl_space *dim) |
2899 | 0 | { |
2900 | 0 | isl_basic_map *bmap; |
2901 | 0 |
|
2902 | 0 | dim = isl_space_add_dims(dim, isl_dim_in, 1); |
2903 | 0 | dim = isl_space_add_dims(dim, isl_dim_out, 1); |
2904 | 0 | bmap = isl_basic_map_universe(dim); |
2905 | 0 | bmap = isl_basic_map_deltas_map(bmap); |
2906 | 0 |
|
2907 | 0 | return isl_union_map_from_map(isl_map_from_basic_map(bmap)); |
2908 | 0 | } |
2909 | | |
2910 | | /* Compute the positive powers of "map", or an overapproximation. |
2911 | | * The result maps the exponent to a nested copy of the corresponding power. |
2912 | | * If the result is exact, then *exact is set to 1. |
2913 | | */ |
2914 | | __isl_give isl_union_map *isl_union_map_power(__isl_take isl_union_map *umap, |
2915 | | int *exact) |
2916 | 0 | { |
2917 | 0 | int n; |
2918 | 0 | isl_union_map *inc; |
2919 | 0 | isl_union_map *dm; |
2920 | 0 |
|
2921 | 0 | if (!umap) |
2922 | 0 | return NULL; |
2923 | 0 | n = isl_union_map_n_map(umap); |
2924 | 0 | if (n == 0) |
2925 | 0 | return umap; |
2926 | 0 | if (n == 1) { |
2927 | 0 | struct isl_union_power up = { NULL, exact }; |
2928 | 0 | isl_union_map_foreach_map(umap, &power, &up); |
2929 | 0 | isl_union_map_free(umap); |
2930 | 0 | return up.pow; |
2931 | 0 | } |
2932 | 0 | inc = increment(isl_union_map_get_space(umap)); |
2933 | 0 | umap = isl_union_map_product(inc, umap); |
2934 | 0 | umap = isl_union_map_transitive_closure(umap, exact); |
2935 | 0 | umap = isl_union_map_zip(umap); |
2936 | 0 | dm = deltas_map(isl_union_map_get_space(umap)); |
2937 | 0 | umap = isl_union_map_apply_domain(umap, dm); |
2938 | 0 | |
2939 | 0 | return umap; |
2940 | 0 | } |
2941 | | |
2942 | | #undef TYPE |
2943 | 1 | #define TYPE isl_map |
2944 | | #include "isl_power_templ.c" |
2945 | | |
2946 | | #undef TYPE |
2947 | 0 | #define TYPE isl_union_map |
2948 | | #include "isl_power_templ.c" |