compmod/90.cpp at master · paavalipopov/compmod · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
#include <iostream>
#include <ctime>
#include <cmath>
#include <vector>
#include <cstdlib>
#include <fstream>
#include <string>

using namespace std;

double stand() {
    return 1.0 * rand()/RAND_MAX;
}


class Chamber {
    int n, n_start;
    int time_av, flag;
    double dispersion, n_av, n_sq_av;
    fstream results;

public:
    Chamber(int n, string fout) {
        this->n = n;
        n_start = n;
        dispersion = n_av = n_sq_av = flag = 0;
        results.open(fout);
        results << "Number of molecules - " << this->n << endl << endl;
        results << "Step \tn" << endl;
    }

    Chamber() {};

    ~Chamber() {
        results << endl << endl;
        results.close();
    };

    void make_step(int step) {
        if(stand() <= 1.0*n/n_start)
			n--;
		else
			n++;

        if(n <= n_start/2 && flag == 0) {
            flag = 1;
            time_av = step;
            n_av = 0;
        }
        if(flag == 1) {
            n_av += n;
            n_sq_av += n*n;
        }
        results << ' ' << step << "\t" << n << endl;
    }


    void average_calculation(int final_step) {
        n_av = n_av / (final_step - time_av);
        n_sq_av = n_sq_av / (final_step - time_av);
        dispersion = n_sq_av - n_av * n_av;
    }

    double get_dispersion() {
        return dispersion;
    }

    double get_n_av() {
        return n_av;
    }

    double get_n_sq_av() {
        return n_sq_av;
    }

    void print_final_state() {
        cout << "Number of molecules - " << n_start << endl;
        cout << "n average - " << n_av << endl;
        cout << "Dispersion - " << dispersion << endl;
        results << endl << "n average" << n_av << "\t dispersion - " << dispersion << endl << endl;
    }
};

int main(void) {
    srand(time(0));

    const int number_of_molecules = 160;
    const int number_of_ensembles = 100;

    double dispersion = 0;
    double n_av = 0;
    double n_sq_av = 0;


	fstream fs("results.txt", ios::out);
	fs.close();


    for(int i = 0; i < number_of_ensembles; i++) {
        cout << endl << endl << "Experiment #" << i <<endl;
        Chamber cmb(number_of_molecules, "results.txt");

        for(int t = 0; t < 10000; t++) {
            cmb.make_step(t);
        }
        cmb.average_calculation(9999);
        cmb.print_final_state();

        n_av += cmb.get_n_av();
        n_sq_av += cmb.get_n_sq_av();

        cout << "dispersion / n average = " << sqrt(cmb.get_dispersion()) / cmb.get_n_av() << endl;
        cout << "1 / sqrt(n) = " << 1.0 / sqrt(number_of_molecules) << endl;

    }


    n_av /= number_of_ensembles;
    n_sq_av /= number_of_ensembles;
    dispersion  = n_sq_av - n_av * n_av;

    cout << endl << endl << "Final results for ensembles" << endl << "n average = " << n_av << endl;
    cout << "n sq average = " << n_sq_av << endl;
    cout << "dispersion = " << dispersion << endl;
    cout << "sqrt(dispersion) / n average = " << sqrt(dispersion)/n_av << endl;
    cout << "1 / sqrt(n) = " << 1.0 / sqrt(number_of_molecules) << endl;

    return 0;
}