/* File CExamples/Recursion/combin3.cpp */
#include <map>
#include <iterator>
#include <stdio.h>

using namespace std;

long int combin(int i, int k)
/* Binomial coefficient of i choose k;
   Note: (i choose k) = 0 for k < 0,
                      = 1 for k == 0,
                      = pow(-1, k) (k-i-1 choose k) for i < 0, k > 0
                      = (i * (i-1) * ... * (i-k+1))/
                        (k * (k-1) * ... * (1))     for i >= 0, k > 0 */
/* This function remembers previous calls */
{long int combinat; /* intermediate function result */
#define MAXARG 50
 typedef map< pair<int, int>, long int > memo_map_t;
 static memo_map_t memo;
 static memo_map_t::iterator memo_p;

#ifdef DEBUG
   /* code within #ifdef DEBUG prints the trace */
   /* compile with -DDEBUG to get trace */
   static int invoke = 0;  /* level of invocation */
   static int invblanks;  /* loop index in debug statement */
   invoke = invoke+1;  /* raise level of invocation */
   for (invblanks = 1; invblanks <= invoke-1; invblanks++)
       printf("%1d   ", invblanks % 10); /* last digit only */
   /* indent properly, writing level numbers also */
   if(memo.find(make_pair(i,k)) != memo.end())// indicate if remembered
   printf("combin(%3d,%3d) found\n",i,k);  /* indicate call */
   else
   printf("combin(%3d,%3d)\n",i,k);  /* indicate call */
#endif

  /* Check arguments */
  if (k < 0)  combinat = 0L;
  else if (k == 0) combinat = 1L;
  /* Now k > 0 */
  else if (i < 0)
         {if (k & 1) /* test lowest order bit of k */
               /* k is odd */  combinat = -combin(k-i-1, k);
          else /* k is even */ combinat = combin(k-i-1, k);
         }
  /* Now k > 0 and i >= 0 */
  else if (i == 0 || k > i) combinat = 0L;
  /* Now i > 0 and 0 < k <= i */
  else if (i == k) combinat = 1L;
  /* Now i > 0 and 0 < k < i */
  else if (k == 1 || k == i-1) combinat = static_cast<long int>(i);
  /* Now i > 0 and 1 < k < i-1 */
  else /* check if already computed */
      {memo_p = memo.find( make_pair(i, k) );
       if(memo_p != memo.end())
         {
          // already computed: do not recur
          combinat = memo_p->second; // value of map
         }
         else /* compute and memorize */
         {
          combinat = memo[make_pair(i, k)] = combin(i-1,k-1)+combin(i-1,k);
         }
      }//

#ifdef DEBUG
    for(invblanks=1; invblanks<=invoke-1; invblanks++)
       printf("%1d   ",invblanks % 10);  
    /* indent properly, writing level numbers also */
    printf("value=%ld\n",combinat); 
    invoke=invoke-1;  /* lower level of invocation before return */
#endif

   return(combinat);
}/*combin*/