]> gitweb.michael.orlitzky.com - sage.d.git/blobdiff - mjo/cone/rearrangement.py
mjo/**/*.py: drop obsolete set_random_seed().
[sage.d.git] / mjo / cone / rearrangement.py
index b2ccfcb58f8aa58b52a8d724ca42725593ed7463..0bbf95bb7603bed819d4b9e91a4b148856dc464c 100644 (file)
@@ -1,6 +1,6 @@
 from sage.all import *
 
-def rearrangement_cone(p,n):
+def rearrangement_cone(p,n,lattice=None):
     r"""
     Return the rearrangement cone of order ``p`` in ``n`` dimensions.
 
@@ -19,14 +19,43 @@ def rearrangement_cone(p,n):
 
     INPUT:
 
-    - ``p`` -- The number of components to "rearrange."
+      - ``p`` -- The number of components to "rearrange."
 
-    - ``n`` -- The dimension of the ambient space for the resulting cone.
+      - ``n`` -- The dimension of the ambient space for the resulting cone.
+
+      - ``lattice`` -- (default: ``None``) an ambient lattice of rank ``n``
+                       to be passed to the :func:`Cone` constructor.
 
     OUTPUT:
 
     A polyhedral closed convex cone object representing a rearrangement
-    cone of order ``p`` in ``n`` dimensions.
+    cone of order ``p`` in ``n`` dimensions. Each generating ray will
+    have the integer ring as its base ring.
+
+    If a ``lattice`` was specified, then the resulting cone will live in
+    that lattice unless its rank is incompatible with the dimension
+    ``n`` (in which case a ``ValueError`` is raised).
+
+    ALGORITHM:
+
+    The generators for the rearrangement cone are given by [Jeong]_
+    Theorem 5.2.3.
+
+    REFERENCES:
+
+    .. [GowdaJeong] Muddappa Seetharama Gowda and Juyoung Jeong.
+       Spectral cones in Euclidean Jordan algebras.
+       Linear Algebra and its Applications, 509, 286-305.
+       doi:10.1016/j.laa.2016.08.004
+
+    .. [HenrionSeeger] Rene Henrion and Alberto Seeger.
+       Inradius and Circumradius of Various Convex Cones Arising in
+       Applications. Set-Valued and Variational Analysis, 18(3-4),
+       483-511, 2010. doi:10.1007/s11228-010-0150-z
+
+    .. [Jeong] Juyoung Jeong.
+       Spectral sets and functions on Euclidean Jordan algebras.
+       University of Maryland, Baltimore County, Ph.D. thesis, 2017.
 
     SETUP::
 
@@ -50,15 +79,17 @@ def rearrangement_cone(p,n):
         sage: rearrangement_cone(5,5).lineality()
         4
 
-    All rearrangement cones are proper::
+    All rearrangement cones are proper when ``p`` is less than ``n`` by
+    [Jeong]_ Proposition 5.2.1::
 
         sage: all( rearrangement_cone(p,n).is_proper()
         ....:              for n in range(10)
-        ....:              for p in range(n) )
+        ....:              for p in range(1, n) )
         True
 
     The Lyapunov rank of the rearrangement cone of order ``p`` in ``n``
-    dimensions is ``n`` for ``p == 1`` or ``p == n`` and one otherwise::
+    dimensions is ``n`` for ``p == 1`` or ``p == n`` and one otherwise,
+    by [Jeong]_ Corollary 5.2.4::
 
         sage: all( rearrangement_cone(p,n).lyapunov_rank() == n
         ....:              for n in range(2, 10)
@@ -71,7 +102,8 @@ def rearrangement_cone(p,n):
 
     TESTS:
 
-    The rearrangement cone is permutation-invariant::
+    All rearrangement cones are permutation-invariant by [Jeong]_
+    Proposition 5.2.1::
 
         sage: n = ZZ.random_element(2,10).abs()
         sage: p = ZZ.random_element(1,n)
@@ -80,68 +112,83 @@ def rearrangement_cone(p,n):
         sage: all( K.contains(P*r) for r in K )
         True
 
-    """
-
-    def d(j):
-        v = [1]*n    # Create the list of all ones...
-        v[j] = 1 - p # Now "fix" the ``j``th entry.
-        return v
-
-    V = VectorSpace(QQ, n)
-    G = V.basis() + [ d(j) for j in range(n) ]
-    return Cone(G)
-
-
-def has_rearrangement_property(v, p):
-    r"""
-    Test if the vector ``v`` has the "rearrangement property."
-
-    The rearrangement cone of order ``p`` in `n` dimensions has its
-    members vectors of length `n`. The "rearrangement property,"
-    satisfied by its elements, is to have its smallest ``p`` components
-    sum to a nonnegative number.
-
-    We believe that we have a description of the extreme vectors of the
-    rearrangement cone: see ``rearrangement_cone()``. This function is
-    used to test that conic combinations of those extreme vectors are in
-    fact elements of the rearrangement cone. We can't test all conic
-    combinations, obviously, but we can test a random one.
+    The smallest ``p`` components of every element of the rearrangement
+    cone should sum to a nonnegative number (this tests that the
+    generators really are what we think they are)::
+
+        sage: def _has_rearrangement_property(v,p):
+        ....:     return sum( sorted(v)[0:p] ) >= 0
+        sage: all( _has_rearrangement_property(
+        ....:      rearrangement_cone(p,n).random_element(),
+        ....:      p
+        ....:    )
+        ....:    for n in range(2, 10)
+        ....:    for p in range(1, n-1)
+        ....: )
+        True
 
-    To become more sure of the result, generate a bunch of vectors with
-    ``random_element()`` and test them with this function.
+    The rearrangenent cone of order ``p`` is contained in the rearrangement
+    cone of order ``p + 1`` by [Jeong]_ Proposition 5.2.1::
 
-    INPUT:
-
-    - ``v`` -- An element of a cone suspected of being the rearrangement
-               cone of order ``p``.
+        sage: n = ZZ.random_element(2,10)
+        sage: p = ZZ.random_element(1,n)
+        sage: K1 = rearrangement_cone(p,n)
+        sage: K2 = rearrangement_cone(p+1,n)
+        sage: all( x in K2 for x in K1 )
+        True
 
-    - ``p`` -- The suspected order of the rearrangement cone.
+    The rearrangement cone of order ``p`` is linearly isomorphic to the
+    rearrangement cone of order ``n - p`` when ``p`` is less than ``n``,
+    by [Jeong]_ Proposition 5.2.1::
 
-    OUTPUT:
+        sage: n = ZZ.random_element(2,10)
+        sage: p = ZZ.random_element(1,n)
+        sage: K1 = rearrangement_cone(p,n)
+        sage: K2 = rearrangement_cone(n-p, n)
+        sage: Mp = (1/p)*matrix.ones(QQ,n) - identity_matrix(QQ,n)
+        sage: Cone( (Mp*K2.rays()).columns() ).is_equivalent(K1)
+        True
 
-    If ``v`` has the rearrangement property (that is, if its smallest ``p``
-    components sum to a nonnegative number), ``True`` is returned. Otherwise
-    ``False`` is returned.
+    The order ``p`` should be between ``1`` and ``n``, inclusive::
 
-    SETUP::
+        sage: rearrangement_cone(0,3)
+        Traceback (most recent call last):
+        ...
+        ValueError: order p=0 should be between 1 and n=3, inclusive
+        sage: rearrangement_cone(5,3)
+        Traceback (most recent call last):
+        ...
+        ValueError: order p=5 should be between 1 and n=3, inclusive
 
-        sage: from mjo.cone.rearrangement import (has_rearrangement_property,
-        ....:                                     rearrangement_cone)
+    If a ``lattice`` was given, it is actually used::
 
-    EXAMPLES:
+        sage: L = ToricLattice(3, 'M')
+        sage: rearrangement_cone(2, 3, lattice=L)
+        3-d cone in 3-d lattice M
 
-    Every element of a rearrangement cone should have the property::
+    Unless the rank of the lattice disagrees with ``n``::
 
-        sage: set_random_seed()
-        sage: all( has_rearrangement_property(
-        ....:        rearrangement_cone(p,n).random_element(),
-        ....:        p
-        ....:      )
-        ....:      for n in range(2, 10)
-        ....:      for p in range(1, n-1)
-        ....: )
-        True
+        sage: L = ToricLattice(1, 'M')
+        sage: rearrangement_cone(2, 3, lattice=L)
+        Traceback (most recent call last):
+        ...
+        ValueError: lattice rank=1 and dimension n=3 are incompatible
 
     """
-    components = sorted(v)[0:p]
-    return sum(components) >= 0
+    if p < 1 or p > n:
+        raise ValueError('order p=%d should be between 1 and n=%d, inclusive'
+                         %
+                         (p,n))
+
+    if lattice is None:
+        lattice = ToricLattice(n)
+
+    if lattice.rank() != n:
+        raise ValueError('lattice rank=%d and dimension n=%d are incompatible'
+                         %
+                         (lattice.rank(), n))
+
+    I = identity_matrix(ZZ,n)
+    M = matrix.ones(ZZ,n) - p*I
+    G = identity_matrix(ZZ,n).rows() + M.rows()
+    return Cone(G, lattice=lattice)