32.19 Dense matrices over the Complex Double Field

Module: sage.matrix.matrix_complex_double_dense

Dense matrices over the Complex Double Field.

IMPLEMENTATION: Specialized matrix operations use GSL and numpy.

sage: b = Mat(CDF,2,3).basis()
sage: b[0]
[1.0   0   0]
[  0   0   0]

TESTS:

sage: a = matrix(CDF,2,range(4), sparse=False)
sage: loads(dumps(a)) == a
True

Class: Matrix_complex_double_dense

class Matrix_complex_double_dense
Class that implements matrices over the complex double field. These are supposed to be fast matrix operations using C doubles. Most operations are implemented using GSl or numpy libraries which will call the underlying BLAS on the system.

Examples:

sage: m = Matrix(CDF, [[1,2*I],[3+I,4]])
sage: m**2
[-1.0 + 6.0*I       10.0*I]
[15.0 + 5.0*I 14.0 + 6.0*I]

sage: n= m^(-1); n
[  0.333333333333 + 0.333333333333*I   0.166666666667 - 0.166666666667*I]
[ -0.166666666667 - 0.333333333333*I 0.0833333333333 + 0.0833333333333*I]

To compute eigenvalues the use the functions left_eigenvectors or right_eigenvectors

sage: p,e = m.right_eigenvectors()

the result of eigen is a pair p,e . p is a list of eigenvalues and the e is a matrix whose columns are the eigenvectors.

To solve a linear system Ax = b for A = [[1.0,2*I] and b = [1,I] [3+I,4]]

sage: b = vector(CDF,[1,I])
sage: m.solve_left(b)
(0.5 + 0.5*I, -0.25 - 0.25*I)

Functions: determinant,$ \,$ eigenspaces,$ \,$ left_eigenvectors,$ \,$ log_determinant,$ \,$ LU,$ \,$ numpy,$ \,$ right_eigenvectors,$ \,$ solve_left,$ \,$ solve_left_LU,$ \,$ SVD,$ \,$ transpose

determinant( )

Return the determinant of self.

ALGORITHM: Use GSL (LU decompositon)

sage: m = matrix(CDF,2,range(4)); m.det()
-2.0
sage: m = matrix(CDF,0,[]); m.det()
1.0
sage: m = matrix(CDF, 2, range(6)); m.det()
Traceback (most recent call last):
...
ValueError: self must be square

eigenspaces( )

Return a list of pairs (e, V) where e runs through all complex eigenvalues of this matrix, and V is the corresponding eigenspace (always a 1-dimensional complex vector space).

sage: m = matrix(CDF, 3, 3, [[0, 0, 1], [1, 0, 1], [0, 1, 0]]); m
[  0   0 1.0]
[1.0   0 1.0]
[  0 1.0   0]
sage: es = m.eigenspaces()
sage: es # random
[(1.32471795724 - 5.55111512313e-17*I, Vector space of degree 3 and
dimension 1 over Complex Double Field
User basis matrix:
[0.413998885523 - 5.55111512313e-17*I 0.548431757932 - 2.22044604925e-16*I 
0.726517398056]), (-0.662358978622 + 0.562279512062*I, Vector space of
degree 3 and dimension 1 over Complex Double Field
User basis matrix:
[                    0.655865618097 -0.434418480916 + 0.368779799722*I
0.0803836554683 - 0.488529222961*I]), (-0.662358978622 - 0.562279512062*I,
Vector space of degree 3 and dimension 1 over Complex Double Field
User basis matrix:
[                    0.655865618097 -0.434418480916 - 0.368779799722*I
0.0803836554683 + 0.488529222961*I])]

sage: e, v = es[0]
sage: v = v.basis()[0]
sage: a = v * m
sage: b = e * v
sage: diff = a.change_ring(CDF) - b
sage: abs(abs(diff)) < 1e-10
True
sage: diff # random -- very small numbers
(6.66133814775e-16 - 1.25526431318e-16*I, 3.24590553708e-16*I,
9.99200722163e-16 - 2.37225939001e-16*I)

left_eigenvectors( )

Computes the eigenvalues and *left* eigenvectors of this matrix m acting *from the right*. I.e., vectors v such that v*m = lambda*v.

Output:

eigenvalues
- as a list
corresponding eigenvectors
- as an CDF matrix whose rows are the eigenvectors.

sage: I = CDF.gen()
sage: m = I*Matrix(CDF, 3, range(9))
sage: vals, vecs = m.left_eigenvectors()
sage: vecs*m      # random precision
[      1.83356748856e-16 + 5.8765683663*I      3.11202346481e-16 +
7.58017348024*I      4.39047944106e-16 + 9.28377859418*I]
[     1.66533453694e-15 - 1.21076184124*I     1.72084568817e-15 -
0.375459730779*I      1.7763568394e-15 + 0.459842379686*I]
[-3.33066907388e-16 - 6.66133814775e-16*I -5.55111512313e-16 -
7.21644966006e-16*I -7.77156117238e-16 - 7.77156117238e-16*I]

sage: vals[0] * vecs.row(0)  # random precision
(1.67272748078e-15 + 5.8765683663*I, 1.82458083857e-15 + 7.58017348024*I,
1.23544530231e-15 + 9.28377859418*I)
sage: vals[1] * vecs.row(1)  # random precision
(-2.00045384696e-16 - 1.21076184124*I, 5.36806149949e-16 -
0.375459730779*I, -5.97719297047e-16 + 0.459842379686*I)
sage: vals[2] * vecs.row(2)  # random precision
(-2.0580841304e-16 + 1.41000047373e-16*I, 4.11616826079e-16 -
2.82000094746e-16*I, -2.0580841304e-16 + 1.41000047373e-16*I)

log_determinant( )

Compute the log of the absolute value of the determinant using GSL (LU decomposition).

NOTE: This is useful if the usual determinant overlows.

sage: m = matrix(CDF,2,2,range(4)); m
[  0 1.0]
[2.0 3.0]
sage: RDF(log(abs(m.determinant())))
0.69314718056
sage: m.log_determinant()
0.69314718056
sage: m = matrix(CDF,0,0,[]); m
[]
sage: m.log_determinant()
0.0

LU( )

Computes the LU decomposition of a matrix.

For and square matrix A we can find matrices P,L, and U. s.t. P*A = L*U where P is a permutation matrix, L is lower triangular and U is upper triangular.

Output:

P, L, U
- as a tuple

sage: m=matrix(CDF,4,range(16))
sage: P,L,U = m.LU()
sage: P*m
[12.0 13.0 14.0 15.0]
[   0  1.0  2.0  3.0]
[ 8.0  9.0 10.0 11.0]
[ 4.0  5.0  6.0  7.0]        
sage: L*U
[12.0 13.0 14.0 15.0]
[   0  1.0  2.0  3.0]
[ 8.0  9.0 10.0 11.0]
[ 4.0  5.0  6.0  7.0]

numpy( )

This method returns a copy of the matrix as a numpy array. It is fast as the copy is done using the numpy C/api.

sage: I = CDF.gen()
sage: m = matrix(CDF,[[1,2],[3,4]])
sage: m = I*m
sage: n = m.numpy()
sage: import numpy
sage: numpy.linalg.eig(n)       # random low bits
(array([  0.00000000e+00-0.37228132j,  -8.03393810e-17+5.37228132j]),
array([[ 0.82456484 +0.00000000e+00j,  0.41597356 +8.77187320e-17j],
       [-0.56576746 -0.00000000e+00j,  0.90937671 +0.00000000e+00j]]))

TESTS:

sage: m = matrix(CDF, 2, range(6)); m
[  0 1.0 2.0]
[3.0 4.0 5.0]
sage: m.numpy()
array([[ 0.+0.j,  1.+0.j,  2.+0.j],
       [ 3.+0.j,  4.+0.j,  5.+0.j]])
sage: m = matrix(CDF,0,5,[]); m
[]
sage: m.numpy()
array([], shape=(0, 5), dtype=complex128)
sage: m = matrix(CDF,5,0,[]); m
[]

right_eigenvectors( )

Computes the eigenvalues and *right* eigenvectors of this matrix m acting *from the left*. I.e., vectors v such that m * v = lambda*v, where v is viewed as a column vector.

Output:

eigenvalues
- as a list
corresponding eigenvectors
- as an CDF matrix whose columns are the eigenvectors.

sage: I = CDF.gen()
sage: m = I*Matrix(CDF, 3, range(9))
sage: vals, vecs = m.right_eigenvectors()
sage: m*vecs      # random precision

[ -2.18763583539e-17 + 2.19934474494*I 2.22220787778e-16 - 1.07837038763*I -1.28830323146e-16] [ -5.79030412929e-16 + 6.75131502804*I 2.46412048752e-16 - 0.140518298155*I -3.67002435386e-16 + 1.66533453694e-16*I] [ -1.1361844675e-15 + 11.3032853111*I 2.70603309725e-16 + 0.797333791317*I -6.05174547627e-16 + 3.33066907388e-16*I]

    sage: vals[0] * vecs.column(0)  # random precision
    (-2.45441964831e-15 + 2.19934474494*I, -1.11280381694e-15 +
6.75131502804*I, -1.37419154761e-15 + 11.3032853111*I)
    sage: vals[1] * vecs.column(1)  # random precision
    (3.94205510736e-18 - 1.07837038763*I, 2.78424120914e-16 -
0.140518298155*I, -2.91698877571e-16 + 0.797333791317*I)
    sage: vals[2] * vecs.column(2)  # random precision
    (-5.48020775807e-17 - 3.30603810196e-17*I, 1.09604155161e-16 +
6.61207620392e-17*I, -5.48020775807e-17 - 3.30603810196e-17*I)
IMPLEMENTATION:
    Uses numpy.

solve_left( )

Solve the equation A*x = b, where

       sage: I = CDF.gen()
       sage: A =I*matrix(CDF, 3,3, [1,2,5,7.6,2.3,1,1,2,-1])
       sage: A   # slightly random output
       [1.0*I             2.0*I                5.0*I]
[7.59999990463*I   2.29999995232*I      1.0*I]
 [1.0*I             2.0*I               -1.0*I]
       sage: b = vector(CDF,[1,2,3])+I*vector(CDF,[1,2,3])
       sage: x = A.solve_left(b); x
        (-0.113695090439 + 0.113695090439*I, 1.39018087855 -
1.39018087855*I, -0.333333333333 + 0.333333333333*I)
       sage: A*x
       (1.0 + 1.0*I, 2.0 + 2.0*I, 3.0 + 3.0*I)

TESTS: We test two degenerate cases:

sage: A = matrix(CDF, 0, 3, [])
sage: A.solve_left(vector(CDF,[]))
()
sage: A = matrix(CDF, 3, 0, [])
sage: A.solve_left(vector(CDF,3, [1,2,3]))
(0, 0, 0)

solve_left_LU( )

Solve the equation A*x = b, where

       sage: I = CDF.gen()
       sage: A = I*matrix(CDF, 3,3, [1,2,5,7.6,2.3,1,1,2,-1])
       sage: A   # slightly random output
       [1.0*I             2.0*I                5.0*I]
[7.59999990463*I   2.29999995232*I      1.0*I]
 [1.0*I             2.0*I               -1.0*I]
       sage: b = vector(CDF,[1,2,3])+I*vector(CDF,[1,2,3])
       sage: x = A.solve_left(b); x
       (-0.113695090439 + 0.113695090439*I, 1.39018087855 -
1.39018087855*I, -0.333333333333 + 0.333333333333*I)
       sage: A*x
       (1.0 + 1.0*I, 2.0 + 2.0*I, 3.0 + 3.0*I)

TESTS: We test two degenerate cases:

sage: A = matrix(CDF, 0, 3, [])
sage: A.solve_left_LU(vector(CDF,[]))
()
sage: A = matrix(CDF, 3, 0, [])
sage: A.solve_left_LU(vector(CDF,3, [1,2,3]))
(0, 0, 0)

This method precomputes and stores the LU decomposition before solving. If many equations of the form Ax=b need to be solved for a singe matrix A, then this method should be used instead of solve. The first time this method is called it will compute the LU decomposition. If the matrix has not changed then subsequent calls will be very fast as the precomputed LU decomposition will be used.

SVD( )

Return the singular value decomposition of this matrix.

Input:

A
- a matrix
Output:
U, S, V
- matrices such that $ A = U S V^t$ , where U and V are orthogonal and S is diagonal.

sage: m = matrix(CDF,4,range(16))
sage: U,S,V = m.SVD()
sage: U*S*V.transpose()    # slightly random output (due to computer architecture)
[3.45569519412e-16               1.0               2.0               3.0]
[4.0               5.0               6.0               7.0]
[8.0               9.0              10.0              11.0]
[12.0              13.0              14.0              15.0]

A non-square example:

sage: m = matrix(CDF, 2, range(6)); m
[  0 1.0 2.0]
[3.0 4.0 5.0]
sage: U, S, V = m.SVD()
sage: U
[-0.274721127897 -0.961523947641]
[-0.961523947641  0.274721127897]
sage: S
[7.34846922835             0             0]
[            0           1.0             0]
sage: V
[-0.392540507864  0.824163383692  0.408248290464]
[-0.560772154092  0.137360563949 -0.816496580928]
[ -0.72900380032 -0.549442255795  0.408248290464]
sage: U*S*V.transpose()           # random low bits
[7.62194127257e-17               1.0               2.0]
[              3.0               4.0               5.0]
sage: m = matrix(CDF,3,2,range(6)); m
[  0 1.0]
[2.0 3.0]
[4.0 5.0]
sage: U,S,V = m.SVD()
sage: U*S*V.transpose()   # random low order bits
[-8.13151629364e-19                1.0]
[               2.0                3.0]
[               4.0                5.0]

TESTS:

sage: m = matrix(CDF, 3, 0, []); m
[]
sage: m.SVD()
([], [], [])
sage: m = matrix(CDF, 0, 3, []); m
[]        
sage: m.SVD()
([], [], [])

transpose( )

Return the transpose of this matrix.

sage: m = matrix(CDF,2,3,range(6)); m
[  0 1.0 2.0]
[3.0 4.0 5.0]
sage: m.transpose()
[  0 3.0]
[1.0 4.0]
[2.0 5.0]
sage: m = matrix(CDF,0,3); m
[]
sage: m.transpose().parent()
Full MatrixSpace of 3 by 0 dense matrices over Complex Double Field

Special Functions: __eigen_numpy,$ \,$ __eq__,$ \,$ __ge__,$ \,$ __gt__,$ \,$ __init__,$ \,$ __invert__,$ \,$ __le__,$ \,$ __lt__,$ \,$ __ne__,$ \,$ __neg__,$ \,$ _multiply_classical,$ \,$ _replace_self_with_numpy,$ \,$ _replace_self_with_numpy32

__eigen_numpy( )

Computes the eigenvalues and eigenvectors of this matrix acting from the left on column vectors.

Output:

eigenvalues
- as a list
corresponding eigenvectors
- as a matrix whose ** ROWS ** are the eigenvectors of self acting from the ** LEFT ** (which is really weird!)

IMPLEMENTATION: Uses numpy.

__invert__( )

__neg__( )

_multiply_classical( )

_replace_self_with_numpy( )

Replaces the entries of self with the elements of a numpy array with dtype = 'complex128'.

sage: import numpy
sage: j = numpy.complex(0,1)
sage: a = numpy.array([[1*j, 2*j],[3*j, 4*j]], 'complex128')
sage: m = matrix(CDF, 2, 2, 0)
sage: m._replace_self_with_numpy(a)
sage: m
[1.0*I 2.0*I]
[3.0*I 4.0*I]

_replace_self_with_numpy32( )

Replaces the entries of self with the elements of a numpy array with dtype = 'complex64'.

sage: import numpy
sage: j = numpy.complex(0,1)
sage: a = numpy.array([[1*j, 2*j],[3*j, 4*j]], 'complex64')
sage: m = matrix(CDF, 2, 2, 0)
sage: m._replace_self_with_numpy32(a)
sage: m
[1.0*I 2.0*I]
[3.0*I 4.0*I]

See About this document... for information on suggesting changes.