// PROGRAM random_walk
// 1次元ランダムウォークのシミュレーション

#include "TrueBASIC.h"
int main();
void initial(double *p, int *N, int xcum[], int x2cum[], int *ntrial);
void walk(double p, int N, int xcum[], int x2cum[]);
void data(int x, int xcum[], int x2cum[], int istep);
void output(int N, int xcum[], int x2cum[], int ntrial);

int main()
{
    int N, itrial, ntrial, x2cum[65], xcum[65];
    double p;

    initial(&p, &N, xcum, x2cum, &ntrial);
    for(itrial = 1; itrial <= ntrial; ++itrial)
        walk(p, N, xcum, x2cum);
    output(N, xcum, x2cum, ntrial);

    return 0;
}

void initial(double *p, int *N, int xcum[], int x2cum[], int *ntrial)
{
    int istep;
    char Stmp1_[_LBUFF_];

    rnd(-1);
    printf("maximum number of steps N = ");
    fgets(Stmp1_, _LBUFF_, stdin);
    sscanf(Stmp1_, "%d", N);
    (*p) = 0.5;
    (*ntrial) = 1000;                       // 試行回数
    for(istep = 1; istep <= (*N); ++istep)
    {
        xcum[istep] = 0;                    // True BASIC では必要ではない
        x2cum[istep] = 0;
    }
}

void walk(double p, int N, int xcum[], int x2cum[])
{
    int istep, x;

    x = 0;                                  // 各試行における酔歩者の初期位置
    for(istep = 1; istep <= N; ++istep)
    {
        if(rnd(0) <= p)
            ++x;
        else
            --x;
        // すべての歩を進めた後にデータを集める
        data(x, xcum, x2cum, istep);
    }
}

void data(int x, int xcum[], int x2cum[], int istep)
{
    xcum[istep] += x;
    x2cum[istep] += x*x;
}

void output(int N, int xcum[], int x2cum[], int ntrial)
{
    int istep;
    double variance, x2bar, xbar;

    printf("     # steps          <x>        <x^2> <x^2> - <x>^2\n");
    printf("\n");
    for(istep = 1; istep <= N; ++istep)
    {
        xbar = (double)xcum[istep]/ntrial;
        x2bar = (double)x2cum[istep]/ntrial;
        variance = x2bar - xbar*xbar;
        printf("%12d %12f %12f %12f\n", istep, xbar, x2bar, variance);
    }
}