isl_basic_map_detect_equalities: only compute recession cone once
[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, 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         hull = isl_basic_set_preimage(hull, Q);
500
501         isl_basic_set_free(cone);
502
503         return hull;
504 error:
505         isl_basic_set_free(bset);
506         isl_basic_set_free(cone);
507         return NULL;
508 }
509
510 /* Look for all equalities satisfied by the integer points in bset,
511  * which is assumed not to have any explicit equalities.
512  *
513  * The equalities are obtained by successively looking for
514  * a point that is affinely independent of the points found so far.
515  * In particular, for each equality satisfied by the points so far,
516  * we check if there is any point on a hyperplane parallel to the
517  * corresponding hyperplane shifted by at least one (in either direction).
518  *
519  * Before looking for any outside points, we first compute the recession
520  * cone.  The directions of this recession cone will always be part
521  * of the affine hull, so there is no need for looking for any points
522  * in these directions.
523  * In particular, if the recession cone is full-dimensional, then
524  * the affine hull is simply the whole universe.
525  */
526 static struct isl_basic_set *uset_affine_hull(struct isl_basic_set *bset)
527 {
528         struct isl_basic_set *cone;
529
530         if (isl_basic_set_fast_is_empty(bset))
531                 return bset;
532
533         cone = isl_basic_set_recession_cone(isl_basic_set_copy(bset));
534         if (!cone)
535                 goto error;
536         if (cone->n_eq == 0) {
537                 struct isl_basic_set *hull;
538                 isl_basic_set_free(cone);
539                 hull = isl_basic_set_universe_like(bset);
540                 isl_basic_set_free(bset);
541                 return hull;
542         }
543
544         if (cone->n_eq < isl_basic_set_total_dim(cone))
545                 return affine_hull_with_cone(bset, cone);
546
547         isl_basic_set_free(cone);
548         return uset_affine_hull_bounded(bset);
549 error:
550         isl_basic_set_free(bset);
551         return NULL;
552 }
553
554 /* Look for all equalities satisfied by the integer points in bmap
555  * that are independent of the equalities already explicitly available
556  * in bmap.
557  *
558  * We first remove all equalities already explicitly available,
559  * then look for additional equalities in the reduced space
560  * and then transform the result to the original space.
561  * The original equalities are _not_ added to this set.  This is
562  * the responsibility of the calling function.
563  * The resulting basic set has all meaning about the dimensions removed.
564  * In particular, dimensions that correspond to existential variables
565  * in bmap and that are found to be fixed are not removed.
566  */
567 static struct isl_basic_set *equalities_in_underlying_set(
568                                                 struct isl_basic_map *bmap)
569 {
570         struct isl_mat *T2 = NULL;
571         struct isl_basic_set *bset = NULL;
572         struct isl_basic_set *hull = NULL;
573
574         bset = isl_basic_map_underlying_set(bmap);
575         bset = isl_basic_set_remove_equalities(bset, NULL, &T2);
576         if (!bset)
577                 goto error;
578
579         hull = uset_affine_hull(bset);
580         if (T2)
581                 hull = isl_basic_set_preimage(hull, T2);
582
583         return hull;
584 error:
585         isl_mat_free(T2);
586         isl_basic_set_free(bset);
587         isl_basic_set_free(hull);
588         return NULL;
589 }
590
591 /* Detect and make explicit all equalities satisfied by the (integer)
592  * points in bmap.
593  */
594 struct isl_basic_map *isl_basic_map_detect_equalities(
595                                                 struct isl_basic_map *bmap)
596 {
597         int i, j;
598         struct isl_basic_set *hull = NULL;
599
600         if (!bmap)
601                 return NULL;
602         if (bmap->n_ineq == 0)
603                 return bmap;
604         if (ISL_F_ISSET(bmap, ISL_BASIC_MAP_EMPTY))
605                 return bmap;
606         if (ISL_F_ISSET(bmap, ISL_BASIC_MAP_ALL_EQUALITIES))
607                 return bmap;
608         if (ISL_F_ISSET(bmap, ISL_BASIC_MAP_RATIONAL))
609                 return isl_basic_map_implicit_equalities(bmap);
610
611         hull = equalities_in_underlying_set(isl_basic_map_copy(bmap));
612         if (!hull)
613                 goto error;
614         if (ISL_F_ISSET(hull, ISL_BASIC_SET_EMPTY)) {
615                 isl_basic_set_free(hull);
616                 return isl_basic_map_set_to_empty(bmap);
617         }
618         bmap = isl_basic_map_extend_dim(bmap, isl_dim_copy(bmap->dim), 0,
619                                         hull->n_eq, 0);
620         for (i = 0; i < hull->n_eq; ++i) {
621                 j = isl_basic_map_alloc_equality(bmap);
622                 if (j < 0)
623                         goto error;
624                 isl_seq_cpy(bmap->eq[j], hull->eq[i],
625                                 1 + isl_basic_set_total_dim(hull));
626         }
627         isl_basic_set_free(hull);
628         ISL_F_SET(bmap, ISL_BASIC_MAP_NO_IMPLICIT | ISL_BASIC_MAP_ALL_EQUALITIES);
629         bmap = isl_basic_map_simplify(bmap);
630         return isl_basic_map_finalize(bmap);
631 error:
632         isl_basic_set_free(hull);
633         isl_basic_map_free(bmap);
634         return NULL;
635 }
636
637 __isl_give isl_basic_set *isl_basic_set_detect_equalities(
638                                                 __isl_take isl_basic_set *bset)
639 {
640         return (isl_basic_set *)
641                 isl_basic_map_detect_equalities((isl_basic_map *)bset);
642 }
643
644 struct isl_map *isl_map_detect_equalities(struct isl_map *map)
645 {
646         struct isl_basic_map *bmap;
647         int i;
648
649         if (!map)
650                 return NULL;
651
652         for (i = 0; i < map->n; ++i) {
653                 bmap = isl_basic_map_copy(map->p[i]);
654                 bmap = isl_basic_map_detect_equalities(bmap);
655                 if (!bmap)
656                         goto error;
657                 isl_basic_map_free(map->p[i]);
658                 map->p[i] = bmap;
659         }
660
661         return map;
662 error:
663         isl_map_free(map);
664         return NULL;
665 }
666
667 __isl_give isl_set *isl_set_detect_equalities(__isl_take isl_set *set)
668 {
669         return (isl_set *)isl_map_detect_equalities((isl_map *)set);
670 }
671
672 /* After computing the rational affine hull (by detecting the implicit
673  * equalities), we compute the additional equalities satisfied by
674  * the integer points (if any) and add the original equalities back in.
675  */
676 struct isl_basic_map *isl_basic_map_affine_hull(struct isl_basic_map *bmap)
677 {
678         bmap = isl_basic_map_detect_equalities(bmap);
679         bmap = isl_basic_map_cow(bmap);
680         isl_basic_map_free_inequality(bmap, bmap->n_ineq);
681         return bmap;
682 }
683
684 struct isl_basic_set *isl_basic_set_affine_hull(struct isl_basic_set *bset)
685 {
686         return (struct isl_basic_set *)
687                 isl_basic_map_affine_hull((struct isl_basic_map *)bset);
688 }
689
690 struct isl_basic_map *isl_map_affine_hull(struct isl_map *map)
691 {
692         int i;
693         struct isl_basic_map *model = NULL;
694         struct isl_basic_map *hull = NULL;
695         struct isl_set *set;
696
697         if (!map)
698                 return NULL;
699
700         if (map->n == 0) {
701                 hull = isl_basic_map_empty_like_map(map);
702                 isl_map_free(map);
703                 return hull;
704         }
705
706         map = isl_map_detect_equalities(map);
707         map = isl_map_align_divs(map);
708         if (!map)
709                 return NULL;
710         model = isl_basic_map_copy(map->p[0]);
711         set = isl_map_underlying_set(map);
712         set = isl_set_cow(set);
713         if (!set)
714                 goto error;
715
716         for (i = 0; i < set->n; ++i) {
717                 set->p[i] = isl_basic_set_cow(set->p[i]);
718                 set->p[i] = isl_basic_set_affine_hull(set->p[i]);
719                 set->p[i] = isl_basic_set_gauss(set->p[i], NULL);
720                 if (!set->p[i])
721                         goto error;
722         }
723         set = isl_set_remove_empty_parts(set);
724         if (set->n == 0) {
725                 hull = isl_basic_map_empty_like(model);
726                 isl_basic_map_free(model);
727         } else {
728                 struct isl_basic_set *bset;
729                 while (set->n > 1) {
730                         set->p[0] = affine_hull(set->p[0], set->p[--set->n]);
731                         if (!set->p[0])
732                                 goto error;
733                 }
734                 bset = isl_basic_set_copy(set->p[0]);
735                 hull = isl_basic_map_overlying_set(bset, model);
736         }
737         isl_set_free(set);
738         hull = isl_basic_map_simplify(hull);
739         return isl_basic_map_finalize(hull);
740 error:
741         isl_basic_map_free(model);
742         isl_set_free(set);
743         return NULL;
744 }
745
746 struct isl_basic_set *isl_set_affine_hull(struct isl_set *set)
747 {
748         return (struct isl_basic_set *)
749                 isl_map_affine_hull((struct isl_map *)set);
750 }