rqmc.cpp


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

#include "low_discrepancy.hpp"
#include <vector>
#include <gsl/gsl_cdf.h>
#include "var_alea.hpp"


double frac_part(double x){
    return x - floor(x);
}

double mean_rqmc(int N, double X) {
    sobol s(1);
    double sum = 0;
    for(int i=0; i<N; i++){
        sum +=gsl_cdf_gaussian_Pinv (frac_part(X+s()[0]), 1);
    }
    return sum/N;
}

int main() {
    init_alea(0);
    int I=100;
    int N= 10000;
    uniform U;
    double m = 0;
    double s = 0;
    double temp;
    
    for(int i=0;i<I;i++){
        temp = mean_rqmc(N,U());
        m+=temp;
        s+=temp*temp;
    }
    m = m/N;
    s = s/N - m*m;
    
    std::cout<<"espérance "<<m<<" taille de l'IC "<<sqrt(s)*1.64/10<<std::endl;

    return 0;
}