X-Git-Url: http://gitweb.michael.orlitzky.com/?a=blobdiff_plain;f=mjo%2Fmatrix_vector.py;h=83bf313a82cde7c689d10337efff8054301acad7;hb=ba04841bc8b5786c75b5e6fa9e767895c6af51d0;hp=e1d1931e6c574053827038a9764bebf9ef0c482e;hpb=971bb0a6f516d838ad3ad35745cb49858fda2d12;p=sage.d.git diff --git a/mjo/matrix_vector.py b/mjo/matrix_vector.py index e1d1931..83bf313 100644 --- a/mjo/matrix_vector.py +++ b/mjo/matrix_vector.py @@ -7,44 +7,220 @@ two vector spaces often. """ from sage.all import * +from sage.matrix.matrix_space import is_MatrixSpace -def isomorphism(matrix_space): +def _mat2vec(m): + return vector(m.base_ring(), m.list()) + +def basis_representation(M): """ - Create isomorphism (i.e. the function) that converts elements - of a matrix space into those of the corresponding finite-dimensional - vector space. + Return the forward (``MatrixSpace`` -> ``VectorSpace``) and + inverse isometries, as a pair, that take elements of the given + ``MatrixSpace`` `M` to their representations as "long vectors," + and vice-versa. + + The argument ``M`` can be either a ``MatrixSpace`` or a basis for + a space of matrices. This function is needed because SageMath does + not know that matrix spaces are vector spaces, and therefore + cannot perform common operations with them -- like computing the + basis representation of an element. + + Moreover, the ability to pass in a basis (rather than a + ``MatrixSpace``) is needed because SageMath has no way to express + that e.g. a (sub)space of symmetric matrices is itself a + ``MatrixSpace``. INPUT: - - matrix_space: A finite-dimensional ``MatrixSpace`` object. + - ``M`` -- Either a ``MatrixSpace``, or a list of matrices that form + a basis for a matrix space. OUTPUT: - - (phi, phi_inverse): If ``matrix_space`` has dimension m*n, then - ``phi`` will map m-by-n matrices to R^(m*n). - The inverse mapping ``phi_inverse`` will go - the other way. + A pair of isometries ``(phi, phi_inv)``. + + If the matrix space associated with `M` has dimension `n`, then + ``phi`` will map its elements to vectors of length `n` over the + same base ring. The inverse map ``phi_inv`` reverses that + operation. + + SETUP:: + + sage: from mjo.matrix_vector import basis_representation EXAMPLES: - sage: M = MatrixSpace(QQ,4,4) - sage: (p, p_inv) = isomorphism(M) - sage: m = M(range(0,16)) - sage: p_inv(p(m)) == m + This function computes the correct coordinate representations (of + length 3) for a basis of the space of two-by-two symmetric + matrices, the the inverse does indeed invert the process:: + + sage: E11 = matrix(QQbar,[ [1,0], + ....: [0,0] ]) + sage: E12 = matrix(QQbar,[ [0, 1/sqrt(2)], + ....: [1/sqrt(2), 0] ]) + sage: E22 = matrix(QQbar,[ [0,0], + ....: [0,1] ]) + sage: basis = [E11, E12, E22] + sage: phi, phi_inv = basis_representation(basis) + sage: phi(E11); phi(E12); phi(E22) + (1, 0, 0) + (0, 1, 0) + (0, 0, 1) + sage: phi_inv(phi(E11)) == E11 + True + sage: phi_inv(phi(E12)) == E12 + True + sage: phi_inv(phi(E22)) == E22 + True + + MatrixSpace arguments work too:: + + sage: M = MatrixSpace(QQ,2) + sage: phi, phi_inv = basis_representation(M) + sage: X = matrix(QQ, [ [1,2], + ....: [3,4] ]) + sage: phi(X) + (1, 2, 3, 4) + sage: phi_inv(phi(X)) == X + True + + TESTS: + + The inverse is generally an inverse:: + + sage: set_random_seed() + sage: n = ZZ.random_element(10) + sage: M = MatrixSpace(QQ,n) + sage: X = M.random_element() + sage: (phi, phi_inv) = basis_representation(M) + sage: phi_inv(phi(X)) == X + True + + """ + if is_MatrixSpace(M): + basis_space = M + basis = list(M.basis()) + else: + basis_space = M[0].matrix_space() + basis = M + + def phi(X): + """ + The isometry sending ``X`` to its representation as a long vector. + """ + if X not in basis_space: + raise ValueError("X does not live in the domain of phi") + + V = VectorSpace(basis_space.base_ring(), X.nrows()*X.ncols()) + W = V.span_of_basis( _mat2vec(s) for s in basis ) + return W.coordinate_vector(_mat2vec(X)) + + def phi_inv(Y): + """ + The isometry sending the long vector `Y` to an element of either + `M` or the span of `M` (depending on whether or not ``M`` + is a ``MatrixSpace`` or a basis). + """ + return basis_space.linear_combination( zip(Y,basis) ) + + return (phi, phi_inv) + + + +def basis_repr_of_operator(M, L): + """ + Return the matrix of the operator `L` with respect to the basis + `M` if `M` is a list of basis vectors for a matrix space; or with + respect to the standard basis of `M` if `M` is a ``MatrixSpace``. + + This function is necessary because SageMath does not know that + matrix spaces are vector spaces, and it moreover it doesn't know + that (for example) the subspace of symmetric matrices is a matrix + space in its own right. + + Use ``linear_transformation().matrix()`` instead if you have a + true ``VectorSpace``. + + INPUT: + + - ``M`` -- Either a ``MatrixSpace``, or a list of matrices that form + a basis for a matrix space. + + OUTPUT: + + If the matrix space associated with `M` has dimension `n`, then an + `n`-by-`n` matrix over the same base ring is returned. + + SETUP:: + + sage: from mjo.matrix_vector import (basis_representation, + ....: basis_repr_of_operator) + + EXAMPLES: + + The matrix of the identity operator on the space of two-by-two + symmetric matrices is the identity matrix, regardless of the basis:: + + sage: E11 = matrix(QQbar,[ [1,0], + ....: [0,0] ]) + sage: E12 = matrix(QQbar,[ [0, 1/sqrt(2)], + ....: [1/sqrt(2), 0] ]) + sage: E22 = matrix(QQbar,[ [0,0], + ....: [0,1] ]) + sage: basis = [E11, E12, E22] + sage: identity = lambda X: X + sage: basis_repr_of_operator(basis, identity) + [1 0 0] + [0 1 0] + [0 0 1] + sage: E11 = matrix(QQ,[[2,0],[0,0]]) + sage: E12 = matrix(QQ,[[0,2],[2,0]]) + sage: basis = [E11, E12, E22] + sage: basis_repr_of_operator(basis, identity) + [1 0 0] + [0 1 0] + [0 0 1] + + A more complicated example confirms that we get a matrix consistent + with our ``matrix_to_vector`` function:: + + sage: M = MatrixSpace(QQ,3,3) + sage: Q = M([[0,1,0],[1,0,0],[0,0,1]]) + sage: def f(x): + ....: return Q*x*Q.inverse() + ....: + sage: F = basis_repr_of_operator(M, f) + sage: F + [0 0 0 0 1 0 0 0 0] + [0 0 0 1 0 0 0 0 0] + [0 0 0 0 0 1 0 0 0] + [0 1 0 0 0 0 0 0 0] + [1 0 0 0 0 0 0 0 0] + [0 0 1 0 0 0 0 0 0] + [0 0 0 0 0 0 0 1 0] + [0 0 0 0 0 0 1 0 0] + [0 0 0 0 0 0 0 0 1] + sage: phi, phi_inv = basis_representation(M) + sage: X = M([[1,2,3],[4,5,6],[7,8,9]]) + sage: F*phi(X) == phi(f(X)) True """ - from sage.matrix.matrix_space import is_MatrixSpace - if not is_MatrixSpace(matrix_space): - raise TypeError('argument must be a matrix space') + if is_MatrixSpace(M): + basis_space = M + basis = list(M.basis()) + else: + basis_space = M[0].matrix_space() + basis = M - base_ring = matrix_space.base_ring() - vector_space = VectorSpace(base_ring, matrix_space.dimension()) + (phi, phi_inv) = basis_representation(M) - def phi(m): - return vector_space(m.list()) + # Get a basis for the image space. Since phi is an isometry, + # it takes one basis to another. + image_basis = [ phi(b) for b in basis ] - def phi_inverse(v): - return matrix_space(v.list()) + # Now construct the image space itself equipped with our custom basis. + W = VectorSpace(basis_space.base_ring(), len(basis)) + W = W.span_of_basis(image_basis) - return (phi, phi_inverse) + return matrix.column( W.coordinates(phi(L(b))) for b in basis )