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
#ifdef HAVE_CONFIG_H
#include "config.h"
#endif
#ifdef HAVE_IPOPT
#undef HAVE_IPOPT
#endif
#include <dune/fufem/assemblers/transferoperatorassembler.hh>
#include <dune/solvers/solvers/solver.hh>
#include "solverfactory.hh"
template <size_t dim, class BlockProblem, class Grid>
SolverFactory<dim, BlockProblem, Grid>::SolverFactory(
Dune::ParameterTree const &parset, Grid const &grid,
Dune::BitSetVector<dim> const &ignoreNodes)
: baseEnergyNorm(linearBaseSolverStep),
linearBaseSolver(&linearBaseSolverStep,
parset.get<size_t>("linear.maxiumumIterations"),
parset.get<double>("linear.tolerance"), &baseEnergyNorm,
Solver::QUIET),
transferOperators(grid.maxLevel()),
multigridStep(
std::make_shared<Step>(linearIterationStep, nonlinearSmoother)) {
// linear iteration step
linearIterationStep.setMGType(parset.get<int>("linear.cycle"),
parset.get<int>("linear.pre"),
parset.get<int>("linear.post"));
linearIterationStep.basesolver_ = &linearBaseSolver;
linearIterationStep.setSmoother(&linearPresmoother, &linearPostsmoother);
// transfer operators
for (auto &&x : transferOperators)
x = new CompressedMultigridTransfer<Vector>;
TransferOperatorAssembler<Grid>(grid)
.assembleOperatorPointerHierarchy(transferOperators);
linearIterationStep.setTransferOperators(transferOperators);
// tnnmg iteration step
multigridStep->setSmoothingSteps(parset.get<int>("main.pre"),
parset.get<int>("main.multi"),
parset.get<int>("main.post"));
multigridStep->ignoreNodes_ = &ignoreNodes;
}
template <size_t dim, class BlockProblem, class Grid>
SolverFactory<dim, BlockProblem, Grid>::~SolverFactory() {
for (auto &&x : transferOperators)
delete x;
}
template <size_t dim, class BlockProblem, class Grid>
auto SolverFactory<dim, BlockProblem, Grid>::getStep()
-> std::shared_ptr<Step> {
return multigridStep;
}
#include "solverfactory_tmpl.cc"