-# Sage doesn't load ~/.sage/init.sage during testing (sage -t), so we
-# have to explicitly mangle our sitedir here so that "mjo.cone"
-# resolves.
-from os.path import abspath
-from site import addsitedir
-addsitedir(abspath('../../'))
-
from sage.all import *
-
-def lyapunov_rank(K):
- r"""
- Compute the Lyapunov (or bilinearity) rank of this cone.
-
- The Lyapunov rank of a cone can be thought of in (mainly) two ways:
-
- 1. The dimension of the Lie algebra of the automorphism group of the
- cone.
-
- 2. The dimension of the linear space of all Lyapunov-like
- transformations on the cone.
-
- INPUT:
-
- A closed, convex polyhedral cone.
-
- OUTPUT:
-
- An integer representing the Lyapunov rank of the cone. If the
- dimension of the ambient vector space is `n`, then the Lyapunov rank
- will be between `1` and `n` inclusive; however a rank of `n-1` is
- not possible (see the first reference).
-
- .. note::
-
- In the references, the cones are always assumed to be proper. We
- do not impose this restriction.
-
- .. seealso::
-
- :meth:`is_proper`
-
- ALGORITHM:
-
- The codimension formula from the second reference is used. We find
- all pairs `(x,s)` in the complementarity set of `K` such that `x`
- and `s` are rays of our cone. It is known that these vectors are
- sufficient to apply the codimension formula. Once we have all such
- pairs, we "brute force" the codimension formula by finding all
- linearly-independent `xs^{T}`.
-
- REFERENCES:
-
- 1. M.S. Gowda and J. Tao. On the bilinearity rank of a proper cone
- and Lyapunov-like transformations, Mathematical Programming, 147
- (2014) 155-170.
-
- 2. G. Rudolf, N. Noyan, D. Papp, and F. Alizadeh, Bilinear
- optimality constraints for the cone of positive polynomials,
- Mathematical Programming, Series B, 129 (2011) 5-31.
-
- EXAMPLES:
-
- The nonnegative orthant in `\mathbb{R}^{n}` always has rank `n`::
-
- sage: positives = Cone([(1,)])
- sage: lyapunov_rank(positives)
- 1
- sage: quadrant = Cone([(1,0), (0,1)])
- sage: lyapunov_rank(quadrant)
- 2
- sage: octant = Cone([(1,0,0), (0,1,0), (0,0,1)])
- sage: lyapunov_rank(octant)
- 3
-
- The `L^{3}_{1}` cone is known to have a Lyapunov rank of one::
-
- sage: L31 = Cone([(1,0,1), (0,-1,1), (-1,0,1), (0,1,1)])
- sage: lyapunov_rank(L31)
- 1
-
- Likewise for the `L^{3}_{\infty}` cone::
-
- sage: L3infty = Cone([(0,1,1), (1,0,1), (0,-1,1), (-1,0,1)])
- sage: lyapunov_rank(L3infty)
- 1
-
- The Lyapunov rank should be additive on a product of cones::
-
- sage: L31 = Cone([(1,0,1), (0,-1,1), (-1,0,1), (0,1,1)])
- sage: octant = Cone([(1,0,0), (0,1,0), (0,0,1)])
- sage: K = L31.cartesian_product(octant)
- sage: lyapunov_rank(K) == lyapunov_rank(L31) + lyapunov_rank(octant)
- True
-
- Two isomorphic cones should have the same Lyapunov rank. The cone
- ``K`` in the following example is isomorphic to the nonnegative
- octant in `\mathbb{R}^{3}`::
-
- sage: K = Cone([(1,2,3), (-1,1,0), (1,0,6)])
- sage: lyapunov_rank(K)
- 3
-
- The dual cone `K^{*}` of ``K`` should have the same Lyapunov rank as ``K``
- itself::
-
- sage: K = Cone([(2,2,4), (-1,9,0), (2,0,6)])
- sage: lyapunov_rank(K) == lyapunov_rank(K.dual())
- True
-
- """
- V = K.lattice().vector_space()
-
- xs = [V(x) for x in K.rays()]
- ss = [V(s) for s in K.dual().rays()]
-
- # WARNING: This isn't really C(K), it only contains the pairs
- # (x,s) in C(K) where x,s are extreme in their respective cones.
- C_of_K = [(x,s) for x in xs for s in ss if x.inner_product(s) == 0]
-
- matrices = [x.column() * s.row() for (x,s) in C_of_K]
-
- # Sage doesn't think matrices are vectors, so we have to convert
- # our matrices to vectors explicitly before we can figure out how
- # many are linearly-indepenedent.
- #
- # The space W has the same base ring as V, but dimension
- # dim(V)^2. So it has the same dimension as the space of linear
- # transformations on V. In other words, it's just the right size
- # to create an isomorphism between it and our matrices.
- W = VectorSpace(V.base_ring(), V.dimension()**2)
-
- def phi(m):
- r"""
- Convert a matrix to a vector isomorphically.
- """
- return W(m.list())
-
- vectors = [phi(m) for m in matrices]
-
- return (W.dimension() - W.span(vectors).rank())
+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):
+ gens = K.Z_operators_gens()
+ L = ToricLattice(K.lattice_dim()**2)
+ return Cone(( g.list() for g in gens ), lattice=L, check=False)
+
+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)