Add positive/Z tests and update code for upstream changes.

[sage.d.git] / mjo / cone / cone.py

diff --git a/mjo/cone/cone.py b/mjo/cone/cone.py
index 21f9862c24a9e9e3d9cffc52a1e789018f59a153..a7c16a520cbe19242398e176635ad3250e5c7671 100644 (file)
--- a/mjo/cone/cone.py
+++ b/mjo/cone/cone.py
@@ -1,14 +1,12 @@
 from sage.all import *
 
 from sage.all import *
 

-def is_lyapunov_like(L,K):
+def is_positive_on(L,K):

     r"""
     r"""

-    Determine whether or not ``L`` is Lyapunov-like on ``K``.
+    Determine whether or not ``L`` is positive on ``K``.

 
 

-    We say that ``L`` is Lyapunov-like on ``K`` if `\left\langle
-    L\left\lparenx\right\rparen,s\right\rangle = 0` for all pairs
-    `\left\langle x,s \right\rangle` in the complementarity set of
-    ``K``. It is known [Orlitzky]_ that this property need only be
-    checked for generators of ``K`` and its dual.
+    We say that ``L`` is positive on ``K`` if `L\left\lparen x
+    \right\rparen` belongs to ``K`` for all `x` in ``K``. This
+    property need only be checked for generators of ``K``.

 
     INPUT:
 
 
     INPUT:
 


@@ -18,31 +16,26 @@ def is_lyapunov_like(L,K):
 
     OUTPUT:
 
 
     OUTPUT:
 

-    ``True`` if it can be proven that ``L`` is Lyapunov-like on ``K``,
+    ``True`` if it can be proven that ``L`` is positive on ``K``,

     and ``False`` otherwise.
 
     .. WARNING::
 
     and ``False`` otherwise.
 
     .. WARNING::
 

-        If this function returns ``True``, then ``L`` is Lyapunov-like
+        If this function returns ``True``, then ``L`` is positive

         on ``K``. However, if ``False`` is returned, that could mean one
         of two things. The first is that ``L`` is definitely not
         on ``K``. However, if ``False`` is returned, that could mean one
         of two things. The first is that ``L`` is definitely not

-        Lyapunov-like on ``K``. The second is more of an "I don't know"
+        positive on ``K``. The second is more of an "I don't know"

         answer, returned (for example) if we cannot prove that an inner
         answer, returned (for example) if we cannot prove that an inner

-        product is zero.
-
-    REFERENCES:
-
-    M. Orlitzky. The Lyapunov rank of an improper cone.
-    http://www.optimization-online.org/DB_HTML/2015/10/5135.html
+        product is nonnegative.

 
     EXAMPLES:
 
 
     EXAMPLES:
 

-    The identity is always Lyapunov-like in a nontrivial space::
+    The identity is always positive in a nontrivial space::

 
         sage: set_random_seed()
         sage: K = random_cone(min_ambient_dim=1, max_ambient_dim=8)
         sage: L = identity_matrix(K.lattice_dim())
 
         sage: set_random_seed()
         sage: K = random_cone(min_ambient_dim=1, max_ambient_dim=8)
         sage: L = identity_matrix(K.lattice_dim())

-        sage: is_lyapunov_like(L,K)
+        sage: is_positive_on(L,K)

         True
 
     As is the "zero" transformation::
         True
 
     As is the "zero" transformation::


@@ -50,683 +43,268 @@ def is_lyapunov_like(L,K):
         sage: K = random_cone(min_ambient_dim=1, max_ambient_dim=8)
         sage: R = K.lattice().vector_space().base_ring()
         sage: L = zero_matrix(R, K.lattice_dim())
         sage: K = random_cone(min_ambient_dim=1, max_ambient_dim=8)
         sage: R = K.lattice().vector_space().base_ring()
         sage: L = zero_matrix(R, K.lattice_dim())

-        sage: is_lyapunov_like(L,K)
+        sage: is_positive_on(L,K)

         True
 
         True
 

-        Everything in ``K.lyapunov_like_basis()`` should be Lyapunov-like
-        on ``K``::
+    TESTS:
+
+    Everything in ``K.positive_operators_gens()`` should be
+    positive on ``K``::

 
         sage: K = random_cone(min_ambient_dim=1, max_ambient_dim=6)
 
         sage: K = random_cone(min_ambient_dim=1, max_ambient_dim=6)

-        sage: all([ is_lyapunov_like(L,K) for L in K.lyapunov_like_basis() ])
+        sage: all([ is_positive_on(L,K)
+        ....:       for L in K.positive_operators_gens() ])
+        True
+        sage: all([ is_positive_on(L.change_ring(SR),K)
+        ....:       for L in K.positive_operators_gens() ])

         True
 
     """
         True
 
     """

-    return all([(L*x).inner_product(s) == 0
-                for (x,s) in K.discrete_complementarity_set()])
+    if L.base_ring().is_exact():
+        # This could potentially be extended to other types of ``K``...
+        return all([ L*x in K for x in K ])
+    elif L.base_ring() is SR:
+        # Fall back to inequality-checking when the entries of ``L``
+        # might be symbolic.
+        return all([ s*(L*x) >= 0 for x in K for s in K ])
+    else:
+        # The only inexact ring that we're willing to work with is SR,
+        # since it can still be exact when working with symbolic
+        # constants like pi and e.
+        raise ValueError('base ring of operator L is neither SR nor exact')
+
+
+def is_cross_positive_on(L,K):
+    r"""
+    Determine whether or not ``L`` is cross-positive on ``K``.

 
 

+    We say that ``L`` is cross-positive on ``K`` if `\left\langle
+    L\left\lparenx\right\rparen,s\right\rangle \ge 0` for all pairs
+    `\left\langle x,s \right\rangle` in the complementarity set of
+    ``K``. This property need only be checked for generators of
+    ``K`` and its dual.

 
 

-def motzkin_decomposition(K):
-    r"""
-    Return the pair of components in the Motzkin decomposition of this cone.
+    INPUT:

 
 

-    Every convex cone is the direct sum of a strictly convex cone and a
-    linear subspace [Stoer-Witzgall]_. Return a pair ``(P,S)`` of cones
-    such that ``P`` is strictly convex, ``S`` is a subspace, and ``K``
-    is the direct sum of ``P`` and ``S``.
+    - ``L`` -- A linear transformation or matrix.
+
+    - ``K`` -- A polyhedral closed convex cone.

 
     OUTPUT:
 
 
     OUTPUT:
 

-    An ordered pair ``(P,S)`` of closed convex polyhedral cones where
-    ``P`` is strictly convex, ``S`` is a subspace, and ``K`` is the
-    direct sum of ``P`` and ``S``.
+    ``True`` if it can be proven that ``L`` is cross-positive on ``K``,
+    and ``False`` otherwise.

 
 

-    REFERENCES:
+    .. WARNING::

 
 

-    .. [Stoer-Witzgall] J. Stoer and C. Witzgall. Convexity and
-       Optimization in Finite Dimensions I. Springer-Verlag, New
-       York, 1970.
+        If this function returns ``True``, then ``L`` is cross-positive
+        on ``K``. However, if ``False`` is returned, that could mean one
+        of two things. The first is that ``L`` is definitely not
+        cross-positive on ``K``. The second is more of an "I don't know"
+        answer, returned (for example) if we cannot prove that an inner
+        product is nonnegative.

 
     EXAMPLES:
 
 
     EXAMPLES:
 

-    The nonnegative orthant is strictly convex, so it is its own
-    strictly convex component and its subspace component is trivial::
+    The identity is always cross-positive in a nontrivial space::

 
 

-        sage: K = Cone([(1,0,0),(0,1,0),(0,0,1)])
-        sage: (P,S) = motzkin_decomposition(K)
-        sage: K.is_equivalent(P)
-        True
-        sage: S.is_trivial()
+        sage: set_random_seed()
+        sage: K = random_cone(min_ambient_dim=1, max_ambient_dim=8)
+        sage: L = identity_matrix(K.lattice_dim())
+        sage: is_cross_positive_on(L,K)

         True
 
         True
 

-    Likewise, full spaces are their own subspace components::
+    As is the "zero" transformation::

 
 

-        sage: K = Cone([(1,0),(-1,0),(0,1),(0,-1)])
-        sage: K.is_full_space()
-        True
-        sage: (P,S) = motzkin_decomposition(K)
-        sage: K.is_equivalent(S)
-        True
-        sage: P.is_trivial()
+        sage: K = random_cone(min_ambient_dim=1, max_ambient_dim=8)
+        sage: R = K.lattice().vector_space().base_ring()
+        sage: L = zero_matrix(R, K.lattice_dim())
+        sage: is_cross_positive_on(L,K)

         True
 
     TESTS:
 
         True
 
     TESTS:
 

-    A random point in the cone should belong to either the strictly
-    convex component or the subspace component. If the point is nonzero,
-    it cannot be in both::
-
-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=8)
-        sage: (P,S) = motzkin_decomposition(K)
-        sage: x = K.random_element(ring=QQ)
-        sage: P.contains(x) or S.contains(x)
-        True
-        sage: x.is_zero() or (P.contains(x) != S.contains(x))
-        True
+    Everything in ``K.cross_positive_operators_gens()`` should be
+    cross-positive on ``K``::

 
 

-    The strictly convex component should always be strictly convex, and
-    the subspace component should always be a subspace::
-
-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=8)
-        sage: (P,S) = motzkin_decomposition(K)
-        sage: P.is_strictly_convex()
+        sage: K = random_cone(min_ambient_dim=1, max_ambient_dim=6)
+        sage: all([ is_cross_positive_on(L,K)
+        ....:       for L in K.cross_positive_operators_gens() ])

         True
         True

-        sage: S.lineality() == S.dim()
+        sage: all([ is_cross_positive_on(L.change_ring(SR),K)
+        ....:       for L in K.cross_positive_operators_gens() ])

         True
 
         True
 

-    A strictly convex cone should be equal to its strictly convex component::
-
-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=8, strictly_convex=True)
-        sage: (P,_) = motzkin_decomposition(K)
-        sage: K.is_equivalent(P)
-        True
-
-    The generators of the components are obtained from orthogonal
-    projections of the original generators [Stoer-Witzgall]_::
-
-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=8)
-        sage: (P,S) = motzkin_decomposition(K)
-        sage: A = S.linear_subspace().complement().matrix()
-        sage: proj_S_perp = A.transpose() * (A*A.transpose()).inverse() * A
-        sage: expected_P = Cone([ proj_S_perp*g for g in K ], K.lattice())
-        sage: P.is_equivalent(expected_P)
-        True
-        sage: A = S.linear_subspace().matrix()
-        sage: proj_S = A.transpose() * (A*A.transpose()).inverse() * A
-        sage: expected_S = Cone([ proj_S*g for g in K ], K.lattice())
-        sage: S.is_equivalent(expected_S)
-        True

     """
     """

-    # The lines() method only returns one generator per line. For a true
-    # line, we also need a generator pointing in the opposite direction.
-    S_gens = [ direction*gen for direction in [1,-1] for gen in K.lines() ]
-    S = Cone(S_gens, K.lattice(), check=False)
+    if L.base_ring().is_exact() or L.base_ring() is SR:
+        return all([ s*(L*x) >= 0
+                     for (x,s) in K.discrete_complementarity_set() ])
+    else:
+        # The only inexact ring that we're willing to work with is SR,
+        # since it can still be exact when working with symbolic
+        # constants like pi and e.
+        raise ValueError('base ring of operator L is neither SR nor exact')

 
 

-    # Since ``S`` is a subspace, the rays of its dual generate its
-    # orthogonal complement.
-    S_perp = Cone(S.dual(), K.lattice(), check=False)
-    P = K.intersection(S_perp)

 
 

-    return (P,S)
-
-
-def positive_operator_gens(K):
+def is_Z_on(L,K):

     r"""
     r"""

-    Compute generators of the cone of positive operators on this cone.
-
-    OUTPUT:
-
-    A list of `n`-by-``n`` matrices where ``n == K.lattice_dim()``.
-    Each matrix ``P`` in the list should have the property that ``P*x``
-    is an element of ``K`` whenever ``x`` is an element of
-    ``K``. Moreover, any nonnegative linear combination of these
-    matrices shares the same property.
-
-    EXAMPLES:
-
-    Positive operators on the nonnegative orthant are nonnegative matrices::
-
-        sage: K = Cone([(1,)])
-        sage: positive_operator_gens(K)
-        [[1]]
-
-        sage: K = Cone([(1,0),(0,1)])
-        sage: positive_operator_gens(K)
-        [
-        [1 0]  [0 1]  [0 0]  [0 0]
-        [0 0], [0 0], [1 0], [0 1]
-        ]
-
-    The trivial cone in a trivial space has no positive operators::
-
-        sage: K = Cone([], ToricLattice(0))
-        sage: positive_operator_gens(K)
-        []
-
-    Every operator is positive on the trivial cone::
-
-        sage: K = Cone([(0,)])
-        sage: positive_operator_gens(K)
-        [[1], [-1]]
-
-        sage: K = Cone([(0,0)])
-        sage: K.is_trivial()
-        True
-        sage: positive_operator_gens(K)
-        [
-        [1 0]  [-1  0]  [0 1]  [ 0 -1]  [0 0]  [ 0  0]  [0 0]  [ 0  0]
-        [0 0], [ 0  0], [0 0], [ 0  0], [1 0], [-1  0], [0 1], [ 0 -1]
-        ]
-
-    Every operator is positive on the ambient vector space::
+    Determine whether or not ``L`` is a Z-operator on ``K``.

 
 

-        sage: K = Cone([(1,),(-1,)])
-        sage: K.is_full_space()
-        True
-        sage: positive_operator_gens(K)
-        [[1], [-1]]
-
-        sage: K = Cone([(1,0),(-1,0),(0,1),(0,-1)])
-        sage: K.is_full_space()
-        True
-        sage: positive_operator_gens(K)
-        [
-        [1 0]  [-1  0]  [0 1]  [ 0 -1]  [0 0]  [ 0  0]  [0 0]  [ 0  0]
-        [0 0], [ 0  0], [0 0], [ 0  0], [1 0], [-1  0], [0 1], [ 0 -1]
-        ]
-
-    A non-obvious application is to find the positive operators on the
-    right half-plane::
-
-        sage: K = Cone([(1,0),(0,1),(0,-1)])
-        sage: positive_operator_gens(K)
-        [
-        [1 0]  [0 0]  [ 0  0]  [0 0]  [ 0  0]
-        [0 0], [1 0], [-1  0], [0 1], [ 0 -1]
-        ]
-
-    TESTS:
-
-    Each positive operator generator should send the generators of the
-    cone into the cone::
-
-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=4)
-        sage: pi_of_K = positive_operator_gens(K)
-        sage: all([ K.contains(P*x) for P in pi_of_K for x in K ])
-        True
-
-    Each positive operator generator should send a random element of the
-    cone into the cone::
+    We say that ``L`` is a Z-operator on ``K`` if `\left\langle
+    L\left\lparenx\right\rparen,s\right\rangle \le 0` for all pairs
+    `\left\langle x,s \right\rangle` in the complementarity set of
+    ``K``. It is known that this property need only be
+    checked for generators of ``K`` and its dual.

 
 

-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=4)
-        sage: pi_of_K = positive_operator_gens(K)
-        sage: all([ K.contains(P*K.random_element(QQ)) for P in pi_of_K ])
-        True
+    A matrix is a Z-operator on ``K`` if and only if its negation is a
+    cross-positive operator on ``K``.

 
 

-    A random element of the positive operator cone should send the
-    generators of the cone into the cone::
+    INPUT:

 
 

-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=4)
-        sage: pi_of_K = positive_operator_gens(K)
-        sage: L = ToricLattice(K.lattice_dim()**2)
-        sage: pi_cone = Cone([ g.list() for g in pi_of_K ],
-        ....:                lattice=L,
-        ....:                check=False)
-        sage: P = matrix(K.lattice_dim(), pi_cone.random_element(QQ).list())
-        sage: all([ K.contains(P*x) for x in K ])
-        True
+    - ``L`` -- A linear transformation or matrix.

 
 

-    A random element of the positive operator cone should send a random
-    element of the cone into the cone::
+    - ``K`` -- A polyhedral closed convex cone.

 
 

-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=4)
-        sage: pi_of_K = positive_operator_gens(K)
-        sage: L = ToricLattice(K.lattice_dim()**2)
-        sage: pi_cone = Cone([ g.list() for g in pi_of_K ],
-        ....:                lattice=L,
-        ....:                check=False)
-        sage: P = matrix(K.lattice_dim(), pi_cone.random_element(QQ).list())
-        sage: K.contains(P*K.random_element(ring=QQ))
-        True
+    OUTPUT:

 
 

-    The lineality space of the dual of the cone of positive operators
-    can be computed from the lineality spaces of the cone and its dual::
+    ``True`` if it can be proven that ``L`` is a Z-operator on ``K``,
+    and ``False`` otherwise.

 
 

-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=4)
-        sage: pi_of_K = positive_operator_gens(K)
-        sage: L = ToricLattice(K.lattice_dim()**2)
-        sage: pi_cone = Cone([ g.list() for g in pi_of_K ],
-        ....:                lattice=L,
-        ....:                check=False)
-        sage: actual = pi_cone.dual().linear_subspace()
-        sage: U1 = [ vector((s.tensor_product(x)).list())
-        ....:        for x in K.lines()
-        ....:        for s in K.dual() ]
-        sage: U2 = [ vector((s.tensor_product(x)).list())
-        ....:        for x in K
-        ....:        for s in K.dual().lines() ]
-        sage: expected = pi_cone.lattice().vector_space().span(U1 + U2)
-        sage: actual == expected
-        True
+    .. WARNING::

 
 

-    The lineality of the dual of the cone of positive operators
-    is known from its lineality space::
+        If this function returns ``True``, then ``L`` is a Z-operator
+        on ``K``. However, if ``False`` is returned, that could mean one
+        of two things. The first is that ``L`` is definitely not
+        a Z-operator on ``K``. The second is more of an "I don't know"
+        answer, returned (for example) if we cannot prove that an inner
+        product is nonnegative.

 
 

-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=4)
-        sage: n = K.lattice_dim()
-        sage: m = K.dim()
-        sage: l = K.lineality()
-        sage: pi_of_K = positive_operator_gens(K)
-        sage: L = ToricLattice(n**2)
-        sage: pi_cone = Cone([p.list() for p in pi_of_K],
-        ....:                 lattice=L,
-        ....:                 check=False)
-        sage: actual = pi_cone.dual().lineality()
-        sage: expected = l*(m - l) + m*(n - m)
-        sage: actual == expected
-        True
+    EXAMPLES:

 
 

-    The dimension of the cone of positive operators is given by the
-    corollary in my paper::
+    The identity is always a Z-operator in a nontrivial space::

 
         sage: set_random_seed()
 
         sage: set_random_seed()

-        sage: K = random_cone(max_ambient_dim=4)
-        sage: n = K.lattice_dim()
-        sage: m = K.dim()
-        sage: l = K.lineality()
-        sage: pi_of_K = positive_operator_gens(K)
-        sage: L = ToricLattice(n**2)
-        sage: pi_cone = Cone([p.list() for p in pi_of_K],
-        ....:                lattice=L,
-        ....:                check=False)
-        sage: actual = pi_cone.dim()
-        sage: expected = n**2 - l*(m - l) - (n - m)*m
-        sage: actual == expected
+        sage: K = random_cone(min_ambient_dim=1, max_ambient_dim=8)
+        sage: L = identity_matrix(K.lattice_dim())
+        sage: is_Z_on(L,K)

         True
 
         True
 

-    The trivial cone, full space, and half-plane all give rise to the
-    expected dimensions::
+    As is the "zero" transformation::

 
 

-        sage: n = ZZ.random_element().abs()
-        sage: K = Cone([[0] * n], ToricLattice(n))
-        sage: K.is_trivial()
-        True
-        sage: L = ToricLattice(n^2)
-        sage: pi_of_K = positive_operator_gens(K)
-        sage: pi_cone = Cone([p.list() for p in pi_of_K],
-        ....:                lattice=L,
-        ....:                check=False)
-        sage: actual = pi_cone.dim()
-        sage: actual == n^2
-        True
-        sage: K = K.dual()
-        sage: K.is_full_space()
-        True
-        sage: pi_of_K = positive_operator_gens(K)
-        sage: pi_cone = Cone([p.list() for p in pi_of_K],
-        ....:                lattice=L,
-        ....:                check=False)
-        sage: actual = pi_cone.dim()
-        sage: actual == n^2
-        True
-        sage: K = Cone([(1,0),(0,1),(0,-1)])
-        sage: pi_of_K = positive_operator_gens(K)
-        sage: actual = Cone([p.list() for p in pi_of_K], check=False).dim()
-        sage: actual == 3
+        sage: K = random_cone(min_ambient_dim=1, max_ambient_dim=8)
+        sage: R = K.lattice().vector_space().base_ring()
+        sage: L = zero_matrix(R, K.lattice_dim())
+        sage: is_Z_on(L,K)

         True
 
         True
 

-    The lineality of the cone of positive operators follows from the
-    description of its generators::
-
-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=4)
-        sage: n = K.lattice_dim()
-        sage: pi_of_K = positive_operator_gens(K)
-        sage: L = ToricLattice(n**2)
-        sage: pi_cone = Cone([p.list() for p in pi_of_K],
-        ....:                lattice=L,
-        ....:                check=False)
-        sage: actual = pi_cone.lineality()
-        sage: expected = n**2 - K.dim()*K.dual().dim()
-        sage: actual == expected
-        True
+    TESTS:

 
 

-    The trivial cone, full space, and half-plane all give rise to the
-    expected linealities::
+    Everything in ``K.Z_operators_gens()`` should be a Z-operator
+    on ``K``::

 
 

-        sage: n = ZZ.random_element().abs()
-        sage: K = Cone([[0] * n], ToricLattice(n))
-        sage: K.is_trivial()
-        True
-        sage: L = ToricLattice(n^2)
-        sage: pi_of_K = positive_operator_gens(K)
-        sage: pi_cone = Cone([p.list() for p in pi_of_K],
-        ....:                lattice=L,
-        ....:                check=False)
-        sage: actual = pi_cone.lineality()
-        sage: actual == n^2
-        True
-        sage: K = K.dual()
-        sage: K.is_full_space()
-        True
-        sage: pi_of_K = positive_operator_gens(K)
-        sage: pi_cone = Cone([p.list() for p in pi_of_K], lattice=L)
-        sage: pi_cone.lineality() == n^2
-        True
-        sage: K = Cone([(1,0),(0,1),(0,-1)])
-        sage: pi_of_K = positive_operator_gens(K)
-        sage: pi_cone = Cone([p.list() for p in pi_of_K], check=False)
-        sage: actual = pi_cone.lineality()
-        sage: actual == 2
+        sage: K = random_cone(min_ambient_dim=1, max_ambient_dim=6)
+        sage: all([ is_Z_on(L,K)
+        ....:       for L in K.Z_operators_gens() ])

         True
         True

-
-    A cone is proper if and only if its cone of positive operators
-    is proper::
-
-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=4)
-        sage: pi_of_K = positive_operator_gens(K)
-        sage: L = ToricLattice(K.lattice_dim()**2)
-        sage: pi_cone = Cone([p.list() for p in pi_of_K],
-        ....:                lattice=L,
-        ....:                check=False)
-        sage: K.is_proper() == pi_cone.is_proper()
+        sage: all([ is_Z_on(L.change_ring(SR),K)
+        ....:       for L in K.Z_operators_gens() ])

         True
 
         True
 

-    The positive operators of a permuted cone can be obtained by
-    conjugation::
-
-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=4)
-        sage: L = ToricLattice(K.lattice_dim()**2)
-        sage: p = SymmetricGroup(K.lattice_dim()).random_element().matrix()
-        sage: pK = Cone([ p*k for k in K ], K.lattice(), check=False)
-        sage: pi_of_pK = positive_operator_gens(pK)
-        sage: actual = Cone([t.list() for t in pi_of_pK],
-        ....:                lattice=L,
-        ....:                check=False)
-        sage: pi_of_K = positive_operator_gens(K)
-        sage: expected = Cone([(p*t*p.inverse()).list() for t in pi_of_K],
-        ....:                   lattice=L,
-        ....:                   check=False)
-        sage: actual.is_equivalent(expected)
-        True

     """
     """

-    # Matrices are not vectors in Sage, so we have to convert them
-    # to vectors explicitly before we can find a basis. We need these
-    # two values to construct the appropriate "long vector" space.
-    F = K.lattice().base_field()
-    n = K.lattice_dim()
-
-    tensor_products = [ s.tensor_product(x) for x in K for s in K.dual() ]
-
-    # Convert those tensor products to long vectors.
-    W = VectorSpace(F, n**2)
-    vectors = [ W(tp.list()) for tp in tensor_products ]
-
-    check = True
-    if K.is_solid() or K.is_strictly_convex():
-        # The lineality space of either ``K`` or ``K.dual()`` is
-        # trivial and it's easy to show that our generating set is
-        # minimal. I would love a proof that this works when ``K`` is
-        # neither pointed nor solid.
-        #
-        # Note that in that case we can get *duplicates*, since the
-        # tensor product of (x,s) is the same as that of (-x,-s).
-        check = False
-
-    # Create the dual cone of the positive operators, expressed as
-    # long vectors.
-    pi_dual = Cone(vectors, ToricLattice(W.dimension()), check=check)
-
-    # Now compute the desired cone from its dual...
-    pi_cone = pi_dual.dual()
-
-    # And finally convert its rays back to matrix representations.
-    M = MatrixSpace(F, n)
-    return [ M(v.list()) for v in pi_cone ]
-
-
-def Z_transformation_gens(K):
-    r"""
-    Compute generators of the cone of Z-transformations on this cone.
-
-    OUTPUT:
+    return is_cross_positive_on(-L,K)

 
 

-    A list of `n`-by-``n`` matrices where ``n == K.lattice_dim()``.
-    Each matrix ``L`` in the list should have the property that
-    ``(L*x).inner_product(s) <= 0`` whenever ``(x,s)`` is an element the
-    discrete complementarity set of ``K``. Moreover, any nonnegative
-    linear combination of these matrices shares the same property.

 
 

-    EXAMPLES:
+def is_lyapunov_like_on(L,K):
+    r"""
+    Determine whether or not ``L`` is Lyapunov-like on ``K``.

 
 

-    Z-transformations on the nonnegative orthant are just Z-matrices.
-    That is, matrices whose off-diagonal elements are nonnegative::
-
-        sage: K = Cone([(1,0),(0,1)])
-        sage: Z_transformation_gens(K)
-        [
-        [ 0 -1]  [ 0  0]  [-1  0]  [1 0]  [ 0  0]  [0 0]
-        [ 0  0], [-1  0], [ 0  0], [0 0], [ 0 -1], [0 1]
-        ]
-        sage: K = Cone([(1,0,0,0),(0,1,0,0),(0,0,1,0),(0,0,0,1)])
-        sage: all([ z[i][j] <= 0 for z in Z_transformation_gens(K)
-        ....:                    for i in range(z.nrows())
-        ....:                    for j in range(z.ncols())
-        ....:                    if i != j ])
-        True
+    We say that ``L`` is Lyapunov-like on ``K`` if `\left\langle
+    L\left\lparenx\right\rparen,s\right\rangle = 0` for all pairs
+    `\left\langle x,s \right\rangle` in the complementarity set of
+    ``K``. This property need only be checked for generators of
+    ``K`` and its dual.

 
 

-    The trivial cone in a trivial space has no Z-transformations::
+    INPUT:

 
 

-        sage: K = Cone([], ToricLattice(0))
-        sage: Z_transformation_gens(K)
-        []
+    - ``L`` -- A linear transformation or matrix.

 
 

-    Every operator is a Z-transformation on the ambient vector space::
+    - ``K`` -- A polyhedral closed convex cone.

 
 

-        sage: K = Cone([(1,),(-1,)])
-        sage: K.is_full_space()
-        True
-        sage: Z_transformation_gens(K)
-        [[-1], [1]]
+    OUTPUT:

 
 

-        sage: K = Cone([(1,0),(-1,0),(0,1),(0,-1)])
-        sage: K.is_full_space()
-        True
-        sage: Z_transformation_gens(K)
-        [
-        [-1  0]  [1 0]  [ 0 -1]  [0 1]  [ 0  0]  [0 0]  [ 0  0]  [0 0]
-        [ 0  0], [0 0], [ 0  0], [0 0], [-1  0], [1 0], [ 0 -1], [0 1]
-        ]
-
-    A non-obvious application is to find the Z-transformations on the
-    right half-plane::
-
-        sage: K = Cone([(1,0),(0,1),(0,-1)])
-        sage: Z_transformation_gens(K)
-        [
-        [-1  0]  [1 0]  [ 0  0]  [0 0]  [ 0  0]  [0 0]
-        [ 0  0], [0 0], [-1  0], [1 0], [ 0 -1], [0 1]
-        ]
-
-    Z-transformations on a subspace are Lyapunov-like and vice-versa::
-
-        sage: K = Cone([(1,0),(-1,0),(0,1),(0,-1)])
-        sage: K.is_full_space()
-        True
-        sage: lls = span([ vector(l.list()) for l in K.lyapunov_like_basis() ])
-        sage: zs  = span([ vector(z.list()) for z in Z_transformation_gens(K) ])
-        sage: zs == lls
-        True
+    ``True`` if it can be proven that ``L`` is Lyapunov-like on ``K``,
+    and ``False`` otherwise.

 
 

-    TESTS:
+    .. WARNING::

 
 

-    The Z-property is possessed by every Z-transformation::
+        If this function returns ``True``, then ``L`` is Lyapunov-like
+        on ``K``. However, if ``False`` is returned, that could mean one
+        of two things. The first is that ``L`` is definitely not
+        Lyapunov-like on ``K``. The second is more of an "I don't know"
+        answer, returned (for example) if we cannot prove that an inner
+        product is zero.

 
 

-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=4)
-        sage: Z_of_K = Z_transformation_gens(K)
-        sage: dcs = K.discrete_complementarity_set()
-        sage: all([(z*x).inner_product(s) <= 0 for z in Z_of_K
-        ....:                                  for (x,s) in dcs])
-        True
+    EXAMPLES:

 
 

-    The lineality space of the cone of Z-transformations is the space of
-    Lyapunov-like transformations::
+    The identity is always Lyapunov-like in a nontrivial space::

 
         sage: set_random_seed()
 
         sage: set_random_seed()

-        sage: K = random_cone(max_ambient_dim=4)
-        sage: L = ToricLattice(K.lattice_dim()**2)
-        sage: Z_cone = Cone([ z.list() for z in Z_transformation_gens(K) ],
-        ....:               lattice=L,
-        ....:               check=False)
-        sage: ll_basis = [ vector(l.list()) for l in K.lyapunov_like_basis() ]
-        sage: lls = L.vector_space().span(ll_basis)
-        sage: Z_cone.linear_subspace() == lls
+        sage: K = random_cone(min_ambient_dim=1, max_ambient_dim=8)
+        sage: L = identity_matrix(K.lattice_dim())
+        sage: is_lyapunov_like_on(L,K)

         True
 
         True
 

-    The lineality of the Z-transformations on a cone is the Lyapunov
-    rank of that cone::
+    As is the "zero" transformation::

 
 

-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=4)
-        sage: Z_of_K = Z_transformation_gens(K)
-        sage: L = ToricLattice(K.lattice_dim()**2)
-        sage: Z_cone  = Cone([ z.list() for z in Z_of_K ],
-        ....:                lattice=L,
-        ....:                check=False)
-        sage: Z_cone.lineality() == K.lyapunov_rank()
+        sage: K = random_cone(min_ambient_dim=1, max_ambient_dim=8)
+        sage: R = K.lattice().vector_space().base_ring()
+        sage: L = zero_matrix(R, K.lattice_dim())
+        sage: is_lyapunov_like_on(L,K)

         True
 
         True
 

-    The lineality spaces of the duals of the positive operator and
-    Z-transformation cones are equal. From this it follows that the
-    dimensions of the Z-transformation cone and positive operator cone
-    are equal::
-
-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=4)
-        sage: pi_of_K = positive_operator_gens(K)
-        sage: Z_of_K = Z_transformation_gens(K)
-        sage: L = ToricLattice(K.lattice_dim()**2)
-        sage: pi_cone = Cone([p.list() for p in pi_of_K],
-        ....:                lattice=L,
-        ....:                check=False)
-        sage: Z_cone = Cone([ z.list() for z in Z_of_K],
-        ....:               lattice=L,
-        ....:               check=False)
-        sage: pi_cone.dim() == Z_cone.dim()
-        True
-        sage: pi_star = pi_cone.dual()
-        sage: z_star = Z_cone.dual()
-        sage: pi_star.linear_subspace() == z_star.linear_subspace()
-        True
+    TESTS:

 
 

-    The trivial cone, full space, and half-plane all give rise to the
-    expected dimensions::
+    Everything in ``K.lyapunov_like_basis()`` should be Lyapunov-like
+    on ``K``::

 
 

-        sage: n = ZZ.random_element().abs()
-        sage: K = Cone([[0] * n], ToricLattice(n))
-        sage: K.is_trivial()
-        True
-        sage: L = ToricLattice(n^2)
-        sage: Z_of_K = Z_transformation_gens(K)
-        sage: Z_cone = Cone([z.list() for z in Z_of_K],
-        ....:               lattice=L,
-        ....:               check=False)
-        sage: actual = Z_cone.dim()
-        sage: actual == n^2
-        True
-        sage: K = K.dual()
-        sage: K.is_full_space()
-        True
-        sage: Z_of_K = Z_transformation_gens(K)
-        sage: Z_cone = Cone([z.list() for z in Z_of_K],
-        ....:                lattice=L,
-        ....:                check=False)
-        sage: actual = Z_cone.dim()
-        sage: actual == n^2
+        sage: K = random_cone(min_ambient_dim=1, max_ambient_dim=6)
+        sage: all([ is_lyapunov_like_on(L,K)
+        ....:       for L in K.lyapunov_like_basis() ])

         True
         True

-        sage: K = Cone([(1,0),(0,1),(0,-1)])
-        sage: Z_of_K = Z_transformation_gens(K)
-        sage: Z_cone = Cone([z.list() for z in Z_of_K], check=False)
-        sage: Z_cone.dim() == 3
+        sage: all([ is_lyapunov_like_on(L.change_ring(SR),K)
+        ....:       for L in K.lyapunov_like_basis() ])

         True
 
         True
 

-    The Z-transformations of a permuted cone can be obtained by
-    conjugation::
-
-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=4)
-        sage: L = ToricLattice(K.lattice_dim()**2)
-        sage: p = SymmetricGroup(K.lattice_dim()).random_element().matrix()
-        sage: pK = Cone([ p*k for k in K ], K.lattice(), check=False)
-        sage: Z_of_pK = Z_transformation_gens(pK)
-        sage: actual = Cone([t.list() for t in Z_of_pK],
-        ....:                lattice=L,
-        ....:                check=False)
-        sage: Z_of_K = Z_transformation_gens(K)
-        sage: expected = Cone([(p*t*p.inverse()).list() for t in Z_of_K],
-        ....:                   lattice=L,
-        ....:                   check=False)
-        sage: actual.is_equivalent(expected)
-        True

     """
     """

-    # Matrices are not vectors in Sage, so we have to convert them
-    # to vectors explicitly before we can find a basis. We need these
-    # two values to construct the appropriate "long vector" space.
-    F = K.lattice().base_field()
-    n = K.lattice_dim()
-
-    # These tensor products contain generators for the dual cone of
-    # the cross-positive transformations.
-    tensor_products = [ s.tensor_product(x)
-                        for (x,s) in K.discrete_complementarity_set() ]
-
-    # Turn our matrices into long vectors...
-    W = VectorSpace(F, n**2)
-    vectors = [ W(m.list()) for m in tensor_products ]
-
-    check = True
-    if K.is_solid() or K.is_strictly_convex():
-        # The lineality space of either ``K`` or ``K.dual()`` is
-        # trivial and it's easy to show that our generating set is
-        # minimal. I would love a proof that this works when ``K`` is
-        # neither pointed nor solid.
-        #
-        # Note that in that case we can get *duplicates*, since the
-        # tensor product of (x,s) is the same as that of (-x,-s).
-        check = False
-
-    # Create the dual cone of the cross-positive operators,
-    # expressed as long vectors.
-    Sigma_dual = Cone(vectors, lattice=ToricLattice(W.dimension()), check=check)
-
-    # Now compute the desired cone from its dual...
-    Sigma_cone = Sigma_dual.dual()
-
-    # And finally convert its rays back to matrix representations.
-    # But first, make them negative, so we get Z-transformations and
-    # not cross-positive ones.
-    M = MatrixSpace(F, n)
-    return [ -M(v.list()) for v in Sigma_cone ]
+    if L.base_ring().is_exact() or L.base_ring() is SR:
+        # The "fast method" of creating a vector space based on a
+        # ``lyapunov_like_basis`` is actually slower than this.
+        return all([ s*(L*x) == 0
+                     for (x,s) in K.discrete_complementarity_set() ])
+    else:
+        # The only inexact ring that we're willing to work with is SR,
+        # since it can still be exact when working with symbolic
+        # constants like pi and e.
+        raise ValueError('base ring of operator L is neither SR nor exact')
+
+def LL_cone(K):
+    gens = K.lyapunov_like_basis()
+    L = ToricLattice(K.lattice_dim()**2)
+    return Cone([ g.list() for g in gens ], lattice=L, check=False)

 
 

+def Sigma_cone(K):
+    gens = K.cross_positive_operators_gens()
+    L = ToricLattice(K.lattice_dim()**2)
+    return Cone([ g.list() for g in gens ], lattice=L, check=False)

 
 def Z_cone(K):
 
 def Z_cone(K):

-    gens = Z_transformation_gens(K)
+    gens = K.Z_operators_gens()

     L = ToricLattice(K.lattice_dim()**2)
     return Cone([ g.list() for g in gens ], lattice=L, check=False)
 
     L = ToricLattice(K.lattice_dim()**2)
     return Cone([ g.list() for g in gens ], lattice=L, check=False)
 

-def pi_cone(K):
-    gens = positive_operator_gens(K)
-    L = ToricLattice(K.lattice_dim()**2)
+def pi_cone(K1, K2=None):
+    if K2 is None:
+        K2 = K1
+    gens = K1.positive_operators_gens(K2)
+    L = ToricLattice(K1.lattice_dim()*K2.lattice_dim())

     return Cone([ g.list() for g in gens ], lattice=L, check=False)
     return Cone([ g.list() for g in gens ], lattice=L, check=False)