from sage.all import *

def rearrangement_cone(p,n):
    r"""
    Return the rearrangement cone of order ``p`` in ``n`` dimensions.

    The rearrangement cone in ``n`` dimensions has as its elements
    vectors of length ``n``. For inclusion in the cone, the smallest
    ``p`` components of a vector must sum to a nonnegative number.

    For example, the rearrangement cone of order ``p == 1`` has its
    single smallest component nonnegative. This implies that all
    components are nonnegative, and that therefore the rearrangement
    cone of order one is the nonnegative orthant.

    When ``p == n``, the sum of all components of a vector must be
    nonnegative for inclusion in the cone. That is, the cone is a
    half-space in ``n`` dimensions.

    INPUT:

    - ``p`` -- The number of components to "rearrange."

    - ``n`` -- The dimension of the ambient space for the resulting cone.

    OUTPUT:

    A polyhedral closed convex cone object representing a rearrangement
    cone of order ``p`` in ``n`` dimensions.

    REFERENCES:

    .. [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

    .. [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

    .. [Jeong] Juyoung Jeong.
       Spectral sets and functions on Euclidean Jordan algebras.

    SETUP::

        sage: from mjo.cone.rearrangement import rearrangement_cone

    EXAMPLES:

    The rearrangement cones of order one are nonnegative orthants::

        sage: rearrangement_cone(1,1) == Cone([(1,)])
        True
        sage: rearrangement_cone(1,2) == Cone([(0,1),(1,0)])
        True
        sage: rearrangement_cone(1,3) == Cone([(0,0,1),(0,1,0),(1,0,0)])
        True

    When ``p == n``, the resulting cone will be a half-space, so we
    expect its lineality to be one less than ``n`` because it will
    contain a hyperplane but is not the entire space::

        sage: rearrangement_cone(5,5).lineality()
        4

    All rearrangement cones are proper::

        sage: all( rearrangement_cone(p,n).is_proper()
        ....:              for n in xrange(10)
        ....:              for p in xrange(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::

        sage: all( rearrangement_cone(p,n).lyapunov_rank() == n
        ....:              for n in xrange(2, 10)
        ....:              for p in [1, n-1] )
        True
        sage: all( rearrangement_cone(p,n).lyapunov_rank() == 1
        ....:              for n in xrange(3, 10)
        ....:              for p in xrange(2, n-1) )
        True

    TESTS:

    The rearrangement cone is permutation-invariant::

        sage: n = ZZ.random_element(2,10).abs()
        sage: p = ZZ.random_element(1,n)
        sage: K = rearrangement_cone(p,n)
        sage: P = SymmetricGroup(n).random_element().matrix()
        sage: all( K.contains(P*r) for r in K )
        True

    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: set_random_seed()
        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 xrange(2, 10)
        ....:    for p in xrange(1, n-1)
        ....: )
        True

    The rearrangenent cone of order ``p`` is contained in the
    rearrangement cone of order ``p + 1``::

        sage: set_random_seed()
        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

    The order ``p`` should be between ``1`` and ``n``, inclusive::

        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

    """
    if p < 1 or p > n:
        raise ValueError('order p=%d should be between 1 and n=%d, inclusive'
                         %
                         (p,n))

    def d(j):
        v = [1]*n    # Create the list of all ones...
        v[j] = 1 - p # Now "fix" the ``j``th entry.
        return v

    G = identity_matrix(ZZ,n).rows() + [ d(j) for j in xrange(n) ]
    return Cone(G)