23#include <opm/common/ErrorMacros.hpp>
24#include <opm/common/TimingMacros.hpp>
26#include <opm/simulators/linalg/PreconditionerFactory.hpp>
28#include <opm/simulators/linalg/DILU.hpp>
29#include <opm/simulators/linalg/ExtraSmoothers.hpp>
30#include <opm/simulators/linalg/FlexibleSolver.hpp>
31#include <opm/simulators/linalg/OwningBlockPreconditioner.hpp>
32#include <opm/simulators/linalg/OwningTwoLevelPreconditioner.hpp>
33#include <opm/simulators/linalg/ParallelOverlappingILU0.hpp>
34#include <opm/simulators/linalg/PressureBhpTransferPolicy.hpp>
35#include <opm/simulators/linalg/PressureTransferPolicy.hpp>
36#include <opm/simulators/linalg/PropertyTree.hpp>
37#include <opm/simulators/linalg/WellOperators.hpp>
39#include <opm/simulators/linalg/ilufirstelement.hh>
40#include <opm/simulators/linalg/matrixblock.hh>
42#include <dune/common/unused.hh>
43#include <dune/istl/owneroverlapcopy.hh>
44#include <dune/istl/paamg/amg.hh>
45#include <dune/istl/paamg/fastamg.hh>
46#include <dune/istl/paamg/kamg.hh>
47#include <dune/istl/preconditioners.hh>
50#include <opm/simulators/linalg/PreconditionerFactoryGPUIncludeWrapper.hpp>
55template <
class Smoother>
60 using SmootherArgs =
typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
66 smootherArgs.relaxationFactor = prm.get<
double>(
"relaxation", 1.0);
71template <
class M,
class V,
class C>
77 using SmootherArgs =
typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
80 const int iluwitdh = prm.get<
int>(
"iluwidth", 0);
82 const MILU_VARIANT milu = convertString2Milu(prm.get<std::string>(
"milutype", std::string(
"ilu")));
87 smootherArgs.relaxationFactor = prm.get<
double>(
"relaxation", 1.0);
99void setUseFixedOrder(C&, ...)
104template <
class Operator,
class Comm,
class Matrix,
class Vector>
105typename AMGHelper<Operator, Comm, Matrix, Vector>::Criterion
106AMGHelper<Operator, Comm, Matrix, Vector>::criterion(
const PropertyTree& prm)
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", 1
e-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));
126 setUseFixedOrder(criterion,
true);
130template <
class Operator,
class Comm,
class Matrix,
class Vector>
131template <
class Smoother>
132typename AMGHelper<Operator, Comm, Matrix, Vector>::PrecPtr
133AMGHelper<Operator, Comm, Matrix, Vector>::makeAmgPreconditioner(
const Operator&
op,
134 const PropertyTree& prm,
137 auto crit = criterion(prm);
138 auto sargs = AMGSmootherArgsHelper<Smoother>::args(prm);
141 return std::make_shared<Type>(
142 op,
crit,
sargs, prm.get<std::size_t>(
"max_krylov", 1), prm.get<
double>(
"min_reduction", 1
e-1));
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) {
161 return createParILU(
op, prm, comm, 0);
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);
177 F::addCreator(
"DILU", [](
const O&
op,
const P& prm,
const std::function<
V()>&, std::size_t,
const C& comm) {
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);
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);
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);
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);
206 if constexpr (std::is_same_v<O, Dune::OverlappingSchwarzOperator<M, V, V, C>> ||
207 std::is_same_v<O, Opm::GhostLastMatrixAdapter<M, V, V, C>>) {
208 F::addCreator(
"amg", [](
const O&
op,
const P& prm,
const std::function<
V()>&, std::size_t,
const C& comm) {
209 using PrecPtr = std::shared_ptr<Dune::PreconditionerWithUpdate<V, V>>;
210 const std::string smoother = prm.get<std::string>(
"smoother",
"ParOverILU0");
212 if (smoother ==
"ILU0" || smoother ==
"ParOverILU0") {
216 PrecPtr prec = std::make_shared<Dune::Amg::AMGCPR<O, V, Smoother, C>>(
op,
crit,
sargs, comm);
218 }
else if (smoother ==
"DILU") {
220 using Smoother = Dune::BlockPreconditioner<V, V, C, SeqSmoother>;
221 using SmootherArgs =
typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
224 PrecPtr prec = std::make_shared<Dune::Amg::AMGCPR<O, V, Smoother, C>>(
op,
crit,
sargs, comm);
226 }
else if (smoother ==
"Jac") {
228 using Smoother = Dune::BlockPreconditioner<V, V, C, SeqSmoother>;
229 using SmootherArgs =
typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
232 PrecPtr prec = std::make_shared<Dune::Amg::AMGCPR<O, V, Smoother, C>>(
op,
crit,
sargs, comm);
234 }
else if (smoother ==
"GS") {
236 using Smoother = Dune::BlockPreconditioner<V, V, C, SeqSmoother>;
237 using SmootherArgs =
typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
240 PrecPtr prec = std::make_shared<Dune::Amg::AMGCPR<O, V, Smoother, C>>(
op,
crit,
sargs, comm);
242 }
else if (smoother ==
"SOR") {
244 using Smoother = Dune::BlockPreconditioner<V, V, C, SeqSmoother>;
245 using SmootherArgs =
typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
248 PrecPtr prec = std::make_shared<Dune::Amg::AMGCPR<O, V, Smoother, C>>(
op,
crit,
sargs, comm);
250 }
else if (smoother ==
"SSOR") {
252 using Smoother = Dune::BlockPreconditioner<V, V, C, SeqSmoother>;
253 using SmootherArgs =
typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
256 PrecPtr prec = std::make_shared<Dune::Amg::AMGCPR<O, V, Smoother, C>>(
op,
crit,
sargs, comm);
258 }
else if (smoother ==
"ILUn") {
260 using Smoother = Dune::BlockPreconditioner<V, V, C, SeqSmoother>;
261 using SmootherArgs =
typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
264 PrecPtr prec = std::make_shared<Dune::Amg::AMGCPR<O, V, Smoother, C>>(
op,
crit,
sargs, comm);
267 OPM_THROW(std::invalid_argument,
"Properties: No smoother with name " + smoother +
".");
276 std::size_t pressureIndex,
279 if (pressureIndex == std::numeric_limits<std::size_t>::max()) {
281 "Pressure index out of bounds. It needs to specified for CPR");
283 using Scalar =
typename V::field_type;
285 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy, Comm>>(
288 F::addCreator(
"cprt",
292 std::size_t pressureIndex,
295 if (pressureIndex == std::numeric_limits<std::size_t>::max()) {
297 "Pressure index out of bounds. It needs to specified for CPR");
299 using Scalar =
typename V::field_type;
301 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy, Comm>>(
305 if constexpr (std::is_same_v<O, WellModelGhostLastMatrixAdapter<M, V, V, true>>) {
306 F::addCreator(
"cprw",
310 std::size_t pressureIndex,
313 if (pressureIndex == std::numeric_limits<std::size_t>::max()) {
315 "Pressure index out of bounds. It needs to specified for CPR");
317 using Scalar =
typename V::field_type;
319 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy, Comm>>(
325 F::addCreator(
"GPUILU0", [](
const O&
op,
const P& prm,
const std::function<
V()>&, std::size_t,
const C& comm) {
326 const double w = prm.get<
double>(
"relaxation", 1.0);
327 using field_type =
typename V::field_type;
328 using GpuILU0 =
typename gpuistl::
329 GpuSeqILU0<M, gpuistl::GpuVector<field_type>, gpuistl::GpuVector<field_type>>;
330 auto gpuILU0 = std::make_shared<GpuILU0>(
op.getmat(),
w);
332 auto adapted = std::make_shared<gpuistl::PreconditionerAdapter<V, V, GpuILU0>>(
gpuILU0);
333 auto wrapped = std::make_shared<gpuistl::GpuBlockPreconditioner<V, V, Comm>>(
adapted, comm);
337 F::addCreator(
"GPUJAC", [](
const O&
op,
const P& prm,
const std::function<
V()>&, std::size_t,
const C& comm) {
338 const double w = prm.get<
double>(
"relaxation", 1.0);
339 using field_type =
typename V::field_type;
342 auto gpuJac = std::make_shared<GpuJac>(
op.getmat(),
w);
344 auto adapted = std::make_shared<gpuistl::PreconditionerAdapter<V, V, GpuJac>>(
gpuJac);
345 auto wrapped = std::make_shared<gpuistl::GpuBlockPreconditioner<V, V, Comm>>(
adapted, comm);
349 F::addCreator(
"GPUDILU", [](
const O&
op, [[
maybe_unused]]
const P& prm,
const std::function<
V()>&, std::size_t,
const C& comm) {
350 const bool split_matrix = prm.get<
bool>(
"split_matrix",
true);
352 using field_type =
typename V::field_type;
356 auto adapted = std::make_shared<gpuistl::PreconditionerAdapter<V, V, GpuDILU>>(
gpuDILU);
357 auto wrapped = std::make_shared<gpuistl::GpuBlockPreconditioner<V, V, Comm>>(
adapted, comm);
361 F::addCreator(
"OPMGPUILU0", [](
const O&
op, [[
maybe_unused]]
const P& prm,
const std::function<
V()>&, std::size_t,
const C& comm) {
362 const bool split_matrix = prm.get<
bool>(
"split_matrix",
true);
365 using field_type =
typename V::field_type;
369 auto adapted = std::make_shared<gpuistl::PreconditionerAdapter<V, V, OpmGpuILU0>>(
gpuilu0);
370 auto wrapped = std::make_shared<gpuistl::GpuBlockPreconditioner<V, V, Comm>>(
adapted, comm);
381 using M =
typename F::Matrix;
382 using V =
typename F::Vector;
384 const double w = prm.get<
double>(
"relaxation", 1.0);
385 const bool redblack = prm.get<
bool>(
"redblack",
false);
390 return std::make_shared<ParallelOverlappingILU0<M, V, V, Comm>>(
393 return std::make_shared<ParallelOverlappingILU0<M, V, V, Comm>>(
406 const auto&
is = comm.indexSet();
407 for (
const auto&
ind :
is) {
408 if (Comm::OwnerSet::contains(
ind.local().attribute())) {
421template <
class Operator>
425 using namespace Dune;
427 using C = Dune::Amg::SequentialInformation;
429 using M =
typename F::Matrix;
430 using V =
typename F::Vector;
432 F::addCreator(
"ILU0", [](
const O&
op,
const P& prm,
const std::function<
V()>&, std::size_t) {
433 const double w = prm.get<
double>(
"relaxation", 1.0);
434 return std::make_shared<ParallelOverlappingILU0<M, V, V, C>>(
437 F::addCreator(
"DuneILU", [](
const O&
op,
const P& prm,
const std::function<
V()>&, std::size_t) {
438 const double w = prm.get<
double>(
"relaxation", 1.0);
439 const int n = prm.get<
int>(
"ilulevel", 0);
440 const bool resort = prm.get<
bool>(
"resort",
false);
443 F::addCreator(
"ParOverILU0", [](
const O&
op,
const P& prm,
const std::function<
V()>&, std::size_t) {
444 const double w = prm.get<
double>(
"relaxation", 1.0);
445 const int n = prm.get<
int>(
"ilulevel", 0);
446 return std::make_shared<ParallelOverlappingILU0<M, V, V, C>>(
449 F::addCreator(
"ILUn", [](
const O&
op,
const P& prm,
const std::function<
V()>&, std::size_t) {
450 const int n = prm.get<
int>(
"ilulevel", 0);
451 const double w = prm.get<
double>(
"relaxation", 1.0);
452 return std::make_shared<ParallelOverlappingILU0<M, V, V, C>>(
455 F::addCreator(
"DILU", [](
const O&
op,
const P& prm,
const std::function<
V()>&, std::size_t) {
457 return std::make_shared<MultithreadDILU<M, V, V>>(
op.getmat());
459 F::addCreator(
"Jac", [](
const O&
op,
const P& prm,
const std::function<
V()>&, std::size_t) {
460 const int n = prm.get<
int>(
"repeats", 1);
461 const double w = prm.get<
double>(
"relaxation", 1.0);
464 F::addCreator(
"GS", [](
const O&
op,
const P& prm,
const std::function<
V()>&, std::size_t) {
465 const int n = prm.get<
int>(
"repeats", 1);
466 const double w = prm.get<
double>(
"relaxation", 1.0);
469 F::addCreator(
"SOR", [](
const O&
op,
const P& prm,
const std::function<
V()>&, std::size_t) {
470 const int n = prm.get<
int>(
"repeats", 1);
471 const double w = prm.get<
double>(
"relaxation", 1.0);
474 F::addCreator(
"SSOR", [](
const O&
op,
const P& prm,
const std::function<
V()>&, std::size_t) {
475 const int n = prm.get<
int>(
"repeats", 1);
476 const double w = prm.get<
double>(
"relaxation", 1.0);
482 if constexpr (std::is_same_v<O, Dune::MatrixAdapter<M, V, V>>) {
483 F::addCreator(
"amg", [](
const O&
op,
const P& prm,
const std::function<
V()>&, std::size_t) {
484 const std::string smoother = prm.get<std::string>(
"smoother",
"ParOverILU0");
485 if (smoother ==
"ILU0" || smoother ==
"ParOverILU0") {
488 }
else if (smoother ==
"Jac") {
491 }
else if (smoother ==
"GS") {
494 }
else if (smoother ==
"DILU") {
497 }
else if (smoother ==
"SOR") {
500 }
else if (smoother ==
"SSOR") {
503 }
else if (smoother ==
"ILUn") {
507 OPM_THROW(std::invalid_argument,
"Properties: No smoother with name " + smoother +
".");
510 F::addCreator(
"kamg", [](
const O&
op,
const P& prm,
const std::function<
V()>&, std::size_t) {
511 const std::string smoother = prm.get<std::string>(
"smoother",
"ParOverILU0");
512 if (smoother ==
"ILU0" || smoother ==
"ParOverILU0") {
515 }
else if (smoother ==
"Jac") {
518 }
else if (smoother ==
"SOR") {
521 }
else if (smoother ==
"GS") {
524 }
else if (smoother ==
"SSOR") {
527 }
else if (smoother ==
"ILUn") {
531 OPM_THROW(std::invalid_argument,
"Properties: No smoother with name " + smoother +
".");
534 F::addCreator(
"famg", [](
const O&
op,
const P& prm,
const std::function<
V()>&, std::size_t) {
535 if constexpr (std::is_same_v<typename V::field_type, float>) {
536 OPM_THROW(std::logic_error,
"famg requires UMFPack which is not available for floats");
540 Dune::Amg::Parameters
parms;
541 parms.setNoPreSmoothSteps(1);
542 parms.setNoPostSmoothSteps(1);
547 if constexpr (std::is_same_v<O, WellModelMatrixAdapter<M, V, V, false>>) {
551 if (pressureIndex == std::numeric_limits<std::size_t>::max()) {
552 OPM_THROW(std::logic_error,
"Pressure index out of bounds. It needs to specified for CPR");
554 using Scalar =
typename V::field_type;
557 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy>>(
565 if (pressureIndex == std::numeric_limits<std::size_t>::max()) {
566 OPM_THROW(std::logic_error,
"Pressure index out of bounds. It needs to specified for CPR");
568 using Scalar =
typename V::field_type;
570 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy>>(
576 if (pressureIndex == std::numeric_limits<std::size_t>::max()) {
577 OPM_THROW(std::logic_error,
"Pressure index out of bounds. It needs to specified for CPR");
579 using Scalar =
typename V::field_type;
581 return std::make_shared<OwningTwoLevelPreconditioner<O, V, LevelTransferPolicy>>(
586 F::addCreator(
"GPUILU0", [](
const O&
op,
const P& prm,
const std::function<
V()>&, std::size_t) {
587 const double w = prm.get<
double>(
"relaxation", 1.0);
588 using field_type =
typename V::field_type;
589 using GpuILU0 =
typename gpuistl::
590 GpuSeqILU0<M, gpuistl::GpuVector<field_type>, gpuistl::GpuVector<field_type>>;
591 return std::make_shared<gpuistl::PreconditionerAdapter<V, V, GpuILU0>>(
592 std::make_shared<GpuILU0>(
op.getmat(),
w));
595 F::addCreator(
"GPUILU0Float", [](
const O&
op,
const P& prm,
const std::function<
V()>&, std::size_t) {
596 const double w = prm.get<
double>(
"relaxation", 1.0);
597 using block_type =
typename V::block_type;
598 using VTo = Dune::BlockVector<Dune::FieldVector<float, block_type::dimension>>;
599 using matrix_type_to =
600 typename Dune::BCRSMatrix<Dune::FieldMatrix<float, block_type::dimension, block_type::dimension>>;
601 using GpuILU0 =
typename gpuistl::
602 GpuSeqILU0<matrix_type_to, gpuistl::GpuVector<float>, gpuistl::GpuVector<float>>;
605 auto converted = std::make_shared<Converter>(
op.getmat());
606 auto adapted = std::make_shared<Adapter>(std::make_shared<GpuILU0>(
converted->getConvertedMatrix(),
w));
611 F::addCreator(
"GPUJAC", [](
const O&
op,
const P& prm,
const std::function<
V()>&, std::size_t) {
612 const double w = prm.get<
double>(
"relaxation", 1.0);
613 using field_type =
typename V::field_type;
616 return std::make_shared<gpuistl::PreconditionerAdapter<V, V, GPUJac>>(
617 std::make_shared<GPUJac>(
op.getmat(),
w));
620 F::addCreator(
"OPMGPUILU0", [](
const O&
op, [[
maybe_unused]]
const P& prm,
const std::function<
V()>&, std::size_t) {
621 const bool split_matrix = prm.get<
bool>(
"split_matrix",
true);
625 using field_type =
typename V::field_type;
631 F::addCreator(
"GPUDILU", [](
const O&
op, [[
maybe_unused]]
const P& prm,
const std::function<
V()>&, std::size_t) {
632 const bool split_matrix = prm.get<
bool>(
"split_matrix",
true);
634 using field_type =
typename V::field_type;
639 F::addCreator(
"GPUDILUFloat", [](
const O&
op, [[
maybe_unused]]
const P& prm,
const std::function<
V()>&, std::size_t) {
640 const bool split_matrix = prm.get<
bool>(
"split_matrix",
true);
643 using block_type =
typename V::block_type;
644 using VTo = Dune::BlockVector<Dune::FieldVector<float, block_type::dimension>>;
645 using matrix_type_to =
typename Dune::BCRSMatrix<Dune::FieldMatrix<float, block_type::dimension, block_type::dimension>>;
649 auto converted = std::make_shared<Converter>(
op.getmat());
658template <
class Operator,
class Comm>
664template <
class Operator,
class Comm>
666PreconditionerFactory<Operator, Comm>::instance()
672template <
class Operator,
class Comm>
674PreconditionerFactory<Operator, Comm>::doCreate(
const Operator&
op,
675 const PropertyTree& prm,
677 std::size_t pressureIndex)
680 StandardPreconditioners<Operator, Comm>::add();
683 const std::string& type = prm.get<std::string>(
"type",
"ParOverILU0");
684 auto it = creators_.find(type);
685 if (it == creators_.end()) {
686 std::ostringstream
msg;
687 msg <<
"Preconditioner type " << type <<
" is not registered in the factory. Available types are: ";
688 for (
const auto& prec : creators_) {
689 msg << prec.first <<
' ';
697template <
class Operator,
class Comm>
699PreconditionerFactory<Operator, Comm>::doCreate(
const Operator&
op,
700 const PropertyTree& prm,
702 std::size_t pressureIndex,
706 StandardPreconditioners<Operator, Comm>::add();
709 const std::string& type = prm.get<std::string>(
"type",
"ParOverILU0");
710 auto it = parallel_creators_.find(type);
711 if (it == parallel_creators_.end()) {
712 std::ostringstream
msg;
713 msg <<
"Parallel preconditioner type " << type <<
" is not registered in the factory. Available types are: ";
714 for (
const auto& prec : parallel_creators_) {
715 msg << prec.first <<
' ';
723template <
class Operator,
class Comm>
725PreconditionerFactory<Operator, Comm>::doAddCreator(
const std::string& type, Creator
c)
730template <
class Operator,
class Comm>
732PreconditionerFactory<Operator, Comm>::doAddCreator(
const std::string& type, ParCreator
c)
734 parallel_creators_[type] =
c;
737template <
class Operator,
class Comm>
742 std::size_t pressureIndex)
747template <
class Operator,
class Comm>
753 std::size_t pressureIndex)
759template <
class Operator,
class Comm>
764 std::size_t pressureIndex)
766 return instance().doCreate(
op, prm, std::function<Vector()>(), pressureIndex, comm);
769template <
class Operator,
class Comm>
773 instance().doAddCreator(type,
creator);
776template <
class Operator,
class Comm>
780 instance().doAddCreator(type,
creator);
783using CommSeq = Dune::Amg::SequentialInformation;
785template<
class Scalar,
int Dim>
786using OpFSeq = Dune::MatrixAdapter<Dune::BCRSMatrix<Dune::FieldMatrix<Scalar, Dim, Dim>>,
787 Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>,
788 Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>>;
789template<
class Scalar,
int Dim>
790using OpBSeq = Dune::MatrixAdapter<Dune::BCRSMatrix<MatrixBlock<Scalar, Dim, Dim>>,
791 Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>,
792 Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>>;
794template<
class Scalar,
int Dim,
bool overlap>
796 Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>,
797 Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>,
800template<
class Scalar,
int Dim,
bool overlap>
802 Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>,
803 Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>,
807using CommPar = Dune::OwnerOverlapCopyCommunication<int, int>;
809template<
class Scalar,
int Dim>
810using OpFPar = Dune::OverlappingSchwarzOperator<Dune::BCRSMatrix<Dune::FieldMatrix<Scalar, Dim, Dim>>,
811 Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>,
812 Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>,
815template<
class Scalar,
int Dim>
816using OpBPar = Dune::OverlappingSchwarzOperator<Dune::BCRSMatrix<MatrixBlock<Scalar, Dim, Dim>>,
817 Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>,
818 Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>,
820template<
class Scalar,
int Dim>
822 Dune::BlockVector<Dune::FieldVector<Scalar,Dim>>,
823 Dune::BlockVector<Dune::FieldVector<Scalar,Dim>>,
826template<
class Scalar,
int Dim>
828 Dune::BlockVector<Dune::FieldVector<Scalar,Dim>>,
829 Dune::BlockVector<Dune::FieldVector<Scalar,Dim>>,
832#define INSTANTIATE_PF_PAR(T,Dim) \
833 template class PreconditionerFactory<OpBSeq<T,Dim>, CommPar>; \
834 template class PreconditionerFactory<OpFPar<T,Dim>, CommPar>; \
835 template class PreconditionerFactory<OpBPar<T,Dim>, CommPar>; \
836 template class PreconditionerFactory<OpGLFPar<T,Dim>, CommPar>; \
837 template class PreconditionerFactory<OpGLBPar<T,Dim>, CommPar>; \
838 template class PreconditionerFactory<OpW<T,Dim, false>, CommPar>; \
839 template class PreconditionerFactory<OpWG<T,Dim, true>, CommPar>; \
840 template class PreconditionerFactory<OpBPar<T,Dim>, CommSeq>; \
841 template class PreconditionerFactory<OpGLBPar<T,Dim>, CommSeq>;
844#define INSTANTIATE_PF_SEQ(T,Dim) \
845 template class PreconditionerFactory<OpFSeq<T,Dim>, CommSeq>; \
846 template class PreconditionerFactory<OpBSeq<T,Dim>, CommSeq>; \
847 template class PreconditionerFactory<OpW<T,Dim, false>, CommSeq>; \
848 template class PreconditionerFactory<OpWG<T,Dim, true>, CommSeq>;
851#define INSTANTIATE_PF(T,Dim) \
852 INSTANTIATE_PF_PAR(T,Dim) \
853 INSTANTIATE_PF_SEQ(T,Dim)
855#define INSTANTIATE_PF(T,Dim) INSTANTIATE_PF_SEQ(T,Dim)
Parallel algebraic multigrid based on agglomeration.
Definition amgcpr.hh:88
Definition PreconditionerWithUpdate.hpp:43
The OpenMP thread parallelized DILU preconditioner.
Definition DILU.hpp:53
Dune linear operator that assumes ghost rows are ordered after interior rows.
Definition WellOperators.hpp:376
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:739
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:771
std::shared_ptr< Dune::PreconditionerWithUpdate< Vector, Vector > > PrecPtr
The type of pointer returned by create().
Definition PreconditionerFactory.hpp:71
Definition PressureBhpTransferPolicy.hpp:99
Definition PressureTransferPolicy.hpp:55
Definition PropertyTree.hpp:37
Adapter to combine a matrix and another linear operator into a combined linear operator.
Definition WellOperators.hpp:273
Adapter to combine a matrix and another linear operator into a combined linear operator.
Definition WellOperators.hpp:179
DILU preconditioner on the GPU.
Definition GpuDILU.hpp:44
Jacobi preconditioner on the GPU.
Definition GpuJac.hpp:47
ILU0 preconditioner on the GPU.
Definition OpmGpuILU0.hpp:46
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 blackoilboundaryratevector.hh:37
MILU_VARIANT
Definition MILU.hpp:34
@ ILU
Do not perform modified ILU.
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
Definition PreconditionerFactory.hpp:43
Definition PreconditionerFactory_impl.hpp:57
Definition PreconditionerFactory_impl.hpp:150
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:402