Line data Source code
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 2955142263 : struct isl_tab *isl_tab_alloc(struct isl_ctx *ctx,
34 : unsigned n_row, unsigned n_var, unsigned M)
35 : {
36 : int i;
37 : struct isl_tab *tab;
38 2955142263 : unsigned off = 2 + M;
39 :
40 2955142263 : tab = isl_calloc_type(ctx, struct isl_tab);
41 2955142263 : if (!tab)
42 0 : return NULL;
43 2955142263 : tab->mat = isl_mat_alloc(ctx, n_row, off + n_var);
44 2955142263 : if (!tab->mat)
45 0 : goto error;
46 2955142263 : tab->var = isl_alloc_array(ctx, struct isl_tab_var, n_var);
47 2955142263 : if (n_var && !tab->var)
48 0 : goto error;
49 2955142263 : tab->con = isl_alloc_array(ctx, struct isl_tab_var, n_row);
50 2955142263 : if (n_row && !tab->con)
51 0 : goto error;
52 2955142263 : tab->col_var = isl_alloc_array(ctx, int, n_var);
53 2955142263 : if (n_var && !tab->col_var)
54 0 : goto error;
55 2955142263 : tab->row_var = isl_alloc_array(ctx, int, n_row);
56 2955142263 : if (n_row && !tab->row_var)
57 0 : goto error;
58 16009677495 : for (i = 0; i < n_var; ++i) {
59 13054535232 : tab->var[i].index = i;
60 13054535232 : tab->var[i].is_row = 0;
61 13054535232 : tab->var[i].is_nonneg = 0;
62 13054535232 : tab->var[i].is_zero = 0;
63 13054535232 : tab->var[i].is_redundant = 0;
64 13054535232 : tab->var[i].frozen = 0;
65 13054535232 : tab->var[i].negated = 0;
66 13054535232 : tab->col_var[i] = i;
67 : }
68 2955142263 : tab->n_row = 0;
69 2955142263 : tab->n_con = 0;
70 2955142263 : tab->n_eq = 0;
71 2955142263 : tab->max_con = n_row;
72 2955142263 : tab->n_col = n_var;
73 2955142263 : tab->n_var = n_var;
74 2955142263 : tab->max_var = n_var;
75 2955142263 : tab->n_param = 0;
76 2955142263 : tab->n_div = 0;
77 2955142263 : tab->n_dead = 0;
78 2955142263 : tab->n_redundant = 0;
79 2955142263 : tab->strict_redundant = 0;
80 2955142263 : tab->need_undo = 0;
81 2955142263 : tab->rational = 0;
82 2955142263 : tab->empty = 0;
83 2955142263 : tab->in_undo = 0;
84 2955142263 : tab->M = M;
85 2955142263 : tab->cone = 0;
86 2955142263 : tab->bottom.type = isl_tab_undo_bottom;
87 2955142263 : tab->bottom.next = NULL;
88 2955142263 : tab->top = &tab->bottom;
89 :
90 2955142263 : tab->n_zero = 0;
91 2955142263 : tab->n_unbounded = 0;
92 2955142263 : tab->basis = NULL;
93 :
94 2955142263 : return tab;
95 : error:
96 0 : isl_tab_free(tab);
97 0 : return NULL;
98 : }
99 :
100 20210094836 : isl_ctx *isl_tab_get_ctx(struct isl_tab *tab)
101 : {
102 20210094836 : return tab ? isl_mat_get_ctx(tab->mat) : NULL;
103 : }
104 :
105 6300884808 : int isl_tab_extend_cons(struct isl_tab *tab, unsigned n_new)
106 : {
107 : unsigned off;
108 :
109 6300884808 : if (!tab)
110 0 : return -1;
111 :
112 6300884808 : off = 2 + tab->M;
113 :
114 6300884808 : if (tab->max_con < tab->n_con + n_new) {
115 : struct isl_tab_var *con;
116 :
117 18742991 : con = isl_realloc_array(tab->mat->ctx, tab->con,
118 : struct isl_tab_var, tab->max_con + n_new);
119 18742991 : if (!con)
120 0 : return -1;
121 18742991 : tab->con = con;
122 18742991 : tab->max_con += n_new;
123 : }
124 6300884808 : if (tab->mat->n_row < tab->n_row + n_new) {
125 : int *row_var;
126 :
127 108135736 : tab->mat = isl_mat_extend(tab->mat,
128 108135736 : tab->n_row + n_new, off + tab->n_col);
129 54067868 : if (!tab->mat)
130 0 : return -1;
131 54067868 : row_var = isl_realloc_array(tab->mat->ctx, tab->row_var,
132 : int, tab->mat->n_row);
133 54067868 : if (!row_var)
134 0 : return -1;
135 54067868 : tab->row_var = row_var;
136 54067868 : if (tab->row_sign) {
137 : enum isl_tab_row_sign *s;
138 0 : s = isl_realloc_array(tab->mat->ctx, tab->row_sign,
139 : enum isl_tab_row_sign, tab->mat->n_row);
140 0 : if (!s)
141 0 : return -1;
142 0 : tab->row_sign = s;
143 : }
144 : }
145 6300884808 : return 0;
146 : }
147 :
148 : /* Make room for at least n_new extra variables.
149 : * Return -1 if anything went wrong.
150 : */
151 0 : int isl_tab_extend_vars(struct isl_tab *tab, unsigned n_new)
152 : {
153 : struct isl_tab_var *var;
154 0 : unsigned off = 2 + tab->M;
155 :
156 0 : if (tab->max_var < tab->n_var + n_new) {
157 0 : var = isl_realloc_array(tab->mat->ctx, tab->var,
158 : struct isl_tab_var, tab->n_var + n_new);
159 0 : if (!var)
160 0 : return -1;
161 0 : tab->var = var;
162 0 : tab->max_var = tab->n_var + n_new;
163 : }
164 :
165 0 : if (tab->mat->n_col < off + tab->n_col + n_new) {
166 : int *p;
167 :
168 0 : tab->mat = isl_mat_extend(tab->mat,
169 0 : tab->mat->n_row, off + tab->n_col + n_new);
170 0 : if (!tab->mat)
171 0 : return -1;
172 0 : p = isl_realloc_array(tab->mat->ctx, tab->col_var,
173 : int, tab->n_col + n_new);
174 0 : if (!p)
175 0 : return -1;
176 0 : tab->col_var = p;
177 : }
178 :
179 0 : return 0;
180 : }
181 :
182 >12568*10^7 : static void free_undo_record(struct isl_tab_undo *undo)
183 : {
184 >12568*10^7 : switch (undo->type) {
185 : case isl_tab_undo_saved_basis:
186 0 : free(undo->u.col_var);
187 0 : break;
188 : default:;
189 : }
190 >12568*10^7 : free(undo);
191 >12568*10^7 : }
192 :
193 2955142559 : static void free_undo(struct isl_tab *tab)
194 : {
195 : struct isl_tab_undo *undo, *next;
196 :
197 3648259254 : for (undo = tab->top; undo && undo != &tab->bottom; undo = next) {
198 693116695 : next = undo->next;
199 693116695 : free_undo_record(undo);
200 : }
201 2955142559 : tab->top = undo;
202 2955142559 : }
203 :
204 2955625229 : void isl_tab_free(struct isl_tab *tab)
205 : {
206 2955625229 : if (!tab)
207 482670 : return;
208 2955142559 : free_undo(tab);
209 2955142559 : isl_mat_free(tab->mat);
210 2955142559 : isl_vec_free(tab->dual);
211 2955142559 : isl_basic_map_free(tab->bmap);
212 2955142559 : free(tab->var);
213 2955142559 : free(tab->con);
214 2955142559 : free(tab->row_var);
215 2955142559 : free(tab->col_var);
216 2955142559 : free(tab->row_sign);
217 2955142559 : isl_mat_free(tab->samples);
218 2955142559 : free(tab->sample_index);
219 2955142559 : isl_mat_free(tab->basis);
220 2955142559 : free(tab);
221 : }
222 :
223 0 : struct isl_tab *isl_tab_dup(struct isl_tab *tab)
224 : {
225 : int i;
226 : struct isl_tab *dup;
227 : unsigned off;
228 :
229 0 : if (!tab)
230 0 : return NULL;
231 :
232 0 : off = 2 + tab->M;
233 0 : dup = isl_calloc_type(tab->mat->ctx, struct isl_tab);
234 0 : if (!dup)
235 0 : return NULL;
236 0 : dup->mat = isl_mat_dup(tab->mat);
237 0 : if (!dup->mat)
238 0 : goto error;
239 0 : dup->var = isl_alloc_array(tab->mat->ctx, struct isl_tab_var, tab->max_var);
240 0 : if (tab->max_var && !dup->var)
241 0 : goto error;
242 0 : for (i = 0; i < tab->n_var; ++i)
243 0 : dup->var[i] = tab->var[i];
244 0 : dup->con = isl_alloc_array(tab->mat->ctx, struct isl_tab_var, tab->max_con);
245 0 : if (tab->max_con && !dup->con)
246 0 : goto error;
247 0 : for (i = 0; i < tab->n_con; ++i)
248 0 : dup->con[i] = tab->con[i];
249 0 : dup->col_var = isl_alloc_array(tab->mat->ctx, int, tab->mat->n_col - off);
250 0 : if ((tab->mat->n_col - off) && !dup->col_var)
251 0 : goto error;
252 0 : for (i = 0; i < tab->n_col; ++i)
253 0 : dup->col_var[i] = tab->col_var[i];
254 0 : dup->row_var = isl_alloc_array(tab->mat->ctx, int, tab->mat->n_row);
255 0 : if (tab->mat->n_row && !dup->row_var)
256 0 : goto error;
257 0 : for (i = 0; i < tab->n_row; ++i)
258 0 : dup->row_var[i] = tab->row_var[i];
259 0 : if (tab->row_sign) {
260 0 : dup->row_sign = isl_alloc_array(tab->mat->ctx, enum isl_tab_row_sign,
261 : tab->mat->n_row);
262 0 : if (tab->mat->n_row && !dup->row_sign)
263 0 : goto error;
264 0 : for (i = 0; i < tab->n_row; ++i)
265 0 : dup->row_sign[i] = tab->row_sign[i];
266 : }
267 0 : 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 : 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 : }
278 0 : dup->n_row = tab->n_row;
279 0 : dup->n_con = tab->n_con;
280 0 : dup->n_eq = tab->n_eq;
281 0 : dup->max_con = tab->max_con;
282 0 : dup->n_col = tab->n_col;
283 0 : dup->n_var = tab->n_var;
284 0 : dup->max_var = tab->max_var;
285 0 : dup->n_param = tab->n_param;
286 0 : dup->n_div = tab->n_div;
287 0 : dup->n_dead = tab->n_dead;
288 0 : dup->n_redundant = tab->n_redundant;
289 0 : dup->rational = tab->rational;
290 0 : dup->empty = tab->empty;
291 0 : dup->strict_redundant = 0;
292 0 : dup->need_undo = 0;
293 0 : dup->in_undo = 0;
294 0 : dup->M = tab->M;
295 0 : tab->cone = tab->cone;
296 0 : dup->bottom.type = isl_tab_undo_bottom;
297 0 : dup->bottom.next = NULL;
298 0 : dup->top = &dup->bottom;
299 :
300 0 : dup->n_zero = tab->n_zero;
301 0 : dup->n_unbounded = tab->n_unbounded;
302 0 : dup->basis = isl_mat_dup(tab->basis);
303 :
304 0 : return dup;
305 : error:
306 0 : isl_tab_free(dup);
307 0 : return NULL;
308 : }
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 296 : 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 : {
328 : int i;
329 : struct isl_mat *prod;
330 : unsigned n;
331 :
332 296 : prod = isl_mat_alloc(mat1->ctx, mat1->n_row + mat2->n_row,
333 296 : off + col1 + col2);
334 296 : if (!prod)
335 0 : return NULL;
336 :
337 296 : n = 0;
338 2211 : for (i = 0; i < r1; ++i) {
339 1915 : isl_seq_cpy(prod->row[n + i], mat1->row[i], off + d1);
340 1915 : isl_seq_clr(prod->row[n + i] + off + d1, d2);
341 3830 : isl_seq_cpy(prod->row[n + i] + off + d1 + d2,
342 1915 : mat1->row[i] + off + d1, col1 - d1);
343 1915 : isl_seq_clr(prod->row[n + i] + off + col1 + d1, col2 - d2);
344 : }
345 :
346 296 : n += r1;
347 2211 : for (i = 0; i < r2; ++i) {
348 1915 : isl_seq_cpy(prod->row[n + i], mat2->row[i], off);
349 1915 : isl_seq_clr(prod->row[n + i] + off, d1);
350 3830 : isl_seq_cpy(prod->row[n + i] + off + d1,
351 3830 : mat2->row[i] + off, d2);
352 1915 : isl_seq_clr(prod->row[n + i] + off + d1 + d2, col1 - d1);
353 3830 : isl_seq_cpy(prod->row[n + i] + off + col1 + d1,
354 1915 : mat2->row[i] + off + d2, col2 - d2);
355 : }
356 :
357 296 : n += r2;
358 1594 : for (i = 0; i < row1 - r1; ++i) {
359 1298 : isl_seq_cpy(prod->row[n + i], mat1->row[r1 + i], off + d1);
360 1298 : isl_seq_clr(prod->row[n + i] + off + d1, d2);
361 2596 : isl_seq_cpy(prod->row[n + i] + off + d1 + d2,
362 1298 : mat1->row[r1 + i] + off + d1, col1 - d1);
363 1298 : isl_seq_clr(prod->row[n + i] + off + col1 + d1, col2 - d2);
364 : }
365 :
366 296 : n += row1 - r1;
367 1594 : for (i = 0; i < row2 - r2; ++i) {
368 1298 : isl_seq_cpy(prod->row[n + i], mat2->row[r2 + i], off);
369 1298 : isl_seq_clr(prod->row[n + i] + off, d1);
370 2596 : isl_seq_cpy(prod->row[n + i] + off + d1,
371 2596 : mat2->row[r2 + i] + off, d2);
372 1298 : isl_seq_clr(prod->row[n + i] + off + d1 + d2, col1 - d1);
373 2596 : isl_seq_cpy(prod->row[n + i] + off + col1 + d1,
374 1298 : mat2->row[r2 + i] + off + d2, col2 - d2);
375 : }
376 :
377 296 : return prod;
378 : }
379 :
380 : /* Update the row or column index of a variable that corresponds
381 : * to a variable in the first input tableau.
382 : */
383 4596 : static void update_index1(struct isl_tab_var *var,
384 : unsigned r1, unsigned r2, unsigned d1, unsigned d2)
385 : {
386 4596 : if (var->index == -1)
387 494 : return;
388 4102 : if (var->is_row && var->index >= r1)
389 1298 : var->index += r2;
390 4102 : if (!var->is_row && var->index >= d1)
391 837 : var->index += d2;
392 : }
393 :
394 : /* Update the row or column index of a variable that corresponds
395 : * to a variable in the second input tableau.
396 : */
397 4596 : static void update_index2(struct isl_tab_var *var,
398 : unsigned row1, unsigned col1,
399 : unsigned r1, unsigned r2, unsigned d1, unsigned d2)
400 : {
401 4596 : if (var->index == -1)
402 494 : return;
403 4102 : if (var->is_row) {
404 3213 : if (var->index < r2)
405 1915 : var->index += r1;
406 : else
407 1298 : var->index += row1;
408 : } else {
409 889 : if (var->index < d2)
410 52 : var->index += d1;
411 : else
412 837 : var->index += col1;
413 : }
414 : }
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 296 : struct isl_tab *isl_tab_product(struct isl_tab *tab1, struct isl_tab *tab2)
435 : {
436 : int i;
437 : struct isl_tab *prod;
438 : unsigned off;
439 : unsigned r1, r2, d1, d2;
440 :
441 296 : if (!tab1 || !tab2)
442 0 : return NULL;
443 :
444 296 : isl_assert(tab1->mat->ctx, tab1->M == tab2->M, return NULL);
445 296 : isl_assert(tab1->mat->ctx, tab1->rational == tab2->rational, return NULL);
446 296 : isl_assert(tab1->mat->ctx, tab1->cone == tab2->cone, return NULL);
447 296 : isl_assert(tab1->mat->ctx, !tab1->row_sign, return NULL);
448 296 : isl_assert(tab1->mat->ctx, !tab2->row_sign, return NULL);
449 296 : isl_assert(tab1->mat->ctx, tab1->n_param == 0, return NULL);
450 296 : isl_assert(tab1->mat->ctx, tab2->n_param == 0, return NULL);
451 296 : isl_assert(tab1->mat->ctx, tab1->n_div == 0, return NULL);
452 296 : isl_assert(tab1->mat->ctx, tab2->n_div == 0, return NULL);
453 :
454 296 : off = 2 + tab1->M;
455 296 : r1 = tab1->n_redundant;
456 296 : r2 = tab2->n_redundant;
457 296 : d1 = tab1->n_dead;
458 296 : d2 = tab2->n_dead;
459 296 : prod = isl_calloc_type(tab1->mat->ctx, struct isl_tab);
460 296 : if (!prod)
461 0 : return NULL;
462 296 : prod->mat = tab_mat_product(tab1->mat, tab2->mat,
463 : tab1->n_row, tab2->n_row,
464 : tab1->n_col, tab2->n_col, off, r1, r2, d1, d2);
465 296 : if (!prod->mat)
466 0 : goto error;
467 296 : prod->var = isl_alloc_array(tab1->mat->ctx, struct isl_tab_var,
468 : tab1->max_var + tab2->max_var);
469 296 : if ((tab1->max_var + tab2->max_var) && !prod->var)
470 0 : goto error;
471 1679 : for (i = 0; i < tab1->n_var; ++i) {
472 1383 : prod->var[i] = tab1->var[i];
473 1383 : update_index1(&prod->var[i], r1, r2, d1, d2);
474 : }
475 1679 : for (i = 0; i < tab2->n_var; ++i) {
476 1383 : prod->var[tab1->n_var + i] = tab2->var[i];
477 1383 : update_index2(&prod->var[tab1->n_var + i],
478 : tab1->n_row, tab1->n_col,
479 : r1, r2, d1, d2);
480 : }
481 296 : prod->con = isl_alloc_array(tab1->mat->ctx, struct isl_tab_var,
482 : tab1->max_con + tab2->max_con);
483 296 : if ((tab1->max_con + tab2->max_con) && !prod->con)
484 0 : goto error;
485 3509 : for (i = 0; i < tab1->n_con; ++i) {
486 3213 : prod->con[i] = tab1->con[i];
487 3213 : update_index1(&prod->con[i], r1, r2, d1, d2);
488 : }
489 3509 : for (i = 0; i < tab2->n_con; ++i) {
490 3213 : prod->con[tab1->n_con + i] = tab2->con[i];
491 3213 : update_index2(&prod->con[tab1->n_con + i],
492 : tab1->n_row, tab1->n_col,
493 : r1, r2, d1, d2);
494 : }
495 296 : prod->col_var = isl_alloc_array(tab1->mat->ctx, int,
496 : tab1->n_col + tab2->n_col);
497 296 : if ((tab1->n_col + tab2->n_col) && !prod->col_var)
498 0 : goto error;
499 1185 : for (i = 0; i < tab1->n_col; ++i) {
500 889 : int pos = i < d1 ? i : i + d2;
501 889 : prod->col_var[pos] = tab1->col_var[i];
502 : }
503 1185 : for (i = 0; i < tab2->n_col; ++i) {
504 889 : int pos = i < d2 ? d1 + i : tab1->n_col + i;
505 889 : int t = tab2->col_var[i];
506 889 : if (t >= 0)
507 0 : t += tab1->n_var;
508 : else
509 889 : t -= tab1->n_con;
510 889 : prod->col_var[pos] = t;
511 : }
512 296 : prod->row_var = isl_alloc_array(tab1->mat->ctx, int,
513 : tab1->mat->n_row + tab2->mat->n_row);
514 296 : if ((tab1->mat->n_row + tab2->mat->n_row) && !prod->row_var)
515 0 : goto error;
516 3509 : for (i = 0; i < tab1->n_row; ++i) {
517 3213 : int pos = i < r1 ? i : i + r2;
518 3213 : prod->row_var[pos] = tab1->row_var[i];
519 : }
520 3509 : for (i = 0; i < tab2->n_row; ++i) {
521 3213 : int pos = i < r2 ? r1 + i : tab1->n_row + i;
522 3213 : int t = tab2->row_var[i];
523 3213 : if (t >= 0)
524 1383 : t += tab1->n_var;
525 : else
526 1830 : t -= tab1->n_con;
527 3213 : prod->row_var[pos] = t;
528 : }
529 296 : prod->samples = NULL;
530 296 : prod->sample_index = NULL;
531 296 : prod->n_row = tab1->n_row + tab2->n_row;
532 296 : prod->n_con = tab1->n_con + tab2->n_con;
533 296 : prod->n_eq = 0;
534 296 : prod->max_con = tab1->max_con + tab2->max_con;
535 296 : prod->n_col = tab1->n_col + tab2->n_col;
536 296 : prod->n_var = tab1->n_var + tab2->n_var;
537 296 : prod->max_var = tab1->max_var + tab2->max_var;
538 296 : prod->n_param = 0;
539 296 : prod->n_div = 0;
540 296 : prod->n_dead = tab1->n_dead + tab2->n_dead;
541 296 : prod->n_redundant = tab1->n_redundant + tab2->n_redundant;
542 296 : prod->rational = tab1->rational;
543 296 : prod->empty = tab1->empty || tab2->empty;
544 296 : prod->strict_redundant = tab1->strict_redundant || tab2->strict_redundant;
545 296 : prod->need_undo = 0;
546 296 : prod->in_undo = 0;
547 296 : prod->M = tab1->M;
548 296 : prod->cone = tab1->cone;
549 296 : prod->bottom.type = isl_tab_undo_bottom;
550 296 : prod->bottom.next = NULL;
551 296 : prod->top = &prod->bottom;
552 :
553 296 : prod->n_zero = 0;
554 296 : prod->n_unbounded = 0;
555 296 : prod->basis = NULL;
556 :
557 296 : return prod;
558 : error:
559 0 : isl_tab_free(prod);
560 0 : return NULL;
561 : }
562 :
563 >77604*10^7 : static struct isl_tab_var *var_from_index(struct isl_tab *tab, int i)
564 : {
565 >77604*10^7 : if (i >= 0)
566 >17726*10^7 : return &tab->var[i];
567 : else
568 >59878*10^7 : return &tab->con[~i];
569 : }
570 :
571 >48445*10^7 : struct isl_tab_var *isl_tab_var_from_row(struct isl_tab *tab, int i)
572 : {
573 >48445*10^7 : return var_from_index(tab, tab->row_var[i]);
574 : }
575 :
576 >18953*10^7 : static struct isl_tab_var *var_from_col(struct isl_tab *tab, int i)
577 : {
578 >18953*10^7 : return var_from_index(tab, tab->col_var[i]);
579 : }
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 11003214679 : static int max_is_manifestly_unbounded(struct isl_tab *tab,
586 : struct isl_tab_var *var)
587 : {
588 : int i;
589 11003214679 : unsigned off = 2 + tab->M;
590 :
591 11003214679 : if (var->is_row)
592 5892966170 : return 0;
593 31587480388 : for (i = tab->n_redundant; i < tab->n_row; ++i) {
594 27698477632 : if (!isl_int_is_neg(tab->mat->row[i][off + var->index]))
595 20972611748 : continue;
596 6725865884 : if (isl_tab_var_from_row(tab, i)->is_nonneg)
597 1221245753 : return 0;
598 : }
599 3889002756 : return 1;
600 : }
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 3228118421 : static int min_is_manifestly_unbounded(struct isl_tab *tab,
607 : struct isl_tab_var *var)
608 : {
609 : int i;
610 3228118421 : unsigned off = 2 + tab->M;
611 :
612 3228118421 : if (var->is_row)
613 384587040 : return 0;
614 13077192055 : for (i = tab->n_redundant; i < tab->n_row; ++i) {
615 12509778105 : if (!isl_int_is_pos(tab->mat->row[i][off + var->index]))
616 7833309506 : continue;
617 4676468599 : if (isl_tab_var_from_row(tab, i)->is_nonneg)
618 2276117431 : return 0;
619 : }
620 567413950 : return 1;
621 : }
622 :
623 16168532245 : static int row_cmp(struct isl_tab *tab, int r1, int r2, int c, isl_int *t)
624 : {
625 16168532245 : unsigned off = 2 + tab->M;
626 :
627 16168532245 : if (tab->M) {
628 : 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 : }
635 16168532245 : isl_int_mul(*t, tab->mat->row[r1][1], tab->mat->row[r2][off + c]);
636 16168532245 : isl_int_submul(*t, tab->mat->row[r2][1], tab->mat->row[r1][off + c]);
637 16168532245 : return isl_int_sgn(*t);
638 : }
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 19954532727 : static int pivot_row(struct isl_tab *tab,
661 : struct isl_tab_var *var, int sgn, int c)
662 : {
663 : int j, r, tsgn;
664 : isl_int t;
665 19954532727 : unsigned off = 2 + tab->M;
666 :
667 19954532727 : isl_int_init(t);
668 19954532727 : r = -1;
669 >16001*10^7 : for (j = tab->n_redundant; j < tab->n_row; ++j) {
670 >14005*10^7 : if (var && j == var->index)
671 17584994312 : continue;
672 >12247*10^7 : if (!isl_tab_var_from_row(tab, j)->is_nonneg)
673 47321186080 : continue;
674 75152630137 : if (sgn * isl_int_sgn(tab->mat->row[j][off + c]) >= 0)
675 45839614514 : continue;
676 29313015623 : if (r < 0) {
677 13144483378 : r = j;
678 13144483378 : continue;
679 : }
680 16168532245 : tsgn = sgn * row_cmp(tab, r, j, c, &t);
681 27689882407 : if (tsgn < 0 || (tsgn == 0 &&
682 11521350162 : tab->row_var[j] < tab->row_var[r]))
683 12850850245 : r = j;
684 : }
685 19954532727 : isl_int_clear(t);
686 19954532727 : return r;
687 : }
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 30577554424 : 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 : {
709 : int j, r, c;
710 : isl_int *tr;
711 :
712 30577554424 : *row = *col = -1;
713 :
714 30577554424 : isl_assert(tab->mat->ctx, var->is_row, return);
715 30577554424 : tr = tab->mat->row[var->index] + 2 + tab->M;
716 :
717 30577554424 : c = -1;
718 >16893*10^7 : for (j = tab->n_dead; j < tab->n_col; ++j) {
719 >13835*10^7 : if (isl_int_is_zero(tr[j]))
720 59545391521 : continue;
721 >13863*10^7 : if (isl_int_sgn(tr[j]) != sgn &&
722 59826031723 : var_from_col(tab, j)->is_nonneg)
723 48405180435 : continue;
724 30408714391 : if (c < 0 || tab->col_var[j] < tab->col_var[c])
725 19079341208 : c = j;
726 : }
727 30577554424 : if (c < 0)
728 12991370298 : return;
729 :
730 17586184126 : sgn *= isl_int_sgn(tr[c]);
731 17586184126 : r = pivot_row(tab, skip_var, sgn, c);
732 17586184126 : *row = r < 0 ? var->index : r;
733 17586184126 : *col = c;
734 : }
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 >14386*10^7 : int isl_tab_row_is_redundant(struct isl_tab *tab, int row)
743 : {
744 : int i;
745 >14386*10^7 : unsigned off = 2 + tab->M;
746 :
747 >14386*10^7 : if (tab->row_var[row] < 0 && !isl_tab_var_from_row(tab, row)->is_nonneg)
748 134926008 : return 0;
749 :
750 >14372*10^7 : if (isl_int_is_neg(tab->mat->row[row][1]))
751 28048518868 : return 0;
752 >11567*10^7 : if (tab->strict_redundant && isl_int_is_zero(tab->mat->row[row][1]))
753 0 : return 0;
754 >11567*10^7 : if (tab->M && isl_int_is_neg(tab->mat->row[row][2]))
755 0 : return 0;
756 :
757 >30956*10^7 : for (i = tab->n_dead; i < tab->n_col; ++i) {
758 >29463*10^7 : if (isl_int_is_zero(tab->mat->row[row][off + i]))
759 >11549*10^7 : continue;
760 >17913*10^7 : if (tab->col_var[i] >= 0)
761 49657920574 : return 0;
762 >12947*10^7 : if (isl_int_is_neg(tab->mat->row[row][off + i]))
763 51065679740 : return 0;
764 78414119565 : if (!var_from_col(tab, i)->is_nonneg)
765 26023260 : return 0;
766 : }
767 14926939047 : return 1;
768 : }
769 :
770 23139250734 : static void swap_rows(struct isl_tab *tab, int row1, int row2)
771 : {
772 : int t;
773 : enum isl_tab_row_sign s;
774 :
775 23139250734 : t = tab->row_var[row1];
776 23139250734 : tab->row_var[row1] = tab->row_var[row2];
777 23139250734 : tab->row_var[row2] = t;
778 23139250734 : isl_tab_var_from_row(tab, row1)->index = row1;
779 23139250734 : isl_tab_var_from_row(tab, row2)->index = row2;
780 23139250734 : tab->mat = isl_mat_swap_rows(tab->mat, row1, row2);
781 :
782 23139250734 : if (!tab->row_sign)
783 23139250734 : return;
784 0 : s = tab->row_sign[row1];
785 0 : tab->row_sign[row1] = tab->row_sign[row2];
786 0 : tab->row_sign[row2] = s;
787 : }
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 >18043*10^7 : static isl_stat push_union(struct isl_tab *tab,
800 : enum isl_tab_undo_type type, union isl_tab_undo_val u)
801 : {
802 : struct isl_tab_undo *undo;
803 :
804 >18043*10^7 : if (!tab)
805 0 : return isl_stat_error;
806 >18043*10^7 : if (!tab->need_undo)
807 54754258721 : return isl_stat_ok;
808 :
809 >12568*10^7 : undo = isl_alloc_type(tab->mat->ctx, struct isl_tab_undo);
810 >12568*10^7 : if (!undo)
811 0 : goto error;
812 >12568*10^7 : undo->type = type;
813 >12568*10^7 : undo->u = u;
814 >12568*10^7 : undo->next = tab->top;
815 >12568*10^7 : tab->top = undo;
816 :
817 >12568*10^7 : return isl_stat_ok;
818 : error:
819 0 : free_undo(tab);
820 0 : tab->top = NULL;
821 0 : return isl_stat_error;
822 : }
823 :
824 >15735*10^7 : isl_stat isl_tab_push_var(struct isl_tab *tab,
825 : enum isl_tab_undo_type type, struct isl_tab_var *var)
826 : {
827 : union isl_tab_undo_val u;
828 >15735*10^7 : if (var->is_row)
829 >12408*10^7 : u.var_index = tab->row_var[var->index];
830 : else
831 33277628123 : u.var_index = tab->col_var[var->index];
832 >15735*10^7 : return push_union(tab, type, u);
833 : }
834 :
835 23078717745 : isl_stat isl_tab_push(struct isl_tab *tab, enum isl_tab_undo_type type)
836 : {
837 23078717745 : union isl_tab_undo_val u = { 0 };
838 23078717745 : return push_union(tab, type, u);
839 : }
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 0 : isl_stat isl_tab_push_basis(struct isl_tab *tab)
845 : {
846 : int i;
847 : union isl_tab_undo_val u;
848 :
849 0 : u.col_var = isl_alloc_array(tab->mat->ctx, int, tab->n_col);
850 0 : if (tab->n_col && !u.col_var)
851 0 : return isl_stat_error;
852 0 : for (i = 0; i < tab->n_col; ++i)
853 0 : u.col_var[i] = tab->col_var[i];
854 0 : return push_union(tab, isl_tab_undo_saved_basis, u);
855 : }
856 :
857 0 : isl_stat isl_tab_push_callback(struct isl_tab *tab,
858 : struct isl_tab_callback *callback)
859 : {
860 : union isl_tab_undo_val u;
861 0 : u.callback = callback;
862 0 : return push_union(tab, isl_tab_undo_callback, u);
863 : }
864 :
865 0 : struct isl_tab *isl_tab_init_samples(struct isl_tab *tab)
866 : {
867 0 : if (!tab)
868 0 : return NULL;
869 :
870 0 : tab->n_sample = 0;
871 0 : tab->n_outside = 0;
872 0 : tab->samples = isl_mat_alloc(tab->mat->ctx, 1, 1 + tab->n_var);
873 0 : if (!tab->samples)
874 0 : goto error;
875 0 : tab->sample_index = isl_alloc_array(tab->mat->ctx, int, 1);
876 0 : if (!tab->sample_index)
877 0 : goto error;
878 0 : return tab;
879 : error:
880 0 : isl_tab_free(tab);
881 0 : return NULL;
882 : }
883 :
884 0 : int isl_tab_add_sample(struct isl_tab *tab, __isl_take isl_vec *sample)
885 : {
886 0 : if (!tab || !sample)
887 : goto error;
888 :
889 0 : if (tab->n_sample + 1 > tab->samples->n_row) {
890 0 : int *t = isl_realloc_array(tab->mat->ctx,
891 : tab->sample_index, int, tab->n_sample + 1);
892 0 : if (!t)
893 0 : goto error;
894 0 : tab->sample_index = t;
895 : }
896 :
897 0 : tab->samples = isl_mat_extend(tab->samples,
898 0 : tab->n_sample + 1, tab->samples->n_col);
899 0 : if (!tab->samples)
900 0 : goto error;
901 :
902 0 : isl_seq_cpy(tab->samples->row[tab->n_sample], sample->el, sample->size);
903 0 : isl_vec_free(sample);
904 0 : tab->sample_index[tab->n_sample] = tab->n_sample;
905 0 : tab->n_sample++;
906 :
907 0 : return 0;
908 : error:
909 0 : isl_vec_free(sample);
910 0 : return -1;
911 : }
912 :
913 0 : struct isl_tab *isl_tab_drop_sample(struct isl_tab *tab, int s)
914 : {
915 0 : if (s != tab->n_outside) {
916 0 : int t = tab->sample_index[tab->n_outside];
917 0 : tab->sample_index[tab->n_outside] = tab->sample_index[s];
918 0 : tab->sample_index[s] = t;
919 0 : isl_mat_swap_rows(tab->samples, tab->n_outside, s);
920 : }
921 0 : tab->n_outside++;
922 0 : if (isl_tab_push(tab, isl_tab_undo_drop_sample) < 0) {
923 0 : isl_tab_free(tab);
924 0 : return NULL;
925 : }
926 :
927 0 : return tab;
928 : }
929 :
930 : /* Record the current number of samples so that we can remove newer
931 : * samples during a rollback.
932 : */
933 0 : isl_stat isl_tab_save_samples(struct isl_tab *tab)
934 : {
935 : union isl_tab_undo_val u;
936 :
937 0 : if (!tab)
938 0 : return isl_stat_error;
939 :
940 0 : u.n = tab->n_sample;
941 0 : return push_union(tab, isl_tab_undo_saved_samples, u);
942 : }
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 19886097797 : int isl_tab_mark_redundant(struct isl_tab *tab, int row)
957 : {
958 19886097797 : struct isl_tab_var *var = isl_tab_var_from_row(tab, row);
959 19886097797 : var->is_redundant = 1;
960 19886097797 : isl_assert(tab->mat->ctx, row >= tab->n_redundant, return -1);
961 19886097797 : if (tab->preserve || tab->need_undo || tab->row_var[row] >= 0) {
962 13289054163 : if (tab->row_var[row] >= 0 && !var->is_nonneg) {
963 5729843691 : var->is_nonneg = 1;
964 5729843691 : if (isl_tab_push_var(tab, isl_tab_undo_nonneg, var) < 0)
965 0 : return -1;
966 : }
967 13289054163 : if (row != tab->n_redundant)
968 11228520054 : swap_rows(tab, row, tab->n_redundant);
969 13289054163 : tab->n_redundant++;
970 13289054163 : return isl_tab_push_var(tab, isl_tab_undo_redundant, var);
971 : } else {
972 6597043634 : if (row != tab->n_row - 1)
973 4913609835 : swap_rows(tab, row, tab->n_row - 1);
974 6597043634 : isl_tab_var_from_row(tab, tab->n_row - 1)->index = -1;
975 6597043634 : tab->n_row--;
976 6597043634 : return 1;
977 : }
978 : }
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 136036 : int isl_tab_mark_rational(struct isl_tab *tab)
986 : {
987 136036 : if (!tab)
988 0 : return -1;
989 136036 : if (!tab->rational && tab->need_undo)
990 136036 : if (isl_tab_push(tab, isl_tab_undo_rational) < 0)
991 0 : return -1;
992 136036 : tab->rational = 1;
993 136036 : return 0;
994 : }
995 :
996 7952065536 : isl_stat isl_tab_mark_empty(struct isl_tab *tab)
997 : {
998 7952065536 : if (!tab)
999 0 : return isl_stat_error;
1000 7952065536 : if (!tab->empty && tab->need_undo)
1001 6088534664 : if (isl_tab_push(tab, isl_tab_undo_empty) < 0)
1002 0 : return isl_stat_error;
1003 7952065536 : tab->empty = 1;
1004 7952065536 : return isl_stat_ok;
1005 : }
1006 :
1007 61076335628 : int isl_tab_freeze_constraint(struct isl_tab *tab, int con)
1008 : {
1009 : struct isl_tab_var *var;
1010 :
1011 61076335628 : if (!tab)
1012 0 : return -1;
1013 :
1014 61076335628 : var = &tab->con[con];
1015 61076335628 : if (var->frozen)
1016 0 : return 0;
1017 61076335628 : if (var->index < 0)
1018 1176353041 : return 0;
1019 59899982587 : var->frozen = 1;
1020 :
1021 59899982587 : if (tab->need_undo)
1022 59899982587 : return isl_tab_push_var(tab, isl_tab_undo_freeze, var);
1023 :
1024 0 : return 0;
1025 : }
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 20210094823 : static void update_row_sign(struct isl_tab *tab, int row, int col, int row_sgn)
1045 : {
1046 : int i;
1047 20210094823 : struct isl_mat *mat = tab->mat;
1048 20210094823 : unsigned off = 2 + tab->M;
1049 :
1050 20210094823 : if (!tab->row_sign)
1051 20210094823 : return;
1052 :
1053 0 : if (tab->row_sign[row] == 0)
1054 0 : return;
1055 0 : isl_assert(mat->ctx, row_sgn > 0, return);
1056 0 : isl_assert(mat->ctx, tab->row_sign[row] == isl_tab_row_neg, return);
1057 0 : tab->row_sign[row] = isl_tab_row_pos;
1058 0 : for (i = 0; i < tab->n_row; ++i) {
1059 : int s;
1060 0 : if (i == row)
1061 0 : continue;
1062 0 : s = isl_int_sgn(mat->row[i][off + col]);
1063 0 : if (!s)
1064 0 : continue;
1065 0 : if (!tab->row_sign[i])
1066 0 : continue;
1067 0 : if (s < 0 && tab->row_sign[i] == isl_tab_row_neg)
1068 0 : continue;
1069 0 : if (s > 0 && tab->row_sign[i] == isl_tab_row_pos)
1070 0 : continue;
1071 0 : tab->row_sign[i] = isl_tab_row_unknown;
1072 : }
1073 : }
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 20210094823 : int isl_tab_pivot(struct isl_tab *tab, int row, int col)
1126 : {
1127 : int i, j;
1128 : int sgn;
1129 : int t;
1130 : isl_ctx *ctx;
1131 20210094823 : struct isl_mat *mat = tab->mat;
1132 : struct isl_tab_var *var;
1133 20210094823 : unsigned off = 2 + tab->M;
1134 :
1135 20210094823 : ctx = isl_tab_get_ctx(tab);
1136 20210094823 : if (isl_ctx_next_operation(ctx) < 0)
1137 0 : return -1;
1138 :
1139 20210094823 : isl_int_swap(mat->row[row][0], mat->row[row][off + col]);
1140 20210094823 : sgn = isl_int_sgn(mat->row[row][0]);
1141 20210094823 : if (sgn < 0) {
1142 10525530446 : isl_int_neg(mat->row[row][0], mat->row[row][0]);
1143 10525530446 : isl_int_neg(mat->row[row][off + col], mat->row[row][off + col]);
1144 : } else
1145 66082324105 : for (j = 0; j < off - 1 + tab->n_col; ++j) {
1146 56397759728 : if (j == off - 1 + col)
1147 9684564377 : continue;
1148 46713195351 : isl_int_neg(mat->row[row][1 + j], mat->row[row][1 + j]);
1149 : }
1150 20210094823 : if (!isl_int_is_one(mat->row[row][0]))
1151 5493001287 : isl_seq_normalize(mat->ctx, mat->row[row], off + tab->n_col);
1152 >18453*10^7 : for (i = 0; i < tab->n_row; ++i) {
1153 >16432*10^7 : if (i == row)
1154 20210094823 : continue;
1155 >14411*10^7 : if (isl_int_is_zero(mat->row[i][off + col]))
1156 53449802035 : continue;
1157 90668621897 : isl_int_mul(mat->row[i][0], mat->row[i][0], mat->row[row][0]);
1158 >66153*10^7 : for (j = 0; j < off - 1 + tab->n_col; ++j) {
1159 >57086*10^7 : if (j == off - 1 + col)
1160 90668621897 : continue;
1161 >48019*10^7 : isl_int_mul(mat->row[i][1 + j],
1162 : mat->row[i][1 + j], mat->row[row][0]);
1163 >48019*10^7 : isl_int_addmul(mat->row[i][1 + j],
1164 : mat->row[i][off + col], mat->row[row][1 + j]);
1165 : }
1166 90668621897 : isl_int_mul(mat->row[i][off + col],
1167 : mat->row[i][off + col], mat->row[row][off + col]);
1168 90668621897 : if (!isl_int_is_one(mat->row[i][0]))
1169 46933256519 : isl_seq_normalize(mat->ctx, mat->row[i], off + tab->n_col);
1170 : }
1171 20210094823 : t = tab->row_var[row];
1172 20210094823 : tab->row_var[row] = tab->col_var[col];
1173 20210094823 : tab->col_var[col] = t;
1174 20210094823 : var = isl_tab_var_from_row(tab, row);
1175 20210094823 : var->is_row = 1;
1176 20210094823 : var->index = row;
1177 20210094823 : var = var_from_col(tab, col);
1178 20210094823 : var->is_row = 0;
1179 20210094823 : var->index = col;
1180 20210094823 : update_row_sign(tab, row, col, sgn);
1181 20210094823 : if (tab->in_undo)
1182 1668417240 : return 0;
1183 >14633*10^7 : for (i = tab->n_redundant; i < tab->n_row; ++i) {
1184 >12779*10^7 : if (isl_int_is_zero(mat->row[i][off + col]))
1185 34568756886 : continue;
1186 >17853*10^7 : if (!isl_tab_var_from_row(tab, i)->frozen &&
1187 85309752098 : isl_tab_row_is_redundant(tab, i)) {
1188 8244494705 : int redo = isl_tab_mark_redundant(tab, i);
1189 8244494705 : if (redo < 0)
1190 0 : return -1;
1191 8244494705 : if (redo)
1192 1611146812 : --i;
1193 : }
1194 : }
1195 18541677583 : return 0;
1196 : }
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 8069566805 : static int to_row(struct isl_tab *tab, struct isl_tab_var *var, int sign)
1206 : {
1207 : int r;
1208 8069566805 : unsigned off = 2 + tab->M;
1209 :
1210 8069566805 : if (var->is_row)
1211 5895565731 : return 0;
1212 :
1213 2174001074 : if (sign == 0) {
1214 609311296 : for (r = tab->n_redundant; r < tab->n_row; ++r)
1215 609311296 : if (!isl_int_is_zero(tab->mat->row[r][off+var->index]))
1216 : break;
1217 191261687 : isl_assert(tab->mat->ctx, r < tab->n_row, return -1);
1218 : } else {
1219 1982739387 : r = pivot_row(tab, NULL, sign, var->index);
1220 1982739387 : isl_assert(tab->mat->ctx, r >= 0, return -1);
1221 : }
1222 :
1223 2174001074 : return isl_tab_pivot(tab, r, var->index);
1224 : }
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 0 : static void check_table(struct isl_tab *tab)
1232 : {
1233 : int i;
1234 :
1235 0 : if (tab->empty)
1236 0 : return;
1237 0 : for (i = tab->n_redundant; i < tab->n_row; ++i) {
1238 : 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 : !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 : }
1248 0 : isl_assert(tab->mat->ctx, !isl_int_is_neg(tab->mat->row[i][1]),
1249 : abort());
1250 : }
1251 : }
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 7141854055 : static int sign_of_max(struct isl_tab *tab, struct isl_tab_var *var)
1265 : {
1266 : int row, col;
1267 :
1268 7141854055 : if (max_is_manifestly_unbounded(tab, var))
1269 825206185 : return 1;
1270 6316647870 : if (to_row(tab, var, 1) < 0)
1271 0 : return -2;
1272 19021485911 : while (!isl_int_is_pos(tab->mat->row[var->index][1])) {
1273 12620530179 : find_pivot(tab, var, var, 1, &row, &col);
1274 12620530179 : if (row == -1)
1275 4946605662 : return isl_int_sgn(tab->mat->row[var->index][1]);
1276 7673924517 : if (isl_tab_pivot(tab, row, col) < 0)
1277 0 : return -2;
1278 7673924517 : if (!var->is_row) /* manifestly unbounded */
1279 1285734346 : return 1;
1280 : }
1281 84307862 : return 1;
1282 : }
1283 :
1284 0 : int isl_tab_sign_of_max(struct isl_tab *tab, int con)
1285 : {
1286 : struct isl_tab_var *var;
1287 :
1288 0 : if (!tab)
1289 0 : return -2;
1290 :
1291 0 : var = &tab->con[con];
1292 0 : isl_assert(tab->mat->ctx, !var->is_redundant, return -2);
1293 0 : isl_assert(tab->mat->ctx, !var->is_zero, return -2);
1294 :
1295 0 : return sign_of_max(tab, var);
1296 : }
1297 :
1298 44597730261 : static int row_is_neg(struct isl_tab *tab, int row)
1299 : {
1300 44597730261 : if (!tab->M)
1301 44597730261 : 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 : }
1308 :
1309 35360563778 : static int row_sgn(struct isl_tab *tab, int row)
1310 : {
1311 35360563778 : if (!tab->M)
1312 35360563778 : 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 : else
1316 0 : return isl_int_sgn(tab->mat->row[row][1]);
1317 : }
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 40601364546 : static int restore_row(struct isl_tab *tab, struct isl_tab_var *var)
1325 : {
1326 : int row, col;
1327 :
1328 85199094807 : while (row_is_neg(tab, var->index)) {
1329 17192412177 : find_pivot(tab, var, var, 1, &row, &col);
1330 17192412177 : if (row == -1)
1331 7955245694 : break;
1332 9237166483 : if (isl_tab_pivot(tab, row, col) < 0)
1333 0 : return -2;
1334 9237166483 : if (!var->is_row) /* manifestly unbounded */
1335 5240800768 : return 1;
1336 : }
1337 35360563778 : return row_sgn(tab, var->index);
1338 : }
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 38629715 : static int at_least_zero(struct isl_tab *tab, struct isl_tab_var *var)
1346 : {
1347 : int row, col;
1348 :
1349 88037026 : while (isl_int_is_neg(tab->mat->row[var->index][1])) {
1350 46980254 : find_pivot(tab, var, var, 1, &row, &col);
1351 46980254 : if (row == -1)
1352 26721360 : break;
1353 20258894 : if (row == var->index) /* manifestly unbounded */
1354 9481298 : return 1;
1355 10777596 : if (isl_tab_pivot(tab, row, col) < 0)
1356 0 : return -1;
1357 : }
1358 29148417 : return !isl_int_is_neg(tab->mat->row[var->index][1]);
1359 : }
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 61013 : static int sign_of_min(struct isl_tab *tab, struct isl_tab_var *var)
1379 : {
1380 : int row, col;
1381 61013 : struct isl_tab_var *pivot_var = NULL;
1382 :
1383 61013 : if (min_is_manifestly_unbounded(tab, var))
1384 0 : return -1;
1385 61013 : if (!var->is_row) {
1386 0 : col = var->index;
1387 0 : row = pivot_row(tab, NULL, -1, col);
1388 0 : pivot_var = var_from_col(tab, col);
1389 0 : if (isl_tab_pivot(tab, row, col) < 0)
1390 0 : return -2;
1391 0 : if (var->is_redundant)
1392 0 : return 0;
1393 0 : if (isl_int_is_neg(tab->mat->row[var->index][1])) {
1394 0 : if (var->is_nonneg) {
1395 0 : if (!pivot_var->is_redundant &&
1396 0 : pivot_var->index == row) {
1397 0 : if (isl_tab_pivot(tab, row, col) < 0)
1398 0 : return -2;
1399 : } else
1400 0 : if (restore_row(tab, var) < -1)
1401 0 : return -2;
1402 : }
1403 0 : return -1;
1404 : }
1405 : }
1406 61013 : if (var->is_redundant)
1407 0 : return 0;
1408 131809 : while (!isl_int_is_neg(tab->mat->row[var->index][1])) {
1409 65416 : find_pivot(tab, var, var, -1, &row, &col);
1410 65416 : if (row == var->index)
1411 49425 : return -1;
1412 15991 : if (row == -1)
1413 6208 : return isl_int_sgn(tab->mat->row[var->index][1]);
1414 9783 : pivot_var = var_from_col(tab, col);
1415 9783 : if (isl_tab_pivot(tab, row, col) < 0)
1416 0 : return -2;
1417 9783 : if (var->is_redundant)
1418 0 : return 0;
1419 : }
1420 5380 : if (pivot_var && var->is_nonneg) {
1421 : /* pivot back to non-negative value */
1422 0 : if (!pivot_var->is_redundant && pivot_var->index == row) {
1423 0 : if (isl_tab_pivot(tab, row, col) < 0)
1424 0 : return -2;
1425 : } else
1426 0 : if (restore_row(tab, var) < -1)
1427 0 : return -2;
1428 : }
1429 5380 : return -1;
1430 : }
1431 :
1432 604663692 : static int row_at_most_neg_one(struct isl_tab *tab, int row)
1433 : {
1434 604663692 : 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 : }
1440 973229170 : return isl_int_is_neg(tab->mat->row[row][1]) &&
1441 368565478 : isl_int_abs_ge(tab->mat->row[row][1],
1442 : tab->mat->row[row][0]);
1443 : }
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 767535680 : int isl_tab_min_at_most_neg_one(struct isl_tab *tab, struct isl_tab_var *var)
1455 : {
1456 : int row, col;
1457 : struct isl_tab_var *pivot_var;
1458 :
1459 767535680 : if (min_is_manifestly_unbounded(tab, var))
1460 0 : return 1;
1461 767535680 : if (!var->is_row) {
1462 385609214 : col = var->index;
1463 385609214 : row = pivot_row(tab, NULL, -1, col);
1464 385609214 : pivot_var = var_from_col(tab, col);
1465 385609214 : if (isl_tab_pivot(tab, row, col) < 0)
1466 0 : return -1;
1467 385609214 : if (var->is_redundant)
1468 37165854 : return 0;
1469 348443360 : if (row_at_most_neg_one(tab, var->index)) {
1470 221718313 : if (var->is_nonneg) {
1471 443436626 : if (!pivot_var->is_redundant &&
1472 221718313 : pivot_var->index == row) {
1473 439508256 : if (isl_tab_pivot(tab, row, col) < 0)
1474 0 : return -1;
1475 : } else
1476 1964185 : if (restore_row(tab, var) < -1)
1477 0 : return -1;
1478 : }
1479 221718313 : return 1;
1480 : }
1481 : }
1482 508651513 : if (var->is_redundant)
1483 0 : return 0;
1484 : do {
1485 686312043 : find_pivot(tab, var, var, -1, &row, &col);
1486 686312043 : if (row == var->index) {
1487 249230552 : if (var->is_nonneg && restore_row(tab, var) < -1)
1488 0 : return -1;
1489 249230552 : return 1;
1490 : }
1491 437081491 : if (row == -1)
1492 61523448 : return 0;
1493 375558043 : pivot_var = var_from_col(tab, col);
1494 375558043 : if (isl_tab_pivot(tab, row, col) < 0)
1495 0 : return -1;
1496 375558043 : if (var->is_redundant)
1497 119337711 : return 0;
1498 256220332 : } while (!row_at_most_neg_one(tab, var->index));
1499 78559802 : if (var->is_nonneg) {
1500 : /* pivot back to non-negative value */
1501 74906909 : if (!pivot_var->is_redundant && pivot_var->index == row)
1502 74100872 : if (isl_tab_pivot(tab, row, col) < 0)
1503 0 : return -1;
1504 74906909 : if (restore_row(tab, var) < -1)
1505 0 : return -1;
1506 : }
1507 78559802 : return 1;
1508 : }
1509 :
1510 : /* Return 1 if "var" can attain values >= 1.
1511 : * Return 0 otherwise.
1512 : */
1513 2195248381 : static int at_least_one(struct isl_tab *tab, struct isl_tab_var *var)
1514 : {
1515 : int row, col;
1516 : isl_int *r;
1517 :
1518 2195248381 : if (max_is_manifestly_unbounded(tab, var))
1519 2110940531 : return 1;
1520 84307850 : if (to_row(tab, var, 1) < 0)
1521 0 : return -1;
1522 84307850 : r = tab->mat->row[var->index];
1523 172535430 : while (isl_int_lt(r[1], r[0])) {
1524 26273186 : find_pivot(tab, var, var, 1, &row, &col);
1525 26273186 : if (row == -1)
1526 14590 : return isl_int_ge(r[1], r[0]);
1527 26258596 : if (row == var->index) /* manifestly unbounded */
1528 22338866 : return 1;
1529 3919730 : if (isl_tab_pivot(tab, row, col) < 0)
1530 0 : return -1;
1531 : }
1532 61954394 : return 1;
1533 : }
1534 :
1535 3692186458 : static void swap_cols(struct isl_tab *tab, int col1, int col2)
1536 : {
1537 : int t;
1538 3692186458 : unsigned off = 2 + tab->M;
1539 3692186458 : t = tab->col_var[col1];
1540 3692186458 : tab->col_var[col1] = tab->col_var[col2];
1541 3692186458 : tab->col_var[col2] = t;
1542 3692186458 : var_from_col(tab, col1)->index = col1;
1543 3692186458 : var_from_col(tab, col2)->index = col2;
1544 3692186458 : tab->mat = isl_mat_swap_cols(tab->mat, off + col1, off + col2);
1545 3692186458 : }
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 5044730281 : int isl_tab_kill_col(struct isl_tab *tab, int col)
1560 : {
1561 5044730281 : var_from_col(tab, col)->is_zero = 1;
1562 5044730281 : if (tab->need_undo) {
1563 348673 : if (isl_tab_push_var(tab, isl_tab_undo_zero,
1564 : var_from_col(tab, col)) < 0)
1565 0 : return -1;
1566 348673 : if (col != tab->n_dead)
1567 179164 : swap_cols(tab, col, tab->n_dead);
1568 348673 : tab->n_dead++;
1569 348673 : return 0;
1570 : } else {
1571 5044381608 : if (col != tab->n_col - 1)
1572 3692007294 : swap_cols(tab, col, tab->n_col - 1);
1573 5044381608 : var_from_col(tab, tab->n_col - 1)->index = -1;
1574 5044381608 : tab->n_col--;
1575 5044381608 : return 1;
1576 : }
1577 : }
1578 :
1579 21622659302 : static int row_is_manifestly_non_integral(struct isl_tab *tab, int row)
1580 : {
1581 21622659302 : unsigned off = 2 + tab->M;
1582 :
1583 21622659302 : if (tab->M && !isl_int_eq(tab->mat->row[row][2],
1584 : tab->mat->row[row][0]))
1585 0 : return 0;
1586 21622659302 : if (isl_seq_first_non_zero(tab->mat->row[row] + off + tab->n_dead,
1587 21622659302 : tab->n_col - tab->n_dead) != -1)
1588 3345233541 : return 0;
1589 :
1590 18277425761 : return !isl_int_is_divisible_by(tab->mat->row[row][1],
1591 : tab->mat->row[row][0]);
1592 : }
1593 :
1594 : /* For integer tableaus, check if any of the coordinates are stuck
1595 : * at a non-integral value.
1596 : */
1597 4946604089 : static int tab_is_manifestly_empty(struct isl_tab *tab)
1598 : {
1599 : int i;
1600 :
1601 4946604089 : if (tab->empty)
1602 783 : return 1;
1603 4946603306 : if (tab->rational)
1604 136977 : return 0;
1605 :
1606 29352138572 : for (i = 0; i < tab->n_var; ++i) {
1607 24405675673 : if (!tab->var[i].is_row)
1608 2783016371 : continue;
1609 21622659302 : if (row_is_manifestly_non_integral(tab, tab->var[i].index))
1610 3430 : return 1;
1611 : }
1612 :
1613 4946462899 : return 0;
1614 : }
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 4946604089 : static isl_stat close_row(struct isl_tab *tab, struct isl_tab_var *var,
1636 : int temp_var)
1637 : {
1638 : int j;
1639 4946604089 : struct isl_mat *mat = tab->mat;
1640 4946604089 : unsigned off = 2 + tab->M;
1641 :
1642 4946604089 : if (!var->is_nonneg)
1643 0 : isl_die(isl_tab_get_ctx(tab), isl_error_internal,
1644 : "expecting non-negative variable",
1645 : return isl_stat_error);
1646 4946604089 : var->is_zero = 1;
1647 4946604089 : if (!temp_var && tab->need_undo)
1648 0 : if (isl_tab_push_var(tab, isl_tab_undo_zero, var) < 0)
1649 0 : return isl_stat_error;
1650 14688882809 : for (j = tab->n_dead; j < tab->n_col; ++j) {
1651 : int recheck;
1652 9742278720 : if (isl_int_is_zero(mat->row[var->index][off + j]))
1653 4751512761 : continue;
1654 4990765959 : if (isl_int_is_pos(mat->row[var->index][off + j]))
1655 0 : isl_die(isl_tab_get_ctx(tab), isl_error_internal,
1656 : "row cannot have positive coefficients",
1657 : return isl_stat_error);
1658 4990765959 : recheck = isl_tab_kill_col(tab, j);
1659 4990765959 : if (recheck < 0)
1660 0 : return isl_stat_error;
1661 4990765959 : if (recheck)
1662 4990452731 : --j;
1663 : }
1664 4946604089 : if (!temp_var && isl_tab_mark_redundant(tab, var->index) < 0)
1665 0 : return isl_stat_error;
1666 4946604089 : if (tab_is_manifestly_empty(tab) && isl_tab_mark_empty(tab) < 0)
1667 0 : return isl_stat_error;
1668 4946604089 : return isl_stat_ok;
1669 : }
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 39297444976 : int isl_tab_allocate_con(struct isl_tab *tab)
1678 : {
1679 : int r;
1680 :
1681 39297444976 : isl_assert(tab->mat->ctx, tab->n_row < tab->mat->n_row, return -1);
1682 39297444976 : isl_assert(tab->mat->ctx, tab->n_con < tab->max_con, return -1);
1683 :
1684 39297444976 : r = tab->n_con;
1685 39297444976 : tab->con[r].index = tab->n_row;
1686 39297444976 : tab->con[r].is_row = 1;
1687 39297444976 : tab->con[r].is_nonneg = 0;
1688 39297444976 : tab->con[r].is_zero = 0;
1689 39297444976 : tab->con[r].is_redundant = 0;
1690 39297444976 : tab->con[r].frozen = 0;
1691 39297444976 : tab->con[r].negated = 0;
1692 39297444976 : tab->row_var[tab->n_row] = ~r;
1693 :
1694 39297444976 : tab->n_row++;
1695 39297444976 : tab->n_con++;
1696 39297444976 : if (isl_tab_push_var(tab, isl_tab_undo_allocate, &tab->con[r]) < 0)
1697 0 : return -1;
1698 :
1699 39297444976 : return r;
1700 : }
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 0 : static int var_insert_entry(struct isl_tab *tab, int first)
1708 : {
1709 : int i;
1710 :
1711 0 : if (tab->n_var >= tab->max_var)
1712 0 : isl_die(isl_tab_get_ctx(tab), isl_error_internal,
1713 : "not enough room for new variable", return -1);
1714 0 : if (first > tab->n_var)
1715 0 : isl_die(isl_tab_get_ctx(tab), isl_error_internal,
1716 : "invalid initial position", return -1);
1717 :
1718 0 : for (i = tab->n_var - 1; i >= first; --i) {
1719 0 : tab->var[i + 1] = tab->var[i];
1720 0 : if (tab->var[i + 1].is_row)
1721 0 : tab->row_var[tab->var[i + 1].index]++;
1722 : else
1723 0 : tab->col_var[tab->var[i + 1].index]++;
1724 : }
1725 :
1726 0 : tab->n_var++;
1727 :
1728 0 : return 0;
1729 : }
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 0 : static int var_drop_entry(struct isl_tab *tab, int first)
1737 : {
1738 : int i;
1739 :
1740 0 : if (first >= tab->n_var)
1741 0 : isl_die(isl_tab_get_ctx(tab), isl_error_internal,
1742 : "invalid initial position", return -1);
1743 :
1744 0 : tab->n_var--;
1745 :
1746 0 : for (i = first; i < tab->n_var; ++i) {
1747 0 : tab->var[i] = tab->var[i + 1];
1748 0 : if (tab->var[i + 1].is_row)
1749 0 : tab->row_var[tab->var[i].index]--;
1750 : else
1751 0 : tab->col_var[tab->var[i].index]--;
1752 : }
1753 :
1754 0 : return 0;
1755 : }
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 0 : int isl_tab_insert_var(struct isl_tab *tab, int r)
1762 : {
1763 : int i;
1764 0 : unsigned off = 2 + tab->M;
1765 :
1766 0 : isl_assert(tab->mat->ctx, tab->n_col < tab->mat->n_col, return -1);
1767 :
1768 0 : if (var_insert_entry(tab, r) < 0)
1769 0 : return -1;
1770 :
1771 0 : tab->var[r].index = tab->n_col;
1772 0 : tab->var[r].is_row = 0;
1773 0 : tab->var[r].is_nonneg = 0;
1774 0 : tab->var[r].is_zero = 0;
1775 0 : tab->var[r].is_redundant = 0;
1776 0 : tab->var[r].frozen = 0;
1777 0 : tab->var[r].negated = 0;
1778 0 : tab->col_var[tab->n_col] = r;
1779 :
1780 0 : for (i = 0; i < tab->n_row; ++i)
1781 0 : isl_int_set_si(tab->mat->row[i][off + tab->n_col], 0);
1782 :
1783 0 : tab->n_col++;
1784 0 : if (isl_tab_push_var(tab, isl_tab_undo_allocate, &tab->var[r]) < 0)
1785 0 : return -1;
1786 :
1787 0 : return r;
1788 : }
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 0 : int isl_tab_allocate_var(struct isl_tab *tab)
1794 : {
1795 0 : if (!tab)
1796 0 : return -1;
1797 :
1798 0 : return isl_tab_insert_var(tab, tab->n_var);
1799 : }
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 39297416180 : int isl_tab_add_row(struct isl_tab *tab, isl_int *line)
1824 : {
1825 : int i;
1826 : int r;
1827 : isl_int *row;
1828 : isl_int a, b;
1829 39297416180 : unsigned off = 2 + tab->M;
1830 :
1831 39297416180 : r = isl_tab_allocate_con(tab);
1832 39297416180 : if (r < 0)
1833 0 : return -1;
1834 :
1835 39297416180 : isl_int_init(a);
1836 39297416180 : isl_int_init(b);
1837 39297416180 : row = tab->mat->row[tab->con[r].index];
1838 39297416180 : isl_int_set_si(row[0], 1);
1839 39297416180 : isl_int_set(row[1], line[0]);
1840 39297416180 : isl_seq_clr(row + 2, tab->M + tab->n_col);
1841 >24813*10^7 : for (i = 0; i < tab->n_var; ++i) {
1842 >20884*10^7 : if (tab->var[i].is_zero)
1843 0 : continue;
1844 >20884*10^7 : if (tab->var[i].is_row) {
1845 >11360*10^7 : isl_int_lcm(a,
1846 : row[0], tab->mat->row[tab->var[i].index][0]);
1847 >11360*10^7 : isl_int_swap(a, row[0]);
1848 >11360*10^7 : isl_int_divexact(a, row[0], a);
1849 >11360*10^7 : isl_int_divexact(b,
1850 : row[0], tab->mat->row[tab->var[i].index][0]);
1851 >11360*10^7 : isl_int_mul(b, b, line[1 + i]);
1852 >22720*10^7 : isl_seq_combine(row + 1, a, row + 1,
1853 >11360*10^7 : b, tab->mat->row[tab->var[i].index] + 1,
1854 >11360*10^7 : 1 + tab->M + tab->n_col);
1855 : } else
1856 95236931886 : isl_int_addmul(row[off + tab->var[i].index],
1857 : line[1 + i], row[0]);
1858 >20884*10^7 : if (tab->M && i >= tab->n_param && i < tab->n_var - tab->n_div)
1859 0 : isl_int_submul(row[2], line[1 + i], row[0]);
1860 : }
1861 39297416180 : isl_seq_normalize(tab->mat->ctx, row, off + tab->n_col);
1862 39297416180 : isl_int_clear(a);
1863 39297416180 : isl_int_clear(b);
1864 :
1865 39297416180 : if (tab->row_sign)
1866 0 : tab->row_sign[tab->con[r].index] = isl_tab_row_unknown;
1867 :
1868 39297416180 : return r;
1869 : }
1870 :
1871 16953225374 : static isl_stat drop_row(struct isl_tab *tab, int row)
1872 : {
1873 16953225374 : isl_assert(tab->mat->ctx, ~tab->row_var[row] == tab->n_con - 1,
1874 : return isl_stat_error);
1875 16953225374 : if (row != tab->n_row - 1)
1876 6997120845 : swap_rows(tab, row, tab->n_row - 1);
1877 16953225374 : tab->n_row--;
1878 16953225374 : tab->n_con--;
1879 16953225374 : return isl_stat_ok;
1880 : }
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 0 : static isl_stat drop_col(struct isl_tab *tab, int col)
1889 : {
1890 : int var;
1891 :
1892 0 : var = tab->col_var[col];
1893 0 : if (col != tab->n_col - 1)
1894 0 : swap_cols(tab, col, tab->n_col - 1);
1895 0 : tab->n_col--;
1896 0 : if (var_drop_entry(tab, var) < 0)
1897 0 : return isl_stat_error;
1898 0 : return isl_stat_ok;
1899 : }
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 39141870197 : isl_stat isl_tab_add_ineq(struct isl_tab *tab, isl_int *ineq)
1908 : {
1909 : int r;
1910 : int sgn;
1911 : isl_int cst;
1912 :
1913 39141870197 : if (!tab)
1914 0 : return isl_stat_error;
1915 39141870197 : if (tab->bmap) {
1916 16989989453 : struct isl_basic_map *bmap = tab->bmap;
1917 :
1918 16989989453 : isl_assert(tab->mat->ctx, tab->n_eq == bmap->n_eq,
1919 : return isl_stat_error);
1920 16989989453 : isl_assert(tab->mat->ctx,
1921 : tab->n_con == bmap->n_eq + bmap->n_ineq,
1922 : return isl_stat_error);
1923 16989989453 : tab->bmap = isl_basic_map_add_ineq(tab->bmap, ineq);
1924 16989989453 : if (isl_tab_push(tab, isl_tab_undo_bmap_ineq) < 0)
1925 0 : return isl_stat_error;
1926 16989989453 : if (!tab->bmap)
1927 0 : return isl_stat_error;
1928 : }
1929 39141870197 : if (tab->cone) {
1930 0 : isl_int_init(cst);
1931 0 : isl_int_set_si(cst, 0);
1932 0 : isl_int_swap(ineq[0], cst);
1933 : }
1934 39141870197 : r = isl_tab_add_row(tab, ineq);
1935 39141870197 : if (tab->cone) {
1936 0 : isl_int_swap(ineq[0], cst);
1937 0 : isl_int_clear(cst);
1938 : }
1939 39141870197 : if (r < 0)
1940 0 : return isl_stat_error;
1941 39141870197 : tab->con[r].is_nonneg = 1;
1942 39141870197 : if (isl_tab_push_var(tab, isl_tab_undo_nonneg, &tab->con[r]) < 0)
1943 0 : return isl_stat_error;
1944 39141870197 : if (isl_tab_row_is_redundant(tab, tab->con[r].index)) {
1945 6682444342 : if (isl_tab_mark_redundant(tab, tab->con[r].index) < 0)
1946 0 : return isl_stat_error;
1947 6682444342 : return isl_stat_ok;
1948 : }
1949 :
1950 32459425855 : sgn = restore_row(tab, &tab->con[r]);
1951 32459425855 : if (sgn < -1)
1952 0 : return isl_stat_error;
1953 32459425855 : if (sgn < 0)
1954 7952059790 : return isl_tab_mark_empty(tab);
1955 24507366065 : if (tab->con[r].is_row && isl_tab_row_is_redundant(tab, tab->con[r].index))
1956 0 : if (isl_tab_mark_redundant(tab, tab->con[r].index) < 0)
1957 0 : return isl_stat_error;
1958 24507366065 : return isl_stat_ok;
1959 : }
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 35445 : static int to_col(struct isl_tab *tab, struct isl_tab_var *var)
1966 : {
1967 : int i;
1968 : int row, col;
1969 35445 : unsigned off = 2 + tab->M;
1970 :
1971 35445 : if (!var->is_row)
1972 0 : return 0;
1973 :
1974 72658 : while (isl_int_is_pos(tab->mat->row[var->index][1])) {
1975 6446 : find_pivot(tab, var, NULL, -1, &row, &col);
1976 6446 : isl_assert(tab->mat->ctx, row != -1, return -1);
1977 6446 : if (isl_tab_pivot(tab, row, col) < 0)
1978 0 : return -1;
1979 6446 : if (!var->is_row)
1980 4678 : return 0;
1981 : }
1982 :
1983 36690 : for (i = tab->n_dead; i < tab->n_col; ++i)
1984 36690 : if (!isl_int_is_zero(tab->mat->row[var->index][off + i]))
1985 : break;
1986 :
1987 30767 : isl_assert(tab->mat->ctx, i < tab->n_col, return -1);
1988 30767 : if (isl_tab_pivot(tab, var->index, i) < 0)
1989 0 : return -1;
1990 :
1991 30767 : return 0;
1992 : }
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 53928877 : static struct isl_tab *add_eq(struct isl_tab *tab, isl_int *eq)
2003 : {
2004 : int i;
2005 : int r;
2006 :
2007 53928877 : if (!tab)
2008 0 : return NULL;
2009 53928877 : r = isl_tab_add_row(tab, eq);
2010 53928877 : if (r < 0)
2011 0 : goto error;
2012 :
2013 53928877 : r = tab->con[r].index;
2014 53928877 : i = isl_seq_first_non_zero(tab->mat->row[r] + 2 + tab->M + tab->n_dead,
2015 53928877 : tab->n_col - tab->n_dead);
2016 53928877 : isl_assert(tab->mat->ctx, i >= 0, goto error);
2017 53928877 : i += tab->n_dead;
2018 53928877 : if (isl_tab_pivot(tab, r, i) < 0)
2019 0 : goto error;
2020 53928877 : if (isl_tab_kill_col(tab, i) < 0)
2021 0 : goto error;
2022 53928877 : tab->n_eq++;
2023 :
2024 53928877 : return tab;
2025 : error:
2026 0 : isl_tab_free(tab);
2027 0 : return NULL;
2028 : }
2029 :
2030 : /* Does the sample value of row "row" of "tab" involve the big parameter,
2031 : * if any?
2032 : */
2033 298041099 : static int row_is_big(struct isl_tab *tab, int row)
2034 : {
2035 298041099 : return tab->M && !isl_int_is_zero(tab->mat->row[row][2]);
2036 : }
2037 :
2038 42184 : static int row_is_manifestly_zero(struct isl_tab *tab, int row)
2039 : {
2040 42184 : unsigned off = 2 + tab->M;
2041 :
2042 42184 : if (!isl_int_is_zero(tab->mat->row[row][1]))
2043 4718 : return 0;
2044 37466 : if (row_is_big(tab, row))
2045 0 : return 0;
2046 74932 : return isl_seq_first_non_zero(tab->mat->row[row] + off + tab->n_dead,
2047 74932 : tab->n_col - tab->n_dead) == -1;
2048 : }
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 11671 : int isl_tab_add_valid_eq(struct isl_tab *tab, isl_int *eq)
2056 : {
2057 : struct isl_tab_var *var;
2058 : int r;
2059 :
2060 11671 : if (!tab)
2061 0 : return -1;
2062 11671 : r = isl_tab_add_row(tab, eq);
2063 11671 : if (r < 0)
2064 0 : return -1;
2065 :
2066 11671 : var = &tab->con[r];
2067 11671 : r = var->index;
2068 11671 : if (row_is_manifestly_zero(tab, r)) {
2069 6508 : var->is_zero = 1;
2070 6508 : if (isl_tab_mark_redundant(tab, r) < 0)
2071 0 : return -1;
2072 6508 : return 0;
2073 : }
2074 :
2075 5163 : if (isl_int_is_neg(tab->mat->row[r][1])) {
2076 600 : isl_seq_neg(tab->mat->row[r] + 1, tab->mat->row[r] + 1,
2077 600 : 1 + tab->n_col);
2078 600 : var->negated = 1;
2079 : }
2080 5163 : var->is_nonneg = 1;
2081 5163 : if (to_col(tab, var) < 0)
2082 0 : return -1;
2083 5163 : var->is_nonneg = 0;
2084 5163 : if (isl_tab_kill_col(tab, var->index) < 0)
2085 0 : return -1;
2086 :
2087 5163 : return 0;
2088 : }
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 28796 : static int add_zero_row(struct isl_tab *tab)
2097 : {
2098 : int r;
2099 : isl_int *row;
2100 :
2101 28796 : r = isl_tab_allocate_con(tab);
2102 28796 : if (r < 0)
2103 0 : return -1;
2104 :
2105 28796 : row = tab->mat->row[tab->con[r].index];
2106 28796 : isl_seq_clr(row + 1, 1 + tab->M + tab->n_col);
2107 28796 : isl_int_set_si(row[0], 1);
2108 :
2109 28796 : return r;
2110 : }
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 30513 : int isl_tab_add_eq(struct isl_tab *tab, isl_int *eq)
2120 : {
2121 30513 : struct isl_tab_undo *snap = NULL;
2122 : struct isl_tab_var *var;
2123 : int r;
2124 : int row;
2125 : int sgn;
2126 : isl_int cst;
2127 :
2128 30513 : if (!tab)
2129 0 : return -1;
2130 30513 : isl_assert(tab->mat->ctx, !tab->M, return -1);
2131 :
2132 30513 : if (tab->need_undo)
2133 30513 : snap = isl_tab_snap(tab);
2134 :
2135 30513 : if (tab->cone) {
2136 0 : isl_int_init(cst);
2137 0 : isl_int_set_si(cst, 0);
2138 0 : isl_int_swap(eq[0], cst);
2139 : }
2140 30513 : r = isl_tab_add_row(tab, eq);
2141 30513 : if (tab->cone) {
2142 0 : isl_int_swap(eq[0], cst);
2143 0 : isl_int_clear(cst);
2144 : }
2145 30513 : if (r < 0)
2146 0 : return -1;
2147 :
2148 30513 : var = &tab->con[r];
2149 30513 : row = var->index;
2150 30513 : if (row_is_manifestly_zero(tab, row)) {
2151 231 : if (snap)
2152 231 : return isl_tab_rollback(tab, snap);
2153 0 : return drop_row(tab, row);
2154 : }
2155 :
2156 30282 : if (tab->bmap) {
2157 28796 : tab->bmap = isl_basic_map_add_ineq(tab->bmap, eq);
2158 28796 : if (isl_tab_push(tab, isl_tab_undo_bmap_ineq) < 0)
2159 0 : return -1;
2160 28796 : isl_seq_neg(eq, eq, 1 + tab->n_var);
2161 28796 : tab->bmap = isl_basic_map_add_ineq(tab->bmap, eq);
2162 28796 : isl_seq_neg(eq, eq, 1 + tab->n_var);
2163 28796 : if (isl_tab_push(tab, isl_tab_undo_bmap_ineq) < 0)
2164 0 : return -1;
2165 28796 : if (!tab->bmap)
2166 0 : return -1;
2167 28796 : if (add_zero_row(tab) < 0)
2168 0 : return -1;
2169 : }
2170 :
2171 30282 : sgn = isl_int_sgn(tab->mat->row[row][1]);
2172 :
2173 30282 : if (sgn > 0) {
2174 21 : isl_seq_neg(tab->mat->row[row] + 1, tab->mat->row[row] + 1,
2175 21 : 1 + tab->n_col);
2176 21 : var->negated = 1;
2177 21 : sgn = -1;
2178 : }
2179 :
2180 30282 : if (sgn < 0) {
2181 52 : sgn = sign_of_max(tab, var);
2182 52 : if (sgn < -1)
2183 0 : return -1;
2184 52 : if (sgn < 0) {
2185 0 : if (isl_tab_mark_empty(tab) < 0)
2186 0 : return -1;
2187 0 : return 0;
2188 : }
2189 : }
2190 :
2191 30282 : var->is_nonneg = 1;
2192 30282 : if (to_col(tab, var) < 0)
2193 0 : return -1;
2194 30282 : var->is_nonneg = 0;
2195 30282 : if (isl_tab_kill_col(tab, var->index) < 0)
2196 0 : return -1;
2197 :
2198 30282 : return 0;
2199 : }
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 0 : static struct isl_vec *ineq_for_div(struct isl_basic_map *bmap, unsigned div)
2212 : {
2213 : unsigned total;
2214 : unsigned div_pos;
2215 : struct isl_vec *ineq;
2216 :
2217 0 : if (!bmap)
2218 0 : return NULL;
2219 :
2220 0 : total = isl_basic_map_total_dim(bmap);
2221 0 : div_pos = 1 + total - bmap->n_div + div;
2222 :
2223 0 : ineq = isl_vec_alloc(bmap->ctx, 1 + total);
2224 0 : if (!ineq)
2225 0 : return NULL;
2226 :
2227 0 : isl_seq_cpy(ineq->el, bmap->div[div] + 1, 1 + total);
2228 0 : isl_int_neg(ineq->el[div_pos], bmap->div[div][0]);
2229 0 : return ineq;
2230 : }
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 0 : static isl_stat add_div_constraints(struct isl_tab *tab, unsigned div,
2248 : isl_stat (*add_ineq)(void *user, isl_int *), void *user)
2249 : {
2250 : unsigned total;
2251 : unsigned div_pos;
2252 : struct isl_vec *ineq;
2253 :
2254 0 : total = isl_basic_map_total_dim(tab->bmap);
2255 0 : div_pos = 1 + total - tab->bmap->n_div + div;
2256 :
2257 0 : ineq = ineq_for_div(tab->bmap, div);
2258 0 : if (!ineq)
2259 0 : goto error;
2260 :
2261 0 : if (add_ineq) {
2262 0 : if (add_ineq(user, ineq->el) < 0)
2263 0 : goto error;
2264 : } else {
2265 0 : if (isl_tab_add_ineq(tab, ineq->el) < 0)
2266 0 : goto error;
2267 : }
2268 :
2269 0 : isl_seq_neg(ineq->el, tab->bmap->div[div] + 1, 1 + total);
2270 0 : isl_int_set(ineq->el[div_pos], tab->bmap->div[div][0]);
2271 0 : isl_int_add(ineq->el[0], ineq->el[0], ineq->el[div_pos]);
2272 0 : isl_int_sub_ui(ineq->el[0], ineq->el[0], 1);
2273 :
2274 0 : if (add_ineq) {
2275 0 : if (add_ineq(user, ineq->el) < 0)
2276 0 : goto error;
2277 : } else {
2278 0 : if (isl_tab_add_ineq(tab, ineq->el) < 0)
2279 0 : goto error;
2280 : }
2281 :
2282 0 : isl_vec_free(ineq);
2283 :
2284 0 : return isl_stat_ok;
2285 : error:
2286 0 : isl_vec_free(ineq);
2287 0 : return isl_stat_error;
2288 : }
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 0 : static int div_is_nonneg(struct isl_tab *tab, __isl_keep isl_vec *div)
2297 : {
2298 : int i;
2299 :
2300 0 : if (tab->M)
2301 0 : return 1;
2302 :
2303 0 : if (isl_int_is_neg(div->el[1]))
2304 0 : return 0;
2305 :
2306 0 : for (i = 0; i < tab->n_var; ++i) {
2307 0 : if (isl_int_is_neg(div->el[2 + i]))
2308 0 : return 0;
2309 0 : if (isl_int_is_zero(div->el[2 + i]))
2310 0 : continue;
2311 0 : if (!tab->var[i].is_nonneg)
2312 0 : return 0;
2313 : }
2314 :
2315 0 : return 1;
2316 : }
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 0 : 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 : {
2332 : int r;
2333 : int nonneg;
2334 : int n_div, o_div;
2335 :
2336 0 : if (!tab || !div)
2337 0 : return -1;
2338 :
2339 0 : if (div->size != 1 + 1 + tab->n_var)
2340 0 : isl_die(isl_tab_get_ctx(tab), isl_error_invalid,
2341 : "unexpected size", return -1);
2342 :
2343 0 : isl_assert(tab->mat->ctx, tab->bmap, return -1);
2344 0 : n_div = isl_basic_map_dim(tab->bmap, isl_dim_div);
2345 0 : o_div = tab->n_var - n_div;
2346 0 : if (pos < o_div || pos > tab->n_var)
2347 0 : isl_die(isl_tab_get_ctx(tab), isl_error_invalid,
2348 : "invalid position", return -1);
2349 :
2350 0 : nonneg = div_is_nonneg(tab, div);
2351 :
2352 0 : if (isl_tab_extend_cons(tab, 3) < 0)
2353 0 : return -1;
2354 0 : if (isl_tab_extend_vars(tab, 1) < 0)
2355 0 : return -1;
2356 0 : r = isl_tab_insert_var(tab, pos);
2357 0 : if (r < 0)
2358 0 : return -1;
2359 :
2360 0 : if (nonneg)
2361 0 : tab->var[r].is_nonneg = 1;
2362 :
2363 0 : tab->bmap = isl_basic_map_insert_div(tab->bmap, pos - o_div, div);
2364 0 : if (!tab->bmap)
2365 0 : return -1;
2366 0 : if (isl_tab_push_var(tab, isl_tab_undo_bmap_div, &tab->var[r]) < 0)
2367 0 : return -1;
2368 :
2369 0 : if (add_div_constraints(tab, pos - o_div, add_ineq, user) < 0)
2370 0 : return -1;
2371 :
2372 0 : return r;
2373 : }
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 0 : int isl_tab_add_div(struct isl_tab *tab, __isl_keep isl_vec *div)
2379 : {
2380 0 : if (!tab)
2381 0 : return -1;
2382 0 : return isl_tab_insert_div(tab, tab->n_var, div, NULL, NULL);
2383 : }
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 2955142263 : __isl_give struct isl_tab *isl_tab_from_basic_map(
2391 : __isl_keep isl_basic_map *bmap, int track)
2392 : {
2393 : int i;
2394 : struct isl_tab *tab;
2395 :
2396 2955142263 : if (!bmap)
2397 0 : return NULL;
2398 5910284526 : tab = isl_tab_alloc(bmap->ctx,
2399 2955142263 : isl_basic_map_total_dim(bmap) + bmap->n_ineq + 1,
2400 : isl_basic_map_total_dim(bmap), 0);
2401 2955142263 : if (!tab)
2402 0 : return NULL;
2403 2955142263 : tab->preserve = track;
2404 2955142263 : tab->rational = ISL_F_ISSET(bmap, ISL_BASIC_MAP_RATIONAL);
2405 2955142263 : if (ISL_F_ISSET(bmap, ISL_BASIC_MAP_EMPTY)) {
2406 0 : if (isl_tab_mark_empty(tab) < 0)
2407 0 : goto error;
2408 0 : goto done;
2409 : }
2410 3009071140 : for (i = 0; i < bmap->n_eq; ++i) {
2411 53928877 : tab = add_eq(tab, bmap->eq[i]);
2412 53928877 : if (!tab)
2413 0 : return tab;
2414 : }
2415 23807754325 : for (i = 0; i < bmap->n_ineq; ++i) {
2416 22149390299 : if (isl_tab_add_ineq(tab, bmap->ineq[i]) < 0)
2417 0 : goto error;
2418 22149390299 : if (tab->empty)
2419 1296778237 : goto done;
2420 : }
2421 : done:
2422 2955142263 : if (track && isl_tab_track_bmap(tab, isl_basic_map_copy(bmap)) < 0)
2423 0 : goto error;
2424 2955142263 : return tab;
2425 : error:
2426 0 : isl_tab_free(tab);
2427 0 : return NULL;
2428 : }
2429 :
2430 1356169760 : __isl_give struct isl_tab *isl_tab_from_basic_set(
2431 : __isl_keep isl_basic_set *bset, int track)
2432 : {
2433 1356169760 : return isl_tab_from_basic_map(bset, track);
2434 : }
2435 :
2436 : /* Construct a tableau corresponding to the recession cone of "bset".
2437 : */
2438 0 : struct isl_tab *isl_tab_from_recession_cone(__isl_keep isl_basic_set *bset,
2439 : int parametric)
2440 : {
2441 : isl_int cst;
2442 : int i;
2443 : struct isl_tab *tab;
2444 0 : unsigned offset = 0;
2445 :
2446 0 : if (!bset)
2447 0 : return NULL;
2448 0 : if (parametric)
2449 0 : offset = isl_basic_set_dim(bset, isl_dim_param);
2450 0 : tab = isl_tab_alloc(bset->ctx, bset->n_eq + bset->n_ineq,
2451 0 : isl_basic_set_total_dim(bset) - offset, 0);
2452 0 : if (!tab)
2453 0 : return NULL;
2454 0 : tab->rational = ISL_F_ISSET(bset, ISL_BASIC_SET_RATIONAL);
2455 0 : tab->cone = 1;
2456 :
2457 0 : isl_int_init(cst);
2458 0 : isl_int_set_si(cst, 0);
2459 0 : for (i = 0; i < bset->n_eq; ++i) {
2460 0 : isl_int_swap(bset->eq[i][offset], cst);
2461 0 : if (offset > 0) {
2462 0 : if (isl_tab_add_eq(tab, bset->eq[i] + offset) < 0)
2463 0 : goto error;
2464 : } else
2465 0 : tab = add_eq(tab, bset->eq[i]);
2466 0 : isl_int_swap(bset->eq[i][offset], cst);
2467 0 : if (!tab)
2468 0 : goto done;
2469 : }
2470 0 : for (i = 0; i < bset->n_ineq; ++i) {
2471 : int r;
2472 0 : isl_int_swap(bset->ineq[i][offset], cst);
2473 0 : r = isl_tab_add_row(tab, bset->ineq[i] + offset);
2474 0 : isl_int_swap(bset->ineq[i][offset], cst);
2475 0 : if (r < 0)
2476 0 : goto error;
2477 0 : tab->con[r].is_nonneg = 1;
2478 0 : if (isl_tab_push_var(tab, isl_tab_undo_nonneg, &tab->con[r]) < 0)
2479 0 : goto error;
2480 : }
2481 : done:
2482 0 : isl_int_clear(cst);
2483 0 : return tab;
2484 : error:
2485 0 : isl_int_clear(cst);
2486 0 : isl_tab_free(tab);
2487 0 : return NULL;
2488 : }
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 0 : isl_bool isl_tab_cone_is_bounded(struct isl_tab *tab)
2494 : {
2495 : int i;
2496 :
2497 0 : if (!tab)
2498 0 : return isl_bool_error;
2499 0 : if (tab->empty)
2500 0 : return isl_bool_true;
2501 0 : if (tab->n_dead == tab->n_col)
2502 0 : return isl_bool_true;
2503 :
2504 : for (;;) {
2505 0 : for (i = tab->n_redundant; i < tab->n_row; ++i) {
2506 : struct isl_tab_var *var;
2507 : int sgn;
2508 0 : var = isl_tab_var_from_row(tab, i);
2509 0 : if (!var->is_nonneg)
2510 0 : continue;
2511 0 : sgn = sign_of_max(tab, var);
2512 0 : if (sgn < -1)
2513 0 : return isl_bool_error;
2514 0 : if (sgn != 0)
2515 0 : return isl_bool_false;
2516 0 : if (close_row(tab, var, 0) < 0)
2517 0 : return isl_bool_error;
2518 0 : break;
2519 : }
2520 0 : if (tab->n_dead == tab->n_col)
2521 0 : return isl_bool_true;
2522 0 : if (i == tab->n_row)
2523 0 : return isl_bool_false;
2524 0 : }
2525 : }
2526 :
2527 1446768798 : int isl_tab_sample_is_integer(struct isl_tab *tab)
2528 : {
2529 : int i;
2530 :
2531 1446768798 : if (!tab)
2532 0 : return -1;
2533 :
2534 8359711171 : for (i = 0; i < tab->n_var; ++i) {
2535 : int row;
2536 6932658308 : if (!tab->var[i].is_row)
2537 649694491 : continue;
2538 6282963817 : row = tab->var[i].index;
2539 6282963817 : if (!isl_int_is_divisible_by(tab->mat->row[row][1],
2540 : tab->mat->row[row][0]))
2541 19715935 : return 0;
2542 : }
2543 1427052863 : return 1;
2544 : }
2545 :
2546 1426627757 : static struct isl_vec *extract_integer_sample(struct isl_tab *tab)
2547 : {
2548 : int i;
2549 : struct isl_vec *vec;
2550 :
2551 1426627757 : vec = isl_vec_alloc(tab->mat->ctx, 1 + tab->n_var);
2552 1426627757 : if (!vec)
2553 0 : return NULL;
2554 :
2555 1426627757 : isl_int_set_si(vec->block.data[0], 1);
2556 8321923361 : for (i = 0; i < tab->n_var; ++i) {
2557 6895295604 : if (!tab->var[i].is_row)
2558 646486075 : isl_int_set_si(vec->block.data[1 + i], 0);
2559 : else {
2560 6248809529 : int row = tab->var[i].index;
2561 6248809529 : isl_int_divexact(vec->block.data[1 + i],
2562 : tab->mat->row[row][1], tab->mat->row[row][0]);
2563 : }
2564 : }
2565 :
2566 1426627757 : return vec;
2567 : }
2568 :
2569 59781867 : struct isl_vec *isl_tab_get_sample_value(struct isl_tab *tab)
2570 : {
2571 : int i;
2572 : struct isl_vec *vec;
2573 : isl_int m;
2574 :
2575 59781867 : if (!tab)
2576 0 : return NULL;
2577 :
2578 59781867 : vec = isl_vec_alloc(tab->mat->ctx, 1 + tab->n_var);
2579 59781867 : if (!vec)
2580 0 : return NULL;
2581 :
2582 59781867 : isl_int_init(m);
2583 :
2584 59781867 : isl_int_set_si(vec->block.data[0], 1);
2585 367677708 : for (i = 0; i < tab->n_var; ++i) {
2586 : int row;
2587 307895841 : if (!tab->var[i].is_row) {
2588 80894247 : isl_int_set_si(vec->block.data[1 + i], 0);
2589 80894247 : continue;
2590 : }
2591 227001594 : row = tab->var[i].index;
2592 227001594 : isl_int_gcd(m, vec->block.data[0], tab->mat->row[row][0]);
2593 227001594 : isl_int_divexact(m, tab->mat->row[row][0], m);
2594 227001594 : isl_seq_scale(vec->block.data, vec->block.data, m, 1 + i);
2595 227001594 : isl_int_divexact(m, vec->block.data[0], tab->mat->row[row][0]);
2596 227001594 : isl_int_mul(vec->block.data[1 + i], m, tab->mat->row[row][1]);
2597 : }
2598 59781867 : vec = isl_vec_normalize(vec);
2599 :
2600 59781867 : isl_int_clear(m);
2601 59781867 : return vec;
2602 : }
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 4164876 : static void get_rounded_sample_value(struct isl_tab *tab,
2608 : struct isl_tab_var *var, int sgn, isl_int *v)
2609 : {
2610 4164876 : if (!var->is_row)
2611 435315 : isl_int_set_si(*v, 0);
2612 3729561 : else if (sgn > 0)
2613 1129987 : isl_int_cdiv_q(*v, tab->mat->row[var->index][1],
2614 : tab->mat->row[var->index][0]);
2615 : else
2616 2599574 : isl_int_fdiv_q(*v, tab->mat->row[var->index][1],
2617 : tab->mat->row[var->index][0]);
2618 4164876 : }
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 1612366175 : struct isl_basic_map *isl_basic_map_update_from_tab(struct isl_basic_map *bmap,
2629 : struct isl_tab *tab)
2630 : {
2631 : int i;
2632 : unsigned n_eq;
2633 :
2634 1612366175 : if (!bmap)
2635 0 : return NULL;
2636 1612366175 : if (!tab)
2637 0 : return bmap;
2638 :
2639 1612366175 : n_eq = tab->n_eq;
2640 1612366175 : if (tab->empty)
2641 410118 : bmap = isl_basic_map_set_to_empty(bmap);
2642 : else
2643 15707184130 : for (i = bmap->n_ineq - 1; i >= 0; --i) {
2644 14095228073 : if (isl_tab_is_equality(tab, n_eq + i))
2645 9935413134 : isl_basic_map_inequality_to_equality(bmap, i);
2646 4159814939 : else if (isl_tab_is_redundant(tab, n_eq + i))
2647 1699749293 : isl_basic_map_drop_inequality(bmap, i);
2648 : }
2649 1612366175 : if (bmap->n_eq != n_eq)
2650 1298401593 : bmap = isl_basic_map_gauss(bmap, NULL);
2651 1612366175 : if (!tab->rational &&
2652 1612366175 : bmap && !bmap->sample && isl_tab_sample_is_integer(tab))
2653 1426627757 : bmap->sample = extract_integer_sample(tab);
2654 1612366175 : return bmap;
2655 : }
2656 :
2657 353387 : struct isl_basic_set *isl_basic_set_update_from_tab(struct isl_basic_set *bset,
2658 : struct isl_tab *tab)
2659 : {
2660 353387 : return bset_from_bmap(isl_basic_map_update_from_tab(bset_to_bmap(bset),
2661 : tab));
2662 : }
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 324681 : static isl_stat drop_last_con_in_row(struct isl_tab *tab, int r)
2668 : {
2669 324681 : if (!tab->con[r].is_row)
2670 0 : isl_die(isl_tab_get_ctx(tab), isl_error_internal,
2671 : "row unexpectedly moved to column",
2672 : return isl_stat_error);
2673 324681 : if (r + 1 != tab->n_con)
2674 0 : isl_die(isl_tab_get_ctx(tab), isl_error_internal,
2675 : "additional constraints added", return isl_stat_error);
2676 324681 : if (drop_row(tab, tab->con[r].index) < 0)
2677 0 : return isl_stat_error;
2678 :
2679 324681 : return isl_stat_ok;
2680 : }
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 324681 : static isl_stat cut_to_hyperplane(struct isl_tab *tab, struct isl_tab_var *var)
2693 : {
2694 : unsigned r;
2695 : isl_int *row;
2696 : int sgn;
2697 324681 : unsigned off = 2 + tab->M;
2698 :
2699 324681 : if (var->is_zero)
2700 0 : return isl_stat_ok;
2701 324681 : if (var->is_redundant || !var->is_nonneg)
2702 0 : isl_die(isl_tab_get_ctx(tab), isl_error_invalid,
2703 : "expecting non-redundant non-negative variable",
2704 : return isl_stat_error);
2705 :
2706 324681 : if (isl_tab_extend_cons(tab, 1) < 0)
2707 0 : return isl_stat_error;
2708 :
2709 324681 : r = tab->n_con;
2710 324681 : tab->con[r].index = tab->n_row;
2711 324681 : tab->con[r].is_row = 1;
2712 324681 : tab->con[r].is_nonneg = 0;
2713 324681 : tab->con[r].is_zero = 0;
2714 324681 : tab->con[r].is_redundant = 0;
2715 324681 : tab->con[r].frozen = 0;
2716 324681 : tab->con[r].negated = 0;
2717 324681 : tab->row_var[tab->n_row] = ~r;
2718 324681 : row = tab->mat->row[tab->n_row];
2719 :
2720 324681 : if (var->is_row) {
2721 69394 : isl_int_set(row[0], tab->mat->row[var->index][0]);
2722 138788 : isl_seq_neg(row + 1,
2723 138788 : tab->mat->row[var->index] + 1, 1 + tab->n_col);
2724 : } else {
2725 255287 : isl_int_set_si(row[0], 1);
2726 255287 : isl_seq_clr(row + 1, 1 + tab->n_col);
2727 255287 : isl_int_set_si(row[off + var->index], -1);
2728 : }
2729 :
2730 324681 : tab->n_row++;
2731 324681 : tab->n_con++;
2732 :
2733 324681 : sgn = sign_of_max(tab, &tab->con[r]);
2734 324681 : if (sgn < -1)
2735 0 : return isl_stat_error;
2736 324681 : if (sgn < 0) {
2737 1533 : if (drop_last_con_in_row(tab, r) < 0)
2738 0 : return isl_stat_error;
2739 1533 : if (isl_tab_mark_empty(tab) < 0)
2740 0 : return isl_stat_error;
2741 1533 : return isl_stat_ok;
2742 : }
2743 323148 : tab->con[r].is_nonneg = 1;
2744 : /* sgn == 0 */
2745 323148 : if (close_row(tab, &tab->con[r], 1) < 0)
2746 0 : return isl_stat_error;
2747 323148 : if (drop_last_con_in_row(tab, r) < 0)
2748 0 : return isl_stat_error;
2749 :
2750 323148 : return isl_stat_ok;
2751 : }
2752 :
2753 : /* Check that "con" is a valid constraint position for "tab".
2754 : */
2755 4961767711 : static isl_stat isl_tab_check_con(struct isl_tab *tab, int con)
2756 : {
2757 4961767711 : if (!tab)
2758 0 : return isl_stat_error;
2759 4961767711 : if (con < 0 || con >= tab->n_con)
2760 0 : isl_die(isl_tab_get_ctx(tab), isl_error_invalid,
2761 : "position out of bounds", return isl_stat_error);
2762 4961767711 : return isl_stat_ok;
2763 : }
2764 :
2765 : /* Given a tableau "tab" and an inequality constraint "con" of the tableau,
2766 : * relax the inequality by one. That is, the inequality r >= 0 is replaced
2767 : * by r' = r + 1 >= 0.
2768 : * If r is a row variable, we simply increase the constant term by one
2769 : * (taking into account the denominator).
2770 : * If r is a column variable, then we need to modify each row that
2771 : * refers to r = r' - 1 by substituting this equality, effectively
2772 : * subtracting the coefficient of the column from the constant.
2773 : * We should only do this if the minimum is manifestly unbounded,
2774 : * however. Otherwise, we may end up with negative sample values
2775 : * for non-negative variables.
2776 : * So, if r is a column variable with a minimum that is not
2777 : * manifestly unbounded, then we need to move it to a row.
2778 : * However, the sample value of this row may be negative,
2779 : * even after the relaxation, so we need to restore it.
2780 : * We therefore prefer to pivot a column up to a row, if possible.
2781 : */
2782 123237 : int isl_tab_relax(struct isl_tab *tab, int con)
2783 : {
2784 : struct isl_tab_var *var;
2785 :
2786 123237 : if (!tab)
2787 0 : return -1;
2788 :
2789 123237 : var = &tab->con[con];
2790 :
2791 123237 : if (var->is_row && (var->index < 0 || var->index < tab->n_redundant))
2792 0 : isl_die(tab->mat->ctx, isl_error_invalid,
2793 : "cannot relax redundant constraint", return -1);
2794 123237 : if (!var->is_row && (var->index < 0 || var->index < tab->n_dead))
2795 0 : isl_die(tab->mat->ctx, isl_error_invalid,
2796 : "cannot relax dead constraint", return -1);
2797 :
2798 123237 : if (!var->is_row && !max_is_manifestly_unbounded(tab, var))
2799 19530 : if (to_row(tab, var, 1) < 0)
2800 0 : return -1;
2801 123237 : if (!var->is_row && !min_is_manifestly_unbounded(tab, var))
2802 38453 : if (to_row(tab, var, -1) < 0)
2803 0 : return -1;
2804 :
2805 123237 : if (var->is_row) {
2806 58603 : isl_int_add(tab->mat->row[var->index][1],
2807 : tab->mat->row[var->index][1], tab->mat->row[var->index][0]);
2808 58603 : if (restore_row(tab, var) < 0)
2809 0 : return -1;
2810 : } else {
2811 : int i;
2812 64634 : unsigned off = 2 + tab->M;
2813 :
2814 365488 : for (i = 0; i < tab->n_row; ++i) {
2815 300854 : if (isl_int_is_zero(tab->mat->row[i][off + var->index]))
2816 152416 : continue;
2817 148438 : isl_int_sub(tab->mat->row[i][1], tab->mat->row[i][1],
2818 : tab->mat->row[i][off + var->index]);
2819 : }
2820 :
2821 : }
2822 :
2823 123237 : if (isl_tab_push_var(tab, isl_tab_undo_relax, var) < 0)
2824 0 : return -1;
2825 :
2826 123237 : return 0;
2827 : }
2828 :
2829 : /* Replace the variable v at position "pos" in the tableau "tab"
2830 : * by v' = v + shift.
2831 : *
2832 : * If the variable is in a column, then we first check if we can
2833 : * simply plug in v = v' - shift. The effect on a row with
2834 : * coefficient f/d for variable v is that the constant term c/d
2835 : * is replaced by (c - f * shift)/d. If shift is positive and
2836 : * f is negative for each row that needs to remain non-negative,
2837 : * then this is clearly safe. In other words, if the minimum of v
2838 : * is manifestly unbounded, then we can keep v in a column position.
2839 : * Otherwise, we can pivot it down to a row.
2840 : * Similarly, if shift is negative, we need to check if the maximum
2841 : * of is manifestly unbounded.
2842 : *
2843 : * If the variable is in a row (from the start or after pivoting),
2844 : * then the constant term c/d is replaced by (c + d * shift)/d.
2845 : */
2846 0 : int isl_tab_shift_var(struct isl_tab *tab, int pos, isl_int shift)
2847 : {
2848 : struct isl_tab_var *var;
2849 :
2850 0 : if (!tab)
2851 0 : return -1;
2852 0 : if (isl_int_is_zero(shift))
2853 0 : return 0;
2854 :
2855 0 : var = &tab->var[pos];
2856 0 : if (!var->is_row) {
2857 0 : if (isl_int_is_neg(shift)) {
2858 0 : if (!max_is_manifestly_unbounded(tab, var))
2859 0 : if (to_row(tab, var, 1) < 0)
2860 0 : return -1;
2861 : } else {
2862 0 : if (!min_is_manifestly_unbounded(tab, var))
2863 0 : if (to_row(tab, var, -1) < 0)
2864 0 : return -1;
2865 : }
2866 : }
2867 :
2868 0 : if (var->is_row) {
2869 0 : isl_int_addmul(tab->mat->row[var->index][1],
2870 : shift, tab->mat->row[var->index][0]);
2871 : } else {
2872 : int i;
2873 0 : unsigned off = 2 + tab->M;
2874 :
2875 0 : for (i = 0; i < tab->n_row; ++i) {
2876 0 : if (isl_int_is_zero(tab->mat->row[i][off + var->index]))
2877 0 : continue;
2878 0 : isl_int_submul(tab->mat->row[i][1],
2879 : shift, tab->mat->row[i][off + var->index]);
2880 : }
2881 :
2882 : }
2883 :
2884 0 : return 0;
2885 : }
2886 :
2887 : /* Remove the sign constraint from constraint "con".
2888 : *
2889 : * If the constraint variable was originally marked non-negative,
2890 : * then we make sure we mark it non-negative again during rollback.
2891 : */
2892 518583 : int isl_tab_unrestrict(struct isl_tab *tab, int con)
2893 : {
2894 : struct isl_tab_var *var;
2895 :
2896 518583 : if (!tab)
2897 0 : return -1;
2898 :
2899 518583 : var = &tab->con[con];
2900 518583 : if (!var->is_nonneg)
2901 0 : return 0;
2902 :
2903 518583 : var->is_nonneg = 0;
2904 518583 : if (isl_tab_push_var(tab, isl_tab_undo_unrestrict, var) < 0)
2905 0 : return -1;
2906 :
2907 518583 : return 0;
2908 : }
2909 :
2910 310091 : int isl_tab_select_facet(struct isl_tab *tab, int con)
2911 : {
2912 310091 : if (!tab)
2913 0 : return -1;
2914 :
2915 310091 : return cut_to_hyperplane(tab, &tab->con[con]);
2916 : }
2917 :
2918 27577425321 : static int may_be_equality(struct isl_tab *tab, int row)
2919 : {
2920 55154850642 : return tab->rational ? isl_int_is_zero(tab->mat->row[row][1])
2921 27577425321 : : isl_int_lt(tab->mat->row[row][1],
2922 : tab->mat->row[row][0]);
2923 : }
2924 :
2925 : /* Return an isl_tab_var that has been marked or NULL if no such
2926 : * variable can be found.
2927 : * The marked field has only been set for variables that
2928 : * appear in non-redundant rows or non-dead columns.
2929 : *
2930 : * Pick the last constraint variable that is marked and
2931 : * that appears in either a non-redundant row or a non-dead columns.
2932 : * Since the returned variable is tested for being a redundant constraint or
2933 : * an implicit equality, there is no need to return any tab variable that
2934 : * corresponds to a variable.
2935 : */
2936 9210084223 : static struct isl_tab_var *select_marked(struct isl_tab *tab)
2937 : {
2938 : int i;
2939 : struct isl_tab_var *var;
2940 :
2941 51232836024 : for (i = tab->n_con - 1; i >= 0; --i) {
2942 49869661657 : var = &tab->con[i];
2943 49869661657 : if (var->index < 0)
2944 34114197872 : continue;
2945 15755463785 : if (var->is_row && var->index < tab->n_redundant)
2946 588346132 : continue;
2947 15167117653 : if (!var->is_row && var->index < tab->n_dead)
2948 70522 : continue;
2949 15167047131 : if (var->marked)
2950 7846909856 : return var;
2951 : }
2952 :
2953 1363174367 : return NULL;
2954 : }
2955 :
2956 : /* Check for (near) equalities among the constraints.
2957 : * A constraint is an equality if it is non-negative and if
2958 : * its maximal value is either
2959 : * - zero (in case of rational tableaus), or
2960 : * - strictly less than 1 (in case of integer tableaus)
2961 : *
2962 : * We first mark all non-redundant and non-dead variables that
2963 : * are not frozen and not obviously not an equality.
2964 : * Then we iterate over all marked variables if they can attain
2965 : * any values larger than zero or at least one.
2966 : * If the maximal value is zero, we mark any column variables
2967 : * that appear in the row as being zero and mark the row as being redundant.
2968 : * Otherwise, if the maximal value is strictly less than one (and the
2969 : * tableau is integer), then we restrict the value to being zero
2970 : * by adding an opposite non-negative variable.
2971 : * The order in which the variables are considered is not important.
2972 : */
2973 1566395957 : int isl_tab_detect_implicit_equalities(struct isl_tab *tab)
2974 : {
2975 : int i;
2976 : unsigned n_marked;
2977 :
2978 1566395957 : if (!tab)
2979 0 : return -1;
2980 1566395957 : if (tab->empty)
2981 410093 : return 0;
2982 1565985864 : if (tab->n_dead == tab->n_col)
2983 1268 : return 0;
2984 :
2985 1565984596 : n_marked = 0;
2986 15027679499 : for (i = tab->n_redundant; i < tab->n_row; ++i) {
2987 13461694903 : struct isl_tab_var *var = isl_tab_var_from_row(tab, i);
2988 26542243981 : var->marked = !var->frozen && var->is_nonneg &&
2989 13080549078 : may_be_equality(tab, i);
2990 13461694903 : if (var->marked)
2991 12892077715 : n_marked++;
2992 : }
2993 9108349332 : for (i = tab->n_dead; i < tab->n_col; ++i) {
2994 7542364736 : struct isl_tab_var *var = var_from_col(tab, i);
2995 7542364736 : var->marked = !var->frozen && var->is_nonneg;
2996 7542364736 : if (var->marked)
2997 522999882 : n_marked++;
2998 : }
2999 10273483924 : while (n_marked) {
3000 : struct isl_tab_var *var;
3001 : int sgn;
3002 8479656316 : var = select_marked(tab);
3003 8479656316 : if (!var)
3004 1338126994 : break;
3005 7141529322 : var->marked = 0;
3006 7141529322 : n_marked--;
3007 7141529322 : sgn = sign_of_max(tab, var);
3008 7141529322 : if (sgn < 0)
3009 0 : return -1;
3010 7141529322 : if (sgn == 0) {
3011 4946280941 : if (close_row(tab, var, 0) < 0)
3012 0 : return -1;
3013 2195248381 : } else if (!tab->rational && !at_least_one(tab, var)) {
3014 14590 : if (cut_to_hyperplane(tab, var) < 0)
3015 0 : return -1;
3016 14590 : return isl_tab_detect_implicit_equalities(tab);
3017 : }
3018 43607922795 : for (i = tab->n_redundant; i < tab->n_row; ++i) {
3019 36466408063 : var = isl_tab_var_from_row(tab, i);
3020 36466408063 : if (!var->marked)
3021 21969531820 : continue;
3022 14496876243 : if (may_be_equality(tab, i))
3023 14488463128 : continue;
3024 8413115 : var->marked = 0;
3025 8413115 : n_marked--;
3026 : }
3027 : }
3028 :
3029 1565970006 : return 0;
3030 : }
3031 :
3032 : /* Update the element of row_var or col_var that corresponds to
3033 : * constraint tab->con[i] to a move from position "old" to position "i".
3034 : */
3035 130101 : static int update_con_after_move(struct isl_tab *tab, int i, int old)
3036 : {
3037 : int *p;
3038 : int index;
3039 :
3040 130101 : index = tab->con[i].index;
3041 130101 : if (index == -1)
3042 73328 : return 0;
3043 56773 : p = tab->con[i].is_row ? tab->row_var : tab->col_var;
3044 56773 : if (p[index] != ~old)
3045 0 : isl_die(tab->mat->ctx, isl_error_internal,
3046 : "broken internal state", return -1);
3047 56773 : p[index] = ~i;
3048 :
3049 56773 : return 0;
3050 : }
3051 :
3052 : /* Interchange constraints "con1" and "con2" in "tab".
3053 : * In particular, interchange the contents of these entries in tab->con.
3054 : * Since tab->col_var and tab->row_var point back into this array,
3055 : * they need to be updated accordingly.
3056 : */
3057 2593 : isl_stat isl_tab_swap_constraints(struct isl_tab *tab, int con1, int con2)
3058 : {
3059 : struct isl_tab_var var;
3060 :
3061 5186 : if (isl_tab_check_con(tab, con1) < 0 ||
3062 2593 : isl_tab_check_con(tab, con2) < 0)
3063 0 : return isl_stat_error;
3064 :
3065 2593 : var = tab->con[con1];
3066 2593 : tab->con[con1] = tab->con[con2];
3067 2593 : if (update_con_after_move(tab, con1, con2) < 0)
3068 0 : return isl_stat_error;
3069 2593 : tab->con[con2] = var;
3070 2593 : if (update_con_after_move(tab, con2, con1) < 0)
3071 0 : return isl_stat_error;
3072 :
3073 2593 : return isl_stat_ok;
3074 : }
3075 :
3076 : /* Rotate the "n" constraints starting at "first" to the right,
3077 : * putting the last constraint in the position of the first constraint.
3078 : */
3079 24480 : static int rotate_constraints(struct isl_tab *tab, int first, int n)
3080 : {
3081 : int i, last;
3082 : struct isl_tab_var var;
3083 :
3084 24480 : if (n <= 1)
3085 5103 : return 0;
3086 :
3087 19377 : last = first + n - 1;
3088 19377 : var = tab->con[last];
3089 106905 : for (i = last; i > first; --i) {
3090 87528 : tab->con[i] = tab->con[i - 1];
3091 87528 : if (update_con_after_move(tab, i, i - 1) < 0)
3092 0 : return -1;
3093 : }
3094 19377 : tab->con[first] = var;
3095 19377 : if (update_con_after_move(tab, first, last) < 0)
3096 0 : return -1;
3097 :
3098 19377 : return 0;
3099 : }
3100 :
3101 : /* Drop the "n" entries starting at position "first" in tab->con, moving all
3102 : * subsequent entries down.
3103 : * Since some of the entries of tab->row_var and tab->col_var contain
3104 : * indices into this array, they have to be updated accordingly.
3105 : */
3106 3387 : static isl_stat con_drop_entries(struct isl_tab *tab,
3107 : unsigned first, unsigned n)
3108 : {
3109 : int i;
3110 :
3111 3387 : if (first + n > tab->n_con || first + n < first)
3112 0 : isl_die(isl_tab_get_ctx(tab), isl_error_internal,
3113 : "invalid range", return isl_stat_error);
3114 :
3115 3387 : tab->n_con -= n;
3116 :
3117 21397 : for (i = first; i < tab->n_con; ++i) {
3118 18010 : tab->con[i] = tab->con[i + n];
3119 18010 : if (update_con_after_move(tab, i, i + n) < 0)
3120 0 : return isl_stat_error;
3121 : }
3122 :
3123 3387 : return isl_stat_ok;
3124 : }
3125 :
3126 : /* isl_basic_map_gauss5 callback that gets called when
3127 : * two (equality) constraints "a" and "b" get interchanged
3128 : * in the basic map. Perform the same interchange in "tab".
3129 : */
3130 2593 : static isl_stat swap_eq(unsigned a, unsigned b, void *user)
3131 : {
3132 2593 : struct isl_tab *tab = user;
3133 :
3134 2593 : return isl_tab_swap_constraints(tab, a, b);
3135 : }
3136 :
3137 : /* isl_basic_map_gauss5 callback that gets called when
3138 : * the final "n" equality constraints get removed.
3139 : * As a special case, if "n" is equal to the total number
3140 : * of equality constraints, then this means the basic map
3141 : * turned out to be empty.
3142 : * Drop the same number of equality constraints from "tab" or
3143 : * mark it empty in the special case.
3144 : */
3145 3387 : static isl_stat drop_eq(unsigned n, void *user)
3146 : {
3147 3387 : struct isl_tab *tab = user;
3148 :
3149 3387 : if (tab->n_eq == n)
3150 0 : return isl_tab_mark_empty(tab);
3151 :
3152 3387 : tab->n_eq -= n;
3153 3387 : return con_drop_entries(tab, tab->n_eq, n);
3154 : }
3155 :
3156 : /* If "bmap" has more than a single reference, then call
3157 : * isl_basic_map_gauss on it, updating "tab" accordingly.
3158 : */
3159 3705 : static __isl_give isl_basic_map *gauss_if_shared(__isl_take isl_basic_map *bmap,
3160 : struct isl_tab *tab)
3161 : {
3162 : isl_bool single;
3163 :
3164 3705 : single = isl_basic_map_has_single_reference(bmap);
3165 3705 : if (single < 0)
3166 0 : return isl_basic_map_free(bmap);
3167 3705 : if (single)
3168 314 : return bmap;
3169 3391 : return isl_basic_map_gauss5(bmap, NULL, &swap_eq, &drop_eq, tab);
3170 : }
3171 :
3172 : /* Make the equalities that are implicit in "bmap" but that have been
3173 : * detected in the corresponding "tab" explicit in "bmap" and update
3174 : * "tab" to reflect the new order of the constraints.
3175 : *
3176 : * In particular, if inequality i is an implicit equality then
3177 : * isl_basic_map_inequality_to_equality will move the inequality
3178 : * in front of the other equality and it will move the last inequality
3179 : * in the position of inequality i.
3180 : * In the tableau, the inequalities of "bmap" are stored after the equalities
3181 : * and so the original order
3182 : *
3183 : * E E E E E A A A I B B B B L
3184 : *
3185 : * is changed into
3186 : *
3187 : * I E E E E E A A A L B B B B
3188 : *
3189 : * where I is the implicit equality, the E are equalities,
3190 : * the A inequalities before I, the B inequalities after I and
3191 : * L the last inequality.
3192 : * We therefore need to rotate to the right two sets of constraints,
3193 : * those up to and including I and those after I.
3194 : *
3195 : * If "tab" contains any constraints that are not in "bmap" then they
3196 : * appear after those in "bmap" and they should be left untouched.
3197 : *
3198 : * Note that this function only calls isl_basic_map_gauss
3199 : * (in case some equality constraints got detected)
3200 : * if "bmap" has more than one reference.
3201 : * If it only has a single reference, then it is left in a temporary state,
3202 : * because the caller may require this state.
3203 : * Calling isl_basic_map_gauss is then the responsibility of the caller.
3204 : */
3205 713039 : __isl_give isl_basic_map *isl_tab_make_equalities_explicit(struct isl_tab *tab,
3206 : __isl_take isl_basic_map *bmap)
3207 : {
3208 : int i;
3209 : unsigned n_eq;
3210 :
3211 713039 : if (!tab || !bmap)
3212 0 : return isl_basic_map_free(bmap);
3213 713039 : if (tab->empty)
3214 3 : return bmap;
3215 :
3216 713036 : n_eq = tab->n_eq;
3217 5162815 : for (i = bmap->n_ineq - 1; i >= 0; --i) {
3218 4449779 : if (!isl_tab_is_equality(tab, bmap->n_eq + i))
3219 4437539 : continue;
3220 12240 : isl_basic_map_inequality_to_equality(bmap, i);
3221 12240 : if (rotate_constraints(tab, 0, tab->n_eq + i + 1) < 0)
3222 0 : return isl_basic_map_free(bmap);
3223 12240 : if (rotate_constraints(tab, tab->n_eq + i + 1,
3224 12240 : bmap->n_ineq - i) < 0)
3225 0 : return isl_basic_map_free(bmap);
3226 12240 : tab->n_eq++;
3227 : }
3228 :
3229 713036 : if (n_eq != tab->n_eq)
3230 3705 : bmap = gauss_if_shared(bmap, tab);
3231 :
3232 713036 : return bmap;
3233 : }
3234 :
3235 767596693 : static int con_is_redundant(struct isl_tab *tab, struct isl_tab_var *var)
3236 : {
3237 767596693 : if (!tab)
3238 0 : return -1;
3239 767596693 : if (tab->rational) {
3240 61013 : int sgn = sign_of_min(tab, var);
3241 61013 : if (sgn < -1)
3242 0 : return -1;
3243 61013 : return sgn >= 0;
3244 : } else {
3245 767535680 : int irred = isl_tab_min_at_most_neg_one(tab, var);
3246 767535680 : if (irred < 0)
3247 0 : return -1;
3248 767535680 : return !irred;
3249 : }
3250 : }
3251 :
3252 : /* Check for (near) redundant constraints.
3253 : * A constraint is redundant if it is non-negative and if
3254 : * its minimal value (temporarily ignoring the non-negativity) is either
3255 : * - zero (in case of rational tableaus), or
3256 : * - strictly larger than -1 (in case of integer tableaus)
3257 : *
3258 : * We first mark all non-redundant and non-dead variables that
3259 : * are not frozen and not obviously negatively unbounded.
3260 : * Then we iterate over all marked variables if they can attain
3261 : * any values smaller than zero or at most negative one.
3262 : * If not, we mark the row as being redundant (assuming it hasn't
3263 : * been detected as being obviously redundant in the mean time).
3264 : */
3265 268320834 : int isl_tab_detect_redundant(struct isl_tab *tab)
3266 : {
3267 : int i;
3268 : unsigned n_marked;
3269 :
3270 268320834 : if (!tab)
3271 0 : return -1;
3272 268320834 : if (tab->empty)
3273 407486 : return 0;
3274 267913348 : if (tab->n_redundant == tab->n_row)
3275 742824 : return 0;
3276 :
3277 267170524 : n_marked = 0;
3278 2045810098 : for (i = tab->n_redundant; i < tab->n_row; ++i) {
3279 1778639574 : struct isl_tab_var *var = isl_tab_var_from_row(tab, i);
3280 1778639574 : var->marked = !var->frozen && var->is_nonneg;
3281 1778639574 : if (var->marked)
3282 451637351 : n_marked++;
3283 : }
3284 1657305695 : for (i = tab->n_dead; i < tab->n_col; ++i) {
3285 1390135171 : struct isl_tab_var *var = var_from_col(tab, i);
3286 2053417832 : var->marked = !var->frozen && var->is_nonneg &&
3287 663282661 : !min_is_manifestly_unbounded(tab, var);
3288 1390135171 : if (var->marked)
3289 368312995 : n_marked++;
3290 : }
3291 1239721582 : while (n_marked) {
3292 : struct isl_tab_var *var;
3293 : int red;
3294 730427907 : var = select_marked(tab);
3295 730427907 : if (!var)
3296 25047373 : break;
3297 705380534 : var->marked = 0;
3298 705380534 : n_marked--;
3299 705380534 : red = con_is_redundant(tab, var);
3300 705380534 : if (red < 0)
3301 0 : return -1;
3302 705380534 : if (red && !var->is_redundant)
3303 12871301 : if (isl_tab_mark_redundant(tab, var->index) < 0)
3304 0 : return -1;
3305 4624080226 : for (i = tab->n_dead; i < tab->n_col; ++i) {
3306 3918699692 : var = var_from_col(tab, i);
3307 3918699692 : if (!var->marked)
3308 3077305397 : continue;
3309 841394295 : if (!min_is_manifestly_unbounded(tab, var))
3310 760701588 : continue;
3311 80692707 : var->marked = 0;
3312 80692707 : n_marked--;
3313 : }
3314 : }
3315 :
3316 267170524 : return 0;
3317 : }
3318 :
3319 14099801089 : int isl_tab_is_equality(struct isl_tab *tab, int con)
3320 : {
3321 : int row;
3322 : unsigned off;
3323 :
3324 14099801089 : if (!tab)
3325 0 : return -1;
3326 14099801089 : if (tab->con[con].is_zero)
3327 9935425374 : return 1;
3328 4164375715 : if (tab->con[con].is_redundant)
3329 1700088571 : return 0;
3330 2464287144 : if (!tab->con[con].is_row)
3331 1972671980 : return tab->con[con].index < tab->n_dead;
3332 :
3333 491615164 : row = tab->con[con].index;
3334 :
3335 491615164 : off = 2 + tab->M;
3336 1763258380 : return isl_int_is_zero(tab->mat->row[row][1]) &&
3337 1068809148 : !row_is_big(tab, row) &&
3338 288596992 : isl_seq_first_non_zero(tab->mat->row[row] + off + tab->n_dead,
3339 288596992 : tab->n_col - tab->n_dead) == -1;
3340 : }
3341 :
3342 : /* Return the minimal value of the affine expression "f" with denominator
3343 : * "denom" in *opt, *opt_denom, assuming the tableau is not empty and
3344 : * the expression cannot attain arbitrarily small values.
3345 : * If opt_denom is NULL, then *opt is rounded up to the nearest integer.
3346 : * The return value reflects the nature of the result (empty, unbounded,
3347 : * minimal value returned in *opt).
3348 : *
3349 : * This function assumes that at least one more row and at least
3350 : * one more element in the constraint array are available in the tableau.
3351 : */
3352 729048 : enum isl_lp_result isl_tab_min(struct isl_tab *tab,
3353 : isl_int *f, isl_int denom, isl_int *opt, isl_int *opt_denom,
3354 : unsigned flags)
3355 : {
3356 : int r;
3357 729048 : enum isl_lp_result res = isl_lp_ok;
3358 : struct isl_tab_var *var;
3359 : struct isl_tab_undo *snap;
3360 :
3361 729048 : if (!tab)
3362 0 : return isl_lp_error;
3363 :
3364 729048 : if (tab->empty)
3365 0 : return isl_lp_empty;
3366 :
3367 729048 : snap = isl_tab_snap(tab);
3368 729048 : r = isl_tab_add_row(tab, f);
3369 729048 : if (r < 0)
3370 0 : return isl_lp_error;
3371 729048 : var = &tab->con[r];
3372 : for (;;) {
3373 : int row, col;
3374 1150085 : find_pivot(tab, var, var, -1, &row, &col);
3375 1150085 : if (row == var->index) {
3376 150868 : res = isl_lp_unbounded;
3377 879916 : break;
3378 : }
3379 999217 : if (row == -1)
3380 578180 : break;
3381 421037 : if (isl_tab_pivot(tab, row, col) < 0)
3382 0 : return isl_lp_error;
3383 421037 : }
3384 729048 : isl_int_mul(tab->mat->row[var->index][0],
3385 : tab->mat->row[var->index][0], denom);
3386 729048 : if (ISL_FL_ISSET(flags, ISL_TAB_SAVE_DUAL)) {
3387 : int i;
3388 :
3389 907 : isl_vec_free(tab->dual);
3390 907 : tab->dual = isl_vec_alloc(tab->mat->ctx, 1 + tab->n_con);
3391 907 : if (!tab->dual)
3392 0 : return isl_lp_error;
3393 907 : isl_int_set(tab->dual->el[0], tab->mat->row[var->index][0]);
3394 18565 : for (i = 0; i < tab->n_con; ++i) {
3395 : int pos;
3396 17658 : if (tab->con[i].is_row) {
3397 9764 : isl_int_set_si(tab->dual->el[1 + i], 0);
3398 9764 : continue;
3399 : }
3400 7894 : pos = 2 + tab->M + tab->con[i].index;
3401 7894 : if (tab->con[i].negated)
3402 1119 : isl_int_neg(tab->dual->el[1 + i],
3403 : tab->mat->row[var->index][pos]);
3404 : else
3405 6775 : isl_int_set(tab->dual->el[1 + i],
3406 : tab->mat->row[var->index][pos]);
3407 : }
3408 : }
3409 729048 : if (opt && res == isl_lp_ok) {
3410 578180 : if (opt_denom) {
3411 123336 : isl_int_set(*opt, tab->mat->row[var->index][1]);
3412 123336 : isl_int_set(*opt_denom, tab->mat->row[var->index][0]);
3413 : } else
3414 454844 : get_rounded_sample_value(tab, var, 1, opt);
3415 : }
3416 729048 : if (isl_tab_rollback(tab, snap) < 0)
3417 0 : return isl_lp_error;
3418 729048 : return res;
3419 : }
3420 :
3421 : /* Is the constraint at position "con" marked as being redundant?
3422 : * If it is marked as representing an equality, then it is not
3423 : * considered to be redundant.
3424 : * Note that isl_tab_mark_redundant marks both the isl_tab_var as
3425 : * redundant and moves the corresponding row into the first
3426 : * tab->n_redundant positions (or removes the row, assigning it index -1),
3427 : * so the final test is actually redundant itself.
3428 : */
3429 4961762525 : int isl_tab_is_redundant(struct isl_tab *tab, int con)
3430 : {
3431 4961762525 : if (isl_tab_check_con(tab, con) < 0)
3432 0 : return -1;
3433 4961762525 : if (tab->con[con].is_zero)
3434 0 : return 0;
3435 4961762525 : if (tab->con[con].is_redundant)
3436 2278791661 : return 1;
3437 2682970864 : return tab->con[con].is_row && tab->con[con].index < tab->n_redundant;
3438 : }
3439 :
3440 : /* Is variable "var" of "tab" fixed to a constant value by its row
3441 : * in the tableau?
3442 : * If so and if "value" is not NULL, then store this constant value
3443 : * in "value".
3444 : *
3445 : * That is, is it a row variable that only has non-zero coefficients
3446 : * for dead columns?
3447 : */
3448 3053126 : static isl_bool is_constant(struct isl_tab *tab, struct isl_tab_var *var,
3449 : isl_int *value)
3450 : {
3451 3053126 : unsigned off = 2 + tab->M;
3452 3053126 : isl_mat *mat = tab->mat;
3453 : int n;
3454 : int row;
3455 : int pos;
3456 :
3457 3053126 : if (!var->is_row)
3458 435315 : return isl_bool_false;
3459 2617811 : row = var->index;
3460 2617811 : if (row_is_big(tab, row))
3461 0 : return isl_bool_false;
3462 2617811 : n = tab->n_col - tab->n_dead;
3463 2617811 : pos = isl_seq_first_non_zero(mat->row[row] + off + tab->n_dead, n);
3464 2617811 : if (pos != -1)
3465 2599561 : return isl_bool_false;
3466 18250 : if (value)
3467 0 : isl_int_divexact(*value, mat->row[row][1], mat->row[row][0]);
3468 18250 : return isl_bool_true;
3469 : }
3470 :
3471 : /* Has the variable "var' of "tab" reached a value that is greater than
3472 : * or equal (if sgn > 0) or smaller than or equal (if sgn < 0) to "target"?
3473 : * "tmp" has been initialized by the caller and can be used
3474 : * to perform local computations.
3475 : *
3476 : * If the sample value involves the big parameter, then any value
3477 : * is reached.
3478 : * Otherwise check if n/d >= t, i.e., n >= d * t (if sgn > 0)
3479 : * or n/d <= t, i.e., n <= d * t (if sgn < 0).
3480 : */
3481 4171019 : static int reached(struct isl_tab *tab, struct isl_tab_var *var, int sgn,
3482 : isl_int target, isl_int *tmp)
3483 : {
3484 4171019 : if (row_is_big(tab, var->index))
3485 0 : return 1;
3486 4171019 : isl_int_mul(*tmp, tab->mat->row[var->index][0], target);
3487 4171019 : if (sgn > 0)
3488 909405 : return isl_int_ge(tab->mat->row[var->index][1], *tmp);
3489 : else
3490 3261614 : return isl_int_le(tab->mat->row[var->index][1], *tmp);
3491 : }
3492 :
3493 : /* Can variable "var" of "tab" attain the value "target" by
3494 : * pivoting up (if sgn > 0) or down (if sgn < 0)?
3495 : * If not, then pivot up [down] to the greatest [smallest]
3496 : * rational value.
3497 : * "tmp" has been initialized by the caller and can be used
3498 : * to perform local computations.
3499 : *
3500 : * If the variable is manifestly unbounded in the desired direction,
3501 : * then it can attain any value.
3502 : * Otherwise, it can be moved to a row.
3503 : * Continue pivoting until the target is reached.
3504 : * If no more pivoting can be performed, the maximal [minimal]
3505 : * rational value has been reached and the target cannot be reached.
3506 : * If the variable would be pivoted into a manifestly unbounded column,
3507 : * then the target can be reached.
3508 : */
3509 3710019 : static isl_bool var_reaches(struct isl_tab *tab, struct isl_tab_var *var,
3510 : int sgn, isl_int target, isl_int *tmp)
3511 : {
3512 : int row, col;
3513 :
3514 3710019 : if (sgn < 0 && min_is_manifestly_unbounded(tab, var))
3515 425256 : return isl_bool_true;
3516 3284763 : if (sgn > 0 && max_is_manifestly_unbounded(tab, var))
3517 0 : return isl_bool_true;
3518 3284763 : if (to_row(tab, var, sgn) < 0)
3519 0 : return isl_bool_error;
3520 7455782 : while (!reached(tab, var, sgn, target, tmp)) {
3521 3824638 : find_pivot(tab, var, var, sgn, &row, &col);
3522 3824638 : if (row == -1)
3523 675156 : return isl_bool_false;
3524 3149482 : if (row == var->index)
3525 2263226 : return isl_bool_true;
3526 886256 : if (isl_tab_pivot(tab, row, col) < 0)
3527 0 : return isl_bool_error;
3528 : }
3529 :
3530 346381 : return isl_bool_true;
3531 : }
3532 :
3533 : /* Check if variable "var" of "tab" can only attain a single (integer)
3534 : * value, and, if so, add an equality constraint to fix the variable
3535 : * to this single value and store the result in "target".
3536 : * "target" and "tmp" have been initialized by the caller.
3537 : *
3538 : * Given the current sample value, round it down and check
3539 : * whether it is possible to attain a strictly smaller integer value.
3540 : * If so, the variable is not restricted to a single integer value.
3541 : * Otherwise, the search stops at the smallest rational value.
3542 : * Round up this value and check whether it is possible to attain
3543 : * a strictly greater integer value.
3544 : * If so, the variable is not restricted to a single integer value.
3545 : * Otherwise, the search stops at the greatest rational value.
3546 : * If rounding down this value yields a value that is different
3547 : * from rounding up the smallest rational value, then the variable
3548 : * cannot attain any integer value. Mark the tableau empty.
3549 : * Otherwise, add an equality constraint that fixes the variable
3550 : * to the single integer value found.
3551 : */
3552 3034876 : static isl_bool detect_constant_with_tmp(struct isl_tab *tab,
3553 : struct isl_tab_var *var, isl_int *target, isl_int *tmp)
3554 : {
3555 : isl_bool reached;
3556 : isl_vec *eq;
3557 : int pos;
3558 : isl_stat r;
3559 :
3560 3034876 : get_rounded_sample_value(tab, var, -1, target);
3561 3034876 : isl_int_sub_ui(*target, *target, 1);
3562 3034876 : reached = var_reaches(tab, var, -1, *target, tmp);
3563 3034876 : if (reached < 0 || reached)
3564 2359733 : return isl_bool_not(reached);
3565 675143 : get_rounded_sample_value(tab, var, 1, target);
3566 675143 : isl_int_add_ui(*target, *target, 1);
3567 675143 : reached = var_reaches(tab, var, 1, *target, tmp);
3568 675143 : if (reached < 0 || reached)
3569 675130 : return isl_bool_not(reached);
3570 13 : get_rounded_sample_value(tab, var, -1, tmp);
3571 13 : isl_int_sub_ui(*target, *target, 1);
3572 13 : if (isl_int_ne(*target, *tmp)) {
3573 0 : if (isl_tab_mark_empty(tab) < 0)
3574 0 : return isl_bool_error;
3575 0 : return isl_bool_false;
3576 : }
3577 :
3578 13 : if (isl_tab_extend_cons(tab, 1) < 0)
3579 0 : return isl_bool_error;
3580 13 : eq = isl_vec_alloc(isl_tab_get_ctx(tab), 1 + tab->n_var);
3581 13 : if (!eq)
3582 0 : return isl_bool_error;
3583 13 : pos = var - tab->var;
3584 13 : isl_seq_clr(eq->el + 1, tab->n_var);
3585 13 : isl_int_set_si(eq->el[1 + pos], -1);
3586 13 : isl_int_set(eq->el[0], *target);
3587 13 : r = isl_tab_add_eq(tab, eq->el);
3588 13 : isl_vec_free(eq);
3589 :
3590 13 : return r < 0 ? isl_bool_error : isl_bool_true;
3591 : }
3592 :
3593 : /* Check if variable "var" of "tab" can only attain a single (integer)
3594 : * value, and, if so, add an equality constraint to fix the variable
3595 : * to this single value and store the result in "value" (if "value"
3596 : * is not NULL).
3597 : *
3598 : * If the current sample value involves the big parameter,
3599 : * then the variable cannot have a fixed integer value.
3600 : * If the variable is already fixed to a single value by its row, then
3601 : * there is no need to add another equality constraint.
3602 : *
3603 : * Otherwise, allocate some temporary variables and continue
3604 : * with detect_constant_with_tmp.
3605 : */
3606 3053126 : static isl_bool get_constant(struct isl_tab *tab, struct isl_tab_var *var,
3607 : isl_int *value)
3608 : {
3609 : isl_int target, tmp;
3610 : isl_bool is_cst;
3611 :
3612 3053126 : if (var->is_row && row_is_big(tab, var->index))
3613 0 : return isl_bool_false;
3614 3053126 : is_cst = is_constant(tab, var, value);
3615 3053126 : if (is_cst < 0 || is_cst)
3616 18250 : return is_cst;
3617 :
3618 3034876 : if (!value)
3619 3034876 : isl_int_init(target);
3620 3034876 : isl_int_init(tmp);
3621 :
3622 3034876 : is_cst = detect_constant_with_tmp(tab, var,
3623 : value ? value : &target, &tmp);
3624 :
3625 3034876 : isl_int_clear(tmp);
3626 3034876 : if (!value)
3627 3034876 : isl_int_clear(target);
3628 :
3629 3034876 : return is_cst;
3630 : }
3631 :
3632 : /* Check if variable "var" of "tab" can only attain a single (integer)
3633 : * value, and, if so, add an equality constraint to fix the variable
3634 : * to this single value and store the result in "value" (if "value"
3635 : * is not NULL).
3636 : *
3637 : * For rational tableaus, nothing needs to be done.
3638 : */
3639 0 : isl_bool isl_tab_is_constant(struct isl_tab *tab, int var, isl_int *value)
3640 : {
3641 0 : if (!tab)
3642 0 : return isl_bool_error;
3643 0 : if (var < 0 || var >= tab->n_var)
3644 0 : isl_die(isl_tab_get_ctx(tab), isl_error_invalid,
3645 : "position out of bounds", return isl_bool_error);
3646 0 : if (tab->rational)
3647 0 : return isl_bool_false;
3648 :
3649 0 : return get_constant(tab, &tab->var[var], value);
3650 : }
3651 :
3652 : /* Check if any of the variables of "tab" can only attain a single (integer)
3653 : * value, and, if so, add equality constraints to fix those variables
3654 : * to these single values.
3655 : *
3656 : * For rational tableaus, nothing needs to be done.
3657 : */
3658 518583 : isl_stat isl_tab_detect_constants(struct isl_tab *tab)
3659 : {
3660 : int i;
3661 :
3662 518583 : if (!tab)
3663 0 : return isl_stat_error;
3664 518583 : if (tab->rational)
3665 0 : return isl_stat_ok;
3666 :
3667 3571709 : for (i = 0; i < tab->n_var; ++i) {
3668 3053126 : if (get_constant(tab, &tab->var[i], NULL) < 0)
3669 0 : return isl_stat_error;
3670 : }
3671 :
3672 518583 : return isl_stat_ok;
3673 : }
3674 :
3675 : /* Take a snapshot of the tableau that can be restored by a call to
3676 : * isl_tab_rollback.
3677 : */
3678 12560611772 : struct isl_tab_undo *isl_tab_snap(struct isl_tab *tab)
3679 : {
3680 12560611772 : if (!tab)
3681 0 : return NULL;
3682 12560611772 : tab->need_undo = 1;
3683 12560611772 : return tab->top;
3684 : }
3685 :
3686 : /* Does "tab" need to keep track of undo information?
3687 : * That is, was a snapshot taken that may need to be restored?
3688 : */
3689 0 : isl_bool isl_tab_need_undo(struct isl_tab *tab)
3690 : {
3691 0 : if (!tab)
3692 0 : return isl_bool_error;
3693 :
3694 0 : return tab->need_undo;
3695 : }
3696 :
3697 : /* Remove all tracking of undo information from "tab", invalidating
3698 : * any snapshots that may have been taken of the tableau.
3699 : * Since all snapshots have been invalidated, there is also
3700 : * no need to start keeping track of undo information again.
3701 : */
3702 0 : void isl_tab_clear_undo(struct isl_tab *tab)
3703 : {
3704 0 : if (!tab)
3705 0 : return;
3706 :
3707 0 : free_undo(tab);
3708 0 : tab->need_undo = 0;
3709 : }
3710 :
3711 : /* Undo the operation performed by isl_tab_relax.
3712 : */
3713 : static isl_stat unrelax(struct isl_tab *tab, struct isl_tab_var *var)
3714 : WARN_UNUSED;
3715 69586 : static isl_stat unrelax(struct isl_tab *tab, struct isl_tab_var *var)
3716 : {
3717 69586 : unsigned off = 2 + tab->M;
3718 :
3719 69586 : if (!var->is_row && !max_is_manifestly_unbounded(tab, var))
3720 19051 : if (to_row(tab, var, 1) < 0)
3721 0 : return isl_stat_error;
3722 :
3723 69586 : if (var->is_row) {
3724 23442 : isl_int_sub(tab->mat->row[var->index][1],
3725 : tab->mat->row[var->index][1], tab->mat->row[var->index][0]);
3726 23442 : if (var->is_nonneg) {
3727 23442 : int sgn = restore_row(tab, var);
3728 23442 : isl_assert(tab->mat->ctx, sgn >= 0,
3729 : return isl_stat_error);
3730 : }
3731 : } else {
3732 : int i;
3733 :
3734 276676 : for (i = 0; i < tab->n_row; ++i) {
3735 230532 : if (isl_int_is_zero(tab->mat->row[i][off + var->index]))
3736 102816 : continue;
3737 127716 : isl_int_add(tab->mat->row[i][1], tab->mat->row[i][1],
3738 : tab->mat->row[i][off + var->index]);
3739 : }
3740 :
3741 : }
3742 :
3743 69586 : return isl_stat_ok;
3744 : }
3745 :
3746 : /* Undo the operation performed by isl_tab_unrestrict.
3747 : *
3748 : * In particular, mark the variable as being non-negative and make
3749 : * sure the sample value respects this constraint.
3750 : */
3751 517606 : static isl_stat ununrestrict(struct isl_tab *tab, struct isl_tab_var *var)
3752 : {
3753 517606 : var->is_nonneg = 1;
3754 :
3755 517606 : if (var->is_row && restore_row(tab, var) < -1)
3756 0 : return isl_stat_error;
3757 :
3758 517606 : return isl_stat_ok;
3759 : }
3760 :
3761 : /* Unmark the last redundant row in "tab" as being redundant.
3762 : * This undoes part of the modifications performed by isl_tab_mark_redundant.
3763 : * In particular, remove the redundant mark and make
3764 : * sure the sample value respects the constraint again.
3765 : * A variable that is marked non-negative by isl_tab_mark_redundant
3766 : * is covered by a separate undo record.
3767 : */
3768 7825121380 : static isl_stat restore_last_redundant(struct isl_tab *tab)
3769 : {
3770 : struct isl_tab_var *var;
3771 :
3772 7825121380 : if (tab->n_redundant < 1)
3773 0 : isl_die(isl_tab_get_ctx(tab), isl_error_internal,
3774 : "no redundant rows", return isl_stat_error);
3775 :
3776 7825121380 : var = isl_tab_var_from_row(tab, tab->n_redundant - 1);
3777 7825121380 : var->is_redundant = 0;
3778 7825121380 : tab->n_redundant--;
3779 7825121380 : restore_row(tab, var);
3780 :
3781 7825121380 : return isl_stat_ok;
3782 : }
3783 :
3784 : static isl_stat perform_undo_var(struct isl_tab *tab, struct isl_tab_undo *undo)
3785 : WARN_UNUSED;
3786 >10205*10^7 : static isl_stat perform_undo_var(struct isl_tab *tab, struct isl_tab_undo *undo)
3787 : {
3788 >10205*10^7 : struct isl_tab_var *var = var_from_index(tab, undo->u.var_index);
3789 >10205*10^7 : switch (undo->type) {
3790 : case isl_tab_undo_nonneg:
3791 17755683853 : var->is_nonneg = 0;
3792 17755683853 : break;
3793 : case isl_tab_undo_redundant:
3794 7573721396 : if (!var->is_row || var->index != tab->n_redundant - 1)
3795 0 : isl_die(isl_tab_get_ctx(tab), isl_error_internal,
3796 : "not undoing last redundant row",
3797 : return isl_stat_error);
3798 7573721396 : return restore_last_redundant(tab);
3799 : case isl_tab_undo_freeze:
3800 59769817461 : var->frozen = 0;
3801 59769817461 : break;
3802 : case isl_tab_undo_zero:
3803 317402 : var->is_zero = 0;
3804 317402 : if (!var->is_row)
3805 317402 : tab->n_dead--;
3806 317402 : break;
3807 : case isl_tab_undo_allocate:
3808 16952900693 : if (undo->u.var_index >= 0) {
3809 0 : isl_assert(tab->mat->ctx, !var->is_row,
3810 : return isl_stat_error);
3811 0 : return drop_col(tab, var->index);
3812 : }
3813 16952900693 : if (!var->is_row) {
3814 1665249288 : if (!max_is_manifestly_unbounded(tab, var)) {
3815 712542479 : if (to_row(tab, var, 1) < 0)
3816 0 : return isl_stat_error;
3817 952706809 : } else if (!min_is_manifestly_unbounded(tab, var)) {
3818 761445122 : if (to_row(tab, var, -1) < 0)
3819 0 : return isl_stat_error;
3820 : } else
3821 191261687 : if (to_row(tab, var, 0) < 0)
3822 0 : return isl_stat_error;
3823 : }
3824 16952900693 : return drop_row(tab, var->index);
3825 : case isl_tab_undo_relax:
3826 69586 : return unrelax(tab, var);
3827 : case isl_tab_undo_unrestrict:
3828 517606 : return ununrestrict(tab, var);
3829 : default:
3830 0 : isl_die(tab->mat->ctx, isl_error_internal,
3831 : "perform_undo_var called on invalid undo record",
3832 : return isl_stat_error);
3833 : }
3834 :
3835 77525818716 : return isl_stat_ok;
3836 : }
3837 :
3838 : /* Restore all rows that have been marked redundant by isl_tab_mark_redundant
3839 : * and that have been preserved in the tableau.
3840 : * Note that isl_tab_mark_redundant may also have marked some variables
3841 : * as being non-negative before marking them redundant. These need
3842 : * to be removed as well as otherwise some constraints could end up
3843 : * getting marked redundant with respect to the variable.
3844 : */
3845 157528804 : isl_stat isl_tab_restore_redundant(struct isl_tab *tab)
3846 : {
3847 157528804 : if (!tab)
3848 0 : return isl_stat_error;
3849 :
3850 157528804 : if (tab->need_undo)
3851 0 : isl_die(isl_tab_get_ctx(tab), isl_error_invalid,
3852 : "manually restoring redundant constraints "
3853 : "interferes with undo history",
3854 : return isl_stat_error);
3855 :
3856 566457592 : while (tab->n_redundant > 0) {
3857 251399984 : if (tab->row_var[tab->n_redundant - 1] >= 0) {
3858 : struct isl_tab_var *var;
3859 :
3860 217675905 : var = isl_tab_var_from_row(tab, tab->n_redundant - 1);
3861 217675905 : var->is_nonneg = 0;
3862 : }
3863 251399984 : restore_last_redundant(tab);
3864 : }
3865 157528804 : return isl_stat_ok;
3866 : }
3867 :
3868 : /* Undo the addition of an integer division to the basic map representation
3869 : * of "tab" in position "pos".
3870 : */
3871 0 : static isl_stat drop_bmap_div(struct isl_tab *tab, int pos)
3872 : {
3873 : int off;
3874 :
3875 0 : off = tab->n_var - isl_basic_map_dim(tab->bmap, isl_dim_div);
3876 0 : if (isl_basic_map_drop_div(tab->bmap, pos - off) < 0)
3877 : return isl_stat_error;
3878 0 : if (tab->samples) {
3879 0 : tab->samples = isl_mat_drop_cols(tab->samples, 1 + pos, 1);
3880 0 : if (!tab->samples)
3881 0 : return isl_stat_error;
3882 : }
3883 :
3884 0 : return isl_stat_ok;
3885 : }
3886 :
3887 : /* Restore the tableau to the state where the basic variables
3888 : * are those in "col_var".
3889 : * We first construct a list of variables that are currently in
3890 : * the basis, but shouldn't. Then we iterate over all variables
3891 : * that should be in the basis and for each one that is currently
3892 : * not in the basis, we exchange it with one of the elements of the
3893 : * list constructed before.
3894 : * We can always find an appropriate variable to pivot with because
3895 : * the current basis is mapped to the old basis by a non-singular
3896 : * matrix and so we can never end up with a zero row.
3897 : */
3898 0 : static int restore_basis(struct isl_tab *tab, int *col_var)
3899 : {
3900 : int i, j;
3901 0 : int n_extra = 0;
3902 0 : int *extra = NULL; /* current columns that contain bad stuff */
3903 0 : unsigned off = 2 + tab->M;
3904 :
3905 0 : extra = isl_alloc_array(tab->mat->ctx, int, tab->n_col);
3906 0 : if (tab->n_col && !extra)
3907 0 : goto error;
3908 0 : for (i = 0; i < tab->n_col; ++i) {
3909 0 : for (j = 0; j < tab->n_col; ++j)
3910 0 : if (tab->col_var[i] == col_var[j])
3911 0 : break;
3912 0 : if (j < tab->n_col)
3913 0 : continue;
3914 0 : extra[n_extra++] = i;
3915 : }
3916 0 : for (i = 0; i < tab->n_col && n_extra > 0; ++i) {
3917 : struct isl_tab_var *var;
3918 : int row;
3919 :
3920 0 : for (j = 0; j < tab->n_col; ++j)
3921 0 : if (col_var[i] == tab->col_var[j])
3922 0 : break;
3923 0 : if (j < tab->n_col)
3924 0 : continue;
3925 0 : var = var_from_index(tab, col_var[i]);
3926 0 : row = var->index;
3927 0 : for (j = 0; j < n_extra; ++j)
3928 0 : if (!isl_int_is_zero(tab->mat->row[row][off+extra[j]]))
3929 : break;
3930 0 : isl_assert(tab->mat->ctx, j < n_extra, goto error);
3931 0 : if (isl_tab_pivot(tab, row, extra[j]) < 0)
3932 0 : goto error;
3933 0 : extra[j] = extra[--n_extra];
3934 : }
3935 :
3936 0 : free(extra);
3937 0 : return 0;
3938 : error:
3939 0 : free(extra);
3940 0 : return -1;
3941 : }
3942 :
3943 : /* Remove all samples with index n or greater, i.e., those samples
3944 : * that were added since we saved this number of samples in
3945 : * isl_tab_save_samples.
3946 : */
3947 0 : static void drop_samples_since(struct isl_tab *tab, int n)
3948 : {
3949 : int i;
3950 :
3951 0 : for (i = tab->n_sample - 1; i >= 0 && tab->n_sample > n; --i) {
3952 0 : if (tab->sample_index[i] < n)
3953 0 : continue;
3954 :
3955 0 : if (i != tab->n_sample - 1) {
3956 0 : int t = tab->sample_index[tab->n_sample-1];
3957 0 : tab->sample_index[tab->n_sample-1] = tab->sample_index[i];
3958 0 : tab->sample_index[i] = t;
3959 0 : isl_mat_swap_rows(tab->samples, tab->n_sample-1, i);
3960 : }
3961 0 : tab->n_sample--;
3962 : }
3963 0 : }
3964 :
3965 : static isl_stat perform_undo(struct isl_tab *tab, struct isl_tab_undo *undo)
3966 : WARN_UNUSED;
3967 >12499*10^7 : static isl_stat perform_undo(struct isl_tab *tab, struct isl_tab_undo *undo)
3968 : {
3969 >12499*10^7 : switch (undo->type) {
3970 : case isl_tab_undo_rational:
3971 135741 : tab->rational = 0;
3972 135741 : break;
3973 : case isl_tab_undo_empty:
3974 6088534501 : tab->empty = 0;
3975 6088534501 : break;
3976 : case isl_tab_undo_nonneg:
3977 : case isl_tab_undo_redundant:
3978 : case isl_tab_undo_freeze:
3979 : case isl_tab_undo_zero:
3980 : case isl_tab_undo_allocate:
3981 : case isl_tab_undo_relax:
3982 : case isl_tab_undo_unrestrict:
3983 >10205*10^7 : return perform_undo_var(tab, undo);
3984 : case isl_tab_undo_bmap_eq:
3985 0 : return isl_basic_map_free_equality(tab->bmap, 1);
3986 : case isl_tab_undo_bmap_ineq:
3987 16848830197 : return isl_basic_map_free_inequality(tab->bmap, 1);
3988 : case isl_tab_undo_bmap_div:
3989 0 : return drop_bmap_div(tab, undo->u.var_index);
3990 : case isl_tab_undo_saved_basis:
3991 0 : if (restore_basis(tab, undo->u.col_var) < 0)
3992 0 : return isl_stat_error;
3993 0 : break;
3994 : case isl_tab_undo_drop_sample:
3995 0 : tab->n_outside--;
3996 0 : break;
3997 : case isl_tab_undo_saved_samples:
3998 0 : drop_samples_since(tab, undo->u.n);
3999 0 : break;
4000 : case isl_tab_undo_callback:
4001 0 : return undo->u.callback->run(undo->u.callback);
4002 : default:
4003 0 : isl_assert(tab->mat->ctx, 0, return isl_stat_error);
4004 : }
4005 6088670242 : return isl_stat_ok;
4006 : }
4007 :
4008 : /* Return the tableau to the state it was in when the snapshot "snap"
4009 : * was taken.
4010 : */
4011 6320919892 : int isl_tab_rollback(struct isl_tab *tab, struct isl_tab_undo *snap)
4012 : {
4013 : struct isl_tab_undo *undo, *next;
4014 :
4015 6320919892 : if (!tab)
4016 0 : return -1;
4017 :
4018 6320919892 : tab->in_undo = 1;
4019 >13131*10^7 : for (undo = tab->top; undo && undo != &tab->bottom; undo = next) {
4020 >12965*10^7 : next = undo->next;
4021 >12965*10^7 : if (undo == snap)
4022 4663248056 : break;
4023 >12499*10^7 : if (perform_undo(tab, undo) < 0) {
4024 0 : tab->top = undo;
4025 0 : free_undo(tab);
4026 0 : tab->in_undo = 0;
4027 0 : return -1;
4028 : }
4029 >12499*10^7 : free_undo_record(undo);
4030 : }
4031 6320919892 : tab->in_undo = 0;
4032 6320919892 : tab->top = undo;
4033 6320919892 : if (!undo)
4034 0 : return -1;
4035 6320919892 : return 0;
4036 : }
4037 :
4038 : /* The given row "row" represents an inequality violated by all
4039 : * points in the tableau. Check for some special cases of such
4040 : * separating constraints.
4041 : * In particular, if the row has been reduced to the constant -1,
4042 : * then we know the inequality is adjacent (but opposite) to
4043 : * an equality in the tableau.
4044 : * If the row has been reduced to r = c*(-1 -r'), with r' an inequality
4045 : * of the tableau and c a positive constant, then the inequality
4046 : * is adjacent (but opposite) to the inequality r'.
4047 : */
4048 26721360 : static enum isl_ineq_type separation_type(struct isl_tab *tab, unsigned row)
4049 : {
4050 : int pos;
4051 26721360 : unsigned off = 2 + tab->M;
4052 :
4053 26721360 : if (tab->rational)
4054 5492 : return isl_ineq_separate;
4055 :
4056 26715868 : if (!isl_int_is_one(tab->mat->row[row][0]))
4057 2009799 : return isl_ineq_separate;
4058 :
4059 24706069 : pos = isl_seq_first_non_zero(tab->mat->row[row] + off + tab->n_dead,
4060 24706069 : tab->n_col - tab->n_dead);
4061 24706069 : if (pos == -1) {
4062 1057421 : if (isl_int_is_negone(tab->mat->row[row][1]))
4063 1048416 : return isl_ineq_adj_eq;
4064 : else
4065 9005 : return isl_ineq_separate;
4066 : }
4067 :
4068 23648648 : if (!isl_int_eq(tab->mat->row[row][1],
4069 : tab->mat->row[row][off + tab->n_dead + pos]))
4070 8532343 : return isl_ineq_separate;
4071 :
4072 30232610 : pos = isl_seq_first_non_zero(
4073 15116305 : tab->mat->row[row] + off + tab->n_dead + pos + 1,
4074 15116305 : tab->n_col - tab->n_dead - pos - 1);
4075 :
4076 15116305 : return pos == -1 ? isl_ineq_adj_ineq : isl_ineq_separate;
4077 : }
4078 :
4079 : /* Check the effect of inequality "ineq" on the tableau "tab".
4080 : * The result may be
4081 : * isl_ineq_redundant: satisfied by all points in the tableau
4082 : * isl_ineq_separate: satisfied by no point in the tableau
4083 : * isl_ineq_cut: satisfied by some by not all points
4084 : * isl_ineq_adj_eq: adjacent to an equality
4085 : * isl_ineq_adj_ineq: adjacent to an inequality.
4086 : */
4087 100845874 : enum isl_ineq_type isl_tab_ineq_type(struct isl_tab *tab, isl_int *ineq)
4088 : {
4089 100845874 : enum isl_ineq_type type = isl_ineq_error;
4090 100845874 : struct isl_tab_undo *snap = NULL;
4091 : int con;
4092 : int row;
4093 :
4094 100845874 : if (!tab)
4095 0 : return isl_ineq_error;
4096 :
4097 100845874 : if (isl_tab_extend_cons(tab, 1) < 0)
4098 0 : return isl_ineq_error;
4099 :
4100 100845874 : snap = isl_tab_snap(tab);
4101 :
4102 100845874 : con = isl_tab_add_row(tab, ineq);
4103 100845874 : if (con < 0)
4104 0 : goto error;
4105 :
4106 100845874 : row = tab->con[con].index;
4107 100845874 : if (isl_tab_row_is_redundant(tab, row))
4108 0 : type = isl_ineq_redundant;
4109 141047034 : else if (isl_int_is_neg(tab->mat->row[row][1]) &&
4110 80321384 : (tab->rational ||
4111 40120224 : isl_int_abs_ge(tab->mat->row[row][1],
4112 38629715 : tab->mat->row[row][0]))) {
4113 38629715 : int nonneg = at_least_zero(tab, &tab->con[con]);
4114 38629715 : if (nonneg < 0)
4115 0 : goto error;
4116 38629715 : if (nonneg)
4117 11908355 : type = isl_ineq_cut;
4118 : else
4119 26721360 : type = separation_type(tab, row);
4120 : } else {
4121 62216159 : int red = con_is_redundant(tab, &tab->con[con]);
4122 62216159 : if (red < 0)
4123 0 : goto error;
4124 62216159 : if (!red)
4125 13557804 : type = isl_ineq_cut;
4126 : else
4127 48658355 : type = isl_ineq_redundant;
4128 : }
4129 :
4130 100845874 : if (isl_tab_rollback(tab, snap))
4131 0 : return isl_ineq_error;
4132 100845874 : return type;
4133 : error:
4134 0 : return isl_ineq_error;
4135 : }
4136 :
4137 1328749450 : isl_stat isl_tab_track_bmap(struct isl_tab *tab, __isl_take isl_basic_map *bmap)
4138 : {
4139 1328749450 : bmap = isl_basic_map_cow(bmap);
4140 1328749450 : if (!tab || !bmap)
4141 : goto error;
4142 :
4143 1328749450 : if (tab->empty) {
4144 1296370751 : bmap = isl_basic_map_set_to_empty(bmap);
4145 1296370751 : if (!bmap)
4146 0 : goto error;
4147 1296370751 : tab->bmap = bmap;
4148 1296370751 : return isl_stat_ok;
4149 : }
4150 :
4151 32378699 : isl_assert(tab->mat->ctx, tab->n_eq == bmap->n_eq, goto error);
4152 32378699 : isl_assert(tab->mat->ctx,
4153 : tab->n_con == bmap->n_eq + bmap->n_ineq, goto error);
4154 :
4155 32378699 : tab->bmap = bmap;
4156 :
4157 32378699 : return isl_stat_ok;
4158 : error:
4159 0 : isl_basic_map_free(bmap);
4160 0 : return isl_stat_error;
4161 : }
4162 :
4163 0 : isl_stat isl_tab_track_bset(struct isl_tab *tab, __isl_take isl_basic_set *bset)
4164 : {
4165 0 : return isl_tab_track_bmap(tab, bset_to_bmap(bset));
4166 : }
4167 :
4168 4955 : __isl_keep isl_basic_set *isl_tab_peek_bset(struct isl_tab *tab)
4169 : {
4170 4955 : if (!tab)
4171 0 : return NULL;
4172 :
4173 4955 : return bset_from_bmap(tab->bmap);
4174 : }
4175 :
4176 0 : static void isl_tab_print_internal(__isl_keep struct isl_tab *tab,
4177 : FILE *out, int indent)
4178 : {
4179 : unsigned r, c;
4180 : int i;
4181 :
4182 0 : if (!tab) {
4183 0 : fprintf(out, "%*snull tab\n", indent, "");
4184 0 : return;
4185 : }
4186 0 : fprintf(out, "%*sn_redundant: %d, n_dead: %d", indent, "",
4187 : tab->n_redundant, tab->n_dead);
4188 0 : if (tab->rational)
4189 0 : fprintf(out, ", rational");
4190 0 : if (tab->empty)
4191 0 : fprintf(out, ", empty");
4192 0 : fprintf(out, "\n");
4193 0 : fprintf(out, "%*s[", indent, "");
4194 0 : for (i = 0; i < tab->n_var; ++i) {
4195 0 : if (i)
4196 0 : fprintf(out, (i == tab->n_param ||
4197 0 : i == tab->n_var - tab->n_div) ? "; "
4198 : : ", ");
4199 0 : fprintf(out, "%c%d%s", tab->var[i].is_row ? 'r' : 'c',
4200 0 : tab->var[i].index,
4201 0 : tab->var[i].is_zero ? " [=0]" :
4202 0 : tab->var[i].is_redundant ? " [R]" : "");
4203 : }
4204 0 : fprintf(out, "]\n");
4205 0 : fprintf(out, "%*s[", indent, "");
4206 0 : for (i = 0; i < tab->n_con; ++i) {
4207 0 : if (i)
4208 0 : fprintf(out, ", ");
4209 0 : fprintf(out, "%c%d%s", tab->con[i].is_row ? 'r' : 'c',
4210 0 : tab->con[i].index,
4211 0 : tab->con[i].is_zero ? " [=0]" :
4212 0 : tab->con[i].is_redundant ? " [R]" : "");
4213 : }
4214 0 : fprintf(out, "]\n");
4215 0 : fprintf(out, "%*s[", indent, "");
4216 0 : for (i = 0; i < tab->n_row; ++i) {
4217 0 : const char *sign = "";
4218 0 : if (i)
4219 0 : fprintf(out, ", ");
4220 0 : if (tab->row_sign) {
4221 0 : if (tab->row_sign[i] == isl_tab_row_unknown)
4222 0 : sign = "?";
4223 0 : else if (tab->row_sign[i] == isl_tab_row_neg)
4224 0 : sign = "-";
4225 0 : else if (tab->row_sign[i] == isl_tab_row_pos)
4226 0 : sign = "+";
4227 : else
4228 0 : sign = "+-";
4229 : }
4230 0 : fprintf(out, "r%d: %d%s%s", i, tab->row_var[i],
4231 0 : isl_tab_var_from_row(tab, i)->is_nonneg ? " [>=0]" : "", sign);
4232 : }
4233 0 : fprintf(out, "]\n");
4234 0 : fprintf(out, "%*s[", indent, "");
4235 0 : for (i = 0; i < tab->n_col; ++i) {
4236 0 : if (i)
4237 0 : fprintf(out, ", ");
4238 0 : fprintf(out, "c%d: %d%s", i, tab->col_var[i],
4239 0 : var_from_col(tab, i)->is_nonneg ? " [>=0]" : "");
4240 : }
4241 0 : fprintf(out, "]\n");
4242 0 : r = tab->mat->n_row;
4243 0 : tab->mat->n_row = tab->n_row;
4244 0 : c = tab->mat->n_col;
4245 0 : tab->mat->n_col = 2 + tab->M + tab->n_col;
4246 0 : isl_mat_print_internal(tab->mat, out, indent);
4247 0 : tab->mat->n_row = r;
4248 0 : tab->mat->n_col = c;
4249 0 : if (tab->bmap)
4250 0 : isl_basic_map_print_internal(tab->bmap, out, indent);
4251 : }
4252 :
4253 0 : void isl_tab_dump(__isl_keep struct isl_tab *tab)
4254 : {
4255 0 : isl_tab_print_internal(tab, stderr, 0);
4256 0 : }
|
