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>
25 #include <opm/simulators/linalg/PressureTransferPolicy.hpp>
26 #include <opm/simulators/linalg/getQuasiImpesWeights.hpp>
29 #include <opm/common/ErrorMacros.hpp>
31 #include <dune/common/fmatrix.hh>
32 #include <dune/istl/bcrsmatrix.hh>
33 #include <dune/istl/paamg/amg.hh>
36 #include <type_traits>
44 template <
class Operator,
class Comm = Dune::Amg::SequentialInformation>
45 class PreconditionerFactory;
53 template <
class MatrixTypeT,
class VectorTypeT>
56 template <
typename T,
typename A,
int i>
57 std::ostream& operator<<(std::ostream& out,
58 const BlockVector< FieldVector< T, i >, A >& vector)
60 Dune::writeMatrixMarket(vector, out);
64 template <
typename T,
typename A,
int i>
65 std::istream& operator>>(std::istream& input,
66 BlockVector< FieldVector< T, i >, A >& vector)
68 Dune::readMatrixMarket(vector, input);
76 template <
class OperatorType,
78 bool transpose =
false,
79 class Communication = Dune::Amg::SequentialInformation>
83 using MatrixType =
typename OperatorType::matrix_type;
87 const std::function<VectorType()> weightsCalculator,
88 std::size_t pressureIndex)
89 : linear_operator_(linearoperator)
91 prm.get_child_optional(
"finesmoother") ?
93 std::function<VectorType()>(), pressureIndex))
95 , weightsCalculator_(weightsCalculator)
96 , weights_(weightsCalculator())
97 , levelTransferPolicy_(dummy_comm_, weights_, pressureIndex)
98 , coarseSolverPolicy_(prm.get_child_optional(
"coarsesolver")? prm.get_child(
"coarsesolver") :
Opm::PropertyTree())
99 , twolevel_method_(linearoperator,
101 levelTransferPolicy_,
103 prm.get<
int>(
"pre_smooth", transpose? 1 : 0),
104 prm.get<
int>(
"post_smooth", transpose? 0 : 1))
107 if (prm.get<
int>(
"verbosity", 0) > 10) {
108 std::string filename = prm.get<std::string>(
"weights_filename",
"impes_weights.txt");
109 std::ofstream outfile(filename);
111 OPM_THROW(std::ofstream::failure,
"Could not write weights to file " << filename <<
".");
113 Dune::writeMatrixMarket(weights_, outfile);
118 const std::function<VectorType()> weightsCalculator,
119 std::size_t pressureIndex,
const Communication& comm)
120 : linear_operator_(linearoperator)
122 prm.get_child_optional(
"finesmoother") ?
124 std::function<VectorType()>(),
125 comm, pressureIndex))
127 , weightsCalculator_(weightsCalculator)
128 , weights_(weightsCalculator())
129 , levelTransferPolicy_(*comm_, weights_, pressureIndex)
130 , coarseSolverPolicy_(prm.get_child_optional(
"coarsesolver")? prm.get_child(
"coarsesolver") :
Opm::PropertyTree())
131 , twolevel_method_(linearoperator,
133 levelTransferPolicy_,
135 prm.get<
int>(
"pre_smooth", transpose? 1 : 0),
136 prm.get<
int>(
"post_smooth", transpose? 0 : 1))
139 if (prm.get<
int>(
"verbosity", 0) > 10 && comm.communicator().rank() == 0) {
140 auto filename = prm.get<std::string>(
"weights_filename",
"impes_weights.txt");
141 std::ofstream outfile(filename);
143 OPM_THROW(std::ofstream::failure,
"Could not write weights to file " << filename <<
".");
145 Dune::writeMatrixMarket(weights_, outfile);
149 virtual void pre(VectorType& x, VectorType& b)
override
151 twolevel_method_.pre(x, b);
154 virtual void apply(VectorType& v,
const VectorType& d)
override
156 twolevel_method_.apply(v, d);
159 virtual void post(VectorType& x)
override
161 twolevel_method_.post(x);
164 virtual void update()
override
166 weights_ = weightsCalculator_();
170 virtual Dune::SolverCategory::Category category()
const override
172 return linear_operator_.category();
176 using PressureMatrixType = Dune::BCRSMatrix<Dune::FieldMatrix<double, 1, 1>>;
177 using PressureVectorType = Dune::BlockVector<Dune::FieldVector<double, 1>>;
178 using SeqCoarseOperatorType = Dune::MatrixAdapter<PressureMatrixType, PressureVectorType, PressureVectorType>;
179 using ParCoarseOperatorType
180 = Dune::OverlappingSchwarzOperator<PressureMatrixType, PressureVectorType, PressureVectorType, Communication>;
181 using CoarseOperatorType = std::conditional_t<std::is_same<Communication, Dune::Amg::SequentialInformation>::value,
182 SeqCoarseOperatorType,
183 ParCoarseOperatorType>;
192 template <
class Comm>
193 void updateImpl(
const Comm*)
196 auto child = prm_.get_child_optional(
"finesmoother");
198 twolevel_method_.updatePreconditioner(finesmoother_, coarseSolverPolicy_);
201 void updateImpl(
const Dune::Amg::SequentialInformation*)
204 auto child = prm_.get_child_optional(
"finesmoother");
206 twolevel_method_.updatePreconditioner(finesmoother_, coarseSolverPolicy_);
209 const OperatorType& linear_operator_;
210 std::shared_ptr<Dune::Preconditioner<VectorType, VectorType>> finesmoother_;
211 const Communication* comm_;
212 std::function<VectorType()> weightsCalculator_;
218 Communication dummy_comm_;
A solver class that encapsulates all needed objects for a linear solver (operator,...
Definition: FlexibleSolver.hpp:40
A version of the two-level preconditioner that is:
Definition: OwningTwoLevelPreconditioner.hpp:81
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:51
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.hpp:71
Definition: PressureTransferPolicy.hpp:33
Definition: PropertyTree.hpp:37
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27
Algebraic twolevel methods.