Arbitrary precision complex ball matrices#

AUTHORS:

  • Clemens Heuberger (2014-10-25): Initial version.

This is an incomplete interface to the acb_mat module of FLINT; it may be useful to refer to its documentation for more details.

class sage.matrix.matrix_complex_ball_dense.Matrix_complex_ball_dense[source]#

Bases: Matrix_dense

Matrix over a complex ball field. Implemented using the acb_mat type of the FLINT library.

EXAMPLES:

sage: MatrixSpace(CBF, 3)(2)
[2.000000000000000                 0                 0]
[                0 2.000000000000000                 0]
[                0                 0 2.000000000000000]
sage: matrix(CBF, 1, 3, [1, 2, -3])
[ 1.000000000000000  2.000000000000000 -3.000000000000000]
>>> from sage.all import *
>>> MatrixSpace(CBF, Integer(3))(Integer(2))
[2.000000000000000                 0                 0]
[                0 2.000000000000000                 0]
[                0                 0 2.000000000000000]
>>> matrix(CBF, Integer(1), Integer(3), [Integer(1), Integer(2), -Integer(3)])
[ 1.000000000000000  2.000000000000000 -3.000000000000000]
charpoly(var='x', algorithm=None)[source]#

Compute the characteristic polynomial of this matrix.

EXAMPLES:

sage: from sage.matrix.benchmark import hilbert_matrix
sage: mat = hilbert_matrix(5).change_ring(ComplexBallField(10))
sage: mat.charpoly()
x^5 + ([-1.8 +/- 0.0258])*x^4 + ([0.3 +/- 0.05...)*x^3 +
([+/- 0.0...])*x^2 + ([+/- 0.0...])*x + [+/- 0.0...]
>>> from sage.all import *
>>> from sage.matrix.benchmark import hilbert_matrix
>>> mat = hilbert_matrix(Integer(5)).change_ring(ComplexBallField(Integer(10)))
>>> mat.charpoly()
x^5 + ([-1.8 +/- 0.0258])*x^4 + ([0.3 +/- 0.05...)*x^3 +
([+/- 0.0...])*x^2 + ([+/- 0.0...])*x + [+/- 0.0...]
contains(other)[source]#

Test if the set of complex matrices represented by self is contained in that represented by other.

EXAMPLES:

sage: b = CBF(0, RBF(0, rad=.1r)); b
[+/- 0.101]*I
sage: matrix(CBF, [0, b]).contains(matrix(CBF, [0, 0]))
True
sage: matrix(CBF, [0, b]).contains(matrix(CBF, [b, 0]))
False
sage: matrix(CBF, [b, b]).contains(matrix(CBF, [b, 0]))
True
>>> from sage.all import *
>>> b = CBF(Integer(0), RBF(Integer(0), rad=.1)); b
[+/- 0.101]*I
>>> matrix(CBF, [Integer(0), b]).contains(matrix(CBF, [Integer(0), Integer(0)]))
True
>>> matrix(CBF, [Integer(0), b]).contains(matrix(CBF, [b, Integer(0)]))
False
>>> matrix(CBF, [b, b]).contains(matrix(CBF, [b, Integer(0)]))
True
determinant()[source]#

Compute the determinant of this matrix.

EXAMPLES:

sage: matrix(CBF, [[1/2, 1/3], [1, 1]]).determinant()
[0.1666666666666667 +/- ...e-17]
sage: matrix(CBF, [[1/2, 1/3], [1, 1]]).det()
[0.1666666666666667 +/- ...e-17]
sage: matrix(CBF, [[1/2, 1/3]]).determinant()
Traceback (most recent call last):
...
ValueError: self must be a square matrix
>>> from sage.all import *
>>> matrix(CBF, [[Integer(1)/Integer(2), Integer(1)/Integer(3)], [Integer(1), Integer(1)]]).determinant()
[0.1666666666666667 +/- ...e-17]
>>> matrix(CBF, [[Integer(1)/Integer(2), Integer(1)/Integer(3)], [Integer(1), Integer(1)]]).det()
[0.1666666666666667 +/- ...e-17]
>>> matrix(CBF, [[Integer(1)/Integer(2), Integer(1)/Integer(3)]]).determinant()
Traceback (most recent call last):
...
ValueError: self must be a square matrix
eigenvalues(other=None, extend=None)[source]#

(Experimental.) Compute rigorous enclosures of the eigenvalues of this matrix.

INPUT:

  • self – an \(n \times n\) matrix

  • other – unsupported (generalized eigenvalue problem), should be None

  • extend – ignored

OUTPUT:

A Sequence of complex balls of length equal to the size of the matrix.

Each element represents one eigenvalue with the correct multiplicities in case of overlap. The output intervals are either disjoint or identical, and identical intervals are guaranteed to be grouped consecutively. Each complete run of \(k\) identical balls thus represents a cluster of exactly \(k\) eigenvalues which could not be separated from each other at the current precision, but which could be isolated from the other eigenvalues.

There is currently no guarantee that the algorithm converges as the working precision is increased.

See the FLINT documentation for more information.

EXAMPLES:

sage: from sage.matrix.benchmark import hilbert_matrix
sage: mat = hilbert_matrix(5).change_ring(CBF)
sage: mat.eigenvalues()
doctest:...: FutureWarning: This class/method/function is marked as experimental.
...
[[1.567050691098...] + [+/- ...]*I, [0.208534218611...] + [+/- ...]*I,
[3.287928...e-6...] + [+/- ...]*I, [0.000305898040...] + [+/- ...]*I,
[0.011407491623...] + [+/- ...]*I]

sage: mat = Permutation([2, 1, 4, 5, 3]).to_matrix().dense_matrix().change_ring(CBF)
sage: mat.eigenvalues()
Traceback (most recent call last):
...
ValueError: unable to certify the eigenvalues
sage: precond = matrix(ZZ, [[-1, -2, 2, 2, -2], [2, -2, -2, -2, 2],
....:     [-2, 2, -1, 2, 1], [2, 1, -1, 0, 2], [-2, 0, 1, -1, 1]])
sage: (~precond*mat*precond).eigenvalues()
[[-0.5000000000000...] + [-0.8660254037844...]*I, [-1.000000000000...] + [+/- ...]*I,
 [-0.5000000000000...] + [0.8660254037844...]*I,
 [1.000000000000...] + [+/- ...]*I, [1.000000000000...] + [+/- ...]*I]
>>> from sage.all import *
>>> from sage.matrix.benchmark import hilbert_matrix
>>> mat = hilbert_matrix(Integer(5)).change_ring(CBF)
>>> mat.eigenvalues()
doctest:...: FutureWarning: This class/method/function is marked as experimental.
...
[[1.567050691098...] + [+/- ...]*I, [0.208534218611...] + [+/- ...]*I,
[3.287928...e-6...] + [+/- ...]*I, [0.000305898040...] + [+/- ...]*I,
[0.011407491623...] + [+/- ...]*I]

>>> mat = Permutation([Integer(2), Integer(1), Integer(4), Integer(5), Integer(3)]).to_matrix().dense_matrix().change_ring(CBF)
>>> mat.eigenvalues()
Traceback (most recent call last):
...
ValueError: unable to certify the eigenvalues
>>> precond = matrix(ZZ, [[-Integer(1), -Integer(2), Integer(2), Integer(2), -Integer(2)], [Integer(2), -Integer(2), -Integer(2), -Integer(2), Integer(2)],
...     [-Integer(2), Integer(2), -Integer(1), Integer(2), Integer(1)], [Integer(2), Integer(1), -Integer(1), Integer(0), Integer(2)], [-Integer(2), Integer(0), Integer(1), -Integer(1), Integer(1)]])
>>> (~precond*mat*precond).eigenvalues()
[[-0.5000000000000...] + [-0.8660254037844...]*I, [-1.000000000000...] + [+/- ...]*I,
 [-0.5000000000000...] + [0.8660254037844...]*I,
 [1.000000000000...] + [+/- ...]*I, [1.000000000000...] + [+/- ...]*I]
eigenvectors_left(other=None, extend=True)[source]#

(Experimental.) Compute rigorous enclosures of the eigenvalues and left eigenvectors of this matrix.

INPUT:

  • self – an \(n \times n\) matrix

  • other – unsupported (generalized eigenvalue problem), should be None

  • extend – ignored

OUTPUT:

A list of triples of the form (eigenvalue, [eigenvector], 1).

Unlike eigenvalues() and eigenvectors_left_approx(), this method currently fails in the presence of multiple eigenvalues.

Additionally, there is currently no guarantee that the algorithm converges as the working precision is increased.

See the FLINT documentation for more information.

EXAMPLES:

sage: mat = matrix(CBF, 3, [2, 3, 5, 7, 11, 13, 17, 19, 23])
sage: eigval, eigvec, _ = mat.eigenvectors_left()[0]
sage: eigval
[1.1052996349...] + [+/- ...]*I
sage: eigvec[0]
([0.69817246751...] + [+/- ...]*I, [-0.67419514369...] + [+/- ...]*I, [0.240865343781...] + [+/- ...]*I)
sage: eigvec[0] * (mat - eigval)
([+/- ...] + [+/- ...]*I, [+/- ...] + [+/- ...]*I, [+/- ...] + [+/- ...]*I)
>>> from sage.all import *
>>> mat = matrix(CBF, Integer(3), [Integer(2), Integer(3), Integer(5), Integer(7), Integer(11), Integer(13), Integer(17), Integer(19), Integer(23)])
>>> eigval, eigvec, _ = mat.eigenvectors_left()[Integer(0)]
>>> eigval
[1.1052996349...] + [+/- ...]*I
>>> eigvec[Integer(0)]
([0.69817246751...] + [+/- ...]*I, [-0.67419514369...] + [+/- ...]*I, [0.240865343781...] + [+/- ...]*I)
>>> eigvec[Integer(0)] * (mat - eigval)
([+/- ...] + [+/- ...]*I, [+/- ...] + [+/- ...]*I, [+/- ...] + [+/- ...]*I)
eigenvectors_left_approx(other=None, extend=None)[source]#

(Experimental.) Compute non-rigorous approximations of the left eigenvalues and eigenvectors of this matrix.

INPUT:

  • self – an \(n \times n\) matrix

  • other – unsupported (generalized eigenvalue problem), should be None

  • extend – ignored

OUTPUT:

A list of triples of the form (eigenvalue, [eigenvector], 1). The eigenvalue and the entries of the eigenvector are complex balls with zero radius.

No guarantees are made about the accuracy of the output.

See the FLINT documentation for more information.

EXAMPLES:

sage: mat = matrix(CBF, 3, [2, 3, 5, 7, 11, 13, 17, 19, 23])
sage: eigval, eigvec, _ = mat.eigenvectors_left_approx()[0]
sage: eigval
[1.1052996349... +/- ...]
sage: eigvec[0]
([0.69817246751...], [-0.67419514369...], [0.240865343781...])
sage: eigvec[0] * (mat - eigval)
([+/- ...], [+/- ...], [+/- ...])
>>> from sage.all import *
>>> mat = matrix(CBF, Integer(3), [Integer(2), Integer(3), Integer(5), Integer(7), Integer(11), Integer(13), Integer(17), Integer(19), Integer(23)])
>>> eigval, eigvec, _ = mat.eigenvectors_left_approx()[Integer(0)]
>>> eigval
[1.1052996349... +/- ...]
>>> eigvec[Integer(0)]
([0.69817246751...], [-0.67419514369...], [0.240865343781...])
>>> eigvec[Integer(0)] * (mat - eigval)
([+/- ...], [+/- ...], [+/- ...])
eigenvectors_right(other=None, extend=None)[source]#

(Experimental.) Compute rigorous enclosures of the eigenvalues and eigenvectors of this matrix.

INPUT:

  • self – an \(n \times n\) matrix

  • other – unsupported (generalized eigenvalue problem), should be None

  • extend – ignored

OUTPUT:

A list of triples of the form (eigenvalue, [eigenvector], 1).

Unlike eigenvalues() and eigenvectors_right_approx(), this method currently fails in the presence of multiple eigenvalues.

Additionally, there is currently no guarantee that the algorithm converges as the working precision is increased.

See the FLINT documentation for more information.

EXAMPLES:

sage: from sage.matrix.benchmark import hilbert_matrix
sage: mat = hilbert_matrix(3).change_ring(CBF)
sage: eigval, eigvec, _ = mat.eigenvectors_right()[0]
doctest:...: FutureWarning: This class/method/function is marked as experimental.
...
sage: eigval
[1.40831892712...] + [+/- ...]*I
sage: eigvec
[([0.82704492697...] + [+/- ...]*I, [0.45986390436...] + [+/- ...]*I, [0.32329843524...] + [+/- ...]*I)]
sage: (mat - eigval)*eigvec[0]
([+/- ...] + [+/- ...]*I, [+/- ...] + [+/- ...]*I, [+/- ...] + [+/- ...]*I)
>>> from sage.all import *
>>> from sage.matrix.benchmark import hilbert_matrix
>>> mat = hilbert_matrix(Integer(3)).change_ring(CBF)
>>> eigval, eigvec, _ = mat.eigenvectors_right()[Integer(0)]
doctest:...: FutureWarning: This class/method/function is marked as experimental.
...
>>> eigval
[1.40831892712...] + [+/- ...]*I
>>> eigvec
[([0.82704492697...] + [+/- ...]*I, [0.45986390436...] + [+/- ...]*I, [0.32329843524...] + [+/- ...]*I)]
>>> (mat - eigval)*eigvec[Integer(0)]
([+/- ...] + [+/- ...]*I, [+/- ...] + [+/- ...]*I, [+/- ...] + [+/- ...]*I)
eigenvectors_right_approx(other=None, extend=None)[source]#

(Experimental.) Compute non-rigorous approximations of the eigenvalues and eigenvectors of this matrix.

INPUT:

  • self – an \(n \times n\) matrix

  • other – unsupported (generalized eigenvalue problem), should be None

  • extend – ignored

OUTPUT:

A list of triples of the form (eigenvalue, [eigenvector], 1). The eigenvalue and the entries of the eigenvector are complex balls with zero radius.

No guarantees are made about the accuracy of the output.

See the FLINT documentation for more information.

EXAMPLES:

sage: from sage.matrix.benchmark import hilbert_matrix
sage: mat = hilbert_matrix(3).change_ring(CBF)
sage: eigval, eigvec, _ = mat.eigenvectors_right_approx()[0]
doctest:...: FutureWarning: This class/method/function is marked as experimental.
...
sage: eigval
[1.40831892712...]
sage: eigval.rad()
0.00000000
sage: eigvec
[([0.8270449269720...], [0.4598639043655...], [0.3232984352444...])]
sage: (mat - eigval)*eigvec[0]
([1e-15 +/- ...], [2e-15 +/- ...], [+/- ...])
>>> from sage.all import *
>>> from sage.matrix.benchmark import hilbert_matrix
>>> mat = hilbert_matrix(Integer(3)).change_ring(CBF)
>>> eigval, eigvec, _ = mat.eigenvectors_right_approx()[Integer(0)]
doctest:...: FutureWarning: This class/method/function is marked as experimental.
...
>>> eigval
[1.40831892712...]
>>> eigval.rad()
0.00000000
>>> eigvec
[([0.8270449269720...], [0.4598639043655...], [0.3232984352444...])]
>>> (mat - eigval)*eigvec[Integer(0)]
([1e-15 +/- ...], [2e-15 +/- ...], [+/- ...])
exp()[source]#

Compute the exponential of this matrix.

EXAMPLES:

sage: matrix(CBF, [[i*pi, 1], [0, i*pi]]).exp()                             # needs sage.symbolic
[[-1.00000000000000 +/- ...e-16] + [+/- ...e-16]*I [-1.00000000000000 +/- ...e-16] + [+/- ...e-16]*I]
[                                                0 [-1.00000000000000 +/- ...e-16] + [+/- ...e-16]*I]
sage: matrix(CBF, [[1/2, 1/3]]).exp()
Traceback (most recent call last):
...
ValueError: self must be a square matrix
>>> from sage.all import *
>>> matrix(CBF, [[i*pi, Integer(1)], [Integer(0), i*pi]]).exp()                             # needs sage.symbolic
[[-1.00000000000000 +/- ...e-16] + [+/- ...e-16]*I [-1.00000000000000 +/- ...e-16] + [+/- ...e-16]*I]
[                                                0 [-1.00000000000000 +/- ...e-16] + [+/- ...e-16]*I]
>>> matrix(CBF, [[Integer(1)/Integer(2), Integer(1)/Integer(3)]]).exp()
Traceback (most recent call last):
...
ValueError: self must be a square matrix
identical(other)[source]#

Test if the corresponding entries of two complex ball matrices represent the same balls.

EXAMPLES:

sage: a = matrix(CBF, [[1/3,2],[3,4]])
sage: b = matrix(CBF, [[1/3,2],[3,4]])
sage: a == b
False
sage: a.identical(b)
True
>>> from sage.all import *
>>> a = matrix(CBF, [[Integer(1)/Integer(3),Integer(2)],[Integer(3),Integer(4)]])
>>> b = matrix(CBF, [[Integer(1)/Integer(3),Integer(2)],[Integer(3),Integer(4)]])
>>> a == b
False
>>> a.identical(b)
True
overlaps(other)[source]#

Test if two matrices with complex ball entries represent overlapping sets of complex matrices.

EXAMPLES:

sage: b = CBF(0, RBF(0, rad=0.1r)); b
[+/- 0.101]*I
sage: matrix(CBF, [0, b]).overlaps(matrix(CBF, [b, 0]))
True
sage: matrix(CBF, [1, 0]).overlaps(matrix(CBF, [b, 0]))
False
>>> from sage.all import *
>>> b = CBF(Integer(0), RBF(Integer(0), rad=0.1)); b
[+/- 0.101]*I
>>> matrix(CBF, [Integer(0), b]).overlaps(matrix(CBF, [b, Integer(0)]))
True
>>> matrix(CBF, [Integer(1), Integer(0)]).overlaps(matrix(CBF, [b, Integer(0)]))
False
trace()[source]#

Compute the trace of this matrix.

EXAMPLES:

sage: matrix(CBF, [[1/3, 1/3], [1, 1]]).trace()
[1.333333333333333 +/- ...e-16]
sage: matrix(CBF, [[1/2, 1/3]]).trace()
Traceback (most recent call last):
...
ValueError: self must be a square matrix
>>> from sage.all import *
>>> matrix(CBF, [[Integer(1)/Integer(3), Integer(1)/Integer(3)], [Integer(1), Integer(1)]]).trace()
[1.333333333333333 +/- ...e-16]
>>> matrix(CBF, [[Integer(1)/Integer(2), Integer(1)/Integer(3)]]).trace()
Traceback (most recent call last):
...
ValueError: self must be a square matrix
transpose()[source]#

Return the transpose of self.

EXAMPLES:

sage: m = matrix(CBF, 2, 3, [1, 2, 3, 4, 5, 6])
sage: m.transpose()
[1.000000000000000 4.000000000000000]
[2.000000000000000 5.000000000000000]
[3.000000000000000 6.000000000000000]
sage: m.transpose().parent()
Full MatrixSpace of 3 by 2 dense matrices over Complex ball field with 53 bits of precision
>>> from sage.all import *
>>> m = matrix(CBF, Integer(2), Integer(3), [Integer(1), Integer(2), Integer(3), Integer(4), Integer(5), Integer(6)])
>>> m.transpose()
[1.000000000000000 4.000000000000000]
[2.000000000000000 5.000000000000000]
[3.000000000000000 6.000000000000000]
>>> m.transpose().parent()
Full MatrixSpace of 3 by 2 dense matrices over Complex ball field with 53 bits of precision