/* Find permutations (of the integers from 0 to SIZE - 1)
   containing a minimum number of monotone subsequences
   of length SEQLEN.

   SIZE and SEQLEN must be defined as macros
   on the compiler command line.  */
/* CVS: $Id: monotone-subseq.c,v 1.1 2002/05/19 15:23:59 jsm28 Exp $ */
/* Copyright 2002 Joseph Samuel Myers.
   All rights reserved.
  
   Redistribution and use in source and binary forms,
   with or without modification, are permitted provided that
   the following conditions are met:
   1. Redistributions of source code must retain the above
      copyright notice, this list of conditions and the
      following disclaimer.
   2. Redistributions in binary form must reproduce the above
      copyright notice, this list of conditions and the
      following disclaimer in the documentation and/or
      other materials provided with the distribution.
   3. The name of the author may not be used to endorse
      or promote products derived from this software without
      specific prior written permission.
  
   THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND
   ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING,
   BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF
   MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
   ARE DISCLAIMED.  IN NO EVENT SHALL THE AUTHOR BE LIABLE
   FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY,
   OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
   PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
   DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED
   AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT,
   STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR
   OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
   SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*/

#include <stdarg.h>
#include <stdbool.h>
#include <stdint.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>

/* Print an error and exit.  */
static void
die(const char *format, ...)
{
  va_list args;
  fprintf(stderr, "monotone-subseq: ");
  va_start(args, format);
  vfprintf(stderr, format, args);
  va_end(args);
  exit(EXIT_FAILURE);
}

/* Compute the binomial coefficient (n choose k).  */
static int
choose(int n, int k)
{
  int r = 1;
  for (int j = 1; j <= k; j++) {
    r *= n;
    r /= j;
    n--;
  }
  return r;
}

/* Compute the expected minimum number of monotone subsequences
   of length seqlen in a permutation of length size.  */
static int
expected_minl(int size, int seqlen)
{
  int segs = seqlen - 1;
  int sseg = size / segs;
  int lseg = sseg + 1;
  int lct = size % segs;
  int sct = segs - lct;
  return (sct * choose(sseg, seqlen)
          + lct * choose(lseg, seqlen));
}

/* Count of number of extremal subsequences found.  */
static uintmax_t count = 0;

/* Count of number of extremal subsequences found
   which have both ascending and descending subsequences
   of length SEQLEN.  */
static uintmax_t scount = 0;

/* True if a count of the sequences found is to be printed.  */
bool count_flag = false;

/* True if the sequences found are to be printed.  */
bool print_flag = false;

/* Print a sequence.  */
static void
print_sequence(const int seq[SIZE])
{
  for (int i = 0; i < SIZE; i++)
    printf(" %d", seq[i]);
  printf("\n");
}

/* Handle an extremal sequence.  */
static void
handle_extremum(const int seq[SIZE], int atotal, int dtotal)
{
  count++;
  if (count == 0)
    die("count overflow");
  if (atotal != 0 && dtotal != 0)
    scount++;
  if (print_flag)
    print_sequence(seq);
}

/* The following command line options are accepted:

   -c   count the number of extrema
   -p   print the extrema found

   Both default to off.  Any extrema found that provide
   counterexamples to the conjectured minimum number of
   monotone subsequences are still printed.

   Apart from these options, two equal-length sequences
   of numbers may be given; these represent the start of
   the first permutation to be considered, and the start
   of the first permutation not to be considered.  This is
   for use in dividing up the search into parts that are
   run in parallel.  If no numbers are given as arguments,
   all permutations are searched.  */

int
main(int argc, char **argv)
{
  /* The permutation under test.  */
  int perm[SIZE];
  /* The numbers of ascending and descending subsequences
     of lengths 1 to SEQLEN ending at each position.  */
  int acounts[SIZE][SEQLEN];
  int dcounts[SIZE][SEQLEN];
  /* next[0] is the least value not in the permutation.
     next[1 + i] is the least value greater than i not
     in the permutation up to and including the place, if
     any, where i is in the permutation.  */
  int next[SIZE + 1];
  /* was[i] is the previous value that was in position i,
     the value in that position having been increased at
     least once while the previous values remained fixed,
     or -1 if the value in position i is the first possible
     value there given the previous values.  */
  int was[SIZE];
  /* The total number of ascending subsequences of length
     SEQLEN ending at or before each position.  */
  int atotals[SIZE];
  /* The total number descending subsequences of length
     SEQLEN ending at or before each position.  */
  int dtotals[SIZE];
  /* The expected minimum number of monotone subsequences
     of length SEQLEN, or the actually found minimum
     if smaller.  */
  int minl = expected_minl(SIZE, SEQLEN);
  /* The current position within the permutation.  */
  int pos = 0;
  /* The permutations at which to start and end.  */
  int chkstart[SIZE];
  int chkend[SIZE];
  /* The number of values in those permutations,
     or -1 if all permutations are to be checked.  */
  int chkct;
  argv++;
  argc--;
  while (argc > 0 && argv[0][0] == '-') {
    if (!strcmp(argv[0], "-c"))
      count_flag = true;
    else if (!strcmp(argv[0], "-p"))
      print_flag = true;
    else
      die("bad argument \"%s\"", argv[0]);
    argv++;
    argc--;
  }
  if (argc % 2 != 0)
    die("bad number of arguments");
  chkct = argc / 2;
  if (chkct >= SIZE)
    die("too many arguments");
  if (chkct > 0) {
    for (int i = 0; i < chkct; i++)
      chkstart[i] = atoi(argv[i]);
    for (int i = 0; i < chkct; i++)
      chkend[i] = atoi(argv[chkct + i]);
    printf("starting at");
    for (int i = 0; i < chkct; i++)
      printf(" %d", chkstart[i]);
    printf("\nending at");
    for (int i = 0; i < chkct; i++)
      printf(" %d", chkend[i]);
    printf("\n");
  } else {
    chkct = -1;
  }
  for (int i = 0; i <= SIZE; i++)
    next[i] = i;
  for (int i = 0; i < SIZE; i++)
    perm[i] = -1;
  /* Fill in some positions with our starting position.
     This uses much the same code as the main loop.  */
  for (pos = 0; pos < chkct; pos++) {
    while (perm[pos] != chkstart[pos]) {
      int j = perm[pos];
      int nx = next[1 + j];
      if (j != -1) {
        int w = was[pos];
        next[1 + w] = j;
      }
      perm[pos] = nx;
      was[pos] = j;
      next[1 + j] = next[1 + nx];
    }
    acounts[pos][0] = 1;
    dcounts[pos][0] = 1;
    for (int k = 1; k < SEQLEN; k++) {
      int jacount = 0, jdcount = 0;
      /* Compute numbers of (k+1)-sequences.  */
      for (int l = k - 1; l < pos; l++) {
        if (perm[l] < perm[pos])
          jacount += acounts[l][k - 1];
        if (perm[l] > perm[pos])
          jdcount += dcounts[l][k - 1];
      }
      acounts[pos][k] = jacount;
      dcounts[pos][k] = jdcount;
    }
    atotals[pos] = ((pos > 0 ? atotals[pos - 1] : 0)
                    + acounts[pos][SEQLEN - 1]);
    dtotals[pos] = ((pos > 0 ? dtotals[pos - 1] : 0)
                    + dcounts[pos][SEQLEN - 1]);
  }
  while (1) {
    /* Fill in position pos.  */
    if (pos == SIZE) {
      int total = atotals[pos - 1] + dtotals[pos - 1];
      if (total <= minl) {
        if (total < minl) {
          count = 0;
          scount = 0;
          minl = total;
          printf("Bettered expected result "
                 "(%d instead of %d): ", minl,
                 expected_minl(SIZE, SEQLEN));
          print_sequence(perm);
        }
        handle_extremum(perm, atotals[pos - 1],
                        dtotals[pos - 1]);
      }
      pos--;
      continue;
    }
    int j = perm[pos];
    int nx = next[1 + j];
    if (j != -1) {
      int w = was[pos];
      next[1 + w] = j;
    }
    if (pos == 0 && nx >= SIZE)
      break;
    if (nx == SIZE) {
      perm[pos] = -1;
      pos--;
      continue;
    }
    perm[pos] = nx;
    was[pos] = j;
    next[1 + j] = next[1 + nx];
    if (pos == chkct
        && !memcmp(chkend, perm, chkct * sizeof(int)))
      break;
    acounts[pos][0] = 1;
    dcounts[pos][0] = 1;
    for (int k = 1; k < SEQLEN; k++) {
      int jacount = 0, jdcount = 0;
      /* Compute numbers of (k+1)-sequences.  */
      for (int l = k - 1; l < pos; l++) {
        if (perm[l] < perm[pos])
          jacount += acounts[l][k - 1];
        if (perm[l] > perm[pos])
          jdcount += dcounts[l][k - 1];
      }
      acounts[pos][k] = jacount;
      dcounts[pos][k] = jdcount;
    }
    atotals[pos] = ((pos > 0 ? atotals[pos - 1] : 0)
                    + acounts[pos][SEQLEN - 1]);
    dtotals[pos] = ((pos > 0 ? dtotals[pos - 1] : 0)
                    + dcounts[pos][SEQLEN - 1]);
    int addv = 0;
    int total = atotals[pos] + dtotals[pos];
    if (pos < SIZE - 1) {
      /* Find the least and greatest values yet to be placed,
         and see if either of them would force the number of
         monotone subsequences to be too large.  */
      int lv, gv;
      lv = next[0];
      int daddv = 0;
      for (int k = 0; k <= pos; k++) {
        if (perm[k] > lv)
          daddv += dcounts[k][SEQLEN - 2];
      }
      if (total + daddv > minl)
        goto no_increase_pos;
      gv = lv;
      while (next[1 + gv] != SIZE)
        gv = next[1 + gv];
      for (int k = 0; k <= pos; k++) {
        if (perm[k] < gv)
          addv += acounts[k][SEQLEN - 2];
      }
    }
    if (total + addv <= minl)
      pos++;
  no_increase_pos: ;
  }
  if (count_flag) {
    printf("Count %ju\n", count);
    printf("Special count %ju\n", scount);
  }
  exit(EXIT_SUCCESS);
}