/* Numerical solution to second Baseball Markov chain model */
/* with double plays allowed when nobody is out yet */

/* Written by Carey Williamson, January 2019 */

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

/* Transition rates for Markov chain model */
#define LAMBDA1 0.20   /* transition rate to get first out */
#define LAMBDA2 0.25   /* transition rate to get second out */
#define LAMBDA3 0.25   /* transition rate to get third out */
#define LAMBDA4 0.50   /* transition rate for TV commercial */

#define DEBUG 1   /* verbose debugging */

int main ()
  {
    float p0, p1, p2, p3;
    float sum, ratio;
    float lambdaD;

    printf("Enter double-play rate (lambdaD): ");
    scanf("%f", &lambdaD);

    /* start with p0 = 1 and compute relative probabilities */
    p0 = 1.0;
#ifdef DEBUG
    printf("p0 is %f\n", p0);
#endif

    ratio = LAMBDA1 / LAMBDA2;
    p1 = ratio * p0;
#ifdef DEBUG
    printf("p1 is %f times p0, so p1 = %f\n", ratio, p1);
#endif

    p2 = (lambdaD * p0 + LAMBDA2 * p1) / LAMBDA3;
#ifdef DEBUG
    printf("p2 is %f\n", p2);
#endif

    ratio = LAMBDA3 / LAMBDA4;
    p3 = ratio * p2;
#ifdef DEBUG
    printf("p3 is %f times p2, so p3 = %f\n", ratio, p3);
#endif

    /* Renormalize the probabilities now */
    sum = p0 + p1 + p2 + p3;
#ifdef DEBUG
    printf("sum is %f\n", sum);
#endif
    p0 /= sum;
    p1 /= sum;
    p2 /= sum;
    p3 /= sum;
    
    printf("\n");
    printf("Analytical Markov Chain Model of Baseball with Double Play\n");
    printf("\n");
    printf("  Prob(nobody out): p0 = %f\n", p0);
    printf("  Prob(one out):    p1 = %f\n", p1);
    printf("  Prob(two out):    p2 = %f\n", p2);
    printf("  Prob(three out):  p3 = %f\n", p3);
    printf("\n");
    printf("Sanity check: p0 + p1 + p2 + p3 = %f\n",
	   p0 + p1 + p2 + p3);
    printf("\n");
  }