diff --git a/doc/source/fmpz_mat.rst b/doc/source/fmpz_mat.rst
index 050feb3c17..7226c487b7 100644
--- a/doc/source/fmpz_mat.rst
+++ b/doc/source/fmpz_mat.rst
@@ -1147,21 +1147,6 @@ Row reduction
``rank``, returns 1 if so and 0 otherwise.
-Modular gaussian elimination
---------------------------------------------------------------------------------
-
-
-.. function:: slong fmpz_mat_rref_mod(slong * perm, fmpz_mat_t A, const fmpz_t p)
-
- Uses fraction-free Gauss-Jordan elimination to set ``A``
- to its reduced row echelon form and returns the rank of ``A``.
- All computations are done modulo p.
-
- Pivot elements are chosen with ``fmpz_mat_find_pivot_any``.
- If ``perm`` is non-``NULL``, the permutation of
- rows in the matrix will also be applied to ``perm``.
-
-
Strong echelon form and Howell form
--------------------------------------------------------------------------------
diff --git a/doc/source/fmpz_mod_mat.rst b/doc/source/fmpz_mod_mat.rst
index b0acb29ec2..598a9c9787 100644
--- a/doc/source/fmpz_mod_mat.rst
+++ b/doc/source/fmpz_mod_mat.rst
@@ -262,13 +262,10 @@ Gaussian elimination
--------------------------------------------------------------------------------
-.. function:: slong fmpz_mod_mat_rref(slong * perm, fmpz_mod_mat_t mat, const fmpz_mod_ctx_t ctx)
+.. function:: slong fmpz_mod_mat_rref(fmpz_mod_mat_t res, const fmpz_mod_mat_t mat, const fmpz_mod_ctx_t ctx)
- Uses Gauss-Jordan elimination to set ``mat`` to its reduced row echelon
- form and returns the rank of ``mat``.
-
- If ``perm`` is non-``NULL``, the permutation of
- rows in the matrix will also be applied to ``perm``.
+ Sets ``res`` to the reduced row echelon form of ``mat``
+ and returns the rank.
The modulus is assumed to be prime.
diff --git a/doc/source/fq_default_mat.rst b/doc/source/fq_default_mat.rst
index 2ce33cb0eb..9fc7ea34fa 100644
--- a/doc/source/fq_default_mat.rst
+++ b/doc/source/fq_default_mat.rst
@@ -338,9 +338,9 @@ Reduced row echelon form
--------------------------------------------------------------------------------
-.. function:: slong fq_default_mat_rref(fq_default_mat_t A, const fq_default_ctx_t ctx)
+.. function:: slong fq_default_mat_rref(fq_default_mat_t B, const fq_default_mat_t A, const fq_default_ctx_t ctx)
- Puts `A` in reduced row echelon form and returns the rank of `A`.
+ Puts `B` in reduced row echelon form and returns the rank of `A`.
The rref is computed by first obtaining an unreduced row echelon
form via LU decomposition and then solving an additional
diff --git a/doc/source/fq_mat.rst b/doc/source/fq_mat.rst
index 67598dff8d..5ccce95d22 100644
--- a/doc/source/fq_mat.rst
+++ b/doc/source/fq_mat.rst
@@ -377,9 +377,9 @@ Reduced row echelon form
--------------------------------------------------------------------------------
-.. function:: slong fq_mat_rref(fq_mat_t A, const fq_ctx_t ctx)
+.. function:: slong fq_mat_rref(fq_mat_t B, const fq_mat_t A, const fq_ctx_t ctx)
- Puts `A` in reduced row echelon form and returns the rank of `A`.
+ Puts `B` in reduced row echelon form and returns the rank of `A`.
The rref is computed by first obtaining an unreduced row echelon
form via LU decomposition and then solving an additional
diff --git a/doc/source/fq_nmod_mat.rst b/doc/source/fq_nmod_mat.rst
index 322cf30757..81d9617663 100644
--- a/doc/source/fq_nmod_mat.rst
+++ b/doc/source/fq_nmod_mat.rst
@@ -376,9 +376,9 @@ Reduced row echelon form
--------------------------------------------------------------------------------
-.. function:: slong fq_nmod_mat_rref(fq_nmod_mat_t A, const fq_nmod_ctx_t ctx)
+.. function:: slong fq_nmod_mat_rref(fq_nmod_mat_t B, const fq_nmod_mat_t A, const fq_nmod_ctx_t ctx)
- Puts `A` in reduced row echelon form and returns the rank of `A`.
+ Puts `B` in reduced row echelon form and returns the rank of `A`.
The rref is computed by first obtaining an unreduced row echelon
form via LU decomposition and then solving an additional
diff --git a/doc/source/fq_zech_mat.rst b/doc/source/fq_zech_mat.rst
index 5d7363aeba..f8723d399e 100644
--- a/doc/source/fq_zech_mat.rst
+++ b/doc/source/fq_zech_mat.rst
@@ -343,9 +343,9 @@ Reduced row echelon form
--------------------------------------------------------------------------------
-.. function:: slong fq_zech_mat_rref(fq_zech_mat_t A, const fq_zech_ctx_t ctx)
+.. function:: slong fq_zech_mat_rref(fq_zech_mat_t B, const fq_zech_mat_t A, const fq_zech_ctx_t ctx)
- Puts `A` in reduced row echelon form and returns the rank of `A`.
+ Puts `B` in reduced row echelon form and returns the rank of `A`.
The rref is computed by first obtaining an unreduced row echelon
form via LU decomposition and then solving an additional
diff --git a/src/fmpz_mat.h b/src/fmpz_mat.h
index da6a30cd4b..b96dad3c6e 100644
--- a/src/fmpz_mat.h
+++ b/src/fmpz_mat.h
@@ -276,10 +276,6 @@ slong fmpz_mat_rref_mul(fmpz_mat_t B, fmpz_t den, const fmpz_mat_t A);
int fmpz_mat_is_in_rref_with_rank(const fmpz_mat_t A, const fmpz_t den,
slong rank);
-/* Modular gaussian elimination *********************************************/
-
-slong fmpz_mat_rref_mod(slong * perm, fmpz_mat_t A, const fmpz_t p);
-
/* Modular Howell and strong echelon form ***********************************/
slong fmpz_mat_howell_form_mod(fmpz_mat_t A, const fmpz_t mod);
diff --git a/src/fmpz_mat/rref_mod.c b/src/fmpz_mat/rref_mod.c
deleted file mode 100644
index 4ef23b3e30..0000000000
--- a/src/fmpz_mat/rref_mod.c
+++ /dev/null
@@ -1,85 +0,0 @@
-/*
- Copyright (C) 2011 Fredrik Johansson
- Copyright (C) 2012 Lina Kulakova
-
- This file is part of FLINT.
-
- FLINT is free software: you can redistribute it and/or modify it under
- the terms of the GNU Lesser General Public License (LGPL) as published
- by the Free Software Foundation; either version 2.1 of the License, or
- (at your option) any later version. See .
-*/
-
-#include "fmpz.h"
-#include "fmpz_mat.h"
-
-#define E(j,k) fmpz_mat_entry(A,j,k)
-
-slong
-fmpz_mat_rref_mod(slong * perm, fmpz_mat_t A, const fmpz_t p)
-{
- fmpz_t t, inv;
- slong m, n, j, k, rank, r, pivot_row, pivot_col;
-
- if (fmpz_mat_is_empty(A))
- {
- return 0;
- }
-
- m = A->r;
- n = A->c;
- rank = pivot_row = pivot_col = 0;
-
- fmpz_init(t);
- fmpz_init(inv);
-
- while (pivot_row < m && pivot_col < n)
- {
- r = fmpz_mat_find_pivot_any(A, pivot_row, m, pivot_col);
-
- if (r == -1)
- {
- pivot_col++;
- continue;
- }
- else if (r != pivot_row)
- {
- fmpz_mat_swap_rows(A, perm, pivot_row, r);
- }
- rank++;
-
- fmpz_invmod(inv, E(pivot_row, pivot_col), p);
-
- /* pivot row */
- for (k = pivot_col + 1; k < n; k++)
- {
- fmpz_mul(E(pivot_row, k), E(pivot_row, k), inv);
- fmpz_mod(E(pivot_row, k), E(pivot_row, k), p);
- }
- fmpz_one(E(pivot_row, pivot_col));
-
- /* other rows */
- for (j = 0; j < m; j++)
- {
- if (j == pivot_row)
- continue;
-
- for (k = pivot_col + 1; k < n; k++)
- {
- fmpz_mul(t, E(j, pivot_col), E(pivot_row, k));
- fmpz_sub(E(j, k), E(j, k), t);
- fmpz_mod(E(j, k), E(j, k), p);
- }
-
- fmpz_zero(E(j, pivot_col));
- }
-
- pivot_row++;
- pivot_col++;
- }
-
- fmpz_clear(inv);
- fmpz_clear(t);
-
- return rank;
-}
diff --git a/src/fmpz_mat/test/main.c b/src/fmpz_mat/test/main.c
index 5649bf7c48..fb265241ef 100644
--- a/src/fmpz_mat/test/main.c
+++ b/src/fmpz_mat/test/main.c
@@ -86,7 +86,6 @@
#include "t-rank.c"
#include "t-rref.c"
#include "t-rref_fflu.c"
-#include "t-rref_mod.c"
#include "t-rref_mul.c"
#include "t-scalar_addmul_fmpz.c"
#include "t-scalar_addmul_nmod_mat_fmpz.c"
@@ -186,7 +185,6 @@ test_struct tests[] =
TEST_FUNCTION(fmpz_mat_rank),
TEST_FUNCTION(fmpz_mat_rref),
TEST_FUNCTION(fmpz_mat_rref_fflu),
- TEST_FUNCTION(fmpz_mat_rref_mod),
TEST_FUNCTION(fmpz_mat_rref_mul),
TEST_FUNCTION(fmpz_mat_scalar_addmul_fmpz),
TEST_FUNCTION(fmpz_mat_scalar_addmul_nmod_mat_fmpz),
diff --git a/src/fmpz_mat/test/t-rref_mod.c b/src/fmpz_mat/test/t-rref_mod.c
deleted file mode 100644
index 10681a5244..0000000000
--- a/src/fmpz_mat/test/t-rref_mod.c
+++ /dev/null
@@ -1,148 +0,0 @@
-/*
- Copyright (C) 2010 Fredrik Johansson
- Copyright (C) 2012 Lina Kulakova
-
- This file is part of FLINT.
-
- FLINT is free software: you can redistribute it and/or modify it under
- the terms of the GNU Lesser General Public License (LGPL) as published
- by the Free Software Foundation; either version 2.1 of the License, or
- (at your option) any later version. See .
-*/
-
-#include "test_helpers.h"
-#include "ulong_extras.h"
-#include "fmpz.h"
-#include "fmpz_mat.h"
-
-static void
-check_rref(fmpz_mat_t A)
-{
- slong i, j, prev_pivot, prev_row_zero;
-
- prev_row_zero = 0;
- prev_pivot = -1;
-
- for (i = 0; i < A->r; i++)
- {
- for (j = 0; j < A->c; j++)
- {
- /* Found nonzero entry */
- if (!fmpz_is_zero(A->rows[i] + j))
- {
- if (prev_row_zero)
- {
- flint_printf("nonzero row after zero row\n");
- fflush(stdout);
- flint_abort();
- }
-
- if (j <= prev_pivot)
- {
- flint_printf("pivot not strictly to the right of previous\n");
- fflush(stdout);
- flint_abort();
- }
-
- prev_pivot = j;
- break;
- }
-
- prev_row_zero = (j + 1 == A->c);
- }
- }
-}
-
-TEST_FUNCTION_START(fmpz_mat_rref_mod, state)
-{
- fmpz_mat_t A;
- fmpz_t p;
- slong i, j, k, m, n, b, d, r, rank;
- slong *perm;
-
- /* Maximally sparse matrices of given rank */
- for (i = 0; i < 10000; i++)
- {
- m = n_randint(state, 10);
- n = n_randint(state, 10);
- perm = flint_malloc(FLINT_MAX(1, m) * sizeof(slong));
-
- fmpz_init_set_ui(p, n_randtest_prime(state, 0));
-
- for (r = 0; r <= FLINT_MIN(m, n); r++)
- {
- b = 1 + n_randint(state, 10) * n_randint(state, 10);
-
- fmpz_mat_init(A, m, n);
-
- fmpz_mat_randrank(A, state, r, b);
-
- for (j = 0; j < m; j++)
- for (k = 0; k < n; k++)
- fmpz_mod(fmpz_mat_entry(A, j, k), fmpz_mat_entry(A, j, k),
- p);
-
- rank = fmpz_mat_rref_mod(perm, A, p);
-
- if (r < rank)
- {
- flint_printf("FAIL:\n");
- flint_printf("wrong rank!\n");
- fflush(stdout);
- flint_abort();
- }
-
- check_rref(A);
-
- fmpz_mat_clear(A);
- }
-
- fmpz_clear(p);
- flint_free(perm);
- }
-
- /* Dense */
- for (i = 0; i < 10000; i++)
- {
- m = n_randint(state, 10);
- n = n_randint(state, 10);
- perm = flint_malloc(FLINT_MAX(1, m) * sizeof(slong));
-
- fmpz_init_set_ui(p, n_randtest_prime(state, 0));
-
- for (r = 0; r <= FLINT_MIN(m, n); r++)
- {
- b = 1 + n_randint(state, 10) * n_randint(state, 10);
- d = n_randint(state, 2 * m * n + 1);
-
- fmpz_mat_init(A, m, n);
- fmpz_mat_randrank(A, state, r, b);
-
- fmpz_mat_randops(A, state, d);
-
- for (j = 0; j < m; j++)
- for (k = 0; k < n; k++)
- fmpz_mod(fmpz_mat_entry(A, j, k), fmpz_mat_entry(A, j, k),
- p);
-
- rank = fmpz_mat_rref_mod(perm, A, p);
-
- if (r < rank)
- {
- flint_printf("FAIL:\n");
- flint_printf("wrong rank!\n");
- fflush(stdout);
- flint_abort();
- }
-
- check_rref(A);
-
- fmpz_mat_clear(A);
- }
-
- fmpz_clear(p);
- flint_free(perm);
- }
-
- TEST_FUNCTION_END(state);
-}
diff --git a/src/fmpz_mod_mat.h b/src/fmpz_mod_mat.h
index 7ed9e32479..22462295b3 100644
--- a/src/fmpz_mod_mat.h
+++ b/src/fmpz_mod_mat.h
@@ -91,10 +91,25 @@ int fmpz_mod_mat_is_square(const fmpz_mod_mat_t mat, const fmpz_mod_ctx_t ctx)
void fmpz_mod_mat_swap(fmpz_mod_mat_t mat1, fmpz_mod_mat_t mat2, const fmpz_mod_ctx_t ctx);
void fmpz_mod_mat_swap_entrywise(fmpz_mod_mat_t mat1, fmpz_mod_mat_t mat2, const fmpz_mod_ctx_t ctx);
-void _fmpz_mod_mat_reduce(fmpz_mod_mat_t mat, const fmpz_mod_ctx_t ctx);
-
void fmpz_mod_mat_set(fmpz_mod_mat_t B, const fmpz_mod_mat_t A, const fmpz_mod_ctx_t ctx);
+/* Conversions */
+
+FMPZ_MOD_MAT_INLINE
+void fmpz_mod_mat_set_nmod_mat(fmpz_mod_mat_t A, const nmod_mat_t B, const fmpz_mod_ctx_t ctx)
+{
+ fmpz_mat_set_nmod_mat_unsigned(A, B);
+}
+
+void fmpz_mod_mat_set_fmpz_mat(fmpz_mod_mat_t A, const fmpz_mat_t B, const fmpz_mod_ctx_t ctx);
+
+void fmpz_mod_mat_get_fmpz_mat(fmpz_mat_t A, const fmpz_mod_mat_t B, const fmpz_mod_ctx_t ctx);
+
+FMPZ_MOD_MAT_INLINE void _fmpz_mod_mat_reduce(fmpz_mod_mat_t mat, const fmpz_mod_ctx_t ctx)
+{
+ fmpz_mod_mat_set_fmpz_mat(mat, mat, ctx);
+}
+
/* Random matrix generation */
void fmpz_mod_mat_randtest(fmpz_mod_mat_t mat, flint_rand_t state, const fmpz_mod_ctx_t ctx);
@@ -158,18 +173,6 @@ void fmpz_mod_mat_transpose(fmpz_mod_mat_t B, const fmpz_mod_mat_t A, const fmpz
fmpz_mat_transpose(B, A);
}
-/* Conversions */
-
-FMPZ_MOD_MAT_INLINE
-void fmpz_mod_mat_set_nmod_mat(fmpz_mod_mat_t A, const nmod_mat_t B, const fmpz_mod_ctx_t ctx)
-{
- fmpz_mat_set_nmod_mat_unsigned(A, B);
-}
-
-void fmpz_mod_mat_set_fmpz_mat(fmpz_mod_mat_t A, const fmpz_mat_t B, const fmpz_mod_ctx_t ctx);
-
-void fmpz_mod_mat_get_fmpz_mat(fmpz_mat_t A, const fmpz_mod_mat_t B, const fmpz_mod_ctx_t ctx);
-
/* Addition and subtraction */
void fmpz_mod_mat_add(fmpz_mod_mat_t C, const fmpz_mod_mat_t A, const fmpz_mod_mat_t B, const fmpz_mod_ctx_t ctx);
@@ -220,7 +223,7 @@ void fmpz_mod_mat_trace(fmpz_t trace, const fmpz_mod_mat_t mat, const fmpz_mod_c
/* Gaussian elimination *********************************************/
-slong fmpz_mod_mat_rref(slong * perm, fmpz_mod_mat_t mat, const fmpz_mod_ctx_t ctx);
+slong fmpz_mod_mat_rref(fmpz_mod_mat_t res, const fmpz_mod_mat_t mat, const fmpz_mod_ctx_t ctx);
slong fmpz_mod_mat_reduce_row(fmpz_mod_mat_t A, slong * P, slong * L, slong m, const fmpz_mod_ctx_t ctx);
diff --git a/src/fmpz_mod_mat/can_solve.c b/src/fmpz_mod_mat/can_solve.c
index edb87a4318..1abcbf2f14 100644
--- a/src/fmpz_mod_mat/can_solve.c
+++ b/src/fmpz_mod_mat/can_solve.c
@@ -15,129 +15,15 @@
#include "fmpz.h"
#include "fmpz_mat.h"
#include "fmpz_mod_mat.h"
+#include "gr.h"
+#include "gr_mat.h"
int fmpz_mod_mat_can_solve(fmpz_mod_mat_t X, const fmpz_mod_mat_t A,
const fmpz_mod_mat_t B, const fmpz_mod_ctx_t ctx)
{
- slong i, j, k, col, *pivots, rank, *perm;
- fmpz_mod_mat_t LU, LU2, PB;
- int result = 1;
-
- if (A->r != B->r || A->c != X->r || X->c != B->c)
- {
- return 0;
- }
-
- if (A->r == 0 || B->c == 0)
- {
- fmpz_mod_mat_zero(X, ctx);
- return 1;
- }
-
- if (A->c == 0)
- {
- fmpz_mod_mat_zero(X, ctx);
- return fmpz_mod_mat_is_zero(B, ctx);
- }
-
- fmpz_mod_mat_init_set(LU, A, ctx);
- perm = flint_malloc(sizeof(slong) * A->r);
- for (i = 0; i < A->r; i++)
- perm[i] = i;
-
- rank = fmpz_mod_mat_lu(perm, LU, 0, ctx);
-
- fmpz_mod_mat_window_init(PB, B, 0, 0, B->r, B->c, ctx);
- for (i = 0; i < B->r; i++)
- PB->rows[i] = B->rows[perm[i]];
-
- fmpz_mod_mat_init(LU2, rank, rank, ctx);
-
- pivots = flint_malloc(sizeof(slong)*rank);
-
- col = 0;
- for (i = 0; i < rank; i++)
- {
- while (fmpz_is_zero(fmpz_mod_mat_entry(LU, i, col)))
- col++;
-
- pivots[i] = col;
-
- for (j = 0; j < rank; j++)
- fmpz_mod_mat_set_entry(LU2, j, i, fmpz_mod_mat_entry(LU, j, col), ctx);
-
- col++;
- }
-
- X->r = rank;
- PB->r = rank;
- LU->r = rank;
- fmpz_mod_mat_solve_tril(X, LU, PB, 1, ctx);
- LU->r = A->r;
-
- if (A->r > rank)
- {
- fmpz_mod_mat_t P;
-
- LU->rows += rank;
- LU->r = A->r - rank;
- X->r = LU->c;
-
- fmpz_mod_mat_init(P, LU->r, B->c, ctx);
-
- fmpz_mod_mat_mul(P, LU, X, ctx);
-
- PB->r = LU->r;
- PB->rows += rank;
-
- result = fmpz_mod_mat_equal(P, PB, ctx);
-
- PB->rows -= rank;
- fmpz_mod_mat_clear(P, ctx);
-
- LU->rows -= rank;
-
- if (!result)
- {
- X->r = A->c;
- fmpz_mod_mat_zero(X, ctx);
- goto cleanup;
- }
- }
-
- fmpz_mod_mat_solve_triu(X, LU2, X, 0, ctx);
-
- X->r = A->c;
-
- k = rank - 1;
- for (i = A->c - 1; i >= 0; i--)
- {
- if (k < 0 || i != pivots[k])
- {
- for (j = 0; j < B->c; j++)
- fmpz_zero(fmpz_mod_mat_entry(X, i, j));
- }
- else
- {
- for (j = 0; j < B->c; j++)
- fmpz_mod_mat_set_entry(X, i, j, fmpz_mod_mat_entry(X, k, j), ctx);
-
- k--;
- }
- }
-
-cleanup:
-
- fmpz_mod_mat_clear(LU2, ctx);
-
- PB->r = B->r;
- fmpz_mod_mat_window_clear(PB, ctx);
-
- LU->r = A->r;
- fmpz_mod_mat_clear(LU, ctx);
- flint_free(perm);
-
- flint_free(pivots);
-
- return result;
+ gr_ctx_t gr_ctx;
+ int status;
+ _gr_ctx_init_fmpz_mod_from_ref(gr_ctx, ctx);
+ status = gr_mat_solve_field((gr_mat_struct *) X, (const gr_mat_struct *) A, (const gr_mat_struct *) B, gr_ctx);
+ return (status == GR_SUCCESS) ? 1 : 0;
}
diff --git a/src/fmpz_mod_mat/nullspace.c b/src/fmpz_mod_mat/nullspace.c
index d2373f6cf0..aa0774037d 100644
--- a/src/fmpz_mod_mat/nullspace.c
+++ b/src/fmpz_mod_mat/nullspace.c
@@ -29,7 +29,7 @@ slong fmpz_mod_mat_nullspace(fmpz_mod_mat_t X, const fmpz_mod_mat_t A, const fmp
p = flint_malloc(sizeof(slong) * FLINT_MAX(m, n));
fmpz_mod_mat_init_set(tmp, A, ctx);
- rank = fmpz_mod_mat_rref(NULL, tmp, ctx);
+ rank = fmpz_mod_mat_rref(tmp, tmp, ctx);
nullity = n - rank;
fmpz_mod_mat_zero(X, ctx);
diff --git a/src/fmpz_mod_mat/reduce.c b/src/fmpz_mod_mat/reduce.c
deleted file mode 100644
index be4965c68f..0000000000
--- a/src/fmpz_mod_mat/reduce.c
+++ /dev/null
@@ -1,99 +0,0 @@
-/*
- Copyright (C) 2021 Daniel Schultz
-
- This file is part of FLINT.
-
- FLINT is free software: you can redistribute it and/or modify it under
- the terms of the GNU Lesser General Public License (LGPL) as published
- by the Free Software Foundation; either version 2.1 of the License, or
- (at your option) any later version. See .
-*/
-
-#include "thread_pool.h"
-#include "thread_support.h"
-#include "fmpz.h"
-#include "fmpz_vec.h"
-#include "fmpz_mod_vec.h"
-#include "fmpz_mod_mat.h"
-
-/* todo: rewrite using parallel_do */
-
-typedef struct {
- slong startrow;
- slong stoprow;
- fmpz_mod_mat_struct * M;
- const fmpz_mod_ctx_struct * ctx;
-} _worker_arg;
-
-static void _red_worker(void * varg)
-{
- _worker_arg * arg = (_worker_arg *) varg;
- slong startrow = arg->startrow;
- slong stoprow = arg->stoprow;
- fmpz_mod_mat_struct * M = arg->M;
- const fmpz_mod_ctx_struct * ctx = arg->ctx;
- slong c = fmpz_mod_mat_ncols(M, ctx);
- slong i;
-
- for (i = startrow; i < stoprow; i++)
- _fmpz_mod_vec_set_fmpz_vec(M->rows[i], M->rows[i], c, ctx);
-}
-
-void _fmpz_mod_mat_reduce(fmpz_mod_mat_t M, const fmpz_mod_ctx_t ctx)
-{
- slong i, r = fmpz_mod_mat_nrows(M, ctx);
- thread_pool_handle * handles;
- slong num_workers;
- _worker_arg mainarg;
- _worker_arg * args;
- slong limit;
-
- /* limit on threads */
- limit = fmpz_size(ctx->n) + r + fmpz_mod_mat_ncols(M, ctx);
- limit = limit < 64 ? 0 : (limit - 64)/64;
- limit = FLINT_MIN(limit, r);
-
- mainarg.startrow = 0;
- mainarg.stoprow = r;
- mainarg.M = M;
- mainarg.ctx = ctx;
-
- if (limit < 2)
- {
-use_one_thread:
- _red_worker(&mainarg);
- return;
- }
-
- num_workers = flint_request_threads(&handles, limit);
- if (num_workers < 1)
- {
- flint_give_back_threads(handles, num_workers);
- goto use_one_thread;
- }
-
- args = FLINT_ARRAY_ALLOC(num_workers, _worker_arg);
-
- for (i = 0; i < num_workers; i++)
- {
- args[i].startrow = (i + 0)*r/(num_workers + 1);
- args[i].stoprow = (i + 1)*r/(num_workers + 1);
- args[i].M = M;
- args[i].ctx = ctx;
- }
-
- i = num_workers;
- mainarg.startrow = (i + 0)*r/(num_workers + 1);
- mainarg.stoprow = (i + 1)*r/(num_workers + 1);
-
- for (i = 0; i < num_workers; i++)
- thread_pool_wake(global_thread_pool, handles[i], 0, _red_worker, &args[i]);
- _red_worker(&mainarg);
- for (i = 0; i < num_workers; i++)
- thread_pool_wait(global_thread_pool, handles[i]);
-
- flint_give_back_threads(handles, num_workers);
- flint_free(args);
-
- return;
-}
diff --git a/src/fmpz_mod_mat/rref.c b/src/fmpz_mod_mat/rref.c
index 1d00fec18e..c64431926c 100644
--- a/src/fmpz_mod_mat/rref.c
+++ b/src/fmpz_mod_mat/rref.c
@@ -1,5 +1,6 @@
/*
Copyright (C) 2019 Tommy Hofmann
+ Copyright (C) 2024 Fredrik Johansson
This file is part of FLINT.
@@ -10,8 +11,14 @@
*/
#include "fmpz_mod_mat.h"
+#include "gr.h"
+#include "gr_mat.h"
-slong fmpz_mod_mat_rref(slong * perm, fmpz_mod_mat_t mat, const fmpz_mod_ctx_t ctx)
+slong fmpz_mod_mat_rref(fmpz_mod_mat_t B, const fmpz_mod_mat_t A, const fmpz_mod_ctx_t ctx)
{
- return fmpz_mat_rref_mod(perm, mat, ctx->n);
+ gr_ctx_t gr_ctx;
+ slong rank;
+ _gr_ctx_init_fmpz_mod_from_ref(gr_ctx, ctx);
+ GR_MUST_SUCCEED(gr_mat_rref_lu(&rank, (gr_mat_struct *) B, (const gr_mat_struct *) A, gr_ctx));
+ return rank;
}
diff --git a/src/fmpz_mod_mat/set_fmpz_mat.c b/src/fmpz_mod_mat/set_fmpz_mat.c
new file mode 100644
index 0000000000..98d4122ff7
--- /dev/null
+++ b/src/fmpz_mod_mat/set_fmpz_mat.c
@@ -0,0 +1,73 @@
+/*
+ Copyright (C) 2021 Daniel Schultz
+ Copyright (C) 2024 Fredrik Johansson
+
+ This file is part of FLINT.
+
+ FLINT is free software: you can redistribute it and/or modify it under
+ the terms of the GNU Lesser General Public License (LGPL) as published
+ by the Free Software Foundation; either version 2.1 of the License, or
+ (at your option) any later version. See .
+*/
+
+#include "thread_pool.h"
+#include "thread_support.h"
+#include "fmpz.h"
+#include "fmpz_vec.h"
+#include "fmpz_mod_vec.h"
+#include "fmpz_mod_mat.h"
+
+typedef struct {
+ fmpz ** Mrows;
+ fmpz * const * Arows;
+ slong c;
+ const fmpz_mod_ctx_struct * ctx;
+} _worker_arg;
+
+static void _red_worker(slong i, void * varg)
+{
+ _worker_arg * arg = (_worker_arg *) varg;
+ _fmpz_mod_vec_set_fmpz_vec(arg->Mrows[i], arg->Arows[i], arg->c, arg->ctx);
+}
+
+static void
+_fmpz_mod_mat_set_fmpz_mat(fmpz_mod_mat_t M, const fmpz_mat_t A, const fmpz_mod_ctx_t ctx)
+{
+ slong i, r, c;
+
+ r = fmpz_mod_mat_nrows(M, ctx);
+ c = fmpz_mod_mat_ncols(M, ctx);
+
+ for (i = 0; i < r; i++)
+ _fmpz_mod_vec_set_fmpz_vec(M->rows[i], A->rows[i], c, ctx);
+}
+
+void fmpz_mod_mat_set_fmpz_mat(fmpz_mod_mat_t M, const fmpz_mat_t A, const fmpz_mod_ctx_t ctx)
+{
+ slong r, c;
+ slong limit;
+
+ r = fmpz_mod_mat_nrows(M, ctx);
+ c = fmpz_mod_mat_ncols(M, ctx);
+
+ /* limit on threads */
+ limit = fmpz_size(ctx->n) + r + c;
+ limit = limit < 64 ? 0 : (limit - 64)/64;
+ limit = FLINT_MIN(limit, r);
+
+ if (limit < 2)
+ {
+ _fmpz_mod_mat_set_fmpz_mat(M, A, ctx);
+ }
+ else
+ {
+ _worker_arg arg;
+
+ arg.Mrows = M->rows;
+ arg.Arows = A->rows;
+ arg.c = c;
+ arg.ctx = ctx;
+
+ flint_parallel_do(_red_worker, &arg, r, limit, FLINT_PARALLEL_UNIFORM);
+ }
+}
diff --git a/src/fmpz_mod_mat/set_get.c b/src/fmpz_mod_mat/set_get.c
index ac457abe84..5e23f6a5fe 100644
--- a/src/fmpz_mod_mat/set_get.c
+++ b/src/fmpz_mod_mat/set_get.c
@@ -19,13 +19,6 @@ void fmpz_mod_mat_set_entry(fmpz_mod_mat_t mat, slong i, slong j, const fmpz_t v
fmpz_set(fmpz_mat_entry(mat, i, j), val);
}
-void fmpz_mod_mat_set_fmpz_mat(fmpz_mod_mat_t A, const fmpz_mat_t B, const fmpz_mod_ctx_t ctx)
-{
- /* todo: an implementation of this function should actually replace _fmpz_mod_mat_reduce */
- fmpz_mat_set(A, B);
- _fmpz_mod_mat_reduce(A, ctx);
-}
-
void fmpz_mod_mat_set(fmpz_mod_mat_t B, const fmpz_mod_mat_t A, const fmpz_mod_ctx_t ctx)
{
fmpz_mat_set(B, A);
diff --git a/src/fmpz_mod_mat/solve.c b/src/fmpz_mod_mat/solve.c
index e5612501c4..4fdcd6da66 100644
--- a/src/fmpz_mod_mat/solve.c
+++ b/src/fmpz_mod_mat/solve.c
@@ -1,6 +1,7 @@
/*
Copyright (C) 2018 Tommy Hofmann
Copyright (C) 2021 Daniel Schultz
+ Copyright (C) 2024 Fredrik Johansson
This file is part of FLINT.
@@ -11,45 +12,15 @@
*/
#include "fmpz_mod_mat.h"
+#include "gr.h"
+#include "gr_mat.h"
int fmpz_mod_mat_solve(fmpz_mod_mat_t X, const fmpz_mod_mat_t A,
const fmpz_mod_mat_t B, const fmpz_mod_ctx_t ctx)
{
- slong i, rank, *perm;
- fmpz_mod_mat_t LU;
- int result;
-
- if (fmpz_mod_mat_is_empty(A, ctx))
- return 1;
-
- fmpz_mod_mat_init_set(LU, A, ctx);
-
- perm = flint_malloc(sizeof(slong) * A->r);
- for (i = 0; i < A->r; i++)
- perm[i] = i;
-
- rank = fmpz_mod_mat_lu(perm, LU, 1, ctx);
-
- if (rank == A->r)
- {
- fmpz_mod_mat_t PB;
- fmpz_mod_mat_window_init(PB, B, 0, 0, B->r, B->c, ctx);
- for (i = 0; i < A->r; i++)
- PB->rows[i] = B->rows[perm[i]];
-
- fmpz_mod_mat_solve_tril(X, LU, PB, 1, ctx);
- fmpz_mod_mat_solve_triu(X, LU, X, 0, ctx);
-
- fmpz_mod_mat_window_clear(PB, ctx);
- result = 1;
- }
- else
- {
- result = 0;
- }
-
- fmpz_mod_mat_clear(LU, ctx);
- flint_free(perm);
-
- return result;
+ gr_ctx_t gr_ctx;
+ int status;
+ _gr_ctx_init_fmpz_mod_from_ref(gr_ctx, ctx);
+ status = gr_mat_nonsingular_solve_lu((gr_mat_struct *) X, (const gr_mat_struct *) A, (const gr_mat_struct *) B, gr_ctx);
+ return (status == GR_SUCCESS) ? 1 : 0;
}
diff --git a/src/fmpz_mod_mat/test/t-rref.c b/src/fmpz_mod_mat/test/t-rref.c
index 6420fd79b4..b6f46dd2f8 100644
--- a/src/fmpz_mod_mat/test/t-rref.c
+++ b/src/fmpz_mod_mat/test/t-rref.c
@@ -59,14 +59,12 @@ TEST_FUNCTION_START(fmpz_mod_mat_rref, state)
fmpz_mod_mat_t A;
fmpz_mod_ctx_t ctx;
slong i, m, n, d, r, rank;
- slong *perm;
/* Maximally sparse matrices of given rank */
for (i = 0; i < 10000; i++)
{
m = n_randint(state, 10);
n = n_randint(state, 10);
- perm = flint_malloc(FLINT_MAX(1, m) * sizeof(slong));
fmpz_mod_ctx_init_rand_bits_prime(ctx, state, 100);
@@ -75,7 +73,7 @@ TEST_FUNCTION_START(fmpz_mod_mat_rref, state)
fmpz_mod_mat_init(A, m, n, ctx);
fmpz_mod_mat_randrank(A, state, r, ctx);
- rank = fmpz_mod_mat_rref(perm, A, ctx);
+ rank = fmpz_mod_mat_rref(A, A, ctx);
if (r < rank)
{
@@ -92,7 +90,6 @@ TEST_FUNCTION_START(fmpz_mod_mat_rref, state)
}
fmpz_mod_ctx_clear(ctx);
- flint_free(perm);
}
/* Dense */
@@ -100,7 +97,6 @@ TEST_FUNCTION_START(fmpz_mod_mat_rref, state)
{
m = n_randint(state, 5);
n = n_randint(state, 4);
- perm = flint_malloc(FLINT_MAX(1, m) * sizeof(slong));
fmpz_mod_ctx_init_rand_bits_prime(ctx, state, 100);
@@ -112,7 +108,7 @@ TEST_FUNCTION_START(fmpz_mod_mat_rref, state)
fmpz_mod_mat_randrank(A, state, r, ctx);
fmpz_mod_mat_randops(A, state, d, ctx);
- rank = fmpz_mod_mat_rref(perm, A, ctx);
+ rank = fmpz_mod_mat_rref(A, A, ctx);
if (r < rank)
{
@@ -128,7 +124,6 @@ TEST_FUNCTION_START(fmpz_mod_mat_rref, state)
}
fmpz_mod_ctx_clear(ctx);
- flint_free(perm);
}
TEST_FUNCTION_END(state);
diff --git a/src/fmpz_mod_mpoly_factor/bpoly_factor_smprime.c b/src/fmpz_mod_mpoly_factor/bpoly_factor_smprime.c
index fd9a1995b9..e2f2de06f6 100644
--- a/src/fmpz_mod_mpoly_factor/bpoly_factor_smprime.c
+++ b/src/fmpz_mod_mpoly_factor/bpoly_factor_smprime.c
@@ -1024,7 +1024,7 @@ static void _lattice(
fmpz_mod_mat_init(T2, fmpz_mod_mat_nrows(T1, ctx), fmpz_mod_mat_ncols(N, ctx), ctx);
fmpz_mod_mat_mul(T2, T1, N, ctx);
fmpz_mod_mat_swap(T2, N, ctx);
- fmpz_mod_mat_rref(NULL, N, ctx);
+ fmpz_mod_mat_rref(N, N, ctx);
fmpz_mod_mat_clear(M, ctx);
fmpz_mod_mat_clear(T1, ctx);
diff --git a/src/fmpz_mod_mpoly_factor/fmpz_mod_mat_extras.c b/src/fmpz_mod_mpoly_factor/fmpz_mod_mat_extras.c
index 60c441718f..8ca494d95d 100644
--- a/src/fmpz_mod_mpoly_factor/fmpz_mod_mat_extras.c
+++ b/src/fmpz_mod_mpoly_factor/fmpz_mod_mat_extras.c
@@ -45,7 +45,7 @@ void fmpz_mod_mat_init_nullspace_tr(fmpz_mod_mat_t X, fmpz_mod_mat_t tmp, const
p = FLINT_ARRAY_ALLOC(FLINT_MAX(m, n), slong);
- rank = fmpz_mod_mat_rref(NULL, tmp, ctx);
+ rank = fmpz_mod_mat_rref(tmp, tmp, ctx);
nullity = n - rank;
diff --git a/src/fmpz_mod_mpoly_factor/gcd_zippel.c b/src/fmpz_mod_mpoly_factor/gcd_zippel.c
index e3fd79fd8c..4234477392 100644
--- a/src/fmpz_mod_mpoly_factor/gcd_zippel.c
+++ b/src/fmpz_mod_mpoly_factor/gcd_zippel.c
@@ -344,7 +344,7 @@ int fmpz_mod_mpolyl_gcds_zippel(
*/
- fmpz_mod_mat_rref(NULL, ML + s, ctx->ffinfo);
+ fmpz_mod_mat_rref(ML + s, ML + s, ctx->ffinfo);
for (i = 0; i < n; i++)
if (!fmpz_is_one(fmpz_mod_mat_entry(ML + s, i, i)))
diff --git a/src/fmpz_mod_poly_factor/factor_berlekamp.c b/src/fmpz_mod_poly_factor/factor_berlekamp.c
index 5bc893c712..677ce96047 100644
--- a/src/fmpz_mod_poly_factor/factor_berlekamp.c
+++ b/src/fmpz_mod_poly_factor/factor_berlekamp.c
@@ -19,6 +19,7 @@
#include "fmpz_mod.h"
#include "fmpz_mod_poly.h"
#include "fmpz_mod_poly_factor.h"
+#include "fmpz_mod_mat.h"
static void
fmpz_mod_poly_to_fmpz_mat_col(fmpz_mat_t mat, slong col, fmpz_mod_poly_t poly)
@@ -69,7 +70,7 @@ __fmpz_mod_poly_factor_berlekamp(fmpz_mod_poly_factor_t factors,
fmpz_mat_t matrix;
fmpz_t coeff, q, mul, pow;
slong i, nullity, col, row;
- slong *shift, *perm;
+ slong *shift;
fmpz_mod_poly_t *basis;
if (f->length <= 2)
@@ -132,9 +133,7 @@ __fmpz_mod_poly_factor_berlekamp(fmpz_mod_poly_factor_t factors,
fmpz_mod_poly_clear(x_pi2, ctx);
/* Row reduce Q - I */
- perm = _perm_init(n);
- nullity = n - fmpz_mat_rref_mod(perm, matrix, p);
- _perm_clear(perm);
+ nullity = n - fmpz_mod_mat_rref(matrix, matrix, ctx);
/* Find a basis for the nullspace */
basis =
diff --git a/src/fq_default_mat.h b/src/fq_default_mat.h
index 32dd684b4e..be6104bb47 100644
--- a/src/fq_default_mat.h
+++ b/src/fq_default_mat.h
@@ -1134,28 +1134,29 @@ FQ_DEFAULT_MAT_INLINE int fq_default_mat_inv(fq_default_mat_t B,
/* Solving *******************************************************************/
-FQ_DEFAULT_MAT_INLINE slong fq_default_mat_rref(fq_default_mat_t A,
+FQ_DEFAULT_MAT_INLINE slong fq_default_mat_rref(fq_default_mat_t B, const fq_default_mat_t A,
const fq_default_ctx_t ctx)
{
if (ctx->type == FQ_DEFAULT_FQ_ZECH)
{
- return fq_zech_mat_rref(A->fq_zech, ctx->ctx.fq_zech);
+ return fq_zech_mat_rref(B->fq_zech, A->fq_zech, ctx->ctx.fq_zech);
}
else if (ctx->type == FQ_DEFAULT_FQ_NMOD)
{
- return fq_nmod_mat_rref(A->fq_nmod, ctx->ctx.fq_nmod);
+ return fq_nmod_mat_rref(B->fq_nmod, A->fq_nmod, ctx->ctx.fq_nmod);
}
else if (ctx->type == FQ_DEFAULT_NMOD)
{
- return nmod_mat_rref(A->nmod);
+ nmod_mat_set(B->nmod, A->nmod);
+ return nmod_mat_rref(B->nmod);
}
else if (ctx->type == FQ_DEFAULT_FMPZ_MOD)
{
- return fmpz_mod_mat_rref(NULL, A->fmpz_mod, ctx->ctx.fmpz_mod.mod);
+ return fmpz_mod_mat_rref(B->fmpz_mod, A->fmpz_mod, ctx->ctx.fmpz_mod.mod);
}
else
{
- return fq_mat_rref(A->fq, ctx->ctx.fq);
+ return fq_mat_rref(B->fq, A->fq, ctx->ctx.fq);
}
}
diff --git a/src/fq_mat_templates.h b/src/fq_mat_templates.h
index 942e7b920b..03d9b072f6 100644
--- a/src/fq_mat_templates.h
+++ b/src/fq_mat_templates.h
@@ -244,7 +244,7 @@ int TEMPLATE(T, mat_inv)(TEMPLATE(T, mat_t) B, TEMPLATE(T, mat_t) A,
/* Solving *******************************************************************/
-slong TEMPLATE(T, mat_rref)(TEMPLATE(T, mat_t) A, const TEMPLATE(T, ctx_t) ctx);
+slong TEMPLATE(T, mat_rref)(TEMPLATE(T, mat_t) B, const TEMPLATE(T, mat_t) A, const TEMPLATE(T, ctx_t) ctx);
slong TEMPLATE(T, mat_reduce_row)(TEMPLATE(T, mat_t) A, slong * P, slong * L,
diff --git a/src/fq_mat_templates/nullspace.c b/src/fq_mat_templates/nullspace.c
index 43cd4efdc7..3b1a8f664e 100644
--- a/src/fq_mat_templates/nullspace.c
+++ b/src/fq_mat_templates/nullspace.c
@@ -30,7 +30,7 @@ TEMPLATE(T, mat_nullspace) (TEMPLATE(T, mat_t) X, const TEMPLATE(T, mat_t) A,
p = flint_malloc(sizeof(slong) * FLINT_MAX(m, n));
TEMPLATE(T, mat_init_set) (tmp, A, ctx);
- rank = TEMPLATE(T, mat_rref) (tmp, ctx);
+ rank = TEMPLATE(T, mat_rref) (tmp, tmp, ctx);
nullity = n - rank;
TEMPLATE(T, mat_zero) (X, ctx);
diff --git a/src/fq_mat_templates/rref.c b/src/fq_mat_templates/rref.c
index 01fba93732..f285ab67a9 100644
--- a/src/fq_mat_templates/rref.c
+++ b/src/fq_mat_templates/rref.c
@@ -18,7 +18,7 @@
#include "perm.h"
slong
-TEMPLATE(T, mat_rref) (TEMPLATE(T, mat_t) A, const TEMPLATE(T, ctx_t) ctx)
+TEMPLATE(T, mat_rref) (TEMPLATE(T, mat_t) B, const TEMPLATE(T, mat_t) A, const TEMPLATE(T, ctx_t) ctx)
{
slong i, j, k, n, rank;
slong *pivots;
@@ -27,18 +27,20 @@ TEMPLATE(T, mat_rref) (TEMPLATE(T, mat_t) A, const TEMPLATE(T, ctx_t) ctx)
TEMPLATE(T, struct) * e;
TEMPLATE(T, mat_t) U, V;
- if (TEMPLATE(T, mat_is_zero)(A, ctx))
+ TEMPLATE(T, mat_set)(B, A, ctx);
+
+ if (TEMPLATE(T, mat_is_zero)(B, ctx))
return 0;
- if (A->r == 1)
+ if (B->r == 1)
{
TEMPLATE(T, struct) * c;
slong i, j;
slong r = 0;
- for (i = 0; i < A->c; i++)
+ for (i = 0; i < B->c; i++)
{
- c = TEMPLATE(T, mat_entry)(A, 0, i);
+ c = TEMPLATE(T, mat_entry)(B, 0, i);
if (!TEMPLATE(T, is_zero)(c, ctx))
{
r = 1;
@@ -46,9 +48,9 @@ TEMPLATE(T, mat_rref) (TEMPLATE(T, mat_t) A, const TEMPLATE(T, ctx_t) ctx)
break;
TEMPLATE(T, inv)(c, c, ctx);
- for (j = i + 1;j < A->c; j++)
+ for (j = i + 1;j < B->c; j++)
{
- TEMPLATE(T, mul)(TEMPLATE(T, mat_entry)(A, 0, j), TEMPLATE(T, mat_entry)(A, 0, j), c, ctx);
+ TEMPLATE(T, mul)(TEMPLATE(T, mat_entry)(B, 0, j), TEMPLATE(T, mat_entry)(B, 0, j), c, ctx);
}
TEMPLATE(T, one)(c, ctx);
break;
@@ -59,17 +61,17 @@ TEMPLATE(T, mat_rref) (TEMPLATE(T, mat_t) A, const TEMPLATE(T, ctx_t) ctx)
n = A->c;
- P = _perm_init(TEMPLATE(T, mat_nrows) (A, ctx));
- rank = TEMPLATE(T, mat_lu) (P, A, 0, ctx);
+ P = _perm_init(TEMPLATE(T, mat_nrows) (B, ctx));
+ rank = TEMPLATE(T, mat_lu) (P, B, 0, ctx);
_perm_clear(P);
if (rank == 0)
return rank;
/* Clear L */
- for (i = 0; i < A->r; i++)
+ for (i = 0; i < B->r; i++)
for (j = 0; j < FLINT_MIN(i, rank); j++)
- TEMPLATE(T, zero) (TEMPLATE(T, mat_entry) (A, i, j), ctx);
+ TEMPLATE(T, zero) (TEMPLATE(T, mat_entry) (B, i, j), ctx);
/* We now reorder U to proper upper triangular form U | V
with U full-rank triangular, set V = U^(-1) V, and then
@@ -86,7 +88,7 @@ TEMPLATE(T, mat_rref) (TEMPLATE(T, mat_t) A, const TEMPLATE(T, ctx_t) ctx)
for (i = j = k = 0; i < rank; i++)
{
- while (TEMPLATE(T, is_zero) (TEMPLATE(T, mat_entry) (A, i, j), ctx))
+ while (TEMPLATE(T, is_zero) (TEMPLATE(T, mat_entry) (B, i, j), ctx))
{
nonpivots[k] = j;
k++;
@@ -106,7 +108,7 @@ TEMPLATE(T, mat_rref) (TEMPLATE(T, mat_t) A, const TEMPLATE(T, ctx_t) ctx)
{
for (j = 0; j <= i; j++)
{
- e = TEMPLATE(T, mat_entry) (A, j, pivots[i]);
+ e = TEMPLATE(T, mat_entry) (B, j, pivots[i]);
TEMPLATE(T, mat_entry_set) (U, j, i, e, ctx);
}
}
@@ -115,7 +117,7 @@ TEMPLATE(T, mat_rref) (TEMPLATE(T, mat_t) A, const TEMPLATE(T, ctx_t) ctx)
{
for (j = 0; j < rank; j++)
{
- e = TEMPLATE(T, mat_entry) (A, j, nonpivots[i]);
+ e = TEMPLATE(T, mat_entry) (B, j, nonpivots[i]);
TEMPLATE(T, mat_entry_set) (V, j, i, e, ctx);
}
}
@@ -129,12 +131,12 @@ TEMPLATE(T, mat_rref) (TEMPLATE(T, mat_t) A, const TEMPLATE(T, ctx_t) ctx)
{
if (i == j)
{
- TEMPLATE(T, one) (TEMPLATE(T, mat_entry) (A, j, pivots[i]),
+ TEMPLATE(T, one) (TEMPLATE(T, mat_entry) (B, j, pivots[i]),
ctx);
}
else
{
- TEMPLATE(T, zero) (TEMPLATE(T, mat_entry) (A, j, pivots[i]),
+ TEMPLATE(T, zero) (TEMPLATE(T, mat_entry) (B, j, pivots[i]),
ctx);
}
}
@@ -144,7 +146,7 @@ TEMPLATE(T, mat_rref) (TEMPLATE(T, mat_t) A, const TEMPLATE(T, ctx_t) ctx)
for (i = 0; i < n - rank; i++)
{
for (j = 0; j < rank; j++)
- TEMPLATE(T, mat_entry_set) (A, j, nonpivots[i],
+ TEMPLATE(T, mat_entry_set) (B, j, nonpivots[i],
TEMPLATE(T, mat_entry) (V, j, i), ctx);
}
diff --git a/src/fq_mat_templates/test/t-rref.c b/src/fq_mat_templates/test/t-rref.c
index 3e1a6d7f37..5041ce3ed5 100644
--- a/src/fq_mat_templates/test/t-rref.c
+++ b/src/fq_mat_templates/test/t-rref.c
@@ -94,7 +94,7 @@ TEST_TEMPLATE_FUNCTION_START(T, mat_rref, state)
TEMPLATE(T, mat_init_set) (B, A, ctx);
TEMPLATE(T, mat_init_set) (C, A, ctx);
- rank1 = TEMPLATE(T, mat_rref) (B, ctx);
+ rank1 = TEMPLATE(T, mat_rref) (B, B, ctx);
if (!check_rref_form(perm, B, rank1, ctx))
{
@@ -132,7 +132,7 @@ TEST_TEMPLATE_FUNCTION_START(T, mat_rref, state)
}
}
- rank2 = TEMPLATE(T, mat_rref) (D, ctx);
+ rank2 = TEMPLATE(T, mat_rref) (D, D, ctx);
equal = (rank1 == rank2);
if (equal)
diff --git a/src/fq_nmod_mpoly/mpolyu_gcdp_zippel.c b/src/fq_nmod_mpoly/mpolyu_gcdp_zippel.c
index 49e9c412ed..d3894b02ef 100644
--- a/src/fq_nmod_mpoly/mpolyu_gcdp_zippel.c
+++ b/src/fq_nmod_mpoly/mpolyu_gcdp_zippel.c
@@ -687,7 +687,7 @@ nmod_gcds_ret_t fq_nmod_mpolyu_gcds_zippel(
}
}
- fq_nmod_mat_rref(ML + s, ctx->fqctx);
+ fq_nmod_mat_rref(ML + s, ML + s, ctx->fqctx);
for (i = 0; i < (f->coeffs + s)->length; i++)
{
@@ -716,7 +716,7 @@ nmod_gcds_ret_t fq_nmod_mpolyu_gcds_zippel(
fq_nmod_mat_window_clear(Mwindow, ctx->fqctx);
}
- nullity = l - fq_nmod_mat_rref(MF, ctx->fqctx);
+ nullity = l - fq_nmod_mat_rref(MF, MF, ctx->fqctx);
if (nullity == 0)
{
diff --git a/src/fq_poly_factor_templates/factor_berlekamp.c b/src/fq_poly_factor_templates/factor_berlekamp.c
index d4055cd126..337bf0c430 100644
--- a/src/fq_poly_factor_templates/factor_berlekamp.c
+++ b/src/fq_poly_factor_templates/factor_berlekamp.c
@@ -148,7 +148,7 @@ __TEMPLATE(T, poly_factor_berlekamp) (TEMPLATE(T, poly_factor_t) factors,
TEMPLATE(T, poly_clear) (x_qi2, ctx);
/* Row reduce Q - I */
- nullity = n - TEMPLATE(T, mat_rref) (matrix, ctx);
+ nullity = n - TEMPLATE(T, mat_rref) (matrix, matrix, ctx);
/* Find a basis for the nullspace */
basis = flint_malloc(nullity * sizeof(TEMPLATE(T, poly_t)));