mjo/cone/rearrangement.py

   1 from sage.all import *
   2
   3 def rearrangement_cone(p,n,lattice=None):
   4     r"""
   5     Return the rearrangement cone of order ``p`` in ``n`` dimensions.
   6
   7     The rearrangement cone in ``n`` dimensions has as its elements
   8     vectors of length ``n``. For inclusion in the cone, the smallest
   9     ``p`` components of a vector must sum to a nonnegative number.
  10
  11     For example, the rearrangement cone of order ``p == 1`` has its
  12     single smallest component nonnegative. This implies that all
  13     components are nonnegative, and that therefore the rearrangement
  14     cone of order one is the nonnegative orthant.
  15
  16     When ``p == n``, the sum of all components of a vector must be
  17     nonnegative for inclusion in the cone. That is, the cone is a
  18     half-space in ``n`` dimensions.
  19
  20     INPUT:
  21
  22       - ``p`` -- The number of components to "rearrange."
  23
  24       - ``n`` -- The dimension of the ambient space for the resulting cone.
  25
  26       - ``lattice`` -- (default: ``None``) an ambient lattice of rank ``n``
  27                        to be passed to the :func:`Cone` constructor.
  28
  29     OUTPUT:
  30
  31     A polyhedral closed convex cone object representing a rearrangement
  32     cone of order ``p`` in ``n`` dimensions. Each generating ray will
  33     have the integer ring as its base ring.
  34
  35     If a ``lattice`` was specified, then the resulting cone will live in
  36     that lattice unless its rank is incompatible with the dimension
  37     ``n`` (in which case a ``ValueError`` is raised).
  38
  39     ALGORITHM:
  40
  41     The generators for the rearrangement cone are given by [Jeong]_
  42     Theorem 5.2.3.
  43
  44     REFERENCES:
  45
  46     .. [GowdaJeong] Muddappa Seetharama Gowda and Juyoung Jeong.
  47        Spectral cones in Euclidean Jordan algebras.
  48        Linear Algebra and its Applications, 509, 286-305.
  49        doi:10.1016/j.laa.2016.08.004
  50
  51     .. [HenrionSeeger] Rene Henrion and Alberto Seeger.
  52        Inradius and Circumradius of Various Convex Cones Arising in
  53        Applications. Set-Valued and Variational Analysis, 18(3-4),
  54        483-511, 2010. doi:10.1007/s11228-010-0150-z
  55
  56     .. [Jeong] Juyoung Jeong.
  57        Spectral sets and functions on Euclidean Jordan algebras.
  58        University of Maryland, Baltimore County, Ph.D. thesis, 2017.
  59
  60     SETUP::
  61
  62         sage: from mjo.cone.rearrangement import rearrangement_cone
  63
  64     EXAMPLES:
  65
  66     The rearrangement cones of order one are nonnegative orthants::
  67
  68         sage: rearrangement_cone(1,1) == Cone([(1,)])
  69         True
  70         sage: rearrangement_cone(1,2) == Cone([(0,1),(1,0)])
  71         True
  72         sage: rearrangement_cone(1,3) == Cone([(0,0,1),(0,1,0),(1,0,0)])
  73         True
  74
  75     When ``p == n``, the resulting cone will be a half-space, so we
  76     expect its lineality to be one less than ``n`` because it will
  77     contain a hyperplane but is not the entire space::
  78
  79         sage: rearrangement_cone(5,5).lineality()
  80         4
  81
  82     All rearrangement cones are proper when ``p`` is less than ``n`` by
  83     [Jeong]_ Proposition 5.2.1::
  84
  85         sage: all( rearrangement_cone(p,n).is_proper()
  86         ....:              for n in range(10)
  87         ....:              for p in range(1, n) )
  88         True
  89
  90     The Lyapunov rank of the rearrangement cone of order ``p`` in ``n``
  91     dimensions is ``n`` for ``p == 1`` or ``p == n`` and one otherwise,
  92     by [Jeong]_ Corollary 5.2.4::
  93
  94         sage: all( rearrangement_cone(p,n).lyapunov_rank() == n
  95         ....:              for n in range(2, 10)
  96         ....:              for p in [1, n-1] )
  97         True
  98         sage: all( rearrangement_cone(p,n).lyapunov_rank() == 1
  99         ....:              for n in range(3, 10)
 100         ....:              for p in range(2, n-1) )
 101         True
 102
 103     TESTS:
 104
 105     All rearrangement cones are permutation-invariant by [Jeong]_
 106     Proposition 5.2.1::
 107
 108         sage: n = ZZ.random_element(2,10).abs()
 109         sage: p = ZZ.random_element(1,n)
 110         sage: K = rearrangement_cone(p,n)
 111         sage: P = SymmetricGroup(n).random_element().matrix()
 112         sage: all( K.contains(P*r) for r in K )
 113         True
 114
 115     The smallest ``p`` components of every element of the rearrangement
 116     cone should sum to a nonnegative number (this tests that the
 117     generators really are what we think they are)::
 118
 119         sage: def _has_rearrangement_property(v,p):
 120         ....:     return sum( sorted(v)[0:p] ) >= 0
 121         sage: all( _has_rearrangement_property(
 122         ....:      rearrangement_cone(p,n).random_element(),
 123         ....:      p
 124         ....:    )
 125         ....:    for n in range(2, 10)
 126         ....:    for p in range(1, n-1)
 127         ....: )
 128         True
 129
 130     The rearrangenent cone of order ``p`` is contained in the rearrangement
 131     cone of order ``p + 1`` by [Jeong]_ Proposition 5.2.1::
 132
 133         sage: n = ZZ.random_element(2,10)
 134         sage: p = ZZ.random_element(1,n)
 135         sage: K1 = rearrangement_cone(p,n)
 136         sage: K2 = rearrangement_cone(p+1,n)
 137         sage: all( x in K2 for x in K1 )
 138         True
 139
 140     The rearrangement cone of order ``p`` is linearly isomorphic to the
 141     rearrangement cone of order ``n - p`` when ``p`` is less than ``n``,
 142     by [Jeong]_ Proposition 5.2.1::
 143
 144         sage: n = ZZ.random_element(2,10)
 145         sage: p = ZZ.random_element(1,n)
 146         sage: K1 = rearrangement_cone(p,n)
 147         sage: K2 = rearrangement_cone(n-p, n)
 148         sage: Mp = (1/p)*matrix.ones(QQ,n) - identity_matrix(QQ,n)
 149         sage: Cone( (Mp*K2.rays()).columns() ).is_equivalent(K1)
 150         True
 151
 152     The order ``p`` should be between ``1`` and ``n``, inclusive::
 153
 154         sage: rearrangement_cone(0,3)
 155         Traceback (most recent call last):
 156         ...
 157         ValueError: order p=0 should be between 1 and n=3, inclusive
 158         sage: rearrangement_cone(5,3)
 159         Traceback (most recent call last):
 160         ...
 161         ValueError: order p=5 should be between 1 and n=3, inclusive
 162
 163     If a ``lattice`` was given, it is actually used::
 164
 165         sage: L = ToricLattice(3, 'M')
 166         sage: rearrangement_cone(2, 3, lattice=L)
 167         3-d cone in 3-d lattice M
 168
 169     Unless the rank of the lattice disagrees with ``n``::
 170
 171         sage: L = ToricLattice(1, 'M')
 172         sage: rearrangement_cone(2, 3, lattice=L)
 173         Traceback (most recent call last):
 174         ...
 175         ValueError: lattice rank=1 and dimension n=3 are incompatible
 176
 177     """
 178     if p < 1 or p > n:
 179         raise ValueError('order p=%d should be between 1 and n=%d, inclusive'
 180                          %
 181                          (p,n))
 182
 183     if lattice is None:
 184         lattice = ToricLattice(n)
 185
 186     if lattice.rank() != n:
 187         raise ValueError('lattice rank=%d and dimension n=%d are incompatible'
 188                          %
 189                          (lattice.rank(), n))
 190
 191     I = identity_matrix(ZZ,n)
 192     M = matrix.ones(ZZ,n) - p*I
 193     G = identity_matrix(ZZ,n).rows() + M.rows()
 194     return Cone(G, lattice=lattice)