Dense matrices over \(\ZZ/n\ZZ\) for \(n < 2^{23}\) using LinBox’s Modular<double>

AUTHORS:

• Burcin Erocal

• Martin Albrecht

class sage.matrix.matrix_modn_dense_double.Matrix_modn_dense_double

Bases: sage.matrix.matrix_modn_dense_double.Matrix_modn_dense_template

Dense matrices over \(\ZZ/n\ZZ\) for \(n < 2^{23}\) using LinBox’s Modular<double>

These are matrices with integer entries mod n represented as floating-point numbers in a 64-bit word for use with LinBox routines. This allows for n up to \(2^{23}\). The analogous Matrix_modn_dense_float class is used for smaller moduli.

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_double.Matrix_modn_dense_template

Bases: sage.matrix.matrix_dense.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)
<type 'sage.matrix.matrix_modn_dense_float.Matrix_modn_dense_float'>
sage: A = random_matrix(Integers(10),1000,1000)
sage: type(A)
<type 'sage.matrix.matrix_modn_dense_float.Matrix_modn_dense_float'>
sage: A = random_matrix(Integers(2^16),1000,1000)
sage: type(A)
<type 'sage.matrix.matrix_modn_dense_double.Matrix_modn_dense_double'>

charpoly(var='x', algorithm='linbox')

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)
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()

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: A = random_matrix(GF(16007), 10, 10)
sage: A.determinant().parent() is GF(16007)
True

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: 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)

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(n^3)\) Gauss elimination implemented in Sage.

  • all - compute using both algorithms and verify that the results are the same.

• **kwds - these are all ignored

OUTPUT:

• self is put in reduced row echelon form.

• the rank of self is computed and cached

• the pivot columns of self are computed and cached.

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

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: 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: 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: 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()

Transforms 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()

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: 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]

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

Returns the minimal polynomial of`` self``.

INPUT:

• var - a variable name

• algorithm - generic 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)
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

randomize(density=1, nonzero=False)

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

INPUT:

• density - Integer; 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(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()

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: A = random_matrix(GF(16007), 100, 100)
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')

Returns a matrix whose rows form a basis for the right kernel of self, where self is a matrix over a (small) finite field.

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)

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 \times 3\) submatrix starting from entry \((1,1)\) in a \(4 \times 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()

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]