Go to the documentation of this file.
20#include <opm/common/ErrorMacros.hpp>
21#include <opm/common/TimingMacros.hpp>
37#include <opm/simulators/linalg/ilufirstelement.hh>
38#include <opm/simulators/linalg/matrixblock.hh>
40#include <dune/common/unused.hh>
41#include <dune/istl/owneroverlapcopy.hh>
42#include <dune/istl/paamg/amg.hh>
43#include <dune/istl/paamg/fastamg.hh>
44#include <dune/istl/paamg/kamg.hh>
45#include <dune/istl/preconditioners.hh>
56template < class Smoother>
60 using SmootherArgs = typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
61 SmootherArgs smootherArgs;
62 smootherArgs.iterations = prm. get< int>( "iterations", 1);
66 smootherArgs.relaxationFactor = prm. get< double>( "relaxation", 1.0);
71template < class M, class V, class C>
76 using SmootherArgs = typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
77 SmootherArgs smootherArgs;
78 smootherArgs.iterations = prm. get< int>( "iterations", 1);
79 const int iluwitdh = prm. get< int>( "iluwidth", 0);
80 smootherArgs.setN(iluwitdh);
82 smootherArgs.setMilu(milu);
86 smootherArgs.relaxationFactor = prm. get< double>( "relaxation", 1.0);
94auto setUseFixedOrder(C criterion, bool booleanValue) -> decltype(criterion.setUseFixedOrder(booleanValue))
96 return criterion.setUseFixedOrder(booleanValue);
104template < class Operator, class Comm, class Matrix, class Vector>
108 Criterion criterion(15, prm. get< int>( "coarsenTarget", 1200));
109 criterion.setDefaultValuesIsotropic(2);
110 criterion.setAlpha(prm. get< double>( "alpha", 0.33));
111 criterion.setBeta(prm. get< double>( "beta", 1e-5));
112 criterion.setMaxLevel(prm. get< int>( "maxlevel", 15));
113 criterion.setSkipIsolated(prm. get< bool>( "skip_isolated", false));
114 criterion.setNoPreSmoothSteps(prm. get< int>( "pre_smooth", 1));
115 criterion.setNoPostSmoothSteps(prm. get< int>( "post_smooth", 1));
116 criterion.setDebugLevel(prm. get< int>( "verbosity", 0));
120 criterion.setAccumulate( static_cast<Dune::Amg::AccumulationMode >(prm. get< int>( "accumulate", 1)));
121 criterion.setProlongationDampingFactor(prm. get< double>( "prolongationdamping", 1.6));
122 criterion.setMaxDistance(prm. get< int>( "maxdistance", 2));
123 criterion.setMaxConnectivity(prm. get< int>( "maxconnectivity", 15));
124 criterion.setMaxAggregateSize(prm. get< int>( "maxaggsize", 6));
125 criterion.setMinAggregateSize(prm. get< int>( "minaggsize", 4));
130template < class Operator, class Comm, class Matrix, class Vector>
131template < class Smoother>
137 auto crit = criterion(prm);
141 return std::make_shared<Type>(
142 op, crit, sargs, prm. get<std::size_t>( "max_krylov", 1), prm. get< double>( "min_reduction", 1e-1));
145 return std::make_shared<Type>(op, crit, sargs);
149template < class Operator, class Comm>
153 using namespace Dune;
157 using M = typename F::Matrix;
158 using V = typename F::Vector;
160 F::addCreator( "ILU0", []( const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
163 F::addCreator( "ParOverILU0",
164 []( const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
165 return createParILU(op, prm, comm, prm.get< int>( "ilulevel", 0));
167 F::addCreator( "ILUn", []( const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
168 return createParILU(op, prm, comm, prm.get< int>( "ilulevel", 0));
170 F::addCreator( "DuneILU", []( const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
171 const int n = prm.get< int>( "ilulevel", 0);
172 const double w = prm.get< double>( "relaxation", 1.0);
173 const bool resort = prm.get< bool>( "resort", false);
174 return wrapBlockPreconditioner<RebuildOnUpdatePreconditioner<Dune::SeqILU<M, V, V>>>(
175 comm, op.getmat(), n, w, resort);
177 F::addCreator( "DILU", []( const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
178 DUNE_UNUSED_PARAMETER(prm);
179 return wrapBlockPreconditioner<MultithreadDILU<M, V, V>>(comm, op.getmat());
181 F::addCreator( "Jac", []( const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
182 const int n = prm.get< int>( "repeats", 1);
183 const double w = prm.get< double>( "relaxation", 1.0);
184 return wrapBlockPreconditioner<DummyUpdatePreconditioner<SeqJac<M, V, V>>>(comm, op.getmat(), n, w);
186 F::addCreator( "GS", []( const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
187 const int n = prm.get< int>( "repeats", 1);
188 const double w = prm.get< double>( "relaxation", 1.0);
189 return wrapBlockPreconditioner<DummyUpdatePreconditioner<SeqGS<M, V, V>>>(comm, op.getmat(), n, w);
191 F::addCreator( "SOR", []( const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
192 const int n = prm.get< int>( "repeats", 1);
193 const double w = prm.get< double>( "relaxation", 1.0);
194 return wrapBlockPreconditioner<DummyUpdatePreconditioner<SeqSOR<M, V, V>>>(comm, op.getmat(), n, w);
196 F::addCreator( "SSOR", []( const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
197 const int n = prm.get< int>( "repeats", 1);
198 const double w = prm.get< double>( "relaxation", 1.0);
199 return wrapBlockPreconditioner<DummyUpdatePreconditioner<SeqSSOR<M, V, V>>>(comm, op.getmat(), n, w);
206 if constexpr (std::is_same_v<O, Dune::OverlappingSchwarzOperator<M, V, V, C>>) {
207 F::addCreator( "amg", []( const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
208 using PrecPtr = std::shared_ptr<Dune::PreconditionerWithUpdate<V, V>>;
209 const std::string smoother = prm.get<std::string>( "smoother", "ParOverILU0");
211 if (smoother == "ILU0" || smoother == "ParOverILU0") {
215 PrecPtr prec = std::make_shared<Dune::Amg::AMGCPR<O, V, Smoother, C>>(op, crit, sargs, comm);
217 } else if (smoother == "DILU") {
219 using Smoother = Dune::BlockPreconditioner<V, V, C, SeqSmoother>;
220 using SmootherArgs = typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
223 PrecPtr prec = std::make_shared<Dune::Amg::AMGCPR<O, V, Smoother, C>>(op, crit, sargs, comm);
225 } else if (smoother == "Jac") {
226 using SeqSmoother = SeqJac<M, V, V>;
227 using Smoother = Dune::BlockPreconditioner<V, V, C, SeqSmoother>;
228 using SmootherArgs = typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
231 PrecPtr prec = std::make_shared<Dune::Amg::AMGCPR<O, V, Smoother, C>>(op, crit, sargs, comm);
233 } else if (smoother == "GS") {
234 using SeqSmoother = SeqGS<M, V, V>;
235 using Smoother = Dune::BlockPreconditioner<V, V, C, SeqSmoother>;
236 using SmootherArgs = typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
239 PrecPtr prec = std::make_shared<Dune::Amg::AMGCPR<O, V, Smoother, C>>(op, crit, sargs, comm);
241 } else if (smoother == "SOR") {
242 using SeqSmoother = SeqSOR<M, V, V>;
243 using Smoother = Dune::BlockPreconditioner<V, V, C, SeqSmoother>;
244 using SmootherArgs = typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
247 PrecPtr prec = std::make_shared<Dune::Amg::AMGCPR<O, V, Smoother, C>>(op, crit, sargs, comm);
249 } else if (smoother == "SSOR") {
250 using SeqSmoother = SeqSSOR<M, V, V>;
251 using Smoother = Dune::BlockPreconditioner<V, V, C, SeqSmoother>;
252 using SmootherArgs = typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
255 PrecPtr prec = std::make_shared<Dune::Amg::AMGCPR<O, V, Smoother, C>>(op, crit, sargs, comm);
257 } else if (smoother == "ILUn") {
258 using SeqSmoother = SeqILU<M, V, V>;
259 using Smoother = Dune::BlockPreconditioner<V, V, C, SeqSmoother>;
260 using SmootherArgs = typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
263 PrecPtr prec = std::make_shared<Dune::Amg::AMGCPR<O, V, Smoother, C>>(op, crit, sargs, comm);
266 OPM_THROW(std::invalid_argument, "Properties: No smoother with name " + smoother + ".");
274 const std::function<V()> weightsCalculator,
275 std::size_t pressureIndex,
277 assert(weightsCalculator);
278 if (pressureIndex == std::numeric_limits<std::size_t>::max()) {
279 OPM_THROW(std::logic_error,
280 "Pressure index out of bounds. It needs to specified for CPR");
283 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy, Comm>>(
284 op, prm, weightsCalculator, pressureIndex, comm);
286 F::addCreator( "cprt",
289 const std::function<V()> weightsCalculator,
290 std::size_t pressureIndex,
292 assert(weightsCalculator);
293 if (pressureIndex == std::numeric_limits<std::size_t>::max()) {
294 OPM_THROW(std::logic_error,
295 "Pressure index out of bounds. It needs to specified for CPR");
298 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy, Comm>>(
299 op, prm, weightsCalculator, pressureIndex, comm);
302 if constexpr (std::is_same_v<O, WellModelGhostLastMatrixAdapter<M, V, V, true>>) {
303 F::addCreator( "cprw",
306 const std::function<V()> weightsCalculator,
307 std::size_t pressureIndex,
309 assert(weightsCalculator);
310 if (pressureIndex == std::numeric_limits<std::size_t>::max()) {
311 OPM_THROW(std::logic_error,
312 "Pressure index out of bounds. It needs to specified for CPR");
315 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy, Comm>>(
316 op, prm, weightsCalculator, pressureIndex, comm);
321 F::addCreator( "CUILU0", []( const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
322 const double w = prm.get< double>( "relaxation", 1.0);
323 using field_type = typename V::field_type;
324 using CuILU0 = typename Opm::cuistl::
325 CuSeqILU0<M, Opm::cuistl::CuVector<field_type>, Opm::cuistl::CuVector<field_type>>;
326 auto cuILU0 = std::make_shared<CuILU0>(op.getmat(), w);
328 auto adapted = std::make_shared<Opm::cuistl::PreconditionerAdapter<V, V, CuILU0>>(cuILU0);
329 auto wrapped = std::make_shared<Opm::cuistl::CuBlockPreconditioner<V, V, Comm>>(adapted, comm);
333 F::addCreator( "CUJac", []( const O& op, const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
334 const double w = prm.get< double>( "relaxation", 1.0);
335 using field_type = typename V::field_type;
338 auto cuJac = std::make_shared<CuJac>(op.getmat(), w);
340 auto adapted = std::make_shared<Opm::cuistl::PreconditionerAdapter<V, V, CuJac>>(cuJac);
341 auto wrapped = std::make_shared<Opm::cuistl::CuBlockPreconditioner<V, V, Comm>>(adapted, comm);
345 F::addCreator( "CUDILU", []( const O& op, [[maybe_unused]] const P& prm, const std::function<V()>&, std::size_t, const C& comm) {
346 using field_type = typename V::field_type;
348 auto cuDILU = std::make_shared<CuDILU>(op.getmat());
350 auto adapted = std::make_shared<Opm::cuistl::PreconditionerAdapter<V, V, CuDILU>>(cuDILU);
351 auto wrapped = std::make_shared<Opm::cuistl::CuBlockPreconditioner<V, V, Comm>>(adapted, comm);
362 using M = typename F::Matrix;
363 using V = typename F::Vector;
365 const double w = prm. get< double>( "relaxation", 1.0);
366 const bool redblack = prm. get< bool>( "redblack", false);
367 const bool reorder_spheres = prm. get< bool>( "reorder_spheres", false);
371 return std::make_shared<Opm::ParallelOverlappingILU0<M, V, V, Comm>>(
374 return std::make_shared<Opm::ParallelOverlappingILU0<M, V, V, Comm>>(
385 std::size_t interior_count = 0;
386 std::size_t highest_interior_index = 0;
387 const auto& is = comm.indexSet();
388 for ( const auto& ind : is) {
389 if (Comm::OwnerSet::contains(ind.local().attribute())) {
391 highest_interior_index = std::max(highest_interior_index, ind.local().local());
394 if (highest_interior_index + 1 == interior_count) {
395 return interior_count;
402template < class Operator>
406 using namespace Dune;
408 using C = Dune::Amg::SequentialInformation;
410 using M = typename F::Matrix;
411 using V = typename F::Vector;
413 F::addCreator( "ILU0", []( const O& op, const P& prm, const std::function<V()>&, std::size_t) {
414 const double w = prm.get< double>( "relaxation", 1.0);
415 return std::make_shared<Opm::ParallelOverlappingILU0<M, V, V, C>>(
418 F::addCreator( "DuneILU", []( const O& op, const P& prm, const std::function<V()>&, std::size_t) {
419 const double w = prm.get< double>( "relaxation", 1.0);
420 const int n = prm.get< int>( "ilulevel", 0);
421 const bool resort = prm.get< bool>( "resort", false);
422 return getRebuildOnUpdateWrapper<Dune::SeqILU<M, V, V>>(op.getmat(), n, w, resort);
424 F::addCreator( "ParOverILU0", []( const O& op, const P& prm, const std::function<V()>&, std::size_t) {
425 const double w = prm.get< double>( "relaxation", 1.0);
426 const int n = prm.get< int>( "ilulevel", 0);
427 return std::make_shared<Opm::ParallelOverlappingILU0<M, V, V, C>>(
430 F::addCreator( "ILUn", []( const O& op, const P& prm, const std::function<V()>&, std::size_t) {
431 const int n = prm.get< int>( "ilulevel", 0);
432 const double w = prm.get< double>( "relaxation", 1.0);
433 return std::make_shared<Opm::ParallelOverlappingILU0<M, V, V, C>>(
436 F::addCreator( "DILU", []( const O& op, const P& prm, const std::function<V()>&, std::size_t) {
437 DUNE_UNUSED_PARAMETER(prm);
438 return std::make_shared<MultithreadDILU<M, V, V>>(op.getmat());
440 F::addCreator( "Jac", []( const O& op, const P& prm, const std::function<V()>&, std::size_t) {
441 const int n = prm.get< int>( "repeats", 1);
442 const double w = prm.get< double>( "relaxation", 1.0);
443 return getDummyUpdateWrapper<SeqJac<M, V, V>>(op.getmat(), n, w);
445 F::addCreator( "GS", []( const O& op, const P& prm, const std::function<V()>&, std::size_t) {
446 const int n = prm.get< int>( "repeats", 1);
447 const double w = prm.get< double>( "relaxation", 1.0);
448 return getDummyUpdateWrapper<SeqGS<M, V, V>>(op.getmat(), n, w);
450 F::addCreator( "SOR", []( const O& op, const P& prm, const std::function<V()>&, std::size_t) {
451 const int n = prm.get< int>( "repeats", 1);
452 const double w = prm.get< double>( "relaxation", 1.0);
453 return getDummyUpdateWrapper<SeqSOR<M, V, V>>(op.getmat(), n, w);
455 F::addCreator( "SSOR", []( const O& op, const P& prm, const std::function<V()>&, std::size_t) {
456 const int n = prm.get< int>( "repeats", 1);
457 const double w = prm.get< double>( "relaxation", 1.0);
458 return getDummyUpdateWrapper<SeqSSOR<M, V, V>>(op.getmat(), n, w);
463 if constexpr (std::is_same_v<O, Dune::MatrixAdapter<M, V, V>>) {
464 F::addCreator( "amg", []( const O& op, const P& prm, const std::function<V()>&, std::size_t) {
465 const std::string smoother = prm.get<std::string>( "smoother", "ParOverILU0");
466 if (smoother == "ILU0" || smoother == "ParOverILU0") {
467 using Smoother = SeqILU<M, V, V>;
469 } else if (smoother == "Jac") {
470 using Smoother = SeqJac<M, V, V>;
472 } else if (smoother == "GS") {
473 using Smoother = SeqGS<M, V, V>;
475 } else if (smoother == "DILU") {
478 } else if (smoother == "SOR") {
479 using Smoother = SeqSOR<M, V, V>;
481 } else if (smoother == "SSOR") {
482 using Smoother = SeqSSOR<M, V, V>;
484 } else if (smoother == "ILUn") {
485 using Smoother = SeqILU<M, V, V>;
488 OPM_THROW(std::invalid_argument, "Properties: No smoother with name " + smoother + ".");
491 F::addCreator( "kamg", []( const O& op, const P& prm, const std::function<V()>&, std::size_t) {
492 const std::string smoother = prm.get<std::string>( "smoother", "ParOverILU0");
493 if (smoother == "ILU0" || smoother == "ParOverILU0") {
494 using Smoother = SeqILU<M, V, V>;
496 } else if (smoother == "Jac") {
497 using Smoother = SeqJac<M, V, V>;
499 } else if (smoother == "SOR") {
500 using Smoother = SeqSOR<M, V, V>;
502 } else if (smoother == "GS") {
503 using Smoother = SeqGS<M, V, V>;
505 } else if (smoother == "SSOR") {
506 using Smoother = SeqSSOR<M, V, V>;
508 } else if (smoother == "ILUn") {
509 using Smoother = SeqILU<M, V, V>;
512 OPM_THROW(std::invalid_argument, "Properties: No smoother with name " + smoother + ".");
515 F::addCreator( "famg", []( const O& op, const P& prm, const std::function<V()>&, std::size_t) {
517 Dune::Amg::Parameters parms;
518 parms.setNoPreSmoothSteps(1);
519 parms.setNoPostSmoothSteps(1);
520 return getRebuildOnUpdateWrapper<Dune::Amg::FastAMG<O, V>>(op, crit, parms);
523 if constexpr (std::is_same_v<O, WellModelMatrixAdapter<M, V, V, false>>) {
526 []( const O& op, const P& prm, const std::function<V()>& weightsCalculator, std::size_t pressureIndex) {
527 if (pressureIndex == std::numeric_limits<std::size_t>::max()) {
528 OPM_THROW(std::logic_error, "Pressure index out of bounds. It needs to specified for CPR");
530 using LevelTransferPolicy
532 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy>>(
533 op, prm, weightsCalculator, pressureIndex);
539 []( const O& op, const P& prm, const std::function<V()>& weightsCalculator, std::size_t pressureIndex) {
540 if (pressureIndex == std::numeric_limits<std::size_t>::max()) {
541 OPM_THROW(std::logic_error, "Pressure index out of bounds. It needs to specified for CPR");
544 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy>>(
545 op, prm, weightsCalculator, pressureIndex);
549 []( const O& op, const P& prm, const std::function<V()>& weightsCalculator, std::size_t pressureIndex) {
550 if (pressureIndex == std::numeric_limits<std::size_t>::max()) {
551 OPM_THROW(std::logic_error, "Pressure index out of bounds. It needs to specified for CPR");
554 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy>>(
555 op, prm, weightsCalculator, pressureIndex);
559 F::addCreator( "CUILU0", []( const O& op, const P& prm, const std::function<V()>&, std::size_t) {
560 const double w = prm.get< double>( "relaxation", 1.0);
561 using field_type = typename V::field_type;
562 using CuILU0 = typename Opm::cuistl::
563 CuSeqILU0<M, Opm::cuistl::CuVector<field_type>, Opm::cuistl::CuVector<field_type>>;
564 return std::make_shared<Opm::cuistl::PreconditionerAdapter<V, V, CuILU0>>(
565 std::make_shared<CuILU0>(op.getmat(), w));
568 F::addCreator( "CUILU0Float", []( const O& op, const P& prm, const std::function<V()>&, std::size_t) {
569 const double w = prm.get< double>( "relaxation", 1.0);
570 using block_type = typename V::block_type;
571 using VTo = Dune::BlockVector<Dune::FieldVector<float, block_type::dimension>>;
572 using matrix_type_to =
573 typename Dune::BCRSMatrix<Dune::FieldMatrix<float, block_type::dimension, block_type::dimension>>;
574 using CuILU0 = typename Opm::cuistl::
575 CuSeqILU0<matrix_type_to, Opm::cuistl::CuVector<float>, Opm::cuistl::CuVector<float>>;
578 auto converted = std::make_shared<Converter>(op.getmat());
579 auto adapted = std::make_shared<Adapter>(std::make_shared<CuILU0>(converted->getConvertedMatrix(), w));
580 converted->setUnderlyingPreconditioner(adapted);
584 F::addCreator( "CUJac", []( const O& op, const P& prm, const std::function<V()>&, std::size_t) {
585 const double w = prm.get< double>( "relaxation", 1.0);
586 using field_type = typename V::field_type;
589 return std::make_shared<Opm::cuistl::PreconditionerAdapter<V, V, CUJac>>(
590 std::make_shared<CUJac>(op.getmat(), w));
593 F::addCreator( "CUDILU", []( const O& op, [[maybe_unused]] const P& prm, const std::function<V()>&, std::size_t) {
594 using field_type = typename V::field_type;
596 return std::make_shared<Opm::cuistl::PreconditionerAdapter<V, V, CUDILU>>(std::make_shared<CUDILU>(op.getmat()));
599 F::addCreator( "CUDILUFloat", []( const O& op, [[maybe_unused]] const P& prm, const std::function<V()>&, std::size_t) {
600 using block_type = typename V::block_type;
601 using VTo = Dune::BlockVector<Dune::FieldVector<float, block_type::dimension>>;
602 using matrix_type_to = typename Dune::BCRSMatrix<Dune::FieldMatrix<float, block_type::dimension, block_type::dimension>>;
606 auto converted = std::make_shared<Converter>(op.getmat());
607 auto adapted = std::make_shared<Adapter>(std::make_shared<CuDILU>(converted->getConvertedMatrix()));
608 converted->setUnderlyingPreconditioner(adapted);
615template < class Operator, class Comm>
616PreconditionerFactory<Operator, Comm>::PreconditionerFactory()
621template < class Operator, class Comm>
622PreconditionerFactory<Operator, Comm>&
623PreconditionerFactory<Operator, Comm>::instance()
625 static PreconditionerFactory singleton;
629template < class Operator, class Comm>
631PreconditionerFactory<Operator, Comm>::doCreate( const Operator& op,
632 const PropertyTree& prm,
633 const std::function<Vector()> weightsCalculator,
634 std::size_t pressureIndex)
640 const std::string& type = prm.get<std::string>( "type", "ParOverILU0");
641 auto it = creators_.find(type);
642 if (it == creators_.end()) {
643 std::ostringstream msg;
644 msg << "Preconditioner type " << type << " is not registered in the factory. Available types are: ";
645 for ( const auto& prec : creators_) {
646 msg << prec.first << ' ';
649 OPM_THROW(std::invalid_argument, msg.str());
651 return it->second(op, prm, weightsCalculator, pressureIndex);
654template < class Operator, class Comm>
656PreconditionerFactory<Operator, Comm>::doCreate( const Operator& op,
657 const PropertyTree& prm,
658 const std::function<Vector()> weightsCalculator,
659 std::size_t pressureIndex,
666 const std::string& type = prm.get<std::string>( "type", "ParOverILU0");
667 auto it = parallel_creators_.find(type);
668 if (it == parallel_creators_.end()) {
669 std::ostringstream msg;
670 msg << "Parallel preconditioner type " << type << " is not registered in the factory. Available types are: ";
671 for ( const auto& prec : parallel_creators_) {
672 msg << prec.first << ' ';
675 OPM_THROW(std::invalid_argument, msg.str());
677 return it->second(op, prm, weightsCalculator, pressureIndex, comm);
680template < class Operator, class Comm>
682PreconditionerFactory<Operator, Comm>::doAddCreator( const std::string& type, Creator c)
687template < class Operator, class Comm>
689PreconditionerFactory<Operator, Comm>::doAddCreator( const std::string& type, ParCreator c)
691 parallel_creators_[type] = c;
694template < class Operator, class Comm>
698 const std::function< Vector()>& weightsCalculator,
699 std::size_t pressureIndex)
701 return instance().doCreate(op, prm, weightsCalculator, pressureIndex);
704template < class Operator, class Comm>
708 const std::function< Vector()>& weightsCalculator,
710 std::size_t pressureIndex)
712 return instance().doCreate(op, prm, weightsCalculator, pressureIndex, comm);
716template < class Operator, class Comm>
721 std::size_t pressureIndex)
723 return instance().doCreate(op, prm, std::function< Vector()>(), pressureIndex, comm);
726template < class Operator, class Comm>
730 instance().doAddCreator(type, creator);
733template < class Operator, class Comm>
737 instance().doAddCreator(type, creator);
740using CommSeq = Dune::Amg::SequentialInformation;
743using OpFSeq = Dune::MatrixAdapter<Dune::BCRSMatrix<Dune::FieldMatrix<double, Dim, Dim>>,
744 Dune::BlockVector<Dune::FieldVector<double, Dim>>,
745 Dune::BlockVector<Dune::FieldVector<double, Dim>>>;
747using OpBSeq = Dune::MatrixAdapter<Dune::BCRSMatrix<Opm::MatrixBlock<double, Dim, Dim>>,
748 Dune::BlockVector<Dune::FieldVector<double, Dim>>,
749 Dune::BlockVector<Dune::FieldVector<double, Dim>>>;
751template < int Dim, bool overlap>
753 Dune::BlockVector<Dune::FieldVector<double, Dim>>,
754 Dune::BlockVector<Dune::FieldVector<double, Dim>>,
757template < int Dim, bool overlap>
759 Dune::BlockVector<Dune::FieldVector<double, Dim>>,
760 Dune::BlockVector<Dune::FieldVector<double, Dim>>,
764using CommPar = Dune::OwnerOverlapCopyCommunication<int, int>;
767using OpFPar = Dune::OverlappingSchwarzOperator<Dune::BCRSMatrix<Dune::FieldMatrix<double, Dim, Dim>>,
768 Dune::BlockVector<Dune::FieldVector<double, Dim>>,
769 Dune::BlockVector<Dune::FieldVector<double, Dim>>,
773using OpBPar = Dune::OverlappingSchwarzOperator<Dune::BCRSMatrix<MatrixBlock<double, Dim, Dim>>,
774 Dune::BlockVector<Dune::FieldVector<double, Dim>>,
775 Dune::BlockVector<Dune::FieldVector<double, Dim>>,
778#define INSTANCE_PF_PAR(Dim) \
779 template class PreconditionerFactory<OpBSeq<Dim>, CommPar>; \
780 template class PreconditionerFactory<OpFPar<Dim>, CommPar>; \
781 template class PreconditionerFactory<OpBPar<Dim>, CommPar>; \
782 template class PreconditionerFactory<OpW<Dim, false>, CommPar>; \
783 template class PreconditionerFactory<OpWG<Dim, true>, CommPar>; \
784 template class PreconditionerFactory<OpBPar<Dim>, CommSeq>;
787#define INSTANCE_PF_SEQ(Dim) \
788 template class PreconditionerFactory<OpFSeq<Dim>, CommSeq>; \
789 template class PreconditionerFactory<OpBSeq<Dim>, CommSeq>; \
790 template class PreconditionerFactory<OpW<Dim, false>, CommSeq>; \
791 template class PreconditionerFactory<OpWG<Dim, true>, CommSeq>;
794#define INSTANCE_PF(Dim) \
795 INSTANCE_PF_PAR(Dim) \
798#define INSTANCE_PF(Dim) INSTANCE_PF_SEQ(Dim)
Dune::OwnerOverlapCopyCommunication< int, int > Comm Definition: FlexibleSolver_impl.hpp:274
Parallel algebraic multigrid based on agglomeration. Definition: amgcpr.hh:88
Definition: PreconditionerWithUpdate.hpp:40
The OpenMP thread parallelized DILU preconditioner. Definition: DILU.hpp:53
A two-step version of an overlapping Schwarz preconditioner using one step ILU0 as. Definition: ParallelOverlappingILU0.hpp:131
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()) Definition: PreconditionerFactory_impl.hpp:696
std::function< PrecPtr(const Operator &, const PropertyTree &, const std::function< Vector()> &, std::size_t, const Comm &)> ParCreator Definition: PreconditionerFactory.hpp:77
typename Operator::domain_type Vector Definition: PreconditionerFactory.hpp:68
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) Definition: PreconditionerFactory_impl.hpp:728
std::shared_ptr< Dune::PreconditionerWithUpdate< Vector, Vector > > PrecPtr The type of pointer returned by create(). Definition: PreconditionerFactory.hpp:71
Definition: PressureBhpTransferPolicy.hpp:93
Definition: PressureTransferPolicy.hpp:53
Definition: PropertyTree.hpp:37
T get(const std::string &key) const
Adapter to combine a matrix and another linear operator into a combined linear operator. Definition: WellOperators.hpp:224
Adapter to combine a matrix and another linear operator into a combined linear operator. Definition: WellOperators.hpp:134
DILU preconditioner on the GPU. Definition: CuDILU.hpp:48
Jacobi preconditioner on the GPU. Definition: CuJac.hpp:47
Makes a CUDA preconditioner available to a CPU simulator. Definition: PreconditionerAdapter.hpp:43
Converts the field type (eg. double to float) to benchmark single precision preconditioners. Definition: PreconditionerConvertFieldTypeAdapter.hpp:86
Definition: SupportsFaceTag.hpp:27
Definition: BlackoilPhases.hpp:27
MILU_VARIANT Definition: MILU.hpp:34
@ ILU Do not perform modified ILU.
Dune::OverlappingSchwarzOperator< Dune::BCRSMatrix< MatrixBlock< double, Dim, Dim > >, Dune::BlockVector< Dune::FieldVector< double, Dim > >, Dune::BlockVector< Dune::FieldVector< double, Dim > >, CommPar > OpBPar Definition: PreconditionerFactory_impl.hpp:776
Dune::MatrixAdapter< Dune::BCRSMatrix< Opm::MatrixBlock< double, Dim, Dim > >, Dune::BlockVector< Dune::FieldVector< double, Dim > >, Dune::BlockVector< Dune::FieldVector< double, Dim > > > OpBSeq Definition: PreconditionerFactory_impl.hpp:749
Dune::OverlappingSchwarzOperator< Dune::BCRSMatrix< Dune::FieldMatrix< double, Dim, Dim > >, Dune::BlockVector< Dune::FieldVector< double, Dim > >, Dune::BlockVector< Dune::FieldVector< double, Dim > >, CommPar > OpFPar Definition: PreconditionerFactory_impl.hpp:770
auto setUseFixedOrder(C criterion, bool booleanValue) -> decltype(criterion.setUseFixedOrder(booleanValue)) Definition: PreconditionerFactory_impl.hpp:94
Dune::Amg::SequentialInformation CommSeq Definition: PreconditionerFactory_impl.hpp:740
Dune::MatrixAdapter< Dune::BCRSMatrix< Dune::FieldMatrix< double, Dim, Dim > >, Dune::BlockVector< Dune::FieldVector< double, Dim > >, Dune::BlockVector< Dune::FieldVector< double, Dim > > > OpFSeq Definition: PreconditionerFactory_impl.hpp:745
MILU_VARIANT convertString2Milu(const std::string &milu)
Dune::OwnerOverlapCopyCommunication< int, int > CommPar Definition: PreconditionerFactory_impl.hpp:764
Definition: PreconditionerFactory.hpp:43
Dune::Amg::CoarsenCriterion< CriterionBase > Criterion Definition: PreconditionerFactory.hpp:47
static Criterion criterion(const PropertyTree &prm) Definition: PreconditionerFactory_impl.hpp:106
std::shared_ptr< Dune::PreconditionerWithUpdate< Vector, Vector > > PrecPtr Definition: PreconditionerFactory.hpp:44
static PrecPtr makeAmgPreconditioner(const Operator &op, const PropertyTree &prm, bool useKamg=false) Definition: PreconditionerFactory_impl.hpp:133
static auto args(const PropertyTree &prm) Definition: PreconditionerFactory_impl.hpp:73
Definition: PreconditionerFactory_impl.hpp:57
static auto args(const PropertyTree &prm) Definition: PreconditionerFactory_impl.hpp:58
Definition: PreconditionerFactory_impl.hpp:150
static std::size_t interiorIfGhostLast(const Comm &comm) Definition: PreconditionerFactory_impl.hpp:383
static PreconditionerFactory< Operator, Comm >::PrecPtr createParILU(const Operator &op, const PropertyTree &prm, const Comm &comm, const int ilulevel) Definition: PreconditionerFactory_impl.hpp:359
static void add() Definition: PreconditionerFactory_impl.hpp:151
|