?  search:
summary | log | graphiclog1 | graphiclog2 | commit | commitdiff | tree | refs | edit | fork
blame | history | raw | HEAD
isl_tab_basic_map_partial_lexopt: copy integer divisions to context
[isl.git] / isl_mat.c
blob5ec99b5a25fccdbebe757bcb5c7979cbebbe33cf
1 /*
2 * Copyright 2008-2009 Katholieke Universiteit Leuven
3 * Copyright 2010 INRIA Saclay
4 * Copyright 2014 Ecole Normale Superieure
5 * Copyright 2017 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 INRIA Saclay - Ile-de-France, Parc Club Orsay Universite,
12 * ZAC des vignes, 4 rue Jacques Monod, 91893 Orsay, France
13 * and Ecole Normale Superieure, 45 rue d'Ulm, 75230 Paris, France
14 */
15
16 #include <isl_ctx_private.h>
17 #include <isl_map_private.h>
18 #include <isl/space.h>
19 #include <isl_seq.h>
20 #include <isl_mat_private.h>
21 #include <isl_vec_private.h>
22 #include <isl_space_private.h>
23 #include <isl_val_private.h>
24
25 isl_ctx *isl_mat_get_ctx(__isl_keep isl_mat *mat)
26 {
27 return mat ? mat->ctx : NULL;
28 }
29
30 /* Return a hash value that digests "mat".
31 */
32 uint32_t isl_mat_get_hash(__isl_keep isl_mat *mat)
33 {
34 int i;
35 uint32_t hash;
36
37 if (!mat)
38 return 0;
39
40 hash = isl_hash_init();
41 isl_hash_byte(hash, mat->n_row & 0xFF);
42 isl_hash_byte(hash, mat->n_col & 0xFF);
43 for (i = 0; i < mat->n_row; ++i) {
44 uint32_t row_hash;
45
46 row_hash = isl_seq_get_hash(mat->row[i], mat->n_col);
47 isl_hash_hash(hash, row_hash);
48 }
49
50 return hash;
51 }
52
53 __isl_give isl_mat *isl_mat_alloc(isl_ctx *ctx,
54 unsigned n_row, unsigned n_col)
55 {
56 int i;
57 struct isl_mat *mat;
58
59 mat = isl_alloc_type(ctx, struct isl_mat);
60 if (!mat)
61 return NULL;
62
63 mat->row = NULL;
64 mat->block = isl_blk_alloc(ctx, n_row * n_col);
65 if (isl_blk_is_error(mat->block))
66 goto error;
67 mat->row = isl_calloc_array(ctx, isl_int *, n_row);
68 if (n_row && !mat->row)
69 goto error;
70
71 if (n_col != 0) {
72 for (i = 0; i < n_row; ++i)
73 mat->row[i] = mat->block.data + i * n_col;
74 }
75
76 mat->ctx = ctx;
77 isl_ctx_ref(ctx);
78 mat->ref = 1;
79 mat->n_row = n_row;
80 mat->n_col = n_col;
81 mat->max_col = n_col;
82 mat->flags = 0;
83
84 return mat;
85 error:
86 isl_blk_free(ctx, mat->block);
87 free(mat);
88 return NULL;
89 }
90
91 __isl_give isl_mat *isl_mat_extend(__isl_take isl_mat *mat,
92 unsigned n_row, unsigned n_col)
93 {
94 int i;
95 isl_int *old;
96 isl_int **row;
97
98 if (!mat)
99 return NULL;
100
101 if (mat->max_col >= n_col && mat->n_row >= n_row) {
102 if (mat->n_col < n_col)
103 mat->n_col = n_col;
104 return mat;
105 }
106
107 if (mat->max_col < n_col) {
108 struct isl_mat *new_mat;
109
110 if (n_row < mat->n_row)
111 n_row = mat->n_row;
112 new_mat = isl_mat_alloc(mat->ctx, n_row, n_col);
113 if (!new_mat)
114 goto error;
115 for (i = 0; i < mat->n_row; ++i)
116 isl_seq_cpy(new_mat->row[i], mat->row[i], mat->n_col);
117 isl_mat_free(mat);
118 return new_mat;
119 }
120
121 mat = isl_mat_cow(mat);
122 if (!mat)
123 goto error;
124
125 old = mat->block.data;
126 mat->block = isl_blk_extend(mat->ctx, mat->block, n_row * mat->max_col);
127 if (isl_blk_is_error(mat->block))
128 goto error;
129 row = isl_realloc_array(mat->ctx, mat->row, isl_int *, n_row);
130 if (n_row && !row)
131 goto error;
132 mat->row = row;
133
134 for (i = 0; i < mat->n_row; ++i)
135 mat->row[i] = mat->block.data + (mat->row[i] - old);
136 for (i = mat->n_row; i < n_row; ++i)
137 mat->row[i] = mat->block.data + i * mat->max_col;
138 mat->n_row = n_row;
139 if (mat->n_col < n_col)
140 mat->n_col = n_col;
141
142 return mat;
143 error:
144 isl_mat_free(mat);
145 return NULL;
146 }
147
148 __isl_give isl_mat *isl_mat_sub_alloc6(isl_ctx *ctx, isl_int **row,
149 unsigned first_row, unsigned n_row, unsigned first_col, unsigned n_col)
150 {
151 int i;
152 struct isl_mat *mat;
153
154 mat = isl_alloc_type(ctx, struct isl_mat);
155 if (!mat)
156 return NULL;
157 mat->row = isl_alloc_array(ctx, isl_int *, n_row);
158 if (n_row && !mat->row)
159 goto error;
160 for (i = 0; i < n_row; ++i)
161 mat->row[i] = row[first_row+i] + first_col;
162 mat->ctx = ctx;
163 isl_ctx_ref(ctx);
164 mat->ref = 1;
165 mat->n_row = n_row;
166 mat->n_col = n_col;
167 mat->block = isl_blk_empty();
168 mat->flags = ISL_MAT_BORROWED;
169 return mat;
170 error:
171 free(mat);
172 return NULL;
173 }
174
175 __isl_give isl_mat *isl_mat_sub_alloc(__isl_keep isl_mat *mat,
176 unsigned first_row, unsigned n_row, unsigned first_col, unsigned n_col)
177 {
178 if (!mat)
179 return NULL;
180 return isl_mat_sub_alloc6(mat->ctx, mat->row, first_row, n_row,
181 first_col, n_col);
182 }
183
184 void isl_mat_sub_copy(struct isl_ctx *ctx, isl_int **dst, isl_int **src,
185 unsigned n_row, unsigned dst_col, unsigned src_col, unsigned n_col)
186 {
187 int i;
188
189 for (i = 0; i < n_row; ++i)
190 isl_seq_cpy(dst[i]+dst_col, src[i]+src_col, n_col);
191 }
192
193 void isl_mat_sub_neg(struct isl_ctx *ctx, isl_int **dst, isl_int **src,
194 unsigned n_row, unsigned dst_col, unsigned src_col, unsigned n_col)
195 {
196 int i;
197
198 for (i = 0; i < n_row; ++i)
199 isl_seq_neg(dst[i]+dst_col, src[i]+src_col, n_col);
200 }
201
202 __isl_give isl_mat *isl_mat_copy(__isl_keep isl_mat *mat)
203 {
204 if (!mat)
205 return NULL;
206
207 mat->ref++;
208 return mat;
209 }
210
211 __isl_give isl_mat *isl_mat_dup(__isl_keep isl_mat *mat)
212 {
213 int i;
214 struct isl_mat *mat2;
215
216 if (!mat)
217 return NULL;
218 mat2 = isl_mat_alloc(mat->ctx, mat->n_row, mat->n_col);
219 if (!mat2)
220 return NULL;
221 for (i = 0; i < mat->n_row; ++i)
222 isl_seq_cpy(mat2->row[i], mat->row[i], mat->n_col);
223 return mat2;
224 }
225
226 __isl_give isl_mat *isl_mat_cow(__isl_take isl_mat *mat)
227 {
228 struct isl_mat *mat2;
229 if (!mat)
230 return NULL;
231
232 if (mat->ref == 1 && !ISL_F_ISSET(mat, ISL_MAT_BORROWED))
233 return mat;
234
235 mat2 = isl_mat_dup(mat);
236 isl_mat_free(mat);
237 return mat2;
238 }
239
240 __isl_null isl_mat *isl_mat_free(__isl_take isl_mat *mat)
241 {
242 if (!mat)
243 return NULL;
244
245 if (--mat->ref > 0)
246 return NULL;
247
248 if (!ISL_F_ISSET(mat, ISL_MAT_BORROWED))
249 isl_blk_free(mat->ctx, mat->block);
250 isl_ctx_deref(mat->ctx);
251 free(mat->row);
252 free(mat);
253
254 return NULL;
255 }
256
257 isl_size isl_mat_rows(__isl_keep isl_mat *mat)
258 {
259 return mat ? mat->n_row : isl_size_error;
260 }
261
262 isl_size isl_mat_cols(__isl_keep isl_mat *mat)
263 {
264 return mat ? mat->n_col : isl_size_error;
265 }
266
267 /* Check that "col" is a valid column position for "mat".
268 */
269 static isl_stat check_col(__isl_keep isl_mat *mat, int col)
270 {
271 if (!mat)
272 return isl_stat_error;
273 if (col < 0 || col >= mat->n_col)
274 isl_die(isl_mat_get_ctx(mat), isl_error_invalid,
275 "column out of range", return isl_stat_error);
276 return isl_stat_ok;
277 }
278
279 /* Check that "row" is a valid row position for "mat".
280 */
281 static isl_stat check_row(__isl_keep isl_mat *mat, int row)
282 {
283 if (!mat)
284 return isl_stat_error;
285 if (row < 0 || row >= mat->n_row)
286 isl_die(isl_mat_get_ctx(mat), isl_error_invalid,
287 "row out of range", return isl_stat_error);
288 return isl_stat_ok;
289 }
290
291 /* Check that there are "n" columns starting at position "first" in "mat".
292 */
293 static isl_stat check_col_range(__isl_keep isl_mat *mat, unsigned first,
294 unsigned n)
295 {
296 if (!mat)
297 return isl_stat_error;
298 if (first + n > mat->n_col || first + n < first)
299 isl_die(isl_mat_get_ctx(mat), isl_error_invalid,
300 "column position or range out of bounds",
301 return isl_stat_error);
302 return isl_stat_ok;
303 }
304
305 /* Check that there are "n" rows starting at position "first" in "mat".
306 */
307 static isl_stat check_row_range(__isl_keep isl_mat *mat, unsigned first,
308 unsigned n)
309 {
310 if (!mat)
311 return isl_stat_error;
312 if (first + n > mat->n_row || first + n < first)
313 isl_die(isl_mat_get_ctx(mat), isl_error_invalid,
314 "row position or range out of bounds",
315 return isl_stat_error);
316 return isl_stat_ok;
317 }
318
319 int isl_mat_get_element(__isl_keep isl_mat *mat, int row, int col, isl_int *v)
320 {
321 if (check_row(mat, row) < 0)
322 return -1;
323 if (check_col(mat, col) < 0)
324 return -1;
325 isl_int_set(*v, mat->row[row][col]);
326 return 0;
327 }
328
329 /* Extract the element at row "row", oolumn "col" of "mat".
330 */
331 __isl_give isl_val *isl_mat_get_element_val(__isl_keep isl_mat *mat,
332 int row, int col)
333 {
334 isl_ctx *ctx;
335
336 if (check_row(mat, row) < 0)
337 return NULL;
338 if (check_col(mat, col) < 0)
339 return NULL;
340 ctx = isl_mat_get_ctx(mat);
341 return isl_val_int_from_isl_int(ctx, mat->row[row][col]);
342 }
343
344 __isl_give isl_mat *isl_mat_set_element(__isl_take isl_mat *mat,
345 int row, int col, isl_int v)
346 {
347 mat = isl_mat_cow(mat);
348 if (check_row(mat, row) < 0)
349 return isl_mat_free(mat);
350 if (check_col(mat, col) < 0)
351 return isl_mat_free(mat);
352 isl_int_set(mat->row[row][col], v);
353 return mat;
354 }
355
356 __isl_give isl_mat *isl_mat_set_element_si(__isl_take isl_mat *mat,
357 int row, int col, int v)
358 {
359 mat = isl_mat_cow(mat);
360 if (check_row(mat, row) < 0)
361 return isl_mat_free(mat);
362 if (check_col(mat, col) < 0)
363 return isl_mat_free(mat);
364 isl_int_set_si(mat->row[row][col], v);
365 return mat;
366 }
367
368 /* Replace the element at row "row", column "col" of "mat" by "v".
369 */
370 __isl_give isl_mat *isl_mat_set_element_val(__isl_take isl_mat *mat,
371 int row, int col, __isl_take isl_val *v)
372 {
373 if (!v)
374 return isl_mat_free(mat);
375 if (!isl_val_is_int(v))
376 isl_die(isl_val_get_ctx(v), isl_error_invalid,
377 "expecting integer value", goto error);
378 mat = isl_mat_set_element(mat, row, col, v->n);
379 isl_val_free(v);
380 return mat;
381 error:
382 isl_val_free(v);
383 return isl_mat_free(mat);
384 }
385
386 __isl_give isl_mat *isl_mat_diag(isl_ctx *ctx, unsigned n_row, isl_int d)
387 {
388 int i;
389 struct isl_mat *mat;
390
391 mat = isl_mat_alloc(ctx, n_row, n_row);
392 if (!mat)
393 return NULL;
394 for (i = 0; i < n_row; ++i) {
395 isl_seq_clr(mat->row[i], i);
396 isl_int_set(mat->row[i][i], d);
397 isl_seq_clr(mat->row[i]+i+1, n_row-(i+1));
398 }
399
400 return mat;
401 }
402
403 /* Create an "n_row" by "n_col" matrix with zero elements.
404 */
405 __isl_give isl_mat *isl_mat_zero(isl_ctx *ctx, unsigned n_row, unsigned n_col)
406 {
407 int i;
408 isl_mat *mat;
409
410 mat = isl_mat_alloc(ctx, n_row, n_col);
411 if (!mat)
412 return NULL;
413 for (i = 0; i < n_row; ++i)
414 isl_seq_clr(mat->row[i], n_col);
415
416 return mat;
417 }
418
419 __isl_give isl_mat *isl_mat_identity(isl_ctx *ctx, unsigned n_row)
420 {
421 if (!ctx)
422 return NULL;
423 return isl_mat_diag(ctx, n_row, ctx->one);
424 }
425
426 /* Is "mat" a (possibly scaled) identity matrix?
427 */
428 isl_bool isl_mat_is_scaled_identity(__isl_keep isl_mat *mat)
429 {
430 int i;
431
432 if (!mat)
433 return isl_bool_error;
434 if (mat->n_row != mat->n_col)
435 return isl_bool_false;
436
437 for (i = 0; i < mat->n_row; ++i) {
438 if (isl_seq_first_non_zero(mat->row[i], i) != -1)
439 return isl_bool_false;
440 if (isl_int_ne(mat->row[0][0], mat->row[i][i]))
441 return isl_bool_false;
442 if (isl_seq_first_non_zero(mat->row[i] + i + 1,
443 mat->n_col - (i + 1)) != -1)
444 return isl_bool_false;
445 }
446
447 return isl_bool_true;
448 }
449
450 __isl_give isl_vec *isl_mat_vec_product(__isl_take isl_mat *mat,
451 __isl_take isl_vec *vec)
452 {
453 int i;
454 struct isl_vec *prod;
455
456 if (!mat || !vec)
457 goto error;
458
459 isl_assert(mat->ctx, mat->n_col == vec->size, goto error);
460
461 prod = isl_vec_alloc(mat->ctx, mat->n_row);
462 if (!prod)
463 goto error;
464
465 for (i = 0; i < prod->size; ++i)
466 isl_seq_inner_product(mat->row[i], vec->el, vec->size,
467 &prod->block.data[i]);
468 isl_mat_free(mat);
469 isl_vec_free(vec);
470 return prod;
471 error:
472 isl_mat_free(mat);
473 isl_vec_free(vec);
474 return NULL;
475 }
476
477 __isl_give isl_vec *isl_mat_vec_inverse_product(__isl_take isl_mat *mat,
478 __isl_take isl_vec *vec)
479 {
480 struct isl_mat *vec_mat;
481 int i;
482
483 if (!mat || !vec)
484 goto error;
485 vec_mat = isl_mat_alloc(vec->ctx, vec->size, 1);
486 if (!vec_mat)
487 goto error;
488 for (i = 0; i < vec->size; ++i)
489 isl_int_set(vec_mat->row[i][0], vec->el[i]);
490 vec_mat = isl_mat_inverse_product(mat, vec_mat);
491 isl_vec_free(vec);
492 if (!vec_mat)
493 return NULL;
494 vec = isl_vec_alloc(vec_mat->ctx, vec_mat->n_row);
495 if (vec)
496 for (i = 0; i < vec->size; ++i)
497 isl_int_set(vec->el[i], vec_mat->row[i][0]);
498 isl_mat_free(vec_mat);
499 return vec;
500 error:
501 isl_mat_free(mat);
502 isl_vec_free(vec);
503 return NULL;
504 }
505
506 __isl_give isl_vec *isl_vec_mat_product(__isl_take isl_vec *vec,
507 __isl_take isl_mat *mat)
508 {
509 int i, j;
510 struct isl_vec *prod;
511
512 if (!mat || !vec)
513 goto error;
514
515 isl_assert(mat->ctx, mat->n_row == vec->size, goto error);
516
517 prod = isl_vec_alloc(mat->ctx, mat->n_col);
518 if (!prod)
519 goto error;
520
521 for (i = 0; i < prod->size; ++i) {
522 isl_int_set_si(prod->el[i], 0);
523 for (j = 0; j < vec->size; ++j)
524 isl_int_addmul(prod->el[i], vec->el[j], mat->row[j][i]);
525 }
526 isl_mat_free(mat);
527 isl_vec_free(vec);
528 return prod;
529 error:
530 isl_mat_free(mat);
531 isl_vec_free(vec);
532 return NULL;
533 }
534
535 __isl_give isl_mat *isl_mat_aff_direct_sum(__isl_take isl_mat *left,
536 __isl_take isl_mat *right)
537 {
538 int i;
539 struct isl_mat *sum;
540
541 if (!left || !right)
542 goto error;
543
544 isl_assert(left->ctx, left->n_row == right->n_row, goto error);
545 isl_assert(left->ctx, left->n_row >= 1, goto error);
546 isl_assert(left->ctx, left->n_col >= 1, goto error);
547 isl_assert(left->ctx, right->n_col >= 1, goto error);
548 isl_assert(left->ctx,
549 isl_seq_first_non_zero(left->row[0]+1, left->n_col-1) == -1,
550 goto error);
551 isl_assert(left->ctx,
552 isl_seq_first_non_zero(right->row[0]+1, right->n_col-1) == -1,
553 goto error);
554
555 sum = isl_mat_alloc(left->ctx, left->n_row, left->n_col + right->n_col - 1);
556 if (!sum)
557 goto error;
558 isl_int_lcm(sum->row[0][0], left->row[0][0], right->row[0][0]);
559 isl_int_divexact(left->row[0][0], sum->row[0][0], left->row[0][0]);
560 isl_int_divexact(right->row[0][0], sum->row[0][0], right->row[0][0]);
561
562 isl_seq_clr(sum->row[0]+1, sum->n_col-1);
563 for (i = 1; i < sum->n_row; ++i) {
564 isl_int_mul(sum->row[i][0], left->row[0][0], left->row[i][0]);
565 isl_int_addmul(sum->row[i][0],
566 right->row[0][0], right->row[i][0]);
567 isl_seq_scale(sum->row[i]+1, left->row[i]+1, left->row[0][0],
568 left->n_col-1);
569 isl_seq_scale(sum->row[i]+left->n_col,
570 right->row[i]+1, right->row[0][0],
571 right->n_col-1);
572 }
573
574 isl_int_divexact(left->row[0][0], sum->row[0][0], left->row[0][0]);
575 isl_int_divexact(right->row[0][0], sum->row[0][0], right->row[0][0]);
576 isl_mat_free(left);
577 isl_mat_free(right);
578 return sum;
579 error:
580 isl_mat_free(left);
581 isl_mat_free(right);
582 return NULL;
583 }
584
585 static void exchange(__isl_keep isl_mat *M, __isl_keep isl_mat **U,
586 __isl_keep isl_mat **Q, unsigned row, unsigned i, unsigned j)
587 {
588 int r;
589 for (r = row; r < M->n_row; ++r)
590 isl_int_swap(M->row[r][i], M->row[r][j]);
591 if (U) {
592 for (r = 0; r < (*U)->n_row; ++r)
593 isl_int_swap((*U)->row[r][i], (*U)->row[r][j]);
594 }
595 if (Q)
596 isl_mat_swap_rows(*Q, i, j);
597 }
598
599 static void subtract(__isl_keep isl_mat *M, __isl_keep isl_mat **U,
600 __isl_keep isl_mat **Q, unsigned row, unsigned i, unsigned j, isl_int m)
601 {
602 int r;
603 for (r = row; r < M->n_row; ++r)
604 isl_int_submul(M->row[r][j], m, M->row[r][i]);
605 if (U) {
606 for (r = 0; r < (*U)->n_row; ++r)
607 isl_int_submul((*U)->row[r][j], m, (*U)->row[r][i]);
608 }
609 if (Q) {
610 for (r = 0; r < (*Q)->n_col; ++r)
611 isl_int_addmul((*Q)->row[i][r], m, (*Q)->row[j][r]);
612 }
613 }
614
615 static void oppose(__isl_keep isl_mat *M, __isl_keep isl_mat **U,
616 __isl_keep isl_mat **Q, unsigned row, unsigned col)
617 {
618 int r;
619 for (r = row; r < M->n_row; ++r)
620 isl_int_neg(M->row[r][col], M->row[r][col]);
621 if (U) {
622 for (r = 0; r < (*U)->n_row; ++r)
623 isl_int_neg((*U)->row[r][col], (*U)->row[r][col]);
624 }
625 if (Q)
626 isl_seq_neg((*Q)->row[col], (*Q)->row[col], (*Q)->n_col);
627 }
628
629 /* Given matrix M, compute
630 *
631 * M U = H
632 * M = H Q
633 *
634 * with U and Q unimodular matrices and H a matrix in column echelon form
635 * such that on each echelon row the entries in the non-echelon column
636 * are non-negative (if neg == 0) or non-positive (if neg == 1)
637 * and strictly smaller (in absolute value) than the entries in the echelon
638 * column.
639 * If U or Q are NULL, then these matrices are not computed.
640 */
641 __isl_give isl_mat *isl_mat_left_hermite(__isl_take isl_mat *M, int neg,
642 __isl_give isl_mat **U, __isl_give isl_mat **Q)
643 {
644 isl_int c;
645 int row, col;
646
647 if (U)
648 *U = NULL;
649 if (Q)
650 *Q = NULL;
651 if (!M)
652 goto error;
653 if (U) {
654 *U = isl_mat_identity(M->ctx, M->n_col);
655 if (!*U)
656 goto error;
657 }
658 if (Q) {
659 *Q = isl_mat_identity(M->ctx, M->n_col);
660 if (!*Q)
661 goto error;
662 }
663
664 if (M->n_col == 0)
665 return M;
666
667 M = isl_mat_cow(M);
668 if (!M)
669 goto error;
670
671 col = 0;
672 isl_int_init(c);
673 for (row = 0; row < M->n_row; ++row) {
674 int first, i, off;
675 first = isl_seq_abs_min_non_zero(M->row[row]+col, M->n_col-col);
676 if (first == -1)
677 continue;
678 first += col;
679 if (first != col)
680 exchange(M, U, Q, row, first, col);
681 if (isl_int_is_neg(M->row[row][col]))
682 oppose(M, U, Q, row, col);
683 first = col+1;
684 while ((off = isl_seq_first_non_zero(M->row[row]+first,
685 M->n_col-first)) != -1) {
686 first += off;
687 isl_int_fdiv_q(c, M->row[row][first], M->row[row][col]);
688 subtract(M, U, Q, row, col, first, c);
689 if (!isl_int_is_zero(M->row[row][first]))
690 exchange(M, U, Q, row, first, col);
691 else
692 ++first;
693 }
694 for (i = 0; i < col; ++i) {
695 if (isl_int_is_zero(M->row[row][i]))
696 continue;
697 if (neg)
698 isl_int_cdiv_q(c, M->row[row][i], M->row[row][col]);
699 else
700 isl_int_fdiv_q(c, M->row[row][i], M->row[row][col]);
701 if (isl_int_is_zero(c))
702 continue;
703 subtract(M, U, Q, row, col, i, c);
704 }
705 ++col;
706 }
707 isl_int_clear(c);
708
709 return M;
710 error:
711 if (Q) {
712 isl_mat_free(*Q);
713 *Q = NULL;
714 }
715 if (U) {
716 isl_mat_free(*U);
717 *U = NULL;
718 }
719 isl_mat_free(M);
720 return NULL;
721 }
722
723 /* Use row "row" of "mat" to eliminate column "col" from all other rows.
724 */
725 static __isl_give isl_mat *eliminate(__isl_take isl_mat *mat, int row, int col)
726 {
727 int k;
728 isl_size nr, nc;
729 isl_ctx *ctx;
730
731 nr = isl_mat_rows(mat);
732 nc = isl_mat_cols(mat);
733 if (nr < 0 || nc < 0)
734 return isl_mat_free(mat);
735
736 ctx = isl_mat_get_ctx(mat);
737
738 for (k = 0; k < nr; ++k) {
739 if (k == row)
740 continue;
741 if (isl_int_is_zero(mat->row[k][col]))
742 continue;
743 mat = isl_mat_cow(mat);
744 if (!mat)
745 return NULL;
746 isl_seq_elim(mat->row[k], mat->row[row], col, nc, NULL);
747 isl_seq_normalize(ctx, mat->row[k], nc);
748 }
749
750 return mat;
751 }
752
753 /* Perform Gaussian elimination on the rows of "mat", but start
754 * from the final row and the final column.
755 * Any zero rows that result from the elimination are removed.
756 *
757 * In particular, for each column from last to first,
758 * look for the last row with a non-zero coefficient in that column,
759 * move it last (but before other rows moved last in previous steps) and
760 * use it to eliminate the column from the other rows.
761 */
762 __isl_give isl_mat *isl_mat_reverse_gauss(__isl_take isl_mat *mat)
763 {
764 int k, row, last;
765 isl_size nr, nc;
766
767 nr = isl_mat_rows(mat);
768 nc = isl_mat_cols(mat);
769 if (nr < 0 || nc < 0)
770 return isl_mat_free(mat);
771
772 last = nc - 1;
773 for (row = nr - 1; row >= 0; --row) {
774 for (; last >= 0; --last) {
775 for (k = row; k >= 0; --k)
776 if (!isl_int_is_zero(mat->row[k][last]))
777 break;
778 if (k >= 0)
779 break;
780 }
781 if (last < 0)
782 break;
783 if (k != row)
784 mat = isl_mat_swap_rows(mat, k, row);
785 if (!mat)
786 return NULL;
787 if (isl_int_is_neg(mat->row[row][last]))
788 mat = isl_mat_row_neg(mat, row);
789 mat = eliminate(mat, row, last);
790 if (!mat)
791 return NULL;
792 }
793 mat = isl_mat_drop_rows(mat, 0, row + 1);
794
795 return mat;
796 }
797
798 /* Negate the lexicographically negative rows of "mat" such that
799 * all rows in the result are lexicographically non-negative.
800 */
801 __isl_give isl_mat *isl_mat_lexnonneg_rows(__isl_take isl_mat *mat)
802 {
803 int i;
804 isl_size nr, nc;
805
806 nr = isl_mat_rows(mat);
807 nc = isl_mat_cols(mat);
808 if (nr < 0 || nc < 0)
809 return isl_mat_free(mat);
810
811 for (i = 0; i < nr; ++i) {
812 int pos;
813
814 pos = isl_seq_first_non_zero(mat->row[i], nc);
815 if (pos < 0)
816 continue;
817 if (isl_int_is_nonneg(mat->row[i][pos]))
818 continue;
819 mat = isl_mat_row_neg(mat, i);
820 if (!mat)
821 return NULL;
822 }
823
824 return mat;
825 }
826
827 /* Given a matrix "H" is column echelon form, what is the first
828 * zero column? That is how many initial columns are non-zero?
829 * Start looking at column "first_col" and only consider
830 * the columns to be of size "n_row".
831 * "H" is assumed to be non-NULL.
832 *
833 * Since "H" is in column echelon form, the first non-zero entry
834 * in a column is always in a later position compared to the previous column.
835 */
836 static int hermite_first_zero_col(__isl_keep isl_mat *H, int first_col,
837 int n_row)
838 {
839 int row, col;
840
841 for (col = first_col, row = 0; col < H->n_col; ++col) {
842 for (; row < n_row; ++row)
843 if (!isl_int_is_zero(H->row[row][col]))
844 break;
845 if (row == n_row)
846 return col;
847 }
848
849 return H->n_col;
850 }
851
852 /* Return the rank of "mat", or isl_size_error in case of error.
853 */
854 isl_size isl_mat_rank(__isl_keep isl_mat *mat)
855 {
856 int rank;
857 isl_mat *H;
858
859 H = isl_mat_left_hermite(isl_mat_copy(mat), 0, NULL, NULL);
860 if (!H)
861 return isl_size_error;
862
863 rank = hermite_first_zero_col(H, 0, H->n_row);
864 isl_mat_free(H);
865
866 return rank;
867 }
868
869 __isl_give isl_mat *isl_mat_right_kernel(__isl_take isl_mat *mat)
870 {
871 int rank;
872 struct isl_mat *U = NULL;
873 struct isl_mat *K;
874
875 mat = isl_mat_left_hermite(mat, 0, &U, NULL);
876 if (!mat || !U)
877 goto error;
878
879 rank = hermite_first_zero_col(mat, 0, mat->n_row);
880 K = isl_mat_alloc(U->ctx, U->n_row, U->n_col - rank);
881 if (!K)
882 goto error;
883 isl_mat_sub_copy(K->ctx, K->row, U->row, U->n_row, 0, rank, U->n_col-rank);
884 isl_mat_free(mat);
885 isl_mat_free(U);
886 return K;
887 error:
888 isl_mat_free(mat);
889 isl_mat_free(U);
890 return NULL;
891 }
892
893 __isl_give isl_mat *isl_mat_lin_to_aff(__isl_take isl_mat *mat)
894 {
895 int i;
896 struct isl_mat *mat2;
897
898 if (!mat)
899 return NULL;
900 mat2 = isl_mat_alloc(mat->ctx, 1+mat->n_row, 1+mat->n_col);
901 if (!mat2)
902 goto error;
903 isl_int_set_si(mat2->row[0][0], 1);
904 isl_seq_clr(mat2->row[0]+1, mat->n_col);
905 for (i = 0; i < mat->n_row; ++i) {
906 isl_int_set_si(mat2->row[1+i][0], 0);
907 isl_seq_cpy(mat2->row[1+i]+1, mat->row[i], mat->n_col);
908 }
909 isl_mat_free(mat);
910 return mat2;
911 error:
912 isl_mat_free(mat);
913 return NULL;
914 }
915
916 /* Given two matrices M1 and M2, return the block matrix
917 *
918 * [ M1 0 ]
919 * [ 0 M2 ]
920 */
921 __isl_give isl_mat *isl_mat_diagonal(__isl_take isl_mat *mat1,
922 __isl_take isl_mat *mat2)
923 {
924 int i;
925 isl_mat *mat;
926
927 if (!mat1 || !mat2)
928 goto error;
929
930 mat = isl_mat_alloc(mat1->ctx, mat1->n_row + mat2->n_row,
931 mat1->n_col + mat2->n_col);
932 if (!mat)
933 goto error;
934 for (i = 0; i < mat1->n_row; ++i) {
935 isl_seq_cpy(mat->row[i], mat1->row[i], mat1->n_col);
936 isl_seq_clr(mat->row[i] + mat1->n_col, mat2->n_col);
937 }
938 for (i = 0; i < mat2->n_row; ++i) {
939 isl_seq_clr(mat->row[mat1->n_row + i], mat1->n_col);
940 isl_seq_cpy(mat->row[mat1->n_row + i] + mat1->n_col,
941 mat2->row[i], mat2->n_col);
942 }
943 isl_mat_free(mat1);
944 isl_mat_free(mat2);
945 return mat;
946 error:
947 isl_mat_free(mat1);
948 isl_mat_free(mat2);
949 return NULL;
950 }
951
952 static int row_first_non_zero(isl_int **row, unsigned n_row, unsigned col)
953 {
954 int i;
955
956 for (i = 0; i < n_row; ++i)
957 if (!isl_int_is_zero(row[i][col]))
958 return i;
959 return -1;
960 }
961
962 static int row_abs_min_non_zero(isl_int **row, unsigned n_row, unsigned col)
963 {
964 int i, min = row_first_non_zero(row, n_row, col);
965 if (min < 0)
966 return -1;
967 for (i = min + 1; i < n_row; ++i) {
968 if (isl_int_is_zero(row[i][col]))
969 continue;
970 if (isl_int_abs_lt(row[i][col], row[min][col]))
971 min = i;
972 }
973 return min;
974 }
975
976 static isl_stat inv_exchange(__isl_keep isl_mat **left,
977 __isl_keep isl_mat **right, unsigned i, unsigned j)
978 {
979 *left = isl_mat_swap_rows(*left, i, j);
980 *right = isl_mat_swap_rows(*right, i, j);
981
982 if (!*left || !*right)
983 return isl_stat_error;
984 return isl_stat_ok;
985 }
986
987 static void inv_oppose(
988 __isl_keep isl_mat *left, __isl_keep isl_mat *right, unsigned row)
989 {
990 isl_seq_neg(left->row[row]+row, left->row[row]+row, left->n_col-row);
991 isl_seq_neg(right->row[row], right->row[row], right->n_col);
992 }
993
994 static void inv_subtract(__isl_keep isl_mat *left, __isl_keep isl_mat *right,
995 unsigned row, unsigned i, isl_int m)
996 {
997 isl_int_neg(m, m);
998 isl_seq_combine(left->row[i]+row,
999 left->ctx->one, left->row[i]+row,
1000 m, left->row[row]+row,
1001 left->n_col-row);
1002 isl_seq_combine(right->row[i], right->ctx->one, right->row[i],
1003 m, right->row[row], right->n_col);
1004 }
1005
1006 /* Compute inv(left)*right
1007 */
1008 __isl_give isl_mat *isl_mat_inverse_product(__isl_take isl_mat *left,
1009 __isl_take isl_mat *right)
1010 {
1011 int row;
1012 isl_int a, b;
1013
1014 if (!left || !right)
1015 goto error;
1016
1017 isl_assert(left->ctx, left->n_row == left->n_col, goto error);
1018 isl_assert(left->ctx, left->n_row == right->n_row, goto error);
1019
1020 if (left->n_row == 0) {
1021 isl_mat_free(left);
1022 return right;
1023 }
1024
1025 left = isl_mat_cow(left);
1026 right = isl_mat_cow(right);
1027 if (!left || !right)
1028 goto error;
1029
1030 isl_int_init(a);
1031 isl_int_init(b);
1032 for (row = 0; row < left->n_row; ++row) {
1033 int pivot, first, i, off;
1034 pivot = row_abs_min_non_zero(left->row+row, left->n_row-row, row);
1035 if (pivot < 0) {
1036 isl_int_clear(a);
1037 isl_int_clear(b);
1038 isl_assert(left->ctx, pivot >= 0, goto error);
1039 }
1040 pivot += row;
1041 if (pivot != row)
1042 if (inv_exchange(&left, &right, pivot, row) < 0)
1043 goto error;
1044 if (isl_int_is_neg(left->row[row][row]))
1045 inv_oppose(left, right, row);
1046 first = row+1;
1047 while ((off = row_first_non_zero(left->row+first,
1048 left->n_row-first, row)) != -1) {
1049 first += off;
1050 isl_int_fdiv_q(a, left->row[first][row],
1051 left->row[row][row]);
1052 inv_subtract(left, right, row, first, a);
1053 if (!isl_int_is_zero(left->row[first][row])) {
1054 if (inv_exchange(&left, &right, row, first) < 0)
1055 goto error;
1056 } else {
1057 ++first;
1058 }
1059 }
1060 for (i = 0; i < row; ++i) {
1061 if (isl_int_is_zero(left->row[i][row]))
1062 continue;
1063 isl_int_gcd(a, left->row[row][row], left->row[i][row]);
1064 isl_int_divexact(b, left->row[i][row], a);
1065 isl_int_divexact(a, left->row[row][row], a);
1066 isl_int_neg(b, b);
1067 isl_seq_combine(left->row[i] + i,
1068 a, left->row[i] + i,
1069 b, left->row[row] + i,
1070 left->n_col - i);
1071 isl_seq_combine(right->row[i], a, right->row[i],
1072 b, right->row[row], right->n_col);
1073 }
1074 }
1075 isl_int_clear(b);
1076
1077 isl_int_set(a, left->row[0][0]);
1078 for (row = 1; row < left->n_row; ++row)
1079 isl_int_lcm(a, a, left->row[row][row]);
1080 if (isl_int_is_zero(a)){
1081 isl_int_clear(a);
1082 isl_assert(left->ctx, 0, goto error);
1083 }
1084 for (row = 0; row < left->n_row; ++row) {
1085 isl_int_divexact(left->row[row][row], a, left->row[row][row]);
1086 if (isl_int_is_one(left->row[row][row]))
1087 continue;
1088 isl_seq_scale(right->row[row], right->row[row],
1089 left->row[row][row], right->n_col);
1090 }
1091 isl_int_clear(a);
1092
1093 isl_mat_free(left);
1094 return right;
1095 error:
1096 isl_mat_free(left);
1097 isl_mat_free(right);
1098 return NULL;
1099 }
1100
1101 void isl_mat_col_scale(__isl_keep isl_mat *mat, unsigned col, isl_int m)
1102 {
1103 int i;
1104
1105 for (i = 0; i < mat->n_row; ++i)
1106 isl_int_mul(mat->row[i][col], mat->row[i][col], m);
1107 }
1108
1109 void isl_mat_col_combine(__isl_keep isl_mat *mat, unsigned dst,
1110 isl_int m1, unsigned src1, isl_int m2, unsigned src2)
1111 {
1112 int i;
1113 isl_int tmp;
1114
1115 isl_int_init(tmp);
1116 for (i = 0; i < mat->n_row; ++i) {
1117 isl_int_mul(tmp, m1, mat->row[i][src1]);
1118 isl_int_addmul(tmp, m2, mat->row[i][src2]);
1119 isl_int_set(mat->row[i][dst], tmp);
1120 }
1121 isl_int_clear(tmp);
1122 }
1123
1124 __isl_give isl_mat *isl_mat_right_inverse(__isl_take isl_mat *mat)
1125 {
1126 struct isl_mat *inv;
1127 int row;
1128 isl_int a, b;
1129
1130 mat = isl_mat_cow(mat);
1131 if (!mat)
1132 return NULL;
1133
1134 inv = isl_mat_identity(mat->ctx, mat->n_col);
1135 inv = isl_mat_cow(inv);
1136 if (!inv)
1137 goto error;
1138
1139 isl_int_init(a);
1140 isl_int_init(b);
1141 for (row = 0; row < mat->n_row; ++row) {
1142 int pivot, first, i, off;
1143 pivot = isl_seq_abs_min_non_zero(mat->row[row]+row, mat->n_col-row);
1144 if (pivot < 0) {
1145 isl_int_clear(a);
1146 isl_int_clear(b);
1147 isl_assert(mat->ctx, pivot >= 0, goto error);
1148 }
1149 pivot += row;
1150 if (pivot != row)
1151 exchange(mat, &inv, NULL, row, pivot, row);
1152 if (isl_int_is_neg(mat->row[row][row]))
1153 oppose(mat, &inv, NULL, row, row);
1154 first = row+1;
1155 while ((off = isl_seq_first_non_zero(mat->row[row]+first,
1156 mat->n_col-first)) != -1) {
1157 first += off;
1158 isl_int_fdiv_q(a, mat->row[row][first],
1159 mat->row[row][row]);
1160 subtract(mat, &inv, NULL, row, row, first, a);
1161 if (!isl_int_is_zero(mat->row[row][first]))
1162 exchange(mat, &inv, NULL, row, row, first);
1163 else
1164 ++first;
1165 }
1166 for (i = 0; i < row; ++i) {
1167 if (isl_int_is_zero(mat->row[row][i]))
1168 continue;
1169 isl_int_gcd(a, mat->row[row][row], mat->row[row][i]);
1170 isl_int_divexact(b, mat->row[row][i], a);
1171 isl_int_divexact(a, mat->row[row][row], a);
1172 isl_int_neg(a, a);
1173 isl_mat_col_combine(mat, i, a, i, b, row);
1174 isl_mat_col_combine(inv, i, a, i, b, row);
1175 }
1176 }
1177 isl_int_clear(b);
1178
1179 isl_int_set(a, mat->row[0][0]);
1180 for (row = 1; row < mat->n_row; ++row)
1181 isl_int_lcm(a, a, mat->row[row][row]);
1182 if (isl_int_is_zero(a)){
1183 isl_int_clear(a);
1184 goto error;
1185 }
1186 for (row = 0; row < mat->n_row; ++row) {
1187 isl_int_divexact(mat->row[row][row], a, mat->row[row][row]);
1188 if (isl_int_is_one(mat->row[row][row]))
1189 continue;
1190 isl_mat_col_scale(inv, row, mat->row[row][row]);
1191 }
1192 isl_int_clear(a);
1193
1194 isl_mat_free(mat);
1195
1196 return inv;
1197 error:
1198 isl_mat_free(mat);
1199 isl_mat_free(inv);
1200 return NULL;
1201 }
1202
1203 __isl_give isl_mat *isl_mat_transpose(__isl_take isl_mat *mat)
1204 {
1205 struct isl_mat *transpose = NULL;
1206 int i, j;
1207
1208 if (!mat)
1209 return NULL;
1210
1211 if (mat->n_col == mat->n_row) {
1212 mat = isl_mat_cow(mat);
1213 if (!mat)
1214 return NULL;
1215 for (i = 0; i < mat->n_row; ++i)
1216 for (j = i + 1; j < mat->n_col; ++j)
1217 isl_int_swap(mat->row[i][j], mat->row[j][i]);
1218 return mat;
1219 }
1220 transpose = isl_mat_alloc(mat->ctx, mat->n_col, mat->n_row);
1221 if (!transpose)
1222 goto error;
1223 for (i = 0; i < mat->n_row; ++i)
1224 for (j = 0; j < mat->n_col; ++j)
1225 isl_int_set(transpose->row[j][i], mat->row[i][j]);
1226 isl_mat_free(mat);
1227 return transpose;
1228 error:
1229 isl_mat_free(mat);
1230 return NULL;
1231 }
1232
1233 __isl_give isl_mat *isl_mat_swap_cols(__isl_take isl_mat *mat,
1234 unsigned i, unsigned j)
1235 {
1236 int r;
1237
1238 mat = isl_mat_cow(mat);
1239 if (check_col_range(mat, i, 1) < 0 ||
1240 check_col_range(mat, j, 1) < 0)
1241 return isl_mat_free(mat);
1242
1243 for (r = 0; r < mat->n_row; ++r)
1244 isl_int_swap(mat->row[r][i], mat->row[r][j]);
1245 return mat;
1246 }
1247
1248 __isl_give isl_mat *isl_mat_swap_rows(__isl_take isl_mat *mat,
1249 unsigned i, unsigned j)
1250 {
1251 isl_int *t;
1252
1253 if (!mat)
1254 return NULL;
1255 mat = isl_mat_cow(mat);
1256 if (check_row_range(mat, i, 1) < 0 ||
1257 check_row_range(mat, j, 1) < 0)
1258 return isl_mat_free(mat);
1259
1260 t = mat->row[i];
1261 mat->row[i] = mat->row[j];
1262 mat->row[j] = t;
1263 return mat;
1264 }
1265
1266 /* Calculate the product of two matrices.
1267 *
1268 * This function is optimized for operand matrices that contain many zeros and
1269 * skips multiplications where we know one of the operands is zero.
1270 */
1271 __isl_give isl_mat *isl_mat_product(__isl_take isl_mat *left,
1272 __isl_take isl_mat *right)
1273 {
1274 int i, j, k;
1275 struct isl_mat *prod;
1276
1277 if (!left || !right)
1278 goto error;
1279 isl_assert(left->ctx, left->n_col == right->n_row, goto error);
1280 prod = isl_mat_alloc(left->ctx, left->n_row, right->n_col);
1281 if (!prod)
1282 goto error;
1283 if (left->n_col == 0) {
1284 for (i = 0; i < prod->n_row; ++i)
1285 isl_seq_clr(prod->row[i], prod->n_col);
1286 isl_mat_free(left);
1287 isl_mat_free(right);
1288 return prod;
1289 }
1290 for (i = 0; i < prod->n_row; ++i) {
1291 for (j = 0; j < prod->n_col; ++j)
1292 isl_int_mul(prod->row[i][j],
1293 left->row[i][0], right->row[0][j]);
1294 for (k = 1; k < left->n_col; ++k) {
1295 if (isl_int_is_zero(left->row[i][k]))
1296 continue;
1297 for (j = 0; j < prod->n_col; ++j)
1298 isl_int_addmul(prod->row[i][j],
1299 left->row[i][k], right->row[k][j]);
1300 }
1301 }
1302 isl_mat_free(left);
1303 isl_mat_free(right);
1304 return prod;
1305 error:
1306 isl_mat_free(left);
1307 isl_mat_free(right);
1308 return NULL;
1309 }
1310
1311 /* Replace the variables x in the rows q by x' given by x = M x',
1312 * with M the matrix mat.
1313 *
1314 * If the number of new variables is greater than the original
1315 * number of variables, then the rows q have already been
1316 * preextended. If the new number is smaller, then the coefficients
1317 * of the divs, which are not changed, need to be shifted down.
1318 * The row q may be the equalities, the inequalities or the
1319 * div expressions. In the latter case, has_div is true and
1320 * we need to take into account the extra denominator column.
1321 */
1322 static int preimage(struct isl_ctx *ctx, isl_int **q, unsigned n,
1323 unsigned n_div, int has_div, struct isl_mat *mat)
1324 {
1325 int i;
1326 struct isl_mat *t;
1327 int e;
1328
1329 if (mat->n_col >= mat->n_row)
1330 e = 0;
1331 else
1332 e = mat->n_row - mat->n_col;
1333 if (has_div)
1334 for (i = 0; i < n; ++i)
1335 isl_int_mul(q[i][0], q[i][0], mat->row[0][0]);
1336 t = isl_mat_sub_alloc6(mat->ctx, q, 0, n, has_div, mat->n_row);
1337 t = isl_mat_product(t, mat);
1338 if (!t)
1339 return -1;
1340 for (i = 0; i < n; ++i) {
1341 isl_seq_swp_or_cpy(q[i] + has_div, t->row[i], t->n_col);
1342 isl_seq_cpy(q[i] + has_div + t->n_col,
1343 q[i] + has_div + t->n_col + e, n_div);
1344 isl_seq_clr(q[i] + has_div + t->n_col + n_div, e);
1345 }
1346 isl_mat_free(t);
1347 return 0;
1348 }
1349
1350 /* Replace the variables x in bset by x' given by x = M x', with
1351 * M the matrix mat.
1352 *
1353 * If there are fewer variables x' then there are x, then we perform
1354 * the transformation in place, which means that, in principle,
1355 * this frees up some extra variables as the number
1356 * of columns remains constant, but we would have to extend
1357 * the div array too as the number of rows in this array is assumed
1358 * to be equal to extra.
1359 */
1360 __isl_give isl_basic_set *isl_basic_set_preimage(
1361 __isl_take isl_basic_set *bset, __isl_take isl_mat *mat)
1362 {
1363 struct isl_ctx *ctx;
1364
1365 if (!bset || !mat)
1366 goto error;
1367
1368 ctx = bset->ctx;
1369 bset = isl_basic_set_cow(bset);
1370 if (isl_basic_set_check_no_params(bset) < 0)
1371 goto error;
1372
1373 isl_assert(ctx, 1+bset->dim->n_out == mat->n_row, goto error);
1374 isl_assert(ctx, mat->n_col > 0, goto error);
1375
1376 if (mat->n_col > mat->n_row) {
1377 bset = isl_basic_set_add_dims(bset, isl_dim_set,
1378 mat->n_col - mat->n_row);
1379 if (!bset)
1380 goto error;
1381 } else if (mat->n_col < mat->n_row) {
1382 bset->dim = isl_space_cow(bset->dim);
1383 if (!bset->dim)
1384 goto error;
1385 bset->dim->n_out -= mat->n_row - mat->n_col;
1386 }
1387
1388 if (preimage(ctx, bset->eq, bset->n_eq, bset->n_div, 0,
1389 isl_mat_copy(mat)) < 0)
1390 goto error;
1391
1392 if (preimage(ctx, bset->ineq, bset->n_ineq, bset->n_div, 0,
1393 isl_mat_copy(mat)) < 0)
1394 goto error;
1395
1396 if (preimage(ctx, bset->div, bset->n_div, bset->n_div, 1, mat) < 0)
1397 goto error2;
1398
1399 ISL_F_CLR(bset, ISL_BASIC_SET_NO_IMPLICIT);
1400 ISL_F_CLR(bset, ISL_BASIC_SET_NO_REDUNDANT);
1401 ISL_F_CLR(bset, ISL_BASIC_SET_SORTED);
1402 ISL_F_CLR(bset, ISL_BASIC_SET_NORMALIZED_DIVS);
1403 ISL_F_CLR(bset, ISL_BASIC_SET_ALL_EQUALITIES);
1404 ISL_F_CLR(bset, ISL_BASIC_SET_REDUCED_COEFFICIENTS);
1405
1406 bset = isl_basic_set_simplify(bset);
1407 bset = isl_basic_set_finalize(bset);
1408
1409 return bset;
1410 error:
1411 isl_mat_free(mat);
1412 error2:
1413 isl_basic_set_free(bset);
1414 return NULL;
1415 }
1416
1417 __isl_give isl_set *isl_set_preimage(
1418 __isl_take isl_set *set, __isl_take isl_mat *mat)
1419 {
1420 int i;
1421
1422 set = isl_set_cow(set);
1423 if (!set)
1424 goto error;
1425
1426 for (i = 0; i < set->n; ++i) {
1427 set->p[i] = isl_basic_set_preimage(set->p[i],
1428 isl_mat_copy(mat));
1429 if (!set->p[i])
1430 goto error;
1431 }
1432 if (mat->n_col != mat->n_row) {
1433 set->dim = isl_space_cow(set->dim);
1434 if (!set->dim)
1435 goto error;
1436 set->dim->n_out += mat->n_col;
1437 set->dim->n_out -= mat->n_row;
1438 }
1439 isl_mat_free(mat);
1440 ISL_F_CLR(set, ISL_SET_NORMALIZED);
1441 return set;
1442 error:
1443 isl_set_free(set);
1444 isl_mat_free(mat);
1445 return NULL;
1446 }
1447
1448 /* Replace the variables x starting at "first_col" in the rows "rows"
1449 * of some coefficient matrix by x' with x = M x' with M the matrix mat.
1450 * That is, replace the corresponding coefficients c by c M.
1451 */
1452 isl_stat isl_mat_sub_transform(isl_int **row, unsigned n_row,
1453 unsigned first_col, __isl_take isl_mat *mat)
1454 {
1455 int i;
1456 isl_ctx *ctx;
1457 isl_mat *t;
1458
1459 if (!mat)
1460 return isl_stat_error;
1461 ctx = isl_mat_get_ctx(mat);
1462 t = isl_mat_sub_alloc6(ctx, row, 0, n_row, first_col, mat->n_row);
1463 t = isl_mat_product(t, mat);
1464 if (!t)
1465 return isl_stat_error;
1466 for (i = 0; i < n_row; ++i)
1467 isl_seq_swp_or_cpy(row[i] + first_col, t->row[i], t->n_col);
1468 isl_mat_free(t);
1469 return isl_stat_ok;
1470 }
1471
1472 void isl_mat_print_internal(__isl_keep isl_mat *mat, FILE *out, int indent)
1473 {
1474 int i, j;
1475
1476 if (!mat) {
1477 fprintf(out, "%*snull mat\n", indent, "");
1478 return;
1479 }
1480
1481 if (mat->n_row == 0)
1482 fprintf(out, "%*s[]\n", indent, "");
1483
1484 for (i = 0; i < mat->n_row; ++i) {
1485 if (!i)
1486 fprintf(out, "%*s[[", indent, "");
1487 else
1488 fprintf(out, "%*s[", indent+1, "");
1489 for (j = 0; j < mat->n_col; ++j) {
1490 if (j)
1491 fprintf(out, ",");
1492 isl_int_print(out, mat->row[i][j], 0);
1493 }
1494 if (i == mat->n_row-1)
1495 fprintf(out, "]]\n");
1496 else
1497 fprintf(out, "]\n");
1498 }
1499 }
1500
1501 void isl_mat_dump(__isl_keep isl_mat *mat)
1502 {
1503 isl_mat_print_internal(mat, stderr, 0);
1504 }
1505
1506 __isl_give isl_mat *isl_mat_drop_cols(__isl_take isl_mat *mat,
1507 unsigned col, unsigned n)
1508 {
1509 int r;
1510
1511 if (n == 0)
1512 return mat;
1513
1514 mat = isl_mat_cow(mat);
1515 if (check_col_range(mat, col, n) < 0)
1516 return isl_mat_free(mat);
1517
1518 if (col != mat->n_col-n) {
1519 for (r = 0; r < mat->n_row; ++r)
1520 isl_seq_cpy(mat->row[r]+col, mat->row[r]+col+n,
1521 mat->n_col - col - n);
1522 }
1523 mat->n_col -= n;
1524 return mat;
1525 }
1526
1527 __isl_give isl_mat *isl_mat_drop_rows(__isl_take isl_mat *mat,
1528 unsigned row, unsigned n)
1529 {
1530 int r;
1531
1532 mat = isl_mat_cow(mat);
1533 if (check_row_range(mat, row, n) < 0)
1534 return isl_mat_free(mat);
1535
1536 for (r = row; r+n < mat->n_row; ++r)
1537 mat->row[r] = mat->row[r+n];
1538
1539 mat->n_row -= n;
1540 return mat;
1541 }
1542
1543 __isl_give isl_mat *isl_mat_insert_cols(__isl_take isl_mat *mat,
1544 unsigned col, unsigned n)
1545 {
1546 isl_mat *ext;
1547
1548 if (check_col_range(mat, col, 0) < 0)
1549 return isl_mat_free(mat);
1550 if (n == 0)
1551 return mat;
1552
1553 ext = isl_mat_alloc(mat->ctx, mat->n_row, mat->n_col + n);
1554 if (!ext)
1555 goto error;
1556
1557 isl_mat_sub_copy(mat->ctx, ext->row, mat->row, mat->n_row, 0, 0, col);
1558 isl_mat_sub_copy(mat->ctx, ext->row, mat->row, mat->n_row,
1559 col + n, col, mat->n_col - col);
1560
1561 isl_mat_free(mat);
1562 return ext;
1563 error:
1564 isl_mat_free(mat);
1565 return NULL;
1566 }
1567
1568 __isl_give isl_mat *isl_mat_insert_zero_cols(__isl_take isl_mat *mat,
1569 unsigned first, unsigned n)
1570 {
1571 int i;
1572
1573 if (!mat)
1574 return NULL;
1575 mat = isl_mat_insert_cols(mat, first, n);
1576 if (!mat)
1577 return NULL;
1578
1579 for (i = 0; i < mat->n_row; ++i)
1580 isl_seq_clr(mat->row[i] + first, n);
1581
1582 return mat;
1583 }
1584
1585 __isl_give isl_mat *isl_mat_add_zero_cols(__isl_take isl_mat *mat, unsigned n)
1586 {
1587 if (!mat)
1588 return NULL;
1589
1590 return isl_mat_insert_zero_cols(mat, mat->n_col, n);
1591 }
1592
1593 __isl_give isl_mat *isl_mat_insert_rows(__isl_take isl_mat *mat,
1594 unsigned row, unsigned n)
1595 {
1596 isl_mat *ext;
1597
1598 if (check_row_range(mat, row, 0) < 0)
1599 return isl_mat_free(mat);
1600 if (n == 0)
1601 return mat;
1602
1603 ext = isl_mat_alloc(mat->ctx, mat->n_row + n, mat->n_col);
1604 if (!ext)
1605 goto error;
1606
1607 isl_mat_sub_copy(mat->ctx, ext->row, mat->row, row, 0, 0, mat->n_col);
1608 isl_mat_sub_copy(mat->ctx, ext->row + row + n, mat->row + row,
1609 mat->n_row - row, 0, 0, mat->n_col);
1610
1611 isl_mat_free(mat);
1612 return ext;
1613 error:
1614 isl_mat_free(mat);
1615 return NULL;
1616 }
1617
1618 __isl_give isl_mat *isl_mat_add_rows(__isl_take isl_mat *mat, unsigned n)
1619 {
1620 if (!mat)
1621 return NULL;
1622
1623 return isl_mat_insert_rows(mat, mat->n_row, n);
1624 }
1625
1626 __isl_give isl_mat *isl_mat_insert_zero_rows(__isl_take isl_mat *mat,
1627 unsigned row, unsigned n)
1628 {
1629 int i;
1630
1631 mat = isl_mat_insert_rows(mat, row, n);
1632 if (!mat)
1633 return NULL;
1634
1635 for (i = 0; i < n; ++i)
1636 isl_seq_clr(mat->row[row + i], mat->n_col);
1637
1638 return mat;
1639 }
1640
1641 __isl_give isl_mat *isl_mat_add_zero_rows(__isl_take isl_mat *mat, unsigned n)
1642 {
1643 if (!mat)
1644 return NULL;
1645
1646 return isl_mat_insert_zero_rows(mat, mat->n_row, n);
1647 }
1648
1649 void isl_mat_col_submul(__isl_keep isl_mat *mat,
1650 int dst_col, isl_int f, int src_col)
1651 {
1652 int i;
1653
1654 for (i = 0; i < mat->n_row; ++i)
1655 isl_int_submul(mat->row[i][dst_col], f, mat->row[i][src_col]);
1656 }
1657
1658 void isl_mat_col_add(__isl_keep isl_mat *mat, int dst_col, int src_col)
1659 {
1660 int i;
1661
1662 if (!mat)
1663 return;
1664
1665 for (i = 0; i < mat->n_row; ++i)
1666 isl_int_add(mat->row[i][dst_col],
1667 mat->row[i][dst_col], mat->row[i][src_col]);
1668 }
1669
1670 void isl_mat_col_mul(__isl_keep isl_mat *mat, int dst_col, isl_int f,
1671 int src_col)
1672 {
1673 int i;
1674
1675 for (i = 0; i < mat->n_row; ++i)
1676 isl_int_mul(mat->row[i][dst_col], f, mat->row[i][src_col]);
1677 }
1678
1679 /* Add "f" times column "src_col" to column "dst_col" of "mat" and
1680 * return the result.
1681 */
1682 __isl_give isl_mat *isl_mat_col_addmul(__isl_take isl_mat *mat, int dst_col,
1683 isl_int f, int src_col)
1684 {
1685 int i;
1686
1687 if (check_col(mat, dst_col) < 0 || check_col(mat, src_col) < 0)
1688 return isl_mat_free(mat);
1689
1690 for (i = 0; i < mat->n_row; ++i) {
1691 if (isl_int_is_zero(mat->row[i][src_col]))
1692 continue;
1693 mat = isl_mat_cow(mat);
1694 if (!mat)
1695 return NULL;
1696 isl_int_addmul(mat->row[i][dst_col], f, mat->row[i][src_col]);
1697 }
1698
1699 return mat;
1700 }
1701
1702 /* Negate column "col" of "mat" and return the result.
1703 */
1704 __isl_give isl_mat *isl_mat_col_neg(__isl_take isl_mat *mat, int col)
1705 {
1706 int i;
1707
1708 if (check_col(mat, col) < 0)
1709 return isl_mat_free(mat);
1710
1711 for (i = 0; i < mat->n_row; ++i) {
1712 if (isl_int_is_zero(mat->row[i][col]))
1713 continue;
1714 mat = isl_mat_cow(mat);
1715 if (!mat)
1716 return NULL;
1717 isl_int_neg(mat->row[i][col], mat->row[i][col]);
1718 }
1719
1720 return mat;
1721 }
1722
1723 /* Negate row "row" of "mat" and return the result.
1724 */
1725 __isl_give isl_mat *isl_mat_row_neg(__isl_take isl_mat *mat, int row)
1726 {
1727 if (check_row(mat, row) < 0)
1728 return isl_mat_free(mat);
1729 if (isl_seq_first_non_zero(mat->row[row], mat->n_col) == -1)
1730 return mat;
1731 mat = isl_mat_cow(mat);
1732 if (!mat)
1733 return NULL;
1734 isl_seq_neg(mat->row[row], mat->row[row], mat->n_col);
1735 return mat;
1736 }
1737
1738 __isl_give isl_mat *isl_mat_unimodular_complete(__isl_take isl_mat *M, int row)
1739 {
1740 int r;
1741 struct isl_mat *H = NULL, *Q = NULL;
1742
1743 if (!M)
1744 return NULL;
1745
1746 isl_assert(M->ctx, M->n_row == M->n_col, goto error);
1747 M->n_row = row;
1748 H = isl_mat_left_hermite(isl_mat_copy(M), 0, NULL, &Q);
1749 M->n_row = M->n_col;
1750 if (!H)
1751 goto error;
1752 for (r = 0; r < row; ++r)
1753 isl_assert(M->ctx, isl_int_is_one(H->row[r][r]), goto error);
1754 for (r = row; r < M->n_row; ++r)
1755 isl_seq_cpy(M->row[r], Q->row[r], M->n_col);
1756 isl_mat_free(H);
1757 isl_mat_free(Q);
1758 return M;
1759 error:
1760 isl_mat_free(H);
1761 isl_mat_free(Q);
1762 isl_mat_free(M);
1763 return NULL;
1764 }
1765
1766 __isl_give isl_mat *isl_mat_concat(__isl_take isl_mat *top,
1767 __isl_take isl_mat *bot)
1768 {
1769 struct isl_mat *mat;
1770
1771 if (!top || !bot)
1772 goto error;
1773
1774 isl_assert(top->ctx, top->n_col == bot->n_col, goto error);
1775 if (top->n_row == 0) {
1776 isl_mat_free(top);
1777 return bot;
1778 }
1779 if (bot->n_row == 0) {
1780 isl_mat_free(bot);
1781 return top;
1782 }
1783
1784 mat = isl_mat_alloc(top->ctx, top->n_row + bot->n_row, top->n_col);
1785 if (!mat)
1786 goto error;
1787 isl_mat_sub_copy(mat->ctx, mat->row, top->row, top->n_row,
1788 0, 0, mat->n_col);
1789 isl_mat_sub_copy(mat->ctx, mat->row + top->n_row, bot->row, bot->n_row,
1790 0, 0, mat->n_col);
1791 isl_mat_free(top);
1792 isl_mat_free(bot);
1793 return mat;
1794 error:
1795 isl_mat_free(top);
1796 isl_mat_free(bot);
1797 return NULL;
1798 }
1799
1800 isl_bool isl_mat_is_equal(__isl_keep isl_mat *mat1, __isl_keep isl_mat *mat2)
1801 {
1802 int i;
1803
1804 if (!mat1 || !mat2)
1805 return isl_bool_error;
1806
1807 if (mat1->n_row != mat2->n_row)
1808 return isl_bool_false;
1809
1810 if (mat1->n_col != mat2->n_col)
1811 return isl_bool_false;
1812
1813 for (i = 0; i < mat1->n_row; ++i)
1814 if (!isl_seq_eq(mat1->row[i], mat2->row[i], mat1->n_col))
1815 return isl_bool_false;
1816
1817 return isl_bool_true;
1818 }
1819
1820 __isl_give isl_mat *isl_mat_from_row_vec(__isl_take isl_vec *vec)
1821 {
1822 struct isl_mat *mat;
1823
1824 if (!vec)
1825 return NULL;
1826 mat = isl_mat_alloc(vec->ctx, 1, vec->size);
1827 if (!mat)
1828 goto error;
1829
1830 isl_seq_cpy(mat->row[0], vec->el, vec->size);
1831
1832 isl_vec_free(vec);
1833 return mat;
1834 error:
1835 isl_vec_free(vec);
1836 return NULL;
1837 }
1838
1839 /* Return a copy of row "row" of "mat" as an isl_vec.
1840 */
1841 __isl_give isl_vec *isl_mat_get_row(__isl_keep isl_mat *mat, unsigned row)
1842 {
1843 isl_vec *v;
1844
1845 if (!mat)
1846 return NULL;
1847 if (row >= mat->n_row)
1848 isl_die(mat->ctx, isl_error_invalid, "row out of range",
1849 return NULL);
1850
1851 v = isl_vec_alloc(isl_mat_get_ctx(mat), mat->n_col);
1852 if (!v)
1853 return NULL;
1854 isl_seq_cpy(v->el, mat->row[row], mat->n_col);
1855
1856 return v;
1857 }
1858
1859 __isl_give isl_mat *isl_mat_vec_concat(__isl_take isl_mat *top,
1860 __isl_take isl_vec *bot)
1861 {
1862 return isl_mat_concat(top, isl_mat_from_row_vec(bot));
1863 }
1864
1865 __isl_give isl_mat *isl_mat_move_cols(__isl_take isl_mat *mat,
1866 unsigned dst_col, unsigned src_col, unsigned n)
1867 {
1868 isl_mat *res;
1869
1870 if (!mat)
1871 return NULL;
1872 if (n == 0 || dst_col == src_col)
1873 return mat;
1874
1875 res = isl_mat_alloc(mat->ctx, mat->n_row, mat->n_col);
1876 if (!res)
1877 goto error;
1878
1879 if (dst_col < src_col) {
1880 isl_mat_sub_copy(res->ctx, res->row, mat->row, mat->n_row,
1881 0, 0, dst_col);
1882 isl_mat_sub_copy(res->ctx, res->row, mat->row, mat->n_row,
1883 dst_col, src_col, n);
1884 isl_mat_sub_copy(res->ctx, res->row, mat->row, mat->n_row,
1885 dst_col + n, dst_col, src_col - dst_col);
1886 isl_mat_sub_copy(res->ctx, res->row, mat->row, mat->n_row,
1887 src_col + n, src_col + n,
1888 res->n_col - src_col - n);
1889 } else {
1890 isl_mat_sub_copy(res->ctx, res->row, mat->row, mat->n_row,
1891 0, 0, src_col);
1892 isl_mat_sub_copy(res->ctx, res->row, mat->row, mat->n_row,
1893 src_col, src_col + n, dst_col - src_col);
1894 isl_mat_sub_copy(res->ctx, res->row, mat->row, mat->n_row,
1895 dst_col, src_col, n);
1896 isl_mat_sub_copy(res->ctx, res->row, mat->row, mat->n_row,
1897 dst_col + n, dst_col + n,
1898 res->n_col - dst_col - n);
1899 }
1900 isl_mat_free(mat);
1901
1902 return res;
1903 error:
1904 isl_mat_free(mat);
1905 return NULL;
1906 }
1907
1908 /* Return the gcd of the elements in row "row" of "mat" in *gcd.
1909 * Return isl_stat_ok on success and isl_stat_error on failure.
1910 */
1911 isl_stat isl_mat_row_gcd(__isl_keep isl_mat *mat, int row, isl_int *gcd)
1912 {
1913 if (check_row(mat, row) < 0)
1914 return isl_stat_error;
1915
1916 isl_seq_gcd(mat->row[row], mat->n_col, gcd);
1917
1918 return isl_stat_ok;
1919 }
1920
1921 void isl_mat_gcd(__isl_keep isl_mat *mat, isl_int *gcd)
1922 {
1923 int i;
1924 isl_int g;
1925
1926 isl_int_set_si(*gcd, 0);
1927 if (!mat)
1928 return;
1929
1930 isl_int_init(g);
1931 for (i = 0; i < mat->n_row; ++i) {
1932 isl_seq_gcd(mat->row[i], mat->n_col, &g);
1933 isl_int_gcd(*gcd, *gcd, g);
1934 }
1935 isl_int_clear(g);
1936 }
1937
1938 /* Return the result of scaling "mat" by a factor of "m".
1939 */
1940 __isl_give isl_mat *isl_mat_scale(__isl_take isl_mat *mat, isl_int m)
1941 {
1942 int i;
1943
1944 if (isl_int_is_one(m))
1945 return mat;
1946
1947 mat = isl_mat_cow(mat);
1948 if (!mat)
1949 return NULL;
1950
1951 for (i = 0; i < mat->n_row; ++i)
1952 isl_seq_scale(mat->row[i], mat->row[i], m, mat->n_col);
1953
1954 return mat;
1955 }
1956
1957 __isl_give isl_mat *isl_mat_scale_down(__isl_take isl_mat *mat, isl_int m)
1958 {
1959 int i;
1960
1961 if (isl_int_is_one(m))
1962 return mat;
1963
1964 mat = isl_mat_cow(mat);
1965 if (!mat)
1966 return NULL;
1967
1968 for (i = 0; i < mat->n_row; ++i)
1969 isl_seq_scale_down(mat->row[i], mat->row[i], m, mat->n_col);
1970
1971 return mat;
1972 }
1973
1974 __isl_give isl_mat *isl_mat_scale_down_row(__isl_take isl_mat *mat, int row,
1975 isl_int m)
1976 {
1977 if (isl_int_is_one(m))
1978 return mat;
1979
1980 mat = isl_mat_cow(mat);
1981 if (!mat)
1982 return NULL;
1983
1984 isl_seq_scale_down(mat->row[row], mat->row[row], m, mat->n_col);
1985
1986 return mat;
1987 }
1988
1989 __isl_give isl_mat *isl_mat_normalize(__isl_take isl_mat *mat)
1990 {
1991 isl_int gcd;
1992
1993 if (!mat)
1994 return NULL;
1995
1996 isl_int_init(gcd);
1997 isl_mat_gcd(mat, &gcd);
1998 mat = isl_mat_scale_down(mat, gcd);
1999 isl_int_clear(gcd);
2000
2001 return mat;
2002 }
2003
2004 __isl_give isl_mat *isl_mat_normalize_row(__isl_take isl_mat *mat, int row)
2005 {
2006 mat = isl_mat_cow(mat);
2007 if (!mat)
2008 return NULL;
2009
2010 isl_seq_normalize(mat->ctx, mat->row[row], mat->n_col);
2011
2012 return mat;
2013 }
2014
2015 /* Number of initial non-zero columns.
2016 */
2017 int isl_mat_initial_non_zero_cols(__isl_keep isl_mat *mat)
2018 {
2019 int i;
2020
2021 if (!mat)
2022 return -1;
2023
2024 for (i = 0; i < mat->n_col; ++i)
2025 if (row_first_non_zero(mat->row, mat->n_row, i) < 0)
2026 break;
2027
2028 return i;
2029 }
2030
2031 /* Return a basis for the space spanned by the rows of "mat".
2032 * Any basis will do, so simply perform Gaussian elimination and
2033 * remove the empty rows.
2034 */
2035 __isl_give isl_mat *isl_mat_row_basis(__isl_take isl_mat *mat)
2036 {
2037 return isl_mat_reverse_gauss(mat);
2038 }
2039
2040 /* Return rows that extend a basis of "mat1" to one
2041 * that covers both "mat1" and "mat2".
2042 * The Hermite normal form of the concatenation of the two matrices is
2043 *
2044 * [ Q1 ]
2045 * [ M1 ] = [ H1 0 0 ] [ Q2 ]
2046 * [ M2 ] = [ H2 H3 0 ] [ Q3 ]
2047 *
2048 * The number of columns in H1 and H3 determine the number of rows
2049 * in Q1 and Q2. Q1 is a basis for M1, while Q2 extends this basis
2050 * to also cover M2.
2051 */
2052 __isl_give isl_mat *isl_mat_row_basis_extension(
2053 __isl_take isl_mat *mat1, __isl_take isl_mat *mat2)
2054 {
2055 isl_size n_row;
2056 int r1, r;
2057 isl_size n1;
2058 isl_mat *H, *Q;
2059
2060 n1 = isl_mat_rows(mat1);
2061 H = isl_mat_concat(mat1, mat2);
2062 H = isl_mat_left_hermite(H, 0, NULL, &Q);
2063 if (n1 < 0 || !H || !Q)
2064 goto error;
2065
2066 r1 = hermite_first_zero_col(H, 0, n1);
2067 r = hermite_first_zero_col(H, r1, H->n_row);
2068 n_row = isl_mat_rows(Q);
2069 if (n_row < 0)
2070 goto error;
2071 Q = isl_mat_drop_rows(Q, r, n_row - r);
2072 Q = isl_mat_drop_rows(Q, 0, r1);
2073
2074 isl_mat_free(H);
2075 return Q;
2076 error:
2077 isl_mat_free(H);
2078 isl_mat_free(Q);
2079 return NULL;
2080 }
2081
2082 /* Are the rows of "mat1" linearly independent of those of "mat2"?
2083 * That is, is there no linear dependence among the combined rows
2084 * that is not already present in either "mat1" or "mat2"?
2085 * In other words, is the rank of "mat1" and "mat2" combined equal
2086 * to the sum of the ranks of "mat1" and "mat2"?
2087 */
2088 isl_bool isl_mat_has_linearly_independent_rows(__isl_keep isl_mat *mat1,
2089 __isl_keep isl_mat *mat2)
2090 {
2091 isl_size r1, r2, r;
2092 isl_mat *mat;
2093
2094 r1 = isl_mat_rank(mat1);
2095 if (r1 < 0)
2096 return isl_bool_error;
2097 if (r1 == 0)
2098 return isl_bool_true;
2099 r2 = isl_mat_rank(mat2);
2100 if (r2 < 0)
2101 return isl_bool_error;
2102 if (r2 == 0)
2103 return isl_bool_true;
2104
2105 mat = isl_mat_concat(isl_mat_copy(mat1), isl_mat_copy(mat2));
2106 r = isl_mat_rank(mat);
2107 isl_mat_free(mat);
2108 if (r < 0)
2109 return isl_bool_error;
2110 return isl_bool_ok(r == r1 + r2);
2111 }
Integer set library