?  search:
summary | log | graphiclog1 | graphiclog2 | commit | commitdiff | tree | refs | edit | fork
blame | history | raw | HEAD
add isl_qpolynomial_domain_reverse
[isl.git] / isl_mat.c
blob37468ee30efae9d2597b94ccadeb25e012c55fe4
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
1405 bset = isl_basic_set_simplify(bset);
1406 bset = isl_basic_set_finalize(bset);
1407
1408 return bset;
1409 error:
1410 isl_mat_free(mat);
1411 error2:
1412 isl_basic_set_free(bset);
1413 return NULL;
1414 }
1415
1416 __isl_give isl_set *isl_set_preimage(
1417 __isl_take isl_set *set, __isl_take isl_mat *mat)
1418 {
1419 int i;
1420
1421 set = isl_set_cow(set);
1422 if (!set)
1423 goto error;
1424
1425 for (i = 0; i < set->n; ++i) {
1426 set->p[i] = isl_basic_set_preimage(set->p[i],
1427 isl_mat_copy(mat));
1428 if (!set->p[i])
1429 goto error;
1430 }
1431 if (mat->n_col != mat->n_row) {
1432 set->dim = isl_space_cow(set->dim);
1433 if (!set->dim)
1434 goto error;
1435 set->dim->n_out += mat->n_col;
1436 set->dim->n_out -= mat->n_row;
1437 }
1438 isl_mat_free(mat);
1439 ISL_F_CLR(set, ISL_SET_NORMALIZED);
1440 return set;
1441 error:
1442 isl_set_free(set);
1443 isl_mat_free(mat);
1444 return NULL;
1445 }
1446
1447 /* Replace the variables x starting at "first_col" in the rows "rows"
1448 * of some coefficient matrix by x' with x = M x' with M the matrix mat.
1449 * That is, replace the corresponding coefficients c by c M.
1450 */
1451 isl_stat isl_mat_sub_transform(isl_int **row, unsigned n_row,
1452 unsigned first_col, __isl_take isl_mat *mat)
1453 {
1454 int i;
1455 isl_ctx *ctx;
1456 isl_mat *t;
1457
1458 if (!mat)
1459 return isl_stat_error;
1460 ctx = isl_mat_get_ctx(mat);
1461 t = isl_mat_sub_alloc6(ctx, row, 0, n_row, first_col, mat->n_row);
1462 t = isl_mat_product(t, mat);
1463 if (!t)
1464 return isl_stat_error;
1465 for (i = 0; i < n_row; ++i)
1466 isl_seq_swp_or_cpy(row[i] + first_col, t->row[i], t->n_col);
1467 isl_mat_free(t);
1468 return isl_stat_ok;
1469 }
1470
1471 void isl_mat_print_internal(__isl_keep isl_mat *mat, FILE *out, int indent)
1472 {
1473 int i, j;
1474
1475 if (!mat) {
1476 fprintf(out, "%*snull mat\n", indent, "");
1477 return;
1478 }
1479
1480 if (mat->n_row == 0)
1481 fprintf(out, "%*s[]\n", indent, "");
1482
1483 for (i = 0; i < mat->n_row; ++i) {
1484 if (!i)
1485 fprintf(out, "%*s[[", indent, "");
1486 else
1487 fprintf(out, "%*s[", indent+1, "");
1488 for (j = 0; j < mat->n_col; ++j) {
1489 if (j)
1490 fprintf(out, ",");
1491 isl_int_print(out, mat->row[i][j], 0);
1492 }
1493 if (i == mat->n_row-1)
1494 fprintf(out, "]]\n");
1495 else
1496 fprintf(out, "]\n");
1497 }
1498 }
1499
1500 void isl_mat_dump(__isl_keep isl_mat *mat)
1501 {
1502 isl_mat_print_internal(mat, stderr, 0);
1503 }
1504
1505 __isl_give isl_mat *isl_mat_drop_cols(__isl_take isl_mat *mat,
1506 unsigned col, unsigned n)
1507 {
1508 int r;
1509
1510 if (n == 0)
1511 return mat;
1512
1513 mat = isl_mat_cow(mat);
1514 if (check_col_range(mat, col, n) < 0)
1515 return isl_mat_free(mat);
1516
1517 if (col != mat->n_col-n) {
1518 for (r = 0; r < mat->n_row; ++r)
1519 isl_seq_cpy(mat->row[r]+col, mat->row[r]+col+n,
1520 mat->n_col - col - n);
1521 }
1522 mat->n_col -= n;
1523 return mat;
1524 }
1525
1526 __isl_give isl_mat *isl_mat_drop_rows(__isl_take isl_mat *mat,
1527 unsigned row, unsigned n)
1528 {
1529 int r;
1530
1531 mat = isl_mat_cow(mat);
1532 if (check_row_range(mat, row, n) < 0)
1533 return isl_mat_free(mat);
1534
1535 for (r = row; r+n < mat->n_row; ++r)
1536 mat->row[r] = mat->row[r+n];
1537
1538 mat->n_row -= n;
1539 return mat;
1540 }
1541
1542 __isl_give isl_mat *isl_mat_insert_cols(__isl_take isl_mat *mat,
1543 unsigned col, unsigned n)
1544 {
1545 isl_mat *ext;
1546
1547 if (check_col_range(mat, col, 0) < 0)
1548 return isl_mat_free(mat);
1549 if (n == 0)
1550 return mat;
1551
1552 ext = isl_mat_alloc(mat->ctx, mat->n_row, mat->n_col + n);
1553 if (!ext)
1554 goto error;
1555
1556 isl_mat_sub_copy(mat->ctx, ext->row, mat->row, mat->n_row, 0, 0, col);
1557 isl_mat_sub_copy(mat->ctx, ext->row, mat->row, mat->n_row,
1558 col + n, col, mat->n_col - col);
1559
1560 isl_mat_free(mat);
1561 return ext;
1562 error:
1563 isl_mat_free(mat);
1564 return NULL;
1565 }
1566
1567 __isl_give isl_mat *isl_mat_insert_zero_cols(__isl_take isl_mat *mat,
1568 unsigned first, unsigned n)
1569 {
1570 int i;
1571
1572 if (!mat)
1573 return NULL;
1574 mat = isl_mat_insert_cols(mat, first, n);
1575 if (!mat)
1576 return NULL;
1577
1578 for (i = 0; i < mat->n_row; ++i)
1579 isl_seq_clr(mat->row[i] + first, n);
1580
1581 return mat;
1582 }
1583
1584 __isl_give isl_mat *isl_mat_add_zero_cols(__isl_take isl_mat *mat, unsigned n)
1585 {
1586 if (!mat)
1587 return NULL;
1588
1589 return isl_mat_insert_zero_cols(mat, mat->n_col, n);
1590 }
1591
1592 __isl_give isl_mat *isl_mat_insert_rows(__isl_take isl_mat *mat,
1593 unsigned row, unsigned n)
1594 {
1595 isl_mat *ext;
1596
1597 if (check_row_range(mat, row, 0) < 0)
1598 return isl_mat_free(mat);
1599 if (n == 0)
1600 return mat;
1601
1602 ext = isl_mat_alloc(mat->ctx, mat->n_row + n, mat->n_col);
1603 if (!ext)
1604 goto error;
1605
1606 isl_mat_sub_copy(mat->ctx, ext->row, mat->row, row, 0, 0, mat->n_col);
1607 isl_mat_sub_copy(mat->ctx, ext->row + row + n, mat->row + row,
1608 mat->n_row - row, 0, 0, mat->n_col);
1609
1610 isl_mat_free(mat);
1611 return ext;
1612 error:
1613 isl_mat_free(mat);
1614 return NULL;
1615 }
1616
1617 __isl_give isl_mat *isl_mat_add_rows(__isl_take isl_mat *mat, unsigned n)
1618 {
1619 if (!mat)
1620 return NULL;
1621
1622 return isl_mat_insert_rows(mat, mat->n_row, n);
1623 }
1624
1625 __isl_give isl_mat *isl_mat_insert_zero_rows(__isl_take isl_mat *mat,
1626 unsigned row, unsigned n)
1627 {
1628 int i;
1629
1630 mat = isl_mat_insert_rows(mat, row, n);
1631 if (!mat)
1632 return NULL;
1633
1634 for (i = 0; i < n; ++i)
1635 isl_seq_clr(mat->row[row + i], mat->n_col);
1636
1637 return mat;
1638 }
1639
1640 __isl_give isl_mat *isl_mat_add_zero_rows(__isl_take isl_mat *mat, unsigned n)
1641 {
1642 if (!mat)
1643 return NULL;
1644
1645 return isl_mat_insert_zero_rows(mat, mat->n_row, n);
1646 }
1647
1648 void isl_mat_col_submul(__isl_keep isl_mat *mat,
1649 int dst_col, isl_int f, int src_col)
1650 {
1651 int i;
1652
1653 for (i = 0; i < mat->n_row; ++i)
1654 isl_int_submul(mat->row[i][dst_col], f, mat->row[i][src_col]);
1655 }
1656
1657 void isl_mat_col_add(__isl_keep isl_mat *mat, int dst_col, int src_col)
1658 {
1659 int i;
1660
1661 if (!mat)
1662 return;
1663
1664 for (i = 0; i < mat->n_row; ++i)
1665 isl_int_add(mat->row[i][dst_col],
1666 mat->row[i][dst_col], mat->row[i][src_col]);
1667 }
1668
1669 void isl_mat_col_mul(__isl_keep isl_mat *mat, int dst_col, isl_int f,
1670 int src_col)
1671 {
1672 int i;
1673
1674 for (i = 0; i < mat->n_row; ++i)
1675 isl_int_mul(mat->row[i][dst_col], f, mat->row[i][src_col]);
1676 }
1677
1678 /* Add "f" times column "src_col" to column "dst_col" of "mat" and
1679 * return the result.
1680 */
1681 __isl_give isl_mat *isl_mat_col_addmul(__isl_take isl_mat *mat, int dst_col,
1682 isl_int f, int src_col)
1683 {
1684 int i;
1685
1686 if (check_col(mat, dst_col) < 0 || check_col(mat, src_col) < 0)
1687 return isl_mat_free(mat);
1688
1689 for (i = 0; i < mat->n_row; ++i) {
1690 if (isl_int_is_zero(mat->row[i][src_col]))
1691 continue;
1692 mat = isl_mat_cow(mat);
1693 if (!mat)
1694 return NULL;
1695 isl_int_addmul(mat->row[i][dst_col], f, mat->row[i][src_col]);
1696 }
1697
1698 return mat;
1699 }
1700
1701 /* Negate column "col" of "mat" and return the result.
1702 */
1703 __isl_give isl_mat *isl_mat_col_neg(__isl_take isl_mat *mat, int col)
1704 {
1705 int i;
1706
1707 if (check_col(mat, col) < 0)
1708 return isl_mat_free(mat);
1709
1710 for (i = 0; i < mat->n_row; ++i) {
1711 if (isl_int_is_zero(mat->row[i][col]))
1712 continue;
1713 mat = isl_mat_cow(mat);
1714 if (!mat)
1715 return NULL;
1716 isl_int_neg(mat->row[i][col], mat->row[i][col]);
1717 }
1718
1719 return mat;
1720 }
1721
1722 /* Negate row "row" of "mat" and return the result.
1723 */
1724 __isl_give isl_mat *isl_mat_row_neg(__isl_take isl_mat *mat, int row)
1725 {
1726 if (check_row(mat, row) < 0)
1727 return isl_mat_free(mat);
1728 if (isl_seq_first_non_zero(mat->row[row], mat->n_col) == -1)
1729 return mat;
1730 mat = isl_mat_cow(mat);
1731 if (!mat)
1732 return NULL;
1733 isl_seq_neg(mat->row[row], mat->row[row], mat->n_col);
1734 return mat;
1735 }
1736
1737 __isl_give isl_mat *isl_mat_unimodular_complete(__isl_take isl_mat *M, int row)
1738 {
1739 int r;
1740 struct isl_mat *H = NULL, *Q = NULL;
1741
1742 if (!M)
1743 return NULL;
1744
1745 isl_assert(M->ctx, M->n_row == M->n_col, goto error);
1746 M->n_row = row;
1747 H = isl_mat_left_hermite(isl_mat_copy(M), 0, NULL, &Q);
1748 M->n_row = M->n_col;
1749 if (!H)
1750 goto error;
1751 for (r = 0; r < row; ++r)
1752 isl_assert(M->ctx, isl_int_is_one(H->row[r][r]), goto error);
1753 for (r = row; r < M->n_row; ++r)
1754 isl_seq_cpy(M->row[r], Q->row[r], M->n_col);
1755 isl_mat_free(H);
1756 isl_mat_free(Q);
1757 return M;
1758 error:
1759 isl_mat_free(H);
1760 isl_mat_free(Q);
1761 isl_mat_free(M);
1762 return NULL;
1763 }
1764
1765 __isl_give isl_mat *isl_mat_concat(__isl_take isl_mat *top,
1766 __isl_take isl_mat *bot)
1767 {
1768 struct isl_mat *mat;
1769
1770 if (!top || !bot)
1771 goto error;
1772
1773 isl_assert(top->ctx, top->n_col == bot->n_col, goto error);
1774 if (top->n_row == 0) {
1775 isl_mat_free(top);
1776 return bot;
1777 }
1778 if (bot->n_row == 0) {
1779 isl_mat_free(bot);
1780 return top;
1781 }
1782
1783 mat = isl_mat_alloc(top->ctx, top->n_row + bot->n_row, top->n_col);
1784 if (!mat)
1785 goto error;
1786 isl_mat_sub_copy(mat->ctx, mat->row, top->row, top->n_row,
1787 0, 0, mat->n_col);
1788 isl_mat_sub_copy(mat->ctx, mat->row + top->n_row, bot->row, bot->n_row,
1789 0, 0, mat->n_col);
1790 isl_mat_free(top);
1791 isl_mat_free(bot);
1792 return mat;
1793 error:
1794 isl_mat_free(top);
1795 isl_mat_free(bot);
1796 return NULL;
1797 }
1798
1799 isl_bool isl_mat_is_equal(__isl_keep isl_mat *mat1, __isl_keep isl_mat *mat2)
1800 {
1801 int i;
1802
1803 if (!mat1 || !mat2)
1804 return isl_bool_error;
1805
1806 if (mat1->n_row != mat2->n_row)
1807 return isl_bool_false;
1808
1809 if (mat1->n_col != mat2->n_col)
1810 return isl_bool_false;
1811
1812 for (i = 0; i < mat1->n_row; ++i)
1813 if (!isl_seq_eq(mat1->row[i], mat2->row[i], mat1->n_col))
1814 return isl_bool_false;
1815
1816 return isl_bool_true;
1817 }
1818
1819 __isl_give isl_mat *isl_mat_from_row_vec(__isl_take isl_vec *vec)
1820 {
1821 struct isl_mat *mat;
1822
1823 if (!vec)
1824 return NULL;
1825 mat = isl_mat_alloc(vec->ctx, 1, vec->size);
1826 if (!mat)
1827 goto error;
1828
1829 isl_seq_cpy(mat->row[0], vec->el, vec->size);
1830
1831 isl_vec_free(vec);
1832 return mat;
1833 error:
1834 isl_vec_free(vec);
1835 return NULL;
1836 }
1837
1838 /* Return a copy of row "row" of "mat" as an isl_vec.
1839 */
1840 __isl_give isl_vec *isl_mat_get_row(__isl_keep isl_mat *mat, unsigned row)
1841 {
1842 isl_vec *v;
1843
1844 if (!mat)
1845 return NULL;
1846 if (row >= mat->n_row)
1847 isl_die(mat->ctx, isl_error_invalid, "row out of range",
1848 return NULL);
1849
1850 v = isl_vec_alloc(isl_mat_get_ctx(mat), mat->n_col);
1851 if (!v)
1852 return NULL;
1853 isl_seq_cpy(v->el, mat->row[row], mat->n_col);
1854
1855 return v;
1856 }
1857
1858 __isl_give isl_mat *isl_mat_vec_concat(__isl_take isl_mat *top,
1859 __isl_take isl_vec *bot)
1860 {
1861 return isl_mat_concat(top, isl_mat_from_row_vec(bot));
1862 }
1863
1864 __isl_give isl_mat *isl_mat_move_cols(__isl_take isl_mat *mat,
1865 unsigned dst_col, unsigned src_col, unsigned n)
1866 {
1867 isl_mat *res;
1868
1869 if (!mat)
1870 return NULL;
1871 if (n == 0 || dst_col == src_col)
1872 return mat;
1873
1874 res = isl_mat_alloc(mat->ctx, mat->n_row, mat->n_col);
1875 if (!res)
1876 goto error;
1877
1878 if (dst_col < src_col) {
1879 isl_mat_sub_copy(res->ctx, res->row, mat->row, mat->n_row,
1880 0, 0, dst_col);
1881 isl_mat_sub_copy(res->ctx, res->row, mat->row, mat->n_row,
1882 dst_col, src_col, n);
1883 isl_mat_sub_copy(res->ctx, res->row, mat->row, mat->n_row,
1884 dst_col + n, dst_col, src_col - dst_col);
1885 isl_mat_sub_copy(res->ctx, res->row, mat->row, mat->n_row,
1886 src_col + n, src_col + n,
1887 res->n_col - src_col - n);
1888 } else {
1889 isl_mat_sub_copy(res->ctx, res->row, mat->row, mat->n_row,
1890 0, 0, src_col);
1891 isl_mat_sub_copy(res->ctx, res->row, mat->row, mat->n_row,
1892 src_col, src_col + n, dst_col - src_col);
1893 isl_mat_sub_copy(res->ctx, res->row, mat->row, mat->n_row,
1894 dst_col, src_col, n);
1895 isl_mat_sub_copy(res->ctx, res->row, mat->row, mat->n_row,
1896 dst_col + n, dst_col + n,
1897 res->n_col - dst_col - n);
1898 }
1899 isl_mat_free(mat);
1900
1901 return res;
1902 error:
1903 isl_mat_free(mat);
1904 return NULL;
1905 }
1906
1907 /* Return the gcd of the elements in row "row" of "mat" in *gcd.
1908 * Return isl_stat_ok on success and isl_stat_error on failure.
1909 */
1910 isl_stat isl_mat_row_gcd(__isl_keep isl_mat *mat, int row, isl_int *gcd)
1911 {
1912 if (check_row(mat, row) < 0)
1913 return isl_stat_error;
1914
1915 isl_seq_gcd(mat->row[row], mat->n_col, gcd);
1916
1917 return isl_stat_ok;
1918 }
1919
1920 void isl_mat_gcd(__isl_keep isl_mat *mat, isl_int *gcd)
1921 {
1922 int i;
1923 isl_int g;
1924
1925 isl_int_set_si(*gcd, 0);
1926 if (!mat)
1927 return;
1928
1929 isl_int_init(g);
1930 for (i = 0; i < mat->n_row; ++i) {
1931 isl_seq_gcd(mat->row[i], mat->n_col, &g);
1932 isl_int_gcd(*gcd, *gcd, g);
1933 }
1934 isl_int_clear(g);
1935 }
1936
1937 /* Return the result of scaling "mat" by a factor of "m".
1938 */
1939 __isl_give isl_mat *isl_mat_scale(__isl_take isl_mat *mat, isl_int m)
1940 {
1941 int i;
1942
1943 if (isl_int_is_one(m))
1944 return mat;
1945
1946 mat = isl_mat_cow(mat);
1947 if (!mat)
1948 return NULL;
1949
1950 for (i = 0; i < mat->n_row; ++i)
1951 isl_seq_scale(mat->row[i], mat->row[i], m, mat->n_col);
1952
1953 return mat;
1954 }
1955
1956 __isl_give isl_mat *isl_mat_scale_down(__isl_take isl_mat *mat, isl_int m)
1957 {
1958 int i;
1959
1960 if (isl_int_is_one(m))
1961 return mat;
1962
1963 mat = isl_mat_cow(mat);
1964 if (!mat)
1965 return NULL;
1966
1967 for (i = 0; i < mat->n_row; ++i)
1968 isl_seq_scale_down(mat->row[i], mat->row[i], m, mat->n_col);
1969
1970 return mat;
1971 }
1972
1973 __isl_give isl_mat *isl_mat_scale_down_row(__isl_take isl_mat *mat, int row,
1974 isl_int m)
1975 {
1976 if (isl_int_is_one(m))
1977 return mat;
1978
1979 mat = isl_mat_cow(mat);
1980 if (!mat)
1981 return NULL;
1982
1983 isl_seq_scale_down(mat->row[row], mat->row[row], m, mat->n_col);
1984
1985 return mat;
1986 }
1987
1988 __isl_give isl_mat *isl_mat_normalize(__isl_take isl_mat *mat)
1989 {
1990 isl_int gcd;
1991
1992 if (!mat)
1993 return NULL;
1994
1995 isl_int_init(gcd);
1996 isl_mat_gcd(mat, &gcd);
1997 mat = isl_mat_scale_down(mat, gcd);
1998 isl_int_clear(gcd);
1999
2000 return mat;
2001 }
2002
2003 __isl_give isl_mat *isl_mat_normalize_row(__isl_take isl_mat *mat, int row)
2004 {
2005 mat = isl_mat_cow(mat);
2006 if (!mat)
2007 return NULL;
2008
2009 isl_seq_normalize(mat->ctx, mat->row[row], mat->n_col);
2010
2011 return mat;
2012 }
2013
2014 /* Number of initial non-zero columns.
2015 */
2016 int isl_mat_initial_non_zero_cols(__isl_keep isl_mat *mat)
2017 {
2018 int i;
2019
2020 if (!mat)
2021 return -1;
2022
2023 for (i = 0; i < mat->n_col; ++i)
2024 if (row_first_non_zero(mat->row, mat->n_row, i) < 0)
2025 break;
2026
2027 return i;
2028 }
2029
2030 /* Return a basis for the space spanned by the rows of "mat".
2031 * Any basis will do, so simply perform Gaussian elimination and
2032 * remove the empty rows.
2033 */
2034 __isl_give isl_mat *isl_mat_row_basis(__isl_take isl_mat *mat)
2035 {
2036 return isl_mat_reverse_gauss(mat);
2037 }
2038
2039 /* Return rows that extend a basis of "mat1" to one
2040 * that covers both "mat1" and "mat2".
2041 * The Hermite normal form of the concatenation of the two matrices is
2042 *
2043 * [ Q1 ]
2044 * [ M1 ] = [ H1 0 0 ] [ Q2 ]
2045 * [ M2 ] = [ H2 H3 0 ] [ Q3 ]
2046 *
2047 * The number of columns in H1 and H3 determine the number of rows
2048 * in Q1 and Q2. Q1 is a basis for M1, while Q2 extends this basis
2049 * to also cover M2.
2050 */
2051 __isl_give isl_mat *isl_mat_row_basis_extension(
2052 __isl_take isl_mat *mat1, __isl_take isl_mat *mat2)
2053 {
2054 isl_size n_row;
2055 int r1, r;
2056 isl_size n1;
2057 isl_mat *H, *Q;
2058
2059 n1 = isl_mat_rows(mat1);
2060 H = isl_mat_concat(mat1, mat2);
2061 H = isl_mat_left_hermite(H, 0, NULL, &Q);
2062 if (n1 < 0 || !H || !Q)
2063 goto error;
2064
2065 r1 = hermite_first_zero_col(H, 0, n1);
2066 r = hermite_first_zero_col(H, r1, H->n_row);
2067 n_row = isl_mat_rows(Q);
2068 if (n_row < 0)
2069 goto error;
2070 Q = isl_mat_drop_rows(Q, r, n_row - r);
2071 Q = isl_mat_drop_rows(Q, 0, r1);
2072
2073 isl_mat_free(H);
2074 return Q;
2075 error:
2076 isl_mat_free(H);
2077 isl_mat_free(Q);
2078 return NULL;
2079 }
2080
2081 /* Are the rows of "mat1" linearly independent of those of "mat2"?
2082 * That is, is there no linear dependence among the combined rows
2083 * that is not already present in either "mat1" or "mat2"?
2084 * In other words, is the rank of "mat1" and "mat2" combined equal
2085 * to the sum of the ranks of "mat1" and "mat2"?
2086 */
2087 isl_bool isl_mat_has_linearly_independent_rows(__isl_keep isl_mat *mat1,
2088 __isl_keep isl_mat *mat2)
2089 {
2090 isl_size r1, r2, r;
2091 isl_mat *mat;
2092
2093 r1 = isl_mat_rank(mat1);
2094 if (r1 < 0)
2095 return isl_bool_error;
2096 if (r1 == 0)
2097 return isl_bool_true;
2098 r2 = isl_mat_rank(mat2);
2099 if (r2 < 0)
2100 return isl_bool_error;
2101 if (r2 == 0)
2102 return isl_bool_true;
2103
2104 mat = isl_mat_concat(isl_mat_copy(mat1), isl_mat_copy(mat2));
2105 r = isl_mat_rank(mat);
2106 isl_mat_free(mat);
2107 if (r < 0)
2108 return isl_bool_error;
2109 return isl_bool_ok(r == r1 + r2);
2110 }
Integer set library