/* civ 3 combat calculator - by vulture (10th Jun 2002)
 *
 * Calculates combat result probabilities between 2 units. If one of the units
 * is able to retreat from combat (fast unit attacking move 1 unit) then
 * that unit is unit 'A'. Change values in the last block of '#define'
 * statements to reflect the values of the units you wish to calculate.
 * Currently set up for veteran cavalry attacking unfortified rifleman on
 * open ground, both with full hit points.
 *
 * You have to work out the terrain modifiers etc. yourself, and use the
 * correct modified attack/defence values for each unit. It doesn't matter
 * whether the attacker is unit A or B, UNLESS ONE OF THE UNITS IS RETREAT
 * CAPABLE - in which case it is unit A regardless of whether it is attacking
 * or defending. For a cavalry being attacked by an infantry on open ground,
 * you would use A_STR 3.3, B_STR 6.0 (A's modified defence, B's attack, in
 * this case).
 *
 * The following assumptions are made (based on memory):
 *
 * When the fast unit is fighting a unit with only 1 hp left it never tries
 * to retreat.
 *
 * The fast unit only gets 1 chance to retreat. E.g. a cavalry attacking
 * a rifleman gets reduced to 1 hp while the rifleman still has 3. The
 * cavalry will try to retreat upon reaching 1 hp. If this fails, the
 * combat continues. If the cavalry wins the next round, then the rifleman
 * is on 2 hp, and the cavalry on 1. I assume that the cavalry doesn't get
 * a second chance to retreat. If it does the calculation gets more complex,
 * and doesn't make a huge difference to the results (no. of wins decreases
 * slightly, no. of retreat increases, no. of losses decreases moderately).
 *
 * Results are hopefully self explanatory */

#include <stdio.h>
#include <math.h>

/* don't change these */
#define NO_RETREAT 0
#define CONSCRIPT 1
#define STANDARD 2
#define VETERAN 3
#define ELITE 4
#define MAXLEVEL 5
#define MAXHP 11 /* actually max. no of hitpoints is 1 less than this */

/* change these if you have altered the retreat chances in the editor */
#define R_CON 0.34
#define R_STD 0.50
#define R_VET 0.58
#define R_ELI 0.66

/* change these for the units you are interested in (or use command line
 * options) */
#define HP_A 4
#define HP_B 3
#define A_STR 6.0
#define B_STR 6.6
#define EXPLEVEL VETERAN

void read_command_line(double *, double *, int *, int *, double *, int *,
		       double [], int, char *[]);
void calculate_probs(double *, double *, double *, double, double, double,
		       double, int, int, double [][]);
double fac(int); /* returns factorial of input value */
void resolve_next_fight(double [], double [], int);

main(int argc, char *argv[])
{
  double a_st = A_STR;
  double b_st = B_STR;
  double a;
  double b;
  int hp_a = HP_A; /* attacker's hit points */
  int hp_b = HP_B; /* defender's hit points */
  int i, j, k; /* loop variables */
  double retreat[MAXLEVEL] = {0.0, R_CON, R_STD, R_VET, R_ELI};
  int level;
  double r = retreat[EXPLEVEL]; /* retreat probability */
  double s; /* probability of not retreating */
  double p_win;
  double p_retreat;
  double p_lose;
  double array[MAXHP][MAXHP];
  double dum_array[MAXHP][MAXHP]; /* not used for anything meaningful */
  double probvec1[MAXHP];
  double probvec2[MAXHP];
  double dievec[MAXHP];
  double p1, p2, p3;
  double av_loss;

  /* this function deals with the command line arguments, and assigns
   * values to some of the variables */
  read_command_line(&a_st, &b_st, &hp_a, &hp_b, &r, &level, retreat, argc,
		    argv);
  a = a_st / (a_st + b_st);  /* probability that A wins 1 round */
  b = 1.0 - a;  /* probability that B wins 1 round */
  s = 1.0 - r; /* probability of not retreating */
  calculate_probs(&p_win, &p_retreat, &p_lose, a, b, r, s, hp_a, hp_b, array);

  printf("\nA:\n  hp: %d  combat strenth: %.2f\n", hp_a, a_st);
  if (r != 0) printf("  retreat chance: %.2f\n", r);
  printf("B:\n  hp: %d  combat strenth: %.2f\n\n", hp_b, b_st);
  printf("Hit points remaining:\nA \\ B");
  for (i = 0; i <= hp_b; i++) printf("%3d  ", i);
  printf("\n");
  for (j = 0; j <= hp_a; j++) {
    printf("%3d  ", j);
    for (i = 0; i <= hp_b; i++) {
      if (array[j][i] != 0.0) printf("%4.1f ", array[j][i] * 100.0);
      else printf("  -  ");
    }
    printf("\n");
  }
  printf("\nA wins:\t\t%4.1f %%\n", p_win * 100.0);
  if (p_retreat != 0.0) printf("A retreats:\t%4.1f %%\n", p_retreat * 100.0);
  printf("B wins:\t\t%4.1f %%\n\n", p_lose * 100.0);

  /* new section - for a given no. of unit A, what are the odds of killing
   * unit B? Calculate this and vice versa, fo no. of units up to 95%
   * confidence of victory */

  /* first for A */
  for (i = 0; i <= MAXHP; i++) {
    probvec1[i] = 0.0;
    probvec2[i] = 0.0;
  }
  for (i = 0; i <= hp_b; i++) {
    for (j = 0; j <= hp_a; j++) probvec1[hp_b - i] += array[j][i];
    probvec2[hp_b - i] = probvec1[hp_b - i];
  }
  av_loss = p_lose;
  printf("Victory chance for n of unit A vs 1 unit B:\n(with expected "
	 " no. of A units killed)\n n   prob    losses\n");
  k = 1;
  while ((p_win < 0.95) && (k < 20)) {
    printf("%2d: %4.1f %%    (%.2f)\n", k, p_win * 100.0, av_loss);
    for (i = 1; i <= hp_b; i++) {
      calculate_probs(&p1, &p2, &p3, a, b, r, s, hp_a, hp_b, dum_array);
      av_loss += probvec1[hp_b - i] * p3;
    }
    resolve_next_fight(probvec1, probvec2, hp_b); /* return val. in probvec1 */
    p_win = probvec1[hp_b];
    k++;
  }
  printf("%2d: %4.1f %%    (%.2f)\n\n", k, p_win * 100.0, av_loss);  
  /* then for B */
  for (i = 0; i <= MAXHP; i++) {
    probvec1[i] = 0.0;
    probvec2[i] = 0.0;
  }
  for (i = 0; i <= hp_a; i++) {
    for (j = 0; j <= hp_b; j++) probvec1[hp_a - i] += array[i][j];
    probvec2[hp_a - i] = probvec1[hp_a - i];
  }
  calculate_probs(&p_win, &p_retreat, &p_lose, a, b, r, s, hp_a, hp_b, array);
  av_loss = p_win;
  printf("Victory chance for n of unit B vs 1 unit A:\n(with expected "
	 " no. of B units killed)\n n   prob    losses\n");
  k = 1;
  while ((p_lose < 0.95) && (k < 20)) {
    printf("%2d: %4.1f %%    (%.2f)\n", k, p_lose * 100.0, av_loss);
    for (i = 1; i <= hp_a; i++) {
      calculate_probs(&p1, &p2, &p3, a, b, r, s, hp_a, hp_b, dum_array);
      av_loss += probvec1[hp_a - i] * p1;
    }
    resolve_next_fight(probvec1, probvec2, hp_a); /* return val in probvec1 */
    p_lose = probvec1[hp_a];
    k++;
  }
  printf("%2d: %4.1f %%    (%.2f)\n\n", k, p_lose * 100.0, av_loss);  

  return(0);
}

void  read_command_line(double *a_st, double *b_st, int *hp_a, int *hp_b,
			double *r, int *level, double retreat[MAXLEVEL],
			int argc, char *argv[])
{
  int i, j;

  /* deal with command line options */
  i = 0;
  while (++i < argc) {
    ++(argv[i]);
    if ((strcmp(argv[i], "A") == 0) || (strcmp(argv[i], "a") == 0)) {
      j = sscanf(argv[++i], "%lf", a_st);
      if (j != 1) {
	fprintf(stderr, "Can't convert a_st (%s)\n", argv[i]);
	exit(-1);
      }
      j = sscanf(argv[++i], "%d", hp_a);
      if (j != 1) {
	fprintf(stderr, "Can't convert hp_a (%s)\n", argv[i]);
	exit(-1);
      }
      if ((*a_st < 0.0) || (*hp_a < 1) || (*hp_a > MAXHP - 1)) {
	fprintf(stderr, "Invalid values for unit A\n");
	exit(-1);
      }
    }
    else if ((strcmp(argv[i], "B") == 0) || (strcmp(argv[i], "b") == 0)) {
      j = sscanf(argv[++i], "%lf", b_st);
      if (j != 1) {
	fprintf(stderr, "Can't convert b_st (%s)\n", argv[i]);
	exit(-1);
      }
      j = sscanf(argv[++i], "%d", hp_b);
      if (j != 1) {
	fprintf(stderr, "Can't convert hp_b (%s)\n", argv[i]);
	exit(-1);
      }
      if ((*b_st < 0.0) || (*hp_b < 1) || (*hp_b > MAXHP - 1)) {
	fprintf(stderr, "Invalid values for unit B\n");
	exit(-1);
      }
    }
    else if ((strcmp(argv[i], "R") == 0) || (strcmp(argv[i], "r") == 0)) {
      j = sscanf(argv[++i], "%lf", r);
      if (j != 1) {
	fprintf(stderr, "Can't convert r (%s)\n", argv[i]);
	exit(-1);
      }
      if ((*r < 0.0) || (*r > 1.0)) {
	fprintf(stderr, "Invalid retreat probability (0.0 - 1.0)\n");
	exit(-1);
      }
    }
    else if ((strcmp(argv[i], "E") == 0) || (strcmp(argv[i], "e") == 0)) {
      j = sscanf(argv[++i], "%d", level);
      if (j != 1) {
	fprintf(stderr, "Can't convert level (%s)\n", argv[i]);
	exit(-1);
      }
      if ((*level < 0) || (*level >= MAXLEVEL)) {
	fprintf(stderr, "Invalid experience level\n");
	exit(-1);
      }
      *r = retreat[*level];
    }
    else if ((strcmp(argv[i], "H") == 0) || (strcmp(argv[i], "h") == 0) ||
	     (strcmp(argv[i], "-help") == 0)) {
      /* help info */
      printf("Civ 3 combat probability calculator\n");
      printf("Usage information\n\ncommand line options:\n");
      printf("\nTo specify combat strengths of units\n");
      printf("-a <float> <int> : Combat strength and hit points for unit"
	     " 'A'\n");
      printf("                   Both values must be given if -a is specified"
	     "\n");
      printf("-b <float> <int> : Combat strength and hit points for unit"
	     " 'B'\n");
      printf("\nThe retreat probability for a fast unit can be specified\n"
	     "in either of two ways\n");
      printf("-r <float>       : retreat probability for unit a\n");
      printf("-e <int>         : experience level\n");
      printf("                   0 - no retreat\n");
      printf("                   1 - conscript\n");
      printf("                   2 - standard\n");
      printf("                   3 - veteran\n");
      printf("                   4 - elite\n\n");
      printf("(Or you can alter values in the program and recompile)\n\n");
      exit(0);
    }
    else {
      fprintf(stderr, "Unrecognised option '%s'\n", argv[i]);
      exit(-1);
    }
  }
  return;
}

void calculate_probs(double *p_win, double *p_retreat, double *p_lose,
		     double a, double b, double r, double s, int hp_a,
		     int hp_b, double array[MAXHP][MAXHP])
{
  double coeff, prob; /* used in calculations to store values */
  int i, j;

  for (i = 0; i < MAXHP; i++) {
    for (j = 0; j < MAXHP; j++) array[i][j] = 0.0;
  }


  *p_win = *p_retreat = *p_lose = 0.0;
  /* wins where retreat isn't important */
  for (i = 0; i <= (hp_a - 2); i++) {
    coeff = fac(hp_b + i - 1) / (fac(hp_b - 1) * fac(i));
    prob = coeff * pow(a, hp_b) * pow(b, i);
    array[hp_a - i][0] += prob;
    *p_win += prob;
  }
  /* wins where no retreat because B is on 1 hp */
  coeff = fac(hp_b + hp_a - 3) / (fac(hp_b - 1) * fac(hp_a - 2));
  prob = coeff * pow(a, hp_b) * pow(b, hp_a - 1);
  array[1][0] += prob;
  *p_win += prob;
  /* wins where attempt to retreat fails */
  for (i = 0; i <= (hp_b - 2); i++) {
    coeff = fac(hp_a + i - 2) / (fac(hp_a - 2) * fac(i));
    prob = coeff * pow(a, hp_b) * pow(b, hp_a - 1) * s;
    array[1][0] += prob;
    *p_win += prob;
  }
  /* succesful retreats */
  for (i = 0; i <= (hp_b - 2); i++) {
    coeff = fac(hp_a + i - 2) / (fac(hp_a - 2) * fac(i));
    prob = coeff * pow(a, i) * pow(b, hp_a - 1) * r;
    array[1][hp_b - i] += prob;
    *p_retreat += prob;
  }
  /* losses where no retreat because B is on 1 hp */
  coeff = fac(hp_b + hp_a - 3) / (fac(hp_b - 1) * fac(hp_a - 2));
  prob = coeff * pow(a, hp_b - 1) * pow(b, hp_a);
  array[0][1] += prob;
  *p_lose += prob;
  /* losses where retreat fails */
  for (i = 0; i <= (hp_b - 2); i++) {
    for (j = 0; j <= (hp_b - i - 1); j++) {
      coeff = fac(hp_a + i - 2) / (fac(hp_a - 2) * fac(i));
      prob = coeff * pow(a, i + j) * pow(b, hp_a) * s;
      array[0][hp_b - i -j] += prob;
      *p_lose += prob;
    }
  }
  return;
}

double fac(int x)
{
  if (x == 0) return(1.0);
  else return(x * fac(x - 1));
}

void resolve_next_fight(double vec1[MAXHP], double vec2[MAXHP], int hp)
{
  /* Takes a vector of probabilities that a unit has taken a certain amount
   * of damage (vec1) and a vector of probabilities of how much damage it
   * takes in a fight against another unit (vec2) and produces a
   * modified vector of probabilities of how much damage it has taken
   * after the fight (vec1 - the return value vector) */
  double tarray[MAXHP][MAXHP];
  int i, j;

  for (i = 0; i <= hp; i++) {
    for (j = 0; j <= hp; j++) tarray[i][j] = vec1[i] * vec2[j];
  }
  for (i = 0; i <= hp; i++) {
    vec1[i] = 0.0;
    for (j = 0; j <= i; j++) vec1[i] += tarray[j][i - j];
  }
  for (i = 1; i <= hp; i++) {
    for (j = hp - i + 1; j <= hp; j++) vec1[hp] += tarray[j][i];
  }
  return;
}