]> gitweb.michael.orlitzky.com - sage.d.git/blobdiff - mjo/basis_repr.py
mjo: rename matrix_vector.py to basis_repr.py.
[sage.d.git] / mjo / basis_repr.py
diff --git a/mjo/basis_repr.py b/mjo/basis_repr.py
new file mode 100644 (file)
index 0000000..5c85998
--- /dev/null
@@ -0,0 +1,230 @@
+r"""
+In an `n`-dimensional vector space, representation with respect to
+a basis is an isometry between that space and `\mathbb{R}^{n}`.
+
+Sage is able to go back/forth relatively easy when you start with a
+``VectorSpace``, but unfortunately, it does not know that a
+``MatrixSpace`` is also a ``VectorSpace``. So, this module exists to
+perform the "basis representation" isometry between a matrix space and
+a vector space of the same dimension.
+
+"""
+
+from sage.all import *
+from sage.matrix.matrix_space import is_MatrixSpace
+
+def _mat2vec(m):
+    return vector(m.base_ring(), m.list())
+
+def basis_repr(M):
+    """
+    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:
+
+    - ``M`` -- Either a ``MatrixSpace``, or a list of matrices that form
+               a basis for a matrix space.
+
+    OUTPUT:
+
+    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.basis_repr import basis_repr
+
+    EXAMPLES:
+
+    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_repr(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_repr(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_repr(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.basis_repr import (basis_repr,
+        ....:                             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_repr(M)
+        sage: X = M([[1,2,3],[4,5,6],[7,8,9]])
+        sage: F*phi(X) == phi(f(X))
+        True
+
+    """
+    if is_MatrixSpace(M):
+        basis_space = M
+        basis = list(M.basis())
+    else:
+        basis_space = M[0].matrix_space()
+        basis = M
+
+    (phi, phi_inv) = basis_repr(M)
+
+    # 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 ]
+
+    # 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 matrix.column( W.coordinates(phi(L(b))) for b in basis )