From 8d1cc85eba46cb1e1ebaa4986312562b893734b8 Mon Sep 17 00:00:00 2001
From: Guillaume Horel <guillaume.horel@gmail.com>
Date: Thu, 28 Apr 2016 20:49:55 -0400
Subject: try to be more generic

---
 src/projet.cpp              | 12 ++++----
 src/rqmc.hpp                | 69 +++++++++++----------------------------------
 src/stratified_sampling.hpp |  4 +--
 3 files changed, 26 insertions(+), 59 deletions(-)

diff --git a/src/projet.cpp b/src/projet.cpp
index 311c5bb..de75fd1 100644
--- a/src/projet.cpp
+++ b/src/projet.cpp
@@ -48,11 +48,12 @@ vector< vector<double> > exemple1_stratified() {
 
 vector< vector<double> > exemple1_rqmc(){
     int I = 100;
-    vector<int> N = {3, 13, 113, 313}; //les N choisis pour que les NI soient égaux aux N de l'exemple 1 stratified_sampling
+    vector<int> N = {3, 13, 113, 313}; //les N choisis pour que les N*I soient égaux aux N de l'exemple 1 stratified_sampling
     first f; //comme quasi_gaussian retourne un vecteur, on doit composer avec f pour avoir le double QG()[0]
     vector< vector<double> > data (4);
     for(size_t i =0; i<N.size(); i++){
-        data[i] = monte_carlo (I, quasi_mean<struct first, sobol> (N[i], 1, f));
+        quasi_gaussian<sobol> QG(1);
+        data[i] = monte_carlo (I, mean(N[i], compose(f, QG)));
     }
     cout<<"moyenne"<<"\t\t"<<"sigma"<<"\t\t"<<"taille IC"<<endl;
     for(int i =0; i<3; i++){
@@ -112,7 +113,7 @@ vector <vector<double> > exemple2_stratified (int d, bool call = true){
 }
 
 vector< vector<double> > exemple2_rqmc(int d, bool call = true) {
-    int N= 10000;
+    int N = 10000;
     double r = 0.05;
     double T = 1.0;
     double S0 = 50;
@@ -122,7 +123,8 @@ vector< vector<double> > exemple2_rqmc(int d, bool call = true) {
     int I=100;//la taille du vrai Monte-Carlo
     for(int i =0; i<3; i++) {
         asian_option A(r, T, S0, V, K[i], call);
-        data[i] = monte_carlo(I, quasi_mean<asian_option, sobol> (N, d, A));
+        quasi_gaussian<sobol> QG(d);
+        data[i] = monte_carlo(I, mean(N, compose(A, QG)));
     }
     return data;
 };
@@ -170,7 +172,7 @@ int main()
     init_alea(0);
     cout<<"comparaison stratified sampling/RQMC sur l'example 1\n"<<endl;
     make_table1(exemple1_stratified(), exemple1_rqmc());
-    
+
     vector< vector<double> > data1 = exemple2_stratified(16);
     vector< vector<double> > data2 = exemple2_rqmc(16);
     vector< vector<double> > data3 = exemple2_stratified(64);
diff --git a/src/rqmc.hpp b/src/rqmc.hpp
index 840e069..a150447 100644
--- a/src/rqmc.hpp
+++ b/src/rqmc.hpp
@@ -9,7 +9,7 @@ double frac_part(double x);
 //fonctions de compositions de monte_carlo.hpp
 
 template <class _Result>
-struct generator 
+struct generator
 {
 	typedef _Result result_type;
 };
@@ -18,9 +18,12 @@ template <typename Fct, typename VA>
 struct compose_t : public generator< typename Fct::result_type >
 {
 	compose_t(Fct f, VA X) : f(f), X(X) {};
-	typename Fct::result_type operator()() { 
-		return f(X()); 
+	typename Fct::result_type operator()() {
+		return f(X());
 	};
+    void reseed() {
+        X.reseed();
+    }
 	private:
 		Fct f; VA X;
 };
@@ -36,6 +39,7 @@ struct mean_t: public generator<double>
 {
     mean_t(int n, VA X): n(n), X(X) {};
     double operator()(){
+        X.reseed();
         double sum = 0;
         for(int i=0; i<n; i++){
             sum+=X();
@@ -45,7 +49,7 @@ struct mean_t: public generator<double>
     private :
         int n; VA X;
     };
-    
+
 template <typename VA>
 inline mean_t<VA>
 mean(int n, VA X){
@@ -74,7 +78,7 @@ std::vector<double> monte_carlo(int n, L X, double p=0.05)
 template <typename L, typename Fct>
 std::vector<double> monte_carlo(int n, Fct f, L X, double p=0.05)
 {
-	return monte_carlo(n, compose(f, X), p); 
+	return monte_carlo(n, compose(f, X), p);
 };
 
 //Les classes de quasi-mean
@@ -84,8 +88,8 @@ struct quasi_gaussian : public var_alea<std::vector<double> >
 {
     quasi_gaussian (int d, double mean = 0, double std =1)
         : d(d), mean(mean), std(std), s(d), U(0,1), seed(d) {
-            for(int i=0; i<d; i++) {seed[i]=U();}
-                };
+        reseed();
+    };
     std::vector<double> operator ()(){
         std::vector<double> result = s();
         for(int i=0;i<d; i++){
@@ -93,54 +97,15 @@ struct quasi_gaussian : public var_alea<std::vector<double> >
         }
         return value = result;
     }
-    private: 
+    void reseed() {
+        for(int i=0; i<d; i++) {
+            seed[i]=U();
+        }
+    }
+    private:
     int d;
     double mean, std;
     LDS s;
     uniform U;
     std::vector<double> seed;
 };
-
-
-
-struct gaussian_d : public var_alea<std::vector<double> >
-{
-    gaussian_d(int d, double mean = 0, double std =1)
-        : d(d), mean(mean), std(std), G(mean, std) {};
-    std::vector<double> operator ()(){
-        std::vector<double> result(d);
-        for(int i=0;i<d; i++){
-            result[i] = G();
-        }
-        return value = result;
-    }
-    private: 
-    int d;
-    double mean, std;
-    gaussian G;
-};
-
-template <typename Fct, typename LDS>
-struct quasi_mean : public generator<typename Fct::result_type>
-{
-    quasi_mean(int n, int d, Fct f) : n(n), d(d), f(f), G(d) {};
-        
-    typename Fct::result_type operator()() {
-        double sum =0;
-        for(int i=0; i<n; i++){
-            sum += f(G());
-            
-        }
-        return sum/n;
-    };
-    
-    private:
-        int n, d;
-        Fct f;
-        quasi_gaussian<LDS> G;
-    };
-
-
-
-
-    
diff --git a/src/stratified_sampling.hpp b/src/stratified_sampling.hpp
index 9175f75..25eef80 100644
--- a/src/stratified_sampling.hpp
+++ b/src/stratified_sampling.hpp
@@ -48,7 +48,7 @@ struct multi_gaussian_truncated : public var_alea<std::vector<double> >
         }
         return X;
     }
-
+    void reseed() {};
     private:
         double a, b;
         uniform V;
@@ -124,7 +124,7 @@ void stratified_sampling<L>::update(int Nk) {
 
 template <typename L>
 void stratified_sampling<L>::draw(int N) {
-		update(N);
+    update(N);
     double m, s, oldmean;
     for(int i=0;i<I;i++){
         m=0;
-- 
cgit v1.2.3-70-g09d2