同期加算で誤差を減らしていくようなやつです

out.csv

synchronized_addition.d

import std.stdio;
import std.random;
import std.math;
import std.range;
import std.array;
import std.algorithm;
import std.conv;

auto boxMuller(T)(double avg, double dist, T uniform)
{
    return (sqrt(dist) * sqrt(-2*log(uniform())) * sin(2*PI*uniform()) + avg);
}


void main()
{    
    auto rnd = Random(unpredictableSeed);
    auto uni = () => uniform(0.0, 1.0, rnd);
    auto nor = () => boxMuller(0, 1, uni);

    const long SAMPLING_N = 128;
    const long AMPLITUDE = 10;
    const long FREQUENCY = 15;

    double[] sin_wave = [];
    foreach (i; 0..SAMPLING_N) {
      double x = AMPLITUDE * sin(2 * cast(double)PI * FREQUENCY * i / SAMPLING_N);
      sin_wave ~= x;
    }

    const long N = 10000; 
    double[][] noised_waves = [];
    double[][] averaged_waves = [];
    foreach (i; 0..N) {
      double[] noise = generate!(nor).takeExactly(SAMPLING_N).array.to!(double[]);
      double[] noised_wave = zip(sin_wave, noise).map!(x => x[0] + x[1])().array;

      noised_waves ~= noised_wave;
      double[] averaged_wave = [];
      foreach (j; 0..SAMPLING_N) {
        double x = 0;
        foreach (n_wave; noised_waves) {
          x += n_wave[j];
        }
        averaged_wave ~= x / noised_waves.length;
      }
      averaged_waves ~= averaged_wave;
    }

    double[] noises = [];
    foreach (i; 0..N) {
      double x = 0;
      foreach (j; 0..SAMPLING_N) {
        x += abs(averaged_waves[i][j] - sin_wave[j]);
      }
      noises ~= x;
    }

    write("A");
    writef(", C%05d", 1);
    writef(", C%05d", 10000);
    writeln(", Nn");

    foreach (i; 0..SAMPLING_N) {
      writef("%.4f", sin_wave[i]);
      writef(", %.4f", averaged_waves[0][i]);
      writef(", %.4f", averaged_waves[10000-1][i]);
      writefln(", %.4f", noises[i]);
    }

    foreach (i; SAMPLING_N..N) {
      write(" , , ,");
      writefln(" %.4f", noises[i]);
    }   
}