A traditional numeric solution to the Cover/Keilers Offensive Earned Run Average equation set is enclosed by the horizontal lines. Save this file as a text file and edit out the html preceeding and following the C++ source code.

// oera_num.c
// Copyright (c) 1997-1998 by John F. Jarvis
// This surce code may be freely copied and distributed as long as
// proper credit is given.
// This software is provided "as is" without any express or 
// implied warranties, including, without limitation, the implied
// warranties of merchantibility and fitness for a particular purpose.
 
#include "matrix.h"
// matrix.h is not a standard header file that defines classes
// mtrx, vctr
 
// Offensive Earned Run Average
// Thomas M. Cover, Carroll W. Keilers
// Operations Reseach, Vol 25, No 5, Sep-Oct 1977 pp 729-740
// numeric solution to equation set
 
// entry points in separately compiled files
double lusolve(mtrx&, vctr&, vctr&);
 
double oera( int outs, int pflag, double bf, double sn, double db,
	double tr, double hr, double bb, double *er)
{ // parameters:
  // outs: number of outs in an inning
  // pflag: use != 0 value to obtain printing from the function
  // bf: batters facing pitching = at bats + walks
  // sn, db, tr, hr: numbers of singles, doubles, triples, home runs
  // bb: walks. can add hit by pitch counts here and to bf if desired
  // er: NULL or array long enough to store solution vector
  // return value is expected runs per game defined 9*e[0]
  
 int states = 8*outs; // outs in an inning
 
 // mtrx and vctr are personal set of classes defined in matrix.h. They can
 // be replaced with one and two dimensional arrays. lusolve() will need
 // to be modified (actually unmodified).
 mtrx q(states,states);
 vctr r(states),e(states);
 
 double p0,p1,p2,p3,p4,pb;
 int i,j,k,m;
 
 p1 = sn/bf; // probablities of each event
 p2 = db/bf;
 p3 = tr/bf;
 p4 = hr/bf;
 pb = bb/bf;
 p0 = 1.0 - (p1+p2+p3+p4+pb); // probability of an out
 if(p0<=0.0){
 	printf("oera p0 = %8.5f, <=0\n",p0);
 	return 0.0;
 }
 
 // build i-q matrix and r vector
 // state definition : 8*outs + 
 //		runner on third|second|first as 3 bit binary number
 // first = 1, second = 2, third = 4
 
 for(i=0;i<states;i++) for(j=0;j<states;j++) q[i][j] = 0.0;
 for(m=0;m<states;m += 8){
 	// initialize basic 8*8 blocks, all possible transistions
 	// without having an out. repeat for the number of outs
	for(i=0;i<8;i++){
 		q[m+i][m] = i==0?(1.0-p4):(-p4);
 		q[m+i][m+2] = i==2?(1.0-p2):(-p2);
 		q[m+i][m+4] = i==4?(1.0-p3):(-p3);
	}
 	q[m][m+1] = -p1 - pb;
 	q[m+2][m+1] = q[m+4][m+1] = q[m+6][m+1] = q[m+1][m+5] = q[m+3][m+5] =
 		 q[m+7][m+5] = -p1;
 	q[m+5][m+5] = 1.0 - p1;
 	q[m+1][m+3] = q[m+2][m+3] = q[m+4][m+5] = q[m+3][m+7] = q[m+5][m+7] =
 		q[m+6][m+7] = -pb;
 	q[m+7][m+7] = 1.0 - pb;
 	q[m+1][m+1] = q[m+3][m+3] = q[m+6][m+6] = 1.0;
 }
 for(i=0;i<states-8;i++)
 	q[i][i+8] = -p0; // prob of an out couples states i to i+8
 r[0] = p4; // probability of runs from state i
 r[1] = p2 + p3 + 2.0*p4;
 r[4] = r[2] = p1 + p2+p3+2.0*p4;
 r[3] = r[5] = p1 + 2.0*(p2+p3) + 3.0*p4;
 r[6] = 2.0*(p1+p2+p3) + 3.0*p4;
 r[7] = 2.0*p1 + 3.0*(p2+p3) + 4.0*p4 + pb;
 for(j=8;j<states;j++) r[j] = r[j-8];
 
 lusolve(q,r,e); // see Numerical Recipes
 
 if(pflag != 0){
 	printf("bf %6.0f sn %6.0f db %6.0f tr %6.0f hr %6.0f bb %6.0f\n",bf,sn,db,tr,hr,bb);
 	printf("p0 %6.4f p1 %6.4f p2 %6.4f p3 %6.4f p4 %6.4f pb %6.4f\n",p0,p1,p2,p3,p4,pb);
 	printf("outs ");
 	for(i=0;i<8;i++) printf("   %1d%1d%1d  ",i&01, (i&02)>>1, (i&04)>>2);
 	for(i=0;i<states;i++){
 		if((i%8)==0) printf("\n %1d  ",i/8);
 		printf(" %7.5f",e[i]);
 	}
 	printf("\n");
 }
 if(er!=NULL) for(i=0;i<states;i++) er[i] = e[i]; // expected runs for game state
 return 9.0*e[0]; // expected runs for 9 innings (game)
}