6 #include "isl_map_private.h"
8 struct isl_basic_map *isl_basic_map_implicit_equalities(struct isl_ctx *ctx,
9 struct isl_basic_map *bmap)
17 if (F_ISSET(bmap, ISL_BASIC_MAP_EMPTY))
19 if (F_ISSET(bmap, ISL_BASIC_MAP_NO_IMPLICIT))
23 for (i = 0; i < bmap->n_ineq; ++i) {
24 enum isl_lp_result res;
25 res = isl_solve_lp(bmap, 1, bmap->ineq[i]+1, ctx->one, &opt);
26 if (res == isl_lp_unbounded)
28 if (res == isl_lp_error)
30 if (res == isl_lp_empty) {
31 bmap = isl_basic_map_set_to_empty(ctx, bmap);
34 isl_int_add(opt, opt, bmap->ineq[i][0]);
35 if (isl_int_is_zero(opt)) {
36 isl_basic_map_inequality_to_equality(ctx, bmap, i);
42 F_SET(bmap, ISL_BASIC_MAP_NO_IMPLICIT);
46 isl_basic_map_free(ctx, bmap);
50 struct isl_basic_map *isl_basic_map_affine_hull(struct isl_ctx *ctx,
51 struct isl_basic_map *bmap)
53 bmap = isl_basic_map_implicit_equalities(ctx, bmap);
56 if (bmap->n_ineq == 0)
58 bmap = isl_basic_map_cow(ctx, bmap);
61 isl_basic_map_free_inequality(ctx, bmap, bmap->n_ineq);
62 return isl_basic_map_finalize(ctx, bmap);
65 struct isl_basic_set *isl_basic_set_affine_hull(struct isl_ctx *ctx,
66 struct isl_basic_set *bset)
68 return (struct isl_basic_set *)
69 isl_basic_map_affine_hull(ctx, (struct isl_basic_map *)bset);
72 /* Make eq[row][col] of both bmaps equal so we can add the row
73 * add the column to the common matrix.
74 * Note that because of the echelon form, the columns of row row
75 * after column col are zero.
77 static void set_common_multiple(
78 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2,
79 unsigned row, unsigned col)
83 if (isl_int_eq(bmap1->eq[row][col], bmap2->eq[row][col]))
88 isl_int_lcm(m, bmap1->eq[row][col], bmap2->eq[row][col]);
89 isl_int_divexact(c, m, bmap1->eq[row][col]);
90 isl_seq_scale(bmap1->eq[row], bmap1->eq[row], c, col+1);
91 isl_int_divexact(c, m, bmap2->eq[row][col]);
92 isl_seq_scale(bmap2->eq[row], bmap2->eq[row], c, col+1);
97 /* Delete a given equality, moving all the following equalities one up.
99 static void delete_row(struct isl_basic_map *bmap, unsigned row)
106 for (r = row; r < bmap->n_eq; ++r)
107 bmap->eq[r] = bmap->eq[r+1];
108 bmap->eq[bmap->n_eq] = t;
111 /* Make first row entries in column col of bmap1 identical to
112 * those of bmap2, using the fact that entry bmap1->eq[row][col]=a
113 * is non-zero. Initially, these elements of bmap1 are all zero.
114 * For each row i < row, we set
115 * A[i] = a * A[i] + B[i][col] * A[row]
118 * A[i][col] = B[i][col] = a * old(B[i][col])
120 static void construct_column(
121 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2,
122 unsigned row, unsigned col)
131 total = 1 + bmap1->nparam + bmap1->n_in + bmap1->n_out + bmap1->n_div;
132 for (r = 0; r < row; ++r) {
133 if (isl_int_is_zero(bmap2->eq[r][col]))
135 isl_int_gcd(b, bmap2->eq[r][col], bmap1->eq[row][col]);
136 isl_int_divexact(a, bmap1->eq[row][col], b);
137 isl_int_divexact(b, bmap2->eq[r][col], b);
138 isl_seq_combine(bmap1->eq[r], a, bmap1->eq[r],
139 b, bmap1->eq[row], total);
140 isl_seq_scale(bmap2->eq[r], bmap2->eq[r], a, total);
144 delete_row(bmap1, row);
147 /* Make first row entries in column col of bmap1 identical to
148 * those of bmap2, using only these entries of the two matrices.
149 * Let t be the last row with different entries.
150 * For each row i < t, we set
151 * A[i] = (A[t][col]-B[t][col]) * A[i] + (B[i][col]-A[i][col) * A[t]
152 * B[i] = (A[t][col]-B[t][col]) * B[i] + (B[i][col]-A[i][col) * B[t]
154 * A[i][col] = B[i][col] = old(A[t][col]*B[i][col]-A[i][col]*B[t][col])
156 static int transform_column(
157 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2,
158 unsigned row, unsigned col)
164 for (t = row-1; t >= 0; --t)
165 if (isl_int_ne(bmap1->eq[t][col], bmap2->eq[t][col]))
170 total = 1 + bmap1->nparam + bmap1->n_in + bmap1->n_out + bmap1->n_div;
174 isl_int_sub(b, bmap1->eq[t][col], bmap2->eq[t][col]);
175 for (i = 0; i < t; ++i) {
176 isl_int_sub(a, bmap2->eq[i][col], bmap1->eq[i][col]);
177 isl_int_gcd(g, a, b);
178 isl_int_divexact(a, a, g);
179 isl_int_divexact(g, b, g);
180 isl_seq_combine(bmap1->eq[i], g, bmap1->eq[i], a, bmap1->eq[t],
182 isl_seq_combine(bmap2->eq[i], g, bmap2->eq[i], a, bmap2->eq[t],
188 delete_row(bmap1, t);
189 delete_row(bmap2, t);
193 /* The implementation is based on Section 5.2 of Michael Karr,
194 * "Affine Relationships Among Variables of a Program",
195 * except that the echelon form we use starts from the last column
196 * and that we are dealing with integer coefficients.
198 static struct isl_basic_map *affine_hull(struct isl_ctx *ctx,
199 struct isl_basic_map *bmap1, struct isl_basic_map *bmap2)
205 total = 1 + bmap1->nparam + bmap1->n_in + bmap1->n_out + bmap1->n_div;
208 for (col = total-1; col >= 0; --col) {
209 int is_zero1 = row >= bmap1->n_eq ||
210 isl_int_is_zero(bmap1->eq[row][col]);
211 int is_zero2 = row >= bmap2->n_eq ||
212 isl_int_is_zero(bmap2->eq[row][col]);
213 if (!is_zero1 && !is_zero2) {
214 set_common_multiple(bmap1, bmap2, row, col);
216 } else if (!is_zero1 && is_zero2) {
217 construct_column(bmap1, bmap2, row, col);
218 } else if (is_zero1 && !is_zero2) {
219 construct_column(bmap2, bmap1, row, col);
221 if (transform_column(bmap1, bmap2, row, col))
225 isl_basic_map_free(ctx, bmap2);
229 struct isl_basic_map *isl_map_affine_hull(struct isl_ctx *ctx,
233 struct isl_basic_map *bmap;
239 bmap = isl_basic_map_empty(ctx,
240 map->nparam, map->n_in, map->n_out);
241 isl_map_free(ctx, map);
245 map = isl_map_align_divs(ctx, map);
246 map = isl_map_cow(ctx, map);
248 for (i = 0; i < map->n; ++i) {
249 map->p[i] = isl_basic_map_cow(ctx, map->p[i]);
250 map->p[i] = isl_basic_map_affine_hull(ctx, map->p[i]);
251 map->p[i] = isl_basic_map_gauss(ctx, map->p[i], NULL);
256 map->p[0] = affine_hull(ctx, map->p[0], map->p[--map->n]);
260 bmap = isl_basic_map_copy(ctx, map->p[0]);
261 isl_map_free(ctx, map);
262 bmap = isl_basic_map_simplify(ctx, bmap);
263 return isl_basic_map_finalize(ctx, bmap);
265 isl_map_free(ctx, map);
269 struct isl_basic_set *isl_set_affine_hull(struct isl_ctx *ctx,
272 return (struct isl_basic_set *)
273 isl_map_affine_hull(ctx, (struct isl_map *)set);