]>
gitweb.michael.orlitzky.com - sage.d.git/blob - mjo/cone/cone.py
f5371d6d803cdbc97fe4ba41cb26f0ee691f1a2a
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('../../'))
11 def drop_dependent(vs
):
13 Return the largest linearly-independent subset of ``vs``.
16 # ...for lazy enough definitions of linearly-independent
20 old_V
= VectorSpace(vs
[0].parent().base_field(), 0)
23 new_V
= span(result
+ [v
])
24 if new_V
.dimension() > old_V
.dimension():
31 def basically_the_same(K1
,K2
):
33 ``True`` if ``K1`` and ``K2`` are basically the same, and ``False``
36 if K1
.lattice_dim() != K2
.lattice_dim():
39 if K1
.nrays() != K2
.nrays():
42 if K1
.dim() != K2
.dim():
45 if lineality(K1
) != lineality(K2
):
48 if K1
.is_solid() != K2
.is_solid():
51 if K1
.is_strictly_convex() != K2
.is_strictly_convex():
54 if len(LL(K1
)) != len(LL(K2
)):
57 C_of_K1
= discrete_complementarity_set(K1
)
58 C_of_K2
= discrete_complementarity_set(K2
)
59 if len(C_of_K1
) != len(C_of_K2
):
62 if len(K1
.facets()) != len(K2
.facets()):
71 Construct the space `W \times W^{\perp}` isomorphic to the ambient space
72 of ``K`` where `W` is equal to the span of ``K``.
74 V
= K
.lattice().vector_space()
76 # Create the space W \times W^{\perp} isomorphic to V.
77 # First we get an orthogonal (but not normal) basis...
78 M
= matrix(V
.base_field(), K
.rays())
79 W_basis
= drop_dependent(K
.rays())
81 W
= V
.subspace_with_basis(W_basis
)
82 W_perp
= W
.complement()
84 return W
.cartesian_product(W_perp
)
89 Construct the IPS isomorphism and its inverse from our paper.
91 Given a cone ``K``, the returned isomorphism will split its ambient
92 vector space `V` into a cartesian product `W \times W^{\perp}` where
93 `W` equals the span of ``K``.
95 V
= K
.lattice().vector_space()
97 (W
, W_perp
) = V_iso
.cartesian_factors()
99 # A space equivalent to V, but using our basis.
100 V_user
= V
.subspace_with_basis( W
.basis() + W_perp
.basis() )
103 # Write v in terms of our custom basis, where the first dim(W)
104 # coordinates are for the W-part of the basis.
105 cs
= V_user
.coordinates(v
)
107 w1
= sum([ V_user
.basis()[idx
]*cs
[idx
]
108 for idx
in range(0, W
.dimension()) ])
109 w2
= sum([ V_user
.basis()[idx
]*cs
[idx
]
110 for idx
in range(W
.dimension(), V
.dimension()) ])
112 return V_iso( (w1
, w2
) )
116 # Crash if the arguments are in the wrong spaces.
119 #w = sum([ sub_w[idx]*W.basis()[idx] for idx in range(0,m) ])
120 #w_prime = sum([ sub_w_prime[idx]*W_perp.basis()[idx]
121 # for idx in range(0,n-m) ])
123 return sum( pair
.cartesian_factors() )
130 def unrestrict_span(K
, K2
=None):
134 _
,phi_inv
= ips_iso(K2
)
135 V_iso
= iso_space(K2
)
136 (W
, W_perp
) = V_iso
.cartesian_factors()
140 w
= sum([ r
[idx
]*W
.basis()[idx
] for idx
in range(0,len(r
)) ])
141 pair
= V_iso( (w
, W_perp
.zero()) )
142 rays
.append( phi_inv(pair
) )
144 L
= ToricLattice(W
.dimension() + W_perp
.dimension())
146 return Cone(rays
, lattice
=L
)
150 def restrict_span(K
, K2
=None):
152 Restrict ``K`` into its own span, or the span of another cone.
156 - ``K2`` -- another cone whose lattice has the same rank as this cone.
160 A new cone in a sublattice.
164 sage: K = Cone([(1,)])
165 sage: restrict_span(K) == K
168 sage: K2 = Cone([(1,0)])
169 sage: restrict_span(K2).rays()
172 sage: K3 = Cone([(1,0,0)])
173 sage: restrict_span(K3).rays()
176 sage: restrict_span(K2) == restrict_span(K3)
181 The projected cone should always be solid::
183 sage: set_random_seed()
184 sage: K = random_cone(max_dim = 8)
185 sage: K_S = restrict_span(K)
189 And the resulting cone should live in a space having the same
190 dimension as the space we restricted it to::
192 sage: set_random_seed()
193 sage: K = random_cone(max_dim = 8)
194 sage: K_S = restrict_span(K, K.dual() )
195 sage: K_S.lattice_dim() == K.dual().dim()
198 This function has ``unrestrict_span()`` as its inverse::
200 sage: set_random_seed()
201 sage: K = random_cone(max_dim = 8, solid=True)
202 sage: J = restrict_span(K)
203 sage: K == unrestrict_span(J,K)
206 This function should not affect the dimension of a cone::
208 sage: set_random_seed()
209 sage: K = random_cone(max_dim = 8)
210 sage: K.dim() == restrict_span(K).dim()
213 Nor should it affect the lineality of a cone::
215 sage: set_random_seed()
216 sage: K = random_cone(max_dim = 8)
217 sage: lineality(K) == lineality(restrict_span(K))
220 No matter which space we restrict to, the lineality should not
223 sage: set_random_seed()
224 sage: K = random_cone(max_dim = 8)
225 sage: lineality(K) >= lineality(restrict_span(K))
227 sage: lineality(K) >= lineality(restrict_span(K, K.dual()))
230 If we do this according to our paper, then the result is proper::
232 sage: set_random_seed()
233 sage: K = random_cone(max_dim = 8, strictly_convex=False, solid=False)
234 sage: K_S = restrict_span(K)
235 sage: P = restrict_span(K_S.dual()).dual()
238 sage: P = restrict_span(K_S, K_S.dual())
244 sage: set_random_seed()
245 sage: K = random_cone(max_dim = 8, strictly_convex=True, solid=False)
246 sage: K_S = restrict_span(K)
247 sage: P = restrict_span(K_S.dual()).dual()
250 sage: P = restrict_span(K_S, K_S.dual())
256 sage: set_random_seed()
257 sage: K = random_cone(max_dim = 8, strictly_convex=False, solid=True)
258 sage: K_S = restrict_span(K)
259 sage: P = restrict_span(K_S.dual()).dual()
262 sage: P = restrict_span(K_S, K_S.dual())
268 sage: set_random_seed()
269 sage: K = random_cone(max_dim = 8, strictly_convex=True, solid=True)
270 sage: K_S = restrict_span(K)
271 sage: P = restrict_span(K_S.dual()).dual()
274 sage: P = restrict_span(K_S, K_S.dual())
278 Test the proposition in our paper concerning the duals, where the
279 subspace `W` is the span of `K^{*}`::
281 sage: set_random_seed()
282 sage: K = random_cone(max_dim = 8, solid=False, strictly_convex=False)
283 sage: K_W = restrict_span(K, K.dual())
284 sage: K_star_W_star = restrict_span(K.dual()).dual()
285 sage: basically_the_same(K_W, K_star_W_star)
290 sage: set_random_seed()
291 sage: K = random_cone(max_dim = 8, solid=True, strictly_convex=False)
292 sage: K_W = restrict_span(K, K.dual())
293 sage: K_star_W_star = restrict_span(K.dual()).dual()
294 sage: basically_the_same(K_W, K_star_W_star)
299 sage: set_random_seed()
300 sage: K = random_cone(max_dim = 8, solid=False, strictly_convex=True)
301 sage: K_W = restrict_span(K, K.dual())
302 sage: K_star_W_star = restrict_span(K.dual()).dual()
303 sage: basically_the_same(K_W, K_star_W_star)
308 sage: set_random_seed()
309 sage: K = random_cone(max_dim = 8, solid=True, strictly_convex=True)
310 sage: K_W = restrict_span(K, K.dual())
311 sage: K_star_W_star = restrict_span(K.dual()).dual()
312 sage: basically_the_same(K_W, K_star_W_star)
320 (W
, W_perp
) = iso_space(K2
).cartesian_factors()
322 ray_pairs
= [ phi(r
) for r
in K
.rays() ]
326 #if any([ w2 != W_perp.zero() for (_, w2) in ray_pairs ]):
327 # msg = 'Cone has nonzero components in W-perp!'
328 # raise ValueError(msg)
330 # Represent the cone in terms of a basis for W, i.e. with smaller
332 ws
= [ W
.coordinate_vector(w1
) for (w1
, _
) in ray_pairs
]
334 L
= ToricLattice(W
.dimension())
336 return Cone(ws
, lattice
=L
)
342 Compute the lineality of this cone.
344 The lineality of a cone is the dimension of the largest linear
345 subspace contained in that cone.
349 A nonnegative integer; the dimension of the largest subspace
350 contained within this cone.
354 .. [Rockafellar] R.T. Rockafellar. Convex Analysis. Princeton
355 University Press, Princeton, 1970.
359 The lineality of the nonnegative orthant is zero, since it clearly
362 sage: K = Cone([(1,0,0), (0,1,0), (0,0,1)])
366 However, if we add another ray so that the entire `x`-axis belongs
367 to the cone, then the resulting cone will have lineality one::
369 sage: K = Cone([(1,0,0), (-1,0,0), (0,1,0), (0,0,1)])
373 If our cone is all of `\mathbb{R}^{2}`, then its lineality is equal
374 to the dimension of the ambient space (i.e. two)::
376 sage: K = Cone([(1,0), (-1,0), (0,1), (0,-1)])
380 Per the definition, the lineality of the trivial cone in a trivial
383 sage: K = Cone([], lattice=ToricLattice(0))
389 The lineality of a cone should be an integer between zero and the
390 dimension of the ambient space, inclusive::
392 sage: set_random_seed()
393 sage: K = random_cone(max_dim = 8)
394 sage: l = lineality(K)
397 sage: (0 <= l) and (l <= K.lattice_dim())
400 A strictly convex cone should have lineality zero::
402 sage: set_random_seed()
403 sage: K = random_cone(max_dim = 8, strictly_convex = True)
408 return K
.linear_subspace().dimension()
411 def discrete_complementarity_set(K
):
413 Compute the discrete complementarity set of this cone.
415 The complementarity set of this cone is the set of all orthogonal
416 pairs `(x,s)` such that `x` is in this cone, and `s` is in its
417 dual. The discrete complementarity set restricts `x` and `s` to be
418 generators of their respective cones.
422 A list of pairs `(x,s)` such that,
424 * `x` is in this cone.
425 * `x` is a generator of this cone.
426 * `s` is in this cone's dual.
427 * `s` is a generator of this cone's dual.
428 * `x` and `s` are orthogonal.
432 The discrete complementarity set of the nonnegative orthant consists
433 of pairs of standard basis vectors::
435 sage: K = Cone([(1,0),(0,1)])
436 sage: discrete_complementarity_set(K)
437 [((1, 0), (0, 1)), ((0, 1), (1, 0))]
439 If the cone consists of a single ray, the second components of the
440 discrete complementarity set should generate the orthogonal
441 complement of that ray::
443 sage: K = Cone([(1,0)])
444 sage: discrete_complementarity_set(K)
445 [((1, 0), (0, 1)), ((1, 0), (0, -1))]
446 sage: K = Cone([(1,0,0)])
447 sage: discrete_complementarity_set(K)
448 [((1, 0, 0), (0, 1, 0)),
449 ((1, 0, 0), (0, -1, 0)),
450 ((1, 0, 0), (0, 0, 1)),
451 ((1, 0, 0), (0, 0, -1))]
453 When the cone is the entire space, its dual is the trivial cone, so
454 the discrete complementarity set is empty::
456 sage: K = Cone([(1,0),(-1,0),(0,1),(0,-1)])
457 sage: discrete_complementarity_set(K)
462 The complementarity set of the dual can be obtained by switching the
463 components of the complementarity set of the original cone::
465 sage: set_random_seed()
466 sage: K1 = random_cone(max_dim=6)
468 sage: expected = [(x,s) for (s,x) in discrete_complementarity_set(K2)]
469 sage: actual = discrete_complementarity_set(K1)
470 sage: sorted(actual) == sorted(expected)
474 V
= K
.lattice().vector_space()
476 # Convert the rays to vectors so that we can compute inner
478 xs
= [V(x
) for x
in K
.rays()]
479 ss
= [V(s
) for s
in K
.dual().rays()]
481 return [(x
,s
) for x
in xs
for s
in ss
if x
.inner_product(s
) == 0]
486 Compute the space `\mathbf{LL}` of all Lyapunov-like transformations
491 A list of matrices forming a basis for the space of all
492 Lyapunov-like transformations on the given cone.
496 The trivial cone has no Lyapunov-like transformations::
498 sage: L = ToricLattice(0)
499 sage: K = Cone([], lattice=L)
503 The Lyapunov-like transformations on the nonnegative orthant are
504 simply diagonal matrices::
506 sage: K = Cone([(1,)])
510 sage: K = Cone([(1,0),(0,1)])
517 sage: K = Cone([(1,0,0),(0,1,0),(0,0,1)])
520 [1 0 0] [0 0 0] [0 0 0]
521 [0 0 0] [0 1 0] [0 0 0]
522 [0 0 0], [0 0 0], [0 0 1]
525 Only the identity matrix is Lyapunov-like on the `L^{3}_{1}` and
526 `L^{3}_{\infty}` cones [Rudolf et al.]_::
528 sage: L31 = Cone([(1,0,1), (0,-1,1), (-1,0,1), (0,1,1)])
536 sage: L3infty = Cone([(0,1,1), (1,0,1), (0,-1,1), (-1,0,1)])
544 If our cone is the entire space, then every transformation on it is
547 sage: K = Cone([(1,0), (-1,0), (0,1), (0,-1)])
548 sage: M = MatrixSpace(QQ,2)
549 sage: M.basis() == LL(K)
554 The inner product `\left< L\left(x\right), s \right>` is zero for
555 every pair `\left( x,s \right)` in the discrete complementarity set
558 sage: set_random_seed()
559 sage: K = random_cone(max_dim=8)
560 sage: C_of_K = discrete_complementarity_set(K)
561 sage: l = [ (L*x).inner_product(s) for (x,s) in C_of_K for L in LL(K) ]
562 sage: sum(map(abs, l))
565 The Lyapunov-like transformations on a cone and its dual are related
566 by transposition, but we're not guaranteed to compute transposed
567 elements of `LL\left( K \right)` as our basis for `LL\left( K^{*}
570 sage: set_random_seed()
571 sage: K = random_cone(max_dim=8)
572 sage: LL2 = [ L.transpose() for L in LL(K.dual()) ]
573 sage: V = VectorSpace( K.lattice().base_field(), K.lattice_dim()^2)
574 sage: LL1_vecs = [ V(m.list()) for m in LL(K) ]
575 sage: LL2_vecs = [ V(m.list()) for m in LL2 ]
576 sage: V.span(LL1_vecs) == V.span(LL2_vecs)
580 V
= K
.lattice().vector_space()
582 C_of_K
= discrete_complementarity_set(K
)
584 tensor_products
= [ s
.tensor_product(x
) for (x
,s
) in C_of_K
]
586 # Sage doesn't think matrices are vectors, so we have to convert
587 # our matrices to vectors explicitly before we can figure out how
588 # many are linearly-indepenedent.
590 # The space W has the same base ring as V, but dimension
591 # dim(V)^2. So it has the same dimension as the space of linear
592 # transformations on V. In other words, it's just the right size
593 # to create an isomorphism between it and our matrices.
594 W
= VectorSpace(V
.base_ring(), V
.dimension()**2)
596 # Turn our matrices into long vectors...
597 vectors
= [ W(m
.list()) for m
in tensor_products
]
599 # Vector space representation of Lyapunov-like matrices
600 # (i.e. vec(L) where L is Luapunov-like).
601 LL_vector
= W
.span(vectors
).complement()
603 # Now construct an ambient MatrixSpace in which to stick our
605 M
= MatrixSpace(V
.base_ring(), V
.dimension())
607 matrix_basis
= [ M(v
.list()) for v
in LL_vector
.basis() ]
613 def lyapunov_rank(K
):
615 Compute the Lyapunov (or bilinearity) rank of this cone.
617 The Lyapunov rank of a cone can be thought of in (mainly) two ways:
619 1. The dimension of the Lie algebra of the automorphism group of the
622 2. The dimension of the linear space of all Lyapunov-like
623 transformations on the cone.
627 A closed, convex polyhedral cone.
631 An integer representing the Lyapunov rank of the cone. If the
632 dimension of the ambient vector space is `n`, then the Lyapunov rank
633 will be between `1` and `n` inclusive; however a rank of `n-1` is
634 not possible (see the first reference).
638 In the references, the cones are always assumed to be proper. We
639 do not impose this restriction.
647 The codimension formula from the second reference is used. We find
648 all pairs `(x,s)` in the complementarity set of `K` such that `x`
649 and `s` are rays of our cone. It is known that these vectors are
650 sufficient to apply the codimension formula. Once we have all such
651 pairs, we "brute force" the codimension formula by finding all
652 linearly-independent `xs^{T}`.
656 .. [Gowda/Tao] M.S. Gowda and J. Tao. On the bilinearity rank of a proper
657 cone and Lyapunov-like transformations, Mathematical Programming, 147
660 .. [Orlitzky/Gowda] M. Orlitzky and M. S. Gowda. The Lyapunov Rank of an
661 Improper Cone. Work in-progress.
663 .. [Rudolf et al.] G. Rudolf, N. Noyan, D. Papp, and F. Alizadeh, Bilinear
664 optimality constraints for the cone of positive polynomials,
665 Mathematical Programming, Series B, 129 (2011) 5-31.
669 The nonnegative orthant in `\mathbb{R}^{n}` always has rank `n`
672 sage: positives = Cone([(1,)])
673 sage: lyapunov_rank(positives)
675 sage: quadrant = Cone([(1,0), (0,1)])
676 sage: lyapunov_rank(quadrant)
678 sage: octant = Cone([(1,0,0), (0,1,0), (0,0,1)])
679 sage: lyapunov_rank(octant)
682 The full space `\mathbb{R}^{n}` has Lyapunov rank `n^{2}`
685 sage: R5 = VectorSpace(QQ, 5)
686 sage: gs = R5.basis() + [ -r for r in R5.basis() ]
688 sage: lyapunov_rank(K)
691 The `L^{3}_{1}` cone is known to have a Lyapunov rank of one
694 sage: L31 = Cone([(1,0,1), (0,-1,1), (-1,0,1), (0,1,1)])
695 sage: lyapunov_rank(L31)
698 Likewise for the `L^{3}_{\infty}` cone [Rudolf et al.]_::
700 sage: L3infty = Cone([(0,1,1), (1,0,1), (0,-1,1), (-1,0,1)])
701 sage: lyapunov_rank(L3infty)
704 A single ray in `n` dimensions should have Lyapunov rank `n^{2} - n
705 + 1` [Orlitzky/Gowda]_::
707 sage: K = Cone([(1,0,0,0,0)])
708 sage: lyapunov_rank(K)
710 sage: K.lattice_dim()**2 - K.lattice_dim() + 1
713 A subspace (of dimension `m`) in `n` dimensions should have a
714 Lyapunov rank of `n^{2} - m\left(n - m)` [Orlitzky/Gowda]_::
716 sage: e1 = (1,0,0,0,0)
717 sage: neg_e1 = (-1,0,0,0,0)
718 sage: e2 = (0,1,0,0,0)
719 sage: neg_e2 = (0,-1,0,0,0)
720 sage: z = (0,0,0,0,0)
721 sage: K = Cone([e1, neg_e1, e2, neg_e2, z, z, z])
722 sage: lyapunov_rank(K)
724 sage: K.lattice_dim()**2 - K.dim()*codim(K)
727 The Lyapunov rank should be additive on a product of proper cones
730 sage: L31 = Cone([(1,0,1), (0,-1,1), (-1,0,1), (0,1,1)])
731 sage: octant = Cone([(1,0,0), (0,1,0), (0,0,1)])
732 sage: K = L31.cartesian_product(octant)
733 sage: lyapunov_rank(K) == lyapunov_rank(L31) + lyapunov_rank(octant)
736 Two isomorphic cones should have the same Lyapunov rank [Rudolf et al.]_.
737 The cone ``K`` in the following example is isomorphic to the nonnegative
738 octant in `\mathbb{R}^{3}`::
740 sage: K = Cone([(1,2,3), (-1,1,0), (1,0,6)])
741 sage: lyapunov_rank(K)
744 The dual cone `K^{*}` of ``K`` should have the same Lyapunov rank as ``K``
745 itself [Rudolf et al.]_::
747 sage: K = Cone([(2,2,4), (-1,9,0), (2,0,6)])
748 sage: lyapunov_rank(K) == lyapunov_rank(K.dual())
753 The Lyapunov rank should be additive on a product of proper cones
756 sage: set_random_seed()
757 sage: K1 = random_cone(max_dim=8, strictly_convex=True, solid=True)
758 sage: K2 = random_cone(max_dim=8, strictly_convex=True, solid=True)
759 sage: K = K1.cartesian_product(K2)
760 sage: lyapunov_rank(K) == lyapunov_rank(K1) + lyapunov_rank(K2)
763 The dual cone `K^{*}` of ``K`` should have the same Lyapunov rank as ``K``
764 itself [Rudolf et al.]_::
766 sage: set_random_seed()
767 sage: K = random_cone(max_dim=8)
768 sage: lyapunov_rank(K) == lyapunov_rank(K.dual())
771 Make sure we exercise the non-strictly-convex/non-solid case::
773 sage: set_random_seed()
774 sage: K = random_cone(max_dim=8, strictly_convex=False, solid=False)
775 sage: lyapunov_rank(K) == lyapunov_rank(K.dual())
778 Let's check the other permutations as well, just to be sure::
780 sage: set_random_seed()
781 sage: K = random_cone(max_dim=8, strictly_convex=False, solid=True)
782 sage: lyapunov_rank(K) == lyapunov_rank(K.dual())
787 sage: set_random_seed()
788 sage: K = random_cone(max_dim=8, strictly_convex=True, solid=False)
789 sage: lyapunov_rank(K) == lyapunov_rank(K.dual())
794 sage: set_random_seed()
795 sage: K = random_cone(max_dim=8, strictly_convex=True, solid=True)
796 sage: lyapunov_rank(K) == lyapunov_rank(K.dual())
799 The Lyapunov rank of a proper polyhedral cone in `n` dimensions can
800 be any number between `1` and `n` inclusive, excluding `n-1`
801 [Gowda/Tao]_. By accident, the `n-1` restriction will hold for the
802 trivial cone in a trivial space as well. However, in zero dimensions,
803 the Lyapunov rank of the trivial cone will be zero::
805 sage: set_random_seed()
806 sage: K = random_cone(max_dim=8, strictly_convex=True, solid=True)
807 sage: b = lyapunov_rank(K)
808 sage: n = K.lattice_dim()
809 sage: (n == 0 or 1 <= b) and b <= n
814 In fact [Orlitzky/Gowda]_, no closed convex polyhedral cone can have
815 Lyapunov rank `n-1` in `n` dimensions::
817 sage: set_random_seed()
818 sage: K = random_cone(max_dim=8)
819 sage: b = lyapunov_rank(K)
820 sage: n = K.lattice_dim()
824 The calculation of the Lyapunov rank of an improper cone can be
825 reduced to that of a proper cone [Orlitzky/Gowda]_::
827 sage: set_random_seed()
828 sage: K = random_cone(max_dim=8)
829 sage: actual = lyapunov_rank(K)
830 sage: K_S = restrict_span(K)
831 sage: P = restrict_span(K_S.dual()).dual()
832 sage: l = lineality(K)
834 sage: expected = lyapunov_rank(P) + K.dim()*(l + c) + c**2
835 sage: actual == expected
838 The Lyapunov rank of a proper cone is just the dimension of ``LL(K)``::
840 sage: set_random_seed()
841 sage: K = random_cone(max_dim=8, strictly_convex=True, solid=True)
842 sage: lyapunov_rank(K) == len(LL(K))
845 In fact the same can be said of any cone. These additional tests
846 just increase our confidence that the reduction scheme works::
848 sage: set_random_seed()
849 sage: K = random_cone(max_dim=8, strictly_convex=True, solid=False)
850 sage: lyapunov_rank(K) == len(LL(K))
855 sage: set_random_seed()
856 sage: K = random_cone(max_dim=8, strictly_convex=False, solid=True)
857 sage: lyapunov_rank(K) == len(LL(K))
862 sage: set_random_seed()
863 sage: K = random_cone(max_dim=8, strictly_convex=False, solid=False)
864 sage: lyapunov_rank(K) == len(LL(K))
876 # K is not solid, project onto its span.
880 beta
+= m
*(n
- m
) + (n
- m
)**2
883 # K is not pointed, project its dual onto its span.
884 # Uses a proposition from our paper, i.e. this is
885 # equivalent to K = restrict_span(K.dual()).dual()
886 #K = restrict_span(intersect_span(K,K.dual()), K.dual())
887 K
= restrict_span(K
, K
.dual())