aboutsummaryrefslogtreecommitdiffstats
path: root/src/rqmc.hpp
diff options
context:
space:
mode:
authorBertrand <bertrand.horel@gmail.com>2016-02-19 15:03:51 +0000
committerBertrand <bertrand.horel@gmail.com>2016-02-19 15:03:51 +0000
commitd2b133901a65244934eb642ec8e20c797efaf650 (patch)
treef8d186f8e8ca0886f8f0a464261ba8747242b4e6 /src/rqmc.hpp
parent355e4567e68a76356714e2e58a42dcd78533cf6c (diff)
downloadprojet_C++-d2b133901a65244934eb642ec8e20c797efaf650.tar.gz
nettoyage du dépôt
Diffstat (limited to 'src/rqmc.hpp')
-rw-r--r--src/rqmc.hpp127
1 files changed, 127 insertions, 0 deletions
diff --git a/src/rqmc.hpp b/src/rqmc.hpp
new file mode 100644
index 0000000..1b33713
--- /dev/null
+++ b/src/rqmc.hpp
@@ -0,0 +1,127 @@
+#include "low_discrepancy.hpp"
+#include <vector>
+#include <gsl/gsl_cdf.h>
+#include "var_alea.hpp"
+
+
+double frac_part(double x);
+
+//fonctions de compositions de monte_carlo.hpp
+
+template <class _Result>
+struct generator
+{
+ typedef _Result result_type;
+};
+
+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());
+ };
+ private:
+ Fct f; VA X;
+};
+
+template <typename Fct, typename VA>
+inline compose_t<Fct, VA>
+compose(Fct f, VA X) {
+ return compose_t<Fct, VA>(f, X);
+};
+
+template <typename VA>
+struct mean_t: public generator<double>
+{
+ mean_t(int n, VA X): n(n), X(X) {};
+ double operator()(){
+ double sum = 0;
+ for(int i=0; i<n; i++){
+ sum+=X();
+ }
+ return sum/n;
+ }
+ private :
+ int n; VA X;
+ };
+
+template <typename VA>
+inline mean_t<VA>
+mean(int n, VA X){
+ return mean_t<VA>(n, X);
+};
+
+
+//Les classes de monte-Carlo
+
+template <typename L>
+std::vector<double> monte_carlo(int n, L X)
+{
+ std::vector<double> result(3,0);
+ double x;
+ for (int j = 0; j < n; j++) {
+ x = X();
+ result[0] += x;
+ result[1] += x*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);
+ return result;
+};
+
+template <typename L, typename Fct>
+std::vector<double> monte_carlo(int n, Fct f, L X)
+{
+ return monte_carlo(n, compose(f, X));
+};
+
+//Les classes de quasi-mean
+
+template <typename LDS>
+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();}
+ };
+ std::vector<double> operator ()(){
+ std::vector<double> result = s();
+ for(int i=0;i<d; i++){
+ result[i] = mean + gsl_cdf_gaussian_Pinv(frac_part(result[i]+seed[i]), std);
+ }
+ return value = result;
+ }
+ 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;
+};
+
+
+
+
+
+