Dense matrices over the integer ring¶
AUTHORS:
- William Stein
- Robert Bradshaw
- Marc Masdeu (August 2014). Implemented using FLINT, see trac ticket #16803.
- Jeroen Demeyer (October 2014): lots of fixes, see trac ticket #17090 and trac ticket #17094.
- Vincent Delecroix (February 2015): make it faster, see trac ticket #17822.
- Vincent Delecroix (May 2017): removed duplication of entries and cleaner linbox interface
EXAMPLES:
sage: a = matrix(ZZ, 3,3, range(9)); a
[0 1 2]
[3 4 5]
[6 7 8]
sage: a.det()
0
sage: a[0,0] = 10; a.det()
-30
sage: a.charpoly()
x^3 - 22*x^2 + 102*x + 30
sage: b = -3*a
sage: a == b
False
sage: b < a
True
-
class
sage.matrix.matrix_integer_dense.
Matrix_integer_dense
¶ Bases:
sage.matrix.matrix_dense.Matrix_dense
Matrix over the integers, implemented using FLINT.
On a 32-bit machine, they can have at most \(2^{32}-1\) rows or columns. On a 64-bit machine, matrices can have at most \(2^{64}-1\) rows or columns.
EXAMPLES:
sage: a = MatrixSpace(ZZ,3)(2); a [2 0 0] [0 2 0] [0 0 2] sage: a = matrix(ZZ,1,3, [1,2,-3]); a [ 1 2 -3] sage: a = MatrixSpace(ZZ,2,4)(2); a Traceback (most recent call last): ... TypeError: nonzero scalar matrix must be square
-
BKZ
(delta=None, algorithm='fpLLL', fp=None, block_size=10, prune=0, use_givens=False, precision=0, proof=None, **kwds)¶ Block Korkin-Zolotarev reduction.
INPUT:
delta
– (default:0.99
) LLL parameteralgorithm
– (default:"fpLLL"
)"fpLLL"
or"NTL"
fp
– floating point number implementationNone
– NTL’s exact reduction or fpLLL’s wrapper (default)'fp'
– double precision: NTL’s FP or fpLLL’s double'ld'
– long doubles (fpLLL only)'qd'
– NTL’s QP
- -
'qd1'
– quad doubles: Usesquad_float
precision to compute Gram-Schmidt, but uses double precision in the search phase of the block reduction algorithm. This seems adequate for most purposes, and is faster than
'qd'
, which uses quad_float precision uniformly throughout (NTL only).
'xd'
– extended exponent: NTL’s XD or fpLLL’s dpe'rr'
– arbitrary precision: NTL’RR or fpLLL’s MPFR
block_size
– (default:10
) Specifies the size of the blocks in the reduction. High values yield shorter vectors, but the running time increases double exponentially withblock_size
.block_size
should be between 2 and the number of rows ofself
.proof
– (default: same asproof.linear_algebra()
) Insist on full BKZ reduction. If disabled and fplll is called, reduction is much faster but the result is not fully BKZ reduced.
NTL SPECIFIC INPUT:
prune
– (default:0
) The optional parameterprune
can be set to any positive number to invoke the Volume Heuristic from [SH1995]. This can significantly reduce the running time, and hence allow much bigger block size, but the quality of the reduction is of course not as good in general. Higher values ofprune
mean better quality, and slower running time. Whenprune
is0
, pruning is disabled. Recommended usage: forblock_size==30
, set10 <= prune <=15
.use_givens
– Use Givens orthogonalization. Only applies to approximate reduction using NTL. This is a bit slower, but generally much more stable, and is really the preferred orthogonalization strategy. For a nice description of this, see Chapter 5 of [GL1996].
fpLLL SPECIFIC INPUT:
precision
– (default:0
for automatic choice) bit precision to use iffp='rr'
is set**kwds
– keywords to be passed tofpylll
. Seefpylll.BKZ.Param
for details.
Also, if the verbose level is at least \(2\), some output is printed during the computation.
EXAMPLES:
sage: A = Matrix(ZZ,3,3,range(1,10)) sage: A.BKZ() [ 0 0 0] [ 2 1 0] [-1 1 3] sage: A = Matrix(ZZ,3,3,range(1,10)) sage: A.BKZ(use_givens=True) [ 0 0 0] [ 2 1 0] [-1 1 3] sage: A = Matrix(ZZ,3,3,range(1,10)) sage: A.BKZ(fp="fp") [ 0 0 0] [ 2 1 0] [-1 1 3]
ALGORITHM:
Calls either NTL or fpLLL.
-
LLL
(delta=None, eta=None, algorithm='fpLLL:wrapper', fp=None, prec=0, early_red=False, use_givens=False, use_siegel=False, **kwds)¶ Return LLL reduced or approximated LLL reduced lattice \(R\) for this matrix interpreted as a lattice.
A lattice \((b_1, b_2, ..., b_d)\) is \((\delta, \eta)\)-LLL-reduced if the two following conditions hold:
- For any \(i > j\), we have \(\lvert \mu_{i,j} \rvert \leq \eta\).
- For any \(i < d\), we have \(\delta \lvert b_i^* \rvert^2 \leq \lvert b_{i + 1}^* + \mu_{i+1, i} b_i^* \rvert^2\),
where \(\mu_{i,j} = \langle b_i, b_j^* \rangle / \langle b_j^*, b_j^* \rangle\) and \(b_i^*\) is the \(i\)-th vector of the Gram-Schmidt orthogonalisation of \((b_1, b_2, ..., b_d)\).
The default reduction parameters are \(\delta = 0.99\) and \(\eta = 0.501\). The parameters \(\delta\) and \(\eta\) must satisfy \(0.25 < \delta \leq 1.0\) and \(0.5 \leq \eta < \sqrt{\delta}\). Polynomial time complexity is only guaranteed for \(\delta < 1\). Not every algorithm admits the case \(\delta = 1\).
The lattice is returned as a matrix. Also the rank (and the determinant) of
self
are cached if those are computed during the reduction. Note that in general this only happens whenself.rank() == self.ncols()
and the exact algorithm is used.INPUT:
delta
– (default:0.99
) \(\delta\) parameter as described aboveeta
– (default:0.501
) \(\eta\) parameter as described above, ignored by NTLalgorithm
– string; one of the algorithms listed below (default:"fpLLL:wrapper"
).fp
– floating point number implementation:None
– NTL’s exact reduction or fpLLL’s wrapper'fp'
– double precision: NTL’s FP or fpLLL’s double'ld'
– long doubles (fpLLL only)'qd'
– NTL’s QP'xd'
– extended exponent: NTL’s XD or fpLLL’s dpe'rr'
– arbitrary precision: NTL’s RR or fpLLL’s MPFR
prec
– (default: auto choose) precision, ignored by NTLearly_red
– (default:False
) perform early reduction, ignored by NTLuse_givens
– (default:False
) use Givens orthogonalization. Only applies to approximate reduction using NTL. This is slower but generally more stable.use_siegel
– (default:False
) use Siegel’s condition instead of Lovász’s condition, ignored by NTL**kwds
– keywords to be passed tofpylll
. Seefpylll.LLL.reduction()
for details.
Also, if the verbose level is at least \(2\), some output is printed during the computation.
AVAILABLE ALGORITHMS:
'NTL:LLL'
- NTL’s LLL + choice offp
.'fpLLL:heuristic'
- fpLLL’s heuristic + choice offp
.'fpLLL:fast'
- fpLLL’s fast + choice offp
.'fpLLL:proved'
- fpLLL’s proved + choice offp
.'fpLLL:wrapper'
- fpLLL’s automatic choice (default).
OUTPUT:
A matrix over the integers.
EXAMPLES:
sage: A = Matrix(ZZ,3,3,range(1,10)) sage: A.LLL() [ 0 0 0] [ 2 1 0] [-1 1 3]
We compute the extended GCD of a list of integers using LLL, this example is from the Magma handbook:
sage: Q = [ 67015143, 248934363018, 109210, 25590011055, 74631449, ....: 10230248, 709487, 68965012139, 972065, 864972271 ] sage: n = len(Q) sage: S = 100 sage: X = Matrix(ZZ, n, n + 1) sage: for i in range(n): ....: X[i, i + 1] = 1 sage: for i in range(n): ....: X[i, 0] = S * Q[i] sage: L = X.LLL() sage: M = L.row(n-1).list()[1:] sage: M [-3, -1, 13, -1, -4, 2, 3, 4, 5, -1] sage: add(Q[i]*M[i] for i in range(n)) -1
The case \(\delta = 1\) is not always supported:
sage: L = X.LLL(delta=2) Traceback (most recent call last): ... TypeError: delta must be <= 1 sage: L = X.LLL(delta=1) # not tested, will eat lots of ram Traceback (most recent call last): ... RuntimeError: infinite loop in LLL sage: L = X.LLL(delta=1, algorithm='NTL:LLL') sage: L[-1] (-100, -3, -1, 13, -1, -4, 2, 3, 4, 5, -1)
Note
See
sage.libs.ntl.ntl_mat_ZZ.ntl_mat_ZZ.LLL
andfpylll.fplll.lll
for details on the algorithms used.Although LLL is a deterministic algorithm, the output for different implementations and CPUs (32-bit vs. 64-bit) may vary, while still being correct.
-
LLL_gram
(flag=0)¶ LLL reduction of the lattice whose gram matrix is
self
, assuming thatself
is positive definite.Warning
The algorithm does not check if
self
is positive definite.INPUT:
self
– a gram matrix of a positive definite quadratic formflag
– an optional flag passed toqflllgram
.According to pari:qflllgram’s documentation the options are:
0
– (default), assume thatself
has either exact- (integral or rational) or real floating point entries.
The matrix is rescaled, converted to integers and the
behavior is then as in
flag = 1
.
1
– assume that G is integral.- Computations involving Gram-Schmidt vectors are approximate, with precision varying as needed.
OUTPUT:
A dense matrix
U
over the integers that represents a unimodular transformation matrix such thatU.T * M * U
is LLL-reduced.ALGORITHM:
Calls PARI’s pari:qflllgram.
EXAMPLES:
sage: M = Matrix(ZZ, 2, 2, [5,3,3,2]) ; M [5 3] [3 2] sage: U = M.LLL_gram(); U [-1 1] [ 1 -2] sage: U.transpose() * M * U [1 0] [0 1]
The algorithm might work for some semidefinite and indefinite forms:
sage: Matrix(ZZ,2,2,[2,6,6,3]).LLL_gram() [-3 -1] [ 1 0] sage: Matrix(ZZ,2,2,[1,0,0,-1]).LLL_gram() [ 0 -1] [ 1 0]
However, it might fail for others by raising a
ValueError
:sage: Matrix(ZZ, 1,1,[0]).LLL_gram() Traceback (most recent call last): ... ValueError: qflllgram did not return a square matrix, perhaps the matrix is not positive definite sage: Matrix(ZZ, 2, 2, [0,1,1,0]).LLL_gram() Traceback (most recent call last): ... ValueError: qflllgram did not return a square matrix, perhaps the matrix is not positive definite
or by running forever:
sage: Matrix(ZZ, 2, 2, [-5, -1, -1, -5]).LLL_gram() # not tested Traceback (most recent call last): ... RuntimeError: infinite loop while calling qflllgram
-
antitranspose
()¶ Returns the antitranspose of self, without changing self.
EXAMPLES:
sage: A = matrix(2,3,range(6)) sage: type(A) <type 'sage.matrix.matrix_integer_dense.Matrix_integer_dense'> sage: A.antitranspose() [5 2] [4 1] [3 0] sage: A [0 1 2] [3 4 5] sage: A.subdivide(1,2); A [0 1|2] [---+-] [3 4|5] sage: A.antitranspose() [5|2] [-+-] [4|1] [3|0]
-
augment
(right, subdivide=False)¶ Returns a new matrix formed by appending the matrix (or vector)
right
on the right side ofself
.INPUT:
right
- a matrix, vector or free module element, whose dimensions are compatible withself
.subdivide
- default:False
- request the resulting matrix to have a new subdivision, separatingself
fromright
.
OUTPUT:
A new matrix formed by appending
right
onto the right side ofself
. Ifright
is a vector (or free module element) then in this context it is appropriate to consider it as a column vector. (The code first converts a vector to a 1-column matrix.)EXAMPLES:
sage: A = matrix(ZZ, 4, 5, range(20)) sage: B = matrix(ZZ, 4, 3, range(12)) sage: A.augment(B) [ 0 1 2 3 4 0 1 2] [ 5 6 7 8 9 3 4 5] [10 11 12 13 14 6 7 8] [15 16 17 18 19 9 10 11]
A vector may be augmented to a matrix.
sage: A = matrix(ZZ, 3, 5, range(15)) sage: v = vector(ZZ, 3, range(3)) sage: A.augment(v) [ 0 1 2 3 4 0] [ 5 6 7 8 9 1] [10 11 12 13 14 2]
The
subdivide
option will add a natural subdivision betweenself
andright
. For more details about how subdivisions are managed when augmenting, seesage.matrix.matrix1.Matrix.augment()
.sage: A = matrix(ZZ, 3, 5, range(15)) sage: B = matrix(ZZ, 3, 3, range(9)) sage: A.augment(B, subdivide=True) [ 0 1 2 3 4| 0 1 2] [ 5 6 7 8 9| 3 4 5] [10 11 12 13 14| 6 7 8]
Errors are raised if the sizes are incompatible.
sage: A = matrix(ZZ, [[1, 2],[3, 4]]) sage: B = matrix(ZZ, [[10, 20], [30, 40], [50, 60]]) sage: A.augment(B) Traceback (most recent call last): ... TypeError: number of rows must be the same, not 2 != 3
-
charpoly
(var='x', algorithm=None)¶ Note
The characteristic polynomial is defined as \(\det(xI-A)\).
INPUT:
var
- a variable namealgorithm
- (optional) either ‘generic’, ‘flint’ or ‘linbox’. Default is set to ‘linbox’.
EXAMPLES:
sage: A = matrix(ZZ,6, range(36)) sage: f = A.charpoly(); f x^6 - 105*x^5 - 630*x^4 sage: f(A) == 0 True sage: g = A.charpoly(algorithm='flint') sage: f == g True sage: n=20; A = Mat(ZZ,n)(range(n^2)) sage: A.charpoly() x^20 - 3990*x^19 - 266000*x^18 sage: A.minpoly() x^3 - 3990*x^2 - 266000*x
On non square matrices, this method raises an ArithmeticError:
sage: matrix(ZZ, 2, 1).charpoly() Traceback (most recent call last): ... ArithmeticError: only valid for square matrix
- {}”.format(
- ….: ans_flint, ans_linbox, ans_generic, m.str()))
-
decomposition
(**kwds)¶ Returns the decomposition of the free module on which this matrix A acts from the right (i.e., the action is x goes to x A), along with whether this matrix acts irreducibly on each factor. The factors are guaranteed to be sorted in the same way as the corresponding factors of the characteristic polynomial, and are saturated as ZZ modules.
INPUT:
self
- a matrix over the integers**kwds
- these are passed onto to the decomposition over QQ command.
EXAMPLES:
sage: t = ModularSymbols(11,sign=1).hecke_matrix(2) sage: w = t.change_ring(ZZ) sage: w [ 3 -2] [ 0 -2] sage: w.charpoly().factor() (x - 3) * (x + 2) sage: w.decomposition() [ (Free module of degree 2 and rank 1 over Integer Ring Echelon basis matrix: [ 5 -2], True), (Free module of degree 2 and rank 1 over Integer Ring Echelon basis matrix: [0 1], True) ]
-
determinant
(algorithm='default', proof=None, stabilize=2)¶ Return the determinant of this matrix.
INPUT:
algorithm
'default'
– useflint
'flint'
– let flint do the determinant'padic'
- uses a p-adic / multimodular algorithm that relies on code in IML and linbox'linbox'
- calls linbox det (you must set proof=False to use this!)'ntl'
- calls NTL’s det function'pari'
- uses PARI
proof
- bool or None; if None use proof.linear_algebra(); only relevant for the padic algorithm.Note
It would be VERY VERY hard for det to fail even with proof=False.
stabilize
- if proof is False, require det to be the same for this many CRT primes in a row. Ignored if proof is True.
ALGORITHM: The p-adic algorithm works by first finding a random vector v, then solving A*x = v and taking the denominator \(d\). This gives a divisor of the determinant. Then we compute \(\det(A)/d\) using a multimodular algorithm and the Hadamard bound, skipping primes that divide \(d\).
EXAMPLES:
sage: A = matrix(ZZ,8,8,[3..66]) sage: A.determinant() 0
sage: A = random_matrix(ZZ,20,20) sage: D1 = A.determinant() sage: A._clear_cache() sage: D2 = A.determinant(algorithm='ntl') sage: D1 == D2 True
We have a special-case algorithm for 4 x 4 determinants:
sage: A = matrix(ZZ,4,[1,2,3,4,4,3,2,1,0,5,0,1,9,1,2,3]) sage: A.determinant() 270
Next we try the Linbox det. Note that we must have proof=False.
sage: A = matrix(ZZ,5,[1,2,3,4,5,4,6,3,2,1,7,9,7,5,2,1,4,6,7,8,3,2,4,6,7]) sage: A.determinant(algorithm='linbox') Traceback (most recent call last): ... RuntimeError: you must pass the proof=False option to the determinant command to use LinBox's det algorithm sage: A.determinant(algorithm='linbox', proof=False) -21 sage: A._clear_cache() sage: A.determinant() -21
Try the other algorithms on the same example:
sage: A._clear_cache(); A.determinant(algorithm='padic') -21 sage: A._clear_cache(); A.determinant(algorithm='pari') -21 sage: A._clear_cache(); A.determinant(algorithm='ntl') -21 sage: A._clear_cache(); A.determinant(algorithm='padic') -21
A bigger example:
sage: A = random_matrix(ZZ,30) sage: d = A.determinant() sage: A._clear_cache() sage: A.determinant(algorithm='linbox',proof=False) == d True
-
echelon_form
(algorithm='default', proof=None, include_zero_rows=True, transformation=False, D=None)¶ Return the echelon form of this matrix over the integers, also known as the hermite normal form (HNF).
INPUT:
algorithm
– String. The algorithm to use. Valid options are:'default'
– Let Sage pick an algorithm (default). Up to 75 rows or columns with no transformation matrix, use pari with flag 0; otherwise, use flint.'flint'
- use flint'ntl'
- use NTL (only works for square matrices of full rank!)'padic'
- an asymptotically fast p-adic modular algorithm, If your matrix has large coefficients and is small, you may also want to try this.'pari'
- use PARI with flag 1'pari0'
- use PARI with flag 0'pari1'
- use PARI with flag 1'pari4'
- use PARI with flag 4 (use heuristic LLL)
proof
- (default: True); if proof=False certain determinants are computed using a randomized hybrid p-adic multimodular strategy until it stabilizes twice (instead of up to the Hadamard bound). It is incredibly unlikely that one would ever get an incorrect result with proof=False.include_zero_rows
- (default: True) if False, don’t include zero rowstransformation
- if given, also compute transformation matrix; only valid for flint and padic algorithmD
- (default: None) if given and the algorithm is ‘ntl’, then D must be a multiple of the determinant and this function will use that fact.
OUTPUT:
The Hermite normal form (=echelon form over \(\ZZ\)) of self as an immutable matrix.
EXAMPLES:
sage: A = MatrixSpace(ZZ,2)([1,2,3,4]) sage: A.echelon_form() [1 0] [0 2] sage: A = MatrixSpace(ZZ,5)(range(25)) sage: A.echelon_form() [ 5 0 -5 -10 -15] [ 0 1 2 3 4] [ 0 0 0 0 0] [ 0 0 0 0 0] [ 0 0 0 0 0]
Getting a transformation matrix in the nonsquare case:
sage: A = matrix(ZZ,5,3,[1..15]) sage: H, U = A.hermite_form(transformation=True, include_zero_rows=False) sage: H [1 2 3] [0 3 6] sage: U [ 0 0 0 4 -3] [ 0 0 0 13 -10] sage: U*A == H True
Note
If ‘ntl’ is chosen for a non square matrix this function raises a ValueError.
Special cases: 0 or 1 rows:
sage: a = matrix(ZZ, 1,2,[0,-1]) sage: a.hermite_form() [0 1] sage: a.pivots() (1,) sage: a = matrix(ZZ, 1,2,[0,0]) sage: a.hermite_form() [0 0] sage: a.pivots() () sage: a = matrix(ZZ,1,3); a [0 0 0] sage: a.echelon_form(include_zero_rows=False) [] sage: a.echelon_form(include_zero_rows=True) [0 0 0]
Illustrate using various algorithms.:
sage: matrix(ZZ,3,[1..9]).hermite_form(algorithm='pari') [1 2 3] [0 3 6] [0 0 0] sage: matrix(ZZ,3,[1..9]).hermite_form(algorithm='pari0') [1 2 3] [0 3 6] [0 0 0] sage: matrix(ZZ,3,[1..9]).hermite_form(algorithm='pari4') [1 2 3] [0 3 6] [0 0 0] sage: matrix(ZZ,3,[1..9]).hermite_form(algorithm='padic') [1 2 3] [0 3 6] [0 0 0] sage: matrix(ZZ,3,[1..9]).hermite_form(algorithm='default') [1 2 3] [0 3 6] [0 0 0]
The ‘ntl’ algorithm doesn’t work on matrices that do not have full rank.:
sage: matrix(ZZ,3,[1..9]).hermite_form(algorithm='ntl') Traceback (most recent call last): ... ValueError: ntl only computes HNF for square matrices of full rank. sage: matrix(ZZ,3,[0] +[2..9]).hermite_form(algorithm='ntl') [1 0 0] [0 1 0] [0 0 3]
-
elementary_divisors
(algorithm='pari')¶ Return the elementary divisors of self, in order.
Warning
This is MUCH faster than the
smith_form()
function.The elementary divisors are the invariants of the finite abelian group that is the cokernel of left multiplication of this matrix. They are ordered in reverse by divisibility.
INPUT:
self
- matrixalgorithm
- (default: ‘pari’)'pari'
: works robustly, but is slower.'linbox'
- use linbox (currently off, broken)
OUTPUT: list of integers
Note
These are the invariants of the cokernel of left multiplication:
sage: M = Matrix([[3,0,1],[0,1,0]]) sage: M [3 0 1] [0 1 0] sage: M.elementary_divisors() [1, 1] sage: M.transpose().elementary_divisors() [1, 1, 0]
EXAMPLES:
sage: matrix(3, range(9)).elementary_divisors() [1, 3, 0] sage: matrix(3, range(9)).elementary_divisors(algorithm='pari') [1, 3, 0] sage: C = MatrixSpace(ZZ,4)([3,4,5,6,7,3,8,10,14,5,6,7,2,2,10,9]) sage: C.elementary_divisors() [1, 1, 1, 687]
sage: M = matrix(ZZ, 3, [1,5,7, 3,6,9, 0,1,2]) sage: M.elementary_divisors() [1, 1, 6]
This returns a copy, which is safe to change:
sage: edivs = M.elementary_divisors() sage: edivs.pop() 6 sage: M.elementary_divisors() [1, 1, 6]
See also
-
frobenius
(flag=0, var='x')¶ Return the Frobenius form (rational canonical form) of this matrix.
INPUT:
flag
– 0 (default), 1 or 2 as follows:0
– (default) return the Frobenius form of this matrix.1
– return only the elementary divisor polynomials, as polynomials in var.2
– return a two-components vector [F,B] where F is the Frobenius form and B is the basis change so that \(M=B^{-1}FB\).
var
– a string (default: ‘x’)
ALGORITHM: uses PARI’s matfrobenius()
EXAMPLES:
sage: A = MatrixSpace(ZZ, 3)(range(9)) sage: A.frobenius(0) [ 0 0 0] [ 1 0 18] [ 0 1 12] sage: A.frobenius(1) [x^3 - 12*x^2 - 18*x] sage: A.frobenius(1, var='y') [y^3 - 12*y^2 - 18*y] sage: F, B = A.frobenius(2) sage: A == B^(-1)*F*B True sage: a=matrix([]) sage: a.frobenius(2) ([], []) sage: a.frobenius(0) [] sage: a.frobenius(1) [] sage: B = random_matrix(ZZ,2,3) sage: B.frobenius() Traceback (most recent call last): ... ArithmeticError: frobenius matrix of non-square matrix not defined.
AUTHORS:
- Martin Albrect (2006-04-02)
TODO: - move this to work for more general matrices than just over Z. This will require fixing how PARI polynomials are coerced to Sage polynomials.
-
gcd
()¶ Return the gcd of all entries of self; very fast.
EXAMPLES:
sage: a = matrix(ZZ,2, [6,15,-6,150]) sage: a.gcd() 3
-
height
()¶ Return the height of this matrix, i.e., the max absolute value of the entries of the matrix.
OUTPUT: A nonnegative integer.
EXAMPLES:
sage: a = Mat(ZZ,3)(range(9)) sage: a.height() 8 sage: a = Mat(ZZ,2,3)([-17,3,-389,15,-1,0]); a [ -17 3 -389] [ 15 -1 0] sage: a.height() 389
-
hermite_form
(algorithm='default', proof=None, include_zero_rows=True, transformation=False, D=None)¶ Return the echelon form of this matrix over the integers, also known as the hermite normal form (HNF).
INPUT:
algorithm
– String. The algorithm to use. Valid options are:'default'
– Let Sage pick an algorithm (default). Up to 75 rows or columns with no transformation matrix, use pari with flag 0; otherwise, use flint.'flint'
- use flint'ntl'
- use NTL (only works for square matrices of full rank!)'padic'
- an asymptotically fast p-adic modular algorithm, If your matrix has large coefficients and is small, you may also want to try this.'pari'
- use PARI with flag 1'pari0'
- use PARI with flag 0'pari1'
- use PARI with flag 1'pari4'
- use PARI with flag 4 (use heuristic LLL)
proof
- (default: True); if proof=False certain determinants are computed using a randomized hybrid p-adic multimodular strategy until it stabilizes twice (instead of up to the Hadamard bound). It is incredibly unlikely that one would ever get an incorrect result with proof=False.include_zero_rows
- (default: True) if False, don’t include zero rowstransformation
- if given, also compute transformation matrix; only valid for flint and padic algorithmD
- (default: None) if given and the algorithm is ‘ntl’, then D must be a multiple of the determinant and this function will use that fact.
OUTPUT:
The Hermite normal form (=echelon form over \(\ZZ\)) of self as an immutable matrix.
EXAMPLES:
sage: A = MatrixSpace(ZZ,2)([1,2,3,4]) sage: A.echelon_form() [1 0] [0 2] sage: A = MatrixSpace(ZZ,5)(range(25)) sage: A.echelon_form() [ 5 0 -5 -10 -15] [ 0 1 2 3 4] [ 0 0 0 0 0] [ 0 0 0 0 0] [ 0 0 0 0 0]
Getting a transformation matrix in the nonsquare case:
sage: A = matrix(ZZ,5,3,[1..15]) sage: H, U = A.hermite_form(transformation=True, include_zero_rows=False) sage: H [1 2 3] [0 3 6] sage: U [ 0 0 0 4 -3] [ 0 0 0 13 -10] sage: U*A == H True
Note
If ‘ntl’ is chosen for a non square matrix this function raises a ValueError.
Special cases: 0 or 1 rows:
sage: a = matrix(ZZ, 1,2,[0,-1]) sage: a.hermite_form() [0 1] sage: a.pivots() (1,) sage: a = matrix(ZZ, 1,2,[0,0]) sage: a.hermite_form() [0 0] sage: a.pivots() () sage: a = matrix(ZZ,1,3); a [0 0 0] sage: a.echelon_form(include_zero_rows=False) [] sage: a.echelon_form(include_zero_rows=True) [0 0 0]
Illustrate using various algorithms.:
sage: matrix(ZZ,3,[1..9]).hermite_form(algorithm='pari') [1 2 3] [0 3 6] [0 0 0] sage: matrix(ZZ,3,[1..9]).hermite_form(algorithm='pari0') [1 2 3] [0 3 6] [0 0 0] sage: matrix(ZZ,3,[1..9]).hermite_form(algorithm='pari4') [1 2 3] [0 3 6] [0 0 0] sage: matrix(ZZ,3,[1..9]).hermite_form(algorithm='padic') [1 2 3] [0 3 6] [0 0 0] sage: matrix(ZZ,3,[1..9]).hermite_form(algorithm='default') [1 2 3] [0 3 6] [0 0 0]
The ‘ntl’ algorithm doesn’t work on matrices that do not have full rank.:
sage: matrix(ZZ,3,[1..9]).hermite_form(algorithm='ntl') Traceback (most recent call last): ... ValueError: ntl only computes HNF for square matrices of full rank. sage: matrix(ZZ,3,[0] +[2..9]).hermite_form(algorithm='ntl') [1 0 0] [0 1 0] [0 0 3]
-
index_in_saturation
(proof=None)¶ Return the index of self in its saturation.
INPUT:
proof
- (default: use proof.linear_algebra()); if False, the determinant calculations are done with proof=False.
OUTPUT:
positive integer
- the index of the row span of this matrix in its saturation
ALGORITHM: Use Hermite normal form twice to find an invertible matrix whose inverse transforms a matrix with the same row span as self to its saturation, then compute the determinant of that matrix.
EXAMPLES:
sage: A = matrix(ZZ, 2,3, [1..6]); A [1 2 3] [4 5 6] sage: A.index_in_saturation() 3 sage: A.saturation() [1 2 3] [1 1 1]
-
insert_row
(index, row)¶ Create a new matrix from self with.
INPUT:
index
- integerrow
- a vector
EXAMPLES:
sage: X = matrix(ZZ,3,range(9)); X [0 1 2] [3 4 5] [6 7 8] sage: X.insert_row(1, [1,5,-10]) [ 0 1 2] [ 1 5 -10] [ 3 4 5] [ 6 7 8] sage: X.insert_row(0, [1,5,-10]) [ 1 5 -10] [ 0 1 2] [ 3 4 5] [ 6 7 8] sage: X.insert_row(3, [1,5,-10]) [ 0 1 2] [ 3 4 5] [ 6 7 8] [ 1 5 -10]
-
integer_valued_polynomials_generators
()¶ Determine the generators of the ring of integer valued polynomials on this matrix.
OUTPUT:
A pair
(mu_B, P)
whereP
is a list of polynomials in \(\QQ[X]\) such that\[\{f \in \QQ[X] \mid f(B) \in M_n(\ZZ)\} = \mu_B \QQ[X] + \sum_{g\in P} g \ZZ[X]\]where \(B\) is this matrix.
EXAMPLES:
sage: B = matrix(ZZ, [[1, 0, 1], [1, -2, -1], [10, 0, 0]]) sage: B.integer_valued_polynomials_generators() (x^3 + x^2 - 12*x - 20, [1, 1/4*x^2 + 3/4*x + 1/2])
-
inverse_of_unit
()¶ If self is a matrix with determinant \(1\) or \(-1\) return the inverse of
self
as a matrix over \(ZZ\).EXAMPLES:
sage: m = matrix(ZZ, 2, [2,1,1,1]).inverse_of_unit() sage: m [ 1 -1] [-1 2] sage: parent(m) Full MatrixSpace of 2 by 2 dense matrices over Integer Ring sage: matrix(2, [2,1,0,1]).inverse_of_unit() Traceback (most recent call last): ... ArithmeticError: non-invertible matrix
-
is_LLL_reduced
(delta=None, eta=None)¶ Return
True
if this lattice is \((\delta, \eta)\)-LLL reduced. Seeself.LLL
for a definition of LLL reduction.INPUT:
delta
– (default: \(0.99\)) parameter \(\delta\) as described aboveeta
– (default: \(0.501\)) parameter \(\eta\) as described above
EXAMPLES:
sage: A = random_matrix(ZZ, 10, 10) sage: L = A.LLL() sage: A.is_LLL_reduced() False sage: L.is_LLL_reduced() True
-
is_one
()¶ Tests whether self is the identity matrix.
EXAMPLES:
sage: matrix(2, [1,0,0,1]).is_one() True sage: matrix(2, [1,1,0,1]).is_one() False sage: matrix(2, 3, [1,0,0,0,1,0]).is_one() False
-
is_primitive
()¶ Test whether the matrix is primitive.
An integral matrix \(A\) is primitive if all its entries are non-negative and for some positive integer \(n\) the matrix \(A^n\) has all its entries positive.
EXAMPLES:
sage: m = matrix(3, [1,1,0,0,0,1,1,0,0]) sage: m.is_primitive() True sage: m**4 [3 2 1] [1 1 1] [2 1 1] sage: m = matrix(4, [[1,1,0,0],[0,0,1,0],[0,0,0,1],[1,0,0,0]]) sage: m.is_primitive() True sage: m**6 [4 3 2 1] [1 1 1 1] [2 1 1 1] [3 2 1 1] sage: m = matrix(4, [[0,1,0,1],[1,0,1,0],[0,1,0,1],[1,0,1,0]]) sage: m.is_primitive() False
Testing extremal matrices:
sage: def matrix1(d): ....: m = matrix(d) ....: m[0,0] = 1 ....: for i in range(d-1): ....: m[i,i+1] = m[i+1,i] = 1 ....: return m sage: all(matrix1(d).is_primitive() for d in range(2,20)) True sage: def matrix2(d): ....: m = matrix(d) ....: for i in range(d-1): ....: m[i,i+1] = 1 ....: m[d-1,0] = m[d-1,1] = 1 ....: return m sage: all(matrix2(d).is_primitive() for d in range(2,20)) True
Non-primitive families:
sage: def matrix3(d): ....: m = matrix(d) ....: m[0,0] = 1 ....: for i in range(d-1): ....: m[i,i+1] = 1 ....: return m sage: any(matrix3(d).is_primitive() for d in range(2,20)) False
-
minpoly
(var='x', algorithm=None)¶ INPUT:
var
- a variable namealgorithm
- (optional) either ‘linbox’ (default) or ‘generic’
EXAMPLES:
sage: A = matrix(ZZ, 6, range(36)) sage: A.minpoly() x^3 - 105*x^2 - 630*x sage: A = Mat(ZZ, 6)([k^2 for k in range(36)]) sage: A.minpoly(algorithm='linbox') x^4 - 2695*x^3 - 257964*x^2 + 1693440*x sage: A.minpoly(algorithm='generic') x^4 - 2695*x^3 - 257964*x^2 + 1693440*x
On non square matrices, this method raises an ArithmeticError:
sage: matrix(ZZ, 2, 1).minpoly() Traceback (most recent call last): ... ArithmeticError: only valid for square matrix
- {}”.format(
- ….: ans_generic, ans_linbox, m.str()))
-
null_ideal
(b=0)¶ Return the \((b)\)-ideal of this matrix.
Let \(B\) be a \(n \times n\) matrix. The null ideal modulo \(b\), or \((b)\)-ideal, is
\[N_{(b)}(B) = \{f \in \ZZ[X] \mid f(B) \in M_n(b\ZZ)\}.\]INPUT:
b
– an element of \(\ZZ\) (default: 0)
OUTPUT:
An ideal in \(\ZZ[X]\).
EXAMPLES:
sage: B = matrix(ZZ, [[1, 0, 1], [1, -2, -1], [10, 0, 0]]) sage: B.null_ideal() Principal ideal (x^3 + x^2 - 12*x - 20) of Univariate Polynomial Ring in x over Integer Ring sage: B.null_ideal(8) Ideal (8, x^3 + x^2 - 12*x - 20, 2*x^2 + 6*x + 4) of Univariate Polynomial Ring in x over Integer Ring sage: B.null_ideal(6) Ideal (6, 2*x^3 + 2*x^2 - 24*x - 40, 3*x^2 + 3*x) of Univariate Polynomial Ring in x over Integer Ring
See also
-
p_minimal_polynomials
(p, s_max=None)¶ Compute \((p^s)\)-minimal polynomials \(\nu_s\) of this matrix.
For \(s \ge 0\), a \((p^s)\)-minimal polynomial of a matrix \(B\) is a monic polynomial \(f \in \ZZ[X]\) of minimal degree such that all entries of \(f(B)\) are divisible by \(p^s\).
Compute a finite subset \(\mathcal{S}\) of the positive integers and \((p^s)\)-minimal polynomials \(\nu_s\) for \(s \in \mathcal{S}\).
For \(0 < t \le \max \mathcal{S}\), a \((p^t)\)-minimal polynomial is given by \(\nu_s\) where \(s = \min\{ r \in \mathcal{S} \mid r\ge t \}\). For \(t > \max\mathcal{S}\), the minimal polynomial of \(B\) is also a \((p^t)\)-minimal polynomial.
INPUT:
p
– a prime in \(\ZZ\)s_max
– a positive integer (default:None
); if set, only \((p^s)\)-minimal polynomials fors <= s_max
are computed (see below for details)
OUTPUT:
A dictionary. Keys are the finite set \(\mathcal{S}\), the values are the associated \((p^s)\)-minimal polynomials \(\nu_s\), \(s\in\mathcal{S}\).
Setting
s_max
only affects the output ifs_max
is at most \(\max\mathcal{S}\) where \(\mathcal{S}\) denotes the full set. In that case, only those \(\nu_s\) withs <= s_max
are returned wheres_max
is always included even if it is not included in the full set \(\mathcal{S}\).EXAMPLES:
sage: B = matrix(ZZ, [[1, 0, 1], [1, -2, -1], [10, 0, 0]]) sage: B.p_minimal_polynomials(2) {2: x^2 + 3*x + 2}
See also
-
pivots
()¶ Return the pivot column positions of this matrix.
OUTPUT: a tuple of Python integers: the position of the first nonzero entry in each row of the echelon form.
EXAMPLES:
sage: n = 3; A = matrix(ZZ,n,range(n^2)); A [0 1 2] [3 4 5] [6 7 8] sage: A.pivots() (0, 1) sage: A.echelon_form() [ 3 0 -3] [ 0 1 2] [ 0 0 0]
-
prod_of_row_sums
(cols)¶ Return the product of the sums of the entries in the submatrix of
self
with given columns.INPUT:
cols
– a list (or set) of integers representing columns ofself
OUTPUT: an integer
EXAMPLES:
sage: a = matrix(ZZ,2,3,[1..6]); a [1 2 3] [4 5 6] sage: a.prod_of_row_sums([0,2]) 40 sage: (1+3)*(4+6) 40 sage: a.prod_of_row_sums(set([0,2])) 40
-
randomize
(density=1, x=None, y=None, distribution=None, nonzero=False)¶ Randomize
density
proportion of the entries of this matrix, leaving the rest unchanged.The parameters are the same as the ones for the integer ring’s
random_element
function.If
x
andy
are given, randomized entries of this matrix have to be betweenx
andy
and have density 1.INPUT:
self
- a mutable matrix over ZZdensity
- a float between 0 and 1x, y
- if notNone
, these are passed to theZZ.random_element
function as the upper and lower endpoints in the uniform distributiondistribution
- would also be passed intoZZ.random_element
if givennonzero
- bool (default:False
); whether the new entries are guaranteed to be zero
OUTPUT:
- None, the matrix is modified in-place
EXAMPLES:
sage: A = matrix(ZZ, 2,3, [1..6]); A [1 2 3] [4 5 6] sage: A.randomize() sage: A [-8 2 0] [ 0 1 -1] sage: A.randomize(x=-30,y=30) sage: A [ 5 -19 24] [ 24 23 -9]
-
rank
(algorithm='modp')¶ Return the rank of this matrix.
INPUT:
algorithm
– either'modp'
(default) or'flint'
or'linbox'
OUTPUT:
- a nonnegative integer – the rank
Note
The rank is cached.
ALGORITHM:
If set to
'modp'
, first check if the matrix has maximum possible rank by working modulo one random prime. If not, call LinBox’s rank function.EXAMPLES:
sage: a = matrix(ZZ,2,3,[1..6]); a [1 2 3] [4 5 6] sage: a.rank() 2 sage: a = matrix(ZZ,3,3,[1..9]); a [1 2 3] [4 5 6] [7 8 9] sage: a.rank() 2
Here is a bigger example - the rank is of course still 2:
sage: a = matrix(ZZ,100,[1..100^2]); a.rank() 2
-
rational_reconstruction
(N)¶ Use rational reconstruction to lift self to a matrix over the rational numbers (if possible), where we view self as a matrix modulo N.
INPUT:
N
- an integer
OUTPUT:
matrix
- over QQ or raise a ValueError
EXAMPLES: We create a random 4x4 matrix over ZZ.
sage: A = matrix(ZZ, 4, [4, -4, 7, 1, -1, 1, -1, -12, -1, -1, 1, -1, -3, 1, 5, -1])
There isn’t a unique rational reconstruction of it:
sage: A.rational_reconstruction(11) Traceback (most recent call last): ... ValueError: rational reconstruction does not exist
We throw in a denominator and reduce the matrix modulo 389 - it does rationally reconstruct.
sage: B = (A/3 % 389).change_ring(ZZ) sage: B.rational_reconstruction(389) == A/3 True
-
saturation
(p=0, proof=None, max_dets=5)¶ Return a saturation matrix of self, which is a matrix whose rows span the saturation of the row span of self. This is not unique.
The saturation of a \(\ZZ\) module \(M\) embedded in \(\ZZ^n\) is a module \(S\) that contains \(M\) with finite index such that \(\ZZ^n/S\) is torsion free. This function takes the row span \(M\) of self, and finds another matrix of full rank with row span the saturation of \(M\).
INPUT:
p
- (default: 0); if nonzero given, saturate only at the prime \(p\), i.e., return a matrix whose row span is a \(\ZZ\)-module \(S\) that contains self and such that the index of \(S\) in its saturation is coprime to \(p\). If \(p\) is None, return full saturation of self.proof
- (default: use proof.linear_algebra()); if False, the determinant calculations are done with proof=False.max_dets
- (default: 5); technical parameter - max number of determinant to compute when bounding prime divisor of self in its saturation.
OUTPUT:
matrix
- a matrix over ZZ
Note
The result is not cached.
ALGORITHM: 1. Replace input by a matrix of full rank got from a subset of the rows. 2. Divide out any common factors from rows. 3. Check max_dets random dets of submatrices to see if their GCD (with p) is 1 - if so matrix is saturated and we’re done. 4. Finally, use that if A is a matrix of full rank, then \(hnf(transpose(A))^{-1}*A\) is a saturation of A.
EXAMPLES:
sage: A = matrix(ZZ, 3, 5, [-51, -1509, -71, -109, -593, -19, -341, 4, 86, 98, 0, -246, -11, 65, 217]) sage: A.echelon_form() [ 1 5 2262 20364 56576] [ 0 6 35653 320873 891313] [ 0 0 42993 386937 1074825] sage: S = A.saturation(); S [ -51 -1509 -71 -109 -593] [ -19 -341 4 86 98] [ 35 994 43 51 347]
Notice that the saturation spans a different module than A.
sage: S.echelon_form() [ 1 2 0 8 32] [ 0 3 0 -2 -6] [ 0 0 1 9 25] sage: V = A.row_space(); W = S.row_space() sage: V.is_submodule(W) True sage: V.index_in(W) 85986 sage: V.index_in_saturation() 85986
We illustrate each option:
sage: S = A.saturation(p=2) sage: S = A.saturation(proof=False) sage: S = A.saturation(max_dets=2)
-
smith_form
(transformation=True, integral=None)¶ Return the smith normal form of this matrix, that is the diagonal matrix \(S\) with diagonal entries the ordered elementary divisors of this matrix.
INPUT:
transformation
– a boolean (default:True
); whether to return the transformation matrices \(U\) and \(V\) such that \(S=U\cdot self\cdot V\).integral
– a subring of the base ring orTrue
(default:None
); ignored for matrices with integer entries.
Note
The
elementary_divisors()
function, which returns the diagonal entries of \(S\), is VASTLY faster than this function.The elementary divisors are the invariants of the finite abelian group that is the cokernel of this matrix. They are ordered in reverse by divisibility.
EXAMPLES:
sage: A = MatrixSpace(IntegerRing(), 3)(range(9)) sage: D, U, V = A.smith_form() sage: D [1 0 0] [0 3 0] [0 0 0] sage: U [ 0 1 0] [ 0 -1 1] [-1 2 -1] sage: V [-1 4 1] [ 1 -3 -2] [ 0 0 1] sage: U*A*V [1 0 0] [0 3 0] [0 0 0]
It also makes sense for nonsquare matrices:
sage: A = Matrix(ZZ,3,2,range(6)) sage: D, U, V = A.smith_form() sage: D [1 0] [0 2] [0 0] sage: U [ 0 1 0] [ 0 -1 1] [-1 2 -1] sage: V [-1 3] [ 1 -2] sage: U * A * V [1 0] [0 2] [0 0]
Empty matrices are handled sensibly (see trac ticket #3068):
sage: m = MatrixSpace(ZZ, 2,0)(0); d,u,v = m.smith_form(); u*m*v == d True sage: m = MatrixSpace(ZZ, 0,2)(0); d,u,v = m.smith_form(); u*m*v == d True sage: m = MatrixSpace(ZZ, 0,0)(0); d,u,v = m.smith_form(); u*m*v == d True
See also
-
symplectic_form
()¶ Find a symplectic basis for self if self is an anti-symmetric, alternating matrix.
Returns a pair (F, C) such that the rows of C form a symplectic basis for self and F = C * self * C.transpose().
Raises a ValueError if self is not anti-symmetric, or self is not alternating.
Anti-symmetric means that \(M = -M^t\). Alternating means that the diagonal of \(M\) is identically zero.
A symplectic basis is a basis of the form \(e_1, \ldots, e_j, f_1, \ldots f_j, z_1, \dots, z_k\) such that
- \(z_i M v^t\) = 0 for all vectors \(v\)
- \(e_i M {e_j}^t = 0\) for all \(i, j\)
- \(f_i M {f_j}^t = 0\) for all \(i, j\)
- \(e_i M {f_i}^t = 1\) for all \(i\)
- \(e_i M {f_j}^t = 0\) for all \(i\) not equal
- \(j\).
The ordering for the factors \(d_{i} | d_{i+1}\) and for the placement of zeroes was chosen to agree with the output of
smith_form()
.See the example for a pictorial description of such a basis.
EXAMPLES:
sage: E = matrix(ZZ, 5, 5, [0, 14, 0, -8, -2, -14, 0, -3, -11, 4, 0, 3, 0, 0, 0, 8, 11, 0, 0, 8, 2, -4, 0, -8, 0]); E [ 0 14 0 -8 -2] [-14 0 -3 -11 4] [ 0 3 0 0 0] [ 8 11 0 0 8] [ 2 -4 0 -8 0] sage: F, C = E.symplectic_form() sage: F [ 0 0 1 0 0] [ 0 0 0 2 0] [-1 0 0 0 0] [ 0 -2 0 0 0] [ 0 0 0 0 0] sage: F == C * E * C.transpose() True sage: E.smith_form()[0] [1 0 0 0 0] [0 1 0 0 0] [0 0 2 0 0] [0 0 0 2 0] [0 0 0 0 0]
-
transpose
()¶ Returns the transpose of self, without changing self.
EXAMPLES:
We create a matrix, compute its transpose, and note that the original matrix is not changed.
sage: A = matrix(ZZ, 2, 3, range(6)) sage: type(A) <type 'sage.matrix.matrix_integer_dense.Matrix_integer_dense'> sage: B = A.transpose() sage: print(B) [0 3] [1 4] [2 5] sage: print(A) [0 1 2] [3 4 5]
.T
is a convenient shortcut for the transpose:sage: A.T [0 3] [1 4] [2 5]
sage: A.subdivide(None, 1); A [0|1 2] [3|4 5] sage: A.transpose() [0 3] [---] [1 4] [2 5]
-