suite nettoyage, écriture exemple normale rqmc, transformation de la fonction monte_carlo en quasi_mean

4 files changed, 47 insertions, 46 deletions

diff --git a/Makefile b/Makefile
index 4db562b..f58a3cc 100644
--- a/Makefile
+++ b/Makefile
@@ -17,16 +17,12 @@ table.tex: make_table
 make_table: make_table.o mt19937.o
 	$(CXX) $^ -o $@
 
-projet: projet.o stratified_sampling.o mt19937.o var_alea.o opti.o
+projet: projet.o stratified_sampling.o mt19937.o var_alea.o opti.o rqmc.o p_adic.o
 	$(CXX) $^ -o $@ $(GSL_FLAGS) $(NLOPT_FLAGS)
 
 stratified_sampling.o: stratified_sampling.hpp
 
-test: test.o mt19937.o
-	$(CXX) $^ -o $@
-
-rqmc: rqmc.o mt19937.o
-	$(CXX) $^ -o $@ $(GSL_FLAGS) 
+rqmc.o: rqmc.hpp 
 
 option: option.o mt19937.o p_adic.o var_alea.o
 	$(CXX) $^ -o $@ $(GSL_FLAGS) 
diff --git a/src/option.cpp b/src/option.cpp
index 918e735..c50a7bf 100644
--- a/src/option.cpp
+++ b/src/option.cpp
@@ -6,31 +6,7 @@ double pos (double x){
 }
 
 
-int main(){
-    init_alea(1);
-    asian_option A(0.05, 1.0, 50.0, 0.1, 16, 45);
-    int N= 10000;
-    
-    int d =16; 
-
-    
-    std::vector<double> result(3);
-    result = monte_carlo(100, quasi_option<asian_option, sobol> (N, d, A));
-    for(int i =0; i<3; i++){
-        std::cout<<result[i]<<std::endl;
-    }
-    
-    std::vector<double> result2(3);
-    result2 = monte_carlo(100, quasi_option<asian_option, halton> (N, d, A));
-    for(int i =0; i<3; i++){
-        std::cout<<result2[i]<<std::endl;
-    }
-
-
-    
-    return 0;
-}
-        
+  
         
         
         
diff --git a/src/projet.cpp b/src/projet.cpp
index 980ef98..57f6b91 100644
--- a/src/projet.cpp
+++ b/src/projet.cpp
@@ -9,8 +9,7 @@
 
 using namespace std;
 
-
-void exemple1() {
+void exemple1_stratified() {
     vector<double> q = quantile_norm(10, 1);
     vector<double> p(10, 0.1);
     vector<gaussian_truncated> rvar;
@@ -35,7 +34,7 @@ void exemple1() {
 };
 
 
- void exemple2 (){
+ void exemple2_stratified (){
      int d= 16;
      std::vector<double> mu(d);
      mu = argmax(0.05, 1.0, 50, 0.1, 45, d);
@@ -62,18 +61,48 @@ void exemple1() {
      stratified_sampling<compose_t <f_mu, multi_gaussian_truncated> > S(p, X);
      S.draw(1000);
      cout<<"l'estimateur de la moyenne est :"<<S.estimator().first<<endl;
-     //~ compose_t <f_mu, multi_gaussian_truncated> X = compose(G,MG);
-     //~ for(int i=0; i<10; i++){
-        //~ std::cout<<X()<<std::endl;
-    //~ }
 }
 
+void exemple2_rqmc() {
+    asian_option A(0.05, 1.0, 50.0, 0.1, 16, 45);
+    int N= 10000;
+    
+    int d =16; 
+
+    
+    std::vector<double> result(3);
+    result = quasi_mean(100, quasi_option<asian_option, sobol> (N, d, A));
+    for(int i =0; i<3; i++){
+        std::cout<<result[i]<<std::endl;
+    }
+    
+    std::vector<double> result2(3);
+    result2 = quasi_mean(100, quasi_option<asian_option, halton> (N, d, A));
+    for(int i =0; i<3; i++){
+        std::cout<<result2[i]<<std::endl;
+    }
+};
+
+struct first:public std::unary_function<std::vector<double>, double>
+{ double operator()(std::vector<double> X){return X[0];}
+};
+
+void exemple1_rqmc(){
+    int N = 100;
+    first f; //comme quasi_gaussian retourne un vecteur, on doit composer avec f pour avoir le double QG()[0]
+    std::vector<double> result(3);
+    result = quasi_mean (100,quasi_option<struct first, sobol> (N, 1, f));
+    for(int i =0; i<3; i++){
+        std::cout<<result[i]<<std::endl;
+    }
+};
 
 int main()
 {   
     init_alea(1);
-    exemple1();
-    //exemple2();
-        
-    return 0; 
+    //exemple2_rqmc();
+    //exemple1_stratified();
+    exemple1_rqmc();  
+    return 0;
+    
 }
diff --git a/src/rqmc.hpp b/src/rqmc.hpp
index 1b33713..03f7e81 100644
--- a/src/rqmc.hpp
+++ b/src/rqmc.hpp
@@ -56,7 +56,7 @@ mean(int n, VA X){
 //Les classes de monte-Carlo
 
 template <typename L>
-std::vector<double> monte_carlo(int n, L X)
+std::vector<double> quasi_mean(int n, L X, double p=0.05)
 {
 	std::vector<double> result(3,0);
 	double x;
@@ -67,14 +67,14 @@ std::vector<double> monte_carlo(int n, L X)
 	}
 	result[0] /= (double) n;
 	result[1] = (result[1] - n*result[0]*result[0])/(double)(n-1);
-	result[2] = 1.96*sqrt(result[1]/(double) n);
+	result[2] = gsl_cdf_gaussian_Pinv(1-p/2,1)*sqrt(result[1]/(double) n);
 	return result;
 };
 
 template <typename L, typename Fct>
-std::vector<double> monte_carlo(int n, Fct f, L X)
+std::vector<double> quasi_mean(int n, Fct f, L X, double p=0.05)
 {
-	return monte_carlo(n, compose(f, X)); 
+	return quasi_mean(n, compose(f, X), p); 
 };
 
 //Les classes de quasi-mean