/* This code is a bit ugly, but it works.  I would write it much cleaner with
	C++ if I were doing it today.
*/


#include <stdio.h>

#define NULL_VAL -1
#define BUST_PROB (sp[s[0]][s[1]][s[2]][0].num / 1296.0)

int i,j;
int d[4];
int s[3];
double temp_val, max_val;
long num_leaves;

typedef struct np {
   int num;
   double val;
} NP;
NP sp[13][13][13][11];

double stop_val[13][13][13];
double exp_val[13][13][13];
double app_exp_val[13][13][13];

double payoff[13];

main()
{
    init_stuff();
    get_probs();
    sort_probs();
    print_all_probs();
/*
    print_bust_probs();
*/
    print_bust_probs_table();
    compute_stop_vals();
/* 
    print_stop_vals();
*/
    print_stop_vals_table();
    compute_exp_vals(); 
/*
    print_exp_vals();
*/
    print_exp_vals_table();
/*
    compute_app_exp_vals();
    print_app_exp_vals();
    print_app_exp_vals_table();
    print_exp_vals_diffs_table();
*/
}

init_stuff()
{
  for( s[0] = 2      ; s[0] <= 10 ; s[0]++ )
  for( s[1] = s[0]+1 ; s[1] <= 11 ; s[1]++ )
  for( s[2] = s[1]+1 ; s[2] <= 12 ; s[2]++ )
    for(i=0 ; i < 11 ; i++)
      sp[s[0]][s[1]][s[2]][i].val = NULL_VAL;

  payoff[2] = 1.0 / 3;
  payoff[3] = 1.0 / 5;
  payoff[4] = 1.0 / 7;
  payoff[5] = 1.0 / 9;
  payoff[6] = 1.0 / 11;
  payoff[7] = 1.0 / 13;
  payoff[8] = 1.0 / 11;
  payoff[9] = 1.0 / 9;
  payoff[10] = 1.0 / 7;
  payoff[11] = 1.0 / 5;
  payoff[12] = 1.0 / 3;
}

get_probs()
{
  for( d[0] = 1 ; d[0] <= 6 ; d[0]++)
  for( d[1] = 1 ; d[1] <= 6 ; d[1]++)
  for( d[2] = 1 ; d[2] <= 6 ; d[2]++)
  for( d[3] = 1 ; d[3] <= 6 ; d[3]++)
  {
    for( s[0] = 2      ; s[0] <= 10 ; s[0]++ )
    for( s[1] = s[0]+1 ; s[1] <= 11 ; s[1]++ )
    for( s[2] = s[1]+1 ; s[2] <= 12 ; s[2]++ )
    {
      max_val = 0;
      check_sums(d[0]+d[1],d[2]+d[3]);
      check_sums(d[0]+d[2],d[3]+d[1]);
      check_sums(d[0]+d[3],d[1]+d[2]);
      update_sp();
    }
  }
}

check_sums(int s0, int s1)
{
  temp_val = 0;
  for(i=0 ; i <= 2 ; i++)
  {
    if(s0 == s[i])
      temp_val += payoff[s[i]];
    if(s1 == s[i])
      temp_val += payoff[s[i]];
  }

  if(temp_val > max_val)
    max_val = temp_val;
}

update_sp()
{
  for(i=0 ; sp[s[0]][s[1]][s[2]][i].val != NULL_VAL &&
            sp[s[0]][s[1]][s[2]][i].val != max_val ; i++) ;

  if(sp[s[0]][s[1]][s[2]][i].val == max_val)
    sp[s[0]][s[1]][s[2]][i].num++;
  else
  {
    sp[s[0]][s[1]][s[2]][i].num = 1;
    sp[s[0]][s[1]][s[2]][i].val = max_val;
  }
}

sort_probs()
{  
    NP temp_sp;

    for( s[0] = 2      ; s[0] <= 10 ; s[0]++ )
    for( s[1] = s[0]+1 ; s[1] <= 11 ; s[1]++ )
    for( s[2] = s[1]+1 ; s[2] <= 12 ; s[2]++ )
    {
      for(i=0 ; sp[s[0]][s[1]][s[2]][i+1].val != NULL_VAL ; i++) 
      for(j=0 ; sp[s[0]][s[1]][s[2]][j+1].val != NULL_VAL ; j++) 
      {
        if(sp[s[0]][s[1]][s[2]][j].val > sp[s[0]][s[1]][s[2]][j+1].val)
        {
          temp_sp = sp[s[0]][s[1]][s[2]][j+1];
  	  sp[s[0]][s[1]][s[2]][j+1] = sp[s[0]][s[1]][s[2]][j];
          sp[s[0]][s[1]][s[2]][j] = temp_sp;
        }
      }
    }
}

print_bust_probs()
{  
    for( s[0] = 2      ; s[0] <= 10 ; s[0]++ )
    for( s[1] = s[0]+1 ; s[1] <= 11 ; s[1]++ )
    for( s[2] = s[1]+1 ; s[2] <= 12 ; s[2]++ )
    {
      printf("Bust for %2d %2d %2d is %f\n",s[0],s[1],s[2], BUST_PROB );
    }
}

print_bust_probs_table()
{  
    printf("Bust Percent Values Table\n\n");

    printf("     ");
    for(i=4 ; i <= 12; i++)
      printf("%4d ",i);
    printf("\n");

    for( s[0] = 2      ; s[0] <= 10 ; s[0]++ )
    for( s[1] = s[0]+1 ; s[1] <= 11 ; s[1]++ )
    {
      printf("%2d %2d ",s[0],s[1]); 
      for(i=0 ; i < 5 * (s[1] - 3) ; i++)  printf(" ");
      for( s[2] = s[1]+1 ; s[2] <= 12 ; s[2]++ )
        printf("%2.1f ", 100.0 * BUST_PROB );
      printf(" %2d %2d ",s[0],s[1]); 
      printf("\n");
    }

    printf("     ");
    for(i=4 ; i <= 12; i++)
      printf("%4d ",i);
    printf("\n");

    printf("\n\nExpected value for a roll = (Stop_value - Current_value) * Bust_probability\n\n");
}

print_all_probs()
{  
    for( s[0] = 2      ; s[0] <= 10 ; s[0]++ )
    for( s[1] = s[0]+1 ; s[1] <= 11 ; s[1]++ )
    for( s[2] = s[1]+1 ; s[2] <= 12 ; s[2]++ )
    {
      printf("Probabilities for %2d %2d %2d\n",s[0],s[1],s[2]);
      for(i=0 ; sp[s[0]][s[1]][s[2]][i].val != NULL_VAL ; i++) 
        printf(" %f occurs %3d times, so probability is %f\n",
		sp[s[0]][s[1]][s[2]][i].val,
		sp[s[0]][s[1]][s[2]][i].num,
		sp[s[0]][s[1]][s[2]][i].num / 1296.0 );
    }
}

compute_stop_vals()
{  
    for( s[0] = 2      ; s[0] <= 10 ; s[0]++ )
    for( s[1] = s[0]+1 ; s[1] <= 11 ; s[1]++ )
    for( s[2] = s[1]+1 ; s[2] <= 12 ; s[2]++ )
    {
      temp_val = 0;
      for(i=0 ; sp[s[0]][s[1]][s[2]][i].val != NULL_VAL ; i++) 
      {
        temp_val += sp[s[0]][s[1]][s[2]][i].val * 
		    sp[s[0]][s[1]][s[2]][i].num / 1296.0;
      }
      
      stop_val[s[0]][s[1]][s[2]] = 
		temp_val / BUST_PROB;
    }
}

print_stop_vals_table()
{  
    printf("Stop Values Table\n\n");

    printf("     ");
    for(i=4 ; i <= 12; i++)
      printf("%4d ",i);
    printf("\n");

    for( s[0] = 2      ; s[0] <= 10 ; s[0]++ )
    for( s[1] = s[0]+1 ; s[1] <= 11 ; s[1]++ )
    {
      printf("%2d %2d ",s[0],s[1]); 
      for(i=0 ; i < 5 * (s[1] - 3) ; i++)  printf(" ");
      for( s[2] = s[1]+1 ; s[2] <= 12 ; s[2]++ )
        printf("%1.2f ",stop_val[s[0]][s[1]][s[2]]);
      printf(" %2d %2d ",s[0],s[1]); 
      printf("\n");
    }

    printf("     ");
    for(i=4 ; i <= 12; i++)
      printf("%4d ",i);
    printf("\n");

    printf("\n\nExpected value for a roll = (Stop_value - Current_value) * Bust_probability\n\n");
}

print_stop_vals()
{  
    for( s[0] = 2      ; s[0] <= 10 ; s[0]++ )
    for( s[1] = s[0]+1 ; s[1] <= 11 ; s[1]++ )
    for( s[2] = s[1]+1 ; s[2] <= 12 ; s[2]++ )
    {
      printf("Stop at %f for %2d %2d %2d\n",
      		stop_val[s[0]][s[1]][s[2]], s[0],s[1],s[2]);
    }
}

expand_cur_val(double cur_prob, double cur_val)
{
  int k;

  if(cur_val > stop_val[s[0]][s[1]][s[2]] )
  {
    exp_val[s[0]][s[1]][s[2]] += cur_val * cur_prob;
    num_leaves++;
  }
  else
    for(k=1 ; sp[s[0]][s[1]][s[2]][k].val != NULL_VAL ; k++) 
      expand_cur_val(cur_prob * sp[s[0]][s[1]][s[2]][k].num / 1296.0,
		     cur_val  + sp[s[0]][s[1]][s[2]][k].val);
}

compute_exp_vals()
{  
    for( s[0] = 2      ; s[0] <= 10 ; s[0]++ )
    for( s[1] = s[0]+1 ; s[1] <= 11 ; s[1]++ )
    for( s[2] = s[1]+1 ; s[2] <= 12 ; s[2]++ )
    {
      num_leaves=0;
      exp_val[s[0]][s[1]][s[2]] = 0;
      expand_cur_val(1.0,payoff[s[0]]+payoff[s[1]]+payoff[s[2]]);

      printf("Expected value is %f for %2d %2d %2d with %d leaves\n",
      		exp_val[s[0]][s[1]][s[2]], s[0],s[1],s[2],num_leaves);
      fflush(stdout);
    }
}

print_exp_vals()
{  
    for( s[0] = 2      ; s[0] <= 10 ; s[0]++ )
    for( s[1] = s[0]+1 ; s[1] <= 11 ; s[1]++ )
    for( s[2] = s[1]+1 ; s[2] <= 12 ; s[2]++ )
    {
      printf("Expected value is %f for %2d %2d %2d\n",
      		exp_val[s[0]][s[1]][s[2]], s[0],s[1],s[2]);
    }
}

print_exp_vals_table()
{  
    printf("Expected Values Table\n\n");

    printf("     ");
    for(i=4 ; i <= 12; i++)
      printf("%4d ",i);
    printf("\n");

    for( s[0] = 2      ; s[0] <= 10 ; s[0]++ )
    for( s[1] = s[0]+1 ; s[1] <= 11 ; s[1]++ )
    {
      printf("%2d %2d ",s[0],s[1]); 
      for(i=0 ; i < 5 * (s[1] - 3) ; i++)  printf(" ");
      for( s[2] = s[1]+1 ; s[2] <= 12 ; s[2]++ )
        printf("%1.2f ",exp_val[s[0]][s[1]][s[2]]);
      printf(" %2d %2d ",s[0],s[1]); 
      printf("\n");
    }

    printf("     ");
    for(i=4 ; i <= 12; i++)
      printf("%4d ",i);
    printf("\n");
}

double power(double x, int n)
{
  int e;
  double product=1.0;

  for(e=1 ; e <= n ; e++ )
    product *= x;

  return(product);
}
  
double app_exp_n_turns(int n)
{ 
    return( power(1.0 - BUST_PROB, n) *
            ( payoff[s[0]]+payoff[s[1]]+payoff[s[2]] +
	      n * stop_val[s[0]][s[1]][s[2]] * 
	          BUST_PROB /
	          (1.0 - BUST_PROB)
	    )
          );
}
		  
compute_app_exp_vals()
{  
    for( s[0] = 2      ; s[0] <= 10 ; s[0]++ )
    for( s[1] = s[0]+1 ; s[1] <= 11 ; s[1]++ )
    for( s[2] = s[1]+1 ; s[2] <= 12 ; s[2]++ )
    {
      for(i=0 ; ; i++) 
      {
        if(app_exp_n_turns(i) > app_exp_n_turns(i+1))
        {
          app_exp_val[s[0]][s[1]][s[2]] = app_exp_n_turns(i);
	  break;
	}
      }
    }
}

print_app_exp_vals()
{  
    for( s[0] = 2      ; s[0] <= 10 ; s[0]++ )
    for( s[1] = s[0]+1 ; s[1] <= 11 ; s[1]++ )
    for( s[2] = s[1]+1 ; s[2] <= 12 ; s[2]++ )
    {
      printf("Approximate Expected value is %f for %2d %2d %2d\n",
      		app_exp_val[s[0]][s[1]][s[2]], s[0],s[1],s[2]);
    }
}

print_app_exp_vals_table()
{  
    printf("Approximate Expected Values Table\n\n");

    printf("     ");
    for(i=4 ; i <= 12; i++)
      printf("%4d ",i);
    printf("\n");

    for( s[0] = 2      ; s[0] <= 10 ; s[0]++ )
    for( s[1] = s[0]+1 ; s[1] <= 11 ; s[1]++ )
    {
      printf("%2d %2d ",s[0],s[1]); 
      for(i=0 ; i < 5 * (s[1] - 3) ; i++)  printf(" ");
      for( s[2] = s[1]+1 ; s[2] <= 12 ; s[2]++ )
        printf("%1.2f ",app_exp_val[s[0]][s[1]][s[2]]);
      printf(" %2d %2d ",s[0],s[1]); 
      printf("\n");
    }

    printf("     ");
    for(i=4 ; i <= 12; i++)
      printf("%4d ",i);
    printf("\n");
}

print_exp_vals_diffs_table()
{  
    printf("Diffs Expected Values Table\n\n");

    printf("     ");
    for(i=4 ; i <= 12; i++)
      printf("%4d ",i);
    printf("\n");

    for( s[0] = 2      ; s[0] <= 10 ; s[0]++ )
    for( s[1] = s[0]+1 ; s[1] <= 11 ; s[1]++ )
    {
      printf("%2d %2d ",s[0],s[1]); 
      for(i=0 ; i < 5 * (s[1] - 3) ; i++)  printf(" ");
      for( s[2] = s[1]+1 ; s[2] <= 12 ; s[2]++ )
        printf("%0.3f ",exp_val[s[0]][s[1]][s[2]] -
                        app_exp_val[s[0]][s[1]][s[2]]);
      printf(" %2d %2d ",s[0],s[1]); 
      printf("\n");
    }

    printf("     ");
    for(i=4 ; i <= 12; i++)
      printf("%4d ",i);
    printf("\n");
}