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
""" Sparse rational matrices.
AUTHORS:
- William Stein (2007-02-21) - Soroosh Yazdani (2007-02-21)
TESTS::
sage: a = matrix(QQ,2,range(4), sparse=True) sage: TestSuite(a).run() sage: matrix(QQ,0,0,sparse=True).inverse() [] """
#***************************************************************************** # Copyright (C) 2007 William Stein <wstein@gmail.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_on, sig_off from cysignals.memory cimport sig_malloc, sig_free
from collections import Iterator, Sequence
from sage.data_structures.binary_search cimport * from sage.modules.vector_integer_sparse cimport * from sage.modules.vector_rational_sparse cimport *
from cpython.sequence cimport *
from sage.rings.rational cimport Rational from sage.rings.integer cimport Integer from .matrix cimport Matrix
from sage.libs.gmp.mpz cimport * from sage.libs.gmp.mpq cimport *
from sage.libs.flint.fmpq cimport fmpq_set_mpq from sage.libs.flint.fmpq_mat cimport fmpq_mat_entry
from sage.rings.integer_ring import ZZ from sage.rings.rational_field import QQ
cimport sage.structure.element
import sage.matrix.matrix_space
from .matrix_integer_sparse cimport Matrix_integer_sparse from .matrix_rational_dense cimport Matrix_rational_dense
from sage.misc.misc import verbose
cdef class Matrix_rational_sparse(Matrix_sparse): def __cinit__(self, parent, entries, copy, coerce): # set the parent, nrows, ncols, etc.
raise MemoryError("error allocating sparse matrix") # initialize the rows
# record that rows have been initialized
def __dealloc__(self):
cdef _dealloc(self): cdef Py_ssize_t i
def __init__(self, parent, entries, copy, coerce): """ Create a sparse matrix over the rational numbers.
INPUT:
- ``parent`` -- a matrix space
- ``entries`` -- can be one of the following:
* a Python dictionary whose items have the form ``(i, j): x``, where ``0 <= i < nrows``, ``0 <= j < ncols``, and ``x`` is coercible to a rational. The ``i,j`` entry of ``self`` is set to ``x``. The ``x``'s can be ``0``. * Alternatively, entries can be a list of *all* the entries of the sparse matrix, read row-by-row from top to bottom (so they would be mostly 0).
- ``copy`` -- ignored
- ``coerce`` -- ignored """ cdef Py_ssize_t i, j, k cdef Rational z cdef PyObject** X
# fill in entries in the dict case # Dense input format -- fill in entries # Get fast access to the entries list. z = R(<object>X[k]) if z != 0: mpq_vector_set_entry(&self._matrix[i], j, z.value) k = k + 1 else: # fill in entries in the scalar case
cdef set_unsafe(self, Py_ssize_t i, Py_ssize_t j, x):
cdef get_unsafe(self, Py_ssize_t i, Py_ssize_t j): cdef Rational x
def add_to_entry(self, Py_ssize_t i, Py_ssize_t j, elt): r""" Add ``elt`` to the entry at position ``(i, j)``.
EXAMPLES::
sage: m = matrix(QQ, 2, 2, sparse=True) sage: m.add_to_entry(0, 0, -1/3) sage: m [-1/3 0] [ 0 0] sage: m.add_to_entry(0, 0, 1/3) sage: m [0 0] [0 0] sage: m.nonzero_positions() [] """ elt = Rational(elt) i += self._nrows raise IndexError("row index out of range") j += self._ncols raise IndexError("column index out of range")
cdef mpq_t z
######################################################################## # LEVEL 2 functionality # * def _pickle # * def _unpickle # * cdef _add_ # * cdef _sub_ # * cdef _mul_ # * cpdef _cmp_ # * __neg__ # * __invert__ # * __copy__ # * _multiply_classical # * _matrix_times_matrix_ # * _list -- list of underlying elements (need not be a copy) # * x _dict -- sparse dictionary of underlying elements (need not be a copy)
cdef sage.structure.element.Matrix _matrix_times_matrix_(self, sage.structure.element.Matrix _right): cdef Matrix_rational_sparse right, ans
cdef mpq_vector* v
# Build a table that gives the nonzero positions in each column of right cdef Py_ssize_t i, j, k
# Now do the multiplication, getting each row completely before filling it in. cdef mpq_t x, y, s
def _matrix_times_matrix_dense(self, sage.structure.element.Matrix _right): """ Do the sparse matrix multiply, but return a dense matrix as the result.
EXAMPLES::
sage: a = matrix(QQ, 2, [1,2,3,4], sparse=True) sage: b = matrix(QQ, 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_rational_dense.Matrix_rational_dense'> """ cdef Matrix_rational_sparse right cdef Matrix_rational_dense ans
cdef mpq_vector* v
# Build a table that gives the nonzero positions in each column of right cdef Py_ssize_t i, j, k
# Now do the multiplication, getting each row completely before filling it in. cdef mpq_t x, y, s
######################################################################## # def _pickle(self): # def _unpickle(self, data, int version): # use version >= 0 # cpdef _add_(self, right): # cdef _mul_(self, Matrix right): # cpdef int _cmp_(self, Matrix right) except -2: # def __neg__(self): # def __invert__(self): # def __copy__(self): # def _multiply_classical(left, matrix.Matrix _right): # def _list(self):
# TODO ## cpdef _lmul_(self, Element right): ## """ ## EXAMPLES: ## sage: a = matrix(QQ,2,range(6)) ## sage: (3/4) * a ## [ 0 3/4 3/2] ## [ 9/4 3 15/4] ## """
def _dict(self): """ 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. """
cdef Py_ssize_t i, j, k
######################################################################## # LEVEL 3 functionality (Optional) # * cdef _sub_ # * __deepcopy__ # * __invert__ # * Matrix windows -- only if you need strassen for that base # * Other functions (list them here): ########################################################################
def _nonzero_positions_by_row(self, copy=True): """ 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).
EXAMPLES::
sage: M = Matrix(QQ, [[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)]
""" cdef Py_ssize_t i, j
def height(self): """ Return the height of this matrix, which is the least common multiple of all numerators and denominators of elements of this matrix.
OUTPUT:
-- Integer
EXAMPLES::
sage: b = matrix(QQ,2,range(6), sparse=True); b[0,0]=-5007/293; b [-5007/293 1 2] [ 3 4 5] sage: b.height() 5007 """
cdef int mpz_height(self, mpz_t height) except -1: cdef mpz_t x, h cdef int i, j
cdef int mpz_denom(self, mpz_t d) except -1: cdef Py_ssize_t i, j cdef mpq_vector *v
def denominator(self): """ Return the denominator of this matrix.
OUTPUT:
-- Sage Integer
EXAMPLES::
sage: b = matrix(QQ,2,range(6)); b[0,0]=-5007/293; b [-5007/293 1 2] [ 3 4 5] sage: b.denominator() 293 """
def _clear_denom(self): """ INPUT:
self -- a matrix
OUTPUT:
D*self, D
The product D*self is a matrix over ZZ
EXAMPLES::
sage: a = matrix(QQ,3,[-2/7, -1/4, -2, 0, 1/7, 1, 0, 1/2, 1/5],sparse=True) sage: a.denominator() 140 sage: a._clear_denom() ( [ -40 -35 -280] [ 0 20 140] [ 0 70 28], 140 ) """ cdef Integer D cdef Py_ssize_t i, j cdef Matrix_integer_sparse A cdef mpz_t t cdef mpq_vector* v
################################################ # Echelon form ################################################ def echelonize(self, height_guess=None, proof=True, **kwds): """ Transform the matrix ``self`` into reduced row echelon form in place.
INPUT:
``height_guess``, ``proof``, ``**kwds`` -- all passed to the multimodular algorithm; ignored by the p-adic algorithm.
OUTPUT:
Nothing. The matrix ``self`` is transformed into reduced row echelon form in place.
ALGORITHM: a multimodular algorithm.
EXAMPLES::
sage: a = matrix(QQ, 4, range(16), sparse=True); a[0,0] = 1/19; a[0,1] = 1/5; a [1/19 1/5 2 3] [ 4 5 6 7] [ 8 9 10 11] [ 12 13 14 15] sage: a.echelonize(); a [ 1 0 0 -76/157] [ 0 1 0 -5/157] [ 0 0 1 238/157] [ 0 0 0 0]
:trac:`10319` has been fixed::
sage: m = Matrix(QQ, [1], sparse=True); m.echelonize() sage: m = Matrix(QQ, [1], sparse=True); m.echelonize(); m [1] """
def echelon_form(self, algorithm='default', height_guess=None, proof=True, **kwds): """ INPUT:
``height_guess``, ``proof``, ``**kwds`` -- all passed to the multimodular algorithm; ignored by the p-adic algorithm.
OUTPUT:
self is no in reduced row echelon form.
EXAMPLES::
sage: a = matrix(QQ, 4, range(16), sparse=True); a[0,0] = 1/19; a[0,1] = 1/5; a [1/19 1/5 2 3] [ 4 5 6 7] [ 8 9 10 11] [ 12 13 14 15] sage: a.echelon_form() [ 1 0 0 -76/157] [ 0 1 0 -5/157] [ 0 0 1 238/157] [ 0 0 0 0] """
# Multimodular echelonization algorithms def _echelonize_multimodular(self, height_guess=None, proof=True, **kwds): cdef Py_ssize_t i, j cdef Matrix_rational_sparse E cdef mpq_vector* v cdef mpq_vector* w # Get rid of self's data # Copy E's data to self's data. raise MemoryError("error allocating sparse matrix")
def _echelon_form_multimodular(self, height_guess=None, proof=True): """ Returns reduced row-echelon form using a multi-modular algorithm. Does not change self.
INPUT:
- height_guess -- integer or None - proof -- boolean (default: True) """ cdef Matrix E
def set_row_to_multiple_of_row(self, i, j, s): """ Set row i equal to s times row j.
EXAMPLES::
sage: a = matrix(QQ,2,3,range(6), sparse=True); a [0 1 2] [3 4 5] sage: a.set_row_to_multiple_of_row(1,0,-3) sage: a [ 0 1 2] [ 0 -3 -6] """ cdef Rational _s
def dense_matrix(self): """ Return dense version of this matrix.
EXAMPLES::
sage: a = matrix(QQ,2,[1..4],sparse=True); type(a) <type 'sage.matrix.matrix_rational_sparse.Matrix_rational_sparse'> sage: type(a.dense_matrix()) <type 'sage.matrix.matrix_rational_dense.Matrix_rational_dense'> sage: a.dense_matrix() [1 2] [3 4]
Check that subdivisions are preserved when converting between dense and sparse matrices::
sage: a.subdivide([1,1], [2]) sage: b = a.dense_matrix().sparse_matrix().dense_matrix() sage: b.subdivisions() == a.subdivisions() True """ cdef Matrix_rational_dense B cdef mpq_vector* v
## def _set_row_to_negative_of_row_of_A_using_subset_of_columns(self, Py_ssize_t i, Matrix A, ## Py_ssize_t r, cols): ## B = self.__copy__() ## B.x_set_row_to_negative_of_row_of_A_using_subset_of_columns(i, A, r, cols) ## cdef Py_ssize_t l ## l = 0 ## for z in range(self.ncols()): ## self[i,z] = 0 ## for k in cols: ## self.set_unsafe(i,l,-A.get_unsafe(r,k)) #self[i,l] = -A[r,k] ## l += 1 ## if self != B: ## print("correct =\n", self.str()) ## print("wrong = \n", B.str()) ## print("diff = \n", (self-B).str())
def _set_row_to_negative_of_row_of_A_using_subset_of_columns(self, Py_ssize_t i, Matrix A, Py_ssize_t r, cols, cols_index=None): """ Set row i of self to -(row r of A), but where we only take the given column positions in that row of A. Note that we *DO* zero out the other entries of self's row i.
INPUT:
- i -- integer, index into the rows of self - A -- a sparse matrix - r -- integer, index into rows of A - cols -- a *sorted* list of integers. - cols_index -- (optional). But set it to this to vastly speed up calls to this function::
dict([(cols[i], i) for i in range(len(cols))])
EXAMPLES::
sage: a = matrix(QQ,2,3,range(6), sparse=True); a [0 1 2] [3 4 5]
Note that the row is zeroed out before being set in the sparse case. ::
sage: a._set_row_to_negative_of_row_of_A_using_subset_of_columns(0,a,1,[1,2]) sage: a [-4 -5 0] [ 3 4 5] """ # this function exists just because it is useful for modular symbols presentations. raise IndexError("invalid row")
A = A.sparse_matrix()
A = A.change_ring(QQ)
cdef Matrix_rational_sparse _A
cdef Py_ssize_t l, n
cdef mpq_vector *v cdef mpq_vector *w
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-iml-rational' that identifies the nature of the basis vectors.
Second item is a matrix whose rows are a basis for the right kernel, over the rationals, as computed by the IML library. Notice that the IML library returns a matrix that is in the 'pivot' format, once the whole matrix is multiplied by -1. So the 'computed' format is very close to the 'pivot' format.
EXAMPLES::
sage: A = matrix(QQ, [ ....: [1, 0, 1, -3, 1], ....: [-5, 1, 0, 7, -3], ....: [0, -1, -4, 6, -2], ....: [4, -1, 0, -6, 2]], ....: sparse=True) sage: result = A._right_kernel_matrix() sage: result[0] 'computed-iml-rational' sage: result[1] [-1 2 -2 -1 0] [ 1 2 0 0 -1] sage: X = result[1].transpose() sage: A*X == zero_matrix(QQ, 4, 2) True
Computed result is the negative of the pivot basis, which is just slightly more efficient to compute. ::
sage: A.right_kernel_matrix(basis='pivot') == -A.right_kernel_matrix(basis='computed') True
TESTS:
We test three trivial cases. ::
sage: A = matrix(QQ, 0, 2, sparse=True) sage: A._right_kernel_matrix()[1] [1 0] [0 1] sage: A = matrix(QQ, 2, 0, sparse=True) sage: A._right_kernel_matrix()[1].parent() Full MatrixSpace of 0 by 0 dense matrices over Rational Field sage: A = zero_matrix(QQ, 4, 3, sparse=True) sage: A._right_kernel_matrix()[1] [1 0 0] [0 1 0] [0 0 1] """
#########################
# used for a function above cdef mpq_t minus_one mpq_init(minus_one) mpq_set_si(minus_one, -1,1)
|