"""
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 isomorphism
+ sage: from mjo.matrix_vector import basis_representation
+
+ 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_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:
- EXAMPLES::
+ The inverse is generally an inverse::
- sage: M = MatrixSpace(QQ,4,4)
- sage: (p, p_inv) = isomorphism(M)
- sage: m = M(range(0,16))
- sage: p_inv(p(m)) == m
+ 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
"""
- 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())
+ 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")
- def phi(m):
- return vector_space(m.list())
+ 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_inverse(v):
- return matrix_space(v.list())
+ 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_inverse)
+ return (phi, phi_inv)
-def matrix_of_transformation(T, V):
+def basis_repr_of_operator(M, L):
"""
- Compute the matrix of a linear transformation ``T``, `$T : V
- \rightarrow V$` with domain/range ``V``. This essentially uses the
- Riesz representation theorem to represent the entries of the matrix
- of ``T`` in terms of inner products.
+ 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``.
- INPUT:
+ 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``.
- - ``T`` -- The linear transformation whose matrix we should
- compute. This should be a callable function that
- takes as its single argument an element of ``V``.
+ INPUT:
- - ``V`` -- The vector or matrix space on which ``T`` is defined.
+ - ``M`` -- Either a ``MatrixSpace``, or a list of matrices that form
+ a basis for a matrix space.
OUTPUT:
- If the dimension of ``V`` is `$n$`, we return an `$n \times n$`
- matrix that represents ``T`` with respect to the standard basis of
- ``V``.
+ 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 isomorphism, matrix_of_transformation
+ sage: from mjo.matrix_vector import (basis_representation,
+ ....: basis_repr_of_operator)
EXAMPLES:
- The matrix of a transformation on a simple vector space should be
- the expected matrix::
-
- sage: V = VectorSpace(QQ, 3)
- sage: def f(x):
- ....: return 3*x
- ....:
- sage: matrix_of_transformation(f, V)
- [3 0 0]
- [0 3 0]
- [0 0 3]
+ 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: def f(x):
....: return Q*x*Q.inverse()
....:
- sage: F = matrix_of_transformation(f, M)
+ 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 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 = isomorphism(M)
+ sage: phi, phi_inv = basis_representation(M)
sage: X = M([[1,2,3],[4,5,6],[7,8,9]])
- sage: F*phi(X)
- (5, 4, 6, 2, 1, 3, 8, 7, 9)
- sage: phi(f(X))
- (5, 4, 6, 2, 1, 3, 8, 7, 9)
sage: F*phi(X) == phi(f(X))
True
"""
- n = V.dimension()
- B = list(V.basis())
-
- def inner_product(v, w):
- # An inner product function that works for both matrices and
- # vectors.
- if callable(getattr(v, 'inner_product', None)):
- return v.inner_product(w)
- elif callable(getattr(v, 'matrix_space', None)):
- # V must be a matrix space?
- return (v*w.transpose()).trace()
- else:
- raise ValueError('inner_product only works on vectors and matrices')
-
- def apply(L, x):
- # A "call" that works for both matrices and functions.
- if callable(getattr(L, 'matrix_space', None)):
- # L is a matrix, and we need to use "multiply" to call it.
- return L*x
- else:
- # If L isn't a matrix, try this. It works for python
- # functions at least.
- return L(x)
-
- entries = []
- for j in range(0,n):
- for i in range(0,n):
- entry = inner_product(apply(T,B[i]), B[j])
- entries.append(entry)
-
- # Construct the matrix space in which our return value will lie.
- W = MatrixSpace(V.base_ring(), n, n)
-
- # And make a matrix out of our list of entries.
- A = W(entries)
-
- return A
+ if is_MatrixSpace(M):
+ basis_space = M
+ basis = list(M.basis())
+ else:
+ basis_space = M[0].matrix_space()
+ basis = M
+
+ (phi, phi_inv) = basis_representation(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 )