Sparse matrices over \ZZ/n\ZZ for n small.

This is a compiled implementation of sparse matrices over \ZZ/n\ZZ for n small.

TODO: - move vectors into a pyrex vector class - add _add_ and _mul_ methods.

EXAMPLES:

sage: a = matrix(Integers(37),3,3,range(9),sparse=True); a
[0 1 2]
[3 4 5]
[6 7 8]
sage: type(a)
<type 'sage.matrix.matrix_modn_sparse.Matrix_modn_sparse'>
sage: parent(a)
Full MatrixSpace of 3 by 3 sparse matrices over Ring of integers modulo 37
sage: a^2
[15 18 21]
[ 5 17 29]
[32 16  0]
sage: a+a
[ 0  2  4]
[ 6  8 10]
[12 14 16]
sage: b = a.new_matrix(2,3,range(6)); b
[0 1 2]
[3 4 5]
sage: a*b
...
TypeError: unsupported operand parent(s) for '*': 'Full MatrixSpace of 3 by 3 sparse matrices over Ring of integers modulo 37' and 'Full MatrixSpace of 2 by 3 sparse matrices over Ring of integers modulo 37'
sage: b*a
[15 18 21]
[ 5 17 29]
sage: a == loads(dumps(a))
True
sage: b == loads(dumps(b))
True
sage: a.echelonize(); a
[ 1  0 36]
[ 0  1  2]
[ 0  0  0]
sage: b.echelonize(); b
[ 1  0 36]
[ 0  1  2]
sage: a.pivots()
[0, 1]
sage: b.pivots()
[0, 1]
sage: a.rank()
2
sage: b.rank()
2
sage: a[2,2] = 5
sage: a.rank()
3
TESTS:
sage: matrix(Integers(37),0,0,sparse=True).inverse() []
class sage.matrix.matrix_modn_sparse.Matrix_modn_sparse
__eq__()
x.__eq__(y) <==> x==y
__ge__()
x.__ge__(y) <==> x>=y
__gt__()
x.__gt__(y) <==> x>y
__hash__()
x.__hash__() <==> hash(x)
__init__()
x.__init__(...) initializes x; see x.__class__.__doc__ for signature
__le__()
x.__le__(y) <==> x<=y
__lt__()
x.__lt__(y) <==> x<y
__ne__()
x.__ne__(y) <==> x!=y
static __new__()
T.__new__(S, ...) -> a new object with type S, a subtype of T
_dict()

Unsafe version of the dict method, mainly for internal use. This may return the dict of elements, but as an unsafe reference to the underlying dict of the object. It might be dangerous if you change entries of the returned dict.

EXAMPLES:
sage: MS = MatrixSpace(GF(13), 50, 50, sparse=True) sage: m = MS.random_element(density=0.002) sage: m._dict() {(7, 43): 11, (29, 44): 10, (35, 4): 11, (16, 26): 2}
TESTS:
sage: parent(m._dict()[7,43]) Finite Field of size 13
_echelon_in_place_classical()

Replace self by its reduction to reduced row echelon form.

ALGORITHM: We use Gauss elimination, in a slightly intelligent way, in that we clear each column using a row with the minimum number of nonzero entries.

TODO: Implement switching to a dense method when the matrix gets dense.

_matrix_times_matrix_dense()

Multiply self by the sparse matrix _right, and return the result as a dense matrix.

EXAMPLES:

sage: a = matrix(GF(10007), 2, [1,2,3,4], sparse=True) sage: b = matrix(GF(10007), 2, 3, [1..6], sparse=True) sage: a * b [ 9 12 15] [19 26 33] sage: c = a._matrix_times_matrix_dense(b); c [ 9 12 15] [19 26 33] sage: type(c) <type ‘sage.matrix.matrix_modn_dense.Matrix_modn_dense’>

sage: a = matrix(GF(2), 20, 20, sparse=True) sage: a*a == a._matrix_times_matrix_dense(a) True sage: type(a._matrix_times_matrix_dense(a)) <type ‘sage.matrix.matrix_mod2_dense.Matrix_mod2_dense’>

_nonzero_positions_by_row()

Returns the list of pairs (i,j) such that self[i,j] != 0.

It is safe to change the resulting list (unless you give the option copy=False).

EXAMPLE::
sage: M = Matrix(GF(7), [[0,0,0,1,0,0,0,0],[0,1,0,0,0,0,1,0]], sparse=True); M [0 0 0 1 0 0 0 0] [0 1 0 0 0 0 1 0] sage: M.nonzero_positions() [(0, 3), (1, 1), (1, 6)]
_pickle()

TESTS:

sage: M = Matrix( GF(2), [[1,1,1,1,0,0,0,0,0,0]], sparse=True )
sage: loads(dumps(M))
[1 1 1 1 0 0 0 0 0 0]
sage: loads(dumps(M)) == M
True
_rank_linbox()
See self.rank().
_solve_right_nonsingular_square()

If self is a matrix A, 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.

Note

In Sage one can also write A  B for A.solve_right(B), i.e., Sage implements the “the MATLAB/Octave backslash operator”.

INPUT:

  • B - a matrix or vector
  • algorithm - one of the following:
  • 'LinBox:BlasElimination' - dense elimination
  • 'LinBox:Blackbox' - LinBox chooses a Blackbox algorithm
  • 'LinBox:Wiedemann' - Wiedemann’s algorithm
  • 'generic' - use generic implementation (inversion)
  • None - LinBox chooses an algorithm (default)
  • check_rank - check rank before attempting to solve (default: True)

OUTPUT: a matrix or vector

EXAMPLES:

sage: A = matrix(GF(127), 3, [1,2,3,-1,2,5,2,3,1], sparse=True)
sage: b = vector(GF(127),[1,2,3])
sage: x = A \ b; x
(73, 76, 10) 
sage: A * x
(1, 2, 3)
_unpickle()
density()

Return the density of self, i.e., the ratio of the number of nonzero entries of self to the total size of self.

EXAMPLES:

sage: A = matrix(QQ,3,3,[0,1,2,3,0,0,6,7,8],sparse=True)
sage: A.density()
2/3

Notice that the density parameter does not ensure the density of a matrix; it is only an upper bound.

sage: A = random_matrix(GF(127),200,200,density=0.3, sparse=True)
sage: A.density()
257/1000
lift()

Return lift of this matrix to a sparse matrix over the integers.

EXAMPLES:
sage: a = matrix(GF(7),2,3,[1..6], sparse=True) sage: a.lift() [1 2 3] [4 5 6] sage: a.lift().parent() Full MatrixSpace of 2 by 3 sparse matrices over Integer Ring

Subdivisions are preserved when lifting:

sage: a.subdivide([], [1,1]); a
[1||2 3]
[4||5 6]
sage: a.lift()
[1||2 3]
[4||5 6]
matrix_from_columns()

Return the matrix constructed from self using columns with indices in the columns list.

EXAMPLES:

sage: M = MatrixSpace(GF(127),3,3,sparse=True)
sage: A = M(range(9)); A
[0 1 2]
[3 4 5]
[6 7 8]
sage: A.matrix_from_columns([2,1])
[2 1]
[5 4]
[8 7]
matrix_from_rows()

Return the matrix constructed from self using rows with indices in the rows list.

INPUT:

  • rows - list or tuple of row indices

EXAMPLE:

sage: M = MatrixSpace(GF(127),3,3,sparse=True)
sage: A = M(range(9)); A
[0 1 2]
[3 4 5]
[6 7 8]
sage: A.matrix_from_rows([2,1])
[6 7 8]
[3 4 5]
p
rank()

Compute the rank of self.

INPUT:

  • gauss - if True LinBox’ Gaussian elimination is used. If False ‘Symbolic Reordering’ as implemented in LinBox is used. If ‘native’ the native Sage implementation is used. (default: False)

EXAMPLE:

sage: A = random_matrix(GF(127),200,200,density=0.01,sparse=True)
sage: r1 = A.rank(gauss=False)
sage: r2 = A.rank(gauss=True)
sage: r3 = A.rank(gauss='native')
sage: r1 == r2 == r3
True
sage: r1
155

ALGORITHM: Uses LinBox or native implementation.

REFERENCES:

Note

For very sparse matrices Gaussian elimination is faster because it barly has anything to do. If the fill in needs to be considered, ‘Symbolic Reordering’ is usually much faster.

swap_rows()
transpose()

Return the transpose of self.

EXAMPLE:

sage: A = matrix(GF(127),3,3,[0,1,0,2,0,0,3,0,0],sparse=True)
sage: A
[0 1 0]
[2 0 0]
[3 0 0]
sage: A.transpose()
[0 2 3]
[1 0 0]
[0 0 0]
visualize_structure()

Write a PNG image to ‘filename’ which visualizes self by putting black pixels in those positions which have nonzero entries.

White pixels are put at positions with zero entries. If ‘maxsize’ is given, then the maximal dimension in either x or y direction is set to ‘maxsize’ depending on which is bigger. If the image is scaled, the darkness of the pixel reflects how many of the represented entries are nonzero. So if e.g. one image pixel actually represents a 2x2 submatrix, the dot is darker the more of the four values are nonzero.

INPUT:

  • filename - either a path or None in which case a filename in the current directory is chosen automatically (default:None)
  • maxsize - maximal dimension in either x or y direction of the resulting image. If None or a maxsize larger than max(self.nrows(),self.ncols()) is given the image will have the same pixelsize as the matrix dimensions (default: 512)

Previous topic

Dense matrices over \ZZ/n\ZZ for n small.

Next topic

Dense matrices over the integer ring.

This Page