MonteCarlo.h 2.83 KB
Newer Older
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
#ifndef MLMC_MC_HPP
#define MLMC_MC_HPP

#include "Elliptic.h"

/*
 * MonteCarlo
 *
 * l:          Level, should be set as the index in the data frames
 * baseLevel:  Bool to indicate if on base level -> Only fine grid calculations
 * dM:         Samples to go
 * M:          Samples so far
 * functional: Type of functional used
 * sum_cost:   Summed cost
 * avg_cost:   Average cost
 * Qf:         Functional on fine level
 * Qc:         Functional on coarse level
 * Q = Qf
 * Y = Qf - Qc
 * sum_Y, sum_Q:       Sum of differences, functionals
 * avg_Y, avg_Q:       Average of differences, functionals
 * sum_Y2, sum_Q2:     Sum of squared differences, functionals
 * E(Y^2), E(Q^2):     Average of squared differences, functionals
 * var_Y = E(Y^2) - (E(Y))^2, V(Q) = E(Q^2) - (E(Q))^2:
 *                     Variance of differences, functionals
 *
 * MCPlotting = 0     for no plots at all
 * MCPlotting = 1     for average plots
 * MCPlotting = 2     for sample and average plots
 * MCPlotting = 3     for sample, average and interpolated fields
 *
 */

class MonteCarlo {
public:
    int verbose = 0, plotting = 0;

    int l = 0;
    bool baseLevel = true;
    string functional = "L2";

    StochasticFields &stochFields;
    EllipticAssemble &assemble;
    MatrixGraphs &graphs;

    Vector u_sum, u_avg;

    int M = 0, dM = 0;

    double sum_cost = 0.0, avg_cost = 0.0;
    double sum_Y = 0.0, avg_Y = 0.0;
    double sum_Y2 = 0.0, avg_Y2 = 0.0, var_Y = 0.0;
    double sum_Y3 = 0.0, avg_Y3 = 0.0;
    double sum_Y4 = 0.0, avg_Y4 = 0.0;
    double sum_Q = 0.0, avg_Q = 0.0;
    double sum_Q2 = 0.0, avg_Q2 = 0.0, var_Q = 0.0;
    double kurtosis = 0.0;

    MonteCarlo(int l, int dM, bool baseLevel, StochasticFields &stochFields,
               EllipticAssemble &assemble, MatrixGraphs &graphs) :
            l(l), dM(dM), baseLevel(baseLevel),
            stochFields(stochFields), assemble(assemble),
            graphs(graphs), u_sum(graphs[l - 1]), u_avg(graphs[l - 1]) {
        ReadConfig(Settings, "MCVerbose", verbose);
        ReadConfig(Settings, "MCPlotting", plotting);
        ReadConfig(Settings, "functional", functional);
    }

    virtual ~MonteCarlo() = default;;

    void updateSum(const double *sums);

    void updateAvg();

    void updateStatistic();

    void upscaleU();

    void updateUSum();

    void updateUAvg();

    virtual void method() = 0;

    virtual void solvePDE(double &Q, double &cost, Vector &u, int m, int l) {};
};

MonteCarlo *getMonteCarlo(int l, int dM, bool baseLevel, StochasticFields &stochFields,
                          EllipticAssemble &assemble, MatrixGraphs &graphs,
                          string MonteCarlo = "");

void estimateExponents(const map<int, MonteCarlo *> &df_mc,
                       double &alpha, double &beta, double &gamma,
                       bool excludeBaseLevel = true);

#endif //MLMC_MC_HPP