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 3033779522 : 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 3033779522 : unsigned off = 2 + M;
39 :
40 3033779522 : tab = isl_calloc_type(ctx, struct isl_tab);
41 3033779522 : if (!tab)
42 0 : return NULL;
43 3033779522 : tab->mat = isl_mat_alloc(ctx, n_row, off + n_var);
44 3033779522 : if (!tab->mat)
45 0 : goto error;
46 3033779522 : tab->var = isl_alloc_array(ctx, struct isl_tab_var, n_var);
47 3033779522 : if (n_var && !tab->var)
48 0 : goto error;
49 3033779522 : tab->con = isl_alloc_array(ctx, struct isl_tab_var, n_row);
50 3033779522 : if (n_row && !tab->con)
51 0 : goto error;
52 3033779522 : tab->col_var = isl_alloc_array(ctx, int, n_var);
53 3033779522 : if (n_var && !tab->col_var)
54 0 : goto error;
55 3033779522 : tab->row_var = isl_alloc_array(ctx, int, n_row);
56 3033779522 : if (n_row && !tab->row_var)
57 0 : goto error;
58 15682690778 : for (i = 0; i < n_var; ++i) {
59 12648911256 : tab->var[i].index = i;
60 12648911256 : tab->var[i].is_row = 0;
61 12648911256 : tab->var[i].is_nonneg = 0;
62 12648911256 : tab->var[i].is_zero = 0;
63 12648911256 : tab->var[i].is_redundant = 0;
64 12648911256 : tab->var[i].frozen = 0;
65 12648911256 : tab->var[i].negated = 0;
66 12648911256 : tab->col_var[i] = i;
67 : }
68 3033779522 : tab->n_row = 0;
69 3033779522 : tab->n_con = 0;
70 3033779522 : tab->n_eq = 0;
71 3033779522 : tab->max_con = n_row;
72 3033779522 : tab->n_col = n_var;
73 3033779522 : tab->n_var = n_var;
74 3033779522 : tab->max_var = n_var;
75 3033779522 : tab->n_param = 0;
76 3033779522 : tab->n_div = 0;
77 3033779522 : tab->n_dead = 0;
78 3033779522 : tab->n_redundant = 0;
79 3033779522 : tab->strict_redundant = 0;
80 3033779522 : tab->need_undo = 0;
81 3033779522 : tab->rational = 0;
82 3033779522 : tab->empty = 0;
83 3033779522 : tab->in_undo = 0;
84 3033779522 : tab->M = M;
85 3033779522 : tab->cone = 0;
86 3033779522 : tab->bottom.type = isl_tab_undo_bottom;
87 3033779522 : tab->bottom.next = NULL;
88 3033779522 : tab->top = &tab->bottom;
89 :
90 3033779522 : tab->n_zero = 0;
91 3033779522 : tab->n_unbounded = 0;
92 3033779522 : tab->basis = NULL;
93 :
94 3033779522 : return tab;
95 : error:
96 0 : isl_tab_free(tab);
97 0 : return NULL;
98 : }
99 :
100 25468168346 : isl_ctx *isl_tab_get_ctx(struct isl_tab *tab)
101 : {
102 25468168346 : return tab ? isl_mat_get_ctx(tab->mat) : NULL;
103 : }
104 :
105 13142973352 : int isl_tab_extend_cons(struct isl_tab *tab, unsigned n_new)
106 : {
107 : unsigned off;
108 :
109 13142973352 : if (!tab)
110 0 : return -1;
111 :
112 13142973352 : off = 2 + tab->M;
113 :
114 13142973352 : if (tab->max_con < tab->n_con + n_new) {
115 : struct isl_tab_var *con;
116 :
117 85413381 : con = isl_realloc_array(tab->mat->ctx, tab->con,
118 : struct isl_tab_var, tab->max_con + n_new);
119 85413381 : if (!con)
120 0 : return -1;
121 85413381 : tab->con = con;
122 85413381 : tab->max_con += n_new;
123 : }
124 13142973352 : if (tab->mat->n_row < tab->n_row + n_new) {
125 : int *row_var;
126 :
127 356611802 : tab->mat = isl_mat_extend(tab->mat,
128 356611802 : tab->n_row + n_new, off + tab->n_col);
129 178305901 : if (!tab->mat)
130 0 : return -1;
131 178305901 : row_var = isl_realloc_array(tab->mat->ctx, tab->row_var,
132 : int, tab->mat->n_row);
133 178305901 : if (!row_var)
134 0 : return -1;
135 178305901 : tab->row_var = row_var;
136 178305901 : if (tab->row_sign) {
137 : enum isl_tab_row_sign *s;
138 7651081 : s = isl_realloc_array(tab->mat->ctx, tab->row_sign,
139 : enum isl_tab_row_sign, tab->mat->n_row);
140 7651081 : if (!s)
141 0 : return -1;
142 7651081 : tab->row_sign = s;
143 : }
144 : }
145 13142973352 : return 0;
146 : }
147 :
148 : /* Make room for at least n_new extra variables.
149 : * Return -1 if anything went wrong.
150 : */
151 903868586 : int isl_tab_extend_vars(struct isl_tab *tab, unsigned n_new)
152 : {
153 : struct isl_tab_var *var;
154 903868586 : unsigned off = 2 + tab->M;
155 :
156 903868586 : if (tab->max_var < tab->n_var + n_new) {
157 47200712 : var = isl_realloc_array(tab->mat->ctx, tab->var,
158 : struct isl_tab_var, tab->n_var + n_new);
159 47200712 : if (!var)
160 0 : return -1;
161 47200712 : tab->var = var;
162 47200712 : tab->max_var = tab->n_var + n_new;
163 : }
164 :
165 903868586 : if (tab->mat->n_col < off + tab->n_col + n_new) {
166 : int *p;
167 :
168 61606826 : tab->mat = isl_mat_extend(tab->mat,
169 61606826 : tab->mat->n_row, off + tab->n_col + n_new);
170 30803413 : if (!tab->mat)
171 0 : return -1;
172 30803413 : p = isl_realloc_array(tab->mat->ctx, tab->col_var,
173 : int, tab->n_col + n_new);
174 30803413 : if (!p)
175 0 : return -1;
176 30803413 : tab->col_var = p;
177 : }
178 :
179 903868586 : return 0;
180 : }
181 :
182 >11595*10^7 : static void free_undo_record(struct isl_tab_undo *undo)
183 : {
184 >11595*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 >11595*10^7 : free(undo);
191 >11595*10^7 : }
192 :
193 3110340287 : static void free_undo(struct isl_tab *tab)
194 : {
195 : struct isl_tab_undo *undo, *next;
196 :
197 4219827190 : for (undo = tab->top; undo && undo != &tab->bottom; undo = next) {
198 1109486903 : next = undo->next;
199 1109486903 : free_undo_record(undo);
200 : }
201 3110340287 : tab->top = undo;
202 3110340287 : }
203 :
204 3116371104 : void isl_tab_free(struct isl_tab *tab)
205 : {
206 3116371104 : if (!tab)
207 9835320 : return;
208 3106535784 : free_undo(tab);
209 3106535784 : isl_mat_free(tab->mat);
210 3106535784 : isl_vec_free(tab->dual);
211 3106535784 : isl_basic_map_free(tab->bmap);
212 3106535784 : free(tab->var);
213 3106535784 : free(tab->con);
214 3106535784 : free(tab->row_var);
215 3106535784 : free(tab->col_var);
216 3106535784 : free(tab->row_sign);
217 3106535784 : isl_mat_free(tab->samples);
218 3106535784 : free(tab->sample_index);
219 3106535784 : isl_mat_free(tab->basis);
220 3106535784 : free(tab);
221 : }
222 :
223 33441720 : struct isl_tab *isl_tab_dup(struct isl_tab *tab)
224 : {
225 : int i;
226 : struct isl_tab *dup;
227 : unsigned off;
228 :
229 33441720 : if (!tab)
230 0 : return NULL;
231 :
232 33441720 : off = 2 + tab->M;
233 33441720 : dup = isl_calloc_type(tab->mat->ctx, struct isl_tab);
234 33441720 : if (!dup)
235 0 : return NULL;
236 33441720 : dup->mat = isl_mat_dup(tab->mat);
237 33441720 : if (!dup->mat)
238 0 : goto error;
239 33441720 : dup->var = isl_alloc_array(tab->mat->ctx, struct isl_tab_var, tab->max_var);
240 33441720 : if (tab->max_var && !dup->var)
241 0 : goto error;
242 200746668 : for (i = 0; i < tab->n_var; ++i)
243 167304948 : dup->var[i] = tab->var[i];
244 33441720 : dup->con = isl_alloc_array(tab->mat->ctx, struct isl_tab_var, tab->max_con);
245 33441720 : if (tab->max_con && !dup->con)
246 0 : goto error;
247 354068292 : for (i = 0; i < tab->n_con; ++i)
248 320626572 : dup->con[i] = tab->con[i];
249 33441720 : dup->col_var = isl_alloc_array(tab->mat->ctx, int, tab->mat->n_col - off);
250 33441720 : if ((tab->mat->n_col - off) && !dup->col_var)
251 0 : goto error;
252 200308015 : for (i = 0; i < tab->n_col; ++i)
253 166866295 : dup->col_var[i] = tab->col_var[i];
254 33441720 : dup->row_var = isl_alloc_array(tab->mat->ctx, int, tab->mat->n_row);
255 33441720 : if (tab->mat->n_row && !dup->row_var)
256 0 : goto error;
257 353986817 : for (i = 0; i < tab->n_row; ++i)
258 320545097 : dup->row_var[i] = tab->row_var[i];
259 33441720 : if (tab->row_sign) {
260 33441720 : dup->row_sign = isl_alloc_array(tab->mat->ctx, enum isl_tab_row_sign,
261 : tab->mat->n_row);
262 33441720 : if (tab->mat->n_row && !dup->row_sign)
263 0 : goto error;
264 353986817 : for (i = 0; i < tab->n_row; ++i)
265 320545097 : dup->row_sign[i] = tab->row_sign[i];
266 : }
267 33441720 : 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 33441720 : dup->n_row = tab->n_row;
279 33441720 : dup->n_con = tab->n_con;
280 33441720 : dup->n_eq = tab->n_eq;
281 33441720 : dup->max_con = tab->max_con;
282 33441720 : dup->n_col = tab->n_col;
283 33441720 : dup->n_var = tab->n_var;
284 33441720 : dup->max_var = tab->max_var;
285 33441720 : dup->n_param = tab->n_param;
286 33441720 : dup->n_div = tab->n_div;
287 33441720 : dup->n_dead = tab->n_dead;
288 33441720 : dup->n_redundant = tab->n_redundant;
289 33441720 : dup->rational = tab->rational;
290 33441720 : dup->empty = tab->empty;
291 33441720 : dup->strict_redundant = 0;
292 33441720 : dup->need_undo = 0;
293 33441720 : dup->in_undo = 0;
294 33441720 : dup->M = tab->M;
295 33441720 : tab->cone = tab->cone;
296 33441720 : dup->bottom.type = isl_tab_undo_bottom;
297 33441720 : dup->bottom.next = NULL;
298 33441720 : dup->top = &dup->bottom;
299 :
300 33441720 : dup->n_zero = tab->n_zero;
301 33441720 : dup->n_unbounded = tab->n_unbounded;
302 33441720 : dup->basis = isl_mat_dup(tab->basis);
303 :
304 33441720 : 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 39314542 : 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 39314542 : prod = isl_mat_alloc(mat1->ctx, mat1->n_row + mat2->n_row,
333 39314542 : off + col1 + col2);
334 39314542 : if (!prod)
335 0 : return NULL;
336 :
337 39314542 : n = 0;
338 188562650 : for (i = 0; i < r1; ++i) {
339 149248108 : isl_seq_cpy(prod->row[n + i], mat1->row[i], off + d1);
340 149248108 : isl_seq_clr(prod->row[n + i] + off + d1, d2);
341 298496216 : isl_seq_cpy(prod->row[n + i] + off + d1 + d2,
342 149248108 : mat1->row[i] + off + d1, col1 - d1);
343 149248108 : isl_seq_clr(prod->row[n + i] + off + col1 + d1, col2 - d2);
344 : }
345 :
346 39314542 : n += r1;
347 188562650 : for (i = 0; i < r2; ++i) {
348 149248108 : isl_seq_cpy(prod->row[n + i], mat2->row[i], off);
349 149248108 : isl_seq_clr(prod->row[n + i] + off, d1);
350 298496216 : isl_seq_cpy(prod->row[n + i] + off + d1,
351 298496216 : mat2->row[i] + off, d2);
352 149248108 : isl_seq_clr(prod->row[n + i] + off + d1 + d2, col1 - d1);
353 298496216 : isl_seq_cpy(prod->row[n + i] + off + col1 + d1,
354 149248108 : mat2->row[i] + off + d2, col2 - d2);
355 : }
356 :
357 39314542 : n += r2;
358 401789152 : for (i = 0; i < row1 - r1; ++i) {
359 362474610 : isl_seq_cpy(prod->row[n + i], mat1->row[r1 + i], off + d1);
360 362474610 : isl_seq_clr(prod->row[n + i] + off + d1, d2);
361 724949220 : isl_seq_cpy(prod->row[n + i] + off + d1 + d2,
362 362474610 : mat1->row[r1 + i] + off + d1, col1 - d1);
363 362474610 : isl_seq_clr(prod->row[n + i] + off + col1 + d1, col2 - d2);
364 : }
365 :
366 39314542 : n += row1 - r1;
367 401789152 : for (i = 0; i < row2 - r2; ++i) {
368 362474610 : isl_seq_cpy(prod->row[n + i], mat2->row[r2 + i], off);
369 362474610 : isl_seq_clr(prod->row[n + i] + off, d1);
370 724949220 : isl_seq_cpy(prod->row[n + i] + off + d1,
371 724949220 : mat2->row[r2 + i] + off, d2);
372 362474610 : isl_seq_clr(prod->row[n + i] + off + d1 + d2, col1 - d1);
373 724949220 : isl_seq_cpy(prod->row[n + i] + off + col1 + d1,
374 362474610 : mat2->row[r2 + i] + off + d2, col2 - d2);
375 : }
376 :
377 39314542 : 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 708608207 : static void update_index1(struct isl_tab_var *var,
384 : unsigned r1, unsigned r2, unsigned d1, unsigned d2)
385 : {
386 708608207 : if (var->index == -1)
387 2807747 : return;
388 705800460 : if (var->is_row && var->index >= r1)
389 362474610 : var->index += r2;
390 705800460 : if (!var->is_row && var->index >= d1)
391 167099285 : 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 708608207 : 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 708608207 : if (var->index == -1)
402 2807747 : return;
403 705800460 : if (var->is_row) {
404 511722718 : if (var->index < r2)
405 149248108 : var->index += r1;
406 : else
407 362474610 : var->index += row1;
408 : } else {
409 194077742 : if (var->index < d2)
410 26978457 : var->index += d1;
411 : else
412 167099285 : 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 39314542 : 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 39314542 : if (!tab1 || !tab2)
442 0 : return NULL;
443 :
444 39314542 : isl_assert(tab1->mat->ctx, tab1->M == tab2->M, return NULL);
445 39314542 : isl_assert(tab1->mat->ctx, tab1->rational == tab2->rational, return NULL);
446 39314542 : isl_assert(tab1->mat->ctx, tab1->cone == tab2->cone, return NULL);
447 39314542 : isl_assert(tab1->mat->ctx, !tab1->row_sign, return NULL);
448 39314542 : isl_assert(tab1->mat->ctx, !tab2->row_sign, return NULL);
449 39314542 : isl_assert(tab1->mat->ctx, tab1->n_param == 0, return NULL);
450 39314542 : isl_assert(tab1->mat->ctx, tab2->n_param == 0, return NULL);
451 39314542 : isl_assert(tab1->mat->ctx, tab1->n_div == 0, return NULL);
452 39314542 : isl_assert(tab1->mat->ctx, tab2->n_div == 0, return NULL);
453 :
454 39314542 : off = 2 + tab1->M;
455 39314542 : r1 = tab1->n_redundant;
456 39314542 : r2 = tab2->n_redundant;
457 39314542 : d1 = tab1->n_dead;
458 39314542 : d2 = tab2->n_dead;
459 39314542 : prod = isl_calloc_type(tab1->mat->ctx, struct isl_tab);
460 39314542 : if (!prod)
461 0 : return NULL;
462 39314542 : 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 39314542 : if (!prod->mat)
466 0 : goto error;
467 39314542 : prod->var = isl_alloc_array(tab1->mat->ctx, struct isl_tab_var,
468 : tab1->max_var + tab2->max_var);
469 39314542 : if ((tab1->max_var + tab2->max_var) && !prod->var)
470 0 : goto error;
471 236200031 : for (i = 0; i < tab1->n_var; ++i) {
472 196885489 : prod->var[i] = tab1->var[i];
473 196885489 : update_index1(&prod->var[i], r1, r2, d1, d2);
474 : }
475 236200031 : for (i = 0; i < tab2->n_var; ++i) {
476 196885489 : prod->var[tab1->n_var + i] = tab2->var[i];
477 196885489 : update_index2(&prod->var[tab1->n_var + i],
478 : tab1->n_row, tab1->n_col,
479 : r1, r2, d1, d2);
480 : }
481 39314542 : prod->con = isl_alloc_array(tab1->mat->ctx, struct isl_tab_var,
482 : tab1->max_con + tab2->max_con);
483 39314542 : if ((tab1->max_con + tab2->max_con) && !prod->con)
484 0 : goto error;
485 551037260 : for (i = 0; i < tab1->n_con; ++i) {
486 511722718 : prod->con[i] = tab1->con[i];
487 511722718 : update_index1(&prod->con[i], r1, r2, d1, d2);
488 : }
489 551037260 : for (i = 0; i < tab2->n_con; ++i) {
490 511722718 : prod->con[tab1->n_con + i] = tab2->con[i];
491 511722718 : update_index2(&prod->con[tab1->n_con + i],
492 : tab1->n_row, tab1->n_col,
493 : r1, r2, d1, d2);
494 : }
495 39314542 : prod->col_var = isl_alloc_array(tab1->mat->ctx, int,
496 : tab1->n_col + tab2->n_col);
497 39314542 : if ((tab1->n_col + tab2->n_col) && !prod->col_var)
498 0 : goto error;
499 233392284 : for (i = 0; i < tab1->n_col; ++i) {
500 194077742 : int pos = i < d1 ? i : i + d2;
501 194077742 : prod->col_var[pos] = tab1->col_var[i];
502 : }
503 233392284 : for (i = 0; i < tab2->n_col; ++i) {
504 194077742 : int pos = i < d2 ? d1 + i : tab1->n_col + i;
505 194077742 : int t = tab2->col_var[i];
506 194077742 : if (t >= 0)
507 328255 : t += tab1->n_var;
508 : else
509 193749487 : t -= tab1->n_con;
510 194077742 : prod->col_var[pos] = t;
511 : }
512 39314542 : prod->row_var = isl_alloc_array(tab1->mat->ctx, int,
513 : tab1->mat->n_row + tab2->mat->n_row);
514 39314542 : if ((tab1->mat->n_row + tab2->mat->n_row) && !prod->row_var)
515 0 : goto error;
516 551037260 : for (i = 0; i < tab1->n_row; ++i) {
517 511722718 : int pos = i < r1 ? i : i + r2;
518 511722718 : prod->row_var[pos] = tab1->row_var[i];
519 : }
520 551037260 : for (i = 0; i < tab2->n_row; ++i) {
521 511722718 : int pos = i < r2 ? r1 + i : tab1->n_row + i;
522 511722718 : int t = tab2->row_var[i];
523 511722718 : if (t >= 0)
524 196557234 : t += tab1->n_var;
525 : else
526 315165484 : t -= tab1->n_con;
527 511722718 : prod->row_var[pos] = t;
528 : }
529 39314542 : prod->samples = NULL;
530 39314542 : prod->sample_index = NULL;
531 39314542 : prod->n_row = tab1->n_row + tab2->n_row;
532 39314542 : prod->n_con = tab1->n_con + tab2->n_con;
533 39314542 : prod->n_eq = 0;
534 39314542 : prod->max_con = tab1->max_con + tab2->max_con;
535 39314542 : prod->n_col = tab1->n_col + tab2->n_col;
536 39314542 : prod->n_var = tab1->n_var + tab2->n_var;
537 39314542 : prod->max_var = tab1->max_var + tab2->max_var;
538 39314542 : prod->n_param = 0;
539 39314542 : prod->n_div = 0;
540 39314542 : prod->n_dead = tab1->n_dead + tab2->n_dead;
541 39314542 : prod->n_redundant = tab1->n_redundant + tab2->n_redundant;
542 39314542 : prod->rational = tab1->rational;
543 39314542 : prod->empty = tab1->empty || tab2->empty;
544 39314542 : prod->strict_redundant = tab1->strict_redundant || tab2->strict_redundant;
545 39314542 : prod->need_undo = 0;
546 39314542 : prod->in_undo = 0;
547 39314542 : prod->M = tab1->M;
548 39314542 : prod->cone = tab1->cone;
549 39314542 : prod->bottom.type = isl_tab_undo_bottom;
550 39314542 : prod->bottom.next = NULL;
551 39314542 : prod->top = &prod->bottom;
552 :
553 39314542 : prod->n_zero = 0;
554 39314542 : prod->n_unbounded = 0;
555 39314542 : prod->basis = NULL;
556 :
557 39314542 : return prod;
558 : error:
559 0 : isl_tab_free(prod);
560 0 : return NULL;
561 : }
562 :
563 >73297*10^7 : static struct isl_tab_var *var_from_index(struct isl_tab *tab, int i)
564 : {
565 >73297*10^7 : if (i >= 0)
566 >18418*10^7 : return &tab->var[i];
567 : else
568 >54879*10^7 : return &tab->con[~i];
569 : }
570 :
571 >47895*10^7 : struct isl_tab_var *isl_tab_var_from_row(struct isl_tab *tab, int i)
572 : {
573 >47895*10^7 : return var_from_index(tab, tab->row_var[i]);
574 : }
575 :
576 >15864*10^7 : static struct isl_tab_var *var_from_col(struct isl_tab *tab, int i)
577 : {
578 >15864*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 4904963797 : static int max_is_manifestly_unbounded(struct isl_tab *tab,
586 : struct isl_tab_var *var)
587 : {
588 : int i;
589 4904963797 : unsigned off = 2 + tab->M;
590 :
591 4904963797 : if (var->is_row)
592 915194585 : return 0;
593 39412040512 : for (i = tab->n_redundant; i < tab->n_row; ++i) {
594 37442407783 : if (!isl_int_is_neg(tab->mat->row[i][off + var->index]))
595 31941303743 : continue;
596 5501104040 : if (isl_tab_var_from_row(tab, i)->is_nonneg)
597 2020136483 : return 0;
598 : }
599 1969632729 : 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 9358863551 : static int min_is_manifestly_unbounded(struct isl_tab *tab,
607 : struct isl_tab_var *var)
608 : {
609 : int i;
610 9358863551 : unsigned off = 2 + tab->M;
611 :
612 9358863551 : if (var->is_row)
613 6484632874 : return 0;
614 27520312314 : for (i = tab->n_redundant; i < tab->n_row; ++i) {
615 26158979355 : if (!isl_int_is_pos(tab->mat->row[i][off + var->index]))
616 22370600761 : continue;
617 3788378594 : if (isl_tab_var_from_row(tab, i)->is_nonneg)
618 1512897718 : return 0;
619 : }
620 1361332959 : return 1;
621 : }
622 :
623 18815811640 : static int row_cmp(struct isl_tab *tab, int r1, int r2, int c, isl_int *t)
624 : {
625 18815811640 : unsigned off = 2 + tab->M;
626 :
627 18815811640 : 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 18815811640 : isl_int_mul(*t, tab->mat->row[r1][1], tab->mat->row[r2][off + c]);
636 18815811640 : isl_int_submul(*t, tab->mat->row[r2][1], tab->mat->row[r1][off + c]);
637 18815811640 : 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 14982198090 : 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 14982198090 : unsigned off = 2 + tab->M;
666 :
667 14982198090 : isl_int_init(t);
668 14982198090 : r = -1;
669 >21112*10^7 : for (j = tab->n_redundant; j < tab->n_row; ++j) {
670 >19614*10^7 : if (var && j == var->index)
671 11288195834 : continue;
672 >18485*10^7 : if (!isl_tab_var_from_row(tab, j)->is_nonneg)
673 73304085104 : continue;
674 >11155*10^7 : if (sgn * isl_int_sgn(tab->mat->row[j][off + c]) >= 0)
675 79935270398 : continue;
676 31614898985 : if (r < 0) {
677 12799087345 : r = j;
678 12799087345 : continue;
679 : }
680 18815811640 : tsgn = sgn * row_cmp(tab, r, j, c, &t);
681 21476389869 : if (tsgn < 0 || (tsgn == 0 &&
682 2660578229 : tab->row_var[j] < tab->row_var[r]))
683 9093484879 : r = j;
684 : }
685 14982198090 : isl_int_clear(t);
686 14982198090 : 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 27440058229 : 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 27440058229 : *row = *col = -1;
713 :
714 27440058229 : isl_assert(tab->mat->ctx, var->is_row, return);
715 27440058229 : tr = tab->mat->row[var->index] + 2 + tab->M;
716 :
717 27440058229 : c = -1;
718 >13623*10^7 : for (j = tab->n_dead; j < tab->n_col; ++j) {
719 >10879*10^7 : if (isl_int_is_zero(tr[j]))
720 53074923658 : continue;
721 90212957128 : if (isl_int_sgn(tr[j]) != sgn &&
722 34496476227 : var_from_col(tab, j)->is_nonneg)
723 30965402021 : continue;
724 24751078880 : if (c < 0 || tab->col_var[j] < tab->col_var[c])
725 16220577733 : c = j;
726 : }
727 27440058229 : if (c < 0)
728 15174496435 : return;
729 :
730 12265561794 : sgn *= isl_int_sgn(tr[c]);
731 12265561794 : r = pivot_row(tab, skip_var, sgn, c);
732 12265561794 : *row = r < 0 ? var->index : r;
733 12265561794 : *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 >12700*10^7 : int isl_tab_row_is_redundant(struct isl_tab *tab, int row)
743 : {
744 : int i;
745 >12700*10^7 : unsigned off = 2 + tab->M;
746 :
747 >12700*10^7 : if (tab->row_var[row] < 0 && !isl_tab_var_from_row(tab, row)->is_nonneg)
748 12373882911 : return 0;
749 :
750 >11463*10^7 : if (isl_int_is_neg(tab->mat->row[row][1]))
751 24224747729 : return 0;
752 90409798803 : if (tab->strict_redundant && isl_int_is_zero(tab->mat->row[row][1]))
753 0 : return 0;
754 90409798803 : if (tab->M && isl_int_is_neg(tab->mat->row[row][2]))
755 3123469 : return 0;
756 :
757 >32352*10^7 : for (i = tab->n_dead; i < tab->n_col; ++i) {
758 >31158*10^7 : if (isl_int_is_zero(tab->mat->row[row][off + i]))
759 >18303*10^7 : continue;
760 >12854*10^7 : if (tab->col_var[i] >= 0)
761 12829867391 : return 0;
762 >11571*10^7 : if (isl_int_is_neg(tab->mat->row[row][off + i]))
763 58315340132 : return 0;
764 57399707290 : if (!var_from_col(tab, i)->is_nonneg)
765 7317070414 : return 0;
766 : }
767 11944397397 : return 1;
768 : }
769 :
770 17225940846 : 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 17225940846 : t = tab->row_var[row1];
776 17225940846 : tab->row_var[row1] = tab->row_var[row2];
777 17225940846 : tab->row_var[row2] = t;
778 17225940846 : isl_tab_var_from_row(tab, row1)->index = row1;
779 17225940846 : isl_tab_var_from_row(tab, row2)->index = row2;
780 17225940846 : tab->mat = isl_mat_swap_rows(tab->mat, row1, row2);
781 :
782 17225940846 : if (!tab->row_sign)
783 17224912318 : return;
784 1028528 : s = tab->row_sign[row1];
785 1028528 : tab->row_sign[row1] = tab->row_sign[row2];
786 1028528 : 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 >13634*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 >13634*10^7 : if (!tab)
805 0 : return isl_stat_error;
806 >13634*10^7 : if (!tab->need_undo)
807 20382088992 : return isl_stat_ok;
808 :
809 >11595*10^7 : undo = isl_alloc_type(tab->mat->ctx, struct isl_tab_undo);
810 >11595*10^7 : if (!undo)
811 0 : goto error;
812 >11595*10^7 : undo->type = type;
813 >11595*10^7 : undo->u = u;
814 >11595*10^7 : undo->next = tab->top;
815 >11595*10^7 : tab->top = undo;
816 :
817 >11595*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 >11755*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 >11755*10^7 : if (var->is_row)
829 99221380462 : u.var_index = tab->row_var[var->index];
830 : else
831 18329598011 : u.var_index = tab->col_var[var->index];
832 >11755*10^7 : return push_union(tab, type, u);
833 : }
834 :
835 18716236112 : isl_stat isl_tab_push(struct isl_tab *tab, enum isl_tab_undo_type type)
836 : {
837 18716236112 : union isl_tab_undo_val u = { 0 };
838 18716236112 : 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 34100935 : 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 34100935 : u.callback = callback;
862 34100935 : return push_union(tab, isl_tab_undo_callback, u);
863 : }
864 :
865 4860613 : struct isl_tab *isl_tab_init_samples(struct isl_tab *tab)
866 : {
867 4860613 : if (!tab)
868 0 : return NULL;
869 :
870 4860613 : tab->n_sample = 0;
871 4860613 : tab->n_outside = 0;
872 4860613 : tab->samples = isl_mat_alloc(tab->mat->ctx, 1, 1 + tab->n_var);
873 4860613 : if (!tab->samples)
874 0 : goto error;
875 4860613 : tab->sample_index = isl_alloc_array(tab->mat->ctx, int, 1);
876 4860613 : if (!tab->sample_index)
877 0 : goto error;
878 4860613 : return tab;
879 : error:
880 0 : isl_tab_free(tab);
881 0 : return NULL;
882 : }
883 :
884 84741085 : int isl_tab_add_sample(struct isl_tab *tab, __isl_take isl_vec *sample)
885 : {
886 84741085 : if (!tab || !sample)
887 : goto error;
888 :
889 84741085 : if (tab->n_sample + 1 > tab->samples->n_row) {
890 43250788 : int *t = isl_realloc_array(tab->mat->ctx,
891 : tab->sample_index, int, tab->n_sample + 1);
892 43250788 : if (!t)
893 0 : goto error;
894 43250788 : tab->sample_index = t;
895 : }
896 :
897 169482170 : tab->samples = isl_mat_extend(tab->samples,
898 169482170 : tab->n_sample + 1, tab->samples->n_col);
899 84741085 : if (!tab->samples)
900 0 : goto error;
901 :
902 84741085 : isl_seq_cpy(tab->samples->row[tab->n_sample], sample->el, sample->size);
903 84741085 : isl_vec_free(sample);
904 84741085 : tab->sample_index[tab->n_sample] = tab->n_sample;
905 84741085 : tab->n_sample++;
906 :
907 84741085 : return 0;
908 : error:
909 0 : isl_vec_free(sample);
910 0 : return -1;
911 : }
912 :
913 130449173 : struct isl_tab *isl_tab_drop_sample(struct isl_tab *tab, int s)
914 : {
915 130449173 : if (s != tab->n_outside) {
916 111606946 : int t = tab->sample_index[tab->n_outside];
917 111606946 : tab->sample_index[tab->n_outside] = tab->sample_index[s];
918 111606946 : tab->sample_index[s] = t;
919 111606946 : isl_mat_swap_rows(tab->samples, tab->n_outside, s);
920 : }
921 130449173 : tab->n_outside++;
922 130449173 : if (isl_tab_push(tab, isl_tab_undo_drop_sample) < 0) {
923 0 : isl_tab_free(tab);
924 0 : return NULL;
925 : }
926 :
927 130449173 : return tab;
928 : }
929 :
930 : /* Record the current number of samples so that we can remove newer
931 : * samples during a rollback.
932 : */
933 38961548 : isl_stat isl_tab_save_samples(struct isl_tab *tab)
934 : {
935 : union isl_tab_undo_val u;
936 :
937 38961548 : if (!tab)
938 0 : return isl_stat_error;
939 :
940 38961548 : u.n = tab->n_sample;
941 38961548 : 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 12377564340 : int isl_tab_mark_redundant(struct isl_tab *tab, int row)
957 : {
958 12377564340 : struct isl_tab_var *var = isl_tab_var_from_row(tab, row);
959 12377564340 : var->is_redundant = 1;
960 12377564340 : isl_assert(tab->mat->ctx, row >= tab->n_redundant, return -1);
961 12377564340 : if (tab->preserve || tab->need_undo || tab->row_var[row] >= 0) {
962 11376720453 : if (tab->row_var[row] >= 0 && !var->is_nonneg) {
963 3623686146 : var->is_nonneg = 1;
964 3623686146 : if (isl_tab_push_var(tab, isl_tab_undo_nonneg, var) < 0)
965 0 : return -1;
966 : }
967 11376720453 : if (row != tab->n_redundant)
968 9906905270 : swap_rows(tab, row, tab->n_redundant);
969 11376720453 : tab->n_redundant++;
970 11376720453 : return isl_tab_push_var(tab, isl_tab_undo_redundant, var);
971 : } else {
972 1000843887 : if (row != tab->n_row - 1)
973 440561609 : swap_rows(tab, row, tab->n_row - 1);
974 1000843887 : isl_tab_var_from_row(tab, tab->n_row - 1)->index = -1;
975 1000843887 : tab->n_row--;
976 1000843887 : 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 31909 : int isl_tab_mark_rational(struct isl_tab *tab)
986 : {
987 31909 : if (!tab)
988 0 : return -1;
989 31909 : if (!tab->rational && tab->need_undo)
990 31880 : if (isl_tab_push(tab, isl_tab_undo_rational) < 0)
991 0 : return -1;
992 31909 : tab->rational = 1;
993 31909 : return 0;
994 : }
995 :
996 3276905924 : isl_stat isl_tab_mark_empty(struct isl_tab *tab)
997 : {
998 3276905924 : if (!tab)
999 0 : return isl_stat_error;
1000 3276905924 : if (!tab->empty && tab->need_undo)
1001 2716683228 : if (isl_tab_push(tab, isl_tab_undo_empty) < 0)
1002 0 : return isl_stat_error;
1003 3276905924 : tab->empty = 1;
1004 3276905924 : return isl_stat_ok;
1005 : }
1006 :
1007 38602522827 : int isl_tab_freeze_constraint(struct isl_tab *tab, int con)
1008 : {
1009 : struct isl_tab_var *var;
1010 :
1011 38602522827 : if (!tab)
1012 0 : return -1;
1013 :
1014 38602522827 : var = &tab->con[con];
1015 38602522827 : if (var->frozen)
1016 0 : return 0;
1017 38602522827 : if (var->index < 0)
1018 1606245165 : return 0;
1019 36996277662 : var->frozen = 1;
1020 :
1021 36996277662 : if (tab->need_undo)
1022 36996277662 : 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 25467432504 : static void update_row_sign(struct isl_tab *tab, int row, int col, int row_sgn)
1045 : {
1046 : int i;
1047 25467432504 : struct isl_mat *mat = tab->mat;
1048 25467432504 : unsigned off = 2 + tab->M;
1049 :
1050 25467432504 : if (!tab->row_sign)
1051 25430779528 : return;
1052 :
1053 36652976 : if (tab->row_sign[row] == 0)
1054 5914381 : return;
1055 30738595 : isl_assert(mat->ctx, row_sgn > 0, return);
1056 30738595 : isl_assert(mat->ctx, tab->row_sign[row] == isl_tab_row_neg, return);
1057 30738595 : tab->row_sign[row] = isl_tab_row_pos;
1058 344970473 : for (i = 0; i < tab->n_row; ++i) {
1059 : int s;
1060 314231878 : if (i == row)
1061 30738595 : continue;
1062 283493283 : s = isl_int_sgn(mat->row[i][off + col]);
1063 283493283 : if (!s)
1064 126558891 : continue;
1065 156934392 : if (!tab->row_sign[i])
1066 20228217 : continue;
1067 136706175 : if (s < 0 && tab->row_sign[i] == isl_tab_row_neg)
1068 0 : continue;
1069 136706175 : if (s > 0 && tab->row_sign[i] == isl_tab_row_pos)
1070 73128129 : continue;
1071 63578046 : 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 25467432504 : 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 25467432504 : struct isl_mat *mat = tab->mat;
1132 : struct isl_tab_var *var;
1133 25467432504 : unsigned off = 2 + tab->M;
1134 :
1135 25467432504 : ctx = isl_tab_get_ctx(tab);
1136 25467432504 : if (isl_ctx_next_operation(ctx) < 0)
1137 0 : return -1;
1138 :
1139 25467432504 : isl_int_swap(mat->row[row][0], mat->row[row][off + col]);
1140 25467432504 : sgn = isl_int_sgn(mat->row[row][0]);
1141 25467432504 : if (sgn < 0) {
1142 13533367746 : isl_int_neg(mat->row[row][0], mat->row[row][0]);
1143 13533367746 : isl_int_neg(mat->row[row][off + col], mat->row[row][off + col]);
1144 : } else
1145 94033631982 : for (j = 0; j < off - 1 + tab->n_col; ++j) {
1146 82099567224 : if (j == off - 1 + col)
1147 11934064758 : continue;
1148 70165502466 : isl_int_neg(mat->row[row][1 + j], mat->row[row][1 + j]);
1149 : }
1150 25467432504 : if (!isl_int_is_one(mat->row[row][0]))
1151 11738609453 : isl_seq_normalize(mat->ctx, mat->row[row], off + tab->n_col);
1152 >34032*10^7 : for (i = 0; i < tab->n_row; ++i) {
1153 >31485*10^7 : if (i == row)
1154 25467432504 : continue;
1155 >28938*10^7 : if (isl_int_is_zero(mat->row[i][off + col]))
1156 >17914*10^7 : continue;
1157 >11023*10^7 : isl_int_mul(mat->row[i][0], mat->row[i][0], mat->row[row][0]);
1158 >11063*10^8 : for (j = 0; j < off - 1 + tab->n_col; ++j) {
1159 >99611*10^7 : if (j == off - 1 + col)
1160 >11023*10^7 : continue;
1161 >88587*10^7 : isl_int_mul(mat->row[i][1 + j],
1162 : mat->row[i][1 + j], mat->row[row][0]);
1163 >88587*10^7 : isl_int_addmul(mat->row[i][1 + j],
1164 : mat->row[i][off + col], mat->row[row][1 + j]);
1165 : }
1166 >11023*10^7 : isl_int_mul(mat->row[i][off + col],
1167 : mat->row[i][off + col], mat->row[row][off + col]);
1168 >11023*10^7 : if (!isl_int_is_one(mat->row[i][0]))
1169 >10126*10^7 : isl_seq_normalize(mat->ctx, mat->row[i], off + tab->n_col);
1170 : }
1171 25467432504 : t = tab->row_var[row];
1172 25467432504 : tab->row_var[row] = tab->col_var[col];
1173 25467432504 : tab->col_var[col] = t;
1174 25467432504 : var = isl_tab_var_from_row(tab, row);
1175 25467432504 : var->is_row = 1;
1176 25467432504 : var->index = row;
1177 25467432504 : var = var_from_col(tab, col);
1178 25467432504 : var->is_row = 0;
1179 25467432504 : var->index = col;
1180 25467432504 : update_row_sign(tab, row, col, sgn);
1181 25467432504 : if (tab->in_undo)
1182 3240021376 : return 0;
1183 >23506*10^7 : for (i = tab->n_redundant; i < tab->n_row; ++i) {
1184 >21283*10^7 : if (isl_int_is_zero(mat->row[i][off + col]))
1185 >11323*10^7 : continue;
1186 >18682*10^7 : if (!isl_tab_var_from_row(tab, i)->frozen &&
1187 87222121821 : isl_tab_row_is_redundant(tab, i)) {
1188 6445472906 : int redo = isl_tab_mark_redundant(tab, i);
1189 6445472906 : if (redo < 0)
1190 0 : return -1;
1191 6445472906 : if (redo)
1192 191995097 : --i;
1193 : }
1194 : }
1195 22227411128 : 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 4447611042 : static int to_row(struct isl_tab *tab, struct isl_tab_var *var, int sign)
1206 : {
1207 : int r;
1208 4447611042 : unsigned off = 2 + tab->M;
1209 :
1210 4447611042 : if (var->is_row)
1211 917256817 : return 0;
1212 :
1213 3530354225 : if (sign == 0) {
1214 8178057132 : for (r = tab->n_redundant; r < tab->n_row; ++r)
1215 8178057132 : if (!isl_int_is_zero(tab->mat->row[r][off+var->index]))
1216 : break;
1217 1091446187 : isl_assert(tab->mat->ctx, r < tab->n_row, return -1);
1218 : } else {
1219 2438908038 : r = pivot_row(tab, NULL, sign, var->index);
1220 2438908038 : isl_assert(tab->mat->ctx, r >= 0, return -1);
1221 : }
1222 :
1223 3530354225 : 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 1054932283 : static int sign_of_max(struct isl_tab *tab, struct isl_tab_var *var)
1265 : {
1266 : int row, col;
1267 :
1268 1054932283 : if (max_is_manifestly_unbounded(tab, var))
1269 136693321 : return 1;
1270 918238962 : if (to_row(tab, var, 1) < 0)
1271 0 : return -2;
1272 2380846554 : while (!isl_int_is_pos(tab->mat->row[var->index][1])) {
1273 1134592144 : find_pivot(tab, var, var, 1, &row, &col);
1274 1134592144 : if (row == -1)
1275 405915703 : return isl_int_sgn(tab->mat->row[var->index][1]);
1276 728676441 : if (isl_tab_pivot(tab, row, col) < 0)
1277 0 : return -2;
1278 728676441 : if (!var->is_row) /* manifestly unbounded */
1279 184307811 : return 1;
1280 : }
1281 328015448 : 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 30292538876 : static int row_is_neg(struct isl_tab *tab, int row)
1299 : {
1300 30292538876 : if (!tab->M)
1301 30292538876 : 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 24746236061 : static int row_sgn(struct isl_tab *tab, int row)
1310 : {
1311 24746236061 : if (!tab->M)
1312 24746236061 : 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 26342020442 : static int restore_row(struct isl_tab *tab, struct isl_tab_var *var)
1325 : {
1326 : int row, col;
1327 :
1328 56634559318 : while (row_is_neg(tab, var->index)) {
1329 8792622548 : find_pivot(tab, var, var, 1, &row, &col);
1330 8792622548 : if (row == -1)
1331 3246319733 : break;
1332 5546302815 : if (isl_tab_pivot(tab, row, col) < 0)
1333 0 : return -2;
1334 5546302815 : if (!var->is_row) /* manifestly unbounded */
1335 1595784381 : return 1;
1336 : }
1337 24746236061 : 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 3403914917 : static int at_least_zero(struct isl_tab *tab, struct isl_tab_var *var)
1346 : {
1347 : int row, col;
1348 :
1349 6818744981 : while (isl_int_is_neg(tab->mat->row[var->index][1])) {
1350 3409302333 : find_pivot(tab, var, var, 1, &row, &col);
1351 3409302333 : if (row == -1)
1352 3376365398 : break;
1353 32936935 : if (row == var->index) /* manifestly unbounded */
1354 22021788 : return 1;
1355 10915147 : if (isl_tab_pivot(tab, row, col) < 0)
1356 0 : return -1;
1357 : }
1358 3381893129 : 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 5272141592 : static int sign_of_min(struct isl_tab *tab, struct isl_tab_var *var)
1379 : {
1380 : int row, col;
1381 5272141592 : struct isl_tab_var *pivot_var = NULL;
1382 :
1383 5272141592 : if (min_is_manifestly_unbounded(tab, var))
1384 0 : return -1;
1385 5272141592 : if (!var->is_row) {
1386 424409 : col = var->index;
1387 424409 : row = pivot_row(tab, NULL, -1, col);
1388 424409 : pivot_var = var_from_col(tab, col);
1389 424409 : if (isl_tab_pivot(tab, row, col) < 0)
1390 0 : return -2;
1391 424409 : if (var->is_redundant)
1392 73823 : return 0;
1393 350586 : if (isl_int_is_neg(tab->mat->row[var->index][1])) {
1394 204700 : if (var->is_nonneg) {
1395 409400 : if (!pivot_var->is_redundant &&
1396 204700 : pivot_var->index == row) {
1397 383520 : if (isl_tab_pivot(tab, row, col) < 0)
1398 0 : return -2;
1399 : } else
1400 12940 : if (restore_row(tab, var) < -1)
1401 0 : return -2;
1402 : }
1403 204700 : return -1;
1404 : }
1405 : }
1406 5271863069 : if (var->is_redundant)
1407 0 : return 0;
1408 10554893787 : while (!isl_int_is_neg(tab->mat->row[var->index][1])) {
1409 5276727583 : find_pivot(tab, var, var, -1, &row, &col);
1410 5276727583 : if (row == var->index)
1411 11899331 : return -1;
1412 5264828252 : if (row == -1)
1413 5253432456 : return isl_int_sgn(tab->mat->row[var->index][1]);
1414 11395796 : pivot_var = var_from_col(tab, col);
1415 11395796 : if (isl_tab_pivot(tab, row, col) < 0)
1416 0 : return -2;
1417 11395796 : if (var->is_redundant)
1418 228147 : return 0;
1419 : }
1420 6303135 : if (pivot_var && var->is_nonneg) {
1421 : /* pivot back to non-negative value */
1422 236429 : if (!pivot_var->is_redundant && pivot_var->index == row) {
1423 450152 : if (isl_tab_pivot(tab, row, col) < 0)
1424 0 : return -2;
1425 : } else
1426 11353 : if (restore_row(tab, var) < -1)
1427 0 : return -2;
1428 : }
1429 6303135 : return -1;
1430 : }
1431 :
1432 986131726 : static int row_at_most_neg_one(struct isl_tab *tab, int row)
1433 : {
1434 986131726 : 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 1560348265 : return isl_int_is_neg(tab->mat->row[row][1]) &&
1441 574216539 : 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 1496149807 : 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 1496149807 : if (min_is_manifestly_unbounded(tab, var))
1460 7992499 : return 1;
1461 1488157308 : if (!var->is_row) {
1462 277303849 : col = var->index;
1463 277303849 : row = pivot_row(tab, NULL, -1, col);
1464 277303849 : pivot_var = var_from_col(tab, col);
1465 277303849 : if (isl_tab_pivot(tab, row, col) < 0)
1466 0 : return -1;
1467 277303849 : if (var->is_redundant)
1468 48603990 : return 0;
1469 228699859 : if (row_at_most_neg_one(tab, var->index)) {
1470 189783226 : if (var->is_nonneg) {
1471 379566452 : if (!pivot_var->is_redundant &&
1472 189783226 : pivot_var->index == row) {
1473 378037506 : if (isl_tab_pivot(tab, row, col) < 0)
1474 0 : return -1;
1475 : } else
1476 764473 : if (restore_row(tab, var) < -1)
1477 0 : return -1;
1478 : }
1479 189783226 : return 1;
1480 : }
1481 : }
1482 1249770092 : if (var->is_redundant)
1483 16434195 : return 0;
1484 : do {
1485 1663372101 : find_pivot(tab, var, var, -1, &row, &col);
1486 1663372101 : if (row == var->index) {
1487 363115104 : if (var->is_nonneg && restore_row(tab, var) < -1)
1488 0 : return -1;
1489 363115104 : return 1;
1490 : }
1491 1300256997 : if (row == -1)
1492 185421440 : return 0;
1493 1114835557 : pivot_var = var_from_col(tab, col);
1494 1114835557 : if (isl_tab_pivot(tab, row, col) < 0)
1495 0 : return -1;
1496 1114835557 : if (var->is_redundant)
1497 357403690 : return 0;
1498 757431867 : } while (!row_at_most_neg_one(tab, var->index));
1499 327395663 : if (var->is_nonneg) {
1500 : /* pivot back to non-negative value */
1501 321963780 : if (!pivot_var->is_redundant && pivot_var->index == row)
1502 320433973 : if (isl_tab_pivot(tab, row, col) < 0)
1503 0 : return -1;
1504 321963780 : if (restore_row(tab, var) < -1)
1505 0 : return -1;
1506 : }
1507 327395663 : return 1;
1508 : }
1509 :
1510 : /* Return 1 if "var" can attain values >= 1.
1511 : * Return 0 otherwise.
1512 : */
1513 626245598 : 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 626245598 : if (max_is_manifestly_unbounded(tab, var))
1519 320736053 : return 1;
1520 305509545 : if (to_row(tab, var, 1) < 0)
1521 0 : return -1;
1522 305509545 : r = tab->mat->row[var->index];
1523 631348405 : while (isl_int_lt(r[1], r[0])) {
1524 32300056 : find_pivot(tab, var, var, 1, &row, &col);
1525 32300056 : if (row == -1)
1526 9376890 : return isl_int_ge(r[1], r[0]);
1527 22923166 : if (row == var->index) /* manifestly unbounded */
1528 2593851 : return 1;
1529 20329315 : if (isl_tab_pivot(tab, row, col) < 0)
1530 0 : return -1;
1531 : }
1532 293538804 : return 1;
1533 : }
1534 :
1535 5214788084 : static void swap_cols(struct isl_tab *tab, int col1, int col2)
1536 : {
1537 : int t;
1538 5214788084 : unsigned off = 2 + tab->M;
1539 5214788084 : t = tab->col_var[col1];
1540 5214788084 : tab->col_var[col1] = tab->col_var[col2];
1541 5214788084 : tab->col_var[col2] = t;
1542 5214788084 : var_from_col(tab, col1)->index = col1;
1543 5214788084 : var_from_col(tab, col2)->index = col2;
1544 5214788084 : tab->mat = isl_mat_swap_cols(tab->mat, off + col1, off + col2);
1545 5214788084 : }
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 10268831957 : int isl_tab_kill_col(struct isl_tab *tab, int col)
1560 : {
1561 10268831957 : var_from_col(tab, col)->is_zero = 1;
1562 10268831957 : if (tab->need_undo) {
1563 737945553 : if (isl_tab_push_var(tab, isl_tab_undo_zero,
1564 : var_from_col(tab, col)) < 0)
1565 0 : return -1;
1566 737945553 : if (col != tab->n_dead)
1567 436766699 : swap_cols(tab, col, tab->n_dead);
1568 737945553 : tab->n_dead++;
1569 737945553 : return 0;
1570 : } else {
1571 9530886404 : if (col != tab->n_col - 1)
1572 4601245666 : swap_cols(tab, col, tab->n_col - 1);
1573 9530886404 : var_from_col(tab, tab->n_col - 1)->index = -1;
1574 9530886404 : tab->n_col--;
1575 9530886404 : return 1;
1576 : }
1577 : }
1578 :
1579 1647523367 : static int row_is_manifestly_non_integral(struct isl_tab *tab, int row)
1580 : {
1581 1647523367 : unsigned off = 2 + tab->M;
1582 :
1583 1647523367 : if (tab->M && !isl_int_eq(tab->mat->row[row][2],
1584 : tab->mat->row[row][0]))
1585 0 : return 0;
1586 1647523367 : if (isl_seq_first_non_zero(tab->mat->row[row] + off + tab->n_dead,
1587 1647523367 : tab->n_col - tab->n_dead) != -1)
1588 206912585 : return 0;
1589 :
1590 1440610782 : 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 404827872 : static int tab_is_manifestly_empty(struct isl_tab *tab)
1598 : {
1599 : int i;
1600 :
1601 404827872 : if (tab->empty)
1602 184474 : return 1;
1603 404643398 : if (tab->rational)
1604 574130 : return 0;
1605 :
1606 2243073029 : for (i = 0; i < tab->n_var; ++i) {
1607 1839543564 : if (!tab->var[i].is_row)
1608 192020197 : continue;
1609 1647523367 : if (row_is_manifestly_non_integral(tab, tab->var[i].index))
1610 539803 : return 1;
1611 : }
1612 :
1613 403529465 : 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 404827872 : static isl_stat close_row(struct isl_tab *tab, struct isl_tab_var *var,
1636 : int temp_var)
1637 : {
1638 : int j;
1639 404827872 : struct isl_mat *mat = tab->mat;
1640 404827872 : unsigned off = 2 + tab->M;
1641 :
1642 404827872 : 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 404827872 : var->is_zero = 1;
1647 404827872 : if (!temp_var && tab->need_undo)
1648 5347714 : if (isl_tab_push_var(tab, isl_tab_undo_zero, var) < 0)
1649 0 : return isl_stat_error;
1650 1089124601 : for (j = tab->n_dead; j < tab->n_col; ++j) {
1651 : int recheck;
1652 684296729 : if (isl_int_is_zero(mat->row[var->index][off + j]))
1653 314915985 : continue;
1654 369380744 : 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 369380744 : recheck = isl_tab_kill_col(tab, j);
1659 369380744 : if (recheck < 0)
1660 0 : return isl_stat_error;
1661 369380744 : if (recheck)
1662 361713036 : --j;
1663 : }
1664 404827872 : if (!temp_var && isl_tab_mark_redundant(tab, var->index) < 0)
1665 0 : return isl_stat_error;
1666 404827872 : if (tab_is_manifestly_empty(tab) && isl_tab_mark_empty(tab) < 0)
1667 0 : return isl_stat_error;
1668 404827872 : 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 42249676737 : int isl_tab_allocate_con(struct isl_tab *tab)
1678 : {
1679 : int r;
1680 :
1681 42249676737 : isl_assert(tab->mat->ctx, tab->n_row < tab->mat->n_row, return -1);
1682 42249676737 : isl_assert(tab->mat->ctx, tab->n_con < tab->max_con, return -1);
1683 :
1684 42249676737 : r = tab->n_con;
1685 42249676737 : tab->con[r].index = tab->n_row;
1686 42249676737 : tab->con[r].is_row = 1;
1687 42249676737 : tab->con[r].is_nonneg = 0;
1688 42249676737 : tab->con[r].is_zero = 0;
1689 42249676737 : tab->con[r].is_redundant = 0;
1690 42249676737 : tab->con[r].frozen = 0;
1691 42249676737 : tab->con[r].negated = 0;
1692 42249676737 : tab->row_var[tab->n_row] = ~r;
1693 :
1694 42249676737 : tab->n_row++;
1695 42249676737 : tab->n_con++;
1696 42249676737 : if (isl_tab_push_var(tab, isl_tab_undo_allocate, &tab->con[r]) < 0)
1697 0 : return -1;
1698 :
1699 42249676737 : 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 904087627 : static int var_insert_entry(struct isl_tab *tab, int first)
1708 : {
1709 : int i;
1710 :
1711 904087627 : 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 904087627 : 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 904177258 : for (i = tab->n_var - 1; i >= first; --i) {
1719 89631 : tab->var[i + 1] = tab->var[i];
1720 89631 : if (tab->var[i + 1].is_row)
1721 46632 : tab->row_var[tab->var[i + 1].index]++;
1722 : else
1723 42999 : tab->col_var[tab->var[i + 1].index]++;
1724 : }
1725 :
1726 904087627 : tab->n_var++;
1727 :
1728 904087627 : 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 883582966 : static int var_drop_entry(struct isl_tab *tab, int first)
1737 : {
1738 : int i;
1739 :
1740 883582966 : 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 883582966 : tab->n_var--;
1745 :
1746 883671065 : for (i = first; i < tab->n_var; ++i) {
1747 88099 : tab->var[i] = tab->var[i + 1];
1748 88099 : if (tab->var[i + 1].is_row)
1749 46249 : tab->row_var[tab->var[i].index]--;
1750 : else
1751 41850 : tab->col_var[tab->var[i].index]--;
1752 : }
1753 :
1754 883582966 : 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 904087627 : int isl_tab_insert_var(struct isl_tab *tab, int r)
1762 : {
1763 : int i;
1764 904087627 : unsigned off = 2 + tab->M;
1765 :
1766 904087627 : isl_assert(tab->mat->ctx, tab->n_col < tab->mat->n_col, return -1);
1767 :
1768 904087627 : if (var_insert_entry(tab, r) < 0)
1769 0 : return -1;
1770 :
1771 904087627 : tab->var[r].index = tab->n_col;
1772 904087627 : tab->var[r].is_row = 0;
1773 904087627 : tab->var[r].is_nonneg = 0;
1774 904087627 : tab->var[r].is_zero = 0;
1775 904087627 : tab->var[r].is_redundant = 0;
1776 904087627 : tab->var[r].frozen = 0;
1777 904087627 : tab->var[r].negated = 0;
1778 904087627 : tab->col_var[tab->n_col] = r;
1779 :
1780 14815113546 : for (i = 0; i < tab->n_row; ++i)
1781 13911025919 : isl_int_set_si(tab->mat->row[i][off + tab->n_col], 0);
1782 :
1783 904087627 : tab->n_col++;
1784 904087627 : if (isl_tab_push_var(tab, isl_tab_undo_allocate, &tab->var[r]) < 0)
1785 0 : return -1;
1786 :
1787 904087627 : 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 42226730212 : 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 42226730212 : unsigned off = 2 + tab->M;
1830 :
1831 42226730212 : r = isl_tab_allocate_con(tab);
1832 42226730212 : if (r < 0)
1833 0 : return -1;
1834 :
1835 42226730212 : isl_int_init(a);
1836 42226730212 : isl_int_init(b);
1837 42226730212 : row = tab->mat->row[tab->con[r].index];
1838 42226730212 : isl_int_set_si(row[0], 1);
1839 42226730212 : isl_int_set(row[1], line[0]);
1840 42226730212 : isl_seq_clr(row + 2, tab->M + tab->n_col);
1841 >34849*10^7 : for (i = 0; i < tab->n_var; ++i) {
1842 >30626*10^7 : if (tab->var[i].is_zero)
1843 0 : continue;
1844 >30626*10^7 : if (tab->var[i].is_row) {
1845 >20952*10^7 : isl_int_lcm(a,
1846 : row[0], tab->mat->row[tab->var[i].index][0]);
1847 >20952*10^7 : isl_int_swap(a, row[0]);
1848 >20952*10^7 : isl_int_divexact(a, row[0], a);
1849 >20952*10^7 : isl_int_divexact(b,
1850 : row[0], tab->mat->row[tab->var[i].index][0]);
1851 >20952*10^7 : isl_int_mul(b, b, line[1 + i]);
1852 >41904*10^7 : isl_seq_combine(row + 1, a, row + 1,
1853 >20952*10^7 : b, tab->mat->row[tab->var[i].index] + 1,
1854 >20952*10^7 : 1 + tab->M + tab->n_col);
1855 : } else
1856 96743460504 : isl_int_addmul(row[off + tab->var[i].index],
1857 : line[1 + i], row[0]);
1858 >30626*10^7 : if (tab->M && i >= tab->n_param && i < tab->n_var - tab->n_div)
1859 60281892 : isl_int_submul(row[2], line[1 + i], row[0]);
1860 : }
1861 42226730212 : isl_seq_normalize(tab->mat->ctx, row, off + tab->n_col);
1862 42226730212 : isl_int_clear(a);
1863 42226730212 : isl_int_clear(b);
1864 :
1865 42226730212 : if (tab->row_sign)
1866 42003299 : tab->row_sign[tab->con[r].index] = isl_tab_row_unknown;
1867 :
1868 42226730212 : return r;
1869 : }
1870 :
1871 27917403684 : static isl_stat drop_row(struct isl_tab *tab, int row)
1872 : {
1873 27917403684 : isl_assert(tab->mat->ctx, ~tab->row_var[row] == tab->n_con - 1,
1874 : return isl_stat_error);
1875 27917403684 : if (row != tab->n_row - 1)
1876 6878473967 : swap_rows(tab, row, tab->n_row - 1);
1877 27917403684 : tab->n_row--;
1878 27917403684 : tab->n_con--;
1879 27917403684 : 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 883582966 : static isl_stat drop_col(struct isl_tab *tab, int col)
1889 : {
1890 : int var;
1891 :
1892 883582966 : var = tab->col_var[col];
1893 883582966 : if (col != tab->n_col - 1)
1894 176775719 : swap_cols(tab, col, tab->n_col - 1);
1895 883582966 : tab->n_col--;
1896 883582966 : if (var_drop_entry(tab, var) < 0)
1897 0 : return isl_stat_error;
1898 883582966 : 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 20597166143 : 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 20597166143 : if (!tab)
1914 0 : return isl_stat_error;
1915 20597166143 : if (tab->bmap) {
1916 15852449897 : struct isl_basic_map *bmap = tab->bmap;
1917 :
1918 15852449897 : isl_assert(tab->mat->ctx, tab->n_eq == bmap->n_eq,
1919 : return isl_stat_error);
1920 15852449897 : isl_assert(tab->mat->ctx,
1921 : tab->n_con == bmap->n_eq + bmap->n_ineq,
1922 : return isl_stat_error);
1923 15852449897 : tab->bmap = isl_basic_map_add_ineq(tab->bmap, ineq);
1924 15852449897 : if (isl_tab_push(tab, isl_tab_undo_bmap_ineq) < 0)
1925 0 : return isl_stat_error;
1926 15852449897 : if (!tab->bmap)
1927 0 : return isl_stat_error;
1928 : }
1929 20597166143 : if (tab->cone) {
1930 14072843 : isl_int_init(cst);
1931 14072843 : isl_int_set_si(cst, 0);
1932 14072843 : isl_int_swap(ineq[0], cst);
1933 : }
1934 20597166143 : r = isl_tab_add_row(tab, ineq);
1935 20597166143 : if (tab->cone) {
1936 14072843 : isl_int_swap(ineq[0], cst);
1937 14072843 : isl_int_clear(cst);
1938 : }
1939 20597166143 : if (r < 0)
1940 0 : return isl_stat_error;
1941 20597166143 : tab->con[r].is_nonneg = 1;
1942 20597166143 : if (isl_tab_push_var(tab, isl_tab_undo_nonneg, &tab->con[r]) < 0)
1943 0 : return isl_stat_error;
1944 20597166143 : if (isl_tab_row_is_redundant(tab, tab->con[r].index)) {
1945 5498924488 : if (isl_tab_mark_redundant(tab, tab->con[r].index) < 0)
1946 0 : return isl_stat_error;
1947 5498924488 : return isl_stat_ok;
1948 : }
1949 :
1950 15098241655 : sgn = restore_row(tab, &tab->con[r]);
1951 15098241655 : if (sgn < -1)
1952 0 : return isl_stat_error;
1953 15098241655 : if (sgn < 0)
1954 3240582310 : return isl_tab_mark_empty(tab);
1955 11857659345 : 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 11857659345 : 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 724706714 : static int to_col(struct isl_tab *tab, struct isl_tab_var *var)
1966 : {
1967 : int i;
1968 : int row, col;
1969 724706714 : unsigned off = 2 + tab->M;
1970 :
1971 724706714 : if (!var->is_row)
1972 265079 : return 0;
1973 :
1974 1782348273 : while (isl_int_is_pos(tab->mat->row[var->index][1])) {
1975 976452514 : find_pivot(tab, var, NULL, -1, &row, &col);
1976 976452514 : isl_assert(tab->mat->ctx, row != -1, return -1);
1977 976452514 : if (isl_tab_pivot(tab, row, col) < 0)
1978 0 : return -1;
1979 976452514 : if (!var->is_row)
1980 642987511 : return 0;
1981 : }
1982 :
1983 144856975 : for (i = tab->n_dead; i < tab->n_col; ++i)
1984 144856975 : if (!isl_int_is_zero(tab->mat->row[var->index][off + i]))
1985 : break;
1986 :
1987 81454124 : isl_assert(tab->mat->ctx, i < tab->n_col, return -1);
1988 81454124 : if (isl_tab_pivot(tab, var->index, i) < 0)
1989 0 : return -1;
1990 :
1991 81454124 : 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 9168830118 : static struct isl_tab *add_eq(struct isl_tab *tab, isl_int *eq)
2003 : {
2004 : int i;
2005 : int r;
2006 :
2007 9168830118 : if (!tab)
2008 0 : return NULL;
2009 9168830118 : r = isl_tab_add_row(tab, eq);
2010 9168830118 : if (r < 0)
2011 0 : goto error;
2012 :
2013 9168830118 : r = tab->con[r].index;
2014 9168830118 : i = isl_seq_first_non_zero(tab->mat->row[r] + 2 + tab->M + tab->n_dead,
2015 9168830118 : tab->n_col - tab->n_dead);
2016 9168830118 : isl_assert(tab->mat->ctx, i >= 0, goto error);
2017 9168830118 : i += tab->n_dead;
2018 9168830118 : if (isl_tab_pivot(tab, r, i) < 0)
2019 0 : goto error;
2020 9168830118 : if (isl_tab_kill_col(tab, i) < 0)
2021 0 : goto error;
2022 9168830118 : tab->n_eq++;
2023 :
2024 9168830118 : 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 483196050 : static int row_is_big(struct isl_tab *tab, int row)
2034 : {
2035 483196050 : return tab->M && !isl_int_is_zero(tab->mat->row[row][2]);
2036 : }
2037 :
2038 744338189 : static int row_is_manifestly_zero(struct isl_tab *tab, int row)
2039 : {
2040 744338189 : unsigned off = 2 + tab->M;
2041 :
2042 744338189 : if (!isl_int_is_zero(tab->mat->row[row][1]))
2043 644193690 : return 0;
2044 100144499 : if (row_is_big(tab, row))
2045 0 : return 0;
2046 200288998 : return isl_seq_first_non_zero(tab->mat->row[row] + off + tab->n_dead,
2047 200288998 : 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 729999256 : 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 729999256 : if (!tab)
2061 0 : return -1;
2062 729999256 : r = isl_tab_add_row(tab, eq);
2063 729999256 : if (r < 0)
2064 0 : return -1;
2065 :
2066 729999256 : var = &tab->con[r];
2067 729999256 : r = var->index;
2068 729999256 : if (row_is_manifestly_zero(tab, r)) {
2069 14119176 : var->is_zero = 1;
2070 14119176 : if (isl_tab_mark_redundant(tab, r) < 0)
2071 0 : return -1;
2072 14119176 : return 0;
2073 : }
2074 :
2075 715880080 : if (isl_int_is_neg(tab->mat->row[r][1])) {
2076 269502176 : isl_seq_neg(tab->mat->row[r] + 1, tab->mat->row[r] + 1,
2077 269502176 : 1 + tab->n_col);
2078 269502176 : var->negated = 1;
2079 : }
2080 715880080 : var->is_nonneg = 1;
2081 715880080 : if (to_col(tab, var) < 0)
2082 0 : return -1;
2083 715880080 : var->is_nonneg = 0;
2084 715880080 : if (isl_tab_kill_col(tab, var->index) < 0)
2085 0 : return -1;
2086 :
2087 715880080 : 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 8310967 : static int add_zero_row(struct isl_tab *tab)
2097 : {
2098 : int r;
2099 : isl_int *row;
2100 :
2101 8310967 : r = isl_tab_allocate_con(tab);
2102 8310967 : if (r < 0)
2103 0 : return -1;
2104 :
2105 8310967 : row = tab->mat->row[tab->con[r].index];
2106 8310967 : isl_seq_clr(row + 1, 1 + tab->M + tab->n_col);
2107 8310967 : isl_int_set_si(row[0], 1);
2108 :
2109 8310967 : 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 14338933 : int isl_tab_add_eq(struct isl_tab *tab, isl_int *eq)
2120 : {
2121 14338933 : 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 14338933 : if (!tab)
2129 0 : return -1;
2130 14338933 : isl_assert(tab->mat->ctx, !tab->M, return -1);
2131 :
2132 14338933 : if (tab->need_undo)
2133 14338933 : snap = isl_tab_snap(tab);
2134 :
2135 14338933 : if (tab->cone) {
2136 36535 : isl_int_init(cst);
2137 36535 : isl_int_set_si(cst, 0);
2138 36535 : isl_int_swap(eq[0], cst);
2139 : }
2140 14338933 : r = isl_tab_add_row(tab, eq);
2141 14338933 : if (tab->cone) {
2142 36535 : isl_int_swap(eq[0], cst);
2143 36535 : isl_int_clear(cst);
2144 : }
2145 14338933 : if (r < 0)
2146 0 : return -1;
2147 :
2148 14338933 : var = &tab->con[r];
2149 14338933 : row = var->index;
2150 14338933 : if (row_is_manifestly_zero(tab, row)) {
2151 5512299 : if (snap)
2152 5512299 : return isl_tab_rollback(tab, snap);
2153 0 : return drop_row(tab, row);
2154 : }
2155 :
2156 8826634 : if (tab->bmap) {
2157 8310967 : tab->bmap = isl_basic_map_add_ineq(tab->bmap, eq);
2158 8310967 : if (isl_tab_push(tab, isl_tab_undo_bmap_ineq) < 0)
2159 0 : return -1;
2160 8310967 : isl_seq_neg(eq, eq, 1 + tab->n_var);
2161 8310967 : tab->bmap = isl_basic_map_add_ineq(tab->bmap, eq);
2162 8310967 : isl_seq_neg(eq, eq, 1 + tab->n_var);
2163 8310967 : if (isl_tab_push(tab, isl_tab_undo_bmap_ineq) < 0)
2164 0 : return -1;
2165 8310967 : if (!tab->bmap)
2166 0 : return -1;
2167 8310967 : if (add_zero_row(tab) < 0)
2168 0 : return -1;
2169 : }
2170 :
2171 8826634 : sgn = isl_int_sgn(tab->mat->row[row][1]);
2172 :
2173 8826634 : if (sgn > 0) {
2174 317570 : isl_seq_neg(tab->mat->row[row] + 1, tab->mat->row[row] + 1,
2175 317570 : 1 + tab->n_col);
2176 317570 : var->negated = 1;
2177 317570 : sgn = -1;
2178 : }
2179 :
2180 8826634 : if (sgn < 0) {
2181 6377629 : sgn = sign_of_max(tab, var);
2182 6377629 : if (sgn < -1)
2183 0 : return -1;
2184 6377629 : if (sgn < 0) {
2185 0 : if (isl_tab_mark_empty(tab) < 0)
2186 0 : return -1;
2187 0 : return 0;
2188 : }
2189 : }
2190 :
2191 8826634 : var->is_nonneg = 1;
2192 8826634 : if (to_col(tab, var) < 0)
2193 0 : return -1;
2194 8826634 : var->is_nonneg = 0;
2195 8826634 : if (isl_tab_kill_col(tab, var->index) < 0)
2196 0 : return -1;
2197 :
2198 8826634 : 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 885050817 : 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 885050817 : if (!bmap)
2218 0 : return NULL;
2219 :
2220 885050817 : total = isl_basic_map_total_dim(bmap);
2221 885050817 : div_pos = 1 + total - bmap->n_div + div;
2222 :
2223 885050817 : ineq = isl_vec_alloc(bmap->ctx, 1 + total);
2224 885050817 : if (!ineq)
2225 0 : return NULL;
2226 :
2227 885050817 : isl_seq_cpy(ineq->el, bmap->div[div] + 1, 1 + total);
2228 885050817 : isl_int_neg(ineq->el[div_pos], bmap->div[div][0]);
2229 885050817 : 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 885050817 : 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 885050817 : total = isl_basic_map_total_dim(tab->bmap);
2255 885050817 : div_pos = 1 + total - tab->bmap->n_div + div;
2256 :
2257 885050817 : ineq = ineq_for_div(tab->bmap, div);
2258 885050817 : if (!ineq)
2259 0 : goto error;
2260 :
2261 885050817 : if (add_ineq) {
2262 13755136 : if (add_ineq(user, ineq->el) < 0)
2263 0 : goto error;
2264 : } else {
2265 871295681 : if (isl_tab_add_ineq(tab, ineq->el) < 0)
2266 0 : goto error;
2267 : }
2268 :
2269 885050817 : isl_seq_neg(ineq->el, tab->bmap->div[div] + 1, 1 + total);
2270 885050817 : isl_int_set(ineq->el[div_pos], tab->bmap->div[div][0]);
2271 885050817 : isl_int_add(ineq->el[0], ineq->el[0], ineq->el[div_pos]);
2272 885050817 : isl_int_sub_ui(ineq->el[0], ineq->el[0], 1);
2273 :
2274 885050817 : if (add_ineq) {
2275 13755136 : if (add_ineq(user, ineq->el) < 0)
2276 0 : goto error;
2277 : } else {
2278 871295681 : if (isl_tab_add_ineq(tab, ineq->el) < 0)
2279 0 : goto error;
2280 : }
2281 :
2282 885050817 : isl_vec_free(ineq);
2283 :
2284 885050817 : 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 885050817 : static int div_is_nonneg(struct isl_tab *tab, __isl_keep isl_vec *div)
2297 : {
2298 : int i;
2299 :
2300 885050817 : if (tab->M)
2301 0 : return 1;
2302 :
2303 885050817 : if (isl_int_is_neg(div->el[1]))
2304 325647541 : return 0;
2305 :
2306 1548123438 : for (i = 0; i < tab->n_var; ++i) {
2307 1491692222 : if (isl_int_is_neg(div->el[2 + i]))
2308 178036966 : return 0;
2309 1313655256 : if (isl_int_is_zero(div->el[2 + i]))
2310 896059033 : continue;
2311 417596223 : if (!tab->var[i].is_nonneg)
2312 324935094 : return 0;
2313 : }
2314 :
2315 56431216 : 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 885050817 : 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 885050817 : if (!tab || !div)
2337 0 : return -1;
2338 :
2339 885050817 : 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 885050817 : isl_assert(tab->mat->ctx, tab->bmap, return -1);
2344 885050817 : n_div = isl_basic_map_dim(tab->bmap, isl_dim_div);
2345 885050817 : o_div = tab->n_var - n_div;
2346 885050817 : 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 885050817 : nonneg = div_is_nonneg(tab, div);
2351 :
2352 885050817 : if (isl_tab_extend_cons(tab, 3) < 0)
2353 0 : return -1;
2354 885050817 : if (isl_tab_extend_vars(tab, 1) < 0)
2355 0 : return -1;
2356 885050817 : r = isl_tab_insert_var(tab, pos);
2357 885050817 : if (r < 0)
2358 0 : return -1;
2359 :
2360 885050817 : if (nonneg)
2361 56431216 : tab->var[r].is_nonneg = 1;
2362 :
2363 885050817 : tab->bmap = isl_basic_map_insert_div(tab->bmap, pos - o_div, div);
2364 885050817 : if (!tab->bmap)
2365 0 : return -1;
2366 885050817 : if (isl_tab_push_var(tab, isl_tab_undo_bmap_div, &tab->var[r]) < 0)
2367 0 : return -1;
2368 :
2369 885050817 : if (add_div_constraints(tab, pos - o_div, add_ineq, user) < 0)
2370 0 : return -1;
2371 :
2372 885050817 : 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 871295681 : int isl_tab_add_div(struct isl_tab *tab, __isl_keep isl_vec *div)
2379 : {
2380 871295681 : if (!tab)
2381 0 : return -1;
2382 871295681 : 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 3011506644 : __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 3011506644 : if (!bmap)
2397 0 : return NULL;
2398 6023013288 : tab = isl_tab_alloc(bmap->ctx,
2399 3011506644 : isl_basic_map_total_dim(bmap) + bmap->n_ineq + 1,
2400 : isl_basic_map_total_dim(bmap), 0);
2401 3011506644 : if (!tab)
2402 0 : return NULL;
2403 3011506644 : tab->preserve = track;
2404 3011506644 : tab->rational = ISL_F_ISSET(bmap, ISL_BASIC_MAP_RATIONAL);
2405 3011506644 : 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 12171162391 : for (i = 0; i < bmap->n_eq; ++i) {
2411 9159655747 : tab = add_eq(tab, bmap->eq[i]);
2412 9159655747 : if (!tab)
2413 0 : return tab;
2414 : }
2415 7704016196 : for (i = 0; i < bmap->n_ineq; ++i) {
2416 4729822087 : if (isl_tab_add_ineq(tab, bmap->ineq[i]) < 0)
2417 0 : goto error;
2418 4729822087 : if (tab->empty)
2419 37312535 : goto done;
2420 : }
2421 : done:
2422 3011506644 : if (track && isl_tab_track_bmap(tab, isl_basic_map_copy(bmap)) < 0)
2423 0 : goto error;
2424 3011506644 : return tab;
2425 : error:
2426 0 : isl_tab_free(tab);
2427 0 : return NULL;
2428 : }
2429 :
2430 196218765 : __isl_give struct isl_tab *isl_tab_from_basic_set(
2431 : __isl_keep isl_basic_set *bset, int track)
2432 : {
2433 196218765 : return isl_tab_from_basic_map(bset, track);
2434 : }
2435 :
2436 : /* Construct a tableau corresponding to the recession cone of "bset".
2437 : */
2438 17412265 : 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 17412265 : unsigned offset = 0;
2445 :
2446 17412265 : if (!bset)
2447 0 : return NULL;
2448 17412265 : if (parametric)
2449 2212793 : offset = isl_basic_set_dim(bset, isl_dim_param);
2450 17412265 : tab = isl_tab_alloc(bset->ctx, bset->n_eq + bset->n_ineq,
2451 17412265 : isl_basic_set_total_dim(bset) - offset, 0);
2452 17412265 : if (!tab)
2453 0 : return NULL;
2454 17412265 : tab->rational = ISL_F_ISSET(bset, ISL_BASIC_SET_RATIONAL);
2455 17412265 : tab->cone = 1;
2456 :
2457 17412265 : isl_int_init(cst);
2458 17412265 : isl_int_set_si(cst, 0);
2459 26586636 : for (i = 0; i < bset->n_eq; ++i) {
2460 9174371 : isl_int_swap(bset->eq[i][offset], cst);
2461 9174371 : if (offset > 0) {
2462 0 : if (isl_tab_add_eq(tab, bset->eq[i] + offset) < 0)
2463 0 : goto error;
2464 : } else
2465 9174371 : tab = add_eq(tab, bset->eq[i]);
2466 9174371 : isl_int_swap(bset->eq[i][offset], cst);
2467 9174371 : if (!tab)
2468 0 : goto done;
2469 : }
2470 137056856 : for (i = 0; i < bset->n_ineq; ++i) {
2471 : int r;
2472 119644591 : isl_int_swap(bset->ineq[i][offset], cst);
2473 119644591 : r = isl_tab_add_row(tab, bset->ineq[i] + offset);
2474 119644591 : isl_int_swap(bset->ineq[i][offset], cst);
2475 119644591 : if (r < 0)
2476 0 : goto error;
2477 119644591 : tab->con[r].is_nonneg = 1;
2478 119644591 : if (isl_tab_push_var(tab, isl_tab_undo_nonneg, &tab->con[r]) < 0)
2479 0 : goto error;
2480 : }
2481 : done:
2482 17412265 : isl_int_clear(cst);
2483 17412265 : 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 2212793 : isl_bool isl_tab_cone_is_bounded(struct isl_tab *tab)
2494 : {
2495 : int i;
2496 :
2497 2212793 : if (!tab)
2498 0 : return isl_bool_error;
2499 2212793 : if (tab->empty)
2500 0 : return isl_bool_true;
2501 2212793 : if (tab->n_dead == tab->n_col)
2502 2146727 : return isl_bool_true;
2503 :
2504 : for (;;) {
2505 482472 : for (i = tab->n_redundant; i < tab->n_row; ++i) {
2506 : struct isl_tab_var *var;
2507 : int sgn;
2508 241236 : var = isl_tab_var_from_row(tab, i);
2509 241236 : if (!var->is_nonneg)
2510 171475 : continue;
2511 69761 : sgn = sign_of_max(tab, var);
2512 69761 : if (sgn < -1)
2513 0 : return isl_bool_error;
2514 69761 : if (sgn != 0)
2515 0 : return isl_bool_false;
2516 69761 : if (close_row(tab, var, 0) < 0)
2517 0 : return isl_bool_error;
2518 69761 : break;
2519 : }
2520 69761 : if (tab->n_dead == tab->n_col)
2521 66066 : return isl_bool_true;
2522 3695 : if (i == tab->n_row)
2523 0 : return isl_bool_false;
2524 3695 : }
2525 : }
2526 :
2527 1021716250 : int isl_tab_sample_is_integer(struct isl_tab *tab)
2528 : {
2529 : int i;
2530 :
2531 1021716250 : if (!tab)
2532 0 : return -1;
2533 :
2534 3199447804 : for (i = 0; i < tab->n_var; ++i) {
2535 : int row;
2536 2967648193 : if (!tab->var[i].is_row)
2537 279998600 : continue;
2538 2687649593 : row = tab->var[i].index;
2539 2687649593 : if (!isl_int_is_divisible_by(tab->mat->row[row][1],
2540 : tab->mat->row[row][0]))
2541 789916639 : return 0;
2542 : }
2543 231799611 : return 1;
2544 : }
2545 :
2546 175239834 : static struct isl_vec *extract_integer_sample(struct isl_tab *tab)
2547 : {
2548 : int i;
2549 : struct isl_vec *vec;
2550 :
2551 175239834 : vec = isl_vec_alloc(tab->mat->ctx, 1 + tab->n_var);
2552 175239834 : if (!vec)
2553 0 : return NULL;
2554 :
2555 175239834 : isl_int_set_si(vec->block.data[0], 1);
2556 1004470746 : for (i = 0; i < tab->n_var; ++i) {
2557 829230912 : if (!tab->var[i].is_row)
2558 204378637 : isl_int_set_si(vec->block.data[1 + i], 0);
2559 : else {
2560 624852275 : int row = tab->var[i].index;
2561 624852275 : isl_int_divexact(vec->block.data[1 + i],
2562 : tab->mat->row[row][1], tab->mat->row[row][0]);
2563 : }
2564 : }
2565 :
2566 175239834 : return vec;
2567 : }
2568 :
2569 340371861 : 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 340371861 : if (!tab)
2576 0 : return NULL;
2577 :
2578 340371861 : vec = isl_vec_alloc(tab->mat->ctx, 1 + tab->n_var);
2579 340371861 : if (!vec)
2580 0 : return NULL;
2581 :
2582 340371861 : isl_int_init(m);
2583 :
2584 340371861 : isl_int_set_si(vec->block.data[0], 1);
2585 2435484845 : for (i = 0; i < tab->n_var; ++i) {
2586 : int row;
2587 2095112984 : if (!tab->var[i].is_row) {
2588 261832384 : isl_int_set_si(vec->block.data[1 + i], 0);
2589 261832384 : continue;
2590 : }
2591 1833280600 : row = tab->var[i].index;
2592 1833280600 : isl_int_gcd(m, vec->block.data[0], tab->mat->row[row][0]);
2593 1833280600 : isl_int_divexact(m, tab->mat->row[row][0], m);
2594 1833280600 : isl_seq_scale(vec->block.data, vec->block.data, m, 1 + i);
2595 1833280600 : isl_int_divexact(m, vec->block.data[0], tab->mat->row[row][0]);
2596 1833280600 : isl_int_mul(vec->block.data[1 + i], m, tab->mat->row[row][1]);
2597 : }
2598 340371861 : vec = isl_vec_normalize(vec);
2599 :
2600 340371861 : isl_int_clear(m);
2601 340371861 : 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 960296052 : static void get_rounded_sample_value(struct isl_tab *tab,
2608 : struct isl_tab_var *var, int sgn, isl_int *v)
2609 : {
2610 960296052 : if (!var->is_row)
2611 1562631 : isl_int_set_si(*v, 0);
2612 958733421 : else if (sgn > 0)
2613 956403940 : isl_int_cdiv_q(*v, tab->mat->row[var->index][1],
2614 : tab->mat->row[var->index][0]);
2615 : else
2616 2329481 : isl_int_fdiv_q(*v, tab->mat->row[var->index][1],
2617 : tab->mat->row[var->index][0]);
2618 960296052 : }
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 1809800654 : 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 1809800654 : if (!bmap)
2635 0 : return NULL;
2636 1809800654 : if (!tab)
2637 0 : return bmap;
2638 :
2639 1809800654 : n_eq = tab->n_eq;
2640 1809800654 : if (tab->empty)
2641 1317917 : bmap = isl_basic_map_set_to_empty(bmap);
2642 : else
2643 5230430358 : for (i = bmap->n_ineq - 1; i >= 0; --i) {
2644 3421947621 : if (isl_tab_is_equality(tab, n_eq + i))
2645 656460188 : isl_basic_map_inequality_to_equality(bmap, i);
2646 2765487433 : else if (isl_tab_is_redundant(tab, n_eq + i))
2647 340347401 : isl_basic_map_drop_inequality(bmap, i);
2648 : }
2649 1809800654 : if (bmap->n_eq != n_eq)
2650 96104534 : bmap = isl_basic_map_gauss(bmap, NULL);
2651 1809800654 : if (!tab->rational &&
2652 419482629 : bmap && !bmap->sample && isl_tab_sample_is_integer(tab))
2653 175239834 : bmap->sample = extract_integer_sample(tab);
2654 1809800654 : return bmap;
2655 : }
2656 :
2657 22440735 : struct isl_basic_set *isl_basic_set_update_from_tab(struct isl_basic_set *bset,
2658 : struct isl_tab *tab)
2659 : {
2660 22440735 : 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 11175886 : static isl_stat drop_last_con_in_row(struct isl_tab *tab, int r)
2668 : {
2669 11175886 : 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 11175886 : 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 11175886 : if (drop_row(tab, tab->con[r].index) < 0)
2677 0 : return isl_stat_error;
2678 :
2679 11175886 : 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 11175886 : 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 11175886 : unsigned off = 2 + tab->M;
2698 :
2699 11175886 : if (var->is_zero)
2700 0 : return isl_stat_ok;
2701 11175886 : 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 11175886 : if (isl_tab_extend_cons(tab, 1) < 0)
2707 0 : return isl_stat_error;
2708 :
2709 11175886 : r = tab->n_con;
2710 11175886 : tab->con[r].index = tab->n_row;
2711 11175886 : tab->con[r].is_row = 1;
2712 11175886 : tab->con[r].is_nonneg = 0;
2713 11175886 : tab->con[r].is_zero = 0;
2714 11175886 : tab->con[r].is_redundant = 0;
2715 11175886 : tab->con[r].frozen = 0;
2716 11175886 : tab->con[r].negated = 0;
2717 11175886 : tab->row_var[tab->n_row] = ~r;
2718 11175886 : row = tab->mat->row[tab->n_row];
2719 :
2720 11175886 : if (var->is_row) {
2721 9866529 : isl_int_set(row[0], tab->mat->row[var->index][0]);
2722 19733058 : isl_seq_neg(row + 1,
2723 19733058 : tab->mat->row[var->index] + 1, 1 + tab->n_col);
2724 : } else {
2725 1309357 : isl_int_set_si(row[0], 1);
2726 1309357 : isl_seq_clr(row + 1, 1 + tab->n_col);
2727 1309357 : isl_int_set_si(row[off + var->index], -1);
2728 : }
2729 :
2730 11175886 : tab->n_row++;
2731 11175886 : tab->n_con++;
2732 :
2733 11175886 : sgn = sign_of_max(tab, &tab->con[r]);
2734 11175886 : if (sgn < -1)
2735 0 : return isl_stat_error;
2736 11175886 : if (sgn < 0) {
2737 146731 : if (drop_last_con_in_row(tab, r) < 0)
2738 0 : return isl_stat_error;
2739 146731 : if (isl_tab_mark_empty(tab) < 0)
2740 0 : return isl_stat_error;
2741 146731 : return isl_stat_ok;
2742 : }
2743 11029155 : tab->con[r].is_nonneg = 1;
2744 : /* sgn == 0 */
2745 11029155 : if (close_row(tab, &tab->con[r], 1) < 0)
2746 0 : return isl_stat_error;
2747 11029155 : if (drop_last_con_in_row(tab, r) < 0)
2748 0 : return isl_stat_error;
2749 :
2750 11029155 : return isl_stat_ok;
2751 : }
2752 :
2753 : /* Check that "con" is a valid constraint position for "tab".
2754 : */
2755 4020607654 : static isl_stat isl_tab_check_con(struct isl_tab *tab, int con)
2756 : {
2757 4020607654 : if (!tab)
2758 0 : return isl_stat_error;
2759 4020607654 : 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 4020607654 : 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 1620703 : int isl_tab_relax(struct isl_tab *tab, int con)
2783 : {
2784 : struct isl_tab_var *var;
2785 :
2786 1620703 : if (!tab)
2787 0 : return -1;
2788 :
2789 1620703 : var = &tab->con[con];
2790 :
2791 1620703 : 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 1620703 : 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 1620703 : if (!var->is_row && !max_is_manifestly_unbounded(tab, var))
2799 231338 : if (to_row(tab, var, 1) < 0)
2800 0 : return -1;
2801 1620703 : if (!var->is_row && !min_is_manifestly_unbounded(tab, var))
2802 22257 : if (to_row(tab, var, -1) < 0)
2803 0 : return -1;
2804 :
2805 1620703 : if (var->is_row) {
2806 296169 : isl_int_add(tab->mat->row[var->index][1],
2807 : tab->mat->row[var->index][1], tab->mat->row[var->index][0]);
2808 296169 : if (restore_row(tab, var) < 0)
2809 0 : return -1;
2810 : } else {
2811 : int i;
2812 1324534 : unsigned off = 2 + tab->M;
2813 :
2814 6255559 : for (i = 0; i < tab->n_row; ++i) {
2815 4931025 : if (isl_int_is_zero(tab->mat->row[i][off + var->index]))
2816 3238131 : continue;
2817 1692894 : 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 1620703 : if (isl_tab_push_var(tab, isl_tab_undo_relax, var) < 0)
2824 0 : return -1;
2825 :
2826 1620703 : 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 46615 : int isl_tab_shift_var(struct isl_tab *tab, int pos, isl_int shift)
2847 : {
2848 : struct isl_tab_var *var;
2849 :
2850 46615 : if (!tab)
2851 0 : return -1;
2852 46615 : if (isl_int_is_zero(shift))
2853 26295 : return 0;
2854 :
2855 20320 : var = &tab->var[pos];
2856 20320 : if (!var->is_row) {
2857 7670 : if (isl_int_is_neg(shift)) {
2858 5840 : if (!max_is_manifestly_unbounded(tab, var))
2859 1023 : if (to_row(tab, var, 1) < 0)
2860 0 : return -1;
2861 : } else {
2862 1830 : if (!min_is_manifestly_unbounded(tab, var))
2863 187 : if (to_row(tab, var, -1) < 0)
2864 0 : return -1;
2865 : }
2866 : }
2867 :
2868 20320 : if (var->is_row) {
2869 13860 : isl_int_addmul(tab->mat->row[var->index][1],
2870 : shift, tab->mat->row[var->index][0]);
2871 : } else {
2872 : int i;
2873 6460 : unsigned off = 2 + tab->M;
2874 :
2875 38765 : for (i = 0; i < tab->n_row; ++i) {
2876 32305 : if (isl_int_is_zero(tab->mat->row[i][off + var->index]))
2877 23518 : continue;
2878 8787 : isl_int_submul(tab->mat->row[i][1],
2879 : shift, tab->mat->row[i][off + var->index]);
2880 : }
2881 :
2882 : }
2883 :
2884 20320 : 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 505574 : int isl_tab_unrestrict(struct isl_tab *tab, int con)
2893 : {
2894 : struct isl_tab_var *var;
2895 :
2896 505574 : if (!tab)
2897 0 : return -1;
2898 :
2899 505574 : var = &tab->con[con];
2900 505574 : if (!var->is_nonneg)
2901 0 : return 0;
2902 :
2903 505574 : var->is_nonneg = 0;
2904 505574 : if (isl_tab_push_var(tab, isl_tab_undo_unrestrict, var) < 0)
2905 0 : return -1;
2906 :
2907 505574 : return 0;
2908 : }
2909 :
2910 1798996 : int isl_tab_select_facet(struct isl_tab *tab, int con)
2911 : {
2912 1798996 : if (!tab)
2913 0 : return -1;
2914 :
2915 1798996 : return cut_to_hyperplane(tab, &tab->con[con]);
2916 : }
2917 :
2918 3192934462 : static int may_be_equality(struct isl_tab *tab, int row)
2919 : {
2920 6412382955 : return tab->rational ? isl_int_is_zero(tab->mat->row[row][1])
2921 3219448493 : : 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 2201740862 : static struct isl_tab_var *select_marked(struct isl_tab *tab)
2937 : {
2938 : int i;
2939 : struct isl_tab_var *var;
2940 :
2941 11126982264 : for (i = tab->n_con - 1; i >= 0; --i) {
2942 10985171246 : var = &tab->con[i];
2943 10985171246 : if (var->index < 0)
2944 3137317851 : continue;
2945 7847853395 : if (var->is_row && var->index < tab->n_redundant)
2946 1808338049 : continue;
2947 6039515346 : if (!var->is_row && var->index < tab->n_dead)
2948 12914808 : continue;
2949 6026600538 : if (var->marked)
2950 2059929844 : return var;
2951 : }
2952 :
2953 141811018 : 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 1737865772 : int isl_tab_detect_implicit_equalities(struct isl_tab *tab)
2974 : {
2975 : int i;
2976 : unsigned n_marked;
2977 :
2978 1737865772 : if (!tab)
2979 0 : return -1;
2980 1737865772 : if (tab->empty)
2981 1559808 : return 0;
2982 1736305964 : if (tab->n_dead == tab->n_col)
2983 1393650194 : return 0;
2984 :
2985 342655770 : n_marked = 0;
2986 2535471149 : for (i = tab->n_redundant; i < tab->n_row; ++i) {
2987 2192815379 : struct isl_tab_var *var = isl_tab_var_from_row(tab, i);
2988 3855809646 : var->marked = !var->frozen && var->is_nonneg &&
2989 1662994267 : may_be_equality(tab, i);
2990 2192815379 : if (var->marked)
2991 1113477557 : n_marked++;
2992 : }
2993 1691445970 : for (i = tab->n_dead; i < tab->n_col; ++i) {
2994 1348790200 : struct isl_tab_var *var = var_from_col(tab, i);
2995 1348790200 : var->marked = !var->frozen && var->is_nonneg;
2996 1348790200 : if (var->marked)
2997 407870454 : n_marked++;
2998 : }
2999 1713243657 : while (n_marked) {
3000 : struct isl_tab_var *var;
3001 : int sgn;
3002 1152492171 : var = select_marked(tab);
3003 1152492171 : if (!var)
3004 115183164 : break;
3005 1037309007 : var->marked = 0;
3006 1037309007 : n_marked--;
3007 1037309007 : sgn = sign_of_max(tab, var);
3008 1037309007 : if (sgn < 0)
3009 0 : return -1;
3010 1037309007 : if (sgn == 0) {
3011 393728956 : if (close_row(tab, var, 0) < 0)
3012 0 : return -1;
3013 643580051 : } else if (!tab->rational && !at_least_one(tab, var)) {
3014 9376890 : if (cut_to_hyperplane(tab, var) < 0)
3015 0 : return -1;
3016 9376890 : return isl_tab_detect_implicit_equalities(tab);
3017 : }
3018 6767239525 : for (i = tab->n_redundant; i < tab->n_row; ++i) {
3019 5739307408 : var = isl_tab_var_from_row(tab, i);
3020 5739307408 : if (!var->marked)
3021 4209367213 : continue;
3022 1529940195 : if (may_be_equality(tab, i))
3023 1509842853 : continue;
3024 20097342 : var->marked = 0;
3025 20097342 : n_marked--;
3026 : }
3027 : }
3028 :
3029 333278880 : 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 595442 : static int update_con_after_move(struct isl_tab *tab, int i, int old)
3036 : {
3037 : int *p;
3038 : int index;
3039 :
3040 595442 : index = tab->con[i].index;
3041 595442 : if (index == -1)
3042 420564 : return 0;
3043 174878 : p = tab->con[i].is_row ? tab->row_var : tab->col_var;
3044 174878 : if (p[index] != ~old)
3045 0 : isl_die(tab->mat->ctx, isl_error_internal,
3046 : "broken internal state", return -1);
3047 174878 : p[index] = ~i;
3048 :
3049 174878 : 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 21766 : isl_stat isl_tab_swap_constraints(struct isl_tab *tab, int con1, int con2)
3058 : {
3059 : struct isl_tab_var var;
3060 :
3061 43532 : if (isl_tab_check_con(tab, con1) < 0 ||
3062 21766 : isl_tab_check_con(tab, con2) < 0)
3063 0 : return isl_stat_error;
3064 :
3065 21766 : var = tab->con[con1];
3066 21766 : tab->con[con1] = tab->con[con2];
3067 21766 : if (update_con_after_move(tab, con1, con2) < 0)
3068 0 : return isl_stat_error;
3069 21766 : tab->con[con2] = var;
3070 21766 : if (update_con_after_move(tab, con2, con1) < 0)
3071 0 : return isl_stat_error;
3072 :
3073 21766 : 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 133424 : 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 133424 : if (n <= 1)
3085 40624 : return 0;
3086 :
3087 92800 : last = first + n - 1;
3088 92800 : var = tab->con[last];
3089 517436 : for (i = last; i > first; --i) {
3090 424636 : tab->con[i] = tab->con[i - 1];
3091 424636 : if (update_con_after_move(tab, i, i - 1) < 0)
3092 0 : return -1;
3093 : }
3094 92800 : tab->con[first] = var;
3095 92800 : if (update_con_after_move(tab, first, last) < 0)
3096 0 : return -1;
3097 :
3098 92800 : 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 13086 : static isl_stat con_drop_entries(struct isl_tab *tab,
3107 : unsigned first, unsigned n)
3108 : {
3109 : int i;
3110 :
3111 13086 : 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 13086 : tab->n_con -= n;
3116 :
3117 47560 : for (i = first; i < tab->n_con; ++i) {
3118 34474 : tab->con[i] = tab->con[i + n];
3119 34474 : if (update_con_after_move(tab, i, i + n) < 0)
3120 0 : return isl_stat_error;
3121 : }
3122 :
3123 13086 : 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 21766 : static isl_stat swap_eq(unsigned a, unsigned b, void *user)
3131 : {
3132 21766 : struct isl_tab *tab = user;
3133 :
3134 21766 : 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 13086 : static isl_stat drop_eq(unsigned n, void *user)
3146 : {
3147 13086 : struct isl_tab *tab = user;
3148 :
3149 13086 : if (tab->n_eq == n)
3150 0 : return isl_tab_mark_empty(tab);
3151 :
3152 13086 : tab->n_eq -= n;
3153 13086 : 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 25902 : 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 25902 : single = isl_basic_map_has_single_reference(bmap);
3165 25902 : if (single < 0)
3166 0 : return isl_basic_map_free(bmap);
3167 25902 : if (single)
3168 6796 : return bmap;
3169 19106 : 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 1378158306 : __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 1378158306 : if (!tab || !bmap)
3212 0 : return isl_basic_map_free(bmap);
3213 1378158306 : if (tab->empty)
3214 2583 : return bmap;
3215 :
3216 1378155723 : n_eq = tab->n_eq;
3217 1396999266 : for (i = bmap->n_ineq - 1; i >= 0; --i) {
3218 18843543 : if (!isl_tab_is_equality(tab, bmap->n_eq + i))
3219 18776831 : continue;
3220 66712 : isl_basic_map_inequality_to_equality(bmap, i);
3221 66712 : if (rotate_constraints(tab, 0, tab->n_eq + i + 1) < 0)
3222 0 : return isl_basic_map_free(bmap);
3223 66712 : if (rotate_constraints(tab, tab->n_eq + i + 1,
3224 66712 : bmap->n_ineq - i) < 0)
3225 0 : return isl_basic_map_free(bmap);
3226 66712 : tab->n_eq++;
3227 : }
3228 :
3229 1378155723 : if (n_eq != tab->n_eq)
3230 25902 : bmap = gauss_if_shared(bmap, tab);
3231 :
3232 1378155723 : return bmap;
3233 : }
3234 :
3235 6475077076 : static int con_is_redundant(struct isl_tab *tab, struct isl_tab_var *var)
3236 : {
3237 6475077076 : if (!tab)
3238 0 : return -1;
3239 6475077076 : if (tab->rational) {
3240 5272141592 : int sgn = sign_of_min(tab, var);
3241 5272141592 : if (sgn < -1)
3242 0 : return -1;
3243 5272141592 : return sgn >= 0;
3244 : } else {
3245 1202935484 : int irred = isl_tab_min_at_most_neg_one(tab, var);
3246 1202935484 : if (irred < 0)
3247 0 : return -1;
3248 1202935484 : 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 1669357177 : int isl_tab_detect_redundant(struct isl_tab *tab)
3266 : {
3267 : int i;
3268 : unsigned n_marked;
3269 :
3270 1669357177 : if (!tab)
3271 0 : return -1;
3272 1669357177 : if (tab->empty)
3273 609861 : return 0;
3274 1668747316 : if (tab->n_redundant == tab->n_row)
3275 496278 : return 0;
3276 :
3277 1668251038 : n_marked = 0;
3278 6771911214 : for (i = tab->n_redundant; i < tab->n_row; ++i) {
3279 5103660176 : struct isl_tab_var *var = isl_tab_var_from_row(tab, i);
3280 5103660176 : var->marked = !var->frozen && var->is_nonneg;
3281 5103660176 : if (var->marked)
3282 887690107 : n_marked++;
3283 : }
3284 3188439283 : for (i = tab->n_dead; i < tab->n_col; ++i) {
3285 1520188245 : struct isl_tab_var *var = var_from_col(tab, i);
3286 1954275808 : var->marked = !var->frozen && var->is_nonneg &&
3287 434087563 : !min_is_manifestly_unbounded(tab, var);
3288 1520188245 : if (var->marked)
3289 222418972 : n_marked++;
3290 : }
3291 4359122913 : while (n_marked) {
3292 : struct isl_tab_var *var;
3293 : int red;
3294 1049248691 : var = select_marked(tab);
3295 1049248691 : if (!var)
3296 26627854 : break;
3297 1022620837 : var->marked = 0;
3298 1022620837 : n_marked--;
3299 1022620837 : red = con_is_redundant(tab, var);
3300 1022620837 : if (red < 0)
3301 0 : return -1;
3302 1022620837 : if (red && !var->is_redundant)
3303 24931177 : if (isl_tab_mark_redundant(tab, var->index) < 0)
3304 0 : return -1;
3305 7064242055 : for (i = tab->n_dead; i < tab->n_col; ++i) {
3306 6041621218 : var = var_from_col(tab, i);
3307 6041621218 : if (!var->marked)
3308 5399633706 : continue;
3309 641987512 : if (!min_is_manifestly_unbounded(tab, var))
3310 593978933 : continue;
3311 48008579 : var->marked = 0;
3312 48008579 : n_marked--;
3313 : }
3314 : }
3315 :
3316 1668251038 : return 0;
3317 : }
3318 :
3319 3453123136 : int isl_tab_is_equality(struct isl_tab *tab, int con)
3320 : {
3321 : int row;
3322 : unsigned off;
3323 :
3324 3453123136 : if (!tab)
3325 0 : return -1;
3326 3453123136 : if (tab->con[con].is_zero)
3327 652210478 : return 1;
3328 2800912658 : if (tab->con[con].is_redundant)
3329 342613299 : return 0;
3330 2458299359 : if (!tab->con[con].is_row)
3331 1304720533 : return tab->con[con].index < tab->n_dead;
3332 :
3333 1153578826 : row = tab->con[con].index;
3334 :
3335 1153578826 : off = 2 + tab->M;
3336 3832613984 : return isl_int_is_zero(tab->mat->row[row][1]) &&
3337 1897703200 : !row_is_big(tab, row) &&
3338 372062187 : isl_seq_first_non_zero(tab->mat->row[row] + off + tab->n_dead,
3339 372062187 : 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 2698376716 : 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 2698376716 : enum isl_lp_result res = isl_lp_ok;
3358 : struct isl_tab_var *var;
3359 : struct isl_tab_undo *snap;
3360 :
3361 2698376716 : if (!tab)
3362 0 : return isl_lp_error;
3363 :
3364 2698376716 : if (tab->empty)
3365 0 : return isl_lp_empty;
3366 :
3367 2698376716 : snap = isl_tab_snap(tab);
3368 2698376716 : r = isl_tab_add_row(tab, f);
3369 2698376716 : if (r < 0)
3370 0 : return isl_lp_error;
3371 2698376716 : var = &tab->con[r];
3372 : for (;;) {
3373 : int row, col;
3374 6150125946 : find_pivot(tab, var, var, -1, &row, &col);
3375 6150125946 : if (row == var->index) {
3376 1964545 : res = isl_lp_unbounded;
3377 2700341261 : break;
3378 : }
3379 6148161401 : if (row == -1)
3380 2696412171 : break;
3381 3451749230 : if (isl_tab_pivot(tab, row, col) < 0)
3382 0 : return isl_lp_error;
3383 3451749230 : }
3384 2698376716 : isl_int_mul(tab->mat->row[var->index][0],
3385 : tab->mat->row[var->index][0], denom);
3386 2698376716 : if (ISL_FL_ISSET(flags, ISL_TAB_SAVE_DUAL)) {
3387 : int i;
3388 :
3389 440289371 : isl_vec_free(tab->dual);
3390 440289371 : tab->dual = isl_vec_alloc(tab->mat->ctx, 1 + tab->n_con);
3391 440289371 : if (!tab->dual)
3392 0 : return isl_lp_error;
3393 440289371 : isl_int_set(tab->dual->el[0], tab->mat->row[var->index][0]);
3394 12731020035 : for (i = 0; i < tab->n_con; ++i) {
3395 : int pos;
3396 12290730664 : if (tab->con[i].is_row) {
3397 7591124743 : isl_int_set_si(tab->dual->el[1 + i], 0);
3398 7591124743 : continue;
3399 : }
3400 4699605921 : pos = 2 + tab->M + tab->con[i].index;
3401 4699605921 : if (tab->con[i].negated)
3402 721053089 : isl_int_neg(tab->dual->el[1 + i],
3403 : tab->mat->row[var->index][pos]);
3404 : else
3405 3978552832 : isl_int_set(tab->dual->el[1 + i],
3406 : tab->mat->row[var->index][pos]);
3407 : }
3408 : }
3409 2698376716 : if (opt && res == isl_lp_ok) {
3410 2696412171 : if (opt_denom) {
3411 1740993626 : isl_int_set(*opt, tab->mat->row[var->index][1]);
3412 1740993626 : isl_int_set(*opt_denom, tab->mat->row[var->index][0]);
3413 : } else
3414 955418545 : get_rounded_sample_value(tab, var, 1, opt);
3415 : }
3416 2698376716 : if (isl_tab_rollback(tab, snap) < 0)
3417 0 : return isl_lp_error;
3418 2698376716 : 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 4020564122 : int isl_tab_is_redundant(struct isl_tab *tab, int con)
3430 : {
3431 4020564122 : if (isl_tab_check_con(tab, con) < 0)
3432 0 : return -1;
3433 4020564122 : if (tab->con[con].is_zero)
3434 259212 : return 0;
3435 4020304910 : if (tab->con[con].is_redundant)
3436 1251501308 : return 1;
3437 2768803602 : 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 4254434 : static isl_bool is_constant(struct isl_tab *tab, struct isl_tab_var *var,
3449 : isl_int *value)
3450 : {
3451 4254434 : unsigned off = 2 + tab->M;
3452 4254434 : isl_mat *mat = tab->mat;
3453 : int n;
3454 : int row;
3455 : int pos;
3456 :
3457 4254434 : if (!var->is_row)
3458 1562631 : return isl_bool_false;
3459 2691803 : row = var->index;
3460 2691803 : if (row_is_big(tab, row))
3461 0 : return isl_bool_false;
3462 2691803 : n = tab->n_col - tab->n_dead;
3463 2691803 : pos = isl_seq_first_non_zero(mat->row[row] + off + tab->n_dead, n);
3464 2691803 : if (pos != -1)
3465 2062232 : return isl_bool_false;
3466 629571 : if (value)
3467 0 : isl_int_divexact(*value, mat->row[row][1], mat->row[row][0]);
3468 629571 : 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 5605758 : static int reached(struct isl_tab *tab, struct isl_tab_var *var, int sgn,
3482 : isl_int target, isl_int *tmp)
3483 : {
3484 5605758 : if (row_is_big(tab, var->index))
3485 0 : return 1;
3486 5605758 : isl_int_mul(*tmp, tab->mat->row[var->index][0], target);
3487 5605758 : if (sgn > 0)
3488 1769064 : return isl_int_ge(tab->mat->row[var->index][1], *tmp);
3489 : else
3490 3836694 : 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 4610258 : 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 4610258 : if (sgn < 0 && min_is_manifestly_unbounded(tab, var))
3515 890926 : return isl_bool_true;
3516 3719332 : if (sgn > 0 && max_is_manifestly_unbounded(tab, var))
3517 0 : return isl_bool_true;
3518 3719332 : if (to_row(tab, var, sgn) < 0)
3519 0 : return isl_bool_error;
3520 9325090 : while (!reached(tab, var, sgn, target, tmp)) {
3521 4563004 : find_pivot(tab, var, var, sgn, &row, &col);
3522 4563004 : if (row == -1)
3523 1252644 : return isl_bool_false;
3524 3310360 : if (row == var->index)
3525 1423934 : return isl_bool_true;
3526 1886426 : if (isl_tab_pivot(tab, row, col) < 0)
3527 0 : return isl_bool_error;
3528 : }
3529 :
3530 1042754 : 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 3624863 : 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 3624863 : get_rounded_sample_value(tab, var, -1, target);
3561 3624863 : isl_int_sub_ui(*target, *target, 1);
3562 3624863 : reached = var_reaches(tab, var, -1, *target, tmp);
3563 3624863 : if (reached < 0 || reached)
3564 2639468 : return isl_bool_not(reached);
3565 985395 : get_rounded_sample_value(tab, var, 1, target);
3566 985395 : isl_int_add_ui(*target, *target, 1);
3567 985395 : reached = var_reaches(tab, var, 1, *target, tmp);
3568 985395 : if (reached < 0 || reached)
3569 718146 : return isl_bool_not(reached);
3570 267249 : get_rounded_sample_value(tab, var, -1, tmp);
3571 267249 : isl_int_sub_ui(*target, *target, 1);
3572 267249 : 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 267249 : if (isl_tab_extend_cons(tab, 1) < 0)
3579 0 : return isl_bool_error;
3580 267249 : eq = isl_vec_alloc(isl_tab_get_ctx(tab), 1 + tab->n_var);
3581 267249 : if (!eq)
3582 0 : return isl_bool_error;
3583 267249 : pos = var - tab->var;
3584 267249 : isl_seq_clr(eq->el + 1, tab->n_var);
3585 267249 : isl_int_set_si(eq->el[1 + pos], -1);
3586 267249 : isl_int_set(eq->el[0], *target);
3587 267249 : r = isl_tab_add_eq(tab, eq->el);
3588 267249 : isl_vec_free(eq);
3589 :
3590 267249 : 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 4254434 : 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 4254434 : if (var->is_row && row_is_big(tab, var->index))
3613 0 : return isl_bool_false;
3614 4254434 : is_cst = is_constant(tab, var, value);
3615 4254434 : if (is_cst < 0 || is_cst)
3616 629571 : return is_cst;
3617 :
3618 3624863 : if (!value)
3619 2376618 : isl_int_init(target);
3620 3624863 : isl_int_init(tmp);
3621 :
3622 3624863 : is_cst = detect_constant_with_tmp(tab, var,
3623 : value ? value : &target, &tmp);
3624 :
3625 3624863 : isl_int_clear(tmp);
3626 3624863 : if (!value)
3627 2376618 : isl_int_clear(target);
3628 :
3629 3624863 : 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 1248245 : isl_bool isl_tab_is_constant(struct isl_tab *tab, int var, isl_int *value)
3640 : {
3641 1248245 : if (!tab)
3642 0 : return isl_bool_error;
3643 1248245 : 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 1248245 : if (tab->rational)
3647 0 : return isl_bool_false;
3648 :
3649 1248245 : 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 505574 : isl_stat isl_tab_detect_constants(struct isl_tab *tab)
3659 : {
3660 : int i;
3661 :
3662 505574 : if (!tab)
3663 0 : return isl_stat_error;
3664 505574 : if (tab->rational)
3665 0 : return isl_stat_ok;
3666 :
3667 3511763 : for (i = 0; i < tab->n_var; ++i) {
3668 3006189 : if (get_constant(tab, &tab->var[i], NULL) < 0)
3669 0 : return isl_stat_error;
3670 : }
3671 :
3672 505574 : 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 18230001677 : struct isl_tab_undo *isl_tab_snap(struct isl_tab *tab)
3679 : {
3680 18230001677 : if (!tab)
3681 0 : return NULL;
3682 18230001677 : tab->need_undo = 1;
3683 18230001677 : 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 3804503 : isl_bool isl_tab_need_undo(struct isl_tab *tab)
3690 : {
3691 3804503 : if (!tab)
3692 0 : return isl_bool_error;
3693 :
3694 3804503 : 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 3804503 : void isl_tab_clear_undo(struct isl_tab *tab)
3703 : {
3704 3804503 : if (!tab)
3705 0 : return;
3706 :
3707 3804503 : free_undo(tab);
3708 3804503 : 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 1570805 : static isl_stat unrelax(struct isl_tab *tab, struct isl_tab_var *var)
3716 : {
3717 1570805 : unsigned off = 2 + tab->M;
3718 :
3719 1570805 : if (!var->is_row && !max_is_manifestly_unbounded(tab, var))
3720 205947 : if (to_row(tab, var, 1) < 0)
3721 0 : return isl_stat_error;
3722 :
3723 1570805 : if (var->is_row) {
3724 242651 : isl_int_sub(tab->mat->row[var->index][1],
3725 : tab->mat->row[var->index][1], tab->mat->row[var->index][0]);
3726 242651 : if (var->is_nonneg) {
3727 242651 : int sgn = restore_row(tab, var);
3728 242651 : isl_assert(tab->mat->ctx, sgn >= 0,
3729 : return isl_stat_error);
3730 : }
3731 : } else {
3732 : int i;
3733 :
3734 6293380 : for (i = 0; i < tab->n_row; ++i) {
3735 4965226 : if (isl_int_is_zero(tab->mat->row[i][off + var->index]))
3736 3232510 : continue;
3737 1732716 : 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 1570805 : 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 502242 : static isl_stat ununrestrict(struct isl_tab *tab, struct isl_tab_var *var)
3752 : {
3753 502242 : var->is_nonneg = 1;
3754 :
3755 502242 : if (var->is_row && restore_row(tab, var) < -1)
3756 0 : return isl_stat_error;
3757 :
3758 502242 : 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 10585139006 : static isl_stat restore_last_redundant(struct isl_tab *tab)
3769 : {
3770 : struct isl_tab_var *var;
3771 :
3772 10585139006 : 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 10585139006 : var = isl_tab_var_from_row(tab, tab->n_redundant - 1);
3777 10585139006 : var->is_redundant = 0;
3778 10585139006 : tab->n_redundant--;
3779 10585139006 : restore_row(tab, var);
3780 :
3781 10585139006 : 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 95379736051 : static isl_stat perform_undo_var(struct isl_tab *tab, struct isl_tab_undo *undo)
3787 : {
3788 95379736051 : struct isl_tab_var *var = var_from_index(tab, undo->u.var_index);
3789 95379736051 : switch (undo->type) {
3790 : case isl_tab_undo_nonneg:
3791 18619936469 : var->is_nonneg = 0;
3792 18619936469 : break;
3793 : case isl_tab_undo_redundant:
3794 10410779973 : 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 10410779973 : return restore_last_redundant(tab);
3799 : case isl_tab_undo_freeze:
3800 36909013347 : var->frozen = 0;
3801 36909013347 : break;
3802 : case isl_tab_undo_zero:
3803 648122451 : var->is_zero = 0;
3804 648122451 : if (!var->is_row)
3805 644229133 : tab->n_dead--;
3806 648122451 : break;
3807 : case isl_tab_undo_allocate:
3808 28789810764 : if (undo->u.var_index >= 0) {
3809 883582966 : isl_assert(tab->mat->ctx, !var->is_row,
3810 : return isl_stat_error);
3811 883582966 : return drop_col(tab, var->index);
3812 : }
3813 27906227798 : if (!var->is_row) {
3814 3219682451 : if (!max_is_manifestly_unbounded(tab, var)) {
3815 1710158858 : if (to_row(tab, var, 1) < 0)
3816 0 : return isl_stat_error;
3817 1509523593 : } else if (!min_is_manifestly_unbounded(tab, var)) {
3818 418077406 : if (to_row(tab, var, -1) < 0)
3819 0 : return isl_stat_error;
3820 : } else
3821 1091446187 : if (to_row(tab, var, 0) < 0)
3822 0 : return isl_stat_error;
3823 : }
3824 27906227798 : return drop_row(tab, var->index);
3825 : case isl_tab_undo_relax:
3826 1570805 : return unrelax(tab, var);
3827 : case isl_tab_undo_unrestrict:
3828 502242 : 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 56177072267 : 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 140040870 : isl_stat isl_tab_restore_redundant(struct isl_tab *tab)
3846 : {
3847 140040870 : if (!tab)
3848 0 : return isl_stat_error;
3849 :
3850 140040870 : 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 454440773 : while (tab->n_redundant > 0) {
3857 174359033 : if (tab->row_var[tab->n_redundant - 1] >= 0) {
3858 : struct isl_tab_var *var;
3859 :
3860 149290217 : var = isl_tab_var_from_row(tab, tab->n_redundant - 1);
3861 149290217 : var->is_nonneg = 0;
3862 : }
3863 174359033 : restore_last_redundant(tab);
3864 : }
3865 140040870 : 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 882100037 : static isl_stat drop_bmap_div(struct isl_tab *tab, int pos)
3872 : {
3873 : int off;
3874 :
3875 882100037 : off = tab->n_var - isl_basic_map_dim(tab->bmap, isl_dim_div);
3876 882100037 : if (isl_basic_map_drop_div(tab->bmap, pos - off) < 0)
3877 : return isl_stat_error;
3878 882100037 : if (tab->samples) {
3879 13372331 : tab->samples = isl_mat_drop_cols(tab->samples, 1 + pos, 1);
3880 13372331 : if (!tab->samples)
3881 0 : return isl_stat_error;
3882 : }
3883 :
3884 882100037 : 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 34292089 : static void drop_samples_since(struct isl_tab *tab, int n)
3948 : {
3949 : int i;
3950 :
3951 116784211 : for (i = tab->n_sample - 1; i >= 0 && tab->n_sample > n; --i) {
3952 82492122 : if (tab->sample_index[i] < n)
3953 14969987 : continue;
3954 :
3955 67522135 : if (i != tab->n_sample - 1) {
3956 31057326 : int t = tab->sample_index[tab->n_sample-1];
3957 31057326 : tab->sample_index[tab->n_sample-1] = tab->sample_index[i];
3958 31057326 : tab->sample_index[i] = t;
3959 31057326 : isl_mat_swap_rows(tab->samples, tab->n_sample-1, i);
3960 : }
3961 67522135 : tab->n_sample--;
3962 : }
3963 34292089 : }
3964 :
3965 : static isl_stat perform_undo(struct isl_tab *tab, struct isl_tab_undo *undo)
3966 : WARN_UNUSED;
3967 >11484*10^7 : static isl_stat perform_undo(struct isl_tab *tab, struct isl_tab_undo *undo)
3968 : {
3969 >11484*10^7 : switch (undo->type) {
3970 : case isl_tab_undo_rational:
3971 31830 : tab->rational = 0;
3972 31830 : break;
3973 : case isl_tab_undo_empty:
3974 2716560881 : tab->empty = 0;
3975 2716560881 : 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 95379736051 : 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 15687782760 : return isl_basic_map_free_inequality(tab->bmap, 1);
3988 : case isl_tab_undo_bmap_div:
3989 882100037 : 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 120663085 : tab->n_outside--;
3996 120663085 : break;
3997 : case isl_tab_undo_saved_samples:
3998 34292089 : drop_samples_since(tab, undo->u.n);
3999 34292089 : break;
4000 : case isl_tab_undo_callback:
4001 27534440 : return undo->u.callback->run(undo->u.callback);
4002 : default:
4003 0 : isl_assert(tab->mat->ctx, 0, return isl_stat_error);
4004 : }
4005 2871547885 : 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 15518787686 : int isl_tab_rollback(struct isl_tab *tab, struct isl_tab_undo *snap)
4012 : {
4013 : struct isl_tab_undo *undo, *next;
4014 :
4015 15518787686 : if (!tab)
4016 0 : return -1;
4017 :
4018 15518787686 : tab->in_undo = 1;
4019 >13036*10^7 : for (undo = tab->top; undo && undo != &tab->bottom; undo = next) {
4020 >11951*10^7 : next = undo->next;
4021 >11951*10^7 : if (undo == snap)
4022 4665004781 : break;
4023 >11484*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 >11484*10^7 : free_undo_record(undo);
4030 : }
4031 15518787686 : tab->in_undo = 0;
4032 15518787686 : tab->top = undo;
4033 15518787686 : if (!undo)
4034 0 : return -1;
4035 15518787686 : 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 3376365398 : static enum isl_ineq_type separation_type(struct isl_tab *tab, unsigned row)
4049 : {
4050 : int pos;
4051 3376365398 : unsigned off = 2 + tab->M;
4052 :
4053 3376365398 : if (tab->rational)
4054 3313470690 : return isl_ineq_separate;
4055 :
4056 62894708 : if (!isl_int_is_one(tab->mat->row[row][0]))
4057 2098916 : return isl_ineq_separate;
4058 :
4059 60795792 : pos = isl_seq_first_non_zero(tab->mat->row[row] + off + tab->n_dead,
4060 60795792 : tab->n_col - tab->n_dead);
4061 60795792 : if (pos == -1) {
4062 44976463 : if (isl_int_is_negone(tab->mat->row[row][1]))
4063 9943929 : return isl_ineq_adj_eq;
4064 : else
4065 35032534 : return isl_ineq_separate;
4066 : }
4067 :
4068 15819329 : if (!isl_int_eq(tab->mat->row[row][1],
4069 : tab->mat->row[row][off + tab->n_dead + pos]))
4070 4972012 : return isl_ineq_separate;
4071 :
4072 21694634 : pos = isl_seq_first_non_zero(
4073 10847317 : tab->mat->row[row] + off + tab->n_dead + pos + 1,
4074 10847317 : tab->n_col - tab->n_dead - pos - 1);
4075 :
4076 10847317 : 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 8856371156 : enum isl_ineq_type isl_tab_ineq_type(struct isl_tab *tab, isl_int *ineq)
4088 : {
4089 8856371156 : enum isl_ineq_type type = isl_ineq_error;
4090 8856371156 : struct isl_tab_undo *snap = NULL;
4091 : int con;
4092 : int row;
4093 :
4094 8856371156 : if (!tab)
4095 0 : return isl_ineq_error;
4096 :
4097 8856371156 : if (isl_tab_extend_cons(tab, 1) < 0)
4098 0 : return isl_ineq_error;
4099 :
4100 8856371156 : snap = isl_tab_snap(tab);
4101 :
4102 8856371156 : con = isl_tab_add_row(tab, ineq);
4103 8856371156 : if (con < 0)
4104 0 : goto error;
4105 :
4106 8856371156 : row = tab->con[con].index;
4107 8856371156 : if (isl_tab_row_is_redundant(tab, row))
4108 0 : type = isl_ineq_redundant;
4109 12260929898 : else if (isl_int_is_neg(tab->mat->row[row][1]) &&
4110 3493782902 : (tab->rational ||
4111 89224160 : isl_int_abs_ge(tab->mat->row[row][1],
4112 3403914917 : tab->mat->row[row][0]))) {
4113 3403914917 : int nonneg = at_least_zero(tab, &tab->con[con]);
4114 3403914917 : if (nonneg < 0)
4115 0 : goto error;
4116 3403914917 : if (nonneg)
4117 27549519 : type = isl_ineq_cut;
4118 : else
4119 3376365398 : type = separation_type(tab, row);
4120 : } else {
4121 5452456239 : int red = con_is_redundant(tab, &tab->con[con]);
4122 5452456239 : if (red < 0)
4123 0 : goto error;
4124 5452456239 : if (!red)
4125 39691959 : type = isl_ineq_cut;
4126 : else
4127 5412764280 : type = isl_ineq_redundant;
4128 : }
4129 :
4130 8856371156 : if (isl_tab_rollback(tab, snap))
4131 0 : return isl_ineq_error;
4132 8856371156 : return type;
4133 : error:
4134 0 : return isl_ineq_error;
4135 : }
4136 :
4137 121974939 : isl_stat isl_tab_track_bmap(struct isl_tab *tab, __isl_take isl_basic_map *bmap)
4138 : {
4139 121974939 : bmap = isl_basic_map_cow(bmap);
4140 121974939 : if (!tab || !bmap)
4141 : goto error;
4142 :
4143 121974939 : if (tab->empty) {
4144 31620957 : bmap = isl_basic_map_set_to_empty(bmap);
4145 31620957 : if (!bmap)
4146 0 : goto error;
4147 31620957 : tab->bmap = bmap;
4148 31620957 : return isl_stat_ok;
4149 : }
4150 :
4151 90353982 : isl_assert(tab->mat->ctx, tab->n_eq == bmap->n_eq, goto error);
4152 90353982 : isl_assert(tab->mat->ctx,
4153 : tab->n_con == bmap->n_eq + bmap->n_ineq, goto error);
4154 :
4155 90353982 : tab->bmap = bmap;
4156 :
4157 90353982 : return isl_stat_ok;
4158 : error:
4159 0 : isl_basic_map_free(bmap);
4160 0 : return isl_stat_error;
4161 : }
4162 :
4163 15199472 : isl_stat isl_tab_track_bset(struct isl_tab *tab, __isl_take isl_basic_set *bset)
4164 : {
4165 15199472 : return isl_tab_track_bmap(tab, bset_to_bmap(bset));
4166 : }
4167 :
4168 73097319 : __isl_keep isl_basic_set *isl_tab_peek_bset(struct isl_tab *tab)
4169 : {
4170 73097319 : if (!tab)
4171 0 : return NULL;
4172 :
4173 73097319 : 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 : }
|
