My Project
PreconditionerFactory.hpp
1
2/*
3 Copyright 2009, 2010 SINTEF ICT, Applied Mathematics.
4 Copyright 2019 SINTEF Digital, Mathematics and Cybernetics.
5
6 This file is part of the Open Porous Media project (OPM).
7
8 OPM is free software: you can redistribute it and/or modify
9 it under the terms of the GNU General Public License as published by
10 the Free Software Foundation, either version 3 of the License, or
11 (at your option) any later version.
12
13 OPM is distributed in the hope that it will be useful,
14 but WITHOUT ANY WARRANTY; without even the implied warranty of
15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 GNU General Public License for more details.
17
18 You should have received a copy of the GNU General Public License
19 along with OPM. If not, see <http://www.gnu.org/licenses/>.
20*/
21
22#ifndef OPM_PRECONDITIONERFACTORY_HEADER
23#define OPM_PRECONDITIONERFACTORY_HEADER
24
25#include <opm/simulators/linalg/PreconditionerWithUpdate.hpp>
26
27#include <dune/common/version.hh>
28#include <dune/istl/paamg/aggregates.hh>
29#if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
30#include <dune/istl/paamg/matrixhierarchy.hh>
31#else
32#include <dune/istl/paamg/hierarchy.hh>
33#endif
34
35#include <cstddef>
36#include <map>
37#include <memory>
38#include <limits>
39#include <string>
40
41namespace Opm
42{
43
44class PropertyTree;
45
46template <class Operator, class Comm, class Matrix, class Vector>
48{
49 using PrecPtr = std::shared_ptr<Dune::PreconditionerWithUpdate<Vector, Vector>>;
50 using CriterionBase
51 = Dune::Amg::AggregationCriterion<Dune::Amg::SymmetricDependency<Matrix, Dune::Amg::FirstDiagonal>>;
52 using Criterion = Dune::Amg::CoarsenCriterion<CriterionBase>;
53
54 static Criterion criterion(const PropertyTree& prm);
55
56 template <class Smoother>
57 static PrecPtr makeAmgPreconditioner(const Operator& op,
58 const PropertyTree& prm,
59 bool useKamg = false);
60};
61
67template <class Operator, class Comm>
69{
70public:
72 using Matrix = typename Operator::matrix_type;
73 using Vector = typename Operator::domain_type; // Assuming symmetry: that domain and range types are the same.
74
76 using PrecPtr = std::shared_ptr<Dune::PreconditionerWithUpdate<Vector, Vector>>;
77
79 using Creator = std::function<PrecPtr(const Operator&, const PropertyTree&,
80 const std::function<Vector()>&, std::size_t)>;
81 using ParCreator = std::function<PrecPtr(const Operator&, const PropertyTree&,
82 const std::function<Vector()>&, std::size_t, const Comm&)>;
83
89 static PrecPtr create(const Operator& op, const PropertyTree& prm,
90 const std::function<Vector()>& weightsCalculator = {},
91 std::size_t pressureIndex = std::numeric_limits<std::size_t>::max());
92
99 static PrecPtr create(const Operator& op, const PropertyTree& prm,
100 const std::function<Vector()>& weightsCalculator, const Comm& comm,
101 std::size_t pressureIndex = std::numeric_limits<std::size_t>::max());
102
108 static PrecPtr create(const Operator& op, const PropertyTree& prm, const Comm& comm,
109 std::size_t pressureIndex = std::numeric_limits<std::size_t>::max());
110
118 static void addCreator(const std::string& type, Creator creator);
119
127 static void addCreator(const std::string& type, ParCreator creator);
128
129 using CriterionBase
130 = Dune::Amg::AggregationCriterion<Dune::Amg::SymmetricDependency<Matrix, Dune::Amg::FirstDiagonal>>;
131 using Criterion = Dune::Amg::CoarsenCriterion<CriterionBase>;
132
133
134private:
135 // The method that implements the singleton pattern,
136 // using the Meyers singleton technique.
137 static PreconditionerFactory& instance();
138
139 // Private constructor, to keep users from creating a PreconditionerFactory.
141
142 // Actually creates the product object.
143 PrecPtr doCreate(const Operator& op, const PropertyTree& prm,
144 const std::function<Vector()> weightsCalculator,
145 std::size_t pressureIndex);
146
147 PrecPtr doCreate(const Operator& op, const PropertyTree& prm,
148 const std::function<Vector()> weightsCalculator,
149 std::size_t pressureIndex, const Comm& comm);
150
151 // Actually adds the creator.
152 void doAddCreator(const std::string& type, Creator c);
153
154 // Actually adds the creator.
155 void doAddCreator(const std::string& type, ParCreator c);
156
157 // This map contains the whole factory, i.e. all the Creators.
158 std::map<std::string, Creator> creators_;
159 std::map<std::string, ParCreator> parallel_creators_;
160 bool defAdded_= false;
161};
162
163} // namespace Dune
164
165#endif // OPM_PRECONDITIONERFACTORY_HEADER
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
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:80
static void addCreator(const std::string &type, Creator creator)
Add a creator for a serial preconditioner to the PreconditionerFactory.
Definition: PreconditionerFactory_impl.hpp:529
std::shared_ptr< Dune::PreconditionerWithUpdate< Vector, Vector > > PrecPtr
The type of pointer returned by create().
Definition: PreconditionerFactory.hpp:76
typename Operator::matrix_type Matrix
Linear algebra types.
Definition: PreconditionerFactory.hpp:72
Definition: PropertyTree.hpp:37
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27
Definition: PreconditionerFactory.hpp:48