#include <stdio.h>
#include <stdlib.h>
#include <math.h>

int x;
int p=2319;
int x=2319;
int a=8005;
int c=1;
int m=65536;
double rpos[100]={};
double rpossq[100]={};
int na;
int nb;
double prob;
int y;
int z;
FILE*finout;

int main(void) {
setvbuf(stdout,NULL,_IONBF,0);
setvbuf(stderr,NULL,_IONBF,0);
finout=fopen("data_output.xls", "w");

for (na=0; na<1000; na++)
{
    p=0;
    y=0;
    z=0;
    for(nb=0; nb<100; nb++)
    {
        x=(a*x+c)%m;
        prob=((double)x/((double)m-1));
        if (prob<0.11595)
        {
            p=p+1;
        }
        if ((0.11595<prob) && (prob<0.2319))
        {
            p=p-1;
        }

        x=(a*x+c)%m;
        prob=((double)x/((double)m-1));
        if (prob<0.11595)
        {
            y=y+1;
        }
        if ((0.11595<prob) && (prob<0.2319))
        {
            y=y-1;
        }

        x=(a*x+c)%m;
        prob=((double)x/((double)m-1));
        if (prob<0.11595)
        {
            z=z+1;
        }
        if ((0.11595<prob) && (prob<0.2319))
        {
            z=z-1;
        }

        rpos[nb]+=sqrt(pow(p,2)+pow(y,2)+pow(z,2));
        rpossq[nb]+=pow(sqrt(pow(p,2)+pow(y,2)+pow(z,2)),2);
    }

}

for(nb=0; nb<100; nb++)
{
    rpos[nb]=rpos[nb]/1000;
    rpossq[nb]=rpossq[nb]/1000;
    fprintf(finout, "%d\t%g\t%g\n",nb,rpos[nb],rpossq[nb]);
}
return 0;
}