/**************************************************************************\ MODULE: mat_ZZ SUMMARY: Defines the class mat_ZZ. \**************************************************************************/ #include #include typedef Mat mat_ZZ; // backward compatibility void add(mat_ZZ& X, const mat_ZZ& A, const mat_ZZ& B); // X = A + B void sub(mat_ZZ& X, const mat_ZZ& A, const mat_ZZ& B); // X = A - B void negate(mat_ZZ& X, const mat_ZZ& A); // X = - A void mul(mat_ZZ& X, const mat_ZZ& A, const mat_ZZ& B); // X = A * B void mul(vec_ZZ& x, const mat_ZZ& A, const vec_ZZ& b); // x = A * b void mul(vec_ZZ& x, const vec_ZZ& a, const mat_ZZ& B); // x = a * B void mul(mat_ZZ& X, const mat_ZZ& A, const ZZ& b); void mul(mat_ZZ& X, const mat_ZZ& A, long b); // X = A * b void mul(mat_ZZ& X, const ZZ& a, const mat_ZZ& B); void mul(mat_ZZ& X, long a, const mat_ZZ& B); // X = a * B void determinant(ZZ& d, const mat_ZZ& A, long deterministic=0); ZZ determinant(const mat_ZZ& a, long deterministic=0); // d = determinant(A). If !deterministic, a randomized strategy may // be used that errs with probability at most 2^{-80}. void solve(ZZ& d, vec_ZZ& x, const mat_ZZ& A, const vec_ZZ& b, long deterministic=0) // computes d = determinant(A) and solves x*A = b*d if d != 0; A must // be a square matrix and have compatible dimensions with b. If // !deterministic, the computation of d may use a randomized strategy // that errs with probability 2^{-80}. void solve1(ZZ& d, vec_ZZ& x, const mat_ZZ& A, const vec_ZZ& b); // A must be a square matrix. // If A is singular, this routine sets d = 0 and returns. // Otherwise, it computes d, x such that x*A == b*d, // such that d > 0 and minimal. // Note that d is a positive divisor of the determinant, // and is not in general equal to the determinant. // The routine is deterministic, and uses a Hensel lifting strategy. // For backward compatability, there is also a routine called // HenselSolve1 that simply calls solve1. void inv(ZZ& d, mat_ZZ& X, const mat_ZZ& A, long deterministic=0); // computes d = determinant(A) and solves X*A = I*d if d != 0; A must // be a square matrix. If !deterministic, the computation of d may // use a randomized strategy that errs with probability 2^{-80}. // NOTE: See LLL.txt for routines that compute the kernel and // image of an integer matrix. // NOTE: See HNF.txt for a routine that computes Hermite Normal Forms. void sqr(mat_ZZ& X, const mat_ZZ& A); mat_ZZ sqr(const mat_ZZ& A); // X = A*A void inv(mat_ZZ& X, const mat_ZZ& A); mat_ZZ inv(const mat_ZZ& A); // X = A^{-1}; error is raised if |det(A)| != 1. void power(mat_ZZ& X, const mat_ZZ& A, const ZZ& e); mat_ZZ power(const mat_ZZ& A, const ZZ& e); void power(mat_ZZ& X, const mat_ZZ& A, long e); mat_ZZ power(const mat_ZZ& A, long e); // X = A^e; e may be negative (in which case A must be nonsingular). void ident(mat_ZZ& X, long n); mat_ZZ ident_mat_ZZ(long n); // X = n x n identity matrix long IsIdent(const mat_ZZ& A, long n); // test if A is the n x n identity matrix void diag(mat_ZZ& X, long n, const ZZ& d); mat_ZZ diag(long n, const ZZ& d); // X = n x n diagonal matrix with d on diagonal long IsDiag(const mat_ZZ& A, long n, const ZZ& d); // test if X is an n x n diagonal matrix with d on diagonal void transpose(mat_ZZ& X, const mat_ZZ& A); mat_ZZ transpose(const mat_ZZ& A); // X = transpose of A long CRT(mat_ZZ& a, ZZ& prod, const mat_zz_p& A); // Incremental Chinese Remaindering: If p is the current zz_p modulus with // (p, prod) = 1; Computes a' such that a' = a mod prod and a' = A mod p, // with coefficients in the interval (-p*prod/2, p*prod/2]; // Sets a := a', prod := p*prod, and returns 1 if a's value changed. // miscellaneous: void clear(mat_ZZ& a); // x = 0 (dimension unchanged) long IsZero(const mat_ZZ& a); // test if a is the zero matrix (any dimension) // operator notation: mat_ZZ operator+(const mat_ZZ& a, const mat_ZZ& b); mat_ZZ operator-(const mat_ZZ& a, const mat_ZZ& b); mat_ZZ operator*(const mat_ZZ& a, const mat_ZZ& b); mat_ZZ operator-(const mat_ZZ& a); // matrix/scalar multiplication: mat_ZZ operator*(const mat_ZZ& a, const ZZ& b); mat_ZZ operator*(const mat_ZZ& a, long b); mat_ZZ operator*(const ZZ& a, const mat_ZZ& b); mat_ZZ operator*(long a, const mat_ZZ& b); // matrix/vector multiplication: vec_ZZ operator*(const mat_ZZ& a, const vec_ZZ& b); vec_ZZ operator*(const vec_ZZ& a, const mat_ZZ& b); // assignment operator notation: mat_ZZ& operator+=(mat_ZZ& x, const mat_ZZ& a); mat_ZZ& operator-=(mat_ZZ& x, const mat_ZZ& a); mat_ZZ& operator*=(mat_ZZ& x, const mat_ZZ& a); mat_ZZ& operator*=(mat_ZZ& x, const ZZ& a); mat_ZZ& operator*=(mat_ZZ& x, long a); vec_ZZ& operator*=(vec_ZZ& x, const mat_ZZ& a);