Dense matrices over the rational field

EXAMPLES:

We create a 3x3 matrix with rational entries and do some operations with it.

sage: a = matrix(QQ, 3,3, [1,2/3, -4/5, 1,1,1, 8,2, -3/19]); a
[    1   2/3  -4/5]
[    1     1     1]
[    8     2 -3/19]
sage: a.det()
2303/285
sage: a.charpoly()
x^3 - 35/19*x^2 + 1259/285*x - 2303/285
sage: b = a^(-1); b
[ -615/2303  -426/2303   418/2303]
[ 2325/2303  1779/2303  -513/2303]
[-1710/2303   950/2303    95/2303]
sage: b.det()
285/2303
sage: a == b
False
sage: a < b
False
sage: b < a
True
sage: a > b
True
sage: a*b
[1 0 0]
[0 1 0]
[0 0 1]
class sage.matrix.matrix_rational_dense.MatrixWindow

Bases: object

class sage.matrix.matrix_rational_dense.Matrix_rational_dense

Bases: sage.matrix.matrix_dense.Matrix_dense

INPUT:

  • parent – a matrix space over QQ

  • entries – see matrix()

  • copy – ignored (for backwards compatibility)

  • coerce – if False, assume without checking that the entries are of type Rational.

LLL(*args, **kwargs)

Return an LLL reduced or approximated LLL reduced lattice for self interpreted as a lattice.

For details on input parameters, see sage.matrix.matrix_integer_dense.Matrix_integer_dense.LLL().

EXAMPLES:

sage: A = Matrix(QQ, 3, 3, [1/n for n in range(1, 10)])
sage: A.LLL()
[ 1/28 -1/40 -1/18]
[ 1/28 -1/40  1/18]
[    0 -3/40     0]
add_to_entry(i, j, elt)

Add elt to the entry at position (i,j)

EXAMPLES:

sage: m = matrix(QQ, 2, 2)
sage: m.add_to_entry(0, 0, -1/3)
sage: m
[-1/3    0]
[   0    0]
antitranspose()

Returns the antitranspose of self, without changing self.

EXAMPLES:

sage: A = matrix(QQ,2,3,range(6))
sage: type(A)
<type 'sage.matrix.matrix_rational_dense.Matrix_rational_dense'>
sage: A.antitranspose()
[5 2]
[4 1]
[3 0]
sage: A
[0 1 2]
[3 4 5]

sage: A.subdivide(1,2); A
[0 1|2]
[---+-]
[3 4|5]
sage: A.antitranspose()
[5|2]
[-+-]
[4|1]
[3|0]
change_ring(R)

Create the matrix over R with entries the entries of self coerced into R.

EXAMPLES:

sage: a = matrix(QQ,2,[1/2,-1,2,3])
sage: a.change_ring(GF(3))
[2 2]
[2 0]
sage: a.change_ring(ZZ)
Traceback (most recent call last):
...
TypeError: matrix has denominators so can...t change to ZZ
sage: b = a.change_ring(QQ['x']); b
[1/2  -1]
[  2   3]
sage: b.parent()
Full MatrixSpace of 2 by 2 dense matrices over Univariate Polynomial Ring in x over Rational Field
charpoly(var='x', algorithm=None)

Return the characteristic polynomial of this matrix.

Note

The characteristic polynomial is defined as \(\det(xI-A)\).

INPUT:

  • var - (optional) name of the variable as a string

  • algorithm – an optional specification of an algorithm. It can be one of:

    • None: (default) will use flint for small dimensions and linbox otherwise

    • 'flint': uses flint library

    • 'linbox': uses linbox library

    • 'generic': uses Sage generic implementation

OUTPUT: a polynomial over the rational numbers.

EXAMPLES:

sage: a = matrix(QQ, 3, [4/3, 2/5, 1/5, 4, -3/2, 0, 0, -2/3, 3/4])
sage: f = a.charpoly(); f
x^3 - 7/12*x^2 - 149/40*x + 97/30
sage: f(a)
[0 0 0]
[0 0 0]
[0 0 0]
column(i, from_list=False)

Return the i-th column of this matrix as a dense vector.

INPUT:

  • i - integer

  • from_list - ignored

EXAMPLES:

sage: m = matrix(QQ, 3, 2, [1/5,-2/3,3/4,4/9,-1,0])
sage: m.column(1)
(-2/3, 4/9, 0)
sage: m.column(1,from_list=True)
(-2/3, 4/9, 0)
sage: m.column(-1)
(-2/3, 4/9, 0)
sage: m.column(-2)
(1/5, 3/4, -1)

sage: m.column(2)
Traceback (most recent call last):
...
IndexError: column index out of range
sage: m.column(-3)
Traceback (most recent call last):
...
IndexError: column index out of range
decomposition(is_diagonalizable=False, dual=False, algorithm=None, height_guess=None, proof=None)

Returns the decomposition of the free module on which this matrix A acts from the right (i.e., the action is x goes to x A), along with whether this matrix acts irreducibly on each factor. The factors are guaranteed to be sorted in the same way as the corresponding factors of the characteristic polynomial.

Let A be the matrix acting from the on the vector space V of column vectors. Assume that A is square. This function computes maximal subspaces W_1, …, W_n corresponding to Galois conjugacy classes of eigenvalues of A. More precisely, let f(X) be the characteristic polynomial of A. This function computes the subspace \(W_i = ker(g_(A)^n)\), where g_i(X) is an irreducible factor of f(X) and g_i(X) exactly divides f(X). If the optional parameter is_diagonalizable is True, then we let W_i = ker(g(A)), since then we know that ker(g(A)) = \(ker(g(A)^n)\).

If dual is True, also returns the corresponding decomposition of V under the action of the transpose of A. The factors are guaranteed to correspond.

INPUT:

  • is_diagonalizable - ignored

  • dual - whether to also return decompositions for the dual

  • algorithm - an optional specification of an algorithm

    • None - (default) use default algorithm for computing Echelon forms

    • ‘multimodular’: much better if the answers factors have small height

  • height_guess - positive integer; only used by the multimodular algorithm

  • proof - bool or None (default: None, see proof.linear_algebra or sage.structure.proof); only used by the multimodular algorithm. Note that the Sage global default is proof=True.

Note

IMPORTANT: If you expect that the subspaces in the answer are spanned by vectors with small height coordinates, use algorithm=’multimodular’ and height_guess=1; this is potentially much faster than the default. If you know for a fact the answer will be very small, use algorithm=’multimodular’, height_guess=bound on height, proof=False.

You can get very very fast decomposition with proof=False.

EXAMPLES:

sage: a = matrix(QQ,3,[1..9])
sage: a.decomposition()
[
(Vector space of degree 3 and dimension 1 over Rational Field
Basis matrix:
[ 1 -2  1], True),
(Vector space of degree 3 and dimension 2 over Rational Field
Basis matrix:
[ 1  0 -1]
[ 0  1  2], True)
]
denominator()

Return the denominator of this matrix.

OUTPUT: a 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

sage: matrix(QQ, 2, [1/2, 1/3, 1/4, 1/5]).denominator()
60
determinant(algorithm=None, proof=None)

Return the determinant of this matrix.

INPUT:

  • algorithm – an optional specification of an algorithm. It can be one of

    • None: (default) uses flint

    • 'flint': uses flint library

    • 'pari': uses PARI library

    • 'integer': removes denominator and call determinant on the corresponding

      integer matrix

    • 'generic': calls the generic Sage implementation

  • proof - bool or None; if None use proof.linear_algebra(); only relevant for the padic algorithm.

Note

It would be VERY VERY hard for det to fail even with proof=False.

EXAMPLES:

sage: m = matrix(QQ,3,[1,2/3,4/5, 2,2,2, 5,3,2/5])
sage: m.determinant()
-34/15
sage: m.charpoly()
x^3 - 17/5*x^2 - 122/15*x + 34/15

sage: m = matrix(QQ, 3, [(1/i)**j for i in range(2,5) for j in range(3)])
sage: m.determinant(algorithm="flint")
-1/288

sage: m = matrix(QQ, 4, [(-1)**n/n for n in range(1,17)])
sage: m.determinant(algorithm="pari")
2/70945875

sage: m = matrix(QQ, 5, [1/(i+j+1) for i in range(5) for j in range(5)])
sage: m.determinant(algorithm="integer")
1/266716800000

On non-square matrices, the method raises a ValueError:

sage: matrix(QQ, 2, 3).determinant(algorithm='flint')
Traceback (most recent call last):
...
ValueError: non square matrix
sage: matrix(QQ, 2, 3).determinant(algorithm='pari')
Traceback (most recent call last):
...
ValueError: non square matrix
sage: matrix(QQ, 2, 3).determinant(algorithm='integer')
Traceback (most recent call last):
...
ValueError: non square matrix
sage: matrix(QQ, 2, 3).determinant(algorithm='generic')
Traceback (most recent call last):
...
ValueError: non square matrix
echelon_form(algorithm=None, height_guess=None, proof=None, **kwds)

Return the echelon form of this matrix.

The (row) echelon form of a matrix, see Wikipedia article Row_echelon_form, is the matrix obtained by performing Gauss elimination on the rows of the matrix.

INPUT: See echelonize() for the options.

EXAMPLES:

sage: a = matrix(QQ, 4, range(16)); 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]
sage: a.echelon_form(algorithm='multimodular')
[      1       0       0 -76/157]
[      0       1       0  -5/157]
[      0       0       1 238/157]
[      0       0       0       0]

The result is an immutable matrix, so if you want to modify the result then you need to make a copy. This checks that trac ticket #10543 is fixed.:

sage: A = matrix(QQ, 2, range(6))
sage: E = A.echelon_form()
sage: E.is_mutable()
False
sage: F = copy(E)
sage: F[0,0] = 50
sage: F
[50  0 -1]
[ 0  1  2]
echelonize(algorithm=None, height_guess=None, proof=None, **kwds)

Transform the matrix self into reduced row echelon form in place.

INPUT:

  • algorithm – an optional specification of an algorithm. One of

  • None: (default) uses flint for small dimension and multimodular otherwise

  • 'flint': use the flint library,

  • 'padic': an algorithm based on the IML p-adic solver,

  • 'multimodular': uses a multimodular algorithm the uses linbox modulo many primes (likely to be faster when coefficients are huge),

  • 'classical': just clear each column using Gauss elimination.

  • height_guess, **kwds - all passed to the multimodular algorithm; ignored by other algorithms.

  • proof - bool or None (default: None, see proof.linear_algebra or sage.structure.proof). Passed to the multimodular algorithm. Note that the Sage global default is proof=True.

EXAMPLES:

sage: a = matrix(QQ, 4, range(16)); 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()
sage: a
[      1       0       0 -76/157]
[      0       1       0  -5/157]
[      0       0       1 238/157]
[      0       0       0       0]
sage: a = matrix(QQ, 4, range(16)); a[0,0] = 1/19; a[0,1] = 1/5
sage: a.echelonize(algorithm='multimodular')
sage: a
[      1       0       0 -76/157]
[      0       1       0  -5/157]
[      0       0       1 238/157]
[      0       0       0       0]
height()

Return the height of this matrix, which is the maximum of the absolute values of all numerators and denominators of entries in this matrix.

OUTPUT: an Integer

EXAMPLES:

sage: b = matrix(QQ,2,range(6)); b[0,0]=-5007/293; b
[-5007/293         1         2]
[        3         4         5]
sage: b.height()
5007
inverse(algorithm=None, check_invertible=True)

Return the inverse of this matrix

INPUT:

  • algorithm – an optional specification of an algorithm. It can be one of

    • None: (default) uses flint

    • 'flint': uses flint library

    • 'pari': uses PARI library

    • 'iml': uses IML library

  • check_invertible - only used when algorithm=iml. Whether to check that matrix is invertible

EXAMPLES:

sage: a = matrix(QQ,3,[1,2,5,3,2,1,1,1,1,])
sage: a.inverse()
[1/2 3/2  -4]
[ -1  -2   7]
[1/2 1/2  -2]

sage: a = matrix(QQ, 2, [1, 5, 17, 3])
sage: a.inverse(algorithm="flint")
[-3/82  5/82]
[17/82 -1/82]
sage: a.inverse(algorithm="flint")  * a
[1 0]
[0 1]

sage: a = matrix(QQ, 2, [-1, 5, 12, -3])
sage: a.inverse(algorithm="iml")
[1/19 5/57]
[4/19 1/57]
sage: a.inverse(algorithm="iml") * a
[1 0]
[0 1]

sage: a = matrix(QQ, 4, primes_first_n(16))
sage: a.inverse(algorithm="pari")
[   3/11  -12/55    -1/5    2/11]
[  -5/11   -2/55    3/10   -3/22]
[ -13/22 307/440   -1/10   -9/88]
[  15/22  -37/88       0    7/88]

On singular matrices this method raises a ZeroDivisionError:

sage: a = matrix(QQ, 2)
sage: a.inverse(algorithm="flint")
Traceback (most recent call last):
...
ZeroDivisionError: input matrix must be nonsingular
sage: a.inverse(algorithm="iml")
Traceback (most recent call last):
...
ZeroDivisionError: input matrix must be nonsingular
sage: a.inverse(algorithm="pari")
Traceback (most recent call last):
...
ZeroDivisionError: input matrix must be nonsingular
matrix_from_columns(columns)

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

EXAMPLES:

sage: A = matrix(QQ, 3, range(9))
sage: A
[0 1 2]
[3 4 5]
[6 7 8]
sage: A.matrix_from_columns([2,1])
[2 1]
[5 4]
[8 7]
sage: A.matrix_from_columns((2,1,0,2))
[2 1 0 2]
[5 4 3 5]
[8 7 6 8]
minpoly(var='x', algorithm=None)

Return the minimal polynomial of this matrix

INPUT:

  • var - (optional) the variable name as a string (default is ‘x’)

  • algorithm - an optional specification of an algorithm. It can be one of

    • None: (default) will use linbox

    • 'linbox': uses the linbox library

    • 'generic': uses the generic Sage implementation

OUTPUT: a polynomial over the rationals

EXAMPLES:

sage: a = matrix(QQ, 3, [4/3, 2/5, 1/5, 4, -3/2, 0, 0, -2/3, 3/4])
sage: f = a.minpoly(); f
x^3 - 7/12*x^2 - 149/40*x + 97/30
sage: a = Mat(ZZ,4)(range(16))
sage: f = a.minpoly(); f.factor()
x * (x^2 - 30*x - 80)
sage: f(a) == 0
True
sage: a = matrix(QQ, 4, [1..4^2])
sage: factor(a.minpoly())
x * (x^2 - 34*x - 80)
sage: factor(a.minpoly('y'))
y * (y^2 - 34*y - 80)
sage: factor(a.charpoly())
x^2 * (x^2 - 34*x - 80)
sage: b = matrix(QQ, 4, [-1, 2, 2, 0, 0, 4, 2, 2, 0, 0, -1, -2, 0, -4, 0, 4])
sage: a = matrix(QQ, 4, [1, 1, 0,0, 0,1,0,0, 0,0,5,0, 0,0,0,5])
sage: c = b^(-1)*a*b
sage: factor(c.minpoly())
(x - 5) * (x - 1)^2
sage: factor(c.charpoly())
(x - 5)^2 * (x - 1)^2

Check consistency:

sage: for _ in range(100):
....:     dim = randint(0, 10)
....:     m = random_matrix(QQ, dim, num_bound=8, den_bound=8)
....:     p_linbox = m.charpoly(algorithm='linbox'); m._clear_cache()
....:     p_generic = m.charpoly(algorithm='generic')
....:     assert p_linbox == p_generic
prod_of_row_sums(cols)
randomize(density=1, num_bound=2, den_bound=2, distribution=None, nonzero=False)

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

If x and y are given, randomized entries of this matrix have numerators and denominators bounded by x and y and have density 1.

INPUT:

  • density - number between 0 and 1 (default: 1)

  • num_bound - numerator bound (default: 2)

  • den_bound - denominator bound (default: 2)

  • distribution - None or ‘1/n’ (default: None); if ‘1/n’ then num_bound, den_bound are ignored and numbers are chosen using the GMP function mpq_randomize_entry_recip_uniform

OUTPUT:

  • None, the matrix is modified in-space

EXAMPLES:

The default distribution:

sage: from collections import defaultdict
sage: total_count = 0
sage: dic = defaultdict(Integer)
sage: def add_samples(distribution=None):
....:     global dic, total_count
....:     for _ in range(100):
....:         A = Matrix(QQ, 2, 4, 0)
....:         A.randomize(distribution=distribution)
....:         for a in A.list():
....:             dic[a] += 1
....:             total_count += 1.0

sage: expected = {-2: 1/9, -1: 3/18, -1/2: 1/18, 0: 3/9,
....:             1/2: 1/18, 1: 3/18, 2: 1/9}
sage: add_samples()
sage: while not all(abs(dic[a]/total_count - expected[a]) < 0.001 for a in dic):
....:     add_samples()

The distribution '1/n':

sage: def mpq_randomize_entry_recip_uniform():
....:     r = 2*random() - 1
....:     if r == 0: r = 1
....:     num = int(4/(5*r))
....:     r = random()
....:     if r == 0: r = 1
....:     den = int(1/random())
....:     return Integer(num)/Integer(den)

sage: total_count = 0
sage: dic = defaultdict(Integer)
sage: dic2 = defaultdict(Integer)
sage: add_samples('1/n')
sage: for _ in range(8):
....:     dic2[mpq_randomize_entry_recip_uniform()] += 1
sage: while not all(abs(dic[a] - dic2[a])/total_count < 0.005 for a in dic):
....:     add_samples('1/n')
....:     for _ in range(800):
....:         dic2[mpq_randomize_entry_recip_uniform()] += 1

The default can be used to obtain matrices of different rank:

sage: ranks = [False]*11
sage: while not all(ranks):
....:     for dens in (0.05, 0.1, 0.2, 0.5):
....:         A = Matrix(QQ, 10, 10, 0)
....:         A.randomize(dens)
....:         ranks[A.rank()] = True

The default density is \(6/9\):

sage: def add_sample(density, num_rows, num_cols):
....:     global density_sum, total_count
....:     total_count += 1.0
....:     A = Matrix(QQ, num_rows, num_cols, 0)
....:     A.randomize(density)
....:     density_sum += float(A.density())

sage: density_sum = 0.0
sage: total_count = 0.0
sage: expected_density = 6/9
sage: add_sample(1.0, 100, 100)
sage: while abs(density_sum/total_count - expected_density) > 0.001:
....:     add_sample(1.0, 100, 100)

The modified density depends on the number of columns:

sage: density_sum = 0.0
sage: total_count = 0.0
sage: expected_density = 6/9*0.5
sage: add_sample(0.5, 100, 2)
sage: while abs(density_sum/total_count - expected_density) > 0.001:
....:     add_sample(0.5, 100, 2)

sage: density_sum = 0.0
sage: total_count = 0.0
sage: expected_density = 6/9*(1.0 - (99/100)^50)
sage: expected_density
0.263...

sage: add_sample(0.5, 100, 100)
sage: while abs(density_sum/total_count - expected_density) > 0.001:
....:     add_sample(0.5, 100, 100)

Modifying the bounds for numerator and denominator:

sage: num_dic = defaultdict(Integer)
sage: den_dic = defaultdict(Integer)
sage: while not (all(num_dic[i] for i in range(-200, 201))
....:            and all(den_dic[i] for i in range(1, 101))):
....:     a = matrix(QQ, 2, 4)
....:     a.randomize(num_bound=200, den_bound=100)
....:     for q in a.list():
....:         num_dic[q.numerator()] += 1
....:         den_dic[q.denominator()] += 1
sage: len(num_dic)
401
sage: len(den_dic)
100
rank(algorithm=None)

Return the rank of this matrix.

INPUT:

  • algorithm - an optional specification of an algorithm. One of

    • None: (default) will use flint

    • 'flint': uses the flint library

    • 'pari': uses the PARI library

    • 'integer': eliminate denominators and calls the rank function on the corresponding integer matrix

EXAMPLES:

sage: matrix(QQ,3,[1..9]).rank()
2
sage: matrix(QQ,100,[1..100^2]).rank()
2
row(i, from_list=False)

Return the i-th row of this matrix as a dense vector.

INPUT:

  • i - integer

  • from_list - ignored

EXAMPLES:

sage: m = matrix(QQ, 2, [1/5, -2/3, 3/4, 4/9])
sage: m.row(0)
(1/5, -2/3)
sage: m.row(1)
(3/4, 4/9)
sage: m.row(1, from_list=True)
(3/4, 4/9)
sage: m.row(-2)
(1/5, -2/3)

sage: m.row(2)
Traceback (most recent call last):
...
IndexError: row index out of range
sage: m.row(-3)
Traceback (most recent call last):
...
IndexError: row index out of range
set_row_to_multiple_of_row(i, j, s)

Set row i equal to s times row j.

EXAMPLES:

sage: a = matrix(QQ,2,3,range(6)); 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]
transpose()

Returns 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: A = matrix(QQ, 2, 3, range(6))
sage: type(A)
<type 'sage.matrix.matrix_rational_dense.Matrix_rational_dense'>
sage: B = A.transpose()
sage: print(B)
[0 3]
[1 4]
[2 5]
sage: print(A)
[0 1 2]
[3 4 5]

.T is a convenient shortcut for the transpose:

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