Minor cleanups and optimizations.

diff --git a/foosrank.cpp b/foosrank.cpp
index cc660764762885eacd62abf9362e7b7cca54e3c3..45c08cc1479264adf49909c28f210c8683fa31f7 100644 (file)
--- a/foosrank.cpp
+++ b/foosrank.cpp
@@ -16,9 +16,9 @@ static const double rating_constant = 455.0;
 
 using namespace std;
 
-double prob_score(int k, double a, double rd);
-double prob_score_real(int k, double a, double binomial, double rd_norm);
-double prodai(int k, double a);
+double prob_score(int k, int a, double rd);
+double prob_score_real(int k, int a, double binomial, double rd_norm);
+double prodai(int k, int a);
 double fac(int x);
 
 
@@ -38,23 +38,40 @@ double fac(int x);
 // Glicko/Bradley-Terry assumption that a player rated 400 points over
 // his/her opponent will win with a probability of 10/11 =~ 0.90909. 
 //
-double prob_score(int k, double a, double rd)
+double prob_score(int k, int a, double rd)
 {
        return prob_score_real(k, a, prodai(k, a) / fac(k-1), rd/rating_constant);
 }
 
+// computes x^a, probably more efficiently than pow(x, a) (but requires that a
+// is n unsigned integer)
+double intpow(double x, unsigned a)
+{
+       double result = 1.0;
+
+       while (a > 0) {
+               if (a & 1) {
+                       result *= x;
+               }
+               a >>= 1;
+               x *= x;
+       }
+
+       return result;
+}
+
 // Same, but takes in binomial(a+k-1, k-1) as an argument in
 // addition to a. Faster if you already have that precomputed, and assumes rd
 // is already divided by 455.
-double prob_score_real(int k, double a, double binomial, double rd_norm)
+double prob_score_real(int k, int a, double binomial, double rd_norm)
 {
-       double nom = binomial * pow(2.0, rd_norm * a); 
-       double denom = pow(1.0 + pow(2.0, rd_norm), k+a);
+       double nom = binomial * intpow(pow(2.0, rd_norm), a); 
+       double denom = intpow(1.0 + pow(2.0, rd_norm), k+a);
        return nom/denom;
 }
 
 // Calculates Product(a+i, i=1..k-1) (see above).
-double prodai(int k, double a)
+double prodai(int k, int a)
 {
        double prod = 1.0;
        for (int i = 1; i < k; ++i)
@@ -70,48 +87,11 @@ double fac(int x)
        return prod;
 }
 
-// 
-// Computes the integral
-//
-//   +inf
-//    /
-//    |
-//    | ProbScore[a] (r1-r2) Gaussian[mu2, sigma2] (r2) dr2
-//    |
-//   /
-// -inf
-//
-// For practical reasons, -inf and +inf are replaced by 0 and 3000, which
-// is reasonable in the this context.
-//
-// The Gaussian is not normalized.
-//
-// Set the last parameter to 1.0 if player 1 won, or -1.0 if player 2 won.
-// In the latter case, ProbScore will be given (r2-r1) instead of (r1-r2).
-//
-class ProbScoreEvaluator {
-private:
-       int k;
-       double a;
-       double binomial_precompute, r1, mu2, sigma2, winfac;
-
-public:
-       ProbScoreEvaluator(int k, double a, double binomial_precompute, double r1, double mu2, double sigma2, double winfac)
-               : k(k), a(a), binomial_precompute(binomial_precompute), r1(r1), mu2(mu2), sigma2(sigma2), winfac(winfac) {}
-       inline double operator() (double x) const
-       {
-               double probscore = prob_score_real(k, a, binomial_precompute, (r1 - x)*winfac);
-               double z = (x - mu2)/sigma2;
-               double gaussian = exp(-(z*z/2.0));
-               return probscore * gaussian;
-       }
-};
-
 void convolve(int size)
 {
 }
 
-void compute_opponent_rating_pdf(int k, double a, double mu2, double sigma2, double winfac, vector<pair<double, double> > &result)
+void compute_opponent_rating_pdf(int k, int a, double mu2, double sigma2, double winfac, vector<pair<double, double> > &result)
 {
        double binomial_precompute = prodai(k, a) / fac(k-1);
        winfac /= rating_constant;