20 #ifndef OPM_OWNINGTWOLEVELPRECONDITIONER_HEADER_INCLUDED
21 #define OPM_OWNINGTWOLEVELPRECONDITIONER_HEADER_INCLUDED
23 #include <opm/simulators/linalg/PreconditionerWithUpdate.hpp>
24 #include <opm/simulators/linalg/PressureSolverPolicy.hpp>
27 #include <opm/common/ErrorMacros.hpp>
29 #include <dune/common/fmatrix.hh>
30 #include <dune/istl/bcrsmatrix.hh>
31 #include <dune/istl/paamg/amg.hh>
34 #include <type_traits>
42 template <
class Operator,
class Comm = Dune::Amg::SequentialInformation>
43 class PreconditionerFactory;
51 template <
class Operator>
54 template <
typename T,
typename A,
int i>
55 std::ostream& operator<<(std::ostream& out,
56 const BlockVector< FieldVector< T, i >, A >& vector)
58 Dune::writeMatrixMarket(vector, out);
62 template <
typename T,
typename A,
int i>
63 std::istream& operator>>(std::istream& input,
64 BlockVector< FieldVector< T, i >, A >& vector)
66 Dune::readMatrixMarket(vector, input);
74 template <
class OperatorType,
76 class LevelTransferPolicy,
77 class Communication = Dune::Amg::SequentialInformation>
81 using MatrixType =
typename OperatorType::matrix_type;
83 using AbstractOperatorType = Dune::AssembledLinearOperator<MatrixType, VectorType, VectorType>;
86 const std::function<VectorType()> weightsCalculator,
87 std::size_t pressureIndex)
88 : linear_operator_(linearoperator)
90 prm.get_child_optional(
"finesmoother") ?
92 std::function<VectorType()>(), pressureIndex))
94 , weightsCalculator_(weightsCalculator)
95 , weights_(weightsCalculator())
96 , levelTransferPolicy_(dummy_comm_, weights_, prm, pressureIndex)
97 , coarseSolverPolicy_(prm.get_child_optional(
"coarsesolver") ? prm.get_child(
"coarsesolver") :
Opm::PropertyTree())
98 , twolevel_method_(linearoperator,
100 levelTransferPolicy_,
102 prm.get<
int>(
"pre_smooth", 0),
103 prm.get<
int>(
"post_smooth", 1))
106 if (prm.get<
int>(
"verbosity", 0) > 10) {
107 std::string filename = prm.get<std::string>(
"weights_filename",
"impes_weights.txt");
108 std::ofstream outfile(filename);
110 OPM_THROW(std::ofstream::failure,
"Could not write weights to file " << filename <<
".");
112 Dune::writeMatrixMarket(weights_, outfile);
117 const std::function<VectorType()> weightsCalculator,
118 std::size_t pressureIndex,
const Communication& comm)
119 : linear_operator_(linearoperator)
121 prm.get_child_optional(
"finesmoother") ?
123 std::function<VectorType()>(),
124 comm, pressureIndex))
126 , weightsCalculator_(weightsCalculator)
127 , weights_(weightsCalculator())
128 , levelTransferPolicy_(*comm_, weights_, prm, pressureIndex)
129 , coarseSolverPolicy_(prm.get_child_optional(
"coarsesolver") ? prm.get_child(
"coarsesolver") :
Opm::PropertyTree())
130 , twolevel_method_(linearoperator,
132 levelTransferPolicy_,
134 prm.get<
int>(
"pre_smooth", 0),
135 prm.get<
int>(
"post_smooth", 1))
138 if (prm.get<
int>(
"verbosity", 0) > 10 && comm.communicator().rank() == 0) {
139 auto filename = prm.get<std::string>(
"weights_filename",
"impes_weights.txt");
140 std::ofstream outfile(filename);
142 OPM_THROW(std::ofstream::failure,
"Could not write weights to file " << filename <<
".");
144 Dune::writeMatrixMarket(weights_, outfile);
148 virtual void pre(VectorType& x, VectorType& b)
override
150 twolevel_method_.pre(x, b);
153 virtual void apply(VectorType& v,
const VectorType& d)
override
155 twolevel_method_.apply(v, d);
158 virtual void post(VectorType& x)
override
160 twolevel_method_.post(x);
163 virtual void update()
override
165 weights_ = weightsCalculator_();
169 virtual Dune::SolverCategory::Category category()
const override
171 return linear_operator_.category();
175 using CoarseOperator =
typename LevelTransferPolicy::CoarseOperator;
178 LevelTransferPolicy>;
183 template <
class Comm>
184 void updateImpl(
const Comm*)
187 auto child = prm_.get_child_optional(
"finesmoother");
189 twolevel_method_.updatePreconditioner(finesmoother_, coarseSolverPolicy_);
192 void updateImpl(
const Dune::Amg::SequentialInformation*)
195 auto child = prm_.get_child_optional(
"finesmoother");
197 twolevel_method_.updatePreconditioner(finesmoother_, coarseSolverPolicy_);
200 const OperatorType& linear_operator_;
201 std::shared_ptr<Dune::Preconditioner<VectorType, VectorType>> finesmoother_;
202 const Communication* comm_;
203 std::function<VectorType()> weightsCalculator_;
205 LevelTransferPolicy levelTransferPolicy_;
209 Communication dummy_comm_;
A solver class that encapsulates all needed objects for a linear solver (operator,...
Definition: FlexibleSolver.hpp:45
A version of the two-level preconditioner that is:
Definition: OwningTwoLevelPreconditioner.hpp:79
Interface class adding the update() method to the preconditioner interface.
Definition: PreconditionerWithUpdate.hpp:32
This is an object factory for creating preconditioners.
Definition: PreconditionerFactory.hpp:69
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:500
Definition: PropertyTree.hpp:37
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27
Algebraic twolevel methods.