Newer
Older
#ifdef HAVE_CONFIG_H
#include "config.h"
#endif
#include <dune/common/exceptions.hh>
#include <dune/solvers/common/arithmetic.hh>
#include <dune/solvers/norms/energynorm.hh>
#include <dune/solvers/solvers/loopsolver.hh>
#include "enums.hh"
#include "enumparser.hh"
#include "fixedpointiterator.hh"
template <class Factory, class StateUpdater, class VelocityUpdater,
class ErrorNorm>
FixedPointIterator<Factory, StateUpdater, VelocityUpdater, ErrorNorm>::
FixedPointIterator(Factory &factory, Dune::ParameterTree const &parset,
std::shared_ptr<Nonlinearity> globalFriction,
ErrorNorm const &errorNorm)
: factory_(factory),
parset_(parset),
globalFriction_(globalFriction),
fixedPointMaxIterations_(parset.get<size_t>("v.fpi.maximumIterations")),
fixedPointTolerance_(parset.get<double>("v.fpi.tolerance")),
lambda_(parset.get<double>("v.fpi.lambda")),
velocityMaxIterations_(parset.get<size_t>("v.solver.maximumIterations")),
velocityTolerance_(parset.get<double>("v.solver.tolerance")),
verbosity_(parset.get<Solver::VerbosityMode>("v.solver.verbosity")),
errorNorm_(errorNorm) {}
template <class Factory, class StateUpdater, class VelocityUpdater,
class ErrorNorm>
int FixedPointIterator<Factory, StateUpdater, VelocityUpdater, ErrorNorm>::run(
std::shared_ptr<StateUpdater> stateUpdater,
std::shared_ptr<VelocityUpdater> velocityUpdater,
Matrix const &velocityMatrix, Vector const &velocityRHS,
Vector &velocityIterate) {
auto multigridStep = factory_.getSolver();
EnergyNorm<Matrix, Vector> energyNorm(velocityMatrix);
LoopSolver<Vector> velocityProblemSolver(
multigridStep, velocityMaxIterations_, velocityTolerance_, &energyNorm,
verbosity_, false); // absolute error
ScalarVector alpha;
stateUpdater->extractAlpha(alpha);
for (fixedPointIteration = 0; fixedPointIteration < fixedPointMaxIterations_;
++fixedPointIteration) {
// solve a velocity problem
globalFriction_->updateAlpha(alpha);
ConvexProblem convexProblem(1.0, velocityMatrix, *globalFriction_,
velocityRHS, velocityIterate);
BlockProblem velocityProblem(parset_, convexProblem);
multigridStep->setProblem(velocityIterate, velocityProblem);
velocityProblemSolver.preprocess();
velocityProblemSolver.solve();
Vector v_m;
velocityUpdater->extractOldVelocity(v_m);
v_m *= 1.0 - lambda_;
Arithmetic::addProduct(v_m, lambda_, velocityIterate);
// solve a state problem
stateUpdater->solve(v_m);
ScalarVector newAlpha;
stateUpdater->extractAlpha(newAlpha);
if (errorNorm_.diff(alpha, newAlpha) < fixedPointTolerance_) {
}
if (fixedPointIteration == fixedPointMaxIterations_)
DUNE_THROW(Dune::Exception, "FPI failed to converge");
velocityUpdater->postProcess(velocityIterate);
return fixedPointIteration;
}
#include "fixedpointiterator_tmpl.cc"