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#include <opm/common/ErrorMacros.hpp>
21#include <opm/common/TimingMacros.hpp>
22
23#include <opm/simulators/linalg/PreconditionerFactory.hpp>
24
25
26#include <opm/simulators/linalg/matrixblock.hh>
28#include <opm/simulators/linalg/FlexibleSolver.hpp>
29#include <opm/simulators/linalg/ilufirstelement.hh>
30#include <opm/simulators/linalg/OwningBlockPreconditioner.hpp>
31#include <opm/simulators/linalg/OwningTwoLevelPreconditioner.hpp>
32#include <opm/simulators/linalg/ParallelOverlappingILU0.hpp>
33#include <opm/simulators/linalg/PressureBhpTransferPolicy.hpp>
34#include <opm/simulators/linalg/PressureTransferPolicy.hpp>
35#include <opm/simulators/linalg/PropertyTree.hpp>
36#include <opm/simulators/linalg/WellOperators.hpp>
37
38#include <dune/istl/owneroverlapcopy.hh>
39#include <dune/istl/preconditioners.hh>
40#include <dune/istl/paamg/amg.hh>
41#include <dune/istl/paamg/kamg.hh>
42#include <dune/istl/paamg/fastamg.hh>
43
44#include <config.h>
45#if HAVE_CUDA
46#include <opm/simulators/linalg/cuistl/CuSeqILU0.hpp>
47#include <opm/simulators/linalg/cuistl/PreconditionerAdapter.hpp>
48#include <opm/simulators/linalg/cuistl/PreconditionerConvertFieldTypeAdapter.hpp>
49#include <opm/simulators/linalg/cuistl/CuBlockPreconditioner.hpp>
50
51#endif
52
53
54namespace Opm {
55
56template<class Smoother>
58{
59 static auto args(const PropertyTree& prm)
60 {
61 using SmootherArgs = typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
62 SmootherArgs smootherArgs;
63 smootherArgs.iterations = prm.get<int>("iterations", 1);
64 // smootherArgs.overlap=SmootherArgs::vertex;
65 // smootherArgs.overlap=SmootherArgs::none;
66 // smootherArgs.overlap=SmootherArgs::aggregate;
67 smootherArgs.relaxationFactor = prm.get<double>("relaxation", 1.0);
68 return smootherArgs;
69 }
70};
71
72template<class M, class V, class C>
74{
75 static auto args(const PropertyTree& prm)
76 {
78 using SmootherArgs = typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
79 SmootherArgs smootherArgs;
80 smootherArgs.iterations = prm.get<int>("iterations", 1);
81 const int iluwitdh = prm.get<int>("iluwidth", 0);
83 const MILU_VARIANT milu = convertString2Milu(prm.get<std::string>("milutype", std::string("ilu")));
84 smootherArgs.setMilu(milu);
85 // smootherArgs.overlap=SmootherArgs::vertex;
86 // smootherArgs.overlap=SmootherArgs::none;
87 // smootherArgs.overlap=SmootherArgs::aggregate;
88 smootherArgs.relaxationFactor = prm.get<double>("relaxation", 1.0);
89 return smootherArgs;
90 }
91};
92
93template <class Operator, class Comm, class Matrix, class Vector>
94typename AMGHelper<Operator, Comm, Matrix, Vector>::Criterion
96{
97 Criterion criterion(15, prm.get<int>("coarsenTarget", 1200));
98 criterion.setDefaultValuesIsotropic(2);
99 criterion.setAlpha(prm.get<double>("alpha", 0.33));
100 criterion.setBeta(prm.get<double>("beta", 1e-5));
101 criterion.setMaxLevel(prm.get<int>("maxlevel", 15));
102 criterion.setSkipIsolated(prm.get<bool>("skip_isolated", false));
103 criterion.setNoPreSmoothSteps(prm.get<int>("pre_smooth", 1));
104 criterion.setNoPostSmoothSteps(prm.get<int>("post_smooth", 1));
105 criterion.setDebugLevel(prm.get<int>("verbosity", 0));
106 // As the default we request to accumulate data to 1 process always as our matrix
107 // graph might be unsymmetric and hence not supported by the PTScotch/ParMetis
108 // calls in DUNE. Accumulating to 1 skips PTScotch/ParMetis
109 criterion.setAccumulate(static_cast<Dune::Amg::AccumulationMode>(prm.get<int>("accumulate", 1)));
110 criterion.setProlongationDampingFactor(prm.get<double>("prolongationdamping", 1.6));
111 criterion.setMaxDistance(prm.get<int>("maxdistance", 2));
112 criterion.setMaxConnectivity(prm.get<int>("maxconnectivity", 15));
113 criterion.setMaxAggregateSize(prm.get<int>("maxaggsize", 6));
114 criterion.setMinAggregateSize(prm.get<int>("minaggsize", 4));
115 return criterion;
116}
117
118template <class Operator, class Comm, class Matrix, class Vector>
119template <class Smoother>
120typename AMGHelper<Operator, Comm, Matrix, Vector>::PrecPtr
121AMGHelper<Operator,Comm,Matrix,Vector>::
122makeAmgPreconditioner(const Operator& op,
123 const PropertyTree& prm,
124 bool useKamg)
125{
126 auto crit = criterion(prm);
127 auto sargs = AMGSmootherArgsHelper<Smoother>::args(prm);
128 if (useKamg) {
130 return std::make_shared<Type>(op, crit, sargs,
131 prm.get<std::size_t>("max_krylov", 1),
132 prm.get<double>("min_reduction", 1e-1));
133 } else {
135 return std::make_shared<Type>(op, crit, sargs);
136 }
137}
138
139template<class Operator, class Comm>
141{
142 static void add()
143 {
144 using namespace Dune;
145 using O = Operator;
146 using C = Comm;
148 using M = typename F::Matrix;
149 using V = typename F::Vector;
150 using P = PropertyTree;
151 F::addCreator("ILU0", [](const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
152 return createParILU(op, prm, comm, 0);
153 });
154 F::addCreator("ParOverILU0", [](const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
155 return createParILU(op, prm, comm, prm.get<int>("ilulevel", 0));
156 });
157 F::addCreator("ILUn", [](const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
158 return createParILU(op, prm, comm, prm.get<int>("ilulevel", 0));
159 });
160 F::addCreator("Jac", [](const O& op, const P& prm, const std::function<V()>&,
161 std::size_t, const C& comm) {
162 const int n = prm.get<int>("repeats", 1);
163 const double w = prm.get<double>("relaxation", 1.0);
165 });
166 F::addCreator("GS", [](const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
167 const int n = prm.get<int>("repeats", 1);
168 const double w = prm.get<double>("relaxation", 1.0);
170 });
171 F::addCreator("SOR", [](const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
172 const int n = prm.get<int>("repeats", 1);
173 const double w = prm.get<double>("relaxation", 1.0);
175 });
176 F::addCreator("SSOR", [](const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
177 const int n = prm.get<int>("repeats", 1);
178 const double w = prm.get<double>("relaxation", 1.0);
180 });
181
182 // Only add AMG preconditioners to the factory if the operator
183 // is the overlapping schwarz operator. This could be extended
184 // later, but at this point no other operators are compatible
185 // with the AMG hierarchy construction.
186 if constexpr (std::is_same_v<O, Dune::OverlappingSchwarzOperator<M, V, V, C>>) {
187 F::addCreator("amg", [](const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
188 const std::string smoother = prm.get<std::string>("smoother", "ParOverILU0");
189 if (smoother == "ILU0" || smoother == "ParOverILU0") {
193 return std::make_shared<Dune::Amg::AMGCPR<O, V, Smoother, C>>(op, crit, sargs, comm);
194 } else {
195 OPM_THROW(std::invalid_argument, "Properties: No smoother with name " + smoother + ".");
196 }
197 });
198 }
199
200 F::addCreator("cpr", [](const O& op, const P& prm, const std::function<V()> weightsCalculator, std::size_t pressureIndex, const C& comm) {
202 if (pressureIndex == std::numeric_limits<std::size_t>::max())
203 {
204 OPM_THROW(std::logic_error, "Pressure index out of bounds. It needs to specified for CPR");
205 }
207 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy, Comm>>(op, prm, weightsCalculator, pressureIndex, comm);
208 });
209 F::addCreator("cprt", [](const O& op, const P& prm, const std::function<V()> weightsCalculator, std::size_t pressureIndex, const C& comm) {
211 if (pressureIndex == std::numeric_limits<std::size_t>::max())
212 {
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 if constexpr (std::is_same_v<O, WellModelGhostLastMatrixAdapter<M, V, V, true>>) {
220 F::addCreator("cprw",
221 [](const O& op, const P& prm, const std::function<V()> weightsCalculator, std::size_t pressureIndex, const C& comm) {
223 if (pressureIndex == std::numeric_limits<std::size_t>::max()) {
224 OPM_THROW(std::logic_error, "Pressure index out of bounds. It needs to specified for CPR");
225 }
227 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy, Comm>>(op, prm, weightsCalculator, pressureIndex, comm);
228 });
229 }
230
231#if HAVE_CUDA
232 F::addCreator("CUILU0", [](const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
233 const double w = prm.get<double>("relaxation", 1.0);
234 using field_type = typename V::field_type;
235 using CuILU0 = typename Opm::cuistl::CuSeqILU0<M, Opm::cuistl::CuVector<field_type>, Opm::cuistl::CuVector<field_type>>;
236 auto cuILU0 = std::make_shared<CuILU0>(op.getmat(), w);
237
238 auto adapted = std::make_shared<Opm::cuistl::PreconditionerAdapter<V, V, CuILU0>>(cuILU0);
239 auto wrapped = std::make_shared<Opm::cuistl::CuBlockPreconditioner<V, V, Comm>>(adapted, comm);
240 return wrapped;
241 });
242#endif
243 }
244
245
247 createParILU(const Operator& op, const PropertyTree& prm, const Comm& comm, const int ilulevel)
248 {
250 using M = typename F::Matrix;
251 using V = typename F::Vector;
252
253 const double w = prm.get<double>("relaxation", 1.0);
254 const bool redblack = prm.get<bool>("redblack", false);
255 const bool reorder_spheres = prm.get<bool>("reorder_spheres", false);
256 // Already a parallel preconditioner. Need to pass comm, but no need to wrap it in a BlockPreconditioner.
257 if (ilulevel == 0) {
258 const std::size_t num_interior = interiorIfGhostLast(comm);
259 return std::make_shared<Opm::ParallelOverlappingILU0<M, V, V, Comm>>(
261 } else {
262 return std::make_shared<Opm::ParallelOverlappingILU0<M, V, V, Comm>>(
264 }
265 }
266
271 static std::size_t interiorIfGhostLast(const Comm& comm)
272 {
273 std::size_t interior_count = 0;
274 std::size_t highest_interior_index = 0;
275 const auto& is = comm.indexSet();
276 for (const auto& ind : is) {
277 if (Comm::OwnerSet::contains(ind.local().attribute())) {
279 highest_interior_index = std::max(highest_interior_index, ind.local().local());
280 }
281 }
283 return interior_count;
284 } else {
285 return is.size();
286 }
287 }
288
289};
290
291template<class Operator>
292struct StandardPreconditioners<Operator,Dune::Amg::SequentialInformation>
293{
294 static void add()
295 {
296 using namespace Dune;
297 using O = Operator;
298 using C = Dune::Amg::SequentialInformation;
300 using M = typename F::Matrix;
301 using V = typename F::Vector;
302 using P = PropertyTree;
303 F::addCreator("ILU0", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
304 const double w = prm.get<double>("relaxation", 1.0);
305 return std::make_shared<Opm::ParallelOverlappingILU0<M, V, V, C>>(
306 op.getmat(), 0, w, Opm::MILU_VARIANT::ILU);
307 });
308 F::addCreator("ParOverILU0", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
309 const double w = prm.get<double>("relaxation", 1.0);
310 const int n = prm.get<int>("ilulevel", 0);
311 return std::make_shared<Opm::ParallelOverlappingILU0<M, V, V, C>>(
312 op.getmat(), n, w, Opm::MILU_VARIANT::ILU);
313 });
314 F::addCreator("ILUn", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
315 const int n = prm.get<int>("ilulevel", 0);
316 const double w = prm.get<double>("relaxation", 1.0);
317 return std::make_shared<Opm::ParallelOverlappingILU0<M, V, V, C>>(
318 op.getmat(), n, w, Opm::MILU_VARIANT::ILU);
319 });
320 F::addCreator("Jac", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
321 const int n = prm.get<int>("repeats", 1);
322 const double w = prm.get<double>("relaxation", 1.0);
323 return wrapPreconditioner<SeqJac<M, V, V>>(op.getmat(), n, w);
324 });
325 F::addCreator("GS", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
326 const int n = prm.get<int>("repeats", 1);
327 const double w = prm.get<double>("relaxation", 1.0);
328 return wrapPreconditioner<SeqGS<M, V, V>>(op.getmat(), n, w);
329 });
330 F::addCreator("SOR", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
331 const int n = prm.get<int>("repeats", 1);
332 const double w = prm.get<double>("relaxation", 1.0);
333 return wrapPreconditioner<SeqSOR<M, V, V>>(op.getmat(), n, w);
334 });
335 F::addCreator("SSOR", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
336 const int n = prm.get<int>("repeats", 1);
337 const double w = prm.get<double>("relaxation", 1.0);
338 return wrapPreconditioner<SeqSSOR<M, V, V>>(op.getmat(), n, w);
339 });
340
341 // Only add AMG preconditioners to the factory if the operator
342 // is an actual matrix operator.
343 if constexpr (std::is_same_v<O, Dune::MatrixAdapter<M, V, V>>) {
344 F::addCreator("amg", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
345 const std::string smoother = prm.get<std::string>("smoother", "ParOverILU0");
346 if (smoother == "ILU0" || smoother == "ParOverILU0") {
347 using Smoother = SeqILU<M, V, V>;
349 } else if (smoother == "Jac") {
350 using Smoother = SeqJac<M, V, V>;
352 } else if (smoother == "SOR") {
353 using Smoother = SeqSOR<M, V, V>;
355 } else if (smoother == "SSOR") {
356 using Smoother = SeqSSOR<M, V, V>;
358 } else if (smoother == "ILUn") {
359 using Smoother = SeqILU<M, V, V>;
361 } else {
362 OPM_THROW(std::invalid_argument,
363 "Properties: No smoother with name " + smoother + ".");
364 }
365 });
366 F::addCreator("kamg", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
367 const std::string smoother = prm.get<std::string>("smoother", "ParOverILU0");
368 if (smoother == "ILU0" || smoother == "ParOverILU0") {
369 using Smoother = SeqILU<M, V, V>;
371 } else if (smoother == "Jac") {
372 using Smoother = SeqJac<M, V, V>;
374 } else if (smoother == "SOR") {
375 using Smoother = SeqSOR<M, V, V>;
377 // } else if (smoother == "GS") {
378 // using Smoother = SeqGS<M, V, V>;
379 // return makeAmgPreconditioner<Smoother>(op, prm, true);
380 } else if (smoother == "SSOR") {
381 using Smoother = SeqSSOR<M, V, V>;
383 } else if (smoother == "ILUn") {
384 using Smoother = SeqILU<M, V, V>;
386 } else {
387 OPM_THROW(std::invalid_argument,
388 "Properties: No smoother with name " + smoother + ".");
389 }
390 });
391 F::addCreator("famg", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
393 Dune::Amg::Parameters parms;
394 parms.setNoPreSmoothSteps(1);
395 parms.setNoPostSmoothSteps(1);
397 });
398 }
399 if constexpr (std::is_same_v<O, WellModelMatrixAdapter<M, V, V, false>>) {
400 F::addCreator("cprw", [](const O& op, const P& prm, const std::function<V()>& weightsCalculator, std::size_t pressureIndex) {
401 if (pressureIndex == std::numeric_limits<std::size_t>::max()) {
402 OPM_THROW(std::logic_error, "Pressure index out of bounds. It needs to specified for CPR");
403 }
405 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy>>(op, prm, weightsCalculator, pressureIndex);
406 });
407 }
408
409 F::addCreator("cpr", [](const O& op, const P& prm, const std::function<V()>& weightsCalculator, std::size_t pressureIndex) {
410 if (pressureIndex == std::numeric_limits<std::size_t>::max())
411 {
412 OPM_THROW(std::logic_error, "Pressure index out of bounds. It needs to specified for CPR");
413 }
415 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy>>(op, prm, weightsCalculator, pressureIndex);
416 });
417 F::addCreator("cprt", [](const O& op, const P& prm, const std::function<V()>& weightsCalculator, std::size_t pressureIndex) {
418 if (pressureIndex == std::numeric_limits<std::size_t>::max())
419 {
420 OPM_THROW(std::logic_error, "Pressure index out of bounds. It needs to specified for CPR");
421 }
423 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy>>(op, prm, weightsCalculator, pressureIndex);
424 });
425
426#if HAVE_CUDA
427 F::addCreator("CUILU0", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
428 const double w = prm.get<double>("relaxation", 1.0);
429 using field_type = typename V::field_type;
430 using CuILU0 = typename Opm::cuistl::CuSeqILU0<M, Opm::cuistl::CuVector<field_type>, Opm::cuistl::CuVector<field_type>>;
431 return std::make_shared<Opm::cuistl::PreconditionerAdapter<V, V, CuILU0>>(std::make_shared<CuILU0>(op.getmat(), w));
432 });
433
434 F::addCreator("CUILU0Float", [](const O& op, const P& prm, const std::function<V()>&, std::size_t) {
435 const double w = prm.get<double>("relaxation", 1.0);
436 using block_type = typename V::block_type;
437 using VTo = Dune::BlockVector<Dune::FieldVector<float, block_type::dimension>>;
438 using matrix_type_to = typename Dune::BCRSMatrix<Dune::FieldMatrix<float, block_type::dimension, block_type::dimension>>;
439 using CuILU0 = typename Opm::cuistl::CuSeqILU0<matrix_type_to, Opm::cuistl::CuVector<float>, Opm::cuistl::CuVector<float>>;
442 auto converted = std::make_shared<Converter>(op.getmat());
443 auto adapted = std::make_shared<Adapter>(std::make_shared<CuILU0>(converted->getConvertedMatrix(), w));
444 converted->setUnderlyingPreconditioner(adapted);
445 return converted;
446
447 });
448#endif
449 }
450};
451
452template <class Operator, class Comm>
454{
455}
456
457
458template <class Operator, class Comm>
459PreconditionerFactory<Operator,Comm>&
460PreconditionerFactory<Operator,Comm>::instance()
461{
462 static PreconditionerFactory singleton;
463 return singleton;
464}
465
466template <class Operator, class Comm>
468PreconditionerFactory<Operator,Comm>::
469doCreate(const Operator& op, const PropertyTree& prm,
470 const std::function<Vector()> weightsCalculator,
471 std::size_t pressureIndex)
472{
473 if (!defAdded_) {
474 StandardPreconditioners<Operator,Comm>::add();
475 defAdded_ = true;
476 }
477 const std::string& type = prm.get<std::string>("type", "ParOverILU0");
478 auto it = creators_.find(type);
479 if (it == creators_.end()) {
480 std::ostringstream msg;
481 msg << "Preconditioner type " << type << " is not registered in the factory. Available types are: ";
482 for (const auto& prec : creators_) {
483 msg << prec.first << ' ';
484 }
485 msg << std::endl;
486 OPM_THROW(std::invalid_argument, msg.str());
487 }
488 return it->second(op, prm, weightsCalculator, pressureIndex);
489}
490
491template <class Operator, class Comm>
493PreconditionerFactory<Operator,Comm>::
494doCreate(const Operator& op, const PropertyTree& prm,
495 const std::function<Vector()> weightsCalculator,
496 std::size_t pressureIndex, const Comm& comm)
497{
498 if (!defAdded_) {
499 StandardPreconditioners<Operator,Comm>::add();
500 defAdded_ = true;
501 }
502 const std::string& type = prm.get<std::string>("type", "ParOverILU0");
503 auto it = parallel_creators_.find(type);
504 if (it == parallel_creators_.end()) {
505 std::ostringstream msg;
506 msg << "Parallel preconditioner type " << type
507 << " is not registered in the factory. Available types are: ";
508 for (const auto& prec : parallel_creators_) {
509 msg << prec.first << ' ';
510 }
511 msg << std::endl;
512 OPM_THROW(std::invalid_argument, msg.str());
513 }
514 return it->second(op, prm, weightsCalculator, pressureIndex, comm);
515}
516
517template <class Operator, class Comm>
518void PreconditionerFactory<Operator,Comm>::
519doAddCreator(const std::string& type, Creator c)
520{
521 creators_[type] = c;
522}
523
524template <class Operator, class Comm>
525void PreconditionerFactory<Operator,Comm>::
526doAddCreator(const std::string& type, ParCreator c)
527{
528 parallel_creators_[type] = c;
529}
530
531template <class Operator, class Comm>
534create(const Operator& op, const PropertyTree& prm,
535 const std::function<Vector()>& weightsCalculator,
536 std::size_t pressureIndex)
537{
538 return instance().doCreate(op, prm, weightsCalculator, pressureIndex);
539}
540
541template <class Operator, class Comm>
544create(const Operator& op, const PropertyTree& prm,
545 const std::function<Vector()>& weightsCalculator, const Comm& comm,
546 std::size_t pressureIndex)
547{
548 return instance().doCreate(op, prm, weightsCalculator, pressureIndex, comm);
549}
550
551
552template <class Operator, class Comm>
555create(const Operator& op, const PropertyTree& prm, const Comm& comm,
556 std::size_t pressureIndex)
557{
558 return instance().doCreate(op, prm, std::function<Vector()>(), pressureIndex, comm);
559}
560
561template <class Operator, class Comm>
563addCreator(const std::string& type, Creator creator)
564{
565 instance().doAddCreator(type, creator);
566}
567
568template <class Operator, class Comm>
570addCreator(const std::string& type, ParCreator creator)
571{
572 instance().doAddCreator(type, creator);
573}
574
575using CommSeq = Dune::Amg::SequentialInformation;
576
577template<int Dim>
578using OpFSeq = Dune::MatrixAdapter<Dune::BCRSMatrix<Dune::FieldMatrix<double,Dim,Dim>>,
579 Dune::BlockVector<Dune::FieldVector<double,Dim>>,
580 Dune::BlockVector<Dune::FieldVector<double,Dim>>>;
581template<int Dim>
582using OpBSeq = Dune::MatrixAdapter<Dune::BCRSMatrix<Opm::MatrixBlock<double,Dim,Dim>>,
583 Dune::BlockVector<Dune::FieldVector<double,Dim>>,
584 Dune::BlockVector<Dune::FieldVector<double,Dim>>>;
585
586template<int Dim, bool overlap>
588 Dune::BlockVector<Dune::FieldVector<double,Dim>>,
589 Dune::BlockVector<Dune::FieldVector<double,Dim>>,
590 overlap>;
591
592template<int Dim, bool overlap>
594 Dune::BlockVector<Dune::FieldVector<double,Dim>>,
595 Dune::BlockVector<Dune::FieldVector<double,Dim>>,
596 overlap>;
597
598#if HAVE_MPI
599using CommPar = Dune::OwnerOverlapCopyCommunication<int,int>;
600
601template<int Dim>
602using OpFPar = Dune::OverlappingSchwarzOperator<Dune::BCRSMatrix<Dune::FieldMatrix<double,Dim,Dim>>,
603 Dune::BlockVector<Dune::FieldVector<double,Dim>>,
604 Dune::BlockVector<Dune::FieldVector<double,Dim>>,
605 CommPar>;
606
607template<int Dim>
608using OpBPar = Dune::OverlappingSchwarzOperator<Dune::BCRSMatrix<MatrixBlock<double,Dim,Dim>>,
609 Dune::BlockVector<Dune::FieldVector<double,Dim>>,
610 Dune::BlockVector<Dune::FieldVector<double,Dim>>,
611 CommPar>;
612
613#define INSTANCE_PF_PAR(Dim) \
614 template class PreconditionerFactory<OpBSeq<Dim>,CommPar>; \
615 template class PreconditionerFactory<OpFPar<Dim>,CommPar>; \
616 template class PreconditionerFactory<OpBPar<Dim>,CommPar>; \
617 template class PreconditionerFactory<OpW<Dim,false>,CommPar>; \
618 template class PreconditionerFactory<OpWG<Dim,true>,CommPar>; \
619 template class PreconditionerFactory<OpBPar<Dim>,CommSeq>;
620#endif
621
622#define INSTANCE_PF_SEQ(Dim) \
623 template class PreconditionerFactory<OpFSeq<Dim>,CommSeq>; \
624 template class PreconditionerFactory<OpBSeq<Dim>,CommSeq>; \
625 template class PreconditionerFactory<OpW<Dim,false>,CommSeq>; \
626 template class PreconditionerFactory<OpWG<Dim,true>,CommSeq>;
627
628#if HAVE_MPI
629#define INSTANCE_PF(Dim) \
630 INSTANCE_PF_PAR(Dim) \
631 INSTANCE_PF_SEQ(Dim)
632#else
633#define INSTANCE_PF(Dim) \
634 INSTANCE_PF_SEQ(Dim)
635#endif
636}
The AMG preconditioner.
Parallel algebraic multigrid based on agglomeration.
Definition amgcpr.hh:88
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:534
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:563
std::shared_ptr< Dune::PreconditionerWithUpdate< Vector, Vector > > PrecPtr
The type of pointer returned by create().
Definition PreconditionerFactory.hpp:71
Definition PropertyTree.hpp:37
Sequential ILU0 preconditioner on the GPU through the CuSparse library.
Definition CuSeqILU0.hpp:50
Makes a CUDA preconditioner available to a CPU simulator.
Definition PreconditionerAdapter.hpp:43
Converts the field type (eg.
Definition PreconditionerConvertFieldTypeAdapter.hpp:86
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:58
Definition PreconditionerFactory_impl.hpp:141
static std::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:271