+++ /dev/null
-def linear_isomorphisms(K1, K2):
- r"""
- Generate invertible linear transformations mapping ``K1``
- onto ``K2``.
-
- In general, if ``K1`` and ``K2`` are isomorphic, there will be an
- infinite number of isomorphisms between them. For example, any
- positive scalar multiple of a cone isomorphism is itself a cone
- isomorphism. As such we do not return *every* isomorphism, only a
- "generating set." For the precise meaning of this, please see the
- description of the algorithm below.
-
- INPUT:
-
- - ``K1`` -- the first cone
-
- - ``K2`` -- :class:`sage.geometry.cone.ConvexRationalPolyhedralCone`;
- the cone to generate isomorphisms with
-
- OUTPUT:
-
- A generator that yields linear isomorphisms between ``K1`` and
- ``K2``, one at a time. The isomorphisms consist of invertible
- matrices ``A`` whose entries are rational numbers and satisfy
- ``Cone(K1.rays()*A).is_equivalent(K2)``. In other words, they map
- the cone ``K1`` onto ``K2``, acting on the rays of ``K1`` from the
- right. The right-action is mainly for compatibility with the ray
- matrices of cones and the basis matrices of vector spaces, which
- in Sage are "lists of rows."
-
- The returned matrices are immutable so that they may be hashed,
- and, for example, de-duplicated using ``set()``.
-
- ALGORITHM:
-
- It is well-known that the linear automorphisms of a cone that
- :meth:`is_strictly_convex` must map extreme rays of the original
- cone (which exist, because it is strictly convex) to extreme rays
- of the target cone (which had better exist, otherwise the cones
- are not isomorphic). If moreover ``K1`` meth:`is_solid`, then it
- has at least `n` extreme rays, where `n` is the dimension of its
- ambient space, and likewise for ``K2`` (otherwise, again, they are
- not isomorphic). As a result, if ``K1`` is both strictly convex
- and solid, any isomorphism must map a basis consisting of its
- extreme rays to a basis of extreme rays for ``K2``. In this
- scenario we generate all such maps, and yield the ones that turn
- out to be isomorphisms when we check them explicitly.
-
- If the cones are both solid and strictly convex, the generated
- isomorphisms are unique up to positive row scalings. When there
- are exactly `n` extreme rays (that is, when the cone
- :meth:`is_simplicial`), each row of an isomorphism can be scaled
- independently without destroying the isomorphism property. But
- when there are more than `n` extreme rays, it is typically only
- safe to scale the entire isomorphism by a single positive amount.
-
- If ``K1`` and ``K2`` are *not* pointed, then in addition to the
- procedure above, we map a basis for the :meth:`linear_subspace` of
- ``K1`` to a basis of the :meth:`linear_subspace` of
- ``K2``. Similarly, if ``K1`` and ``K2`` are not solid, we map
- bases for their orthogonal complements to one another. There are
- many isomorphisms (corresponding to many choices of bases) between
- two vector subspaces of the same dimension, so in this case, the
- cone isomorphisms we generate should be considered unique only up
- to vector-space isomorphisms of the lineality spaces and othogonal
- complements.
-
- REFERENCES:
-
- .. [GowdaTrott2014] \M. Seetharama Gowda and D. Trott. On the
- irreducibility, Lyapunov rank, and automorphisms of special
- Bishop–Phelps cones. Journal of Mathematical Analysis and
- Applications, 419(1):172-184, 2014.
-
- SETUP::
-
- sage: from mjo.cone.isomorphism import linear_isomorphisms
-
- EXAMPLES:
-
- In this example, the lower half-space is obtained from the upper
- half space by flipping the y-coordinate. Note that the x-axis is
- mapped to the x-axis via the canonical basis (the ``1`` in the
- top-left corner), but any other non-zero number would suffice
- there. This is discussed in the description of the algorithm::
-
- sage: K1 = Cone([(0,1),(1,0),(-1,0)])
- sage: K2 = Cone([(0,-1),(1,0),(-1,0)])
- sage: set(linear_isomorphisms(K1,K2))
- {[ 1 0]
- [ 0 -1]}
- sage: A = matrix(QQ, [[-6, 0],
- ....: [ 0,-1]])
- sage: Cone(K1.rays()*A).is_equivalent(K2)
- True
-
- A randomly-generated example constructed to be isomorphic::
-
- sage: K1 = Cone( [(17, 55, -13, 0),
- ....: (-9, -32, 914, 0),
- ....: ( 1, 16, -2, 7),
- ....: (-1, -16, 2, -7)])
- sage: A = matrix(QQ, [[ 1, 4, 7, 43],
- ....: [ 0, 1, 3, 16],
- ....: [ 2, 6, 9, 58],
- ....: [ -1, -7, -16, -90]])
- sage: K2 = Cone(K1.rays()*A, K1.lattice())
- sage: g = next(linear_isomorphisms(K1,K2))
- sage: Cone(K1.rays()*g).is_equivalent(K2)
- True
-
- Automorphisms can be obtained by passing ``K2 == K1``. There are
- often many duplicates, so we use ``set()`` to obtain distinct
- transformations. Gowda and Trott [GowdaTrott2014]_ have computed
- the automorphism groups of these cone, and we recover them all up
- to a positive scalar::
-
- sage: K1 = Cone([(1,0,1), (-1,0,1), (0,1,1), (0,-1,1)])
- sage: G = set(linear_isomorphisms(K1,K1)); G
- {[-1 0 0]
- [ 0 -1 0]
- [ 0 0 1],
- [-1 0 0]
- [ 0 1 0]
- [ 0 0 1],
- [ 0 -1 0]
- [-1 0 0]
- [ 0 0 1],
- [ 0 -1 0]
- [ 1 0 0]
- [ 0 0 1],
- [ 0 1 0]
- [-1 0 0]
- [ 0 0 1],
- [0 1 0]
- [1 0 0]
- [0 0 1],
- [ 1 0 0]
- [ 0 -1 0]
- [ 0 0 1],
- [1 0 0]
- [0 1 0]
- [0 0 1]}
- sage: K2 = Cone([(1,1,1), (-1,1,1), (-1,-1,1), (1,-1,1)])
- sage: H = set(linear_isomorphisms(K2,K2))
- sage: G == H
- True
-
- Up to a positive row-scalings, the automorphism group of the
- nonnegative orthant is the set of all permutation matrices. Only
- the permutations (not the scalar factors) are returned; this is
- discussed in the description of the algorithm::
-
- sage: K1 = cones.nonnegative_orthant(3)
- sage: set(linear_isomorphisms(K1,K1))
- {[0 0 1]
- [0 1 0]
- [1 0 0],
- [0 0 1]
- [1 0 0]
- [0 1 0],
- [0 1 0]
- [0 0 1]
- [1 0 0],
- [0 1 0]
- [1 0 0]
- [0 0 1],
- [1 0 0]
- [0 0 1]
- [0 1 0],
- [1 0 0]
- [0 1 0]
- [0 0 1]}
-
- TESTS:
-
- Every cone should be isomorphic to itself under (at least) the
- identity transformation::
-
- sage: K = random_cone(max_ambient_dim=5)
- sage: I = matrix.identity(QQ, K.lattice_dim())
- sage: I in linear_isomorphisms(K,K)
- True
-
- Check the properties of the generated isomorphisms. We build
- ``K2`` from ``K1`` using an invertible rational matrix, so we know
- that there is at least one isomorphism between them::
-
- sage: K1 = random_cone(max_ambient_dim=5)
- sage: R = K1.lattice().vector_space().base_ring()
- sage: n = K1.lattice_dim()
- sage: A = matrix.random(R, n, algorithm='unimodular')
- sage: q = QQ.random_element()
- sage: while q.is_zero():
- ....: q = QQ.random_element()
- sage: K2 = Cone(K1.rays()*A, K1.lattice())
- sage: all(
- ....: g.is_invertible()
- ....: and
- ....: g.is_immutable()
- ....: and
- ....: Cone(K1.rays()*g, K1.lattice()).is_equivalent(K2)
- ....: and
- ....: Cone(K2.rays()*g.inverse(), K2.lattice()).is_equivalent(K1)
- ....: for g in linear_isomorphisms(K1,K2)
- ....: )
- True
- """
- # There are no invertible maps between lattices of different
- # dimensions.
- if K1.lattice_dim() != K2.lattice_dim():
- return
-
- # Create two vector spaces V and W (resp.) to hold K1 and K2,
- # and find their lineality / perp subspaces.
- V = K1.lattice().vector_space()
- L = K1.linear_subspace()
- K1_perp = K1.span().vector_space().complement()
- W = K2.lattice().vector_space()
- M = K2.linear_subspace()
- K2_perp = K2.span().vector_space().complement()
-
- if L.dimension() != M.dimension():
- # Linealities don't match, not isomorphic.
- return
-
- if K1_perp.dimension() != K2_perp.dimension():
- # Cone dimensions don't match, not isomorphic.
- return
-
- # Standard trick for orthogonal projections onto L and M.
- A = L.basis_matrix()
- L_proj = A.T * (A*A.T).inverse() * A
- B = M.basis_matrix()
- M_proj = B.T * (B*B.T).inverse() * B
-
- # Project onto the lineality space and subtract to obtain the
- # "pointed parts"...
- L_perps = [ r - r*L_proj for r in K1 ]
- M_perps = [ r - r*M_proj for r in K2 ]
-
- # ...which need to be immutable to be put into sets.
- for i in range(len(L_perps)):
- L_perps[i].set_immutable()
- for j in range(len(M_perps)):
- M_perps[j].set_immutable()
-
- # Our generators are assumed to be minimal, so none should be
- # eliminated by the Cone(...) below, although some may be
- # reordered (which is why we still use L_perps/M_perps directly
- # even after constructing these cones). In any case, if we wind up
- # with different numbers of rays here, we aren't going to get any
- # isomorphisms. This is slow, but it's the most common case (cones
- # not isomorphic), so probably worth doing before proceeding.
- from sage.geometry.cone import Cone
- K1_p = Cone(L_perps, K1.lattice())
- K2_p = Cone(M_perps, K2.lattice())
- if not K1_p.nrays() == K2_p.nrays():
- return
-
- # The dimension of the non-lineal part of our cones. Needs to be a
- # python int so that we can pass it to itertools.permutations.
- n = int(K1.dim() - K1.lineality())
-
- # The base ring of the matrices we'll construct later. We pass it
- # explicitly to avoid the rare situation where our rational
- # matrices all contain integral entries, making sage think that we
- # want an integral matrix. This becomes important when passing
- # ``extend=False`` to solve_right().
- R = K1.lattice().vector_space().base_ring() # QQ
-
- # If the non-lineal parts of our cones are dimension zero, both
- # cones are subspaces of the same dimension, and the identity map
- # between them should work. This special case is needed to avoid
- # solving linear equations with the empty set below.
- from sage.matrix.constructor import matrix
- if n == 0:
- X = matrix.identity(R, V.dimension())
- X.set_immutable()
- yield X
- return
-
- # Compute these outside of the loop. These are extra rows that we
- # append the the system of equations to ensure that the lineality
- # and perp spaces of K1 get mapped to those of K2.
- K1_extra_rows = A.rows() + K1_perp.basis_matrix().rows()
- K2_extra_rows = B.rows() + K2_perp.basis_matrix().rows()
-
- # Now try to find isomorphisms in the first component (non-lineal
- # part), and extend them to the whole thing. Use lists to keep
- # everything in the right order.
- from itertools import permutations
- for rs1 in permutations(L_perps, n):
- for rs2 in permutations(M_perps, n):
- # Try mapping every subset of the right size in the domain
- # to every subset of the right size in the codomain. Extend
- # to the lineal part, and then check that the result is
- # actually an isomorphism.
- A = matrix(R, list(rs1) + K1_extra_rows)
- B = matrix(R, list(rs2) + K2_extra_rows)
- try:
- X = A.solve_right(B, extend=False)
- except ValueError:
- continue
-
- if Cone(K1.rays()*X).is_equivalent(K2):
- X.set_immutable()
- yield X
-
-
-def is_linearly_isomorphic(K1,K2):
- r"""
- Return ``True`` if ``K1`` and ``K2`` are linearly isomorphic
- over the rationals, and ``False`` otherwise.
-
- INPUT:
-
- - ``K1`` -- the first cone
-
- - ``K2`` -- :class:`sage.geometry.cone.ConvexRationalPolyhedralCone`;
- the cone to check for isomorphism with ``K1``
-
- OUTPUT:
-
- ``True`` if there exists an invertible linear matrix with rational
- entries mapping ``K1`` to ``K2``, and ``False`` otherwise.
-
- ALGORITHM:
-
- We attempt to construct an explicit isomorphism using
- :meth:`linear_isomorphisms`. For more information, see the
- documentation of that method.
-
- SETUP::
-
- sage: from mjo.cone.isomorphism import is_linearly_isomorphic
-
- EXAMPLES:
-
- Two different descriptions of the plane should be isomorphic::
-
- sage: K1 = Cone([(1,0), (-1,0), (0,1), (0,-1)])
- sage: K2 = Cone([(1,1), (-1,1), (0,-1)])
- sage: K1.is_equivalent(K2)
- True
- sage: is_linearly_isomorphic(K1,K2)
- True
-
- All simplicial cones with the same number of rays are isomorphic::
-
- sage: K1 = random_cone(max_ambient_dim=5)
- sage: while not K1.is_simplicial():
- ....: K1 = random_cone(max_ambient_dim=5)
- sage: n = K1.lattice_dim()
- sage: r = K1.nrays()
- sage: K2 = random_cone(min_ambient_dim=n,
- ....: max_ambient_dim=n,
- ....: min_rays=r,
- ....: max_rays=r,
- ....: strictly_convex=True)
- sage: K2.is_simplicial()
- True
- sage: is_linearly_isomorphic(K1,K2)
- True
-
- TESTS:
-
- Every cone is isomorphic to itself::
-
- sage: K = random_cone(max_ambient_dim=5)
- sage: is_linearly_isomorphic(K,K)
- True
-
- Isomorphism is symmetric::
-
- sage: K = random_cone(max_ambient_dim=5)
- sage: J = random_cone(min_ambient_dim=K.lattice_dim(),
- ....: max_ambient_dim=K.lattice_dim(),
- ....: min_rays=K.nrays(),
- ....: max_rays=K.nrays())
- sage: is_linearly_isomorphic(K,J) == is_linearly_isomorphic(K,J)
- True
-
- """
- try:
- next(linear_isomorphisms(K1,K2))
- return True
- except StopIteration:
- return False