MultilevelEstimator.cpp

#include "MultilevelEstimator.hpp"


void MultilevelEstimator::initializeCollocationMap() {
}

void MultilevelEstimator::initializeMonteCarloMap() {
  estimatorMap.clear();
  for (unsigned long i = 0; i < initLevels.size(); i++) {
    bool onlyFine = (i == 0);
    int M = initSampleAmount[i];
    estimatorMap.insert(
      {initLevels[i], new MonteCarlo(initLevels[i], M, onlyFine, true)}
    );
  }
}

void MultilevelEstimator::Method() {
  mout.StartBlock("MLMC Method");

  if (epsilon == 0.0) {
    vout(1) << "Non adaptive method" << endl;
    method();
  } else {
    vout(1) << "eps=" << to_string(epsilon) << endl;
    adaptiveMethod();
  }

  data.ExtractDataFrom(estimatorMap);
  exponents.ComputeExponents(data);

  mout.EndBlock(verbose == 0);
}

void MultilevelEstimator::method() {
  for (auto &[level, estimator] : estimatorMap)
    estimator->Method();
}

void MultilevelEstimator::adaptiveMethod() {
  bool converged = false;
  while (!converged) {
    for (auto &[level, estimator] : estimatorMap)
      if (estimator->aggregate.ctr.dM > 0)
        estimator->Method();

    estimatorMap.UpdateSampleCounter(epsilon);
    data.ExtractDataFrom(estimatorMap);
    exponents.ComputeExponents(data);

    vout(1) << exponents;
    vout(1) << "dM=" << vec2str(data.ctrs.GetdMVector()) << endl;

    if (data.ctrs.NoSamplesLeft()) {
      errors.Estimate(data);
      vout(1) << "numErr=" << errors.numeric << endl;

      if (errors.numeric > epsilon / 2) {
        estimatorMap.AppendLevel(epsilon, exponents, estimatorMap.rbegin()->first + 1);
      } else {
        errors.Estimate(data);
        vout(1) << "statErr=" << errors.stochastic << endl;
        converged = true;
      }
    }
  }
}

void MultilevelEstimator::PrintInfo() const {
  mout.PrintInfo("MultilevelEstimator", verbose,
                 PrintInfoEntry("initLevels", vec2str(initLevels)),
                 PrintInfoEntry("initSample", vec2str(initSampleAmount)));
}

void MultilevelEstimator::MultilevelResults() const {
  data.MultilevelResults();
}

void MultilevelEstimator::EstimatorResults() const {
  mout.PrintInfo("MLMC Results", verbose,
                 PrintInfoEntry("Value", data.means.Q),
                 PrintInfoEntry("Cost", data.means.Cost),
                 PrintInfoEntry("Epsilon", epsilon),
                 PrintInfoEntry("Total Error", errors.total),
                 PrintInfoEntry("Statistical Error", errors.stochastic),
                 PrintInfoEntry("Numerical Error", errors.numeric));
}

void MultilevelEstimator::ExponentResults() const {
  mout.PrintInfo("Exponents", verbose,
                 PrintInfoEntry("Final alpha", exponents.alpha),
                 PrintInfoEntry("Final beta", exponents.beta),
                 PrintInfoEntry("Final gamma", exponents.gamma));
}