21#ifndef OPM_FLEXIBLE_SOLVER_IMPL_HEADER_INCLUDED
22#define OPM_FLEXIBLE_SOLVER_IMPL_HEADER_INCLUDED
24#include <opm/common/ErrorMacros.hpp>
26#include <opm/simulators/linalg/matrixblock.hh>
27#include <opm/simulators/linalg/ilufirstelement.hh>
28#include <opm/simulators/linalg/FlexibleSolver.hpp>
29#include <opm/simulators/linalg/PreconditionerFactory.hpp>
30#include <opm/simulators/linalg/PropertyTree.hpp>
31#include <opm/simulators/linalg/WellOperators.hpp>
33#include <dune/common/fmatrix.hh>
34#include <dune/istl/bcrsmatrix.hh>
35#include <dune/istl/solvers.hh>
36#include <dune/istl/umfpack.hh>
37#include <dune/istl/owneroverlapcopy.hh>
38#include <dune/istl/paamg/pinfo.hh>
43 template <
class Operator>
47 const std::function<VectorType()>& weightsCalculator,
48 std::size_t pressureIndex)
50 init(op, Dune::Amg::SequentialInformation(), prm, weightsCalculator,
55 template <
class Operator>
61 const std::function<VectorType()>& weightsCalculator,
62 std::size_t pressureIndex)
64 init(op, comm, prm, weightsCalculator, pressureIndex);
67 template <
class Operator>
70 apply(VectorType& x, VectorType& rhs, Dune::InverseOperatorResult& res)
72 linsolver_->apply(x, rhs, res);
75 template <
class Operator>
77 FlexibleSolver<Operator>::
78 apply(VectorType& x, VectorType& rhs,
double reduction, Dune::InverseOperatorResult& res)
80 linsolver_->apply(x, rhs, reduction, res);
84 template <
class Operator>
89 return *preconditioner_;
92 template <
class Operator>
93 Dune::SolverCategory::Category
97 return linearoperator_for_solver_->category();
101 template <
class Operator>
102 template <
class Comm>
104 FlexibleSolver<Operator>::
105 initOpPrecSp(Operator& op,
107 const std::function<VectorType()> weightsCalculator,
109 std::size_t pressureIndex)
112 linearoperator_for_solver_ = &op;
113 auto child = prm.get_child_optional(
"preconditioner");
119 scalarproduct_ = Dune::createScalarProduct<VectorType, Comm>(comm, op.category());
122 template <
class Operator>
124 FlexibleSolver<Operator>::
125 initOpPrecSp(Operator& op,
127 const std::function<VectorType()> weightsCalculator,
128 const Dune::Amg::SequentialInformation&,
129 std::size_t pressureIndex)
132 linearoperator_for_solver_ = &op;
133 auto child = prm.get_child_optional(
"preconditioner");
138 scalarproduct_ = std::make_shared<Dune::SeqScalarProduct<VectorType>>();
141 template <
class Operator>
143 FlexibleSolver<Operator>::
146 const double tol = prm.get<
double>(
"tol", 1e-2);
147 const int maxiter = prm.get<
int>(
"maxiter", 200);
148 const int verbosity = is_iorank ? prm.get<
int>(
"verbosity", 0) : 0;
149 const std::string solver_type = prm.get<std::string>(
"solver",
"bicgstab");
150 if (solver_type ==
"bicgstab") {
151 linsolver_ = std::make_shared<Dune::BiCGSTABSolver<VectorType>>(*linearoperator_for_solver_,
157 }
else if (solver_type ==
"loopsolver") {
158 linsolver_ = std::make_shared<Dune::LoopSolver<VectorType>>(*linearoperator_for_solver_,
164 }
else if (solver_type ==
"gmres") {
165 int restart = prm.get<
int>(
"restart", 15);
166 linsolver_ = std::make_shared<Dune::RestartedGMResSolver<VectorType>>(*linearoperator_for_solver_,
173 }
else if (solver_type ==
"flexgmres") {
174 int restart = prm.get<
int>(
"restart", 15);
175 linsolver_ = std::make_shared<Dune::RestartedFlexibleGMResSolver<VectorType>>(*linearoperator_for_solver_,
182#if HAVE_SUITESPARSE_UMFPACK
183 }
else if (solver_type ==
"umfpack") {
185 using MatrixType = std::remove_const_t<std::remove_reference_t<
decltype(linearoperator_for_solver_->getmat())>>;
186 linsolver_ = std::make_shared<Dune::UMFPack<MatrixType>>(linearoperator_for_solver_->getmat(), verbosity, dummy);
189 OPM_THROW(std::invalid_argument,
190 "Properties: Solver " + solver_type +
" not known.");
197 template <
class Operator>
198 template <
class Comm>
200 FlexibleSolver<Operator>::
204 const std::function<VectorType()> weightsCalculator,
205 std::size_t pressureIndex)
207 initOpPrecSp(op, prm, weightsCalculator, comm, pressureIndex);
208 initSolver(prm, comm.communicator().rank() == 0);
218using BV = Dune::BlockVector<Dune::FieldVector<double, N>>;
220using OBM = Dune::BCRSMatrix<Opm::MatrixBlock<double, N, N>>;
224using SeqOpM = Dune::MatrixAdapter<OBM<N>, BV<N>, BV<N>>;
231using Comm = Dune::OwnerOverlapCopyCommunication<int, int>;
233using ParOpM = Dune::OverlappingSchwarzOperator<OBM<N>, BV<N>, BV<N>, Comm>;
241#define INSTANTIATE_FLEXIBLESOLVER_OP(Operator) \
242template class Dune::FlexibleSolver<Operator>; \
243template Dune::FlexibleSolver<Operator>::FlexibleSolver(Operator& op, \
245 const Opm::PropertyTree& prm, \
246 const std::function<typename Operator::domain_type()>& weightsCalculator, \
247 std::size_t pressureIndex);
248#define INSTANTIATE_FLEXIBLESOLVER(N) \
249INSTANTIATE_FLEXIBLESOLVER_OP(SeqOpM<N>); \
250INSTANTIATE_FLEXIBLESOLVER_OP(SeqOpW<N>); \
251INSTANTIATE_FLEXIBLESOLVER_OP(ParOpM<N>); \
252INSTANTIATE_FLEXIBLESOLVER_OP(ParOpW<N>);
256#define INSTANTIATE_FLEXIBLESOLVER_OP(Operator) \
257template class Dune::FlexibleSolver<Operator>;
259#define INSTANTIATE_FLEXIBLESOLVER(N) \
260INSTANTIATE_FLEXIBLESOLVER_OP(SeqOpM<N>); \
261INSTANTIATE_FLEXIBLESOLVER_OP(SeqOpW<N>);
A solver class that encapsulates all needed objects for a linear solver (operator,...
Definition: FlexibleSolver.hpp:45
FlexibleSolver(Operator &op, const Opm::PropertyTree &prm, const std::function< VectorType()> &weightsCalculator, std::size_t pressureIndex)
Create a sequential solver.
Definition: FlexibleSolver_impl.hpp:45
AbstractPrecondType & preconditioner()
Access the contained preconditioner.
Definition: FlexibleSolver_impl.hpp:87
static PrecPtr create(const Operator &op, const PropertyTree &prm, const std::function< Vector()> &weightsCalculator={}, std::size_t pressureIndex=std::numeric_limits< std::size_t >::max())
Create a new serial preconditioner and return a pointer to it.
Definition: PreconditionerFactory_impl.hpp:486
Definition: PropertyTree.hpp:37
Adapter to combine a matrix and another linear operator into a combined linear operator.
Definition: WellOperators.hpp:209
Adapter to combine a matrix and another linear operator into a combined linear operator.
Definition: WellOperators.hpp:123