My Project
Loading...
Searching...
No Matches
PreconditionerFactory_impl.hpp
1/*
2 Copyright 2009, 2010 SINTEF ICT, Applied Mathematics.
3 Copyright 2019 SINTEF Digital, Mathematics and Cybernetics.
4
5 This file is part of the Open Porous Media project (OPM).
6
7 OPM is free software: you can redistribute it and/or modify
8 it under the terms of the GNU General Public License as published by
9 the Free Software Foundation, either version 3 of the License, or
10 (at your option) any later version.
11
12 OPM is distributed in the hope that it will be useful,
13 but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 GNU General Public License for more details.
16
17 You should have received a copy of the GNU General Public License
18 along with OPM. If not, see <http://www.gnu.org/licenses/>.
19*/
20
21#include <opm/simulators/linalg/PreconditionerFactory.hpp>
22
23#include <opm/common/ErrorMacros.hpp>
24
25#include <opm/simulators/linalg/matrixblock.hh>
27#include <opm/simulators/linalg/FlexibleSolver.hpp>
28#include <opm/simulators/linalg/ilufirstelement.hh>
29#include <opm/simulators/linalg/OwningBlockPreconditioner.hpp>
30#include <opm/simulators/linalg/OwningTwoLevelPreconditioner.hpp>
31#include <opm/simulators/linalg/ParallelOverlappingILU0.hpp>
32#include <opm/simulators/linalg/PressureBhpTransferPolicy.hpp>
33#include <opm/simulators/linalg/PressureTransferPolicy.hpp>
34#include <opm/simulators/linalg/PropertyTree.hpp>
35#include <opm/simulators/linalg/WellOperators.hpp>
36
37#include <dune/istl/owneroverlapcopy.hh>
38#include <dune/istl/preconditioners.hh>
39#include <dune/istl/paamg/amg.hh>
40#include <dune/istl/paamg/kamg.hh>
41#include <dune/istl/paamg/fastamg.hh>
42
43namespace Opm {
44
45template<class Smoother>
47{
48 static auto args(const PropertyTree& prm)
49 {
50 using SmootherArgs = typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
51 SmootherArgs smootherArgs;
52 smootherArgs.iterations = prm.get<int>("iterations", 1);
53 // smootherArgs.overlap=SmootherArgs::vertex;
54 // smootherArgs.overlap=SmootherArgs::none;
55 // smootherArgs.overlap=SmootherArgs::aggregate;
56 smootherArgs.relaxationFactor = prm.get<double>("relaxation", 1.0);
57 return smootherArgs;
58 }
59};
60
61template<class M, class V, class C>
63{
64 static auto args(const PropertyTree& prm)
65 {
67 using SmootherArgs = typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
68 SmootherArgs smootherArgs;
69 smootherArgs.iterations = prm.get<int>("iterations", 1);
70 const int iluwitdh = prm.get<int>("iluwidth", 0);
72 const MILU_VARIANT milu = convertString2Milu(prm.get<std::string>("milutype", std::string("ilu")));
73 smootherArgs.setMilu(milu);
74 // smootherArgs.overlap=SmootherArgs::vertex;
75 // smootherArgs.overlap=SmootherArgs::none;
76 // smootherArgs.overlap=SmootherArgs::aggregate;
77 smootherArgs.relaxationFactor = prm.get<double>("relaxation", 1.0);
78 return smootherArgs;
79 }
80};
81
82template <class Operator, class Comm, class Matrix, class Vector>
83typename AMGHelper<Operator, Comm, Matrix, Vector>::Criterion
85{
86 Criterion criterion(15, prm.get<int>("coarsenTarget", 1200));
87 criterion.setDefaultValuesIsotropic(2);
88 criterion.setAlpha(prm.get<double>("alpha", 0.33));
89 criterion.setBeta(prm.get<double>("beta", 1e-5));
90 criterion.setMaxLevel(prm.get<int>("maxlevel", 15));
91 criterion.setSkipIsolated(prm.get<bool>("skip_isolated", false));
92 criterion.setNoPreSmoothSteps(prm.get<int>("pre_smooth", 1));
93 criterion.setNoPostSmoothSteps(prm.get<int>("post_smooth", 1));
94 criterion.setDebugLevel(prm.get<int>("verbosity", 0));
95 // As the default we request to accumulate data to 1 process always as our matrix
96 // graph might be unsymmetric and hence not supported by the PTScotch/ParMetis
97 // calls in DUNE. Accumulating to 1 skips PTScotch/ParMetis
98 criterion.setAccumulate(static_cast<Dune::Amg::AccumulationMode>(prm.get<int>("accumulate", 1)));
99 criterion.setProlongationDampingFactor(prm.get<double>("prolongationdamping", 1.6));
100 criterion.setMaxDistance(prm.get<int>("maxdistance", 2));
101 criterion.setMaxConnectivity(prm.get<int>("maxconnectivity", 15));
102 criterion.setMaxAggregateSize(prm.get<int>("maxaggsize", 6));
103 criterion.setMinAggregateSize(prm.get<int>("minaggsize", 4));
104 return criterion;
105}
106
107template <class Operator, class Comm, class Matrix, class Vector>
108template <class Smoother>
109typename AMGHelper<Operator, Comm, Matrix, Vector>::PrecPtr
110AMGHelper<Operator,Comm,Matrix,Vector>::
111makeAmgPreconditioner(const Operator& op,
112 const PropertyTree& prm,
113 bool useKamg)
114{
115 auto crit = criterion(prm);
116 auto sargs = AMGSmootherArgsHelper<Smoother>::args(prm);
117 if (useKamg) {
119 return std::make_shared<Type>(op, crit, sargs,
120 prm.get<size_t>("max_krylov", 1),
121 prm.get<double>("min_reduction", 1e-1));
122 } else {
124 return std::make_shared<Type>(op, crit, sargs);
125 }
126}
127
128template<class Operator, class Comm>
130{
131 static void add()
132 {
133 using namespace Dune;
134 using O = Operator;
135 using C = Comm;
137 using M = typename F::Matrix;
138 using V = typename F::Vector;
139 using P = PropertyTree;
140 F::addCreator("ILU0", [](const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
141 return createParILU(op, prm, comm, 0);
142 });
143 F::addCreator("ParOverILU0", [](const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
144 return createParILU(op, prm, comm, prm.get<int>("ilulevel", 0));
145 });
146 F::addCreator("ILUn", [](const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
147 return createParILU(op, prm, comm, prm.get<int>("ilulevel", 0));
148 });
149 F::addCreator("Jac", [](const O& op, const P& prm, const std::function<V()>&,
150 std::size_t, const C& comm) {
151 const int n = prm.get<int>("repeats", 1);
152 const double w = prm.get<double>("relaxation", 1.0);
154 });
155 F::addCreator("GS", [](const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
156 const int n = prm.get<int>("repeats", 1);
157 const double w = prm.get<double>("relaxation", 1.0);
159 });
160 F::addCreator("SOR", [](const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
161 const int n = prm.get<int>("repeats", 1);
162 const double w = prm.get<double>("relaxation", 1.0);
164 });
165 F::addCreator("SSOR", [](const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
166 const int n = prm.get<int>("repeats", 1);
167 const double w = prm.get<double>("relaxation", 1.0);
169 });
170
171 // Only add AMG preconditioners to the factory if the operator
172 // is the overlapping schwarz operator. This could be extended
173 // later, but at this point no other operators are compatible
174 // with the AMG hierarchy construction.
175 if constexpr (std::is_same_v<O, Dune::OverlappingSchwarzOperator<M, V, V, C>>) {
176 F::addCreator("amg", [](const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
177 const std::string smoother = prm.get<std::string>("smoother", "ParOverILU0");
178 if (smoother == "ILU0" || smoother == "ParOverILU0") {
182 return std::make_shared<Dune::Amg::AMGCPR<O, V, Smoother, C>>(op, crit, sargs, comm);
183 } else {
184 OPM_THROW(std::invalid_argument, "Properties: No smoother with name " + smoother + ".");
185 }
186 });
187 }
188
189 F::addCreator("cpr", [](const O& op, const P& prm, const std::function<V()> weightsCalculator, std::size_t pressureIndex, const C& comm) {
191 if (pressureIndex == std::numeric_limits<std::size_t>::max())
192 {
193 OPM_THROW(std::logic_error, "Pressure index out of bounds. It needs to specified for CPR");
194 }
196 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy, Comm>>(op, prm, weightsCalculator, pressureIndex, comm);
197 });
198 F::addCreator("cprt", [](const O& op, const P& prm, const std::function<V()> weightsCalculator, std::size_t pressureIndex, const C& comm) {
200 if (pressureIndex == std::numeric_limits<std::size_t>::max())
201 {
202 OPM_THROW(std::logic_error, "Pressure index out of bounds. It needs to specified for CPR");
203 }
205 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy, Comm>>(op, prm, weightsCalculator, pressureIndex, comm);
206 });
207
208 if constexpr (std::is_same_v<O, WellModelGhostLastMatrixAdapter<M, V, V, true>>) {
209 F::addCreator("cprw",
210 [](const O& op, const P& prm, const std::function<V()> weightsCalculator, std::size_t pressureIndex, const C& comm) {
212 if (pressureIndex == std::numeric_limits<std::size_t>::max()) {
213 OPM_THROW(std::logic_error, "Pressure index out of bounds. It needs to specified for CPR");
214 }
216 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy, Comm>>(op, prm, weightsCalculator, pressureIndex, comm);
217 });
218 }
219 }
220
221
223 createParILU(const Operator& op, const PropertyTree& prm, const Comm& comm, const int ilulevel)
224 {
226 using M = typename F::Matrix;
227 using V = typename F::Vector;
228
229 const double w = prm.get<double>("relaxation", 1.0);
230 const bool redblack = prm.get<bool>("redblack", false);
231 const bool reorder_spheres = prm.get<bool>("reorder_spheres", false);
232 // Already a parallel preconditioner. Need to pass comm, but no need to wrap it in a BlockPreconditioner.
233 if (ilulevel == 0) {
234 const size_t num_interior = interiorIfGhostLast(comm);
235 return std::make_shared<Opm::ParallelOverlappingILU0<M, V, V, Comm>>(
237 } else {
238 return std::make_shared<Opm::ParallelOverlappingILU0<M, V, V, Comm>>(
240 }
241 }
242
247 static size_t interiorIfGhostLast(const Comm& comm)
248 {
249 size_t interior_count = 0;
250 size_t highest_interior_index = 0;
251 const auto& is = comm.indexSet();
252 for (const auto& ind : is) {
253 if (Comm::OwnerSet::contains(ind.local().attribute())) {
255 highest_interior_index = std::max(highest_interior_index, ind.local().local());
256 }
257 }
259 return interior_count;
260 } else {
261 return is.size();
262 }
263 }
264
265};
266
267template<class Operator>
268struct StandardPreconditioners<Operator,Dune::Amg::SequentialInformation>
269{
270 static void add()
271 {
272 using namespace Dune;
273 using O = Operator;
274 using C = Dune::Amg::SequentialInformation;
276 using M = typename F::Matrix;
277 using V = typename F::Vector;
278 using P = PropertyTree;
279 F::addCreator("ILU0", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
280 const double w = prm.get<double>("relaxation", 1.0);
281 return std::make_shared<Opm::ParallelOverlappingILU0<M, V, V, C>>(
282 op.getmat(), 0, w, Opm::MILU_VARIANT::ILU);
283 });
284 F::addCreator("ParOverILU0", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
285 const double w = prm.get<double>("relaxation", 1.0);
286 const int n = prm.get<int>("ilulevel", 0);
287 return std::make_shared<Opm::ParallelOverlappingILU0<M, V, V, C>>(
288 op.getmat(), n, w, Opm::MILU_VARIANT::ILU);
289 });
290 F::addCreator("ILUn", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
291 const int n = prm.get<int>("ilulevel", 0);
292 const double w = prm.get<double>("relaxation", 1.0);
293 return std::make_shared<Opm::ParallelOverlappingILU0<M, V, V, C>>(
294 op.getmat(), n, w, Opm::MILU_VARIANT::ILU);
295 });
296 F::addCreator("Jac", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
297 const int n = prm.get<int>("repeats", 1);
298 const double w = prm.get<double>("relaxation", 1.0);
299 return wrapPreconditioner<SeqJac<M, V, V>>(op.getmat(), n, w);
300 });
301 F::addCreator("GS", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
302 const int n = prm.get<int>("repeats", 1);
303 const double w = prm.get<double>("relaxation", 1.0);
304 return wrapPreconditioner<SeqGS<M, V, V>>(op.getmat(), n, w);
305 });
306 F::addCreator("SOR", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
307 const int n = prm.get<int>("repeats", 1);
308 const double w = prm.get<double>("relaxation", 1.0);
309 return wrapPreconditioner<SeqSOR<M, V, V>>(op.getmat(), n, w);
310 });
311 F::addCreator("SSOR", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
312 const int n = prm.get<int>("repeats", 1);
313 const double w = prm.get<double>("relaxation", 1.0);
314 return wrapPreconditioner<SeqSSOR<M, V, V>>(op.getmat(), n, w);
315 });
316
317 // Only add AMG preconditioners to the factory if the operator
318 // is an actual matrix operator.
319 if constexpr (std::is_same_v<O, Dune::MatrixAdapter<M, V, V>>) {
320 F::addCreator("amg", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
321 const std::string smoother = prm.get<std::string>("smoother", "ParOverILU0");
322 if (smoother == "ILU0" || smoother == "ParOverILU0") {
323 using Smoother = SeqILU<M, V, V>;
325 } else if (smoother == "Jac") {
326 using Smoother = SeqJac<M, V, V>;
328 } else if (smoother == "SOR") {
329 using Smoother = SeqSOR<M, V, V>;
331 } else if (smoother == "SSOR") {
332 using Smoother = SeqSSOR<M, V, V>;
334 } else if (smoother == "ILUn") {
335 using Smoother = SeqILU<M, V, V>;
337 } else {
338 OPM_THROW(std::invalid_argument,
339 "Properties: No smoother with name " + smoother + ".");
340 }
341 });
342 F::addCreator("kamg", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
343 const std::string smoother = prm.get<std::string>("smoother", "ParOverILU0");
344 if (smoother == "ILU0" || smoother == "ParOverILU0") {
345 using Smoother = SeqILU<M, V, V>;
347 } else if (smoother == "Jac") {
348 using Smoother = SeqJac<M, V, V>;
350 } else if (smoother == "SOR") {
351 using Smoother = SeqSOR<M, V, V>;
353 // } else if (smoother == "GS") {
354 // using Smoother = SeqGS<M, V, V>;
355 // return makeAmgPreconditioner<Smoother>(op, prm, true);
356 } else if (smoother == "SSOR") {
357 using Smoother = SeqSSOR<M, V, V>;
359 } else if (smoother == "ILUn") {
360 using Smoother = SeqILU<M, V, V>;
362 } else {
363 OPM_THROW(std::invalid_argument,
364 "Properties: No smoother with name " + smoother + ".");
365 }
366 });
367 F::addCreator("famg", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
369 Dune::Amg::Parameters parms;
370 parms.setNoPreSmoothSteps(1);
371 parms.setNoPostSmoothSteps(1);
373 });
374 }
375 if constexpr (std::is_same_v<O, WellModelMatrixAdapter<M, V, V, false>>) {
376 F::addCreator("cprw", [](const O& op, const P& prm, const std::function<V()>& weightsCalculator, std::size_t pressureIndex) {
377 if (pressureIndex == std::numeric_limits<std::size_t>::max()) {
378 OPM_THROW(std::logic_error, "Pressure index out of bounds. It needs to specified for CPR");
379 }
381 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy>>(op, prm, weightsCalculator, pressureIndex);
382 });
383 }
384
385 F::addCreator("cpr", [](const O& op, const P& prm, const std::function<V()>& weightsCalculator, std::size_t pressureIndex) {
386 if (pressureIndex == std::numeric_limits<std::size_t>::max())
387 {
388 OPM_THROW(std::logic_error, "Pressure index out of bounds. It needs to specified for CPR");
389 }
391 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy>>(op, prm, weightsCalculator, pressureIndex);
392 });
393 F::addCreator("cprt", [](const O& op, const P& prm, const std::function<V()>& weightsCalculator, std::size_t pressureIndex) {
394 if (pressureIndex == std::numeric_limits<std::size_t>::max())
395 {
396 OPM_THROW(std::logic_error, "Pressure index out of bounds. It needs to specified for CPR");
397 }
399 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy>>(op, prm, weightsCalculator, pressureIndex);
400 });
401 }
402};
403
404template <class Operator, class Comm>
406{
407}
408
409
410template <class Operator, class Comm>
411PreconditionerFactory<Operator,Comm>&
412PreconditionerFactory<Operator,Comm>::instance()
413{
414 static PreconditionerFactory singleton;
415 return singleton;
416}
417
418template <class Operator, class Comm>
420PreconditionerFactory<Operator,Comm>::
421doCreate(const Operator& op, const PropertyTree& prm,
422 const std::function<Vector()> weightsCalculator,
423 std::size_t pressureIndex)
424{
425 if (!defAdded_) {
426 StandardPreconditioners<Operator,Comm>::add();
427 defAdded_ = true;
428 }
429 const std::string& type = prm.get<std::string>("type", "ParOverILU0");
430 auto it = creators_.find(type);
431 if (it == creators_.end()) {
432 std::ostringstream msg;
433 msg << "Preconditioner type " << type << " is not registered in the factory. Available types are: ";
434 for (const auto& prec : creators_) {
435 msg << prec.first << ' ';
436 }
437 msg << std::endl;
438 OPM_THROW(std::invalid_argument, msg.str());
439 }
440 return it->second(op, prm, weightsCalculator, pressureIndex);
441}
442
443template <class Operator, class Comm>
445PreconditionerFactory<Operator,Comm>::
446doCreate(const Operator& op, const PropertyTree& prm,
447 const std::function<Vector()> weightsCalculator,
448 std::size_t pressureIndex, const Comm& comm)
449{
450 if (!defAdded_) {
451 StandardPreconditioners<Operator,Comm>::add();
452 defAdded_ = true;
453 }
454 const std::string& type = prm.get<std::string>("type", "ParOverILU0");
455 auto it = parallel_creators_.find(type);
456 if (it == parallel_creators_.end()) {
457 std::ostringstream msg;
458 msg << "Parallel preconditioner type " << type
459 << " is not registered in the factory. Available types are: ";
460 for (const auto& prec : parallel_creators_) {
461 msg << prec.first << ' ';
462 }
463 msg << std::endl;
464 OPM_THROW(std::invalid_argument, msg.str());
465 }
466 return it->second(op, prm, weightsCalculator, pressureIndex, comm);
467}
468
469template <class Operator, class Comm>
470void PreconditionerFactory<Operator,Comm>::
471doAddCreator(const std::string& type, Creator c)
472{
473 creators_[type] = c;
474}
475
476template <class Operator, class Comm>
477void PreconditionerFactory<Operator,Comm>::
478doAddCreator(const std::string& type, ParCreator c)
479{
480 parallel_creators_[type] = c;
481}
482
483template <class Operator, class Comm>
486create(const Operator& op, const PropertyTree& prm,
487 const std::function<Vector()>& weightsCalculator,
488 std::size_t pressureIndex)
489{
490 return instance().doCreate(op, prm, weightsCalculator, pressureIndex);
491}
492
493template <class Operator, class Comm>
496create(const Operator& op, const PropertyTree& prm,
497 const std::function<Vector()>& weightsCalculator, const Comm& comm,
498 std::size_t pressureIndex)
499{
500 return instance().doCreate(op, prm, weightsCalculator, pressureIndex, comm);
501}
502
503
504template <class Operator, class Comm>
507create(const Operator& op, const PropertyTree& prm, const Comm& comm,
508 std::size_t pressureIndex)
509{
510 return instance().doCreate(op, prm, std::function<Vector()>(), pressureIndex, comm);
511}
512
513template <class Operator, class Comm>
515addCreator(const std::string& type, Creator creator)
516{
517 instance().doAddCreator(type, creator);
518}
519
520template <class Operator, class Comm>
522addCreator(const std::string& type, ParCreator creator)
523{
524 instance().doAddCreator(type, creator);
525}
526
527using CommSeq = Dune::Amg::SequentialInformation;
528
529template<int Dim>
530using OpFSeq = Dune::MatrixAdapter<Dune::BCRSMatrix<Dune::FieldMatrix<double,Dim,Dim>>,
531 Dune::BlockVector<Dune::FieldVector<double,Dim>>,
532 Dune::BlockVector<Dune::FieldVector<double,Dim>>>;
533template<int Dim>
534using OpBSeq = Dune::MatrixAdapter<Dune::BCRSMatrix<Opm::MatrixBlock<double,Dim,Dim>>,
535 Dune::BlockVector<Dune::FieldVector<double,Dim>>,
536 Dune::BlockVector<Dune::FieldVector<double,Dim>>>;
537
538template<int Dim, bool overlap>
540 Dune::BlockVector<Dune::FieldVector<double,Dim>>,
541 Dune::BlockVector<Dune::FieldVector<double,Dim>>,
542 overlap>;
543
544template<int Dim, bool overlap>
546 Dune::BlockVector<Dune::FieldVector<double,Dim>>,
547 Dune::BlockVector<Dune::FieldVector<double,Dim>>,
548 overlap>;
549
550#if HAVE_MPI
551using CommPar = Dune::OwnerOverlapCopyCommunication<int,int>;
552
553template<int Dim>
554using OpFPar = Dune::OverlappingSchwarzOperator<Dune::BCRSMatrix<Dune::FieldMatrix<double,Dim,Dim>>,
555 Dune::BlockVector<Dune::FieldVector<double,Dim>>,
556 Dune::BlockVector<Dune::FieldVector<double,Dim>>,
557 CommPar>;
558
559template<int Dim>
560using OpBPar = Dune::OverlappingSchwarzOperator<Dune::BCRSMatrix<MatrixBlock<double,Dim,Dim>>,
561 Dune::BlockVector<Dune::FieldVector<double,Dim>>,
562 Dune::BlockVector<Dune::FieldVector<double,Dim>>,
563 CommPar>;
564
565#define INSTANCE_PF_PAR(Dim) \
566 template class PreconditionerFactory<OpBSeq<Dim>,CommPar>; \
567 template class PreconditionerFactory<OpFPar<Dim>,CommPar>; \
568 template class PreconditionerFactory<OpBPar<Dim>,CommPar>; \
569 template class PreconditionerFactory<OpW<Dim,false>,CommPar>; \
570 template class PreconditionerFactory<OpWG<Dim,true>,CommPar>; \
571 template class PreconditionerFactory<OpBPar<Dim>,CommSeq>;
572#endif
573
574#define INSTANCE_PF_SEQ(Dim) \
575 template class PreconditionerFactory<OpFSeq<Dim>,CommSeq>; \
576 template class PreconditionerFactory<OpBSeq<Dim>,CommSeq>; \
577 template class PreconditionerFactory<OpW<Dim,false>,CommSeq>; \
578 template class PreconditionerFactory<OpWG<Dim,true>,CommSeq>;
579
580#if HAVE_MPI
581#define INSTANCE_PF(Dim) \
582 INSTANCE_PF_PAR(Dim) \
583 INSTANCE_PF_SEQ(Dim)
584#else
585#define INSTANCE_PF(Dim) \
586 INSTANCE_PF_SEQ(Dim)
587#endif
588}
The AMG preconditioner.
Parallel algebraic multigrid based on agglomeration.
Definition amgcpr.hh:87
Definition PreconditionerWithUpdate.hpp:40
Definition AquiferInterface.hpp:35
A two-step version of an overlapping Schwarz preconditioner using one step ILU0 as.
Definition ParallelOverlappingILU0.hpp:131
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:486
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:515
std::shared_ptr< Dune::PreconditionerWithUpdate< Vector, Vector > > PrecPtr
The type of pointer returned by create().
Definition PreconditionerFactory.hpp:71
Definition PropertyTree.hpp:37
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition BlackoilPhases.hpp:27
MILU_VARIANT
Definition MILU.hpp:34
@ ILU
Do not perform modified ILU.
Definition PreconditionerFactory.hpp:43
Definition PreconditionerFactory_impl.hpp:47
Definition PreconditionerFactory_impl.hpp:130
static size_t interiorIfGhostLast(const Comm &comm)
Helper method to determine if the local partitioning has the K interior cells from [0,...
Definition PreconditionerFactory_impl.hpp:247