Dense matrices over Z/nZ for n<28 using LinBox’s Modular<float>

AUTHORS: - Burcin Erocal - Martin Albrecht

class sage.matrix.matrix_modn_dense_float.Matrix_modn_dense_float[source]

Bases: Matrix_modn_dense_template

Dense matrices over Z/nZ for n<28 using LinBox’s Modular<float>.

These are matrices with integer entries mod n represented as floating-point numbers in a 32-bit word for use with LinBox routines. This could allow for n up to 211, but for performance reasons this is limited to n up to 28, and for larger moduli the Matrix_modn_dense_double class is used.

Routines here are for the most basic access, see the matrix_modn_dense_template.pxi file for higher-level routines.

class sage.matrix.matrix_modn_dense_float.Matrix_modn_dense_template[source]

Bases: Matrix_dense

Create a new matrix.

INPUT:

  • parent – a matrix space

  • entries – see matrix()

  • copy – ignored (for backwards compatibility)

  • coerce – perform modular reduction first?

EXAMPLES:

sage: A = random_matrix(GF(3),1000,1000)
sage: type(A)
<class 'sage.matrix.matrix_modn_dense_float.Matrix_modn_dense_float'>
sage: A = random_matrix(Integers(10),1000,1000)
sage: type(A)
<class 'sage.matrix.matrix_modn_dense_float.Matrix_modn_dense_float'>
sage: A = random_matrix(Integers(2^16),1000,1000)
sage: type(A)
<class 'sage.matrix.matrix_modn_dense_double.Matrix_modn_dense_double'>
charpoly(var='x', algorithm='linbox')[source]

Return the characteristic polynomial of self.

INPUT:

  • var – a variable name

  • algorithm – ‘generic’, ‘linbox’ or ‘all’ (default: linbox)

EXAMPLES:

sage: A = random_matrix(GF(19), 10, 10)
sage: B = copy(A)
sage: char_p = A.characteristic_polynomial()
sage: char_p(A) == 0
True
sage: B == A              # A is not modified
True

sage: min_p = A.minimal_polynomial(proof=True)
sage: min_p.divides(char_p)
True
sage: A = random_matrix(GF(2916337), 7, 7)                                  # needs sage.rings.finite_rings
sage: B = copy(A)
sage: char_p = A.characteristic_polynomial()
sage: char_p(A) == 0
True
sage: B == A               # A is not modified
True

sage: min_p = A.minimal_polynomial(proof=True)
sage: min_p.divides(char_p)
True

sage: A = Mat(Integers(6),3,3)(range(9))
sage: A.charpoly()
x^3

ALGORITHM: Uses LinBox if self.base_ring() is a field, otherwise use Hessenberg form algorithm.

determinant()[source]

Return the determinant of this matrix.

EXAMPLES:

sage: s = set()
sage: while s != set(GF(7)):
....:     A = random_matrix(GF(7), 10, 10)
....:     s.add(A.determinant())
sage: A = random_matrix(GF(7), 100, 100)
sage: A.determinant() == A.transpose().determinant()
True

sage: B = random_matrix(GF(7), 100, 100)
sage: (A*B).determinant() == A.determinant() * B.determinant()
True
sage: # needs sage.rings.finite_rings
sage: A = random_matrix(GF(16007), 10, 10)
sage: A.determinant().parent() is GF(16007)
True
sage: # needs sage.rings.finite_rings
sage: A = random_matrix(GF(16007), 100, 100)
sage: A.determinant().parent() is GF(16007)
True
sage: A.determinant() == A.transpose().determinant()
True
sage: B = random_matrix(GF(16007), 100, 100)
sage: (A*B).determinant() == A.determinant() * B.determinant()
True

Parallel computation:

sage: # needs sage.rings.finite_rings
sage: A = random_matrix(GF(65521),200)
sage: B = copy(A)
sage: Parallelism().set('linbox', nproc=2)
sage: d = A.determinant()
sage: Parallelism().set('linbox', nproc=1) # switch off parallelization
sage: e = B.determinant()
sage: d==e
True
echelonize(algorithm='linbox_noefd', **kwds)[source]

Put self in reduced row echelon form.

INPUT:

  • self – a mutable matrix

  • algorithm

    • linbox – uses the LinBox library (wrapping fflas-ffpack)

    • linbox_noefd – uses the FFPACK directly, less memory and faster (default)

    • gauss – uses a custom slower O(n3) Gauss elimination implemented in Sage

    • all – compute using all algorithms and verify that the results are the same

  • **kwds – these are all ignored

OUTPUT: if self is known to be echelonized (the information is cached), nothing is done; else, self is put in reduced row echelon form and

  • the fact that self is now in echelon form is cached so future calls to echelonize return immediately

  • the rank of self is computed and cached

  • the pivot columns (a.k.a. column rank profile) of self are computed and cached

  • only if using algorithm linbox_noefd: the pivot rows of self before echelonization (a.k.a. its row rank profile) are computed and cached

EXAMPLES:

sage: A = random_matrix(GF(7), 10, 20)
sage: E = A.echelon_form()
sage: A.row_space() == E.row_space()
True
sage: all(r[r.nonzero_positions()[0]] == 1 for r in E.rows() if r)
True
sage: A = random_matrix(GF(13), 10, 10)
sage: while A.rank() != 10:
....:     A = random_matrix(GF(13), 10, 10)
sage: MS = parent(A)
sage: B = A.augment(MS(1))
sage: B.echelonize()
sage: A.rank()
10
sage: C = B.submatrix(0,10,10,10)
sage: ~A == C
True

sage: B = A.augment(MS(1))
sage: B.echelonize(algorithm='linbox')
sage: A.rank()
10
sage: C = B.submatrix(0,10,10,10)
sage: ~A == C
True
sage: A = random_matrix(Integers(10), 10, 20)
sage: A.echelon_form()
Traceback (most recent call last):
...
NotImplementedError: Echelon form not implemented over 'Ring of integers modulo 10'.
sage: # needs sage.rings.finite_rings
sage: A = random_matrix(GF(16007), 10, 20)
sage: E = A.echelon_form()
sage: A.row_space() == E.row_space()
True
sage: all(r[r.nonzero_positions()[0]] == 1 for r in E.rows() if r)
True
sage: A = random_matrix(Integers(10000), 10, 20)
sage: A.echelon_form()
Traceback (most recent call last):
...
NotImplementedError: Echelon form not implemented over 'Ring of integers modulo 10000'.

Parallel computation:

sage: # needs sage.rings.finite_rings
sage: A = random_matrix(GF(65521),100,200)
sage: Parallelism().set('linbox', nproc=2)
sage: E = A.echelon_form()
sage: Parallelism().set('linbox', nproc=1) # switch off parallelization
sage: F = A.echelon_form()
sage: E==F
True
hessenbergize()[source]

Transform self in place to its Hessenberg form.

EXAMPLES:

sage: A = random_matrix(GF(17), 10, 10, density=0.1)
sage: B = copy(A)
sage: A.hessenbergize()
sage: all(A[i,j] == 0 for j in range(10) for i in range(j+2, 10))
True
sage: A.charpoly() == B.charpoly()
True
lift()[source]

Return the lift of this matrix to the integers.

EXAMPLES:

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

sage: # needs sage.rings.finite_rings
sage: A = matrix(GF(16007),2,3,[1..6])
sage: A.lift()
[1 2 3]
[4 5 6]
sage: A.lift().parent()
Full MatrixSpace of 2 by 3 dense 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(columns)[source]

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

EXAMPLES:

sage: M = MatrixSpace(Integers(8),3,3)
sage: A = M(range(9)); A
[0 1 2]
[3 4 5]
[6 7 0]
sage: A.matrix_from_columns([2,1])
[2 1]
[5 4]
[0 7]
matrix_from_rows(rows)[source]

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

EXAMPLES:

sage: M = MatrixSpace(Integers(8),3,3)
sage: A = M(range(9)); A
[0 1 2]
[3 4 5]
[6 7 0]
sage: A.matrix_from_rows([2,1])
[6 7 0]
[3 4 5]
matrix_from_rows_and_columns(rows, columns)[source]

Return the matrix constructed from self from the given rows and columns.

EXAMPLES:

sage: M = MatrixSpace(Integers(8),3,3)
sage: A = M(range(9)); A
[0 1 2]
[3 4 5]
[6 7 0]
sage: A.matrix_from_rows_and_columns([1], [0,2])
[3 5]
sage: A.matrix_from_rows_and_columns([1,2], [1,2])
[4 5]
[7 0]

Note that row and column indices can be reordered or repeated:

sage: A.matrix_from_rows_and_columns([2,1], [2,1])
[0 7]
[5 4]

For example here we take from row 1 columns 2 then 0 twice, and do this 3 times:

sage: A.matrix_from_rows_and_columns([1,1,1],[2,0,0])
[5 3 3]
[5 3 3]
[5 3 3]

AUTHORS:

  • Jaap Spies (2006-02-18)

  • Didier Deshommes: some Pyrex speedups implemented

minpoly(var='x', algorithm='linbox', proof=None)[source]

Return the minimal polynomial of self.

INPUT:

  • var – a variable name

  • algorithmgeneric or linbox (default: linbox)

  • proof – (default: True) whether to provably return the true minimal polynomial; if False, we only guarantee to return a divisor of the minimal polynomial. There are also certainly cases where the computed results is frequently not exactly equal to the minimal polynomial (but is instead merely a divisor of it).

Warning

If proof=True, minpoly is insanely slow compared to proof=False. This matters since proof=True is the default, unless you first type proof.linear_algebra(False).

EXAMPLES:

sage: A = random_matrix(GF(17), 10, 10)
sage: B = copy(A)
sage: min_p = A.minimal_polynomial(proof=True)
sage: min_p(A) == 0
True
sage: B == A
True

sage: char_p = A.characteristic_polynomial()
sage: min_p.divides(char_p)
True
sage: A = random_matrix(GF(1214471), 10, 10)                                # needs sage.rings.finite_rings
sage: B = copy(A)
sage: min_p = A.minimal_polynomial(proof=True)
sage: min_p(A) == 0
True
sage: B == A
True

sage: char_p = A.characteristic_polynomial()
sage: min_p.divides(char_p)
True

EXAMPLES:

sage: R.<x>=GF(3)[]
sage: A = matrix(GF(3),2,[0,0,1,2])
sage: A.minpoly()
x^2 + x

sage: A.minpoly(proof=False) in [x, x+1, x^2+x]
True
pivot_rows()[source]

Return the row pivot positions for this matrix, which form the topmost subset of the rows that span the row space and are linearly independent. This coincides with the position of the first nonzero entry in each column of the reduced column echelon form of self, and is also known as the row rank profile of self. The returned tuple is ordered increasingly.

If this has already been computed and cached, this returns the cached value. Otherwise, this computes an echelon form (using algorithm linbox_noefd), and deduces the pivot row positions to be cached and returned. In the latter case, this also caches other attributes at the same time: the reduced row echelon form of self and the rank of self, as well as its pivot indices (also known as column rank profile).

OUTPUT: a tuple of r integers where r is the rank of self

See also

The method Matrix_modn_dense_template.pivots() computes the column pivot positions, also known as column rank profile.

EXAMPLES:

sage: A = matrix(GF(3), [[1,0,1,0],[1,0,0,0],[1,0,0,0],[0,1,0,0]])
sage: A
[1 0 1 0]
[1 0 0 0]
[1 0 0 0]
[0 1 0 0]
sage: A.pivot_rows() == (0, 1, 3)
True

sage: A = matrix(GF(65537), 3, 5,
....:            [  223, 13996, 39034, 14660, 12016,
....:               669, 41988, 51565, 43980, 36048,
....:             21130, 21387, 40500,  3899,  9308])
sage: A[1,:] == 3 * A[0,:]
True
sage: A.pivot_rows() == (0, 2)
True
sage: A = matrix(GF(7), 5, 3, [2, 4, 5,
....:                          2, 4, 2,
....:                          4, 1, 5,
....:                          1, 4, 6,
....:                          4, 6, 6])
sage: A[0,:] == 2 * A[1,:] + 3 * A[2,:]
True
sage: A.pivot_rows() == (0, 1, 3)
True
pivots()[source]

Return the column pivot positions for this matrix, which form the leftmost subset of the columns that span the column space and are linearly independent. This coincides with the position of the first nonzero entry in each row of the reduced row echelon form of self, and is also known as the column rank profile of self. The returned tuple is ordered increasingly.

If this has already been computed and cached, this returns the cached value. Otherwise, this computes an echelon form (using algorithm linbox_noefd), and deduces the pivot positions to be cached and returned. In the latter case, this also caches other attributes at the same time: the reduced row echelon form of self and the rank of self, as well as its row pivot positions (also known as row rank profile).

OUTPUT: a tuple of r integers where r is the rank of self

See also

The method Matrix_modn_dense_template.pivot_rows() computes the row pivot positions, also known as row rank profile.

EXAMPLES:

sage: A = matrix(GF(7), 2, 2, range(4))
sage: A.pivots()
(0, 1)

sage: A = matrix(GF(3), [[1,1,1,0],[0,0,0,1],[1,0,0,0]])
sage: A
[1 1 1 0]
[0 0 0 1]
[1 0 0 0]
sage: A.pivots() == (0, 1, 3)
True

sage: A = matrix(GF(65537), 5, 3,
....:            [  223,   669, 21130,
....:             13996, 41988, 21387,
....:             39034, 51565, 40500,
....:             14660, 43980,  3899,
....:             12016, 36048,  9308])
sage: A[:,1] == 3 * A[:,0]
True
sage: A.pivots() == (0, 2)
True
sage: A = matrix(GF(7), 3, 5, [2, 2, 4, 1, 4,
....:                          4, 4, 1, 4, 6,
....:                          5, 2, 5, 6, 6])
sage: A[:,0] == 2 * A[:,1] + 3 * A[:,2]
True
sage: A.pivots() == (0, 1, 3)
True
randomize(density=1, nonzero=False)[source]

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

INPUT:

  • density – integer; proportion (roughly) to be considered for changes

  • nonzero – boolean (default: False); whether the new entries are forced to be nonzero

OUTPUT: none, the matrix is modified in-space

EXAMPLES:

sage: A = matrix(GF(5), 5, 5, 0)
sage: total_count = 0
sage: from collections import defaultdict
sage: dic = defaultdict(Integer)
sage: def add_samples(density):
....:     global dic, total_count
....:     for _ in range(100):
....:         A = Matrix(GF(5), 5, 5, 0)
....:         A.randomize(density)
....:         for a in A.list():
....:             dic[a] += 1
....:             total_count += 1.0

sage: add_samples(1.0)
sage: while not all(abs(dic[a]/total_count - 1/5) < 0.01 for a in dic):
....:     add_samples(1.0)

sage: def add_sample(density):
....:     global density_sum, total_count
....:     total_count += 1.0
....:     density_sum += random_matrix(GF(5), 1000, 1000, density=density).density()

sage: density_sum = 0.0
sage: total_count = 0.0
sage: add_sample(0.5)
sage: expected_density = 1.0 - (999/1000)^500
sage: expected_density
0.3936...
sage: while abs(density_sum/total_count - expected_density) > 0.001:
....:     add_sample(0.5)

The matrix is updated instead of overwritten:

sage: def add_sample(density):
....:     global density_sum, total_count
....:     total_count += 1.0
....:     A = random_matrix(GF(5), 1000, 1000, density=density)
....:     A.randomize(density=density, nonzero=True)
....:     density_sum += A.density()

sage: density_sum = 0.0
sage: total_count = 0.0
sage: add_sample(0.5)
sage: expected_density = 1.0 - (999/1000)^1000
sage: expected_density
0.6323...
sage: while abs(density_sum/total_count - expected_density) > 0.001:
....:     add_sample(0.5)

sage: density_sum = 0.0
sage: total_count = 0.0
sage: add_sample(0.1)
sage: expected_density = 1.0 - (999/1000)^200
sage: expected_density
0.1813...
sage: while abs(density_sum/total_count - expected_density) > 0.001:
....:     add_sample(0.1)
rank()[source]

Return the rank of this matrix.

EXAMPLES:

sage: A = random_matrix(GF(3), 100, 100)
sage: B = copy(A)
sage: _ = A.rank()
sage: B == A
True

sage: A = random_matrix(GF(3), 100, 100, density=0.01)
sage: A.transpose().rank() == A.rank()
True

sage: A = matrix(GF(3), 100, 100)
sage: A.rank()
0

Rank is not implemented over the integers modulo a composite yet.:

sage: M = matrix(Integers(4), 2, [2,2,2,2])
sage: M.rank()
Traceback (most recent call last):
...
NotImplementedError: Echelon form not implemented over 'Ring of integers modulo 4'.
sage: # needs sage.rings.finite_rings
sage: while True:
....:     A = random_matrix(GF(16007), 100, 100)
....:     if A.rank() == 100:
....:         break
sage: B = copy(A)
sage: A.rank()
100
sage: B == A
True
sage: MS = A.parent()
sage: MS(1) == ~A*A
True
right_kernel_matrix(algorithm='linbox', basis='echelon')[source]

Return a matrix whose rows form a basis for the right kernel of self.

If the base ring is the ring of integers modulo a composite, the keyword arguments are ignored and the computation is delegated to Matrix_dense.right_kernel_matrix().

INPUT:

  • algorithm – (default: 'linbox') a parameter that is passed on to self.echelon_form, if computation of an echelon form is required; see that routine for allowable values

  • basis – (default: 'echelon') a keyword that describes the format of the basis returned, allowable values are:

    • 'echelon': the basis matrix is in echelon form

    • 'pivot': the basis matrix is such that the submatrix obtained

      by taking the columns that in self contain no pivots, is the identity matrix

    • 'computed': no work is done to transform the basis; in

      the current implementation the result is the negative of that returned by 'pivot'

OUTPUT:

A matrix X whose rows are a basis for the right kernel of self. This means that self * X.transpose() is a zero matrix.

The result is not cached, but the routine benefits when self is known to be in echelon form already.

EXAMPLES:

sage: M = matrix(GF(5),6,6,range(36))
sage: M.right_kernel_matrix(basis='computed')
[4 2 4 0 0 0]
[3 3 0 4 0 0]
[2 4 0 0 4 0]
[1 0 0 0 0 4]
sage: M.right_kernel_matrix(basis='pivot')
[1 3 1 0 0 0]
[2 2 0 1 0 0]
[3 1 0 0 1 0]
[4 0 0 0 0 1]
sage: M.right_kernel_matrix()
[1 0 0 0 0 4]
[0 1 0 0 1 3]
[0 0 1 0 2 2]
[0 0 0 1 3 1]
sage: M * M.right_kernel_matrix().transpose()
[0 0 0 0]
[0 0 0 0]
[0 0 0 0]
[0 0 0 0]
[0 0 0 0]
[0 0 0 0]
submatrix(row=0, col=0, nrows=-1, ncols=-1)[source]

Return the matrix constructed from self using the specified range of rows and columns.

INPUT:

  • row, col – index of the starting row and column; indices start at zero

  • nrows, ncols – (optional) number of rows and columns to take; if not provided, take all rows below and all columns to the right of the starting entry

See also

The functions matrix_from_rows(), matrix_from_columns(), and matrix_from_rows_and_columns() allow one to select arbitrary subsets of rows and/or columns.

EXAMPLES:

Take the 3×3 submatrix starting from entry (1,1) in a 4×4 matrix:

sage: m = matrix(GF(17),4, [1..16])
sage: m.submatrix(1, 1)
[ 6  7  8]
[10 11 12]
[14 15 16]

Same thing, except take only two rows:

sage: m.submatrix(1, 1, 2)
[ 6  7  8]
[10 11 12]

And now take only one column:

sage: m.submatrix(1, 1, 2, 1)
[ 6]
[10]

You can take zero rows or columns if you want:

sage: m.submatrix(0, 0, 0)
[]
sage: parent(m.submatrix(0, 0, 0))
Full MatrixSpace of 0 by 4 dense matrices over Finite Field of size 17
transpose()[source]

Return the transpose of self, without changing self.

EXAMPLES:

We create a matrix, compute its transpose, and note that the original matrix is not changed.

sage: M = MatrixSpace(GF(41),  2)
sage: A = M([1,2,3,4])
sage: B = A.transpose()
sage: B
[1 3]
[2 4]
sage: A
[1 2]
[3 4]

.T is a convenient shortcut for the transpose:

sage: A.T
[1 3]
[2 4]
sage: A.subdivide(None, 1); A
[1|2]
[3|4]
sage: A.transpose()
[1 3]
[---]
[2 4]