Module: sage.matrix.matrix_real_double_dense
Dense matrices over the real double field.
Matrix operations use GSl and numpy.
sage: b=Mat(RDF,2,3).basis() sage: b[0] [1.0 0.0 0.0] [0.0 0.0 0.0]
We deal with the case of zero rows or zero columns:
sage: m = MatrixSpace(RDF,0,3) sage: m.zero_matrix() []
TESTS:
sage: a = matrix(RDF,2,range(4), sparse=False) sage: loads(dumps(a)) == a True
Author Log:
Class: Matrix_real_double_dense
sage: m = Matrix(RDF, [[1,2],[3,4]]) sage: m**2 [ 7.0 10.0] [15.0 22.0] sage: n= m^(-1); n [-2.0 1.0] [ 1.5 -0.5]
To compute eigenvalues the use the functions left eigen_vectors or right_eigenvectors
sage: p,e = m.right_eigenvectors()
the result of eigen is a pair p,e, where 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,2] and b = [5,6] [3,4]]
sage: b = vector(RDF,[5,6]) sage: m.solve_left(b) (-4.0, 4.5)
See the commands QR,LU,SVD for QR, LU, and singular value decomposition.
Functions: cholesky,
determinant,
eigenspaces,
is_symmetric,
left_eigenvectors,
log_determinant,
LU,
LU_valid,
numpy,
QR,
right_eigenvectors,
solve_left,
solve_left_LU,
SVD,
transpose
) |
Return the cholesky factorization of this matrix. The input matrix must be symmetric and positive definite or an exception will be raised.
sage: M = MatrixSpace(RDF,5) sage: r = matrix(RDF,[[ 0., 0., 0., 0., 1.],[ 1., 1., 1., 1., 1.],[ 16., 8., 4., 2., 1.],[ 81., 27., 9., 3., 1.],[ 256., 64., 16., 4., 1.]])
sage: m = r*M(1)*r.transpose() sage: c = m.cholesky() sage: c*c.transpose() [ 1.0 1.0 1.0 1.0 1.0] [ 1.0 5.0 31.0 121.0 341.0] [ 1.0 31.0 341.0 1555.0 4681.0] [ 1.0 121.0 1555.0 7381.0 22621.0] [ 1.0 341.0 4681.0 22621.0 69905.0]
) |
Return the determinant of self.
ALGORITHM: Use GSL (LU decompositon)
sage: m = matrix(RDF,2,range(4)); m.det() -2.0 sage: m = matrix(RDF,0,[]); m.det() 1.0 sage: m = matrix(RDF, 2, range(6)); m.det() Traceback (most recent call last): ... ValueError: self must be square
) |
Return a list of pairs (e, V) where e runs through all complex eigenvalues of this matrix, and V is the corresponding left eigenspace (always a 1-dimensional complex vector space).
sage: m = matrix(RDF, 3, range(9)); m [0.0 1.0 2.0] [3.0 4.0 5.0] [6.0 7.0 8.0] sage: es = m.eigenspaces() sage: es # random [(13.3484692283, Vector space of degree 3 and dimension 1 over Real Double Field User basis matrix: [-0.440242867236 -0.567868371314 -0.695493875393]), (-1.34846922835, Vector space of degree 3 and dimension 1 over Real Double Field User basis matrix: [-0.897878732262 -0.278434036822 0.341010658618]), (-9.10854412047e-16, Vector space of degree 3 and dimension 1 over Real Double Field User basis matrix: [ 0.408248290464 -0.816496580928 0.408248290464])]
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 (-2.6645352591e-15, -7.1054273576e-15, -3.5527136788e-15)
) |
Return whether this matrix is symmetric, to the given tolerance.
sage: m = matrix(RDF,2,2,range(4)); m [0.0 1.0] [2.0 3.0] sage: m.is_symmetric() False sage: m[1,0] = 1.1; m [0.0 1.0] [1.1 3.0] sage: m.is_symmetric() False
The tolerance inequality is strict:
sage: m.is_symmetric(tol=0.1) False sage: m.is_symmetric(tol=0.11) True
) |
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:
sage: m = Matrix(RDF, 3, range(9)) sage: es = m.left_eigenvectors() sage: es # random-ish platform-dependent output (low order digits) ([13.3484692283, -1.34846922835, -9.10854412047e-16], [-0.440242867236 -0.567868371314 -0.695493875393] [-0.897878732262 -0.278434036822 0.341010658618] [ 0.408248290464 -0.816496580928 0.408248290464]) sage: e, v = es; e = e[0]; v = v[0] sage: abs(abs(e*v - v*m)) < 1e-10 True
IMPLEMENTATION: Uses numpy.
) |
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(RDF,2,2,range(4)); m [0.0 1.0] [2.0 3.0] sage: RDF(log(abs(m.determinant()))) 0.69314718056 sage: m.log_determinant() 0.69314718056 sage: m = matrix(RDF,0,0,[]); m [] sage: m.log_determinant() 0.0
) |
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:
sage: m = matrix(RDF,4,range(16)) sage: P,L,U = m.LU() sage: P*m [12.0 13.0 14.0 15.0] [ 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.0 1.0 2.0 3.0] [ 8.0 9.0 10.0 11.0] [ 4.0 5.0 6.0 7.0]
) |
Returns True
if the LU form of this matrix has
already been computed.
sage: A= random_matrix(RDF,3) ; A.LU_valid() False sage: L,U,P = A.LU() sage: A.LU_valid() True
) |
This method returns a copy of the matrix as a numpy array. It uses the numpy C/api so is very fast.
sage: m = matrix(RDF,[[1,2],[3,4]]) sage: n = m.numpy() sage: import numpy sage: numpy.linalg.eig(n) (array([-0.37228132, 5.37228132]), array([[-0.82456484, -0.41597356], [ 0.56576746, -0.90937671]])) sage: m = matrix(RDF, 2, range(6)); m [0.0 1.0 2.0] [3.0 4.0 5.0] sage: m.numpy() array([[ 0., 1., 2.], [ 3., 4., 5.]])
TESTS:
sage: m = matrix(RDF,0,5,[]); m [] sage: m.numpy() array([], shape=(0, 5), dtype=float64) sage: m = matrix(RDF,5,0,[]); m [] sage: m.numpy() array([], shape=(5, 0), dtype=float64)
) |
Return the Q,R factorization of a real matrix A.
Input:
sage: m = matrix(RDF,3,range(12)); m [ 0.0 1.0 2.0 3.0] [ 4.0 5.0 6.0 7.0] [ 8.0 9.0 10.0 11.0] sage: Q,R = m.QR() sage: Q*R [ 0.0 1.0 2.0 3.0] [ 4.0 5.0 6.0 7.0] [ 8.0 9.0 10.0 11.0]
Note that the columns of Q will be an orthogonal
sage: Q*Q.transpose() # slightly random output. [1.0 5.55111512313e-17 -1.11022302463e-16] [ 5.55111512313e-17 1.0 -5.55111512313e-17] [-1.11022302463e-16 -5.55111512313e-17 1.0]
) |
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:
sage: m = Matrix(RDF, 3, range(9)) sage: m.right_eigenvectors() # random-ish platform-dependent output (low order digits) ([13.3484692283, -1.34846922835, -6.43352295265e-16], [ 0.164763817282 0.799699663112 0.408248290464] [ 0.505774475901 0.104205787719 -0.816496580928] [ 0.846785134519 -0.591288087674 0.408248290464])
IMPLEMENTATION: Uses numpy.
) |
Solve the equation A*x = b, where
sage: A = matrix(RDF, 3,3, [1,2,5,7.6,2.3,1,1,2,-1]); A [ 1.0 2.0 5.0] [ 7.6 2.3 1.0] [ 1.0 2.0 -1.0] sage: b = vector(RDF,[1,2,3]) sage: x = A.solve_left(b); x (-0.113695090439, 1.39018087855, -0.333333333333) sage: A*x (1.0, 2.0, 3.0)
TESTS: We test two degenerate cases:
sage: A = matrix(RDF, 0, 3, []) sage: A.solve_left(vector(RDF,[])) () sage: A = matrix(RDF, 3, 0, []) sage: A.solve_left(vector(RDF,3, [1,2,3])) (0.0, 0.0, 0.0)
) |
Solve the equation A*x = b.
Input:
NOTES: 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 reused.
sage: A = matrix(RDF, 3,3, [1,2,5,7.6,2.3,1,1,2,-1]); A [ 1.0 2.0 5.0] [ 7.6 2.3 1.0] [ 1.0 2.0 -1.0] sage: b = vector(RDF,[1,2,3]) sage: x = A.solve_left(b); x (-0.113695090439, 1.39018087855, -0.333333333333) sage: A*x (1.0, 2.0, 3.0)
TESTS: We test two degenerate cases:
sage: A = matrix(RDF, 0, 3, []) sage: A.solve_left_LU(vector(RDF,[])) (0.0, 0.0, 0.0) sage: A = matrix(RDF, 3, 0, []) sage: A.solve_left_LU(vector(RDF,3, [1,2,3])) ()
) |
Return the singular value decomposition of this matrix.
Input:
sage: m = matrix(RDF,4,range(16)) sage: U,S,V = m.SVD(algorithm='gsl') 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(RDF, 2, range(6)); m [0.0 1.0 2.0] [3.0 4.0 5.0] sage: U, S, V = m.SVD(algorithm='numpy') sage: U*S*V.transpose() # random low bits [7.62194127257e-17 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.0 1.0] sage: V [-0.392540507864 0.824163383692] [-0.560772154092 0.137360563949] [ -0.72900380032 -0.549442255795] sage: U*S*V.transpose() # random low bits [7.62194127257e-17 1.0 2.0] [ 3.0 4.0 5.0] sage: m = matrix(RDF,3,2,range(6)); m [0.0 1.0] [2.0 3.0] [4.0 5.0] sage: U,S,V = m.SVD(algorithm='numpy') sage: U*S*V.transpose() # random low order bits [-8.13151629364e-19 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(RDF, 3, 0, []); m [] sage: m.SVD() ([], [], []) sage: m = matrix(RDF, 0, 3, []); m [] sage: m.SVD() ([], [], [])
) |
Return the transpose of this matrix.
sage: m = matrix(RDF,2,3,range(6)); m [0.0 1.0 2.0] [3.0 4.0 5.0] sage: m.transpose() [0.0 3.0] [1.0 4.0] [2.0 5.0] sage: m = matrix(RDF,0,3); m [] sage: m.transpose() [] sage: m.transpose().parent() Full MatrixSpace of 3 by 0 dense matrices over Real Double Field
Special Functions: __eq__,
__ge__,
__gt__,
__init__,
__invert__,
__le__,
__lt__,
__ne__,
__neg__,
_cholesky_gsl,
_hadamard_row_bound,
_replace_self_with_numpy,
_replace_self_with_numpy32,
_SVD_gsl,
_SVD_numpy
) |
Invert this matrix.
sage: A = Matrix(RDF, [[10, 0], [0, 100]]) sage: (~A).det() 0.001
sage: A = matrix(RDF,3,[2,3,5,7,8,9,11,13,17]);A [ 2.0 3.0 5.0] [ 7.0 8.0 9.0] [11.0 13.0 17.0] sage: ~A [ -2.71428571429 -2.0 1.85714285714] [ 2.85714285714 3.0 -2.42857142857] [-0.428571428571 -1.0 0.714285714286]
Note that if this matrix is (nearly) singular, finding its inverse will not help much and will give slightly different answers on similar platforms depending on the hardware and tuning options given to ATLAS:
sage: A = Matrix(RDF, [[1, 0], [0, 0]]) sage: A.inverse().det() nan sage: A = matrix(RDF,3,range(1,10));A [1.0 2.0 3.0] [4.0 5.0 6.0] [7.0 8.0 9.0]
sage: A.determinant() < 10e-12 True sage: ~A # slightly random [-4.50359962737e+15 9.00719925474e+15 -4.50359962737e+15] [ 9.00719925474e+15 -1.80143985095e+16 9.00719925474e+15] [-4.50359962737e+15 9.00719925474e+15 -4.50359962737e+15]
) |
Negate this matrix
sage: A = matrix(RDF,3,range(1,10)) sage: -A [-1.0 -2.0 -3.0] [-4.0 -5.0 -6.0] [-7.0 -8.0 -9.0] sage: B = -A ; (A+B).is_zero() True
) |
) |
Return an integer n such that the absolute value of the
determinant of this matrix is at most
.
sage: a = matrix(RDF, 3, [1,2,5,7,-3,4,2,1,123]) sage: a._hadamard_row_bound() 4 sage: a.det() -2014.0 sage: 10^4 10000
) |
sage: import numpy sage: a = numpy.array([[1,2],[3,4]], 'float64') sage: m = matrix(RDF,2,2,0) sage: m._replace_self_with_numpy(a) sage: m [1.0 2.0] [3.0 4.0]
) |
sage: import numpy sage: a = numpy.array([[1,2],[3,4]], 'float32') sage: m = matrix(RDF,2,2,0) sage: m._replace_self_with_numpy32(a) sage: m [1.0 2.0] [3.0 4.0]
) |
Return the singular value decomposition of this matrix using GSL. Note that the matrices the dimensions of the matrices that this returns are (m,p), (p,p), and (n, p) where p = min(m,n).
sage: def shape(x): return (x.nrows(), x.ncols()) sage: m = matrix(RDF, 2, 3, range(6)) sage: map(shape, m._SVD_gsl()) [(2, 2), (2, 2), (3, 2)]
) |
Return the singular value decomposition of this matrix using GSL. Note that the matrices the dimensions of the matrices that this returns are (m,m), (m,n), and (n, n).
sage: def shape(x): return (x.nrows(), x.ncols()) sage: m = matrix(RDF, 2, 3, range(6)) sage: map(shape, m._SVD_numpy()) [(2, 2), (2, 3), (3, 3)]
See About this document... for information on suggesting changes.