Up One Level SAGE Reference Manual Contents Module Index Index
Previous: Up: 31. Matrices and Spaces Next:

31.16 Dense matrices over the integer ring

Module: sage.matrix.matrix_integer_dense

File: sage/matrix/matrix_integer_dense.pyx (starting at line 1)

Dense matrices over the integer ring.

Author Log: 

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

TESTS: 

sage: a = matrix(ZZ,2,range(4), sparse=False)
sage: loads(dumps(a)) == a
True

Module-level Functions

clear_mpz_globals( )

convert_parimatrix( )

gmp_randrange( )

init_mpz_globals( )

tune_multiplication( )

Class: Matrix_integer_dense

class Matrix_integer_dense
File: sage/matrix/matrix_integer_dense.pyx (starting at line 112)

Matrix over the integers.

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. 

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

Functions: charpoly,$ $ decomposition,$ $ determinant,$ $ echelon_form,$ $ elementary_divisors,$ $ frobenius,$ $ height,$ $ kernel,$ $ LLL,$ $ LLL_gram,$ $ minpoly,$ $ pivots,$ $ prod_of_row_sums,$ $ randomize,$ $ rank,$ $ rational_reconstruction,$ $ smith_form,$ $ solve_right

charpoly( )

INPUT:
    var -- a variable name
    algorithm -- 'linbox' (default)
                 'generic'

NOTE: Linbox charpoly disabled on 64-bit machines, since it hangs in many cases. 

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: n=20; A = Mat(ZZ,n)(range(n^2))
sage: A.charpoly()
x^20 - 3990*x^19 - 266000*x^18
sage: A.minpoly()                # optional -- os x only right now
x^3 - 3990*x^2 - 266000*x

decomposition( )

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.


sage: t = ModularSymbols(11,sign=1).hecke_matrix(2)
sage: w = t.change_ring(ZZ)
sage: w.list()
[3, -1, 0, -2]

determinant( )

Return the determinant of this matrix. 

INPUT:
    algorithm -- 'linbox', 'ntl' or None (default: None)

ALGORITHM: Uses LinBox or NTL. 

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

echelon_form( )

Return the echelon form of this matrix over the integers also known as the hermit normal form (HNF). 

INPUT:
    algorithm -- 'pari', 'ntl' or 'default';
          the default is ntl if self is square and of full rank;
          otherwise, it is ntl.   
    cutoff -- ignored currently
    include_zero_rows -- (default: True)
                         if False, don't include zero rows


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]

TESTS: Make sure the zero matrices are handled correctly: 

sage: m = matrix(ZZ,3,3,[0]*9)
sage: m.echelon_form() 
[0 0 0]
[0 0 0]
[0 0 0]
sage: m = matrix(ZZ,3,1,[0]*3)
sage: m.echelon_form()
[0]
[0]
[0]
sage: m = matrix(ZZ,1,3,[0]*3)
sage: m.echelon_form()
[0 0 0]

The ultimate border case! 

sage: m = matrix(ZZ,0,0,[])        
sage: m.echelon_form()
[]

NOTE: If 'ntl' is chosen for a non square matrix we silently fall back to the 'pari' implementation. Also, if the matrix's rank is not sufficient for 'ntl' we silently fall back to 'pari'.

elementary_divisors( )

Return the elementary divisors of self, in order.

IMPLEMENTATION: Uses linbox, except sometimes linbox doesn't work (errors about pre-conditioning), in which case PARI is used.

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 this matrix. They are ordered in reverse by divisibility. 

INPUT:
    self -- matrix
    algorithm -- (default: 'pari')
         'pari': works robustless, but is slower.
         'linbox' -- use linbox (currently off, broken)
    
OUTPUT:
    list of int's


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: smith_form

frobenius( )

Return the Frobenius form (rational canonical form) of this matrix. 

INPUT:
    flag --an integer:
        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')


INPUT:
   flag -- 0 (default), 1 or 2 as described above

ALGORITHM: uses pari's matfrobenius() 

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: A.frobenius(2)
([ 0  0  0]
[ 1  0 18]
[ 0  1 12],
[    -1      2     -1]
[     0  23/15 -14/15]
[     0  -2/15   1/15])

Author: 2006-04-02: Martin Albrecht

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.

height( )

Return the height of this matrix, i.e., the max absolute value of the entries of the matrix.

OUTPUT: A nonnegative integer. 

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

kernel( )

Return the kernel of this matrix, as a module over the integers. 

INPUT:
   LLL -- bool (default: False); if True the basis is an LLL
          reduced basis; otherwise, it is an echelon basis.

By convention if self has 0 rows, the kernel is of dimension 0, whereas the kernel is whole domain if self has 0 columns. 

sage: M = MatrixSpace(ZZ,4,2)(range(8))
sage: M.kernel()
Free module of degree 4 and rank 2 over Integer Ring
Echelon basis matrix:
[ 1  0 -3  2]
[ 0  1 -2  1]

LLL( )

Returns LLL reduced or approximated LLL reduced lattice R for self.

The lattice is returned as a matrix. Also the rank (and the determinant) of self are cached if those are computed during the reduction.

More specifically, elementary row transformations are performed on a copy of self so that the non-zero rows of R form an LLL-reduced basis for the lattice spanned by the rows of self. The default reduction parameters are $ \delta=3/4$ and $ eta=0.501$.

For a basis reduced with parameter $ \delta$, the squared length of the first non-zero basis vector is no more than $ 1/(delta-1/4)^{r-1}$ times that of the shortest vector in the lattice.

If we can compute the determinant of self using this method, we also cache it. Note that in general this only happens when self.rank() == self.ncols() and the exact algorithm is used. 

INPUT:
    delta -- parameter as described above (default: 3/4)
    eta -- parameter as described above (default: 0.501), ignored
           by NTL
    algorithm -- string, one of the algorithms mentioned below
                or None (default: None)
    fp -- None -- NTL's exact reduction or fpLLL's wrapper
       'fp' -- double precision: NTL's FP or fpLLL's double 
       'qd' -- quad doubles: NTL's QP
       'xd' -- extended exponent: NTL's XD or fpLLL's dpe
       'rr' -- arbitrary precision: NTL'RR or fpLLL's MPFR
    prec -- precision, ignored by NTL (default: auto choose)
    early_red -- perform early reduction, ignored by NTL (default: False)
    use_givens -- use Givens orthogonalization (default: False)
                  only applicable to approximate reductions and NTL.
                  This is more stable but slower.
                  

Also, if the verbose level is >= 2, some more verbose output
is printed during the calculation if NTL is used.

AVAILABLE ALGORITHMS:
    NTL:LLL -- NTL's LLL + fp
    fpLLL:heuristic -- fpLLL's heuristic + fp
    fpLLL:fast -- fpLLL's fast
    fpLLL:wrapper -- fpLLL's automatic choice (default)

OUTPUT:
    a matrix over the integers


sage: A = Matrix(ZZ,3,3,range(1,10))
sage: A.LLL()
[ 0  0  0]
[ 2  1  0]
[-1  1  3]

ALGORITHM: Uses NTL or fpLLL.

REFERENCES: ntl.mat_ZZ or sage.libs.fplll.fplll for details on the used algorithms.

LLL_gram( )

LLL reduction of the lattice whose gram matrix is self. 

INPUT:
    M -- gram matrix of a definite quadratic form

OUTPUT:
    U -- unimodular transformation matrix such that

U.transpose() * M * U

is LLL-reduced

ALGORITHM: Use PARI 

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]

Semidefinite and indefinite forms raise a ValueError: 

sage: Matrix(ZZ,2,2,[2,6,6,3]).LLL_gram()
Traceback (most recent call last):
...
ValueError: not a definite matrix
sage: Matrix(ZZ,2,2,[1,0,0,-1]).LLL_gram()
Traceback (most recent call last):
...
ValueError: not a definite matrix

BUGS: should work for semidefinite forms (PARI is ok)

minpoly( )

INPUT:
    var -- a variable name
    algorithm -- 'linbox' (default)
                 'generic'

NOTE: Linbox charpoly disabled on 64-bit machines, since it hangs in many cases. 

sage: A = matrix(ZZ,6, range(36))
sage: A.minpoly()           # optional -- os x only right now
x^3 - 105*x^2 - 630*x
sage: n=6; A = Mat(ZZ,n)([k^2 for k in range(n^2)])
sage: A.minpoly()           # optional -- os x only right now
x^4 - 2695*x^3 - 257964*x^2 + 1693440*x

pivots( )

Return the pivot column positions of this matrix as a list of Python integers.

This returns a list, of the position of the first nonzero entry in each row of the echelon form.

OUTPUT: list - a list of Python ints 

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]

randomize( )

Randomize density proportion of the entries of this matrix, leaving the rest unchanged.

The parameters are the same as the integer ring's random_element function.

If x and y are given, randomized entries of this matrix to be between x and y and have density 1.

rational_reconstruction( )

Use rational reconstruction to lift self to a matrix over the rational numbers (if possible), where we view self as a matrix modulo N.

smith_form( )

Returns matrices S, U, and V such that S = U*self*V, and S is in Smith normal form. Thus S is diagonal with diagonal entries the ordered elementary divisors of S.

WARNING: 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. 

sage: A = MatrixSpace(IntegerRing(), 3)(range(9))
sage: D, U, V = A.smith_form()
sage: D
[0 0 0]
[0 3 0]
[0 0 1]
sage: U
[-1  2 -1]
[ 0 -1  1]
[ 0  1  0]
sage: V
[ 1  4 -1]
[-2 -3  1]
[ 1  0  0]
sage: U*A*V
[0 0 0]
[0 3 0]
[0 0 1]

It also makes sense for nonsquare matrices: 

sage: A = Matrix(ZZ,3,2,range(6))
sage: D, U, V = A.smith_form()
sage: D
[0 0]
[2 0]
[0 1]
sage: U
[-1  2 -1]
[ 0 -1  1]
[ 0  1  0]
sage: V
[ 3 -1]
[-2  1]
sage: U * A * V
[0 0]
[2 0]
[0 1]

SEE ALSO: elementary_divisors

solve_right( )

If self is a matrix $ A$ of full rank, then this function returns a vector or matrix $ X$ such that $ A X = B$. If $ B$ is a vector then $ X$ is a vector and if $ B$ is a matrix, then $ X$ is a matrix. The base ring of $ X$ is the integers unless a denominator is needed in which case the base ring is the rational numbers.

NOTE: In SAGE one can also write A B for A.solve_right(B), i.e., SAGE implements the ``the MATLAB/Octave backslash operator''.

NOTE: This is currently only implemented when A is square. 

INPUT:
    B -- a matrix or vector
OUTPUT:
    a matrix or vector


sage: a = matrix(ZZ, 2, [0, -1, 1, 0])
sage: v = vector(ZZ, [2, 3])
sage: a \ v
(3, -2)
sage: parent(a\v)
Ambient free module of rank 2 over the principal ideal domain Integer Ring

We solve a bigger system where the answer is over the rationals. 

sage: a = matrix(ZZ, 3, 3, [1,2,3,4, 5, 6, 8, -2, 3])
sage: v = vector(ZZ, [1,2,3])
sage: w = a \ v; w
(2/15, -4/15, 7/15)
sage: parent(w)
Vector space of dimension 3 over Rational Field
sage: a * w
(1, 2, 3)

We solve a system where the right hand matrix has multiple columns. 

sage: a = matrix(ZZ, 3, 3, [1,2,3,4, 5, 6, 8, -2, 3])
sage: b = matrix(ZZ, 3, 2, [1,5, 2, -3, 3, 0])
sage: w = a \ b; w
[ 2/15 -19/5]
[-4/15 -27/5]
[ 7/15 98/15]
sage: a * w
[ 1  5]
[ 2 -3]
[ 3  0]

TESTS: We create a random 100x100 matrix and solve the corresponding system, then verify that the result is correct. (NOTE: This test is very risky without having a seeded random number generator!) 

sage: n = 100
sage: a = random_matrix(ZZ,n)
sage: v = vector(ZZ,n,range(n))
sage: x = a \ v
sage: a * x == v
True

Special Functions: __copy__,$ $ __eq__,$ $ __ge__,$ $ __gt__,$ $ __le__,$ $ __lt__,$ $ __ne__,$ $ _adjoint,$ $ _charpoly_linbox,$ $ _det_linbox,$ $ _det_ntl,$ $ _echelon_in_place_classical,$ $ _echelon_strassen,$ $ _elementary_divisors_linbox,$ $ _invert_iml,$ $ _linbox_dense,$ $ _linbox_sparse,$ $ _minpoly_linbox,$ $ _mod_int,$ $ _multiply_classical,$ $ _multiply_linbox,$ $ _multiply_multi_modular,$ $ _ntl_,$ $ _pickle,$ $ _poly_linbox,$ $ _rank_linbox,$ $ _rank_modp,$ $ _rational_echelon_via_solve,$ $ _rational_kernel_iml,$ $ _reduce,$ $ _solve_iml,$ $ _unpickle

__copy__( )

_adjoint( )

Return the adjoint of this matrix.

Assumes self is a square matrix (checked in adjoint).

_charpoly_linbox( )

_det_linbox( )

Compute the determinant of this matrix using Linbox.

_det_ntl( )

Compute the determinant of this matrix using NTL.

_echelon_in_place_classical( )

_echelon_strassen( )

_elementary_divisors_linbox( )

_invert_iml( )

Invert this matrix using IML. The output matrix is an integer matrix and a denominator. 

INPUT:
   self -- an invertible matrix
   use_nullspace -- (default: False): whether to use nullspace
             algorithm, which is slower, but doesn't require
             checking that the matrix is invertible as a
             precondition.
   check_invertible -- (default: True) whether to check that
             the matrix is invertible.

OUTPUT: A, d such that A*self = d
   A -- a matrix over ZZ
   d -- an integer

ALGORITHM: Uses IML's p-adic nullspace function. 

sage: a = matrix(ZZ,3,[1,2,5, 3,7,8, 2,2,1])
sage: b, d = a._invert_iml(); b,d
([  9  -8  19]
[-13   9  -7]
[  8  -2  -1], 23)
sage: a*b
[23  0  0]
[ 0 23  0]
[ 0  0 23]

Author: William Stein

_linbox_dense( )

_linbox_sparse( )

_minpoly_linbox( )

_mod_int( )

_multiply_classical( )

sage: n = 3
sage: a = MatrixSpace(ZZ,n,n)(range(n^2))
sage: a*a
[ 15  18  21]
[ 42  54  66]
[ 69  90 111]

_multiply_linbox( )

Multiply matrices over ZZ using linbox.

WARNING: This is very slow right now, i.e., linbox is very slow. 

sage: A = matrix(ZZ,2,3,range(6))
sage: A*A.transpose()
[ 5 14]
[14 50]
sage: A._multiply_linbox(A.transpose())
[ 5 14]
[14 50]

_multiply_multi_modular( )

_ntl_( )

ntl.mat_ZZ representation of self. 

sage: a = MatrixSpace(ZZ,200).random_element(x=-2, y=2)    # -2 to 2
sage: A = a._ntl_()

Note: NTL only knows dense matrices, so if you provide a sparse matrix NTL will allocate memory for every zero entry.

_pickle( )

_poly_linbox( )

INPUT:
    var -- 'x'
    typ -- 'minpoly' or 'charpoly'

_rank_linbox( )

Compute the rank of this matrix using Linbox.

_rank_modp( )

_rational_echelon_via_solve( )

Computes information that gives the reduced row echelon form (over QQ!) of a matrix with integer entries. 

INPUT:
    self -- a matrix over the integers.
    
OUTPUT:
    pivots -- ordered list of integers that give the pivot column positions
    nonpivots -- ordered list of the nonpivot column positions
    X -- matrix with integer entries
    d -- integer

If you put standard basis vectors in order at the pivot columns, and put the matrix (1/d)*X everywhere else, then you get the reduced row echelon form of self, without zero rows at the bottom.

NOTE: IML is the actual underlying $ p$-adic solver that we use.

Author: William Stein

ALGORITHM: I came up with this algorithm from scratch. As far as I know it is new. It's pretty simple, but it is ... (fast?!).

Let A be the input matrix.

1
Compute r = rank(A).

2
Compute the pivot columns of the transpose $ A^t$ of $ A$. One way to do this is to choose a random prime $ p$ and compute the row echelon form of $ A^t$ modulo $ p$ (if the number of pivot columns is not equal to $ r$, pick another prime).

3
Let $ B$ be the submatrix of $ A$ consisting of the rows corresponding to the pivot columns found in the previous step. Note that, aside from zero rows at the bottom, $ B$ and $ A$ have the same reduced row echelon form.

4
Compute the pivot columns of $ B$, again possibly by choosing a random prime $ p$ as in [2] and computing the Echelon form modulo $ p$. If the number of pivot columns is not $ r$, repeat with a different prime. Note - in this step it is possible that we mistakenly choose a bad prime $ p$ such that there are the right number of pivot columns modulo $ p$, but they are at the wrong positions - e.g., imagine the augmented matrix $ [pI\vert I]$ - modulo $ p$ we would miss all the pivot columns. This is OK, since in the next step we would detect this, as the matrix we would obtain would not be in echelon form.

5
Let $ C$ be the submatrix of $ B$ of pivot columns. Let $ D$ be the complementary submatrix of $ B$ of all all non-pivot columns. Use a $ p$-adic solver to find the matrix $ X$ and integer $ d$ such that $ C (1/d) X=D$. I.e., solve a bunch of linear systems of the form $ Cx = v$, where the columns of $ X$ are the solutions.

6
Verify that we had chosen the correct pivot columns. Inspect the matrix $ X$ obtained in step 5. If when used to construct the echelon form of $ B$, $ X$ indeed gives a matrix in reduced row echelon form, then we are done - output the pivot columns, $ X$, and $ d$. To tell if $ X$ is correct, do the following: For each column of $ X$ (corresponding to non-pivot column $ i$ of $ B$), read up the column of $ X$ until finding the first nonzero position $ j$; then verify that $ j$ is strictly less than the number of pivot columns of $ B$ that are strictly less than $ i$. Otherwise, we got the pivots of $ B$ wrong - try again starting at step 4, but with a different random prime.

_rational_kernel_iml( )

IML: Return the rational kernel of this matrix (acting from the left), considered as a matrix over QQ. I.e., returns a matrix K such that self*K = 0, and the number of columns of K equals the nullity of self.

Author: William Stein

_reduce( )

_solve_iml( )

Let A equal self be a square matrix. Given B return an integer matrix C and an integer d such that self C*A == d*B if right is False or A*C == d*B if right is True.

OUTPUT: C - integer matrix d - integer denominator 

sage: A = matrix(ZZ,4,4,[0, 1, -2, -1, -1, 1, 0, 2, 2, 2, 2, -1, 0, 2, 2, 1])
sage: B = matrix(ZZ,3,4, [-1, 1, 1, 0, 2, 0, -2, -1, 0, -2, -2, -2])
sage: C,d = A._solve_iml(B,right=False); C
[  6 -18 -15  27]
[  0  24  24 -36]
[  4 -12  -6  -2]


sage: d
12


sage: C*A == d*B
True


sage: A = matrix(ZZ,4,4,[0, 1, -2, -1, -1, 1, 0, 2, 2, 2, 2, -1, 0, 2, 2, 1])
sage: B = matrix(ZZ,4,3, [-1, 1, 1, 0, 2, 0, -2, -1, 0, -2, -2, -2])
sage: C,d = A._solve_iml(B)
sage: C
[ 12  40  28]
[-12  -4  -4]
[ -6 -25 -16]
[ 12  34  16]


sage: d
12


sage: A*C == d*B
True

ALGORITHM: Uses IML.

Author: Martin Albrecht

_unpickle( )

Up One Level SAGE Reference Manual Contents Module Index Index
Previous: Up: 31. Matrices and Spaces Next:
Release 2007.10.28, documentation updated on October 28, 2007.
See About this document... for information on suggesting changes.