]>
gitweb.michael.orlitzky.com - sage.d.git/blob - mjo/cone/rearrangement.py
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 rearrangement_cone(p
,n
):
12 Return the rearrangement cone of order ``p`` in ``n`` dimensions.
14 The rearrangement cone in ``n`` dimensions has as its elements
15 vectors of length ``n``. For inclusion in the cone, the smallest
16 ``p`` components of a vector must sum to a nonnegative number.
18 For example, the rearrangement cone of order ``p == 1`` has its
19 single smallest component nonnegative. This implies that all
20 components are nonnegative, and that therefore the rearrangement
21 cone of order one is the nonnegative orthant.
23 When ``p == n``, the sum of all components of a vector must be
24 nonnegative for inclusion in the cone. That is, the cone is a
25 half-space in ``n`` dimensions.
29 - ``p`` -- The number of components to "rearrange."
31 - ``n`` -- The dimension of the ambient space for the resulting cone.
35 A polyhedral closed convex cone object representing a rearrangement
36 cone of order ``p`` in ``n`` dimensions.
40 The rearrangement cones of order one are nonnegative orthants::
42 sage: rearrangement_cone(1,1) == Cone([(1,)])
44 sage: rearrangement_cone(1,2) == Cone([(0,1),(1,0)])
46 sage: rearrangement_cone(1,3) == Cone([(0,0,1),(0,1,0),(1,0,0)])
49 When ``p == n``, the resulting cone will be a half-space, so we
50 expect its lineality to be one less than ``n`` because it will
51 contain a hyperplane but is not the entire space::
53 sage: rearrangement_cone(5,5).lineality()
56 All rearrangement cones are proper::
58 sage: all([ rearrangement_cone(p,n).is_proper()
59 ....: for n in range(10)
60 ....: for p in range(n) ])
63 The Lyapunov rank of the rearrangement cone of order ``p`` in ``n``
64 dimensions is ``n`` for ``p == 1`` or ``p == n`` and one otherwise::
66 sage: all([ rearrangement_cone(p,n).lyapunov_rank() == n
67 ....: for n in range(2, 10)
68 ....: for p in [1, n-1] ])
70 sage: all([ rearrangement_cone(p,n).lyapunov_rank() == 1
71 ....: for n in range(3, 10)
72 ....: for p in range(2, n-1) ])
77 The rearrangement cone is permutation-invariant::
79 sage: n = ZZ.random_element(2,10).abs()
80 sage: p = ZZ.random_element(1,n)
81 sage: K = rearrangement_cone(p,n)
82 sage: P = SymmetricGroup(n).random_element().matrix()
83 sage: all([ K.contains(P*r) for r in K.rays() ])
89 v
= [1]*n
# Create the list of all ones...
90 v
[j
] = 1 - p
# Now "fix" the ``j``th entry.
93 V
= VectorSpace(QQ
, n
)
94 G
= V
.basis() + [ d(j
) for j
in range(n
) ]
98 def has_rearrangement_property(v
, p
):
100 Test if the vector ``v`` has the "rearrangement property."
102 The rearrangement cone of order ``p`` in `n` dimensions has its
103 members vectors of length `n`. The "rearrangement property,"
104 satisfied by its elements, is to have its smallest ``p`` components
105 sum to a nonnegative number.
107 We believe that we have a description of the extreme vectors of the
108 rearrangement cone: see ``rearrangement_cone()``. This function is
109 used to test that conic combinations of those extreme vectors are in
110 fact elements of the rearrangement cone. We can't test all conic
111 combinations, obviously, but we can test a random one.
113 To become more sure of the result, generate a bunch of vectors with
114 ``random_element()`` and test them with this function.
118 - ``v`` -- An element of a cone suspected of being the rearrangement
121 - ``p`` -- The suspected order of the rearrangement cone.
125 If ``v`` has the rearrangement property (that is, if its smallest ``p``
126 components sum to a nonnegative number), ``True`` is returned. Otherwise
127 ``False`` is returned.
131 Every element of a rearrangement cone should have the property::
133 sage: for n in range(2,10):
134 ....: for p in range(1, n-1):
135 ....: K = rearrangement_cone(p,n)
136 ....: v = K.random_element()
137 ....: if not has_rearrangement_property(v,p): print v
140 components
= sorted(v
)[0:p
]
141 return sum(components
) >= 0