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^e}` for `2 <= e <= 10` using the M4RIE library.
The M4RIE library offers two matrix representations:
1) ``mzed_t``
m x n matrices over `\GF{2^e}` are internally represented roughly as m x (en) matrices over `\GF{2}`. Several elements are packed into words such that each element is filled with zeroes until the next power of two. Thus, for example, elements of `\GF{2^3}` are represented as ``[0xxx|0xxx|0xxx|0xxx|...]``. This representation is wrapped as :class:`Matrix_gf2e_dense` in Sage.
Multiplication and elimination both use "Newton-John" tables. These tables are simply all possible multiples of a given row in a matrix such that a scale+add operation is reduced to a table lookup + add. On top of Newton-John multiplication M4RIE implements asymptotically fast Strassen-Winograd multiplication. Elimination uses simple Gaussian elimination which requires `O(n^3)` additions but only `O(n^2 * 2^e)` multiplications.
2) ``mzd_slice_t``
m x n matrices over `\GF{2^e}` are internally represented as slices of m x n matrices over `\GF{2}`. This representation allows for very fast matrix times matrix products using Karatsuba's polynomial multiplication for polynomials over matrices. However, it is not feature complete yet and hence not wrapped in Sage for now.
See http://m4ri.sagemath.org for more details on the M4RIE library.
EXAMPLES::
sage: K.<a> = GF(2^8) sage: A = random_matrix(K, 3,4) sage: A [ a^6 + a^5 + a^4 + a^2 a^6 + a^3 + a + 1 a^5 + a^3 + a^2 + a + 1 a^7 + a^6 + a + 1] [ a^7 + a^6 + a^3 a^7 + a^6 + a^5 + 1 a^5 + a^4 + a^3 + a + 1 a^6 + a^5 + a^4 + a^3 + a^2 + 1] [ a^6 + a^5 + a + 1 a^7 + a^3 + 1 a^7 + a^3 + a + 1 a^7 + a^6 + a^3 + a^2 + a + 1]
sage: A.echelon_form() [ 1 0 0 a^6 + a^5 + a^4 + a^2] [ 0 1 0 a^7 + a^5 + a^3 + a + 1] [ 0 0 1 a^6 + a^4 + a^3 + a^2 + 1]
AUTHOR:
* Martin Albrecht <martinralbrecht@googlemail.com>
TESTS::
sage: TestSuite(sage.matrix.matrix_gf2e_dense.Matrix_gf2e_dense).run(verbose=True) running ._test_new() . . . pass running ._test_pickling() . . . pass
Test hashing::
sage: K.<a> = GF(2^4) sage: A = random_matrix(K, 1000, 1000) sage: A.set_immutable() sage: {A:1} {1000 x 1000 dense matrix over Finite Field in a of size 2^4: 1}
.. TODO::
Wrap ``mzd_slice_t``.
REFERENCES:
- [BB2009]_
"""
#***************************************************************************** # Copyright (C) 2010 Martin Albrecht <martinralbrecht@googlemail.com> # # 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.signals cimport sig_check, sig_on, sig_off
cimport sage.matrix.matrix_dense as matrix_dense from sage.structure.element cimport Matrix, Vector from sage.structure.element cimport ModuleElement, Element, RingElement
from sage.rings.all import FiniteField as GF from sage.misc.randstate cimport randstate, current_randstate
from sage.matrix.matrix_mod2_dense cimport Matrix_mod2_dense
from sage.libs.m4ri cimport m4ri_word, mzd_copy, mzd_init from sage.libs.m4rie cimport * from sage.libs.m4rie cimport mzed_t
# we must keep a copy of the internal finite field representation # around to avoid re-creating it over and over again. Furthermore, # M4RIE assumes pointer equivalence of identical fields.
_m4rie_finite_field_cache = {}
cdef class M4RIE_finite_field: """ A thin wrapper around the M4RIE finite field class such that we can put it in a hash table. This class is not meant for public consumption. """ cdef gf2e *ff
def __cinit__(self): """ EXAMPLES::
sage: from sage.matrix.matrix_gf2e_dense import M4RIE_finite_field sage: K = M4RIE_finite_field(); K <sage.matrix.matrix_gf2e_dense.M4RIE_finite_field object at 0x...> """ pass
def __dealloc__(self): """ EXAMPLES::
sage: from sage.matrix.matrix_gf2e_dense import M4RIE_finite_field sage: K = M4RIE_finite_field() sage: del K """ gf2e_free(self.ff)
cdef m4ri_word poly_to_word(f):
cdef object word_to_poly(w, F):
cdef class Matrix_gf2e_dense(matrix_dense.Matrix_dense): def __cinit__(self, parent, entries, copy, coerce, alloc=True): """ Create new matrix over `GF(2^e)` for 2<=e<=10.
INPUT:
- ``parent`` - a :class:`MatrixSpace`. - ``entries`` - may be list or a finite field element. - ``copy`` - ignored, elements are always copied - ``coerce`` - ignored, elements are always coerced - ``alloc`` - if ``True`` the matrix is allocated first (default: ``True``)
EXAMPLES::
sage: K.<a> = GF(2^4) sage: A = Matrix(K, 3, 4); A [0 0 0 0] [0 0 0 0] [0 0 0 0]
sage: A.randomize(); A [ a^2 a^3 + a + 1 a^3 + a^2 + a + 1 a + 1] [ a^3 1 a^3 + a + 1 a^3 + a^2 + 1] [ a + 1 a^3 + 1 a^3 + a + 1 a^3 + a^2 + a + 1]
sage: K.<a> = GF(2^3) sage: A = Matrix(K,3,4); A [0 0 0 0] [0 0 0 0] [0 0 0 0]
sage: A.randomize(); A [ a^2 + a a^2 + 1 a^2 + a a^2 + a] [ a^2 + 1 a^2 + a + 1 a^2 + 1 a^2] [ a^2 + a a^2 + 1 a^2 + a + 1 a + 1] """
cdef M4RIE_finite_field FF
else:
# cache elements
def __dealloc__(self): """ TESTS::
sage: K.<a> = GF(2^4) sage: A = Matrix(K, 1000, 1000) sage: del A sage: A = Matrix(K, 1000, 1000) sage: del A """
def __init__(self, parent, entries, copy, coerce): """ Create new matrix over `GF(2^e)` for 2<=e<=10.
INPUT:
- ``parent`` - a :class:`MatrixSpace`. - ``entries`` - may be list or a finite field element. - ``copy`` - ignored, elements are always copied - ``coerce`` - ignored, elements are always coerced
EXAMPLES::
sage: K.<a> = GF(2^4) sage: l = [K.random_element() for _ in range(3*4)]; l [a^2 + 1, a^3 + 1, 0, 0, a, a^3 + a + 1, a + 1, a + 1, a^2, a^3 + a + 1, a^3 + a, a^3 + a]
sage: A = Matrix(K, 3, 4, l); A [ a^2 + 1 a^3 + 1 0 0] [ a a^3 + a + 1 a + 1 a + 1] [ a^2 a^3 + a + 1 a^3 + a a^3 + a]
sage: A.list() [a^2 + 1, a^3 + 1, 0, 0, a, a^3 + a + 1, a + 1, a + 1, a^2, a^3 + a + 1, a^3 + a, a^3 + a]
sage: l[0], A[0,0] (a^2 + 1, a^2 + 1)
sage: A = Matrix(K, 3, 3, a); A [a 0 0] [0 a 0] [0 0 a] """ cdef int i,j
# scalar ?
# all entries are given as a long list raise IndexError("The vector of entries has the wrong length.")
cdef set_unsafe(self, Py_ssize_t i, Py_ssize_t j, value): """ A[i,j] = value without bound checks
INPUT: - ``i`` - row index - ``j`` - column index - ``value`` - a finite field element (not checked but assumed)
EXAMPLES::
sage: K.<a> = GF(2^4) sage: A = Matrix(K,3,4,[K.random_element() for _ in range(3*4)]); A [ a^2 + 1 a^3 + 1 0 0] [ a a^3 + a + 1 a + 1 a + 1] [ a^2 a^3 + a + 1 a^3 + a a^3 + a]
sage: A[0,0] = a # indirect doctest sage: A [ a a^3 + 1 0 0] [ a a^3 + a + 1 a + 1 a + 1] [ a^2 a^3 + a + 1 a^3 + a a^3 + a] """
cdef get_unsafe(self, Py_ssize_t i, Py_ssize_t j): """ Get A[i,j] without bound checks.
INPUT: - ``i`` - row index - ``j`` - column index
EXAMPLES::
sage: K.<a> = GF(2^4) sage: A = random_matrix(K,3,4) sage: A[2,3] # indirect doctest a^3 + a^2 + a + 1 sage: K.<a> = GF(2^3) sage: m,n = 3, 4 sage: A = random_matrix(K,3,4); A [ a^2 + a a^2 + 1 a^2 + a a^2 + a] [ a^2 + 1 a^2 + a + 1 a^2 + 1 a^2] [ a^2 + a a^2 + 1 a^2 + a + 1 a + 1] """
cpdef _add_(self, right): """ Return A+B
INPUT:
- ``right`` - a matrix
EXAMPLES::
sage: K.<a> = GF(2^4) sage: A = random_matrix(K,3,4); A [ a^2 a^3 + a + 1 a^3 + a^2 + a + 1 a + 1] [ a^3 1 a^3 + a + 1 a^3 + a^2 + 1] [ a + 1 a^3 + 1 a^3 + a + 1 a^3 + a^2 + a + 1]
sage: B = random_matrix(K,3,4); B [ a^2 + a a^2 + 1 a^2 + a a^3 + a^2 + a] [ a^2 + 1 a^3 + a^2 + a + 1 a^2 + 1 a^2] [ a^3 + a^2 + a a^2 + 1 a^2 + a + 1 a^3 + a + 1]
sage: C = A + B; C # indirect doctest [ a a^3 + a^2 + a a^3 + 1 a^3 + a^2 + 1] [a^3 + a^2 + 1 a^3 + a^2 + a a^3 + a^2 + a a^3 + 1] [a^3 + a^2 + 1 a^3 + a^2 a^3 + a^2 a^2] """ cdef Matrix_gf2e_dense A return A
cpdef _sub_(self, right): """ EXAMPLES::
sage: from sage.matrix.matrix_gf2e_dense import Matrix_gf2e_dense sage: K.<a> = GF(2^4) sage: m,n = 3, 4 sage: MS = MatrixSpace(K,m,n) sage: A = random_matrix(K,3,4); A [ a^2 a^3 + a + 1 a^3 + a^2 + a + 1 a + 1] [ a^3 1 a^3 + a + 1 a^3 + a^2 + 1] [ a + 1 a^3 + 1 a^3 + a + 1 a^3 + a^2 + a + 1]
sage: B = random_matrix(K,3,4); B [ a^2 + a a^2 + 1 a^2 + a a^3 + a^2 + a] [ a^2 + 1 a^3 + a^2 + a + 1 a^2 + 1 a^2] [ a^3 + a^2 + a a^2 + 1 a^2 + a + 1 a^3 + a + 1]
sage: C = A - B; C # indirect doctest [ a a^3 + a^2 + a a^3 + 1 a^3 + a^2 + 1] [a^3 + a^2 + 1 a^3 + a^2 + a a^3 + a^2 + a a^3 + 1] [a^3 + a^2 + 1 a^3 + a^2 a^3 + a^2 a^2] """
def _multiply_classical(self, Matrix right): """ Classical cubic matrix multiplication.
EXAMPLES::
sage: K.<a> = GF(2^2) sage: A = random_matrix(K, 50, 50) sage: B = random_matrix(K, 50, 50) sage: A*B == A._multiply_classical(B) True
sage: K.<a> = GF(2^4) sage: A = random_matrix(K, 50, 50) sage: B = random_matrix(K, 50, 50) sage: A*B == A._multiply_classical(B) True
sage: K.<a> = GF(2^8) sage: A = random_matrix(K, 50, 50) sage: B = random_matrix(K, 50, 50) sage: A*B == A._multiply_classical(B) True
sage: K.<a> = GF(2^10) sage: A = random_matrix(K, 50, 50) sage: B = random_matrix(K, 50, 50) sage: A*B == A._multiply_classical(B) True
.. NOTE::
This function is very slow. Use ``*`` operator instead. """ raise ArithmeticError("left ncols must match right nrows")
cdef Matrix_gf2e_dense ans
return ans
cdef _matrix_times_matrix_(self, Matrix right): """ Return A*B
Uses the M4RIE machinery to decide which function to call.
INPUT:
- ``right`` - a matrix
EXAMPLES::
sage: K.<a> = GF(2^3) sage: A = random_matrix(K, 51, 50) sage: B = random_matrix(K, 50, 52) sage: A*B == A._multiply_newton_john(B) # indirect doctest True
sage: K.<a> = GF(2^5) sage: A = random_matrix(K, 10, 50) sage: B = random_matrix(K, 50, 12) sage: A*B == A._multiply_newton_john(B) True
sage: K.<a> = GF(2^7) sage: A = random_matrix(K,100, 50) sage: B = random_matrix(K, 50, 17) sage: A*B == A._multiply_classical(B) True """ raise ArithmeticError("left ncols must match right nrows")
cdef Matrix_gf2e_dense ans
return ans
cpdef Matrix_gf2e_dense _multiply_newton_john(Matrix_gf2e_dense self, Matrix_gf2e_dense right): """ Return A*B using Newton-John tables.
We can write classical cubic multiplication (``C=A*B``) as::
for i in range(A.ncols()): for j in range(A.nrows()): C[j] += A[j,i] * B[j]
Hence, in the inner-most loop we compute multiples of ``B[j]`` by the values ``A[j,i]``. If the matrix ``A`` is big and the finite field is small, there is a very high chance that ``e * B[j]`` is computed more than once for any ``e`` in the finite field. Instead, we compute all possible multiples of ``B[j]`` and re-use this data in the inner loop. This is what is called a "Newton-John" table in M4RIE.
INPUT:
- ``right`` - a matrix
EXAMPLES::
sage: K.<a> = GF(2^2) sage: A = random_matrix(K, 50, 50) sage: B = random_matrix(K, 50, 50) sage: A._multiply_newton_john(B) == A._multiply_classical(B) # indirect doctest True
sage: K.<a> = GF(2^4) sage: A = random_matrix(K, 50, 50) sage: B = random_matrix(K, 50, 50) sage: A._multiply_newton_john(B) == A._multiply_classical(B) True
sage: K.<a> = GF(2^8) sage: A = random_matrix(K, 50, 50) sage: B = random_matrix(K, 50, 50) sage: A._multiply_newton_john(B) == A._multiply_classical(B) True
sage: K.<a> = GF(2^10) sage: A = random_matrix(K, 50, 50) sage: B = random_matrix(K, 50, 50) sage: A._multiply_newton_john(B) == A._multiply_classical(B) True """ raise ArithmeticError("left ncols must match right nrows")
cdef Matrix_gf2e_dense ans
return ans
cpdef Matrix_gf2e_dense _multiply_karatsuba(Matrix_gf2e_dense self, Matrix_gf2e_dense right): """ Matrix multiplication using Karatsuba over polynomials with matrix coefficients over GF(2).
The idea behind Karatsuba multiplication for matrices over `\GF{p^n}` is to treat these matrices as polynomials with coefficients of matrices over `\GF{p}`. Then, Karatsuba-style formulas can be used to perform multiplication, cf. [BB2009]_.
INPUT:
- ``right`` - a matrix
EXAMPLES::
sage: K.<a> = GF(2^2) sage: A = random_matrix(K, 50, 50) sage: B = random_matrix(K, 50, 50) sage: A._multiply_karatsuba(B) == A._multiply_classical(B) # indirect doctest True
sage: K.<a> = GF(2^2) sage: A = random_matrix(K, 137, 11) sage: B = random_matrix(K, 11, 23) sage: A._multiply_karatsuba(B) == A._multiply_classical(B) True
sage: K.<a> = GF(2^10) sage: A = random_matrix(K, 50, 50) sage: B = random_matrix(K, 50, 50) sage: A._multiply_karatsuba(B) == A._multiply_classical(B) True """ raise ArithmeticError("left ncols must match right nrows")
cdef Matrix_gf2e_dense ans
return ans
cpdef Matrix_gf2e_dense _multiply_strassen(Matrix_gf2e_dense self, Matrix_gf2e_dense right, cutoff=0): """ Winograd-Strassen matrix multiplication with Newton-John multiplication as base case.
INPUT:
- ``right`` - a matrix - ``cutoff`` - row or column dimension to switch over to Newton-John multiplication (default: 64)
EXAMPLES::
sage: K.<a> = GF(2^2) sage: A = random_matrix(K, 50, 50) sage: B = random_matrix(K, 50, 50) sage: A._multiply_strassen(B) == A._multiply_classical(B) # indirect doctest True
sage: K.<a> = GF(2^4) sage: A = random_matrix(K, 50, 50) sage: B = random_matrix(K, 50, 50) sage: A._multiply_strassen(B) == A._multiply_classical(B) True
sage: K.<a> = GF(2^8) sage: A = random_matrix(K, 50, 50) sage: B = random_matrix(K, 50, 50) sage: A._multiply_strassen(B) == A._multiply_classical(B) True
sage: K.<a> = GF(2^10) sage: A = random_matrix(K, 50, 50) sage: B = random_matrix(K, 50, 50) sage: A._multiply_strassen(B) == A._multiply_classical(B) True """ raise ArithmeticError("left ncols must match right nrows")
cdef Matrix_gf2e_dense ans
return ans
cpdef _lmul_(self, Element right): """ Return ``a*B`` for ``a`` an element of the base field.
INPUT:
- ``right`` - an element of the base field
EXAMPLES::
sage: K.<a> = GF(4) sage: A = random_matrix(K,10,10) sage: A [ 0 a + 1 a + 1 a + 1 0 1 a + 1 1 a + 1 1] [a + 1 a + 1 a 1 a a 1 a + 1 1 0] [ a 1 a + 1 a + 1 0 a + 1 a 1 a a] [a + 1 a 0 0 1 a + 1 a + 1 0 a + 1 1] [ a 0 a + 1 a a 0 a + 1 a 1 a + 1] [ a 0 a a + 1 a 1 a + 1 a a a] [a + 1 a 0 1 0 a + 1 a + 1 a 0 a + 1] [a + 1 a + 1 0 a + 1 a 1 a + 1 a + 1 a + 1 0] [ 0 0 0 a + 1 1 a + 1 0 a + 1 1 0] [ 1 a + 1 0 1 a 0 0 a a + 1 a]
sage: a*A # indirect doctest [ 0 1 1 1 0 a 1 a 1 a] [ 1 1 a + 1 a a + 1 a + 1 a 1 a 0] [a + 1 a 1 1 0 1 a + 1 a a + 1 a + 1] [ 1 a + 1 0 0 a 1 1 0 1 a] [a + 1 0 1 a + 1 a + 1 0 1 a + 1 a 1] [a + 1 0 a + 1 1 a + 1 a 1 a + 1 a + 1 a + 1] [ 1 a + 1 0 a 0 1 1 a + 1 0 1] [ 1 1 0 1 a + 1 a 1 1 1 0] [ 0 0 0 1 a 1 0 1 a 0] [ a 1 0 a a + 1 0 0 a + 1 1 a + 1] """
def __neg__(self): """ EXAMPLES::
sage: K.<a> = GF(2^4) sage: A = random_matrix(K, 3, 4); A [ a^2 a^3 + a + 1 a^3 + a^2 + a + 1 a + 1] [ a^3 1 a^3 + a + 1 a^3 + a^2 + 1] [ a + 1 a^3 + 1 a^3 + a + 1 a^3 + a^2 + a + 1]
sage: -A [ a^2 a^3 + a + 1 a^3 + a^2 + a + 1 a + 1] [ a^3 1 a^3 + a + 1 a^3 + a^2 + 1] [ a + 1 a^3 + 1 a^3 + a + 1 a^3 + a^2 + a + 1] """
cpdef int _cmp_(self, right) except -2: """ EXAMPLES::
sage: K.<a> = GF(2^4) sage: A = random_matrix(K,3,4) sage: B = copy(A) sage: A == B True sage: A[0,0] = a sage: A == B False """ if self._nrows == 0 or self._ncols == 0: return 0 return mzed_cmp(self._entries, (<Matrix_gf2e_dense>right)._entries)
def __copy__(self): """ EXAMPLES::
sage: K.<a> = GF(2^4) sage: m,n = 3, 4 sage: A = random_matrix(K,3,4); A [ a^2 a^3 + a + 1 a^3 + a^2 + a + 1 a + 1] [ a^3 1 a^3 + a + 1 a^3 + a^2 + 1] [ a + 1 a^3 + 1 a^3 + a + 1 a^3 + a^2 + a + 1]
sage: A2 = copy(A); A2 [ a^2 a^3 + a + 1 a^3 + a^2 + a + 1 a + 1] [ a^3 1 a^3 + a + 1 a^3 + a^2 + 1] [ a + 1 a^3 + 1 a^3 + a + 1 a^3 + a^2 + a + 1]
sage: A[0,0] = 1 sage: A2[0,0] a^2 """ cdef Matrix_gf2e_dense A
def _list(self): """ EXAMPLES::
sage: K.<a> = GF(2^4) sage: m,n = 3, 4 sage: A = random_matrix(K,3,4); A [ a^2 a^3 + a + 1 a^3 + a^2 + a + 1 a + 1] [ a^3 1 a^3 + a + 1 a^3 + a^2 + 1] [ a + 1 a^3 + 1 a^3 + a + 1 a^3 + a^2 + a + 1]
sage: A.list() # indirect doctest [a^2, a^3 + a + 1, a^3 + a^2 + a + 1, a + 1, a^3, 1, a^3 + a + 1, a^3 + a^2 + 1, a + 1, a^3 + 1, a^3 + a + 1, a^3 + a^2 + a + 1] """ cdef int i,j
def randomize(self, density=1, nonzero=False, *args, **kwds): """ 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-place
EXAMPLES::
sage: K.<a> = GF(2^4) sage: A = Matrix(K,3,3)
sage: A.randomize(); A [ a^2 a^3 + a + 1 a^3 + a^2 + a + 1] [ a + 1 a^3 1] [ a^3 + a + 1 a^3 + a^2 + 1 a + 1]
sage: K.<a> = GF(2^4) sage: A = random_matrix(K,1000,1000,density=0.1) sage: float(A.density()) 0.0999...
sage: A = random_matrix(K,1000,1000,density=1.0) sage: float(A.density()) 1.0
sage: A = random_matrix(K,1000,1000,density=0.5) sage: float(A.density()) 0.4996...
Note, that the matrix is updated and not zero-ed out before being randomized::
sage: A = matrix(K,1000,1000) sage: A.randomize(nonzero=False, density=0.1) sage: float(A.density()) 0.0936...
sage: A.randomize(nonzero=False, density=0.05) sage: float(A.density()) 0.135854 """
cdef Py_ssize_t i,j
return
_density = 1.0
else: else: else:
def echelonize(self, algorithm='heuristic', reduced=True, **kwds): """ Compute the row echelon form of ``self`` in place.
INPUT:
- ``algorithm`` - one of the following - ``heuristic`` - let M4RIE decide (default) - ``newton_john`` - use newton_john table based algorithm - ``ple`` - use PLE decomposition - ``naive`` - use naive cubic Gaussian elimination (M4RIE implementation) - ``builtin`` - use naive cubic Gaussian elimination (Sage implementation) - ``reduced`` - if ``True`` return reduced echelon form. No guarantee is given that the matrix is *not* reduced if ``False`` (default: ``True``)
EXAMPLES::
sage: K.<a> = GF(2^4) sage: m,n = 3, 5 sage: A = random_matrix(K, 3, 5); A [ a^2 a^3 + a + 1 a^3 + a^2 + a + 1 a + 1 a^3] [ 1 a^3 + a + 1 a^3 + a^2 + 1 a + 1 a^3 + 1] [ a^3 + a + 1 a^3 + a^2 + a + 1 a^2 + a a^2 + 1 a^2 + a]
sage: A.echelonize(); A [ 1 0 0 a + 1 a^2 + 1] [ 0 1 0 a^2 a + 1] [ 0 0 1 a^3 + a^2 + a a^3]
sage: K.<a> = GF(2^3) sage: m,n = 3, 5 sage: MS = MatrixSpace(K,m,n) sage: A = random_matrix(K, 3, 5)
sage: copy(A).echelon_form('newton_john') [ 1 0 a + 1 0 a^2 + 1] [ 0 1 a^2 + a + 1 0 a] [ 0 0 0 1 a^2 + a + 1]
sage: copy(A).echelon_form('naive'); [ 1 0 a + 1 0 a^2 + 1] [ 0 1 a^2 + a + 1 0 a] [ 0 0 0 1 a^2 + a + 1]
sage: copy(A).echelon_form('builtin'); [ 1 0 a + 1 0 a^2 + 1] [ 0 1 a^2 + a + 1 0 a] [ 0 0 0 1 a^2 + a + 1] """ self.cache('in_echelon_form',True) self.cache('rank', 0) self.cache('pivots', []) return self
cdef int k, n, full
sig_on() r = mzed_echelonize_ple(self._entries, full) sig_off()
else: raise ValueError("No algorithm '%s'."%algorithm)
def _pivots(self): """ EXAMPLES::
sage: K.<a> = GF(2^8) sage: A = random_matrix(K, 15, 15) sage: A.pivots() # indirect doctest (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14) """ raise ValueError("self must be in reduced row echelon form first.") cdef Py_ssize_t i, j, nc
def __invert__(self): """ EXAMPLES::
sage: K.<a> = GF(2^3) sage: A = random_matrix(K,3,3); A [ a^2 a + 1 a^2 + a + 1] [ a + 1 0 1] [ a + 1 a^2 + 1 a + 1]
sage: B = ~A; B [a^2 + a + 1 a^2 a^2] [ a + 1 a^2 + a + 1 a + 1] [ a a^2 + a a^2 + a + 1]
sage: A*B [1 0 0] [0 1 0] [0 0 1] """ cdef Matrix_gf2e_dense A
cdef rescale_row_c(self, Py_ssize_t row, multiple, Py_ssize_t start_col): """ Return ``multiple * self[row][start_col:]``
INPUT:
- ``row`` - row index for row to rescale - ``multiple`` - finite field element to scale by - ``start_col`` - only start at this column index.
EXAMPLES::
sage: K.<a> = GF(2^3) sage: A = random_matrix(K,3,3); A [ a^2 a + 1 a^2 + a + 1] [ a + 1 0 1] [ a + 1 a^2 + 1 a + 1]
sage: A.rescale_row(0, a , 0); A [ a + 1 a^2 + a a^2 + 1] [ a + 1 0 1] [ a + 1 a^2 + 1 a + 1]
sage: A.rescale_row(0,0,0); A [ 0 0 0] [ a + 1 0 1] [ a + 1 a^2 + 1 a + 1] """
cdef add_multiple_of_row_c(self, Py_ssize_t row_to, Py_ssize_t row_from, multiple, Py_ssize_t start_col): """ Compute ``self[row_to][start_col:] += multiple * self[row_from][start_col:]``.
INPUT:
- ``row_to`` - row index of source - ``row_from`` - row index of destination - ``multiple`` - finite field element - ``start_col`` - only start at this column index
EXAMPLES::
sage: K.<a> = GF(2^3) sage: A = random_matrix(K,3,3); A [ a^2 a + 1 a^2 + a + 1] [ a + 1 0 1] [ a + 1 a^2 + 1 a + 1]
sage: A.add_multiple_of_row(0,1,a,0); A [ a a + 1 a^2 + 1] [ a + 1 0 1] [ a + 1 a^2 + 1 a + 1] """
cdef swap_rows_c(self, Py_ssize_t row1, Py_ssize_t row2): """ Swap rows ``row1`` and ``row2``.
INPUT:
- ``row1`` - row index - ``row2`` - row index
EXAMPLES::
sage: K.<a> = GF(2^3) sage: A = random_matrix(K,3,3) sage: A [ a^2 a + 1 a^2 + a + 1] [ a + 1 0 1] [ a + 1 a^2 + 1 a + 1]
sage: A.swap_rows(0,1); A [ a + 1 0 1] [ a^2 a + 1 a^2 + a + 1] [ a + 1 a^2 + 1 a + 1]
"""
cdef swap_columns_c(self, Py_ssize_t col1, Py_ssize_t col2): """ Swap columns ``col1`` and ``col2``.
INPUT:
- ``col1`` - column index - ``col2`` - column index
EXAMPLES::
sage: K.<a> = GF(2^3) sage: A = random_matrix(K,3,3) sage: A [ a^2 a + 1 a^2 + a + 1] [ a + 1 0 1] [ a + 1 a^2 + 1 a + 1]
sage: A.swap_columns(0,1); A [ a + 1 a^2 a^2 + a + 1] [ 0 a + 1 1] [ a^2 + 1 a + 1 a + 1]
sage: A = random_matrix(K,4,16)
sage: B = copy(A) sage: B.swap_columns(0,1) sage: B.swap_columns(0,1) sage: A == B True
sage: A.swap_columns(0,15) sage: A.column(0) == B.column(15) True sage: A.swap_columns(14,15) sage: A.column(14) == B.column(0) True """
def augment(self, Matrix_gf2e_dense right): """ Augments ``self`` with ``right``.
INPUT:
- ``right`` - a matrix
EXAMPLES::
sage: K.<a> = GF(2^4) sage: MS = MatrixSpace(K,3,3) sage: A = random_matrix(K,3,3) sage: B = A.augment(MS(1)); B [ a^2 a^3 + a + 1 a^3 + a^2 + a + 1 1 0 0] [ a + 1 a^3 1 0 1 0] [ a^3 + a + 1 a^3 + a^2 + 1 a + 1 0 0 1]
sage: B.echelonize(); B [ 1 0 0 a^2 + a a^3 + 1 a^3 + a] [ 0 1 0 a^3 + a^2 + a a^3 + a^2 + a + 1 a^2 + a] [ 0 0 1 a + 1 a^3 a^3]
sage: C = B.matrix_from_columns([3,4,5]); C [ a^2 + a a^3 + 1 a^3 + a] [ a^3 + a^2 + a a^3 + a^2 + a + 1 a^2 + a] [ a + 1 a^3 a^3]
sage: C == ~A True
sage: C*A == MS(1) True
TESTS::
sage: K.<a> = GF(2^4) sage: A = random_matrix(K,2,3) sage: B = random_matrix(K,2,0) sage: A.augment(B) [ a^3 + 1 a^3 + a + 1 a^3 + a^2 + a + 1] [ a^2 + a a^2 + 1 a^2 + a]
sage: B.augment(A) [ a^3 + 1 a^3 + a + 1 a^3 + a^2 + a + 1] [ a^2 + a a^2 + 1 a^2 + a]
sage: M = Matrix(K, 0, 0, 0) sage: N = Matrix(K, 0, 19, 0) sage: W = M.augment(N) sage: W.ncols() 19
sage: M = Matrix(K, 0, 1, 0) sage: N = Matrix(K, 0, 1, 0) sage: M.augment(N) [] """ cdef Matrix_gf2e_dense A
raise TypeError("Both numbers of rows must match.")
def stack(self, Matrix_gf2e_dense other): """ Stack ``self`` on top of ``other``.
INPUT:
- ``other`` - a matrix
EXAMPLES::
sage: K.<a> = GF(2^4) sage: A = random_matrix(K,2,2); A [ a^2 a^3 + a + 1] [a^3 + a^2 + a + 1 a + 1]
sage: B = random_matrix(K,2,2); B [ a^3 1] [ a^3 + a + 1 a^3 + a^2 + 1]
sage: A.stack(B) [ a^2 a^3 + a + 1] [a^3 + a^2 + a + 1 a + 1] [ a^3 1] [ a^3 + a + 1 a^3 + a^2 + 1]
sage: B.stack(A) [ a^3 1] [ a^3 + a + 1 a^3 + a^2 + 1] [ a^2 a^3 + a + 1] [a^3 + a^2 + a + 1 a + 1]
TESTS::
sage: A = random_matrix(K,0,3) sage: B = random_matrix(K,3,3) sage: A.stack(B) [ a + 1 a^3 + 1 a^3 + a + 1] [a^3 + a^2 + a + 1 a^2 + a a^2 + 1] [ a^2 + a a^3 + a^2 + a a^2 + 1]
sage: B.stack(A) [ a + 1 a^3 + 1 a^3 + a + 1] [a^3 + a^2 + a + 1 a^2 + a a^2 + 1] [ a^2 + a a^3 + a^2 + a a^2 + 1]
sage: M = Matrix(K, 0, 0, 0) sage: N = Matrix(K, 19, 0, 0) sage: W = M.stack(N) sage: W.nrows() 19 sage: M = Matrix(K, 1, 0, 0) sage: N = Matrix(K, 1, 0, 0) sage: M.stack(N) [] """ raise TypeError("Both numbers of columns must match.")
cdef Matrix_gf2e_dense A
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: K.<a> = GF(2^10) sage: A = random_matrix(K,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 """
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 rank(self): """ Return the rank of this matrix (cached).
EXAMPLES::
sage: K.<a> = GF(2^4) sage: A = random_matrix(K, 1000, 1000) sage: A.rank() 1000
sage: A = matrix(K, 10, 0) sage: A.rank() 0 """
def __reduce__(self): """ EXAMPLES::
sage: K.<a> = GF(2^8) sage: A = random_matrix(K,70,70) sage: f, s= A.__reduce__() sage: from sage.matrix.matrix_gf2e_dense import unpickle_matrix_gf2e_dense_v0 sage: f == unpickle_matrix_gf2e_dense_v0 True sage: f(*s) == A True
See :trac:`21669`::
sage: all(f(*s) == B ....: for r,c in [(0,0),(0,1),(1,0)] ....: for B in [Matrix(GF(4, 'a'), r,c)] ....: for f,s in [B.__reduce__()]) True """
cdef Matrix_mod2_dense A cdef int r,c
def slice(self): r""" Unpack this matrix into matrices over `\GF{2}`.
Elements in `\GF{2^e}` can be represented as `\sum c_i a^i` where `a` is a root the minimal polynomial. This function returns a tuple of matrices `C` whose entry `C_i[x,y]` is the coefficient of `c_i` in `A[x,y]` if this matrix is `A`.
EXAMPLES::
sage: K.<a> = GF(2^2) sage: A = random_matrix(K, 5, 5); A [ 0 a + 1 a + 1 a + 1 0] [ 1 a + 1 1 a + 1 1] [a + 1 a + 1 a 1 a] [ a 1 a + 1 1 0] [ a 1 a + 1 a + 1 0]
sage: A1,A0 = A.slice() sage: A0 [0 1 1 1 0] [0 1 0 1 0] [1 1 1 0 1] [1 0 1 0 0] [1 0 1 1 0]
sage: A1 [0 1 1 1 0] [1 1 1 1 1] [1 1 0 1 0] [0 1 1 1 0] [0 1 1 1 0]
sage: A0[2,4]*a + A1[2,4], A[2,4] (a, a)
sage: K.<a> = GF(2^3) sage: A = random_matrix(K, 5, 5); A [ a + 1 a^2 + a 1 a a^2 + a] [ a + 1 a^2 + a a^2 a^2 a^2 + 1] [a^2 + a + 1 a^2 + a + 1 0 a^2 + a + 1 a^2 + 1] [ a^2 + a 0 a^2 + a + 1 a a] [ a^2 a + 1 a a^2 + 1 a^2 + a + 1]
sage: A0,A1,A2 = A.slice() sage: A0 [1 0 1 0 0] [1 0 0 0 1] [1 1 0 1 1] [0 0 1 0 0] [0 1 0 1 1]
Slicing and clinging are inverse operations::
sage: B = matrix(K, 5, 5) sage: B.cling(A0,A1,A2) sage: B == A True """ raise NotImplementedError("Slicing is only implemented for degree <= 4.")
cdef Matrix_mod2_dense A0, A1, A2, A3 A3 = Matrix_mod2_dense.__new__(Matrix_mod2_dense, MS, 0, 0, 0, alloc = True) mzd_copy(A3._entries, a.x[3])
elif self._entries.finite_field.degree == 3: elif self._entries.finite_field.degree == 4: return A0,A1,A2,A3
def cling(self, *C): r""" Pack the matrices over `\GF{2}` into this matrix over `\GF{2^e}`.
Elements in `\GF{2^e}` can be represented as `\sum c_i a^i` where `a` is a root the minimal polynomial. If this matrix is `A` then this function writes `c_i a^i` to the entry `A[x,y]` where `c_i` is the entry `C_i[x,y]`.
INPUT:
- ``C`` - a list of matrices over GF(2)
EXAMPLES::
sage: K.<a> = GF(2^2) sage: A = matrix(K, 5, 5) sage: A0 = random_matrix(GF(2), 5, 5); A0 [0 1 0 1 1] [0 1 1 1 0] [0 0 0 1 0] [0 1 1 0 0] [0 0 0 1 1]
sage: A1 = random_matrix(GF(2), 5, 5); A1 [0 0 1 1 1] [1 1 1 1 0] [0 0 0 1 1] [1 0 0 0 1] [1 0 0 1 1]
sage: A.cling(A1, A0); A [ 0 a 1 a + 1 a + 1] [ 1 a + 1 a + 1 a + 1 0] [ 0 0 0 a + 1 1] [ 1 a a 0 1] [ 1 0 0 a + 1 a + 1]
sage: A0[0,3]*a + A1[0,3], A[0,3] (a + 1, a + 1)
Slicing and clinging are inverse operations::
sage: B1, B0 = A.slice() sage: B0 == A0 and B1 == A1 True
TESTS::
sage: K.<a> = GF(2^2) sage: A = matrix(K, 5, 5) sage: A0 = random_matrix(GF(2), 5, 5) sage: A1 = random_matrix(GF(2), 5, 5) sage: A.cling(A0, A1) sage: B = copy(A) sage: A.cling(A0, A1) sage: A == B True
sage: A.cling(A0) Traceback (most recent call last): ... ValueError: The number of input matrices must be equal to the degree of the base field.
sage: K.<a> = GF(2^5) sage: A = matrix(K, 5, 5) sage: A0 = random_matrix(GF(2), 5, 5) sage: A1 = random_matrix(GF(2), 5, 5) sage: A2 = random_matrix(GF(2), 5, 5) sage: A3 = random_matrix(GF(2), 5, 5) sage: A4 = random_matrix(GF(2), 5, 5) sage: A.cling(A0, A1, A2, A3, A4) Traceback (most recent call last): ... NotImplementedError: Cling is only implemented for degree <= 4. """ cdef Py_ssize_t i
raise TypeError("Immutable matrices cannot be modified.")
mzd_slice_free(v) raise TypeError("All input matrices must be over GF(2).")
def unpickle_matrix_gf2e_dense_v0(Matrix_mod2_dense a, base_ring, nrows, ncols): r""" EXAMPLES::
sage: K.<a> = GF(2^2) sage: A = random_matrix(K,10,10) sage: f, s= A.__reduce__() sage: from sage.matrix.matrix_gf2e_dense import unpickle_matrix_gf2e_dense_v0 sage: f == unpickle_matrix_gf2e_dense_v0 True sage: f(*s) == A True
We can still unpickle pickles from before :trac:`19240`::
sage: old_pickle = 'x\x9c\x85RKo\xd3@\x10\xae\xdd$$\xdb&\xe5U\x1e-\x8f\xc2\xc9\x12RD#$\xce\xa0\xb4\x80\x07\xa2\xca\xc2\x07\x0e\xd5\xe2:\x1b\xdb\x8acg\x1c\xa7J\x85*!\xa4\x90\xe6\x07p\xe0\xc4\x01q\xe5\xc4\x19\xf5\xd0?\xc1\x81\xdf\x80\xb8q\x0b\xb3\x8eMS\xa1\x82V;;\xb3\xdf\xce\xf7\xcd\x8e\xe6\xb5j\xf7,GT;V\x1cy\x83\xf4\xe0\x9d\xb0Y\x13\xbc)\x82\x9e`\xfd\xa0\xeb\xd9m_\xf0\xbf1\xbe{\x97\xa1\xa2\x9d\xc6\xf0\x0f\x82,\x7f\x9d\xa1\xaa\x81\n\xb9m\x9c\xd7\xf4\xf1d2\x81-h\xc0#(\x03\x83\x15\xdas\xc9*\xc3\x13x\x0cu0\xd28\x97\x9e*(0\x9f\xfa\x1b\xd0\xd2\x7fH\x82\xb5\xf4\xa2@TO\xe19\x01I\xac\x136\x991\x9f\xa4\xf9&\xcd\x07i\xbe\xcb\xd4ib\t\xba\xa4\xf6\x02zIT\xd1\x8f2(u\x15\xfd\x9d<\xee@\x05V\xd3\x94E*\xb0\x0e\x0fH\xad\xa8\xbf\x97\xa0\r\x03\xfd\xf0\xb8\x1aU\xff\x92\x90\xe8?\xa5\xd6\x814_\xa5\xf9(\xcd\xafc\xe99\xe2\xd9\xa0\x06\xd4\xf5\xcf\xf2\xf2!\xbc\xd4\xdf\x90#\xc0\x8f\r\xccM\x1b\xdd\x8b\xa3\xbe\x1d\xf7#QmYv\x1cF{\xcc\x11\x81\x88<\x9b\xa71\xcf:\xce0\xaf\x9d\x96\xe3\x87a\xbb\xdf\xe5\x8e\x1f\xeeX>\xc3\x82\xb9\xb0\xe9\x05^,6=\xe17\xf1\xcc\xd0\xc0"u\xb0d\xe6wDl\xdd\x1fa)e\x8a\xbc\xc0\xe9U\xbd \x16\x8e\x88X\xc7j\x0b\x9e\x05\xc8L\xe5\x1e%.\x98\x8a5\xc4\xc5\xd9\xf7\xdd\xd0\xdf\x0b\xc2\x8eg\xf93.wZ\xb5\xc1\x94B\xf8\xa2#\x82\x98a\xf9\xffY\x12\xe3v\x18L\xff\x14Fl\xeb\x0ff\x10\xc4\xb0\xa2\xb9y\xcd-\xba%\xcd\xa5\x8ajT\xd1\x92\xa9\x0c\x86x\xb6a\xe6h\xf8\x02<g\xaa\xaf\xf6\xdd%\x89\xae\x13z\xfe \xc6\x0b\xfb1^4p\x99\x1e6\xc6\xd4\xebK\xdbx\xf9\xc4\x8f[Iw\xf8\x89\xef\xcbQf\xcfh\xe3\x95\x8c\xebj&\xb9\xe2.\x8f\x0c\\ui\x89\xf1x\xf4\xd6\xc0kf\xc1\xf1v\xad(\xc4\xeb\x89~\xfa\xf0\x06\xa8\xa4\x7f\x93\xf4\xd7\x0c\xbcE#\xad\x92\xfc\xed\xeao\xefX\\\x03' sage: loads(old_pickle) [ 0 a] [a + 1 1] """
|