[gnumeric] sstest: add lots of fractile tests.

[Morten Welinder <mortenw src gnome org>]

commit 07fc11b67df46e7f1818eb7f738caa5061fe7be5
Author: Morten Welinder <terra gnome org>
Date:   Fri Mar 20 15:55:02 2015 -0400

    sstest: add lots of fractile tests.

 src/sstest.c |  140 ++++++++++++++++++++++++++++++++++++++++++++--------------
 1 files changed, 107 insertions(+), 33 deletions(-)
---
diff --git a/src/sstest.c b/src/sstest.c
index 021d4d9..eb47fe7 100644
--- a/src/sstest.c
+++ b/src/sstest.c
@@ -530,6 +530,43 @@ define_cell (Sheet *sheet, int c, int r, const char *expr)
        sheet_cell_set_text (cell, expr, NULL);
 }
 
+static gboolean
+rand_fractile_test (gnm_float const *vals, int N, int nf, gnm_float const *fractiles)
+{
+       gnm_float f = 1.0 / nf, T = f * N;
+       int *fractilecount = g_new (int, nf + 1);
+       int i;
+       gboolean ok = TRUE;
+
+       for (i = 0; i <= nf; i++)
+               fractilecount[i] = 0;
+
+       for (i = 0; i < N; i++) {
+               gnm_float r = vals[i];
+               int j;
+               for (j = 1; j < nf; j++)
+                       if (r < fractiles[j])
+                               break;
+               fractilecount[j]++;
+       }
+       g_printerr ("Fractile counts:");
+       for (i = 1; i <= nf; i++)
+               g_printerr ("%s%d", (i == 1) ? " " : ", ", fractilecount[i]);
+       g_printerr ("\n");
+
+       g_printerr ("Expected count in each fractile: %.10" GNM_FORMAT_g "\n", T);
+       for (i = 1; i <= nf; i++) {
+               if (!(gnm_abs (fractilecount[i] - T) < 3 * gnm_sqrt (f * N))) {
+                       g_printerr ("Fractile test failure.\n");
+                       ok = FALSE;
+               }
+       }
+
+       g_free (fractilecount);
+
+       return ok;
+}
+
 static gnm_float *
 test_random_1 (int N, const char *expr,
               gnm_float *mean, gnm_float *var,
@@ -667,6 +704,8 @@ test_random_rand (int N)
        int i;
        gboolean ok;
        gnm_float T;
+       gnm_float fractiles[10];
+       const int nf = G_N_ELEMENTS (fractiles);
 
        vals = test_random_1 (N, "=RAND()", &mean, &var, &skew, &kurt);
        ok = TRUE;
@@ -678,7 +717,6 @@ test_random_rand (int N)
                        break;
                }
        }
-       g_free (vals);
 
        T = 0.5;
        if (gnm_abs (mean - T) > 0.01) {
@@ -700,9 +738,18 @@ test_random_rand (int N)
                g_printerr ("Kurt failure [%.10" GNM_FORMAT_g "]\n", T);
                ok = FALSE;
        }
+
+       /* Fractile test */
+       for (i = 1; i < nf; i++)
+               fractiles[i] = i / (double)nf;
+       if (!rand_fractile_test (vals, N, nf, fractiles))
+               ok = FALSE;
+
        if (ok)
                g_printerr ("OK\n");
        g_printerr ("\n");
+
+       g_free (vals);
 }
 
 static void
@@ -722,6 +769,8 @@ test_random_randuniform (int N)
        char *expr;
        gnm_float T;
        int i;
+       gnm_float fractiles[10];
+       const int nf = G_N_ELEMENTS (fractiles);
 
        expr = g_strdup_printf ("=RANDUNIFORM(%.10" GNM_FORMAT_g ",%.10" GNM_FORMAT_g ")", param_l, param_h);
        vals = test_random_1 (N, expr, &mean, &var, &skew, &kurt);
@@ -736,7 +785,6 @@ test_random_randuniform (int N)
                        break;
                }
        }
-       g_free (vals);
 
        T = mean_target;
        g_printerr ("Expected mean: %.10" GNM_FORMAT_g "\n", T);
@@ -769,9 +817,17 @@ test_random_randuniform (int N)
                ok = FALSE;
        }
 
+       /* Fractile test */
+       for (i = 1; i < nf; i++)
+               fractiles[i] = param_l + n * i / (double)nf;
+       if (!rand_fractile_test (vals, N, nf, fractiles))
+               ok = FALSE;
+
        if (ok)
                g_printerr ("OK\n");
        g_printerr ("\n");
+
+       g_free (vals);
 }
 
 static void
@@ -1022,6 +1078,8 @@ test_random_randexp (int N)
        char *expr;
        gnm_float T;
        int i;
+       gnm_float fractiles[10];
+       const int nf = G_N_ELEMENTS (fractiles);
 
        expr = g_strdup_printf ("=RANDEXP(%.10" GNM_FORMAT_g ")", param_l);
        vals = test_random_1 (N, expr, &mean, &var, &skew, &kurt);
@@ -1036,7 +1094,6 @@ test_random_randexp (int N)
                        break;
                }
        }
-       g_free (vals);
 
        T = mean_target;
        g_printerr ("Expected mean: %.10" GNM_FORMAT_g "\n", T);
@@ -1069,9 +1126,17 @@ test_random_randexp (int N)
                ok = FALSE;
        }
 
+       /* Fractile test */
+       for (i = 1; i < nf; i++)
+               fractiles[i] = qexp (i / (double)nf, param_l, TRUE, FALSE);
+       if (!rand_fractile_test (vals, N, nf, fractiles))
+               ok = FALSE;
+
        if (ok)
                g_printerr ("OK\n");
        g_printerr ("\n");
+
+       g_free (vals);
 }
 
 static void
@@ -1089,6 +1154,8 @@ test_random_randgamma (int N)
        char *expr;
        gnm_float T;
        int i;
+       gnm_float fractiles[10];
+       const int nf = G_N_ELEMENTS (fractiles);
 
        expr = g_strdup_printf ("=RANDGAMMA(%.0" GNM_FORMAT_f ",%.10" GNM_FORMAT_g ")", param_shape, 
param_scale);
        vals = test_random_1 (N, expr, &mean, &var, &skew, &kurt);
@@ -1103,7 +1170,6 @@ test_random_randgamma (int N)
                        break;
                }
        }
-       g_free (vals);
 
        T = mean_target;
        g_printerr ("Expected mean: %.10" GNM_FORMAT_g "\n", T);
@@ -1136,9 +1202,17 @@ test_random_randgamma (int N)
                ok = FALSE;
        }
 
+       /* Fractile test */
+       for (i = 1; i < nf; i++)
+               fractiles[i] = qgamma (i / (double)nf, param_shape, param_scale, TRUE, FALSE);
+       if (!rand_fractile_test (vals, N, nf, fractiles))
+               ok = FALSE;
+
        if (ok)
                g_printerr ("OK\n");
        g_printerr ("\n");
+
+       g_free (vals);
 }
 
 static void
@@ -1155,6 +1229,8 @@ test_random_randtdist (int N)
        char *expr;
        gnm_float T;
        int i;
+       gnm_float fractiles[10];
+       const int nf = G_N_ELEMENTS (fractiles);
 
        expr = g_strdup_printf ("=RANDTDIST(%.0" GNM_FORMAT_f ")", param_df);
        vals = test_random_1 (N, expr, &mean, &var, &skew, &kurt);
@@ -1169,7 +1245,6 @@ test_random_randtdist (int N)
                        break;
                }
        }
-       g_free (vals);
 
        T = mean_target;
        g_printerr ("Expected mean: %.10" GNM_FORMAT_g "\n", T);
@@ -1202,9 +1277,17 @@ test_random_randtdist (int N)
                ok = FALSE;
        }
 
+       /* Fractile test */
+       for (i = 1; i < nf; i++)
+               fractiles[i] = qt (i / (double)nf, param_df, TRUE, FALSE);
+       if (!rand_fractile_test (vals, N, nf, fractiles))
+               ok = FALSE;
+
        if (ok)
                g_printerr ("OK\n");
        g_printerr ("\n");
+
+       g_free (vals);
 }
 
 static void
@@ -1225,6 +1308,8 @@ test_random_randfdist (int N)
        char *expr;
        gnm_float T;
        int i;
+       gnm_float fractiles[10];
+       const int nf = G_N_ELEMENTS (fractiles);
 
        expr = g_strdup_printf ("=RANDFDIST(%.0" GNM_FORMAT_f ",%.0" GNM_FORMAT_f ")", param_df1, param_df2);
        vals = test_random_1 (N, expr, &mean, &var, &skew, &kurt);
@@ -1239,7 +1324,6 @@ test_random_randfdist (int N)
                        break;
                }
        }
-       g_free (vals);
 
        T = mean_target;
        g_printerr ("Expected mean: %.10" GNM_FORMAT_g "\n", T);
@@ -1272,9 +1356,17 @@ test_random_randfdist (int N)
                ok = FALSE;
        }
 
+       /* Fractile test */
+       for (i = 1; i < nf; i++)
+               fractiles[i] = qf (i / (double)nf, param_df1, param_df2, TRUE, FALSE);
+       if (!rand_fractile_test (vals, N, nf, fractiles))
+               ok = FALSE;
+
        if (ok)
                g_printerr ("OK\n");
        g_printerr ("\n");
+
+       g_free (vals);
 }
 
 static void
@@ -1290,21 +1382,14 @@ test_random_randcauchy (int N)
        gnm_float kurt_target = gnm_nan;
        char *expr;
        gnm_float T;
-       gnm_float fractiles[10 + 1];
-       int i, fractilecount[G_N_ELEMENTS (fractiles)];
-       const int nf = G_N_ELEMENTS (fractiles) - 1;
-       gnm_float f = 1.0 / nf;
+       int i;
+       gnm_float fractiles[10];
+       const int nf = G_N_ELEMENTS (fractiles);
 
        /*
         * The distribution has no mean, no variance, no skew, and no kurtosis.
         * The support is all reals.
         */
-       fractiles[0] = gnm_ninf;
-       for (i = 1; i < nf; i++)
-               fractiles[i] = qcauchy (i * f, 0.0, param_scale, TRUE, FALSE);
-       fractiles[nf] = gnm_pinf;
-       for (i = 0; i <= nf; i++)
-               fractilecount[i] = 0;
 
        expr = g_strdup_printf ("=RANDCAUCHY(%.10" GNM_FORMAT_g ")", param_scale);
        vals = test_random_1 (N, expr, &mean, &var, &skew, &kurt);
@@ -1312,22 +1397,13 @@ test_random_randcauchy (int N)
 
        ok = TRUE;
        for (i = 0; i < N; i++) {
-               int j;
                gnm_float r = vals[i];
                if (!(gnm_finite (r))) {
                        g_printerr ("Range failure.\n");
                        ok = FALSE;
                        break;
                }
-               for (j = 1; j <= nf; j++)
-                       if (r < fractiles[j])
-                               break;
-               fractilecount[j]++;
        }
-       g_free (vals);
-
-       for (i = 1; i <= nf; i++)
-               g_printerr ("Fractile %d count: %d\n", i, fractilecount[i]);
 
        T = mean_target;
        g_printerr ("Expected mean: %.10" GNM_FORMAT_g "\n", T);
@@ -1361,18 +1437,16 @@ test_random_randcauchy (int N)
        }
 
        /* Fractile test */
-       T = f * N;
-       g_printerr ("Expected count in each fractile: %.10" GNM_FORMAT_g "\n", T);
-       for (i = 1; i <= nf; i++) {
-               if (!(gnm_abs (fractilecount[i] - T) < 3 * gnm_sqrt (d * N))) {
-                       g_printerr ("Fractile test failure.\n");
-                       ok = FALSE;
-               }
-       }
+       for (i = 1; i < nf; i++)
+               fractiles[i] = qcauchy (i / (double)nf, 0.0, param_scale, TRUE, FALSE);
+       if (!rand_fractile_test (vals, N, nf, fractiles))
+               ok = FALSE;
 
        if (ok)
                g_printerr ("OK\n");
        g_printerr ("\n");
+
+       g_free (vals);
 }
 
 static void