Hot-keys on this page
r m x p toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
""" Dense matrices over GF(2) using the M4RI library.
AUTHOR: Martin Albrecht <malb@informatik.uni-bremen.de>
EXAMPLES::
sage: a = matrix(GF(2),3,range(9),sparse=False); a [0 1 0] [1 0 1] [0 1 0] sage: a.rank() 2 sage: type(a) <type 'sage.matrix.matrix_mod2_dense.Matrix_mod2_dense'> sage: a[0,0] = 1 sage: a.rank() 3 sage: parent(a) Full MatrixSpace of 3 by 3 dense matrices over Finite Field of size 2
sage: a^2 [0 1 1] [1 0 0] [1 0 1] sage: a+a [0 0 0] [0 0 0] [0 0 0]
sage: b = a.new_matrix(2,3,range(6)); b [0 1 0] [1 0 1]
sage: a*b Traceback (most recent call last): ... TypeError: unsupported operand parent(s) for *: 'Full MatrixSpace of 3 by 3 dense matrices over Finite Field of size 2' and 'Full MatrixSpace of 2 by 3 dense matrices over Finite Field of size 2' sage: b*a [1 0 1] [1 0 0]
sage: TestSuite(a).run() sage: TestSuite(b).run()
sage: a.echelonize(); a [1 0 0] [0 1 0] [0 0 1] sage: b.echelonize(); b [1 0 1] [0 1 0]
TESTS::
sage: FF = FiniteField(2) sage: V = VectorSpace(FF,2) sage: v = V([0,1]); v (0, 1) sage: W = V.subspace([v]) sage: W Vector space of degree 2 and dimension 1 over Finite Field of size 2 Basis matrix: [0 1] sage: v in W True
sage: M = Matrix(GF(2), [[1,1,0],[0,1,0]]) sage: M.row_space() Vector space of degree 3 and dimension 2 over Finite Field of size 2 Basis matrix: [1 0 0] [0 1 0]
sage: M = Matrix(GF(2), [[1,1,0],[0,0,1]]) sage: M.row_space() Vector space of degree 3 and dimension 2 over Finite Field of size 2 Basis matrix: [1 1 0] [0 0 1]
.. TODO::
- make LinBox frontend and use it
- charpoly ? - minpoly ?
- make Matrix_modn_frontend and use it (?) """
#***************************************************************************** # Copyright (C) 2004,2005,2006 William Stein <wstein@gmail.com> # Copyright (C) 2007,2008,2009 Martin Albrecht <M.R.Albrecht@rhul.ac.uk> # # This program 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. # http://www.gnu.org/licenses/ #***************************************************************************** from __future__ import absolute_import
from cysignals.memory cimport check_malloc, sig_free from cysignals.signals cimport sig_check, sig_on, sig_str, sig_off
from collections import Iterator, Sequence cimport sage.matrix.matrix_dense as matrix_dense from libc.stdio cimport * from sage.structure.element cimport (Matrix, Vector, ModuleElement, Element) from sage.modules.free_module_element cimport FreeModuleElement from sage.libs.gmp.random cimport * from sage.misc.randstate cimport randstate, current_randstate from sage.misc.misc import verbose, get_verbose, cputime from sage.modules.free_module import VectorSpace from sage.modules.vector_mod2_dense cimport Vector_mod2_dense
cdef extern from "gd.h": ctypedef struct gdImagePtr "gdImagePtr": pass
gdImagePtr gdImageCreateFromPng(FILE *f) gdImagePtr gdImageCreateFromPngPtr(int size, void *s) gdImagePtr gdImageCreate(int x, int y) void gdImagePng(gdImagePtr im, FILE *out) void *gdImagePngPtr(gdImagePtr im, int *size) void gdImageDestroy(gdImagePtr im) int gdImageSX(gdImagePtr im) int gdImageSY(gdImagePtr im) int gdImageGetPixel(gdImagePtr im, int x, int y) void gdImageSetPixel(gdImagePtr im, int x, int y, int value) int gdImageColorAllocate(gdImagePtr im, int r, int g, int b) void gdImageFilledRectangle(gdImagePtr im, int x1, int y1, int x2, int y2, int color) void gdFree(void *m)
# Construct global Gray code tables m4ri_build_all_codes()
def free_m4ri(): """ Free global Gray code tables. """ m4ri_destroy_all_codes()
cdef class Matrix_mod2_dense(matrix_dense.Matrix_dense): # dense or sparse """ Dense matrix over GF(2). """ def __cinit__(self, parent, entries, copy, coerce, alloc=True): """ Dense matrix over GF(2) constructor.
INPUT:
- ``parent`` - MatrixSpace. - ``entries`` - may be list or 0 or 1 - ``copy`` - ignored, elements are always copied - ``coerce`` - ignored, elements are always coerced to ints % 2
EXAMPLES::
sage: type(random_matrix(GF(2),2,2)) <type 'sage.matrix.matrix_mod2_dense.Matrix_mod2_dense'>
sage: Matrix(GF(2),3,3,1) # indirect doctest [1 0 0] [0 1 0] [0 0 1]
TESTS:
See :trac:`10858`::
sage: matrix(GF(2),0,[]) * vector(GF(2),0,[]) ()
Large matrices fail gracefully::
sage: import resource sage: if resource.RLIMIT_AS == getattr(resource, 'RLIMIT_RSS', None): ....: # Skip this test if RLIMIT_AS is not properly ....: # supported like on OS X, see Trac #24190 ....: raise RuntimeError("matrix allocation failed") ....: else: # Real test ....: MatrixSpace(GF(2), 2^30)(1) Traceback (most recent call last): ... RuntimeError: matrix allocation failed sage: MatrixSpace(GF(2), 1, 2^40).zero() Traceback (most recent call last): ... OverflowError: ... sage: MatrixSpace(GF(2), 2^40, 1).zero() Traceback (most recent call last): ... OverflowError: ... """
# m4ri assumes that nrows and ncols are of type rci_t: # check for overflow
# cache elements
def __dealloc__(self):
def __init__(self, parent, entries, copy, coerce): """ Dense matrix over GF(2) constructor.
INPUT:
- ``parent`` - MatrixSpace. - ``entries`` - may be list or 0 or 1 - ``copy`` - ignored, elements are always copied - ``coerce`` - ignored, elements are always coerced to ints % 2
EXAMPLES::
sage: type(random_matrix(GF(2),2,2)) <type 'sage.matrix.matrix_mod2_dense.Matrix_mod2_dense'>
sage: Matrix(GF(2),3,3,1) [1 0 0] [0 1 0] [0 0 1]
sage: Matrix(GF(2),2,2,[1,1,1,0]) [1 1] [1 0]
sage: Matrix(GF(2),2,2,4) [0 0] [0 0]
sage: Matrix(GF(2),1,1, 1/3) [1] sage: Matrix(GF(2),1,1, [1/3]) [1]
TESTS::
sage: Matrix(GF(2),0,0) [] sage: Matrix(GF(2),2,0) [] sage: Matrix(GF(2),0,2) [] """ cdef int i,j,e
# scalar ?
# all entries are given as a long iterable
cdef long _hash_(self) except -1: r""" EXAMPLES::
sage: B = random_matrix(GF(2),3,3) sage: B.set_immutable() sage: {B:0} # indirect doctest {[0 1 0] [0 1 1] [0 0 0]: 0} sage: M = random_matrix(GF(2), 123, 321) sage: M.set_immutable() sage: MZ = M.change_ring(ZZ) sage: MZ.set_immutable() sage: hash(M) == hash(MZ) True sage: MS = M.sparse_matrix() sage: MS.set_immutable() sage: hash(M) == hash(MS) True
TESTS::
sage: A = matrix(GF(2),2,0) sage: hash(A) Traceback (most recent call last): ... TypeError: mutable matrices are unhashable sage: A.set_immutable() sage: hash(A) 0
Check that there are no collisions for all matrices up to 4x4, except for the zero matrix and the scalar matrix 1::
sage: L = [] sage: for nr in [1, 2, 3, 4]: # long time ....: for nc in [1, 2, 3, 4]: ....: MS = MatrixSpace(GF(2), nr, nc) ....: for M in MS: ....: if (M == 0) or (M == 1): continue ....: M.set_immutable() ....: L.append(hash(M)) sage: len(L) # long time 74934 sage: len(set(L)) # long time 74934 """ cdef long C[5]
cdef Py_ssize_t i, j
return -2
# this exists for compatibility with matrix_modn_dense cdef set_unsafe_int(self, Py_ssize_t i, Py_ssize_t j, int value): """ Set the (i,j) entry of self to the int value. """
cdef set_unsafe(self, Py_ssize_t i, Py_ssize_t j, value):
cdef get_unsafe(self, Py_ssize_t i, Py_ssize_t j): else:
def str(self, rep_mapping=None, zero=None, plus_one=None, minus_one=None): r""" Return a nice string representation of the matrix.
INPUT:
- ``rep_mapping`` - a dictionary or callable used to override the usual representation of elements. For a dictionary, keys should be elements of the base ring and values the desired string representation.
- ``zero`` - string (default: ``None``); if not ``None`` use the value of ``zero`` as the representation of the zero element.
- ``plus_one`` - string (default: ``None``); if not ``None`` use the value of ``plus_one`` as the representation of the one element.
- ``minus_one`` - Ignored. Only for compatibility with generic matrices.
EXAMPLES::
sage: B = random_matrix(GF(2),3,3) sage: B # indirect doctest [0 1 0] [0 1 1] [0 0 0] sage: block_matrix([[B, 1], [0, B]]) [0 1 0|1 0 0] [0 1 1|0 1 0] [0 0 0|0 0 1] [-----+-----] [0 0 0|0 1 0] [0 0 0|0 1 1] [0 0 0|0 0 0] sage: B.str(zero='.') '[. 1 .]\n[. 1 1]\n[. . .]' """
cdef int i,j, last_i cdef char *row_s cdef char *div_s
# Set the mapping based on keyword arguments # We ignore minus_one (it's only there for compatibility with Matrix) # Shunt mappings off to the generic code since they might not be single characters
cdef list row_div, col_div # Adjacent column divisions messy, use generic code return matrix_dense.Matrix_dense.str(self, rep_mapping)
def row(self, Py_ssize_t i, from_list=False): """ Return the ``i``'th row of this matrix as a vector.
This row is a dense vector if and only if the matrix is a dense matrix.
INPUT:
- ``i`` - integer
- ``from_list`` - bool (default: ``False``); if ``True``, returns the ``i``'th element of ``self.rows()`` (see :func:`rows`), which may be faster, but requires building a list of all rows the first time it is called after an entry of the matrix is changed.
EXAMPLES::
sage: A = random_matrix(GF(2),10,10); A [0 1 0 1 1 0 0 0 1 1] [0 1 1 1 0 1 1 0 0 1] [0 0 0 1 0 1 0 0 1 0] [0 1 1 0 0 1 0 1 1 0] [0 0 0 1 1 1 1 0 1 1] [0 0 1 1 1 1 0 0 0 0] [1 1 1 1 0 1 0 1 1 0] [0 0 0 1 1 0 0 0 1 1] [1 0 0 0 1 1 1 0 1 1] [1 0 0 1 1 0 1 0 0 0]
sage: A.row(0) (0, 1, 0, 1, 1, 0, 0, 0, 1, 1)
sage: A.row(-1) (1, 0, 0, 1, 1, 0, 1, 0, 0, 0)
sage: A.row(2,from_list=True) (0, 0, 0, 1, 0, 1, 0, 0, 1, 0)
sage: A = Matrix(GF(2),1,0) sage: A.row(0) () """ raise IndexError("matrix has no rows") raise IndexError("row index out of range") cdef Py_ssize_t j
######################################################################## # LEVEL 2 functionality # * def _pickle # * def _unpickle # * cdef _mul_ # * cpdef _cmp_ # * _list -- list of underlying elements (need not be a copy) # * _dict -- sparse dictionary of underlying elements (need not be a copy) ######################################################################## # def _pickle(self): # def _unpickle(self, data, int version): # use version >= 0
cpdef _add_(self, right): """ Matrix addition.
INPUT:
- right -- matrix of dimension self.nrows() x self.ncols()
EXAMPLES::
sage: A = random_matrix(GF(2),10,10) sage: A + A == Matrix(GF(2),10,10,0) True
sage: A = random_matrix(GF(2),257,253) sage: A + A == Matrix(GF(2),257,253,0) # indirect doctest True
TESTS::
sage: A = matrix(GF(2),2,0) sage: A+A [] sage: A = matrix(GF(2),0,2) sage: A+A [] sage: A = matrix(GF(2),0,0) sage: A+A [] """ cdef Matrix_mod2_dense A
cpdef _sub_(self, right): """ Matrix addition.
INPUT:
- right -- matrix of dimension self.nrows() x self.ncols()
EXAMPLES::
sage: A = random_matrix(GF(2),10,10) sage: A - A == Matrix(GF(2),10,10,0) # indirect doctest True """
cdef _matrix_times_vector_(self, Vector v): """ EXAMPLES::
sage: A = random_matrix(GF(2),10^4,10^4) sage: v0 = random_matrix(GF(2),10^4,1) sage: v1 = v0.column(0) sage: r0 = A*v0 sage: r1 = A*v1 sage: r0.column(0) == r1 True
TESTS:
Check that :trac:`19378` is fixed::
sage: m = matrix(GF(2), 11, 0) sage: v = vector(GF(2), 0) sage: m * v (0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0) """ cdef mzd_t *tmp M = VectorSpace(self._base_ring, self._nrows) v = M(v) raise ArithmeticError("number of columns of matrix must equal degree of vector")
# If the vector is 0-dimensional, the result will be the 0-vector
cdef _matrix_times_matrix_(self, Matrix right): """ Matrix multiplication.
ALGORITHM: Uses the 'Method of the Four Russians Multiplication', see :func:`_multiply_m4rm`. """ verbose('matrix multiply of %s x %s matrix by %s x %s matrix'%( self._nrows, self._ncols, right._nrows, right._ncols))
cpdef Matrix_mod2_dense _multiply_m4rm(Matrix_mod2_dense self, Matrix_mod2_dense right, int k): """ Multiply matrices using the 'Method of the Four Russians Multiplication' (M4RM) or Konrod's method.
The algorithm is based on an algorithm by Arlazarov, Dinic, Kronrod, and Faradzev [ADKF1970]_ and appeared in [AHU1974]_. This implementation is based on a description given in Gregory Bard's 'Method of the Four Russians Inversion' paper [Bar06]_.
INPUT:
- right -- Matrix - k -- parameter `k` for the Gray Code table size. If `k=0` a suitable value is chosen by the function. (`0<= k <= 16`, default: 0)
EXAMPLES::
sage: A = Matrix(GF(2), 4, 3, [0, 0, 0, 0, 1, 0, 0, 1, 1, 0, 0, 1] ) sage: B = Matrix(GF(2), 3, 4, [0, 0, 1, 0, 1, 0, 0, 1, 1, 1, 0, 0] ) sage: A [0 0 0] [0 1 0] [0 1 1] [0 0 1] sage: B [0 0 1 0] [1 0 0 1] [1 1 0 0] sage: A._multiply_m4rm(B, 0) [0 0 0 0] [1 0 0 1] [0 1 0 1] [1 1 0 0]
TESTS::
sage: A = random_matrix(GF(2),0,0) sage: B = random_matrix(GF(2),0,0) sage: A._multiply_m4rm(B, 0) [] sage: A = random_matrix(GF(2),3,0) sage: B = random_matrix(GF(2),0,3) sage: A._multiply_m4rm(B, 0) [0 0 0] [0 0 0] [0 0 0] sage: A = random_matrix(GF(2),0,3) sage: B = random_matrix(GF(2),3,0) sage: A._multiply_m4rm(B, 0) []
ALGORITHM: Uses the 'Method of the Four Russians' multiplication as implemented in the M4RI library.
REFERENCES:
- [AHU1974]_ - [AKF1970]_ - [Bar2006]_ """ raise ArithmeticError("left ncols must match right nrows")
verbose('m4rm multiply of %s x %s matrix by %s x %s matrix'%( self._nrows, self._ncols, right._nrows, right._ncols))
cdef Matrix_mod2_dense ans
def _multiply_classical(Matrix_mod2_dense self, Matrix_mod2_dense right): """ Classical `O(n^3)` multiplication.
This can be quite fast for matrix vector multiplication but the other routines fall back to this implementation in that case anyway.
EXAMPLES::
sage: A = Matrix(GF(2), 4, 3, [0, 0, 0, 0, 1, 0, 0, 1, 1, 0, 0, 1] ) sage: B = Matrix(GF(2), 3, 4, [0, 0, 1, 0, 1, 0, 0, 1, 1, 1, 0, 0] ) sage: A [0 0 0] [0 1 0] [0 1 1] [0 0 1] sage: B [0 0 1 0] [1 0 0 1] [1 1 0 0] sage: A._multiply_classical(B) [0 0 0 0] [1 0 0 1] [0 1 0 1] [1 1 0 0]
TESTS::
sage: A = random_matrix(GF(2),0,0) sage: B = random_matrix(GF(2),0,0) sage: A._multiply_classical(B) [] sage: A = random_matrix(GF(2),3,0) sage: B = random_matrix(GF(2),0,3) sage: A._multiply_classical(B) [0 0 0] [0 0 0] [0 0 0] sage: A = random_matrix(GF(2),0,3) sage: B = random_matrix(GF(2),3,0) sage: A._multiply_classical(B) [] """ cdef Matrix_mod2_dense A
cpdef Matrix_mod2_dense _multiply_strassen(Matrix_mod2_dense self, Matrix_mod2_dense right, int cutoff): r""" Strassen-Winograd `O(n^{2.807})` multiplication [Str1969]_.
This implementation in M4RI is inspired by Sage's generic Strassen implementation [BHS2008]_ but uses a more memory efficient operation schedule [DP2008]_.
The performance of this routine depends on the parameter cutoff. On many modern machines 2048 should give acceptable performance, a good rule of thumb for calculating the optimal cutoff would that two matrices of the cutoff size should fit in L2 cache, so: `cutoff = \sqrt{L2 * 8 * 1024^2 / 2}` where `L2` is the size of the L2 cache in MB.
INPUT:
- ``right`` - a matrix of matching dimensions. - ``cutoff`` - matrix dimension where M4RM should be used instead of Strassen (default: let M4RI decide)
EXAMPLES::
sage: A = Matrix(GF(2), 4, 3, [0, 0, 0, 0, 1, 0, 0, 1, 1, 0, 0, 1] ) sage: B = Matrix(GF(2), 3, 4, [0, 0, 1, 0, 1, 0, 0, 1, 1, 1, 0, 0] ) sage: A [0 0 0] [0 1 0] [0 1 1] [0 0 1] sage: B [0 0 1 0] [1 0 0 1] [1 1 0 0] sage: A._multiply_strassen(B, 0) [0 0 0 0] [1 0 0 1] [0 1 0 1] [1 1 0 0] sage: A = random_matrix(GF(2),2701,3000) sage: B = random_matrix(GF(2),3000,3172) sage: A._multiply_strassen(B, 256) == A._multiply_m4rm(B, 0) # long time True
TESTS::
sage: A = random_matrix(GF(2),0,0) sage: B = random_matrix(GF(2),0,0) sage: A._multiply_strassen(B, 0) [] sage: A = random_matrix(GF(2),3,0) sage: B = random_matrix(GF(2),0,3) sage: A._multiply_strassen(B, 0) [0 0 0] [0 0 0] [0 0 0] sage: A = random_matrix(GF(2),0,3) sage: B = random_matrix(GF(2),3,0) sage: A._multiply_strassen(B, 0) []
ALGORITHM: Uses Strassen-Winograd matrix multiplication with M4RM as base case as implemented in the M4RI library. """ raise ArithmeticError("left ncols must match right nrows")
cdef Matrix_mod2_dense ans #ans = self.new_matrix(nrows = self._nrows, ncols = right._ncols) # The following is a little faster:
def __neg__(self): """ EXAMPLES::
sage: A = random_matrix(GF(2),100,100) sage: A - A == A - -A True """
def __invert__(self): """ Inverts self using the 'Method of the Four Russians' inversion.
If ``self`` is not invertible a ``ZeroDivisionError`` is raised.
EXAMPLES::
sage: A = Matrix(GF(2),3,3, [0, 0, 1, 0, 1, 1, 1, 0, 1]) sage: MS = A.parent() sage: A [0 0 1] [0 1 1] [1 0 1] sage: ~A [1 0 1] [1 1 0] [1 0 0] sage: A * ~A == ~A * A == MS(1) True
TESTS::
sage: A = matrix(GF(2),0,0) sage: A^(-1) [] """ cdef mzd_t *I cdef Matrix_mod2_dense A
raise ZeroDivisionError("input matrix must be nonsingular") else:
def __copy__(self): """ Returns a copy of ``self``.
EXAMPLES::
sage: MS = MatrixSpace(GF(2),3,3) sage: A = MS(1) sage: A.__copy__() == A True sage: A.__copy__() is A False
sage: A = random_matrix(GF(2),100,100) sage: A.__copy__() == A True sage: A.__copy__() is A False
sage: A.echelonize() sage: A.__copy__() == A True
""" cdef Matrix_mod2_dense A
def _list(self): """ Returns list of the elements of ``self`` in row major ordering.
EXAMPLES::
sage: A = Matrix(GF(2),2,2,[1,0,1,1]) sage: A [1 0] [1 1] sage: A.list() #indirect doctest [1, 0, 1, 1]
TESTS::
sage: A = Matrix(GF(2),3,0) sage: A.list() [] """ cdef int i,j else:
# def _dict(self):
######################################################################## # LEVEL 3 functionality (Optional) # * __deepcopy__ # * Matrix windows -- only if you need strassen for that base ########################################################################
def echelonize(self, algorithm='heuristic', cutoff=0, reduced=True, **kwds): """ Puts self in (reduced) row echelon form.
INPUT:
- self -- a mutable matrix - algorithm
- 'heuristic' -- uses M4RI and PLUQ (default) - 'm4ri' -- uses M4RI - 'pluq' -- uses PLUQ factorization - 'classical' -- uses classical Gaussian elimination
- k -- the parameter 'k' of the M4RI algorithm. It MUST be between 1 and 16 (inclusive). If it is not specified it will be calculated as 3/4 * log_2( min(nrows, ncols) ) as suggested in the M4RI paper. - reduced -- return reduced row echelon form (default:True)
EXAMPLES::
sage: A = random_matrix(GF(2), 10, 10) sage: B = A.__copy__(); B.echelonize() # fastest sage: C = A.__copy__(); C.echelonize(k=2) # force k sage: E = A.__copy__(); E.echelonize(algorithm='classical') # force Gaussian elimination sage: B == C == E True
TESTS::
sage: VF2 = VectorSpace(GF(2),2) sage: WF2 = VF2.submodule([VF2([1,1])]) sage: WF2 Vector space of degree 2 and dimension 1 over Finite Field of size 2 Basis matrix: [1 1]
sage: A2 = matrix(GF(2),2,[1,0,0,1]) sage: A2.kernel() Vector space of degree 2 and dimension 0 over Finite Field of size 2 Basis matrix: []
ALGORITHM:
Uses M4RI library
REFERENCES:
- [Bar2006]_ """ cdef int k, n, full
self.check_mutability() self.clear_cache()
if 'k' in kwds: k = int(kwds['k'])
if k<1 or k>16: raise RuntimeError("k must be between 1 and 16") k = round(k) else: k = 0
sig_on() r = mzd_echelonize_m4ri(self._entries, full, k) sig_off()
self.cache('in_echelon_form',True) self.cache('rank', r) self.cache('pivots', tuple(self._pivots()))
self.check_mutability() self.clear_cache()
sig_on() r = mzd_echelonize_pluq(self._entries, full) sig_off()
self.cache('in_echelon_form',True) self.cache('rank', r) self.cache('pivots', tuple(self._pivots()))
#self._echelonize_linbox() raise NotImplementedError
# for debugging purposes only, it is slow else: raise ValueError("no algorithm '%s'"%algorithm)
def _pivots(self): """ Returns the pivot columns of ``self`` if ``self`` is in row echelon form.
EXAMPLES::
sage: A = matrix(GF(2),5,5,[0,1,0,1,0,0,1,0,1,1,0,1,0,1,0,0,0,0,1,0,0,0,1,0,1]) sage: E = A.echelon_form() sage: E [0 1 0 0 0] [0 0 1 0 0] [0 0 0 1 0] [0 0 0 0 1] [0 0 0 0 0] sage: E._pivots() [1, 2, 3, 4] """ raise RuntimeError("self must be in reduced row echelon form first.") cdef Py_ssize_t i, j, nc
def randomize(self, density=1, nonzero=False): """ Randomize ``density`` proportion of the entries of this matrix, leaving the rest unchanged.
INPUT:
- ``density`` - float; proportion (roughly) to be considered for changes - ``nonzero`` - Bool (default: ``False``); whether the new entries are forced to be non-zero
OUTPUT:
- None, the matrix is modified in-space
EXAMPLES::
sage: A = matrix(GF(2), 5, 5, 0) sage: A.randomize(0.5); A [0 0 0 1 1] [0 1 0 0 1] [1 0 0 0 0] [0 1 0 0 0] [0 0 0 1 0] sage: A.randomize(); A [0 0 1 1 0] [1 1 0 0 1] [1 1 1 1 0] [1 1 1 1 1] [0 0 1 1 0]
TESTS:
With the libc random number generator random(), we had problems where the ranks of all of these matrices would be the same (and they would all be far too low). This verifies that the problem is gone, with Mersenne Twister::
sage: MS2 = MatrixSpace(GF(2), 1000) sage: [MS2.random_element().rank() for i in range(5)] [999, 998, 1000, 999, 999]
Testing corner case::
sage: A = random_matrix(GF(2),3,0) sage: A [] """
return density = float(1)
cdef int i, j, k cdef int nc cdef unsigned int low, high
# Original code, before adding the ``nonzero`` option. # for portability we get 32-bit twice rather than 64-bit once else:
# New code for the case when ``nonzero`` is ``True``. else: sig_on() for i from 0 <= i < self._nrows: for j from 0 <= j < self._ncols: if rstate.c_rand_double() <= density: mzd_write_bit(self._entries, i, j, 1) sig_off()
cdef rescale_row_c(self, Py_ssize_t row, multiple, Py_ssize_t start_col): """ EXAMPLES::
sage: A = random_matrix(GF(2),3,3); A [0 1 0] [0 1 1] [0 0 0]
sage: A.rescale_row(0,0,0); A [0 0 0] [0 1 1] [0 0 0] """
cdef add_multiple_of_row_c(self, Py_ssize_t row_to, Py_ssize_t row_from, multiple, Py_ssize_t start_col): """ EXAMPLES::
sage: A = random_matrix(GF(2),3,3); A [0 1 0] [0 1 1] [0 0 0] sage: A.add_multiple_of_row(0,1,1,0); A [0 0 1] [0 1 1] [0 0 0] """
cdef swap_rows_c(self, Py_ssize_t row1, Py_ssize_t row2): """ EXAMPLES::
sage: A = random_matrix(GF(2),3,3); A [0 1 0] [0 1 1] [0 0 0] sage: A.swap_rows(0,1); A [0 1 1] [0 1 0] [0 0 0] """
cdef swap_columns_c(self, Py_ssize_t col1, Py_ssize_t col2): """ EXAMPLES::
sage: A = random_matrix(GF(2),3,3); A [0 1 0] [0 1 1] [0 0 0]
sage: A.swap_columns(0,1); A [1 0 0] [1 0 1] [0 0 0]
sage: A = random_matrix(GF(2),3,65)
sage: B = A.__copy__() sage: B.swap_columns(0,1) sage: B.swap_columns(0,1) sage: A == B True
sage: A.swap_columns(0,64) sage: A.column(0) == B.column(64) True sage: A.swap_columns(63,64) sage: A.column(63) == B.column(0) True """
def _magma_init_(self, magma): """ Returns a string of self in ``Magma`` form. Does not return ``Magma`` object but string.
EXAMPLES::
sage: A = random_matrix(GF(2),3,3) sage: A._magma_init_(magma) # optional - magma 'Matrix(GF(2),3,3,StringToIntegerSequence("0 1 0 0 1 1 0 0 0"))' sage: A = random_matrix(GF(2),100,100) sage: B = random_matrix(GF(2),100,100) sage: magma(A*B) == magma(A) * magma(B) # optional - magma True
TESTS::
sage: A = random_matrix(GF(2),0,3) sage: magma(A) # optional - magma Matrix with 0 rows and 3 columns sage: A = matrix(GF(2),2,3,[0,1,1,1,0,0]) sage: A._magma_init_(magma) # optional - magma 'Matrix(GF(2),2,3,StringToIntegerSequence("0 1 1 1 0 0"))' sage: magma(A) # optional - magma [0 1 1] [1 0 0] """ s = self.base_ring()._magma_init_(magma) return 'Matrix(%s,%s,%s,StringToIntegerSequence("%s"))'%( s, self._nrows, self._ncols, self._export_as_string())
def determinant(self): """ Return the determinant of this matrix over GF(2).
EXAMPLES::
sage: matrix(GF(2),2,[1,1,0,1]).determinant() 1 sage: matrix(GF(2),2,[1,1,1,1]).determinant() 0 """
def transpose(self): """ Returns transpose of self and leaves self untouched.
EXAMPLES::
sage: A = Matrix(GF(2),3,5,[1, 0, 1, 0, 0, 0, 1, 1, 0, 0, 1, 1, 0, 1, 0]) sage: A [1 0 1 0 0] [0 1 1 0 0] [1 1 0 1 0] sage: B = A.transpose(); B [1 0 1] [0 1 1] [1 1 0] [0 0 1] [0 0 0] sage: B.transpose() == A True
``.T`` is a convenient shortcut for the transpose::
sage: A.T [1 0 1] [0 1 1] [1 1 0] [0 0 1] [0 0 0]
TESTS::
sage: A = random_matrix(GF(2),0,40) sage: A.transpose() 40 x 0 dense matrix over Finite Field of size 2 (use the '.str()' method to see the entries)
sage: A = Matrix(GF(2), [1,0]) sage: B = A.transpose() sage: A[0,0] = 0 sage: B[0,0] 1 """
cpdef int _cmp_(self, right) except -2: """ Compares ``self`` with ``right``. While equality and inequality are clearly defined, ``<`` and ``>`` are not. For those first the matrix dimensions of ``self`` and ``right`` are compared. If these match then ``<`` means that there is a position smallest (i,j) in ``self`` where ``self[i,j]`` is zero but ``right[i,j]`` is one. This (i,j) is smaller than the (i,j) if ``self`` and ``right`` are exchanged for the comparison.
EXAMPLES::
sage: A = MatrixSpace(GF(2),3,3).one() sage: B = copy(MatrixSpace(GF(2),3,3).one()) sage: B[0,1] = 1 sage: A < B True
TESTS::
sage: A = matrix(GF(2),2,0) sage: B = matrix(GF(2),2,0) sage: A < B False """ if self._nrows == 0 or self._ncols == 0: return 0 return mzd_cmp(self._entries, (<Matrix_mod2_dense>right)._entries)
def augment(self, right, subdivide=False): r""" Augments ``self`` with ``right``.
EXAMPLES::
sage: MS = MatrixSpace(GF(2),3,3) sage: A = MS([0, 1, 0, 1, 1, 0, 1, 1, 1]); A [0 1 0] [1 1 0] [1 1 1] sage: B = A.augment(MS(1)); B [0 1 0 1 0 0] [1 1 0 0 1 0] [1 1 1 0 0 1] sage: B.echelonize(); B [1 0 0 1 1 0] [0 1 0 1 0 0] [0 0 1 0 1 1] sage: C = B.matrix_from_columns([3,4,5]); C [1 1 0] [1 0 0] [0 1 1] sage: C == ~A True sage: C*A == MS(1) True
A vector may be augmented to a matrix. ::
sage: A = matrix(GF(2), 3, 4, range(12)) sage: v = vector(GF(2), 3, range(3)) sage: A.augment(v) [0 1 0 1 0] [0 1 0 1 1] [0 1 0 1 0]
The ``subdivide`` option will add a natural subdivision between ``self`` and ``right``. For more details about how subdivisions are managed when augmenting, see :meth:`sage.matrix.matrix1.Matrix.augment`. ::
sage: A = matrix(GF(2), 3, 5, range(15)) sage: B = matrix(GF(2), 3, 3, range(9)) sage: A.augment(B, subdivide=True) [0 1 0 1 0|0 1 0] [1 0 1 0 1|1 0 1] [0 1 0 1 0|0 1 0]
TESTS::
sage: A = random_matrix(GF(2),2,3) sage: B = random_matrix(GF(2),2,0) sage: A.augment(B) [0 1 0] [0 1 1]
sage: B.augment(A) [0 1 0] [0 1 1]
sage: M = Matrix(GF(2), 0, 0, 0) sage: N = Matrix(GF(2), 0, 19, 0) sage: W = M.augment(N) sage: W.ncols() 19 sage: M = Matrix(GF(2), 0, 1, 0) sage: N = Matrix(GF(2), 0, 1, 0) sage: M.augment(N) []
Check that :trac:`19165` is solved::
sage: m = matrix(GF(2), 2, range(4)) sage: m.augment(matrix(GF(2), 2, range(4), sparse=True)) [0 1 0 1] [0 1 0 1]
sage: m.augment(1) Traceback (most recent call last): ... TypeError: right must either be a matrix or a vector. Not <type 'sage.rings.integer.Integer'> """ cdef Matrix_mod2_dense other
else:
raise TypeError("Both numbers of rows must match.")
cdef Matrix_mod2_dense Z
cdef _stack_impl(self, bottom): r""" Stack ``self`` on top of ``bottom``.
EXAMPLES::
sage: A = matrix(GF(2),2,2,[1,0,0,1]) sage: B = matrix(GF(2),2,2,[0,1,1,0]) sage: A.stack(B) [1 0] [0 1] [0 1] [1 0] sage: B.stack(A) [0 1] [1 0] [1 0] [0 1]
A vector may be stacked below a matrix. ::
sage: A = matrix(GF(2), 2, 5, range(10)) sage: v = vector(GF(2), 5, range(5)) sage: A.stack(v) [0 1 0 1 0] [1 0 1 0 1] [0 1 0 1 0]
The ``subdivide`` option will add a natural subdivision between ``self`` and ``bottom``. For more details about how subdivisions are managed when stacking, see :meth:`sage.matrix.matrix1.Matrix.stack`. ::
sage: A = matrix(GF(2), 3, 5, range(15)) sage: B = matrix(GF(2), 1, 5, range(5)) sage: A.stack(B, subdivide=True) [0 1 0 1 0] [1 0 1 0 1] [0 1 0 1 0] [---------] [0 1 0 1 0]
TESTS::
sage: A = random_matrix(GF(2),0,3) sage: B = random_matrix(GF(2),3,3) sage: A.stack(B) [0 1 0] [0 1 1] [0 0 0]
sage: B.stack(A) [0 1 0] [0 1 1] [0 0 0]
sage: M = Matrix(GF(2), 0, 0, 0) sage: N = Matrix(GF(2), 19, 0, 0) sage: W = M.stack(N) sage: W.nrows() 19 sage: M = Matrix(GF(2), 1, 0, 0) sage: N = Matrix(GF(2), 1, 0, 0) sage: M.stack(N) [] """ cdef Matrix_mod2_dense Z
def submatrix(self, Py_ssize_t row=0, Py_ssize_t col=0, Py_ssize_t nrows=-1, Py_ssize_t ncols=-1): """ Return submatrix from the index row, col (inclusive) with dimension nrows x ncols.
INPUT:
- row -- index of start row - col -- index of start column - nrows -- number of rows of submatrix - ncols -- number of columns of submatrix
EXAMPLES::
sage: A = random_matrix(GF(2),200,200) sage: A[0:2,0:2] == A.submatrix(0,0,2,2) True sage: A[0:100,0:100] == A.submatrix(0,0,100,100) True sage: A == A.submatrix(0,0,200,200) True
sage: A[1:3,1:3] == A.submatrix(1,1,2,2) True sage: A[1:100,1:100] == A.submatrix(1,1,99,99) True sage: A[1:200,1:200] == A.submatrix(1,1,199,199) True
TESTS for handling of default arguments (:trac:`18761`)::
sage: A.submatrix(17,15) == A.submatrix(17,15,183,185) True sage: A.submatrix(row=100,col=37,nrows=1,ncols=3) == A.submatrix(100,37,1,3) True """ cdef Matrix_mod2_dense A
cdef int highr, highc
nrows = 0
ncols = 0
raise TypeError("Expected row >= 0, but got %d instead."%row)
raise TypeError("Expected col >= 0, but got %d instead."%col)
raise TypeError("Expected highc <= self.ncols(), but got %d > %d instead."%(highc, self._entries.ncols))
raise TypeError("Expected highr <= self.nrows(), but got %d > %d instead."%(highr, self._entries.nrows))
return A
def __reduce__(self): """ Serialize ``self``.
EXAMPLES::
sage: A = random_matrix(GF(2),10,10) sage: f,s = A.__reduce__() sage: f(*s) == A True """ cdef int i,j, r,c, size
cpdef _export_as_string(self): """ Return space separated string of the entries in this matrix.
EXAMPLES::
sage: w = matrix(GF(2),2,3,[1,0,1,1,1,0]) sage: w._export_as_string() '1 0 1 1 1 0' """ cdef Py_ssize_t i, j, k, n cdef char *s cdef char *t
data = '' else:
def density(self, approx=False): """ Return the density of this matrix.
By density we understand the ration of the number of nonzero positions and the self.nrows() * self.ncols(), i.e. the number of possible nonzero positions.
INPUT:
- approx -- return floating point approximation (default: False)
EXAMPLES::
sage: A = random_matrix(GF(2),1000,1000) sage: d = A.density(); d 62483/125000
sage: float(d) 0.499864
sage: A.density(approx=True) 0.499864000...
sage: float(len(A.nonzero_positions())/1000^2) 0.499864 """ else:
def rank(self, algorithm='ple'): """ Return the rank of this matrix.
On average 'ple' should be faster than 'm4ri' and hence it is the default choice. However, for small - i.e. quite few thousand rows & columns - and sparse matrices 'm4ri' might be a better choice.
INPUT:
- ``algorithm`` - either "ple" or "m4ri"
EXAMPLES::
sage: A = random_matrix(GF(2), 1000, 1000) sage: A.rank() 999
sage: A = matrix(GF(2),10, 0) sage: A.rank() 0 """ cdef mzp_t *P cdef mzp_t *Q
elif algorithm == 'm4ri': r = mzd_echelonize_m4ri(A, 0, 0) else: raise ValueError("Algorithm '%s' unknown."%algorithm)
def _right_kernel_matrix(self, **kwds): r""" Returns a pair that includes a matrix of basis vectors for the right kernel of ``self``.
INPUT:
- ``kwds`` - these are provided for consistency with other versions of this method. Here they are ignored as there is no optional behavior available.
OUTPUT:
Returns a pair. First item is the string 'computed-pluq' that identifies the nature of the basis vectors.
Second item is a matrix whose rows are a basis for the right kernel, over the integers mod 2, as computed by the M4RI library using PLUQ matrix decomposition.
EXAMPLES::
sage: A = matrix(GF(2), [ ....: [0, 1, 0, 0, 1, 0, 1, 1], ....: [1, 0, 1, 0, 0, 1, 1, 0], ....: [0, 0, 1, 0, 0, 1, 0, 1], ....: [1, 0, 1, 1, 0, 1, 1, 0], ....: [0, 0, 1, 0, 0, 1, 0, 1], ....: [1, 1, 0, 1, 1, 0, 0, 0]]) sage: A [0 1 0 0 1 0 1 1] [1 0 1 0 0 1 1 0] [0 0 1 0 0 1 0 1] [1 0 1 1 0 1 1 0] [0 0 1 0 0 1 0 1] [1 1 0 1 1 0 0 0] sage: result = A._right_kernel_matrix() sage: result[0] 'computed-pluq' sage: result[1] [0 1 0 0 1 0 0 0] [0 0 1 0 0 1 0 0] [1 1 0 0 0 0 1 0] [1 1 1 0 0 0 0 1] sage: X = result[1].transpose() sage: A*X == zero_matrix(GF(2), 6, 4) True
TESTS:
We test the three trivial cases. Matrices with no rows or columns will cause segfaults in the M4RI code, so we protect against that instance. Computing a kernel or a right kernel matrix should never pass these problem matrices here. ::
sage: A = matrix(GF(2), 0, 2) sage: A._right_kernel_matrix() Traceback (most recent call last): ... ValueError: kernels of matrices mod 2 with zero rows or zero columns cannot be computed sage: A = matrix(GF(2), 2, 0) sage: A._right_kernel_matrix() Traceback (most recent call last): ... ValueError: kernels of matrices mod 2 with zero rows or zero columns cannot be computed sage: A = zero_matrix(GF(2), 4, 3) sage: A._right_kernel_matrix()[1] [1 0 0] [0 1 0] [0 0 1] """ cdef Matrix_mod2_dense M # Despite the name, this next call returns X such that M*X = 0 else:
# Used for hashing cdef int i, k cdef unsigned long parity_table[256] for i from 0 <= i < 256: parity_table[i] = 1 & ((i) ^ (i>>1) ^ (i>>2) ^ (i>>3) ^ (i>>4) ^ (i>>5) ^ (i>>6) ^ (i>>7))
# gmp's ULONG_PARITY may use special # assembly instructions, could be faster cpdef inline unsigned long parity(m4ri_word a): """ Returns the parity of the number of bits in a.
EXAMPLES::
sage: from sage.matrix.matrix_mod2_dense import parity sage: parity(1) 1L sage: parity(3) 0L sage: parity(0x10000101011) 1L """
cdef inline unsigned long parity_mask(m4ri_word a): return -parity(a)
def unpickle_matrix_mod2_dense_v1(r, c, data, size): """ Deserialize a matrix encoded in the string ``s``.
INPUT:
- r -- number of rows of matrix - c -- number of columns of matrix - s -- a string - size -- length of the string s
EXAMPLES::
sage: A = random_matrix(GF(2),100,101) sage: _,(r,c,s,s2) = A.__reduce__() sage: from sage.matrix.matrix_mod2_dense import unpickle_matrix_mod2_dense_v1 sage: unpickle_matrix_mod2_dense_v1(r,c,s,s2) == A True sage: loads(dumps(A)) == A True """
cdef int i, j cdef Matrix_mod2_dense A
raise TypeError("Pickled data dimension doesn't match.")
def from_png(filename): """ Returns a dense matrix over GF(2) from a 1-bit PNG image read from ``filename``. No attempt is made to verify that the filename string actually points to a PNG image.
INPUT:
- filename -- a string
EXAMPLES::
sage: from sage.matrix.matrix_mod2_dense import from_png, to_png sage: A = random_matrix(GF(2),10,10) sage: fn = tmp_filename() sage: to_png(A, fn) sage: B = from_png(fn) sage: A == B True """
cdef int i,j,r,c cdef Matrix_mod2_dense A
def to_png(Matrix_mod2_dense A, filename): """ Saves the matrix ``A`` to filename as a 1-bit PNG image.
INPUT:
- ``A`` - a matrix over GF(2) - ``filename`` - a string for a file in a writable position
EXAMPLES::
sage: from sage.matrix.matrix_mod2_dense import from_png, to_png sage: A = random_matrix(GF(2),10,10) sage: fn = tmp_filename() sage: to_png(A, fn) sage: B = from_png(fn) sage: A == B True """ cdef int i,j, r,c raise TypeError("Cannot write image with dimensions %d x %d"%(c,r))
def pluq(Matrix_mod2_dense A, algorithm="standard", int param=0): """ Return PLUQ factorization of A.
INPUT:
- A -- matrix - algorithm
* 'standard' asymptotically fast (default) * 'mmpf' M4RI inspired * 'naive' naive cubic
- param -- either k for 'mmpf' is chosen or matrix multiplication cutoff for 'standard' (default: 0)
EXAMPLES::
sage: from sage.matrix.matrix_mod2_dense import pluq sage: A = random_matrix(GF(2),4,4); A [0 1 0 1] [0 1 1 1] [0 0 0 1] [0 1 1 0]
sage: LU, P, Q = pluq(A) sage: LU [1 0 1 0] [1 1 0 0] [0 0 1 0] [1 1 1 0]
sage: P [0, 1, 2, 3]
sage: Q [1, 2, 3, 3] """
elif algorithm == "mmpf": sig_on() _mzd_pluq_russian(B._entries, p, q, param) sig_off() elif algorithm == "naive": sig_on() _mzd_pluq_naive(B._entries, p, q) sig_off() else: raise ValueError("Algorithm '%s' unknown."%algorithm)
def ple(Matrix_mod2_dense A, algorithm="standard", int param=0): """ Return PLE factorization of A.
INPUT:
- A -- matrix - algorithm
- 'standard' asymptotically fast (default) - 'russian' M4RI inspired - 'naive' naive cubic
- param -- either k for 'mmpf' is chosen or matrix multiplication cutoff for 'standard' (default: 0)
EXAMPLES::
sage: from sage.matrix.matrix_mod2_dense import ple sage: A = random_matrix(GF(2),4,4); A [0 1 0 1] [0 1 1 1] [0 0 0 1] [0 1 1 0]
sage: LU, P, Q = ple(A) sage: LU [1 0 0 1] [1 1 0 0] [0 0 1 0] [1 1 1 0]
sage: P [0, 1, 2, 3]
sage: Q [1, 2, 3, 3]
sage: A = random_matrix(GF(2),1000,1000) sage: ple(A) == ple(A,'russian') == ple(A,'naive') True """
else: raise ValueError("Algorithm '%s' unknown."%algorithm)
|