]>
gitweb.michael.orlitzky.com - sage.d.git/blob - mjo/cone/cone.py
32e13861fdbb0cb2f806f5892cd9c90c1df051f4
1 # Sage doesn't load ~/.sage/init.sage during testing (sage -t), so we
2 # have to explicitly mangle our sitedir here so that "mjo.cone"
4 from os
.path
import abspath
5 from site
import addsitedir
6 addsitedir(abspath('../../'))
10 def is_lyapunov_like(L
,K
):
12 Determine whether or not ``L`` is Lyapunov-like on ``K``.
14 We say that ``L`` is Lyapunov-like on ``K`` if `\left\langle
15 L\left\lparenx\right\rparen,s\right\rangle = 0` for all pairs
16 `\left\langle x,s \right\rangle` in the complementarity set of
17 ``K``. It is known [Orlitzky]_ that this property need only be
18 checked for generators of ``K`` and its dual.
22 - ``L`` -- A linear transformation or matrix.
24 - ``K`` -- A polyhedral closed convex cone.
28 ``True`` if it can be proven that ``L`` is Lyapunov-like on ``K``,
29 and ``False`` otherwise.
33 If this function returns ``True``, then ``L`` is Lyapunov-like
34 on ``K``. However, if ``False`` is returned, that could mean one
35 of two things. The first is that ``L`` is definitely not
36 Lyapunov-like on ``K``. The second is more of an "I don't know"
37 answer, returned (for example) if we cannot prove that an inner
42 M. Orlitzky. The Lyapunov rank of an improper cone.
43 http://www.optimization-online.org/DB_HTML/2015/10/5135.html
47 The identity is always Lyapunov-like in a nontrivial space::
49 sage: set_random_seed()
50 sage: K = random_cone(min_ambient_dim = 1, max_rays = 8)
51 sage: L = identity_matrix(K.lattice_dim())
52 sage: is_lyapunov_like(L,K)
55 As is the "zero" transformation::
57 sage: K = random_cone(min_ambient_dim = 1, max_rays = 5)
58 sage: R = K.lattice().vector_space().base_ring()
59 sage: L = zero_matrix(R, K.lattice_dim())
60 sage: is_lyapunov_like(L,K)
63 Everything in ``K.lyapunov_like_basis()`` should be Lyapunov-like
66 sage: K = random_cone(min_ambient_dim = 1, max_rays = 5)
67 sage: all([ is_lyapunov_like(L,K) for L in K.lyapunov_like_basis() ])
71 return all([(L
*x
).inner_product(s
) == 0
72 for (x
,s
) in K
.discrete_complementarity_set()])
75 def random_element(K
):
77 Return a random element of ``K`` from its ambient vector space.
81 The cone ``K`` is specified in terms of its generators, so that
82 ``K`` is equal to the convex conic combination of those generators.
83 To choose a random element of ``K``, we assign random nonnegative
84 coefficients to each generator of ``K`` and construct a new vector
87 A vector, rather than a ray, is returned so that the element may
88 have non-integer coordinates. Thus the element may have an
89 arbitrarily small norm.
93 A random element of the trivial cone is zero::
95 sage: set_random_seed()
96 sage: K = Cone([], ToricLattice(0))
97 sage: random_element(K)
99 sage: K = Cone([(0,)])
100 sage: random_element(K)
102 sage: K = Cone([(0,0)])
103 sage: random_element(K)
105 sage: K = Cone([(0,0,0)])
106 sage: random_element(K)
111 Any cone should contain an element of itself::
113 sage: set_random_seed()
114 sage: K = random_cone(max_rays = 8)
115 sage: K.contains(random_element(K))
119 V
= K
.lattice().vector_space()
121 coefficients
= [ F
.random_element().abs() for i
in range(K
.nrays()) ]
122 vector_gens
= map(V
, K
.rays())
123 scaled_gens
= [ coefficients
[i
]*vector_gens
[i
]
124 for i
in range(len(vector_gens
)) ]
126 # Make sure we return a vector. Without the coercion, we might
127 # return ``0`` when ``K`` has no rays.
128 v
= V(sum(scaled_gens
))
132 def positive_operators(K
):
134 Compute generators of the cone of positive operators on this cone.
138 A list of `n`-by-``n`` matrices where ``n == K.lattice_dim()``.
139 Each matrix ``P`` in the list should have the property that ``P*x``
140 is an element of ``K`` whenever ``x`` is an element of
141 ``K``. Moreover, any nonnegative linear combination of these
142 matrices shares the same property.
146 The trivial cone in a trivial space has no positive operators::
148 sage: K = Cone([], ToricLattice(0))
149 sage: positive_operators(K)
152 Positive operators on the nonnegative orthant are nonnegative matrices::
154 sage: K = Cone([(1,)])
155 sage: positive_operators(K)
158 sage: K = Cone([(1,0),(0,1)])
159 sage: positive_operators(K)
161 [1 0] [0 1] [0 0] [0 0]
162 [0 0], [0 0], [1 0], [0 1]
165 Every operator is positive on the ambient vector space::
167 sage: K = Cone([(1,),(-1,)])
168 sage: K.is_full_space()
170 sage: positive_operators(K)
173 sage: K = Cone([(1,0),(-1,0),(0,1),(0,-1)])
174 sage: K.is_full_space()
176 sage: positive_operators(K)
178 [1 0] [-1 0] [0 1] [ 0 -1] [0 0] [ 0 0] [0 0] [ 0 0]
179 [0 0], [ 0 0], [0 0], [ 0 0], [1 0], [-1 0], [0 1], [ 0 -1]
184 A positive operator on a cone should send its generators into the cone::
186 sage: K = random_cone(max_ambient_dim = 6)
187 sage: pi_of_K = positive_operators(K)
188 sage: all([K.contains(p*x) for p in pi_of_K for x in K.rays()])
192 # Sage doesn't think matrices are vectors, so we have to convert
193 # our matrices to vectors explicitly before we can figure out how
194 # many are linearly-indepenedent.
196 # The space W has the same base ring as V, but dimension
197 # dim(V)^2. So it has the same dimension as the space of linear
198 # transformations on V. In other words, it's just the right size
199 # to create an isomorphism between it and our matrices.
200 V
= K
.lattice().vector_space()
201 W
= VectorSpace(V
.base_ring(), V
.dimension()**2)
203 tensor_products
= [ s
.tensor_product(x
) for x
in K
for s
in K
.dual() ]
205 # Turn our matrices into long vectors...
206 vectors
= [ W(m
.list()) for m
in tensor_products
]
208 # Create the *dual* cone of the positive operators, expressed as
210 L
= ToricLattice(W
.dimension())
211 pi_dual
= Cone(vectors
, lattice
=L
)
213 # Now compute the desired cone from its dual...
214 pi_cone
= pi_dual
.dual()
216 # And finally convert its rays back to matrix representations.
217 M
= MatrixSpace(V
.base_ring(), V
.dimension())
219 return [ M(v
.list()) for v
in pi_cone
.rays() ]
222 def Z_transformations(K
):
224 Compute generators of the cone of Z-transformations on this cone.
228 A list of `n`-by-``n`` matrices where ``n == K.lattice_dim()``.
229 Each matrix ``L`` in the list should have the property that
230 ``(L*x).inner_product(s) <= 0`` whenever ``(x,s)`` is an element the
231 discrete complementarity set of ``K``. Moreover, any nonnegative
232 linear combination of these matrices shares the same property.
236 Z-transformations on the nonnegative orthant are just Z-matrices.
237 That is, matrices whose off-diagonal elements are nonnegative::
239 sage: K = Cone([(1,0),(0,1)])
240 sage: Z_transformations(K)
242 [ 0 -1] [ 0 0] [-1 0] [1 0] [ 0 0] [0 0]
243 [ 0 0], [-1 0], [ 0 0], [0 0], [ 0 -1], [0 1]
245 sage: K = Cone([(1,0,0,0),(0,1,0,0),(0,0,1,0),(0,0,0,1)])
246 sage: all([ z[i][j] <= 0 for z in Z_transformations(K)
247 ....: for i in range(z.nrows())
248 ....: for j in range(z.ncols())
252 The trivial cone in a trivial space has no Z-transformations::
254 sage: K = Cone([], ToricLattice(0))
255 sage: Z_transformations(K)
258 Z-transformations on a subspace are Lyapunov-like and vice-versa::
260 sage: K = Cone([(1,0),(-1,0),(0,1),(0,-1)])
261 sage: K.is_full_space()
263 sage: lls = span([ vector(l.list()) for l in K.lyapunov_like_basis() ])
264 sage: zs = span([ vector(z.list()) for z in Z_transformations(K) ])
270 The Z-property is possessed by every Z-transformation::
272 sage: set_random_seed()
273 sage: K = random_cone(max_ambient_dim = 6)
274 sage: Z_of_K = Z_transformations(K)
275 sage: dcs = K.discrete_complementarity_set()
276 sage: all([(z*x).inner_product(s) <= 0 for z in Z_of_K
277 ....: for (x,s) in dcs])
280 The lineality space of Z is LL::
282 sage: set_random_seed()
283 sage: K = random_cone(min_ambient_dim = 1, max_ambient_dim = 6)
284 sage: lls = span([ vector(l.list()) for l in K.lyapunov_like_basis() ])
285 sage: z_cone = Cone([ z.list() for z in Z_transformations(K) ])
286 sage: z_cone.linear_subspace() == lls
290 # Sage doesn't think matrices are vectors, so we have to convert
291 # our matrices to vectors explicitly before we can figure out how
292 # many are linearly-indepenedent.
294 # The space W has the same base ring as V, but dimension
295 # dim(V)^2. So it has the same dimension as the space of linear
296 # transformations on V. In other words, it's just the right size
297 # to create an isomorphism between it and our matrices.
298 V
= K
.lattice().vector_space()
299 W
= VectorSpace(V
.base_ring(), V
.dimension()**2)
301 C_of_K
= K
.discrete_complementarity_set()
302 tensor_products
= [ s
.tensor_product(x
) for (x
,s
) in C_of_K
]
304 # Turn our matrices into long vectors...
305 vectors
= [ W(m
.list()) for m
in tensor_products
]
307 # Create the *dual* cone of the cross-positive operators,
308 # expressed as long vectors..
309 L
= ToricLattice(W
.dimension())
310 Sigma_dual
= Cone(vectors
, lattice
=L
)
312 # Now compute the desired cone from its dual...
313 Sigma_cone
= Sigma_dual
.dual()
315 # And finally convert its rays back to matrix representations.
316 # But first, make them negative, so we get Z-transformations and
317 # not cross-positive ones.
318 M
= MatrixSpace(V
.base_ring(), V
.dimension())
320 return [ -M(v
.list()) for v
in Sigma_cone
.rays() ]