// Example solution to Assignment 7.
// CS 509 - Advanced Programming II, Fall 2002

// A little bit of fiddling with a paper and pencil is enough to show that
// repeated halving produces an additive up-sequence of O(log n) elements,
// which is pretty small.  However, it's not at all clear (at least it wasn't
// to me) that it produces a minimal up-sequence.  A little bit more fiddling
// with a pencil and paper shows that it doesn't.
//
// Repeated halving adds the largest number of numbers to an up-sequence when
// the result of halving is odd, in which case n/2 and (n/2) + 1 are added.
// For example, when n = 15, repeated halving gives the up-sequence 15, 8, 7,
// 4, 3, 2, 1.  Unfortunately, the up-sequence 15, 10, 5, 3, 2, 1 has one less
// element, and so the halving up-sequence isn't minimal.
//
// Now we have two questions: What is the pattern behind the smaller
// up-sequence?, and Is the smaller up-sequence minimal, or is there an even
// smaller one?  The pattern appears to be based on doubling, which makes sense
// because doubling adds the smallest number of elements possible to the
// up-sequence.  Unfortunately, the doubling pattern isn't clear (at least not
// to me it isn't), and if the pattern isn't clear, then neither is the agument
// that it leads to a minimal up-sequence.
//
// Fortunately, even though intelligence may fail, brute force always works.
// The up-sequence for any valid n ends with 2, 1 (why 2?)  Given the partial
// up-sequence 2, 1, the next number in the up-sequence must be the sum of two
// of the existing numbers: either 4 (= 2 + 2) or 3 (= 2 + 1) (notice that 2 =
// 1 + 1) isn't a valid candidate because the next number in the up-sequence
// must be strictly larger than any of the existing numbers).
//
// These observations give us the general outline of the brute-force search:
// take a partial up-sequence and find all possible next numbers for the
// up-sequence by summing all existing-number pairs.  If the possible next
// numbers contains n, then the up-sequence is done.  Otherwise, add each
// possible next number to the partial up-sequence to form a new partial
// up-sequence, and recurse the search with the new partial up-sequence.
//
// This is a ridiculously expensive exponential search, but it's guaranteed to
// find a minimal up-sequence because it looks at all possible up-sequences.
// There are also some opportunities for pruning the search - such as limiting
// the next possible numbers to fall between n and the largest number in the
// partial up-sequence - but at heart the search will still be exponential.

#include <algorithm>
#include <cstdio>
#include <iostream>
#include <string>
#include <sys/time.h>
#include <unistd.h>
#include <vector>
#include <strings.h>
#include <cstdlib>
#include <iterator>
#include <cassert>
#include <set>
#include <iomanip>
#include <sstream>
#include <cmath>
#include "mau.h"


// Non-negative integers.

   typedef unsigned nat;

// An up-sequence is a vector of naturals.

   typedef std::vector<nat>       natvec;
   typedef natvec::iterator       nvec_iter;
   typedef natvec::const_iterator nvec_citer;

// The set of next-possible numbers is a set of naturals.

   typedef std::set<nat>                  natset;
   typedef natset::iterator               nset_iter;
   typedef natset::const_iterator         nset_citer;
   typedef natset::const_reverse_iterator nset_criter;

// The function that finds minimal additive upsequences.

   typedef natvec (* mau_fun)(nat);

// Deja vu, c++ style.

   static void check_times(time_t, time_t);
   static void run_test(bool, nat);
   static void print_time(time_t);
   static void print_times(const char * const, time_t, time_t);
   static natvec timed_run(mau_fun, nat, time_t &);
   std::ostream & operator << (std::ostream &, const natvec &);
   static void print_emsg(const char * const);
   static natset next_numbers(const natvec &, nat);


static void
check_answer(const natvec & mau, nat n) {

  // Determine if mau is a valid up-sequence to n; die with an error message if
  // it isn't.

  const nat s = 100;
  char emsg[s];

  if ((mau.size() == 0) and (n < 2))
    return;

  if (n < 2) {
    snprintf(emsg, s, "should contain no numbers, not %u", mau.size());
    print_emsg(emsg);
    }

  if (mau.size() < 2) {
    snprintf(emsg, s,
	     "should contain at least two numbers, not %u", mau.size());
    print_emsg(emsg);
    }

  if (mau.front() != 1)
    print_emsg("doesn't start with 1");

  if (mau.back() != n) {
    snprintf(emsg, s, "doesn't end with %u", n);
    print_emsg(emsg);
    }

  nat i;
  for (i = mau.size() - 1; i > 0; i--)
    if (mau[i - 1] >= mau[i])
      print_emsg("isn't in strictly increasing order");

  for (nvec_citer i = mau.end() - 1; i != mau.begin(); --i) {
    bool found;
    nvec_citer j = i - 1;
    const nvec_citer b = mau.begin();

    do found = std::find(b, j + 1, *i - *j) != (j + 1);
    while (!found and (j-- != b));

    if (!found) {
      snprintf(emsg, s, 
	       "value %u isn't the sum of two previous up-sequence values",
	       *i);
      print_emsg(emsg);
      }
    }
  }


static nat
check_mau(natvec candidate, natvec & answer, nat n, nat max_cnt) {

  // Extend the partial up-sequence candiate by one element.  answer contains
  // the current best up-sequence, which has max_cnt elements; n is the largest
  // element in the up-sequence.

  // Find all the possible next elements for candidate.

     natset nexts = next_numbers(candidate, n);

  // Make space for the new sum.  For this call to be useful, candidate should
  // still be smaller than max-cnt after the new element is added.

     candidate.push_back(0);
     const nat candidate_size = candidate.size();
     assert(candidate_size < max_cnt);

  if (nexts.count(n)) {

    // The possible next elements contains the largest element in the
    // up-sequence, which makes candidate the new best answer.

       candidate.back() = n;
       std::swap(answer, candidate);

    }
  else

    // Try each possible next element as long as the partial up-sequence
    // remains small enough to form a new best answer.  Notice the next
    // elements are taken in decreasing order, which forces failures to occur
    // earlier than they would if the next elements were taken in increasing
    // order.

       for (nset_criter i = nexts.rbegin();
	    (i != nexts.rend()) and (candidate_size + 1 < max_cnt); 
	    i++) {
	 candidate.back() = *i;
	 max_cnt = check_mau(candidate, answer, n, max_cnt);
	 }

  return answer.size();
  }


static natvec
check_mau(nat n) {
  
  // Return the minimal additive up-sequence having n as its largest element.

  natvec answer;

  if (n > 1) {
    const nat m[] = { 1, 2 };
    std::copy(m, m + 2, back_inserter(answer));

    if (n > 2) {
      nat max_cnt = 
	2*static_cast<nat>(ceil(log(static_cast<double>(n))/log(2.0)));
      check_mau(answer, answer, n, max_cnt);
      assert(answer.size() <= max_cnt);
      }

    assert(answer.front() == 1);
    assert(answer.back() == n);
    }

  return answer;
  }


static void
check_times(time_t answert, time_t checkt) { 

  // Compare an answer time to a check time and print a message if the answer
  // time's too slow compared to the check time.  Return true iff the time's
  // are good.

  if (answert > 2*checkt) {
    print_times("\n!!! ", answert, checkt);
    exit(EXIT_FAILURE);
    }
  }


int
main(int argc, char * argv[]) {
  
  bool telltime = false;
  nat seed = getpid();
  int n = -1;

  int c;
  extern char *optarg;
  extern int optind;

  while ((c = getopt(argc, argv, "ts:n:")) != EOF)
     switch (c) {
     case 't':
       telltime = true;
       break;

     case 's':
       seed = atoi(optarg);
       break;

     case 'n': {
       int i = atoi(optarg);
       if (i < 0) {
	 std::cerr << "The -n option should be positive.\n";
	 exit(EXIT_FAILURE);
	 }
       n = static_cast<nat>(i);
       break;
       }

     case '?':
       std::cerr << "Command format is \"" << argv[0] 
		 << " [-t] [-sn | -nn]\".\n";
       return EXIT_FAILURE;

     default:
       assert(!"getopt() went berserk");
     }

  if (optind != argc) {
    std::cerr << "Command format is \"" << argv[0] 
	      << " [-t] [-sn | -nn]\".\n";
    return EXIT_FAILURE;
    }

  if (n < 0) {
    srandom(seed);

    check_mau(16);
    for (nat i = 0; i < 5; i++)
      run_test(telltime,  random() % 300);
    }
  else 
    run_test(telltime, static_cast<nat>(n));
  }


std::ostream & 
operator << (std::ostream & os, const natvec & iv) {

  // Print iv to the output stream os, returning os.

  for (nvec_citer i = iv.begin(); i != iv.end(); i++)
    os << std::setw(4) << *i;

  return os;
  }


static natset
next_numbers(const natvec & mau, nat n) {

  // Given the up-sequence mau and the largest possible value of the
  // up-sequence n, return a set of possible next elements for the up-sequence.

  const nat min_sum = mau.back();
  natset nexts;

  assert(mau.back() < n);

  for (nvec_citer i = mau.begin(); i != mau.end(); i++)
    for (nvec_citer j = i; j != mau.end(); j++) {
      const nat sum = *i + *j;
      if ((min_sum < sum) and (sum <= n))
	nexts.insert(sum);
      }

  return nexts;
  }


static void
print_emsg(const char * const emsg) {

  // Print to std-err the error message given in emsg, then die.

  std::cerr << "\n!!! Computed up-sequence " << emsg << ".\n";
  exit(EXIT_FAILURE);
  }


static void 
print_time(time_t t) {

  // Convert the micosecond time interval t into a conveniently small number
  // and print it to std-err.

  if (t < 1000)
    fprintf(stderr, "%5u usec", static_cast<nat>(t));
  else if (t < 100000)
    fprintf(stderr, "%5.1f msec", static_cast<double>(t)/1000.0);
  else
    fprintf(stderr, "%5.1f  sec", static_cast<double>(t)/1000000.0);
  }


static void
print_times(const char * const prefix, time_t answert, time_t checkt) { 

  // Compare an answer time to a check time and print a message if the answer
  // time's too slow compared to the check time.  Return true iff the time's
  // are good.

  fprintf(stderr, "%sComputed mau() time is", prefix);
  print_time(answert);
  fprintf(stderr, " vs. check mau() time ");
  print_time(checkt);
  fprintf(stderr, " (x%3.1f).\n",
	  static_cast<double>(answert)/static_cast<double>(checkt));
  }


static void
run_test(bool telltime, nat n) {

  // Feed n to the two mau procedures, and compare the results.  If telltime is
  // true, print the execution times for each mau procedure.

  time_t answert, check_answert;
  const natvec
    check_ans = timed_run(check_mau, n, check_answert),
    ans = timed_run(mau, n, answert);

  if (telltime)
    print_times("", answert, check_answert);

  check_answer(ans, n);

  if (ans.size() > check_ans.size()) {
    const nat s = 100;
    char emsg[s];
    snprintf(emsg, s, "has %u numbers, should have no more than %u numbers.\n",
	     ans.size(), check_ans.size());

    std::ostringstream oss;
    oss << "ans   = " << ans << "\n";
    oss << "check = " << check_ans;

    print_emsg((std::string(emsg) + oss.str()).c_str());
    }

  check_times(answert, check_answert);
  }


static natvec
timed_run(mau_fun mf, nat n, time_t & et) {

  // Run the mau procedure mf against n, returning the result returned by mf.
  // Also return in et the execution time in usecs.

  struct timeval s, e;

  gettimeofday(&s, NULL);
  const natvec a = (* mf)(n);
  gettimeofday(&e, NULL);

  if (s.tv_sec < e.tv_sec) { 
    e.tv_sec--;
    e.tv_usec += 1000000; 
    }
  et = (e.tv_sec - s.tv_sec)*1000000 + (e.tv_usec - s.tv_usec);

  return a;
  }


// $Log: main.cc,v $
// Revision 1.6  2002/12/19 17:47:20  rclayton
// Fix the documentation.
//
// Revision 1.5  2002/12/19 00:05:23  rclayton
// Implement the -n option; add brute-force searching; document.
//
// Revision 1.4  2002/12/02 16:29:56  rclayton
// Explictly include <cassert>.
//
// Revision 1.3  2002/12/02 16:04:24  rclayton
// Implement the -t option.
//
// Revision 1.2  2002/11/30 21:54:18  rclayton
// Implement a log(n) mau().
//
// Revision 1.1  2002/11/30 20:26:19  rclayton
// Initial revision
//


syntax highlighted by Code2HTML, v. 0.9