librcd normal distribution

gistfile1.txt

static double gaussian_noise(double mu, double sigma) {
    double u1;
    do {
        u1 = ((double) lwt_rdrand64()) * (1.0 / (double) UINT64_MAX);
    } while (u1 <= (DBL_EPSILON));
    double u2 = ((double) lwt_rdrand64()) * (1.0 / (double) UINT64_MAX);
    double z0 = sqrt(-2.0 * log(u1)) * cos((2.0 *  M_PI) * u2);
    return z0 * sigma + mu;
}

void rcd_main(list(fstr_t)* main_args, list(fstr_t)* main_env) {
    const int64_t stdd_width = 10;
    const int64_t n_stdd = 4;
    const size_t total = 100000;
    size_t hits[stdd_width * n_stdd * 2] = {0};
    for (size_t i = 0; i < total; i++) {
        double x = gaussian_noise(0, 1);
        int64_t offs = (stdd_width * n_stdd) + floor(x * stdd_width);
        offs = MIN(MAX(offs, 0), LENGTHOF(hits) - 1);
        hits[offs]++;
    }
    for (size_t i = 0; i < LENGTHOF(hits); i++) {
        rio_debug(concs("[", i, "]: ", hits[i], "\n"));
    }
    for (int64_t dev = 1; dev <= n_stdd; dev++) {
        size_t dev_total = 0;
        for (int64_t i = -(dev * stdd_width); i < (dev * stdd_width); i++) {
            int64_t offs = (stdd_width * n_stdd) + i;
            // DBGFN(dev, " ", offs);
            dev_total += hits[offs];
        }
        rio_debug(concs("[", dev, "]: ", dev_total, " [", (100.0 * ((double) dev_total / (double) total)), "] %\n"));
    }
}