#! /usr/bin/env python
#
def mono_rank_grlex ( m, x ):

#*****************************************************************************80
#
## MONO_RANK_GRLEX computes the graded lexicographic rank of a monomial.
#
#  Discussion:
#
#    The graded lexicographic ordering is used, over all monomials in 
#    M dimensions, for total degree = 0, 1, 2, ...
#
#    For example, if M = 3, the ranking begins:
#
#    Rank  Sum    1  2  3
#    ----  ---   -- -- --
#       1    0    0  0  0
#
#       2    1    0  0  1
#       3    1    0  1  0
#       4    1    1  0  1
#
#       5    2    0  0  2
#       6    2    0  1  1
#       7    2    0  2  0
#       8    2    1  0  1
#       9    2    1  1  0
#      10    2    2  0  0
#
#      11    3    0  0  3
#      12    3    0  1  2
#      13    3    0  2  1
#      14    3    0  3  0
#      15    3    1  0  2
#      16    3    1  1  1
#      17    3    1  2  0
#      18    3    2  0  1
#      19    3    2  1  0
#      20    3    3  0  0
#
#      21    4    0  0  4
#      ..   ..   .. .. ..
#
#  Licensing:
#
#   This code is distributed under the GNU LGPL license.
#
#  Modified:
#
#    31 October 2015
#
#  Author:
#
#    John Burkardt
#
#  Parameters:
#
#    Input, integer M, the spatial dimension.
#    1 <= M.
#
#    Input, integer X[M], the composition.
#    For each 1 <= I <= M, we have 0 <= X(I).
#
#    Output, integer RANK, the rank.
#
  import numpy as np
  from i4_choose import i4_choose
  from i4vec_sum import i4vec_sum
  from sys import exit
#
#  Ensure that 1 <= M.
#
  if ( m < 1 ):
    print ( '' )
    print ( 'MONO_RANK_GRLEX - Fatal error!' )
    print ( '  M < 1' )
    exit ( 'MONO_RANK_GRLEX - Fatal error!' )
#
#  Ensure that 0 <= X(I).
#
  for i in range ( 0, m ):
    if ( x[i] < 0 ):
      print ( '' )
      print ( 'MONO_RANK_GRLEX - Fatal error!' )
      print ( '  X[I] < 0' )
      exit ( 'MONO_RANK_GRLEX - Fatal error!' )
#
#  NM = sum ( X )
#
  nm = i4vec_sum ( m, x )
#
#  Convert to KSUBSET format.
#
  ns = nm + m - 1
  ks = m - 1
  if ( 0 < ks ):
    xs = np.zeros ( ks, dtype = np.int32 )
    xs[0] = x[0] + 1
    for i in range ( 2, m ):
      xs[i-1] = xs[i-2] + x[i-1] + 1
#
#  Compute the rank.
#
  rank = 1

  for i in range ( 1, ks + 1 ):
    if ( i == 1 ):
      tim1 = 0
    else:
      tim1 = xs[i-2]

    if ( tim1 + 1 <= xs[i-1] - 1 ):
      for j in range ( tim1 + 1, xs[i-1] ):
        rank = rank + i4_choose ( ns - j, ks - i )

  for n in range ( 0, nm ):
    rank = rank + i4_choose ( n + m - 1, n )

  return rank

def mono_rank_grlex_test ( ):

#******************************************************************************/
#
## MONO_RANK_GRLEX_TEST tests MONO_RANK_GRLEX.
#
#  Licensing:
#
#    This code is distributed under the GNU LGPL license.
#
#  Modified:
#
#    25 October 2014
#
#  Author:
#
#    John Burkardt
#
  import numpy as np
  import platform
  from mono_upto_next_grlex import mono_upto_next_grlex

  m = 3
  test_num = 8
  x_test = np.array (  [ \
    0, 0, 0, \
    1, 0, 0, \
    0, 0, 1, \
    0, 2, 0, \
    1, 0, 2, \
    0, 3, 1, \
    3, 2, 1, \
    5, 2, 1 ], dtype = np.int32 )

  print ( '' )
  print ( 'MONO_RANK_GRLEX_TEST' )
  print ( '  Python version: %s' % ( platform.python_version ( ) ) )
  print ( '  MONO_RANK_GRLEX returns the rank of a monomial in the sequence' )
  print ( '  of all monomials in M dimensions, in grlex order.' )

  print ( '' )
  print ( '  Print a monomial sequence with ranks assigned.' )

  n = 4

  print ( '' )
  print ( '  Let M = %d' % ( m ) )
  print ( '      N = %d' % ( n ) )
  print ( '' )

  x = np.zeros ( m, dtype = np.int32 )

  x[0] = 0
  x[1] = 0
  x[2] = 0

  i = 1

  while ( True ):
    print ( '  %2d    ' % ( i ), end = '' )
    for j in range ( 0, m ):
      print ( '%2d' % ( x[j] ), end = '' )
    print ( '' )

    if ( x[0] == n ):
      break

    mono_upto_next_grlex ( m, n, x )
    i = i + 1

  print ( '' )
  print ( '  Now, given a monomial, retrieve its rank in the sequence:' )
  print ( '' )

  for test in range ( 0, test_num ):
    for j in range ( 0, m ):
      x[j] = x_test[j+test*m]
    rank = mono_rank_grlex ( m, x )

    print ( '  %3d    ' % ( rank ), end = '' )
    for j in range ( 0, m ):
      print ( '%2d' % ( x[j] ), end = '' )
    print ( '' )
#
#  Terminate.
#
  print ( '' )
  print ( 'MONO_RANK_GRLEX_TEST' )
  print ( '  Normal end of execution.' )
  return

if ( __name__ == '__main__' ):
  from timestamp import timestamp
  timestamp ( )
  mono_rank_grlex_test ( )
  timestamp ( )