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 124bf5415f8841d36343e11c78ed5cd61b7ed276..a7c16a520cbe19242398e176635ad3250e5c7671 100644 (file)
--- a/mjo/cone/cone.py
+++ b/mjo/cone/cone.py
@@ -1,20 +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``.
-
-    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.
+    Determine whether or not ``L`` is positive on ``K``.

 
 

-    There are faster ways of checking this property. For example, we
-    could compute a `lyapunov_like_basis` of the cone, and then test
-    whether or not the given matrix is contained in the span of that
-    basis. The value of this function is that it works on symbolic
-    matrices.
+    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:
 


@@ -24,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::


@@ -56,757 +43,249 @@ 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()])
-
-
-def positive_operator_gens(K1, K2 = None):
+    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"""
     r"""

-    Compute generators of the cone of positive operators on this cone. A
-    linear operator on a cone is positive if the image of the cone under
-    the operator is a subset of the cone. This concept can be extended
-    to two cones, where the image of the first cone under a positive
-    operator is a subset of the second cone.
+    Determine whether or not ``L`` is cross-positive on ``K``.

 
 

-    INPUT:
-
-    - ``K2`` -- (default: ``K1``) the codomain cone; the image of this
-                cone under the returned operators is a subset of ``K2``.
+    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.

 
 

-    OUTPUT:
+    INPUT:

 
 

-    A list of `m`-by-``n`` matrices where ``m == K2.lattice_dim()`` and
-    ``n == K1.lattice_dim()``. Each matrix ``P`` in the list should have
-    the property that ``P*x`` is an element of ``K2`` whenever ``x`` is
-    an element of ``K1``. Moreover, any nonnegative linear combination of
-    these matrices shares the same property.
+    - ``L`` -- A linear transformation or matrix.

 
 

-    .. SEEALSO::
+    - ``K`` -- A polyhedral closed convex cone.

 
 

-           :meth:`cross_positive_operator_gens`, :meth:`Z_operator_gens`,
+    OUTPUT:

 
 

-    REFERENCES:
+    ``True`` if it can be proven that ``L`` is cross-positive on ``K``,
+    and ``False`` otherwise.

 
 

-    .. [Orlitzky-Pi-Z]
-       M. Orlitzky.
-       Positive and Z-operators on closed convex cones.
+    .. WARNING::

 
 

-    .. [Tam]
-       B.-S. Tam.
-       Some results of polyhedral cones and simplicial cones.
-       Linear and Multilinear Algebra, 4:4 (1977) 281--284.
+        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:
 

-    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::
+    The identity is always cross-positive in a nontrivial space::

 
 

-        sage: K = Cone([(0,)])
-        sage: positive_operator_gens(K)
-        [[1], [-1]]
-
-        sage: K = Cone([(0,0)])
-        sage: K.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

-        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::

 
 

-        sage: K = Cone([(1,),(-1,)])
-        sage: K.is_full_space()
-        True
-        sage: positive_operator_gens(K)
-        [[1], [-1]]
+    As is the "zero" transformation::

 
 

-        sage: K = Cone([(1,0),(-1,0),(0,1),(0,-1)])
-        sage: K.is_full_space()
+        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
         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:
 
 
     TESTS:
 

-    Each positive operator generator should send the generators of one
-    cone into the other cone::
-
-        sage: set_random_seed()
-        sage: K1 = random_cone(max_ambient_dim=4)
-        sage: K2 = random_cone(max_ambient_dim=4)
-        sage: pi_K1_K2 = positive_operator_gens(K1,K2)
-        sage: all([ K2.contains(P*x) for P in pi_K1_K2 for x in K1 ])
-        True
-
-    Each positive operator generator should send a random element of one
-    cone into the other cone::
+    Everything in ``K.cross_positive_operators_gens()`` should be
+    cross-positive on ``K``::

 
 

-        sage: set_random_seed()
-        sage: K1 = random_cone(max_ambient_dim=4)
-        sage: K2 = random_cone(max_ambient_dim=4)
-        sage: pi_K1_K2 = positive_operator_gens(K1,K2)
-        sage: all([ K2.contains(P*K1.random_element(QQ)) for P in pi_K1_K2 ])
+        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

-
-    A random element of the positive operator cone should send the
-    generators of one cone into the other cone::
-
-        sage: set_random_seed()
-        sage: K1 = random_cone(max_ambient_dim=4)
-        sage: K2 = random_cone(max_ambient_dim=4)
-        sage: pi_K1_K2 = positive_operator_gens(K1,K2)
-        sage: L = ToricLattice(K1.lattice_dim() * K2.lattice_dim())
-        sage: pi_cone = Cone([ g.list() for g in pi_K1_K2 ],
-        ....:                lattice=L,
-        ....:                check=False)
-        sage: P = matrix(K2.lattice_dim(),
-        ....:            K1.lattice_dim(),
-        ....:            pi_cone.random_element(QQ).list())
-        sage: all([ K2.contains(P*x) for x in K1 ])
+        sage: all([ is_cross_positive_on(L.change_ring(SR),K)
+        ....:       for L in K.cross_positive_operators_gens() ])

         True
 
         True
 

-    A random element of the positive operator cone should send a random
-    element of one cone into the other cone::
-
-        sage: set_random_seed()
-        sage: K1 = random_cone(max_ambient_dim=4)
-        sage: K2 = random_cone(max_ambient_dim=4)
-        sage: pi_K1_K2 = positive_operator_gens(K1,K2)
-        sage: L = ToricLattice(K1.lattice_dim() * K2.lattice_dim())
-        sage: pi_cone = Cone([ g.list() for g in pi_K1_K2 ],
-        ....:                lattice=L,
-        ....:                check=False)
-        sage: P = matrix(K2.lattice_dim(),
-        ....:            K1.lattice_dim(),
-        ....:            pi_cone.random_element(QQ).list())
-        sage: K2.contains(P*K1.random_element(ring=QQ))
-        True
+    """
+    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')

 
 

-    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::

 
 

-        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
+def is_Z_on(L,K):
+    r"""
+    Determine whether or not ``L`` is a Z-operator on ``K``.

 
 

-    The lineality of the dual of the cone of positive operators
-    is known from its lineality space::
+    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: 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
+    A matrix is a Z-operator on ``K`` if and only if its negation is a
+    cross-positive operator on ``K``.

 
 

-    The dimension of the cone of positive operators is given by the
-    corollary in my paper::
+    INPUT:

 
 

-        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
-        True
+    - ``L`` -- A linear transformation or matrix.

 
 

-    The trivial cone, full space, and half-plane all give rise to the
-    expected dimensions::
+    - ``K`` -- A polyhedral closed convex cone.

 
 

-        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
-        True
+    OUTPUT:

 
 

-    The lineality of the cone of positive operators follows from the
-    description of its generators::
+    ``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: 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
+    .. WARNING::

 
 

-    The trivial cone, full space, and half-plane all give rise to the
-    expected linealities::
+        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: 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
-        True
+    EXAMPLES:

 
 

-    A cone is proper if and only if its cone of positive operators
-    is proper::
+    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: 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: 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 positive operators of a permuted cone can be obtained by
-    conjugation::
+    As is the "zero" transformation::

 
 

-        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)
+        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
 

-    A transformation is positive on a cone if and only if its adjoint is
-    positive on the dual of that cone::
-
-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=4)
-        sage: F = K.lattice().vector_space().base_field()
-        sage: n = K.lattice_dim()
-        sage: L = ToricLattice(n**2)
-        sage: W = VectorSpace(F, n**2)
-        sage: pi_of_K = positive_operator_gens(K)
-        sage: pi_of_K_star = positive_operator_gens(K.dual())
-        sage: pi_cone = Cone([p.list() for p in pi_of_K],
-        ....:                lattice=L,
-        ....:                check=False)
-        sage: pi_star = Cone([p.list() for p in pi_of_K_star],
-        ....:                lattice=L,
-        ....:                check=False)
-        sage: M = MatrixSpace(F, n)
-        sage: L = M(pi_cone.random_element(ring=QQ).list())
-        sage: pi_star.contains(W(L.transpose().list()))
-        True
+    TESTS:

 
 

-        sage: L = W.random_element()
-        sage: L_star = W(M(L.list()).transpose().list())
-        sage: pi_cone.contains(L) ==  pi_star.contains(L_star)
-        True
+    Everything in ``K.Z_operators_gens()`` should be a Z-operator
+    on ``K``::

 
 

-    The Lyapunov rank of the positive operator cone is the product of
-    the Lyapunov ranks of the associated cones if they're all proper::
-
-        sage: K1 = random_cone(max_ambient_dim=4,
-        ....:                  strictly_convex=True,
-        ....:                  solid=True)
-        sage: K2 = random_cone(max_ambient_dim=4,
-        ....:                  strictly_convex=True,
-        ....:                  solid=True)
-        sage: pi_K1_K2 = positive_operator_gens(K1,K2)
-        sage: L = ToricLattice(K1.lattice_dim() * K2.lattice_dim())
-        sage: pi_cone = Cone([ g.list() for g in pi_K1_K2 ],
-        ....:                lattice=L,
-        ....:                check=False)
-        sage: beta1 = K1.lyapunov_rank()
-        sage: beta2 = K2.lyapunov_rank()
-        sage: pi_cone.lyapunov_rank() == beta1*beta2
+        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

-
-    The Lyapunov-like operators on a proper polyhedral positive operator
-    cone can be computed from the Lyapunov-like operators on the cones
-    with respect to which the operators are positive::
-
-        sage: K1 = random_cone(max_ambient_dim=4,
-        ....:                  strictly_convex=True,
-        ....:                  solid=True)
-        sage: K2 = random_cone(max_ambient_dim=4,
-        ....:                  strictly_convex=True,
-        ....:                  solid=True)
-        sage: pi_K1_K2 = positive_operator_gens(K1,K2)
-        sage: F = K1.lattice().base_field()
-        sage: m = K1.lattice_dim()
-        sage: n = K2.lattice_dim()
-        sage: L = ToricLattice(m*n)
-        sage: M1 = MatrixSpace(F, m, m)
-        sage: M2 = MatrixSpace(F, n, n)
-        sage: LL_K1 = [ M1(x.list()) for x in K1.dual().lyapunov_like_basis() ]
-        sage: LL_K2 = [ M2(x.list()) for x in K2.lyapunov_like_basis() ]
-        sage: tps = [ s.tensor_product(x) for x in LL_K1 for s in LL_K2 ]
-        sage: W = VectorSpace(F, (m**2)*(n**2))
-        sage: expected = span(F, [ W(x.list()) for x in tps ])
-        sage: pi_cone = Cone([p.list() for p in pi_K1_K2],
-        ....:                 lattice=L,
-        ....:                 check=False)
-        sage: LL_pi = pi_cone.lyapunov_like_basis()
-        sage: actual = span(F, [ W(x.list()) for x in LL_pi ])
-        sage: actual == expected
+        sage: all([ is_Z_on(L.change_ring(SR),K)
+        ....:       for L in K.Z_operators_gens() ])

         True
 
     """
         True
 
     """

-    if K2 is None:
-        K2 = K1
-
-    # 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 = K1.lattice().base_field()
-    n = K1.lattice_dim()
-    m = K2.lattice_dim()
-
-    tensor_products = [ s.tensor_product(x) for x in K1 for s in K2.dual() ]
+    return is_cross_positive_on(-L,K)

 
 

-    # Convert those tensor products to long vectors.
-    W = VectorSpace(F, n*m)
-    vectors = [ W(tp.list()) for tp in tensor_products ]

 
 

-    check = True
-    if K1.is_proper() and K2.is_proper():
-        # All of the generators involved are extreme vectors and
-        # therefore minimal [Tam]_. If this cone is neither solid nor
-        # strictly convex, then the tensor product of ``s`` and ``x``
-        # is the same as that of ``-s`` and ``-x``. However, as a
-        # /set/, ``tensor_products`` may still be minimal.
-        check = False
-
-    # Create the dual cone of the positive operators, expressed as
-    # long vectors.
-    pi_dual = Cone(vectors, ToricLattice(W.dimension()), check=check)
+def is_lyapunov_like_on(L,K):
+    r"""
+    Determine whether or not ``L`` is Lyapunov-like on ``K``.

 
 

-    # Now compute the desired cone from its dual...
-    pi_cone = pi_dual.dual()
+    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.

 
 

-    # And finally convert its rays back to matrix representations.
-    M = MatrixSpace(F, m, n)
-    return [ M(v.list()) for v in pi_cone ]
+    INPUT:

 
 

+    - ``L`` -- A linear transformation or matrix.

 
 

-def cross_positive_operator_gens(K):
-    r"""
-    Compute generators of the cone of cross-positive operators on this
-    cone.
+    - ``K`` -- A polyhedral closed convex cone.

 
     OUTPUT:
 
 
     OUTPUT:
 

-    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 of
-    this cone's :meth:`discrete_complementarity_set`. Moreover, any
-    conic (nonnegative linear) combination of these matrices shares the
-    same property.
-
-    .. SEEALSO::
-
-       :meth:`positive_operator_gens`, :meth:`Z_operator_gens`,
+    ``True`` if it can be proven that ``L`` is Lyapunov-like on ``K``,
+    and ``False`` otherwise.

 
 

-    REFERENCES:
+    .. WARNING::

 
 

-    M. Orlitzky.
-    Positive and Z-operators on closed convex cones.
+        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.

 
     EXAMPLES:
 
 
     EXAMPLES:
 

-    Cross-positive operators on the nonnegative orthant are negations
-    of Z-matrices; that is, matrices whose off-diagonal elements are
-    nonnegative::
-
-        sage: K = Cone([(1,0),(0,1)])
-        sage: cross_positive_operator_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([ c[i][j] >= 0 for c in cross_positive_operator_gens(K)
-        ....:                    for i in range(c.nrows())
-        ....:                    for j in range(c.ncols())
-        ....:                    if i != j ])
-        True
-
-    The trivial cone in a trivial space has no cross-positive operators::
-
-        sage: K = Cone([], ToricLattice(0))
-        sage: cross_positive_operator_gens(K)
-        []
-
-    Every operator is a cross-positive operator on the ambient vector
-    space::
-
-        sage: K = Cone([(1,),(-1,)])
-        sage: K.is_full_space()
-        True
-        sage: cross_positive_operator_gens(K)
-        [[1], [-1]]
-
-        sage: K = Cone([(1,0),(-1,0),(0,1),(0,-1)])
-        sage: K.is_full_space()
-        True
-        sage: cross_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 cross-positive operators
-    on the right half-plane::
-
-        sage: K = Cone([(1,0),(0,1),(0,-1)])
-        sage: cross_positive_operator_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]
-        ]
-
-    Cross-positive operators 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: cs  = span([ vector(c.list()) for c in cross_positive_operator_gens(K) ])
-        sage: cs == lls
-        True
-
-    TESTS:
-
-    The cross-positive property is possessed by every cross-positive
-    operator::
-
-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=4)
-        sage: Sigma_of_K = cross_positive_operator_gens(K)
-        sage: dcs = K.discrete_complementarity_set()
-        sage: all([(c*x).inner_product(s) >= 0 for c in Sigma_of_K
-        ....:                                  for (x,s) in dcs])
-        True
-
-    The lineality space of the cone of cross-positive operators is the
-    space of Lyapunov-like operators::
-
-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=4)
-        sage: L = ToricLattice(K.lattice_dim()**2)
-        sage: Sigma_cone = Cone([ c.list() for c in cross_positive_operator_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: Sigma_cone.linear_subspace() == lls
-        True
-
-    The lineality of the cross-positive operators on a cone is the
-    Lyapunov rank of that cone::
-
-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=4)
-        sage: Sigma_of_K = cross_positive_operator_gens(K)
-        sage: L = ToricLattice(K.lattice_dim()**2)
-        sage: Sigma_cone  = Cone([ c.list() for c in Sigma_of_K ],
-        ....:                      lattice=L,
-        ....:                      check=False)
-        sage: Sigma_cone.lineality() == K.lyapunov_rank()
-        True
-
-    The lineality spaces of the duals of the positive and cross-positive
-    operator cones are equal. From this it follows that the dimensions of
-    the cross-positive operator cone and positive operator cone are equal::
+    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: pi_of_K = positive_operator_gens(K)
-        sage: Sigma_of_K = cross_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: Sigma_cone = Cone([ c.list() for c in Sigma_of_K],
-        ....:                     lattice=L,
-        ....:                     check=False)
-        sage: pi_cone.dim() == Sigma_cone.dim()
-        True
-        sage: pi_star = pi_cone.dual()
-        sage: sigma_star = Sigma_cone.dual()
-        sage: pi_star.linear_subspace() == sigma_star.linear_subspace()
+        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 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: Sigma_of_K = cross_positive_operator_gens(K)
-        sage: Sigma_cone = Cone([c.list() for c in Sigma_of_K],
-        ....:                    lattice=L,
-        ....:                    check=False)
-        sage: actual = Sigma_cone.dim()
-        sage: actual == n^2
-        True
-        sage: K = K.dual()
-        sage: K.is_full_space()
-        True
-        sage: Sigma_of_K = cross_positive_operator_gens(K)
-        sage: Sigma_cone = Cone([ c.list() for c in Sigma_of_K ],
-        ....:                     lattice=L,
-        ....:                     check=False)
-        sage: actual = Sigma_cone.dim()
-        sage: actual == n^2
-        True
-        sage: K = Cone([(1,0),(0,1),(0,-1)])
-        sage: Sigma_of_K = cross_positive_operator_gens(K)
-        sage: Sigma_cone = Cone([ c.list() for c in Sigma_of_K ], check=False)
-        sage: Sigma_cone.dim() == 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_lyapunov_like_on(L,K)

         True
 
         True
 

-    The cross-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: Sigma_of_pK = cross_positive_operator_gens(pK)
-        sage: actual = Cone([t.list() for t in Sigma_of_pK],
-        ....:                lattice=L,
-        ....:                check=False)
-        sage: Sigma_of_K = cross_positive_operator_gens(K)
-        sage: expected = Cone([ (p*t*p.inverse()).list() for t in Sigma_of_K ],
-        ....:                   lattice=L,
-        ....:                   check=False)
-        sage: actual.is_equivalent(expected)
-        True
+    TESTS:

 
 

-    An operator is cross-positive on a cone if and only if its
-    adjoint is cross-positive on the dual of that cone::
+    Everything in ``K.lyapunov_like_basis()`` should be Lyapunov-like
+    on ``K``::

 
 

-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=4)
-        sage: F = K.lattice().vector_space().base_field()
-        sage: n = K.lattice_dim()
-        sage: L = ToricLattice(n**2)
-        sage: W = VectorSpace(F, n**2)
-        sage: Sigma_of_K = cross_positive_operator_gens(K)
-        sage: Sigma_of_K_star = cross_positive_operator_gens(K.dual())
-        sage: Sigma_cone = Cone([ p.list() for p in Sigma_of_K ],
-        ....:                     lattice=L,
-        ....:                     check=False)
-        sage: Sigma_star = Cone([ p.list() for p in Sigma_of_K_star ],
-        ....:                     lattice=L,
-        ....:                     check=False)
-        sage: M = MatrixSpace(F, n)
-        sage: L = M(Sigma_cone.random_element(ring=QQ).list())
-        sage: Sigma_star.contains(W(L.transpose().list()))
+        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: L = W.random_element()
-        sage: L_star = W(M(L.list()).transpose().list())
-        sage: Sigma_cone.contains(L) ==  Sigma_star.contains(L_star)
+        sage: all([ is_lyapunov_like_on(L.change_ring(SR),K)
+        ....:       for L in K.lyapunov_like_basis() ])

         True
         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 operators.
-    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_proper():
-        # All of the generators involved are extreme vectors and
-        # therefore minimal. If this cone is neither solid nor
-        # strictly convex, then the tensor product of ``s`` and ``x``
-        # is the same as that of ``-s`` and ``-x``. However, as a
-        # /set/, ``tensor_products`` may still be minimal.
-        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.
-    M = MatrixSpace(F, n)
-    return [ M(v.list()) for v in Sigma_cone ]
-
-
-def Z_operator_gens(K):
-    r"""
-    Compute generators of the cone of Z-operators on this cone.
-
-    The Z-operators on a cone generalize the Z-matrices over the
-    nonnegative orthant. They are simply negations of the
-    :meth:`cross_positive_operators`.
-
-    OUTPUT:
-
-    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 of
-    this cone's :meth:`discrete_complementarity_set`. Moreover, any
-    conic (nonnegative linear) combination of these matrices shares the
-    same property.
-
-    .. SEEALSO::

 
 

-       :meth:`positive_operator_gens`, :meth:`cross_positive_operator_gens`,
-
-    REFERENCES:
-
-    M. Orlitzky.
-    Positive and Z-operators on closed convex cones.
-
-    TESTS:
-
-    The Z-property is possessed by every Z-operator::
-
-        sage: set_random_seed()
-        sage: K = random_cone(max_ambient_dim=4)
-        sage: Z_of_K = Z_operator_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

     """
     """

-    return [ -cp for cp in cross_positive_operator_gens(K) ]
-
+    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()
 
 def LL_cone(K):
     gens = K.lyapunov_like_basis()


@@ -814,16 +293,18 @@ def LL_cone(K):
     return Cone([ g.list() for g in gens ], lattice=L, check=False)
 
 def Sigma_cone(K):
     return Cone([ g.list() for g in gens ], lattice=L, check=False)
 
 def Sigma_cone(K):

-    gens = cross_positive_operator_gens(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):
     L = ToricLattice(K.lattice_dim()**2)
     return Cone([ g.list() for g in gens ], lattice=L, check=False)
 
 def Z_cone(K):

-    gens = Z_operator_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)