From ba04841bc8b5786c75b5e6fa9e767895c6af51d0 Mon Sep 17 00:00:00 2001 From: Michael Orlitzky Date: Sun, 13 Sep 2020 16:33:06 -0400 Subject: [PATCH] mjo/matrix_vector: replace matrix_of_transformation with basis_repr_of_operator. --- mjo/matrix_vector.py | 125 ++++++++++++++++++++----------------------- 1 file changed, 58 insertions(+), 67 deletions(-) diff --git a/mjo/matrix_vector.py b/mjo/matrix_vector.py index 62b75fd..83bf313 100644 --- a/mjo/matrix_vector.py +++ b/mjo/matrix_vector.py @@ -127,45 +127,59 @@ def basis_representation(M): -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. - - ``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``. + Use ``linear_transformation().matrix()`` instead if you have a + true ``VectorSpace``. - - ``V`` -- The vector or matrix space on which ``T`` is defined. + INPUT: + + - ``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 (basis_representation, - ....: matrix_of_transformation) + ....: basis_repr_of_operator) EXAMPLES: - The matrix of a transformation on a simple vector space should be - the expected matrix:: + The matrix of the identity operator on the space of two-by-two + symmetric matrices is the identity matrix, regardless of the basis:: - 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] + 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:: @@ -175,7 +189,7 @@ def matrix_of_transformation(T, V): 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] @@ -186,50 +200,27 @@ def matrix_of_transformation(T, V): [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(n): - for i in range(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 ) -- 2.44.2