Coverage Report

Created: 2017-08-21 19:50

/Users/buildslave/jenkins/sharedspace/clang-stage2-coverage-R@2/llvm/tools/polly/lib/External/isl/isl_convex_hull.c
Line
Count
Source (jump to first uncovered line)
1
/*
2
 * Copyright 2008-2009 Katholieke Universiteit Leuven
3
 * Copyright 2014      INRIA Rocquencourt
4
 *
5
 * Use of this software is governed by the MIT license
6
 *
7
 * Written by Sven Verdoolaege, K.U.Leuven, Departement
8
 * Computerwetenschappen, Celestijnenlaan 200A, B-3001 Leuven, Belgium
9
 * and Inria Paris - Rocquencourt, Domaine de Voluceau - Rocquencourt,
10
 * B.P. 105 - 78153 Le Chesnay, France
11
 */
12
13
#include <isl_ctx_private.h>
14
#include <isl_map_private.h>
15
#include <isl_lp_private.h>
16
#include <isl/map.h>
17
#include <isl_mat_private.h>
18
#include <isl_vec_private.h>
19
#include <isl/set.h>
20
#include <isl_seq.h>
21
#include <isl_options_private.h>
22
#include "isl_equalities.h"
23
#include "isl_tab.h"
24
#include <isl_sort.h>
25
26
#include <bset_to_bmap.c>
27
#include <bset_from_bmap.c>
28
#include <set_to_map.c>
29
30
static __isl_give isl_basic_set *uset_convex_hull_wrap_bounded(
31
  __isl_take isl_set *set);
32
33
/* Remove redundant
34
 * constraints.  If the minimal value along the normal of a constraint
35
 * is the same if the constraint is removed, then the constraint is redundant.
36
 *
37
 * Since some constraints may be mutually redundant, sort the constraints
38
 * first such that constraints that involve existentially quantified
39
 * variables are considered for removal before those that do not.
40
 * The sorting is also needed for the use in map_simple_hull.
41
 *
42
 * Note that isl_tab_detect_implicit_equalities may also end up
43
 * marking some constraints as redundant.  Make sure the constraints
44
 * are preserved and undo those marking such that isl_tab_detect_redundant
45
 * can consider the constraints in the sorted order.
46
 *
47
 * Alternatively, we could have intersected the basic map with the
48
 * corresponding equality and then checked if the dimension was that
49
 * of a facet.
50
 */
51
__isl_give isl_basic_map *isl_basic_map_remove_redundancies(
52
  __isl_take isl_basic_map *bmap)
53
276k
{
54
276k
  struct isl_tab *tab;
55
276k
56
276k
  if (!bmap)
57
0
    return NULL;
58
276k
59
276k
  bmap = isl_basic_map_gauss(bmap, NULL);
60
276k
  if (ISL_F_ISSET(bmap, ISL_BASIC_MAP_EMPTY))
61
1.10k
    return bmap;
62
275k
  
if (275k
ISL_F_ISSET275k
(bmap, ISL_BASIC_MAP_NO_REDUNDANT))
63
88.7k
    return bmap;
64
187k
  
if (187k
bmap->n_ineq <= 1187k
)
65
91.1k
    return bmap;
66
187k
67
187k
  bmap = isl_basic_map_sort_constraints(bmap);
68
95.9k
  tab = isl_tab_from_basic_map(bmap, 0);
69
95.9k
  if (!tab)
70
0
    goto error;
71
95.9k
  tab->preserve = 1;
72
95.9k
  if (isl_tab_detect_implicit_equalities(tab) < 0)
73
0
    goto error;
74
95.9k
  
if (95.9k
isl_tab_restore_redundant(tab) < 095.9k
)
75
0
    goto error;
76
95.9k
  tab->preserve = 0;
77
95.9k
  if (isl_tab_detect_redundant(tab) < 0)
78
0
    goto error;
79
95.9k
  bmap = isl_basic_map_update_from_tab(bmap, tab);
80
95.9k
  isl_tab_free(tab);
81
95.9k
  if (!bmap)
82
0
    return NULL;
83
95.9k
  
ISL_F_SET95.9k
(bmap, ISL_BASIC_MAP_NO_IMPLICIT);95.9k
84
95.9k
  ISL_F_SET(bmap, ISL_BASIC_MAP_NO_REDUNDANT);
85
95.9k
  return bmap;
86
95.9k
error:
87
0
  isl_tab_free(tab);
88
0
  isl_basic_map_free(bmap);
89
95.9k
  return NULL;
90
276k
}
91
92
__isl_give isl_basic_set *isl_basic_set_remove_redundancies(
93
  __isl_take isl_basic_set *bset)
94
2.92k
{
95
2.92k
  return bset_from_bmap(
96
2.92k
    isl_basic_map_remove_redundancies(bset_to_bmap(bset)));
97
2.92k
}
98
99
/* Remove redundant constraints in each of the basic maps.
100
 */
101
__isl_give isl_map *isl_map_remove_redundancies(__isl_take isl_map *map)
102
863
{
103
863
  return isl_map_inline_foreach_basic_map(map,
104
863
              &isl_basic_map_remove_redundancies);
105
863
}
106
107
__isl_give isl_set *isl_set_remove_redundancies(__isl_take isl_set *set)
108
852
{
109
852
  return isl_map_remove_redundancies(set);
110
852
}
111
112
/* Check if the set set is bound in the direction of the affine
113
 * constraint c and if so, set the constant term such that the
114
 * resulting constraint is a bounding constraint for the set.
115
 */
116
static int uset_is_bound(__isl_keep isl_set *set, isl_int *c, unsigned len)
117
3.88k
{
118
3.88k
  int first;
119
3.88k
  int j;
120
3.88k
  isl_int opt;
121
3.88k
  isl_int opt_denom;
122
3.88k
123
3.88k
  isl_int_init(opt);
124
3.88k
  isl_int_init(opt_denom);
125
3.88k
  first = 1;
126
11.6k
  for (j = 0; 
j < set->n11.6k
;
++j7.76k
)
{7.76k
127
7.76k
    enum isl_lp_result res;
128
7.76k
129
7.76k
    if (ISL_F_ISSET(set->p[j], ISL_BASIC_SET_EMPTY))
130
0
      continue;
131
7.76k
132
7.76k
    res = isl_basic_set_solve_lp(set->p[j],
133
7.76k
        0, c, set->ctx->one, &opt, &opt_denom, NULL);
134
7.76k
    if (res == isl_lp_unbounded)
135
0
      break;
136
7.76k
    
if (7.76k
res == isl_lp_error7.76k
)
137
0
      goto error;
138
7.76k
    
if (7.76k
res == isl_lp_empty7.76k
)
{0
139
0
      set->p[j] = isl_basic_set_set_to_empty(set->p[j]);
140
0
      if (!set->p[j])
141
0
        goto error;
142
0
      continue;
143
7.76k
    }
144
7.76k
    
if (7.76k
first || 7.76k
isl_int_is_neg3.88k
(opt))
{4.75k
145
4.75k
      if (
!4.75k
isl_int_is_one4.75k
(opt_denom))
146
4.54k
        isl_seq_scale(c, c, opt_denom, len);
147
4.75k
      isl_int_sub(c[0], c[0], opt);
148
7.76k
    }
149
7.76k
    first = 0;
150
7.76k
  }
151
3.88k
  
isl_int_clear3.88k
(opt);3.88k
152
3.88k
  isl_int_clear(opt_denom);
153
3.88k
  return j >= set->n;
154
3.88k
error:
155
0
  isl_int_clear(opt);
156
0
  isl_int_clear(opt_denom);
157
3.88k
  return -1;
158
3.88k
}
159
160
static struct isl_basic_set *isl_basic_set_add_equality(
161
  struct isl_basic_set *bset, isl_int *c)
162
42.3k
{
163
42.3k
  int i;
164
42.3k
  unsigned dim;
165
42.3k
166
42.3k
  if (!bset)
167
0
    return NULL;
168
42.3k
169
42.3k
  
if (42.3k
ISL_F_ISSET42.3k
(bset, ISL_BASIC_SET_EMPTY))
170
0
    return bset;
171
42.3k
172
42.3k
  
isl_assert42.3k
(bset->ctx, isl_basic_set_n_param(bset) == 0, goto error);42.3k
173
42.3k
  
isl_assert42.3k
(bset->ctx, bset->n_div == 0, goto error);42.3k
174
42.3k
  dim = isl_basic_set_n_dim(bset);
175
42.3k
  bset = isl_basic_set_cow(bset);
176
42.3k
  bset = isl_basic_set_extend(bset, 0, dim, 0, 1, 0);
177
42.3k
  i = isl_basic_set_alloc_equality(bset);
178
42.3k
  if (i < 0)
179
0
    goto error;
180
42.3k
  isl_seq_cpy(bset->eq[i], c, 1 + dim);
181
42.3k
  return bset;
182
42.3k
error:
183
0
  isl_basic_set_free(bset);
184
42.3k
  return NULL;
185
42.3k
}
186
187
static __isl_give isl_set *isl_set_add_basic_set_equality(
188
  __isl_take isl_set *set, isl_int *c)
189
7.59k
{
190
7.59k
  int i;
191
7.59k
192
7.59k
  set = isl_set_cow(set);
193
7.59k
  if (!set)
194
0
    return NULL;
195
22.7k
  
for (i = 0; 7.59k
i < set->n22.7k
;
++i15.1k
)
{15.1k
196
15.1k
    set->p[i] = isl_basic_set_add_equality(set->p[i], c);
197
15.1k
    if (!set->p[i])
198
0
      goto error;
199
15.1k
  }
200
7.59k
  return set;
201
7.59k
error:
202
0
  isl_set_free(set);
203
7.59k
  return NULL;
204
7.59k
}
205
206
/* Given a union of basic sets, construct the constraints for wrapping
207
 * a facet around one of its ridges.
208
 * In particular, if each of n the d-dimensional basic sets i in "set"
209
 * contains the origin, satisfies the constraints x_1 >= 0 and x_2 >= 0
210
 * and is defined by the constraints
211
 *            [ 1 ]
212
 *        A_i [ x ]  >= 0
213
 *
214
 * then the resulting set is of dimension n*(1+d) and has as constraints
215
 *
216
 *            [ a_i ]
217
 *        A_i [ x_i ] >= 0
218
 *
219
 *              a_i   >= 0
220
 *
221
 *      \sum_i x_{i,1} = 1
222
 */
223
static __isl_give isl_basic_set *wrap_constraints(__isl_keep isl_set *set)
224
15.8k
{
225
15.8k
  struct isl_basic_set *lp;
226
15.8k
  unsigned n_eq;
227
15.8k
  unsigned n_ineq;
228
15.8k
  int i, j, k;
229
15.8k
  unsigned dim, lp_dim;
230
15.8k
231
15.8k
  if (!set)
232
0
    return NULL;
233
15.8k
234
15.8k
  dim = 1 + isl_set_n_dim(set);
235
15.8k
  n_eq = 1;
236
15.8k
  n_ineq = set->n;
237
44.9k
  for (i = 0; 
i < set->n44.9k
;
++i29.1k
)
{29.1k
238
29.1k
    n_eq += set->p[i]->n_eq;
239
29.1k
    n_ineq += set->p[i]->n_ineq;
240
29.1k
  }
241
15.8k
  lp = isl_basic_set_alloc(set->ctx, 0, dim * set->n, 0, n_eq, n_ineq);
242
15.8k
  lp = isl_basic_set_set_rational(lp);
243
15.8k
  if (!lp)
244
0
    return NULL;
245
15.8k
  lp_dim = isl_basic_set_n_dim(lp);
246
15.8k
  k = isl_basic_set_alloc_equality(lp);
247
15.8k
  isl_int_set_si(lp->eq[k][0], -1);
248
44.9k
  for (i = 0; 
i < set->n44.9k
;
++i29.1k
)
{29.1k
249
29.1k
    isl_int_set_si(lp->eq[k][1+dim*i], 0);
250
29.1k
    isl_int_set_si(lp->eq[k][1+dim*i+1], 1);
251
29.1k
    isl_seq_clr(lp->eq[k]+1+dim*i+2, dim-2);
252
29.1k
  }
253
44.9k
  for (i = 0; 
i < set->n44.9k
;
++i29.1k
)
{29.1k
254
29.1k
    k = isl_basic_set_alloc_inequality(lp);
255
29.1k
    isl_seq_clr(lp->ineq[k], 1+lp_dim);
256
29.1k
    isl_int_set_si(lp->ineq[k][1+dim*i], 1);
257
29.1k
258
65.9k
    for (j = 0; 
j < set->p[i]->n_eq65.9k
;
++j36.8k
)
{36.8k
259
36.8k
      k = isl_basic_set_alloc_equality(lp);
260
36.8k
      isl_seq_clr(lp->eq[k], 1+dim*i);
261
36.8k
      isl_seq_cpy(lp->eq[k]+1+dim*i, set->p[i]->eq[j], dim);
262
36.8k
      isl_seq_clr(lp->eq[k]+1+dim*(i+1), dim*(set->n-i-1));
263
36.8k
    }
264
29.1k
265
98.3k
    for (j = 0; 
j < set->p[i]->n_ineq98.3k
;
++j69.1k
)
{69.1k
266
69.1k
      k = isl_basic_set_alloc_inequality(lp);
267
69.1k
      isl_seq_clr(lp->ineq[k], 1+dim*i);
268
69.1k
      isl_seq_cpy(lp->ineq[k]+1+dim*i, set->p[i]->ineq[j], dim);
269
69.1k
      isl_seq_clr(lp->ineq[k]+1+dim*(i+1), dim*(set->n-i-1));
270
69.1k
    }
271
29.1k
  }
272
15.8k
  return lp;
273
15.8k
}
274
275
/* Given a facet "facet" of the convex hull of "set" and a facet "ridge"
276
 * of that facet, compute the other facet of the convex hull that contains
277
 * the ridge.
278
 *
279
 * We first transform the set such that the facet constraint becomes
280
 *
281
 *      x_1 >= 0
282
 *
283
 * I.e., the facet lies in
284
 *
285
 *      x_1 = 0
286
 *
287
 * and on that facet, the constraint that defines the ridge is
288
 *
289
 *      x_2 >= 0
290
 *
291
 * (This transformation is not strictly needed, all that is needed is
292
 * that the ridge contains the origin.)
293
 *
294
 * Since the ridge contains the origin, the cone of the convex hull
295
 * will be of the form
296
 *
297
 *      x_1 >= 0
298
 *      x_2 >= a x_1
299
 *
300
 * with this second constraint defining the new facet.
301
 * The constant a is obtained by settting x_1 in the cone of the
302
 * convex hull to 1 and minimizing x_2.
303
 * Now, each element in the cone of the convex hull is the sum
304
 * of elements in the cones of the basic sets.
305
 * If a_i is the dilation factor of basic set i, then the problem
306
 * we need to solve is
307
 *
308
 *      min \sum_i x_{i,2}
309
 *      st
310
 *        \sum_i x_{i,1} = 1
311
 *            a_i   >= 0
312
 *          [ a_i ]
313
 *        A [ x_i ] >= 0
314
 *
315
 * with
316
 *            [  1  ]
317
 *        A_i [ x_i ] >= 0
318
 *
319
 * the constraints of each (transformed) basic set.
320
 * If a = n/d, then the constraint defining the new facet (in the transformed
321
 * space) is
322
 *
323
 *      -n x_1 + d x_2 >= 0
324
 *
325
 * In the original space, we need to take the same combination of the
326
 * corresponding constraints "facet" and "ridge".
327
 *
328
 * If a = -infty = "-1/0", then we just return the original facet constraint.
329
 * This means that the facet is unbounded, but has a bounded intersection
330
 * with the union of sets.
331
 */
332
isl_int *isl_set_wrap_facet(__isl_keep isl_set *set,
333
  isl_int *facet, isl_int *ridge)
334
15.8k
{
335
15.8k
  int i;
336
15.8k
  isl_ctx *ctx;
337
15.8k
  struct isl_mat *T = NULL;
338
15.8k
  struct isl_basic_set *lp = NULL;
339
15.8k
  struct isl_vec *obj;
340
15.8k
  enum isl_lp_result res;
341
15.8k
  isl_int num, den;
342
15.8k
  unsigned dim;
343
15.8k
344
15.8k
  if (!set)
345
0
    return NULL;
346
15.8k
  ctx = set->ctx;
347
15.8k
  set = isl_set_copy(set);
348
15.8k
  set = isl_set_set_rational(set);
349
15.8k
350
15.8k
  dim = 1 + isl_set_n_dim(set);
351
15.8k
  T = isl_mat_alloc(ctx, 3, dim);
352
15.8k
  if (!T)
353
0
    goto error;
354
15.8k
  
isl_int_set_si15.8k
(T->row[0][0], 1);15.8k
355
15.8k
  isl_seq_clr(T->row[0]+1, dim - 1);
356
15.8k
  isl_seq_cpy(T->row[1], facet, dim);
357
15.8k
  isl_seq_cpy(T->row[2], ridge, dim);
358
15.8k
  T = isl_mat_right_inverse(T);
359
15.8k
  set = isl_set_preimage(set, T);
360
15.8k
  T = NULL;
361
15.8k
  if (!set)
362
0
    goto error;
363
15.8k
  lp = wrap_constraints(set);
364
15.8k
  obj = isl_vec_alloc(ctx, 1 + dim*set->n);
365
15.8k
  if (!obj)
366
0
    goto error;
367
15.8k
  
isl_int_set_si15.8k
(obj->block.data[0], 0);15.8k
368
44.9k
  for (i = 0; 
i < set->n44.9k
;
++i29.1k
)
{29.1k
369
29.1k
    isl_seq_clr(obj->block.data + 1 + dim*i, 2);
370
29.1k
    isl_int_set_si(obj->block.data[1 + dim*i+2], 1);
371
29.1k
    isl_seq_clr(obj->block.data + 1 + dim*i+3, dim-3);
372
29.1k
  }
373
15.8k
  isl_int_init(num);
374
15.8k
  isl_int_init(den);
375
15.8k
  res = isl_basic_set_solve_lp(lp, 0,
376
15.8k
          obj->block.data, ctx->one, &num, &den, NULL);
377
15.8k
  if (
res == isl_lp_ok15.8k
)
{14.5k
378
14.5k
    isl_int_neg(num, num);
379
14.5k
    isl_seq_combine(facet, num, facet, den, ridge, dim);
380
14.5k
    isl_seq_normalize(ctx, facet, dim);
381
15.8k
  }
382
15.8k
  isl_int_clear(num);
383
15.8k
  isl_int_clear(den);
384
15.8k
  isl_vec_free(obj);
385
15.8k
  isl_basic_set_free(lp);
386
15.8k
  isl_set_free(set);
387
15.8k
  if (res == isl_lp_error)
388
0
    return NULL;
389
15.8k
  
isl_assert15.8k
(ctx, res == isl_lp_ok || res == isl_lp_unbounded, 15.8k
390
15.8k
       return NULL);
391
15.8k
  return facet;
392
15.8k
error:
393
0
  isl_basic_set_free(lp);
394
0
  isl_mat_free(T);
395
0
  isl_set_free(set);
396
15.8k
  return NULL;
397
15.8k
}
398
399
/* Compute the constraint of a facet of "set".
400
 *
401
 * We first compute the intersection with a bounding constraint
402
 * that is orthogonal to one of the coordinate axes.
403
 * If the affine hull of this intersection has only one equality,
404
 * we have found a facet.
405
 * Otherwise, we wrap the current bounding constraint around
406
 * one of the equalities of the face (one that is not equal to
407
 * the current bounding constraint).
408
 * This process continues until we have found a facet.
409
 * The dimension of the intersection increases by at least
410
 * one on each iteration, so termination is guaranteed.
411
 */
412
static __isl_give isl_mat *initial_facet_constraint(__isl_keep isl_set *set)
413
3.88k
{
414
3.88k
  struct isl_set *slice = NULL;
415
3.88k
  struct isl_basic_set *face = NULL;
416
3.88k
  int i;
417
3.88k
  unsigned dim = isl_set_n_dim(set);
418
3.88k
  int is_bound;
419
3.88k
  isl_mat *bounds = NULL;
420
3.88k
421
3.88k
  isl_assert(set->ctx, set->n > 0, goto error);
422
3.88k
  bounds = isl_mat_alloc(set->ctx, 1, 1 + dim);
423
3.88k
  if (!bounds)
424
0
    return NULL;
425
3.88k
426
3.88k
  isl_seq_clr(bounds->row[0], dim);
427
3.88k
  isl_int_set_si(bounds->row[0][1 + dim - 1], 1);
428
3.88k
  is_bound = uset_is_bound(set, bounds->row[0], 1 + dim);
429
3.88k
  if (is_bound < 0)
430
0
    goto error;
431
3.88k
  
isl_assert3.88k
(set->ctx, is_bound, goto error);3.88k
432
3.88k
  isl_seq_normalize(set->ctx, bounds->row[0], 1 + dim);
433
3.88k
  bounds->n_row = 1;
434
3.88k
435
7.59k
  for (;;) {
436
7.59k
    slice = isl_set_copy(set);
437
7.59k
    slice = isl_set_add_basic_set_equality(slice, bounds->row[0]);
438
7.59k
    face = isl_set_affine_hull(slice);
439
7.59k
    if (!face)
440
0
      goto error;
441
7.59k
    
if (7.59k
face->n_eq == 17.59k
)
{3.88k
442
3.88k
      isl_basic_set_free(face);
443
3.88k
      break;
444
7.59k
    }
445
6.98k
    
for (i = 0; 3.71k
i < face->n_eq6.98k
;
++i3.27k
)
446
6.98k
      
if (6.98k
!isl_seq_eq(bounds->row[0], face->eq[i], 1 + dim) &&6.98k
447
6.98k
          !isl_seq_is_neg(bounds->row[0],
448
6.98k
            face->eq[i], 1 + dim))
449
3.71k
        break;
450
3.71k
    isl_assert(set->ctx, i < face->n_eq, goto error);
451
3.71k
    
if (3.71k
!isl_set_wrap_facet(set, bounds->row[0], face->eq[i])3.71k
)
452
0
      goto error;
453
3.71k
    isl_seq_normalize(set->ctx, bounds->row[0], bounds->n_col);
454
3.71k
    isl_basic_set_free(face);
455
3.88k
  }
456
3.88k
457
3.88k
  return bounds;
458
3.88k
error:
459
0
  isl_basic_set_free(face);
460
0
  isl_mat_free(bounds);
461
3.88k
  return NULL;
462
3.88k
}
463
464
/* Given the bounding constraint "c" of a facet of the convex hull of "set",
465
 * compute a hyperplane description of the facet, i.e., compute the facets
466
 * of the facet.
467
 *
468
 * We compute an affine transformation that transforms the constraint
469
 *
470
 *        [ 1 ]
471
 *      c [ x ] = 0
472
 *
473
 * to the constraint
474
 *
475
 *         z_1  = 0
476
 *
477
 * by computing the right inverse U of a matrix that starts with the rows
478
 *
479
 *      [ 1 0 ]
480
 *      [  c  ]
481
 *
482
 * Then
483
 *      [ 1 ]     [ 1 ]
484
 *      [ x ] = U [ z ]
485
 * and
486
 *      [ 1 ]     [ 1 ]
487
 *      [ z ] = Q [ x ]
488
 *
489
 * with Q = U^{-1}
490
 * Since z_1 is zero, we can drop this variable as well as the corresponding
491
 * column of U to obtain
492
 *
493
 *      [ 1 ]      [ 1  ]
494
 *      [ x ] = U' [ z' ]
495
 * and
496
 *      [ 1  ]      [ 1 ]
497
 *      [ z' ] = Q' [ x ]
498
 *
499
 * with Q' equal to Q, but without the corresponding row.
500
 * After computing the facets of the facet in the z' space,
501
 * we convert them back to the x space through Q.
502
 */
503
static __isl_give isl_basic_set *compute_facet(__isl_keep isl_set *set,
504
  isl_int *c)
505
13.5k
{
506
13.5k
  struct isl_mat *m, *U, *Q;
507
13.5k
  struct isl_basic_set *facet = NULL;
508
13.5k
  struct isl_ctx *ctx;
509
13.5k
  unsigned dim;
510
13.5k
511
13.5k
  ctx = set->ctx;
512
13.5k
  set = isl_set_copy(set);
513
13.5k
  dim = isl_set_n_dim(set);
514
13.5k
  m = isl_mat_alloc(set->ctx, 2, 1 + dim);
515
13.5k
  if (!m)
516
0
    goto error;
517
13.5k
  
isl_int_set_si13.5k
(m->row[0][0], 1);13.5k
518
13.5k
  isl_seq_clr(m->row[0]+1, dim);
519
13.5k
  isl_seq_cpy(m->row[1], c, 1+dim);
520
13.5k
  U = isl_mat_right_inverse(m);
521
13.5k
  Q = isl_mat_right_inverse(isl_mat_copy(U));
522
13.5k
  U = isl_mat_drop_cols(U, 1, 1);
523
13.5k
  Q = isl_mat_drop_rows(Q, 1, 1);
524
13.5k
  set = isl_set_preimage(set, U);
525
13.5k
  facet = uset_convex_hull_wrap_bounded(set);
526
13.5k
  facet = isl_basic_set_preimage(facet, Q);
527
13.5k
  if (
facet && 13.5k
facet->n_eq != 013.5k
)
528
0
    isl_die(ctx, isl_error_internal, "unexpected equality",
529
13.5k
      return isl_basic_set_free(facet));
530
13.5k
  return facet;
531
13.5k
error:
532
0
  isl_basic_set_free(facet);
533
0
  isl_set_free(set);
534
13.5k
  return NULL;
535
13.5k
}
536
537
/* Given an initial facet constraint, compute the remaining facets.
538
 * We do this by running through all facets found so far and computing
539
 * the adjacent facets through wrapping, adding those facets that we
540
 * hadn't already found before.
541
 *
542
 * For each facet we have found so far, we first compute its facets
543
 * in the resulting convex hull.  That is, we compute the ridges
544
 * of the resulting convex hull contained in the facet.
545
 * We also compute the corresponding facet in the current approximation
546
 * of the convex hull.  There is no need to wrap around the ridges
547
 * in this facet since that would result in a facet that is already
548
 * present in the current approximation.
549
 *
550
 * This function can still be significantly optimized by checking which of
551
 * the facets of the basic sets are also facets of the convex hull and
552
 * using all the facets so far to help in constructing the facets of the
553
 * facets
554
 * and/or
555
 * using the technique in section "3.1 Ridge Generation" of
556
 * "Extended Convex Hull" by Fukuda et al.
557
 */
558
static __isl_give isl_basic_set *extend(__isl_take isl_basic_set *hull,
559
  __isl_keep isl_set *set)
560
3.90k
{
561
3.90k
  int i, j, f;
562
3.90k
  int k;
563
3.90k
  struct isl_basic_set *facet = NULL;
564
3.90k
  struct isl_basic_set *hull_facet = NULL;
565
3.90k
  unsigned dim;
566
3.90k
567
3.90k
  if (!hull)
568
0
    return NULL;
569
3.90k
570
3.90k
  
isl_assert3.90k
(set->ctx, set->n > 0, goto error);3.90k
571
3.90k
572
3.90k
  dim = isl_set_n_dim(set);
573
3.90k
574
17.4k
  for (i = 0; 
i < hull->n_ineq17.4k
;
++i13.5k
)
{13.5k
575
13.5k
    facet = compute_facet(set, hull->ineq[i]);
576
13.5k
    facet = isl_basic_set_add_equality(facet, hull->ineq[i]);
577
13.5k
    facet = isl_basic_set_gauss(facet, NULL);
578
13.5k
    facet = isl_basic_set_normalize_constraints(facet);
579
13.5k
    hull_facet = isl_basic_set_copy(hull);
580
13.5k
    hull_facet = isl_basic_set_add_equality(hull_facet, hull->ineq[i]);
581
13.5k
    hull_facet = isl_basic_set_gauss(hull_facet, NULL);
582
13.5k
    hull_facet = isl_basic_set_normalize_constraints(hull_facet);
583
13.5k
    if (
!facet || 13.5k
!hull_facet13.5k
)
584
0
      goto error;
585
13.5k
    hull = isl_basic_set_cow(hull);
586
13.5k
    hull = isl_basic_set_extend_space(hull,
587
13.5k
      isl_space_copy(hull->dim), 0, 0, facet->n_ineq);
588
13.5k
    if (!hull)
589
0
      goto error;
590
48.2k
    
for (j = 0; 13.5k
j < facet->n_ineq48.2k
;
++j34.7k
)
{34.7k
591
67.1k
      for (f = 0; 
f < hull_facet->n_ineq67.1k
;
++f32.4k
)
592
57.4k
        
if (57.4k
isl_seq_eq(facet->ineq[j],57.4k
593
57.4k
            hull_facet->ineq[f], 1 + dim))
594
25.0k
          break;
595
34.7k
      if (f < hull_facet->n_ineq)
596
25.0k
        continue;
597
34.7k
      k = isl_basic_set_alloc_inequality(hull);
598
9.65k
      if (k < 0)
599
0
        goto error;
600
9.65k
      isl_seq_cpy(hull->ineq[k], hull->ineq[i], 1+dim);
601
9.65k
      if (!isl_set_wrap_facet(set, hull->ineq[k], facet->ineq[j]))
602
0
        goto error;
603
13.5k
    }
604
13.5k
    isl_basic_set_free(hull_facet);
605
13.5k
    isl_basic_set_free(facet);
606
13.5k
  }
607
3.90k
  hull = isl_basic_set_simplify(hull);
608
3.90k
  hull = isl_basic_set_finalize(hull);
609
3.90k
  return hull;
610
3.90k
error:
611
0
  isl_basic_set_free(hull_facet);
612
0
  isl_basic_set_free(facet);
613
0
  isl_basic_set_free(hull);
614
3.90k
  return NULL;
615
3.90k
}
616
617
/* Special case for computing the convex hull of a one dimensional set.
618
 * We simply collect the lower and upper bounds of each basic set
619
 * and the biggest of those.
620
 */
621
static __isl_give isl_basic_set *convex_hull_1d(__isl_take isl_set *set)
622
4.06k
{
623
4.06k
  struct isl_mat *c = NULL;
624
4.06k
  isl_int *lower = NULL;
625
4.06k
  isl_int *upper = NULL;
626
4.06k
  int i, j, k;
627
4.06k
  isl_int a, b;
628
4.06k
  struct isl_basic_set *hull;
629
4.06k
630
12.2k
  for (i = 0; 
i < set->n12.2k
;
++i8.13k
)
{8.13k
631
8.13k
    set->p[i] = isl_basic_set_simplify(set->p[i]);
632
8.13k
    if (!set->p[i])
633
0
      goto error;
634
8.13k
  }
635
4.06k
  set = isl_set_remove_empty_parts(set);
636
4.06k
  if (!set)
637
0
    goto error;
638
4.06k
  
isl_assert4.06k
(set->ctx, set->n > 0, goto error);4.06k
639
4.06k
  c = isl_mat_alloc(set->ctx, 2, 2);
640
4.06k
  if (!c)
641
0
    goto error;
642
4.06k
643
4.06k
  
if (4.06k
set->p[0]->n_eq > 04.06k
)
{4.06k
644
4.06k
    isl_assert(set->ctx, set->p[0]->n_eq == 1, goto error);
645
4.06k
    lower = c->row[0];
646
4.06k
    upper = c->row[1];
647
4.06k
    if (
isl_int_is_pos4.06k
(set->p[0]->eq[0][1]))
{4.06k
648
4.06k
      isl_seq_cpy(lower, set->p[0]->eq[0], 2);
649
4.06k
      isl_seq_neg(upper, set->p[0]->eq[0], 2);
650
4.06k
    } else {
651
0
      isl_seq_neg(lower, set->p[0]->eq[0], 2);
652
0
      isl_seq_cpy(upper, set->p[0]->eq[0], 2);
653
4.06k
    }
654
4.06k
  } else {
655
16
    for (j = 0; 
j < set->p[0]->n_ineq16
;
++j10
)
{10
656
10
      if (
isl_int_is_pos10
(set->p[0]->ineq[j][1]))
{6
657
6
        lower = c->row[0];
658
6
        isl_seq_cpy(lower, set->p[0]->ineq[j], 2);
659
10
      } else {
660
4
        upper = c->row[1];
661
4
        isl_seq_cpy(upper, set->p[0]->ineq[j], 2);
662
10
      }
663
10
    }
664
4.06k
  }
665
4.06k
666
4.06k
  
isl_int_init4.06k
(a);4.06k
667
4.06k
  isl_int_init(b);
668
12.2k
  for (i = 0; 
i < set->n12.2k
;
++i8.13k
)
{8.13k
669
8.13k
    struct isl_basic_set *bset = set->p[i];
670
8.13k
    int has_lower = 0;
671
8.13k
    int has_upper = 0;
672
8.13k
673
15.7k
    for (j = 0; 
j < bset->n_eq15.7k
;
++j7.59k
)
{7.59k
674
7.59k
      has_lower = 1;
675
7.59k
      has_upper = 1;
676
7.59k
      if (
lower7.59k
)
{7.59k
677
7.59k
        isl_int_mul(a, lower[0], bset->eq[j][1]);
678
7.59k
        isl_int_mul(b, lower[1], bset->eq[j][0]);
679
7.59k
        if (
isl_int_lt7.59k
(a, b) && 7.59k
isl_int_is_pos1.64k
(bset->eq[j][1]))
680
1.64k
          isl_seq_cpy(lower, bset->eq[j], 2);
681
7.59k
        if (
isl_int_gt7.59k
(a, b) && 7.59k
isl_int_is_neg1.89k
(bset->eq[j][1]))
682
0
          isl_seq_neg(lower, bset->eq[j], 2);
683
7.59k
      }
684
7.59k
      if (
upper7.59k
)
{7.59k
685
7.59k
        isl_int_mul(a, upper[0], bset->eq[j][1]);
686
7.59k
        isl_int_mul(b, upper[1], bset->eq[j][0]);
687
7.59k
        if (
isl_int_lt7.59k
(a, b) && 7.59k
isl_int_is_pos1.89k
(bset->eq[j][1]))
688
1.89k
          isl_seq_neg(upper, bset->eq[j], 2);
689
7.59k
        if (
isl_int_gt7.59k
(a, b) && 7.59k
isl_int_is_neg1.64k
(bset->eq[j][1]))
690
0
          isl_seq_cpy(upper, bset->eq[j], 2);
691
7.59k
      }
692
8.13k
    }
693
9.20k
    for (j = 0; 
j < bset->n_ineq9.20k
;
++j1.07k
)
{1.07k
694
1.07k
      if (isl_int_is_pos(bset->ineq[j][1]))
695
537
        has_lower = 1;
696
1.07k
      if (isl_int_is_neg(bset->ineq[j][1]))
697
535
        has_upper = 1;
698
1.07k
      if (
lower && 1.07k
isl_int_is_pos1.07k
(bset->ineq[j][1]))
{537
699
537
        isl_int_mul(a, lower[0], bset->ineq[j][1]);
700
537
        isl_int_mul(b, lower[1], bset->ineq[j][0]);
701
537
        if (isl_int_lt(a, b))
702
286
          isl_seq_cpy(lower, bset->ineq[j], 2);
703
1.07k
      }
704
1.07k
      if (
upper && 1.07k
isl_int_is_neg1.07k
(bset->ineq[j][1]))
{535
705
535
        isl_int_mul(a, upper[0], bset->ineq[j][1]);
706
535
        isl_int_mul(b, upper[1], bset->ineq[j][0]);
707
535
        if (isl_int_gt(a, b))
708
245
          isl_seq_cpy(upper, bset->ineq[j], 2);
709
1.07k
      }
710
8.13k
    }
711
8.13k
    if (!has_lower)
712
0
      lower = NULL;
713
8.13k
    if (!has_upper)
714
2
      upper = NULL;
715
8.13k
  }
716
4.06k
  isl_int_clear(a);
717
4.06k
  isl_int_clear(b);
718
4.06k
719
4.06k
  hull = isl_basic_set_alloc(set->ctx, 0, 1, 0, 0, 2);
720
4.06k
  hull = isl_basic_set_set_rational(hull);
721
4.06k
  if (!hull)
722
0
    goto error;
723
4.06k
  
if (4.06k
lower4.06k
)
{4.06k
724
4.06k
    k = isl_basic_set_alloc_inequality(hull);
725
4.06k
    isl_seq_cpy(hull->ineq[k], lower, 2);
726
4.06k
  }
727
4.06k
  if (
upper4.06k
)
{4.06k
728
4.06k
    k = isl_basic_set_alloc_inequality(hull);
729
4.06k
    isl_seq_cpy(hull->ineq[k], upper, 2);
730
4.06k
  }
731
4.06k
  hull = isl_basic_set_finalize(hull);
732
4.06k
  isl_set_free(set);
733
4.06k
  isl_mat_free(c);
734
4.06k
  return hull;
735
4.06k
error:
736
0
  isl_set_free(set);
737
0
  isl_mat_free(c);
738
4.06k
  return NULL;
739
4.06k
}
740
741
static __isl_give isl_basic_set *convex_hull_0d(__isl_take isl_set *set)
742
0
{
743
0
  struct isl_basic_set *convex_hull;
744
0
745
0
  if (!set)
746
0
    return NULL;
747
0
748
0
  
if (0
isl_set_is_empty(set)0
)
749
0
    convex_hull = isl_basic_set_empty(isl_space_copy(set->dim));
750
0
  else
751
0
    convex_hull = isl_basic_set_universe(isl_space_copy(set->dim));
752
0
  isl_set_free(set);
753
0
  return convex_hull;
754
0
}
755
756
/* Compute the convex hull of a pair of basic sets without any parameters or
757
 * integer divisions using Fourier-Motzkin elimination.
758
 * The convex hull is the set of all points that can be written as
759
 * the sum of points from both basic sets (in homogeneous coordinates).
760
 * We set up the constraints in a space with dimensions for each of
761
 * the three sets and then project out the dimensions corresponding
762
 * to the two original basic sets, retaining only those corresponding
763
 * to the convex hull.
764
 */
765
static __isl_give isl_basic_set *convex_hull_pair_elim(
766
  __isl_take isl_basic_set *bset1, __isl_take isl_basic_set *bset2)
767
14
{
768
14
  int i, j, k;
769
14
  struct isl_basic_set *bset[2];
770
14
  struct isl_basic_set *hull = NULL;
771
14
  unsigned dim;
772
14
773
14
  if (
!bset1 || 14
!bset214
)
774
0
    goto error;
775
14
776
14
  dim = isl_basic_set_n_dim(bset1);
777
14
  hull = isl_basic_set_alloc(bset1->ctx, 0, 2 + 3 * dim, 0,
778
14
        1 + dim + bset1->n_eq + bset2->n_eq,
779
14
        2 + bset1->n_ineq + bset2->n_ineq);
780
14
  bset[0] = bset1;
781
14
  bset[1] = bset2;
782
42
  for (i = 0; 
i < 242
;
++i28
)
{28
783
41
    for (j = 0; 
j < bset[i]->n_eq41
;
++j13
)
{13
784
13
      k = isl_basic_set_alloc_equality(hull);
785
13
      if (k < 0)
786
0
        goto error;
787
13
      isl_seq_clr(hull->eq[k], (i+1) * (1+dim));
788
13
      isl_seq_clr(hull->eq[k]+(i+2)*(1+dim), (1-i)*(1+dim));
789
13
      isl_seq_cpy(hull->eq[k]+(i+1)*(1+dim), bset[i]->eq[j],
790
13
          1+dim);
791
28
    }
792
126
    
for (j = 0; 28
j < bset[i]->n_ineq126
;
++j98
)
{98
793
98
      k = isl_basic_set_alloc_inequality(hull);
794
98
      if (k < 0)
795
0
        goto error;
796
98
      isl_seq_clr(hull->ineq[k], (i+1) * (1+dim));
797
98
      isl_seq_clr(hull->ineq[k]+(i+2)*(1+dim), (1-i)*(1+dim));
798
98
      isl_seq_cpy(hull->ineq[k]+(i+1)*(1+dim),
799
98
          bset[i]->ineq[j], 1+dim);
800
98
    }
801
28
    k = isl_basic_set_alloc_inequality(hull);
802
28
    if (k < 0)
803
0
      goto error;
804
28
    isl_seq_clr(hull->ineq[k], 1+2+3*dim);
805
28
    isl_int_set_si(hull->ineq[k][(i+1)*(1+dim)], 1);
806
28
  }
807
65
  
for (j = 0; 14
j < 1+dim65
;
++j51
)
{51
808
51
    k = isl_basic_set_alloc_equality(hull);
809
51
    if (k < 0)
810
0
      goto error;
811
51
    isl_seq_clr(hull->eq[k], 1+2+3*dim);
812
51
    isl_int_set_si(hull->eq[k][j], -1);
813
51
    isl_int_set_si(hull->eq[k][1+dim+j], 1);
814
51
    isl_int_set_si(hull->eq[k][2*(1+dim)+j], 1);
815
51
  }
816
14
  hull = isl_basic_set_set_rational(hull);
817
14
  hull = isl_basic_set_remove_dims(hull, isl_dim_set, dim, 2*(1+dim));
818
14
  hull = isl_basic_set_remove_redundancies(hull);
819
14
  isl_basic_set_free(bset1);
820
14
  isl_basic_set_free(bset2);
821
14
  return hull;
822
14
error:
823
0
  isl_basic_set_free(bset1);
824
0
  isl_basic_set_free(bset2);
825
0
  isl_basic_set_free(hull);
826
14
  return NULL;
827
14
}
828
829
/* Is the set bounded for each value of the parameters?
830
 */
831
isl_bool isl_basic_set_is_bounded(__isl_keep isl_basic_set *bset)
832
2.77k
{
833
2.77k
  struct isl_tab *tab;
834
2.77k
  isl_bool bounded;
835
2.77k
836
2.77k
  if (!bset)
837
0
    return isl_bool_error;
838
2.77k
  
if (2.77k
isl_basic_set_plain_is_empty(bset)2.77k
)
839
0
    return isl_bool_true;
840
2.77k
841
2.77k
  tab = isl_tab_from_recession_cone(bset, 1);
842
2.77k
  bounded = isl_tab_cone_is_bounded(tab);
843
2.77k
  isl_tab_free(tab);
844
2.77k
  return bounded;
845
2.77k
}
846
847
/* Is the image bounded for each value of the parameters and
848
 * the domain variables?
849
 */
850
isl_bool isl_basic_map_image_is_bounded(__isl_keep isl_basic_map *bmap)
851
0
{
852
0
  unsigned nparam = isl_basic_map_dim(bmap, isl_dim_param);
853
0
  unsigned n_in = isl_basic_map_dim(bmap, isl_dim_in);
854
0
  isl_bool bounded;
855
0
856
0
  bmap = isl_basic_map_copy(bmap);
857
0
  bmap = isl_basic_map_cow(bmap);
858
0
  bmap = isl_basic_map_move_dims(bmap, isl_dim_param, nparam,
859
0
          isl_dim_in, 0, n_in);
860
0
  bounded = isl_basic_set_is_bounded(bset_from_bmap(bmap));
861
0
  isl_basic_map_free(bmap);
862
0
863
0
  return bounded;
864
0
}
865
866
/* Is the set bounded for each value of the parameters?
867
 */
868
isl_bool isl_set_is_bounded(__isl_keep isl_set *set)
869
234
{
870
234
  int i;
871
234
872
234
  if (!set)
873
0
    return isl_bool_error;
874
234
875
712
  
for (i = 0; 234
i < set->n712
;
++i478
)
{510
876
510
    isl_bool bounded = isl_basic_set_is_bounded(set->p[i]);
877
510
    if (
!bounded || 510
bounded < 0478
)
878
32
      return bounded;
879
510
  }
880
202
  return isl_bool_true;
881
234
}
882
883
/* Compute the lineality space of the convex hull of bset1 and bset2.
884
 *
885
 * We first compute the intersection of the recession cone of bset1
886
 * with the negative of the recession cone of bset2 and then compute
887
 * the linear hull of the resulting cone.
888
 */
889
static __isl_give isl_basic_set *induced_lineality_space(
890
  __isl_take isl_basic_set *bset1, __isl_take isl_basic_set *bset2)
891
4
{
892
4
  int i, k;
893
4
  struct isl_basic_set *lin = NULL;
894
4
  unsigned dim;
895
4
896
4
  if (
!bset1 || 4
!bset24
)
897
0
    goto error;
898
4
899
4
  dim = isl_basic_set_total_dim(bset1);
900
4
  lin = isl_basic_set_alloc_space(isl_basic_set_get_space(bset1), 0,
901
4
          bset1->n_eq + bset2->n_eq,
902
4
          bset1->n_ineq + bset2->n_ineq);
903
4
  lin = isl_basic_set_set_rational(lin);
904
4
  if (!lin)
905
0
    goto error;
906
4
  
for (i = 0; 4
i < bset1->n_eq4
;
++i0
)
{0
907
0
    k = isl_basic_set_alloc_equality(lin);
908
0
    if (k < 0)
909
0
      goto error;
910
0
    
isl_int_set_si0
(lin->eq[k][0], 0);0
911
0
    isl_seq_cpy(lin->eq[k] + 1, bset1->eq[i] + 1, dim);
912
4
  }
913
25
  
for (i = 0; 4
i < bset1->n_ineq25
;
++i21
)
{21
914
21
    k = isl_basic_set_alloc_inequality(lin);
915
21
    if (k < 0)
916
0
      goto error;
917
21
    
isl_int_set_si21
(lin->ineq[k][0], 0);21
918
21
    isl_seq_cpy(lin->ineq[k] + 1, bset1->ineq[i] + 1, dim);
919
21
  }
920
4
  
for (i = 0; 4
i < bset2->n_eq4
;
++i0
)
{0
921
0
    k = isl_basic_set_alloc_equality(lin);
922
0
    if (k < 0)
923
0
      goto error;
924
0
    
isl_int_set_si0
(lin->eq[k][0], 0);0
925
0
    isl_seq_neg(lin->eq[k] + 1, bset2->eq[i] + 1, dim);
926
4
  }
927
29
  
for (i = 0; 4
i < bset2->n_ineq29
;
++i25
)
{25
928
25
    k = isl_basic_set_alloc_inequality(lin);
929
25
    if (k < 0)
930
0
      goto error;
931
25
    
isl_int_set_si25
(lin->ineq[k][0], 0);25
932
25
    isl_seq_neg(lin->ineq[k] + 1, bset2->ineq[i] + 1, dim);
933
25
  }
934
4
935
4
  isl_basic_set_free(bset1);
936
4
  isl_basic_set_free(bset2);
937
4
  return isl_basic_set_affine_hull(lin);
938
4
error:
939
0
  isl_basic_set_free(lin);
940
0
  isl_basic_set_free(bset1);
941
0
  isl_basic_set_free(bset2);
942
4
  return NULL;
943
4
}
944
945
static __isl_give isl_basic_set *uset_convex_hull(__isl_take isl_set *set);
946
947
/* Given a set and a linear space "lin" of dimension n > 0,
948
 * project the linear space from the set, compute the convex hull
949
 * and then map the set back to the original space.
950
 *
951
 * Let
952
 *
953
 *  M x = 0
954
 *
955
 * describe the linear space.  We first compute the Hermite normal
956
 * form H = M U of M = H Q, to obtain
957
 *
958
 *  H Q x = 0
959
 *
960
 * The last n rows of H will be zero, so the last n variables of x' = Q x
961
 * are the one we want to project out.  We do this by transforming each
962
 * basic set A x >= b to A U x' >= b and then removing the last n dimensions.
963
 * After computing the convex hull in x'_1, i.e., A' x'_1 >= b',
964
 * we transform the hull back to the original space as A' Q_1 x >= b',
965
 * with Q_1 all but the last n rows of Q.
966
 */
967
static __isl_give isl_basic_set *modulo_lineality(__isl_take isl_set *set,
968
  __isl_take isl_basic_set *lin)
969
12
{
970
12
  unsigned total = isl_basic_set_total_dim(lin);
971
12
  unsigned lin_dim;
972
12
  struct isl_basic_set *hull;
973
12
  struct isl_mat *M, *U, *Q;
974
12
975
12
  if (
!set || 12
!lin12
)
976
0
    goto error;
977
12
  lin_dim = total - lin->n_eq;
978
12
  M = isl_mat_sub_alloc6(set->ctx, lin->eq, 0, lin->n_eq, 1, total);
979
12
  M = isl_mat_left_hermite(M, 0, &U, &Q);
980
12
  if (!M)
981
0
    goto error;
982
12
  isl_mat_free(M);
983
12
  isl_basic_set_free(lin);
984
12
985
12
  Q = isl_mat_drop_rows(Q, Q->n_row - lin_dim, lin_dim);
986
12
987
12
  U = isl_mat_lin_to_aff(U);
988
12
  Q = isl_mat_lin_to_aff(Q);
989
12
990
12
  set = isl_set_preimage(set, U);
991
12
  set = isl_set_remove_dims(set, isl_dim_set, total - lin_dim, lin_dim);
992
12
  hull = uset_convex_hull(set);
993
12
  hull = isl_basic_set_preimage(hull, Q);
994
12
995
12
  return hull;
996
12
error:
997
0
  isl_basic_set_free(lin);
998
0
  isl_set_free(set);
999
12
  return NULL;
1000
12
}
1001
1002
/* Given two polyhedra with as constraints h_{ij} x >= 0 in homegeneous space,
1003
 * set up an LP for solving
1004
 *
1005
 *  \sum_j \alpha_{1j} h_{1j} = \sum_j \alpha_{2j} h_{2j}
1006
 *
1007
 * \alpha{i0} corresponds to the (implicit) positivity constraint 1 >= 0
1008
 * The next \alpha{ij} correspond to the equalities and come in pairs.
1009
 * The final \alpha{ij} correspond to the inequalities.
1010
 */
1011
static __isl_give isl_basic_set *valid_direction_lp(
1012
  __isl_take isl_basic_set *bset1, __isl_take isl_basic_set *bset2)
1013
8
{
1014
8
  isl_space *dim;
1015
8
  struct isl_basic_set *lp;
1016
8
  unsigned d;
1017
8
  int n;
1018
8
  int i, j, k;
1019
8
1020
8
  if (
!bset1 || 8
!bset28
)
1021
0
    goto error;
1022
8
  d = 1 + isl_basic_set_total_dim(bset1);
1023
8
  n = 2 +
1024
8
      2 * bset1->n_eq + bset1->n_ineq + 2 * bset2->n_eq + bset2->n_ineq;
1025
8
  dim = isl_space_set_alloc(bset1->ctx, 0, n);
1026
8
  lp = isl_basic_set_alloc_space(dim, 0, d, n);
1027
8
  if (!lp)
1028
0
    goto error;
1029
105
  
for (i = 0; 8
i < n105
;
++i97
)
{97
1030
97
    k = isl_basic_set_alloc_inequality(lp);
1031
97
    if (k < 0)
1032
0
      goto error;
1033
97
    isl_seq_clr(lp->ineq[k] + 1, n);
1034
97
    isl_int_set_si(lp->ineq[k][0], -1);
1035
97
    isl_int_set_si(lp->ineq[k][1 + i], 1);
1036
97
  }
1037
38
  
for (i = 0; 8
i < d38
;
++i30
)
{30
1038
30
    k = isl_basic_set_alloc_equality(lp);
1039
30
    if (k < 0)
1040
0
      goto error;
1041
30
    n = 0;
1042
30
    isl_int_set_si(lp->eq[k][n], 0); n++;
1043
30
    /* positivity constraint 1 >= 0 */
1044
30
    isl_int_set_si(lp->eq[k][n], i == 0); n++;
1045
44
    for (j = 0; 
j < bset1->n_eq44
;
++j14
)
{14
1046
14
      isl_int_set(lp->eq[k][n], bset1->eq[j][i]); n++;
1047
14
      isl_int_neg(lp->eq[k][n], bset1->eq[j][i]); n++;
1048
30
    }
1049
180
    for (j = 0; 
j < bset1->n_ineq180
;
++j150
)
{150
1050
150
      isl_int_set(lp->eq[k][n], bset1->ineq[j][i]); n++;
1051
150
    }
1052
30
    /* positivity constraint 1 >= 0 */
1053
30
    isl_int_set_si(lp->eq[k][n], -(i == 0)); n++;
1054
44
    for (j = 0; 
j < bset2->n_eq44
;
++j14
)
{14
1055
14
      isl_int_neg(lp->eq[k][n], bset2->eq[j][i]); n++;
1056
14
      isl_int_set(lp->eq[k][n], bset2->eq[j][i]); n++;
1057
30
    }
1058
204
    for (j = 0; 
j < bset2->n_ineq204
;
++j174
)
{174
1059
174
      isl_int_neg(lp->eq[k][n], bset2->ineq[j][i]); n++;
1060
174
    }
1061
30
  }
1062
8
  lp = isl_basic_set_gauss(lp, NULL);
1063
8
  isl_basic_set_free(bset1);
1064
8
  isl_basic_set_free(bset2);
1065
8
  return lp;
1066
8
error:
1067
0
  isl_basic_set_free(bset1);
1068
0
  isl_basic_set_free(bset2);
1069
8
  return NULL;
1070
8
}
1071
1072
/* Compute a vector s in the homogeneous space such that <s, r> > 0
1073
 * for all rays in the homogeneous space of the two cones that correspond
1074
 * to the input polyhedra bset1 and bset2.
1075
 *
1076
 * We compute s as a vector that satisfies
1077
 *
1078
 *  s = \sum_j \alpha_{ij} h_{ij} for i = 1,2     (*)
1079
 *
1080
 * with h_{ij} the normals of the facets of polyhedron i
1081
 * (including the "positivity constraint" 1 >= 0) and \alpha_{ij}
1082
 * strictly positive numbers.  For simplicity we impose \alpha_{ij} >= 1.
1083
 * We first set up an LP with as variables the \alpha{ij}.
1084
 * In this formulation, for each polyhedron i,
1085
 * the first constraint is the positivity constraint, followed by pairs
1086
 * of variables for the equalities, followed by variables for the inequalities.
1087
 * We then simply pick a feasible solution and compute s using (*).
1088
 *
1089
 * Note that we simply pick any valid direction and make no attempt
1090
 * to pick a "good" or even the "best" valid direction.
1091
 */
1092
static __isl_give isl_vec *valid_direction(
1093
  __isl_take isl_basic_set *bset1, __isl_take isl_basic_set *bset2)
1094
8
{
1095
8
  struct isl_basic_set *lp;
1096
8
  struct isl_tab *tab;
1097
8
  struct isl_vec *sample = NULL;
1098
8
  struct isl_vec *dir;
1099
8
  unsigned d;
1100
8
  int i;
1101
8
  int n;
1102
8
1103
8
  if (
!bset1 || 8
!bset28
)
1104
0
    goto error;
1105
8
  lp = valid_direction_lp(isl_basic_set_copy(bset1),
1106
8
        isl_basic_set_copy(bset2));
1107
8
  tab = isl_tab_from_basic_set(lp, 0);
1108
8
  sample = isl_tab_get_sample_value(tab);
1109
8
  isl_tab_free(tab);
1110
8
  isl_basic_set_free(lp);
1111
8
  if (!sample)
1112
0
    goto error;
1113
8
  d = isl_basic_set_total_dim(bset1);
1114
8
  dir = isl_vec_alloc(bset1->ctx, 1 + d);
1115
8
  if (!dir)
1116
0
    goto error;
1117
8
  isl_seq_clr(dir->block.data + 1, dir->size - 1);
1118
8
  n = 1;
1119
8
  /* positivity constraint 1 >= 0 */
1120
8
  isl_int_set(dir->block.data[0], sample->block.data[n]); n++;
1121
12
  for (i = 0; 
i < bset1->n_eq12
;
++i4
)
{4
1122
4
    isl_int_sub(sample->block.data[n],
1123
4
          sample->block.data[n], sample->block.data[n+1]);
1124
4
    isl_seq_combine(dir->block.data,
1125
4
        bset1->ctx->one, dir->block.data,
1126
4
        sample->block.data[n], bset1->eq[i], 1 + d);
1127
4
1128
4
    n += 2;
1129
8
  }
1130
39
  for (i = 0; 
i < bset1->n_ineq39
;
++i31
)
1131
31
    isl_seq_combine(dir->block.data,
1132
31
        bset1->ctx->one, dir->block.data,
1133
31
        sample->block.data[n++], bset1->ineq[i], 1 + d);
1134
8
  isl_vec_free(sample);
1135
8
  isl_seq_normalize(bset1->ctx, dir->el, dir->size);
1136
8
  isl_basic_set_free(bset1);
1137
8
  isl_basic_set_free(bset2);
1138
8
  return dir;
1139
8
error:
1140
0
  isl_vec_free(sample);
1141
0
  isl_basic_set_free(bset1);
1142
0
  isl_basic_set_free(bset2);
1143
8
  return NULL;
1144
8
}
1145
1146
/* Given a polyhedron b_i + A_i x >= 0 and a map T = S^{-1},
1147
 * compute b_i' + A_i' x' >= 0, with
1148
 *
1149
 *  [ b_i A_i ]        [ y' ]                 [ y' ]
1150
 *  [  1   0  ] S^{-1} [ x' ] >= 0  or  [ b_i' A_i' ] [ x' ] >= 0
1151
 *
1152
 * In particular, add the "positivity constraint" and then perform
1153
 * the mapping.
1154
 */
1155
static __isl_give isl_basic_set *homogeneous_map(__isl_take isl_basic_set *bset,
1156
  __isl_take isl_mat *T)
1157
16
{
1158
16
  int k;
1159
16
1160
16
  if (!bset)
1161
0
    goto error;
1162
16
  bset = isl_basic_set_extend_constraints(bset, 0, 1);
1163
16
  k = isl_basic_set_alloc_inequality(bset);
1164
16
  if (k < 0)
1165
0
    goto error;
1166
16
  isl_seq_clr(bset->ineq[k] + 1, isl_basic_set_total_dim(bset));
1167
16
  isl_int_set_si(bset->ineq[k][0], 1);
1168
16
  bset = isl_basic_set_preimage(bset, T);
1169
16
  return bset;
1170
16
error:
1171
0
  isl_mat_free(T);
1172
0
  isl_basic_set_free(bset);
1173
16
  return NULL;
1174
16
}
1175
1176
/* Compute the convex hull of a pair of basic sets without any parameters or
1177
 * integer divisions, where the convex hull is known to be pointed,
1178
 * but the basic sets may be unbounded.
1179
 *
1180
 * We turn this problem into the computation of a convex hull of a pair
1181
 * _bounded_ polyhedra by "changing the direction of the homogeneous
1182
 * dimension".  This idea is due to Matthias Koeppe.
1183
 *
1184
 * Consider the cones in homogeneous space that correspond to the
1185
 * input polyhedra.  The rays of these cones are also rays of the
1186
 * polyhedra if the coordinate that corresponds to the homogeneous
1187
 * dimension is zero.  That is, if the inner product of the rays
1188
 * with the homogeneous direction is zero.
1189
 * The cones in the homogeneous space can also be considered to
1190
 * correspond to other pairs of polyhedra by chosing a different
1191
 * homogeneous direction.  To ensure that both of these polyhedra
1192
 * are bounded, we need to make sure that all rays of the cones
1193
 * correspond to vertices and not to rays.
1194
 * Let s be a direction such that <s, r> > 0 for all rays r of both cones.
1195
 * Then using s as a homogeneous direction, we obtain a pair of polytopes.
1196
 * The vector s is computed in valid_direction.
1197
 *
1198
 * Note that we need to consider _all_ rays of the cones and not just
1199
 * the rays that correspond to rays in the polyhedra.  If we were to
1200
 * only consider those rays and turn them into vertices, then we
1201
 * may inadvertently turn some vertices into rays.
1202
 *
1203
 * The standard homogeneous direction is the unit vector in the 0th coordinate.
1204
 * We therefore transform the two polyhedra such that the selected
1205
 * direction is mapped onto this standard direction and then proceed
1206
 * with the normal computation.
1207
 * Let S be a non-singular square matrix with s as its first row,
1208
 * then we want to map the polyhedra to the space
1209
 *
1210
 *  [ y' ]     [ y ]    [ y ]          [ y' ]
1211
 *  [ x' ] = S [ x ]  i.e., [ x ] = S^{-1} [ x' ]
1212
 *
1213
 * We take S to be the unimodular completion of s to limit the growth
1214
 * of the coefficients in the following computations.
1215
 *
1216
 * Let b_i + A_i x >= 0 be the constraints of polyhedron i.
1217
 * We first move to the homogeneous dimension
1218
 *
1219
 *  b_i y + A_i x >= 0    [ b_i A_i ] [ y ]    [ 0 ]
1220
 *      y         >= 0  or  [  1   0  ] [ x ] >= [ 0 ]
1221
 *
1222
 * Then we change directoin
1223
 *
1224
 *  [ b_i A_i ]        [ y' ]                 [ y' ]
1225
 *  [  1   0  ] S^{-1} [ x' ] >= 0  or  [ b_i' A_i' ] [ x' ] >= 0
1226
 *
1227
 * Then we compute the convex hull of the polytopes b_i' + A_i' x' >= 0
1228
 * resulting in b' + A' x' >= 0, which we then convert back
1229
 *
1230
 *              [ y ]           [ y ]
1231
 *  [ b' A' ] S [ x ] >= 0  or  [ b A ] [ x ] >= 0
1232
 *
1233
 * The polyhedron b + A x >= 0 is then the convex hull of the input polyhedra.
1234
 */
1235
static __isl_give isl_basic_set *convex_hull_pair_pointed(
1236
  __isl_take isl_basic_set *bset1, __isl_take isl_basic_set *bset2)
1237
8
{
1238
8
  struct isl_ctx *ctx = NULL;
1239
8
  struct isl_vec *dir = NULL;
1240
8
  struct isl_mat *T = NULL;
1241
8
  struct isl_mat *T2 = NULL;
1242
8
  struct isl_basic_set *hull;
1243
8
  struct isl_set *set;
1244
8
1245
8
  if (
!bset1 || 8
!bset28
)
1246
0
    goto error;
1247
8
  ctx = isl_basic_set_get_ctx(bset1);
1248
8
  dir = valid_direction(isl_basic_set_copy(bset1),
1249
8
        isl_basic_set_copy(bset2));
1250
8
  if (!dir)
1251
0
    goto error;
1252
8
  T = isl_mat_alloc(ctx, dir->size, dir->size);
1253
8
  if (!T)
1254
0
    goto error;
1255
8
  isl_seq_cpy(T->row[0], dir->block.data, dir->size);
1256
8
  T = isl_mat_unimodular_complete(T, 1);
1257
8
  T2 = isl_mat_right_inverse(isl_mat_copy(T));
1258
8
1259
8
  bset1 = homogeneous_map(bset1, isl_mat_copy(T2));
1260
8
  bset2 = homogeneous_map(bset2, T2);
1261
8
  set = isl_set_alloc_space(isl_basic_set_get_space(bset1), 2, 0);
1262
8
  set = isl_set_add_basic_set(set, bset1);
1263
8
  set = isl_set_add_basic_set(set, bset2);
1264
8
  hull = uset_convex_hull(set);
1265
8
  hull = isl_basic_set_preimage(hull, T);
1266
8
   
1267
8
  isl_vec_free(dir);
1268
8
1269
8
  return hull;
1270
8
error:
1271
0
  isl_vec_free(dir);
1272
0
  isl_basic_set_free(bset1);
1273
0
  isl_basic_set_free(bset2);
1274
8
  return NULL;
1275
8
}
1276
1277
static __isl_give isl_basic_set *uset_convex_hull_wrap(__isl_take isl_set *set);
1278
static __isl_give isl_basic_set *modulo_affine_hull(
1279
  __isl_take isl_set *set, __isl_take isl_basic_set *affine_hull);
1280
1281
/* Compute the convex hull of a pair of basic sets without any parameters or
1282
 * integer divisions.
1283
 *
1284
 * This function is called from uset_convex_hull_unbounded, which
1285
 * means that the complete convex hull is unbounded.  Some pairs
1286
 * of basic sets may still be bounded, though.
1287
 * They may even lie inside a lower dimensional space, in which
1288
 * case they need to be handled inside their affine hull since
1289
 * the main algorithm assumes that the result is full-dimensional.
1290
 *
1291
 * If the convex hull of the two basic sets would have a non-trivial
1292
 * lineality space, we first project out this lineality space.
1293
 */
1294
static __isl_give isl_basic_set *convex_hull_pair(
1295
  __isl_take isl_basic_set *bset1, __isl_take isl_basic_set *bset2)
1296
26
{
1297
26
  isl_basic_set *lin, *aff;
1298
26
  int bounded1, bounded2;
1299
26
1300
26
  if (
bset1->ctx->opt->convex == 26
ISL_CONVEX_HULL_FM26
)
1301
14
    return convex_hull_pair_elim(bset1, bset2);
1302
26
1303
26
  aff = isl_set_affine_hull(isl_basic_set_union(isl_basic_set_copy(bset1),
1304
12
                isl_basic_set_copy(bset2)));
1305
12
  if (!aff)
1306
0
    goto error;
1307
12
  
if (12
aff->n_eq != 012
)
1308
1
    return modulo_affine_hull(isl_basic_set_union(bset1, bset2), aff);
1309
12
  isl_basic_set_free(aff);
1310
11
1311
11
  bounded1 = isl_basic_set_is_bounded(bset1);
1312
11
  bounded2 = isl_basic_set_is_bounded(bset2);
1313
11
1314
11
  if (
bounded1 < 0 || 11
bounded2 < 011
)
1315
0
    goto error;
1316
11
1317
11
  
if (11
bounded1 && 11
bounded24
)
1318
1
    return uset_convex_hull_wrap(isl_basic_set_union(bset1, bset2));
1319
11
1320
10
  
if (10
bounded1 || 10
bounded27
)
1321
6
    return convex_hull_pair_pointed(bset1, bset2);
1322
10
1323
10
  lin = induced_lineality_space(isl_basic_set_copy(bset1),
1324
4
              isl_basic_set_copy(bset2));
1325
4
  if (!lin)
1326
0
    goto error;
1327
4
  
if (4
isl_basic_set_plain_is_universe(lin)4
)
{0
1328
0
    isl_basic_set_free(bset1);
1329
0
    isl_basic_set_free(bset2);
1330
0
    return lin;
1331
4
  }
1332
4
  
if (4
lin->n_eq < isl_basic_set_total_dim(lin)4
)
{2
1333
2
    struct isl_set *set;
1334
2
    set = isl_set_alloc_space(isl_basic_set_get_space(bset1), 2, 0);
1335
2
    set = isl_set_add_basic_set(set, bset1);
1336
2
    set = isl_set_add_basic_set(set, bset2);
1337
2
    return modulo_lineality(set, lin);
1338
4
  }
1339
4
  isl_basic_set_free(lin);
1340
2
1341
4
  return convex_hull_pair_pointed(bset1, bset2);
1342
4
error:
1343
0
  isl_basic_set_free(bset1);
1344
0
  isl_basic_set_free(bset2);
1345
4
  return NULL;
1346
26
}
1347
1348
/* Compute the lineality space of a basic set.
1349
 * We basically just drop the constants and turn every inequality
1350
 * into an equality.
1351
 * Any explicit representations of local variables are removed
1352
 * because they may no longer be valid representations
1353
 * in the lineality space.
1354
 */
1355
__isl_give isl_basic_set *isl_basic_set_lineality_space(
1356
  __isl_take isl_basic_set *bset)
1357
119
{
1358
119
  int i, k;
1359
119
  struct isl_basic_set *lin = NULL;
1360
119
  unsigned n_div, dim;
1361
119
1362
119
  if (!bset)
1363
0
    goto error;
1364
119
  n_div = isl_basic_set_dim(bset, isl_dim_div);
1365
119
  dim = isl_basic_set_total_dim(bset);
1366
119
1367
119
  lin = isl_basic_set_alloc_space(isl_basic_set_get_space(bset),
1368
119
          n_div, dim, 0);
1369
119
  for (i = 0; 
i < n_div119
;
++i0
)
1370
0
    
if (0
isl_basic_set_alloc_div(lin) < 00
)
1371
0
      goto error;
1372
119
  
if (119
!lin119
)
1373
0
    goto error;
1374
237
  
for (i = 0; 119
i < bset->n_eq237
;
++i118
)
{118
1375
118
    k = isl_basic_set_alloc_equality(lin);
1376
118
    if (k < 0)
1377
0
      goto error;
1378
118
    
isl_int_set_si118
(lin->eq[k][0], 0);118
1379
118
    isl_seq_cpy(lin->eq[k] + 1, bset->eq[i] + 1, dim);
1380
119
  }
1381
119
  lin = isl_basic_set_gauss(lin, NULL);
1382
119
  if (!lin)
1383
0
    goto error;
1384
304
  
for (i = 0; 119
i < bset->n_ineq && 304
lin->n_eq < dim222
;
++i185
)
{185
1385
185
    k = isl_basic_set_alloc_equality(lin);
1386
185
    if (k < 0)
1387
0
      goto error;
1388
185
    
isl_int_set_si185
(lin->eq[k][0], 0);185
1389
185
    isl_seq_cpy(lin->eq[k] + 1, bset->ineq[i] + 1, dim);
1390
185
    lin = isl_basic_set_gauss(lin, NULL);
1391
185
    if (!lin)
1392
0
      goto error;
1393
185
  }
1394
119
  isl_basic_set_free(bset);
1395
119
  return lin;
1396
119
error:
1397
0
  isl_basic_set_free(lin);
1398
0
  isl_basic_set_free(bset);
1399
119
  return NULL;
1400
119
}
1401
1402
/* Compute the (linear) hull of the lineality spaces of the basic sets in the
1403
 * set "set".
1404
 */
1405
__isl_give isl_basic_set *isl_set_combined_lineality_space(
1406
  __isl_take isl_set *set)
1407
55
{
1408
55
  int i;
1409
55
  struct isl_set *lin = NULL;
1410
55
1411
55
  if (!set)
1412
0
    return NULL;
1413
55
  
if (55
set->n == 055
)
{0
1414
0
    isl_space *space = isl_set_get_space(set);
1415
0
    isl_set_free(set);
1416
0
    return isl_basic_set_empty(space);
1417
55
  }
1418
55
1419
55
  lin = isl_set_alloc_space(isl_set_get_space(set), set->n, 0);
1420
170
  for (i = 0; 
i < set->n170
;
++i115
)
1421
115
    lin = isl_set_add_basic_set(lin,
1422
115
        isl_basic_set_lineality_space(isl_basic_set_copy(set->p[i])));
1423
55
  isl_set_free(set);
1424
55
  return isl_set_affine_hull(lin);
1425
55
}
1426
1427
/* Compute the convex hull of a set without any parameters or
1428
 * integer divisions.
1429
 * In each step, we combined two basic sets until only one
1430
 * basic set is left.
1431
 * The input basic sets are assumed not to have a non-trivial
1432
 * lineality space.  If any of the intermediate results has
1433
 * a non-trivial lineality space, it is projected out.
1434
 */
1435
static __isl_give isl_basic_set *uset_convex_hull_unbounded(
1436
  __isl_take isl_set *set)
1437
24
{
1438
24
  isl_basic_set_list *list;
1439
24
1440
24
  list = isl_set_get_basic_set_list(set);
1441
24
  isl_set_free(set);
1442
24
1443
26
  while (
list26
)
{26
1444
26
    int n;
1445
26
    struct isl_basic_set *t;
1446
26
    isl_basic_set *bset1, *bset2;
1447
26
1448
26
    n = isl_basic_set_list_n_basic_set(list);
1449
26
    if (n < 2)
1450
0
      isl_die(isl_basic_set_list_get_ctx(list),
1451
26
        isl_error_internal,
1452
26
        "expecting at least two elements", goto error);
1453
26
    bset1 = isl_basic_set_list_get_basic_set(list, n - 1);
1454
26
    bset2 = isl_basic_set_list_get_basic_set(list, n - 2);
1455
26
    bset1 = convex_hull_pair(bset1, bset2);
1456
26
    if (
n == 226
)
{22
1457
22
      isl_basic_set_list_free(list);
1458
22
      return bset1;
1459
26
    }
1460
26
    bset1 = isl_basic_set_underlying_set(bset1);
1461
4
    list = isl_basic_set_list_drop(list, n - 2, 2);
1462
4
    list = isl_basic_set_list_add(list, bset1);
1463
4
1464
4
    t = isl_basic_set_list_get_basic_set(list, n - 2);
1465
4
    t = isl_basic_set_lineality_space(t);
1466
4
    if (!t)
1467
0
      goto error;
1468
4
    
if (4
isl_basic_set_plain_is_universe(t)4
)
{0
1469
0
      isl_basic_set_list_free(list);
1470
0
      return t;
1471
4
    }
1472
4
    
if (4
t->n_eq < isl_basic_set_total_dim(t)4
)
{2
1473
2
      set = isl_basic_set_list_union(list);
1474
2
      return modulo_lineality(set, t);
1475
4
    }
1476
4
    isl_basic_set_free(t);
1477
24
  }
1478
24
1479
0
  return NULL;
1480
24
error:
1481
0
  isl_basic_set_list_free(list);
1482
24
  return NULL;
1483
24
}
1484
1485
/* Compute an initial hull for wrapping containing a single initial
1486
 * facet.
1487
 * This function assumes that the given set is bounded.
1488
 */
1489
static __isl_give isl_basic_set *initial_hull(__isl_take isl_basic_set *hull,
1490
  __isl_keep isl_set *set)
1491
3.88k
{
1492
3.88k
  struct isl_mat *bounds = NULL;
1493
3.88k
  unsigned dim;
1494
3.88k
  int k;
1495
3.88k
1496
3.88k
  if (!hull)
1497
0
    goto error;
1498
3.88k
  bounds = initial_facet_constraint(set);
1499
3.88k
  if (!bounds)
1500
0
    goto error;
1501
3.88k
  k = isl_basic_set_alloc_inequality(hull);
1502
3.88k
  if (k < 0)
1503
0
    goto error;
1504
3.88k
  dim = isl_set_n_dim(set);
1505
3.88k
  isl_assert(set->ctx, 1 + dim == bounds->n_col, goto error);
1506
3.88k
  isl_seq_cpy(hull->ineq[k], bounds->row[0], bounds->n_col);
1507
3.88k
  isl_mat_free(bounds);
1508
3.88k
1509
3.88k
  return hull;
1510
3.88k
error:
1511
0
  isl_basic_set_free(hull);
1512
0
  isl_mat_free(bounds);
1513
3.88k
  return NULL;
1514
3.88k
}
1515
1516
struct max_constraint {
1517
  struct isl_mat *c;
1518
  int   count;
1519
  int   ineq;
1520
};
1521
1522
static int max_constraint_equal(const void *entry, const void *val)
1523
37
{
1524
37
  struct max_constraint *a = (struct max_constraint *)entry;
1525
37
  isl_int *b = (isl_int *)val;
1526
37
1527
37
  return isl_seq_eq(a->c->row[0] + 1, b, a->c->n_col - 1);
1528
37
}
1529
1530
static void update_constraint(struct isl_ctx *ctx, struct isl_hash_table *table,
1531
  isl_int *con, unsigned len, int n, int ineq)
1532
3.08k
{
1533
3.08k
  struct isl_hash_table_entry *entry;
1534
3.08k
  struct max_constraint *c;
1535
3.08k
  uint32_t c_hash;
1536
3.08k
1537
3.08k
  c_hash = isl_seq_get_hash(con + 1, len);
1538
3.08k
  entry = isl_hash_table_find(ctx, table, c_hash, max_constraint_equal,
1539
3.08k
      con + 1, 0);
1540
3.08k
  if (!entry)
1541
3.05k
    return;
1542
3.08k
  c = entry->data;
1543
37
  if (
c->count < n37
)
{0
1544
0
    isl_hash_table_remove(ctx, table, entry);
1545
0
    return;
1546
37
  }
1547
37
  c->count++;
1548
37
  if (isl_int_gt(c->c->row[0][0], con[0]))
1549
3
    return;
1550
34
  
if (34
isl_int_eq34
(c->c->row[0][0], con[0]))
{34
1551
34
    if (ineq)
1552
16
      c->ineq = ineq;
1553
34
    return;
1554
34
  }
1555
34
  c->c = isl_mat_cow(c->c);
1556
0
  isl_int_set(c->c->row[0][0], con[0]);
1557
0
  c->ineq = ineq;
1558
3.08k
}
1559
1560
/* Check whether the constraint hash table "table" constains the constraint
1561
 * "con".
1562
 */
1563
static int has_constraint(struct isl_ctx *ctx, struct isl_hash_table *table,
1564
  isl_int *con, unsigned len, int n)
1565
0
{
1566
0
  struct isl_hash_table_entry *entry;
1567
0
  struct max_constraint *c;
1568
0
  uint32_t c_hash;
1569
0
1570
0
  c_hash = isl_seq_get_hash(con + 1, len);
1571
0
  entry = isl_hash_table_find(ctx, table, c_hash, max_constraint_equal,
1572
0
      con + 1, 0);
1573
0
  if (!entry)
1574
0
    return 0;
1575
0
  c = entry->data;
1576
0
  if (c->count < n)
1577
0
    return 0;
1578
0
  
return 0
isl_int_eq0
(c->c->row[0][0], con[0]);
1579
0
}
1580
1581
/* Check for inequality constraints of a basic set without equalities
1582
 * such that the same or more stringent copies of the constraint appear
1583
 * in all of the basic sets.  Such constraints are necessarily facet
1584
 * constraints of the convex hull.
1585
 *
1586
 * If the resulting basic set is by chance identical to one of
1587
 * the basic sets in "set", then we know that this basic set contains
1588
 * all other basic sets and is therefore the convex hull of set.
1589
 * In this case we set *is_hull to 1.
1590
 */
1591
static __isl_give isl_basic_set *common_constraints(
1592
  __isl_take isl_basic_set *hull, __isl_keep isl_set *set, int *is_hull)
1593
3.90k
{
1594
3.90k
  int i, j, s, n;
1595
3.90k
  int min_constraints;
1596
3.90k
  int best;
1597
3.90k
  struct max_constraint *constraints = NULL;
1598
3.90k
  struct isl_hash_table *table = NULL;
1599
3.90k
  unsigned total;
1600
3.90k
1601
3.90k
  *is_hull = 0;
1602
3.90k
1603
11.0k
  for (i = 0; 
i < set->n11.0k
;
++i7.10k
)
1604
7.77k
    
if (7.77k
set->p[i]->n_eq == 07.77k
)
1605
674
      break;
1606
3.90k
  if (i >= set->n)
1607
3.23k
    return hull;
1608
3.90k
  min_constraints = set->p[i]->n_ineq;
1609
674
  best = i;
1610
710
  for (i = best + 1; 
i < set->n710
;
++i36
)
{36
1611
36
    if (set->p[i]->n_eq != 0)
1612
8
      continue;
1613
28
    
if (28
set->p[i]->n_ineq >= min_constraints28
)
1614
27
      continue;
1615
28
    min_constraints = set->p[i]->n_ineq;
1616
1
    best = i;
1617
674
  }
1618
674
  constraints = isl_calloc_array(hull->ctx, struct max_constraint,
1619
674
          min_constraints);
1620
674
  if (!constraints)
1621
0
    return hull;
1622
674
  
table = 674
isl_alloc_type674
(hull->ctx, struct isl_hash_table);
1623
674
  if (isl_hash_table_init(hull->ctx, table, min_constraints))
1624
0
    goto error;
1625
674
1626
674
  total = isl_space_dim(set->dim, isl_dim_all);
1627
3.15k
  for (i = 0; 
i < set->p[best]->n_ineq3.15k
;
++i2.48k
)
{2.48k
1628
2.48k
    constraints[i].c = isl_mat_sub_alloc6(hull->ctx,
1629
2.48k
      set->p[best]->ineq + i, 0, 1, 0, 1 + total);
1630
2.48k
    if (!constraints[i].c)
1631
0
      goto error;
1632
2.48k
    constraints[i].ineq = 1;
1633
2.48k
  }
1634
3.15k
  
for (i = 0; 674
i < min_constraints3.15k
;
++i2.48k
)
{2.48k
1635
2.48k
    struct isl_hash_table_entry *entry;
1636
2.48k
    uint32_t c_hash;
1637
2.48k
    c_hash = isl_seq_get_hash(constraints[i].c->row[0] + 1, total);
1638
2.48k
    entry = isl_hash_table_find(hull->ctx, table, c_hash,
1639
2.48k
      max_constraint_equal, constraints[i].c->row[0] + 1, 1);
1640
2.48k
    if (!entry)
1641
0
      goto error;
1642
2.48k
    
isl_assert2.48k
(hull->ctx, !entry->data, goto error);2.48k
1643
2.48k
    entry->data = &constraints[i];
1644
2.48k
  }
1645
674
1646
674
  n = 0;
1647
2.02k
  for (s = 0; 
s < set->n2.02k
;
++s1.34k
)
{1.34k
1648
1.34k
    if (s == best)
1649
674
      continue;
1650
1.34k
1651
1.37k
    
for (i = 0; 674
i < set->p[s]->n_eq1.37k
;
++i696
)
{696
1652
696
      isl_int *eq = set->p[s]->eq[i];
1653
2.08k
      for (j = 0; 
j < 22.08k
;
++j1.39k
)
{1.39k
1654
1.39k
        isl_seq_neg(eq, eq, 1 + total);
1655
1.39k
        update_constraint(hull->ctx, table,
1656
1.39k
                  eq, total, n, 0);
1657
1.39k
      }
1658
696
    }
1659
2.37k
    for (i = 0; 
i < set->p[s]->n_ineq2.37k
;
++i1.69k
)
{1.69k
1660
1.69k
      isl_int *ineq = set->p[s]->ineq[i];
1661
1.69k
      update_constraint(hull->ctx, table, ineq, total, n,
1662
1.69k
        set->p[s]->n_eq == 0);
1663
1.69k
    }
1664
674
    ++n;
1665
674
  }
1666
674
1667
3.15k
  for (i = 0; 
i < min_constraints3.15k
;
++i2.48k
)
{2.48k
1668
2.48k
    if (constraints[i].count < n)
1669
2.44k
      continue;
1670
37
    
if (37
!constraints[i].ineq37
)
1671
0
      continue;
1672
37
    j = isl_basic_set_alloc_inequality(hull);
1673
37
    if (j < 0)
1674
0
      goto error;
1675
37
    isl_seq_cpy(hull->ineq[j], constraints[i].c->row[0], 1 + total);
1676
674
  }
1677
674
1678
2.02k
  
for (s = 0; 674
s < set->n2.02k
;
++s1.34k
)
{1.34k
1679
1.34k
    if (set->p[s]->n_eq)
1680
646
      continue;
1681
702
    
if (702
set->p[s]->n_ineq != hull->n_ineq702
)
1682
702
      continue;
1683
0
    
for (i = 0; 0
i < set->p[s]->n_ineq0
;
++i0
)
{0
1684
0
      isl_int *ineq = set->p[s]->ineq[i];
1685
0
      if (!has_constraint(hull->ctx, table, ineq, total, n))
1686
0
        break;
1687
0
    }
1688
0
    if (i == set->p[s]->n_ineq)
1689
0
      *is_hull = 1;
1690
674
  }
1691
674
1692
674
  isl_hash_table_clear(table);
1693
3.15k
  for (i = 0; 
i < min_constraints3.15k
;
++i2.48k
)
1694
2.48k
    isl_mat_free(constraints[i].c);
1695
674
  free(constraints);
1696
674
  free(table);
1697
674
  return hull;
1698
674
error:
1699
0
  isl_hash_table_clear(table);
1700
0
  free(table);
1701
0
  if (constraints)
1702
0
    
for (i = 0; 0
i < min_constraints0
;
++i0
)
1703
0
      isl_mat_free(constraints[i].c);
1704
0
  free(constraints);
1705
674
  return hull;
1706
3.90k
}
1707
1708
/* Create a template for the convex hull of "set" and fill it up
1709
 * obvious facet constraints, if any.  If the result happens to
1710
 * be the convex hull of "set" then *is_hull is set to 1.
1711
 */
1712
static __isl_give isl_basic_set *proto_hull(__isl_keep isl_set *set,
1713
  int *is_hull)
1714
3.90k
{
1715
3.90k
  struct isl_basic_set *hull;
1716
3.90k
  unsigned n_ineq;
1717
3.90k
  int i;
1718
3.90k
1719
3.90k
  n_ineq = 1;
1720
11.7k
  for (i = 0; 
i < set->n11.7k
;
++i7.81k
)
{7.81k
1721
7.81k
    n_ineq += set->p[i]->n_eq;
1722
7.81k
    n_ineq += set->p[i]->n_ineq;
1723
7.81k
  }
1724
3.90k
  hull = isl_basic_set_alloc_space(isl_space_copy(set->dim), 0, 0, n_ineq);
1725
3.90k
  hull = isl_basic_set_set_rational(hull);
1726
3.90k
  if (!hull)
1727
0
    return NULL;
1728
3.90k
  return common_constraints(hull, set, is_hull);
1729
3.90k
}
1730
1731
static __isl_give isl_basic_set *uset_convex_hull_wrap(__isl_take isl_set *set)
1732
3.90k
{
1733
3.90k
  struct isl_basic_set *hull;
1734
3.90k
  int is_hull;
1735
3.90k
1736
3.90k
  hull = proto_hull(set, &is_hull);
1737
3.90k
  if (
hull && 3.90k
!is_hull3.90k
)
{3.90k
1738
3.90k
    if (hull->n_ineq == 0)
1739
3.88k
      hull = initial_hull(hull, set);
1740
3.90k
    hull = extend(hull, set);
1741
3.90k
  }
1742
3.90k
  isl_set_free(set);
1743
3.90k
1744
3.90k
  return hull;
1745
3.90k
}
1746
1747
/* Compute the convex hull of a set without any parameters or
1748
 * integer divisions.  Depending on whether the set is bounded,
1749
 * we pass control to the wrapping based convex hull or
1750
 * the Fourier-Motzkin elimination based convex hull.
1751
 * We also handle a few special cases before checking the boundedness.
1752
 */
1753
static __isl_give isl_basic_set *uset_convex_hull(__isl_take isl_set *set)
1754
59
{
1755
59
  isl_bool bounded;
1756
59
  struct isl_basic_set *convex_hull = NULL;
1757
59
  struct isl_basic_set *lin;
1758
59
1759
59
  if (isl_set_n_dim(set) == 0)
1760
0
    return convex_hull_0d(set);
1761
59
1762
59
  set = isl_set_coalesce(set);
1763
59
  set = isl_set_set_rational(set);
1764
59
1765
59
  if (!set)
1766
0
    return NULL;
1767
59
  
if (59
set->n == 159
)
{14
1768
14
    convex_hull = isl_basic_set_copy(set->p[0]);
1769
14
    isl_set_free(set);
1770
14
    return convex_hull;
1771
59
  }
1772
45
  
if (45
isl_set_n_dim(set) == 145
)
1773
2
    return convex_hull_1d(set);
1774
45
1775
45
  bounded = isl_set_is_bounded(set);
1776
43
  if (bounded < 0)
1777
0
    goto error;
1778
43
  
if (43
bounded && 43
set->ctx->opt->convex == 14
ISL_CONVEX_HULL_WRAP14
)
1779
11
    return uset_convex_hull_wrap(set);
1780
43
1781
43
  lin = isl_set_combined_lineality_space(isl_set_copy(set));
1782
32
  if (!lin)
1783
0
    goto error;
1784
32
  
if (32
isl_basic_set_plain_is_universe(lin)32
)
{0
1785
0
    isl_set_free(set);
1786
0
    return lin;
1787
32
  }
1788
32
  
if (32
lin->n_eq < isl_basic_set_total_dim(lin)32
)
1789
8
    return modulo_lineality(set, lin);
1790
32
  isl_basic_set_free(lin);
1791
24
1792
32
  return uset_convex_hull_unbounded(set);
1793
32
error:
1794
0
  isl_set_free(set);
1795
0
  isl_basic_set_free(convex_hull);
1796
32
  return NULL;
1797
59
}
1798
1799
/* This is the core procedure, where "set" is a "pure" set, i.e.,
1800
 * without parameters or divs and where the convex hull of set is
1801
 * known to be full-dimensional.
1802
 */
1803
static __isl_give isl_basic_set *uset_convex_hull_wrap_bounded(
1804
  __isl_take isl_set *set)
1805
13.5k
{
1806
13.5k
  struct isl_basic_set *convex_hull = NULL;
1807
13.5k
1808
13.5k
  if (!set)
1809
0
    goto error;
1810
13.5k
1811
13.5k
  
if (13.5k
isl_set_n_dim(set) == 013.5k
)
{0
1812
0
    convex_hull = isl_basic_set_universe(isl_space_copy(set->dim));
1813
0
    isl_set_free(set);
1814
0
    convex_hull = isl_basic_set_set_rational(convex_hull);
1815
0
    return convex_hull;
1816
13.5k
  }
1817
13.5k
1818
13.5k
  set = isl_set_set_rational(set);
1819
13.5k
  set = isl_set_coalesce(set);
1820
13.5k
  if (!set)
1821
0
    goto error;
1822
13.5k
  
if (13.5k
set->n == 113.5k
)
{5.61k
1823
5.61k
    convex_hull = isl_basic_set_copy(set->p[0]);
1824
5.61k
    isl_set_free(set);
1825
5.61k
    convex_hull = isl_basic_map_remove_redundancies(convex_hull);
1826
5.61k
    return convex_hull;
1827
13.5k
  }
1828
7.96k
  
if (7.96k
isl_set_n_dim(set) == 17.96k
)
1829
4.06k
    return convex_hull_1d(set);
1830
7.96k
1831
3.89k
  return uset_convex_hull_wrap(set);
1832
7.96k
error:
1833
0
  isl_set_free(set);
1834
7.96k
  return NULL;
1835
13.5k
}
1836
1837
/* Compute the convex hull of set "set" with affine hull "affine_hull",
1838
 * We first remove the equalities (transforming the set), compute the
1839
 * convex hull of the transformed set and then add the equalities back
1840
 * (after performing the inverse transformation.
1841
 */
1842
static __isl_give isl_basic_set *modulo_affine_hull(
1843
  __isl_take isl_set *set, __isl_take isl_basic_set *affine_hull)
1844
3
{
1845
3
  struct isl_mat *T;
1846
3
  struct isl_mat *T2;
1847
3
  struct isl_basic_set *dummy;
1848
3
  struct isl_basic_set *convex_hull;
1849
3
1850
3
  dummy = isl_basic_set_remove_equalities(
1851
3
      isl_basic_set_copy(affine_hull), &T, &T2);
1852
3
  if (!dummy)
1853
0
    goto error;
1854
3
  isl_basic_set_free(dummy);
1855
3
  set = isl_set_preimage(set, T);
1856
3
  convex_hull = uset_convex_hull(set);
1857
3
  convex_hull = isl_basic_set_preimage(convex_hull, T2);
1858
3
  convex_hull = isl_basic_set_intersect(convex_hull, affine_hull);
1859
3
  return convex_hull;
1860
3
error:
1861
0
  isl_mat_free(T);
1862
0
  isl_mat_free(T2);
1863
0
  isl_basic_set_free(affine_hull);
1864
0
  isl_set_free(set);
1865
3
  return NULL;
1866
3
}
1867
1868
/* Return an empty basic map living in the same space as "map".
1869
 */
1870
static __isl_give isl_basic_map *replace_map_by_empty_basic_map(
1871
  __isl_take isl_map *map)
1872
794
{
1873
794
  isl_space *space;
1874
794
1875
794
  space = isl_map_get_space(map);
1876
794
  isl_map_free(map);
1877
794
  return isl_basic_map_empty(space);
1878
794
}
1879
1880
/* Compute the convex hull of a map.
1881
 *
1882
 * The implementation was inspired by "Extended Convex Hull" by Fukuda et al.,
1883
 * specifically, the wrapping of facets to obtain new facets.
1884
 */
1885
__isl_give isl_basic_map *isl_map_convex_hull(__isl_take isl_map *map)
1886
40
{
1887
40
  struct isl_basic_set *bset;
1888
40
  struct isl_basic_map *model = NULL;
1889
40
  struct isl_basic_set *affine_hull = NULL;
1890
40
  struct isl_basic_map *convex_hull = NULL;
1891
40
  struct isl_set *set = NULL;
1892
40
1893
40
  map = isl_map_detect_equalities(map);
1894
40
  map = isl_map_align_divs_internal(map);
1895
40
  if (!map)
1896
0
    goto error;
1897
40
1898
40
  
if (40
map->n == 040
)
1899
2
    return replace_map_by_empty_basic_map(map);
1900
40
1901
40
  model = isl_basic_map_copy(map->p[0]);
1902
38
  set = isl_map_underlying_set(map);
1903
38
  if (!set)
1904
0
    goto error;
1905
38
1906
38
  affine_hull = isl_set_affine_hull(isl_set_copy(set));
1907
38
  if (!affine_hull)
1908
0
    goto error;
1909
38
  
if (38
affine_hull->n_eq != 038
)
1910
2
    bset = modulo_affine_hull(set, affine_hull);
1911
38
  else {
1912
36
    isl_basic_set_free(affine_hull);
1913
36
    bset = uset_convex_hull(set);
1914
38
  }
1915
38
1916
38
  convex_hull = isl_basic_map_overlying_set(bset, model);
1917
38
  if (!convex_hull)
1918
0
    return NULL;
1919
38
1920
38
  
ISL_F_SET38
(convex_hull, ISL_BASIC_MAP_NO_IMPLICIT);38
1921
38
  ISL_F_SET(convex_hull, ISL_BASIC_MAP_ALL_EQUALITIES);
1922
38
  ISL_F_CLR(convex_hull, ISL_BASIC_MAP_RATIONAL);
1923
38
  return convex_hull;
1924
38
error:
1925
0
  isl_set_free(set);
1926
0
  isl_basic_map_free(model);
1927
38
  return NULL;
1928
40
}
1929
1930
struct isl_basic_set *isl_set_convex_hull(struct isl_set *set)
1931
40
{
1932
40
  return bset_from_bmap(isl_map_convex_hull(set_to_map(set)));
1933
40
}
1934
1935
__isl_give isl_basic_map *isl_map_polyhedral_hull(__isl_take isl_map *map)
1936
0
{
1937
0
  isl_basic_map *hull;
1938
0
1939
0
  hull = isl_map_convex_hull(map);
1940
0
  return isl_basic_map_remove_divs(hull);
1941
0
}
1942
1943
__isl_give isl_basic_set *isl_set_polyhedral_hull(__isl_take isl_set *set)
1944
0
{
1945
0
  return bset_from_bmap(isl_map_polyhedral_hull(set_to_map(set)));
1946
0
}
1947
1948
struct sh_data_entry {
1949
  struct isl_hash_table *table;
1950
  struct isl_tab    *tab;
1951
};
1952
1953
/* Holds the data needed during the simple hull computation.
1954
 * In particular,
1955
 *  n   the number of basic sets in the original set
1956
 *  hull_table  a hash table of already computed constraints
1957
 *      in the simple hull
1958
 *  p   for each basic set,
1959
 *    table   a hash table of the constraints
1960
 *    tab   the tableau corresponding to the basic set
1961
 */
1962
struct sh_data {
1963
  struct isl_ctx    *ctx;
1964
  unsigned    n;
1965
  struct isl_hash_table *hull_table;
1966
  struct sh_data_entry  p[1];
1967
};
1968
1969
static void sh_data_free(struct sh_data *data)
1970
2.79k
{
1971
2.79k
  int i;
1972
2.79k
1973
2.79k
  if (!data)
1974
0
    return;
1975
2.79k
  isl_hash_table_free(data->ctx, data->hull_table);
1976
9.26k
  for (i = 0; 
i < data->n9.26k
;
++i6.46k
)
{6.46k
1977
6.46k
    isl_hash_table_free(data->ctx, data->p[i].table);
1978
6.46k
    isl_tab_free(data->p[i].tab);
1979
6.46k
  }
1980
2.79k
  free(data);
1981
2.79k
}
1982
1983
struct ineq_cmp_data {
1984
  unsigned  len;
1985
  isl_int   *p;
1986
};
1987
1988
static int has_ineq(const void *entry, const void *val)
1989
68.2k
{
1990
68.2k
  isl_int *row = (isl_int *)entry;
1991
68.2k
  struct ineq_cmp_data *v = (struct ineq_cmp_data *)val;
1992
68.2k
1993
68.2k
  return isl_seq_eq(row + 1, v->p + 1, v->len) ||
1994
2.47k
         isl_seq_is_neg(row + 1, v->p + 1, v->len);
1995
68.2k
}
1996
1997
static int hash_ineq(struct isl_ctx *ctx, struct isl_hash_table *table,
1998
      isl_int *ineq, unsigned len)
1999
56.9k
{
2000
56.9k
  uint32_t c_hash;
2001
56.9k
  struct ineq_cmp_data v;
2002
56.9k
  struct isl_hash_table_entry *entry;
2003
56.9k
2004
56.9k
  v.len = len;
2005
56.9k
  v.p = ineq;
2006
56.9k
  c_hash = isl_seq_get_hash(ineq + 1, len);
2007
56.9k
  entry = isl_hash_table_find(ctx, table, c_hash, has_ineq, &v, 1);
2008
56.9k
  if (!entry)
2009
0
    return - 1;
2010
56.9k
  entry->data = ineq;
2011
56.9k
  return 0;
2012
56.9k
}
2013
2014
/* Fill hash table "table" with the constraints of "bset".
2015
 * Equalities are added as two inequalities.
2016
 * The value in the hash table is a pointer to the (in)equality of "bset".
2017
 */
2018
static int hash_basic_set(struct isl_hash_table *table,
2019
  __isl_keep isl_basic_set *bset)
2020
6.46k
{
2021
6.46k
  int i, j;
2022
6.46k
  unsigned dim = isl_basic_set_total_dim(bset);
2023
6.46k
2024
9.01k
  for (i = 0; 
i < bset->n_eq9.01k
;
++i2.54k
)
{2.54k
2025
7.63k
    for (j = 0; 
j < 27.63k
;
++j5.09k
)
{5.09k
2026
5.09k
      isl_seq_neg(bset->eq[i], bset->eq[i], 1 + dim);
2027
5.09k
      if (hash_ineq(bset->ctx, table, bset->eq[i], dim) < 0)
2028
0
        return -1;
2029
5.09k
    }
2030
6.46k
  }
2031
58.3k
  
for (i = 0; 6.46k
i < bset->n_ineq58.3k
;
++i51.8k
)
{51.8k
2032
51.8k
    if (hash_ineq(bset->ctx, table, bset->ineq[i], dim) < 0)
2033
0
      return -1;
2034
51.8k
  }
2035
6.46k
  return 0;
2036
6.46k
}
2037
2038
static struct sh_data *sh_data_alloc(__isl_keep isl_set *set, unsigned n_ineq)
2039
2.79k
{
2040
2.79k
  struct sh_data *data;
2041
2.79k
  int i;
2042
2.79k
2043
2.79k
  data = isl_calloc(set->ctx, struct sh_data,
2044
2.79k
    sizeof(struct sh_data) +
2045
2.79k
    (set->n - 1) * sizeof(struct sh_data_entry));
2046
2.79k
  if (!data)
2047
0
    return NULL;
2048
2.79k
  data->ctx = set->ctx;
2049
2.79k
  data->n = set->n;
2050
2.79k
  data->hull_table = isl_hash_table_alloc(set->ctx, n_ineq);
2051
2.79k
  if (!data->hull_table)
2052
0
    goto error;
2053
9.26k
  
for (i = 0; 2.79k
i < set->n9.26k
;
++i6.46k
)
{6.46k
2054
6.46k
    data->p[i].table = isl_hash_table_alloc(set->ctx,
2055
6.46k
            2 * set->p[i]->n_eq + set->p[i]->n_ineq);
2056
6.46k
    if (!data->p[i].table)
2057
0
      goto error;
2058
6.46k
    
if (6.46k
hash_basic_set(data->p[i].table, set->p[i]) < 06.46k
)
2059
0
      goto error;
2060
6.46k
  }
2061
2.79k
  return data;
2062
2.79k
error:
2063
0
  sh_data_free(data);
2064
2.79k
  return NULL;
2065
2.79k
}
2066
2067
/* Check if inequality "ineq" is a bound for basic set "j" or if
2068
 * it can be relaxed (by increasing the constant term) to become
2069
 * a bound for that basic set.  In the latter case, the constant
2070
 * term is updated.
2071
 * Relaxation of the constant term is only allowed if "shift" is set.
2072
 *
2073
 * Return 1 if "ineq" is a bound
2074
 *    0 if "ineq" may attain arbitrarily small values on basic set "j"
2075
 *   -1 if some error occurred
2076
 */
2077
static int is_bound(struct sh_data *data, __isl_keep isl_set *set, int j,
2078
  isl_int *ineq, int shift)
2079
4.39k
{
2080
4.39k
  enum isl_lp_result res;
2081
4.39k
  isl_int opt;
2082
4.39k
2083
4.39k
  if (
!data->p[j].tab4.39k
)
{2.62k
2084
2.62k
    data->p[j].tab = isl_tab_from_basic_set(set->p[j], 0);
2085
2.62k
    if (!data->p[j].tab)
2086
0
      return -1;
2087
4.39k
  }
2088
4.39k
2089
4.39k
  
isl_int_init4.39k
(opt);4.39k
2090
4.39k
2091
4.39k
  res = isl_tab_min(data->p[j].tab, ineq, data->ctx->one,
2092
4.39k
        &opt, NULL, 0);
2093
4.39k
  if (
res == isl_lp_ok && 4.39k
isl_int_is_neg802
(opt))
{268
2094
268
    if (shift)
2095
71
      isl_int_sub(ineq[0], ineq[0], opt);
2096
268
    else
2097
197
      res = isl_lp_unbounded;
2098
4.39k
  }
2099
4.39k
2100
4.39k
  isl_int_clear(opt);
2101
4.39k
2102
4.39k
  return (res == isl_lp_ok || 
res == isl_lp_empty3.78k
) ?
1605
:
2103
3.78k
         
res == isl_lp_unbounded ? 3.78k
03.78k
:
-10
;
2104
4.39k
}
2105
2106
/* Set the constant term of "ineq" to the maximum of those of the constraints
2107
 * in the basic sets of "set" following "i" that are parallel to "ineq".
2108
 * That is, if any of the basic sets of "set" following "i" have a more
2109
 * relaxed copy of "ineq", then replace "ineq" by the most relaxed copy.
2110
 * "c_hash" is the hash value of the linear part of "ineq".
2111
 * "v" has been set up for use by has_ineq.
2112
 *
2113
 * Note that the two inequality constraints corresponding to an equality are
2114
 * represented by the same inequality constraint in data->p[j].table
2115
 * (but with different hash values).  This means the constraint (or at
2116
 * least its constant term) may need to be temporarily negated to get
2117
 * the actually hashed constraint.
2118
 */
2119
static void set_max_constant_term(struct sh_data *data, __isl_keep isl_set *set,
2120
  int i, isl_int *ineq, uint32_t c_hash, struct ineq_cmp_data *v)
2121
13.4k
{
2122
13.4k
  int j;
2123
13.4k
  isl_ctx *ctx;
2124
13.4k
  struct isl_hash_table_entry *entry;
2125
13.4k
2126
13.4k
  ctx = isl_set_get_ctx(set);
2127
26.7k
  for (j = i + 1; 
j < set->n26.7k
;
++j13.2k
)
{13.2k
2128
13.2k
    int neg;
2129
13.2k
    isl_int *ineq_j;
2130
13.2k
2131
13.2k
    entry = isl_hash_table_find(ctx, data->p[j].table,
2132
13.2k
            c_hash, &has_ineq, v, 0);
2133
13.2k
    if (!entry)
2134
1.31k
      continue;
2135
13.2k
2136
13.2k
    ineq_j = entry->data;
2137
11.9k
    neg = isl_seq_is_neg(ineq_j + 1, ineq + 1, v->len);
2138
11.9k
    if (neg)
2139
1.22k
      isl_int_neg(ineq_j[0], ineq_j[0]);
2140
11.9k
    if (isl_int_gt(ineq_j[0], ineq[0]))
2141
362
      isl_int_set(ineq[0], ineq_j[0]);
2142
11.9k
    if (neg)
2143
1.22k
      isl_int_neg(ineq_j[0], ineq_j[0]);
2144
13.4k
  }
2145
13.4k
}
2146
2147
/* Check if inequality "ineq" from basic set "i" is or can be relaxed to
2148
 * become a bound on the whole set.  If so, add the (relaxed) inequality
2149
 * to "hull".  Relaxation is only allowed if "shift" is set.
2150
 *
2151
 * We first check if "hull" already contains a translate of the inequality.
2152
 * If so, we are done.
2153
 * Then, we check if any of the previous basic sets contains a translate
2154
 * of the inequality.  If so, then we have already considered this
2155
 * inequality and we are done.
2156
 * Otherwise, for each basic set other than "i", we check if the inequality
2157
 * is a bound on the basic set, but first replace the constant term
2158
 * by the maximal value of any translate of the inequality in any
2159
 * of the following basic sets.
2160
 * For previous basic sets, we know that they do not contain a translate
2161
 * of the inequality, so we directly call is_bound.
2162
 * For following basic sets, we first check if a translate of the
2163
 * inequality appears in its description.  If so, the constant term
2164
 * of the inequality has already been updated with respect to this
2165
 * translate and the inequality is therefore known to be a bound
2166
 * of this basic set.
2167
 */
2168
static __isl_give isl_basic_set *add_bound(__isl_take isl_basic_set *hull,
2169
  struct sh_data *data, __isl_keep isl_set *set, int i, isl_int *ineq,
2170
  int shift)
2171
25.3k
{
2172
25.3k
  uint32_t c_hash;
2173
25.3k
  struct ineq_cmp_data v;
2174
25.3k
  struct isl_hash_table_entry *entry;
2175
25.3k
  int j, k;
2176
25.3k
2177
25.3k
  if (!hull)
2178
0
    return NULL;
2179
25.3k
2180
25.3k
  v.len = isl_basic_set_total_dim(hull);
2181
25.3k
  v.p = ineq;
2182
25.3k
  c_hash = isl_seq_get_hash(ineq + 1, v.len);
2183
25.3k
2184
25.3k
  entry = isl_hash_table_find(hull->ctx, data->hull_table, c_hash,
2185
25.3k
          has_ineq, &v, 0);
2186
25.3k
  if (entry)
2187
11.8k
    return hull;
2188
25.3k
2189
14.9k
  
for (j = 0; 13.5k
j < i14.9k
;
++j1.33k
)
{1.43k
2190
1.43k
    entry = isl_hash_table_find(hull->ctx, data->p[j].table,
2191
1.43k
            c_hash, has_ineq, &v, 0);
2192
1.43k
    if (entry)
2193
95
      break;
2194
13.5k
  }
2195
13.5k
  if (j < i)
2196
95
    return hull;
2197
13.5k
2198
13.5k
  k = isl_basic_set_alloc_inequality(hull);
2199
13.4k
  if (k < 0)
2200
0
    goto error;
2201
13.4k
  isl_seq_cpy(hull->ineq[k], ineq, 1 + v.len);
2202
13.4k
2203
13.4k
  set_max_constant_term(data, set, i, hull->ineq[k], c_hash, &v);
2204
13.6k
  for (j = 0; 
j < i13.6k
;
++j167
)
{1.16k
2205
1.16k
    int bound;
2206
1.16k
    bound = is_bound(data, set, j, hull->ineq[k], shift);
2207
1.16k
    if (bound < 0)
2208
0
      goto error;
2209
1.16k
    
if (1.16k
!bound1.16k
)
2210
993
      break;
2211
13.4k
  }
2212
13.4k
  
if (13.4k
j < i13.4k
)
{993
2213
993
    isl_basic_set_free_inequality(hull, 1);
2214
993
    return hull;
2215
13.4k
  }
2216
13.4k
2217
24.5k
  
for (j = i + 1; 12.4k
j < set->n24.5k
;
++j12.0k
)
{12.9k
2218
12.9k
    int bound;
2219
12.9k
    entry = isl_hash_table_find(hull->ctx, data->p[j].table,
2220
12.9k
            c_hash, has_ineq, &v, 0);
2221
12.9k
    if (entry)
2222
11.8k
      continue;
2223
12.9k
    bound = is_bound(data, set, j, hull->ineq[k], shift);
2224
1.17k
    if (bound < 0)
2225
0
      goto error;
2226
1.17k
    
if (1.17k
!bound1.17k
)
2227
960
      break;
2228
12.4k
  }
2229
12.4k
  
if (12.4k
j < set->n12.4k
)
{960
2230
960
    isl_basic_set_free_inequality(hull, 1);
2231
960
    return hull;
2232
12.4k
  }
2233
12.4k
2234
12.4k
  entry = isl_hash_table_find(hull->ctx, data->hull_table, c_hash,
2235
11.5k
          has_ineq, &v, 1);
2236
11.5k
  if (!entry)
2237
0
    goto error;
2238
11.5k
  entry->data = hull->ineq[k];
2239
11.5k
2240
11.5k
  return hull;
2241
11.5k
error:
2242
0
  isl_basic_set_free(hull);
2243
11.5k
  return NULL;
2244
25.3k
}
2245
2246
/* Check if any inequality from basic set "i" is or can be relaxed to
2247
 * become a bound on the whole set.  If so, add the (relaxed) inequality
2248
 * to "hull".  Relaxation is only allowed if "shift" is set.
2249
 */
2250
static __isl_give isl_basic_set *add_bounds(__isl_take isl_basic_set *bset,
2251
  struct sh_data *data, __isl_keep isl_set *set, int i, int shift)
2252
4.77k
{
2253
4.77k
  int j, k;
2254
4.77k
  unsigned dim = isl_basic_set_total_dim(bset);
2255
4.77k
2256
7.14k
  for (j = 0; 
j < set->p[i]->n_eq7.14k
;
++j2.36k
)
{2.36k
2257
7.09k
    for (k = 0; 
k < 27.09k
;
++k4.73k
)
{4.73k
2258
4.73k
      isl_seq_neg(set->p[i]->eq[j], set->p[i]->eq[j], 1+dim);
2259
4.73k
      bset = add_bound(bset, data, set, i, set->p[i]->eq[j],
2260
4.73k
              shift);
2261
4.73k
    }
2262
4.77k
  }
2263
25.4k
  for (j = 0; 
j < set->p[i]->n_ineq25.4k
;
++j20.6k
)
2264
20.6k
    bset = add_bound(bset, data, set, i, set->p[i]->ineq[j], shift);
2265
4.77k
  return bset;
2266
4.77k
}
2267
2268
/* Compute a superset of the convex hull of set that is described
2269
 * by only (translates of) the constraints in the constituents of set.
2270
 * Translation is only allowed if "shift" is set.
2271
 */
2272
static __isl_give isl_basic_set *uset_simple_hull(__isl_take isl_set *set,
2273
  int shift)
2274
2.32k
{
2275
2.32k
  struct sh_data *data = NULL;
2276
2.32k
  struct isl_basic_set *hull = NULL;
2277
2.32k
  unsigned n_ineq;
2278
2.32k
  int i;
2279
2.32k
2280
2.32k
  if (!set)
2281
0
    return NULL;
2282
2.32k
2283
2.32k
  n_ineq = 0;
2284
7.10k
  for (i = 0; 
i < set->n7.10k
;
++i4.77k
)
{4.77k
2285
4.77k
    if (!set->p[i])
2286
0
      goto error;
2287
4.77k
    n_ineq += 2 * set->p[i]->n_eq + set->p[i]->n_ineq;
2288
4.77k
  }
2289
2.32k
2290
2.32k
  hull = isl_basic_set_alloc_space(isl_space_copy(set->dim), 0, 0, n_ineq);
2291
2.32k
  if (!hull)
2292
0
    goto error;
2293
2.32k
2294
2.32k
  data = sh_data_alloc(set, n_ineq);
2295
2.32k
  if (!data)
2296
0
    goto error;
2297
2.32k
2298
7.10k
  
for (i = 0; 2.32k
i < set->n7.10k
;
++i4.77k
)
2299
4.77k
    hull = add_bounds(hull, data, set, i, shift);
2300
2.32k
2301
2.32k
  sh_data_free(data);
2302
2.32k
  isl_set_free(set);
2303
2.32k
2304
2.32k
  return hull;
2305
2.32k
error:
2306
0
  sh_data_free(data);
2307
0
  isl_basic_set_free(hull);
2308
0
  isl_set_free(set);
2309
2.32k
  return NULL;
2310
2.32k
}
2311
2312
/* Compute a superset of the convex hull of map that is described
2313
 * by only (translates of) the constraints in the constituents of map.
2314
 * Handle trivial cases where map is NULL or contains at most one disjunct.
2315
 */
2316
static __isl_give isl_basic_map *map_simple_hull_trivial(
2317
  __isl_take isl_map *map)
2318
35.9k
{
2319
35.9k
  isl_basic_map *hull;
2320
35.9k
2321
35.9k
  if (!map)
2322
0
    return NULL;
2323
35.9k
  
if (35.9k
map->n == 035.9k
)
2324
792
    return replace_map_by_empty_basic_map(map);
2325
35.9k
2326
35.9k
  hull = isl_basic_map_copy(map->p[0]);
2327
35.1k
  isl_map_free(map);
2328
35.9k
  return hull;
2329
35.9k
}
2330
2331
/* Return a copy of the simple hull cached inside "map".
2332
 * "shift" determines whether to return the cached unshifted or shifted
2333
 * simple hull.
2334
 */
2335
static __isl_give isl_basic_map *cached_simple_hull(__isl_take isl_map *map,
2336
  int shift)
2337
45
{
2338
45
  isl_basic_map *hull;
2339
45
2340
45
  hull = isl_basic_map_copy(map->cached_simple_hull[shift]);
2341
45
  isl_map_free(map);
2342
45
2343
45
  return hull;
2344
45
}
2345
2346
/* Compute a superset of the convex hull of map that is described
2347
 * by only (translates of) the constraints in the constituents of map.
2348
 * Translation is only allowed if "shift" is set.
2349
 *
2350
 * The constraints are sorted while removing redundant constraints
2351
 * in order to indicate a preference of which constraints should
2352
 * be preserved.  In particular, pairs of constraints that are
2353
 * sorted together are preferred to either both be preserved
2354
 * or both be removed.  The sorting is performed inside
2355
 * isl_basic_map_remove_redundancies.
2356
 *
2357
 * The result of the computation is stored in map->cached_simple_hull[shift]
2358
 * such that it can be reused in subsequent calls.  The cache is cleared
2359
 * whenever the map is modified (in isl_map_cow).
2360
 * Note that the results need to be stored in the input map for there
2361
 * to be any chance that they may get reused.  In particular, they
2362
 * are stored in a copy of the input map that is saved before
2363
 * the integer division alignment.
2364
 */
2365
static __isl_give isl_basic_map *map_simple_hull(__isl_take isl_map *map,
2366
  int shift)
2367
38.3k
{
2368
38.3k
  struct isl_set *set = NULL;
2369
38.3k
  struct isl_basic_map *model = NULL;
2370
38.3k
  struct isl_basic_map *hull;
2371
38.3k
  struct isl_basic_map *affine_hull;
2372
38.3k
  struct isl_basic_set *bset = NULL;
2373
38.3k
  isl_map *input;
2374
38.3k
2375
38.3k
  if (
!map || 38.3k
map->n <= 138.3k
)
2376
35.1k
    return map_simple_hull_trivial(map);
2377
38.3k
2378
3.16k
  
if (3.16k
map->cached_simple_hull[shift]3.16k
)
2379
45
    return cached_simple_hull(map, shift);
2380
3.16k
2381
3.16k
  map = isl_map_detect_equalities(map);
2382
3.11k
  if (
!map || 3.11k
map->n <= 13.11k
)
2383
792
    return map_simple_hull_trivial(map);
2384
3.11k
  affine_hull = isl_map_affine_hull(isl_map_copy(map));
2385
2.32k
  input = isl_map_copy(map);
2386
2.32k
  map = isl_map_align_divs_internal(map);
2387
2.32k
  model = map ? isl_basic_map_copy(map->p[0]) : NULL;
2388
2.32k
2389
2.32k
  set = isl_map_underlying_set(map);
2390
2.32k
2391
2.32k
  bset = uset_simple_hull(set, shift);
2392
2.32k
2393
2.32k
  hull = isl_basic_map_overlying_set(bset, model);
2394
2.32k
2395
2.32k
  hull = isl_basic_map_intersect(hull, affine_hull);
2396
2.32k
  hull = isl_basic_map_remove_redundancies(hull);
2397
2.32k
2398
2.32k
  if (
hull2.32k
)
{2.32k
2399
2.32k
    ISL_F_SET(hull, ISL_BASIC_MAP_NO_IMPLICIT);
2400
2.32k
    ISL_F_SET(hull, ISL_BASIC_MAP_ALL_EQUALITIES);
2401
2.32k
  }
2402
2.32k
2403
2.32k
  hull = isl_basic_map_finalize(hull);
2404
2.32k
  if (input)
2405
2.32k
    input->cached_simple_hull[shift] = isl_basic_map_copy(hull);
2406
2.32k
  isl_map_free(input);
2407
2.32k
2408
3.11k
  return hull;
2409
38.3k
}
2410
2411
/* Compute a superset of the convex hull of map that is described
2412
 * by only translates of the constraints in the constituents of map.
2413
 */
2414
__isl_give isl_basic_map *isl_map_simple_hull(__isl_take isl_map *map)
2415
25.7k
{
2416
25.7k
  return map_simple_hull(map, 1);
2417
25.7k
}
2418
2419
struct isl_basic_set *isl_set_simple_hull(struct isl_set *set)
2420
11.0k
{
2421
11.0k
  return bset_from_bmap(isl_map_simple_hull(set_to_map(set)));
2422
11.0k
}
2423
2424
/* Compute a superset of the convex hull of map that is described
2425
 * by only the constraints in the constituents of map.
2426
 */
2427
__isl_give isl_basic_map *isl_map_unshifted_simple_hull(
2428
  __isl_take isl_map *map)
2429
12.6k
{
2430
12.6k
  return map_simple_hull(map, 0);
2431
12.6k
}
2432
2433
__isl_give isl_basic_set *isl_set_unshifted_simple_hull(
2434
  __isl_take isl_set *set)
2435
7.47k
{
2436
7.47k
  return isl_map_unshifted_simple_hull(set);
2437
7.47k
}
2438
2439
/* Drop all inequalities from "bmap1" that do not also appear in "bmap2".
2440
 * A constraint that appears with different constant terms
2441
 * in "bmap1" and "bmap2" is also kept, with the least restrictive
2442
 * (i.e., greatest) constant term.
2443
 * "bmap1" and "bmap2" are assumed to have the same (known)
2444
 * integer divisions.
2445
 * The constraints of both "bmap1" and "bmap2" are assumed
2446
 * to have been sorted using isl_basic_map_sort_constraints.
2447
 *
2448
 * Run through the inequality constraints of "bmap1" and "bmap2"
2449
 * in sorted order.
2450
 * Each constraint of "bmap1" without a matching constraint in "bmap2"
2451
 * is removed.
2452
 * If a match is found, the constraint is kept.  If needed, the constant
2453
 * term of the constraint is adjusted.
2454
 */
2455
static __isl_give isl_basic_map *select_shared_inequalities(
2456
  __isl_take isl_basic_map *bmap1, __isl_keep isl_basic_map *bmap2)
2457
1.44k
{
2458
1.44k
  int i1, i2;
2459
1.44k
2460
1.44k
  bmap1 = isl_basic_map_cow(bmap1);
2461
1.44k
  if (
!bmap1 || 1.44k
!bmap21.44k
)
2462
0
    return isl_basic_map_free(bmap1);
2463
1.44k
2464
1.44k
  i1 = bmap1->n_ineq - 1;
2465
1.44k
  i2 = bmap2->n_ineq - 1;
2466
3.01k
  while (
bmap1 && 3.01k
i1 >= 03.01k
&&
i2 >= 01.72k
)
{1.56k
2467
1.56k
    int cmp;
2468
1.56k
2469
1.56k
    cmp = isl_basic_map_constraint_cmp(bmap1, bmap1->ineq[i1],
2470
1.56k
              bmap2->ineq[i2]);
2471
1.56k
    if (
cmp < 01.56k
)
{364
2472
364
      --i2;
2473
364
      continue;
2474
1.56k
    }
2475
1.20k
    
if (1.20k
cmp > 01.20k
)
{453
2476
453
      if (isl_basic_map_drop_inequality(bmap1, i1) < 0)
2477
0
        bmap1 = isl_basic_map_free(bmap1);
2478
453
      --i1;
2479
453
      continue;
2480
1.20k
    }
2481
748
    
if (748
isl_int_lt748
(bmap1->ineq[i1][0], bmap2->ineq[i2][0]))
2482
29
      isl_int_set(bmap1->ineq[i1][0], bmap2->ineq[i2][0]);
2483
748
    --i1;
2484
748
    --i2;
2485
1.44k
  }
2486
1.62k
  for (; 
i1 >= 01.62k
;
--i1181
)
2487
181
    
if (181
isl_basic_map_drop_inequality(bmap1, i1) < 0181
)
2488
0
      bmap1 = isl_basic_map_free(bmap1);
2489
1.44k
2490
1.44k
  return bmap1;
2491
1.44k
}
2492
2493
/* Drop all equalities from "bmap1" that do not also appear in "bmap2".
2494
 * "bmap1" and "bmap2" are assumed to have the same (known)
2495
 * integer divisions.
2496
 *
2497
 * Run through the equality constraints of "bmap1" and "bmap2".
2498
 * Each constraint of "bmap1" without a matching constraint in "bmap2"
2499
 * is removed.
2500
 */
2501
static __isl_give isl_basic_map *select_shared_equalities(
2502
  __isl_take isl_basic_map *bmap1, __isl_keep isl_basic_map *bmap2)
2503
1.44k
{
2504
1.44k
  int i1, i2;
2505
1.44k
  unsigned total;
2506
1.44k
2507
1.44k
  bmap1 = isl_basic_map_cow(bmap1);
2508
1.44k
  if (
!bmap1 || 1.44k
!bmap21.44k
)
2509
0
    return isl_basic_map_free(bmap1);
2510
1.44k
2511
1.44k
  total = isl_basic_map_total_dim(bmap1);
2512
1.44k
2513
1.44k
  i1 = bmap1->n_eq - 1;
2514
1.44k
  i2 = bmap2->n_eq - 1;
2515
2.21k
  while (
bmap1 && 2.21k
i1 >= 02.21k
&&
i2 >= 0774
)
{767
2516
767
    int last1, last2;
2517
767
2518
767
    last1 = isl_seq_last_non_zero(bmap1->eq[i1] + 1, total);
2519
767
    last2 = isl_seq_last_non_zero(bmap2->eq[i2] + 1, total);
2520
767
    if (
last1 > last2767
)
{3
2521
3
      --i2;
2522
3
      continue;
2523
767
    }
2524
764
    
if (764
last1 < last2764
)
{3
2525
3
      if (isl_basic_map_drop_equality(bmap1, i1) < 0)
2526
0
        bmap1 = isl_basic_map_free(bmap1);
2527
3
      --i1;
2528
3
      continue;
2529
764
    }
2530
761
    
if (761
!isl_seq_eq(bmap1->eq[i1], bmap2->eq[i2], 1 + total)761
)
{2
2531
2
      if (isl_basic_map_drop_equality(bmap1, i1) < 0)
2532
0
        bmap1 = isl_basic_map_free(bmap1);
2533
761
    }
2534
761
    --i1;
2535
761
    --i2;
2536
1.44k
  }
2537
1.45k
  for (; 
i1 >= 01.45k
;
--i17
)
2538
7
    
if (7
isl_basic_map_drop_equality(bmap1, i1) < 07
)
2539
0
      bmap1 = isl_basic_map_free(bmap1);
2540
1.44k
2541
1.44k
  return bmap1;
2542
1.44k
}
2543
2544
/* Compute a superset of "bmap1" and "bmap2" that is described
2545
 * by only the constraints that appear in both "bmap1" and "bmap2".
2546
 *
2547
 * First drop constraints that involve unknown integer divisions
2548
 * since it is not trivial to check whether two such integer divisions
2549
 * in different basic maps are the same.
2550
 * Then align the remaining (known) divs and sort the constraints.
2551
 * Finally drop all inequalities and equalities from "bmap1" that
2552
 * do not also appear in "bmap2".
2553
 */
2554
__isl_give isl_basic_map *isl_basic_map_plain_unshifted_simple_hull(
2555
  __isl_take isl_basic_map *bmap1, __isl_take isl_basic_map *bmap2)
2556
1.44k
{
2557
1.44k
  bmap1 = isl_basic_map_drop_constraint_involving_unknown_divs(bmap1);
2558
1.44k
  bmap2 = isl_basic_map_drop_constraint_involving_unknown_divs(bmap2);
2559
1.44k
  bmap2 = isl_basic_map_align_divs(bmap2, bmap1);
2560
1.44k
  bmap1 = isl_basic_map_align_divs(bmap1, bmap2);
2561
1.44k
  bmap1 = isl_basic_map_gauss(bmap1, NULL);
2562
1.44k
  bmap2 = isl_basic_map_gauss(bmap2, NULL);
2563
1.44k
  bmap1 = isl_basic_map_sort_constraints(bmap1);
2564
1.44k
  bmap2 = isl_basic_map_sort_constraints(bmap2);
2565
1.44k
2566
1.44k
  bmap1 = select_shared_inequalities(bmap1, bmap2);
2567
1.44k
  bmap1 = select_shared_equalities(bmap1, bmap2);
2568
1.44k
2569
1.44k
  isl_basic_map_free(bmap2);
2570
1.44k
  bmap1 = isl_basic_map_finalize(bmap1);
2571
1.44k
  return bmap1;
2572
1.44k
}
2573
2574
/* Compute a superset of the convex hull of "map" that is described
2575
 * by only the constraints in the constituents of "map".
2576
 * In particular, the result is composed of constraints that appear
2577
 * in each of the basic maps of "map"
2578
 *
2579
 * Constraints that involve unknown integer divisions are dropped
2580
 * since it is not trivial to check whether two such integer divisions
2581
 * in different basic maps are the same.
2582
 *
2583
 * The hull is initialized from the first basic map and then
2584
 * updated with respect to the other basic maps in turn.
2585
 */
2586
__isl_give isl_basic_map *isl_map_plain_unshifted_simple_hull(
2587
  __isl_take isl_map *map)
2588
1.28k
{
2589
1.28k
  int i;
2590
1.28k
  isl_basic_map *hull;
2591
1.28k
2592
1.28k
  if (!map)
2593
0
    return NULL;
2594
1.28k
  
if (1.28k
map->n <= 11.28k
)
2595
0
    return map_simple_hull_trivial(map);
2596
1.28k
  map = isl_map_drop_constraint_involving_unknown_divs(map);
2597
1.28k
  hull = isl_basic_map_copy(map->p[0]);
2598
2.73k
  for (i = 1; 
i < map->n2.73k
;
++i1.44k
)
{1.44k
2599
1.44k
    isl_basic_map *bmap_i;
2600
1.44k
2601
1.44k
    bmap_i = isl_basic_map_copy(map->p[i]);
2602
1.44k
    hull = isl_basic_map_plain_unshifted_simple_hull(hull, bmap_i);
2603
1.44k
  }
2604
1.28k
2605
1.28k
  isl_map_free(map);
2606
1.28k
  return hull;
2607
1.28k
}
2608
2609
/* Compute a superset of the convex hull of "set" that is described
2610
 * by only the constraints in the constituents of "set".
2611
 * In particular, the result is composed of constraints that appear
2612
 * in each of the basic sets of "set"
2613
 */
2614
__isl_give isl_basic_set *isl_set_plain_unshifted_simple_hull(
2615
  __isl_take isl_set *set)
2616
836
{
2617
836
  return isl_map_plain_unshifted_simple_hull(set);
2618
836
}
2619
2620
/* Check if "ineq" is a bound on "set" and, if so, add it to "hull".
2621
 *
2622
 * For each basic set in "set", we first check if the basic set
2623
 * contains a translate of "ineq".  If this translate is more relaxed,
2624
 * then we assume that "ineq" is not a bound on this basic set.
2625
 * Otherwise, we know that it is a bound.
2626
 * If the basic set does not contain a translate of "ineq", then
2627
 * we call is_bound to perform the test.
2628
 */
2629
static __isl_give isl_basic_set *add_bound_from_constraint(
2630
  __isl_take isl_basic_set *hull, struct sh_data *data,
2631
  __isl_keep isl_set *set, isl_int *ineq)
2632
7.38k
{
2633
7.38k
  int i, k;
2634
7.38k
  isl_ctx *ctx;
2635
7.38k
  uint32_t c_hash;
2636
7.38k
  struct ineq_cmp_data v;
2637
7.38k
2638
7.38k
  if (
!hull || 7.38k
!set7.38k
)
2639
0
    return isl_basic_set_free(hull);
2640
7.38k
2641
7.38k
  v.len = isl_basic_set_total_dim(hull);
2642
7.38k
  v.p = ineq;
2643
7.38k
  c_hash = isl_seq_get_hash(ineq + 1, v.len);
2644
7.38k
2645
7.38k
  ctx = isl_basic_set_get_ctx(hull);
2646
39.2k
  for (i = 0; 
i < set->n39.2k
;
++i31.8k
)
{34.7k
2647
34.7k
    int bound;
2648
34.7k
    struct isl_hash_table_entry *entry;
2649
34.7k
2650
34.7k
    entry = isl_hash_table_find(ctx, data->p[i].table,
2651
34.7k
            c_hash, &has_ineq, &v, 0);
2652
34.7k
    if (
entry34.7k
)
{32.6k
2653
32.6k
      isl_int *ineq_i = entry->data;
2654
32.6k
      int neg, more_relaxed;
2655
32.6k
2656
32.6k
      neg = isl_seq_is_neg(ineq_i + 1, ineq + 1, v.len);
2657
32.6k
      if (neg)
2658
60
        isl_int_neg(ineq_i[0], ineq_i[0]);
2659
32.6k
      more_relaxed = isl_int_gt(ineq_i[0], ineq[0]);
2660
32.6k
      if (neg)
2661
60
        isl_int_neg(ineq_i[0], ineq_i[0]);
2662
32.6k
      if (more_relaxed)
2663
1.05k
        break;
2664
32.6k
      else
2665
31.6k
        continue;
2666
34.7k
    }
2667
34.7k
    bound = is_bound(data, set, i, ineq, 0);
2668
2.05k
    if (bound < 0)
2669
0
      return isl_basic_set_free(hull);
2670
2.05k
    
if (2.05k
!bound2.05k
)
2671
1.83k
      break;
2672
7.38k
  }
2673
7.38k
  
if (7.38k
i < set->n7.38k
)
2674
2.88k
    return hull;
2675
7.38k
2676
7.38k
  k = isl_basic_set_alloc_inequality(hull);
2677
4.49k
  if (k < 0)
2678
0
    return isl_basic_set_free(hull);
2679
4.49k
  isl_seq_cpy(hull->ineq[k], ineq, 1 + v.len);
2680
4.49k
2681
4.49k
  return hull;
2682
7.38k
}
2683
2684
/* Compute a superset of the convex hull of "set" that is described
2685
 * by only some of the "n_ineq" constraints in the list "ineq", where "set"
2686
 * has no parameters or integer divisions.
2687
 *
2688
 * The inequalities in "ineq" are assumed to have been sorted such
2689
 * that constraints with the same linear part appear together and
2690
 * that among constraints with the same linear part, those with
2691
 * smaller constant term appear first.
2692
 *
2693
 * We reuse the same data structure that is used by uset_simple_hull,
2694
 * but we do not need the hull table since we will not consider the
2695
 * same constraint more than once.  We therefore allocate it with zero size.
2696
 *
2697
 * We run through the constraints and try to add them one by one,
2698
 * skipping identical constraints.  If we have added a constraint and
2699
 * the next constraint is a more relaxed translate, then we skip this
2700
 * next constraint as well.
2701
 */
2702
static __isl_give isl_basic_set *uset_unshifted_simple_hull_from_constraints(
2703
  __isl_take isl_set *set, int n_ineq, isl_int **ineq)
2704
472
{
2705
472
  int i;
2706
472
  int last_added = 0;
2707
472
  struct sh_data *data = NULL;
2708
472
  isl_basic_set *hull = NULL;
2709
472
  unsigned dim;
2710
472
2711
472
  hull = isl_basic_set_alloc_space(isl_set_get_space(set), 0, 0, n_ineq);
2712
472
  if (!hull)
2713
0
    goto error;
2714
472
2715
472
  data = sh_data_alloc(set, 0);
2716
472
  if (!data)
2717
0
    goto error;
2718
472
2719
472
  dim = isl_set_dim(set, isl_dim_set);
2720
33.9k
  for (i = 0; 
i < n_ineq33.9k
;
++i33.4k
)
{33.4k
2721
33.4k
    int hull_n_ineq = hull->n_ineq;
2722
33.4k
    int parallel;
2723
33.4k
2724
33.4k
    parallel = i > 0 && isl_seq_eq(ineq[i - 1] + 1, ineq[i] + 1,
2725
33.4k
            dim);
2726
33.4k
    if (parallel &&
2727
27.2k
        
(last_added || 27.2k
isl_int_eq3.92k
(ineq[i - 1][0], ineq[i][0])))
2728
26.1k
      continue;
2729
33.4k
    hull = add_bound_from_constraint(hull, data, set, ineq[i]);
2730
7.38k
    if (!hull)
2731
0
      goto error;
2732
7.38k
    last_added = hull->n_ineq > hull_n_ineq;
2733
7.38k
  }
2734
472
2735
472
  sh_data_free(data);
2736
472
  isl_set_free(set);
2737
472
  return hull;
2738
472
error:
2739
0
  sh_data_free(data);
2740
0
  isl_set_free(set);
2741
0
  isl_basic_set_free(hull);
2742
472
  return NULL;
2743
472
}
2744
2745
/* Collect pointers to all the inequalities in the elements of "list"
2746
 * in "ineq".  For equalities, store both a pointer to the equality and
2747
 * a pointer to its opposite, which is first copied to "mat".
2748
 * "ineq" and "mat" are assumed to have been preallocated to the right size
2749
 * (the number of inequalities + 2 times the number of equalites and
2750
 * the number of equalities, respectively).
2751
 */
2752
static __isl_give isl_mat *collect_inequalities(__isl_take isl_mat *mat,
2753
  __isl_keep isl_basic_set_list *list, isl_int **ineq)
2754
472
{
2755
472
  int i, j, n, n_eq, n_ineq;
2756
472
2757
472
  if (!mat)
2758
0
    return NULL;
2759
472
2760
472
  n_eq = 0;
2761
472
  n_ineq = 0;
2762
472
  n = isl_basic_set_list_n_basic_set(list);
2763
2.89k
  for (i = 0; 
i < n2.89k
;
++i2.42k
)
{2.42k
2764
2.42k
    isl_basic_set *bset;
2765
2.42k
    bset = isl_basic_set_list_get_basic_set(list, i);
2766
2.42k
    if (!bset)
2767
0
      return isl_mat_free(mat);
2768
3.17k
    
for (j = 0; 2.42k
j < bset->n_eq3.17k
;
++j759
)
{759
2769
759
      ineq[n_ineq++] = mat->row[n_eq];
2770
759
      ineq[n_ineq++] = bset->eq[j];
2771
759
      isl_seq_neg(mat->row[n_eq++], bset->eq[j], mat->n_col);
2772
2.42k
    }
2773
34.3k
    for (j = 0; 
j < bset->n_ineq34.3k
;
++j31.9k
)
2774
31.9k
      ineq[n_ineq++] = bset->ineq[j];
2775
2.42k
    isl_basic_set_free(bset);
2776
2.42k
  }
2777
472
2778
472
  return mat;
2779
472
}
2780
2781
/* Comparison routine for use as an isl_sort callback.
2782
 *
2783
 * Constraints with the same linear part are sorted together and
2784
 * among constraints with the same linear part, those with smaller
2785
 * constant term are sorted first.
2786
 */
2787
static int cmp_ineq(const void *a, const void *b, void *arg)
2788
208k
{
2789
208k
  unsigned dim = *(unsigned *) arg;
2790
208k
  isl_int * const *ineq1 = a;
2791
208k
  isl_int * const *ineq2 = b;
2792
208k
  int cmp;
2793
208k
2794
208k
  cmp = isl_seq_cmp((*ineq1) + 1, (*ineq2) + 1, dim);
2795
208k
  if (cmp != 0)
2796
160k
    return cmp;
2797
47.8k
  
return 47.8k
isl_int_cmp47.8k
((*ineq1)[0], (*ineq2)[0]);
2798
208k
}
2799
2800
/* Compute a superset of the convex hull of "set" that is described
2801
 * by only constraints in the elements of "list", where "set" has
2802
 * no parameters or integer divisions.
2803
 *
2804
 * We collect all the constraints in those elements and then
2805
 * sort the constraints such that constraints with the same linear part
2806
 * are sorted together and that those with smaller constant term are
2807
 * sorted first.
2808
 */
2809
static __isl_give isl_basic_set *uset_unshifted_simple_hull_from_basic_set_list(
2810
  __isl_take isl_set *set, __isl_take isl_basic_set_list *list)
2811
472
{
2812
472
  int i, n, n_eq, n_ineq;
2813
472
  unsigned dim;
2814
472
  isl_ctx *ctx;
2815
472
  isl_mat *mat = NULL;
2816
472
  isl_int **ineq = NULL;
2817
472
  isl_basic_set *hull;
2818
472
2819
472
  if (!set)
2820
0
    goto error;
2821
472
  ctx = isl_set_get_ctx(set);
2822
472
2823
472
  n_eq = 0;
2824
472
  n_ineq = 0;
2825
472
  n = isl_basic_set_list_n_basic_set(list);
2826
2.89k
  for (i = 0; 
i < n2.89k
;
++i2.42k
)
{2.42k
2827
2.42k
    isl_basic_set *bset;
2828
2.42k
    bset = isl_basic_set_list_get_basic_set(list, i);
2829
2.42k
    if (!bset)
2830
0
      goto error;
2831
2.42k
    n_eq += bset->n_eq;
2832
2.42k
    n_ineq += 2 * bset->n_eq + bset->n_ineq;
2833
2.42k
    isl_basic_set_free(bset);
2834
2.42k
  }
2835
472
2836
472
  
ineq = 472
isl_alloc_array472
(ctx, isl_int *, n_ineq);
2837
472
  if (
n_ineq > 0 && 472
!ineq472
)
2838
0
    goto error;
2839
472
2840
472
  dim = isl_set_dim(set, isl_dim_set);
2841
472
  mat = isl_mat_alloc(ctx, n_eq, 1 + dim);
2842
472
  mat = collect_inequalities(mat, list, ineq);
2843
472
  if (!mat)
2844
0
    goto error;
2845
472
2846
472
  
if (472
isl_sort(ineq, n_ineq, sizeof(ineq[0]), &cmp_ineq, &dim) < 0472
)
2847
0
    goto error;
2848
472
2849
472
  hull = uset_unshifted_simple_hull_from_constraints(set, n_ineq, ineq);
2850
472
2851
472
  isl_mat_free(mat);
2852
472
  free(ineq);
2853
472
  isl_basic_set_list_free(list);
2854
472
  return hull;
2855
472
error:
2856
0
  isl_mat_free(mat);
2857
0
  free(ineq);
2858
0
  isl_set_free(set);
2859
0
  isl_basic_set_list_free(list);
2860
472
  return NULL;
2861
472
}
2862
2863
/* Compute a superset of the convex hull of "map" that is described
2864
 * by only constraints in the elements of "list".
2865
 *
2866
 * If the list is empty, then we can only describe the universe set.
2867
 * If the input map is empty, then all constraints are valid, so
2868
 * we return the intersection of the elements in "list".
2869
 *
2870
 * Otherwise, we align all divs and temporarily treat them
2871
 * as regular variables, computing the unshifted simple hull in
2872
 * uset_unshifted_simple_hull_from_basic_set_list.
2873
 */
2874
static __isl_give isl_basic_map *map_unshifted_simple_hull_from_basic_map_list(
2875
  __isl_take isl_map *map, __isl_take isl_basic_map_list *list)
2876
472
{
2877
472
  isl_basic_map *model;
2878
472
  isl_basic_map *hull;
2879
472
  isl_set *set;
2880
472
  isl_basic_set_list *bset_list;
2881
472
2882
472
  if (
!map || 472
!list472
)
2883
0
    goto error;
2884
472
2885
472
  
if (472
isl_basic_map_list_n_basic_map(list) == 0472
)
{0
2886
0
    isl_space *space;
2887
0
2888
0
    space = isl_map_get_space(map);
2889
0
    isl_map_free(map);
2890
0
    isl_basic_map_list_free(list);
2891
0
    return isl_basic_map_universe(space);
2892
472
  }
2893
472
  
if (472
isl_map_plain_is_empty(map)472
)
{0
2894
0
    isl_map_free(map);
2895
0
    return isl_basic_map_list_intersect(list);
2896
472
  }
2897
472
2898
472
  map = isl_map_align_divs_to_basic_map_list(map, list);
2899
472
  if (!map)
2900
0
    goto error;
2901
472
  list = isl_basic_map_list_align_divs_to_basic_map(list, map->p[0]);
2902
472
2903
472
  model = isl_basic_map_list_get_basic_map(list, 0);
2904
472
2905
472
  set = isl_map_underlying_set(map);
2906
472
  bset_list = isl_basic_map_list_underlying_set(list);
2907
472
2908
472
  hull = uset_unshifted_simple_hull_from_basic_set_list(set, bset_list);
2909
472
  hull = isl_basic_map_overlying_set(hull, model);
2910
472
2911
472
  return hull;
2912
472
error:
2913
0
  isl_map_free(map);
2914
0
  isl_basic_map_list_free(list);
2915
472
  return NULL;
2916
472
}
2917
2918
/* Return a sequence of the basic maps that make up the maps in "list".
2919
 */
2920
static __isl_give isl_basic_map_list *collect_basic_maps(
2921
  __isl_take isl_map_list *list)
2922
472
{
2923
472
  int i, n;
2924
472
  isl_ctx *ctx;
2925
472
  isl_basic_map_list *bmap_list;
2926
472
2927
472
  if (!list)
2928
0
    return NULL;
2929
472
  n = isl_map_list_n_map(list);
2930
472
  ctx = isl_map_list_get_ctx(list);
2931
472
  bmap_list = isl_basic_map_list_alloc(ctx, 0);
2932
472
2933
1.56k
  for (i = 0; 
i < n1.56k
;
++i1.09k
)
{1.09k
2934
1.09k
    isl_map *map;
2935
1.09k
    isl_basic_map_list *list_i;
2936
1.09k
2937
1.09k
    map = isl_map_list_get_map(list, i);
2938
1.09k
    map = isl_map_compute_divs(map);
2939
1.09k
    list_i = isl_map_get_basic_map_list(map);
2940
1.09k
    isl_map_free(map);
2941
1.09k
    bmap_list = isl_basic_map_list_concat(bmap_list, list_i);
2942
1.09k
  }
2943
472
2944
472
  isl_map_list_free(list);
2945
472
  return bmap_list;
2946
472
}
2947
2948
/* Compute a superset of the convex hull of "map" that is described
2949
 * by only constraints in the elements of "list".
2950
 *
2951
 * If "map" is the universe, then the convex hull (and therefore
2952
 * any superset of the convexhull) is the universe as well.
2953
 *
2954
 * Otherwise, we collect all the basic maps in the map list and
2955
 * continue with map_unshifted_simple_hull_from_basic_map_list.
2956
 */
2957
__isl_give isl_basic_map *isl_map_unshifted_simple_hull_from_map_list(
2958
  __isl_take isl_map *map, __isl_take isl_map_list *list)
2959
483
{
2960
483
  isl_basic_map_list *bmap_list;
2961
483
  int is_universe;
2962
483
2963
483
  is_universe = isl_map_plain_is_universe(map);
2964
483
  if (is_universe < 0)
2965
0
    map = isl_map_free(map);
2966
483
  if (
is_universe < 0 || 483
is_universe483
)
{11
2967
11
    isl_map_list_free(list);
2968
11
    return isl_map_unshifted_simple_hull(map);
2969
483
  }
2970
483
2971
483
  bmap_list = collect_basic_maps(list);
2972
483
  return map_unshifted_simple_hull_from_basic_map_list(map, bmap_list);
2973
483
}
2974
2975
/* Compute a superset of the convex hull of "set" that is described
2976
 * by only constraints in the elements of "list".
2977
 */
2978
__isl_give isl_basic_set *isl_set_unshifted_simple_hull_from_set_list(
2979
  __isl_take isl_set *set, __isl_take isl_set_list *list)
2980
50
{
2981
50
  return isl_map_unshifted_simple_hull_from_map_list(set, list);
2982
50
}
2983
2984
/* Given a set "set", return parametric bounds on the dimension "dim".
2985
 */
2986
static struct isl_basic_set *set_bounds(struct isl_set *set, int dim)
2987
0
{
2988
0
  unsigned set_dim = isl_set_dim(set, isl_dim_set);
2989
0
  set = isl_set_copy(set);
2990
0
  set = isl_set_eliminate_dims(set, dim + 1, set_dim - (dim + 1));
2991
0
  set = isl_set_eliminate_dims(set, 0, dim);
2992
0
  return isl_set_convex_hull(set);
2993
0
}
2994
2995
/* Computes a "simple hull" and then check if each dimension in the
2996
 * resulting hull is bounded by a symbolic constant.  If not, the
2997
 * hull is intersected with the corresponding bounds on the whole set.
2998
 */
2999
__isl_give isl_basic_set *isl_set_bounded_simple_hull(__isl_take isl_set *set)
3000
0
{
3001
0
  int i, j;
3002
0
  struct isl_basic_set *hull;
3003
0
  unsigned nparam, left;
3004
0
  int removed_divs = 0;
3005
0
3006
0
  hull = isl_set_simple_hull(isl_set_copy(set));
3007
0
  if (!hull)
3008
0
    goto error;
3009
0
3010
0
  nparam = isl_basic_set_dim(hull, isl_dim_param);
3011
0
  for (i = 0; 
i < isl_basic_set_dim(hull, isl_dim_set)0
;
++i0
)
{0
3012
0
    int lower = 0, upper = 0;
3013
0
    struct isl_basic_set *bounds;
3014
0
3015
0
    left = isl_basic_set_total_dim(hull) - nparam - i - 1;
3016
0
    for (j = 0; 
j < hull->n_eq0
;
++j0
)
{0
3017
0
      if (isl_int_is_zero(hull->eq[j][1 + nparam + i]))
3018
0
        continue;
3019
0
      
if (0
isl_seq_first_non_zero(hull->eq[j]+1+nparam+i+1,0
3020
0
                left) == -1)
3021
0
        break;
3022
0
    }
3023
0
    if (j < hull->n_eq)
3024
0
      continue;
3025
0
3026
0
    
for (j = 0; 0
j < hull->n_ineq0
;
++j0
)
{0
3027
0
      if (isl_int_is_zero(hull->ineq[j][1 + nparam + i]))
3028
0
        continue;
3029
0
      
if (0
isl_seq_first_non_zero(hull->ineq[j]+1+nparam+i+1,0
3030
0
                left) != -1 ||
3031
0
          isl_seq_first_non_zero(hull->ineq[j]+1+nparam,
3032
0
                i) != -1)
3033
0
        continue;
3034
0
      
if (0
isl_int_is_pos0
(hull->ineq[j][1 + nparam + i]))
3035
0
        lower = 1;
3036
0
      else
3037
0
        upper = 1;
3038
0
      if (
lower && 0
upper0
)
3039
0
        break;
3040
0
    }
3041
0
3042
0
    if (
lower && 0
upper0
)
3043
0
      continue;
3044
0
3045
0
    
if (0
!removed_divs0
)
{0
3046
0
      set = isl_set_remove_divs(set);
3047
0
      if (!set)
3048
0
        goto error;
3049
0
      removed_divs = 1;
3050
0
    }
3051
0
    bounds = set_bounds(set, i);
3052
0
    hull = isl_basic_set_intersect(hull, bounds);
3053
0
    if (!hull)
3054
0
      goto error;
3055
0
  }
3056
0
3057
0
  isl_set_free(set);
3058
0
  return hull;
3059
0
error:
3060
0
  isl_set_free(set);
3061
0
  isl_basic_set_free(hull);
3062
0
  return NULL;
3063
0
}