3303d3b508a7e0076b78f3c3bdc2d30ad10db573
[platform/upstream/isl.git] / isl_affine_hull.c
1 #include "isl_ctx.h"
2 #include "isl_seq.h"
3 #include "isl_set.h"
4 #include "isl_lp.h"
5 #include "isl_map.h"
6 #include "isl_map_private.h"
7 #include "isl_equalities.h"
8 #include "isl_sample.h"
9 #include "isl_tab.h"
10
11 struct isl_basic_map *isl_basic_map_implicit_equalities(
12                                                 struct isl_basic_map *bmap)
13 {
14         struct isl_tab *tab;
15
16         if (!bmap)
17                 return bmap;
18
19         bmap = isl_basic_map_gauss(bmap, NULL);
20         if (ISL_F_ISSET(bmap, ISL_BASIC_MAP_EMPTY))
21                 return bmap;
22         if (ISL_F_ISSET(bmap, ISL_BASIC_MAP_NO_IMPLICIT))
23                 return bmap;
24         if (bmap->n_ineq <= 1)
25                 return bmap;
26
27         tab = isl_tab_from_basic_map(bmap);
28         tab = isl_tab_detect_implicit_equalities(tab);
29         bmap = isl_basic_map_update_from_tab(bmap, tab);
30         isl_tab_free(tab);
31         bmap = isl_basic_map_gauss(bmap, NULL);
32         ISL_F_SET(bmap, ISL_BASIC_MAP_NO_IMPLICIT);
33         return bmap;
34 }
35
36 struct isl_basic_set *isl_basic_set_implicit_equalities(
37                                                 struct isl_basic_set *bset)
38 {
39         return (struct isl_basic_set *)
40                 isl_basic_map_implicit_equalities((struct isl_basic_map*)bset);
41 }
42
43 struct isl_map *isl_map_implicit_equalities(struct isl_map *map)
44 {
45         int i;
46
47         if (!map)
48                 return map;
49
50         for (i = 0; i < map->n; ++i) {
51                 map->p[i] = isl_basic_map_implicit_equalities(map->p[i]);
52                 if (!map->p[i])
53                         goto error;
54         }
55
56         return map;
57 error:
58         isl_map_free(map);
59         return NULL;
60 }
61
62 /* Make eq[row][col] of both bmaps equal so we can add the row
63  * add the column to the common matrix.
64  * Note that because of the echelon form, the columns of row row
65  * after column col are zero.
66  */
67 static void set_common_multiple(
68         struct isl_basic_set *bset1, struct isl_basic_set *bset2,
69         unsigned row, unsigned col)
70 {
71         isl_int m, c;
72
73         if (isl_int_eq(bset1->eq[row][col], bset2->eq[row][col]))
74                 return;
75
76         isl_int_init(c);
77         isl_int_init(m);
78         isl_int_lcm(m, bset1->eq[row][col], bset2->eq[row][col]);
79         isl_int_divexact(c, m, bset1->eq[row][col]);
80         isl_seq_scale(bset1->eq[row], bset1->eq[row], c, col+1);
81         isl_int_divexact(c, m, bset2->eq[row][col]);
82         isl_seq_scale(bset2->eq[row], bset2->eq[row], c, col+1);
83         isl_int_clear(c);
84         isl_int_clear(m);
85 }
86
87 /* Delete a given equality, moving all the following equalities one up.
88  */
89 static void delete_row(struct isl_basic_set *bset, unsigned row)
90 {
91         isl_int *t;
92         int r;
93
94         t = bset->eq[row];
95         bset->n_eq--;
96         for (r = row; r < bset->n_eq; ++r)
97                 bset->eq[r] = bset->eq[r+1];
98         bset->eq[bset->n_eq] = t;
99 }
100
101 /* Make first row entries in column col of bset1 identical to
102  * those of bset2, using the fact that entry bset1->eq[row][col]=a
103  * is non-zero.  Initially, these elements of bset1 are all zero.
104  * For each row i < row, we set
105  *              A[i] = a * A[i] + B[i][col] * A[row]
106  *              B[i] = a * B[i]
107  * so that
108  *              A[i][col] = B[i][col] = a * old(B[i][col])
109  */
110 static void construct_column(
111         struct isl_basic_set *bset1, struct isl_basic_set *bset2,
112         unsigned row, unsigned col)
113 {
114         int r;
115         isl_int a;
116         isl_int b;
117         unsigned total;
118
119         isl_int_init(a);
120         isl_int_init(b);
121         total = 1 + isl_basic_set_n_dim(bset1);
122         for (r = 0; r < row; ++r) {
123                 if (isl_int_is_zero(bset2->eq[r][col]))
124                         continue;
125                 isl_int_gcd(b, bset2->eq[r][col], bset1->eq[row][col]);
126                 isl_int_divexact(a, bset1->eq[row][col], b);
127                 isl_int_divexact(b, bset2->eq[r][col], b);
128                 isl_seq_combine(bset1->eq[r], a, bset1->eq[r],
129                                               b, bset1->eq[row], total);
130                 isl_seq_scale(bset2->eq[r], bset2->eq[r], a, total);
131         }
132         isl_int_clear(a);
133         isl_int_clear(b);
134         delete_row(bset1, row);
135 }
136
137 /* Make first row entries in column col of bset1 identical to
138  * those of bset2, using only these entries of the two matrices.
139  * Let t be the last row with different entries.
140  * For each row i < t, we set
141  *      A[i] = (A[t][col]-B[t][col]) * A[i] + (B[i][col]-A[i][col) * A[t]
142  *      B[i] = (A[t][col]-B[t][col]) * B[i] + (B[i][col]-A[i][col) * B[t]
143  * so that
144  *      A[i][col] = B[i][col] = old(A[t][col]*B[i][col]-A[i][col]*B[t][col])
145  */
146 static int transform_column(
147         struct isl_basic_set *bset1, struct isl_basic_set *bset2,
148         unsigned row, unsigned col)
149 {
150         int i, t;
151         isl_int a, b, g;
152         unsigned total;
153
154         for (t = row-1; t >= 0; --t)
155                 if (isl_int_ne(bset1->eq[t][col], bset2->eq[t][col]))
156                         break;
157         if (t < 0)
158                 return 0;
159
160         total = 1 + isl_basic_set_n_dim(bset1);
161         isl_int_init(a);
162         isl_int_init(b);
163         isl_int_init(g);
164         isl_int_sub(b, bset1->eq[t][col], bset2->eq[t][col]);
165         for (i = 0; i < t; ++i) {
166                 isl_int_sub(a, bset2->eq[i][col], bset1->eq[i][col]);
167                 isl_int_gcd(g, a, b);
168                 isl_int_divexact(a, a, g);
169                 isl_int_divexact(g, b, g);
170                 isl_seq_combine(bset1->eq[i], g, bset1->eq[i], a, bset1->eq[t],
171                                 total);
172                 isl_seq_combine(bset2->eq[i], g, bset2->eq[i], a, bset2->eq[t],
173                                 total);
174         }
175         isl_int_clear(a);
176         isl_int_clear(b);
177         isl_int_clear(g);
178         delete_row(bset1, t);
179         delete_row(bset2, t);
180         return 1;
181 }
182
183 /* The implementation is based on Section 5.2 of Michael Karr,
184  * "Affine Relationships Among Variables of a Program",
185  * except that the echelon form we use starts from the last column
186  * and that we are dealing with integer coefficients.
187  */
188 static struct isl_basic_set *affine_hull(
189         struct isl_basic_set *bset1, struct isl_basic_set *bset2)
190 {
191         unsigned total;
192         int col;
193         int row;
194
195         total = 1 + isl_basic_set_n_dim(bset1);
196
197         row = 0;
198         for (col = total-1; col >= 0; --col) {
199                 int is_zero1 = row >= bset1->n_eq ||
200                         isl_int_is_zero(bset1->eq[row][col]);
201                 int is_zero2 = row >= bset2->n_eq ||
202                         isl_int_is_zero(bset2->eq[row][col]);
203                 if (!is_zero1 && !is_zero2) {
204                         set_common_multiple(bset1, bset2, row, col);
205                         ++row;
206                 } else if (!is_zero1 && is_zero2) {
207                         construct_column(bset1, bset2, row, col);
208                 } else if (is_zero1 && !is_zero2) {
209                         construct_column(bset2, bset1, row, col);
210                 } else {
211                         if (transform_column(bset1, bset2, row, col))
212                                 --row;
213                 }
214         }
215         isl_basic_set_free(bset2);
216         isl_assert(bset1->ctx, row == bset1->n_eq, goto error);
217         bset1 = isl_basic_set_normalize_constraints(bset1);
218         return bset1;
219 error:
220         isl_basic_set_free(bset1);
221         return NULL;
222 }
223
224 /* Find an integer point in "bset" that lies outside of the equality
225  * "eq" e(x) = 0.
226  * If "up" is true, look for a point satisfying e(x) - 1 >= 0.
227  * Otherwise, look for a point satisfying -e(x) - 1 >= 0 (i.e., e(x) <= -1).
228  * The point, if found, is returned as a singleton set.
229  * If no point can be found, the empty set is returned.
230  *
231  * Before solving an ILP problem, we first check if simply
232  * adding the normal of the constraint to one of the known
233  * integer points in the basic set yields another point
234  * inside the basic set.
235  *
236  * The caller of this function ensures that "bset" is bounded.
237  */
238 static struct isl_basic_set *outside_point(struct isl_ctx *ctx,
239         struct isl_basic_set *bset, isl_int *eq, int up)
240 {
241         struct isl_basic_set *slice = NULL;
242         struct isl_vec *sample;
243         struct isl_basic_set *point;
244         unsigned dim;
245         int k;
246
247         dim = isl_basic_set_n_dim(bset);
248         sample = isl_vec_alloc(ctx, 1 + dim);
249         if (!sample)
250                 return NULL;
251         isl_int_set_si(sample->block.data[0], 1);
252         isl_seq_combine(sample->block.data + 1,
253                 ctx->one, bset->sample->block.data + 1,
254                 up ? ctx->one : ctx->negone, eq + 1, dim);
255         if (isl_basic_set_contains(bset, sample))
256                 return isl_basic_set_from_vec(sample);
257         isl_vec_free(sample);
258         sample = NULL;
259
260         slice = isl_basic_set_copy(bset);
261         if (!slice)
262                 goto error;
263         slice = isl_basic_set_cow(slice);
264         slice = isl_basic_set_extend(slice, 0, dim, 0, 0, 1);
265         k = isl_basic_set_alloc_inequality(slice);
266         if (k < 0)
267                 goto error;
268         if (up)
269                 isl_seq_cpy(slice->ineq[k], eq, 1 + dim);
270         else
271                 isl_seq_neg(slice->ineq[k], eq, 1 + dim);
272         isl_int_sub_ui(slice->ineq[k][0], slice->ineq[k][0], 1);
273
274         sample = isl_basic_set_sample_bounded(slice);
275         if (!sample)
276                 goto error;
277         if (sample->size == 0) {
278                 isl_vec_free(sample);
279                 point = isl_basic_set_empty_like(bset);
280         } else
281                 point = isl_basic_set_from_vec(sample);
282
283         return point;
284 error:
285         isl_basic_set_free(slice);
286         return NULL;
287 }
288
289 struct isl_basic_set *isl_basic_set_recession_cone(struct isl_basic_set *bset)
290 {
291         int i;
292
293         bset = isl_basic_set_cow(bset);
294         if (!bset)
295                 return NULL;
296         isl_assert(bset->ctx, bset->n_div == 0, goto error);
297
298         for (i = 0; i < bset->n_eq; ++i)
299                 isl_int_set_si(bset->eq[i][0], 0);
300
301         for (i = 0; i < bset->n_ineq; ++i)
302                 isl_int_set_si(bset->ineq[i][0], 0);
303
304         ISL_F_CLR(bset, ISL_BASIC_SET_NO_IMPLICIT);
305         return isl_basic_set_implicit_equalities(bset);
306 error:
307         isl_basic_set_free(bset);
308         return NULL;
309 }
310
311 /* Extend an initial (under-)approximation of the affine hull of "bset"
312  * by looking for points that do not satisfy one of the equalities
313  * in the current approximation and adding them to that approximation
314  * until no such points can be found any more.
315  *
316  * The caller of this function ensures that "bset" is bounded.
317  */
318 static struct isl_basic_set *extend_affine_hull(struct isl_basic_set *bset,
319         struct isl_basic_set *hull)
320 {
321         int i, j, k;
322         struct isl_ctx *ctx;
323         unsigned dim;
324
325         ctx = bset->ctx;
326         dim = isl_basic_set_n_dim(bset);
327         for (i = 0; i < dim; ++i) {
328                 struct isl_basic_set *point;
329                 for (j = 0; j < hull->n_eq; ++j) {
330                         point = outside_point(ctx, bset, hull->eq[j], 1);
331                         if (!point)
332                                 goto error;
333                         if (!ISL_F_ISSET(point, ISL_BASIC_SET_EMPTY))
334                                 break;
335                         isl_basic_set_free(point);
336                         point = outside_point(ctx, bset, hull->eq[j], 0);
337                         if (!point)
338                                 goto error;
339                         if (!ISL_F_ISSET(point, ISL_BASIC_SET_EMPTY))
340                                 break;
341                         isl_basic_set_free(point);
342
343                         bset = isl_basic_set_extend_constraints(bset, 1, 0);
344                         k = isl_basic_set_alloc_equality(bset);
345                         if (k < 0)
346                                 goto error;
347                         isl_seq_cpy(bset->eq[k], hull->eq[j],
348                                         1 + isl_basic_set_total_dim(hull));
349                         bset = isl_basic_set_gauss(bset, NULL);
350                         if (!bset)
351                                 goto error;
352                 }
353                 if (j == hull->n_eq)
354                         break;
355                 hull = affine_hull(hull, point);
356         }
357         isl_basic_set_free(bset);
358
359         return hull;
360 error:
361         isl_basic_set_free(bset);
362         isl_basic_set_free(hull);
363         return NULL;
364 }
365
366 /* Drop all constraints in bset that involve any of the dimensions
367  * first to first+n-1.
368  */
369 static struct isl_basic_set *drop_constraints_involving
370         (struct isl_basic_set *bset, unsigned first, unsigned n)
371 {
372         int i;
373
374         if (!bset)
375                 return NULL;
376
377         bset = isl_basic_set_cow(bset);
378
379         for (i = bset->n_eq - 1; i >= 0; --i) {
380                 if (isl_seq_first_non_zero(bset->eq[i] + 1 + first, n) == -1)
381                         continue;
382                 isl_basic_set_drop_equality(bset, i);
383         }
384
385         for (i = bset->n_ineq - 1; i >= 0; --i) {
386                 if (isl_seq_first_non_zero(bset->ineq[i] + 1 + first, n) == -1)
387                         continue;
388                 isl_basic_set_drop_inequality(bset, i);
389         }
390
391         return bset;
392 }
393
394 /* Look for all equalities satisfied by the integer points in bset,
395  * which is assumed to be bounded.
396  *
397  * The equalities are obtained by successively looking for
398  * a point that is affinely independent of the points found so far.
399  * In particular, for each equality satisfied by the points so far,
400  * we check if there is any point on a hyperplane parallel to the
401  * corresponding hyperplane shifted by at least one (in either direction).
402  */
403 static struct isl_basic_set *uset_affine_hull_bounded(struct isl_basic_set *bset)
404 {
405         struct isl_vec *sample = NULL;
406         struct isl_basic_set *hull;
407
408         if (isl_basic_set_is_empty(bset))
409                 return bset;
410
411         sample = isl_basic_set_sample_vec(isl_basic_set_copy(bset));
412         if (!sample)
413                 goto error;
414         if (sample->size == 0) {
415                 struct isl_basic_set *hull;
416                 isl_vec_free(sample);
417                 hull = isl_basic_set_empty_like(bset);
418                 isl_basic_set_free(bset);
419                 return hull;
420         }
421         if (sample->size == 1) {
422                 isl_vec_free(sample);
423                 return bset;
424         }
425
426         hull = isl_basic_set_from_vec(sample);
427
428         return extend_affine_hull(bset, hull);
429 error:
430         isl_basic_set_free(bset);
431         return NULL;
432 }
433
434 /* Compute the affine hull of "bset", where "cone" is the recession cone
435  * of "bset".
436  *
437  * We first compute a unimodular transformation that puts the unbounded
438  * directions in the last dimensions.  In particular, we take a transformation
439  * that maps all equalities to equalities (in HNF) on the first dimensions.
440  * Let x be the original dimensions and y the transformed, with y_1 bounded
441  * and y_2 unbounded.
442  *
443  *             [ y_1 ]                  [ y_1 ]   [ Q_1 ]
444  *      x = U  [ y_2 ]                  [ y_2 ] = [ Q_2 ] x
445  *
446  * Let's call the input basic set S.  We compute S' = preimage(S, U)
447  * and drop the final dimensions including any constraints involving them.
448  * This results in set S''.
449  * Then we compute the affine hull A'' of S''.
450  * Let F y_1 >= g be the constraint system of A''.  In the transformed
451  * space the y_2 are unbounded, so we can add them back without any constraints,
452  * resulting in
453  *
454  *                      [ y_1 ]
455  *              [ F 0 ] [ y_2 ] >= g
456  * or
457  *                      [ Q_1 ]
458  *              [ F 0 ] [ Q_2 ] x >= g
459  * or
460  *              F Q_1 x >= g
461  *
462  * The affine hull in the original space is then obtained as
463  * A = preimage(A'', Q_1).
464  */
465 static struct isl_basic_set *affine_hull_with_cone(struct isl_basic_set *bset,
466         struct isl_basic_set *cone)
467 {
468         unsigned total;
469         unsigned cone_dim;
470         struct isl_basic_set *hull;
471         struct isl_mat *M, *U, *Q;
472
473         if (!bset || !cone)
474                 goto error;
475
476         total = isl_basic_set_total_dim(cone);
477         cone_dim = total - cone->n_eq;
478
479         M = isl_mat_sub_alloc(bset->ctx, cone->eq, 0, cone->n_eq, 1, total);
480         M = isl_mat_left_hermite(M, 0, &U, &Q);
481         if (!M)
482                 goto error;
483         isl_mat_free(M);
484
485         U = isl_mat_lin_to_aff(U);
486         bset = isl_basic_set_preimage(bset, isl_mat_copy(U));
487
488         bset = drop_constraints_involving(bset, total - cone_dim, cone_dim);
489         bset = isl_basic_set_drop_dims(bset, total - cone_dim, cone_dim);
490
491         Q = isl_mat_lin_to_aff(Q);
492         Q = isl_mat_drop_rows(Q, 1 + total - cone_dim, cone_dim);
493
494         if (bset && bset->sample && bset->sample->size == 1 + total)
495                 bset->sample = isl_mat_vec_product(isl_mat_copy(Q), bset->sample);
496
497         hull = uset_affine_hull_bounded(bset);
498
499         if (!hull)
500                 isl_mat_free(U);
501         else {
502                 struct isl_vec *sample = isl_vec_copy(hull->sample);
503                 U = isl_mat_drop_cols(U, 1 + total - cone_dim, cone_dim);
504                 if (sample && sample->size > 0)
505                         sample = isl_mat_vec_product(U, sample);
506                 else
507                         isl_mat_free(U);
508                 hull = isl_basic_set_preimage(hull, Q);
509                 isl_vec_free(hull->sample);
510                 hull->sample = sample;
511         }
512
513         isl_basic_set_free(cone);
514
515         return hull;
516 error:
517         isl_basic_set_free(bset);
518         isl_basic_set_free(cone);
519         return NULL;
520 }
521
522 /* Look for all equalities satisfied by the integer points in bset,
523  * which is assumed not to have any explicit equalities.
524  *
525  * The equalities are obtained by successively looking for
526  * a point that is affinely independent of the points found so far.
527  * In particular, for each equality satisfied by the points so far,
528  * we check if there is any point on a hyperplane parallel to the
529  * corresponding hyperplane shifted by at least one (in either direction).
530  *
531  * Before looking for any outside points, we first compute the recession
532  * cone.  The directions of this recession cone will always be part
533  * of the affine hull, so there is no need for looking for any points
534  * in these directions.
535  * In particular, if the recession cone is full-dimensional, then
536  * the affine hull is simply the whole universe.
537  */
538 static struct isl_basic_set *uset_affine_hull(struct isl_basic_set *bset)
539 {
540         struct isl_basic_set *cone;
541
542         if (isl_basic_set_fast_is_empty(bset))
543                 return bset;
544
545         cone = isl_basic_set_recession_cone(isl_basic_set_copy(bset));
546         if (!cone)
547                 goto error;
548         if (cone->n_eq == 0) {
549                 struct isl_basic_set *hull;
550                 isl_basic_set_free(cone);
551                 hull = isl_basic_set_universe_like(bset);
552                 isl_basic_set_free(bset);
553                 return hull;
554         }
555
556         if (cone->n_eq < isl_basic_set_total_dim(cone))
557                 return affine_hull_with_cone(bset, cone);
558
559         isl_basic_set_free(cone);
560         return uset_affine_hull_bounded(bset);
561 error:
562         isl_basic_set_free(bset);
563         return NULL;
564 }
565
566 /* Look for all equalities satisfied by the integer points in bmap
567  * that are independent of the equalities already explicitly available
568  * in bmap.
569  *
570  * We first remove all equalities already explicitly available,
571  * then look for additional equalities in the reduced space
572  * and then transform the result to the original space.
573  * The original equalities are _not_ added to this set.  This is
574  * the responsibility of the calling function.
575  * The resulting basic set has all meaning about the dimensions removed.
576  * In particular, dimensions that correspond to existential variables
577  * in bmap and that are found to be fixed are not removed.
578  */
579 static struct isl_basic_set *equalities_in_underlying_set(
580                                                 struct isl_basic_map *bmap)
581 {
582         struct isl_mat *T1 = NULL;
583         struct isl_mat *T2 = NULL;
584         struct isl_basic_set *bset = NULL;
585         struct isl_basic_set *hull = NULL;
586
587         bset = isl_basic_map_underlying_set(bmap);
588         if (!bset)
589                 return NULL;
590         if (bset->n_eq)
591                 bset = isl_basic_set_remove_equalities(bset, &T1, &T2);
592         if (!bset)
593                 goto error;
594
595         hull = uset_affine_hull(bset);
596         if (!T2)
597                 return hull;
598
599         if (!hull)
600                 isl_mat_free(T1);
601         else {
602                 struct isl_vec *sample = isl_vec_copy(hull->sample);
603                 if (sample && sample->size > 0)
604                         sample = isl_mat_vec_product(T1, sample);
605                 else
606                         isl_mat_free(T1);
607                 hull = isl_basic_set_preimage(hull, T2);
608                 isl_vec_free(hull->sample);
609                 hull->sample = sample;
610         }
611
612         return hull;
613 error:
614         isl_mat_free(T2);
615         isl_basic_set_free(bset);
616         isl_basic_set_free(hull);
617         return NULL;
618 }
619
620 /* Detect and make explicit all equalities satisfied by the (integer)
621  * points in bmap.
622  */
623 struct isl_basic_map *isl_basic_map_detect_equalities(
624                                                 struct isl_basic_map *bmap)
625 {
626         int i, j;
627         struct isl_basic_set *hull = NULL;
628
629         if (!bmap)
630                 return NULL;
631         if (bmap->n_ineq == 0)
632                 return bmap;
633         if (ISL_F_ISSET(bmap, ISL_BASIC_MAP_EMPTY))
634                 return bmap;
635         if (ISL_F_ISSET(bmap, ISL_BASIC_MAP_ALL_EQUALITIES))
636                 return bmap;
637         if (ISL_F_ISSET(bmap, ISL_BASIC_MAP_RATIONAL))
638                 return isl_basic_map_implicit_equalities(bmap);
639
640         hull = equalities_in_underlying_set(isl_basic_map_copy(bmap));
641         if (!hull)
642                 goto error;
643         if (ISL_F_ISSET(hull, ISL_BASIC_SET_EMPTY)) {
644                 isl_basic_set_free(hull);
645                 return isl_basic_map_set_to_empty(bmap);
646         }
647         bmap = isl_basic_map_extend_dim(bmap, isl_dim_copy(bmap->dim), 0,
648                                         hull->n_eq, 0);
649         for (i = 0; i < hull->n_eq; ++i) {
650                 j = isl_basic_map_alloc_equality(bmap);
651                 if (j < 0)
652                         goto error;
653                 isl_seq_cpy(bmap->eq[j], hull->eq[i],
654                                 1 + isl_basic_set_total_dim(hull));
655         }
656         isl_vec_free(bmap->sample);
657         bmap->sample = isl_vec_copy(hull->sample);
658         isl_basic_set_free(hull);
659         ISL_F_SET(bmap, ISL_BASIC_MAP_NO_IMPLICIT | ISL_BASIC_MAP_ALL_EQUALITIES);
660         bmap = isl_basic_map_simplify(bmap);
661         return isl_basic_map_finalize(bmap);
662 error:
663         isl_basic_set_free(hull);
664         isl_basic_map_free(bmap);
665         return NULL;
666 }
667
668 __isl_give isl_basic_set *isl_basic_set_detect_equalities(
669                                                 __isl_take isl_basic_set *bset)
670 {
671         return (isl_basic_set *)
672                 isl_basic_map_detect_equalities((isl_basic_map *)bset);
673 }
674
675 struct isl_map *isl_map_detect_equalities(struct isl_map *map)
676 {
677         struct isl_basic_map *bmap;
678         int i;
679
680         if (!map)
681                 return NULL;
682
683         for (i = 0; i < map->n; ++i) {
684                 bmap = isl_basic_map_copy(map->p[i]);
685                 bmap = isl_basic_map_detect_equalities(bmap);
686                 if (!bmap)
687                         goto error;
688                 isl_basic_map_free(map->p[i]);
689                 map->p[i] = bmap;
690         }
691
692         return map;
693 error:
694         isl_map_free(map);
695         return NULL;
696 }
697
698 __isl_give isl_set *isl_set_detect_equalities(__isl_take isl_set *set)
699 {
700         return (isl_set *)isl_map_detect_equalities((isl_map *)set);
701 }
702
703 /* After computing the rational affine hull (by detecting the implicit
704  * equalities), we compute the additional equalities satisfied by
705  * the integer points (if any) and add the original equalities back in.
706  */
707 struct isl_basic_map *isl_basic_map_affine_hull(struct isl_basic_map *bmap)
708 {
709         bmap = isl_basic_map_detect_equalities(bmap);
710         bmap = isl_basic_map_cow(bmap);
711         isl_basic_map_free_inequality(bmap, bmap->n_ineq);
712         return bmap;
713 }
714
715 struct isl_basic_set *isl_basic_set_affine_hull(struct isl_basic_set *bset)
716 {
717         return (struct isl_basic_set *)
718                 isl_basic_map_affine_hull((struct isl_basic_map *)bset);
719 }
720
721 struct isl_basic_map *isl_map_affine_hull(struct isl_map *map)
722 {
723         int i;
724         struct isl_basic_map *model = NULL;
725         struct isl_basic_map *hull = NULL;
726         struct isl_set *set;
727
728         if (!map)
729                 return NULL;
730
731         if (map->n == 0) {
732                 hull = isl_basic_map_empty_like_map(map);
733                 isl_map_free(map);
734                 return hull;
735         }
736
737         map = isl_map_detect_equalities(map);
738         map = isl_map_align_divs(map);
739         if (!map)
740                 return NULL;
741         model = isl_basic_map_copy(map->p[0]);
742         set = isl_map_underlying_set(map);
743         set = isl_set_cow(set);
744         if (!set)
745                 goto error;
746
747         for (i = 0; i < set->n; ++i) {
748                 set->p[i] = isl_basic_set_cow(set->p[i]);
749                 set->p[i] = isl_basic_set_affine_hull(set->p[i]);
750                 set->p[i] = isl_basic_set_gauss(set->p[i], NULL);
751                 if (!set->p[i])
752                         goto error;
753         }
754         set = isl_set_remove_empty_parts(set);
755         if (set->n == 0) {
756                 hull = isl_basic_map_empty_like(model);
757                 isl_basic_map_free(model);
758         } else {
759                 struct isl_basic_set *bset;
760                 while (set->n > 1) {
761                         set->p[0] = affine_hull(set->p[0], set->p[--set->n]);
762                         if (!set->p[0])
763                                 goto error;
764                 }
765                 bset = isl_basic_set_copy(set->p[0]);
766                 hull = isl_basic_map_overlying_set(bset, model);
767         }
768         isl_set_free(set);
769         hull = isl_basic_map_simplify(hull);
770         return isl_basic_map_finalize(hull);
771 error:
772         isl_basic_map_free(model);
773         isl_set_free(set);
774         return NULL;
775 }
776
777 struct isl_basic_set *isl_set_affine_hull(struct isl_set *set)
778 {
779         return (struct isl_basic_set *)
780                 isl_map_affine_hull((struct isl_map *)set);
781 }