isl_mat.c

   1 #include "isl_seq.h"
   2 #include "isl_mat.h"
   3 #include "isl_map_private.h"
   4
   5 struct isl_mat *isl_mat_alloc(struct isl_ctx *ctx,
   6         unsigned n_row, unsigned n_col)
   7 {
   8         int i;
   9         struct isl_mat *mat;
  10
  11         mat = isl_alloc_type(ctx, struct isl_mat);
  12         if (!mat)
  13                 return NULL;
  14
  15         mat->block.size = 0;
  16         mat->row = NULL;
  17         mat->block = isl_blk_alloc(ctx, n_row * n_col);
  18         if (!mat->block.data)
  19                 goto error;
  20         mat->row = isl_alloc_array(ctx, isl_int *, n_row);
  21         if (!mat->row)
  22                 goto error;
  23
  24         for (i = 0; i < n_row; ++i)
  25                 mat->row[i] = mat->block.data + i * n_col;
  26
  27         mat->ref = 1;
  28         mat->n_row = n_row;
  29         mat->n_col = n_col;
  30         mat->flags = 0;
  31
  32         return mat;
  33 error:
  34         isl_blk_free(ctx, mat->block);
  35         free(mat);
  36         return NULL;
  37 }
  38
  39 struct isl_mat *isl_mat_sub_alloc(struct isl_ctx *ctx, isl_int **row,
  40         unsigned first_row, unsigned n_row, unsigned first_col, unsigned n_col)
  41 {
  42         int i;
  43         struct isl_mat *mat;
  44
  45         mat = isl_alloc_type(ctx, struct isl_mat);
  46         if (!mat)
  47                 return NULL;
  48         mat->row = isl_alloc_array(ctx, isl_int *, n_row);
  49         if (!mat->row)
  50                 goto error;
  51         for (i = 0; i < n_row; ++i)
  52                 mat->row[i] = row[first_row+i] + first_col;
  53         mat->ref = 1;
  54         mat->n_row = n_row;
  55         mat->n_col = n_col;
  56         mat->block = isl_blk_empty();
  57         mat->flags = ISL_MAT_BORROWED;
  58         return mat;
  59 error:
  60         free(mat);
  61         return NULL;
  62 }
  63
  64 void isl_mat_sub_copy(struct isl_ctx *ctx, isl_int **dst, isl_int **src,
  65         unsigned n_row, unsigned dst_col, unsigned src_col, unsigned n_col)
  66 {
  67         int i;
  68
  69         for (i = 0; i < n_row; ++i)
  70                 isl_seq_cpy(dst[i]+dst_col, src[i]+src_col, n_col);
  71 }
  72
  73 void isl_mat_sub_neg(struct isl_ctx *ctx, isl_int **dst, isl_int **src,
  74         unsigned n_row, unsigned dst_col, unsigned src_col, unsigned n_col)
  75 {
  76         int i;
  77
  78         for (i = 0; i < n_row; ++i)
  79                 isl_seq_neg(dst[i]+dst_col, src[i]+src_col, n_col);
  80 }
  81
  82 struct isl_mat *isl_mat_copy(struct isl_ctx *ctx, struct isl_mat *mat)
  83 {
  84         if (!mat)
  85                 return NULL;
  86
  87         mat->ref++;
  88         return mat;
  89 }
  90
  91 struct isl_mat *isl_mat_dup(struct isl_ctx *ctx, struct isl_mat *mat)
  92 {
  93         int i;
  94         struct isl_mat *mat2;
  95
  96         if (!mat)
  97                 return NULL;
  98         mat2 = isl_mat_alloc(ctx, mat->n_row, mat->n_col);
  99         if (!mat2)
 100                 return NULL;
 101         for (i = 0; i < mat->n_row; ++i)
 102                 isl_seq_cpy(mat2->row[i], mat->row[i], mat->n_col);
 103         return mat2;
 104 }
 105
 106 struct isl_mat *isl_mat_cow(struct isl_ctx *ctx, struct isl_mat *mat)
 107 {
 108         struct isl_mat *mat2;
 109         if (!mat)
 110                 return NULL;
 111
 112         if (mat->ref == 1 && !F_ISSET(mat, ISL_MAT_BORROWED))
 113                 return mat;
 114
 115         mat2 = isl_mat_dup(ctx, mat);
 116         isl_mat_free(ctx, mat);
 117         return mat2;
 118 }
 119
 120 void isl_mat_free(struct isl_ctx *ctx, struct isl_mat *mat)
 121 {
 122         if (!mat)
 123                 return;
 124
 125         if (--mat->ref > 0)
 126                 return;
 127
 128         if (!F_ISSET(mat, ISL_MAT_BORROWED))
 129                 isl_blk_free(ctx, mat->block);
 130         free(mat->row);
 131         free(mat);
 132 }
 133
 134 struct isl_mat *isl_mat_identity(struct isl_ctx *ctx, unsigned n_row)
 135 {
 136         int i;
 137         struct isl_mat *mat;
 138
 139         mat = isl_mat_alloc(ctx, n_row, n_row);
 140         if (!mat)
 141                 return NULL;
 142         for (i = 0; i < n_row; ++i) {
 143                 isl_seq_clr(mat->row[i], i);
 144                 isl_int_set_si(mat->row[i][i], 1);
 145                 isl_seq_clr(mat->row[i]+i+1, n_row-(i+1));
 146         }
 147
 148         return mat;
 149 }
 150
 151 struct isl_vec *isl_mat_vec_product(struct isl_ctx *ctx,
 152         struct isl_mat *mat, struct isl_vec *vec)
 153 {
 154         int i;
 155         struct isl_vec *prod;
 156
 157         if (!mat || !vec)
 158                 goto error;
 159
 160         isl_assert(ctx, mat->n_col == vec->size, goto error);
 161
 162         prod = isl_vec_alloc(ctx, mat->n_row);
 163         if (!prod)
 164                 goto error;
 165
 166         for (i = 0; i < prod->size; ++i)
 167                 isl_seq_inner_product(mat->row[i], vec->block.data, vec->size,
 168                                         &prod->block.data[i]);
 169         isl_mat_free(ctx, mat);
 170         isl_vec_free(ctx, vec);
 171         return prod;
 172 error:
 173         isl_mat_free(ctx, mat);
 174         isl_vec_free(ctx, vec);
 175         return NULL;
 176 }
 177
 178 struct isl_mat *isl_mat_aff_direct_sum(struct isl_ctx *ctx,
 179         struct isl_mat *left, struct isl_mat *right)
 180 {
 181         int i;
 182         struct isl_mat *sum;
 183
 184         if (!left || !right)
 185                 goto error;
 186
 187         isl_assert(ctx, left->n_row == right->n_row, goto error);
 188         isl_assert(ctx, left->n_row >= 1, goto error);
 189         isl_assert(ctx, left->n_col >= 1, goto error);
 190         isl_assert(ctx, right->n_col >= 1, goto error);
 191         isl_assert(ctx,
 192             isl_seq_first_non_zero(left->row[0]+1, left->n_col-1) == -1,
 193             goto error);
 194         isl_assert(ctx,
 195             isl_seq_first_non_zero(right->row[0]+1, right->n_col-1) == -1,
 196             goto error);
 197
 198         sum = isl_mat_alloc(ctx, left->n_row, left->n_col + right->n_col - 1);
 199         if (!sum)
 200                 goto error;
 201         isl_int_lcm(sum->row[0][0], left->row[0][0], right->row[0][0]);
 202         isl_int_divexact(left->row[0][0], sum->row[0][0], left->row[0][0]);
 203         isl_int_divexact(right->row[0][0], sum->row[0][0], right->row[0][0]);
 204
 205         isl_seq_clr(sum->row[0]+1, sum->n_col-1);
 206         for (i = 1; i < sum->n_row; ++i) {
 207                 isl_int_mul(sum->row[i][0], left->row[0][0], left->row[i][0]);
 208                 isl_int_addmul(sum->row[i][0],
 209                                 right->row[0][0], right->row[i][0]);
 210                 isl_seq_scale(sum->row[i]+1, left->row[i]+1, left->row[0][0],
 211                                 left->n_col-1);
 212                 isl_seq_scale(sum->row[i]+left->n_col,
 213                                 right->row[i]+1, right->row[0][0],
 214                                 right->n_col-1);
 215         }
 216
 217         isl_int_divexact(left->row[0][0], sum->row[0][0], left->row[0][0]);
 218         isl_int_divexact(right->row[0][0], sum->row[0][0], right->row[0][0]);
 219         isl_mat_free(ctx, left);
 220         isl_mat_free(ctx, right);
 221         return sum;
 222 error:
 223         isl_mat_free(ctx, left);
 224         isl_mat_free(ctx, right);
 225         return NULL;
 226 }
 227
 228 static void exchange(struct isl_ctx *ctx, struct isl_mat *M, struct isl_mat **U,
 229         struct isl_mat **Q, unsigned row, unsigned i, unsigned j)
 230 {
 231         int r;
 232         for (r = row; r < M->n_row; ++r)
 233                 isl_int_swap(M->row[r][i], M->row[r][j]);
 234         if (U) {
 235                 for (r = 0; r < (*U)->n_row; ++r)
 236                         isl_int_swap((*U)->row[r][i], (*U)->row[r][j]);
 237         }
 238         if (Q)
 239                 isl_mat_swap_rows(ctx, *Q, i, j);
 240 }
 241
 242 static void subtract(struct isl_ctx *ctx, struct isl_mat *M, struct isl_mat **U,
 243         struct isl_mat **Q, unsigned row, unsigned i, unsigned j, isl_int m)
 244 {
 245         int r, c;
 246         for (r = row; r < M->n_row; ++r)
 247                 isl_int_submul(M->row[r][j], m, M->row[r][i]);
 248         if (U) {
 249                 for (r = 0; r < (*U)->n_row; ++r)
 250                         isl_int_submul((*U)->row[r][j], m, (*U)->row[r][i]);
 251         }
 252         if (Q) {
 253                 for (r = 0; r < (*Q)->n_col; ++r)
 254                         isl_int_addmul((*Q)->row[i][r], m, (*Q)->row[j][r]);
 255         }
 256 }
 257
 258 static void oppose(struct isl_ctx *ctx, struct isl_mat *M, struct isl_mat **U,
 259         struct isl_mat **Q, unsigned row, unsigned col)
 260 {
 261         int r;
 262         for (r = row; r < M->n_row; ++r)
 263                 isl_int_neg(M->row[r][col], M->row[r][col]);
 264         if (U) {
 265                 for (r = 0; r < (*U)->n_row; ++r)
 266                         isl_int_neg((*U)->row[r][col], (*U)->row[r][col]);
 267         }
 268         if (Q)
 269                 isl_seq_neg((*Q)->row[col], (*Q)->row[col], (*Q)->n_col);
 270 }
 271
 272 /* Given matrix M, compute
 273  *
 274  *              M U = H
 275  *              M   = H Q
 276  *
 277  * with U and Q unimodular matrices and H a matrix in column echelon form
 278  * such that on each echelon row the entries in the non-echelon column
 279  * are non-negative and stricly smaller than the entries in the echelon
 280  * column.  If U or Q are NULL, then these matrices are not computed.
 281  */
 282 struct isl_mat *isl_mat_left_hermite(struct isl_ctx *ctx,
 283         struct isl_mat *M, struct isl_mat **U, struct isl_mat **Q)
 284 {
 285         isl_int c;
 286         int row, col;
 287
 288         if (U)
 289                 *U = NULL;
 290         if (Q)
 291                 *Q = NULL;
 292         if (!M)
 293                 goto error;
 294         M = isl_mat_cow(ctx, M);
 295         if (!M)
 296                 goto error;
 297         if (U) {
 298                 *U = isl_mat_identity(ctx, M->n_col);
 299                 if (!*U)
 300                         goto error;
 301         }
 302         if (Q) {
 303                 *Q = isl_mat_identity(ctx, M->n_col);
 304                 if (!*Q)
 305                         goto error;
 306         }
 307
 308         col = 0;
 309         isl_int_init(c);
 310         for (row = 0; row < M->n_row; ++row) {
 311                 int first, i, off;
 312                 first = isl_seq_abs_min_non_zero(M->row[row]+col, M->n_col-col);
 313                 if (first == -1)
 314                         continue;
 315                 first += col;
 316                 if (first != col)
 317                         exchange(ctx, M, U, Q, row, first, col);
 318                 if (isl_int_is_neg(M->row[row][col]))
 319                         oppose(ctx, M, U, Q, row, col);
 320                 first = col+1;
 321                 while ((off = isl_seq_first_non_zero(M->row[row]+first,
 322                                                        M->n_col-first)) != -1) {
 323                         first += off;
 324                         isl_int_fdiv_q(c, M->row[row][first], M->row[row][col]);
 325                         subtract(ctx, M, U, Q, row, col, first, c);
 326                         if (!isl_int_is_zero(M->row[row][first]))
 327                                 exchange(ctx, M, U, Q, row, first, col);
 328                         else
 329                                 ++first;
 330                 }
 331                 for (i = 0; i < col; ++i) {
 332                         if (isl_int_is_zero(M->row[row][i]))
 333                                 continue;
 334                         isl_int_fdiv_q(c, M->row[row][i], M->row[row][col]);
 335                         if (isl_int_is_zero(c))
 336                                 continue;
 337                         subtract(ctx, M, U, Q, row, col, i, c);
 338                 }
 339                 ++col;
 340         }
 341         isl_int_clear(c);
 342
 343         return M;
 344 error:
 345         if (Q) {
 346                 isl_mat_free(ctx, *Q);
 347                 *Q = NULL;
 348         }
 349         if (U) {
 350                 isl_mat_free(ctx, *U);
 351                 *U = NULL;
 352         }
 353         return NULL;
 354 }
 355
 356 struct isl_mat *isl_mat_lin_to_aff(struct isl_ctx *ctx, struct isl_mat *mat)
 357 {
 358         int i;
 359         struct isl_mat *mat2;
 360
 361         if (!mat)
 362                 return NULL;
 363         mat2 = isl_mat_alloc(ctx, 1+mat->n_row, 1+mat->n_col);
 364         if (!mat2)
 365                 return NULL;
 366         isl_int_set_si(mat2->row[0][0], 1);
 367         isl_seq_clr(mat2->row[0]+1, mat->n_col);
 368         for (i = 0; i < mat->n_row; ++i) {
 369                 isl_int_set_si(mat2->row[1+i][0], 0);
 370                 isl_seq_cpy(mat2->row[1+i]+1, mat->row[i], mat->n_col);
 371         }
 372         isl_mat_free(ctx, mat);
 373         return mat2;
 374 }
 375
 376 static int row_first_non_zero(isl_int **row, unsigned n_row, unsigned col)
 377 {
 378         int i;
 379
 380         for (i = 0; i < n_row; ++i)
 381                 if (!isl_int_is_zero(row[i][col]))
 382                         return i;
 383         return -1;
 384 }
 385
 386 static int row_abs_min_non_zero(isl_int **row, unsigned n_row, unsigned col)
 387 {
 388         int i, min = row_first_non_zero(row, n_row, col);
 389         if (min < 0)
 390                 return -1;
 391         for (i = min + 1; i < n_row; ++i) {
 392                 if (isl_int_is_zero(row[i][col]))
 393                         continue;
 394                 if (isl_int_abs_lt(row[i][col], row[min][col]))
 395                         min = i;
 396         }
 397         return min;
 398 }
 399
 400 static void inv_exchange(struct isl_ctx *ctx,
 401         struct isl_mat *left, struct isl_mat *right,
 402         unsigned i, unsigned j)
 403 {
 404         left = isl_mat_swap_rows(ctx, left, i, j);
 405         right = isl_mat_swap_rows(ctx, right, i, j);
 406 }
 407
 408 static void inv_oppose(struct isl_ctx *ctx,
 409         struct isl_mat *left, struct isl_mat *right, unsigned row)
 410 {
 411         isl_seq_neg(left->row[row]+row, left->row[row]+row, left->n_col-row);
 412         isl_seq_neg(right->row[row], right->row[row], right->n_col);
 413 }
 414
 415 static void inv_subtract(struct isl_ctx *ctx,
 416         struct isl_mat *left, struct isl_mat *right,
 417         unsigned row, unsigned i, isl_int m)
 418 {
 419         isl_int_neg(m, m);
 420         isl_seq_combine(left->row[i]+row,
 421                         ctx->one, left->row[i]+row,
 422                         m, left->row[row]+row,
 423                         left->n_col-row);
 424         isl_seq_combine(right->row[i], ctx->one, right->row[i],
 425                         m, right->row[row], right->n_col);
 426 }
 427
 428 /* Compute inv(left)*right
 429  */
 430 struct isl_mat *isl_mat_inverse_product(struct isl_ctx *ctx,
 431         struct isl_mat *left, struct isl_mat *right)
 432 {
 433         int row;
 434         isl_int a, b;
 435
 436         if (!left || !right)
 437                 goto error;
 438
 439         isl_assert(ctx, left->n_row == left->n_col, goto error);
 440         isl_assert(ctx, left->n_row == right->n_row, goto error);
 441
 442         if (left->n_row == 0) {
 443                 isl_mat_free(ctx, left);
 444                 return right;
 445         }
 446
 447         left = isl_mat_cow(ctx, left);
 448         right = isl_mat_cow(ctx, right);
 449         if (!left || !right)
 450                 goto error;
 451
 452         isl_int_init(a);
 453         isl_int_init(b);
 454         for (row = 0; row < left->n_row; ++row) {
 455                 int pivot, first, i, off;
 456                 pivot = row_abs_min_non_zero(left->row+row, left->n_row-row, row);
 457                 if (pivot < 0) {
 458                         isl_int_clear(a);
 459                         isl_int_clear(b);
 460                         isl_assert(ctx, pivot >= 0, goto error);
 461                 }
 462                 pivot += row;
 463                 if (pivot != row)
 464                         inv_exchange(ctx, left, right, pivot, row);
 465                 if (isl_int_is_neg(left->row[row][row]))
 466                         inv_oppose(ctx, left, right, row);
 467                 first = row+1;
 468                 while ((off = row_first_non_zero(left->row+first,
 469                                         left->n_row-first, row)) != -1) {
 470                         first += off;
 471                         isl_int_fdiv_q(a, left->row[first][row],
 472                                         left->row[row][row]);
 473                         inv_subtract(ctx, left, right, row, first, a);
 474                         if (!isl_int_is_zero(left->row[first][row]))
 475                                 inv_exchange(ctx, left, right, row, first);
 476                         else
 477                                 ++first;
 478                 }
 479                 for (i = 0; i < row; ++i) {
 480                         if (isl_int_is_zero(left->row[i][row]))
 481                                 continue;
 482                         isl_int_gcd(a, left->row[row][row], left->row[i][row]);
 483                         isl_int_divexact(b, left->row[i][row], a);
 484                         isl_int_divexact(a, left->row[row][row], a);
 485                         isl_int_neg(a, a);
 486                         isl_seq_combine(left->row[i]+row,
 487                                         a, left->row[i]+row,
 488                                         b, left->row[row]+row,
 489                                         left->n_col-row);
 490                         isl_seq_combine(right->row[i], a, right->row[i],
 491                                         b, right->row[row], right->n_col);
 492                 }
 493         }
 494         isl_int_clear(b);
 495
 496         isl_int_set(a, left->row[0][0]);
 497         for (row = 1; row < left->n_row; ++row)
 498                 isl_int_lcm(a, a, left->row[row][row]);
 499         if (isl_int_is_zero(a)){
 500                 isl_int_clear(a);
 501                 isl_assert(ctx, 0, goto error);
 502         }
 503         for (row = 0; row < left->n_row; ++row) {
 504                 isl_int_divexact(left->row[row][row], a, left->row[row][row]);
 505                 if (isl_int_is_one(left->row[row][row]))
 506                         continue;
 507                 isl_seq_scale(right->row[row], right->row[row],
 508                                 left->row[row][row], right->n_col);
 509         }
 510         isl_int_clear(a);
 511
 512         isl_mat_free(ctx, left);
 513         return right;
 514 error:
 515         isl_mat_free(ctx, left);
 516         isl_mat_free(ctx, right);
 517         return NULL;
 518 }
 519
 520 struct isl_mat *isl_mat_swap_rows(struct isl_ctx *ctx,
 521         struct isl_mat *mat, unsigned i, unsigned j)
 522 {
 523         isl_int *t;
 524
 525         if (!mat)
 526                 return NULL;
 527         mat = isl_mat_cow(ctx, mat);
 528         if (!mat)
 529                 return NULL;
 530         t = mat->row[i];
 531         mat->row[i] = mat->row[j];
 532         mat->row[j] = t;
 533         return mat;
 534 }
 535
 536 struct isl_mat *isl_mat_product(struct isl_ctx *ctx,
 537         struct isl_mat *left, struct isl_mat *right)
 538 {
 539         int i, j, k;
 540         struct isl_mat *prod;
 541
 542         if (!left || !right)
 543                 goto error;
 544         isl_assert(ctx, left->n_col == right->n_row, goto error);
 545         prod = isl_mat_alloc(ctx, left->n_row, right->n_col);
 546         if (!prod)
 547                 goto error;
 548         if (left->n_col == 0) {
 549                 for (i = 0; i < prod->n_row; ++i)
 550                         isl_seq_clr(prod->row[i], prod->n_col);
 551                 return prod;
 552         }
 553         for (i = 0; i < prod->n_row; ++i) {
 554                 for (j = 0; j < prod->n_col; ++j) {
 555                         isl_int_mul(prod->row[i][j],
 556                                     left->row[i][0], right->row[0][j]);
 557                         for (k = 1; k < left->n_col; ++k)
 558                                 isl_int_addmul(prod->row[i][j],
 559                                             left->row[i][k], right->row[k][j]);
 560                 }
 561         }
 562         isl_mat_free(ctx, left);
 563         isl_mat_free(ctx, right);
 564         return prod;
 565 error:
 566         isl_mat_free(ctx, left);
 567         isl_mat_free(ctx, right);
 568         return NULL;
 569 }
 570
 571 struct isl_basic_set *isl_basic_set_preimage(struct isl_ctx *ctx,
 572         struct isl_basic_set *bset, struct isl_mat *mat)
 573 {
 574         struct isl_mat *t;
 575         int i;
 576
 577         if (!bset || !mat)
 578                 goto error;
 579
 580         bset = isl_basic_set_cow(ctx, bset);
 581         if (!bset)
 582                 goto error;
 583
 584         isl_assert(ctx, bset->nparam == 0, goto error);
 585         isl_assert(ctx, bset->n_div == 0, goto error);
 586         isl_assert(ctx, 1+bset->dim == mat->n_row, goto error);
 587
 588         t = isl_mat_sub_alloc(ctx, bset->eq, 0, bset->n_eq, 0, mat->n_row);
 589         t = isl_mat_product(ctx, t, isl_mat_copy(ctx, mat));
 590         if (!t)
 591                 goto error;
 592         for (i = 0; i < bset->n_eq; ++i)
 593                 isl_seq_swp_or_cpy(bset->eq[i], t->row[i], t->n_col);
 594         isl_mat_free(ctx, t);
 595
 596         t = isl_mat_sub_alloc(ctx, bset->ineq, 0, bset->n_ineq, 0, mat->n_row);
 597         t = isl_mat_product(ctx, t, isl_mat_copy(ctx, mat));
 598         if (!t)
 599                 goto error;
 600         for (i = 0; i < bset->n_ineq; ++i)
 601                 isl_seq_swp_or_cpy(bset->ineq[i], t->row[i], t->n_col);
 602         isl_mat_free(ctx, t);
 603
 604         bset->dim -= mat->n_row - mat->n_col;
 605         bset = isl_basic_set_simplify(ctx, bset);
 606         bset = isl_basic_set_finalize(ctx, bset);
 607
 608         isl_mat_free(ctx, mat);
 609         return bset;
 610 error:
 611         isl_mat_free(ctx, mat);
 612         isl_basic_set_free(ctx, bset);
 613         return NULL;
 614 }
 615
 616 void isl_mat_dump(struct isl_ctx *ctx, struct isl_mat *mat,
 617                                 FILE *out, int indent)
 618 {
 619         int i, j;
 620
 621         if (mat->n_row == 0)
 622                 fprintf(out, "%*s[]\n", indent, "");
 623
 624         for (i = 0; i < mat->n_row; ++i) {
 625                 if (!i)
 626                         fprintf(out, "%*s[[", indent, "");
 627                 else
 628                         fprintf(out, "%*s[", indent+1, "");
 629                 for (j = 0; j < mat->n_col; ++j) {
 630                         if (j)
 631                             fprintf(out, ",");
 632                         isl_int_print(out, mat->row[i][j]);
 633                 }
 634                 if (i == mat->n_row-1)
 635                         fprintf(out, "]]\n");
 636                 else
 637                         fprintf(out, "]\n");
 638         }
 639 }
 640
 641 struct isl_mat *isl_mat_drop_col(struct isl_ctx *ctx, struct isl_mat *mat,
 642                                 unsigned col)
 643 {
 644         int r;
 645
 646         if (!mat)
 647                 return NULL;
 648
 649         if (col != mat->n_col-1) {
 650                 for (r = 0; r < mat->n_row; ++r)
 651                         isl_seq_cpy(mat->row[r]+col, mat->row[r]+col+1,
 652                                         mat->n_col - col - 1);
 653         }
 654         mat->n_col--;
 655         return mat;
 656 }
 657
 658 struct isl_mat *isl_mat_drop_rows(struct isl_ctx *ctx, struct isl_mat *mat,
 659                                 unsigned row, unsigned n)
 660 {
 661         int r;
 662
 663         if (!mat)
 664                 return NULL;
 665
 666         for (r = row; r+n < mat->n_row; ++r)
 667                 mat->row[r] = mat->row[r+n];
 668
 669         mat->n_row -= n;
 670         return mat;
 671 }