pqc/external/flint-2.4.3/fmpz_mat/doc/fmpz_mat.txt

/*=============================================================================

    This file is part of FLINT.

    FLINT is free software; you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
    the Free Software Foundation; either version 2 of the License, or
    (at your option) any later version.

    FLINT is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU General Public License for more details.

    You should have received a copy of the GNU General Public License
    along with FLINT; if not, write to the Free Software
    Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301 USA

=============================================================================*/
/******************************************************************************

    Copyright (C) 2010 William Hart
    Copyright (C) 2010-2011 Andy Novocin
    Copyright (C) 2010-2011 Fredrik Johansson

******************************************************************************/


*******************************************************************************

    Memory management

*******************************************************************************

void fmpz_mat_init(fmpz_mat_t mat, slong rows, slong cols)

    Initialises a matrix with the given number of rows and columns for use. 

void fmpz_mat_clear(fmpz_mat_t mat)

    Clears the given matrix.

*******************************************************************************

    Basic assignment and manipulation

*******************************************************************************

void fmpz_mat_set(fmpz_mat_t mat1, const fmpz_mat_t mat2)

    Sets \code{mat1} to a copy of \code{mat2}. The dimensions of 
    \code{mat1} and \code{mat2} must be the same.

void fmpz_mat_init_set(fmpz_mat_t mat, const fmpz_mat_t src)

    Initialises the matrix \code{mat} to the same size as \code{src} and 
    sets it to a copy of \code{src}.

void fmpz_mat_swap(fmpz_mat_t mat1, fmpz_mat_t mat2)

    Swaps two matrices. The dimensions of \code{mat1} and \code{mat2} 
    are allowed to be different.

fmpz * fmpz_mat_entry(fmpz_mat_t mat, slong i, slong j)

    Returns a reference to the entry of \code{mat} at row $i$ and column $j$.
    This reference can be passed as an input or output variable to any
    function in the \code{fmpz} module for direct manipulation.

    Both $i$ and $j$ must not exceed the dimensions of the matrix.

    This function is implemented as a macro.

void fmpz_mat_zero(fmpz_mat_t mat)

    Sets all entries of \code{mat} to 0.

void fmpz_mat_one(fmpz_mat_t mat)

    Sets \code{mat} to the unit matrix, having ones on the main diagonal
    and zeroes elsewhere. If \code{mat} is nonsquare, it is set to the
    truncation of a unit matrix.


*******************************************************************************

    Random matrix generation

*******************************************************************************

void fmpz_mat_randbits(fmpz_mat_t mat, flint_rand_t state, mp_bitcnt_t bits)

    Sets the entries of \code{mat} to random signed integers whose absolute 
    values have the given number of binary bits.

void fmpz_mat_randtest(fmpz_mat_t mat, flint_rand_t state, mp_bitcnt_t bits)

    Sets the entries of \code{mat} to random signed integers whose 
    absolute values have a random number of bits up to the given number 
    of bits inclusive.

void fmpz_mat_randintrel(fmpz_mat_t mat, flint_rand_t state, mp_bitcnt_t bits)

    Sets \code{mat} to be a random \emph{integer relations} matrix, with 
    signed entries up to the given number of bits.

    The number of columns of \code{mat} must be equal to one more than 
    the number of rows. The format of the matrix is a set of random integers 
    in the left hand column and an identity matrix in the remaining square 
    submatrix.

void fmpz_mat_randsimdioph(fmpz_mat_t mat, flint_rand_t state, 
                                           mp_bitcnt_t bits, mp_bitcnt_t bits2)

    Sets \code{mat} to a random \emph{simultaneous diophantine} matrix.

    The matrix must be square. The top left entry is set to \code{2^bits2}. 
    The remainder of that row is then set to signed random integers of the 
    given number of binary bits. The remainder of the first column is zero. 
    Running down the rest of the diagonal are the values \code{2^bits} with 
    all remaining entries zero.

void fmpz_mat_randntrulike(fmpz_mat_t mat, flint_rand_t state, 
                                                     mp_bitcnt_t bits, ulong q)

    Sets a square matrix \code{mat} of even dimension to a random 
    \emph{NTRU like} matrix.

    The matrix is broken into four square submatrices. The top left submatrix
    is set to the identity. The bottom left submatrix is set to the zero 
    matrix. The bottom right submatrix is set to $q$ times the identity matrix.
    Finally the top right submatrix has the following format. A random vector
    $h$ of length $r/2$ is created, with random signed entries of the given 
    number of bits. Then entry $(i, j)$ of the submatrix is set to 
    $h[i + j \bmod{r/2}]$. 

void fmpz_mat_randntrulike2(fmpz_mat_t mat, flint_rand_t state, 
                                                     mp_bitcnt_t bits, ulong q)

    Sets a square matrix \code{mat} of even dimension to a random 
    \emph{NTRU like} matrix.

    The matrix is broken into four square submatrices. The top left submatrix
    is set to $q$ times the identity matrix. The top right submatrix is set to 
    the zero matrix. The bottom right submatrix is set to the identity matrix.
    Finally the bottom left submatrix has the following format. A random vector
    $h$ of length $r/2$ is created, with random signed entries of the given 
    number of bits. Then entry $(i, j)$ of the submatrix is set to 
    $h[i + j \bmod{r/2}]$.

void fmpz_mat_randajtai(fmpz_mat_t mat, flint_rand_t state, double alpha)

    Sets a square matrix \code{mat} to a random \emph{ajtai} matrix. 
    The diagonal entries $(i, i)$ are set to a random entry in the range 
    $[1, 2^{b-1}]$ inclusive where $b = \floor{(2 r - i)^\alpha}$ for some 
    double parameter~$\alpha$. The entries below the diagonal in column~$i$ 
    are set to a random entry in the range $(-2^b + 1, 2^b - 1)$ whilst the 
    entries to the right of the diagonal in row~$i$ are set to zero. 

int fmpz_mat_randpermdiag(fmpz_mat_t mat, flint_rand_t state, 
                                                   const fmpz * diag, slong n)

    Sets \code{mat} to a random permutation of the rows and columns of a
    given diagonal matrix. The diagonal matrix is specified in the form of
    an array of the $n$ initial entries on the main diagonal.

    The return value is $0$ or $1$ depending on whether the permutation is
    even or odd.

void fmpz_mat_randrank(fmpz_mat_t mat, flint_rand_t state, slong rank, 
                                                              mp_bitcnt_t bits)

    Sets \code{mat} to a random sparse matrix with the given rank, 
    having exactly as many non-zero elements as the rank, with the 
    nonzero elements being random integers of the given bit size.

    The matrix can be transformed into a dense matrix with unchanged
    rank by subsequently calling \code{fmpz_mat_randops()}.

void fmpz_mat_randdet(fmpz_mat_t mat, flint_rand_t state, const fmpz_t det)

    Sets \code{mat} to a random sparse matrix with minimal number of
    nonzero entries such that its determinant has the given value.

    Note that the matrix will be zero if \code{det} is zero.  
    In order to generate a non-zero singular matrix, the function 
    \code{fmpz_mat_randrank()} can be used.

    The matrix can be transformed into a dense matrix with unchanged
    determinant by subsequently calling \code{fmpz_mat_randops()}.

void fmpz_mat_randops(fmpz_mat_t mat, flint_rand_t state, slong count)

    Randomises \code{mat} by performing elementary row or column operations.
    More precisely, at most \code{count} random additions or subtractions of
    distinct rows and columns will be performed. This leaves the rank
    (and for square matrices, the determinant) unchanged.


*******************************************************************************

    Input and output

*******************************************************************************

int fmpz_mat_fprint(FILE * file, const fmpz_mat_t mat)

    Prints the given matrix to the stream \code{file}.  The format is 
    the number of rows, a space, the number of columns, two spaces, then 
    a space separated list of coefficients, one row after the other.

    In case of success, returns a positive value;  otherwise, returns 
    a non-positive value.

int fmpz_mat_fprint_pretty(FILE * file, const fmpz_mat_t mat)

    Prints the given matrix to the stream \code{file}.  The format is an 
    opening square bracket then on each line a row of the matrix, followed 
    by a closing square bracket. Each row is written as an opening square 
    bracket followed by a space separated list of coefficients followed 
    by a closing square bracket.

    In case of success, returns a positive value;  otherwise, returns 
    a non-positive value.

int fmpz_mat_print(const fmpz_mat_t mat)

    Prints the given matrix to the stream \code{stdout}.  For further 
    details, see \code{fmpz_mat_fprint()}.

int fmpz_mat_print_pretty(const fmpz_mat_t mat)

    Prints the given matrix to \code{stdout}.  For further details, 
    see \code{fmpz_mat_fprint_pretty()}.

int fmpz_mat_fread(FILE* file, fmpz_mat_t mat)

    Reads a matrix from the stream \code{file}, storing the result 
    in \code{mat}.  The expected format is the number of rows, a 
    space, the number of columns, two spaces, then a space separated
    list of coefficients, one row after the other.

    In case of success, returns a positive number.  In case of failure, 
    returns a non-positive value.

int fmpz_mat_read(fmpz_mat_t mat)

    Reads a matrix from \code{stdin}, storing the result 
    in \code{mat}.

    In case of success, returns a positive number.  In case of failure, 
    returns a non-positive value.

*******************************************************************************

    Comparison

*******************************************************************************

int fmpz_mat_equal(const fmpz_mat_t mat1, const fmpz_mat_t mat2)

    Returns a non-zero value if \code{mat1} and \code{mat2} have 
    the same dimensions and entries, and zero otherwise.

int fmpz_mat_is_zero(const fmpz_mat_t mat)

    Returns a non-zero value if all entries \code{mat} are zero, and
    otherwise returns zero.

int fmpz_mat_is_empty(const fmpz_mat_t mat)

    Returns a non-zero value if the number of rows or the number of
    columns in \code{mat} is zero, and otherwise returns
    zero.

int fmpz_mat_is_square(const fmpz_mat_t mat)

    Returns a non-zero value if the number of rows is equal to the
    number of columns in \code{mat}, and otherwise returns zero.

*******************************************************************************

    Transpose

*******************************************************************************

void fmpz_mat_transpose(fmpz_mat_t B, const fmpz_mat_t A)

    Sets $B$ to $A^T$, the transpose of $A$. Dimensions must be compatible.
    $A$ and $B$ are allowed to be the same object if $A$ is a square matrix.


*******************************************************************************

    Modular reduction and reconstruction

*******************************************************************************

void fmpz_mat_get_nmod_mat(nmod_mat_t Amod, const fmpz_mat_t A)

    Sets the entries of \code{Amod} to the entries of \code{A} reduced
    by the modulus of \code{Amod}.

void fmpz_mat_set_nmod_mat(fmpz_mat_t A, const nmod_mat_t Amod)

    Sets the entries of \code{Amod} to the residues in \code{Amod},
    normalised to the interval $-m/2 <= r < m/2$ where $m$ is the modulus.

void fmpz_mat_set_nmod_mat_unsigned(fmpz_mat_t A, const nmod_mat_t Amod)

    Sets the entries of \code{Amod} to the residues in \code{Amod},
    normalised to the interval $0 <= r < m$ where $m$ is the modulus.

void fmpz_mat_CRT_ui(fmpz_mat_t res, const fmpz_mat_t mat1,
                        const fmpz_t m1, const nmod_mat_t mat2, int sign)

    Given \code{mat1} with entries modulo \code{m} and \code{mat2}
    with modulus $n$, sets \code{res} to the CRT reconstruction modulo $mn$
    with entries satisfying $-mn/2 <= c < mn/2$ (if sign = 1)
    or $0 <= c < mn$ (if sign = 0).

void fmpz_mat_multi_mod_ui_precomp(nmod_mat_t * residues, slong nres, 
    const fmpz_mat_t mat, fmpz_comb_t comb, fmpz_comb_temp_t temp)

    Sets each of the \code{nres} matrices in \code{residues} to \code{mat}
    reduced modulo the modulus of the respective matrix, given
    precomputed \code{comb} and \code{comb_temp} structures.

void fmpz_mat_multi_mod_ui(nmod_mat_t * residues, slong nres,
        const fmpz_mat_t mat)

    Sets each of the \code{nres} matrices in \code{residues} to \code{mat}
    reduced modulo the modulus of the respective matrix.

    This function is provided for convenience purposes.
    For reducing or reconstructing multiple integer matrices over the same
    set of moduli, it is faster to use\\ \code{fmpz_mat_multi_mod_precomp}.

void fmpz_mat_multi_CRT_ui_precomp(fmpz_mat_t mat, nmod_mat_t * const residues,
    slong nres, fmpz_comb_t comb, fmpz_comb_temp_t temp, int sign)

    Reconstructs \code{mat} from its images modulo the \code{nres} matrices
    in \code{residues}, given precomputed \code{comb} and \code{comb_temp}
    structures.

void fmpz_mat_multi_CRT_ui(fmpz_mat_t mat, nmod_mat_t * const residues,
    slong nres, int sign)

    Reconstructs \code{mat} from its images modulo the \code{nres} matrices
    in \code{residues}.

    This function is provided for convenience purposes.
    For reducing or reconstructing multiple integer matrices over the same
    set of moduli, it is faster to use\\ \code{fmpz_mat_multi_CRT_ui_precomp}.

*******************************************************************************

    Addition and subtraction

*******************************************************************************

void fmpz_mat_add(fmpz_mat_t C, const fmpz_mat_t A, const fmpz_mat_t B)

    Sets \code{C} to the elementwise sum $A + B$. All inputs must
    be of the same size. Aliasing is allowed.

void fmpz_mat_sub(fmpz_mat_t C, const fmpz_mat_t A, const fmpz_mat_t B)

    Sets \code{C} to the elementwise difference $A - B$. All inputs must
    be of the same size. Aliasing is allowed.

void fmpz_mat_neg(fmpz_mat_t B, const fmpz_mat_t A)

    Sets \code{B} to the elementwise negation of \code{A}. Both inputs
    must be of the same size. Aliasing is allowed.

*******************************************************************************

    Matrix-scalar arithmetic

*******************************************************************************

void fmpz_mat_scalar_mul_si(fmpz_mat_t B, const fmpz_mat_t A, slong c)

void fmpz_mat_scalar_mul_ui(fmpz_mat_t B, const fmpz_mat_t A, ulong c)

void fmpz_mat_scalar_mul_fmpz(fmpz_mat_t B, const fmpz_mat_t A, const fmpz_t c)

    Set \code{A = B*c} where \code{B} is an \code{fmpz_mat_t} and \code{c}
    is a scalar respectively of type \code{slong}, \code{ulong},
    or \code{fmpz_t}. The dimensions of \code{A} and \code{B} must
    be compatible.

void fmpz_mat_scalar_addmul_si(fmpz_mat_t B, const fmpz_mat_t A, slong c)

void fmpz_mat_scalar_addmul_ui(fmpz_mat_t B, const fmpz_mat_t A, ulong c)

void fmpz_mat_scalar_addmul_fmpz(fmpz_mat_t B, const fmpz_mat_t A,
                                    const fmpz_t c)

    Set \code{A = A + B*c} where \code{B} is an \code{fmpz_mat_t} and \code{c}
    is a scalar respectively of type \code{slong}, \code{ulong},
    or \code{fmpz_t}. The dimensions of \code{A} and \code{B} must
    be compatible.

void fmpz_mat_scalar_submul_si(fmpz_mat_t B, const fmpz_mat_t A, slong c)

void fmpz_mat_scalar_submul_ui(fmpz_mat_t B, const fmpz_mat_t A, ulong c)

void fmpz_mat_scalar_submul_fmpz(fmpz_mat_t B, const fmpz_mat_t A,
                                    const fmpz_t c)

    Set \code{A = A - B*c} where \code{B} is an \code{fmpz_mat_t} and \code{c}
    is a scalar respectively of type \code{slong}, \code{ulong},
    or \code{fmpz_t}. The dimensions of \code{A} and \code{B} must
    be compatible.

void fmpz_mat_scalar_addmul_nmod_mat_ui(fmpz_mat_t B, const nmod_mat_t A,
                                    ulong c)
void fmpz_mat_scalar_addmul_nmod_mat_fmpz(fmpz_mat_t B, const nmod_mat_t A,
                                    const fmpz_t c)

    Set \code{A = A + B*c} where \code{B} is an \code{nmod_mat_t} and \code{c}
    is a scalar respectively of type \code{ulong} or \code{fmpz_t}.
    The dimensions of \code{A} and \code{B} must be compatible.

void fmpz_mat_scalar_divexact_si(fmpz_mat_t B, const fmpz_mat_t A, slong c)

void fmpz_mat_scalar_divexact_ui(fmpz_mat_t B, const fmpz_mat_t A, ulong c)

void fmpz_mat_scalar_divexact_fmpz(fmpz_mat_t B, const fmpz_mat_t A,
                                    const fmpz_t c)

    Set \code{A = B / c}, where \code{B} is an \code{fmpz_mat_t} and \code{c}
    is a scalar respectively of type \code{slong}, \code{ulong},
    or \code{fmpz_t}, which is assumed to divide all elements of
    \code{B} exactly.


*******************************************************************************

    Matrix multiplication

*******************************************************************************

void fmpz_mat_mul(fmpz_mat_t C, const fmpz_mat_t A, const fmpz_mat_t B)

    Sets \code{C} to the matrix product $C = A B$. The matrices must have
    compatible dimensions for matrix multiplication. Aliasing
    is allowed.

    This function automatically switches between classical and
    multimodular multiplication, based on a heuristic comparison of
    the dimensions and entry sizes.

void fmpz_mat_mul_classical(fmpz_mat_t C, 
                                        const fmpz_mat_t A, const fmpz_mat_t B)

    Sets \code{C} to the matrix product $C = A B$ computed using
    classical matrix algorithm.

    The matrices must have compatible dimensions for matrix multiplication.
    No aliasing is allowed.

void _fmpz_mat_mul_multi_mod(fmpz_mat_t C,
            const fmpz_mat_t A, const fmpz_mat_t B, mp_bitcnt_t bits)

void fmpz_mat_mul_multi_mod(fmpz_mat_t C,
            const fmpz_mat_t A, const fmpz_mat_t B)

    Sets \code{C} to the matrix product $C = AB$ computed using a multimodular 
    algorithm. $C$ is computed modulo several small prime numbers
    and reconstructed using the Chinese Remainder Theorem. This generally
    becomes more efficient than classical multiplication for large matrices.

    The \code{bits} parameter is a bound for the bit size of largest 
    element of $C$, or twice the absolute value of the largest element 
    if any elements of $C$ are negative. The function
    \code{fmpz_mat_mul_multi_mod} calculates a rigorous bound automatically.
    If the default bound is too pessimistic, \code{_fmpz_mat_mul_multi_mod}
    can be used with a custom bound.

    The matrices must have compatible dimensions for matrix multiplication.
    No aliasing is allowed.

void fmpz_mat_sqr(fmpz_mat_t B, const fmpz_mat_t A)

    Sets \code{B} to the square of the matrix \code{A}, which must be
    a square matrix. Aliasing is allowed.

void fmpz_mat_pow(fmpz_mat_t B, const fmpz_mat_t A, ulong e)

    Sets \code{B} to the matrix \code{A} raised to the power \code{e},
    where \code{A} must be a square matrix. Aliasing is allowed.


*******************************************************************************

    Inverse

*******************************************************************************

int fmpz_mat_inv(fmpz_mat_t Ainv, fmpz_t den, const fmpz_mat_t A)

    Sets (\code{Ainv}, \code{den}) to the inverse matrix of \code{A}.
    Returns 1 if \code{A} is nonsingular and 0 if \code{A} is singular.
    Aliasing of \code{Ainv} and \code{A} is allowed.

    The denominator is not guaranteed to be minimal, but is guaranteed
    to be a divisor of the determinant of \code{A}.

    This function uses a direct formula for matrices of size two or less,
    and otherwise solves for the identity matrix using
    fraction-free LU decomposition.


*******************************************************************************

    Trace

*******************************************************************************

void fmpz_mat_trace(fmpz_t trace, const fmpz_mat_t mat)

    Computes the trace of the matrix, i.e. the sum of the entries on
    the main diagonal. The matrix is required to be square.


*******************************************************************************

    Determinant

*******************************************************************************

void fmpz_mat_det(fmpz_t det, const fmpz_mat_t A)

    Sets \code{det} to the determinant of the square matrix $A$.
    The matrix of dimension $0 \times 0$ is defined to have determinant 1.

    This function automatically chooses between \code{fmpz_mat_det_cofactor},\\
    \code{fmpz_mat_det_bareiss}, \code{fmpz_mat_det_modular} and\\
    \code{fmpz_mat_det_modular_accelerated}
    (with \code{proved} = 1), depending on the size of the matrix
    and its entries.

void fmpz_mat_det_cofactor(fmpz_t det, const fmpz_mat_t A)

    Sets \code{det} to the determinant of the square matrix $A$
    computed using direct cofactor expansion. This function only
    supports matrices up to size $4 \times 4$.

void fmpz_mat_det_bareiss(fmpz_t det, const fmpz_mat_t A)

    Sets \code{det} to the determinant of the square matrix $A$
    computed using the Bareiss algorithm. A copy of the input matrix is
    row reduced using fraction-free Gaussian elimination, and the
    determinant is read off from the last element on the main
    diagonal.

void fmpz_mat_det_modular(fmpz_t det, const fmpz_mat_t A, int proved)

    Sets \code{det} to the determinant of the square matrix $A$
    (if \code{proved} = 1), or a probabilistic value for the
    determinant (\code{proved} = 0), computed using a multimodular
    algorithm.

    The determinant is computed modulo several small primes and
    reconstructed using the Chinese Remainder Theorem.
    With \code{proved} = 1, sufficiently many primes are chosen
    to satisfy the bound computed by \code{fmpz_mat_det_bound}.
    With \code{proved} = 0, the determinant is considered determined
    if it remains unchanged modulo several consecutive primes
    (currently if their product exceeds $2^{100}$).

void fmpz_mat_det_modular_accelerated(fmpz_t det,
        const fmpz_mat_t A, int proved)

    Sets \code{det} to the determinant of the square matrix $A$
    (if \code{proved} = 1), or a probabilistic value for the
    determinant (\code{proved} = 0), computed using a multimodular
    algorithm.

    This function uses the same basic algorithm as \code{fmpz_mat_det_modular},
    but instead of computing $\det(A)$ directly, it generates a divisor $d$
    of $\det(A)$ and then computes $x = \det(A) / d$ modulo several
    small primes not dividing $d$. This typically accelerates the
    computation by requiring fewer primes for large matrices, since $d$
    with high probability will be nearly as large as the determinant.
    This trick is described in \citep{AbbottBronsteinMulders1999}.

void fmpz_mat_det_modular_given_divisor(fmpz_t det, const fmpz_mat_t A,
        const fmpz_t d, int proved)

    Given a positive divisor $d$ of $\det(A)$, sets \code{det} to the
    determinant of the square matrix $A$ (if \code{proved} = 1), or a
    probabilistic value for the determinant (\code{proved} = 0), computed
    using a multimodular algorithm.

void fmpz_mat_det_bound(fmpz_t bound, const fmpz_mat_t A)

    Sets \code{bound} to a nonnegative integer $B$ such that
    $|\det(A)| \le B$. Assumes $A$ to be a square matrix.
    The bound is computed from the Hadamard inequality
    $|\det(A)| \le \prod \|a_i\|_2$ where the product is taken
    over the rows $a_i$ of $A$.

void fmpz_mat_det_divisor(fmpz_t d, const fmpz_mat_t A)

    Sets $d$ to some positive divisor of the determinant of the given
    square matrix $A$, if the determinant is nonzero. If $|\det(A)| = 0$,
    $d$ will always be set to zero.

    A divisor is obtained by solving $Ax = b$ for an arbitrarily chosen
    right-hand side $b$ using Dixon's algorithm and computing the least
    common multiple of the denominators in $x$. This yields a divisor $d$
    such that $|\det(A)| / d$ is tiny with very high probability.


*******************************************************************************

    Charactersitic polynomial

*******************************************************************************

void _fmpz_mat_charpoly(fmpz * cp, const fmpz_mat_t mat)

    Sets \code{(cp, n+1)} to the characteristic polynomial of 
    an $n \times n$ square matrix.

void fmpz_mat_charpoly(fmpz_poly_t cp, const fmpz_mat_t mat)

    Computes the characteristic polynomial of length $n + 1$ of 
    an $n \times n$ square matrix.

*******************************************************************************

    Rank

*******************************************************************************

slong fmpz_mat_rank(const fmpz_mat_t A)

    Returns the rank, that is, the number of linearly independent columns
    (equivalently, rows), of $A$. The rank is computed by row reducing
    a copy of $A$.


*******************************************************************************

    Nonsingular solving

    The following functions allow solving matrix-matrix equations $AX = B$
    where the system matrix $A$ is square and has full rank. The solving
    is implicitly done over the field of rational numbers: except
    where otherwise noted, an integer matrix $\hat X$ and a separate
    denominator $d$ (\code{den}) are computed such that $A(\hat X/d) = b$,
    equivalently such that $A\hat X = bd$ holds over the integers.

    No guarantee is made that the numerators and denominator
    are reduced to lowest terms, but the denominator is always guaranteed
    to be a divisor of the determinant of $A$. If $A$ is singular,
    \code{den} will be set to zero and the elements of the solution
    vector or matrix will have undefined values. No aliasing is
    allowed between arguments.

*******************************************************************************

int fmpz_mat_solve(fmpz_mat_t X, fmpz_t den,
                    const fmpz_mat_t A, const fmpz_mat_t B)

    Solves the equation $AX = B$ for nonsingular $A$. More precisely, computes
    (\code{X}, \code{den}) such that $AX = B \times \operatorname{den}$.
    Returns 1 if $A$ is nonsingular and 0 if $A$ is singular.
    The computed denominator will not generally be minimal.

    This function uses Cramer's rule for small systems and
    fraction-free LU decomposition followed by fraction-free forward
    and back substitution for larger systems.

    Note that for very large systems, it is faster to compute a modular
    solution using \code{fmpz_mat_solve_dixon}.

int fmpz_mat_solve_fflu(fmpz_mat_t X, fmpz_t den,
                            const fmpz_mat_t A, const fmpz_mat_t B)

    Solves the equation $AX = B$ for nonsingular $A$. More precisely, computes
    (\code{X}, \code{den}) such that $AX = B \times \operatorname{den}$.
    Returns 1 if $A$ is nonsingular and 0 if $A$ is singular.
    The computed denominator will not generally be minimal.

    Uses fraction-free LU decomposition followed by fraction-free
    forward and back substitution.

void fmpz_mat_solve_fflu_precomp(fmpz_mat_t X,
                    const slong * perm,
                    const fmpz_mat_t FFLU, const fmpz_mat_t B)

    Performs fraction-free forward and back substitution given a precomputed
    fraction-free LU decomposition and corresponding permutation.

int fmpz_mat_solve_cramer(fmpz_mat_t X, fmpz_t den,
                            const fmpz_mat_t A, const fmpz_mat_t B)

    Solves the equation $AX = B$ for nonsingular $A$. More precisely, computes
    (\code{X}, \code{den}) such that $AX = B \times \operatorname{den}$.
    Returns 1 if $A$ is nonsingular and 0 if $A$ is singular.

    Uses Cramer's rule. Only systems of size up to $3 \times 3$ are allowed.

void fmpz_mat_solve_bound(fmpz_t N, fmpz_t D, const fmpz_mat_t A,
    const fmpz_mat_t B)

    Assuming that $A$ is nonsingular, computes integers $N$ and $D$
    such that the reduced numerators and denominators $n/d$ in
    $A^{-1} B$ satisfy the bounds $0 \le |n| \le N$ and $0 \le d \le D$.

int fmpz_mat_solve_dixon(fmpz_mat_t X, fmpz_t M, const fmpz_mat_t A,
        const fmpz_mat_t B)

    Solves $AX = B$ given a nonsingular square matrix $A$ and a matrix $B$ of
    compatible dimensions, using a modular algorithm. In particular,
    Dixon's p-adic lifting algorithm is used (currently a non-adaptive version).
    This is generally the preferred method for large dimensions.

    More precisely, this function computes an integer $M$ and an integer
    matrix $X$ such that $AX = B \bmod M$ and such that all the reduced
    numerators and denominators of the elements $x = p/q$ in the full
    solution satisfy $2|p|q < B$. As such, the explicit rational solution
    matrix can be recovered uniquely by passing the output of this
    function to \code{fmpq_mat_set_fmpz_mat_mod}.

    A nonzero value is returned if $A$ is nonsingular. If $A$ is singular,
    zero is returned and the values of the output variables will be
    undefined.

    Aliasing between input and output matrices is allowed.

*******************************************************************************

    Row reduction

*******************************************************************************

slong fmpz_mat_find_pivot_any(const fmpz_mat_t mat,
                                    slong start_row, slong end_row, slong c)

    Attempts to find a pivot entry for row reduction.
    Returns a row index $r$ between \code{start_row} (inclusive) and
    \code{stop_row} (exclusive) such that column $c$ in \code{mat} has
    a nonzero entry on row $r$, or returns -1 if no such entry exists.

    This implementation simply chooses the first nonzero entry from
    it encounters. This is likely to be a nearly optimal choice if all
    entries in the matrix have roughly the same size, but can lead to
    unnecessary coefficient growth if the entries vary in size.

slong fmpz_mat_fflu(fmpz_mat_t B, fmpz_poly_t den, slong * perm,
                            const fmpz_mat_t A, int rank_check)

    Uses fraction-free Gaussian elimination to set (\code{B}, \code{den}) to a
    fraction-free LU decomposition of \code{A} and returns the
    rank of \code{A}. Aliasing of \code{A} and \code{B} is allowed.

    Pivot elements are chosen with \code{fmpz_mat_find_pivot_any}.
    If \code{perm} is non-\code{NULL}, the permutation of
    rows in the matrix will also be applied to \code{perm}.

    If \code{rank_check} is set, the function aborts and returns 0 if the
    matrix is detected not to have full rank without completing the
    elimination.

    The denominator \code{den} is set to $\pm \operatorname{det}(S)$ where
    $S$ is an appropriate submatrix of $A$ ($S = A$ if $A$ is square)
    and the sign is decided by the parity of the permutation. Note that the
    determinant is not generally the minimal denominator.

    The fraction-free LU decomposition is defined in \citep{NakTurWil1997}.

slong fmpz_mat_rref(fmpz_mat_t B, fmpz_t den, const fmpz_mat_t A)

    Sets (\code{B}, \code{den}) to the reduced row echelon form of \code{A}
    and returns the rank of \code{A}. Aliasing of \code{A} and \code{B}
    is allowed.

    The algorithm proceeds by first computing a row echelon form using
    \code{fmpz_mat_fflu}. Letting the upper part of this matrix be
    $(U | V) P$ where $U$ is full rank upper triangular and $P$ is a
    permutation matrix, we obtain the rref by setting $V$ to $U^{-1} V$
    using back substitution. Scaling each completed row in the back
    substitution to the denominator \code{den}, we avoid introducing
    new fractions. This strategy is equivalent to the fraction-free
    Gauss-Jordan elimination in \citep{NakTurWil1997}, but faster since
    only the part $V$ corresponding to the null space has to be updated.

    The denominator \code{den} is set to $\pm \operatorname{det}(S)$ where
    $S$ is an appropriate submatrix of $A$ ($S = A$ if $A$ is square).
    Note that the determinant is not generally the minimal denominator.


*******************************************************************************

    Modular gaussian elimination

*******************************************************************************

slong fmpz_mat_rref_mod(slong * perm, fmpz_mat_t A, const fmpz_t p)

    Uses fraction-free Gauss-Jordan elimination to set \code{A}
    to its reduced row echelon form and returns the rank of \code{A}.
    All computations are done modulo p.

    Pivot elements are chosen with \code{fmpz_mat_find_pivot_any}.
    If \code{perm} is non-\code{NULL}, the permutation of
    rows in the matrix will also be applied to \code{perm}.

*******************************************************************************

    Nullspace

*******************************************************************************

slong fmpz_mat_nullspace(fmpz_mat_t B, const fmpz_mat_t A)

    Computes a basis for the right rational nullspace of $A$ and returns
    the dimension of the nullspace (or nullity). $B$ is set to a matrix with
    linearly independent columns and maximal rank such that $AB = 0$
    (i.e. $Ab = 0$ for each column $b$ in $B$), and the rank of $B$ is
    returned.

    In general, the entries in $B$ will not be minimal: in particular,
    the pivot entries in $B$ will generally differ from unity.
    $B$ must be allocated with sufficient space to represent the result
    (at most $n \times n$ where $n$ is the number of column of $A$).


*******************************************************************************

    Echelon form

*******************************************************************************

slong fmpz_mat_rref_fraction_free(slong * perm, fmpz_mat_t B, fmpz_t den,
        const fmpz_mat_t A)

    Computes an integer matrix \code{B} and an integer \code{den} such that
    \code{B / den} is the unique row reduced echelon form (RREF) of \code{A}
    and returns the rank, i.e. the number of nonzero rows in \code{B}.

    Aliasing of \code{B} and \code{A} is allowed, with an in-place
    computation being more efficient. The size of \code{B} must be
    the same as that of \code{A}.

    The permutation order will be written to \code{perm} unless this
    argument is \code{NULL}. That is, row \code{i} of the output matrix will
    correspond to row \code{perm[i]} of the input matrix.

    The denominator will always be a divisor of the determinant of (some
    submatrix of) $A$, but is not guaranteed to be minimal or canonical in
    any other sense.