22#ifndef OPM_PRECONDITIONERFACTORY_HEADER
23#define OPM_PRECONDITIONERFACTORY_HEADER
24#include <opm/common/TimingMacros.hpp>
25#include <opm/simulators/linalg/PreconditionerWithUpdate.hpp>
27#include <dune/istl/paamg/aggregates.hh>
28#include <dune/istl/paamg/matrixhierarchy.hh>
41template <
class Operator,
class Comm,
class Matrix,
class Vector>
44 using PrecPtr = std::shared_ptr<Dune::PreconditionerWithUpdate<Vector, Vector>>;
46 = Dune::Amg::AggregationCriterion<Dune::Amg::SymmetricDependency<Matrix, Dune::Amg::FirstDiagonal>>;
47 using Criterion = Dune::Amg::CoarsenCriterion<CriterionBase>;
51 template <
class Smoother>
52 static PrecPtr makeAmgPreconditioner(
const Operator&
op,
62template <
class Operator,
class Comm>
67 using Matrix =
typename Operator::matrix_type;
68 using Vector =
typename Operator::domain_type;
71 using PrecPtr = std::shared_ptr<Dune::PreconditionerWithUpdate<Vector, Vector>>;
75 const std::function<Vector()>&, std::size_t)>;
77 const std::function<Vector()>&, std::size_t,
const Comm&)>;
86 std::size_t pressureIndex = std::numeric_limits<std::size_t>::max());
96 std::size_t pressureIndex = std::numeric_limits<std::size_t>::max());
104 std::size_t pressureIndex = std::numeric_limits<std::size_t>::max());
125 = Dune::Amg::AggregationCriterion<Dune::Amg::SymmetricDependency<Matrix, Dune::Amg::FirstDiagonal>>;
126 using Criterion = Dune::Amg::CoarsenCriterion<CriterionBase>;
140 std::size_t pressureIndex);
144 std::size_t pressureIndex,
const Comm& comm);
147 void doAddCreator(
const std::string& type,
Creator c);
150 void doAddCreator(
const std::string& type, ParCreator
c);
153 std::map<std::string, Creator> creators_;
154 std::map<std::string, ParCreator> parallel_creators_;
155 bool defAdded_=
false;
This is an object factory for creating preconditioners.
Definition PreconditionerFactory.hpp:64
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:774
std::function< PrecPtr(const Operator &, const PropertyTree &, const std::function< Vector()> &, std::size_t)> Creator
The type of creator functions passed to addCreator().
Definition PreconditionerFactory.hpp:75
static void addCreator(const std::string &type, Creator creator)
Add a creator for a serial preconditioner to the PreconditionerFactory.
Definition PreconditionerFactory_impl.hpp:806
std::shared_ptr< Dune::PreconditionerWithUpdate< Vector, Vector > > PrecPtr
The type of pointer returned by create().
Definition PreconditionerFactory.hpp:71
typename Operator::matrix_type Matrix
Linear algebra types.
Definition PreconditionerFactory.hpp:67
Hierarchical collection of key/value pairs.
Definition PropertyTree.hpp:39
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
Definition PreconditionerFactory.hpp:43