2 * Copyright 2008-2009 Katholieke Universiteit Leuven
3 * Copyright 2014 Ecole Normale Superieure
5 * Use of this software is governed by the MIT license
7 * Written by Sven Verdoolaege, K.U.Leuven, Departement
8 * Computerwetenschappen, Celestijnenlaan 200A, B-3001 Leuven, Belgium
9 * and Ecole Normale Superieure, 45 rue d'Ulm, 75230 Paris, France
12 #include <isl_ctx_private.h>
13 #include <isl_map_private.h>
14 #include <isl/space.h>
16 #include <isl_mat_private.h>
17 #include <isl_vec_private.h>
18 #include <isl_space_private.h>
19 #include <isl_val_private.h>
20 #include <isl/deprecated/mat_int.h>
22 isl_ctx
*isl_mat_get_ctx(__isl_keep isl_mat
*mat
)
24 return mat
? mat
->ctx
: NULL
;
27 /* Return a hash value that digests "mat".
29 uint32_t isl_mat_get_hash(__isl_keep isl_mat
*mat
)
37 hash
= isl_hash_init();
38 isl_hash_byte(hash
, mat
->n_row
& 0xFF);
39 isl_hash_byte(hash
, mat
->n_col
& 0xFF);
40 for (i
= 0; i
< mat
->n_row
; ++i
) {
43 row_hash
= isl_seq_get_hash(mat
->row
[i
], mat
->n_col
);
44 isl_hash_hash(hash
, row_hash
);
50 struct isl_mat
*isl_mat_alloc(struct isl_ctx
*ctx
,
51 unsigned n_row
, unsigned n_col
)
56 mat
= isl_alloc_type(ctx
, struct isl_mat
);
61 mat
->block
= isl_blk_alloc(ctx
, n_row
* n_col
);
62 if (isl_blk_is_error(mat
->block
))
64 mat
->row
= isl_alloc_array(ctx
, isl_int
*, n_row
);
65 if (n_row
&& !mat
->row
)
68 for (i
= 0; i
< n_row
; ++i
)
69 mat
->row
[i
] = mat
->block
.data
+ i
* n_col
;
81 isl_blk_free(ctx
, mat
->block
);
86 struct isl_mat
*isl_mat_extend(struct isl_mat
*mat
,
87 unsigned n_row
, unsigned n_col
)
96 if (mat
->max_col
>= n_col
&& mat
->n_row
>= n_row
) {
97 if (mat
->n_col
< n_col
)
102 if (mat
->max_col
< n_col
) {
103 struct isl_mat
*new_mat
;
105 if (n_row
< mat
->n_row
)
107 new_mat
= isl_mat_alloc(mat
->ctx
, n_row
, n_col
);
110 for (i
= 0; i
< mat
->n_row
; ++i
)
111 isl_seq_cpy(new_mat
->row
[i
], mat
->row
[i
], mat
->n_col
);
116 mat
= isl_mat_cow(mat
);
120 old
= mat
->block
.data
;
121 mat
->block
= isl_blk_extend(mat
->ctx
, mat
->block
, n_row
* mat
->max_col
);
122 if (isl_blk_is_error(mat
->block
))
124 row
= isl_realloc_array(mat
->ctx
, mat
->row
, isl_int
*, n_row
);
129 for (i
= 0; i
< mat
->n_row
; ++i
)
130 mat
->row
[i
] = mat
->block
.data
+ (mat
->row
[i
] - old
);
131 for (i
= mat
->n_row
; i
< n_row
; ++i
)
132 mat
->row
[i
] = mat
->block
.data
+ i
* mat
->max_col
;
134 if (mat
->n_col
< n_col
)
143 __isl_give isl_mat
*isl_mat_sub_alloc6(isl_ctx
*ctx
, isl_int
**row
,
144 unsigned first_row
, unsigned n_row
, unsigned first_col
, unsigned n_col
)
149 mat
= isl_alloc_type(ctx
, struct isl_mat
);
152 mat
->row
= isl_alloc_array(ctx
, isl_int
*, n_row
);
153 if (n_row
&& !mat
->row
)
155 for (i
= 0; i
< n_row
; ++i
)
156 mat
->row
[i
] = row
[first_row
+i
] + first_col
;
162 mat
->block
= isl_blk_empty();
163 mat
->flags
= ISL_MAT_BORROWED
;
170 __isl_give isl_mat
*isl_mat_sub_alloc(__isl_keep isl_mat
*mat
,
171 unsigned first_row
, unsigned n_row
, unsigned first_col
, unsigned n_col
)
175 return isl_mat_sub_alloc6(mat
->ctx
, mat
->row
, first_row
, n_row
,
179 void isl_mat_sub_copy(struct isl_ctx
*ctx
, isl_int
**dst
, isl_int
**src
,
180 unsigned n_row
, unsigned dst_col
, unsigned src_col
, unsigned n_col
)
184 for (i
= 0; i
< n_row
; ++i
)
185 isl_seq_cpy(dst
[i
]+dst_col
, src
[i
]+src_col
, n_col
);
188 void isl_mat_sub_neg(struct isl_ctx
*ctx
, isl_int
**dst
, isl_int
**src
,
189 unsigned n_row
, unsigned dst_col
, unsigned src_col
, unsigned n_col
)
193 for (i
= 0; i
< n_row
; ++i
)
194 isl_seq_neg(dst
[i
]+dst_col
, src
[i
]+src_col
, n_col
);
197 __isl_give isl_mat
*isl_mat_copy(__isl_keep isl_mat
*mat
)
206 struct isl_mat
*isl_mat_dup(struct isl_mat
*mat
)
209 struct isl_mat
*mat2
;
213 mat2
= isl_mat_alloc(mat
->ctx
, mat
->n_row
, mat
->n_col
);
216 for (i
= 0; i
< mat
->n_row
; ++i
)
217 isl_seq_cpy(mat2
->row
[i
], mat
->row
[i
], mat
->n_col
);
221 struct isl_mat
*isl_mat_cow(struct isl_mat
*mat
)
223 struct isl_mat
*mat2
;
227 if (mat
->ref
== 1 && !ISL_F_ISSET(mat
, ISL_MAT_BORROWED
))
230 mat2
= isl_mat_dup(mat
);
235 __isl_null isl_mat
*isl_mat_free(__isl_take isl_mat
*mat
)
243 if (!ISL_F_ISSET(mat
, ISL_MAT_BORROWED
))
244 isl_blk_free(mat
->ctx
, mat
->block
);
245 isl_ctx_deref(mat
->ctx
);
252 int isl_mat_rows(__isl_keep isl_mat
*mat
)
254 return mat
? mat
->n_row
: -1;
257 int isl_mat_cols(__isl_keep isl_mat
*mat
)
259 return mat
? mat
->n_col
: -1;
262 /* Check that "col" is a valid column position for "mat".
264 static isl_stat
check_col(__isl_keep isl_mat
*mat
, int col
)
267 return isl_stat_error
;
268 if (col
< 0 || col
>= mat
->n_col
)
269 isl_die(isl_mat_get_ctx(mat
), isl_error_invalid
,
270 "column out of range", return isl_stat_error
);
274 int isl_mat_get_element(__isl_keep isl_mat
*mat
, int row
, int col
, isl_int
*v
)
278 if (row
< 0 || row
>= mat
->n_row
)
279 isl_die(mat
->ctx
, isl_error_invalid
, "row out of range",
281 if (check_col(mat
, col
) < 0)
283 isl_int_set(*v
, mat
->row
[row
][col
]);
287 /* Extract the element at row "row", oolumn "col" of "mat".
289 __isl_give isl_val
*isl_mat_get_element_val(__isl_keep isl_mat
*mat
,
296 ctx
= isl_mat_get_ctx(mat
);
297 if (row
< 0 || row
>= mat
->n_row
)
298 isl_die(ctx
, isl_error_invalid
, "row out of range",
300 if (check_col(mat
, col
) < 0)
302 return isl_val_int_from_isl_int(ctx
, mat
->row
[row
][col
]);
305 __isl_give isl_mat
*isl_mat_set_element(__isl_take isl_mat
*mat
,
306 int row
, int col
, isl_int v
)
308 mat
= isl_mat_cow(mat
);
311 if (row
< 0 || row
>= mat
->n_row
)
312 isl_die(mat
->ctx
, isl_error_invalid
, "row out of range",
314 if (check_col(mat
, col
) < 0)
315 return isl_mat_free(mat
);
316 isl_int_set(mat
->row
[row
][col
], v
);
323 __isl_give isl_mat
*isl_mat_set_element_si(__isl_take isl_mat
*mat
,
324 int row
, int col
, int v
)
326 mat
= isl_mat_cow(mat
);
329 if (row
< 0 || row
>= mat
->n_row
)
330 isl_die(mat
->ctx
, isl_error_invalid
, "row out of range",
332 if (check_col(mat
, col
) < 0)
333 return isl_mat_free(mat
);
334 isl_int_set_si(mat
->row
[row
][col
], v
);
341 /* Replace the element at row "row", column "col" of "mat" by "v".
343 __isl_give isl_mat
*isl_mat_set_element_val(__isl_take isl_mat
*mat
,
344 int row
, int col
, __isl_take isl_val
*v
)
347 return isl_mat_free(mat
);
348 if (!isl_val_is_int(v
))
349 isl_die(isl_val_get_ctx(v
), isl_error_invalid
,
350 "expecting integer value", goto error
);
351 mat
= isl_mat_set_element(mat
, row
, col
, v
->n
);
356 return isl_mat_free(mat
);
359 __isl_give isl_mat
*isl_mat_diag(isl_ctx
*ctx
, unsigned n_row
, isl_int d
)
364 mat
= isl_mat_alloc(ctx
, n_row
, n_row
);
367 for (i
= 0; i
< n_row
; ++i
) {
368 isl_seq_clr(mat
->row
[i
], i
);
369 isl_int_set(mat
->row
[i
][i
], d
);
370 isl_seq_clr(mat
->row
[i
]+i
+1, n_row
-(i
+1));
376 /* Create an "n_row" by "n_col" matrix with zero elements.
378 __isl_give isl_mat
*isl_mat_zero(isl_ctx
*ctx
, unsigned n_row
, unsigned n_col
)
383 mat
= isl_mat_alloc(ctx
, n_row
, n_col
);
386 for (i
= 0; i
< n_row
; ++i
)
387 isl_seq_clr(mat
->row
[i
], n_col
);
392 __isl_give isl_mat
*isl_mat_identity(isl_ctx
*ctx
, unsigned n_row
)
396 return isl_mat_diag(ctx
, n_row
, ctx
->one
);
399 /* Is "mat" a (possibly scaled) identity matrix?
401 int isl_mat_is_scaled_identity(__isl_keep isl_mat
*mat
)
407 if (mat
->n_row
!= mat
->n_col
)
410 for (i
= 0; i
< mat
->n_row
; ++i
) {
411 if (isl_seq_first_non_zero(mat
->row
[i
], i
) != -1)
413 if (isl_int_ne(mat
->row
[0][0], mat
->row
[i
][i
]))
415 if (isl_seq_first_non_zero(mat
->row
[i
] + i
+ 1,
416 mat
->n_col
- (i
+ 1)) != -1)
423 __isl_give isl_vec
*isl_mat_vec_product(__isl_take isl_mat
*mat
,
424 __isl_take isl_vec
*vec
)
427 struct isl_vec
*prod
;
432 isl_assert(mat
->ctx
, mat
->n_col
== vec
->size
, goto error
);
434 prod
= isl_vec_alloc(mat
->ctx
, mat
->n_row
);
438 for (i
= 0; i
< prod
->size
; ++i
)
439 isl_seq_inner_product(mat
->row
[i
], vec
->el
, vec
->size
,
440 &prod
->block
.data
[i
]);
450 __isl_give isl_vec
*isl_mat_vec_inverse_product(__isl_take isl_mat
*mat
,
451 __isl_take isl_vec
*vec
)
453 struct isl_mat
*vec_mat
;
458 vec_mat
= isl_mat_alloc(vec
->ctx
, vec
->size
, 1);
461 for (i
= 0; i
< vec
->size
; ++i
)
462 isl_int_set(vec_mat
->row
[i
][0], vec
->el
[i
]);
463 vec_mat
= isl_mat_inverse_product(mat
, vec_mat
);
467 vec
= isl_vec_alloc(vec_mat
->ctx
, vec_mat
->n_row
);
469 for (i
= 0; i
< vec
->size
; ++i
)
470 isl_int_set(vec
->el
[i
], vec_mat
->row
[i
][0]);
471 isl_mat_free(vec_mat
);
479 __isl_give isl_vec
*isl_vec_mat_product(__isl_take isl_vec
*vec
,
480 __isl_take isl_mat
*mat
)
483 struct isl_vec
*prod
;
488 isl_assert(mat
->ctx
, mat
->n_row
== vec
->size
, goto error
);
490 prod
= isl_vec_alloc(mat
->ctx
, mat
->n_col
);
494 for (i
= 0; i
< prod
->size
; ++i
) {
495 isl_int_set_si(prod
->el
[i
], 0);
496 for (j
= 0; j
< vec
->size
; ++j
)
497 isl_int_addmul(prod
->el
[i
], vec
->el
[j
], mat
->row
[j
][i
]);
508 __isl_give isl_mat
*isl_mat_aff_direct_sum(__isl_take isl_mat
*left
,
509 __isl_take isl_mat
*right
)
517 isl_assert(left
->ctx
, left
->n_row
== right
->n_row
, goto error
);
518 isl_assert(left
->ctx
, left
->n_row
>= 1, goto error
);
519 isl_assert(left
->ctx
, left
->n_col
>= 1, goto error
);
520 isl_assert(left
->ctx
, right
->n_col
>= 1, goto error
);
521 isl_assert(left
->ctx
,
522 isl_seq_first_non_zero(left
->row
[0]+1, left
->n_col
-1) == -1,
524 isl_assert(left
->ctx
,
525 isl_seq_first_non_zero(right
->row
[0]+1, right
->n_col
-1) == -1,
528 sum
= isl_mat_alloc(left
->ctx
, left
->n_row
, left
->n_col
+ right
->n_col
- 1);
531 isl_int_lcm(sum
->row
[0][0], left
->row
[0][0], right
->row
[0][0]);
532 isl_int_divexact(left
->row
[0][0], sum
->row
[0][0], left
->row
[0][0]);
533 isl_int_divexact(right
->row
[0][0], sum
->row
[0][0], right
->row
[0][0]);
535 isl_seq_clr(sum
->row
[0]+1, sum
->n_col
-1);
536 for (i
= 1; i
< sum
->n_row
; ++i
) {
537 isl_int_mul(sum
->row
[i
][0], left
->row
[0][0], left
->row
[i
][0]);
538 isl_int_addmul(sum
->row
[i
][0],
539 right
->row
[0][0], right
->row
[i
][0]);
540 isl_seq_scale(sum
->row
[i
]+1, left
->row
[i
]+1, left
->row
[0][0],
542 isl_seq_scale(sum
->row
[i
]+left
->n_col
,
543 right
->row
[i
]+1, right
->row
[0][0],
547 isl_int_divexact(left
->row
[0][0], sum
->row
[0][0], left
->row
[0][0]);
548 isl_int_divexact(right
->row
[0][0], sum
->row
[0][0], right
->row
[0][0]);
558 static void exchange(struct isl_mat
*M
, struct isl_mat
**U
,
559 struct isl_mat
**Q
, unsigned row
, unsigned i
, unsigned j
)
562 for (r
= row
; r
< M
->n_row
; ++r
)
563 isl_int_swap(M
->row
[r
][i
], M
->row
[r
][j
]);
565 for (r
= 0; r
< (*U
)->n_row
; ++r
)
566 isl_int_swap((*U
)->row
[r
][i
], (*U
)->row
[r
][j
]);
569 isl_mat_swap_rows(*Q
, i
, j
);
572 static void subtract(struct isl_mat
*M
, struct isl_mat
**U
,
573 struct isl_mat
**Q
, unsigned row
, unsigned i
, unsigned j
, isl_int m
)
576 for (r
= row
; r
< M
->n_row
; ++r
)
577 isl_int_submul(M
->row
[r
][j
], m
, M
->row
[r
][i
]);
579 for (r
= 0; r
< (*U
)->n_row
; ++r
)
580 isl_int_submul((*U
)->row
[r
][j
], m
, (*U
)->row
[r
][i
]);
583 for (r
= 0; r
< (*Q
)->n_col
; ++r
)
584 isl_int_addmul((*Q
)->row
[i
][r
], m
, (*Q
)->row
[j
][r
]);
588 static void oppose(struct isl_mat
*M
, struct isl_mat
**U
,
589 struct isl_mat
**Q
, unsigned row
, unsigned col
)
592 for (r
= row
; r
< M
->n_row
; ++r
)
593 isl_int_neg(M
->row
[r
][col
], M
->row
[r
][col
]);
595 for (r
= 0; r
< (*U
)->n_row
; ++r
)
596 isl_int_neg((*U
)->row
[r
][col
], (*U
)->row
[r
][col
]);
599 isl_seq_neg((*Q
)->row
[col
], (*Q
)->row
[col
], (*Q
)->n_col
);
602 /* Given matrix M, compute
607 * with U and Q unimodular matrices and H a matrix in column echelon form
608 * such that on each echelon row the entries in the non-echelon column
609 * are non-negative (if neg == 0) or non-positive (if neg == 1)
610 * and strictly smaller (in absolute value) than the entries in the echelon
612 * If U or Q are NULL, then these matrices are not computed.
614 __isl_give isl_mat
*isl_mat_left_hermite(__isl_take isl_mat
*M
, int neg
,
615 __isl_give isl_mat
**U
, __isl_give isl_mat
**Q
)
630 *U
= isl_mat_identity(M
->ctx
, M
->n_col
);
635 *Q
= isl_mat_identity(M
->ctx
, M
->n_col
);
642 for (row
= 0; row
< M
->n_row
; ++row
) {
644 first
= isl_seq_abs_min_non_zero(M
->row
[row
]+col
, M
->n_col
-col
);
649 exchange(M
, U
, Q
, row
, first
, col
);
650 if (isl_int_is_neg(M
->row
[row
][col
]))
651 oppose(M
, U
, Q
, row
, col
);
653 while ((off
= isl_seq_first_non_zero(M
->row
[row
]+first
,
654 M
->n_col
-first
)) != -1) {
656 isl_int_fdiv_q(c
, M
->row
[row
][first
], M
->row
[row
][col
]);
657 subtract(M
, U
, Q
, row
, col
, first
, c
);
658 if (!isl_int_is_zero(M
->row
[row
][first
]))
659 exchange(M
, U
, Q
, row
, first
, col
);
663 for (i
= 0; i
< col
; ++i
) {
664 if (isl_int_is_zero(M
->row
[row
][i
]))
667 isl_int_cdiv_q(c
, M
->row
[row
][i
], M
->row
[row
][col
]);
669 isl_int_fdiv_q(c
, M
->row
[row
][i
], M
->row
[row
][col
]);
670 if (isl_int_is_zero(c
))
672 subtract(M
, U
, Q
, row
, col
, i
, c
);
692 struct isl_mat
*isl_mat_right_kernel(struct isl_mat
*mat
)
695 struct isl_mat
*U
= NULL
;
698 mat
= isl_mat_left_hermite(mat
, 0, &U
, NULL
);
702 for (i
= 0, rank
= 0; rank
< mat
->n_col
; ++rank
) {
703 while (i
< mat
->n_row
&& isl_int_is_zero(mat
->row
[i
][rank
]))
708 K
= isl_mat_alloc(U
->ctx
, U
->n_row
, U
->n_col
- rank
);
711 isl_mat_sub_copy(K
->ctx
, K
->row
, U
->row
, U
->n_row
, 0, rank
, U
->n_col
-rank
);
721 __isl_give isl_mat
*isl_mat_lin_to_aff(__isl_take isl_mat
*mat
)
724 struct isl_mat
*mat2
;
728 mat2
= isl_mat_alloc(mat
->ctx
, 1+mat
->n_row
, 1+mat
->n_col
);
731 isl_int_set_si(mat2
->row
[0][0], 1);
732 isl_seq_clr(mat2
->row
[0]+1, mat
->n_col
);
733 for (i
= 0; i
< mat
->n_row
; ++i
) {
734 isl_int_set_si(mat2
->row
[1+i
][0], 0);
735 isl_seq_cpy(mat2
->row
[1+i
]+1, mat
->row
[i
], mat
->n_col
);
744 /* Given two matrices M1 and M2, return the block matrix
749 __isl_give isl_mat
*isl_mat_diagonal(__isl_take isl_mat
*mat1
,
750 __isl_take isl_mat
*mat2
)
758 mat
= isl_mat_alloc(mat1
->ctx
, mat1
->n_row
+ mat2
->n_row
,
759 mat1
->n_col
+ mat2
->n_col
);
762 for (i
= 0; i
< mat1
->n_row
; ++i
) {
763 isl_seq_cpy(mat
->row
[i
], mat1
->row
[i
], mat1
->n_col
);
764 isl_seq_clr(mat
->row
[i
] + mat1
->n_col
, mat2
->n_col
);
766 for (i
= 0; i
< mat2
->n_row
; ++i
) {
767 isl_seq_clr(mat
->row
[mat1
->n_row
+ i
], mat1
->n_col
);
768 isl_seq_cpy(mat
->row
[mat1
->n_row
+ i
] + mat1
->n_col
,
769 mat2
->row
[i
], mat2
->n_col
);
780 static int row_first_non_zero(isl_int
**row
, unsigned n_row
, unsigned col
)
784 for (i
= 0; i
< n_row
; ++i
)
785 if (!isl_int_is_zero(row
[i
][col
]))
790 static int row_abs_min_non_zero(isl_int
**row
, unsigned n_row
, unsigned col
)
792 int i
, min
= row_first_non_zero(row
, n_row
, col
);
795 for (i
= min
+ 1; i
< n_row
; ++i
) {
796 if (isl_int_is_zero(row
[i
][col
]))
798 if (isl_int_abs_lt(row
[i
][col
], row
[min
][col
]))
804 static isl_stat
inv_exchange(__isl_keep isl_mat
**left
,
805 __isl_keep isl_mat
**right
, unsigned i
, unsigned j
)
807 *left
= isl_mat_swap_rows(*left
, i
, j
);
808 *right
= isl_mat_swap_rows(*right
, i
, j
);
810 if (!*left
|| !*right
)
811 return isl_stat_error
;
815 static void inv_oppose(
816 struct isl_mat
*left
, struct isl_mat
*right
, unsigned row
)
818 isl_seq_neg(left
->row
[row
]+row
, left
->row
[row
]+row
, left
->n_col
-row
);
819 isl_seq_neg(right
->row
[row
], right
->row
[row
], right
->n_col
);
822 static void inv_subtract(struct isl_mat
*left
, struct isl_mat
*right
,
823 unsigned row
, unsigned i
, isl_int m
)
826 isl_seq_combine(left
->row
[i
]+row
,
827 left
->ctx
->one
, left
->row
[i
]+row
,
828 m
, left
->row
[row
]+row
,
830 isl_seq_combine(right
->row
[i
], right
->ctx
->one
, right
->row
[i
],
831 m
, right
->row
[row
], right
->n_col
);
834 /* Compute inv(left)*right
836 __isl_give isl_mat
*isl_mat_inverse_product(__isl_take isl_mat
*left
,
837 __isl_take isl_mat
*right
)
845 isl_assert(left
->ctx
, left
->n_row
== left
->n_col
, goto error
);
846 isl_assert(left
->ctx
, left
->n_row
== right
->n_row
, goto error
);
848 if (left
->n_row
== 0) {
853 left
= isl_mat_cow(left
);
854 right
= isl_mat_cow(right
);
860 for (row
= 0; row
< left
->n_row
; ++row
) {
861 int pivot
, first
, i
, off
;
862 pivot
= row_abs_min_non_zero(left
->row
+row
, left
->n_row
-row
, row
);
866 isl_assert(left
->ctx
, pivot
>= 0, goto error
);
870 if (inv_exchange(&left
, &right
, pivot
, row
) < 0)
872 if (isl_int_is_neg(left
->row
[row
][row
]))
873 inv_oppose(left
, right
, row
);
875 while ((off
= row_first_non_zero(left
->row
+first
,
876 left
->n_row
-first
, row
)) != -1) {
878 isl_int_fdiv_q(a
, left
->row
[first
][row
],
879 left
->row
[row
][row
]);
880 inv_subtract(left
, right
, row
, first
, a
);
881 if (!isl_int_is_zero(left
->row
[first
][row
])) {
882 if (inv_exchange(&left
, &right
, row
, first
) < 0)
888 for (i
= 0; i
< row
; ++i
) {
889 if (isl_int_is_zero(left
->row
[i
][row
]))
891 isl_int_gcd(a
, left
->row
[row
][row
], left
->row
[i
][row
]);
892 isl_int_divexact(b
, left
->row
[i
][row
], a
);
893 isl_int_divexact(a
, left
->row
[row
][row
], a
);
895 isl_seq_combine(left
->row
[i
] + i
,
897 b
, left
->row
[row
] + i
,
899 isl_seq_combine(right
->row
[i
], a
, right
->row
[i
],
900 b
, right
->row
[row
], right
->n_col
);
905 isl_int_set(a
, left
->row
[0][0]);
906 for (row
= 1; row
< left
->n_row
; ++row
)
907 isl_int_lcm(a
, a
, left
->row
[row
][row
]);
908 if (isl_int_is_zero(a
)){
910 isl_assert(left
->ctx
, 0, goto error
);
912 for (row
= 0; row
< left
->n_row
; ++row
) {
913 isl_int_divexact(left
->row
[row
][row
], a
, left
->row
[row
][row
]);
914 if (isl_int_is_one(left
->row
[row
][row
]))
916 isl_seq_scale(right
->row
[row
], right
->row
[row
],
917 left
->row
[row
][row
], right
->n_col
);
929 void isl_mat_col_scale(struct isl_mat
*mat
, unsigned col
, isl_int m
)
933 for (i
= 0; i
< mat
->n_row
; ++i
)
934 isl_int_mul(mat
->row
[i
][col
], mat
->row
[i
][col
], m
);
937 void isl_mat_col_combine(struct isl_mat
*mat
, unsigned dst
,
938 isl_int m1
, unsigned src1
, isl_int m2
, unsigned src2
)
944 for (i
= 0; i
< mat
->n_row
; ++i
) {
945 isl_int_mul(tmp
, m1
, mat
->row
[i
][src1
]);
946 isl_int_addmul(tmp
, m2
, mat
->row
[i
][src2
]);
947 isl_int_set(mat
->row
[i
][dst
], tmp
);
952 __isl_give isl_mat
*isl_mat_right_inverse(__isl_take isl_mat
*mat
)
958 mat
= isl_mat_cow(mat
);
962 inv
= isl_mat_identity(mat
->ctx
, mat
->n_col
);
963 inv
= isl_mat_cow(inv
);
969 for (row
= 0; row
< mat
->n_row
; ++row
) {
970 int pivot
, first
, i
, off
;
971 pivot
= isl_seq_abs_min_non_zero(mat
->row
[row
]+row
, mat
->n_col
-row
);
975 isl_assert(mat
->ctx
, pivot
>= 0, goto error
);
979 exchange(mat
, &inv
, NULL
, row
, pivot
, row
);
980 if (isl_int_is_neg(mat
->row
[row
][row
]))
981 oppose(mat
, &inv
, NULL
, row
, row
);
983 while ((off
= isl_seq_first_non_zero(mat
->row
[row
]+first
,
984 mat
->n_col
-first
)) != -1) {
986 isl_int_fdiv_q(a
, mat
->row
[row
][first
],
988 subtract(mat
, &inv
, NULL
, row
, row
, first
, a
);
989 if (!isl_int_is_zero(mat
->row
[row
][first
]))
990 exchange(mat
, &inv
, NULL
, row
, row
, first
);
994 for (i
= 0; i
< row
; ++i
) {
995 if (isl_int_is_zero(mat
->row
[row
][i
]))
997 isl_int_gcd(a
, mat
->row
[row
][row
], mat
->row
[row
][i
]);
998 isl_int_divexact(b
, mat
->row
[row
][i
], a
);
999 isl_int_divexact(a
, mat
->row
[row
][row
], a
);
1001 isl_mat_col_combine(mat
, i
, a
, i
, b
, row
);
1002 isl_mat_col_combine(inv
, i
, a
, i
, b
, row
);
1007 isl_int_set(a
, mat
->row
[0][0]);
1008 for (row
= 1; row
< mat
->n_row
; ++row
)
1009 isl_int_lcm(a
, a
, mat
->row
[row
][row
]);
1010 if (isl_int_is_zero(a
)){
1014 for (row
= 0; row
< mat
->n_row
; ++row
) {
1015 isl_int_divexact(mat
->row
[row
][row
], a
, mat
->row
[row
][row
]);
1016 if (isl_int_is_one(mat
->row
[row
][row
]))
1018 isl_mat_col_scale(inv
, row
, mat
->row
[row
][row
]);
1031 __isl_give isl_mat
*isl_mat_transpose(__isl_take isl_mat
*mat
)
1033 struct isl_mat
*transpose
= NULL
;
1039 if (mat
->n_col
== mat
->n_row
) {
1040 mat
= isl_mat_cow(mat
);
1043 for (i
= 0; i
< mat
->n_row
; ++i
)
1044 for (j
= i
+ 1; j
< mat
->n_col
; ++j
)
1045 isl_int_swap(mat
->row
[i
][j
], mat
->row
[j
][i
]);
1048 transpose
= isl_mat_alloc(mat
->ctx
, mat
->n_col
, mat
->n_row
);
1051 for (i
= 0; i
< mat
->n_row
; ++i
)
1052 for (j
= 0; j
< mat
->n_col
; ++j
)
1053 isl_int_set(transpose
->row
[j
][i
], mat
->row
[i
][j
]);
1061 __isl_give isl_mat
*isl_mat_swap_cols(__isl_take isl_mat
*mat
,
1062 unsigned i
, unsigned j
)
1066 mat
= isl_mat_cow(mat
);
1069 isl_assert(mat
->ctx
, i
< mat
->n_col
, goto error
);
1070 isl_assert(mat
->ctx
, j
< mat
->n_col
, goto error
);
1072 for (r
= 0; r
< mat
->n_row
; ++r
)
1073 isl_int_swap(mat
->row
[r
][i
], mat
->row
[r
][j
]);
1080 __isl_give isl_mat
*isl_mat_swap_rows(__isl_take isl_mat
*mat
,
1081 unsigned i
, unsigned j
)
1087 mat
= isl_mat_cow(mat
);
1091 mat
->row
[i
] = mat
->row
[j
];
1096 /* Calculate the product of two matrices.
1098 * This function is optimized for operand matrices that contain many zeros and
1099 * skips multiplications where we know one of the operands is zero.
1101 __isl_give isl_mat
*isl_mat_product(__isl_take isl_mat
*left
,
1102 __isl_take isl_mat
*right
)
1105 struct isl_mat
*prod
;
1107 if (!left
|| !right
)
1109 isl_assert(left
->ctx
, left
->n_col
== right
->n_row
, goto error
);
1110 prod
= isl_mat_alloc(left
->ctx
, left
->n_row
, right
->n_col
);
1113 if (left
->n_col
== 0) {
1114 for (i
= 0; i
< prod
->n_row
; ++i
)
1115 isl_seq_clr(prod
->row
[i
], prod
->n_col
);
1117 isl_mat_free(right
);
1120 for (i
= 0; i
< prod
->n_row
; ++i
) {
1121 for (j
= 0; j
< prod
->n_col
; ++j
)
1122 isl_int_mul(prod
->row
[i
][j
],
1123 left
->row
[i
][0], right
->row
[0][j
]);
1124 for (k
= 1; k
< left
->n_col
; ++k
) {
1125 if (isl_int_is_zero(left
->row
[i
][k
]))
1127 for (j
= 0; j
< prod
->n_col
; ++j
)
1128 isl_int_addmul(prod
->row
[i
][j
],
1129 left
->row
[i
][k
], right
->row
[k
][j
]);
1133 isl_mat_free(right
);
1137 isl_mat_free(right
);
1141 /* Replace the variables x in the rows q by x' given by x = M x',
1142 * with M the matrix mat.
1144 * If the number of new variables is greater than the original
1145 * number of variables, then the rows q have already been
1146 * preextended. If the new number is smaller, then the coefficients
1147 * of the divs, which are not changed, need to be shifted down.
1148 * The row q may be the equalities, the inequalities or the
1149 * div expressions. In the latter case, has_div is true and
1150 * we need to take into account the extra denominator column.
1152 static int preimage(struct isl_ctx
*ctx
, isl_int
**q
, unsigned n
,
1153 unsigned n_div
, int has_div
, struct isl_mat
*mat
)
1159 if (mat
->n_col
>= mat
->n_row
)
1162 e
= mat
->n_row
- mat
->n_col
;
1164 for (i
= 0; i
< n
; ++i
)
1165 isl_int_mul(q
[i
][0], q
[i
][0], mat
->row
[0][0]);
1166 t
= isl_mat_sub_alloc6(mat
->ctx
, q
, 0, n
, has_div
, mat
->n_row
);
1167 t
= isl_mat_product(t
, mat
);
1170 for (i
= 0; i
< n
; ++i
) {
1171 isl_seq_swp_or_cpy(q
[i
] + has_div
, t
->row
[i
], t
->n_col
);
1172 isl_seq_cpy(q
[i
] + has_div
+ t
->n_col
,
1173 q
[i
] + has_div
+ t
->n_col
+ e
, n_div
);
1174 isl_seq_clr(q
[i
] + has_div
+ t
->n_col
+ n_div
, e
);
1180 /* Replace the variables x in bset by x' given by x = M x', with
1183 * If there are fewer variables x' then there are x, then we perform
1184 * the transformation in place, which means that, in principle,
1185 * this frees up some extra variables as the number
1186 * of columns remains constant, but we would have to extend
1187 * the div array too as the number of rows in this array is assumed
1188 * to be equal to extra.
1190 struct isl_basic_set
*isl_basic_set_preimage(struct isl_basic_set
*bset
,
1191 struct isl_mat
*mat
)
1193 struct isl_ctx
*ctx
;
1199 bset
= isl_basic_set_cow(bset
);
1203 isl_assert(ctx
, bset
->dim
->nparam
== 0, goto error
);
1204 isl_assert(ctx
, 1+bset
->dim
->n_out
== mat
->n_row
, goto error
);
1205 isl_assert(ctx
, mat
->n_col
> 0, goto error
);
1207 if (mat
->n_col
> mat
->n_row
) {
1208 bset
= isl_basic_set_extend(bset
, 0, mat
->n_col
-1, 0, 0, 0);
1211 } else if (mat
->n_col
< mat
->n_row
) {
1212 bset
->dim
= isl_space_cow(bset
->dim
);
1215 bset
->dim
->n_out
-= mat
->n_row
- mat
->n_col
;
1218 if (preimage(ctx
, bset
->eq
, bset
->n_eq
, bset
->n_div
, 0,
1219 isl_mat_copy(mat
)) < 0)
1222 if (preimage(ctx
, bset
->ineq
, bset
->n_ineq
, bset
->n_div
, 0,
1223 isl_mat_copy(mat
)) < 0)
1226 if (preimage(ctx
, bset
->div
, bset
->n_div
, bset
->n_div
, 1, mat
) < 0)
1229 ISL_F_CLR(bset
, ISL_BASIC_SET_NO_IMPLICIT
);
1230 ISL_F_CLR(bset
, ISL_BASIC_SET_NO_REDUNDANT
);
1231 ISL_F_CLR(bset
, ISL_BASIC_SET_NORMALIZED
);
1232 ISL_F_CLR(bset
, ISL_BASIC_SET_NORMALIZED_DIVS
);
1233 ISL_F_CLR(bset
, ISL_BASIC_SET_ALL_EQUALITIES
);
1235 bset
= isl_basic_set_simplify(bset
);
1236 bset
= isl_basic_set_finalize(bset
);
1242 isl_basic_set_free(bset
);
1246 struct isl_set
*isl_set_preimage(struct isl_set
*set
, struct isl_mat
*mat
)
1250 set
= isl_set_cow(set
);
1254 for (i
= 0; i
< set
->n
; ++i
) {
1255 set
->p
[i
] = isl_basic_set_preimage(set
->p
[i
],
1260 if (mat
->n_col
!= mat
->n_row
) {
1261 set
->dim
= isl_space_cow(set
->dim
);
1264 set
->dim
->n_out
+= mat
->n_col
;
1265 set
->dim
->n_out
-= mat
->n_row
;
1268 ISL_F_CLR(set
, ISL_SET_NORMALIZED
);
1276 /* Replace the variables x starting at "first_col" in the rows "rows"
1277 * of some coefficient matrix by x' with x = M x' with M the matrix mat.
1278 * That is, replace the corresponding coefficients c by c M.
1280 isl_stat
isl_mat_sub_transform(isl_int
**row
, unsigned n_row
,
1281 unsigned first_col
, __isl_take isl_mat
*mat
)
1288 return isl_stat_error
;
1289 ctx
= isl_mat_get_ctx(mat
);
1290 t
= isl_mat_sub_alloc6(ctx
, row
, 0, n_row
, first_col
, mat
->n_row
);
1291 t
= isl_mat_product(t
, mat
);
1293 return isl_stat_error
;
1294 for (i
= 0; i
< n_row
; ++i
)
1295 isl_seq_swp_or_cpy(row
[i
] + first_col
, t
->row
[i
], t
->n_col
);
1300 void isl_mat_print_internal(__isl_keep isl_mat
*mat
, FILE *out
, int indent
)
1305 fprintf(out
, "%*snull mat\n", indent
, "");
1309 if (mat
->n_row
== 0)
1310 fprintf(out
, "%*s[]\n", indent
, "");
1312 for (i
= 0; i
< mat
->n_row
; ++i
) {
1314 fprintf(out
, "%*s[[", indent
, "");
1316 fprintf(out
, "%*s[", indent
+1, "");
1317 for (j
= 0; j
< mat
->n_col
; ++j
) {
1320 isl_int_print(out
, mat
->row
[i
][j
], 0);
1322 if (i
== mat
->n_row
-1)
1323 fprintf(out
, "]]\n");
1325 fprintf(out
, "]\n");
1329 void isl_mat_dump(__isl_keep isl_mat
*mat
)
1331 isl_mat_print_internal(mat
, stderr
, 0);
1334 struct isl_mat
*isl_mat_drop_cols(struct isl_mat
*mat
, unsigned col
, unsigned n
)
1341 mat
= isl_mat_cow(mat
);
1345 if (col
!= mat
->n_col
-n
) {
1346 for (r
= 0; r
< mat
->n_row
; ++r
)
1347 isl_seq_cpy(mat
->row
[r
]+col
, mat
->row
[r
]+col
+n
,
1348 mat
->n_col
- col
- n
);
1354 struct isl_mat
*isl_mat_drop_rows(struct isl_mat
*mat
, unsigned row
, unsigned n
)
1358 mat
= isl_mat_cow(mat
);
1362 for (r
= row
; r
+n
< mat
->n_row
; ++r
)
1363 mat
->row
[r
] = mat
->row
[r
+n
];
1369 __isl_give isl_mat
*isl_mat_insert_cols(__isl_take isl_mat
*mat
,
1370 unsigned col
, unsigned n
)
1379 ext
= isl_mat_alloc(mat
->ctx
, mat
->n_row
, mat
->n_col
+ n
);
1383 isl_mat_sub_copy(mat
->ctx
, ext
->row
, mat
->row
, mat
->n_row
, 0, 0, col
);
1384 isl_mat_sub_copy(mat
->ctx
, ext
->row
, mat
->row
, mat
->n_row
,
1385 col
+ n
, col
, mat
->n_col
- col
);
1394 __isl_give isl_mat
*isl_mat_insert_zero_cols(__isl_take isl_mat
*mat
,
1395 unsigned first
, unsigned n
)
1401 mat
= isl_mat_insert_cols(mat
, first
, n
);
1405 for (i
= 0; i
< mat
->n_row
; ++i
)
1406 isl_seq_clr(mat
->row
[i
] + first
, n
);
1411 __isl_give isl_mat
*isl_mat_add_zero_cols(__isl_take isl_mat
*mat
, unsigned n
)
1416 return isl_mat_insert_zero_cols(mat
, mat
->n_col
, n
);
1419 __isl_give isl_mat
*isl_mat_insert_rows(__isl_take isl_mat
*mat
,
1420 unsigned row
, unsigned n
)
1429 ext
= isl_mat_alloc(mat
->ctx
, mat
->n_row
+ n
, mat
->n_col
);
1433 isl_mat_sub_copy(mat
->ctx
, ext
->row
, mat
->row
, row
, 0, 0, mat
->n_col
);
1434 isl_mat_sub_copy(mat
->ctx
, ext
->row
+ row
+ n
, mat
->row
+ row
,
1435 mat
->n_row
- row
, 0, 0, mat
->n_col
);
1444 __isl_give isl_mat
*isl_mat_add_rows(__isl_take isl_mat
*mat
, unsigned n
)
1449 return isl_mat_insert_rows(mat
, mat
->n_row
, n
);
1452 __isl_give isl_mat
*isl_mat_insert_zero_rows(__isl_take isl_mat
*mat
,
1453 unsigned row
, unsigned n
)
1457 mat
= isl_mat_insert_rows(mat
, row
, n
);
1461 for (i
= 0; i
< n
; ++i
)
1462 isl_seq_clr(mat
->row
[row
+ i
], mat
->n_col
);
1467 __isl_give isl_mat
*isl_mat_add_zero_rows(__isl_take isl_mat
*mat
, unsigned n
)
1472 return isl_mat_insert_zero_rows(mat
, mat
->n_row
, n
);
1475 void isl_mat_col_submul(struct isl_mat
*mat
,
1476 int dst_col
, isl_int f
, int src_col
)
1480 for (i
= 0; i
< mat
->n_row
; ++i
)
1481 isl_int_submul(mat
->row
[i
][dst_col
], f
, mat
->row
[i
][src_col
]);
1484 void isl_mat_col_add(__isl_keep isl_mat
*mat
, int dst_col
, int src_col
)
1491 for (i
= 0; i
< mat
->n_row
; ++i
)
1492 isl_int_add(mat
->row
[i
][dst_col
],
1493 mat
->row
[i
][dst_col
], mat
->row
[i
][src_col
]);
1496 void isl_mat_col_mul(struct isl_mat
*mat
, int dst_col
, isl_int f
, int src_col
)
1500 for (i
= 0; i
< mat
->n_row
; ++i
)
1501 isl_int_mul(mat
->row
[i
][dst_col
], f
, mat
->row
[i
][src_col
]);
1504 /* Add "f" times column "src_col" to column "dst_col" of "mat" and
1505 * return the result.
1507 __isl_give isl_mat
*isl_mat_col_addmul(__isl_take isl_mat
*mat
, int dst_col
,
1508 isl_int f
, int src_col
)
1512 if (check_col(mat
, dst_col
) < 0 || check_col(mat
, src_col
) < 0)
1513 return isl_mat_free(mat
);
1515 for (i
= 0; i
< mat
->n_row
; ++i
) {
1516 if (isl_int_is_zero(mat
->row
[i
][src_col
]))
1518 mat
= isl_mat_cow(mat
);
1521 isl_int_addmul(mat
->row
[i
][dst_col
], f
, mat
->row
[i
][src_col
]);
1527 /* Negate column "col" of "mat" and return the result.
1529 __isl_give isl_mat
*isl_mat_col_neg(__isl_take isl_mat
*mat
, int col
)
1533 if (check_col(mat
, col
) < 0)
1534 return isl_mat_free(mat
);
1536 for (i
= 0; i
< mat
->n_row
; ++i
) {
1537 if (isl_int_is_zero(mat
->row
[i
][col
]))
1539 mat
= isl_mat_cow(mat
);
1542 isl_int_neg(mat
->row
[i
][col
], mat
->row
[i
][col
]);
1548 __isl_give isl_mat
*isl_mat_unimodular_complete(__isl_take isl_mat
*M
, int row
)
1551 struct isl_mat
*H
= NULL
, *Q
= NULL
;
1556 isl_assert(M
->ctx
, M
->n_row
== M
->n_col
, goto error
);
1558 H
= isl_mat_left_hermite(isl_mat_copy(M
), 0, NULL
, &Q
);
1559 M
->n_row
= M
->n_col
;
1562 for (r
= 0; r
< row
; ++r
)
1563 isl_assert(M
->ctx
, isl_int_is_one(H
->row
[r
][r
]), goto error
);
1564 for (r
= row
; r
< M
->n_row
; ++r
)
1565 isl_seq_cpy(M
->row
[r
], Q
->row
[r
], M
->n_col
);
1576 __isl_give isl_mat
*isl_mat_concat(__isl_take isl_mat
*top
,
1577 __isl_take isl_mat
*bot
)
1579 struct isl_mat
*mat
;
1584 isl_assert(top
->ctx
, top
->n_col
== bot
->n_col
, goto error
);
1585 if (top
->n_row
== 0) {
1589 if (bot
->n_row
== 0) {
1594 mat
= isl_mat_alloc(top
->ctx
, top
->n_row
+ bot
->n_row
, top
->n_col
);
1597 isl_mat_sub_copy(mat
->ctx
, mat
->row
, top
->row
, top
->n_row
,
1599 isl_mat_sub_copy(mat
->ctx
, mat
->row
+ top
->n_row
, bot
->row
, bot
->n_row
,
1610 isl_bool
isl_mat_is_equal(__isl_keep isl_mat
*mat1
, __isl_keep isl_mat
*mat2
)
1615 return isl_bool_error
;
1617 if (mat1
->n_row
!= mat2
->n_row
)
1618 return isl_bool_false
;
1620 if (mat1
->n_col
!= mat2
->n_col
)
1621 return isl_bool_false
;
1623 for (i
= 0; i
< mat1
->n_row
; ++i
)
1624 if (!isl_seq_eq(mat1
->row
[i
], mat2
->row
[i
], mat1
->n_col
))
1625 return isl_bool_false
;
1627 return isl_bool_true
;
1630 __isl_give isl_mat
*isl_mat_from_row_vec(__isl_take isl_vec
*vec
)
1632 struct isl_mat
*mat
;
1636 mat
= isl_mat_alloc(vec
->ctx
, 1, vec
->size
);
1640 isl_seq_cpy(mat
->row
[0], vec
->el
, vec
->size
);
1649 /* Return a copy of row "row" of "mat" as an isl_vec.
1651 __isl_give isl_vec
*isl_mat_get_row(__isl_keep isl_mat
*mat
, unsigned row
)
1657 if (row
>= mat
->n_row
)
1658 isl_die(mat
->ctx
, isl_error_invalid
, "row out of range",
1661 v
= isl_vec_alloc(isl_mat_get_ctx(mat
), mat
->n_col
);
1664 isl_seq_cpy(v
->el
, mat
->row
[row
], mat
->n_col
);
1669 __isl_give isl_mat
*isl_mat_vec_concat(__isl_take isl_mat
*top
,
1670 __isl_take isl_vec
*bot
)
1672 return isl_mat_concat(top
, isl_mat_from_row_vec(bot
));
1675 __isl_give isl_mat
*isl_mat_move_cols(__isl_take isl_mat
*mat
,
1676 unsigned dst_col
, unsigned src_col
, unsigned n
)
1682 if (n
== 0 || dst_col
== src_col
)
1685 res
= isl_mat_alloc(mat
->ctx
, mat
->n_row
, mat
->n_col
);
1689 if (dst_col
< src_col
) {
1690 isl_mat_sub_copy(res
->ctx
, res
->row
, mat
->row
, mat
->n_row
,
1692 isl_mat_sub_copy(res
->ctx
, res
->row
, mat
->row
, mat
->n_row
,
1693 dst_col
, src_col
, n
);
1694 isl_mat_sub_copy(res
->ctx
, res
->row
, mat
->row
, mat
->n_row
,
1695 dst_col
+ n
, dst_col
, src_col
- dst_col
);
1696 isl_mat_sub_copy(res
->ctx
, res
->row
, mat
->row
, mat
->n_row
,
1697 src_col
+ n
, src_col
+ n
,
1698 res
->n_col
- src_col
- n
);
1700 isl_mat_sub_copy(res
->ctx
, res
->row
, mat
->row
, mat
->n_row
,
1702 isl_mat_sub_copy(res
->ctx
, res
->row
, mat
->row
, mat
->n_row
,
1703 src_col
, src_col
+ n
, dst_col
- src_col
);
1704 isl_mat_sub_copy(res
->ctx
, res
->row
, mat
->row
, mat
->n_row
,
1705 dst_col
, src_col
, n
);
1706 isl_mat_sub_copy(res
->ctx
, res
->row
, mat
->row
, mat
->n_row
,
1707 dst_col
+ n
, dst_col
+ n
,
1708 res
->n_col
- dst_col
- n
);
1718 /* Return the gcd of the elements in row "row" of "mat" in *gcd.
1719 * Return isl_stat_ok on success and isl_stat_error on failure.
1721 isl_stat
isl_mat_row_gcd(__isl_keep isl_mat
*mat
, int row
, isl_int
*gcd
)
1724 return isl_stat_error
;
1726 if (row
< 0 || row
>= mat
->n_row
)
1727 isl_die(isl_mat_get_ctx(mat
), isl_error_invalid
,
1728 "row out of range", return isl_stat_error
);
1729 isl_seq_gcd(mat
->row
[row
], mat
->n_col
, gcd
);
1734 void isl_mat_gcd(__isl_keep isl_mat
*mat
, isl_int
*gcd
)
1739 isl_int_set_si(*gcd
, 0);
1744 for (i
= 0; i
< mat
->n_row
; ++i
) {
1745 isl_seq_gcd(mat
->row
[i
], mat
->n_col
, &g
);
1746 isl_int_gcd(*gcd
, *gcd
, g
);
1751 /* Return the result of scaling "mat" by a factor of "m".
1753 __isl_give isl_mat
*isl_mat_scale(__isl_take isl_mat
*mat
, isl_int m
)
1757 if (isl_int_is_one(m
))
1760 mat
= isl_mat_cow(mat
);
1764 for (i
= 0; i
< mat
->n_row
; ++i
)
1765 isl_seq_scale(mat
->row
[i
], mat
->row
[i
], m
, mat
->n_col
);
1770 __isl_give isl_mat
*isl_mat_scale_down(__isl_take isl_mat
*mat
, isl_int m
)
1774 if (isl_int_is_one(m
))
1777 mat
= isl_mat_cow(mat
);
1781 for (i
= 0; i
< mat
->n_row
; ++i
)
1782 isl_seq_scale_down(mat
->row
[i
], mat
->row
[i
], m
, mat
->n_col
);
1787 __isl_give isl_mat
*isl_mat_scale_down_row(__isl_take isl_mat
*mat
, int row
,
1790 if (isl_int_is_one(m
))
1793 mat
= isl_mat_cow(mat
);
1797 isl_seq_scale_down(mat
->row
[row
], mat
->row
[row
], m
, mat
->n_col
);
1802 __isl_give isl_mat
*isl_mat_normalize(__isl_take isl_mat
*mat
)
1810 isl_mat_gcd(mat
, &gcd
);
1811 mat
= isl_mat_scale_down(mat
, gcd
);
1817 __isl_give isl_mat
*isl_mat_normalize_row(__isl_take isl_mat
*mat
, int row
)
1819 mat
= isl_mat_cow(mat
);
1823 isl_seq_normalize(mat
->ctx
, mat
->row
[row
], mat
->n_col
);
1828 /* Number of initial non-zero columns.
1830 int isl_mat_initial_non_zero_cols(__isl_keep isl_mat
*mat
)
1837 for (i
= 0; i
< mat
->n_col
; ++i
)
1838 if (row_first_non_zero(mat
->row
, mat
->n_row
, i
) < 0)