Coverage Report

Created: 2018-02-20 23:11

/Users/buildslave/jenkins/workspace/clang-stage2-coverage-R/llvm/tools/polly/lib/External/isl/isl_tab.c
Line
Count
Source (jump to first uncovered line)
1
/*
2
 * Copyright 2008-2009 Katholieke Universiteit Leuven
3
 * Copyright 2013      Ecole Normale Superieure
4
 * Copyright 2014      INRIA Rocquencourt
5
 * Copyright 2016      Sven Verdoolaege
6
 *
7
 * Use of this software is governed by the MIT license
8
 *
9
 * Written by Sven Verdoolaege, K.U.Leuven, Departement
10
 * Computerwetenschappen, Celestijnenlaan 200A, B-3001 Leuven, Belgium
11
 * and Ecole Normale Superieure, 45 rue d'Ulm, 75230 Paris, France
12
 * and Inria Paris - Rocquencourt, Domaine de Voluceau - Rocquencourt,
13
 * B.P. 105 - 78153 Le Chesnay, France
14
 */
15
16
#include <isl_ctx_private.h>
17
#include <isl_mat_private.h>
18
#include <isl_vec_private.h>
19
#include "isl_map_private.h"
20
#include "isl_tab.h"
21
#include <isl_seq.h>
22
#include <isl_config.h>
23
24
#include <bset_to_bmap.c>
25
#include <bset_from_bmap.c>
26
27
/*
28
 * The implementation of tableaus in this file was inspired by Section 8
29
 * of David Detlefs, Greg Nelson and James B. Saxe, "Simplify: a theorem
30
 * prover for program checking".
31
 */
32
33
struct isl_tab *isl_tab_alloc(struct isl_ctx *ctx,
34
  unsigned n_row, unsigned n_var, unsigned M)
35
733k
{
36
733k
  int i;
37
733k
  struct isl_tab *tab;
38
733k
  unsigned off = 2 + M;
39
733k
40
733k
  tab = isl_calloc_type(ctx, struct isl_tab);
41
733k
  if (!tab)
42
0
    return NULL;
43
733k
  tab->mat = isl_mat_alloc(ctx, n_row, off + n_var);
44
733k
  if (!tab->mat)
45
0
    goto error;
46
733k
  tab->var = isl_alloc_array(ctx, struct isl_tab_var, n_var);
47
733k
  if (n_var && 
!tab->var732k
)
48
0
    goto error;
49
733k
  tab->con = isl_alloc_array(ctx, struct isl_tab_var, n_row);
50
733k
  if (n_row && 
!tab->con733k
)
51
0
    goto error;
52
733k
  tab->col_var = isl_alloc_array(ctx, int, n_var);
53
733k
  if (n_var && 
!tab->col_var732k
)
54
0
    goto error;
55
733k
  tab->row_var = isl_alloc_array(ctx, int, n_row);
56
733k
  if (n_row && 
!tab->row_var733k
)
57
0
    goto error;
58
3.96M
  
for (i = 0; 733k
i < n_var;
++i3.23M
) {
59
3.23M
    tab->var[i].index = i;
60
3.23M
    tab->var[i].is_row = 0;
61
3.23M
    tab->var[i].is_nonneg = 0;
62
3.23M
    tab->var[i].is_zero = 0;
63
3.23M
    tab->var[i].is_redundant = 0;
64
3.23M
    tab->var[i].frozen = 0;
65
3.23M
    tab->var[i].negated = 0;
66
3.23M
    tab->col_var[i] = i;
67
3.23M
  }
68
733k
  tab->n_row = 0;
69
733k
  tab->n_con = 0;
70
733k
  tab->n_eq = 0;
71
733k
  tab->max_con = n_row;
72
733k
  tab->n_col = n_var;
73
733k
  tab->n_var = n_var;
74
733k
  tab->max_var = n_var;
75
733k
  tab->n_param = 0;
76
733k
  tab->n_div = 0;
77
733k
  tab->n_dead = 0;
78
733k
  tab->n_redundant = 0;
79
733k
  tab->strict_redundant = 0;
80
733k
  tab->need_undo = 0;
81
733k
  tab->rational = 0;
82
733k
  tab->empty = 0;
83
733k
  tab->in_undo = 0;
84
733k
  tab->M = M;
85
733k
  tab->cone = 0;
86
733k
  tab->bottom.type = isl_tab_undo_bottom;
87
733k
  tab->bottom.next = NULL;
88
733k
  tab->top = &tab->bottom;
89
733k
90
733k
  tab->n_zero = 0;
91
733k
  tab->n_unbounded = 0;
92
733k
  tab->basis = NULL;
93
733k
94
733k
  return tab;
95
0
error:
96
0
  isl_tab_free(tab);
97
0
  return NULL;
98
733k
}
99
100
isl_ctx *isl_tab_get_ctx(struct isl_tab *tab)
101
3.22M
{
102
3.22M
  return tab ? isl_mat_get_ctx(tab->mat) : NULL;
103
3.22M
}
104
105
int isl_tab_extend_cons(struct isl_tab *tab, unsigned n_new)
106
808k
{
107
808k
  unsigned off;
108
808k
109
808k
  if (!tab)
110
0
    return -1;
111
808k
112
808k
  off = 2 + tab->M;
113
808k
114
808k
  if (tab->max_con < tab->n_con + n_new) {
115
91.7k
    struct isl_tab_var *con;
116
91.7k
117
91.7k
    con = isl_realloc_array(tab->mat->ctx, tab->con,
118
91.7k
            struct isl_tab_var, tab->max_con + n_new);
119
91.7k
    if (!con)
120
0
      return -1;
121
91.7k
    tab->con = con;
122
91.7k
    tab->max_con += n_new;
123
91.7k
  }
124
808k
  if (tab->mat->n_row < tab->n_row + n_new) {
125
94.5k
    int *row_var;
126
94.5k
127
94.5k
    tab->mat = isl_mat_extend(tab->mat,
128
94.5k
          tab->n_row + n_new, off + tab->n_col);
129
94.5k
    if (!tab->mat)
130
0
      return -1;
131
94.5k
    row_var = isl_realloc_array(tab->mat->ctx, tab->row_var,
132
94.5k
              int, tab->mat->n_row);
133
94.5k
    if (!row_var)
134
0
      return -1;
135
94.5k
    tab->row_var = row_var;
136
94.5k
    if (tab->row_sign) {
137
219
      enum isl_tab_row_sign *s;
138
219
      s = isl_realloc_array(tab->mat->ctx, tab->row_sign,
139
219
          enum isl_tab_row_sign, tab->mat->n_row);
140
219
      if (!s)
141
0
        return -1;
142
219
      tab->row_sign = s;
143
219
    }
144
94.5k
  }
145
808k
  return 0;
146
808k
}
147
148
/* Make room for at least n_new extra variables.
149
 * Return -1 if anything went wrong.
150
 */
151
int isl_tab_extend_vars(struct isl_tab *tab, unsigned n_new)
152
5.77k
{
153
5.77k
  struct isl_tab_var *var;
154
5.77k
  unsigned off = 2 + tab->M;
155
5.77k
156
5.77k
  if (tab->max_var < tab->n_var + n_new) {
157
4.48k
    var = isl_realloc_array(tab->mat->ctx, tab->var,
158
4.48k
            struct isl_tab_var, tab->n_var + n_new);
159
4.48k
    if (!var)
160
0
      return -1;
161
4.48k
    tab->var = var;
162
4.48k
    tab->max_var = tab->n_var + n_new;
163
4.48k
  }
164
5.77k
165
5.77k
  if (tab->mat->n_col < off + tab->n_col + n_new) {
166
3.15k
    int *p;
167
3.15k
168
3.15k
    tab->mat = isl_mat_extend(tab->mat,
169
3.15k
            tab->mat->n_row, off + tab->n_col + n_new);
170
3.15k
    if (!tab->mat)
171
0
      return -1;
172
3.15k
    p = isl_realloc_array(tab->mat->ctx, tab->col_var,
173
3.15k
              int, tab->n_col + n_new);
174
3.15k
    if (!p)
175
0
      return -1;
176
3.15k
    tab->col_var = p;
177
3.15k
  }
178
5.77k
179
5.77k
  return 0;
180
5.77k
}
181
182
static void free_undo_record(struct isl_tab_undo *undo)
183
3.35M
{
184
3.35M
  switch (undo->type) {
185
3.35M
  case isl_tab_undo_saved_basis:
186
781
    free(undo->u.col_var);
187
781
    break;
188
3.35M
  
default:;3.35M
189
3.35M
  }
190
3.35M
  free(undo);
191
3.35M
}
192
193
static void free_undo(struct isl_tab *tab)
194
739k
{
195
739k
  struct isl_tab_undo *undo, *next;
196
739k
197
1.26M
  for (undo = tab->top; undo && undo != &tab->bottom; 
undo = next530k
) {
198
530k
    next = undo->next;
199
530k
    free_undo_record(undo);
200
530k
  }
201
739k
  tab->top = undo;
202
739k
}
203
204
void isl_tab_free(struct isl_tab *tab)
205
773k
{
206
773k
  if (!tab)
207
33.2k
    return;
208
739k
  free_undo(tab);
209
739k
  isl_mat_free(tab->mat);
210
739k
  isl_vec_free(tab->dual);
211
739k
  isl_basic_map_free(tab->bmap);
212
739k
  free(tab->var);
213
739k
  free(tab->con);
214
739k
  free(tab->row_var);
215
739k
  free(tab->col_var);
216
739k
  free(tab->row_sign);
217
739k
  isl_mat_free(tab->samples);
218
739k
  free(tab->sample_index);
219
739k
  isl_mat_free(tab->basis);
220
739k
  free(tab);
221
739k
}
222
223
struct isl_tab *isl_tab_dup(struct isl_tab *tab)
224
2.83k
{
225
2.83k
  int i;
226
2.83k
  struct isl_tab *dup;
227
2.83k
  unsigned off;
228
2.83k
229
2.83k
  if (!tab)
230
0
    return NULL;
231
2.83k
232
2.83k
  off = 2 + tab->M;
233
2.83k
  dup = isl_calloc_type(tab->mat->ctx, struct isl_tab);
234
2.83k
  if (!dup)
235
0
    return NULL;
236
2.83k
  dup->mat = isl_mat_dup(tab->mat);
237
2.83k
  if (!dup->mat)
238
0
    goto error;
239
2.83k
  dup->var = isl_alloc_array(tab->mat->ctx, struct isl_tab_var, tab->max_var);
240
2.83k
  if (tab->max_var && !dup->var)
241
0
    goto error;
242
26.6k
  
for (i = 0; 2.83k
i < tab->n_var;
++i23.8k
)
243
23.8k
    dup->var[i] = tab->var[i];
244
2.83k
  dup->con = isl_alloc_array(tab->mat->ctx, struct isl_tab_var, tab->max_con);
245
2.83k
  if (tab->max_con && !dup->con)
246
0
    goto error;
247
29.4k
  
for (i = 0; 2.83k
i < tab->n_con;
++i26.5k
)
248
26.5k
    dup->con[i] = tab->con[i];
249
2.83k
  dup->col_var = isl_alloc_array(tab->mat->ctx, int, tab->mat->n_col - off);
250
2.83k
  if ((tab->mat->n_col - off) && !dup->col_var)
251
0
    goto error;
252
13.4k
  
for (i = 0; 2.83k
i < tab->n_col;
++i10.6k
)
253
10.6k
    dup->col_var[i] = tab->col_var[i];
254
2.83k
  dup->row_var = isl_alloc_array(tab->mat->ctx, int, tab->mat->n_row);
255
2.83k
  if (tab->mat->n_row && !dup->row_var)
256
0
    goto error;
257
29.2k
  
for (i = 0; 2.83k
i < tab->n_row;
++i26.4k
)
258
26.4k
    dup->row_var[i] = tab->row_var[i];
259
2.83k
  if (tab->row_sign) {
260
2.82k
    dup->row_sign = isl_alloc_array(tab->mat->ctx, enum isl_tab_row_sign,
261
2.82k
            tab->mat->n_row);
262
2.82k
    if (tab->mat->n_row && !dup->row_sign)
263
0
      goto error;
264
29.2k
    
for (i = 0; 2.82k
i < tab->n_row;
++i26.3k
)
265
26.3k
      dup->row_sign[i] = tab->row_sign[i];
266
2.82k
  }
267
2.83k
  if (tab->samples) {
268
0
    dup->samples = isl_mat_dup(tab->samples);
269
0
    if (!dup->samples)
270
0
      goto error;
271
0
    dup->sample_index = isl_alloc_array(tab->mat->ctx, int,
272
0
              tab->samples->n_row);
273
0
    if (tab->samples->n_row && !dup->sample_index)
274
0
      goto error;
275
0
    dup->n_sample = tab->n_sample;
276
0
    dup->n_outside = tab->n_outside;
277
0
  }
278
2.83k
  dup->n_row = tab->n_row;
279
2.83k
  dup->n_con = tab->n_con;
280
2.83k
  dup->n_eq = tab->n_eq;
281
2.83k
  dup->max_con = tab->max_con;
282
2.83k
  dup->n_col = tab->n_col;
283
2.83k
  dup->n_var = tab->n_var;
284
2.83k
  dup->max_var = tab->max_var;
285
2.83k
  dup->n_param = tab->n_param;
286
2.83k
  dup->n_div = tab->n_div;
287
2.83k
  dup->n_dead = tab->n_dead;
288
2.83k
  dup->n_redundant = tab->n_redundant;
289
2.83k
  dup->rational = tab->rational;
290
2.83k
  dup->empty = tab->empty;
291
2.83k
  dup->strict_redundant = 0;
292
2.83k
  dup->need_undo = 0;
293
2.83k
  dup->in_undo = 0;
294
2.83k
  dup->M = tab->M;
295
2.83k
  tab->cone = tab->cone;
296
2.83k
  dup->bottom.type = isl_tab_undo_bottom;
297
2.83k
  dup->bottom.next = NULL;
298
2.83k
  dup->top = &dup->bottom;
299
2.83k
300
2.83k
  dup->n_zero = tab->n_zero;
301
2.83k
  dup->n_unbounded = tab->n_unbounded;
302
2.83k
  dup->basis = isl_mat_dup(tab->basis);
303
2.83k
304
2.83k
  return dup;
305
0
error:
306
0
  isl_tab_free(dup);
307
0
  return NULL;
308
2.83k
}
309
310
/* Construct the coefficient matrix of the product tableau
311
 * of two tableaus.
312
 * mat{1,2} is the coefficient matrix of tableau {1,2}
313
 * row{1,2} is the number of rows in tableau {1,2}
314
 * col{1,2} is the number of columns in tableau {1,2}
315
 * off is the offset to the coefficient column (skipping the
316
 *  denominator, the constant term and the big parameter if any)
317
 * r{1,2} is the number of redundant rows in tableau {1,2}
318
 * d{1,2} is the number of dead columns in tableau {1,2}
319
 *
320
 * The order of the rows and columns in the result is as explained
321
 * in isl_tab_product.
322
 */
323
static struct isl_mat *tab_mat_product(struct isl_mat *mat1,
324
  struct isl_mat *mat2, unsigned row1, unsigned row2,
325
  unsigned col1, unsigned col2,
326
  unsigned off, unsigned r1, unsigned r2, unsigned d1, unsigned d2)
327
3.12k
{
328
3.12k
  int i;
329
3.12k
  struct isl_mat *prod;
330
3.12k
  unsigned n;
331
3.12k
332
3.12k
  prod = isl_mat_alloc(mat1->ctx, mat1->n_row + mat2->n_row,
333
3.12k
          off + col1 + col2);
334
3.12k
  if (!prod)
335
0
    return NULL;
336
3.12k
337
3.12k
  n = 0;
338
11.0k
  for (i = 0; i < r1; 
++i7.89k
) {
339
7.89k
    isl_seq_cpy(prod->row[n + i], mat1->row[i], off + d1);
340
7.89k
    isl_seq_clr(prod->row[n + i] + off + d1, d2);
341
7.89k
    isl_seq_cpy(prod->row[n + i] + off + d1 + d2,
342
7.89k
        mat1->row[i] + off + d1, col1 - d1);
343
7.89k
    isl_seq_clr(prod->row[n + i] + off + col1 + d1, col2 - d2);
344
7.89k
  }
345
3.12k
346
3.12k
  n += r1;
347
11.0k
  for (i = 0; i < r2; 
++i7.89k
) {
348
7.89k
    isl_seq_cpy(prod->row[n + i], mat2->row[i], off);
349
7.89k
    isl_seq_clr(prod->row[n + i] + off, d1);
350
7.89k
    isl_seq_cpy(prod->row[n + i] + off + d1,
351
7.89k
          mat2->row[i] + off, d2);
352
7.89k
    isl_seq_clr(prod->row[n + i] + off + d1 + d2, col1 - d1);
353
7.89k
    isl_seq_cpy(prod->row[n + i] + off + col1 + d1,
354
7.89k
          mat2->row[i] + off + d2, col2 - d2);
355
7.89k
  }
356
3.12k
357
3.12k
  n += r2;
358
36.4k
  for (i = 0; i < row1 - r1; 
++i33.3k
) {
359
33.3k
    isl_seq_cpy(prod->row[n + i], mat1->row[r1 + i], off + d1);
360
33.3k
    isl_seq_clr(prod->row[n + i] + off + d1, d2);
361
33.3k
    isl_seq_cpy(prod->row[n + i] + off + d1 + d2,
362
33.3k
        mat1->row[r1 + i] + off + d1, col1 - d1);
363
33.3k
    isl_seq_clr(prod->row[n + i] + off + col1 + d1, col2 - d2);
364
33.3k
  }
365
3.12k
366
3.12k
  n += row1 - r1;
367
36.4k
  for (i = 0; i < row2 - r2; 
++i33.3k
) {
368
33.3k
    isl_seq_cpy(prod->row[n + i], mat2->row[r2 + i], off);
369
33.3k
    isl_seq_clr(prod->row[n + i] + off, d1);
370
33.3k
    isl_seq_cpy(prod->row[n + i] + off + d1,
371
33.3k
          mat2->row[r2 + i] + off, d2);
372
33.3k
    isl_seq_clr(prod->row[n + i] + off + d1 + d2, col1 - d1);
373
33.3k
    isl_seq_cpy(prod->row[n + i] + off + col1 + d1,
374
33.3k
          mat2->row[r2 + i] + off + d2, col2 - d2);
375
33.3k
  }
376
3.12k
377
3.12k
  return prod;
378
3.12k
}
379
380
/* Update the row or column index of a variable that corresponds
381
 * to a variable in the first input tableau.
382
 */
383
static void update_index1(struct isl_tab_var *var,
384
  unsigned r1, unsigned r2, unsigned d1, unsigned d2)
385
58.3k
{
386
58.3k
  if (var->index == -1)
387
288
    return;
388
58.1k
  if (var->is_row && 
var->index >= r141.2k
)
389
33.3k
    var->index += r2;
390
58.1k
  if (!var->is_row && 
var->index >= d116.8k
)
391
15.7k
    var->index += d2;
392
58.1k
}
393
394
/* Update the row or column index of a variable that corresponds
395
 * to a variable in the second input tableau.
396
 */
397
static void update_index2(struct isl_tab_var *var,
398
  unsigned row1, unsigned col1,
399
  unsigned r1, unsigned r2, unsigned d1, unsigned d2)
400
58.3k
{
401
58.3k
  if (var->index == -1)
402
288
    return;
403
58.1k
  if (var->is_row) {
404
41.2k
    if (var->index < r2)
405
7.89k
      var->index += r1;
406
33.3k
    else
407
33.3k
      var->index += row1;
408
41.2k
  } else {
409
16.8k
    if (var->index < d2)
410
1.19k
      var->index += d1;
411
15.7k
    else
412
15.7k
      var->index += col1;
413
16.8k
  }
414
58.1k
}
415
416
/* Create a tableau that represents the Cartesian product of the sets
417
 * represented by tableaus tab1 and tab2.
418
 * The order of the rows in the product is
419
 *  - redundant rows of tab1
420
 *  - redundant rows of tab2
421
 *  - non-redundant rows of tab1
422
 *  - non-redundant rows of tab2
423
 * The order of the columns is
424
 *  - denominator
425
 *  - constant term
426
 *  - coefficient of big parameter, if any
427
 *  - dead columns of tab1
428
 *  - dead columns of tab2
429
 *  - live columns of tab1
430
 *  - live columns of tab2
431
 * The order of the variables and the constraints is a concatenation
432
 * of order in the two input tableaus.
433
 */
434
struct isl_tab *isl_tab_product(struct isl_tab *tab1, struct isl_tab *tab2)
435
3.12k
{
436
3.12k
  int i;
437
3.12k
  struct isl_tab *prod;
438
3.12k
  unsigned off;
439
3.12k
  unsigned r1, r2, d1, d2;
440
3.12k
441
3.12k
  if (!tab1 || !tab2)
442
0
    return NULL;
443
3.12k
444
3.12k
  isl_assert(tab1->mat->ctx, tab1->M == tab2->M, return NULL);
445
3.12k
  isl_assert(tab1->mat->ctx, tab1->rational == tab2->rational, return NULL);
446
3.12k
  isl_assert(tab1->mat->ctx, tab1->cone == tab2->cone, return NULL);
447
3.12k
  isl_assert(tab1->mat->ctx, !tab1->row_sign, return NULL);
448
3.12k
  isl_assert(tab1->mat->ctx, !tab2->row_sign, return NULL);
449
3.12k
  isl_assert(tab1->mat->ctx, tab1->n_param == 0, return NULL);
450
3.12k
  isl_assert(tab1->mat->ctx, tab2->n_param == 0, return NULL);
451
3.12k
  isl_assert(tab1->mat->ctx, tab1->n_div == 0, return NULL);
452
3.12k
  isl_assert(tab1->mat->ctx, tab2->n_div == 0, return NULL);
453
3.12k
454
3.12k
  off = 2 + tab1->M;
455
3.12k
  r1 = tab1->n_redundant;
456
3.12k
  r2 = tab2->n_redundant;
457
3.12k
  d1 = tab1->n_dead;
458
3.12k
  d2 = tab2->n_dead;
459
3.12k
  prod = isl_calloc_type(tab1->mat->ctx, struct isl_tab);
460
3.12k
  if (!prod)
461
0
    return NULL;
462
3.12k
  prod->mat = tab_mat_product(tab1->mat, tab2->mat,
463
3.12k
        tab1->n_row, tab2->n_row,
464
3.12k
        tab1->n_col, tab2->n_col, off, r1, r2, d1, d2);
465
3.12k
  if (!prod->mat)
466
0
    goto error;
467
3.12k
  prod->var = isl_alloc_array(tab1->mat->ctx, struct isl_tab_var,
468
3.12k
          tab1->max_var + tab2->max_var);
469
3.12k
  if ((tab1->max_var + tab2->max_var) && !prod->var)
470
0
    goto error;
471
20.3k
  
for (i = 0; 3.12k
i < tab1->n_var;
++i17.1k
) {
472
17.1k
    prod->var[i] = tab1->var[i];
473
17.1k
    update_index1(&prod->var[i], r1, r2, d1, d2);
474
17.1k
  }
475
20.3k
  for (i = 0; i < tab2->n_var; 
++i17.1k
) {
476
17.1k
    prod->var[tab1->n_var + i] = tab2->var[i];
477
17.1k
    update_index2(&prod->var[tab1->n_var + i],
478
17.1k
        tab1->n_row, tab1->n_col,
479
17.1k
        r1, r2, d1, d2);
480
17.1k
  }
481
3.12k
  prod->con = isl_alloc_array(tab1->mat->ctx, struct isl_tab_var,
482
3.12k
          tab1->max_con +  tab2->max_con);
483
3.12k
  if ((tab1->max_con + tab2->max_con) && !prod->con)
484
0
    goto error;
485
44.3k
  
for (i = 0; 3.12k
i < tab1->n_con;
++i41.2k
) {
486
41.2k
    prod->con[i] = tab1->con[i];
487
41.2k
    update_index1(&prod->con[i], r1, r2, d1, d2);
488
41.2k
  }
489
44.3k
  for (i = 0; i < tab2->n_con; 
++i41.2k
) {
490
41.2k
    prod->con[tab1->n_con + i] = tab2->con[i];
491
41.2k
    update_index2(&prod->con[tab1->n_con + i],
492
41.2k
        tab1->n_row, tab1->n_col,
493
41.2k
        r1, r2, d1, d2);
494
41.2k
  }
495
3.12k
  prod->col_var = isl_alloc_array(tab1->mat->ctx, int,
496
3.12k
          tab1->n_col + tab2->n_col);
497
3.12k
  if ((tab1->n_col + tab2->n_col) && !prod->col_var)
498
0
    goto error;
499
20.0k
  
for (i = 0; 3.12k
i < tab1->n_col;
++i16.8k
) {
500
16.8k
    int pos = i < d1 ? 
i1.19k
:
i + d215.7k
;
501
16.8k
    prod->col_var[pos] = tab1->col_var[i];
502
16.8k
  }
503
20.0k
  for (i = 0; i < tab2->n_col; 
++i16.8k
) {
504
16.8k
    int pos = i < d2 ? 
d1 + i1.19k
:
tab1->n_col + i15.7k
;
505
16.8k
    int t = tab2->col_var[i];
506
16.8k
    if (t >= 0)
507
181
      t += tab1->n_var;
508
16.7k
    else
509
16.7k
      t -= tab1->n_con;
510
16.8k
    prod->col_var[pos] = t;
511
16.8k
  }
512
3.12k
  prod->row_var = isl_alloc_array(tab1->mat->ctx, int,
513
3.12k
          tab1->mat->n_row + tab2->mat->n_row);
514
3.12k
  if ((tab1->mat->n_row + tab2->mat->n_row) && !prod->row_var)
515
0
    goto error;
516
44.3k
  
for (i = 0; 3.12k
i < tab1->n_row;
++i41.2k
) {
517
41.2k
    int pos = i < r1 ? 
i7.89k
:
i + r233.3k
;
518
41.2k
    prod->row_var[pos] = tab1->row_var[i];
519
41.2k
  }
520
44.3k
  for (i = 0; i < tab2->n_row; 
++i41.2k
) {
521
41.2k
    int pos = i < r2 ? 
r1 + i7.89k
:
tab1->n_row + i33.3k
;
522
41.2k
    int t = tab2->row_var[i];
523
41.2k
    if (t >= 0)
524
17.0k
      t += tab1->n_var;
525
24.2k
    else
526
24.2k
      t -= tab1->n_con;
527
41.2k
    prod->row_var[pos] = t;
528
41.2k
  }
529
3.12k
  prod->samples = NULL;
530
3.12k
  prod->sample_index = NULL;
531
3.12k
  prod->n_row = tab1->n_row + tab2->n_row;
532
3.12k
  prod->n_con = tab1->n_con + tab2->n_con;
533
3.12k
  prod->n_eq = 0;
534
3.12k
  prod->max_con = tab1->max_con + tab2->max_con;
535
3.12k
  prod->n_col = tab1->n_col + tab2->n_col;
536
3.12k
  prod->n_var = tab1->n_var + tab2->n_var;
537
3.12k
  prod->max_var = tab1->max_var + tab2->max_var;
538
3.12k
  prod->n_param = 0;
539
3.12k
  prod->n_div = 0;
540
3.12k
  prod->n_dead = tab1->n_dead + tab2->n_dead;
541
3.12k
  prod->n_redundant = tab1->n_redundant + tab2->n_redundant;
542
3.12k
  prod->rational = tab1->rational;
543
3.12k
  prod->empty = tab1->empty || tab2->empty;
544
3.12k
  prod->strict_redundant = tab1->strict_redundant || tab2->strict_redundant;
545
3.12k
  prod->need_undo = 0;
546
3.12k
  prod->in_undo = 0;
547
3.12k
  prod->M = tab1->M;
548
3.12k
  prod->cone = tab1->cone;
549
3.12k
  prod->bottom.type = isl_tab_undo_bottom;
550
3.12k
  prod->bottom.next = NULL;
551
3.12k
  prod->top = &prod->bottom;
552
3.12k
553
3.12k
  prod->n_zero = 0;
554
3.12k
  prod->n_unbounded = 0;
555
3.12k
  prod->basis = NULL;
556
3.12k
557
3.12k
  return prod;
558
0
error:
559
0
  isl_tab_free(prod);
560
0
  return NULL;
561
3.12k
}
562
563
static struct isl_tab_var *var_from_index(struct isl_tab *tab, int i)
564
116M
{
565
116M
  if (i >= 0)
566
30.9M
    return &tab->var[i];
567
85.0M
  else
568
85.0M
    return &tab->con[~i];
569
116M
}
570
571
struct isl_tab_var *isl_tab_var_from_row(struct isl_tab *tab, int i)
572
87.7M
{
573
87.7M
  return var_from_index(tab, tab->row_var[i]);
574
87.7M
}
575
576
static struct isl_tab_var *var_from_col(struct isl_tab *tab, int i)
577
25.8M
{
578
25.8M
  return var_from_index(tab, tab->col_var[i]);
579
25.8M
}
580
581
/* Check if there are any upper bounds on column variable "var",
582
 * i.e., non-negative rows where var appears with a negative coefficient.
583
 * Return 1 if there are no such bounds.
584
 */
585
static int max_is_manifestly_unbounded(struct isl_tab *tab,
586
  struct isl_tab_var *var)
587
2.09M
{
588
2.09M
  int i;
589
2.09M
  unsigned off = 2 + tab->M;
590
2.09M
591
2.09M
  if (var->is_row)
592
1.43M
    return 0;
593
3.73M
  
for (i = tab->n_redundant; 657k
i < tab->n_row;
++i3.07M
) {
594
3.28M
    if (!isl_int_is_neg(tab->mat->row[i][off + var->index]))
595
3.28M
      
continue2.84M
;
596
443k
    if (isl_tab_var_from_row(tab, i)->is_nonneg)
597
208k
      return 0;
598
443k
  }
599
657k
  
return 1449k
;
600
657k
}
601
602
/* Check if there are any lower bounds on column variable "var",
603
 * i.e., non-negative rows where var appears with a positive coefficient.
604
 * Return 1 if there are no such bounds.
605
 */
606
static int min_is_manifestly_unbounded(struct isl_tab *tab,
607
  struct isl_tab_var *var)
608
1.62M
{
609
1.62M
  int i;
610
1.62M
  unsigned off = 2 + tab->M;
611
1.62M
612
1.62M
  if (var->is_row)
613
775k
    return 0;
614
6.36M
  
for (i = tab->n_redundant; 852k
i < tab->n_row;
++i5.50M
) {
615
5.99M
    if (!isl_int_is_pos(tab->mat->row[i][off + var->index]))
616
5.99M
      
continue5.21M
;
617
783k
    if (isl_tab_var_from_row(tab, i)->is_nonneg)
618
491k
      return 0;
619
783k
  }
620
852k
  
return 1361k
;
621
852k
}
622
623
static int row_cmp(struct isl_tab *tab, int r1, int r2, int c, isl_int *t)
624
1.55M
{
625
1.55M
  unsigned off = 2 + tab->M;
626
1.55M
627
1.55M
  if (tab->M) {
628
0
    int s;
629
0
    isl_int_mul(*t, tab->mat->row[r1][2], tab->mat->row[r2][off+c]);
630
0
    isl_int_submul(*t, tab->mat->row[r2][2], tab->mat->row[r1][off+c]);
631
0
    s = isl_int_sgn(*t);
632
0
    if (s)
633
0
      return s;
634
1.55M
  }
635
1.55M
  isl_int_mul(*t, tab->mat->row[r1][1], tab->mat->row[r2][off + c]);
636
1.55M
  isl_int_submul(*t, tab->mat->row[r2][1], tab->mat->row[r1][off + c]);
637
1.55M
  return isl_int_sgn(*t);
638
1.55M
}
639
640
/* Given the index of a column "c", return the index of a row
641
 * that can be used to pivot the column in, with either an increase
642
 * (sgn > 0) or a decrease (sgn < 0) of the corresponding variable.
643
 * If "var" is not NULL, then the row returned will be different from
644
 * the one associated with "var".
645
 *
646
 * Each row in the tableau is of the form
647
 *
648
 *  x_r = a_r0 + \sum_i a_ri x_i
649
 *
650
 * Only rows with x_r >= 0 and with the sign of a_ri opposite to "sgn"
651
 * impose any limit on the increase or decrease in the value of x_c
652
 * and this bound is equal to a_r0 / |a_rc|.  We are therefore looking
653
 * for the row with the smallest (most stringent) such bound.
654
 * Note that the common denominator of each row drops out of the fraction.
655
 * To check if row j has a smaller bound than row r, i.e.,
656
 * a_j0 / |a_jc| < a_r0 / |a_rc| or a_j0 |a_rc| < a_r0 |a_jc|,
657
 * we check if -sign(a_jc) (a_j0 a_rc - a_r0 a_jc) < 0,
658
 * where -sign(a_jc) is equal to "sgn".
659
 */
660
static int pivot_row(struct isl_tab *tab,
661
  struct isl_tab_var *var, int sgn, int c)
662
3.25M
{
663
3.25M
  int j, r, tsgn;
664
3.25M
  isl_int t;
665
3.25M
  unsigned off = 2 + tab->M;
666
3.25M
667
3.25M
  isl_int_init(t);
668
3.25M
  r = -1;
669
42.1M
  for (j = tab->n_redundant; j < tab->n_row; 
++j38.8M
) {
670
38.8M
    if (var && 
j == var->index33.2M
)
671
2.82M
      continue;
672
36.0M
    if (!isl_tab_var_from_row(tab, j)->is_nonneg)
673
7.49M
      continue;
674
28.5M
    if (sgn * isl_int_sgn(tab->mat->row[j][off + c]) >= 0)
675
24.8M
      continue;
676
3.76M
    if (r < 0) {
677
2.20M
      r = j;
678
2.20M
      continue;
679
2.20M
    }
680
1.55M
    tsgn = sgn * row_cmp(tab, r, j, c, &t);
681
1.55M
    if (tsgn < 0 || 
(1.18M
tsgn == 01.18M
&&
682
1.18M
              
tab->row_var[j] < tab->row_var[r]526k
))
683
781k
      r = j;
684
1.55M
  }
685
3.25M
  isl_int_clear(t);
686
3.25M
  return r;
687
3.25M
}
688
689
/* Find a pivot (row and col) that will increase (sgn > 0) or decrease
690
 * (sgn < 0) the value of row variable var.
691
 * If not NULL, then skip_var is a row variable that should be ignored
692
 * while looking for a pivot row.  It is usually equal to var.
693
 *
694
 * As the given row in the tableau is of the form
695
 *
696
 *  x_r = a_r0 + \sum_i a_ri x_i
697
 *
698
 * we need to find a column such that the sign of a_ri is equal to "sgn"
699
 * (such that an increase in x_i will have the desired effect) or a
700
 * column with a variable that may attain negative values.
701
 * If a_ri is positive, then we need to move x_i in the same direction
702
 * to obtain the desired effect.  Otherwise, x_i has to move in the
703
 * opposite direction.
704
 */
705
static void find_pivot(struct isl_tab *tab,
706
  struct isl_tab_var *var, struct isl_tab_var *skip_var,
707
  int sgn, int *row, int *col)
708
4.17M
{
709
4.17M
  int j, r, c;
710
4.17M
  isl_int *tr;
711
4.17M
712
4.17M
  *row = *col = -1;
713
4.17M
714
4.17M
  isl_assert(tab->mat->ctx, var->is_row, return);
715
4.17M
  tr = tab->mat->row[var->index] + 2 + tab->M;
716
4.17M
717
4.17M
  c = -1;
718
32.7M
  for (j = tab->n_dead; j < tab->n_col; 
++j28.5M
) {
719
28.5M
    if (isl_int_is_zero(tr[j]))
720
28.5M
      
continue22.3M
;
721
6.17M
    if (isl_int_sgn(tr[j]) != sgn &&
722
6.17M
        
var_from_col(tab, j)->is_nonneg3.46M
)
723
2.30M
      continue;
724
3.87M
    if (c < 0 || 
tab->col_var[j] < tab->col_var[c]940k
)
725
3.15M
      c = j;
726
3.87M
  }
727
4.17M
  if (c < 0)
728
1.23M
    return;
729
2.93M
730
2.93M
  sgn *= isl_int_sgn(tr[c]);
731
2.93M
  r = pivot_row(tab, skip_var, sgn, c);
732
2.93M
  *row = r < 0 ? 
var->index1.04M
:
r1.89M
;
733
2.93M
  *col = c;
734
2.93M
}
735
736
/* Return 1 if row "row" represents an obviously redundant inequality.
737
 * This means
738
 *  - it represents an inequality or a variable
739
 *  - that is the sum of a non-negative sample value and a positive
740
 *    combination of zero or more non-negative constraints.
741
 */
742
int isl_tab_row_is_redundant(struct isl_tab *tab, int row)
743
18.5M
{
744
18.5M
  int i;
745
18.5M
  unsigned off = 2 + tab->M;
746
18.5M
747
18.5M
  if (tab->row_var[row] < 0 && 
!isl_tab_var_from_row(tab, row)->is_nonneg14.5M
)
748
894k
    return 0;
749
17.6M
750
17.6M
  if (isl_int_is_neg(tab->mat->row[row][1]))
751
17.6M
    
return 01.49M
;
752
16.1M
  if (tab->strict_redundant && 
isl_int_is_zero43
(tab->mat->row[row][1]))
753
16.1M
    
return 042
;
754
16.1M
  if (tab->M && 
isl_int_is_neg41.6k
(tab->mat->row[row][2]))
755
16.1M
    
return 02.24k
;
756
16.1M
757
65.1M
  
for (i = tab->n_dead; 16.1M
i < tab->n_col;
++i49.0M
) {
758
63.4M
    if (isl_int_is_zero(tab->mat->row[row][off + i]))
759
63.4M
      
continue44.0M
;
760
19.4M
    if (tab->col_var[i] >= 0)
761
7.64M
      return 0;
762
11.8M
    if (isl_int_is_neg(tab->mat->row[row][off + i]))
763
11.8M
      
return 06.60M
;
764
5.21M
    if (!var_from_col(tab, i)->is_nonneg)
765
162k
      return 0;
766
5.21M
  }
767
16.1M
  
return 11.71M
;
768
16.1M
}
769
770
static void swap_rows(struct isl_tab *tab, int row1, int row2)
771
1.70M
{
772
1.70M
  int t;
773
1.70M
  enum isl_tab_row_sign s;
774
1.70M
775
1.70M
  t = tab->row_var[row1];
776
1.70M
  tab->row_var[row1] = tab->row_var[row2];
777
1.70M
  tab->row_var[row2] = t;
778
1.70M
  isl_tab_var_from_row(tab, row1)->index = row1;
779
1.70M
  isl_tab_var_from_row(tab, row2)->index = row2;
780
1.70M
  tab->mat = isl_mat_swap_rows(tab->mat, row1, row2);
781
1.70M
782
1.70M
  if (!tab->row_sign)
783
1.69M
    return;
784
9.58k
  s = tab->row_sign[row1];
785
9.58k
  tab->row_sign[row1] = tab->row_sign[row2];
786
9.58k
  tab->row_sign[row2] = s;
787
9.58k
}
788
789
static isl_stat push_union(struct isl_tab *tab,
790
  enum isl_tab_undo_type type, union isl_tab_undo_val u) WARN_UNUSED;
791
792
/* Push record "u" onto the undo stack of "tab", provided "tab"
793
 * keeps track of undo information.
794
 *
795
 * If the record cannot be pushed, then mark the undo stack as invalid
796
 * such that a later rollback attempt will not try to undo earlier
797
 * records without having been able to undo the current record.
798
 */
799
static isl_stat push_union(struct isl_tab *tab,
800
  enum isl_tab_undo_type type, union isl_tab_undo_val u)
801
12.9M
{
802
12.9M
  struct isl_tab_undo *undo;
803
12.9M
804
12.9M
  if (!tab)
805
0
    return isl_stat_error;
806
12.9M
  if (!tab->need_undo)
807
9.55M
    return isl_stat_ok;
808
3.35M
809
3.35M
  undo = isl_alloc_type(tab->mat->ctx, struct isl_tab_undo);
810
3.35M
  if (!undo)
811
0
    goto error;
812
3.35M
  undo->type = type;
813
3.35M
  undo->u = u;
814
3.35M
  undo->next = tab->top;
815
3.35M
  tab->top = undo;
816
3.35M
817
3.35M
  return isl_stat_ok;
818
0
error:
819
0
  free_undo(tab);
820
0
  tab->top = NULL;
821
0
  return isl_stat_error;
822
3.35M
}
823
824
isl_stat isl_tab_push_var(struct isl_tab *tab,
825
  enum isl_tab_undo_type type, struct isl_tab_var *var)
826
12.4M
{
827
12.4M
  union isl_tab_undo_val u;
828
12.4M
  if (var->is_row)
829
12.1M
    u.var_index = tab->row_var[var->index];
830
243k
  else
831
243k
    u.var_index = tab->col_var[var->index];
832
12.4M
  return push_union(tab, type, u);
833
12.4M
}
834
835
isl_stat isl_tab_push(struct isl_tab *tab, enum isl_tab_undo_type type)
836
424k
{
837
424k
  union isl_tab_undo_val u = { 0 };
838
424k
  return push_union(tab, type, u);
839
424k
}
840
841
/* Push a record on the undo stack describing the current basic
842
 * variables, so that the this state can be restored during rollback.
843
 */
844
isl_stat isl_tab_push_basis(struct isl_tab *tab)
845
781
{
846
781
  int i;
847
781
  union isl_tab_undo_val u;
848
781
849
781
  u.col_var = isl_alloc_array(tab->mat->ctx, int, tab->n_col);
850
781
  if (tab->n_col && !u.col_var)
851
0
    return isl_stat_error;
852
8.86k
  
for (i = 0; 781
i < tab->n_col;
++i8.08k
)
853
8.08k
    u.col_var[i] = tab->col_var[i];
854
781
  return push_union(tab, isl_tab_undo_saved_basis, u);
855
781
}
856
857
isl_stat isl_tab_push_callback(struct isl_tab *tab,
858
  struct isl_tab_callback *callback)
859
21.8k
{
860
21.8k
  union isl_tab_undo_val u;
861
21.8k
  u.callback = callback;
862
21.8k
  return push_union(tab, isl_tab_undo_callback, u);
863
21.8k
}
864
865
struct isl_tab *isl_tab_init_samples(struct isl_tab *tab)
866
7.83k
{
867
7.83k
  if (!tab)
868
0
    return NULL;
869
7.83k
870
7.83k
  tab->n_sample = 0;
871
7.83k
  tab->n_outside = 0;
872
7.83k
  tab->samples = isl_mat_alloc(tab->mat->ctx, 1, 1 + tab->n_var);
873
7.83k
  if (!tab->samples)
874
0
    goto error;
875
7.83k
  tab->sample_index = isl_alloc_array(tab->mat->ctx, int, 1);
876
7.83k
  if (!tab->sample_index)
877
0
    goto error;
878
7.83k
  return tab;
879
0
error:
880
0
  isl_tab_free(tab);
881
0
  return NULL;
882
7.83k
}
883
884
int isl_tab_add_sample(struct isl_tab *tab, __isl_take isl_vec *sample)
885
11.8k
{
886
11.8k
  if (!tab || !sample)
887
0
    goto error;
888
11.8k
889
11.8k
  if (tab->n_sample + 1 > tab->samples->n_row) {
890
3.87k
    int *t = isl_realloc_array(tab->mat->ctx,
891
3.87k
          tab->sample_index, int, tab->n_sample + 1);
892
3.87k
    if (!t)
893
0
      goto error;
894
3.87k
    tab->sample_index = t;
895
3.87k
  }
896
11.8k
897
11.8k
  tab->samples = isl_mat_extend(tab->samples,
898
11.8k
        tab->n_sample + 1, tab->samples->n_col);
899
11.8k
  if (!tab->samples)
900
0
    goto error;
901
11.8k
902
11.8k
  isl_seq_cpy(tab->samples->row[tab->n_sample], sample->el, sample->size);
903
11.8k
  isl_vec_free(sample);
904
11.8k
  tab->sample_index[tab->n_sample] = tab->n_sample;
905
11.8k
  tab->n_sample++;
906
11.8k
907
11.8k
  return 0;
908
0
error:
909
0
  isl_vec_free(sample);
910
0
  return -1;
911
11.8k
}
912
913
struct isl_tab *isl_tab_drop_sample(struct isl_tab *tab, int s)
914
6.41k
{
915
6.41k
  if (s != tab->n_outside) {
916
4.00k
    int t = tab->sample_index[tab->n_outside];
917
4.00k
    tab->sample_index[tab->n_outside] = tab->sample_index[s];
918
4.00k
    tab->sample_index[s] = t;
919
4.00k
    isl_mat_swap_rows(tab->samples, tab->n_outside, s);
920
4.00k
  }
921
6.41k
  tab->n_outside++;
922
6.41k
  if (isl_tab_push(tab, isl_tab_undo_drop_sample) < 0) {
923
0
    isl_tab_free(tab);
924
0
    return NULL;
925
0
  }
926
6.41k
927
6.41k
  return tab;
928
6.41k
}
929
930
/* Record the current number of samples so that we can remove newer
931
 * samples during a rollback.
932
 */
933
isl_stat isl_tab_save_samples(struct isl_tab *tab)
934
28.8k
{
935
28.8k
  union isl_tab_undo_val u;
936
28.8k
937
28.8k
  if (!tab)
938
0
    return isl_stat_error;
939
28.8k
940
28.8k
  u.n = tab->n_sample;
941
28.8k
  return push_union(tab, isl_tab_undo_saved_samples, u);
942
28.8k
}
943
944
/* Mark row with index "row" as being redundant.
945
 * If we may need to undo the operation or if the row represents
946
 * a variable of the original problem, the row is kept,
947
 * but no longer considered when looking for a pivot row.
948
 * Otherwise, the row is simply removed.
949
 *
950
 * The row may be interchanged with some other row.  If it
951
 * is interchanged with a later row, return 1.  Otherwise return 0.
952
 * If the rows are checked in order in the calling function,
953
 * then a return value of 1 means that the row with the given
954
 * row number may now contain a different row that hasn't been checked yet.
955
 */
956
int isl_tab_mark_redundant(struct isl_tab *tab, int row)
957
2.30M
{
958
2.30M
  struct isl_tab_var *var = isl_tab_var_from_row(tab, row);
959
2.30M
  var->is_redundant = 1;
960
2.30M
  isl_assert(tab->mat->ctx, row >= tab->n_redundant, return -1);
961
2.30M
  if (tab->preserve || 
tab->need_undo1.63M
||
tab->row_var[row] >= 01.51M
) {
962
1.54M
    if (tab->row_var[row] >= 0 && 
!var->is_nonneg1.16M
) {
963
1.15M
      var->is_nonneg = 1;
964
1.15M
      if (isl_tab_push_var(tab, isl_tab_undo_nonneg, var) < 0)
965
0
        return -1;
966
1.54M
    }
967
1.54M
    if (row != tab->n_redundant)
968
1.07M
      swap_rows(tab, row, tab->n_redundant);
969
1.54M
    tab->n_redundant++;
970
1.54M
    return isl_tab_push_var(tab, isl_tab_undo_redundant, var);
971
1.54M
  } else {
972
759k
    if (row != tab->n_row - 1)
973
509k
      swap_rows(tab, row, tab->n_row - 1);
974
759k
    isl_tab_var_from_row(tab, tab->n_row - 1)->index = -1;
975
759k
    tab->n_row--;
976
759k
    return 1;
977
759k
  }
978
2.30M
}
979
980
/* Mark "tab" as a rational tableau.
981
 * If it wasn't marked as a rational tableau already and if we may
982
 * need to undo changes, then arrange for the marking to be undone
983
 * during the undo.
984
 */
985
int isl_tab_mark_rational(struct isl_tab *tab)
986
10.4k
{
987
10.4k
  if (!tab)
988
0
    return -1;
989
10.4k
  if (!tab->rational && 
tab->need_undo10.4k
)
990
10.4k
    if (isl_tab_push(tab, isl_tab_undo_rational) < 0)
991
0
      return -1;
992
10.4k
  tab->rational = 1;
993
10.4k
  return 0;
994
10.4k
}
995
996
isl_stat isl_tab_mark_empty(struct isl_tab *tab)
997
74.0k
{
998
74.0k
  if (!tab)
999
0
    return isl_stat_error;
1000
74.0k
  if (!tab->empty && 
tab->need_undo73.2k
)
1001
63.0k
    if (isl_tab_push(tab, isl_tab_undo_empty) < 0)
1002
0
      return isl_stat_error;
1003
74.0k
  tab->empty = 1;
1004
74.0k
  return isl_stat_ok;
1005
74.0k
}
1006
1007
int isl_tab_freeze_constraint(struct isl_tab *tab, int con)
1008
377k
{
1009
377k
  struct isl_tab_var *var;
1010
377k
1011
377k
  if (!tab)
1012
0
    return -1;
1013
377k
1014
377k
  var = &tab->con[con];
1015
377k
  if (var->frozen)
1016
0
    return 0;
1017
377k
  if (var->index < 0)
1018
25.3k
    return 0;
1019
352k
  var->frozen = 1;
1020
352k
1021
352k
  if (tab->need_undo)
1022
318k
    return isl_tab_push_var(tab, isl_tab_undo_freeze, var);
1023
33.4k
1024
33.4k
  return 0;
1025
33.4k
}
1026
1027
/* Update the rows signs after a pivot of "row" and "col", with "row_sgn"
1028
 * the original sign of the pivot element.
1029
 * We only keep track of row signs during PILP solving and in this case
1030
 * we only pivot a row with negative sign (meaning the value is always
1031
 * non-positive) using a positive pivot element.
1032
 *
1033
 * For each row j, the new value of the parametric constant is equal to
1034
 *
1035
 *  a_j0 - a_jc a_r0/a_rc
1036
 *
1037
 * where a_j0 is the original parametric constant, a_rc is the pivot element,
1038
 * a_r0 is the parametric constant of the pivot row and a_jc is the
1039
 * pivot column entry of the row j.
1040
 * Since a_r0 is non-positive and a_rc is positive, the sign of row j
1041
 * remains the same if a_jc has the same sign as the row j or if
1042
 * a_jc is zero.  In all other cases, we reset the sign to "unknown".
1043
 */
1044
static void update_row_sign(struct isl_tab *tab, int row, int col, int row_sgn)
1045
3.22M
{
1046
3.22M
  int i;
1047
3.22M
  struct isl_mat *mat = tab->mat;
1048
3.22M
  unsigned off = 2 + tab->M;
1049
3.22M
1050
3.22M
  if (!tab->row_sign)
1051
3.19M
    return;
1052
28.1k
1053
28.1k
  if (tab->row_sign[row] == 0)
1054
22.0k
    return;
1055
6.04k
  isl_assert(mat->ctx, row_sgn > 0, return);
1056
6.04k
  isl_assert(mat->ctx, tab->row_sign[row] == isl_tab_row_neg, return);
1057
6.04k
  tab->row_sign[row] = isl_tab_row_pos;
1058
54.4k
  for (i = 0; i < tab->n_row; 
++i48.4k
) {
1059
48.4k
    int s;
1060
48.4k
    if (i == row)
1061
6.04k
      continue;
1062
42.3k
    s = isl_int_sgn(mat->row[i][off + col]);
1063
42.3k
    if (!s)
1064
24.7k
      continue;
1065
17.5k
    if (!tab->row_sign[i])
1066
6.93k
      continue;
1067
10.6k
    if (s < 0 && 
tab->row_sign[i] == isl_tab_row_neg5.94k
)
1068
0
      continue;
1069
10.6k
    if (s > 0 && 
tab->row_sign[i] == isl_tab_row_pos4.68k
)
1070
4.68k
      continue;
1071
5.94k
    tab->row_sign[i] = isl_tab_row_unknown;
1072
5.94k
  }
1073
6.04k
}
1074
1075
/* Given a row number "row" and a column number "col", pivot the tableau
1076
 * such that the associated variables are interchanged.
1077
 * The given row in the tableau expresses
1078
 *
1079
 *  x_r = a_r0 + \sum_i a_ri x_i
1080
 *
1081
 * or
1082
 *
1083
 *  x_c = 1/a_rc x_r - a_r0/a_rc + sum_{i \ne r} -a_ri/a_rc
1084
 *
1085
 * Substituting this equality into the other rows
1086
 *
1087
 *  x_j = a_j0 + \sum_i a_ji x_i
1088
 *
1089
 * with a_jc \ne 0, we obtain
1090
 *
1091
 *  x_j = a_jc/a_rc x_r + a_j0 - a_jc a_r0/a_rc + sum a_ji - a_jc a_ri/a_rc 
1092
 *
1093
 * The tableau
1094
 *
1095
 *  n_rc/d_r    n_ri/d_r
1096
 *  n_jc/d_j    n_ji/d_j
1097
 *
1098
 * where i is any other column and j is any other row,
1099
 * is therefore transformed into
1100
 *
1101
 * s(n_rc)d_r/|n_rc|    -s(n_rc)n_ri/|n_rc|
1102
 * s(n_rc)d_r n_jc/(|n_rc| d_j) (n_ji |n_rc| - s(n_rc)n_jc n_ri)/(|n_rc| d_j)
1103
 *
1104
 * The transformation is performed along the following steps
1105
 *
1106
 *  d_r/n_rc    n_ri/n_rc
1107
 *  n_jc/d_j    n_ji/d_j
1108
 *
1109
 *  s(n_rc)d_r/|n_rc| -s(n_rc)n_ri/|n_rc|
1110
 *  n_jc/d_j    n_ji/d_j
1111
 *
1112
 *  s(n_rc)d_r/|n_rc| -s(n_rc)n_ri/|n_rc|
1113
 *  n_jc/(|n_rc| d_j) n_ji/(|n_rc| d_j)
1114
 *
1115
 *  s(n_rc)d_r/|n_rc| -s(n_rc)n_ri/|n_rc|
1116
 *  n_jc/(|n_rc| d_j) (n_ji |n_rc|)/(|n_rc| d_j)
1117
 *
1118
 *  s(n_rc)d_r/|n_rc| -s(n_rc)n_ri/|n_rc|
1119
 *  n_jc/(|n_rc| d_j) (n_ji |n_rc| - s(n_rc)n_jc n_ri)/(|n_rc| d_j)
1120
 *
1121
 * s(n_rc)d_r/|n_rc|    -s(n_rc)n_ri/|n_rc|
1122
 * s(n_rc)d_r n_jc/(|n_rc| d_j) (n_ji |n_rc| - s(n_rc)n_jc n_ri)/(|n_rc| d_j)
1123
 *
1124
 */
1125
int isl_tab_pivot(struct isl_tab *tab, int row, int col)
1126
3.22M
{
1127
3.22M
  int i, j;
1128
3.22M
  int sgn;
1129
3.22M
  int t;
1130
3.22M
  isl_ctx *ctx;
1131
3.22M
  struct isl_mat *mat = tab->mat;
1132
3.22M
  struct isl_tab_var *var;
1133
3.22M
  unsigned off = 2 + tab->M;
1134
3.22M
1135
3.22M
  ctx = isl_tab_get_ctx(tab);
1136
3.22M
  if (isl_ctx_next_operation(ctx) < 0)
1137
0
    return -1;
1138
3.22M
1139
3.22M
  isl_int_swap(mat->row[row][0], mat->row[row][off + col]);
1140
3.22M
  sgn = isl_int_sgn(mat->row[row][0]);
1141
3.22M
  if (sgn < 0) {
1142
1.84M
    isl_int_neg(mat->row[row][0], mat->row[row][0]);
1143
1.84M
    isl_int_neg(mat->row[row][off + col], mat->row[row][off + col]);
1144
1.84M
  } else
1145
11.4M
    
for (j = 0; 1.37M
j < off - 1 + tab->n_col;
++j10.0M
) {
1146
10.0M
      if (j == off - 1 + col)
1147
1.37M
        continue;
1148
8.65M
      isl_int_neg(mat->row[row][1 + j], mat->row[row][1 + j]);
1149
8.65M
    }
1150
3.22M
  if (!isl_int_is_one(mat->row[row][0]))
1151
3.22M
    
isl_seq_normalize(mat->ctx, mat->row[row], off + tab->n_col)570k
;
1152
40.4M
  for (i = 0; i < tab->n_row; 
++i37.2M
) {
1153
37.2M
    if (i == row)
1154
3.22M
      continue;
1155
33.9M
    if (isl_int_is_zero(mat->row[i][off + col]))
1156
33.9M
      
continue25.0M
;
1157
8.91M
    isl_int_mul(mat->row[i][0], mat->row[i][0], mat->row[row][0]);
1158
94.4M
    for (j = 0; j < off - 1 + tab->n_col; 
++j85.5M
) {
1159
85.5M
      if (j == off - 1 + col)
1160
8.91M
        continue;
1161
76.6M
      isl_int_mul(mat->row[i][1 + j],
1162
76.6M
            mat->row[i][1 + j], mat->row[row][0]);
1163
76.6M
      isl_int_addmul(mat->row[i][1 + j],
1164
76.6M
            mat->row[i][off + col], mat->row[row][1 + j]);
1165
76.6M
    }
1166
8.91M
    isl_int_mul(mat->row[i][off + col],
1167
8.91M
          mat->row[i][off + col], mat->row[row][off + col]);
1168
8.91M
    if (!isl_int_is_one(mat->row[i][0]))
1169
8.91M
      
isl_seq_normalize(mat->ctx, mat->row[i], off + tab->n_col)4.08M
;
1170
8.91M
  }
1171
3.22M
  t = tab->row_var[row];
1172
3.22M
  tab->row_var[row] = tab->col_var[col];
1173
3.22M
  tab->col_var[col] = t;
1174
3.22M
  var = isl_tab_var_from_row(tab, row);
1175
3.22M
  var->is_row = 1;
1176
3.22M
  var->index = row;
1177
3.22M
  var = var_from_col(tab, col);
1178
3.22M
  var->is_row = 0;
1179
3.22M
  var->index = col;
1180
3.22M
  update_row_sign(tab, row, col, sgn);
1181
3.22M
  if (tab->in_undo)
1182
117k
    return 0;
1183
31.9M
  
for (i = tab->n_redundant; 3.10M
i < tab->n_row;
++i28.8M
) {
1184
28.8M
    if (isl_int_is_zero(mat->row[i][off + col]))
1185
28.8M
      
continue18.0M
;
1186
10.7M
    if (!isl_tab_var_from_row(tab, i)->frozen &&
1187
10.7M
        
isl_tab_row_is_redundant(tab, i)10.4M
) {
1188
1.60M
      int redo = isl_tab_mark_redundant(tab, i);
1189
1.60M
      if (redo < 0)
1190
0
        return -1;
1191
1.60M
      if (redo)
1192
149k
        --i;
1193
1.60M
    }
1194
10.7M
  }
1195
3.10M
  return 0;
1196
3.10M
}
1197
1198
/* If "var" represents a column variable, then pivot is up (sgn > 0)
1199
 * or down (sgn < 0) to a row.  The variable is assumed not to be
1200
 * unbounded in the specified direction.
1201
 * If sgn = 0, then the variable is unbounded in both directions,
1202
 * and we pivot with any row we can find.
1203
 */
1204
static int to_row(struct isl_tab *tab, struct isl_tab_var *var, int sign) WARN_UNUSED;
1205
static int to_row(struct isl_tab *tab, struct isl_tab_var *var, int sign)
1206
1.67M
{
1207
1.67M
  int r;
1208
1.67M
  unsigned off = 2 + tab->M;
1209
1.67M
1210
1.67M
  if (var->is_row)
1211
1.43M
    return 0;
1212
240k
1213
240k
  if (sign == 0) {
1214
46.9k
    for (r = tab->n_redundant; r < tab->n_row; 
++r29.0k
)
1215
46.9k
      if (!isl_int_is_zero(tab->mat->row[r][off+var->index]))
1216
46.9k
        
break17.9k
;
1217
17.9k
    isl_assert(tab->mat->ctx, r < tab->n_row, return -1);
1218
222k
  } else {
1219
222k
    r = pivot_row(tab, NULL, sign, var->index);
1220
222k
    isl_assert(tab->mat->ctx, r >= 0, return -1);
1221
222k
  }
1222
240k
1223
240k
  return isl_tab_pivot(tab, r, var->index);
1224
240k
}
1225
1226
/* Check whether all variables that are marked as non-negative
1227
 * also have a non-negative sample value.  This function is not
1228
 * called from the current code but is useful during debugging.
1229
 */
1230
static void check_table(struct isl_tab *tab) __attribute__ ((unused));
1231
static void check_table(struct isl_tab *tab)
1232
0
{
1233
0
  int i;
1234
0
1235
0
  if (tab->empty)
1236
0
    return;
1237
0
  for (i = tab->n_redundant; i < tab->n_row; ++i) {
1238
0
    struct isl_tab_var *var;
1239
0
    var = isl_tab_var_from_row(tab, i);
1240
0
    if (!var->is_nonneg)
1241
0
      continue;
1242
0
    if (tab->M) {
1243
0
      isl_assert(tab->mat->ctx,
1244
0
        !isl_int_is_neg(tab->mat->row[i][2]), abort());
1245
0
      if (isl_int_is_pos(tab->mat->row[i][2]))
1246
0
        continue;
1247
0
    }
1248
0
    isl_assert(tab->mat->ctx, !isl_int_is_neg(tab->mat->row[i][1]),
1249
0
        abort());
1250
0
  }
1251
0
}
1252
1253
/* Return the sign of the maximal value of "var".
1254
 * If the sign is not negative, then on return from this function,
1255
 * the sample value will also be non-negative.
1256
 *
1257
 * If "var" is manifestly unbounded wrt positive values, we are done.
1258
 * Otherwise, we pivot the variable up to a row if needed
1259
 * Then we continue pivoting down until either
1260
 *  - no more down pivots can be performed
1261
 *  - the sample value is positive
1262
 *  - the variable is pivoted into a manifestly unbounded column
1263
 */
1264
static int sign_of_max(struct isl_tab *tab, struct isl_tab_var *var)
1265
1.28M
{
1266
1.28M
  int row, col;
1267
1.28M
1268
1.28M
  if (max_is_manifestly_unbounded(tab, var))
1269
95.0k
    return 1;
1270
1.19M
  if (to_row(tab, var, 1) < 0)
1271
0
    return -2;
1272
2.01M
  
while (1.19M
!isl_int_is_pos(tab->mat->row[var->index][1])) {
1273
1.62M
    find_pivot(tab, var, var, 1, &row, &col);
1274
1.62M
    if (row == -1)
1275
574k
      return isl_int_sgn(tab->mat->row[var->index][1]);
1276
1.05M
    if (isl_tab_pivot(tab, row, col) < 0)
1277
0
      return -2;
1278
1.05M
    if (!var->is_row) /* manifestly unbounded */
1279
226k
      return 1;
1280
1.05M
  }
1281
1.19M
  
return 1389k
;
1282
1.19M
}
1283
1284
int isl_tab_sign_of_max(struct isl_tab *tab, int con)
1285
130
{
1286
130
  struct isl_tab_var *var;
1287
130
1288
130
  if (!tab)
1289
0
    return -2;
1290
130
1291
130
  var = &tab->con[con];
1292
130
  isl_assert(tab->mat->ctx, !var->is_redundant, return -2);
1293
130
  isl_assert(tab->mat->ctx, !var->is_zero, return -2);
1294
130
1295
130
  return sign_of_max(tab, var);
1296
130
}
1297
1298
static int row_is_neg(struct isl_tab *tab, int row)
1299
5.21M
{
1300
5.21M
  if (!tab->M)
1301
5.21M
    return isl_int_is_neg(tab->mat->row[row][1]);
1302
0
  if (isl_int_is_pos(tab->mat->row[row][2]))
1303
0
    return 0;
1304
0
  if (isl_int_is_neg(tab->mat->row[row][2]))
1305
0
    return 1;
1306
0
  return isl_int_is_neg(tab->mat->row[row][1]);
1307
0
}
1308
1309
static int row_sgn(struct isl_tab *tab, int row)
1310
4.63M
{
1311
4.63M
  if (!tab->M)
1312
4.63M
    return isl_int_sgn(tab->mat->row[row][1]);
1313
0
  if (!isl_int_is_zero(tab->mat->row[row][2]))
1314
0
    return isl_int_sgn(tab->mat->row[row][2]);
1315
0
  else
1316
0
    return isl_int_sgn(tab->mat->row[row][1]);
1317
0
}
1318
1319
/* Perform pivots until the row variable "var" has a non-negative
1320
 * sample value or until no more upward pivots can be performed.
1321
 * Return the sign of the sample value after the pivots have been
1322
 * performed.
1323
 */
1324
static int restore_row(struct isl_tab *tab, struct isl_tab_var *var)
1325
4.99M
{
1326
4.99M
  int row, col;
1327
4.99M
1328
5.21M
  while (row_is_neg(tab, var->index)) {
1329
649k
    find_pivot(tab, var, var, 1, &row, &col);
1330
649k
    if (row == -1)
1331
70.6k
      break;
1332
578k
    if (isl_tab_pivot(tab, row, col) < 0)
1333
0
      return -2;
1334
578k
    if (!var->is_row) /* manifestly unbounded */
1335
359k
      return 1;
1336
578k
  }
1337
4.99M
  
return row_sgn(tab, var->index)4.63M
;
1338
4.99M
}
1339
1340
/* Perform pivots until we are sure that the row variable "var"
1341
 * can attain non-negative values.  After return from this
1342
 * function, "var" is still a row variable, but its sample
1343
 * value may not be non-negative, even if the function returns 1.
1344
 */
1345
static int at_least_zero(struct isl_tab *tab, struct isl_tab_var *var)
1346
156k
{
1347
156k
  int row, col;
1348
156k
1349
174k
  while (isl_int_is_neg(tab->mat->row[var->index][1])) {
1350
166k
    find_pivot(tab, var, var, 1, &row, &col);
1351
166k
    if (row == -1)
1352
77.0k
      break;
1353
89.7k
    if (row == var->index) /* manifestly unbounded */
1354
71.3k
      return 1;
1355
18.3k
    if (isl_tab_pivot(tab, row, col) < 0)
1356
0
      return -1;
1357
18.3k
  }
1358
156k
  
return !85.1k
isl_int_is_neg85.1k
(tab->mat->row[var->index][1]);
1359
156k
}
1360
1361
/* Return a negative value if "var" can attain negative values.
1362
 * Return a non-negative value otherwise.
1363
 *
1364
 * If "var" is manifestly unbounded wrt negative values, we are done.
1365
 * Otherwise, if var is in a column, we can pivot it down to a row.
1366
 * Then we continue pivoting down until either
1367
 *  - the pivot would result in a manifestly unbounded column
1368
 *    => we don't perform the pivot, but simply return -1
1369
 *  - no more down pivots can be performed
1370
 *  - the sample value is negative
1371
 * If the sample value becomes negative and the variable is supposed
1372
 * to be nonnegative, then we undo the last pivot.
1373
 * However, if the last pivot has made the pivoting variable
1374
 * obviously redundant, then it may have moved to another row.
1375
 * In that case we look for upward pivots until we reach a non-negative
1376
 * value again.
1377
 */
1378
static int sign_of_min(struct isl_tab *tab, struct isl_tab_var *var)
1379
70.5k
{
1380
70.5k
  int row, col;
1381
70.5k
  struct isl_tab_var *pivot_var = NULL;
1382
70.5k
1383
70.5k
  if (min_is_manifestly_unbounded(tab, var))
1384
0
    return -1;
1385
70.5k
  if (!var->is_row) {
1386
1.86k
    col = var->index;
1387
1.86k
    row = pivot_row(tab, NULL, -1, col);
1388
1.86k
    pivot_var = var_from_col(tab, col);
1389
1.86k
    if (isl_tab_pivot(tab, row, col) < 0)
1390
0
      return -2;
1391
1.86k
    if (var->is_redundant)
1392
317
      return 0;
1393
1.54k
    if (isl_int_is_neg(tab->mat->row[var->index][1])) {
1394
693
      if (var->is_nonneg) {
1395
693
        if (!pivot_var->is_redundant &&
1396
693
            pivot_var->index == row) {
1397
680
          if (isl_tab_pivot(tab, row, col) < 0)
1398
0
            return -2;
1399
13
        } else
1400
13
          if (restore_row(tab, var) < -1)
1401
0
            return -2;
1402
693
      }
1403
693
      return -1;
1404
693
    }
1405
1.54k
  }
1406
69.5k
  if (var->is_redundant)
1407
0
    return 0;
1408
107k
  
while (69.5k
!isl_int_is_neg(tab->mat->row[var->index][1])) {
1409
93.6k
    find_pivot(tab, var, var, -1, &row, &col);
1410
93.6k
    if (row == var->index)
1411
16.2k
      return -1;
1412
77.4k
    if (row == -1)
1413
38.7k
      return isl_int_sgn(tab->mat->row[var->index][1]);
1414
38.6k
    pivot_var = var_from_col(tab, col);
1415
38.6k
    if (isl_tab_pivot(tab, row, col) < 0)
1416
0
      return -2;
1417
38.6k
    if (var->is_redundant)
1418
1.14k
      return 0;
1419
38.6k
  }
1420
69.5k
  
if (13.3k
pivot_var13.3k
&&
var->is_nonneg13.3k
) {
1421
968
    /* pivot back to non-negative value */
1422
968
    if (!pivot_var->is_redundant && pivot_var->index == row) {
1423
955
      if (isl_tab_pivot(tab, row, col) < 0)
1424
0
        return -2;
1425
13
    } else
1426
13
      if (restore_row(tab, var) < -1)
1427
0
        return -2;
1428
13.3k
  }
1429
13.3k
  return -1;
1430
13.3k
}
1431
1432
static int row_at_most_neg_one(struct isl_tab *tab, int row)
1433
277k
{
1434
277k
  if (tab->M) {
1435
0
    if (isl_int_is_pos(tab->mat->row[row][2]))
1436
0
      return 0;
1437
0
    if (isl_int_is_neg(tab->mat->row[row][2]))
1438
0
      return 1;
1439
277k
  }
1440
277k
  return isl_int_is_neg(tab->mat->row[row][1]) &&
1441
277k
         
isl_int_abs_ge159k
(tab->mat->row[row][1],
1442
277k
            tab->mat->row[row][0]);
1443
277k
}
1444
1445
/* Return 1 if "var" can attain values <= -1.
1446
 * Return 0 otherwise.
1447
 *
1448
 * If the variable "var" is supposed to be non-negative (is_nonneg is set),
1449
 * then the sample value of "var" is assumed to be non-negative when the
1450
 * the function is called.  If 1 is returned then the constraint
1451
 * is not redundant and the sample value is made non-negative again before
1452
 * the function returns.
1453
 */
1454
int isl_tab_min_at_most_neg_one(struct isl_tab *tab, struct isl_tab_var *var)
1455
795k
{
1456
795k
  int row, col;
1457
795k
  struct isl_tab_var *pivot_var;
1458
795k
1459
795k
  if (min_is_manifestly_unbounded(tab, var))
1460
188
    return 1;
1461
794k
  if (!var->is_row) {
1462
89.8k
    col = var->index;
1463
89.8k
    row = pivot_row(tab, NULL, -1, col);
1464
89.8k
    pivot_var = var_from_col(tab, col);
1465
89.8k
    if (isl_tab_pivot(tab, row, col) < 0)
1466
0
      return -1;
1467
89.8k
    if (var->is_redundant)
1468
12.3k
      return 0;
1469
77.4k
    if (row_at_most_neg_one(tab, var->index)) {
1470
59.0k
      if (var->is_nonneg) {
1471
59.0k
        if (!pivot_var->is_redundant &&
1472
59.0k
            pivot_var->index == row) {
1473
55.5k
          if (isl_tab_pivot(tab, row, col) < 0)
1474
0
            return -1;
1475
3.53k
        } else
1476
3.53k
          if (restore_row(tab, var) < -1)
1477
0
            return -1;
1478
59.0k
      }
1479
59.0k
      return 1;
1480
59.0k
    }
1481
77.4k
  }
1482
723k
  if (var->is_redundant)
1483
11.5k
    return 0;
1484
837k
  
do 711k
{
1485
837k
    find_pivot(tab, var, var, -1, &row, &col);
1486
837k
    if (row == var->index) {
1487
359k
      if (var->is_nonneg && 
restore_row(tab, var) < -1318k
)
1488
0
        return -1;
1489
359k
      return 1;
1490
359k
    }
1491
477k
    if (row == -1)
1492
179k
      return 0;
1493
297k
    pivot_var = var_from_col(tab, col);
1494
297k
    if (isl_tab_pivot(tab, row, col) < 0)
1495
0
      return -1;
1496
297k
    if (var->is_redundant)
1497
97.5k
      return 0;
1498
200k
  } while (!row_at_most_neg_one(tab, var->index));
1499
711k
  
if (75.3k
var->is_nonneg75.3k
) {
1500
64.1k
    /* pivot back to non-negative value */
1501
64.1k
    if (!pivot_var->is_redundant && pivot_var->index == row)
1502
59.4k
      if (isl_tab_pivot(tab, row, col) < 0)
1503
0
        return -1;
1504
64.1k
    if (restore_row(tab, var) < -1)
1505
0
      return -1;
1506
75.3k
  }
1507
75.3k
  return 1;
1508
75.3k
}
1509
1510
/* Return 1 if "var" can attain values >= 1.
1511
 * Return 0 otherwise.
1512
 */
1513
static int at_least_one(struct isl_tab *tab, struct isl_tab_var *var)
1514
685k
{
1515
685k
  int row, col;
1516
685k
  isl_int *r;
1517
685k
1518
685k
  if (max_is_manifestly_unbounded(tab, var))
1519
320k
    return 1;
1520
365k
  if (to_row(tab, var, 1) < 0)
1521
0
    return -1;
1522
365k
  r = tab->mat->row[var->index];
1523
379k
  while (isl_int_lt(r[1], r[0])) {
1524
14.7k
    find_pivot(tab, var, var, 1, &row, &col);
1525
14.7k
    if (row == -1)
1526
537
      return isl_int_ge(r[1], r[0]);
1527
14.2k
    if (row == var->index) /* manifestly unbounded */
1528
71
      return 1;
1529
14.1k
    if (isl_tab_pivot(tab, row, col) < 0)
1530
0
      return -1;
1531
14.1k
  }
1532
365k
  
return 1364k
;
1533
365k
}
1534
1535
static void swap_cols(struct isl_tab *tab, int col1, int col2)
1536
615k
{
1537
615k
  int t;
1538
615k
  unsigned off = 2 + tab->M;
1539
615k
  t = tab->col_var[col1];
1540
615k
  tab->col_var[col1] = tab->col_var[col2];
1541
615k
  tab->col_var[col2] = t;
1542
615k
  var_from_col(tab, col1)->index = col1;
1543
615k
  var_from_col(tab, col2)->index = col2;
1544
615k
  tab->mat = isl_mat_swap_cols(tab->mat, off + col1, off + col2);
1545
615k
}
1546
1547
/* Mark column with index "col" as representing a zero variable.
1548
 * If we may need to undo the operation the column is kept,
1549
 * but no longer considered.
1550
 * Otherwise, the column is simply removed.
1551
 *
1552
 * The column may be interchanged with some other column.  If it
1553
 * is interchanged with a later column, return 1.  Otherwise return 0.
1554
 * If the columns are checked in order in the calling function,
1555
 * then a return value of 1 means that the column with the given
1556
 * column number may now contain a different column that
1557
 * hasn't been checked yet.
1558
 */
1559
int isl_tab_kill_col(struct isl_tab *tab, int col)
1560
879k
{
1561
879k
  var_from_col(tab, col)->is_zero = 1;
1562
879k
  if (tab->need_undo) {
1563
103k
    if (isl_tab_push_var(tab, isl_tab_undo_zero,
1564
103k
              var_from_col(tab, col)) < 0)
1565
0
      return -1;
1566
103k
    if (col != tab->n_dead)
1567
42.7k
      swap_cols(tab, col, tab->n_dead);
1568
103k
    tab->n_dead++;
1569
103k
    return 0;
1570
776k
  } else {
1571
776k
    if (col != tab->n_col - 1)
1572
572k
      swap_cols(tab, col, tab->n_col - 1);
1573
776k
    var_from_col(tab, tab->n_col - 1)->index = -1;
1574
776k
    tab->n_col--;
1575
776k
    return 1;
1576
776k
  }
1577
879k
}
1578
1579
static int row_is_manifestly_non_integral(struct isl_tab *tab, int row)
1580
3.24M
{
1581
3.24M
  unsigned off = 2 + tab->M;
1582
3.24M
1583
3.24M
  if (tab->M && 
!0
isl_int_eq0
(tab->mat->row[row][2],
1584
3.24M
          tab->mat->row[row][0]))
1585
3.24M
    
return 00
;
1586
3.24M
  if (isl_seq_first_non_zero(tab->mat->row[row] + off + tab->n_dead,
1587
3.24M
            tab->n_col - tab->n_dead) != -1)
1588
257k
    return 0;
1589
2.99M
1590
2.99M
  return !isl_int_is_divisible_by(tab->mat->row[row][1],
1591
2.99M
          tab->mat->row[row][0]);
1592
2.99M
}
1593
1594
/* For integer tableaus, check if any of the coordinates are stuck
1595
 * at a non-integral value.
1596
 */
1597
static int tab_is_manifestly_empty(struct isl_tab *tab)
1598
574k
{
1599
574k
  int i;
1600
574k
1601
574k
  if (tab->empty)
1602
0
    return 1;
1603
574k
  if (tab->rational)
1604
12.4k
    return 0;
1605
561k
1606
7.00M
  
for (i = 0; 561k
i < tab->n_var;
++i6.44M
) {
1607
6.44M
    if (!tab->var[i].is_row)
1608
3.19M
      continue;
1609
3.24M
    if (row_is_manifestly_non_integral(tab, tab->var[i].index))
1610
28
      return 1;
1611
3.24M
  }
1612
561k
1613
561k
  
return 0561k
;
1614
561k
}
1615
1616
/* Row variable "var" is non-negative and cannot attain any values
1617
 * larger than zero.  This means that the coefficients of the unrestricted
1618
 * column variables are zero and that the coefficients of the non-negative
1619
 * column variables are zero or negative.
1620
 * Each of the non-negative variables with a negative coefficient can
1621
 * then also be written as the negative sum of non-negative variables
1622
 * and must therefore also be zero.
1623
 *
1624
 * If "temp_var" is set, then "var" is a temporary variable that
1625
 * will be removed after this function returns and for which
1626
 * no information is recorded on the undo stack.
1627
 * Do not add any undo records involving this variable in this case
1628
 * since the variable will have been removed before any future undo
1629
 * operations.  Also avoid marking the variable as redundant,
1630
 * since that either adds an undo record or needlessly removes the row
1631
 * (the caller will take care of removing the row).
1632
 */
1633
static isl_stat close_row(struct isl_tab *tab, struct isl_tab_var *var,
1634
  int temp_var) WARN_UNUSED;
1635
static isl_stat close_row(struct isl_tab *tab, struct isl_tab_var *var,
1636
  int temp_var)
1637
574k
{
1638
574k
  int j;
1639
574k
  struct isl_mat *mat = tab->mat;
1640
574k
  unsigned off = 2 + tab->M;
1641
574k
1642
574k
  if (!var->is_nonneg)
1643
574k
    
isl_die0
(isl_tab_get_ctx(tab), isl_error_internal,
1644
574k
      "expecting non-negative variable",
1645
574k
      return isl_stat_error);
1646
574k
  var->is_zero = 1;
1647
574k
  if (!temp_var && 
tab->need_undo561k
)
1648
424
    if (isl_tab_push_var(tab, isl_tab_undo_zero, var) < 0)
1649
0
      return isl_stat_error;
1650
4.47M
  
for (j = tab->n_dead; 574k
j < tab->n_col;
++j3.90M
) {
1651
3.90M
    int recheck;
1652
3.90M
    if (isl_int_is_zero(mat->row[var->index][off + j]))
1653
3.90M
      
continue3.39M
;
1654
503k
    if (isl_int_is_pos(mat->row[var->index][off + j]))
1655
503k
      
isl_die0
(isl_tab_get_ctx(tab), isl_error_internal,
1656
503k
        "row cannot have positive coefficients",
1657
503k
        return isl_stat_error);
1658
503k
    recheck = isl_tab_kill_col(tab, j);
1659
503k
    if (recheck < 0)
1660
0
      return isl_stat_error;
1661
503k
    if (recheck)
1662
490k
      --j;
1663
503k
  }
1664
574k
  if (!temp_var && 
isl_tab_mark_redundant(tab, var->index) < 0561k
)
1665
0
    return isl_stat_error;
1666
574k
  if (tab_is_manifestly_empty(tab) && 
isl_tab_mark_empty(tab) < 028
)
1667
0
    return isl_stat_error;
1668
574k
  return isl_stat_ok;
1669
574k
}
1670
1671
/* Add a constraint to the tableau and allocate a row for it.
1672
 * Return the index into the constraint array "con".
1673
 *
1674
 * This function assumes that at least one more row and at least
1675
 * one more element in the constraint array are available in the tableau.
1676
 */
1677
int isl_tab_allocate_con(struct isl_tab *tab)
1678
5.20M
{
1679
5.20M
  int r;
1680
5.20M
1681
5.20M
  isl_assert(tab->mat->ctx, tab->n_row < tab->mat->n_row, return -1);
1682
5.20M
  isl_assert(tab->mat->ctx, tab->n_con < tab->max_con, return -1);
1683
5.20M
1684
5.20M
  r = tab->n_con;
1685
5.20M
  tab->con[r].index = tab->n_row;
1686
5.20M
  tab->con[r].is_row = 1;
1687
5.20M
  tab->con[r].is_nonneg = 0;
1688
5.20M
  tab->con[r].is_zero = 0;
1689
5.20M
  tab->con[r].is_redundant = 0;
1690
5.20M
  tab->con[r].frozen = 0;
1691
5.20M
  tab->con[r].negated = 0;
1692
5.20M
  tab->row_var[tab->n_row] = ~r;
1693
5.20M
1694
5.20M
  tab->n_row++;
1695
5.20M
  tab->n_con++;
1696
5.20M
  if (isl_tab_push_var(tab, isl_tab_undo_allocate, &tab->con[r]) < 0)
1697
0
    return -1;
1698
5.20M
1699
5.20M
  return r;
1700
5.20M
}
1701
1702
/* Move the entries in tab->var up one position, starting at "first",
1703
 * creating room for an extra entry at position "first".
1704
 * Since some of the entries of tab->row_var and tab->col_var contain
1705
 * indices into this array, they have to be updated accordingly.
1706
 */
1707
static int var_insert_entry(struct isl_tab *tab, int first)
1708
5.84k
{
1709
5.84k
  int i;
1710
5.84k
1711
5.84k
  if (tab->n_var >= tab->max_var)
1712
5.84k
    
isl_die0
(isl_tab_get_ctx(tab), isl_error_internal,
1713
5.84k
      "not enough room for new variable", return -1);
1714
5.84k
  if (first > tab->n_var)
1715
5.84k
    
isl_die0
(isl_tab_get_ctx(tab), isl_error_internal,
1716
5.84k
      "invalid initial position", return -1);
1717
5.84k
1718
6.65k
  
for (i = tab->n_var - 1; 5.84k
i >= first;
--i817
) {
1719
817
    tab->var[i + 1] = tab->var[i];
1720
817
    if (tab->var[i + 1].is_row)
1721
547
      tab->row_var[tab->var[i + 1].index]++;
1722
270
    else
1723
270
      tab->col_var[tab->var[i + 1].index]++;
1724
817
  }
1725
5.84k
1726
5.84k
  tab->n_var++;
1727
5.84k
1728
5.84k
  return 0;
1729
5.84k
}
1730
1731
/* Drop the entry at position "first" in tab->var, moving all
1732
 * subsequent entries down.
1733
 * Since some of the entries of tab->row_var and tab->col_var contain
1734
 * indices into this array, they have to be updated accordingly.
1735
 */
1736
static int var_drop_entry(struct isl_tab *tab, int first)
1737
3.43k
{
1738
3.43k
  int i;
1739
3.43k
1740
3.43k
  if (first >= tab->n_var)
1741
3.43k
    
isl_die0
(isl_tab_get_ctx(tab), isl_error_internal,
1742
3.43k
      "invalid initial position", return -1);
1743
3.43k
1744
3.43k
  tab->n_var--;
1745
3.43k
1746
3.97k
  for (i = first; i < tab->n_var; 
++i533
) {
1747
533
    tab->var[i] = tab->var[i + 1];
1748
533
    if (tab->var[i + 1].is_row)
1749
527
      tab->row_var[tab->var[i].index]--;
1750
6
    else
1751
6
      tab->col_var[tab->var[i].index]--;
1752
533
  }
1753
3.43k
1754
3.43k
  return 0;
1755
3.43k
}
1756
1757
/* Add a variable to the tableau at position "r" and allocate a column for it.
1758
 * Return the index into the variable array "var", i.e., "r",
1759
 * or -1 on error.
1760
 */
1761
int isl_tab_insert_var(struct isl_tab *tab, int r)
1762
5.84k
{
1763
5.84k
  int i;
1764
5.84k
  unsigned off = 2 + tab->M;
1765
5.84k
1766
5.84k
  isl_assert(tab->mat->ctx, tab->n_col < tab->mat->n_col, return -1);
1767
5.84k
1768
5.84k
  if (var_insert_entry(tab, r) < 0)
1769
0
    return -1;
1770
5.84k
1771
5.84k
  tab->var[r].index = tab->n_col;
1772
5.84k
  tab->var[r].is_row = 0;
1773
5.84k
  tab->var[r].is_nonneg = 0;
1774
5.84k
  tab->var[r].is_zero = 0;
1775
5.84k
  tab->var[r].is_redundant = 0;
1776
5.84k
  tab->var[r].frozen = 0;
1777
5.84k
  tab->var[r].negated = 0;
1778
5.84k
  tab->col_var[tab->n_col] = r;
1779
5.84k
1780
27.6k
  for (i = 0; i < tab->n_row; 
++i21.8k
)
1781
21.8k
    isl_int_set_si(tab->mat->row[i][off + tab->n_col], 0);
1782
5.84k
1783
5.84k
  tab->n_col++;
1784
5.84k
  if (isl_tab_push_var(tab, isl_tab_undo_allocate, &tab->var[r]) < 0)
1785
0
    return -1;
1786
5.84k
1787
5.84k
  return r;
1788
5.84k
}
1789
1790
/* Add a variable to the tableau and allocate a column for it.
1791
 * Return the index into the variable array "var".
1792
 */
1793
int isl_tab_allocate_var(struct isl_tab *tab)
1794
0
{
1795
0
  if (!tab)
1796
0
    return -1;
1797
0
1798
0
  return isl_tab_insert_var(tab, tab->n_var);
1799
0
}
1800
1801
/* Add a row to the tableau.  The row is given as an affine combination
1802
 * of the original variables and needs to be expressed in terms of the
1803
 * column variables.
1804
 *
1805
 * This function assumes that at least one more row and at least
1806
 * one more element in the constraint array are available in the tableau.
1807
 *
1808
 * We add each term in turn.
1809
 * If r = n/d_r is the current sum and we need to add k x, then
1810
 *  if x is a column variable, we increase the numerator of
1811
 *    this column by k d_r
1812
 *  if x = f/d_x is a row variable, then the new representation of r is
1813
 *
1814
 *     n    k f   d_x/g n + d_r/g k f   m/d_r n + m/d_g k f
1815
 *    --- + --- = ------------------- = -------------------
1816
 *    d_r   d_r        d_r d_x/g                m
1817
 *
1818
 *  with g the gcd of d_r and d_x and m the lcm of d_r and d_x.
1819
 *
1820
 * If tab->M is set, then, internally, each variable x is represented
1821
 * as x' - M.  We then also need no subtract k d_r from the coefficient of M.
1822
 */
1823
int isl_tab_add_row(struct isl_tab *tab, isl_int *line)
1824
5.19M
{
1825
5.19M
  int i;
1826
5.19M
  int r;
1827
5.19M
  isl_int *row;
1828
5.19M
  isl_int a, b;
1829
5.19M
  unsigned off = 2 + tab->M;
1830
5.19M
1831
5.19M
  r = isl_tab_allocate_con(tab);
1832
5.19M
  if (r < 0)
1833
0
    return -1;
1834
5.19M
1835
5.19M
  isl_int_init(a);
1836
5.19M
  isl_int_init(b);
1837
5.19M
  row = tab->mat->row[tab->con[r].index];
1838
5.19M
  isl_int_set_si(row[0], 1);
1839
5.19M
  isl_int_set(row[1], line[0]);
1840
5.19M
  isl_seq_clr(row + 2, tab->M + tab->n_col);
1841
47.4M
  for (i = 0; i < tab->n_var; 
++i42.2M
) {
1842
42.2M
    if (tab->var[i].is_zero)
1843
0
      continue;
1844
42.2M
    if (tab->var[i].is_row) {
1845
7.73M
      isl_int_lcm(a,
1846
7.73M
        row[0], tab->mat->row[tab->var[i].index][0]);
1847
7.73M
      isl_int_swap(a, row[0]);
1848
7.73M
      isl_int_divexact(a, row[0], a);
1849
7.73M
      isl_int_divexact(b,
1850
7.73M
        row[0], tab->mat->row[tab->var[i].index][0]);
1851
7.73M
      isl_int_mul(b, b, line[1 + i]);
1852
7.73M
      isl_seq_combine(row + 1, a, row + 1,
1853
7.73M
          b, tab->mat->row[tab->var[i].index] + 1,
1854
7.73M
          1 + tab->M + tab->n_col);
1855
7.73M
    } else
1856
42.2M
      
isl_int_addmul34.4M
(row[off + tab->var[i].index],
1857
42.2M
              line[1 + i], row[0]);
1858
42.2M
    if (tab->M && 
i >= tab->n_param334k
&&
i < tab->n_var - tab->n_div145k
)
1859
42.2M
      
isl_int_submul143k
(row[2], line[1 + i], row[0]);
1860
42.2M
  }
1861
5.19M
  isl_seq_normalize(tab->mat->ctx, row, off + tab->n_col);
1862
5.19M
  isl_int_clear(a);
1863
5.19M
  isl_int_clear(b);
1864
5.19M
1865
5.19M
  if (tab->row_sign)
1866
40.9k
    tab->row_sign[tab->con[r].index] = isl_tab_row_unknown;
1867
5.19M
1868
5.19M
  return r;
1869
5.19M
}
1870
1871
static isl_stat drop_row(struct isl_tab *tab, int row)
1872
1.07M
{
1873
1.07M
  isl_assert(tab->mat->ctx, ~tab->row_var[row] == tab->n_con - 1,
1874
1.07M
    return isl_stat_error);
1875
1.07M
  if (row != tab->n_row - 1)
1876
124k
    swap_rows(tab, row, tab->n_row - 1);
1877
1.07M
  tab->n_row--;
1878
1.07M
  tab->n_con--;
1879
1.07M
  return isl_stat_ok;
1880
1.07M
}
1881
1882
/* Drop the variable in column "col" along with the column.
1883
 * The column is removed first because it may need to be moved
1884
 * into the last position and this process requires
1885
 * the contents of the col_var array in a state
1886
 * before the removal of the variable.
1887
 */
1888
static isl_stat drop_col(struct isl_tab *tab, int col)
1889
3.43k
{
1890
3.43k
  int var;
1891
3.43k
1892
3.43k
  var = tab->col_var[col];
1893
3.43k
  if (col != tab->n_col - 1)
1894
641
    swap_cols(tab, col, tab->n_col - 1);
1895
3.43k
  tab->n_col--;
1896
3.43k
  if (var_drop_entry(tab, var) < 0)
1897
0
    return isl_stat_error;
1898
3.43k
  return isl_stat_ok;
1899
3.43k
}
1900
1901
/* Add inequality "ineq" and check if it conflicts with the
1902
 * previously added constraints or if it is obviously redundant.
1903
 *
1904
 * This function assumes that at least one more row and at least
1905
 * one more element in the constraint array are available in the tableau.
1906
 */
1907
isl_stat isl_tab_add_ineq(struct isl_tab *tab, isl_int *ineq)
1908
4.05M
{
1909
4.05M
  int r;
1910
4.05M
  int sgn;
1911
4.05M
  isl_int cst;
1912
4.05M
1913
4.05M
  if (!tab)
1914
0
    return isl_stat_error;
1915
4.05M
  if (tab->bmap) {
1916
339k
    struct isl_basic_map *bmap = tab->bmap;
1917
339k
1918
339k
    isl_assert(tab->mat->ctx, tab->n_eq == bmap->n_eq,
1919
339k
      return isl_stat_error);
1920
339k
    isl_assert(tab->mat->ctx,
1921
339k
          tab->n_con == bmap->n_eq + bmap->n_ineq,
1922
339k
          return isl_stat_error);
1923
339k
    tab->bmap = isl_basic_map_add_ineq(tab->bmap, ineq);
1924
339k
    if (isl_tab_push(tab, isl_tab_undo_bmap_ineq) < 0)
1925
0
      return isl_stat_error;
1926
339k
    if (!tab->bmap)
1927
0
      return isl_stat_error;
1928
4.05M
  }
1929
4.05M
  if (tab->cone) {
1930
3.50k
    isl_int_init(cst);
1931
3.50k
    isl_int_set_si(cst, 0);
1932
3.50k
    isl_int_swap(ineq[0], cst);
1933
3.50k
  }
1934
4.05M
  r = isl_tab_add_row(tab, ineq);
1935
4.05M
  if (tab->cone) {
1936
3.50k
    isl_int_swap(ineq[0], cst);
1937
3.50k
    isl_int_clear(cst);
1938
3.50k
  }
1939
4.05M
  if (r < 0)
1940
0
    return isl_stat_error;
1941
4.05M
  tab->con[r].is_nonneg = 1;
1942
4.05M
  if (isl_tab_push_var(tab, isl_tab_undo_nonneg, &tab->con[r]) < 0)
1943
0
    return isl_stat_error;
1944
4.05M
  if (isl_tab_row_is_redundant(tab, tab->con[r].index)) {
1945
106k
    if (isl_tab_mark_redundant(tab, tab->con[r].index) < 0)
1946
0
      return isl_stat_error;
1947
106k
    return isl_stat_ok;
1948
106k
  }
1949
3.94M
1950
3.94M
  sgn = restore_row(tab, &tab->con[r]);
1951
3.94M
  if (sgn < -1)
1952
0
    return isl_stat_error;
1953
3.94M
  if (sgn < 0)
1954
70.6k
    return isl_tab_mark_empty(tab);
1955
3.87M
  if (tab->con[r].is_row && 
isl_tab_row_is_redundant(tab, tab->con[r].index)3.51M
)
1956
0
    if (isl_tab_mark_redundant(tab, tab->con[r].index) < 0)
1957
0
      return isl_stat_error;
1958
3.87M
  return isl_stat_ok;
1959
3.87M
}
1960
1961
/* Pivot a non-negative variable down until it reaches the value zero
1962
 * and then pivot the variable into a column position.
1963
 */
1964
static int to_col(struct isl_tab *tab, struct isl_tab_var *var) WARN_UNUSED;
1965
static int to_col(struct isl_tab *tab, struct isl_tab_var *var)
1966
90.5k
{
1967
90.5k
  int i;
1968
90.5k
  int row, col;
1969
90.5k
  unsigned off = 2 + tab->M;
1970
90.5k
1971
90.5k
  if (!var->is_row)
1972
108
    return 0;
1973
90.3k
1974
111k
  
while (90.3k
isl_int_is_pos(tab->mat->row[var->index][1])) {
1975
107k
    find_pivot(tab, var, NULL, -1, &row, &col);
1976
107k
    isl_assert(tab->mat->ctx, row != -1, return -1);
1977
107k
    if (isl_tab_pivot(tab, row, col) < 0)
1978
0
      return -1;
1979
107k
    if (!var->is_row)
1980
86.6k
      return 0;
1981
107k
  }
1982
90.3k
1983
90.3k
  
for (i = tab->n_dead; 3.72k
i < tab->n_col6.15k
;
++i2.43k
)
1984
6.15k
    if (!isl_int_is_zero(tab->mat->row[var->index][off + i]))
1985
6.15k
      
break3.72k
;
1986
3.72k
1987
3.72k
  isl_assert(tab->mat->ctx, i < tab->n_col, return -1);
1988
3.72k
  if (isl_tab_pivot(tab, var->index, i) < 0)
1989
0
    return -1;
1990
3.72k
1991
3.72k
  return 0;
1992
3.72k
}
1993
1994
/* We assume Gaussian elimination has been performed on the equalities.
1995
 * The equalities can therefore never conflict.
1996
 * Adding the equalities is currently only really useful for a later call
1997
 * to isl_tab_ineq_type.
1998
 *
1999
 * This function assumes that at least one more row and at least
2000
 * one more element in the constraint array are available in the tableau.
2001
 */
2002
static struct isl_tab *add_eq(struct isl_tab *tab, isl_int *eq)
2003
261k
{
2004
261k
  int i;
2005
261k
  int r;
2006
261k
2007
261k
  if (!tab)
2008
0
    return NULL;
2009
261k
  r = isl_tab_add_row(tab, eq);
2010
261k
  if (r < 0)
2011
0
    goto error;
2012
261k
2013
261k
  r = tab->con[r].index;
2014
261k
  i = isl_seq_first_non_zero(tab->mat->row[r] + 2 + tab->M + tab->n_dead,
2015
261k
          tab->n_col - tab->n_dead);
2016
261k
  isl_assert(tab->mat->ctx, i >= 0, goto error);
2017
261k
  i += tab->n_dead;
2018
261k
  if (isl_tab_pivot(tab, r, i) < 0)
2019
0
    goto error;
2020
261k
  if (isl_tab_kill_col(tab, i) < 0)
2021
0
    goto error;
2022
261k
  tab->n_eq++;
2023
261k
2024
261k
  return tab;
2025
0
error:
2026
0
  isl_tab_free(tab);
2027
0
  return NULL;
2028
261k
}
2029
2030
/* Does the sample value of row "row" of "tab" involve the big parameter,
2031
 * if any?
2032
 */
2033
static int row_is_big(struct isl_tab *tab, int row)
2034
79.7k
{
2035
79.7k
  return tab->M && 
!0
isl_int_is_zero0
(tab->mat->row[row][2]);
2036
79.7k
}
2037
2038
static int row_is_manifestly_zero(struct isl_tab *tab, int row)
2039
101k
{
2040
101k
  unsigned off = 2 + tab->M;
2041
101k
2042
101k
  if (!isl_int_is_zero(tab->mat->row[row][1]))
2043
101k
    
return 087.0k
;
2044
14.9k
  if (row_is_big(tab, row))
2045
0
    return 0;
2046
14.9k
  return isl_seq_first_non_zero(tab->mat->row[row] + off + tab->n_dead,
2047
14.9k
          tab->n_col - tab->n_dead) == -1;
2048
14.9k
}
2049
2050
/* Add an equality that is known to be valid for the given tableau.
2051
 *
2052
 * This function assumes that at least one more row and at least
2053
 * one more element in the constraint array are available in the tableau.
2054
 */
2055
int isl_tab_add_valid_eq(struct isl_tab *tab, isl_int *eq)
2056
87.5k
{
2057
87.5k
  struct isl_tab_var *var;
2058
87.5k
  int r;
2059
87.5k
2060
87.5k
  if (!tab)
2061
0
    return -1;
2062
87.5k
  r = isl_tab_add_row(tab, eq);
2063
87.5k
  if (r < 0)
2064
0
    return -1;
2065
87.5k
2066
87.5k
  var = &tab->con[r];
2067
87.5k
  r = var->index;
2068
87.5k
  if (row_is_manifestly_zero(tab, r)) {
2069
798
    var->is_zero = 1;
2070
798
    if (isl_tab_mark_redundant(tab, r) < 0)
2071
0
      return -1;
2072
798
    return 0;
2073
798
  }
2074
86.7k
2075
86.7k
  if (isl_int_is_neg(tab->mat->row[r][1])) {
2076
33.8k
    isl_seq_neg(tab->mat->row[r] + 1, tab->mat->row[r] + 1,
2077
33.8k
          1 + tab->n_col);
2078
33.8k
    var->negated = 1;
2079
33.8k
  }
2080
86.7k
  var->is_nonneg = 1;
2081
86.7k
  if (to_col(tab, var) < 0)
2082
0
    return -1;
2083
86.7k
  var->is_nonneg = 0;
2084
86.7k
  if (isl_tab_kill_col(tab, var->index) < 0)
2085
0
    return -1;
2086
86.7k
2087
86.7k
  return 0;
2088
86.7k
}
2089
2090
/* Add a zero row to "tab" and return the corresponding index
2091
 * in the constraint array.
2092
 *
2093
 * This function assumes that at least one more row and at least
2094
 * one more element in the constraint array are available in the tableau.
2095
 */
2096
static int add_zero_row(struct isl_tab *tab)
2097
2.79k
{
2098
2.79k
  int r;
2099
2.79k
  isl_int *row;
2100
2.79k
2101
2.79k
  r = isl_tab_allocate_con(tab);
2102
2.79k
  if (r < 0)
2103
0
    return -1;
2104
2.79k
2105
2.79k
  row = tab->mat->row[tab->con[r].index];
2106
2.79k
  isl_seq_clr(row + 1, 1 + tab->M + tab->n_col);
2107
2.79k
  isl_int_set_si(row[0], 1);
2108
2.79k
2109
2.79k
  return r;
2110
2.79k
}
2111
2112
/* Add equality "eq" and check if it conflicts with the
2113
 * previously added constraints or if it is obviously redundant.
2114
 *
2115
 * This function assumes that at least one more row and at least
2116
 * one more element in the constraint array are available in the tableau.
2117
 * If tab->bmap is set, then two rows are needed instead of one.
2118
 */
2119
int isl_tab_add_eq(struct isl_tab *tab, isl_int *eq)
2120
14.3k
{
2121
14.3k
  struct isl_tab_undo *snap = NULL;
2122
14.3k
  struct isl_tab_var *var;
2123
14.3k
  int r;
2124
14.3k
  int row;
2125
14.3k
  int sgn;
2126
14.3k
  isl_int cst;
2127
14.3k
2128
14.3k
  if (!tab)
2129
0
    return -1;
2130
14.3k
  isl_assert(tab->mat->ctx, !tab->M, return -1);
2131
14.3k
2132
14.3k
  if (tab->need_undo)
2133
13.6k
    snap = isl_tab_snap(tab);
2134
14.3k
2135
14.3k
  if (tab->cone) {
2136
1.20k
    isl_int_init(cst);
2137
1.20k
    isl_int_set_si(cst, 0);
2138
1.20k
    isl_int_swap(eq[0], cst);
2139
1.20k
  }
2140
14.3k
  r = isl_tab_add_row(tab, eq);
2141
14.3k
  if (tab->cone) {
2142
1.20k
    isl_int_swap(eq[0], cst);
2143
1.20k
    isl_int_clear(cst);
2144
1.20k
  }
2145
14.3k
  if (r < 0)
2146
0
    return -1;
2147
14.3k
2148
14.3k
  var = &tab->con[r];
2149
14.3k
  row = var->index;
2150
14.3k
  if (row_is_manifestly_zero(tab, row)) {
2151
10.6k
    if (snap)
2152
10.5k
      return isl_tab_rollback(tab, snap);
2153
83
    return drop_row(tab, row);
2154
83
  }
2155
3.72k
2156
3.72k
  if (tab->bmap) {
2157
2.79k
    tab->bmap = isl_basic_map_add_ineq(tab->bmap, eq);
2158
2.79k
    if (isl_tab_push(tab, isl_tab_undo_bmap_ineq) < 0)
2159
0
      return -1;
2160
2.79k
    isl_seq_neg(eq, eq, 1 + tab->n_var);
2161
2.79k
    tab->bmap = isl_basic_map_add_ineq(tab->bmap, eq);
2162
2.79k
    isl_seq_neg(eq, eq, 1 + tab->n_var);
2163
2.79k
    if (isl_tab_push(tab, isl_tab_undo_bmap_ineq) < 0)
2164
0
      return -1;
2165
2.79k
    if (!tab->bmap)
2166
0
      return -1;
2167
2.79k
    if (add_zero_row(tab) < 0)
2168
0
      return -1;
2169
3.72k
  }
2170
3.72k
2171
3.72k
  sgn = isl_int_sgn(tab->mat->row[row][1]);
2172
3.72k
2173
3.72k
  if (sgn > 0) {
2174
115
    isl_seq_neg(tab->mat->row[row] + 1, tab->mat->row[row] + 1,
2175
115
          1 + tab->n_col);
2176
115
    var->negated = 1;
2177
115
    sgn = -1;
2178
115
  }
2179
3.72k
2180
3.72k
  if (sgn < 0) {
2181
2.43k
    sgn = sign_of_max(tab, var);
2182
2.43k
    if (sgn < -1)
2183
0
      return -1;
2184
2.43k
    if (sgn < 0) {
2185
0
      if (isl_tab_mark_empty(tab) < 0)
2186
0
        return -1;
2187
0
      return 0;
2188
0
    }
2189
2.43k
  }
2190
3.72k
2191
3.72k
  var->is_nonneg = 1;
2192
3.72k
  if (to_col(tab, var) < 0)
2193
0
    return -1;
2194
3.72k
  var->is_nonneg = 0;
2195
3.72k
  if (isl_tab_kill_col(tab, var->index) < 0)
2196
0
    return -1;
2197
3.72k
2198
3.72k
  return 0;
2199
3.72k
}
2200
2201
/* Construct and return an inequality that expresses an upper bound
2202
 * on the given div.
2203
 * In particular, if the div is given by
2204
 *
2205
 *  d = floor(e/m)
2206
 *
2207
 * then the inequality expresses
2208
 *
2209
 *  m d <= e
2210
 */
2211
static struct isl_vec *ineq_for_div(struct isl_basic_map *bmap, unsigned div)
2212
1.97k
{
2213
1.97k
  unsigned total;
2214
1.97k
  unsigned div_pos;
2215
1.97k
  struct isl_vec *ineq;
2216
1.97k
2217
1.97k
  if (!bmap)
2218
0
    return NULL;
2219
1.97k
2220
1.97k
  total = isl_basic_map_total_dim(bmap);
2221
1.97k
  div_pos = 1 + total - bmap->n_div + div;
2222
1.97k
2223
1.97k
  ineq = isl_vec_alloc(bmap->ctx, 1 + total);
2224
1.97k
  if (!ineq)
2225
0
    return NULL;
2226
1.97k
2227
1.97k
  isl_seq_cpy(ineq->el, bmap->div[div] + 1, 1 + total);
2228
1.97k
  isl_int_neg(ineq->el[div_pos], bmap->div[div][0]);
2229
1.97k
  return ineq;
2230
1.97k
}
2231
2232
/* For a div d = floor(f/m), add the constraints
2233
 *
2234
 *    f - m d >= 0
2235
 *    -(f-(m-1)) + m d >= 0
2236
 *
2237
 * Note that the second constraint is the negation of
2238
 *
2239
 *    f - m d >= m
2240
 *
2241
 * If add_ineq is not NULL, then this function is used
2242
 * instead of isl_tab_add_ineq to effectively add the inequalities.
2243
 *
2244
 * This function assumes that at least two more rows and at least
2245
 * two more elements in the constraint array are available in the tableau.
2246
 */
2247
static isl_stat add_div_constraints(struct isl_tab *tab, unsigned div,
2248
  isl_stat (*add_ineq)(void *user, isl_int *), void *user)
2249
1.97k
{
2250
1.97k
  unsigned total;
2251
1.97k
  unsigned div_pos;
2252
1.97k
  struct isl_vec *ineq;
2253
1.97k
2254
1.97k
  total = isl_basic_map_total_dim(tab->bmap);
2255
1.97k
  div_pos = 1 + total - tab->bmap->n_div + div;
2256
1.97k
2257
1.97k
  ineq = ineq_for_div(tab->bmap, div);
2258
1.97k
  if (!ineq)
2259
0
    goto error;
2260
1.97k
2261
1.97k
  if (add_ineq) {
2262
601
    if (add_ineq(user, ineq->el) < 0)
2263
0
      goto error;
2264
1.37k
  } else {
2265
1.37k
    if (isl_tab_add_ineq(tab, ineq->el) < 0)
2266
0
      goto error;
2267
1.97k
  }
2268
1.97k
2269
1.97k
  isl_seq_neg(ineq->el, tab->bmap->div[div] + 1, 1 + total);
2270
1.97k
  isl_int_set(ineq->el[div_pos], tab->bmap->div[div][0]);
2271
1.97k
  isl_int_add(ineq->el[0], ineq->el[0], ineq->el[div_pos]);
2272
1.97k
  isl_int_sub_ui(ineq->el[0], ineq->el[0], 1);
2273
1.97k
2274
1.97k
  if (add_ineq) {
2275
601
    if (add_ineq(user, ineq->el) < 0)
2276
0
      goto error;
2277
1.37k
  } else {
2278
1.37k
    if (isl_tab_add_ineq(tab, ineq->el) < 0)
2279
0
      goto error;
2280
1.97k
  }
2281
1.97k
2282
1.97k
  isl_vec_free(ineq);
2283
1.97k
2284
1.97k
  return 0;
2285
0
error:
2286
0
  isl_vec_free(ineq);
2287
0
  return -1;
2288
1.97k
}
2289
2290
/* Check whether the div described by "div" is obviously non-negative.
2291
 * If we are using a big parameter, then we will encode the div
2292
 * as div' = M + div, which is always non-negative.
2293
 * Otherwise, we check whether div is a non-negative affine combination
2294
 * of non-negative variables.
2295
 */
2296
static int div_is_nonneg(struct isl_tab *tab, __isl_keep isl_vec *div)
2297
1.97k
{
2298
1.97k
  int i;
2299
1.97k
2300
1.97k
  if (tab->M)
2301
0
    return 1;
2302
1.97k
2303
1.97k
  if (isl_int_is_neg(div->el[1]))
2304
1.97k
    
return 0225
;
2305
1.74k
2306
3.94k
  
for (i = 0; 1.74k
i < tab->n_var;
++i2.19k
) {
2307
3.66k
    if (isl_int_is_neg(div->el[2 + i]))
2308
3.66k
      
return 0330
;
2309
3.33k
    if (isl_int_is_zero(div->el[2 + i]))
2310
3.33k
      
continue1.87k
;
2311
1.45k
    if (!tab->var[i].is_nonneg)
2312
1.13k
      return 0;
2313
1.45k
  }
2314
1.74k
2315
1.74k
  
return 1281
;
2316
1.74k
}
2317
2318
/* Insert an extra div, prescribed by "div", to the tableau and
2319
 * the associated bmap (which is assumed to be non-NULL).
2320
 * The extra integer division is inserted at (tableau) position "pos".
2321
 * Return "pos" or -1 if an error occurred.
2322
 *
2323
 * If add_ineq is not NULL, then this function is used instead
2324
 * of isl_tab_add_ineq to add the div constraints.
2325
 * This complication is needed because the code in isl_tab_pip
2326
 * wants to perform some extra processing when an inequality
2327
 * is added to the tableau.
2328
 */
2329
int isl_tab_insert_div(struct isl_tab *tab, int pos, __isl_keep isl_vec *div,
2330
  isl_stat (*add_ineq)(void *user, isl_int *), void *user)
2331
1.97k
{
2332
1.97k
  int r;
2333
1.97k
  int nonneg;
2334
1.97k
  int n_div, o_div;
2335
1.97k
2336
1.97k
  if (!tab || !div)
2337
0
    return -1;
2338
1.97k
2339
1.97k
  if (div->size != 1 + 1 + tab->n_var)
2340
1.97k
    
isl_die0
(isl_tab_get_ctx(tab), isl_error_invalid,
2341
1.97k
      "unexpected size", return -1);
2342
1.97k
2343
1.97k
  isl_assert(tab->mat->ctx, tab->bmap, return -1);
2344
1.97k
  n_div = isl_basic_map_dim(tab->bmap, isl_dim_div);
2345
1.97k
  o_div = tab->n_var - n_div;
2346
1.97k
  if (pos < o_div || pos > tab->n_var)
2347
1.97k
    
isl_die0
(isl_tab_get_ctx(tab), isl_error_invalid,
2348
1.97k
      "invalid position", return -1);
2349
1.97k
2350
1.97k
  nonneg = div_is_nonneg(tab, div);
2351
1.97k
2352
1.97k
  if (isl_tab_extend_cons(tab, 3) < 0)
2353
0
    return -1;
2354
1.97k
  if (isl_tab_extend_vars(tab, 1) < 0)
2355
0
    return -1;
2356
1.97k
  r = isl_tab_insert_var(tab, pos);
2357
1.97k
  if (r < 0)
2358
0
    return -1;
2359
1.97k
2360
1.97k
  if (nonneg)
2361
281
    tab->var[r].is_nonneg = 1;
2362
1.97k
2363
1.97k
  tab->bmap = isl_basic_map_insert_div(tab->bmap, pos - o_div, div);
2364
1.97k
  if (!tab->bmap)
2365
0
    return -1;
2366
1.97k
  if (isl_tab_push_var(tab, isl_tab_undo_bmap_div, &tab->var[r]) < 0)
2367
0
    return -1;
2368
1.97k
2369
1.97k
  if (add_div_constraints(tab, pos - o_div, add_ineq, user) < 0)
2370
0
    return -1;
2371
1.97k
2372
1.97k
  return r;
2373
1.97k
}
2374
2375
/* Add an extra div, prescribed by "div", to the tableau and
2376
 * the associated bmap (which is assumed to be non-NULL).
2377
 */
2378
int isl_tab_add_div(struct isl_tab *tab, __isl_keep isl_vec *div)
2379
1.37k
{
2380
1.37k
  if (!tab)
2381
0
    return -1;
2382
1.37k
  return isl_tab_insert_div(tab, tab->n_var, div, NULL, NULL);
2383
1.37k
}
2384
2385
/* If "track" is set, then we want to keep track of all constraints in tab
2386
 * in its bmap field.  This field is initialized from a copy of "bmap",
2387
 * so we need to make sure that all constraints in "bmap" also appear
2388
 * in the constructed tab.
2389
 */
2390
__isl_give struct isl_tab *isl_tab_from_basic_map(
2391
  __isl_keep isl_basic_map *bmap, int track)
2392
722k
{
2393
722k
  int i;
2394
722k
  struct isl_tab *tab;
2395
722k
2396
722k
  if (!bmap)
2397
0
    return NULL;
2398
722k
  tab = isl_tab_alloc(bmap->ctx,
2399
722k
          isl_basic_map_total_dim(bmap) + bmap->n_ineq + 1,
2400
722k
          isl_basic_map_total_dim(bmap), 0);
2401
722k
  if (!tab)
2402
0
    return NULL;
2403
722k
  tab->preserve = track;
2404
722k
  tab->rational = ISL_F_ISSET(bmap, ISL_BASIC_MAP_RATIONAL);
2405
722k
  if (ISL_F_ISSET(bmap, ISL_BASIC_MAP_EMPTY)) {
2406
19
    if (isl_tab_mark_empty(tab) < 0)
2407
0
      goto error;
2408
19
    goto done;
2409
19
  }
2410
982k
  
for (i = 0; 722k
i < bmap->n_eq;
++i259k
) {
2411
259k
    tab = add_eq(tab, bmap->eq[i]);
2412
259k
    if (!tab)
2413
0
      return tab;
2414
259k
  }
2415
4.40M
  
for (i = 0; 722k
i < bmap->n_ineq;
++i3.67M
) {
2416
3.68M
    if (isl_tab_add_ineq(tab, bmap->ineq[i]) < 0)
2417
0
      goto error;
2418
3.68M
    if (tab->empty)
2419
5.03k
      goto done;
2420
3.68M
  }
2421
722k
done:
2422
722k
  if (track && 
isl_tab_track_bmap(tab, isl_basic_map_copy(bmap)) < 0189k
)
2423
0
    goto error;
2424
722k
  return tab;
2425
0
error:
2426
0
  isl_tab_free(tab);
2427
0
  return NULL;
2428
722k
}
2429
2430
__isl_give struct isl_tab *isl_tab_from_basic_set(
2431
  __isl_keep isl_basic_set *bset, int track)
2432
275k
{
2433
275k
  return isl_tab_from_basic_map(bset, track);
2434
275k
}
2435
2436
/* Construct a tableau corresponding to the recession cone of "bset".
2437
 */
2438
struct isl_tab *isl_tab_from_recession_cone(__isl_keep isl_basic_set *bset,
2439
  int parametric)
2440
3.82k
{
2441
3.82k
  isl_int cst;
2442
3.82k
  int i;
2443
3.82k
  struct isl_tab *tab;
2444
3.82k
  unsigned offset = 0;
2445
3.82k
2446
3.82k
  if (!bset)
2447
0
    return NULL;
2448
3.82k
  if (parametric)
2449
2.93k
    offset = isl_basic_set_dim(bset, isl_dim_param);
2450
3.82k
  tab = isl_tab_alloc(bset->ctx, bset->n_eq + bset->n_ineq,
2451
3.82k
        isl_basic_set_total_dim(bset) - offset, 0);
2452
3.82k
  if (!tab)
2453
0
    return NULL;
2454
3.82k
  tab->rational = ISL_F_ISSET(bset, ISL_BASIC_SET_RATIONAL);
2455
3.82k
  tab->cone = 1;
2456
3.82k
2457
3.82k
  isl_int_init(cst);
2458
3.82k
  isl_int_set_si(cst, 0);
2459
6.19k
  for (i = 0; i < bset->n_eq; 
++i2.37k
) {
2460
2.37k
    isl_int_swap(bset->eq[i][offset], cst);
2461
2.37k
    if (offset > 0) {
2462
697
      if (isl_tab_add_eq(tab, bset->eq[i] + offset) < 0)
2463
0
        goto error;
2464
1.67k
    } else
2465
1.67k
      tab = add_eq(tab, bset->eq[i]);
2466
2.37k
    isl_int_swap(bset->eq[i][offset], cst);
2467
2.37k
    if (!tab)
2468
0
      goto done;
2469
2.37k
  }
2470
15.8k
  
for (i = 0; 3.82k
i < bset->n_ineq;
++i12.0k
) {
2471
12.0k
    int r;
2472
12.0k
    isl_int_swap(bset->ineq[i][offset], cst);
2473
12.0k
    r = isl_tab_add_row(tab, bset->ineq[i] + offset);
2474
12.0k
    isl_int_swap(bset->ineq[i][offset], cst);
2475
12.0k
    if (r < 0)
2476
0
      goto error;
2477
12.0k
    tab->con[r].is_nonneg = 1;
2478
12.0k
    if (isl_tab_push_var(tab, isl_tab_undo_nonneg, &tab->con[r]) < 0)
2479
0
      goto error;
2480
12.0k
  }
2481
3.82k
done:
2482
3.82k
  isl_int_clear(cst);
2483
3.82k
  return tab;
2484
0
error:
2485
0
  isl_int_clear(cst);
2486
0
  isl_tab_free(tab);
2487
0
  return NULL;
2488
3.82k
}
2489
2490
/* Assuming "tab" is the tableau of a cone, check if the cone is
2491
 * bounded, i.e., if it is empty or only contains the origin.
2492
 */
2493
isl_bool isl_tab_cone_is_bounded(struct isl_tab *tab)
2494
2.93k
{
2495
2.93k
  int i;
2496
2.93k
2497
2.93k
  if (!tab)
2498
0
    return isl_bool_error;
2499
2.93k
  if (tab->empty)
2500
0
    return isl_bool_true;
2501
2.93k
  if (tab->n_dead == tab->n_col)
2502
789
    return isl_bool_true;
2503
2.14k
2504
3.12k
  
for (;;)2.14k
{
2505
3.40k
    for (i = tab->n_redundant; i < tab->n_row; 
++i279
) {
2506
3.40k
      struct isl_tab_var *var;
2507
3.40k
      int sgn;
2508
3.40k
      var = isl_tab_var_from_row(tab, i);
2509
3.40k
      if (!var->is_nonneg)
2510
279
        continue;
2511
3.12k
      sgn = sign_of_max(tab, var);
2512
3.12k
      if (sgn < -1)
2513
0
        return isl_bool_error;
2514
3.12k
      if (sgn != 0)
2515
311
        return isl_bool_false;
2516
2.81k
      if (close_row(tab, var, 0) < 0)
2517
0
        return isl_bool_error;
2518
2.81k
      break;
2519
2.81k
    }
2520
3.12k
    
if (2.81k
tab->n_dead == tab->n_col2.81k
)
2521
1.83k
      return isl_bool_true;
2522
979
    if (i == tab->n_row)
2523
3
      return isl_bool_false;
2524
979
  }
2525
2.14k
}
2526
2527
int isl_tab_sample_is_integer(struct isl_tab *tab)
2528
427k
{
2529
427k
  int i;
2530
427k
2531
427k
  if (!tab)
2532
0
    return -1;
2533
427k
2534
2.07M
  
for (i = 0; 427k
i < tab->n_var;
++i1.64M
) {
2535
1.74M
    int row;
2536
1.74M
    if (!tab->var[i].is_row)
2537
499k
      continue;
2538
1.24M
    row = tab->var[i].index;
2539
1.24M
    if (!isl_int_is_divisible_by(tab->mat->row[row][1],
2540
1.24M
            tab->mat->row[row][0]))
2541
1.24M
      
return 096.1k
;
2542
1.24M
  }
2543
427k
  
return 1331k
;
2544
427k
}
2545
2546
static struct isl_vec *extract_integer_sample(struct isl_tab *tab)
2547
201k
{
2548
201k
  int i;
2549
201k
  struct isl_vec *vec;
2550
201k
2551
201k
  vec = isl_vec_alloc(tab->mat->ctx, 1 + tab->n_var);
2552
201k
  if (!vec)
2553
0
    return NULL;
2554
201k
2555
201k
  isl_int_set_si(vec->block.data[0], 1);
2556
1.27M
  for (i = 0; i < tab->n_var; 
++i1.07M
) {
2557
1.07M
    if (!tab->var[i].is_row)
2558
1.07M
      
isl_int_set_si388k
(vec->block.data[1 + i], 0);
2559
1.07M
    else {
2560
689k
      int row = tab->var[i].index;
2561
689k
      isl_int_divexact(vec->block.data[1 + i],
2562
689k
        tab->mat->row[row][1], tab->mat->row[row][0]);
2563
689k
    }
2564
1.07M
  }
2565
201k
2566
201k
  return vec;
2567
201k
}
2568
2569
struct isl_vec *isl_tab_get_sample_value(struct isl_tab *tab)
2570
255k
{
2571
255k
  int i;
2572
255k
  struct isl_vec *vec;
2573
255k
  isl_int m;
2574
255k
2575
255k
  if (!tab)
2576
0
    return NULL;
2577
255k
2578
255k
  vec = isl_vec_alloc(tab->mat->ctx, 1 + tab->n_var);
2579
255k
  if (!vec)
2580
0
    return NULL;
2581
255k
2582
255k
  isl_int_init(m);
2583
255k
2584
255k
  isl_int_set_si(vec->block.data[0], 1);
2585
1.12M
  for (i = 0; i < tab->n_var; 
++i868k
) {
2586
868k
    int row;
2587
868k
    if (!tab->var[i].is_row) {
2588
412k
      isl_int_set_si(vec->block.data[1 + i], 0);
2589
412k
      continue;
2590
412k
    }
2591
455k
    row = tab->var[i].index;
2592
455k
    isl_int_gcd(m, vec->block.data[0], tab->mat->row[row][0]);
2593
455k
    isl_int_divexact(m, tab->mat->row[row][0], m);
2594
455k
    isl_seq_scale(vec->block.data, vec->block.data, m, 1 + i);
2595
455k
    isl_int_divexact(m, vec->block.data[0], tab->mat->row[row][0]);
2596
455k
    isl_int_mul(vec->block.data[1 + i], m, tab->mat->row[row][1]);
2597
455k
  }
2598
255k
  vec = isl_vec_normalize(vec);
2599
255k
2600
255k
  isl_int_clear(m);
2601
255k
  return vec;
2602
255k
}
2603
2604
/* Store the sample value of "var" of "tab" rounded up (if sgn > 0)
2605
 * or down (if sgn < 0) to the nearest integer in *v.
2606
 */
2607
static void get_rounded_sample_value(struct isl_tab *tab,
2608
  struct isl_tab_var *var, int sgn, isl_int *v)
2609
242k
{
2610
242k
  if (!var->is_row)
2611
242k
    
isl_int_set_si3.02k
(*v, 0);
2612
242k
  else 
if (239k
sgn > 0239k
)
2613
239k
    
isl_int_cdiv_q237k
(*v, tab->mat->row[var->index][1],
2614
239k
           tab->mat->row[var->index][0]);
2615
239k
  else
2616
239k
    
isl_int_fdiv_q1.87k
(*v, tab->mat->row[var->index][1],
2617
242k
           tab->mat->row[var->index][0]);
2618
242k
}
2619
2620
/* Update "bmap" based on the results of the tableau "tab".
2621
 * In particular, implicit equalities are made explicit, redundant constraints
2622
 * are removed and if the sample value happens to be integer, it is stored
2623
 * in "bmap" (unless "bmap" already had an integer sample).
2624
 *
2625
 * The tableau is assumed to have been created from "bmap" using
2626
 * isl_tab_from_basic_map.
2627
 */
2628
struct isl_basic_map *isl_basic_map_update_from_tab(struct isl_basic_map *bmap,
2629
  struct isl_tab *tab)
2630
429k
{
2631
429k
  int i;
2632
429k
  unsigned n_eq;
2633
429k
2634
429k
  if (!bmap)
2635
0
    return NULL;
2636
429k
  if (!tab)
2637
0
    return bmap;
2638
429k
2639
429k
  n_eq = tab->n_eq;
2640
429k
  if (tab->empty)
2641
6.59k
    bmap = isl_basic_map_set_to_empty(bmap);
2642
423k
  else
2643
2.94M
    
for (i = bmap->n_ineq - 1; 423k
i >= 0;
--i2.52M
) {
2644
2.52M
      if (isl_tab_is_equality(tab, n_eq + i))
2645
1.04M
        isl_basic_map_inequality_to_equality(bmap, i);
2646
1.47M
      else if (isl_tab_is_redundant(tab, n_eq + i))
2647
191k
        isl_basic_map_drop_inequality(bmap, i);
2648
2.52M
    }
2649
429k
  if (bmap->n_eq != n_eq)
2650
156k
    bmap = isl_basic_map_gauss(bmap, NULL);
2651
429k
  if (!tab->rational &&
2652
429k
      
bmap395k
&&
!bmap->sample395k
&&
isl_tab_sample_is_integer(tab)217k
)
2653
201k
    bmap->sample = extract_integer_sample(tab);
2654
429k
  return bmap;
2655
429k
}
2656
2657
struct isl_basic_set *isl_basic_set_update_from_tab(struct isl_basic_set *bset,
2658
  struct isl_tab *tab)
2659
38.9k
{
2660
38.9k
  return bset_from_bmap(isl_basic_map_update_from_tab(bset_to_bmap(bset),
2661
38.9k
                tab));
2662
38.9k
}
2663
2664
/* Drop the last constraint added to "tab" in position "r".
2665
 * The constraint is expected to have remained in a row.
2666
 */
2667
static isl_stat drop_last_con_in_row(struct isl_tab *tab, int r)
2668
13.0k
{
2669
13.0k
  if (!tab->con[r].is_row)
2670
13.0k
    
isl_die0
(isl_tab_get_ctx(tab), isl_error_internal,
2671
13.0k
      "row unexpectedly moved to column",
2672
13.0k
      return isl_stat_error);
2673
13.0k
  if (r + 1 != tab->n_con)
2674
13.0k
    
isl_die0
(isl_tab_get_ctx(tab), isl_error_internal,
2675
13.0k
      "additional constraints added", return isl_stat_error);
2676
13.0k
  if (drop_row(tab, tab->con[r].index) < 0)
2677
0
    return isl_stat_error;
2678
13.0k
2679
13.0k
  return isl_stat_ok;
2680
13.0k
}
2681
2682
/* Given a non-negative variable "var", temporarily add a new non-negative
2683
 * variable that is the opposite of "var", ensuring that "var" can only attain
2684
 * the value zero.  The new variable is removed again before this function
2685
 * returns.  However, the effect of forcing "var" to be zero remains.
2686
 * If var = n/d is a row variable, then the new variable = -n/d.
2687
 * If var is a column variables, then the new variable = -var.
2688
 * If the new variable cannot attain non-negative values, then
2689
 * the resulting tableau is empty.
2690
 * Otherwise, we know the value will be zero and we close the row.
2691
 */
2692
static isl_stat cut_to_hyperplane(struct isl_tab *tab, struct isl_tab_var *var)
2693
13.0k
{
2694
13.0k
  unsigned r;
2695
13.0k
  isl_int *row;
2696
13.0k
  int sgn;
2697
13.0k
  unsigned off = 2 + tab->M;
2698
13.0k
2699
13.0k
  if (var->is_zero)
2700
0
    return isl_stat_ok;
2701
13.0k
  if (var->is_redundant || !var->is_nonneg)
2702
13.0k
    
isl_die0
(isl_tab_get_ctx(tab), isl_error_invalid,
2703
13.0k
      "expecting non-redundant non-negative variable",
2704
13.0k
      return isl_stat_error);
2705
13.0k
2706
13.0k
  if (isl_tab_extend_cons(tab, 1) < 0)
2707
0
    return isl_stat_error;
2708
13.0k
2709
13.0k
  r = tab->n_con;
2710
13.0k
  tab->con[r].index = tab->n_row;
2711
13.0k
  tab->con[r].is_row = 1;
2712
13.0k
  tab->con[r].is_nonneg = 0;
2713
13.0k
  tab->con[r].is_zero = 0;
2714
13.0k
  tab->con[r].is_redundant = 0;
2715
13.0k
  tab->con[r].frozen = 0;
2716
13.0k
  tab->con[r].negated = 0;
2717
13.0k
  tab->row_var[tab->n_row] = ~r;
2718
13.0k
  row = tab->mat->row[tab->n_row];
2719
13.0k
2720
13.0k
  if (var->is_row) {
2721
1.22k
    isl_int_set(row[0], tab->mat->row[var->index][0]);
2722
1.22k
    isl_seq_neg(row + 1,
2723
1.22k
          tab->mat->row[var->index] + 1, 1 + tab->n_col);
2724
11.8k
  } else {
2725
11.8k
    isl_int_set_si(row[0], 1);
2726
11.8k
    isl_seq_clr(row + 1, 1 + tab->n_col);
2727
11.8k
    isl_int_set_si(row[off + var->index], -1);
2728
11.8k
  }
2729
13.0k
2730
13.0k
  tab->n_row++;
2731
13.0k
  tab->n_con++;
2732
13.0k
2733
13.0k
  sgn = sign_of_max(tab, &tab->con[r]);
2734
13.0k
  if (sgn < -1)
2735
0
    return isl_stat_error;
2736
13.0k
  if (sgn < 0) {
2737
51
    if (drop_last_con_in_row(tab, r) < 0)
2738
0
      return isl_stat_error;
2739
51
    if (isl_tab_mark_empty(tab) < 0)
2740
0
      return isl_stat_error;
2741
51
    return isl_stat_ok;
2742
51
  }
2743
13.0k
  tab->con[r].is_nonneg = 1;
2744
13.0k
  /* sgn == 0 */
2745
13.0k
  if (close_row(tab, &tab->con[r], 1) < 0)
2746
0
    return isl_stat_error;
2747
13.0k
  if (drop_last_con_in_row(tab, r) < 0)
2748
0
    return isl_stat_error;
2749
13.0k
2750
13.0k
  return isl_stat_ok;
2751
13.0k
}
2752
2753
/* Given a tableau "tab" and an inequality constraint "con" of the tableau,
2754
 * relax the inequality by one.  That is, the inequality r >= 0 is replaced
2755
 * by r' = r + 1 >= 0.
2756
 * If r is a row variable, we simply increase the constant term by one
2757
 * (taking into account the denominator).
2758
 * If r is a column variable, then we need to modify each row that
2759
 * refers to r = r' - 1 by substituting this equality, effectively
2760
 * subtracting the coefficient of the column from the constant.
2761
 * We should only do this if the minimum is manifestly unbounded,
2762
 * however.  Otherwise, we may end up with negative sample values
2763
 * for non-negative variables.
2764
 * So, if r is a column variable with a minimum that is not
2765
 * manifestly unbounded, then we need to move it to a row.
2766
 * However, the sample value of this row may be negative,
2767
 * even after the relaxation, so we need to restore it.
2768
 * We therefore prefer to pivot a column up to a row, if possible.
2769
 */
2770
int isl_tab_relax(struct isl_tab *tab, int con)
2771
1.67k
{
2772
1.67k
  struct isl_tab_var *var;
2773
1.67k
2774
1.67k
  if (!tab)
2775
0
    return -1;
2776
1.67k
2777
1.67k
  var = &tab->con[con];
2778
1.67k
2779
1.67k
  if (var->is_row && 
(53
var->index < 053
||
var->index < tab->n_redundant53
))
2780
1.67k
    
isl_die0
(tab->mat->ctx, isl_error_invalid,
2781
1.67k
      "cannot relax redundant constraint", return -1);
2782
1.67k
  if (!var->is_row && 
(1.62k
var->index < 01.62k
||
var->index < tab->n_dead1.62k
))
2783
1.67k
    
isl_die0
(tab->mat->ctx, isl_error_invalid,
2784
1.67k
      "cannot relax dead constraint", return -1);
2785
1.67k
2786
1.67k
  if (!var->is_row && 
!max_is_manifestly_unbounded(tab, var)1.62k
)
2787
405
    if (to_row(tab, var, 1) < 0)
2788
0
      return -1;
2789
1.67k
  if (!var->is_row && 
!min_is_manifestly_unbounded(tab, var)1.22k
)
2790
17
    if (to_row(tab, var, -1) < 0)
2791
0
      return -1;
2792
1.67k
2793
1.67k
  if (var->is_row) {
2794
475
    isl_int_add(tab->mat->row[var->index][1],
2795
475
        tab->mat->row[var->index][1], tab->mat->row[var->index][0]);
2796
475
    if (restore_row(tab, var) < 0)
2797
0
      return -1;
2798
1.20k
  } else {
2799
1.20k
    int i;
2800
1.20k
    unsigned off = 2 + tab->M;
2801
1.20k
2802
6.66k
    for (i = 0; i < tab->n_row; 
++i5.45k
) {
2803
5.45k
      if (isl_int_is_zero(tab->mat->row[i][off + var->index]))
2804
5.45k
        
continue4.13k
;
2805
1.32k
      isl_int_sub(tab->mat->row[i][1], tab->mat->row[i][1],
2806
1.32k
          tab->mat->row[i][off + var->index]);
2807
1.32k
    }
2808
1.20k
2809
1.20k
  }
2810
1.67k
2811
1.67k
  if (isl_tab_push_var(tab, isl_tab_undo_relax, var) < 0)
2812
0
    return -1;
2813
1.67k
2814
1.67k
  return 0;
2815
1.67k
}
2816
2817
/* Replace the variable v at position "pos" in the tableau "tab"
2818
 * by v' = v + shift.
2819
 *
2820
 * If the variable is in a column, then we first check if we can
2821
 * simply plug in v = v' - shift.  The effect on a row with
2822
 * coefficient f/d for variable v is that the constant term c/d
2823
 * is replaced by (c - f * shift)/d.  If shift is positive and
2824
 * f is negative for each row that needs to remain non-negative,
2825
 * then this is clearly safe.  In other words, if the minimum of v
2826
 * is manifestly unbounded, then we can keep v in a column position.
2827
 * Otherwise, we can pivot it down to a row.
2828
 * Similarly, if shift is negative, we need to check if the maximum
2829
 * of is manifestly unbounded.
2830
 *
2831
 * If the variable is in a row (from the start or after pivoting),
2832
 * then the constant term c/d is replaced by (c + d * shift)/d.
2833
 */
2834
int isl_tab_shift_var(struct isl_tab *tab, int pos, isl_int shift)
2835
142
{
2836
142
  struct isl_tab_var *var;
2837
142
2838
142
  if (!tab)
2839
0
    return -1;
2840
142
  if (isl_int_is_zero(shift))
2841
142
    
return 073
;
2842
69
2843
69
  var = &tab->var[pos];
2844
69
  if (!var->is_row) {
2845
9
    if (isl_int_is_neg(shift)) {
2846
6
      if (!max_is_manifestly_unbounded(tab, var))
2847
4
        if (to_row(tab, var, 1) < 0)
2848
0
          return -1;
2849
3
    } else {
2850
3
      if (!min_is_manifestly_unbounded(tab, var))
2851
0
        if (to_row(tab, var, -1) < 0)
2852
0
          return -1;
2853
69
    }
2854
9
  }
2855
69
2856
69
  if (var->is_row) {
2857
64
    isl_int_addmul(tab->mat->row[var->index][1],
2858
64
        shift, tab->mat->row[var->index][0]);
2859
64
  } else {
2860
5
    int i;
2861
5
    unsigned off = 2 + tab->M;
2862
5
2863
26
    for (i = 0; i < tab->n_row; 
++i21
) {
2864
21
      if (isl_int_is_zero(tab->mat->row[i][off + var->index]))
2865
21
        
continue15
;
2866
6
      isl_int_submul(tab->mat->row[i][1],
2867
6
            shift, tab->mat->row[i][off + var->index]);
2868
6
    }
2869
5
2870
5
  }
2871
69
2872
69
  return 0;
2873
69
}
2874
2875
/* Remove the sign constraint from constraint "con".
2876
 *
2877
 * If the constraint variable was originally marked non-negative,
2878
 * then we make sure we mark it non-negative again during rollback.
2879
 */
2880
int isl_tab_unrestrict(struct isl_tab *tab, int con)
2881
436
{
2882
436
  struct isl_tab_var *var;
2883
436
2884
436
  if (!tab)
2885
0
    return -1;
2886
436
2887
436
  var = &tab->con[con];
2888
436
  if (!var->is_nonneg)
2889
0
    return 0;
2890
436
2891
436
  var->is_nonneg = 0;
2892
436
  if (isl_tab_push_var(tab, isl_tab_undo_unrestrict, var) < 0)
2893
0
    return -1;
2894
436
2895
436
  return 0;
2896
436
}
2897
2898
int isl_tab_select_facet(struct isl_tab *tab, int con)
2899
12.5k
{
2900
12.5k
  if (!tab)
2901
0
    return -1;
2902
12.5k
2903
12.5k
  return cut_to_hyperplane(tab, &tab->con[con]);
2904
12.5k
}
2905
2906
static int may_be_equality(struct isl_tab *tab, int row)
2907
8.18M
{
2908
8.18M
  return tab->rational ? 
isl_int_is_zero36.3k
(tab->mat->row[row][1])
2909
8.18M
           : 
isl_int_lt8.15M
(tab->mat->row[row][1],
2910
8.18M
              tab->mat->row[row][0]);
2911
8.18M
}
2912
2913
/* Return an isl_tab_var that has been marked or NULL if no such
2914
 * variable can be found.
2915
 * The marked field has only been set for variables that
2916
 * appear in non-redundant rows or non-dead columns.
2917
 *
2918
 * Pick the last constraint variable that is marked and
2919
 * that appears in either a non-redundant row or a non-dead columns.
2920
 * Since the returned variable is tested for being a redundant constraint or
2921
 * an implicit equality, there is no need to return any tab variable that
2922
 * corresponds to a variable.
2923
 */
2924
static struct isl_tab_var *select_marked(struct isl_tab *tab)
2925
2.01M
{
2926
2.01M
  int i;
2927
2.01M
  struct isl_tab_var *var;
2928
2.01M
2929
16.8M
  for (i = tab->n_con - 1; i >= 0; 
--i14.8M
) {
2930
16.7M
    var = &tab->con[i];
2931
16.7M
    if (var->index < 0)
2932
6.49M
      continue;
2933
10.2M
    if (var->is_row && 
var->index < tab->n_redundant8.05M
)
2934
606k
      continue;
2935
9.60M
    if (!var->is_row && 
var->index < tab->n_dead2.15M
)
2936
1.71k
      continue;
2937
9.60M
    if (var->marked)
2938
1.85M
      return var;
2939
9.60M
  }
2940
2.01M
2941
2.01M
  
return NULL162k
;
2942
2.01M
}
2943
2944
/* Check for (near) equalities among the constraints.
2945
 * A constraint is an equality if it is non-negative and if
2946
 * its maximal value is either
2947
 *  - zero (in case of rational tableaus), or
2948
 *  - strictly less than 1 (in case of integer tableaus)
2949
 *
2950
 * We first mark all non-redundant and non-dead variables that
2951
 * are not frozen and not obviously not an equality.
2952
 * Then we iterate over all marked variables if they can attain
2953
 * any values larger than zero or at least one.
2954
 * If the maximal value is zero, we mark any column variables
2955
 * that appear in the row as being zero and mark the row as being redundant.
2956
 * Otherwise, if the maximal value is strictly less than one (and the
2957
 * tableau is integer), then we restrict the value to being zero
2958
 * by adding an opposite non-negative variable.
2959
 * The order in which the variables are considered is not important.
2960
 */
2961
int isl_tab_detect_implicit_equalities(struct isl_tab *tab)
2962
425k
{
2963
425k
  int i;
2964
425k
  unsigned n_marked;
2965
425k
2966
425k
  if (!tab)
2967
0
    return -1;
2968
425k
  if (tab->empty)
2969
3.69k
    return 0;
2970
421k
  if (tab->n_dead == tab->n_col)
2971
18.1k
    return 0;
2972
403k
2973
403k
  n_marked = 0;
2974
3.00M
  for (i = tab->n_redundant; i < tab->n_row; 
++i2.60M
) {
2975
2.60M
    struct isl_tab_var *var = isl_tab_var_from_row(tab, i);
2976
2.60M
    var->marked = !var->frozen && 
var->is_nonneg2.58M
&&
2977
2.60M
      
may_be_equality(tab, i)2.36M
;
2978
2.60M
    if (var->marked)
2979
1.74M
      n_marked++;
2980
2.60M
  }
2981
2.29M
  for (i = tab->n_dead; i < tab->n_col; 
++i1.89M
) {
2982
1.89M
    struct isl_tab_var *var = var_from_col(tab, i);
2983
1.89M
    var->marked = !var->frozen && 
var->is_nonneg1.88M
;
2984
1.89M
    if (var->marked)
2985
171k
      n_marked++;
2986
1.89M
  }
2987
1.67M
  while (n_marked) {
2988
1.42M
    struct isl_tab_var *var;
2989
1.42M
    int sgn;
2990
1.42M
    var = select_marked(tab);
2991
1.42M
    if (!var)
2992
156k
      break;
2993
1.26M
    var->marked = 0;
2994
1.26M
    n_marked--;
2995
1.26M
    sgn = sign_of_max(tab, var);
2996
1.26M
    if (sgn < 0)
2997
0
      return -1;
2998
1.26M
    if (sgn == 0) {
2999
558k
      if (close_row(tab, var, 0) < 0)
3000
0
        return -1;
3001
708k
    } else if (!tab->rational && 
!at_least_one(tab, var)685k
) {
3002
537
      if (cut_to_hyperplane(tab, var) < 0)
3003
0
        return -1;
3004
537
      return isl_tab_detect_implicit_equalities(tab);
3005
537
    }
3006
11.4M
    
for (i = tab->n_redundant; 1.26M
i < tab->n_row;
++i10.2M
) {
3007
10.2M
      var = isl_tab_var_from_row(tab, i);
3008
10.2M
      if (!var->marked)
3009
4.40M
        continue;
3010
5.81M
      if (may_be_equality(tab, i))
3011
5.73M
        continue;
3012
79.8k
      var->marked = 0;
3013
79.8k
      n_marked--;
3014
79.8k
    }
3015
1.26M
  }
3016
403k
3017
403k
  
return 0403k
;
3018
403k
}
3019
3020
/* Update the element of row_var or col_var that corresponds to
3021
 * constraint tab->con[i] to a move from position "old" to position "i".
3022
 */
3023
static int update_con_after_move(struct isl_tab *tab, int i, int old)
3024
6.35k
{
3025
6.35k
  int *p;
3026
6.35k
  int index;
3027
6.35k
3028
6.35k
  index = tab->con[i].index;
3029
6.35k
  if (index == -1)
3030
4.20k
    return 0;
3031
2.15k
  p = tab->con[i].is_row ? 
tab->row_var1.47k
:
tab->col_var672
;
3032
2.15k
  if (p[index] != ~old)
3033
2.15k
    
isl_die0
(tab->mat->ctx, isl_error_internal,
3034
2.15k
      "broken internal state", return -1);
3035
2.15k
  p[index] = ~i;
3036
2.15k
3037
2.15k
  return 0;
3038
2.15k
}
3039
3040
/* Rotate the "n" constraints starting at "first" to the right,
3041
 * putting the last constraint in the position of the first constraint.
3042
 */
3043
static int rotate_constraints(struct isl_tab *tab, int first, int n)
3044
1.83k
{
3045
1.83k
  int i, last;
3046
1.83k
  struct isl_tab_var var;
3047
1.83k
3048
1.83k
  if (n <= 1)
3049
596
    return 0;
3050
1.24k
3051
1.24k
  last = first + n - 1;
3052
1.24k
  var = tab->con[last];
3053
6.35k
  for (i = last; i > first; 
--i5.11k
) {
3054
5.11k
    tab->con[i] = tab->con[i - 1];
3055
5.11k
    if (update_con_after_move(tab, i, i - 1) < 0)
3056
0
      return -1;
3057
5.11k
  }
3058
1.24k
  tab->con[first] = var;
3059
1.24k
  if (update_con_after_move(tab, first, last) < 0)
3060
0
    return -1;
3061
1.24k
3062
1.24k
  return 0;
3063
1.24k
}
3064
3065
/* Make the equalities that are implicit in "bmap" but that have been
3066
 * detected in the corresponding "tab" explicit in "bmap" and update
3067
 * "tab" to reflect the new order of the constraints.
3068
 *
3069
 * In particular, if inequality i is an implicit equality then
3070
 * isl_basic_map_inequality_to_equality will move the inequality
3071
 * in front of the other equality and it will move the last inequality
3072
 * in the position of inequality i.
3073
 * In the tableau, the inequalities of "bmap" are stored after the equalities
3074
 * and so the original order
3075
 *
3076
 *    E E E E E A A A I B B B B L
3077
 *
3078
 * is changed into
3079
 *
3080
 *    I E E E E E A A A L B B B B
3081
 *
3082
 * where I is the implicit equality, the E are equalities,
3083
 * the A inequalities before I, the B inequalities after I and
3084
 * L the last inequality.
3085
 * We therefore need to rotate to the right two sets of constraints,
3086
 * those up to and including I and those after I.
3087
 *
3088
 * If "tab" contains any constraints that are not in "bmap" then they
3089
 * appear after those in "bmap" and they should be left untouched.
3090
 *
3091
 * Note that this function leaves "bmap" in a temporary state
3092
 * as it does not call isl_basic_map_gauss.  Calling this function
3093
 * is the responsibility of the caller.
3094
 */
3095
__isl_give isl_basic_map *isl_tab_make_equalities_explicit(struct isl_tab *tab,
3096
  __isl_take isl_basic_map *bmap)
3097
57.8k
{
3098
57.8k
  int i;
3099
57.8k
3100
57.8k
  if (!tab || !bmap)
3101
0
    return isl_basic_map_free(bmap);
3102
57.8k
  if (tab->empty)
3103
75
    return bmap;
3104
57.7k
3105
169k
  
for (i = bmap->n_ineq - 1; 57.7k
i >= 0;
--i111k
) {
3106
111k
    if (!isl_tab_is_equality(tab, bmap->n_eq + i))
3107
110k
      continue;
3108
918
    isl_basic_map_inequality_to_equality(bmap, i);
3109
918
    if (rotate_constraints(tab, 0, tab->n_eq + i + 1) < 0)
3110
0
      return isl_basic_map_free(bmap);
3111
918
    if (rotate_constraints(tab, tab->n_eq + i + 1,
3112
918
          bmap->n_ineq - i) < 0)
3113
0
      return isl_basic_map_free(bmap);
3114
918
    tab->n_eq++;
3115
918
  }
3116
57.7k
3117
57.7k
  return bmap;
3118
57.7k
}
3119
3120
static int con_is_redundant(struct isl_tab *tab, struct isl_tab_var *var)
3121
844k
{
3122
844k
  if (!tab)
3123
0
    return -1;
3124
844k
  if (tab->rational) {
3125
70.5k
    int sgn = sign_of_min(tab, var);
3126
70.5k
    if (sgn < -1)
3127
0
      return -1;
3128
70.5k
    return sgn >= 0;
3129
773k
  } else {
3130
773k
    int irred = isl_tab_min_at_most_neg_one(tab, var);
3131
773k
    if (irred < 0)
3132
0
      return -1;
3133
773k
    return !irred;
3134
773k
  }
3135
844k
}
3136
3137
/* Check for (near) redundant constraints.
3138
 * A constraint is redundant if it is non-negative and if
3139
 * its minimal value (temporarily ignoring the non-negativity) is either
3140
 *  - zero (in case of rational tableaus), or
3141
 *  - strictly larger than -1 (in case of integer tableaus)
3142
 *
3143
 * We first mark all non-redundant and non-dead variables that
3144
 * are not frozen and not obviously negatively unbounded.
3145
 * Then we iterate over all marked variables if they can attain
3146
 * any values smaller than zero or at most negative one.
3147
 * If not, we mark the row as being redundant (assuming it hasn't
3148
 * been detected as being obviously redundant in the mean time).
3149
 */
3150
int isl_tab_detect_redundant(struct isl_tab *tab)
3151
202k
{
3152
202k
  int i;
3153
202k
  unsigned n_marked;
3154
202k
3155
202k
  if (!tab)
3156
0
    return -1;
3157
202k
  if (tab->empty)
3158
2.69k
    return 0;
3159
199k
  if (tab->n_redundant == tab->n_row)
3160
3.52k
    return 0;
3161
195k
3162
195k
  n_marked = 0;
3163
1.43M
  for (i = tab->n_redundant; i < tab->n_row; 
++i1.24M
) {
3164
1.24M
    struct isl_tab_var *var = isl_tab_var_from_row(tab, i);
3165
1.24M
    var->marked = !var->frozen && 
var->is_nonneg1.09M
;
3166
1.24M
    if (var->marked)
3167
563k
      n_marked++;
3168
1.24M
  }
3169
1.07M
  for (i = tab->n_dead; i < tab->n_col; 
++i877k
) {
3170
877k
    struct isl_tab_var *var = var_from_col(tab, i);
3171
877k
    var->marked = !var->frozen && 
var->is_nonneg820k
&&
3172
877k
      
!min_is_manifestly_unbounded(tab, var)369k
;
3173
877k
    if (var->marked)
3174
117k
      n_marked++;
3175
877k
  }
3176
778k
  while (n_marked) {
3177
589k
    struct isl_tab_var *var;
3178
589k
    int red;
3179
589k
    var = select_marked(tab);
3180
589k
    if (!var)
3181
6.96k
      break;
3182
582k
    var->marked = 0;
3183
582k
    n_marked--;
3184
582k
    red = con_is_redundant(tab, var);
3185
582k
    if (red < 0)
3186
0
      return -1;
3187
582k
    if (red && 
!var->is_redundant130k
)
3188
23.3k
      if (isl_tab_mark_redundant(tab, var->index) < 0)
3189
0
        return -1;
3190
8.29M
    
for (i = tab->n_dead; 582k
i < tab->n_col;
++i7.71M
) {
3191
7.71M
      var = var_from_col(tab, i);
3192
7.71M
      if (!var->marked)
3193
7.36M
        continue;
3194
355k
      if (!min_is_manifestly_unbounded(tab, var))
3195
267k
        continue;
3196
88.2k
      var->marked = 0;
3197
88.2k
      n_marked--;
3198
88.2k
    }
3199
582k
  }
3200
195k
3201
195k
  return 0;
3202
195k
}
3203
3204
int isl_tab_is_equality(struct isl_tab *tab, int con)
3205
2.65M
{
3206
2.65M
  int row;
3207
2.65M
  unsigned off;
3208
2.65M
3209
2.65M
  if (!tab)
3210
0
    return -1;
3211
2.65M
  if (tab->con[con].is_zero)
3212
1.04M
    return 1;
3213
1.60M
  if (tab->con[con].is_redundant)
3214
193k
    return 0;
3215
1.41M
  if (!tab->con[con].is_row)
3216
784k
    return tab->con[con].index < tab->n_dead;
3217
625k
3218
625k
  row = tab->con[con].index;
3219
625k
3220
625k
  off = 2 + tab->M;
3221
625k
  return isl_int_is_zero(tab->mat->row[row][1]) &&
3222
625k
    
!row_is_big(tab, row)55.2k
&&
3223
625k
    isl_seq_first_non_zero(tab->mat->row[row] + off + tab->n_dead,
3224
55.2k
          tab->n_col - tab->n_dead) == -1;
3225
625k
}
3226
3227
/* Return the minimal value of the affine expression "f" with denominator
3228
 * "denom" in *opt, *opt_denom, assuming the tableau is not empty and
3229
 * the expression cannot attain arbitrarily small values.
3230
 * If opt_denom is NULL, then *opt is rounded up to the nearest integer.
3231
 * The return value reflects the nature of the result (empty, unbounded,
3232
 * minimal value returned in *opt).
3233
 *
3234
 * This function assumes that at least one more row and at least
3235
 * one more element in the constraint array are available in the tableau.
3236
 */
3237
enum isl_lp_result isl_tab_min(struct isl_tab *tab,
3238
  isl_int *f, isl_int denom, isl_int *opt, isl_int *opt_denom,
3239
  unsigned flags)
3240
306k
{
3241
306k
  int r;
3242
306k
  enum isl_lp_result res = isl_lp_ok;
3243
306k
  struct isl_tab_var *var;
3244
306k
  struct isl_tab_undo *snap;
3245
306k
3246
306k
  if (!tab)
3247
0
    return isl_lp_error;
3248
306k
3249
306k
  if (tab->empty)
3250
27
    return isl_lp_empty;
3251
306k
3252
306k
  snap = isl_tab_snap(tab);
3253
306k
  r = isl_tab_add_row(tab, f);
3254
306k
  if (r < 0)
3255
0
    return isl_lp_error;
3256
306k
  var = &tab->con[r];
3257
670k
  for (;;) {
3258
670k
    int row, col;
3259
670k
    find_pivot(tab, var, var, -1, &row, &col);
3260
670k
    if (row == var->index) {
3261
9.78k
      res = isl_lp_unbounded;
3262
9.78k
      break;
3263
9.78k
    }
3264
660k
    if (row == -1)
3265
296k
      break;
3266
364k
    if (isl_tab_pivot(tab, row, col) < 0)
3267
0
      return isl_lp_error;
3268
364k
  }
3269
306k
  isl_int_mul(tab->mat->row[var->index][0],
3270
306k
        tab->mat->row[var->index][0], denom);
3271
306k
  if (ISL_FL_ISSET(flags, ISL_TAB_SAVE_DUAL)) {
3272
26.7k
    int i;
3273
26.7k
3274
26.7k
    isl_vec_free(tab->dual);
3275
26.7k
    tab->dual = isl_vec_alloc(tab->mat->ctx, 1 + tab->n_con);
3276
26.7k
    if (!tab->dual)
3277
0
      return isl_lp_error;
3278
26.7k
    isl_int_set(tab->dual->el[0], tab->mat->row[var->index][0]);
3279
1.05M
    for (i = 0; i < tab->n_con; 
++i1.02M
) {
3280
1.02M
      int pos;
3281
1.02M
      if (tab->con[i].is_row) {
3282
632k
        isl_int_set_si(tab->dual->el[1 + i], 0);
3283
632k
        continue;
3284
632k
      }
3285
392k
      pos = 2 + tab->M + tab->con[i].index;
3286
392k
      if (tab->con[i].negated)
3287
392k
        
isl_int_neg70.0k
(tab->dual->el[1 + i],
3288
392k
              tab->mat->row[var->index][pos]);
3289
392k
      else
3290
392k
        
isl_int_set322k
(tab->dual->el[1 + i],
3291
392k
              tab->mat->row[var->index][pos]);
3292
392k
    }
3293
26.7k
  }
3294
306k
  if (opt && 
res == isl_lp_ok305k
) {
3295
296k
    if (opt_denom) {
3296
59.8k
      isl_int_set(*opt, tab->mat->row[var->index][1]);
3297
59.8k
      isl_int_set(*opt_denom, tab->mat->row[var->index][0]);
3298
59.8k
    } else
3299
236k
      get_rounded_sample_value(tab, var, 1, opt);
3300
296k
  }
3301
306k
  if (isl_tab_rollback(tab, snap) < 0)
3302
0
    return isl_lp_error;
3303
306k
  return res;
3304
306k
}
3305
3306
/* Is the constraint at position "con" marked as being redundant?
3307
 * If it is marked as representing an equality, then it is not
3308
 * considered to be redundant.
3309
 * Note that isl_tab_mark_redundant marks both the isl_tab_var as
3310
 * redundant and moves the corresponding row into the first
3311
 * tab->n_redundant positions (or removes the row, assigning it index -1),
3312
 * so the final test is actually redundant itself.
3313
 */
3314
int isl_tab_is_redundant(struct isl_tab *tab, int con)
3315
1.89M
{
3316
1.89M
  if (!tab)
3317
0
    return -1;
3318
1.89M
  if (con < 0 || con >= tab->n_con)
3319
1.89M
    
isl_die0
(isl_tab_get_ctx(tab), isl_error_invalid,
3320
1.89M
      "position out of bounds", return -1);
3321
1.89M
  if (tab->con[con].is_zero)
3322
132
    return 0;
3323
1.89M
  if (tab->con[con].is_redundant)
3324
307k
    return 1;
3325
1.58M
  return tab->con[con].is_row && 
tab->con[con].index < tab->n_redundant658k
;
3326
1.58M
}
3327
3328
/* Is variable "var" of "tab" fixed to a constant value by its row
3329
 * in the tableau?
3330
 * If so and if "value" is not NULL, then store this constant value
3331
 * in "value".
3332
 *
3333
 * That is, is it a row variable that only has non-zero coefficients
3334
 * for dead columns?
3335
 */
3336
static isl_bool is_constant(struct isl_tab *tab, struct isl_tab_var *var,
3337
  isl_int *value)
3338
4.91k
{
3339
4.91k
  unsigned off = 2 + tab->M;
3340
4.91k
  isl_mat *mat = tab->mat;
3341
4.91k
  int n;
3342
4.91k
  int row;
3343
4.91k
  int pos;
3344
4.91k
3345
4.91k
  if (!var->is_row)
3346
3.02k
    return isl_bool_false;
3347
1.89k
  row = var->index;
3348
1.89k
  if (row_is_big(tab, row))
3349
0
    return isl_bool_false;
3350
1.89k
  n = tab->n_col - tab->n_dead;
3351
1.89k
  pos = isl_seq_first_non_zero(mat->row[row] + off + tab->n_dead, n);
3352
1.89k
  if (pos != -1)
3353
1.71k
    return isl_bool_false;
3354
172
  if (value)
3355
172
    
isl_int_divexact0
(*value, mat->row[row][1], mat->row[row][0]);
3356
172
  return isl_bool_true;
3357
172
}
3358
3359
/* Has the variable "var' of "tab" reached a value that is greater than
3360
 * or equal (if sgn > 0) or smaller than or equal (if sgn < 0) to "target"?
3361
 * "tmp" has been initialized by the caller and can be used
3362
 * to perform local computations.
3363
 *
3364
 * If the sample value involves the big parameter, then any value
3365
 * is reached.
3366
 * Otherwise check if n/d >= t, i.e., n >= d * t (if sgn > 0)
3367
 * or n/d <= t, i.e., n <= d * t (if sgn < 0).
3368
 */
3369
static int reached(struct isl_tab *tab, struct isl_tab_var *var, int sgn,
3370
  isl_int target, isl_int *tmp)
3371
5.86k
{
3372
5.86k
  if (row_is_big(tab, var->index))
3373
0
    return 1;
3374
5.86k
  isl_int_mul(*tmp, tab->mat->row[var->index][0], target);
3375
5.86k
  if (sgn > 0)
3376
2.30k
    return isl_int_ge(tab->mat->row[var->index][1], *tmp);
3377
5.86k
  else
3378
5.86k
    
return 3.56k
isl_int_le3.56k
(tab->mat->row[var->index][1], *tmp);
3379
5.86k
}
3380
3381
/* Can variable "var" of "tab" attain the value "target" by
3382
 * pivoting up (if sgn > 0) or down (if sgn < 0)?
3383
 * If not, then pivot up [down] to the greatest [smallest]
3384
 * rational value.
3385
 * "tmp" has been initialized by the caller and can be used
3386
 * to perform local computations.
3387
 *
3388
 * If the variable is manifestly unbounded in the desired direction,
3389
 * then it can attain any value.
3390
 * Otherwise, it can be moved to a row.
3391
 * Continue pivoting until the target is reached.
3392
 * If no more pivoting can be performed, the maximal [minimal]
3393
 * rational value has been reached and the target cannot be reached.
3394
 * If the variable would be pivoted into a manifestly unbounded column,
3395
 * then the target can be reached.
3396
 */
3397
static isl_bool var_reaches(struct isl_tab *tab, struct isl_tab_var *var,
3398
  int sgn, isl_int target, isl_int *tmp)
3399
5.97k
{
3400
5.97k
  int row, col;
3401
5.97k
3402
5.97k
  if (sgn < 0 && 
min_is_manifestly_unbounded(tab, var)4.74k
)
3403
2.16k
    return isl_bool_true;
3404
3.80k
  if (sgn > 0 && 
max_is_manifestly_unbounded(tab, var)1.23k
)
3405
0
    return isl_bool_true;
3406
3.80k
  if (to_row(tab, var, sgn) < 0)
3407
0
    return isl_bool_error;
3408
5.86k
  
while (3.80k
!reached(tab, var, sgn, target, tmp)) {
3409
5.26k
    find_pivot(tab, var, var, sgn, &row, &col);
3410
5.26k
    if (row == -1)
3411
1.38k
      return isl_bool_false;
3412
3.88k
    if (row == var->index)
3413
1.82k
      return isl_bool_true;
3414
2.06k
    if (isl_tab_pivot(tab, row, col) < 0)
3415
0
      return isl_bool_error;
3416
2.06k
  }
3417
3.80k
3418
3.80k
  
return isl_bool_true597
;
3419
3.80k
}
3420
3421
/* Check if variable "var" of "tab" can only attain a single (integer)
3422
 * value, and, if so, add an equality constraint to fix the variable
3423
 * to this single value and store the result in "target".
3424
 * "target" and "tmp" have been initialized by the caller.
3425
 *
3426
 * Given the current sample value, round it down and check
3427
 * whether it is possible to attain a strictly smaller integer value.
3428
 * If so, the variable is not restricted to a single integer value.
3429
 * Otherwise, the search stops at the smallest rational value.
3430
 * Round up this value and check whether it is possible to attain
3431
 * a strictly greater integer value.
3432
 * If so, the variable is not restricted to a single integer value.
3433
 * Otherwise, the search stops at the greatest rational value.
3434
 * If rounding down this value yields a value that is different
3435
 * from rounding up the smallest rational value, then the variable
3436
 * cannot attain any integer value.  Mark the tableau empty.
3437
 * Otherwise, add an equality constraint that fixes the variable
3438
 * to the single integer value found.
3439
 */
3440
static isl_bool detect_constant_with_tmp(struct isl_tab *tab,
3441
  struct isl_tab_var *var, isl_int *target, isl_int *tmp)
3442
4.74k
{
3443
4.74k
  isl_bool reached;
3444
4.74k
  isl_vec *eq;
3445
4.74k
  int pos;
3446
4.74k
  isl_stat r;
3447
4.74k
3448
4.74k
  get_rounded_sample_value(tab, var, -1, target);
3449
4.74k
  isl_int_sub_ui(*target, *target, 1);
3450
4.74k
  reached = var_reaches(tab, var, -1, *target, tmp);
3451
4.74k
  if (reached < 0 || reached)
3452
3.51k
    return isl_bool_not(reached);
3453
1.23k
  get_rounded_sample_value(tab, var, 1, target);
3454
1.23k
  isl_int_add_ui(*target, *target, 1);
3455
1.23k
  reached = var_reaches(tab, var, 1, *target, tmp);
3456
1.23k
  if (reached < 0 || reached)
3457
1.07k
    return isl_bool_not(reached);
3458
155
  get_rounded_sample_value(tab, var, -1, tmp);
3459
155
  isl_int_sub_ui(*target, *target, 1);
3460
155
  if (isl_int_ne(*target, *tmp)) {
3461
0
    if (isl_tab_mark_empty(tab) < 0)
3462
0
      return isl_bool_error;
3463
0
    return isl_bool_false;
3464
0
  }
3465
155
3466
155
  if (isl_tab_extend_cons(tab, 1) < 0)
3467
0
    return isl_bool_error;
3468
155
  eq = isl_vec_alloc(isl_tab_get_ctx(tab), 1 + tab->n_var);
3469
155
  if (!eq)
3470
0
    return isl_bool_error;
3471
155
  pos = var - tab->var;
3472
155
  isl_seq_clr(eq->el + 1, tab->n_var);
3473
155
  isl_int_set_si(eq->el[1 + pos], -1);
3474
155
  isl_int_set(eq->el[0], *target);
3475
155
  r = isl_tab_add_eq(tab, eq->el);
3476
155
  isl_vec_free(eq);
3477
155
3478
155
  return r < 0 ? 
isl_bool_error0
: isl_bool_true;
3479
155
}
3480
3481
/* Check if variable "var" of "tab" can only attain a single (integer)
3482
 * value, and, if so, add an equality constraint to fix the variable
3483
 * to this single value and store the result in "value" (if "value"
3484
 * is not NULL).
3485
 *
3486
 * If the current sample value involves the big parameter,
3487
 * then the variable cannot have a fixed integer value.
3488
 * If the variable is already fixed to a single value by its row, then
3489
 * there is no need to add another equality constraint.
3490
 *
3491
 * Otherwise, allocate some temporary variables and continue
3492
 * with detect_constant_with_tmp.
3493
 */
3494
static isl_bool get_constant(struct isl_tab *tab, struct isl_tab_var *var,
3495
  isl_int *value)
3496
4.91k
{
3497
4.91k
  isl_int target, tmp;
3498
4.91k
  isl_bool is_cst;
3499
4.91k
3500
4.91k
  if (var->is_row && 
row_is_big(tab, var->index)1.89k
)
3501
0
    return isl_bool_false;
3502
4.91k
  is_cst = is_constant(tab, var, value);
3503
4.91k
  if (is_cst < 0 || is_cst)
3504
172
    return is_cst;
3505
4.74k
3506
4.74k
  if (!value)
3507
4.74k
    
isl_int_init1.85k
(target);
3508
4.74k
  isl_int_init(tmp);
3509
4.74k
3510
4.74k
  is_cst = detect_constant_with_tmp(tab, var,
3511
4.74k
              value ? 
value2.88k
:
&target1.85k
, &tmp);
3512
4.74k
3513
4.74k
  isl_int_clear(tmp);
3514
4.74k
  if (!value)
3515
4.74k
    
isl_int_clear1.85k
(target);
3516
4.74k
3517
4.74k
  return is_cst;
3518
4.74k
}
3519
3520
/* Check if variable "var" of "tab" can only attain a single (integer)
3521
 * value, and, if so, add an equality constraint to fix the variable
3522
 * to this single value and store the result in "value" (if "value"
3523
 * is not NULL).
3524
 *
3525
 * For rational tableaus, nothing needs to be done.
3526
 */
3527
isl_bool isl_tab_is_constant(struct isl_tab *tab, int var, isl_int *value)
3528
2.88k
{
3529
2.88k
  if (!tab)
3530
0
    return isl_bool_error;
3531
2.88k
  if (var < 0 || var >= tab->n_var)
3532
2.88k
    
isl_die0
(isl_tab_get_ctx(tab), isl_error_invalid,
3533
2.88k
      "position out of bounds", return isl_bool_error);
3534
2.88k
  if (tab->rational)
3535
0
    return isl_bool_false;
3536
2.88k
3537
2.88k
  return get_constant(tab, &tab->var[var], value);
3538
2.88k
}
3539
3540
/* Check if any of the variables of "tab" can only attain a single (integer)
3541
 * value, and, if so, add equality constraints to fix those variables
3542
 * to these single values.
3543
 *
3544
 * For rational tableaus, nothing needs to be done.
3545
 */
3546
isl_stat isl_tab_detect_constants(struct isl_tab *tab)
3547
436
{
3548
436
  int i;
3549
436
3550
436
  if (!tab)
3551
0
    return isl_stat_error;
3552
436
  if (tab->rational)
3553
0
    return isl_stat_ok;
3554
436
3555
2.46k
  
for (i = 0; 436
i < tab->n_var;
++i2.03k
) {
3556
2.03k
    if (get_constant(tab, &tab->var[i], NULL) < 0)
3557
0
      return isl_stat_error;
3558
2.03k
  }
3559
436
3560
436
  return isl_stat_ok;
3561
436
}
3562
3563
/* Take a snapshot of the tableau that can be restored by a call to
3564
 * isl_tab_rollback.
3565
 */
3566
struct isl_tab_undo *isl_tab_snap(struct isl_tab *tab)
3567
1.09M
{
3568
1.09M
  if (!tab)
3569
0
    return NULL;
3570
1.09M
  tab->need_undo = 1;
3571
1.09M
  return tab->top;
3572
1.09M
}
3573
3574
/* Does "tab" need to keep track of undo information?
3575
 * That is, was a snapshot taken that may need to be restored?
3576
 */
3577
isl_bool isl_tab_need_undo(struct isl_tab *tab)
3578
115
{
3579
115
  if (!tab)
3580
0
    return isl_bool_error;
3581
115
3582
115
  return tab->need_undo;
3583
115
}
3584
3585
/* Remove all tracking of undo information from "tab", invalidating
3586
 * any snapshots that may have been taken of the tableau.
3587
 * Since all snapshots have been invalidated, there is also
3588
 * no need to start keeping track of undo information again.
3589
 */
3590
void isl_tab_clear_undo(struct isl_tab *tab)
3591
115
{
3592
115
  if (!tab)
3593
0
    return;
3594
115
3595
115
  free_undo(tab);
3596
115
  tab->need_undo = 0;
3597
115
}
3598
3599
/* Undo the operation performed by isl_tab_relax.
3600
 */
3601
static isl_stat unrelax(struct isl_tab *tab, struct isl_tab_var *var)
3602
  WARN_UNUSED;
3603
static isl_stat unrelax(struct isl_tab *tab, struct isl_tab_var *var)
3604
1.09k
{
3605
1.09k
  unsigned off = 2 + tab->M;
3606
1.09k
3607
1.09k
  if (!var->is_row && 
!max_is_manifestly_unbounded(tab, var)1.08k
)
3608
85
    if (to_row(tab, var, 1) < 0)
3609
0
      return isl_stat_error;
3610
1.09k
3611
1.09k
  if (var->is_row) {
3612
90
    isl_int_sub(tab->mat->row[var->index][1],
3613
90
        tab->mat->row[var->index][1], tab->mat->row[var->index][0]);
3614
90
    if (var->is_nonneg) {
3615
90
      int sgn = restore_row(tab, var);
3616
90
      isl_assert(tab->mat->ctx, sgn >= 0,
3617
90
        return isl_stat_error);
3618
90
    }
3619
1.00k
  } else {
3620
1.00k
    int i;
3621
1.00k
3622
5.70k
    for (i = 0; i < tab->n_row; 
++i4.70k
) {
3623
4.70k
      if (isl_int_is_zero(tab->mat->row[i][off + var->index]))
3624
4.70k
        
continue3.64k
;
3625
1.05k
      isl_int_add(tab->mat->row[i][1], tab->mat->row[i][1],
3626
1.05k
          tab->mat->row[i][off + var->index]);
3627
1.05k
    }
3628
1.00k
3629
1.00k
  }
3630
1.09k
3631
1.09k
  return isl_stat_ok;
3632
1.09k
}
3633
3634
/* Undo the operation performed by isl_tab_unrestrict.
3635
 *
3636
 * In particular, mark the variable as being non-negative and make
3637
 * sure the sample value respects this constraint.
3638
 */
3639
static isl_stat ununrestrict(struct isl_tab *tab, struct isl_tab_var *var)
3640
413
{
3641
413
  var->is_nonneg = 1;
3642
413
3643
413
  if (var->is_row && 
restore_row(tab, var) < -1338
)
3644
0
    return isl_stat_error;
3645
413
3646
413
  return isl_stat_ok;
3647
413
}
3648
3649
/* Unmark the last redundant row in "tab" as being redundant.
3650
 * This undoes part of the modifications performed by isl_tab_mark_redundant.
3651
 * In particular, remove the redundant mark and make
3652
 * sure the sample value respects the constraint again.
3653
 * A variable that is marked non-negative by isl_tab_mark_redundant
3654
 * is covered by a separate undo record.
3655
 */
3656
static isl_stat restore_last_redundant(struct isl_tab *tab)
3657
665k
{
3658
665k
  struct isl_tab_var *var;
3659
665k
3660
665k
  if (tab->n_redundant < 1)
3661
665k
    
isl_die0
(isl_tab_get_ctx(tab), isl_error_internal,
3662
665k
      "no redundant rows", return isl_stat_error);
3663
665k
3664
665k
  var = isl_tab_var_from_row(tab, tab->n_redundant - 1);
3665
665k
  var->is_redundant = 0;
3666
665k
  tab->n_redundant--;
3667
665k
  restore_row(tab, var);
3668
665k
3669
665k
  return isl_stat_ok;
3670
665k
}
3671
3672
static isl_stat perform_undo_var(struct isl_tab *tab, struct isl_tab_undo *undo)
3673
  WARN_UNUSED;
3674
static isl_stat perform_undo_var(struct isl_tab *tab, struct isl_tab_undo *undo)
3675
2.46M
{
3676
2.46M
  struct isl_tab_var *var = var_from_index(tab, undo->u.var_index);
3677
2.46M
  switch (undo->type) {
3678
2.46M
  case isl_tab_undo_nonneg:
3679
519k
    var->is_nonneg = 0;
3680
519k
    break;
3681
2.46M
  case isl_tab_undo_redundant:
3682
527k
    if (!var->is_row || var->index != tab->n_redundant - 1)
3683
527k
      
isl_die0
(isl_tab_get_ctx(tab), isl_error_internal,
3684
527k
        "not undoing last redundant row", return -1);
3685
527k
    return restore_last_redundant(tab);
3686
527k
  case isl_tab_undo_freeze:
3687
259k
    var->frozen = 0;
3688
259k
    break;
3689
527k
  case isl_tab_undo_zero:
3690
82.7k
    var->is_zero = 0;
3691
82.7k
    if (!var->is_row)
3692
82.4k
      tab->n_dead--;
3693
82.7k
    break;
3694
1.07M
  case isl_tab_undo_allocate:
3695
1.07M
    if (undo->u.var_index >= 0) {
3696
3.43k
      isl_assert(tab->mat->ctx, !var->is_row,
3697
3.43k
        return isl_stat_error);
3698
3.43k
      return drop_col(tab, var->index);
3699
1.06M
    }
3700
1.06M
    if (!var->is_row) {
3701
116k
      if (!max_is_manifestly_unbounded(tab, var)) {
3702
84.3k
        if (to_row(tab, var, 1) < 0)
3703
0
          return isl_stat_error;
3704
31.6k
      } else if (!min_is_manifestly_unbounded(tab, var)) {
3705
13.7k
        if (to_row(tab, var, -1) < 0)
3706
0
          return isl_stat_error;
3707
17.9k
      } else
3708
17.9k
        if (to_row(tab, var, 0) < 0)
3709
0
          return isl_stat_error;
3710
1.06M
    }
3711
1.06M
    return drop_row(tab, var->index);
3712
1.06M
  case isl_tab_undo_relax:
3713
1.09k
    return unrelax(tab, var);
3714
1.06M
  case isl_tab_undo_unrestrict:
3715
413
    return ununrestrict(tab, var);
3716
1.06M
  default:
3717
0
    isl_die(tab->mat->ctx, isl_error_internal,
3718
2.46M
      "perform_undo_var called on invalid undo record",
3719
2.46M
      return isl_stat_error);
3720
2.46M
  }
3721
2.46M
3722
2.46M
  
return isl_stat_ok862k
;
3723
2.46M
}
3724
3725
/* Restore all rows that have been marked redundant by isl_tab_mark_redundant
3726
 * and that have been preserved in the tableau.
3727
 * Note that isl_tab_mark_redundant may also have marked some variables
3728
 * as being non-negative before marking them redundant.  These need
3729
 * to be removed as well as otherwise some constraints could end up
3730
 * getting marked redundant with respect to the variable.
3731
 */
3732
isl_stat isl_tab_restore_redundant(struct isl_tab *tab)
3733
112k
{
3734
112k
  if (!tab)
3735
0
    return isl_stat_error;
3736
112k
3737
112k
  if (tab->need_undo)
3738
112k
    
isl_die0
(isl_tab_get_ctx(tab), isl_error_invalid,
3739
112k
      "manually restoring redundant constraints "
3740
112k
      "interferes with undo history",
3741
112k
      return isl_stat_error);
3742
112k
3743
250k
  
while (112k
tab->n_redundant > 0) {
3744
138k
    if (tab->row_var[tab->n_redundant - 1] >= 0) {
3745
121k
      struct isl_tab_var *var;
3746
121k
3747
121k
      var = isl_tab_var_from_row(tab, tab->n_redundant - 1);
3748
121k
      var->is_nonneg = 0;
3749
121k
    }
3750
138k
    restore_last_redundant(tab);
3751
138k
  }
3752
112k
  return isl_stat_ok;
3753
112k
}
3754
3755
/* Undo the addition of an integer division to the basic map representation
3756
 * of "tab" in position "pos".
3757
 */
3758
static isl_stat drop_bmap_div(struct isl_tab *tab, int pos)
3759
1.18k
{
3760
1.18k
  int off;
3761
1.18k
3762
1.18k
  off = tab->n_var - isl_basic_map_dim(tab->bmap, isl_dim_div);
3763
1.18k
  if (isl_basic_map_drop_div(tab->bmap, pos - off) < 0)
3764
0
    return isl_stat_error;
3765
1.18k
  if (tab->samples) {
3766
468
    tab->samples = isl_mat_drop_cols(tab->samples, 1 + pos, 1);
3767
468
    if (!tab->samples)
3768
0
      return isl_stat_error;
3769
1.18k
  }
3770
1.18k
3771
1.18k
  return isl_stat_ok;
3772
1.18k
}
3773
3774
/* Restore the tableau to the state where the basic variables
3775
 * are those in "col_var".
3776
 * We first construct a list of variables that are currently in
3777
 * the basis, but shouldn't.  Then we iterate over all variables
3778
 * that should be in the basis and for each one that is currently
3779
 * not in the basis, we exchange it with one of the elements of the
3780
 * list constructed before.
3781
 * We can always find an appropriate variable to pivot with because
3782
 * the current basis is mapped to the old basis by a non-singular
3783
 * matrix and so we can never end up with a zero row.
3784
 */
3785
static int restore_basis(struct isl_tab *tab, int *col_var)
3786
587
{
3787
587
  int i, j;
3788
587
  int n_extra = 0;
3789
587
  int *extra = NULL;  /* current columns that contain bad stuff */
3790
587
  unsigned off = 2 + tab->M;
3791
587
3792
587
  extra = isl_alloc_array(tab->mat->ctx, int, tab->n_col);
3793
587
  if (tab->n_col && !extra)
3794
0
    goto error;
3795
6.92k
  
for (i = 0; 587
i < tab->n_col;
++i6.33k
) {
3796
51.6k
    for (j = 0; j < tab->n_col; 
++j45.3k
)
3797
50.2k
      if (tab->col_var[i] == col_var[j])
3798
4.84k
        break;
3799
6.33k
    if (j < tab->n_col)
3800
4.84k
      continue;
3801
1.48k
    extra[n_extra++] = i;
3802
1.48k
  }
3803
5.41k
  for (i = 0; i < tab->n_col && 
n_extra > 05.30k
;
++i4.83k
) {
3804
4.83k
    struct isl_tab_var *var;
3805
4.83k
    int row;
3806
4.83k
3807
37.0k
    for (j = 0; j < tab->n_col; 
++j32.2k
)
3808
35.5k
      if (col_var[i] == tab->col_var[j])
3809
3.34k
        break;
3810
4.83k
    if (j < tab->n_col)
3811
3.34k
      continue;
3812
1.48k
    var = var_from_index(tab, col_var[i]);
3813
1.48k
    row = var->index;
3814
1.88k
    for (j = 0; j < n_extra; 
++j400
)
3815
1.88k
      if (!isl_int_is_zero(tab->mat->row[row][off+extra[j]]))
3816
1.88k
        
break1.48k
;
3817
1.48k
    isl_assert(tab->mat->ctx, j < n_extra, goto error);
3818
1.48k
    if (isl_tab_pivot(tab, row, extra[j]) < 0)
3819
0
      goto error;
3820
1.48k
    extra[j] = extra[--n_extra];
3821
1.48k
  }
3822
587
3823
587
  free(extra);
3824
587
  return 0;
3825
0
error:
3826
0
  free(extra);
3827
0
  return -1;
3828
587
}
3829
3830
/* Remove all samples with index n or greater, i.e., those samples
3831
 * that were added since we saved this number of samples in
3832
 * isl_tab_save_samples.
3833
 */
3834
static void drop_samples_since(struct isl_tab *tab, int n)
3835
23.7k
{
3836
23.7k
  int i;
3837
23.7k
3838
28.9k
  for (i = tab->n_sample - 1; i >= 0 && 
tab->n_sample > n27.2k
;
--i5.19k
) {
3839
5.19k
    if (tab->sample_index[i] < n)
3840
2.00k
      continue;
3841
3.18k
3842
3.18k
    if (i != tab->n_sample - 1) {
3843
2.20k
      int t = tab->sample_index[tab->n_sample-1];
3844
2.20k
      tab->sample_index[tab->n_sample-1] = tab->sample_index[i];
3845
2.20k
      tab->sample_index[i] = t;
3846
2.20k
      isl_mat_swap_rows(tab->samples, tab->n_sample-1, i);
3847
2.20k
    }
3848
3.18k
    tab->n_sample--;
3849
3.18k
  }
3850
23.7k
}
3851
3852
static isl_stat perform_undo(struct isl_tab *tab, struct isl_tab_undo *undo)
3853
  WARN_UNUSED;
3854
static isl_stat perform_undo(struct isl_tab *tab, struct isl_tab_undo *undo)
3855
2.82M
{
3856
2.82M
  switch (undo->type) {
3857
2.82M
  case isl_tab_undo_rational:
3858
10.3k
    tab->rational = 0;
3859
10.3k
    break;
3860
2.82M
  case isl_tab_undo_empty:
3861
63.0k
    tab->empty = 0;
3862
63.0k
    break;
3863
2.82M
  case isl_tab_undo_nonneg:
3864
2.46M
  case isl_tab_undo_redundant:
3865
2.46M
  case isl_tab_undo_freeze:
3866
2.46M
  case isl_tab_undo_zero:
3867
2.46M
  case isl_tab_undo_allocate:
3868
2.46M
  case isl_tab_undo_relax:
3869
2.46M
  case isl_tab_undo_unrestrict:
3870
2.46M
    return perform_undo_var(tab, undo);
3871
2.46M
  case isl_tab_undo_bmap_eq:
3872
0
    return isl_basic_map_free_equality(tab->bmap, 1);
3873
2.46M
  case isl_tab_undo_bmap_ineq:
3874
255k
    return isl_basic_map_free_inequality(tab->bmap, 1);
3875
2.46M
  case isl_tab_undo_bmap_div:
3876
1.18k
    return drop_bmap_div(tab, undo->u.var_index);
3877
2.46M
  case isl_tab_undo_saved_basis:
3878
587
    if (restore_basis(tab, undo->u.col_var) < 0)
3879
0
      return isl_stat_error;
3880
587
    break;
3881
5.69k
  case isl_tab_undo_drop_sample:
3882
5.69k
    tab->n_outside--;
3883
5.69k
    break;
3884
23.7k
  case isl_tab_undo_saved_samples:
3885
23.7k
    drop_samples_since(tab, undo->u.n);
3886
23.7k
    break;
3887
2.39k
  case isl_tab_undo_callback:
3888
2.39k
    return undo->u.callback->run(undo->u.callback);
3889
587
  default:
3890
0
    isl_assert(tab->mat->ctx, 0, return isl_stat_error);
3891
2.82M
  }
3892
2.82M
  
return isl_stat_ok103k
;
3893
2.82M
}
3894
3895
/* Return the tableau to the state it was in when the snapshot "snap"
3896
 * was taken.
3897
 */
3898
int isl_tab_rollback(struct isl_tab *tab, struct isl_tab_undo *snap)
3899
991k
{
3900
991k
  struct isl_tab_undo *undo, *next;
3901
991k
3902
991k
  if (!tab)
3903
0
    return -1;
3904
991k
3905
991k
  tab->in_undo = 1;
3906
3.81M
  for (undo = tab->top; undo && undo != &tab->bottom; 
undo = next2.82M
) {
3907
3.10M
    next = undo->next;
3908
3.10M
    if (undo == snap)
3909
284k
      break;
3910
2.82M
    if (perform_undo(tab, undo) < 0) {
3911
0
      tab->top = undo;
3912
0
      free_undo(tab);
3913
0
      tab->in_undo = 0;
3914
0
      return -1;
3915
0
    }
3916
2.82M
    free_undo_record(undo);
3917
2.82M
  }
3918
991k
  tab->in_undo = 0;
3919
991k
  tab->top = undo;
3920
991k
  if (!undo)
3921
0
    return -1;
3922
991k
  return 0;
3923
991k
}
3924
3925
/* The given row "row" represents an inequality violated by all
3926
 * points in the tableau.  Check for some special cases of such
3927
 * separating constraints.
3928
 * In particular, if the row has been reduced to the constant -1,
3929
 * then we know the inequality is adjacent (but opposite) to
3930
 * an equality in the tableau.
3931
 * If the row has been reduced to r = c*(-1 -r'), with r' an inequality
3932
 * of the tableau and c a positive constant, then the inequality
3933
 * is adjacent (but opposite) to the inequality r'.
3934
 */
3935
static enum isl_ineq_type separation_type(struct isl_tab *tab, unsigned row)
3936
77.0k
{
3937
77.0k
  int pos;
3938
77.0k
  unsigned off = 2 + tab->M;
3939
77.0k
3940
77.0k
  if (tab->rational)
3941
9.20k
    return isl_ineq_separate;
3942
67.8k
3943
67.8k
  if (!isl_int_is_one(tab->mat->row[row][0]))
3944
67.8k
    
return isl_ineq_separate173
;
3945
67.6k
3946
67.6k
  pos = isl_seq_first_non_zero(tab->mat->row[row] + off + tab->n_dead,
3947
67.6k
          tab->n_col - tab->n_dead);
3948
67.6k
  if (pos == -1) {
3949
12.4k
    if (isl_int_is_negone(tab->mat->row[row][1]))
3950
12.4k
      
return isl_ineq_adj_eq10.7k
;
3951
1.72k
    else
3952
1.72k
      return isl_ineq_separate;
3953
55.1k
  }
3954
55.1k
3955
55.1k
  if (!isl_int_eq(tab->mat->row[row][1],
3956
55.1k
      tab->mat->row[row][off + tab->n_dead + pos]))
3957
55.1k
    
return isl_ineq_separate20.7k
;
3958
34.4k
3959
34.4k
  pos = isl_seq_first_non_zero(
3960
34.4k
      tab->mat->row[row] + off + tab->n_dead + pos + 1,
3961
34.4k
      tab->n_col - tab->n_dead - pos - 1);
3962
34.4k
3963
34.4k
  return pos == -1 ? 
isl_ineq_adj_ineq32.2k
:
isl_ineq_separate2.19k
;
3964
34.4k
}
3965
3966
/* Check the effect of inequality "ineq" on the tableau "tab".
3967
 * The result may be
3968
 *  isl_ineq_redundant: satisfied by all points in the tableau
3969
 *  isl_ineq_separate:  satisfied by no point in the tableau
3970
 *  isl_ineq_cut:   satisfied by some by not all points
3971
 *  isl_ineq_adj_eq:  adjacent to an equality
3972
 *  isl_ineq_adj_ineq:  adjacent to an inequality.
3973
 */
3974
enum isl_ineq_type isl_tab_ineq_type(struct isl_tab *tab, isl_int *ineq)
3975
417k
{
3976
417k
  enum isl_ineq_type type = isl_ineq_error;
3977
417k
  struct isl_tab_undo *snap = NULL;
3978
417k
  int con;
3979
417k
  int row;
3980
417k
3981
417k
  if (!tab)
3982
0
    return isl_ineq_error;
3983
417k
3984
417k
  if (isl_tab_extend_cons(tab, 1) < 0)
3985
0
    return isl_ineq_error;
3986
417k
3987
417k
  snap = isl_tab_snap(tab);
3988
417k
3989
417k
  con = isl_tab_add_row(tab, ineq);
3990
417k
  if (con < 0)
3991
0
    goto error;
3992
417k
3993
417k
  row = tab->con[con].index;
3994
417k
  if (isl_tab_row_is_redundant(tab, row))
3995
0
    type = isl_ineq_redundant;
3996
417k
  else if (isl_int_is_neg(tab->mat->row[row][1]) &&
3997
417k
     
(157k
tab->rational157k
||
3998
157k
        
isl_int_abs_ge134k
(tab->mat->row[row][1],
3999
157k
           tab->mat->row[row][0]))) {
4000
156k
    int nonneg = at_least_zero(tab, &tab->con[con]);
4001
156k
    if (nonneg < 0)
4002
0
      goto error;
4003
156k
    if (nonneg)
4004
79.5k
      type = isl_ineq_cut;
4005
77.0k
    else
4006
77.0k
      type = separation_type(tab, row);
4007
261k
  } else {
4008
261k
    int red = con_is_redundant(tab, &tab->con[con]);
4009
261k
    if (red < 0)
4010
0
      goto error;
4011
261k
    if (!red)
4012
66.8k
      type = isl_ineq_cut;
4013
194k
    else
4014
194k
      type = isl_ineq_redundant;
4015
261k
  }
4016
417k
4017
417k
  if (isl_tab_rollback(tab, snap))
4018
0
    return isl_ineq_error;
4019
417k
  return type;
4020
0
error:
4021
0
  return isl_ineq_error;
4022
0
}
4023
4024
isl_stat isl_tab_track_bmap(struct isl_tab *tab, __isl_take isl_basic_map *bmap)
4025
190k
{
4026
190k
  bmap = isl_basic_map_cow(bmap);
4027
190k
  if (!tab || !bmap)
4028
0
    goto error;
4029
190k
4030
190k
  if (tab->empty) {
4031
3.77k
    bmap = isl_basic_map_set_to_empty(bmap);
4032
3.77k
    if (!bmap)
4033
0
      goto error;
4034
3.77k
    tab->bmap = bmap;
4035
3.77k
    return isl_stat_ok;
4036
3.77k
  }
4037
186k
4038
186k
  isl_assert(tab->mat->ctx, tab->n_eq == bmap->n_eq, goto error);
4039
186k
  isl_assert(tab->mat->ctx,
4040
186k
        tab->n_con == bmap->n_eq + bmap->n_ineq, goto error);
4041
186k
4042
186k
  tab->bmap = bmap;
4043
186k
4044
186k
  return isl_stat_ok;
4045
0
error:
4046
0
  isl_basic_map_free(bmap);
4047
0
  return isl_stat_error;
4048
186k
}
4049
4050
isl_stat isl_tab_track_bset(struct isl_tab *tab, __isl_take isl_basic_set *bset)
4051
885
{
4052
885
  return isl_tab_track_bmap(tab, bset_to_bmap(bset));
4053
885
}
4054
4055
__isl_keep isl_basic_set *isl_tab_peek_bset(struct isl_tab *tab)
4056
29.9k
{
4057
29.9k
  if (!tab)
4058
0
    return NULL;
4059
29.9k
4060
29.9k
  return bset_from_bmap(tab->bmap);
4061
29.9k
}
4062
4063
static void isl_tab_print_internal(__isl_keep struct isl_tab *tab,
4064
  FILE *out, int indent)
4065
0
{
4066
0
  unsigned r, c;
4067
0
  int i;
4068
0
4069
0
  if (!tab) {
4070
0
    fprintf(out, "%*snull tab\n", indent, "");
4071
0
    return;
4072
0
  }
4073
0
  fprintf(out, "%*sn_redundant: %d, n_dead: %d", indent, "",
4074
0
    tab->n_redundant, tab->n_dead);
4075
0
  if (tab->rational)
4076
0
    fprintf(out, ", rational");
4077
0
  if (tab->empty)
4078
0
    fprintf(out, ", empty");
4079
0
  fprintf(out, "\n");
4080
0
  fprintf(out, "%*s[", indent, "");
4081
0
  for (i = 0; i < tab->n_var; ++i) {
4082
0
    if (i)
4083
0
      fprintf(out, (i == tab->n_param ||
4084
0
              i == tab->n_var - tab->n_div) ? "; "
4085
0
                    : ", ");
4086
0
    fprintf(out, "%c%d%s", tab->var[i].is_row ? 'r' : 'c',
4087
0
          tab->var[i].index,
4088
0
          tab->var[i].is_zero ? " [=0]" :
4089
0
          tab->var[i].is_redundant ? " [R]" : "");
4090
0
  }
4091
0
  fprintf(out, "]\n");
4092
0
  fprintf(out, "%*s[", indent, "");
4093
0
  for (i = 0; i < tab->n_con; ++i) {
4094
0
    if (i)
4095
0
      fprintf(out, ", ");
4096
0
    fprintf(out, "%c%d%s", tab->con[i].is_row ? 'r' : 'c',
4097
0
          tab->con[i].index,
4098
0
          tab->con[i].is_zero ? " [=0]" :
4099
0
          tab->con[i].is_redundant ? " [R]" : "");
4100
0
  }
4101
0
  fprintf(out, "]\n");
4102
0
  fprintf(out, "%*s[", indent, "");
4103
0
  for (i = 0; i < tab->n_row; ++i) {
4104
0
    const char *sign = "";
4105
0
    if (i)
4106
0
      fprintf(out, ", ");
4107
0
    if (tab->row_sign) {
4108
0
      if (tab->row_sign[i] == isl_tab_row_unknown)
4109
0
        sign = "?";
4110
0
      else if (tab->row_sign[i] == isl_tab_row_neg)
4111
0
        sign = "-";
4112
0
      else if (tab->row_sign[i] == isl_tab_row_pos)
4113
0
        sign = "+";
4114
0
      else
4115
0
        sign = "+-";
4116
0
    }
4117
0
    fprintf(out, "r%d: %d%s%s", i, tab->row_var[i],
4118
0
        isl_tab_var_from_row(tab, i)->is_nonneg ? " [>=0]" : "", sign);
4119
0
  }
4120
0
  fprintf(out, "]\n");
4121
0
  fprintf(out, "%*s[", indent, "");
4122
0
  for (i = 0; i < tab->n_col; ++i) {
4123
0
    if (i)
4124
0
      fprintf(out, ", ");
4125
0
    fprintf(out, "c%d: %d%s", i, tab->col_var[i],
4126
0
        var_from_col(tab, i)->is_nonneg ? " [>=0]" : "");
4127
0
  }
4128
0
  fprintf(out, "]\n");
4129
0
  r = tab->mat->n_row;
4130
0
  tab->mat->n_row = tab->n_row;
4131
0
  c = tab->mat->n_col;
4132
0
  tab->mat->n_col = 2 + tab->M + tab->n_col;
4133
0
  isl_mat_print_internal(tab->mat, out, indent);
4134
0
  tab->mat->n_row = r;
4135
0
  tab->mat->n_col = c;
4136
0
  if (tab->bmap)
4137
0
    isl_basic_map_print_internal(tab->bmap, out, indent);
4138
0
}
4139
4140
void isl_tab_dump(__isl_keep struct isl_tab *tab)
4141
0
{
4142
0
  isl_tab_print_internal(tab, stderr, 0);
4143
0
}