Improve dudect test with cropped time analysis

dudect/fixture.c

@@ -43,7 +43,11 @@
 #define ENOUGH_MEASURE 10000
 #define TEST_TRIES 10
 
-static t_context_t *t;
+/* Number of percentiles to calculate */
+#define NUM_PERCENTILES (100)
+#define DUDECT_TESTS (NUM_PERCENTILES + 1)
+
+static t_context_t *ctxs[DUDECT_TESTS];
 
 /* threshold values for Welch's t-test */
 enum {
@@ -56,6 +60,37 @@ static void __attribute__((noreturn)) die(void)
     exit(111);
 }
 
+static int64_t percentile(const int64_t *a_sorted, double which, size_t size)
+{
+    assert(which >= 0 && which <= 1.0);
+    size_t pos = (size_t) (which * size);
+    return a_sorted[pos];
+}
+
+/* leverages the fact that comparison expressions return 1 or 0. */
+static int cmp(const void *aa, const void *bb)
+{
+    int64_t a = *(const int64_t *) aa, b = *(const int64_t *) bb;
+    return (a > b) - (a < b);
+}
+
+/* This function is used to set different thresholds for cropping measurements.
+ * To filter out slow measurements, we keep only the fastest ones by a
+ * complementary exponential decay scale as thresholds for cropping
+ * measurements: threshold(x) = 1 - 0.5^(10 * x / N_MEASURES), where x is the
+ * counter of the measurement.
+ */
+static void prepare_percentiles(int64_t *exec_times, int64_t *percentiles)
+{
+    qsort(exec_times, N_MEASURES, sizeof(int64_t), cmp);
+
+    for (size_t i = 0; i < NUM_PERCENTILES; i++) {
+        percentiles[i] = percentile(
+            exec_times, 1 - (pow(0.5, 10 * (double) (i + 1) / NUM_PERCENTILES)),
+            N_MEASURES);
+    }
+}
+
 static void differentiate(int64_t *exec_times,
                           const int64_t *before_ticks,
                           const int64_t *after_ticks)
@@ -64,7 +99,9 @@ static void differentiate(int64_t *exec_times,
         exec_times[i] = after_ticks[i] - before_ticks[i];
 }
 
-static void update_statistics(const int64_t *exec_times, uint8_t *classes)
+static void update_statistics(const int64_t *exec_times,
+                              uint8_t *classes,
+                              int64_t *percentiles)
 {
     for (size_t i = 0; i < N_MEASURES; i++) {
         int64_t difference = exec_times[i];
@@ -73,15 +110,36 @@ static void update_statistics(const int64_t *exec_times, uint8_t *classes)
             continue;
 
         /* do a t-test on the execution time */
-        t_push(t, difference, classes[i]);
+        t_push(ctxs[0], difference, classes[i]);
+
+        /* t-test on cropped execution times, for several cropping thresholds.
+         */
+        for (size_t j = 0; j < NUM_PERCENTILES; j++) {
+            if (difference < percentiles[j]) {
+                t_push(ctxs[j + 1], difference, classes[i]);
+            }
+        }
     }
 }
 
+static t_context_t *max_test()
+{
+    size_t max_idx = 0;
+    double max_t = 0.0f;
+    for (size_t i = 0; i < NUM_PERCENTILES + 1; i++) {
+        double t = fabs(t_compute(ctxs[i]));
+        if (t > max_t) {
+            max_t = t;
+            max_idx = i;
+        }
+    }
+    return ctxs[max_idx];
+}
+
 static bool report(void)
 {
-    double max_t = fabs(t_compute(t));
+    t_context_t *t = max_test();
     double number_traces_max_t = t->n[0] + t->n[1];
-    double max_tau = max_t / sqrt(number_traces_max_t);
 
     printf("\033[A\033[2K");
     printf("measure: %7.2lf M, ", (number_traces_max_t / 1e6));
@@ -91,6 +149,9 @@ static bool report(void)
         return false;
     }
 
+    double max_t = fabs(t_compute(t));
+    double max_tau = max_t / sqrt(number_traces_max_t);
+
     /* max_t: the t statistic value
      * max_tau: a t value normalized by sqrt(number of measurements).
      *          this way we can compare max_tau taken with different
@@ -123,6 +184,7 @@ static bool doit(int mode)
     int64_t *exec_times = calloc(N_MEASURES, sizeof(int64_t));
     uint8_t *classes = calloc(N_MEASURES, sizeof(uint8_t));
     uint8_t *input_data = calloc(N_MEASURES * CHUNK_SIZE, sizeof(uint8_t));
+    int64_t *percentiles = calloc(NUM_PERCENTILES, sizeof(int64_t));
 
     if (!before_ticks || !after_ticks || !exec_times || !classes ||
         !input_data) {
@@ -133,28 +195,32 @@ static bool doit(int mode)
 
     bool ret = measure(before_ticks, after_ticks, input_data, mode);
     differentiate(exec_times, before_ticks, after_ticks);
-    update_statistics(exec_times, classes);
+    prepare_percentiles(exec_times, percentiles);
+    update_statistics(exec_times, classes, percentiles);
     ret &= report();
 
     free(before_ticks);
     free(after_ticks);
     free(exec_times);
     free(classes);
     free(input_data);
+    free(percentiles);
 
     return ret;
 }
 
 static void init_once(void)
 {
     init_dut();
-    t_init(t);
+    for (size_t i = 0; i < DUDECT_TESTS; i++) {
+        ctxs[i] = malloc(sizeof(t_context_t));
+        t_init(ctxs[i]);
+    }
 }
 
 static bool test_const(char *text, int mode)
 {
     bool result = false;
-    t = malloc(sizeof(t_context_t));
 
     for (int cnt = 0; cnt < TEST_TRIES; ++cnt) {
         printf("Testing %s...(%d/%d)\n\n", text, cnt, TEST_TRIES);
@@ -166,7 +232,11 @@ static bool test_const(char *text, int mode)
         if (result)
             break;
     }
-    free(t);
+
+    for (size_t i = 0; i < DUDECT_TESTS; i++) {
+        free(ctxs[i]);
+    }
+
     return result;
 }