Commit f3e9bb10 authored by niklas.baumgarten's avatar niklas.baumgarten
Browse files

minor refactoring

parent c393f8ed
......@@ -8,7 +8,7 @@ void MultilevelMonteCarlo::initializeMapMonteCarlo() {
for (unsigned long i = 0; i < initLevels.size(); i++) {
bool onlyFine = (i == 0) || (mcOnly);
int l = initLevels[i], M = initSampleAmount[i];
mapMonteCarlo[l] = createMonteCarlo(l, M, onlyFine);
mcMap[l] = createMonteCarlo(l, M, onlyFine);
}
}
......@@ -19,15 +19,15 @@ MonteCarlo *MultilevelMonteCarlo::createMonteCarlo(int l, int dM, bool onlyFine)
}
void MultilevelMonteCarlo::clearMapMonteCarlo() {
for (auto &iter: mapMonteCarlo) {
for (auto &iter: mcMap) {
delete iter.second;
}
mapMonteCarlo.clear();
mcMap.clear();
}
void MultilevelMonteCarlo::Method() {
mout.StartBlock("MLMC Method");
for (auto &mc : mapMonteCarlo)
for (auto &mc : mcMap)
mc.second->Method();
estimateExponents();
wrapUpResults();
......@@ -40,7 +40,7 @@ void MultilevelMonteCarlo::Method(const double eps) {
vout(1) << "eps=" << to_string(eps) << endl;
bool converged = false;
while (!converged) {
for (auto &mc : mapMonteCarlo)
for (auto &mc : mcMap)
if (mc.second->ctr.dM > 0)
mc.second->Method();
estimateExponents();
......@@ -62,11 +62,11 @@ void MultilevelMonteCarlo::Method(const double eps) {
void MultilevelMonteCarlo::estimateExponents(bool excludeBaseLevel) {
vector<double> levelList, avgsY, varsY, costs;
auto iter = mapMonteCarlo.begin();
auto iter = mcMap.begin();
if (excludeBaseLevel)
iter++;
for (; iter != mapMonteCarlo.end(); iter++) {
for (; iter != mcMap.end(); iter++) {
levelList.push_back((double) iter->first);
avgsY.push_back(-log2(iter->second->avgs.Y));
varsY.push_back(-log2(iter->second->vars.Y));
......@@ -86,23 +86,23 @@ void MultilevelMonteCarlo::estimateExponents(bool excludeBaseLevel) {
void MultilevelMonteCarlo::updateSampleAmount(const double &eps) {
int optimalM;
double factor = 0.0;
for (auto &mc : mapMonteCarlo)
for (auto &mc : mcMap)
factor += sqrt(mc.second->vars.Y * mc.second->avgs.Cost);
for (auto &mc : mapMonteCarlo) {
for (auto &mc : mcMap) {
optimalM = (int) (ceil(2 * pow(eps, -2) * factor *
sqrt(mc.second->vars.Y / mc.second->avgs.Cost)));
if (optimalM == 1) optimalM++; // Hack
mapMonteCarlo[mc.first]->ctr.dM = optimalM - mapMonteCarlo[mc.first]->ctr.M;
mcMap[mc.first]->ctr.dM = optimalM - mcMap[mc.first]->ctr.M;
}
vout(1) << "dM=[";
for (auto &mc:mapMonteCarlo) vout(1) << " " << mc.second->ctr.dM;
for (auto &mc:mcMap) vout(1) << " " << mc.second->ctr.dM;
vout(1) << "]" << endl;
}
bool MultilevelMonteCarlo::checkForNewLevel() {
bool newLevel = true;
for (auto &mc : mapMonteCarlo) {
for (auto &mc : mcMap) {
if (mc.second->ctr.dM > 0) {
newLevel = false;
break;
......@@ -113,9 +113,9 @@ bool MultilevelMonteCarlo::checkForNewLevel() {
void MultilevelMonteCarlo::estimateNumericalError() {
numErr = 0.0;
unsigned long exponent = mapMonteCarlo.size() - 2;
for (auto &mc : mapMonteCarlo) {
if (mc.first == mapMonteCarlo.begin()->first) continue;
unsigned long exponent = mcMap.size() - 2;
for (auto &mc : mcMap) {
if (mc.first == mcMap.begin()->first) continue;
numErr = max(mc.second->avgs.Y / (pow(2.0, alpha) - 1) / (pow(2.0, exponent)),
numErr);
exponent--;
......@@ -125,7 +125,7 @@ void MultilevelMonteCarlo::estimateNumericalError() {
void MultilevelMonteCarlo::estimateStatisticalError() {
statErr = 0.0;
for (auto &mc : mapMonteCarlo) {
for (auto &mc : mcMap) {
statErr += mc.second->vars.Y / mc.second->ctr.M;
}
statErr = sqrt(statErr);
......@@ -133,20 +133,20 @@ void MultilevelMonteCarlo::estimateStatisticalError() {
}
void MultilevelMonteCarlo::appendLevel(double eps) {
int L = mapMonteCarlo.rbegin()->first + 1;
int L = mcMap.rbegin()->first + 1;
if (L > maxLevel) {
Exit ("Maximum level has been reached.")
} else {
vout(1) << "L_new=" << to_string(L) << endl;
mapMonteCarlo[L] = createMonteCarlo(L, 0, false);
mapMonteCarlo[L]->vars.Y = mapMonteCarlo[L - 1]->vars.Y / pow(2.0, beta);
mapMonteCarlo[L]->avgs.Cost = mapMonteCarlo[L - 1]->avgs.Cost * pow(2.0, gamma);
mcMap[L] = createMonteCarlo(L, 0, false);
mcMap[L]->vars.Y = mcMap[L - 1]->vars.Y / pow(2.0, beta);
mcMap[L]->avgs.Cost = mcMap[L - 1]->avgs.Cost * pow(2.0, gamma);
updateSampleAmount(eps);
}
}
void MultilevelMonteCarlo::wrapUpResults() {
for (auto &mc : mapMonteCarlo) {
for (auto &mc : mcMap) {
value += mc.second->avgs.Y;
cost += mc.second->sums.Cost;
......@@ -168,7 +168,7 @@ void MultilevelMonteCarlo::PrintInfo() const {
PrintInfoEntry("maxLevel", maxLevel));
}
void MultilevelMonteCarlo::PrintMCResults() {
void MultilevelMonteCarlo::PrintMCResults() const {
mout.PrintInfo("MC Results", verbose,
PrintInfoEntry("E(Qf-Qc)", vec2str(avgY)),
PrintInfoEntry("E(Qf)", vec2str(avgQ)),
......@@ -180,7 +180,7 @@ void MultilevelMonteCarlo::PrintMCResults() {
PrintInfoEntry("Used Samples", vec2str(samples)));
}
void MultilevelMonteCarlo::PrintMLMCResults() {
void MultilevelMonteCarlo::PrintMLMCResults() const {
mout.PrintInfo("MLMC Results", verbose,
PrintInfoEntry("Value", value),
PrintInfoEntry("Cost", cost),
......@@ -190,7 +190,7 @@ void MultilevelMonteCarlo::PrintMLMCResults() {
PrintInfoEntry("Numerical Error", numErr));
}
void MultilevelMonteCarlo::PrintExponentResults() {
void MultilevelMonteCarlo::PrintExponentResults() const {
mout.PrintInfo("Exponents", verbose,
PrintInfoEntry("Final alpha", alpha),
PrintInfoEntry("Final beta", beta),
......
......@@ -5,7 +5,7 @@
#include "main/Utils.hpp"
typedef std::map<int, MonteCarlo *> MapMonteCarlo;
typedef std::map<int, MonteCarlo *> MonteCarloMap;
class MultilevelMonteCarlo {
private:
......@@ -15,7 +15,7 @@ public:
int maxLevel, pLevel;
vector<int> &initLevels, &initSampleAmount;
MapMonteCarlo mapMonteCarlo;
MonteCarloMap mcMap;
Meshes *meshes;
IStochasticAssemble *assemble;
......@@ -67,17 +67,15 @@ public:
void PrintInfo() const;
void PrintMCResults();
void PrintMCResults() const;
void PrintMLMCResults();
void PrintMLMCResults() const;
void PrintExponentResults();
void PrintExponentResults() const;
private:
void initializeMapMonteCarlo();
void initializeValues();
MonteCarlo *createMonteCarlo(int l, int dM, bool onlyFine);
void clearMapMonteCarlo();
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment